x = imread('flower.jpg');
%image is 700 by 940
x=x(1:640, 201:840,:); %crop the image so that its dimentsion is a multiple of 8
x=double(x);
r=x(:,:,1); 
g=x(:,:,2);
b=x(:,:,3);
xgray=0.2126*r+0.7152*g+0.0722*b;
xgray=uint8(xgray);
%figure(1);imagesc(xgray);colormap(gray);
s=0;
image= zeros(640); % 640 by 640 matrix
for i = 1:8:640
    for j = 1:8:640
        X=xgray(i:i+7, j:j+7); %extract an 8 by 8 pixel block
        Xd=double(X);
        Xc=Xd-128;
        C=dct(Xc);
        Y=C*Xc*C';
        p=1; % loss parameter
        Q=p*8./hilb(8); %linear quantization matrix
        %Recommended Matrix for JPEG
%         Q=p*[16 11 10 16 24 40 51 61;
%             12 12 14 19 26 58 60 55 ;
%             14 13 16 24 40 57 69 56 ;
%             14 17 22 29 51 87 80 62 ;
%             18 22 37 56 68 109 103 77; 
%             24 35 55 64 81 104 113 92;
%             49 64 78 87 103 121 120 101;
%             72 92 95 98 112 100 103 99; ]
        Yq=round(Y./Q);
            
        %to recover the image back
        Ydq=Yq.*Q; 
        Xdq=C'*Ydq*C;
        Xe=Xdq+128;
        s=s+(sum(sum((Xd-Xe).^2)));
              
        Xf=uint8(Xe);
        image(i:i+7, j:j+7)=Xf;
    end
end
rms=s/(640*640)
figure(1);imagesc(image);colormap(gray);