%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Author: Jingting Zeng (Nikki)
% Homework 5: Find the Gear
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function Newmain

% read the image
a=imread('objects.JPG');
a=im2double(a);

% convert the image format to prepare tracing_boundary
b=rgb2gray(a);
s = size(b);
for i=1:s(1)
    b(i,:) = 1-b(i,:);
end
b(find(b>0.5))=1;
b(find(b<0.5))=0;

% trace boundary and extract the boundary coordinate
c = imdilate(b,strel('line',2,0));
c = imerode(c,strel('disk',1));
[test, num] =  bwlabel(c,8);
imshow(test);
s = size(test);
flag=false;
figure;
k=1;
for i=1:num
    for row = 1:s(1)
        for col=1:s(2)
            if test(row,col)==i
                flag=true;
                break;
            end
        end

        if flag==true
            temp = bwtraceboundary(test, [row, col], 'N', 8);
            % exclude the noise
            if size(temp,1)>60
                contour{k} = temp;
                subplot(4,4,k);
                plot(contour{k}(:,2),-contour{k}(:,1),'g','LineWidth',2);
                title(strcat('Object ',num2str(k)));
                k=k+1;
            end

            flag=false;
            break;
        end
    end
end

figure;
s = size(contour);
for i=1:s(2)
    xcenter = [];
    ycenter = [];
    leng = [];
    fourier = [];
    
    % calculate the length between central point and boundary
     pointsize = size(contour{1,i},1);
%     xpoint = max(contour{1,i}(:,1))/2;
%     ypoint = max(contour{1,i}(:,2))/2;
%     xcenter(1:pointsize) = xpoint;
%     ycenter(1:pointsize) = ypoint;
%     leng = sqrt((contour{1,i}(:,1)-xcenter').*(contour{1,i}(:,1)-xcenter') + (contour{1,i}(:,2)-ycenter').*(contour{1,i}(:,2)-ycenter'));

    p=contour{i};
    center=sum(p)/length(p);
    p(:,1)=p(:,1)-center(1);
    p(:,2)=p(:,2)-center(2);
    leng=sqrt(sum(p.^2,2));
    
    
    % use FFT to calculate 
    fourier=abs(fft(leng,pointsize));
    cdf=fourier(2: floor(pointsize/2))./ fourier(1); %start from frequency 1
    
    % find the maximum number between 5 and pointsize/2 
    result(i) = find(cdf(:)==max(cdf(5:pointsize/2-1)));
   
    subplot(4,4,i);
    stairs(cdf);
    title(strcat('Object ',num2str(i)));
    
    % show the similar graphs as in "Fourier Shape Descriptors"
        cdf=fourier(2: 50)./ fourier(1); %show frequency 1-49
        figure;
        peri=perimeter(contour{1,i});
        plotcentroiddistance(leng,peri);
        title(strcat('Object ',num2str(i)));
        figure;
        stairs(cdf);
        title(strcat('Object ',num2str(i)));
    
end

% if the maximum number between 2 and pointsize/2 is larger than 10, it
% means that this object is a gear
gear_No = find(result>15)
wheel_No = result(gear_No)