1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
|
## Copyright (C) 2015-2019 - Juan Pablo Carbajal
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 3 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program. If not, see <http://www.gnu.org/licenses/>.
## Author: Juan Pablo Carbajal <ajuanpi+dev@gmail.com>
## Updated: 2019-05-14
## -*- texinfo -*-
## @deftypefn {} {[@var{pline2} @var{idx}] = } simplifyPolyline_geometry (@var{pline})
## @deftypefnx {} {@dots{} = } simplifyPolyline_geometry (@dots{},@var{property},@var{value},@dots{})
## Simplify or subsample a polyline using the Ramer-Douglas-Peucker algorithm,
## a.k.a. the iterative end-point fit algorithm or the split-and-merge algorithm.
##
## The @var{pline} as a N-by-2 matrix. Rows correspond to the
## verices (compatible with @code{polygons2d}). The vector @var{idx} constains
## the indexes on vetices in @var{pline} that generates @var{pline2}, i.e.
## @code{pline2 = pline(idx,:)}.
##
## @strong{Parameters}
## @table @samp
## @item 'Nmax'
## Maximum number of vertices. Default value @code{1e3}.
## @item 'Tol'
## Tolerance for the error criteria. Default value @code{1e-4}.
## @item 'MaxIter'
## Maximum number of iterations. Default value @code{10}.
## @item 'Method'
## Not implemented.
## @end table
##
## Run @code{demo simplifyPolyline_geometry} to see an example.
##
## @seealso{curve2polyline, curveval, simplifyPolygon_geometry}
## @end deftypefn
function [pline idx] = simplifyPolyline_geometry (pline_o, varargin)
## TODO do not print warnings if user provided Nmax or MaxIter.
# --- Parse arguments --- #
parser = inputParser ();
parser.FunctionName = "simplifyPolyline_geometry";
toldef = 1e-4;#max(sumsq(diff(pline_o),2))*2;
parser.addParamValue ('Tol', toldef, @(x)x>0);
parser.addParamValue ('Nmax', 1e3, @(x)x>0);
parser.addParamValue ('MaxIter', 100, @(x)x>0);
parser.parse(varargin{:});
Nmax = parser.Results.Nmax;
tol = parser.Results.Tol;
MaxIter = parser.Results.MaxIter;
clear parser toldef
msg = ["simplifyPolyline_geometry: Maximum number of points reached with maximum error %g." ...
" Increase %s if the result is not satisfactory."];
# ------ #
[N dim] = size (pline_o);
idx = [1 N];
for iter = 1:MaxIter
# Find the point with the maximum distance.
[dist ii] = maxdistance (pline_o, idx);
tf = dist > tol;
n = sum (tf);
if ( all (!tf) );
break;
end
idx(end+1:end+n) = ii(tf);
idx = sort (idx);
if length(idx) >= Nmax
## TODO remove extra points
warning ('geometry:MayBeWrongOutput', sprintf (msg, max (dist),'Nmax'));
break;
end
end
if iter == MaxIter
warning ('geometry:MayBeWrongOutput', sprintf (msg ,max (dist), 'MaxIter'));
end
pline = pline_o(idx,:);
endfunction
function [dist ii] = maxdistance (p, idx)
## Separate the groups of points according to the edge they can divide.
idxc = arrayfun (@colon, idx(1:end-1), idx(2:end), "UniformOutput",false);
points = cellfun (@(x)p(x,:), idxc, "UniformOutput",false);
## Build the edges
edges = [p(idx(1:end-1),:) p(idx(2:end),:)];
edges = mat2cell (edges, ones (1,size (edges,1)), 4)';
## Calculate distance between the points and the corresponding edge
[dist ii] = cellfun (@dd, points, edges, idxc);
endfunction
function [dist ii] = dd (p,e,idx)
[d pos] = distancePointEdge (p, e);
[dist ii] = max (d);
ii = idx(ii);
endfunction
%!demo
%! t = linspace(0,1,100).';
%! y = polyval([1 -1.5 0.5 0],t);
%! pline = [t y];
%!
%! figure(1)
%! clf
%! plot (t,y,'-r;Original;','linewidth',2);
%! hold on
%!
%! tol = [8 2 1 0.5]*1e-2;
%! colors = jet(4);
%!
%! for i=1:4
%! pline_ = simplifyPolyline_geometry(pline,'tol',tol(i));
%! msg = sprintf('-;%g;',tol(i));
%! h = plot (pline_(:,1),pline_(:,2),msg);
%! set(h,'color',colors(i,:),'linewidth',2,'markersize',4);
%! end
%! hold off
%!
%! # ---------------------------------------------------------
%! # Four approximations of the initial polyline with decreasing tolerances.
%!demo
%! P = [0 0; 3 1; 3 4; 1 3; 2 2; 1 1];
%! func = @(x,y) linspace(x,y,5);
%! P2(:,1) = cell2mat( ...
%! arrayfun (func, P(1:end-1,1),P(2:end,1), ...
%! 'uniformoutput',false))'(:);
%! P2(:,2) = cell2mat( ...
%! arrayfun (func, P(1:end-1,2),P(2:end,2), ...
%! 'uniformoutput',false))'(:);
%!
%! P2s = simplifyPolyline_geometry (P2);
%!
%! plot(P(:,1),P(:,2),'s',P2(:,1),P2(:,2),'o',P2s(:,1),P2s(:,2),'-ok');
%!
%! # ---------------------------------------------------------
%! # Simplification of a polyline in the plane.
|
