多聚焦图像融合源代码演示教学

多聚焦图像融合源代

码

针对经典的最大系数法不准确和方差法计算量大的问题，本文给出了一种混合多级式多聚焦图像融合方法。对于三层小波分解的多聚焦图像融合，每幅图像被分解为三层十个频带。对这十个频带本文分别采用三种方法进行融合。对于低频系数，本文仍然采用求平均法；对于高频系数本文采用方差法和最大系数法进行融合。它们的计算量比最大系数法大一些，但是融合结果更接近于原始清晰图像，而相比于方差法，它们的计算量小的多，但是融合质量稍差一些，应用者可以根据不同的需要进行选择。

本文还给出了一种基于Canny算子边缘检测的小波变换多聚焦图像融合方法。首先对图像进行三层小波分解，然后用Canny算子进行边缘检测，得到各层分辨率下的边缘图像；对相应分辨率的高频小波系数根据其是否为图像的边缘点采用最大系数法或方差法分别进行融合。仿真实验证明该方法效果良好，计算量可以灵活调节。

关键词:小波变换；多尺度几何分析；多聚焦图像融合；边缘检测

主要程序：

clear all;

close all;

leo1=imread('a1.bmp');%读入图片

leo2=imread('a2.bmp')

T=0.4;k1=0.5;k2=0.5;w='db4';m='edge'; tic;

outdoor1=leo1;

outdoor2=leo2;

%三层小波分解

[ca11,chd11,cvd11,cdd11]=dwt2(outdoo r1,w);

[ca12,chd12,cvd12,cdd12]=dwt2(ca11,w );

[ca13,chd13,cvd13,cdd13]=dwt2(ca12,w );

[ca21,chd21,cvd21,cdd21]=dwt2(outdoo r2,w);

[ca22,chd22,cvd22,cdd22]=dwt2(ca21,w );

[ca23,chd23,cvd23,cdd23]=dwt2(ca22,w );

%求边缘图像

e11=edge(ca11,'canny',T);

e12=edge(ca12,'canny',T);

e13=edge(ca13,'canny',T);

e21=edge(ca21,'canny',T);

e22=edge(ca22,'canny',T);

e23=edge(ca23,'canny',T);

%矩阵融合

chd3=matfusion(chd13,chd23,e13,e23); cvd3=matfusion(cvd13,cvd23,e13,e23); cdd3=matfusion(cdd13,cdd23,e13,e23); chd2=matfusion(chd12,chd22,e12,e22); cvd2=matfusion(cvd12,cvd22,e12,e22); cdd2=matfusion(cdd12,cdd22,e12,e22); chd1=matfusion(chd11,chd21,e11,e21); cvd1=matfusion(cvd11,cvd21,e11,e21); cdd1=matfusion(cdd11,cdd21,e11,e21); ca3=k1*ca13+k2*ca23;

%反小波变换

L2=size(chd2);L1=size(chd1);

ca2=idwt2(ca3,chd3,cvd3,cdd3,w);

ca1=idwt2(ca2(1:L2(1),1:L2(2)),chd2,cv d2,cdd2,w);

I=idwt2(ca1(1:L1(1),1:L1(2)),chd1,cvd1 ,cdd1,w);

toc;

%显示图像

subplot(2,2,1);imshow(leo1); colormap(gray);

title('input1');

axis square

subplot(2,2,2);imshow(leo2); colormap(gray);

title('input2');

axis square

subplot(2,2,3);imshow(uint8(I)); colormap(gray);

title('最大系数法与方差法融合'); axis square;

%标准差

f_stdvar=sta_variance1(I)

%清晰度，又称为平均梯度

f_avggrad=avg_gradient(I)

f=PSNR(I)%%%%% 信噪比

%其值越大，说明融合图像和标准参考图像%的差异越小，融合效果越好se=RMSE(I)% %RMSE(均方根误差,也称%为标准偏差，融合图像与标准图像差异%%程度，差异越小说明融合图像与标准图%像越接近%%%%%%%%%%%%

%I3 = imread('clockB.jpg');

[sI,cI,bI] = size(uint8(I));

I3 = double(uint8(I))+1;

s(256) = 0;

for i = 1:sI

for j = 1:cI

va = I3(i,j);

s(va) = s(va)+1;

end

end

p = s./(sI*cI);

en = 0.;

for i = 1:256

if p(i)~=0

en = en+p(i)*log2(p(i));

end

end

entropy = -en