基于小波变换的图像融合算法详解

图像融合的笔记【1】——小波融合

参考的以下大佬的文章：
小波融合的代码实例【1】
小波融合的代码实例【2】
小波融合代码详解

先是最基本的基于小波变换的图像融合，融合规则采用简单的低频部分取平均值，高频部分取基于区域能量。
我们使用的是matalb中自带的小波分解函数

- wavedec2函数

首先来了解matlab中小波分解函数中各个参数的作用：`

[C,S] = wavedec2(A,N,'wname');

其中A表示的是待处理的图像；
N表示分解的层数；
‘wname’表示的是使用的小波基函数；（具体小波基函数的选择可以查看help文档）
C为各层的分解系数
S表示的是各层分解系数长度的矩阵
C其实是一个**1*size（A）**的行向量；具体表现为：
C=[A(N)H(N)V(N)D(N)H(N-1)V(N-1)D(N-1)…H(1)V(1)D(1)]
其中A(N)是第N层的低频系数；H(N)V(N)D(N)则是第N层的高频系数，分别代表水平、数值以及对角线方向的高频系数。

盗两张图表示为：
生成矩阵C为：

S作为储存各层分解系数长度的矩阵，第一行是A(N)的长度（即A(N)原矩阵的行数和列数），第二行是H(N)V(N)D(N)的长度，以此类推最后一行为原图A的长度，所以S是一个N+2×2的矩阵。

搞清楚小波分解函数的参数意义，再去做图像融合就比较简单了。

function y=wtfusion(I1,I2,N,wname)

%函数功能：
%     函数x=wtfusion(I1,I2,N,wname)将两幅原图像I1,I2进行基于小波变换的图像融合，得到融合后的图像y
%     近似分量采用加权平均的融合规则，各细节分量采用基于区域特性量测的融合规则
%输入参数：
%     I1----原图像1
%     I2----原图像2
%     N----小波分解的层数
%     wname----小波基函数
%输出参数：
%     y----原图像融合后得到的图像
%-----------------------------------------------------------------%

I1=double(I1);                   %将uint8的图像数据类型转换成double型进行数据处理,这一步考虑的是计算时的精度问题。
I2=double(I2);

 %将原图像I1,I2分别进行N层小波分解，wname为小波基函数，
 %C为各层分解系数,S为各层分解系数长度,也就是图像的大小.
 %S(N+2行2列)的结构是储存各层分解系数长度的,即第一行是A(N)的长度（其实是A(N)的原矩阵的行数和列数）,
 %第二行是H(N)|V(N)|D(N)|的长度,第三行是H(N-1)|V(N-1)|D(N-1)的长度,
 %倒数第二行是H(1)|V(1)|D(1)长度,最后一行是X的长度(大小)
 %对于256×256的原图像而言，进行3层分解得到的S=[32 32；32 32；64 64； 128 128；256 256];
 %得到的C大小为1×（256×256）即1×65536；
[C1,S1]=wavedec2(I1,N,wname); 
[C2,S2]=wavedec2(I2,N,wname);  


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%小波分解系数的融合处理规则：采取低频系数取平均值；高频部分基于区域能量的方法；
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


KK=size（C1）;%得到C1矩阵的行、列数值；
coef_fuse = zeros(kk(1),KK(2));%根据建立一个与C1矩阵相同大小的融合系数矩阵；
Temp = zeros(1,2);%创建一个1×2的全0矩阵；


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%进行低频小波系数的融合：取平均值
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

A1=appcoef2(C1,S1,wname,N);            %提取出小波分解的近似分量
A2=appcoef2(C2,S2,wname,N);
A = (A1+A2)/2;%对小波分解得到的低频分量取平均值；
A = reshape(A,1,S1(1,1)*S1(1,2));

coef_fuse(1,1:(S1(1,1)*S1(1,2))) = A;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%进行高频小波系数的融合：基于区域能量
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

for i=N:-1:1
 [H1,V1,D1]=detcoef2('all',C1,S1,i);  %提取出小波分解的各层细节分量
 [H2,V2,D2]=detcoef2('all',C2,S2,i);
 H=f(H1,H2); %各层细节分量的融合规则采用基于区域特性量测的融合规则
 V=f(V1,V2);
 D=f(D1,D2);
 h=reshape(H,1,S1(N+2-i,1)*S1(N+2-i,2));%分别将融合后的细节分量转成行向量，并存入行向量coef_fuse中
 v=reshape(V,1,S1(N+2-i,1)*S1(N+2-i,2));
 d=reshape(D,1,S1(N+2-i,1)*S1(N+2-i,2)); 
 coef_fuse = [A,h,v,d];
end
%进行小波重构得到融合之后的图像

S=S1;
y=waverec2(coef_fuse,S,wname);      %重构原图像
figure(1);imshow(uint8(y));title('基于小波变换的融合图像')
end
function y=f(x1,x2)

%函数功能：
%       y=f(x1,x2)将两幅原图像x1和x2基于区域特性量测的融合规则进行融合,得到融合后的图像y
%       首先计算两幅图像的匹配度，若匹配度大于阈值，说明两幅图像对应局部能量较接近，
%       因此采用加权平均的融合方法；若匹配度小于阈值，说明两幅图像对应局部能量相差较大，
%       因此选取局部区域能量较大的小波系数作为融合图像的小波系数
%输入参数：
%      x1----输入原图像1
%      x2----输入原图像2
%输出参数：
%      y----融合后的图像
%------------------------------------------------------------%

w=1/16*[1 2 1;2 4 2;1 2 1];   %权系数
E1=conv2(x1.^2,w,'same');     %分别计算两幅图像相应分解层上对应局部区域的“能量”
E2=conv2(x2.^2,w,'same');
M=2*conv2(x1.*x2,w,'same')./(E1+E2);%计算两幅图像对应局部区域的匹配度
T=0.7;                              %定义匹配阈值,可以自己改变
Wmin=1/2-1/2*((1-M)/(1-T));
Wmax=1-Wmin;
[m,n]=size(M);

for i=1:m
    for j=1:n
        if M(i,j)<T                %如果匹配度小于匹配阈值，说明两幅图像对应局部区域能量距离较远；
            if E1(i,j)>=E2(i,j)    %那么就直接选取区域能量较大的小波系数
                y(i,j)=x1(i,j);
            else
                y(i,j)=x2(i,j);
            end
        else                       %如果匹配度大于匹配阈值，说明两幅图像对应局部区域能量比较接近；
            if E1(i,j)>=E2(i,j)    %那么就采用加权的融合算法
                y(i,j)=Wmax(i,j)*x1(i,j)+Wmin(i,j)*x2(i,j);
            else
                y(i,j)=Wmin(i,j)*x1(i,j)+Wmax(i,j)*x2(i,j);
            end
        end
    end
end

end