function Wavelet_Vector_Comprssion()
%X=imread('j:\數字圖像處理－活用Matlab程序\CH7\lena.bmp','bmp');
%[LL1,LH1,HL1,HH1]=dwtper2(X,'bior3.3');
%load lena
X=imread('E:\matlab1\LENA256.bmp','bmp');
%figure(100);imshow(uint8(X));title('original image');
%利用正交小波函數對X進行多尺度分解
%[C,S]=wavedec2(double(X),N,Wavelet_Name)是使用指定的小波基‘Wavelet_Name’對圖像X進行N層分解
Wavelet_Name='db1';
%Wavelet_Scale=1;
[c,s]=wavedec2(double(X),2,Wavelet_Name);
%appcoef2用來提取二維信號小波分解的近似分量
%   提取一維小波系數
ca1=appcoef2(c,s,Wavelet_Name,1);

%detcoef2用來提取二維信號小波分解的細節分量
%D=detcoef2(O,C,S,N)，O指定細節分量的類型，
%當O=‘h’時，重構水平細節分量，當O=‘v’時，重構垂直細節分量，當O=‘d’時，重構細節細節分量
cv1=detcoef2('v',c,s,1);
ch1=detcoef2('h',c,s,1);
cd1=detcoef2('d',c,s,1);


%   提取二維小波系數
ca2=appcoef2(c,s,Wavelet_Name,2);
cv2=detcoef2('v',c,s,2);
ch2=detcoef2('h',c,s,2);
cd2=detcoef2('d',c,s,2);
%figure;imshow(uint8(ca1))

%訓練LL1頻帶碼本，將碼本大小設為128，維度為4 %
th=0.05;
nc=64;nd=4;
[cbLL2]=LBG(ca2,nc,nd,th); %LBG.m在L3.3節中




%利用一層小波變換加上向量量化進行Lena圖像壓縮%
%輸入一幅Lena圖像%
X=imread('E:\matlab1\LENA256.bmp','bmp');
%load lena
[M,N]=size(X);
%使用Matlab中Wavelet工具箱的周期性小波變換函數dwtper2.m%
%并使用雙正交小波bior3.3進行一層小波變換%
%[LL1,LH1,HL1,HH1]=dwtper2(X,'bior3.3');
[C,S]=wavedec2(double(X),2,Wavelet_Name);
%[C,S]=wavedec2(double(X),Wavelet_Scale,Wavelet_Name);
LL1=appcoef2(C,S,Wavelet_Name,1);
LH1=detcoef2('v',C,S,1);
HL1=detcoef2('h',C,S,1);
HH1=detcoef2('d',C,S,1);

LL2=appcoef2(C,S,Wavelet_Name,2);
LH2=detcoef2('v',C,S,2);
HL2=detcoef2('h',C,S,2);
HH2=detcoef2('d',C,S,2);

wt_X(1:M/4,1:N/4)=LL2;
wt_X(1:M/4,N/4+1:N/2)=LH2;
wt_X(M/4+1:M/2,1:N/4)=HL2;
wt_X(M/4+1:M/2,N/4+1:N/2)=HH2;
%分別針對不同的頻帶進行向量量化，其中VQ.m置于L3.4節中%
[R_LL2,bitLL2]=VQ(LL2,cbLL2);
wt_Y(1:M/4,1:N/4)=R_LL2;
u=mean2(LH2);
sum1=0;
[m,n]=size(LH2)
for i=1:1:m
    for j=1:1:n
   if LH2(i,j)<u
        LH2(i,j)=0;
       sum1=sum1+1; 
   end
    end
end
sum1
u=mean2(HL2);
sum1=0;
[m,n]=size(HL2)
for i=1:1:m
    for j=1:1:n
   if HL2(i,j)<u
        HL2(i,j)=0;
       sum1=sum1+1; 
   end
    end
end
sum1
u=mean2(HH2);
sum1=0;
[m,n]=size(HH2)
for i=1:1:m
    for j=1:1:n
   if HH2(i,j)<u
        HH2(i,j)=0;
       sum1=sum1+1; 
   end
    end
end
sum1
u=mean2(HL1);
sum1=0;
[m,n]=size(HL1)
for i=1:1:m
    for j=1:1:n
   if HL1(i,j)<u
        HL1(i,j)=0;
       sum1=sum1+1; 
   end
    end
end
sum1
u=mean2(LH1);
sum1=0;
[m,n]=size(LH1)
for i=1:1:m
    for j=1:1:n
   if LH1(i,j)<u
       LH1(i,j)=0;
       sum1=sum1+1; 
   end
    end
end
sum1
u=mean2(HH1);
sum1=0;
[m,n]=size(HH1)
for i=1:1:m
    for j=1:1:n
   if HH1(i,j)<u
        HH1(i,j)=0;
       sum1=sum1+1; 
   end
    end
end
sum1



%計算圖像經過向量量化后的壓縮率CR及失真PSNR%
totalbit=bitLL2        %花費的總位數  =4416
%Y=idwtper2(R_LL1,R_LH1,R_HL1,R_HH1,'bior3.3'); %反小波變換

C2= [R_LL2 HL2 LH2 HH2];
s2=s;
%a2=wrcoef2('a',R_LL2,s2,'db1',2);
Y1=idwt2(R_LL2, HL2 ,LH2, HH2,'db1',size(LL1));
%figure(66);imshow(uint8(Y2));title('compressed image');

Y2=idwt2(Y1,HL1 ,LH1, HH1,'db1',size(X));
%size(a2)
figure(100);imshow(uint8(X));title('original image');
figure(200);imshow(uint8(Y2));title('compressed image');

%Y2=waverec2(C2,s2,'db1',1);
%size(a2)size(Y)=[128 128]
%figure(300);imshow(uint8(Y2));
X=double(X);
Y2=double(Y2);

MSE=(sum(sum((X-Y2).^2)))/(M*N)
PSNR=20*log10(255/sqrt(MSE))
CR=M*N*8/totalbit



%figure;imshow(uint8(Y));

%figure;imshow(uint8(c1));
%figure;imshow(uint8(c2));
%figure;imshow(uint8(ca2));