close all;
clear all;

data=imread('cameraman.tif');
figure(1);
colormap(gray);
subplot(2,2,1);
imagesc(data);
axis image;
title('original image');
a=-1.586134342;
b=-0.052980118;
r=0.882911075;
e=0.443506852;
p=1.149604398;
% di=zeros(256,128);
% si=zeros(256,128);
for i=1:256
    for j=1:128
        d(1,j)=data(i,j*2-1);
        s(1,j)=data(i,j*2);
    end
    d1(1,1)=d(1,1)+fix(a*(s(1,1)+s(1,2))+0.5);
    d1=d;
    s1=s;
    for j=2:127
        d1(1,j)=d(1,j)+fix(a*(s(1,j)+s(1,j+1))+0.5);
    end
    for j=2:127
        s1(1,j)=s(1,j)+fix(b*(d1(1,j-1)+d1(1,j))+0.5);
    end
    d1(1,128)=d(1,128)+fix(a*(s(1,128)+s(1,1))+0.5);
    s1(1,1)=s(1,1)+fix(b*(d1(1,128)+d1(1,1))+0.5);
    s1(1,128)=s(1,128)+fix(b*(d1(1,127)+d1(1,128))+0.5);
    d2(1,1)=d1(1,1)+fix(r*(s1(1,1)+s1(1,2))+0.5);
    d2=d1;
    s2=s1;
    for j=2:127
        d2(1,j)=d1(1,j)+fix(r*(s1(1,j)+s1(1,j+1))+0.5);
    end   
    for j=2:127
        s2(1,j)=s1(1,j)+fix(e*(d2(1,j-1)+d2(1,j))+0.5);
    end
    d2(1,128)=d1(1,128)+fix(r*(s1(1,128)+s1(1,1))+0.5);
    s2(1,1)=s1(1,1)+fix(e*(d2(1,128)+d2(1,1))+0.5);
    s2(1,128)=s1(1,128)+fix(e*(d2(1,127)+d2(1,128))+0.5);
    for j=1:128
        di(i,j)=d2(1,j)*p;
        si(i,j)=s2(1,j)/p;
    end
end
D=zeros(128,128);
S=zeros(128,128);
for i=1:128
    for j=1:128
        d(1,j)=di(i*2-1,j);
        s(1,j)=data(i*2,j);
    end
    d1(1,1)=d(1,1)+fix(a*(s(1,1)+s(1,2))+0.5);
    d1=d;
    s1=s;
    for j=2:127
        d1(1,j)=d(1,j)+fix(a*(s(1,j)+s(1,j+1))+0.5);
    end
    for j=2:127
        s1(1,j)=s(1,j)+fix(b*(d1(1,j-1)+d1(1,j))+0.5);
    end
    d1(1,128)=d(1,128)+fix(a*(s(1,128)+s(1,1))+0.5);
    s1(1,1)=s(1,1)+fix(b*(d1(1,128)+d1(1,1))+0.5);
    s1(1,128)=s(1,128)+fix(b*(d1(1,127)+d1(1,128))+0.5);
    d2(1,1)=d1(1,1)+fix(r*(s1(1,1)+s1(1,2))+0.5);
    d2=d1;
    s2=s1;
    for j=2:127
        d2(1,j)=d1(1,j)+fix(r*(s1(1,j)+s1(1,j+1))+0.5);
    end   
    for j=2:127
        s2(1,j)=s1(1,j)+fix(e*(d2(1,j-1)+d2(1,j))+0.5);
    end
    d2(1,128)=d1(1,128)+fix(r*(s1(1,128)+s1(1,1))+0.5);
    s2(1,1)=s1(1,1)+fix(e*(d2(1,128)+d2(1,1))+0.5);
    s2(1,128)=s1(1,128)+fix(e*(d2(1,127)+d2(1,128))+0.5);
    for j=1:128
        D(i,j)=d2(1,j)*p;
        S(i,j)=s2(1,j)/p;
    end
end
subplot(2,2,2);
imagesc(D);
axis equal;
subplot(2