CN115452167A - 基于不变像元的卫星遥感器交叉定标方法和装置 - Google Patents

Abstract

本发明提供一种基于不变像元的卫星遥感器交叉定标方法和装置。方法包括：确定图像对序列；获取每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得每个图像对中的不变像元；获取待定标遥感器和参考遥感器的光谱匹配因子进行光谱匹配，确定每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；将每个图像对中的不变像元的参考遥感器光谱矫正表观反射率与对应的待定标遥感器表观反射率进行正交回归，确定交叉定标系数。本发明通过自动检测场景中的不变像元目标，计算光谱匹配因子以修正两传感器间的光谱响应差异，建立不变像元样本表观反射率的线性拟合关系，得到交叉定标结果及其长时间序列交叉定标系数。

Description

基于不变像元的卫星遥感器交叉定标方法和装置

技术领域

本发明涉及遥感器定标技术领域，尤其涉及一种基于不变像元的卫星遥感器交叉定标方法和装置。

背景技术

随着信息化和全球化的快速发展以及航天遥感探测技术的不断进步，遥感应用逐步深入到人类活动的各个领域，遥感器的精确定标是定量遥感应用的重要前提。

遥感器的辐射定标方法包括很多种：在绝对辐射定标中，需要预先获取观测目标的绝对辐射值，实施过程难度高且难以应用于历史卫星数据再定标；在轨星上定标方法依赖于星上定标设备，许多卫星缺少星上定标设备或受限于星上定标设备的研制水平，定标精度低或者无法进行在轨星上定标。

交叉定标是一种有效易行的在轨替代定标方法，其以具有高绝对定标精度的传感器为基准，通过对同一目标同时观测，对待定标的传感器进行交叉定标。然而，现有的交叉定标方法中存在不变目标选取困难，观测目标地表特性单一，反射率动态范围小以及定标频次低的缺陷。

发明内容

本发明提供一种基于不变像元的卫星遥感器交叉定标方法和装置，用以解决现有技术中对遥感器进行交叉定标时，不变目标选取困难、目标反射率动态范围小、定标频次低的缺陷。本发明无须通过人工选取不变目标，能有效提高传感器间的交叉定标频次，同时，大量的不变像元有效增加了反射率的覆盖范围，可适用于大场景数据。

本发明提供一种基于不变像元的卫星遥感器交叉定标方法，包括：

确定图像对序列，其中，每个所述图像对包括待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得所述每个图像对中的不变像元；

获取待定标遥感器和参考遥感器的光谱匹配因子，基于所述光谱匹配因子对所述每个图像对中的不变像元的参考遥感器表观反射率进行光谱匹配，确定所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；

将所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率与对应的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述确定图像对序列，包括：

采集卫星同时过境时待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

对所述待定标遥感器采集到的多光谱图像和所述参考遥感器采集到的多光谱图像进行预处理，得到目标图像对。

基于多天采集的所述目标图像对，确定所述图像对序列。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述预处理包括分辨率重采样、栅格化和剔除无效像元。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述剔除无效像元包括：

剔除云污染像元、水体目标像元和卫星观测天顶角大于等于30°的像元。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述不变像元检测模型是基于迭代加权多元变化检测IR-MAD方法建立的。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，输出所述每个图像对的不变像元，包括：

获取所述每个图像对中单个像元的探测结果并计算像元表观反射率；

将所述每个图像对中单个像元的像元表观反射率输入不变像元检测模型，构建MAD变量并基于所述MAD变量构建每个图像对中单个像元的观测值；

基于所述每个图像对中单个像元的观测值和不变概率决策阈值，确定所述每个图像对的不变像元。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述获取待定标遥感器和参考遥感器的光谱匹配因子，包括：

获取所述待定标遥感器的入瞳辐亮度和所述参考遥感器的入瞳辐亮度；

基于所述待定标遥感器的入瞳辐亮度和所述参考遥感器的入瞳辐亮度，确定所述待定标遥感器和参考遥感器的光谱匹配因子。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述将所述每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数，包括：

对所述每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率建立线性拟合关系进行正交回归；

基于正交回归的斜率和截距，确定获得交叉定标系数。

根据本发明提供的一种基于不变像元的卫星遥感器交叉定标方法，所述确定交叉定标系数，之后还包括：

对所述待定标遥感器长时间序列数据各通道进行长时间序列的交叉定标；

基于所述长时间序列的交叉定标的结果，确定长时间序列的交叉定标系数。

本发明还提供一种基于不变像元的卫星遥感器交叉定标装置，包括：

采集模块，用于确定图像对序列，每个所述图像对包括待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

不变像元检测模块，用于获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得所述每个图像对中的不变像元；

光谱匹配模块，用于获取待定标遥感器和参考遥感器的光谱匹配因子，基于所述光谱匹配因子对所述每个图像对中的不变像元的参考遥感器表观反射率进行光谱匹配，确定所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；

回归模块，用于将所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率与对应的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数。

本发明提供的基于不变像元的卫星遥感器交叉定标方法及装置，有别于传统的交叉定标方法，不依赖卫星过境人工选取的不变定标场，也不受限于地面同步观测的苛刻条件，能有效提高传感器间的交叉定标频次。对于不变像元目标的检测是自动进行的，对于处理长时间序列数据和历史卫星数据再定标有积极的作用。获取的不变像元目标在空间上是不连续的，虽然不能完整分析其地表的空间特性，但大量的不变像元样本有效增加了反射率的覆盖范围，且可适用于大场景数据，改善了传统交叉定标方法观测目标场地反射率范围小、样本单一的问题。定标结果可以用于遥感数据后续定量化反演产品的研究应用以及遥感器辐射响应特性的长期监测。

附图说明

为了更清楚地说明本发明或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明实施例提供的基于不变像元的卫星遥感器交叉定标方法的流程示意图之一；

图2是本发明实施例提供的基于不变像元的卫星遥感器交叉定标方法的流程示意图之二；

图3是本发明实施例提供的FY-3B/VIRR、MERSI通道光谱响应函数及SCIAMACHY高光谱样本示意图；

图4是本发明实施例提供的不变像元表观反射率正交回归结果示意图；

图5是本发明实施例提供的VIRR相对定标斜率长时间序列示意图；

图6是本发明实施例提供的基于不变像元的卫星遥感器交叉定标装置的结构示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

现有的交叉定标方法主要分为两类，一种是SNO交叉定标方法，其主要步骤为：

1)通过轨道预报获取两颗卫星轨道交叉点的地理位置和观测时间；

2)基于轨道预报交叉点的时间和地理位置，对观测数据进行像元匹配已选取用于交叉定标的数据，包括时间、观测几何、空间位置匹配；

3)根据两个遥感器相似通道的光谱响应差异进行光谱匹配；

4)回归分析计算定标系数。

另一种是场地交叉定标方法，其主要步骤为：

1)选取合适的不变场地作为观测目标；

2)获取两遥感器观测条件匹配的不变场地观测数据并投影到相同的地理栅格；

3)利用辐射传输模型和光谱匹配获取参考传感器的模拟表观反射率或表观辐亮度；

4)回归分析计算定标系数。

基于上述两种现有交叉定标方法，可以得出以下三个交叉定标存在的难题：

不变目标选取困难。SNO交叉定标需要利用轨道预报模型选取卫星的轨道交叉点，在有限的交叉点中经过严格的像元匹配才能获取可用的交叉定标数据，场地交叉定标同样需要人工选取划定不变场地，耗费人力物力。

观测目标地表特性单一，反射率动态范围小。如利用沙漠不变场进行场地交叉定标，沙漠场地的可见光波段反射率大约在0.2-0.3之间，但遥感器的辐射性能存在一定的目标依赖性，即定标系数会随着目标辐射的动态发生变化，因此需要大量的目标样本尽可能覆盖遥感器的动态范围，才能获得较高精度的定标结果。

定标频次低。无论是SNO交叉定标还是场地交叉定标，均由于人工选取不变目标的实施难度造成交叉定标的频次较低，无法展开连续的辐射定标以达到对遥感器辐射性能长期检测的目标。

基于上述问题，本发明实施例提供了一种基于不变像元的卫星遥感器交叉定标方法，实现对遥感器的交叉定标，解决对遥感器进行交叉定标时，不变目标选取困难、目标反射率动态范围小、定标频次低的缺陷。

下面结合图1-图5描述本发明的基于不变像元的卫星遥感器交叉定标方法，如图1所示，方法至少包括如下步骤：

步骤101、确定图像对序列，其中，每个图像对包括待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

步骤102、获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得所述每个图像对中的不变像元；

步骤103、获取待定标遥感器和参考遥感器的光谱匹配因子，基于所述光谱匹配因子对所述每个图像对中的不变像元的参考遥感器表观反射率进行光谱匹配，确定所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；

步骤104、将所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率与对应的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数。

针对步骤101，需要说明的是，图像对序列中包括多个图像对，而每组图像对都是选取卫星同时过境时两遥感器采集的同一目标场景的同时相多光谱图像。两个图像对可能是相邻两天的在同一时刻对于同一观测场景采集的两组多光谱图像，也可以是在同一天的不同时刻对于同一观测场景采集的两组多光谱图像，再或者也可以是日期、采集时刻、和观测场景均不同的两组多光谱图像。图像对序列需要采集连续多日的卫星数据，连续采集的时间要在卫星轨道回归周期内。

另外，每张多光谱图像均包含多个通道波段的卫星数据，在选取图像对时，还需要对两张同时相并且为同一个目标场景的多光谱图像进行通道匹配，将相似波段的通道一一对应。

针对步骤102，需要说明的是，每个图像对中单个像元的探测结果一般为遥感遥感器辐射测量值，即DN计数值。基于探测结果可以计算每个像元的表观反射率。不变像元检测模型是预先建立好的可以从输入的图像对中依据表观反射率筛选出不变像元的模型。对于每一个单日图像对，经过不变像元检测模型后，都可以输出该图像对的不变像元，将图像对序列依次输入不变像元检测模型后，可以得到不同图像对对应的不变像元的序列。

针对步骤103，需要说明的是，由于待定标遥感器与参考遥感器对应通道的光谱响应存在差异，即对于相同的入瞳辐射量，两遥感器会得到不同的测量值，因此需要对遥感器对应通道进行光谱匹配，以修正由于光谱响应差异造成的定标误差。

在本实施例中，待定标遥感器和参考遥感器的光谱匹配因子，需要使用高光谱仪器观测数据作为高光谱样本，分别与待定标遥感器与参考遥感器的通道光谱响应函数卷积，得到相应遥感器的入瞳辐亮度，再通过建立待定标遥感器与参考遥感器匹配通道入瞳幅亮度间的关系式即可获得。

针对步骤104，需要说明的是，任意一个图像对经不变像元检测模型获得的不变像元，都可通过步骤103和步骤104的处理后，获得交叉定标系数。然而，由于用于统计分析的场景广阔，单日图像对结果仅包含场景中符合观测天顶角限制条件的部分区域信息，为尽可能在单次回归中包含更丰富的不变像元样本，本实施例采集了图像对序列将连续多日检测得到的不变像元合并，进行拟合处理。基于遥感器辐射响应的短期稳定性，这种合并处理是可行的，并将有助于扩大像元样本的反射率值范围，有效提升回归效果与定标精度。

本发明实施例的基于不变像元的卫星遥感器交叉定标方法，针对性现有交叉定标方法存在的不变目标选取问题、观测目标的反射率动态范围限制问题以及长时间序列的遥感器连续定标监测问题，提出了一种行之有效的解决方案。该方法有别于传统的交叉定标方法，不依赖卫星过境人工选取的不变定标场，也不受限于地面同步观测的苛刻条件，能有效提高传感器间的交叉定标频次。对于不变像元目标的检测是自动进行的，对于处理长时间序列数据和历史卫星数据再定标有积极的作用。获取的不变像元目标在空间上是不连续的，虽然不能完整分析其地表的空间特性，但大量的不变像元样本有效增加了反射率的覆盖范围，且可适用于大场景数据，改善了传统交叉定标方法观测目标场地反射率范围小、样本单一的问题。定标结果可以用于遥感数据后续定量化反演产品的研究应用以及遥感器辐射响应特性的长期监测。

可以理解的是，确定图像对序列，包括：

采集卫星同时过境时待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

对待定标遥感器采集到的多光谱图像和参考遥感器采集到的多光谱图像进行预处理，得到目标图像对。

基于多天采集的目标图像对，确定图像对序列。

需要说明的是，由于不同遥感器之间以及不同通道之间的空间分辨率不同，数据维数不同则无法展开后续的统计学分析，因此需要对待定标遥感器采集到的多光谱图像和参考遥感器采集到的多光谱图像进行预处理，预处理后的图像对即为目标图像对。采集连续多日内的目标图像对，构成图像对序列。

可以理解的是，预处理包括分辨率重采样、栅格化和剔除无效像元。

需要说明的是，预处理时，需要对同一图像对中的两张多光谱图像的相同的通道波段进行一一匹配，随后步骤处理的基础即为同一通道波段下的图像对，并获得每个图像对的两张多光谱图像在各个通道波段下的DN值。在本发明实施例中，通过每个图像对的卫星数据进行分辨率重采样，即可实现通道波段的匹配。将多光谱图像对应的卫星遥感数据投影到同一地理栅格，将栅格化后的卫星遥感数据组建为一个卫星数据集。其中。地理栅格就是将地理空间分割成有规律的网格，每一个网格称为一个单元，并在各单元上赋予相应的属性值来表示卫星遥感数据的一种数据形式。

两遥感器获取的不同通道的DN值中，均存在失效数据。由于不变像元检测模型基于统计学原理展开，要求输入的两组数据维度一致，故需要将每个图像对中两张图像中的无效像元对应剔除。

可以理解的是，剔除无效像元包括：

剔除云污染像元、水体目标像元和卫星观测天顶角大于等于30°的像元。

需要说明的是，不变像元检测模型用于检测场景中的不变像元时，变化像元会在算法迭代过程中逐步剔除，但考虑到算法的运行效率，对于场景中显著的非稳定像元，又称云污染像元，例如云、沙尘等目标应在不变像元检测模型对应的算法运行前予以剔除。而水体及海洋等区域虽然有较为稳定的反射率特性，但其反射光谱特性差异较大，对部分通道波段信号微弱，对于水体目标像元利用数据中的陆地海洋掩膜数据进行剔除。

另外，图像对序列的数据经过无效点、云、沙尘和海洋水体像元去除后，仍有数百万数量级的像元样本，其中存在大量卫星观测天顶角过大的像元样本。卫星观测天顶角过大会导致像元空间分辨率和探测辐***度下降，为了提高不变像元检测模型对不变像元的识别精度，使用卫星观测天顶角小于30°的像元样本用于后续的统计分析。

可以理解的是，不变像元检测模型是基于迭代加权多元变化检测IR-MAD方法建立的。

需要说明的是，随着计算机信息技术的进步以及统计分析方法在科学研究中的广泛应用，遥感科研工作者逐渐将一些数学方法应用于卫星数据的处理分析中。现有技术中，差异图像比较法以及比值法适用于单通道图像的分析，不适用于多通道的卫星遥感器数据。主成分分析(PCA)法可以综合各个通道的变化信息，但不能消除不同遥感器通道之间的相关性。而一种具有线性尺度不变性的多变量变化检测技术(MAD)，可以消除同时相传感器通道内部的相关性以及传感器不同通道间的相关性。基于此，本发明实施例利用IR-MAD方法，构建了一种应用于不同遥感器同时相观测场景的不变像元检测模型，并将其应用于遥感器的交叉定标。

可以理解的是，本发明实施例中通过IR-MAD方法构建不变像元检测模型，获取每个图像对中的不变像元。其中，对同一时间待定标遥感器和参考遥感器采集的同一场景两幅n通道多光谱图像使用不变像元检测模型进行分析处理，至少包括如下步骤：

步骤201、获取每个图像对中单个像元的探测结果并计算像元表观反射率；

步骤202、将每个图像对中单个像元的像元表观反射率输入不变像元检测模型，构建MAD变量并基于MAD变量构建每个图像对中单个像元的观测值；

步骤203、基于每个图像对中单个像元的观测值和不变概率决策阈值，确定每个图像对的不变像元。

针对步骤201、需要说明的是，表观反射率(TOA)是指大气层顶的反射率，其值为地表反射率与大气透射率贡献之和。利用卫星运行初期的固定定标系数，由遥感器通道DN值计算像元表观反射率：

其中slope_i、bias_i分别是任一遥感器第i通道对应的固定定标斜率与截距，DN_i为任一遥感器获取的DN值，d²为日地距离修正因子，θ_s为太阳天顶角。

针对步骤202、需要说明的是，将每个图像对中单个像元的像元表观反射率输入不变像元检测模型后，首先需要构建MAD变量。构建MAD变量的过程包括：

步骤2021、分别用向量F＝(F₁…F_n)^T、G＝(G₁…G_n)^T表示一个图像对中两幅多光谱图像中DN值转换后的表观反射率，对所有谱带进行线性组合以获得典型变量U、V，如式1所示：

其中，n当前图像对的总匹配通道数，m和n是常向量，可通过求解式2的耦合广义特征值方程得出：

其中，ρ为典型变量U、V的相关系数，∑_ff、∑_gg、∑_fg、∑_gf为图像向量F、G的协方差矩阵。

步骤2022、基于典型变量U、V的差分得到MAD变量，如式3所示：

其中，MAD_i表示遥感器第i通道的MAD变量。

步骤2023、基于MAD变量的线性尺度不变性，得到观测值Z，如式4所示：

其中，

为遥感器第i通道的MAD变量的标准差，Z为标准化MAD变量的平方和，近似为卡方分布(χ²)，有n个自由度，Z值越小，像元不变概率越高。

针对步骤203、需要说明的是，对任一个图像对中两个遥感器同一时间采集的同一场景两幅n通道多光谱图像，理想情况下它们之间的差异是由仪器响应本身以及噪声、大气波动等随机效应引起的，从中心极限定理来看，MAD变量符合理想的正态分布。然而与变化观测相关的MAD变量或多或少会偏离这种多变量正态分布，在存在变化的情况下，把重点放在逐次迭代以建立一个越来越好的无变化背景上，那么MAD变换的敏感性将得到提高。

由于单次的MAD变换难以达到理想的不变像元检测效果，因此，当估计样本均值和协方差矩阵时，本发明实施例通过单次迭代确定的不变概率对样本数据进行加权，赋予更高不变概率的像元更大的权重，通过典型相关分析确定下一次迭代的MAD变量，多次迭代以获得更好的不变像元检测效果。

具体的，步骤203包括如下子步骤：

步骤2031、基于观测值，计算不变概率权重值；

需要说明的是，对于每次迭代，不变概率权重值Pr可由观测值Z的卡方分布检验结果式5确定：

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

步骤2032、设定不变概率决策阈值，基于每个图像对中单个像元的观测值和不变概率决策阈值，确定每个图像对的不变像元。

需要说明的是，由于本实施例寻求找到更大概率的不变像元，因此当匹配通道组数n确定时，对多个的不变概率置信度进行测试分析，综合考虑算法运行效率和不变像元检测效果，选择合适的不变概率决策阈值t，当观测值Z＜t时，认为该像元样本为不变像元的置信度高于

可用于后续交叉定标分析。

因此，本实施例中选择像元样本中满足式6条件的样本作为不变像元(NCPs)：

其中，t为具有n个自由度的卡方分布上分位数，即不变像元的概率决策阈值，

P为卡方检验值小于等于t的概率。

另外，在没有地表先验知识的情况下由数据统计特性选择的NCPs，其空间位置与图像对之间的不变特征相对应，但NCPs的位置会随不同图像对辐射信息的差异而改变。

可以理解的是，获取待定标遥感器和参考遥感器的光谱匹配因子，包括：

获取待定标遥感器的入瞳辐亮度和参考遥感器的入瞳辐亮度；

基于待定标遥感器的入瞳辐亮度和参考遥感器的入瞳辐亮度，确定待定标遥感器和参考遥感器的光谱匹配因子。

需要说明的是，本发明实施例使用高光谱仪器观测数据作为高光谱样本，分别与待定标遥感器与参考遥感器的通道光谱响应函数卷积，得到相应遥感器的入瞳辐亮度，如式7所示：

其中，R_h为高光谱样本辐亮度，f_{i_sensor}为高光谱仪器i通道的光谱响应函数，R_{i_sensor}为卷积后的待定标遥感器或参考遥感器入瞳辐亮度。

在相同的地表、大气和观测时间几何条件下，待定标遥感器和参考遥感器匹配通道的入同辐射量之比称为光谱匹配因子(SBAF)。建立待定标遥感器与参考遥感器匹配通道入瞳幅亮度间的关系式如式8所示，获得两遥感器对应通道的SBAF系数：

R_{i_CAL}＝A_i，j×R_{j_REF}+B_i，j 式8

其中R表示辐亮度，A、B为光谱匹配因子，i、j分别为待定标遥感器(CAL)与参考遥感器(REF)的匹配通道序号，A_i，j和B_i，j表示待定标遥感器i通道和参考遥感器j通道的光谱匹配因子，通过最小二乘回归拟合计算各匹配通道的SBAF。

可以理解的是，将每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数，包括：

对每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率建立线性拟合关系进行正交回归；

基于正交回归的斜率和截距，确定获得交叉定标系数。

需要说明的是，对IR-MAD检测获取的不变像元，将参考遥感器j通道的表观反射率探测结果ρ_{j_REF}利用光谱匹配因子进行转换，获得i通道光谱响应下的光谱矫正表观反射率光谱矫正表观反射率ρ_{i_REF}，如式9所示：

ρ_{i_REF}＝A_i，j×ρ_{j_REF}+B_i，j 式9

对任一图像对经IR-MAD检测获取的不变像元，将经光谱匹配得到的i通道光谱矫正表观反射率光谱矫正表观反射率ρ_{i_REF}与待定标遥感器的表观反射率探测结果ρ_{i_CAL}建立线性拟合关系进行正交回归如式10所示，从而获得交叉定标系数：

ρ_{i_REF}＝a×ρ_{i_CAL}+b 式10

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

可以理解的是，确定交叉定标系数，之后还包括：

对所述待定标遥感器长时间序列数据各通道进行长时间序列的交叉定标；

基于所述长时间序列的交叉定标的结果，确定长时间序列的交叉定标系数。

需要说明的是，对待定标遥感器各通道进行长时间序列的交叉定标，用于监测遥感器通道的长期辐射响应情况，确定长时间序列交叉定标系数趋势。

下面给出本发明方法的具体实施例：

风云三号B星(FY-3B)是中国第二代极轨气象卫星中的第二颗卫星，携带有11个遥感仪器，可见光红外扫描辐射计(VIRR)和中分辨率成像光谱仪(MERSI)是其中的两个主要仪器。FY-3B/VIRR共有10个通道波段，光谱范围0.43-12.5μm，其中有7个可见光近红外波段和3个红外发射波段，星下点空间分辨率为1.1km。VIRR主要用于探测识别云信息、探测地表植被覆盖、获取地表海面相关信息、监测大气水汽总量等。FY3B/MERSI具有19个反射太阳波段(0.41～2.13μm)和1个热红外波段(11.25μm)，星下点空间分辨率为250m和1000m。VIRR的1，2，6，7，8，9通道与MERSI的3，4，6，1，10_11，2通道设置相似(VIRR的8通道对应MERSI的10和11两个通道)，具体指标见表1，故将MERSI作为参考遥感器，将VIRR作为待定标遥感器，对VIRR共计6组对应通道进行交叉定标。

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

表1 FY3B/VIRR、MERSI光谱波段指标

可以理解的是，针对步骤101，本实施例选取的图像对序列的观测场景为中国西北地区，中心坐标为东经91°、北纬39°。地形以高原、盆地、山地为主，海拔1km以上。地貌景观主要为黄土高原、戈壁滩、荒漠草原和戈壁沙漠。此区域属于大陆性气候，受喜马拉雅山脉南部的影响，冬季寒冷干燥，夏季炎热，降水稀少，干旱是该地区的主要自然特征。研究区域地表反射率季节性变化小，符合理想的低气溶胶量、低水汽含量和高晴空概率的交叉定标观测区域标准。

如表1所示，VIRR与MERSI各通道分辨率存在差异，有250m，1km以及1.1km三种空间分辨率，对各通道探测数据重采样为1km空间分辨率，并投影到同一地理栅格，单幅场景数据大小为1400*3200像素，时间范围从2011年1月21日至2018年11月14日。

可以理解的是，针对步骤101，VIRR和MERSI采集得到的同时相多光谱图像DN值中，各通道分别存在一些无效点，须确保用于统计分析的每个像元各通道数据均为有效值，因此对存在无效值的像元进行剔除。

由于观测场景究区域显著变化场景主要是云，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)的辐射比值可以较好的区分地上云层与地面高反射率物体。对于海洋水体像元目标，利用VIRR与MERSI数据中的陆地海洋掩膜(landseamask)数据进行剔除。VIRR和MERSI遥感器视场扫描范围均为±55.4°，选择卫星观测天顶角小于30°的像元用于后续的统计分析。

可以理解的是，针对步骤102，利用FY-3B卫星运行初期遥感器的固定定标系数，由通道DN值计算像元表观反射率，获得研究场景同时相的表观反射率图像。

可以理解的是，针对步骤102，对同一时间VIRR和MERSI采集的中国西北地区两幅多光谱图像，经数据筛选和预处理后，保留用于IR-MAD分析的样本像元的6组匹配通道表观反射率数据，即n＝6。

在本实施例中，对多组图像对数据对进行IR-MAD分析，用于测试选取合适的迭代停止阈值参数，在迭代次数0-10次时，不变像元数量迅速下降，10-20次以后下降趋势逐渐平缓，20-25次时算法趋于收敛，此时两次迭代间的典型变量的相关系数变化小于0.001，因此可设定两次典型变量相关系数的变化差值小于0.001时即可视为算法收敛并停止迭代，并设定算法最大迭代次数为30次。同时，为了保证有足够的不变像元样本用于后续的交叉定标回归分析，设定最小的NCPs数量为400，当某次迭代后NCPs数量小于此数值，迭代停止。

在本实施例中，需要获取不变概率为90％以上的不变像元，当自由度n＝6(6组匹配通道)时，对90％、92.5％、95％、97.5％的不变概率置信度进行测试分析，综合考虑算法运行效率和不变像元检测效果，选择

作为不变概率决策阈值，当观测值Z＜t＝1.635时，认为该像元样本为不变像元的置信度高于95％，可用于后续交叉定标分析。

可以理解的是，针对步骤203，建立VIRR与MERSI匹配通道入瞳幅亮度间的关系式，获得两遥感器对应通道的SBAF系数，具体包括：

R_{i_VIRR}＝A_i，j×R_{j_MERSI}+B_i，j 式11

R_{8_VIRR}＝A_8，10×R_{10_MERSI}+A_8，11×R_{11_MERSI}+B_i，j 式12

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

表2 FY-3B VIRR与MERSI光谱匹配因子(SBAF)

可以理解的是，针对步骤203，对IR-MAD检测获取的不变像元，将MERSI j通道的表观反射率探测结果ρ_{j_MERSI}利用光谱匹配因子进行转换，获得i通道光谱响应下的光谱矫正表观反射率ρ_{i_MERSI}。

可以理解的是，针对步骤204，对任一图像对经IR-MAD检测获取的不变像元，将经光谱匹配得到的i通道光谱矫正表观反射率ρ_{i_MERSI}与VIRR的表观反射率探测结果ρ_{i_VIRR}建立线性拟合关系如式12所示，从而获得交叉定标系数。

ρ_{i_MERSI}＝a×ρ_{i_VIRR}+b 式12

其中，a为定标斜率，b为定标截距。

需要说明的是，风云三号卫星标称轨道回归周期为5.5天，对连续5天检测获取的不变像元合并进行正交回归分析。图4为2011年4月8日-12日不变像元各组匹配通道正交回归结果。各通道不变像元线性拟合效果良好，TOA动态范围大。

图5为2011年1月21日-2018年11月14日VIRR各通道相对定标斜率的长时间序列趋势，用于监测遥感器通道的长期辐射响应情况。本发明方法基于数学统计学原理展开，利用单日图像对进行分析时，不变像元的识别结果会受到云量大小、水汽气溶胶含量、极端恶劣天气情况(如暴雨暴雪、沙尘等)的影响，进而导致不变像元的TOA回归效果较差，因此可利用各通道回归结果的相关系数、残差等回归质量评价指标对长时间序列中不符合期望效果的数据对结果进行剔除。对某些数据对来说，仅有其中部分通道受到干扰时，考虑到方法消除了同时相传感器通道内部的相关性以及不同传感器不同通道间的相关性，且对通道的定标与分析具有独立性，因此可保留此类数据对有效的通道数据结果。长时间序列结果表明本发明的方法可进行近乎不间断的自动交叉定标实施。

下面对本发明提供的基于不变像元的卫星遥感器交叉定标装置进行描述，下文描述的基于不变像元的卫星遥感器交叉定标装置与上文描述的基于不变像元的卫星遥感器交叉定标方法可相互对应参照。

如图6所示，本发明实施例公开了一种基于不变像元的卫星遥感器交叉定标装置，包括：

采集模块601，用于确定图像对序列，每个图像对包括待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

不变像元检测模块602，用于获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得所述每个图像对中的不变像元；

光谱匹配模块603，用于获取待定标遥感器和参考遥感器的光谱匹配因子，基于所述光谱匹配因子对所述每个图像对中的不变像元的参考遥感器表观反射率进行光谱匹配，确定所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；

回归模块604，用于将所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率与对应的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数。

本发明实施例的基于不变像元的卫星遥感器交叉定标装置，针对光学遥感成像仪的定标问题和历史数据再定标需求，提出了基于智能选取不变像元的遥感器交叉定标方法。该方法利用两个遥感器获取的同时相卫星场景图像，经数据筛选与前处理后，通过不变像元检测模型自动检测场景中的不变像元目标，计算光谱匹配因子以修正两传感器间的光谱响应差异，建立不变像元样本表观反射率的线性拟合关系，得到交叉定标结果及其长时间序列定标趋势。

可以理解的是，采集模块601中确定图像对序列，包括：

采集卫星同时过境时待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

对待定标遥感器采集到的多光谱图像和参考遥感器采集到的多光谱图像进行预处理，得到目标图像对。

基于多天采集的目标图像对，确定图像对序列。

可以理解的是，预处理包括分辨率重采样、栅格化和剔除无效像元。

可以理解的是，剔除无效像元包括：

剔除云污染像元、水体目标像元和卫星观测天顶角大于等于30°的像元。

可以理解的是，不变像元检测模型是基于迭代加权多元变化检测IR-MAD方法建立的。

可以理解的是，不变像元检测模块602包括：

获取每个图像对中单个像元的探测结果并计算像元表观反射率；

将每个图像对中单个像元的像元表观反射率输入不变像元检测模型，构建MAD变量并基于MAD变量构建每个图像对中单个像元的观测值；

基于每个图像对中单个像元的观测值和不变概率决策阈值，确定每个图像对中的不变像元。

可以理解的是，光谱匹配模块603中获取待定标遥感器和参考遥感器的光谱匹配因子，包括：

获取待定标遥感器的入瞳辐亮度和参考遥感器的入瞳辐亮度；

基于待定标遥感器的入瞳辐亮度和参考遥感器的入瞳辐亮度，确定待定标遥感器和参考遥感器的光谱匹配因子。

可以理解的是，回归模块604包括：

对每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率建立线性拟合关系进行正交回归；

基于正交回归的斜率和截距，确定获得交叉定标系数。

可以理解的是，回归模块604还包括：

对待定标遥感器长时间序列数据各通道进行长时间序列的交叉定标；

基于长时间序列的交叉定标的结果，确定长时间序列的交叉定标系数。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (10)

1.一种基于不变像元的卫星遥感器交叉定标方法，其特征在于，包括：

确定图像对序列，其中，每个所述图像对包括待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得所述每个图像对中的不变像元；

获取待定标遥感器和参考遥感器的光谱匹配因子，基于所述光谱匹配因子对所述每个图像对中的不变像元的参考遥感器表观反射率进行光谱匹配，确定所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；

将所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率与对应的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数。

2.根据权利要求1所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述确定图像对序列，包括：

采集卫星同时过境时待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

对所述待定标遥感器采集到的多光谱图像和所述参考遥感器采集到的多光谱图像进行预处理，得到目标图像对；

基于多天采集的所述目标图像对，确定所述图像对序列。

3.根据权利要求2所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述预处理包括分辨率重采样、栅格化和剔除无效像元。

4.根据权利要求3所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述剔除无效像元包括：

剔除云污染像元、水体目标像元和卫星观测天顶角大于等于30°的像元。

5.根据权利要求1所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述不变像元检测模型是基于迭代加权多元变化检测IR-MAD方法建立的。

6.根据权利要求5所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，输出所述每个图像对中的不变像元，包括：

获取所述每个图像对中单个像元的探测结果并计算像元表观反射率；

将所述每个图像对中单个像元的像元表观反射率输入不变像元检测模型，构建MAD变量并基于所述MAD变量构建每个图像对中单个像元的观测值；

基于所述每个图像对中单个像元的观测值和不变概率决策阈值，确定所述每个图像对中的不变像元。

7.根据权利要求1所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述获取待定标遥感器和参考遥感器的光谱匹配因子，包括：

获取所述待定标遥感器的入瞳辐亮度和所述参考遥感器的入瞳辐亮度；

基于所述待定标遥感器的入瞳辐亮度和所述参考遥感器的入瞳辐亮度，确定所述待定标遥感器和参考遥感器的光谱匹配因子。

8.根据权利要求1所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述将所述每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数，包括：

对所述每个图像对的不变像元的参考遥感器的光谱矫正表观反射率与对应图像对的不变像元的待定标遥感器的表观反射率建立线性拟合关系进行正交回归；

基于正交回归的斜率和截距，确定获得交叉定标系数。

9.根据权利要求1所述的基于不变像元的卫星遥感器交叉定标方法，其特征在于，所述确定交叉定标系数，之后还包括：

对所述待定标遥感器长时间序列数据各通道进行长时间序列的交叉定标；

基于所述长时间序列的交叉定标的结果，确定长时间序列的交叉定标系数。

10.一种基于不变像元的卫星遥感器交叉定标装置，其特征在于，包括：

采集模块，用于确定图像对序列，每个所述图像对包括待定标遥感器和参考遥感器在同一天同一时刻下对同一观测场景分别观测得到的多光谱图像；

不变像元检测模块，用于获取所述每个图像对中单个像元的表观反射率并输入不变像元检测模型，获得所述每个图像对中的不变像元；

光谱匹配模块，用于获取待定标遥感器和参考遥感器的光谱匹配因子，基于所述光谱匹配因子对所述每个图像对中的不变像元的参考遥感器表观反射率进行光谱匹配，确定所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率；

回归模块，用于将所述每个图像对中的不变像元的参考遥感器的光谱矫正表观反射率与对应的待定标遥感器的表观反射率进行正交回归，确定交叉定标系数。

Images