%%%%% 基于小波的图像融合
clear
clc
I1=imread('F:\else\融合图像/2.jpg');
I2=imread('F:\else\融合图像/3.jpg');
I1=double(rgb2gray(I1));
[m,n]=size(I1);
I2=double(rgb2gray(I2));
figure;
subplot(1,3,1);
imshow(uint8(I1))
subplot(1,3,2);
imshow(uint8(I2))
[c1,s1]=wavedec2(I1,2,'haar'); %小波分解
%提取第二幅图的分解系数
ca11=appcoef2(c1,s1,'haar',2);
ch12=detcoef2('h',c1,s1,2);
cv12=detcoef2('v',c1,s1,2);
cd12=detcoef2('d',c1,s1,2);
ch11=detcoef2('h',c1,s1,1);
cv11=detcoef2('v',c1,s1,1);
cd11=detcoef2('d',c1,s1,1);
[c2,s2]=wavedec2(I2,2,'haar'); %小波分解
%提取第二幅图的分解系数
ca21=appcoef2(c2,s2,'haar',2);
ch22=detcoef2('h',c2,s2,2);
cv22=detcoef2('v',c2,s2,2);
cd22=detcoef2('d',c2,s2,2);
ch21=detcoef2('h',c2,s2,1);
cv21=detcoef2('v',c2,s2,1);
cd21=detcoef2('d',c2,s2,1);
%构造第一幅图的小波分解系数
facI1(1:s1(1,1),1:s1(1,2))=ca11;
facI1((s1(1,1)+1):(s1(1,1)*2),1:s1(1,2))=cv12;
facI1(1:s1(1,1),(s1(1,2)+1):(s1(1,2)*2))=ch12;
facI1((s1(1,1)+1):(s1(1,1)*2),(s1(1,2)+1):(s1(1,2)*2))=cd12;
facI1((s1(3,1)+1):(s1(3,1)*2),1:s1(3,2))=cv11;
facI1(1:s1(3,1),(s1(3,2)+1):(s1(3,2)*2))=ch11;
facI1((s1(3,1)+1):(s1(3,1)*2),(s1(3,2)+1):(s1(3,2)*2))=cd11;
%构造第二幅图的小波分解系数
facI2(1:s2(1,1),1:s2(1,2))=ca21;
facI2((s2(1,1)+1):(s2(1,1)*2),1:s2(1,2))=cv22;
facI2(1:s2(1,1),(s2(1,2)+1):(s2(1,2)*2))=ch22;
facI2((s2(1,1)+1):(s2(1,1)*2),(s2(1,2)+1):(s2(1,2)*2))=cd22;
facI2((s2(3,1)+1):(s2(3,1)*2),1:s2(3,2))=cv21;
facI2(1:s2(3,1),(s2(3,2)+1):(s2(3,2)*2))=ch21;
facI2((s2(3,1)+1):(s2(3,1)*2),(s2(3,2)+1):(s2(3,2)*2))=cd21;
%计算第一幅图的显著性
S1=zeros(m,n);
aveX1=sum(I1(:))/(m*n);
detX1=zeros(m,n);
detY1=zeros(m,n);
M=3;N=3;
for i=1:m-1
for j=1:n-1
detX1(i,j)=abs(I1(i,j)-I1(i,j+1));
detY1(i,j)=abs(I1(i,j)-I1(i+1,j));
end
end
for iImg=2:(m-1)
for jImg=2:(n-1)
firS=(1/(M*N))^2;
secS=0;
for i=(iImg-1):(iImg+1)
for j=(jImg-1):(jImg+1)
secS=secS+I1(i,j)-aveX1;
end
end
thiS=0;
for i=(iImg-1):(iImg+1)
for j=(jImg-1):(jImg+1)
thiS=thiS+sqrt(((I1(i,j)-aveX1)^2)*(abs(detX1(i,j)+detY1(i,j)))*(I1(i,j)^2));
end
end
S1(iImg,jImg)=firS*secS*thiS;
end
end
S1=imresize(S1,[s1(1,1),s1(1,2)]);
%计算第二幅图的显著性
S2=zeros(m,n);
aveX2=sum(I2(:))/(m*n);
detX2=zeros(m,n);
detY2=zeros(m,n);
M=3;N=3;
for i=1:m-1
for j=1:n-1
detX2(i,j)=abs(I2(i,j)-I2(i,j+1));
detY2(i,j)=abs(I2(i,j)-I2(i+1,j));
end
end
for iImg=2:(m-1)
for jImg=2:(n-1)
firS=(1/(M*N))^2;
secS=0;
for i=(iImg-1):(iImg+1)
for j=(jImg-1):(jImg+1)
secS=secS+I2(i,j)-aveX2;
end
end
thiS=0;
for i=(iImg-1):(iImg+1)
for j=(jImg-1):(jImg+1)
thiS=thiS+sqrt(((I2(i,j)-aveX2)^2)*(abs(detX2(i,j)+detY2(i,j)))*(I2(i,j)^2));
end
end
S2(iImg,jImg)=firS*secS*thiS;
end
end
S2=imresize(S2,[s2(1,1),s2(1,2)]);
%利用显著性进行融合，构造新的小波系数
facNew=zeros(m,n);
for i=1:s1(1,1)
for j=1:s1(1,2)
if S1(i,j) >= S2(i,j)
facNew(i,j)=facI1(i,j);
facNew((i+s1(1,1)),j)=facI1((i+s1(1,1)),j);
facNew(i,(j+s1(1,2)))=facI1(i,(j+s1(1,2)));
facNew((i+s1(1,1)),(j+s1(1,2)))=facI1((i+s1(1,1)),(j+s1(1,2)));
facNew(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((j*2-1):(j*2)))=facI1(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((j*2-1):(j*2)));
facNew(((i*2-1):(i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)))=facI1(((i*2-1):(i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)));
facNew(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)))=facI1(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)));
else
facNew(i,j)=facI2(i,j);
facNew((i+s1(1,1)),j)=facI2((i+s1(1,1)),j);
facNew(i,(j+s1(1,2)))=facI2(i,(j+s1(1,2)));
facNew((i+s1(1,1)),(j+s1(1,2)))=facI2((i+s1(1,1)),(j+s1(1,2)));
facNew(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((j*2-1):(j*2)))=facI2(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((j*2-1):(j*2)));
facNew(((i*2-1):(i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)))=facI2(((i*2-1):(i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)));
facNew(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)))=facI2(((s1(3,1)+i*2-1):(s1(3,1)+i*2)),((s1(3,2)+j*2-1):(s1(3,2)+j*2)));
end
end
end
%获取融合后的小波系数
caNew=facNew(1:s1(1,1),1:s1(1,2));
ch2New=facNew(1:s1(1,1),(s1(1,2)+1):(s1(1,2)*2));
cv2New=facNew((s1(1,1)+1):(s1(1,1)*2),1:s1(1,2));
cd2New=facNew((s1(1,1)+1):(s1(1,1)*2),(s1(1,2)+1):(s1(1,2)*2));
ch1New=facNew(1:s1(3,1),(s1(3,2)+1):(s1(3,2)*2));
cv1New=facNew((s1(3,1)+1):(s1(3,1)*2),1:s1(3,2));
cd1New=facNew((s1(3,1)+1):(s1(3,1)*2),(s1(3,2)+1):(s1(3,2)*2));
%重构小波系数，然后重构图像
for i=1:s1(1,1)
for j=1:s1(1,2)
caNewGen((j-1)*s1(1,1)+i)=caNew(i,j);
end
end
for i=1:s1(2,1)
for j=1:s1(2,2)
ch2NewGen((j-1)*s1(2,1)+i)=ch2New(i,j);
end
end
for i=1:s1(2,1)
for j=1:s1(2,2)
cv2NewGen((j-1)*s1(2,1)+i)=cv2New(i,j);
end
end
for i=1:s1(2,1)
for j=1:s1(2,2)
cd2NewGen((j-1)*s1(2,1)+i)=cd2New(i,j);
end
end
for i=1:s1(3,1)
for j=1:s1(3,2)
ch1NewGen((j-1)*s1(3,1)+i)=ch1New(i,j);
end
end
for i=1:s1(3,1)
for j=1:s1(3,2)
cv1NewGen((j-1)*s1(3,1)+i)=cv1New(i,j);
end
end
for i=1:s1(3,1)
for j=1:s1(3,2)
cd1NewGen((j-1)*s1(3,1)+i)=cd1New(i,j);
end
end
CNew=[caNewGen,ch2NewGen,cv2NewGen,cd2NewGen,ch1NewGen,cv1NewGen,cd1NewGen]; %排列处理后生成的系数
imgComRes=waverec2(CNew,s1,'haar'); %重构
subplot(1,3,3);
imshow(uint8(imgComRes))
