☰
FathomDEM中国30米分省高程数据:优化原理、裁切流程与水文应用
2026/10/21 4:27:21 网站建设 项目流程

如果你经常倒腾DEM数据,看到“全国分省打包”这几个字,大概能猜出这东西的含金量。过去两年里,我为了给一个流域模拟项目找一套没有条纹噪声的30米高程数据,先后试过SRTM、ALOS、原始哥白尼DEM,每个版本都有让人头大的地方。直到最近用上FathomDEM——这套由哥白尼DEM优化而来的全球高精度高程数据,我才觉得终于能把精力从“反复清洗数据”转移到“真正分析数据”上。这篇内容就把中国30米分省包的底细、整理思路、切分流程和实际使用变化一次性讲清楚,给做地质、水利、规划、遥感和GIS开发的同学做个参考。

另外多说一句,这次分享无套路,不需要先关注公众号、不需要转发三个人、不需要后台回复什么口令。文章里该给的目录结构、裁切命令、容易踩的坑,全部放在明面上。你唯一要做的是拿到数据后先花十分钟确认坐标系统和NoData设置——这两个问题比数据本身更能影响后续分析。

1. FathomDEM到底优化了什么:先搞懂你下载的是哪一版数据

1.1 哥白尼DEM的“血统”与痛点

哥白尼DEM(Copernicus DEM)来自欧盟哥白尼计划的哨兵一号雷达卫星干涉测量,全球覆盖,公开分发,分辨率分为1弧秒(约30米)和3弧秒(约90米)两个档位。相比老一代SRTM,它在空间连续性和全球覆盖完整性上已经强了不少,尤其是纬度覆盖范围和内陆水体的处理,比SRTM要干净。

但原始哥白尼DEM并不完美。雷达干涉测量受地表散射、大气延迟和轨道误差影响,在地形平缓的区域容易出现一种横向条带纹理,业内俗称“条纹噪声”。在平原耕地区域,这种条纹会让地形渲染图看起来像被拉过丝一样,导致等高线出现大量锯齿状的小波动。更麻烦的是,在高山峡谷、陡峭崖壁区域,雷达信号会产生阴影和叠掩,原始的哥白尼DEM在某些山谷里会有明显的高程异常空洞,直接影响到后续水文分析中河网的提取。

我最早用原始Copernicus DEM跑一个小流域的汇水区划分,结果在农田区提取出了一堆平行伪河道,跟实际水系对不上。这个现象在很多GIS论坛里都能看到,反正视觉上就是“平原上的河流像被挠过的头发一样乱糟糟”。

1.2 “优化版”的再加工到底做了什么

FathomDEM这个产品,本质上不是用一颗新卫星重新测了一遍全球,而是对哥白尼DEM进行的深度后处理。简单说,它把原始DEM当成一张“脏图”,用一整套自动化的地形增强流程去把噪声、空洞和细节问题修掉。

具体来说,核心处理包括三块:

  • 雷达斑点与条带噪声抑制。通过自适应滤波和深度学习模型识别条纹、斑点这类非地形信号,把它们从高程场中分离出去。这一步比传统的中值滤波或低通滤波高明在它能保留真正的微地形细节,而不是把所有起伏都抹平。
  • 空洞与异常值修复。针对雷达阴影区、叠掩区、水域边缘的空值,利用周边地形趋势、坡度坡向以及辅助数据完成插补。陡峭山区原来那种“缺一块”的现象会明显减少。
  • 地形细节增强。在去噪的同时对山脊线、山谷线、河流下切等微地形进行适度的锐化处理。说白了,就是让河谷和山脊看起来更“有骨头”,而不是软绵绵的一片。

这套优化由一家长期做海洋与地理空间分析的国际团队完成,后来也被整合进某主流GIS软件的全球底图服务里,可以直接在线预览。不过在线预览始终不方便做本地分析,这也是我整理分省包的原因之一。

1.3 同一片地形的直观差异:原始版与优化版对照

很多读者可能关心,FathomDEM和原始Copernicus DEM放在一起,到底能不能看出区别。我直接给一个对比表格,都是我自己在实际项目中观察到的典型差异:

对比维度原始Copernicus DEMFathomDEM
平缓平原的纹理常见横向条带、斑点噪声纹理干净,等高线平滑
河谷微地形细小沟谷容易被噪声淹没河流下切趋势清晰可见
陡峭山地雷达阴影区域偶见空洞空洞大幅减少,连续性更好
山脊线局部破碎、锯齿感明显山脊连续,形态更接近真实地形
垂直精度观感中低坡度区域残留噪声明显中低坡度区域噪声显著收窄
插值倾向原始测量值,噪声保留经过模型修正,细节更干净

如果你手上有原始版的30米数据,最直观的检验方法是生成一张山体阴影图(Hillshade),在平原和丘陵交界处放大看。原始版会看到很多不该存在的颗粒感,而FathomDEM则更像一张“理过发”的地形,该平滑的地方平滑,该锋利的地方锋利。

需要提醒的是,“优化”不等于“绝对准确”。它属于对原始雷达高程的系统性修正,修正结果在整体趋势上更符合真实地形,但在极陡的悬崖、茂密森林覆盖区、采矿区等场景下,仍会存在偏差。

2. 中国区30米分省包的数据底细:目录、坐标、NoData这些容易被忽略的地方

2.1 目录长什么样,文件怎么命名

这套分省包是按省级行政区划切割整理好的,整体目录结构类似下面这样:

China_FathomDEM_30m/ ├── README.md ├── boundary/ │ ├── china_province.shp │ └── china_province.dbf ├── hb/ │ ├── hb_30m_dem.tif │ └── hb_30m_dem.sha256 ├── ah/ │ ├── ah_30m_dem.tif │ └── ah_30m_dem.sha256 └── xj/ ├── xj_30m_dem.tif └── xj_30m_dem.sha256

文件命名采用省份拼音首字母缩写加分辨率后缀。我用了两个字符的缩写,比如“hb”代表某中部省份,“ah”代表某东部省份,“xj”代表面积较大的西部省份。直辖市也按省级单位单独建目录,方便直接按省级行政单元调用。面积特别大的省份,原始栅格裁出来后文件体积比较大,我做了分块存储,每一个分块仍然是GeoTIFF,可以用GIS软件直接打开。

2.2 坐标系选择与理由

FathomDEM官方分发版常见的坐标系是Web墨卡托(EPSG:3857),因为发布方主要面向Web地图服务场景,直接叠在线底图上非常方便。但对于做地形分析的人来说,直接用3857有个隐患:它的纬度方向被拉伸,越往高纬度变形越大,直接用它计算坡度、坡向会有明显误差。

我整理分省包时做了个决定:统一重投影成WGS84经纬度坐标(EPSG:4326),高程值保持不变。这样做有几个好处:

  • 各省之间的坐标系统一致,跨省拼接时不会出现“一个省是3857、另一个省是4326”的尴尬。
  • 经纬度坐标系是GIS软件通用的默认坐标系,加载时不会出现软件莫名报错。
  • 后续需要做面积、坡度等分析时,再按项目需求投影到对应的本地投影坐标系即可。

当然,这样做也有代价:在QGIS或ArcGIS里打开时,软件会默认按WGS84显示,如果你叠的是在线墨卡托底图,会出现坐标转换的轻微重采样,视觉上略微发虚。不过对分析而言,这个影响不大。真要做坡度和面积统计,我建议在项目里新建一个合适的投影坐标系后再做。

2.3 NoData、高程单位与负值:三个最隐蔽的坑

这三个问题我特别想展开说,因为几乎每个拿到DEM的人都会在这里栽跟头。

高程单位是米,基准是EGM2008大地水准面。这不是椭球高,而是近似海拔高。在中国区域,它与常用的“1985国家高程基准”之间存在一定系统差,不同地区差值方向不完全一样。如果你只是做流域分析、坡度分析、视域分析,这个差值影响不大;但如果要跟实测高程点严格对比,建议先做区域性的基准偏差矫正。

NoData值统一设为-9999。很多人习惯用0作为空值,这在平原地区问题不大,但在一些低海拔洼地,真实地形高程本身就可能低于0米,比如典型的极低海拔内陆盆地、某些干涸湖泊底部。如果用0当空值,等于把这些真实存在的负高程全部吃掉。这一点务必检查。

有些省份的边缘存在真实空洞。我裁切时对部分沿海省份进行了近岸岛屿的完整保留,但极端远海的小岛礁在原始全球数据中本身是缺失的,我在README里做了说明。遇到这些空洞,建议把它视作“无数据区域”而不是“海平面0米”,否则做淹没分析时会得到完全错误的结果。

拿到任何一省数据后,第一件事建议执行一下gdalinfo命令,看看NoData、坐标范围、像素类型是否符合预期。不要偷懒,这一步能帮你省下后面半天排查时间。

2.4 分省包的范围怎么来的

范围边界采用的是公开渠道获取的省级行政区划矢量数据,裁剪时以矢量边界为准。需要明确的是,行政区划数据存在不同的统计口径和发布版本,边界线在细节上可能有几十米到几百米的差异。这个差异对宏观分析影响不大,但如果涉及到具体地类面积统计或边界接边,建议你自己再叠加最新的影像底图校对一次。

对于跨省项目,我的建议是不要只拿一个省的包去做,而是把相邻省份一起加载,拼接后再分析。这也是下一章要讲的内容。

3. 从全国瓦片到分省包:我的裁切流程与翻车点

3.1 工具组合:为什么用GDAL而不是可视化工具

很多人习惯在可视化GIS软件里用“按掩膜提取”工具切省界,小范围研究区完全没问题。但全国级别的原始FathomDEM数据量很大,可视化软件一次性加载全国图层非常吃力,切到一半可能就卡死了。

我用的是GDAL命令行加Python脚本的组合。GDAL对于栅格裁切、重投影、压缩、构建概览这些操作非常成熟,关键是一套命令写好后,可以批量跑完所有省份,中间不用人工干预。整个过程可控、可复现,后续有人问你“这个省怎么切的”,你直接把命令发过去就行。

我用到的完整流程大致是:先获取省界矢量,再统一坐标参考,接着用gdalwarp裁切,裁完做压缩优化,最后生成校验文件。

3.2 裁切代码与关键参数

这里给一套最核心的裁切命令,命令本身不复杂,但参数里全是细节:

#!/bin/bash # 以某省边界裁剪FathomDEM全国栅格 INPUT="/data/fathomdem/world_30m_4326.tif" BOUNDARY="/data/boundary/hb.shp" OUTPUT="/data/dem_province/hb_30m_dem.tif" # 确保边界矢量与栅格坐标系一致 ogr2ogr -t_srs EPSG:4326 hb_4326.shp "$BOUNDARY" # 用cutline裁切,同时设置NoData和压缩参数 gdalwarp -cutline hb_4326.shp \ -crop_to_cutline \ -dstnodata -9999 \ -tr 0.0002777777778 0.0002777777778 \ -r bilinear \ -co COMPRESS=DEFLATE \ -co TILED=YES \ -co BLOCKXSIZE=256 \ -co BLOCKYSIZE=256 \ "$INPUT" "$OUTPUT"

逐条解释一下核心参数:

  • -crop_to_cutline:让输出范围严格贴着省界,裁出来的栅格是“凸多边形”而不是一个巨大的矩形。不加这个参数,一个省的文件体积可能膨胀好几倍,里面全是NoData空地。
  • -dstnodata -9999:设置输出NoData值,保证负海拔区域不被吃掉。
  • -tr 0.0002777777778 0.0002777777778:这是30米分辨率的经纬度网格约值。严格来说,1弧秒对应约30.87米,转换成十进制度数就是0.0002777777778度。这里用了约等于,实际项目中够用。
  • -r bilinear:重采样方法选双线性插值。对于高程连续表面,双线性比最邻近法平滑,比三次卷积法速度快很多。如果你后续做严格水文分析,建议改用-r cubicspline,保留的地形连续性更好,但耗时更长。
  • -co COMPRESS=DEFLATE:压缩成DEFLATE格式,体积减少非常明显,尤其对DEM这种连续缓变数据,压缩率通常能到8到12倍。
  • -co TILED=YES:输出瓦片化存储,便于局部读取和快速渲染。

跑完上述命令后,建议再做一步,给每个省的GeoTIFF构建概览(overview)。构建概览后,在GIS软件中缩放浏览时会明显感觉流畅很多:

gdaladdo -r average "$OUTPUT" 2 4 8 16 32

3.3 三个翻车点一次说清楚

翻车点一:边界矢量本身有缝隙或重叠。部分公开的省级面矢量在省与省交界处会有微小重叠或者缝隙,这会导致裁出来的栅格在接边处出现一些三角形孔洞或重叠带。处理方式是先对边界做一次拓扑修复,用QGIS的修复工具或GDAL的缓冲运算,把边界整体外扩1米再内缩1米,这样能消除绝大部分微缝隙。

翻车点二:省界边缘的NoData锯齿。直接按省界裁切,边缘像素和矢量边界会产生锯齿感,这在分辨率较低的数据上尤其明显。如果你只是看整体,锯齿可以忽略;但如果你做全省栅格转矢量,边界上的“城墙锯齿”会非常碍眼。解决方案是在裁切前把边界外扩一小段缓冲区,裁完以后再按边界裁剪一次。不过我实测下来,对30米DEM而言,这个锯齿的影响远小于河流分析和坡度分析的误差,如果不是出精美图件,可以不做缓冲区。

翻车点三:投影坐标系的隐性错位。我一开始直接拿原始3857栅格和WGS84省界做裁切,没转坐标系,结果切出来的地形明显偏移,误差在几百米量级。后来统一先把边界转成和栅格一致的坐标系再裁。这类错误在GIS里经常发生,整体表现是“地形和边界对不上”,排查时可以先用点图层检查边界点和栅格的坐标是否一致。

3.4 压缩、校验与发布细节

裁切完成后,我用sha256sum给每个文件生成了校验码文件,便于下载后核对完整性。这一步对大数据包特别重要,网络上传输过程中偶尔会丢字节,没有校验文件很难定位问题。

此外,README里我还记录了每个省的坐标范围、像素尺寸、NoData设置方式,以及数据来源和处理时间。别小看这份说明,三个月后你自己回来看包时,就靠它回忆“当时这个NoData到底怎么设置的”。

4. 换上FathomDEM后,水文分析和坡度计算的结果变化

4.1 河网提取:伪洼地明显变少,但填洼仍要设置阈值

用原始Copernicus DEM跑水文分析的朋友,大概率都有过被“伪洼地”支配的恐惧。雷达噪声在平原区域制造出大量不存在的凹陷,用GIS软件填洼时,一次性填掉的像元数量可能多到离谱,生成的河网在农田区常常出现一堆平行的小河道。

换成FathomDEM之后,最直观的变化是:需要填充的虚假洼地数量显著减少,提取出来的河网更接近真实水系形态。这是它最值得称道的地方。它保留了真实地形中的小沟谷和下切痕迹,但不会再像原始版那样把平地弄得坑坑洼洼。

不过这不代表你可以完全不填洼。某些喀斯特地区、采空区、人工开挖区仍然存在真实的洼地。我的建议是:不要用“直接填到0”这种粗暴方式,而是采用阈值填洼法,比如在SAGA GIS里用Fill Sinks(Wang & Liu)算法,设置合理的最小坡度阈值(比如0.01度),只填充那些小于阈值的伪洼地,保留真实负地形。具体阈值可以根据你的研究区平均坡度做几次试验。

4.2 坡度和山体阴影:细节变多,噪声也会被放大

FathomDEM细节增强带来的另一个问题是:如果你直接用原始像素计算坡度,会发现坡面在细节上有点“过活跃”,局部会出现细碎的坡度跳变。这不是说数据变得更差,而是细节增强以后,微地形信号被放大了,坡度算法又把噪声一并放大了。

我的处理顺序是:先对DEM做一个轻度的低通滤波(比如3x3均值滤波或高斯滤波),再做坡度、坡向计算。这样既能保留宏观地形趋势,又不会让坡度图布满细点。

山体阴影同样如此。直接用FathomDEM生成的Hillshade,在放大到1:1时细节非常漂亮,山脊、冲沟看得一清二楚;但在高倍率下,部分陡峭区域的阴影块会比原始DEM更“尖锐”。这属于正常现象,做可视化时可以适当调整对比度或透明度。

4.3 与实测数据的预期差:别拿DEM当天量高程

很多人问:这套30米数据的垂直精度到底多少米?恕我直言,任何单一数字放在全国尺度上都没有太大意义。区域不同、地形不同、地表覆盖不同,误差特征完全不同。

从公开的评估情况来看,FathomDEM在中低坡度、低植被覆盖区域的提升明显,垂直误差分布比原始Copernicus DEM更加集中,整体上更适合做流域划分、视域分析、淹没模拟的大范围初筛。但在高海拔陡峭山区和茂密森林区,雷达信号本身穿透深度有限,修正模型也难完全消除误差,这时候不要指望它达到精确测量级别的精度。

我的建议是:如果你手上有局部实测高程点,把分省包中的对应位置提取出来做一次差值直方图,看看偏差趋势。多数情况下你会发现差值呈正态分布,均值不大,但标准差可能比你想象的大。这个“标准差”才是你后续做误差传播分析时要用的关键参数。

5. 跨省拼接与可视化:让分省包真正派上用场

5.1 用虚拟栅格拼接多个省

分省包最大的便捷是“按需拿”,但如果你研究区恰好跨了两个省,你不需要把整个文件再重新裁剪一遍。最优雅的方案是用GDAL构建虚拟栅格(VRT),只生成一个轻量的XML文件,指向若干个省的GeoTIFF,加载和分析时直接以虚拟栅格作为输入。

gdalbuildvrt cross_province.vrt hb_30m_dem.tif ah_30m_dem.tif

构建完成后,在GIS软件里打开这个VRT文件,感觉就像打开了一个完整的跨省DEM,不需要等待巨大TIF的加载。后续做填洼、坡度、河网分析,都可以直接把VRT作为输入数据。

要注意的是,VRT只是逻辑拼接,不会把数据真正合并。如果分析工具不支持VRT,或者你要生成高程分层文件,可以再执行一次gdal_translate把它转成真正的GeoTIFF,这时建议保留与分省包相同的NoData和压缩设置。

5.2 快速出一张能直接放进报告的地形渲染图

分省包配合GDAL自带的山体阴影算法,可以快速做一张观感不错的地形渲染图。我最常用的两段命令是:

# 生成山体阴影 gdaldem hillshade hb_30m_dem.tif hb_hillshade.tif \ -z 1.3 -az 315 -alt 45 -compute_edges # 生成彩色高程渲染 gdaldem color-relief hb_30m_dem.tif color_ramp.txt hb_color.tif

参数说明:-z 1.3是垂直夸张系数,适合中纬度平缓地形,山区可以调低到1.0或0.8;-az 315是光源方位角,来自西北方向的光源比较符合大多数人看地形的习惯;-alt 45是太阳高度角,45度最能体现地形立体感,如果地形非常平坦,可以调到30度左右增强阴影对比。

彩色高程渲染的color_ramp.txt是一个文本颜色表,格式很简单,网上有很多现成的色带可以直接用。把彩色图放到图层上,山体阴影图放在下面,用“叠加”混合模式叠加,透明度调到70%左右,一张能直接放进报告里的地形底图就出来了。这套方法不挑软件,QGIS和ArcGIS里都用得上。

5.3 两个值得继续折腾的扩展方向

如果你愿意再花一点时间,接下来可以尝试两个方向。

一是把分省DEM转成等高线并做抽稀。用GDAL的gdaldem contour可以快速生成等高线,但直接生成的结果在小比例尺下会非常杂乱,需要用简化算法抽稀,或者按等高距分级导出。得到的等高线可以直接用于制图,也可以转成CAD格式给非GIS背景的同事使用。

二是把DEM与土地覆盖数据叠加,做地形位分类。把海拔、坡度、坡向和最新的土地覆盖分类放在一起做叠加分析,可以快速识别“高山森林带”“平原农田带”等地形-生态功能区。这类分析对规划、林业、生态评估项目非常有用,也是把单纯的DEM数据转化为决策信息的典型路径。


这套中国30米分省包我自己整理了好几天,随后在流域模拟里连续用了两三个月,最明显的感受是过去需要反复实验的预处理脚本基本可以删掉了。真要说缺点,也不是没有:它对平坦水面的处理仍然和原版一样会显得有些扁平,如果你做湖泊库容这类研究,还是得结合实测水深再做调整。最后提醒一句,不管从哪个渠道拿到DEM,开工前一定先看坐标和NoData,这两个问题最隐蔽,也最容易让前面的分析全部白做。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询