CN105160647A - 一种全色多光谱影像融合方法 - Google Patents

Abstract

本发明提供了一种全色多光谱影像融合方法，对全色影像Pan降质，将原始多光谱影像MS重采样，确定模拟全色影像Pan^*中各像素的基础融合系数，根据多光谱影像判断模拟全色影像Pan^*中各像素的地物类型，根据地物类型对基础融合系数进行区别性拉伸，得到最终融合系数，最后根据最终融合系数进行融合。本发明的一种全色多光谱影像融合方法，能够减小融合影像颜色失真程度；消除云雾、图幅大小等因素对融合质量的影响；平衡不同地物对色彩保真度和清晰度的要求。

Description

一种全色多光谱影像融合方法

技术领域

本发明涉及一种全色多光谱影像融合方法，属于遥感影像信息处理技术领域。

背景技术

当前，很多卫星遥感***都能够同时获取地面同一场景的全色影像和多光谱影像。全色影像的空间分辨率高但光谱分辨率低，多光谱影像虽然具有多个波段，但空间分辨率往往不能满足应用需求。因此，利用影像融合技术可以充分利用多源影像中互补信息，获取一个高空间分辨率的多光谱影像，具有重要的研究意义和应用价值。

当前的影像融合方法主要分为如下几类：第一类是基于算术运算的方法，利用内插后的多光谱影像与全色影像进行某种算术运算得到融合影像，如相加法、相乘法、Brovey方法等；第二类是基于成分替换的方法，通过对多光谱影像进行某种空间变换，然后利用直方图匹配后的全色影像对某成分进行替换，再进行反变换得到融合影像，其中HIS方法和PCA方法是最为常用的两种方法，另外还有快速光谱函数法、UNB方法、Gram-Schmidt方法等；第三类方法是基于多分辨率分析的方法，把影像转换到不同的分辨率层次，通过多分辨率分析实现融合处理，如Laplacian和小波金字塔分析；第四类方法是基于滤波的方法，如高通滤波融合法是先对全色影像进行高通滤波，然后把滤波后的高频信息叠加在多光谱影像上，提高其空间分辨率；第五类方法是混合分析方法，将光谱混合分解技术引入到了融合框架之中，典型的方法包括线性混合模型和统计混合模型方法。

上述融合方法都可以在一定程度上提高多光谱遥感影像的空间分辨率。然而，现有算法的融合准则存在不足，导致融合结果中依旧存在有一定程度的光谱畸变，同时现有算法没有考虑图幅大小及地物类型等因素对于拟合系数产生的影响，导致部分融合结果存在偏色，部分细节信息丢失等情况，对高分辨率影像的处理效果不佳，在云雾的影响下容易色彩失真；同时基于统计的融合方法对大图幅的处理会在局部区域无法做到最佳。

发明内容

为了解决现有技术的不足，本发明提供了一种全色多光谱影像融合方法，应用均值滤波来产生模拟全色影像，可以在现有算法的基础上更好地抑制光谱畸变；同时对系数进行拉伸等方法来增强融合影像的细节信息；此外通过判断地物类型来自适应地设置权重以平衡色彩保真度与清晰度。

本发明为解决其技术问题所采用的技术方案是：提供了一种全色多光谱影像融合方法，包括以下步骤：

(1)对全色影像Pan降质，得到模拟全色影像Pan^*；

(2)将原始多光谱影像MS重采样，得到多光谱影像MS′，使得多光谱影像MS′在长和宽具有和全色影像Pan相同数量的采样点；

(3)确定模拟全色影像Pan^*中各像素的基础融合系数；

(4)根据多光谱影像判断模拟全色影像Pan^*中各像素的地物类型；

(5)根据地物类型对基础融合系数进行区别性拉伸，得到最终融合系数；

(6)根据最终融合系数进行融合。

步骤(1)中，通过以下公式对全色影像Pan降质，得到模拟全色影像Pan^*：

其中，m和n分别表示全色影像Pan的长度和宽度，表示经过均值滤波后位于模拟全色影像Pan^*的(x,y)处的灰度值，f(i,j)表示全色影像Pan的(i,j)处的灰度值。

步骤(2)中，使用双线性内插法对原始多光谱影像MS进行重采样，以得到多光谱影像MS′：

首先，通过以下公式计算重采样后位于多光谱影像MS′中(dstY,dstX)处的像素在原始多光谱影像MS中的位置(srcY,srcX)：

srcX＝dstX*(srcWidth/dstWidth)

srcY＝dstY*(srcHeight/dstHeight)

其中srcHeight和srcWidth分别表示原始多光谱影像MS的长与宽，dstHeight和dstWidth分别表示重采样后的多光谱影像MS′的长与宽；

然后，将srcY分为整数部分i₂与小数部分u，将srcX分为整数部分j₂与小数部分v，通过以下公式对重采样进行描述：

$\begin{array}{c}\hat{f}(dstY,dstX)=f(srcY,srcX)=f(i_2+u,j_2+v)\\ =(1-u)(1-v)f(i_2,j_2)+(1-u)vf(i_2,j_2+1)+u(1-v)f(i_2+1,j_2)+uvf(i_2+1,j_2+1)\end{array}$

其中，表示重采样得到的多光谱影像MS′在(dstY,dstX)处的像素灰度值，f(srcY,srcX)表示原始多光谱影像MS中(dstY,dstX)的位置。

步骤(3)中，通过以下公式确定基础融合系数：

${\mathrm{\λ}}_{basic}(i,j)={\left\{\frac{pan(i,j)+{\mathrm{\μ}}_1\[pan(i,j)-{\mathrm{pan}}^{*}(i,j)\]}{{\mathrm{pan}}^{*}(i,j)}\right\}}^{{\mathrm{\μ}}_2}$

其中，μ₁和μ₂为设置的系数，0<μ₁<3，1<μ₂<3，λ_basic(i,j)表示融合影像在(i,j)处的基础融合系数，pan(i,j)表示全色影像Pan在(i,j)处的灰度值，pan^*(i,j)表示模拟全色影像Pan^*在(i,j)处的灰度值。

μ₁和μ₂均设置为2。

步骤(4)中，地物分为植被、水体和建筑用地三类；利用归一化植被指数NDVI检测植被：

$NDVI=\frac{G-NIR}{G+NIR}>{\mathrm{\&delta;}}_1$

其中，G是绿光波段，NIR是近红外波段，δ₁是设置的阈值，当NDVI大于δ₁则判断为植被；

利用归一化水体指数NDWI检测水体：

$NDWI=\frac{NIR-R}{NIR+R}>{\mathrm{\&delta;}}_2$

其中，R是红光波段，δ₂是设置的阈值，当NDWI大于δ₂则判断为水体；

植被与水体之外的地物判断为建筑用地。

δ₁设置为0.3，δ₂设置为0.05。

步骤(5)中，通过以下公式对基础融合系数进行区别性拉伸，得到最终融合系数：

λ(i,j)＝ωλ_basic(i,j)+1-ω

其中，λ(i,j)表示融合影像上(i,j)处的最终融合系数，ω为根据地物设置的权重值，若地物类型为建筑用地则ω设置为ω_building，若地物类型为植被则ω设置为ω_vegetation，若地物类型为水体则ω设置为ω_water，ω_building＞ω_vegetation＞ω_water。

ω_building＝0.8，ω_vegetation＝0.64，ω_water＝0.4。

步骤(6)中，通过以下公式进行融合：

${\mathrm{MS}}_{band}^{fusion}(i,j)={\mathrm{MS}}_{band}^{\&prime;}(i,j)*\mathrm{\&lambda;}(i,j)$

其中，表示融合影像在(i,j)处的灰度值，MS′_band(i,j)表示重采样得到的多光谱影像MS′在(i,j)处的灰度值。

本发明基于其技术方案所具有的有益效果在于：本文算法应用均值滤波来产生模拟全色影像，可以在现有算法的基础上更好地抑制光谱畸变；同时对系数进行拉伸等方法来增强融合影像的细节信息；此外通过判断地物类型来自适应地设置权重以平衡色彩保真度与清晰度。据多次实验，该方法实现的融合影像的相关系数平均可达95.20％，平均梯度可达13.84。作为对比，应用PCIGeomatics软件中的Pansharp融合模块所做的融合影像其相关系数平均为94.57％，平均梯度为9.02；应用ENVI软件中的Gram-Schmidt融合模块所做的融合影像其相关系数平均为94.44％，平均梯度为7.57。注：相关系数为评价色彩保真度的指标，系数越大说明色彩保持能力越好；平均梯度为评价清晰度指标，系数越大说明图像清晰度越好。

附图说明

图1是本发明的流程框图。

具体实施方式

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

本发明提供了一种全色多光谱影像融合方法，包括以下步骤：

(1)对全色影像Pan降质，得到模拟全色影像Pan^*：

通过以下公式对全色影像Pan降质，得到模拟全色影像Pan^*：

其中，m和n分别表示全色影像Pan的长度和宽度，表示经过均值滤波后位于模拟全色影像Pan^*的(x,y)处的灰度值，f(i,j)表示全色影像Pan的(i,j)处的灰度值。

(2)将原始多光谱影像MS重采样，得到多光谱影像MS′，使得多光谱影像MS′在长和宽具有和全色影像Pan相同数量的采样点：

使用双线性内插法对原始多光谱影像MS进行重采样，以得到多光谱影像MS′：

首先，通过以下公式计算重采样后位于多光谱影像MS′中(dstY,dstX)处的像素在原始多光谱影像MS中的位置(srcY,srcX)：

srcX＝dstX*(srcWidth/dstWidth)

srcY＝dstY*(srcHeight/dstHeight)

其中srcHeight和srcWidth分别表示原始多光谱影像MS的长与宽，dstHeight和dstWidth分别表示重采样后的多光谱影像MS′的长与宽；

然后，将srcY分为整数部分i₂与小数部分u，将srcX分为整数部分j₂与小数部分v，通过以下公式对重采样进行描述：

$\begin{array}{c}\hat{f}(dstY,dstX)=f(srcY,srcX)=f(i_2+u,j_2+v)\\ =(1-u)(1-v)f(i_2,j_2)+(1-u)vf(i_2,j_2+1)+u(1-v)f(i_2+1,j_2)+uvf(i_2+1,j_2+1)\end{array}$

其中，表示重采样得到的多光谱影像MS′在(dstY,dstX)处的像素灰度值，f(srcY,srcX)表示原始多光谱影像MS中(dstY,dstX)的位置。

(3)确定模拟全色影像Pan^*中各像素的基础融合系数：

通过以下公式确定基础融合系数：

${\mathrm{\&lambda;}}_{basic}(i,j)={\left\{\frac{pan(i,j)+{\mathrm{\&mu;}}_1\&lsqb;pan(i,j)-{\mathrm{pan}}^{*}(i,j)\&rsqb;}{{\mathrm{pan}}^{*}(i,j)}\right\}}^{{\mathrm{\&mu;}}_2}$

其中，μ₁和μ₂为设置的系数，0<μ₁<3，1<μ₂<3，μ₁和μ₂可以均设置为2，λ_basic(i,j)表示融合影像在(i,j)处的基础融合系数，pan(i,j)表示全色影像Pan在(i,j)处的灰度值，pan^*(i,j)表示模拟全色影像Pan^*在(i,j)处的灰度值。

(4)根据多光谱影像判断模拟全色影像Pan^*中各像素的地物类型：

地物分为植被、水体和建筑用地三类；利用归一化植被指数NDVI检测植被：

$NDVI=\frac{G-NIR}{G+NIR}>{\mathrm{\&delta;}}_1$

其中，G是绿光波段，NIR是近红外波段，δ₁是设置的阈值，可以设置为0.3，当NDVI大于δ₁则判断为植被；

利用归一化水体指数NDWI检测水体：

$NDWI=\frac{NIR-R}{NIR+R}>{\mathrm{\&delta;}}_2$

其中，R是红光波段，δ₂是设置的阈值，可以设置为0.05，当NDWI大于δ₂则判断为水体；

植被与水体之外的地物判断为建筑用地。

(5)根据地物类型对基础融合系数进行区别性拉伸，得到最终融合系数：

通过以下公式对基础融合系数进行区别性拉伸，得到最终融合系数：

λ(i,j)＝ωλ_basic(i,j)+1-ω

其中，λ(i,j)表示融合影像上(i,j)处的最终融合系数，ω为根据地物设置的权重值，若地物类型为建筑用地则ω设置为ω_building，若地物类型为植被则ω设置为ω_vegetation，若地物类型为水体则ω设置为ω_water，ω_building＞ω_vegetation＞ω_water，ω_building＝0.8，ω_vegetation＝0.64，ω_water＝0.4。

(6)根据最终融合系数进行融合：

通过以下公式进行融合：

${\mathrm{MS}}_{band}^{fusion}(i,j)={\mathrm{MS}}_{band}^{\&prime;}(i,j)*\mathrm{\&lambda;}(i,j)$

其中，表示融合影像在(i,j)处的灰度值，MS′_band(i,j)表示重采样得到的多光谱影像MS′在(i,j)处的灰度值。

本发明的一种全色多光谱影像融合方法，能够减小融合影像颜色失真程度；消除云雾、图幅大小等因素对融合质量的影响；平衡不同地物对色彩保真度和清晰度的要求。

Claims (10)

1.一种全色多光谱影像融合方法，其特征在于包括以下步骤：

(1)对全色影像Pan降质，得到模拟全色影像Pan^*；

(2)将原始多光谱影像MS重采样，得到多光谱影像MS′，使得多光谱影像MS′在长和宽具有和全色影像Pan相同数量的采样点；

(3)确定模拟全色影像Pan^*中各像素的基础融合系数；

(4)根据多光谱影像判断模拟全色影像Pan^*中各像素的地物类型；

(5)根据地物类型对基础融合系数进行区别性拉伸，得到最终融合系数；

(6)根据最终融合系数进行融合。

2.根据权利要求1所述的全色多光谱影像融合方法，其特征在于：步骤(1)中，通过以下公式对全色影像Pan降质，得到模拟全色影像Pan^*：

其中，m和n分别表示全色影像Pan的长度和宽度，表示经过均值滤波后位于模拟全色影像Pan^*的(x,y)处的灰度值，f(i,j)表示全色影像Pan的(i,j)处的灰度值。

3.根据权利要求1所述的全色多光谱影像融合方法，其特征在于：步骤(2)中，使用双线性内插法对原始多光谱影像MS进行重采样，以得到多光谱影像MS′：

首先，通过以下公式计算重采样后位于多光谱影像MS′中(dstY,dstX)处的像素在原始多光谱影像MS中的位置(srcY,srcX)：

srcX＝dstX*(srcWidth/dstWidth)

srcY＝dstY*(srcHeight/dstHeight)

其中srcHeight和srcWidth分别表示原始多光谱影像MS的长与宽，dstHeight和dstWidth分别表示重采样后的多光谱影像MS′的长与宽；

然后，将srcY分为整数部分i₂与小数部分u，将srcX分为整数部分j₂与小数部分v，通过以下公式对重采样进行描述：

$\begin{array}{c}\hat{f}(dstY,dstX)=f(srcY,srcX)=f(i_2+u,j_2+v)\\ =(1-u)(1-v)f(i_2,j_2)+(1-u)vf(i_2,j_2+1)+u(1-v)f(i_2+1,j_2)+uvf(i_2+1,j_2+1)\end{array}$

其中，表示重采样得到的多光谱影像MS′在(dstY,dstX)处的像素灰度值，f(srcY,srcX)表示原始多光谱影像MS中(dstY,dstX)的位置。

4.根据权利要求1所述的全色多光谱影像融合方法，其特征在于：步骤(3)中，通过以下公式确定基础融合系数：

${\mathrm{\&lambda;}}_{basic}(i,j)={\left\{\frac{pan(i,j)+{\mathrm{\&mu;}}_1\&lsqb;pan(i,j)-{\mathrm{pan}}^{*}(i,j)\&rsqb;}{{\mathrm{pan}}^{*}(i,j)}\right\}}^{{\mathrm{\&mu;}}_2}$

其中，μ₁和μ₂为设置的系数，0<μ₁<3，1<μ₂<3，λ_basic(i,j)表示融合影像在(i,j)处的基础融合系数，pan(i,j)表示全色影像Pan在(i,j)处的灰度值，pan^*(i,j)表示模拟全色影像Pan^*在(i,j)处的灰度值。

5.根据权利要求4所述的全色多光谱影像融合方法，其特征在于：μ₁和μ₂均设置为2。

6.根据权利要求1所述的全色多光谱影像融合方法，其特征在于：步骤(4)中，地物分为植被、水体和建筑用地三类；利用归一化植被指数NDVI检测植被：

$NDVI=\frac{G-NIR}{G+NIR}>{\mathrm{\&delta;}}_1$

其中，G是绿光波段，NIR是近红外波段，δ₁是设置的阈值，当NDVI大于δ₁则判断为植被；

利用归一化水体指数NDWI检测水体：

$NDWI=\frac{NIR-R}{NIR+R}>{\mathrm{\&delta;}}_2$

其中，R是红光波段，δ₂是设置的阈值，当NDWI大于δ₂则判断为水体；

植被与水体之外的地物判断为建筑用地。

7.根据权利要求6所述的全色多光谱影像融合方法，其特征在于：δ₁设置为0.3，δ₂设置为0.05。

8.根据权利要求1所述的全色多光谱影像融合方法，其特征在于：步骤(5)中，通过以下公式对基础融合系数进行区别性拉伸，得到最终融合系数：

λ(i,j)＝ωλ_basic(i,j)+1-ω

其中，λ(i,j)表示融合影像上(i,j)处的最终融合系数，ω为根据地物设置的权重值，若地物类型为建筑用地则ω设置为ω_building，若地物类型为植被则ω设置为ω_vegetation，若地物类型为水体则ω设置为ω_water，ω_building＞ω_vegetation＞ω_water。

9.根据权利要求8所述的全色多光谱影像融合方法，其特征在于：ω_building＝0.8，ω_vegetation＝0.64，ω_water＝0.4。

10.根据权利要求1所述的全色多光谱影像融合方法，其特征在于：步骤(6)中，通过以下公式进行融合：

${\mathrm{MS}}_{band}^{fusion}(i,j)={\mathrm{MS}}_{band}^{\&prime;}(i,j)*\mathrm{\&lambda;}(i,j)$

其中，表示融合影像在(i,j)处的灰度值，MS′_band(i,j)表示重采样得到的多光谱影像MS′在(i,j)处的灰度值。