grey/BOSL2
1
0
Fork 0
mirror of https://github.com/BelfrySCAD/BOSL2.git synced 2025-01-01 09:49:45 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

Corrections in _cleave_connected_region and polygon_triangulate and some minor tweaks,

Browse source
This commit is contained in:
RonaldoCMP 2021-11-03 20:57:35 +00:00
parent c5da41ed8c
commit 8c4022fa3c
5 changed files with 112 additions and 86 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

86
geometry.scad
View file

@ -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) );
@ -1695,7 +1696,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.
// .
@ -1704,47 +1705,50 @@ 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`.
// error = If false, reports triangulation errors returning `undef` otherwise report them by an assert error. Default: true
// 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);
@ -1756,18 +1760,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)<len(poly) ),
assert(is_undef(ind) || (is_vector(ind) && min(ind)>=0 && max(ind)<len(poly) ),
"Improper or out of bounds list of indices")
(! is_undef(ind) ) && len(ind) == 0 ? [] :
let( ind = is_undef(ind) ? count(len(poly)) : ind )
len(ind) <=2 ? [] :
len(ind) == 3
? _degenerate_tri([poly[ind[0]], poly[ind[1]], poly[ind[2]]], eps) ? [] :
// non zero area
assert( norm(scalar_vec3(cross(poly[ind[1]]-poly[ind[0]], poly[ind[2]]-poly[ind[0]]))) > 2*eps,
"The polygon vertices are collinear.")
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
@ -1779,44 +1783,47 @@ 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=[]) =
// is homeomoph to the interior of a circle
function _triangulate(poly, ind, error, eps=EPSILON, tris=[]) =
//echo(eps=eps)
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:
@ -1825,6 +1832,7 @@ is_undef(ear) ? tris :
// 2. or a degenerate triangle where two vertices are coincident
// the returned ear is specified by the index of `ind` of its first vertex
function _get_ear(poly, ind, eps, _i=0) =
//let( x = ! is_num(eps) ? echo(ind=ind,eps=eps) 0:0 )
let( lind = len(ind) )
lind==3 ? 0 :
let( // the _i-th ear candidate
@ -1841,9 +1849,7 @@ function _get_ear(poly, ind, eps, _i=0) =
// otherwise check the next ear candidate
_i<lind-1 ? _get_ear(poly, ind, eps, _i=_i+1) :
// poly has no ears, look for wiskers
let(
wiskers = [for(j=idx(ind)) if(norm(poly[ind[j]]-poly[ind[(j+2)%lind]])<eps) j ]
)
let( wiskers = [for(j=idx(ind)) if(norm(poly[ind[j]]-poly[ind[(j+2)%lind]])<eps) j ] )
wiskers==[] ? undef : [wiskers[0]];
@ -1867,8 +1873,8 @@ function _none_inside(idxs,poly,p0,p1,p2,eps,i=0) =
_tri_class([p2,p0,vert],eps)>=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

25
lists.scad
View file

@ -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()
@ -952,7 +953,7 @@ function enumerate(l,idx=undef) =
// Example:
// l = ["A","B","C","D"];
// echo([for (p=pair(l)) str(p.y,p.x)]); // Outputs: ["BA", "CB", "DC"]
function pair(list, wrap=false) =
function _old_pair(list, wrap=false) =
assert(is_list(list)||is_string(list), "Invalid input." )
assert(is_bool(wrap))
let(
@ -961,6 +962,12 @@ function pair(list, wrap=false) =
? [for (i=[0:1:ll-1]) [list[i], list[(i+1) % ll]]]
: [for (i=[0:1:ll-2]) [list[i], list[i+1]]];
function pair(list, wrap=false) =
assert(is_list(list)||is_string(list), "Invalid input." )
assert(is_bool(wrap))
let( ll = len(list)-1 )
[for (i=[0:1:ll-1]) [list[i], list[i+1]], if(wrap) [list[ll], list[0]] ];
// Function: triplet()
// Usage:
@ -981,7 +988,7 @@ function pair(list, wrap=false) =
// translate(b) rot(from=FWD,to=v) anchor_arrow2d();
// }
// stroke(path);
function triplet(list, wrap=false) =
function _old_triplet(list, wrap=false) =
assert(is_list(list)||is_string(list), "Invalid input." )
assert(is_bool(wrap))
let(
@ -991,6 +998,14 @@ function triplet(list, wrap=false) =
: [for (i=[0:1:ll-3]) [ list[i], list[i+1], list[i+2] ]];
function triplet(list, wrap=false) =
assert(is_list(list)||is_string(list), "Invalid input." )
assert(is_bool(wrap))
let( ll = len(list) )
[ for (i=[0:1:ll-3]) [ list[i], list[i+1], list[i+2] ],
each if(wrap) [ [ list[ll-2], list[ll-1], list[0] ], [ list[ll-1], list[0], list[1] ] ] ];
// Function: combinations()
// Usage:
// list = combinations(l, [n]);

6
paths.scad
View file

@ -274,9 +274,7 @@ 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)])
@ -286,10 +284,8 @@ function _path_self_intersections(path, closed=true, eps=EPSILON) =
isect = _general_line_intersection([a1,a2],[b1,b2],eps=eps)
)
if (isect
// && isect[1]> (i==0 && !closed? -eps: 0) // Apparently too strict
&& isect[1]>=-eps
&& isect[1]<= 1+eps
// && isect[2]> 0
&& isect[2]>= -eps
&& isect[2]<= 1+eps)
[isect[0], i, isect[1], j, isect[2]]

6
regions.scad
View file

@ -216,9 +216,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
@ -260,7 +260,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

75
vnf.scad
View file

@ -414,15 +414,15 @@ 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.
/// It expect that all region paths be simple closed paths with arbitrary widings.
/// 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(
@ -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,18 @@ 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")
// 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 +473,12 @@ 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), 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)<eps // accepts touching vertices
|| _tri_class([pt, edge[0], edge[1]], eps)>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)<eps ? edge[1] :
@ -487,22 +487,25 @@ function _bridge(pt, outer,eps) =
]
]
)
assert(crxs!=[], "Error: check input polygon restrictions")
// assert(crxs!=[], "Error: check input polygon restrictions")
crxs == [] ? undef :
let(
// the intersection point nearest to pt
// the intersection point of the nearest edge to pt with minimum slope
minX = min([for(p=crxs) p[1].x]),
crxcand = [for(crx=crxs) if(crx[1].x < minX+eps) crx ],
nearest = min_index([for(crx=crxcand) outer[crx[0]].y]),
crxcand = [for(crx=crxs) if(crx[1].x < minX+eps) crx ], // nearest edges
nearest = min_index([for(crx=crxcand)
(outer[crx[0]].x - pt.x) / (outer[crx[0]].y - pt.y) ]), // minimum slope
proj = crxcand[nearest],
vert0 = outer[proj[0]], // the two vertices of the nearest crossing edge
vert1 = outer[(proj[0]+1)%l],
isect = proj[1] // the intersection point
)
// if pt touches the middle of an outer edge -> 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 :
// as vert0.y > pt.y then pt!=vert0
norm(pt-vert1) < eps ? (proj[0]+1)%l : // if pt touches an outer vertex, return its index
// assert( ! approx(pt,isect,eps) , // if pt touches the middle of an outer edge -> error
// "A hole vertex is in the middle of an edge" )
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
@ -555,7 +558,10 @@ function vnf_from_region(region, transform, reverse=false) =
regions = region_parts(force_region(region)),
vnfs = [
for (rgn = regions) let(
cleaved = path3d(_cleave_connected_region(rgn)),
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]]
@ -668,8 +674,11 @@ function vnf_triangulate(vnf) =
let(
verts = vnf[0],
faces = [for (face=vnf[1]) each len(face)==3 ? [face] :
polygon_triangulate(verts, face)]
) [verts, faces];
polygon_triangulate(verts, face)],
invalid = [for(face=faces) if(face==undef) 1 ]
)
assert( invalid==[], "Some `vnf` face cannot be triangulated.")
[verts, faces];