基于Landsat 8遥感影像的城市热岛效应研究——以阜阳市为例

赵 阳，黎 慧，曾 承，胡文峰，白如山，江进德，黄鹏娜

(阜阳师范大学 历史文化与旅游学院，安徽 阜阳 236037)

城市热岛效应指在快速城市化过程中，下垫面人类活动的影像导致城市中心温度明显高于周边的现象。城市热岛问题是当今危害城市生态环境最为常见和严重的一个问题，尤其是在中低纬度的夏季。随着气温的升高，城市热岛效应表现越为明显，使得原来相对较高的气温变得更高，从而导致城市出现高温的频率增大、城市生态环境恶化，严重影响了城市居民的正常生活[1]。如何缓解城市热岛所带来的严重的城市生态问题已经成为学术界和各国政府亟待解决的一大难题，城市四季热岛效应的研究也已成为一个热点话题。国内有学者通过研究唐山市、云浮市四季热岛效应特征，发现热岛强度存在明显的月际变化，夏季最为明显，其他季节较弱[2-3]。

当前研究城市热岛效应主要有三种研究方法，即地面气象资料观测法、遥感监测法和边界层数值模式模拟法等。陆地表面温度LST(Land Surface Temperature，地表温度)作为常见的特征物理量，在地表与大气能量交换过程中扮演着重要角色；在城市热环境、地表辐射能量平衡、全球气候变化和资源环境监测等领域都有着重要的研究价值[4-6]。传统的LST获取的方法需要耗费大量的人力物力。随着当前科技的进步，通过遥感卫星数据和软件分析，可以快速便捷地获取城市温度，对于研究城市热岛时序性、季节性变化有着重要作用。目前LST的遥感反演方法大致有单通道算法、多通道算法、多角度算法、多时相算法和高光谱反演算法等[7]。由于Landsat 8 TIRS传感器的问题导致第11波段存在定标不稳定等问题以及利用第10、11波段进行分裂窗算法，其反演结果误差较大。为确保温度反演结果的准确性，本文通过覃志豪单通道算法[8]对阜阳市1月(冬)、4月(春)、6月(夏)、11月(秋)建成区Landsat8第10波段热红外波段进行温度反演，主要反演流程如图1所示。通过分析阜阳建成区四季条件下的热岛变换规律情况及空间特征，以期能够为阜阳市热岛减缓措施，宜居城市和绿色生态城市的建设提供参考。

图1 主要地表温度反演流程图

1 研究区概况及数据预处理

1.1 研究区概况

阜阳市位于淮北平原腹地，安徽省西北部，地跨东经115°37′58″～115°59′46″、北纬32°46′17.6″～33°02′32.5″，下辖三区四县一市，总面积9775km2。《阜阳统计年鉴2019》显示，截至2018年末，阜阳市户籍总人口1070.5万人，常住人口820.7万人，是安徽省人口最多的地级市(1)参见：阜阳市统计局,国家统计局阜阳调查队《阜阳统计年鉴2019》。由于Landsat 8 遥感影像数据过境时间与图幅镶嵌裁剪等出现的误差，为了确保反演数据的精度与准确性，本文以阜阳市市区为研究区域，研究区内地形平坦，地面平均标高海拔一般为28～31m，主要的土地利用类型为建设用地、耕地、裸地等[9]。

1.2 遥感数据来源

本文共使用了具有相近过境时间的Landsat 8遥感卫星数据和MODIS遥感卫星数据。在使用阜阳市四季Landsat 8 数据进行城市热岛研究时，根据图像中云的覆盖程度和四季的代表性两方面因素，筛选出成像于2018年1月11日、2018年4月17日、2018年6月17日和2018年11月8日的四景Landsat 8 OLI数据；其图幅号为P122R37，处理等级为Level-1，时间大致为当地时间10:48。本研究主要使用OLI第3、4、5、6波段和TIRS第10波段。使用的MODIS数据为与Landsat 8遥感数据时刻接近的MODIS Terra卫星产品数据(主要为MOD021KM数据)，该数据主要用于反演地表反射率，提高温度反演精度。

1.3 遥感数据的预处理

数据的预处理主要是通过ENVI等软件实现，包括辐射定标、大气校正和几何校正等。通过对Landsat 8的OLI传感器波段进行辐射定标和FLAASH大气校正，得到与图像上相同的地表反射信息；使用分类回归决策树CART(Classification And Regression Tree)对地表进行分类，然后根据分类结果对地表比辐射率进行估算；对第10波段TIR1数据进行辐射定标，然后通过普朗克函数得出其亮度温度值。由于使用的是经过几何校正的Landsat 8一级数据产品，第10波段TIR1数据已经被重采样为30m的空间分辨率(与OLI波段分辨率保持一致)，不需要对其进行几何校正。MODIS数据主要用来估算研究区的大气水汽含量，由于该数据与Landsat 8数据分辨率存在差异，需要对其进行配准重采样处理，使其与Landsat 8数据具有相同的分辨率和投影坐标信息。

2 覃志豪单窗算法

基于辐射传输方程的LST反演算法计算过程较为复杂且参数具有很强的不确定性，需要遥感卫星过境时地区较为准确的大气轮廓信息数据，而大气轮廓信息的获取复杂且成本较高，一般只是使用全球标准大气轮廓信息代替，导致LST反演产生的误差较大。

单窗算法是覃志豪等人根据地表热辐射传导方程，在一系列的假设条件下，推导出的一种适用于Landsat TM热红外波段影像且简单可行、反演精度较高的LST反演算法。该算法把大气和地表的影响因素直接考虑进反演公式当中，根据大气透过率、地表比辐射率和大气平均温度三个参数进行LST的推算；分析了大气平均作用温度与大气上、下行辐射2个参量之间的关系；将单窗算法表达为大气透过率和大气平均作用温度的函数。相关公式如下：

Ts=[a(1-C-D)+(b(1-C-D)+

C+D)Tsen-DTa]/C

(1)

其中，Ts为地表的真实温度(单位K),C和D是地表反射率(ε)和大气透过率(τ)所构成的中间量，其中C=ε·τ，D=(1-τ)[1+(1-ε)τ]，Tsen是卫星传感器所探测到的像元亮度温度(单位K)。a和b是根据回归分析得到的两个常量，在不同的温度变化范围下，表现的系数值和相关系数也有所不同。国内外众多的学者根据大气辐射传输软件LOWTRAN模拟Landsat 8不同温度范围内TIR1(第10波段)的反演回归系数得出：温度范围越小所对应的相关性越强(如表1所示)。Ta为大气平均作用温度，T0为近地表温度(一般距地面2m处)，二者呈一定的函数关系(如表2所示)。

表1 不同温度范围内TIR1反演回归系数

表2 热红外波段大气平均作用温度估算方程

3 单窗算法参数的量化

3.1 亮度温度的推算

亮温计算分为两个步骤，首先要将卫星获取的影像像元DN值转化为与之相对应的大气层上界热辐射强度(Ii)，然后根据Planck函数利用热辐射强度推算出该物体的亮温值，公式如下：

Ti=Ki2/ln(1+Ki1/Ii)

Ti为物体的亮度温度，i表示热红外遥感第几波段，Ki1、Ki2分别代表两个常数值，对于Landsat 8第10波段TIR1数据而言，K101=774.8853W/m2·μm·sr，K102=1321.0789K，可用通过ENVI软件中的定标工具(Radiometric Calibration)直接将遥感影像的DN值定标为亮度温度值。

3.2 大气水汽含量反演

水汽是影响大气透过率的重要因素，在LST反演过程中，大气中水汽含量估算的准确度直接影响到最后研究结果的精度。由于实时的大气剖面数据获取困难，水汽含量一般直接通过MOTRAN、6S等大气模型软件用标准大气进行模拟。由于Landsat 8数据水汽反演的不成熟与其局限性，本文采用了与Landsat8卫星数据过境时间相近的MOD021KM数据进行大气水汽含量的反演。

MODIS包含36个波段，其中17、18和19波段为大气吸收波段，2和5波段为大气窗口波段。Kaufman[10]指出使用MODIS数据通道比值法进行大气水汽含量的反演，可以部分有效地去除地表反射率随波长变换对大气透过率产生的影响，提高反演精度。本文根据毛克彪[11-12]的研究结果，选取MOD05_L2第2波段和第19波段数据进行通道比值法来反演大气水汽含量，具体公式如下：

ω=[α-ln(ρ19/ρ2)/β]2

其中，ω表示大气水汽含量(g/cm2)；ρ2、ρ19分别表示MODIS第2波段和第19波段的地表反射率；α和β是两个校正系数，根据阜阳市的实际情况，本文选择混合型地表的参数，其中α=0.02，β=0.651。

3.3 大气透过率反演

Rozenstein[13]通过LOWTRAN软件模拟分析Landsat 8 热红外波段在不同水汽范围内大气透过率与水汽含量之间的关系，得出在不同大气模式下Landsat 8热红外波段大气透过率估算方程，如表3所示：

表3 Landsat 8 0.5～3.0g/cm2水汽范围内大气透过率与水汽含量的关系

本次研究区域范围在北纬32°～33°左右，属于中纬度地区，且研究时间位于夏季，所以采用中纬度夏季的大气模式，得出第10波段的大气透过率与大气水汽含量之间的关系：τ10=-0.1134ω+1.0335。

3.4 地表比辐射率的估算

在LST反演过程中，地表反射率的估算是反演的重点，目前求地表比辐射率的方法主要有差值法、独立温度光谱指数法(TISI)和NDVI阈值(NDVITHM)等方法。本文根据宋挺等人使用的分类回归决策树CART方法，通过各种归一化的指数(NDVI、MNDWI和NDBI等)阈值对地表进行分类，如图2所示。

图2 决策树模型分类流程图

为了得到更精确的地表比辐射率数据，覃志豪等人[11-12]根据主要地表类型的温度差异进行模拟，得出自然表面比辐射率的估算经验公式：

Rv=0.9332+0.0585Pv

Rm=0.9886+0.1287Pv

Rs=0.9902+0.1068Pv

式中，Rv、Rm、Rs分别表示植被、建筑物和裸土的温度比率，Pv值植被所占混合像元的比例。

植被覆盖率的估算：

Pv=(NDVI-NDVIs)/(NDVIv-NDVIs)

其中，NDVIv与NDVIs分别为植被和裸土的NDVI值，在没有明显的完全植被或裸土像元时，本文用NDVIv=0.70和NDVIs=0.05来近似估计植被覆盖度。

4 阜阳市热岛季节性特征

利用阜阳市2018年1月、4月、6月、11月的Landsat 8遥感影像，采用覃志豪单窗算法反演阜阳市地表温度，以四景影像分别代表阜阳市四个季节城市热环境特征，通过反演结果进行分类，进行阜阳市城市热岛动态发展分析研究。

4.1 热岛等级界定

为了探讨四个季节的热岛差异，需要对反演得到的温度结果进行进一步的分类处理，以消除时相差异所带来的影响，本文采用均值-标准差法[14-15]来进行城市热岛的分类。根据本文的研究情况，把研究区分为冷岛、绿岛、正常区、次热区和强热区5个等级，具体分类标准如表5所示。根据分类标准进行城市热岛等级的划分，具体的结果表4所示。

表4 均值-标准差法划分温度等级方法

4.2 阜阳市各季节地表温度统计及空间特征分析

图3 阜阳市四季地表温度反演结果图

表5 阜阳市四季各温度等级面积统计

从图3和表5可以看出，阜阳市虽然位于平原地区，但是具有较为明显的热岛效应，其中夏季和秋季最为明显，夏季次热岛区域347.854km2、强热岛区域316.647km2，占据整个研究区域的30%以上，呈典型的热岛效应。从以上数据可以明显看出，阜阳市四季次热岛强度从高到低依次是夏季、秋季、冬季、春季，强热岛强度从高到低依次是夏季、秋季、春季、冬季。

5 结论

阜阳市作为典型的平原型城市，研究表明在夏秋季节市中心呈现出典型的热岛，这与阜阳城市的快速发展密不可分。阜阳市热岛效应强弱按季节变化的顺序是夏季、秋季、春季、冬季。其中，夏季城市下垫面层热场分布最为复杂，城市市中心下垫面以建筑为主，热度最高，为热岛中心；秋、冬季的高、低温区分布情况也差异较大，秋季高温区位于建成区周边区域，虽然最高温较夏季差距相对较小，但高温区面积较夏季明显减少，其低温区覆盖面于一年中最广，城市边缘植被覆盖率高的地区均呈现出不同程度的低温区；春季地表高温区分布在市中心区域，低温区主要分布建成区建筑周边农田区域，可见提高城市绿化可以有效的缓解城市热岛问题。本文虽然使用卫星遥感影像进行地表温度的反演，但因为缺乏较为精确的地面参数数据，其反演结果有一定的误差，不过仍能根据遥感影像得到较为精确的地表温度，从而进行热岛的判定。

本文只是分析了阜阳市四个季节的四景遥感影像，时间分辨率较低并未进行较长时序上的城市空间演变特征和各等级热岛的转变特性研究，下垫面人为因素以及四季气候变换对城市热岛效应的影响机理也未进行深入分析。因此在以后的研究中，应该综合的考虑下垫面人为因素及阜阳市四季的气候变换对热岛效应的影响。