From 99e40a5679ab29db044aeb969a4f7910e421966f Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Wed, 7 Apr 2021 07:13:03 +0100 Subject: [PATCH 01/21] Test convexity of 3d polygons --- geometry.scad | 34 ++++++++++++++++++---------------- tests/test_geometry.scad | 1 + 2 files changed, 19 insertions(+), 16 deletions(-) diff --git a/geometry.scad b/geometry.scad index 73abe90..21b4f29 100644 --- a/geometry.scad +++ b/geometry.scad @@ -1074,6 +1074,7 @@ function distance_from_plane(plane, point) = let( plane = normalize_plane(plane) ) point3d(plane)* point - plane[3]; + // Returns [POINT, U] if line intersects plane at one point. // Returns [LINE, undef] if the line is on the plane. // Returns undef if line is parallel to, but not on the given plane. @@ -1596,7 +1597,6 @@ function circle_circle_tangents(c1,r1,c2,r2,d1,d2) = ]; - // Function: circle_line_intersection() // Usage: // isect = circle_line_intersection(c,r,line,,); @@ -1629,7 +1629,6 @@ function circle_line_intersection(c,r,line,d,bounded=false,eps=EPSILON) = let( offset = sqrt(r*r-d*d), uvec=unit(line[1]-line[0]) ) [closest-offset*uvec, closest+offset*uvec] - ) [for(p=isect) if ((!bounded[0] || (p-line[0])*(line[1]-line[0])>=0) @@ -1637,7 +1636,6 @@ function circle_line_intersection(c,r,line,d,bounded=false,eps=EPSILON) = - // Section: Pointlists @@ -1758,27 +1756,31 @@ function polygon_area(poly, signed=false) = // Usage: // is_convex_polygon(poly); // Description: -// Returns true if the given 2D polygon is convex. The result is meaningless if the polygon is not simple (self-intersecting). -// If the points are collinear the result is true. +// Returns true if the given 2D or 3D polygon is convex. +// The result is meaningless if the polygon is not simple (self-intersecting) or non coplanar. +// If the points are collinear an error is generated. // Arguments: // poly = Polygon to check. // Example: // is_convex_polygon(circle(d=50)); // Returns: true +// is_convex_polygon(rot([50,120,30], p=path3d(circle(1,$fn=50)))); // Returns: true // Example: // spiral = [for (i=[0:36]) let(a=-i*10) (10+i)*[cos(a),sin(a)]]; // is_convex_polygon(spiral); // Returns: false function is_convex_polygon(poly) = - assert(is_path(poly,dim=2), "The input should be a 2D polygon." ) - let( l = len(poly) ) - len([for( i = l-1, - c = cross(poly[(i+1)%l]-poly[i], poly[(i+2)%l]-poly[(i+1)%l]), - s = sign(c); - i>=0 && sign(c)==s; - i = i-1, - c = i<0? 0: cross(poly[(i+1)%l]-poly[i],poly[(i+2)%l]-poly[(i+1)%l]), - s = s==0 ? sign(c) : s - ) i - ])== l; + assert(is_path(poly), "The input should be a 2D or 3D polygon." ) + let( lp = len(poly), + p0 = poly[0] ) + assert( lp>=3 , "A polygon must have at least 3 points" ) + let( crosses = [for(i=[0:1:lp-1]) cross(poly[(i+1)%lp]-poly[i], poly[(i+2)%lp]-poly[(i+1)%lp]) ] ) + len(p0)==2 + ? assert( !approx(max(crosses)) && !approx(min(crosses)), "The points are collinear" ) + min(crosses) >=0 || max(crosses)<=0 + : let( prod = crosses*sum(crosses), + minc = min(prod), + maxc = max(prod) ) + assert( !approx(maxc-minc), "The points are collinear" ) + minc>=0 || maxc<=0; // Function: polygon_shift() diff --git a/tests/test_geometry.scad b/tests/test_geometry.scad index 0e70f31..492deb1 100644 --- a/tests/test_geometry.scad +++ b/tests/test_geometry.scad @@ -844,6 +844,7 @@ module test_polygon_area() { module test_is_convex_polygon() { assert(is_convex_polygon([[1,1],[-1,1],[-1,-1],[1,-1]])); assert(is_convex_polygon(circle(r=50,$fn=1000))); + assert(is_convex_polygon(rot([50,120,30], p=path3d(circle(1,$fn=50))))); assert(!is_convex_polygon([[1,1],[0,0],[-1,1],[-1,-1],[1,-1]])); } *test_is_convex_polygon(); From ca29e2ef8707b9c419b31dcf7be2d4418a9c66d2 Mon Sep 17 00:00:00 2001 From: github-actions Date: Wed, 7 Apr 2021 07:27:34 +0000 Subject: [PATCH 02/21] Bump release version. --- version.scad | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/version.scad b/version.scad index 89b3852..457f177 100644 --- a/version.scad +++ b/version.scad @@ -6,7 +6,7 @@ ////////////////////////////////////////////////////////////////////// -BOSL_VERSION = [2,0,601]; +BOSL_VERSION = [2,0,602]; // Section: BOSL Library Version Functions From ed2daac8c2a48ec1247e510d731540d79238717a Mon Sep 17 00:00:00 2001 From: Adrian Mariano Date: Wed, 7 Apr 2021 16:20:54 -0400 Subject: [PATCH 03/21] Fix problem with face vertex order. --- polyhedra.scad | 44 ++++++++++++++++++++++++++++++-------------- 1 file changed, 30 insertions(+), 14 deletions(-) diff --git a/polyhedra.scad b/polyhedra.scad index ac0d575..5a30172 100644 --- a/polyhedra.scad +++ b/polyhedra.scad @@ -147,7 +147,7 @@ function _unique_groups(m) = [ // Side Effects: // `$faceindex` - Index number of the face // `$face` - Coordinates of the face (2d if rotate_children==true, 3d if not) -// `$center` - Polyhedron center in the child coordinate system +// `$center` - Face center in the child coordinate system // // Examples: All of the available polyhedra by name in their native orientation // regular_polyhedron("tetrahedron", facedown=false); @@ -275,6 +275,10 @@ function _unique_groups(m) = [ // %sphere(r=.98); // regular_polyhedron("pentagonal hexecontahedron", or=1,facedown=false); // } +// Example: Stellate an Archimedian solid, which has mixed faces +// regular_polyhedron("truncated icosahedron",stellate=1.5,or=1); +// Example: Stellate a Catalan solid where faces are not regular +// regular_polyhedron("triakis tetrahedron",stellate=0.5,or=1); module regular_polyhedron( name=undef, index=undef, @@ -330,20 +334,19 @@ module regular_polyhedron( } } } + translate(translation) if ($children>0) { maxrange = repeat ? len(faces)-1 : $children-1; for(i=[0:1:maxrange]) { // Would like to orient so an edge (longest edge?) is parallel to x axis - facepts = move(translation, p=select(scaled_points, faces[i])); - center = mean(facepts); - rotatedface = rot(from=face_normals[i], to=[0,0,1], p=move(-center, p=facepts)); - clockwise = sortidx([for(pt=rotatedface) -atan2(pt.y,pt.x)]); - $face = rotate_children? - path2d(select(rotatedface,clockwise)) : - select(move(-center,p=facepts), clockwise); + facepts = select(scaled_points, faces[i]); + $center = -mean(facepts); + cfacepts = move($center, p=facepts); + $face = rotate_children + ? path2d(rot(from=face_normals[i], to=[0,0,1], p=cfacepts)) + : cfacepts; $faceindex = i; - $center = -translation-center; - translate(center) + translate(-$center) if (rotate_children) { rot(from=[0,0,1], to=face_normals[i]) children(i % $children); @@ -699,7 +702,7 @@ function regular_polyhedron_info( face_normals, radius_scale*entry[in_radius] ] : - info == "vnf" ? [move(translation,p=scaled_points), stellate ? faces : face_triangles] : + info == "vnf" ? [move(translation,p=scaled_points), faces] : info == "vertices" ? move(translation,p=scaled_points) : info == "faces" ? faces : info == "face normals" ? face_normals : @@ -770,15 +773,28 @@ function _facenormal(pts, face) = unit(cross(pts[face[2]]-pts[face[0]], pts[face // hull() function returns triangulated faces. This function identifies the vertices that belong to each face // by grouping together the face triangles that share normal vectors. The output gives the face polygon -// point indices in arbitrary order (not usable as input to a polygon call) and a normal vector. +// point indices in arbitrary order (not usable as input to a polygon call) and a normal vector. Finally +// the faces are ordered based on angle with their center (will always give a valid order for convex polygons). +// Final return is [ordered_faces, facenormals] where the first is a list of indices into the point list +// and the second is a list of vectors. function _full_faces(pts,faces) = let( normals = [for(face=faces) quant(_facenormal(pts,face),1e-12)], groups = _unique_groups(normals), faces = [for(entry=groups) unique(flatten(select(faces, entry)))], - facenormals = [for(entry=groups) normals[entry[0]]] - ) [faces, facenormals]; + facenormals = [for(entry=groups) normals[entry[0]]], + ordered_faces = [ + for(i=idx(faces)) + let( + facepts = select(pts, faces[i]), + center = mean(facepts), + rotatedface = rot(from=facenormals[i], to=[0,0,1], p=move(-center, p=facepts)), + clockwise = sortidx([for(pt=rotatedface) -atan2(pt.y,pt.x)]) + ) + select(faces[i],clockwise) + ] + ) [ordered_faces, facenormals]; // vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap From 218246e080d5abc4764dc4e61c731a5578e4b700 Mon Sep 17 00:00:00 2001 From: Adrian Mariano Date: Wed, 7 Apr 2021 22:54:46 -0400 Subject: [PATCH 04/21] added small triambic icosahedron --- polyhedra.scad | 3 +++ 1 file changed, 3 insertions(+) diff --git a/polyhedra.scad b/polyhedra.scad index 5a30172..0914a60 100644 --- a/polyhedra.scad +++ b/polyhedra.scad @@ -118,6 +118,7 @@ function _unique_groups(m) = [ // * `"great dodecahedron"` // * `"small stellated dodecahedron"` // * `"great stellated dodecahedron"` +// * `"small triambic icosahedron"` // // Arguments: // name = Name of polyhedron to create. @@ -185,6 +186,7 @@ function _unique_groups(m) = [ // regular_polyhedron("great dodecahedron"); // regular_polyhedron("small stellated dodecahedron"); // regular_polyhedron("great stellated dodecahedron"); +// regular_polyhedron("small triambic icosahedron"); // Example: Third Archimedean solid // regular_polyhedron(type="archimedean", index=2); // Example(Med): Solids that have 8 or 10 vertex faces @@ -540,6 +542,7 @@ _stellated_polyhedra_ = [ ["great dodecahedron", "icosahedron", -sqrt(5/3-PHI)], ["small stellated dodecahedron", "dodecahedron", sqrt((5+2*sqrt(5))/5)], ["great stellated dodecahedron", "icosahedron", sqrt(2/3+PHI)], + ["small triambic icosahedron", "icosahedron", sqrt(3/5) - 1/sqrt(3)] ]; From 52ddecbafe80c6c352684fa648e4149b16c67420 Mon Sep 17 00:00:00 2001 From: github-actions Date: Thu, 8 Apr 2021 03:02:23 +0000 Subject: [PATCH 05/21] Bump release version. --- version.scad | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/version.scad b/version.scad index 457f177..7b6c730 100644 --- a/version.scad +++ b/version.scad @@ -6,7 +6,7 @@ ////////////////////////////////////////////////////////////////////// -BOSL_VERSION = [2,0,602]; +BOSL_VERSION = [2,0,603]; // Section: BOSL Library Version Functions From bdf30a56233d7193a7bac10e0a0bcda3d1441ca1 Mon Sep 17 00:00:00 2001 From: github-actions Date: Fri, 9 Apr 2021 02:55:52 +0000 Subject: [PATCH 06/21] Bump release version. --- version.scad | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/version.scad b/version.scad index 7b6c730..efc6eff 100644 --- a/version.scad +++ b/version.scad @@ -6,7 +6,7 @@ ////////////////////////////////////////////////////////////////////// -BOSL_VERSION = [2,0,603]; +BOSL_VERSION = [2,0,604]; // Section: BOSL Library Version Functions From 1451b827ecdffdf07da69a7052a96c0c2328914d Mon Sep 17 00:00:00 2001 From: Adrian Mariano Date: Fri, 9 Apr 2021 21:43:41 -0400 Subject: [PATCH 07/21] Moved vnf_tri_array to vnf.scad and added docs --- vnf.scad | 81 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 81 insertions(+) diff --git a/vnf.scad b/vnf.scad index 9c05e0f..55095c8 100644 --- a/vnf.scad +++ b/vnf.scad @@ -367,6 +367,87 @@ function vnf_vertex_array( ]); +// Function: vnf_tri_array() +// Usage: +// vnf = vnf_tri_array(points, , ) +// Description: +// Produces a vnf from an array of points where each row length can differ from the adjacent rows by up to 2 in length. This enables +// the construction of triangular VNF patches. The resulting VNF can be wrapped along the rows by setting `row_wrap` to true. +// Arguments: +// points = List of point lists for each row +// row_wrap = If true then add faces connecting the first row and last row. These rows must differ by at most 2 in length. +// reverse = Set this to reverse the direction of the faces +// Examples: Each row has one more point than the preceeding one. +// pts = [for(y=[1:1:10]) [for(x=[0:y-1]) [x,y,y]]]; +// vnf = vnf_tri_array(pts); +// vnf_wireframe(vnf,d=.1); +// color("red")move_copies(flatten(pts)) sphere(r=.15,$fn=9); +// Examples: Each row has one more point than the preceeding one. +// pts = [for(y=[0:2:10]) [for(x=[-y/2:y/2]) [x,y,y]]]; +// vnf = vnf_tri_array(pts); +// vnf_wireframe(vnf,d=.1); +// color("red")move_copies(flatten(pts)) sphere(r=.15,$fn=9); +// Example: Chaining two VNFs to construct a cone with one point length change between rows. +// pts1 = [for(z=[0:10]) path3d(arc(3+z,r=z/2+1, angle=[0,180]),10-z)]; +// pts2 = [for(z=[0:10]) path3d(arc(3+z,r=z/2+1, angle=[180,360]),10-z)]; +// vnf = vnf_tri_array(pts1, +// vnf=vnf_tri_array(pts2)); +// color("green")vnf_wireframe(vnf,d=.1); +// vnf_polyhedron(vnf); +// Example: Cone with length change two between rows +// pts1 = [for(z=[0:1:10]) path3d(arc(3+2*z,r=z/2+1, angle=[0,180]),10-z)]; +// pts2 = [for(z=[0:1:10]) path3d(arc(3+2*z,r=z/2+1, angle=[180,360]),10-z)]; +// vnf = vnf_tri_array(pts1, +// vnf=vnf_tri_array(pts2)); +// color("green")vnf_wireframe(vnf,d=.1); +// vnf_polyhedron(vnf); +// Example: Point count can change irregularly +// lens = [10,9,7,5,6,8,8,10]; +// pts = [for(y=idx(lens)) lerpn([-lens[y],y,y],[lens[y],y,y],lens[y])]; +// vnf = vnf_tri_array(pts); +// vnf_wireframe(vnf,d=.1); +// color("red")move_copies(flatten(pts)) sphere(r=.15,$fn=9); +function vnf_tri_array(points, row_wrap=false, reverse=false, vnf=EMPTY_VNF) = + let( + lens = [for(row=points) len(row)], + rowstarts = [0,each cumsum(lens)], + faces = + [for(i=[0:1:len(points) - 1 - (row_wrap ? 0 : 1)]) each + let( + rowstart = rowstarts[i], + nextrow = select(rowstarts,i+1), + delta = select(lens,i+1)-lens[i] + ) + delta == 0 ? + [for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow] : [j+rowstart, j+rowstart+1, j+nextrow], + for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+nextrow, j+nextrow+1] : [j+rowstart+1, j+nextrow+1, j+nextrow]] : + delta == 1 ? + [for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+1] : [j+rowstart, j+rowstart+1, j+nextrow+1], + for(j=[0:1:lens[i]-1]) reverse ? [j+rowstart, j+nextrow, j+nextrow+1] : [j+rowstart, j+nextrow+1, j+nextrow]] : + delta == -1 ? + [for(j=[0:1:lens[i]-3]) reverse ? [j+rowstart+1, j+nextrow, j+nextrow+1]: [j+rowstart+1, j+nextrow+1, j+nextrow], + for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow] : [j+rowstart, j+rowstart+1, j+nextrow]] : + let(count = floor((lens[i]-1)/2)) + delta == 2 ? + [ + for(j=[0:1:count-1]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+1] : [j+rowstart, j+rowstart+1, j+nextrow+1], // top triangles left + for(j=[count:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+2] : [j+rowstart, j+rowstart+1, j+nextrow+2], // top triangles right + for(j=[0:1:count]) reverse ? [j+rowstart, j+nextrow, j+nextrow+1] : [j+rowstart, j+nextrow+1, j+nextrow], // bot triangles left + for(j=[count+1:1:select(lens,i+1)-2]) reverse ? [j+rowstart-1, j+nextrow, j+nextrow+1] : [j+rowstart-1, j+nextrow+1, j+nextrow], // bot triangles right + ] : + delta == -2 ? + [ + for(j=[0:1:count-2]) reverse ? [j+nextrow, j+nextrow+1, j+rowstart+1] : [j+nextrow, j+rowstart+1, j+nextrow+1], + for(j=[count-1:1:lens[i]-4]) reverse ? [j+nextrow,j+nextrow+1,j+rowstart+2] : [j+nextrow,j+rowstart+2, j+nextrow+1], + for(j=[0:1:count-1]) reverse ? [j+nextrow, j+rowstart+1, j+rowstart] : [j+nextrow, j+rowstart, j+rowstart+1], + for(j=[count:1:select(lens,i+1)]) reverse ? [ j+nextrow-1, j+rowstart+1, j+rowstart]: [ j+nextrow-1, j+rowstart, j+rowstart+1], + ] : + assert(false,str("Unsupported row length difference of ",delta, " between row ",i," and ",(i+1)%len(points))) + ]) + vnf_merge(cleanup=true, [vnf, [flatten(points), faces]]); + + + // Module: vnf_polyhedron() // Usage: // vnf_polyhedron(vnf); From 49a3a166eb153105c8c667b4b356e3904f88d237 Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sat, 10 Apr 2021 20:09:03 +0100 Subject: [PATCH 08/21] Review of geometry.scad for speed --- arrays.scad | 16 +++--- common.scad | 11 ++-- geometry.scad | 112 ++++++++++++++++++++++++++------------- tests/test_geometry.scad | 8 +-- 4 files changed, 95 insertions(+), 52 deletions(-) diff --git a/arrays.scad b/arrays.scad index ac97444..604186b 100644 --- a/arrays.scad +++ b/arrays.scad @@ -41,6 +41,7 @@ function is_homogeneous(l, depth=10) = [] == [for(i=[1:len(l)-1]) if( ! _same_type(l[i],l0, depth+1) ) 0 ]; function is_homogenous(l, depth=10) = is_homogeneous(l, depth); + function _same_type(a,b, depth) = (depth==0) || @@ -50,7 +51,7 @@ function _same_type(a,b, depth) = (is_string(a) && is_string(b)) || (is_list(a) && is_list(b) && len(a)==len(b) && []==[for(i=idx(a)) if( ! _same_type(a[i],b[i],depth-1) ) 0] ); - + // Function: select() // Topics: List Handling @@ -97,7 +98,6 @@ function select(list, start, end) = // Function: slice() -// Topics: List Handling // Usage: // list = slice(list,s,e); // Description: @@ -476,7 +476,7 @@ function reverse(x) = // l9 = list_rotate([1,2,3,4,5],6); // Returns: [2,3,4,5,1] function list_rotate(list,n=1) = assert(is_list(list)||is_string(list), "Invalid list or string.") - assert(is_finite(n), "Invalid number") + assert(is_int(n), "The rotation number should be integer") let ( ll = len(list), n = ((n % ll) + ll) % ll, @@ -990,7 +990,7 @@ function _sort_vectors(arr, idxlist, _i=0) = _sort_vectors(equal, idxlist, _i+1), _sort_vectors(greater, idxlist, _i ) ); - + // sorting using compare_vals(); returns indexed list when `indexed==true` function _sort_general(arr, idx=undef, indexed=false) = (len(arr)<=1) ? arr : @@ -1332,6 +1332,8 @@ function permutations(l,n=2) = // pairs = zip(a,b); // triples = zip(a,b,c); // quads = zip([LIST1,LIST2,LIST3,LIST4]); +// Topics: List Handling, Iteration +// See Also: zip_long() // Description: // Zips together two or more lists into a single list. For example, if you have two // lists [3,4,5], and [8,7,6], and zip them together, you get [[3,8],[4,7],[5,6]]. @@ -1357,6 +1359,8 @@ function zip(a,b,c) = // pairs = zip_long(a,b); // triples = zip_long(a,b,c); // quads = zip_long([LIST1,LIST2,LIST3,LIST4]); +// Topics: List Handling, Iteration +// See Also: zip() // Description: // Zips together two or more lists into a single list. For example, if you have two // lists [3,4,5], and [8,7,6], and zip them together, you get [[3,8],[4,7],[5,6]]. @@ -1526,7 +1530,6 @@ function subindex(M, idx) = // [[4,2], 91, false], // [6, [3,4], undef]]; // submatrix(A,[0,2],[1,2]); // Returns [[17, "test"], [[3, 4], undef]] - function submatrix(M,idx1,idx2) = [for(i=idx1) [for(j=idx2) M[i][j] ] ]; @@ -1629,7 +1632,6 @@ function block_matrix(M) = assert(badrows==[], "Inconsistent or invalid input") bigM; - // Function: diagonal_matrix() // Usage: // mat = diagonal_matrix(diag, ); @@ -1855,7 +1857,7 @@ function transpose(arr, reverse=false) = // A = matrix to test // eps = epsilon for comparing equality. Default: 1e-12 function is_matrix_symmetric(A,eps=1e-12) = - approx(A,transpose(A)); + approx(A,transpose(A), eps); // vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap diff --git a/common.scad b/common.scad index e24c744..6d2d44f 100644 --- a/common.scad +++ b/common.scad @@ -205,7 +205,8 @@ function is_func(x) = version_num()>20210000 && is_function(x); // Description: // Tests whether input is a list of entries which all have the same list structure // and are filled with finite numerical data. You can optionally specify a required -// list structure with the pattern argument. It returns `true` for the empty list. +// list structure with the pattern argument. +// It returns `true` for the empty list regardless the value of the `pattern`. // Arguments: // list = list to check // pattern = optional pattern required to match @@ -293,7 +294,7 @@ function default(v,dflt=undef) = is_undef(v)? dflt : v; // v = The list whose items are being checked. // recursive = If true, sublists are checked recursively for defined values. The first sublist that has a defined item is returned. // Examples: -// val = first_defined([undef,7,undef,true]); // Returns: 1 +// val = first_defined([undef,7,undef,true]); // Returns: 7 function first_defined(v,recursive=false,_i=0) = _i, ); // Description: // Given a list of 3 or more coplanar 3D points, returns the coefficients of the normalized cartesian equation of a plane, // that is [A,B,C,D] where Ax+By+Cz=D is the equation of the plane where norm([A,B,C])=1. -// If `fast` is false and the points in the list are collinear or not coplanar, then `undef` is returned. -// if `fast` is true, then the coplanarity test is skipped and a plane passing through 3 non-collinear arbitrary points is returned. +// If the points in the list are collinear or not coplanar, then `undef` is returned. // Arguments: // points = The list of points to find the plane of. -// fast = If true, don't verify that all points in the list are coplanar. Default: false // eps = Tolerance in geometric comparisons. Default: `EPSILON` (1e-9) // Example(3D): // xyzpath = rot(45, v=[-0.3,1,0], p=path3d(star(n=6,id=70,d=100), 70)); @@ -911,20 +947,21 @@ function plane_from_normal(normal, pt=[0,0,0]) = // #stroke(xyzpath,closed=true); // cp = centroid(xyzpath); // move(cp) rot(from=UP,to=plane_normal(plane)) anchor_arrow(); -function plane_from_points(points, fast=false, eps=EPSILON) = +function plane_from_points(points,fast=false, eps=EPSILON) = assert( is_path(points,dim=3), "Improper 3d point list." ) assert( is_finite(eps) && (eps>=0), "The tolerance should be a non-negative value." ) - let( - indices = noncollinear_triple(points,error=false) - ) - indices==[] ? undef : - let( - p1 = points[indices[0]], - p2 = points[indices[1]], - p3 = points[indices[2]], - plane = plane3pt(p1,p2,p3) - ) - fast || points_on_plane(points,plane,eps=eps) ? plane : undef; + len(points) == 3 + ? let( plane = plane3pt(points[0],points[1],points[2]) ) + plane==[] ? undef : plane + : let( + cov_evals = _covariance_evals(points), + pm = cov_evals[0], + evals = cov_evals[1], + M = cov_evals[2], + evec = _eigenvec_symm_3(M,evals,i=2) ) +// echo(error_points_plane= abs(max(points*evec)-pm*evec), limit=eps) + !fast && abs(max(points*evec)-pm*evec)>eps*evals[0] ? undef : + [ each evec, pm*evec] ; // Function: plane_from_polygon() @@ -945,17 +982,16 @@ function plane_from_points(points, fast=false, eps=EPSILON) = // #stroke(xyzpath,closed=true); // cp = centroid(xyzpath); // move(cp) rot(from=UP,to=plane_normal(plane)) anchor_arrow(); -function plane_from_polygon(poly, fast=false, eps=EPSILON) = +function plane_from_polygon(poly, fast=false, eps=EPSILON) = assert( is_path(poly,dim=3), "Invalid polygon." ) assert( is_finite(eps) && (eps>=0), "The tolerance should be a non-negative value." ) - let( - poly = deduplicate(poly), - n = polygon_normal(poly), - plane = [n.x, n.y, n.z, n*poly[0]] - ) - fast? plane: coplanar(poly,eps=eps)? plane: []; - - + len(poly)==3 ? plane3pt(poly[0],poly[1],poly[2]) : + let( triple = sort(noncollinear_triple(poly,error=false)) ) + triple==[] ? [] : + let( plane = plane3pt(poly[triple[0]],poly[triple[1]],poly[triple[2]])) + fast? plane: points_on_plane(poly, plane, eps=eps)? plane: []; + + // Function: plane_normal() // Usage: // plane_normal(plane); @@ -1252,9 +1288,11 @@ function coplanar(points, eps=EPSILON) = len(points)<=2 ? false : let( ip = noncollinear_triple(points,error=false,eps=eps) ) ip == [] ? false : - let( plane = plane3pt(points[ip[0]],points[ip[1]],points[ip[2]]), - normal = point3d(plane) ) - max( points*normal ) - plane[3]< eps*norm(normal); + let( + plane = plane3pt(points[ip[0]],points[ip[1]],points[ip[2]]), + normal = point3d(plane), + pt_nrm = points*normal ) + abs(max(max(pt_nrm)-plane[3], -min(pt_nrm)+plane[3])) < eps; // Function: points_on_plane() @@ -1665,11 +1703,11 @@ function noncollinear_triple(points,error=true,eps=EPSILON) = n = (pb-pa)/nrm, distlist = [for(i=[0:len(points)-1]) _dist2line(points[i]-pa, n)] ) - max(distlist) Date: Sat, 10 Apr 2021 20:14:40 +0100 Subject: [PATCH 09/21] Revert "Review of geometry.scad for speed" This reverts commit 49a3a166eb153105c8c667b4b356e3904f88d237. --- arrays.scad | 16 +++--- common.scad | 11 ++-- geometry.scad | 112 +++++++++++++-------------------------- tests/test_geometry.scad | 8 ++- 4 files changed, 52 insertions(+), 95 deletions(-) diff --git a/arrays.scad b/arrays.scad index 604186b..ac97444 100644 --- a/arrays.scad +++ b/arrays.scad @@ -41,7 +41,6 @@ function is_homogeneous(l, depth=10) = [] == [for(i=[1:len(l)-1]) if( ! _same_type(l[i],l0, depth+1) ) 0 ]; function is_homogenous(l, depth=10) = is_homogeneous(l, depth); - function _same_type(a,b, depth) = (depth==0) || @@ -51,7 +50,7 @@ function _same_type(a,b, depth) = (is_string(a) && is_string(b)) || (is_list(a) && is_list(b) && len(a)==len(b) && []==[for(i=idx(a)) if( ! _same_type(a[i],b[i],depth-1) ) 0] ); - + // Function: select() // Topics: List Handling @@ -98,6 +97,7 @@ function select(list, start, end) = // Function: slice() +// Topics: List Handling // Usage: // list = slice(list,s,e); // Description: @@ -476,7 +476,7 @@ function reverse(x) = // l9 = list_rotate([1,2,3,4,5],6); // Returns: [2,3,4,5,1] function list_rotate(list,n=1) = assert(is_list(list)||is_string(list), "Invalid list or string.") - assert(is_int(n), "The rotation number should be integer") + assert(is_finite(n), "Invalid number") let ( ll = len(list), n = ((n % ll) + ll) % ll, @@ -990,7 +990,7 @@ function _sort_vectors(arr, idxlist, _i=0) = _sort_vectors(equal, idxlist, _i+1), _sort_vectors(greater, idxlist, _i ) ); - + // sorting using compare_vals(); returns indexed list when `indexed==true` function _sort_general(arr, idx=undef, indexed=false) = (len(arr)<=1) ? arr : @@ -1332,8 +1332,6 @@ function permutations(l,n=2) = // pairs = zip(a,b); // triples = zip(a,b,c); // quads = zip([LIST1,LIST2,LIST3,LIST4]); -// Topics: List Handling, Iteration -// See Also: zip_long() // Description: // Zips together two or more lists into a single list. For example, if you have two // lists [3,4,5], and [8,7,6], and zip them together, you get [[3,8],[4,7],[5,6]]. @@ -1359,8 +1357,6 @@ function zip(a,b,c) = // pairs = zip_long(a,b); // triples = zip_long(a,b,c); // quads = zip_long([LIST1,LIST2,LIST3,LIST4]); -// Topics: List Handling, Iteration -// See Also: zip() // Description: // Zips together two or more lists into a single list. For example, if you have two // lists [3,4,5], and [8,7,6], and zip them together, you get [[3,8],[4,7],[5,6]]. @@ -1530,6 +1526,7 @@ function subindex(M, idx) = // [[4,2], 91, false], // [6, [3,4], undef]]; // submatrix(A,[0,2],[1,2]); // Returns [[17, "test"], [[3, 4], undef]] + function submatrix(M,idx1,idx2) = [for(i=idx1) [for(j=idx2) M[i][j] ] ]; @@ -1632,6 +1629,7 @@ function block_matrix(M) = assert(badrows==[], "Inconsistent or invalid input") bigM; + // Function: diagonal_matrix() // Usage: // mat = diagonal_matrix(diag, ); @@ -1857,7 +1855,7 @@ function transpose(arr, reverse=false) = // A = matrix to test // eps = epsilon for comparing equality. Default: 1e-12 function is_matrix_symmetric(A,eps=1e-12) = - approx(A,transpose(A), eps); + approx(A,transpose(A)); // vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap diff --git a/common.scad b/common.scad index 6d2d44f..e24c744 100644 --- a/common.scad +++ b/common.scad @@ -205,8 +205,7 @@ function is_func(x) = version_num()>20210000 && is_function(x); // Description: // Tests whether input is a list of entries which all have the same list structure // and are filled with finite numerical data. You can optionally specify a required -// list structure with the pattern argument. -// It returns `true` for the empty list regardless the value of the `pattern`. +// list structure with the pattern argument. It returns `true` for the empty list. // Arguments: // list = list to check // pattern = optional pattern required to match @@ -294,7 +293,7 @@ function default(v,dflt=undef) = is_undef(v)? dflt : v; // v = The list whose items are being checked. // recursive = If true, sublists are checked recursively for defined values. The first sublist that has a defined item is returned. // Examples: -// val = first_defined([undef,7,undef,true]); // Returns: 7 +// val = first_defined([undef,7,undef,true]); // Returns: 1 function first_defined(v,recursive=false,_i=0) = _i, ); // Description: // Given a list of 3 or more coplanar 3D points, returns the coefficients of the normalized cartesian equation of a plane, // that is [A,B,C,D] where Ax+By+Cz=D is the equation of the plane where norm([A,B,C])=1. -// If the points in the list are collinear or not coplanar, then `undef` is returned. +// If `fast` is false and the points in the list are collinear or not coplanar, then `undef` is returned. +// if `fast` is true, then the coplanarity test is skipped and a plane passing through 3 non-collinear arbitrary points is returned. // Arguments: // points = The list of points to find the plane of. +// fast = If true, don't verify that all points in the list are coplanar. Default: false // eps = Tolerance in geometric comparisons. Default: `EPSILON` (1e-9) // Example(3D): // xyzpath = rot(45, v=[-0.3,1,0], p=path3d(star(n=6,id=70,d=100), 70)); @@ -947,21 +911,20 @@ function _covariance_evals(points) = // #stroke(xyzpath,closed=true); // cp = centroid(xyzpath); // move(cp) rot(from=UP,to=plane_normal(plane)) anchor_arrow(); -function plane_from_points(points,fast=false, eps=EPSILON) = +function plane_from_points(points, fast=false, eps=EPSILON) = assert( is_path(points,dim=3), "Improper 3d point list." ) assert( is_finite(eps) && (eps>=0), "The tolerance should be a non-negative value." ) - len(points) == 3 - ? let( plane = plane3pt(points[0],points[1],points[2]) ) - plane==[] ? undef : plane - : let( - cov_evals = _covariance_evals(points), - pm = cov_evals[0], - evals = cov_evals[1], - M = cov_evals[2], - evec = _eigenvec_symm_3(M,evals,i=2) ) -// echo(error_points_plane= abs(max(points*evec)-pm*evec), limit=eps) - !fast && abs(max(points*evec)-pm*evec)>eps*evals[0] ? undef : - [ each evec, pm*evec] ; + let( + indices = noncollinear_triple(points,error=false) + ) + indices==[] ? undef : + let( + p1 = points[indices[0]], + p2 = points[indices[1]], + p3 = points[indices[2]], + plane = plane3pt(p1,p2,p3) + ) + fast || points_on_plane(points,plane,eps=eps) ? plane : undef; // Function: plane_from_polygon() @@ -982,16 +945,17 @@ function plane_from_points(points,fast=false, eps=EPSILON) = // #stroke(xyzpath,closed=true); // cp = centroid(xyzpath); // move(cp) rot(from=UP,to=plane_normal(plane)) anchor_arrow(); -function plane_from_polygon(poly, fast=false, eps=EPSILON) = +function plane_from_polygon(poly, fast=false, eps=EPSILON) = assert( is_path(poly,dim=3), "Invalid polygon." ) assert( is_finite(eps) && (eps>=0), "The tolerance should be a non-negative value." ) - len(poly)==3 ? plane3pt(poly[0],poly[1],poly[2]) : - let( triple = sort(noncollinear_triple(poly,error=false)) ) - triple==[] ? [] : - let( plane = plane3pt(poly[triple[0]],poly[triple[1]],poly[triple[2]])) - fast? plane: points_on_plane(poly, plane, eps=eps)? plane: []; - - + let( + poly = deduplicate(poly), + n = polygon_normal(poly), + plane = [n.x, n.y, n.z, n*poly[0]] + ) + fast? plane: coplanar(poly,eps=eps)? plane: []; + + // Function: plane_normal() // Usage: // plane_normal(plane); @@ -1288,11 +1252,9 @@ function coplanar(points, eps=EPSILON) = len(points)<=2 ? false : let( ip = noncollinear_triple(points,error=false,eps=eps) ) ip == [] ? false : - let( - plane = plane3pt(points[ip[0]],points[ip[1]],points[ip[2]]), - normal = point3d(plane), - pt_nrm = points*normal ) - abs(max(max(pt_nrm)-plane[3], -min(pt_nrm)+plane[3])) < eps; + let( plane = plane3pt(points[ip[0]],points[ip[1]],points[ip[2]]), + normal = point3d(plane) ) + max( points*normal ) - plane[3]< eps*norm(normal); // Function: points_on_plane() @@ -1703,11 +1665,11 @@ function noncollinear_triple(points,error=true,eps=EPSILON) = n = (pb-pa)/nrm, distlist = [for(i=[0:len(points)-1]) _dist2line(points[i]-pa, n)] ) - max(distlist) Date: Sat, 10 Apr 2021 21:07:23 +0100 Subject: [PATCH 10/21] Review of geometry.scad for speed --- geometry.scad | 107 +++++++++++++++++++++++++++------------ tests/test_geometry.scad | 49 +++++++++--------- 2 files changed, 99 insertions(+), 57 deletions(-) diff --git a/geometry.scad b/geometry.scad index 21b4f29..80f6756 100644 --- a/geometry.scad +++ b/geometry.scad @@ -888,11 +888,49 @@ function plane3pt_indexed(points, i1, i2, i3) = // Example: // plane_from_normal([0,0,1], [2,2,2]); // Returns the xy plane passing through the point (2,2,2) function plane_from_normal(normal, pt=[0,0,0]) = - assert( is_matrix([normal,pt],2,3) && !approx(norm(normal),0), - "Inputs `normal` and `pt` should 3d vectors/points and `normal` cannot be zero." ) - concat(normal, normal*pt) / norm(normal); + assert( is_matrix([normal,pt],2,3) && !approx(norm(normal),0), + "Inputs `normal` and `pt` should be 3d vectors/points and `normal` cannot be zero." ) + concat(normal, normal*pt) / norm(normal); +// Eigenvalues for a 3x3 symmetrical matrix in decreasing order +// Based on: https://en.wikipedia.org/wiki/Eigenvalue_algorithm +function _eigenvals_symm_3(M) = + let( p1 = pow(M[0][1],2) + pow(M[0][2],2) + pow(M[1][2],2) ) + (p1, ); @@ -900,7 +938,7 @@ function plane_from_normal(normal, pt=[0,0,0]) = // Given a list of 3 or more coplanar 3D points, returns the coefficients of the normalized cartesian equation of a plane, // that is [A,B,C,D] where Ax+By+Cz=D is the equation of the plane where norm([A,B,C])=1. // If `fast` is false and the points in the list are collinear or not coplanar, then `undef` is returned. -// if `fast` is true, then the coplanarity test is skipped and a plane passing through 3 non-collinear arbitrary points is returned. +// If `fast` is true, the polygon coplanarity check is skipped and a best fitted plane is returned. // Arguments: // points = The list of points to find the plane of. // fast = If true, don't verify that all points in the list are coplanar. Default: false @@ -914,17 +952,18 @@ function plane_from_normal(normal, pt=[0,0,0]) = function plane_from_points(points, fast=false, eps=EPSILON) = assert( is_path(points,dim=3), "Improper 3d point list." ) assert( is_finite(eps) && (eps>=0), "The tolerance should be a non-negative value." ) - let( - indices = noncollinear_triple(points,error=false) - ) - indices==[] ? undef : - let( - p1 = points[indices[0]], - p2 = points[indices[1]], - p3 = points[indices[2]], - plane = plane3pt(p1,p2,p3) - ) - fast || points_on_plane(points,plane,eps=eps) ? plane : undef; + len(points) == 3 + ? let( plane = plane3pt(points[0],points[1],points[2]) ) + plane==[] ? undef : plane + : let( + cov_evals = _covariance_evals(points), + pm = cov_evals[0], + evals = cov_evals[1], + M = cov_evals[2], + evec = _eigenvec_symm_3(M,evals,i=2) ) +// echo(error_points_plane= abs(max(points*evec)-pm*evec), limit=eps) + !fast && abs(max(points*evec)-pm*evec)>eps*evals[0] ? undef : + [ each evec, pm*evec] ; // Function: plane_from_polygon() @@ -934,7 +973,8 @@ function plane_from_points(points, fast=false, eps=EPSILON) = // Given a 3D planar polygon, returns the normalized cartesian equation of its plane. // Returns [A,B,C,D] where Ax+By+Cz=D is the equation of the plane where norm([A,B,C])=1. // If not all the points in the polygon are coplanar, then [] is returned. -// If `fast` is true, the polygon coplanarity check is skipped and the plane may not contain all polygon points. +// If `fast` is false and the points in the list are collinear or not coplanar, then `undef` is returned. +// if `fast` is true, then the coplanarity test is skipped and a plane passing through 3 non-collinear arbitrary points is returned. // Arguments: // poly = The planar 3D polygon to find the plane of. // fast = If true, doesn't verify that all points in the polygon are coplanar. Default: false @@ -948,14 +988,13 @@ function plane_from_points(points, fast=false, eps=EPSILON) = function plane_from_polygon(poly, fast=false, eps=EPSILON) = assert( is_path(poly,dim=3), "Invalid polygon." ) assert( is_finite(eps) && (eps>=0), "The tolerance should be a non-negative value." ) - let( - poly = deduplicate(poly), - n = polygon_normal(poly), - plane = [n.x, n.y, n.z, n*poly[0]] - ) - fast? plane: coplanar(poly,eps=eps)? plane: []; - - + len(poly)==3 ? plane3pt(poly[0],poly[1],poly[2]) : + let( triple = sort(noncollinear_triple(poly,error=false)) ) + triple==[] ? [] : + let( plane = plane3pt(poly[triple[0]],poly[triple[1]],poly[triple[2]])) + fast? plane: points_on_plane(poly, plane, eps=eps)? plane: []; + + // Function: plane_normal() // Usage: // plane_normal(plane); @@ -1252,9 +1291,11 @@ function coplanar(points, eps=EPSILON) = len(points)<=2 ? false : let( ip = noncollinear_triple(points,error=false,eps=eps) ) ip == [] ? false : - let( plane = plane3pt(points[ip[0]],points[ip[1]],points[ip[2]]), - normal = point3d(plane) ) - max( points*normal ) - plane[3]< eps*norm(normal); + let( + plane = plane3pt(points[ip[0]],points[ip[1]],points[ip[2]]), + normal = point3d(plane), + pt_nrm = points*normal ) + abs(max(max(pt_nrm)-plane[3], -min(pt_nrm)+plane[3])) < eps; // Function: points_on_plane() @@ -1665,11 +1706,11 @@ function noncollinear_triple(points,error=true,eps=EPSILON) = n = (pb-pa)/nrm, distlist = [for(i=[0:len(points)-1]) _dist2line(points[i]-pa, n)] ) - max(distlist)0, "Invalid pointlist." ) // using is_matrix it accepts len(points)==1 assert( is_finite(eps) && eps>=0, "The tolerance should be a positive number." ) - let( normal = point3d(plane), - pt_nrm = points*normal ) - abs(max( max(pt_nrm) - plane[3], -min(pt_nrm)+plane[3]))< eps*norm(normal); + _pointlist_greatest_distance(points,plane) < eps; // Function: in_front_of_plane() From 160e3f3edd390a91c28f1a373d413ca41ab2525e Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 11:35:38 +0100 Subject: [PATCH 12/21] Correction of polygon_area and a better polygon_normal --- geometry.scad | 15 ++++++--------- 1 file changed, 6 insertions(+), 9 deletions(-) diff --git a/geometry.scad b/geometry.scad index 160ecd5..7f06614 100644 --- a/geometry.scad +++ b/geometry.scad @@ -1784,7 +1784,7 @@ function polygon_area(poly, signed=false) = len(poly[0])==2 ? let( total = sum([for(i=[1:1:len(poly)-2]) cross(poly[i]-poly[0],poly[i+1]-poly[0]) ])/2 ) signed ? total : abs(total) - : let( plane = plane_from_points(poly) ) + : let( plane = plane_from_polygon(poly) ) plane==undef? undef : let( n = plane_normal(plane), @@ -2101,18 +2101,15 @@ function reverse_polygon(poly) = // Description: // Given a 3D planar polygon, returns a unit-length normal vector for the // clockwise orientation of the polygon. If the polygon points are collinear, returns `undef`. +// It doesn't check for coplanarity. // Arguments: // poly = The list of 3D path points for the perimeter of the polygon. function polygon_normal(poly) = assert(is_path(poly,dim=3), "Invalid 3D polygon." ) - let( - poly = cleanup_path(poly), - p0 = poly[0], - n = sum([ - for (i=[1:1:len(poly)-2]) - cross(poly[i+1]-p0, poly[i]-p0) - ]) - ) unit(n,undef); + len(poly)==3 ? point3d(plane3pt(poly[0],poly[1],poly[2])) : + let( triple = sort(noncollinear_triple(poly,error=false)) ) + triple==[] ? undef : + point3d(plane3pt(poly[triple[0]],poly[triple[1]],poly[triple[2]])) ; function _split_polygon_at_x(poly, x) = From 80feb93c98fc03bf6bb0b16c03c5628e11637c2c Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 12:32:49 +0100 Subject: [PATCH 13/21] Change undef to [] as return of polygon functions --- geometry.scad | 39 +++++++++++++++++++-------------------- 1 file changed, 19 insertions(+), 20 deletions(-) diff --git a/geometry.scad b/geometry.scad index 7f06614..d0e98df 100644 --- a/geometry.scad +++ b/geometry.scad @@ -453,7 +453,7 @@ function segment_closest_point(seg,pt) = // Usage: // line_from_points(points, [fast], [eps]); // Description: -// Given a list of 2 or more colinear points, returns a line containing them. +// Given a list of 2 or more collinear points, returns a line containing them. // If `fast` is false and the points are coincident, then `undef` is returned. // if `fast` is true, then the collinearity test is skipped and a line passing through 2 distinct arbitrary points is returned. // Arguments: @@ -916,12 +916,12 @@ function _eigenvals_symm_3(M) = // using Cayley–Hamilton theorem according to: // https://en.wikipedia.org/wiki/Eigenvalue_algorithm function _eigenvec_symm_3(M,evals,i=0) = - let( - I = ident(3), - A = (M - evals[(i+1)%3]*I) * (M - evals[(i+2)%3]*I) , - k = max_index( [for(i=[0:2]) norm(A[i]) ]) - ) - norm(A[k])=0), "The tolerance should be a non-negative value." ) len(points) == 3 ? let( plane = plane3pt(points[0],points[1],points[2]) ) - plane==[] ? undef : plane + plane==[] ? [] : plane : let( covmix = _covariance_evec_eval(points), pm = covmix[0], @@ -976,7 +976,7 @@ function plane_from_points(points, fast=false, eps=EPSILON) = // Given a 3D planar polygon, returns the normalized cartesian equation of its plane. // Returns [A,B,C,D] where Ax+By+Cz=D is the equation of the plane where norm([A,B,C])=1. // If not all the points in the polygon are coplanar, then [] is returned. -// If `fast` is false and the points in the list are collinear or not coplanar, then `undef` is returned. +// If `fast` is false and the points in the list are collinear or not coplanar, then [] is returned. // if `fast` is true, then the coplanarity test is skipped and a plane passing through 3 non-collinear arbitrary points is returned. // Arguments: // poly = The planar 3D polygon to find the plane of. @@ -1301,10 +1301,10 @@ function coplanar(points, eps=EPSILON) = // the maximum distance from points to the plane function _pointlist_greatest_distance(points,plane) = let( - normal = point3d(plane), - pt_nrm = points*normal + normal = point3d(plane), + pt_nrm = points*normal ) - abs(max( max(pt_nrm) - plane[3], -min(pt_nrm)+plane[3])) / norm(normal); + abs(max( max(pt_nrm) - plane[3], -min(pt_nrm)+plane[3])) / norm(normal); // Function: points_on_plane() @@ -1647,20 +1647,19 @@ function circle_circle_tangents(c1,r1,c2,r2,d1,d2) = // Function: circle_line_intersection() // Usage: -// isect = circle_line_intersection(c,r,line,,); -// isect = circle_line_intersection(c,d,line,,); +// isect = circle_line_intersection(c,,,,); // Description: // Find intersection points between a 2d circle and a line, ray or segment specified by two points. // By default the line is unbounded. // Arguments: // c = center of circle // r = radius of circle +// --- +// d = diameter of circle // line = two points defining the unbounded line // bounded = false for unbounded line, true for a segment, or a vector [false,true] or [true,false] to specify a ray with the first or second end unbounded. Default: false // eps = epsilon used for identifying the case with one solution. Default: 1e-9 -// --- -// d = diameter of circle -function circle_line_intersection(c,r,line,d,bounded=false,eps=EPSILON) = +function circle_line_intersection(c,r,d,line,bounded=false,eps=EPSILON) = let(r=get_radius(r=r,d=d,dflt=undef)) assert(_valid_line(line,2), "Input 'line' is not a valid 2d line.") assert(is_vector(c,2), "Circle center must be a 2-vector") @@ -2100,7 +2099,7 @@ function reverse_polygon(poly) = // n = polygon_normal(poly); // Description: // Given a 3D planar polygon, returns a unit-length normal vector for the -// clockwise orientation of the polygon. If the polygon points are collinear, returns `undef`. +// clockwise orientation of the polygon. If the polygon points are collinear, returns []. // It doesn't check for coplanarity. // Arguments: // poly = The list of 3D path points for the perimeter of the polygon. @@ -2108,7 +2107,7 @@ function polygon_normal(poly) = assert(is_path(poly,dim=3), "Invalid 3D polygon." ) len(poly)==3 ? point3d(plane3pt(poly[0],poly[1],poly[2])) : let( triple = sort(noncollinear_triple(poly,error=false)) ) - triple==[] ? undef : + triple==[] ? [] : point3d(plane3pt(poly[triple[0]],poly[triple[1]],poly[triple[2]])) ; From f9e7baa5157c7970eb5fbd36c46fbe68821875f1 Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 12:33:40 +0100 Subject: [PATCH 14/21] Review of is_matrix --- math.scad | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/math.scad b/math.scad index 43b1347..1f6cb92 100644 --- a/math.scad +++ b/math.scad @@ -985,10 +985,8 @@ function determinant(M) = function is_matrix(A,m,n,square=false) = is_list(A) && (( is_undef(m) && len(A) ) || len(A)==m) - && is_list(A[0]) - && (( is_undef(n) && len(A[0]) ) || len(A[0])==n) && (!square || len(A) == len(A[0])) - && is_vector(A[0]) + && is_vector(A[0],n) && is_consistent(A); From 376609d52f1ad754d2aadc1538dcd42bca825461 Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 12:34:59 +0100 Subject: [PATCH 15/21] Cosmetics --- common.scad | 11 +-- paths.scad | 216 ++++++++++++++++++++++++++-------------------------- 2 files changed, 114 insertions(+), 113 deletions(-) diff --git a/common.scad b/common.scad index e24c744..6d2d44f 100644 --- a/common.scad +++ b/common.scad @@ -205,7 +205,8 @@ function is_func(x) = version_num()>20210000 && is_function(x); // Description: // Tests whether input is a list of entries which all have the same list structure // and are filled with finite numerical data. You can optionally specify a required -// list structure with the pattern argument. It returns `true` for the empty list. +// list structure with the pattern argument. +// It returns `true` for the empty list regardless the value of the `pattern`. // Arguments: // list = list to check // pattern = optional pattern required to match @@ -293,7 +294,7 @@ function default(v,dflt=undef) = is_undef(v)? dflt : v; // v = The list whose items are being checked. // recursive = If true, sublists are checked recursively for defined values. The first sublist that has a defined item is returned. // Examples: -// val = first_defined([undef,7,undef,true]); // Returns: 1 +// val = first_defined([undef,7,undef,true]); // Returns: 7 function first_defined(v,recursive=false,_i=0) = _i 0? _corner_chamfer_path(p1, p2, p3, side=chamfer[i]) : - rounding[i] > 0? _corner_roundover_path(p1, p2, p3, r=rounding[i]) : - [p2] - ], - out = [ - if (!closed) path[0], - for (i=(closed? [0:1:lp-1] : [1:1:lp-2])) let( - p1 = select(path,i-1), - p2 = select(path,i), - crn1 = select(corner_paths,i-1), - crn2 = corner_paths[i], - l1 = norm(last(crn1)-p1), - l2 = norm(crn2[0]-p2), - needed = l1 + l2, - seglen = norm(p2-p1), - check = assert(seglen >= needed, str("Path segment ",i," is too short to fulfill rounding/chamfering for the adjacent corners.")) - ) each crn2, - if (!closed) last(path) - ] - ) deduplicate(out); + let ( + path = deduplicate(path,closed=true), + lp = len(path), + chamfer = is_undef(chamfer)? repeat(0,lp) : + is_vector(chamfer)? list_pad(chamfer,lp,0) : + is_num(chamfer)? repeat(chamfer,lp) : + assert(false, "Bad chamfer value."), + rounding = is_undef(rounding)? repeat(0,lp) : + is_vector(rounding)? list_pad(rounding,lp,0) : + is_num(rounding)? repeat(rounding,lp) : + assert(false, "Bad rounding value."), + corner_paths = [ + for (i=(closed? [0:1:lp-1] : [1:1:lp-2])) let( + p1 = select(path,i-1), + p2 = select(path,i), + p3 = select(path,i+1) + ) + chamfer[i] > 0? _corner_chamfer_path(p1, p2, p3, side=chamfer[i]) : + rounding[i] > 0? _corner_roundover_path(p1, p2, p3, r=rounding[i]) : + [p2] + ], + out = [ + if (!closed) path[0], + for (i=(closed? [0:1:lp-1] : [1:1:lp-2])) let( + p1 = select(path,i-1), + p2 = select(path,i), + crn1 = select(corner_paths,i-1), + crn2 = corner_paths[i], + l1 = norm(last(crn1)-p1), + l2 = norm(crn2[0]-p2), + needed = l1 + l2, + seglen = norm(p2-p1), + check = assert(seglen >= needed, str("Path segment ",i," is too short to fulfill rounding/chamfering for the adjacent corners.")) + ) each crn2, + if (!closed) last(path) + ] + ) deduplicate(out); function _corner_chamfer_path(p1, p2, p3, dist1, dist2, side, angle) = - let( - v1 = unit(p1 - p2), - v2 = unit(p3 - p2), - n = vector_axis(v1,v2), - ang = vector_angle(v1,v2), - path = (is_num(dist1) && is_undef(dist2) && is_undef(side))? ( - // dist1 & optional angle - assert(dist1 > 0) - let(angle = default(angle,(180-ang)/2)) - assert(is_num(angle)) - assert(angle > 0 && angle < 180) - let( - pta = p2 + dist1*v1, - a3 = 180 - angle - ang - ) assert(a3>0, "Angle too extreme.") - let( - side = sin(angle) * dist1/sin(a3), - ptb = p2 + side*v2 - ) [pta, ptb] - ) : (is_undef(dist1) && is_num(dist2) && is_undef(side))? ( - // dist2 & optional angle - assert(dist2 > 0) - let(angle = default(angle,(180-ang)/2)) - assert(is_num(angle)) - assert(angle > 0 && angle < 180) - let( - ptb = p2 + dist2*v2, - a3 = 180 - angle - ang - ) assert(a3>0, "Angle too extreme.") - let( - side = sin(angle) * dist2/sin(a3), - pta = p2 + side*v1 - ) [pta, ptb] - ) : (is_undef(dist1) && is_undef(dist2) && is_num(side))? ( - // side & optional angle - assert(side > 0) - let(angle = default(angle,(180-ang)/2)) - assert(is_num(angle)) - assert(angle > 0 && angle < 180) - let( - a3 = 180 - angle - ang - ) assert(a3>0, "Angle too extreme.") - let( - dist1 = sin(a3) * side/sin(ang), - dist2 = sin(angle) * side/sin(ang), - pta = p2 + dist1*v1, - ptb = p2 + dist2*v2 - ) [pta, ptb] - ) : (is_num(dist1) && is_num(dist2) && is_undef(side) && is_undef(side))? ( - // dist1 & dist2 - assert(dist1 > 0) - assert(dist2 > 0) - let( - pta = p2 + dist1*v1, - ptb = p2 + dist2*v2 - ) [pta, ptb] - ) : ( - assert(false,"Bad arguments.") - ) - ) path; + let( + v1 = unit(p1 - p2), + v2 = unit(p3 - p2), + n = vector_axis(v1,v2), + ang = vector_angle(v1,v2), + path = (is_num(dist1) && is_undef(dist2) && is_undef(side))? ( + // dist1 & optional angle + assert(dist1 > 0) + let(angle = default(angle,(180-ang)/2)) + assert(is_num(angle)) + assert(angle > 0 && angle < 180) + let( + pta = p2 + dist1*v1, + a3 = 180 - angle - ang + ) assert(a3>0, "Angle too extreme.") + let( + side = sin(angle) * dist1/sin(a3), + ptb = p2 + side*v2 + ) [pta, ptb] + ) : (is_undef(dist1) && is_num(dist2) && is_undef(side))? ( + // dist2 & optional angle + assert(dist2 > 0) + let(angle = default(angle,(180-ang)/2)) + assert(is_num(angle)) + assert(angle > 0 && angle < 180) + let( + ptb = p2 + dist2*v2, + a3 = 180 - angle - ang + ) assert(a3>0, "Angle too extreme.") + let( + side = sin(angle) * dist2/sin(a3), + pta = p2 + side*v1 + ) [pta, ptb] + ) : (is_undef(dist1) && is_undef(dist2) && is_num(side))? ( + // side & optional angle + assert(side > 0) + let(angle = default(angle,(180-ang)/2)) + assert(is_num(angle)) + assert(angle > 0 && angle < 180) + let( + a3 = 180 - angle - ang + ) assert(a3>0, "Angle too extreme.") + let( + dist1 = sin(a3) * side/sin(ang), + dist2 = sin(angle) * side/sin(ang), + pta = p2 + dist1*v1, + ptb = p2 + dist2*v2 + ) [pta, ptb] + ) : (is_num(dist1) && is_num(dist2) && is_undef(side) && is_undef(side))? ( + // dist1 & dist2 + assert(dist1 > 0) + assert(dist2 > 0) + let( + pta = p2 + dist1*v1, + ptb = p2 + dist2*v2 + ) [pta, ptb] + ) : ( + assert(false,"Bad arguments.") + ) + ) path; function _corner_roundover_path(p1, p2, p3, r, d) = - let( - r = get_radius(r=r,d=d,dflt=undef), - res = circle_2tangents(p1, p2, p3, r=r, tangents=true), - cp = res[0], - n = res[1], - tp1 = res[2], - ang = res[4]+res[5], - steps = floor(segs(r)*ang/360+0.5), - step = ang / steps, - path = [for (i=[0:1:steps]) move(cp, p=rot(a=-i*step, v=n, p=tp1-cp))] - ) path; + let( + r = get_radius(r=r,d=d,dflt=undef), + res = circle_2tangents(p1, p2, p3, r=r, tangents=true), + cp = res[0], + n = res[1], + tp1 = res[2], + ang = res[4]+res[5], + steps = floor(segs(r)*ang/360+0.5), + step = ang / steps, + path = [for (i=[0:1:steps]) move(cp, p=rot(a=-i*step, v=n, p=tp1-cp))] + ) path; From c5799d65396a1cb711363f38b6b9f40d2877a2f4 Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 12:36:05 +0100 Subject: [PATCH 16/21] Review of vnf_centroid() --- tests/test_vnf.scad | 3 ++- vnf.scad | 16 ++++++---------- 2 files changed, 8 insertions(+), 11 deletions(-) diff --git a/tests/test_vnf.scad b/tests/test_vnf.scad index bb171d4..b83775e 100644 --- a/tests/test_vnf.scad +++ b/tests/test_vnf.scad @@ -74,7 +74,8 @@ module test_vnf_centroid() { assert_approx(vnf_centroid(cube(100, anchor=TOP)), [0,0,-50]); assert_approx(vnf_centroid(sphere(d=100, anchor=CENTER, $fn=36)), [0,0,0]); assert_approx(vnf_centroid(sphere(d=100, anchor=BOT, $fn=36)), [0,0,50]); -} + ellipse = xscale(2, p=circle($fn=24, r=3)); + assert_approx(vnf_centroid(path_sweep(pentagon(r=1), path3d(ellipse), closed=true)),[0,0,0]);} test_vnf_centroid(); diff --git a/vnf.scad b/vnf.scad index 9c05e0f..7d1c021 100644 --- a/vnf.scad +++ b/vnf.scad @@ -450,18 +450,13 @@ function vnf_volume(vnf) = // Returns the centroid of the given manifold VNF. The VNF must describe a valid polyhedron with consistent face direction and // no holes; otherwise the results are undefined. -// Divide the solid up into tetrahedra with the origin as one vertex. The centroid of a tetrahedron is the average of its vertices. +// Divide the solid up into tetrahedra with the origin as one vertex. +// The centroid of a tetrahedron is the average of its vertices. // The centroid of the total is the volume weighted average. function vnf_centroid(vnf) = + assert(is_vnf(vnf) && len(vnf[0])!=0 ) let( verts = vnf[0], - vol = sum([ - for(face=vnf[1], j=[1:1:len(face)-2]) let( - v0 = verts[face[0]], - v1 = verts[face[j]], - v2 = verts[face[j+1]] - ) cross(v2,v1)*v0 - ]), pos = sum([ for(face=vnf[1], j=[1:1:len(face)-2]) let( v0 = verts[face[0]], @@ -469,10 +464,11 @@ function vnf_centroid(vnf) = v2 = verts[face[j+1]], vol = cross(v2,v1)*v0 ) - (v0+v1+v2)*vol + [ vol, (v0+v1+v2)*vol ] ]) ) - pos/vol/4; + assert(!approx(pos[0],0, EPSILON), "The vnf has self-intersections.") + pos[1]/pos[0]/4; function _triangulate_planar_convex_polygons(polys) = From 5062ee1605db3827d703b3242467e8109a9aa44b Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 12:37:49 +0100 Subject: [PATCH 17/21] Minor change in is_matrix_symmetric --- arrays.scad | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/arrays.scad b/arrays.scad index ac97444..604186b 100644 --- a/arrays.scad +++ b/arrays.scad @@ -41,6 +41,7 @@ function is_homogeneous(l, depth=10) = [] == [for(i=[1:len(l)-1]) if( ! _same_type(l[i],l0, depth+1) ) 0 ]; function is_homogenous(l, depth=10) = is_homogeneous(l, depth); + function _same_type(a,b, depth) = (depth==0) || @@ -50,7 +51,7 @@ function _same_type(a,b, depth) = (is_string(a) && is_string(b)) || (is_list(a) && is_list(b) && len(a)==len(b) && []==[for(i=idx(a)) if( ! _same_type(a[i],b[i],depth-1) ) 0] ); - + // Function: select() // Topics: List Handling @@ -97,7 +98,6 @@ function select(list, start, end) = // Function: slice() -// Topics: List Handling // Usage: // list = slice(list,s,e); // Description: @@ -476,7 +476,7 @@ function reverse(x) = // l9 = list_rotate([1,2,3,4,5],6); // Returns: [2,3,4,5,1] function list_rotate(list,n=1) = assert(is_list(list)||is_string(list), "Invalid list or string.") - assert(is_finite(n), "Invalid number") + assert(is_int(n), "The rotation number should be integer") let ( ll = len(list), n = ((n % ll) + ll) % ll, @@ -990,7 +990,7 @@ function _sort_vectors(arr, idxlist, _i=0) = _sort_vectors(equal, idxlist, _i+1), _sort_vectors(greater, idxlist, _i ) ); - + // sorting using compare_vals(); returns indexed list when `indexed==true` function _sort_general(arr, idx=undef, indexed=false) = (len(arr)<=1) ? arr : @@ -1332,6 +1332,8 @@ function permutations(l,n=2) = // pairs = zip(a,b); // triples = zip(a,b,c); // quads = zip([LIST1,LIST2,LIST3,LIST4]); +// Topics: List Handling, Iteration +// See Also: zip_long() // Description: // Zips together two or more lists into a single list. For example, if you have two // lists [3,4,5], and [8,7,6], and zip them together, you get [[3,8],[4,7],[5,6]]. @@ -1357,6 +1359,8 @@ function zip(a,b,c) = // pairs = zip_long(a,b); // triples = zip_long(a,b,c); // quads = zip_long([LIST1,LIST2,LIST3,LIST4]); +// Topics: List Handling, Iteration +// See Also: zip() // Description: // Zips together two or more lists into a single list. For example, if you have two // lists [3,4,5], and [8,7,6], and zip them together, you get [[3,8],[4,7],[5,6]]. @@ -1526,7 +1530,6 @@ function subindex(M, idx) = // [[4,2], 91, false], // [6, [3,4], undef]]; // submatrix(A,[0,2],[1,2]); // Returns [[17, "test"], [[3, 4], undef]] - function submatrix(M,idx1,idx2) = [for(i=idx1) [for(j=idx2) M[i][j] ] ]; @@ -1629,7 +1632,6 @@ function block_matrix(M) = assert(badrows==[], "Inconsistent or invalid input") bigM; - // Function: diagonal_matrix() // Usage: // mat = diagonal_matrix(diag, ); @@ -1855,7 +1857,7 @@ function transpose(arr, reverse=false) = // A = matrix to test // eps = epsilon for comparing equality. Default: 1e-12 function is_matrix_symmetric(A,eps=1e-12) = - approx(A,transpose(A)); + approx(A,transpose(A), eps); // vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap From dab6407b0f0bcf610e81ffd87fba0e36cdad065c Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 13:01:06 +0100 Subject: [PATCH 18/21] Minor tweak --- geometry.scad | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/geometry.scad b/geometry.scad index d0e98df..26fa9ad 100644 --- a/geometry.scad +++ b/geometry.scad @@ -1772,7 +1772,7 @@ function furthest_point(pt, points) = // area = polygon_area(poly); // Description: // Given a 2D or 3D planar polygon, returns the area of that polygon. -// If the polygon is self-crossing, the results are undefined. For non-planar 3D polygon the result is undef. +// If the polygon is self-crossing, the results are undefined. For non-planar 3D polygon the result is []. // When `signed` is true, a signed area is returned; a positive area indicates a clockwise polygon. // Arguments: // poly = Polygon to compute the area of. @@ -1784,9 +1784,9 @@ function polygon_area(poly, signed=false) = ? let( total = sum([for(i=[1:1:len(poly)-2]) cross(poly[i]-poly[0],poly[i+1]-poly[0]) ])/2 ) signed ? total : abs(total) : let( plane = plane_from_polygon(poly) ) - plane==undef? undef : + plane==[]? [] : let( - n = plane_normal(plane), + n = plane_normal(plane), total = sum([ for(i=[1:1:len(poly)-2]) let( From 41616872fe0185627ce014b2cc791642c824d726 Mon Sep 17 00:00:00 2001 From: RonaldoCMP Date: Sun, 11 Apr 2021 23:45:21 +0100 Subject: [PATCH 19/21] Correction in is_convex_polygon --- geometry.scad | 11 ++++++----- tests/test_geometry.scad | 1 + 2 files changed, 7 insertions(+), 5 deletions(-) diff --git a/geometry.scad b/geometry.scad index 26fa9ad..6995f7f 100644 --- a/geometry.scad +++ b/geometry.scad @@ -1304,7 +1304,7 @@ function _pointlist_greatest_distance(points,plane) = normal = point3d(plane), pt_nrm = points*normal ) - abs(max( max(pt_nrm) - plane[3], -min(pt_nrm)+plane[3])) / norm(normal); + abs(max( max(pt_nrm) - plane[3], -min(pt_nrm) + plane[3])) / norm(normal); // Function: points_on_plane() @@ -1786,7 +1786,7 @@ function polygon_area(poly, signed=false) = : let( plane = plane_from_polygon(poly) ) plane==[]? [] : let( - n = plane_normal(plane), + n = plane_normal(plane), total = sum([ for(i=[1:1:len(poly)-2]) let( @@ -1808,25 +1808,26 @@ function polygon_area(poly, signed=false) = // If the points are collinear an error is generated. // Arguments: // poly = Polygon to check. +// eps = Tolerance for the collinearity test. Default: EPSILON. // Example: // is_convex_polygon(circle(d=50)); // Returns: true // is_convex_polygon(rot([50,120,30], p=path3d(circle(1,$fn=50)))); // Returns: true // Example: // spiral = [for (i=[0:36]) let(a=-i*10) (10+i)*[cos(a),sin(a)]]; // is_convex_polygon(spiral); // Returns: false -function is_convex_polygon(poly) = +function is_convex_polygon(poly,eps=EPSILON) = assert(is_path(poly), "The input should be a 2D or 3D polygon." ) let( lp = len(poly), p0 = poly[0] ) assert( lp>=3 , "A polygon must have at least 3 points" ) let( crosses = [for(i=[0:1:lp-1]) cross(poly[(i+1)%lp]-poly[i], poly[(i+2)%lp]-poly[(i+1)%lp]) ] ) len(p0)==2 - ? assert( !approx(max(crosses)) && !approx(min(crosses)), "The points are collinear" ) + ? assert( !approx(sqrt(max(max(crosses),-min(crosses))),eps), "The points are collinear" ) min(crosses) >=0 || max(crosses)<=0 : let( prod = crosses*sum(crosses), minc = min(prod), maxc = max(prod) ) - assert( !approx(maxc-minc), "The points are collinear" ) + assert( !approx(sqrt(max(maxc,-minc)),eps), "The points are collinear" ) minc>=0 || maxc<=0; diff --git a/tests/test_geometry.scad b/tests/test_geometry.scad index 69d4705..ccb4ac7 100644 --- a/tests/test_geometry.scad +++ b/tests/test_geometry.scad @@ -847,6 +847,7 @@ module test_is_convex_polygon() { assert(is_convex_polygon(circle(r=50,$fn=1000))); assert(is_convex_polygon(rot([50,120,30], p=path3d(circle(1,$fn=50))))); assert(!is_convex_polygon([[1,1],[0,0],[-1,1],[-1,-1],[1,-1]])); + assert(!is_convex_polygon([for (i=[0:36]) let(a=-i*10) (10+i)*[cos(a),sin(a)]])); // spiral } *test_is_convex_polygon(); From d4d5794ef22e6421ceebdb8cdbe9d24628492360 Mon Sep 17 00:00:00 2001 From: Adrian Mariano Date: Sun, 11 Apr 2021 22:29:25 -0400 Subject: [PATCH 20/21] fixed sum to handle larger input added all_equal added vnf_tri_array as public added bezier_patch_degenerate as public --- beziers.scad | 187 ++++++++++++++++++++++++++++++++++++++++++++++++++ geometry.scad | 6 +- math.scad | 14 +++- rounding.scad | 132 +++-------------------------------- 4 files changed, 214 insertions(+), 125 deletions(-) diff --git a/beziers.scad b/beziers.scad index bd6a8f0..7cda75b 100644 --- a/beziers.scad +++ b/beziers.scad @@ -1222,6 +1222,193 @@ function bezier_patch(patch, splinesteps=16, vnf=EMPTY_VNF, style="default") = ) vnf; + +// Function: bezier_patch_degenerate() +// Usage: +// vnf = bezier_patch_degenerate(patch, , ); +// vnf_edges = bezier_patch_degenerate(patch, , , return_edges=true); +// Description: +// Returns a VNF for a degenerate rectangular bezier patch where some of the corners of the patch are +// equal. If the resulting patch has no faces then returns an empty VNF. Note that due to the degeneracy, +// the shape of the patch can be triangular even though the actual underlying patch is a rectangle. This is +// a different method for creating triangular bezier patches than the triangular patch. +// If you specify return_edges then the return is a list whose first element is the vnf and whose second +// element lists the edges in the order [left, right, top, bottom], where each list is a list of the actual +// point values, but possibly only a single point if that edge is degenerate. +// The method checks for various types of degeneracy and uses a triangular or partly triangular array of sample points. +// See examples below for the types of degeneracy detected and how the patch is sampled for those cases. +// Note that splinesteps is the same for both directions of the patch, so it cannot be an array. +// Arguments: +// patch = Patch to process +// splinesteps = Number of segments to produce on each side. Default: 16 +// reverse = reverse direction of faces. Default: false +// return_edges = if true return the points on the four edges: [left, right, top, bottom]. Default: false +// Example: This quartic patch is degenerate at one corner, where a row of control points are equal. Processing this degenerate patch normally produces excess triangles near the degenerate point. +// splinesteps=8; +// patch=[ +// repeat([-12.5, 12.5, 15],5), +// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]], +// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]], +// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]], +// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]] +// ]; +// vnf_wireframe((bezier_patch(patch, splinesteps)),d=0.1); +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// Example: With bezier_patch_degenerate the degenerate point does not have excess triangles. The top half of the patch decreases the number of sampled points by 2 for each row. +// splinesteps=8; +// patch=[ +// repeat([-12.5, 12.5, 15],5), +// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]], +// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]], +// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]], +// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]] +// ]; +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// Example: With splinesteps odd you get one "odd" row where the point count decreases by 1 instead of 2. You may prefer even values for splinesteps to avoid this. +// splinesteps=7; +// patch=[ +// repeat([-12.5, 12.5, 15],5), +// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]], +// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]], +// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]], +// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]] +// ]; +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// Example: A more extreme degeneracy occurs when the top half of a patch is degenerate to a line. (For odd length patches the middle row must be degenerate to trigger this style.) In this case the number of points in each row decreases by 1 for every row. It doesn't matter of splinesteps is odd or even. +// splinesteps=8; +// patch = [[[10, 0, 0], [10, -10.4, 0], [10, -20.8, 0], [1.876, -14.30, 0], [-6.24, -7.8, 0]], +// [[5, 0, 0], [5, -5.2, 0], [5, -10.4, 0], [0.938, -7.15, 0], [-3.12, -3.9, 0]], +// repeat([0,0,0],5), +// repeat([0,0,5],5), +// repeat([0,0,10],5) +// ]; +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// Example: Here is a degenerate cubic patch. +// splinesteps=8; +// patch = [ [ [-20,0,0], [-10,0,0],[0,10,0],[0,20,0] ], +// [ [-20,0,10], [-10,0,10],[0,10,10],[0,20,10]], +// [ [-10,0,20], [-5,0,20], [0,5,20], [0,10,20]], +// repeat([0,0,30],4) +// ]; +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// Example: A more extreme degenerate cubic patch, where two rows are equal. +// splinesteps=8; +// patch = [ [ [-20,0,0], [-10,0,0],[0,10,0],[0,20,0] ], +// [ [-20,0,10], [-10,0,10],[0,10,10],[0,20,10] ], +// repeat([-10,10,20],4), +// repeat([-10,10,30],4) +// ]; +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// Example: Quadratic patch degenerate at the right side: +// splinesteps=8; +// patch = [[[0, -10, 0],[10, -5, 0],[20, 0, 0]], +// [[0, 0, 0], [10, 0, 0], [20, 0, 0]], +// [[0, 0, 10], [10, 0, 5], [20, 0, 0]]]; +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +// Example: Cubic patch degenerate at both ends. In this case the point count changes by 2 at every row. +// splinesteps=8; +// patch = [ +// repeat([10,-10,0],4), +// [ [-20,0,0], [-1,0,0],[0,10,0],[0,20,0] ], +// [ [-20,0,10], [-10,0,10],[0,10,10],[0,20,10] ], +// repeat([-10,10,20],4), +// ]; +// vnf_wireframe(bezier_patch_degenerate(patch, splinesteps),d=0.1); +// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9); +function bezier_patch_degenerate(patch, splinesteps=16, reverse=false, return_edges=false) = + !return_edges ? bezier_patch_degenerate(patch, splinesteps, reverse, true)[0] : + assert(is_rectpatch(patch), "Must supply rectangular bezier patch") + assert(is_int(splinesteps) && splinesteps>=3, "splinesteps must be an integer 3 or larger") + let( + row_degen = [for(row=patch) all_equal(row)], + col_degen = [for(col=transpose(patch)) all_equal(col)], + + top_degen = row_degen[0], + bot_degen = last(row_degen), + left_degen = col_degen[0], + right_degen = last(col_degen), + samplepts = lerpn(0,1,splinesteps+1) + ) + all(row_degen) && all(col_degen) ? // fully degenerate case + [EMPTY_VNF, repeat([patch[0][0]],4)] : + all(row_degen) ? // degenerate to a line (top to bottom) + let(pts = bezier_points(subindex(patch,0), samplepts)) + [EMPTY_VNF, [pts,pts,[pts[0]],[last(pts)]]] : + all(col_degen) ? // degenerate to a line (left to right) + let(pts = bezier_points(patch[0], samplepts)) + [EMPTY_VNF, [[pts[0]], [last(pts)], pts, pts]] : + !top_degen && !bot_degen && !left_degen && !right_degen ? // non-degenerate case + let(pts = bezier_patch_points(patch, samplepts, samplepts)) + [ + vnf_vertex_array(pts, reverse=!reverse), + [subindex(pts,0), subindex(pts,len(pts)-1), pts[0], last(pts)] + ] : + top_degen && bot_degen ? + let( + rowcount = [ + each list([3:2:splinesteps]), + if (splinesteps%2==0) splinesteps+1, + each reverse(list([3:2:splinesteps])) + ], + bpatch = [for(i=[0:1:len(patch[0])-1]) bezier_points(subindex(patch,i), samplepts)], + pts = [ + [bpatch[0][0]], + for(j=[0:splinesteps-2]) bezier_points(subindex(bpatch,j+1), lerpn(0,1,rowcount[j])), + [last(bpatch[0])] + ], + vnf = vnf_tri_array(pts, reverse=reverse) + ) [ + vnf, + [ + subindex(pts,0), + [for(row=pts) last(row)], + pts[0], + last(pts), + ] + ] : + bot_degen ? // only bottom is degenerate + let( + result = bezier_patch_degenerate(reverse(patch), splinesteps=splinesteps, reverse=!reverse, return_edges=true) + ) + [ + result[0], + [reverse(result[1][0]), reverse(result[1][1]), (result[1][3]), (result[1][2])] + ] : + top_degen ? // only top is degenerate + let( + full_degen = len(patch)>=4 && all(select(row_degen,1,ceil(len(patch)/2-1))), + rowmax = full_degen ? count(splinesteps+1) : + [for(j=[0:splinesteps]) j<=splinesteps/2 ? 2*j : splinesteps], + bpatch = [for(i=[0:1:len(patch[0])-1]) bezier_points(subindex(patch,i), samplepts)], + pts = [ + [bpatch[0][0]], + for(j=[1:splinesteps]) bezier_points(subindex(bpatch,j), lerpn(0,1,rowmax[j]+1)) + ], + vnf = vnf_tri_array(pts, reverse=reverse) + ) [ + vnf, + [ + subindex(pts,0), + [for(row=pts) last(row)], + pts[0], + last(pts), + ] + ] : + // must have left or right degeneracy, so transpose and recurse + let( + result = bezier_patch_degenerate(transpose(patch), splinesteps=splinesteps, reverse=!reverse, return_edges=true) + ) + [result[0], + select(result[1],[2,3,0,1]) + ]; + + function _tri_count(n) = (n*(1+n))/2; diff --git a/geometry.scad b/geometry.scad index 21b4f29..2be3719 100644 --- a/geometry.scad +++ b/geometry.scad @@ -1856,10 +1856,10 @@ function reindex_polygon(reference, poly, return_error=false) = polygon_is_clockwise(reference) ? clockwise_polygon(poly) : ccw_polygon(poly), - I = [for(i=[0:N-1]) 1], + I = [for(i=reference) 1], val = [ for(k=[0:N-1]) -           [for(i=[0:N-1]) -              (reference[i]*poly[(i+k)%N]) ] ]*I, + [for(i=[0:N-1]) + (reference[i]*poly[(i+k)%N]) ] ]*I, optimal_poly = polygon_shift(fixpoly, max_index(val)) ) return_error? [optimal_poly, min(poly*(I*poly)-2*val)] : diff --git a/math.scad b/math.scad index 43b1347..052f6e0 100644 --- a/math.scad +++ b/math.scad @@ -583,7 +583,7 @@ function lcm(a,b=[]) = function sum(v, dflt=0) = v==[]? dflt : assert(is_consistent(v), "Input to sum is non-numeric or inconsistent") - is_vector(v) || is_matrix(v) ? [for(i=[1:len(v)]) 1]*v : + is_vector(v) || is_matrix(v) ? [for(i=v) 1]*v : _sum(v,v[0]*0); function _sum(v,_total,_i=0) = _i>=len(v) ? _total : _sum(v,_total+v[_i], _i+1); @@ -1155,6 +1155,18 @@ function all_nonnegative(x) = false; +// Function all_equal() +// Usage: +// b = all_equal(vec,); +// Description: +// Returns true if all of the entries in vec are equal to each other, or approximately equal to each other if eps is set. +// Arguments: +// vec = vector to check +// eps = Set to tolerance for approximate equality. Default: 0 +function all_equal(vec,eps=0) = + eps==0 ? [for(v=vec) if (v!=vec[0]) v] == [] + : [for(v=vec) if (!approx(v,vec[0])) v] == []; + // Function: approx() // Usage: // b = approx(a,b,) diff --git a/rounding.scad b/rounding.scad index b85a70c..d2493c6 100644 --- a/rounding.scad +++ b/rounding.scad @@ -1779,8 +1779,8 @@ function _rp_compute_patches(top, bot, rtop, rsides, ktop, ksides, concave) = // M = path3d(turtle(["left", 180, "length",3,"move", "left", "move", 3, "right", "move", "right", "move", 4, "right", "move", 3, "right", "move", 2])); // rounded_prism(M, apply(right(1)*scale(.75)*up(3),M), joint_top=0.5, joint_bot=0.2, joint_sides=[.2,1,1,0.5,1.5,.5,2], splinesteps=32); // Example: this example shows most of the different types of patches that rounded_prism creates. Note that some of the patches are close to interfering with each other across the top of the polyhedron, which would create an invalid result. -// N = apply(rot(180)*yscale(.8),turtle(["length",3,"left", "move", 2, "right", 135, "move", sqrt(2), "left", "move", sqrt(2), "right", 135, "move", 2])); -// rounded_prism(N, height=3, joint_bot=0.5, joint_top=1.25, joint_sides=[[1,1.75],0,.5,.5,2], debug=true); + N = apply(rot(180)*yscale(.8),turtle(["length",3,"left", "move", 2, "right", 135, "move", sqrt(2), "left", "move", sqrt(2), "right", 135, "move", 2])); + rounded_prism(N, height=3, joint_bot=0.5, joint_top=1.25, joint_sides=[[1,1.75],0,.5,.5,2], debug=true); // Example: This object has different scales on its different axies. Here is the largest symmetric rounding that fits. Note that the rounding is slightly smaller than the object dimensions because of roundoff error. // rounded_prism(square([100.1,30.1]), height=8.1, joint_top=4, joint_bot=4, joint_sides=15, k_sides=0.3, splinesteps=32); // Example: Using asymetric rounding enables a much more rounded form: @@ -1886,8 +1886,8 @@ function rounded_prism(bottom, top, joint_bot=0, joint_top=0, joint_sides=0, k_b let( // Entries in the next two lists have the form [edges, vnf] where // edges is a list [leftedge, rightedge, topedge, botedge] - top_samples = [for(patch=top_patch) bezier_patch_degenerate(patch,splinesteps,reverse=true) ], - bot_samples = [for(patch=bot_patch) bezier_patch_degenerate(patch,splinesteps,reverse=false) ], + top_samples = [for(patch=top_patch) bezier_patch_degenerate(patch,splinesteps,reverse=true,return_edges=true) ], + bot_samples = [for(patch=bot_patch) bezier_patch_degenerate(patch,splinesteps,reverse=false,return_edges=true) ], leftidx=0, rightidx=1, topidx=2, @@ -1895,14 +1895,14 @@ function rounded_prism(bottom, top, joint_bot=0, joint_top=0, joint_sides=0, k_b edge_points = [for(i=[0:N-1]) let( - top_edge = [ top_samples[i][0][rightidx], select(top_samples, i+1)[0][leftidx]], - bot_edge = [ select(bot_samples, i+1)[0][leftidx], bot_samples[i][0][rightidx]], - vert_edge = [ bot_samples[i][0][botidx], top_samples[i][0][botidx]] + top_edge = [ top_samples[i][1][rightidx], select(top_samples, i+1)[1][leftidx]], + bot_edge = [ select(bot_samples, i+1)[1][leftidx], bot_samples[i][1][rightidx]], + vert_edge = [ bot_samples[i][1][botidx], top_samples[i][1][botidx]] ) each [top_edge, bot_edge, vert_edge] ], faces = [ - [for(i=[0:N-1]) each top_samples[i][0][topidx]], - [for(i=[N-1:-1:0]) each reverse(bot_samples[i][0][topidx])], + [for(i=[0:N-1]) each top_samples[i][1][topidx]], + [for(i=[N-1:-1:0]) each reverse(bot_samples[i][1][topidx])], for(i=[0:N-1]) [ bot_patch[i][4][4], select(bot_patch,i+1)[4][0], @@ -1934,8 +1934,8 @@ function rounded_prism(bottom, top, joint_bot=0, joint_top=0, joint_sides=0, k_b "Roundovers interfere with each other on bottom face: either input is self intersecting or top joint length is too large") assert(debug || (verify_vert==[] && verify_horiz==[]), "Curvature continuity failed") let( - vnf = vnf_merge([ each subindex(top_samples,1), - each subindex(bot_samples,1), + vnf = vnf_merge([ each subindex(top_samples,0), + each subindex(bot_samples,0), for(pts=edge_points) vnf_vertex_array(pts), vnf_triangulate(vnf_add_faces(EMPTY_VNF,faces)) ]) @@ -1943,116 +1943,6 @@ function rounded_prism(bottom, top, joint_bot=0, joint_top=0, joint_sides=0, k_b debug ? [concat(top_patch, bot_patch), vnf] : vnf; -// This function takes a bezier patch as input and returns [edges, vnf], where -// edges = [leftedge, rightedge, topedge, bottomedge] -// gives the points along the edges of the patch, and the vnf is the patch vnf. -// It checks for various types of degeneracy and uses half or full triangular -// sampling on degenerate patches. - -function bezier_patch_degenerate(patch, splinesteps=16, reverse=false) = - assert(is_num(splinesteps), "splinesteps must be a number") - let( - top_degen = patch[0][0] == last(patch[0]), - bot_degen = last(patch)[0] == last(last(patch)), - left_degen = patch[0][0] == last(patch)[0], - right_degen = last(patch[0]) == last(last(patch)), - samplepts = count(splinesteps+1)/splinesteps - ) - top_degen && bot_degen && left_degen && right_degen ? // fully degenerate case - [repeat([patch[0][0]],4), EMPTY_VNF] : - top_degen && bot_degen ? // double degenerate case (top/bot) - let( - pts = bezier_points(subindex(patch,0), samplepts) - ) - [[pts,pts,[pts[0]],[last(pts)]], EMPTY_VNF] : - left_degen && right_degen ? // double degenerate case (sides) - let( - pts = bezier_points(patch[0], samplepts) - ) - [[[pts[0]], [last(pts)], pts, pts], EMPTY_VNF] : - !top_degen && !bot_degen ? // non-degenerate case - let( - k=echo("non-degenerate case"), - pts = bezier_patch_points(patch, samplepts, samplepts) - ) - [ - [subindex(pts,0), subindex(pts,len(pts)-1), pts[0], last(pts)], - vnf_vertex_array(pts, reverse=reverse) - ] : - bot_degen ? // only bottom is degenerate - let( - result = bezier_patch_degenerate(reverse(patch), splinesteps=splinesteps, reverse=!reverse) - ) - [ - [reverse(result[0][0]), reverse(result[0][1]), (result[0][3]), (result[0][2])], - result[1] - ] : - // at this point top_degen is true // only top is degenerate - let( - full_degen = patch[1][0] == last(patch[1]), - rowmax = full_degen ? count(splinesteps+1) : - [for(j=[0:splinesteps]) j<=splinesteps/2 ? 2*j : splinesteps], - vbb=echo("single degenerate case"), - bpatch = [for(i=[0:1:len(patch[0])-1]) bezier_points(subindex(patch,i), samplepts)], - pts = [ - [bpatch[0][0]], - for(j=[1:splinesteps]) bezier_points(subindex(bpatch,j), lerpn(0,1,rowmax[j]+1)) - ], - vnf = vnf_tri_array(pts, reverse=reverse) - ) [ - [ - subindex(pts,0), - [for(row=pts) last(row)], - pts[0], - last(pts), - ], - vnf - ]; - - -// This function produces a vnf with a triangulation for a list of rows -// where the number of points between rows differs by at most 2. -// It's a generalization of vnf_vertex_array. -function vnf_tri_array(points, row_wrap=false, reverse=false) = - let( - lens = [for(row=points) len(row)], - rowstarts = [0,each cumsum(lens)], - faces = - [for(i=[0:1:len(points) - 1 - (row_wrap ? 0 : 1)]) each - let( - rowstart = rowstarts[i], - nextrow = select(rowstarts,i+1), - delta = select(lens,i+1)-lens[i] - ) - delta == 0 ? - [for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow] : [j+rowstart, j+rowstart+1, j+nextrow], - for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+nextrow, j+nextrow+1] : [j+rowstart+1, j+nextrow+1, j+nextrow]] : - delta == 1 ? - [for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+1] : [j+rowstart, j+rowstart+1, j+nextrow+1], - for(j=[0:1:lens[i]-1]) reverse ? [j+rowstart, j+nextrow, j+nextrow+1] : [j+rowstart, j+nextrow+1, j+nextrow]] : - delta == -1 ? - [for(j=[0:1:lens[i]-3]) reverse ? [j+rowstart+1, j+nextrow, j+nextrow+1]: [j+rowstart+1, j+nextrow+1, j+nextrow], - for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow] : [j+rowstart, j+rowstart+1, j+nextrow]] : - let(count = floor((lens[i]-1)/2)) - delta == 2 ? - [ - for(j=[0:1:count-1]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+1] : [j+rowstart, j+rowstart+1, j+nextrow+1], // top triangles left - for(j=[count:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+2] : [j+rowstart, j+rowstart+1, j+nextrow+2], // top triangles right - for(j=[0:1:count]) reverse ? [j+rowstart, j+nextrow, j+nextrow+1] : [j+rowstart, j+nextrow+1, j+nextrow], // bot triangles left - for(j=[count+1:1:select(lens,i+1)-2]) reverse ? [j+rowstart-1, j+nextrow, j+nextrow+1] : [j+rowstart-1, j+nextrow+1, j+nextrow], // bot triangles right - ] : - delta == -2 ? - [ - for(j=[0:1:count-2]) reverse ? [j+nextrow, j+nextrow+1, j+rowstart+1] : [j+nextrow, j+rowstart+1, j+nextrow+1], - for(j=[count-1:1:lens[i]-4]) reverse ? [j+nextrow,j+nextrow+1,j+rowstart+2] : [j+nextrow,j+rowstart+2, j+nextrow+1], - for(j=[0:1:count-1]) reverse ? [j+nextrow, j+rowstart+1, j+rowstart] : [j+nextrow, j+rowstart, j+rowstart+1], - for(j=[count:1:select(lens,i+1)]) reverse ? [ j+nextrow-1, j+rowstart+1, j+rowstart]: [ j+nextrow-1, j+rowstart, j+rowstart+1], - ] : - assert(false,str("Unsupported row length difference of ",delta, " between row ",i," and ",i+1)) - ]) - [flatten(points), faces]; - - // Converts a 2d path to a path on a cylinder at radius r function _cyl_hole(r, path) = From 0dbd2f7c6907f42c57dd60777c67467e41cda796 Mon Sep 17 00:00:00 2001 From: github-actions Date: Mon, 12 Apr 2021 07:37:27 +0000 Subject: [PATCH 21/21] Bump release version. --- version.scad | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/version.scad b/version.scad index efc6eff..25bda53 100644 --- a/version.scad +++ b/version.scad @@ -6,7 +6,7 @@ ////////////////////////////////////////////////////////////////////// -BOSL_VERSION = [2,0,604]; +BOSL_VERSION = [2,0,605]; // Section: BOSL Library Version Functions