医学图像处理案例（十四）——基于小波变换的图像融合
医学图像处理案例（十四）——基于小波变换的图像融合
今天将简单介绍使用小波变换来对多模态图像进行融合。
1、图像融合概述
图像融合（Image Fusion）是指将多源信道所采集到的关于同一目标的图像数据经过图像处理和计算机技术等，最大限度的提取各自信道中的有利信息，最后综合成高质量的图像，以提高图像信息的利用率、改善计算机解译精度和可靠性、提升原始图像的空间分辨率和光谱分辨率，利于监测。
2、小波变换特点介绍
小波变换的固有特性使其在图像处理中有如下优点：完善的重构能力，保证信号在分解过程中没有信息损失和冗余信息；把图像分解成低频图像和细节（高频）图像的组合，分别代表了图像的不同结构，因此容易提取原始图像的结构信息和细节信息；小波分析提供了与人类视觉系统方向相吻合的选择性图像。
一般图像融合的小波分解采用离散小波变换(Discrete Wavelet Transform, DWT)。DWT的函数基由一个称为母小波或分析小波的单一函数通过膨胀和平移获得。因而，DWT同时具有时域和频域分析能力，与一般的金字塔分解相比，DWT图像分解具有以下优势：
1）具有方向性，在提取图像低频信息的同时，还可获得了水平、垂直和对角三个方向的高频信息；
2）通过合理的选择母小波，可使DWT在压缩噪声的同时更有效的提取纹理、边缘等显著信息；
3）金字塔分解各尺度之间具有信息的相关性，而DWT在不同尺度上具有更高的独立性。
3、基于小波变换的图像融合
DWT 融合算法基本思想：首先对源图像进行小波变换，然后按照一定规则对变换系数进行合并；最后对合并后的系数进行小波逆变换得到融合图像。
3.1、小波分解原理简介
LL：水平低频，垂直低频
LH：水平低频，垂直高频
HL：水平高频，垂直低频
HH：水平高频，垂直高频
其中，L表示低频，H表示高频，下标1、2表示一级或二级分解。在每一分解层上，图像均被分解为LL，LH，HH和HL四个频带，下一层的分解仅对低频分量LL进行分解。这四个子图像中的每一个都是由原图与一个小波基函数的内积后，再经过在x和y方向都进行2倍的间隔采样而生成的，这是正变换，也就是图像的分解；逆变换，也就是图像的重建，是通过图像的增频采样和卷积来实现的。
3.2、融合规则
规则一：系数绝对值较大法
该融合规则适合高频成分比较丰富，亮度、对比度比较高的源图像，否则在融合图像中只保留一幅源图像的特征，其他的特征被覆盖。小波变换的实际作用是对信号解相关，并将信号的全部信息集中到一部分具有大幅值的小波系数中。这些大的小波系数含有的能量远比小系数含有的能量大，从而在信号的重构中，大的系数比小的系数更重要。
规则二：加权平均法
权重系数可调，适用范围广，可消除部分噪声，源图像信息损失较少，但会造成图像对比度的下降，需要增强图像灰度。
4、基于小波变换的图像融合代码实现
我将分享matlab和python版本代码来融合红外和可见光图像，融合策略是低频图像采用平均值法，高频图像采用最大值法。python版本中需要用到PyWavelets库，可以使用下面命令来安装，具体可以见原文链接。
pip install PyWaveletspython版本代码：
import pywt
import cv2
import numpy as np
# This function does the coefficient fusing according to the fusion method
def fuseCoeff(cooef1, cooef2, method):
if (method == 'mean'):
cooef = (cooef1 + cooef2) / 2
elif (method == 'min'):
cooef = np.minimum(cooef1, cooef2)
elif (method == 'max'):
cooef = np.maximum(cooef1, cooef2)
return cooef
# Params
FUSION_METHOD = 'mean' # Can be 'min' || 'max || anything you choose according theory
FUSION_METHOD1 = 'max'
# Read the two image
I1 = cv2.imread('IR1.png', 0)
I2 = cv2.imread('VIS1.png', 0)
# First: Do wavelet transform on each image
wavelet = 'db2'
cooef1 = pywt.wavedec2(I1[:, :], wavelet, level=1)
cooef2 = pywt.wavedec2(I2[:, :], wavelet, level=1)
# Second: for each level in both image do the fusion according to the desire option
fusedCooef = []
for i in range(len(cooef1)):
# The first values in each decomposition is the apprximation values of the top level
if (i == 0):
fusedCooef.append(fuseCoeff(cooef1[0], cooef2[0], FUSION_METHOD))
else:
# For the rest of the levels we have tupels with 3 coeeficents
c1 = fuseCoeff(cooef1[i][0], cooef2[i][0], FUSION_METHOD1)
c2 = fuseCoeff(cooef1[i][1], cooef2[i][1], FUSION_METHOD1)
c3 = fuseCoeff(cooef1[i][2], cooef2[i][2], FUSION_METHOD1)
fusedCooef.append((c1, c2, c3))
# Third: After we fused the cooefficent we nned to transfor back to get the image
fusedImage = pywt.waverec2(fusedCooef, wavelet)
# Forth: normmalize values to be in uint8
fusedImage1 = np.multiply(np.divide(fusedImage - np.min(fusedImage), (np.max(fusedImage) - np.min(fusedImage))), 255)
fusedImage1 = fusedImage1.astype(np.uint8)
# Fith: Show image
cv2.imwrite("win.bmp", fusedImage1)matlab版本代码：
close all
clear
clc
%% Read Images
% the size of images must be equal
a=imread('IR1.png');
b=imread('VIS1.png');
%% Wavelet Transform
[a1,b1,c1,d1]=dwt2(a,'db2');
[a2,b2,c2,d2]=dwt2(b,'db2');
[k1,k2]=size(a1);
%% Fusion Rules
%% Average Rule
for i=1:k1
for j=1:k2
a3(i,j)=(a1(i,j)+a2(i,j))/2;
end
end
%% Max Rule
for i=1:k1
for j=1:k2
b3(i,j)=max(b1(i,j),b2(i,j));
c3(i,j)=max(c1(i,j),c2(i,j));
d3(i,j)=max(d1(i,j),d2(i,j));
end
end
%% Inverse Wavelet Transform
c=idwt2(a3,b3,c3,d3,'db2');
imshow(a)
title('First Image')
figure,imshow(b)
title('Second Image')
figure,imshow(c,[])
title('Fused Image')
fusedimage=(c-min(c(:)))/(max(c(:))-min(c(:)))*255;
fusedimage1=uint8(fusedimage);
imwrite(fusedimage1,'fusedimage.bmp')5、融合结果
下图是红外图像和可见光图像。
matlab融合结果
python融合结果
如果碰到任何问题，随时留言，我会尽量去回答的。