夏安全+齊建國+姜振飛+馬津

摘要：當前地表溫度反演的遙感數據源多為Landsat TM/ETM+、MODIS數據，Landsat 8熱紅外數據的使用還不是很多，許多針對Landsat 8熱紅外數據的地表溫度反演算法雖然被提出，但是否能滿足不同試驗區(qū)的精度要求還有待考究，同時，經驗模型的使用可能會對求解地表比輻射率和大氣透過率等地表溫度反演參數造成不同程度的影響。因此，本研究以山東省濟南市為研究區(qū)，Landsat 8數據為數據源，分別利用大氣校正法、JM_SC10算法、TIRS10_SC算法，結合分類回歸樹（classification and regression tree，簡稱CART）算法與中等光譜分辨率大氣透過率算法計算機模型（moderate spectral resolution atmospheric transmittance algorithm and computer model，簡稱MODTRAN），構建適合濟南地區(qū)的溫度反演參數，實現(xiàn)地表溫度的反演，并以濟南市16個氣象站的溫度數據為基準進行精度驗證。結果表明，3種算法反演溫度平均誤差為1.78 ℃，TIRS10_SC算法反演精度最高，其次是大氣校正法、JM_SC10算法。在一定誤差要求下，3種算法均可應用于濟南地區(qū)的地表溫度反演。

關鍵詞：Landsat 8；地表溫度反演；單通道算法；大氣透過率；地表比輻射率

中圖分類號： S127 文獻標志碼： A 文章編號：1002-1302（2017）20-0254-05

陸地表面溫度（land surface temperature，簡稱LST）是地表能量平衡和資源環(huán)境變化的重要參數。地面監(jiān)測站獲取溫度信息準確，但不能獲取大范圍溫度值及其時空分布，熱紅外遙感的出現(xiàn)解決了這個難題。熱紅外遙感技術能獲取熱紅外波段的輻射能量，根據地表物體的發(fā)射率特性反演得到其溫度，實現(xiàn)大范圍的溫度信息獲取，因而熱紅外遙感在地表溫度反演方面有著重要的作用[1]。

目前地表溫度的反演算法主要有單通道算法、劈窗算法、多通道多角度算法。單通道算法是指只利用1個熱紅外通道反演地表溫度的方法，適用于幾乎所有的熱紅外波段；劈窗算法是指利用相鄰的2個熱紅外通道來進行地表溫度反演的方法，是目前為止發(fā)展最為成熟的地表溫度反演算法；多通道算法還在發(fā)展之中，目前還沒有一個簡便可行的多通道算法可以用來進行地表溫度反演。Landsat 8有2個傳感器：陸地成像儀（operational land imager，簡稱OLI）和熱紅外傳感器（thermal infrared sensor，簡稱TRIS），其中TRIS有2個熱紅外波段，這雖然為劈窗算法的使用提供了條件，但一直以來美國地質調查局（United States Geological Survey，簡稱USGS）對TIRS11波段的定標準確性把握不準，所以不鼓勵使用劈窗算法來反演Landsat 8數據的LST[2]。因而本研究基于Landsat 8第10波段使用3種單通道算法，包括大氣校正法[3]、JM_SC10 算法[4]、TIRS10_SC算法[5]進行LST反演。

本研究以山東省濟南市為研究區(qū)，Landsat 8數據為數據源，根據濟南市實時的大氣水汽含量，選用合適的大氣透過率模型進行大氣透過率的估算；計算歸一化植被指數（normalized difference vegetation index，簡稱NDVI）、改進型歸一化水體指數（modified normalized difference water index，簡稱MNDWI）、壓縮數據維建筑用地指數（index-based build-up index，簡稱IBI），采用基于分類回歸樹（classification and regression tree，簡稱CART）算法的數據挖掘技術，獲取用于區(qū)分影像中不同類別地物的規(guī)則閾值，實現(xiàn)影像地物的精確分類；獲取每一類別純凈植被像元與純凈裸土像元的NDVI，摒棄植被覆蓋度計算的經驗模型，實現(xiàn)研究區(qū)植被覆蓋度及地表比輻射率的精確計算。在此基礎上，對3種算法反演的地表溫度進行比較分析，確定精度最高的單通道地表溫度反演算法，為后期應用Landsat 8數據進行地表溫度的反演提供參考。

1 材料與方法

1.1 研究區(qū)概況

濟南市為山東省省會，位于山東省中西部，介于36°01′～37°32′N、116°11′～117°44′E。地處中緯度地帶，季風明顯，四季分明，年平均氣溫為13.8 ℃，夏季平均氣溫為26.7 ℃，冬季平均氣溫在1.0 ℃左右。濟南三面環(huán)山，南依泰山，北跨黃河，地勢南高北低，落差達500 m，平均海拔高度為118 m，這種地勢構造令水汽和熱空氣回流聚集且不易擴散，一旦發(fā)生熱污染將很難消除，因此實時監(jiān)測大范圍的濟南市地表溫度是很有必要的。

1.2 數據處理

遙感數據選用2015年4月25號覆蓋研究區(qū)的Landsat 8 OLI、TIRS數據以及MODIS L1B 1KM數據。其中OLI數據用來計算研究區(qū)的地表比輻射率；TIRS數據用來反演LST；MODIS L1B 1KM數據用來計算研究區(qū)大氣水汽含量；地面實測數據由濟南市氣象站提供。

數據處理包括MODIS L1B 1KM數據、Landsat 8多光譜數據輻射定標與大氣校正、Landsat 8熱紅外數據輻射定標等，總體技術路線如圖1所示。

2 地表溫度反演算法

2.1 大氣校正法

衛(wèi)星傳感器接收到的熱紅外輻射亮度值Lλ由3部分組成：大氣向上輻射亮度L↑、地面的真實輻射亮度經過大氣層之后到達衛(wèi)星傳感器的能量、大氣向下輻射能量L↓到達地面后反射的能量。衛(wèi)星傳感器接收到的熱紅外輻射亮度值Lλ的表達式（輻射傳輸方程）可寫為：endprint

式中：λ為波長，μm；ελ為地表比輻射率；Ts為地表真實溫度，K；B（Ts）為黑體熱輻射亮度，W/（m2·sr·μm）；τλ為大氣在波長λ處的透過率。變換（1）式得溫度為Ts的黑體在熱紅外波段的輻射亮度B（Ts）為：

Ts可以用普朗克公式的函數獲?。?/p>

對于TIRS 10，K1=774.89 W/（m2·sr·μm），K2=1 321.08 K，將式（3）中B（Ts）替換為Lλ即為亮度溫度Tλ的計算表達式。

2.2 Jiménez-Muoz改進單通道算法

JM_SC10算法是由Jiménez-Muoz等于2014年在其原有單通道算法（single channel method，簡稱SC）的基礎上提出來的，增加了針對Landsat 8的大氣參數，JM_SC10算法如下：

式中：γ、δ是基于Planck函數的2個參數；bγ為常數，在TIRS10波段為1 324，Lλ表示熱輻射亮度值，W/（m2·sr·μm）；Tλ表示亮度溫度，K。ψ1、ψ2、ψ3與大氣水汽含量（g/cm2）ω有關，當ω>3 g/cm2時，

2.3 TIRS10_SC算法

胡德勇等于2015年提出了TIRS10_SC算法[5]，該算法是在總結輻射傳輸方程和覃志豪單窗算法的基礎上，專門針對Landsat 8 TIRS傳感器開發(fā)的。

T10為TIRS10的亮溫，K；K2=1 321.08 K；Ta為大氣平均作用溫度，K；根據影像的獲取時間，采用中緯度夏季的大氣平均作用溫度估算方程進行Ta的計算，估算方程如下：

式中：T0為近地表溫度（K），可根據當地氣象資料獲取。

3 大氣透過率τ與地表比輻射率ε計算

3.1 大氣透過率τ的計算

熱紅外波段大氣透過率主要取決于大氣水汽含量，可通過中等光譜分辨率大氣透過率算法計算機模型（moderate spectral resolution atmospheric transmittance algorithm and computer model，簡稱MODTRAN）模擬計算出大氣透過率與大氣水汽含量的關系，根據獲取的水汽含量信息來計算各像元大氣透過率。Landsat 8數據很難進行大氣水汽含量反演，但是MODIS 數據卻可以。MODIS數據包含36個波段，其中17、18、19為大氣吸收波段，第2、第5波段為大氣窗口波段。毛克彪等研究發(fā)現(xiàn)，利用MODIS數據的第2、19波段可反演大氣水汽含量[6]，大氣水汽含量ω為：

式中：ρ19和ρ2 分別為MODIS數據的第19、第2波段的地表反射率；α、β為常數，分別為0.020 0、0.632 1。

不同研究區(qū)、不同水汽含量下，大氣透過率的估算模型是不同的。當前針對Landsat 8數據的通用估算模型要求水汽含量在0.5～3.0 g/cm2，而研究區(qū)平均的水汽含量為 4.086 g/cm2，因此本研究區(qū)不能使用Landsat 8數據通用模型。鑒于Landsat 8數據的TIRS10與MODIS L1B數據的31波段中心波長和波寬相似（圖2），那么應用于31波段的大氣透過率估算方程，基本也適用于Landsat 8數據的TIRS10。根據中緯度大氣剖面數據進行模擬[7]，得到Band31在各種水汽含量下的通用大氣透過率估算方程（相關系數為 0.997 48）

通過式（12）實現(xiàn)TIRS10大氣透過率的計算，并得到該波段大氣透過率影像（圖3）。

3.2 地表比輻射率ε計算

地表比輻射率的計算一直是地表溫度反演的重點和難點。地表比輻射率別稱地表發(fā)射率，與地表的物質結構有關。要實現(xiàn)對地表比輻射率的精確計算，前提就是對影像地物進行精確分類?；贑ART算法的數據挖掘技術可以獲取不同地物的規(guī)則閾值，大大提高決策樹分類效率和精度。本研究將影像地物分為4類，分別是水體、自然表面（耕地、林地、草地）、城鎮(zhèn)表面（建筑物、道路）、裸土，利用CART算法獲取的各類規(guī)則閾值，確定決策樹如圖4所示，分類結果影像如圖5所示。

根據ASTER光譜庫和Nichol的研究成果[8]，獲取各地類純凈像元在TIRS10波段的地表比輻射率，水體為0.996 83、建筑為0.964 885、植被為0.986 72、裸土為0.967 67。影像中雖然存在各地類的純凈像元，但也存在大量的混合像元，其中自然表面和裸土區(qū)域可以看作是植被與裸土的混合，城鎮(zhèn)表面可以看作建筑物與植被的混合，考慮到混合像元這種情況，覃志豪等提出以下模型[9]計算地表比輻射率：

式中：ε自、ε城、ε裸分別為自然地表比輻射率、城鎮(zhèn)地表比輻射率以及裸土地表比輻射率；RV、RB、RS分別為植被、建筑、裸土的溫度比率；εV、εS、εB分別為純凈植被、裸土、建筑的地表輻射率；鑒于研究區(qū)南北高差相距500 m，可依據植被的構成比例簡單估計dε，經驗模型[10]為：

當0≤PV≤0.5時，dε=0.003 796PV；

當0.5< PV≤1時，dε=0.003 796（1- PV）

需要注意的是，如果應用公式計算出的ε自大于εV，取ε自=εV；ε城>εB，取ε城=εB；ε裸>εS，取ε裸=εS。PV為植被覆蓋度，由改進像元二分模型計算：

NDVIS為純凈裸土或者建筑像元的NDVI值；NDVIV為純凈植被像元的NDVI值。溫度比率Ri與植被覆蓋度PV相關，覃志豪等根據各地表類型的溫度差異進行模擬，確定出植被、建筑、裸土表面的溫度比率[9]：

計算地表比輻射率需要參數多、運算復雜，總體波段運算表達式如下：

式中：b1、b6、b7、b9分別為植被、裸地、建筑、水體的掩膜影像；b2為植被覆蓋度影像；b3、b4、b8分別為RV、RS、RB影像；b5為dε影像。地表比輻射率影像如圖6所示。endprint

4 地表溫度反演

4.1 地表溫度反演結果

逐一計算各個參數之后，進行3種算法的地表溫度反演，反演結果（取置信區(qū)間99.9%以去除異常值）如圖7～9所示。

4.2 精度驗證

根據濟南市氣象網統(tǒng)計，2015年4月25日濟南市溫度大致在16～31 ℃之間，部分地區(qū)溫度達到36～40 ℃，與影像反演結果大致相同。結合反演影像得出高溫天氣主要發(fā)生在章丘市以及濟南市區(qū)部分地區(qū)。究其原因是濟南市重工業(yè)多集中在該地區(qū)，如三一重工、濟南重工、梟龍重工等，重工業(yè)產生的熱污染使章丘市與濟南城區(qū)的部分地區(qū)溫度居高不下，嚴重影響到人們的生活。

對3種算法反演的溫度數據進行分類統(tǒng)計發(fā)現(xiàn)，耕地植被區(qū)域溫度相對適中，溫度在22～37 ℃之間，平均溫度為 29 ℃；水體區(qū)域溫度最低，溫度在16～33 ℃之間，其中黃河水域的溫度主要集中在17～21 ℃之間，各大水庫如鵲山水庫、臥虎山水庫、玉清湖水庫的溫度在16～19 ℃之間，溫度較高的水體區(qū)域多為零星分布的池塘以及海岸地帶；裸地和建筑區(qū)域溫度最高，在28～40 ℃之間，平均溫度為33 ℃。整體來看，3種算法對濟南市地表溫度的反演比較合理。

為精確驗證反演結果的精度，本研究選取濟南市氣象站16個站點溫度數據作為實測數據，對3種算法反演的地表溫度進行精度驗證。各站點實測溫度與應用3種算法反演的溫度數據對比結果如表1所示。

由表1可看出，3種單通道算法中TIRS10_SC算法反演地表溫度誤差最小，平均誤差為1.10 ℃，絕對平均誤差為 1.24 ℃；其次是大氣校正法，平均誤差為1.40 ℃，絕對平均誤差為1.45 ℃；最大誤差來自JM_SC10算法反演的地表溫度，平均誤差為2.85 ℃，絕對平均誤差為2.85 ℃。算法反演溫度與溫度站實測溫度不盡相同，其中大部分反演溫度略高于實測站溫度，分析原因主要有：（1）溫度站實測溫度的獲取時間為10：30，而影像的獲取時間為10：48，在晴朗天氣下，溫度是逐漸上升的。（2）對于算法中應用到的大氣向上和向下的輻射亮度L↑、L↓，采用美國航空航天局（National Aeronautics and Space Administration，簡稱NASA）提供的計算模型。要提高模型計算的精度，輸入參數必須要有實時的氣壓、風速、海拔高度、相對濕度等參數，本研究沒有這些數據，使用了模型的缺省值。（3）對JM_SC10算法，大氣水汽含量是影響該算法精度的一個關鍵參數，當大氣水汽含量降低到2 g/cm2時，地表溫度的反演誤差會降低1.5～3.0 ℃[4]。試驗中研究區(qū)平均大氣水汽含量為4.086 g/cm2，因此可以推斷，一旦大氣水汽含量降低到2 g/cm2，應用JM_SC10算法反演溫度的誤差有可能至少降低1.5 ℃， 使該算法反演的地表溫度誤差大大降低。

盡管上述原因造成了一定的誤差，但3種算法反演的溫度誤差平均值為1.78 ℃，絕對誤差的平均值為1.85 ℃，反演結果相對理想。

5 討論

本研究選用3種地表溫度反演算法進行地表溫度的反演，根據濟南地區(qū)的實際情況，構建適合濟南地區(qū)的大氣透過率和地表比輻射率參數。本研究創(chuàng)新點有：（1）對大氣透過率的計算，考慮到各研究區(qū)不同大氣狀況、不同水汽含量的影響，選用與Landsat 8數據第10波段性質相近的MODIS第31波段的大氣透過率估算模型，消除了含水量的影響，摒棄先前整景影像的單一大氣透過率計算，實現(xiàn)影像上每像元的大氣透過率求解。（2）基于CART算法的數據挖掘技術，結合多源影像，如NDVI影像、MNDWI影像、IBI影像，獲取用于區(qū)分影像中不同類別地物的規(guī)則閾值，實現(xiàn)影像地物的精確分類。在此基礎上，獲取每一類別純凈植被像元與純凈裸土像元的NDVI值，摒棄植被覆蓋度計算的經驗模型，求解溫度比率，實現(xiàn)研究區(qū)植被覆蓋度以及地表比輻射率的精確計算。

研究展望：（1）考慮如果能根據濟南地區(qū)的實際情況利用MODTRAN模型模擬出Landsat 8數據第10波段大氣水汽含量與大氣透過率的關系模型，對反演精度的提高會大有幫助；（2）獲取濟南市實時的氣壓、風速、海拔高度、相對濕度等參數，實現(xiàn)L↑、L↓精確計算；（3）JM_SC10算法提出的前提是ω>3 g/cm2，考慮對JM_SC10算法的使用進行改進，在原有算法的基礎上 ±1.5 ℃，反演2幅溫度影像，根據部分實測溫度數據，確定出真實的溫度反演影像，提高算法的精度。

6 結束語

不管是通用的大氣校正法還是最近提出的地表溫度反演算法，在實現(xiàn)濟南市大氣透過率和地表比輻射率的精確計算后，各算法均對濟南市地表溫度的反演有很好的適用性，尤其是TIRS10_SC算法精度相對較高，而當大氣水汽含量低于 2 g/cm2 時，JM_SC10算法精度也可大大提高。

參考文獻：

[1]宋 挺，段 崢，劉軍志，等. Landsat 8數據地表溫度反演算法對比[J]. 遙感學報，2015，19（3）：451-464.

[2]徐涵秋. 新型Landsat 8衛(wèi)星影像的反射率和地表溫度反演[J]. 地球物理學報，2015，58（3）：741-747.

[3]Qin Z H，Karnieli A，Berliner P. A mono-window algorithm for retrieving land surface temperature from Landsat TM data and its application to the Israel-Egypt border region[J]. lnternational Journal of Remote Sensing，2001，22（18）：3719-3746.

[4]Jiménez-Muoz J C，Sobrino J A，Skokovic D，et al. Land surface temperature retrieval methods from Landsat-8 thermal infrared sensor data[J]. IEEE Geoscience and Remote Sensing Letters，2014，11（10）：1840-1843.

[5]胡德勇，喬 琨，王興玲，等. 單窗算法結合Landsat 8熱紅外數據反演地表溫度[J]. 遙感學報，2015，19（6）：964-976.

[6]毛克彪，覃志豪，王建明，等. 針對MODIS數據的大氣水汽含量反演及31和32波段透過率計算[J]. 國土資源遙感，2005（1）：26-29.

[7]Mao K B，Qin Z H，Shi J. A practical split-window algorithm for retrieving land surface temperature from MODIS data[J]. International Journal of Remote Sensing，2005，26（15）：3181-3204.

[8]Nichol J. An emissivity modulation method for spatial enhancement of thermal satellite images in urban heat island analysis[J]. Photogrammetric Engineering and Remote Sensing，2009，75（5）：547-556.

[9]覃志豪，李文娟，徐 斌，等. 陸地衛(wèi)星TM 6波段范圍內地表比輻射率的估計[J]. 國土資源遙感， 2004，16（3）：27-32.

[10]Sobrino J A，Jiménez-Munoz J C，Paolini L. Land surface temperature retrieval from Landsat TM 5[J]. Remote Sensing of Environment，2004，90（4）：434-440.endprint