function mydwt2 %主程序用于測試圖像,
I=imread('hao','bmp');
f=I(1:400,1:400,1);
f=double(f);
g=mydwt(3,f);%分解后圖,3為層數(shù),可以任意大小的整數(shù)
;%任意層數(shù)重構(gòu)，這里用3層
r=imydwt(3,g);
t=f-r
f=uint8(f);

imshow(f);
figure%重構(gòu)后圖
g=uint8(g);
imshow(g);
figure
r=uint8(r);
imshow(r);

function B=mydwt(k,I)    %分解子程序
%分解層數(shù)為k，采用純二維小波分解，二維提升.可原址運算.
B=I;
[N,M]=size(I);
while k~=0
    [L,J]=size(zeros(N/2,M/2));
    T=ones(2,2);
for t=1:2:N-1
    for tt=1:2:M-1
        T=I(t:t+1,tt:tt+1);
       T=1/4*T;
       e=sum(T(:));
       right=e-I(t,tt+1);
       down=e-I(t+1,tt);
       diag=e-I(t+1,tt+1);
       B((t+1)/2,(tt+1)/2)=e;
       B((t+1)/2,J+(tt+1)/2)=right;
       B(L+(t+1)/2,(tt+1)/2)=down;
       B(L+(t+1)/2,(tt+1)/2+J)=diag;
end
end
I=B;
k=k-1;
N=N/2;
M=M/2;
end
%重構(gòu)為分解的反向運算
function B=imydwt(k,I) %重構(gòu)子程序
B=I;
[N,M]=size(I);
N=N/2^(k-1);
M=M/2^(k-1);
while k~=0
    [L,J]=size(zeros(N/2,M/2));
    for i=1:2:N-1
    for j=1:2:M-1
        B(i,j+1)=I((i+1)/2,(j+1)/2)-I((i+1)/2,(j+1)/2+J);
        B(i+1,j)=I((i+1)/2,(j+1)/2)-I(L+(i+1)/2,(j+1)/2);
        B(i+1,j+1)=I((i+1)/2,(j+1)/2)-I(L+(i+1)/2,J+(j+1)/2);
        B(i,j)=4*I((i+1)/2,(j+1)/2)-B(i,j+1)-B(i+1,j)-B(i+1,j+1);
    end
end
I(1:N,1:M)=B(1:N,1:M);
 k=k-1;
 N=N*2;
 M=M*2;
end