grey/BOSL2
1
0
Fork 0
mirror of https://github.com/BelfrySCAD/BOSL2.git synced 2024-12-29 00:09:41 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

Correction in _cleave_connected_region and polygon_triangulate and some few tweaks

Browse source
This commit is contained in:
RonaldoCMP 2021-11-04 12:09:29 +00:00
parent 10604cd20b
commit 6bd1dd918f
6 changed files with 105 additions and 98 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

4
comparisons.scad
View file

@ -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) =

88
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) );
@ -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)<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.")
[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<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]];
@ -1873,12 +1875,12 @@ 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
: _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()

7
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()

5
paths.scad
View file

@ -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)])

6
regions.scad
View file

@ -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

93
vnf.scad
View file

@ -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)<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,21 @@ function _bridge(pt, outer,eps) =
]
]
)
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 :
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];