小波变换mallat

1、-. z.实验目的：通过编程实现离散快速小波变换Mallat算法，从而加深理解二维小波变换的分解与合成，同时，提高编程能力和matlab的应用，为以后的学习打下根底。实验原理：1、Mallat快速算法本实验使用离散快速小波变换快速算法Mallat算法，算法原理如下 12重构算法： 3对于1、2等效于经过冲击响应为和的数字滤波器，然后再分别进展二抽取，Mallat分解算法的滤波器表示形式如下列图用滤波器表示如下列图2、 分别表示原始图像和重建后的图像，。3、边界延拓方法有零延拓、周期延拓、对称周期延拓、常数连续延拓等，本实验采用以上四种方法进展原图像的1/8延拓，并进展重构，各种延拓方法所对应的

2、函数为yan0(*)、yancir*、yan(*)、yanc(*),在主程序中，需要*种延拓，便调用*种函数。实验编程思路：为使程序易于理解，在不考虑算法复杂度的情况下，分解程序采用简洁的循环计算出下一级的分解系数，程序采用的编程思想如下以上矩阵等式左面是进展二抽样的结果，是分解的低频局部。同理，对于分解的高频局部有如下矩阵形式：分解程序： len*=size(*,2);%*为一维向量lenh=size(h,2);h=h,zeros(1,(len*-lenh);g=g,zeros(1,(len*-lenh);r1(1)=sum(h.*);r2(1)=sum(g.*);for k=1:1:(le

3、n*/2-1) %循环求出下一级低频和高频分量h=h(end-1:end),h(1:(end-2);r1(k+1)=sum(h.*);g=g(end-1:end),g(1:1:(end-2);r2(k+1)=sum(g.*);endy=r1,r2;对于重构算法，其等效形式为上式等号右边局部实质上是对变量的数字卷积运算,程序采用频域相乘代替卷积，重建程序为 y=ifft(fft(c3,len*).*fft(h,len*)+ ifft(fft(d3,len*).*fft(g,len*);实验结果及分析：1、多尺度分解与重构图像二维小波变换采用小波采用db3，其峰值信噪比PNSR=230.13db，

4、并对三级分解图像进展归一化，求出0的个数为37626，其所占的百分比为57.41%。2、延拓重建图像延拓方法周期延拓对称周期延拓零延拓常数连续延拓PNSR230.18db230.48db229.80db230.13db从PNSR结果可知，在各种延拓中，对称周期延拓的重建图像结果最好，相比之下零延拓图像效果不如其他方法延拓。3、不同小波下重构图像的性质 用不同小波进展图像重构，所得的重构图像能量分布如下用各种小波进展重构后的图像的均值方差如下表。小波db1db2db3db4均值124.0309124.0509124.0509124.0509方差2272227222722272附录：1、主函数程序

5、clc;clear;*=imread(LENA.bmp);%路径*=double(*);% S=yancir(*);A=mallatdec2(*,db3, 3);image(abs(A);colormap(gray(255);title(3级多尺度分解图像);Y=mallatrec2(A,db3,3);Y=real(Y);figure(2);subplot(1,2,1);image(*);colormap(gray(255);title(原始图像);subplot(1,2,2);image(Y);colormap(gray(255);title(重建图像);zerosn=numberzeros

6、(A);% Y=Y(33:288,33:288); %当调用延拓图像时，从延拓的重建图像进展截取csize=size(*);sr=csize(1);sc=csize(2);mse=sum(sum( (Y-*).2,1)/(sr*sc);psnr=10*log(255*255/mse)/log(10)2、分解程序function Y=mallatdec2(*,wname,level)%输入:* 载入的二维图像像数值;% level 小波分解次(级)数设定值(如果设定值超过最高可分解次数,按最高分解次数分% wname 小波名字wavelet name%输出：Y 多极小波分解后的小波系数矩阵h,g

7、=wfilters(wname,d); %h,g 分别为低通和高通滤波器*=double(*);t=1;hh,ll=size(*);while t=level%先进展行小波变换for row=1:hh Y(row,1:ll)=mdec1(*(row,1:ll),h,g) ;end%再进展列小波变换for col=1:ll temp=mdec1( Y(1:hh,col),h,g); Y(1:hh,col)=temp;end t=t+1; hh=hh/2; ll=ll/2; *=Y;endfunction y=mdec1(*,h,g)%输入：* 行数组% h 为低通滤波器% g 为高通滤波器%输出

8、： y 进展一级小波分解后的系数len*=size(*,2);lenh=size(h,2);h=h,zeros(1,(len*-lenh);g=g,zeros(1,(len*-lenh);r1(1)=sum(h.*);r2(1)=sum(g.*);for k=1:1:(len*/2-1) h=h(end-1:end),h(1:(end-2); r1(k+1)=sum(h.*); g=g(end-1:end),g(1:1:(end-2); r2(k+1)=sum(g.*);endy=r1,r2;3、重建程序function Y=mallatrec2(*,wname,level)%输入:* 载入的

9、小波系数矩阵;% level 小波分解次(级)数设定值(如果设定值超过最高可分解次数,按最高分解次数% wname 小波名字wavelet name%输出：Y 重构图像矩阵h,g=wfilters(wname,d); %h,g 分别为重构低通滤波器和重构高通滤波器hz=size(*,2);h1=hz/(2(level-1);while h1=hz% 对列变换for col=1:h1 temp=mrec1(*(1:h1,col),h,g); *(1:h1,col)=temp;end%再对行变换for row=1:h1 temp=mrec1(*(row,1:h1),h,g); *(row,1:h1

10、)=temp;endh1=h1*2;endY=*;function y=mrec1(*,h,g)%输入：* 行数组% h 为低通滤波器% g 为高通滤波器%输出： y 进展一级小波重构后值len*=size(*,2);lenh=size(h,2);h=h,zeros(1,(len*-lenh);g=g,zeros(1,(len*-lenh);c3=dyadup(*(1,1:len*0.5),0,len*); %内插零use para 0d3=dyadup(*(1,(len*0.5+1):len*),0,len*); %use para 0y=ifft(fft(c3,len*).*fft(h,l

11、en*)+ ifft(fft(d3,len*).*fft(g,len*);4、延拓程序%零延拓程序function T=yan0(*) r0,c0=size(*);for s=1:r0 a(s,:)=zeros(1,c0/8),*(s,:),zeros(1,c0/8);endc0=c0+c0/4;for t=1:c0 T(:,t)=zeros(r0/8,1),a(:,t),zeros(r0/8,1);End%连续常数延拓程序function T=yanc(*) r0,c0=size(*);for s=1:r0 r1(1:c0/8)=*(s,1); r2(1:c0/8)=*(s,c0); a(s,:)=r1,*(s,:),r2;endc0=c0+c0/4;for t=1:c0 c1(1:r0/8)=a(1,t); c2(1:r0/8)=a(r0,t); T(:,t)=c1,a(:,t),c2;end%对称周期延拓程序function T=yancir(*)r0,c0=size(*);for s=1:r0 a(s,:)=*(s,end-(c0/8-1):1:end),*(s,:),*(s,1:1:c0/8);endc0=c0+c0/4;for t=1:c0 T(:,t)=a(end-(r0/8-1):1:end,t),a(:,t),a(1:1:r0/8,t)