简介:这套数据面向地理信息、水土保持与生态规划相关的科研人员和工程技术人员,提供黄土高原DEM高程标准成图以及按国家相关部门划分原则构建的四大生态分区矢量数据。分区依据气候、地貌、水土流失治理模式与生态恢复工程差异,在县级行政边界基础上通过空间分析与分区统计划定,涵盖黄土高塬沟壑区、黄土丘陵沟壑区、沙地与农灌区、土石山区及河谷平原区,可用于水土流失治理、生态修复和土地利用规划等场景。压缩包共12个文件,约4.78MB,核心内容包括标准TIF格式的DEM成图、Shapefile格式的生态分区矢量层,以及可编辑的MXD制图工程文件;配套的dbf、prj、sbn等文件保证了属性表、投影信息与空间索引的完整性,在ArcMap中可直接加载出图或进一步做流域分析、地貌分区等操作。目前已有135人学习,适合需要标准基础底图与生态分区边界开展研究、教学或项目制图的人员直接取用。
1. 黄土高原DEM高程与生态分区:一份能直接打开的mxd+shape+TIF数据包
黄土高原大半个省跨度的地形数据,自己拼DEM、自己画生态分区、自己配图层符号,整套搞下来少说一周。这份资源把三个最常见的痛点在同一个包里解决掉:DEM高程栅格、4大类生态分区矢量,以及把二者揉在一起的可编辑mxd工程和标准成图TIF。换句话说,它不是单个DEM文件,而是一套“打开就能看、要改能动手”的黄土高原地理底图。适合做区域生态评估、水土保持分区,也适合课程设计和横向项目里被地形数据折腾过的从业者。至于刚接触ArcMap、想知道mxd和shape怎么配合用的新手,按后文流程走一遍,能少踩不少坑。
2. 先拆包:DEM、生态分区shape、mxd工程和TIF成图各是什么角色
打开压缩包,看到的通常是四类文件混在一起:tif栅格、shp矢量、mxd工程文档,以及配套的图层样式。很多第一次拿到这类资源的人第一反应是双击mxd,然后被红色感叹号吓退。先别急着打开,把四类数据的角色分清,后面所有操作就有了判断依据。
2.1 黄土高原DEM栅格:高程基准和投影,先统一再谈分析
DEM(Digital Elevation Model)记录的是地表高程。黄土高原这种地形破碎、沟壑密集的区域,DEM的细节直接决定后续坡度、坡向、地形起伏度计算的可靠性,也直接关系到生态分区边界的实际投放效果。
拿到DEM后,第一件事不是看颜色,而是右键图层属性,切到“源”选项卡,检查三样东西:像元大小(Cell Size)、空间参考(Spatial Reference)、NoData值。如果空间参考那里写的是Unknown,后面所有空间分析都站不住。常见的一个误区是把DEM当普通图片直接叠加,等发现和生态分区shape对不上时,已经做完半轮分析。
这里多说一句:数据包里给的是DEM,不是DSM。DSM是数字表面模型,地表建筑、树冠都在里面;DEM是滤掉这些地物之后的裸地高程,所以它适合做水文分析、沟壑提取和生态分区底图,而不是做城市建筑高度之类的研究。如果你的原始数据来源是第三方下载工具,拿到的是带建筑信息的DSM,那要先做滤波生成DEM,再进入这套流程。
2.2 4大类生态分区shape:属性表才是真正的分区依据
生态分区的矢量文件通常是面要素(Polygon)。标题里写的4大类分区,常见对应风沙区、丘陵沟壑区、高原沟壑区、土石山区这类水土保持类型区划分,但不同数据来源对类目的命名和编码不一样。拿到shape后先打开属性表,找TYPE、CLASS或NAME字段,以里面的枚举值为准,不要凭图形颜色猜分区。
shape文件是一个多文件格式,至少需要.shp(几何)、.shx(索引)、.dbf(属性)、.prj(投影)四个文件同时存在。只拷一个.shp给别人,对方打开后属性表是空的、图形位置是错的,因为投影和属性全留在了你电脑上。所谓“标准shape文件”,最基本的一条就是.prj必须存在。缺少.prj的shape只能叫图形,不能叫空间数据。
4大类分区在分析里的角色有两种:一是作为掩膜去裁剪DEM,把研究区外的地形切掉;二是作为分区统计的区划单元,统计每个区的高程均值、起伏度、坡度分布。同一个面图层,两种用途在ArcMap里的处理路径完全不同,一个走Clip或Extract by Mask,一个走Zonal Statistics,具体见第4章和第6章。
另外,属性表字段名不建议用中文。ArcMap虽然支持中文字段,一旦导出到GeoJSON或转给其他GIS平台,编码问题会非常头疼。如果这份资源的属性表里已经有中文备注字段,给它补一个英文短字段做主键,后面写代码循环出图会方便很多。
2.3 mxd工程和成图TIF:一个装流程,一个装结果
mxd是ArcMap的工程文件,本质上记录的是一堆图层路径、符号化规则、标注和整饰要素,不是数据本身。这个特性决定了两件事:第一,mxd换一台电脑、换一个目录,图层路径就会失效,出现满屏红色感叹号;第二,mxd里调好的配色、图例、比例尺、图廓,才是可编辑的“生产车间”。
成图TIF则是从这个车间里导出的成品,已经带图例、比例尺、图廓,任何看图软件都能打开,交付给不装GIS的同事最省事。代价是图例和文字已经像素化,想改一个分区名字只能回到mxd里改完重新导出。标准成图TIF不是分析数据,是汇报材料。
这份资源和单独一份DEM数据的价值差别就在这里:DEM给你算,shape给你圈,mxd给你改,TIF给你交差。四件套的分工,先用一张表记住。
| 资源组件 | 格式 | 在项目里的角色 |
|---|---|---|
| DEM高程 | TIF栅格 | 提供连续地表高程,用于坡度、起伏度、等高线提取 |
| 生态分区 | Shapefile面 | 提供4大类分区空间范围,用于裁剪和分区统计 |
| mxd工程 | ArcMap文档 | 保存图层路径、符号化、标注和整饰,可编辑调整 |
| 成图TIF | TIF栅格 | 已经配好图廓图例的标准成图,直接交付汇报 |
3. 打开mxd工程:数据源修复、图层叠放和坐标系“三板斧”
mxd打开后能不能正常工作,取决于三个问题:数据路径对不对、图层顺序对不对、坐标系统一没统一。这三个问题按顺序处理,比在ArcMap里乱点要快得多。
3.1 mxd停用路径失效:用arcpy批量修复数据源
mxd默认记录的是绝对路径。比如作者本机路径是C盘某个目录,你解压到了D盘的huangtu文件夹,mxd里的每个图层都会找不到原始数据,显示红色感叹号。修复方式有两种:一是把mxd和shape、TIF按压缩包里的目录结构完整解压,不改变相对位置;二是写脚本批量重定向数据源。第二种在数据图层多的时候更实用。
import arcpy # mxd是ArcMap工程文件,用arcpy.mapping读取 mxd = arcpy.mapping.MapDocument(r"E:\huangtu\huangtu_eco.mxd") for lyr in arcpy.mapping.ListLayers(mxd): if lyr.isBroken: print("丢失图层:", lyr.name) # 第一个参数是数据所在新目录,第二个参数是workspace类型 # SHAPEFILE对应shapefile,栅格请换成 RASTER lyr.replaceDataSource( r"E:\huangtu\shape", # 新数据所在目录 "SHAPEFILE", # 数据类型 lyr.name + ".shp", # 数据文件名 False # 不区分大小写 ) mxd.save() print("修复完成")逻辑很简单:ListLayers遍历mxd中全部图层,isBroken判断是否有感叹号,replaceDataSource只改路径,不重新加载符号。参数里最容易错的是第二项workspace类型,英文必须用SHAPEFILE、RASTER或FILEGDB这几个固定值,写错会直接报Invalid workspace type。如果你的mxd里还有TIF图层,单独加一个分支,把workspace类型换成RASTER。如果你用的是ArcGIS Pro,mxd会通过导入工程的方式转换,但接口变成arcpy.mp,思路一致。
3.2 图层顺序和透明叠加:让地形和分区叠而不乱
ArcMap默认新拖入的图层在最上层。合理的顺序是:DEM栅格放最底层,山体阴影或高程渲染层放中间,生态分区面放最上面。分区面的填充色尽量用“无填充”或低透明度,只保留边界线,用带对比度的线色勾勒,这样既能看清4大类分区的界线,又不会把地形细节整个盖住。
具体做法:右键分区图层,属性→显示,透明度根据DEM色带深度调,通常20%到40%之间合适。如果要保留分区内部属性的视觉区分,用“按属性渲染”里的“唯一值”,给4大类各配一个柔和色,透明度统一设置。资源包里的mxd如果已经调好了这套顺序,不要随便拖拽图层,否则后续导出成图TIF会乱。
我一般会在DEM上叠加一层山体阴影(Hillshade),再在阴影层上放一个带高程配色的DEM,把阴影设成按0到100拉伸的灰度,透明度放60%左右,就做出等高线设色加阴影的地形效果。这一步看着是渲染,实际影响出图TIF的信息量,做生态分区叠加之前非常值。
3.3 坐标系不一致时怎么投影:先查再定义
“图能打开但位置不对”是空间数据最常见的翻车现场。黄土高原区域数据常出现三种情况:原始DEM用WGS84经纬度,生态分区用CGCS2000高斯投影;或者两个数据都用WGS84但单位不同;还有的DEM压根没写坐标系。ArcMap的动态投影让显示上看起来是套在一起的,但一旦做裁剪、统计这类空间分析,动态投影不参与计算,结果仍然按各自的原始坐标系算。
统一坐标系的稳妥做法:以生态分区shape的坐标系为准,因为专题矢量数据的坐标系口径通常更权威。先右键分区图层查空间参考,把它的prj文件导出保存;再把DEM投影过去。
import arcpy # 1. 如果DEM缺少投影信息,先用DefineProjection补上 arcpy.management.DefineProjection( in_dataset=r"E:\huangtu\dem_raw.tif", coor_system=r"E:\huangtu\WGS1984.prj" ) # 2. 投影到分区shape所在坐标系,高程连续变量用BILINEAR arcpy.management.ProjectRaster( in_raster=r"E:\huangtu\dem_raw.tif", out_raster=r"E:\huangtu\dem_proj.tif", out_coor_system=r"E:\huangtu\CGCS2000_Albers.prj", resampling_type="BILINEAR", cell_size="30 30" )这里有个先后问题:DefineProjection只是给栅格补一个投影元数据,不做任何重采样。如果原始数据是有坐标系的,你乱定义了一个“正确但不匹配”的坐标系,后面所有结果比“未定义”还难救。所以第一步要先确认原始数据到底是什么坐标系,最常见的检查方式是用GDAL或GlobalMapper看栅格头文件。
ProjectRaster把栅格从原坐标系重采样到目标坐标系,高程这类连续变量用BILINEAR双线性内插,土地利用分类这种离散变量必须用NEAREST最邻近。cell_size参数决定输出分辨率,如果你手里的DEM是30米,就写“30 30”;如果是第三方下载的5米级数据,就得改成“5 5”,照抄会白白损失精度。
4. 用分区面去裁剪DEM:Clip、按掩膜提取和从DEM提shp的实操
生态分区面图层最常见的用途就是裁剪DEM。ArcToolbox里有两个长得像、名字像、结果也像的工具:Clip和Extract by Mask,很多人在这里反复试错,主要因为不知道它们的底层差异。
4.1 Clip和Extract by Mask到底差在哪
“用面图层裁剪DEM栅格TIF”对应的是数据管理工具集里的Clip;“用面图层作为掩膜提取”对应的是Spatial Analyst里的Extract by Mask,中文界面下叫“按掩膜提取”,也有人叫“掩码提取”。两个工具结果看起来都像“把研究区外的DEM切掉”,但处理逻辑不一样。
Clip是把面要素当成一把剪刀,按面几何直接切栅格。注意ArcMap的Clip默认输出范围是面图层的整个外接矩形,只是把面外的像元设为NoData;如果勾选了“使用输入要素裁剪几何”选项,输出才会裁成不规则形状。Extract by Mask则是先把面转成掩膜栅格,再判断每个像元中心是否落在面内,输出保留的仍然是原DEM的像元网格,只是面外全部设成NoData。
实际影响主要在边界和网格对齐上。Clip的像元边缘可能和原DEM网格不完全对齐;Extract by Mask如果不在环境设置里指定snapRaster,输出网格也可能和原DEM有偏移,边缘像元会出现瑕疵。两种工具的适用场景,我用一张表列出来。
| 对比项 | Clip(数据管理) | Extract by Mask(Spatial Analyst) |
|---|---|---|
| 边界对齐 | 默认按面外接矩形,勾选裁剪几何后沿面边界 | 保持原像元网格,面外像元置NoData |
| 关键环境设置 | 输出范围可手动指定 | snapRaster必须与原始DEM一致 |
| 适用场景 | 快速切出非规则范围,出图前裁边 | 做坡度、分区统计前的标准栅格提取 |
| 输出NoData | 可自定义NoData值 | 继承输入栅格的NoData设定 |
我的习惯是:做分区统计和坡度分析前用Extract by Mask,并且强制把snapRaster设成原始DEM;只是出图裁个范围,直接用Clip,速度快、参数少。
4.2 掩膜提取落地:环境设置和输出位深
用arcpy实现按掩膜提取,代码本身很短,坑全在环境设置里。
import arcpy from arcpy.sa import * arcpy.env.workspace = r"E:\huangtu\result" arcpy.env.snapRaster = r"E:\huangtu\dem_proj.tif" arcpy.env.extent = r"E:\huangtu\shape\eco_zones.shp" # 按掩膜提取:生态分区面作为掩膜 out_mask = ExtractByMask( in_raster=r"E:\huangtu\dem_proj.tif", in_mask_data=r"E:\huangtu\shape\eco_zones.shp" ) out_mask.save(r"E:\huangtu\result\dem_by_mask.tif")env.snapRaster的作用是让输出像元网格和原始DEM严格对齐,避免边缘出现半像元缝隙。env.extent如果不设,脚本会按掩膜面范围计算,但有些版本会默认按输入栅格全范围算,输出文件会大得离谱。save之前结果只存在内存里,脚本里必须显式调用save才会写盘。
如果你只想保留4大类分区里的某一类地形,比如只提取风沙区的DEM,常见做法是先按属性筛选面要素,再拿筛选结果做掩膜,不要用整张分区图裁完再手动删。筛选用arcpy.analysis.Select或属性查询都行,关键是掩膜数据源的几何要素要落在目标区域内。
Clip版本的代码更短,但参数细节不同:
arcpy.management.Clip( in_raster=r"E:\huangtu\dem_proj.tif", in_template_dataset=r"E:\huangtu\shape\eco_zones.shp", out_raster=r"E:\huangtu\result\dem_clip.tif", no_data_value="255" )这里的no_data_value是输出NoData的设定值,默认是“保留输入NoData”,但8位无符号栅格常见把255当NoData。如果你的输入DEM是浮点型,NoData是-3.4e38那种数字,写255会让输出栅格凭空多出一圈异常值,正确做法是留空,让工具继承输入值。
4.3 从DEM提等值线shp:等高距、基准高程和NoData陷阱
有人会问“怎么从DEM里提取shp”,最常用的就是把DEM转成等高线shp,也就是ArcToolbox里的“等值线”工具。这个工具在3D Analyst工具集下,ArcMap里叫Raster Surface→Contour。资源包里如果有已经裁好的DEM,直接拿来提取是最顺的路径。
import arcpy arcpy.management.Contour( in_raster=r"E:\huangtu\result\dem_by_mask.tif", out_polyline_features=r"E:\huangtu\result\contour_100m.shp", contour_interval=100, base_contour=0, z_factor=1.0, contour_field_name="ELEV" )contour_interval是等高距。黄土高原地形起伏大,100米等高距出的图可读性较好,图面简洁;想表现沟壑细节可以改成20米,但输出线要素量会成倍上涨,打开和渲染都会变慢。base_contour是基准等高线,设0表示从能被等高距整除的第一个高程值开始绘制,不是从0米开始画。比如区域实际海拔在800到2000米之间,100米等高距下第一根线会落在800米而不是0米,这是正常的。
z_factor是垂直比例因子,数据单位是米就写1,如果原始DEM单位是英尺或厘米,必须换算,否则等高线高程整体偏掉。contour_field_name是写入属性表的字段名,后面做等高线标注直接引这一个字段即可。
一个必须提醒的坑:如果裁剪后的DEM存在NoData缝隙,Contour会在缝隙边缘生成异常短线,甚至拉出横穿NoData的直线。提取前先检查有没有NoData区域,有的话用IsNull和Con工具把NoData填成一个比区域最小高程还低的值,或者干脆再用一次掩膜提取把范围卡死。这个坑我早年踩过,一条等值线白提取几百米。
5. 避坑与常见问题:路径失效、坐标系漂移、NoData黑边
这一章把处理黄土高原DEM和生态分区数据时最常见的5个问题按“现象→原因→解决”写清楚。每一条都是真实复现过的,对照着查能省半天时间。
5.1 换电脑打开mxd全是红色感叹号
现象:mxd文件拷到另一台电脑,打开后每个图层前面一个红色感叹号,地图窗口一片空白,图层右键属性也看不到数据源。
原因:mxd保存的是绝对路径。从压缩包解压到不同的盘符或目录后,原路径指向的文件不存在了,图层自然处于Broken状态。这是mxd机制决定的,不是文件损坏。
解决:按3.1的repairDataSource方法修复。更快的方式是在ArcMap里右键图层→数据源→更换数据源,手动指向新目录,但图层多的时候太慢。从那以后我拿到mxd的第一件事就是先看数据框属性里的路径存储方式,如果数据包整体解压后目录结构不变,可以设成相对路径,换机器时一次改对。
5.2 图层叠加错位几十米:坐标系漂移
现象:DEM和生态分区面在视图中叠加后,边界错开几十米甚至几百米,放大后更明显,做裁剪时输出范围也对不上。
原因:两个图层坐标系不一致。ArcMap的动态投影只负责显示时把不同坐标系换算到同一画面,所有空间分析仍然基于各层自身坐标系计算。经纬度数据和高斯投影数据混用,错位是必然的。
解决:右键图层→属性→源选项卡,把两个图层的空间参考都读出来。确认哪一层坐标系是权威的,一般生态分区shape更可靠。然后用ProjectRaster把DEM投影到分区所在坐标系。重点提醒:不要用DefineProjection去“猜坐标系”,来源不明的DEM先用GDAL或GlobalMapper查看头文件里的投影信息,确认无误再定义。
5.3 裁剪完一圈黑边:NoData没识别
现象:Clip或Extract by Mask之后,输出TIF在ArcMap里显示大块黑色区域,或者研究区外围不是白色透明而是纯黑,出图TIF上很难看。
原因:输入DEM的NoData值和输出设定不一致。浮点型DEM的NoData通常是-3.4028235e+38,但Clip默认输出可能把NoData写成255;8位栅格里255被当作有效值显示成黑色,于是一圈黑边就出现了。
解决:Clip时不要乱填no_data_value,让工具继承输入的NoData设定;如果已经生成了黑边,可以在栅格图层属性里重新指定NoData值,或者用“条件函数”把小于某个阈值的像元直接置NoData。导出成图TIF前,用环境设置里的“将NoData显示为透明”检查一遍。
5.4 shape只拷了一个文件:四件套缺一不可
现象:把“生态分区.shp”单独拷给队友,对方打开后图形变形、属性表空白、投影全偏。
原因:shapefile是复合格式,.shp存几何,.shx存索引,.dbf存属性,.prj存投影。四个文件是一个整体,只拷一个等于把其他三个都留在原机器上。
解决:拷贝shape时整文件夹复制,或者用ArcMap目录窗口的复制粘贴功能,它会自动带全所有伴生文件。检查资源包时,看每个shape是否配齐了.shp、.shx、.dbf、.prj;缺.prj的优先用原始mxd修复后重新导出,不要手工拼一个坐标系进去。
5.5 导出成图TIF发灰发糊:导出参数的锅
现象:mxd里看的好好的,导出成TIF后颜色发灰,图例文字边缘锯齿,放大了模糊。
原因:ArcMap导出对话框里几个参数没设对。8位颜色深度在色带上本来就容易丢失渐变;分辨率低于150dpi时文字边缘会明显发虚;还有的看图软件不识别ArcMap写入的ICC颜色配置文件,导致颜色整体变灰。
解决:导出时分辨率至少设300dpi,颜色模式选RGB,压缩选LZW,背景色设白色;尝试取消“将颜色配置文件写入文档”选项。我自己更常用的路径是先从mxd输出PDF,再用PDF转TIF,因为PDF保留矢量,图例和标注不会变糊,转出来的TIF交付效果明显好于直接导出。
6. 把4大类分区用起来:分区统计高差,再按分区批量出图
6.1 分区统计高差:Zonal Statistics as Table
最实用的产出:按4大类分区统计DEM的高程特征,输出一张属性表,直接贴进报告。工具是Spatial Analyst下的Zonal Statistics as Table。
import arcpy from arcpy.sa import * arcpy.env.snapRaster = r"E:\huangtu\result\dem_by_mask.tif" tbl = ZonalStatisticsAsTable( in_zone_data=r"E:\huangtu\shape\eco_zones.shp", zone_field="TYPE", in_value_raster=r"E:\huangtu\result\dem_by_mask.tif", out_table=r"E:\huangtu\result\zonal_stats.dbf", statistics_type="ALL", ignore_nodata="DATA" )zone_field要填分区shape属性表里表示4大类的那个字段,代码里写的是TYPE,如果你的资源包里字段叫CLASS或NAME,替换一下就行。statistics_type设ALL,会输出MIN、MAX、MEAN、STD、SUM等字段,分区高差(MAX减MIN)、平均高程、地形离散程度一次算齐。输出是dbf,Excel能直接打开,但字段名如果是中文可能乱码,所以前面才建议保留英文短字段做主键。
6.2 按TYPE循环出图:批量给每个生态分区出图
如果报告要求每个生态分区单独配一张图,手动切4次再出图太慢。用definitionQuery循环切换显示范围,再逐个导出,是效率最高的做法。
import arcpy mxd = arcpy.mapping.MapDocument(r"E:\huangtu\huangtu_eco.mxd") df = arcpy.mapping.ListDataFrames(mxd)[0] lyr = arcpy.mapping.ListLayers(mxd, "生态分区", df)[0] # 读取分区字段的所有唯一值 codes = [r[0] for r in arcpy.da.SearchCursor(lyr, "TYPE")] for code in codes: lyr.definitionQuery = "TYPE = '{}'".format(code) df.zoomToSelectedFeatures() arcpy.mapping.ExportToJPEG( mxd, r"E:\huangtu\result\zone_{}.jpg".format(code), df, 150 ) lyr.definitionQuery = ""SearchCursor自动把4大类分区的TYPE字段值全部读出来,比手动写死4个名字通用得多。ExportToJPEG的最后一个参数是分辨率,150表示150dpi;如果交付需要TIF格式,换成ExportToTIFF,参数结构基本一样。注意:脚本只负责切换和导出,图层符号、图例、比例尺这些必须提前在mxd里调好,脚本不会帮你美化。
6.3 一个我固定的收尾习惯
每次拿到这类数据包,我都强制走一遍“三查三算”:查mxd路径,查shape的prj,查DEM的NoData;算范围、算像元大小、算分区面积。三项检查都通过,再做后续分析和出图。以前我吃过一次亏,拿到一份DEM没查NoData就直接提等值线,边界处拉出一条贯穿NoData的长直线,整条等值线白提。从那以后,我每次做任何栅格处理前都强制检查环境设置里的snapRaster和NoData,再把结果保存为带NoData信息的TIF,确认无误才交给下一步。拿到这份黄土高原DEM和生态分区数据包,建议你也把这一遍流程走完再动手,希望帮到你。
本文还有配套的精品资源,点击获取