一种遥感图像相对辐射校正方法及系统与流程

本发明涉及声像技术领域，特别是涉及一种遥感图像相对辐射校正方法及系统。

背景技术：

随着遥感技术飞速发展，利用多传感器、多时相遥感数据进行土地利用和土地覆盖变化监测、全球资源环境分析、气候变化监测等应用越来越广泛，拓宽了遥感数据在时间维度和空间维度的利用，这种大区域大维度的数据应用，往往对数据集一致性要求非常高，甚至高于对单幅图像的精度要求。

对于多景数据处理，由于每景数据的处理误差、获取时间的差异以及传感器差异等问题，造成拼接后的在数据拼接线附近出现“边界”现象，这种“边界”现象增大了后续工程应用的难度。因此，现有多数采用相对辐射校正技术进行辐射处理，能最大限度保证图像间的辐射一致性，不仅能够纠正大气状况变化带来的差异，而且能够削弱传感器等其它原因产生的噪音。

虽然现有的相对辐射校正方法能较好的保证图像之间辐射一致性，但对图像时相有严格的约束条件：要求参考图像与目标图像的成像时间相近，即要求可以是不同年份的，但月份要求在一个月以内，成像时刻在1小时以内，这一约束大大限制了全国、甚至全球范围内遥感数据进行土地利用、资源分析等应用，特别是针对时相差别较大、传感器成像几何差异等原因造成的图像本身辐射差异较大，很难保证图像之间辐射一致性，因此大大降低了遥信数据的利用率。

技术实现要素：

本发明的目的是提供一种遥感图像相对辐射校正方法及系统，以保证图像间的辐射一致性，以提高遥信数据的利用率。

为实现上述目的，本发明提供了一种遥感图像相对辐射校正方法，该包括以下步骤：

获取目标图像；

根据所述目标图像确定与所述目标图像相对应的参考图像；

分别对所述参考图像、所述目标图像进行预处理后，获得对应所述参考图像的参考图像像素点集和对应所述目标图像的目标图像像素点集；

基于核典型相关分析方法根据所述参考图像像素点集和所述目标图像像素点集确定非线性辐射控制点；

根据所述非线性辐射控制点对所述目标图像各阶段像素点进行校正，获得地表反射率值。

可选的，所述根据所述目标图像确定与所述目标图像相对应的参考图像的具体步骤包括：根据与所述目标图像具有相同重叠区域、近似成像几何、坐标范围的图像确定为与所述目标图像相对应的参考图像。

可选的，所述基于核典型相关分析方法根据所述参考图像像素点集和所述目标图像像素点集确定非线性辐射控制点的具体步骤包括：

分别将所述参考图像像素点集、所述目标图像像素点集映射到核空间，利用核典型相关分析方法确定参考图像核典型变量、目标图像核典型变量；

基于mad中的权重公式对所述参考图像核典型变量、所述目标图像核典型变量进行筛选，确定非线性辐射控制点。

可选的，所述权重公式为：

其中：u为映射到高维特征空间的目标图像核典型变量集合，ui为u的第i个元素，v为映射到高维特征空间的参考图像核典型变量集合，vi是v的第i个元素，var()为方差，z为标准化后的方差，ρλi为第i个特征根，为卡方检验，pr(nochange)为不变像素点的概率值，n为重叠区域的点的个数。

可选的，所述根据所述非线性辐射控制点对所述目标图像各阶段像素点进行校正，获得地表反射率值的具体步骤包括：

对所述非线性辐射控制点进行非线性关系拟合获得非线性拟合公式；

根据所述非线性拟合公式对所述目标图像各波段像素点进行校正，获得地表反射率。

可选的，所述对所述非线性辐射控制点进行非线性关系拟合获得非线性拟合公式具体步骤包括：

根据偏差平方和最小公式对所述非线性辐射控制点进行非线性关系拟合获得非线性拟合公式；

所述偏差平方和最小公式为：

其中，δi为第i个非线性辐射控制点的偏差，ri为第i个非线性辐射控制点的原始目标图像像素值，f(ri)为将第i个非线性辐射控制点的原始目标图像像素值代入拟合的非线性关系表达式后，计算得到的像素值，si为第i个非线性辐射控制点对应的参考图像像素值，r为非线性辐射控制点的总个数。

本发明还提供了一种遥感图像辐射校正系统，该系统包括：

获取模块，用于获取目标图像；

匹配模块，用于根据所述目标图像确定与所述目标图像相对应的参考图像；

处理模块，用于分别对所述参考图像、所述目标图像进行预处理后，获得对应所述参考图像的参考图像像素点集和对应所述目标图像的目标图像像素点集；

非线性辐射控制点确定模块，用于基于核典型相关分析方法根据所述参考图像像素点集和所述目标图像像素点集确定非线性辐射控制点；

校正模块，用于根据所述非线性辐射控制点对所述目标图像各阶段像素点进行校正，获得地表反射率值。

可选的，所述非线性辐射控制点确定模块的具体包括：

第一确定单元，用于分别将所述参考图像像素点集、所述目标图像像素点集映射到核空间，利用核典型相关分析方法确定参考图像核典型变量、目标图像核典型变量；

第二确定单元，用于基于mad中的权重公式对所述参考图像核典型变量、所述目标图像核典型变量进行筛选，确定非线性辐射控制点。

可选的，所述权重公式为：

其中：u为映射到高维特征空间的目标图像核典型变量集合，ui为u的第i个元素，v为映射到高维特征空间的参考图像核典型变量集合，vi是v的第i个元素，var()为方差，z为标准化后的方差，ρλi为第i个特征根，为卡方检验，pr(nochange)为不变像素点的概率值，n为重叠区域的点的个数。

可选的，所述校正模块具体包括：

拟合单元，用于对所述非线性辐射控制点进行非线性关系拟合获得非线性拟合公式；

校正单元，用于根据所述非线性拟合公式对所述目标图像各波段像素点进行校正，获得地表反射率。

根据本发明提供的具体实施例，本发明公开了以下技术效果：

1)、本发明基于核典型相关分析方法自动确定非线性辐射控制点，保证目标图像和参考图像间光谱一致，消除图像之间的辐射差异，提高图像间的辐射一致性；特别是针对时相差别较大、传感器成像几何差异等原因造成的图像本身辐射差异较大，也能保证图像之间辐射一致性，大大提高了遥信数据的利用率。

2)、本发明采用核典型相关分析的办法自动确定非线性辐射控制点，避免了选择参考图像时对时相的约束和受传感器的限制、通用性强、无需人工交互、处理流程单、运算速度快、算法自动性高、适用性强等特点。

附图说明

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

图1为本发明实施例遥感图像相对辐射校正方法的具体流程图；

图2为本发明实施例遥感图像相对辐射校正系统的结构示意图。

具体实施方式

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

本发明的目的是提供一种遥感图像相对辐射校正方法及系统，以保证图像间的辐射一致性，以提高遥信数据的利用率。

为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。

图1为本发明实施例遥感图像相对辐射校正方法的具体流程图；详见图1。其具体的步骤包括：

步骤100：根据所述目标图像确定与所述目标图像相对应的参考图像。

根据与所述目标图像具有相同重叠区域、近似成像几何、坐标范围的图像确定为与所述目标图像相对应的参考图像，所述参考图像与所述目标图像无需时相相近，大大提高辐射处理图像间的辐射一致性，提高了遥信数据的利用率。

步骤200：分别对所述参考图像、所述目标图像进行预处理后，分别获得参考图像像素点集和目标图像像素点集。

所述预处理包括坐标重投影、几何配准、计算地理重叠区域、各波段像素点对应、并使用掩膜技术排除极值点、过饱和点、云、云下阴影及水体的像素点。对所述参考图像、所述目标图像进一步处理，为了获得较为精确的参考图像像素点集和目标图像像素点集。

坐标重投影与几何配准属通常的遥感图像处理步骤。

计算地理重叠区域：对于单个目标图像，有可能找到多景参考图像，此时就需要进行重叠区域内数据的裁切与拼接整合。

比各波段像素点对应：对于目标图像与参考图像为同种传感器数据，分辨率相同，当坐标重投影、几何配准精度足够高时，像素点本身就是对应的。对于不同传感器的遥感图像，分辨率差异超过3倍以上时，就需要插值。

掩膜过滤异常点：该步先通过统计目标图像与参考图像重叠区内的直方图，过滤掉过亮和过暗的像素点，然后检测参考图像中的云、云下阴影、水体并将其掩膜掉。具体的检测算法依据不同的传感器数据而不同。

对掩膜过滤后得到的数据点集进行典型相关分析。

对掩膜过滤后得到的数据点集进行典型相关分析，进一步从这些点集中选出典型相关点集。所选用的典型相关分析算法应当具有不受总体的大气状况和传感器标定导致的线性变换的影响的特性。

典型相关分析算法为了遮蔽两个不同时相图像中的变化像素，首先形成两幅图像n个波段内像素值的线性组合。

步骤300：基于核典型相关分析方法根据所述参考图像像素点集和所述目标图像像素点集确定非线性辐射控制点。

步骤301：分别将所述参考图像像素点集、所述目标图像像素点集映射到核空间，利用核典型相关分析方法确定参考图像核典型变量、目标图像核典型变量。

核典型相关分析方法是一种常用的非线性典型相关方法，将核函数的思想引入典型相关分析当中，把低维数据映射到高维的特征空间中，为在高维特征中解决复杂关联分析奠定基础。

常用的核函数种类非常多，本次实施中选取的核函数为anova核函数该核函数属于径向基的一种，非常适合多维回归的问题。其具体规则为：设x,z∈x，x属于rⁿ空间，非线性函数实现输入空间x到特征空间f的映射，其中f属于r^m，n<<m。则核函数公式为：

其中<,>为内积，k(x,z)即为核函数。

x为目标图像像素点集，用目标图像像素点矩阵表示，此处作为样本矩阵输入，n为重叠区域的点的个数，n为每个像素点的数据维度，即为各波段的灰度值，则x可为n个n维向量，xn×n＝(x1,x2...,xn)，xi为第i个n维向量，y为参考图像像素点集，用参考图像像素点矩阵表示，yi为第i个n维向量，此处作为样本矩阵输入，n为重叠区域的点的个数，n为每个像素点的数据维度，即为各波段的灰度值，y可为n个n维向量，yn×n＝(y1,y2...,yn)。

依据核典型相关分析方法，将经过预处理的目标图像和参考图像的n个波段内像素值分别映射为高维空间的变量，获取其在高维特征空间的线性组合，具体公式为：

其中，c和d分别为高维空间(m维)的常向量矩阵，c^t为c的转置，d^t为d的转置，将目标图像的向量x由低维空间(rⁿ)映射到高维空间(r^m)的非线性变换函数(n<<m)，为将参考图像的向量y由低维空间(rⁿ)映射到高维空间(r^m)的非线性变换函数，变换后的和均为m×n维矩阵，u为映射到高维特征空间的目标图像核典型变量集合，v为映射到高维特征空间的参考图像核典型变量集合。

常向量矩阵c和d可以用n维空间的常向量矩阵α和β与和为：

其中，α＝(α1,α2,,...,αn)，αi是α的第i个元素，每一个αi是一个n维向量；同理βi是β的第i个元素。此时将求解高维空间常向量c和d的问题转化为求解低维空间常向量α和β的问题。将上述公式(5)代入公式(4)，则可得到综合变量集合u和v为：

步骤302：基于mad中的权重公式对所述参考图像核典型变量、所述目标图像核典型变量进行筛选，确定非线性辐射控制点，其中、mad为多元变化检测算法(multivariatealterationdetection)。

根据核函数公式(3)确定在x和y上定义核矩阵kx和ky为：

其中，kx(i,j)表示核矩阵kx的第i行j列元素，ky(i,j)表示核矩阵ky的第i行j列元素。

通过综合变量u和v变换公式(6)和核矩阵公式(7)，我们可以计算出u和v的方差和协方差为公式为：

其中，var()为方差，cov()为协方差。

结合上述公式(8)，可得到相关系数的表达公式(9)：

其中，ρ为pearson相关系数。

令相关系数ρ取最大值时，求解时固定分母来求分子，即该问题简化转换为一个优化问题，即优化公式为：

这个问题为有等式约束条件的优化问题，可以采用拉格朗日乘数法(lagrangemultiplier)求解。因此根据公式(10)构建的拉格朗日公式为：

其中，ρα和ρβ都是待求特征根。

分别对拉格朗日公式中的α、β求偏导，且令偏导分别等于0，具体详见公式(12)：

令ρλ＝ρα＝ρβ，同时防止矩阵不可逆对其正则化，上述问题可转换为求解公式(12)的特征根问题：

其中，ε1和ε2是正则化参数，取值范围[0，1]。

根据公式(13)求解不同的ρλ以及对应的典型向量α和β。

本发明依据nielsen和canty在低维空间中提出的概念，在var(u)＝var(v)＝1条件下，使var(u-v)最大的情况下，ui-vi中最后几个变量集中包含着两组数据之间的线性不变关系，将上述规则映射到高位空间进行求解，以选取图像间的非线性不变点，采用迭代的办法选取不变信息，因为mad变量是测量变量的线性组合，由中心极限定理可得它们呈现的高斯分布。每一个独立的mad变量互相正交，标准化后的像素方差符合自由度为n的分布，其具体权重公式为：

其中，z为标准化后的方差，ui为u中的第i个元素，vi是v中的第i个元素，ρλi为第i个特征根，为卡方检验，pr(nochange)为不变像素点的概率值，用来确定哪些像素可以最终被选为不变特征的，可以设定一个固定的域值t(比如t＝95％)，当pr(nochange)＞t时，判断为特征不变像素点，以此获取非线性辐射控制点。

步骤400：所述根据所述非线性辐射控制点对所述目标图像各阶段像素点进行校正，获得地表反射率值。

步骤401：采用最小二乘法对所述非线性辐射控制点进行非线性关系拟合获得非线性拟合公式。

拟合公式为：

f(r)＝a0+a1r+...+apr^p(14)；

其中，f(r)为待拟合的表达式，将辐射控制点的原始像素值输入该式可获取校正后的像素值，a0、a1...ap为待求解的常系数，r为非线性辐射控制点的原始像素值。

偏差平方最小公式为：

其中，δi为第i个非线性辐射控制点的偏差，ri为第i个非线性辐射控制点的原始目标图像像素值，f(ri)为将第i个非线性辐射控制点的原始目标图像像素值代入拟合的非线性关系表达式后，计算得到的像素值，si为第i个非线性辐射控制点的参考图像像素值，r为非线性辐射控制点的总个数。

令p为一个常数，将第i个非线性辐射控制点的原始目标图像像素值ri和第i个非线性辐射控制点对应的参考图像目标像素值si带入拟合公式(14)便可确定待求解的常系数a0、a1...ap，进而求解出具体拟合公式。

在实验中取不同p值，可以得到不同的非线性关系表达式，依据校正结果进行调整，选取校正结果最好时的p值。

将第i个非线性辐射控制点的原始目标图像像素值ri带入多个不同的拟合公式便可得出多个不同的f(ri)，再将多个不同的f(ri)和多个第i个非线性辐射控制点对应的参考图像目标像素值si带入偏差平方最小公式2便可求取所有非线性辐射控制点总体最小偏差。

根据最小偏差便可确定该拟合公式为非线性拟合公式。

步骤402：根据所述非线性拟合公式对所述目标图像各波段像素点进行校正，获得地表反射率，后续可作为地表反射率产品的应用。

根据上述方法获取的非线性关系表达式，对所述目标图像的各波段像素点带入非线性拟合公式便可求出校正后的像素点，所述校正后的像素点即为地表反射率。

本发明还提供了一种遥感图像相对辐射校正系统，以保证图像间的辐射一致性，以提高遥信数据的利用率。

图2为本发明实施例遥感图像相对辐射校正系统的结构示意图，如图2所示，所述系统包括：获取模块1、匹配模块2、处理模块3、非线性辐射控制点确定模块4、校正模块5。

获取模块1，用于获取目标图像。

匹配模块2，用于根据所述目标图像确定与所述目标图像相对应的参考图像。

根据与所述目标图像具有相同重叠区域、近似成像几何、坐标范围的图像确定为与所述目标图像相对应的参考图像，所述参考图像与所述目标图像无需时相相近，大大提高辐射处理图像间的辐射一致性，提高了遥信数据的利用率。

处理模块3，用于分别对所述参考图像、所述目标图像进行预处理后，获得对应所述参考图像的参考图像像素点集和对应所述目标图像的目标图像像素点集。

所述预处理包括坐标重投影、几何配准、计算地理重叠区域、各波段像素点对应、并使用掩膜技术排除极值点、过饱和点、云、云下阴影及水体的像素点。对所述参考图像、所述目标图像进一步处理，为了获得较为精确的参考图像像素点集和目标图像像素点集。

非线性辐射控制点确定模块4，用于基于核典型相关分析方法根据所述参考图像像素点集和所述目标图像像素点集确定非线性辐射控制点。

所述非线性辐射控制点确定模块4，具体包括：第一确定单元、第二确定单元。

第一确定单元分别将所述参考图像像素点集、所述目标图像像素点集映射到核空间，利用核典型相关分析方法确定参考图像核典型变量、目标图像核典型变量；

第二确定单元基于mad中的权重公式对所述参考图像核典型变量、所述目标图像核典型变量进行筛选，确定非线性辐射控制点。

所述权重公式为：

其中：u为映射到高维特征空间的目标图像核典型变量集合，ui为u的第i个元素，v为映射到高维特征空间的参考图像核典型变量集合，vi是v的第i个元素，var()为方差，z为标准化后的方差，ρλi为第i个特征根，为卡方检验，pr(nochange)为不变像素点的概率值，n为重叠区域的点的个数。

校正模块5，用于根据所述非线性辐射控制点对所述目标图像各阶段像素点进行校正，获得地表反射率值；所述校正模块具体包括：拟合单元、校正单元。

拟合单元根据偏差平方和最小公式对所述非线性辐射控制点进行非线性关系拟合获得非线性拟合公式；

所述偏差平方和最小公式为：

其中，δi为第i个非线性辐射控制点的偏差，ri为第i个非线性辐射控制点的原始目标图像像素值，f(ri)为将第i个非线性辐射控制点的原始目标图像像素值代入拟合的非线性关系表达式后，计算得到的像素值，si为第i个非线性辐射控制点对应的参考图像像素值，r为非线性辐射控制点的总个数。

校正单元根据所述非线性拟合公式对所述目标图像各波段像素点进行校正，获得地表反射率。

本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的系统而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。

本文中应用了具体个例对本发明的规则及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处。综上所述，本说明书内容不应理解为对本发明的限制。