CN102609615A - 基于修正通用水土流失方程的水土流失遥感动态监测方法 - Google Patents

Abstract

本发明公开了一种基于修正通用水土流失方程的水土流失遥感动态监测方法，实施步骤如下：1)获取降雨侵蚀力因子、土壤可蚀性因子、地形因子；2)获取待监测区域的遥感图像；根据遥感图像进行土地利用/覆盖类别划分得到多个地类，分别估算各个地类的植被覆盖-管理因子并进行综合得到植被覆盖-管理因子；3)从遥感图像中提取影响土壤水土保持的特征信息，根据所述影响土壤水土保持的特征信息估算水土保持措施因子；4)根据修正通用水土流失方程计算待监测区域的年均土壤流失量估算值并输出。本发明具有估算结果合理、通用性好、适合区域水土流失的长期动态监测、适用范围广的优点。

Description

基于修正通用水土流失方程的水土流失遥感动态监测方法

技术领域

本发明涉及水土流失风险评估领域，具体涉及一种基于修正通用水土流失方程的水土流失遥感动态监测方法。

背景技术

水土流失是指由于水的侵蚀或者风力的作用使得土壤迁出土体，导致地力下降、严重的甚至完全失去地力。古巴比伦国的消亡就是因为水土流失、环境破坏而导致的。美国在1930年代初，也发生了严重的水土流失。然而，政府、社会各阶层和民众都十分重视，采取强有力的水土保持措施，已经很好地控制了水土流失。我国的水土流失现象十分严重。黄河流域是最明显的例子，巨大的泥沙从戈壁、草原、黄土高原随着河水流入黄海。南方长江流域的水土流失同样极为严重。水土流失不仅破坏地力，而且污染环境。目前我国大多水体存在的富营养化在很大程度上就是水土流失的造成的另一恶果。因此，我国的大多数地方都急需进行水土流失治理。然而，要治理水土流失，我们必须明确水土流失的范围和程度。这样，才能有针对性地和开展有效地水土流失治理。

美国在二十世纪三十年代初发生严重的水土流失之后，成立了土壤保持局和水土保持国家实验室，开展了大量、细致的研究工作。经过30多年的努力，在1965年，得出了著名的通用水土流失方程(Universal Soil Loss Equation，USLE)。在此基础上，又经过近三十年的努力，获得了修订通用水土流失方程(Modified Universal Soil Loss Equation，MUSLE)和修正通用水土流失方程(Revised Universal Soil Loss Equation，RUSLE)。目前，国际上大多采用RUSLE来计算年均土壤流失量。其由六个因子决定：A＝R×K×LS×C×P，其中，A为年均土壤流失量估算值(t ha^-1 yr^-1)，R为降雨侵蚀力(MJ mm ha^-1 h^-1 yr^-1)，K为土壤可蚀性因子(t ha h ha^-1 MJ^-1 mm^-1)，LS为坡长与坡度结合量(无量纲)，C为植被覆盖-管理因子(无量纲)，P是水土保持措施因子(无量纲)。

根据以上公式，土壤侵蚀量由气候(降雨)、土壤、地形、植被和土地利用/覆盖等因子共同作用。其中，降雨、土壤和地形因子受自然条件影响，不同时期变化不大。而植被覆盖和土地利用等因子受人类活动影响，是易于变化但又是具有重要作用、影响极大的因子。

目前，在小范围或者局部田块，可以通过田间实测获得水土流失量。大范围的则通过遥感信息来获取，主要涉及两个方面：

1)直接利用遥感影像对水土流失或风险进行分级分类，包括目视解译和计算机自动解译。

目视解译一直是我国水土保持部门进行水土流失遥感调查采取的主要手段。我国于1985年、1999年、2001年先后开展了全国第一、二、三次土壤侵蚀遥感调查均采用该方法。该方法的优点在于可以将人的经验和知识与遥感技术结合起来，充分利用专家的先验知识，避免了单纯的光谱分析可能带来的误差。然而解译没有明确的标准，解译过程中主观性极强，使得其结果难以在空间区域和时间序列上进行对比，而且其需要投入大量的人力、资金和时间，使得成本和时效不能兼顾。

计算机自动解译是为了避免采用目视解译技术所带来的大量人力、资金和时间的投入，从卫星图像中自动提取土壤侵蚀信息的另一种可行技术。但该方法基本上是依靠单纯的光谱信息进行分类，由于土壤侵蚀本身并不是以特定的土地覆盖等地表特征出现，而且指示土壤侵蚀的土壤属性光谱信息往往被植被覆盖、田间管理和耕种方式等等这样的土壤表层信息所掩盖，理论上只利用单独利用遥感信息是难以准确判断土壤侵蚀状况的。影像分类法在土壤侵蚀研究中的应用往往局限在某些特定的半干旱地区，而对于我国南方多云雨、地形复杂区域，限制了该方法的使用。

2)利用遥感数据提取影响水土流失的某些因子，然后利用一定的数据集成模式将因子综合得到水土流失或风险等级或分类图。

土壤侵蚀影响因子主要有降雨侵蚀力因子R，土壤可蚀性因子K，地形因子LS，植被覆盖-管理因子C和水土保持措施因子P。土壤侵蚀是这些因子共同作用的结果。在这些因子中，降水一般由气象数据获得，地形因子与土壤因子可以从遥感数据获取，但由于数据与技术的限制，目前还未见从遥感数据获取降雨侵蚀力因子的报道，地形数据更多是由等高线或DEM数据产生，土壤因子更多从土壤图获取。所以，从遥感影像获得的主要是植被覆盖-管理因子和水土保持措施因子。也正是由于目前数据和技术的这种限制，利用遥感信息进行水土流失监测时，除了遥感数据，往往还要求有大量翔实的其它来源的数据辅助完成。但实际上，这些辅助数据的获取或更新存在一定难度，尤其在欠发达地区，这些数据在大面积上不易获得，要获取与遥感图像时相相匹配的数据更是难上加难。因此，当在一个区域或国家尺度时，数据的可获取性已成为水土流失动态监测的瓶颈之一。

由于数据的难获取性，现有方法要么只能在小范围通过田间实测获得土壤侵蚀量，若是在大范围，获得土壤侵蚀量不客观，和实际情况有较大差别，目前人们无法获得大范围的、有明确科学依据的水土流失信息，更无法进行区域水土流失的长期动态监测。

发明内容

本发明要解决的技术问题是提供一种基于修正通用水土流失方程的水土流失遥感动态监测方法，其具有估算结果合理、通用性好、适合区域水土流失的长期动态监测、适用范围广的特点。

为解决上述技术问题，本发明采用的技术方案为：

一种基于修正通用水土流失方程的水土流失遥感动态监测方法，其实施步骤如下：

1)获取待监测区域指定期间的月降雨数据，将所获取的月降雨数据利用包括海拔高程、经度和纬度在内的三维坐标进行空间回归计算得到降雨侵蚀力因子；根据待监测区域的大比例尺土壤空间数据库中的土壤有机碳含量和土壤颗粒分布计算得到土壤可蚀性因子；根据待监测区域的数字地形模型数据计算得到地形因子；

2)获取待监测区域的遥感图像；根据遥感图像进行土地利用/覆盖类别划分得到多个地类，分别估算各个地类的植被覆盖-管理因子，并将各个地类的植被覆盖-管理因子进行综合得到遥感图像的植被覆盖-管理因子；

3)从所述遥感图像中提取影响土壤水土保持的特征信息，根据所述影响土壤水土保持的特征信息估算水土保持措施因子；

4)根据修正通用水土流失方程，A＝R×K×LS×C×P，计算待监测区域的年均土壤流失量估算值A并输出。其中，R为降雨侵蚀力因子，K为土壤可蚀性因子，LS为地形因子，C为植被覆盖-管理因子，P为水土保持措施因子。

作为本发明上述技术方案的进一步改进：

所述步骤1)中将所述月降雨数据利用包括海拔高程、经度和纬度在内的三维坐标进行空间回归计算得到降雨侵蚀力因子的详细步骤包括：

A1)根据月降雨数据获取月均降雨量P_i和年均降雨量P；

A2)根据

获取校正后的Fournier系数MFI，其中P_i为各月月均降雨量，P为年均降雨量；

A3)对所述Fournier系数MFI进行判断，如果Fournier系数MFI≥55mm则根据R₁＝95.77-6.081MFI+0.4770MFI²计算降雨侵蚀力R₁；如果Fournier系数MFI＜55mm则根据R₁＝0.07397MFI^1.847计算降雨侵蚀力R₁；

A4)将降雨侵蚀力R₁作为因变量，将对应气象站的海拔高度、经度和纬度作为自变量，建立回归方程

其中h为海拔高度，φ为纬度，λ为经度，R即为最终得到的降雨侵蚀力因子。

所述步骤1)中根据所述大比例尺土壤空间数据库中的土壤有机碳含量和土壤颗粒分布计算得到土壤可蚀性因子具体是指根据

$K=\{0.2+0.3\exp [-0.0256\mathrm{SAN}(1-\frac{\mathrm{SIL}}{100})\left]\right\}{\left(\frac{\mathrm{SIL}}{\mathrm{CLA}+\mathrm{SIL}}\right)}^{0.3}(1.0-\frac{0.25C}{C+\exp (3.72-2.95C)})(1.0-\frac{0.7\mathrm{SN}1}{\mathrm{SN}1+\exp (-5.51+22.9\mathrm{SN}1)})$

计算得到土壤可蚀性因子K，其中，SAN、SIL、CLA分别为土壤中砂粒的比例、粉粒的比例、粘粒的比例，C为土壤中的有机碳含量，SN1等于1-SAN/100。

所述步骤1)中根据待监测区域的数字地形模型数据计算得到地形因子具体是指根据

$\mathrm{LS}={\left(\frac{\mathrm{\λ}}{22.13}\right)}^n(65.41{\sin }^2\mathrm{\θ}+4.56\sin \mathrm{\θ}+0.065)$

计算得到地形因子LS，其中，λ为坡长，θ为坡度角，n为基于坡度的常量。

所述步骤2)的详细步骤包括：

B1)将遥感图像进行土地利用/覆盖类别划分得到多种地类；将所述地类进行分组，得到空间异质性小、季节变化小、年际变化小的第一组地类，空间异质性小、季节变化大、年际变化小的第二组地类，空间异质性大、季节变化小、年际变化大的第三组地类；

B2)采用遥感分类法获取第一组地类和第二组地类中每一个地类的植被覆盖-管理因子；

B3)将遥感图像进行线性光谱混合像元分解，得到绿色植被分量图、土壤分量图和阴影/水分量图；通过C＝F_土壤/(1+F_绿色植被+F_阴影)计算得到第三组地类的植被覆盖-管理因子，其中，F_土壤为土壤分量，F_绿色植被为绿色植被分量，F_阴影为阴影/水分量；

B4)将各个地类的植被覆盖-管理因子进行综合得到遥感图像的植被覆盖-管理因子。

所述步骤B4)中将遥感图像进行线性光谱混合像元分解的详细步骤包括：

B41)采用主成分分析、最大噪声比变换或者奇异值分解对遥感图像进行数据降维；

B42)选择端元模型，根据所述端元模型确定端元的数量和光谱；

B43)根据所述端元的数量和光谱对数据降维后的遥感图像进行反演得到绿色植被分量图、土壤分量图和阴影/水分量图。

所述步骤B43)中对数据降维后的遥感图像采用全约束最小二乘法进行反演，所述全约束最小二乘法的模型为

其中，R_ib为第b波段第i像元的光谱反射率，N为端元数目，f_ki为对应于i像元的第k个端元所占分量值，C_ib为第b波段端元K的光谱反射率，ε_ib为该像素第b波段的分解残差；所述模型的约束条件为

且所述模型通过误差ε_i来评价模型的优劣，误差ε_i的表达式为

其中m为光谱波段数目。

所述多种地类包括密林地、建设用地、水体、水田、旱地、疏林地、园地、裸土地八种地类；所述第一组地类包括密林地、建设用地、水体；所述第二组地类包括水田、旱地；所述第三组地类包括疏林地、园地、裸土地。

所述步骤3)中从所述遥感图像中提取影响土壤水土保持的特征包括等高线种植、条带种植和修筑梯田中的至少一种。

本发明具有下述优点：

1、本发明在修正型通用水土流失方程的基础上，提出水土流失监测的背景因子和动态因子概念，将降雨侵蚀力因子、土壤可蚀性因子、地形因子作为背景因子，根据背景因子随时间变化较小或人为因素无法控制的特性，在监测期内可视为背景值以体现水土流失的空间异质性，背景因子只通过获取监测期内的一期数据即可实现估算，大大降低了数据更新的难度；将植被覆盖-管理因子、水土保持措施因子作为动态因子，根据动态因子随时间变化大且受人为影响大、是造成区域水土流失发生变化的最关键因子的特性，直接从遥感图像获取，利用监测期内目标年份的影像获取相应年份的植被覆盖-管理因子、水土保持措施因子，不存在数据获取的限制，在监测期内在保证水土流失空间异质性得以体现的同时，又能够反映水土流失随时间变化的规律，实现区域水土流失的长期动态监测。

2、本发明借助遥感和地理信息系统获取水土流失的背景因子和动态因子，克服了由于数据缺乏无法进行水土流失的长期监测，既能充分体现水土流失的空间特征，又能充分体现水土流失随时间变化的性状，适合区域水土流失的长期动态监测，其监测结果所揭示的水土流失热点区域可为管理者和决策者提供重要信息，具有估算结果合理、通用性好、适合区域水土流失的长期动态监测、适用范围广的优点。

3、本发明通过获得待监测区域指定期间的月降雨数据，将月降雨数据与海拔高程、经度和纬度三维坐标进行空间回归计算得到降雨侵蚀力因子，充分体现了区域内降雨侵蚀力因子的空间异质性，避免了常用方法中区域内采用单一降雨侵蚀力因子值空间异质性无法体现的问题。

4、本发明中对植被覆盖-管理因子的估算综合了土地利用分类赋值和线性像元分解的方法，既克服了由于地物季变对植被覆盖-管理因子带来的影响，又能充分体现地物植被覆盖-管理因子的空间异质性，使得本发明的水土流失监测技术既不受遥感获取时相的影像，又能客观估算区域的水土流失量，具有估算结果准确、通用性好的优点。适用于多云雨地区，特别适合较难获得光学遥感图像信息的我国南方多云雨地区，可以广泛地应用于农业、环保、水利、国土等部门。

附图说明

图1为本发明实施例的基本流程示意图。

图2为本发明实施例的详细流程示意图。

图3为本发明实施例得到的降雨侵蚀力因子图。

图4为本发明实施例得到的土壤可蚀性因子图。

图5为本发明实施例得到的地形因子图。

图6为本发明实施例输入的遥感图像。

图7为本发明实施例得到的土地利用/覆盖分类图。

图8为本发明实施例得到的植被覆盖-管理因子图。

图9为本发明实施例得到的水土保持措施因子图。

图10为本发明实施例得到的年均土壤流失量估算值图。

具体实施方式

如图1和图2所示，本实施例基于修正通用水土流失方程的水土流失遥感动态监测方法的实施步骤如下：

1)获取待监测区域指定期间的月降雨数据，将月降雨数据与海拔高程、经度和纬度三维坐标进行空间回归计算得到降雨侵蚀力因子(R)；

根据待监测区域的大比例尺土壤空间数据库中的土壤有机碳含量和土壤颗粒分布计算得到土壤可蚀性因子(K)；

根据待监测区域的数字地形模型数据计算得到地形因子(LS)；

2)获取待监测区域的遥感图像；根据遥感图像进行土地利用/覆盖类别划分得到多个地类，分别估算各个地类的植被覆盖-管理因子(C)，并将各个地类的植被覆盖-管理因子进行综合得到遥感图像的植被覆盖-管理因子(C)；

3)从遥感图像中提取影响土壤水土保持的特征，根据影响土壤水土保持的特征估算水土保持措施因子(P)；

4)根据修正通用水土流失方程，A＝R×K×LS×C×P(RUSLE模型)，计算出待监测区域的年均土壤流失量估算值A。

本实施例中，步骤1)中将月降雨数据与海拔高程、经度和纬度三维坐标进行空间回归计算得到降雨侵蚀力因子的详细步骤包括：

A1)根据月降雨数据获取月均降雨量P_i和年均降雨量P。

A2)根据

获取校正后的Fournier系数MFI，其中P_i(单位：mm)为各月月均降雨量，P(单位：mm)为年均降雨量；本实施例中，待监测区域的采样气象站数量为11个。

A3)对Fournier系数MFI进行判断，如果Fournier系数MFI≥55mm则根据R₁＝95.77-6.081MFI+0.4770MFI²计算降雨侵蚀力R₁；如果Fournier系数MFI＜55mm则根据R₁＝0.07397MFI^1.847计算降雨侵蚀力R₁；

A4)将降雨侵蚀力R₁作为因变量，将对应气象站的海拔高度、经度和纬度作为自变量，建立回归方程

其中h(单位：m)为海拔高度，φ(单位：°)为纬度，λ(单位：°)为经度，R即为最终得到的降雨侵蚀力因子。本实施例中，自变量是采用浙江省安吉县30m×30m高程图和经纬度，经过建立回归方程最终生成的生成30m×30m的降雨侵蚀力因子(R)栅格图如图3所示。

本实施例中，步骤1)中根据大比例尺土壤空间数据库中的土壤有机碳含量和土壤颗粒分布计算得到土壤可蚀性因子具体是指根据

$K=\{0.2+0.3\exp [-0.0256\mathrm{SAN}(1-\frac{\mathrm{SIL}}{100})\left]\right\}{\left(\frac{\mathrm{SIL}}{\mathrm{CLA}+\mathrm{SIL}}\right)}^{0.3}(1.0-\frac{0.25C}{C+\exp (3.72-2.95C)})(1.0-\frac{0.7\mathrm{SN}1}{\mathrm{SN}1+\exp (-5.51+22.9\mathrm{SN}1)})$

计算得到土壤可蚀性因子K(单位：t ha h ha^-1 MJ^-1 mm^-1)，其中，SAN、SIL、CLA分别为土壤中砂粒(粒径0.05-2.00mm)的百分比例、粉粒(粒径0.002-0.05mm)的百分比例、粘粒(粒径＜0.002mm)的百分比例，C为土壤中的有机碳含量，SN1等于1-SAN/100。本实施例中，最终得到浙江省安吉县的土壤可蚀性因子(K)图层如图4所示。

本实施例中，步骤1)中根据待监测区域的数字地形模型数据计算得到地形因子具体是指根据利用Wischmeier和Smith设计的经验模型：

$\mathrm{LS}={\left(\frac{\mathrm{\λ}}{22.13}\right)}^n(65.41{\sin }^2\mathrm{\θ}+4.56\sin \mathrm{\θ}+0.065)$

计算得到地形因子LS，其中λ为坡长，θ为坡度角，n为基于坡度的常量。本实施例中，如坡度大于2.86°则n取值0.5；如坡度位于1.72°～2.85°之间，则n取值0.4；如坡度位于0.57°～1.72°之间，则n取0.3；如坡度小于0.7°，则n取0.2。本实施例的数字地形模型数据源于安吉县30m的数字高程图，全球大部分地方的30m分辨率数字高程图都可以从网站上获得，最终得到浙江省安吉县的地形因子(LS)的图形如图5所示。

本实施例中，步骤2)和步骤3)均分别基于图6所示的遥感图像进行估算植被覆盖-管理因子(C)和水土保持措施因子(P)。

步骤2)的详细步骤包括：

B1)将遥感图像进行土地利用/覆盖类别划分得到多种地类；将地类进行分类，得到空间异质性小、季节变化小、年际变化小的第一组地类，空间异质性小、季节变化大、年际变化小的第二组地类，空间异质性大、季节变化小、年际变化大的第三组地类；

B2)采用遥感分类法获取第一组地类和第二组地类中每一个地类的植被覆盖-管理因子。

B3)将遥感图像进行线性光谱混合像元分解，得到绿色植被分量图、土壤分量图和阴影/水分量图；通过C＝F_土壤/(1+F_绿色植被+F_阴影)计算得到第三组地类的植被覆盖-管理因子，其中，F_土壤为土壤分量，F_绿色植被为绿色植被分量，F_阴影为阴影/水分量。

B4)将各个地类的植被覆盖-管理因子进行综合得到遥感图像的植被覆盖-管理因子(C)，最终得到植被覆盖-管理因子(C)图层如图8所示。

本实施例中，步骤B4)中将遥感图像进行线性光谱混合像元分解的详细步骤包括：

B41)采用主成分分析、最大噪声比变换或者奇异值分解对遥感图像进行数据降维；本实施例中，步骤A)中采用主成分分析对遥感图像进行数据降维。确定端元数量以及端元的光谱是成功的混合像元分解的关键。对于城市地区，最常用的端元选取方式是由Ridd(1995)提出的植被-不透水层-土壤端元模型(Vegetation-Impervious surface-Soil，V-I-S)。而在非城市地区，一般采用植被-土壤-阴影(或干植被)端元模型。本实施例步骤B)中端元模型选择的是植被-土壤-阴影/水端元模型，利用像元纯度指数(Pixel Purity Index，PPI)进行交互式提取端元的数量和光谱。

B42)选择端元模型，根据端元模型确定端元的数量和光谱；

B43)根据端元的数量和光谱对数据降维后的遥感图像进行反演得到绿色植被分量图、土壤分量图和阴影/水分量图。

步骤B43)中对数据降维后的遥感图像采用全约束最小二乘法进行反演，全约束最小二乘法的模型为

其中，R_ib为第b波段第i像元的光谱反射率，N为端元数目，f_ki为对应于i像元的第k个端元所占分量值，C_ib为第b波段端元K的光谱反射率，ε_ib为该像素第b波段的分解残差；模型的约束条件为

且0≤f_ki≤1，且模型通过误差ε_i来评价模型的优劣，误差ε_i的表达式为

其中m为光谱波段数目。

本实施例估算植被覆盖-管理因子(C)的具体实施步骤如下：

I)、将如图6所示的浙江省安吉县2003年和2008年获得的遥感图像进行土地利用/覆盖类别划分得到包括密林地、建设用地、水体、水田、旱地、疏林地、园地、裸土地八种地类，得到的土地利用/覆盖分类图如图7所示。

II)、将前述的八种地类进行分类，得到：

空间异质性小、季节变化小、年际变化小的第一组地类，该组地类包括密林地、建设用地、水体；

空间异质性小、季节变化大、年际变化小的第二组地类，本组地类包括水田、旱地；

空间异质性大、季节变化小、年际变化大的第三组地类，这组地类包括疏林地、园地、裸土地。

III)、采用遥感分类法获取第一组地类和第二组地类中每一个地类的植被覆盖-管理因子，即根据相关文献报道并结合当地条件对第一组地类和第二组地类中的每一个地类指定一个植被覆盖-管理因子值，以消除遥感图像的季相变化带来的不确定性影响；本实施例针对浙江省安吉县，根据相关文献报道并结合当地条件，第一组地类和第二组地类中的每一个地类的植被覆盖-管理因子值指定为：密林地，0.0；水体，0.0；建设用地，0.2；水田，0.05；旱地，0.44。

IV)、将如图3所示的遥感图像进行线性光谱混合像元分解，得到绿色植被分量图、土壤分量图和阴影/水分量图；通过C＝F土壤/(1+F绿色植被+F阴影)计算得到第三组地类中每一个地类的植被覆盖-管理因子，其中，F土壤为土壤分量，F绿色植被为绿色植被分量，F阴影为阴影/水分量，F绿色植被、F土壤和F阴影的值均为0至1之间的值，且F绿色植被、F土壤和F阴影三者之和为1。

V)、将各个地类的植被覆盖-管理因子进行综合得到遥感图像的植被覆盖-管理因子(如图8所示)。

本实施例中，步骤3)中从图6所示的遥感图像中提取影响土壤水土保持的特征包括等高线种植、条带种植和修筑梯田中的至少一种。根据影响土壤水土保持的特征估算水土保持措施因子(P)的值范围从0到1，值为0时表示有良好的人为防侵蚀措施，为1时表示无措施。在对研究区实际考察的基础上，结合经验值，确定水田的P值为0.1，其它地类基本上没有采用水保措施，P值为1，最终得到水土保持措施因子图层如图9所示。

本实施例中，步骤4)根据修正通用水土流失方程A＝R×K×LS×C×P(RUSLE模型)计算待监测区域的年均土壤流失量估算值A的图像如图10所示。

以上所述仅为本发明的优选实施方式，本发明的保护范围并不仅限于上述实施方式，凡是属于本发明原理的技术方案均属于本发明的保护范围。对于本领域的技术人员而言，在不脱离本发明的原理的前提下进行的若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (9)

1.一种基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于其实施步骤如下：

1)获取待监测区域指定期间的月降雨数据，将所述月降雨数据与海拔高程、经度和纬度三维坐标进行空间回归计算得到降雨侵蚀力因子；根据待监测区域的大比例尺土壤空间数据库中的土壤有机碳含量和土壤颗粒分布计算得到土壤可蚀性因子；根据待监测区域的数字地形模型数据计算得到地形因子；

2)获取待监测区域的遥感图像；根据遥感图像进行土地利用/覆盖类别划分得到多个地类，分别估算各个地类的植被覆盖-管理因子，并将各个地类的植被覆盖-管理因子进行综合得到植被覆盖-管理因子；

3)从所述遥感图像中提取影响土壤水土保持的特征信息，根据所述影响土壤水土保持的特征信息估算水土保持措施因子；

4)根据修正通用水土流失方程A＝R×K×LS×C×P计算待监测区域的年均土壤流失量估算值A并输出，其中，R为降雨侵蚀力因子，K为土壤可蚀性因子，LS为地形因子，C为植被覆盖-管理因子，P为水土保持措施因子。

2.根据权利要求1所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于，所述步骤1)中将所述月降雨数据与海拔高程、经度和纬度三维坐标进行空间回归计算得到降雨侵蚀力因子的详细步骤包括：

A1)根据月降雨数据获取月均降雨量P_i和年均降雨量P；

A2)根据获取校正后的Fournier系数MFI，其中P_i为各月月均降雨量，P为年均降雨量；

A3)对所述Fournier系数MFI进行判断，如果Fournier系数MFI≥55mm则根据R₁＝95.77-6.081MFI+0.4770MFI²计算降雨侵蚀力R₁；如果Fournier系数MFI＜55mm则根据R₁＝0.07397MFI^1.847计算降雨侵蚀力R₁；

A4)将降雨侵蚀力R₁作为因变量，将对应气象站的海拔高度、经度和纬度作为自变量，建立回归方程

其中h为海拔高度，φ为纬度，λ为经度，R即为最终得到的降雨侵蚀力因子。

3.根据权利要求1所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于：所述步骤1)中根据所述大比例尺土壤空间数据库中的土壤有机碳含量和土壤颗粒分布计算得到土壤可蚀性因子具体是指根据

$K=\{0.2+0.3\exp [-0.0256\mathrm{SAN}(1-\frac{\mathrm{SIL}}{100})\left]\right\}{\left(\frac{\mathrm{SIL}}{\mathrm{CLA}+\mathrm{SIL}}\right)}^{0.3}(1.0-\frac{0.25C}{C+\exp (3.72-2.95C)})(1.0-\frac{0.7\mathrm{SN}1}{\mathrm{SN}1+\exp (-5.51+22.9\mathrm{SN}1)})$

计算得到土壤可蚀性因子K，其中，SAN、SIL、CLA分别为土壤中砂粒的比例、粉粒的比例、粘粒的比例，C为土壤中的有机碳含量，SN1等于1-SAN/100。

4.根据权利要求1所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于，所述步骤1)中根据监测待监测区域的数字地形模型数据计算得到地形因子具体是指根据

$\mathrm{LS}={\left(\frac{\mathrm{\λ}}{22.13}\right)}^n(65.41{\sin }^2\mathrm{\θ}+4.56\sin \mathrm{\θ}+0.065)$

计算得到地形因子LS，其中，λ为坡长，θ为坡度角，n为基于坡度的常量。

5.根据权利要求1所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于所述步骤2)的详细步骤包括：

B1)将遥感图像进行土地利用/覆盖类别划分得到多种地类；将所述地类进行分类，得到空间异质性小、季节变化小、年际变化小的第一组地类，空间异质性小、季节变化大、年际变化小的第二组地类，空间异质性大、季节变化小、年际变化大的第三组地类；

B2)采用遥感分类法获取第一组地类和第二组地类中每一个地类的植被覆盖-管理因子；

B3)将遥感图像进行线性光谱混合像元分解，得到绿色植被分量图、土壤分量图和阴影/水分量图；通过C＝F_土壤/(1+F_绿色植被+F_阴影)计算得到第三组地类的植被覆盖-管理因子，其中，F_土壤为土壤分量，F_绿色植被为绿色植被分量，F_阴影为阴影/水分量；

B4)将各个地类的植被覆盖-管理因子进行综合得到遥感图像的植被覆盖-管理因子。

6.根据权利要求5所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于，所述步骤B4)中将遥感图像进行线性光谱混合像元分解的详细步骤包括：

B41)采用主成分分析、最大噪声比变换或者奇异值分解对遥感图像进行数据降维；

B42)选择端元模型，根据所述端元模型确定端元的数量和光谱；

B43)根据所述端元的数量和光谱对数据降维后的遥感图像进行反演得到绿色植被分量图、土壤分量图和阴影/水分量图。

7.根据权利要求6所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于，所述步骤B43)中对数据降维后的遥感图像采用全约束最小二乘法进行反演，所述全约束最小二乘法的模型为

其中，R_ib为第b波段第i像元的光谱反射率，N为端元数目，f_ki为对应于i像元的第k个端元所占分量值，C_ib为第b波段端元K的光谱反射率，ε_ib为该像素第b波段的分解残差；所述模型的约束条件为

且0≤f_ki≤1，且所述模型通过误差ε_i来评价模型的优劣，误差ε_i的表达式为

其中m为光谱波段数目。

8.根据权利要求7所述的土壤侵蚀植被覆盖-管理因子的遥感估算方法，其特征在于：所述多种地类包括密林地、建设用地、水体、水田、旱地、疏林地、园地、裸土地八种地类；所述第一组地类包括密林地、建设用地、水体；所述第二组地类包括水田、旱地；所述第三组地类包括疏林地、园地、裸土地。

9.根据权利要求1所述的基于修正通用水土流失方程的水土流失遥感动态监测方法，其特征在于：所述步骤3)中从所述遥感图像中提取影响土壤水土保持的特征包括等高线种植、条带种植和修筑梯田中的至少一种。

Images