☰
ArcMap滑坡易发性评价:数据预处理与栅格计算器全流程
2026/10/15 4:23:31 网站建设 项目流程

1. 从一张乱糟糟的DEM开始说起

做滑坡预测的人都知道,ArcMap里最磨人的不是模型本身,而是数据预处理那一关。前面三篇笔记我写过关于图层配准、矢量化和基础地形因子提取的操作,这篇直接进入正题:怎么把多源数据整理成一张能用的栅格集,然后交给信息量模型去跑。

先说个真实场景。某次做某山区滑坡易发性评价,手上数据包括12.5米分辨率的DEM、1:50000的地质图、多年平均降雨量站点数据、Landsat影像反演的NDVI,还有一份土地利用现状分类图。这些东西各有各的坐标系,各有各的栅格尺寸,有的还是矢量格式。如果直接扔进栅格计算器,出来的结果要么是一片空值,要么是范围对不齐导致的错位。

这篇笔记就是围绕这个场景来写的。核心内容是把ArcMap里最常用的几个空间分析工具串成一条完整的生产流程:统一环境、重采样、重分类、栅格计算器加权叠加,再加上中间最容易翻车的几个细节处理。适合正在做滑坡易发性区划、地质灾害危险性评价这类工作的同行参考,也适合选修地信方向课程、正在为毕业设计攒数据的同学按步骤复现。

2. 动手之前先理解为什么要统一坐标和像元

很多教程上来就教操作,但不说为什么。操作本身不难,难的是出了问题之后不知道怎么排查。所以先花点篇幅把原理铺开。

2.1 多源数据为什么不能直接叠加

做滑坡预测需要用到地形因子、地质因子、植被因子、降雨因子,这些因子来自完全不同的数据源。DEM是栅格,地质图是矢量面,降雨量是站点插值出来的栅格,NDVI是遥感反演栅格。它们各自的坐标系可能不同,有的用地理坐标系(经纬度),有的用投影坐标系(平面坐标)。

ArcMap里栅格计算器要求所有参与计算的图层必须严格对齐,不仅是坐标系一致,还包括像元大小、像元位置、栅格范围都要一致。否则计算器会在边缘位置产生大量NoData,而NoData参与四则运算的结果永远是NoData。最后你得到一张绝大多数区域都是空值的"成果图",根本没法用。

这一点用生活类比来说,就像把一份使用公制单位标注尺寸的零件图,和另一份使用英制单位的装配图强行拼在一起加工,尺寸对不上是必然的。

2.2 坐标系的选择逻辑

滑坡预测项目里最常见的坐标系设定是:原始DEM如果已经是投影坐标系,保持不动;如果是地理坐标系(比如WGS84经纬度),就需要投影到适合研究区的投影坐标系统。

选择投影坐标系有个基本原则:研究区所在纬度越高,越不宜使用等距圆柱投影,会产生明显的面积和距离变形。国内区域通常使用高斯-克吕格投影的CGS2000或西安80、北京54版本。实际操作中我的习惯是:先把所有矢量图层和栅格图层全部投影到同一投影坐标系,然后统一转成栅格再做因子计算。

这里补一个细节,ArcMap中的"投影栅格"工具与"定义投影"工具是两个完全不同的功能。很多新手把"定义投影"当作转换坐标系来用,结果只是给数据打了一个错误的坐标标签,后面全乱套。定义投影是告诉软件"这个数据是什么坐标系",投影栅格才是真正把数据从A坐标系换算到B坐标系。

3. 因子图层的生产流程与参数设定

这部分是全文的实操核心。我会按实际项目推进的顺序来写,每一步都给参数和理由。

3.1 第一步:统一投影与裁剪范围

打开ArcToolbox,依次执行"数据管理工具 → 投影和变换 → 栅格 → 投影栅格"。输入原始DEM,输出坐标系选择你提前建好的投影坐标系,重采样技术默认选择双线性插值(Bilinear)。对于DEM这类连续表面数据,双线性插值比最邻近法更平滑,不会出现明显的锯齿边界。

研究区范围裁剪也是在这一步完成。在环境设置里找到"处理范围",选择"与图层相同",然后点选你准备好的研究区边界矢量层。这一步的目的不是让栅格变好看,而是让所有后处理都在同一范围内计算,避免边缘效应。从这一步开始,后面的每个工具执行前我都建议先在环境设置里统一四件事:当前工作空间、处理范围、像元大小、坐标系。

3.2 第二步:基础地形因子提取

坡度(Slope)、坡向(Aspect)、地形曲率(Curvature)这三个因子是滑坡预测里面出场频率最高的地形指标,直接从DEM上提取。

坡度工具在"3D分析工具 → 栅格表面 → 坡度",关键参数是输出测量单位。默认是度(Degree),但在某些模型里你需要的是坡度的百分比,两种单位在数值分布上差异很大,会影响后面重分类的阈值选择。百分比的公式是tan(坡度角)×100%,平地为0,45度坡为100%,90度接近无穷。我习惯保存为度,后续重分类时直接按度数区间划分,比较好解释。

坡向工具的参数相对简单,但要注意输出栅格中平地的值为-1,这一部分在后续重分类时要单独处理,不能和其他坡向类别混在一起。坡向的环形分类通常分成9类:平地、北、东北、东、东南、南、西南、西、西北。

地形曲率工具的输出包含剖面曲率和平面曲率两个波段。剖面曲率沿最大坡度方向计算,正值表示坡面下凸,负值表示上凹,和滑坡的汇水、启滑位置有直接关联。实际操作中我会把平面曲率也保留下来,因为它在后续做沟道提取时经常用到,一次算完就不必返工。

3.3 第三步:地质岩性与植被因子处理

地质图是矢量数据,需要转成栅格。这里最需要注意的是字段选择。属性表里通常有地层代号、岩性名称、地质年代等多个字段,转栅格时要选最符合因子分级目标的那个字段。滑坡预测中常用的是岩性分组,不是具体地层代号,因为代号太细会导致类别数过多,给后续重分类和权重打分带来麻烦。

转栅格工具在"转换工具 → 从矢量到栅格 → 多边形转栅格",字段选择岩性分组字段,像元大小保持和DEM一致。如果原始地质图的比例尺较粗,转出来的栅格会有大量锯齿边界,这是正常现象,不必为这个去过度平滑。

NDVI因子直接从Landsat影像用栅格计算器算出来:NDVI = (近红外 - 红光) / (近红外 + 红光)。计算前先确认影像波段顺序,避免把红波段和近红外波段搞反。波段反了,NDVI值域会变成负值为主,这样的植被因子拿来直接参与滑坡预测,可能得出"植被越稀少越安全"这种背离常理的结论。

3.4 第四步:降雨因子的插值细节

降雨量站点数据通常是点状矢量,要变成连续栅格面,需要插值。ArcMap空间分析工具里提供多种插值方法,滑坡预测最常用的两个是反距离权重法(IDW)和普通克里金法。

对于站点数据较少且分布不均的研究区,我的经验是IDW比克里金更容易控制参数,结果也不会差太多。IDW的核心参数是幂指数(Power),默认值是2,代表距离衰减的强度。幂指数越大,距离越远的站点对插值点的影响越小,结果越"尖锐"。站点密度低时建议用幂指数2,站点密度高时可以用3,但不宜再大,否则会出现明显的以站点为中心的"靶心"效应。

插值的另一个关键参数是输出像元大小和搜索半径。输出像元大小必须与其他因子图层一致,否则后续叠加还是需要重采样。搜索半径设置为变量,让软件自动根据站点分布决定参与插值的站点数量,比固定半径更可靠。

4. 重分类:把所有因子放在同一把尺子上

4.1 为什么非重分类不可

滑坡预测的信息量模型要求每个因子图层都是整数型的类别栅格。为什么?因为信息量计算本质上是统计每个类别内滑坡点的比例,栅格值必须是整数才能对应到具体的类别编号。如果你直接把坡度这类连续值栅格扔进去,模型程序无法区分"坡度15.3"和"坡度15.7"这两个像元属于同一类别,统计结果就会混乱。

重分类的另外一个作用是去量纲。坡度是度数,NDVI是-1到1之间的无量纲数,降雨量是毫米,这些数值量级完全不同,直接加权求和的物理意义不明确。通过重分类,把每个因子的原始值映射到统一的等级序列上,之后的计算才有意义。

4.2 重分类操作与分级标准

利用重分类工具,将坡度按易发性经验阈值划分为5级。分级标准没有唯一答案,但有一条原则:每个级别内滑坡点数量不能为零,否则该级别的信息量值为负无穷,计算会报错。如果发现某个分级区间内没有滑坡点,处理方式有两种:一是调整分级边界,让区间内至少包含一定数量的滑坡点;二是把该区间并入相邻区间。

以坡度因子为例,我会先对研究区内的滑坡点做叠加统计,查看每个坡度分级的滑坡点占比,再据此调整阈值。某次项目中,原始分级把0到5度作为第一级,但统计结果显示该区间内滑坡点为零,最后我把它和5到10度合并。这个细节看起来小,实际操作中经常因为偷懒忽略掉,导致模型输出结果异常。

NDVI分级相对灵活,通常按植被覆盖度从低到高分为5级。土地利用分类图如果有现成的分类编码,可以不重分类,但要检查编码是否连续,跳号会导致某个类别在后续叠加时被当成NoData。

4.3 关于重分类表的一个易错点

重分类工具的界面上有一张重分类表,左边是旧值,右边是新值。最容易犯的错误是忘记修改"NoData"那一行的设置。默认情况下NoData保持为NoData,如果你希望NoData区域参与运算并赋值为0,需要手动把新值改成0。但在滑坡预测中,我更推荐保持NoData为NoData,因为参与运算的像元必须是研究区范围内的,NoData区域本来就不该有值,如果强行赋0,会人为压低该区域的易发性等级,造成结果失真。

5. 栅格计算器实现信息量模型的全流程

5.1 信息量法的数学本质

信息量模型的计算公式在教科书里有一长串推导,落到ArcMap里其实就两个步骤:先算单因子的信息量值,再对所有因子做加权求和。

单因子信息量值的核心思想是:比较某个因子类别内滑坡发生频率与研究区整体滑坡发生频率的比值。比值大于1,说明该因子类别有利于滑坡发生,信息量为正;比值小于1,说明不利,信息量为负。

落到实操中,这个比值用滑坡点数占比来近似。在重分类后的因子图层上叠加滑坡点,统计每个类别内的滑坡点数,再统计研究区总滑坡点数,两者相除后取对数,得到该类别的信息量值。

5.2 用连接表给类别赋权重

具体操作路径是这样的:先用"制图工具 → 要素转栅格"把滑坡点转成栅格,值为1,像元大小与其他因子一致。然后用"空间分析工具 → 叠加分析 → 栅格代数"或直接栅格计算器,对每个因子图层执行Con条件判断,实现"类别编号→信息量值"的映射。

这里我补充一个高效率的替代方案:不用逐个归类写Con嵌套,而是把信息量值表做成一张独立属性表,用"连接"功能挂到因子栅格上。听起来有点绕,但实际操作中能节省大量时间,特别是因子数量较多(比如8个因子以上)时,嵌套式Con的表达式长到让人崩溃。

具体做法是,基于重分类后的栅格,以栅格值作为连接字段,匹配一张包含"类别与信息量值对应关系"的表格。连接完成后,利用"复制栅格"工具或转换到整型栅格,把信息量值固化到栅格属性表里。

5.3 最终加权叠加的计算

所有因子都转成信息量值栅格后,最后一步是栅格计算器加权求和,表达式形如:

"InfoResult" = 0.25 * "InfoSlope" + 0.20 * "InfoAspect" + 0.15 * "InfoCurvature" + ...

权重系数从哪来?最稳妥的做法是层次分析法(AHP)构造判断矩阵。实际操作中判断矩阵的构造有一定主观性,我的建议是:权重不必过度追求精度,但方向不能违背常识。比如在降雨量大的湿润地区,降雨因子的权重不应低于地形因子;在构造活动强烈的区域,地质因子的权重可以适当提高。

权重计算出来后,建议先做一次一致性检验,CR值小于0.1才说明判断矩阵可接受。ArcMap本身没有内置AHP工具,你可以用Excel手算,或者找现成计算脚本。我在项目里通常是先在Excel里算好权重,再拿到ArcMap里写表达式,避免在GIS里反复试错。

5.4 输出结果的栅格渲染与分级

加权叠加的结果是一个连续值的栅格,数值范围可能从负值到正值都有,直接出图很难看,也不符合易发性评价的表达习惯。我会用"重分类工具"把结果划分为四级:低易发、中易发、高易发、极高易发。

分级方法推荐用自然间断点分级法(Jenks),它会在数据分布的自然断点处切分等级,比等间距分级更符合数据本身的分布特征。在ArcMap图层属性里设置分类时选择分类方法为"自然间断点分级法",类别数设为4,即可得到统计上较合理的分区效果。

6. 常见问题与排查技巧实录

6.1 栅格计算器结果全为NoData

这个问题出现频率最高。原因通常是参与运算的图层像元位置没有对齐。排查步骤如下:

  • 检查所有图层的属性 → 源 → 栅格信息,确认像元大小一致、范围一致。
  • 如果不同,返回重采样工具,将像元大小统一(比如全部设为DEM相同分辨率),并在环境设置里固定处理范围。
  • 确认坐标系一致,否则先投影转换再重采样。

6.2 重分类中某类别没有滑坡点

前面提过,这会导致该类别信息量值为负无穷。直接报错还好,最怕的是不报错,但结果里出现异常极值。处理办法是检查统计表,手动合并无滑坡点的类别,重新重分类。不要试图给该类别的信息量值强行设一个0或极小值,这样会影响后续加权叠加的数值平衡。

6.3 插值结果出现负值或异常值

IDW插值的结果在理论上不会超出样本数据范围,但仍可能出现局部异常值。常见原因是站点数据本身包含异常记录,比如气象站的降雨总量在某个年份因维护异常出现跳变。处理方式是在插值前对站点数据进行探索性数据分析,绘制直方图和散点图,剔除明显异常点后再插值。

6.4 滑坡点与因子栅格位置错位

这个坑最隐蔽,因为表面上坐标系一致、范围一致,但错位是亚像元级别的。导致错位的常见原因是最初的配准误差累积。滑坡点在野外用手持GPS采集时本身就有几米到十几米的误差,在高分辨率DEM下就会体现为点位轻微偏移。解决办法是在进行滑坡点缓冲分析时,生成一个合理半径的缓冲区(例如10米),再做叠加统计,而不是直接让点和栅格像元做点对值提取。

6.5 工具执行报错"ERROR 010240"

这是ArcGIS空间分析工具在处理大范围高分辨率栅格时的内存错误。处理办法是在地理处理选项里启用后台处理并增加临时文件空间,或者将研究区分块处理后再拼接。某次处理某跨平台系统数据时,整个研究区DEM有6000×8000个像元,直接处理总是报错,最后只能分区计算再合并。

7. 我个人在实际操作中的几点体会

做了四年滑坡易发性相关的项目,中间踩过的坑远不止上面列出的这些。有几条经验是我现在每次做项目都会提醒自己的,放在这里作为补充。

重采样方法的选择要因数据类型而异。坡度、降雨量这类连续表面用双线性或三次卷积;岩性、土地利用这类离散类别数据必须用最邻近法,否则会让类别边界出现"过渡像元"。这个细节不处理好,后面重分类时会出现莫名其妙的中间类别值。

权重方案不要追求一次算对。AHP判断矩阵构造时,如果一致性检验不通过,不要反复微调一个格子的数值,更好的做法是重新审视整个矩阵的层级逻辑。经常出现的情况是:某两个因子之间的重要性关系本身就设反了,微调数值治标不治本。

出图时注意设置适当的符号化方案。滑坡易发性图的配色要遵循"冷色代表低易发、暖色代表高易发"的视觉直觉,如果颜色顺序反了,评审时很容易被指出问题。这是个小细节,但在项目汇报中直接影响成果的专业感。

最后一条建议关于数据管理。每个中间结果都要保留清晰的命名规范,建议统一格式为"因子名_处理阶段_分辨率",例如"Slope_reclass_12p5m"。这个习惯在我自己同时处理多个研究区时帮了大忙,少了无数重复劳动。

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

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

立即咨询