基于非负稀疏矩阵的WorldView‑2遥感全色与多光谱图像融合方法

摘要

本发明提供了一种基于非负稀疏矩阵的WorldView‑2遥感全色与多光谱图像融合方法，涉及图像融合领域，通过使用一种非负矩阵分解算法对多光谱图像进行亮度分量提取，然后使用HCS变换对图像进行融合，得到融合图像，在细节注入和光谱保持方面都得到了一定的提升，最终得到了高质量的融合结果，由于对I分量的提取采用了NMF方法，提高了亮度分量的提取精度，较对比算法更合理，使得WV‑2卫星融合图像整体质量较高，在细节信息融入和光谱保持方面都有一定的提高，主观评价和客观分析结果能够达到一致，得到的融合图像可视性更好，图片更清晰。较传统的遥感全色与多光谱图像融合方法更有优势。

说明书

技术领域

本发明涉及图像融合领域，尤其是遥感全色与多光谱图像融合的一种方法。

背景技术

近年来，新型卫星遥感多光谱图像的波段数在不断增加，图像的分辨率也在快速地提高。以WorldView-2(WV-2)卫星为例，WV-2卫星于2009发射升空，能够提供8波段1.84米分辨率的多光谱图像和单波段0.46米分辨率的全色图像。与传统卫星图像相比具有以下特点：波段增加，光谱划分更细；全色图像的光谱覆盖范围变窄，使之与多光谱波段的光谱匹配发生较大的变化。在遥感应用上，往往需要同时具有高空间和高光谱分辨率的图像。图像融合技术就是利用高空间分辨率的全色图像去提高多光谱图像的空间分辨率，同时尽量保持多光谱图像的光谱特性不变。WV-2卫星图像代表了超高分辨率遥感图像的发展趋势，同时也为遥感图像的融合提出了更高的要求。正是由于这些变化，使现有的融合方法效果不佳。

目前很多遥感图像融合算法都有提取亮度分量的过程，例如，Brovey变换、HCS变换和MSFIM方法，亮度分量提取的好坏将直接影响融合结果，对算法的性能有很大的影响。其中，MSFIM算法是一种基于亮度平滑滤波调制(SFIM)算法的改进算法，但此改进算法虽然改善了细节信息的融入，却产生了较SFIM算法更差的光谱畸变。为降低MSFIM方法中的光谱畸变，需要改变亮度分量与全色图像的比值，使比值更接近1，即应该使亮度分量与全色图像的光谱响应特性更相似。其他对亮度分量的提取方法有均值法、加权均值法和各波段算术平均法。其实这些方法都是传统方法的演变，提取的亮度分量不够精确。

发明内容

现有的融合方法针对的是传统卫星遥感多光谱图像，对于超高分辨率遥感图像而言，并不是最优的图像融合方法，且现有的亮度分量提取方法并不能较好的解决细节融入和光谱畸变问题，使得提取的亮度分量不够精确。

为了克服现有技术的不足，本发明通过使用一种非负矩阵分解(Non-negative Matrix Factorization，NMF)的算法对多光谱图像进行亮度分量提取，然后使用HCS变换对图像进行融合，得到融合图像，本发明的融合图像在细节注入和光谱保持方面都得到了一定的提升，最终得到了高质量的融合结果。

本发明解决其技术问题所采用的技术方案包括如下步骤：

步骤1.使用非负矩阵分解法提取I分量

首先，将全色图像Pan和八波段的多光谱图像X₁,X₂,…,X₈按行拉直，得到P,M₁,M₂,…,M₈向量，再按公式(1)组成待分解矩阵V，即

V＝[P,M₁,M₂,...,M₈]>

其中，P,M₁,M₂,…,M₈分别为全色图像Pan及八波段多光谱图像X₁,X₂,…,X₈运算时将图像矩阵拉成相应的列向量；

其次，再令

[P,M₁,M₂,...,M₈]＝WH>

其中W为n*r矩阵，n为矩阵W的行数，r为矩阵W的列数，H为r*9矩阵，分解后W为一个列向量，将W恢复成图像矩阵即所得到的I分量；

步骤2.使用Pan分量匹配I分量

令

P^'2＝(Pan)²>

其中Pan为全色图像，即用P'变量代替Pan变量，然后，以下公式即使用P^'2分量匹配步骤1所得的I²分量：

$P^{\prime \prime 2}=\frac{{\sigma }_0}{{\sigma }_1}(P^{\prime 2}-{\mu }_1+{\sigma }_1)+{\mu }_0-{\sigma }_0---\left(4\right)$

其中μ₀、σ₀分别为I²的均值和标准方差，μ₁、σ₁分别为P^'2的均值和标准方差，P^”2为匹配后的I²分量；

$I_{adj}=\sqrt{P^{\prime \prime 2}}---\left(5\right)$

用I_adj分量代替P”分量来表示匹配后的I分量；

步骤3.使用超球面彩色空间分辨率融合HCS变换最终得到新的八波段分量

首先，对八波段的多光谱图像X₁,X₂,...,X₈进行HCS正变换得到相应的角度分量φ1,φ2,...,φ7，HCS正变换如下：

${\phi }_1=\arctan \left(\frac{\sqrt{{X_8}^2+{X_7}^2+\mathrm{...}+{X_2}^2}}{X_1}\right)$

${\phi }_6=\arctan \left(\frac{\sqrt{{X_8}^2+{X_7}^2}}{X_6}\right)$

${\phi }_7=\arctan \left(\frac{X_8}{X_7}\right)---\left(6\right)$

其次，对步骤2求得的匹配后的I分量I_adj和角度分量φ1,φ2,...,φ7做HCS反变换得到新的八波段分量X'₁,X'₂,…,X'₈，HCS反变换如下：

X₁'＝I_adjcosφ₁

X₂'＝I_adjsinφ₁cosφ₂

X₇'＝I_adjsinφ₁sinφ₂...sinφ₆cosφ₇

X₈'＝I_adjsinφ₁sinφ₂...sinφ₆sinφ₇>

步骤4.各波段融合

选取步骤3中X'₁,X'₂,…,X'₈中任意三个波段进行融合，即将三波段图像直接放入RGB三通道中，即可得到融合图像。

本发明的有益效果是由于对I分量的提取采用了NMF方法，提高了亮度分量的提取精度，较对比算法SFIM的改进算法MSFIM算法中I分量的提取算法更合理，使得WV-2卫星融合图像整体质量较高，在细节信息融入和光谱保持方面都有一定的提高，主观评价和客观分析结果能够达到一致，另外本发明方法较对比算法MSFIM算法得到的融合图像可视性更好，图片更清晰。较传统的遥感全色与多光谱图像融合方法更有优势。

附图说明

图1为本发明的技术路线框图。

图2为本发明多光谱图像X₂,X₃,X₅融合结果。

具体实施方式

下面结合附图和实施例对本发明进一步说明。

步骤1.使用非负矩阵分解(NMF)法提取I分量

非负矩阵分解作为一个相对成熟的数据分析手段，在图像分析、文本聚类、数据挖掘、语音处理等方面得到了广泛应用。考虑到从多光谱图像中提取的I分量要尽可能的与全色图像的光谱响应相一致，因此在NMF提取I分量的方法中用全色图像和八波段的多光谱图像组成待分解矩阵V。

首先，将全色图像Pan和八波段的多光谱图像X₁,X₂,…,X₈按行拉直，得到P,M₁,M₂,…,M₈向量，再按公式(1)组成待分解矩阵V，即

V＝[P,M₁,M₂,...,M₈]>

其中，P,M₁,M₂,…,M₈分别为全色图像Pan及八波段多光谱图像X₁,X₂,…,X₈运算时将图像矩阵拉成相应的列向量；

其次，再令

[P,M₁,M₂,...,M₈]＝WH>

其中W为n*r矩阵，n为矩阵W的行数，r为矩阵W的列数，H为r*9矩阵。分解后的W为一个列向量，将W恢复成图像矩阵即所得到的I分量。

其中，W为n*r矩阵，n为矩阵W的行数，r为矩阵W的列数，H为r*9矩阵，考虑到提取的I分量要和全色图像相匹配，此时应提取全色图像和多光谱图像的全局特征W，将分解后的W恢复成图像矩阵即所得到的I分量；

步骤2.使用Pan分量匹配I分量

令

P^'2＝(Pan)²>

其中Pan为全色图像，即用P'变量代替Pan变量，然后，以下公式即使用P^'2分量匹配步骤1所得的I²分量：

$P^{\prime \prime 2}=\frac{{\sigma }_0}{{\sigma }_1}(P^{\prime 2}-{\mu }_1+{\sigma }_1)+{\mu }_0-{\sigma }_0---\left(4\right)$

其中μ₀、σ₀分别为I²的均值和标准方差，μ₁、σ₁分别为P^'2的均值和标准方差，P^”2为匹配后的I²分量；

$I_{adj}=\sqrt{P^{\prime \prime 2}}---\left(5\right)$

用I_adj分量代替P”分量来表示匹配后的I分量；

步骤3.使用超球面彩色空间分辨率融合HCS变换最终得到新的八波段分量

面向WorldView-2卫星图像的融合方法HCS(Hyperspherical Color Space Resolution Merge)变换对波段数没有限制，因此适合多波段遥感图像融合。在HCS变换中，角度变量决定图像的光谱信息，而I决定图像的亮度信息，I变量的变化对图像的光谱信息没有影响，这是HCS变换的关键点。

首先，对八波段的多光谱图像X₁,X₂,...,X₈进行HCS正变换得到相应的角度分量φ1,φ2,...,φ7，HCS正变换如下：

${\phi }_1=\arctan \left(\frac{\sqrt{{X_8}^2+{X_7}^2+\mathrm{...}+{X_2}^2}}{X_1}\right)$

${\phi }_6=\arctan \left(\frac{\sqrt{{X_8}^2+{X_7}^2}}{X_6}\right)$

${\phi }_7=\arctan \left(\frac{X_8}{X_7}\right)---\left(6\right)$

其次，对步骤2求得的匹配后的I分量I_adj和角度分量φ1,φ2,...,φ7做HCS反变换得到新的八波段分量X'₁,X'₂,…,X'₈，HCS反变换如下：

X₁'＝I_adjcosφ₁

X₂'＝I_adjsinφ₁cosφ₂

X₇'＝I_adjsinφ₁sinφ₂...sinφ₆cosφ₇

X₈'＝I_adjsinφ₁sinφ₂...sinφ₆sinφ₇>

步骤4.各波段融合

选取步骤3中X'₁,X'₂,…,X'₈中任意三个波段进行融合，即将三波段图像直接放入RGB三通道中，即可得到融合图像。

为验证基于非负稀疏矩阵的WorldView-2遥感全色与多光谱图像融合新方法的正确性和优良性，进行以下实验。实验中所用的融合方法为HCS变换。两种融合方法的I分量提取方法分别采用MSFIM法和NMF法。实验数据为一组真实的WV-2图像，为2011年4月3日拍摄的Sydney Australia，为验证算法并能清晰的展示图像用于主观评价，实验数据采用图像的一部分，大小为300*300pixels，灰度级256。为获取评价结果时的参考图像，首先把多光谱图像下采样到原来的四分之一，再上采样到原来大小，这样原始多光谱图像可以作为参考图像。

图2中(a)是原始的全色图像下采样到原来图像的四分之一后的图像；图2(b)是原始的多光谱图像，用于实验结果对比。图2(c)是HCS变换融合方法中原始提取I分量方法的融合结果；图2(d)是本文用NMF方法提取I分量的融合结果。观察图2，图2(d)的光谱特征与图2(b)更为相近，即本文方法融合结果在光谱保持上比原始方法要好，由此可见，本发明所得到的融合结果在保持光谱特性不变的情况，细节信息的融入较传统的融合方法要好，且本发明的融合结果在光谱保持上也有一定的提高。

本发明选取空间相关系数(spatial correlate coefficient，sCC)、相关系数(correlate coeffcient，CC)、平均梯度(average gradient，AG)、信息熵(information entropy，IE)等常用客观评价指标对融合结果进行客观评价，这些指标值越大表明融合结果越好。

表1客观评价结果

Band1…Band8表示多光谱图像的八个波段。表1中可以看出本发明方法的各项评价指标均较大，这说明本发明在遥感图像融合中表现出了较传统方法更好的性能。