基于IR-MAD不变像元的风云三号光学成像仪同平台交叉定标

李国荣，何玉青，胡秀清，王俊伟

1.北京理工大学 光电学院 光电成像技术与系统教育部重点实验室,北京 100081;2.国家卫星气象中心（国家空间天气监测预警中心）中国气象局,北京 100081

1 引言

可见光红外扫描辐射计（VIRR）是中国风云系列极轨卫星搭载的主要仪器之一，其积累了20年10 通道对地观测数据，这些数据对天气分析、气候预测和环境灾害业务应用发挥了重要作用。卫星传感器在轨运行期间受其运行环境和器件老化等因素的影响，辐射探测性能和稳定性会发生改变（Xu等，2014），严重影响其所获取辐射观测信号的准确性，因此需要对卫星传感器进行有效的辐射定标，为其数据产品的有效性和可用性提供保障。VIRR 没有配备相应的星载定标设备，主要依靠每年夏季在敦煌辐射校正场进行的场地替代定标方法进行定标校正（Hu等，2010）。然而一年一次的场地定标实验无法支持高频次和稳定的在轨校准系数更新，缺少对VIRR 辐射响应衰变的连续性监测和评估。此外，敦煌辐射校正场地表脆弱，晴空天气少，容易受到TSP不稳定和沙尘天气因素影响（王敏 等，2015），限制了VIRR 场地定标的精度。

交叉定标方法以具有高绝对定标精度的传感器为基准，通过对同一目标同时观测，对待定标的传感器进行相对定标，具有不受场地同步观测条件限制、有效提高定标频次的优势。Tong 等（2005）基于MODIS 敦煌辐射校正场的同步观测数据，对FY-1D/VIRR 的4 个VNIR 波段进行校准，部分通道校准精度达到5%；Zhong等（2018）利用MODIS反演了巴丹吉林沙漠场地BRDF和气溶胶光学厚度（AOD），得到2012 年内FY-3A/VIRR 仪器46 次交叉校准结果，误差约为5%。探测器辐射响应存在一定目标依赖性，基于单个场地展开的场地定标或交叉定标，均存在场地辐射特性单一、反射率动态范围小的问题。为了更好的表征仪器宽动态范围内的辐射响应情况，韦玮（2017）利用全球定标场网提高定标频次和目标辐射范围，得到的FY-3B/VIRR 全球定标场网定标结果与官方定标结果相对偏差为5.15%；Wang 等（2018）使用多站点（multisite，MST）交叉校准方法，将不同卫星平台搭载的VIRR 大气顶反射率差异从5%—10%降低至3%。考虑到定标频次和定标成本；何灵莉等（2020）使用自动化辐射定标方法对FY-3C/VIRR 进行了长时间序列辐射定标，3 年时间尺度上可定标观测100 天，在VNIR 具有3%的定标精度，在SWIR 具有5%的定标精度。探索一种自动化、高频次、多目标特点的交叉定标方法将有助于解决目前VIRR的现实定标问题。

随着计算机信息技术的进步以及统计分析方法在科学研究中的广泛应用，遥感科研工作者逐渐将一些数学方法应用于卫星数据的处理分析中。研究者们提出实现多光谱图像相对辐射标准化的方法是基于一个理想假设下进行的：在两个不同时间获取的来自恒定反射区域的传感器辐射之间的关系在空间上是均匀的，并且可用一个线性函数来近似（Moran 等，1992；Du 等，2002）。Nielsen等（1998）利用典型相关性分析（CCA）对两组卫星通道数据进行线性组合，提出了一种具有线性尺度不变性的多变量变化检测技术（MAD），可以消除同时相传感器通道内部的相关性以及双时相传感器不同通道间的相关性。Canty 等（2004）基于MAD 方法建立了一种自动化的变化检测算法用于获取卫星观测场景中的不变特征。Schmidt等（2008）应用MAD 方法分析了AVHRR 传感器长时间序列的响应变化，结果证明了此方法的有效性。随后Canty 和Nielsen（2008）提出了基于迭代重加权的多变量变化检测方法（IR-MAD），对于变化像元占主要部分的场景，例如季节性更替明显的耕地牧场等，IR-MAD方法在检测不变特征方面表现优异。王俊伟等（2019）利用IR-MAD方法识别同一场景不同时相卫星图像中的不变像元，并基于这些不变像元的辐射变化对FY-3A/MERSI 进行仪器响应衰变分析，其结果与传统方法获取的仪器衰变结果具有很好的一致性。IRMAD 方法对同一传感器不同时相场景中不变像元的检测具有优异性，对于具有相似通道的两个传感器，其相似通道在恒定反射区域获取的辐射量也具有明显的线性关系，因此此方法可以扩展到不同传感器同时相场景的不变像元检测中，以用于交叉定标分析。

本文针对FY-3B 上搭载的可见光红外扫描辐射计（VIRR）定标精度差、定标频次低的问题，以同平台的中分辨率光谱成像仪（MERSI）为基准，基于IR-MAD 方法对两传感器同时相观测场景中的不变像元进行自动化检测，利用不变像元的辐射信息展开交叉定标，并对交叉定标结果和长时间序列监测情况进行验证分析。

2 研究数据与区域

2.1 FY-3B/VIRR与MERSI

FY-3B 是中国第二代极轨气象卫星中的第二颗卫星，可见光红外扫描辐射计（VIRR）和中分辨率成像光谱仪（MERSI）是其中的两个主要仪器。FY-3B/VIRR 共有10 个通道波段，光谱范围0.43—12.5 µm，其中有7 个可见光近红外波段和3 个红外发射波段，星下点空间分辨率为1.1 km。FY-3B/MERSI 具有19 个反射太阳波段（0.41—2.13 µm）和1 个热红外波段（11.25 µm）。VIRR的1，2，6，7，8，9 通道与MERSI 的3，4，6，1，10_11，2 通道设置相似（VIRR 的8 通道对应MERSI 的10 和11 两个通道），具体指标见表1，故将MERSI 作为参考数据，对VIRR 共计6 组对应通道进行交叉定标。本研究使用的数据为FY-3B/VIRR 和MERSI 2011 年1 月21 日 至2018 年11 月14日近8年的1 km空间分辨率的L1数据。

表1 FY-3B/VIRR、MERSI 光谱波段指标Table 1 Spectral band specifications of FY-3B/VIRR、MERSI

一般的交叉定标研究需考虑不同卫星传感器由于观测时间、观测几何条件、辐射光谱响应、大气特性差异以及目标辐射特性带来的不确定性（Lacherade 等，2013）。本文研究仪器VIRR 及MERSI 位于同一卫星平台，且两探测器视场扫描范围均为±55.4°，在视场大小、观测时间、观测几何、大气条件等方面一致，这有利于IR-MAD 算法对场景内不变像元的检测，同时可近似认为两探测器的辐射标准化差异仅由仪器本身响应特性和光谱差异引起（Obata等，2017）。

2.2 研究目标区域

本文以中国西北地区为主要研究区域（图1（a）），中心坐标为（91°E、39°N）；主要包括内蒙古自治区，青海省，甘肃省和新疆维吾尔自治区；地形以高原、盆地、山地为主，海拔1 km 以上；地貌景观主要为黄土高原、戈壁滩、荒漠草原和戈壁沙漠。此区域属于大陆性气候，受喜马拉雅山脉南部的影响，冬季寒冷干燥，夏季炎热，降水稀少，干旱是该地区的主要自然特征。研究区域地表反射率季节性变化小，符合理想的低气溶胶量、低水汽含量和高晴空概率的交叉定标观测区域标准。同时为验证所提出方法的通用性，选用同时期北非地区（图1（b））数据进行比对验证。北非研究区域中心坐标为（15°E、23°N），以沙漠荒漠地貌为主，气候干旱，包含有地球观测卫星委员会和可见光遥感器定标和验证工作小组选定的多个国际校准目标场地，其地表反射率的稳定性和均匀性良好（Lacherade等，2013）。

3 基于IR-MAD 不变像元的交叉定标方法

基于IR-MAD 不变像元进行VIRR 与MERSI 交叉定标的处理过程如图2 所示。首先对同时相的VIRR 和MERSI 图像数据对进行筛选与预处理，将两探测器的匹配通道数据进行线性组合以构造MAD 变量，经IR-MAD 不变像元检测模型分析得到场景内的不变像元样本。基于获得的不变像元样本，以MERSI 探测结果为基准，经光谱匹配修正后，与VIRR 探测结果进行正交回归以得到交叉定标系数。

图2 基于IR-MAD不变像元的FY-3B/VIRR与MERSI交叉定标流程图Fig.2 Cross calibration flow chart of FY-3B/VIRR and MERSI based on IR-MAD No-change pixels

3.1 数据筛选和前处理

FY-3B/VIRR 与MERSI同平台观测，且可实现每天一次的全球扫描覆盖，因此以每日获取的研究区域同时相观测场景为一组匹配图像数据对。两传感器获取的通道DN 值数据中，均存在失效数据。IR-MAD方法基于统计学原理展开，输入模型的两组数据维度一致，故需要将两图像中的无效像元对应剔除。

IR-MAD 算法用于检测场景中的不变像元，变化像元会在算法迭代过程中逐步剔除，考虑到算法的运行效率，对于场景中显著的非稳定像元，如云、沙尘等目标应在算法运行前予以剔除。研究区域显著变化场景主要是云，VIRR 和MERSI 的L1 数据没有云掩膜数据，因此本文使用VIRR 的通道数据，采用阈值判别法进行初步的云识别。首先利用VIRR 的1 通道（0.58—0.68 µm）和2 通道（0.84—0.89 µm）两个波段的辐射比值区分云与晴空区域，再使用2 通道（0.84—0.89 µm）和6通道（1.55—1.64 µm）的辐射比值区分地上云层与地面高反射率物体（齐修东和岳继博，2016）。由于去云处理的主要目的是提高算法的效率，因此云识别结果存在部分的云漏判不会影响到不变像元的检测，这种简单的阈值法去云效果是可接受的。

水体及海洋等区域虽然有较为稳定的反射率特性，但其反射特性光谱差异较大，且对于研究选取的部分通道波段信号微弱。因此对于海洋水体像元目标，利用VIRR 与MERSI数据中的陆地海洋掩膜（landseamask）数据进行剔除。

研究使用的西北地区数据单幅场景大小为1400×3200 像素，北非地区单幅场景大小为1600×2600 像素，经过无效点、云和海洋水体像元去除后，仍有数百万数量级的像元样本，其中存在大量卫星观测天顶角过大的像元样本。卫星观测天顶角过大会导致像元空间分辨率和探测辐射精度下降，为了提高IR-MAD 对不变像元的识别精度，使用卫星观测天顶角小于30°的像元样本用于后续的统计分析。

FY-3B/VIRR 和MERSI同时观测同一地物目标时，具有相同的观测几何和大气路径，避免了水汽气溶胶、臭氧等大气因素的影响，使用卫星运行初期的固定定标系数，由VIRR 探测器通道DN值计算像元表观反射率：

式中，slopei、biasi分别是探测器第i通道对应的固定定标斜率与截距，DNi为探测器获取的DN 值，d2为日地距离修正因子，θs为太阳天顶角。

FY-3B/MERSI在轨定标分为两个阶段：2013年3 月6 日前采用线性定标模型，基于敦煌辐射校正场替代定标结果进行在轨校正，定标精度在波长小于1 µm 的窗区波段优于3%，波长大于1 µm 的波段优于5%（除2.1 µm 波段）（孙凌 等，2012）；2013年3月6日之后MERSI采用基于全球多目标辐射定标跟踪建立的日更新模型定标，其1—4、10、11波段与Aqua MODIS的TOA反射率的相对平均偏差在2%以内，6通道在5%以内（Sun和Li，2014）。为进一步分析MERSI数据的稳定性与同期VIRR 数据质量，对两仪器测得的敦煌辐射校正场TOA 反射率进行时间序列追踪，并计算各通道的年际相对平均偏差（RSD）（表2）。图3（a）所示为MERSI采用日更新定标模型后第一个自然年的时间序列结果，散点序列较为稳定，其1，2，3，10，11通道的RSD 值小于4%，4，6 通道约为5%；图3（b）为同时期VIRR 观测结果，时间序列存在较为明显的震荡情况，研究选取的各通道RSD 值在6%—8%，数据稳定性较差。MERSI 在定标精度及数据稳定性上较VIRR 有显著优势，虽然整体精度较MODIS仍有一定差距，但考虑到MERSI与VIRR同平台观测的一致性优势以及应用到IR-MAD 匹配场景的便利性，将其作为同期VIRR 的交叉定标参考基准是可行的。

图3 2013年3月—2014年3月FY-3B/VIRR与MERSI敦煌场地TOA反射率时间序列Fig.3 TOA reflectance time series of FY-3B/VIRR and MERSI in Dunhuang site from March 2013 to March 2014

表2 2013年3月—2014年3月VIRR与MERSI敦煌场TOA反射率序列相对标准偏差Table 2 RSD of TOA series between VIRR and MERSI in Dunhuang Site from March 2013 to March 2014

3.2 IR-MAD不变像元检测模型

将IR-MAD不变像元检测应用于同一时间VIRR和MERSI采集的同一场景的两幅n通道多光谱图像。

分别用向量F=(F1，…，Fn)T、G=(G1，…，Gn)T表示VIRR、MERSI 多光谱图像中DN 值转换后的表观反射率，对所有谱带进行线性组合以获得典型变量U、V：

式中，m和n是常向量，可通过求解耦合广义特征值方程得出：

式中，ρ为典型变量U、V的相关系数，∑ff、∑gg、∑fg、∑gf为图像向量F、G的协方差矩阵。

MAD变量由典型变量U、V的差分得到，即：

由于MAD 变量具有线性尺度不变性，因此我们选择像元样本中满足

单次的MAD 变换难以达到理想的不变像元检测效果。通过单次迭代确定的不变概率对样本数据进行加权，赋予更高不变概率的像元更大的权重，使用典型相关分析确定下一次迭代的MAD 变量，多次迭代以获得更好的不变像元检测效果。对于每次迭代，不变概率权重值Pr可由观测值Z的卡方分布检验结果确定：

为寻求找到更大概率（90%以上）的不变像元，当自由度n=6（6 组匹配通道）时，对90%、92.5%、95%、97.5%的不变概率置信度进行测试分析，综合考虑算法运行效率和不变像元检测效果，选择t==1.635 作为不变概率决策阈值，当观测值Z＜t=1.635时，即可认为该像元样本为不变像元的置信度高于95%，可用于后续交叉定标分析。

对于IR-MAD算法，需要设定典型变量相关系数变化阈值、最大迭代次数和最小NCPs 数量3 个迭代停止阈值（王俊伟 等，2019），使其在算法优化收敛以及小概率检测失效时停止对某一图像对的迭代处理，以保证算法的自动运行。本研究设定两次典型变量相关系数的变化差值小于0.001 时即可视为算法收敛并停止迭代，并设定算法最大迭代次数为30 次。同时，为保证有足够的不变像元样本用于后续的交叉定标回归分析，设定最小的NCPs 数量为400，当某次迭代后NCPs 数量小于此数值，迭代停止。

IR-MAD 算法从像元样本的辐射信息出发，构造通道间线性组合，消除单个传感器不同通道以及不同传感器非匹配通道间的相关性，经过多次不变概率加权迭代，获得两幅观测场景中的不变像元，不变像元辐射量具有对VIRR 与MERSI匹配通道间的最大线性相关性。对不变像元的空间分布特征进行分析（图4），图4（a）为2013 年6 月3 日获取的不变像元检测结果，宽幅场景图像包含两个FY-3B轨道，算法检测范围包含被分割的左右两块观测天顶角＜30°的区域，检测得到的不变像元（红色标记点）主要分布于新疆塔克拉玛干沙漠、内蒙巴丹吉林沙漠及腾格里沙漠区域。图4（b）为研究区域内的像元不变概率分布图，通过对获取的长时间序列不变像元结果进行概率频次的叠加得到，红色表示像元的不变概率更高。时间序列结果显示，高概率不变像元主要聚集在塔克拉玛干沙漠，敦煌沙漠，巴丹吉林沙漠，腾格里沙漠和一些小戈壁地区，这与使用VIRR 以及MODIS等单个传感器不同时相场景获取的分布情况相一致（Hu等，2020）。相较于场景内其他目标，不变像元具有更好的时间稳定性，基于稳定定标场的交叉定标也常选取沙漠目标进行，这体现了IRMAD方法在目标选用上的合理性。

图4 中国西北地区不变像元分布结果Fig.4 NCPs distribution in Northwest China region

3.3 光谱匹配

VIRR与MERSI对应通道的光谱响应（图5）存在差异，即对于相同的入瞳辐射量，两传感器会得到不同的测量值。因此需要对传感器对应通道进行光谱匹配，以修正由于光谱响应差异造成的定标误差。同平台观测条件下（相同的地表、大气和观测时间几何条件），两传感器匹配通道的入同辐射量之比称为光谱匹配因子（SBAF）。使用SCIAMACHY仪器星下点观测模式的数据作为高光谱样本（图5），分别与VIRR 和MERSI的光谱响应函数卷积，得到相应传感器的入瞳辐亮度，关系式为

图5 FY-3B/VIRR、MERSI 通道光谱响应函数及SCIAMACHY高光谱样本Fig.5 Channel SRF of FY-3B/VIRR and MERSI，IAMACHY hyperspectral samples

式中，Rscia为SCIAMACHY 光谱样本辐亮度，fi_sensor为仪器i通道的光谱响应函数，Ri_sensor为卷积后的仪器入瞳辐亮度。

建立VIRR 与MERSI匹配通道入瞳幅亮度间的关系式：

式中，R表示辐亮度，A、B为光谱匹配因子，i、j分别为VIRR与MERSI的匹配通道序号，对于VIRR与MERSI 的单对单匹配通道使用式（8），VIRR 的8 通道和MERSI 的10、11 通道使用式（9），通过最小二乘回归拟合计算SBAF，结果如表3。

表3 FY-3B VIRR与MERSI光谱匹配因子（SBAF）Table 3 Spectral band adjustment factor of FY-3B VIRR and MERSI

3.4 回归分析与定标计算

对IR-MAD 检测获取的不变像元，以MERSI探测结果作为辐射基准，利用光谱匹配因子计算VIRR 表观反射率准真值，将其与VIRR 探测结果建立线性拟合关系，从而获得交叉定标系数。

式中，a为定标斜率，b为定标截距。ρVIRR为基于发射初期固定定标斜率计算得到的表观反射率，因此定标系数a，b是相较于发射初期固定定标系数的相对定标结果。

本研究区域广阔，单日图像对结果仅包含场景中符合观测天顶角限制条件的部分区域信息，为尽可能在单次回归中包含更丰富的不变像元样本，将连续多日检测得到的不变像元合并处理。基于传感器辐射响应的短期稳定性，这种合并处理是可行的，并将有助于扩大像元样本的反射率值范围，有效提升回归效果与定标精度。考虑到风云三号卫星标称轨道回归周期为5.5 天，同时对于识别出的不变像元样本，Canty 等（2004）分析得出正交回归是更好的回归处理方法，对连续5天检测获取的不变像元合并进行正交回归分析。图6为2011 年4 月8 日—12 日不变像元各组匹配通道正交回归结果。

图6 2011年4月8日—12日不变像元表观反射率正交回归结果Fig.6 Orthogonal regression results of TOA reflectance of NCPs from April 8 to 12，2011

4 交叉定标结果与分析

4.1 交叉定标质量精度分析

基于上述的IR-MAD 不变像元检测和交叉定标回归分析，对2011 年1 月21 日至2018 年11 月14 日数据进行分析处理，共得到中国西北地区、北非地区的5 日合并正交回归结果数量分别为523、551 组。为保证回归结果的质量，对每个回归结果的NCPs 样本表观反射率相关系数进行考察，图7 给出了所有回归结果8-10_11 匹配通道的NCPs 样本反射率相关系数分布情况，相关系数值大于0.99，对其余匹配通道相关系数小于0.99 的结果进行剔除。

图7 8—10_11匹配通道的图像对NCPs样本相关系数Fig.7 NCPs sample correlation coefficients of image pairs with 8-10_11 matching channels

为验证基于IR-MAD 不变像元交叉定标方法所得到的VIRR 交叉定标结果精度，对2012 年1 月15 日—2 月15 日中国西北、北非地区IR-MAD 交叉定标结果与2012 年2 月1 日VIRR 业务定标值进行对比，如图8所示。需指出的是，本文后续表中所列斜率结果均为相对定标系数与发射初期固定定标系数的乘积结果。如表4所列具体数值，中国西北地区IR-MAD 自然月内交叉定标斜率的相对标准偏差值小于1.5%，表明了定标结果的稳定性。除VIRR 7 通道外，其余各通道的IR-MAD 中国西北与北非的定标结果一致性优于2%，说明所述方法基本不受地理区域选择限制，有较好的通用性。VIRR 7 通道为短波蓝色光波段，受大气和气溶胶分子散射影响明显，虽然用于IR-MAD 分析的图像对大气状况几乎一致，但同期中国西北和北非两区域大气状况存在显著差异，这造成了两区域所得VIRR 7 通道IR-MAD 定标结果相较于其余通道的较大偏差（3.73%）。

图8 中国西北、北非地区IR-MAD定标结果与2012年2月1日业务定标结果对比Fig.8 Comparison of calibration results between IR-MAD in Northwest China and North Africa with business calibration results on February 1，2012

表4 中国西北、北非地区IR-MAD交叉定标结果与VIRR业务定标对比Table 4 Comparison of cross calibration results between IR-MAD in northwest China and north Africa with VIRR business calibration

同时，IR-MAD 西北地区各通道的交叉定标结果与业务定标结果的相对偏差均小于4%，1 通道（红色波段）和2 通道（近红外波段）的结果最为接近，符合这两个通道的较为稳定的特征。

VIRR 业务定标频次很低，为进一步验证定标结果精度，将IR-MAD 结果与基于敦煌辐射校正场的FY-3B/VIRR 定标结果（吕佳彦 等，2017）进行比较，结果如表5。两方法定标结果显示，VIRR 1，2通道的相对偏差基本在2%以内，7，8，9 通道的相对偏差低于5%。VIRR 7，8，9 这3 个通道是蓝绿光波段，相较于红色、近红外和短波红外波段，VIRR 和MERSI 在这几个波段有较大的衰变情况，MERSI 的业务定标结果对这种衰变的辐射校正精度在一定程度上影响了IR-MAD 方法的交叉定标精度。此外，两方法结果显示的各通道定标斜率的年际变化趋势一致，其中6通道（短波红外）表现出与其他5个通道不一致的趋势，具体为先增大后回落，此情况发生的原因可能是短波红外通道电子增益出现随机跃动，导致定标系数存在波动，胡秀清等（2013）对风云三号A 星MERSI 的两个短波红外通道响应追踪中也发现了类似现象。

表5 IR-MAD不变像元交叉定标与敦煌场地定标结果对比Table 5 Comparison of cross calibration results between IR-MAD NCPs and Dunhuang site calibration

表6 给出了2012 年不变像元样本数量及VIRR定标后的不变像元反射率较MERSI 观测的相对平均偏差。可以明显看出，冬季月份的不变像元样本数量高于夏季月份，这是由于冬季晴空概率高且水汽含量低，场景中可用于统计分析的像元数量更多。定标校正后，两仪器表观反射率的相对平均偏差小于2%。

表6 2012年不变像元样本VIRR定标反射率较MERSI观测的相对平均偏差Table 6 Relative mean deviation of VIRR calibration reflectivity of NCPs compared with MERSI observation in 2012

4.2 长时间序列监测分析

图9考察了2011年1月21日—2018年11月14日VIRR 各通道相对定标斜率的长时间序列趋势，分别基于中国西北地区和北非地区得到的IR-MAD交叉定标结果显现出高度的一致性。从图9中可以明显发现，北非地区有更好的散点拟合线聚合情况，这与北非地区优异的地表反射率稳定性和空间均匀性相关。

图9 2011年1月21日—2018年11月14日VIRR相对定标斜率长时间序列Fig.9 Long time series of VIRR relative calibration slope from January 21，2011 to November 14，2018

除VIRR 6 通道外，其余各通道自研究数据的初始时间（2011 年1 月21 日）开始，距相对值1存在较大的差异。这是因为回归分析所用的VIRR像元表观反射率通过卫星发射初期使用的固定定标系数（发射前定标结果）计算获取，而在轨运行一段时间后，仪器经历发射阶段和在轨运行期间的复杂环境，在轨辐射性能与发射前定标结果存在较大差异。VIRR 的通道响应衰减趋势与波长表现出一定的相关性，整体上呈现波长越短，衰减越大的趋势。VIRR 的7 和8 通道在轨运行前两年期间，响应衰减幅度较小，第3—5 年出现大幅度的衰减，之后趋于相对平缓，这与吕佳彦等（2017）基于敦煌场地定标的FY-3B/VIRR 通道衰减规律一致。

本文提出的方法可进行近乎不间断的自动交叉定标实施，在VIRR 2 和6 通道的中国西北地区定标结果中发现了明显的季节性波动现象；2通道的北非地区结果也存在类似现象，但幅度较小；两地区的1 通道结果也存在小幅度的季节性波动。首先，本文IR-MAD 算法使用的同时相单日图像对大气条件一致，因此未对TOA 反射率进行大气校正，大气的季节性变化会带来相应的序列跟踪不确定度；IR-MAD 算法从像素级展开统计分析，不需要预先获取地表的先验知识，不变像元的BRDF 季节性变化以及对地表均匀性产生影响的季节性降雪、沙尘等天气活动也是定标序列波动的重要影响因素。此外，孙凌等（2013）对FY-3A/MERSI 的多场地定标跟踪以及何灵莉等（2020）对FY-3C/VIRR 的自动化定标序列中均发现了季节性波动的情况。除地表和大气因素外，仪器本身的温度响应、通道光谱响应特性也会造成不同季节的响应差异，如VIRR 6 波段对云雪较为敏感，中国西北地区相较于北非地区，存在季节性的降雪活动，由于在预处理中未进行积雪识别和去除，这一定程度上导致了两地区在6通道序列结果的不同表现。北非地区主要为低纬度区域，而中国西北地区为中纬度区域，其太阳活动和高度角的季节性变化更为明显，因此近红外波段的响应衰减季节性变化也可能与近红外波段太阳辐射和地表反射率的分光谱特征季节相关性有关。

5 结论

本文针对风云三号卫星光学成像仪的定标问题和历史数据再定标需求，以FY-3B 卫星的VIRR（待定标）和MERSI（参考基准）为研究对象，提出了基于IR-MAD 不变像元的交叉定标方法。该方法利用同平台传感器获取的同时相卫星场景图像，经数据筛选与预处理后，通过IR-MAD 不变像元检测模型自动检测场景中的不变像元目标，计算光谱匹配因子以修正两传感器间的光谱响应差异，建立不变像元样本表观反射率的线性拟合关系，得到接近8年的时间序列定标结果。该方法有别于传统的交叉定标方法，不依赖卫星过境人工选取的伪不变定标场，也不受限于地面同步观测的苛刻条件，能有效提高传感器间的交叉定标频次。研究结果表明，IR-MAD 方法获取的交叉定标结果与业务定标结果的差异小于4%；同时，定标精度和通道波段有一定的相关性，与敦煌沙漠场定标结果相比，VIRR 红色和近红外通道的定标结果优于3个蓝绿短波通道。此外，中国西北和北非地区获得的定标结果长时间序列表现出良好的一致性。

IR-MAD 方法对于不变像元目标的检测是自动进行的，对于处理长时间序列数据和历史卫星数据再定标有积极的作用。同时算法获取的不变像元目标在空间上是不连续的，虽然不能完整分析其地表的空间特性，但大量的不变像元样本有效增加了反射率的覆盖范围，可适用于大场景数据，改善了传统交叉定标方法观测目标场地反射率范围小、样本单一的问题。

目前基于全球定标场网或自动化辐射计的连续定标方法，对于数量庞大的卫星仪器，每个仪器均实施连续定标的成本代价是巨大的。本文研究的IR-MAD 方法，在具有一个连续定标精准仪器的情况下，可实现对同平台其他仪器的高效相对响应监测，解决宽动态、全视场以及数据处理效率方面的问题，有效降低定标成本。IR-MAD方法已实现单个传感器响应衰变分析以及同平台不同传感器间相对响应的交叉定标，相比于以上两种情况，异平台传感器间存在观测几何、时间、大气以及光谱差异等误差来源。在宽视角条件下，异平台传感器的观测几何条件匹配并不困难，光谱差异可通过光谱匹配因子加以校正，在辐射来源稳定的情况下，时间因素的主要误差来源于大气变化，因此大气差异的校正是实现异平台仪器间IR-MAD 交叉定标的关键，借助精细的大气参量和辐射传输模型，我们对异平台传感器间的IRMAD算法表现有较好的期望。

本文方法在仪器的非线性响应和场景的空间特征分析方面仍有不足。一些研究者们提出基于核典型相关分析（KCCA）的IR-MAD 方法和对MAD 变量进行最大自相关因子（MAF）分析将有助于解决非线性和空间特性问题，这将是我们优化核心算法的工作方向。