1//////////////////////////////////////////////////////////////////////
   2// LibFile: skin.scad
   3//   This file provides functions and modules that construct shapes from a list of cross sections.
   4//   In the case of skin() you specify each cross sectional shape yourself, and the number of
   5//   points can vary.  The various forms of sweep use a fixed shape, which may follow a path, or
   6//   be transformed in other ways to produce the list of cross sections.  In all cases it is the
   7//   user's responsibility to avoid creating a self-intersecting shape, which will produce
   8//   cryptic CGAL errors.  This file was inspired by list-comprehension-demos skin():
   9//   - https://github.com/openscad/list-comprehension-demos/blob/master/skin.scad
  10// Includes:
  11//   include <BOSL2/std.scad>
  12// FileGroup: Advanced Modeling
  13// FileSummary: Construct 3D shapes from 2D cross sections of the desired shape.
  14// FileFootnotes: STD=Included in std.scad
  15//////////////////////////////////////////////////////////////////////
  16
  17
  18// Section: Skin and sweep
  19
  20// Function&Module: skin()
  21// Synopsis: Connect a sequence of arbitrary polygons into a 3D object. 
  22// SynTags: VNF, Geom
  23// Topics: Extrusion, Skin
  24// See Also: sweep(), linear_sweep(), rotate_sweep(), spiral_sweep(), path_sweep(), offset_sweep()
  25// Usage: As module:
  26//   skin(profiles, slices, [z=], [refine=], [method=], [sampling=], [caps=], [closed=], [style=], [convexity=], [anchor=],[cp=],[spin=],[orient=],[atype=]) [ATTACHMENTS];
  27// Usage: As function:
  28//   vnf = skin(profiles, slices, [z=], [refine=], [method=], [sampling=], [caps=], [closed=], [style=], [anchor=],[cp=],[spin=],[orient=],[atype=]);
  29// Description:
  30//   Given a list of two or more path `profiles` in 3d space, produces faces to skin a surface between
  31//   the profiles.  Optionally the first and last profiles can have endcaps, or the first and last profiles
  32//   can be connected together.  Each profile should be roughly planar, but some variation is allowed.
  33//   Each profile must rotate in the same clockwise direction.  If called as a function, returns a
  34//   [VNF structure](vnf.scad) `[VERTICES, FACES]`.  If called as a module, creates a polyhedron
  35//    of the skinned profiles.
  36//   .
  37//   The profiles can be specified either as a list of 3d curves or they can be specified as
  38//   2d curves with heights given in the `z` parameter.  It is your responsibility to ensure
  39//   that the resulting polyhedron is free from self-intersections, which would make it invalid
  40//   and can result in cryptic CGAL errors upon rendering with a second object present, even though the polyhedron appears
  41//   OK during preview or when rendered by itself.
  42//   .
  43//   For this operation to be well-defined, the profiles must all have the same vertex count and
  44//   we must assume that profiles are aligned so that vertex `i` links to vertex `i` on all polygons.
  45//   Many interesting cases do not comply with this restriction.  Two basic methods can handle
  46//   these cases: either subdivide edges (insert additional points along edges)
  47//   or duplicate vertcies (insert edges of length 0) so that both polygons have
  48//   the same number of points.
  49//   Duplicating vertices allows two distinct points in one polygon to connect to a single point
  50//   in the other one, creating
  51//   triangular faces.  You can adjust non-matching polygons yourself
  52//   either by resampling them using {{subdivide_path()}} or by duplicating vertices using
  53//   `repeat_entries`.  It is OK to pass a polygon that has the same vertex repeated, such as
  54//   a square with 5 points (two of which are identical), so that it can match up to a pentagon.
  55//   Such a combination would create a triangular face at the location of the duplicated vertex.
  56//   Alternatively, `skin` provides methods (described below) for inserting additional vertices
  57//   automatically to make incompatible paths match.
  58//   .
  59//   In order for skinned surfaces to look good it is usually necessary to use a fine sampling of
  60//   points on all of the profiles, and a large number of extra interpolated slices between the
  61//   profiles that you specify.  It is generally best if the triangles forming your polyhedron
  62//   are approximately equilateral.  The `slices` parameter specifies the number of slices to insert
  63//   between each pair of profiles, either a scalar to insert the same number everywhere, or a vector
  64//   to insert a different number between each pair.
  65//   .
  66//   Resampling may occur, depending on the `method` parameter, to make profiles compatible.
  67//   To force (possibly additional) resampling of the profiles to increase the point density you can set `refine=N`, which
  68//   will multiply the number of points on your profile by `N`.  You can choose between two resampling
  69//   schemes using the `sampling` option, which you can set to `"length"` or `"segment"`.
  70//   The length resampling method resamples proportional to length.
  71//   The segment method divides each segment of a profile into the same number of points.
  72//   This means that if you refine a profile with the "segment" method you will get N points
  73//   on each edge, but if you refine a profile with the "length" method you will get new points
  74//   distributed around the profile based on length, so small segments will get fewer new points than longer ones.
  75//   A uniform division may be impossible, in which case the code computes an approximation, which may result
  76//   in arbitrary distribution of extra points.  See {{subdivide_path()}} for more details.
  77//   Note that when dealing with continuous curves it is always better to adjust the
  78//   sampling in your code to generate the desired sampling rather than using the `refine` argument.
  79//   .
  80//   You can choose from five methods for specifying alignment for incommensurate profiles.
  81//   The available methods are `"distance"`, `"fast_distance"`, `"tangent"`, `"direct"` and `"reindex"`.
  82//   It is useful to distinguish between continuous curves like a circle and discrete profiles
  83//   like a hexagon or star, because the algorithms' suitability depend on this distinction.
  84//   .
  85//   The default method for aligning profiles is `method="direct"`.
  86//   If you simply supply a list of compatible profiles it will link them up
  87//   exactly as you have provided them.  You may find that profiles you want to connect define the
  88//   right shapes but the point lists don't start from points that you want aligned in your skinned
  89//   polyhedron.  You can correct this yourself using `reindex_polygon`, or you can use the "reindex"
  90//   method which will look for the index choice that will minimize the length of all of the edges
  91//   in the polyhedron&mdash;it will produce the least twisted possible result.  This algorithm has quadratic
  92//   run time so it can be slow with very large profiles.
  93//   .
  94//   When the profiles are incommensurate, the "direct" and "reindex" resample them to match.  As noted above,
  95//   for continuous input curves, it is better to generate your curves directly at the desired sample size,
  96//   but for mapping between a discrete profile like a hexagon and a circle, the hexagon must be resampled
  97//   to match the circle.  When you use "direct" or "reindex" the default `sampling` value is
  98//   of `sampling="length"` to approximate a uniform length sampling of the profile.  This will generally
  99//   produce the natural result for connecting two continuously sampled profiles or a continuous
 100//   profile and a polygonal one.  However depending on your particular case,
 101//   `sampling="segment"` may produce a more pleasing result.  These two approaches differ only when
 102//   the segments of your input profiles have unequal length.
 103//   .
 104//   The "distance", "fast_distance" and "tangent" methods work by duplicating vertices to create
 105//   triangular faces.  In the skined object created by two polygons, every vertex of a polygon must
 106//   have an edge that connects to some vertex on the other one.  If you connect two squares this can be
 107//   accomplished with four edges, but if you want to connect a square to a pentagon you must add a
 108//   fifth edge for the "extra" vertex on the pentagon.  You must now decide which vertex on the square to
 109//   connect the "extra" edge to.  How do you decide where to put that fifth edge?  The "distance" method answers this
 110//   question by using an optimization: it minimizes the total length of all the edges connecting
 111//   the two polygons.   This algorithm generally produces a good result when both profiles are discrete ones with
 112//   a small number of vertices.  It is computationally intensive (O(N^3)) and may be
 113//   slow on large inputs.  The resulting surfaces generally have curved faces, so be
 114//   sure to select a sufficiently large value for `slices` and `refine`.  Note that for
 115//   this method, `sampling` must be set to `"segment"`, and hence this is the default setting.
 116//   Using sampling by length would ignore the repeated vertices and ruin the alignment.
 117//   The "fast_distance" method restricts the optimization by assuming that an edge should connect
 118//   vertex 0 of the two polygons.  This reduces the run time to O(N^2) and makes
 119//   the method usable on profiles with more points if you take care to index the inputs to match.
 120//   .
 121//   The `"tangent"` method generally produces good results when
 122//   connecting a discrete polygon to a convex, finely sampled curve.  Given a polygon and a curve, consider one edge
 123//   on the polygon.  Find a plane passing through the edge that is tangent to the curve.  The endpoints of the edge and
 124//   the point of tangency define a triangular face in the output polyhedron.  If you work your way around the polygon
 125//   edges, you can establish a series of triangular faces in this way, with edges linking the polygon to the curve.
 126//   You can then complete the edge assignment by connecting all the edges in between the triangular faces together,
 127//   with many edges meeting at each polygon vertex.  The result is an alternation of flat triangular faces with conical
 128//   curves joining them.  Another way to think about it is that it splits the points on the curve up into groups and
 129//   connects all the points in one group to the same vertex on the polygon.
 130//   .
 131//   The "tangent" method may fail if the curved profile is non-convex, or doesn't have enough points to distinguish
 132//   all of the tangent points from each other.    The algorithm treats whichever input profile has fewer points as the polygon
 133//   and the other one as the curve.  Using `refine` with this method will have little effect on the model, so
 134//   you should do it only for agreement with other profiles, and these models are linear, so extra slices also
 135//   have no effect.  For best efficiency set `refine=1` and `slices=0`.  As with the "distance" method, refinement
 136//   must be done using the "segment" sampling scheme to preserve alignment across duplicated points.
 137//   Note that the "tangent" method produces similar results to the "distance" method on curved inputs.  If this
 138//   method fails due to concavity, "fast_distance" may be a good option.
 139//   .
 140//   It is possible to specify `method` and `refine` as arrays, but it is important to observe
 141//   matching rules when you do this.  If a pair of profiles is connected using "tangent" or "distance"
 142//   then the `refine` values for those two profiles must be equal.  If a profile is connected by
 143//   a vertex duplicating method on one side and a resampling method on the other side, then
 144//   `refine` must be set so that the resulting number of vertices matches the number that is
 145//   used for the resampled profiles.  The best way to avoid confusion is to ensure that the
 146//   profiles connected by "direct" or "reindex" all have the same number of points and at the
 147//   transition, the refined number of points matches.
 148//   .
 149// Arguments:
 150//   profiles = list of 2d or 3d profiles to be skinned.  (If 2d must also give `z`.)
 151//   slices = scalar or vector number of slices to insert between each pair of profiles.  Set to zero to use only the profiles you provided.  Recommend starting with a value around 10.
 152//   ---
 153//   refine = resample profiles to this number of points per edge.  Can be a list to give a refinement for each profile.  Recommend using a value above 10 when using the "distance" or "fast_distance" methods.  Default: 1.
 154//   sampling = sampling method to use with "direct" and "reindex" methods.  Can be "length" or "segment".  Ignored if any profile pair uses either the "distance", "fast_distance", or "tangent" methods.  Default: "length".
 155//   closed = set to true to connect first and last profile (to make a torus).  Default: false
 156//   caps = true to create endcap faces when closed is false.  Can be a length 2 boolean array.  Default is true if closed is false.
 157//   method = method for connecting profiles, one of "distance", "fast_distance", "tangent", "direct" or "reindex".  Default: "direct".
 158//   z = array of height values for each profile if the profiles are 2d
 159//   convexity = convexity setting for use with polyhedron.  (module only) Default: 10
 160//   anchor = Translate so anchor point is at the origin.  Default: "origin"
 161//   spin = Rotate this many degrees around Z axis after anchor.  Default: 0
 162//   orient = Vector to rotate top towards after spin
 163//   atype = Select "hull" or "intersect" anchor types. Default: "hull"
 164//   cp = Centerpoint for determining "intersect" anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: "centroid"
 165//   style = vnf_vertex_array style.  Default: "min_edge"
 166// Anchor Types:
 167//   "hull" = Anchors to the virtual convex hull of the shape.
 168//   "intersect" = Anchors to the surface of the shape.
 169// Example:
 170//   skin([octagon(4), circle($fn=70,r=2)], z=[0,3], slices=10);
 171// Example: Rotating the pentagon place the zero index at different locations, giving a twist
 172//   skin([rot(90,p=pentagon(4)), circle($fn=80,r=2)], z=[0,3], slices=10);
 173// Example: You can untwist it with the "reindex" method
 174//   skin([rot(90,p=pentagon(4)), circle($fn=80,r=2)], z=[0,3], slices=10, method="reindex");
 175// Example: Offsetting the starting edge connects to circles in an interesting way:
 176//   circ = circle($fn=80, r=3);
 177//   skin([circ, rot(110,p=circ)], z=[0,5], slices=20);
 178// Example(FlatSpin,VPD=20):
 179//   skin([ yrot(37,p=path3d(circle($fn=128, r=4))), path3d(square(3),3)], method="reindex",slices=10);
 180// Example(FlatSpin,VPD=16): Ellipses connected with twist
 181//   ellipse = xscale(2.5,p=circle($fn=80));
 182//   skin([ellipse, rot(45,p=ellipse)], z=[0,1.5], slices=10);
 183// Example(FlatSpin,VPD=16): Ellipses connected without a twist.  (Note ellipses stay in the same position: just the connecting edges are different.)
 184//   ellipse = xscale(2.5,p=circle($fn=80));
 185//   skin([ellipse, rot(45,p=ellipse)], z=[0,1.5], slices=10, method="reindex");
 186// Example(FlatSpin,VPD=500):
 187//   $fn=24;
 188//   skin([
 189//         yrot(0, p=yscale(2,p=path3d(circle(d=75)))),
 190//         [[40,0,100], [35,-15,100], [20,-30,100],[0,-40,100],[-40,0,100],[0,40,100],[20,30,100], [35,15,100]]
 191//   ],slices=10);
 192// Example(FlatSpin,VPD=600):
 193//   $fn=48;
 194//   skin([
 195//       for (b=[0,90]) [
 196//           for (a=[360:-360/$fn:0.01])
 197//               point3d(polar_to_xy((100+50*cos((a+b)*2))/2,a),b/90*100)
 198//       ]
 199//   ], slices=20);
 200// Example: Vaccum connector example from list-comprehension-demos
 201//   include <BOSL2/rounding.scad>
 202//   $fn=32;
 203//   base = round_corners(square([2,4],center=true), radius=0.5);
 204//   skin([
 205//       path3d(base,0),
 206//       path3d(base,2),
 207//       path3d(circle(r=0.5),3),
 208//       path3d(circle(r=0.5),4),
 209//       for(i=[0:2]) each [path3d(circle(r=0.6), i+4),
 210//                          path3d(circle(r=0.5), i+5)]
 211//   ],slices=0);
 212// Example: Vaccum nozzle example from list-comprehension-demos, using "length" sampling (the default)
 213//   xrot(90)down(1.5)
 214//   difference() {
 215//       skin(
 216//           [square([2,.2],center=true),
 217//            circle($fn=64,r=0.5)], z=[0,3],
 218//           slices=40,sampling="length",method="reindex");
 219//       skin(
 220//           [square([1.9,.1],center=true),
 221//            circle($fn=64,r=0.45)], z=[-.01,3.01],
 222//           slices=40,sampling="length",method="reindex");
 223//   }
 224// Example: Same thing with "segment" sampling
 225//   xrot(90)down(1.5)
 226//   difference() {
 227//       skin(
 228//           [square([2,.2],center=true),
 229//            circle($fn=64,r=0.5)], z=[0,3],
 230//           slices=40,sampling="segment",method="reindex");
 231//       skin(
 232//           [square([1.9,.1],center=true),
 233//            circle($fn=64,r=0.45)], z=[-.01,3.01],
 234//           slices=40,sampling="segment",method="reindex");
 235//   }
 236// Example: Forma Candle Holder (from list-comprehension-demos)
 237//   r = 50;
 238//   height = 140;
 239//   layers = 10;
 240//   wallthickness = 5;
 241//   holeradius = r - wallthickness;
 242//   difference() {
 243//       skin([for (i=[0:layers-1]) zrot(-30*i,p=path3d(hexagon(ir=r),i*height/layers))],slices=0);
 244//       up(height/layers) cylinder(r=holeradius, h=height);
 245//   }
 246// Example(FlatSpin,VPD=300): A box that is octagonal on the outside and circular on the inside
 247//   height = 45;
 248//   sub_base = octagon(d=71, rounding=2, $fn=128);
 249//   base = octagon(d=75, rounding=2, $fn=128);
 250//   interior = regular_ngon(n=len(base), d=60);
 251//   right_half()
 252//     skin([ sub_base, base, base, sub_base, interior], z=[0,2,height, height, 2], slices=0, refine=1, method="reindex");
 253// Example: Connecting a pentagon and circle with the "tangent" method produces large triangular faces and cone shaped corners.
 254//   skin([pentagon(4), circle($fn=80,r=2)], z=[0,3], slices=10, method="tangent");
 255// Example: rounding corners of a square.  Note that `$fn` makes the number of points constant, and avoiding the `rounding=0` case keeps everything simple.  In this case, the connections between profiles are linear, so there is no benefit to setting `slices` bigger than zero.
 256//   shapes = [for(i=[.01:.045:2])zrot(-i*180/2,cp=[-8,0,0],p=xrot(90,p=path3d(regular_ngon(n=4, side=4, rounding=i, $fn=64))))];
 257//   rotate(180) skin( shapes, slices=0);
 258// Example: Here's a simplified version of the above, with `i=0` included.  That first layer doesn't look good.
 259//   shapes = [for(i=[0:.2:1]) path3d(regular_ngon(n=4, side=4, rounding=i, $fn=32),i*5)];
 260//   skin(shapes, slices=0);
 261// Example: You can fix it by specifying "tangent" for the first method, but you still need "direct" for the rest.
 262//   shapes = [for(i=[0:.2:1]) path3d(regular_ngon(n=4, side=4, rounding=i, $fn=32),i*5)];
 263//   skin(shapes, slices=0, method=concat(["tangent"],repeat("direct",len(shapes)-2)));
 264// Example(FlatSpin,VPD=35): Connecting square to pentagon using "direct" method.
 265//   skin([regular_ngon(n=4, r=4), regular_ngon(n=5,r=5)], z=[0,4], refine=10, slices=10);
 266// Example(FlatSpin,VPD=35): Connecting square to shifted pentagon using "direct" method.
 267//   skin([regular_ngon(n=4, r=4), right(4,p=regular_ngon(n=5,r=5))], z=[0,4], refine=10, slices=10);
 268// Example(FlatSpin,VPD=185): In this example reindexing does not fix the orientation of the triangle because it happens in 3d within skin(), so we have to reverse the triangle manually
 269//   ellipse = yscale(3,circle(r=10, $fn=32));
 270//   tri = move([-50/3,-9],[[0,0], [50,0], [0,27]]);
 271//   skin([ellipse, reverse(tri)], z=[0,20], slices=20, method="reindex");
 272// Example(FlatSpin,VPD=185): You can get a nicer transition by rotating the polygons for better alignment.  You have to resample yourself before calling `align_polygon`. The orientation is fixed so we do not need to reverse.
 273//   ellipse = yscale(3,circle(r=10, $fn=32));
 274//   tri = move([-50/3,-9],
 275//              subdivide_path([[0,0], [50,0], [0,27]], 32));
 276//   aligned = align_polygon(ellipse,tri, [0:5:180]);
 277//   skin([ellipse, aligned], z=[0,20], slices=20);
 278// Example(FlatSpin,VPD=35): The "distance" method is a completely different approach.
 279//   skin([regular_ngon(n=4, r=4), regular_ngon(n=5,r=5)], z=[0,4], refine=10, slices=10, method="distance");
 280// Example(FlatSpin,VPD=35,VPT=[0,0,4]): Connecting pentagon to heptagon inserts two triangular faces on each side
 281//   small = path3d(circle(r=3, $fn=5));
 282//   big = up(2,p=yrot( 0,p=path3d(circle(r=3, $fn=7), 6)));
 283//   skin([small,big],method="distance", slices=10, refine=10);
 284// Example(FlatSpin,VPD=35,VPT=[0,0,4]): But just a slight rotation of the top profile moves the two triangles to one end
 285//   small = path3d(circle(r=3, $fn=5));
 286//   big = up(2,p=yrot(14,p=path3d(circle(r=3, $fn=7), 6)));
 287//   skin([small,big],method="distance", slices=10, refine=10);
 288// Example(FlatSpin,VPD=32,VPT=[1.2,4.3,2]): Another "distance" example:
 289//   off = [0,2];
 290//   shape = turtle(["right",45,"move", "left",45,"move", "left",45, "move", "jump", [.5+sqrt(2)/2,8]]);
 291//   rshape = rot(180,cp=centroid(shape)+off, p=shape);
 292//   skin([shape,rshape],z=[0,4], method="distance",slices=10,refine=15);
 293// Example(FlatSpin,VPD=32,VPT=[1.2,4.3,2]): Slightly shifting the profile changes the optimal linkage
 294//   off = [0,1];
 295//   shape = turtle(["right",45,"move", "left",45,"move", "left",45, "move", "jump", [.5+sqrt(2)/2,8]]);
 296//   rshape = rot(180,cp=centroid(shape)+off, p=shape);
 297//   skin([shape,rshape],z=[0,4], method="distance",slices=10,refine=15);
 298// Example(FlatSpin,VPD=444,VPT=[0,0,50]): This optimal solution doesn't look terrible:
 299//   prof1 = path3d([[-50,-50], [-50,50], [50,50], [25,25], [50,0], [25,-25], [50,-50]]);
 300//   prof2 = path3d(regular_ngon(n=7, r=50),100);
 301//   skin([prof1, prof2], method="distance", slices=10, refine=10);
 302// Example(FlatSpin,VPD=444,VPT=[0,0,50]): But this one looks better.  The "distance" method doesn't find it because it uses two more edges, so it clearly has a higher total edge distance.  We force it by doubling the first two vertices of one of the profiles.
 303//   prof1 = path3d([[-50,-50], [-50,50], [50,50], [25,25], [50,0], [25,-25], [50,-50]]);
 304//   prof2 = path3d(regular_ngon(n=7, r=50),100);
 305//   skin([repeat_entries(prof1,[2,2,1,1,1,1,1]),
 306//         prof2],
 307//        method="distance", slices=10, refine=10);
 308// Example(FlatSpin,VPD=80,VPT=[0,0,7]): The "distance" method will often produces results similar to the "tangent" method if you use it with a polygon and a curve, but the results can also look like this:
 309//   skin([path3d(circle($fn=128, r=10)), xrot(39, p=path3d(square([8,10]),10))],  method="distance", slices=0);
 310// Example(FlatSpin,VPD=80,VPT=[0,0,7]): Using the "tangent" method produces:
 311//   skin([path3d(circle($fn=128, r=10)), xrot(39, p=path3d(square([8,10]),10))],  method="tangent", slices=0);
 312// Example(FlatSpin,VPD=74): Torus using hexagons and pentagons, where `closed=true`
 313//   hex = right(7,p=path3d(hexagon(r=3)));
 314//   pent = right(7,p=path3d(pentagon(r=3)));
 315//   N=5;
 316//   skin(
 317//        [for(i=[0:2*N-1]) yrot(360*i/2/N, p=(i%2==0 ? hex : pent))],
 318//        refine=1,slices=0,method="distance",closed=true);
 319// Example: A smooth morph is achieved when you can calculate all the slices yourself.  Since you provide all the slices, set `slices=0`.
 320//   skin([for(n=[.1:.02:.5])
 321//            yrot(n*60-.5*60,p=path3d(supershape(step=360/128,m1=5,n1=n, n2=1.7),5-10*n))],
 322//        slices=0);
 323// Example: Another smooth supershape morph:
 324//   skin([for(alpha=[-.2:.05:1.5])
 325//            path3d(supershape(step=360/256,m1=7, n1=lerp(2,3,alpha),
 326//                              n2=lerp(8,4,alpha), n3=lerp(4,17,alpha)),alpha*5)],
 327//        slices=0);
 328// Example: Several polygons connected using "distance"
 329//   skin([regular_ngon(n=4, r=3),
 330//         regular_ngon(n=6, r=3),
 331//         regular_ngon(n=9, r=4),
 332//         rot(17,p=regular_ngon(n=6, r=3)),
 333//         rot(37,p=regular_ngon(n=4, r=3))],
 334//        z=[0,2,4,6,9], method="distance", slices=10, refine=10);
 335// Example(FlatSpin,VPD=935,VPT=[75,0,123]): Vertex count of the polygon changes at every profile
 336//   skin([
 337//       for (ang = [0:10:90])
 338//       rot([0,ang,0], cp=[200,0,0], p=path3d(circle(d=100,$fn=12-(ang/10))))
 339//   ],method="distance",slices=10,refine=10);
 340// Example: Möbius Strip.  This is a tricky model because when you work your way around to the connection, the direction of the profiles is flipped, so how can the proper geometry be created?  The trick is to duplicate the first profile and turn the caps off.  The model closes up and forms a valid polyhedron.
 341//   skin([
 342//     for (ang = [0:5:360])
 343//     rot([0,ang,0], cp=[100,0,0], p=rot(ang/2, p=path3d(square([1,30],center=true))))
 344//   ], caps=false, slices=0, refine=20);
 345// Example: This model of two scutoids packed together is based on https://www.thingiverse.com/thing:3024272 by mathgrrl
 346//   sidelen = 10;  // Side length of scutoid
 347//   height = 25;   // Height of scutoid
 348//   angle = -15;   // Angle (twists the entire form)
 349//   push = -5;     // Push (translates the base away from the top)
 350//   flare = 1;     // Flare (the two pieces will be different unless this is 1)
 351//   midpoint = .5; // Height of the extra vertex (as a fraction of total height); the two pieces will be different unless this is .5)
 352//   pushvec = rot(angle/2,p=push*RIGHT);  // Push direction is the average of the top and bottom mating edges
 353//   pent = path3d(apply(move(pushvec)*rot(angle),pentagon(side=sidelen,align_side=RIGHT,anchor="side0")));
 354//   hex = path3d(hexagon(side=flare*sidelen, align_side=RIGHT, anchor="side0"),height);
 355//   pentmate = path3d(pentagon(side=flare*sidelen,align_side=LEFT,anchor="side0"),height);
 356//             // Native index would require mapping first and last vertices together, which is not allowed, so shift
 357//   hexmate = list_rotate(
 358//                           path3d(apply(move(pushvec)*rot(angle),hexagon(side=sidelen,align_side=LEFT,anchor="side0"))),
 359//                           -1);
 360//   join_vertex = lerp(
 361//                       mean(select(hex,1,2)),     // midpoint of "extra" hex edge
 362//                       mean(select(hexmate,0,1)), // midpoint of "extra" hexmate edge
 363//                       midpoint);
 364//   augpent = repeat_entries(pent, [1,2,1,1,1]);         // Vertex 1 will split at the top forming a triangular face with the hexagon
 365//   augpent_mate = repeat_entries(pentmate,[2,1,1,1,1]); // For mating pentagon it is vertex 0 that splits
 366//              // Middle is the interpolation between top and bottom except for the join vertex, which is doubled because it splits
 367//   middle = list_set(lerp(augpent,hex,midpoint),[1,2],[join_vertex,join_vertex]);
 368//   middle_mate = list_set(lerp(hexmate,augpent_mate,midpoint), [0,1], [join_vertex,join_vertex]);
 369//   skin([augpent,middle,hex],  slices=10, refine=10, sampling="segment");
 370//   color("green")skin([augpent_mate,middle_mate,hexmate],  slices=10,refine=10, sampling="segment");
 371// Example: If you create a self-intersecting polyhedron the result is invalid.  In some cases self-intersection may be obvous.  Here is a more subtle example.
 372//   skin([
 373//          for (a = [0:30:180]) let(
 374//              pos  = [-60*sin(a),     0, a    ],
 375//              pos2 = [-60*sin(a+0.1), 0, a+0.1]
 376//          ) move(pos,
 377//              p=rot(from=UP, to=pos2-pos,
 378//                  p=path3d(circle(d=150))
 379//              )
 380//          )
 381//      ],refine=1,slices=0);
 382//      color("red") {
 383//          zrot(25) fwd(130) xrot(75) {
 384//              linear_extrude(height=0.1) {
 385//                  ydistribute(25) {
 386//                      text(text="BAD POLYHEDRONS!", size=20, halign="center", valign="center");
 387//                      text(text="CREASES MAKE", size=20, halign="center", valign="center");
 388//                  }
 389//              }
 390//          }
 391//          up(160) zrot(25) fwd(130) xrot(75) {
 392//              stroke(zrot(30, p=yscale(0.5, p=circle(d=120))),width=10,closed=true);
 393//          }
 394//      }
 395module skin(profiles, slices, refine=1, method="direct", sampling, caps, closed=false, z, style="min_edge", convexity=10,
 396            anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull")
 397{
 398    vnf = skin(profiles, slices, refine, method, sampling, caps, closed, z, style=style);
 399    vnf_polyhedron(vnf,convexity=convexity,spin=spin,anchor=anchor,orient=orient,atype=atype,cp=cp)
 400        children();
 401}
 402
 403
 404function skin(profiles, slices, refine=1, method="direct", sampling, caps, closed=false, z, style="min_edge",
 405              anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull") =
 406  assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"")
 407  assert(is_def(slices),"The slices argument must be specified.")
 408  assert(is_list(profiles) && len(profiles)>1, "Must provide at least two profiles")
 409  let(
 410       profiles = [for(p=profiles) if (is_region(p) && len(p)==1) p[0] else p]
 411  )
 412  let( bad = [for(i=idx(profiles)) if (!(is_path(profiles[i]) && len(profiles[i])>2)) i])
 413  assert(len(bad)==0, str("Profiles ",bad," are not a paths or have length less than 3"))
 414  let(
 415    profcount = len(profiles) - (closed?0:1),
 416    legal_methods = ["direct","reindex","distance","fast_distance","tangent"],
 417    caps = is_def(caps) ? caps :
 418           closed ? false : true,
 419    capsOK = is_bool(caps) || is_bool_list(caps,2),
 420    fullcaps = is_bool(caps) ? [caps,caps] : caps,
 421    refine = is_list(refine) ? refine : repeat(refine, len(profiles)),
 422    slices = is_list(slices) ? slices : repeat(slices, profcount),
 423    refineOK = [for(i=idx(refine)) if (refine[i]<=0 || !is_integer(refine[i])) i],
 424    slicesOK = [for(i=idx(slices)) if (!is_integer(slices[i]) || slices[i]<0) i],
 425    maxsize = max_length(profiles),
 426    methodok = is_list(method) || in_list(method, legal_methods),
 427    methodlistok = is_list(method) ? [for(i=idx(method)) if (!in_list(method[i], legal_methods)) i] : [],
 428    method = is_string(method) ? repeat(method, profcount) : method,
 429    // Define to be zero where a resampling method is used and 1 where a vertex duplicator is used
 430    RESAMPLING = 0,
 431    DUPLICATOR = 1,
 432    method_type = [for(m = method) m=="direct" || m=="reindex" ? 0 : 1],
 433    sampling = is_def(sampling) ? sampling :
 434               in_list(DUPLICATOR,method_type) ? "segment" : "length"
 435  )
 436  assert(len(refine)==len(profiles), "refine list is the wrong length")
 437  assert(len(slices)==profcount, str("slices list must have length ",profcount))
 438  assert(slicesOK==[],str("slices must be nonnegative integers"))
 439  assert(refineOK==[],str("refine must be postive integer"))
 440  assert(methodok,str("method must be one of ",legal_methods,". Got ",method))
 441  assert(methodlistok==[], str("method list contains invalid method at ",methodlistok))
 442  assert(len(method) == profcount,"Method list is the wrong length")
 443  assert(in_list(sampling,["length","segment"]), "sampling must be set to \"length\" or \"segment\"")
 444  assert(sampling=="segment" || (!in_list("distance",method) && !in_list("fast_distance",method) && !in_list("tangent",method)), "sampling is set to \"length\" which is only allowed with methods \"direct\" and \"reindex\"")
 445  assert(capsOK, "caps must be boolean or a list of two booleans")
 446  assert(!closed || !caps, "Cannot make closed shape with caps")
 447  let(
 448    profile_dim=list_shape(profiles,2),
 449    profiles_zcheck = (profile_dim != 2) || (profile_dim==2 && is_list(z) && len(z)==len(profiles)),
 450    profiles_ok = (profile_dim==2 && is_list(z) && len(z)==len(profiles)) || profile_dim==3
 451  )
 452  assert(profiles_zcheck, "z parameter is invalid or has the wrong length.")
 453  assert(profiles_ok,"Profiles must all be 3d or must all be 2d, with matching length z parameter.")
 454  assert(is_undef(z) || profile_dim==2, "Do not specify z with 3d profiles")
 455  assert(profile_dim==3 || len(z)==len(profiles),"Length of z does not match length of profiles.")
 456  let(
 457    // Adjoin Z coordinates to 2d profiles
 458    profiles = profile_dim==3 ? profiles :
 459               [for(i=idx(profiles)) path3d(profiles[i], z[i])],
 460    // True length (not counting repeated vertices) of profiles after refinement
 461    refined_len = [for(i=idx(profiles)) refine[i]*len(profiles[i])],
 462    // Define this to be 1 if a profile is used on either side by a resampling method, zero otherwise.
 463    profile_resampled = [for(i=idx(profiles))
 464      1-(
 465           i==0 ?  method_type[0] * (closed? last(method_type) : 1) :
 466           i==len(profiles)-1 ? last(method_type) * (closed ? select(method_type,-2) : 1) :
 467         method_type[i] * method_type[i-1])],
 468    parts = search(1,[1,for(i=[0:1:len(profile_resampled)-2]) profile_resampled[i]!=profile_resampled[i+1] ? 1 : 0],0),
 469    plen = [for(i=idx(parts)) (i== len(parts)-1? len(refined_len) : parts[i+1]) - parts[i]],
 470    max_list = [for(i=idx(parts)) each repeat(max(select(refined_len, parts[i], parts[i]+plen[i]-1)), plen[i])],
 471    transition_profiles = [for(i=[(closed?0:1):1:profcount-1]) if (select(method_type,i-1) != method_type[i]) i],
 472    badind = [for(tranprof=transition_profiles) if (refined_len[tranprof] != max_list[tranprof]) tranprof]
 473  )
 474  assert(badind==[],str("Profile length mismatch at method transition at indices ",badind," in skin()"))
 475  let(
 476    full_list =    // If there are no duplicators then use more efficient where the whole input is treated together
 477      !in_list(DUPLICATOR,method_type) ?
 478         let(
 479             resampled = [for(i=idx(profiles)) subdivide_path(profiles[i], max_list[i], method=sampling)],
 480             fixedprof = [for(i=idx(profiles))
 481                             i==0 || method[i-1]=="direct" ? resampled[i]
 482                                                         : reindex_polygon(resampled[i-1],resampled[i])],
 483             sliced = slice_profiles(fixedprof, slices, closed)
 484            )
 485            [!closed ? sliced : concat(sliced,[sliced[0]])]
 486      :  // There are duplicators, so use approach where each pair is treated separately
 487      [for(i=[0:profcount-1])
 488        let(
 489          pair =
 490            method[i]=="distance" ? _skin_distance_match(profiles[i],select(profiles,i+1)) :
 491            method[i]=="fast_distance" ? _skin_aligned_distance_match(profiles[i], select(profiles,i+1)) :
 492            method[i]=="tangent" ? _skin_tangent_match(profiles[i],select(profiles,i+1)) :
 493            /*method[i]=="reindex" || method[i]=="direct" ?*/
 494               let( p1 = subdivide_path(profiles[i],max_list[i], method=sampling),
 495                    p2 = subdivide_path(select(profiles,i+1),max_list[i], method=sampling)
 496               ) (method[i]=="direct" ? [p1,p2] : [p1, reindex_polygon(p1, p2)]),
 497            nsamples =  method_type[i]==RESAMPLING ? len(pair[0]) :
 498               assert(refine[i]==select(refine,i+1),str("Refine value mismatch at indices ",[i,(i+1)%len(refine)],
 499                                                        ".  Method ",method[i]," requires equal values"))
 500               refine[i] * len(pair[0])
 501          )
 502          subdivide_and_slice(pair,slices[i], nsamples, method=sampling)],
 503      vnf=vnf_join(
 504          [for(i=idx(full_list))
 505              vnf_vertex_array(full_list[i], cap1=i==0 && fullcaps[0], cap2=i==len(full_list)-1 && fullcaps[1],
 506                               col_wrap=true, style=style)])
 507  )
 508  reorient(anchor,spin,orient,vnf=vnf,p=vnf,extent=atype=="hull",cp=cp);
 509
 510
 511
 512// Function&Module: linear_sweep()
 513// Synopsis: Create a linear extrusion from a path, with optional texturing. 
 514// SynTags: VNF, Geom
 515// Topics: Extrusion, Textures, Sweep
 516// See Also: rotate_sweep(), sweep(), spiral_sweep(), path_sweep(), offset_sweep()
 517// Usage: As Module
 518//   linear_sweep(region, [height], [center=], [slices=], [twist=], [scale=], [style=], [caps=], [convexity=]) [ATTACHMENTS];
 519// Usage: With Texturing
 520//   linear_sweep(region, [height], [center=], texture=, [tex_size=]|[tex_counts=], [tex_scale=], [style=], [tex_samples=], ...) [ATTACHMENTS];
 521// Usage: As Function
 522//   vnf = linear_sweep(region, [height], [center=], [slices=], [twist=], [scale=], [style=], [caps=]);
 523//   vnf = linear_sweep(region, [height], [center=], texture=, [tex_size=]|[tex_counts=], [tex_scale=], [style=], [tex_samples=], ...);
 524// Description:
 525//   If called as a module, creates a polyhedron that is the linear extrusion of the given 2D region or polygon.
 526//   If called as a function, returns a VNF that can be used to generate a polyhedron of the linear extrusion
 527//   of the given 2D region or polygon.  The benefit of using this, over using `linear_extrude region(rgn)` is
 528//   that it supports `anchor`, `spin`, `orient` and attachments.  You can also make more refined
 529//   twisted extrusions by using `maxseg` to subsample flat faces.
 530//   .
 531//   Anchoring for linear_sweep is based on the anchors for the swept region rather than from the polyhedron that is created.  This can produce more
 532//   predictable anchors for LEFT, RIGHT, FWD and BACK in many cases, but the anchors may only
 533//   be aproximately correct for twisted objects, and corner anchors may point in unexpected directions in some cases.
 534//   If you need anchors directly computed from the surface you can pass the vnf from linear_sweep
 535//   to {{vnf_polyhedron()}}, which will compute anchors directly from the full VNF.  
 536// Arguments:
 537//   region = The 2D [Region](regions.scad) or polygon that is to be extruded.
 538//   h / height / l / length = The height to extrude the region.  Default: 1
 539//   center = If true, the created polyhedron will be vertically centered.  If false, it will be extruded upwards from the XY plane.  Default: `false`
 540//   ---
 541//   twist = The number of degrees to rotate the top of the shape, clockwise around the Z axis, relative to the bottom.  Default: 0
 542//   scale = The amount to scale the top of the shape, in the X and Y directions, relative to the size of the bottom.  Default: 1
 543//   shift = The amount to shift the top of the shape, in the X and Y directions, relative to the position of the bottom.  Default: [0,0]
 544//   slices = The number of slices to divide the shape into along the Z axis, to allow refinement of detail, especially when working with a twist.  Default: `twist/5`
 545//   maxseg = If given, then any long segments of the region will be subdivided to be shorter than this length.  This can refine twisting flat faces a lot.  Default: `undef` (no subsampling)
 546//   texture = A texture name string, or a rectangular array of scalar height values (0.0 to 1.0), or a VNF tile that defines the texture to apply to vertical surfaces.  See {{texture()}} for what named textures are supported.
 547//   tex_size = An optional 2D target size for the textures.  Actual texture sizes will be scaled somewhat to evenly fit the available surface. Default: `[5,5]`
 548//   tex_counts = If given instead of tex_size, gives the tile repetition counts for textures over the surface length and height.
 549//   tex_inset = If numeric, lowers the texture into the surface by that amount, before the tex_scale multiplier is applied.  If `true`, insets by exactly `1`.  Default: `false`
 550//   tex_rot = If true, rotates the texture 90º.
 551//   tex_scale = Scaling multiplier for the texture depth.
 552//   tex_samples = Minimum number of "bend points" to have in VNF texture tiles.  Default: 8
 553//   style = The style to use when triangulating the surface of the object.  Valid values are `"default"`, `"alt"`, or `"quincunx"`.
 554//   caps = If false do not create end caps.  Can be a boolean vector.  Default: true
 555//   convexity = Max number of surfaces any single ray could pass through.  Module use only.
 556//   cp = Centerpoint for determining intersection anchors or centering the shape.  Determines the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: `"centroid"`
 557//   atype = Set to "hull" or "intersect" to select anchor type.  Default: "hull"
 558//   anchor = Translate so anchor point is at origin (0,0,0).  See [anchor](attachments.scad#subsection-anchor).  Default: `"origin"`
 559//   spin = Rotate this many degrees around the Z axis after anchor.  See [spin](attachments.scad#subsection-spin).  Default: `0`
 560//   orient = Vector to rotate top towards, after spin.  See [orient](attachments.scad#subsection-orient).  Default: `UP`
 561// Anchor Types:
 562//   "hull" = Anchors to the virtual convex hull of the shape.
 563//   "intersect" = Anchors to the surface of the shape.
 564//   "bbox" = Anchors to the bounding box of the extruded shape.
 565// Extra Anchors:
 566//   "origin" = Centers the extruded shape vertically only, but keeps the original path positions in the X and Y.  Oriented UP.
 567//   "original_base" = Keeps the original path positions in the X and Y, but at the bottom of the extrusion.  Oriented UP.
 568// Example: Extruding a Compound Region.
 569//   rgn1 = [for (d=[10:10:60]) circle(d=d,$fn=8)];
 570//   rgn2 = [square(30,center=false)];
 571//   rgn3 = [for (size=[10:10:20]) move([15,15],p=square(size=size, center=true))];
 572//   mrgn = union(rgn1,rgn2);
 573//   orgn = difference(mrgn,rgn3);
 574//   linear_sweep(orgn,height=20,convexity=16);
 575// Example: With Twist, Scale, Shift, Slices and Maxseg.
 576//   rgn1 = [for (d=[10:10:60]) circle(d=d,$fn=8)];
 577//   rgn2 = [square(30,center=false)];
 578//   rgn3 = [
 579//       for (size=[10:10:20])
 580//       apply(
 581//          move([15,15]),
 582//          square(size=size, center=true)
 583//       )
 584//   ];
 585//   mrgn = union(rgn1,rgn2);
 586//   orgn = difference(mrgn,rgn3);
 587//   linear_sweep(
 588//       orgn, height=50, maxseg=2, slices=40,
 589//       twist=90, scale=0.5, shift=[10,5],
 590//       convexity=16
 591//   );
 592// Example: Anchors on an Extruded Region
 593//   rgn1 = [for (d=[10:10:60]) circle(d=d,$fn=8)];
 594//   rgn2 = [square(30,center=false)];
 595//   rgn3 = [
 596//       for (size=[10:10:20])
 597//       apply(
 598//           move([15,15]),
 599//           rect(size=size)
 600//       )
 601//   ];
 602//   mrgn = union(rgn1,rgn2);
 603//   orgn = difference(mrgn,rgn3);
 604//   linear_sweep(orgn,height=20,convexity=16)
 605//       show_anchors();
 606// Example: "diamonds" texture.
 607//   path = glued_circles(r=15, spread=40, tangent=45);
 608//   linear_sweep(
 609//       path, texture="diamonds", tex_size=[5,10],
 610//       h=40, style="concave");
 611// Example: "pyramids" texture.
 612//   linear_sweep(
 613//       rect(50), texture="pyramids", tex_size=[10,10],
 614//       h=40, style="convex");
 615// Example: "bricks_vnf" texture.
 616//   path = glued_circles(r=15, spread=40, tangent=45);
 617//   linear_sweep(
 618//       path, texture="bricks_vnf", tex_size=[10,10],
 619//       tex_scale=0.25, h=40);
 620// Example: User defined heightfield texture.
 621//   path = ellipse(r=[20,10]);
 622//   texture = [for (i=[0:9])
 623//       [for (j=[0:9])
 624//           1/max(0.5,norm([i,j]-[5,5])) ]];
 625//   linear_sweep(
 626//       path, texture=texture, tex_size=[5,5],
 627//       h=40, style="min_edge", anchor=BOT);
 628// Example: User defined VNF tile texture.
 629//   path = ellipse(r=[20,10]);
 630//   tex = let(n=16,m=0.25) [
 631//        [
 632//            each resample_path(path3d(square(1)),n),
 633//            each move([0.5,0.5],
 634//                p=path3d(circle(d=0.5,$fn=n),m)),
 635//            [1/2,1/2,0],
 636//        ], [
 637//            for (i=[0:1:n-1]) each [
 638//                [i,(i+1)%n,(i+3)%n+n],
 639//                [i,(i+3)%n+n,(i+2)%n+n],
 640//                [2*n,n+i,n+(i+1)%n],
 641//            ]
 642//        ]
 643//   ];
 644//   linear_sweep(path, texture=tex, tex_size=[5,5], h=40);
 645// Example: Textured with twist and scale.
 646//   linear_sweep(regular_ngon(n=3, d=50),
 647//       texture="rough", h=100, tex_scale=2,
 648//       tex_size=[20,20], style="min_edge",
 649//       convexity=10, scale=0.2, twist=120);
 650// Example: As Function
 651//   path = glued_circles(r=15, spread=40, tangent=45);
 652//   vnf = linear_sweep(
 653//       path, h=40, texture="trunc_pyramids", tex_size=[5,5],
 654//       tex_scale=1, style="convex");
 655//   vnf_polyhedron(vnf, convexity=10);
 656// Example: VNF tile that has no top/bottom edges and produces a disconnected result
 657//   shape = skin([rect(2/5),
 658//                 rect(2/3),
 659//                 rect(2/5)],
 660//                z=[0,1/2,1],
 661//                slices=0,
 662//                caps=false);
 663//   tile = move([0,1/2,2/3],yrot(90,shape));
 664//   linear_sweep(circle(20), texture=tile,
 665//                tex_size=[10,10],tex_scale=5,
 666//                h=40,convexity=4);
 667// Example: The same tile from above, turned 90 degrees, creates problems at the ends, because the end cap is not a connected polygon.  When the ends are disconnected you may find that some parts of the end cap are missing and spurious polygons included.  
 668//  shape = skin([rect(2/5),
 669//                rect(2/3),
 670//                rect(2/5)],
 671//               z=[0,1/2,1],
 672//               slices=0,
 673//               caps=false);
 674//  tile = move([1/2,1,2/3],xrot(90,shape));
 675//  linear_sweep(circle(20), texture=tile,
 676//               tex_size=[30,20],tex_scale=15,
 677//               h=40,convexity=4);
 678// Example: This example shows some endcap polygons missing and a spurious triangle
 679//   shape = skin([rect(2/5),
 680//                 rect(2/3),
 681//                 rect(2/5)],
 682//                z=[0,1/2,1],
 683//                slices=0,
 684//                caps=false);
 685//   tile = xscale(.5,move([1/2,1,2/3],xrot(90,shape)));
 686//   doubletile = vnf_join([tile, right(.5,tile)]);
 687//   linear_sweep(circle(20), texture=doubletile,
 688//                tex_size=[45,45],tex_scale=15, h=40);
 689// Example: You can fix ends for disconnected cases using {{top_half()}} and {{bottom_half()}}
 690//   shape = skin([rect(2/5),
 691//                 rect(2/3),
 692//                 rect(2/5)],
 693//                z=[0,1/2,1],
 694//                slices=0,
 695//                caps=false);
 696//   tile = move([1/2,1,2/3],xrot(90,shape));
 697//   vnf_polyhedron(
 698//     top_half(
 699//       bottom_half(
 700//         linear_sweep(circle(20), texture=tile,
 701//                     tex_size=[30,20],tex_scale=15,
 702//                     h=40.2,caps=false),
 703//       z=20),
 704//     z=-20)); 
 705
 706module linear_sweep(
 707    region, height, center,
 708    twist=0, scale=1, shift=[0,0],
 709    slices, maxseg, style="default", convexity, caps=true, 
 710    texture, tex_size=[5,5], tex_counts,
 711    tex_inset=false, tex_rot=false,
 712    tex_scale=1, tex_samples,
 713    cp, atype="hull", h,l,length,
 714    anchor, spin=0, orient=UP
 715) {
 716    h = one_defined([h, height,l,length],"h,height,l,length",dflt=1);
 717    region = force_region(region);
 718    check = assert(is_region(region),"Input is not a region");
 719    anchor = center==true? "origin" :
 720        center == false? "original_base" :
 721        default(anchor, "original_base");
 722    vnf = linear_sweep(
 723        region, height=h, style=style, caps=caps, 
 724        twist=twist, scale=scale, shift=shift,
 725        texture=texture,
 726        tex_size=tex_size,
 727        tex_counts=tex_counts,
 728        tex_inset=tex_inset,
 729        tex_rot=tex_rot,
 730        tex_scale=tex_scale,
 731        tex_samples=tex_samples,
 732        slices=slices,
 733        maxseg=maxseg,
 734        anchor="origin"
 735    );
 736    anchors = [
 737        named_anchor("original_base", [0,0,-h/2], UP)
 738    ];
 739    cp = default(cp, "centroid");
 740    geom = atype=="hull"?  attach_geom(cp=cp, region=region, h=h, extent=true, shift=shift, scale=scale, twist=twist, anchors=anchors) :
 741        atype=="intersect"?  attach_geom(cp=cp, region=region, h=h, extent=false, shift=shift, scale=scale, twist=twist, anchors=anchors) :
 742        atype=="bbox"?
 743            let(
 744                bounds = pointlist_bounds(flatten(region)),
 745                size = bounds[1] - bounds[0],
 746                midpt = (bounds[0] + bounds[1])/2
 747            )
 748            attach_geom(cp=[0,0,0], size=point3d(size,h), offset=point3d(midpt), shift=shift, scale=scale, twist=twist, anchors=anchors) :
 749        assert(in_list(atype, ["hull","intersect","bbox"]), "Anchor type must be \"hull\", \"intersect\", or \"bbox\".");
 750    attachable(anchor,spin,orient, geom=geom) {
 751        vnf_polyhedron(vnf, convexity=convexity);
 752        children();
 753    }
 754}
 755
 756
 757function linear_sweep(
 758    region, height, center,
 759    twist=0, scale=1, shift=[0,0],
 760    slices, maxseg, style="default", caps=true, 
 761    cp, atype="hull", h,
 762    texture, tex_size=[5,5], tex_counts,
 763    tex_inset=false, tex_rot=false,
 764    tex_scale=1, tex_samples, h, l, length, 
 765    anchor, spin=0, orient=UP
 766) =
 767    let( region = force_region(region) )
 768    assert(is_region(region), "Input is not a region or polygon.")
 769    assert(is_num(scale) || is_vector(scale))
 770    assert(is_vector(shift, 2), str(shift))
 771    assert(is_bool(caps) || is_bool_list(caps,2), "caps must be boolean or a list of two booleans")
 772    let(
 773        h = one_defined([h, height,l,length],"h,height,l,length",dflt=1)
 774    )
 775    !is_undef(texture)? _textured_linear_sweep(
 776        region, h=h, caps=caps, 
 777        texture=texture, tex_size=tex_size,
 778        counts=tex_counts, inset=tex_inset,
 779        rot=tex_rot, tex_scale=tex_scale,
 780        twist=twist, scale=scale, shift=shift,
 781        style=style, samples=tex_samples,
 782        anchor=anchor, spin=spin, orient=orient
 783    ) :
 784    let(
 785        caps = is_bool(caps) ? [caps,caps] : caps, 
 786        anchor = center==true? "origin" :
 787            center == false? "original_base" :
 788            default(anchor, "original_base"),
 789        regions = region_parts(region),
 790        slices = default(slices, max(1,ceil(abs(twist)/5))),
 791        scale = is_num(scale)? [scale,scale] : point2d(scale),
 792        topmat = move(shift) * scale(scale) * rot(-twist),
 793        trgns = [
 794            for (rgn = regions) [
 795                for (path = rgn) let(
 796                    p = list_unwrap(path),
 797                    path = is_undef(maxseg)? p : [
 798                        for (seg = pair(p,true)) each
 799                        let( steps = ceil(norm(seg.y - seg.x) / maxseg) )
 800                        lerpn(seg.x, seg.y, steps, false)
 801                    ]
 802                ) apply(topmat, path)
 803            ]
 804        ],
 805        vnf = vnf_join([
 806            for (rgn = regions)
 807            for (pathnum = idx(rgn)) let(
 808                p = list_unwrap(rgn[pathnum]),
 809                path = is_undef(maxseg)? p : [
 810                    for (seg=pair(p,true)) each
 811                    let(steps=ceil(norm(seg.y-seg.x)/maxseg))
 812                    lerpn(seg.x, seg.y, steps, false)
 813                ],
 814                verts = [
 815                    for (i=[0:1:slices]) let(
 816                        u = i / slices,
 817                        scl = lerp([1,1], scale, u),
 818                        ang = lerp(0, -twist, u),
 819                        off = lerp([0,0,-h/2], point3d(shift,h/2), u),
 820                        m = move(off) * scale(scl) * rot(ang)
 821                    ) apply(m, path3d(path))
 822                ]
 823            ) vnf_vertex_array(verts, caps=false, col_wrap=true, style=style),
 824            if (caps[0]) for (rgn = regions) vnf_from_region(rgn, down(h/2), reverse=true),
 825            if (caps[1]) for (rgn = trgns) vnf_from_region(rgn, up(h/2), reverse=false)
 826        ]),
 827        anchors = [
 828            named_anchor("original_base", [0,0,-h/2], UP)
 829        ],
 830        cp = default(cp, "centroid"),
 831        geom = atype=="hull"?  attach_geom(cp=cp, region=region, h=h, extent=true, shift=shift, scale=scale, twist=twist, anchors=anchors) :
 832            atype=="intersect"?  attach_geom(cp=cp, region=region, h=h, extent=false, shift=shift, scale=scale, twist=twist, anchors=anchors) :
 833            atype=="bbox"?
 834                let(
 835                    bounds = pointlist_bounds(flatten(region)),
 836                    size = bounds[1] - bounds[0],
 837                    midpt = (bounds[0] + bounds[1])/2
 838                )
 839                attach_geom(cp=[0,0,0], size=point3d(size,h), offset=point3d(midpt), shift=shift, scale=scale, twist=twist, anchors=anchors) :
 840            assert(in_list(atype, ["hull","intersect","bbox"]), "Anchor type must be \"hull\", \"intersect\", or \"bbox\".")
 841    ) reorient(anchor,spin,orient, geom=geom, p=vnf);
 842
 843
 844// Function&Module: rotate_sweep()
 845// Synopsis: Create a surface of revolution from a path with optional texturing. 
 846// SynTags: VNF, Geom
 847// Topics: Extrusion, Sweep, Revolution, Textures
 848// See Also: linear_sweep(), sweep(), spiral_sweep(), path_sweep(), offset_sweep()
 849// Usage: As Function
 850//   vnf = rotate_sweep(shape, [angle], ...);
 851// Usage: As Module
 852//   rotate_sweep(shape, [angle], ...) [ATTACHMENTS];
 853// Usage: With Texturing
 854//   rotate_sweep(shape, texture=, [tex_size=]|[tex_counts=], [tex_scale=], [tex_samples=], [tex_rot=], [tex_inset=], ...) [ATTACHMENTS];
 855// Description:
 856//   Takes a polygon or [region](regions.scad) and sweeps it in a rotation around the Z axis, with optional texturing.
 857//   When called as a function, returns a [VNF](vnf.scad).
 858//   When called as a module, creates the sweep as geometry.
 859// Arguments:
 860//   shape = The polygon or [region](regions.scad) to sweep around the Z axis.
 861//   angle = If given, specifies the number of degrees to sweep the shape around the Z axis, counterclockwise from the X+ axis.  Default: 360 (full rotation)
 862//   ---
 863//   texture = A texture name string, or a rectangular array of scalar height values (0.0 to 1.0), or a VNF tile that defines the texture to apply to vertical surfaces.  See {{texture()}} for what named textures are supported.
 864//   tex_size = An optional 2D target size for the textures.  Actual texture sizes will be scaled somewhat to evenly fit the available surface. Default: `[5,5]`
 865//   tex_counts = If given instead of tex_size, gives the tile repetition counts for textures over the surface length and height.
 866//   tex_inset = If numeric, lowers the texture into the surface by that amount, before the tex_scale multiplier is applied.  If `true`, insets by exactly `1`.  Default: `false`
 867//   tex_rot = If true, rotates the texture 90º.
 868//   tex_scale = Scaling multiplier for the texture depth.
 869//   tex_samples = Minimum number of "bend points" to have in VNF texture tiles.  Default: 8
 870//   style = {{vnf_vertex_array()}} style.  Default: "min_edge"
 871//   closed = If false, and shape is given as a path, then the revolved path will be sealed to the axis of rotation with untextured caps.  Default: `true`
 872//   convexity = (Module only) Convexity setting for use with polyhedron.  Default: 10
 873//   cp = Centerpoint for determining "intersect" anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: "centroid"
 874//   atype = Select "hull" or "intersect" anchor types.  Default: "hull"
 875//   anchor = Translate so anchor point is at the origin. Default: "origin"
 876//   spin = Rotate this many degrees around Z axis after anchor. Default: 0
 877//   orient = Vector to rotate top towards after spin  (module only)
 878// Anchor Types:
 879//   "hull" = Anchors to the virtual convex hull of the shape.
 880//   "intersect" = Anchors to the surface of the shape.
 881// Example:
 882//   rgn = [
 883//       for (a = [0, 120, 240]) let(
 884//           cp = polar_to_xy(15, a) + [30,0]
 885//       ) each [
 886//           move(cp, p=circle(r=10)),
 887//           move(cp, p=hexagon(d=15)),
 888//       ]
 889//   ];
 890//   rotate_sweep(rgn, angle=240);
 891// Example:
 892//   rgn = right(30, p=union([for (a = [0, 90]) rot(a, p=rect([15,5]))]));
 893//   rotate_sweep(rgn);
 894// Example:
 895//   path = right(50, p=circle(d=40));
 896//   rotate_sweep(path, texture="bricks_vnf", tex_size=[10,10], tex_scale=0.5, style="concave");
 897// Example:
 898//   tex = [
 899//       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
 900//       [0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1],
 901//       [0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1],
 902//       [0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1],
 903//       [0, 0, 1, 0, 0, 1, 1, 1, 0, 0, 1],
 904//       [0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1],
 905//       [0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1],
 906//       [0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1],
 907//       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1],
 908//       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1],
 909//       [0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1],
 910//       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
 911//   ];
 912//   path = arc(cp=[0,0], r=40, start=60, angle=-120);
 913//   rotate_sweep(
 914//       path, closed=false,
 915//       texture=tex, tex_size=[20,20],
 916//       tex_scale=1, style="concave");
 917// Example:
 918//   include <BOSL2/beziers.scad>
 919//   bezpath = [
 920//       [15, 30], [10,15],
 921//       [10,  0], [20, 10], [30,12],
 922//       [30,-12], [20,-10], [10, 0],
 923//       [10,-15], [15,-30]
 924//   ];
 925//   path = bezpath_curve(bezpath, splinesteps=32);
 926//   rotate_sweep(
 927//       path, closed=false,
 928//       texture="diamonds", tex_size=[10,10],
 929//       tex_scale=1, style="concave");
 930// Example:
 931//   path = [
 932//       [20, 30], [20, 20],
 933//       each arc(r=20, corner=[[20,20],[10,0],[20,-20]]),
 934//       [20,-20], [20,-30],
 935//   ];
 936//   vnf = rotate_sweep(
 937//       path, closed=false,
 938//       texture="trunc_pyramids",
 939//       tex_size=[5,5], tex_scale=1,
 940//       style="convex");
 941//   vnf_polyhedron(vnf, convexity=10);
 942// Example:
 943//   rgn = [
 944//       right(40, p=circle(d=50)),
 945//       right(40, p=circle(d=40,$fn=6)),
 946//   ];
 947//   rotate_sweep(
 948//       rgn, texture="diamonds",
 949//       tex_size=[10,10], tex_scale=1,
 950//       angle=240, style="concave");
 951
 952function rotate_sweep(
 953    shape, angle=360,
 954    texture, tex_size=[5,5], tex_counts,
 955    tex_inset=false, tex_rot=false,
 956    tex_scale=1, tex_samples,
 957    tex_taper, shift=[0,0], closed=true,
 958    style="min_edge", cp="centroid",
 959    atype="hull", anchor="origin",
 960    spin=0, orient=UP
 961) =
 962    let( region = force_region(shape) )
 963    assert(is_region(region), "Input is not a region or polygon.")
 964    let(
 965        bounds = pointlist_bounds(flatten(region)),
 966        min_x = bounds[0].x,
 967        max_x = bounds[1].x,
 968        min_y = bounds[0].y,
 969        max_y = bounds[1].y,
 970        h = max_y - min_y
 971    )
 972    assert(min_x>=0, "Input region must exist entirely in the X+ half-plane.")
 973    !is_undef(texture)? _textured_revolution(
 974        shape,
 975        texture=texture,
 976        tex_size=tex_size,
 977        counts=tex_counts,
 978        tex_scale=tex_scale,
 979        inset=tex_inset,
 980        rot=tex_rot,
 981        samples=tex_samples,
 982        taper=tex_taper,
 983        shift=shift,
 984        closed=closed,
 985        angle=angle,
 986        style=style
 987    ) :
 988    let(
 989        steps = segs(max_x),
 990        skmat = down(min_y) * skew(sxz=shift.x/h, syz=shift.y/h) * up(min_y),
 991        transforms = [
 992            if (angle==360) for (i=[0:1:steps-1]) skmat * rot([90,0,360-i*360/steps]),
 993            if (angle<360) for (i=[0:1:steps-1]) skmat * rot([90,0,angle-i*angle/(steps-1)]),
 994        ],
 995        vnf = sweep(
 996            region, transforms,
 997            closed=angle==360,
 998            caps=angle!=360,
 999            style=style, cp=cp,
1000            atype=atype, anchor=anchor,
1001            spin=spin, orient=orient
1002        )
1003    ) vnf;
1004
1005
1006module rotate_sweep(
1007    shape, angle=360,
1008    texture, tex_size=[5,5], tex_counts,
1009    tex_inset=false, tex_rot=false,
1010    tex_scale=1, tex_samples,
1011    tex_taper, shift=[0,0],
1012    style="min_edge",
1013    closed=true,
1014    cp="centroid",
1015    convexity=10,
1016    atype="hull",
1017    anchor="origin",
1018    spin=0,
1019    orient=UP
1020) {
1021    region = force_region(shape);
1022    check = assert(is_region(region), "Input is not a region or polygon.");
1023    bounds = pointlist_bounds(flatten(region));
1024    min_x = bounds[0].x;
1025    max_x = bounds[1].x;
1026    min_y = bounds[0].y;
1027    max_y = bounds[1].y;
1028    h = max_y - min_y;
1029    check2 = assert(min_x>=0, "Input region must exist entirely in the X+ half-plane.");
1030    steps = segs(max_x);
1031    if (!is_undef(texture)) {
1032        _textured_revolution(
1033            shape,
1034            texture=texture,
1035            tex_size=tex_size,
1036            counts=tex_counts,
1037            tex_scale=tex_scale,
1038            inset=tex_inset,
1039            rot=tex_rot,
1040            samples=tex_samples,
1041            taper=tex_taper,
1042            shift=shift,
1043            closed=closed,
1044            angle=angle,
1045            style=style,
1046            atype=atype, anchor=anchor,
1047            spin=spin, orient=orient
1048        ) children();
1049    } else {
1050        skmat = down(min_y) * skew(sxz=shift.x/h, syz=shift.y/h) * up(min_y);
1051        transforms = [
1052            if (angle==360) for (i=[0:1:steps-1]) skmat * rot([90,0,360-i*360/steps]),
1053            if (angle<360) for (i=[0:1:steps-1]) skmat * rot([90,0,angle-i*angle/(steps-1)]),
1054        ];
1055        sweep(
1056            region, transforms,
1057            closed=angle==360,
1058            caps=angle!=360,
1059            style=style, cp=cp,
1060            convexity=convexity,
1061            atype=atype, anchor=anchor,
1062            spin=spin, orient=orient
1063        ) children();
1064    }
1065}
1066
1067
1068// Function&Module: spiral_sweep()
1069// Synopsis: Sweep a path along a helix.
1070// SynTags: VNF, Geom
1071// Topics: Extrusion, Sweep, Spiral
1072// See Also: thread_helix(), linear_sweep(), rotate_sweep(), sweep(), path_sweep(), offset_sweep()
1073// Usage: As Module
1074//   spiral_sweep(poly, h, r|d=, turns, [taper=], [center=], [taper1=], [taper2=], [internal=], ...)[ATTACHMENTS];
1075//   spiral_sweep(poly, h, r1=|d1=, r2=|d2=, turns, [taper=], [center=], [taper1=], [taper2=], [internal=], ...)[ATTACHMENTS];
1076// Usage: As Function
1077//   vnf = spiral_sweep(poly, h, r|d=, turns, ...);
1078//   vnf = spiral_sweep(poly, h, r1=|d1=, r1=|d2=, turns, ...);
1079// Description:
1080//   Takes a closed 2D polygon path, centered on the XY plane, and sweeps/extrudes it along a 3D spiral path
1081//   of a given radius, height and degrees of rotation.  The origin in the profile traces out the helix of the specified radius.
1082//   If turns is positive the path will be right-handed;  if turns is negative the path will be left-handed.
1083//   Such an extrusion can be used to make screw threads.  
1084//   .
1085//   The lead_in options specify a lead-in setiton where the ends of the spiral scale down to avoid a sharp cut face at the ends.
1086//   You can specify the length of this scaling directly with the lead_in parameters or as an angle using the lead_in_ang parameters.
1087//   If you give a positive value, the extrusion is lengthenend by the specified distance or angle; if you give a negative
1088//   value then the scaled end is included in the extrusion length specified by `turns`.  If the value is zero then no scaled ends
1089//   are produced.  The shape of the scaled ends can be controlled with the lead_in_shape parameter.  Supported options are "sqrt", "linear"
1090//   "smooth" and "cut".  
1091//   .
1092//   The inside argument changes how the extrusion lead-in sections are formed.  If it is true then they scale
1093//   towards the outside, like would be needed for internal threading.  If internal is fale then the lead-in sections scale
1094//   towards the inside, like would be appropriate for external threads.  
1095// Arguments:
1096//   poly = Array of points of a polygon path, to be extruded.
1097//   h = height of the spiral extrusion path
1098//   r = Radius of the spiral extrusion path
1099//   turns = number of revolutions to include in the spiral
1100//   ---
1101//   d = Diameter of the spiral extrusion path.
1102//   d1/r1 = Bottom inside diameter or radius of spiral to extrude along.
1103//   d2/r2 = Top inside diameter or radius of spiral to extrude along.
1104//   lead_in = Specify linear length of the lead-in scaled section of the spiral.  Default: 0
1105//   lead_in1 = Specify linear length of the lead-in scaled section of the spiral at the bottom
1106//   lead_in2 = Specify linear length of the lead-in scaled section of the spiral at the top
1107//   lead_in_ang = Specify angular  length of the lead-in scaled section of the spiral
1108//   lead_in_ang1 = Specify angular length of the lead-in scaled section of the spiral at the bottom
1109//   lead_in_ang2 = Specify angular length of the lead-in scaled section of the spiral at the top
1110//   lead_in_shape = Specify the shape of the thread lead in by giving a text string or function.  Default: "sqrt"
1111//   lead_in_shape1 = Specify the shape of the thread lead-in at the bottom by giving a text string or function.  
1112//   lead_in_shape2 = Specify the shape of the thread lead-in at the top by giving a text string or function.
1113//   lead_in_sample = Factor to increase sample rate in the lead-in section.  Default: 10
1114//   internal = if true make internal threads.  The only effect this has is to change how the extrusion lead-in section are formed. When true, the extrusion scales towards the outside; when false, it scales towards the inside.  Default: false
1115//   anchor = Translate so anchor point is at origin (0,0,0).  See [anchor](attachments.scad#subsection-anchor).  Default: `CENTER`
1116//   spin = Rotate this many degrees around the Z axis after anchor.  See [spin](attachments.scad#subsection-spin).  Default: `0`
1117//   orient = Vector to rotate top towards, after spin.  See [orient](attachments.scad#subsection-orient).  Default: `UP`
1118// Example:
1119//   poly = [[-10,0], [-3,-5], [3,-5], [10,0], [0,-30]];
1120//   spiral_sweep(poly, h=200, r=50, turns=3, $fn=36);
1121_leadin_ogive=function (x,L) 
1122     let( minscale = .05,
1123          r=(L^2+(1-minscale^2))/2/(1-minscale),
1124          scale = sqrt(r^2-(L*(1-x))^2) -(r-1)
1125     )
1126     x>1 ? [1,1]
1127   : x<0 ? [lerp(minscale,1,.25),0] 
1128   : [lerp(scale,1,.25),scale];     
1129
1130_leadin_cut = function(x,L) x>0 ? [1,1] : [1,0];
1131
1132_leadin_sqrt = function(x,L)
1133     let(end=0.05)   // Smallest scale at the end
1134     x>1 ? [1,1]
1135   : x<0 ? [lerp(end,1,.25),0]   
1136   : let(  
1137          s = sqrt(x + end^2 * (1-x))
1138     )
1139     [lerp(s,1,.25),s];    // thread width scale, thread height scale
1140
1141_leadin_linear = function(x,L)
1142     let(minscale=.1)
1143     x>1 ? [1,1]
1144   : x<0 ? [lerp(minscale,1,.25),0]
1145   : let(scale = lerp(minscale,1,x))
1146     [lerp(scale,1,.25),scale];
1147
1148_lead_in_table = [
1149     ["default", _leadin_sqrt],
1150     ["sqrt", _leadin_sqrt],
1151     ["cut", _leadin_cut],
1152     ["smooth", _leadin_ogive],
1153     ["linear", _leadin_linear]
1154];
1155
1156      
1157function _ss_polygon_r(N,theta) =
1158        let( alpha = 360/N )
1159        cos(alpha/2)/(cos(posmod(theta,alpha)-alpha/2));
1160function spiral_sweep(poly, h, r, turns=1, taper, r1, r2, d, d1, d2, internal=false,
1161                      lead_in_shape,lead_in_shape1, lead_in_shape2,
1162                      lead_in, lead_in1, lead_in2,
1163                      lead_in_ang, lead_in_ang1, lead_in_ang2,
1164                      height,l,length,
1165                      lead_in_sample = 10,
1166                      anchor=CENTER, spin=0, orient=UP) =
1167    assert(is_num(turns) && turns != 0, "turns must be a nonzero number")
1168    assert(all_positive([h]), "Spiral height must be a positive number")
1169    let(
1170        dir = sign(turns),
1171        r1 = get_radius(r1=r1, r=r, d1=d1, d=d),
1172        r2 = get_radius(r1=r2, r=r, d1=d2, d=d),
1173        bounds = pointlist_bounds(poly),
1174        yctr = (bounds[0].y+bounds[1].y)/2,
1175        xmin = bounds[0].x,
1176        xmax = bounds[1].x,
1177        poly = path3d(clockwise_polygon(poly)),
1178        sides = segs(max(r1,r2)),
1179        ang_step = 360/sides,
1180        turns = abs(turns),
1181        lead_in1 = first_defined([lead_in1, lead_in]),
1182        lead_in2 = first_defined([lead_in1, lead_in]),        
1183        lead_in_ang1 =
1184                      let(
1185                           user_ang = first_defined([lead_in_ang1,lead_in_ang])
1186                      )
1187                      assert(is_undef(user_ang) || is_undef(lead_in1), "Cannot define lead_in/lead_in1 by both length and angle")
1188                      is_def(user_ang) ? user_ang : default(lead_in1,0)*360/(2*PI*r1),
1189        lead_in_ang2 =
1190                      let(
1191                           user_ang = first_defined([lead_in_ang2,lead_in_ang])
1192                      )
1193                      assert(is_undef(user_ang) || is_undef(lead_in2), "Cannot define lead_in/lead_in2 by both length and angle")
1194                      is_def(user_ang) ? user_ang : default(lead_in2,0)*360/(2*PI*r2),
1195        minang = -max(0,lead_in_ang1),
1196        maxang = 360*turns + max(0,lead_in_ang2),
1197        cut_ang1 = minang+abs(lead_in_ang1),
1198        cut_ang2 = maxang-abs(lead_in_ang1),        
1199        lead_in_shape1 = first_defined([lead_in_shape1, lead_in_shape, "default"]),
1200        lead_in_shape2 = first_defined([lead_in_shape2, lead_in_shape, "default"]),             
1201        lead_in_func1 = is_func(lead_in_shape1) ? lead_in_shape1
1202                      : assert(is_string(lead_in_shape1),"lead_in_shape/lead_in_shape1 must be a function or string")
1203                        let(ind = search([lead_in_shape1], _lead_in_table,0)[0])
1204                        assert(ind!=[],str("Unknown lead_in_shape, \"",lead_in_shape1,"\""))
1205                        _lead_in_table[ind[0]][1],
1206        lead_in_func2 = is_func(lead_in_shape2) ? lead_in_shape2
1207                      : assert(is_string(lead_in_shape2),"lead_in_shape/lead_in_shape2 must be a function or string")
1208                        let(ind = search([lead_in_shape2], _lead_in_table,0)[0])
1209                        assert(ind!=[],str("Unknown lead_in_shape, \"",lead_in_shape2,"\""))
1210                        _lead_in_table[ind[0]][1]
1211    )
1212    assert( cut_ang1<cut_ang2, "Tapers are too long to fit")
1213    assert( all_positive([r1,r2]), "Diameter/radius must be positive")
1214    let(
1215  
1216        // This complicated sampling scheme is designed to ensure that faceting always starts at angle zero
1217        // for alignment with cylinders, and there is always a facet boundary at the $fn specified locations, 
1218        // regardless of what kind of subsampling occurs for tapers.
1219        orig_anglist = [
1220            if (minang<0) minang,
1221            each reverse([for(ang = [-ang_step:-ang_step:minang+EPSILON]) ang]),
1222            for(ang = [0:ang_step:maxang-EPSILON]) ang,
1223            maxang
1224        ],
1225        anglist = [
1226           for(a=orig_anglist) if (a<cut_ang1-EPSILON) a,
1227           cut_ang1,
1228           for(a=orig_anglist) if (a>cut_ang1+EPSILON && a<cut_ang2-EPSILON) a,
1229           cut_ang2,
1230           for(a=orig_anglist) if (a>cut_ang2+EPSILON) a
1231        ],
1232        interp_ang = [
1233                      for(i=idx(anglist,e=-2)) 
1234                          each lerpn(anglist[i],anglist[i+1],
1235                                         (lead_in_ang1!=0 && anglist[i+1]<=cut_ang1) || (lead_in_ang2!=0 && anglist[i]>=cut_ang2)
1236                                            ? ceil((anglist[i+1]-anglist[i])/ang_step*lead_in_sample)
1237                                            : 1,
1238                                     endpoint=false),
1239                      last(anglist)
1240                     ],
1241        skewmat = affine3d_skew_xz(xa=atan2(r2-r1,h)),
1242        points = [
1243            for (a = interp_ang) let (
1244                hsc = a<cut_ang1 ? lead_in_func1((a-minang)/abs(lead_in_ang1),abs(lead_in_ang1)*2*PI*r1/360)
1245                    : a>cut_ang2 ? lead_in_func2((maxang-a)/abs(lead_in_ang2),abs(lead_in_ang2)*2*PI*r2/360)
1246                    : [1,1],
1247                u = a/(360*turns), 
1248                r = lerp(r1,r2,u),
1249                mat = affine3d_zrot(dir*a)
1250                    * affine3d_translate([_ss_polygon_r(sides,dir*a)*r, 0, h * (u-0.5)])
1251                    * affine3d_xrot(90)
1252                    * skewmat
1253                    * scale([hsc.y,hsc.x,1], cp=[internal ? xmax : xmin, yctr, 0]),
1254                pts = apply(mat, poly)
1255            ) pts
1256        ],
1257        vnf = vnf_vertex_array(
1258            points, col_wrap=true, caps=true, reverse=dir>0,
1259        //    style=higbee1>0 || higbee2>0 ? "quincunx" : "alt"
1260            style="convex"
1261        )
1262    )
1263    reorient(anchor,spin,orient, vnf=vnf, r1=r1, r2=r2, l=h, p=vnf);
1264
1265
1266
1267module spiral_sweep(poly, h, r, turns=1, taper, r1, r2, d, d1, d2, internal=false,
1268                    lead_in_shape,lead_in_shape1, lead_in_shape2,
1269                    lead_in, lead_in1, lead_in2,
1270                    lead_in_ang, lead_in_ang1, lead_in_ang2,
1271                    height,l,length,
1272                    lead_in_sample=10,
1273                    anchor=CENTER, spin=0, orient=UP)
1274{
1275    vnf = spiral_sweep(poly=poly, h=h, r=r, turns=turns, r1=r1, r2=r2, d=d, d1=d1, d2=d2, internal=internal,
1276                       lead_in_shape=lead_in_shape,lead_in_shape1=lead_in_shape1, lead_in_shape2=lead_in_shape2,
1277                       lead_in=lead_in, lead_in1=lead_in1, lead_in2=lead_in2,
1278                       lead_in_ang=lead_in_ang, lead_in_ang1=lead_in_ang1, lead_in_ang2=lead_in_ang2,
1279                       height=height,l=length,length=length,
1280                       lead_in_sample=lead_in_sample);
1281    h = one_defined([h,height,length,l],"h,height,length,l");
1282    r1 = get_radius(r1=r1, r=r, d1=d1, d=d);
1283    r2 = get_radius(r1=r2, r=r, d1=d2, d=d);
1284    lead_in1 = u_mul(first_defined([lead_in1,lead_in]),1/(2*PI*r1));
1285    lead_in2 = u_mul(first_defined([lead_in2,lead_in]),1/(2*PI*r2));
1286    lead_in_ang1 = first_defined([lead_in_ang1,lead_in_ang]);
1287    lead_in_ang2 = first_defined([lead_in_ang2,lead_in_ang]);
1288    extra_turns = max(0,first_defined([lead_in1,lead_in_ang1,0]))+max(0,first_defined([lead_in2,lead_in_ang2,0]));
1289    attachable(anchor,spin,orient, r1=r1, r2=r2, l=h) {
1290        vnf_polyhedron(vnf, convexity=ceil(2*(abs(turns)+extra_turns)));
1291        children();
1292    }
1293}
1294
1295
1296
1297// Function&Module: path_sweep()
1298// Synopsis: Sweep a 2d polygon path along a 2d or 3d path. 
1299// SynTags: VNF, Geom
1300// Topics: Extrusion, Sweep, Paths
1301// See Also: linear_sweep(), rotate_sweep(), sweep(), spiral_sweep(), path_sweep2d(), offset_sweep()
1302// Usage: As module
1303//   path_sweep(shape, path, [method], [normal=], [closed=], [twist=], [twist_by_length=], [symmetry=], [scale=], [scale_by_length=], [last_normal=], [tangent=], [uniform=], [relaxed=], [caps=], [style=], [convexity=], [anchor=], [cp=], [spin=], [orient=], [atype=]) [ATTACHMENTS];
1304// Usage: As function
1305//   vnf = path_sweep(shape, path, [method], [normal=], [closed=], [twist=], [twist_by_length=], [symmetry=], [scale=], [scale_by_length=], [last_normal=], [tangent=], [uniform=], [relaxed=], [caps=], [style=], [transforms=], [anchor=], [cp=], [spin=], [orient=], [atype=]);
1306// Description:
1307//   Takes as input `shape`, a 2D polygon path (list of points), and `path`, a 2d or 3d path (also a list of points)
1308//   and constructs a polyhedron by sweeping the shape along the path. When run as a module returns the polyhedron geometry.
1309//   When run as a function returns a VNF by default or if you set `transforms=true` then it returns a list of transformations suitable as input to `sweep`.
1310//   .
1311//   The sweeping process places one copy of the shape for each point in the path.  The origin in `shape` is translated to
1312//   the point in `path`.  The normal vector of the shape, which points in the Z direction, is aligned with the tangent
1313//   vector for the path, so this process is constructing a shape whose normal cross sections are equal to your specified shape.
1314//   If you do not supply a list of tangent vectors then an approximate tangent vector is computed
1315//   based on the path points you supply using {{path_tangents()}}.
1316// Figure(3D,Big,VPR=[70,0,345],VPD=20,VPT=[5.5,10.8,-2.7],NoScales): This example shows how the shape, in this case the quadrilateral defined by `[[0, 0], [0, 1], [0.25, 1], [1, 0]]`, appears as the cross section of the swept polyhedron.  The blue line shows the path.  The normal vector to the shape is shown in black; it is based at the origin and points upwards in the Z direction.  The sweep aligns this normal vector with the blue path tangent, which in this case, flips the shape around.  Note that for a 2D path like this one, the Y direction in the shape is mapped to the Z direction in the sweep.
1317//   tri= [[0, 0], [0, 1], [.25,1], [1, 0]];
1318//   path = arc(r=5,n=81,angle=[-20,65]);
1319//   % path_sweep(tri,path);
1320//   T = path_sweep(tri,path,transforms=true);
1321//   color("red")for(i=[0:20:80]) stroke(apply(T[i],path3d(tri)),width=.1,closed=true);
1322//   color("blue")stroke(path3d(arc(r=5,n=101,angle=[-20,80])),width=.1,endcap2="arrow2");
1323//   color("red")stroke([path3d(tri)],width=.1);
1324//   stroke([CENTER,UP], width=.07,endcap2="arrow2",color="black");
1325// Continues:
1326//   In the figure you can see that the swept polyhedron, shown in transparent gray, has the quadrilateral as its cross
1327//   section.  The quadrilateral is positioned perpendicular to the path, which is shown in blue, so that the normal
1328//   vector for the quadrilateral is parallel to the tangent vector for the path.  The origin for the shape is the point
1329//   which follows the path.  For a 2D path, the Y axis of the shape is mapped to the Z axis and in this case,
1330//   pointing the quadrilateral's normal vector (in black) along the tangent line of
1331//   the path, which is going in the direction of the blue arrow, requires that the quadrilateral be "turned around".  If we
1332//   reverse the order of points in the path we get a different result:
1333// Figure(3D,Big,VPR=[70,0,20],VPD=20,VPT=[1.25,9.25,-2.65],NoScales): The same sweep operation with the path traveling in the opposite direction.  Note that in order to line up the normal correctly, the shape is reversed compared to Figure 1, so the resulting sweep looks quite different.
1334//   tri= [[0, 0], [0, 1], [.25,1], [1, 0]];
1335//   path = reverse(arc(r=5,n=81,angle=[-20,65]));
1336//   % path_sweep(tri,path);
1337//   T = path_sweep(tri,path,transforms=true);
1338//   color("red")for(i=[0:20:80]) stroke(apply(T[i],path3d(tri)),width=.1,closed=true);
1339//   color("blue")stroke(reverse(path3d(arc(r=5,n=101,angle=[-20-15,65]))),width=.1,endcap2="arrow2");
1340//   color("red")stroke([path3d(tri)],width=.1);
1341//   stroke([CENTER,UP], width=.07,endcap2="arrow2",color="black");
1342// Continues:
1343//   If your shape is too large for the curves in the path you can create a situation where the shapes cross each
1344//   other.  This results in an invalid polyhedron, which may appear OK when previewed or rendered alone, but will give rise
1345//   to cryptic CGAL errors when rendered with a second object in your model.  You may be able to use {{path_sweep2d()}}
1346//   to produce a valid model in cases like this.  You can debug models like this using the `profiles=true` option which will show all
1347//   the cross sections in your polyhedron.  If any of them intersect, the polyhedron will be invalid.
1348// Figure(3D,Big,VPR=[47,0,325],VPD=23,VPT=[6.8,4,-3.8],NoScales): We have scaled the path to an ellipse and show a large triangle as the shape.  The triangle is sometimes bigger than the local radius of the path, leading to an invalid polyhedron, which you can identify because the red lines cross in the middle.
1349//   tri= scale([4.5,2.5],[[0, 0], [0, 1], [1, 0]]);
1350//   path = xscale(1.5,arc(r=5,n=81,angle=[-70,70]));
1351//   % path_sweep(tri,path);
1352//   T = path_sweep(tri,path,transforms=true);
1353//   color("red")for(i=[0:20:80]) stroke(apply(T[i],path3d(tri)),width=.1,closed=true);
1354//   color("blue")stroke(path3d(xscale(1.5,arc(r=5,n=81,angle=[-70,80]))),width=.1,endcap2="arrow2");
1355// Continues:
1356//   During the sweep operation the shape's normal vector aligns with the tangent vector of the path.  Note that
1357//   this leaves an ambiguity about how the shape is rotated as it sweeps along the path.
1358//   For 2D paths, this ambiguity is resolved by aligning the Y axis of the shape to the Z axis of the swept polyhedron.
1359//   You can force the  shape to twist as it sweeps along the path using the `twist` parameter, which specifies the total
1360//   number of degrees to twist along the whole swept polyhedron.  This produces a result like the one shown below.
1361// Figure(3D,Big,VPR=[66,0,14],VPD=20,VPT=[3.4,4.5,-0.8]): The shape twists as we sweep.  Note that it still aligns the origin in the shape with the path, and still aligns the normal vector with the path tangent vector.
1362//   tri= [[0, 0], [0, 1], [.25,1],[1, 0]];
1363//   path = arc(r=5,n=81,angle=[-20,65]);
1364//   % path_sweep(tri,path,twist=-60);
1365//   T = path_sweep(tri,path,transforms=true,twist=-60);
1366//   color("red")for(i=[0:20:80]) stroke(apply(T[i],path3d(tri)),width=.1,closed=true);
1367//   color("blue")stroke(path3d(arc(r=5,n=101,angle=[-20,80])),width=.1,endcap2="arrow2");
1368// Continues:
1369//   The `twist` argument adds the specified number of degrees of twist into the model, and it may be positive or
1370//   negative.  When `closed=true` the starting shape and ending shape must match to avoid a sudden extreme twist at the
1371//   joint.  By default `twist` is therefore required to be a multiple of 360.  However, if your shape has rotational
1372//   symmetry, this requirement is overly strict.  You can specify the symmetry using the `symmetry` argument, and then
1373//   you can choose smaller twists consistent with the specified symmetry.  The symmetry argument gives the number of
1374//   rotations that map the shape exactly onto itself, so a pentagon has 5-fold symmetry.  This argument is only valid
1375//   for closed sweeps.  When you specify symmetry, the twist must be a multiple of 360/symmetry.
1376//   .
1377//   The twist is normally spread uniformly along your shape based on the path length.  If you set `twist_by_length` to
1378//   false then the twist will be uniform based on the point count of your path.  Twisted shapes will produce twisted
1379//   faces, so if you want them to look good you should use lots of points on your path and also lots of points on the
1380//   shape.  If your shape is a simple polygon, use {{subdivide_path()}} to increase
1381//   the number of points.
1382//   .
1383//   As noted above, the sweep process has an ambiguity regarding the twist.  For 2D paths it is easy to resolve this
1384//   ambiguity by aligning the Y axis in the shape to the Z axis in the swept polyhedron.  When the path is
1385//   three-dimensional, things become more complex.  It is no longer possible to use a simple alignment rule like the
1386//   one we use in 2D.  You may find that the shape rotates unexpectedly around its axis as it traverses the path.  The
1387//   `method` parameter allows you to specify how the shapes are aligned, resulting in different twist in the resulting
1388//   polyhedron.  You can choose from three different methods for selecting the rotation of your shape.  None of these
1389//   methods will produce good, or even valid, results on all inputs, so it is important to select a suitable method.
1390//   .
1391//   The three methods you can choose using the `method` parameter are:
1392//   .
1393//   The "incremental" method (the default) works by adjusting the shape at each step by the minimal rotation that makes the shape normal to the tangent
1394//   at the next point.  This method is robust in that it always produces a valid result for well-behaved paths with sufficiently high
1395//   sampling.  Unfortunately, it can produce a large amount of undesirable twist.  When constructing a closed shape this algorithm in
1396//   its basic form provides no guarantee that the start and end shapes match up.  To prevent a sudden twist at the last segment,
1397//   the method calculates the required twist for a good match and distributes it over the whole model (as if you had specified a
1398//   twist amount).  If you specify `symmetry` this may allow the algorithm to choose a smaller twist for this alignment.
1399//   To start the algorithm, we need an initial condition.  This is supplied by
1400//   using the `normal` argument to give a direction to align the Y axis of your shape.  By default the normal points UP if the path
1401//   makes an angle of 45 deg or less with the xy plane and it points BACK if the path makes a higher angle with the XY plane.  You
1402//   can also supply `last_normal` which provides an ending orientation constraint.  Be aware that the curve may still exhibit
1403//   twisting in the middle.  This method is the default because it is the most robust, not because it generally produces the best result.
1404//   .
1405//   The "natural" method works by computing the Frenet frame at each point on the path.  This is defined by the tangent to the curve and
1406//   the normal which lies in the plane defined by the curve at each point.  This normal points in the direction of curvature of the curve.
1407//   The result is a very well behaved set of shape positions without any unexpected twisting&mdash;as long as the curvature never falls to zero.  At a
1408//   point of zero curvature (a flat point), the curve does not define a plane and the natural normal is not defined.  Furthermore, even if
1409//   you skip over this troublesome point so the normal is defined, it can change direction abruptly when the curvature is zero, leading to
1410//   a nasty twist and an invalid model.  A simple example is a circular arc joined to another arc that curves the other direction.  Note
1411//   that the X axis of the shape is aligned with the normal from the Frenet frame.
1412//   .
1413//   The "manual" method allows you to specify your desired normal either globally with a single vector, or locally with
1414//   a list of normal vectors for every path point.  The normal you supply is projected to be orthogonal to the tangent to the
1415//   path and the Y direction of your shape will be aligned with the projected normal.  (Note this is different from the "natural" method.)
1416//   Careless choice of a normal may result in a twist in the shape, or an error if your normal is parallel to the path tangent.
1417//   If you set `relax=true` then the condition that the cross sections are orthogonal to the path is relaxed and the swept object
1418//   uses the actual specified normal.  In this case, the tangent is projected to be orthogonal to your supplied normal to define
1419//   the cross section orientation.  Specifying a list of normal vectors gives you complete control over the orientation of your
1420//   cross sections and can be useful if you want to position your model to be on the surface of some solid.
1421//   .
1422//   You can also apply scaling to the profile along the path.  You can give a list of scalar scale factors or a list of 2-vector scale. 
1423//   In the latter scale the x and y scales of the profile are scaled separately before the profile is placed onto the path.  For non-closed
1424//   paths you can also give a single scale value or a 2-vector which is treated as the final scale.  The intermediate sections
1425//   are then scaled by linear interpolation either relative to length (if scale_by_length is true) or by point count otherwise.  
1426//   .
1427//   You can use set `transforms` to true to return a list of transformation matrices instead of the swept shape.  In this case, you can
1428//   often omit shape entirely.  The exception is when `closed=true` and you are using the "incremental" method.  In this case, `path_sweep`
1429//   uses the shape to correct for twist when the shape closes on itself, so you must include a valid shape.
1430// Arguments:
1431//   shape = A 2D polygon path or region describing the shape to be swept.
1432//   path = 2D or 3D path giving the path to sweep over
1433//   method = one of "incremental", "natural" or "manual".  Default: "incremental"
1434//   ---
1435//   normal = normal vector for initializing the incremental method, or for setting normals with method="manual".  Default: UP if the path makes an angle lower than 45 degrees to the xy plane, BACK otherwise.
1436//   closed = path is a closed loop.  Default: false
1437//   twist = amount of twist to add in degrees.  For closed sweeps must be a multiple of 360/symmetry.  Default: 0
1438//   twist_by_length = if true then interpolate twist based on the path length of the path. If false interoplate based on point count.  Default: true
1439//   symmetry = symmetry of the shape when closed=true.  Allows the shape to join with a 360/symmetry rotation instead of a full 360 rotation.  Default: 1
1440//   scale = Amount to scale the profiles.  If you give a scalar the scale starts at 1 and ends at your specified value. The same is true for a 2-vector, but x and y are scaled separately.   You can also give a vector of values, one for each path point, and you can give a list of 2-vectors that give the x and y scales of your profile for every point on the path (a Nx2 matrix for a path of length N.  Default: 1 (no scaling)
1441//   scale_by_length = if true then interpolate scale based on the path length of the path. If false interoplate based on point count.  Default: true
1442//   last_normal = normal to last point in the path for the "incremental" method.  Constrains the orientation of the last cross section if you supply it.
1443//   uniform = if set to false then compute tangents using the uniform=false argument, which may give better results when your path is non-uniformly sampled.  This argument is passed to {{path_tangents()}}.  Default: true
1444//   tangent = a list of tangent vectors in case you need more accuracy (particularly at the end points of your curve)
1445//   relaxed = set to true with the "manual" method to relax the orthogonality requirement of cross sections to the path tangent.  Default: false
1446//   caps = Can be a boolean or vector of two booleans.  Set to false to disable caps at the two ends.  Default: true
1447//   style = vnf_vertex_array style.  Default: "min_edge"
1448//   profiles = if true then display all the cross section profiles instead of the solid shape.  Can help debug a sweep.  (module only) Default: false
1449//   width = the width of lines used for profile display.  (module only) Default: 1
1450//   transforms = set to true to return transforms instead of a VNF.  These transforms can be manipulated and passed to sweep().  (function only)  Default: false.
1451//   convexity = convexity parameter for polyhedron().  (module only)  Default: 10
1452//   anchor = Translate so anchor point is at the origin. Default: "origin"
1453//   spin = Rotate this many degrees around Z axis after anchor. Default: 0
1454//   orient = Vector to rotate top towards after spin
1455//   atype  = Select "hull" or "intersect" anchor types.  Default: "hull"
1456//   cp = Centerpoint for determining "intersect" anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: "centroid"
1457// Anchor Types:
1458//   "hull" = Anchors to the virtual convex hull of the shape.
1459//   "intersect" = Anchors to the surface of the shape.
1460// Example(NoScales): A simple sweep of a square along a sine wave:
1461//   path = [for(theta=[-180:5:180]) [theta/10, 10*sin(theta)]];
1462//   sq = square(6,center=true);
1463//   path_sweep(sq,path);
1464// Example(NoScales): If the square is not centered, then we get a different result because the shape is in a different place relative to the origin:
1465//   path = [for(theta=[-180:5:180]) [theta/10, 10*sin(theta)]];
1466//   sq = square(6);
1467//   path_sweep(sq,path);
1468// Example(Med,VPR=[34,0,8],NoScales): It may not be obvious, but the polyhedron in the previous example is invalid.  It will eventually give CGAL errors when you combine it with other shapes.  To see this, set profiles to true and look at the left side.  The profiles cross each other and intersect.  Any time this happens, your polyhedron is invalid, even if it seems to be working at first.  Another observation from the profile display is that we have more profiles than needed over a lot of the shape, so if the model is slow, using fewer profiles in the flat portion of the curve might speed up the calculation.
1469//   path = [for(theta=[-180:5:180]) [theta/10, 10*sin(theta)]];
1470//   sq = square(6);
1471//   path_sweep(sq,path,profiles=true,width=.1,$fn=8);
1472// Example(2D): We'll use this shape in several examples
1473//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1474//   polygon(ushape);
1475// Example(NoScales): Sweep along a clockwise elliptical arc, using default "incremental" method.
1476//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1477//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[180,00], r=30));  // Clockwise
1478//   path_sweep(ushape, path3d(elliptic_arc));
1479// Example(NoScales): Sweep along a counter-clockwise elliptical arc.  Note that the orientation of the shape flips.
1480//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1481//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[0,180], r=30));   // Counter-clockwise
1482//   path_sweep(ushape, path3d(elliptic_arc));
1483// Example(NoScales): Sweep along a clockwise elliptical arc, using "natural" method, which lines up the X axis of the shape with the direction of curvature.  This means the X axis will point inward, so a counterclockwise arc gives:
1484//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1485//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[0,180], r=30));  // Counter-clockwise
1486//   path_sweep(ushape, elliptic_arc, method="natural");
1487// Example(NoScales): Sweep along a clockwise elliptical arc, using "natural" method.  If the curve is clockwise then the shape flips upside-down to align the X axis.
1488//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1489//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[180,0], r=30));  // Clockwise
1490//   path_sweep(ushape, path3d(elliptic_arc), method="natural");
1491// Example(NoScales): Sweep along a clockwise elliptical arc, using "manual" method.  You can orient the shape in a direction you choose (subject to the constraint that the profiles remain normal to the path):
1492//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1493//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[180,0], r=30));  // Clockwise
1494//   path_sweep(ushape, path3d(elliptic_arc), method="manual", normal=UP+RIGHT);
1495// Example(NoScales): Here we changed the ellipse to be more pointy, and with the same results as above we get a shape with an irregularity in the middle where it maintains the specified direction around the point of the ellipse.  If the ellipse were more pointy, this would result in a bad polyhedron:
1496//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1497//   elliptic_arc = yscale(2, p=arc($fn=64,angle=[180,0], r=30));  // Clockwise
1498//   path_sweep(ushape, path3d(elliptic_arc), method="manual", normal=UP+RIGHT);
1499// Example(NoScales): It is easy to produce an invalid shape when your path has a smaller radius of curvature than the width of your shape.  The exact threshold where the shape becomes invalid depends on the density of points on your path.  The error may not be immediately obvious, as the swept shape appears fine when alone in your model, but adding a cube to the model reveals the problem.  In this case the pentagon is turned so its longest direction points inward to create the singularity.
1500//   qpath = [for(x=[-3:.01:3]) [x,x*x/1.8,0]];
1501//   // Prints 0.9, but we use pentagon with radius of 1.0 > 0.9
1502//   echo(radius_of_curvature = 1/max(path_curvature(qpath)));
1503//   path_sweep(apply(rot(90),pentagon(r=1)), qpath, normal=BACK, method="manual");
1504//   cube(0.5);    // Adding a small cube forces a CGAL computation which reveals
1505//                 // the error by displaying nothing or giving a cryptic message
1506// Example(NoScales): Using the `relax` option we allow the profiles to deviate from orthogonality to the path.  This eliminates the crease that broke the previous example because the sections are all parallel to each other.
1507//   qpath = [for(x=[-3:.01:3]) [x,x*x/1.8,0]];
1508//   path_sweep(apply(rot(90),pentagon(r=1)), qpath, normal=BACK, method="manual", relaxed=true);
1509//   cube(0.5);    // Adding a small cube is not a problem with this valid model
1510// Example(Med,VPR=[16,0,100],VPT=[0.05,0.6,0.6],VPD=25,NoScales): Using the `profiles=true` option can help debug bad polyhedra such as this one.  If any of the profiles intersect or cross each other, the polyhedron will be invalid.  In this case, you can see these intersections in the middle of the shape, which may give insight into how to fix your shape.   The profiles may also help you identify cases with a valid polyhedron where you have more profiles than needed to adequately define the shape.
1511//   tri= scale([4.5,2.5],[[0, 0], [0, 1], [1, 0]]);
1512//   path = left(4,xscale(1.5,arc(r=5,n=25,angle=[-70,70])));
1513//   path_sweep(tri,path,profiles=true,width=.1);
1514// Example(NoScales):  This 3d arc produces a result that twists to an undefined angle.  By default the incremental method sets the starting normal to UP, but the ending normal is unconstrained.
1515//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1516//   arc = yrot(37, p=path3d(arc($fn=64, r=30, angle=[0,180])));
1517//   path_sweep(ushape, arc, method="incremental");
1518// Example(NoScales): You can constrain the last normal as well.  Here we point it right, which produces a nice result.
1519//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1520//   arc = yrot(37, p=path3d(arc($fn=64, r=30, angle=[0,180])));
1521//   path_sweep(ushape, arc, method="incremental", last_normal=RIGHT);
1522// Example(NoScales): Here we constrain the last normal to UP.  Be aware that the behavior in the middle is unconstrained.
1523//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1524//   arc = yrot(37, p=path3d(arc($fn=64, r=30, angle=[0,180])));
1525//   path_sweep(ushape, arc, method="incremental", last_normal=UP);
1526// Example(NoScales): The "natural" method produces a very different result
1527//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1528//   arc = yrot(37, p=path3d(arc($fn=64, r=30, angle=[0,180])));
1529//   path_sweep(ushape, arc, method="natural");
1530// Example(NoScales): When the path starts at an angle of more that 45 deg to the xy plane the initial normal for "incremental" is BACK.  This produces the effect of the shape rising up out of the xy plane.  (Using UP for a vertical path is invalid, hence the need for a split in the defaults.)
1531//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1532//   arc = xrot(75, p=path3d(arc($fn=64, r=30, angle=[0,180])));
1533//   path_sweep(ushape, arc, method="incremental");
1534// Example(NoScales): Adding twist
1535//   // Counter-clockwise
1536//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[0,180], r=3));
1537//   path_sweep(pentagon(r=1), path3d(elliptic_arc), twist=72);
1538// Example(NoScales): Closed shape
1539//   ellipse = xscale(2, p=circle($fn=64, r=3));
1540//   path_sweep(pentagon(r=1), path3d(ellipse), closed=true);
1541// Example(NoScales): Closed shape with added twist
1542//   ellipse = xscale(2, p=circle($fn=64, r=3));
1543//   // Looks better with finer sampling
1544//   pentagon = subdivide_path(pentagon(r=1), 30);
1545//   path_sweep(pentagon, path3d(ellipse),
1546//              closed=true, twist=360);
1547// Example(NoScales): The last example was a lot of twist.  In order to use less twist you have to tell `path_sweep` that your shape has symmetry, in this case 5-fold.  Mobius strip with pentagon cross section:
1548//   ellipse = xscale(2, p=circle($fn=64, r=3));
1549//   // Looks better with finer sampling
1550//   pentagon = subdivide_path(pentagon(r=1), 30);
1551//   path_sweep(pentagon, path3d(ellipse), closed=true,
1552//              symmetry = 5, twist=2*360/5);
1553// Example(Med,NoScales): A helical path reveals the big problem with the "incremental" method: it can introduce unexpected and extreme twisting.  (Note helix example came from list-comprehension-demos)
1554//   function helix(t) = [(t / 1.5 + 0.5) * 30 * cos(6 * 360 * t),
1555//                        (t / 1.5 + 0.5) * 30 * sin(6 * 360 * t),
1556//                         200 * (1 - t)];
1557//   helix_steps = 200;
1558//   helix = [for (i=[0:helix_steps]) helix(i/helix_steps)];
1559//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1560//   path_sweep(ushape, helix);
1561// Example(Med,NoScales): You can constrain both ends, but still the twist remains:
1562//   function helix(t) = [(t / 1.5 + 0.5) * 30 * cos(6 * 360 * t),
1563//                        (t / 1.5 + 0.5) * 30 * sin(6 * 360 * t),
1564//                         200 * (1 - t)];
1565//   helix_steps = 200;
1566//   helix = [for (i=[0:helix_steps]) helix(i/helix_steps)];
1567//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1568//   path_sweep(ushape, helix, normal=UP, last_normal=UP);
1569// Example(Med,NoScales): Even if you manually guess the amount of twist and remove it, the result twists one way and then the other:
1570//   function helix(t) = [(t / 1.5 + 0.5) * 30 * cos(6 * 360 * t),
1571//                        (t / 1.5 + 0.5) * 30 * sin(6 * 360 * t),
1572//                         200 * (1 - t)];
1573//   helix_steps = 200;
1574//   helix = [for (i=[0:helix_steps]) helix(i/helix_steps)];
1575//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1576//   path_sweep(ushape, helix, normal=UP, last_normal=UP, twist=360);
1577// Example(Med,NoScales): To get a good result you must use a different method.
1578//   function helix(t) = [(t / 1.5 + 0.5) * 30 * cos(6 * 360 * t),
1579//                        (t / 1.5 + 0.5) * 30 * sin(6 * 360 * t),
1580//                         200 * (1 - t)];
1581//   helix_steps = 200;
1582//   helix = [for (i=[0:helix_steps]) helix(i/helix_steps)];
1583//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1584//   path_sweep(ushape, helix, method="natural");
1585// Example(Med,NoScales): Note that it may look like the shape above is flat, but the profiles are very slightly tilted due to the nonzero torsion of the curve.  If you want as flat as possible, specify it so with the "manual" method:
1586//   function helix(t) = [(t / 1.5 + 0.5) * 30 * cos(6 * 360 * t),
1587//                        (t / 1.5 + 0.5) * 30 * sin(6 * 360 * t),
1588//                         200 * (1 - t)];
1589//   helix_steps = 200;
1590//   helix = [for (i=[0:helix_steps]) helix(i/helix_steps)];
1591//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1592//   path_sweep(ushape, helix, method="manual", normal=UP);
1593// Example(Med,NoScales): What if you want to angle the shape inward?  This requires a different normal at every point in the path:
1594//   function helix(t) = [(t / 1.5 + 0.5) * 30 * cos(6 * 360 * t),
1595//                        (t / 1.5 + 0.5) * 30 * sin(6 * 360 * t),
1596//                         200 * (1 - t)];
1597//   helix_steps = 200;
1598//   helix = [for (i=[0:helix_steps]) helix(i/helix_steps)];
1599//   normals = [for(i=[0:helix_steps]) [-cos(6*360*i/helix_steps), -sin(6*360*i/helix_steps), 2.5]];
1600//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1601//   path_sweep(ushape, helix, method="manual", normal=normals);
1602// Example(NoScales): When using "manual" it is important to choose a normal that works for the whole path, producing a consistent result.  Here we have specified an upward normal, and indeed the shape is pointed up everywhere, but two abrupt transitional twists render the model invalid.
1603//   yzcircle = yrot(90,p=path3d(circle($fn=64, r=30)));
1604//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1605//   path_sweep(ushape, yzcircle, method="manual", normal=UP, closed=true);
1606// Example(NoScales): The "natural" method will introduce twists when the curvature changes direction.  A warning is displayed.
1607//   arc1 = path3d(arc(angle=90, r=30));
1608//   arc2 = xrot(-90, cp=[0,30],p=path3d(arc(angle=[90,180], r=30)));
1609//   two_arcs = path_merge_collinear(concat(arc1,arc2));
1610//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1611//   path_sweep(ushape, two_arcs, method="natural");
1612// Example(NoScales): The only simple way to get a good result is the "incremental" method:
1613//   arc1 = path3d(arc(angle=90, r=30));
1614//   arc2 = xrot(-90, cp=[0,30],p=path3d(arc(angle=[90,180], r=30)));
1615//   arc3 = apply( translate([-30,60,30])*yrot(90), path3d(arc(angle=[270,180], r=30)));
1616//   three_arcs = path_merge_collinear(concat(arc1,arc2,arc3));
1617//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1618//   path_sweep(ushape, three_arcs, method="incremental");
1619// Example(Med,NoScales): knot example from list-comprehension-demos, "incremental" method
1620//   function knot(a,b,t) =   // rolling knot
1621//        [ a * cos (3 * t) / (1 - b* sin (2 *t)),
1622//          a * sin( 3 * t) / (1 - b* sin (2 *t)),
1623//        1.8 * b * cos (2 * t) /(1 - b* sin (2 *t))];
1624//   a = 0.8; b = sqrt (1 - a * a);
1625//   ksteps = 400;
1626//   knot_path = [for (i=[0:ksteps-1]) 50 * knot(a,b,(i/ksteps)*360)];
1627//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1628//   path_sweep(ushape, knot_path, closed=true, method="incremental");
1629// Example(Med,NoScales): knot example from list-comprehension-demos, "natural" method.  Which one do you like better?
1630//   function knot(a,b,t) =   // rolling knot
1631//        [ a * cos (3 * t) / (1 - b* sin (2 *t)),
1632//          a * sin( 3 * t) / (1 - b* sin (2 *t)),
1633//        1.8 * b * cos (2 * t) /(1 - b* sin (2 *t))];
1634//   a = 0.8; b = sqrt (1 - a * a);
1635//   ksteps = 400;
1636//   knot_path = [for (i=[0:ksteps-1]) 50 * knot(a,b,(i/ksteps)*360)];
1637//   ushape = [[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1638//   path_sweep(ushape, knot_path, closed=true, method="natural");
1639// Example(Med,NoScales): knot with twist.  Note if you twist it the other direction the center section untwists because of the natural twist there.  Also compare to the "incremental" method which has less twist in the center.
1640//   function knot(a,b,t) =   // rolling knot
1641//        [ a * cos (3 * t) / (1 - b* sin (2 *t)),
1642//          a * sin( 3 * t) / (1 - b* sin (2 *t)),
1643//        1.8 * b * cos (2 * t) /(1 - b* sin (2 *t))];
1644//   a = 0.8; b = sqrt (1 - a * a);
1645//   ksteps = 400;
1646//   knot_path = [for (i=[0:ksteps-1]) 50 * knot(a,b,(i/ksteps)*360)];
1647//   path_sweep(subdivide_path(pentagon(r=12),30), knot_path, closed=true,
1648//              twist=-360*8, symmetry=5, method="natural");
1649// Example(Med,NoScales): twisted knot with twist distributed by path sample points instead of by length using `twist_by_length=false`
1650//   function knot(a,b,t) =   // rolling knot
1651//           [ a * cos (3 * t) / (1 - b* sin (2 *t)),
1652//             a * sin( 3 * t) / (1 - b* sin (2 *t)),
1653//           1.8 * b * cos (2 * t) /(1 - b* sin (2 *t))];
1654//   a = 0.8; b = sqrt (1 - a * a);
1655//   ksteps = 400;
1656//   knot_path = [for (i=[0:ksteps-1]) 50 * knot(a,b,(i/ksteps)*360)];
1657//   path_sweep(subdivide_path(pentagon(r=12),30), knot_path, closed=true,
1658//              twist=-360*8, symmetry=5, method="natural", twist_by_length=false);
1659// Example(Big,NoScales): This torus knot example comes from list-comprehension-demos.  The knot lies on the surface of a torus.  When we use the "natural" method the swept figure is angled compared to the surface of the torus because the curve doesn't follow geodesics of the torus.
1660//   function knot(phi,R,r,p,q) =
1661//       [ (r * cos(q * phi) + R) * cos(p * phi),
1662//         (r * cos(q * phi) + R) * sin(p * phi),
1663//          r * sin(q * phi) ];
1664//   ushape = 3*[[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1665//   points = 50;       // points per loop
1666//   R = 400; r = 150;  // Torus size
1667//   p = 2;  q = 5;     // Knot parameters
1668//   %torus(r_maj=R,r_min=r);
1669//   k = max(p,q) / gcd(p,q) * points;
1670//   knot_path   = [ for (i=[0:k-1]) knot(360*i/k/gcd(p,q),R,r,p,q) ];
1671//   path_sweep(rot(90,p=ushape),knot_path,  method="natural", closed=true);
1672// Example(Big,NoScales): By computing the normal to the torus at the path we can orient the path to lie on the surface of the torus:
1673//   function knot(phi,R,r,p,q) =
1674//       [ (r * cos(q * phi) + R) * cos(p * phi),
1675//         (r * cos(q * phi) + R) * sin(p * phi),
1676//          r * sin(q * phi) ];
1677//   function knot_normal(phi,R,r,p,q) =
1678//       knot(phi,R,r,p,q)
1679//           - R*unit(knot(phi,R,r,p,q)
1680//               - [0,0, knot(phi,R,r,p,q)[2]]) ;
1681//   ushape = 3*[[-10, 0],[-10, 10],[ -7, 10],[ -7, 2],[  7, 2],[  7, 7],[ 10, 7],[ 10, 0]];
1682//   points = 50;       // points per loop
1683//   R = 400; r = 150;  // Torus size
1684//   p = 2;  q = 5;     // Knot parameters
1685//   %torus(r_maj=R,r_min=r);
1686//   k = max(p,q) / gcd(p,q) * points;
1687//   knot_path   = [ for (i=[0:k-1]) knot(360*i/k/gcd(p,q),R,r,p,q) ];
1688//   normals = [ for (i=[0:k-1]) knot_normal(360*i/k/gcd(p,q),R,r,p,q) ];
1689//   path_sweep(ushape,knot_path,normal=normals, method="manual", closed=true);
1690// Example(NoScales): You can request the transformations and manipulate them before passing them on to sweep.  Here we construct a tube that changes scale by first generating the transforms and then applying the scale factor and connecting the inside and outside.  Note that the wall thickness varies because it is produced by scaling.
1691//   shape = star(n=5, r=10, ir=5);
1692//   rpath = arc(25, points=[[29,6,-4], [3,4,6], [1,1,7]]);
1693//   trans = path_sweep(shape, rpath, transforms=true);
1694//   outside = [for(i=[0:len(trans)-1]) trans[i]*scale(lerp(1,1.5,i/(len(trans)-1)))];
1695//   inside = [for(i=[len(trans)-1:-1:0]) trans[i]*scale(lerp(1.1,1.4,i/(len(trans)-1)))];
1696//   sweep(shape, concat(outside,inside),closed=true);
1697// Example(NoScales): An easier way to scale your model is to use the scale parameter.
1698//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[0,180], r=3));
1699//   path_sweep(pentagon(r=1), path3d(elliptic_arc), scale=2);
1700// Example(NoScales): Scaling only in the y direction of the profile (z direction in the model in this case)
1701//   elliptic_arc = xscale(2, p=arc($fn=64,angle=[0,180], r=3));
1702//   path_sweep(rect(2), path3d(elliptic_arc), scale=[1,2]);
1703// Example(NoScales): Specifying scale at every point for a closed path
1704//   N=64;
1705//   path = circle(r=5, $fn=64);
1706//   theta = lerpn(0,360,N,endpoint=false);
1707//   scale = [for(t=theta) sin(6*t)/5+1];
1708//   path_sweep(rect(2), path3d(path), closed=true, scale=scale);
1709// Example(Med,NoScales): Using path_sweep on a region
1710//   rgn1 = [for (d=[10:10:60]) circle(d=d,$fn=8)];
1711//   rgn2 = [square(30,center=false)];
1712//   rgn3 = [for (size=[10:10:20]) move([15,15],p=square(size=size, center=true))];
1713//   mrgn = union(rgn1,rgn2);
1714//   orgn = difference(mrgn,rgn3);
1715//   path_sweep(orgn,arc(r=40,angle=180));
1716// Example(Med,NoScales): A region with a twist
1717//   region = [for(i=pentagon(5)) move(i,p=circle(r=2,$fn=25))];
1718//   path_sweep(region,
1719//              circle(r=16,$fn=75),closed=true,
1720//              twist=360/5*2,symmetry=5);
1721// Example(Med,NoScales): Cutting a cylinder with a curved path.  Note that in this case, the incremental method produces just a slight twist but the natural method produces an extreme twist.  But manual specification produces no twist, as desired:
1722//   $fn=90;
1723//   r=8;
1724//   thickness=1;
1725//   len=21;
1726//   curve = [for(theta=[0:4:359])
1727//              [r*cos(theta), r*sin(theta), 10+sin(6*theta)]];
1728//   difference(){
1729//     cylinder(r=r, h=len);
1730//     down(.5)cylinder(r=r-thickness, h=len+1);
1731//     path_sweep(left(.05,square([1.1,1])), curve, closed=true,
1732//                method="manual", normal=UP);
1733//   }
1734
1735module path_sweep(shape, path, method="incremental", normal, closed, twist=0, twist_by_length=true, scale=1, scale_by_length=true,
1736                    symmetry=1, last_normal, tangent, uniform=true, relaxed=false, caps, style="min_edge", convexity=10,
1737                    anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull",profiles=false,width=1)
1738{
1739    dummy = assert(is_region(shape) || is_path(shape,2), "shape must be a 2D path or region");
1740    vnf = path_sweep(shape, path, method, normal, closed, twist, twist_by_length, scale, scale_by_length,
1741                    symmetry, last_normal, tangent, uniform, relaxed, caps, style);
1742
1743    if (profiles){
1744        assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"");
1745        tran = path_sweep(shape, path, method, normal, closed, twist, twist_by_length, scale, scale_by_length, 
1746                          symmetry, last_normal, tangent, uniform, relaxed,transforms=true);
1747        rshape = is_path(shape) ? [path3d(shape)]
1748                                : [for(s=shape) path3d(s)];
1749        attachable(anchor,spin,orient, vnf=vnf, extent=atype=="hull", cp=cp) {
1750            for(T=tran) stroke([for(part=rshape)apply(T,part)],width=width);
1751            children();
1752        }
1753    }
1754    else
1755      vnf_polyhedron(vnf,convexity=convexity,anchor=anchor, spin=spin, orient=orient, atype=atype, cp=cp)
1756          children();
1757}
1758
1759
1760function path_sweep(shape, path, method="incremental", normal, closed, twist=0, twist_by_length=true, scale=1, scale_by_length=true, 
1761                    symmetry=1, last_normal, tangent, uniform=true, relaxed=false, caps, style="min_edge", transforms=false,
1762                    anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull") =
1763  is_1region(path) ? path_sweep(shape=shape,path=path[0], method=method, normal=normal, closed=default(closed,true), 
1764                                twist=twist, scale=scale, scale_by_length=scale_by_length, twist_by_length=twist_by_length, symmetry=symmetry, last_normal=last_normal,
1765                                tangent=tangent, uniform=uniform, relaxed=relaxed, caps=caps, style=style, transforms=transforms,
1766                                anchor=anchor, cp=cp, spin=spin, orient=orient, atype=atype) :
1767  let(closed=default(closed,false))
1768  assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"")
1769  assert(!closed || twist % (360/symmetry)==0, str("For a closed sweep, twist must be a multiple of 360/symmetry = ",360/symmetry))
1770  assert(closed || symmetry==1, "symmetry must be 1 when closed is false")
1771  assert(is_integer(symmetry) && symmetry>0, "symmetry must be a positive integer")
1772  let(path = force_path(path))
1773  assert(is_path(path,[2,3]), "input path is not a 2D or 3D path")
1774  assert(!closed || !approx(path[0],last(path)), "Closed path includes start point at the end")
1775  assert((is_region(shape) || is_path(shape,2)) || (transforms && !(closed && method=="incremental")),"shape must be a 2d path or region")
1776  let(
1777    path = path3d(path),
1778    caps = is_def(caps) ? caps :
1779           closed ? false : true,
1780    capsOK = is_bool(caps) || is_bool_list(caps,2),
1781    fullcaps = is_bool(caps) ? [caps,caps] : caps,
1782    normalOK = is_undef(normal) || (method!="natural" && is_vector(normal,3))
1783                                || (method=="manual" && same_shape(normal,path)),
1784    scaleOK = scale==1 || ((is_num(scale) || is_vector(scale,2)) && !closed) || is_vector(scale,len(path)) || is_matrix(scale,len(path),2)
1785    
1786  )
1787  assert(normalOK,  method=="natural" ? "Cannot specify normal with the \"natural\" method"
1788                  : method=="incremental" ? "Normal with \"incremental\" method must be a 3-vector"
1789                  : str("Incompatible normal given.  Must be a 3-vector or a list of ",len(path)," 3-vectors"))
1790  assert(capsOK, "caps must be boolean or a list of two booleans")
1791  assert(!closed || !caps, "Cannot make closed shape with caps")
1792  assert(is_undef(normal) || (is_vector(normal) && len(normal)==3) || (is_path(normal) && len(normal)==len(path) && len(normal[0])==3), "Invalid normal specified")
1793  assert(is_undef(tangent) || (is_path(tangent) && len(tangent)==len(path) && len(tangent[0])==3), "Invalid tangent specified")
1794  assert(scaleOK,str("Incompatible or invalid scale",closed?" for closed path":"",": must be ", closed?"":"a scalar, a 2-vector, ",
1795                     "a vector of length ",len(path)," or a ",len(path),"x2 matrix of scales"))
1796  let(
1797    scale = !(is_num(scale) || is_vector(scale,2)) ? scale
1798          : let(s=is_num(scale) ? [scale,scale] : scale)
1799            !scale_by_length ? lerpn([1,1],s,len(path))
1800          : lerp([1,1],s, path_length_fractions(path,false)),
1801    scale_list = [for(s=scale) scale(s),if (closed) scale(scale[0])],
1802    tangents = is_undef(tangent) ? path_tangents(path,uniform=uniform,closed=closed) : [for(t=tangent) unit(t)],
1803    normal = is_path(normal) ? [for(n=normal) unit(n)] :
1804             is_def(normal) ? unit(normal) :
1805             method =="incremental" && abs(tangents[0].z) > 1/sqrt(2) ? BACK : UP,
1806    normals = is_path(normal) ? normal : repeat(normal,len(path)),
1807    tpathfrac = twist_by_length ? path_length_fractions(path, closed) : [for(i=[0:1:len(path)]) i / (len(path)-(closed?0:1))],
1808    spathfrac = scale_by_length ? path_length_fractions(path, closed) : [for(i=[0:1:len(path)]) i / (len(path)-(closed?0:1))],    
1809    L = len(path),
1810    unscaled_transform_list =
1811        method=="incremental" ?
1812          let(rotations =
1813                 [for( i  = 0,
1814                       ynormal = normal - (normal * tangents[0])*tangents[0],
1815                       rotation = frame_map(y=ynormal, z=tangents[0])
1816                         ;
1817                       i < len(tangents) + (closed?1:0) ;
1818                       rotation = i<len(tangents)-1+(closed?1:0)? rot(from=tangents[i],to=tangents[(i+1)%L])*rotation : undef,
1819                       i=i+1
1820                      )
1821                   rotation],
1822              // The mismatch is the inverse of the last transform times the first one for the closed case, or the inverse of the
1823              // desired final transform times the realized final transform in the open case.  Note that when closed==true the last transform
1824              // is a actually looped around and applies to the first point position, so if we got back exactly where we started
1825              // then it will be the identity, but we might have accumulated some twist which will show up as a rotation around the
1826              // X axis.  Similarly, in the closed==false case the desired and actual transformations can only differ in the twist,
1827              // so we can need to calculate the twist angle so we can apply a correction, which we distribute uniformly over the whole path.
1828              reference_rot = closed ? rotations[0] :
1829                           is_undef(last_normal) ? last(rotations) :
1830                             let(
1831                                 last_tangent = last(tangents),
1832                                 lastynormal = last_normal - (last_normal * last_tangent) * last_tangent
1833                             )
1834                           frame_map(y=lastynormal, z=last_tangent),
1835              mismatch = transpose(last(rotations)) * reference_rot,
1836              correction_twist = atan2(mismatch[1][0], mismatch[0][0]),
1837              // Spread out this extra twist over the whole sweep so that it doesn't occur
1838              // abruptly as an artifact at the last step.
1839              twistfix = correction_twist%(360/symmetry),
1840              adjusted_final = !closed ? undef :
1841                            translate(path[0]) * rotations[0] * zrot(-correction_twist+correction_twist%(360/symmetry)-twist)
1842          )  [for(i=idx(path)) translate(path[i]) * rotations[i] * zrot((twistfix-twist)*tpathfrac[i]), if(closed) adjusted_final] 
1843      : method=="manual" ?
1844              [for(i=[0:L-(closed?0:1)]) let(
1845                       ynormal = relaxed ? normals[i%L] : normals[i%L] - (normals[i%L] * tangents[i%L])*tangents[i%L],
1846                       znormal = relaxed ? tangents[i%L] - (normals[i%L] * tangents[i%L])*normals[i%L] : tangents[i%L],
1847                       rotation = frame_map(y=ynormal, z=znormal)
1848                   )
1849                   assert(approx(ynormal*znormal,0),str("Supplied normal is parallel to the path tangent at point ",i))
1850                   translate(path[i%L])*rotation*zrot(-twist*tpathfrac[i])
1851              ]
1852      : method=="cross"?
1853              let(
1854                  crossnormal_mid = [for(i=[(closed?0:1):L-(closed?1:2)])
1855                                       let(v=    cross(  select(path,i+1)-path[i], path[i]-select(path,i-1)),
1856                                           f=assert(norm(v)>EPSILON)
1857                                       )
1858                                       v
1859                                    ],
1860                  crossnormal = closed ? crossnormal_mid : [crossnormal_mid[0], each crossnormal_mid, last(crossnormal_mid)]
1861              )
1862              [for(i=[0:L-(closed?0:1)]) let(
1863                       rotation = frame_map(x=crossnormal[i%L], z=tangents[i%L])
1864                   )
1865                   translate(path[i%L])*rotation*zrot(-twist*tpathfrac[i])
1866                 ] 
1867      : method=="natural" ?   // map x axis of shape to the path normal, which points in direction of curvature
1868              let (pathnormal = path_normals(path, tangents, closed))
1869              assert(all_defined(pathnormal),"Natural normal vanishes on your curve, select a different method")
1870              let( testnormals = [for(i=[0:len(pathnormal)-1-(closed?1:2)]) pathnormal[i]*select(pathnormal,i+2)],
1871                   a=[for(i=idx(testnormals)) testnormals[i]<.5 ? echo(str("Big change at index ",i," pn=",pathnormal[i]," pn2= ",select(pathnormal,i+2))):0],
1872                   dummy = min(testnormals) < .5 ? echo("WARNING: ***** Abrupt change in normal direction.  Consider a different method in path_sweep() *****") :0
1873                 )
1874              [for(i=[0:L-(closed?0:1)]) let(
1875                       rotation = frame_map(x=pathnormal[i%L], z=tangents[i%L])
1876                   )
1877                   translate(path[i%L])*rotation*zrot(-twist*tpathfrac[i])
1878                 ] 
1879      : assert(false,"Unknown method or no method given"), // unknown method
1880    transform_list = v_mul(unscaled_transform_list, scale_list),
1881    ends_match = !closed ? true
1882                 : let( rshape = is_path(shape) ? [path3d(shape)]
1883                                                : [for(s=shape) path3d(s)]
1884                   )
1885                   are_regions_equal(apply(transform_list[0], rshape),
1886                                     apply(transform_list[L], rshape)),
1887    dummy = ends_match ? 0 : echo("WARNING: ***** The points do not match when closing the model in path_sweep() *****")
1888  )
1889  transforms ? transform_list
1890             : sweep(is_path(shape)?clockwise_polygon(shape):shape, transform_list, closed=false, caps=fullcaps,style=style,
1891                       anchor=anchor,cp=cp,spin=spin,orient=orient,atype=atype);
1892
1893
1894// Function&Module: path_sweep2d()
1895// Synopsis: Sweep a 2d polygon path along a 2d path allowing self-intersection. 
1896// SynTags: VNF, Geom
1897// Topics: Extrusion, Sweep, Paths
1898// See Also: linear_sweep(), rotate_sweep(), sweep(), spiral_sweep(), path_sweep(), offset_sweep()
1899// Usage: as module
1900//   path_sweep2d(shape, path, [closed], [caps], [quality], [style], [convexity=], [anchor=], [spin=], [orient=], [atype=], [cp=]) [ATTACHMENTS];
1901// Usage: as function
1902//   vnf = path_sweep2d(shape, path, [closed], [caps], [quality], [style], [anchor=], [spin=], [orient=], [atype=], [cp=]);
1903// Description:
1904//   Takes an input 2D polygon (the shape) and a 2d path, and constructs a polyhedron by sweeping the shape along the path.
1905//   When run as a module returns the polyhedron geometry.  When run as a function returns a VNF.
1906//   .
1907//   See {{path_sweep()}} for more details on how the sweep operation works and for introductory examples.
1908//   This 2d version is different because local self-intersections (creases in the output) are allowed and do not produce CGAL errors.
1909//   This is accomplished by using offset() calculations, which are more expensive than simply copying the shape along
1910//   the path, so if you do not have local self-intersections, use {{path_sweep()}} instead.  If xmax is the largest x value (in absolute value)
1911//   of the shape, then path_sweep2d() will work as long as the offset of `path` exists at `delta=xmax`.  If the offset vanishes, as in the
1912//   case of a circle offset by more than its radius, then you will get an error about a degenerate offset.
1913//   Note that global self-intersections will still give rise to CGAL errors.  You should be able to handle these by partitioning your model.  The y axis of the
1914//   shape is mapped to the z axis in the swept polyhedron, and no twisting can occur.
1915//   The quality parameter is passed to offset to determine the offset quality.
1916// Arguments:
1917//   shape = a 2D polygon describing the shape to be swept
1918//   path = a 2D path giving the path to sweep over
1919//   closed = path is a closed loop.  Default: false
1920//   caps = true to create endcap faces when closed is false.  Can be a length 2 boolean array.  Default is true if closed is false.
1921//   quality = quality of offset used in calculation.  Default: 1
1922//   style = vnf_vertex_array style.  Default: "min_edge"
1923//   ---
1924//   convexity = convexity parameter for polyhedron (module only)  Default: 10
1925//   anchor = Translate so anchor point is at the origin.  Default: "origin"
1926//   spin = Rotate this many degrees around Z axis after anchor.  Default: 0
1927//   orient = Vector to rotate top towards after spin
1928//   atype = Select "hull" or "intersect" anchor types.  Default: "hull"
1929//   cp = Centerpoint for determining "intersect" anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: "centroid"
1930// Anchor Types:
1931//   "hull" = Anchors to the virtual convex hull of the shape.
1932//   "intersect" = Anchors to the surface of the shape.
1933// Example: Sine wave example with self-intersections at each peak.  This would fail with path_sweep().
1934//   sinewave = [for(i=[-30:10:360*2+30]) [i/40,3*sin(i)]];
1935//   path_sweep2d(circle(r=3,$fn=15), sinewave);
1936// Example: The ends can look weird if they are in a place where self intersection occurs.  This is a natural result of how offset behaves at ends of a path.
1937//   coswave = [for(i=[0:10:360*1.5]) [i/40,3*cos(i)]];
1938//   zrot(-20)
1939//     path_sweep2d( circle(r=3,$fn=15), coswave);
1940// Example: This closed path example works ok as long as the hole in the center remains open.
1941//   ellipse = yscale(3,p=circle(r=3,$fn=120));
1942//   path_sweep2d(circle(r=2.5,$fn=32), reverse(ellipse), closed=true);
1943// Example: When the hole is closed a global intersection renders the model invalid.  You can fix this by taking the union of the two (valid) halves.
1944//   ellipse = yscale(3,p=circle(r=3,$fn=120));
1945//   L = len(ellipse);
1946//   path_sweep2d(circle(r=3.25, $fn=32), select(ellipse,floor(L*.2),ceil(L*.8)),closed=false);
1947//   path_sweep2d(circle(r=3.25, $fn=32), select(ellipse,floor(L*.7),ceil(L*.3)),closed=false);
1948
1949function path_sweep2d(shape, path, closed=false, caps, quality=1, style="min_edge",
1950                      anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull") =
1951   let(
1952        caps = is_def(caps) ? caps
1953             : closed ? false : true,
1954        capsOK = is_bool(caps) || is_bool_list(caps,2),
1955        fullcaps = is_bool(caps) ? [caps,caps] : caps,
1956        shape = force_path(shape,"shape"),
1957        path = force_path(path)
1958   )
1959   assert(is_path(shape,2), "shape must be a 2D path")
1960   assert(is_path(path,2), "path must be a 2D path")
1961   assert(capsOK, "caps must be boolean or a list of two booleans")
1962   assert(!closed || !caps, "Cannot make closed shape with caps")
1963   let(
1964        profile = ccw_polygon(shape),
1965        flip = closed && is_polygon_clockwise(path) ? -1 : 1,
1966        path = flip ? reverse(path) : path,
1967        proflist= transpose(
1968                     [for(pt = profile)
1969                        let(
1970                            ofs = offset(path, delta=-flip*pt.x, return_faces=true,closed=closed, quality=quality),
1971                            map = column(_ofs_vmap(ofs,closed=closed),1)
1972                        )
1973                        select(path3d(ofs[0],pt.y),map)
1974                      ]
1975                  ),
1976        vnf = vnf_vertex_array([
1977                         each proflist,
1978                         if (closed) proflist[0]
1979                        ],cap1=fullcaps[0],cap2=fullcaps[1],col_wrap=true,style=style)
1980   )
1981   reorient(anchor,spin,orient,vnf=vnf,p=vnf,extent=atype=="hull",cp=cp);
1982
1983
1984module path_sweep2d(profile, path, closed=false, caps, quality=1, style="min_edge", convexity=10,
1985                    anchor="origin", cp="centroid", spin=0, orient=UP, atype="hull")
1986{
1987   vnf = path_sweep2d(profile, path, closed, caps, quality, style);
1988   vnf_polyhedron(vnf,convexity=convexity,anchor=anchor, spin=spin, orient=orient, atype=atype, cp=cp)
1989        children();
1990}
1991
1992// Extract vertex mapping from offset face list.  The output of this function
1993// is a list of pairs [i,j] where i is an index into the parent curve and j is
1994// an index into the offset curve.  It would probably make sense to rewrite
1995// offset() to return this instead of the face list and have offset_sweep
1996// use this input to assemble the faces it needs.
1997
1998function _ofs_vmap(ofs,closed=false) =
1999    let(   // Caclulate length of the first (parent) curve
2000        firstlen = max(flatten(ofs[1]))+1-len(ofs[0])
2001    )
2002    [
2003     for(entry=ofs[1]) _ofs_face_edge(entry,firstlen),
2004     if (!closed) _ofs_face_edge(last(ofs[1]),firstlen,second=true)
2005    ];
2006
2007
2008// Extract first (default) or second edge that connects the parent curve to its offset.  The first input
2009// face is a list of 3 or 4 vertices as indices into the two curves where the parent curve vertices are
2010// numbered from 0 to firstlen-1 and the offset from firstlen and up.  The firstlen pararameter is used
2011// to determine which curve the vertices belong to and to remove the offset so that the return gives
2012// the index into each curve with a 0 base.
2013function _ofs_face_edge(face,firstlen,second=false) =
2014   let(
2015       itry = min_index(face),
2016       i = select(face,itry-1)<firstlen ? itry-1:itry,
2017       edge1 = select(face,[i,i-1]),
2018       edge2 = select(face,i+1)<firstlen ? select(face,[i+1,i+2])
2019                                         : select(face,[i,i+1])
2020   )
2021   (second ? edge2 : edge1)-[0,firstlen];
2022
2023
2024
2025// Function&Module: sweep()
2026// Synopsis: Construct a 3d object from arbitrary transformations of a 2d polygon path.
2027// SynTags: VNF, Geom
2028// Topics: Extrusion, Sweep, Paths
2029// See Also: linear_sweep(), rotate_sweep(), spiral_sweep(), path_sweep(), path_sweep2d(), offset_sweep()
2030// Usage: As Module
2031//   sweep(shape, transforms, [closed], [caps], [style], [convexity=], [anchor=], [spin=], [orient=], [atype=]) [ATTACHMENTS];
2032// Usage: As Function
2033//   vnf = sweep(shape, transforms, [closed], [caps], [style], [anchor=], [spin=], [orient=], [atype=]);
2034// Description:
2035//   The input `shape` must be a non-self-intersecting 2D polygon or region, and `transforms`
2036//   is a list of 4x4 transformation matrices.  The sweep algorithm applies each transformation in sequence
2037//   to the shape input and links the resulting polygons together to form a polyhedron.
2038//   If `closed=true` then the first and last transformation are linked together.
2039//   The `caps` parameter controls whether the ends of the shape are closed.
2040//   As a function, returns the VNF for the polyhedron.  As a module, computes the polyhedron.
2041//   .
2042//   Note that this is a very powerful, general framework for producing polyhedra.  It is important
2043//   to ensure that your resulting polyhedron does not include any self-intersections, or it will
2044//   be invalid and will generate CGAL errors.  If you get such errors, most likely you have an
2045//   overlooked self-intersection.  Note also that the errors will not occur when your shape is alone
2046//   in your model, but will arise if you add a second object to the model.  This may mislead you into
2047//   thinking the second object caused a problem.  Even adding a simple cube to the model will reveal the problem.
2048// Arguments:
2049//   shape = 2d path or region, describing the shape to be swept.
2050//   transforms = list of 4x4 matrices to apply
2051//   closed = set to true to form a closed (torus) model.  Default: false
2052//   caps = true to create endcap faces when closed is false.  Can be a singe boolean to specify endcaps at both ends, or a length 2 boolean array.  Default is true if closed is false.
2053//   style = vnf_vertex_array style.  Default: "min_edge"
2054//   ---
2055//   convexity = convexity setting for use with polyhedron.  (module only) Default: 10
2056//   cp = Centerpoint for determining "intersect" anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 3D point.  Default: "centroid"
2057//   atype = Select "hull" or "intersect" anchor types.  Default: "hull"
2058//   anchor = Translate so anchor point is at the origin. Default: "origin"
2059//   spin = Rotate this many degrees around Z axis after anchor. Default: 0
2060//   orient = Vector to rotate top towards after spin  (module only)
2061// Anchor Types:
2062//   "hull" = Anchors to the virtual convex hull of the shape.
2063//   "intersect" = Anchors to the surface of the shape.
2064// Example(VPR=[45,0,74],VPD=175,VPT=[-3.8,12.4,19]): A bent object that also changes shape along its length.
2065//   radius = 75;
2066//   angle = 40;
2067//   shape = circle(r=5,$fn=32);
2068//   T = [for(i=[0:25]) xrot(-angle*i/25,cp=[0,radius,0])*scale([1+i/25, 2-i/25,1])];
2069//   sweep(shape,T);
2070// Example: This is the "sweep-drop" example from list-comprehension-demos.
2071//   function drop(t) = 100 * 0.5 * (1 - cos(180 * t)) * sin(180 * t) + 1;
2072//   function path(t) = [0, 0, 80 + 80 * cos(180 * t)];
2073//   function rotate(t) = 180 * pow((1 - t), 3);
2074//   step = 0.01;
2075//   path_transforms = [for (t=[0:step:1-step]) translate(path(t)) * zrot(rotate(t)) * scale([drop(t), drop(t), 1])];
2076//   sweep(circle(1, $fn=12), path_transforms);
2077// Example: Another example from list-comprehension-demos
2078//   function f(x) = 3 - 2.5 * x;
2079//   function r(x) = 2 * 180 * x * x * x;
2080//   pathstep = 1;
2081//   height = 100;
2082//   shape_points = subdivide_path(square(10),40,closed=true);
2083//   path_transforms = [for (i=[0:pathstep:height]) let(t=i/height) up(i) * scale([f(t),f(t),i]) * zrot(r(t))];
2084//   sweep(shape_points, path_transforms);
2085// Example: Twisted container.  Note that this technique doesn't create a fixed container wall thickness.
2086//   shape = subdivide_path(square(30,center=true), 40, closed=true);
2087//   outside = [for(i=[0:24]) up(i)*rot(i)*scale(1.25*i/24+1)];
2088//   inside = [for(i=[24:-1:2]) up(i)*rot(i)*scale(1.2*i/24+1)];
2089//   sweep(shape, concat(outside,inside));
2090
2091function sweep(shape, transforms, closed=false, caps, style="min_edge",
2092               anchor="origin", cp="centroid", spin=0, orient=UP, atype="hull") =
2093    assert(is_consistent(transforms, ident(4)), "Input transforms must be a list of numeric 4x4 matrices in sweep")
2094    assert(is_path(shape,2) || is_region(shape), "Input shape must be a 2d path or a region.")
2095    let(
2096        caps = is_def(caps) ? caps :
2097            closed ? false : true,
2098        capsOK = is_bool(caps) || is_bool_list(caps,2),
2099        fullcaps = is_bool(caps) ? [caps,caps] : caps
2100    )
2101    assert(len(transforms), "transformation must be length 2 or more")
2102    assert(capsOK, "caps must be boolean or a list of two booleans")
2103    assert(!closed || !caps, "Cannot make closed shape with caps")
2104    is_region(shape)? let(
2105        regions = region_parts(shape),
2106        rtrans = reverse(transforms),
2107        vnfs = [
2108            for (rgn=regions) each [
2109                for (path=rgn)
2110                    sweep(path, transforms, closed=closed, caps=false),
2111                if (fullcaps[0]) vnf_from_region(rgn, transform=transforms[0], reverse=true),
2112                if (fullcaps[1]) vnf_from_region(rgn, transform=last(transforms)),
2113            ],
2114        ],
2115        vnf = vnf_join(vnfs)
2116    ) vnf :
2117    assert(len(shape)>=3, "shape must be a path of at least 3 non-colinear points")
2118    vnf_vertex_array([for(i=[0:len(transforms)-(closed?0:1)]) apply(transforms[i%len(transforms)],path3d(shape))],
2119                     cap1=fullcaps[0],cap2=fullcaps[1],col_wrap=true,style=style);
2120
2121
2122module sweep(shape, transforms, closed=false, caps, style="min_edge", convexity=10,
2123             anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull")
2124{
2125    vnf = sweep(shape, transforms, closed, caps, style);
2126    vnf_polyhedron(vnf, convexity=convexity, anchor=anchor, spin=spin, orient=orient, atype=atype, cp=cp)
2127        children();
2128}
2129
2130
2131
2132// Section: Functions for resampling and slicing profile lists
2133
2134// Function: subdivide_and_slice()
2135// Synopsis: Resample list of paths to have the same point count and interpolate additional paths. 
2136// SynTags: PathList
2137// Topics: Paths, Path Subdivision
2138// See Also: slice_profiles()
2139// Usage:
2140//   newprof = subdivide_and_slice(profiles, slices, [numpoints], [method], [closed]);
2141// Description:
2142//   Subdivides the input profiles to have length `numpoints` where `numpoints` must be at least as
2143//   big as the largest input profile.  By default `numpoints` is set equal to the length of the
2144//   largest profile.  You can set `numpoints="lcm"` to sample to the least common multiple of all
2145//   curves, which will avoid sampling artifacts but may produce a huge output.  After subdivision,
2146//   profiles are sliced.
2147// Arguments:
2148//   profiles = profiles to operate on
2149//   slices = number of slices to insert between each pair of profiles.  May be a vector
2150//   numpoints = number of points after sampling.
2151//   method = method used for calling {{subdivide_path()}}, either `"length"` or `"segment"`.  Default: `"length"`
2152//   closed = the first and last profile are connected.  Default: false
2153function subdivide_and_slice(profiles, slices, numpoints, method="length", closed=false) =
2154  let(
2155    maxsize = max_length(profiles),
2156    numpoints = is_undef(numpoints) ? maxsize :
2157                numpoints == "lcm" ? lcmlist([for(p=profiles) len(p)]) :
2158                is_num(numpoints) ? round(numpoints) : undef
2159  )
2160  assert(is_def(numpoints), "Parameter numpoints must be \"max\", \"lcm\" or a positive number")
2161  assert(numpoints>=maxsize, "Number of points requested is smaller than largest profile")
2162  let(fixpoly = [for(poly=profiles) subdivide_path(poly, numpoints,method=method)])
2163  slice_profiles(fixpoly, slices, closed);
2164
2165
2166
2167// Function: slice_profiles()
2168// Synopsis: Linearly interpolates between path profiles.
2169// SynTags: PathList
2170// Topics: Paths, Path Subdivision
2171// See Also: subdivide_and_slice()
2172// Usage:
2173//   profs = slice_profiles(profiles, slices, [closed]);
2174// Description:
2175//   Given an input list of profiles, linearly interpolate between each pair to produce a
2176//   more finely sampled list.  The parameters `slices` specifies the number of slices to
2177//   be inserted between each pair of profiles and can be a number or a list.
2178// Arguments:
2179//   profiles = list of paths to operate on.  They must be lists of the same shape and length.
2180//   slices = number of slices to insert between each pair, or a list to vary the number inserted.
2181//   closed = set to true if last profile connects to first one.  Default: false
2182function slice_profiles(profiles,slices,closed=false) =
2183  assert(is_num(slices) || is_list(slices))
2184  let(listok = !is_list(slices) || len(slices)==len(profiles)-(closed?0:1))
2185  assert(listok, "Input slices to slice_profiles is a list with the wrong length")
2186  let(
2187    count = is_num(slices) ? repeat(slices,len(profiles)-(closed?0:1)) : slices,
2188    slicelist = [for (i=[0:len(profiles)-(closed?1:2)])
2189      each lerpn(profiles[i], select(profiles,i+1), count[i]+1, false)
2190    ]
2191  )
2192  concat(slicelist, closed?[]:[profiles[len(profiles)-1]]);
2193
2194
2195
2196function _closest_angle(alpha,beta) =
2197    is_vector(beta) ? [for(entry=beta) _closest_angle(alpha,entry)]
2198  : beta-alpha > 180 ? beta - ceil((beta-alpha-180)/360) * 360
2199  : beta-alpha < -180 ? beta + ceil((alpha-beta-180)/360) * 360
2200  : beta;
2201
2202
2203// Smooth data with N point moving average.  If angle=true handles data as angles.
2204// If closed=true assumes last point is adjacent to the first one.
2205// If closed=false pads data with left/right value (probably wrong behavior...should do linear interp)
2206function _smooth(data,len,closed=false,angle=false) =
2207  let(  halfwidth = floor(len/2),
2208        result = closed ? [for(i=idx(data))
2209                           let(
2210                             window = angle ? _closest_angle(data[i],select(data,i-halfwidth,i+halfwidth))
2211                                            : select(data,i-halfwidth,i+halfwidth)
2212                           )
2213                           mean(window)]
2214               : [for(i=idx(data))
2215                   let(
2216                       window = select(data,max(i-halfwidth,0),min(i+halfwidth,len(data)-1)),
2217                       left = i-halfwidth<0,
2218                       pad = left ? data[0] : last(data)
2219                   )
2220                   sum(window)+pad*(len-len(window))] / len
2221   )
2222   result;
2223
2224
2225// Function: rot_resample()
2226// Synopsis: Resample a list of rotation operators. 
2227// SynTags: MatList
2228// Topics: Matrices, Interpolation, Rotation
2229// See Also: subdivide_and_slice(), slice_profiles()
2230// Usage:
2231//   rlist = rot_resample(rotlist, n, [method=], [twist=], [scale=], [smoothlen=], [long=], [turns=], [closed=])
2232// Description:
2233//   Takes as input a list of rotation matrices in 3d.  Produces as output a resampled
2234//   list of rotation operators (4x4 matrixes) suitable for use with sweep().  You can optionally apply twist to
2235//   the output with the twist parameter, which is either a scalar to apply a uniform
2236//   overall twist, or a vector to apply twist non-uniformly.  Similarly you can apply
2237//   scaling either overall or with a vector.  The smoothlen parameter applies smoothing
2238//   to the twist and scaling to prevent abrupt changes.  This is done by a moving average
2239//   of the smoothing or scaling values.  The default of 1 means no smoothing.  The long parameter causes
2240//   the interpolation to be done the "long" way around the rotation instead of the short way.
2241//   Note that the rotation matrix cannot distinguish which way you rotate, only the place you
2242//   end after rotation.  Another ambiguity arises if your rotation is more than 360 degrees.
2243//   You can add turns with the turns parameter, so giving turns=1 will add 360 degrees to the
2244//   rotation so it completes one full turn plus the additional rotation given my the transform.
2245//   You can give long as a scalar or as a vector.  Finally if closed is true then the
2246//   resampling will connect back to the beginning.
2247//   .
2248//   The default is to resample based on the length of the arc defined by each rotation operator.  This produces
2249//   uniform sampling over all of the transformations.  It requires that each rotation has nonzero length.
2250//   In this case n specifies the total number of samples.  If you set method to "count" then you get
2251//   n samples for each transform.  You can set n to a vector to vary the samples at each step.
2252// Arguments:
2253//   rotlist = list of rotation operators in 3d to resample
2254//   n = Number of rotations to produce as output when method is "length" or number for each transformation if method is "count".  Can be a vector when method is "count"
2255//   ---
2256//   method = sampling method, either "length" or "count"
2257//   twist = scalar or vector giving twist to add overall or at each rotation.  Default: none
2258//   scale = scalar or vector giving scale factor to add overall or at each rotation.  Default: none
2259//   smoothlen = amount of smoothing to apply to scaling and twist.  Should be an odd integer.  Default: 1
2260//   long = resample the "long way" around the rotation, a boolean or list of booleans.  Default: false
2261//   turns = add extra turns.  If a scalar adds the turns to every rotation, or give a vector.  Default: 0
2262//   closed = if true then the rotation list is treated as closed.  Default: false
2263// Example(3D): Resampling the arc from a compound rotation with translations thrown in.
2264//   tran = rot_resample([ident(4), back(5)*up(4)*xrot(-10)*zrot(-20)*yrot(117,cp=[10,0,0])], n=25);
2265//   sweep(circle(r=1,$fn=3), tran);
2266// Example(3D): Applying a scale factor
2267//   tran = rot_resample([ident(4), back(5)*up(4)*xrot(-10)*zrot(-20)*yrot(117,cp=[10,0,0])], n=25, scale=2);
2268//   sweep(circle(r=1,$fn=3), tran);
2269// Example(3D): Applying twist
2270//   tran = rot_resample([ident(4), back(5)*up(4)*xrot(-10)*zrot(-20)*yrot(117,cp=[10,0,0])], n=25, twist=60);
2271//   sweep(circle(r=1,$fn=3), tran);
2272// Example(3D): Going the long way
2273//   tran = rot_resample([ident(4), back(5)*up(4)*xrot(-10)*zrot(-20)*yrot(117,cp=[10,0,0])], n=25, long=true);
2274//   sweep(circle(r=1,$fn=3), tran);
2275// Example(3D): Getting transformations from turtle3d
2276//   include<BOSL2/turtle3d.scad>
2277//   tran=turtle3d(["arcsteps",1,"up", 10, "arczrot", 10,170],transforms=true);
2278//   sweep(circle(r=1,$fn=3),rot_resample(tran, n=40));
2279// Example(3D): If you specify a larger angle in turtle you need to use the long argument
2280//   include<BOSL2/turtle3d.scad>
2281//   tran=turtle3d(["arcsteps",1,"up", 10, "arczrot", 10,270],transforms=true);
2282//   sweep(circle(r=1,$fn=3),rot_resample(tran, n=40,long=true));
2283// Example(3D): And if the angle is over 360 you need to add turns to get the right result.  Note long is false when the remaining angle after subtracting full turns is below 180:
2284//   include<BOSL2/turtle3d.scad>
2285//   tran=turtle3d(["arcsteps",1,"up", 10, "arczrot", 10,90+360],transforms=true);
2286//   sweep(circle(r=1,$fn=3),rot_resample(tran, n=40,long=false,turns=1));
2287// Example(3D): Here the remaining angle is 270, so long must be set to true
2288//   include<BOSL2/turtle3d.scad>
2289//   tran=turtle3d(["arcsteps",1,"up", 10, "arczrot", 10,270+360],transforms=true);
2290//   sweep(circle(r=1,$fn=3),rot_resample(tran, n=40,long=true,turns=1));
2291// Example(3D): Note the visible line at the scale transition
2292//   include<BOSL2/turtle3d.scad>
2293//   tran = turtle3d(["arcsteps",1,"arcup", 10, 90, "arcdown", 10, 90], transforms=true);
2294//   rtran = rot_resample(tran,200,scale=[1,6]);
2295//   sweep(circle(1,$fn=32),rtran);
2296// Example(3D): Observe how using a large smoothlen value eases that transition
2297//   include<BOSL2/turtle3d.scad>
2298//   tran = turtle3d(["arcsteps",1,"arcup", 10, 90, "arcdown", 10, 90], transforms=true);
2299//   rtran = rot_resample(tran,200,scale=[1,6],smoothlen=17);
2300//   sweep(circle(1,$fn=32),rtran);
2301// Example(3D): A similar issues can arise with twist, where a "line" is visible at the transition
2302//   include<BOSL2/turtle3d.scad>
2303//   tran = turtle3d(["arcsteps", 1, "arcup", 10, 90, "move", 10], transforms=true,state=[1,-.5,0]);
2304//   rtran = rot_resample(tran,100,twist=[0,60],smoothlen=1);
2305//   sweep(subdivide_path(rect([3,3]),40),rtran);
2306// Example(3D): Here's the smoothed twist transition
2307//   include<BOSL2/turtle3d.scad>
2308//   tran = turtle3d(["arcsteps", 1, "arcup", 10, 90, "move", 10], transforms=true,state=[1,-.5,0]);
2309//   rtran = rot_resample(tran,100,twist=[0,60],smoothlen=17);
2310//   sweep(subdivide_path(rect([3,3]),40),rtran);
2311// Example(3D): Toothed belt based on a list-comprehension-demos example.  This version has a smoothed twist transition.  Try changing smoothlen to 1 to see the more abrupt transition that occurs without smoothing.
2312//   include<BOSL2/turtle3d.scad>
2313//   r_small = 19;       // radius of small curve
2314//   r_large = 46;       // radius of large curve
2315//   flat_length = 100;  // length of flat belt section
2316//   teeth=42;           // number of teeth
2317//   belt_width = 12;
2318//   tooth_height = 9;
2319//   belt_thickness = 3;
2320//   angle = 180 - 2*atan((r_large-r_small)/flat_length);
2321//   beltprofile = path3d(subdivide_path(
2322//                   square([belt_width, belt_thickness],anchor=FWD),
2323//                   20));
2324//   beltrots =
2325//     turtle3d(["arcsteps",1,
2326//               "move", flat_length,
2327//               "arcleft", r_small, angle,
2328//               "move", flat_length,
2329//     // Closing path will be interpolated
2330//     //        "arcleft", r_large, 360-angle
2331//              ],transforms=true);
2332//   beltpath = rot_resample(beltrots,teeth*4,
2333//                           twist=[180,0,-180,0],
2334//                           long=[false,false,false,true],
2335//                           smoothlen=15,closed=true);
2336//   belt = [for(i=idx(beltpath))
2337//             let(tooth = floor((i+$t*4)/2)%2)
2338//             apply(beltpath[i]*
2339//                     yscale(tooth
2340//                            ? tooth_height/belt_thickness
2341//                            : 1),
2342//                   beltprofile)
2343//          ];
2344//   skin(belt,slices=0,closed=true);
2345function rot_resample(rotlist,n,twist,scale,smoothlen=1,long=false,turns=0,closed=false,method="length") =
2346    assert(is_int(smoothlen) && smoothlen>0 && smoothlen%2==1, "smoothlen must be a positive odd integer")
2347    assert(method=="length" || method=="count")
2348    let(tcount = len(rotlist) + (closed?0:-1))
2349    assert(method=="count" || is_int(n), "n must be an integer when method is \"length\"")
2350    assert(is_int(n) || is_vector(n,tcount), str("n must be scalar or vector with length ",tcount))
2351    let(
2352          count = method=="length" ? (closed ? n+1 : n)
2353                                   : (is_vector(n) ? sum(n) : tcount*n)+1  //(closed?0:1)
2354    )
2355    assert(is_bool(long) || len(long)==tcount,str("Input long must be a scalar or have length ",tcount))
2356    let(
2357        long = force_list(long,tcount),
2358        turns = force_list(turns,tcount),
2359        T = [for(i=[0:1:tcount-1]) rot_inverse(rotlist[i])*select(rotlist,i+1)],
2360        parms = [for(i=idx(T))
2361                    let(tparm = rot_decode(T[i],long[i]))
2362                    [tparm[0]+turns[i]*360,tparm[1],tparm[2],tparm[3]]
2363                ],
2364        radius = [for(i=idx(parms)) norm(parms[i][2])],
2365        length = [for(i=idx(parms)) norm([norm(parms[i][3]), parms[i][0]/360*2*PI*radius[i]])]
2366    )
2367    assert(method=="count" || all_positive(length),
2368           "Rotation list includes a repeated entry or a rotation around the origin, not allowed when method=\"length\"")
2369    let(
2370        cumlen = [0, each cumsum(length)],
2371        totlen = last(cumlen),
2372        stepsize = totlen/(count-1),
2373        samples = method=="count"
2374                  ? let( n = force_list(n,tcount))
2375                    [for(N=n) lerpn(0,1,N,endpoint=false)]
2376                  :[for(i=idx(parms))
2377                    let(
2378                        remainder = cumlen[i] % stepsize,
2379                        offset = remainder==0 ? 0
2380                                              : stepsize-remainder,
2381                        num = ceil((length[i]-offset)/stepsize)
2382                    )
2383                    count(num,offset,stepsize)/length[i]],
2384         twist = first_defined([twist,0]),
2385         scale = first_defined([scale,1]),
2386         needlast = !approx(last(last(samples)),1),
2387         sampletwist = is_num(twist) ? lerpn(0,twist,count)
2388                     : let(
2389                          cumtwist = [0,each cumsum(twist)]
2390                      )
2391                      [for(i=idx(parms)) each lerp(cumtwist[i],cumtwist[i+1],samples[i]),
2392                      if (needlast) last(cumtwist)
2393                      ],
2394         samplescale = is_num(scale) ? lerp(1,scale,lerpn(0,1,count))
2395                     : let(
2396                          cumscale = [1,each cumprod(scale)]
2397                      )
2398                      [for(i=idx(parms)) each lerp(cumscale[i],cumscale[i+1],samples[i]),
2399                       if (needlast) last(cumscale)],
2400         smoothtwist = _smooth(closed?select(sampletwist,0,-2):sampletwist,smoothlen,closed=closed,angle=true),
2401         smoothscale = _smooth(samplescale,smoothlen,closed=closed),
2402         interpolated = [
2403           for(i=idx(parms))
2404             each [for(u=samples[i]) rotlist[i] * move(u*parms[i][3]) * rot(a=u*parms[i][0],v=parms[i][1],cp=parms[i][2])],
2405           if (needlast) last(rotlist)
2406         ]
2407     )
2408     [for(i=idx(interpolated,e=closed?-2:-1)) interpolated[i]*zrot(smoothtwist[i])*scale(smoothscale[i])];
2409
2410
2411
2412
2413
2414//////////////////////////////////////////////////////////////////
2415//
2416// Minimum Distance Mapping using Dynamic Programming
2417//
2418// Given inputs of a two polygons, computes a mapping between their vertices that minimizes the sum the sum of
2419// the distances between every matched pair of vertices.  The algorithm uses dynamic programming to calculate
2420// the optimal mapping under the assumption that poly1[0] <-> poly2[0].  We then rotate through all the
2421// possible indexings of the longer polygon.  The theoretical run time is quadratic in the longer polygon and
2422// linear in the shorter one.
2423//
2424// The top level function, _skin_distance_match(), cycles through all the of the indexings of the larger
2425// polygon, computes the optimal value for each indexing, and chooses the overall best result.  It uses
2426// _dp_extract_map() to thread back through the dynamic programming array to determine the actual mapping, and
2427// then converts the result to an index repetition count list, which is passed to repeat_entries().
2428//
2429// The function _dp_distance_array builds up the rows of the dynamic programming matrix with reference
2430// to the previous rows, where `tdist` holds the total distance for a given mapping, and `map`
2431// holds the information about which path was optimal for each position.
2432//
2433// The function _dp_distance_row constructs each row of the dynamic programming matrix in the usual
2434// way where entries fill in based on the three entries above and to the left.  Note that we duplicate
2435// entry zero so account for wrap-around at the ends, and we initialize the distance to zero to avoid
2436// double counting the length of the 0-0 pair.
2437//
2438// This function builds up the dynamic programming distance array where each entry in the
2439// array gives the optimal distance for aligning the corresponding subparts of the two inputs.
2440// When the array is fully populated, the bottom right corner gives the minimum distance
2441// for matching the full input lists.  The `map` array contains a the three key values for the three
2442// directions, where _MAP_DIAG means you map the next vertex of `big` to the next vertex of `small`,
2443// _MAP_LEFT means you map the next vertex of `big` to the current vertex of `small`, and _MAP_UP
2444// means you map the next vertex of `small` to the current vertex of `big`.
2445//
2446// Return value is [min_distance, map], where map is the array that is used to extract the actual
2447// vertex map.
2448
2449_MAP_DIAG = 0;
2450_MAP_LEFT = 1;
2451_MAP_UP = 2;
2452
2453/*
2454function _dp_distance_array(small, big, abort_thresh=1/0, small_ind=0, tdist=[], map=[]) =
2455   small_ind == len(small)+1 ? [tdist[len(tdist)-1][len(big)-1], map] :
2456   let( newrow = _dp_distance_row(small, big, small_ind, tdist) )
2457   min(newrow[0]) > abort_thresh ? [tdist[len(tdist)-1][len(big)-1],map] :
2458   _dp_distance_array(small, big, abort_thresh, small_ind+1, concat(tdist, [newrow[0]]), concat(map, [newrow[1]]));
2459*/
2460
2461
2462function _dp_distance_array(small, big, abort_thresh=1/0) =
2463   [for(
2464        small_ind = 0,
2465        tdist = [],
2466        map = []
2467           ;
2468        small_ind<=len(small)+1
2469           ;
2470        newrow =small_ind==len(small)+1 ? [0,0,0] :  // dummy end case
2471                           _dp_distance_row(small,big,small_ind,tdist),
2472        tdist = concat(tdist, [newrow[0]]),
2473        map = concat(map, [newrow[1]]),
2474        small_ind = min(newrow[0])>abort_thresh ? len(small)+1 : small_ind+1
2475       )
2476     if (small_ind==len(small)+1) each [tdist[len(tdist)-1][len(big)], map]];
2477                                     //[tdist,map]];
2478
2479
2480function _dp_distance_row(small, big, small_ind, tdist) =
2481                    // Top left corner is zero because it gets counted at the end in bottom right corner
2482   small_ind == 0 ? [cumsum([0,for(i=[1:len(big)]) norm(big[i%len(big)]-small[0])]), repeat(_MAP_LEFT,len(big)+1)] :
2483   [for(big_ind=1,
2484       newrow=[ norm(big[0] - small[small_ind%len(small)]) + tdist[small_ind-1][0] ],
2485       newmap = [_MAP_UP]
2486         ;
2487       big_ind<=len(big)+1
2488         ;
2489       costs = big_ind == len(big)+1 ? [0] :    // handle extra iteration
2490                             [tdist[small_ind-1][big_ind-1],  // diag
2491                              newrow[big_ind-1],              // left
2492                              tdist[small_ind-1][big_ind]],   // up
2493       newrow = concat(newrow, [min(costs)+norm(big[big_ind%len(big)]-small[small_ind%len(small)])]),
2494       newmap = concat(newmap, [min_index(costs)]),
2495       big_ind = big_ind+1
2496   ) if (big_ind==len(big)+1) each [newrow,newmap]];
2497
2498
2499function _dp_extract_map(map) =
2500      [for(
2501           i=len(map)-1,
2502           j=len(map[0])-1,
2503           smallmap=[],
2504           bigmap = []
2505              ;
2506           j >= 0
2507              ;
2508           advance_i = map[i][j]==_MAP_UP || map[i][j]==_MAP_DIAG,
2509           advance_j = map[i][j]==_MAP_LEFT || map[i][j]==_MAP_DIAG,
2510           i = i - (advance_i ? 1 : 0),
2511           j = j - (advance_j ? 1 : 0),
2512           bigmap = concat( [j%(len(map[0])-1)] ,  bigmap),
2513           smallmap = concat( [i%(len(map)-1)]  , smallmap)
2514          )
2515        if (i==0 && j==0) each [smallmap,bigmap]];
2516
2517
2518/// Internal Function: _skin_distance_match(poly1,poly2)
2519/// Usage:
2520///   polys = _skin_distance_match(poly1,poly2);
2521/// Description:
2522///   Find a way of associating the vertices of poly1 and vertices of poly2
2523///   that minimizes the sum of the length of the edges that connect the two polygons.
2524///   Polygons can be in 2d or 3d.  The algorithm has cubic run time, so it can be
2525///   slow if you pass large polygons.  The output is a pair of polygons with vertices
2526///   duplicated as appropriate to be used as input to `skin()`.
2527/// Arguments:
2528///   poly1 = first polygon to match
2529///   poly2 = second polygon to match
2530function _skin_distance_match(poly1,poly2) =
2531   let(
2532      swap = len(poly1)>len(poly2),
2533      big = swap ? poly1 : poly2,
2534      small = swap ? poly2 : poly1,
2535      map_poly = [ for(
2536              i=0,
2537              bestcost = 1/0,
2538              bestmap = -1,
2539              bestpoly = -1
2540              ;
2541              i<=len(big)
2542              ;
2543              shifted = list_rotate(big,i),
2544              result =_dp_distance_array(small, shifted, abort_thresh = bestcost),
2545              bestmap = result[0]<bestcost ? result[1] : bestmap,
2546              bestpoly = result[0]<bestcost ? shifted : bestpoly,
2547              best_i = result[0]<bestcost ? i : best_i,
2548              bestcost = min(result[0], bestcost),
2549              i=i+1
2550              )
2551              if (i==len(big)) each [bestmap,bestpoly,best_i]],
2552      map = _dp_extract_map(map_poly[0]),
2553      smallmap = map[0],
2554      bigmap = map[1],
2555      // These shifts are needed to handle the case when points from both ends of one curve map to a single point on the other
2556      bigshift =  len(bigmap) - max(max_index(bigmap,all=true))-1,
2557      smallshift = len(smallmap) - max(max_index(smallmap,all=true))-1,
2558      newsmall = list_rotate(repeat_entries(small,unique_count(smallmap)[1]),smallshift),
2559      newbig = list_rotate(repeat_entries(map_poly[1],unique_count(bigmap)[1]),bigshift)
2560      )
2561      swap ? [newbig, newsmall] : [newsmall,newbig];
2562
2563
2564// This function associates vertices but with the assumption that index 0 is associated between the
2565// two inputs.  This gives only quadratic run time.  As above, output is pair of polygons with
2566// vertices duplicated as suited to use as input to skin().
2567
2568function _skin_aligned_distance_match(poly1, poly2) =
2569    let(
2570      result = _dp_distance_array(poly1, poly2, abort_thresh=1/0),
2571      map = _dp_extract_map(result[1]),
2572      shift0 = len(map[0]) - max(max_index(map[0],all=true))-1,
2573      shift1 = len(map[1]) - max(max_index(map[1],all=true))-1,
2574      new0 = list_rotate(repeat_entries(poly1,unique_count(map[0])[1]),shift0),
2575      new1 = list_rotate(repeat_entries(poly2,unique_count(map[1])[1]),shift1)
2576  )
2577  [new0,new1];
2578
2579
2580//////////////////////////////////////////////////////////////////////////////////////////////////////////////
2581/// Internal Function: _skin_tangent_match()
2582/// Usage:
2583///   x = _skin_tangent_match(poly1, poly2)
2584/// Description:
2585///   Finds a mapping of the vertices of the larger polygon onto the smaller one.  Whichever input is the
2586///   shorter path is the polygon, and the longer input is the curve.  For every edge of the polygon, the algorithm seeks a plane that contains that
2587///   edge and is tangent to the curve.  There will be more than one such point.  To choose one, the algorithm centers the polygon and curve on their centroids
2588///   and chooses the closer tangent point.  The algorithm works its way around the polygon, computing a series of tangent points and then maps all of the
2589///   points on the curve between two tangent points into one vertex of the polygon.  This algorithm can fail if the curve has too few points or if it is concave.
2590/// Arguments:
2591///   poly1 = input polygon
2592///   poly2 = input polygon
2593function _skin_tangent_match(poly1, poly2) =
2594    let(
2595        swap = len(poly1)>len(poly2),
2596        big = swap ? poly1 : poly2,
2597        small = swap ? poly2 : poly1,
2598        curve_offset = centroid(small)-centroid(big),
2599        cutpts = [for(i=[0:len(small)-1]) _find_one_tangent(big, select(small,i,i+1),curve_offset=curve_offset)],
2600        shift = last(cutpts)+1,
2601        newbig = list_rotate(big, shift),
2602        repeat_counts = [for(i=[0:len(small)-1]) posmod(cutpts[i]-select(cutpts,i-1),len(big))],
2603        newsmall = repeat_entries(small,repeat_counts)
2604    )
2605    assert(len(newsmall)==len(newbig), "Tangent alignment failed, probably because of insufficient points or a concave curve")
2606    swap ? [newbig, newsmall] : [newsmall, newbig];
2607
2608
2609function _find_one_tangent(curve, edge, curve_offset=[0,0,0], closed=true) =
2610    let(
2611        angles = [
2612            for (i = [0:len(curve)-(closed?1:2)])
2613            let(
2614                plane = plane3pt( edge[0], edge[1], curve[i]),
2615                tangent = [curve[i], select(curve,i+1)]
2616            ) plane_line_angle(plane,tangent)
2617        ],
2618        zero_cross = [
2619            for (i = [0:len(curve)-(closed?1:2)])
2620            if (sign(angles[i]) != sign(select(angles,i+1)))
2621            i
2622        ],
2623        d = [
2624            for (i = zero_cross)
2625            point_line_distance(curve[i]+curve_offset, edge)
2626        ]
2627    ) zero_cross[min_index(d)];
2628
2629
2630// Function: associate_vertices()
2631// Synopsis: Create vertex association to control how {{skin()}} links vertices. 
2632// SynTags: PathList
2633// Topics: Extrusion, Skinning, Paths
2634// See Also: skin()
2635// Usage:
2636//   newpoly = associate_vertices(polygons, split);
2637// Description:
2638//   Takes as input a list of polygons and duplicates specified vertices in each polygon in the list through the series so
2639//   that the input can be passed to `skin()`.  This allows you to decide how the vertices are linked up rather than accepting
2640//   the automatically computed minimal distance linkage.  However, the number of vertices in the polygons must not decrease in the list.
2641//   The output is a list of polygons that all have the same number of vertices with some duplicates.  You specify the vertex splitting
2642//   using the `split` which is a list where each entry corresponds to a polygon: split[i] is a value or list specifying which vertices in polygon i to split.
2643//   Give the empty list if you don't want a split for a particular polygon.  If you list a vertex once then it will be split and mapped to
2644//   two vertices in the next polygon.  If you list it N times then N copies will be created to map to N+1 vertices in the next polygon.
2645//   You must ensure that each mapping produces the correct number of vertices to exactly map onto every vertex of the next polygon.
2646//   Note that if you split (only) vertex i of a polygon that means it will map to vertices i and i+1 of the next polygon.  Vertex 0 will always
2647//   map to vertex 0 and the last vertices will always map to each other, so if you want something different than that you'll need to reindex
2648//   your polygons.
2649// Arguments:
2650//   polygons = list of polygons to split
2651//   split = list of lists of split vertices
2652// Example(FlatSpin,VPD=17,VPT=[0,0,2]):  If you skin together a square and hexagon using the optimal distance method you get two triangular faces on opposite sides:
2653//   sq = regular_ngon(4,side=2);
2654//   hex = apply(rot(15),hexagon(side=2));
2655//   skin([sq,hex], slices=10, refine=10, method="distance", z=[0,4]);
2656// Example(FlatSpin,VPD=17,VPT=[0,0,2]):  Using associate_vertices you can change the location of the triangular faces.  Here they are connect to two adjacent vertices of the square:
2657//   sq = regular_ngon(4,side=2);
2658//   hex = apply(rot(15),hexagon(side=2));
2659//   skin(associate_vertices([sq,hex],[[1,2]]), slices=10, refine=10, sampling="segment", z=[0,4]);
2660// Example(FlatSpin,VPD=17,VPT=[0,0,2]): Here the two triangular faces connect to a single vertex on the square.  Note that we had to rotate the hexagon to line them up because the vertices match counting forward, so in this case vertex 0 of the square matches to vertices 0, 1, and 2 of the hexagon.
2661//   sq = regular_ngon(4,side=2);
2662//   hex = apply(rot(60),hexagon(side=2));
2663//   skin(associate_vertices([sq,hex],[[0,0]]), slices=10, refine=10, sampling="segment", z=[0,4]);
2664// Example(3D): This example shows several polygons, with only a single vertex split at each step:
2665//   sq = regular_ngon(4,side=2);
2666//   pent = pentagon(side=2);
2667//   hex = hexagon(side=2);
2668//   sep = regular_ngon(7,side=2);
2669//   profiles = associate_vertices([sq,pent,hex,sep], [1,3,4]);
2670//   skin(profiles ,slices=10, refine=10, method="distance", z=[0,2,4,6]);
2671// Example(3D): The polygons cannot shrink, so if you want to have decreasing polygons you'll need to concatenate multiple results.  Note that it is perfectly ok to duplicate a profile as shown here, where the pentagon is duplicated:
2672//   sq = regular_ngon(4,side=2);
2673//   pent = pentagon(side=2);
2674//   grow = associate_vertices([sq,pent], [1]);
2675//   shrink = associate_vertices([sq,pent], [2]);
2676//   skin(concat(grow, reverse(shrink)), slices=10, refine=10, method="distance", z=[0,2,2,4]);
2677function associate_vertices(polygons, split, curpoly=0) =
2678   curpoly==len(polygons)-1 ? polygons :
2679   let(
2680      polylen = len(polygons[curpoly]),
2681      cursplit = force_list(split[curpoly])
2682   )
2683    assert(len(split)==len(polygons)-1,str(split,"Split list length mismatch: it has length ", len(split)," but must have length ",len(polygons)-1))
2684    assert(polylen<=len(polygons[curpoly+1]),str("Polygon ",curpoly," has more vertices than the next one."))
2685    assert(len(cursplit)+polylen == len(polygons[curpoly+1]),
2686           str("Polygon ", curpoly, " has ", polylen, " vertices.  Next polygon has ", len(polygons[curpoly+1]),
2687                  " vertices.  Split list has length ", len(cursplit), " but must have length ", len(polygons[curpoly+1])-polylen))
2688    assert(len(cursplit) == 0 || max(cursplit)<polylen && min(curpoly)>=0,
2689           str("Split ",cursplit," at polygon ",curpoly," has invalid vertices.  Must be in [0:",polylen-1,"]"))
2690    len(cursplit)==0 ? associate_vertices(polygons,split,curpoly+1) :
2691    let(
2692      splitindex = sort(concat(count(polylen), cursplit)),
2693      newpoly = [for(i=[0:len(polygons)-1]) i<=curpoly ? select(polygons[i],splitindex) : polygons[i]]
2694    )
2695   associate_vertices(newpoly, split, curpoly+1);
2696
2697
2698
2699// DefineHeader(Table;Headers=Texture Name|Type|Description): Texture Values
2700
2701// Section: Texturing
2702//   Some operations are able to add texture to the objects they create.  A texture can be any regularly repeated variation in the height of the surface.
2703//   To define a texture you need to specify how the height should vary over a rectangular block that will be repeated to tile the object.  Because textures
2704//   are based on rectangular tiling, this means adding textures to curved shapes may result in distortion of the basic texture unit.  For example, if you
2705//   texture a cone, the scale of the texture will be larger at the wide end of the cone and smaller at the narrower end of the cone.
2706//   .
2707//   You can specify a texture using to method: a height field or a VNF.  For each method you also must specify the scale of the texture, which
2708//   gives the size of the rectangular unit in your object that will correspond to one texture tile.  Note that this scale does not preserve
2709//   aspect ratio: you can stretch the texture as desired.  
2710// Subsection: Height Field Texture Maps
2711//   The simplest way to specify a texture map is to give a 2d array of
2712//   height values which specify the height of the texture on a grid.
2713//   Values in the height field should range from 0 to 1.  A zero height
2714//   in the height field corresponds to the height of the surface and 1
2715//   the heighest point in the texture.
2716// Figure(2D,Big,NoScales): Here is a 2d texture described by a "grid" that just contains a single row.  Such a texture can be used to create ribbing. The texture is `[[0, 1, 1, 0]]`, and the fixture shows three repetitions of the basic texture unit.
2717//   ftex1 = [0,1,1,0,0];
2718//   stroke( transpose([count(5),ftex1]), dots=true, dots_width=3,width=.05);
2719//   right(4)stroke( transpose([count(5),ftex1]), dots=true, width=.05,dots_color="red",color="blue",dots_width=3);
2720//   right(8)stroke( transpose([count(5),ftex1]), dots=true, dots_width=3,width=.05);
2721//   stroke([[4,-.3],[8,-.3]],width=.05,endcaps="arrow2",color="black");
2722//   move([6,-.4])color("black")text("Texture Size", size=0.3,anchor=BACK);
2723// Continues:
2724//   Line segments connect the dots within the texture and also the dots between adjacent texture tiles.
2725//   The size of the texture (specified with `tex_size`) includes the segment that connects the tile to the next one.
2726//   Note that the grid is always uniformly spaced.  If you want to keep the texture the same size but make the slope
2727//   steeper you need to add more points.  
2728// Figure(2D,Big,NoScales):  
2729//   ftex2 = xscale(4/11,transpose([count(12),[0,1,1,1,1,1,1,1,1,1,0,0]]));
2730//   stroke( ftex2, dots=true, dots_width=3,width=.05);
2731//   right(4)stroke( ftex2, dots=true, width=.05,dots_color="red",color="blue",dots_width=3);
2732//   right(8)stroke( ftex2, dots=true, dots_width=3,width=.05);
2733//   stroke([[4,-.3],[8,-.3]],width=.05,endcaps="arrow2",color="black");
2734//   move([6,-.4])color("black")text("Texture Size", size=0.3,anchor=BACK);
2735// Continues:
2736//   A more serious limitation of height field textures is that some shapes, such as hexagons or circles, cannot be accurately represented because
2737//   their points don't fall on a grid.  Trying to create such shapes is difficult and will require many points to approximate the
2738//   true point positions for the desired shape.  This will make the texture slow to compute.  
2739//   Another serious limitation is more subtle.  In the 2D examples above, it is obvious how to connect the
2740//   dots together.  But in 3D example we need to triangulate the points on a grid, and this triangulation is not unique.
2741//   The `style` argument lets you specify how the points are triangulated using the styles supported by {{vnf_vertex_array()}}.
2742//   In the example below we have expanded the 2D example into 3D:
2743//   ```openscad
2744//       [[0,0,0,0],
2745//        [0,1,1,0],
2746//        [0,1,1,0],
2747//        [0,0,0,0]]
2748//   ```
2749//   and we show the 3D triangulations produced by the different styles:
2750// Figure(3D,Big,NoAxes,VPR=[39.2,0,13.3],VPT=[3.76242,-5.50969,4.51854],VPD=32.0275):
2751//   tex = [
2752//          [0,0,0,0,0],
2753//          [0,1,1,0,0],
2754//          [0,1,1,0,0],
2755//          [0,0,0,0,0],
2756//          [0,0,0,0,0]       
2757//         ];
2758//   hm = [for(i=[0:4]) [for(j=[0:4]) [i,-j,tex[i][j]]]];      
2759//   types = ["quincunx", "convex", "concave","default","alt","min_edge"]; 
2760//   grid2d(spacing=5, n=[3,2]){
2761//     let(s = types[$row*3+$col]){
2762//       vnf_polyhedron(vnf_vertex_array(hm,style=s));
2763//       if ($row==1)
2764//         back(.8)right(2)rotate($vpr)color("black")text(s,size=.5,anchor=CENTER);
2765//       else
2766//         fwd(4.7)right(2)rotate($vpr)color("black")text(s,size=.5,anchor=CENTER);    
2767//     }
2768//   }  
2769// Continues:
2770//   Note that of the six available styles, five produce a different result.  There may exist some concave shape where none of the styles
2771//   produce the right result everywhere on the shape.  If this happens it would be another limitation of height field textures.  (If you have an
2772//   example of such a texture and shape please let us know!)
2773// Subsection: VNF Textures
2774//   VNF textures overcome all of the limitations of height field textures, but with two costs.  They can be more difficult to construct than
2775//   a simple array of height values, and they are significantly slower to compute for a tile with the same number of points.  Note, however, for
2776//   textures that don't neatly lie on a grid, a VNF tile will be more efficient than a finely sampled height field.  With VNF textures you can create
2777//   textures that have disconnected components, or concavities that cannot be expressed with a single valued height map.  However, you can also
2778//   create invalid textures that fail to close at the ends, so care is required to ensure that your resulting shape is valid.  
2779//   .
2780//   A VNF texture is defined by defining the texture tile with a VNF whose projection onto the XY plane is contained in the unit square [0,1] x [0,1] so
2781//   that the VNF can be tiled.   The VNF is tiled without a gap, matching the edges, so the vertices along corresponding edges must match to make a
2782//   consistent triangulation possible.  The VNF cannot have any X or Y values outside the interval [0,1].  If you want a valid polyhedron
2783//   that OpenSCAD will render then you need to take care with edges of the tiles that correspond to endcap faces in the textured object.
2784//   So for example, in a linear sweep, the top and bottom edges of tiles end abruptly to form the end cap of the object.  You can make a valid object
2785//   in two ways.  One way is to create a tile with a single, complete edge along Y=0, and of course a corresponding edges along Y=1.  The second way
2786//   to make a valid object is to have no points at all on the Y=0 line, and of course none on Y=1.  In this case, the resulting texture produces
2787//   a collection of disconnected objects.  Note that the Z coordinates of your tile can be anything, but for the dimensional settings on textures
2788//   to work intuitively, you should construct your tile so that Z ranges from 0 to 1.
2789// Figure(3D): This is the "hexgrid" VNF tile, which creates a hexagonal grid texture, something which doesn't work well with a height field because the edges of the hexagon don't align with the grid.  Note how the tile ranges between 0 and 1 in both X, Y and Z.  In fact, to get a proper aspect ratio in your final texture you need to use the `tex_size` parameter to introduct a sqrt(3) scale factor.  
2790//   tex = texture("hex_grid");
2791//   vnf_polyhedron(tex);
2792// Figure(3D): This is an example of a tile that has no edges at the top or bottom, so it creates disconnected rings.  See {{linear_sweep()}} for examples showing this tile in use.
2793//   shape = skin([
2794//                 rect(2/5),
2795//                 rect(2/3),
2796//                 rect(2/5)
2797//                ],
2798//                z=[0,1/2,1],
2799//                slices=0,
2800//                caps=false);
2801//   tile = move([0,1/2,2/3],yrot(90,shape));
2802//   vnf_polyhedron(tile);
2803
2804
2805
2806// Function: texture()
2807// Topics: Textures, Knurling
2808// Synopsis: Produce a standard texture. 
2809// Topics: Extrusion, Textures
2810// See Also: linear_sweep(), rotate_sweep(), heightfield(), cylindrical_heightfield()
2811// Usage:
2812//   tx = texture(tex, [n=], [inset=], [gap=], [roughness=]);
2813// Description:
2814//   Given a texture name, returns a texture.  Textures can come in two varieties:
2815//   - Heightfield textures which are 2D arrays of scalars.  These are usually faster to render, but can be less precise and prone to triangulation errors.  The table below gives the recommended style for the best triangulation.  If results are still incorrect, switch to the similar VNF tile by adding the "_vnf" suffix.
2816//   - VNF Tile textures, which are VNFs that cover the unit square [0,0] x [1,1].  These tend to be slower to render, but allow greater flexibility and precision for shapes that don't align with a grid.
2817//   In the descriptions below, imagine the textures positioned on the XY plane, so "horizontal" refers to the "sideways" dimensions of the texture and
2818//   "up" and "down" refer to the depth dimension.  If a texture is placed on a cylinder the "depth" will become the radial direction and the "horizontal"
2819//   direction will be the vertical and tangential directions on the cylindrical surface.  All horizontal dimensions for VNF textures are relative to the unit square
2820//   on which the textures are defined, so a value of 0.25 for a gap or inset will refer to 1/4 of the texture's full length and/or width.  All supported textures appear below in the examples.  
2821// Arguments:
2822//   tex = The name of the texture to get.
2823//   ---
2824//   n = The general number of vertices to use to refine the resolution of the texture.
2825//   inset = The amount to inset part of a VNF tile texture.  Generally between 0 and 0.5.
2826//   gap = The gap between logically distinct parts of a VNF tile.  (ie: gap between bricks, gap between truncated ribs, etc.)
2827//   roughness = The amount of roughness used on the surface of some heightfield textures.  Generally between 0 and 0.5.
2828// Example(3D): **"bricks"** (Heightfield) = A brick-wall pattern.  Giving `n=` sets the number of heightfield samples to `n x n`.  Default: 24.  Giving `roughness=` adds a level of height randomization to add roughness to the texture.  Default: 0.05.  Use `style="convex"`.
2829//   tex = texture("bricks");
2830//   linear_sweep(
2831//       rect(30), texture=tex, h=30,
2832//       tex_size=[10,10]
2833//   );
2834// Example(3D): **"bricks_vnf"** (VNF) = VNF version of "bricks".  Giving `gap=` sets the "mortar" gap between adjacent bricks, default 0.05.  Giving `inset=` specifies that the top face of the brick is smaller than the bottom of the brick by `inset` on each of the four sides.  If `gap` is zero then an `inset` value close to 0.5 will cause bricks to come to a sharp pointed edge, with just a tiny flat top surface.  Note that `gap+inset` must be strictly smaller than 0.5.   Default is `inset=0.05`.  
2835//   tex = texture("bricks_vnf");
2836//   linear_sweep(
2837//       rect(30), texture=tex, h=30,
2838//       tex_size=[10,10]
2839//   );
2840// Example(3D): "bricks_vnf" texture with large inset. 
2841//   tex = texture("bricks_vnf",inset=0.25);
2842//   linear_sweep(
2843//       rect(30), texture=tex, h=30,
2844//       tex_size=[10,10]
2845//   );
2846// Example(3D,VPR=[84.4,0,4.7],VPT=[2.44496,6.53317,14.6135],VPD = 126): **"checkers"** (VNF) = A pattern of alternating checkerboard squares.  Giving `inset=` specifies that the top face of the checker surface is smaller than the bottom by `inset` on each of the four sides.  As `inset` approaches 0.5 the tops come to sharp corners.  You must set `inset` strictly between 0 and 0.5.  Default: 0.05.
2847//   tex = texture("checkers");
2848//   linear_sweep(
2849//       rect(30), texture=tex, h=30,
2850//       tex_size=[10,10]
2851//   );
2852// Example(3D,VPR=[84.4,0,4.7],VPT=[2.44496,6.53317,14.6135],VPD = 126): "checkers" texture with large inset.  
2853//   tex = texture("checkers",inset=0.25);
2854//   linear_sweep(
2855//       rect(30), texture=tex, h=30,
2856//       tex_size=[10,10]
2857//   );
2858// Example(3D): **"cones"** (VNF) = Raised conical spikes.  Giving `n=` sets `$fn` for the cone (will be rounded to a multiple of 4).  Default: 16.  Giving `inset=` specifies the horizontal inset of the base of the cone, relative to the tile edges.  The `inset` value must be nonnegative and smaller than 0.5.  Default: 0.
2859//   tex = texture("cones");
2860//   linear_sweep(
2861//       rect(30), texture=tex, h=30, tex_scale=3,
2862//       tex_size=[10,10]
2863//   );
2864// Example(3D): **"cubes"** (VNF) = Corner-cubes texture.  This texture needs to be scaled in vertically by sqrt(3) to have its correct aspect
2865//   tex = texture("cubes");
2866//   linear_sweep(
2867//       rect(30), texture=tex, h=30,
2868//       tex_size=[10,10]
2869//   );
2870// Example(3D): "cubes" texture at the correct scale.  
2871//   tex = texture("cubes");
2872//   linear_sweep(
2873//       rect(30), texture=tex, h=20*sqrt(3), tex_scale=3,
2874//       tex_size=[10,10*sqrt(3)]
2875//   );
2876// Example(3D): **"diamonds"** (Heightfield) = Four-sided pyramid with the corners of the base aligned with the axes.  Compare to "pyramids".  Useful for knurling.  Giving `n=` sets the number of heightfield samples to `n x n`. Default: 2.  Use `style="concave"` for pointed bumps, or `style="default"` or `style="alt"` for a diagonal ribs.  
2877//   tex = texture("diamonds");
2878//   linear_sweep(
2879//       rect(30), texture=tex, h=30,
2880//       tex_size=[10,10], style="concave"
2881//   );
2882// Example(3D): "diamonds" texture can give diagonal ribbing with "default" style. 
2883//   tex = texture("diamonds");
2884//   linear_sweep(
2885//       rect(30), texture=tex, h=30,
2886//       tex_size=[10,10], style="default"
2887//   );
2888// Example(3D): "diamonds" texture gives diagonal ribbing the other direction with "alt" style.  
2889//   tex = texture("diamonds");
2890//   linear_sweep(
2891//       rect(30), texture=tex, h=30,
2892//       tex_size=[10,10], style="alt"
2893//   );
2894// Example(3D): **"diamonds_vnf"** (VNF) = VNF version of "diamonds".
2895//   tex = texture("diamonds_vnf");
2896//   linear_sweep(
2897//       rect(30), texture=tex, h=30,
2898//       tex_size=[10,10]
2899//   );
2900// Example(3D): **"dimples"** (VNF) = Round divots.  Giving `n=` sets `$fn` for the curve (will be rounded to a multiple of 4).  Default: 16.  Giving `inset=` specifies the horizontal distance of the flat region around the dimple relative to the edge of the tile.  Must be nonnegative and strictly less than 0.5.  Default: 0.05.  
2901//   tex = texture("dimples");
2902//   linear_sweep(
2903//       rect(30), texture=tex, h=30, 
2904//       tex_size=[10,10]
2905//   );
2906// Example(3D): **"dots"** (VNF) = Raised round bumps.  Giving `n=` sets `$fn` for the curve (will be rounded to a multiple of 4).  Default: 16.   Giving `inset=` specifies the horizontal inset of the dots, relative to the edge of the tile.  Must be nonnegative and strictly less than 0.5.  Default: 0.05.
2907//   tex = texture("dots");
2908//   linear_sweep(
2909//       rect(30), texture=tex, h=30,
2910//       tex_size=[10,10]
2911//   );
2912// Example(3D): **"hex_grid"** (VNF) = A hexagonal grid defined by V-grove borders.  Giving `inset=` specifies the horizontal inset of the left and right edges of the hexagonal tops, relative to their bottoms.  This means the V-groove top width for grooves running parallel to the Y axis will be double the inset value.  If the texture is scaled in the Y direction by sqrt(3) then the groove will be uniform on all six sides of the hexagon.  Inset must be strictly between 0 and 0.5, default: 0.1.
2913//   tex = texture("hex_grid");
2914//   linear_sweep(
2915//       rect(30), texture=tex, h=30,
2916//       tex_size=[10,10]
2917//   );
2918// Example(3D): "hex_grid" texture with large inset
2919//   tex = texture("hex_grid", inset=0.4);
2920//   linear_sweep(
2921//       rect(30), texture=tex, h=30,
2922//       tex_size=[10,10]
2923//   );
2924// Example(3D): "hex_grid" scaled in Y by sqrt(3) so hexagons are regular and grooves are all the same width.  Note height of cube is also scaled so tile fits without being automatically adjusted to fit, ruining our choice of scale.
2925//   tex = texture("hex_grid",inset=.07);
2926//   linear_sweep(
2927//       rect(30), texture=tex, h=quantup(30,10*sqrt(3)),
2928//       tex_size=[10,10*sqrt(3)], tex_scale=3
2929//   );
2930// Example(3D): "hex_grid" texture, with approximate scaling because 17 is close to sqrt(3) times 10.
2931//   tex = texture("hex_grid");
2932//   linear_sweep(
2933//       rect(30), texture=tex, h=34,
2934//       tex_size=[10,17]
2935//   );
2936// Example(3D): **"hills"** (Heightfield) = Wavy sine-wave hills and valleys,  Giving `n=` sets the number of heightfield samples to `n` x `n`.  Default: 12.  Set `style="quincunx"`.
2937//   tex = texture("hills");
2938//   linear_sweep(
2939//       rect(30), texture=tex, h=30,
2940//       tex_size=[10,10], style="quincunx"
2941//   );
2942// Example(3D): **"pyramids"** (Heightfield) = Four-sided pyramid with the edges of the base aligned with the axess.  Compare to "diamonds". Useful for knurling.  Giving `n=` sets the number of heightfield samples to `n` by `n`. Default: 2. Set style to "convex".  Note that style="concave" or style="min_edge" produce mini-diamonds with flat squares in between.
2943//   tex = texture("pyramids");
2944//   linear_sweep(
2945//       rect(30), texture=tex, h=30,
2946//       tex_size=[10,10], style="convex"
2947//   );
2948// Example(3D): "pyramids" texture, with "concave" produces a mini-diamond texture.  Note that "min_edge" also gives this result.
2949//   tex = texture("pyramids");
2950//   linear_sweep(
2951//       rect(30), texture=tex, h=30,
2952//       tex_size=[10,10], style="concave"
2953//   );
2954// Example(3D): **"pyramids_vnf"** (VNF) = VNF version of "pyramids".
2955//   tex = texture("pyramids_vnf");
2956//   linear_sweep(
2957//       rect(30), texture=tex, h=30,
2958//       tex_size=[10,10]
2959//   );
2960// Example(3D): **"ribs"** (Heightfield) = Vertically aligned triangular ribs.  Giving `n=` sets the number of heightfield samples to `n` by 1.  Default: 2.  The choice of style does not matter.
2961//   tex = texture("ribs");
2962//   linear_sweep(
2963//       rect(30), texture=tex, h=30, tex_scale=3,
2964//       tex_size=[10,10], style="concave"
2965//   );
2966// Example(3D): **"rough"** (Heightfield) = A pseudo-randomized rough texture.  Giving `n=` sets the number of heightfield samples to `n` by `n`.  Default: 32.  The `roughness=` parameter specifies the height of the random texture.  Default: 0.2.
2967//   tex = texture("rough");
2968//   linear_sweep(
2969//       rect(30), texture=tex, h=30,
2970//       tex_size=[10,10], style="min_edge"
2971//   );
2972// Example(3D): **"tri_grid"** (VNF) = A triangular grid defined by V-groove borders  Giving `inset=` specifies the horizontal inset of the triangular tops, relative to their bottoms, along the horizontal edges (parallel to the X axis) of the triangles.  This means the V-groove top width of the grooves parallel to the X axis will be double the inset value.  (The other grooves are wider.) If the tile is scaled in the Y direction by sqrt(3) then the groove will be uniform on the three sides of the triangle.  The inset must be strictly between 0 and 1/6, default: 0.05.
2973//   tex = texture("tri_grid");
2974//   linear_sweep(
2975//       rect(30), texture=tex, h=30,
2976//       tex_size=[10,10]
2977//   );
2978// Example(3D): "tri_grid" texture with large inset.  (Max inset for tri_grid is 1/6.)  
2979//   tex = texture("tri_grid",inset=.12);
2980//   linear_sweep(
2981//       rect(30), texture=tex, h=30,
2982//       tex_size=[10,10]
2983//   );
2984// Example(3D): "tri_grid" texture scaled in Y by sqrt(3) so triangles are equilateral and grooves are all the same width.  Note we have to ensure the height evenly fits the scaled texture tiles.
2985//   tex = texture("tri_grid",inset=.04);
2986//   linear_sweep(
2987//       rect(30), texture=tex, h=quantup(30,10*sqrt(3)),
2988//       tex_size=[10,10*sqrt(3)], tex_scale=3
2989//   );
2990// Example(3D): "tri_grid" texture.  Here scale makes Y approximately sqrt(3) larger than X so triangles are close to equilateral.
2991//   tex = texture("tri_grid");
2992//   linear_sweep(
2993//       rect(30), texture=tex, h=34,
2994//       tex_size=[10,17]
2995//   );
2996// Example(3D): **"trunc_diamonds"** (VNF) = Truncated diamonds, four-sided pyramids with the base corners aligned with the axes and the top cut off.  Or you can interpret it as V-groove lines at 45º angles.  Giving `inset=` specifies the horizontal inset of the square top face compared to the bottom face along all four edges.  This means the V-groove top width will be double the inset value.  The inset must be strictly between 0 and sqrt(2)/4, which is about 0.35.  Default: 0.1.
2997//   tex = texture("trunc_diamonds");
2998//   linear_sweep(
2999//       rect(30), texture=tex, h=30,
3000//       tex_size=[10,10]
3001//   );
3002// Example(3D): "trunc_diamonds" texture with large inset. 
3003//   tex = texture("trunc_diamonds",inset=.25);
3004//   linear_sweep(
3005//       rect(30), texture=tex, h=30,
3006//       tex_size=[10,10]
3007//   );
3008// Example(3D): **"trunc_pyramids"** (Heightfield) = Truncated pyramids, four sided pyramids with the base edges aligned to the axes and the top cut off.  Giving `n=` sets the number of heightfield samples to `n` by `n`.  Default: 6.  Set `style="convex"`.
3009//   tex = texture("trunc_pyramids");
3010//   linear_sweep(
3011//       rect(30), texture=tex, h=30,
3012//       tex_size=[10,10], style="convex"
3013//   );
3014// Example(3D): **"trunc_pyramids_vnf"** (VNF) = Truncated pyramids, four sided pyramids with the base edges aligned to the axes and the top cut off.  You can also regard this as a grid of V-grooves.  Giving `inset=` specifies the horizontal inset of the flat square tops on all four sides relative to their bottoms.  This means the V-groove top width will be double the inset value.  The inset must be strictly between 0 and 0.5.  Default: 0.1.
3015//   tex = texture("trunc_pyramids_vnf");
3016//   linear_sweep(
3017//       rect(30), texture=tex, h=30,
3018//       tex_size=[10,10]
3019//   );
3020// Example(3D): "trunc_pyramids_vnf" texture with large inset
3021//   tex = texture("trunc_pyramids_vnf", inset=.4);
3022//   linear_sweep(
3023//       rect(30), texture=tex, h=30,
3024//       tex_size=[10,10]
3025//   );
3026// Example(3D): **"trunc_ribs"** (Heightfield) = Truncated ribs.  Vertically aligned triangular ribs with the tops cut off, and with rib separation equal to the width of the flat tops.  Giving `n=` sets the number of heightfield samples to `n` by `1`.  Default: 4.  The style does not matter.
3027//   tex = texture("trunc_ribs");
3028//   linear_sweep(
3029//       rect(30), h=30, texture=tex,
3030//       tex_scale=3, tex_size=[10,10],
3031//       style="concave"
3032//   );
3033// Example(3D): **"trunc_ribs_vnf"** (VNF) = Vertically aligned triangular ribs with the tops cut off.  Giving `gap=` sets the bottom gap between ribs.  Giving `inset=` specifies the horizontal inset of the rib top face, relative to the bottom on both sides.  In order to fit, gap+2*inset must be less than 1.  (This is because the gap is counted once but the inset counts on both sides.)  Defaults: gap=1/4, inset=1/4.
3034//   tex = texture("trunc_ribs_vnf", gap=0.25, inset=1/6);
3035//   linear_sweep(
3036//       rect(30), h=30, texture=tex,
3037//       tex_scale=3, tex_size=[10,10]
3038//   );
3039// Example(3D): **"wave_ribs"** (Heightfield) = Vertically aligned wavy ribs.  Giving `n=` sets the number of heightfield samples to `n` by `1`.  Default: 8.  The style does not matter.  
3040//   tex = texture("wave_ribs");
3041//   linear_sweep(
3042//       rect(30), h=30, texture=tex, 
3043//       tex_size=[10,10], tex_scale=3, style="concave"
3044//   );
3045
3046function texture(tex, n, inset, gap, roughness) =
3047    assert(is_undef(n) || all_positive([n]), "n must be a positive value if given")
3048    assert(is_undef(inset) || is_finite(inset), "inset must be a number if given")
3049    assert(is_undef(gap) || is_finite(gap), "gap must be a number if given")
3050    assert(is_undef(roughness) || all_nonnegative([roughness]), "roughness must be a nonnegative value if given")  
3051    tex=="ribs"?
3052        assert(num_defined([gap, inset, roughness])==0, "ribs texture does not accept gap, inset or roughness")
3053
3054        let(
3055            n = quantup(default(n,2),2)
3056        ) [[
3057            each lerpn(1,0,n/2,endpoint=false),
3058            each lerpn(0,1,n/2,endpoint=false),
3059        ]] :
3060    tex=="trunc_ribs"?
3061        assert(num_defined([gap, inset, roughness])==0, "trunc_ribs texture does not accept gap, inset or roughness")
3062        let(
3063            n = quantup(default(n,4),4)
3064        ) [[
3065            each repeat(0,n/4),
3066            each lerpn(0,1,n/4,endpoint=false),
3067            each repeat(1,n/4),
3068            each lerpn(1,0,n/4,endpoint=false),
3069        ]] :
3070    tex=="trunc_ribs_vnf"?
3071        assert(is_undef(n), "trunc_ribs_vnf texture does not accept n")
3072        let(
3073            inset = default(inset,1/4)*2,
3074            gap = default(gap,1/4)
3075        )
3076        assert(all_nonnegative([inset,gap]), "trunc_ribs_vnf texture requires gap>=0 and inset>=0")
3077        assert(gap+inset > 0, "trunc_ribs_vnf texture requires that gap+inset>0")
3078        assert(gap+inset <= 1, "trunc_ribs_vnf texture requires that gap+2*inset<=1")
3079        [
3080            [
3081               each move([0.5,0.5], p=path3d(rect([1-gap,1]),0)),
3082               each move([0.5,0.5], p=path3d(rect([1-gap-inset,1]),1)),
3083               each path3d(square(1)),
3084            ], [
3085                [1,2,6], [1,6,5], [0,4,3], [3,4,7],
3086                if (gap+inset < 1-EPSILON) each [[4,5,6], [4,6,7]],
3087                if (gap > EPSILON) each [[1,9,10], [1,10,2], [0,3,8], [3,11,8]],
3088            ]
3089        ] :
3090    tex=="wave_ribs"?
3091        assert(num_defined([gap, inset, roughness])==0, "wave_ribs texture does not accept gap, inset or roughness")  
3092        let(
3093            n = max(6,default(n,8))
3094        ) [[
3095            for(a=[0:360/n:360-EPSILON])
3096            (cos(a)+1)/2
3097        ]] :
3098    tex=="diamonds"?
3099        assert(num_defined([gap, inset, roughness])==0, "diamonds texture does not accept gap, inset or roughness")  
3100        let(
3101            n = quantup(default(n,2),2)
3102        ) [
3103            let(
3104                path = [
3105                    each lerpn(0,1,n/2,endpoint=false),
3106                    each lerpn(1,0,n/2,endpoint=false),
3107                ]
3108            )
3109            for (i=[0:1:n-1]) [
3110                for (j=[0:1:n-1]) min(
3111                    select(path,i+j),
3112                    select(path,i-j)
3113                )
3114            ],
3115        ] :
3116    tex=="diamonds_vnf"?
3117        assert(num_defined([n,gap, inset, roughness])==0, "diamonds_vnf texture does not accept n, gap, inset or roughness")
3118        [
3119            [
3120                [0,   1, 1], [1/2,   1, 0], [1,   1, 1],
3121                [0, 1/2, 0], [1/2, 1/2, 1], [1, 1/2, 0],
3122                [0,   0, 1], [1/2,   0, 0], [1,   0, 1],
3123            ], [
3124                [0,1,3], [2,5,1], [8,7,5], [6,3,7],
3125                [1,5,4], [5,7,4], [7,3,4], [4,3,1],
3126            ]
3127        ] :
3128    tex=="pyramids"?
3129        assert(num_defined([gap, inset, roughness])==0, "pyramids texture does not accept gap, inset or roughness")
3130        let(
3131            n = quantup(default(n,2),2)
3132        ) [
3133            for (i = [0:1:n-1]) [
3134                for (j = [0:1:n-1])
3135                1 - (max(abs(i-n/2), abs(j-n/2)) / (n/2))
3136            ]
3137        ] :
3138    tex=="pyramids_vnf"?
3139        assert(num_defined([n,gap, inset, roughness])==0, "pyramids_Vnf texture does not accept n, gap, inset or roughness")  
3140        [
3141            [ [0,1,0], [1,1,0], [1/2,1/2,1], [0,0,0], [1,0,0] ],
3142            [ [2,0,1], [2,1,4], [2,4,3], [2,3,0] ]
3143        ] :
3144    tex=="trunc_pyramids"?
3145        assert(num_defined([gap, inset, roughness])==0, "trunc_pyramids texture does not accept gap, inset or roughness")  
3146        let(
3147            n = quantup(default(n,6),3)
3148        ) [
3149            for (i = [0:1:n-1]) [
3150                for (j = [0:1:n-1])
3151                (1 - (max(n/6, abs(i-n/2), abs(j-n/2)) / (n/2))) * 1.5
3152            ]
3153        ] :
3154    tex=="trunc_pyramids_vnf"?
3155        assert(num_defined([gap, n, roughness])==0, "trunc_pyramids_vnf texture does not accept gap, n or roughness")
3156        let(
3157            inset = default(inset,0.1)
3158        )
3159        assert(inset>0 && inset<.5, "trunc_pyramids_vnf texture requires inset in (0,0.5)")
3160        [
3161            [
3162                each path3d(square(1)),
3163                each move([1/2,1/2,1], p=path3d(rect(1-2*inset))),
3164            ], [
3165                for (i=[0:3]) each [
3166                    [i, (i+1)%4, i+4],
3167                    [(i+1)%4, (i+1)%4+4, i+4],
3168                ],
3169                [4,5,6], [4,6,7],
3170            ]
3171        ] :
3172    tex=="hills"?
3173        assert(num_defined([gap, inset, roughness])==0, "hills texture does not accept gap, inset or roughness")  
3174        let(
3175            n = default(n,12)
3176        ) [
3177            for (a=[0:360/n:359.999]) [
3178                for (b=[0:360/n:359.999])
3179                (cos(a)*cos(b)+1)/2
3180            ]
3181        ] :
3182    tex=="bricks"?
3183        assert(num_defined([gap,inset])==0, "bricks texture does not accept gap or inset")  
3184        let(
3185            n = quantup(default(n,24),2),
3186            rough = default(roughness,0.05)
3187        ) [
3188            for (y = [0:1:n-1])
3189            rands(-rough/2, rough/2, n, seed=12345+y*678) + [
3190                for (x = [0:1:n-1])
3191                (y%(n/2) <= max(1,n/16))? 0 :
3192                let( even = floor(y/(n/2))%2? n/2 : 0 )
3193                (x+even) % n <= max(1,n/16)? 0 : 0.5
3194            ]
3195        ] :
3196    tex=="bricks_vnf"?
3197        assert(num_defined([n,roughness])==0, "bricks_vnf texture does not accept n or roughness")
3198        let(
3199            inset = default(inset,0.05),
3200            gap = default(gap,0.05)
3201        )
3202        assert(inset>=0,"bricks_vnf texture requires nonnegative inset")
3203        assert(gap>0, "bricks_vnf requires gap greater than 0")
3204        assert(gap+inset<0.5, "bricks_vnf requires gap+inset<0.5")
3205          [
3206            [
3207                each path3d(square(1)),
3208                each move([gap/2, gap/2, 0], p=path3d(square([1-gap, 0.5-gap]))),
3209                each move([gap/2+inset/2, gap/2+inset/2, 1], p=path3d(square([1-gap-inset, 0.5-gap-inset]))),
3210                each move([0, 0.5+gap/2, 0], p=path3d(square([0.5-gap/2, 0.5-gap]))),
3211                each move([0, 0.5+gap/2+inset/2, 1], p=path3d(square([0.5-gap/2-inset/2, 0.5-gap-inset]))),
3212                each move([0.5+gap/2, 0.5+gap/2, 0], p=path3d(square([0.5-gap/2, 0.5-gap]))),
3213                each move([0.5+gap/2+inset/2, 0.5+gap/2+inset/2, 1], p=path3d(square([0.5-gap/2-inset/2, 0.5-gap-inset]))),
3214            ], [
3215                [ 8, 9,10], [ 8,10,11], [16,17,18], [16,18,19], [24,25,26],
3216                [24,26,27], [ 0, 1, 5], [ 0, 5, 4], [ 1,13, 6], [ 1, 6, 5],
3217                [ 6,13,12], [ 6,12,21], [ 7,21,20], [ 6,21, 7], [ 0, 4, 7],
3218                [ 0, 7,20], [21,12,15], [21,15,22], [ 3,23,22], [ 3,22,15],
3219                [ 2,15,14], [ 2, 3,15], [23,27,26], [23,26,22], [21,22,26],
3220                [21,26,25], [21,25,24], [21,24,20], [12,16,19], [12,19,15],
3221                [14,15,19], [14,19,18], [13,17,16], [13,16,12], [ 6,10, 9],
3222                [ 6, 9, 5], [ 5, 9, 8], [ 5, 8, 4], [ 4, 8,11], [ 4,11, 7],
3223                [ 7,11,10], [ 7,10, 6],
3224            ]
3225        ] :
3226    tex=="checkers"?
3227        assert(num_defined([gap, n, roughness])==0, "checkers texture does not accept gap, n or roughness")
3228        let(
3229            inset = default(inset,0.05)
3230        )
3231        assert(inset>0 && inset<.5, "checkers texture requires inset in (0,0.5)")
3232          [
3233            [
3234                each move([0,0], p=path3d(square(0.5-inset),1)),
3235                each move([0,0.5], p=path3d(square(0.5-inset))),
3236                each move([0.5,0], p=path3d(square(0.5-inset))),
3237                each move([0.5,0.5], p=path3d(square(0.5-inset),1)),
3238                [1/2-inset/2,1/2-inset/2,1/2], [0,1,1], [1/2-inset,1,1],
3239                [1/2,1,0], [1-inset,1,0], [1,0,1], [1,1/2-inset,1],
3240                [1,1/2,0], [1,1-inset,0], [1,1,1], [1/2-inset/2,1-inset/2,1/2],
3241                [1-inset/2,1-inset/2,1/2], [1-inset/2,1/2-inset/2,1/2],
3242            ], [
3243                for (i=[0:4:12]) each [[i,i+1,i+2], [i, i+2, i+3]],
3244                [10,13,11], [13,12,11], [2,5,4], [4,3,2],
3245                [0,3,10], [10,9,0], [4,7,14], [4,14,13],
3246                [4,13,16], [10,16,13], [10,3,16], [3,4,16],
3247                [7,6,17], [7,17,18], [14,19,20], [14,20,15],
3248                [8,11,22], [8,22,21], [12,15,24], [12,24,23],
3249                [7,18,26], [7,26,14], [14,26,19], [18,19,26],
3250                [15,20,27], [20,25,27], [24,27,25], [15,27,24],
3251                [11,12,28], [12,23,28], [11,28,22], [23,22,28],
3252            ]
3253        ] :
3254    tex=="cones"?
3255        assert(num_defined([gap,roughness])==0, "cones texture does not accept gap or roughness")  
3256        let(
3257            n = quant(default(n,16),4),
3258            inset = default(inset,0)
3259        )
3260        assert(inset>=0 && inset<0.5)
3261        [
3262            [
3263                each move([1/2,1/2], p=path3d(circle(d=1-2*inset,$fn=n))),
3264                [1/2,1/2,1],
3265                each path3d(square(1)),
3266            ], [
3267                for (i=[0:1:n-1]) [i, (i+1)%n, n],
3268                for (i=[0:1:3], j=[0:1:n/4-1]) [n+1+i, (i*n/4+j+1)%n, i*n/4+j],
3269                if (inset > 0) for (i = [0:1:3]) [i+n+1, (i+1)%4+n+1, ((i+1)*n/4)%n],
3270            ]
3271        ] :
3272    tex=="cubes"?
3273        assert(num_defined([n, gap, inset, roughness])==0, "cubes texture does not accept n, gap, inset or roughness")  
3274        [
3275            [
3276                [0,1,1/2], [1,1,1/2], [1/2,5/6,1], [0,4/6,0], [1,4/6,0],
3277                [1/2,3/6,1/2], [0,2/6,1], [1,2/6,1], [1/2,1/6,0], [0,0,1/2],
3278                [1,0,1/2],
3279            ], [
3280                [0,1,2], [0,2,3], [1,4,2], [2,5,3], [2,4,5],
3281                [6,3,5], [4,7,5], [7,8,5], [6,5,8], [10,8,7],
3282                [9,6,8], [10,9,8],
3283            ]
3284        ] :
3285    tex=="trunc_diamonds"?
3286        assert(num_defined([gap, n, roughness])==0, "trunc_diamonds texture does not accept gap, n or roughness")
3287        let(
3288            inset = default(inset,0.1)/sqrt(2)*2
3289        )
3290        assert(inset>0 && inset<0.5)
3291        [
3292            [
3293                each move([1/2,1/2,0], p=path3d(circle(d=1,$fn=4))),
3294                each move([1/2,1/2,1], p=path3d(circle(d=1-inset*2,$fn=4))),
3295                for (a=[0:90:359]) each move([1/2,1/2], p=zrot(-a, p=[[1/2,inset,1], [inset,1/2,1], [1/2,1/2,1]]))
3296            ], [
3297                for (i=[0:3]) each let(j=i*3+8) [
3298                    [i,(i+1)%4,(i+1)%4+4], [i,(i+1)%4+4,i+4],
3299                    [j,j+1,j+2], [i, (i+3)%4, j], [(i+3)%4, j+1, j],
3300                ],
3301                [4,5,6], [4,6,7],
3302            ]
3303        ] :
3304    tex=="dimples" || tex=="dots" ?
3305        assert(num_defined([gap,roughness])==0, str(tex," texture does not accept gap or roughness"))
3306        let(
3307            n = quant(default(n,16),4),
3308            inset = default(inset,0.05)
3309        )
3310        assert(inset>=0 && inset < 0.5)
3311        let(
3312            rows=ceil(n/4),
3313            r=adj_ang_to_hyp(1/2-inset,45),
3314            dots = tex=="dots",
3315            cp = [1/2, 1/2, r*sin(45)*(dots?-1:1)],
3316            sc = 1 / (r - abs(cp.z)),
3317            uverts = [
3318                each path3d(square(1)),
3319                for (p=[0:1:rows-1], t=[0:360/n:359.999])
3320                    cp + (
3321                        dots? spherical_to_xyz(r, -t, 45-45*p/rows) :
3322                        spherical_to_xyz(r, -t, 135+45*p/rows)
3323                    ),
3324                cp + r * (dots?UP:DOWN),
3325            ],
3326            verts = zscale(sc, p=uverts),
3327            faces = [
3328                for (i=[0:1:3], j=[0:1:n/4-1]) [i, 4+(i*n/4+j+1)%n, 4+i*n/4+j],
3329                for (i=[0:1:rows-2], j=[0:1:n-1]) each [
3330                    [4+i*n+j, 4+(i+1)*n+(j+1)%n, 4+(i+1)*n+j],
3331                    [4+i*n+j, 4+i*n+(j+1)%n, 4+(i+1)*n+(j+1)%n],
3332                ],
3333                for (i=[0:1:n-1]) [4+(rows-1)*n+i, 4+(rows-1)*n+(i+1)%n, 4+rows*n],
3334                if (inset>0) for (i=[0:3]) [i, (i+1)%4, 4+(i+1)%4*n/4]
3335            ]
3336        ) [verts, faces] :
3337    tex=="tri_grid"?
3338        assert(num_defined([gap, n, roughness])==0, str(tex," texture does not accept gap, n or roughness"))  
3339        let(
3340            inset = default(inset,0.05)*sqrt(3)
3341        )
3342        assert(inset>0 && inset<sqrt(3)/6, "tri_grid texture requires inset in (0,1/6)")
3343        let(
3344            adj = opp_ang_to_adj(inset, 30),
3345            y1 = inset / adj_ang_to_opp(1,60),     // i/sqrt(3)
3346            y2 = 2*y1,            // 2*i/sqrt(3)
3347            y3 = 0.5 - y1,
3348            y4 = 0.5 + y1,
3349            y5 = 1 - y2,
3350            y6 = 1 - y1
3351        )
3352        [
3353            [
3354                [0,0,0], [1,0,0],
3355                [adj,y1,1], [1-adj,y1,1],
3356                [0,y2,1], [1,y2,1],
3357                [0.5,0.5-y2,1],
3358                [0,y3,1], [0.5-adj,y3,1], [0.5+adj,y3,1], [1,y3,1],
3359                [0,0.5,0], [0.5,0.5,0], [1,0.5,0],
3360                [0,y4,1], [0.5-adj,y4,1], [0.5+adj,y4,1], [1,y4,1],
3361                [0.5,0.5+y2,1],
3362                [0,y5,1], [1,y5,1],
3363                [adj,y6,1], [1-adj,y6,1],
3364                [0,1,0], [1,1,0],
3365            ], [
3366               [0,2,3], [0,3,1], [2,6,3], [0,12,2], [2,12,6], [3,6,12], [3,12,1],
3367               [0,4,8], [0,8,12], [4,7,8], [7,11,12], [7,12,8],
3368               [1,12,9], [1,9,5], [5,9,10], [9,12,13], [9,13,10],
3369               [11,14,15], [11,15,12], [19,15,14], [19,23,12], [19,12,15],
3370               [12,16,13], [16,17,13], [16,20,17], [12,24,20], [12,20,16],
3371               [21,22,18], [21,23,24], [21,24,22], [12,23,21], [12,21,18],
3372               [12,18,22], [12,22,24],
3373            ]
3374        ] :
3375    tex=="hex_grid"?
3376        assert(num_defined([gap, n, roughness])==0, str(tex," texture does not accept gap, n or roughness"))
3377        let(
3378            inset=default(inset,0.1)
3379        )
3380        assert(inset>0 && inset<0.5)
3381        let(
3382            diag=opp_ang_to_hyp(inset,60),
3383            side=adj_ang_to_opp(1,30),
3384            hyp=adj_ang_to_hyp(0.5,30),
3385            sc = 1/3/hyp,
3386            hex=[ [1,2/6,0], [1/2,1/6,0], [0,2/6,0], [0,4/6,0], [1/2,5/6,0], [1,4/6,0] ]
3387        ) [
3388            [
3389                each hex,
3390                each move([0.5,0.5], p=yscale(sc, p=path3d(ellipse(d=1-2*inset, circum=true, spin=-30,$fn=6),1))),
3391                hex[0]-[0,diag*sc,-1],
3392                for (ang=[270+60,270-60]) hex[1]+yscale(sc, p=cylindrical_to_xyz(diag,ang,1)),
3393                hex[2]-[0,diag*sc,-1],
3394                [0,0,1], [0.5-inset,0,1], [0.5,0,0], [0.5+inset,0,1], [1,0,1],
3395                hex[3]+[0,diag*sc,1],
3396                for (ang=[90+60,90-60]) hex[4]+yscale(sc, p=cylindrical_to_xyz(diag,ang,1)),
3397                hex[5]+[0,diag*sc,1],
3398                [0,1,1], [0.5-inset,1,1], [0.5,1,0], [0.5+inset,1,1], [1,1,1],
3399            ], [
3400                for (i=[0:2:5]) let(b=6) [b+i, b+(i+1)%6, b+(i+2)%6], [6,8,10],
3401                for (i=[0:1:5]) each [ [i, (i+1)%6, (i+1)%6+6], [i, (i+1)%6+6, i+6] ],
3402                [19,13,12], [19,12,20], [17,16,15], [17,15,14],
3403                [21,25,26], [21,26,22], [23,28,29], [23,29,24],
3404                [0,12,13], [0,13,1], [1,14,15], [1,15,2],
3405                [3,21,22], [3,22,4], [4,23,24], [4,24,5],
3406                [1,13,19], [1,19,18], [1,18,17], [1,17,14],
3407                [4,22,26], [4,26,27], [4,27,28], [4,28,23],
3408            ]
3409        ] :
3410    tex=="rough"?
3411        assert(num_defined([gap,inset])==0, str(tex," texture does not accept gap or inset"))
3412        let(
3413            n = default(n,32),
3414            rough = default(roughness, 0.2)
3415        ) [
3416            for (y = [0:1:n-1])
3417            rands(0, rough, n, seed=123456+29*y)
3418        ] :
3419    assert(false, str("Unrecognized texture name: ", tex));
3420
3421
3422/// Function&Module: _textured_linear_sweep()
3423/// Usage: As Function
3424///   vnf = _textured_linear_sweep(region, texture, tex_size, h, ...);
3425///   vnf = _textured_linear_sweep(region, texture, counts=, h=, ...);
3426/// Usage: As Module
3427///   _textured_linear_sweep(region, texture, tex_size, h, ...) [ATTACHMENTS];
3428///   _textured_linear_sweep(region, texture, counts=, h=, ...) [ATTACHMENTS];
3429/// Topics: Sweep, Extrusion, Textures, Knurling
3430/// See Also: heightfield(), cylindrical_heightfield(), texture()
3431/// Description:
3432///   Given a [[Region|regions.scad]], creates a linear extrusion of it vertically, optionally twisted, scaled, and/or shifted,
3433///   with a given texture tiled evenly over the side surfaces.  The texture can be given in one of three ways:
3434///   - As a texture name string. (See {{texture()}} for supported named textures.)
3435///   - As a 2D array of evenly spread height values. (AKA a heightfield.)
3436///   - As a VNF texture tile.  A VNF tile exactly defines a surface from `[0,0]` to `[1,1]`, with the Z coordinates
3437///     being the height of the texture point from the surface.  VNF tiles MUST be able to tile in both X and Y
3438///     directions with no gaps, with the front and back edges aligned exactly, and the left and right edges as well.
3439///   One script to convert a grayscale image to a texture heightfield array in a .scad file can be found at:
3440///   https://raw.githubusercontent.com/BelfrySCAD/BOSL2/master/scripts/img2scad.py
3441/// Arguments:
3442///   region = The [[Region|regions.scad]] to sweep/extrude.
3443///   texture = A texture name string, or a rectangular array of scalar height values (0.0 to 1.0), or a VNF tile that defines the texture to apply to vertical surfaces.  See {{texture()}} for what named textures are supported.
3444///   tex_size = An optional 2D target size for the textures.  Actual texture sizes will be scaled somewhat to evenly fit the available surface. Default: `[5,5]`
3445///   h / l = The height to extrude/sweep the path.
3446///   ---
3447///   counts = If given instead of tex_size, gives the tile repetition counts for textures over the surface length and height.
3448///   inset = If numeric, lowers the texture into the surface by that amount, before the tex_scale multiplier is applied.  If `true`, insets by exactly `1`.  Default: `false`
3449///   rot = If true, rotates the texture 90º.
3450///   tex_scale = Scaling multiplier for the texture depth.
3451///   twist = Degrees of twist for the top of the extrustion/sweep, compared to the bottom.  Default: 0
3452///   scale = Scaling multiplier for the top of the extrustion/sweep, compared to the bottom.  Default: 1
3453///   shift = [X,Y] amount to translate the top, relative to the bottom.  Default: [0,0]
3454///   style = The triangulation style used.  See {{vnf_vertex_array()}} for valid styles.  Used only with heightfield type textures. Default: `"min_edge"`
3455///   samples = Minimum number of "bend points" to have in VNF texture tiles.  Default: 8
3456///   anchor = Translate so anchor point is at origin (0,0,0).  See [anchor](attachments.scad#subsection-anchor).  Default: `CENTER`
3457///   spin = Rotate this many degrees around the Z axis after anchor.  See [spin](attachments.scad#subsection-spin).  Default: `0`
3458///   orient = Vector to rotate top towards, after spin.  See [orient](attachments.scad#subsection-orient).  Default: `UP`
3459/// Extra Anchors:
3460///   centroid_top = The centroid of the top of the shape, oriented UP.
3461///   centroid = The centroid of the center of the shape, oriented UP.
3462///   centroid_bot = The centroid of the bottom of the shape, oriented DOWN.
3463
3464function _get_vnf_tile_edges(texture) =
3465    let(
3466        verts = texture[0],
3467        faces = texture[1],
3468        everts = [for (v = verts) (v.x==0 || v.y==0 || v.x==1 || v.y==1)],
3469        uc = unique_count([
3470            for (face = faces, i = idx(face))
3471            let(edge = select(face,i,i+1), i1 = min(edge), i2 = max(edge))
3472            if (everts[i1] && everts[i2])
3473            [i1, i2]
3474        ]),
3475        edges = uc[0], counts = uc[1],
3476        uedges = [for (i = idx(edges)) if (counts[i] == 1) edges[i] ]
3477    ) uedges;
3478
3479
3480function _validate_texture(texture) =
3481    is_vnf(texture)
3482      ? let( // Validate VNF tile texture
3483            bounds = pointlist_bounds(texture[0]),
3484            min_xy = point2d(bounds[0]),
3485            max_xy = point2d(bounds[1])
3486        )
3487        //assert(min_xy==[0,0] && max_xy==[1,1],"VNF tiles must span exactly from [0,0] to [1,1] in the X and Y components."))
3488        assert(all_nonnegative(concat(min_xy,[1,1]-max_xy)), "VNF tile X and Y components must be between 0 and 1.")
3489        let(
3490            verts = texture[0],
3491            uedges = _get_vnf_tile_edges(texture),
3492            edge_verts = [for (i = unique(flatten(uedges))) verts[i] ],
3493            hverts = [for(v = edge_verts) if(v.x==0 || v.x==1) v],
3494            vverts = [for(v = edge_verts) if(v.y==0 || v.y==1) v],
3495            allgoodx = all(hverts, function(v) any(hverts, function(w) approx(w,[1-v.x, v.y, v.z]))),
3496            allgoody = all(vverts, function(v) any(vverts, function(w) approx(w,[v.x, 1-v.y, v.z])))
3497        )
3498        assert(allgoodx && allgoody, "All VNF tile edge vertices must line up with a vertex on the opposite side of the tile.")
3499        true
3500      : // Validate heightfield texture.
3501        assert(is_matrix(texture), "Malformed texture.")
3502        let( tex_dim = list_shape(texture) )
3503        assert(len(tex_dim) == 2, "Heightfield texture must be a 2D square array of scalar heights.")
3504        assert(all_defined(tex_dim), "Heightfield texture must be a 2D square array of scalar heights.")
3505        true;
3506
3507
3508function _textured_linear_sweep(
3509    region, texture, tex_size=[5,5],
3510    h, counts, inset=false, rot=false,
3511    tex_scale=1, twist, scale, shift,
3512    style="min_edge", l, caps=true, 
3513    height, length, samples,
3514    anchor=CENTER, spin=0, orient=UP
3515) =
3516    assert(is_path(region,[2]) || is_region(region))
3517    assert(is_undef(samples) || is_int(samples))
3518    assert(counts==undef || is_vector(counts,2))
3519    assert(tex_size==undef || is_vector(tex_size,2))
3520    assert(is_bool(rot) || in_list(rot,[0,90,180,270]))
3521    assert(is_bool(caps) || is_bool_list(caps,2))
3522    let(
3523        caps = is_bool(caps) ? [caps,caps] : caps,
3524        regions = is_path(region,2)? [[region]] : region_parts(region),
3525        tex = is_string(texture)? texture(texture) : texture,
3526        texture = !rot? tex :
3527            is_vnf(tex)? zrot(is_num(rot)?rot:90, cp=[1/2,1/2], p=tex) :
3528            rot==180? reverse([for (row=tex) reverse(row)]) :
3529            rot==270? [for (row=transpose(tex)) reverse(row)] :
3530            reverse(transpose(tex)),
3531        h = first_defined([h, l, height, length, 1]),
3532        inset = is_num(inset)? inset : inset? 1 : 0,
3533        twist = default(twist, 0),
3534        shift = default(shift, [0,0]),
3535        scale = scale==undef? [1,1,1] :
3536            is_num(scale)? [scale,scale,1] : scale,
3537        samples = !is_vnf(texture)? len(texture[0]) :
3538            is_num(samples)? samples : 8,
3539        check_tex = _validate_texture(texture),
3540        sorted_tile =
3541            !is_vnf(texture)? texture :
3542            let(
3543                s = 1 / max(1, samples),
3544                vnf = samples<=1? texture :
3545                    let(
3546                        vnft = vnf_slice(texture, "X", list([s:s:1-s/2])),
3547                        zvnf = [
3548                            [
3549                                for (p=vnft[0]) [
3550                                    approx(p.x,0)? 0 : approx(p.x,1)? 1 : p.x,
3551                                    approx(p.y,0)? 0 : approx(p.y,1)? 1 : p.y,
3552                                    p.z
3553                                ]
3554                            ],
3555                            vnft[1]
3556                        ]
3557                    ) zvnf
3558            ) _vnf_sort_vertices(vnf, idx=[1,0]),
3559        vertzs = !is_vnf(sorted_tile)? undef :
3560            group_sort(sorted_tile[0], idx=1),
3561        tpath = is_vnf(sorted_tile)
3562            ? _find_vnf_tile_edge_path(sorted_tile,0)
3563            : let(
3564                  row = sorted_tile[0],
3565                  rlen = len(row)
3566              ) [for (i = [0:1:rlen]) [i/rlen, row[i%rlen]]],
3567        tmat = scale(scale) * zrot(twist) * up(h/2),
3568        pre_skew_vnf = vnf_join([
3569            for (rgn = regions) let(
3570                walls_vnf = vnf_join([
3571                    for (path = rgn) let(
3572                        path = reverse(path),
3573                        plen = path_length(path, closed=true),
3574                        counts = is_vector(counts,2)? counts :
3575                            is_vector(tex_size,2)
3576                              ? [round(plen/tex_size.x), max(1,round(h/tex_size.y)), ]
3577                              : [ceil(6*plen/h), 6],
3578                        obases = resample_path(path, n=counts.x * samples, closed=true),
3579                        onorms = path_normals(obases, closed=true),
3580                        bases = list_wrap(obases),
3581                        norms = list_wrap(onorms),
3582                        vnf = is_vnf(texture)
3583                          ? let( // VNF tile texture
3584                                row_vnf = vnf_join([
3585                                    for (j = [0:1:counts.x-1]) [
3586                                        [
3587                                            for (group = vertzs)
3588                                            each [
3589                                                for (vert = group) let(
3590                                                    u = floor((j + vert.x) * samples),
3591                                                    uu = ((j + vert.x) * samples) - u,
3592                                                    texh = (vert.z - inset) * tex_scale,
3593                                                    base = lerp(bases[u], select(bases,u+1), uu),
3594                                                    norm = unit(lerp(norms[u], select(norms,u+1), uu)),
3595                                                    xy = base + norm * texh
3596                                                ) point3d(xy,vert.y)
3597                                            ]
3598                                        ],
3599                                        sorted_tile[1]
3600                                    ]
3601                                ]),
3602                                sorted_row = _vnf_sort_vertices(row_vnf, idx=[1,0]),
3603                                rvertzs = group_sort(sorted_row[0], idx=1),
3604                                vnf1 = vnf_join([
3605                                    for (i = [0:1:counts.y-1]) [
3606                                        [
3607                                            for (group = rvertzs) let(
3608                                                v = (i + group[0].z) / counts.y,
3609                                                sc = lerp([1,1,1], scale, v),
3610                                                mat = scale(sc) *
3611                                                    zrot(twist*v) *
3612                                                    up(((i/counts.y)-0.5)*h) *
3613                                                    zscale(h/counts.y)
3614                                            ) each apply(mat, group)
3615                                        ],
3616                                        sorted_row[1]
3617                                    ]
3618                                ])
3619                            ) vnf1
3620                          : let( // Heightfield texture
3621                                texcnt = [len(texture[0]), len(texture)],
3622                                tile_rows = [
3623                                    for (ti = [0:1:texcnt.y-1])
3624                                    path3d([
3625                                        for (j = [0:1:counts.x])
3626                                        for (tj = [0:1:texcnt.x-1])
3627                                        if (j != counts.x || tj == 0)
3628                                        let(
3629                                            part = (j + (tj/texcnt.x)) * samples,
3630                                            u = floor(part),
3631                                            uu = part - u,
3632                                            texh = (texture[ti][tj] - inset) * tex_scale,
3633                                            base = lerp(bases[u], select(bases,u+1), uu),
3634                                            norm = unit(lerp(norms[u], select(norms,u+1), uu)),
3635                                            xy = base + norm * texh
3636                                        ) xy
3637                                    ])
3638                                ],
3639                                tiles = [
3640                                    for (i = [0:1:counts.y], ti = [0:1:texcnt.y-1])
3641                                    if (i != counts.y || ti == 0)
3642                                    let(
3643                                        v = (i + (ti/texcnt.y)) / counts.y,
3644                                        sc = lerp([1, 1, 1], scale, v),
3645                                        mat = up((v-0.5)*h) *
3646                                              scale(sc) *
3647                                              zrot(twist*v)
3648                                    ) apply(mat, tile_rows[(texcnt.y-ti)%texcnt.y])
3649                                ]
3650                            ) vnf_vertex_array(
3651                                tiles, caps=false, style=style,
3652                                col_wrap=true, row_wrap=false,
3653                                reverse=true
3654                            )
3655                    ) vnf
3656                ]),
3657                brgn = [
3658                    for (path = rgn) let(
3659                        path = reverse(path),
3660                        plen = path_length(path, closed=true),
3661                        counts = is_vector(counts,2)? counts :
3662                            is_vector(tex_size,2)
3663                              ? [round(plen/tex_size.x), max(1,round(h/tex_size.y)), ]
3664                              : [ceil(6*plen/h), 6],
3665                        obases = resample_path(path, n=counts.x * samples, closed=true),
3666                        onorms = path_normals(obases, closed=true),
3667                        bases = list_wrap(obases),
3668                        norms = list_wrap(onorms),
3669                        nupath = [
3670                            for (j = [0:1:counts.x-1], vert = tpath) let(
3671                                part = (j + vert.x) * samples,
3672                                u = floor(part),
3673                                uu = part - u,
3674                                texh = (vert.y - inset) * tex_scale,
3675                                base = lerp(bases[u], select(bases,u+1), uu),
3676                                norm = unit(lerp(norms[u], select(norms,u+1), uu)),
3677                                xy = base + norm * texh
3678                            ) xy
3679                        ]
3680                    ) nupath
3681                ],
3682                bot_vnf = !caps[0] || brgn==[[]] ? EMPTY_VNF
3683                    : vnf_from_region(brgn, down(h/2), reverse=true),
3684                top_vnf = !caps[1] || brgn==[[]] ? EMPTY_VNF
3685                    : vnf_from_region(brgn, tmat, reverse=false)
3686            ) vnf_join([walls_vnf, bot_vnf, top_vnf])
3687        ]),
3688        skmat = down(h/2) * skew(sxz=shift.x/h, syz=shift.y/h) * up(h/2),
3689        final_vnf = apply(skmat, pre_skew_vnf),
3690        cent = centroid(region),
3691        anchors = [
3692            named_anchor("centroid_top", point3d(cent, h/2), UP),
3693            named_anchor("centroid",     point3d(cent),      UP),
3694            named_anchor("centroid_bot", point3d(cent,-h/2), DOWN)
3695        ]
3696    ) reorient(anchor,spin,orient, vnf=final_vnf, extent=true, anchors=anchors, p=final_vnf);
3697
3698
3699
3700function _find_vnf_tile_edge_path(vnf, val) =
3701    let(
3702        verts = vnf[0],
3703        fragments = [
3704            for(edge = _get_vnf_tile_edges(vnf))
3705            let(v0 = verts[edge[0]], v1 = verts[edge[1]])
3706            if (approx(v0.y, val) && approx(v1.y, val))
3707            v0.x <= v1.x? [[v0.x,v0.z], [v1.x,v1.z]] :
3708            [[v1.x,v1.z], [v0.x,v0.z]]
3709        ],
3710        sfrags = sort(fragments, idx=[0,1]),
3711        rpath = _assemble_a_path_from_fragments(sfrags)[0],
3712        opath = rpath==[]? []
3713              : rpath[0].x > last(rpath).x ? reverse(rpath)
3714              : rpath
3715    ) opath;
3716
3717
3718/// Function&Module: _textured_revolution()
3719/// Usage: As Function
3720///   vnf = _textured_revolution(shape, texture, tex_size, [tex_scale=], ...);
3721///   vnf = _textured_revolution(shape, texture, counts=, [tex_scale=], ...);
3722/// Usage: As Module
3723///   _textured_revolution(shape, texture, tex_size, [tex_scale=], ...) [ATTACHMENTS];
3724///   _textured_revolution(shape, texture, counts=, [tex_scale=], ...) [ATTACHMENTS];
3725/// Topics: Sweep, Extrusion, Textures, Knurling
3726/// See Also: heightfield(), cylindrical_heightfield(), texture()
3727/// Description:
3728///   Given a 2D region or path, fully in the X+ half-plane, revolves that shape around the Z axis (after rotating its Y+ to Z+).
3729///   This creates a solid from that surface of revolution, possibly capped top and bottom, with the sides covered in a given tiled texture.
3730///   The texture can be given in one of three ways:
3731///   - As a texture name string. (See {{texture()}} for supported named textures.)
3732///   - As a 2D array of evenly spread height values. (AKA a heightfield.)
3733///   - As a VNF texture tile.  A VNF tile exactly defines a surface from `[0,0]` to `[1,1]`, with the Z coordinates
3734///     being the height of the texture point from the surface.  VNF tiles MUST be able to tile in both X and Y
3735///     directions with no gaps, with the front and back edges aligned exactly, and the left and right edges as well.
3736///   One script to convert a grayscale image to a texture heightfield array in a .scad file can be found at:
3737///   https://raw.githubusercontent.com/BelfrySCAD/BOSL2/master/scripts/img2scad.py
3738/// Arguments:
3739///   shape = The path or region to sweep/extrude.
3740///   texture = A texture name string, or a rectangular array of scalar height values (0.0 to 1.0), or a VNF tile that defines the texture to apply to the revolution surface.  See {{texture()}} for what named textures are supported.
3741///   tex_size = An optional 2D target size for the textures.  Actual texture sizes will be scaled somewhat to evenly fit the available surface. Default: `[5,5]`
3742///   tex_scale = Scaling multiplier for the texture depth.
3743///   ---
3744///   inset = If numeric, lowers the texture into the surface by that amount, before the tex_scale multiplier is applied.  If `true`, insets by exactly `1`.  Default: `false`
3745///   rot = If true, rotates the texture 90º.
3746///   shift = [X,Y] amount to translate the top, relative to the bottom.  Default: [0,0]
3747///   closed = If false, and shape is given as a path, then the revolved path will be sealed to the axis of rotation with untextured caps.  Default: `true`
3748///   taper = If given, and `closed=false`, tapers the texture height to zero over the first and last given percentage of the path.  If given as a lookup table with indices between 0 and 100, uses the percentage lookup table to ramp the texture heights.  Default: `undef` (no taper)
3749///   angle = The number of degrees counter-clockwise from X+ to revolve around the Z axis.  Default: `360`
3750///   style = The triangulation style used.  See {{vnf_vertex_array()}} for valid styles.  Used only with heightfield type textures. Default: `"min_edge"`
3751///   counts = If given instead of tex_size, gives the tile repetition counts for textures over the surface length and height.
3752///   samples = Minimum number of "bend points" to have in VNF texture tiles.  Default: 8
3753///   anchor = Translate so anchor point is at origin (0,0,0).  See [anchor](attachments.scad#subsection-anchor).  Default: `CENTER`
3754///   spin = Rotate this many degrees around the Z axis after anchor.  See [spin](attachments.scad#subsection-spin).  Default: `0`
3755///   orient = Vector to rotate top towards, after spin.  See [orient](attachments.scad#subsection-orient).  Default: `UP`
3756/// Anchor Types:
3757///   "hull" = Anchors to the virtual convex hull of the shape.
3758///   "intersect" = Anchors to the surface of the shape.
3759
3760function _textured_revolution(
3761    shape, texture, tex_size, tex_scale=1,
3762    inset=false, rot=false, shift=[0,0],
3763    taper, closed=true, angle=360,
3764    counts, samples,
3765    style="min_edge", atype="intersect",
3766    anchor=CENTER, spin=0, orient=UP
3767) =
3768    assert(angle>0 && angle<=360)
3769    assert(is_path(shape,[2]) || is_region(shape))
3770    assert(is_undef(samples) || is_int(samples))
3771    assert(is_bool(closed))
3772    assert(counts==undef || is_vector(counts,2))
3773    assert(tex_size==undef || is_vector(tex_size,2))
3774    assert(is_bool(rot) || in_list(rot,[0,90,180,270]))
3775    let( taper_is_ok = is_undef(taper) || (is_finite(taper) && taper>=0 && taper<50) || is_path(taper,2) )
3776    assert(taper_is_ok, "Bad taper= value.")
3777    assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"")
3778    let(
3779        regions = !is_path(shape,2)? region_parts(shape) :
3780            closed? region_parts([shape]) :
3781            let(
3782                clpoly = [[0,shape[0].y], each shape, [0,last(shape).y]],
3783                dpoly = deduplicate(clpoly),
3784                cwpoly = is_polygon_clockwise(dpoly) ? dpoly : reverse(dpoly)
3785            )
3786            [[ select(cwpoly,1,-2) ]],
3787        checks = [
3788            for (rgn=regions, path=rgn)
3789            assert(all(path, function(pt) pt.x>=0))
3790        ]
3791    )
3792    assert(closed || is_path(shape,2))
3793    let(
3794        tex = is_string(texture)? texture(texture) : texture,
3795        texture = !rot? tex :
3796            is_vnf(tex)? zrot(is_num(rot)?rot:90, cp=[1/2,1/2], p=tex) :
3797            rot==180? reverse([for (row=tex) reverse(row)]) :
3798            rot==270? [for (row=transpose(tex)) reverse(row)] :
3799            reverse(transpose(tex)),
3800        check_tex = _validate_texture(texture),
3801        inset = is_num(inset)? inset : inset? 1 : 0,
3802        samples = !is_vnf(texture)? len(texture) :
3803            is_num(samples)? samples : 8,
3804        bounds = pointlist_bounds(flatten(flatten(regions))),
3805        maxx = bounds[1].x,
3806        miny = bounds[0].y,
3807        maxy = bounds[1].y,
3808        h = maxy - miny,
3809        circumf = 2 * PI * maxx,
3810        tile = !is_vnf(texture)? texture :
3811            let(
3812                utex = samples<=1? texture :
3813                    let(
3814                        s = 1 / samples,
3815                        slices = list([s : s : 1-s/2]),
3816                        vnfx = vnf_slice(texture, "X", slices),
3817                        vnfy = vnf_slice(vnfx, "Y", slices),
3818                        vnft = vnf_triangulate(vnfy),
3819                        zvnf = [
3820                            [
3821                                for (p=vnft[0]) [
3822                                    approx(p.x,0)? 0 : approx(p.x,1)? 1 : p.x,
3823                                    approx(p.y,0)? 0 : approx(p.y,1)? 1 : p.y,
3824                                    p.z
3825                                ]
3826                            ],
3827                            vnft[1]
3828                        ]
3829                    ) zvnf
3830            ) _vnf_sort_vertices(utex, idx=[0,1]),
3831        vertzs = is_vnf(texture)? group_sort(tile[0], idx=0) : undef,
3832        bpath = is_vnf(tile)
3833            ? _find_vnf_tile_edge_path(tile,1)
3834            : let(
3835                  row = tile[0],
3836                  rlen = len(row)
3837              ) [for (i = [0:1:rlen]) [i/rlen, row[i%rlen]]],
3838        counts_x = is_vector(counts,2)? counts.x :
3839            is_vector(tex_size,2)
3840              ? max(1,round(angle/360*circumf/tex_size.x))
3841              : ceil(6*angle/360*circumf/h),
3842        taper_lup = closed || is_undef(taper)? [[-1,1],[2,1]] :
3843            is_num(taper)? [[-1,0], [0,0], [taper/100+EPSILON,1], [1-taper/100-EPSILON,1], [1,0], [2,0]] :
3844            is_path(taper,2)? let(
3845                retaper = [
3846                    for (t=taper)
3847                    assert(t[0]>=0 && t[0]<=100, "taper lookup indices must be between 0 and 100 inclusive.")
3848                    [t[0]/100, t[1]]
3849                ],
3850                taperout = [[-1,retaper[0][1]], each retaper, [2,last(retaper)[1]]]
3851            ) taperout :
3852            assert(false, "Bad taper= argument value."),
3853        full_vnf = vnf_join([
3854            for (rgn = regions) let(
3855                rgn_wall_vnf = vnf_join([
3856                    for (path = rgn) let(
3857                        plen = path_length(path, closed=closed),
3858                        counts_y = is_vector(counts,2)? counts.y :
3859                            is_vector(tex_size,2)? max(1,round(plen/tex_size.y)) : 6,
3860                        obases = resample_path(path, n=counts_y * samples + (closed?0:1), closed=closed),
3861                        onorms = path_normals(obases, closed=closed),
3862                        rbases = closed? list_wrap(obases) : obases,
3863                        rnorms = closed? list_wrap(onorms) : onorms,
3864                        bases = xrot(90, p=path3d(rbases)),
3865                        norms = xrot(90, p=path3d(rnorms)),
3866                        vnf = is_vnf(texture)
3867                          ? vnf_join([ // VNF tile texture
3868                                for (j = [0:1:counts_y-1])
3869                                [
3870                                    [
3871                                        for (group = vertzs) each [
3872                                            for (vert = group) let(
3873                                                part = (j + (1-vert.y)) * samples,
3874                                                u = floor(part),
3875                                                uu = part - u,
3876                                                base = lerp(select(bases,u), select(bases,u+1), uu),
3877                                                norm = unit(lerp(select(norms,u), select(norms,u+1), uu)),
3878                                                tex_scale = tex_scale * lookup(part/samples/counts_y, taper_lup),
3879                                                texh = (vert.z - inset) * tex_scale * (base.x / maxx),
3880                                                xyz = base - norm * texh
3881                                            ) zrot(vert.x*angle/counts_x, p=xyz)
3882                                        ]
3883                                    ],
3884                                    tile[1]
3885                                ]
3886                            ])
3887                          : let( // Heightfield texture
3888                                texcnt = [len(texture[0]), len(texture)],
3889                                tiles = transpose([
3890                                    for (j = [0,1], tj = [0:1:texcnt.x-1])
3891                                    if (j == 0 || tj == 0)
3892                                    let(
3893                                        v = (j + (tj/texcnt.x)) / counts_x,
3894                                        mat = zrot(v*angle)
3895                                    ) apply(mat, [
3896                                        for (i = [0:1:counts_y-(closed?1:0)], ti = [0:1:texcnt.y-1])
3897                                        if (i != counts_y || ti == 0)
3898                                        let(
3899                                            part = (i + (ti/texcnt.y)) * samples,
3900                                            u = floor(part),
3901                                            uu = part - u,
3902                                            base = lerp(bases[u], select(bases,u+1), uu),
3903                                            norm = unit(lerp(norms[u], select(norms,u+1), uu)),
3904                                            tex_scale = tex_scale * lookup(part/samples/counts_y, taper_lup),
3905                                            texh = (texture[ti][tj] - inset) * tex_scale * (base.x / maxx),
3906                                            xyz = base - norm * texh
3907                                        ) xyz
3908                                    ])
3909                                ])
3910                            ) vnf_vertex_array(
3911                                tiles, caps=false, style=style,
3912                                col_wrap=false, row_wrap=closed
3913                            )
3914                    ) vnf
3915                ]),
3916                walls_vnf = vnf_join([
3917                    for (i = [0:1:counts_x-1])
3918                    zrot(i*angle/counts_x, rgn_wall_vnf)
3919                ]),
3920                endcap_vnf = angle == 360? EMPTY_VNF :
3921                    let(
3922                        cap_rgn = [
3923                            for (path = rgn) let(
3924                                plen = path_length(path, closed=closed),
3925                                counts_y = is_vector(counts,2)? counts.y :
3926                                    is_vector(tex_size,2)? max(1,round(plen/tex_size.y)) : 6,
3927                                obases = resample_path(path, n=counts_y * samples + (closed?0:1), closed=closed),
3928                                onorms = path_normals(obases, closed=closed),
3929                                bases = closed? list_wrap(obases) : obases,
3930                                norms = closed? list_wrap(onorms) : onorms,
3931                                ppath = is_vnf(texture)
3932                                  ? [ // VNF tile texture
3933                                        for (j = [0:1:counts_y-1])
3934                                        for (group = vertzs, vert = reverse(group))
3935                                        if (approx(vert.x, 0)) let(
3936                                            part = (j + (1 - vert.y)) * samples,
3937                                            u = floor(part),
3938                                            uu = part - u,
3939                                            base = lerp(select(bases,u), select(bases,u+1), uu),
3940                                            norm = unit(lerp(select(norms,u), select(norms,u+1), uu)),
3941                                            tex_scale = tex_scale * lookup(part/samples/counts_y, taper_lup),
3942                                            texh = (vert.z - inset) * tex_scale * (base.x / maxx),
3943                                            xyz = base - norm * texh
3944                                        ) xyz
3945                                    ]
3946                                  : let( // Heightfield texture
3947                                        texcnt = [len(texture[0]), len(texture)]
3948                                    ) [
3949                                        for (i = [0:1:counts_y-(closed?1:0)], ti = [0:1:texcnt.y-1])
3950                                        if (i != counts_y || ti == 0)
3951                                        let(
3952                                            part = (i + (ti/texcnt.y)) * samples,
3953                                            u = floor(part),
3954                                            uu = part - u,
3955                                            base = lerp(bases[u], select(bases,u+1), uu),
3956                                            norm = unit(lerp(norms[u], select(norms,u+1), uu)),
3957                                            tex_scale = tex_scale * lookup(part/samples/counts_y, taper_lup),
3958                                            texh = (texture[ti][0] - inset) * tex_scale * (base.x / maxx),
3959                                            xyz = base - norm * texh
3960                                        ) xyz
3961                                    ],
3962                                path = closed? ppath : [
3963                                    [0, ppath[0].y],
3964                                    each ppath,
3965                                    [0, last(ppath).y],
3966                                ]
3967                            ) deduplicate(path, closed=closed)
3968                        ],
3969                        vnf2 = vnf_from_region(cap_rgn, xrot(90), reverse=false),
3970                        vnf3 = vnf_from_region(cap_rgn, rot([90,0,angle]), reverse=true)
3971                    ) vnf_join([vnf2, vnf3]),
3972                allcaps_vnf = closed? EMPTY_VNF :
3973                    let(
3974                        plen = path_length(rgn[0], closed=closed),
3975                        counts_y = is_vector(counts,2)? counts.y :
3976                            is_vector(tex_size,2)? max(1,round(plen/tex_size.y)) : 6,
3977                        obases = resample_path(rgn[0], n=counts_y * samples + (closed?0:1), closed=closed),
3978                        onorms = path_normals(obases, closed=closed),
3979                        rbases = closed? list_wrap(obases) : obases,
3980                        rnorms = closed? list_wrap(onorms) : onorms,
3981                        bases = xrot(90, p=path3d(rbases)),
3982                        norms = xrot(90, p=path3d(rnorms)),
3983                        caps_vnf = vnf_join([
3984                            for (j = [-1,0]) let(
3985                                base = select(bases,j),
3986                                norm = unit(select(norms,j)),
3987                                ppath = [
3988                                    for (vert = bpath) let(
3989                                        uang = vert.x / counts_x,
3990                                        tex_scale = tex_scale * lookup([0,1][j+1], taper_lup),
3991                                        texh = (vert.y - inset) * tex_scale * (base.x / maxx),
3992                                        xyz = base - norm * texh
3993                                    ) zrot(angle*uang, p=xyz)
3994                                ],
3995                                pplen = len(ppath),
3996                                zed = j<0? max(column(ppath,2)) :
3997                                    min(column(ppath,2)),
3998                                slice_vnf = [
3999                                    [
4000                                        each ppath,
4001                                        [0, 0, zed],
4002                                    ], [
4003                                        for (i = [0:1:pplen-2])
4004                                            j<0? [pplen, i, (i+1)%pplen] :
4005                                            [pplen, (i+1)%pplen, i]
4006                                    ]
4007                                ],
4008                                cap_vnf = vnf_join([
4009                                    for (i = [0:1:counts_x-1])
4010                                        zrot(i*angle/counts_x, p=slice_vnf)
4011                                ])
4012                            ) cap_vnf
4013                        ])
4014                    ) caps_vnf
4015            ) vnf_join([walls_vnf, endcap_vnf, allcaps_vnf])
4016        ]),
4017        skmat = down(-miny) * skew(sxz=shift.x/h, syz=shift.y/h) * up(-miny),
4018        skvnf = apply(skmat, full_vnf),
4019        geom = atype=="intersect"
4020              ? attach_geom(vnf=skvnf, extent=false)
4021              : attach_geom(vnf=skvnf, extent=true)
4022    ) reorient(anchor,spin,orient, geom=geom, p=skvnf);
4023
4024
4025module _textured_revolution(
4026    shape, texture, tex_size, tex_scale=1,
4027    inset=false, rot=false, shift=[0,0],
4028    taper, closed=true, angle=360,
4029    style="min_edge", atype="intersect",
4030    convexity=10, counts, samples,
4031    anchor=CENTER, spin=0, orient=UP
4032) {
4033    dummy = assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"");
4034    vnf = _textured_revolution(
4035        shape, texture, tex_size=tex_size,
4036        tex_scale=tex_scale, inset=inset, rot=rot,
4037        taper=taper, closed=closed, style=style,
4038        shift=shift, angle=angle,
4039        samples=samples, counts=counts
4040    );
4041    geom = atype=="intersect"
4042          ? attach_geom(vnf=vnf, extent=false)
4043          : attach_geom(vnf=vnf, extent=true);
4044    attachable(anchor,spin,orient, geom=geom) {
4045        vnf_polyhedron(vnf, convexity=convexity);
4046        children();
4047    }
4048}
4049
4050
4051
4052// vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap