ARTICLE DETAIL

资讯详情

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

全国地形地貌地质数据集处理全流程:从DEM到统一底图的避坑指南

全国地形地貌地质数据集处理全流程:从DEM到统一底图的避坑指南 简介《全国地形地貌地质数据集》是一份面向GIS专业人员、科研工作者及地理爱好者的综合性地理数据包可满足地质勘探、灾害评估、农业规划、生态保护等多场景的空间分析需求。压缩包共609个文件大小约772.86MB包含adf、nit、tif、shp等栅格与矢量格式以及配套的lyr图层文件、prj投影文件和dbf属性表方便在ArcGIS、ENVI中直接加载与制图dat、log等辅助文件可用于数据转换与过程记录。已有4129人学习下载热度较高。内容覆盖全国地质图、地质矿产、各省耕地面积、水系、地貌、地形DEM、森林分布、土地利用、土壤类型和植被分布十大类主题可为学术研究、政策制定与工程应用提供完整的数据底座。通过该数据包用户既能提取高程、坡度、坡向等地形因子也可结合属性表完成专题制图与统计分析大幅减少数据采集与预处理时间适合地理信息相关课程设计、区域规划项目及科研前期探索使用。1. 全国地形地貌地质数据集是什么值不值得为它搭一条处理流水线在水利、交通、风电场选址还有跨区域环评这类项目里我拿到任务后做的第一件事几乎都一样先把项目区的地形骨架和地质底子摸清楚。所谓“全国地形地貌地质数据集”不是某一个固定网盘压缩包而是一条可复现的组合路线——把全国范围的高程网格DEM、地貌类型面、岩性/地层分布和断层构造线统一到一个坐标系和裁剪边界里形成一张能直接进 GIS、能套进建模软件的底图。它回答的核心问题就是这里高程多少、属于什么地貌、地下以什么岩性为主、有没有断层穿过选址红线。适合做勘察布点、区域规划、灾害评估和行业数据服务的工程师快速建底图。2. 拆开再看高程、地貌、地质三部分各是什么规格2.1 高程层DEM分辨率、垂直基准与“这里的海拔不可直接当真”高程数据是整个数据集里最直观也最常用的一层常见来源是 SRTM、ASTER GDEM、ALOS World 3D 这类全球覆盖的栅格产品。普通项目里30 m 和 90 m 两种分辨率是绝对主流全流域或省级范围用 90 m 能显著减少处理时间几十平方公里的关键区域再换 30 m 做细活。市面上也有 12.5 m 甚至 1 m 的商业或机载数据但那些不属于“全国数据集”的默认范畴需要单独立项采购。选 DEM 时比分辨率更需要先确认的是高程基准。SRTM 系列的公开产品通常基于 EGM96 大地水准面也就是我们日常理解的海拔高度但部分 GDEM 产品或者经过后处理的栅格会直接使用 WGS84 椭球面高程两者之间的差值在山地可能达到十几到几十米。这个数值不是误差而是“口径不一致”。新手最常见的翻车是直接把椭球高当成海拔用叠加到工程测量成果上时发现整体偏移还以为是坐标系坏了。我把三套常用 DEM 的选型口径整理成一张对比表方便你在拿到原始文件后先对号入座数据系列常用分辨率垂直基准适用范围使用注意SRTM 系列30 m / 90 mEGM96 大地水准面流域分析、坡度坡向、区域高程统计河谷和陡崖处容易有空洞需要填洼ASTER GDEM30 m局部拼接噪声明显需查看元数据部分区域按椭球高处理大范围快速浏览、定性判断建筑和树木密集区会出现高程拉高ALOS World 3D30 m局部有 12.5 m 版本与大地水准面相关山区细节要求较高的项目文件体积偏大全国范围不宜一次性全量加载拿到任何 DEM 文件第一件事不是急着切片而是先用元数据文件确认垂直基准。国内不少数据处理平台会对原始 DEM 做二次投影和换带这一步最容易把基准信息弄丢。如果元数据里没写基准我通常会拿项目区内已知的等级点或连续运行参考站高程做一次对比偏差超过 10 m 就要怀疑是椭球高而不是真实海拔。2.2 地貌层面状分类数据边界是“大概齐”而不是精确红线地貌层解决的是“这里长什么样”的问题。项目里用到的地貌数据往往来自全国或区域尺度地貌图属性字段会给出较大类型比如平原、台地、丘陵、山地、高原、盆地再细一点还有中地貌类型比如冲积扇、冲积平原、喀斯特地貌、黄土塬等。它是一份面状矢量数据每个多边形附带类型代码和名称方便做面积统计。使用这层数据时要记住一个原则地貌图的边界不等于实地突变线。我从野外核图经验看平原和丘陵的交界往往是一段过渡带而不是图上画出的某一条固定折线。所以地貌层适合做分区规划、评价单元划分、景观格局分析不适合用于征地边界这类精度要求极高的用途。拿到数据后优先检查属性表里是否有两个关键字段地貌大类代码和中地貌代码。有些老数据只有图形没有属性这种情况要么补字典要么放弃因为画在地图上分辨不出意义。2.3 地质层岩性、地层时代与断层三条线一起看地质层是数据集中专业门槛最高的一部分。通常包含面状的地层岩性分布图以及线状或点状的构造要素比如断层、褶皱轴、产状点。属性表里常见字段包括地层代号、岩性名称、时代、成岩环境等。全国尺度的公开地质图比例尺多为 1:100 万或 1:50 万能看出大的岩性分区和区域性断裂但不足以指导单体工程的定点和验算。我习惯把地质信息拆成三层来用第一层是岩性大类判断场地是土质还是岩质第二层是地层时代判断地层老新关系和可能的承压条件第三层是断层缓冲判断工程布点离活动构造是否有安全距离。这三层不需要同时在一个图上可视化但在做综合分析时缺一不可。较麻烦的是不同来源的地质图图例和字段命名并不统一有的用中文名称有的用代码。先把字段字典统一好后面所有叠加分析才走得通。2.4 三套数据的文件格式选型决定后期会不会返工高程层默认保存为 GeoTIFF保留原始位深和无数据值。如果担心后期要反复裁剪建议保留一份未做拉伸的原始整型 DEM不要直接保存成带有色彩映射的 8 位渲染图否则高程信息和精度都会丢失这是很多项目里“文件还在但数字不可用”的根源。地貌层和地质层我用 GeoPackage 比较多它把属性表、坐标系和图层都存在一个文件里避免 Shapefile 那种 .shp、.dbf、.prj 丢一个文件就开不了的情况。全国范围的数据量大但按项目边界裁剪后通常只有几十到几百兆GeoPackage 完全扛得住。老数据还是有大量 Shapefile 和 MapGIS 格式需要统一转一遍。转换时一定要保留属性字段的编码用 UTF-8 输出不然中文属性在 QGIS 里会变成乱码这在处理国产 GIS 软件导出的数据时经常碰到。3. 从公开数据源拼起下载、解压与原数据检查的固定动作3.1 高程数据的下载策略别一次性拉全国要按项目范围圈常见做法是去公开科学数据平台搜索“全国 DEM”或“SRTM 30 m”分幅数据国内平台通常按经纬度分幅组织文件。下载时不要抱着“存全天下”的心态而是按项目区外扩 510 公里后涉及的图幅范围选择文件。一次全量下载全国 30 m DEM 不只耗时还会让你在后续整理时因为文件数量太多而无从下手。我一般按两级目录管理原始下载包放downloads/解压后按数据类别放进raw/。下载完成后先做一个最基础的可读性检查文件是否损坏、能否正常打开栅格、是否带有地理范围信息。有些下载包文件名看着正常解压后里面套了一层子目录直接通配符读取会漏掉一半有的 TIF 文件缺失投影信息后期叠加会错位得莫名其妙。3.2 地貌与地质数据的搜索关键词以及“要矢量不要扫描图”搜索地貌与地质数据时用对关键词能省很多时间。地貌数据搜“全国地貌类型分布”、“1:100 万地貌图 Shapefile”地质数据搜“全国地质图矢量”、“区域地质图 GeoPackage”。判断一份数据能不能用第一看格式优先 Shapefile、GeoPackage、GeoJSON、DWG 转出的矢量看到 JPG、TIF 加配准文件的扫描图直接跳过因为后期矢量化代价太高且精度无法保证。地质数据还常见一个坑网页上显示是矢量下载下来却是栅格底图带一堆冗余图层。我一般先在本地用软体扫一遍属性表确认地层代号字段真的存在而不是只有颜色字段。只有图层、没有属性字段的“地质图”本质上就是一张彩色画做不了任何空间查询和统计分析。3.3 用一条统一命令把下载包变成规整的 raw 目录下载完一堆压缩包后先建目录结构再统一解压。不要让数据裸奔在桌面或系统下载目录里否则项目还没开始文件路径就已经失控了。下面这段脚本我每次都用逻辑很直白mkdir -p raw/dem raw/landform raw/geology for zip in downloads/*.zip; do base$(basename $zip .zip) mkdir -p raw/${base} unzip -o $zip -d raw/${base} done find raw -name *.tif | head -20 find raw -name *.shp | wc -l这段脚本先创建三个原始数据目录然后遍历downloads/下所有 zip 包按文件名建子目录并解压。最后两行分别查看已解压的 TIF 文件样本和统计 Shapefile 数量用数字确认数据真的落地了。注意此处按压缩包名建子目录解压嵌套的问题就不会把多个图幅的文件混在一起。如果你下载的包是 rar 或 7z 格式把unzip换成对应解压工具即可参数思路一致。3.4 解压后的体检用 gdalinfo 快速判断文件能不能用拿到解压后的 TIF 文件我习惯在 QGIS 里做正式加载之前先用命令行工具给这批数据做个体检。GDAL 自带的 gdalinfo 能输出一个栅格文件的基础信息行数、列数、像素尺寸、坐标系、无数据值、最值统计。这些字段决定后续裁剪和重投影的参数对不对在命令行里批量检查比在 QGIS 里逐个打开图层快得多。for f in raw/dem/*/*.tif; do echo ${f} gdalinfo -stats $f | grep -E Size is|Pixel Size|Coordinate System|NoData|STATISTICS_MINIMUM|STATISTICS_MAXIMUM done这段脚本对raw/dem下所有 TIF 逐个打印关键信息。重点看三个地方第一Size is后面两个数字是否和预期分辨率对得上第二NoData值是否一致第三STATISTICS_MINIMUM和MAXIMUM有没有出现超出高程合理区间的怪值比如海沟深度或珠峰高度出现这种值说明原始数据可能带了粗差或坏像元。凡是体检不合格的文件我会单独移到一个quarantine/目录不给它混进工作区的机会。4. 统一投影与空间裁剪把三套数据揉成一张工作底图4.1 全国底图为什么不用 UTM而选兰伯特等角圆锥投影高程、地貌、地质三套数据来自不同生产单位原始坐标系往往不一致DEM 栅格大多是 WGS 84 地理坐标区域地质图可能是 CGCS2000 或北京 54 那个年代的历史投影。把这些数据直接叠在一个画布里错位几十到几百米是家常便饭。所以工作的第一步是先确定一套统一输出投影让所有图层“说同一种语言”。全国尺度下我不倾向于用 UTM。UTM 分带之间接缝明显一旦项目范围跨两个中央经线在带上边界的变形会让面积统计和边长测量出现明显误差。做全国或大区级综合分析时我一般用兰伯特等角圆锥投影两条标准纬线取 25°N 和 47°N中央经线取 105°E这也是一张大区域底图里常见的投影参数套路。这种投影在东西方向上变形小适合中国东西跨度大的版图。处理时用下面这串参数作为统一目标坐标系SRC_PROJprojaea lat_125 lat_247 lat_00 lon_0105 x_00 y_00 ellpsWGS84 unitsm no_defs4.2 用 gdalbuildvrt 拼瓦片再交给 gdalwarp 重投影全国 DEM 是按图幅分散存储的直接一个个打开再另存太慢。更高效的做法是先做一个 VRT 虚拟拼接让 GDAL 在逻辑上把这些分幅当成一个连续文件然后做一次代真正的重投影输出成统一坐标系、统一分辨率的 GeoTIFF。VRT 本身不复制像素数据创建速度快到可以忽略真正的耗时时在 warp 阶段。gdalbuildvrt dem_raw.vrt raw/dem/*/*.tif gdalwarp -t_srs $SRC_PROJ \ -tr 30 30 \ -r cubic \ -srcnodata -32768 \ -dstnodata -32768 \ dem_raw.vrt dem_aea_30m.tif这里-tr 30 30指定输出分辨率为 30 米-r cubic用三次卷积重采样比最近邻更平滑适合连续高程表面。-srcnodata和-dstnodata同时指定把原始数据中的无效像素值在整个处理链里固定成-32768避免重投影后出现黑块。全部参数的作用是把多幅 DEM 拼接为一个覆盖项目区的高程底图并把无数据值统一固化下来方便后续在做坡度或者流域分析时排除。4.3 用项目边界裁剪 DEM并顺手生成一张山体阴影图DEM 重投影完成后项目连边界还不需要整个全国数据参与后续分析直接裁剪到外扩范围即可。裁剪推荐用掩膜方式让 gdalwarp 读取边界矢量作为裁切线输出栅格的范围和边界完全吻合不需要二次对齐。gdalwarp -cutline work/boundary_aea.shp \ -crop_to_cutline \ -of GTiff \ -srcnodata -32768 \ -dstnodata -32768 \ dem_aea_30m.tif work/dem_cut.tif gdaldem hillshade -z 2 -az 315 -alt 45 work/dem_cut.tif work/hillshade.tif第一条命令把全国 DEM 裁剪到boundary_aea.shp的范围并保留-32768无数据值。第二条用 gdaldem 生成山体阴影图-z 2是垂直方向夸张系数能让微地形更明显-az 315和-alt 45模拟西北方向的光照和 45 度太阳高度角。山体阴影不参与定量计算但它是后期做图时最实用的立体底图能快速看出地形起伏是否合理。4.4 地貌与地质矢量的重投影和裁剪一条 Geopandas 命令走完矢量数据处理比栅格简单但要注意重投影和裁剪的顺序。先重投影到目标坐标系再做空间相交能避免在经纬度坐标系下算错面积。用 Geopandas 做这套操作代码量很短而且每一步都有清晰的中间结果可以检查。import geopandas as gpd target_crs projaea lat_125 lat_247 lat_00 lon_0105 x_00 y_00 ellpsWGS84 unitsm no_defs # 读取项目边界并统一到目标投影 boundary gpd.read_file(raw/landform/boundary.shp) boundary boundary.to_crs(target_crs) # 读取地貌图做投影转换 landform gpd.read_file(raw/landform/landform_1m.shp) landform landform.to_crs(target_crs) # 按边界做空间裁剪保留与边界相关的所有要素 clipped gpd.overlay(landform, boundary, howintersection) # 输出为 GeoPackage图层名保持可读 clipped.to_file(work/landform_clipped.gpkg, layerlandform, driverGPKG) print(裁剪后要素数:, len(clipped))这段代码的处理逻辑是先把边界和地貌图都投影到目标坐标系再用gpd.overlay做面面相交最后输出为 GeoPackage 单文件。需要注意to_crs参数用的是和前期 DEM 完全相同的投影字符串保持一致才能让矢量与栅格在出图时严丝合缝。如果你是重复跑不同项目建议把投影字符串单独定义一次放到脚本头部减少手滑出错的风险。5. 避坑坐标系错位、黑格、高程基准混用与精度翻车5.1 图层叠加错位几百米方向对了但位置不对现象地貌面和 DEM 的山脊线在 QGIS 里叠加看起来“大致方向对”但边界位移明显有时几百米有时一两公里。原因只有两类一类是原始数据的坐标系不一致两个文件各自带了一套不同的投影参数另一类是某个 Shapefile 丢了.prj投影文件GIS 软件默认按 WGS84 经纬度猜测于是整体偏到海里去。解决加载图层后立即查看属性里的坐标系信息不匹配就先用to_crs统一再做任何分析。如果遇到没有.prj的文件先不要乱猜去数据来源页面找投影说明找不回来就对已有点位或者等高线做空间配准用 3 个以上控制点验证。铁律是裁剪和分析之前先统一坐标系不要到最后做验收时才去对坐标那是给自己留后悔药。5.2 DEM 文件加载后一片黑砸掉一整个处理流程现象裁剪后的 DEM 在 QGIS 里显示全黑或者某个区域明显凹陷成一个黑洞用查看工具点下去全是同一个极端值。原因原始 DEM 的无数据值没有被正确识别。常见无数据值有-32768、-9999、0这三种如果软件把无数据值当成真实高程参与渲染就会占据整个数据范围让正常地形全部被压缩到尾段。解决在前期体检时就要确认 NoData 值处理链里贯穿使用-srcnodata和-dstnodata。如果文件已经裁坏用gdal_fillnodata对黑块做插值补洞更稳妥的做法是回到重投影那一步重新输出。以后所有下游产品都会继承坏值所以黑块一定要在源头处理掉而不是等最后出图时再打补丁。5.3 同一片场地DEM 高程和实测高程差了 30 米现象现场用 GNSS 测了几个点拿来和 DEM 提取的高程对比发现系统性偏差 2030 米且方向一致。原因这不是 DEM 误差而是高程基准口径不一致。SRTM 和 ALOS 默认使用 EGM96 大地水准面成果而部分后处理 DEM 或 3D 建模提取的高程是 WGS84 椭球高两者在全国多数地区有十几到几十米的系统差。解决做高程对比前先查数据源元数据里的垂直基准说明统一换算到同一个高程基准再比较。方法是在项目区找一个等级点或者连续运行参考站的已知海拔算出一个常数改正量应用到整个 DEM 上。这一步看起来不起眼但直接影响后续汇水分析和挖填方估算属于那种项目做一半才发现就得全重跑的黑匣子问题。5.4 小比例尺地质图的断层线你把它当现场鉴定用了现象1:100 万地质图上的断层从项目场地正中间穿过现场拉剖面却找不到任何断层迹象最后审查时被质疑数据不可靠。原因小比例尺地质图是区域构造的概化表达断层线通常经过制图综合位置精度在本就有限。它表达的是“存在反向断层”这个地质认识而不是“断层就在这条线上”。解决小比例尺地质图只用来做区域稳定性初判、选址红线避让和勘察方案布置。进入详细勘察阶段必须去收集项目区 1:5 万或更大比例尺的区域地质图再配合野外踏勘验证。凡是大断裂通过的位置现场至少要安排物探或钻探控制点用实据去证明或否定图上认识。5.5 全国 DEM 一次性加载机器开始原地转圈现象把整块全国 DEM 加进工程后视图缩放一卡就是 30 秒内存占用一路飙升。原因没有对栅格做金字塔优化GIS 软件被迫读取全分辨率像素来做缩略显示文件越大越明显。这不代表电脑配置不够而是数据组织方式不合适。解决用gdaladdo给重投影后的 DEM 生成金字塔概览命令是gdaladdo -r average dem_aea_30m.tif 2 4 8 16 32让软件在缩小时直接读取低分辨率概览。同时将工作区内的中间成果统一放在 SSD 上避免通过慢速网络盘直接拖拽大数据层。这套组合拳基本能解决绝大多数“大 DEM 卡到爆”的现状。提示只要你用的是全国尺度数据处理前先建金字塔总是没错的。这比后期反复调机器性能省心得多属于投入一分钟、后期少挨半个小时的稳赚操作。6. 进阶用现场实测点位反向验证这套全国数据集靠不靠得住数据集整理完不等于工作结束用真实点位验证才是让成果站住脚的关键。做法很简单在项目区选一批有代表性的现场点位记录坐标和实测高程再从处理好的 DEM 里提取对应坐标的高程值做对比。顺带把点位所在位置的地貌类型和地质属性也一并取出和现场观察结果对照看数据集在项目区到底可用到什么程度。下面这段脚本完成核心动作读取实测点矢量批量从 DEM 提取高程并把地貌层的类型属性关联到点上。import rasterio import geopandas as gpd target_crs projaea lat_125 lat_247 lat_00 lon_0105 x_00 y_00 ellpsWGS84 unitsm no_defs # 读取实测点并投影到与 DEM 一致的坐标系 pts gpd.read_file(work/field_points.shp) pts pts.to_crs(target_crs) # 读取裁剪后的地貌面数据 landform gpd.read_file(work/landform_clipped.gpkg) landform landform.to_crs(target_crs) # 先做空间连接把所在面要素属性带回到点 pts gpd.sjoin(pts, landform, howleft, predicatewithin) # 再从 DEM 逐点采样高程 elevations [] with rasterio.open(work/dem_cut.tif) as src: for _, row in pts.iterrows(): try: value list(src.sample([(row.geometry.x, row.geometry.y)], 1))[0][0] elevations.append(value if value ! -32768 else None) except Exception: elevations.append(None) pts[提取高程] elevations # 输出对比表供后续统计中误差 pts[[点号, 实测高程, 提取高程, 地貌类型]].to_csv(work/validate_result.csv, indexFalse) print(pts.head())这段脚本先做空间连接把地貌类型匹配到点上再逐个点位从 DEM 提取高程。运行完后输出一个 CSV里面同时有实测高程、提取高程和地貌属性。把提取高程减去实测高程按点位数统计出平均差和中误差如果偏差在项目允许范围内这张底图就可以放心交出去如果系统性偏差过大就该回去检查高程基准和原始 DEM 分辨率是否真的适配。多年做区域数据整合下来我最大的感受是全国数据集不怕信息量少就怕你没搞清楚每一层的口径就叠在一起用。宁可多花半天做基准统一点位验证也别让一张错位的底图带着整个项目翻车。希望这套方法帮你在下次开工时少踩几个坑。本文还有配套的精品资源点击获取
返回列表