clear;

InitialImage=imread('lena.bmp');
%
InitialImage=double(InitialImage)/255;
%  
figure(1);
imshow(InitialImage);

T=hadamard(8);
DCTCoe=blkproc(InitialImage,[8 8],'P1*x*P2',T,T);
%
CoeVar=im2col(DCTCoe,[8 8],'distinct');
Coe=CoeVar;
%
[Y,Ind]=sort(CoeVar);  
%
[m,n]=size(CoeVar);
Snum=64-16;
for i=1:n
   Coe(Ind(1:Snum),i)=0;
end
%
B2=col2im(Coe,[8 8],[256 256],'distinct');
%
I2=blkproc(B2,[8 8],'P1*x*P2',T, T);
%
figure(4);
imagesc(I2);

% 計(jì)算歸一化圖像的均方誤差
error=InitialImage.^2-I2.^2;
MSE=sum(error(:))/prod(size(I2))