function [vertices] = estimate_contours_vertices(label, max_error, max_delta, plot_images)
% -------------------------------------------------------------------------
% MSc Project - MSc in Communications and Signal Processing
% Department of Electrical and Electronic Engineering
% Imperial College London, 2011
%
% Supervisor  : Dr Pier Luigi Dragotti
% Contributor : Andriy Gelman
% Author      : Jon Onativia Bravo
%
% File        : estimate_contours_vertices.m
% -------------------------------------------------------------------------
% Estimate contour vertices for a given binary shape allowing a
% max_error distance between approximating polygon and original shape.
%
% [vertices] = estimate_contours_vertices(label, max_error[, max_delta, plot_images])
%
% INPUT:
%  - label       : Binary shape image
%  - max_error   : Maximum allowed distance between real contour and
%                  approximated polygon.
%  - max_delta   : Optional argument.
%                  Maximum offset in pixels for X and Y coordinates
%                  between two consecutive vertices. If not transmitted,
%                  no offset limit.
%  - plot_images : Optional boolean argument. Plot the images or not 
%                  (false by default).
%
% OUTPUT:
%  - vertices     : Array with the coordinates of the vertices of the
%                   approximating polygon.


%% Check input arguments
if nargin < 4
    plot_images = false;
end
if nargin < 3
    max_delta = max(size(label));
end


%% Output parameter
vertices = [];


%% Binarize label image
label(label ~= 0) = 1;


%% Open the label image
if plot_images == true
    figure
    imagesc(label)
end


%% Find the contour
figure
[c, ~] = contour(label, [0 0]);
close gcf;
c = round(c);


%% Determine the number of segments.
% The c vector is 2xN matrix, with first column [0; len_of_segment] and
% then len_of_segment columns  with first row contains x coordinate 
% and second row contains y coordinate of the points of the contour. If
% there are more contours another [0; len_of_segment] element will appear.
segment_inds = find(c(1,:) == 0);
seg_lengths  = c(2,segment_inds);
num_of_segs  = length(segment_inds);


%% Remove duplicated values
c1 = [];
for i_seg = 1 : num_of_segs
    % Retrieve the segment
    seg_len   = seg_lengths(i_seg);
    seg_start = segment_inds(i_seg) + 1;
    seg_end   = seg_start + seg_len - 1;
    seg       = c(:,seg_start:seg_end);
    
    % Remove duplicated values
%     [~, l]    = unique(seg', 'rows');
%     seg       = seg(:,sort(l));
    new_seg = seg(:,1);
    for i = 2 : seg_len-1
        if seg(1,i) ~= new_seg(1,end) || seg(2,i) ~= new_seg(2,end)
            new_seg = [new_seg, seg(:,i)];
        end
    end
    % Insert the last vertex only if it's different from previous and first
    if (seg(1,end) ~= new_seg(1,end) || seg(2,end) ~= new_seg(2,end)) ...
            && (seg(1,end) ~= new_seg(1,1)  || seg(2,end) ~= new_seg(2,1))
    	new_seg = [new_seg, seg(:,end)];
    end
    % Remove last vertex if it is equal to the first one
    if size(new_seg,2) > 1 ...
            && new_seg(1,end) == new_seg(1,1) && new_seg(2,end) == new_seg(2,1)
        new_seg = new_seg(:,1:end-1);
    end
    seg = new_seg;
    
    % Insert segment in vertices array
    c1  = [c1, [0; length(seg(1,:))]];
    c1  = [c1, seg];
end


%% Get the different new segments
segment_inds = find(c1(1,:) == 0);
seg_lengths  = c1(2,segment_inds);
num_of_segs  = length(segment_inds);


%% Estimate vertices

if max_error == 0
    vertices = c1;

else
    
    for i_seg = 1 : num_of_segs
        
        % Retrieve the segment
        seg_len   = seg_lengths(i_seg);
        seg_start = segment_inds(i_seg) + 1;
        seg_end   = seg_start + seg_len - 1;
        seg       = c1(:,seg_start:seg_end);
        
        cur_vert  = 1;
        next_vert = 2;
        
        new_seg = [];
        while cur_vert < seg_len
            
            new_seg = [ new_seg, [seg(1,cur_vert); seg(2,cur_vert)] ];
            
            if next_vert < seg_len
                d = get_max_distance([seg(1,cur_vert); seg(2,cur_vert)], ...
                                         [seg(1,next_vert+1); seg(2,next_vert+1)], ...
                                         seg(:,cur_vert+1:next_vert));
            end
            
            while d < max_error ...
                    && next_vert+1 < seg_len ...
                    && abs(seg(1,next_vert+1) - seg(1,cur_vert)) <= max_delta ...
                    && abs(seg(2,next_vert+1) - seg(2,cur_vert)) <= max_delta
                next_vert = next_vert + 1;
                d = get_max_distance([seg(1,cur_vert); seg(2,cur_vert)], ...
                                         [seg(1,next_vert+1); seg(2,next_vert+1)], ...
                                         seg(:,cur_vert+1:next_vert));
            end
            
            cur_vert  = next_vert;
            next_vert = next_vert + 1;
            
        end
        
        % If the original contour is not closed (the contour ends in the
        % border of the image), force the insertion of the last vertex so
        % the contour contains the vertices of the image border.
        if abs(seg(1,end) - seg(1,1)) > 1 || abs(seg(2,end) - seg(2,1)) > 1
            % Check if the last vertex of the original contour has already
            % been inserted
            if seg(1,end) ~= new_seg(1,end) || seg(2,end) ~= new_seg(2,end)
                new_seg = [ new_seg, [seg(1,end); seg(2,end)] ];
            end
        end
        
        % Insert the reduced segment in the vertices vector
        if ~isempty(new_seg)
            vertices = [vertices, [0 ; length(new_seg(1,:))], new_seg ];
        end
        
    end
    
end

segment_inds = find(vertices(1,:) == 0);
seg_lengths  = vertices(2,segment_inds);
num_of_segs  = length(segment_inds);

%% Check that the max delta between vertices satisfies the restriction
for i_seg = 1 : num_of_segs
    % Retrieve the segment
    seg_len   = seg_lengths(i_seg);
    seg_start = segment_inds(i_seg) + 1;
    seg_end   = seg_start + seg_len - 1;
    seg       = vertices(:,seg_start:seg_end);
        
    max_cur_delta = get_max_offset([ [0; size(seg, 2)] , seg ]);
    if max_cur_delta > max_delta
        disp('The maximum delta is bigger than the maximum allowed')
    end
end

%% Plot the vertices
if plot_images == true
    hold on
    for i_seg = 1 : num_of_segs

        % Retrieve the segment
        seg_len   = seg_lengths(i_seg);
        seg_start = segment_inds(i_seg) + 1;
        seg_end   = seg_start + seg_len - 1;
        seg       = vertices(:,seg_start:seg_end);

        % Plot it
        plot(seg(1,:), seg(2,:), 'LineWidth', 2, 'color', [0 0 0], 'Marker', 'x', 'MarkerSize', 8, 'MarkerEdgeColor', 'b');
    end
end
%%

end