-
高空间分辨率的地表温度(land surface temperature,LST)数据对于城市热岛监测、土壤水分估算、地气水热交换等有着重要的作用[1-3]。目前主要搭载热红外波段传感器的遥感平台中,可见光/近红外波段空间分辨率都比热红外波段空间分辨率高[4],使得从热红外波段反演得到的地表温度空间分辨率也不高。因此,如何综合可见光/近红外波段的信息提高地表温度的空间分辨率一直是热红外遥感应用研究的重要方向[5-6]。
目前常用的分解方法是文献[7]提出的Distrad算法,该算法通过拟合地表温度与归一化植被指数(normalized difference vegetation index, NDVI)的线性统计关系LST-NDVI实现地表温度的降尺度。文献[8]基于LST-NDVI统计模型中干、湿边组成的三角形关系,提出了热红外波段像元分解的三角算法。TsHARP方法[9]提出用植被覆盖度来代替NDVI。但上述模型仅考虑了单一的植被指数,因此,其他不同的遥感指数也逐渐被引入。文献[10]分析了地表温度与归一化裸土指数的统计关系。文献[11]根据不同的地表类型,分别统计了地表温度与15种遥感指数之间的线性关系,提出了E-DisTrad算法,并分析得出在城市地表中不透水层覆盖度与地表温度的拟合关系最佳。上述模型利用地表温度与遥感指数的统计关系实现了地表温度像元分解,但是整个分解过程中没有考虑像元分解前后辐射能量值的一致性。
本文以北京市Landsat TM影像为数据源,首先基于地表类型的线性统计模型(E-DisTrad)获得高空间分辨率的初次分解地表温度,计算得到初次分解子像元的辐亮度。然后,采用面向对象的图像分割方法对初次分解子像元辐亮度进行4×4规则分割,得到初次分解子像元的辐亮度与其分割后每个对象的辐亮度之间的权重。最后,根据上述权重实现对地表温度父像元的二次分解。由于该方法涉及地表温度两次分解,因此,本文将此方法称之为地表温度二次像元分解方法。
-
本文使用Landsat TM数据,成像地区位于北京,过境时间是2009-07-20T10:42:26。对TM数据进行辐射定标、大气校正等预处理后,选取了北京市一个矩形区域作为本文的实验区,实验区TM假彩色合成影像(RGB432)如图 1所示。
-
Landsat TM可见光/近红外波段的空间分辨率为30 m,热红外波段的空间分辨率为120 m,反演得到的地表温度的空间分辨率也是120 m,从而极大地限制了地表温度的应用。本文的研究目标是以高空间分辨的可见光/近红外波段为辅助数据,提出一种新的地表温度像元分解算法,以期将地表温度的空间分辨率提高到30 m。本文提出的二次像元分解方法的原理如图 2所示。
利用120 m空间分辨率的地表温度数据Ts,对于任一个父像元,可以计算得到其分解前的辐亮度值(I)。以高空间分辨率的可见光/近红外波段数据为基础,以每个分解后子像元占其父像元辐射能量值的权重(Wi)为辅助数据,从而可以确定每个分解后子像元的辐亮度值(Ii),获得16个30 m空间分辨率的子像元(Tsi)地表温度。因此,如何获取子像元的权重是本文方法的重点。
本文首先利用E-DisTrad方法分别统计地表温度与不同地表类型的遥感指数之间的线性关系,从而获得高空间分辨率的初次分解地表温度,计算得到初次分解子像元的辐亮度值,为二次分解子像元的权重的计算提供基础。因此,如何获取初次分解地表温度是本文方法的第一步。
-
初次分解地表温度的计算是通过E-DisTrad方法进行的。水体的地表温度趋于稳定,故可以在分解时赋予一个合理的统计值,对植被和建筑用地两类地表类型,分别统计NDVI和不透水层(impervious surface area, ISA)覆盖度与地表温度的线性统计关系。
利用可见光/近红外波段计算得到30 m空间分辨率的归一化植被指数NDVI30和不透水层指数ISA30。由于从热红外波段反演得到的地表温度空间分辨率为120 m,因此,需要将NDVI30和ISA30重采样到120 m,从而分别统计地表温度LST120与归一化植被指数NDVI120和不透水层指数ISA120的线性关系,得:
$$ \begin{array}{l} {\rm{植被: }}{\rm{LS}}{{\rm{T}}_{120}} = a{\rm{NDV}}{{\rm{I}}_{120}} + b\\ 建筑用地:{\rm{LS}}{{\rm{T}}_{120}} = c{\rm{IS}}{{\rm{A}}_{120}} + d \end{array} $$ (1) 式中,a、b、c、d分别为地表温度与遥感指数之间的线性拟合参数。LST120采用单窗算法反演,ISA120采用线性光谱分解模型计算,NDVI120采用TM波段运算获得,上述参数的计算详见参考文献[12-14]。
E-DisTrad方法假设地表温度与遥感指数的拟合关系在不同空间分辨率上是一致的,因此,可以将上述线性关系直接应用于NDVI30和ISA30遥感指数,从而获得30 m空间分辨率的初次分解地表温度LST30,得:
$$ \begin{array}{l} 植被:{\rm{LS}}{{\rm{T}}_{30}} = a{\rm{NDV}}{{\rm{I}}_{30}} + b\\ 建筑用地:{\rm{LS}}{{\rm{T}}_{30}} = c{\rm{IS}}{{\rm{A}}_{30}} + d \end{array} $$ (2) 通过上述分解方法获取的初次分解地表温度并不能保证能量守恒,因此,以初次分解地表温度为基础,获取每个二次分解子像元占其父像元辐亮度的权重,从而得到每个二次分解子像元的辐亮度值,实现对地表温度的二次分解。
-
根据120 m空间分辨率的地表温度,对于任一个父像元,根据Stenfan-Boltzman定律,可以计算得到辐亮度值I:
$$ I = \varepsilon \sigma {T^4} $$ (3) 式中,I为父像元的辐射亮度值(W·m-2);ε为父像元的比辐射率;T为父像元的地表温度(K);σ为Stenfan-Boltzman常数,σ=5.67×10-8W·m-2·K-4。二次分解子像元权重计算方法见图 3。
如图 3所示,对于任一个30 m空间分辨率的初次分解子像元地表温度计算得到辐亮度值:
$$ {I_k} = {\varepsilon _k}\sigma T_k^4 $$ (4) 式中,Ik为第k个初次分解子像元的辐亮度值(W·m-2);εk为第k个初次分解子像元的比辐射率;Tk为第k个初次分解子像元的地表温度(K)。
然后,根据面向对象的图像分割方法,对初次分解子像元的辐亮度进行4×4规则分割,使分割得到的每个对象与其分解前父像元是相对应的。则每个分割后对象的辐亮度值为:
$$ I\prime = \sum\limits_{k = 1}^{16} {{I_k}{S_2}/{S_1}} $$ (5) 式中,I′为分割后对象的辐亮度值(W·m-2);S1和S2分别为分割后对象和初次分解子像元的面积(m2)。从而可以确定初次分解子像元的辐射能量与其分割后每个对象的辐射能量之间的权重,得到每个二次分解子像元的权重Wk:
$$ {W_k} = \frac{{{I_k}{S_2}}}{{I\prime {S_1}}} = {I_k}/\sum\limits_{k = 1}^{16} {{I_k}} $$ (6) 最后,利用上述权重值实现对其分解前120 m空间分辨率的父像元的二次分解,得到每个二次分解子像元的辐亮度值(I′k):
$$ {{I'}_k} = {W_k}I{S_1}/{S_2} = 16{W_k}I $$ (7) 因此,二次分解子像元的地表温温度Tsk为:
$$ {T_{{s_k}}} = {({{I'}_k}/\sigma {\varepsilon _k})^{0.25}} $$ (8) -
由于地表温度每天变化较快,很难找到实时的温度场来验证TM影像的分解结果,故本文借鉴文献[15]提出的“先升尺度再分解”验证方法,评价方法的优劣性。
首先,利用单窗算法反演得到120 m空间分辨率的地表温度,可以认为该地表温度能够反映出实验区真实的地表温度空间分布,即为“真值”。然后,采用均值法将上述地表温度重采样到480 m。最后,分别采用E-DisTrad方法和二次像元分解方法将480 m空间分辨率的地表温度数据降尺度到120 m,与“真值”进行比较分析,评价两种分解方法的差异性。
-
利用单窗算法反演得到的地表温度图像的空间分辨率是120 m,可以认为该图像代表真实的地表温度空间分布 Ts (见图 4(a))。图 4(b)、4(c)分别表示E-DisTrad和二次像元分解两种地表温度图像分解模型的结果Te和 Tw。
由图 4可知,实验区内建筑用地的地表温度最高,植被次之,水体的地表温度最低,3种算法计算得到的地表温度具有相似的整体空间格局。高温区(≥310 K)主要集中于北京市城区的南部,次高温区(305~310 K)主要集中于中部,低温区(≤305 K)主要集中于水体和植被覆盖度较高的郊区。这与北京中心城区快速的城市化进程相吻合。与Ts(图 4(a))相比,Te和Tw(图 4(b)、4(c))地表温度的空间分辨率都为30 m,相较于热红外波段的空间分辨率(120 m)提高了3倍。对于E-DisTrad模型,高温区和低温区的像元明显增多,地表温度趋于305~307 K的像元几乎没有,与分解前实验区地表温度的空间分布一致性较差。对于二次像元分解模型,地表温度的整体变化趋势更为平缓,但是,地表温度较低(302~305 K)的像元增多,多发现于郊区植被中。
-
为了验证二次像元分解方法的合理性,本文采用“先升尺度再分解”的方法获得高空间分辨率的地表温度数据(见图 5)。其中,图 5(a)是TM6波段采用单窗算法反演得到的120 m空间分辨率的地表温度(Ts),作为两种方法对比的真值;图 5(b)、5(c)分别是将120 m的热红外波段数据升尺度到480 m后,经过E-DisTrad方法和二次像元分解方法得到的120 m空间分辨率的地表温度结果(T′e和T′w)。
与Ts相比较,T′e地表温度的统计值明显偏高,中心区域的地表温度普遍高于320 K,多被发现在植被和城市用地的混合地区,地表温度趋于305~310 K的像元几乎没有。T′w地表温度的统计值与观测值较为接近,地表温度的空间分布与观测结果更为一致,地表温度的变化趋势更为平缓。这与地表温度的频数统计曲线(见图 6(b))相一致。
从图 6(b)中可以看出,二次像元分解方法结果T′w(紫色)的地表温度频数统计曲线与单窗算法反演结果Ts(红色)更为接近,不同地表温度区间的像元数目也较为接近。与单窗算法反演结果相比,在302~310 K地表温度区间内,E-DisTrad方法降尺度得到的地表温度像元数目明显少于TM6反演结果,而在314 K和300 K处(图 6(b)中A、B)地表温度像元频数出现两个峰值,这是由于表征建筑用地的ISA遥感指数在0.7~0.9区间内的像元频数最大(图 6(a)中A),从而导致城市地区地表温度初次分解结果多集中在314 K左右,高估了城市地表的温度。与此相同的是,表征自然表面的NDVI遥感指数在0.6~0.8区间内的像元频数也比较多(图 6(a)中B),使得初次分解温度在300 K处出现另一个峰值。这两个峰值的出现,使得302~310 K地表温度区间内的像元数目远远低于实际真值。但是,通过二次像元分解方法分解后,上述情况得到明显改善,二次像元分解方法得到的不同地表温度区间内像元频数与单窗算法反演结果更为接近。需要特别指出的是,在312~314 K地表温度区间内,二次像元分解方法结果的像元数目增多,多发现于建筑用地和林草地等的混合区域。图 7给出了T′e、T′w与Ts之间地表温度的散点图。从图 7中可以看出,T′w与Ts之间地表温度的散点图拟合线的斜率更接近1,R2值更大,更加印证了基于二次像元分解模型降尺度得到的地表温度与观测结果更为接近。
-
本文提出了一种能量平衡角度的地表温度二次像元分解方法,并将该方法应用到具有复杂地表覆盖类型的北京市Landsat TM遥感影像中,实现了地表温度的分解,将空间分辨率提高到30 m,相较于热红外波段的空间分辨率(120 m)提高了3倍。
同TM6波段直接利用单窗算法反演得到的地表温度相比,本文采用升尺度再降尺度的验证方法,讨论了E-DisTrad模型和二次像元分解模型分解后获得的地表温度的精度。分解后地表温度的空间分布、频数统计曲线和散点图分析表明,利用二次像元分解模型降尺度得到的地表温度不仅能准确地反映出实验区不同地表类型地表温度空间变化的差异性,而且有效地保证了分解前后影像的总体及局部能量的守恒,非常适合用于复杂地表覆盖地区的热红外波段影像数据的降尺度处理。
需要指出的是,本文中初次分解地表温度的计算以E-DisTrad模型为基础,其他不同的地表温度像元分解方法也可以进行引入实验;子像元权重计算过程中采用面向对象的4×4规则分割,但是对象的定义不单单是规则的矩形,其他分割方法仍需进行进一步的研究。
An Efficient Approach for Pixel Decomposition of Land Surface Temperature from Landsat TM Data
-
摘要: 提出了一种基于Landsat TM的地表温度二次像元分解方法,将地表温度的空间分辨率从120 m提高到30 m。首先,利用地表类型的线性统计模型(E-DisTrad)获取初次分解子像元的地表温度,计算得到初次分解子像元的辐亮度;然后,利用面向对象的图像分割方法获取二次分解子像元的权重,实现对地表温度的二次分解;最后,采用升尺度再分解的验证方法进行精度分析,并选取了北京市TM影像进行实例分析。实验结果表明,二次像元分解模型不仅能有效地提高地表温度的空间分辨率,反映出不同地表类型地表温度的空间差异性,而且保证了像元分解前后能量值的一致性,非常适合于复杂地表覆盖地区的热红外波段遥感影像数据的降尺度处理。
-
关键词:
- 地表温度 /
- 像元分解 /
- Landsat TM
Abstract: The spatial resolution of land surface temperature (LST, 120 m) retrieved from thermal infrared (TIR) band is lower than its visible/near-infrared (VNIR) bands (30 m). LST image with high spatial resolution compatible with VNIR bands of Landsat TM is very important for application of the LST image to many studies such as environmental monitoring. The objective of this study is to decompose the coarse pixels of the LST image data into the same pixel scale of the auxiliary VNIR data. Firstly, the E-DisTrad method is used to divide LST of the parent pixels into the sub-pixels to obtain the first decomposition temperature and the theoretical radiance at-sensor of each sub-pixel. Then, a chess-segmentation is done to the temperature of the sub-pixels on the basis of the object-oriented image segmentation method to compute the weight for each sub-pixel, which is consequently used to allocate the thermal radiance for each sub-pixel to generate the decomposed LST image. At last, a method of increase the spatial resolution firstly and then used for double-step pixel decomposition is executed in the study to validate the accuracy and efficiency of the decomposition method for the Landsat TM image of Beijing. The result showed that the double-step pixel decomposition method is applicable for decomposition of LST images for high spatial resolution, reflecting the spatial differences of different land use and land cover and ensure the conservation of energy before and after pixel decomposition. Therefore we can conclude that it is ideally suited for complex covered surface area of downscaling thermal infrared band remote sensing data.-
Key words:
- land surface temperature /
- pixel decomposition /
- Landsat TM
-
[1] Lu D, Weng Q. SpectralMixture Analysis of ASTER Images for Examining the Relationship Between Urban Thermal Features and Biophysical Descriptors in Indianapolis, Indiana, USA[J]. Remote Sensing of Environment, 2006, 104(2):157-167 doi: 10.1016/j.rse.2005.11.015 [2] Li Zhaoliang, Tang Bohui. Wu Hua, et al. Satellite-Derived Land Surface Temperature:Current Status and Perspective[J]. Remote Sensing of Environment, 2013, 131:14-37 doi: 10.1016/j.rse.2012.12.008 [3] 曹丽琴, 张良培, 李平湘, 等.城市下垫面覆盖类型变化对热岛效应影响的模拟研究[J].武汉大学学报·信息科学版, 2008, 33(12):1229-1232 Simulation Study of Influence of Change of Land Surface Types on Urban Heat Island[J]. Geomatics and Information Science of Wuhan University, 2008, 33(12):1229-1232 [4] Agam N, Kustas W P, Anderson M C, et al. A Vegetation Index Based Technique for Spatial Sharpening of Thermal Imagery[J]. Remote Sensing of Environment, 2007, 107(4):545-558 doi: 10.1016/j.rse.2006.10.006 [5] 杨贵军, 柳钦火, 刘强, 等.基于遗传自组织神经元网络的可见光与热红外遥感数据融合方法[J].武汉大学学报·信息科学版, 2007, 32(9):786-790 http://ch.whu.edu.cn/CN/abstract/abstract1986.shtml Yang Guijun, Liu Qinhuo, Liu Qiang, et al. Fusion of Visible and Thermal Infrared Remote Sensing Data Based on GA-SOFM Netural Network[J]. Geomatics and Information Science of Wuhan University, 2007, 32(9):786-790 http://ch.whu.edu.cn/CN/abstract/abstract1986.shtml [6] Friedl M. Forward and Inverse Modeling of Land Surface Energy Balance Using Surface Temperature Measurements[J].Remote Sensing of Environment, 2002, 79(2):344-354 http://citeseerx.ist.psu.edu/showciting?cid=5811350 [7] Kustas W P, Norman J M, Anderson M C, et al. Estimating Subpixel Surface Temperatures and Energy Fluxes from the Vegetation Index-Radiometric Temperature Relationship[J]. Remote Sensing of Environment, 2003, 85(4):429-440 doi: 10.1016/S0034-4257(03)00036-1 [8] Yang H J, Cong Z T, Liu Z W, et al.Estimating Subpixel Temperatures Using the Triangle Algorithm[J]. International Journal of Remote Sensing, 2010, 31:6047-6060 doi: 10.1080/01431160903376373 [9] Sandholt I. A Simple Interpretation of the Surface Temperature/Vegetation Index Space for Assessment of Surface Moisture Status[J]. Remote Sensing of Environment, 2002, 79(2):213-224 http://www.docin.com/p-1155858698.html [10] Chen X L. Remote Sensing Image-Based Analysis of the Relationship Between Urban Heat Island and Land Use/Cover Changes[J]. Remote Sensing of Environment, 2006, 104(2):133-146 doi: 10.1016/j.rse.2005.11.016 [11] Essa W, Verbeiren B, Kwast J, et al.Evaluation of the DisTrad Thermal Sharpening Methodology for Urban Areas[J]. International Journal of Applied Earth Observation and Geoinformation, 2012, 19:163-172 doi: 10.1016/j.jag.2012.05.010 [12] 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]. International Journal of Remote Sensing, 2001, 22(18):3719-3746 doi: 10.1080/01431160010006971 [13] 覃志豪, Zhang Minghua, Karnieli A, 等.用陆地卫星TM6数据演算地表温度的单窗算法[J].地理学报, 2001, 56(4):456-466 http://www.cnki.com.cn/Article/CJFDTOTAL-DLXB200104008.htm Qin Zhihao, Zhang Minghua, Karnieli A, et al. Mono-Window Algorithm for Retriving Land Surface Temperature from Landsat TM 6 Data[J]. Acta Geographica Sinica, 2001. 56(4):456-466 http://www.cnki.com.cn/Article/CJFDTOTAL-DLXB200104008.htm [14] 王斐, 覃志豪, 王倩倩.基于地表类型的TM6波段像元分解方法[J].国土资源遥感, 2012, 94(3):54-59 http://www.cnki.com.cn/Article/CJFDTOTAL-GTYG201203012.htm Wang Fei, Qin Zhihao, Wang Qianqian. A Method of TM6 Band Pixel Decomposition Based on the Earth Surface Types[J]. Remote Sensing for Land & Resource, 2012, 94(3):54-59 http://www.cnki.com.cn/Article/CJFDTOTAL-GTYG201203012.htm [15] Zhu S, Guan H, Millington A C, et al. Disaggregation of Land Surface Temperature over a Heterogeneous Urban and Surrounding Suburban Area:A Case Study in Shanghai, China[J]. International Journal of Remote Sensing, 2013, 34(5):1707-1723 doi: 10.1080/01431161.2012.725957 -