李彥彬,梁 坤,張文鴿
(1.華北水利水電大學(xué),河南 鄭州 450046;2.黃河水利科學(xué)研究院,河南 鄭州 450003)
蒸散發(fā)(Evapotranspiration,ET)包括水體、土壤、植被表面截留水的蒸發(fā)和植物的蒸騰[1],是地表水循環(huán)與能量循環(huán)的重要組成部分,也是連接生態(tài)過程和水文過程的關(guān)鍵環(huán)節(jié)。蒸散發(fā)的準(zhǔn)確估算,對于區(qū)域水資源的合理配置、干旱監(jiān)測、農(nóng)業(yè)精準(zhǔn)灌溉具有重要意義[2]。
傳統(tǒng)的蒸散發(fā)估算方法多以站點觀測為主,無法滿足區(qū)域尺度的研究要求。隨著近幾十年遙感技術(shù)的不斷發(fā)展,具有多種時間和空間分辨率的衛(wèi)星傳感器已經(jīng)成為獲取區(qū)域蒸散發(fā)的理想工具[3]。根據(jù)不同的理論基礎(chǔ),基于遙感的蒸散發(fā)模型大致可以分為三類[4]。①經(jīng)驗統(tǒng)計模型[5],利用遙感數(shù)據(jù)與蒸散發(fā)之間的經(jīng)驗關(guān)系進(jìn)行估算。這類模型數(shù)據(jù)量少、使用方便,但模型的統(tǒng)計參數(shù)在不同的區(qū)域都需重新校正,具有明顯的區(qū)域局限性。②特征空間模型,研究學(xué)者們發(fā)現(xiàn)反演得到的地表溫度與植被指數(shù)的散點圖多呈現(xiàn)出梯形或三角形的形狀,如Jiang等[6]提出的三角特征空間模型。這類模型的原理為有植被覆蓋地區(qū)的地表溫度要遠(yuǎn)低于裸土地區(qū),由此構(gòu)建的特征空間形狀可以確定出區(qū)域蒸散發(fā)的上下邊界,然后根據(jù)像元在特征空間的位置以及離上下邊界的距離系數(shù),反演出該像元的實際蒸散發(fā)。③能量平衡余項模型,基于地表能量平衡方程,將潛熱通量作為剩余項計算得到,該類模型又可分為單源模型和雙源模型[7]。雙源模型將土壤參數(shù)和植被參數(shù)分離,能更好反映植被特征對蒸散發(fā)的影響,但模型復(fù)雜、參數(shù)眾多,模型的不確定性也會相應(yīng)增加[8];相比之下,單源模型結(jié)構(gòu)簡單,模型參數(shù)更容易獲得,因而應(yīng)用較為廣泛。
陸面地表能量平衡(Surface Energy Balance Algorithm for Fand,SEBAL)模型作為典型的單源模型,具有適用性廣、反演精度高的特點,被用于不同地區(qū)的蒸散發(fā)研究。曾麗紅[9]、高海東等[10]利用不同數(shù)據(jù)源結(jié)合SEBAL模型,在流域和平原地區(qū)進(jìn)行了日蒸散發(fā)的反演;劉蓉等[11]利用SEBAL模型結(jié)合風(fēng)云氣象衛(wèi)星,反演了藏北高原地區(qū)各地表通量的空間日變化過程;在干旱半干旱區(qū),喬成龍等[12]利用SEBAL模型對寧夏鹽池縣人工灌叢與草地的蒸散量做了對照研究,評估了人工灌叢植被對陸地蒸散的影響;在灌區(qū)研究方面,SEBAL模型多用于農(nóng)田蒸散反演和用水管理[13]。如:Wei等[14]將SEBAL模型應(yīng)用于濕潤地帶的贛撫平原灌區(qū),分析了不同水文年早、中、晚熟水稻的蒸散發(fā)時空變化特征。
河套灌區(qū)屬于干旱大陸性氣候,氣候干燥,水資源短缺,農(nóng)業(yè)生產(chǎn)主要依靠引黃河水灌溉。蒸散發(fā)是該地區(qū)水分的主要消耗項,準(zhǔn)確估算出灌區(qū)尺度的蒸散發(fā)量就顯得十分必要。目前,針對河套灌區(qū)的蒸散發(fā)研究已取得部分成果,但難以準(zhǔn)確反映出蒸散發(fā)的空間異質(zhì)性。土地利用類型對蒸散發(fā)有著顯著影響,不同土地類型下蒸散發(fā)的分布亦有所差異。因此,明確蒸散發(fā)的空間格局變化特征,對于灌區(qū)水資源的高效利用具有重要意義。選擇河套灌區(qū)永濟(jì)灌域作為研究區(qū)域,利用Landsat衛(wèi)星數(shù)據(jù),基于SEBAL遙感蒸散發(fā)模型,估算永濟(jì)灌域2019年作物生長季(4—10月)的蒸散發(fā)量,并與該地區(qū)的土地利用數(shù)據(jù)結(jié)合,分析不同土地利用類型下蒸散發(fā)的時空變化特征,為區(qū)域農(nóng)田用水量提供參考和水資源的合理高效配置提供決策依據(jù)。
以內(nèi)蒙古河套灌區(qū)永濟(jì)灌域為研究區(qū),永濟(jì)灌域位于河套平原的中部,東與義長灌域為鄰,西與解放閘灌域毗鄰,南鄰黃河,北抵陰山,地理坐標(biāo)40°36′~41°13′N,107°13′~107°42′E。灌域南北長60 km,東西寬約40 km,涉及臨河區(qū)、五原縣及烏拉特中旗。該地區(qū)屬于干旱大陸性氣候,氣候干燥,降雨量少,日照時間長,蒸發(fā)強(qiáng)烈,年均氣溫8.9℃,年均降雨量135.2 mm,其中夏季降水量占全年降水總量的60%左右。地面高程在974~1 080 m,地形較為平緩,地勢呈現(xiàn)西南高、東北低的特點,見圖1。灌域總土地面積約1.836×105hm2,土地利用類型多以農(nóng)田為主,主要種植小麥、玉米和葵花等作物。
圖1 研究區(qū)高程及地理位置
1.2.1遙感數(shù)據(jù)及預(yù)處理
Landsat系列衛(wèi)星影像具有空間分辨率高,時間覆蓋連續(xù),獲取免費(fèi)等特點,是進(jìn)行遙感蒸散發(fā)反演的理想數(shù)據(jù)源之一。本文采用分辨率為30 m,過境周期為16 d的Landsat8 OIL/TIRS及Landsat7 ETM+遙感影像,數(shù)據(jù)獲取于USGS官方(https://earthexplorer.usgs.gov/)。選取2019年4—10月共7期質(zhì)量良好(晴空或少量云覆蓋)的遙感影像作為反演對象,根據(jù)研究區(qū)控制面積,需要影像軌道號為129-31和129-32的兩景影像鑲嵌得到,見表1。
表1 Landsat影像數(shù)據(jù)
遙感影像反演前,需要進(jìn)行預(yù)處理。遙感影像預(yù)處理包括輻射定標(biāo)、大氣校正、圖像裁剪及鑲嵌。對獲取的原始遙感影像,首先進(jìn)行輻射定標(biāo),以確定影像準(zhǔn)確輻射值。接著利用ENVI軟件下的FLAASH模型進(jìn)行大氣校正,這一步驟是為了將定標(biāo)后的輻射亮度值轉(zhuǎn)化為實際地表反射率,以消除大氣反射、散射、吸收等引起的誤差。然后將處理后的多時相影像鑲嵌并裁剪,得到永濟(jì)灌域2019年7景遙感影像。
1.2.2DEM數(shù)據(jù)
數(shù)字高程數(shù)據(jù)DEM來自地理空間數(shù)據(jù)云網(wǎng)站(http://www.gscloud.cn/),選取分辨率為30 m的ASTER GDEM V2全球數(shù)字高程數(shù)據(jù)。根據(jù)研究區(qū)跨度,在ArcGIS中將兩幅DEM數(shù)據(jù)拼接在一起,并利用矢量文件裁剪出研究區(qū)高程影像。
1.2.3氣象數(shù)據(jù)
本文的氣象數(shù)據(jù)均在中國數(shù)據(jù)氣象網(wǎng)(http://data.cma.cn)下載。收集研究區(qū)臨河氣象站點2019年4—10月的氣象數(shù)據(jù),包括最高氣溫、最低氣溫、平均氣溫、平均風(fēng)速、相對濕度和日照時數(shù),見表2。
表2 研究區(qū)各氣象數(shù)據(jù)
1.2.4土地利用數(shù)據(jù)
土地利用數(shù)據(jù)來源于武漢大學(xué)團(tuán)隊制作的2019年中國土地覆蓋30 m分辨率數(shù)據(jù)集[15],總體分類精度為79.31%。結(jié)合永濟(jì)灌域的實際情況,將土地利用類型分為五類,即耕地、草地、水體、荒地、城鄉(xiāng)用地。
SEBAL(Surface Energy Balance Algorithm for Land)模型由荷蘭學(xué)者Bastiaanssen開發(fā),是一種基于能量平衡的遙感蒸散發(fā)模型。該模型利用遙感數(shù)據(jù)和少量氣象數(shù)據(jù)估算出凈輻射通量、土壤熱通量和顯熱通量,將潛熱通量作為能量平衡的余項求出,計算見式(1)。
λET=Rn-G-H
(1)
式中Rn——地表凈輻射通量,W/m2;G——土壤熱通量,W/m2;H——感熱通量,W/m2;λET——潛熱通量,W/m2。
2.1.1地表凈輻射通量
地表凈輻射通量是單位面積內(nèi)的入射輻射能量與出輻射能量的差值,代表著地表可用的實際輻射能量,見式(2)—(5)。
Rn=(1-α)Rs↓+RL↓-RL↑-(1-ε)RL↓
(2)
(3)
RL↓=εaσTa4
(4)
RL↑=εσTs4
(5)
式中Rs↓——入射短波輻射,W/m2;RL↓——入射長波輻射,W/m2;RL↑——出射長波輻射,W/m2;α——地表反照率,無量綱;ε——地表比輻射率,無量綱;Gsc——太陽常數(shù),取1 367 W/m2;θ——太陽天頂角;dr——日地距離,天文單位;τsw——大氣單向透射率,無量綱;εa——大氣比輻射率(無量綱);σ——斯蒂芬-玻爾茲曼常數(shù),取5.67×10-8W/(m2·K4);Ta——空氣溫度,K;Ts——地表溫度,K。
2.1.2土壤熱通量
土壤熱通量是指由于傳導(dǎo)作用而儲存在土壤和植被中的那部分能量,是一個相對較小的量,由式(6)計算得到[16]。
(6)
式中,NDVI為歸一化植被指數(shù)。
2.1.3感熱通量
感熱通量是指由于對流和傳導(dǎo)作用而從地球表面散失到大氣中的能量,表征下墊面與大氣間湍流形式的熱交換,見式(7)—(11)。
(7)
(8)
(9)
(10)
dT=aTs+b
(11)
式中ρa(bǔ)ir——空氣密度,kg/m3;Cp——空氣熱量常數(shù),1 004 J/kg/K;dT——溫差,K;rah——空氣動力學(xué)阻力,s/m;U*——摩擦風(fēng)速,m/s;Z0——地面粗糙度,m;U200——距離地面200 m處的風(fēng)速[17];U2——距離地面2 m處的風(fēng)速,由氣象數(shù)據(jù)獲得;Z1——近地表高度,Z2——參考高度,Z1=0.1 m,Z2=2 m;k——波爾茨曼常數(shù),取值0.41;a、b——回歸系數(shù)。
計算a、b需要引入2個極端像元點,即“冷點”和“熱點”?!袄潼c”處水分充足,植被狀況良好,地表溫度較低,感熱通量幾乎為0;“熱點”處植被覆蓋很少,地表溫度較高,潛熱通量幾乎為0。同時,由于近地表大氣的不穩(wěn)定性,還需要引入Monin-Obukhov長度通過循環(huán)迭代運(yùn)算對進(jìn)行校正,直到取到穩(wěn)定的H值為止,計算方法見參考文獻(xiàn)[18]。
2.1.4日蒸散發(fā)
以上通量均為衛(wèi)星過境時刻的瞬時能量通量,而研究要用到每日蒸散發(fā),因此采用蒸發(fā)比不變法將瞬時蒸散發(fā)量轉(zhuǎn)換為日蒸散發(fā)量,見式(12)—(14)。
(12)
(13)
λ=[2.501-0.002361×(Ts-273.15) ]×106
(14)
式中 ET24——日蒸散發(fā)量,mm/d;Λ——蒸發(fā)比,無量綱;Rn24——一天內(nèi)的凈輻射通量,W/m2;λ——水的汽化潛熱,J/kg,可以由與溫度的關(guān)系式得到。
P-M公式是根據(jù)能量平衡與水分輸送原理研究非飽和下墊面蒸散發(fā)的方法[19],被認(rèn)為是世界上應(yīng)用最為廣泛、精度較高的公式之一。1948年P(guān)enman提出了在無水汽平流輸送的情況下,濕潤下墊面計算蒸發(fā)的公式[20]。之后,Monteith在Penman研究的基礎(chǔ)上將公式進(jìn)一步完善[21]。1998年,經(jīng)聯(lián)合國糧農(nóng)組織(FAO)修改驗證,發(fā)展為廣泛使用的FAO P-M公式,見式(15)—(16)。
(15)
ET=ET0×Kc
(16)
式中 ET0——參考作物蒸散發(fā)量,mm/d;γ——干濕表常數(shù),kPa/℃;Δ——飽和水汽壓曲線斜率,kPa/℃;G——土壤熱通量,W/m2;Rn——凈輻射通量,W/m2;u2——2 m處的風(fēng)速,m/s;T——日平均氣溫,℃;ea——實際水汽壓,kPa;es——飽和水汽壓,kPa;ET——實際蒸散發(fā)量;Kc——作物系數(shù)。
圖2所示,為評價SEBAL模型反演ET的精度,選取決定系數(shù)(R2),均方根誤差(RMSE)和相對誤差(MRE)這3個指標(biāo)對模型反演結(jié)果進(jìn)行分析,計算過程見式(17)—(19)。對于一個模型而言,評價指標(biāo)R2值越高,RMSE和MRE值越小,證明模型越可靠。
圖2 SEBAL模型計算流程
(17)
(18)
(19)
式中n——總樣品數(shù)量;ETPi——用FAO P-M公式計算ET;ETSi——SEBAL模型估算ET。
采用世界糧食組織(FAO)推薦的P-M公式計算蒸散發(fā),與SEBAL模型反演結(jié)果進(jìn)行對比驗證。利用FAO P-M公式結(jié)合永濟(jì)灌域氣象站點觀測數(shù)據(jù),計算得到參考作物蒸散量ET0,將ET0與作物系數(shù)相乘得到實際蒸散發(fā)量,其中作物系數(shù)參考《內(nèi)蒙古河套灌區(qū)不同作物騰發(fā)量及作物系數(shù)的研究》[22]中的取值。將氣象站點的FAO P-M計算值與對應(yīng)氣象站經(jīng)緯度的遙感像元蒸散值進(jìn)行對比分析,畫出散點圖,見圖3。為減少衛(wèi)星漂移以及大氣狀況的影響,在提取反演的ET值時,以站點中心3×3范圍的像元平均值作為模型的估算值。
圖3 日蒸散發(fā)精度驗證
圖3所示,SEBAL模型估算結(jié)果與FAO P-M公式計算結(jié)果比較:R2為0.94,說明二者之間具有較好的相關(guān)性和變化趨勢;RMSE為0.43 mm/d,MRE為8.62%,表明模型估算值與FAO P-M計算值接近,SEBAL模型的估算結(jié)果可信?;谏鲜鲵炞C結(jié)果,可以說明SEBAL模型具有良好的反演精度,模型在永濟(jì)灌域的適應(yīng)性和穩(wěn)定性較好。
本文SEBAL模型反演誤差為8.62%,與劉馨井雨等[23]利用SEBAL模型估算義長灌域ET的誤差結(jié)果為5%左右,YANG等[24]采用MODIS數(shù)據(jù)反演河套灌區(qū)ET的誤差結(jié)果為5.6%相比略高。影響SEBAL模型反演精度的因素有:①SEBAL模型冷熱像元點的挑選是手動進(jìn)行的,具有主觀性,對ET估算的精度影響較大;②在計算感熱通量時,即使是薄薄的一層云也會導(dǎo)致較大的誤差。本文選取的研究區(qū)遙感影像有少量云層,可能是導(dǎo)致模型反演結(jié)果誤差有所增加的原因。
選取2019年4月21日、5月15日、6月8日、7月10日、8月27日、9月20日、10月6日永濟(jì)灌域作物生長季內(nèi)晴空或云量較少的Landsat影像,分析永濟(jì)灌域2019年生長季日蒸散估算結(jié)果,見圖4。由于云量等因素的影響,部分像素點的處理結(jié)果為空白。
a)04-21
反演結(jié)果可以看出:2019年永濟(jì)灌域4—10月的日蒸散量在0~6.17 mm,日最大蒸散量由6月8日的6.17 mm,逐漸降低到10月6日的2.61 mm;日最小蒸散量為0 mm。從空間分布上可以看出,蒸散量的最大值主要位于研究區(qū)的南部邊界,這一區(qū)域緊鄰黃河,水分充足,蒸發(fā)強(qiáng)烈;而中部的荒地和草地,由于植被稀疏,蒸散量要明顯低于其他地區(qū)。
圖5統(tǒng)計了SEBAL模型反演的日均蒸散量。從時間序列上來看,永濟(jì)灌域4—6月的日均蒸散量逐漸增加,分別達(dá)到3.30、3.24、4.38 mm/d;7月日均蒸散量達(dá)到最大,為4.56 mm/d;8—10月日均蒸散逐漸下降,但8、9月的日均蒸散量還保持在相對較高的水平,分別為3.43、3.39 mm/d;10月的日均蒸散量下降明顯,為1.87 mm/d,達(dá)到生長季內(nèi)最低值??偟膩碚f,研究區(qū)生長季日均蒸散量呈現(xiàn)先增加后減少的變化趨勢,最大值出現(xiàn)在7月,最小值出現(xiàn)在10月,生長季日均蒸散量在1.87~4.56 mm內(nèi)變化,存在明顯的時間分布差異。
圖5 永濟(jì)灌域日均蒸散量
4月份氣溫逐漸回暖,大部分作物開始播種,蒸散量逐漸增大;6、7月氣溫升高,降雨量和灌溉水量增多,植被生長旺盛,蒸散量增加速率較快,在7月達(dá)到最大;8月末、9月氣溫開始下降,降雨量也逐漸減少,蒸散量呈現(xiàn)下降趨勢;10月初,作物成熟進(jìn)入收割階段,蒸發(fā)量降到最小。因此,模型反演的蒸散量在時間分布上符合永濟(jì)灌域作物生長規(guī)律。
根據(jù)研究區(qū)土地利用和覆蓋特點將研究區(qū)地形分為耕地、草地、水體、荒地和城鄉(xiāng)用地5種類型,見圖6。不同土地類型由于下墊面理化特性的不同,如植被狀況,土壤含水量等,在蒸散發(fā)的表現(xiàn)上也會有差異。結(jié)合土地利用分布圖,在ArcGIS軟件中利用空間分析方法對永濟(jì)灌域2019年4月21日、5月15日、6月8日、7月10日、8月27日、9月20日、10月6日的蒸散發(fā)分布圖進(jìn)行分區(qū)統(tǒng)計,分別提取耕地、草地、水體、荒地和城鄉(xiāng)用地5種土地利用類型的蒸散量。
從表3和圖7中可以看出,各土地利用類型的日均蒸散量,在時間上都是先升高而后逐漸降低,呈單峰型分布。水體的日均蒸散量最大,其次是耕地,草地和城鄉(xiāng)用地比較接近,荒地的日均蒸散量最小,不同土地利用類型的蒸散量存在較大差異。
表3 各土地利用類型的日均蒸散量 單位:mm/d
分析各種土地利用類型的蒸散特性,可以得到以下幾點認(rèn)識。
a)耕地是永濟(jì)灌域面積最多的土地利用類型,日蒸散量由4月21日開始增加,在7月10日達(dá)到峰值,而后逐漸減少,在10月6日達(dá)到最小值。在整個生長季內(nèi),耕地的日均蒸散量分布與整個永濟(jì)灌域的日均蒸散量十分接近,表明耕地的蒸散耗水是永濟(jì)灌域作物生長季內(nèi)水分消耗的主要組成部分,很大程度上體現(xiàn)了研究區(qū)的蒸散特征。
b)水體的日均蒸散量在2.39~6.01 mm,整體變化幅度較大。由于水體與其他土地類型的下墊面不同,有充足的水分蒸發(fā),相當(dāng)于蒸散過程中的 “冷點”區(qū)域,將接收的凈輻射通量幾乎都轉(zhuǎn)化為潛熱通量,因此水體的蒸散量是所有土地利用類型中最大的。
c)城鄉(xiāng)用地的蒸散量與草地接近,城鄉(xiāng)用地蒸散量略高于草地0.3 mm/d。
d)荒地日均蒸散量在1.38~2.68 mm,變化幅度不大。由于研究區(qū)氣候干旱,降水稀少,加上荒地的植被覆蓋度低,嚴(yán)重限制了地表蒸散活動的進(jìn)行,因此荒地的蒸散量最低。
綜上所述,永濟(jì)灌域各土地利用類型的蒸散量大小依次為水體>耕地>城鄉(xiāng)用地>草地>荒地。土地利用類型決定了蒸散量的高低,其空間分布也決定了蒸散發(fā)在空間格局上的差異,而同時蒸散發(fā)也能夠很好地反映下墊面的物理特性,兩者之間有著緊密的聯(lián)系。
將SEBAL模型應(yīng)用于干旱區(qū)的永濟(jì)灌域,評估其日尺度ET估算的準(zhǔn)確性和適用性,分析研究區(qū)生長季蒸散發(fā)的時空變化特征,主要結(jié)論如下。
a)在日尺度上,SEBAL模型估算ET與FAO P-M公式計算結(jié)果一致,決定系數(shù)為0.94,均方根誤差RMSE為0.43 mm/d,相對誤差MRE為8.62%。說明可以利用SEBAL模型和Landsat數(shù)據(jù)對永濟(jì)灌域的ET分布進(jìn)行估計。
b)從蒸散量的時間特征來看,研究區(qū)日均蒸散量在4月開始逐漸增加,到7月日均蒸散量達(dá)到最大,為4.56 mm/d。8—10月逐漸減少,10月日均蒸散量最小,為1.87 mm/d,蒸散量變化總體呈現(xiàn)單峰型分布。
c)蒸散量受土地利用類型的影響很大,呈現(xiàn)出明顯的差異性,不同土地利用類型生長季蒸散量大小依次為水體>耕地>城鄉(xiāng)用地>草地>荒地。
在利用SEBAL模型進(jìn)行蒸散發(fā)反演時,所需要的地表溫度、地表反照率、植被覆蓋度等地表參數(shù)也是經(jīng)遙感反演得到的,并未對這些參數(shù)進(jìn)行精度評價。所以可能會存在誤差,進(jìn)而導(dǎo)致誤差傳遞使模型中各通量的反演誤差增大,從而使反演得到的蒸散發(fā)量誤差也增大。因此,在以后的遙感反演研究中要注重提高地表參數(shù)的精度,從源頭上控制蒸散發(fā)反演的精度。