本发明属于遥感图像融合领域，涉及一种基于空间-光谱全变分正则化(spatial-spectraltotalvariation，sstv)的高光谱、多光谱和全色数据融合的方法，适用于由低分辨率的高光谱图像和低分辨率的多光谱图像与高分辨率的全色图像融合的得到高分辨率高光谱图像的应用场景。
背景技术：
在遥感领域，全色(panchromatic，pan)，多光谱(multispectral，ms)和高光谱(hyperspectral，hs)数据具有广泛的应用，包括土地分类，变化检测等。由于光学传感器的成本限制和数据存储与传输带宽的限制，对光学传感器的设计需要在空间分辨率和光谱分辨率之间做一个权衡。例如波段数多达数百个的高光谱数据提供了较高的光谱分辨率，但其空间画面一般较为模糊。具有单个波段的全色图像的光谱信息相对单一，但其空间信息较为丰富。出于对兼具高空间分辨率和光谱分辨率的遥感数据的应用需求，已经研究出多种图像融合技术。
多光谱融合技术将低分辨率的多光谱图像(lowresolutionmultispectralimage，lr-msi)与高分辨率的全色图像(highresolutionpanchromaticimage，hr-pan)融合，以提升多光谱图像的空间分辨率，这种融合技术当前发展较为成熟，主要分为成分替换的方法，多分辨率分析的方法，基于变分法和机器学习的方法。高光谱融合技术则通常将低分辨率的高光谱图像(lowresolutionhyperspectralimage，lr-hsi)与高分辨率的多光谱图像(highresolutionmultispectralimage，hr-msi)融合，以提升高光谱图像的空间分辨率，与多光谱图像融合技术相比，高光谱图像融合过程中使用的数据维度更高，需要估计的变量也更多。高光谱图像融合技术主要包括多光谱图像融合方法的扩展，基于贝叶斯的方法，和基于线性光谱分解的方法。
考虑到高光谱，多光谱和全色图像包含的信息为互补的，因此将拍摄同一画面的三种数据进行融合可以更好的挖掘其中包含的空间和光谱信息。这类数据融合方法的研究较少，为了更好的利用高光谱图像，多光谱图像和全色图像总的信息，本发明提出一种新的将空间-光谱正则化与高光谱、多光谱和全色数据融合相结合的基于模型的方法。
技术实现要素：
本发明针对现有遥感图像融合方法的不足，提供一种新的将空间-光谱正则化与高光谱、多光谱和全色数据融合相结合的基于模型的方法。
本发明所采用的技术方案是：一种将空间-光谱正则化与高光谱、多光谱和全色图像融合相结合的基于模型的方法。我们将hs图像利用图像分解的原理在低维子空间进行表示，把对融合目标的求解转为对低维张量系数矩阵的求解，而对系数矩阵的求解表述为包含三个二次数据拟合项和一个空间-光谱全变分正则化项的最小化目标函数。通过交替方向乘子法(alternativedirectionmultipliermethod，admm)制定并求解该模型，从而获得具有高分辨率的高光谱图像。
步骤1：选取同一时间段拍摄相同画面的经过配准等预处理的低分辨率高光谱图像、多光谱图像和全色图像。低分辨率的高光谱图像记为张量
步骤2：考虑到高光谱图像不同波段之间的具有较高的相关性，且在低维空间中比在融合目标
其中
(1)式等效于：
f(3)＝ex(3)(2)
其中
步骤3：构建基于子空间的关于全色图像，多光谱图像和高光谱图像三者融合的模型，得到由三个二次数据拟合项和一个正则化项组成的关于系数矩阵的目标能量函数。
步骤4：基于步骤3中的到的目标能量函数，通过admm来求解优化问题，从而得到系数矩阵x(3)，进而得到融合图像f(3)。
作为优选，步骤3的具体实现包括以下子步骤：
步骤3.1：高光谱图像h(3)可以看做f(3)的空间下采样版本：
h(3)＝f(3)bhsh(3)
其中矩阵
步骤3.2：多光谱图像m(3)可以看做f(3)经过空间下采样以及光谱下采样后的版本：
m(3)＝rmf(3)bmsm(4)
其中
步骤3.3：全色图像p(3)可以看做f(3)的光谱下采样版本：
p(3)＝rpf(3)(5)
其中
作为优选，步骤3.1中的bh矩阵，步骤3.2中的rm，bm矩阵和步骤3.3中的rp矩阵可由高光谱图像h(3)，多光谱图像m(3)，全色图像p(3)估算得到，我们使用的是名为hysure的高光谱图像融合方法中提出的估算方法。
步骤3.4：使用空间-光谱全变分正则化建立关于3d张量
其中d为对张量的梯度运算，
步骤3.5：基于式(1)-式(6)建立目标能量函数：
其中
作为优选，步骤4的具体实现包括以下子步骤：
步骤4.1：通过引入4个辅助变量o，u，v和w，(7)式中的优化问题表示为：
s.t.o＝x(3)bh
u＝x(3)bm
v＝x(3)
(8)式的增广拉格朗日函数表示为：
其中，y1，y2，y3和y4为尺度对偶变量(scaleddualvariable)，μ为惩罚参数，(9)式的优化问题可以被分解为6个子问题，我们将分别进行求解。
步骤4.2：求解x子问题：
其中右上角的t表示迭代次数。上述问题的解为：
步骤4.3：求解o子问题：
我们将o分为osh和
步骤4.4：求解u子问题：
我们将u分为usm和
步骤4.5：求解v子问题：
上述问题的解为：
步骤4.6：求解w子问题：
上述问题可通过软阈值直接求解。
步骤4.7：求解y子问题。该子问题可由以下式子直接更新：
步骤4.8：经过对上述子问题的迭代求解，得到系数张量
本发明主要针对由低分辨率高光谱图像与低分辨率多光谱和高分辨率全色图像融合得到高分辨率高光谱图像的应用需求。我们的融合方法将图像的子空间正则化与图像融合相结合，从而将对融合目标的求解转为对低维系数矩阵的求解。我们提出的模型包含三个二次数据拟合项和一个空间-光谱全变分正则化项，并通过admm算法进行求解，从而获得定性和定量方面均优于对比方法的高光谱融合图像。
附图说明
图1是本发明所提出方法与其他对比方法的融合图像的视觉质量比较。
图2是本发明的实验流程图。
具体实施方式
为了便于本领域普通技术人员理解和实施本发明，下面结合附图及实施例对本发明做进一步的详细描述，应当理解的是，此处所描述的实施施例仅用于解释本发明，并不用于限定本发明。
我们所使用的数据来自德国的rosis航空成像光谱仪拍摄的pavia大学市区。我们从该数据的115波段中选取了29个波段得到hr-hsi(320×320×29)，将其作为参考图像ref用于融合后进行质量评价。为了生成lr-hsi(80×80×29)，我们首先使用7×7高斯模糊(均值为零，标准差为2)对hr-hsi进行模糊处理，然后在hr-hsi的宽度和高度模式下每四个像素下采样一次。为了模拟同一场景的lr-msi(160×160×4)，使用类似于ikonos的反射光谱响应滤波器沿光谱模式下采样hr-hsi，然后使用7×7高斯模糊(均值为零，标准差为2)对hr-hsi进行模糊处理，最后在hr-hsi的宽度和高度模式下每两个像素下采样一次。生成同一场景的hr-pan(320×320)的方法为，先使用类似于ikonos的反射光谱响应滤波器沿光谱模式下采样hr-hsi为四个波段的图像，然后将前三个波段分别按照[0.114，0.587，0.299]权重进行线性组合。该融合实验的流程图如图2所示。
本实施例提供一种将空间-光谱正则化与高光谱、多光谱和全色图像融合相结合的基于模型的方法，包括以下步骤：
步骤1：选取同一时间段拍摄相同画面的经过配准等预处理的低分辨率高光谱图像、多光谱图像和全色图像。低分辨率的高光谱图像记为张量
步骤2：考虑到高光谱图像不同波段之间的具有较高的相关性，且在低维空间中比在
其中
(1)式等效于：
f(3)＝ex(3)(2)
其中
步骤3：构建基于子空间的关于全色图像，多光谱图像和高光谱图像的模型，得到由三个二次数据拟合项和一个正则化项组成的关于系数矩阵的目标能量函数。步骤3的具体实现包括以下子步骤：
步骤3.1：h(3)可以看做f(3)的空间下采样版本：
h(3)＝f(3)bhsh(3)
其中矩阵
步骤3.2：m(3)可以看做f(3)经过空间下采样以及光谱下采样后的版本：
m(3)＝rmf(3)bmsm(4)
其中
步骤3.3：p(3)可以看做f(3)的光谱下采样版本：
p(3)＝rpf(3)(5)
其中
作为优选，步骤3.1中的bh矩阵，步骤3.2中的rm，bm矩阵和步骤3.3中的rp矩阵可由高光谱图像h(3)，多光谱图像m(3)，全色图像p(3)估算得到，我们使用的估算方法来自于为名为hysure的高光谱图像融合方法。
步骤3.4：使用空间-光谱全变分正则化建立关于3d张量
其中d为对张量的梯度运算，
步骤3.5：基于式(1)-式(6)建立目标能量函数：
其中
步骤4：基于步骤3中的到的目标能量函数，通过admm来求解优化问题，从而得到系数矩阵x(3)，进而得到融合图像f(3)。步骤4的具体实现包括以下子步骤：
步骤4.1：通过引入4个辅助变量o，u，v和w，(7)式中的优化问题表示为：
s.t.o＝x(3)bh
u＝x(3)bm
v＝x(3)
(8)式的增广拉格朗日函数表示为：
y1，y2，y3和y4为尺度对偶变量(scaleddualvariable)。μ为惩罚参数。(9)式的优化问题可以被分解为6个子问题，我们将分别进行求解。
步骤4.2：求解x子问题：
其中右上角的t表示迭代次数。上述问题的解为：
步骤4.3：求解o子问题：
我们将o分为osh和
步骤4.4：求解u子问题：
我们将u分为usm和
步骤4.5：求解v子问题：
上述问题的解为：
步骤4.6：求解w子问题：
上述问题可通过软阈值直接求解。
步骤4.7：求解y子问题。该子问题可由以下式子直接更新：
步骤4.8：经过对上述子问题的迭代求解，得到系数张量
优化求解的流程如表1所示。
表1
步骤5：基于上述步骤我们将由高光谱图像
可以看出，由于我们的融合方法利用了全色图像中的信息，与仅使用hs和ms数据进行高光谱融合的其他融合结果相比，我们的融合结果清晰度更高，各项指标也最接近于理想值，在定性和定量方面均优于对比方法。
应当理解的是，本说明未详细阐述的部分均属于现有技术。
应当理解的是，上述针对实施例的描述较为详细，并不能因此而认为是对本发明专利保护范围的限制，本领域的普通技术人员在本发明的启示下，在不脱离本发明权利要求所保护的范围情况下，还可以做出替换或变形，均落入本发明的保护范围之内，本发明的请求保护范围应以所附权利要求为准。
1.一种基于空间-光谱全变分正则化的高光谱融合方法，其特征在于，包括如下步骤：
步骤1，选取同一时间段拍摄相同画面的经过配准预处理后的低分辨率高光谱图像、多光谱图像和全色图像；其中低分辨率的高光谱图像记为张量
步骤2，将融合目标分解为子空间矩阵和系数矩阵的乘积，将对融合目标
其中
(1)式等效于：
f(3)＝ex(3)(2)
其中
步骤3，构建基于子空间的关于全色图像，多光谱图像和高光谱图像三者融合的模型，得到由三个二次数据拟合项和一个正则化项组成的关于系数矩阵的目标能量函数；
步骤4，求解基于步骤3中得到的目标能量函数，从而得到系数矩阵x(3)，进而得到融合图像f(3)。
2.如权利要求1所述的一种基于空间-光谱全变分正则化的高光谱融合方法，其特征在于：步骤3中步骤3的具体实现包括以下子步骤：
步骤3.1，高光谱图像h(3)可以看做f(3)的空间下采样版本：
h(3)＝f(3)bhsh(3)
其中矩阵
步骤3.2，多光谱图像m(3)可以看做f(3)经过空间下采样以及光谱下采样后的版本：
m(3)＝rmf(3)bmsm(4)
其中
步骤3.3，全色图像p(3)可以看做f(3)的光谱下采样版本：
p(3)＝rpf(3)(5)
其中
步骤3.4，使用空间-光谱全变分正则化建立关于3d张量
其中d为对张量的梯度运算，
步骤3.5，基于式(1)-式(6)建立目标能量函数：
其中
3.如权利要求2所述的一种基于空间-光谱全变分正则化的高光谱融合方法，其特征在于：步骤3.1中的bh矩阵，步骤3.2中的rm，bm矩阵和步骤3.3中的rp矩阵可由高光谱图像h(3)，多光谱图像m(3)，全色图像p(3)估算得到，使用的是名为hysure的高光谱图像融合方法中提出的估算方法。
4.如权利要求1所述的一种基于空间-光谱全变分正则化的高光谱融合方法，其特征在于：步骤4中通过admm基于步骤3中得到的目标能量函数。
5.如权利要求2所述的一种基于空间-光谱全变分正则化的高光谱融合方法，其特征在于：步骤4的具体实现包括以下子步骤：
步骤4.1，通过引入4个辅助变量o，u，v和w，(7)式中的优化问题表示为：
s.t.0＝x(3)bh
u＝x(3)bm
v＝x(3)
(8)式的增广拉格朗日函数表示为：
其中，y1，y2，y3和y4为尺度对偶变量，μ为惩罚参数，(9)式的优化问题被分解为以下6个子问题；
步骤4.2，求解x子问题：
其中右上角的t表示迭代次数，上述问题的解为：
步骤4.3，求解o子问题：
将o分为osh和
步骤4.4，求解u子问题：
将u分为usm和
步骤4.5，求解v子问题：
上述问题的解为：
步骤4.6，求解w子问题：
上述问题通过软阈值直接求解；
步骤4.7，求解y子问题，该子问题可由以下式子直接更新：
步骤4.8，经过对上述子问题的迭代求解，得到系数张量
6.如权利要求1所述的一种基于空间-光谱全变分正则化的高光谱融合方法，其特征在于：步骤1中，全色图像与高光谱图像之间的分辨率之比为sh，全色图像与多光谱图像之间的分辨率之比为sm，即
