CN107644284A - 一种农田蒸散估算方法及系统 - Google Patents

Abstract

本发明提供一种农田蒸散估算方法及系统，所述方法包括：S1、获取农田的遥感影像数据和气象数据；S2、基于所述农田的遥感影像数据和气象数据，计算潜热通量、净辐射通量和土壤热通量；S3、基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。本发明提供的农田蒸散估算方法及系统，通过充分利用气象数据与光学遥感和热红外遥感相结合，物理基础明确，空间分辨率高，估算结果具有更强的可靠性。

Description

一种农田蒸散估算方法及系统

技术领域

本发明涉及水文气象技术领域，更具体地，涉及一种农田蒸散估算方法及系统。

背景技术

蒸散包括地表水分蒸发、植被截留蒸发和植被蒸腾。蒸散联系着水圈、大气圈和生物圈，是地表水量平衡和能量平衡的重要组成部分，在农业水资源管理、作物估产和干旱监测等领域具有重要作用，而农田蒸散是发生在农田地表的蒸散过程，因此，实时地监测农田蒸散，特别是开展基于遥感信息的农田蒸散算法研究，对于合理管理水资源和防旱抗旱有重要的科学意义和现实指导价值。

早期的蒸散估算主要是通过波文比能量平衡和涡度相关的方法获得准确的点测量数据。而由于地表能量传输的复杂性以及空间分布的异质性，点观测数据代表性很差。而大范围密集观测成本很高，基于单点的蒸散结果很难满足各种应用需求。遥感技术能够提供地表的多源多维多时相信息，具有大面积、宏观、实时和动态等优势，为农田蒸散估算开辟了新的途径。遥感反演农田蒸散的方法很多，主要包括可见光-近红外遥感蒸散算法和热红外遥感蒸散算法，可见光-近红外的遥感蒸散算法主要通过可见光-近红外数据反演的植被指数叶面积指数等相关关系对模型的部分参数进行参数化并结合气象数据对蒸散进行估算。

但现有技术的农田蒸散算法受气象数据和植被指数的影响较大，因此在低植被区域或非生长季，蒸散估算结果会比较平滑，降低了蒸散估算结果的准确性和可靠性。

发明内容

本发明提供一种克服上述问题或者至少部分地解决上述问题的一种农田蒸散估算方法及系统。

根据本发明的一个方面，提供一种农田蒸散估算方法，其特征在于，包括：

S1、获取农田的遥感影像数据和气象数据；

S2、基于所述农田的遥感影像数据和气象数据，计算潜热通量、净辐射通量和土壤热通量；

S3、基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。

其中，步骤S2包括：

S21、基于所述农田的遥感影像数据，计算农田的植被指数；

S22、基于所述农田的植被指数和所述气象数据，计算地表温度；

S23、基于所述农田的植被指数、所述地表温度以及所述农田的遥感影像数据和气象数据，计算所述潜热通量、所述净辐射通量和所述土壤热通量。

其中，步骤S23包括：

S231、基于所述农田的遥感影像数据和气象数据，计算净辐射通量；

S232、基于所述净辐射通量、所述地表温度，计算所述土壤热通量；

S233、基于所述农田的遥感影像数据和气象数据、所述净辐射通量以及所述土壤热通量，计算所述潜热通量。

其中，步骤S233包括：

基于所述农田的遥感影像数据和气象数据，计算显热通量；

基于能量余项法、所述显热通量、所述净辐射通量以及所述土壤热通量，计算所述潜热通量。

其中，所述基于所述农田的遥感影像数据和气象数据，计算显热通量，包括：

基于所述农田的植被指数和地表粗糙度的估算模型，估算地表粗糙度；

基于所述地表粗糙度和所述农田的遥感影像数据和气象数据，计算所述显热通量。

其中，步骤S21具体包括：

基于所述农田的遥感影像数据，获取近红外波段的地表反射率ρ_nir、可见光红波段的地表反射率ρ_red和蓝光波段的地表反射率ρ_blue；

基于公式计算农田的植被指数；

其中，EVI为所述农田的植被指数。

其中，步骤S231具体包括：

基于所述农田的植被指数和所述气象数据，获取地表反照率α、下行的太阳短波辐射R_s↓、大气长波辐射R_L↓、地表长波辐射R_L↑以及地表比辐射率ε；

基于公式R_n＝(1-α)R_s↓+R_L↓-R_L↑-(1-ε)R_L↓，计算净辐射通量R_n。

其中，步骤S232具体包括：

基于所述净辐射通量、所述农田的植被指数以及所述地表温度，获取地表反照率α、地表温度T_s、农田的植被指数EVI；

基于公式计算所述土壤热通量G。

其中，所述基于所述农田的遥感影像数据和气象数据，计算显热通量，具体包括：

基于所述农田的遥感影像数据和气象数据，获取空气密度ρ_air、空气定压比热C_p、空气动力学阻抗r_ah，温差dT；

基于公式计算显热通量H。

根据本发明提供的另一方面，本发明提供一种农田蒸散估算系统，包括：

数据模块，用于获取农田的遥感影像数据和气象数据；

土壤热通量模块，用于基于所述农田的遥感影像数据和气象数据，计算所述土壤热通量；

净辐射通量模块，用于基于所述农田的遥感影像数据和气象数据，计算所述净辐射通量；

潜热通量模块，用于基于所述农田的遥感影像数据和气象数据，计算所述潜热通量；

农田蒸散估算模块，用于基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。

本发明提供的农田蒸散估算方法及系统，通过充分利用气象数据与光学遥感和热红外遥感相结合，物理基础明确，空间分辨率高，估算结果具有更强的可靠性。

附图说明

图1是本发明实施例提供一种农田蒸散估算方法流程图；

图2是本发明实施例提供的各个观测点的测量的蒸散值分布图；

图3是本发明实施例提供的一种农田蒸散估算系统结构图。

具体实施方式

下面结合附图和实施例，对本发明的具体实施方式作进一步详细描述。以下实施例用于说明本发明，但不用来限制本发明的范围。

图1是本发明实施例提供一种农田蒸散估算方法流程图，如图1所示，包括：

S1、获取农田的遥感影像数据和气象数据；

S2、基于所述农田的遥感影像数据和气象数据，计算潜热通量、净辐射通量和土壤热通量；

S3、基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。

S1中，所述农田的遥感影像数据，是卫星的遥感数据，优选的，本发明实施例提供的是Landsat卫星热红外遥感数据，但本发明不对具体的卫星数据做具体限定。

S1中，所述气象数据为农田周边的气象观测战获得的气象数据，例如空气温度、风速等气象参数。可以理解的是，所述气象数据包括了本发明实施例计算过程中所需要的大量气象参数值，本发明实施例在此不一一列举。

S2中，所述潜热通量为地表水分和土壤水分相变过程所需的能量，所述净辐射通量为地表接收到的总输入能量，所述土壤热通量为土壤吸收的能量。

可以理解的是，净辐射通量最终分成三部分，一部分能量以对流和热传导加热空气，这部分能量为显热通量，一部分作为地表水分和土壤水分相变过程所需的能量，为潜热通量，还有一部分是被土壤吸收的能量，为土壤热通量。

那么通过能量余项法，可以由公式LE＝R_n-G-H表示，其中，R_n为净辐射通量(W/m²)，G是土壤热通量(W/m²)，H是显热通量(W/m²)，LE是潜热通量(W/m²)。

S3中，所述蒸散值一般指代为日蒸散值，但例如月蒸散值、年蒸散值，均可根据本发明实施例提供的农田蒸散方法进行估算。

具体的，首先计算出瞬时净辐射通量R_n,inst、瞬时土壤热通量G_inst、瞬时潜热通量LE_inst，然后根据公式计算蒸发比EF、日净辐射通量R_n,24以及全天的土壤热通量G_24，根据蒸发比不变的假设，进行尺度转换得到日蒸散量ET₂₄，具体公式为：

其中，LE_inst和ET₂₄分别是瞬时潜热通量值和日蒸散值，(R_n-G)_inst和(R_n-G)₂₄是瞬时可利用能量值和日可利用能量值，λ为汽化潜热。

可以理解的是，由于土壤白天吸收的能量和晚上辐射出的能量几乎相等，故全天的土壤热通量约为零，晴空条件下，辐射平衡的各分量日曲线在白天都呈正弦曲线变化，具体的计算全天净辐射R_n,24的公式如下：

其中，t为卫星获取影像的瞬间离日出时的时间长度(h)，N为日照时数，R_n,inst为瞬时净辐射通量(W/m²)。

本发明实施例提供的农田蒸散估算方法，通过充分利用气象数据与光学遥感和热红外遥感相结合，物理基础明确，空间分辨率高，估算结果具有更强的可靠性。

在上述实施例的基础上，步骤S2包括：

S21、基于所述农田的遥感影像数据，计算农田的植被指数；

S22、基于所述农田的植被指数和所述气象数据，计算地表温度；

S23、基于所述农田的植被指数、所述地表温度以及所述农田的遥感影像数据和气象数据，计算所述潜热通量、所述净辐射通量和所述土壤热通量。

其中，所述农田的植被指数，由近红外波段、可见光红色波段以及蓝色波段的反射率计算得到的，具体的，遥感数据来源于Landsat 5数据的第3、第4和第1波段或Landsat 8数据的第4、第5和第2波段。

其中，所述地表温度采用的是单窗算法，并且考虑地表温度与高程的相关关系和对反演的地表温度进行高程校正。

具体的，计算表达式为：

T_s＝{a(1-C-D)+[(b-1)(1-C-D)+1]T_b-DT_air}/C，

T_{s_DEM}＝T_s-0.006ΔZ，

其中，T_s是地表温度(K)，a和b是经验系数，C和D是中间变量，T_b是亮度温度(K)，T_air是大气平均作用温度(K)，通过近地表气温(T₀)的经验关系计算，T_{s_DEM}是校正后的地表温度(K)，ΔZ为像元点高程与参考高程之差(m)，参考高程一般取气象站高程。

可以理解的是，比辐射率根据步骤计算的NDVI计算所得，大气透射率、空气温度由气象数据中获取。

在上述实施例的基础上，步骤S23包括：

S231、基于所述农田的遥感影像数据和气象数据，计算净辐射通量；

S232、基于所述净辐射通量、所述地表温度，计算所述土壤热通量；

S233、基于所述农田的遥感影像数据和气象数据、所述净辐射通量以及所述土壤热通量，计算所述潜热通量。

具体的，S231中，所述净辐射通量是通过所有输入能量与所有输出能量作差所得，其表达式为：

R_n＝(1-α)R_s↓+R_L↓-R_L↑-(1-ε)R_L↓，

其中，α为地表反照率，R_s↓为下行的太阳短波辐射(W/m²)；R_L↓为大气长波辐射(W/m²)；R_L↑为地表长波辐射(W/m²)；ε为地表比辐射率。(1-ε)R_L↓为下行长波辐射被反射的那部分能量(W/m²)。

需要说明的是，计算净辐射通量中，地表反照率通过遥感影像各个波段的反射率加权平均得到反照率。

S232中，所述土壤热通量是根据地表反照率和地表温度之间的经验关系来计算，计算式为：

其中，R_n为净辐射通量(W/m²)，α是地表反照率，T_s是地表温度(K)。

在上述实施例的基础上，步骤S233包括：

基于所述农田的遥感影像数据和气象数据，计算显热通量；

基于能量余项法、所述显热通量、所述净辐射通量以及所述土壤热通量，计算所述潜热通量。

可以理解的是，所述潜热通量是由能量余项法的差值计算出的，计算式如下：

LE＝R_n-G-H，

其中，R_n为净辐射通量(W/m²)，G是土壤热通量(W/m²)，H是显热通量(W/m²)，LE是潜热通量(W/m²)。

所述显热通量是利用迭代运算来对影像的像元进行内部校准从而得到稳定的显热通量，通过地气之间的温度梯度和空气动力学阻抗估算所得，其表达式为：

其中，ρ_air为空气密度(kg/m³)，C_p为空气定压比热(J/(kg·K))，r_ah为空气动力学阻抗(s/m)，dT(K)是z₁和z₂高度之间的温差，z₁和z₂通常分别取0.01m和2m。

可以理解的是，计算显热通量时，需通过选取影像中的冷热像元，构建温差与地表温度直接的经验关系计算地气温差，具体公式为：

dT＝a+bT_{s_DEM}，

其中，T_{s_DEM}为经过数字高程校正后的地表温度(K)，a和b为经验系数，由遥感图像上的极端点(冷热像元)的相关参数计算的到。

在上述实施例的基础上，所述基于所述农田的遥感影像数据和气象数据，计算显热通量，包括：

基于所述农田的植被指数和地表粗糙度的估算模型，估算地表粗糙度；

基于所述地表粗糙度和所述农田的遥感影像数据和气象数据，计算所述显热通量。

可以理解的是，显热通量的估算依赖于空气动力学阻抗的稳定，而阻抗的计算中如何可靠估算地表粗糙度尤为重要。

现有技术采用地表粗糙度与叶面积指数的线性经验关系计算，但在实践中发现，这种线性经验关系仍然存在不稳定性。

针对上述现有技术存在的问题，本发明实施例根据地表粗糙度与光谱反射的相关性并结合相关研究，采用地表粗糙度与EVI的经验关系进行估算，其中，农田的植被指数和地表粗糙度的估算模型，即为该经验关系的估算表达式，所述估算表达式为：

z_om＝exp(a₁+b₁EVI)，

其中，z_om为地表粗糙度，a₁、b₁为回归参数。

本发明实施例采用在研究区选取一定的数量样点，测出样点的株高h，根据公式z_om＝0.123h以算出样点的粗糙度z_om，然后在EVI影像上找出对应样点的像元，得到样本点的EVI值，最后利用非线性拟合的方法算出回归参数a₁、b₁。

本发明实施例拟合出的关系式为：

z_om＝exp(-6.57+7.33EVI)，

可以理解的是，本发明实施例不对拟合的具体数值做限定。

在上述实施例的基础上，步骤S21具体包括：

基于所述农田的遥感影像数据，获取近红外波段的地表反射率ρ_nir、可见光红波段的地表反射率ρ_red和蓝光波段的地表反射率ρ_blue；

基于公式计算农田的植被指数；

其中，EVI为所述农田的植被指数。

可以理解的是，本发明实施例构建了一个EVI计算模块，根据获取的农田的遥感影像数据，自适应的计算农田的植被指数。

在上述实施例的基础上，步骤S231具体包括：

基于所述农田的植被指数和所述气象数据，获取地表反照率α、下行的太阳短波辐射R_s↓、大气长波辐射R_L↓、地表长波辐射R_L↑以及地表比辐射率ε；

基于公式R_n＝(1-α)R_s↓+R_L↓-R_L↑-(1-ε)R_L↓，计算净辐射通量R_n。

可以理解的是，本发明实施例构建了一个净辐射通量计算模块，通过获取的农田的植被指数和所述气象数据，自适应的计算净辐射通量。

在上述实施例的基础上，步骤S232具体包括：

基于所述净辐射通量、所述农田的植被指数以及所述地表温度，获取地表反照率α、地表温度T_s；

基于公式计算所述土壤热通量G。

可以理解的是，本发明实施例构建了一个土壤热通量计算模块，当获取到净辐射通量、所述农田的植被指数以及所述地表温度时，自适应的计算土壤热通量。

在上述实施例的基础上，所述基于所述农田的遥感影像数据和气象数据，计算显热通量，具体包括：

基于所述农田的遥感影像数据和气象数据，获取空气密度ρ_air、空气定压比热C_p、空气动力学阻抗r_ah，温差dT；

基于公式计算显热通量H。

本发明实施例通过具体的物理模型构建，系统准确的对各个通量进行计算，提高了农田蒸散估算的准确率。

图2是本发明实施例提供的各个观测点的测量的蒸散值散点图，具体的，本发明实施例选用中国海河流域地区作为实验场地，研究样带区域主要在海河流域，处于东经112°-120°，北纬35°-43°，地貌类型主要为高原、山地以及平原，地势西北高东南低，北部区域为蒙古高原区域，东部地区为华北平原区，而西北区域则为黄土高原和太行山区。属于温带东亚季风气候区，属半湿润半干旱地带，春季干旱多沙，夏季高温多雨，冬季寒冷少雪，年平均降水量540mm，年蒸散量470mm。

以海河流域的4个通量观测站点的周边为研究区，该研究区数据主要是2008年到2010年的Landsat5数据2013年到2014年的Landsat8数据以及对应时间段海河流域自动气象站观测数据集。

如图2所示，本发明实施例提供的蒸散估算方法平均误差为-0.35mm，相关系数平方R²为0.81，均方差为0.88mm。

由上述结果可以证明，本发明实施例提供的蒸散估算算法，估算结果较好，对测量蒸散量很有效，具有极高的可靠性。

图3是本发明实施例提供的一种农田蒸散估算系统结构图，如图3所示，包括：数据模块1、土壤热通量模块2、净辐射通量模块3、潜热通量模块4和农田蒸散估算模块5，其中，

数据模块1用于获取农田的遥感影像数据和气象数据；

土壤热通量模块2用于基于所述农田的遥感影像数据和气象数据，计算所述土壤热通量；

净辐射通量模块3用于基于所述农田的遥感影像数据和气象数据，计算所述净辐射通量；

潜热通量模块4用于基于所述农田的遥感影像数据和气象数据，计算所述潜热通量；

农田蒸散估算模块5用于基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。

具体的农田蒸散方法可参见上述实施例，本发明实施例在此不再赘述。

本发明提供的农田蒸散估算系统，通过充分利用气象数据与光学遥感和热红外遥感相结合，物理基础明确，空间分辨率高，估算结果具有更强的可靠性。

最后，本申请的方法仅为较佳的实施方案，并非用于限定本发明的保护范围。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种农田蒸散估算方法，其特征在于，包括：

S1、获取农田的遥感影像数据和气象数据；

S2、基于所述农田的遥感影像数据和气象数据，计算潜热通量、净辐射通量和土壤热通量；

S3、基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。

2.根据权利要求1所述的方法，其特征在于，步骤S2包括：

S21、基于所述农田的遥感影像数据，计算农田的植被指数；

S22、基于所述农田的植被指数和所述气象数据，计算地表温度；

S23、基于所述农田的植被指数、所述地表温度以及所述农田的遥感影像数据和气象数据，计算所述潜热通量、所述净辐射通量和所述土壤热通量。

3.根据权利要求2所述的方法，其特征在于，步骤S23包括：

S231、基于所述农田的遥感影像数据和气象数据，计算净辐射通量；

S232、基于所述净辐射通量、所述农田的植被指数以及所述地表温度，计算所述土壤热通量；

S233、基于所述农田的遥感影像数据和气象数据、所述净辐射通量以及所述土壤热通量，计算所述潜热通量。

4.根据权利要求3所述的方法，其特征在于，步骤S233包括：

基于所述农田的遥感影像数据和气象数据，计算显热通量；

基于能量余项法、所述显热通量、所述净辐射通量以及所述土壤热通量，计算所述潜热通量。

5.根据权利要求4所述的方法，其特征在于，所述基于所述农田的遥感影像数据和气象数据，计算显热通量，包括：

基于所述农田的植被指数和地表粗糙度的估算模型，估算地表粗糙度；

基于所述地表粗糙度和所述农田的遥感影像数据和气象数据，计算所述显热通量。

6.根据权利要求2所述的方法，其特征在于，步骤S21具体包括：

基于所述农田的遥感影像数据，获取近红外波段的地表反射率ρ_nir、可见光红波段的地表反射率ρ_red和蓝光波段的地表反射率ρ_blue；

基于公式计算农田的植被指数；

其中，EVI为所述农田的植被指数。

7.根据权利3所述的方法，其特征在于，步骤S231具体包括：

基于所述农田的植被指数和所述气象数据，获取地表反照率α、下行的太阳短波辐射R_s↓、大气长波辐射R_L↓、地表长波辐射R_L↑以及地表比辐射率ε；

基于公式R_n＝(1-α)R_s↓+R_L↓-R_L↑-(1-ε)R_L↓，计算净辐射通量R_n。

8.根据权利要求3所述的方法，其特征在于，步骤S232具体包括：

基于所述净辐射通量、所述农田的植被指数以及所述地表温度，获取地表反照率α、地表温度T_s；

基于公式计算所述土壤热通量G。

9.根据权利要求4所述的方法，其特征在于，所述基于所述农田的遥感影像数据和气象数据，计算显热通量，具体包括：

基于所述农田的遥感影像数据和气象数据，获取空气密度ρ_air、空气定压比热C_p、空气动力学阻抗r_ah，温差dT；

基于公式计算显热通量H。

10.一种农田蒸散估算系统，其特征在于，包括：

数据模块，用于获取农田的遥感影像数据和气象数据；

土壤热通量模块，用于基于所述农田的遥感影像数据和气象数据，计算所述土壤热通量；

净辐射通量模块，用于基于所述农田的遥感影像数据和气象数据，计算所述净辐射通量；

潜热通量模块，用于基于所述农田的遥感影像数据和气象数据，计算所述潜热通量；

农田蒸散估算模块，用于基于所述潜热通量、所述净辐射通量和所述土壤热通量，估算蒸散值。