1/////////////////////////////////////////////////////////////////////
   2// LibFile: rounding.scad
   3//   Routines to create rounded corners, with either circular rounding,
   4//   or continuous curvature rounding with no sudden curvature transitions.
   5//   Provides rounding of corners or rounding that preserves corner points and curves the edges.
   6//   Also provides some 3D rounding functions, and a powerful function for joining
   7//   two prisms together with a rounded fillet at the joint.  
   8// Includes:
   9//   include <BOSL2/std.scad>
  10//   include <BOSL2/rounding.scad>
  11// FileGroup: Advanced Modeling
  12// FileSummary: Round path corners, rounded prisms, rounded cutouts in tubes, filleted prism joints
  13//////////////////////////////////////////////////////////////////////
  14include <beziers.scad>
  15include <structs.scad>
  16
  17// Section: Types of Roundovers
  18//   The functions and modules in this file support two different types of roundovers and some different mechanisms for specifying
  19//   the size of the roundover.  The usual circular roundover can produce a tactile "bump" where the curvature changes from flat to
  20//   circular.  See https://hackernoon.com/apples-icons-have-that-shape-for-a-very-good-reason-720d4e7c8a14 for details.
  21//   We compute continuous curvature rounding using 4th order Bezier curves.  This type of rounding, which we call "smooth" rounding,
  22//   does not have a "radius" so we need different ways to specify the size of the roundover.  We introduce the `cut` and `joint`
  23//   parameters for this purpose.  They can specify dimensions of circular roundovers, continuous curvature "smooth" roundovers, and even chamfers.  
  24//   .
  25//   The `cut` parameter specifies the distance from the unrounded corner to the rounded tip, so how
  26//   much of the corner to "cut" off.  This can be easier to understand than setting a circular radius, which can be
  27//   unexpectedly extreme when the corner is very sharp.  It also allows a systematic specification of
  28//   corner treatments that are the same size for all corner treatments.
  29//   .
  30//   The `joint` parameter specifies the distance
  31//   away from the corner along the path where the roundover or chamfer should start.  This parameter is good for ensuring that
  32//   your roundover will fit on the polygon or polyhedron, since you can easily tell whether you have enough space, and whether
  33//   adjacent corner treatments will interfere.
  34//   .
  35//   For circular rounding you can use the `radius` or `r` parameter to set the rounding radius.
  36//   .
  37//   For chamfers you can use `width` to set the width of the chamfer.  
  38//   .
  39//   The "smooth" rounding method also has a parameter that specifies how smooth the curvature match is.  This parameter, `k`,
  40//   ranges from 0 to 1, with a default of 0.5.  Larger values gives a more 
  41//   abrupt transition and smaller ones a more gradual transition.  If you set the value much higher
  42//   than 0.8 the curvature changes abruptly enough that though it is theoretically continuous, it may
  43//   not be continuous in practice.  If you set it very small then the transition is so gradual that
  44//   the length of the roundover may be extremely long, and the actual rounded part of the curve may be very small.  
  45// Figure(2D,Med,NoAxes):  Parameters of a "circle" roundover
  46//   h = 18;
  47//   w = 12.6;
  48//   strokewidth = .3;
  49//   example = [[0,0],[w,h],[2*w,0]];
  50//   stroke(example, width=strokewidth*1.5);
  51//   textangle = 90-vector_angle(example)/2;
  52//   theta = vector_angle(example)/2;
  53//   color("green"){ stroke([[w,h], [w,h-18*(1-sin(theta))/cos(theta)]], width=strokewidth, endcaps="arrow2");
  54//                   translate([w-1.75,h-7])scale(.1)rotate(textangle)text("cut",size=14); }
  55//   ll=lerp([w,h], [0,0],18/norm([w,h]-[0,0]) );
  56//   color("blue"){ stroke(_shift_segment([[w,h], ll], -.7), width=strokewidth,endcaps="arrow2");
  57//                  translate([w/2-1.3,h/2+.6])  scale(.1)rotate(textangle)text("joint",size=14);}
  58//   color("red")stroke(
  59//         select(round_corners(example, joint=18, method="circle",$fn=64,closed=false),1,-2),
  60//         width=strokewidth);
  61//   r=18*tan(theta);
  62//   color("black"){
  63//     stroke([ll, [w,h-r-18*(1-sin(theta))/cos(theta)]], width=strokewidth, endcaps="arrow2");
  64//     translate([w/1.6,0])text("radius", size=1.4);
  65//   }
  66// Figure(2D,Med,NoAxes):  Parameters of a "smooth" roundover with the default of `k=0.5`.  Note the long, slow transition from flat to round.  
  67//   h = 18;
  68//   w = 12.6;
  69//   strokewidth = .3;
  70//   example = [[0,0],[w,h],[2*w,0]];
  71//   stroke(example, width=strokewidth*1.5);
  72//   textangle = 90-vector_angle(example)/2;
  73//   color("green"){ stroke([[w,h], [w,h-cos(vector_angle(example)/2) *3/8*h]], width=strokewidth, endcaps="arrow2");
  74//                   translate([w-1.75,h-5.5])scale(.1)rotate(textangle)text("cut",size=14); }
  75//   ll=lerp([w,h], [0,0],18/norm([w,h]-[0,0]) );
  76//   color("blue"){ stroke(_shift_segment([[w,h], ll], -.7), width=strokewidth,endcaps="arrow2");
  77//                  translate([w/2-1.3,h/2+.6])  scale(.1)rotate(textangle)text("joint",size=14);}
  78//   color("red")stroke(
  79//         select(round_corners(example, joint=18, method="smooth",closed=false),1,-2),
  80//         width=strokewidth);
  81// Figure(2D,Med,NoAxes):  Parameters of a "smooth" roundover, with `k=0.75`.  The transition into the roundover is shorter, and faster.  The cut length is bigger for the same joint length.
  82//   h = 18;
  83//   w = 12.6;
  84//   strokewidth = .3;
  85//   example = [[0,0],[w,h],[2*w,0]];
  86//   stroke(example, width=strokewidth*1.5);
  87//   textangle = 90-vector_angle(example)/2;
  88//   color("green"){ stroke([[w,h], [w,h-cos(vector_angle(example)/2) *4/8*h]], width=strokewidth, endcaps="arrow2");
  89//                   translate([w-1.75,h-5.5])scale(.1)rotate(textangle)text("cut",size=14); }
  90//   ll=lerp([w,h], [0,0],18/norm([w,h]-[0,0]) );
  91//   color("blue"){ stroke(_shift_segment([[w,h], ll], -.7), width=strokewidth,endcaps="arrow2");
  92//                  translate([w/2-1.3,h/2+.6])  scale(.1)rotate(textangle)text("joint",size=14);}
  93//   color("red")stroke(
  94//         select(round_corners(example, joint=18, method="smooth",closed=false,k=.75),1,-2),
  95//         width=strokewidth);
  96// Figure(2D,Med,NoAxes):  Parameters of a "smooth" roundover, with `k=0.15`.  The transition is so gradual that it appears that the roundover is much smaller than specified.  The cut length is much smaller for the same joint length.  
  97//   h = 18;
  98//   w = 12.6;
  99//   strokewidth = .3;
 100//   example = [[0,0],[w,h],[2*w,0]];
 101//   stroke(example, width=strokewidth*1.5);
 102//   textangle = 90-vector_angle(example)/2;
 103//   color("green"){ stroke([[w,h], [w,h-cos(vector_angle(example)/2) *1.6/8*h]], width=strokewidth, endcaps="arrow2");
 104//                   translate([w+.3,h])text("cut",size=1.4); }
 105//   ll=lerp([w,h], [0,0],18/norm([w,h]-[0,0]) );
 106//   color("blue"){ stroke(_shift_segment([[w,h], ll], -.7), width=strokewidth,endcaps="arrow2");
 107//                  translate([w/2-1.3,h/2+.6])  scale(.1)rotate(textangle)text("joint",size=14);}
 108//   color("red")stroke(
 109//         select(round_corners(example, joint=18, method="smooth",closed=false,k=.15),1,-2),
 110//         width=strokewidth);
 111// Figure(2D,Med,NoAxes):  Parameters of a symmetric "chamfer".
 112//   h = 18;
 113//   w = 12.6;
 114//   strokewidth = .3;
 115//   example = [[0,0],[w,h],[2*w,0]];
 116//   stroke(example, width=strokewidth*1.5);
 117//   textangle = 90-vector_angle(example)/2;
 118//   color("black"){
 119//        stroke(fwd(1,
 120//         select(round_corners(example, joint=18, method="chamfer",closed=false),1,-2)),
 121//         width=strokewidth,endcaps="arrow2");
 122//        translate([w,.3])text("width", size=1.4,halign="center");
 123//   }     
 124//   color("green"){ stroke([[w,h], [w,h-18*cos(vector_angle(example)/2)]], width=strokewidth, endcaps="arrow2");
 125//                   translate([w-1.75,h-5.5])scale(.1)rotate(textangle)text("cut",size=14); }
 126//   ll=lerp([w,h], [0,0],18/norm([w,h]-[0,0]) );
 127//   color("blue"){ stroke(_shift_segment([[w,h], ll], -.7), width=strokewidth,endcaps="arrow2");
 128//                  translate([w/2-1.3,h/2+.6]) rotate(textangle)text("joint",size=1.4);}
 129//   color("red")stroke(
 130//         select(round_corners(example, joint=18, method="chamfer",closed=false),1,-2),
 131//         width=strokewidth);
 132
 133
 134// Section: Rounding Paths
 135
 136// Function: round_corners()
 137//
 138// Usage:
 139//   rounded_path = round_corners(path, [method], [radius=], [cut=], [joint=], [closed=], [verbose=]);
 140//
 141// Description:
 142//   Takes a 2D or 3D path as input and rounds each corner
 143//   by a specified amount.  The rounding at each point can be different and some points can have zero
 144//   rounding.  The `round_corners()` function supports three types of corner treatment: chamfers, circular rounding,
 145//   and continuous curvature rounding using 4th order bezier curves.  See
 146//   [Types of Roundover](rounding.scad#subsection-types-of-roundover) for details on rounding types.  
 147//   .
 148//   You select the type of rounding using the `method` parameter, which should be `"smooth"` to
 149//   get continuous curvature rounding, `"circle"` to get circular rounding, or `"chamfer"` to get chamfers.  The default is circle
 150//   rounding.  Each method accepts multiple options to specify the amount of rounding.  See
 151//   [Types of Roundover](rounding.scad#subsection-types-of-roundover) for example diagrams.
 152//   .
 153//   * The `cut` parameter specifies the distance from the unrounded corner to the rounded tip, so how
 154//   much of the corner to "cut" off.  
 155//   * The `joint` parameter specifies the distance
 156//   away from the corner along the path where the roundover or chamfer should start.  This makes it easy to ensure your roundover will fit,
 157//   so use it if you want the largest possible roundover.  
 158//   * For circular rounding you can use the `radius` or `r` parameter to set the rounding radius.
 159//   * For chamfers you can use the `width` parameter, which sets the width of the chamfer edge.  
 160//   .
 161//   As explained in [Types of Roundover](rounding.scad#subsection-types-of-roundover), the continuous curvature "smooth"
 162//   type of rounding also accepts the `k` parameter, between 0 and 1, which specifies how fast the curvature changes at
 163//   the joint.  The default is `k=0.5`.  
 164//   .
 165//   If you select curves that are too large to fit the function will fail with an error.  You can set `verbose=true` to
 166//   get a message showing a list of scale factors you can apply to your rounding parameters so that the
 167//   roundovers will fit on the curve.  If the scale factors are larger than one
 168//   then they indicate how much you can increase the curve sizes before collisions will occur.
 169//   .
 170//   The parameters `radius`, `cut`, `joint` and `k` can be numbers, which round every corner using the same parameters, or you
 171//   can specify a list to round each corner with different parameters.  If the curve is not closed then the first and last points
 172//   of the curve are not rounded.  In this case you can specify a full list of points anyway, and the endpoint values are ignored,
 173//   or you can specify a list that has length len(path)-2, omitting the two dummy values.
 174//   .
 175//   If your input path includes collinear points you must use a cut or radius value of zero for those "corners".  You can
 176//   choose a nonzero joint parameter when the collinear points form a 180 degree angle.  This will cause extra points to be inserted. 
 177//   If the collinear points form a spike (0 degree angle) then round_corners will fail. 
 178//   .
 179//   Examples:
 180//   * `method="circle", radius=2`:
 181//       Rounds every point with circular, radius 2 roundover
 182//   * `method="smooth", cut=2`:
 183//       Rounds every point with continuous curvature rounding with a cut of 2, and a default 0.5 smoothing parameter
 184//   * `method="smooth", cut=2, k=0.3`:
 185//       Rounds every point with continuous curvature rounding with a cut of 2, and a very gentle 0.3 smoothness setting
 186//   .
 187//   The number of segments used for roundovers is determined by `$fa`, `$fs` and `$fn` as usual for
 188//   circular roundovers.  For continuous curvature roundovers `$fs` and `$fn` are used and `$fa` is
 189//   ignored.  Note that $fn is interpreted as the number of points on the roundover curve, which is
 190//   not equivalent to its meaning for rounding circles because roundovers are usually small fractions
 191//   of a circular arc.  As usual, $fn overrides $fs.  When doing continuous curvature rounding be sure to use lots of segments or the effect
 192//   will be hidden by the discretization.  Note that if you use $fn with "smooth" then $fn points are added at each corner.
 193//   This guarantees a specific output length.  It also means that if
 194//   you set `joint` nonzero on a flat "corner", with collinear points, you will get $fn points at that "corner."
 195//   If you have two roundovers that fully consume a segment then they share a point where they meet in the segment, which means the output
 196//   point count will be decreased by one.  
 197// Arguments:
 198//   path = list of 2d or 3d points defining the path to be rounded.
 199//   method = rounding method to use.  Set to "chamfer" for chamfers, "circle" for circular rounding and "smooth" for continuous curvature 4th order bezier rounding.  Default: "circle"
 200//   ---
 201//   radius/r = rounding radius, only compatible with `method="circle"`. Can be a number or vector.
 202//   cut = rounding cut distance, compatible with all methods.  Can be a number or vector.
 203//   joint = rounding joint distance, compatible with `method="chamfer"` and `method="smooth"`.  Can be a number or vector.
 204//   width = width of the flat edge created by chamfering, compatible with `method="chamfer"`.  Can be a number or vector. 
 205//   k = continuous curvature smoothness parameter for `method="smooth"`.  Can be a number or vector.  Default: 0.5
 206//   closed = if true treat the path as a closed polygon, otherwise treat it as open.  Default: true.
 207//   verbose = if true display rounding scale factors that show how close roundovers are to overlapping.  Default: false
 208//
 209// Example(2D,Med): Standard circular roundover with radius the same at every point. Compare results at the different corners.
 210//   $fn=36;
 211//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 212//   polygon(round_corners(shape, radius=1));
 213//   color("red") down(.1) polygon(shape);
 214// Example(2D,Med): Circular roundover using the "cut" specification, the same at every corner.
 215//   $fn=36;
 216//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 217//   polygon(round_corners(shape, cut=1));
 218//   color("red") down(.1) polygon(shape);
 219// Example(2D,Med): Continous curvature roundover using "cut", still the same at every corner.  The default smoothness parameter of 0.5 was too gradual for these roundovers to fit, but 0.7 works.
 220//   $fn=36;
 221//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 222//   polygon(round_corners(shape, method="smooth", cut=1, k=0.7));
 223//   color("red") down(.1) polygon(shape);
 224// Example(2D,Med): Continuous curvature roundover using "joint", for the last time the same at every corner.  Notice how small the roundovers are.
 225//   $fn=36;
 226//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 227//   polygon(round_corners(shape, method="smooth", joint=1, k=0.7));
 228//   color("red") down(.1) polygon(shape);
 229// Example(2D,Med): Circular rounding, different at every corner, some corners left unrounded
 230//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 231//   radii = [1.8, 0, 2, 0.3, 1.2, 0];
 232//   polygon(round_corners(shape, radius = radii,$fn=64));
 233//   color("red") down(.1) polygon(shape);
 234// Example(2D,Med): Continuous curvature rounding, different at every corner, with varying smoothness parameters as well, and `$fs` set very small.  Note that `$fa` is ignored here with method set to "smooth".
 235//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 236//   cuts = [1.5,0,2,0.3, 1.2, 0];
 237//   k = [0.6, 0.5, 0.5, 0.7, 0.3, 0.5];
 238//   polygon(round_corners(shape, method="smooth", cut=cuts, k=k, $fs=0.1));
 239//   color("red") down(.1) polygon(shape);
 240// Example(2D,Med): Chamfers
 241//   $fn=36;
 242//   shape = [[0,0], [10,0], [15,12], [6,6], [6, 12], [-3,7]];
 243//   polygon(round_corners(shape, method="chamfer", cut=1));
 244//   color("red") down(.1) polygon(shape);
 245// Example(Med3D): 3D printing test pieces to display different curvature shapes.  You can see the discontinuity in the curvature on the "C" piece in the rendered image.
 246//   ten = square(50);
 247//   cut = 5;
 248//   linear_extrude(height=14) {
 249//     translate([25,25,0])text("C",size=30, valign="center", halign="center");
 250//     translate([85,25,0])text("5",size=30, valign="center", halign="center");
 251//     translate([85,85,0])text("3",size=30, valign="center", halign="center");
 252//     translate([25,85,0])text("7",size=30, valign="center", halign="center");
 253//   }
 254//   linear_extrude(height=13) {
 255//     polygon(round_corners(ten, cut=cut, $fn=96*4));
 256//     translate([60,0,0])polygon(round_corners(ten,  method="smooth", cut=cut, $fn=96));
 257//     translate([60,60,0])polygon(round_corners(ten, method="smooth", cut=cut, k=0.32, $fn=96));
 258//     translate([0,60,0])polygon(round_corners(ten, method="smooth", cut=cut, k=0.7, $fn=96));
 259//   }
 260// Example(2D,Med): Rounding a path that is not closed in a three different ways.
 261//   $fs=.1;
 262//   $fa=1;
 263//   zigzagx = [-10, 0, 10, 20, 29, 38, 46, 52, 59, 66, 72, 78, 83, 88, 92, 96, 99, 102, 112];
 264//   zigzagy = concat([0], flatten(repeat([-10,10],8)), [-10,0]);
 265//   zig = hstack(zigzagx,zigzagy);
 266//   stroke(zig,width=1);   // Original shape
 267//   fwd(20)            // Smooth size corners with a cut of 4 and curvature parameter 0.6
 268//     stroke(round_corners(zig,cut=4, k=0.6, method="smooth", closed=false),width=1);
 269//   fwd(40)            // Smooth size corners with circular arcs and a cut of 4
 270//     stroke(round_corners(zig,cut=4,closed=false, method="circle"),width=1);
 271//                      // Smooth size corners with a circular arc and radius 1.5 (close to maximum possible)
 272//   fwd(60)            // Note how the different points are cut back by different amounts
 273//     stroke(round_corners(zig,radius=1.5,closed=false),width=1);
 274// Example(FlatSpin,VPD=42,VPT=[7.75,6.69,5.22]): Rounding some random 3D paths
 275//   $fn=36;
 276//   list1= [
 277//     [2.887360, 4.03497, 6.372090],
 278//     [5.682210, 9.37103, 0.783548],
 279//     [7.808460, 4.39414, 1.843770],
 280//     [0.941085, 5.30548, 4.467530],
 281//     [1.860540, 9.81574, 6.497530],
 282//     [6.938180, 7.21163, 5.794530]
 283//   ];
 284//   list2= [
 285//     [1.079070, 4.74091, 6.900390],
 286//     [8.775850, 4.42248, 6.651850],
 287//     [5.947140, 9.17137, 6.156420],
 288//     [0.662660, 6.95630, 5.884230],
 289//     [6.564540, 8.86334, 9.953110],
 290//     [5.420150, 4.91874, 3.866960]
 291//   ];
 292//   path_sweep(regular_ngon(n=36,or=.1),round_corners(list1,closed=false, method="smooth", cut = 0.65));
 293//   right(6)
 294//     path_sweep(regular_ngon(n=36,or=.1),round_corners(list2,closed=false, method="circle", cut = 0.75));
 295// Example(3D,Med):  Rounding a spiral with increased rounding along the length
 296//   // Construct a square spiral path in 3D
 297//   $fn=36;
 298//   square = [[0,0],[1,0],[1,1],[0,1]];
 299//   spiral = flatten(repeat(concat(square,reverse(square)),5));  // Squares repeat 10x, forward and backward
 300//   squareind = [for(i=[0:9]) each [i,i,i,i]];                   // Index of the square for each point
 301//   z = count(40)*.2+squareind;
 302//   path3d = hstack(spiral,z);                                   // 3D spiral
 303//   rounding = squareind/20;
 304//       // Setting k=1 means curvature won't be continuous, but curves are as round as possible
 305//       // Try changing the value to see the effect.
 306//   rpath = round_corners(path3d, joint=rounding, k=1, method="smooth", closed=false);
 307//   path_sweep( regular_ngon(n=36, or=.1), rpath);
 308// Example(2D): The rounding invocation that is commented out gives an error because the rounding parameters interfere with each other.  The error message gives a list of factors that can help you fix this: [0.852094, 0.852094, 1.85457, 10.1529]
 309//   $fn=64;
 310//   path = [[0, 0],[10, 0],[20, 20],[30, -10]];
 311//   debug_polygon(path);
 312//   //polygon(round_corners(path,cut = [1,3,1,1],
 313//   //        method="circle"));
 314// Example(2D): The list of factors shows that the problem is in the first two rounding values, because the factors are smaller than one.  If we multiply the first two parameters by 0.85 then the roundings fit.  The verbose option gives us the same fit factors.  
 315//   $fn=64;
 316//   path = [[0, 0],[10, 0],[20, 20],[30, -10]];
 317//   polygon(round_corners(path,cut = [0.85,3*0.85,1,1],
 318//                         method="circle", verbose=true));
 319// Example(2D): From the fit factors we can see that rounding at vertices 2 and 3 could be increased a lot.  Applying those factors we get this more rounded shape.  The new fit factors show that we can still further increase the rounding parameters if we wish.  
 320//   $fn=64;
 321//   path = [[0, 0],[10, 0],[20, 20],[30, -10]];
 322//   polygon(round_corners(path,cut = [0.85,3*0.85,2.13, 10.15],
 323//                         method="circle",verbose=true));
 324// Example(2D): Using the `joint` parameter it's easier to understand whether your roundvers will fit.  We can guarantee a fairly large roundover on any path by picking each one to use up half the segment distance along the shorter of its two segments:
 325//   $fn=64;
 326//   path = [[0, 0],[10, 0],[20, 20],[30, -10]];
 327//   path_len = path_segment_lengths(path,closed=true);
 328//   halflen = [for(i=idx(path)) min(select(path_len,i-1,i))/2];
 329//   polygon(round_corners(path,joint = halflen,
 330//                         method="circle",verbose=true));
 331// Example(2D): Chamfering, specifying the chamfer width
 332//   path = star(5, step=2, d=100);
 333//   path2 = round_corners(path, method="chamfer", width=5);
 334//   polygon(path2);
 335// Example(2D): Chamfering, specifying the cut
 336//   path = star(5, step=2, d=100);
 337//   path2 = round_corners(path, method="chamfer", cut=5);
 338//   polygon(path2);
 339// Example(2D): Chamfering, specifying joint length
 340//   path = star(5, step=2, d=100);
 341//   path2 = round_corners(path, method="chamfer", joint=5);
 342//   polygon(path2);
 343// Example(2D): Two passes to apply chamfers first, and then round the unchamfered corners.  Chamfers always add one point, so it's not hard to keep track of the vertices
 344//   $fn=32;
 345//   shape = square(10);
 346//   chamfered = round_corners(shape, method="chamfer",
 347//                             cut=[2,0,2,0]);
 348//   rounded = round_corners(chamfered, 
 349//              cut = [0, 0,  // 1st original vertex, chamfered
 350//                     1.5,   // 2nd original vertex
 351//                     0, 0,  // 3rd original vertex, chamfered
 352//                     2.5]); // 4th original vertex
 353//   polygon(rounded);
 354// Example(2D): Another example of mixing chamfers and roundings with two passes
 355//   path = star(5, step=2, d=100);
 356//   chamfcut = [for (i=[0:4]) each [7,0]];
 357//   radii = [for (i=[0:4]) each [0,0,10]];
 358//   path2=round_corners(
 359//           round_corners(path,
 360//                         method="chamfer",
 361//                         cut=chamfcut),
 362//           radius=radii);
 363//   stroke(path2, closed=true);
 364// Example(2D,Med,NoAxes): Specifying by corner index.  Use {{list_set()}} to construct the full chamfer cut list. 
 365//   path = star(47, ir=25, or=50);  // long path, lots of corners
 366//   chamfind = [8, 28, 60];         // But only want 3 chamfers
 367//   chamfcut = list_set([],chamfind,[10,13,15],minlen=len(path)-1);
 368//   rpath = round_corners(path, cut=chamfcut, method="chamfer");   
 369//   polygon(rpath);
 370// Example(2D,Med,NoAxes): Two-pass to chamfer and round by index.  Use {{repeat_entries()}} to correct for first pass chamfers.
 371//   $fn=32;
 372//   path = star(47, ir=32, or=65);  // long path, lots of corners
 373//   chamfind = [8, 28, 60];         // But only want 3 chamfers
 374//   roundind = [7,9,27,29,59,61];   // And 6 roundovers
 375//   chamfcut = list_set([],chamfind,[10,13,15],minlen=len(path)-1);
 376//   roundcut = list_set([],roundind,repeat(8,6),minlen=len(path)-1);
 377//   dups = list_set([], chamfind, repeat(2,len(chamfind)), dflt=1, minlen=len(path)-1);
 378//   rpath1 = round_corners(path, cut=chamfcut, method="chamfer");
 379//   rpath2 = round_corners(rpath1, cut=repeat_entries(roundcut,dups));
 380//   polygon(rpath2);
 381module round_corners(path, method="circle", radius, r, cut, joint, width, k, closed=true, verbose=false) {no_module();}
 382function round_corners(path, method="circle", radius, r, cut, joint, width, k, closed=true, verbose=false) =
 383    assert(in_list(method,["circle", "smooth", "chamfer"]), "method must be one of \"circle\", \"smooth\" or \"chamfer\"")
 384    let(
 385        default_k = 0.5,
 386        size=one_defined([radius, r, cut, joint, width], "radius,r,cut,joint,width"),
 387        path = force_path(path), 
 388        size_ok = is_num(size) || len(size)==len(path) || (!closed && len(size)==len(path)-2),
 389        k_ok = is_undef(k) || (method=="smooth" && (is_num(k) || len(k)==len(path) || (!closed && len(k)==len(path)-2))),
 390        measure = is_def(radius) ? "radius"
 391                : is_def(r) ? "radius"
 392                : is_def(cut) ? "cut" 
 393                : is_def(joint) ? "joint"
 394                : "width"
 395    )
 396    assert(is_path(path,[2,3]), "input path must be a 2d or 3d path")
 397    assert(len(path)>2,str("Path has length ",len(path),".  Length must be 3 or more."))
 398    assert(size_ok,str("Input ",measure," must be a number or list with length ",len(path), closed?"":str(" or ",len(path)-2)))
 399    assert(k_ok,method=="smooth" ? str("Input k must be a number or list with length ",len(path), closed?"":str(" or ",len(path)-2)) :
 400                                   "Input k is only allowed with method=\"smooth\"")
 401    assert(method=="circle" || measure!="radius", "radius parameter allowed only with method=\"circle\"")
 402    assert(method=="chamfer" || measure!="width", "width parameter  allowed only with method=\"chamfer\"")
 403    let(
 404        parm = is_num(size) ? repeat(size, len(path)) :
 405               len(size)<len(path) ? [0, each size, 0] :
 406                                     size,
 407        k = is_undef(k) ? repeat(default_k,len(path)) :
 408            is_num(k) ? repeat(k, len(path)) :
 409            len(k)<len(path) ? [0, each k, 0] :
 410                               k,
 411        badparm = [for(i=idx(parm)) if(parm[i]<0)i],
 412        badk = [for(i=idx(k)) if(k[i]<0 || k[i]>1)i]
 413     )
 414     assert(is_vector(parm) && badparm==[], str(measure," must be nonnegative"))
 415     assert(is_vector(k) && badk==[], "k parameter must be in the interval [0,1]")
 416     let(
 417        // dk is a list of parameters, where distance is the joint length to move away from the corner
 418        //     "smooth" method: [distance, curvature]
 419        //     "circle" method: [distance, radius]
 420        //     "chamfer" method: [distance]
 421        dk = [
 422              for(i=[0:1:len(path)-1])
 423                  let(
 424                      pathbit = select(path,i-1,i+1),
 425                      // This is the half-angle at the corner
 426                      angle = approx(pathbit[0],pathbit[1]) || approx(pathbit[1],pathbit[2]) ? undef
 427                            : vector_angle(select(path,i-1,i+1))/2
 428                  )
 429                  (!closed && (i==0 || i==len(path)-1))  ? [0] :          // Force zeros at ends for non-closed
 430                  parm[i]==0 ? [0]    : // If no rounding requested then don't try to compute parameters
 431                  assert(is_def(angle), str("Repeated point in path at index ",i," with nonzero rounding"))
 432                  assert(!approx(angle,0), closed && i==0 ? "Closing the path causes it to turn back on itself at the end" :
 433                                                            str("Path turns back on itself at index ",i," with nonzero rounding"))
 434                  (method=="chamfer" && measure=="joint")? [parm[i]] :
 435                  (method=="chamfer" && measure=="cut")  ? [parm[i]/cos(angle)] :
 436                  (method=="chamfer" && measure=="width") ? [parm[i]/sin(angle)/2] :
 437                  (method=="smooth" && measure=="joint") ? [parm[i],k[i]] :
 438                  (method=="smooth" && measure=="cut")   ? [8*parm[i]/cos(angle)/(1+4*k[i]),k[i]] :
 439                  (method=="circle" && measure=="radius")? [parm[i]/tan(angle), parm[i]] :
 440                  (method=="circle" && measure=="joint") ? [parm[i], parm[i]*tan(angle)] : 
 441                /*(method=="circle" && measure=="cut")*/   approx(angle,90) ? [INF] : 
 442                                                           let( circ_radius = parm[i] / (1/sin(angle) - 1))
 443                                                           [circ_radius/tan(angle), circ_radius],
 444        ],
 445        lengths = [for(i=[0:1:len(path)]) norm(select(path,i)-select(path,i-1))],
 446        scalefactors = [
 447            for(i=[0:1:len(path)-1])
 448                if (closed || (i!=0 && i!=len(path)-1))
 449                 min(
 450                    lengths[i]/(select(dk,i-1)[0]+dk[i][0]),
 451                    lengths[i+1]/(dk[i][0]+select(dk,i+1)[0])
 452                 )
 453        ],
 454        dummy = verbose ? echo("Roundover scale factors:",scalefactors) : 0
 455    )
 456    assert(min(scalefactors)>=1,str("Roundovers are too big for the path.  If you multitply them by this vector they should fit: ",scalefactors))
 457    // duplicates are introduced when roundings fully consume a segment, so remove them
 458    deduplicate([
 459        for(i=[0:1:len(path)-1]) each
 460            (dk[i][0] == 0)? [path[i]] :
 461            (method=="smooth")? _bezcorner(select(path,i-1,i+1), dk[i]) :
 462            (method=="chamfer") ? _chamfcorner(select(path,i-1,i+1), dk[i]) :
 463            _circlecorner(select(path,i-1,i+1), dk[i])
 464    ]);
 465
 466// Computes the continuous curvature control points for a corner when given as
 467// input three points in a list defining the corner.  The points must be
 468// equidistant from each other to produce the continuous curvature result.
 469// The output control points will include the 3 input points plus two
 470// interpolated points.
 471//
 472// k is the curvature parameter, ranging from 0 for very slow transition
 473// up to 1 for a sharp transition that doesn't have continuous curvature any more
 474function _smooth_bez_fill(points,k) = [
 475        points[0],
 476        lerp(points[1],points[0],k),
 477        points[1],
 478        lerp(points[1],points[2],k),
 479        points[2],
 480];
 481
 482// Computes the points of a continuous curvature roundover given as input
 483// the list of 3 points defining the corner and a parameter specification
 484//
 485// If parm is a scalar then it is treated as the curvature and the control
 486// points are calculated using _smooth_bez_fill.  Otherwise, parm is assumed
 487// to be a pair [d,k] where d is the length of the curve.  The length is
 488// calculated from the input point list and the control point list will not
 489// necessarily include points[0] or points[2] on its output.
 490//
 491// The number of points output is $fn if it is set.  Otherwise $fs is used
 492// to calculate the point count.
 493
 494function _bezcorner(points, parm) =
 495        let(
 496                P = is_list(parm)?
 497                        let(
 498                                d = parm[0],
 499                                k = parm[1],
 500                                prev = unit(points[0]-points[1]),
 501                                next = unit(points[2]-points[1])
 502                        ) [
 503                                points[1]+d*prev,
 504                                points[1]+k*d*prev,
 505                                points[1],
 506                                points[1]+k*d*next,
 507                                points[1]+d*next
 508                        ] : _smooth_bez_fill(points,parm),
 509                N = max(3,$fn>0 ?$fn : ceil(bezier_length(P)/$fs))
 510        )
 511        bezier_curve(P,N,endpoint=true);
 512
 513function _chamfcorner(points, parm) =
 514        let(
 515                d = parm[0],
 516                prev = unit(points[0]-points[1]),
 517                next = unit(points[2]-points[1])
 518          )
 519       [points[1]+prev*d, points[1]+next*d];
 520
 521function _circlecorner(points, parm) =
 522        let(
 523            angle = vector_angle(points)/2,
 524            d = parm[0],
 525            r = parm[1],
 526            prev = unit(points[0]-points[1]),
 527            next = unit(points[2]-points[1])
 528        )
 529        approx(angle,90) ? [points[1]+prev*d, points[1]+next*d] :
 530        let(
 531            center = r/sin(angle) * unit(prev+next)+points[1],
 532                    start = points[1]+prev*d,
 533                    end = points[1]+next*d
 534        )     // 90-angle is half the angle of the circular arc
 535        arc(max(3,ceil((90-angle)/180*segs(r))), cp=center, points=[start,end]);
 536
 537
 538// Used by offset_sweep and convex_offset_extrude.
 539// Produce edge profile curve from the edge specification
 540// z_dir is the direction multiplier (1 to build up, -1 to build down)
 541function _rounding_offsets(edgespec,z_dir=1) =
 542        let(
 543                edgetype = struct_val(edgespec, "type"),
 544                extra = struct_val(edgespec,"extra"),
 545                N = struct_val(edgespec, "steps"),
 546                r = struct_val(edgespec,"r"),
 547                cut = struct_val(edgespec,"cut"),
 548                k = struct_val(edgespec,"k"),
 549                radius = in_list(edgetype,["circle","teardrop"])
 550                            ? (is_def(cut) ? cut/(sqrt(2)-1) : r)
 551                         :edgetype=="chamfer"
 552                            ? (is_def(cut) ? sqrt(2)*cut : r)
 553                         : undef,
 554                chamf_angle = struct_val(edgespec, "angle"),
 555                cheight = struct_val(edgespec, "chamfer_height"),
 556                cwidth = struct_val(edgespec, "chamfer_width"),
 557                chamf_width = first_defined([!all_defined([cut,chamf_angle]) ? undef : cut/cos(chamf_angle),
 558                                             cwidth,
 559                                             !all_defined([cheight,chamf_angle]) ? undef : cheight*tan(chamf_angle)]),
 560                chamf_height = first_defined([
 561                                              !all_defined([cut,chamf_angle]) ? undef : cut/sin(chamf_angle),
 562                                              cheight,
 563                                              !all_defined([cwidth, chamf_angle]) ? undef : cwidth/tan(chamf_angle)]),
 564                joint = first_defined([
 565                        struct_val(edgespec,"joint"),
 566                        all_defined([cut,k]) ? 16*cut/sqrt(2)/(1+4*k) : undef
 567                ]),
 568                points = struct_val(edgespec, "points"),
 569                argsOK = in_list(edgetype,["circle","teardrop"])? is_def(radius) :
 570                        edgetype == "chamfer"? chamf_angle>0 && chamf_angle<90 && num_defined([chamf_height,chamf_width])==2 :
 571                        edgetype == "smooth"? num_defined([k,joint])==2 :
 572                        edgetype == "profile"? points[0]==[0,0] :
 573                        false
 574        )
 575        assert(argsOK,str("Invalid specification with type ",edgetype))
 576        let(
 577                offsets =
 578                        edgetype == "profile"? scale([-1,z_dir], p=list_tail(points)) :
 579                        edgetype == "chamfer"?  chamf_width==0 && chamf_height==0? [] : [[-chamf_width,z_dir*abs(chamf_height)]] :
 580                        edgetype == "teardrop"? (
 581                                radius==0? [] : concat(
 582                                        [for(i=[1:N]) [radius*(cos(i*45/N)-1),z_dir*abs(radius)* sin(i*45/N)]],
 583                                        [[-2*radius*(1-sqrt(2)/2), z_dir*abs(radius)]]
 584                                )
 585                        ) :
 586                        edgetype == "circle"? radius==0? [] : [for(i=[1:N]) [radius*(cos(i*90/N)-1), z_dir*abs(radius)*sin(i*90/N)]] :
 587                        /* smooth */ joint==0 ? [] :
 588                        list_tail(
 589                                _bezcorner([[0,0],[0,z_dir*abs(joint)],[-joint,z_dir*abs(joint)]], k, $fn=N+2)
 590                        )
 591        )
 592  
 593        quant(extra > 0? concat(offsets, [last(offsets)+[0,z_dir*extra]]) : offsets, 1/1024);
 594
 595
 596
 597// Function: smooth_path()
 598// Usage:
 599//   smoothed = smooth_path(path, [tangents], [size=|relsize=], [splinesteps=], [closed=], [uniform=]);
 600// Description:
 601//   Smooths the input path using a cubic spline.  Every segment of the path will be replaced by a cubic curve
 602//   with `splinesteps` points.  The cubic interpolation will pass through every input point on the path
 603//   and will match the tangents at every point.  If you do not specify tangents they will be computed using
 604//   path_tangents with uniform=false by default.  Note that setting uniform to true with non-uniform
 605//   sampling may be desirable in some cases but tends to produces curves that overshoot the point on the path.  
 606//   .
 607//   The size or relsize parameter determines how far the curve can bend away from
 608//   the input path.  In the case where the curve has a single hump, the size specifies the exact distance
 609//   between the specified path and the curve.  If you give relsize then it is relative to the segment
 610//   length (e.g. 0.05 means 5% of the segment length).  In 2d when the spline may make an S-curve,
 611//   in which case the size parameter specifies the sum of the deviations of the two peaks of the curve.  In 3-space
 612//   the bezier curve may have three extrema: two maxima and one minimum.  In this case the size specifies
 613//   the sum of the maxima minus the minimum.  At a given segment there is a maximum size: if your size
 614//   value is too large it will be rounded down.  See also path_to_bezpath().
 615// Arguments:
 616//   path = path to smooth
 617//   tangents = tangents constraining curve direction at each point.  Default: computed automatically
 618//   ---
 619//   relsize = relative size specification for the curve, a number or vector.  Default: 0.1
 620//   size = absolute size specification for the curve, a number or vector
 621//   uniform = set to true to compute tangents with uniform=true.  Default: false
 622//   closed = true if the curve is closed.  Default: false. 
 623// Example(2D): Original path in green, smoothed path in yellow:
 624//   color("green")stroke(square(4), width=0.1);
 625//   stroke(smooth_path(square(4),size=0.4), width=0.1);
 626// Example(2D): Closing the path changes the end tangents
 627//   polygon(smooth_path(square(4),size=0.4,closed=true));
 628// Example(2D): Turning on uniform tangent calculation also changes the end derivatives:
 629//   color("green")stroke(square(4), width=0.1);
 630//   stroke(smooth_path(square(4),size=0.4,uniform=true),
 631//          width=0.1);
 632// Example(2D): Here's a wide rectangle.  Using size means all edges bulge the same amount, regardless of their length. 
 633//   color("green")
 634//     stroke(square([10,4]), closed=true, width=0.1);
 635//   stroke(smooth_path(square([10,4]),size=1,closed=true),
 636//          width=0.1);
 637// Example(2D): With relsize the bulge is proportional to the side length. 
 638//   color("green")stroke(square([10,4]), closed=true, width=0.1);
 639//   stroke(smooth_path(square([10,4]),relsize=0.1,closed=true),
 640//          width=0.1);
 641// Example(2D): Settting uniform to true biases the tangents to aline more with the line sides
 642//   color("green")
 643//     stroke(square([10,4]), closed=true, width=0.1);
 644//   stroke(smooth_path(square([10,4]),uniform=true,
 645//                      relsize=0.1,closed=true),
 646//          width=0.1);
 647// Example(2D): A more interesting shape:
 648//   path = [[0,0], [4,0], [7,14], [-3,12]];
 649//   polygon(smooth_path(path,size=1,closed=true));
 650// Example(2D): Here's the square again with less smoothing.
 651//   polygon(smooth_path(square(4), size=.25,closed=true));
 652// Example(2D): Here's the square with a size that's too big to achieve, so you get the maximum possible curve:
 653//   color("green")stroke(square(4), width=0.1,closed=true);
 654//   stroke(smooth_path(square(4), size=4, closed=true),
 655//          closed=true,width=.1);
 656// Example(2D): You can alter the shape of the curve by specifying your own arbitrary tangent values
 657//   polygon(smooth_path(square(4),
 658//           tangents=1.25*[[-2,-1], [-4,1], [1,2], [6,-1]],
 659//           size=0.4,closed=true));
 660// Example(2D): Or you can give a different size for each segment
 661//   polygon(smooth_path(square(4),size = [.4, .05, 1, .3],
 662//                       closed=true));
 663// Example(FlatSpin,VPD=35,VPT=[4.5,4.5,1]):  Works on 3d paths as well
 664//   path = [[0,0,0],[3,3,2],[6,0,1],[9,9,0]];
 665//   stroke(smooth_path(path,relsize=.1),width=.3);
 666// Example(2D): This shows the type of overshoot that can occur with uniform=true.  You can produce overshoots like this if you supply a tangent that is difficult to connect to the adjacent points  
 667//   pts = [[-3.3, 1.7], [-3.7, -2.2], [3.8, -4.8], [-0.9, -2.4]];
 668//   stroke(smooth_path(pts, uniform=true, relsize=0.1),width=.1);
 669//   color("red")move_copies(pts)circle(r=.15,$fn=12);
 670// Example(2D): With the default of uniform false no overshoot occurs.  Note that the shape of the curve is quite different.  
 671//   pts = [[-3.3, 1.7], [-3.7, -2.2], [3.8, -4.8], [-0.9, -2.4]];
 672//   stroke(smooth_path(pts, uniform=false, relsize=0.1),width=.1);
 673//   color("red")move_copies(pts)circle(r=.15,$fn=12);
 674module smooth_path(path, tangents, size, relsize, splinesteps=10, uniform=false, closed=false) {no_module();}
 675function smooth_path(path, tangents, size, relsize, splinesteps=10, uniform=false, closed) =
 676  is_1region(path) ? smooth_path(path[0], tangents, size, relsize, splinesteps, uniform, default(closed,true)) :
 677  let (
 678     bez = path_to_bezpath(path, tangents=tangents, size=size, relsize=relsize, uniform=uniform, closed=default(closed,false)),
 679     smoothed = bezpath_curve(bez,splinesteps=splinesteps)
 680  )
 681  closed ? cleanup_path(smoothed) : smoothed;
 682
 683
 684function _scalar_to_vector(value,length,varname) = 
 685  is_vector(value)
 686    ? assert(len(value)==length, str(varname," must be length ",length))
 687      value
 688    : assert(is_num(value), str(varname, " must be a numerical value"))
 689      repeat(value, length);
 690
 691
 692// Function: path_join()
 693// Usage:
 694//   joined_path = path_join(paths, [joint], [k=], [relocate=], [closed=]);
 695// Description:
 696//   Connect a sequence of paths together into a single path with optional continuous curvature rounding
 697//   applied at the joints.  By default the first path is taken as specified and subsequent paths are
 698//   translated into position so that each path starts where the previous path ended.
 699//   If you set relocate to false then this relocation is skipped.
 700//   You specify rounding using the `joint` parameter, which specifies the distance away from the corner
 701//   where the roundover should start.  The path_join function may remove many path points to cut the path 
 702//   back by the joint length.  Rounding is using continous curvature 4th order bezier splines and
 703//   the parameter `k` specifies how smooth the curvature match is.  This parameter ranges from 0 to 1 with
 704//   a default of 0.5.  Use a larger k value to get a curve that is bigger for the same joint value.  When
 705//   k=1 the curve may be similar to a circle if your curves are symmetric.  As the path is built up, the joint
 706//   parameter applies to the growing path, so if you pick a large joint parameter it may interact with the
 707//   previous path sections.  See [Types of Roundover](rounding.scad#subsection-types-of-roundover) for more details
 708//   on continuous curvature rounding. 
 709//   .
 710//   The rounding is created by extending the two clipped paths to define a corner point.  If the extensions of
 711//   the paths do not intersect, the function issues an error.  When closed=true the final path should actually close
 712//   the shape, repeating the starting point of the shape.  If it does not, then the rounding will fill the gap.
 713//   .
 714//   The number of segments in the roundovers is set based on $fn and $fs.  If you use $fn it specifies the number of
 715//   segments in the roundover, regardless of its angular extent.
 716// Arguments:
 717//   paths = list of paths to join
 718//   joint = joint distance, either a number, a pair (giving the previous and next joint distance) or a list of numbers and pairs.  Default: 0
 719//   ---
 720//   k = curvature parameter, either a number or vector.  Default: 0.5
 721//   relocate = set to false to prevent paths from being arranged tail to head.  Default: true
 722//   closed = set to true to round the junction between the last and first paths.  Default: false
 723// Example(2D): Connection of 3 simple paths.  
 724//   horiz = [[0,0],[10,0]];
 725//   vert = [[0,0],[0,10]];
 726//   stroke(path_join([horiz, vert, -horiz]));
 727// Example(2D): Adding curvature with joint of 3
 728//   horiz = [[0,0],[10,0]];
 729//   vert = [[0,0],[0,10]];
 730//   stroke(path_join([horiz, vert, -horiz],joint=3,$fn=16));
 731// Example(2D): Setting k=1 increases the amount of curvature
 732//   horiz = [[0,0],[10,0]];
 733//   vert = [[0,0],[0,10]];
 734//   stroke(path_join([horiz, vert, -horiz],joint=3,k=1,$fn=16));
 735// Example(2D): Specifying pairs of joint values at a path joint creates an asymmetric curve
 736//   horiz = [[0,0],[10,0]];
 737//   vert = [[0,0],[0,10]];
 738//   stroke(path_join([horiz, vert, -horiz],
 739//                    joint=[[4,1],[1,4]],$fn=16),width=.3);
 740// Example(2D): A closed square
 741//   horiz = [[0,0],[10,0]];
 742//   vert = [[0,0],[0,10]];
 743//   stroke(path_join([horiz, vert, -horiz, -vert],
 744//                    joint=3,k=1,closed=true,$fn=16),closed=true);
 745// Example(2D): Different curve at each corner by changing the joint size
 746//   horiz = [[0,0],[10,0]];
 747//   vert = [[0,0],[0,10]];
 748//   stroke(path_join([horiz, vert, -horiz, -vert],
 749//                    joint=[3,0,1,2],k=1,closed=true,$fn=16),
 750//          closed=true,width=0.4);
 751// Example(2D): Different curve at each corner by changing the curvature parameter.  Note that k=0 still gives a small curve, unlike joint=0 which gives a sharp corner.
 752//   horiz = [[0,0],[10,0]];
 753//   vert = [[0,0],[0,10]];
 754//   stroke(path_join([horiz, vert, -horiz, -vert],joint=3,
 755//                    k=[1,.5,0,.7],closed=true,$fn=16),
 756//          closed=true,width=0.4);
 757// Example(2D): Joint value of 7 is larger than half the square so curves interfere with each other, which breaks symmetry because they are computed sequentially
 758//   horiz = [[0,0],[10,0]];
 759//   vert = [[0,0],[0,10]];
 760//   stroke(path_join([horiz, vert, -horiz, -vert],joint=7,
 761//                     k=.4,closed=true,$fn=16),
 762//          closed=true);
 763// Example(2D): Unlike round_corners, we can add curves onto curves.
 764//   $fn=64;
 765//   myarc = arc(width=20, thickness=5 );
 766//   stroke(path_join(repeat(myarc,3), joint=4));
 767// Example(2D): Here we make a closed shape from two arcs and round the sharp tips
 768//   arc1 = arc(width=20, thickness=4,$fn=75);
 769//   arc2 = reverse(arc(width=20, thickness=2,$fn=75));
 770//   // Without rounding
 771//   stroke(path_join([arc1,arc2]),width=.3);
 772//   // With rounding
 773//   color("red")stroke(path_join([arc1,arc2], 3,k=1,closed=true),
 774//                      width=.3,closed=true,$fn=12); 
 775// Example(2D): Combining arcs with segments
 776//   arc1 = arc(width=20, thickness=4,$fn=75);
 777//   arc2 = reverse(arc(width=20, thickness=2,$fn=75));
 778//   vpath = [[0,0],[0,-5]];
 779//   stroke(path_join([arc1,vpath,arc2,reverse(vpath)]),width=.2);
 780//   color("red")stroke(path_join([arc1,vpath,arc2,reverse(vpath)],
 781//                                [1,2,2,1],k=1,closed=true),
 782//                      width=.2,closed=true,$fn=12);
 783// Example(2D): Here relocation is off.  We have three segments (in yellow) and add the curves to the segments.  Notice that joint zero still produces a curve because it refers to the endpoints of the supplied paths.  
 784//   p1 = [[0,0],[2,0]];
 785//   p2 = [[3,1],[1,3]];
 786//   p3 = [[0,3],[-1,1]];
 787//   color("red")stroke(
 788//     path_join([p1,p2,p3], joint=0, relocate=false,
 789//               closed=true),
 790//     width=.3,$fn=12);
 791//   for(x=[p1,p2,p3]) stroke(x,width=.3);
 792// Example(2D): If you specify closed=true when the last path doesn't meet the first one then it is similar to using relocate=false: the function tries to close the path using a curve.  In the example below, this results in a long curve to the left, when given the unclosed three segments as input.  Note that if the segments are parallel the function fails with an error.  The extension of the curves must intersect in a corner for the rounding to be well-defined.  To get a normal rounding of the closed shape, you must include a fourth path, the last segment that closes the shape.
 793//   horiz = [[0,0],[10,0]];
 794//   vert = [[0,0],[0,10]];
 795//   h2 = [[0,-3],[10,0]];
 796//   color("red")stroke(
 797//     path_join([horiz, vert, -h2],closed=true,
 798//               joint=3,$fn=25),
 799//     closed=true,width=.5);
 800//   stroke(path_join([horiz, vert, -h2]),width=.3);
 801// Example(2D): With a single path with closed=true the start and end junction is rounded.
 802//   tri = regular_ngon(n=3, r=7);
 803//   stroke(path_join([tri], joint=3,closed=true,$fn=12),
 804//          closed=true,width=.5);
 805module path_join(paths,joint=0,k=0.5,relocate=true,closed=false) { no_module();}
 806function path_join(paths,joint=0,k=0.5,relocate=true,closed=false)=
 807  assert(is_list(paths),"Input paths must be a list of paths")
 808  let(
 809      paths = [for(i=idx(paths)) force_path(paths[i],str("paths[",i,"]"))],
 810      badpath = [for(j=idx(paths)) if (!is_path(paths[j])) j]
 811  )
 812  assert(badpath==[], str("Entries in paths are not valid paths: ",badpath))
 813  len(paths)==0 ? [] :
 814  len(paths)==1 && !closed ? paths[0] :
 815  let(
 816      paths = !closed || len(paths)>1
 817            ? paths
 818            : [close_path(paths[0])],
 819      N = len(paths) + (closed?0:-1),
 820      k = _scalar_to_vector(k,N),
 821      repjoint = is_num(joint) || (is_vector(joint,2) && len(paths)!=3),
 822      joint = repjoint ? repeat(joint,N) : joint
 823  )
 824  assert(all_nonnegative(k), "k must be nonnegative")
 825  assert(len(joint)==N,str("Input joint must be scalar or length ",N))
 826  let(
 827      bad_j = [for(j=idx(joint)) if (!is_num(joint[j]) && !is_vector(joint[j],2)) j]
 828  )
 829  assert(bad_j==[], str("Invalid joint values at indices ",bad_j))
 830  let(result=_path_join(paths,joint,k, relocate=relocate, closed=closed))
 831  closed ? cleanup_path(result) : result;
 832
 833function _path_join(paths,joint,k=0.5,i=0,result=[],relocate=true,closed=false) =
 834  let( 
 835      result = result==[] ? paths[0] : result,
 836      loop = i==len(paths)-1,
 837      revresult = reverse(result),
 838      nextpath = loop     ? result
 839               : relocate ? move(revresult[0]-paths[i+1][0], p=paths[i+1])
 840               : paths[i+1],
 841      d_first = is_vector(joint[i]) ? joint[i][0] : joint[i],
 842      d_next = is_vector(joint[i]) ? joint[i][1] : joint[i]
 843  )
 844  assert(d_first>=0 && d_next>=0, str("Joint value negative when adding path ",i+1))
 845  assert(d_first<path_length(revresult),str("Path ",i," is too short for specified cut distance ",d_first))
 846  assert(d_next<path_length(nextpath), str("Path ",i+1," is too short for specified cut distance ",d_next))
 847  let(
 848      firstcut = path_cut_points(revresult, d_first, direction=true),
 849      nextcut = path_cut_points(nextpath, d_next, direction=true)
 850  )
 851  assert(!loop || nextcut[1] < len(revresult)-1-firstcut[1], "Path is too short to close the loop")
 852  let(
 853     first_dir=firstcut[2],
 854     next_dir=nextcut[2],
 855     corner = line_intersection([firstcut[0], firstcut[0]-first_dir], [nextcut[0], nextcut[0]-next_dir],RAY,RAY)
 856  )
 857  assert(is_def(corner), str("Curve directions at cut points don't intersect in a corner when ",
 858                             loop?"closing the path":str("adding path ",i+1)))
 859  let(
 860      bezpts = _smooth_bez_fill([firstcut[0], corner, nextcut[0]],k[i]),
 861      N = max(3,$fn>0 ?$fn : ceil(bezier_length(bezpts)/$fs)),
 862      bezpath = approx(firstcut[0],corner) && approx(corner,nextcut[0])
 863                  ? []
 864                  : bezier_curve(bezpts,N),
 865      new_result = [each select(result,loop?nextcut[1]:0,len(revresult)-1-firstcut[1]),
 866                    each bezpath,
 867                    nextcut[0],
 868                    if (!loop) each list_tail(nextpath,nextcut[1])
 869                   ]
 870  )
 871  i==len(paths)-(closed?1:2)
 872     ? new_result
 873     : _path_join(paths,joint,k,i+1,new_result, relocate,closed);
 874
 875
 876
 877// Function&Module: offset_stroke()
 878// Usage: as module
 879//   offset_stroke(path, [width], [rounded=], [chamfer=], [start=], [end=], [check_valid=], [quality=], [closed=],...) [ATTACHMENTS];
 880// Usage: as function
 881//   path = offset_stroke(path, [width], closed=false, [rounded=], [chamfer=], [start=], [end=], [check_valid=], [quality=],...);
 882//   region = offset_stroke(path, [width], closed=true, [rounded=], [chamfer=], [start=], [end=], [check_valid=], [quality=],...);
 883// Description:
 884//   Uses `offset()` to compute a stroke for the input path.  Unlike `stroke`, the result does not need to be
 885//   centered on the input path.  The corners can be rounded, pointed, or chamfered, and you can make the ends
 886//   rounded, flat or pointed with the `start` and `end` parameters.
 887//   .
 888//   The `check_valid` and `quality`  parameters are passed through to `offset()`
 889//   .
 890//   If `width` is a scalar then the output will be a centered stroke of the specified width.  If width
 891//   is a list of two values then those two values will define the stroke side positions relative to the center line, where
 892//   as with offset(), the shift is to the left for open paths and outward for closed paths.  For example,
 893//   setting `width` to `[0,1]` will create a stroke of width 1 that extends entirely to the left of the input, and and [-4,-6]
 894//   will create a stroke of width 2 offset 4 units to the right of the input path.
 895//   .
 896//   If closed==false then the function form will return a path.  If closed==true then it will return a region.  The `start` and
 897//   `end` parameters are forbidden for closed paths.
 898//   .
 899//   Three simple end treatments are supported, "flat" (the default), "round" and "pointed".  The "flat" treatment
 900//   cuts off the ends perpendicular to the path and the "round" treatment applies a semicircle to the end.  The
 901//   "pointed" end treatment caps the stroke with a centered triangle that has 45 degree angles on each side.
 902//   .
 903//   More complex end treatments are available through parameter lists with helper functions to ease parameter passing.  The parameter list
 904//   keywords are
 905//      - "for" : must appear first in the list and have the value "offset_stroke"
 906//      - "type": the type of end treatment, one of "shifted_point", "roundover", or "flat"
 907//      - "angle": relative angle (relative to the path)
 908//      - "abs_angle": absolute angle (angle relative to x-axis)
 909//      - "cut": cut distance for roundovers, a single value to round both corners identically or a list of two values for the two corners.  Negative values round outward.
 910//      - "k": curvature smoothness parameter for roundovers, default 0.75
 911//   .
 912//   Function helpers for defining ends, prefixed by "os" for offset_stroke, are:
 913//      - os_flat(angle|absangle): specify a flat end either relative to the path or relative to the x-axis
 914//      - os_pointed(dist, [loc]): specify a pointed tip where the point is distance `loc` from the centerline (positive is the left direction as for offset), and `dist` is the distance from the path end to the point tip.  The default value for `loc` is zero (the center).  You must specify `dist` when using this option.
 915//      - os_round(cut, [angle|absangle], [k]).  Rounded ends with the specified cut distance, based on the specified angle or absolute angle.  The `k` parameter is the smoothness parameter for continuous curvature rounding.  See [Types of Roundover](rounding.scad#subsection-types-of-roundover) for more details on
 916//        continuous curvature rounding.  
 917//   .
 918//   Note that `offset_stroke()` will attempt to apply roundovers and angles at the ends even when it means deleting segments of the stroke, unlike round_corners which only works on a segment adjacent to a corner.  If you specify an overly extreme angle it will fail to find an intersection with the stroke and display an error.  When you specify an angle the end segment is rotated around the center of the stroke and the last segment of the stroke one one side is extended to the corner.
 919//   .
 920//   The `$fn` and `$fs` variables are used in the usual way to determine the number of segments for roundings produced by the offset
 921//   invocations and roundings produced by the semi-circular "round" end treatment.  The os_round() end treatment
 922//   uses a bezier curve, and will produce segments of approximate length `$fs` or it will produce `$fn` segments.
 923//   (This means that even a quarter circle will have `$fn` segments, unlike the usual case where it would have `$fn/4` segments.)
 924// Arguments:
 925//   path = 2d path that defines the stroke
 926//   width = width of the stroke, a scalar or a vector of 2 values giving the offset from the path.  Default: 1
 927//   ---
 928//   rounded = set to true to use rounded offsets, false to use sharp (delta) offsets.  Default: true
 929//   chamfer = set to true to use chamfers when `rounded=false`.  Default: false
 930//   start = end treatment for the start of the stroke when closed=false.  See above for details.  Default: "flat"
 931//   end = end treatment for the end of the stroke when closed=false.  See above for details.  Default: "flat"
 932//   check_valid = passed to offset().  Default: true
 933//   quality = passed to offset().  Default: 1
 934//   closed = true if the curve is closed, false otherwise.  Default: false
 935//   anchor = Translate so anchor point is at origin (0,0,0).  See [anchor](attachments.scad#subsection-anchor).  Default: `"origin"`
 936//   spin = Rotate this many degrees after anchor.  See [spin](attachments.scad#subsection-spin).  Default: `0`
 937//   cp = Centerpoint for determining intersection anchors or centering the shape.  Determintes the base of the anchor vector.  Can be "centroid", "mean", "box" or a 2D point.  Default: "centroid"
 938//   atype = Set to "hull" or "intersect" to select anchor type.  Default: "hull"
 939// Example(2D):  Basic examples illustrating flat, round, and pointed ends, on a finely sampled arc and a path made from 3 segments.
 940//   arc = arc(points=[[1,1],[3,4],[6,3]],n=50);
 941//   path = [[0,0],[6,2],[9,7],[8,10]];
 942//   xdistribute(spacing=10){
 943//     offset_stroke(path, width = 2);
 944//     offset_stroke(path, start="round", end="round", width = 2);
 945//     offset_stroke(path, start="pointed", end="pointed", width = 2);
 946//   }
 947//   fwd(10) xdistribute(spacing=10){
 948//     offset_stroke(arc, width = 2);
 949//     offset_stroke(arc, start="round", end="round", width = 2);
 950//     offset_stroke(arc, start="pointed", end="pointed", width = 2);
 951//   }
 952// Example(2D):  The effect of the `rounded` and `chamfer` options is most evident at sharp corners.  This only affects the middle of the path, not the ends.
 953//   sharppath = [[0,0], [1.5,5], [3,0]];
 954//   xdistribute(spacing=5){
 955//     offset_stroke(sharppath, $fn=16);
 956//     offset_stroke(sharppath, rounded=false);
 957//     offset_stroke(sharppath, rounded=false, chamfer=true);
 958//   }
 959// Example(2D):  When closed is enabled all the corners are affected by those options.
 960//   sharppath = [[0,0], [1.5,5], [3,0]];
 961//   xdistribute(spacing=5){
 962//     offset_stroke(sharppath,closed=true, $fn=16);
 963//     offset_stroke(sharppath, rounded=false, closed=true);
 964//     offset_stroke(sharppath, rounded=false, chamfer=true,
 965//                   closed=true);
 966//   }
 967// Example(2D):  The left stroke uses flat ends with a relative angle of zero.  The right hand one uses flat ends with an absolute angle of zero, so the ends are parallel to the x-axis.
 968//   path = [[0,0],[6,2],[9,7],[8,10]];
 969//   offset_stroke(path, start=os_flat(angle=0), end=os_flat(angle=0));
 970//   right(5)
 971//     offset_stroke(path, start=os_flat(abs_angle=0), end=os_flat(abs_angle=0));
 972// Example(2D):  With continuous sampling the end treatment can remove segments or extend the last segment linearly, as shown here.  Again the left side uses relative angle flat ends and the right hand example uses absolute angle.
 973//   arc = arc(points=[[4,0],[3,4],[6,3]],n=50);
 974//   offset_stroke(arc, start=os_flat(angle=45), end=os_flat(angle=45));
 975//   right(5)
 976//     offset_stroke(arc, start=os_flat(abs_angle=45), end=os_flat(abs_angle=45));
 977// Example(2D):  The os_pointed() end treatment allows adjustment of the point tip, as shown here.  The width is 2 so a location of 1 is at the edge.
 978//   arc = arc(points=[[1,1],[3,4],[6,3]],n=50);
 979//   offset_stroke(arc, width=2, start=os_pointed(loc=1,dist=3),end=os_pointed(loc=1,dist=3));
 980//   right(10)
 981//     offset_stroke(arc, width=2, start=os_pointed(dist=4),end=os_pointed(dist=-1));
 982//   fwd(7)
 983//     offset_stroke(arc, width=2, start=os_pointed(loc=2,dist=2),end=os_pointed(loc=.5,dist=-1));
 984// Example(2D):  The os_round() end treatment adds roundovers to the end corners by specifying the `cut` parameter.  In the first example, the cut parameter is the same at each corner.  The bezier smoothness parameter `k` is given to allow a larger cut.  In the second example, each corner is given a different roundover, including zero for no rounding at all.  The red shows the same strokes without the roundover.
 985//   $fn=36;
 986//   arc = arc(points=[[1,1],[3,4],[6,3]],n=50);
 987//   path = [[0,0],[6,2],[9,7],[8,10]];
 988//   offset_stroke(path, width=2, rounded=false,start=os_round(angle=-20, cut=0.4,k=.9),
 989//                                              end=os_round(angle=-35, cut=0.4,k=.9));
 990//   color("red")down(.1)offset_stroke(path, width=2, rounded=false,start=os_flat(-20),
 991//                                                                  end=os_flat(-35));
 992//   right(9){
 993//     offset_stroke(arc, width=2, rounded=false, start=os_round(cut=[.3,.6],angle=-45),
 994//                                                end=os_round(angle=20,cut=[.6,0]));
 995//     color("red")down(.1)offset_stroke(arc, width=2, rounded=false, start=os_flat(-45),
 996//                                                                    end=os_flat(20));
 997//   }
 998// Example(2D):  Negative cut values produce a flaring end.  Note how the absolute angle aligns the ends of the first example withi the axes.  In the second example positive and negative cut values are combined.  Note also that very different cuts are needed at the start end to produce a similar looking flare.
 999//   arc = arc(points=[[1,1],[3,4],[6,3]],n=50);
1000//   path = [[0,0],[6,2],[9,7],[8,10]];
1001//   offset_stroke(path, width=2, rounded=false,start=os_round(cut=-1, abs_angle=90),
1002//                                              end=os_round(cut=-0.5, abs_angle=0),$fn=36);
1003//   right(10)
1004//      offset_stroke(arc, width=2, rounded=false, start=os_round(cut=[-.75,-.2], angle=-45),
1005//                                                 end=os_round(cut=[-.2,.2], angle=20),$fn=36);
1006// Example(2D):  Setting the width to a vector allows you to offset the stroke.  Here with successive increasing offsets we create a set of parallel strokes
1007//   path = [[0,0],[4,4],[8,4],[2,9],[10,10]];
1008//   for(i=[0:.25:2])
1009//     offset_stroke(path, rounded=false,width = [i,i+.08]);
1010// Example(2D):  Setting rounded=true in the above example makes a very big difference in the result.  
1011//   path = [[0,0],[4,4],[8,4],[2,9],[10,10]];
1012//   for(i=[0:.25:2])
1013//     offset_stroke(path, rounded=true,width = [i,i+.08], $fn=36);
1014// Example(2D):  In this example a spurious triangle appears.  This results from overly enthusiastic validity checking.  Turning validity checking off fixes it in this case.
1015//   path = [[0,0],[4,4],[8,4],[2,9],[10,10]];
1016//   offset_stroke(path, check_valid=true,rounded=false,
1017//                 width = [1.4, 1.5]);
1018//   right(2)
1019//     offset_stroke(path, check_valid=false,rounded=false,
1020//                   width = [1.4, 1.5]);
1021// Example(2D):  But in this case, disabling the validity check produces an invalid result.
1022//   path = [[0,0],[4,4],[8,4],[2,9],[10,10]];
1023//   offset_stroke(path, check_valid=true,rounded=false,
1024//                 width = [1.9, 2]);
1025//   translate([1,-0.25])
1026//     offset_stroke(path, check_valid=false,rounded=false,
1027//                   width = [1.9, 2]);
1028// Example(2D): Self-intersecting paths are handled differently than with the `stroke()` module.
1029//   $fn=16;
1030//   path = turtle(["move",10,"left",144], repeat=4);
1031//   stroke(path, closed=true);
1032//   right(12)
1033//     offset_stroke(path, width=1, closed=true);
1034function offset_stroke(path, width=1, rounded=true, start, end, check_valid=true, quality=1, chamfer=false, closed=false,
1035                       atype="hull", anchor, spin, cp="centroid") =
1036        let(path = force_path(path))
1037        assert(is_path(path,2),"path is not a 2d path")
1038        let(
1039            closedok = !closed || (is_undef(start) && is_undef(end)),
1040            start = default(start,"flat"),
1041            end = default(end,"flat")
1042        )
1043        assert(closedok, "Parameters `start` and `end` not allowed with closed path")
1044        let(
1045            start = closed? [] : _parse_stroke_end(default(start,"flat"),"start"),
1046            end = closed? [] : _parse_stroke_end(default(end,"flat"),"end"),
1047            width = is_list(width)? reverse(sort(width)) : [1,-1]*width/2,
1048            left_r = !rounded? undef : width[0],
1049            left_delta = rounded? undef : width[0],
1050            right_r = !rounded? undef : width[1],
1051            right_delta = rounded? undef : width[1],
1052            left_path = offset(
1053                    path, delta=left_delta, r=left_r, closed=closed,
1054                    check_valid=check_valid, quality=quality,
1055                    chamfer=chamfer 
1056            ),
1057            right_path = offset(
1058                    path, delta=right_delta, r=right_r, closed=closed,
1059                    check_valid=check_valid, quality=quality,
1060                    chamfer=chamfer 
1061            )
1062         )
1063         closed? let(pts = [left_path, right_path])
1064                 reorient(anchor=anchor, spin=spin, two_d=true, region=pts, extent=atype=="hull", cp=cp, p=pts)
1065         :
1066         let(
1067             startpath = _stroke_end(width,left_path, right_path, start),
1068             endpath = _stroke_end(reverse(width),reverse(right_path), reverse(left_path),end),
1069             clipping_ok = startpath[1]+endpath[2]<=len(left_path) && startpath[2]+endpath[1]<=len(right_path)
1070         )
1071         assert(clipping_ok, "End treatment removed the whole stroke")
1072         let(
1073             pts = concat(
1074                          slice(left_path,startpath[1],-1-endpath[2]),
1075                          endpath[0],
1076                          reverse(slice(right_path,startpath[2],-1-endpath[1])),
1077                          startpath[0]
1078                  )
1079         )
1080         reorient(anchor=anchor, spin=spin, two_d=true, path=pts, extent=atype=="hull", cp=cp, p=pts);
1081
1082function os_pointed(dist,loc=0) =
1083        assert(is_def(dist), "Must specify `dist`")
1084        [
1085                "for", "offset_stroke",
1086                "type", "shifted_point",
1087                "loc",loc,
1088                "dist",dist
1089        ];
1090
1091function os_round(cut, angle, abs_angle, k, r) =
1092        assert(is_undef(r), "Radius not supported for os_round with offset_stroke.  (Did you mean os_circle for offset_sweep?)")
1093        let(
1094                acount = num_defined([angle,abs_angle]),
1095                use_angle = first_defined([angle,abs_angle,0])
1096        )
1097        assert(acount<2, "You must define only one of `angle` and `abs_angle`")
1098        assert(is_def(cut), "Parameter `cut` not defined.")
1099        [
1100                "for", "offset_stroke",
1101                "type", "roundover",
1102                "angle", use_angle,
1103                "absolute", is_def(abs_angle),
1104                "cut", is_vector(cut)? point2d(cut) : [cut,cut],
1105                "k", first_defined([k, 0.75])
1106        ];
1107
1108
1109function os_flat(angle, abs_angle) =
1110        let(
1111                acount = num_defined([angle,abs_angle]),
1112                use_angle = first_defined([angle,abs_angle,0])
1113        )
1114        assert(acount<2, "You must define only one of `angle` and `abs_angle`")
1115        [
1116                "for", "offset_stroke",
1117                "type", "flat",
1118                "angle", use_angle,
1119                "absolute", is_def(abs_angle)
1120        ];
1121
1122
1123
1124// Return angle in (-90,90] required to map line1 onto line2 (lines specified as lists of two points)
1125function angle_between_lines(line1,line2) =
1126        let(angle = atan2(det2([line1,line2]),line1*line2))
1127        angle > 90 ? angle-180 :
1128        angle <= -90 ? angle+180 :
1129        angle;
1130
1131
1132function _parse_stroke_end(spec,name) =
1133        is_string(spec)?
1134            assert(
1135                    in_list(spec,["flat","round","pointed"]),
1136                    str("Unknown \"",name,"\" string specification \"", spec,"\".  Must be \"flat\", \"round\", or \"pointed\"")
1137            )
1138            [["type", spec]]
1139        : let(
1140              dummy = _struct_valid(spec,"offset_stroke",name)
1141          )
1142          struct_set([], spec);
1143
1144
1145function _stroke_end(width,left, right, spec) =
1146        let(
1147                type = struct_val(spec, "type"),
1148                user_angle = default(struct_val(spec, "angle"), 0),
1149                normal_seg = _normal_segment(right[0], left[0]),
1150                normal_pt = normal_seg[1],
1151                center = normal_seg[0],
1152                parallel_dir = unit(left[0]-right[0]),
1153                normal_dir = unit(normal_seg[1]-normal_seg[0]),
1154                width_dir = sign(width[0]-width[1])
1155        )
1156        type == "round"? [arc(points=[right[0],normal_pt,left[0]],n=ceil(segs(width/2)/2)),1,1]  :
1157        type == "pointed"? [[normal_pt],0,0] :
1158        type == "shifted_point"? (
1159                let(shiftedcenter = center + width_dir * parallel_dir * struct_val(spec, "loc"))
1160                [[shiftedcenter+normal_dir*struct_val(spec, "dist")],0,0]
1161        ) :
1162        // Remaining types all support angled cutoff, so compute that
1163        assert(abs(user_angle)<=90, "End angle must be in [-90,90]")
1164        let(
1165                angle = struct_val(spec,"absolute")?
1166                        angle_between_lines(left[0]-right[0],[cos(user_angle),sin(user_angle)]) :
1167                        user_angle,
1168                endseg = [center, rot(p=[left[0]], angle, cp=center)[0]],
1169                intright = angle>0,
1170                pathclip = _path_line_intersection(intright? right : left, endseg),
1171                pathextend = line_intersection(endseg, select(intright? left:right,0,1))
1172        )
1173        type == "flat"? (
1174                intright?
1175                        [[pathclip[0], pathextend], 1, pathclip[1]] :
1176                        [[pathextend, pathclip[0]], pathclip[1],1]
1177        ) :
1178        type == "roundover"? (
1179                let(
1180                        bez_k = struct_val(spec,"k"),
1181                        cut = struct_val(spec,"cut"),
1182                        cutleft = cut[0],
1183                        cutright = cut[1],
1184                        // Create updated paths taking into account clipping for end rotation
1185                        newright = intright?
1186                                concat([pathclip[0]],list_tail(right,pathclip[1])) :
1187                                concat([pathextend],list_tail(right)),
1188                        newleft = !intright?
1189                                concat([pathclip[0]],list_tail(left,pathclip[1])) :
1190                                concat([pathextend],list_tail(left)),
1191                        // calculate corner angles, which are different when the cut is negative (outside corner)
1192                        leftangle = cutleft>=0?
1193                                vector_angle([newleft[1],newleft[0],newright[0]])/2 :
1194                                90-vector_angle([newleft[1],newleft[0],newright[0]])/2,
1195                        rightangle = cutright>=0?
1196                                vector_angle([newright[1],newright[0],newleft[0]])/2 :
1197                                90-vector_angle([newright[1],newright[0],newleft[0]])/2,
1198                        jointleft = 8*cutleft/cos(leftangle)/(1+4*bez_k),
1199                        jointright = 8*cutright/cos(rightangle)/(1+4*bez_k),
1200                        pathcutleft = path_cut_points(newleft,abs(jointleft)),
1201                        pathcutright = path_cut_points(newright,abs(jointright)),
1202                        leftdelete = intright? pathcutleft[1] : pathcutleft[1] + pathclip[1] -1,
1203                        rightdelete = intright? pathcutright[1] + pathclip[1] -1 : pathcutright[1],
1204                        leftcorner = line_intersection([pathcutleft[0], newleft[pathcutleft[1]]], [newright[0],newleft[0]]),
1205                        rightcorner = line_intersection([pathcutright[0], newright[pathcutright[1]]], [newright[0],newleft[0]]),
1206                        roundover_fits = jointleft+jointright < norm(rightcorner-leftcorner)
1207                )
1208                assert(roundover_fits,"Roundover too large to fit")
1209                let(
1210                        angled_dir = unit(newleft[0]-newright[0]),
1211                        nPleft = [
1212                                leftcorner - jointleft*angled_dir,
1213                                leftcorner,
1214                                pathcutleft[0]
1215                        ],
1216                        nPright = [
1217                                pathcutright[0],
1218                                rightcorner,
1219                                rightcorner + jointright*angled_dir
1220                        ],
1221                        leftcurve = _bezcorner(nPleft, bez_k),
1222                        rightcurve = _bezcorner(nPright, bez_k)
1223                )
1224                [concat(rightcurve, leftcurve), leftdelete, rightdelete]
1225        ) : [[],0,0];  // This case shouldn't occur
1226
1227// returns [intersection_pt, index of first point in path after the intersection]
1228function _path_line_intersection(path, line, ind=0) =
1229        ind==len(path)-1 ? undef :
1230        let(intersect=line_intersection(line, select(path,ind,ind+1),LINE,SEGMENT))
1231        // If it intersects the segment excluding it's final point, then we're done
1232        // The final point is treated as part of the next segment
1233        is_def(intersect) && intersect != path[ind+1]?
1234                [intersect, ind+1] :
1235                _path_line_intersection(path, line, ind+1);
1236
1237module offset_stroke(path, width=1, rounded=true, start, end, check_valid=true, quality=1, chamfer=false, closed=false,
1238                     atype="hull", anchor, spin, cp="centroid")
1239{
1240        result = offset_stroke(
1241                path, width=width, rounded=rounded,
1242                start=start, end=end,
1243                check_valid=check_valid, quality=quality,
1244                chamfer=chamfer,
1245                closed=closed
1246        );
1247        region(result,atype=atype, anchor=anchor, spin=spin, cp=cp) children();
1248}
1249
1250
1251// Section: Three-Dimensional Rounding
1252
1253// Function&Module: offset_sweep()
1254// Usage: most common module arguments.  See Arguments list below for more.
1255//    offset_sweep(path, [height|length|h|l|], [bottom], [top], [offset=], [convexity=],...) [ATTACHMENTS];
1256// Usage: most common function arguments.  See Arguments list below for more.
1257//    vnf = offset_sweep(path, [height|h|l|length], [bottom], [top], [offset=], ...);
1258// Description:
1259//   Takes a 2d path as input and extrudes it upwards and/or downward.  Each layer in the extrusion is produced using `offset()` to expand or shrink the previous layer.  When invoked as a function returns a VNF; when invoked as a module produces geometry.  
1260//   Using the `top` and/or `bottom` arguments you can specify a sequence of offsets values, or you can use several built-in offset profiles that
1261//   provide end treatments such as roundovers.
1262//   The height of the resulting object can be specified using the `height` argument, in which case `height` must be larger than the combined height
1263//   of the end treatments.  If you omit `height` then the object height will be the height of just the top and bottom end treatments.  
1264//   .
1265//   The path is shifted by `offset()` multiple times in sequence
1266//   to produce the final shape (not multiple shifts from one parent), so coarse definition of the input path will degrade
1267//   from the successive shifts.  If the result seems rough or strange try increasing the number of points you use for
1268//   your input.  If you get unexpected corners in your result you may have forgotten to set `$fn` or `$fa` and `$fs`.  
1269//   Be aware that large numbers of points (especially when check_valid is true) can lead to lengthy run times.  If your
1270//   shape doesn't develop new corners from the offsetting you may be able to save a lot of time by setting `check_valid=false`.  Be aware that
1271//   disabling the validity check when it is needed can generate invalid polyhedra that will produce CGAL errors upon
1272//   rendering.  Such validity errors will also occur if you specify a self-intersecting shape.
1273//   The offset profile is quantized to 1/1024 steps to avoid failures in offset() that can occur with very tiny offsets.
1274//   .
1275//   The build-in profiles are: circular rounding, teardrop rounding, continuous curvature rounding, and chamfer.
1276//   Also note that when a rounding radius is negative the rounding will flare outwards.  The easiest way to specify
1277//   the profile is by using the profile helper functions.  These functions take profile parameters, as well as some
1278//   general settings and translate them into a profile specification, with error checking on your input.  The description below
1279//   describes the helper functions and the parameters specific to each function.  Below that is a description of the generic
1280//   settings that you can optionally use with all of the helper functions.  For more details on the "cut" and "joint" rounding parameters, and
1281//   on continuous curvature rounding, see [Types of Roundover](rounding.scad#subsection-types-of-roundover). 
1282//   .
1283//   - profile: os_profile(points)
1284//     Define the offset profile with a list of points.  The first point must be [0,0] and the roundover should rise in the positive y direction, with positive x values for inward motion (standard roundover) and negative x values for flaring outward.  If the y value ever decreases then you might create a self-intersecting polyhedron, which is invalid.  Such invalid polyhedra will create cryptic assertion errors when you render your model and it is your responsibility to avoid creating them.  Note that the starting point of the profile is the center of the extrusion.  If you use a profile as the top it will rise upwards.  If you use it as the bottom it will be inverted, and will go downward.
1285//   - circle: os_circle(r|cut).  Define circular rounding either by specifying the radius or cut distance.
1286//   - smooth: os_smooth(cut|joint, [k]).  Define continuous curvature rounding, with `cut` and `joint` as for round_corners. The k parameter controls how fast the curvature changes and should be between 0 and 1.  
1287//   - teardrop: os_teardrop(r|cut).  Rounding using a 1/8 circle that then changes to a 45 degree chamfer.  The chamfer is at the end, and enables the object to be 3d printed without support.  The radius gives the radius of the circular part.
1288//   - chamfer: os_chamfer([height], [width], [cut], [angle]).  Chamfer the edge at desired angle or with desired height and width.  You can specify height and width together and the angle will be ignored, or specify just one of height and width and the angle is used to determine the shape.  Alternatively, specify "cut" along with angle to specify the cut back distance of the chamfer.
1289//   - mask: os_mask(mask, [out]).  Create a profile from one of the [2d masking shapes](shapes2d.scad#5-2d-masking-shapes).  The `out` parameter specifies that the mask should flare outward (like crown molding or baseboard).  This is set false by default.  
1290//   .
1291//   The general settings that you can use with all of the helper functions are mostly used to control how offset_sweep() calls the offset() function.
1292//   - extra: Add an extra vertical step of the specified height, to be used for intersections or differences.  This extra step will extend the resulting object beyond the height you specify.  Default: 0
1293//   - check_valid: passed to offset().  Default: true
1294//   - quality: passed to offset().  Default: 1
1295//   - steps: Number of vertical steps to use for the profile.  (Not used by os_profile).  Default: 16
1296//   - offset: Select "round" (r=) or "delta" (delta=) offset types for offset. You can also choose "chamfer" but this leads to exponential growth in the number of vertices with the steps parameter.  Default: "round"
1297//   .
1298//   Many of the arguments are described as setting "default" values because they establish settings which may be overridden by
1299//   the top and bottom profile specifications.
1300//   .
1301//   You will generally want to use the above helper functions to generate the profiles.
1302//   The profile specification is a list of pairs of keywords and values, e.g. ["for","offset_sweep","r",12, type, "circle"]. The keywords are
1303//   - "for" - must appear first in the list and have the value "offset_sweep"
1304//   - "type" - type of rounding to apply, one of "circle", "teardrop", "chamfer", "smooth", or "profile" (Default: "circle")
1305//   - "r" - the radius of the roundover, which may be zero for no roundover, or negative to round or flare outward.  Default: 0
1306//   - "cut" - the cut distance for the roundover or chamfer, which may be negative for flares
1307//   - "chamfer_width" - the width of a chamfer
1308//   - "chamfer_height" - the height of a chamfer
1309//   - "angle" - the chamfer angle, measured from the vertical (so zero is vertical, 90 is horizontal).  Default: 45
1310//   - "joint" - the joint distance for a "smooth" roundover
1311//   - "k" - the curvature smoothness parameter for "smooth" roundovers, a value in [0,1].  Default: 0.75
1312//   - "points" - point list for use with the "profile" type
1313//   - "extra" - extra height added for unions/differences.  This makes the shape taller than the requested height.  (Default: 0)
1314//   - "check_valid" - passed to offset.  Default: true.
1315//   - "quality" - passed to offset.  Default: 1.
1316//   - "steps" - number of vertical steps to use for the roundover.  Default: 16.
1317//   - "offset" - select "round" (r=), "delta" (delta=), or "chamfer" offset type for offset.  Default: "round"
1318//   .
1319//   Note that if you set the "offset" parameter to "chamfer" then every exterior corner turns from one vertex into two vertices with
1320//   each offset operation.  Since the offsets are done one after another, each on the output of the previous one, this leads to
1321//   exponential growth in the number of vertices.  This can lead to long run times or yield models that
1322//   run out of recursion depth and give a cryptic error.  Furthermore, the generated vertices are distributed non-uniformly.  Generally you
1323//   will get a similar or better looking model with fewer vertices using "round" instead of
1324//   "chamfer".  Use the "chamfer" style offset only in cases where the number of steps is very small or just one (such as when using
1325//   the `os_chamfer` profile type).
1326//
1327// Arguments:
1328//   path = 2d path (list of points) to extrude
1329//   height / l / h = total height (including rounded portions, but not extra sections) of the output.  Default: combined height of top and bottom end treatments.
1330//   bottom = rounding spec for the bottom end
1331//   top = rounding spec for the top end.
1332//   ---
1333//   offset = default offset, `"round"` or `"delta"`.  Default: `"round"`
1334//   steps = default step count.  Default: 16
1335//   quality = default quality.  Default: 1
1336//   check_valid = default check_valid.  Default: true.
1337//   extra = default extra height.  Default: 0
1338//   cut = default cut value.
1339//   chamfer_width = default width value for chamfers.
1340//   chamfer_height = default height value for chamfers.
1341//   angle = default angle for chamfers.  Default: 45
1342//   joint = default joint value for smooth roundover.
1343//   k = default curvature parameter value for "smooth" roundover
1344//   convexity = convexity setting for use with polyhedron.  (module only) Default: 10
1345//   anchor = Translate so anchor point is at the origin.  (module only) Default: "origin"
1346//   spin = Rotate this many degrees around Z axis after anchor.  (module only) Default: 0
1347//   orient = Vector to rotate top towards after spin  (module only)
1348//   atype = Select "hull" or "intersect" anchor types.  Default: "hull"
1349//   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"
1350// Example: Rounding a star shaped prism with postive radius values
1351//   star = star(5, r=22, ir=13);
1352//   rounded_star = round_corners(star, cut=flatten(repeat([.5,0],5)), $fn=24);
1353//   offset_sweep(rounded_star, height=20, bottom=os_circle(r=4), top=os_circle(r=1), steps=15);
1354// Example: Rounding a star shaped prism with negative radius values.  The starting shape has no corners, so the value of `$fn` does not matter.
1355//   star = star(5, r=22, ir=13); 
1356//   rounded_star = round_corners(star, cut=flatten(repeat([.5,0],5)), $fn=36);
1357//   offset_sweep(rounded_star, height=20, bottom=os_circle(r=-4), top=os_circle(r=-1), steps=15);
1358// Example: If the shape has sharp corners, make sure to set `$fn/$fs/$fa`.  The corners of this triangle are not round, even though `offset="round"` (the default) because the number of segments is small.
1359//   triangle = [[0,0],[10,0],[5,10]];
1360//   offset_sweep(triangle, height=6, bottom = os_circle(r=-2),steps=4);
1361// Example: Can improve the result by increasing $fn
1362//   $fn=12;
1363//   triangle = [[0,0],[10,0],[5,10]];
1364//   offset_sweep(triangle, height=6, bottom = os_circle(r=-2),steps=4);
1365// Example: Using $fa and $fs works too; it produces a different looking triangulation of the rounded corner
1366//   $fa=1;$fs=0.3;
1367//   triangle = [[0,0],[10,0],[5,10]];
1368//   offset_sweep(triangle, height=6, bottom = os_circle(r=-2),steps=4);
1369// Example: Here is the star chamfered at the top with a teardrop rounding at the bottom. Check out the rounded corners on the chamfer.  The large $fn value ensures a smooth curve on the concave corners of the chamfer.  It has no effect anywhere else on the model.  Observe how the rounded star points vanish at the bottom in the teardrop: the number of vertices does not remain constant from layer to layer.
1370//    star = star(5, r=22, ir=13);
1371//    rounded_star = round_corners(star, cut=flatten(repeat([.5,0],5)), $fn=24);
1372//    offset_sweep(rounded_star, height=20, bottom=os_teardrop(r=4), top=os_chamfer(width=4),$fn=64);
1373// Example: We round a cube using the continous curvature rounding profile.  But note that the corners are not smooth because the curved square collapses into a square with corners.    When a collapse like this occurs, we cannot turn `check_valid` off.  For a better result use `rounded_prism()` instead.
1374//   square = square(1);
1375//   rsquare = round_corners(square, method="smooth", cut=0.1, k=0.7, $fn=36);
1376//   end_spec = os_smooth(cut=0.1, k=0.7, steps=22);
1377//   offset_sweep(rsquare, height=1, bottom=end_spec, top=end_spec);
1378// Example(3D,Med): A nice rounded box, with a teardrop base and circular rounded interior and top
1379//   box = square([255,50]);
1380//   rbox = round_corners(box, method="smooth", cut=4, $fn=12);
1381//   thickness = 2;
1382//   difference(){
1383//     offset_sweep(rbox, height=50, check_valid=false, steps=22,
1384//                  bottom=os_teardrop(r=2), top=os_circle(r=1));
1385//     up(thickness)
1386//       offset_sweep(offset(rbox, r=-thickness, closed=true,check_valid=false),
1387//                    height=48, steps=22, check_valid=false,
1388//                    bottom=os_circle(r=4), top=os_circle(r=-1,extra=1));
1389//   }
1390// Example: This box is much thicker, and cut in half to show the profiles.  Note also that we can turn `check_valid` off for the outside and for the top inside, but not for the bottom inside.  This example shows use of the direct keyword syntax without the helper functions.
1391//   smallbox = square([75,50]);
1392//   roundbox = round_corners(smallbox, method="smooth", cut=4, $fn=12);
1393//   thickness=4;
1394//   height=50;
1395//   back_half(y=25, s=200)
1396//     difference(){
1397//       offset_sweep(roundbox, height=height, bottom=["for","offset_sweep","r",10,"type","teardrop"],
1398//                                             top=["for","offset_sweep","r",2], steps = 22, check_valid=false);
1399//       up(thickness)
1400//         offset_sweep(offset(roundbox, r=-thickness, closed=true),
1401//                       height=height-thickness, steps=22,
1402//                       bottom=["for","offset_sweep","r",6],
1403//                       top=["for","offset_sweep","type","chamfer","angle",30,
1404//                            "chamfer_height",-3,"extra",1,"check_valid",false]);
1405//   }
1406// Example(3D,Med): A box with multiple sections and rounded dividers
1407//   thickness = 2;
1408//   box = square([255,50]);
1409//   cutpoints = [0, 125, 190, 255];
1410//   rbox = round_corners(box, method="smooth", cut=4, $fn=12);
1411//   back_half(y=25, s=700)
1412//     difference(){
1413//       offset_sweep(rbox, height=50, check_valid=false, steps=22,
1414//                    bottom=os_teardrop(r=2), top=os_circle(r=1));
1415//       up(thickness)
1416//         for(i=[0:2]){
1417//           ofs = i==1 ? 2 : 0;
1418//           hole = round_corners([[cutpoints[i]-ofs,0], [cutpoints[i]-ofs,50],
1419//                                 [cutpoints[i+1]+ofs, 50], [cutpoints[i+1]+ofs,0]],
1420//                                method="smooth", cut=4, $fn=36);
1421//           offset_sweep(offset(hole, r=-thickness, closed=true,check_valid=false),
1422//                         height=48, steps=22, check_valid=false,
1423//                         bottom=os_circle(r=4), top=os_circle(r=-1,extra=1));
1424//         }
1425//     }
1426// Example(3D,Med): Star shaped box
1427//   star = star(5, r=22, ir=13);
1428//   rounded_star = round_corners(star, cut=flatten(repeat([.5,0],5)), $fn=24);
1429//   thickness = 2;
1430//   ht=20;
1431//   difference(){
1432//     offset_sweep(rounded_star, height=ht, bottom=["for","offset_sweep","r",4],
1433//                                           top=["for","offset_sweep","r",1], steps=15);
1434//     up(thickness)
1435//         offset_sweep(offset(rounded_star,r=-thickness,closed=true),
1436//                       height=ht-thickness, check_valid=false,
1437//                       bottom=os_circle(r=7), top=os_circle(r=-1, extra=1),$fn=40);
1438//     }
1439// Example: A profile defined by an arbitrary sequence of points.
1440//   star = star(5, r=22, ir=13);
1441//   rounded_star = round_corners(star, cut=flatten(repeat([.5,0],5)), $fn=24);
1442//   profile = os_profile(points=[[0,0],[.3,.1],[.6,.3],[.9,.9], [1.2, 2.7],[.8,2.7],[.8,3]]);
1443//   offset_sweep(reverse(rounded_star), height=20, top=profile, bottom=profile, $fn=32);
1444// Example: Parabolic rounding
1445//   star = star(5, r=22, ir=13);
1446//   rounded_star = round_corners(star, cut=flatten(repeat([.5,0],5)), $fn=24);
1447//   offset_sweep(rounded_star, height=20, top=os_profile(points=[for(r=[0:.1:2])[sqr(r),r]]),
1448//                                          bottom=os_profile(points=[for(r=[0:.2:5])[-sqrt(r),r]]),$fn=32);
1449// Example: This example uses a sine wave offset profile.  Note that we give no specification for the bottom, so it is straight.
1450//   sq = [[0,0],[20,0],[20,20],[0,20]];
1451//   sinwave = os_profile(points=[for(theta=[0:5:720]) [4*sin(theta), theta/700*15]]);
1452//   offset_sweep(sq, height=20, top=sinwave, $fn=32);
1453// Example: The same as the previous example but `offset="delta"`
1454//   sq = [[0,0],[20,0],[20,20],[0,20]];
1455//   sinwave = os_profile(points=[for(theta=[0:5:720]) [4*sin(theta), theta/700*15]]);
1456//   offset_sweep(sq, height=20, top=sinwave, offset="delta");
1457// Example: a box with a flared top.  A nice roundover on the top requires a profile edge, but we can use "extra" to create a small chamfer.
1458//   rhex = round_corners(hexagon(side=10), method="smooth", joint=2, $fs=0.2);
1459//   back_half()
1460//     difference(){
1461//       offset_sweep(rhex, height=10, bottom=os_teardrop(r=2), top=os_teardrop(r=-4, extra=0.2));
1462//       up(1)
1463//         offset_sweep(offset(rhex,r=-1), height=9.5, bottom=os_circle(r=2), top=os_teardrop(r=-4));
1464//     }
1465// Example: Using os_mask to create ogee profiles:
1466//   ogee = mask2d_ogee([
1467//       "xstep",1,  "ystep",1,  // Starting shoulder.
1468//       "fillet",5, "round",5,  // S-curve.
1469//       "ystep",1,              // Ending shoulder.
1470//   ]);
1471//   star = star(5, r=220, ir=130);
1472//   rounded_star = round_corners(star, cut=flatten(repeat([5,0],5)), $fn=24);
1473//   offset_sweep(rounded_star, height=100, top=os_mask(ogee), bottom=os_mask(ogee,out=true));
1474
1475
1476// This function does the actual work of repeatedly calling offset() and concatenating the resulting face and vertex lists to produce
1477// the inputs for the polyhedron module.
1478function _make_offset_polyhedron(path,offsets, offset_type, flip_faces, quality, check_valid, offsetind=0,
1479                                 vertexcount=0, vertices=[], faces=[] )=
1480        offsetind==len(offsets)? (
1481                let(
1482                        bottom = count(len(path),vertexcount),
1483                        oriented_bottom = !flip_faces? bottom : reverse(bottom)
1484                ) [vertices, concat(faces,[oriented_bottom])]
1485        ) : (
1486                let(
1487                        this_offset = offsetind==0? offsets[0][0] : offsets[offsetind][0] - offsets[offsetind-1][0],
1488                        delta = offset_type=="delta" || offset_type=="chamfer" ? this_offset : undef,
1489                        r = offset_type=="round"? this_offset : undef,
1490                        do_chamfer = offset_type == "chamfer"
1491                )
1492                let(
1493                        vertices_faces = offset(
1494                                path, r=r, delta=delta, chamfer = do_chamfer, closed=true,
1495                                check_valid=check_valid, quality=quality,
1496                                return_faces=true,
1497                                firstface_index=vertexcount,
1498                                flip_faces=flip_faces
1499                        )
1500                )
1501                _make_offset_polyhedron(
1502                        vertices_faces[0], offsets, offset_type,
1503                        flip_faces, quality, check_valid, 
1504                        offsetind+1, vertexcount+len(path),
1505                        vertices=concat(
1506                                vertices,
1507                                path3d(vertices_faces[0],offsets[offsetind][1])
1508                        ),
1509                        faces=concat(faces, vertices_faces[1])
1510                )
1511        );
1512
1513
1514function _struct_valid(spec, func, name) =
1515        spec==[] ? true :
1516        assert(is_list(spec) && len(spec)>=2 && spec[0]=="for",str("Specification for \"", name, "\" is an invalid structure"))
1517        assert(spec[1]==func, str("Specification for \"",name,"\" is for a different function (",func,")"));
1518
1519function offset_sweep(
1520                       path, height, 
1521                       bottom=[], top=[], 
1522                       h, l, length, 
1523                       offset="round", r=0, steps=16,
1524                       quality=1, check_valid=true,
1525                       extra=0,
1526                       cut=undef, chamfer_width=undef, chamfer_height=undef,
1527                       joint=undef, k=0.75, angle=45
1528                      ) =
1529    let(
1530        argspec = [
1531                   ["for",""],
1532                   ["r",r],
1533                   ["extra",extra],
1534                   ["type","circle"],
1535                   ["check_valid",check_valid],
1536                   ["quality",quality],
1537                   ["steps",steps],
1538                   ["offset",offset],
1539                   ["chamfer_width",chamfer_width],
1540                   ["chamfer_height",chamfer_height],
1541                   ["angle",angle],
1542                   ["cut",cut],
1543                   ["joint",joint],
1544                   ["k", k],
1545                   ["points", []],
1546        ],
1547        path = force_path(path)
1548    )
1549    assert(is_path(path,2), "Input path must be a 2D path")
1550    let(
1551        clockwise = is_polygon_clockwise(path),
1552        dummy1 = _struct_valid(top,"offset_sweep","top"),
1553        dummy2 = _struct_valid(bottom,"offset_sweep","bottom"),
1554        top = struct_set(argspec, top, grow=false),
1555        bottom = struct_set(argspec, bottom, grow=false),
1556
1557        //  This code does not work.  It hits the error in _make_offset_polyhedron from offset being wrong
1558        //  before this code executes.  Had to move the test into _make_offset_polyhedron, which is ugly since it's in the loop
1559        offsetsok = in_list(struct_val(top, "offset"),["round","delta","chamfer"])
1560                    && in_list(struct_val(bottom, "offset"),["round","delta","chamfer"])
1561    )
1562    assert(offsetsok,"Offsets must be one of \"round\", \"delta\", or \"chamfer\"")
1563    let(
1564        offsets_bot = _rounding_offsets(bottom, -1),
1565        offsets_top = _rounding_offsets(top, 1),
1566        dummy = (struct_val(top,"offset")=="chamfer" && len(offsets_top)>5)
1567                        || (struct_val(bottom,"offset")=="chamfer" && len(offsets_bot)>5)
1568                ? echo("WARNING: You have selected offset=\"chamfer\", which leads to exponential growth in the vertex count and requested more than 5 layers.  This can be slow or run out of recursion depth.")
1569                : 0,
1570
1571        // "Extra" height enlarges the result beyond the requested height, so subtract it
1572        bottom_height = len(offsets_bot)==0 ? 0 : abs(last(offsets_bot)[1]) - struct_val(bottom,"extra"),
1573        top_height = len(offsets_top)==0 ? 0 : abs(last(offsets_top)[1]) - struct_val(top,"extra"),
1574
1575        height = one_defined([l,h,height,length], "l,h,height,length", dflt=u_add(bottom_height,top_height)),
1576        middle = height-bottom_height-top_height
1577    )
1578    assert(height>0, "Height must be positive") 
1579    assert(middle>=0, str("Specified end treatments (bottom height = ",bottom_height,
1580                          " top_height = ",top_height,") are too large for extrusion height (",height,")"
1581                         )
1582    )
1583    let(
1584        initial_vertices_bot = path3d(path),
1585
1586        vertices_faces_bot = _make_offset_polyhedron(
1587                path, offsets_bot, struct_val(bottom,"offset"), clockwise,
1588                struct_val(bottom,"quality"),
1589                struct_val(bottom,"check_valid"),
1590                vertices=initial_vertices_bot
1591        ),
1592
1593        top_start_ind = len(vertices_faces_bot[0]),
1594        initial_vertices_top = path3d(path, middle),
1595        vertices_faces_top = _make_offset_polyhedron(
1596                path, move(p=offsets_top,[0,middle]),
1597                struct_val(top,"offset"), !clockwise,
1598                struct_val(top,"quality"),
1599                struct_val(top,"check_valid"),
1600                vertexcount=top_start_ind,
1601                vertices=initial_vertices_top
1602        ),
1603        middle_faces = middle==0 ? [] : [
1604                for(i=[0:len(path)-1]) let(
1605                        oneface=[i, (i+1)%len(path), top_start_ind+(i+1)%len(path), top_start_ind+i]
1606                ) !clockwise ? reverse(oneface) : oneface
1607        ]
1608    )
1609    [up(bottom_height, concat(vertices_faces_bot[0],vertices_faces_top[0])),  // Vertices
1610     concat(vertices_faces_bot[1], vertices_faces_top[1], middle_faces)];     // Faces
1611
1612
1613module offset_sweep(path, height, 
1614                    bottom=[], top=[], 
1615                    h, l,
1616                    offset="round", r=0, steps=16,
1617                    quality=1, check_valid=true,
1618                    extra=0,
1619                    cut=undef, chamfer_width=undef, chamfer_height=undef,
1620                    joint=undef, k=0.75, angle=45,
1621                    convexity=10,anchor="origin",cp="centroid",
1622                    spin=0, orient=UP, atype="hull")
1623{
1624    assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"");
1625    vnf = offset_sweep(path=path, height=height, h=h, l=l, top=top, bottom=bottom, offset=offset, r=r, steps=steps,
1626                       quality=quality, check_valid=check_valid, extra=extra, cut=cut, chamfer_width=chamfer_width,
1627                       chamfer_height=chamfer_height, joint=joint, k=k, angle=angle);
1628    vnf_polyhedron(vnf,convexity=convexity,anchor=anchor, spin=spin, orient=orient, atype=atype, cp=cp)
1629        children();
1630}   
1631
1632
1633
1634function os_circle(r,cut,extra,check_valid, quality,steps, offset) =
1635        assert(num_defined([r,cut])==1, "Must define exactly one of `r` and `cut`")
1636        _remove_undefined_vals([
1637                "for", "offset_sweep",
1638                "type", "circle",
1639                "r",r,
1640                "cut",cut,
1641                "extra",extra,
1642                "check_valid",check_valid,
1643                "quality", quality,
1644                "steps", steps,
1645                "offset", offset
1646        ]);
1647
1648function os_teardrop(r,cut,extra,check_valid, quality,steps, offset) =
1649        assert(num_defined([r,cut])==1, "Must define exactly one of `r` and `cut`")
1650        _remove_undefined_vals([
1651                "for", "offset_sweep",
1652                "type", "teardrop",
1653                "r",r,
1654                "cut",cut,
1655                "extra",extra,
1656                "check_valid",check_valid,
1657                "quality", quality,
1658                "steps", steps,
1659                "offset", offset
1660        ]);
1661
1662function os_chamfer(height, width, cut, angle, extra,check_valid, quality,steps, offset) =
1663        let(ok = (is_def(cut) && num_defined([height,width])==0) || num_defined([height,width])>0)
1664        assert(ok, "Must define `cut`, or one or both of `width` and `height`")
1665        _remove_undefined_vals([
1666                "for", "offset_sweep",
1667                "type", "chamfer",
1668                "chamfer_width",width,
1669                "chamfer_height",height,
1670                "cut",cut,
1671                "angle",angle,
1672                "extra",extra,
1673                "check_valid",check_valid,
1674                "quality", quality,
1675                "steps", steps,
1676                "offset", offset
1677        ]);
1678
1679function os_smooth(cut, joint, k, extra,check_valid, quality,steps, offset) =
1680        assert(num_defined([joint,cut])==1, "Must define exactly one of `joint` and `cut`")
1681        _remove_undefined_vals([
1682                "for", "offset_sweep",
1683                "type", "smooth",
1684                "joint",joint,
1685                "k",k,
1686                "cut",cut,
1687                "extra",extra,
1688                "check_valid",check_valid,
1689                "quality", quality,
1690                "steps", steps,
1691                "offset", offset
1692        ]);
1693
1694function os_profile(points, extra,check_valid, quality, offset) =
1695        assert(is_path(points),"Profile point list is not valid")
1696        _remove_undefined_vals([
1697                "for", "offset_sweep",
1698                "type", "profile",
1699                "points", points,
1700                "extra",extra,
1701                "check_valid",check_valid,
1702                "quality", quality,
1703                "offset", offset
1704        ]);
1705
1706
1707function os_mask(mask, out=false, extra,check_valid, quality, offset) =
1708  let(
1709      origin_index = [for(i=idx(mask)) if (mask[i].x<0 && mask[i].y<0) i],
1710      xfactor = out ? -1 : 1
1711  )
1712  assert(len(origin_index)==1,"Cannot find origin in the mask")
1713  let(
1714      points = ([for(pt=list_rotate(mask,origin_index[0])) [xfactor*max(pt.x,0),-max(pt.y,0)]])
1715  )
1716  os_profile(deduplicate(move(-points[1],p=list_tail(points))), extra,check_valid,quality,offset);
1717
1718
1719// Module: convex_offset_extrude()
1720// Usage: Basic usage.  See below for full options
1721//   convex_offset_extrude(height, [bottom], [top], ...) 2D-CHILDREN;
1722// Description:
1723//   Extrudes 2d children with layers formed from the convex hull of the offset of each child according to a sequence of offset values.
1724//   Like `offset_sweep` this module can use built-in offset profiles to provide treatments such as roundovers or chamfers but unlike `offset_sweep()` it
1725//   operates on 2d children rather than a point list.  Each offset is computed using
1726//   the native `offset()` module from the input geometry.
1727//   If your shape has corners that you want rounded by offset be sure to set `$fn` or `$fs` appropriately.
1728//   If your geometry has internal holes or is too small for the specified offset then you may get
1729//   unexpected results.
1730//   .
1731//   The build-in profiles are: circular rounding, teardrop rounding, continuous curvature rounding, and chamfer.
1732//   Also note that when a rounding radius is negative the rounding will flare outwards.  The easiest way to specify
1733//   the profile is by using the profile helper functions.  These functions take profile parameters, as well as some
1734//   general settings and translate them into a profile specification, with error checking on your input.  The description below
1735//   describes the helper functions and the parameters specific to each function.  Below that is a description of the generic
1736//   settings that you can optionally use with all of the helper functions.
1737//   For more details on the "cut" and "joint" rounding parameters, and
1738//   on continuous curvature rounding, see [Types of Roundover](rounding.scad#subsection-types-of-roundover). 
1739//   .
1740//   The final shape is created by combining convex hulls of small extrusions.  The thickness of these small extrusions may result
1741//   your model being slightly too long (if the curvature at the end is flaring outward), so if the exact length is very important
1742//   you may need to intersect with a bounding cube.  (Note that extra length can also be intentionally added with the `extra` argument.)
1743//   .
1744//   - profile: os_profile(points)
1745//     Define the offset profile with a list of points.  The first point must be [0,0] and the roundover should rise in the positive y direction, with positive x values for inward motion (standard roundover) and negative x values for flaring outward.  If the y value ever decreases then you might create a self-intersecting polyhedron, which is invalid.  Such invalid polyhedra will create cryptic assertion errors when you render your model and it is your responsibility to avoid creating them.  Note that the starting point of the profile is the center of the extrusion.  If you use a profile as the top it will rise upwards.  If you use it as the bottom it will be inverted, and will go downward.
1746//   - circle: os_circle(r|cut).  Define circular rounding either by specifying the radius or cut distance.
1747//   - smooth: os_smooth(cut|joint, [k]).  Define continuous curvature rounding, with `cut` and `joint` as for round_corners.  The k parameter controls how fast the curvature changes and should be between 0 and 1.
1748//   - teardrop: os_teardrop(r|cut).  Rounding using a 1/8 circle that then changes to a 45 degree chamfer.  The chamfer is at the end, and enables the object to be 3d printed without support.  The radius gives the radius of the circular part.
1749//   - chamfer: os_chamfer([height], [width], [cut], [angle]).  Chamfer the edge at desired angle or with desired height and width.  You can specify height and width together and the angle will be ignored, or specify just one of height and width and the angle is used to determine the shape.  Alternatively, specify "cut" along with angle to specify the cut back distance of the chamfer.
1750//   .
1751//   The general settings that you can use with all of the helper functions are mostly used to control how offset_sweep() calls the offset() function.
1752//   - extra: Add an extra vertical step of the specified height, to be used for intersections or differences.  This extra step will extend the resulting object beyond the height you specify.  Default: 0
1753//   - steps: Number of vertical steps to use for the profile.  (Not used by os_profile).  Default: 16
1754//   - offset: Select "round" (r=), "delta" (delta=), or "chamfer" offset types for offset.  Default: "round"
1755//   .
1756//   Many of the arguments are described as setting "default" values because they establish settings which may be overridden by
1757//   the top and bottom profile specifications.
1758//   .
1759//   You will generally want to use the above helper functions to generate the profiles.
1760//   The profile specification is a list of pairs of keywords and values, e.g. ["r",12, type, "circle"]. The keywords are
1761//   - "type" - type of rounding to apply, one of "circle", "teardrop", "chamfer", "smooth", or "profile" (Default: "circle")
1762//   - "r" - the radius of the roundover, which may be zero for no roundover, or negative to round or flare outward.  Default: 0
1763//   - "cut" - the cut distance for the roundover or chamfer, which may be negative for flares
1764//   - "chamfer_width" - the width of a chamfer
1765//   - "chamfer_height" - the height of a chamfer
1766//   - "angle" - the chamfer angle, measured from the vertical (so zero is vertical, 90 is horizontal).  Default: 45
1767//   - "joint" - the joint distance for a "smooth" roundover
1768//   - "k" - the curvature smoothness parameter for "smooth" roundovers, a value in [0,1].  Default: 0.75
1769//   - "points" - point list for use with the "profile" type
1770//   - "extra" - extra height added for unions/differences.  This makes the shape taller than the requested height.  (Default: 0)
1771//   - "steps" - number of vertical steps to use for the roundover.  Default: 16.
1772//   - "offset" - select "round" (r=) or "delta" (delta=) offset type for offset.  Default: "round"
1773//   .
1774//   Note that unlike `offset_sweep`, because the offset operation is always performed from the base shape, using chamfered offsets does not increase the
1775//   number of vertices or lead to any special complications.
1776//
1777// Arguments:
1778//   height / length / l / h = total height (including rounded portions, but not extra sections) of the output.  Default: combined height of top and bottom end treatments.
1779//   bottom = rounding spec for the bottom end
1780//   top = rounding spec for the top end.
1781//   ---
1782//   offset = default offset, `"round"`, `"delta"`, or `"chamfer"`.  Default: `"round"`
1783//   steps = default step count.  Default: 16
1784//   extra = default extra height.  Default: 0
1785//   cut = default cut value.
1786//   chamfer_width = default width value for chamfers.
1787//   chamfer_height = default height value for chamfers.
1788//   angle = default angle for chamfers.  Default: 45
1789//   joint = default joint value for smooth roundover.
1790//   k = default curvature parameter value for "smooth" roundover
1791//   convexity = convexity setting for use with polyhedron.  Default: 10
1792// Example: Chamfered elliptical prism.  If you stretch a chamfered cylinder the chamfer will be uneven.
1793//   convex_offset_extrude(bottom = os_chamfer(height=-2),
1794//                         top=os_chamfer(height=1), height=7)
1795//     xscale(4)circle(r=6,$fn=64);
1796// Example: Elliptical prism with circular roundovers.
1797//   convex_offset_extrude(bottom=os_circle(r=-2),
1798//                         top=os_circle(r=1), height=7,steps=10)
1799//     xscale(4)circle(r=6,$fn=64);
1800// Example: If you give a non-convex input you get a convex hull output
1801//   right(50) linear_extrude(height=7) star(5,r=22,ir=13);
1802//   convex_offset_extrude(bottom = os_chamfer(height=-2),
1803//                         top=os_chamfer(height=1), height=7, $fn=32)
1804//     star(5,r=22,ir=13);
1805function convex_offset_extrude(
1806        height, 
1807        bottom=[], top=[], 
1808        h, l, length,
1809        offset="round", r=0, steps=16,
1810        extra=0,
1811        cut=undef, chamfer_width=undef, chamfer_height=undef,
1812        joint=undef, k=0.75, angle=45,
1813        convexity=10, thickness = 1/1024
1814) = no_function("convex_offset_extrude");
1815module convex_offset_extrude(
1816        height,
1817        bottom=[],
1818        top=[], 
1819        h, l, length,
1820        offset="round", r=0, steps=16,
1821        extra=0,
1822        cut=undef, chamfer_width=undef, chamfer_height=undef,
1823        joint=undef, k=0.75, angle=45,
1824        convexity=10, thickness = 1/1024
1825) {
1826        req_children($children);  
1827        argspec = [
1828                ["for", ""],
1829                ["r",r],
1830                ["extra",extra],
1831                ["type","circle"],
1832                ["steps",steps],
1833                ["offset",offset],
1834                ["chamfer_width",chamfer_width],
1835                ["chamfer_height",chamfer_height],
1836                ["angle",angle],
1837                ["cut",cut],
1838                ["joint",joint],
1839                ["k", k],
1840                ["points", []],
1841        ];
1842        top = struct_set(argspec, top, grow=false);
1843        bottom = struct_set(argspec, bottom, grow=false);
1844
1845        offsets_bot = _rounding_offsets(bottom, -1);
1846        offsets_top = _rounding_offsets(top, 1);
1847
1848        // "Extra" height enlarges the result beyond the requested height, so subtract it
1849        bottom_height = len(offsets_bot)==0 ? 0 : abs(last(offsets_bot)[1]) - struct_val(bottom,"extra");
1850        top_height = len(offsets_top)==0 ? 0 : abs(last(offsets_top)[1]) - struct_val(top,"extra");
1851
1852        height = one_defined([l,h,height,length], "l,h,height,length", dflt=u_add(bottom_height,top_height));
1853        middle = height-bottom_height-top_height;
1854        check =
1855          assert(height>=0, "Height must be nonnegative")
1856          assert(middle>=0, str(
1857                                "Specified end treatments (bottom height = ",bottom_height,
1858                                " top_height = ",top_height,") are too large for extrusion height (",height,")"
1859                            )
1860          );
1861        // The entry r[i] is [radius,z] for a given layer
1862        r = move([0,bottom_height],p=concat(
1863                          reverse(offsets_bot), [[0,0], [0,middle]], move([0,middle], p=offsets_top)));
1864        delta = [for(val=deltas(column(r,0))) sign(val)];
1865        below=[-thickness,0];
1866        above=[0,thickness];
1867           // layers is a list of pairs of the relative positions for each layer, e.g. [0,thickness]
1868           // puts the layer above the polygon, and [-thickness,0] puts it below.
1869        layers = [for (i=[0:len(r)-1])
1870          i==0 ? (delta[0]<0 ? below : above) :
1871          i==len(r)-1 ? (delta[len(delta)-1] < 0 ? below : above) :
1872          delta[i]==0 ? above :
1873          delta[i+1]==0 ? below :
1874          delta[i]==delta[i-1] ? [-thickness/2, thickness/2] :
1875          delta[i] == 1 ? above :
1876          /* delta[i] == -1 ? */ below];
1877        dochamfer = offset=="chamfer";
1878        for(i=[0:len(r)-2])
1879          for(j=[0:$children-1])
1880           hull(){
1881             up(r[i][1]+layers[i][0])
1882               linear_extrude(convexity=convexity,height=layers[i][1]-layers[i][0])
1883                 if (offset=="round")
1884                   offset(r=r[i][0])
1885                     children(j);
1886                 else
1887                   offset(delta=r[i][0],chamfer = dochamfer)
1888                     children(j);
1889             up(r[i+1][1]+layers[i+1][0])
1890               linear_extrude(convexity=convexity,height=layers[i+1][1]-layers[i+1][0])
1891                 if (offset=="round")
1892                   offset(r=r[i+1][0])
1893                     children(j);
1894                 else
1895                   offset(delta=r[i+1][0],chamfer=dochamfer)
1896                     children(j);
1897           }
1898}
1899
1900
1901
1902function _remove_undefined_vals(list) =
1903        let(ind=search([undef],list,0)[0])
1904        list_remove(list, concat(ind, add_scalar(ind,-1)));
1905
1906
1907
1908function _rp_compute_patches(top, bot, rtop, rsides, ktop, ksides, concave) =
1909   let(
1910     N = len(top),
1911     plane = plane3pt_indexed(top,0,1,2),
1912     rtop_in = is_list(rtop) ? rtop[0] : rtop,
1913     rtop_down = is_list(rtop) ? rtop[1] : abs(rtop)
1914   )
1915  [for(i=[0:N-1])
1916           let(
1917               rside_prev = is_list(rsides[i])? rsides[i][0] : rsides[i],
1918               rside_next = is_list(rsides[i])? rsides[i][1] : rsides[i],
1919               concave_sign = (concave[i] ? -1 : 1) * (rtop_in>=0 ? 1 : -1),  // Negative if normals need to go "out"
1920               prev = select(top,i-1) - top[i],
1921               next = select(top,i+1) - top[i],
1922               prev_offset = top[i] + rside_prev * unit(prev) / sin(vector_angle(prev,bot[i]-top[i])),
1923               next_offset = top[i] + rside_next * unit(next) / sin(vector_angle(next,bot[i]-top[i])),
1924               down = rtop_down * unit(bot[i]-top[i]) / sin(abs(plane_line_angle(plane, [bot[i],top[i]]))),
1925               row2 = [prev_offset,     top[i],     next_offset     ],
1926               row4 = [prev_offset+down,top[i]+down,next_offset+down],
1927               in_prev = concave_sign * unit(next-(next*prev)*prev/(prev*prev)),
1928               in_next = concave_sign * unit(prev-(prev*next)*next/(next*next)),
1929               far_corner = top[i]+ concave_sign*unit(unit(prev)+unit(next))* abs(rtop_in) / sin(vector_angle(prev,next)/2),
1930               row0 =
1931                 concave_sign<0 ?
1932                    [prev_offset+abs(rtop_in)*in_prev, far_corner, next_offset+abs(rtop_in)*in_next]
1933                 :
1934                    let(
1935                       prev_corner = prev_offset + abs(rtop_in)*in_prev,
1936                       next_corner = next_offset + abs(rtop_in)*in_next,
1937                       prev_degenerate = is_undef(line_intersection(path2d([far_corner, far_corner+prev]),
1938                                                                   path2d([prev_offset, prev_offset+in_prev]),RAY,RAY)),
1939                       next_degenerate = is_undef(line_intersection(path2d([far_corner, far_corner+next]),
1940                                                                   path2d([next_offset, next_offset+in_next]),RAY,RAY))
1941                    )
1942                    [ prev_degenerate ? far_corner : prev_corner,
1943                      far_corner,
1944                      next_degenerate ? far_corner : next_corner]
1945            ) _smooth_bez_fill(
1946                      [for(row=[row0, row2, row4]) _smooth_bez_fill(row,ksides[i])],
1947                      ktop)];
1948
1949
1950// Function&Module: rounded_prism()
1951// Usage: as a module
1952//   rounded_prism(bottom, [top], [height=|h=|length=|l=], [joint_top=], [joint_bot=], [joint_sides=], [k=], [k_top=], [k_bot=], [k_sides=], [splinesteps=], [debug=], [convexity=],...) [ATTACHMENTS];
1953// Usage: as a function
1954//   vnf = rounded_prism(bottom, [top], [height=|h=|length=|l=], [joint_top=], [joint_bot=], [joint_sides=], [k=], [k_top=], [k_bot=], [k_sides=], [splinesteps=], [debug=]);
1955// Description:
1956//   Construct a generalized prism with continuous curvature rounding.  You supply the polygons for the top and bottom of the prism.  The only
1957//   limitation is that joining the edges must produce a valid polyhedron with coplanar side faces.  You specify the rounding by giving
1958//   the joint distance away from the corner for the rounding curve.  The k parameter ranges from 0 to 1 with a default of 0.5.  Larger
1959//   values give a more abrupt transition and smaller ones a more gradual transition.  If you set the value much higher
1960//   than 0.8 the curvature changes abruptly enough that though it is theoretically continuous, it may
1961//   not be continuous in practice.  A value of 0.92 is a good approximation to a circle.  If you set it very small then the transition
1962//   is so gradual that the roundover may be very small.  If you want a very smooth roundover, set the joint parameter as large as possible and
1963//   then adjust the k value down as low as gives a sufficiently large roundover.  See
1964//   [Types of Roundover](rounding.scad#subsection-types-of-roundover) for more information on continuous curvature rounding.  
1965//   .
1966//   You can specify the bottom and top polygons by giving two compatible 3d paths.  You can also give 2d paths and a height/length and the
1967//   two shapes will be offset in the z direction from each other.  The final option is to specify just the bottom along with a height/length;
1968//   in this case the top will be a copy of the bottom, offset in the z direction by the specified height.
1969//   .
1970//   You define rounding for all of the top edges, all of the bottom edges, and independently for each of the connecting side edges.
1971//   You specify rounding the rounding by giving the joint distance for where the curved section should start.  If the joint distance is 1 then
1972//   it means the curved section begins 1 unit away from the edge (in the perpendicular direction).  Typically each joint distance is a scalar
1973//   value and the rounding is symmetric around each edge.  However, you can specify a 2-vector for the joint distance to produce asymmetric
1974//   rounding which is different on the two sides of the edge.  This may be useful when one one edge in your polygon is much larger than another.
1975//   For the top and bottom you can specify negative joint distances.  If you give a scalar negative value then the roundover will flare
1976//   outward.  If you give a vector value then a negative value then if joint_top[0] is negative the shape will flare outward, but if
1977//   joint_top[1] is negative the shape will flare upward.  At least one value must be non-negative.  The same rules apply for joint_bot.
1978//   The joint_sides parameter must be entirely nonnegative.
1979//   .
1980//   If you set `debug` to true the module version will display the polyhedron even when it is invalid and it will show the bezier patches at the corners.
1981//   This can help troubleshoot problems with your parameters.  With the function form setting debug to true causes it to return [patches,vnf] where
1982//   patches is a list of the bezier control points for the corner patches.
1983//   .
1984//   Note that rounded_prism() is not well suited to rounding shapes that have already been rounded, or that have many points.
1985//   It works best when the top and bottom are polygons with well-defined corners.  When the polygons have been rounded already,
1986//   further rounding generates tiny bezier patches patches that can more easily
1987//   interfere, giving rise to an invalid polyhedron.  It's also slow because you get bezier patches for every corner in the model.  
1988//   .
1989// Arguments:
1990//   bottom = 2d or 3d path describing bottom polygon
1991//   top = 2d or 3d path describing top polygon (must be the same dimension as bottom)
1992//   ---
1993//   height/length/h/l = height of the shape when you give 2d bottom
1994//   joint_top = rounding length for top (number or 2-vector).  Default: 0
1995//   joint_bot = rounding length for bottom (number or 2-vector).  Default: 0
1996//   joint_sides = rounding length for side edges, a number/2-vector or list of them.  Default: 0
1997//   k = continuous curvature rounding parameter for all edges.  Default: 0.5
1998//   k_top = continuous curvature rounding parameter for top
1999//   k_bot = continuous curvature rounding parameter for bottom
2000//   k_bot = continuous curvature rounding parameter for bottom
2001//   splinesteps = number of segments to use for curved patches.  Default: 16
2002//   debug = turn on debug mode which displays illegal polyhedra and shows the bezier corner patches for troubleshooting purposes.  Default: False
2003//   convexity = convexity parameter for polyhedron(), only for module version.  Default: 10
2004//   anchor = Translate so anchor point is at the origin.  (module only) Default: "origin"
2005//   spin = Rotate this many degrees around Z axis after anchor.  (module only) Default: 0
2006//   orient = Vector to rotate top towards after spin  (module only)
2007//   atype = Select "hull" or "intersect" anchor types.  (module only) Default: "hull"
2008//   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.  (module only) Default: "centroid"
2009// Example: Uniformly rounded pentagonal prism
2010//   rounded_prism(pentagon(3), height=3,
2011//                 joint_top=0.5, joint_bot=0.5, joint_sides=0.5);
2012// Example: Maximum possible rounding.
2013//   rounded_prism(pentagon(3), height=3,
2014//                 joint_top=1.5, joint_bot=1.5, joint_sides=1.5);
2015// Example: Decreasing k from the default of 0.5 to 0.3 gives a smoother round over which takes up more space, so it appears less rounded.
2016//   rounded_prism(pentagon(3), height=3, joint_top=1.5, joint_bot=1.5,
2017//                 joint_sides=1.5, k=0.3, splinesteps=32);
2018// Example: Increasing k from the default of 0.5 to 0.92 approximates a circular roundover, which does not have continuous curvature.  Notice the visible "edges" at the boundary of the corner and edge patches.  
2019//   rounded_prism(pentagon(3), height=3, joint_top=0.5,
2020//                 joint_bot=0.5, joint_sides=0.5, k=0.92);
2021// Example: rounding just one edge
2022//   rounded_prism(pentagon(side=3), height=3, joint_top=0.5, joint_bot=0.5,
2023//                 joint_sides=[0,0,0.5,0,0], splinesteps=16);
2024// Example: rounding all the edges differently
2025//   rounded_prism(pentagon(side=3), height=3, joint_top=0.25, joint_bot=0.5,
2026//                 joint_sides=[1.7,.5,.7,1.2,.4], splinesteps=32);
2027// Example: different k values for top, bottom and sides
2028//   rounded_prism(pentagon(side=3.0), height=3.0, joint_top=1.4, joint_bot=1.4,
2029//                 joint_sides=0.7, k_top=0.7, k_bot=0.3, k_sides=0.5, splinesteps=48);
2030// Example: flared bottom
2031//   rounded_prism(pentagon(3), height=3, joint_top=1.0,
2032//                 joint_bot=-0.5, joint_sides=0.5);
2033// Example: truncated pyramid
2034//   rounded_prism(pentagon(3), apply(scale(.7),pentagon(3)),
2035//                 height=3, joint_top=0.5, joint_bot=0.5, joint_sides=0.5);
2036// Example: top translated
2037//   rounded_prism(pentagon(3), apply(right(2),pentagon(3)),
2038//                 height=3, joint_top=0.5, joint_bot=0.5, joint_sides=0.5);
2039// Example(NORENDER): top rotated: fails due to non-coplanar side faces
2040//   rounded_prism(pentagon(3), apply(rot(45),pentagon(3)), height=3,
2041//                 joint_top=0.5, joint_bot=0.5, joint_sides=0.5);
2042// Example: skew top
2043//   rounded_prism(path3d(pentagon(3)), apply(affine3d_skew_yz(0,-20),path3d(pentagon(3),3)),
2044//                 joint_top=0.5, joint_bot=0.5, joint_sides=0.5);
2045// Example: this rotation gives coplanar sides
2046//   rounded_prism(path3d(square(4)), apply(yrot(-100)*right(2),path3d(square(4),3)),
2047//                 joint_top=0.5, joint_bot=0.5, joint_sides=0.5);
2048// Example: a shape with concave corners
2049//   M = path3d(turtle(["left", 180, "length",3,"move", "left", "move", 3, "right",
2050//                      "move", "right", "move", 4, "right", "move", 3, "right", "move", 2]));
2051//   rounded_prism(M, apply(up(3),M), joint_top=0.75, joint_bot=0.2,
2052//                 joint_sides=[.2,2.5,2,0.5,1.5,.5,2.5], splinesteps=32);
2053// Example: using debug mode to see the corner patch sizes, which may help figure out problems with interfering corners or invalid polyhedra.  The corner patches must not intersect each other.
2054//   M = path3d(turtle(["left", 180, "length",3,"move", "left", "move", 3, "right",
2055//                      "move", "right", "move", 4, "right", "move", 3, "right", "move", 2]));
2056//   rounded_prism(M, apply(up(3),M), joint_top=0.75, joint_bot=0.2,
2057//                 joint_sides=[.2,2.5,2,0.5,1.5,.5,2.5], splinesteps=16,debug=true);
2058// Example: applying transformation to the previous example
2059//   M = path3d(turtle(["left", 180, "length",3,"move", "left", "move", 3, "right",
2060//                      "move", "right", "move", 4, "right", "move", 3, "right", "move", 2]));
2061//   rounded_prism(M, apply(right(1)*scale(.75)*up(3),M), joint_top=0.5, joint_bot=0.2,
2062//                 joint_sides=[.2,1,1,0.5,1.5,.5,2], splinesteps=32);
2063// Example: this example shows most of the different types of patches that rounded_prism creates.  Note that some of the patches are close to interfering with each other across the top of the polyhedron, which would create an invalid result.
2064//   N = apply(rot(180)*yscale(.8),turtle(["length",3,"left", "move", 2, "right", 135,"move", sqrt(2), 
2065//                                         "left", "move", sqrt(2), "right", 135, "move", 2]));
2066//   rounded_prism(N, height=3, joint_bot=0.5, joint_top=1.25, joint_sides=[[1,1.75],0,.5,.5,2], debug=true);
2067// Example: This object has different scales on its different axies.  Here is the largest symmetric rounding that fits.  Note that the rounding is slightly smaller than the object dimensions because of roundoff error.
2068//   rounded_prism(square([100.1,30.1]), height=8.1, joint_top=4, joint_bot=4,
2069//                 joint_sides=15, k_sides=0.3, splinesteps=32);
2070// Example: Using asymetric rounding enables a much more rounded form:
2071//   rounded_prism(square([100.1,30.1]), height=8.1, joint_top=[15,4], joint_bot=[15,4],
2072//                 joint_sides=[[15,50],[50,15],[15,50],[50,15]], k_sides=0.3, splinesteps=32);
2073// Example: Flaring the top upward instead of outward.  The bottom has an asymmetric rounding with a small flare but a large rounding up the side.
2074//   rounded_prism(pentagon(3), height=3, joint_top=[1,-1],
2075//                 joint_bot=[-0.5,2], joint_sides=0.5);
2076// Example: Sideways polygons:
2077//   rounded_prism(apply(yrot(95),path3d(hexagon(3))), apply(yrot(95), path3d(hexagon(3),3)),
2078//                 joint_top=2, joint_bot=1, joint_sides=1);
2079// Example: Chamfer a polyhedron by setting splinesteps to 1
2080//   N = apply(rot(180)*yscale(.8),turtle(["length",3,"left", "move", 2, "right", 135,"move", sqrt(2), 
2081//                                         "left", "move", sqrt(2), "right", 135, "move", 2]));
2082//   rounded_prism(N, height=3, joint_bot=-0.3, joint_top=.4, joint_sides=[.75,0,.2,.2,.7], splinesteps=1);
2083
2084
2085module rounded_prism(bottom, top, joint_bot=0, joint_top=0, joint_sides=0, k_bot, k_top, k_sides,
2086                     k=0.5, splinesteps=16, h, length, l, height, convexity=10, debug=false,
2087                     anchor="origin",cp="centroid",spin=0, orient=UP, atype="hull")
2088{
2089  assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"");
2090  result = rounded_prism(bottom=bottom, top=top, joint_bot=joint_bot, joint_top=joint_top, joint_sides=joint_sides,
2091                         k_bot=k_bot, k_top=k_top, k_sides=k_sides, k=k, splinesteps=splinesteps, h=h, length=length, height=height, l=l,debug=debug);
2092  vnf = debug ? result[1] : result;
2093  attachable(anchor=anchor, spin=spin, orient=orient, vnf=vnf, extent=atype=="hull", cp=cp)
2094  {
2095    if (debug){
2096        vnf_polyhedron(vnf, convexity=convexity);
2097        debug_bezier_patches(result[0], showcps=true, splinesteps=splinesteps, $fn=16, showdots=false, showpatch=false);
2098    }
2099    else vnf_polyhedron(vnf,convexity=convexity);
2100    children();
2101  }
2102}
2103
2104
2105function rounded_prism(bottom, top, joint_bot=0, joint_top=0, joint_sides=0, k_bot, k_top, k_sides, k=0.5, splinesteps=16,
2106                       h, length, l, height, debug=false) =
2107   let(
2108       bottom = force_path(bottom,"bottom"),
2109       top = force_path(top,"top")
2110   )
2111   assert(is_path(bottom,[2,3]) && len(bottom)>=3, "bottom must be a 2D or 3D path")
2112   assert(is_num(k) && k>=0 && k<=1, "Curvature parameter k must be in interval [0,1]")
2113   let(
2114     N=len(bottom),
2115     k_top = default(k_top, k),
2116     k_bot = default(k_bot, k),
2117     k_sides = default(k_sides, k),
2118     height = one_defined([h,l,height,length],"height,length,l,h", dflt=undef),
2119     shapedimok = (len(bottom[0])==3 && is_path(top,3)) || (len(bottom[0])==2 && (is_undef(top) || is_path(top,2)))
2120   )
2121   assert(is_num(k_top) && k_top>=0 && k_top<=1, "Curvature parameter k_top must be in interval [0,1]")
2122   assert(is_num(k_bot) && k_bot>=0 && k_bot<=1, "Curvature parameter k_bot must be in interval [0,1]")
2123   assert(shapedimok,"bottom and top must be 2d or 3d paths with the same dimension")
2124   assert(len(bottom[0])==3 || is_num(height),"Must give height/length with 2d polygon input")
2125   let(
2126     // Determine which points are concave by making bottom 2d if necessary
2127     bot_proj = len(bottom[0])==2 ? bottom :  project_plane(select(bottom,0,2),bottom),
2128     bottom_sign = is_polygon_clockwise(bot_proj) ? 1 : -1,
2129     concave = [for(i=[0:N-1]) bottom_sign*sign(_point_left_of_line2d(select(bot_proj,i+1), select(bot_proj, i-1,i)))>0],
2130     top = is_undef(top) ? path3d(bottom,height/2) :
2131           len(top[0])==2 ? path3d(top,height/2) :
2132           top,
2133     bottom = len(bottom[0])==2 ? path3d(bottom,-height/2) : bottom,
2134     jssingleok = (is_num(joint_sides) && joint_sides >= 0) || (is_vector(joint_sides,2) && joint_sides[0]>=0 && joint_sides[1]>=0),
2135     jsvecok = is_list(joint_sides) && len(joint_sides)==N && []==[for(entry=joint_sides) if (!(is_num(entry) || is_vector(entry,2))) entry]
2136   )
2137   assert(is_num(joint_top) || is_vector(joint_top,2))
2138   assert(is_num(joint_bot) || is_vector(joint_bot,2))
2139   assert(is_num(joint_top) || (joint_top[0]>=0 ||joint_top[1]>=0), "Both entries in joint_top cannot be negative")
2140   assert(is_num(joint_bot) || (joint_bot[0]>=0 ||joint_bot[1]>=0), "Both entries in joint_bot cannot be negative")
2141   assert(jsvecok || jssingleok,
2142          str("Argument joint_sides is invalid.  All entries must be nonnegative, and it must be a number, 2-vector, or a length ",N," list those."))
2143   assert(is_num(k_sides) || is_vector(k_sides,N), str("Curvature parameter k_sides must be a number or length ",N," vector"))
2144   assert(is_coplanar(bottom))
2145   assert(is_coplanar(top))
2146   assert(!is_num(k_sides) || (k_sides>=0 && k_sides<=1), "Curvature parameter k_sides must be in interval [0,1]")
2147   let(
2148     non_coplanar=[for(i=[0:N-1]) if (!is_coplanar(concat(select(top,i,i+1), select(bottom,i,i+1)))) [i,(i+1)%N]],
2149     k_sides_vec = is_num(k_sides) ? repeat(k_sides, N) : k_sides,
2150     kbad = [for(i=[0:N-1]) if (k_sides_vec[i]<0 || k_sides_vec[i]>1) i],
2151     joint_sides_vec = jssingleok ? repeat(joint_sides,N) : joint_sides,
2152     top_collinear = [for(i=[0:N-1]) if (is_collinear(select(top,i-1,i+1))) i],
2153     bot_collinear = [for(i=[0:N-1]) if (is_collinear(select(bottom,i-1,i+1))) i]
2154   )
2155   assert(non_coplanar==[], str("Side faces are non-coplanar at edges: ",non_coplanar))
2156   assert(top_collinear==[], str("Top has collinear or duplicated points at indices: ",top_collinear))
2157   assert(bot_collinear==[], str("Bottom has collinear or duplicated points at indices: ",bot_collinear))
2158   assert(kbad==[], str("k_sides parameter outside interval [0,1] at indices: ",kbad))
2159   let(
2160     top_patch = _rp_compute_patches(top, bottom, joint_top, joint_sides_vec, k_top, k_sides_vec, concave),
2161     bot_patch = _rp_compute_patches(bottom, top, joint_bot, joint_sides_vec, k_bot, k_sides_vec, concave),
2162
2163     vertbad = [for(i=[0:N-1])
2164                   if (norm(top[i]-top_patch[i][4][2]) + norm(bottom[i]-bot_patch[i][4][2]) > norm(bottom[i]-top[i])) i],
2165     topbad = [for(i=[0:N-1])
2166                   if (norm(top_patch[i][2][4]-top_patch[i][2][2]) + norm(select(top_patch,i+1)[2][0]-select(top_patch,i+1)[2][2])
2167                  > norm(top_patch[i][2][2] - select(top_patch,i+1)[2][2]))   [i,(i+1)%N]],
2168     botbad = [for(i=[0:N-1])
2169                   if (norm(bot_patch[i][2][4]-bot_patch[i][2][2]) + norm(select(bot_patch,i+1)[2][0]-select(bot_patch,i+1)[2][2])
2170                  > norm(bot_patch[i][2][2] - select(bot_patch,i+1)[2][2]))   [i,(i+1)%N]],
2171     topinbad = [for(i=[0:N-1])
2172                   if (norm(top_patch[i][0][2]-top_patch[i][0][4]) + norm(select(top_patch,i+1)[0][0]-select(top_patch,i+1)[0][2])
2173                          > norm(top_patch[i][0][2]-select(top_patch,i+1)[0][2])) [i,(i+1)%N]],
2174     botinbad = [for(i=[0:N-1])
2175                   if (norm(bot_patch[i][0][2]-bot_patch[i][0][4]) + norm(select(bot_patch,i+1)[0][0]-select(bot_patch,i+1)[0][2])
2176                          > norm(bot_patch[i][0][2]-select(bot_patch,i+1)[0][2])) [i,(i+1)%N]]
2177   )
2178   assert(debug || vertbad==[], str("Top and bottom joint lengths are too large; they interfere with each other at vertices: ",vertbad))
2179   assert(debug || topbad==[], str("Joint lengths too large at top edges: ",topbad))
2180   assert(debug || botbad==[], str("Joint lengths too large at bottom edges: ",botbad))
2181   assert(debug || topinbad==[], str("Joint length too large on the top face at edges: ", topinbad))
2182   assert(debug || botinbad==[], str("Joint length too large on the bottom face at edges: ", botinbad))
2183   let(
2184     // Entries in the next two lists have the form [edges, vnf] where
2185     // edges is a list [leftedge, rightedge, topedge, botedge]
2186     top_samples = [for(patch=top_patch) bezier_vnf_degenerate_patch(patch,splinesteps,reverse=false,return_edges=true) ],
2187     bot_samples = [for(patch=bot_patch) bezier_vnf_degenerate_patch(patch,splinesteps,reverse=true,return_edges=true) ],
2188     leftidx=0,
2189     rightidx=1,
2190     topidx=2,
2191     botidx=3,
2192     edge_points =
2193       [for(i=[0:N-1])
2194            let(
2195               top_edge  = [ top_samples[i][1][rightidx], select(top_samples, i+1)[1][leftidx]],
2196               bot_edge  = [ select(bot_samples, i+1)[1][leftidx], bot_samples[i][1][rightidx]],
2197               vert_edge = [ bot_samples[i][1][botidx], top_samples[i][1][botidx]]
2198               )
2199               each [top_edge, bot_edge, vert_edge] ],
2200     faces = [
2201              [for(i=[0:N-1]) each top_samples[i][1][topidx]],
2202              [for(i=[N-1:-1:0]) each reverse(bot_samples[i][1][topidx])],
2203              for(i=[0:N-1]) [
2204                                 bot_patch[i][4][4],
2205                                 select(bot_patch,i+1)[4][0],
2206                                 select(top_patch,i+1)[4][0],
2207                                 top_patch[i][4][4]
2208                             ]
2209             ],
2210     top_simple = is_path_simple(project_plane(faces[0],faces[0]),closed=true),
2211     bot_simple = is_path_simple(project_plane(faces[1],faces[1]),closed=true),
2212     // verify vertical edges
2213     verify_vert =
2214       [for(i=[0:N-1],j=[0:4])
2215         let(
2216               vline = concat(select(column(top_patch[i],j),2,4),
2217                              select(column(bot_patch[i],j),2,4))
2218             )
2219         if (!is_collinear(vline)) [i,j]],
2220     //verify horiz edges
2221     verify_horiz=[for(i=[0:N-1], j=[0:4])
2222         let(
2223             hline_top = concat(select(top_patch[i][j],2,4), select(select(top_patch, i+1)[j],0,2)),
2224             hline_bot = concat(select(bot_patch[i][j],2,4), select(select(bot_patch, i+1)[j],0,2))
2225         )
2226         if (!is_collinear(hline_top) || !is_collinear(hline_bot)) [i,j]]
2227    )
2228    assert(debug || top_simple,
2229          "Roundovers interfere with each other on top face: either input is self intersecting or top joint length is too large")
2230    assert(debug || bot_simple,
2231          "Roundovers interfere with each other on bottom face: either input is self intersecting or top joint length is too large")
2232    assert(debug || (verify_vert==[] && verify_horiz==[]), "Curvature continuity failed")
2233    let( 
2234        vnf = vnf_join([ each column(top_samples,0),
2235                          each column(bot_samples,0),
2236                          for(pts=edge_points) vnf_vertex_array(pts),
2237                          debug ? vnf_from_polygons(faces) 
2238                                : vnf_triangulate(vnf_from_polygons(faces))
2239                       ])
2240    )
2241    debug ? [concat(top_patch, bot_patch), vnf] : vnf;
2242
2243
2244
2245// Converts a 2d path to a path on a cylinder at radius r
2246function _cyl_hole(r, path) =
2247    [for(point=path) cylindrical_to_xyz(concat([r],xscale(360/(2*PI*r),p=point)))];
2248
2249// Mask profile of 180 deg of a circle to round an edge
2250function _circle_mask(r) =
2251   let(eps=0.1)
2252
2253   fwd(r+.01,p=
2254   [
2255    [r+eps,0],
2256    each arc(r=r, angle=[0, 180]),
2257    [-r-eps,0],
2258    [-r-eps, r+3*eps],
2259    [r+eps, r+3*eps]
2260   ]);
2261
2262
2263// Module: bent_cutout_mask()
2264// Usage:
2265//   bent_cutout_mask(r|radius, thickness, path);
2266// Description:
2267//   Creates a mask for cutting a round-edged hole out of a vertical cylindrical shell.  The specified radius
2268//   is the center radius of the cylindrical shell.  The path needs to be sampled finely enough
2269//   so that it can follow the curve of the cylinder.  The thickness may need to be slighly oversized to
2270//   handle the faceting of the cylinder.  The path is wrapped around a cylinder, keeping the
2271//   same dimensions that is has on the plane, with y axis mapping to the z axis and the x axis bending
2272//   around the curve of the cylinder.  The angular span of the path on the cylinder must be somewhat
2273//   less than 180 degrees, and the path shouldn't have closely spaced points at concave points of high curvature because
2274//   this will cause self-intersection in the mask polyhedron, resulting in CGAL failures.
2275// Arguments:
2276//   r / radius = center radius of the cylindrical shell to cut a hole in
2277//   thickness = thickness of cylindrical shell (may need to be slighly oversized)
2278//   path = 2d path that defines the hole to cut
2279// Example: The mask as long pointed ends because this was the most efficient way to close off those ends.
2280//   bent_cutout_mask(10, 1, apply(xscale(3),circle(r=3)),$fn=64);
2281// Example: An elliptical hole.  Note the thickness is slightly increased to 1.05 compared to the actual thickness of 1.
2282//   rot(-90) {
2283//     $fn=128;
2284//     difference(){
2285//       cyl(r=10.5, h=10);
2286//       cyl(r=9.5, h=11);
2287//       bent_cutout_mask(10, 1.05, apply(xscale(3),circle(r=3)),
2288//                        $fn=64);
2289//     }
2290//   }
2291// Example: An elliptical hole in a thick cylinder
2292//   rot(-90) {
2293//     $fn=128;
2294//     difference(){
2295//       cyl(r=14.5, h=15);
2296//       cyl(r=9.5, h=16);
2297//       bent_cutout_mask(12, 5.1, apply(xscale(3),circle(r=3)));
2298//     }
2299//   }
2300// Example: Complex shape example
2301//   rot(-90) {
2302//     $fn=128;
2303//     difference(){
2304//       cyl(r=10.5, h=10, $fn=128);
2305//       cyl(r=9.5, h=11, $fn=128);
2306//       bent_cutout_mask(10, 1.05,
2307//         apply(scale(3),
2308//           supershape(step=2,m1=5, n1=0.3,n2=1.7)),$fn=32);
2309//     }
2310//   }
2311// Example: this shape is invalid due to self-intersections at the inner corners
2312//   rot(-90) {
2313//     $fn=128;
2314//     difference(){
2315//       cylinder(r=10.5, h=10,center=true);
2316//       cylinder(r=9.5, h=11,center=true);
2317//       bent_cutout_mask(10, 1.05,
2318//         apply(scale(3),
2319//           supershape(step=2,m1=5, n1=0.1,n2=1.7)),$fn=32);
2320//     }
2321//   }
2322// Example: increasing the step gives a valid shape, but the shape looks terrible with so few points.
2323//   rot(-90) {
2324//     $fn=128;
2325//     difference(){
2326//       cylinder(r=10.5, h=10,center=true);
2327//       cylinder(r=9.5, h=11,center=true);
2328//       bent_cutout_mask(10, 1.05,
2329//         apply(scale(3),
2330//           supershape(step=12,m1=5, n1=0.1,n2=1.7)),$fn=32);
2331//     }
2332//   }
2333// Example: uniform resampling produces a somewhat better result, but room remains for improvement.  The lesson is that concave corners in your cutout cause trouble.  To get a very good result we need to non-uniformly sample the supershape with more points at the star tips and few points at the inner corners.
2334//   rot(-90) {
2335//     $fn=128;
2336//     difference(){
2337//       cylinder(r=10.5, h=10,center=true);
2338//       cylinder(r=9.5, h=11,center=true);
2339//       bent_cutout_mask(10, 1.05,
2340//         apply(scale(3), resample_path(
2341//              supershape(step=1,m1=5, n1=0.10,n2=1.7),
2342//              60,closed=true)),
2343//         $fn=32);
2344//     }
2345//   }
2346// Example: The cutout spans 177 degrees.  If you decrease the tube radius to 2.5 the cutout spans over 180 degrees and the model fails.
2347//   r=2.6;     // Don't make this much smaller or it will fail
2348//   rot(-90) {
2349//     $fn=128;
2350//     difference(){
2351//       tube(or=r, wall=1, h=10, anchor=CENTER);
2352//       bent_cutout_mask(r-0.5, 1.05,
2353//         apply(scale(3),
2354//           supershape(step=1,m1=5, n1=0.15,n2=1.7)),$fn=32);
2355//     }
2356//   }
2357// Example: A square hole is not as simple as it seems.  The model valid, but wrong, because the square didn't have enough samples to follow the curvature of the cylinder.
2358//   r=25;
2359//   rot(-90) {
2360//     $fn=128;
2361//     difference(){
2362//       tube(or=r, wall=2, h=35, anchor=BOT);
2363//       bent_cutout_mask(r-1, 2.1, back(5,p=square([18,18])));
2364//     }
2365//   }
2366// Example: Adding additional points fixed this problem
2367//   r=25;
2368//   rot(-90) {
2369//     $fn=128;
2370//     difference(){
2371//       tube(or=r, wall=2, h=35, anchor=BOT);
2372//       bent_cutout_mask(r-1, 2.1,
2373//         subdivide_path(back(5,p=square([18,18])),64,closed=true));
2374//     }
2375//   }
2376// Example: Rounding just the exterior corners of this star avoids the problems we had above with concave corners of the supershape, as long as we don't oversample the star.
2377//   r=25;
2378//   rot(-90) {
2379//     $fn=128;
2380//     difference(){
2381//       tube(or=r, wall=2, h=35, anchor=BOT);
2382//       bent_cutout_mask(r-1, 2.1,
2383//         apply(back(15),
2384//           subdivide_path(
2385//             round_corners(star(n=7,ir=5,or=10),
2386//                           cut=flatten(repeat([0.5,0],7)),$fn=32),
2387//             14*15,closed=true)));
2388//     }
2389//   }
2390// Example(2D): Cutting a slot in a cylinder is tricky if you want rounded corners at the top.  This slot profile has slightly angled top edges to blend into the top edge of the cylinder.
2391//   function slot(slotwidth, slotheight, slotradius) = let(
2392//       angle = 85,
2393//       slot = round_corners(
2394//           turtle([
2395//               "right",
2396//               "move", slotwidth,
2397//               "left", angle,
2398//               "move", 2*slotwidth,
2399//               "right", angle,
2400//               "move", slotheight,
2401//               "left",
2402//               "move", slotwidth,
2403//               "left",
2404//               "move", slotheight,
2405//               "right", angle,
2406//               "move", 2*slotwidth,
2407//               "left", angle,
2408//               "move", slotwidth
2409//           ]),
2410//           radius = [0,0,each repeat(slotradius,4),0,0], closed=false
2411//       )
2412//   ) apply(left(max(column(slot,0))/2)*fwd(min(column(slot,1))), slot);
2413//   stroke(slot(15,29,7));
2414// Example: A cylindrical container with rounded edges and a rounded finger slot.
2415//   function slot(slotwidth, slotheight, slotradius) = let(
2416//       angle = 85,
2417//       slot = round_corners(
2418//           turtle([
2419//               "right",
2420//               "move", slotwidth,
2421//               "left", angle,
2422//               "move", 2*slotwidth,
2423//               "right", angle,
2424//               "move", slotheight,
2425//               "left",
2426//               "move", slotwidth,
2427//               "left",
2428//               "move", slotheight,
2429//               "right", angle,
2430//               "move", 2*slotwidth,
2431//               "left", angle,
2432//               "move", slotwidth
2433//           ]),
2434//           radius = [0,0,each repeat(slotradius,4),0,0], closed=false
2435//       )
2436//   ) apply(left(max(column(slot,0))/2)*fwd(min(column(slot,1))), slot);
2437//   diam = 80;
2438//   wall = 4;
2439//   height = 40;
2440//   rot(-90) {
2441//       $fn=128;
2442//       difference(){
2443//           cyl(d=diam, rounding=wall/2, h=height, anchor=BOTTOM);
2444//           up(wall)cyl(d=diam-2*wall, rounding1=wall, rounding2=-wall/2, h=height-wall+.01, anchor=BOTTOM);
2445//           bent_cutout_mask(diam/2-wall/2, wall+.1, subdivide_path(apply(back(10),slot(15, 29, 7)),250));
2446//       }
2447//   }
2448function bent_cutout_mask(r, thickness, path, radius, convexity=10) = no_function("bent_cutout_mask");
2449module bent_cutout_mask(r, thickness, path, radius, convexity=10)
2450{
2451  no_children($children);
2452  r = get_radius(r1=r, r2=radius);
2453  dummy1=assert(is_def(r) && r>0,"Radius of the cylinder to bend around must be positive");
2454  path2 = force_path(path);
2455  dummy2=assert(is_path(path2,2),"Input path must be a 2D path")
2456         assert(r-thickness>0, "Thickness too large for radius")
2457         assert(thickness>0, "Thickness must be positive");
2458  fixpath = clockwise_polygon(path2);
2459  curvepoints = arc(d=thickness, angle = [-180,0]);
2460  profiles = [for(pt=curvepoints) _cyl_hole(r+pt.x,apply(xscale((r+pt.x)/r), offset(fixpath,delta=thickness/2+pt.y,check_valid=false,closed=true)))];
2461  pathx = column(fixpath,0);
2462  minangle = (min(pathx)-thickness/2)*360/(2*PI*r);
2463  maxangle = (max(pathx)+thickness/2)*360/(2*PI*r);
2464  mindist = (r+thickness/2)/cos((maxangle-minangle)/2);
2465  dummy3 = assert(maxangle-minangle<180,"Cutout angle span is too large.  Must be smaller than 180.");
2466  zmean = mean(column(fixpath,1));
2467  innerzero = repeat([0,0,zmean], len(fixpath));
2468  outerpt = repeat( [1.5*mindist*cos((maxangle+minangle)/2),1.5*mindist*sin((maxangle+minangle)/2),zmean], len(fixpath));
2469  default_tag("remove")
2470    vnf_polyhedron(vnf_vertex_array([innerzero, each profiles, outerpt],col_wrap=true),convexity=convexity);
2471}
2472
2473
2474
2475/*
2476
2477join_prism To Do List:
2478
2479special handling for planar joins?
2480   offset method
2481   cut, radius?
2482Access to the derivative smoothing parameter?   
2483
2484*/
2485
2486
2487
2488// Function&Module: join_prism()
2489// Usage: The two main forms with most common options
2490//   join_prism(polygon, base, length=|height=|l=|h=, fillet=, [base_T=], [scale=], [prism_end_T=], [short=], ...) [ATTACHMENTS];
2491//   join_prism(polygon, base, aux=, fillet=, [base_T=], [aux_T=], [scale=], [prism_end_T=], [short=], ...) [ATTACHMENTS];
2492// Usage: As function
2493//   vnf = join_prism( ... );
2494// Description:
2495//   This function creates a smooth fillet between one or both ends of an arbitrary prism and various other shapes: a plane, a sphere, a cylinder,
2496//   or another arbitrary prism.  The fillet is a continuous curvature rounding with a specified width/height.  This module is very general
2497//   and hence has a complex interface.  The examples below form a tutorial on how to use `join_prism` that steps
2498//   through the various options and how they affect the results.  Be sure to check the examples for help understanding how the various options work.
2499//   .
2500//   When joining between planes this function produces similar results to {{rounded_prism()}}.  This function works best when the prism
2501//   cross section is a continuous shape with a high sampling rate and without sharp corners.  If you have sharp corners you should consider
2502//   giving them a small rounding first.  When the prism cross section has concavities the fillet size will be limited by the curvature of those concavities.
2503//   In contrast, {{rounded_prism()}} works best on a prism that has fewer points.  A high sampling rate can lead to problems, and rounding
2504//   over sharp corners leads to poor results.  
2505//   .
2506//   You specify the prism by giving its cross section as a 2D path.  The cross section will always be the orthogonal cross
2507//   section of the prism.  Depending on end conditions, the ends may not be perpendicular to the
2508//   axis of the prism, but the cross section you give *is* always perpendicular to that cross section.
2509// Figure(3D,Big,NoScales,VPR=[74.6, 0, 329.7], VPT=[28.5524, 35.3006, 22.522], VPD=325.228): The layout and terminology used by `join_prism`.  The "base object" is centered on the origin.  The "auxiliary object" (if present) is some distance away so there is room for the "joiner prism" to connect the two objects.  The blue line is the axis of the jointer prism.  It will be at the origin of the shape you supply for defining that prism.  The "root" point of the joiner prism is the point where the prism axis intersects the base.  The prism end point is where the prism axis intersects the auxiliary object.  If you don't give an auxiliary object then the prism end point is distance `length` along the axis from the root.  
2510//   aT = right(-10)*back(0)*up(75)*xrot(-35)*zrot(75);
2511//   br = 17;
2512//   ar = 15;
2513//   xcyl(r=br, l=50, circum=true, $fn=64);
2514//   multmatrix(aT)right(15)xcyl(r=ar,circum=true,l=50,$fn=64);
2515//   %join_prism(circle(r=10), base = "cyl", base_r=br, aux="cyl", aux_r=ar, aux_T=aT,fillet=3);
2516//   root = [-2.26667, 0, 17];
2517//   rback = [15,0,25];
2518//   endpt =  [-7.55915, 0, 56.6937];
2519//   endback = [10,0,55];
2520//   stroke([root,endpt],
2521//          width=1,endcap_width=3,endcaps="dot",endcap_color="red",color="blue",$fn=16);
2522//   stroke(move(3*unit(rback-root), [rback,root]), endcap2="arrow2",width=1/2,$fn=16,color="black");
2523//   down(0)right(4)color("black")move(rback)rot($vpr)text("prism root point",size=4);
2524//   stroke(move(3*unit(endback-endpt), [endback,endpt]), endcap2="arrow2", width=1/2, $fn=16, color="black");
2525//   down(2)right(4)color("black")move(endback)rot($vpr)text("prism end point",size=4);
2526//   right(4)move(-20*[1,1])color("black")rot($vpr)text("base",size=8);
2527//   up(83)right(-10)move(-20*[1,1])color("black")rot($vpr)text("aux",size=8);
2528//   aend=[-13,13,30];
2529//   ast=aend+10*[-1,1,0];
2530//   stroke([ast,aend],endcap2="arrow2", width=1/2, color="black");
2531//   left(2)move(ast)rot($vpr)color("black")text("joiner prism",size=5,anchor=RIGHT);
2532// Continues:
2533//   You must include a base ("plane", "sphere", "cylinder", "cyl"), or a polygon describing the cross section of a base prism.  If you specify a
2534//   sphere or cylinder you must give `base_r` or `base_d` to specify the radius or diameter of the base object.  If you choose a cylinder or a polygonal
2535//   prism then the base object appears aligned with the X axis.  In the case of the planar base, the
2536//   joining prism will have one end of its axis at the origin.  As shown above, the point where the joining prism attaches to its base is the "root" of the prism.
2537//   If you use some other base shape, the root will be adjusted so that it is on the boundary of your shape.  This happens by finding the intersection
2538//   of the joiner prisms's axis and using that as the root.  By default the prism axis is parallel to the Z axis.  
2539//   .
2540//   You may give `base_T`, a rotation operator that will be applied to the base.  This is
2541//   useful to tilt a planar or cylindrical base.  The `base_T` operator must be an origin-centered rotation like yrot(25).  
2542//   .
2543//   You may optionally specify an auxiliary shape.  When you do this, the joining prism connects the base to the auxiliary shape,
2544//   which must be one of "none", "plane", "sphere", "cyl", or "cylinder".  You can also set it to a polygon to create an arbitrary
2545//   prism for the auxiliary shape.  As is the case for the base, auxiliary cylinders and prisms appear oriented along the X axis.  
2546//   For a cylinder or sphere you must use `aux_r` or `aux_d` to specify the radius or diameter.
2547//   The auxiliary shape appears centered on the origin and will most likely be invalid as an end location unless you translate it to a position
2548//   away from the base object.  The `aux_T` operator operates on the auxiliary object, and unlike `base_T` can be a rotation that includes translation
2549//   operations (or is a non-centered rotation).
2550//   .
2551//   When you specify an auxiliary object, the joiner prism axis is initially the line connecting the origin (the base center point) to the auxiliary
2552//   object center point.  The joiner prism end point is determined analogously to how the root is determined, by intersecting the joiner
2553//   prism axis with the auxiliary object.  Note that this means that if `aux_T` is a rotation it will change the joiner prism root, because
2554//   the rotated prism axis will intersect the base in a different location.  If you do not give an auxiliary object then you must give
2555//   the length/height parameter to specify the prism length.  This gives the length of the prism measured from the root to the end point.
2556//   Note that the joint with a curved base may significantly extend the length of the joiner prism: it's total length will often be larger than
2557//   the length you request.  
2558//   .
2559//   For the cylinder and spherical objects you may with to joint a prism to the concave surface.  You can do this by setting a negative
2560//   radius for the base or auxiliary object.  When `base_r` is negative, and the joiner prism axis is vertical, the prism root will be **below** the
2561//   XY plane.  In this case it is actually possible to use the same object for base and aux and you can get a joiner prism that crosses a cylindrical
2562//   or spherical hole.
2563//   .
2564//   When placing prisms inside a hole, an ambiguity can arise about how to identify the root and end of the joiner prism.  The prism axis will have
2565//   two intersections with a cylinder and both are potentially valid roots.  When the auxiliary object is entirely inside the hole, or the auxiliary
2566//   object is a sphere or cylinder with negative radius that intersections the base, both prism directions produce a valid
2567//   joiner prism that meets the hole's concave surface, so two valid interpretations exist.  By default, the longer prism will be returned.
2568//   You can select the shorter prism by setting `short=true`.  If you specify `short=true` when the base has a negative radius, but only one valid
2569//   prism exists, you'll get an error, but it won't clearly identify that a bogus `short=true` was the real cause.  
2570//   .
2571//   You can also alter your prism by using the `prism_end_T` operator which applies to the end point of the prism.  It does not effect
2572//   the root  of the prism.  The `prism_end_T` operator is applied in a coordinate system where the root of the
2573//   prism is the origin, so if you set it to a rotation the prism base will stay rooted at the same location and the prism will rotate 
2574//   in the specified fashion.  After `prism_end_T` is applied, the prism axis will probably be different and the resulting new end point will
2575//   probably not be on the auxiliary object, or it will have changed the length of the prism.  Therefore, the end point is recalculated
2576//   to achieve the specified length (if aux is "none") or to contact the auxiliary object, if you have specified one.  This means, for example,
2577//   that setting `prism_end_T` to a scale operation won't change the result because it doesn't alter the prism axis.  
2578//   .
2579//   The size of the fillets is determined by the fillet, `fillet_base`, and `fillet_aux` parameters.  The fillet parameter will control both
2580//   ends of the prism, or you can set the ends independently.  The fillets must be nonnegative except when the prism joints a plane.
2581//   In this case a negative fillet gives a roundover.  In the case of no auxiliary object you can use `round_end` to round over the planar
2582//   far end of the joiner prism.  By default, the fillet is constructed using a method that produces a fillet with a uniform height along
2583//   the joiner prism.  This can be limiting when connectijng to objects with high curvature, so you can turn it off using the `uniform` option.
2584//   See the figures below for an explanation of the uniform and non-uniform filleting methods.  
2585//   .
2586//   The overlap is a potentially tricky parameter.  It specifies how much extra material to
2587//   create underneath the filleted prism so it overlaps the object that it joins to, ensuring valid unions.
2588//   For joins to convex objects you can choose a small value, but when joining to a concave object the overlap may need to be
2589//   very large to ensure that the base of the joiner prism is well-behaved.  In such cases you may need to use an intersection
2590//   remove excess base.
2591// Figure(2D,Med,NoAxes): Uniform fillet method.  This image shows how the fillet we construct a uniform fillet.  The pictures shows the cross section that is perpendicular to the prism.  The blue curve represents the base object surface.  The vertical line is the side of the prism.  To construct a fillet we travel along the surface of the base, following the curve, until we have moved the fillet length, `a`.  This defines the point `u`.  We then construct a tangent line to the base and find its intersection, `v`, with the prism.  Note that if the base is steeply curved, this tangent may fail to intersect, and the algorithm will fail with an error because `v` does not exist.  Finally we locate `w` to be distance `a` above the point where the prism intersects the base object.  The fillet is defined by the `[u,v,w]` triple and is shown in red.  Note that with this method, the fillet is always height `a` above the base, so it makes a uniform curve parallel to the base object.  However, when the base curvature is more extreme, point `v` may end up above point `w`, resulting in an invalid configuration.  It also happens that point `v`, while below `w`, is very close to `w`, so the resulting fillet has an abrupt angle near `w` instead of a smooth transition.  
2592//   R=60;
2593//   base = R*[cos(70),sin(70)];
2594//   end = R*[cos(45),sin(45)];
2595//   tang = [-sin(45),cos(45)];
2596//   isect = line_intersection([base,back(1,base)], [end,end+tang]);
2597//   toppt = base+[0,2*PI*R*25/360];
2598//   bez = _smooth_bez_fill([toppt, isect,end], 0.8);
2599//   color("red")
2600//     stroke(bezier_curve(bez,30,endpoint=true), width=.5);
2601//   color("blue"){
2602//      stroke(arc(n=50,angle=[35,80], r=R), width=1);
2603//      stroke([base, back(40,base)]);
2604//      move(R*[cos(35),sin(35)])text("Base", size=5,anchor=BACK);
2605//      back(1)move(base+[0,40]) text("Prism", size=5, anchor=FWD);
2606//   }
2607//   color([.3,1,.3]){
2608//     right(2)move(toppt)text("w",size=5);
2609//     right(2)move(end)text("u",size=5);
2610//     stroke([isect+[1,1/4], isect+[16,4]], width=.5, endcap1="arrow2");
2611//     move([16.5,3])move(isect)text("v",size=5);
2612//     stroke([end,isect],dots=true);
2613//     stroke([isect,toppt], dots=true);
2614//   }
2615//   color("black")  {
2616//      stroke(arc(n=50, angle=[45,70], r=R-3), color="black", width=.6, endcaps="arrow2");
2617//       move( (R-10)*[cos(57.5),sin(57.5)]) text("a",size=5);
2618//      left(3)move( base+[0,PI*R*25/360]) text("a", size=5,anchor=RIGHT);
2619//      left(2)stroke( [base, toppt],endcaps="arrow2",width=.6);
2620//   }
2621// Figure(2D,Med,NoAxes): Non-Uniform fillet method.  This method differs because point `w` is found by moving the fillet distance `a` starting at the intersection point `v` instead of at the base surface.  This means that the `[u,v,w]` triple is always in the correct order to produce a valid fillet.  However, the height of the fillet above the surface will vary.  When the base concave, point `v` is below the surface of the base, which in more extreme cases can produce a fillet that goes below the base surface.  The uniform method is less likely to produce this kind of result.  When the base surface is a plane, the uniform and non-uniform methods are identical.
2622//   R=60;
2623//   base = R*[cos(70),sin(70)];
2624//   end = R*[cos(45),sin(45)];
2625//   tang = [-sin(45),cos(45)];
2626//   isect = line_intersection([base,back(1,base)], [end,end+tang]);
2627//   toppt = isect+[0,2*PI*R*25/360];
2628//   bez = _smooth_bez_fill([toppt, isect,end], 0.8);
2629//   color("red")stroke(bezier_curve(bez,30,endpoint=true), width=.5);
2630//   color("blue"){
2631//      stroke(arc(n=50,angle=[35,80], r=R), width=1);
2632//      stroke([base, back(40,base)]);
2633//      move(R*[cos(35),sin(35)])text("Base", size=5,anchor=BACK);
2634//      back(1)move(base+[0,40]) text("Prism", size=5, anchor=FWD);
2635//   }
2636//   color([.3,1,.3]){
2637//     right(2)move(toppt)text("w",size=5);
2638//     right(2)move(end)text("u",size=5);
2639//     stroke([isect+[1,1/4], isect+[16,4]], width=.5, endcap1="arrow2");
2640//     move([16.5,3])move(isect)text("v",size=5);
2641//     stroke([end,isect],dots=true);
2642//     stroke([isect,toppt], dots=true);
2643//   }
2644//   color("black")  {
2645//      stroke(arc(n=50, angle=[45,70], r=R-3), width=.6, endcaps="arrow2");
2646//      move( (R-10)*[cos(57.5),sin(57.5)]) text("a",size=5);
2647//      left(3)move( (isect+toppt)/2) text("a", size=5,anchor=RIGHT);
2648//      left(2)stroke( [isect, toppt],endcaps="arrow2",width=.6);
2649//   }
2650// Arguments:
2651//   polygon = polygon giving prism cross section
2652//   base = string specifying base object to join to ("plane","cyl","cylinder", "sphere") or a point list to use an arbitrary prism as the base.
2653//   ---
2654//   length / height / l / h = length/height of prism if aux=="none"
2655//   scale = scale factor for prism far end.  Default: 1
2656//   prism_end_T = root-centered arbitrary transform to apply to the prism's far point.  Default: IDENT
2657//   short = flip prism direction for concave sphere or cylinder base, when there are two valid prisms.  Default: false
2658//   base_T = origin-centered rotation operator to apply to the base
2659//   base_r / base_d = base radius or diameter if you picked sphere or cylinder
2660//   aux = string specifying auxilary object to connect to ("none", "plane", "cyl", "cylinder", or "sphere") or a point list to use an arbitrary prism.  Default: "none"
2661//   aux_T = rotation operator that may include translation when aux is not "none" to apply to aux
2662//   aux_r / aux_d = radius or diameter of auxiliary object if you picked sphere or cylinder
2663//   n = number of segments in the fillet at both ends.  Default: 15
2664//   base_n = number of segments to use in fillet at the base
2665//   aux_n = number of segments to use in fillet at the aux object
2666//   end_n = number of segments to use in roundover at the end of prism with no aux object
2667//   fillet = fillet for both ends of the prism (if applicable)  Must be nonnegative except for joiner prisms with planar ends
2668//   base_fillet = fillet for base end of prism 
2669//   aux_fillet = fillet for joint with aux object
2670//   end_round = roundover of end of prism with no aux object 
2671//   overlap = amount of overlap of prism fillet into objects at both ends.  Default: 1 for normal fillets, 0 for negative fillets and roundovers
2672//   base_overlap = amount of overlap of prism fillet into the base object
2673//   aux_overlap = amount of overlap of the prism fillet into aux object
2674//   k = fillet curvature parameter for both ends of prism
2675//   base_k = fillet curvature parameter for base end of prism
2676//   end_k / aux_k = fillet curvature parameter for end of prism where the aux object is
2677//   uniform = set to false to get non-uniform filleting at both ends (see Figures 2-3).  Default: true
2678//   base_uniform = set to false to get non-uniform filleting at the base
2679//   aux_uniform = set to false to get non-uniform filleting at the auxiliary object
2680//   debug = set to true to allow return of various cases where self-intersection was detected
2681//   anchor = Translate so anchor point is at the origin.  (module only) Default: "origin"
2682//   spin = Rotate this many degrees around Z axis after anchor.  (module only) Default: 0
2683//   orient = Vector to rotate top towards after spin  (module only)
2684//   atype = Select "hull" or "intersect" anchor types.  (module only) Default: "hull"
2685//   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.  (module only) Default: "centroid"
2686// Extra Anchors:
2687//   "root" = Root point of the joiner prism, pointing out in the direction of the prism axis
2688//   "end" = End point of the joiner prism, pointing out in the direction of the prism axis
2689// Example(3D,NoScales): Here is the simplest case, a circular prism with a specified length standing vertically on a plane.  
2690//   join_prism(circle(r=15,$fn=60),base="plane",
2691//              length=18, fillet=3, n=12);
2692//   cube([50,50,5],anchor=TOP);
2693// Example(3D,NoScales): Here we substitute an abitrary prism. 
2694//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2695//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2696//   join_prism(flower,base="plane",length=18, fillet=3, n=12);
2697//   cube([50,50,5],anchor=TOP);
2698// Example(3D,NoScales): Here we apply a rotation of the prism, using prism_end_T, which rotates around the prism root.  Note that aux_T will rotate around the origin, which is the same when the prism is joined to a plane.  
2699//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2700//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2701//   join_prism(flower,base="plane",length=18, fillet=3,
2702//              n=12, prism_end_T=yrot(25));
2703//   cube([50,50,5],anchor=TOP);
2704// Example(3D,NoScales): We can use `end_round` to get a roundover
2705//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2706//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2707//   join_prism(flower,base="plane",length=18, fillet=3,
2708//              n=12, prism_end_T=yrot(25), end_round=4);
2709//   cube([50,50,5],anchor=TOP);
2710// Example(3D,NoScales): We can tilt the base plane by applying a base rotation.  Note that because we did not tilt the prism, it still points upwards.  
2711//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2712//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2713//   join_prism(flower,base="plane",length=18, fillet=3,
2714//              n=12, base_T=yrot(25));
2715//   yrot(25)cube([50,50,5],anchor=TOP);
2716// Example(3D,NoScales): Next consider attaching the prism to a sphere.  You must use a circumscribed sphere to avoid a lip or gap between the sphere and prism.  Note that the prism is attached to the sphere's boundary above the origin and projects by the specified length away from the attachment point.  
2717//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2718//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2719//   join_prism(flower,base="sphere",base_r=30, length=18,
2720//              fillet=3, n=12);
2721//   spheroid(r=30,circum=true,$fn=64);
2722// Example(3D,NoScales): Rotating using the prism_end_T option rotates around the attachment point.  Note that if you rotate too far, some points of the prism will miss the sphere, which is an error.  
2723//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2724//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2725//   join_prism(flower,base="sphere",base_r=30, length=18,
2726//              fillet=3, n=12, prism_end_T=yrot(-15));
2727//   spheroid(r=30,circum=true,$fn=64);
2728// Example(3D,NoScales): Rotating using the aux_T option rotates around the origin.  You could get the same result in this case by rotating the whole model.  
2729//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2730//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2731//   join_prism(flower,base="sphere",base_r=30, length=18,
2732//              fillet=3, n=12, aux_T=yrot(-45));
2733//   spheroid(r=30,circum=true,$fn=64);
2734// Example(3D,NoScales): The origin in the prism cross section always aligns with the origin of the object you attach to.  If you want to attach off center, then shift your prism cross section.  If you shift too far so that parts of the prism miss the base object then you will get an error.  
2735//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2736//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2737//   join_prism(right(10,flower),base="sphere",base_r=30,
2738//              length=18, fillet=3, n=12);
2739//   spheroid(r=30,circum=true,$fn=64);
2740// Example(3D,NoScales): The third available base shape is the cylinder.  
2741//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2742//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2743//   join_prism(flower,base="cylinder",base_r=30,
2744//              length=18, fillet=4, n=12); 
2745//   xcyl(r=30,l=75,circum=true,$fn=64);
2746// Example(3D,NoScales): You can rotate the cylinder the same way we rotated the plane.
2747//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2748//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2749//   join_prism(flower,base="cylinder",base_r=30, length=18,
2750//              fillet=4, n=12, base_T=zrot(33)); 
2751//   zrot(33)xcyl(r=30,l=75,circum=true,$fn=64);
2752// Example(3D,NoScales): And you can rotate the prism around its attachment point with prism_end_T
2753//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2754//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2755//   join_prism(flower,base="cylinder",base_r=30, length=18,
2756//              fillet=4, n=12, prism_end_T=yrot(22));
2757//   xcyl(r=30,l=75,circum=true,$fn=64);
2758// Example(3D,NoScales): Or you can rotate the prism around the origin with aux_T
2759//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2760//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2761//   join_prism(flower,base="cylinder",base_r=30, length=18,
2762//              fillet=4, n=12, aux_T=xrot(22));
2763//   xcyl(r=30,l=75,circum=true,$fn=64);
2764// Example(3D,NoScales): Here's a prism where the scale changes
2765//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2766//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2767//   join_prism(flower,base="cylinder",base_r=30, length=18,
2768//              fillet=4, n=12,scale=.5);
2769//   xcyl(r=30,l=75,circum=true,$fn=64);
2770// Example(3D,NoScales,VPD=190,VPR=[61.3,0,69.1],VPT=[41.8956,-9.49649,4.896]): Giving a negative radius attaches to the inside of a sphere or cylinder.  Note you want the inscribed cylinder for the inner wall.  
2771//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2772//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2773//   join_prism(flower,base="cylinder",base_r=-30, length=18,
2774//              fillet=4, n=12);
2775//   bottom_half(z=-10)
2776//     tube(ir=30,wall=3,l=74,$fn=64,orient=RIGHT,anchor=CENTER);
2777// Example(3D,NoScales,VPD=140,VPR=[72.5,0,73.3],VPT=[40.961,-19.8319,-3.03302]): A hidden problem lurks with concave attachments.  The bottom of the prism does not follow the curvature of the base.  Here you can see a gap.  In some cases you can create a self-intersection in the prism.  
2778//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2779//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2780//   left_half(){
2781//     join_prism(flower,base="cylinder",base_r=-30, length=18,
2782//                fillet=4, n=12);
2783//     bottom_half(z=-10)
2784//       tube(ir=30,wall=3,l=74,$fn=64,orient=RIGHT,anchor=CENTER);
2785//   }
2786// Example(3D,NoScales,VPD=140,VPR=[72.5,0,73.3],VPT=[40.961,-19.8319,-3.03302]): The solution to both problems is to increase the overlap parameter, but you may then have excess base that must be differenced or intersected away.  In this case, an overlap of 2 is sufficient to eliminate the hole.  
2787//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2788//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2789//   left_half(){
2790//     join_prism(flower,base="cylinder",base_r=-30, length=18,
2791//                fillet=4, n=12, overlap=2);     
2792//     bottom_half(z=-10)
2793//       tube(ir=30,wall=3,l=74,$fn=64,orient=RIGHT,anchor=CENTER);
2794//   }
2795// Example(3D,NoScales,VPD=126,VPR=[76.7,0,111.1],VPT=[6.99093,2.52831,-14.8461]): Here is an example with a spherical base.  This overlap is near the minimum required to eliminate the gap, but it creates a large excess structure around the base of the prism.  
2796//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2797//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2798//   left_half(){
2799//     join_prism(flower,base="sphere",base_r=-30, length=18,
2800//                fillet=4, n=12, overlap=7);
2801//     bottom_half(z=-10) difference(){
2802//       sphere(r=33,$fn=16);
2803//       sphere(r=30,$fn=64);
2804//     }
2805//   }
2806// Example(3D,NoScales,VPD=126,VPR=[55,0,25],VPT=[1.23541,-1.80334,-16.9789]): Here is an example with a spherical base.  This overlap is near the minimum required to eliminate the gap, but it creates a large excess structure around the base of the prism.  
2807//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2808//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2809//   intersection(){
2810//     union(){
2811//       join_prism(flower,base="sphere",base_r=-30, length=18, 
2812//                  fillet=4, n=12, overlap=7);
2813//       difference(){
2814//         down(18)cuboid([68,68,30],anchor=TOP);
2815//         sphere(r=30,$fn=64);
2816//       }
2817//     }
2818//     sphere(r=33,$fn=16);
2819//   }
2820// Example(3D,NoScales,VPD=126,VPR=[55,0,25],VPT=[1.23541,-1.80334,-16.9789]): As before, rotating with aux_T rotates around the origin. 
2821//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2822//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2823//   intersection(){
2824//     union(){
2825//       join_prism(flower,base="sphere",base_r=-30, length=18,
2826//                  fillet=4, n=12, overlap=7, aux_T=yrot(13));
2827//       difference(){
2828//         down(18)cuboid([68,68,30],anchor=TOP);
2829//         sphere(r=30,$fn=64);
2830//       }
2831//     }
2832//     sphere(r=33,$fn=16);
2833//   }
2834// Example(3D,NoScales,VPD=102.06,VPR=[55,0,25],VPT=[3.96744,-2.80884,-19.9293]): Rotating with prism_end_T rotates around the attachment point.  We shrank the prism to allow a significant rotation.
2835//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2836//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2837//   intersection(){
2838//     union(){
2839//       join_prism(scale(.5,flower),base="sphere",base_r=-30,
2840//                  length=18, fillet=2, n=12, overlap=7,
2841//                  prism_end_T=yrot(25));
2842//       difference(){
2843//         down(23)cuboid([68,68,30],anchor=TOP);
2844//         sphere(r=30,$fn=64);
2845//       }
2846//     }
2847//     sphere(r=33,$fn=16);
2848//   }
2849// Example(3D,NoScales,VPR=[65.5,0,105.3],VPT=[8.36329,13.0211,9.98397],VPD=237.091): You can create a prism that crosses the inside of a cylinder or sphere by giving the same negative radius twice and leaving both objects with the same center, as shown here.  
2850//   left_half(x=7){
2851//     join_prism(circle(r=15),base="cylinder",base_r=-30, n=12,
2852//                aux="cylinder", aux_r=-30, fillet=8, overlap=3);
2853//     tube(ir=30,wall=5,l=74,$fn=64,orient=RIGHT,anchor=CENTER);     
2854//   }
2855// Example(3D,NoScales,VPR=[65.5,0,105.3],VPT=[8.36329,13.0211,9.98397],VPD=237.091): Here's a similar example with a plane for the auxiliary object.  Note that we observe the 1 unit overlap on the top surface.  
2856//   left_half(x=7){
2857//     join_prism(circle(r=15),base="cylinder",base_r=-30,
2858//                aux="plane", fillet=8, n=12, overlap=3);
2859//     tube(ir=30,wall=5,l=74,$fn=64,orient=RIGHT,anchor=CENTER);     
2860//   }
2861// Example(3D,NoScales,VPR=[65.5,0,105.3],VPT=[8.36329,13.0211,9.98397],VPD=237.091): We have tweaked the previous example just slightly by lowering the height of the plane.  The result is a bit of a surprise:  the prism flips upside down!  This happens because there is an ambiguity in creating a prism between a plane and the inside of the cylinder.  By default, this ambiguity is resolved by choosing the longer prism.  
2862//   left_half(x=7){
2863//     join_prism(circle(r=15),base="cylinder",base_r=-30, n=12,
2864//                aux="plane", aux_T=down(5), fillet=8, overlap=3);
2865//     tube(ir=30,wall=5,l=74,$fn=64,orient=RIGHT,anchor=CENTER);     
2866//   }
2867// Example(3D,NoScales,VPR=[65.5,0,105.3],VPT=[8.36329,13.0211,9.98397],VPD=237.091): Adding `short=true` resolves the ambiguity of which prism to construct in the other way, by choosing the shorter option.  
2868//   left_half(x=7){
2869//     join_prism(circle(r=15),base="cylinder",base_r=-30,
2870//                aux="plane", aux_T=down(5), fillet=8,
2871//                n=12, overlap=3, short=true);
2872//     tube(ir=30,wall=5,l=74,$fn=64,orient=RIGHT,anchor=CENTER);
2873//   }
2874// Example(3D,NoScales,VPR=[85.1,0,107.4],VPT=[8.36329,13.0211,9.98397],VPD=237.091): The problem does not arise in this case because the auxiliary object only allows one possible way to make the connection. 
2875//   left_half(x=7){
2876//     join_prism(circle(r=15),base="cylinder",base_r=-30,
2877//                aux="cylinder", aux_r=30, aux_T=up(20),
2878//                fillet=8, n=12, overlap=3);
2879//     tube(ir=30,wall=5,l=74,$fn=64,orient=RIGHT,anchor=CENTER);
2880//     up(20)xcyl(r=30,l=74,$fn=64);
2881//   }
2882// Example(3D,NoScales,VPT=[-1.23129,-3.61202,-0.249883],VPR=[87.9,0,295.7],VPD=213.382): When the aux cylinder is inside the base cylinder we can select the two options, shown here as red for the default and blue for the `short=true` case. 
2883//   color("red")
2884//     join_prism(circle(r=5),base="cylinder",base_r=-30, 
2885//                aux="cyl",aux_r=10, aux_T=up(12), fillet=4,
2886//                 n=12, overlap=3, short=false);
2887//   color("blue")
2888//     join_prism(circle(r=5),base="cylinder",base_r=-30, 
2889//                aux="cyl",aux_r=10, aux_T=up(12), fillet=4,
2890//                n=12, overlap=3, short=true);
2891//   tube(ir=30,wall=5,$fn=64,l=18,orient=RIGHT,anchor=CENTER);
2892//   up(12)xcyl(r=10, circum=true, l=18);
2893// Example(3D,NoScales,VPR=[94.9,0,106.7],VPT=[4.34503,1.48579,-2.32228],VPD=237.091): The same thing is true when you use a negative radius for the aux cylinder. This is the default long case.  
2894//   join_prism(circle(r=5,$fn=64),base="cylinder",base_r=-30, 
2895//              aux="cyl",aux_r=-10, aux_T=up(12), fillet=4,
2896//              n=12, overlap=3, short=false);
2897//   tube(ir=30,wall=5,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2898//   up(12) top_half()
2899//      tube(ir=10,wall=4,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2900// Example(3D,NoScales,VPR=[94.9,0,106.7],VPT=[4.34503,1.48579,-2.32228],VPD=237.091): And here is the short case:
2901//   join_prism(circle(r=5,$fn=64),base="cylinder",base_r=-30, 
2902//              aux="cyl",aux_r=-10, aux_T=up(12), fillet=4,
2903//              n=12, overlap=3, short=true);
2904//   tube(ir=30,l=24,wall=5,$fn=64,orient=RIGHT,anchor=CENTER);
2905//   up(12) bottom_half()
2906//     tube(ir=10,wall=4,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2907// Example(3D,NoScales,VPR=[94.9,0,106.7],VPT=[0.138465,6.78002,24.2731],VPD=325.228): Another example where the cylinders overlap, with the long case here:
2908//   auxT=up(40);
2909//   join_prism(circle(r=5,$fn=64),base="cylinder",base_r=-30, 
2910//              aux="cyl",aux_r=-40, aux_T=auxT, fillet=4,
2911//              n=12, overlap=3, short=false);
2912//   tube(ir=30,wall=4,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2913//   multmatrix(auxT)
2914//     tube(ir=40,wall=4,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2915// Example(3D,NoScales,VPR=[94.9,0,106.7],VPT=[0.138465,6.78002,24.2731],VPD=325.228): And the short case:
2916//   auxT=up(40);
2917//   join_prism(circle(r=5,$fn=64),base="cylinder",base_r=-30, 
2918//              aux="cyl",aux_r=-40, aux_T=auxT, fillet=4,
2919//              n=12, overlap=3, short=true);
2920//   tube(ir=30,wall=4,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2921//   multmatrix(auxT)
2922//     tube(ir=40,wall=4,l=24,$fn=64,orient=RIGHT,anchor=CENTER);
2923// Example(3D,NoScales): Many of the preceeding examples feature a prism with a concave shape cross section.  Concave regions can limit the amount of rounding that is possible.  This occurs because the algorithm is not able to handle a fillet that intersects itself.  Fillets on a convex prism always grow larger as they move away from the prism, so they cannot self intersect.  This means that you can make the fillet as big as will fit on the base shape.  The fillet will fail to fit if the tangent plane to the base at the fillet distance from the prism fails to intersect the prism.  Here is an extreme example, almost the largest possible fillet to the convex elliptical convex prism.  
2924//   ellipse = ellipse([17,10],$fn=164);  
2925//   join_prism(ellipse,base="sphere",base_r=30, length=18,
2926//              fillet=18, n=25, overlap=1);
2927//   sphere(r=30,circum=true, $fn=96);
2928// Example(3D,NoScales): This example shows a failed rounding attempt where the result is self-intersecting.  Using the `debug=true` option makes it possible to view the result to understand what went wrong.  Note that the concave corners have a crease where the fillet crosses itself.  The error message will advise you to decrease the size of the fillet.  You can also fix the problem by making your concave curves shallower.  
2929//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2930//             (15+2.5*sin(6*theta))*[cos(theta),sin(theta)]];
2931//   join_prism(flower,base="cylinder",base_r=30, length=18,
2932//              fillet=6, n=12, debug=true); 
2933// Example(3D,NoScales): Your prism needs to be finely sampled enough to follow the contour of the base you are attaching it to.  If it is not, you get a result like this.  The fillet joints the prism smoothly, but makes a poor transition to the sphere. 
2934//   sq = rect(15);
2935//   join_prism(sq, base="sphere", base_r=25,
2936//              length=18, fillet=4, n=12);
2937//   sphere(r=25, circum=true, $fn=96);
2938// Example(3D,NoScales): To fix the problem, you must subdivide the polygon that defines the prism.  But note that the join_prism method works poorly at sharp corners.
2939//   sq = subdivide_path(rect(15),n=64);
2940//   join_prism(sq, base="sphere", base_r=25,
2941//              length=18, fillet=4, n=12);
2942//   sphere(r=25, circum=true,$fn=96);
2943// Example(3D,NoScales): In the previous example, a small rounding of the prism corners produces a nicer result.
2944//   sq = subdivide_path(
2945//          round_corners(rect(15),cut=.5,$fn=32),
2946//          n=128);
2947//   join_prism(sq, base="sphere", base_r=25,
2948//              length=18, fillet=4, n=12);
2949//   sphere(r=25, circum=true,$fn=96);
2950// Example(3D,NoScales): The final option for specifying the base is to use an arbitrary prism, specified by a polygon.  Note that the base prism is oriented to the RIGHT, so the attached prism remains Z oriented.  
2951//   ellipse = ellipse([17,10],$fn=164);  
2952//   join_prism(zrot(90,ellipse), base=2*ellipse, length=19,
2953//              fillet=4, n=12);
2954//   linear_sweep(2*ellipse,height=60, center=true, orient=RIGHT);
2955// Example(3D,NoScales): As usual, you can rotate around the attachment point using prism_end_T. 
2956//   ellipse = ellipse([17,10],$fn=164);  
2957//   join_prism(zrot(90,ellipse), base=2*ellipse, length=19,
2958//              fillet=4, n=12, prism_end_T=yrot(22));
2959//   linear_sweep(2*ellipse,height=60, center=true, orient=RIGHT);
2960// Example(3D,NoScales): And you can rotate around the origin with aux_T.
2961//   ellipse = ellipse([17,10],$fn=164);  
2962//   join_prism(zrot(90,ellipse), base=2*ellipse, length=19,
2963//              fillet=4, n=12, aux_T=yrot(22));
2964//   linear_sweep(2*ellipse,height=60, center=true, orient=RIGHT);
2965// Example(3D,NoScales): The base prism can be a more complicated shape.
2966//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2967//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2968//   join_prism(flower,base=1.4*flower, fillet=3,
2969//              n=15, length=20);
2970//   linear_sweep(1.4*flower,height=60,center=true,
2971//                convexity=10,orient=RIGHT);
2972// Example(3D,NoScales): Here's an example with both prism_end_T and aux_T 
2973//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2974//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2975//   join_prism(flower,base=1.4*flower, length=20,
2976//              prism_end_T=yrot(20),aux_T=xrot(10),
2977//              fillet=3, n=25);
2978//   linear_sweep(1.4*flower,height=60,center=true,
2979//                convexity=10,orient=RIGHT);
2980// Example(3D,NoScales,VPR=[78,0,42],VPT=[12.45,-12.45,10.4],VPD=130): Instead of terminating your prism in a flat face perpendicular to its axis you can attach it to a second object.  The simplest case is to connect to planar attachments.  When connecting to a second object you must position and orient the second object using aux_T, which is now allowed to be a rotation and translation operator.  The `length` parameter is no longer allowed.  
2981//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2982//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2983//   join_prism(flower,base="plane", fillet=4, n=12,
2984//              aux="plane", aux_T=up(12));
2985//   %up(12)cuboid([40,40,4],anchor=BOT); 
2986//   cuboid([40,40,4],anchor=TOP);
2987// Example(3D,NoScales,VPR=[78,0,42],VPT=[12.45,-12.45,10.4],VPD=130): Here's an example where the second object is rotated.  Note that the prism will go from the origin to the origin point of the object.  In this case because the rotation is applied first, the prism is vertical.  
2988//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2989//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2990//   aux_T = up(12)*xrot(-22);
2991//   join_prism(flower,base="plane",fillet=4, n=12,
2992//              aux="plane", aux_T=aux_T); 
2993//   multmatrix(aux_T)cuboid([42,42,4],anchor=BOT);
2994//   cuboid([40,40,4],anchor=TOP);
2995// Example(3D,NoScales,VPR=[78,0,42],VPT=[12.45,-12.45,10.4],VPD=130): In this example, the aux_T transform moves the centerpoint (origin) of the aux object, and the resulting prism connects centerpoints, so it is no longer vertical. 
2996//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
2997//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
2998//   aux_T = xrot(-22)*up(12);
2999//   join_prism(flower,base="plane",fillet=4, n=12,
3000//              aux="plane", aux_T=aux_T);
3001//   multmatrix(aux_T)cuboid([42,42,4],anchor=BOT);
3002//   cuboid([43,43,4],anchor=TOP);
3003// Example(3D,NoScales,VPR=[78,0,42],VPT=[9.95,-9.98,13.0],VPD=142]): You can combine with base_T
3004//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3005//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3006//   aux_T = xrot(-22)*up(22);
3007//   base_T = xrot(5)*yrot(-12);
3008//   join_prism(flower,base="plane",base_T=base_T, 
3009//              aux="plane",aux_T=aux_T, fillet=4, n=12);
3010//   multmatrix(aux_T)cuboid([42,42,4],anchor=BOT);
3011//   multmatrix(base_T)cuboid([45,45,4],anchor=TOP);
3012// Example(3D,NoScales,VPR=[76.6,0,29.4],VPT=[11.4009,-8.43978,16.1934],VPD=157.778): Using prism_end_T shifts the prism's end without tilting the plane, so the prism ends are not perpendicular to the prism axis.  
3013//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3014//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3015//   join_prism(flower,base="plane", prism_end_T=right(14),
3016//              aux="plane",aux_T=up(24), fillet=4, n=12);
3017//   right(7){
3018//     %up(24)cuboid([65,42,4],anchor=BOT);
3019//     cuboid([65,42,4],anchor=TOP);
3020//   }
3021// Example(3D,NoAxes,NoScales,VPR=[101.9, 0, 205.6], VPT=[5.62846, -5.13283, 12.0751], VPD=102.06): Negative fillets give roundovers and are pemitted only for joints to planes.  Note that overlap defaults to zero for negative fillets.  
3022//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3023//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3024//   aux_T = xrot(-22)*up(22);
3025//   base_T = xrot(5)*yrot(-12);
3026//   join_prism(flower,base="plane",base_T=base_T,
3027//              aux="plane", aux_T=aux_T, fillet=-4,n=12);
3028// Example(3D,NoScales,VPR=[84,0,21],VPT=[13.6,-1,46.8],VPD=446): It works the same way with the other shapes, but make sure you move the shapes far enough apart that there is room for a prism.  
3029//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3030//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3031//   aux_T = up(85);
3032//   base_T = xrot(5)*yrot(-12);
3033//   join_prism(flower,base="cylinder",base_r=25, fillet=4, n=12,
3034//              aux="sphere",aux_r=35,base_T=base_T, aux_T=aux_T);
3035//   multmatrix(aux_T)sphere(35,circum=true);
3036//   multmatrix(base_T)xcyl(l=75,r=25,circum=true);
3037// Example(3D,NoScales,VPR=[84,0,21],VPT=[13.6,-1,46.8],VPD=446): Here we translate the sphere to the right and the prism goes with it
3038//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3039//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3040//   aux_T = right(40)*up(85);
3041//   join_prism(flower,base="cylinder",base_r=25, n=12,
3042//              aux="sphere",aux_r=35, aux_T=aux_T, fillet=4);
3043//   multmatrix(aux_T)sphere(35,circum=true);
3044//   xcyl(l=75,r=25,circum=true);
3045// Example(3D,NoScales,VPR=[84,0,21],VPT=[13.6,-1,46.8],VPD=446): This is the previous example with the prism_end_T transformation used to shift the far end of the prism away from the sphere center.  Note that prism_end_T can be any transformation, but it just acts on the location of the prism endpoint to shift the direction the prism points.  
3046//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3047//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3048//   aux_T = right(40)*up(85);
3049//   join_prism(flower,base="cylinder",base_r=25,
3050//              prism_end_T=left(4), fillet=3, n=12, 
3051//              aux="sphere",aux_r=35, aux_T=aux_T); 
3052//   multmatrix(aux_T)sphere(35,circum=true);
3053//   xcyl(l=75,r=25,circum=true);
3054// Example(3D,NoScales,VPR=[96.9,0,157.5],VPT=[-7.77616,-2.272,37.9424],VPD=366.527): Here the base is a cylinder but the auxilary object is a generic prism, and the joiner prism has a scale factor.  
3055//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3056//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3057//   aux_T = up(85)*zrot(-75);
3058//   ellipse = ellipse([17,10],$fn=164);  
3059//   join_prism(flower,base="cylinder",base_r=25,
3060//              fillet=4, n=12,
3061//              aux=ellipse, aux_T=aux_T,scale=.5);
3062//   multmatrix(aux_T)
3063//     linear_sweep(ellipse,orient=RIGHT,height=75,center=true);
3064//   xcyl(l=75,r=25,circum=true,$fn=100);
3065// Example(3D,NoAxes,VPT=[10.0389,1.71153,26.4635],VPR=[89.3,0,39],VPD=237.091): Base and aux are both a general prism in this case.
3066//   ellipse = ellipse([10,17]/2,$fn=96);  
3067//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3068//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3069//   aux_T=up(50);   
3070//   join_prism(ellipse,base=flower,aux_T=aux_T,aux=flower,
3071//              fillet=3, n=12, prism_end_T=right(9));
3072//   multmatrix(aux_T)
3073//     linear_sweep(flower,height=60,center=true,orient=RIGHT);
3074//   linear_sweep(flower,height=60,center=true,orient=RIGHT);
3075// Example(3D,NoAxes,VPT=[8.57543,0.531762,26.8046],VPR=[89.3,0,39],VPD=172.84): Shifting the joiner prism forward brings it close to a steeply curved edge of the auxiliary prism at the top.  Note that a funny looking bump with a sharp corner has appeared in the fillet.  This bump/corner is a result of the uniform filleting method running out of space.  If we move the joiner prism farther forward, the algorithm fails completely.  
3076//   ellipse = ellipse([10,17]/2,$fn=96);  
3077//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3078//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3079//   aux_T=up(50);   
3080//   join_prism(ellipse,base=flower,aux_T=aux_T,aux=flower,
3081//              fillet=3, n=12, prism_end_T=fwd(1.6));
3082//   multmatrix(aux_T)
3083//     linear_sweep(flower,height=60,center=true,orient=RIGHT);
3084//   linear_sweep(flower,height=60,center=true,orient=RIGHT);
3085// Example(3D,NoAxes,VPT=[8.57543,0.531762,26.8046],VPR=[89.3,0,39],VPD=172.84): This is the same example as above but with uniform turned off.  Note how the line the fillet makes on the joiner prism is not uniform, but the overall curved shape is more pleasing than the previous result, and we can bring the joiner prism a little farther forward and still construct a model. 
3086//   ellipse = ellipse([10,17]/2,$fn=96);  
3087//   flower = [for(theta=lerpn(0,360,180,endpoint=false))
3088//             (15+1.3*sin(6*theta))*[cos(theta),sin(theta)]];
3089//   aux_T=up(50);   
3090//   join_prism(ellipse,base=flower,aux_T=aux_T,aux=flower,
3091//              fillet=3, n=12, prism_end_T=fwd(1.7),
3092//              uniform=false);
3093//   multmatrix(aux_T)
3094//     linear_sweep(flower,height=60,center=true,orient=RIGHT);
3095//   linear_sweep(flower,height=60,center=true,orient=RIGHT);
3096// Example(3D): Positioning a joiner prism as an attachment
3097//   cuboid([20,30,40])
3098//     attach(RIGHT,"root")
3099//       join_prism(circle(r=8,$fn=32),
3100//                  l=10, base="plane", fillet=4);
3101module join_prism(polygon, base, base_r, base_d, base_T=IDENT,
3102                    scale=1, prism_end_T=IDENT, short=false, 
3103                    length, l, height, h,
3104                    aux="none", aux_T=IDENT, aux_r, aux_d,
3105                    overlap, base_overlap,aux_overlap,
3106                    n=15, base_n, end_n, aux_n,
3107                    fillet, base_fillet,aux_fillet,end_round,
3108                    k=0.7, base_k,aux_k,end_k,
3109                    uniform=true, base_uniform, aux_uniform, 
3110                    debug=false, anchor="origin", extent=true, cp="centroid", atype="hull", orient=UP, spin=0,
3111                    convexity=10)
3112{
3113    assert(in_list(atype, _ANCHOR_TYPES), "Anchor type must be \"hull\" or \"intersect\"");
3114    vnf_start_end = join_prism(polygon,base, base_r=base_r, base_d=base_d, base_T=base_T,
3115                   scale=scale, prism_end_T=prism_end_T, short=short,
3116                   length=length, l=l, height=height, h=h,
3117                   aux=aux, aux_T=aux_T, aux_r=aux_r, aux_d=aux_d,
3118                   overlap=overlap, base_overlap=base_overlap, aux_overlap=aux_overlap,
3119                   n=n,base_n=base_n, end_n=end_n, aux_n=aux_n,
3120                   fillet=fillet, base_fillet=base_fillet, aux_fillet=aux_fillet, end_round=end_round,
3121                   k=k, base_k=base_k, aux_k=aux_k, end_k=end_k,
3122                   uniform=uniform, base_uniform=base_uniform, aux_uniform=aux_uniform, 
3123                   debug=debug,
3124                   return_axis=true
3125    );
3126    axis = vnf_start_end[2] - vnf_start_end[1];
3127    anchors = [
3128               named_anchor("root",vnf_start_end[1], -axis),
3129               named_anchor("end",vnf_start_end[2], axis)
3130              ];
3131    attachable(anchor,spin,orient,vnf=vnf_start_end[0], extent=atype=="hull", cp=cp, anchors=anchors) {
3132      vnf_polyhedron(vnf_start_end[0],convexity=convexity);
3133      children();
3134    }
3135}
3136
3137
3138
3139function join_prism(polygon, base, base_r, base_d, base_T=IDENT,
3140                    scale=1, prism_end_T=IDENT, short=false, 
3141                    length, l, height, h,
3142                    aux="none", aux_T=IDENT, aux_r, aux_d,
3143                    overlap, base_overlap,aux_overlap,
3144                    n=15, base_n, aux_n, end_n, 
3145                    fillet, base_fillet,aux_fillet,end_round,
3146                    k=0.7, base_k,aux_k,end_k,
3147                    uniform=true, base_uniform, aux_uniform, 
3148                    debug=false, return_axis=false) =
3149  let(
3150      objects=["cyl","cylinder","plane","sphere"],
3151      length = one_defined([h,height,l,length], "h,height,l,length", dflt=undef)
3152  )
3153  assert(is_path(polygon,2),"Prism polygon must be a 2d path")
3154  assert(is_rotation(base_T,3,centered=true),"Base transformation must be a rotation around the origin")
3155  assert(is_rotation(aux_T,3),"Aux transformation must be a rotation")
3156  assert(aux!="none" || is_rotation(aux_T,centered=true), "With no aux, aux_T must be a rotation centered on the origin")
3157  assert(is_matrix(prism_end_T,4), "Prism endpoint transformation is invalid")
3158  assert(aux!="none" || (is_num(length) && length>0),"With no aux must give positive length")
3159  assert(aux=="none" || is_undef(length), "length parameter allowed only when aux is \"none\"")
3160  assert(aux=="none" || is_path(aux,2) || in_list(aux,objects), "Unknown aux type")
3161  assert(is_path(base,2) || in_list(base,objects), "Unknown base type")
3162  assert(is_undef(length) || (is_num(length) && length>0), "Prism length must be positive")
3163  assert(is_num(scale) && scale>=0, "Prism scale must be non-negative")
3164  assert(num_defined([end_k,aux_k])<2, "Cannot define both end_k and aux_k")
3165  assert(num_defined([end_n,aux_n])<2, "Cannot define both end_n and aux_n")
3166  let(
3167      base_r = get_radius(r=base_r,d=base_d),
3168      aux_r = get_radius(r=aux_r,d=aux_d),
3169      base_k= first_defined([base_k,k]),
3170      aux_k = first_defined([end_k,aux_k,k]),
3171      aux_n = first_defined([end_n,aux_n,n]),
3172      base_n = first_defined([base_n,n]),
3173      base_fillet = one_defined([fillet,base_fillet],"fillet,base_fillet"),
3174      aux_fillet = aux=="none" ? one_defined([aux_fillet,u_mul(-1,end_round)],"aux_fillet,end_round",0)
3175              : one_defined([fillet,aux_fillet],"fillet,aux_fillet"),
3176      base_overlap = one_defined([base_overlap,overlap],"base_overlap,overlap",base_fillet>0?1:0),
3177      aux_overlap = one_defined([aux_overlap,overlap],"aux_overlap,overlap",aux_fillet>0?1:0),
3178      base_uniform = first_defined([base_uniform, uniform]),
3179      aux_uniform = first_defined([aux_uniform, uniform])
3180  )
3181  assert(is_num(base_fillet),"Must give a numeric fillet or base_fillet value")
3182  assert(base=="plane" || base_fillet>=0, "Fillet for non-planar base object must be nonnegative")
3183  assert(is_num(aux_fillet), "Must give numeric fillet or aux_fillet")
3184  assert(in_list(aux,["none","plane"]) || aux_fillet>=0, "Fillet for aux object must be nonnegative")
3185  assert(!in_list(base,["sphere","cyl","cylinder"]) || (is_num(base_r) && !approx(base_r,0)), str("Must give nonzero base_r with base ",base))
3186  assert(!in_list(aux,["sphere","cyl","cylinder"]) || (is_num(aux_r) && !approx(aux_r,0)), str("Must give nonzero aux_r with base ",base))
3187  assert(!short || (in_list(base,["sphere","cyl","cylinder"]) && base_r<0), "You can only set short to true if the base is a sphere or cylinder with radius<0")
3188  let(
3189      base_r=default(base_r,0),
3190      polygon=clockwise_polygon(polygon),
3191      start_center = CENTER,
3192      dir = aux=="none" ? apply(aux_T,UP)
3193          : apply(aux_T,CENTER) == CENTER ? apply(aux_T,UP)
3194          : apply(aux_T,CENTER),
3195      flip = short ? -1 : 1,
3196      start = base=="sphere" ?
3197                let( answer = _sphere_line_isect_best(abs(base_r),[CENTER,flip*dir], sign(base_r)*flip*dir))
3198                assert(answer,"Prism center doesn't intersect sphere (base)")
3199                answer
3200            : base=="cyl" || base=="cylinder" ?
3201                let(
3202                     mapped = apply(yrot(90),[CENTER,flip*dir]),
3203                     answer = _cyl_line_intersection(abs(base_r),mapped,sign(base_r)*mapped[1])
3204                 )
3205                 assert(answer,"Prism center doesn't intersect cylinder (base)")
3206                 apply(yrot(-90),answer)
3207            : is_path(base) ?
3208                let( 
3209                     mapped = apply(yrot(90),[CENTER,flip*dir]),
3210                     answer = _prism_line_isect(pair(base,wrap=true),mapped,mapped[1])[0]
3211                 )
3212                 assert(answer,"Prism center doesn't intersect prism (base)")
3213                 apply(yrot(-90),answer)
3214            : start_center,
3215      aux_T = aux=="none" ? move(start)*prism_end_T*move(-start)*move(length*dir)*move(start)
3216              : aux_T,
3217      prism_end_T = aux=="none" ? IDENT : prism_end_T,
3218      aux = aux=="none" && aux_fillet!=0 ? "plane" : aux, 
3219      end_center = apply(aux_T,CENTER), 
3220      ndir = base_r<0 ? unit(start_center-start) : unit(end_center-start_center,UP),
3221      end_prelim = apply(move(start)*prism_end_T*move(-start),
3222            aux=="sphere" ?
3223                let( answer = _sphere_line_isect_best(abs(aux_r), [start,start+ndir], -sign(aux_r)*ndir))
3224                assert(answer,"Prism center doesn't intersect sphere (aux)")
3225                apply(aux_T,answer)
3226          : aux=="cyl" || aux=="cylinder" ? 
3227                let(
3228                     mapped = apply(yrot(90)*rot_inverse(aux_T),[start,start+ndir]),
3229                     answer = _cyl_line_intersection(abs(aux_r),mapped, -sign(aux_r)*(mapped[1]-mapped[0]))
3230                 )
3231                 assert(answer,"Prism center doesn't intersect cylinder (aux)")
3232                 apply(aux_T*yrot(-90),answer)
3233          : is_path(aux) ?
3234                let( 
3235                     mapped = apply(yrot(90),[start,start+ndir]),
3236                     answer = _prism_line_isect(pair(aux,wrap=true),mapped,mapped[0]-mapped[1])[0]
3237                 )
3238                 assert(answer,"Prism center doesn't intersect prism (aux)")
3239                 apply(aux_T*yrot(-90),answer)
3240          : end_center
3241      ),
3242      end = prism_end_T == IDENT ? end_prelim
3243          : aux=="sphere" ?
3244                let( answer = _sphere_line_isect_best(abs(aux_r), move(-end_center,[start,end_prelim]), -sign(aux_r)*(end_prelim-start)))
3245                assert(answer,"Prism center doesn't intersect sphere (aux)")
3246                answer+end_center
3247          : aux=="cyl" || aux=="cylinder" ? 
3248                let(
3249                     mapped = apply(yrot(90)*move(-end_center),[start,end_prelim]),
3250                     answer = _cyl_line_intersection(abs(aux_r),mapped, -sign(aux_r)*(mapped[1]-mapped[0]))
3251                 )
3252                 assert(answer,"Prism center doesn't intersect cylinder (aux)")
3253                 apply(move(end_center)*yrot(-90),answer)
3254          : is_path(aux) ?
3255                let( 
3256                     mapped = apply(yrot(90)*move(-end_center),[start,end_prelim]),
3257                     answer = _prism_line_isect(pair(aux,wrap=true),mapped,mapped[0]-mapped[1])[0]
3258                 )
3259                 assert(answer,"Prism center doesn't intersect prism (aux)")
3260                 apply(move(end_center)*yrot(-90),answer)
3261          : plane_line_intersection( plane_from_normal(apply(aux_T,UP), end_prelim),[start,end_prelim]),
3262      pangle = rot(from=UP, to=end-start),
3263      truetop = apply(move(start)*pangle,path3d(scale(scale,polygon),norm(start-end))),      
3264      truebot = apply(move(start)*pangle,path3d(polygon)),
3265      base_trans = rot_inverse(base_T),
3266      base_top = apply(base_trans, truetop),
3267      base_bot = apply(base_trans, truebot),
3268      botmesh = apply(base_T,_prism_fillet("base", base, base_r, base_bot, base_top, base_fillet, base_k, n, base_overlap,base_uniform,debug)),
3269      aux_trans = rot_inverse(aux_T),
3270      aux_top = apply(aux_trans, reverse_polygon(truetop)),
3271      aux_bot = apply(aux_trans, reverse_polygon(truebot)),
3272      topmesh_reversed = _prism_fillet("aux",aux, aux_r, aux_top, aux_bot, aux_fillet, aux_k, n, aux_overlap,aux_uniform,debug),
3273      topmesh = apply(aux_T,[for(i=[len(topmesh_reversed)-1:-1:0]) reverse_polygon(topmesh_reversed[i])]),
3274      round_dir = select(topmesh,-1)-botmesh[0],
3275      roundings_cross = [for(i=idx(topmesh)) if (round_dir[i]*(truetop[i]-truebot[i])<0) i],
3276      vnf = vnf_vertex_array(concat(topmesh,botmesh),col_wrap=true, caps=true, reverse=true)
3277  )
3278  assert(debug || roundings_cross==[],"Roundings from the two ends cross on the prism: decrease size of roundings")
3279  return_axis ? [vnf,start,end] : vnf;
3280
3281function _fix_angle_list(list,ind=0, result=[]) =
3282    ind==0 ? _fix_angle_list(list,1,[list[0]])
3283  : ind==len(list) ? result 
3284  : list[ind]-result[ind-1]>90 ? _fix_angle_list(list,ind+1,concat(result,[list[ind]-360]))
3285  : list[ind]-result[ind-1]<-90 ? _fix_angle_list(list,ind+1,concat(result,[list[ind]+360]))
3286  : _fix_angle_list(list,ind+1,concat(result,[list[ind]]));
3287                 
3288
3289
3290// intersection with cylinder of radius R oriented on Z axis, with infinite extent
3291// if ref is given, return point with larger inner product with ref.  
3292function _cyl_line_intersection(R, line, ref) =
3293   let(
3294       line2d = path2d(line),
3295       cisect = circle_line_intersection(r=R, cp=[0,0], line=line2d)
3296   )
3297   len(cisect)<2 ? [] :
3298   let(
3299       linevec = line2d[1]-line2d[0],
3300       dz = line[1].z-line[0].z,
3301       pts = [for(pt=cisect)
3302          let(t = (pt-line2d[0])*linevec/(linevec*linevec))  // position parameter for line
3303          [pt.x,pt.y,dz * t + line[0].z]]
3304   )
3305   is_undef(ref) ? pts :
3306   let(   
3307      dist = [for(pt=pts) ref*pt]
3308   )
3309   dist[0]>dist[1] ? pts[0] : pts[1];
3310
3311
3312function _sphere_line_isect_best(R, line, ref) =
3313   let(
3314        pts = sphere_line_intersection(abs(R), [0,0,0], line=line)
3315   )
3316   len(pts)<2 ? [] :
3317   let(  
3318        dist = [for(pt=pts) ref*pt]
3319   )
3320   dist[0]>dist[1] ? pts[0] : pts[1];
3321
3322// First input is all the pairs of the polygon, e.g. pair(poly,wrap=true)
3323// Unlike the others this returns [point, ind, u], where point is the actual intersection
3324// point, ind ind and u are the segment index and u value.  Prism is z-aligned.  
3325function _prism_line_isect(poly_pairs, line, ref) =
3326   let(
3327       line2d = path2d(line),
3328       ref=point2d(ref),
3329       ilist = [for(j=idx(poly_pairs)) 
3330                 let(segisect = _general_line_intersection(poly_pairs[j],line2d))
3331                 if (segisect && segisect[1]>=-EPSILON && segisect[1]<=1+EPSILON)
3332                    [segisect[0],j,segisect[1],segisect[0]*ref]]
3333   )
3334   len(ilist)==0 ? [] :
3335   let (
3336       ind = max_index(column(ilist,3)),
3337       isect2d = ilist[ind][0],
3338       isect_ind = ilist[ind][1],
3339       isect_u = ilist[ind][2],
3340       slope = (line[1].z-line[0].z)/norm(line[1]-line[0]),
3341       z = slope * norm(line2d[0]-isect2d) + line[0].z
3342   )
3343   [point3d(isect2d,z),isect_ind, isect_u];
3344
3345  
3346function _prism_fillet(name, base, R, bot, top, d, k, N, overlap,uniform,debug) =
3347    base=="none" ? [bot] 
3348  : base=="plane" ? _prism_fillet_plane(name,bot, top, d, k, N, overlap,debug)
3349  : base=="cyl" || base=="cylinder" ? _prism_fillet_cyl(name, R, bot, top, d, k, N, overlap,uniform,debug)
3350  : base=="sphere" ? _prism_fillet_sphere(name, R, bot, top, d, k, N, overlap,uniform,debug)
3351  : is_path(base,2) ? _prism_fillet_prism(name, base, bot, top, d, k, N, overlap,uniform,debug)
3352  : assert(false,"Unknown base type");
3353
3354function _prism_fillet_plane(name, bot, top, d, k, N, overlap,debug) = 
3355    let(
3356        dir = sign(top[0].z-bot[0].z),
3357        isect = [for (i=idx(top)) plane_line_intersection([0,0,1,0], [top[i],bot[i]])],
3358        base_normal = -path3d(path_normals(path2d(isect), closed=true)),
3359        mesh = transpose([for(i=idx(top))
3360          let(
3361              
3362              base_angle = vector_angle(top[i],isect[i],isect[i]+sign(d)*base_normal[i]),
3363              // joint length
3364              // d = r,
3365              r=abs(d)*tan(base_angle/2),
3366              // radius
3367              //d = r/tan(base_angle/2),
3368              // cut
3369              //r = r / (1/sin(base_angle/2) - 1),
3370              //d = r/tan(base_angle/2),
3371              prev = unit(top[i]-isect[i]),
3372              next = sign(d)*dir*base_normal[i],
3373              center = r/sin(base_angle/2) * unit(prev+next) + isect[i]
3374          )
3375          [
3376            each arc(N, cp=center, points = [isect[i]+prev*abs(d), isect[i]+next*d]),
3377            isect[i]+next*d+[0,0,-overlap*dir]
3378          ]
3379        ])
3380    )
3381    assert(debug || is_path_simple(path2d(select(mesh,-2)),closed=true),"Fillet doesn't fit: it intersects itself")
3382    mesh;
3383
3384function _prism_fillet_plane(name, bot, top, d, k, N, overlap,debug) = 
3385    let(
3386        dir = sign(top[0].z-bot[0].z),    // Negative if we are upside down, with "top" below "bot"
3387        isect = [for (i=idx(top)) plane_line_intersection([0,0,1,0], [top[i],bot[i]])]
3388    )
3389    d==0 ? [isect, if (overlap!=0) isect + overlap*dir*DOWN] :
3390    let(
3391        base_normal = -path3d(path_normals(path2d(isect), closed=true)),
3392        mesh = transpose([for(i=idx(top))
3393          assert(norm(top[i]-isect[i])>=d,"Prism is too short for fillet to fit")
3394          let(
3395              d_step = isect[i]+abs(d)*unit(top[i]-isect[i]),
3396              edgepoint = isect[i]+d*dir*base_normal[i],
3397              bez = _smooth_bez_fill([d_step, isect[i], edgepoint],k)
3398          )
3399          [
3400            each bezier_curve(bez,N,endpoint=true),
3401            if (overlap!=0) edgepoint + overlap*dir*DOWN
3402          ]
3403        ])
3404    )
3405    assert(debug || is_path_simple(path2d(select(mesh,-2)),closed=true),"Fillet doesn't fit: it intersects itself")
3406    mesh;
3407
3408
3409// This function was written for a z-aligned cylinder but the actual
3410// upstream assumption is an x-aligned cylinder, so input is rotated and
3411// output is un-rotated.  
3412function _prism_fillet_cyl(name, R, bot, top, d, k, N, overlap, uniform, debug) =
3413    let(
3414        top = yrot(-90,top),
3415        bot = yrot(-90,bot),
3416        isect = [for (i=idx(top))
3417                   let (cisect = _cyl_line_intersection(abs(R), [top[i],bot[i]], sign(R)*(top[i]-bot[i])))
3418                   assert(cisect, str("Prism doesn't fully intersect cylinder (",name,")"))
3419                   cisect
3420                ]
3421    )
3422    d==0 ? [ 
3423             isect,
3424             if (overlap!=0) [for(p=isect) point3d(unit(point2d(p))*(norm(point2d(p))-sign(R)*overlap),p.z)]
3425           ] :
3426    let(
3427        tangent = path_tangents(isect,closed=true),
3428        mesh = transpose([for(i=idx(top))
3429           assert(norm(top[i]-isect[i])>=d,str("Prism is too short for fillet to fit (",name,")"))
3430           let(
3431               dir = sign(R)*unit(cross([isect[i].x,isect[i].y,0],tangent[i])),
3432               zpart = d*dir.z,
3433               curvepart = d*norm(point2d(dir)),
3434               curveang = sign(cross(point2d(isect[i]),point2d(dir))) * curvepart * 180 / PI / abs(R), 
3435               edgepoint = apply(up(zpart)*zrot(curveang), isect[i]),
3436               corner = plane_line_intersection(plane_from_normal([edgepoint.x,edgepoint.y,0], edgepoint),
3437                                                [isect[i],top[i]],
3438                                                bounded=false/*[R>0,true]*/),
3439               d_step = abs(d)*unit(top[i]-isect[i])+(uniform?isect[i]:corner)
3440           )
3441           assert(is_vector(corner,3),str("Fillet does not fit.  Decrease size of fillet (",name,")."))
3442           assert(debug || R<0 || (d_step-corner)*(corner-isect[i])>=0,
3443                 str("Unable to fit fillet, probably due to steep curvature of the cylinder (",name,")."))
3444           let(
3445                bez = _smooth_bez_fill([d_step,corner,edgepoint], k)
3446           )
3447           [ 
3448             each bezier_curve(bez, N, endpoint=true),
3449             if (overlap!=0) point3d(unit(point2d(edgepoint))*(norm(point2d(edgepoint))-sign(R)*overlap),edgepoint.z)
3450           ]
3451        ]),
3452        angle_list = _fix_angle_list([for(pt=select(mesh,-2)) atan2(pt.y,pt.x)]),
3453        z_list = [for(pt=select(mesh,-2)) pt.z],
3454        is_simple = debug || is_path_simple(hstack([angle_list,z_list]), closed=true)
3455    )
3456    assert(is_simple, str("Fillet doesn't fit: its edge is self-intersecting.  Decrease size of roundover. (",name,")"))
3457    yrot(90,mesh);
3458
3459
3460
3461function _prism_fillet_sphere(name, R,bot, top, d, k, N, overlap, uniform, debug) = 
3462    let(
3463        isect = [for (i=idx(top))
3464                    let( isect_pt = _sphere_line_isect_best(abs(R), [top[i],bot[i]],sign(R)*(top[i]-bot[i])))
3465                    assert(isect_pt, str("Prism doesn't fully intersect sphere (",name,")"))
3466                    isect_pt
3467                ]
3468    )
3469    d==0 ? [isect,
3470            if (overlap!=0) [for(p=isect) p - overlap*sign(R)*unit(p)]
3471           ] :
3472    let(          
3473        tangent = path_tangents(isect,closed=true),
3474        mesh = transpose([for(i=idx(top))
3475           assert(norm(top[i]-isect[i])>=d,str("Prism is too short for fillet to fit (",name,")"))
3476           let(   
3477               dir = sign(R)*unit(cross(isect[i],tangent[i])),
3478               curveang = d * 180 / PI / R,
3479               edgepoint = rot(-curveang,v=tangent[i],p=isect[i]),
3480               corner = plane_line_intersection(plane_from_normal(edgepoint, edgepoint),
3481                                                [isect[i],top[i]],
3482                                                bounded=[R>0,true]),
3483               d_step = d*unit(top[i]-isect[i])+(uniform?isect[i]:corner)
3484           ) 
3485           assert(is_vector(corner,3),str("Fillet does not fit (",name,")"))
3486           assert(debug || R<0 || (d_step-corner)*(corner-isect[i])>0, 
3487                  str("Unable to fit fillet, probably due to steep curvature of the sphere (",name,")."))
3488           let(
3489               bez = _smooth_bez_fill([d_step,corner,edgepoint], k)         
3490           ) 
3491           [ 
3492             each bezier_curve(bez, N, endpoint=true),
3493             if (overlap!=0) edgepoint - overlap*sign(R)*unit(edgepoint)
3494           ]
3495        ])
3496      )
3497      // this test will fail if the prism isn't "vertical".  Project along prism direction?  
3498      assert(debug || is_path_simple(path2d(select(mesh,-2)),closed=true),str("Fillet doesn't fit: it intersects itself (",name,")"))
3499      mesh;
3500
3501
3502
3503// Return an interpolated normal to the polygon at segment i, fraction u along the segment.
3504
3505function _getnormal(polygon,index,u,) =
3506  let(
3507      //flat=1/3,
3508      flat=1/8,
3509//     flat=0,
3510      edge = (1-flat)/2,
3511      L=len(polygon),
3512      next_ind = posmod(index+1,L),
3513      prev_ind = posmod(index-1,L),
3514      this_normal = line_normal(select(polygon,index,index+1))
3515  )
3516    u > 1-edge ? lerp(this_normal,line_normal(select(polygon,index+1,index+2)), (u-edge-flat)/edge/2)
3517  : u < edge ? lerp(line_normal(select(polygon,index-1,index)),this_normal, 0.5+u/edge/2)
3518  : this_normal;
3519
3520
3521// Start at segment ind, position u on the polygon and find a point length units
3522// from that starting point.  If dir<0 goes backwards through polygon segments
3523// and if dir>0 goes forwards through polygon segments.
3524// Returns [ point, ind, u] where point is the actual point desired.  
3525function _polygon_step(poly, ind, u, dir, length) =
3526    let(ind = posmod(ind,len(poly)))
3527    u==0 && dir<0 ? _polygon_step(poly, ind-1, 1, dir, length)
3528  : u==1 && dir>0 ? _polygon_step(poly, ind+1, 0, dir, length)
3529  : let(
3530        seg = select(poly,ind,ind+1),
3531        seglen = norm(seg[1]-seg[0]),
3532        frac_needed = length / seglen
3533    )
3534    dir>0 ?
3535            ( (1-u) < frac_needed ? _polygon_step(poly,ind+1,0,dir,length-(1-u)*seglen)
3536                                 : [lerp(seg[0],seg[1],u+frac_needed),ind,u+frac_needed]
3537            )
3538          :
3539            ( u < frac_needed ? _polygon_step(poly,ind-1,1,dir,length-u*seglen)
3540                                 : [lerp(seg[0],seg[1],u-frac_needed),ind,u-frac_needed]
3541            );
3542
3543
3544// This function needs more error checking?
3545// Needs check for zero overlap case and zero joint case
3546function _prism_fillet_prism(name, basepoly, bot, top, d, k, N, overlap, uniform, debug)=
3547    let(
3548         top = yrot(-90,top),
3549         bot = yrot(-90,bot),
3550         basepoly = clockwise_polygon(basepoly),
3551         segpairs = pair(basepoly,wrap=true),
3552         isect_ind = [for (i=idx(top))
3553                         let(isect = _prism_line_isect(segpairs, [top[i], bot[i]], top[i]))
3554                         assert(isect, str("Prism doesn't fully intersect prism (",name,")"))
3555                         isect
3556                     ],
3557         isect=column(isect_ind,0),
3558         index = column(isect_ind,1),
3559         uval = column(isect_ind,2),
3560         tangent = path_tangents(isect,closed=true),
3561         mesh = transpose([for(i=idx(top))
3562           let(
3563               normal = point3d(_getnormal(basepoly,index[i],uval[i])),
3564               dir = unit(cross(normal,tangent[i])),
3565               zpart = d*dir.z,
3566               length_needed = d*norm(point2d(dir)),
3567               edgept2d = _polygon_step(basepoly, index[i], uval[i], sign(cross(point2d(dir),point2d(normal))), length_needed),
3568               edgepoint = point3d(edgept2d[0],isect[i].z+zpart),
3569               corner = plane_line_intersection(plane_from_normal(point3d(_getnormal(basepoly, edgept2d[1],edgept2d[2])),edgepoint),
3570                                                [top[i],isect[i]],
3571                                                bounded=false), // should be true!!!  But fails to intersect if given true.
3572               d_step = abs(d)*unit(top[i]-isect[i])+(uniform?isect[i]:corner)
3573           )
3574           assert(is_vector(corner,3),str("Fillet does not fit.  Decrease size of fillet (",name,")."))
3575           assert(debug  || (top[i]-d_step)*(d_step-corner)>=0,
3576                   str("Unable to fit fillet, probably due to steep curvature of the prism (",name,").",
3577                     d_step," ",corner," ", edgepoint," ", isect[i]
3578                     ))
3579           let(
3580                bez = _smooth_bez_fill([d_step,corner,edgepoint], k)
3581           )
3582           [ 
3583             each bezier_curve(bez, N, endpoint=true),
3584             if (overlap!=0) edgepoint-point3d(normal)*overlap
3585           ]
3586          ])
3587         )
3588        yrot(90,mesh);
3589
3590
3591// vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap