基于地面观测和Sentinel-2数据的玉米实际蒸散发估算

蒋 磊 蔡甲冰 张宝忠 许 迪 魏 征

(1.天津农学院水利工程学院，天津 300392； 2.中国水利水电科学研究院流域水循环模拟与调控国家重点实验室，北京 100038；3.国家节水灌溉北京工程技术研究中心，北京 100048)

0 引言

作物蒸散发(Evapotranspiration，ET)是农业水管理的重要组成部分，合理估算区域作物蒸散发对作物用水效率评价和精量灌溉决策具有重要的意义[1-2]。随着遥感技术的发展，基于遥感数据和地面观测数据建立的蒸散发模型为科学、准确获取区域作物ET信息提供了一种有效的途径和方法[3-4]。

在区域作物遥感蒸散发估算方法中，发展较为成熟的是基于地表垂直方向能量平衡的思路，即将蒸散发(潜热)作为能量平衡中的余项求出[5-6]。根据不同的辐射通量和湍流交换方式，能量平衡余项法可分为单源模型和双源模型。单源模型又称“大叶模型”，即不区分土壤蒸发和植被蒸腾，将土壤和植被视为一个均一的整体来考虑潜热通量，常见的单源模型包括SEBAL(Surface energy balance algorithm for land)模型[7]和SEBS(Surface energy balance system)模型[8]等。与单源模型不同，双源模型分别计算土壤和冠层两个组分的净辐射、显热和潜热，进而能够区分土壤蒸发和植被蒸腾。NORMAN等[9]开发的TSEB(Two-source energy balance)模型是双源模型的代表。YANG等[10]在前人研究的基础上提出了混合双源蒸散发模型(Hybrid dual-source scheme and trapezoid framework-based evapotranspiration model，HTEM)，该模型已在国内外多个地区进行验证，并取得较高的精度[3,11]，与TSEB、MOD16-ET产品的对比表明，其计算精度高于其他两类模型[12]。

作为遥感蒸散发模型的主要输入信息，遥感数据对模型的估算精度有重要影响。目前，常用于估算区域作物蒸散发的遥感数据主要有MODIS(Moderate resolution imaging spectroradiometer)数据[13-14]、Landsat ETM+数据[15]、Sentinel-2数据[16]等。国内灌区常见地形、田块和种植结构复杂，往往需要高时间分辨率和高空间分辨率的遥感数据，以实现农田作物蒸散发的精准、连续监测和反演[17]。尽管可以利用一些数据融合算法调和MODIS数据和Landsat ETM+数据在空间分辨率和时间分辨率上的矛盾[18]，但复杂的算法和步骤限制了其大范围实际应用和推广。Sentinel-2数据具备高时间分辨率(5 d)和高空间分辨率(10～60 m)，弥补了MODIS数据和Landsat ETM+数据在“时、空”分辨率上的缺陷。近年来，Sentinel-2数据在作物识别、植被指数反演等方面得到了广泛应用[19-22]。但是基于Sentinel-2数据对作物蒸散发量进行估算的研究还较为鲜见。

地表温度(Land surface temperature, LST)作为双源模型中区分植被温度和土壤温度的关键因子，是作物需水状况和土壤水分信息的综合反映。目前，对于遥感地面温度数据同样存在着“时、空”分辨率矛盾的问题，因此，高时空分辨率的LST数据对提高遥感作物蒸散发量估算精度具有重要价值。基于上述问题，本研究利用混合双源蒸散发模型(HTEM)，结合Sentinel-2数据和地面农田实际观测数据，对宁夏回族自治区中卫市玉米主要生育期(5—8月)的实际蒸散发量进行计算，以探究利用Sentinel-2数据和地表作物冠层温度实测值进行区域作物实际蒸散发量估算的可行性，为区域玉米精量灌溉决策提供理论基础。

1 材料和方法

1.1 试验区概况

试验区位于宁夏回族自治区中卫市沙坡头区常乐镇农业示范区(105°7′6.76″E, 37°28′5.7″N)(图1)。研究区属于温带大陆性季风气候，海拔1 225 m。受沙漠影响，日照充足，昼夜温差大，年内平均气温为7.3～9.5℃，年均降水量180～367 mm，年蒸发量为1 829 mm，夏季降水量占全年降水量的61%，全年日照时长2 800 h左右[23]。当地土壤主要为灌淤土。研究区种植的农作物以玉米为主，土地面积约为566.95 hm2。于2019年在玉米的主要生育期(5—8月)进行试验观测，2019年5月初播种，9月中旬收青贮玉米。研究区内自然降雨时段与玉米整个生育期吻合；玉米田块共灌溉5次，总灌溉量为380.81 mm，灌溉日期及灌水量如图2所示。

1.2 地面数据监测与采集

根据研究区2019年玉米分布情况，选取2个代表性地块，在每个地块中间布置1套CTMS-On line型作物冠层温度及环境因子测量系统(分别为H1和H2，图1)。观测系统田间实物如图3所示，数据采集时间间隔均为30 min。本系统同步连续监测的农田信息包括：空气温度/湿度、风速、太阳辐射、光合有效辐射、大气压强、作物冠层红外温度、根区土壤墒情等[24]。每7～10 d人工进行观测和采集玉米株高、干物质量、叶面积等作物信息。地下水位观测值通过宁夏水利科学研究院在示范区安装的观测井获得。其他气象数据(如降水量等)从国家气象科学数据中心下载得到。

1.3 Sentinel-2数据下载及处理

Sentienl-2数据来源包括Sentienl-2A和Sentienl-2B卫星，各波段信息如表1所示。空间分辨率为10～60 m，时间分辨率为5 d。Sentinel-2数据下载及处理流程主要包括原始数据下载(https:∥scihub.copernicus.eu)和校正、辐射定标和裁剪、反演相关参数等3个关键步骤。数据下载完成后利用SNAP平台对原始数据进行处理，包括几何校正、辐射定标、重采样和裁剪等处理。并进一步反演获得遥感蒸散发模型所需的反照率、归一化植被指数(NDVI)、植被覆盖度(fc)等关键参数。通过去除生育期云遮盖等的影响，共选取10景数据进行研究，如表2所示。

表1 Sentinel-2数据各波段信息Tab.1 Statistical band information of Sentinel-2 data

表2 选用的Sentinel-2数据Tab.2 Selections of Sentinel-2 remote sensed data

1.4 混合双源蒸散发模型(HTEM)

基于双源模型的原理，忽略光合作用耗能和水平方向能量交换的地表能量平衡方程为

Rn=H+LE+G

(1)

式中Rn——净辐射通量，W/m2

H——感热通量，W/m2

LE——潜热通量，W/m2

G——土壤热通量，W/m2

潜热通量(LE)作为能量余项，可通过求得地表能量平衡方程中的Rn、H和G求得。

净辐射通量采用ALLEN等[25]提出的基于遥感信息的地表净辐射计算方法获得。土壤热通量采用BASTIAANSSEN[26]提出的半经验模型进行估算，其表达式为

G=Rn(Lst-273.16)(0.003 8+0.007 4α)·

(1-0.98NDVI4)

(2)

其中

α=0.356ρ2+0.130ρ4+0.373ρ8+
0.085ρ11+0.072ρ12-0.001 8

(3)

(4)

式中Lst——地表温度，K

α——地表宽波段反照率[27]

NDVI——归一化植被指数

其中，ρi为各窄波段反射率。地表宽波段反照率需用Band 2、Band 4、Bnad 8、Band 11和Band 12波段的反射率；归一化植被指数需用Band 4和Band 8波段的反射率。

HTEM模型采用层状模型对净辐射在植被和土壤组分间进行分配，采用块状模型计算地表的显热通量、潜热通量及土壤热通量。植被和土壤潜热通量计算式为

(5)

(6)

式中LEc、LEs——植被、土壤潜热通量，W/m2

Rnc、Rns——植被、土壤净辐射，W/m2，可通过消光系数计算得到

ρair——空气密度，kg/m3

Cp——空气比热容，J/(kg·K)

Tc、Ts——植被、土壤表面温度，K，可由植被指数-地表温度梯形特征空间来确定

Ta——空气温度，K

参数具体计算方法和步骤见文献[10]。

1.5 蒸散发时间尺度扩展

由遥感蒸散发模型计算得到的是卫星过境瞬时的潜热通量，需要通过一定的方法将瞬时潜热通量进行时间尺度扩展来获得日蒸散发量及更长时段内的蒸散发量。采用修正蒸发比法对瞬时潜热通量进行时间尺度扩展得到卫星过境日的蒸散发量，计算公式为

EF=λETi/(Rn-G)i

(7)

λETd=mEF(Rn-G)d

(8)

式中EF——蒸发比

m——修正系数，取1.1

λ——水的汽化潜热，取2.45 MJ/kg

下标i、d表示瞬时值和日值。

对于缺乏遥感数据的时段，蒸散发量采用日参考蒸发比三次样条插值法得到。

1.6 精度评价

采用水量平衡方法对遥感反演作物蒸散发量进行验证和评价，水量平衡公式为

ET=Pr+U+I-D-R-ΔW

(9)

式中ET——作物蒸散发量，mm

Pr——有效降雨量，mm

U——地下水补给量，mm

I——灌水量，mm

D——深层渗漏量，mm

R——径流量，mm

ΔW——试验初期和试验末期土壤储水量的变化量，mm

因试验区属于干旱半干旱地区且地势平坦，R和D可忽略不计，则水量平衡公式可简化为

ET=Pr+U+I-ΔW

(10)

有效降雨量采用有效降水系数和实际降水量的乘积得到，试验区5—8月有效降水系数分别为69.87%、78.17%、81.15%、88.16%[26]。地下水补给量采用地下水给水度和研究时段内地下水位变化量乘积得到，本研究中地下水给水度取0.045[28]。土壤储水量的变化量采用试验末期和试验初期土壤储水量差值求得，其中土壤储水量计算公式为

(11)

式中W——土壤储水量，mm

θi——第i层土壤体积含水率，%

hi——第i层土壤厚度，cm

n——土层的层数，本试验设0～20 cm、20～40 cm和40～60 cm共计3个土层

选用相对误差(Relative error，RE)、绝对误差(Absolute error，AE)和均方根误差(Root mean square error，RMSE)作为精度评价指标，其中相对误差、绝对误差和均方根误差越小，说明估算的精度越高。采用OriginPro 9.1软件制图，统计参数分析在Excel 2017中进行。

2 结果与分析

遥感蒸散发计算的关键输入数据包括宽波段反照率(α)、归一化植被指数(NDVI)、株高、地表冠层温度、气象数据等。其中α和NDVI通过遥感数据反演得到，株高通过NDVI进行估算；地表温度(玉米冠层温度)和气象数据通过田间实测获得。

2.1 遥感参数反演结果对比

2.1.1反照率和NDVI反演结果

Sentinel-2数据反演研究区反照率和NDVI结果如图4所示。以2019年6月24日遥感数据为例，研究区农田宽波段反照率主要集中在0.1～0.2之间，农田NDVI主要集中在0.5～0.8之间，其中玉米观测站点南部有部分田块NDVI明显低于周围田块，主要是由于2019年该部分田块种植作物为枸杞。枸杞属灌木，NDVI小于玉米等农作物，反演结果空间分布与实际种植情况相吻合。另外，从反照率和NDVI反演结果可以看出，Sentinel-2数据具有较高的空间分辨率(反照率空间分辨率20 m，NDVI空间分辨率10 m)，适合用于研究区农田地块复杂、种植规模小的情况，降低了混合像元的存在，提高了反演精度。

2.1.2玉米NDVI时间序列拟合

研究区玉米试验田块NDVI时间序列拟合曲线选取DoseResp“S”形曲线，迭代算法采用Levenberg-Marquardt优化算法。DoseResp曲线表达式为

(12)

式中Z——儒略日(DOY)

a、b、c、e——待拟合参数，其中a为曲线下渐近线，b为曲线上渐近线

对研究区玉米NDVI时间序列进行拟合，结果如图5和表3所示。可以看出，“S”形曲线能够很好地描述研究区玉米NDVI时间序列，决定系数R2达到0.99以上。从玉米NDVI时间序列变化曲线来看，玉米播种后至出苗前期，NDVI约为0.2，表现为与裸土NDVI接近。随着玉米生育期的增加，在第160天到第180天，NDVI由0.2迅速增至0.8，这一阶段为玉米快速生长阶段。此后至收获青贮玉米之前，NDVI一直保持在0.8左右。由拟合参数可以看出，a和b表示拟合曲线理论上的最大值和最小值，玉米NDVI时间序列理论最大值和最小值分别为0.797和0.183。

表3 2019年研究区玉米NDVI时间序列曲线拟合结果Tab.3 Maize NDVI series of study area and logistic curve fitting results in 2019

2.1.3玉米株高与NDVI关系

遥感蒸散发模型中，空气动力学阻抗的计算需要冠层高度作为输入数据，2019年共进行10次玉米株高观测，分别为5月19日、6月4日、6月14日、6月24日、7月4日、7月18日、7月30日、8月9日、8月19日和8月29日。玉米株高与NDVI关系如图6所示，并利用指数曲线对株高与NDVI关系进行拟合。由拟合结果可以看出，玉米株高与NDVI呈指数关系，决定系数为0.94。

2.2 计算结果与验证

2.2.1遥感净辐射和蒸散发计算结果

2019年研究区玉米站点(H1和H2)的净辐射观测值和卫星过境瞬时净辐射计算结果如图7所示。可见，研究时段内卫星过境瞬时净辐射变化范围为590～710 W/m2，利用遥感信息获取地表净辐射具有较高的精度，计算值与观测值决定系数为0.786，均方根误差(RMSE)为36.256 W/m2。

利用HTEM模型计算得到的研究区内2个玉米试验田块瞬时潜热通量和采用修正的蒸发比法得到的日蒸散发量结果如表4(E为土壤蒸发量，T为植被蒸腾量)和图8所示。可见，2个观测站点(H1和H2)在玉米生育期内潜热通量和日蒸散发量变化趋势相似，生育初期净辐射量约为600 W/m2，潜热通量约为350 W/m2，潜热通量主要为土壤潜热，占总潜热通量的95%以上。时间尺度扩展得到的日蒸散发量变化范围为3.742～4.676 mm/d，水分消耗主要以土壤蒸发为主，范围为3.470～4.303 mm/d，植被蒸腾很小，变化范围为0.112～0.373 mm/d。在玉米快速生长阶段，净辐射逐渐增加，约为700 W/m2，潜热通量也随之增加，变化范围为495.481～613.873 W/m2。日实际蒸散发量在这一阶段也有明显增加，变化范围为5.573～5.999 mm/d，主要的水分消耗形式从土壤蒸发逐渐变为植被蒸腾，在第175天时，植被蒸腾量占日蒸散发量的95%以上，并且这一比例一直保持到青贮玉米收获。玉米日实际蒸散发量在第185天达到峰值，2个玉米观测站点的日实际蒸散发量分别为6.311 mm/d和6.547 mm/d。随后日实际蒸散发量逐渐减小至4.106 mm/d和4.617 mm/d。

采用日参考蒸发比插值的方法获取无卫星过境日的蒸散发量，进而得到生育期内连续的玉米日实际蒸散发量，采用三次样条插值，并对插值后的参考蒸发比异常值进行修正，得到玉米生育期内连续的日实际蒸散发量(图9)。结果表明，逐日实际蒸散发量同样呈现出单峰的变化趋势，即生育初期日实际蒸散发量在4 mm/d左右浮动。随着生育期的增加，日实际蒸散发量逐渐增大，在第180天到第220天达到峰值，期间个别时段日实际蒸散发量超过8 mm/d，这是由于这一阶段降水较为集中，其中年内两次较大的降水均发生在该时段，并且进行了两次灌溉，土壤含水率较高，水分胁迫较小，加之这一时段作物潜在蒸散发量较大，导致这一阶段实际蒸散发量达到了峰值。至生育期末，日实际蒸散发量逐渐减小为2～4 mm。从统计结果可以看出(表5)，研究区2个玉米试验田块5—8月实际蒸散发总量分别为525.114 mm和533.690 mm，日均实际蒸散发量为4.269 mm/d和4.339 mm/d，其中植被蒸腾总量分别为363.483 mm和358.196 mm，从整个生育期来看，植被蒸腾是作物水分消耗的主要形式。

表5 2019年研究区玉米实际蒸散发量Tab.5 Statistical results of maize evapotranspiration of study area in 2019

2.2.2ET计算结果验证

采用水量平衡法对计算得到的玉米生育期内蒸散发总量进行验证和评价(表6)。结果表明，基于冠层温度观测系统和Sentinel-2数据反演的研究区玉米蒸散发总量与水量平衡法计算得到的蒸散发总量相比，精度较高，2个观测站点的相对误差分别为-2.512%和-1.469%，绝对误差分别为13.533、7.774 mm。可见，利用当地田间多参数观测系统和Sentinel-2数据并结合混合双源蒸散发模型获取作物蒸散发量比较可靠。

表6 2019年研究区玉米蒸散发量验证结果Tab.6 Verification of maize evapotranspiration in study area in 2019

2.3 讨论

Sentinel-2数据在可见光、近红外和短波红外工作谱段范围内具有高时空分辨率，能够很好地描述下垫面的高异质性，大大降低了混合像元的存在。但是地表温度作为遥感蒸散发模型的关键输入数据，目前仍存在着“时、空”上的矛盾，在很大程度上增加了蒸散发估算的不确定性。随着遥感技术的发展以及地表温度解译方法的成熟，建立高时空分辨率的地表温度算法和产品是提高区域蒸散发估算精度的重要课题。

在瞬时蒸散发到日蒸散发量的时间尺度扩展方面，不同方法均有其适用条件，如何根据不同的环境因素、卫星过境时刻、下垫面类型等选取适宜的时间尺度扩展方法仍值得深入研究。此外，利用扩展得到的日蒸散发量获取连续时段的蒸散发量的插值方法仍需进一步研究，特别是针对有降水前后由于云遮蔽的影响缺乏遥感数据，同时又由于降水导致的土壤含水率突变的情况下，如何提高插值的精度也是下一步需继续开展和完善的工作。

3 结论

(1)Sentinel-2数据具有高时空分辨率，能够与研究区复杂的种植地块相匹配，降低了混合像元的存在，提高了玉米蒸散发量估算的精度。利用“S”形曲线拟合研究区玉米NDVI时间序列具有较高精度，决定系数达到0.99以上。同时，研究区玉米株高与NDVI的关系可以用指数曲线拟合，决定系数达到0.94。

(2)基于Sentinel-2数据反演得到的卫星过境瞬时净辐射具有较高精度，其变化范围为590～710 W/m2，与地面观测值间的决定系数为0.786，均方根误差为36.256 W/m2。

(3)卫星过境日实际蒸散发量呈单峰变化，在第185天左右达到峰值，日蒸散发量达到6 mm/d以上。生育初期水分主要以土壤蒸发形式消耗，从第150天之后逐渐过渡到以植被蒸腾为主，到第175天时，植被蒸腾量占日实际蒸散发量的95%以上，并且一直保持到青贮玉米收获。由参考蒸发比三次样条插值得到的逐日实际蒸散发量结果表明，研究区玉米5—8月实际蒸散发总量为529.402 mm，日均实际蒸散发量为4.304 mm/d，其中植被蒸腾总量为360.840 mm，从整个生育期来看，植被蒸腾是作物水分消耗的主要形式。采用水量平衡法对玉米生育期蒸散发总量进行验证，两个观测站点的平均相对误差和绝对误差分别为-1.990 5%和10.653 5 mm。

(4)将地面冠层温度观测数据和高时空分辨率遥感数据相结合，采用混合双源遥感蒸散发模型(HTEM)获取研究区玉米生育期蒸散发量(蒸发量、蒸腾量)具有较高的精度，对精量灌溉决策的制定和区域农业水管理具有重要的科学意义和应用价值。