1//////////////////////////////////////////////////////////////////////
   2// LibFile: vnf.scad
   3//   The Vertices'N'Faces structure (VNF) holds the data used by polyhedron() to construct objects: a vertex
   4//   list and a list of faces.  This library makes it easier to construct polyhedra by providing
   5//   functions to construct, merge, and modify VNF data, while avoiding common pitfalls such as
   6//   reversed faces.  It can find faults in your polyhedrons.  Note that this file is for low level manipulation
   7//   of lists of vertices and faces: it can perform some simple transformations on VNF structures
   8//   but cannot perform boolean operations on the polyhedrons represented by VNFs.
   9// Includes:
  10//   include <BOSL2/std.scad>
  11// FileGroup: Advanced Modeling
  12// FileSummary: Vertices 'n' Faces structure.  Makes polyhedron() easier to use.
  13// FileFootnotes: STD=Included in std.scad
  14//////////////////////////////////////////////////////////////////////
  15
  16
  17// Section: Creating Polyhedrons with VNF Structures
  18//   VNF stands for "Vertices'N'Faces".  VNF structures are 2-item lists, `[VERTICES,FACES]` where the
  19//   first item is a list of vertex points, and the second is a list of face indices into the vertex
  20//   list.  Each VNF is self contained, with face indices referring only to its own vertex list.
  21//   You can construct a `polyhedron()` in parts by describing each part in a self-contained VNF, then
  22//   merge the various VNFs to get the completed polyhedron vertex list and faces.
  23
  24/// Constant: EMPTY_VNF
  25/// Description:
  26///   The empty VNF data structure.  Equal to `[[],[]]`.
  27EMPTY_VNF = [[],[]];  // The standard empty VNF with no vertices or faces.
  28
  29
  30// Function: vnf_vertex_array()
  31// Synopsis: Returns a VNF structure from a rectangular vertex list.
  32// SynTags: VNF
  33// Topics: VNF Generators, Lists
  34// See Also: vnf_tri_array(), vnf_join(), vnf_from_polygons(), vnf_from_region()
  35// Usage:
  36//   vnf = vnf_vertex_array(points, [caps=], [cap1=], [cap2=], [style=], [reverse=], [col_wrap=], [row_wrap=]);
  37// Description:
  38//   Creates a VNF structure from a rectangular vertex list, by dividing the vertices into columns and rows,
  39//   adding faces to tile the surface.  You can optionally have faces added to wrap the last column
  40//   back to the first column, or wrap the last row to the first.  Endcaps can be added to either
  41//   the first and/or last rows.  The style parameter determines how the quadrilaterals are divided into
  42//   triangles.  The default style is an arbitrary, systematic subdivision in the same direction.  The "alt" style
  43//   is the uniform subdivision in the other (alternate) direction.  The "min_edge" style picks the shorter edge to
  44//   subdivide for each quadrilateral, so the division may not be uniform across the shape.  The "quincunx" style
  45//   adds a vertex in the center of each quadrilateral and creates four triangles, and the "convex" and "concave" styles
  46//   chooses the locally convex/concave subdivision.  Degenerate faces
  47//   are not included in the output, but if this results in unused vertices they will still appear in the output.
  48// Arguments:
  49//   points = A list of vertices to divide into columns and rows.
  50//   ---
  51//   caps = If true, add endcap faces to the first AND last rows.
  52//   cap1 = If true, add an endcap face to the first row.
  53//   cap2 = If true, add an endcap face to the last row.
  54//   col_wrap = If true, add faces to connect the last column to the first.
  55//   row_wrap = If true, add faces to connect the last row to the first.
  56//   reverse = If true, reverse all face normals.
  57//   style = The style of subdividing the quads into faces.  Valid options are "default", "alt", "min_edge", "quincunx", "convex" and "concave".
  58// Example(3D):
  59//   vnf = vnf_vertex_array(
  60//       points=[
  61//           for (h = [0:5:180-EPSILON]) [
  62//               for (t = [0:5:360-EPSILON])
  63//                   cylindrical_to_xyz(100 + 12 * cos((h/2 + t)*6), t, h)
  64//           ]
  65//       ],
  66//       col_wrap=true, caps=true, reverse=true, style="alt"
  67//   );
  68//   vnf_polyhedron(vnf);
  69// Example(3D): Both `col_wrap` and `row_wrap` are true to make a torus.
  70//   vnf = vnf_vertex_array(
  71//       points=[
  72//           for (a=[0:5:360-EPSILON])
  73//               apply(
  74//                   zrot(a) * right(30) * xrot(90),
  75//                   path3d(circle(d=20))
  76//               )
  77//       ],
  78//       col_wrap=true, row_wrap=true, reverse=true
  79//   );
  80//   vnf_polyhedron(vnf);
  81// Example(3D): Möbius Strip.  Note that `row_wrap` is not used, and the first and last profile copies are the same.
  82//   vnf = vnf_vertex_array(
  83//       points=[
  84//           for (a=[0:5:360]) apply(
  85//               zrot(a) * right(30) * xrot(90) * zrot(a/2+60),
  86//               path3d(square([1,10], center=true))
  87//           )
  88//       ],
  89//       col_wrap=true, reverse=true
  90//   );
  91//   vnf_polyhedron(vnf);
  92// Example(3D): Assembling a Polyhedron from Multiple Parts
  93//   wall_points = [
  94//       for (a = [-90:2:90]) apply(
  95//           up(a) * scale([1-0.1*cos(a*6),1-0.1*cos((a+90)*6),1]),
  96//           path3d(circle(d=100))
  97//       )
  98//   ];
  99//   cap = [
 100//       for (a = [0:0.01:1+EPSILON]) apply(
 101//           up(90-5*sin(a*360*2)) * scale([a,a,1]),
 102//           wall_points[0]
 103//       )
 104//   ];
 105//   cap1 = [for (p=cap) down(90, p=zscale(-1, p=p))];
 106//   cap2 = [for (p=cap) up(90, p=p)];
 107//   vnf1 = vnf_vertex_array(points=wall_points, col_wrap=true);
 108//   vnf2 = vnf_vertex_array(points=cap1, col_wrap=true);
 109//   vnf3 = vnf_vertex_array(points=cap2, col_wrap=true, reverse=true);
 110//   vnf_polyhedron([vnf1, vnf2, vnf3]);
 111// Example(3D): Building a Multi-Stage Cylindrical Ramp
 112//   include <BOSL2/rounding.scad>
 113//   major_r = 50;
 114//   groove_profile = [
 115//       [-10,0], each arc(points=[[-7,0],[0,-3],[7,0]]), [10,0]
 116//   ];
 117//   ramp_profile = [ [-10,25], [90,25], [180,5], [190,5] ];
 118//   rgroove = apply(right(major_r) * xrot(90), path3d(groove_profile));
 119//   rprofile = round_corners(ramp_profile, radius=20, closed=false, $fn=72);
 120//   vnf = vnf_vertex_array([
 121//       for (a = [ramp_profile[0].x : 1 : last(ramp_profile).x]) let(
 122//           z = lookup(a,rprofile),
 123//           m = zrot(a) * up(z)
 124//       )
 125//       apply(m, [ [rgroove[0].x,0,-z], each rgroove, [last(rgroove).x,0,-z] ])
 126//   ], caps=true, col_wrap=true, reverse=true);
 127//   vnf_polyhedron(vnf, convexity=8);
 128function vnf_vertex_array(
 129    points,
 130    caps, cap1, cap2,
 131    col_wrap=false,
 132    row_wrap=false,
 133    reverse=false,
 134    style="default"
 135) =
 136    assert(!(any([caps,cap1,cap2]) && !col_wrap), "col_wrap must be true if caps are requested")
 137    assert(!(any([caps,cap1,cap2]) && row_wrap), "Cannot combine caps with row_wrap")
 138    assert(in_list(style,["default","alt","quincunx", "convex","concave", "min_edge","min_area"]))
 139    assert(is_matrix(points[0], n=3),"Point array has the wrong shape or points are not 3d")
 140    assert(is_consistent(points), "Non-rectangular or invalid point array")
 141    let(
 142        pts = flatten(points),
 143        pcnt = len(pts),
 144        rows = len(points),
 145        cols = len(points[0])
 146    )
 147    rows<=1 || cols<=1 ? EMPTY_VNF :
 148    let(
 149        cap1 = first_defined([cap1,caps,false]),
 150        cap2 = first_defined([cap2,caps,false]),
 151        colcnt = cols - (col_wrap?0:1),
 152        rowcnt = rows - (row_wrap?0:1),
 153        verts = [
 154            each pts,
 155            if (style=="quincunx")
 156                for (r = [0:1:rowcnt-1], c = [0:1:colcnt-1])
 157                   let(
 158                       i1 = ((r+0)%rows)*cols + ((c+0)%cols),
 159                       i2 = ((r+1)%rows)*cols + ((c+0)%cols),
 160                       i3 = ((r+1)%rows)*cols + ((c+1)%cols),
 161                       i4 = ((r+0)%rows)*cols + ((c+1)%cols)
 162                   )
 163                   mean([pts[i1], pts[i2], pts[i3], pts[i4]])
 164        ],
 165        allfaces = [
 166            if (cap1) count(cols,reverse=!reverse),
 167            if (cap2) count(cols,(rows-1)*cols, reverse=reverse),
 168            for (r = [0:1:rowcnt-1], c=[0:1:colcnt-1])
 169               each
 170               let(
 171                   i1 = ((r+0)%rows)*cols + ((c+0)%cols),
 172                   i2 = ((r+1)%rows)*cols + ((c+0)%cols),
 173                   i3 = ((r+1)%rows)*cols + ((c+1)%cols),
 174                   i4 = ((r+0)%rows)*cols + ((c+1)%cols),
 175                   faces =
 176                        style=="quincunx"?
 177                          let(i5 = pcnt + r*colcnt + c)
 178                          [[i1,i5,i2],[i2,i5,i3],[i3,i5,i4],[i4,i5,i1]]
 179                      : style=="alt"?
 180                          [[i1,i4,i2],[i2,i4,i3]]
 181                      : style=="min_area"?
 182                          let(
 183                               area42 = norm(cross(pts[i2]-pts[i1], pts[i4]-pts[i1]))+norm(cross(pts[i4]-pts[i3], pts[i2]-pts[i3])),
 184                               area13 = norm(cross(pts[i1]-pts[i4], pts[i3]-pts[i4]))+norm(cross(pts[i3]-pts[i2], pts[i1]-pts[i2])),
 185                               minarea_edge = area42 < area13 + EPSILON
 186                                 ? [[i1,i4,i2],[i2,i4,i3]]
 187                                 : [[i1,i3,i2],[i1,i4,i3]]
 188                          )
 189                          minarea_edge
 190                      : style=="min_edge"?
 191                          let(
 192                               d42=norm(pts[i4]-pts[i2]),
 193                               d13=norm(pts[i1]-pts[i3]),
 194                               shortedge = d42<d13+EPSILON
 195                                 ? [[i1,i4,i2],[i2,i4,i3]]
 196                                 : [[i1,i3,i2],[i1,i4,i3]]
 197                          )
 198                          shortedge
 199                      : style=="convex"?
 200                          let(   // Find normal for 3 of the points.  Is the other point above or below?
 201                              n = (reverse?-1:1)*cross(pts[i2]-pts[i1],pts[i3]-pts[i1]),
 202                              convexfaces = n==0
 203                                ? [[i1,i4,i3]]
 204                                : n*pts[i4] > n*pts[i1]
 205                                    ? [[i1,i4,i2],[i2,i4,i3]]
 206                                    : [[i1,i3,i2],[i1,i4,i3]]
 207                          )
 208                          convexfaces
 209                      : style=="concave"?
 210                          let(   // Find normal for 3 of the points.  Is the other point above or below?
 211                              n = (reverse?-1:1)*cross(pts[i2]-pts[i1],pts[i3]-pts[i1]),
 212                              concavefaces = n==0
 213                                ? [[i1,i4,i3]]
 214                                : n*pts[i4] <= n*pts[i1]
 215                                    ? [[i1,i4,i2],[i2,i4,i3]]
 216                                    : [[i1,i3,i2],[i1,i4,i3]]
 217                          )
 218                          concavefaces
 219                      : [[i1,i3,i2],[i1,i4,i3]],
 220                   // remove degenerate faces
 221                   culled_faces= [for(face=faces)
 222                       if (norm(cross(verts[face[1]]-verts[face[0]],
 223                                      verts[face[2]]-verts[face[0]]))>EPSILON)
 224                           face
 225                   ],
 226                   rfaces = reverse? [for (face=culled_faces) reverse(face)] : culled_faces
 227               )
 228               rfaces,
 229        ]
 230    )
 231    [verts,allfaces];
 232
 233
 234// Function: vnf_tri_array()
 235// Synopsis: Returns a VNF from an array of points.
 236// SynTags: VNF
 237// Topics: VNF Generators, Lists
 238// See Also: vnf_vertex_array(), vnf_join(), vnf_from_polygons(), vnf_from_region()
 239// Usage:
 240//   vnf = vnf_tri_array(points, [row_wrap], [reverse])
 241// Description:
 242//   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
 243//   the construction of triangular VNF patches.  The resulting VNF can be wrapped along the rows by setting `row_wrap` to true.
 244//   You cannot wrap columns: if you need to do that you'll need to merge two VNF arrays that share edges.  Degenerate faces
 245//   are not included in the output, but if this results in unused vertices they will still appear in the output.
 246// Arguments:
 247//   points = List of point lists for each row
 248//   row_wrap = If true then add faces connecting the first row and last row.  These rows must differ by at most 2 in length.
 249//   reverse = Set this to reverse the direction of the faces
 250// Example(3D,NoAxes): Each row has one more point than the preceeding one.
 251//   pts = [for(y=[1:1:10]) [for(x=[0:y-1]) [x,y,y]]];
 252//   vnf = vnf_tri_array(pts);
 253//   vnf_wireframe(vnf,width=0.1);
 254//   color("red")move_copies(flatten(pts)) sphere(r=.15,$fn=9);
 255// Example(3D,NoAxes): Each row has two more points than the preceeding one.
 256//   pts = [for(y=[0:2:10]) [for(x=[-y/2:y/2]) [x,y,y]]];
 257//   vnf = vnf_tri_array(pts);
 258//   vnf_wireframe(vnf,width=0.1);
 259//   color("red")move_copies(flatten(pts)) sphere(r=.15,$fn=9);
 260// Example(3D): Merging two VNFs to construct a cone with one point length change between rows.
 261//   pts1 = [for(z=[0:10]) path3d(arc(3+z,r=z/2+1, angle=[0,180]),10-z)];
 262//   pts2 = [for(z=[0:10]) path3d(arc(3+z,r=z/2+1, angle=[180,360]),10-z)];
 263//   vnf = vnf_join([vnf_tri_array(pts1),
 264//                     vnf_tri_array(pts2)]);
 265//   color("green")vnf_wireframe(vnf,width=0.1);
 266//   vnf_polyhedron(vnf);
 267// Example(3D): Cone with length change two between rows
 268//   pts1 = [for(z=[0:1:10]) path3d(arc(3+2*z,r=z/2+1, angle=[0,180]),10-z)];
 269//   pts2 = [for(z=[0:1:10]) path3d(arc(3+2*z,r=z/2+1, angle=[180,360]),10-z)];
 270//   vnf = vnf_join([vnf_tri_array(pts1),
 271//                    vnf_tri_array(pts2)]);
 272//   color("green")vnf_wireframe(vnf,width=0.1);
 273//   vnf_polyhedron(vnf);
 274// Example(3D,NoAxes): Point count can change irregularly
 275//   lens = [10,9,7,5,6,8,8,10];
 276//   pts = [for(y=idx(lens)) lerpn([-lens[y],y,y],[lens[y],y,y],lens[y])];
 277//   vnf = vnf_tri_array(pts);
 278//   vnf_wireframe(vnf,width=0.1);
 279//   color("red")move_copies(flatten(pts)) sphere(r=.15,$fn=9);
 280function vnf_tri_array(points, row_wrap=false, reverse=false) =
 281    let(
 282       lens = [for(row=points) len(row)],
 283       rowstarts = [0,each cumsum(lens)],
 284       faces =
 285          [for(i=[0:1:len(points) - 1 - (row_wrap ? 0 : 1)]) each
 286            let(
 287                rowstart = rowstarts[i],
 288                nextrow = select(rowstarts,i+1),
 289                delta = select(lens,i+1)-lens[i]
 290            )
 291            delta == 0 ?
 292              [for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow] : [j+rowstart, j+rowstart+1, j+nextrow],
 293               for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+nextrow, j+nextrow+1] : [j+rowstart+1, j+nextrow+1, j+nextrow]] :
 294            delta == 1 ?
 295              [for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow+1] : [j+rowstart, j+rowstart+1, j+nextrow+1],
 296               for(j=[0:1:lens[i]-1]) reverse ? [j+rowstart, j+nextrow, j+nextrow+1] : [j+rowstart, j+nextrow+1, j+nextrow]] :
 297            delta == -1 ?
 298              [for(j=[0:1:lens[i]-3]) reverse ? [j+rowstart+1, j+nextrow, j+nextrow+1]: [j+rowstart+1, j+nextrow+1, j+nextrow],
 299               for(j=[0:1:lens[i]-2]) reverse ? [j+rowstart+1, j+rowstart, j+nextrow] : [j+rowstart, j+rowstart+1, j+nextrow]] :
 300            let(count = floor((lens[i]-1)/2))
 301            delta == 2 ?
 302              [
 303               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
 304               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
 305               for(j=[0:1:count]) reverse ? [j+rowstart, j+nextrow, j+nextrow+1] : [j+rowstart, j+nextrow+1, j+nextrow],                        // bot triangles left
 306               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
 307              ] :
 308            delta == -2 ?
 309              [
 310               for(j=[0:1:count-2]) reverse ? [j+nextrow, j+nextrow+1, j+rowstart+1] : [j+nextrow, j+rowstart+1, j+nextrow+1],
 311               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],
 312               for(j=[0:1:count-1]) reverse ? [j+nextrow, j+rowstart+1, j+rowstart] : [j+nextrow, j+rowstart, j+rowstart+1],
 313               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],
 314              ] :
 315            assert(false,str("Unsupported row length difference of ",delta, " between row ",i," and ",(i+1)%len(points)))
 316          ],
 317       verts = flatten(points),
 318       culled_faces=
 319           [for(face=faces)
 320               if (norm(verts[face[0]]-verts[face[1]])>EPSILON &&
 321                   norm(verts[face[1]]-verts[face[2]])>EPSILON &&
 322                   norm(verts[face[2]]-verts[face[0]])>EPSILON)
 323                   face
 324           ]
 325    )
 326    [flatten(points), culled_faces];
 327
 328
 329
 330// Function: vnf_join()
 331// Synopsis: Returns a single VNF structure from a list of VNF structures.
 332// SynTags: VNF
 333// Topics: VNF Generators, Lists
 334// See Also: vnf_tri_array(), vnf_vertex_array(), vnf_from_polygons(), vnf_from_region()
 335// Usage:
 336//   vnf = vnf_join([VNF, VNF, VNF, ...]);
 337// Description:
 338//   Given a list of VNF structures, merges them all into a single VNF structure.
 339//   Combines all the points of the input VNFs and labels the faces appropriately.
 340//   All the points in the input VNFs will appear in the output, even if they are
 341//   duplicates of each other.  It is valid to repeat points in a VNF, but if you
 342//   with to remove the duplicates that will occur along joined edges, use {{vnf_merge_points()}}.
 343//   .
 344//   Note that this is a tool for manipulating polyhedron data.  It is for
 345//   building up a full polyhedron from partial polyhedra.
 346//   It is *not* a union operator for VNFs.  The VNFs to be joined must not intersect each other,
 347//   except at edges, or the result will be an invalid polyhedron.  Similarly the
 348//   result must not have any other illegal polyhedron characteristics, such as creating
 349//   more than two faces sharing the same edge.
 350//   If you want a valid result it is your responsibility to ensure that the polyhedron
 351//   has no holes, no intersecting faces or edges, and obeys all the requirements
 352//   that CGAL expects.
 353//   .
 354//   For example, if you combine two pyramids to try to make an octahedron, the result will
 355//   be invalid because of the two internal faces created by the pyramid bases.  A valid
 356//   use would be to build a cube missing one face and a pyramid missing its base and
 357//   then join them into a cube with a point.
 358// Arguments:
 359//   vnfs = a list of the VNFs to joint into one VNF.
 360// Example(3D,VPR=[60,0,26],VPD=55,VPT=[5.6,-5.3,9.8]): Here is a VNF where the top face is missing.  It is not a valid polyhedron like this, but we can use it as a building block to make a polyhedron.
 361//   bottom = vnf_vertex_array([path3d(rect(8)), path3d(rect(5),4)],col_wrap=true,cap1=true);
 362//   vnf_polyhedron(bottom);
 363// Example(3D,VPR=[60,0,26],VPD=55,VPT=[5.6,-5.3,9.8]): Here is a VNF that also has a missing face.
 364//   triangle = yrot(-90,path3d(regular_ngon(n=3,side=5,anchor=LEFT)));
 365//   top = up(4,vnf_vertex_array([list_set(right(2.5,triangle),0,[0,0,7]),
 366//                               right(6,triangle)
 367//                               ], col_wrap=true, cap2=true));
 368//   vnf_polyhedron(zrot(90,top));
 369// Example(3D,VPR=[60,0,26],VPD=55,VPT=[5.6,-5.3,9.8]): Using vnf_join combines the two VNFs into a single VNF.  Note that they share an edge.  But the result still isn't closed, so it is not yet a valid polyhedron.
 370//   bottom = vnf_vertex_array([path3d(rect(8)), path3d(rect(5),4)],col_wrap=true,cap1=true);
 371//   triangle = yrot(-90,path3d(regular_ngon(n=3,side=5,anchor=LEFT)));
 372//   top = up(4,vnf_vertex_array([list_set(right(2.5,triangle),0,[0,0,7]),
 373//                                right(6,triangle)
 374//                               ], col_wrap=true, cap2=true));
 375//   full = vnf_join([bottom,zrot(90,top)]);
 376//   vnf_polyhedron(full);
 377// Example(3D,VPR=[60,0,26],VPD=55,VPT=[5.6,-5.3,9.8]): If we add enough pieces, and the pieces are all consistent with each other, then we can arrive at a valid polyhedron like this one.  To be valid you need to meet all the CGAL requirements: every edge has exactly two faces, all faces are in clockwise order, no intersections of edges.
 378//   bottom = vnf_vertex_array([path3d(rect(8)), path3d(rect(5),4)],col_wrap=true,cap1=true);
 379//   triangle = yrot(-90,path3d(regular_ngon(n=3,side=5,anchor=LEFT)));
 380//   top = up(4,vnf_vertex_array([list_set(right(2.5,triangle),0,[0,0,7]),
 381//                                right(6,triangle)
 382//                               ], col_wrap=true, cap2=true));
 383//   full = vnf_join([bottom,
 384//                     for(theta=[0:90:359]) zrot(theta,top)
 385//                    ]);
 386//   vnf_polyhedron(full);
 387// Example(3D): The vnf_join function is not a union operator for polyhedra.  If any faces intersect, like they do in this example where we combine the faces of two cubes, the result is invalid and will give rise to CGAL errors when you add more objects into the model.
 388//   cube1 = cube(5);
 389//   cube2 = move([2,2,2],cube1);
 390//   badvnf = vnf_join([cube1,cube2]);
 391//   vnf_polyhedron(badvnf);
 392//   right(2.5)up(3)color("red")
 393//         text3d("Invalid",size=1,anchor=CENTER,
 394//         orient=FRONT,h=.1);
 395function vnf_join(vnfs) =
 396    assert(is_vnf_list(vnfs) , "Input must be a list of VNFs")
 397    len(vnfs)==1 ? vnfs[0]
 398    :
 399    let (
 400        offs  = cumsum([ 0, for (vnf = vnfs) len(vnf[0]) ]),
 401        verts = [for (vnf=vnfs) each vnf[0]],
 402        faces =
 403            [ for (i = idx(vnfs))
 404                let( faces = vnfs[i][1] )
 405                for (face = faces)
 406                    if ( len(face) >= 3 )
 407                        [ for (j = face)
 408                            assert( j>=0 && j<len(vnfs[i][0]),
 409                                    str("VNF number ", i, " has a face indexing an nonexistent vertex") )
 410                            offs[i] + j ]
 411            ]
 412    )
 413    [verts,faces];
 414
 415
 416
 417// Function: vnf_from_polygons()
 418// Synopsis: Returns a VNF from a list of 3D polygons.
 419// SynTags: VNF
 420// Topics: VNF Generators, Lists
 421// See Also: vnf_tri_array(), vnf_join(), vnf_vertex_array(), vnf_from_region()
 422// Usage:
 423//   vnf = vnf_from_polygons(polygons);
 424// Description:
 425//   Given a list of 3D polygons, produces a VNF containing those polygons.
 426//   It is up to the caller to make sure that the points are in the correct order to make the face
 427//   normals point outwards.  No checking for duplicate vertices is done.  If you want to
 428//   remove duplicate vertices use {{vnf_merge_points()}}.
 429// Arguments:
 430//   polygons = The list of 3D polygons to turn into a VNF
 431function vnf_from_polygons(polygons) =
 432   assert(is_list(polygons) && is_path(polygons[0]),"Input should be a list of polygons")
 433   let(
 434       offs = cumsum([0, for(p=polygons) len(p)]),
 435       faces = [for(i=idx(polygons))
 436                  [for (j=idx(polygons[i])) offs[i]+j]
 437               ]
 438   )
 439   [flatten(polygons), faces];
 440
 441
 442
 443
 444function _path_path_closest_vertices(path1,path2) =
 445    let(
 446        dists = [for (i=idx(path1)) let(j=closest_point(path1[i],path2)) [j,norm(path2[j]-path1[i])]],
 447        i1 = min_index(column(dists,1)),
 448        i2 = dists[i1][0]
 449    ) [dists[i1][1], i1, i2];
 450
 451
 452function _join_paths_at_vertices(path1,path2,v1,v2) =
 453    let(
 454        repeat_start = !approx(path1[v1],path2[v2]),
 455        path1 = clockwise_polygon(list_rotate(path1,v1)),
 456        path2 = ccw_polygon(list_rotate(path2,v2))
 457    )
 458    [
 459        each path1,
 460        if (repeat_start) path1[0],
 461        each path2,
 462        if (repeat_start) path2[0],
 463    ];
 464
 465
 466/// Internal Function: _cleave_connected_region(region, eps)
 467/// Description:
 468///   Given a region that is connected and has its outer border in region[0],
 469///   produces a overlapping connected path to join internal holes to
 470///   the outer border without adding points. Output is a single non-simple polygon.
 471/// Requirements:
 472///   It expects that all region paths be simple closed paths, with region[0] CW and
 473///   the other paths CCW and encircled by region[0]. The input region paths are also
 474///   supposed to be disjoint except for common vertices and common edges but with
 475///   no crossings. It may return `undef` if these conditions are not met.
 476///   This function implements an extension of the algorithm discussed in:
 477///   https://www.geometrictools.com/Documentation/TriangulationByEarClipping.pdf
 478function _cleave_connected_region(region, eps=EPSILON) =
 479    len(region)==1 ? region[0] :
 480    let(
 481        outer   = deduplicate(region[0]),             //
 482        holes   = [for(i=[1:1:len(region)-1])         // deduplication possibly unneeded
 483                      deduplicate( region[i] ) ],     //
 484        extridx = [for(li=holes) max_index(column(li,0)) ],
 485        // the right extreme vertex for each hole sorted by decreasing x values
 486        extremes = sort( [for(i=idx(holes)) [ i, extridx[i], -holes[i][extridx[i]].x] ], idx=2 )
 487    )
 488    _polyHoles(outer, holes, extremes, eps, 0);
 489
 490
 491// connect the hole paths one at a time to the outer path.
 492// 'extremes' is the list of the right extreme vertex of each hole sorted by decreasing abscissas
 493// see: _cleave_connected_region(region, eps)
 494function _polyHoles(outer, holes, extremes, eps=EPSILON, n=0) =
 495    let(
 496        extr = extremes[n],    //
 497        hole = holes[extr[0]], // hole path to bridge to the outer path
 498        ipt  = extr[1],        // index of the hole point with maximum abscissa
 499        brdg = _bridge(hole[ipt], outer, eps)  // the index of a point in outer to bridge hole[ipt] to
 500    )
 501    brdg == undef ? undef :
 502    let(
 503        l  = len(outer),
 504        lh = len(hole),
 505        // the new outer polygon bridging the hole to the old outer
 506        npoly =
 507            approx(outer[brdg], hole[ipt], eps)
 508            ?   [ for(i=[brdg:  1: brdg+l])   outer[i%l] ,
 509                  for(i=[ipt+1: 1: ipt+lh-1]) hole[i%lh] ]
 510            :   [ for(i=[brdg:  1: brdg+l])   outer[i%l] ,
 511                  for(i=[ipt:   1: ipt+lh])   hole[i%lh] ]
 512    )
 513    n==len(holes)-1 ?  npoly :
 514    _polyHoles(npoly, holes, extremes, eps, n+1);
 515
 516// find a point in outer to be connected to pt in the interior of outer
 517// by a segment that not cross or touch any non adjacente edge of outer.
 518// return the index of a vertex in the outer path where the bridge should end
 519// see _polyHoles(outer, holes, extremes, eps)
 520function _bridge(pt, outer,eps) =
 521    // find the intersection of a ray from pt to the right
 522    // with the boundary of the outer cycle
 523    let(
 524        l    = len(outer),
 525        crxs =
 526            let( edges = pair(outer,wrap=true) )
 527            [for( i = idx(edges) )
 528                let( edge = edges[i] )
 529                // consider just descending outer edges at right of pt crossing ordinate pt.y
 530                if(    (edge[0].y >  pt.y) //+eps)
 531                    && (edge[1].y <= pt.y)
 532                    && _is_at_left(pt, [edge[1], edge[0]], eps) )
 533                    [ i,
 534                      // the point of edge with ordinate pt.y
 535                      abs(pt.y-edge[1].y)<eps ? edge[1] :
 536                      let( u = (pt-edge[1]).y / (edge[0]-edge[1]).y )
 537                      (1-u)*edge[1] + u*edge[0]
 538                    ]
 539             ]
 540    )
 541    crxs == [] ? undef :
 542    let(
 543        // the intersection point of the nearest edge to pt with minimum slope
 544        minX    = min([for(p=crxs) p[1].x]),
 545        crxcand = [for(crx=crxs) if(crx[1].x < minX+eps) crx ], // nearest edges
 546        nearest = min_index([for(crx=crxcand)
 547                                (outer[crx[0]].x - pt.x) / (outer[crx[0]].y - pt.y) ]), // minimum slope
 548        proj    = crxcand[nearest],
 549        vert0   = outer[proj[0]],    // the two vertices of the nearest crossing edge
 550        vert1   = outer[(proj[0]+1)%l],
 551        isect   = proj[1]            // the intersection point
 552    )
 553    norm(pt-vert1) < eps ? (proj[0]+1)%l : // if pt touches an outer vertex, return its index
 554    // as vert0.y > pt.y then pt!=vert0
 555    norm(pt-isect) < eps ? undef :         // if pt touches the middle of an outer edge -> error
 556    let(
 557        // the edge [vert0, vert1] necessarily satisfies vert0.y > vert1.y
 558        // indices of candidates to an outer bridge point
 559        cand  =
 560            (vert0.x > pt.x)
 561            ?   [ proj[0],
 562                  // select reflex vertices inside of the triangle [pt, vert0, isect]
 563                  for(i=idx(outer))
 564                      if( _tri_class(select(outer,i-1,i+1),eps) <= 0
 565                          && _pt_in_tri(outer[i], [pt, vert0, isect], eps)>=0 )
 566                        i
 567                ]
 568            :   [ (proj[0]+1)%l,
 569                  // select reflex vertices inside of the triangle [pt, isect, vert1]
 570                  for(i=idx(outer))
 571                      if( _tri_class(select(outer,i-1,i+1),eps) <= 0
 572                          &&  _pt_in_tri(outer[i], [pt, isect, vert1], eps)>=0 )
 573                        i
 574                ],
 575        // choose the candidate outer[i] such that the line [pt, outer[i]] has minimum slope
 576        // among those with minimum slope choose the nearest to pt
 577        slopes  = [for(i=cand) 1-abs(outer[i].x-pt.x)/norm(outer[i]-pt) ],
 578        min_slp = min(slopes),
 579        cand2   = [for(i=idx(cand)) if(slopes[i]<=min_slp+eps) cand[i] ],
 580        nearest = min_index([for(i=cand2) norm(pt-outer[i]) ])
 581    )
 582    cand2[nearest];
 583
 584
 585// Function: vnf_from_region()
 586// Synopsis: Returns a 3D VNF given a 2D region.
 587// SynTags: VNF
 588// Topics: VNF Generators, Lists
 589// See Also: vnf_vertex_array(), vnf_tri_array(), vnf_join(), vnf_from_polygons()
 590// Usage:
 591//   vnf = vnf_from_region(region, [transform], [reverse]);
 592// Description:
 593//   Given a (two-dimensional) region, applies the given transformation matrix to it and makes a (three-dimensional) triangulated VNF of
 594//   faces for that region, reversed if desired.
 595// Arguments:
 596//   region = The region to convert to a vnf.
 597//   transform = If given, a transformation matrix to apply to the faces generated from the region.  Default: No transformation applied.
 598//   reverse = If true, reverse the normals of the faces generated from the region.  An untransformed region will have face normals pointing `UP`.  Default: false
 599// Example(3D):
 600//   region = [square([20,10],center=true),
 601//             right(5,square(4,center=true)),
 602//             left(5,square(6,center=true))];
 603//   vnf = vnf_from_region(region);
 604//   color("gray")down(.125)
 605//        linear_extrude(height=.125)region(region);
 606//   vnf_wireframe(vnf,width=.25);
 607function vnf_from_region(region, transform, reverse=false) =
 608    let (
 609        region = [for (path = region) deduplicate(path, closed=true)],
 610        regions = region_parts(force_region(region)),
 611        vnfs =
 612            [
 613                for (rgn = regions)
 614                let(
 615                    cleaved = path3d(_cleave_connected_region(rgn))
 616                )
 617                assert( cleaved, "The region is invalid")
 618                let(
 619                    face = is_undef(transform)? cleaved : apply(transform,cleaved),
 620                    faceidxs = reverse? [for (i=[len(face)-1:-1:0]) i] : [for (i=[0:1:len(face)-1]) i]
 621                ) [face, [faceidxs]]
 622            ],
 623        outvnf = vnf_join(vnfs)
 624    )
 625    vnf_triangulate(outvnf);
 626
 627
 628
 629// Section: VNF Testing and Access
 630
 631
 632// Function: is_vnf()
 633// Synopsis: Returns true given a VNF-like structure.
 634// Topics: VNF Manipulation
 635// See Also: is_vnf_list(), vnf_vertices(), vnf_faces()
 636// Usage:
 637//   bool = is_vnf(x);
 638// Description:
 639//   Returns true if the given value looks like a VNF structure.
 640function is_vnf(x) =
 641    is_list(x) &&
 642    len(x)==2 &&
 643    is_list(x[0]) &&
 644    is_list(x[1]) &&
 645    (x[0]==[] || (len(x[0])>=3 && is_vector(x[0][0],3))) &&
 646    (x[1]==[] || is_vector(x[1][0]));
 647
 648
 649// Function: is_vnf_list()
 650// Synopsis: Returns true given a list of VNF-like structures.
 651// Topics: VNF Manipulation
 652// See Also: is_vnf(), vnf_vertices(), vnf_faces()
 653//
 654// Description: Returns true if the given value looks passingly like a list of VNF structures.
 655function is_vnf_list(x) = is_list(x) && all([for (v=x) is_vnf(v)]);
 656
 657
 658// Function: vnf_vertices()
 659// Synopsis: Returns the list of vertex points from a VNF.
 660// Topics: VNF Manipulation
 661// See Also: is_vnf(), is_vnf_list(), vnf_faces()
 662// Description: Given a VNF structure, returns the list of vertex points.
 663function vnf_vertices(vnf) = vnf[0];
 664
 665
 666// Function: vnf_faces()
 667// Synopsis: Returns the list of faces from a VNF.
 668// Topics: VNF Manipulation
 669// See Also: is_vnf(), is_vnf_list(), vnf_vertices()
 670// Description: Given a VNF structure, returns the list of faces, where each face is a list of indices into the VNF vertex list.
 671function vnf_faces(vnf) = vnf[1];
 672
 673
 674
 675// Section: Altering the VNF Internals
 676
 677
 678// Function: vnf_reverse_faces()
 679// Synopsis: Reverses the faces of a VNF.
 680// SynTags: VNF
 681// Topics: VNF Manipulation
 682// See Also: vnf_quantize(), vnf_merge_points(), vnf_drop_unused_points(), vnf_triangulate(), vnf_slice() 
 683// Usage:
 684//   rvnf = vnf_reverse_faces(vnf);
 685// Description:
 686//   Reverses the orientation of all the faces in the given VNF.
 687function vnf_reverse_faces(vnf) =
 688    [vnf[0], [for (face=vnf[1]) reverse(face)]];
 689
 690
 691// Function: vnf_quantize()
 692// Synopsis: Quantizes the vertex coordinates of a VNF.
 693// SynTags: VNF
 694// Topics: VNF Manipulation
 695// See Also: vnf_reverse_faces(), vnf_merge_points(), vnf_drop_unused_points(), vnf_triangulate(), vnf_slice() 
 696// Usage:
 697//   vnf2 = vnf_quantize(vnf,[q]);
 698// Description:
 699//   Quantizes the vertex coordinates of the VNF to the given quanta `q`.
 700// Arguments:
 701//   vnf = The VNF to quantize.
 702//   q = The quanta to quantize the VNF coordinates to.
 703function vnf_quantize(vnf,q=pow(2,-12)) =
 704    [[for (pt = vnf[0]) quant(pt,q)], vnf[1]];
 705
 706
 707
 708// Function: vnf_merge_points()
 709// Synopsis: Consolidates duplicate vertices of a VNF.
 710// SynTags: VNF
 711// Topics: VNF Manipulation
 712// See Also: vnf_reverse_faces(), vnf_quantize(), vnf_drop_unused_points(), vnf_triangulate(), vnf_slice() 
 713// Usage:
 714//   new_vnf = vnf_merge_points(vnf, [eps]);
 715// Description:
 716//   Given a VNF, consolidates all duplicate vertices with a tolerance `eps`, relabeling the faces as necessary,
 717//   and eliminating any face with fewer than 3 vertices.  Unreferenced vertices of the input VNF are not dropped.
 718//   To remove such vertices uses {{vnf_drop_unused_points()}}.
 719// Arguments:
 720//   vnf = a VNF to consolidate
 721//   eps = the tolerance in finding duplicates. Default: EPSILON
 722function vnf_merge_points(vnf,eps=EPSILON) =
 723    let(
 724        verts = vnf[0],
 725        dedup  = vector_search(verts,eps,verts),                 // collect vertex duplicates
 726        map    = [for(i=idx(verts)) min(dedup[i]) ],             // remap duplic vertices
 727        offset = cumsum([for(i=idx(verts)) map[i]==i ? 0 : 1 ]), // remaping face vertex offsets
 728        map2   = list(idx(verts))-offset,                        // map old vertex indices to new indices
 729        nverts = [for(i=idx(verts)) if(map[i]==i) verts[i] ],    // this doesn't eliminate unreferenced vertices
 730        nfaces =
 731            [ for(face=vnf[1])
 732                let(
 733                    nface = [ for(vi=face) map2[map[vi]] ],
 734                    dface = [for (i=idx(nface))
 735                                if( nface[i]!=nface[(i+1)%len(nface)])
 736                                    nface[i] ]
 737                )
 738                if(len(dface) >= 3) dface
 739            ]
 740    )
 741    [nverts, nfaces];
 742
 743
 744// Function: vnf_drop_unused_points()
 745// Synopsis: Removes unreferenced vertices from a VNF.
 746// SynTags: VNF
 747// Topics: VNF Manipulation
 748// See Also: vnf_reverse_faces(), vnf_quantize(), vnf_merge_points(), vnf_triangulate(), vnf_slice() 
 749// Usage:
 750//   clean_vnf = vnf_drop_unused_points(vnf);
 751// Description:
 752//   Remove all unreferenced vertices from a VNF.  Note that in most
 753//   cases unreferenced vertices cause no harm, and this function may
 754//   be slow on large VNFs.
 755function vnf_drop_unused_points(vnf) =
 756    let(
 757        flat = flatten(vnf[1]),
 758        ind  = _link_indicator(flat,0,len(vnf[0])-1),
 759        verts = [for(i=idx(vnf[0])) if(ind[i]==1) vnf[0][i] ],
 760        map   = cumsum(ind)
 761    )
 762    [ verts, [for(face=vnf[1]) [for(v=face) map[v]-1 ] ] ];
 763
 764function _link_indicator(l,imin,imax) =
 765    len(l) == 0  ? repeat(imax-imin+1,0) :
 766    imax-imin<100 || len(l)<400 ? [for(si=search(list([imin:1:imax]),l,1)) si!=[] ? 1: 0 ] :
 767    let(
 768        pivot   = floor((imax+imin)/2),
 769        lesser  = [ for(li=l) if( li< pivot) li ],
 770        greater = [ for(li=l) if( li> pivot) li ]
 771    )
 772    concat( _link_indicator(lesser ,imin,pivot-1),
 773            search(pivot,l,1) ? 1 : 0 ,
 774            _link_indicator(greater,pivot+1,imax) ) ;
 775
 776// Function: vnf_triangulate()
 777// Synopsis: Triangulates the faces of a VNF.
 778// SynTags: VNF
 779// Topics: VNF Manipulation
 780// See Also: vnf_reverse_faces(), vnf_quantize(), vnf_merge_points(), vnf_drop_unused_points(), vnf_slice() 
 781// Usage:
 782//   vnf2 = vnf_triangulate(vnf);
 783// Description:
 784//   Triangulates faces in the VNF that have more than 3 vertices.
 785// Arguments:
 786//   vnf = VNF to triangulate
 787// Example(3D):
 788//   include <BOSL2/polyhedra.scad>
 789//   vnf = zrot(33,regular_polyhedron_info("vnf", "dodecahedron", side=12));
 790//   vnf_polyhedron(vnf);
 791//   triangulated = vnf_triangulate(vnf);
 792//   color("red")vnf_wireframe(triangulated,width=.3);
 793function vnf_triangulate(vnf) =
 794    let(
 795        verts = vnf[0],
 796        faces = [for (face=vnf[1])
 797                    each (len(face)==3 ? [face] :
 798                    let( tris = polygon_triangulate(verts, face) )
 799                    assert( tris!=undef, "Some `vnf` face cannot be triangulated.")
 800                    tris ) ]
 801    )
 802    [verts, faces];
 803
 804
 805
 806function _vnf_sort_vertices(vnf, idx=[2,1,0]) =
 807    let(
 808        verts = vnf[0],
 809        faces = vnf[1],
 810        vidx = sortidx(verts, idx=idx),
 811        rvidx = sortidx(vidx),
 812        sorted_vnf = [
 813            [ for (i = vidx) verts[i] ],
 814            [ for (face = faces) [ for (i = face) rvidx[i] ] ],
 815        ]
 816    ) sorted_vnf;
 817
 818
 819
 820// Function: vnf_slice()
 821// Synopsis: Slice the faces of a VNF along an axis.
 822// SynTags: VNF
 823// Topics: VNF Manipulation
 824// See Also: vnf_reverse_faces(), vnf_quantize(), vnf_merge_points(), vnf_drop_unused_points(), vnf_triangulate()
 825// Usage:
 826//   sliced = vnf_slice(vnf, dir, cuts);
 827// Description:
 828//   Slice the faces of a VNF along a specified axis direction at a given list of cut points.
 829//   The cut points can appear in any order.  You can use this to refine the faces of a VNF before
 830//   applying a nonlinear transformation to its vertex set.
 831// Arguments:
 832//   vnf = VNF to slice
 833//   dir = normal direction to the slices, either "X", "Y" or "Z"
 834//   cuts = X, Y or Z values where cuts occur
 835// Example(3D):
 836//   include <BOSL2/polyhedra.scad>
 837//   vnf = regular_polyhedron_info("vnf", "dodecahedron", side=12);
 838//   vnf_polyhedron(vnf);
 839//   sliced = vnf_slice(vnf, "X", [-6,-1,10]);
 840//   color("red")vnf_wireframe(sliced,width=.3);
 841function vnf_slice(vnf,dir,cuts) =
 842    let(
 843        cuts = [for (cut=cuts) _shift_cut_plane(vnf,dir,cut)],
 844        vert = vnf[0],
 845        faces = [for(face=vnf[1]) select(vert,face)],
 846        poly_list = _slice_3dpolygons(faces, dir, cuts)
 847    )
 848    vnf_merge_points(vnf_from_polygons(poly_list));
 849
 850
 851function _shift_cut_plane(vnf,dir,cut,off=0.001) =
 852    let(
 853        I = ident(3),
 854        dir_ind = ord(dir)-ord("X"),
 855        verts = vnf[0],
 856        on_cut = [for (x = verts * I[dir_ind]) if(approx(x,cut,eps=1e-4)) 1] != []
 857    ) !on_cut? cut :
 858    _shift_cut_plane(vnf,dir,cut+off);
 859
 860
 861function _split_polygon_at_x(poly, x) =
 862    let(
 863        xs = column(poly,0)
 864    ) (min(xs) >= x || max(xs) <= x)? [poly] :
 865    let(
 866        poly2 = [
 867            for (p = pair(poly,true)) each [
 868                p[0],
 869                if(
 870                    (p[0].x < x && p[1].x > x) ||
 871                    (p[1].x < x && p[0].x > x)
 872                ) let(
 873                    u = (x - p[0].x) / (p[1].x - p[0].x)
 874                ) [
 875                    x,  // Important for later exact match tests
 876                    u*(p[1].y-p[0].y)+p[0].y
 877                ]
 878            ]
 879        ],
 880        out1 = [for (p = poly2) if(p.x <= x) p],
 881        out2 = [for (p = poly2) if(p.x >= x) p],
 882        out3 = [
 883            if (len(out1)>=3) each split_path_at_self_crossings(out1),
 884            if (len(out2)>=3) each split_path_at_self_crossings(out2),
 885        ],
 886        out = [for (p=out3) if (len(p) > 2) list_unwrap(p)]
 887    ) out;
 888
 889
 890function _split_2dpolygons_at_each_x(polys, xs, _i=0) =
 891    _i>=len(xs)? polys :
 892    _split_2dpolygons_at_each_x(
 893        [
 894            for (poly = polys)
 895            each _split_polygon_at_x(poly, xs[_i])
 896        ], xs, _i=_i+1
 897    );
 898
 899/// Internal Function: _slice_3dpolygons()
 900/// Usage:
 901///   splitpolys = _slice_3dpolygons(polys, dir, cuts);
 902/// Topics: Geometry, Polygons, Intersections
 903/// Description:
 904///   Given a list of 3D polygons, a choice of X, Y, or Z, and a cut list, `cuts`, splits all of the polygons where they cross
 905///   X/Y/Z at any value given in cuts.
 906/// Arguments:
 907///   polys = A list of 3D polygons to split.
 908///   dir_ind = slice direction, 0=X, 1=Y, or 2=Z
 909///   cuts = A list of scalar values for locating the cuts
 910function _slice_3dpolygons(polys, dir, cuts) =
 911    assert( [for (poly=polys) if (!is_path(poly,3)) 1] == [], "Expects list of 3D paths.")
 912    assert( is_vector(cuts), "The split list must be a vector.")
 913    assert( in_list(dir, ["X", "Y", "Z"]))
 914    let(
 915        I = ident(3),
 916        dir_ind = ord(dir)-ord("X")
 917    )
 918    flatten([
 919        for (poly = polys)
 920        let( plane = plane_from_polygon(poly))
 921        assert(plane,"Found non-coplanar face.")
 922        let(
 923            normal = point3d(plane),
 924            pnormal = normal - (normal*I[dir_ind])*I[dir_ind]
 925        )
 926        approx(pnormal,[0,0,0]) ? [poly] :
 927        let (
 928            pind = max_index(v_abs(pnormal)),  // project along this direction
 929            otherind = 3-pind-dir_ind,         // keep dir_ind and this direction
 930            keep = [I[dir_ind], I[otherind]],  // dir ind becomes the x dir
 931            poly2d = poly*transpose(keep),     // project to 2d, putting selected direction in the X position
 932            poly_list = [for(p=_split_2dpolygons_at_each_x([poly2d], cuts))
 933                            let(
 934                                a = p*keep,    // unproject, but pind dimension data is missing
 935                                ofs = outer_product((repeat(plane[3], len(a))-a*normal)/plane[pind],I[pind])
 936                             )
 937                             a+ofs]    // ofs computes the missing pind dimension data and adds it back in
 938        )
 939        poly_list
 940    ]);
 941
 942
 943
 944
 945
 946// Section: Turning a VNF into geometry
 947
 948
 949// Module: vnf_polyhedron()
 950// Synopsis: Returns a polyhedron from a VNF or list of VNFs.
 951// SynTags: Geom
 952// Topics: VNF Manipulation
 953// See Also: vnf_wireframe()
 954// Usage:
 955//   vnf_polyhedron(vnf) [ATTACHMENTS];
 956//   vnf_polyhedron([VNF, VNF, VNF, ...]) [ATTACHMENTS];
 957// Description:
 958//   Given a VNF structure, or a list of VNF structures, creates a polyhedron from them.
 959// Arguments:
 960//   vnf = A VNF structure, or list of VNF structures.
 961//   convexity = Max number of times a line could intersect a wall of the shape.
 962//   cp = Centerpoint for determining intersection anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: "centroid"
 963//   anchor = Translate so anchor point is at origin (0,0,0).  See [anchor](attachments.scad#subsection-anchor).  Default: `"origin"`
 964//   spin = Rotate this many degrees around the Z axis after anchor.  See [spin](attachments.scad#subsection-spin).  Default: `0`
 965//   orient = Vector to rotate top towards, after spin.  See [orient](attachments.scad#subsection-orient).  Default: `UP`
 966//   atype = Select "hull" or "intersect" anchor type.  Default: "hull"
 967// Anchor Types:
 968//   "hull" = Anchors to the virtual convex hull of the shape.
 969//   "intersect" = Anchors to the surface of the shape.
 970// Extra Anchors:
 971//   "origin" = Anchor at the origin, oriented UP.
 972module vnf_polyhedron(vnf, convexity=2, cp="centroid", anchor="origin", spin=0, orient=UP, atype="hull") {
 973    vnf = is_vnf_list(vnf)? vnf_join(vnf) : vnf;
 974    assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"");
 975    attachable(anchor,spin,orient, vnf=vnf, extent=atype=="hull", cp=cp) {
 976        polyhedron(vnf[0], vnf[1], convexity=convexity);
 977        children();
 978    }
 979}
 980
 981
 982// Module: vnf_wireframe()
 983// Synopsis: Creates a wireframe model from a VNF.
 984// SynTags: VNF
 985// Topics: VNF Manipulation
 986// See Also: vnf_polyhedron()
 987// Usage:
 988//   vnf_wireframe(vnf, [width]);
 989// Description:
 990//   Given a VNF, creates a wire frame ball-and-stick model of the polyhedron with a cylinder for
 991//   each edge and a sphere at each vertex.  The width parameter specifies the width of the sticks
 992//   that form the wire frame and the diameter of the balls.
 993// Arguments:
 994//   vnf = A vnf structure
 995//   width = width of the cylinders forming the wire frame.  Default: 1
 996// Example:
 997//   $fn=32;
 998//   ball = sphere(r=20, $fn=6);
 999//   vnf_wireframe(ball,width=1);
1000// Example:
1001//   include <BOSL2/polyhedra.scad>
1002//   $fn=32;
1003//   cube_oct = regular_polyhedron_info("vnf",
1004//                      name="cuboctahedron", or=20);
1005//   vnf_wireframe(cube_oct);
1006// Example: The spheres at the vertex are imperfect at aligning with the cylinders, so especially at low $fn things look prety ugly.  This is normal.
1007//   include <BOSL2/polyhedra.scad>
1008//   $fn=8;
1009//   octahedron = regular_polyhedron_info("vnf",
1010//                         name="octahedron", or=20);
1011//   vnf_wireframe(octahedron,width=5);
1012module vnf_wireframe(vnf, width=1)
1013{
1014  no_children($children);
1015  vertex = vnf[0];
1016  edges = unique([for (face=vnf[1], i=idx(face))
1017                    sort([face[i], select(face,i+1)])
1018                 ]);
1019  attachable()
1020  {
1021    union(){
1022      for (e=edges) extrude_from_to(vertex[e[0]],vertex[e[1]]) circle(d=width);
1023      // Identify vertices actually used and draw them
1024      vertused = search(count(len(vertex)), flatten(edges), 1);
1025      for(i=idx(vertex)) if(vertused[i]!=[]) move(vertex[i]) sphere(d=width);
1026    }
1027    union();
1028  }
1029}
1030
1031
1032// Section: Operations on VNFs
1033
1034// Function: vnf_volume()
1035// Synopsis: Returns the volume of a VNF.
1036// Topics: VNF Manipulation
1037// See Also: vnf_area(), vnf_halfspace(), vnf_bend() 
1038// Usage:
1039//   vol = vnf_volume(vnf);
1040// Description:
1041//   Returns the volume enclosed by the given manifold VNF.   The VNF must describe a valid polyhedron with consistent face direction and
1042//   no holes; otherwise the results are undefined.  Returns a positive volume if face direction is clockwise and a negative volume
1043//   if face direction is counter-clockwise.
1044
1045// Divide the polyhedron into tetrahedra with the origin as one vertex and sum up the signed volume.
1046function vnf_volume(vnf) =
1047    let(verts = vnf[0])
1048    sum([
1049         for(face=vnf[1], j=[1:1:len(face)-2])
1050             cross(verts[face[j+1]], verts[face[j]]) * verts[face[0]]
1051    ])/6;
1052
1053
1054// Function: vnf_area()
1055// Synopsis: Returns the surface area of a VNF.
1056// Topics: VNF Manipulation
1057// See Also: vnf_volume(), vnf_halfspace(), vnf_bend() 
1058// Usage:
1059//   area = vnf_area(vnf);
1060// Description:
1061//   Returns the surface area in any VNF by adding up the area of all its faces.  The VNF need not be a manifold.
1062function vnf_area(vnf) =
1063    let(verts=vnf[0])
1064    sum([for(face=vnf[1]) polygon_area(select(verts,face))]);
1065
1066
1067/// Internal Function: _vnf_centroid()
1068/// Usage:
1069///   vol = _vnf_centroid(vnf);
1070/// Description:
1071///   Returns the centroid of the given manifold VNF.  The VNF must describe a valid polyhedron with consistent face direction and
1072///   no holes; otherwise the results are undefined.
1073
1074/// Divide the solid up into tetrahedra with the origin as one vertex.
1075/// The centroid of a tetrahedron is the average of its vertices.
1076/// The centroid of the total is the volume weighted average.
1077function _vnf_centroid(vnf,eps=EPSILON) =
1078    assert(is_vnf(vnf) && len(vnf[0])!=0 && len(vnf[1])!=0,"Invalid or empty VNF given to centroid")
1079    let(
1080        verts = vnf[0],
1081        pos = sum([
1082            for(face=vnf[1], j=[1:1:len(face)-2]) let(
1083                v0  = verts[face[0]],
1084                v1  = verts[face[j]],
1085                v2  = verts[face[j+1]],
1086                vol = cross(v2,v1)*v0
1087            )
1088            [ vol, (v0+v1+v2)*vol ]
1089        ])
1090    )
1091    assert(!approx(pos[0],0, eps), "The vnf has self-intersections.")
1092    pos[1]/pos[0]/4;
1093
1094
1095// Function: vnf_halfspace()
1096// Synopsis: Returns the intersection of the vnf with a half space.
1097// SynTags: VNF
1098// Topics: VNF Manipulation
1099// See Also: vnf_volume(), vnf_area(), vnf_bend() 
1100// Usage:
1101//   newvnf = vnf_halfspace(plane, vnf, [closed], [boundary]);
1102// Description:
1103//   Returns the intersection of the vnf with a half space.  The half space is defined by
1104//   plane = [A,B,C,D], taking the side where the normal [A,B,C] points: Ax+By+Cz≥D.
1105//   If closed is set to false then the cut face is not included in the vnf.  This could
1106//   allow further extension of the vnf by join with other vnfs using {{vnf_join()}}.
1107//   Note that if your given VNF has holes (missing faces) or is not a complete polyhedron
1108//   then closed=true is may produce invalid results when it tries to construct closing faces
1109//   on the cut plane.  Set closed=false for such inputs.
1110//   .
1111//   If you set boundary to true then the return will be the pair [vnf,boundary] where vnf is the
1112//   vnf as usual (with closed=false) and boundary is a list giving each connected component of the cut
1113//   boundary surface.  Each entry in boundary is a list of index values that index into the vnf vertex list (vnf[0]).
1114//   This makes it possible to construct mating shapes, e.g. with {{skin()}} or {{vnf_vertex_array()}} that
1115//   can be combined using {{vnf_join()}} to make a valid polyhedron.
1116//   .
1117//   Note that the input to vnf_halfspace() does not need to be a closed, manifold polyhedron.
1118//   Because it adds the faces on the cut surface, you can use vnf_halfspace() to cap off an open shape if you
1119//   slice through a region that excludes all of the gaps in the input VNF.  
1120// Arguments:
1121//   plane = plane defining the boundary of the half space
1122//   vnf = VNF to cut
1123//   closed = if false do not return the cut face(s) in the returned VNF.  Default: true
1124//   boundary = if true return a pair [vnf,boundary] where boundary is a list of paths on the cut boundary indexed into the VNF vertex list.  If boundary is true, then closed is set to false.  Default: false
1125// Example(3D):
1126//   vnf = cube(10,center=true);
1127//   cutvnf = vnf_halfspace([-1,1,-1,0], vnf);
1128//   vnf_polyhedron(cutvnf);
1129// Example(3D):  Cut face has 2 components
1130//   vnf = path_sweep(circle(r=4, $fn=16),
1131//                    circle(r=20, $fn=64),closed=true);
1132//   cutvnf = vnf_halfspace([-1,1,-4,0], vnf);
1133//   vnf_polyhedron(cutvnf);
1134// Example(3D): Cut face is not simply connected
1135//   vnf = path_sweep(circle(r=4, $fn=16),
1136//                    circle(r=20, $fn=64),closed=true);
1137//   cutvnf = vnf_halfspace([0,0.7,-4,0], vnf);
1138//   vnf_polyhedron(cutvnf);
1139// Example(3D): Cut object has multiple components
1140//   function knot(a,b,t) =   // rolling knot
1141//        [ a * cos (3 * t) / (1 - b* sin (2 *t)),
1142//          a * sin( 3 * t) / (1 - b* sin (2 *t)),
1143//        1.8 * b * cos (2 * t) /(1 - b* sin (2 *t))];
1144//   a = 0.8; b = sqrt (1 - a * a);
1145//   ksteps = 400;
1146//   knot_path = [for (i=[0:ksteps-1]) 50 * knot(a,b,(i/ksteps)*360)];
1147//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1148//   knot=path_sweep(ushape, knot_path, closed=true, method="incremental");
1149//   cut_knot = vnf_halfspace([1,0,0,0], knot);
1150//   vnf_polyhedron(cut_knot);
1151// Example(VPR=[80,0,15]): Cut a sphere with an arbitrary plane
1152//   vnf1=sphere(r=50, style="icosa", $fn=16);
1153//   vnf2=vnf_halfspace([.8,1,-1.5,0], vnf1);
1154//   vnf_polyhedron(vnf2);
1155// Example(VPR=[80,0,15]): Cut it again, but with closed=false to leave an open boundary. 
1156//   vnf1=sphere(r=50, style="icosa", $fn=16);
1157//   vnf2=vnf_halfspace([.8,1,-1.5,0], vnf1);
1158//   vnf3=vnf_halfspace([0,0,-1,0], vnf2, closed=false);
1159//   vnf_polyhedron(vnf3);
1160// Example(VPR=[80,0,15]): Use {vnf_join()} to combine with a mating vnf, in this case a reflection of the part we made. 
1161//   vnf1=sphere(r=50, style="icosa", $fn=16);
1162//   vnf2=vnf_halfspace([.8,1,-1.5,0], vnf1);
1163//   vnf3=vnf_halfspace([0,0,-1,0], vnf2, closed=false);
1164//   vnf4=vnf_join([vnf3, zflip(vnf3,1)]);
1165//   vnf_polyhedron(vnf4);
1166// Example: When the input VNF is a surface with a boundary, if you use the default setting closed=true, then vnf_halfspace() tries to construct closing faces from the edges created by the cut.  These faces may be invalid, for example if the cut points are collinear.  In this example the constructed face is a valid face.
1167//   include <BOSL2/beziers.scad>
1168//   patch=[
1169//          [[10,-10,0],[1,-1,0],[-1,-1,0],[-10,-10,0]],
1170//          [[10,-10,20],[1,-1,20],[-1,-1,20],[-10,-10,20]]
1171//         ];
1172//   vnf=bezier_vnf(patch);
1173//   vnfcut = vnf_halfspace([-.8,0,-1,-14],vnf);
1174//   vnf_polyhedron(vnfcut);
1175// Example: Setting closed to false eliminates this (possibly invalid) face:
1176//   include <BOSL2/beziers.scad>
1177//   patch=[
1178//          [[10,-10,0],[1,-1,0],[-1,-1,0],[-10,-10,0]],
1179//          [[10,-10,20],[1,-1,20],[-1,-1,20],[-10,-10,20]]
1180//         ];
1181//   vnf=bezier_vnf(patch);
1182//   vnfcut = vnf_halfspace([-.8,0,-1,-14],vnf,closed=false);
1183//   vnf_polyhedron(vnfcut);
1184// Example: Here is a VNF that has holes, so it is not a valid manifold. 
1185//   outside = linear_sweep(circle(r=30), h=100, caps=false);
1186//   inside = yrot(7,linear_sweep(circle(r=10), h=120, caps=false));
1187//   open_vnf=vnf_join([outside, vnf_reverse_faces(inside)]);
1188//   vnf_polyhedron(open_vnf);
1189// Example: By cutting it at each end we can create closing faces, resulting in a valid manifold without holes.
1190//   outside = linear_sweep(circle(r=30), h=100, caps=false);
1191//   inside = yrot(11,linear_sweep(circle(r=10), h=120, caps=false));
1192//   open_vnf=vnf_join([outside, vnf_reverse_faces(inside)]);
1193//   vnf = vnf_halfspace([0,0,1,5], vnf_halfspace([0,.7,-1,-75], open_vnf));
1194//   vnf_polyhedron(vnf);
1195// Example: If boundary=true then the return is a list with the VNF and boundary data.  
1196//   vnf = path_sweep(circle(r=4, $fn=16),
1197//                    circle(r=20, $fn=64),closed=true);
1198//   cut_bnd = vnf_halfspace([-1,1,-4,0], vnf, boundary=true);*/
1199//   cutvnf = cut_bnd[0];
1200//   boundary = [for(b=cut_bnd[1]) select(cutvnf[0],b)];
1201//   vnf_polyhedron(cutvnf);
1202//   stroke(boundary,color="red");
1203function vnf_halfspace(plane, vnf, closed=true, boundary=false) =
1204    assert(_valid_plane(plane), "Invalid plane")
1205    assert(is_vnf(vnf), "Invalid vnf")
1206    let(
1207         inside = [for(x=vnf[0]) plane*[each x,-1] >= -EPSILON ? 1 : 0],
1208         vertexmap = [0,each cumsum(inside)],
1209         faces_edges_vertices = _vnfcut(plane, vnf[0],vertexmap,inside, vnf[1], last(vertexmap)),
1210         newvert = concat(bselect(vnf[0],inside), faces_edges_vertices[2])
1211    )
1212    closed==false && !boundary ? [newvert, faces_edges_vertices[0]]
1213  : let(
1214        allpaths = _assemble_paths(newvert, faces_edges_vertices[1]),
1215        newpaths = [for(p=allpaths) if (len(p)>=3) p
1216                                    else assert(approx(p[0],p[1]),"Orphan edge found when assembling cut edges.")
1217           ]
1218    )
1219    boundary ? [[newvert, faces_edges_vertices[0]], newpaths]
1220  : len(newpaths)<=1 ? [newvert, concat(faces_edges_vertices[0], newpaths)]
1221  : let(
1222           M = project_plane(plane),
1223           faceregion = [for(path=newpaths) path2d(apply(M,select(newvert,path)))],
1224           facevnf = vnf_from_region(faceregion,transform=rot_inverse(M),reverse=true)
1225      )
1226      vnf_join([[newvert, faces_edges_vertices[0]], facevnf]);
1227
1228function _assemble_paths(vertices, edges, paths=[],i=0) =
1229     i==len(edges) ? paths :
1230     norm(vertices[edges[i][0]]-vertices[edges[i][1]])<EPSILON ? _assemble_paths(vertices,edges,paths,i+1) :
1231     let(    // Find paths that connects on left side and right side of the edges (if one exists)
1232         left = [for(j=idx(paths)) if (approx(vertices[last(paths[j])],vertices[edges[i][0]])) j],
1233         right = [for(j=idx(paths)) if (approx(vertices[edges[i][1]],vertices[paths[j][0]])) j]
1234     )
1235     assert(len(left)<=1 && len(right)<=1)
1236     let(
1237          keep_path = list_remove(paths,concat(left,right)),
1238          update_path = left==[] && right==[] ? edges[i]
1239                      : left==[] ? concat([edges[i][0]],paths[right[0]])
1240                      : right==[] ?  concat(paths[left[0]],[edges[i][1]])
1241                      : left != right ? concat(paths[left[0]], paths[right[0]])
1242                      : paths[left[0]]
1243     )
1244     _assemble_paths(vertices, edges, concat(keep_path, [update_path]), i+1);
1245
1246
1247function _vnfcut(plane, vertices, vertexmap, inside, faces, vertcount, newfaces=[], newedges=[], newvertices=[], i=0) =
1248   i==len(faces) ? [newfaces, newedges, newvertices] :
1249   let(
1250        pts_inside = select(inside,faces[i])
1251   )
1252   all(pts_inside) ? _vnfcut(plane, vertices, vertexmap, inside, faces, vertcount,
1253                             concat(newfaces, [select(vertexmap,faces[i])]), newedges, newvertices, i+1):
1254   !any(pts_inside) ? _vnfcut(plane, vertices, vertexmap,inside, faces, vertcount, newfaces, newedges, newvertices, i+1):
1255   let(
1256        first = search([[1,0]],pair(pts_inside,wrap=true),0)[0],
1257        second = search([[0,1]],pair(pts_inside,wrap=true),0)[0]
1258   )
1259   assert(len(first)==1 && len(second)==1, "Found concave face in VNF.  Run vnf_triangulate first to ensure convex faces.")
1260   let(
1261        newface = [each select(vertexmap,select(faces[i],second[0]+1,first[0])),vertcount, vertcount+1],
1262        newvert = [plane_line_intersection(plane, select(vertices,select(faces[i],first[0],first[0]+1)),eps=0),
1263                   plane_line_intersection(plane, select(vertices,select(faces[i],second[0],second[0]+1)),eps=0)]
1264   )
1265   true //!approx(newvert[0],newvert[1])
1266       ? _vnfcut(plane, vertices, vertexmap, inside, faces, vertcount+2,
1267                 concat(newfaces, [newface]), concat(newedges,[[vertcount+1,vertcount]]),concat(newvertices,newvert),i+1)
1268   :len(newface)>3
1269       ? _vnfcut(plane, vertices, vertexmap, inside, faces, vertcount+1,
1270                 concat(newfaces, [list_head(newface)]), newedges,concat(newvertices,[newvert[0]]),i+1)
1271   :
1272   _vnfcut(plane, vertices, vertexmap, inside, faces, vertcount,newfaces, newedges, newvert, i+1);
1273
1274
1275
1276
1277function _triangulate_planar_convex_polygons(polys) =
1278    polys==[]? [] :
1279    let(
1280        tris = [for (poly=polys) if (len(poly)==3) poly],
1281        bigs = [for (poly=polys) if (len(poly)>3) poly],
1282        newtris = [for (poly=bigs) select(poly,-2,0)],
1283        newbigs = [for (poly=bigs) select(poly,0,-2)],
1284        newtris2 = _triangulate_planar_convex_polygons(newbigs),
1285        outtris = concat(tris, newtris, newtris2)
1286    ) outtris;
1287
1288//**
1289// this function may produce degenerate triangles:
1290//    _triangulate_planar_convex_polygons([ [for(i=[0:1]) [i,i],
1291//                                           [1,-1], [-1,-1],
1292//                                           for(i=[-1:0]) [i,i] ] ] )
1293//    == [[[-1, -1], [ 0,  0], [0,  0]]
1294//        [[-1, -1], [-1, -1], [0,  0]]
1295//        [[ 1, -1], [-1, -1], [0,  0]]
1296//        [[ 0,  0], [ 1,  1], [1, -1]] ]
1297//
1298
1299// Function: vnf_bend()
1300// Synopsis: Bends a VNF around an axis.
1301// SynTags: VNF
1302// Topics: VNF Manipulation
1303// See Also: vnf_volume(), vnf_area(), vnf_halfspace() 
1304// Usage:
1305//   bentvnf = vnf_bend(vnf,r|d=,[axis=]);
1306// Description:
1307//   Bend a VNF around the X, Y or Z axis, splitting up faces as necessary.  Returns the bent
1308//   VNF.  For bending around the Z axis the input VNF must not cross the Y=0 plane.  For bending
1309//   around the X or Y axes the VNF must not cross the Z=0 plane.  Note that if you wrap a VNF all the way around
1310//   it may intersect itself, which produces an invalid polyhedron.  It is your responsibility to
1311//   avoid this situation.  The 1:1
1312//   radius is where the curved length of the bent VNF matches the length of the original VNF.  If the
1313//   `r` or `d` arguments are given, then they will specify the 1:1 radius or diameter.  If they are
1314//   not given, then the 1:1 radius will be defined by the distance of the furthest vertex in the
1315//   original VNF from the Z=0 plane.  You can adjust the granularity of the bend using the standard
1316//   `$fa`, `$fs`, and `$fn` variables.
1317// Arguments:
1318//   vnf = The original VNF to bend.
1319//   r = If given, the radius where the size of the original shape is the same as in the original.
1320//   ---
1321//   d = If given, the diameter where the size of the original shape is the same as in the original.
1322//   axis = The axis to wrap around.  "X", "Y", or "Z".  Default: "Z"
1323// Example(3D):
1324//   vnf0 = cube([100,40,10], center=true);
1325//   vnf1 = up(50, p=vnf0);
1326//   vnf2 = down(50, p=vnf0);
1327//   bent1 = vnf_bend(vnf1, axis="Y");
1328//   bent2 = vnf_bend(vnf2, axis="Y");
1329//   vnf_polyhedron([bent1,bent2]);
1330// Example(3D):
1331//   vnf0 = linear_sweep(star(n=5,step=2,d=100), height=10);
1332//   vnf1 = up(50, p=vnf0);
1333//   vnf2 = down(50, p=vnf0);
1334//   bent1 = vnf_bend(vnf1, axis="Y");
1335//   bent2 = vnf_bend(vnf2, axis="Y");
1336//   vnf_polyhedron([bent1,bent2]);
1337// Example(3D):
1338//   rgn = union(rect([100,20]),
1339//               rect([20,100]));
1340//   vnf0 = linear_sweep(zrot(45,p=rgn), height=10);
1341//   vnf1 = up(50, p=vnf0);
1342//   vnf2 = down(50, p=vnf0);
1343//   bent1 = vnf_bend(vnf1, axis="Y");
1344//   bent2 = vnf_bend(vnf2, axis="Y");
1345//   vnf_polyhedron([bent1,bent2]);
1346// Example(3D): Bending Around X Axis.
1347//   rgnr = union(
1348//       rect([20,100]),
1349//       back(50, p=trapezoid(w1=40, w2=0, h=20, anchor=FRONT))
1350//   );
1351//   vnf0 = xrot(00,p=linear_sweep(rgnr, height=10));
1352//   vnf1 = up(50, p=vnf0);
1353//   #vnf_polyhedron(vnf1);
1354//   bent1 = vnf_bend(vnf1, axis="X");
1355//   vnf_polyhedron([bent1]);
1356// Example(3D): Bending Around Y Axis.
1357//   rgn = union(
1358//       rect([20,100]),
1359//       back(50, p=trapezoid(w1=40, w2=0, h=20, anchor=FRONT))
1360//   );
1361//   rgnr = zrot(-90, p=rgn);
1362//   vnf0 = xrot(00,p=linear_sweep(rgnr, height=10));
1363//   vnf1 = up(50, p=vnf0);
1364//   #vnf_polyhedron(vnf1);
1365//   bent1 = vnf_bend(vnf1, axis="Y");
1366//   vnf_polyhedron([bent1]);
1367// Example(3D): Bending Around Z Axis.
1368//   rgn = union(
1369//       rect([20,100]),
1370//       back(50, p=trapezoid(w1=40, w2=0, h=20, anchor=FRONT))
1371//   );
1372//   rgnr = zrot(90, p=rgn);
1373//   vnf0 = xrot(90,p=linear_sweep(rgnr, height=10));
1374//   vnf1 = fwd(50, p=vnf0);
1375//   #vnf_polyhedron(vnf1);
1376//   bent1 = vnf_bend(vnf1, axis="Z");
1377//   vnf_polyhedron([bent1]);
1378// Example(3D): Bending more than once around the cylinder
1379//   $fn=32;
1380//   vnf = apply(fwd(5)*yrot(30),cube([100,2,5],center=true));
1381//   bent = vnf_bend(vnf, axis="Z");
1382//   vnf_polyhedron(bent);
1383function vnf_bend(vnf,r,d,axis="Z") =
1384    let(
1385        chk_axis = assert(in_list(axis,["X","Y","Z"])),
1386        verts = vnf[0],
1387        bounds = pointlist_bounds(verts),
1388        bmin = bounds[0],
1389        bmax = bounds[1],
1390        dflt = axis=="Z"?
1391            max(abs(bmax.y), abs(bmin.y)) :
1392            max(abs(bmax.z), abs(bmin.z)),
1393        r = get_radius(r=r,d=d,dflt=dflt),
1394        extent = axis=="X" ? [bmin.y, bmax.y] : [bmin.x, bmax.x]
1395    )
1396    let(
1397        span_chk = axis=="Z"?
1398            assert(bmin.y > 0 || bmax.y < 0, "Entire shape MUST be completely in front of or behind y=0.") :
1399            assert(bmin.z > 0 || bmax.z < 0, "Entire shape MUST be completely above or below z=0."),
1400        steps = 1+ceil(segs(r) * (extent[1]-extent[0])/(2*PI*r)),
1401        step = (extent[1]-extent[0]) / steps,
1402        bend_at = [for(i = [1:1:steps-1]) i*step+extent[0]],
1403        slicedir = axis=="X"? "Y" : "X",   // slice in y dir for X axis case, and x dir otherwise
1404        sliced = vnf_slice(vnf, slicedir, bend_at),
1405        coord = axis=="X" ? [0,sign(bmax.z),0] : axis=="Y" ? [sign(bmax.z),0,0] : [sign(bmax.y),0,0],
1406        new_vert = [for(p=sliced[0])
1407                       let(a=coord*p*180/(PI*r))
1408                       axis=="X"? [p.x, p.z*sin(a), p.z*cos(a)] :
1409                       axis=="Y"? [p.z*sin(a), p.y, p.z*cos(a)] :
1410                       [p.y*sin(a), p.y*cos(a), p.z]]
1411   ) [new_vert,sliced[1]];
1412
1413
1414// Section: Debugging Polyhedrons
1415
1416/// Internal Module: _show_vertices()
1417/// Usage:
1418///   _show_vertices(vertices, [size], [filter=])
1419/// Description:
1420///   Draws all the vertices in an array, at their 3D position, numbered by their
1421///   position in the vertex array.  Also draws any children of this module with
1422///   transparency.
1423/// Arguments:
1424///   vertices = Array of point vertices.
1425///   size = The size of the text used to label the vertices.  Default: 1
1426/// Example:
1427///   verts = [for (z=[-10,10], y=[-10,10], x=[-10,10]) [x,y,z]];
1428///   faces = [[0,1,2], [1,3,2], [0,4,5], [0,5,1], [1,5,7], [1,7,3], [3,7,6], [3,6,2], [2,6,4], [2,4,0], [4,6,7], [4,7,5]];
1429///   _show_vertices(vertices=verts, size=2) {
1430///       polyhedron(points=verts, faces=faces);
1431///   }
1432module _show_vertices(vertices, size=1, filter) {
1433    color("blue") {
1434        dups = vector_search(vertices, EPSILON, vertices);
1435        for (ind = dups) {
1436            if (is_undef(filter) || any(ind, filter)) {
1437                numstr = str_join([for(i=ind) str(i)],",");
1438                v = vertices[ind[0]];
1439                translate(v) {
1440                    rot($vpr) back(size/8){
1441                       linear_extrude(height=size/10, center=true, convexity=10) {
1442                          text(text=numstr, size=size, halign="center");
1443                       }
1444                    }
1445                    sphere(size/10);
1446                }
1447            }
1448        }
1449    }
1450}
1451
1452
1453/// Internal Module: _show_faces()
1454/// Usage:
1455///   _show_faces(vertices, faces, [size=], [filter=]);
1456/// Description:
1457///   Draws all the vertices at their 3D position, numbered in blue by their
1458///   position in the vertex array.  Each face will have their face number drawn
1459///   in red, aligned with the center of face.  All children of this module are drawn
1460///   with transparency.
1461/// Arguments:
1462///   vertices = Array of point vertices.
1463///   faces = Array of faces by vertex numbers.
1464///   size = The size of the text used to label the faces and vertices.  Default: 1
1465/// Example(EdgesMed):
1466///   verts = [for (z=[-10,10], y=[-10,10], x=[-10,10]) [x,y,z]];
1467///   faces = [[0,1,2], [1,3,2], [0,4,5], [0,5,1], [1,5,7], [1,7,3], [3,7,6], [3,6,2], [2,6,4], [2,4,0], [4,6,7], [4,7,5]];
1468///   _show_faces(vertices=verts, faces=faces, size=2) {
1469///       polyhedron(points=verts, faces=faces);
1470///   }
1471module _show_faces(vertices, faces, size=1, filter) {
1472    vlen = len(vertices);
1473    color("red") {
1474        for (i = [0:1:len(faces)-1]) {
1475            face = faces[i];
1476            if (face[0] < 0 || face[1] < 0 || face[2] < 0 || face[0] >= vlen || face[1] >= vlen || face[2] >= vlen) {
1477                echo("BAD FACE: ", vlen=vlen, face=face);
1478            } else if (is_undef(filter) || any(face,filter)) {
1479                verts = select(vertices,face);
1480                c = mean(verts);
1481                v0 = verts[0];
1482                v1 = verts[1];
1483                v2 = verts[2];
1484                dv0 = unit(v1 - v0);
1485                dv1 = unit(v2 - v0);
1486                nrm0 = cross(dv0, dv1);
1487                nrm1 = UP;
1488                axis = vector_axis(nrm0, nrm1);
1489                ang = vector_angle(nrm0, nrm1);
1490                theta = atan2(nrm0[1], nrm0[0]);
1491                translate(c) {
1492                    rotate(a=180-ang, v=axis) {
1493                        zrot(theta-90)
1494                        linear_extrude(height=size/10, center=true, convexity=10) {
1495                            union() {
1496                                text(text=str(i), size=size, halign="center");
1497                                text(text=str("_"), size=size, halign="center");
1498                            }
1499                        }
1500                    }
1501                }
1502            }
1503        }
1504    }
1505}
1506
1507
1508
1509// Module: debug_vnf()
1510// Synopsis: A replacement for `vnf_polyhedron()` to help with debugging.
1511// SynTags: VNF
1512// Topics: VNF Manipulation, Debugging
1513// See Also: vnf_validate()
1514// Usage:
1515//   debug_vnf(vnfs, [faces=], [vertices=], [opacity=], [size=], [convexity=], [filter=]);
1516// Description:
1517//   A drop-in module to replace `vnf_polyhedron()` to help debug vertices and faces.
1518//   Draws all the vertices at their 3D position, numbered in blue by their
1519//   position in the vertex array.  Each face will have its face number drawn
1520//   in red, aligned with the center of face.  All given faces are drawn with
1521//   transparency. All children of this module are drawn with transparency.
1522//   Works best with Thrown-Together preview mode, to see reversed faces.
1523//   You can set opacity to 0 if you want to supress the display of the polyhedron faces.
1524//   .
1525//   The vertex numbers are shown rotated to face you.  As you rotate your polyhedron you
1526//   can rerun the preview to display them oriented for viewing from a different viewpoint.
1527// Topics: Polyhedra, Debugging
1528// Arguments:
1529//   vnf = VNF to display
1530//   ---
1531//   faces = if true display face numbers.  Default: true
1532//   vertices = if true display vertex numbers.  Default: true
1533//   opacity = Opacity of the polyhedron faces.  Default: 0.5
1534//   convexity = The max number of walls a ray can pass through the given polygon paths.
1535//   size = The size of the text used to label the faces and vertices.  Default: 1
1536//   filter = If given a function literal of signature `function(i)`, will only show labels for vertices and faces that have a vertex index that gets a true result from that function.  Default: no filter.
1537// Example(EdgesMed):
1538//   verts = [for (z=[-10,10], a=[0:120:359.9]) [10*cos(a),10*sin(a),z]];
1539//   faces = [[0,1,2], [5,4,3], [0,3,4], [0,4,1], [1,4,5], [1,5,2], [2,5,3], [2,3,0]];
1540//   debug_vnf([verts,faces], size=2);
1541module debug_vnf(vnf, faces=true, vertices=true, opacity=0.5, size=1, convexity=6, filter ) {
1542    no_children($children);
1543    if (faces)
1544      _show_faces(vertices=vnf[0], faces=vnf[1], size=size, filter=filter);
1545    if (vertices)
1546      _show_vertices(vertices=vnf[0], size=size, filter=filter);
1547    if (opacity > 0)
1548      color([0.2, 1.0, 0, opacity])
1549        vnf_polyhedron(vnf,convexity=convexity);
1550}
1551
1552
1553// Module: vnf_validate()
1554// Synopsis: Echos non-manifold VNF errors to the console.
1555// SynTags: VNF
1556// Topics: VNF Manipulation, Debugging
1557// See Also: debug_vnf()
1558// 
1559// Usage: 
1560//   vnf_validate(vnf, [size], [show_warns=], [check_isects=], [opacity=], [adjacent=], [label_verts=], [label_faces=], [wireframe=]);
1561// Description:
1562//   When called as a module, echoes the non-manifold errors to the console, and color hilites the
1563//   bad edges and vertices, overlaid on a transparent gray polyhedron of the VNF.
1564//   .
1565//   Currently checks for these problems:
1566//   .
1567//   Type    | Color    | Code         | Message
1568//   ------- | -------- | ------------ | ---------------------------------
1569//   WARNING | Yellow   | BIG_FACE     | Face has more than 3 vertices, and may confuse CGAL.
1570//   WARNING | Blue     | NULL_FACE    | Face has zero area.
1571//   ERROR   | Cyan     | NONPLANAR    | Face vertices are not coplanar.
1572//   ERROR   | Brown    | DUP_FACE     | Multiple instances of the same face.
1573//   ERROR   | Orange   | MULTCONN     | Multiply Connected Geometry. Too many faces attached at Edge.
1574//   ERROR   | Violet   | REVERSAL     | Faces reverse across edge.
1575//   ERROR   | Red      | T_JUNCTION   | Vertex is mid-edge on another Face.
1576//   ERROR   | Brown    | FACE_ISECT   | Faces intersect.
1577//   ERROR   | Magenta  | HOLE_EDGE    | Edge bounds Hole.
1578//   .
1579//   Still to implement:
1580//   - Overlapping coplanar faces.
1581// Arguments:
1582//   vnf = The VNF to validate.
1583//   size = The width of the lines and diameter of points used to highlight edges and vertices.  Module only.  Default: 1
1584//   ---
1585//   show_warns = If true show warnings for non-triangular faces.  Default: true
1586//   check_isects = If true, performs slow checks for intersecting faces.  Default: false
1587//   opacity = The opacity level to show the polyhedron itself with.    Default: 0.67
1588//   label_verts = If true, shows labels at each vertex that show the vertex number.    Default: false
1589//   label_faces = If true, shows labels at the center of each face that show the face number.    Default: false
1590//   wireframe = If true, shows edges more clearly so you can see them in Thrown Together mode.    Default: false
1591//   adjacent = If true, only display faces adjacent to a vertex listed in the errors.    Default: false
1592// Example(3D,Edges): BIG_FACE Warnings; Faces with More Than 3 Vertices.  CGAL often will fail to accept that a face is planar after a rotation, if it has more than 3 vertices.
1593//   vnf = skin([
1594//       path3d(regular_ngon(n=3, d=100),0),
1595//       path3d(regular_ngon(n=5, d=100),100)
1596//   ], slices=0, caps=true, method="tangent");
1597//   vnf_validate(vnf);
1598// Example(3D,Edges): NONPLANAR Errors; Face Vertices are Not Coplanar
1599//   a = [  0,  0,-50];
1600//   b = [-50,-50, 50];
1601//   c = [-50, 50, 50];
1602//   d = [ 50, 50, 60];
1603//   e = [ 50,-50, 50];
1604//   vnf = vnf_from_polygons([
1605//       [a, b, e], [a, c, b], [a, d, c], [a, e, d], [b, c, d, e]
1606//   ]);
1607//   vnf_validate(vnf);
1608// Example(3D,Edges): MULTCONN Errors; More Than Two Faces Attached to the Same Edge.  This confuses CGAL, and can lead to failed renders.
1609//   vnf = vnf_triangulate(linear_sweep(union(square(50), square(50,anchor=BACK+RIGHT)), height=50));
1610//   vnf_validate(vnf);
1611// Example(3D,Edges): REVERSAL Errors; Faces Reversed Across Edge
1612//   vnf1 = skin([
1613//       path3d(square(100,center=true),0),
1614//       path3d(square(100,center=true),100),
1615//   ], slices=0, caps=false);
1616//   vnf = vnf_join([vnf1, vnf_from_polygons([
1617//       [[-50,-50,  0], [ 50, 50,  0], [-50, 50,  0]],
1618//       [[-50,-50,  0], [ 50,-50,  0], [ 50, 50,  0]],
1619//       [[-50,-50,100], [-50, 50,100], [ 50, 50,100]],
1620//       [[-50,-50,100], [ 50,-50,100], [ 50, 50,100]],
1621//   ])]);
1622//   vnf_validate(vnf);
1623// Example(3D,Edges): T_JUNCTION Errors; Vertex is Mid-Edge on Another Face.
1624//   vnf = [
1625//       [
1626//           each path3d(square(100,center=true),0),
1627//           each path3d(square(100,center=true),100),
1628//           [0,-50,100],
1629//       ], [
1630//          [0,2,1], [0,3,2], [0,8,4], [0,1,8], [1,5,8],
1631//          [0,4,3], [4,7,3], [1,2,5], [2,6,5], [3,7,6],
1632//          [3,6,2], [4,5,6], [4,6,7],
1633//       ]
1634//   ];
1635//   vnf_validate(vnf);
1636// Example(3D,Edges): FACE_ISECT Errors; Faces Intersect
1637//   vnf = vnf_join([
1638//       vnf_triangulate(linear_sweep(square(100,center=true), height=100)),
1639//       move([75,35,30],p=vnf_triangulate(linear_sweep(square(100,center=true), height=100)))
1640//   ]);
1641//   vnf_validate(vnf,size=2,check_isects=true);
1642// Example(3D,Edges): HOLE_EDGE Errors; Edges Adjacent to Holes.
1643//   vnf = skin([
1644//       path3d(regular_ngon(n=4, d=100),0),
1645//       path3d(regular_ngon(n=5, d=100),100)
1646//   ], slices=0, caps=false);
1647//   vnf_validate(vnf,size=2);
1648
1649
1650//   Returns a list of non-manifold errors with the given VNF.
1651//   Each error has the format `[ERR_OR_WARN,CODE,MESG,POINTS,COLOR]`.
1652function _vnf_validate(vnf, show_warns=true, check_isects=false) =
1653    assert(is_vnf(vnf), "Invalid VNF")
1654    let(
1655        varr = vnf[0],
1656        faces = vnf[1],
1657        lvarr = len(varr),
1658        edges = sort([
1659            for (face=faces, edge=pair(face,true))
1660            edge[0]<edge[1]? edge : [edge[1],edge[0]]
1661        ]),
1662        dfaces = [
1663            for (face=faces) let(
1664                face=deduplicate_indexed(varr,face,closed=true)
1665            ) if(len(face)>=3)
1666            face
1667        ],
1668        face_areas = [
1669            for (face = faces)
1670            len(face) < 3? 0 :
1671            polygon_area([for (k=face) varr[k]])
1672        ],
1673        edgecnts = unique_count(edges),
1674        uniq_edges = edgecnts[0],
1675        issues = []
1676    )
1677    let(
1678        big_faces = !show_warns? [] : [
1679            for (face = faces)
1680            if (len(face) > 3)
1681            _vnf_validate_err("BIG_FACE", face)
1682        ],
1683        null_faces = !show_warns? [] : [
1684            for (i = idx(faces)) let(
1685                face = faces[i],
1686                area = face_areas[i]
1687            )
1688            if (is_num(area) && abs(area) < EPSILON)
1689            _vnf_validate_err("NULL_FACE", face)
1690        ],
1691        issues = concat(big_faces, null_faces)
1692    )
1693    let(
1694        bad_indices = [
1695            for (face = faces, idx = face)
1696            if (idx < 0 || idx >= lvarr)
1697            _vnf_validate_err("BAD_INDEX", [idx])
1698        ],
1699        issues = concat(issues, bad_indices)
1700    ) bad_indices? issues :
1701    let(
1702        repeated_faces = [
1703            for (i=idx(dfaces), j=idx(dfaces))
1704            if (i!=j) let(
1705                face1 = dfaces[i],
1706                face2 = dfaces[j]
1707            ) if (min(face1) == min(face2)) let(
1708                min1 = min_index(face1),
1709                min2 = min_index(face2)
1710            ) if (min1 == min2) let(
1711                sface1 = list_rotate(face1,min1),
1712                sface2 = list_rotate(face2,min2)
1713            ) if (sface1 == sface2)
1714            _vnf_validate_err("DUP_FACE", sface1)
1715        ],
1716        issues = concat(issues, repeated_faces)
1717    ) repeated_faces? issues :
1718    let(
1719        multconn_edges = unique([
1720            for (i = idx(uniq_edges))
1721            if (edgecnts[1][i]>2)
1722            _vnf_validate_err("MULTCONN", uniq_edges[i])
1723        ]),
1724        issues = concat(issues, multconn_edges)
1725    ) multconn_edges? issues :
1726    let(
1727        reversals = unique([
1728            for(i = idx(dfaces), j = idx(dfaces)) if(i != j)
1729            for(edge1 = pair(faces[i],true))
1730            for(edge2 = pair(faces[j],true))
1731            if(edge1 == edge2)  // Valid adjacent faces will never have the same vertex ordering.
1732            if(_edge_not_reported(edge1, varr, multconn_edges))
1733            _vnf_validate_err("REVERSAL", edge1)
1734        ]),
1735        issues = concat(issues, reversals)
1736    ) reversals? issues :
1737    let(
1738        t_juncts = unique([
1739            for (v=idx(varr), edge=uniq_edges) let(
1740                ia = edge[0],
1741                ib = v,
1742                ic = edge[1]
1743            )
1744            if (ia!=ib && ib!=ic && ia!=ic) let(
1745                a = varr[ia],
1746                b = varr[ib],
1747                c = varr[ic]
1748            )
1749            if (!approx(a,b) && !approx(b,c) && !approx(a,c)) let(
1750                pt = line_closest_point([a,c],b,SEGMENT)
1751            )
1752            if (approx(pt,b))
1753            _vnf_validate_err("T_JUNCTION", [ib])
1754        ]),
1755        issues = concat(issues, t_juncts)
1756    ) t_juncts? issues :
1757    let(
1758        isect_faces = !check_isects? [] : unique([
1759            for (i = [0:1:len(faces)-2]) let(
1760                f1 = faces[i],
1761                poly1   = select(varr, faces[i]),
1762                plane1  = plane3pt(poly1[0], poly1[1], poly1[2]),
1763                normal1 = [plane1[0], plane1[1], plane1[2]]
1764            )
1765            for (j = [i+1:1:len(faces)-1]) let(
1766                f2 = faces[j],
1767                poly2 = select(varr, f2),
1768                val = poly2 * normal1
1769            )
1770            if( min(val)<=plane1[3] && max(val)>=plane1[3] ) let(
1771                plane2  = plane_from_polygon(poly2),
1772                normal2 = [plane2[0], plane2[1], plane2[2]],
1773                val = poly1 * normal2
1774            )
1775            if( min(val)<=plane2[3] && max(val)>=plane2[3] ) let(
1776                shared_edges = [
1777                    for (edge1 = pair(f1, true), edge2 = pair(f2, true))
1778                    if (edge1 == [edge2[1], edge2[0]]) 1
1779                ]
1780            )
1781            if (!shared_edges) let(
1782                line = plane_intersection(plane1, plane2)
1783            )
1784            if (!is_undef(line)) let(
1785                isects = polygon_line_intersection(poly1, line)
1786            )
1787            if (!is_undef(isects))
1788            for (isect = isects)
1789            if (len(isect) > 1) let(
1790                isects2 = polygon_line_intersection(poly2, isect, bounded=true)
1791            )
1792            if (!is_undef(isects2))
1793            for (seg = isects2)
1794            if (seg[0] != seg[1])
1795            _vnf_validate_err("FACE_ISECT", seg)
1796        ]),
1797        issues = concat(issues, isect_faces)
1798    ) isect_faces? issues :
1799    let(
1800        hole_edges = unique([
1801            for (i=idx(uniq_edges))
1802            if (edgecnts[1][i]<2)
1803            if (_pts_not_reported(uniq_edges[i], varr, t_juncts))
1804            if (_pts_not_reported(uniq_edges[i], varr, isect_faces))
1805            _vnf_validate_err("HOLE_EDGE", uniq_edges[i])
1806        ]),
1807        issues = concat(issues, hole_edges)
1808    ) hole_edges? issues :
1809    let(
1810        nonplanars = unique([
1811            for (i = idx(faces)) let(
1812                face = faces[i],
1813                area = face_areas[i],
1814                faceverts = [for (k=face) varr[k]]
1815            )
1816            if (is_num(area) && abs(area) > EPSILON)
1817            if (!is_coplanar(faceverts))
1818            _vnf_validate_err("NONPLANAR", face)
1819        ]),
1820        issues = concat(issues, nonplanars)
1821    ) issues;
1822
1823
1824_vnf_validate_errs = [
1825    ["BIG_FACE",    "WARNING", "cyan",    "Face has more than 3 vertices, and may confuse CGAL"],
1826    ["NULL_FACE",   "WARNING", "blue",    "Face has zero area."],
1827    ["BAD_INDEX",   "ERROR",   "cyan",    "Invalid face vertex index."],
1828    ["NONPLANAR",   "ERROR",   "yellow",  "Face vertices are not coplanar"],
1829    ["DUP_FACE",    "ERROR",   "brown",   "Multiple instances of the same face."],
1830    ["MULTCONN",    "ERROR",   "orange",  "Multiply Connected Geometry. Too many faces attached at Edge"],
1831    ["REVERSAL",    "ERROR",   "violet",  "Faces Reverse Across Edge"],
1832    ["T_JUNCTION",  "ERROR",   "magenta", "Vertex is mid-edge on another Face"],
1833    ["FACE_ISECT",  "ERROR",   "brown",   "Faces intersect"],
1834    ["HOLE_EDGE",   "ERROR",   "red",     "Edge bounds Hole"]
1835];
1836
1837
1838function _vnf_validate_err(name, extra) =
1839    let(
1840        info = [for (x = _vnf_validate_errs) if (x[0] == name) x][0]
1841    ) concat(info, [extra]);
1842
1843
1844function _pts_not_reported(pts, varr, reports) =
1845    [
1846        for (i = pts, report = reports, pt = report[3])
1847        if (varr[i] == pt) 1
1848    ] == [];
1849
1850
1851function _edge_not_reported(edge, varr, reports) =
1852    let(
1853        edge = sort([for (i=edge) varr[i]])
1854    ) [
1855        for (report = reports) let(
1856            pts = sort(report[3])
1857        ) if (len(pts)==2 && edge == pts) 1
1858    ] == [];
1859
1860
1861module vnf_validate(vnf, size=1, show_warns=true, check_isects=false, opacity=0.67, adjacent=false, label_verts=false, label_faces=false, wireframe=false) {
1862    no_children($children);
1863    vcount = len(vnf[0]);
1864    fcount = len(vnf[1]);
1865    vnf = vnf_merge_points(vnf);
1866    faults = _vnf_validate(
1867        vnf, show_warns=show_warns,
1868        check_isects=check_isects
1869    );
1870    verts = vnf[0];
1871    vnf_changed = len(verts)!=vcount || len(vnf[1])!=fcount;
1872    if (!faults) {
1873        echo("VNF appears valid.");
1874    }
1875    if (vnf_changed) echo("VNF changed when merging points; unable to display indices");
1876    for (fault = faults) {
1877        err = fault[0];
1878        typ = fault[1];
1879        clr = fault[2];
1880        msg = fault[3];
1881        idxs = fault[4];
1882        pts = err=="FACE_ISECT" ? idxs : [for (i=idxs) if(is_finite(i) && i>=0 && i<len(verts)) verts[i]];
1883        if (vnf_changed || err=="FACE_ISECT")
1884          echo(str(typ, " ", err, " (", clr ,"): ", msg, " at ", pts));
1885        else
1886          echo(str(typ, " ", err, " (", clr ,"): ", msg, " at ", pts, " indices: ", idxs));
1887        color(clr) {
1888            if (is_vector(pts[0])) {
1889                if (len(pts)==2) {
1890                    stroke(pts, width=size, endcaps="butt", $fn=8);
1891                } else if (len(pts)>2) {
1892                    stroke(pts, width=size, closed=true, $fn=8);
1893                    polyhedron(pts,[[for (i=idx(pts)) i]]);
1894                } else {
1895                    move_copies(pts) sphere(d=size*3, $fn=18);
1896                }
1897            }
1898        }
1899    }
1900    badverts = unique([for (fault=faults) each fault[4]]);
1901    badverts2 = unique([for (j=idx(verts), i=badverts) if (i!=j && verts[i]==verts[j]) j]);
1902    all_badverts = unique(concat(badverts, badverts2));
1903    adjacent = !faults? false : adjacent;
1904    filter_fn = !adjacent? undef : function(i) in_list(i,all_badverts);
1905    adj_vnf = !adjacent? vnf : [
1906        verts, [for (face=vnf[1]) if (any(face,filter_fn)) face]
1907    ];
1908    if (wireframe) {
1909        vnf_wireframe(adj_vnf, width=size*0.25);
1910    }
1911    if (label_verts) {
1912        debug_vnf(adj_vnf, size=size*3, opacity=0, faces=false, vertices=true, filter=filter_fn);
1913    }
1914    if (label_faces) {
1915        debug_vnf(vnf, size=size*3, opacity=0, faces=true, vertices=false, filter=filter_fn);
1916    }
1917    if (opacity > 0) {
1918        color([0.5,1,0.5,opacity]) vnf_polyhedron(adj_vnf);
1919    }
1920}
1921
1922
1923
1924// vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap