毕业设计（论文）-基于ETM+数据的福州市地表温度反演.doc

基于etm+数据的福州市地表温度反演1 引言随着城市化进程的不断加速，城市热岛效应表现的尤为突出。地表温度作为城市热岛效应、地球环境分析的重要指标，而遥感技术作为现代重要的对地观测手段，使得基于遥感图像的地表温度反演的研究越来越多。地表温度反演方法有很多，如：大气校正法，单窗算法，单通道法等等。本实验主要采用landsat etm+遥感影像，运用辐射传输方程法（大气校正法）对地表温度进行反演。该方法需要进行大气校正，消除大气层对地表辐射能量的影响，这就需要从卫星观测得到的热辐射能量中扣除大气层的辐射分量，并利用热红外波段（band61、62）范围内的地表发射率作为参数，反演出地表的真实温度。本实验详细介绍了利用landsat etm+图像反演福州市地表温度的完整流程，实验涉及landsat etm+的数据读取、辐射定标、反演模型、遥感反演过程、反演结果验证等内容。主要采用的技术：envi 主模块中的landsat etm+数据读取、热红外波段数据的辐射定标、bandmath使用、最终的结果统计分析等。在envi4.1版本完成整个操作过程。2实验区概况福州市是福建省的省会，位于闽东南沿海，台湾海峡西岸，地处东经11808-12031，北纬2515-2639，土地利用类型主要有林地、建设用地、耕地、水体和未利用地。属暖热湿润的亚热带季风气候，年平均温度1722,平均降雨量1400mm2000mm，盛行东北风。福州市属河口盆地地貌，地势由西北内陆向东南沿海倾斜，主城区位于盆地的盆心，东部、北部为山地丘陵，南部为平原。闽江自西北向东南穿流而过，由于向海的开口很小，地形相对封闭，加之建筑物的阻挡,不利于气流流动,造成城区温度较高，夏季中午气温高达36以上。特别是近年来，福州市经济的迅速发展，加快了城市化进程的步伐和建成区扩展的速度，热岛效应严重。3数据本文数据来源于地理空间数据(/)上下载的福州市地区landsat7etm+影像。影像获取时间为北京时间2002年5月26日，14:23；条带号为119，行号为42；中心经度118.8428，中心维度25.9892；平均云量0.34；太阳高度角66.5530太阳方位角96.8191。4 实验步骤4.1数据读取和定标数据定标采用的是envi 下自带的landsat 定标工具“landsat calibration”，在读取文件数据后，会根据元数据中的信息进行增益值与偏移值的自动读取，直接根据功能对话框点击“确定”即可。4.2 可见光-近红外大气校正消除大气的影响以获取精度更高的ndvi 及植被覆盖度，用于计算地表比辐射量。这里采用的是快速大气校正。4.3 反演1)辅助数据的获取本专题采用的是大气校正法进行的地表温度反演，其中涉及到大气上行辐射、下行辐射以及大气的透过率数据，这些数据可以在nasa 官网(/)中输入成像时间以及其中心经纬度，则会提供上式中所需要的参数。2)地表比辐射率的计算地表比辐射率值的计算是地表温度反演的一个难点，本专题主要采用先验公式的方法利用植被覆盖度进行地表比辐射率值的计算。3)地表真实温度的获取专题中反演部分是采用反演公式在envi 主模块下的band math 进行运算实现的。4.4 路线图landsat etm+原始数据辐射定标热红外数据可见光数据快速大气定标ndvi值计算植被盖度热红外波段辐射亮度值地表真实温度分布不规则裁剪结果输出、分析结果反演数据预处理5 详细步骤5.1 数据读取与辐射定标1)数据读取在envi 主菜单中，选择file- open external file- landsat-geotiff with metadata，打开l71119042_04220020526_mtl.txt。包含了三种数据：可将光波段数据（hrf），热红外波段数据（htm）和全色波段数据（b80）。2)辐射定标主菜单-basic tools-preprocessing-calibration utilities-landsat calibration，选择需要定标的数据，包括可见光波段数据和热红外波段数据，该过程功能可以从读取文件中直接获取实验数据的元数据（成像时间，定标参数等等）。定标后值可见光数据辐射定标参数设置5.2 快速大气纠正选择spectral -quick atmospheric correction，对定标后的可见光数据执行快速大气校正。快速大气校正可见光波段校正后结果5.3 地表比辐射率值的计算物体的比辐射率是物体向外辐射电磁波的能力表征。它不仅依赖于地表物体的组成，而且与物体的表面状态(表面粗糙度等)及物理性质(介电常数、含水量等)有关，并随着所测定的波长和观测角度等因素有关。在大尺度上对比辐射率精确测量的难度很大，目前只是基于某些假设获得比辐射率的相对值，本实验主要根据可见光和近红外光谱信息来估计比辐射率。1)归一化植被指数(ndvi) 利用etm3、4波段的象元dn值利用下面的公式（1）求得归一化植被指数ndvi：ndvi = (nir-r)/(nir+r) （1）这里，nir 和r 分别是etm 的近红外波段(波段4)和红光波段(波段3)的dn 值。采用envi的直接采用envi 主模块下的ndvi 模块获取：主菜单-transform-ndvi计算结果如下图：2)植被覆盖度本实验计算植被覆盖度pv采用的是混合像元分解法，将整景影像的地类大致分为水体、植被和建筑，具体的计算公式如下： pv= (ndvi- ndvis)(ndviv - ndvis)2 （2）其中，ndvi 为归一化差异植被指数，取ndviv = 0.70 和ndvis = 0.00，且有，当某个像元的ndvi 大于0.70 时，pv 取值为1；当ndvi 小于0.00，pv 取值为0。利用envi主菜单-basic tools-band math，在公式输入栏中输入：（b1 gt 0.7）*1+（b1 lt 0.0）*0+（b1 ge 0 and b1 le 0.7）*（(b1-0.0)/(0.7-0.0)）其中，b1为获取的ndvi值。输入公式b1为ndvi植被覆盖影像图以及统计结果如下：植被覆盖度影像数据统计结果：3)地表比辐射率根据前人的研究，将遥感影像分为水体、城镇和自然表面3 种类型。本专题采取以下方法计算研究区地表比辐射率：水体像元的比辐射率赋值为0.995，自然表面和城镇像元的比辐射率估算则分别根据下式(3)(4)进行计算： surface = 0.9625 + 0.0614pv- 0.0461pv2 （3） building= 0.9589 + 0.086pv- 0.0671 pv2 （4）式中，surface和building分别代表自然表面像元和城镇像元的比辐射率。利用envi主菜单-basic tools-band math，在公式输入栏中输入：（b1 le 0）*0.995+（b1 gt 0 and b1 lt 0.7）*(0.9589 + 0.086*b2 - 0.0671*b22)+(b1 ge0.7)*(0.9625 + 0.0614*b2 - 0.0461*b22)式中，b1表示归一化植被指数值；b2表示植被盖度值。b1为ndvi,b2为植被盖度值输入公式下图为30米分辨率的地表比辐射率值。上述公式计算的结果分辨率为30m，由于热红外波段数据分辨率为60m，故而需要将其分辨率重采样为60m，目的是为下面反演地表温度提供基础数据。具体的操作步骤为利用envi 主菜单-basic tools-resize data（spatial/spectral），出现对话框设置：图像数据重采样6辅助参数获取卫星传感器接收到的热红外辐射亮度值l由三部分组成：大气向上辐射亮度l，以及地面的真实辐射亮度经过大气层之后到达卫星传感器的能量。地面的真实辐射亮度为同温度黑体的辐射亮度值lt与地物发射率的乘积lt。即，大气校正法的表达式可写为：l = lt + (1-)l + l （5）这里，t为地表真实温度，为大气在热红外波段的透过率.则温度为t的黑体在热红外波段的辐射亮度lt为：lt = l - l- (1-)l/（6）在nasa官网(/)中输入成影时间以及中心经纬度，则会提供上式中所需要的参数。本专题输入的数据是福州市地区2002年9月2日北京时间10：30成像的landsat7 etm+影像，影像中心的经纬度为：25.9892n, 118.8428 e。得到下图参数图: 大气在热红外波段的透过率为0.7，大气向上辐射亮度l为2.17 w/(m2srm)，大气向下辐射亮辐射亮度l为3.47w/(m2srm)。如能获取气象站当日气象信息则输入，否则不填（采用默认值）nasa大气校正网站页面输入待校正影像信息下行辐射亮度上行辐射亮度大气透过率2002年5月26日landsat etm+数据的大气辅助参数利用envi主菜单-basic tools-band math，在公式输入栏中输入：（b2-2.17-0.7*（1-b1）*3.47）/（0.7*b1）式中，b1为60m 分辨率的地表比辐射率值；b2表示热红外波段大气校正后的辐射定标值。获取的黑体在热红外波段的辐射亮度结果。7地表温度反演在获取热红外波段辐射亮度值以后根据普朗克公式的反函数，求得地表真实温度t：t = k2/ln(k1/lt + 1)对于etm+，k1 =666.09w/(m2srm)，k2 =1282.71k。利用envi主菜单-basic tools-band math，在公式输入栏中输入：（1282.71）/alog（666.09/b1 +1）式中，b1为热红外辐射亮度值。计算结果如下：地表真实温度结果8 裁剪工程区（1）显示最终反演后的数据t.img，在image 窗口选择overlay-vectors，打开vector parameters 面板，选择file-open vector file，打开fz.shp；（2）在available vectors list 面板中选择该矢量文件，点击load selected，选择显示图像的display，单击ok，矢量叠加在影像上；（3）在available vectors list 面板中，选择file-export layers to roi，在select data fileto associate with new roi 面板中，选择t.img，单击ok；在export evf layers to roi 中，选择convert all records of an evf layer to one roi，单击ok，将矢量转为一个roi；（4）在图像窗口，选择overlay-region of interest，打开roi 面板，太湖区roi 显示在图像上，在roi tool 面板中，选择file-subset data via rois，在select input file tosubset via roi 面板中，选择t.img，单击ok；（5）在spatial subset via roi parameters 面板中选择工程区roi，mask pixels outside ofroi 选择yes。（6）设置输出路径及文件名t_subset.img，单击ok。裁剪后结果9 反演结果输出在display 中，选择tools-color mapping-density slice，将结果以进行密度分割，分割为4 类，具体的分割阈值以及参数设置如下图：反演结果10 反演结果分析利用envi 软件，采用单波段彩色变换方法对地表真实温度t 的灰度图进行密度分割，得到地表真实温度的分布情况，如上图所示。其中，红色代表的温度范围是：39以上（其中40以上的象元数很少，为异常点）；黄色代表的温度范围是：35至39；绿色代表的温度范围是：30至35；蓝色代表的温度范围是：低于30。对于2002 年5月26日天气温度分析，反演结果符合当时的天气温度情况。11 参考文献1 张友水,冯学智,周成虎.多时相tm影像相对辐射校正研究j.测绘学报,2006,35(2):122-127.2 季青,贺伶俐,余明等.基于landsat etm+数据的福州市土地利用/覆被与城市热岛研究j.福建师范大学学报(自然科学版),2009,25(6):106-1132.3 覃志豪，张明华，等用陆地卫