基于HJ－1B 衛(wèi)星數(shù)據(jù)的地表溫度反演方法研究

李 雪，鐘仕全

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

1 引言

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

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

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

2.1 研究區(qū)概況

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

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

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

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

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

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

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

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

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

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

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

3.1 普適性單通道算法

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

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

3.2 覃志豪單窗算法

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

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

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

3.3 基于影像的Artis 反演算法

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

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

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

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

4.1 地表比輻射率ε

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

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

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

4.2 大氣水汽含量ω

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

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

5 結(jié)果分析

5.1 反演結(jié)果分析

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

圖1 3 種算法反演的地表溫度影像圖(注:白色區(qū)域?yàn)樵茀^(qū))

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

5.2 反演結(jié)果驗(yàn)證

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

表2 列出了各檢驗(yàn)樣區(qū)對(duì)應(yīng)的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)品差異最大，這是因?yàn)樵摲椒ㄖ豢紤]了地表比輻射率的影響而忽略了大氣的影響。覃志豪單窗算法所得結(jié)果與MODIS 溫度產(chǎn)品平均溫差為2.39 K，主要原因?yàn)轳竞绬未八惴ㄊ歉鶕?jù)MODTRAN 軟件4 種標(biāo)準(zhǔn)大氣廓線(xiàn)(美國(guó)1976 標(biāo)準(zhǔn)大氣廓線(xiàn)、低緯度標(biāo)準(zhǔn)大氣廓線(xiàn)、中緯度夏季和中緯度冬季標(biāo)準(zhǔn)大氣廓線(xiàn))進(jìn)行參數(shù)擬合，并建立大氣平均溫度和近地層空氣溫度的經(jīng)驗(yàn)關(guān)系式;本文采用中緯度夏季標(biāo)準(zhǔn)大氣廓線(xiàn)建立關(guān)系式，計(jì)算的大氣平均作用溫度比實(shí)際大氣平均溫度高，導(dǎo)致反演得到的溫度比實(shí)際溫度低。

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

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

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

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

對(duì)于耕地檢驗(yàn)區(qū)有:

對(duì)于城市檢驗(yàn)區(qū)有:

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

6 結(jié)論

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

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

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

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

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

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

［6］劉三超，柳欽火，高懋芳，等.波譜響應(yīng)函數(shù)和波寬對(duì)地表溫度反演的影響［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］段四波，閻廣建，錢(qián)永剛，等.利用HJ－1B 模擬數(shù)據(jù)反演地表溫度的兩種單通道算法［J］.自然科學(xué)進(jìn)展，2008，18(9):1001－1008.

［10］覃志豪，李文娟，徐斌，等.陸地衛(wèi)星TM6 波段范圍內(nèi)地表比輻射率的估計(jì)［J］.國(guó)土資源遙感，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.