李 雪,钟仕全
(1.贵州省黔西南自治州气象局,贵州 兴义 562400;2.广西壮族自治区气象减灾研究所/国家卫星气象中心遥感应用试验基地,广西 南宁 530022)
2008年9月6日,我国环境与灾害监测预报小卫星星座A、B 两颗卫星(简称HJ-1A/B 星)发射升空。该星具有重访能力强、分辨率高等特点,在许多方面都得到了广泛的应用。此外,HJ-1B 星搭载了一台红外相机(IRS),其第4 通道光谱范围为10.5~12.5 μm,该通道对热特性敏感,可用来记录地表的发热特性,星下点空间分辨率为300 m,重访周期为4 d,这为地表温度的反演提供新的数据源[3]。
目前已经发展了多种地表温度的遥感反演方法,考虑到HJ-1B 热红外波段范围与TM 非常相近,因此,本文借鉴TM的地表温度反演算法来对环境卫星热红外通道进行反演,并根据HJ-1B 热红外波段的光谱响应特性来修订算法中的一些经验关系,最后利用MODIS 温度产品(MOD11_L2)进行反演精度的验证,以获得适用于岩溶地貌区的HJ卫星遥感数据反演地表温度的算法。
岩溶区多分布于广西、云南、贵州等省,本文根据广西的地质地貌特征,选择广西中部干旱农业区中的典型岩溶地貌区[4]。该区为亚热带季风气候,农业人口众多,占总人口75.12%,以农业耕作为主,耕地面积为601 620 hm2,其中旱地占耕地面积的46.57%。因此,干旱是制约该地区农业和社会经济发展的最为重大的问题之一。
本文使用的遥感数据见表1,其中,环境减灾卫星的CCD、热红外遥感影像数据以及MODIS 遥感数据用于对地表温度的反演,MOD11_L2 为MODIS 温度产品数据,用于对反演结果精度的验证。
表1 所用遥感影像数据一览表
采用经过精校正后的TM 影像对HJ-1B 数据进行几何精校正,将几何误差控制在1 个像元内,并对HJ-CCD 数据重采样到300m;利用专业遥感软件对MOD 02 和MOD11_L2 数据进行几何校正,并重采样到300 m,并根据MOD11_L2 数据头文件提供的计算公式获得地表温度影像。最后利用矢量图裁剪出研究区,以进行下一步的处理和应用。
2.2.1 传感器接收到的辐射强度LλHJ-1B 热红外影像是用DN 值来表示的,根据影像所带头文件提供的信息,传感器所接收到的辐射强度和DN值之间有以下关系:
⑱Rose,R.,“What is lesson - drawing?”Journal of Public Policy,1991,11(1),pp.3 ~30.
式中,Lλ为热红外传感器所接收到的辐射强度,DN 为影像的灰度值,b 和g 为常数。
2.2.2 亮温T 亮温是遥感器在卫星高度所观测到的热辐射强度相对应的黑体的温度,这一温度包含有大气和地表对热辐射传导的影响,因而不是真正意义上的地表温度[5]。根据Plank 辐射函数,由上述得到的辐射强度值Lλ就可以反算出影像的亮温,公式如下:
式中,T 为亮温(K);Lλ为辐射强度值,由公式(1)获得;h 为普朗克常量,h=6.626×10-34J·s-1;c 为光速,c=2.998 ×108m·s-1;k 为玻尔兹曼常量,k=1.38 ×10-23J·K-1;λ 为有效波长,参考刘三超[6]文献,HJ-1B的热红外通道有效波长为11.511 μm。
其中,Lsensor是传感器接收的辐射强度,ε 是地表比辐射率,参数γ 和δ 可由公式(5)计算得到,Ψ1、Ψ2和Ψ3为大气函数,可由公式(6)计算得到。
其中,Tsensor是传感器亮温;λ 是波段的有效波长;c1=1.191 043 56×108Wμm4m-2sr-1,c2=14 387.685 μmK。大气函数(Ψ1,Ψ2,Ψ3)由大气中水汽含量(ω)的函数而获得,针对HJ-1B 热红外波段,参考段四波等[9]文献,可由以下公式计算:
该算法是覃志豪等[5](2001)根据地表热辐射传导方程,通过一系列假设,建立的适用于Landsat TM 第6 波段的反演方法。其表达式为:
其中,a 和b 为常数,针对HJ-1B 热红外波段,根据图像的温度变化范围得到,根据文献[9]可知,对于0~30℃,a=-60.896 9 和b =0.439 078,对于20~50℃,a =-68.330 1 和b =0.464 012;C 和D 为参数,由公式(8)计算得到;Ta为大气平均作用温度,针对HJ-1B 热红外波段,对于中纬度夏季大气,得到Ta的估算方法[9]。
式中ε 是地表比辐射率;τ 是大气透过率,研究表明,大气透过率的变化主要取决于大气水汽含量ω的动态变化,因此,可以根据大气总水汽含量来估算,针对HJ-1B 热红外波段进行修订的大气透过率估算方程(公式(10));T0为近地表温度,可利用研究区域的站点平均气温数据差值得到。
Artis 等人[8](1982)认为,辐射亮温仅仅代表了黑体的温度,而自然界的大部分物体并非黑体,所以应该用比辐射率对其进行校正,获得绝对表面温度,计算如下:
式中,T(K)为传感器的亮温;λ 为有效波长;ε为地表比辐射率;h 为普朗克常量,其值为6.626×10-34J·s;c 为光速,值为2.998 ×108m·s-1;σ 为玻尔兹曼常量,值为1.38 ×10-23J·K-1。
由上述一系列公式,可以看出需要计算的参数有:地表比辐射率ε 和大气水汽含量ω。
地表比辐射率主要取决于地表的物质结构和遥感器的波段区间。计算地表比辐射率的方法很多,本文采用覃志豪[10]等(2004)根据Sobrino 提出的NDVITEM方法结合研究区地表类型的多样性,得出的不同地表类型的比辐射率计算方法。认为地球表面不同区域的地表结构虽然很复杂,但从卫星像元的尺度来看,可以大体视作由3 种类型构成:水体、城镇和自然表面。水体在热波段范围内的比辐射率很高,接近于黑体,可以用εwater=0.995 来进行估计。对于自然表面εsurface和城市地表εbuilt-up,与植被构成比例Pv的关系式如下:
式中,Pv为植被构成比例,其估计方法为[11]:
式中,ρ3、ρ4分别为HJ-1B CCD 数据的第3 和第4 波段的反射率。一般情况下,如果没有详细的区域植被和土壤光谱或图幅上没有明显的完全植被或裸土像元,则用NDVIv=0.70 和NDVIs=0.05来进行Pv的近似估计。根据上述的方法,将CCD遥感数据分为水体、自然表面和城镇3 种类型,然后按公式(12)~(14)计算出地表比辐射率。
在缺少实时大气水汽含量数据时,可以通过同步MODIS 数据反演大气水汽含量。根据GRIEND等[11]文献,得到大气水汽含量的计算公式:
式中,ω 为大气水汽含量;Tω(19/2)为MODIS第19 波段和第2 波段表观反射率之比;α 和β分别为公式参数,对于混合型地表,α =0.02,β =0.651,其他类型地表参考文献[12]。
根据上述3 种算法反演得到研究区的地表温度影像图(图1),能明显看出三者温度变化趋势是相同的,表现在直方图上其走势相似(图2)。其中,研究区的东部和南部区域地表温度较高,温度较低区域主要分布在西北部的河池、环江和融水等县市。
图1 3 种算法反演的地表温度影像图(注:白色区域为云区)
图2 3 种算法反演的研究区地表温度直方图
为了定量地分析各算法的反演结果,由于实时地面测量温度很难获取,验证3 种算法的精度选用MODIS 地表温度产品进行相对验证。选取同天成像时间为11 时30分覆盖研究区的MOD11_L2 温度产品作为验证数据。在低温、中温和高温区中,分别选取对应的森林、耕地和城市3 种地表类型样区进行结果验证。
表2 列出了各检验样区对应的HJ-1B 反演温度值和MODIS 地表温度产品的温度值。从表中可以看出,普适性单通道算法反演结果与MODIS 地表温度产品最为接近,平均温差为0.36 K(覃志豪单窗算法为2.39 K,基于影像的Artis 反演算法为4.69 K),森林温差最大,平均为0.42 K,耕地温差最小,平均为0.26 K。基于影像的Artis 反演算法与MODIS温度产品差异最大,这是因为该方法只考虑了地表比辐射率的影响而忽略了大气的影响。覃志豪单窗算法所得结果与MODIS 温度产品平均温差为2.39 K,主要原因为覃志豪单窗算法是根据MODTRAN 软件4 种标准大气廓线(美国1976 标准大气廓线、低纬度标准大气廓线、中纬度夏季和中纬度冬季标准大气廓线)进行参数拟合,并建立大气平均温度和近地层空气温度的经验关系式;本文采用中纬度夏季标准大气廓线建立关系式,计算的大气平均作用温度比实际大气平均温度高,导致反演得到的温度比实际温度低。
表2 3 种算法的反演结果与MODIS 地表温度产品的温度值统计 (单位:k)
为了便于观察,对每一个检验区分别建立3 种算法反演的地表温度与MODIS 地表温度产品之间的零截距线性方程,见图3:
图3 森林检验区3 种算法反演的温度与MODIS 地表温度产品之间的零截距线性关系
从方程式可以初步看到,对于低温区,普适性单通道算法反演结果优于覃志豪单窗算法和Artis反演算法。为了比较各算法的反演结果,同样对中温区(耕地)和高温区(城市)进行对比,以便更全面地比较各算法的反演结果。得到如下线性方程式:
对于耕地检验区有:
对于城市检验区有:
从方程式(16)~(17)可以看到,普适性单通道算法反演得到的结果优于其他算法,反演得到的地表温度更接近于MODIS 地表温度产品的温度。同时,建立的线性方程显示了温度变化的一个规律,即森林<耕地<城市,这个结果与实际情况相符合。
本文针对HJ-1B 卫星热红外数据的特点,分别用修订后的普适性单通道算法、覃志豪单窗算法和基于影像的Artis 反演算法反演典型岩溶地貌区的地表温度,并与MODIS 温度产品相比较,以验证算法精度。结果表明,修订后的普适性单通道算法反演结果优于其它两种算法,其与MODIS 温度产品平均温差相差0.36 K,保证误差精度在1 K 之内,说明该算法经过修订后适用于反演岩溶地貌区的地表温度,这为利用HJ-1B 遥感数据开展岩溶地貌区的干旱监测提供技术借鉴。但由于研究中大气参数和地表比辐射率只能依靠估算方程计算,这对反演精度会有一定的影响,需要在今后的研究中进行下一步验证。
[1]夏日元,朱远峰,李兆林.广西岩溶区农业发展的资源及地质环境特征[J].广西科学,1997,4(2):192-195.
[2]徐希孺,柳钦火,陈家宜.遥感陆面温度[J].北京大学学报(自然科学版),1998,34(2):248-253.
[3]李艳芳,李小娟,孟丹.环境减灾卫星热红外数据的地表温度反演及LST分布分析——以北京市城八区为例[J].首都师范大学学报(自然科学版),2010,31(3):70-75.
[4]周游游,蒋忠诚,韦珍莲.广西中部喀斯特干旱农业区的干旱程度及干旱成因分析[J].中国岩溶,2003,22(2):144-149.
[5]覃志豪,Zhang M H,Arnon K,et al.用陆地卫星TM6 数据演算地表温度的单窗算法[J].地理学报,2001,56(4):456-466.
[6]刘三超,柳钦火,高懋芳,等.波谱响应函数和波宽对地表温度反演的影响[J].遥感信息,2007(05):3-6.
[7]Jiménez-Mu?oz,J C.,Sobrino,J A.A Generalized Single Channel Method for Retrieving Land Surface Temperature From Remote Sensing Data [J].Journal of Geophysical Research,2003,108(D22):4688.
[8]Artis,D A,Carnahan,W H.Survey of Emissivity Variability in Thermography of Urban Areas[J]Remote Sensing of Environment,1982,12(4):313-329.
[9]段四波,阎广建,钱永刚,等.利用HJ-1B 模拟数据反演地表温度的两种单通道算法[J].自然科学进展,2008,18(9):1001-1008.
[10]覃志豪,李文娟,徐斌,等.陆地卫星TM6 波段范围内地表比辐射率的估计[J].国土资源遥感,2004(03):28-32.
[11]Carlson T N,Ripley D A.On the relation between NDVI,fractional vegetation cover,and leaf area index[J].Remote Sens.Environ.,1997,62:241-252.
[12]GRIEND A,OWE M.On the Relationship Between Thermal Emissivity and the Normalized Difference Vegetation Index for Nature Surfaces[J].International Journal of Remote Sensing,1993,14(6):1119-1131.