function X = find_edges(I, sigma, Th, Tl)

[F D Fx Fy] = filtered_gradient(I, sigma);

%D = D + pi/2;
D(isnan(D)) = 0;
D = abs(D) ./ (pi/2);
N = zeros(size(D));
N(D >= .25) = 1;
N(D >= .50) = 2;
N(D >= .75) = 3;

M = zeros(size(F));

for i=2:size(N,1)-1;
    for j=2:size(N,2)-1;
           
        if N(i,j) == 3;
            left = F(i+1,j);
            right = F(i-1,j);
            
        elseif N(i,j) == 2;
            left = F(i+1, j+1);
            right = F(i-1,j-1);
            
        elseif N(i,j) == 0;
            left = F(i, j-1);
            right = F(i, j+1);
            
        elseif N(i, j) == 1;
            left = F(i+1,j-1);
            right = F(i-1,j+1);
        end;
        
        if( F(i,j) < left || F(i,j) < right )
            M(i,j) = 0;
        else
            if(F(i,j) >= Th)
                M(i,j) = 1;
            elseif(F(i,j) < Tl)
                 M(i,j) = 0;
                 
            elseif( any(any( F(i-1:i+1,j-1:j+1) > Th )) )
                 M(i,j) = 1;               
                
            end;
        end;
            
            
%         if N(i,j) == 3; %OK 3
%             if F(i,j) < F(i+1,j) || F(i,j) < F(i-1,j);
%                 M(i,j) = 0;
%             end;
%         end;
%         
%          if N(i,j) == 2; %I think this should be 2
%              if F(i,j) < F(i+1,j+1) || F(i,j) < F(i-1,j-1);
%                  M(i,j) = 0;
%             end;
%         end;
%          
%          if N(i,j) == 0; %OK
%              if F(i,j) < F(i,j+1) || F(i,j) < F(i,j-1);
%                  M(i,j) = 0;
%              end;        
%          end;
%          
%          if N(i,j) == 1;
%              if F(i,j) < F(i+1,j-1) || F(i,j) < F(i-1,j+1);
%                  M(i,j) = 0;
%              end;
%          end;    
        
    end;
end;

figure;
subplot(2,2, 1)
title('Gradient in x-direction');
colormap(gray);
imagesc(Fx);

subplot(2,2, 2)
title('Gradient in y-direction');
colormap(gray);
imagesc(Fy);

subplot(2,2, 3);
title('Gradient Magnitude');
colormap(gray);
imagesc(F);

subplot(2,2, 4);
title('Final Thinned and Thresholded Image');
colormap(gray);
imagesc(M);

figure;
colormap(gray);
imagesc(M);

function follow_chain(Mt)