1//////////////////////////////////////////////////////////////////////
2// LibFile: paths.scad
3// A `path` is a list of points of the same dimensions, usually 2D or 3D, that can
4// be connected together to form a sequence of line segments or a polygon.
5// A `region` is a list of paths that represent polygons, and the functions
6// in this file work on paths and also 1-regions, which are regions
7// that include exactly one path. When you pass a 1-region to a function, the default
8// value for `closed` is always `true` because regions represent polygons.
9// Capabilities include computing length of paths, computing
10// path tangents and normals, resampling of paths, and cutting paths up into smaller paths.
11// Includes:
12// include <BOSL2/std.scad>
13// FileGroup: Advanced Modeling
14// FileSummary: Operations on paths: length, resampling, tangents, splitting into subpaths
15// FileFootnotes: STD=Included in std.scad
16//////////////////////////////////////////////////////////////////////
17
18// Section: Utility Functions
19
20// Function: is_path()
21// Usage:
22// is_path(list, [dim], [fast])
23// Description:
24// Returns true if `list` is a path. A path is a list of two or more numeric vectors (AKA points).
25// All vectors must of the same size, and may only contain numbers that are not inf or nan.
26// By default the vectors in a path must be 2d or 3d. Set the `dim` parameter to specify a list
27// of allowed dimensions, or set it to `undef` to allow any dimension. (Note that this function
28// returns `false` on 1-regions.)
29// Example:
30// bool1 = is_path([[3,4],[5,6]]); // Returns true
31// bool2 = is_path([[3,4]]); // Returns false
32// bool3 = is_path([[3,4],[4,5]],2); // Returns true
33// bool4 = is_path([[3,4,3],[5,4,5]],2); // Returns false
34// bool5 = is_path([[3,4,3],[5,4,5]],2); // Returns false
35// bool6 = is_path([[3,4,5],undef,[4,5,6]]); // Returns false
36// bool7 = is_path([[3,5],[undef,undef],[4,5]]); // Returns false
37// bool8 = is_path([[3,4],[5,6],[5,3]]); // Returns true
38// bool9 = is_path([3,4,5,6,7,8]); // Returns false
39// bool10 = is_path([[3,4],[5,6]], dim=[2,3]);// Returns true
40// bool11 = is_path([[3,4],[5,6]], dim=[1,3]);// Returns false
41// bool12 = is_path([[3,4],"hello"], fast=true); // Returns true
42// bool13 = is_path([[3,4],[3,4,5]]); // Returns false
43// bool14 = is_path([[1,2,3,4],[2,3,4,5]]); // Returns false
44// bool15 = is_path([[1,2,3,4],[2,3,4,5]],undef);// Returns true
45// Arguments:
46// list = list to check
47// dim = list of allowed dimensions of the vectors in the path. Default: [2,3]
48// fast = set to true for fast check that only looks at first entry. Default: false
49function is_path(list, dim=[2,3], fast=false) =
50 fast
51 ? is_list(list) && is_vector(list[0])
52 : is_matrix(list)
53 && len(list)>1
54 && len(list[0])>0
55 && (is_undef(dim) || in_list(len(list[0]), force_list(dim)));
56
57// Function: is_1region()
58// Usage:
59// bool = is_1region(path, [name])
60// Description:
61// If `path` is a region with one component (a 1-region) then return true. If path is a region with more components
62// then display an error message about the parameter `name` requiring a path or a single component region. If the input
63// is not a region then return false. This function helps path functions accept 1-regions.
64// Arguments:
65// path = input to process
66// name = name of parameter to use in error message. Default: "path"
67function is_1region(path, name="path") =
68 !is_region(path)? false
69 :assert(len(path)==1,str("Parameter \"",name,"\" must be a path or singleton region, but is a multicomponent region"))
70 true;
71
72
73// Function: force_path()
74// Usage:
75// outpath = force_path(path, [name])
76// Description:
77// If `path` is a region with one component (a 1-region) then return that component as a path. If path is a region with more components
78// then display an error message about the parameter `name` requiring a path or a single component region. If the input
79// is not a region then return the input without any checks. This function helps path functions accept 1-regions.
80// Arguments:
81// path = input to process
82// name = name of parameter to use in error message. Default: "path"
83function force_path(path, name="path") =
84 is_region(path) ?
85 assert(len(path)==1, str("Parameter \"",name,"\" must be a path or singleton region, but is a multicomponent region"))
86 path[0]
87 : path;
88
89
90// Function: is_closed_path()
91// Usage:
92// is_closed_path(path, [eps]);
93// Description:
94// Returns true if the first and last points in the given path are coincident.
95function is_closed_path(path, eps=EPSILON) = approx(path[0], path[len(path)-1], eps=eps);
96
97
98// Function: close_path()
99// Usage:
100// close_path(path);
101// Description:
102// If a path's last point does not coincide with its first point, closes the path so it does.
103function close_path(path, eps=EPSILON) =
104 is_closed_path(path,eps=eps)? path : concat(path,[path[0]]);
105
106
107// Function: cleanup_path()
108// Usage:
109// cleanup_path(path);
110// Description:
111// If a path's last point coincides with its first point, deletes the last point in the path.
112function cleanup_path(path, eps=EPSILON) =
113 is_closed_path(path,eps=eps)? [for (i=[0:1:len(path)-2]) path[i]] : path;
114
115
116/// Internal Function: _path_select()
117/// Usage:
118/// _path_select(path,s1,u1,s2,u2,[closed]):
119/// Description:
120/// Returns a portion of a path, from between the `u1` part of segment `s1`, to the `u2` part of
121/// segment `s2`. Both `u1` and `u2` are values between 0.0 and 1.0, inclusive, where 0 is the start
122/// of the segment, and 1 is the end. Both `s1` and `s2` are integers, where 0 is the first segment.
123/// Arguments:
124/// path = The path to get a section of.
125/// s1 = The number of the starting segment.
126/// u1 = The proportion along the starting segment, between 0.0 and 1.0, inclusive.
127/// s2 = The number of the ending segment.
128/// u2 = The proportion along the ending segment, between 0.0 and 1.0, inclusive.
129/// closed = If true, treat path as a closed polygon.
130function _path_select(path, s1, u1, s2, u2, closed=false) =
131 let(
132 lp = len(path),
133 l = lp-(closed?0:1),
134 u1 = s1<0? 0 : s1>l? 1 : u1,
135 u2 = s2<0? 0 : s2>l? 1 : u2,
136 s1 = constrain(s1,0,l),
137 s2 = constrain(s2,0,l),
138 pathout = concat(
139 (s1<l && u1<1)? [lerp(path[s1],path[(s1+1)%lp],u1)] : [],
140 [for (i=[s1+1:1:s2]) path[i]],
141 (s2<l && u2>0)? [lerp(path[s2],path[(s2+1)%lp],u2)] : []
142 )
143 ) pathout;
144
145
146
147// Function: path_merge_collinear()
148// Description:
149// Takes a path and removes unnecessary sequential collinear points.
150// Usage:
151// path_merge_collinear(path, [eps])
152// Arguments:
153// path = A path of any dimension or a 1-region
154// closed = treat as closed polygon. Default: false
155// eps = Largest positional variance allowed. Default: `EPSILON` (1-e9)
156function path_merge_collinear(path, closed, eps=EPSILON) =
157 is_1region(path) ? path_merge_collinear(path[0], default(closed,true), eps) :
158 let(closed=default(closed,false))
159 assert(is_bool(closed))
160 assert( is_path(path), "Invalid path in path_merge_collinear." )
161 assert( is_undef(eps) || (is_finite(eps) && (eps>=0) ), "Invalid tolerance." )
162 len(path)<=2 ? path :
163 let(
164 indices = [
165 0,
166 for (i=[1:1:len(path)-(closed?1:2)])
167 if (!is_collinear(path[i-1], path[i], select(path,i+1), eps=eps)) i,
168 if (!closed) len(path)-1
169 ]
170 ) [for (i=indices) path[i]];
171
172
173
174// Section: Path length calculation
175
176
177// Function: path_length()
178// Usage:
179// path_length(path,[closed])
180// Description:
181// Returns the length of the path.
182// Arguments:
183// path = Path of any dimension or 1-region.
184// closed = true if the path is closed. Default: false
185// Example:
186// path = [[0,0], [5,35], [60,-25], [80,0]];
187// echo(path_length(path));
188function path_length(path,closed) =
189 is_1region(path) ? path_length(path[0], default(closed,true)) :
190 assert(is_path(path), "Invalid path in path_length")
191 let(closed=default(closed,false))
192 assert(is_bool(closed))
193 len(path)<2? 0 :
194 sum([for (i = [0:1:len(path)-2]) norm(path[i+1]-path[i])])+(closed?norm(path[len(path)-1]-path[0]):0);
195
196
197// Function: path_segment_lengths()
198// Usage:
199// path_segment_lengths(path,[closed])
200// Description:
201// Returns list of the length of each segment in a path
202// Arguments:
203// path = path in any dimension or 1-region
204// closed = true if the path is closed. Default: false
205function path_segment_lengths(path, closed) =
206 is_1region(path) ? path_segment_lengths(path[0], default(closed,true)) :
207 let(closed=default(closed,false))
208 assert(is_path(path),"Invalid path in path_segment_lengths.")
209 assert(is_bool(closed))
210 [
211 for (i=[0:1:len(path)-2]) norm(path[i+1]-path[i]),
212 if (closed) norm(path[0]-last(path))
213 ];
214
215
216// Function: path_length_fractions()
217// Usage:
218// fracs = path_length_fractions(path, [closed]);
219// Description:
220// Returns the distance fraction of each point in the path along the path, so the first
221// point is zero and the final point is 1. If the path is closed the length of the output
222// will have one extra point because of the final connecting segment that connects the last
223// point of the path to the first point.
224// Arguments:
225// path = path in any dimension or a 1-region
226// closed = set to true if path is closed. Default: false
227function path_length_fractions(path, closed) =
228 is_1region(path) ? path_length_fractions(path[0], default(closed,true)):
229 let(closed=default(closed, false))
230 assert(is_path(path))
231 assert(is_bool(closed))
232 let(
233 lengths = [
234 0,
235 each path_segment_lengths(path,closed)
236 ],
237 partial_len = cumsum(lengths),
238 total_len = last(partial_len)
239 )
240 partial_len / total_len;
241
242
243
244/// Internal Function: _path_self_intersections()
245/// Usage:
246/// isects = _path_self_intersections(path, [closed], [eps]);
247/// Description:
248/// Locates all self intersection points of the given path. Returns a list of intersections, where
249/// each intersection is a list like [POINT, SEGNUM1, PROPORTION1, SEGNUM2, PROPORTION2] where
250/// POINT is the coordinates of the intersection point, SEGNUMs are the integer indices of the
251/// intersecting segments along the path, and the PROPORTIONS are the 0.0 to 1.0 proportions
252/// of how far along those segments they intersect at. A proportion of 0.0 indicates the start
253/// of the segment, and a proportion of 1.0 indicates the end of the segment.
254/// .
255/// Note that this function does not return self-intersecting segments, only the points
256/// where non-parallel segments intersect.
257/// Arguments:
258/// path = The path to find self intersections of.
259/// closed = If true, treat path like a closed polygon. Default: true
260/// eps = The epsilon error value to determine whether two points coincide. Default: `EPSILON` (1e-9)
261/// Example(2D):
262/// path = [
263/// [-100,100], [0,-50], [100,100], [100,-100], [0,50], [-100,-100]
264/// ];
265/// isects = _path_self_intersections(path, closed=true);
266/// // isects == [[[-33.3333, 0], 0, 0.666667, 4, 0.333333], [[33.3333, 0], 1, 0.333333, 3, 0.666667]]
267/// stroke(path, closed=true, width=1);
268/// for (isect=isects) translate(isect[0]) color("blue") sphere(d=10);
269function _path_self_intersections(path, closed=true, eps=EPSILON) =
270 let(
271 path = closed ? close_path(path,eps=eps) : path,
272 plen = len(path)
273 )
274 [ for (i = [0:1:plen-3]) let(
275 a1 = path[i],
276 a2 = path[i+1],
277 seg_normal = unit([-(a2-a1).y, (a2-a1).x],[0,0]),
278 vals = path*seg_normal,
279 ref = a1*seg_normal,
280 // The value of vals[j]-ref is positive if vertex j is one one side of the
281 // line [a1,a2] and negative on the other side. Only a segment with opposite
282 // signs at its two vertices can have an intersection with segment
283 // [a1,a2]. The variable signals is zero when abs(vals[j]-ref) is less than
284 // eps and the sign of vals[j]-ref otherwise.
285 signals = [for(j=[i+2:1:plen-(i==0 && closed? 2: 1)])
286 abs(vals[j]-ref) < eps ? 0 : sign(vals[j]-ref) ]
287 )
288 if(max(signals)>=0 && min(signals)<=0 ) // some remaining edge intersects line [a1,a2]
289 for(j=[i+2:1:plen-(i==0 && closed? 3: 2)])
290 if( signals[j-i-2]*signals[j-i-1]<=0 ) let( // segm [b1,b2] intersects line [a1,a2]
291 b1 = path[j],
292 b2 = path[j+1],
293 isect = _general_line_intersection([a1,a2],[b1,b2],eps=eps)
294 )
295 if (isect
296 && isect[1]>=-eps
297 && isect[1]<= 1+eps
298 && isect[2]>= -eps
299 && isect[2]<= 1+eps)
300 [isect[0], i, isect[1], j, isect[2]]
301 ];
302
303// Section: Resampling—changing the number of points in a path
304
305
306// Input `data` is a list that sums to an integer.
307// Returns rounded version of input data so that every
308// entry is rounded to an integer and the sum is the same as
309// that of the input. Works by rounding an entry in the list
310// and passing the rounding error forward to the next entry.
311// This will generally distribute the error in a uniform manner.
312function _sum_preserving_round(data, index=0) =
313 index == len(data)-1 ? list_set(data, len(data)-1, round(data[len(data)-1])) :
314 let(
315 newval = round(data[index]),
316 error = newval - data[index]
317 ) _sum_preserving_round(
318 list_set(data, [index,index+1], [newval, data[index+1]-error]),
319 index+1
320 );
321
322
323// Function: subdivide_path()
324// See Also: subdivide_and_slice(), resample_path(), jittered_poly()
325// Usage:
326// newpath = subdivide_path(path, n|refine=|maxlen=, [method=], [closed=], [exact=]);
327// Description:
328// Takes a path as input (closed or open) and subdivides the path to produce a more
329// finely sampled path. You control the subdivision process by using the `maxlen` arg
330// to specify a maximum segment length, or by specifying `n` or `refine`, which request
331// a certain point count in the output.
332// .
333// You can specify the point count using the `n` option, where
334// you give the number of points you want in the output, or you can use
335// the `refine` option, where you specify a resampling factor. If `refine=3` then
336// the number of points would increase by a factor of three, so a four point square would
337// have 12 points after subdivision. With point-count subdivision, the new points can be distributed
338// proportional to length (`method="length"`), which is the default, or they can be divided up evenly among all the path segments
339// (`method="segment"`). If the extra points don't fit evenly on the path then the
340// algorithm attempts to distribute them as uniformly as possible, but the result may be uneven.
341// The `exact` option, which is true by default, requires that the final point count is
342// exactly as requested. For example, if you subdivide a four point square and request `n=13` then one edge will have
343// an extra point compared to the others.
344// If you set `exact=false` then the
345// algorithm will favor uniformity and the output path may have a different number of
346// points than you requested, but the sampling will be uniform. In our example of the
347// square with `n=13`, you will get only 12 points output, with the same number of points on each edge.
348// .
349// The points are always distributed uniformly on each segment. The `method="length"` option does
350// means that the number of points on a segment is based on its length, but the points are still
351// distributed uniformly on each segment, independent of the other segments.
352// With the `"segment"` method you can also give `n` as a vector of counts. This
353// specifies the desired point count on each segment: with vector valued `n` the `subdivide_path`
354// function places `n[i]-1` points on segment `i`. The reason for the -1 is to avoid
355// double counting the endpoints, which are shared by pairs of segments, so that for
356// a closed polygon the total number of points will be sum(n). Note that with an open
357// path there is an extra point at the end, so the number of points will be sum(n)+1.
358// .
359// If you use the `maxlen` option then you specify the maximum length segment allowed in the output.
360// Each segment is subdivided into the largest number of segments meeting your requirement. As above,
361// the sampling is uniform on each segment, independent of the other segments. With the `maxlen` option
362// you cannot specify `method` or `exact`.
363// Arguments:
364// path = path in any dimension or a 1-region
365// n = scalar total number of points desired or with `method="segment"` can be a vector requesting `n[i]-1` new points added to segment i.
366// ---
367// refine = increase total number of points by this factor (Specify only one of n, refine and maxlen)
368// maxlen = maximum length segment in the output (Specify only one of n, refine and maxlen)
369// closed = set to false if the path is open. Default: True
370// exact = if true return exactly the requested number of points, possibly sacrificing uniformity. If false, return uniform point sample that may not match the number of points requested. (Not allowed with maxlen.) Default: true
371// method = One of `"length"` or `"segment"`. If `"length"`, adds vertices in proportion to segment length, so short segments get fewer points. If `"segment"`, add points evenly among the segments, so all segments get the same number of points. (Not allowed with maxlen.) Default: `"length"`
372// Example(2D):
373// mypath = subdivide_path(square([2,2],center=true), 12);
374// move_copies(mypath)circle(r=.1,$fn=32);
375// Example(2D):
376// mypath = subdivide_path(square([8,2],center=true), 12);
377// move_copies(mypath)circle(r=.2,$fn=32);
378// Example(2D):
379// mypath = subdivide_path(square([8,2],center=true), 12, method="segment");
380// move_copies(mypath)circle(r=.2,$fn=32);
381// Example(2D):
382// mypath = subdivide_path(square([2,2],center=true), 17, closed=false);
383// move_copies(mypath)circle(r=.1,$fn=32);
384// Example(2D): Specifying different numbers of points on each segment
385// mypath = subdivide_path(hexagon(side=2), [2,3,4,5,6,7], method="segment");
386// move_copies(mypath)circle(r=.1,$fn=32);
387// Example(2D): Requested point total is 14 but 15 points output due to extra end point
388// mypath = subdivide_path(pentagon(side=2), [3,4,3,4], method="segment", closed=false);
389// move_copies(mypath)circle(r=.1,$fn=32);
390// Example(2D): Since 17 is not divisible by 5, a completely uniform distribution is not possible.
391// mypath = subdivide_path(pentagon(side=2), 17);
392// move_copies(mypath)circle(r=.1,$fn=32);
393// Example(2D): With `exact=false` a uniform distribution, but only 15 points
394// mypath = subdivide_path(pentagon(side=2), 17, exact=false);
395// move_copies(mypath)circle(r=.1,$fn=32);
396// Example(2D): With `exact=false` you can also get extra points, here 20 instead of requested 18
397// mypath = subdivide_path(pentagon(side=2), 18, exact=false);
398// move_copies(mypath)circle(r=.1,$fn=32);
399// Example(2D): Using refine in this example multiplies the point count by 3 by adding 2 points to each edge
400// mypath = subdivide_path(pentagon(side=2), refine=3);
401// move_copies(mypath)circle(r=.1,$fn=32);
402// Example(2D): But note that refine doesn't distribute evenly by segment unless you change the method. with the default method set to `"length"`, the points are distributed with more on the long segments in this example using refine.
403// mypath = subdivide_path(square([8,2],center=true), refine=3);
404// move_copies(mypath)circle(r=.2,$fn=32);
405// Example(2D): In this example with maxlen, every side gets a different number of new points
406// path = [[0,0],[0,4],[10,6],[10,0]];
407// spath = subdivide_path(path, maxlen=2, closed=true);
408// move_copies(spath) circle(r=.25,$fn=12);
409// Example(FlatSpin,VPD=15,VPT=[0,0,1.5]): Three-dimensional paths also work
410// mypath = subdivide_path([[0,0,0],[2,0,1],[2,3,2]], 12);
411// move_copies(mypath)sphere(r=.1,$fn=32);
412function subdivide_path(path, n, refine, maxlen, closed=true, exact, method) =
413 let(path = force_path(path))
414 assert(is_path(path))
415 assert(num_defined([n,refine,maxlen]),"Must give exactly one of n, refine, and maxlen")
416 refine==1 || n==len(path) ? path :
417 is_def(maxlen) ?
418 assert(is_undef(method), "Cannot give method with maxlen")
419 assert(is_undef(exact), "Cannot give exact with maxlen")
420 [
421 for (p=pair(path,closed))
422 let(steps = ceil(norm(p[1]-p[0])/maxlen))
423 each lerpn(p[0], p[1], steps, false),
424 if (!closed) last(path)
425 ]
426 :
427 let(
428 exact = default(exact, true),
429 method = default(method, "length")
430 )
431 assert(method=="length" || method=="segment")
432 let(
433 n = !is_undef(n)? n :
434 !is_undef(refine)? len(path) * refine :
435 undef
436 )
437 assert((is_num(n) && n>0) || is_vector(n),"Parameter n to subdivide_path must be postive number or vector")
438 let(
439 count = len(path) - (closed?0:1),
440 add_guess = method=="segment"?
441 (
442 is_list(n)
443 ? assert(len(n)==count,"Vector parameter n to subdivide_path has the wrong length")
444 add_scalar(n,-1)
445 : repeat((n-len(path)) / count, count)
446 )
447 : // method=="length"
448 assert(is_num(n),"Parameter n to subdivide path must be a number when method=\"length\"")
449 let(
450 path_lens = path_segment_lengths(path,closed),
451 add_density = (n - len(path)) / sum(path_lens)
452 )
453 path_lens * add_density,
454 add = exact? _sum_preserving_round(add_guess)
455 : [for (val=add_guess) round(val)]
456 )
457 [
458 for (i=[0:1:count-1])
459 each lerpn(path[i],select(path,i+1), 1+add[i],endpoint=false),
460 if (!closed) last(path)
461 ];
462
463
464
465
466// Function: resample_path()
467// Usage:
468// newpath = resample_path(path, n|spacing=, [closed=]);
469// Description:
470// Compute a uniform resampling of the input path. If you specify `n` then the output path will have n
471// points spaced uniformly (by linear interpolation along the input path segments). The only points of the
472// input path that are guaranteed to appear in the output path are the starting and ending points.
473// If you specify `spacing` then the length you give will be rounded to the nearest spacing that gives
474// a uniform sampling of the path and the resulting uniformly sampled path is returned.
475// Note that because this function operates on a discrete input path the quality of the output depends on
476// the sampling of the input. If you want very accurate output, use a lot of points for the input.
477// Arguments:
478// path = path in any dimension or a 1-region
479// n = Number of points in output
480// ---
481// spacing = Approximate spacing desired
482// closed = set to true if path is closed. Default: true
483// Example(2D): Subsampling lots of points from a smooth curve
484// path = xscale(2,circle($fn=250, r=10));
485// sampled = resample_path(path, 16);
486// stroke(path);
487// color("red")move_copies(sampled) circle($fn=16);
488// Example(2D): Specified spacing is rounded to make a uniform sampling
489// path = xscale(2,circle($fn=250, r=10));
490// sampled = resample_path(path, spacing=17);
491// stroke(path);
492// color("red")move_copies(sampled) circle($fn=16);
493// Example(2D): Notice that the corners are excluded
494// path = square(20);
495// sampled = resample_path(path, spacing=6);
496// stroke(path,closed=true);
497// color("red")move_copies(sampled) circle($fn=16);
498// Example(2D): Closed set to false
499// path = square(20);
500// sampled = resample_path(path, spacing=6,closed=false);
501// stroke(path);
502// color("red")move_copies(sampled) circle($fn=16);
503
504
505function resample_path(path, n, spacing, closed=true) =
506 let(path = force_path(path))
507 assert(is_path(path))
508 assert(num_defined([n,spacing])==1,"Must define exactly one of n and spacing")
509 assert(is_bool(closed))
510 let(
511 length = path_length(path,closed),
512 // In the open path case decrease n by 1 so that we don't try to get
513 // path_cut to return the endpoint (which might fail due to rounding)
514 // Add last point later
515 n = is_def(n) ? n-(closed?0:1) : round(length/spacing),
516 distlist = lerpn(0,length,n,false),
517 cuts = path_cut_points(path, distlist, closed=closed)
518 )
519 [ each column(cuts,0),
520 if (!closed) last(path) // Then add last point here
521 ];
522
523
524// Section: Path Geometry
525
526// Function: is_path_simple()
527// Usage:
528// bool = is_path_simple(path, [closed], [eps]);
529// Description:
530// Returns true if the given 2D path is simple, meaning that it has no self-intersections.
531// Repeated points are not considered self-intersections: a path with such points can
532// still be simple.
533// If closed is set to true then treat the path as a polygon.
534// Arguments:
535// path = 2D path or 1-region
536// closed = set to true to treat path as a polygon. Default: false
537// eps = Epsilon error value used for determine if points coincide. Default: `EPSILON` (1e-9)
538function is_path_simple(path, closed, eps=EPSILON) =
539 is_1region(path) ? is_path_simple(path[0], default(closed,true), eps) :
540 let(closed=default(closed,false))
541 assert(is_path(path, 2),"Must give a 2D path")
542 assert(is_bool(closed))
543 // check for path reversals
544 [for(i=[0:1:len(path)-(closed?2:3)])
545 let(v1=path[i+1]-path[i],
546 v2=select(path,i+2)-path[i+1],
547 normv1 = norm(v1),
548 normv2 = norm(v2)
549 )
550 if (approx(v1*v2/normv1/normv2,-1)) 1
551 ] == []
552 &&
553 _path_self_intersections(path,closed=closed,eps=eps) == [];
554
555
556// Function: path_closest_point()
557// Usage:
558// index_pt = path_closest_point(path, pt);
559// Description:
560// Finds the closest path segment, and point on that segment to the given point.
561// Returns `[SEGNUM, POINT]`
562// Arguments:
563// path = path of any dimension or a 1-region
564// pt = the point to find the closest point to
565// closed =
566// Example(2D):
567// path = circle(d=100,$fn=6);
568// pt = [20,10];
569// closest = path_closest_point(path, pt);
570// stroke(path, closed=true);
571// color("blue") translate(pt) circle(d=3, $fn=12);
572// color("red") translate(closest[1]) circle(d=3, $fn=12);
573function path_closest_point(path, pt, closed=true) =
574 let(path = force_path(path))
575 assert(is_path(path), "Input must be a path")
576 assert(is_vector(pt, len(path[0])), "Input pt must be a compatible vector")
577 assert(is_bool(closed))
578 let(
579 pts = [for (seg=pair(path,closed)) line_closest_point(seg,pt,SEGMENT)],
580 dists = [for (p=pts) norm(p-pt)],
581 min_seg = min_index(dists)
582 ) [min_seg, pts[min_seg]];
583
584
585// Function: path_tangents()
586// Usage:
587// tangs = path_tangents(path, [closed], [uniform]);
588// Description:
589// Compute the tangent vector to the input path. The derivative approximation is described in deriv().
590// The returns vectors will be normalized to length 1. If any derivatives are zero then
591// the function fails with an error. If you set `uniform` to false then the sampling is
592// assumed to be non-uniform and the derivative is computed with adjustments to produce corrected
593// values.
594// Arguments:
595// path = path of any dimension or a 1-region
596// closed = set to true of the path is closed. Default: false
597// uniform = set to false to correct for non-uniform sampling. Default: true
598// Example(2D): A shape with non-uniform sampling gives distorted derivatives that may be undesirable. Note that derivatives tilt towards the long edges of the rectangle.
599// rect = square([10,3]);
600// tangents = path_tangents(rect,closed=true);
601// stroke(rect,closed=true, width=0.25);
602// color("purple")
603// for(i=[0:len(tangents)-1])
604// stroke([rect[i]-tangents[i], rect[i]+tangents[i]],width=.25, endcap2="arrow2");
605// Example(2D): Setting uniform to false corrects the distorted derivatives for this example:
606// rect = square([10,3]);
607// tangents = path_tangents(rect,closed=true,uniform=false);
608// stroke(rect,closed=true, width=0.25);
609// color("purple")
610// for(i=[0:len(tangents)-1])
611// stroke([rect[i]-tangents[i], rect[i]+tangents[i]],width=.25, endcap2="arrow2");
612function path_tangents(path, closed, uniform=true) =
613 is_1region(path) ? path_tangents(path[0], default(closed,true), uniform) :
614 let(closed=default(closed,false))
615 assert(is_bool(closed))
616 assert(is_path(path))
617 !uniform ? [for(t=deriv(path,closed=closed, h=path_segment_lengths(path,closed))) unit(t)]
618 : [for(t=deriv(path,closed=closed)) unit(t)];
619
620
621// Function: path_normals()
622// Usage:
623// norms = path_normals(path, [tangents], [closed]);
624// Description:
625// Compute the normal vector to the input path. This vector is perpendicular to the
626// path tangent and lies in the plane of the curve. For 3d paths we define the plane of the curve
627// at path point i to be the plane defined by point i and its two neighbors. At the endpoints of open paths
628// we use the three end points. For 3d paths the computed normal is the one lying in this plane that points
629// towards the center of curvature at that path point. For 2d paths, which lie in the xy plane, the normal
630// is the path pointing to the right of the direction the path is traveling. If points are collinear then
631// a 3d path has no center of curvature, and hence the
632// normal is not uniquely defined. In this case the function issues an error.
633// For 2d paths the plane is always defined so the normal fails to exist only
634// when the derivative is zero (in the case of repeated points).
635// Arguments:
636// path = 2D or 3D path or a 1-region
637// tangents = path tangents optionally supplied
638// closed = if true path is treated as a polygon. Default: false
639function path_normals(path, tangents, closed) =
640 is_1region(path) ? path_normals(path[0], tangents, default(closed,true)) :
641 let(closed=default(closed,false))
642 assert(is_path(path,[2,3]))
643 assert(is_bool(closed))
644 let(
645 tangents = default(tangents, path_tangents(path,closed)),
646 dim=len(path[0])
647 )
648 assert(is_path(tangents) && len(tangents[0])==dim,"Dimensions of path and tangents must match")
649 [
650 for(i=idx(path))
651 let(
652 pts = i==0 ? (closed? select(path,-1,1) : select(path,0,2))
653 : i==len(path)-1 ? (closed? select(path,i-1,i+1) : select(path,i-2,i))
654 : select(path,i-1,i+1)
655 )
656 dim == 2 ? [tangents[i].y,-tangents[i].x]
657 : let( v=cross(cross(pts[1]-pts[0], pts[2]-pts[0]),tangents[i]))
658 assert(norm(v)>EPSILON, "3D path contains collinear points")
659 unit(v)
660 ];
661
662
663// Function: path_curvature()
664// Usage:
665// curvs = path_curvature(path, [closed]);
666// Description:
667// Numerically estimate the curvature of the path (in any dimension).
668// Arguments:
669// path = path in any dimension or a 1-region
670// closed = if true then treat the path as a polygon. Default: false
671function path_curvature(path, closed) =
672 is_1region(path) ? path_curvature(path[0], default(closed,true)) :
673 let(closed=default(closed,false))
674 assert(is_bool(closed))
675 assert(is_path(path))
676 let(
677 d1 = deriv(path, closed=closed),
678 d2 = deriv2(path, closed=closed)
679 ) [
680 for(i=idx(path))
681 sqrt(
682 sqr(norm(d1[i])*norm(d2[i])) -
683 sqr(d1[i]*d2[i])
684 ) / pow(norm(d1[i]),3)
685 ];
686
687
688// Function: path_torsion()
689// Usage:
690// torsions = path_torsion(path, [closed]);
691// Description:
692// Numerically estimate the torsion of a 3d path.
693// Arguments:
694// path = 3D path
695// closed = if true then treat path as a polygon. Default: false
696function path_torsion(path, closed=false) =
697 assert(is_path(path,3), "Input path must be a 3d path")
698 assert(is_bool(closed))
699 let(
700 d1 = deriv(path,closed=closed),
701 d2 = deriv2(path,closed=closed),
702 d3 = deriv3(path,closed=closed)
703 ) [
704 for (i=idx(path)) let(
705 crossterm = cross(d1[i],d2[i])
706 ) crossterm * d3[i] / sqr(norm(crossterm))
707 ];
708
709
710// Section: Breaking paths up into subpaths
711
712
713
714// Function: path_cut()
715// Topics: Paths
716// See Also: split_path_at_self_crossings()
717// Usage:
718// path_list = path_cut(path, cutdist, [closed]);
719// Description:
720// Given a list of distances in `cutdist`, cut the path into
721// subpaths at those lengths, returning a list of paths.
722// If the input path is closed then the final path will include the
723// original starting point. The list of cut distances must be
724// in ascending order and should not include the endpoints: 0
725// or len(path). If you repeat a distance you will get an
726// empty list in that position in the output. If you give an
727// empty cutdist array you will get the input path as output
728// (without the final vertex doubled in the case of a closed path).
729// Arguments:
730// path = path of any dimension or a 1-region
731// cutdist = Distance or list of distances where path is cut
732// closed = If true, treat the path as a closed polygon. Default: false
733// Example(2D,NoAxes):
734// path = circle(d=100);
735// segs = path_cut(path, [50, 200], closed=true);
736// rainbow(segs) stroke($item, endcaps="butt", width=3);
737function path_cut(path,cutdist,closed) =
738 is_num(cutdist) ? path_cut(path,[cutdist],closed) :
739 is_1region(path) ? path_cut(path[0], cutdist, default(closed,true)):
740 let(closed=default(closed,false))
741 assert(is_bool(closed))
742 assert(is_vector(cutdist))
743 assert(last(cutdist)<path_length(path,closed=closed),"Cut distances must be smaller than the path length")
744 assert(cutdist[0]>0, "Cut distances must be strictly positive")
745 let(
746 cutlist = path_cut_points(path,cutdist,closed=closed)
747 )
748 _path_cut_getpaths(path, cutlist, closed);
749
750
751function _path_cut_getpaths(path, cutlist, closed) =
752 let(
753 cuts = len(cutlist)
754 )
755 [
756 [ each list_head(path,cutlist[0][1]-1),
757 if (!approx(cutlist[0][0], path[cutlist[0][1]-1])) cutlist[0][0]
758 ],
759 for(i=[0:1:cuts-2])
760 cutlist[i][0]==cutlist[i+1][0] && cutlist[i][1]==cutlist[i+1][1] ? []
761 :
762 [ if (!approx(cutlist[i][0], select(path,cutlist[i][1]))) cutlist[i][0],
763 each slice(path, cutlist[i][1], cutlist[i+1][1]-1),
764 if (!approx(cutlist[i+1][0], select(path,cutlist[i+1][1]-1))) cutlist[i+1][0],
765 ],
766 [
767 if (!approx(cutlist[cuts-1][0], select(path,cutlist[cuts-1][1]))) cutlist[cuts-1][0],
768 each select(path,cutlist[cuts-1][1],closed ? 0 : -1)
769 ]
770 ];
771
772
773
774// Function: path_cut_points()
775//
776// Usage:
777// cuts = path_cut_points(path, cutdist, [closed=], [direction=]);
778//
779// Description:
780// Cuts a path at a list of distances from the first point in the path. Returns a list of the cut
781// points and indices of the next point in the path after that point. So for example, a return
782// value entry of [[2,3], 5] means that the cut point was [2,3] and the next point on the path after
783// this point is path[5]. If the path is too short then path_cut_points returns undef. If you set
784// `direction` to true then `path_cut_points` will also return the tangent vector to the path and a normal
785// vector to the path. It tries to find a normal vector that is coplanar to the path near the cut
786// point. If this fails it will return a normal vector parallel to the xy plane. The output with
787// direction vectors will be `[point, next_index, tangent, normal]`.
788// .
789// If you give the very last point of the path as a cut point then the returned index will be
790// one larger than the last index (so it will not be a valid index). If you use the closed
791// option then the returned index will be equal to the path length for cuts along the closing
792// path segment, and if you give a point equal to the path length you will get an
793// index of len(path)+1 for the index.
794//
795// Arguments:
796// path = path to cut
797// cutdist = distances where the path should be cut (a list) or a scalar single distance
798// ---
799// closed = set to true if the curve is closed. Default: false
800// direction = set to true to return direction vectors. Default: false
801//
802// Example(NORENDER):
803// square=[[0,0],[1,0],[1,1],[0,1]];
804// path_cut_points(square, [.5,1.5,2.5]); // Returns [[[0.5, 0], 1], [[1, 0.5], 2], [[0.5, 1], 3]]
805// path_cut_points(square, [0,1,2,3]); // Returns [[[0, 0], 1], [[1, 0], 2], [[1, 1], 3], [[0, 1], 4]]
806// path_cut_points(square, [0,0.8,1.6,2.4,3.2], closed=true); // Returns [[[0, 0], 1], [[0.8, 0], 1], [[1, 0.6], 2], [[0.6, 1], 3], [[0, 0.8], 4]]
807// path_cut_points(square, [0,0.8,1.6,2.4,3.2]); // Returns [[[0, 0], 1], [[0.8, 0], 1], [[1, 0.6], 2], [[0.6, 1], 3], undef]
808function path_cut_points(path, cutdist, closed=false, direction=false) =
809 let(long_enough = len(path) >= (closed ? 3 : 2))
810 assert(long_enough,len(path)<2 ? "Two points needed to define a path" : "Closed path must include three points")
811 is_num(cutdist) ? path_cut_points(path, [cutdist],closed, direction)[0] :
812 assert(is_vector(cutdist))
813 assert(is_increasing(cutdist), "Cut distances must be an increasing list")
814 let(cuts = path_cut_points_recurse(path,cutdist,closed))
815 !direction
816 ? cuts
817 : let(
818 dir = _path_cuts_dir(path, cuts, closed),
819 normals = _path_cuts_normals(path, cuts, dir, closed)
820 )
821 hstack(cuts, list_to_matrix(dir,1), list_to_matrix(normals,1));
822
823// Main recursive path cut function
824function path_cut_points_recurse(path, dists, closed=false, pind=0, dtotal=0, dind=0, result=[]) =
825 dind == len(dists) ? result :
826 let(
827 lastpt = len(result)==0? [] : last(result)[0], // location of last cut point
828 dpartial = len(result)==0? 0 : norm(lastpt-select(path,pind)), // remaining length in segment
829 nextpoint = dists[dind] < dpartial+dtotal // Do we have enough length left on the current segment?
830 ? [lerp(lastpt,select(path,pind),(dists[dind]-dtotal)/dpartial),pind]
831 : _path_cut_single(path, dists[dind]-dtotal-dpartial, closed, pind)
832 )
833 path_cut_points_recurse(path, dists, closed, nextpoint[1], dists[dind],dind+1, concat(result, [nextpoint]));
834
835
836// Search for a single cut point in the path
837function _path_cut_single(path, dist, closed=false, ind=0, eps=1e-7) =
838 // If we get to the very end of the path (ind is last point or wraparound for closed case) then
839 // check if we are within epsilon of the final path point. If not we're out of path, so we fail
840 ind==len(path)-(closed?0:1) ?
841 assert(dist<eps,"Path is too short for specified cut distance")
842 [select(path,ind),ind+1]
843 :let(d = norm(path[ind]-select(path,ind+1))) d > dist ?
844 [lerp(path[ind],select(path,ind+1),dist/d), ind+1] :
845 _path_cut_single(path, dist-d,closed, ind+1, eps);
846
847// Find normal directions to the path, coplanar to local part of the path
848// Or return a vector parallel to the x-y plane if the above fails
849function _path_cuts_normals(path, cuts, dirs, closed=false) =
850 [for(i=[0:len(cuts)-1])
851 len(path[0])==2? [-dirs[i].y, dirs[i].x]
852 :
853 let(
854 plane = len(path)<3 ? undef :
855 let(start = max(min(cuts[i][1],len(path)-1),2)) _path_plane(path, start, start-2)
856 )
857 plane==undef?
858 ( dirs[i].x==0 && dirs[i].y==0 ? [1,0,0] // If it's z direction return x vector
859 : unit([-dirs[i].y, dirs[i].x,0])) // otherwise perpendicular to projection
860 : unit(cross(dirs[i],cross(plane[0],plane[1])))
861 ];
862
863// Scan from the specified point (ind) to find a noncoplanar triple to use
864// to define the plane of the path.
865function _path_plane(path, ind, i,closed) =
866 i<(closed?-1:0) ? undef :
867 !is_collinear(path[ind],path[ind-1], select(path,i))?
868 [select(path,i)-path[ind-1],path[ind]-path[ind-1]] :
869 _path_plane(path, ind, i-1);
870
871// Find the direction of the path at the cut points
872function _path_cuts_dir(path, cuts, closed=false, eps=1e-2) =
873 [for(ind=[0:len(cuts)-1])
874 let(
875 zeros = path[0]*0,
876 nextind = cuts[ind][1],
877 nextpath = unit(select(path, nextind+1)-select(path, nextind),zeros),
878 thispath = unit(select(path, nextind) - select(path,nextind-1),zeros),
879 lastpath = unit(select(path,nextind-1) - select(path, nextind-2),zeros),
880 nextdir =
881 nextind==len(path) && !closed? lastpath :
882 (nextind<=len(path)-2 || closed) && approx(cuts[ind][0], path[nextind],eps)
883 ? unit(nextpath+thispath)
884 : (nextind>1 || closed) && approx(cuts[ind][0],select(path,nextind-1),eps)
885 ? unit(thispath+lastpath)
886 : thispath
887 ) nextdir
888 ];
889
890
891// internal function
892// converts pathcut output form to a [segment, u]
893// form list that works withi path_select
894function _cut_to_seg_u_form(pathcut, path, closed) =
895 let(lastind = len(path) - (closed?0:1))
896 [for(entry=pathcut)
897 entry[1] > lastind ? [lastind,0] :
898 let(
899 a = path[entry[1]-1],
900 b = path[entry[1]],
901 c = entry[0],
902 i = max_index(v_abs(b-a)),
903 factor = (c[i]-a[i])/(b[i]-a[i])
904 )
905 [entry[1]-1,factor]
906 ];
907
908
909
910// Function: split_path_at_self_crossings()
911// Usage:
912// paths = split_path_at_self_crossings(path, [closed], [eps]);
913// Description:
914// Splits a 2D path into sub-paths wherever the original path crosses itself.
915// Splits may occur mid-segment, so new vertices will be created at the intersection points.
916// Returns a list of the resulting subpaths.
917// Arguments:
918// path = A 2D path or a 1-region.
919// closed = If true, treat path as a closed polygon. Default: true
920// eps = Acceptable variance. Default: `EPSILON` (1e-9)
921// Example(2D,NoAxes):
922// path = [ [-100,100], [0,-50], [100,100], [100,-100], [0,50], [-100,-100] ];
923// paths = split_path_at_self_crossings(path);
924// rainbow(paths) stroke($item, closed=false, width=3);
925function split_path_at_self_crossings(path, closed=true, eps=EPSILON) =
926 let(path = force_path(path))
927 assert(is_path(path,2), "Must give a 2D path")
928 assert(is_bool(closed))
929 let(
930 path = cleanup_path(path, eps=eps),
931 isects = deduplicate(
932 eps=eps,
933 concat(
934 [[0, 0]],
935 sort([
936 for (
937 a = _path_self_intersections(path, closed=closed, eps=eps),
938 ss = [ [a[1],a[2]], [a[3],a[4]] ]
939 ) if (ss[0] != undef) ss
940 ]),
941 [[len(path)-(closed?1:2), 1]]
942 )
943 )
944 ) [
945 for (p = pair(isects))
946 let(
947 s1 = p[0][0],
948 u1 = p[0][1],
949 s2 = p[1][0],
950 u2 = p[1][1],
951 section = _path_select(path, s1, u1, s2, u2, closed=closed),
952 outpath = deduplicate(eps=eps, section)
953 )
954 if (len(outpath)>1) outpath
955 ];
956
957
958function _tag_self_crossing_subpaths(path, nonzero, closed=true, eps=EPSILON) =
959 let(
960 subpaths = split_path_at_self_crossings(
961 path, closed=true, eps=eps
962 )
963 ) [
964 for (subpath = subpaths) let(
965 seg = select(subpath,0,1),
966 mp = mean(seg),
967 n = line_normal(seg) / 2048,
968 p1 = mp + n,
969 p2 = mp - n,
970 p1in = point_in_polygon(p1, path, nonzero=nonzero) >= 0,
971 p2in = point_in_polygon(p2, path, nonzero=nonzero) >= 0,
972 tag = (p1in && p2in)? "I" : "O"
973 ) [tag, subpath]
974 ];
975
976
977// Function: polygon_parts()
978// Usage:
979// splitpolys = polygon_parts(poly, [nonzero], [eps]);
980// Description:
981// Given a possibly self-intersecting 2d polygon, constructs a representation of the original polygon as a list of
982// non-intersecting simple polygons. If nonzero is set to true then it uses the nonzero method for defining polygon membership.
983// For simple cases, such as the pentagram, this will produce the outer perimeter of a self-intersecting polygon.
984// Arguments:
985// poly = a 2D polygon or 1-region
986// nonzero = If true use the nonzero method for checking if a point is in a polygon. Otherwise use the even-odd method. Default: false
987// eps = The epsilon error value to determine whether two points coincide. Default: `EPSILON` (1e-9)
988// Example(2D,NoAxes): This cross-crossing polygon breaks up into its 3 components (regardless of the value of nonzero).
989// poly = [
990// [-100,100], [0,-50], [100,100],
991// [100,-100], [0,50], [-100,-100]
992// ];
993// splitpolys = polygon_parts(poly);
994// rainbow(splitpolys) stroke($item, closed=true, width=3);
995// Example(2D,NoAxes): With nonzero=false you get even-odd mode which matches OpenSCAD, so the pentagram breaks apart into its five points.
996// pentagram = turtle(["move",100,"left",144], repeat=4);
997// left(100)polygon(pentagram);
998// rainbow(polygon_parts(pentagram,nonzero=false))
999// stroke($item,closed=true,width=2.5);
1000// Example(2D,NoAxes): With nonzero=true you get only the outer perimeter. You can use this to create the polygon using the nonzero method, which is not supported by OpenSCAD.
1001// pentagram = turtle(["move",100,"left",144], repeat=4);
1002// outside = polygon_parts(pentagram,nonzero=true);
1003// left(100)region(outside);
1004// rainbow(outside)
1005// stroke($item,closed=true,width=2.5);
1006// Example(2D,NoAxes):
1007// N=12;
1008// ang=360/N;
1009// sr=10;
1010// poly = turtle(["angle", 90+ang/2,
1011// "move", sr, "left",
1012// "move", 2*sr*sin(ang/2), "left",
1013// "repeat", 4,
1014// ["move", 2*sr, "left",
1015// "move", 2*sr*sin(ang/2), "left"],
1016// "move", sr]);
1017// stroke(poly, width=.3);
1018// right(20)rainbow(polygon_parts(poly)) polygon($item);
1019// Example(2D,NoAxes): overlapping poly segments disappear
1020// poly = [[0,0], [10,0], [10,10], [0,10],[0,20], [20,10],[10,10], [0,10],[0,0]];
1021// stroke(poly,width=0.3);
1022// right(22)stroke(polygon_parts(poly)[0], width=0.3, closed=true);
1023// Example(2D,NoAxes): Poly segments disappear outside as well
1024// poly = turtle(["repeat", 3, ["move", 17, "left", "move", 10, "left", "move", 7, "left", "move", 10, "left"]]);
1025// back(2)stroke(poly,width=.5);
1026// fwd(12)rainbow(polygon_parts(poly)) stroke($item, closed=true, width=0.5);
1027// Example(2D,NoAxes): This shape has six components
1028// poly = turtle(["repeat", 3, ["move", 15, "left", "move", 7, "left", "move", 10, "left", "move", 17, "left"]]);
1029// polygon(poly);
1030// right(22)rainbow(polygon_parts(poly)) polygon($item);
1031// Example(2D,NoAxes): When the loops of the shape overlap then nonzero gives a different result than the even-odd method.
1032// poly = turtle(["repeat", 3, ["move", 15, "left", "move", 7, "left", "move", 10, "left", "move", 10, "left"]]);
1033// polygon(poly);
1034// right(27)rainbow(polygon_parts(poly)) polygon($item);
1035// move([16,-14])rainbow(polygon_parts(poly,nonzero=true)) polygon($item);
1036function polygon_parts(poly, nonzero=false, eps=EPSILON) =
1037 let(poly = force_path(poly))
1038 assert(is_path(poly,2), "Must give 2D polygon")
1039 assert(is_bool(nonzero))
1040 let(
1041 poly = cleanup_path(poly, eps=eps),
1042 tagged = _tag_self_crossing_subpaths(poly, nonzero=nonzero, closed=true, eps=eps),
1043 kept = [for (sub = tagged) if(sub[0] == "O") sub[1]],
1044 outregion = _assemble_path_fragments(kept, eps=eps)
1045 ) outregion;
1046
1047
1048function _extreme_angle_fragment(seg, fragments, rightmost=true, eps=EPSILON) =
1049 !fragments? [undef, []] :
1050 let(
1051 delta = seg[1] - seg[0],
1052 segang = atan2(delta.y,delta.x),
1053 frags = [
1054 for (i = idx(fragments)) let(
1055 fragment = fragments[i],
1056 fwdmatch = approx(seg[1], fragment[0], eps=eps),
1057 bakmatch = approx(seg[1], last(fragment), eps=eps)
1058 ) [
1059 fwdmatch,
1060 bakmatch,
1061 bakmatch? reverse(fragment) : fragment
1062 ]
1063 ],
1064 angs = [
1065 for (frag = frags)
1066 (frag[0] || frag[1])? let(
1067 delta2 = frag[2][1] - frag[2][0],
1068 segang2 = atan2(delta2.y, delta2.x)
1069 ) modang(segang2 - segang) : (
1070 rightmost? 999 : -999
1071 )
1072 ],
1073 fi = rightmost? min_index(angs) : max_index(angs)
1074 ) abs(angs[fi]) > 360? [undef, fragments] : let(
1075 remainder = [for (i=idx(fragments)) if (i!=fi) fragments[i]],
1076 frag = frags[fi],
1077 foundfrag = frag[2]
1078 ) [foundfrag, remainder];
1079
1080
1081/// Internal Function: _assemble_a_path_from_fragments()
1082/// Usage:
1083/// _assemble_a_path_from_fragments(subpaths);
1084/// Description:
1085/// Given a list of paths, assembles them together into one complete closed polygon path, and
1086/// remainder fragments. Returns [PATH, FRAGMENTS] where FRAGMENTS is the list of remaining
1087/// unused path fragments.
1088/// Arguments:
1089/// fragments = List of paths to be assembled into complete polygons.
1090/// rightmost = If true, assemble paths using rightmost turns. Leftmost if false.
1091/// startfrag = The fragment to start with. Default: 0
1092/// eps = The epsilon error value to determine whether two points coincide. Default: `EPSILON` (1e-9)
1093function _assemble_a_path_from_fragments(fragments, rightmost=true, startfrag=0, eps=EPSILON) =
1094 len(fragments)==0? [[],[]] :
1095 len(fragments)==1? [fragments[0],[]] :
1096 let(
1097 path = fragments[startfrag],
1098 newfrags = [for (i=idx(fragments)) if (i!=startfrag) fragments[i]]
1099 ) is_closed_path(path, eps=eps)? (
1100 // starting fragment is already closed
1101 [path, newfrags]
1102 ) : let(
1103 // Find rightmost/leftmost continuation fragment
1104 seg = select(path,-2,-1),
1105 extrema = _extreme_angle_fragment(seg=seg, fragments=newfrags, rightmost=rightmost, eps=eps),
1106 foundfrag = extrema[0],
1107 remainder = extrema[1]
1108 ) is_undef(foundfrag)? (
1109 // No remaining fragments connect! INCOMPLETE PATH!
1110 // Treat it as complete.
1111 [path, remainder]
1112 ) : is_closed_path(foundfrag, eps=eps)? (
1113 // Found fragment is already closed
1114 [foundfrag, concat([path], remainder)]
1115 ) : let(
1116 fragend = last(foundfrag),
1117 hits = [for (i = idx(path,e=-2)) if(approx(path[i],fragend,eps=eps)) i]
1118 ) hits? (
1119 let(
1120 // Found fragment intersects with initial path
1121 hitidx = last(hits),
1122 newpath = list_head(path,hitidx),
1123 newfrags = concat(len(newpath)>1? [newpath] : [], remainder),
1124 outpath = concat(slice(path,hitidx,-2), foundfrag)
1125 )
1126 [outpath, newfrags]
1127 ) : let(
1128 // Path still incomplete. Continue building it.
1129 newpath = concat(path, list_tail(foundfrag)),
1130 newfrags = concat([newpath], remainder)
1131 )
1132 _assemble_a_path_from_fragments(
1133 fragments=newfrags,
1134 rightmost=rightmost,
1135 eps=eps
1136 );
1137
1138
1139/// Internal Function: _assemble_path_fragments()
1140/// Usage:
1141/// _assemble_path_fragments(subpaths);
1142/// Description:
1143/// Given a list of paths, assembles them together into complete closed polygon paths if it can.
1144/// Polygons with area < eps will be discarded and not returned.
1145/// Arguments:
1146/// fragments = List of paths to be assembled into complete polygons.
1147/// eps = The epsilon error value to determine whether two points coincide. Default: `EPSILON` (1e-9)
1148function _assemble_path_fragments(fragments, eps=EPSILON, _finished=[]) =
1149 len(fragments)==0? _finished :
1150 let(
1151 minxidx = min_index([
1152 for (frag=fragments) min(column(frag,0))
1153 ]),
1154 result_l = _assemble_a_path_from_fragments(
1155 fragments=fragments,
1156 startfrag=minxidx,
1157 rightmost=false,
1158 eps=eps
1159 ),
1160 result_r = _assemble_a_path_from_fragments(
1161 fragments=fragments,
1162 startfrag=minxidx,
1163 rightmost=true,
1164 eps=eps
1165 ),
1166 l_area = abs(polygon_area(result_l[0])),
1167 r_area = abs(polygon_area(result_r[0])),
1168 result = l_area < r_area? result_l : result_r,
1169 newpath = cleanup_path(result[0]),
1170 remainder = result[1],
1171 finished = min(l_area,r_area)<eps ? _finished : concat(_finished, [newpath])
1172 ) _assemble_path_fragments(
1173 fragments=remainder,
1174 eps=eps,
1175 _finished=finished
1176 );
1177
1178
1179
1180// vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap