基于比值法的农田土壤水分反演方法及系统

摘要

本发明公开一种基于比值法的农田土壤水分反演方法及系统，首先获取SAR影像和光学影像，并进行预处理；然后提取后向散射系数以及植被描述参数；构建土壤水分估算数据集，并确定训练集和测试集；在训练集上构建比值方程；基于最小代价函数策略估算土壤水分和估算精度；有效粗糙度参数条件下构建的比值方程为最终的土壤水分估算模型；最后进行精度验证。本发明提供的方法利用比值方法分离植被散射贡献，简便且计算高效，减少了对于实测粗糙度参数的依赖，通过SAR提取的植被描述参数，以及借助光学数据提取植被参数，用于表征比值方程中植被的散射特征，发挥了光学数据能较好的表征植被生长状况的优势，在一定程度上提高了土壤水分的反演精度。

说明书

技术领域

本发明涉及遥感图像数据处理技术领域，特别是一种基于比值法的农田土壤水分反演方法。

背景技术

土壤水分是生态水循环中的关键一环，控制着地表径流、地表水渗透以及地表水蒸散发等过程。在农业生产中，土壤水分含量的降低会导致土壤干裂，引发农业干旱，进而造成农作物大面积减产、威胁到社会的粮食安全。因此实时准确获取大面积农田的土壤水分信息，对于农作物估产、旱涝预测以及农业生产管理有着重要的指导意义。

传统土壤水分依赖于野外实测，费时费力且只能获取少量有限点的信息，遥感技术的出现为大范围土壤水分信息的实时准确监测提供了一种有效途径。微波遥感具有较强的穿透能力，具备全天时监测的能力，不受云雾等天气条件影响，并且相较于光学遥感，地表微波散射辐射对于土壤水分含量的变化加敏感，因此微波遥感是土壤水分遥感监测中的主要手段。目前基于主动微波遥感尤其是SAR的土壤水分反演已经提出了不少模型，并且取得了较好的反演结果。但是大多数模型的建立均是针对裸土区域，无法直接应用于植被覆盖区域，这是因为在植被覆盖区域，SAR传感器发出的微波信号在植被与地表间会产生复杂的相互作用，使得最终接收到的后向散射中包含了多种形式的散射分量：来自植被、来自地表或来自于植被与地表相互作用。估算植被覆盖区域如农田种植区地表的土壤水分，关键是需要将传感器接收的信号中来自植被的散射贡献分离出去，获得来自地表的后向散射。

基于SAR的土壤水分反演模型受多种参数影响，包括土壤水分、植被冠层描述参数、传感器配置参数、地表粗糙度参数等，但是在植被茂密的区域，部分参数是难以测量和确定的，如粗糙度参数和植被冠层描述参数。对于缺失部分实测参数的区域，土壤水分反演模型的使用将受到限制。

发明内容

有鉴于此，本发明的目的在于提供一种基于比值法的农田土壤水分反演方法，该方法利用解决了基于SAR反演土壤水分可能存在的上述问题，利用比值方法从传感器观测到的总的SAR信号中分离植被散射分量，从而消除植被的散射贡献，进而估算土壤含水量。

为达到上述目的，本发明提供如下技术方案：

本发明提供的基于比值法的农田土壤水分反演方法，包括以下步骤：

获取SAR影像和光学影像，并进行预处理，所述影像数据为覆盖了种植植物和土壤图像；

从预处理后的SAR影像和光学影像中提取后向散射系数以及植被描述参数；

构建土壤水分估算数据集，并确定训练集和测试集；

在训练集上，构建不同粗糙度条件的地表散射模型，通过地表散射模型模拟裸土后向散射系数；结合观测的总后向散射系数、提取的植被描述参数、模拟的裸土后向散射系数构建比值方程；计算得到训练集上估算的裸土后向散射系数，基于最小代价函数策略，从当前粗糙度条件下地表散射模型构建的查找表中反演土壤水分；计算估算土壤水分与实测土壤水分的均方根误差作为训练集土壤水分的估算精度；

遍历给定粗糙度范围内的粗糙度参数，重复上述步骤，取估算精度最高时的粗糙度参数作为研究区内的有效粗糙度参数；有效粗糙度参数条件下构建的比值方程为最终的土壤水分估算模型；

基于有效粗糙度参数构建地表散射模型查找表，然后使用最小代价函数反演策略，从查找表中反演测试集的土壤水分，然后与测试集实测土壤含水量进行比较，进行精度验证。

进一步，所述比值方程按照以下公式建立：

其中，

F＝aV+bV

其中，V代表植被描述参数，a、b、c为比值函数的待定系数；

使用训练集的数据结合最小二乘法，求出比值函数的待定系数。

进一步，所述后向散射系数的提取按照以下步骤进行：

根据采样点的经纬度信息，从预处理后的SAR影像中提取出采样点的四种极化方式的后向散射系数，所述四种极化方式包括HH/VV/HV/VH；

对提取出的后向散射系数进行入射角的归一化处理，计算方法如下：

其中，θ为局部入射角，θ

从多光谱数据中提取出红光波段和近红外波段的反射率信息，然后计算出采样点的两种植被指数，按照以下公式计算RVI和NDVI，将RVI和NDVI作为比值方程中的植被描述参数：

其中，

其中，NIR和R分别表示近红外波段和红光波段的地表反射率。

进一步，所述训练集和测试集是按照以下方式进行的：

将提取出的总后向散射系数、植被参数与实测的土壤水分一同构建土壤水分估算数据集，并按照预设比例分成训练集与测试集。

进一步，所述裸土后向散射系数是使用CIEM地表散射模型来计算，所述CIEM模型中相关长度与均方根高度经验关系按照以下公式计算：

L(s,θ,pp)＝a(sinθ)

其中，L为地表相关长度，s为均方根高度，pp表示HH或VV极化方式，θ表示局部入射角；a、b、c、d为待定系数。

进一步，所述最小代价函数反演策略按照以下公式计算：

其中，

进一步，所述测试集的土壤水分的计算按照以下步骤进行：

使用构建好的最终给的比值方程计算得到估算的裸土后向散射系数，然后基于有效粗糙度参数下CIEM模型模拟的裸土后向散射系数，结合最小代价函数策略，反演得到测试集上的土壤水分，最后结合测试集上的实测土壤水分完成精度评估。

进一步，还包括土壤含水量区域性估算与制图，具体步骤如下：

对于每一个采样日期，使用有效粗糙度参数构建CIEM模型的后向散射系数查找表，然后遍历遥感数据中的种植区的所有像元，通过土壤含水量估算方法计算每个像元的土壤含水量信息，最终完成土壤含水量估算并进行土壤含水量区域制图。

本发明提供的基于比值法的农田土壤水分反演系统，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现上述方法的步骤。

本发明的有益效果在于：

本发明提供的基于SAR和比值方法的土壤含水量方法，该方法利用比值方法分离植被散射贡献，简便且计算高效。使用有效粗糙度参数来参数化地表散射模型，减少了对于实测粗糙度参数的依赖。同时除了使用SAR提取的植被描述参数外，还借助光学数据提取植被参数，用于表征比值方程中植被的散射特征，发挥了光学数据能较好的表征植被生长状况的优势，在一定程度上提高了土壤水分的反演精度。

本发明的其他优点、目标和特征在某种程度上将在随后的说明书中进行阐述，并且在某种程度上，基于对下文的考察研究对本领域技术人员而言将是显而易见的，或者可以从本发明的实践中得到教导。本发明的目标和其他优点可以通过下面的说明书来实现和获得。

附图说明

为了使本发明的目的、技术方案和有益效果更加清楚，本发明提供如下附图进行说明：

图1为本发明的具体实施方案中小麦研究区示意图，所选研究区为加拿大安大略省东南部某小麦种植区。

图2为本发明的整体方案流程图。

图3为本发明具体实施方案中使用两种不同植被参数来参数化比值方程在测试集上的估算精度。

图4为本发明具体案例实施方案中获得的不同采样日期小麦种植区的土壤含水量区域性制图。

具体实施方式

下面结合附图和具体实施例对本发明作进一步说明，以使本领域的技术人员可以更好的理解本发明并能予以实施，但所举实施例不作为对本发明的限定。

实施例1

如图1所示，本实施例提供的基于比值法的农田土壤水分反演方法，包括以下步骤：

获取SAR影像和光学影像，并进行预处理，所述影像数据为覆盖了种植植物和土壤图像。

从预处理后的SAR影像和光学影像中提取后向散射系数(VV、HH、HV和VH)以及植被描述参数(RVI以及NDVI)，其中RVI从SAR影像中获取，NDVI从光学影像中获取。

构建土壤水分估算数据集，并确定训练集和测试集。

在训练集上，构建不同粗糙度条件的地表散射模型，通过地表散射模型模拟裸土后向散射系数。结合观测的总后向散射系数、提取的植被描述参数、模拟的裸土后向散射系数构建比值方程，使用最小二乘法求出比值方程的待定系数。求出比值方程的待定系数后，计算得到训练集上估算的裸土后向散射系数，基于最小代价函数策略，从当前粗糙度条件下地表散射模型构建的查找表中反演土壤水分。计算估算土壤水分与实测土壤水分的均方根误差(RMSE)作为训练集土壤水分的估算精度。

遍历给定粗糙度范围内的粗糙度参数，重复上述步骤，取估算精度最高时的粗糙度参数作为研究区内的有效粗糙度参数。有效粗糙度参数条件下构建的比值方程为最终的土壤水分估算模型。

基于有效粗糙度参数构建地表散射模型查找表，然后使用最小代价函数反演策略，从查找表中反演测试集的土壤水分，然后与测试集实测土壤含水量进行比较，进行精度验证。

实施例2

下面结合具体实施例和说明书附图对本发明提供的基于SAR和比值方法的土壤含水量估算方法作进一步说明：

数据获取：本方法所采用的遥感数据为RADARSAT-2全极化SAR数据和Sentinel-2多光谱数据，实测数据包括小麦种植区部分采样点的土壤含水量信息。研究区位置如图1所示。SAR数据与实测数据获取日期一致，涵盖了小麦生长期内的八个日期，具体日期如图4所注。光学数据获取日期与采样日期尽可能接近，并且保证影像内研究区范围内无云覆盖，最终有6个采样日期有满足要求的光学影像(不包含5月9日与5月20日)。

如图2所示，该估算土壤含水量的方法，包括以下步骤：

步骤一：遥感影像的预处理，主要包括辐射定标、斑点滤波、地形校正与投影转换。

步骤二：根据采样点的经纬度信息，从预处理后的SAR影像中提取出采样点的四种极化方式(HH/VV/HV/VH)的后向散射系数，

其中，HH表示水平发射水平接收电磁波；VV表示垂直发射垂直接收电磁波；HV表示水平发射垂直接收电磁波；VH表示垂直发射水平接收电磁波。

SAR的局部入射角通常会对地表后向散射系数产生影响，这在土壤水分估算中是不可忽略的。为了限制不同局部入射角的影响，对提取出的后向散射系数进行入射角的归一化处理，计算方法如下：

其中，θ为局部入射角，θ

从Sentinel-2影像中提取出红光波段(Red)和近红外波段(NIR)的反射率信息，然后计算出采样点的两种植被指数，RVI和NDVI作为比值方程中的植被描述参数。RVI和NDVI的计算方式如下：

其中，

其中，NIR和R分别表示近红外波段和红光波段的地表反射率。

将提取出的总后向散射系数、植被参数与实测的土壤水分一同构建土壤水分估算数据集，并按照7:3的比例分成训练集与测试集。

步骤三：使用CIEM地表散射模型来模拟裸露土壤的后向散射系数。IEM模型构建了同极化后向散射系数与地表粗糙度参数(均方根高度、相关长度)、入射角、入射波长、土壤介电常数间的函数关系，可以表示为公式4：

其中：

其中，s表示裸土表面均方根高度；k表示雷达波数，其计算方式为k＝2π/λ，λ为入射波长；θ表示入射角；

Baghdadi等根据野外调查的大量粗糙度数据与C波段SAR数据，结合IEM模型的模拟值，构建了均方根高度与相关长度间的函数表达式，提出了CIEM模型，减少了对于相关长度或均方根高度的依赖。针对C波段SAR建立的CIEM模型中关于相关长度与均方根高度经验关系可以表示为：

L(s,θ,pp)＝a(sinθ)

其中，L为地表相关长度，θ为入射角，s为均方根高度，pp表示HH或VV极化方式，a、b、c、d为待定系数，根据极化方式的不同，取值如下：

a

b

其中，a

给定一个粗糙度参数范围，对于本实例，均方根高度s的范围设置为1mm-30mm，步长为0.1mm。然后使用给定范围内的粗糙度参数构建地表散射模型，模拟小麦种植区裸土的后向散射系数，然后结合观测的总后向散射系数、从遥感影像中提取的植被描述参数，构建比值方程，比值方程的表达方式如下：

其中，pp表示HH或VV极化方式，

F＝aV+bV

其中，V代表植被描述参数，本实例中的V使用RVI和NDVI进行参数化，a、b、c等参数为比值函数的待定系数。

使用训练集的数据结合最小二乘法，求出比值函数的待定系数。然后得到训练集上估算的裸土后向散射系数。进而结合当前粗糙度条件下CIEM模型构建的查找表，基于最小代价函数策略(公式12)，反演得到训练集上估算的土壤水分。进一步与训练集上实测的土壤水分进行比较，求出估算值与实测值的RMSE作为训练集上的估算精度。

其中，

步骤四：遍历给定粗糙度参数区间内的所有均方根高度(s)，重复步骤三，直到s达到给定范围的最大值。比较不同s下训练集上土壤水分的估算精度，取估算精度最高时的s作为该时相下研究区的有效粗糙度参数，有效粗糙度参数下构建的比值方程作为最终的土壤水分反演模型。

步骤五：在测试集上，使用构建好的最终给的比值方程计算得到估算的裸土后向散射系数，然后基于有效粗糙度参数下CIEM模型模拟的裸土后向散射系数，结合最小代价函数策略，反演得到测试集上的土壤水分，最后结合测试集上的实测土壤水分完成精度评估(如图3所示)。

步骤六：小麦农田土壤含水量区域性估算与制图。对于每一个采样日期，使用有效粗糙度参数构建CIEM模型的后向散射系数查找表，然后遍历RADARSAT-2影像中的小麦种植区的所有像元，使用上述步骤中的土壤含水量估算方法，估算每个像元的土壤含水量信息，最终完成小麦农田区域的土壤含水量估算并进行土壤含水量区域制图。基于RVI的研究区土壤水分反演制图如图4所示。

以上所述实施例仅是为充分说明本发明而所举的较佳的实施例，本发明的保护范围不限于此。本技术领域的技术人员在本发明基础上所作的等同替代或变换，均在本发明的保护范围之内。本发明的保护范围以权利要求书为准。