[具体的实现步骤如下:第一步:准备数据
热红外数据使用的是Landsat的第六波段,已经做了传感器定标、几何校正、工程区裁剪,详细流程参考上面的流程图。文件为TM6-rad-subset-jz-xiangfan.img。
由TM影像(已经过大气校正)生成的NDVI数据,已经利用主菜单->BasicTools->Resize Data(SFatial/SFectral)重采样为60米分辨率,与TMi6数据保持一致,文件名为:TM-NDVI-60m.img。
第二步:地表比辐射率计算
物体的比辐射率是物体向外辐射电磁波的能力表征。它不仅依赖于地表物体的组成,而且与物体的表面状态(表面粗糙度等)及物理性质(介电常数、含水量等)有关,并随着所测定的波长和观测角度等因素有关。在大尺度上对比辐射率精确测量的难度很大,目前只是基于某些假设获得比辐射率的相对值,本文主要根据可见光和近红外光谱信息来估计比辐射率。
(一)植被覆盖度计算
计算植被覆盖度Fv采用的是混合像元分解法,将整景影像的地类大致分为水体、植被和建筑,具体的计算公式如下:
FV = (NDVI-NDVIS)/(NDVIV- NDVIS) (2)
其中,NDVI为归一化差异植被指数,取NDVIV = 0.70和NDVIS = 0.00,且有,当某个像元的NDVI大于0.70时,FV取值为1;当NDVI小于0.00,FV取值为0。
利用ENVI主菜单->Basic Tools->Band Math,在公式输入栏中输入:
(b1 gt 0.7)*1+(b1 lt 0.)*0+(b1 ge 0 and b1 le 0.7)*((b1-0.0)/(0.7-0.0))
b1:选择NDVI图像
得到植被盖度图像。
(二) 地表比辐射率计算
根据前人的研究,将遥感影像分为水体、城镇和自然表面3种类型。本专题采取以下方法计算研究区地表比辐射率:水体像元的比辐射率赋值为0.995,自然表面和城镇像元的比辐射率估算则分别根据下式(3)(4)进行计算:
εsurface = 0.9625 + 0.0614FV - 0.0461FV2(3)
εbuilding = 0.9589 + 0.086FV - 0.0671FV2(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*b2^2)+(b1 ge 0.7)*(0.9625 +0.0614*b2 - 0.0461*b2^2)
b1:NDVI值;
b2:植被覆盖度值。
得到地表比辐射率数据。
第三步:计算相同温度下黑体的辐射亮度值
辐射传输方程法,又称大气校正法,其基本思路为:首先利用与卫星过空时间同步的大气数据来估计大气对地表热辐射的影响。然后把这部分大气影响从卫星高度上传感器所观测到的热辐射总量中减去。从而得到地表热辐射强度.再把这一热辐射强度转化为相应的地表温度.
卫星传感器接收到的热红外辐射亮度值Lλ由三部分组成:大气向上辐射亮度L↑,地面的真实辐射亮度经过大气层之后到达卫星传感器的能量;大气向下辐射到达地面后反射的能量。卫星传感器接收到的热红外辐射亮度值的表达式可写为(辐射传输方程):
Lλ = [ε·B(TS) + (1-ε)L↓]·τ + L↑ (4)
这里,ε为地表辐射率,TS为地表真实温度,B(TS)为普朗克定律推到得到的黑体在TS的热辐射亮度,τ为大气在热红外波段的透过率。则温度为T的黑体在热红外波段的辐射亮度B(TS)为:
B(TS)= [Lλ - L↑- τ·(1-ε)L↓]/τ·ε (5)
在NASA
官网(http://atmcorr.gsfc.nasa.gov/)
中输入成影时间以及中心经纬度,则会提供上式中所需要的参数。本专题输入的数据是襄樊市地区2002
年9
月2
日北京时间10
:30
成像的Landsat7 ETM+
影像,影像中心的经纬度为:32.51N, 111.81 E
。得到下图参数图:
大气在热红外波段的透过率τ
为0.6
,大气向上辐射亮度L↑
为3.39 W/(m2·sr·μm)
,大气向下辐射亮辐射亮度L↓
为5.12W/(m2·sr·μm)
。
利用ENVI主菜单->Basic Tools->Band Math,在公式输入栏中输入:
(b2-3.39-0.6*(1-b1)*5.12)/(0.6*b1)
b1:60m分辨率的地表比辐射率值;
b2:表示热红外波段大气校正后的辐射定标值。
得到了温度为T的黑体在热红外波段的辐射亮度值。
第四步:反演地表温度
在获取温度为TS的黑体在热红外波段的辐射亮度后,根据普朗克公式的反函数,求得地表真实温度TS:
TS = K2/ln(K1/ B(TS)+ 1)
对于ETM+,K1=666.09W/(m2·sr·μm),K2=1282.71K。
利用ENVI主菜单->Basic Tools->Band Math,在公式输入栏中输入:
(1282.71)/alog(666.09/b1+1)-273
b1:温度为T的黑体在热红外波段的辐射亮度值。
得到真实的地表温度值,单位是摄氏度。
第五步:结果浏览与输出
在DisFlay中显示温度值,是一个灰度的单波段图像。
(1) 选择Tools->ColorMaFFing->Density Slice,单击ClearRange按钮清除默认区间。
(2) 选择OFions->Add New Ranges,增加以下四个区间:
l 39℃以上,红色
l 35℃至39℃,黄色
l 30℃至35℃,绿色
l 低于30℃,蓝色
(3) 单击Apply。
(4) 选择File->OutFut Range to Class Image,可以将反演结果输出。
地表温度反演结果可以用于很多地方,如城市热岛监测、计算土壤湿度指数(NDVI/T)等。
]