企业官网建设流程全解析

做水土流失模拟这几年，被问得最多的问题不是“模型怎么跑”，而是“数据怎么凑”“因子怎么算”。很多人拿到ArcGIS、InVEST、RUSLE这三个词就开始套模型，结果导出来的图和实测站点对不上，要么是侵蚀量离谱，要么是栅格一片空白。这篇就把我从数据预处理到InVEST SDR结果解读的完整流程盘一遍，重点说那些实际操作中绕不开的坑和对应解法，给正在做区域土壤侵蚀评估、生态服务功能评价或者国土空间规划的同行一个能直接参考的路线。

1. 水土流失模拟的两种主流路线：RUSLE与InVEST SDR怎么选

1.1 RUSLE模型的三层理解

RUSLE（Revised Universal Soil Loss Equation）本质上是个经验统计模型，算的不是某一次暴雨的侵蚀量，而是多年平均的年土壤侵蚀量。它的数学形式很简单：

A = R × K × LS × C × P

A是年平均土壤流失量（t/(hm²·a)），R是降雨侵蚀力因子，K是土壤可蚀性因子，LS是坡度坡长因子，C是植被覆盖与管理因子，P是水土保持措施因子。六个变量乘在一起，每个因子都对应一个空间栅格，最后在ArcGIS的栅格计算器里一乘，就得到一张侵蚀强度分布图。

这个模型之所以在科研和实际项目中用得广，是因为它把复杂的水土流失过程压缩成了六个相对容易获取的参数。但代价是它不考虑泥沙在坡面、沟道中的输移过程。也就是说，RUSLE算出来的是“某个地块一年平均被冲走多少土”，至于这些土有没有进到河道、在哪个位置沉积下来，它管不了。

1.2 InVEST SDR模块与RUSLE的关系

InVEST的SDR（Sediment Delivery Ratio，泥沙输移比）模块，底层用的还是RUSLE的思路，但它往前多走了一步——在计算每个像元的侵蚀量之后，通过一个基于地形连通性的输移比模型，把坡面侵蚀量换算成进入河道的产沙量，还能算出泥沙在空间上的拦截与保留量。

SDR模块的核心逻辑是这样：先用USLE类的方程计算出每个栅格的潜在侵蚀量，然后计算一个连通性指数（IC），这个指数反映了坡面物质到达河道的能力。地形越陡、坡面越连续、径流路径越短，IC值越高，泥沙越容易进河。最后用SDR值把侵蚀量和产沙量联系起来。SDR的取值不是固定常数，而是根据IC动态计算出来的。

所以你可以这么理解：RUSLE是算“哪里在掉土”，InVEST SDR是算“掉了多少土被送进了河”。做区域水土流失防治规划，两个都有用，但回答的问题不一样。

1.3 我的选型建议

结合我自己的项目经验，选型的判断标准可以简化成两句话：如果只是想画一张土壤侵蚀强度分级图，给水土保持区划提供底图，那就老老实实算RUSLE，工作量小、解释起来也容易；如果项目要求算入河泥沙量、评估不同土地利用情景下的产沙变化，或者要对比生态修复前后的泥沙拦截效果，那就必须上InVEST SDR。

另外提醒一句：InVEST SDR和RUSLE在数据要求上有大量重叠，DEM、降雨、土壤、土地利用都是共用的。哪怕最终目标是跑InVEST，也建议先把RUSLE的六个因子完整算一遍。因为SDR模块里的不少参数需要你先理解RUSLE因子的空间分布特征，跑完RUSLE再看InVEST的输出，才分得清哪些变化是泥沙输移造成，哪些是侵蚀量本身在变化。

2. 环境配置与数据准备：别在起跑线上耗时间

2.1 ArcGIS版本选用与常见启动问题

水土流失模拟涉及的栅格计算、水文分析、插值，用ArcMap 10.8或者ArcGIS Pro 3.x都能完成。我个人习惯是ArcMap 10.8做因子计算，InVEST模型单独跑，最后用ArcGIS Pro做制图和数据管理。

很多人在ArcGIS安装环节就把时间耗完了。启动许可没反应，十有八九是License Manager服务没起来，或者被杀毒软件拦了。Windows 10系统下，打开服务管理器，找到ArcGIS License Service，确认状态是“正在运行”。如果启动失败，把安装目录下的license文件删掉重新激活一次。还有更常见的情况是.NET Framework 3.5没装全，ArcMap 10.x系列对这个依赖很强，装完系统补丁之后务必在“启用或关闭Windows功能”里把.NET Framework 3.5勾上。这些步骤做完，启动问题基本能解决百分之八十。

2.2 数据清单与来源

跑一个完整的水土流失模拟，至少需要四类数据：

• DEM数字高程模型：分辨率建议10~30米。国内可以用地理空间数据云的ASTER GDEM或SRTM数据，如果做县域级别，建议用12.5米的ALOS PALSAR数据，地形细节远好于30米。裁剪到研究区范围，统一投影坐标系。
• 降雨数据：至少找个三到五年的各站点月降雨量，用于推算降雨侵蚀力R值。站点数量少于十个的话，后面插值会很难受。
• 土壤数据：需要土壤类型图或者至少能查到土壤的砂粒、粉粒、黏粒、有机碳含量。国内有不少省级土壤普查成果数据，实在没有就退一步用HWSD（世界土壤数据库），但精度要打折，表述上一定要写清楚数据来源。
• 土地利用/覆盖数据：GlobeLand30、FROM-GLC或者自己解译的遥感分类结果都行。关键是要能对应上C因子和P因子的赋值标准。

2.3 统一投影与分辨率：一切计算的前提

这是整个流程里最容易被轻视的一步，但也是出错率最高的一步。RUSLE的LS因子涉及坡长和坡度运算，如果用经纬度坐标系直接算，哪怕是WGS84也不行，因为面积、长度单位在不同纬度下全乱了。所有栅格数据必须统一到一个投影坐标系下。国内项目我一般用CGCS2000 / 3-degree Gauss-Kruger zone或者UTM对应分带，统一到中央经线附近，变形的幅度最小。

统一到投影坐标系之后，还要统一分辨率。常用的做法是设定成30米，用最邻近法重采样土地利用数据，用双线性插值重采样DEM和降雨栅格。土壤数据按矢量转栅格时的像元大小对齐。我见过有人把土地利用是30米、DEM是90米的数据直接扔进模型，InVEST跑完结果奇形怪状，排查半天才发现是分辨率不一致导致的。

3. 最难算对的LS因子：ArcGIS实操全流程

3.1 DEM预处理：填洼与流向

LS因子是整个RUSLE模型中对空间计算要求最高的部分，也是最容易算错的。第一步是填洼（Fill）。原始DEM里存在大量凹陷区域，如果不填，水流会在洼地中断，流向计算和坡长累积都会出错。ArcGIS里的工具路径是Spatial Analyst → Hydrology → Fill。

但填洼不能无脑填。有些区域的洼地是真实地形，比如喀斯特地貌、采石坑，这种真实洼地一旦填平，坡长计算就会严重偏高。实操中我的习惯是先看一眼填洼前后的DEM差值，如果某个区域的填高超过几十米，十有八九是数据本身有错误或噪声，需要用Focal Statistics对DEM做一次中值滤波，或者用原始等高线重新插值DEM来修正。

3.2 流向与坡长累积计算

填洼完成后，按顺序执行Flow Direction、Flow Accumulation。Flow Direction用的是D8单流向算法，每个像元的流向指向周围八个像元中坡度下降最陡的那个。ArcGIS里生成的方向编码是1、2、4、8、16、32、64、128，不是常规角度的东南西北，这个编码后续不会直接用到，但检查结果时要知道它代表什么，避免把方向栅格当普通栅格处理。

Flow Accumulation得到的是每个像元上游有多少个像元汇入。这个值乘以像元面积，就是水文意义上的汇流面积。LS公式里的“坡长”并不是从山脊到坡脚的实际坡面长度，而是用汇流面积来近似单位宽度上的累积径流量。这是RUSLE在实际空间计算中一个很重要的替代处理，理解这一点，你才能明白为什么需要用Flow Accumulation来算坡长。

3.3 栅格计算器里LS公式的写法

LS因子在ArcGIS里的常用计算式是Moore和Burch提出的简化形式。公式长这样：

LS = Pow(（FlowAcc × CellSize）/ 22.13, 0.4) × Pow(Sin(SlopeRad) / 0.0896, 1.3) × 1.4

其中FlowAcc是汇流累积量，CellSize是像元尺寸（米），22.13是标准径流小区坡长（米），0.0896是标准坡度（9%对应的弧度值），1.4是前人修正系数。

实际操作中我拆成多步来做，便于检查每一级的输出：

第一步，用Spatial Analyst的Slope工具提取坡度。注意输出类型选“DEGREE”（度）。但栅格计算器里的Sin函数要求弧度，所以必须先把角度值乘以π/180转成弧度，Scratch的写法就是Slope × 0.017453。

第二步，计算坡长栅格。栅格计算器写：

LangthFactor = Pow（FlowAcc * 30 / 22.13, 0.4）

第三步，计算坡度因子：

SlopeFactor = Pow（Sin(Slope * 0.017453) / 0.0896, 1.3）

第四步，乘在一起得到LS：

LS_Final = LangthFactor * SlopeFactor * 1.4

这里最容易翻车的点有三个。一个是Slope工具输出的是百分数还是度数，默认是度数，但很多人设成Percent就忘改了，这时候Sin（Slope）算出来完全不是那么回事。二是分辨率。FlowAcc乘30的前提是分辨率正好30米，如果用的是12.5米数据，就要乘12.5。三是栅格计算器里英文字段名不能有中文和空格，多重命名搞混的情况我见过太多次。

如果只是做粗略分析，也可以用贴块查表法代替连续LS计算。方法是把坡度分好几级，每一级给一个经验LS值，再用重分类赋值。这个办法在数据精度不足时反而更稳健，但会丢失地形细节，做精细化评估时我不建议用。

4. R因子与K因子：插值与查表的正确姿势

4.1 降雨侵蚀力R的计算与空间插值

R因子代表降雨对土壤的剥离和搬运能力，在年尺度上它与降雨量、降雨强度密切相关。完整的EI30算法需要逐场降雨过程数据，绝大多数项目拿不到，最常用的是用经验公式从年降雨量或月降雨量推算。国内不同省区各自发表了不少率定公式，使用时要搞清楚自己的研究区适合哪一组参数，避免跨区域硬套。

比如有一种常见形式是R = 0.1834 × P^1.4269，P是多年平均降雨量（mm）。这类公式是用特定地区的降雨资料率定出来的，换个气候区就可能失真严重。

空间插值这一步是R因子的重头戏。几十个气象站点的点数据要变成覆盖整个研究区的连续栅格，多数人第一反应是克里金。如果你的站点少于十五个，克里金的变差函数拟合基本就是摆设，出来的结果不稳定。实测下来，站点少的时候IDW（反距离权重）和样条函数反而更稳，误差可控。站点够多而且分布均匀，克里金能给出带置信区间的、更平滑的结果。

插值之后记得做一步交叉验证，看看实测值和插值结果的均方根误差有多大。不要只盯着好看的成图，一幅平滑得像绸缎一样的降雨分布图往往是过度平滑的产物，真实降雨的空间异质性比你看到的大得多。

4.2 K因子查表法与EPIC公式

K因子反映土壤抵抗侵蚀的能力，它的取值取决于土壤质地、有机质含量、渗透性这些理化性质。最省事的办法是查表：拿到土壤类型图，把土种对应的K值列表，国内常见土壤类型的K值在各类文献里都能找到。红壤大概在0.22到0.35之间，黑土偏高，砂质土偏低。

如果你手里的土壤数据包含机械组成信息，可以用EPIC公式精细计算：

K = 0.1317 × [0.2 + 0.3exp(-0.0256 × m_s × (1 - m_silt/100))] × [m_silt/(m_clay + m_silt)]^0.3 × [1 - 0.25 × orgC/(orgC + exp(3.72 - 2.95 × orgC))] × [1 - 0.7 × (1 - m_s/100)/((1 - m_s/100) + exp(-5.51 + 22.9 × (1 - m_s/100)))]

其中m_s是砂粒含量百分比，m_silt是粉粒，m_clay是黏粒，orgC是有机碳含量百分比。公式算出来之后乘以0.1317转换成国际单位。

实际操作中，土壤数据常常是矢量多边形。用Polygon to Raster转栅格时，要选K值字段，像元大小和其他栅格对齐。如果一张土壤图里同一个图斑套了好几个土壤亚类的属性，先用Dissolve合并碎斑，再赋值，否则转栅格会很碎，输出结果和底图完全对不上。

5. C因子和P因子：土地利用数据里藏着最容易忽略的细节

5.1 C因子：从土地利用类型到覆盖管理系数

C因子是植被覆盖和田间管理对侵蚀的抑制效应，取值在0到1之间。裸土是1，意味着没有任何保护；茂密森林可以低到0.001~0.01。给土地利用类型赋值这个是常识，但细节在于同一类地类内部差异极大。

比如耕地，玉米地、小麦地、果园的C因子完全不同，刚翻耕过的地和封行后的地也差好几倍。如果项目精度要求高，光靠土地利用类型给固定值是不够的，更可靠的办法是用NDVI反演植被覆盖度，再建立植被覆盖度与C因子的对应关系。常用公式是：

C = exp(-a × FVC / (b - FVC))

FVC是植被覆盖度，a、b是经验参数。算出来C的连续栅格在空间上比单纯查表赋值平滑得多，更接近真实情况。但要注意，这个公式对水体、建设用地、裸岩的处理很别扭，需要先把非植被区单独掩膜赋值。

5.2 P因子：这套流程里最主观的一环

P因子是水土保持措施对侵蚀的削减作用，传统上按措施类型赋值。等高耕作约0.6，等高带状种植约0.4到0.5，梯田约0.15到0.3，无措施为1.0。问题在于，P因子的空间分布很难从公开数据中获得。大范围项目里几乎拿不到每一块地用了什么措施的矢量图，基本都是按土地利用类型粗给一个固定系数，甚至全部取1.0。

这里我说句实在话：如果做的是几万平方公里的区域评估，P因子统一取1.0是可以接受的，因为区域尺度上坡面措施类型数据缺失是客观现实。但论文里或者评审阶段一定会有人问，必须如实交代这个假设。如果你有一小块实验区的梯田分布图，可以单独提出来做敏感性分析，看看P因子取值变化对总侵蚀量的影响幅度有多大。

5.3 重分类与栅格计算器的整合

把土地利用数据转成C因子和P因子栅格，最顺手的方法是建一张属性映射表，然后用Reclassify或者Lookup工具批量赋值，而不是在栅格计算器里写一堆if嵌套。

注意土地利用类型的字段值是用整数代码还是文本名称。GlobeLand30这类数据是数字代码，先去出一个代码清单，确认每个代码代表什么类型，防止代码错位。赋值完成后看一眼直方图和属性表，C因子栅格的像元值应该都落在0到1之间，如果出现负值或大于1的异常像元，八成是原始分类数据里有未识别的类别被公式代入了错误值。

6. InVEST SDR实操：参数设置与结果解读

6.1 SDR模块需要准备哪些输入

InVEST SDR跑起来之前，要先把输入数据盘一遍。打开模型的SDR选项卡之后，实际操作中必填的项目有：

• DEM高程模型。注意这里要用填洼之后的DEM，InVEST自己不会帮你做水文预处理。
• LULC土地利用栅格。
• Rainfall Erosivity Index，也就是R因子栅格。
• Soil Erodibility，也就是K因子栅格。
• Biophysical Table。这是一个CSV格式的表，里面至少包含lulc_veg（是否属于植被覆盖）、usle_c（C因子）、usle_p（P因子）三列，以土地利用代码为链接键。

很多人第一次跑就卡在Biophysical Table这步。CSV文件的表头必须是英文，urf代码不能有重复，最后一列别带多余空行。inVEST对表格格式的要求很严格，我见过因为CSV里多了个空格导致整个运算报错的，十分浪费排查时间。

6.2 阈值参数和SDRmax的取值逻辑

SDR模块的几个模型参数值得多说几句。ic_param是连通性指数公式里的校准系数，模型给的默认值为0.5。这个值控制着IC对产沙量的敏感程度，调高会让地形连通性差异对产沙的影响更明显。科研文献里常见做法是先用默认值跑一版，再结合实测径流小区的泥沙观测数据做灵敏度率定。如果你手里没有任何实测数据，就别随意大调，默认值在大多数区域都说得过去。

kb参数默认是2，它控制着坡面SDR的衰减速度，最合适的取值与气候和岩性有关。SDRmax默认0.8，表示一个像元最大可能输送到河道的泥沙比例。这些参数在InVEST用户手册里都解释得很清楚，但项目报告中你必须写清楚“模型参数采用默认值，未做本地率定”这句话，才能堵住评审的追问。

6.3 结果栅格与RUSLE输出的对照解读

跑完之后在output文件夹里重点看三个栅格：sed_export（进入河道的泥沙量），sed_retention（被植被和地形拦截保留的泥沙量），SDR（每一像元的泥沙输移比）。

我习惯的做法是把sed_export和RUSLE的侵蚀总量放一起对照。侵蚀量大的地方产沙量不一定大，因为那些区域可能是洼地、湿地或者植被茂密的山谷，泥沙在到达河道前就被拦住了。反过来，陡坡直连河道的地方，哪怕侵蚀量中等，产沙量也可能很高。这种空间错位正是InVEST SDR相比RUSLE最有增值的地方。

如果一个地块的SDR值全是0.8附近，没有空间差异，说明连通性指数算出来的结果太平滑，可能DEM分辨率偏粗或者ic_param取值不合适。这时候回去看看IC那一步的中间结果，确认是不是填洼后地形细节被过度磨平了。

7. 我踩过的那些坑：ArcGIS水土流失模拟避坑实录

7.1 投影坐标系不统一的连锁反应

我之前做过一个跨两个县的区域评估，从省里拿到的土壤数据是WGS84经纬度，自己提的DEM是UTM投影，土地利用是CGCS2000高斯投影。三种坐标系直接叠加之后，表面看着区域轮廓差不多，实际上错位了百来米。LS因子里的坡长计算一旦基于错位数据，栅格乘出来的结果完全是错的，而且这种错误不会报警告，只有把几个因子的栅格叠加对比对齐关系时才看得出来。

现在我的流程固定为：拿到所有数据的第一个动作就是统一投影，用Project Raster和Project工具处理栅格与矢量，全部转换之后再进下一步。宁可多花十分钟，不要在跑完模型之后推倒重来。

7.2 栅格计算器里NoData的传染效应

栅格计算器有个特别坑人的特性：只要参与运算的任何一个栅格在某位置是NoData，结果栅格在那个位置也是NoData。在用C因子、P因子稀疏区域乘起来的时候，一个NoData会像瘟疫一样传染给成片的结果区域。

解决的办法是：在相乘之前，先对每个因子栅格做一次NoData清洗。用Con工具写：

Con(IsNull(R_g), 0, R_g)

把所有栅格统一清洗之后，最后的侵蚀量结果图就不会出现诡异的空洞。有人担心把NoData补零会拉低侵蚀量，这其实要分情况，C因子和P因子在无数据区补0是偏保守的处理，但在最终报告中一定要明确说明NoData的处理方式，否则数据审查这关过不去。

7.3 批量出图与成果表达

模拟结果最终要变成图件，这步才是很多人的最后一公里。ArcGIS里做批量出图，比较痛苦的是如果每个子流域要出一张专题图，还要在图上插入对应的数据表。检索里也老有人问批量出图怎么插Excel表格。

我的做法是：用ArcMap的数据驱动页面（Data Driven Pages）按子流域索引切分图幅，然后利用动态文本和图表元素，通过XML或地图样式引用外部表格数据。如果嫌配置麻烦，更省事的路线是脚本出图：用arcpy.mapping或者ArcGIS Pro的arcpy.mp遍历要素类，替换标题和比例尺，把统计表写到图片下方。商业项目里我经常把所有图件用Python一次性批量导出到PDF，速度和排版远胜手工切图。

7.4 关于ArcGIS版本迁移的提醒

如果是从ArcMap向ArcGIS Pro迁移，同样的栅格计算器语法基本兼容，但要注意字段计算器和Python环境的关系。ArcMap 10.8自带的是Python 2.7，ArcGIS Pro 3.x用的是Python 3，很多老脚本直接跑会报语法错误。InVEST和ArcGIS之间是通过文件交换的，不存在需要互相调API的场景，所以把版本差异处理在脚本适配这一步就好。

还有一个冷门但真实的问题：ArcGIS Pro 3.7版本装完打开时要登录ArcGIS Online账号，做离线项目的同学很容易被卡住。这个解决路径是配置Portal连接为离线模式，或者在设置里关闭在线登录检测，别被弹窗误导去注册在线账号。

做水土流失模拟这件事，模型公式是全世界通用的，但真正考验人的是把公式落到具体数据上那一步。RUSLE和InVEST SDR的输出都只是初步成果，它们告诉你的是“在现有数据和参数设置下，区域的空间格局大致如此”，而不是精确的数学答案。我个人的体会是，每一张侵蚀分布图背后都藏着一串数据来源的假设，把这些假设如实写清楚，图件的价值才站得住。你把这六个因子和泥沙输移流程完整走一遍之后，再回头看那些高高低低的栅格值，心里对哪些区域该优先治理、哪些地方适合做什么类型的措施，基本就有数了。