1、引　言

地表温度LST(Land Surface Temperature)作为区域和全球尺度地表物理过程的关键参数，同时是表征陆地表面能量平衡和气候变化的重要指标(Wan和Dozier，1996；Zhou 等，2010)。地表温度的准确获取对地气能量、地球水汽平衡、城市热岛效应以及生态系统研究具有重要的意义(Gillespie 等，1998；Huang 等，2008；Islam 等，2017)。地表发射率LSE(Land Surface Emissivity)是地物固有的物理特性，是地—气能量平衡和矿物填图的关键参数(Li 等，2013；Sobrino 等，2008)；同时地表发射率也是地表温度反演中一个重要输入变量，在不同的大气条件下，1%的地表发射率误差能够引起0.3—0.7 K的地表温度误差，因此地表发射率的准确反演也是高精度地表温度反演的必要条件(Dash 等，2002)。地表的非均质特征导致地表温度随着时空变化具有较大的差异，因此利用地面站点无法获取大区域地表温度。随着热红外遥感的不断发展，以及LST反演算法的不断成熟，卫星遥感成为获取区域和全球尺度LST时空分布唯一可能的手段。不考虑夜晚和多角度观测，地表温度反演主要采用3类算法，分别是单通道算法SC(Single-Channel methods) (Jiménez-Muñoz和Sobrino，2003；Qin 等，2001；Li 等，2010)、劈窗算法SW(Split-Window algorithm) (McMillin，1975；Sun和Pinker，2007；Wan和Dozier，1996；Du 等，2015；Liu 等，2017)和温度与发射率分离算法TES(Temperature Emissivity Separation method) (Gillespie 等，1998；Wang 等，2011)。

中国高分五号卫星(GF-5)装载了6 台遥感仪器，能够有效的进行大气污染气体、温室气体、环境质量以及生态环境等遥感监测(范斌 等，2017)。全谱段光谱成像仪VIMI(Visual and Infrared Mutispectral Imager)是搭载在GF-5卫星的传感器之一，光谱范围覆盖可见光到热红外区域，共12个通道；其中可见光多光谱和短波红外通道6个，空间分辨率为20 m；中红外通道2个，空间分辨率40 m；热红外通道4个，空间分辨率40 m，中心波长分别为8.20 (8.01—8.39 μm)，8.63 (8.42—8.83 μm)，10.80 (10.30—11.30 μm)和11.95 (11.40—12.50 μm) (Ren 等，2018)。鉴于GF-5 VIMI具有多个热红外波段，可以采用SC算法、SW算法以及TES算法进行温度反演。SC算法是在LSE和大气廓线已知的情况下，基于辐射传输方程反演LST(Hook 等，1992；Meng 等，2016)，该算法在高水汽含量下精度较差(Sobrino 等，2004)。SW算法是利用两个相邻且具有不同的水汽吸收特性的热红外通道来消除大气的影响。Ye等(2017)提出了针对GF-5数据的SW算法，但利用SW算法进行LST反演仍需要LSE，并且SW算法对LSE误差较为敏感(Li 等，2014)。与SC和SW算法不同，TES算法能够同时反演地表温度和发射率，该方法通过建立发射率光谱对比MMD(Maximum-Minimum Difference)与最小发射率的经验关系增加约束条件，实现LST与LSE的同时反演(Gillespie 等，1998)。在精确的大气校正前提下，TES算法反演LST和LSE的理论精度分别为±1.5 K和±1.5%(Gillespie 等，1999)。然而，一些研究指出对于低发射率光谱对比的地物(水体、冰雪、植被)及高水汽含量的大气条件下，TES算法反演的LST&E具有较大的误差(Coll 等，2007；Gillespie 等，2011；Hulley 等，2009；Hulley和Hook，2011)。

大气校正的残差是TES算法反演LST&E的主要误差来源，在高水汽含量的大气条件下能够对ASTER和MODIS LST分别造成3 K和4.5 K的误差(Hulley和Hook，2012；Hulley 等，2012)。为了改进热红外波段大气校正的精度，Tonooka(2001)提出了WVS (Water Vapor Scaling)大气校正方法，该方法通过水汽缩放因子γ对大气参数进行修正从而提高大气校正的精度。Tonooka(2005a)以ASTER数据为例，对比了标准大气校正算法和WVS算法的精度，实验结果表明WVS大气校正方法能够提供更高精度的LST&E反演结果，LST的RMSE期望值约为0.6 K。之后WVS方法被进一步应用到了MODIS数据(Hulley和Hook，2011；Tonooka，2005b)。Malakar和Hulley(2016)针对MODIS数据提出了考虑传感器角度的改进WVS方法，该方法与TES算法结合被用来生产MOD21地表温度和发射率产品，验证结果表明改进的WVS方法有望使LST&E反演精度达到1 K和1%。同样该方法也被应用到VIIRS数据进行生产LST&E产品(VNP21)(Islam 等，2017)。

本文提出了一个适用于GF-5 VIMI数据的地表温度和发射率分离算法，为了提高大气校正精度，首先基于Seebor V5.0全球大气廓线数据库提出适用于GF-5 VIMI数据的WVS大气校正方法，再利用TES算法进行温度和发射率分离(简称：WVSTES算法)。首先利用模拟数据对算法精度进行验证，然后利用黑河流域生态—水文过程综合遥感观测联合试验中获取的ASTER数据作为替代数据和同步的实测地表温度和发射率对算法进行精度验证。

2、温度与发射率分离算法     (2.1) 热红外辐射传输方程

在热红外波段(8—14 μm)传感器接收到的地表辐射主要包括3部分：地表发射辐射、地表反射的大气下行辐射以及大气上行辐射。因此，传感器入瞳处的辐射亮度值 ${L_i}$ 可以表示为

${L_i} = {\varepsilon _i}{B_i}({T_{\rm{s}}}){\tau _i} + (1 - {\varepsilon _i})R_i^ \downarrow {\tau _i} + R_i^ \uparrow $

(1)

式中， ${\varepsilon _i}$ 为第 $i$ 波段的地表发射率； ${B_i}({T_{\rm{s}}})$ 为Planck公式计算的黑体辐射， ${T_{\rm{s}}}$ 为地表温度； ${\tau _i}$ 为 $i$ 波段大气透过率； $R_i^ \downarrow $ 和 $R_i^ \uparrow $ 分别为 $i$ 波段大气下行和上行辐射。式(1)中与地表相关的辐射项即离地辐射可以表示为

${L_{{\rm{s}},i}} = {\varepsilon _i}{B_i}({T_{\rm{s}}}) + (1 - {\varepsilon _i})R_i^ \downarrow = \frac{{{L_i} - R_i^ \uparrow }}{{{\tau _i}}}$

(2)

因此，利用式(1)和式(2)进行LST&E反演时，首先需要利用辐射传输模型计算大气参数( ${\tau _i}$ ， $R_i^ \downarrow $ 和 $R_i^ \uparrow $ )，进行大气校正。

    (2.2) 大气校正

地表辐射经过与大气的相互作用(大气反射，衰减以及散射)最终到达传感器。因此，求解地表辐射时需要进行大气校正。传统基于辐射传输模型的大气校正方法通过输入与卫星数据过境时间匹配的大气廓线进行大气参数的计算。目前，常用的辐射传输模型有MODTRAN(Kneizys 等，1996)和RTTOV(Hocking 等，2015；Matricardi 等，2004)；Meng等(2016)利用TIGR3大气廓线对比了MODTRAN和RTTOV两种辐射传输模型的精度，结果表明二者计算的大气参数精度基本一致，但RTTOV具有更快的运行速度。常用的大气廓线数据主要通过两种方式获取：(1)利用再分析数据，通过大气再分析数据进行时间和空间插值模拟实时大气状况，目前常用的再分析数据主要有NOAA的NCEP(National Centers for Environmental Prediction)、NASA的MERRA(Modern Era Retrospective-analysis for Research and Applications)和欧洲中期天气预报中心ECMWF(European Centre for Medium-Range Weather Forecasts)再分析数据；(2)利用遥感反演的大气廓线产品，例如MODIS MOD07和AIRS大气廓线产品(Li 等，2013)。由于再分析数据空间上连续并且时间分辨率较高，本文采用MERRA再分析数据来计算大气透过率和大气上、下行辐射，同时为了提高运行效率，对图像数据的大气校正采用快速辐射传输模型RTTOV。

    (2.3) WVS方法

热红外波段大气校正的残差是TES算法反演LST&E误差的主要来源，而传统的大气校正算法在高水汽含量时具有较大误差，这种不确定性会导致更大的LST&E反演误差(Hulley 等，2009；Li 等，2014；Sabol 等，2009)。因此，本文将WVS大气校正方法应用于GF-5 VIMI数据，以提高大气校正精度。

WVS大气校正方法是由WVD(Water Vapor Dependent)算法发展而来(Francois和Ottle，1996)，地表亮温可以表示为水汽和大气层顶亮温的二次多项式(Francois和Ottle，1996；Malakar和Hulley，2016；Tonooka，2005a)

${T_{{\rm{g}},i}} = {a_{i,0}} + \sum\limits_{k = 1}^n {{a_{i,k}}{T_k}} $

(3)

${a_{i,k}}{\rm{ = }}{p_{i,k}} + {q_{i,k}}W + {r_{i,k}}{W^2}$

(4)

式中， ${T_{{\rm{g}},i}}$ 为第 $i$ 波段的地表亮温； $p$ 、 $q$ 和 $r$ 为系数； $W$ 为大气水汽含量，单位为g/cm²； ${T_k}$ 为第 $k$ 波段大气层顶亮温。

WVS方法假设大气透过率 ${\tau _i}$ 可以表示为Pierluissi双指数波段模型(Pierluissi 等，1985)，进而推导出水汽缩放因子 $\gamma $ 。首先利用MODTRAN 5.2模拟不同水汽缩放因子( ${\gamma_1}$ =1.0、 ${\gamma _2}$ =0.7)条件下的大气透过率 ${\tau _i}({\gamma _1})$ 和 ${\tau _i}({\gamma _2})$ ，然后利用 ${T_{{\rm{g}},i}}$ 、 ${\tau _i}({\gamma _1})$ 和 ${\tau _i}({\gamma _2})$ 计算灰体像元水汽缩放因子 $\gamma $

${\gamma ^i}{\rm{ = }}\frac{{\ln \left(\frac{{{\tau _i}{{(\theta,{\gamma _2})}^{\gamma _1^{{\beta _{^i}}}}}}}{{{\tau _i}{{(\theta,{\gamma _1})}^{\gamma _2^{{\beta _{^i}}}}}}} \cdot {{\left(\frac{{{B_i}({T_{{\rm{g}},i}}) - L_i^ \uparrow (\theta,{\gamma _1})/(1 - {\tau _i}(\theta,{\gamma _1}))}}{{{L_i} - L_i^ \uparrow (\theta,{\gamma _1})/(1 - {\tau _i}(\theta,{\gamma _1})}}\right)}^{\gamma _1^{{\beta _{^i}}} - \gamma _2^{{\beta _{^i}}}}}\right)}}{{\ln ({\tau _i}(\theta,{\gamma _2})/{\tau _i}(\theta,{\gamma _1}))}}$

(5)

式中， ${\beta _i}$ 为波段模型参数； ${\tau _i}(\theta,{\gamma _1})$ 和 ${\tau _i}(\theta,{\gamma _2})$ 分别为 ${\gamma _1}$ 和 ${\gamma _2}$ 水汽缩放条件下的大气透过率； $L_i^ \uparrow (\theta,{\gamma _1})$ 为大气上行辐射； ${B_i}({T_{{\rm{g}},i}})$ 和 ${L_i}$ 分别为地表和大气层顶辐亮度，θ为卫星观测天顶角，对于GF-5 VIMI数据设置为0。通过式(5)可以获取灰体像元的水汽缩放因子，对于非灰体像元采用反距离加权法进行插值得到。

本文利用全球大气廓线库Seebor V5.0模拟WVS系数，其中包括(1) WVS回归系数p、q和r；(2) β_i波段模型参数BMP(Band Model Parameters)；(3) 大气下行辐射回归系数。Seebor V5.0大气廓线库共包含了15704条全球范围均匀分布的晴空条件下的大气廓线数据；每条廓线记录了101个大气压强层的温度，湿度和臭氧含量等(Borbas 等，2005)。由于无线探空仪可能向水体区域偏移，因此，以陆地面积>50%为条件，共得到9136条陆地上空的大气廓线；此外又利用日出和日落时间，对廓线数据进行二次筛选，最终分别获取了4948条白天和4188条夜晚的陆地上空的廓线数据。

(1) WVS系数模拟。本文利用MODTRAN 5.2辐射传输模型，以筛选的大气廓线作为输入，模拟了GF-5传感器大气层顶亮温 ${T_k}$ 和地表辐射亮温 ${T_{\rm{g}}}$ 。然后利用式(3)和式(4)基于最小二乘方法对WVS系数进行模拟。表1为白天WVS系数的模拟结果。

(2) BMP模拟。与WVS系数模拟相似，BMP模拟时利用MODTRAN 5.2分别获取了不同水汽缩放因子γ(分别为0.7、0.9和1.0)下的大气透过率，利用经验模型对参数进行最小二乘回归

${\tau _i}(\theta,\gamma) = {\tau _i}{(\theta,{\gamma _1})^{\frac{{{\gamma ^{{\beta _i}}} - {\gamma _2}^{{\beta _i}}}}{{{\gamma _1}^{{\beta _i}} - {\gamma _2}^{{\beta _i}}}}}} \cdot {\tau _i}{(\theta,{\gamma _2})^{\frac{{{\gamma _1}^{{\beta _i}} - {\gamma ^{{\beta _i}}}}}{{{\gamma _1}^{{\beta _i}} - {\gamma _2}^{{\beta _i}}}}}}$

(6)

式中， $\gamma $ 、 ${\gamma _1}$ 和 ${\gamma _2}$ 分别为0.9、1.0和0.7； ${\tau _i}(\theta,\gamma)$ 、 ${\tau _i}(\theta,{\gamma _1})$ 和 ${\tau _i}(\theta,{\gamma _2})$ 分别为第 $i$ 波段不同 $\gamma $ 时的大气透过率； ${\beta _i}$ 为第 $i$ 波段的BMP系数。模拟结果如表2所示。

(3) 大气下行辐射回归系数模拟。为了减少辐射传输模型的运行次数和时间，大气下行辐射可以表示为天底观测方向大气上行辐射的非线性函数

$L_i^ \downarrow (\gamma) = {a_i} + {b_i} \cdot L_{_i}^ \uparrow (0,\gamma) + {c_i} \cdot L_{_i}^ \uparrow {(0,\gamma)^2}$

(7)

式中， ${a_i}$ 、 ${b_i}$ 和 ${c_i}$ 为回归系数； $L_i^ \downarrow (\gamma)$ 为大气下行辐射； $L_{_i}^ \uparrow (0,\gamma)$ 为天底观测方向的大气上行辐射。

$L_i^ \uparrow (0,\gamma) = L_i^ \uparrow (\theta,\gamma) \cdot \frac{{1 - {\tau _i}{{(\theta,\gamma)}^{\cos \theta }}}}{{1 - {\tau _i}(\theta,\gamma)}}$

(8)

式中， $L_i^ \uparrow (\theta,\gamma)$ 为观测天顶角 $\theta $ 下的大气上行辐射。为了提高模拟的效率，本文利用MODTRAN 5.2模拟大气下行辐射时直接采用53度来替代半球积分的结果(Kondratyev，1969)。

    (2.4) 大气参数调整

对于非灰体像元 $\gamma $ 值采用反距离加权法进行插值，并对结果进行平滑处理。然后利用公式(6)对大气透过率进行调整，而大气上行辐射与经过 $\gamma $ 校正的透过率有以下关系

$L_i^ \uparrow (\theta,\gamma) = L_i^ \uparrow (\theta,{\gamma _1}) \cdot \frac{{1 - {\tau _i}(\theta,\gamma)}}{{1 - {\tau _i}(\theta,{\gamma _1})}}$

(9)

式中， ${\tau _i}(\theta,\gamma)$ 为经过 $\gamma $ 校正的大气透过率； $L_i^ \uparrow (\theta,{\gamma _1})$ 和 ${\tau _i}(\theta,{\gamma _1})$ 为水汽缩放因子为 ${\gamma _1}$ 时大气上行辐射和透过率。

由于GF-5 VIMI幅宽较小，可近似看作垂直观测；因此在进行大气下行辐射调整时可以直接利用式(7)进行计算，表3为回归系数。

    (2.5) TES算法

Gillespie等(1998) 提出了针对 ASTER 数据的温度与发射率分离算法。该算法综合了NEM (Normalized Emissivity Method)、RATIO (Ratio algorithm)以及MMD算法的优点。NEM算法主要通过式(10)—式(12)迭代计算地表温度T和发射率的初始值

$R_i'{\rm{ = }}L_{{\rm{s}},i}' - (1 - {\varepsilon _{\max }})L_i^ \downarrow $

(10)

${T_{\rm{NEM}}} = \max (B_i^{ - 1}(R_i'))$

(11)

${\varepsilon '} = R_i'/{B_i}({T_{\rm{NEM}}})$

(12)

式中， $R_i'$ 为通道i的地表发射辐射； $L_i^ \downarrow $ 为大气下行辐射； ${B_i}$ 和 $B_i^{ - 1}$ 分别为正向和逆向Planck公式； ${\varepsilon _{\max }}$ 为给定的最大发射率，对于植被、冰雪和水体 ${\varepsilon _{\max }}$ 约为0.99，而沙石相对于植被等物体具有较低的 ${\varepsilon _{\max }}$ ，约为0.96。由于TES算法迭代计算时可以对 ${\varepsilon _{\max }}$ 进行调整，因此一般将 ${\varepsilon _{\max }}$ 初始化为0.99； $L_{{\rm{s}},i}'$ 为离地辐亮度，可以通过下式计算得到

$L_{s,i}' = (L - L_i^ \uparrow)/{\tau _i}$

(13)

式中， $L_i^ \uparrow $ 和 ${\tau _i}$ 分别表示大气上行辐射和大气透过率，可以通过大气校正获得； $L$ 为大气层顶辐射亮度。

RATIO算法利用各通道光谱发射率及其平均值进行比值计算，作为发射率光谱波形的无偏估计，即

${\beta _i}{\rm{ = }}\frac{{{\varepsilon _i}}}{{\bar \varepsilon }}$

(14)

最后利用MMD算法即最大最小发射率之差与最小发射率之间的经验关系确定发射率

${\rm{MMD}} = \max ({\beta _i}) - \min ({\beta _i})$

(15)

${\varepsilon _{{\rm{min}}}}{\rm{ = }}{\alpha _1}{\rm{ - }}{\alpha _2}{{\rm{MMD}}^{{\alpha _3}}}$

(16)

$\varepsilon _i^{\rm{TES}} = {\beta _i}\left(\frac{{{\varepsilon _{{\rm{min}}}}}}{{\min ({\beta _i})}}\right)$

(17)

式中， ${\alpha _i}$ 为MMD系数，可利用发射率波谱库数据模拟得到。

本文利用GF-5 VIMI热红外通道光谱响应函数对81条ASTER发射率波谱数据(包括水体、植被、冰雪、人造地物、土壤和岩石等)进行卷积运算，获取了各热红外波段的通道发射率，然后基于最小发射率回归模型对系数进行模拟，式(18)为GF-5 VIMI热红外通道的MMD与最小发射率的拟合结果；图1为MMD拟合时散点分布图。

${\varepsilon _{{\rm{min}}}}{\rm{ = }}0.9865 - 0.7451 \times {{\rm{MMD}}^{0.8455}}$

(18)

3、试验数据与处理     (3.1) GF-5数据模拟

ASTER传感器与GF-5 VIMI具有相似的波谱分布(图2)，为了评价WVSTES算法的精度，本文选取了2012年5月30日—2012年9月19日黑河流域11景无云ASTER作为替代数据(表4)，利用数值模拟的方法将ASTER 热红外波段的数据转换成GF-5 VIMI热红外波段的数据。

ASTER传感器搭载在Terra卫星，由日本和美国合作研制，于1999年12月18日发射。ASTER波谱覆盖了0.52—11.65 μm，共14个波段，其中包含3个可见光近红外波段，空间分辨率15 m；6个短波红外波段，空间分布率为30 m；以及5个热红外波段，空间分辨率为90 m(Yamaguchi 等，1998)。

由于VIMI和ASTER波段设置具有一定的差异(图2)，本文以大气辐射传输模型MODTRAN 5.2为基础，利用9136条陆地上空的Seebor大气廓线，同时结合ASTER发射率数据，利用辐射传输方程分别模拟了ASTER和VIMI热红外波段对应的大气层顶辐射亮度，建立了VIMI与ASTER各波段的转换关系(Ye 等，2017)，其模拟系数见表5，其中VIMI第11波段采用了ASTER 14和15两个波段。

$L_{_i}^{\rm{GF-5}} = a \times L_j^{\rm{ASTER}} + b \times L_k^{\rm{ASTER}} + c$

(19)

$L_i^{\rm{GF-5}} = a \times L_j^{\rm{ASTER}} + b$

(20)

式中， $L_i^{\rm{GF-5}}$ 为VIMI第 $i$ 波段辐亮度； $L_j^{\rm{ASTER}}$ 为ASTER第 $j$ 波段辐亮度；a、b 和 c 为回归系数。

    (3.2) 地面验证数据          3.2.1. 地表温度数据

“黑河流域生态—水文过程综合遥感观测联合试验”(简称“黑河生态水文遥感试验”，英文名称Heihe Watershed Allied Telemetry Experimental Research，简称HiWATER)是2012年—2015年在黑河流域开展的卫星和航空遥感及地面观测互相配合的多尺度综合观测试验(Li 等，2013)。本文利用了HiWATER 3种地表温度的观测网络数据进行验证，(1) WATERNET无线传感器网络(Jin 等，2012)，(2) MUSOEXE自动气象观测站矩阵(Xu 等，2013)，以及(3) 定标站。共包括52个地表温度测量站点(图3)。

WATERNET是在黑河中游盈科/大满灌区5.5 km×5.5 km观测矩阵(图3黄色图框)布设的获取生态水文关键参数的无线传感器网络，该传感器网络共有50个节点，其中29个下垫面为玉米的节点安装了SI-111红外辐射计(图3右图)；数据观测频率为10 min，当卫星或飞机过境时每分钟观测一次(Jin 等，2012)。MUSOEXE是在黑河流域中游张掖地区构建的两个嵌套的自动气象站(AMS)观测矩阵，利用CNR1/CNR4四分量辐射仪进行观测，观测频率为1 min。张掖绿洲周边的神沙窝(SSW)、花寨子(HZZ)、戈壁(GB)与湿地(SD)4个普通站，观测矩阵大小为30 km×30 km，下垫面涉及绿洲农田、沙漠、荒漠、戈壁和湿地。在绿洲内5.5 km×5.5 km观测区有17个观测站点，利用四分量辐射仪进行测量(Xu 等，2013)。两个定标站主要用于试验期间机载热红外传感器的绝对辐射定标，安装有SI-111红外辐射计；其中低温目标场地位于张掖市正东的二坝水库(EB)，高温目标场地位于张掖机场附近的荒漠(JCHM)。

遥感影像像元代表了一定尺度范围内地表信息的综合，而利用地面站点对遥感数据进行验证时，容易受地表异质性结构的影响。为了对验证站点进行空间均一性评价，本文利用11景ASTER影像统计了各站点周围3×3像元大气层顶亮温的标准差，以标准差小于1 K作为衡量标准(Wan，2014；杨贵军 等，2015)，最终选取了32个验证站点，其中包括13个WATERNET站点，17个MUSOEXE站点以及2个定标站。

对于SI-111观测，使用以下公式计算地表温度

$B({T_{\rm{s}}}) = (B({T_{\rm{r}}}) - (1 - \varepsilon){L_{{\rm{sky}}}})/\varepsilon $

(21)

式中，B为Planck公式， ${T_{\rm{s}}}$ 为地表温度， ${T_{\rm{r}}}$ 为SI-111观测的温度， ${L_{{\rm{sky}}}}$ 为大气下行辐射，利用天顶角为55°观测的SI-111计算得到， $\varepsilon $ 是SI-111通道发射率，利用BOMEM MR304红外光谱仪实测的组分发射率结合基于植被覆盖度的方法计算得到(Li 等，2014；Valor和Caselles, 1996)。

对于CNR1/CNR4观测的地表长波辐射数据，使用式(2)计算地表温度

${T_{\rm{s}}} = {\left(\frac{{{F^ \uparrow } - (1 - {\varepsilon _{\rm{b}}}) \cdot {F^ \downarrow }}}{{{\varepsilon _{\rm{b}}} \cdot \sigma }}\right)^{1/4}}$

(22)

式中， ${T_{\rm{s}}}$ 为地表温度， ${F^ \uparrow }$ 为上行长波辐射， ${F^ \downarrow }$ 为下行长波辐射， ${\varepsilon _{\rm{b}}}$ 为宽波段发射率， $\sigma $ 为Stefan Boltzmann常数(5.67×10⁻⁸ W·m⁻²·K⁻⁴)。本文采用8—13.5 μm的宽波段发射率(Cheng 等，2013)，利用BOMEM MR304红外光谱仪实测的组分发射率结合基于植被覆盖度的方法计算得到(Valor和Caselles, 1996)，其中玉米的宽波段发射率为0.985，裸土的宽波段发射率为0.958，植被覆盖度利用ASTER可见光近红外波段计算得到(Mu 等，2015)。

         3.2.2. 地表发射率数据

野外试验期间利用CIMEL CE312-2热红外辐射仪和漫反射金板同步测量了地表辐射和大气下行辐射，并使用TES算法反演了地表发射率，测量的样点包括EB站的水体、SD站的芦苇和MUSOEXE AMS15站点的玉米。CE312-2包括5个窄波段和1个宽波段，光谱范围覆盖了8—13 μm。由于CE312-2窄波段与ASTER TIR波段基本一致，为了能够利用CE312-2发射率对VIMI反演的发射率进行验证，本文分别利用CE312-2和VIMI TIR波谱响应函数对81条ASTER发射率波谱进行卷积计算，获得与 CE312-2和VIMI各TIR通道对应的通道发射率，然后利用最小二乘方法对发射率进行回归，最后利用回归系数将CE312-2的通道发射率转换成与VIMI通道对应的发射率，回归系数如表6所示。

    (3.3) MERRA数据

MERRA再分析数据是美国航空航天局通过地球观测系统的资料同化系统提供的再分析资料，致力于天气和气候时间尺度与水循环相关研究。本文使用的MERRA数据空间分辨率为2/3°×1/2°，垂直方向分为42个压强层，共记录了包括地面压强，地面温度以及臭氧混合比等多个数据集(Rienecker 等，2011)，其时间分辨率为6 h，本文根据影像成像时间对MEERA数据进行线性插值，得到卫星过境时刻的大气廓线数据。

4、结果与分析

为了验证WVSTES算法的精度，本文利用两种方法对算法进行了验证。(1)利用含误差的Seebor廓线对算法进行验证，由于TES算法对低水汽含量敏感性较低，因此本文将20%水汽误差加入廓线模拟大气校正的残差(Jiménez-Muñoz等，2014)。(2)利用黑河流域11景ASTER数据作为替代数据，分别利用TES和WVSTES算法反演了LST&E，同时将反演结果与地表温度和发射率实测数据进行对比。

    (4.1) 利用模拟数据验证

为评价WVSTES算法对VIMI的适应性，首先利用无水汽误差的大气廓线及81条ASTER发射率数据基于辐射传输方程模拟了大气层顶辐射亮度；然后利用加入20%水汽误差的大气廓线获取含误差的大气参数，用于模拟LST&E反演时大气校正的残差；最后分别利用TES和WVSTES算法进行LST&E反演，验证时计算了偏差(Bias)、标准差(STD)以及均方根误差RMSE(Root-Mean-Square-Error)，表7为TES和WVSTES算法的验证结果。

从表7可以看出，TES算法对地表温度高估了1.45 K，严重低估了地表发射率，特别是前3个波段；相比于TES算法，WVSTES算法LST&E总体精度有较大的提高，LST偏差由1.45 K降低到−0.52 K；STD由2.15 K降低到1.45 K；均方根误差由2.59 K减小到1.54 K；同时，各波段LSE误差也得到了不同程度的降低(0.01—0.08)，特别是前3个波段RMSE分别由0.122、0.120和0.102降低到0.042、0.040和0.028。图4(a)和图4(b)分别为TES和WVSTES算法LST误差(ΔLST)的直方图，从图中可以看出，当增加20%水汽误差时，TES算法对地表温度存在明显高估；这是因为当水汽含量增大时，大气透过率降低，为了保持大气层顶辐射不变，离地辐射亮度增加，进而导致地表温度升高。而对于WVSTES算法而言(图4(b))，ΔLST直方图符合均值为0的正态分布且误差基本小于2 K，说明该算法能够降低大气校正的误差，提高LST&E反演精度。

图5为两种算法4个热红外波段LSE的误差直方图。从图5中可以看出TES算法受大气校正误差的影响更为严重，由于高估温度，导致发射率被低估，前3个波段的绝对Bias都大于0.04，RMSE除第4波段外均大于0.1。而WVSTES误差直方图呈对称分布，无明显高估或低估现象，并且Bias较小，均小于0.01。与TES算法相似，WVSTES算法前两个波段RMSE比后两个波段更大，主要是由于前两个波段位于8—9 μm，更容易受大气的影响，且前两个波段的地表发射率动态范围也更大。

本文还利用模拟数据对不同地表类型下TES和WVSTES算法的精度进行了验证，选取的地物类型包括灰体(水体、冰雪、植被)、人造地物以及土壤。表8为TES算法地表温度和各通道发射率的Bias和RMSE，表9为WVSTES算法地表温度和各通道发射率的Bias和RMSE。对比表8和表9可以看出WVSTES算法反演的地表温度具有更高的精度，尤其是灰体区域，其RMSE从TES算法的2.53 K下降到WVSTES算法的0.87 K；不同地表类型下LSE RMSE精度约提高了0.01—0.09，特别是前3个波段提高了0.08左右。此外，从灰体模拟结果来看，WVSTES算法能够较好的改进TES算法在高水汽含量下灰体区域精度较差的情况，反演得到更高精度的LST&E。

    (4.2) 地面站点验证          4.2.1. 地表温度验证

本文通过将卫星过境时间与地表实测数据时间相匹配，提取了HiWATER各站点对应的LST反演结果，共284个验证数据，同时分别计算了Bias、STD以及RMSE对反演结果进行定量化描述。对比发现部分反演结果与实测值之间存在较大的误差，因此为了获取可靠的验证结果，在验证时利用阈值对异常值进行剔除，本文采用3倍算法理论误差作为异常值去除的条件(Sobrino和Julien，2016)。对于TES算法，地表温度的理论误差是1.5 K，因此本文用4.5 K作为异常值的剔除条件，最终有277个点用于LST结果的验证。图6为反演LST与地面实测数据的散点图，其中红色为WVSTES算法反演结果，蓝色为TES算法反演结果。由整体统计结果可知，TES算法高估了LST，Bias为1.08 K，STD和RMSE分别为1.89 K和2.17 K；而WVSTES算法具有更高的LST反演精度，Bias、STD和RMSE都更小，其中Bias为0.47 K，与TES算法相比降低了0.61 K；STD和RMSE分别为1.64 K和1.7 K，分别减小了0.25 K和0.47 K。WVSTES算法反演结果与实测值的决定系数 ${R^2}$ 达到了0.96，而TES算法为0.95。

表10是两种算法11天数据的验证结果。从表中可以看出TES算法11天LST结果的Bias介于−1.06—3.01 K，其中5天的Bias大于1 K；而WVSTES算法的Bias介于−1.01—1.07 K，基本都在±1 K以内，说明WVSTES具有更高的精度。TES算法11天的RMSE介于1.2—3.5 K之间，WVSTES算法11天的RMSE介于1.1—2.1 K之间。TES算法Bias和RMSE的最大值出现在8月11号，分别为3.01 K和3.47 K，这是由于高水汽含量(3.05 g/cm²)导致的结果，而WVSTES算法的Bias和RMSE分别为−0.12 K和1.73 K，这说明WVSTES算法能在高水汽含量下具有更高的反演精度，这与Tonooka (2005a)和Hulley等(2012)的结果一致。

         4.2.2. 地表发射率验证

为了验证GF-5 WVSTES算法LSE反演结果，将WVSTES和TES算法反演的LSE与CE312-2测量的发射率进行对比。选取了6月24日水体(EB)、芦苇(SD)和玉米(AMS15)作为LSE验证站点，图7为WVSTES, TES以及实测LSE验证结果。从图中可以看出，TES算法的LSE结果低估严重，特别是前两个通道，这与TES算法高估LST的结果一致。WVSTES算法反演的LSE发射率波谱形状与实测值基本保持一致，二者相差结果较小，而TES算法的发射率波谱有较大偏差，与实际情况不符。图7中直方图代表以上验证点分别利用WVSTES和TES算法反演LSE与实测值之间的差异，其中黑色为WVSTES算法结果，灰色为TES算法结果。从图中可以看出，TES算法受大气校正残差的影响较大，LSE反演结果与实测值之间误差均大于0.01，其中各波段最小偏差分别为0.048、0.030、0.011和0.016；而WVSTES算法能够提高LSE反演的精度，误差能够达到1%的理论精度。

5、结　论

本文提出了一个适用于GF-5 VIMI数据的地表温度和发射率分离算法，为了提高大气校正精度，该算法加入了WVS大气校正方法。该算法首先利用快速辐射传输模型RTTOV和MERRA再分析数据进行大气校正，再利用WVS方法改进大气校正结果，最后利用TES算法同时反演地表温度和发射率。为了评价算法精度，选取了Seebor V5.0大气廓线库中9136条陆地上空的晴空大气廓线进行了模拟验证。首先利用真实的大气廓线模拟了大气层顶辐射亮度，然后对大气廓线数据加入20%的水汽误差以模拟大气校正的残差，再利用模拟数据分别对TES算法和WVSTES算法进行精度验证。验证结果表明加入WVS方法后TES算法反演地表温度的RMSE由2.59 K降低到1.54 K，4个波段地表发射率RMSE分别由0.122、0.12、0.10和0.037降低到0.042、0.04、0.028和0.026；此外，不同地表类型下TES和WVSTES算法精度也进行了模拟，验证结果表明加入WVS算法后不同地表类型地表温度RMSE降低了约1.5 K，各通道发射率RMSE降低了0.01—0.09，特别是前3个通道约0.08；模拟结果说明WVSTES算法能够精确进行大气校正，反演更加准确的地表温度和发射率。

此外，本文还利用黑河流域的ASTER数据模拟了GF-5 VIMI影像，并利用同步的实测地表温度和发射率数据对算法精度进行了验证。结果表明TES算法高估了地表温度，Bias和RMSE分别为1.08 K和2.17 K；WVSTES算法反演结果与实测LST有更好的一致性，Bias和RMSE分别为0.47 K和1.7 K，相比TES算法分别降低了0.61 K和0.47 K。同时利用CE312-2实测的地表发射率数据对LSE反演结果进行了验证，结果表明WVSTES算法能够明显降低大气校正误差对结果的影响，各波段发射率RMSE均小于0.01。以上结果表明，WVSTES算法能够取得较高的地表温度和发射率反演精度，适用于GF-5 VIMI的地表温度和发射率反演。本文在利用ASTER数据模拟GF-5 VIMI的数据时只进行了光谱转换，并没有考虑两种数据空间分辨率的差异，因此今后将采用更多真实GF-5数据对算法进行测试和验证。

本文在进行地表温度和发射率反演时仅计算了灰体区域的水汽调节因子，而对于非灰体区域的水汽调节因子通过反距离加权法插值和平滑得到可能会产生一定的误差；下一步工作将对非灰体区域水汽调节因子的计算进行改进，实现更高精度的大气校正。

志 谢 此次实验的MERRA数据从NASA GES DISC下载，ASTER数据从Earthdata(https://earthdata.nasa.gov/)下载，地面验证数据来自于HiWATER实验，GF-5 RTTOV系数由James hocking提供，在此表示衷心的感谢！