ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

临沧市30米DEM数据处理全流程:裁剪、坡度水文分析与可视化避坑

临沧市30米DEM数据处理全流程:裁剪、坡度水文分析与可视化避坑 简介云南省临沧市DEM数字高程数据为30米分辨率、含区域范围Shapefile的GIS数据集面向地理信息初学者与研究人员适用于地形分析、坡度计算、洪水淹没模拟、城市规划等场景。压缩包共12个文件核心为TIFF格式的高程栅格并配套Shapefile矢量边界及投影、索引、属性等文件另有XML元数据说明整体大小约87.69MB。目前已有330人学习/下载。数据覆盖临沧市及周边部分区域每个像元代表30米×30米的地面高程信息。用户可将高程栅格与边界结合进行裁剪、重分类、坡向提取、可视域分析等也可借此熟悉DEM栅格组织、Shapefile文件结构以及投影坐标定义是学习GIS数据格式与地形分析的实用素材。1. 临沧DEM 30m数据包一份能直接开工的地形底图而不是一堆瓦片接手临沧市地形分析项目的时候手头最缺的就是一份现成的数字高程数据。自己从数据平台下载分幅DEM再拼接裁剪折腾一天是常有的事。这份“云南省临沧市DEM数字高程数据30m含区域范围shp文件.zip”不一样解压出来就是一个覆盖临沧市范围的30米分辨率DEM一个tif文件加配套的区域范围shp。30米分辨率意味着每个像元对应地面上约30米×30米做坡度分析、水文分析、选址评价这类市级尺度的活儿完全够用。不再需要自己拼瓦片、对坐标系、找边界适合GIS从业者、规划设计师、水利和林业方向的人也适合需要临沧地形底图的学生项目。2. 拿到zip先做数据体检坐标系、分辨率、空值一个都不能漏做数据工作最忌拿到文件就直接往软件里拖。这个zip包到手我习惯先花十分钟做一遍体检摸清元数据再动手能省掉后面一半的返工。体检的内容无非三件事坐标系是否明确、分辨率是否真实、NoData设置是否合理。2.1 解压后先检查文件清单shp是一组文件不是单一文件先用命令解压文件名带了中文和括号注意引号unzip 云南省临沧市DEM数字高程数据30m含区域范围shp文件.zip -d lincang_dem命令行处理带空格的参数时引号是必须的。-d lincang_dem指定解压到新建的lincang_dem目录避免文件散落一地。解压后先别急着拖进ArcGIS对照文件清单检查一遍。shp格式不是单一文件它至少需要.shp、.shx、.dbf三个核心文件配套外加一个.prj记录坐标系。缺了.shxArcGIS会提示无法打开缺了.prj坐标系会显示未知后面所有投影相关操作都会出问题。文件后缀作用缺失后表现.tifDEM主体栅格数据无数据.shp临沧市范围面要素几何无边界可用.shx要素几何索引打不开shp.dbf属性表属性缺失.prj坐标系定义坐标系显示未知2.2 元数据体检用rasterio确认坐标系和分辨率不能迷信文件名里的“30m”三个字。有些标称30米的数据实际是1弧秒采样换算过来约30.8米有些数据本身是WGS84经纬度坐标需要转投影才能算面积和坡度。用一小段Python确认import rasterio dem_path lincang_dem_30m.tif with rasterio.open(dem_path) as src: # 读取六个关键元数据先确认数据底细再做分析 print(坐标系:, src.crs) print(分辨率:, src.resolution) print(边界范围:, src.bounds) print(无效值:, src.nodata) print(波段数:, src.count) print(宽度x高度:, src.width, src.height)crs决定后续计算能不能直接用地理坐标resolution如果输出(0.0002777778, 0.0002777778)这类度单位说明数据是经纬度格式一个像元大约对应30.8米bounds可以判断范围是否真的覆盖临沧全境nodata标记无效像元的取值后续裁剪和分析识别空值都靠它。如果zip里包含多个分幅tif不想一个个手查可以写个循环批量输出import rasterio import pathlib for tif_path in pathlib.Path().glob(*.tif): with rasterio.open(tif_path) as src: # 批量输出所有tif的元数据快速找出坐标系异常项 print(tif_path.name) print( crs:, src.crs) print( res:, src.resolution) print( nodata:, src.nodata) print( bounds:, src.bounds)注意crs如果显示None说明tif丢失了空间参考信息。这种情况在ArcGIS里加载不报错但后续所有测量都会变成无坐标单位的无效计算需要手动用“定义投影”补回坐标系信息。2.3 在ArcGIS里做一次目视检查脚本只能看数字实际地形细节还得用眼睛过一遍。加载tif后右键图层属性在Source页签查看像元大小和坐标系再用符号系统的拉伸显示如果图像边缘有大片黑色或白色往往是没有正确识别NoData。这里有一个关键区分NoData是不参与分析的无效值0值是真实的海拔高度两者混在一起后面的坡度、填洼分析都会出现虚假的断层和深谷。看到DEM边缘发黑先检查图层属性里的NoData设置别急着用“复制栅格”工具去重写。3. 用shp裁剪DEMClip和掩膜提取差在哪参数怎么设不翻车拿这份数据的人第一步基本都是把DEM裁到临沧市范围内去掉周边区域。但ArcGIS里“裁剪”和“按掩膜提取”两个工具名字相近逻辑完全不同用错的人不在少数。3.1 “裁剪”与“掩膜提取”的实质差别“裁剪”工具位于数据管理工具→栅格→栅格处理→裁剪默认行为是生成一个矩形输出范围等于裁剪面要素的外接矩形。勾选“使用输入要素裁剪几何”之后才会按照面要素的实际形状裁剪。如果不勾裁出来是矩形框框外区域被填为NoData勾选后输出才是贴合边界的栅格。“按掩膜提取”位于Spatial Analyst工具→提取分析→按掩膜提取它直接按掩膜面的范围提取栅格像元不会生成矩形包围盒里的空白区域边缘按照面边界切割。两者对边缘像元的处理也有差别Clip按几何边界切割掩膜提取按像元中心是否落在面内来取舍。在30米分辨率下这个差异肉眼很难察觉但做面积统计和边界制图时会体现出来。对比项裁剪 Clip按掩膜提取 Extract by Mask工具位置数据管理工具→栅格→栅格处理Spatial Analyst→提取分析默认输出范围面要素外接矩形面要素实际范围是否依赖Spatial Analyst许可不需要需要边缘像元取舍几何切割像元中心判定常用场景快速裁出矩形研究区精确按行政区提取如果本机没有Spatial Analyst扩展模块用Clip勾选“使用输入要素裁剪几何”也能达到接近的效果。不依赖ArcGIS图形界面的话用Python的rasterio.mask也能完成同样的操作参数显式写出来出问题更好排查import geopandas as gpd import rasterio from rasterio.mask import mask boundary gpd.read_file(lincang_boundary.shp) boundary boundary.to_crs(EPSG:4326) # 与DEM坐标系对齐 geom boundary.geometry.values[0] with rasterio.open(lincang_dem_30m.tif) as src: out_image, out_transform mask( src, [geom], cropTrue, # 只保留面范围内的栅格 nodata-9999 # 面外像元统一标记为无效值 ) out_meta src.meta.copy() out_meta.update({ driver: GTiff, height: out_image.shape[1], width: out_image.shape[2], transform: out_transform, nodata: -9999 }) with rasterio.open(lincang_dem_cropped.tif, w, **out_meta) as dst: dst.write(out_image)逻辑说明先用geopandas读取shp把几何统一到与DEM一致的坐标系mask函数用临沧市面要素当作掩膜cropTrue只保留面内栅格nodata-9999把面外区域标记为无效防止后续处理把边界外的区域当作海拔0的真值。参数方面如果shp里包含多个面要素建议先用unary_union合并成一个多边形再传入避免几何对象冲突。3.2 先投影还是先裁剪这个顺序问题最容易翻车。如果DEM是WGS84经纬度shp是CGCS2000投影坐标系直接做掩膜提取会提示空间参考不一致或者输出结果偏移半个像元。我一般习惯先用“投影栅格”把DEM转到与shp一致的投影坐标系再做掩膜提取。投影栅格工具里有重采样方法选项最近邻、双线性、三次卷积。做地形分析用双线性即可它在保持平滑的同时不会像三次卷积那样在陡峭山地区域产生过冲高值和低值后者会让坡度分析出现虚假极值。重采样后检查输出像元大小UTM投影下应保持在30米左右如果变成29.9或30.1说明累计误差已出现长距离剖面分析时要留意。顺带说一下投影带临沧大部分区域位于UTM 47N分区具体以项目要求为准不要只凭数据范围去猜。3.3 批量裁剪多个乡镇如果shp里有多个乡镇面要素想每个乡镇单独输出一份DEM手动一个个跑效率太低。可以直接循环裁剪import geopandas as gpd import rasterio from rasterio.mask import mask boundary gpd.read_file(townships.shp) for idx, row in boundary.iterrows(): name row[name] geom row.geometry with rasterio.open(lincang_dem_30m.tif) as src: out_image, out_transform mask(src, [geom], cropTrue, nodata-9999) out_meta src.meta.copy() out_meta.update({ height: out_image.shape[1], width: out_image.shape[2], transform: out_transform }) with rasterio.open(f{name}_dem.tif, w, **out_meta) as dst: dst.write(out_image) print(f{name} 处理完成)逻辑说明iterrows逐行读取shp属性表row.geometry取出每个乡镇的面几何每次循环打开一次DEM并裁剪输出文件名直接用属性表的name字段省掉手动重命名。参数注意name字段如果包含斜杠、星号等非法字符Windows下写入会报错建议先做一次字符串清洗把非法字符替换成下划线再拼文件名。4. 从30米DEM提取地形信息坡度、坡向与河网提取的实用参数DEM不是拿来看颜色的它的价值在于能提取出真正用于决策的地形因子。临沧这类山地城市坡度、坡向决定建设用地选择水文分析决定流域和河网边界。这一章把关键参数讲透。4.1 坡度与坡向单位选择和Z因子坡度工具在Spatial Analyst→表面分析→坡度。输出单位有两种度和百分比。做坡度分级图、建设用地适宜性评价用度工程排水、坡面流计算常用百分比。临沧大部分像元的坡度集中在15度到40度之间如果算出来全区平均坡度超过50度先别怀疑数据检查坐标系和Z因子。Z因子是这地方最容易被忽略的参数。当DEM是经纬度坐标系时x和y方向单位是度z方向单位是米直接算坡度等于把水平距离和垂直距离当同单位处理结果完全失真。解决办法是给Z因子填一个换算系数经纬度DEM大约填0.000009它把高程米数换算到与坐标单位一致的尺度。DEM若是UTM或高斯投影x、y、z单位都是米Z因子保持1即可。场景坐标系Z因子WGS84经纬度DEMEPSG:43260.000009附近UTM投影DEMEPSG:32647等1CGCS2000高斯投影EPSG:4490等1ArcGIS的坡度工具不会自动识别DEM的坐标系并调整Z因子每换一个数据源都要确认一次这也是所有人都会踩一遍的坑。坡向工具相对简单输出0到360度的坡面朝向0度和360度都代表正北-1表示平地。制图时把坡向按8方位重分类做农林选址和建筑朝向分析很实用。4.2 水文分析填洼参数、流向与累积阈值水文分析是流程式工具链填洼→流向→流量累积→栅格计算器→河流链接→河网分级。第一步填洼决定后续所有结果它把DEM里的洼地填充到周边最低溢出口高度保证水流能连续流向流域出口。填洼工具的Z限制参数默认是100米意思是深度超过100米的洼地不会被填。临沧有喀斯特地貌天然溶蚀洼地不少有些真实洼地深度超过50米。Z限制设太大这些真实地形会被填成平地提取的河网偏离实际汇流路径设太小又填不掉DEM噪声造成的伪洼地。我一般先用默认值跑一遍再改用50米限值对比看河网结果差异大不大再决定最终参数不盲目固执于单个数值。流向工具用D8算法把每个像元的水流方向指定到八个相邻像元中坡度最陡的方向。D8算法决定了平坦区域容易提取出平行排列的伪河道这是30米分辨率DEM的通病。流量累积之后栅格计算器提取河道阈值1000代表累积流量超过1000个像元的区域被识别为河道。换算成面积是1000 × 30 × 30 900000平方米约0.9平方公里。想提取集水面积大于1平方公里的河道阈值设1200左右只要主要干流设5000以上。计算器表达式为Con(FlowAcc 1000, 1)Con是条件函数满足条件的像元输出1不满足输出NoDataFlowAcc是流量累积栅格图层的名称路径含中文时注意引号。5. 避坑记录临沧DEM使用中的五个数据坑与排查办法数据本身扎实但使用过程中翻车的情况太多了。梳理五条实际踩过的坑每条按现象、原因、解决的顺序写清楚动手前扫一眼能省不少排查时间。5.1 投影与坐标系的两类典型坑坑一面积统计明显偏小或偏大。 现象裁剪出的临沧DEM计算面积结果和官方公布的行政区面积差很多有时差一个量级。 原因DEM仍是WGS84经纬度坐标面积计算工具没有动态投影统计结果直接在度平面上算单位错误。 解决先用投影栅格工具把DEM转换到UTM或CGCS2000投影坐标系再计算面积单位从度变成米结果才正确。坑二shp和DEM加载后不在一个位置。 现象ArcGIS里同时加载shp和DEM两个图层一个在临沧一个偏移了几十公里甚至到邻国境内。 原因两个文件的坐标系定义不一致。常见是shp为CGCS2000DEM为WGS84两者都是经纬度表示但椭球体参数有差异大范围数据叠加时偏移被放大。 解决先看两个图层的属性定义确认坐标系差再用“投影”工具把shp转到与DEM一致的坐标系。注意不要用“定义投影”去强行改定义投影只改标注不改坐标数值会造成更大偏移。5.2 数据值与形态的三类坑坑三裁剪后图像边缘出现大片黑色或白色区域。 现象掩膜提取后边缘有一圈明显的无效区域坡度分析结果里也多出一圈异常高值。 原因裁剪时没正确指定NoData值或原始DEM本身没定义NoData面外像元被当成真实低海拔数据参与计算。 解决先用rasterio检查nodata设置没有定义就统一设为-9999或-32768裁剪后再用“按属性提取”把该值重新定义为NoData。坑四水文分析在平坝区提取出平行排列的假河道。 现象河网shp在河谷平地、山间坝区出现大量平行短线段和实际水系走向完全不符。 原因30米DEM在平坦区域无法分辨细微坡度差D8流向算法把水流强制分配到固定方向网格产生人为平行流。 解决填洼前对DEM用“焦点统计”做一次3×3窗口平滑再走完整水文工具链提取河网后还要对照遥感影像目视检查。临沧河谷坝区面积不大但坝区边缘特别容易出现这类假河道。坑五坡度分析结果出现刺状突变的假高坡。 现象坡度图上看到很多细碎的极陡区域分布呈刺状和山体走向不吻合。 原因重采样到投影坐标系时选了三次卷积方法在陡峭山谷边缘产生过冲值局部高差被算法夸大。 解决投影栅格重采样时选双线性插值如果已经生成错误结果用焦点统计平滑后重新计算坡度。30米数据本身细节有限不要用强锐化算法人为制造虚假地形特征。6. 收官技巧DEM转等高线叠加晕渲做一张能汇报的地形图分析做完最后要出图。直接扔一张DEM拉伸渲染图不懂GIS的人很难读懂但等高线叠加山体阴影和彩色高程不需要解释也能看出哪里是山、哪里是谷。做法分三步。第一步生成等高线。用Contour工具市域范围出图选100米等距合适因为临沧从河谷到山脊高差经常超过2000米等距设太细图面会密成一片。县区或乡镇尺度的重点区域用20米到50米。起始等高线建议取DEM最小值的整百值让首条线落在整齐的数值上。第二步生成山体阴影。用Hillshade工具方位角默认315度太阳高度角默认45度这套组合在大多数地形图上观感都不错。想突出沟谷细节把高度角降到30到35度阴影会更重山谷走向更清楚。第三步组合图层参数可以按自己的出图偏好微调图层顺序图层类型推荐参数最底层DEM彩色渲染黄绿棕渐变色带拉伸方式选最值中层山体阴影透明度约40%灰度显示最上层等高线间隔100米线宽0.4浅灰色透明度设在40%左右山地明暗变化能透出来又不抢DEM颜色的信息。等高线用浅灰而不是纯白避免在图上喧宾夺主。QGIS里同样可以按这个三层结构叠加效果基本一致。出图前我会再叠一层水系shp做交叉验证确认等高线走向没有和真实河谷冲突。从那以后我每拿到一个新区划的DEM数据包都先强制走一遍元数据检查、投影确认、空值和裁剪四个固定步骤再进入分析流程这套顺序帮我少返工很多次。这份临沧市数据包里包含DEM和shp范围解压后按流程从头跑一遍半小时就能出第一张可汇报的地形图希望帮到你。本文还有配套的精品资源点击获取
返回列表