图像的Mallat算法分解(Matlab代码)

Mallat 算法的分析与综合框架参考书上的资料很多，这里就不多说了。

下面是我写的关于图像的程序，分别是：一维分解，二维分解；一维合成，二维合成。最后是测试主程序。

谢谢参考，错了请反馈一下！

%内部子函数，对一行(row)矢量进行一次小波变换，利用fft实现
function y=mdec1(x,h,g)
%x 行数组
%h为低通滤波器
%g为高通滤波器
%输出： y 进行一级小波分解后的系数
lenx=size(x,2);                 %求行的长度
lenh=size(h,2);                 %求低通滤波器的长度
rh=h(end:-1:1);                 %rh是h的逆序  h0(n)=h(-n)
rrh=[zeros(1,(lenx-lenh)),rh];  %rrg在rg前插入(lenx-lenh))个零
rrh=circshift(rrh',1)';         %循环位移

rg=g(end:-1:1);                 %rg为高通滤波器的逆序g0(n)=g(-n)
rrg=[zeros(1,(lenx-lenh)),rg];  %rrg在rg前插入(lenx-lenh))个零
rrg=circshift(rrg',1)';         %循环位移

r1=dyaddown(ifft(fft(x).*fft(rrh,lenx)),1); %use para 1，下2奇采样
r2=dyaddown(ifft(fft(x).*fft(rrg,lenx)),1); %use para 1，下2奇采样
y=[r1,r2];                     %相加得到小波系数

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function Y=mallatdec2(X,wname,level)
%输入:X????? 载入的二维图像像数值;
%???? level? 小波分解次(级)数设定值(如果设定值超过最高可分解次数,按最高分解次数分解)
%????? wname? 小波名字wavelet name
%输出：Y???? 多极小波分解后的小波系数矩阵
%wfilters：小波滤波器

[h,g]=wfilters(wname,'d');    %h,g分别为分解的低通和高通滤波器
X=double(X);
t=1;
hh=size(X,2);
Y=zeros(1,hh);
while t<=level   %循环
    %先进行行小波变换
    for row=1:hh
        Y(row,1:hh)=mdec1(X(row,1:hh),h,g) ;
    end
    %再进行列小波变换
    for col=1:hh
        temp=mdec1( Y(1:hh,col)',h,g); %Y(1:hh,col)为矩阵中的一列，转置后为行向量，得到的结果temp也是行向量
        Y(1:hh,col)=temp';             %将行向量转置为列向量存储在Y中
    end
    t=t+1;                             %循环次数加1，大于level，循环终止
    hh=hh/2;                           %分解的矩阵的大小为原来的四分之一
    X=Y;                               %将Y值赋给X等待下一次分解                     
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%内部子函数，对一行小波系数进行重构
function y=mrec1(x,h,g)
%输入：x 行数组
%???? h为低通滤波器
%???? g为高通滤波器
%输出： y 进行一级小波重构后值
lenx=size(x,2);

r3=dyadup(x(1,1:lenx*0.5),0);        %内插零use para 0，上采样
r4=dyadup(x(1,(lenx*0.5+1):lenx),0); %内插零use para 0，上采样
y=ifft(fft(r3,lenx).*fft(h,lenx))+ ifft(fft(r4,lenx).*fft(g,lenx));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%二维小波重构函数
function Y=mallatrec2(X,wname,level)
%输入:X????? 载入的小波系数矩阵;
%?? ??level? 小波分解次(级)数设定值(如果设定值超过最高可分解次数,按最高分解次数分解)
%????? wname? 小波名字wavelet name
%输出：Y???? 重构图像矩阵

[h,g]=wfilters(wname,'d');%h为重构的低通滤波器；g为重构高通滤波器

hz=size(X,2);
h1=hz/(2^(level-1));
while h1<=hz
    % 对列变换
    for col=1:h1
        temp=mrec1(X(1:h1,col)',h,g)';
        X(1:h1,col)=temp;
    end
    %再对行变换
    for row=1:h1
        temp=mrec1(X(row,1:h1),h,g);
        X(row,1:h1)=temp;
    end
    h1=h1*2;
end
Y=X;

 

%%%%%%%%%%%%%%%%%%%%%%%%%%

%测试函数（主函数）
clc;
clear;
close all;

X=imread('lena.bmp');;%路径
figure(1);imshow(X);

X=double(X);
level=2;
A = mallatdec2(X,'db4',level);

figure(2);
imshow(uint8(A));title('多尺度分解图像');

Y= mallatrec2(A,'db4',level);

figure(3);
imshow(uint8(Y));title('重构图像');