%=====================================================
% Demo program for chain-coding
% Image processing course, Ghent University, TELIN
%
% Aleksandra Pizurica, Ghent University, TELIN, 2003.
%======================================================

clear all;
x=imread('Contour.tif');
C=bwperim(x,8);
%preprocessing
C(26,56)=0;
subplot(2,1,1),imagesc(C)
title('Input Contour');
grid on;



CONTOUR(1,:)=input('Coordinates of the first point ([27 57])');
CONTOUR(2,:)=input('Coordinates of the second point ([27 58])');
Dir=input('Direction between the first two points: (0) ');

Count=1;
DIR(1)=Dir
istart=CONTOUR(1,1);
jstart=CONTOUR(1,2);

i=CONTOUR(2,1);
j=CONTOUR(2,2);
while(i~=istart | j~=jstart)
    
ForbiddenDir=mod(Dir+2,4);

if(ForbiddenDir~=0 & C(i,j+1)==1)
        Dir=0;
        j=j+1;
        %break;
elseif(ForbiddenDir~=1 & C(i-1,j)==1)
        Dir=1;
        i=i-1;
        %break;
elseif(ForbiddenDir~=2 & C(i,j-1)==1)
        Dir=2;
        j=j-1;
        %break;
elseif(ForbiddenDir~=3 & C(i+1,j)==1)
        Dir=3;
        i=i+1;
        %break;
end;
Count=Count+1
CONTOUR(Count+1,1)=i;
CONTOUR(Count+1,2)=j;
DIR(Count)=Dir;
end



%Plot the calculated contour points
[M,N]=size(x);
A=zeros(M,N);
for i=1:(Count+1)
    A(CONTOUR(i,1),CONTOUR(i,2))=1;
end
subplot(2,1,2)
imagesc(A);
title('Stored contour points');


% Manipulate the Fourier descriptors

z=CONTOUR(:,1)+sqrt(-1)*CONTOUR(:,2);
s=fft(z);
z1=z/5;
s1=fft(z1);
C1=ifft(s1);
u=round(real(C1));
v=round(imag(C1));
T=zeros(M,N);
for i=1:length(u)
  T(u(i),v(i))=1;
end
figure
subplot(2,2,1),plot(abs(s(10:650)))
axis([10 650 0 40])
title('(a) Fourier magnitude spectrum');
subplot(2,2,2),plot(abs(s1(10:650)))
title('(b) Magnitude spectrum divided by 5');
axis([10 650 0 40])
subplot(2,2,3),imagesc(A)
title('Reconstruction from (a)');
subplot(2,2,4),imagesc(T)
title('Reconstruction from (b)');

