全球DEM数据校验、格式选型与GIS地形分析实操 简介面向GIS从业者与地理信息学习者的全球DEM地形数据获取文档聚焦全球高程数据下载渠道分散、免费工具与学习资料难配套的问题。压缩包仅含1个docx文件约12KB篇幅精炼集中整理全球范围DEM数据的百度网盘与天翼云盘入口并附免费GIS软件及GIS教程资源链接便于按需选用。文档兼顾数据获取、软件处理与入门学习三类需求涉及QGIS、GRASS GIS、MapInfo等工具以及视频教程、操作指南和案例分析可帮助读者完成地形加载、坡度与排水路径分析、空间数据处理和制图等任务也可与ArcGIS的地形切割、坡度计算、流线分析等功能衔接。已有6002人学习适合预算有限或希望快速补齐数据与工具链的初中级GIS用户参考。1. 全球DEM地形数据拿到手之后先跑一遍校验再谈分析全球DEM地形数据这东西找起来比用起来麻烦十倍。网上搜一圈要么是几TB的付费商业数据集要么是残缺不全的测试瓦片下完发现缺了半个大陆。这份资源包把三样东西拢到了一起全球范围的DEM高程瓦片、能直接打开这些数据的免费GIS工具、以及一套从加载到分析的教程。走的是云盘分发两个渠道选快的那个就行。它适合的人很具体GIS刚入门需要一套完整数据练手的新手、做方案验证时需要快速拿到全球地形底图的工程师、以及预算有限但不想在数据源上反复试错的小团队。后面几章按“格式判断→解压验证→分析实操→踩坑排查→批量预处理”的顺序走每一步对应一个能直接抄的操作不绕弯子。2. DEM格式与分辨率选型从瓦片命名到空洞分布拿到一份DEM数据第一件事不是拖进软件看效果而是搞清楚它的格式、分辨率和数据结构。这一步决定了后面所有分析结果靠不靠谱。DEM本质上是一个栅格矩阵每个像元存一个海拔值但存法五花八门——有的带完整地理标签有的只有裸数据靠文件名猜范围有的分辨率30米有的90米有的海洋区域是NoData有的直接填0。不先把这些搞清楚后面算坡度能算出全零做水文分析能给你画出一条穿过山脊的河流。2.1 常见格式与瓦片命名规则全球DEM数据从分发形态上通常分三类。GeoTIFF是最常见的一个文件包含栅格数据加内嵌的地理坐标信息拖进ArcGIS或QGIS直接能定位。HGT是航天飞机雷达地形测绘任务的原始格式纯二进制高程值没有头文件坐标范围靠文件名推出来——比如N45W120表示北纬45度到46度、西经120度到119度这一个方块。BIL/BSQ/BIP是原始二进制格式需要配一个.hdr头文件说明行列数和数据类型缺了头文件就只剩一堆字节。这三种格式之间可以互转GDAL一条命令的事但前提是你得先认出它是什么。全球DEM瓦片一般按1°×1°切分每个瓦片3601×3601像素比标准1度对应3600像素多一圈多出来的边缘那一圈用于相邻瓦片拼接时的重叠匹配。30米分辨率下单个GeoTIFF大约25MBHGT约20MB90米分辨率的瓦片则小得多。如果你解压后看到文件名是N45W120、S30E135这种格式基本可以确定是按经纬度分块的全球数据。文件名就是坐标不需要打开文件就能知道它覆盖哪块区域。格式是否自带坐标信息典型文件大小30m1°×1°常见来源GeoTIFF是内嵌约20-30MB多数处理后分发数据HGT否靠文件名约20MBSRTM原始数据BIL需配.hdr约20MB部分科研数据集提示HGT文件拖进软件之前先确认软件是否自动按文件名赋坐标系。ArcGIS对HGT的识别有时会默认给个错误投影加载后手动设WGS84更稳。2.2 三个检查步骤判断数据能不能直接用第一个检查点是文件大小。一个30米分辨率的全球瓦片GeoTIFF正常在20-30MB之间。如果某个文件只有几百KB大概率是下载不完整或者本身就是低分辨率缩略图。如果几个瓦片大小差异很大也要留意——可能是数据源混用了不同分辨率。第二个检查点用GDAL读元数据。命令行一行命令就能把关键信息全拉出来# 查看单个DEM瓦片的元数据与统计值 gdalinfo -mm -stats N45W120.tif-mm计算最小/最大值-stats额外算均值标准差。输出里重点看三行Size is 3601, 3601确认瓦片尺寸是否正常Coordinate System is ...确认坐标系是否WGS 84NoData Value看无数据值设成了多少。SRTM不同版本里这个值可能是-32768也可能是-9999还有的版本根本没设NoData。没设NoData的后果是海洋区域以0米参与计算坡度分析时海岸线附近会出现异常平坦或异常陡峭的假象。第三个检查点是看空洞分布。SRTM数据在雷达反射率低的区域水面、雪盖、城市高楼密集区存在数据空洞表现为条带状或块状的NoData。快速办法是在QGIS里加载后把渲染方式改成单波段伪彩色NoData区域会显示为透明或统一颜色一眼能看出空洞在哪。大面积的空洞需要后续插补零散的小坑对区域级分析影响有限。# 批量快速检查目录下所有DEM文件的尺寸和极值 for f in *.tif; do echo $f gdalinfo -mm $f | grep -E Size is|Min/Max|NoData done这段循环只提取关键行适合在拿到一批瓦片后快速过一遍确认没有明显异常的文件。grep -E里的三个关键词对应尺寸、高程范围和无数据值。如果某个文件的Min和Max相等说明整个瓦片是平的要么是水域要么是坏文件。2.3 选型30米够不够用什么时候该找更高分辨率分辨率直接决定分析结果的精度。坡度计算中分辨率越粗陡坡被相邻像元平均化的程度越高算出来的坡度偏缓。有实际对比测试表明90米DEM算出的坡度比30米DEM系统性偏低10%-20%在地形起伏剧烈的山区差距更大。水文分析中分辨率决定了提取出的河网密度——30米数据能识别出较细的支流90米只能看到主干。这份资源包提供的全球数据以30米级为主对以下几类工作完全够用区域尺度的坡度坡向制图、流域边界提取和汇流分析、风电场和光伏电站选址初筛、洪水淹没范围模拟的底图。做城市内部精细土方计算、单体建筑选址、或者需要区分田埂和沟渠的农业排水设计30米就不够了得找12米甚至更高精度的数据那份资源的覆盖范围不一定满足需要另找来源。判断标准很简单如果你的研究区在一个县或一个流域的尺度30米数据是合适的起点。如果你要在一个平方公里内做精细地形分析先确认这份数据的分辨率再动手别花半天跑完流程才发现像元比地块还大。3. 从云盘压缩包到可加载图层解压、验证与投影处理下载完的压缩包解压之后目录结构往往比想象中乱——有的按大洲分文件夹有的把所有瓦片平铺在一个目录里还有的夹着说明文档和无关文件。先摸清结构再动手比直接全选拖进软件省时间。这一章走完解压、完整性验证、坐标系处理三步数据才算真正可用。3.1 解压后的目录结构怎么读拿到压缩包先别急着解压用列表命令看一眼里面有什么# 查看压缩包内容清单不实际解压 unzip -l global_dem.zip | head -50 # 如果是7z格式 7z l global_dem.zip | head -50重点看目录层级和文件命名。如果瓦片平铺在一个目录里、文件名全是经纬度格式那说明是标准化分块数据可以直接进下一步。如果按大洲或国家分了子目录先确认你要的研究区对应哪个目录只解压需要的部分省磁盘也省时间。解压时注意路径不要带中文或空格GDAL和ArcGIS对中文路径的兼容性时好时坏血泪经验是养成全英文路径的习惯# 解压到指定目录-d 指定输出路径-q 静默模式 unzip -q global_dem.zip -d /data/dem_global/ # 解压后快速统计tif文件数量 find /data/dem_global/ -name *.tif | wc -l-d指定解压目标目录-q抑制解压过程的逐文件输出。解压后用find统计tif文件数量和压缩包清单里的数量对一下少了就是解压不完整。3.2 用GDAL验证文件完整性解压完逐个文件跑gdalinfo太慢写个脚本批量过一遍只输出异常项import subprocess import os dem_dir /data/dem_global/ issues [] for f in sorted(os.listdir(dem_dir)): if not f.endswith(.tif): continue path os.path.join(dem_dir, f) # 只读元数据不加载全图速度快 result subprocess.run( [gdalinfo, -mm, path], capture_outputTrue, textTrue ) out result.stdout # 检查是否包含坐标系信息 if Coordinate System is not in out: issues.append(f{f}: 缺少坐标系信息) # 检查是否有NoData定义 if NoData Value not in out: issues.append(f{f}: 未定义NoData值) # 检查文件大小是否异常小于5MB视为可疑 size_mb os.path.getsize(path) / (1024 * 1024) if size_mb 5: issues.append(f{f}: 文件仅{size_mb:.1f}MB疑似不完整) print(f检查完成发现 {len(issues)} 个问题) for i in issues: print( i)这段脚本做了三件事调gdalinfo读元数据、检查坐标系和NoData是否存在、用文件大小做粗筛。capture_outputTrue把命令输出捕获到变量里而不是打到终端textTrue让输出以字符串而非字节形式返回。实际使用中可以把issues列表写到CSV里备查。问题文件不多的话手动补元数据或重新下载就行问题文件多说明下载环节出了系统性故障得回去检查云盘下载的完整性。注意有些DEM数据在制作时故意不写NoData值比如全球陆地数据把海洋裁掉了整个文件都是有数值的陆地这种情况不算异常跳过即可。关键看数据来源说明里有没有提到。3.3 坐标系识别与统一投影WGS84地理坐标系下水平坐标单位是度垂直高程单位是米两者量纲不一致。直接在这个坐标系下算坡度软件会把1度当成1米来除结果完全错误。ArcGIS的坡度工具在地理坐标系输入时会弹警告QGIS的坡度算法则可能静默算出错误值。两种处理路径手动设Z因子做换算或者先投影到米制坐标系再算。Z因子的思路是把经度方向的实际地面距离算出来在纬度φ处1度经度约等于111320×cos(φ)米1度纬度约等于110574米。Z因子取1/(111320×cos(φ))量级在10⁻⁵左右。但这个方法在跨纬度较大的区域误差明显因为cos(φ)随纬度变化。更稳妥的做法是先重投影到UTM或Albers等积投影让水平垂直单位统一为米# 将WGS84地理坐标系的DEM重投影到UTM 50NEPSG:32650 gdalwarp -t_srs EPSG:32650 \ -tr 30 30 \ -r bilinear \ -co COMPRESSLZW \ N45W120.tif N45W120_utm50.tif-t_srs指定目标坐标系EPSG:32650对应北半球UTM第50带覆盖东经117°到123°。-tr 30 30把输出分辨率固定为30米×30米避免重投影后像元大小变化。-r bilinear用双线性插值重采样比最近邻平滑但会轻微改变极值做坡度分析够用如果后续要做精确的水文分析用-r cubicspline更平滑。-co COMPRESSLZW输出时做无损压缩文件体积能小一半左右。UTM分带的选法北半球EPSG从32601开始南半球从32701开始带号 floor((经度 180) / 6) 1。研究区跨多个带时选中心经度所在的带或者改用Albers等积投影适合大区域。投影后的DEM再进坡度计算Z因子保持默认1即可不用再做换算。4. ArcGIS与QGIS地形分析坡度、水文与镶嵌实操数据准备好了进分析环节。ArcGIS和QGIS都能做地形分析操作路径不同但底层逻辑一致坡度看的是像元与邻域的高程变化率水文分析看的是水流方向累积。这一章把坡度、山体阴影和水文分析三个高频操作的参数设法和注意点讲清楚外加大区域拼接的VRT方案。4.1 坡度与山体阴影Z因子和单位怎么设ArcGIS里的路径是ArcToolbox → 3D Analyst Tools → Raster Surface → Slope。输入投影后的DEM输出测量单位选DEGREE度或PERCENT_RISE百分比。如果输入是地理坐标系度单位且没有投影Z因子必须手动设成约1/(111320×cos(中心纬度))否则结果不可用。如果输入已经是UTM投影米单位Z因子保持1。QGIS里的路径是Processing Toolbox → Raster terrain analysis → Slope参数界面更简洁单位选度或百分比Z因子在高级参数里。QGIS不会在输入地理坐标系时弹出警告需要自己检查。参数ArcGISQGIS说明输入图层投影后DEM投影后DEM地理坐标系需设Z因子测量单位DEGREE / PERCENT_RISE度 / 百分比坡度制图常用度Z因子1投影后1投影后地理坐标系需换算输出像元大小默认与输入一致默认与输入一致重采样后要核对山体阴影Hillshade用来做地形可视化底图方位角设315度西北方向、高度角设45度是通用默认值。方位角315能让山脊线和沟谷产生明暗对比比正北方向更能突出地形起伏。这两个参数在ArcGIS的Hillshade工具和QGIS的Hillshade工具里位置相同效果也基本一致。4.2 水文分析链路填洼、流向、流量累积ArcGIS的水文分析工具集在Spatial Analyst Tools → Hydrology下。标准链路是三步填洼Fill→ 流向Flow Direction→ 流量累积Flow Accumulation。填洼的目的是消除DEM中的伪洼地——真实洼地湖泊和雷达测高的噪声小坑分不清不填的话水流会在小坑里断掉提取出的河网是碎的。填洼参数里有个Z limit填洼阈值不设的话所有洼地不管大小全部填平。对全球DEM来说大型湖泊会被填成平地后续流向计算会把湖泊当作一个方向流出的点影响下游汇流。稳妥做法是先设一个阈值比如10米只填平深度小于10米的小坑保留真实地形洼地。阈值怎么定取决于研究区地形起伏——平原区可以设小一点山区可以设大一点。流向计算用D8算法八方向输出栅格值1到255对应八个方向。ArcGIS的Flow Direction工具有个选项“Force edge cells to flow outward”选是的话边缘像元强制向外流避免边缘积水导致累积量异常。流量累积的输入是流向栅格权重默认1输出每个像元的上游汇水像元数。做全球或大洲尺度的水文分析之前务必先按研究区裁剪。一个30米分辨率的全球拼接DEM动辄几十GBArcGIS加载后做填洼运算内存直接爆掉。用gdalwarp的-te参数按经纬度范围裁出研究区# 按经纬度范围裁剪DEM-te 后跟 xmin ymin xmax ymax gdalwarp -te 100 25 110 35 \ -t_srs EPSG:32649 \ input_dem.tif clipped_dem.tif-te指定裁剪范围顺序是西经、南纬、东经、北纬。-t_srs可同时做投影。裁剪后再进ArcGIS做水文分析运算量小几个数量级。4.3 大区域镶嵌VRT虚拟拼接与重采样研究区跨多个瓦片时直接拼成一个大TIFF有两个问题文件巨大全球拼接轻松超过100GB、拼接过程耗时长。GDAL的VRT虚拟栅格方案解决了这个问题——它不生成新文件只写一个XML索引记录瓦片的位置关系读取时按需拼合。# 第一步将所有瓦片建成VRT索引-resolution highest 保留最高分辨率 gdalbuildvrt -resolution highest -input_file_list tile_list.txt mosaic.vrt # 第二步将VRT转换为实际GeoTIFF按需指定输出分辨率和压缩 gdal_translate -tr 30 30 \ -r bilinear \ -co COMPRESSLZW \ -co BIGTIFFYES \ mosaic.vrt mosaic_30m.tif-input_file_list指定一个文本文件里面每行一个瓦片路径适合瓦片数量多时批量传入。-resolution highest让VRT以最高分辨率瓦片为基准避免混合分辨率时自动取平均。gdal_translate那一步才真正生成物理文件-tr 30 30统一输出分辨率-co BIGTIFFYES允许输出超过4GB的文件普通TIFF上限4GB全球30米数据远超这个值COMPRESSLZW做压缩。如果只是临时查看不打算长期保存直接加载mosaic.vrt到QGIS即可不用跑第二步。ArcGIS对VRT的支持不如QGIS流畅大区域建议在QGIS里用VRT浏览和裁切导出小范围结果后再进ArcGIS做分析。5. 避坑排查五条踩坑记录与解决路径下面这五条是处理全球DEM数据时最常撞上的问题每条按“现象→原因→解决”拆开写。有些坑我踩过不止一次写下来省得下次再翻。5.1 解压报错CRC错误还是分卷缺失现象解压到一半报“unexpected end of file”或“CRC failed”部分文件解压失败。原因从云盘下载大文件时网络中断导致文件不完整这是最常见的情况。另一种是分卷压缩包.z01、.z02、.zip少下了一个分卷也会报错。还有一种少见情况是解压工具版本太老不支持压缩包使用的算法。解决先对照云盘上显示的文件大小和本地文件大小不一致就重新下载。分卷压缩必须从第一个文件.zip或.z01开始解压缺任何一个分卷都不行。如果大小一致但报CRC错误换个解压工具7-Zip或Bandizip再试。下载大文件时尽量用支持断点续传的方式一次下完不中断。5.2 加载后位置偏移坐标系缺失还是标注错误现象DEM拖进ArcGIS后和矢量边界对不上整体偏移几百米到几公里。原因HGT格式没有内嵌坐标系软件按默认值处理导致位置偏移。另一种是数据本身标注了坐标系但标错了比如实际是WGS84却标成了CGCS2000几十米的偏移在低精度底图上看不出来和精确矢量叠加时暴露。解决先用gdalinfo确认坐标系字段。HGT文件手动设WGS84。如果确认坐标系标注有误用gdal_edit.py修改元数据而不动栅格数据本身# 强制将坐标系设为WGS84不重采样不改变像元值 gdal_edit.py -a_srs EPSG:4326 N45W120.tif-a_srs只改元数据里的坐标系标签不触碰栅格数据。改完再用gdalinfo确认。这个操作不可逆改之前备份原文件。5.3 坡度计算大面积空值NoData没设对现象原始DEM看起来有数据但坡度计算结果大块区域为空和原始数据的空洞区域对不上。原因DEM的NoData值设置不对。比如原始数据海洋区域填的是-9999但文件元数据里没声明NoData软件把-9999当成真实海拔参与坡度计算导致海岸线附近坡度异常后续处理时又被当成异常值剔除。反过来有数据区域被误标为NoData也会造成大面积空值。解决gdalinfo确认NoData值用gdal_translate重新指定# 将NoData值重设为-9999或数据实际的无数据值 gdal_translate -a_nodata -9999 input.tif output.tif-a_nodata修改输出文件的NoData声明。如果原数据没有NoData但实际有填充值先确认填充值是多少用gdalinfo -mm看Min值如果Min是一个异常大的负数就是填充值再用-a_nodata声明它。5.4 拼接后接缝裂缝重叠区没处理现象两块DEM拼接后接缝处出现细线状裂缝或高程跳变山体阴影图上尤其明显。原因相邻瓦片边缘有1像素重叠直接拼接时重叠区被重复计算或丢弃导致接缝处出现一条异常带。不同来源的数据比如SRTM和ASTER混用高程基准不一致也会在接缝处产生跳变。解决用gdalbuildvrt的-resolution highest或-resolution average统一处理重叠区。拼接前先用gdalwarp把所有瓦片重采样到完全一致的网格相同分辨率、相同范围对齐再建VRT。不同来源的数据不要混拼先统一到同一高程基准再合。5.5 ArcGIS打开大文件崩溃内存不足还是金字塔没建现象加载拼接后的大DEM几十GB时ArcGIS卡死或报内存错误。原因32位ArcGIS有约2-4GB的内存上限大TIFF加载后直接爆内存。即使64位版本没有金字塔的TIFF每次缩放都要重读全图速度极慢。超过4GB的普通TIFF格式本身也会出问题。解决用gdaladdo建金字塔用BIGTIFF格式存储# 建立金字塔2 4 8 16 表示逐级抽稀的倍率 gdaladdo -r average mosaic_30m.tif 2 4 8 16-r average用平均值法生成缩略层适合连续高程数据。金字塔建好后ArcGIS加载时只读当前缩放级别对应的层级速度大幅提升。超过4GB的文件用gdal_translate -co BIGTIFFYES转格式或者干脆用VRT方式在QGIS里浏览只把裁切后的小范围导出到ArcGIS。6. 进阶用法Python批量验证与预处理脚本骨架单个文件处理手动操作没问题一批瓦片几十上百个的时候就得靠脚本。我一般用Python调GDAL做两件事批量体检和自动拼接预处理。6.1 批量检查空洞率与极值分布from osgeo import gdal import numpy as np import csv, os dem_dir /data/dem_global/ report [] for f in sorted(os.listdir(dem_dir)): if not f.endswith(.tif): continue ds gdal.Open(os.path.join(dem_dir, f)) band ds.GetRasterBand(1) arr band.ReadAsArray() nodata band.GetNoDataValue() # 计算空洞比例 if nodata is not None: hole_ratio np.sum(arr nodata) / arr.size valid arr[arr ! nodata] else: hole_ratio 0.0 valid arr.ravel() report.append({ file: f, min: float(valid.min()), max: float(valid.max()), mean: float(valid.mean()), hole_pct: round(hole_ratio * 100, 2) }) ds None # 关闭数据集释放内存 with open(dem_report.csv, w, newline) as fp: writer csv.DictWriter(fp, fieldnamesreport[0].keys()) writer.writeheader() writer.writerows(report) print(f完成 {len(report)} 个文件的检查结果写入 dem_report.csv)ReadAsArray()把整个波段读进numpy数组30米瓦片约3601×3601像元内存占用可控。GetNoDataValue()取元数据里的无数据值为None时说明文件没有声明NoData跳过空洞计算。输出CSV里hole_pct超过5%的瓦片单独标记这些瓦片后续做水文分析前需要插补处理。ds None那行不能省GDAL的数据集不手动释放的话循环几十个文件后内存会持续上涨。6.2 自动化拼接与重采样的脚本骨架批量预处理的核心是“建VRT→重投影→导出”三步自动化脚本把瓦片列表和参数传进去就行from osgeo import gdal # 第一步建VRT tile_list [N45W120.tif, N45W121.tif, N44W120.tif] vrt gdal.BuildVRT(mosaic.vrt, tile_list, resolutionhighest) vrt None # 写完关闭 # 第二步重投影重采样压缩输出 gdal.Warp(mosaic_utm.tif, mosaic.vrt, dstSRSEPSG:32650, xRes30, yRes30, resampleAlgbilinear, creationOptions[COMPRESSLZW, BIGTIFFYES]) print(拼接与重投影完成)BuildVRT的resolutionhighest与命令行参数对应。Warp函数把-t_srs、-tr、-r等参数统一收进关键字参数creationOptions对应-co。这种写法比subprocess调命令行更容易在循环里批量替换参数比如按研究区换EPSG编码。脚本跑完后用gdalinfo确认输出文件的尺寸和坐标系确认没问题再进ArcGIS做后续分析。从那以后我每次拿到新的DEM数据都强制走一遍“gdalinfo看元数据→脚本扫空洞→确认投影再进软件”三步不再直接拖进ArcGIS碰运气。这套流程写进脚本里跑一次不到两分钟省掉的是后面反复排查“为什么坡度算出来不对”的几个小时。希望帮到你。本文还有配套的精品资源点击获取