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

    基于MODIS的GF-4/PMS遥感器交叉定标以巴丹吉林沙漠为参考目标    张浩，刘涛，2，闫东川，阎跃观，崔珍珍，51.中国科学院空天信息创新研究院,北京 100094;2.中科星图空间技术有限公司,西安 710100;3.中国冶金地质总局矿产资源研究院,北京 101300;4.中国矿业大学（北京）,北京 100083;5.河南理工大学 测绘与国土信息工程学院,焦作 4540031 引言高分四号（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 calibration2.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 site2.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：// 中，计算出两个传感器模拟的表观辐亮度，得到两个传感器表观辐亮度的光谱匹配因子为评估地表二向性的影响，也同时利用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 对。