ARTICLE DETAIL

资讯详情

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

青海省30米DEM处理全流程:从数据源选型到地形分析实战

青海省30米DEM处理全流程:从数据源选型到地形分析实战 简介青海省30米分辨率DEM数据基于ASTER GDEM V3生成面向GIS从业者、地理科研人员及环境规划相关学习者可用于地形特征提取、坡度坡向分析、流域与地质灾害研究等场景。压缩包共10个文件核心为GeoTIFF格式的DEM栅格配套坐标参考、元数据及属性表文件另含青海省边界矢量数据Shapefile涵盖shp、dbf、shx、sbn、sbx、prj等类型便于划定研究区并进行叠加分析。资源包约944.59MB精度适合省级尺度的地形研究WGS84坐标系可兼容ArcGIS、QGIS、GlobalMapper等主流软件。已有494人学习下载。使用者可结合遥感影像分析地形对植被、水文的影响也可服务于城市规划、灾害风险评估与交通路线规划等实际项目为青藏高原环境变化研究提供基础数据支撑。1. 青海省30米DEM高原工程的第一张地形底图在青海做光伏选址、输电线路路径规划或者矿山复垦方案第一个要喂给建模软件的不是卫星图而是DEM。这份30米分辨率的青海省数字高程模型覆盖全省约72万平方公里的连续高程信息从祁连山北麓到唐古拉山口都有完整的栅格表达。它解决的问题非常直接哪片坡地超过20度不适合布置场坪、哪段河谷的填挖方量会超出概算、哪条山脊线两侧的高差会影响线路弧垂设计。适合做省级规划、厂址初筛、水文分析预判这一类前期工作。如果你要出施工图级别的精细设计单靠它还撑不住需要再叠加密测量数据但它一定是你在青海落项目的第一张底图。2. DEM数据源与选型SRTM、ASTER、ALOS在青海怎么选2.1 三种公开30米数据源在青海高原的差异青海省的地形特征决定了没有哪一款DEM是绝对安全的。省内既有祁连山的深切割地形又有柴达木盆地的大面积平坦盐碱地还有长江源头的冰碛湖群和冻胀丘不同数据源的获取方式在这些地物上的表现差异很大。SRTM是2000年航天飞机雷达干涉测量的产物C波段雷达信号在柴达木盆地这样的平坦干沙区干涉质量稳定在祁连山区的峡谷里却容易产生雷达阴影和叠掩表现为山体陡坡上的条带状异常值。ASTER GDEM是光学立体像对匹配生成的在青海湖周边雪线以上、裸岩区、盐湖反光区会出现成片的匹配失败典型症状是地形表面有成片的凹陷坑和尖锥凸起局部高差能差到三五十米。ALOS AW3D30是日本JAXA发布的L波段SAR数据波长更长对冰雪和干沙的穿透性更好在高原地区的空洞率明显低于前两者但公开版本在部分区域的平坦地形上会看到轻微的水印状痕迹。我拿到一份青海省DEM时第一件事就是看它的来源字段。GDAL命令行里一条指令就能查gdalinfo -proj4 qinghai_dem.tif注意看Metadata里是否有SRTMGL3、ASTGTMV3、AW3D30这样的标识。如果没有来源信息用剖面工具在已知高程点上拉两条线一条穿祁连山陡坡一条穿柴达木盆地对比ASTER和SRTM在陡坡处的数值跳动幅度。实际操作中ALOS在青海的山区表现最稳定SRTM在盆地和平原区最干净ASTER只在局部小范围内可用。如果条件允许推荐用ALOS或SRTM作为基础用ASTER只做局部插补。2.2 坐标系与投影跨带项目怎么定参数青海省东西跨度约13.5度经度从东经89度35分到103度04分。按高斯-克吕格3度带划分西端属于第30带中央经线90°E东端湟水谷地属于第34带中央经线102°E。如果项目区正好跨越中央经线比如在格尔木到都兰一带做输变电线路直接用其中某一带投影会导致远离中央经线的一侧边长变形迅速累积。省级尺度的分析我一般用Albers等积圆锥投影双标准纬线设在36°N和38°N中央经线设在96°E。这样全省范围内面积量算不会出现系统偏差。如果项目区是单县域或单流域则用所在区域的3度带高斯投影保证局部形变精度。原始DEM数据通常是WGS84经纬度坐标计算坡度坡向或面积之前必须先重投影到米制坐标系。这里有一个经常被忽略的细节30米分辨率在赤道附近对应约0.00027度但青海处于北纬31度到39度之间同样0.00027度的经度间隔在39度纬度上的地面距离只有约23米。如果直接用经纬度栅格计算坡度又不做比例修正得到的坡度值普遍偏小。要统一到30米×30米的真实地面网格必须用重投影和重采样解决。2.3 拿到数据先做三项检查第一项检查是边界完整性。用gdalinfo查看栅格的四至范围与青海省的省界矢量叠加确认没有缺角或者大范围空洞区。三十米分辨率的全省数据覆盖面积大经常出现某块分幅瓦片缺失的情况。第二项检查是空值分布。用QGIS打开数据用无数据值渲染找出空洞集中的区域。青海的冰川覆盖区、盐湖周边、高山峡谷带往往是空洞高发区。统计空洞面积占全省的比例如果超过5%后续需要做填充或重新找源。第三项是剖面验证。选择几个高程已知的检查点比如西宁市区约2260米、格尔木市区约2780米、青海湖面约3196米在QGIS里用Profile Tool拉剖面看DEM读出的数值与实际高程是否吻合。这一步能快速判断数据的高程基准是不是标准产品也能发现是否存在整体偏移。做完这三项检查数据能不能投入生产心里就有了底再进入预处理流程。3. 预处理实操拼接、裁剪、空洞填补与重采样3.1 用GDAL批量拼接并按边界裁剪青海省DEM如果是从分幅下载源获取的每个文件是一度或半度的一个小格子。推荐用GDAL的VRT机制先拼接不实际合并文件让磁盘压力和内存占用降下来。VRT是一个虚栅格文件记录各分幅的路径和位置关系后续任何处理都可以直接把它当作单一文件对待。# 把青海省目录下所有tif分幅按文件名字符串排列构建一个VRT gdalbuildvrt qinghai_raw.vrt ./青海省分幅/*.tif # 按省界矢量裁剪同时完成到CGCS2000投影坐标系的转换 gdalwarp -cutline qinghai_boundary.shp -crop_to_cutline \ -t_srs EPSG:XXXX -tr 30 30 -r bilinear \ qinghai_raw.vrt qinghai_dem_proj.tif构建VRT这一步如果各分幅的分辨率、坐标系不统一gdalbuildvrt会报错或拼出畸形网格。最常见的分幅下载数据都是WGS84经纬度坐标分辨率也一致这一步通常不会出问题。如果某几个分幅是从不同渠道拿的先单独用gdalinfo核对不要把坐标系混乱的分幅混进同一批。gdalwarp的参数里-cutline指定省界矢量-crop_to_cutline是让输出栅格范围严格贴合边界边界外的像元不输出。-t_srs后面的EPSG:XXXX需要替换成目标投影代码比如项目区所在的3度带高斯投影代码或者省级Albers投影代码。-tr 30 30表示输出分辨率是30米×30米-r bilinear表示重采样方法用双线性内插适合地形数据。如果对地形细节有更高要求可以用cubic三次卷积但计算量会大不少且在空洞区域容易产生过冲波纹。3.2 高原空洞区的填补策略青海的空洞区域有明显聚集性冰川作用区的陡峭岩壁、盐湖湖面、宽河谷的河漫滩。这些区域在SRTM和ASTER中经常是NoData。直接留空会影响后续坡度计算和水文分析因为填洼算法会把空洞当作绝对深坑。用GDAL自带的FillNodata算法处理是最快的方式from osgeo import gdal src_ds gdal.Open(qinghai_dem_proj.tif, gdal.GA_Update) band src_ds.GetRasterBand(1) # 先确保NoData值被正确识别SRTM的NoData通常是-32768 band.SetNoDataValue(-32768) # maxSearchDist200搜索半径200个像素也就是6公里smoothingIterations0不做额外平滑 gdal.FillNodata(band, None, maxSearchDist200, smoothingIterations0) # 关闭文件使写入生效 src_ds None print(空洞填补完成)maxSearchDist的取值直接决定填补质量。设小了大面积的空洞中间部分填不进去设大了填补结果趋向于一个局部均值地形细节被抹平。200个像素对于30米数据来说是一个折中值既覆盖了省内绝大多数空洞尺度又不至于把雅丹地貌的纹理都抹掉。smoothingIterations参数我保持为0因为填补算法本身已经做了一次插值再加平滑会拖垮边缘的锐度。填补完之后必须复查让QGIS把填补区域的边界渲染出来确认填补值没有在湖盆、冰川槽谷里留下突兀的平台或尖锥。如果发现高原面上的填补结果呈现锅盖状说明搜索半径偏大改小一半再跑一次。3.3 重投影与分辨率重采样参数选择依据预处理流程中重投影和重采样可能前后各出现一次。第一次是把经纬度坐标的原始数据转成米制投影供坡度、坡向计算使用第二次是在数据拼接完成后按项目要求把分辨率调整到目标值比如省级分析统一到50米县级工程加密到10米内插。# 从30米降到50米分辨率用于省级景观格局分析 gdal_translate -outsize 50 50 \ -r average \ qinghai_dem_proj.tif qinghai_dem_50m.tif注意这里用的是average重采样而不是bilinear。降分辨率时用平均值聚合像元能最大限度保留地形体积的真实性如果用bilinear会把山脊和沟谷的极值拉低导致后续坡度分析失真。反过来如果因为项目需要从30米加密到15米用cubic三次卷积插值不要用bilinear后者在加密时会产生阶梯状伪影。重采样之后观察直方图。正常地形的DEM直方图应该接近单峰右偏分布峰值在高原平均海拔附近。如果重采样后直方图出现明显的双峰说明插值方法把大片低洼区顶起来了需要检查原始数据中是否存在系统性异常斑块。4. 地形分析实战坡度、坡向、山体阴影与等高线参数4.1 坡度与坡向提取GDAL自带模块的参数细节坡度坡向计算是整个DEM应用里用得最多的环节。GDAL提供了专门的DEM工具集一条命令就能出结果。关键是要选对算法和处理参数。# 在投影坐标系米制下计算坡度输出单位为度 gdaldem slope qinghai_dem_proj.tif qinghai_slope.tif # 计算坡向输出为0-360度正北为0顺时针 gdaldem aspect qinghai_dem_proj.tif qinghai_aspect.tifgdaldem slope默认使用Horn算法该算法在3×3窗口内对中心像元的八个邻域做加权差分。Horn算法对噪声相对保守适合祁连山地这种地形起伏大、且原始数据有残余噪声的场景。如果你的数据来源是平滑的水准测量成果可以改用-ZevenbergenThorne算法它对表面细节更敏感但对单个像元的异常值几乎没有抵抗力。在青海的应用场景中坡度的分级阈值通常和工程规范挂钩光伏场地要求坡度通常小于15度到20度输电线路塔位要求小于30度泥石流沟道识别则关注大于35度的陡坡段。导出坡度栅格后用gdalwarp做一个条件分类按阈值把可建设区域直接矢量化导出。# 提取坡度小于15度的区域 gdal_calc.py -A qinghai_slope.tif --outfileslope_lt15.tif \ --calcA15 --NoDataValue0gdal_calc.py是GDAL自带的栅格计算器A代表输入的坡度栅格计算结果中True自动变为1False变为0。后面接的NoDataValue0把所有不可建区置为NoData方便后续转为矢量。4.2 山体阴影制图太阳方位角和高度的设定逻辑山体阴影是DEM最直观的展示形式也是纸质图件最常用的地形基底。GDAL的hillshade模块参数设定直接影响图面可读性。gdaldem hillshade qinghai_dem_proj.tif qinghai_hillshade.tif \ -z 3.0 -az 315 -alt 45-z 3.0是垂直放大系数。青海高原面大体平坦绝大多数区域坡度不超过5度直接生成的阴影整体偏灰、层次弱放大三倍后地形脉络才清晰。-az 315表示光源方位角315度即从西北方向打光阴影投向东南。这个方向在制图规范中最符合人眼的读图习惯。-alt 45是太阳高度角45度在高原地区既不会让阴影过浓也不会让地形过于扁平。山体阴影文件通常还需要与真实色彩或高程分层设色合成使用。做法是用QGIS的栅格混合模式把阴影栅格作为Multiply图层垫底DEM的颜色渲染作为Normal图层叠在上面透明度调整到70%到80%。这种阴影浮雕效果在一张图里既能表达绝对高程又能表达地形的起伏纹理。4.3 等高线生成与抽稀制图比例尺的匹配逻辑从DEM提取等高线有两种工具一条命令解决# 生成间距100米的等高线高程值写入属性elev gdal_contour -a elev -i 100 qinghai_dem_proj.tif qinghai_contour_100m.shp-i 100的等高线间距选择跟目标图件的比例尺直接相关。省级挂图用200米间距州县级工作底图用100米县级工程设计用10米到20米。间距越小生成的Shapefile文件越大线条越密。青海省地形高差大如果强行在工程图里用5米间距提取等高线峡谷地段会出现线线粘连图面一塌糊涂。等高线提取后通常会做一次抽稀和平滑因为DEM栅格是离散的直接提线在陡坡处会呈现锯齿。用QGIS的Simplify工具容差设为30米算法选Douglas-Peucker保留山脊线的基本骨架。注意抽稀时不要把闭合等高线在鞍部截断否则后期填色会产生错误。5. 避坑指南青海DEM处理中的五个真实翻车现场5.1 湖面高程大面积异常现象在青海湖、扎陵湖、鄂陵湖的湖面范围内DEM高程值不是平滑的平面而是呈现密集的随机凹凸局部高差达20米以上湖岸线处甚至出现夸张的断层。原因雷达和光学数据在水体区域都会失效。水面几乎没有雷达回波立体像对匹配也会因为纹理缺失而失败数据生产方用内插补出的湖面高程自然不可信。解决用已有的湖泊水面边界矢量把湖面范围单独提取出来直接赋统一高程值。青海湖实测湖面高程约为3196米具体取值以最新水利普查或实测数据为准。先用gdal_rasterize把湖泊矢量栅格化再用gdal_calc.py把该范围的高程替换成常量输出新的DEM文件作为工作版本。5.2 深切峡谷区的条带状噪声现象玉树、果洛一带的深切峡谷坡面上DEM高程呈条带状起伏条带方向与山体走向一致每条带宽约100到300米像搓衣板纹理。原因这是雷达阴影效应的典型表现。SRTM的C波段在陡峭峡谷里接收不到有效回波数据生产方用周边像元插补时留下了方向性痕迹。ASTER在同样位置的立体匹配也经常失败表现为横向撕裂。解决先用坡度数据把噪声区识别出来坡度大于40度的区域作掩膜然后用中值滤波窗口5×5或7×7只对掩膜区域做平滑。这样能压制条带又不会破坏缓坡地的真实地形。如果平滑后依然明显换ALOS AW3D30重新提取该区域数据做局部镶嵌。5.3 与GPS实测高程偏差达到20米以上现象拿着RTK实测的工程控制点高程与DEM读取的高程对比发现系统性偏低或偏高15到25米且误差方向一致。原因高程基准不一致。DEM产品通常采用EGM96或EGM2008大地水准面模型高程值是基于大地水准面的正高而国内工程测量使用1985国家高程基准两者在青藏高原地区的差异可以达到几十米。另有部分商业DEM使用WGS84椭球高那是把地形表达成相对椭球面的距离与正常高之间隔着一个大地水准面差距。解决首先查数据的官方元数据确认高程基准类型。如果是EGM96则需要用高精度似大地水准面模型如省级似大地水准面精化成果做格网改正把正高转为1985正常高。快速修正办法是在项目区均匀采集5到10个已知高程点计算平均差AfterDEM和实测差再用gdal_calc对整个DEM做常数平移。这个办法适用于局部工程但不适合全省统一成果。5.4 分幅接边处出现错位跳变现象拼接后的DEM图面上沿着原分幅边界有一条清晰的高程台阶山脊线在边界两侧错开20到50米但每幅内部看起来都正常。原因相邻分幅来自不同版本的数据源或同一版本但不同期处理几何配准和垂直基准有细微差异。尤其常见于跨省边界青海和西藏、新疆的交界分幅来自不同数据中心。解决如果所有分幅是同一个源同一版本先检查各分幅的仿射变换参数是否完全一致不一致的用gdal_translate修正到统一网格。不同版本混拼的最简单方案是不拼接、各幅单独处理后再用加权融合。GDAL有个VRT的edge羽化功能可以在VRT构建时加blend-dist参数让接边在数百米范围内渐晕过渡视觉效果会好很多但严格定量分析仍然不推荐。5.5 重采样后地形细节消失现象把30米DEM重采样到90米后原本清晰的山脊线变钝小型冲沟消失计算出的流域面积和实际对不上。原因重采样方法选错。降分辨率时用bilinear相当于对局部地形做了平均山脊和沟谷这些高频信息被当成噪声滤掉了。直接后果就是水文分析结果偏差大源头的细小河道会被合并甚至消失。解决降分辨率用average聚合它先把30米像元划分成3×3块再对每个块求均值能保留总量不丢。如果特别在意地形骨架更好的做法是用一个名为terrain position index的中间层先提取地形特征再对特征层做聚合重采样。从那以后我每次做重采样都会在输出文件命名里强加一个后缀标明采样方式避免半年后自己都不记得当时用了什么参数。6. 进阶玩法从DEM到水文分析与三维场景6.1 河网提取填洼阈值怎么定才不翻车水文分析是DEM应用里最容易出玄学问题的模块。标准流程是填洼、算流向、算汇流累积、按阈值提取河网。填洼这一步的阈值选择直接决定提取出的河网密度和形态。# 使用WhiteboxTools填洼修复所有凹陷 whitebox_tools --runFillsDepressions -idem_proj.tif -odem_filled.tif # 计算D8流向 whitebox_tools --runD8Pointer -idem_filled.tif -od8.tif # 计算汇流累积量 whitebox_tools --runD8FlowAccumulation -id8.tif -oacc.tif在青海的高原面上天然浅洼地特别多比如长江源区的草甸冻胀坑、湖盆边缘的浅水塘。完全填洼会把所有洼地都填平后续提取的河网会变成沈阳棋盘格一样的平行线。对于这类区域我不改默认填洼而是在汇流累积量这一步把河网提取阈值调大一般用积累量大于5000像元才定义为河道对应约450万平方米的汇水面积。如果提取结果过于稀疏逐步降低阈值到3000每次跑完对比一下与真实河流叠加的吻合度不追求一刀切。6.2 三维场景快速构建QGIS加Qgis2threejs地形汇报场景里一张山体阴影叠加彩色高程的静态图往往不够直观。用QGIS的Qgis2threejs插件可以把DEM直接导出成一个HTML三维场景任何人用浏览器就能打开旋转查看不需要安装专业GIS软件。操作流程QGIS里加载填洼后的DEM设置好高程渲染色带叠加山体阴影图层。打开Qgis2threejs面板将DEM拖入作为高程层设定垂直拉伸系数。青海的地形平坦区域垂直拉伸系数建议4到5五千米级的高山区用2到3否则峡谷会显得过于狰狞或者高原面显得一片平坦。地形纹理选择叠加山体阴影的合成图层导出HTML后即可交付。这个流程不生成磁盘上成百GB的三维模型却能在方案评审现场直接拉视角看线路走廊和场平范围比对着等高线图解释高效得多。我把这套流程整理成了一个固定模板每次新增项目区换DEM、换边界十分钟出场景拿来汇报和复核都够用。希望这些处理习惯能帮到你。本文还有配套的精品资源点击获取
返回列表