diff --git a/comparisons.scad b/comparisons.scad index bd2a25c..300b468 100644 --- a/comparisons.scad +++ b/comparisons.scad @@ -290,12 +290,12 @@ function compare_lists(a, b) = // idx = find_approx(val, list, [start=], [eps=]); // indices = find_approx(val, list, all=true, [start=], [eps=]); // Description: -// Finds the first item in `list` that matches `val`, returning the index. +// Finds the first item in `list` that matches `val`, returning the index. Returns `undef` if there is no match. // Arguments: // val = The value to search for. If given a function literal of signature `function (x)`, uses that function to check list items. Returns true for a match. // list = The list to search through. // --- -// start = The index to start searching from. +// start = The index to start searching from. Default: 0 // all = If true, returns a list of all matching item indices. // eps = The maximum allowed floating point rounding error for numeric comparisons. function find_approx(val, list, start=0, all=false, eps=EPSILON) = diff --git a/geometry.scad b/geometry.scad index c94d5ca..77930bb 100644 --- a/geometry.scad +++ b/geometry.scad @@ -104,15 +104,16 @@ function _tri_class(tri, eps=EPSILON) = /// class = _pt_in_tri(point, tri); /// Topics: Geometry, Points, Triangles /// Description: -/// Return 1 if point is inside the triangle interion. -/// Return =0 if point is on the triangle border. -/// Return -1 if point is outside the triangle. +// For CW triangles `tri` : +/// return 1 if point is inside the triangle interior. +/// return =0 if point is on the triangle border. +/// return -1 if point is outside the triangle. /// Arguments: /// point = The point to check position of. /// tri = A list of the three 2d vertices of a triangle. /// eps = Tolerance in the geometrical tests. function _pt_in_tri(point, tri, eps=EPSILON) = - min( _tri_class([tri[0],tri[1],point],eps), + min( _tri_class([tri[0],tri[1],point],eps), _tri_class([tri[1],tri[2],point],eps), _tri_class([tri[2],tri[0],point],eps) ); @@ -1701,7 +1702,7 @@ function point_in_polygon(point, poly, nonzero=false, eps=EPSILON) = // Description: // Given a simple polygon in 2D or 3D, triangulates it and returns a list // of triples indexing into the polygon vertices. When the optional argument `ind` is -// given, it is used as an index list into `poly` to define the polygon. In that case, +// given, it is used as an index list into `poly` to define the polygon vertices. In that case, // `poly` may have a length greater than `ind`. When `ind` is undefined, all points in `poly` // are considered as vertices of the polygon. // . @@ -1710,47 +1711,49 @@ function point_in_polygon(point, poly, nonzero=false, eps=EPSILON) = // vector with the same direction of the polygon normal. // . // The function produce correct triangulations for some non-twisted non-simple polygons. -// A polygon is non-twisted iff it is simple or there is a partition of it in +// A polygon is non-twisted iff it is simple or it has a partition in // simple polygons with the same winding such that the intersection of any two partitions is -// made of full edges of both partitions. These polygons may have "touching" vertices +// made of full edges and/or vertices of both partitions. These polygons may have "touching" vertices // (two vertices having the same coordinates, but distinct adjacencies) and "contact" edges // (edges whose vertex pairs have the same pairwise coordinates but are in reversed order) but has // no self-crossing. See examples bellow. If all polygon edges are contact edges (polygons with -// zero area), it returns an empty list for 2d polygons and issues an error for 3d polygons. +// zero area), it returns an empty list for 2d polygons and reports an error for 3d polygons. +// Triangulation errors are reported either by an assert error (when `error=true`) or by returning +// `undef` (when `error=false`). Invalid arguments always produce an assert error. // . -// Twisted polygons have no consistent winding and when input to this function usually produce -// an error but when an error is not issued the outputs are not correct triangulations. The function +// Twisted polygons have no consistent winding and when input to this function usually reports +// an error but when an error is not reported the outputs are not correct triangulations. The function // can work for 3d non-planar polygons if they are close enough to planar but may otherwise -// issue an error for this case. +// report an error for this case. // Arguments: // poly = Array of the polygon vertices. // ind = A list indexing the vertices of the polygon in `poly`. // eps = A maximum tolerance in geometrical tests. Default: EPSILON -// Example(2D,NoAxes): +// Example(2D,NoAxes): a simple polygon; see from above // poly = star(id=10, od=15,n=11); // tris = polygon_triangulate(poly); // color("lightblue") for(tri=tris) polygon(select(poly,tri)); // color("blue") up(1) for(tri=tris) { stroke(select(poly,tri),.15,closed=true); } // color("magenta") up(2) stroke(poly,.25,closed=true); // color("black") up(3) vnf_debug([path3d(poly),[]],faces=false,size=1); -// Example(2D,NoAxes): a polygon with a hole and one "contact" edge +// Example(2D,NoAxes): a polygon with a hole and one "contact" edge; see from above // poly = [ [-10,0], [10,0], [0,10], [-10,0], [-4,4], [4,4], [0,2], [-4,4] ]; // tris = polygon_triangulate(poly); // color("lightblue") for(tri=tris) polygon(select(poly,tri)); // color("blue") up(1) for(tri=tris) { stroke(select(poly,tri),.15,closed=true); } // color("magenta") up(2) stroke(poly,.25,closed=true); // color("black") up(3) vnf_debug([path3d(poly),[]],faces=false,size=1); -// Example(2D,NoAxes): a polygon with "touching" vertices and no holes +// Example(2D,NoAxes): a polygon with "touching" vertices and no holes; see from above // poly = [ [0,0], [5,5], [-5,5], [0,0], [-5,-5], [5,-5] ]; // tris = polygon_triangulate(poly); // color("lightblue") for(tri=tris) polygon(select(poly,tri)); // color("blue") up(1) for(tri=tris) { stroke(select(poly,tri),.15,closed=true); } // color("magenta") up(2) stroke(poly,.25,closed=true); // color("black") up(3) vnf_debug([path3d(poly),[]],faces=false,size=1); -// Example(2D,NoAxes): a polygon with "contact" edges and no holes +// Example(2D,NoAxes): a polygon with "contact" edges and no holes; see from above // poly = [ [0,0], [10,0], [10,10], [0,10], [0,0], [3,3], [7,3], // [7,7], [7,3], [3,3] ]; -// tris = polygon_triangulate(poly); // see from above +// tris = polygon_triangulate(poly); // color("lightblue") for(tri=tris) polygon(select(poly,tri)); // color("blue") up(1) for(tri=tris) { stroke(select(poly,tri),.15,closed=true); } // color("magenta") up(2) stroke(poly,.25,closed=true); @@ -1762,19 +1765,18 @@ function point_in_polygon(point, poly, nonzero=false, eps=EPSILON) = // vnf_tri = [vnf[0], [for(face=vnf[1]) each polygon_triangulate(vnf[0], face) ] ]; // color("blue") // vnf_wireframe(vnf_tri, width=.15); -function polygon_triangulate(poly, ind, eps=EPSILON) = +function polygon_triangulate(poly, ind, error=true, eps=EPSILON) = assert(is_path(poly) && len(poly)>=3, "Polygon `poly` should be a list of at least three 2d or 3d points") - assert(is_undef(ind) - || (is_vector(ind) && min(ind)>=0 && max(ind)=0 && max(ind) 2*eps, - "The polygon vertices are collinear.") - [ind] + let( degen = norm(scalar_vec3(cross(poly[ind[1]]-poly[ind[0]], poly[ind[2]]-poly[ind[0]]))) < 2*eps ) + assert( ! error || ! degen, "The polygon vertices are collinear.") + degen ? undef : [ind] : len(poly[ind[0]]) == 3 ? // find a representation of the polygon as a 2d polygon by projecting it on its own plane let( @@ -1785,44 +1787,46 @@ function polygon_triangulate(poly, ind, eps=EPSILON) = pts = select(poly,ind), nrm = -polygon_normal(pts) ) - assert( nrm!=undef, + assert( ! error || (nrm != undef), "The polygon has self-intersections or zero area or its vertices are collinear or non coplanar.") + nrm == undef ? undef : let( imax = max_index([for(p=pts) norm(p-pts[0]) ]), v1 = unit( pts[imax] - pts[0] ), v2 = cross(v1,nrm), prpts = pts*transpose([v1,v2]) // the 2d projection of pts on the polygon plane ) - [for(tri=_triangulate(prpts, count(len(ind)), eps)) select(ind,tri) ] + let( tris = _triangulate(prpts, count(len(ind)), error, eps) ) + tris == undef ? undef : + [for(tri=tris) select(ind,tri) ] : is_polygon_clockwise(select(poly, ind)) - ? _triangulate( poly, ind, eps ) - : [for(tri=_triangulate( poly, reverse(ind), eps )) reverse(tri) ]; + ? _triangulate( poly, ind, error, eps ) + : let( tris = _triangulate( poly, reverse(ind), error, eps ) ) + tris == undef ? undef : + [for(tri=tris) reverse(tri) ]; // poly is supposed to be a 2d cw polygon // implements a modified version of ear cut method for non-twisted polygons -// the polygons accepted by this function are (tecnically) the ones whose interior -// is homeomoph to the interior of a simple polygon -function _triangulate(poly, ind, eps=EPSILON, tris=[]) = +// the polygons accepted by this function are those decomposable in simple +// CW polygons. +function _triangulate(poly, ind, error, eps=EPSILON, tris=[]) = len(ind)==3 ? _degenerate_tri(select(poly,ind),eps) ? tris // if last 3 pts perform a degenerate triangle, ignore it : concat(tris,[ind]) // otherwise, include it : let( ear = _get_ear(poly,ind,eps) ) -/* -let( x= [if(is_undef(ear)) echo(ind=ind) 0] ) -is_undef(ear) ? tris : -*/ - assert( ear!=undef, + assert( ! error || (ear != undef), "The polygon has twists or all its vertices are collinear or non coplanar.") + ear == undef ? undef : is_list(ear) // is it a degenerate ear ? ? len(ind) <= 4 ? tris : - _triangulate(poly, select(ind,ear[0]+3, ear[0]), eps, tris) // discard it + _triangulate(poly, select(ind,ear[0]+3, ear[0]), error, eps, tris) // discard it : let( ear_tri = select(ind,ear,ear+2), indr = select(ind,ear+2, ear) // indices of the remaining path ) - _triangulate(poly, indr, eps, concat(tris,[ear_tri])); + _triangulate(poly, indr, error, eps, concat(tris,[ear_tri])); // a returned ear will be: @@ -1847,9 +1851,7 @@ function _get_ear(poly, ind, eps, _i=0) = // otherwise check the next ear candidate _i=0 ) // or it is equal to p1 and some of its adjacent edges cross the open segment (p0,p2) || ( norm(vert-p1) < eps - && ( _is_at_left(p0,[prev_vert,p1],eps) - && _is_at_left(p2,[p1,next_vert],eps) ) + && _is_at_left(p0,[prev_vert,p1],eps) && _is_at_left(p2,[p1,prev_vert],eps) + && _is_at_left(p2,[p1,next_vert],eps) && _is_at_left(p0,[next_vert,p1],eps) ) ) ? false - : _none_inside(idxs,poly,p0,p1,p2,eps,i=i+1); + : _none_inside(idxs,poly,p0,p1,p2,eps,i=i+1); // Function: is_polygon_clockwise() diff --git a/lists.scad b/lists.scad index 9da1972..0720545 100644 --- a/lists.scad +++ b/lists.scad @@ -230,14 +230,15 @@ function select(list, start, end) = : end==undef ? is_num(start) ? list[ (start%l+l)%l ] - : assert( is_list(start) || is_range(start), "Invalid start parameter") + : assert( start==[] || is_vector(start) || is_range(start), "Invalid start parameter") [for (i=start) list[ (i%l+l)%l ] ] : assert(is_finite(start), "When `end` is given, `start` parameter should be a number.") assert(is_finite(end), "Invalid end parameter.") let( s = (start%l+l)%l, e = (end%l+l)%l ) (s <= e) - ? [for (i = [s:1:e]) list[i]] - : concat([for (i = [s:1:l-1]) list[i]], [for (i = [0:1:e]) list[i]]) ; + ? [ for (i = [s:1:e]) list[i] ] + : [ for (i = [s:1:l-1]) list[i], + for (i = [0:1:e]) list[i] ] ; // Function: slice() diff --git a/paths.scad b/paths.scad index 15349dd..1c76251 100644 --- a/paths.scad +++ b/paths.scad @@ -278,9 +278,8 @@ function _path_self_intersections(path, closed=true, eps=EPSILON) = // signs at its two vertices can have an intersection with segment // [a1,a2]. The variable signals is zero when abs(vals[j]-ref) is less than // eps and the sign of vals[j]-ref otherwise. - signals = [for(j=[i+2:1:plen-(i==0 && closed? 2: 1)]) vals[j]-ref > eps ? 1 - : vals[j]-ref < -eps ? -1 - : 0] + signals = [for(j=[i+2:1:plen-(i==0 && closed? 2: 1)]) + abs(vals[j]-ref) < eps ? 0 : sign(vals[j]-ref) ] ) if(max(signals)>=0 && min(signals)<=0 ) // some remaining edge intersects line [a1,a2] for(j=[i+2:1:plen-(i==0 && closed? 3: 2)]) diff --git a/regions.scad b/regions.scad index 36549b3..05bdd57 100644 --- a/regions.scad +++ b/regions.scad @@ -282,9 +282,9 @@ function _region_region_intersections(region1, region2, closed1=true,closed2=tru for(p2=idx(region2)) let( poly = closed2?close_path(region2[p2]):region2[p2], - signs = [for(v=poly*seg_normal) v-ref> eps ? 1 : v-ref<-eps ? -1 : 0] + signs = [for(v=poly*seg_normal) abs(v-ref) < eps ? 0 : sign(v-ref) ] ) - if(max(signs)>=0 && min(signs)<=0) // some edge edge intersects line [a1,a2] + if(max(signs)>=0 && min(signs)<=0) // some edge intersects line [a1,a2] for(j=[0:1:len(poly)-2]) if(signs[j]!=signs[j+1]) let( // exclude non-crossing and collinear segments @@ -329,7 +329,7 @@ function _region_region_intersections(region1, region2, closed1=true,closed2=tru // where region1 intersections region2. Split region2 similarly with respect to region1. // The return is a pair of results of the form [split1, split2] where split1=[frags1,frags2,...] // and frags1 is a list of path pieces (in order) from the first path of the region. -// You can pass a single path in for either region, but the output will be a singleton list, as ify +// You can pass a single path in for either region, but the output will be a singleton list, as if // you passed in a singleton region. // Arguments: // region1 = first region diff --git a/vnf.scad b/vnf.scad index 3457a45..3351bd9 100644 --- a/vnf.scad +++ b/vnf.scad @@ -414,22 +414,22 @@ function _old_cleave_connected_region(region) = /// Internal Function: _cleave_connected_region(region, eps) /// Description: -/// Given a region that is connected and has its outer border in region[0], -/// produces a polygon with the same points that has overlapping connected paths -/// to join internal holes to the outer border. Output is a single path. -/// It expect that region[0] be a simple closed CW path and that each hole, -/// region[i] for i>0, be a simple closed CCW path. -/// The paths are also supposed to be disjoint except for common vertices and -/// common edges but no crossing. -/// This function implements an extension of the algorithm discussed in: -/// https://www.geometrictools.com/Documentation/TriangulationByEarClipping.pdf +/// Given a region that is connected and has its outer border in region[0], +/// produces a overlapping connected path to join internal holes to +/// the outer border without adding points. Output is a single non-simple polygon. +/// Requirements: +/// It expects that all region paths be simple closed paths, with region[0] CW and +/// the other paths CCW and encircled by region[0]. The input region paths are also +/// supposed to be disjoint except for common vertices and common edges but with +/// no crossings. It may return `undef` if these conditions are not met. +/// This function implements an extension of the algorithm discussed in: +/// https://www.geometrictools.com/Documentation/TriangulationByEarClipping.pdf function _cleave_connected_region(region, eps=EPSILON) = len(region)==1 ? region[0] : let( - outer = deduplicate(clockwise_polygon(region[0])), // - holes = [for(i=[1:1:len(region)-1]) // possibly unneeded - let(poly=region[i]) // - deduplicate( ccw_polygon(poly) ) ], // + outer = deduplicate(region[0]), // + holes = [for(i=[1:1:len(region)-1]) // deduplication possibly unneeded + deduplicate( region[i] ) ], // extridx = [for(li=holes) max_index(column(li,0)) ], // the right extreme vertex for each hole sorted by decreasing x values extremes = sort( [for(i=idx(holes)) [ i, extridx[i], -holes[i][extridx[i]].x] ], idx=2 ) @@ -439,7 +439,7 @@ function _cleave_connected_region(region, eps=EPSILON) = // connect the hole paths one at a time to the outer path. // 'extremes' is the list of the right extreme vertex of each hole sorted by decreasing abscissas -// see _cleave_connected_region(region, eps) +// see: _cleave_connected_region(region, eps) function _polyHoles(outer, holes, extremes, eps=EPSILON, n=0) = let( extr = extremes[n], // @@ -447,17 +447,17 @@ function _polyHoles(outer, holes, extremes, eps=EPSILON, n=0) = ipt = extr[1], // index of the hole point with maximum abscissa brdg = _bridge(hole[ipt], outer, eps) // the index of a point in outer to bridge hole[ipt] to ) - assert(brdg!=undef, "Error: check input polygon restrictions") + brdg == undef ? undef : let( l = len(outer), lh = len(hole), // the new outer polygon bridging the hole to the old outer npoly = approx(outer[brdg], hole[ipt], eps) - ? [ for(i=[brdg: 1: brdg+l]) outer[i%l] , - for(i=[ipt+1:1: ipt+lh-1]) hole[i%lh] ] - : [ for(i=[brdg: 1: brdg+l]) outer[i%l] , - for(i=[ipt:1: ipt+lh]) hole[i%lh] ] + ? [ for(i=[brdg: 1: brdg+l]) outer[i%l] , + for(i=[ipt+1: 1: ipt+lh-1]) hole[i%lh] ] + : [ for(i=[brdg: 1: brdg+l]) outer[i%l] , + for(i=[ipt: 1: ipt+lh]) hole[i%lh] ] ) n==len(holes)-1 ? npoly : _polyHoles(npoly, holes, extremes, eps, n+1); @@ -472,13 +472,13 @@ function _bridge(pt, outer,eps) = let( l = len(outer), crxs = - [for( i=idx(outer) ) - let( edge = select(outer,i,i+1) ) + let( edges = pair(outer,wrap=true) ) + [for( i = idx(edges) ) + let( edge = edges[i] ) // consider just descending outer edges at right of pt crossing ordinate pt.y - if( (edge[0].y> pt.y) - && (edge[1].y<=pt.y) - && ( norm(edge[1]-pt)0 ) ) + if( (edge[0].y > pt.y+eps) + && (edge[1].y <= pt.y) + && _is_at_left(pt, [edge[1], edge[0]], eps) ) [ i, // the point of edge with ordinate pt.y abs(pt.y-edge[1].y) error - assert( ! approx(pt,isect,eps) || approx(pt,vert0,eps) || approx(pt,vert1,eps), - "There is a forbidden self_intersection" ) - norm(pt-vert0) < eps ? proj[0] : // if pt touches an outer vertex, return its index - norm(pt-vert1) < eps ? (proj[0]+1)%l : + norm(pt-vert1) < eps ? (proj[0]+1)%l : // if pt touches an outer vertex, return its index + // as vert0.y > pt.y then pt!=vert0 + norm(pt-isect) < eps ? undef : // if pt touches the middle of an outer edge -> error let( // the edge [vert0, vert1] necessarily satisfies vert0.y > vert1.y // indices of candidates to an outer bridge point @@ -553,13 +552,15 @@ function _bridge(pt, outer,eps) = function vnf_from_region(region, transform, reverse=false) = let ( regions = region_parts(force_region(region)), - vnfs = [ - for (rgn = regions) let( - cleaved = path3d(_cleave_connected_region(rgn)), - face = is_undef(transform)? cleaved : apply(transform,cleaved), - faceidxs = reverse? [for (i=[len(face)-1:-1:0]) i] : [for (i=[0:1:len(face)-1]) i] - ) [face, [faceidxs]] - ], + vnfs = + [ for (rgn = regions) + let( cleaved = path3d(_cleave_connected_region(rgn)) ) + assert( cleaved, "The region is invalid") + let( + face = is_undef(transform)? cleaved : apply(transform,cleaved), + faceidxs = reverse? [for (i=[len(face)-1:-1:0]) i] : [for (i=[0:1:len(face)-1]) i] + ) [face, [faceidxs]] + ], outvnf = vnf_merge(vnfs) ) vnf_triangulate(outvnf); @@ -667,9 +668,13 @@ function _link_indicator(l,imin,imax) = function vnf_triangulate(vnf) = let( verts = vnf[0], - faces = [for (face=vnf[1]) each len(face)==3 ? [face] : - polygon_triangulate(verts, face)] - ) [verts, faces]; + faces = [for (face=vnf[1]) + each (len(face)==3 ? [face] : + let( tris = polygon_triangulate(verts, face) ) + assert( tris!=undef, "Some `vnf` face cannot be triangulated.") + tris ) ] + ) + [verts, faces];