简介压缩包提供黄土高原地区DEM高程数据及四大类生态分区成果面向地理信息、水土保持与生态规划方向的从业者和研究者适用于ArcMap环境下的流域分析、地貌制图与生态修复区划。包内共12个文件压缩包约4.78MB含tif格式高程栅格、可编辑mxd工程文档、docx数据说明以及配套shp/shx/dbf/prj等标准shapefile组件可直接加载或二次制图。已有136人下载学习。生态分区依据气候、地形及治理模式差异将黄土高原划分为高塬沟壑区、丘陵沟壑区、沙地与农灌区、土石山区及河谷平原区并细分子区附面积与涉及县域信息。DEM数据可支撑水文分析、坡面提取等研究分区成果可作为水土流失治理、土地整治与生态保护规划的基础依据。1. 从一份真实需求说起黄土高原DEM高程与4大类生态分区为什么必须三件套一起交付一个做生态修复立项的朋友来找我说评审专家要求提供「带地形背景的生态分区底图」。他手上有DEM高程栅格、一张cad格式的分区草图、还有一份打印出来的分区图jpg——三样东西互相对不上坐标不统一、格式不一样根本没法在一张图里叠起来看。我帮他整理成了三件套一个可编辑的mxd工程文件、一份标准shapefile、一份标准成图TIF。所谓黄土高原DEM高程及4大类生态分区核心就是这个数据组织问题高程地形做底生态分区做面工程文件做状态保存矢量格式做入库交换栅格成图做汇报交付。这篇文章面向GIS数据处理、水土保持、国土空间规划方向的从业者讲清楚这三样东西分别是什么、怎么做出来、哪里最容易翻车。2. 拆开看数据DEM高程、生态分区与三种交付格式的本质2.1 DEM高程数字化高程栅格放进TIF前必须确认的四个属性DEM高程本质上是一张单波段栅格每个像元的值代表该处的地面海拔单位通常是米。它本身没有任何颜色你在ArcMap里看到的绿色到棕色的渐变是符号化渲染的结果不是数据自带的属性。这一点很多人一开始会弄混拿到一个TIF打开发现是灰色或者黑色就以为数据坏了其实只是没有做拉伸或者色带渲染。我一般在处理DEM时会确认四个属性像素类型、NoData背景值、坐标系、范围。像素类型常见是浮点型或整型浮点型保留小数高程整型则会有取舍NoData背景值决定范围外区域是否参与渲染和计算坐标系决定它能和什么数据叠合范围则决定了栅格的边界和面积统计是否准确。这四个属性只要有一个不对后续做坡度分析、填挖方、面积统计都会出问题而且问题往往藏得很深表面看不出来。拿到一套黄土高原的DEM数据时我通常先做一次「体检」用gdalinfo看坐标系和像素信息再叠加一幅自己信得过的边界矢量看一眼范围对不对。很多公开DEM的数据范围比你实际需要的范围大得多直接拿去出图会有大片无数据区图面难看不说面积统计也是错的。2.2 标准shapefile4大类生态分区在属性表里如何组织shapefile是GIS里最通用的矢量交换格式所谓「标准」主要体现在属性表字段设计和几何类型上。标题里的4大类生态分区常见口径是水土保持功能类、植被覆盖类、地形地貌类、土地利用类每一类下面还可以细分二级类型。你可以根据项目实际口径替换但字段结构是通用的。我建分区属性表一般按这组字段来设计字段名类型说明FID整型要素编号系统自带type_code文本或整型分区代码建议用数字排序如101、102type_name文本分区名称如「黄土塬水土保持区」area_km2双精度分区面积单位平方公里出图时可以直接用data_src文本数据来源或划分依据便于追溯为什么要把type_code和type_name分开因为shapefile的属性表格式DBF对字段名长度有限制而且有些老旧软件对中文支持不好。用一个短代码做主键中文名称只在显示层面用这是行业里最常见的做法后面讲坑的时候会细说。几何类型上生态分区用面要素Polygon每个分区的边界要闭合不能有自相交。我一般会在出图前先把area_km2字段算好用正确的投影坐标系算面积然后在属性表里校验一次所有分区面积之和应该等于整个研究区总面积误差不超过0.5%。这步不用脚本ArcMap里就能做但它能挡掉很大一部分几何错误。2.3 可编辑mxd的真正职责它不存数据只存「图纸状态」mxd是ArcMap的工程文件很多交付方和接收方都对它有误解。mxd本身不包含任何数据它只记录了三件事数据源的路径、图层的符号化设置、图版要素的位置和样式。你可以把mxd理解成一份「菜谱」数据是食材菜谱只写怎么搭配和烹调菜谱丢了食材还在但做出来的菜就未必是你想要的样子了。所以当标题写「可编辑mxd文件」时它的价值在于接收方打开mxd看到的是一张已经配好符号、修好图例、摆好比例尺的图纸不需要自己再花半天去调。DEM的色带渐变、4大类和二级分区的配色、标注的字号、图例的排列顺序全部在这个mxd里保存着。这三种交付物对照关系大概是这样的文件类型用途可编辑性mxd工程文件工程文件保存符号化与图版状态高打开即改shapefile矢量数据数据入库、空间分析、数据交换属性可编辑TIF栅格栅格成图汇报、打印、做底图低改起来麻烦一个容易踩的误区是把mxd当数据备份觉得mxd在就万事大吉。实际上mxd里的路径一旦失效打开就是一片红色感叹号。所以我在交付时一定会连同整套数据一起打包并且把mxd的路径设置为相对路径这部分在避坑章节细讲。3. 从散数据到三件套用脚本把坐标系、字段与图层状态一次调齐3.1 先搭目录规范一个工程文件夹撑起所有数据动手处理数据之前我建议先建一个标准目录结构因为mxd的路径管理直接和目录结构绑定后面省下的时间远大于现在多花的两分钟。我常用的结构是project文件夹下分data、mxd、output三个子目录data里再分dem和shp两类。# 建立标准工程目录 mkdir -p loess_project/data/dem mkdir -p loess_project/data/shp mkdir -p loess_project/mxd mkdir -p loess_project/output/tif # 目录结构说明 # data/dem 放原始DEM与重投影后的DEM # data/shp 放生态分区shapefile及配套文件 # mxd 放工程文件 # output/tif 放最终导出的标准成图这段命令没有技术含量但值得你养成习惯。我见过太多人把shp、mxd、tif散落在桌面和各个下载文件夹里等交付的时候才发现某个关键文件被移动过mxd整个失效。把目录先立住后面所有脚本都基于这个相对结构去写mxd的相对路径才能稳定。3.2 统一坐标系与矢栅对齐GDAL脚本一次跑完接下来处理最核心的一步把DEM和生态分区shp统一到同一个坐标系下。黄土高原横跨多个投影带如果使用低精度地理坐标或默认的Web墨卡托面积统计和坡度分析会严重失真。我一般使用适合全国中纬度地区的Albers等积投影或CGCS2000下按区域就近分带的高斯投影具体选哪个取决于你要不要算面积——算面积就用等积投影。下面这个Python脚本使用GDAL库完成重投影、重采样和按矢量范围裁剪from osgeo import gdal, ogr, osr # 目标投影Albers等积投影适合面积统计和区域制图 target_srs osr.SpatialReference() target_srs.ImportFromEPSG(9822) # EPSG 9822为Albers等积投影 # 打开矢量分区获取其范围作为裁剪边界 shp_path data/shp/ecotype.shp ds_vec ogr.Open(shp_path) lyr ds_vec.GetLayer(0) extent lyr.GetExtent() # (minX, maxX, minY, maxY) # 打开原始DEM src_path data/dem/raw_dem.tif src_ds gdal.Open(src_path) # 重投影并裁剪目标像元大小设为90米接近SRTM/ALOS常用分辨率 gdal.Warp( data/dem/dem_albers.tif, src_ds, dstSRStarget_srs, outputBounds(extent[0], extent[2], extent[1], extent[3]), # 用分区范围裁剪DEM xRes90, yRes90, resampleAlgbilinear, # 高程用双线性内插保留平滑过渡 srcNodata-32768, # 原始DEM的NoData值 dstNodata-32768, creationOptions[COMPRESSLZW], # 无损压缩缩小文件体积 ) print(DEM已重投影并裁剪完成)逻辑说明脚本先拿到分区shp的范围用这个范围去裁剪DEM这样DEM和shp在空间范围上完全一致后续做任何叠合分析都不会出界。重采样方法选了bilinear这是因为高程是连续表面用最近邻会产生明显的锯齿边界用双线性更平滑。如果你要拿DEM做汇水分析可以换用cubic但普通出图和叠合显示bilinear就够了。两个参数值得你注意target_srs直接指定了Albers投影这不是「顺便选的」等积投影能保证图斑面积与实际相符做生态分区面积统计时不会被投影变形干扰dstNodata必须要和源数据一致否则裁剪后的DEM会出现一圈假的小值渲染出来就是一片不合理的深色。3.3 构建可编辑mxd挂图层、配符号、摆图版坐标系统一之后下一步是构建mxd工程文件。用arcpy有两种思路一是用ArcMap界面手动搭建一次存成mxd模板后续用脚本替换数据源二是完全用arcpy从零创建。我推荐前者因为界面调的符号更直观但为了让读者能抄作业这里给出一个用arcpy创建mxd并挂接图层、设置符号化的可运行脚本import arcpy # 设置工作空间与数据路径 mxd_path mxd/loess_eco.mxd dem_path data/dem/dem_albers.tif shp_path data/shp/ecotype.shp # 创建空白地图文档添加图层 mxd arcpy.mapping.MapDocument(CURRENT) df arcpy.mapping.ListDataFrames(mxd, Layers)[0] # 先加分区矢量再加DEM作为底图保证图层顺序为DEM在下、分区在上 arcpy.mapping.AddLayer(df, arcpy.mapping.Layer(shp_path), BOTTOM) arcpy.mapping.AddLayer(df, arcpy.mapping.Layer(dem_path), AUTO) # 给DEM渲染渐变高程色带 dem_lyr arcpy.mapping.ListLayers(mxd, dem_albers, df)[0] sym dem_lyr.symbology if sym.classes: # 按海拔分8级颜色由低到高渐变 sym.updateClassification(EQUAL_INTERVAL, 8) sym.classBreakValues [0, 500, 800, 1000, 1200, 1500, 2000, 2500, 3000] dem_lyr.symbology sym # 给分区图层按type_code做唯一值渲染 eco_lyr arcpy.mapping.ListLayers(mxd, ecotype, df)[0] eco_lyr.symbology.renderer UNIQUE_VALUES eco_lyr.symbology.valueField type_code eco_lyr.symbology.addAllValues() eco_lyr.symbology.sortValues(type_code, ascendingTrue) # 按代码排序而非拼音或字典序 mxd.saveACopy(mxd_path) print(mxd已生成 mxd_path)逻辑说明脚本先从当前ArcMap会话创建地图文档然后按顺序添加图层。注意添加顺序用了AddLayer的placement参数确保DEM在下面做底分区在上面做主题。符号化部分DEM用等间隔分成8级高程范围按黄土高原常见海拔来分段分区图层按type_code用唯一值渲染这部分决定了图例的显示顺序。这个脚本我建议在ArcMap的Python窗口里运行而不是直接跑在外部IDE。因为arcpy.mapping.MapDocument(CURRENT)需要一个打开的ArcMap会话。如果你是完全批量化场景可以改成不依赖会话的方式但对于绝大多数项目一次搭建好mxd模板比每次都跑脚本要实用得多。脚本跑完后你还需要在图版视图中手动加上图例、指北针、比例尺、标题和来源注记。这些要素用脚本也能加但调整位置很费劲推荐在ArcMap界面里用Layout视图摆一次然后把mxd存档。摆图版消耗的是审美和耐心脚本帮不了你。3.4 导出标准成图TIF栅格化的是图纸不是地形标题里的「标准成图TIF」和原始DEM数据TIF是两回事。成图TIF是把整个地图版面栅格化之后得到的图片它包含的是渲染后的地图内容——色带DEM、分区边界、图例、比例尺、标题、指北针全都在一个文件里。别人拿到这个TIF可以直接放进PPT也可以用看图软件打开就得到完整信息不需要装任何GIS软件。用arcpy从mxd导出成图TIF的标准流程import arcpy mxd_path mxd/loess_eco.mxd out_tif output/tif/黄土高原生态分区成图.tif # 打开mxd按页面尺寸导出 mxd arcpy.mapping.MapDocument(mxd_path) # 导出参数分辨率300dpi背景白色LZW压缩 arcpy.mapping.ExportToTIFF( mxd, out_tif, resolution300, geotiffFalse, # 成图不需要地理配准信息关掉避免文件过大 world_fileFalse, color_modeRGB, compressionLZW, background_color(255, 255, 255) # 白色背景避免透明底变黑 ) print(成图TIF已导出 out_tif)参数说明resolution设为300dpi是为了保证打印和PPT放大都清晰低于150dpi时图上的小字会虚。geotiff参数这里特意设成了False原因和很多人直觉相反——成图是给人看的不是用来做空间分析的带地理坐标只会让文件体积变大而且有些看图软件对带坐标的tif另眼相看。背景色设成白色是防止导出时背景透明在部分软件里显示成黑底。实际操作中我建议在导出前先打开mxd确认一遍图版没有要素出界。你可以在Layout视图里切换到预览模式看看整页效果这一步虽然传统但比任何脚本都可靠。4. 排查避坑五处最容易翻车的位置按这个顺序查4.1 叠合不上坐标系「黑匣子」问题现象把DEM和生态分区shp放进同一张图两者位置差了几百米甚至几公里zoom to layer也救不回来。原因最常见是DEM是WGS84地理坐标而shp是CGCS2000投影坐标或者两个数据用了不同的投影带。ArcMap的「动态投影」能在显示层面临时对齐但它不会改动你的源数据一旦导出或做分析错位还会回来。解决在动手做三件套之前先做一次坐标系体检。用gdalinfo读取两个数据源的坐标系再统一执行第3章的GDAL重投影脚本。记住一个原则源数据不做原地修改全部通过脚本输出到data/dem和data/shp目标文件夹这样原始数据保留排查问题时还能溯源。4.2 mxd打开一片红感叹号所有图层链接全部断裂现象把交付包解压后打开mxd图层面板里全是红底白色感叹号一个图层都加载不出来。原因mxd记录的是数据源的绝对路径。如果对方把mxd和数据放在不同盘符或者解压目录结构和压缩包内不一致绝对路径就失效了。很多初学者打包交付时只拷了mxd和几个关键文件却漏掉了同级的gdb或栅格数据集。解决在保存mxd之前打开文档属性勾选「存储数据源的相对路径」然后重新保存。相对路径让mxd按自身所在位置找数据只要mxd和数据文件夹的相对位置不变整个包挪到任何机器都能打开。交付前再做一个测试把文件夹复制到另一个位置打开一遍mxd确认没有感叹号再发出去。4.3 DEM背景在成图TIF里显示为黑色方块现象DEM数据本身没问题在ArcMap里显示正常但导出成图TIF后研究区边界外围是一大片黑块整个图面很难看。原因DEM范围外是NoData值在ArcMap里你设置了NoData显示为透明但导出成图时背景颜色处理逻辑不同NoData区域被渲染成了黑色。如果你用的是透明背景导出这部分就会变成黑。这属于符号化设置没有传递到导出环节。解决有两个办法推荐组合使用。第一步在DEM图层符号化里把NoData颜色显式设为白色或无色感第二步在导出时按第3章的参数明确background_color为白色。这样即使NoData区域被渲染出的图也是干净的白色底。4.4 shapefile字段名被截断或中文乱码现象在ArcMap属性表里建好了一个字段叫「水土保持功能区」保存关闭再打开字段名变成「水土保持」或者打开时出现问号乱码。原因shapefile的属性文件DBF对字段名长度限制很严格一般按字节算中文占用较多而且不同编码环境下中文兼容性很差换个软件打开就乱码。这不是单个软件的问题是shapefile这种老格式的固有限制。解决字段设计时全部用短英文代码如type_code、type_name、area_km2中文含义放到单独的分区代码对照表里交付时附一份Excel或CSV即可。属性表里的type_name字段存中文没问题只要字段名是英文短代码就行。这是所有shapefile数据都该遵守的纪律和本次项目无关也要坚持。4.5 图例顺序乱套唯一值渲染的顺序是玄学现象mxd里图例的顺序不是按业务逻辑排列的比如「水土保持区」排在了「植被覆盖区」前面或者二级分区顺序错乱怎么拖都拖不对。原因唯一值符号化默认按代码排序如果你的type_code设计成了中文如「一类」「二类」或者没有考虑业务顺序图例就会按字典序、拼音序或者入库顺序排。ArcMap的图例排序逻辑有时不符合直觉手动调整后一旦刷新又变回去。解决从根本设计上解决——type_code设计成数字编码如101、102、201、202第一位是大类后两位是二级类。数字排序天然符合业务顺序图例和符号化都会按这个顺序走。如果mxd已经做完了可以通过arcpy脚本对图层的valueField排序或者在图例属性里手动上移下移。但我建议回到字段设计去改一劳永逸。5. 交付前做这三件事数据核验、审图清单与「相对路径」收尾习惯5.1 跑一遍数据一致性核验脚本交付之前我建议用一段简单的Python脚本对三个交付物做最后体检。它不会修复任何问题但能在你打包之前发现「坐标系不一致」「范围对不上」「字段缺失」这类硬伤避免交付后再被打回。from osgeo import gdal, ogr import json report {} # 1. 检查DEM是否已重投影、是否有NoData dem gdal.Open(data/dem/dem_albers.tif) dem_srs dem.GetProjection() dem_extent dem.GetGeoTransform() # 栅格范围 report[dem_srs] dem_srs.split()[3] if dem_srs else NONE report[dem_nodata] dem.GetRasterBand(1).GetNoDataValue() # 2. 检查shp是否与DEM在同一坐标系、范围能否对上 vec ogr.Open(data/shp/ecotype.shp) lyr vec.GetLayer(0) vec_srs lyr.GetSpatialRef().ExportToWkt().split()[3] vec_extent lyr.GetExtent() # 3. 检查mxd文件是否存在、大小是否正常 import os mxd_size os.path.getsize(mxd/loess_eco.mxd) if os.path.exists(mxd/loess_eco.mxd) else -1 report[mxd_exists] mxd_size 0 # 4. 范围对齐判断矢量范围是否落在DEM范围内 dem_min_x, dem_min_y, dem_max_x, dem_max_y \ dem_extent[0], dem_extent[3] dem_extent[5]*dem.RasterYSize, \ dem_extent[0] dem_extent[1]*dem.RasterXSize, dem_extent[3] range_ok (vec_extent[0] dem_min_x - 100 and vec_extent[2] dem_max_x 100) report[range_aligned] bool(range_ok) print(json.dumps(report, ensure_asciiFalse, indent2))逻辑说明脚本读取DEM的坐标系、NoData设置和四至范围读取shp的坐标系和四至范围比对两者在100米容差内是否对齐最后确认mxd文件存在且非空。输出一份JSON结果一目了然。100米容差的来源是我们重投影时允许的边界偏移如果shp比DEM范围大太多说明裁剪步骤没做对。5.2 对照审图清单看一遍图纸数据检查过了不代表图能直接交付。我出图前会按下面这张清单过一遍图例是否包含全部分区类型顺序是否与业务口径一致比例尺数值是否合理放在图面不被其他要素遮挡指北针方向是否与地图实际方向一致标题是否包含区域名称、专题名称和日期DEM色带有没有被压缩得层次不清彩色打印是否友好分区边界线和标注在成图TIF里是否清晰分辨率300dpi下不虚图外空白处是否干净有没有残留的十字丝或编辑辅助点这张清单也是mxd「可编辑」价值的一种体现——审图发现问题直接改mxd的图版要素重新导出TIF三分钟搞定。如果没有mxd只有一张成品图改一个错别字都要找源文件从头来。5.3 收尾习惯交付前把mxd路径「相对化」一次最后分享一个血泪换来的习惯。有次我给别人交了一套数据文件夹结构是这样的根目录下放mxddata子目录放shp和dem我自认为没有问题结果对方解压后打开红感叹号一片。后来排查发现我在存mxd的时候没有勾选相对路径所有图层记录的都是本机的绝对路径换一台机器必然失效。从那以后我给自己定了一条规矩无论是不是自己用交付前都重设一遍数据源路径。在ArcMap里打开文档属性勾选「存储数据源的相对路径」保存然后把这个文件夹完整复制到另一个位置重新打开测试一次。测试通过才打包交付。这条规矩看着简单但能挡掉绝大多数「我明明都放进去了你为什么打不开」的扯皮。三件套的价值最后都落在「别人拿到就能用」这五个字上。mxd保存的是图纸状态shp是可交换的数据本体TIF是零门槛可查看的成果表现三样东西缺一样接收方就会多花时间自己补。补的东西可能和你的原意有偏差来回沟通的精力远大于整理时的半小时。希望这篇笔记能帮你一次就把三件套做齐希望帮到你。本文还有配套的精品资源点击获取