技术领域
[0001] 本
发明涉及遥感图像的
可视化增强处理技术，尤其是针对遥感多
光谱图像 和全色图像的图像融合方法
背景技术
[0002] 由于受到
传感器物理特性等技术的限制，目前的遥感图像无法同时获得高空 间
分辨率与高光谱分辨率，若获得较高的光谱分辨率，则常常需要牺牲一定的高 空间分辨率，反之亦然。为了解决这一问题，遥感平台经常需要搭载不同特性的 传感器来分别获取具有高空间分辨率的全色图像和高光谱分辨率的多光谱图像， 组成具有互补性的图像数据集。通过对互补数据集的处理即可得到同时具有高空 间分辨率和高光谱分辨率的遥感图像。然而，互补数据集的处理并不只是简单的 相加运算，由于互补数据集之间的相关性，互补数据集的处理更加强调信息的优 化，以突出游泳的专题信息，消除或抑制存在于互补数据集中的冗余信息，以提 高单幅图像的信息量。
[0003] 在过去几十年中，已出现许多遥感图像融合的方法，主要包括成分替换法和 多尺度分析法。经典的成分替换方法有HIS(Intensity-Hue-Saturation)变换、 GS(Gram Schmidt)变换、PCA(Principal Component Analysis)变换等，由于全色图 像和多光谱图像间复杂的光谱响应，成分替换法会伴随着严重的光谱扭曲。经典 的多尺度分析方法有ATWT(àtrous wavelet transform)和拉普拉斯金字塔 (Laplacian Pyramid)等，这类方法提取的空间细节可能与多光谱图像不匹配，因 此，融合图像存在一定的空间扭曲。近来，有学者提出了细节注入模型，并将主 成分替换以及多尺度分析都归在细节注入模型内，两者的区别主要在于细节信息 获取方式不同，前者是通过计算全色图像与多光谱图像强度分量之间的差别来获 取细节信息，而后者则是通过计算全色图像与降
采样的全色图像来获取细节信息。
发明内容
[0004] 为了克服
现有技术的不足，本发明提供一种结合引导滤波和IHS变换的图 像融合方法，预处理阶段对原始多光谱图像进行上采样和进行图像间的直方图匹 配，设计多尺度引导滤波提取全色图像的空间信息，并将得到的空间信息通过细 节注入模型插入到上采样的多光谱图像中，由于保留了原多光谱图像的光谱信息， 融合图像将同时具有高空艰难分辨率和高光谱分辨率。
[0005] 本发明解决其技术问题所采用的技术方案包括以下步骤：
[0006] 第一步，运用双三次插值法对包含N波段的原始多光谱图像MSi进行上采样， 上采样后的图像记为MSi；
[0007] 第二步，根据优化模型计算得到多光谱影像的强度分量
[0008] 第三步，使用步骤102得到的强度分量I对
全色影像进行直方图匹配，匹配 后的全色影像记为PAN；
[0009] 第四步，构建引导
滤波器，并将引导滤波
算法记为G(·)，对匹配后的全色影 像PAN进行多尺度滤波，得到全色影像的空间细节信息PD；
[0010] 第五步，计算细节插入模型的插入系数
[0011] 第六步，将步骤102得到的空间细节信息插入到上采样的多光谱影像中，得 到最终的融合影像MSi＝MSi+giPD。
[0012] 进一步的，步骤二具体为：根据公式 计算多光谱的
亮度分量， 并通过求解以下优化模型得到加权系数αi：
[0013]
[0014] 进一步的，步骤四具体如下所示：
[0015] 4a)构建引导滤波器 用G(·)表示引导滤波，其中ωk是 以像素k为中心的局部窗口，系数ak和bk通过最小化代价函数求解，ε为正则化 系数：
[0016]
[0017] 4b)对全色影像进行两层滤波分解：将直方图匹配后的全色影像PAN作为初 始引导滤波的输入影像，I分量作为引导影像，输出影像为全色影像的低频分量， 全色影像的空间信息则通过计算输入影像与输出影像的差值得到：
[0018] Pl＝G(Pl-1,I)
[0019] PDl＝Pl-1-Pl
[0020] 其中，Pl是输入影像l次滤波的输出图像，Pl-1是输入影像l-1次滤波的输出 图像，当l为1时，Pl-1为直方图匹配后的全色影像，PDl是l次滤波得到的全色 影像高频信息。
[0021] 进一步的，步骤六具体为通过公式 得 到最终融合影像。其中，m表示多光谱影像的波段数MSj表示多光谱影像的第j个波段，MSj表示融合影像的第j个波段，PD为空间细节信息。
[0022] 本发明的有益效果是：引导滤波可以充分挖掘图像之间的相关性，利用多光 谱的强度分量作为引导滤波对全色图像进行引导滤波，可以利用图像间的相关性， 减少细节信息的冗余性，从而避免信息的过度插入现象。采用细节注入模型可以 最大程度上的保留原有多光谱图像的光谱信息，从而在增强空间分辨率的同时保 持较高的光谱分辨率，得到高
质量的融合图像，是一种适用于高分辨率多光谱图 像与全色图像融合的有效融合方法。
附图说明
具体实施方式
[0024] 下面结合附图和
实施例对本发明进一步说明，本发明包括但不限于下述实施 例。
[0025] 本发明包括以下步骤：
[0026] 假设原始多光谱图像(MSi)包含N波段，下标i＝1,2,...,N，表示对应多波段 图像的第i波段。
[0027] 第一步、上采样多光谱图像
[0028] 运用双三次插值法对原始多光谱图像进行上采样，得到与全色图像匹配的上 采样多光谱图像，并记为MSi。
[0029] 第二步、计算多光谱图像的强度分量I
[0030] 计算强度分量I，I可以由合适的上采样多光谱图像MSi的线性组合模型估 计得到：
[0031]
[0032] 其中线性组合系数可由求解以下优化模型得到：
[0033]
[0034] 第三步、直方图匹配
[0035] 将全色图像(PAN)与第二步得到的强度分量I进行直方图匹配，匹配后的全 色图像记为PAN，直方图匹配具体计算过程如下：
[0036]
[0037] 式(3)中，μp、μI分别是全色图像PAN和强度分量I的均值，而δp、δI分别 是全色图像PAN和强度分量I的方差。
[0038] 第四步、引导滤波
[0039] 引导滤波器包含一幅输入图像Q，一幅引导图像P。滤波输出图像记为N， 用G(·)表示引导滤波。在以
像素k为中心的窗口ωk中，本方法采用正方形窗口， 窗口半径为r，窗口大小为(2r+1)×(2r+1)。N由P经过线性变换得到：
[0040]
[0041] 其中，系数ak和bk通过最小化代价函数求解，ε为正则化系数：
[0042]
[0043] 对全色影像进行两层滤波分解：将直方图匹配后的全色影像PAN作为初始 引导滤波的输入影像，I分量作为引导影像，输出影像为全色影像的低频分量， 全色影像的空间信息则通过计算输入影像与输出影像的差值得到：
[0044] Pl＝G(Pl-1,I) (6)
[0045] PDl＝Pl-1-Pl (7)
[0046] 其中，Pl是输入影像l次滤波的输出图像，Pl-1是输入影像l-1次滤波的输出 图像，当l为1时，Pl-1为直方图匹配后的全色影像，PDl是l次滤波得到的全色 影像高频信息。
[0047] 第五步、细节插入模型
[0048] 由于细节插入模型是将细节信息插入到上采样的多光谱图像的每个波段中， 模型公式如下所示：
[0049] MSi＝MSi+giPD (8)
[0050] 其中，插入系数应根据多光谱图像每个波段信息分别计算，即：
[0051]
[0052] 第六步、融合影像
[0053] 由于第四步已得到细节信息，第五步得到了插入系数gi，故而融合影像MSi可由式(8)获得。
[0054] 方法实施例：
[0055] 采用资源三号卫星遥感多光谱图像和全色图像，多光谱图像包含红、蓝、绿 以及
近红外4个波段，全色图像为单波段。多光谱图像的空间分辨率为5.8m， 大小为128×128，全色图像的空间分辨率为2m，大小为512×512。实施本发明 包括以下步骤：
[0056] 第一步、上采样多光谱图像
[0057] 由于全色图像的空间分辨率不同于多光谱
图像空间分辨率，故采用的多光谱 图像和全色图像大小也不相同，需要将多光谱图像上采样得到和全色图像相同的 图像。
[0058] 这里采用双三次插值法对原始多光谱图像进行上采样，上采样后的多光谱图 像记为MSi，上采样后的图像大小为512×512。
[0059] 第二步、计算多光谱强度分量
[0060] 计算强度分量I，I可以由合适的上采样多光谱图像MSi的线性组合模型估 计得到：
[0061]
[0062] 其中线性组合系数可由求解以下优化模型得到：
[0063]
[0064] 第三步、直方图匹配
[0065] 利用公式(3)，将全色图像与强度分量I进行直方图匹配，匹配后的图像记为 PAN。
[0066] 第四步、引导滤波
[0067] 用G(·)表示引导滤波，其中输入图像和引导图像均为单波段图像，使用引 导滤波对PAN进行多尺度滤波：
[0068] P1＝G(PAN,I)
[0069] P2＝G(P1,I)
[0070] 其中，强度分量I作为引导滤波的引导图像，第一层滤波的输入影像为PAN， 第二层输入图像为第一层滤波结果P1。在该两层上，分别用输入图像减去多尺度 引导滤波输出即可得到细节信息：
[0071] PD1＝PAN-P1
[0072] PD2＝P1-P2
[0073] 故细节信息为PD＝PD1+PD2
[0074] 第五步、计算插入系数
[0075] 利用公式(9)，得到多光谱图像各波段的插入系数gi
[0076] 第六步、融合图像
[0077] 利用公式(8)，得到融合图像的各个波段信息。
