Code covered by the BSD License  

Highlights from
distributePointsOnPolygon

image thumbnail

distributePointsOnPolygon

by

 

14 May 2013 (Updated )

This function distributes polygon points evenly across the polygon.

distributePointsOnPolygon( varargin )
function varargout = distributePointsOnPolygon( varargin )
% DISTRIBUTEPOINTSONPOLYGON this function uses faces and 2D/3D vertices as
% input. The vertices are redistributed over the polygon to decrease the
% variability in distance between the points.
% This is done by shifting the points iteratively over the original contour
% described by the polygon.
%
% Notes:
% - the polygon should be connected; i.e. the first and last point must be
%   connected by faces.
% - the vertices may not be connected > 2 times.
% - the resulting polygon might be slightly smoohted, depending on the
%   number of original points and the distance between those points; you can
%   avoid this by using more input points or adding points at the longest faces
%
% Usage:
% [ v ]    = distributePointsOnPolygon ( v )
% [ f, v ] = distributePointsOnPolygon ( f, v )
% [ f, v ] = distributePointsOnPolygon ( f, v, debug )
%
% use debug = 1 to see the progress.
% tweak the PARAMETRS on line 70-76 to speed up the process/ improve accuracy
%
%
% Demo:
% distributePointsOnPolygon()
%
%
% --- Author information
% Wouter Potters
% Academic Medical Center, Amsterdam, The Netherlands
% w.v.potters@amc.nl
% Date: 4-June-2013


% check input and output arguments
narginchk(0,3);
nargoutchk(0,2);

switch nargin
    case 0
        % demo time!!!
        demoFunction();
        return;
    case 1 % assume faces [1:end] and [2:end 1]
        v = varargin{1};
        f = [ (1:size(v,1))' [2:size(v,1) 1]' ];
        debug = 0;
    case 2
        f = varargin{1};
        v = varargin{2};
        debug = 0;
    case 3
        f = varargin{1};
        v = varargin{2};
        debug = varargin{3};
end

if debug == 1
    disp('Debugging on; please note that the function takes longer when debugging is on')
end

if ( ( size(v,2) ~= 2 ) && ( size(v,2) ~= 3 ) ) || ( size(f,2) ~= 2 )
    error('Vertices and faces should describe a 2D polygon (but may vertices may be 3D).')
end


v2 = v;  f2 = f; % save original arguments

% add a point to connect the endpoints if required
if (f2(end,end) ~= f2(1,1)) && (f2(end,end) ~= f2(1,end)) && (f2(end,1) ~= f2(1,1)) && (f2(end,1) ~= f2(1,end))
    f2(end+1,:) = [f2(end,end) f2(1,1)]; % add 1 connecting face
end

%%% PARAMETERS START
threshold          = 0.10; % initial threshold
threshold_steps    = 0.04; % steps for threshold
final_threshold    = 0.90; % final threshold
move_scaling       = 0.8;  % this scales the movement of a point to its next position
maximum_iterations = 100;
%%% PARAMETERS END

continue_updates = 1;   % while loop condition
iteration_number = 1;

while continue_updates % is aborted when dividing 2 distances
    % calculate length of faces = distance between points
    dist_f2 = (sqrt(sum((v2(f2(:,1),:) - v2(f2(:,2),:)).^2,2)));
    
    % calculate distances to neighboring points
    diff_v = zeros(length(v2),2);
    row_in_f = zeros(length(v2),2);
    for iv = 1:length(v2)
        [row,~] = find(f2 == iv);
        if length(row) ~= 2
            error('error length r ~= 2')
        end
        row_in_f(iv,:) = row; clear temp;
        diff_v(iv,:) = [dist_f2(row_in_f(iv,1)),dist_f2(row_in_f(iv,2))];
    end
    
    % Calculate relative difference between connected faces - force value
    % between 0 and 1
    divided_diff_v = diff_v(:,1) ./ diff_v(:,2);
    divided_diff_v( divided_diff_v > 1 ) = 1 ./ divided_diff_v( divided_diff_v > 1 );
    
    % define 'badly' distributed points based on threshold
    change_pts = find( divided_diff_v < threshold );
    
    % plot progress - show which points will change
    if debug == 1
        debugPlotFunction( v, v2, divided_diff_v, change_pts, threshold );
    end
    
    % calculate new position for current (ic) vertex
    for ic = change_pts'
        current_total_distance = sum( diff_v(ic,:) );
        smallest_distance      = min( diff_v(ic,:) );
        biggest_distance       = max( diff_v(ic,:) );
        
        selection = row_in_f( ic, ...
            diff_v(ic,:) == biggest_distance ); % select target positon based on largest face
        
        id_of_target = f2( selection, ...
            f2( selection , : ) ~= ic );            % get index of vertex from facelist
        
        
        distance_to_go = ( (current_total_distance / 2) - smallest_distance ) / biggest_distance;
        direction_to_go = ( v2( id_of_target, : ) - v2( ic, : ) );
        vector = move_scaling * distance_to_go * direction_to_go ; % direction and size of movement
        
        v2( ic, : ) = v2( ic, : ) + vector; % move v2(ic,:) according to vector
    end
    
    if threshold > final_threshold
        if ~any(divided_diff_v < threshold)
            continue_updates = 0; % abort loop
        end
    else
        threshold = threshold + threshold_steps; % increase threshold
    end
    
    if iteration_number > maximum_iterations
        error('Maximum iterations reached')
    else
        iteration_number = iteration_number + 1;
        if debug
            disp(['Now performing iteration #' num2str(iteration_number) ' - current threshold: ' num2str(threshold)])
        end
    end
    
end

if debug == 1
    debugPlotFunction( v, v2, divided_diff_v, [], threshold );
end

% set outputs
switch nargout
    case 1
        varargout{1} = v2;
    case 2
        varargout{1} = f2;
        varargout{2} = v2;
end

end % function distributePointsOnPolygon


%% ADDITIONAL FUNCTIONS
%- DEBUGPLOTFUNCTION
function fi = debugPlotFunction( v0, v, divided_diff_v, change_pts, threshold)


fi = figure(1); set(fi,'color','w')

sp1 = subplot(1,2,1); cla;
switch size(v,2)
    case 2
        plot(v0(:,1),v0(:,2),'g.'); hold on
        plot(v(:,1),v(:,2),'ro')
    case 3
        plot3(v0(:,1),v0(:,2),v0(:,3),'g.')
        plot3( v(:,1), v(:,2), v(:,3),'ro')
end

axis equal tight off
view(2)
hold on;

if ~isempty(change_pts)
    switch size(v,2)
        case 2
            plot(v(change_pts,1),v(change_pts,2),'b+')
        case 3
            plot3(v(change_pts,1),v(change_pts,2),v(change_pts,3),'b+')
    end
    lg = legend(sp1,'original points','current points','point subject to change in this iteration');
    set(lg,'Location','NorthOutside')
end

subplot(1,2,2); cla;
hist(divided_diff_v,10);
c = hist(divided_diff_v,10);
axis([0 1 0 1.5*max(c)])
hold on
plot([threshold threshold],[0 max(c)],'r--')
title('Histogram: length of 1st and 2st face divided')
drawnow; % update graph

end % function debugPlotFunction

%- DEMOFUNCTION
function demoFunction()
answer = questdlg('Do you wish to draw your own polygon points?','DEMOFUNCTION FOR DISTRIBUTEPOINTSONPOLYGON','Yes','No','No');
switch answer
    case 'Yes'
        fh = figure(1313); axis([0 10 0 10]);
        title('Please make a closed polygon of some points')
        [v] = impoly();
        try v = v.getPosition; catch, error('Poly failed; try again.'); end
        
        if length(v) < 5
            error('Please select at least 5 points')
        end
        
        % now increase the number of points to avoid smoothing of the
        % contour
        v = [ v; v(1,:)]; % add first point before interpolation
        v2 = zeros(length(v)*5,2);
        v2(:,1) = interp1(1:length(v),v(:,1),linspace(1,length(v),length(v)*5));
        v2(:,2) = interp1(1:length(v),v(:,2),linspace(1,length(v),length(v)*5));
        v = v2(1:end-1,:); % remove last point after interpolation to avoid duplicate points
        
        try close(fh); end % try to close the figure window.
        
        % create faces belonging to the vertices
        f = [(1:size(v,1))' [2:size(v,1) 1]' ];
        
    case 'No'
        v = [-21.2341   18.9143   -0.4923; -21.1714   19.3139   -0.4760
            -21.1404   18.6266   -0.5040;  -21.1404   19.3731   -0.4736
            -21.1046   18.5146   -0.5086
            -21.0156   18.1149   -0.5249
            -21.0111   19.7136   -0.4597
            -20.9628   17.7153   -0.5412
            -20.8586   17.3156   -0.5575
            -20.8479   20.1133   -0.4434
            -20.7404   17.1215   -0.5654
            -20.7404   20.2266   -0.4388
            -20.5938   16.9159   -0.5738
            -20.5716   20.5129   -0.4271
            -20.3404   16.5325   -0.5894
            -20.3404   20.7760   -0.4164
            -20.3235   16.5163   -0.5901
            -20.2054   20.9126   -0.4108
            -19.9404   16.2182   -0.6022
            -19.9404   21.1424   -0.4014
            -19.8198   16.1166   -0.6063
            -19.7167   21.3123   -0.3945
            -19.5404   15.9387   -0.6136
            -19.5404   21.4561   -0.3886
            -19.2952   15.7169   -0.6226
            -19.1404   15.5690   -0.6287
            -19.1404   21.5608   -0.3844
            -18.7894   15.3173   -0.6389
            -18.7404   15.2796   -0.6405
            -18.7404   21.5791   -0.3836
            -18.3404   15.1024   -0.6477
            -18.3404   21.6371   -0.3813
            -18.0091   21.7119   -0.3782
            -17.9404   15.0633   -0.6493
            -17.9404   21.7330   -0.3774
            -17.9078   21.7119   -0.3782
            -17.5404   15.0953   -0.6480
            -17.5404   21.5781   -0.3837
            -17.1404   15.1217   -0.6469
            -17.1404   21.5629   -0.3843
            -16.7404   15.1836   -0.6444
            -16.7404   21.5248   -0.3858
            -16.5346   15.3173   -0.6389
            -16.3404   15.4492   -0.6336
            -16.3404   21.4559   -0.3887
            -16.0155   21.3123   -0.3945
            -15.9404   15.6387   -0.6258
            -15.9404   21.2698   -0.3962
            -15.8272   15.7169   -0.6226
            -15.5404   15.9218   -0.6143
            -15.5404   21.0013   -0.4072
            -15.4566   20.9126   -0.4108
            -15.3008   16.1166   -0.6063
            -15.2483   20.5129   -0.4271
            -15.1404   16.2758   -0.5999
            -15.1404   20.3836   -0.4324
            -14.9339   16.5163   -0.5901
            -14.9029   20.1133   -0.4434
            -14.7521   16.9159   -0.5738
            -14.7404   16.9339   -0.5730
            -14.7404   19.9332   -0.4507
            -14.5494   19.7136   -0.4597
            -14.5325   17.3156   -0.5575
            -14.4762   17.7153   -0.5412
            -14.4442   18.1149   -0.5249
            -14.3404   18.4358   -0.5118
            -14.3404   19.4301   -0.4713
            -14.3124   18.5146   -0.5086
            -14.2728   19.3139   -0.4760
            -14.2667   18.9143   -0.4923];
        f = [14    16; 12    14
            20    22
            27    30
            32    33
            35    36
            38    40
            42    45
            46    48
            51    52
            54    56
            61    62
            69    70
            68    66
            66    65
            64    63
            63    60
            59    57
            55    53
            50    49
            44    43
            43    41
            39    37
            31    29
            29    28
            26    25
            23    21
            19    17
            13    11
            11     9
            3     1
            2     4
            1     2
            5     3
            6     5
            8     6
            9     8
            15    13
            17    15
            21    19
            25    23
            28    26
            34    31
            37    34
            41    39
            47    44
            49    47
            53    50
            57    55
            60    59
            65    64
            70    68
            67    69
            62    67
            58    61
            56    58
            52    54
            48    51
            45    46
            40    42
            36    38
            33    35
            30    32
            24    27
            22    24
            18    20
            16    18
            10    12
            4     7
            7    10];
end
distributePointsOnPolygon(f,v,1);
end % function demoFunction

Contact us