基于PAN调制和多元线性回归的MS与PAN图像融合方法与流程

本发明涉及遥感图像处理技术领域，特别是一种基于PAN调制和多元线性回归的MS与PAN图像融合方法。

背景技术：

目前，大量在轨的卫星(例如Landsat 7ETM+、IKONOS、QuickBird、SPOT-5和WorldView-2/3等)，同时提供高空间分辨率的全色波段(PAN)和低空间分辨率的多光谱波段(MS)遥感图像。由于大量的应用需要使用高空间分辨率的多光谱图像，因此破切地需要将MS与PAN图像进行融合处理，以获得空间分辨率增强的MS图像，以应用于遥感图像解译、地表覆盖分类、目标检测等应用中。近年来国内外研究人员发展了大量遥感融合技术来融合MS和PAN图像以得到高空间分辨率的MS图像。

基于PAN调制的融合方法(例如PCI Geomatica软件中的PANSHARP方法)是目前使用较多的融合方法之一。基于PAN调制的MS与PAN融合基于假设：融合后MS波段与原始MS波段的比值等于原始PAN图像与合成PAN图像的比值。由于该计算公式暗含了光谱失真最小化模型，因此融合图像具有较好的目视效果和较小的光谱失真。基于PAN调制技术的融合算法主要有Brovey变换、Pradines’、合成变化比率(Synthetic Variable Ratio)、Smoothing Filter-based Intensity Modulation、PANSHARP(PS)、Haze-and Ratio-based (HR)等。基于PAN调制的融合方法，由于具有计算简单、鲁棒性强等优点，被广泛地应用于卫星数据的融合。待融合图像的配准精度、大气辐射影响、混合像元问题等是影响融合图像质量的几个关键问题。

Jing and Cheng(2009)提出的Haze-and Ratio-based(HR)方法采用原始PAN图像进行融合；且HR方法所使用的P_S的生成方式为：首先对原始PAN分辨率图像采用平均法降空间分辨率，然后再上采样到原始PAN空间分辨率。该方式对MS和PAN图像存在少量空间错位的情况无法改善效果。

针对相关技术中的问题，目前尚未提出有效的解决方案。

技术实现要素：

针对相关技术中的上述技术问题，本发明提出了一种基于PAN调制和多元线性回归的MS与PAN图像融合方法，可减少MS和PAN图像的少量空间错位对融合图像质量的影响。

本发明的目的通过以下技术方案实现：

一种基于PAN调制和多元线性回归的MS与PAN图像融合方法，包括以下步骤：

S1：将原始低分辨率MS图像I^L采用立方卷积的方式上采样到原始PAN图像P的分辨率，得到上采样MS图像I；

S2：将原始PAN图像采用平均法下采样到MS图像空间分辨率，得到MS空间分辨率的PAN图像P_L；

S3：采用各个MS波段的加权和以合成MS空间分辨率的PAN图像P_S；

S4：计算PAN波段和MS图像第i波段的雾气值H_P和H_i，其中，PAN波段的雾气值为H_p＝min(P)；MS图像第i个波段的雾气值为H_i＝min(MS_i)：

S5：根据空间信息增强程度系数k，得到对比度增强PAN图像P_F；

S6：根据空间信息增强程度系数k的不同取值，对上采样MS图像I、合成PAN图像P_S和对比度增强PAN图像P_F进行融合。

进一步的，所述步骤S3中，以P_L基为因变量，I^L为自变量，采用最小二乘法公式(1)求解各个MS波段的回归系数a_i和b：

利用公式(2)计算上采样MS各个波段的加权和以合成PAN图像P_S：

公式(1)和公式(2)中，N为MS图像的波段数量，为原始MS图像的第i波段，I_i为上采样MS图像的第i波段，P_S为合成PAN图像，a_i和b分别为第i波段的系数和常数项。

进一步的，所述步骤S5中，对比度增强图像P_F由以下公式获得，阈值T的值设置为P_F的方差；

公式(3)中P_E为对原始PAN图像进行Laplacian(拉普拉斯)滤波器滤波后的边缘细节图像。

进一步的，生成边缘细节图像P_E采用以下Laplacian滤波器：

公式(4)中，g为Laplacian滤波器。

进一步的，所述步骤S6中：

对于PAN图像中灰度大于或等于阈值T的像元(m，n)，利用以下公式(5)计算其融合光谱F_i；

对于PAN图像中灰度小于阈值T的像元(m，n)，根据以下公式(6)、(7)计算相对低的雾气值和

公式(6)和公式(7)中，0.5≤p≤1，和分别为图像中暗像元的PAN波段和第i波段的雾气值；

利用以下公式(8)计算其融合光谱F_i(m，n)：

公式(5)和公式(8)中，k为空间信息增强程度参数，P为原始PAN波段对比度增强后的图像，P_F为PAN波段对比度增强后的图像。

本发明的有益效果：显著降低融合图像的光谱失真，并能提供不同程度空间细节增强的融合图像，特别是在MS和PAN像元存在配准误差的情况下，该方法融合图像具有较少的光谱失真和较好的目视效果。

附图说明

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

图1是根据本发明实施例所述的一种基于PAN调制和多元线性回归的MS与PAN图像融合方法的流程图。

具体实施方式

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

如图1所示，本发明实施例所述的一种基于PAN调制和多元线性回归的MS与PAN图像融合方法，包括以下步骤：

S1：将原始低分辨率MS图像(I^L)采用立方卷积方式上采样到原始PAN图像(P)的分辨率，得到上采样MS图像I；

S2：将P采用平均法下采样到MS图像空间分辨率，得到MS空间分辨率的PAN图像P_L；

S3：以P_L基为因变量，I^L为自变量，采用最小二乘法公式(1)求解各个MS波段的回归系数a_i和b：

求解回归系数a_i和b，利用以下公式生成合成PAN图像P_S：

公式(1)和公式(2)中，N为MS图像的波段数量，为原始MS图像的第i波段，I_i为上采样MS图像的第i波段，P_S为合成PAN图像，a_i和b分别为第i波段的系数和常数项。

S4：计算PAN波段和MS图像第i波段的雾气值H_P和H_i，分别由PAN波段和MS图像第i波段的最小值确定，即PAN波段的雾气值为H_P＝min(P)，MS图像第i个波段的雾气值为H_i＝min(MS_i)。

S5：根据空间信息增强程度系数k，得到对比度增强PAN图像P_F。对比度增强图像P_F由以下公式获得，阈值T的值设置为P_F的方差；

公式(3)中P_E为对原始PAN图像进行Laplacian(拉普拉斯)滤波后的边缘细节图像，生成边缘细节图像P_E采用以下Laplacian滤波器：

公式(4)中，g为Laplacian(拉普拉斯)滤波器。

S6：根据空间信息增强程度系数k的不同取值，对上采样MS图像I、合成PAN图像P_S和对比度增强PAN图像P_F进行融合。

对于PAN图像中灰度大于或等于阈值T的像元(m，n)，利用以下公式计算其融合光谱F_i；

对于PAN图像中灰度小于阈值T的像元(m，n)，根据以下公式计算相对低的雾气值和

公式(6)和公式(7)中，0.5≤p≤1，和分别为图像中暗像元的PAN波段和第i波段的雾气值，利用以下公式计算其融合光谱F_i(m，n)：

公式(5)和公式(8)中，k为空间信息增强程度参数，P为原始PAN波段对比度增强后的图像，P_F为PAN波段对比度增强后的图像。

为了方便理解本发明的上述技术方案，以下通过具体使用方式对本发明的上述技术方案进行详细说明。

以下内容涉及了2个对比实验。第一个实验是在MS和PAN图像无配准误差的情况下进行，而第二个实验是在MS和PAN图像存在不同程度配准误差的情况下进行。

第一个对比实验中，实验数据包括来3个分别自WorldView-2、Pleiades、IKONOS等传感器的高分辨率遥感图像；对比方法选择了5种目前公认比较优秀的融合算法，具体包括Gram-Schmidt(GS)方法的Mode 1(以下简称为GS1)和Model 2(以下简称为GS2)、Adaptive GS(GSA)、Generalized Gaussian Generalized(GLP)、“à trous”wavelet transform(ATWT)和Additive Wavelet Luminance Proportional(AWLP)；融合图像质量评价指标选用了相对全局维数综合指标(ERGAS)、光谱角(SAM)、综合质量指数Q4/Q8、和空间相关系数(SCC)。其中，EASE反映了融合图像跟参考图像的偏差，值越小融合效果越好；ERGAS反映了融合图像跟参考图像的全局光谱辐射变形误差，越小越好；SAM反映融合图像光谱跟参考图像光谱的差异，值越小融合效果越好；Q4与Q8是同时考虑了融合图像跟参考图像的局部均值偏差、对比度变化以及相关性丢失情况的综合质量指标，值越大越好；SCC是考虑了融合图像跟PAN图像空间细节相关性的指标，值越大越好。3个实验图像的融合图像质量评价指标的统计见表1-表3。

表1 WorldVIew-2卫星图像的融合图像的质量评价指标统计

表2 Pleiades卫星图像的融合图像的质量评价指标统计

表3 IKONOS卫星图像的融合图像的质量评价指标统计

从表1-3中的统计指标可见，从光谱质量评价指标(RASE、ERGAS、SAM和Q4/Q8等)和空间质量评价指标(SCC)来看，本发明方法在空间细节调节参数k为零时，都优于其它融合方法。这表明，在不刻意进行空间细节增强的情况下(即k＝0)，本发明方法的融合图像在光谱指标和空间指标上都显著优于GS1、GS2、GSA、GLP、ATWT和ATWP等方法的融合图像。尽管在本发明方法中k值取1或2时，融合图像的RASE、ERGAS、SAM、Q4/Q8、SCC等质量评价指标并不完全优于其它方法，但从目视比较上看，融合图像没有明显的光谱失真，且更好地增强了空间细节。这表明，当融合图像用于遥感制图、目视解译等应用时，选用本发明方法是较好的选择。

在第二个对比实验中，实验数据为第一个对比实验中的IKONOS数据；实验方法为将上采样MS图像分别进行不同程度的偏移后，采用本发明方法、GSA、GLP、ATWT、AWLP等方法进行融合，最后对融合图像进行质量评价；融合图像质量评价选用了ERGAS、SAM、Q4和SCC。在本实验中，偏移量选择了(0，1)、(1，1)、(2，2)、(3，2)、(3，3)、(4，3)、(4，4)共七种情况。本实验中生成的融合图像质量评价指标的统计见表4。

从表4中的统计指标可见，在不同偏移量情况下，本发明方法的融合图像从ERGAS、SAM和Q4光谱质量评价指标和空间质量评价指标SCC来看均优于GSA、GLP、ATWT和ATWP等方法的融合图像。这表明，在MS和PAN图像存在空间错位的情况下，本发明方法跟其它方法相比具有显著优势。对各个融合图像进行目视比较上看，本发明方法的融合图像具有较小的光谱失真和较清晰的边界。

由上述两个实验可见，本发明方法跟同类方法相比，进一步降低了融合图像光谱失真并增强了空间细节，特别是在MS和PAN图像存在空间错位的情况下，融合图像从光谱、空间质量评价指标和目视效果上，都优于对比方法。

表4 IKONOS卫星图像在不同偏移情况下的融合图像的质量评价指标统计

本发明提出的融合方法所使用的MS空间分辨率的PAN图像P_S是采用各个MS波段的加权和得到，其中，权重系数为对PAN波段和各MS波段进行多元线性回归时的最优系数。在MS和PAN存在轻微空间错位的情况，最优系数较好地反应了PAN和各MS波段之间的关系，因此能减少错位对融合图像质量的影响。

本发明方法是一种基于PAN调制的、考虑雾气影响的融合方法。基于PAN调制的融合方法在空间细节注入的过程利用了光谱失真最小化模型。因此能限制融合图像光谱失真，保证融合图像的目视效果。由于大气程辐射的影响，图像中每个波段的雾气值会影响融合像元光谱矢量的方向，因此影响了融合图像的光谱失真程度。因此在PAN调制之前去除雾气影响能降低融合图像光谱失真。本方法在实施过程中对图像中相对暗的像元，采用了稍微低的雾气值，因此能改善水体、阴影等融合像元的光谱失真现象。

本发明方法可通过对空间细节增强系数k的设置，对原始PAN图像进行不同程度的对比度增强，然后用对比度增强后的PAN图像进行融合。该方式有利于用户根据应用需求，自主控制融合图像的空间细节增强程度，获得满足不用应用需求的融合图像。

综上所述，借助于本发明的上述技术方案，考虑了雾气影响的PAN调制方法，显著降低融合图像的光谱失真，并能提供不同程度空间细节增强的融合图像，特别是在MS和PAN像元存在配准误差的情况下，该方法融合图像具有较少的光谱失真和较好的目视效果。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。