基于TM影像的河套灌區區域蒸散發研究

李 亮，翟 進，史海濱

（1.水利部牧區水利科學研究所，呼和浩特010020；2.內蒙古自治區水利水電勘測設計院，呼和浩特010020；3.內蒙古農業大學 水利與土木建筑工程學院，呼和浩特010018）

目前區域蒸散研究具有相當完善的理論基礎和較成熟的模型算法。但是所有模型算法都不是通用和萬能的，必須針對研究區的實際情況，建立適合本研究區的具體蒸散模型，只有這樣才能使估算結果更加準確有效。蒸發散計算生態需水量是較為普遍的方法之一。由于下墊面因素、水文參數等空間的變異性和不均勻性，使得由蒸發散計算的生態需水在向大尺度的轉換過程中產生了很大的誤差，影響了計算結果的準確性。遙感技術的實時性、區域性，為監測大范圍陸面地表能量和水分狀況提供了方便。同時隨著遙感空間分辨率、時間分辨率和光譜分辨率的提高，利用遙感技術定量反演地表參數和地表通量，立足于地表能量平衡方程，進而推算陸面蒸散量已成為區域蒸散估算的發展方向。基于遙感技術的區域蒸散研究在水資源缺乏的西北干旱、半干旱地區得到了廣泛的應用，而內蒙古河套灌區在區域蒸散研究方面尚屬空白。

1 研究區概況

內蒙古河套灌區位于中國西部，北緯40°19′—41°18′，東經106°20′—109°19′，是全國三大灌區之一［1］。東西長270km，南北寬40～75km。灌區地形平坦，西南高，東北低，海拔1 007～1 050m，坡度0.125‰～0.2‰。灌區總土地面積約為1.12×106hm2，現有灌溉面積約5.74×105hm2，占總土地面積的51.2%左右［2］。年降水量136.8～213.5mm，年蒸發量1 993～2 372mm，年平均氣溫6～8℃，自東向西升高，平均相對濕度40%～50%。全年封凍期5～6個月，最大凍結深度1.0～1.3m。封凍期為每年11月下旬至翌年4月，無霜期135～150d，全年日照期3 100～3 300h［3－4］。

2 數據選擇及研究方法

Landsat TM 5影像的精度高，重復訪問周期為16d，本研究采用30m分辨率數據，可以滿足灌區尺度的應用，選取2005年7月21日河套灌區1排干溝至7排干溝的灌區中、上游區域為研究對象，1景TM 5影像可覆蓋。運用遙感技術，使用Landsat TM 5、IKONOS多種分辨率影像數據對灌區基本特征參數進行提取。運用Erdas 9.1軟件進行非監督分類提取耕地的種植結構，用SEBAL模型反演灌區蒸散量。

3 數據解譯

3.1 數據預處理

原始遙感數據在使用前應進行圖像校正，本研究所用的數據為已經經過輻射校正和幾何粗校正的1B產品，因此，只需進行幾何精校正，校正流程如圖1所示。

3.2 波段選擇

對影像各波段數據進行統計特征分析、主成份分析、相關分析的結果表明，TM影像的7個波段數據中，TM5的信息量最大，其次為TM1、TM3和TM4，但TM1、TM3的波段相關性大，數據疊加多，因此本文選擇TM5、3、4偽彩色波段組合，疊加結果如附圖14所示。

圖1 TM影像校正流程

3.3 分類結果

3.3.1 分類精度 在ERDAS IMAGINE 9.1軟件中，采用非監督分類法，將影像分為40類，根據經驗以及地面實際數據逐一進行判別，最后合并為7個大類。由于實測資料不足，本文采用IKONOS影像對分類結果進行精度評價。IKONOS影像全色波段精度為1m，它能夠很好地反映地表覆蓋，目視即可判別大部分地表覆蓋。結果顯示，除林地和荒地外，其他各個類別都達到了較好的分類精度，總體精度為71.27%，Kappa系數為0.681 7。水體因其顯著的光譜差異性而具有很高的識別精度。由此可見，影像的分類精度較高，所得數據可靠性好，分類結果如附圖14所示。

3.3.2 分類結果 各類別名稱及在研究區域所占面積及比例見表1及附圖15。由表1可見，耕地所占的比例為50.8%，與實地調查所得的51.2%十分接近。

表1 研究區遙感土地利用分布

4 蒸散發反演

4.1 基本原理

太陽輻射是地表能量交換的基礎，當輻射能量經過大氣衰減到達地表后，其能量主要被用于加熱空氣與土壤以及促進水分蒸發，SEBAL模型就是利用了地表能量平衡原理來計算蒸騰量，其表達式如下：

式中：λET——潛熱通量，其中λ為汽化潛熱；ET——蒸騰量；Rn——凈輻射通量；G——土壤熱通量；H——感熱通量；PH——用于植物光合作用的能量（其值很小可以忽略）［5－7］。SEBAL模型根據 Landsat 5數據以及相關氣象資料逐像元地計算出研究區地面反照率，植被指數，比輻射率和地表溫度資料，并依據反演參數逐步計算出衛星過境時刻的Rn，G，H值，求出瞬時ET值，最終通過計算蒸發比分的方法推求出時段的ET量。

4.2 參數推求

（1）大氣外光譜反射率rb。Landsat TM／ETM波段1～5和7的波長為0.45～2.35μm，接收的主要是地面物體反射的太陽輻射，因此可以計算地面物體在大氣外光譜反射率rb。

式中：Lb——地面物體在波段b處的大氣外光譜輻射亮度［W／（m2·μm·sr）］；d——日地天文單位距離；Eb——波段b處的大氣頂層太陽光譜照射度［W／（m2·μm·sr）］；θ——太陽天頂角。

計算每個波段的反射率，需將灰度值（DN）轉化為輻射量。

式中：Gain——增益［W／（m2·μm·sr）］；Bias——偏置［W／（m2·μm·sr）］；QCAL——經過定標和量子化的比輻射率，無量綱；Lmin——QCAL＝0（或1）時的波譜輻射率，Lmax——QCALmax時的波譜輻射率；QCALmax——新的比輻射率（Rescaled Radiance）范圍，對于所有的TM值，QCALmax＝255。Lmin和Lmax的值可以從Landsat技術手冊得到（Landsat 7Science Data Users Handbook，2002）。

（2）地面反照率α。各波段反射率計算大氣頂層反照率αtoair的公式為：

式中：wb——波段b的權重系數，各波段權重參見表2；rb——波段b的反射率。

表2 各波段權重

對大氣外反照率作簡單的大氣輻射校正，得到地面的反照率α：

式中：αtoair——大氣外反射率；αpath——考慮了大氣影響的程輻射。本文中原始影像已進行了大氣校正，故忽略了大氣對地表反射率計算的影響，取αpath＝0。晴空單向大氣透射率的值一般可以由經驗公式估算［8］：

式中：z——地面高程（m）。

（3）歸一化差值植被數NDVI。

式中：CH3，CH4——式（9）計算的波段3和4的反射率。NDVI主要在－1～＋1之間，水體、建筑物、沙地及裸地的值接近于0；對應于高覆蓋度植被NDVI的值越大。

（4）比輻射率ε。比輻射率是一個無量綱值，取值在0～1之間。假定研究目標對熱輻射是不透明的，取值為0［9－11］。SEBAL中采用經驗公式計算比輻射率ε：

式中：INDV＞0，否則假設ε為0。研究區地表比輻射率均值為0.915，植被好的區域比輻射率較高，一般在0.96以上；鹽堿地及荒地的比輻射率最低，基本在0.90以下。

（5）地面溫度Ts。Landsat TM／ETM波段6的波譜范圍是10.4～12.5μm，主要接收地面長波輻射，可以用來計算地面溫度。Stafan—Boltzman定律反映物體溫度與輻射之間的關系，單波段6的波譜范圍太窄，因此利用Plank公式計算地面物體的亮度溫度Ts：

式中：L6——地面物體在波段6處的大氣頂層光譜輻射亮度，K1和K2為計算常數，見表3。

表3 計算常數K1和K2

4.3 模型參數確定

（1）地表凈輻射Rn。由所有的入射能量減去出射能量來計算，如圖2所示：

圖2 地表輻射平衡示意圖

式中：Sin——入射短波輻射；Lin——入射長波輻射；Lout——出射長波輻射；（1－ε）Lin——經地表反射的入射長波輻射項；α——地表反射率。研究區凈輻射量（Rn）值集中在500～600W／m2，占研究區域面積的90%以上，均值為581.4W／m2。

（2）土壤熱通量G。土壤熱通量取決于地表特征和土壤含水率等因素，在本文中，通過對多種計算土壤熱通量經驗公式的比較，我們采用Bastiaanssen［12］提出的經驗公式來估算地表（包括植被覆蓋地區和裸地）的G：

G＝（Ts－273.16）（0.0038＋0.0074α）

式中：Ts——地表溫度（K）。

研究區土壤熱通量（G）集中在40～80W／m2之間的值占研究區域面積的90%以上，其分布和地表溫度大致相同，而與地表反射率相反，即地表反射率越高，相應土壤吸收熱量的能力就越弱。

（3）感熱通量H。感熱通量的計算，假定研究區內Ts為線性關系，通過在地表溫度分布圖上選擇“冷點”與“熱點”，采用 Monin－Obukor迭代方法；通過對摩擦風速u＊和空氣動力學阻力rah經過多次循環遞歸最終求出穩定的H 值，計算公式如下：

式中：H——感熱通量（W／m2）；ρair——空氣密度；Cpair＝1004［J／（kg·K）］；rah——空氣動力學阻力（s／m）。

式中：P——大氣壓。H、aTS＋b和rah都是未知量，且彼此直接相關。計算Monin－Obukov長度時還需用到感熱通量H，因此只能進行迭代求解確定感熱通量H。

通過計算得到a＝0.499，b＝－143.66，rah＝21.15s／m。計算得到感熱通量分布在0～599.8W／m2，均值為282.81W／m2。低值主要分布在水域，高值主要分布在沙地，植被條件好的感熱通量為200～300W／m2，占區域面積的50%左右。

（4）潛熱通量。潛熱通量是下墊面與大氣之間交換的水汽通量，是水分循環和能量平衡的重要組成部分，潛熱通量與顯熱通量正好相反。研究區潛熱通量在0.022～704.77W／m2，均值為212.12W／m2。其分布規律與感熱通量相反。遙感反演區域潛熱通量見圖3—4。

（5）時段蒸騰量ET。根據上式計算出的能量平衡方程的各項結果，即可求得潛熱通量λET。但這一結果僅為衛片拍攝時的瞬時蒸散值，可以通過蒸發比在一天之中為常數的特性，通過計算蒸發比率Λ將瞬時的蒸散值ET延伸為全天蒸散值ET24。計算公式如下：

式中：ETinst——區域瞬時蒸散值；cos（s）——地表坡度余弦；λET——潛熱通量；ET——蒸騰量；λ——汽化潛熱；Λ——蒸發比率；G——土壤熱通量；G24——全天土壤熱通量；Rn——凈輻射通量；Rn24——全天太陽凈輻射。

圖3 遙感反演區域潛熱通量Histogram圖

圖4 遙感反演研究區蒸散量Histogram圖

4.4 區域蒸散量反演

利用蒸發比率推算出的日蒸散量結果如表4所示，根據區域內5個試驗點田間微氣象站及蒸散量反演值，日蒸散量均值為4.81mm／d，實測均值為5.09 mm／d，相對誤差平均為5.8%，結果合理可靠。由于外界熱量平流輸入會破壞自我穩定狀態，使能量構成比例發生變化，因而在使用瞬時蒸發比計算日蒸散量時存在一定誤差。對于更長時間的蒸散量計算，需要長序列的遙感影像，鑒于本地區遙感影像數據的不足，本文暫未考慮。

表4 反演值與實測值對比

5 結論

通過遙感數據分析，確定了TM影像的最佳波段組合，經過幾何校正采用非監督分類對影像進行分類，最終得到灌區土地分類及耕地種植結構。利用高分辨率IKONOS影像對分類結果作了精度評估，精度分析證明多時相中分辨率TM影像用于灌區尺度土地利用分類有很高的精度。

利用基于地表能量平衡原理的SEBAL模型，根據Landsat TM 5數據以及相關氣象資料逐像元地計算出研究區地面反照率，植被指數，比輻射率和地表溫度資料，并依據反演參數逐步計算出衛星過境時刻的Rn，G，H值，求出瞬時ET值，最終通過計算蒸發比分的方法推求出時段的ET量。經過反演得到日蒸散量均值為4.81mm／d，實測均值為5.09mm／d，相對誤差平均為5.8%，結果合理。利用地表能量平衡原理的SEBAL模型對蒸散發進行反演，反演結果與實測值誤差在允許范圍內，為河套灌區區域用水量研究提供新的方法。

［1］ 李亮.史海濱，賈錦鳳，等.內蒙古河套灌區荒地水鹽運移規律模擬［J］.農業工程學報，2010，26（1）：31－35.

［2］ 李亮.內蒙古河套灌區秋澆荒地水鹽運移規律的研究［J］.中國農村水利水電；2012（4）：41－44.

［3］ 李瑞平，史海濱，赤江剛夫，等.季節性凍融土壤水鹽動態預測BP網絡模型研究［J］.農業工程學報，2007，23（11）：125－128.

［4］ Shi Haibin，Akae Takeo.Simulation of leaching requirement for Hetao irrigation district considering salt redistribution after irrigation［J］.Transactions of the CSAE，2002，18（5）：67－72.

［5］ 李亮，史海濱，赤江剛夫，等.內蒙古河套灌區耕地與荒地間水鹽補排規律的研究［J］.灌溉排水學報，2010，29（5）：73－77.

［6］ Jackson R D，Idso S B，Reginato R J，et al.Canopy temperature as a crop water stress indication［J］.Water Resour.Res.，1981，17（4）：1133－1138.

［7］ Medina J L，Camacho E，Reca J，et al.Determination and analysis of regional evapotranspiration in southern Spain based on remote sensing and GIS［J］.Phys.Chem.Earth，1998，23（4）：427－432.

［8］ Moran M S，Clarke T R，Inoue Y，et al.Estimating crop water deficit using the relation between surface－air temperature and spectral vegetation index［J］.Remote Sens.Environ.，1994，49（3）：246－263.

［9］ 黃妙芬，劉紹民.地表溫度和地表輻射溫度差值分析［J］.地球科學進展，2005，10，20（10）：1075－1081.

［10］ 藺文靜，董華.基于SEBAL模型的區域蒸發蒸騰遙感估算［J］.遙感信息，2008（5）：50－54.

［11］ 曾麗紅.宋開山，張柏，等.應用Landsat數據和SEBAL模型反演區域蒸散發及其參數估算［J］.遙感技術與應用，2008，23（3）：255－263.

［12］ Bastiaanssen W G M.SEBAL－based sensible and latent heat fluxes in the irrigated Gediz Basin，Turkey［J］.Journal of Hydrology，2000，229（1／2）：87－100.