%%%%%%%%%讀入圖像信號
I=imread('girl.bmp');
figure;
imshow(I);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%第一次圖像壓縮
[cA1,cH1,cV1,cD1] = dwt2(I,'db1');
%[cA1,cH1,cV1,cD1] = dwt2(I,'harr');采用不同的小波函數
%[cA1,cH1,cV1,cD1] = dwt2(I,'bior3.7');
cA2=uint8(cA1);cH2=uint8(cH1);cV2=uint8(cV1);cD2=uint8(cD1);
C1=[cA2 cH2;cV2 cD2];
figure;
subplot(2,2,1);
imshow(cA2);
title('低頻近似cA1');
subplot(2,2,2);
imshow(cH2);
title('水平細節cH1');
subplot(2,2,3);
imshow(cV2);
title('垂直細節cV1');
subplot(2,2,4);
imshow(cD2);
title('對角線細節cD1');
%高頻置零
[M,N]=size(cA1);
%cH1=zeros(M,N);cV1=zeros(M,N);cD1=zeros(M,N);
% [cH1,k1]=thre(cH1,0);
% [cV1,k2]=thre(cV1,0);
% [cD1,k3]=thre(cD1,0);
[cH1,cV1,cD1,cpr]=thre2(cH1,cV1,cD1);
Isyn=idwt2(cA1,cH1,cV1,cD1,'db1');
%Isyn=idwt2(cA1,cH1,cV1,cD1,'harr');
%Isyn=idwt2(cA1,cH1,cV1,cD1,'bior3.7');

%imwrite(Isyn,'girl2.bmp');

% 第一次壓縮后的信噪比
I1=double(I);
r1=snr(I1,Isyn);%必需轉換成需要的double格式，否則無法使用log10函數

Isyn=uint8(Isyn);
figure;
imshow(Isyn);
title('第一次壓縮后的圖像');



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%第2次壓縮
[cA3,cH3,cV3,cD3] = dwt2(cA1,'db1');
%[cA3,cH3,cV3,cD3] = dwt2(cA1,'harr');
%[cA3,cH3,cV3,cD3] = dwt2(cA1,'bior3.7');

cA4=uint8(cA3);cH4=uint8(cH3);cV4=uint8(cV3);cD4=uint8(cD3);
C2=[cA4 cH4;cV4 cD4];
figure;
imshow(cA4);
title('第二次小波變換后的低頻近似cA3');
%高頻置零
[M1,N1]=size(cA3);
%cH3=zeros(M,N);cV3=zeros(M,N);cD3=zeros(M,N);
[cH3,k1]=thre(cH3,0);
[cV3,k2]=thre(cV3,0);
[cD3,k3]=thre(cD3,0);
Isyn2=idwt2(cA3,cH3,cV3,cD3,'db1');
%imwrite(Isyn2,'girl2.bmp');
Isyn2=uint8(Isyn2);
figure;
imshow(Isyn2);
title('第2次壓縮后的低頻近似部分');


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%從圖像上顯示壓縮過程
C1=[C2 cH2;cV2 cD2];
figure;
imshow(C1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%計算信噪比

ca=idwt2(cA3,cH3,cV3,cD3,'db1');
Isyn3=idwt2(ca,cH1,cV1,cD1,'db1');
r2=snr(I1,Isyn3);% 第二次壓縮后的信噪比

Isyn3=uint8(Isyn3);
imshow(Isyn3);
title('恢復后的圖像');