基于MODIS的GF-4/PMS遥感器交叉定标——以巴丹吉林沙漠为参考目标

张浩，刘涛，2，闫东川，阎跃观，崔珍珍，5

1.中国科学院空天信息创新研究院,北京 100094;

2.中科星图空间技术有限公司,西安 710100;

3.中国冶金地质总局矿产资源研究院,北京 101300;

4.中国矿业大学（北京）,北京 100083;

5.河南理工大学 测绘与国土信息工程学院,焦作 454003

1 引言

高分四号（GF−4）卫星于2015 年12 月29 日在西昌卫星发射中心发射成功，是中国第一颗地球同步轨道遥感卫星，搭载了一台可见光50 m/中波红外400 m 分辨率、大于400 km 幅宽的凝视相机，采用面阵凝视方式成像，具备可见光、多光谱和红外成像能力，设计寿命8 a，重访周期达20 s，通过指向控制，实现对中国及周边地区的观测。

GF−4 为减灾、林业、地震和气象等各种应用提供快速、可靠、稳定的光学遥感数据，为灾害风险预报、森林火灾监测、地震构造信息增添新的技术方法，在灾害监测等方面具有巨大的潜力和广阔的应用空间（聂娟 等，2018；张磊，2018；吴玮，2019）。遗憾的是GF−4 没有星上定标系统，传感器的定标完全依赖于场地定标，而场地定标受限于场地、设备、成本等因素，无法满足高频率定标需求，这限制了它的应用（高海亮 等，2010；Zhang 等，2018）。交叉定标是利用定标精度较高的传感器来标定待标定的卫星传感器的方法，该方法不需要高成本的野外同步试验，也不需要很精确的大气参数测量；同时交叉定标方法可以对历史数据进行标定，是目前最具有广泛应用前景的定标方法之一（赵维宁 等，2015）。

基于交叉定标的优势，国内外学者先后将交叉定标方法应用于NOAA、BJ−1、HJ−1、GF−1等卫星遥感器。Teillet等（1990）基于美国White Sands场，针 对NOAA−10/AVHRR，利 用Landsat 5 和SPOT/HVR 影像，进行了AVHRR 的交叉定标。Teillet 等（2007）分析了不同地物光谱对交叉定标光谱匹配因子的影响，得出利用Railroad 试验场交叉定标精度优于利用草地交叉定标。Rao 等（2001）基于Sonoran 沙漠试验场，以NOAA−14/AVHRR 为参考传感器，对GOES−8 传感器进行交叉定标。Cao 等（2005）提出SNO（Simultaneous Nadir Overpass）交叉定标方法，该方法利用极地地区的卫星影像对完成了NOAA 系列极轨卫星15、16、17 号搭载的HIRS 传感器之间的相互交叉定标。杨忠东等（2004）使用Landsat 7/ETM+为参考传感器，通过辐射传输模拟和统计分析，实现了对CBERS−01（中巴地球资源卫星）CCD 相机的交叉定标。陈正超等（2008）等在缺少光谱响应函数的情况下，利用SPOT4/HRVIR2、Landst 5/TM和Terra/MODIS 这3 种传感对北京一号小卫星进行交叉辐射定标，得到比较可信的定标系数。马晓红（2011）选取阿拉伯沙漠作为定标试验场，以MODIS 为参考传感器对HJ−1/CCD 相机进行交叉辐射定标，结果表明交叉辐射定标的系数具有较高的精度。Zhong等（2014）和Yang 等（2015，2017）选取巴丹吉林沙漠，基于分辨率较高的Landsat 7/ETM+、Landsat 8/OL 和DEM 数据构建地表BRDF模型，实现对HJ−1/CCD、GF−1/WFV 和GF−4/PMS的交叉辐射定标。张玉环等（2016）基于反射率低、中、高（植被区、沙漠和雪）场景，使用MODIS作为参考传感器对GOCI 进行交叉辐射定标。Chen等（2017）将交叉辐射定标作为最优逼近问题，以Landsat 8/OLI作为参考传感器，采用SCE−UA 算法，通过迭代方程找到最优定标系数和BRDF调整因子，结果表明交叉辐射定标系数计算的地表反射率误差小于5%。

为减少观测角度差异造成的定标误差，通常交叉定标选取两个传感器观测角较小且较为接近的影像进行交叉定标。GF−4/PMS 具有较宽的视角，由于纬度原因在中国的大部分地区都有较大的观测角，这对GF−4/PMS 交叉定标带来较大困难。本文通过影像搜索巴丹吉林沙漠均匀场地作为交叉定标试验场，以MODIS 传感器为参考，利用MODTRAN 辐射传输模型，在约束MODIS 与GF−4 成像角度差异和二者过境时间差异条件下，对GF−4/PMS 时序数据进行交叉定标，获取时序的交叉定标系数，实现对GF−4/PMS 传感器性能的时序监测和评估。

2 研究方法

本文通过较高分辨率影像（Landsat 8/OLI）数据搜索均匀区域对大角度观测的中高分辨率数据（GF−4/PMS）与Terra（Aqua）/MODIS 数据进行交叉定标。首先，使用几何和辐射校准较好的Landsat 8/OLI数据寻找合适的均匀区域作为定标试验场；然后，根据选取的试验场选取MODIS 和GF−4/PMS 时间序列数据；接着，时序数据筛选，除去云和成像质量较差的数据后，利用如下3个条件进一步筛选：（1）550 nm 气溶胶光学厚度小于0.3；（2）MODIS 成像散射角度（观测方向和太阳入射方向夹角）和GF−4 成像散射角度差异小于20°；（3）MODIS 和GF−4 成像时刻小于2 h；最后，对时间序列数据进行预处理，提取试验场时间序列数据的观测几何等信息，使用辐射传输模型计算MODIS和GF−4/PMS模拟的表观辐亮度，计算光谱匹配因子；最后对GF−4/PMS 进行交叉定标，总体流程如图1所示。

图1 交叉定标流程Fig.1 The calculation flow of cross calibration

2.1 定标原理

交叉定标是利用定标精度较高的传感器作为参考，对待定标的传感器进行定标；其原理是选取对同一目标成像的同步或近似同步的影像对，在分析两个传感器光谱响应、观测几何、大气参数等匹配的基础上，建立两个传感器图像数字计数值之间的关系，利用参考传感器已知的辐射定标系数求解待标定传感器的定标系数（高海亮 等，2010；吕文博，2014）。通常情况下，卫星遥感器的DN值与入瞳辐亮度存在线性关系：

式中，gain和offset分别是定标系数的增益和截距，L表示传感器的入瞳辐亮度，DN 表示图像的数字计数值。

通过两个传感器的光谱匹配因子，计算得到待标定传感器的表观辐亮度的模拟值：

2.2 试验场

为了选取合适的实验场，本文使用巴丹吉林沙漠地区无云的Landsat 8/OLI影像，位置如图2所示，搜索DN 值的均值和方差最小的500 m×500 m大小区域作为试验场，最终选取了巴丹吉林沙漠南部边缘的500 m×500 m的均匀区域（图2），图2（a）为试验场的位置，其中右侧部分为Landsat 8/OLI真彩色图像；图2（b）为GF−4/PMS影像上试验场的位置；图2（c）为Landsat 8/OLI 影像上试验场的位置。本文选取目标试验场主要考虑以下几点：（1）巴丹吉林沙漠随时间在亮度、空间均匀性、季节变化、长时期的稳定，且地表主要地物为沙子（Yang等，2017），符合作为交叉定标试验场的条件（Scott 等，1996）；（2）马晓红等使用影像从阿拉伯沙漠选取均匀场地实现了对环境星的交叉定标（马晓红，2011），结果表明使用影像搜索沙漠地区均匀区域作为交叉辐射定标的试验场是可行的且有较高的精度；（3）在此区域能够获取一定数量连续的影像对数据，可实现时序交叉定标。

图2 试验场的位置及试验场影像Fig.2 The location and images of calibration site

2.3 数据

交叉辐射定标的精度依赖于参考传感器的精度，搭载在EOS−Terra/Aqua 上的MODIS 传感器，配备星上太阳辐射校正系统，绝对定标系数不确定度在3%左右（Chang 等，2017）；另外MODIS传感器重访周期短，覆盖范围广、拥有丰富的数据的特点使得其经常被用来作为参考传感器对其他卫星传感器进行交叉定标研究。HJ−1/CCD、CBERS−02/CCD、北京一号小卫星、NOAA−16、GOCI、ETM+等都以MODIS 为参考传感器进行过交叉定标研究（Hu 等，2001；李小英 等，2005；Vermote和Saleous，2006；陈正超 等，2008；马晓红，2011；张玉环 等，2016）。考虑稳定的辐射性能、较高的辐射定标精度和丰富的数据，因此本文选择MODIS作为参考传感器对GF−4/PMS进行交叉辐射定标研究。

2.4 光谱匹配

在对GF−4/PMS 进行交叉定前需要考虑与参考传感器（MODIS）光谱响应函数（SRF）的差异，图3 绘制了GF−4/PMS 和MODIS 可见光和近红外波段的光谱响应函数。其中，MODIS 的光谱响应函数和波段大气层顶太阳辐照度数据来自NASA 官网MCST 团队（https：//mcst.gsfc.nasa.gov/calibration/parameters［2021−12−03］），GF−4/PMS 的光谱响应函数和波段大气层顶太阳辐照度数据来自中国资源卫星中心（http：//www.cresda.com/CN/Downloads/dbcs/［2021−12−03］），二者在可见光和近红外波段的波段信息如表1所示。将影像中提取的试验场观测几何和地表反射光谱等信息输入到辐射传输模型MODTRAN 中，计算出两个传感器模拟的表观辐亮度，得到两个传感器表观辐亮度的光谱匹配因子。为评估地表二向性的影响，也同时利用Ross−Li 模型进行计算了光谱匹配因子（Li 和Strahler，1992），相关输入参数来自于MCD43 产品（Schaaf等，2011），对应如下形式：

表1 GF-4/PMS和MODIS波段范围Table 1 The spectral range of GF-4/PMS and MODIS

图3 GF−4/PMS和MODIS光谱响应曲线Fig.3 The Spectral response function for the corresponding channels of the GF−4/PM and MODIS

式中，fiso(λ)、fvol(λ)与fgeo(λ)分别为各向同性散射、体散射与几何光学散射对应的权重系数，Kvol与Kgeo分别为体散射核与几何光学核。MCD43 A1产品提供了1—7 个波段的fiso(λ)、fvol(λ)与fgeo(λ)，按照MODIS 3、4、1、2 波段与GF−4/PMS 的2—5 波段对应关系从MCD43 A1 提取权重因子，然后按照MODTRAN 运行要求输入对应参数，利用式（3）计算光谱匹配因子。

其中，光谱匹配因子计算使用的地面光谱数据来自（Yang等，2017）2012年7月在巴丹吉林沙漠地区实测，光谱曲线如图4所示；根据（Lacherade等，2013）的研究，在缺少可靠数据的情况下沙漠地区550 nm 气溶胶光学厚度可以采用默认值0.2，其余相关的大气参数使用MODTRAN 模型的默认参数，由参数假设带来光谱配因子的计算误差分析见4.1 节。为尽可能减少误差，这里气溶胶光学厚度数据采用了MAIAC 算法生产的MCD19 产品（Lyapustin等，2021）。

图4 巴丹吉林沙漠地面光谱曲线（Yang等，2017）Fig.4 The spectral curve of Badain Jaran desert（Yang et al.，2017）

3 数据处理

（1）数据。为获取尽量多的数据，选择MODIS 和GF−4 在试验区域同日过境的无云、清晰的影像对，在大气、成像角度和成像时刻约束下选取了2016 年5 月至2018 年9 月过境的13 对GF−4和MODIS 影像作为本次交叉定标的数据，各影像对的成像时间、时间差、试验场位置的成像几何和成像散射角度差异信息如表2所示。

表2 选择的GF-4和MODIS影像对的信息Table 2 Information of selected GF-4/PMS and MODIS image pairs

（2）影像对处理。在选取定标试验场之后GF−4 影像处理步骤：（1）使用ENVI 进行正射校正，裁剪出实验区域2397×1955 大小（对应Landsat 8/OLI：3000×2000 大 小）；（2）使 用Landsat 8/OLI 影像为基准对GF−4/PMS 影像进行几何配准，配准精度优于1 个像元；（3）提取GF−4/PMS 试验场的观测几何信息（观测天顶角、观测方位角、太阳天顶角和太阳方位角）及各波段的DN 值。MODIS 影像的处理步骤：（1）使用ENVI进行几何校正；（2）提取MODIS和GF−4在试验场位置的观测几何以及各波段的DN 值；（3）根据交叉定标区域位置提取MCD19A2 中气溶胶光学厚度、MCD43A1 中BRDF 数值，其中MCD43A1 数值只采用MCD43A2 中质量标记为0 或者1 的数值。将提取的观测及大气参数等信息输入MODTRAN辐射传输模型，计算出GF−4/PMS和MODIS各波段的模拟的表观辐亮度（TOA），求得参考传感器与待标定传感器的光谱匹配因子。利用光谱匹配因子和MODIS 观测得到的表观辐亮度计算出GF−4/PMS 模拟的入瞳辐亮度；将GF−4/PMS 的DN 值与模拟出的入瞳辐亮度进行线性拟合，得到交叉辐射定标系数。为进行BRDF效应的对比分析，分别基于实测光谱的朗伯体假设情况和基于MCD43 产品的BRDF模型情况两种方式计算光谱匹配因子。

4 结果

4.1 实验结果与分析

本文得到的朗伯体假设和BRDF 模型计算的GF−4/PMS 交叉辐射定标系数如表3 所示，中国资源卫星应用中心公布的2016 年—2018 年场地定标系数见表4。为了更直观的分析交叉定标结果，本文绘制了朗伯体假设和BRDF模型两种情况下计算的交叉定标系数和中国资源卫星应用对比散点图（图5），其中虚线代表交叉定标系数随时间的变化趋势。结果表明：（1）由于限定了过境时间差和成像角度差，朗伯体假设和BRDF 模型假设条件下计算的交叉定标系数具有较高一致性，绿、蓝、红和近红外波段二者的相对差异分别为3.86±2.36%、3.36±1.87%、3.88±2.21%和3.89±2.00%；（2）2016 年—2018 年间GF−4/PMS4 各个波段的定标系数整体呈现缓慢增加趋势，表明传感器辐射性能出现下降趋势，朗伯体假设和BRDF模型假设条件下4 个波段的年衰减率分别为1.00%、0.69%、0.43%、0.28%和0.90%、0.39%、0.18%、0.06%；（3）每年8—10月附近交叉定标结果平均值与中国资源卫星应用中心场地定标公布的定标系数的相对误差，2017 年偏差明显（超过9%），这一方面可能与遥感器的辐射稳定性有关，另一方面与单次定标出现的随机误差有关（Chen 等，2014）；（4）相比朗伯体假定条件下计算的交叉定标系数，BRDF 模型假定下计算结果与场地定标系数偏差略高，这可能与后者采用朗伯体假定计算的定标方法有关，该问题需要后续更多场地实测BRDF支持情况下进一步验证。

表3 两种情况的交叉定标结果Table 3 Cross calibration results for two cases

表4 资源卫星中心公布的定标系数Table 4 Radiometric calibration coefficients from CRESDA

图5 交叉定标系数及与中国资源卫星应用中心系数对比（Lambertian代表在朗伯体假定条件下的计算结果；BRDF代表采用Ross−Li模型的计算结果；CRESDA代表中国资源卫星应用中心的结果）Fig.5 The comparison between the cross calibration coefficients and the calibration coefficient of China Center for Resources Satellite Data and Application（Lambertian denotes the results of cross calibration methods；BRDF denotes the results of Ross−Li BRDF model；and CRESDA denotes the results of China Center for Resources Satellite Data and Application）

从图5 中进一步发现2018 年5 月30 日两幅GF−4 分别与MODIS Terra 和Aqua 进行交叉定标，二者在4 个波段差异小于6%，说明本文方法具有较高精度。为进一步检验该结果，我们放宽角度和成像时间差的差异，发现同一天GF−4 对MODIS不同时刻交叉定标的结果保持了较好的一致性，相对误差均小于6%（如表5）。

表5 同一天多组交叉定标系数对比Table 5 Cross calibration coefficients comparison in the same day

4.2 不确定性分析

通过对GF−4/PMS 与Terra（Aqua）/MODIS 的交叉定标过程分析，影响定标精度的因素主要有以下几方面：（a）MODIS 自身定标的不确定性；（b）影像空间配准带来的影响；（c）大气参数误差的影响；（d）地表反射光谱的误差；（e）大气辐射传输模型本身计算的不确定性；（f）地表二向性的不确定性；（g）其他因素引起的误差。

（1）MODIS 定标系数的不确定性为±3%（Chang 等，2017）。巩慧等（2010）用2007 年二连浩特对MODIS 的同步实验数据，采用反射率基法对MODIS 进行了辐射定标和真实性检验研究，证明了MODIS具有较高的定标精度。

（2）为了分析空间像元匹配误差对交叉定标结果的影响，本文将所有选用影像的试验场扩大一倍（1050 m×1050 m）的平均DN 值与试验场的平均DN值做比较最大差异分别为：蓝波段0.28%、绿波段0.24%、红波段0.4%和近红外波段0.19%。因此，即使几何配准误差在1个像元甚至更多像元情况下，DN值的差异也很小。

（3）大气条件对光谱匹配因子的计算的影响主要为：气溶胶类型、气溶胶光学厚度和水汽含量。由于本文交叉定标缺乏同步测量数据，我们使用MODTRAN 中所提供的标准大气成分组成，550 nm 气溶胶光学厚度设置为0.2。在其他变量不变的情况下，分别改变气溶胶类型（沙漠型、乡村型和城市型）、550 nm 气溶胶光学厚度（0.05 至0.4）和水汽含量（0.4 至4.0）的值计算光谱匹配因子的差异，计算结果如表5 所示；图6 列出了2016 年8 月25日不同条件对光谱匹配因子的影响；从表5 和图6（b）中均可看出水汽含量对近红外波段光谱匹配因子的影响较大。

图6 不同条件对光谱匹配因子的影响Fig.6 The impact of difference factors on spectral matching factor

（4）本文所使用的地面光谱为巴丹吉林沙漠中心地带的光谱曲线，并非试验场实测的光谱，为了分析误差我们使用敦煌定标场的地面光谱来计算光谱匹配因子，分析在没有实测光谱的情况下，使用同类地物光谱替代对交叉定标结果所带来的影响；4 个波段计算的光谱匹配因子差异分别为：5.04%、4.02%、0.02%和2.94%，说明使用同类地物光谱计算光谱匹配因子进行交叉定标在可接受的范围之内；本文所使用的验证地物的光谱曲线如图7所示。

图7 替代地物光谱曲线Fig.7 Alternative spectral curves

（5）前面对比了采用朗伯体假定和Ross−Li BRDF 模型假定得到交叉定标差异，蓝、绿、红、近红外波段的相对差异为3.86±2.36%、3.36±1.87%、3.88±2.21%和3.89±2.00%。

综合上述分析，通过对GF−4/PMS与MODIS交叉定标误差来源分析，得到此次交叉定标总的不确定度如表6 所示，总体误差在7.4%以内，高于垂直观测卫星的交叉定标误差估算结果（Chen 等，2014）。

表6 交叉定标结果不确定性分析Table 6 Uncertainty analysis for cross calibration results

5 结论

受到太空环境、遥感器自身光电性能下降等因素的影响，遥感器辐射性能不可避免发生变化，因此，一年一次的场地定标频次难以监测遥感器性能的变化情况，无法满足卫星传感器数据长期定量化的需求。本文以辐射定标精度较高的MODIS 为参考传感器，对2016 年—2018 年GF−4/PMS 可见光近红外波段进行时序定标，为GF−4/PMS 长期辐射性能评估和监测提供一种技术手段。结果表明：

（1）GF−4 在中国区域具有较大观测角度，通过约束与MODIS 成像角度差异和成像时间差异，能较好减少因为大气变化和BRDF 差异带来的影响，忽略和考虑沙地BRDF 效应带来的差异约4%，两种情况下得到的交叉定标系数具有较高的一致性，定标不确定度小于7.4%。

（2）与中国资源卫星中心的场地定标系数相比，2016 年和2018 年交叉定标结果与场地定标系数较为接近，但2017 年偏离较大，这可能与场地定标朗伯体假定、交叉定标自身不确定度有关。

（3）高频次的交叉定标能够有效检查单次定标可能存在的误差，通过交叉定标和场地定标结果相互检查能够剔除存在的粗差点，提高结果的可靠性。

（4）通过时序的交叉定标，我们发现GF−4/PMS 各个波段的辐射性能从2016 年—2018 年呈现缓慢下降趋势，年衰减率小于1%。

由于GF−4 卫星静止轨道特征，与极轨卫星有更多交叉成像机会，在约束成像角度和时间情况下，继续监测GF−4的辐射性能变化情况。

志 谢本文所使用的巴丹吉林沙漠的地表实测光谱曲线由中国科学院空天信息研究院的杨爱霞博士提供，在此表示衷心的感谢；感谢中国资源卫星中心提供了大量的GF-4数据；感谢匿名审稿人细致的评审。