蔣 磊 蔡甲冰 張寶忠 許 迪 魏 征
(1.天津農學院水利工程學院,天津 300392; 2.中國水利水電科學研究院流域水循環模擬與調控國家重點實驗室,北京 100038;3.國家節水灌溉北京工程技術研究中心,北京 100048)
作物蒸散發(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數據和地表作物冠層溫度實測值進行區域作物實際蒸散發量估算的可行性,為區域玉米精量灌溉決策提供理論基礎。
試驗區位于寧夏回族自治區中衛市沙坡頭區常樂鎮農業示范區(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所示。
根據研究區2019年玉米分布情況,選取2個代表性地塊,在每個地塊中間布置1套CTMS-On line型作物冠層溫度及環境因子測量系統(分別為H1和H2,圖1)。觀測系統田間實物如圖3所示,數據采集時間間隔均為30 min。本系統同步連續監測的農田信息包括:空氣溫度/濕度、風速、太陽輻射、光合有效輻射、大氣壓強、作物冠層紅外溫度、根區土壤墑情等[24]。每7~10 d人工進行觀測和采集玉米株高、干物質量、葉面積等作物信息。地下水位觀測值通過寧夏水利科學研究院在示范區安裝的觀測井獲得。其他氣象數據(如降水量等)從國家氣象科學數據中心下載得到。
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
基于雙源模型的原理,忽略光合作用耗能和水平方向能量交換的地表能量平衡方程為
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]。
由遙感蒸散發模型計算得到的是衛星過境瞬時的潛熱通量,需要通過一定的方法將瞬時潛熱通量進行時間尺度擴展來獲得日蒸散發量及更長時段內的蒸散發量。采用修正蒸發比法對瞬時潛熱通量進行時間尺度擴展得到衛星過境日的蒸散發量,計算公式為
EF=λETi/(Rn-G)i
(7)
λETd=mEF(Rn-G)d
(8)
式中EF——蒸發比
m——修正系數,取1.1
λ——水的汽化潛熱,取2.45 MJ/kg
下標i、d表示瞬時值和日值。
對于缺乏遙感數據的時段,蒸散發量采用日參考蒸發比三次樣條插值法得到。
采用水量平衡方法對遙感反演作物蒸散發量進行驗證和評價,水量平衡公式為
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中進行。
遙感蒸散發計算的關鍵輸入數據包括寬波段反照率(α)、歸一化植被指數(NDVI)、株高、地表冠層溫度、氣象數據等。其中α和NDVI通過遙感數據反演得到,株高通過NDVI進行估算;地表溫度(玉米冠層溫度)和氣象數據通過田間實測獲得。
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.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
Sentinel-2數據在可見光、近紅外和短波紅外工作譜段范圍內具有高時空分辨率,能夠很好地描述下墊面的高異質性,大大降低了混合像元的存在。但是地表溫度作為遙感蒸散發模型的關鍵輸入數據,目前仍存在著“時、空”上的矛盾,在很大程度上增加了蒸散發估算的不確定性。隨著遙感技術的發展以及地表溫度解譯方法的成熟,建立高時空分辨率的地表溫度算法和產品是提高區域蒸散發估算精度的重要課題。
在瞬時蒸散發到日蒸散發量的時間尺度擴展方面,不同方法均有其適用條件,如何根據不同的環境因素、衛星過境時刻、下墊面類型等選取適宜的時間尺度擴展方法仍值得深入研究。此外,利用擴展得到的日蒸散發量獲取連續時段的蒸散發量的插值方法仍需進一步研究,特別是針對有降水前后由于云遮蔽的影響缺乏遙感數據,同時又由于降水導致的土壤含水率突變的情況下,如何提高插值的精度也是下一步需繼續開展和完善的工作。
(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)獲取研究區玉米生育期蒸散發量(蒸發量、蒸騰量)具有較高的精度,對精量灌溉決策的制定和區域農業水管理具有重要的科學意義和應用價值。