李 雪，鐘仕全

(1.貴州省黔西南自治州氣象局，貴州 興義 562400;2.廣西壯族自治區(qū)氣象減災研究所/國家衛(wèi)星氣象中心遙感應用試驗基地，廣西 南寧 530022)

1 引言

2008年9月6日，我國環(huán)境與災害監(jiān)測預報小衛(wèi)星星座A、B 兩顆衛(wèi)星(簡稱HJ－1A/B 星)發(fā)射升空。該星具有重訪能力強、分辨率高等特點，在許多方面都得到了廣泛的應用。此外，HJ－1B 星搭載了一臺紅外相機(IRS)，其第4 通道光譜范圍為10.5～12.5 μm，該通道對熱特性敏感，可用來記錄地表的發(fā)熱特性，星下點空間分辨率為300 m，重訪周期為4 d，這為地表溫度的反演提供新的數(shù)據(jù)源［3］。

目前已經(jīng)發(fā)展了多種地表溫度的遙感反演方法，考慮到HJ－1B 熱紅外波段范圍與TM 非常相近，因此，本文借鑒TM的地表溫度反演算法來對環(huán)境衛(wèi)星熱紅外通道進行反演，并根據(jù)HJ－1B 熱紅外波段的光譜響應特性來修訂算法中的一些經(jīng)驗關(guān)系，最后利用MODIS 溫度產(chǎn)品(MOD11_L2)進行反演精度的驗證，以獲得適用于巖溶地貌區(qū)的HJ衛(wèi)星遙感數(shù)據(jù)反演地表溫度的算法。

2 研究區(qū)概況及數(shù)據(jù)預處理

2.1 研究區(qū)概況

巖溶區(qū)多分布于廣西、云南、貴州等省，本文根據(jù)廣西的地質(zhì)地貌特征，選擇廣西中部干旱農(nóng)業(yè)區(qū)中的典型巖溶地貌區(qū)［4］。該區(qū)為亞熱帶季風氣候，農(nóng)業(yè)人口眾多，占總?cè)丝?5.12%，以農(nóng)業(yè)耕作為主，耕地面積為601 620 hm2，其中旱地占耕地面積的46.57%。因此，干旱是制約該地區(qū)農(nóng)業(yè)和社會經(jīng)濟發(fā)展的最為重大的問題之一。

2.2 數(shù)據(jù)預處理

本文使用的遙感數(shù)據(jù)見表1，其中，環(huán)境減災衛(wèi)星的CCD、熱紅外遙感影像數(shù)據(jù)以及MODIS 遙感數(shù)據(jù)用于對地表溫度的反演，MOD11_L2 為MODIS 溫度產(chǎn)品數(shù)據(jù)，用于對反演結(jié)果精度的驗證。

表1 所用遙感影像數(shù)據(jù)一覽表

采用經(jīng)過精校正后的TM 影像對HJ－1B 數(shù)據(jù)進行幾何精校正，將幾何誤差控制在1 個像元內(nèi)，并對HJ－CCD 數(shù)據(jù)重采樣到300m;利用專業(yè)遙感軟件對MOD 02 和MOD11_L2 數(shù)據(jù)進行幾何校正，并重采樣到300 m，并根據(jù)MOD11_L2 數(shù)據(jù)頭文件提供的計算公式獲得地表溫度影像。最后利用矢量圖裁剪出研究區(qū)，以進行下一步的處理和應用。

2.2.1 傳感器接收到的輻射強度LλHJ－1B 熱紅外影像是用DN 值來表示的，根據(jù)影像所帶頭文件提供的信息，傳感器所接收到的輻射強度和DN值之間有以下關(guān)系:

?Rose，R．，“What is lesson － drawing?”Journal of Public Policy，1991，11(1)，pp．3 ～30．

式中，Lλ為熱紅外傳感器所接收到的輻射強度，DN 為影像的灰度值，b 和g 為常數(shù)。

2.2.2 亮溫T 亮溫是遙感器在衛(wèi)星高度所觀測到的熱輻射強度相對應的黑體的溫度，這一溫度包含有大氣和地表對熱輻射傳導的影響，因而不是真正意義上的地表溫度［5］。根據(jù)Plank 輻射函數(shù)，由上述得到的輻射強度值Lλ就可以反算出影像的亮溫，公式如下:

式中，T 為亮溫(K);Lλ為輻射強度值，由公式(1)獲得;h 為普朗克常量，h=6.626×10－34J·s－1;c 為光速，c=2.998 ×108m·s－1;k 為玻爾茲曼常量，k=1.38 ×10－23J·K－1;λ 為有效波長，參考劉三超［6］文獻，HJ－1B的熱紅外通道有效波長為11.511 μm。

3 地表溫度反演算法及經(jīng)驗公式的修訂

3.1 普適性單通道算法

其中，Lsensor是傳感器接收的輻射強度，ε 是地表比輻射率，參數(shù)γ 和δ 可由公式(5)計算得到，Ψ1、Ψ2和Ψ3為大氣函數(shù)，可由公式(6)計算得到。

其中，Tsensor是傳感器亮溫;λ 是波段的有效波長;c1=1.191 043 56×108Wμm4m－2sr－1，c2=14 387.685 μmK。大氣函數(shù)(Ψ1，Ψ2，Ψ3)由大氣中水汽含量(ω)的函數(shù)而獲得，針對HJ－1B 熱紅外波段，參考段四波等［9］文獻，可由以下公式計算:

3.2 覃志豪單窗算法

該算法是覃志豪等［5］(2001)根據(jù)地表熱輻射傳導方程，通過一系列假設(shè)，建立的適用于Landsat TM 第6 波段的反演方法。其表達式為:

其中，a 和b 為常數(shù)，針對HJ－1B 熱紅外波段，根據(jù)圖像的溫度變化范圍得到，根據(jù)文獻［9］可知，對于0～30℃，a=－60.896 9 和b =0.439 078，對于20～50℃，a =－68.330 1 和b =0.464 012;C 和D 為參數(shù)，由公式(8)計算得到;Ta為大氣平均作用溫度，針對HJ－1B 熱紅外波段，對于中緯度夏季大氣，得到Ta的估算方法［9］。

式中ε 是地表比輻射率;τ 是大氣透過率，研究表明，大氣透過率的變化主要取決于大氣水汽含量ω的動態(tài)變化，因此，可以根據(jù)大氣總水汽含量來估算，針對HJ－1B 熱紅外波段進行修訂的大氣透過率估算方程(公式(10));T0為近地表溫度，可利用研究區(qū)域的站點平均氣溫數(shù)據(jù)差值得到。

3.3 基于影像的Artis 反演算法

Artis 等人［8］(1982)認為，輻射亮溫僅僅代表了黑體的溫度，而自然界的大部分物體并非黑體，所以應該用比輻射率對其進行校正，獲得絕對表面溫度，計算如下:

式中，T(K)為傳感器的亮溫;λ 為有效波長;ε為地表比輻射率;h 為普朗克常量，其值為6.626×10－34J·s;c 為光速，值為2.998 ×108m·s－1;σ 為玻爾茲曼常量，值為1.38 ×10－23J·K－1。

4 反演參數(shù)的計算

由上述一系列公式，可以看出需要計算的參數(shù)有:地表比輻射率ε 和大氣水汽含量ω。

4.1 地表比輻射率ε

地表比輻射率主要取決于地表的物質(zhì)結(jié)構(gòu)和遙感器的波段區(qū)間。計算地表比輻射率的方法很多，本文采用覃志豪［10］等(2004)根據(jù)Sobrino 提出的NDVITEM方法結(jié)合研究區(qū)地表類型的多樣性，得出的不同地表類型的比輻射率計算方法。認為地球表面不同區(qū)域的地表結(jié)構(gòu)雖然很復雜，但從衛(wèi)星像元的尺度來看，可以大體視作由3 種類型構(gòu)成:水體、城鎮(zhèn)和自然表面。水體在熱波段范圍內(nèi)的比輻射率很高，接近于黑體，可以用εwater=0.995 來進行估計。對于自然表面εsurface和城市地表εbuilt－up，與植被構(gòu)成比例Pv的關(guān)系式如下:

式中，Pv為植被構(gòu)成比例，其估計方法為［11］:

式中，ρ3、ρ4分別為HJ－1B CCD 數(shù)據(jù)的第3 和第4 波段的反射率。一般情況下，如果沒有詳細的區(qū)域植被和土壤光譜或圖幅上沒有明顯的完全植被或裸土像元，則用NDVIv=0.70 和NDVIs=0.05來進行Pv的近似估計。根據(jù)上述的方法，將CCD遙感數(shù)據(jù)分為水體、自然表面和城鎮(zhèn)3 種類型，然后按公式(12)～(14)計算出地表比輻射率。

4.2 大氣水汽含量ω

在缺少實時大氣水汽含量數(shù)據(jù)時，可以通過同步MODIS 數(shù)據(jù)反演大氣水汽含量。根據(jù)GRIEND等［11］文獻，得到大氣水汽含量的計算公式:

式中，ω 為大氣水汽含量;Tω(19/2)為MODIS第19 波段和第2 波段表觀反射率之比;α 和β分別為公式參數(shù)，對于混合型地表，α =0.02，β =0.651，其他類型地表參考文獻［12］。

5 結(jié)果分析

5.1 反演結(jié)果分析

根據(jù)上述3 種算法反演得到研究區(qū)的地表溫度影像圖(圖1)，能明顯看出三者溫度變化趨勢是相同的，表現(xiàn)在直方圖上其走勢相似(圖2)。其中，研究區(qū)的東部和南部區(qū)域地表溫度較高，溫度較低區(qū)域主要分布在西北部的河池、環(huán)江和融水等縣市。

圖1 3 種算法反演的地表溫度影像圖(注:白色區(qū)域為云區(qū))

圖2 3 種算法反演的研究區(qū)地表溫度直方圖

5.2 反演結(jié)果驗證

為了定量地分析各算法的反演結(jié)果，由于實時地面測量溫度很難獲取，驗證3 種算法的精度選用MODIS 地表溫度產(chǎn)品進行相對驗證。選取同天成像時間為11 時30分覆蓋研究區(qū)的MOD11_L2 溫度產(chǎn)品作為驗證數(shù)據(jù)。在低溫、中溫和高溫區(qū)中，分別選取對應的森林、耕地和城市3 種地表類型樣區(qū)進行結(jié)果驗證。

表2 列出了各檢驗樣區(qū)對應的HJ－1B 反演溫度值和MODIS 地表溫度產(chǎn)品的溫度值。從表中可以看出，普適性單通道算法反演結(jié)果與MODIS 地表溫度產(chǎn)品最為接近，平均溫差為0.36 K(覃志豪單窗算法為2.39 K，基于影像的Artis 反演算法為4.69 K)，森林溫差最大，平均為0.42 K，耕地溫差最小，平均為0.26 K?；谟跋竦腁rtis 反演算法與MODIS溫度產(chǎn)品差異最大，這是因為該方法只考慮了地表比輻射率的影響而忽略了大氣的影響。覃志豪單窗算法所得結(jié)果與MODIS 溫度產(chǎn)品平均溫差為2.39 K，主要原因為覃志豪單窗算法是根據(jù)MODTRAN 軟件4 種標準大氣廓線(美國1976 標準大氣廓線、低緯度標準大氣廓線、中緯度夏季和中緯度冬季標準大氣廓線)進行參數(shù)擬合，并建立大氣平均溫度和近地層空氣溫度的經(jīng)驗關(guān)系式;本文采用中緯度夏季標準大氣廓線建立關(guān)系式，計算的大氣平均作用溫度比實際大氣平均溫度高，導致反演得到的溫度比實際溫度低。

表2 3 種算法的反演結(jié)果與MODIS 地表溫度產(chǎn)品的溫度值統(tǒng)計 (單位:k)

為了便于觀察，對每一個檢驗區(qū)分別建立3 種算法反演的地表溫度與MODIS 地表溫度產(chǎn)品之間的零截距線性方程，見圖3:

圖3 森林檢驗區(qū)3 種算法反演的溫度與MODIS 地表溫度產(chǎn)品之間的零截距線性關(guān)系

從方程式可以初步看到，對于低溫區(qū)，普適性單通道算法反演結(jié)果優(yōu)于覃志豪單窗算法和Artis反演算法。為了比較各算法的反演結(jié)果，同樣對中溫區(qū)(耕地)和高溫區(qū)(城市)進行對比，以便更全面地比較各算法的反演結(jié)果。得到如下線性方程式:

對于耕地檢驗區(qū)有:

對于城市檢驗區(qū)有:

從方程式(16)～(17)可以看到，普適性單通道算法反演得到的結(jié)果優(yōu)于其他算法，反演得到的地表溫度更接近于MODIS 地表溫度產(chǎn)品的溫度。同時，建立的線性方程顯示了溫度變化的一個規(guī)律，即森林＜耕地＜城市，這個結(jié)果與實際情況相符合。

6 結(jié)論

本文針對HJ－1B 衛(wèi)星熱紅外數(shù)據(jù)的特點，分別用修訂后的普適性單通道算法、覃志豪單窗算法和基于影像的Artis 反演算法反演典型巖溶地貌區(qū)的地表溫度，并與MODIS 溫度產(chǎn)品相比較，以驗證算法精度。結(jié)果表明，修訂后的普適性單通道算法反演結(jié)果優(yōu)于其它兩種算法，其與MODIS 溫度產(chǎn)品平均溫差相差0.36 K，保證誤差精度在1 K 之內(nèi)，說明該算法經(jīng)過修訂后適用于反演巖溶地貌區(qū)的地表溫度，這為利用HJ－1B 遙感數(shù)據(jù)開展巖溶地貌區(qū)的干旱監(jiān)測提供技術(shù)借鑒。但由于研究中大氣參數(shù)和地表比輻射率只能依靠估算方程計算，這對反演精度會有一定的影響，需要在今后的研究中進行下一步驗證。

［1］夏日元，朱遠峰，李兆林.廣西巖溶區(qū)農(nóng)業(yè)發(fā)展的資源及地質(zhì)環(huán)境特征［J］.廣西科學，1997，4(2):192－195.

［2］徐希孺，柳欽火，陳家宜.遙感陸面溫度［J］.北京大學學報(自然科學版)，1998，34(2):248－253.

［3］李艷芳，李小娟，孟丹.環(huán)境減災衛(wèi)星熱紅外數(shù)據(jù)的地表溫度反演及LST分布分析——以北京市城八區(qū)為例［J］.首都師范大學學報(自然科學版)，2010，31(3):70－75.

［4］周游游，蔣忠誠，韋珍蓮.廣西中部喀斯特干旱農(nóng)業(yè)區(qū)的干旱程度及干旱成因分析［J］.中國巖溶，2003，22(2):144－149.

［5］覃志豪，Zhang M H，Arnon K，et al.用陸地衛(wèi)星TM6 數(shù)據(jù)演算地表溫度的單窗算法［J］.地理學報，2001，56(4):456－466.

［6］劉三超，柳欽火，高懋芳，等.波譜響應函數(shù)和波寬對地表溫度反演的影響［J］.遙感信息，2007(05):3－6.

［7］Jiménez－Mu?oz，J C.，Sobrino，J A.A Generalized Single Channel Method for Retrieving Land Surface Temperature From Remote Sensing Data ［J］.Journal of Geophysical Research，2003，108(D22):4688.

［8］Artis，D A，Carnahan，W H.Survey of Emissivity Variability in Thermography of Urban Areas［J］Remote Sensing of Environment，1982，12(4):313－329.

［9］段四波，閻廣建，錢永剛，等.利用HJ－1B 模擬數(shù)據(jù)反演地表溫度的兩種單通道算法［J］.自然科學進展，2008，18(9):1001－1008.

［10］覃志豪，李文娟，徐斌，等.陸地衛(wèi)星TM6 波段范圍內(nèi)地表比輻射率的估計［J］.國土資源遙感，2004(03):28－32.

［11］Carlson T N，Ripley D A.On the relation between NDVI，fractional vegetation cover，and leaf area index［J］.Remote Sens.Environ.，1997，62:241－252.

［12］GRIEND A，OWE M.On the Relationship Between Thermal Emissivity and the Normalized Difference Vegetation Index for Nature Surfaces［J］.International Journal of Remote Sensing，1993，14(6):1119－1131.