地表温度反演.doc

地表温度反演

目录

一：单窗算法 3
1.1实验原理 3
1.1.1TM/ETM波段的热辐射传导方程： 3
1.1.2化简后最终的单窗体算法模型为： 3
1.1.3大气平均作用温度Ta的近似估计 3
1.1.4大气透射率t6的估计 3
1.1.5地表比辐射率的估计 4
1.1.6像元亮度温度计算 4
1.1.7遥感器接收的辐射强度计算 4
1.2操作步骤 5
1.2.1研究区示意图 5
1.3实验结果 7
1.3.1灰度图像 7
1.3.2密度分割后图像 7
二：单通道算法 8
2.1实验原理 8
2.1.1单通道算法模型为： 8
2.1.2大气平均作用温度Ta的近似估计 8
2.1.3大气透射率t6的估计 8
2.1.5像元亮度温度计算 8
2.1.6遥感器接收的辐射强度计算 9
2.2操作步骤 9
研究区示意图 9
2.2.1计算L6 10
2.2.2T6e6的求算 10
2.2.3计算R 10
2.2.4计算t 10
2.3实验结果 11
2.3.1温度反演灰度图像 11
2.3.2密度分割后的图像 11
三：辐射方程 12
3.1实验过程 12
3.1.1数据准备 12
3.1.2地表比辐射率的估计 12
3.1.3计算同温度下黑体的辐射亮度值 12
3.1.4反演地表温度 13
3.2温度反演结果 13

一：单窗算法
1.1实验原理
1.1.1TM/ETM波段的热辐射传导方程：
B6(T6)=t6(q)[ ε6B6(Ts)+(1-ε6)I6~]+I6_
Ts是地表温度；
T6是TM6的亮度温度；
t6是大气透射率；
ε6是地表辐射率。
B6(T6)表示TM6遥感器所接收到的热辐射强度；
B6(Ts)是地表在TM6波段区间内的实际热辐射强度，直接决取于地表温度；
I6~和I6_分别是大气在TM6波段区间内的向上和向下热辐射强度。

1.1.2化简后最终的单窗体算法模型计算Ts（地表温度）：
Ts={a(1-C-D)+[b(1-C-D)+C+D]T6-DTa}/C
式中
C6=τ6ε6(ε6为比辐射率，τ6为透射率）
D6=(1-τ6)[1+t6(1-ε6)]
a =-67.355351，b=0.458606
1.1.3大气平均作用温度Ta的近似估计
温度换算：T＝t＋273.15
本图为9月份拍摄，对于中纬度夏季平均大气Ta=16.0110+0.92621T0
取平均气温为25摄氏度时Ta = 312.15753
1.1.4大气透射率τ6的估计
τ6=0.974290-0.08007w，0.4≤w≤1.6。
w为水分含量，单位(g/cm2)，这里，取w=1.0，计算得到τ6=0.89422
1.1.5地表比辐射率ε6的估计
典型地表类型的比辐射率
Ndvi<=0 水体， εw=0.995
NDVI>0 陆地（surface） εsurface=0.9625+0.0614Pv-0.0461Pv^2
Pv=[(NDVI- NDVIs)/(NDVIv- NDVIs)]2=NDVI^2/0.49
(NDVI>0 and NDVI <0.70)
式中NDVI 为归一化植被指数,
取NDVIv=0.70 和 NDVIs=0（分别取自5%及95%数据）

ew算法：0.995*(b1 le 0)
Pv算法：((b1*b1)/0.49)*(b1 gt 0 and b1 le 0.7)
esurface 算法：(0.9625+0.0614*b1-0.0461*b1^2)*(b2 gt 0)
e6算法：ew+esurface
1.1.6像元亮度温度T6计算
T6=K2/ln(1+K1/L6)
其中：K1=607.76，K2=1260.56（覃志豪，用陆地卫星TM6数据演算地表温度的单窗算法）；
L6为遥感器接收的辐射强度。
1.1.7遥感器接收的辐射强度L6计算
L6=(15.303-1.238)*b1/255.0+1.238 （覃志豪，用陆地卫星TM6数据演算地表温度的单窗算法）； b1为第六波段像元灰度值（DN值）。
LMAX_BAND6 = 15.303
LMIN_BAND6 = 1.238

1.2操作步骤
1.2.1研究区示意图

1.2.2计算L6
ENVI->Basic Tools->Band Math

在Enter an expression中输入(15.303-1.238)*b1/255.0+1.238，点击Add to List.
选中刚才输入的公式，点击OK。

选中B1变量，在Available Bands List中选Band6，保存文件名为L6，点击OK。
1.2.3计算T6（同计算L6,只是计算公式不同）。
1.2.4计算C6（同计算L6,只是计算公式不同）。
1.2.5计算D6（同计算L6,只是计算公式不同）。
1.2.36计算Ts
同计算C6相同，在Enter an expression中输入：
(((-67.355351*(1-b1-b2)+(0.458606*(1-b1-b2)+b1+b2)*b3-b2*292.15753)/b1-273.15))，点击Add to List.
选中b1变量，在Available Bands List中选择C6，
选中b2变量，在Available Bands List中选择D6，
选中b3变量，在Available Bands List中选择T6，保存文件名为Ts ，点击OK。
1.3实验结果
1.3.1灰度图像

1.3.2密度分割后图像

二：单通道算法
2.1实验原理
2.1.1单通道算法模型为：
Ts=r*[(y1*L6+y2)/e6+y3]+t

2.1.2大气平均作用温度Ta的近似估计
温度换算：T＝t＋273.15
本图为7月份拍摄，对于中纬度夏季平均大气Ta=16.0110+0.92621T0
取乌鲁木齐市平均气温为25摄氏度时Ta = 312.15753

2.1.3大气透射率t6的估计
t6=0.974290-0.08007w，0.4≤w≤1.6。
w为水分含量，单位(g/cm2)，这里，取w=1.0，计算得到t6=0.89422

2.1.4地表比辐射率的估计
典型地表类型的比辐射率
ew=0.995
ev=0.986
em=0.970
Pv=[(NDVI- NDVIs)/(NDVIv- NDVIs)]2
NDVI 为归一化植被指数, 取NDVIv=0.70 和 NDVIs=0（分别取自5%及95%数据）
e surface=0.9625+0.0614Pv-0.0461Pv^2

2.1.5像元亮度温度计算
T6=K2/ln(1+K1/L6)
其中：K1=607.76，K2=1260.56（覃志豪，用陆地卫星TM6数据演算地表温度的单窗算法）；
L6为遥感器接收的辐射强度。
2.1.6遥感器接收的辐射强度计算
L6=(15.303-1.238)*b1/255.0+1.238 （覃志豪，用陆地卫星TM6数据演算地表温度的单窗算法）； b1为第六波段像元灰度值（DN值）。
LMAX_BAND6 = 15.303
LMIN_BAND6 = 1.238
单通道算法模型表达式中, Ts是陆地表面温度, L6是卫星高度上遥感传感器测得的辐射强度(W·m- 2·sr- 1·m- 1),e6是地表发射率; r, t,y1,y2,y3 是中间变量, 分别由以下公式计算:
R=1/[(c2*L6/T6^2)*(λ^4*L6/c1+1/λ)]
y1=0.14714w^2- 0.15583w+1.1234=1.11471
y2=-1.1836w^2- 0.37607w- 0.52894=0.27859
y3=-0.04554w^2+1.8719w- 0.39071=1.43565
t=T6-r*L6
其中, C1 和C2 是Plank 函数的常量, C1=1.19104×10^8
W·μm4·m- 2·sr- 1, C2=14387.7 μm·K; T6是卫星高
度上传感器所探测到的像元亮度温度, 单位为K; λ
是有效作用波长, TM6 有效作用波长为11.457μm;
w为大气水分含量
在此我们也取w=取w=1.0

2.2操作步骤
研究区示意图

2.2.1计算L6
在单窗算法中我们已经具体说明，在此请参阅基于TM单窗算法反演地表温度文本。
2.2.2T6e6的求算
T6同理。
e6同理。
2.2.3计算R
在Band Math中编写1/[(14387.7*L6/T6^2)*(11.457^4*L6/1.19104×10^8+1/11.457)]
选中b1变量，在Available Bands List中选择L6，
选中b2变量，在Available Bands List中选择T6，
保存文件为R。
2.2.4计算t
在Band Math中编写(b1-b2*b3)
选中b1变量，在Available Bands List中选择T6，
选中b2变量，在Available Bands List中选择R，
选中b3变量，在Available Bands List中选择L6，保存文件名为Ts ，点击OK。

2.3实验结果
2.3.1温度反演灰度图像

操作过程中许多未知参数大都取参考值，所以可能不能很准确的表示真实地表温度，在日后的学习过程中，需要收集更多的资料做好温度反演。
2.3.2密度分割后的图像

三：辐射方程
3.1实验过程
3.1.1数据准备
热红外数据使用的是TM6，我们需要已经经过大气校正，几何校正，传感器定标的剪裁好的TM影像，以及由次影像生成的NDVI图像。

3.1.2地表比辐射率的估计
典型地表类型的比辐射率
ew=0.995
ev=0.986
em=0.970
Pv=[(NDVI- NDVIs)/(NDVIv- NDVIs)]2
NDVI 为归一化植被指数, 取NDVIv=0.70 和 NDVIs=0（分别取自5%及95%数据）
e surface=0.9625+0.0614Pv-0.0461Pv^2
具体操作步骤与其他两种算法一直，在Band Math中编写即可。

3.1.3计算同温度下黑体的辐射亮度值
Lx=[e6B(Ts)+(1-e6)I6~]t6+I6_
其中e6是地表辐射率。
B(Ts)为普朗克定律推导得到的黑体在Ts的热辐射亮度。
I6~和I6_分别是大气在TM6波段区间内的向上和向下热辐射强度。
Ts是地表真是温度。
t6为透射度。

则B(Ts)=[Lx-I6~-t6(1-e6)I6_]/t6e6

上式中所需要的参数可以通过各种数学模型进行计算，也可以在NASA官网(http://atmcorr.gsfc.nasa.gov/)进行查询，本次试验中，我们取t6=0.89422，I6~=3.39W/(m2*sr*um)，I6_=5.12W/(m2*sr*um)。

在Band Math中编写
(b2-3.39-0.89422*(1-b1)*5.12)/(0.89422*b1)
其中b1为比辐射率e6。
b2表示红外波段大气校正后的辐射定标值。

注:辐射定标
1. 由于ENVI 4.5 中有专门进行辐射定标的模块，因此实际的操作十分简单。将原始TM 影像打开以后，选择
Basic Tools–Preprocessing–Calibration Utilities–Landsat TM
2. 进入下一步参数选择：根据传感器类型选择Landsat 4,5 或者7。从遥感影像的头文件中获取Data Acquisition 的时间，Sun elevation。如果你是用File–Open External File–Landsat–Fast 的方法打开header.dat（头文件）
的话，sun elevation 就已经填好了。这里Calibration Type 注意选择为Radiance。输出文件，定标就完成了。

3.1.4反演地表温度
在获取温度为Ts的黑体在热红外波段的辐射亮度后，根据普朗克公式的反函数，求得地表真是温度Ts。

Ts=K2/ln(K1/B(Ts)+1)
其中：K1=607.76，K2=1260.56（覃志豪，用陆地卫星TM6数据演算地表温度的单窗算法）；
在Band Math下进行编写。
1260.56/alog(607.76/b1+1)-273
其中b1为温度Ts下的黑体辐射亮度值。

3.2温度反演结果