1//////////////////////////////////////////////////////////////////////
2// LibFile: beziers.scad
3// Bezier curves and surfaces are way to represent smooth curves and smoothly curving
4// surfaces with a set of control points. The curve or surface is defined by
5// the control points, but usually only passes through the first and last control point (the endpoints).
6// This file provides some
7// aids to constructing the control points, and highly optimized functions for
8// computing the Bezier curves and surfaces given by the control points,
9// Includes:
10// include <BOSL2/std.scad>
11// include <BOSL2/beziers.scad>
12// FileGroup: Advanced Modeling
13// FileSummary: Bezier curves and surfaces.
14//////////////////////////////////////////////////////////////////////
15
16// Terminology:
17// Path = A series of points joined by straight line segements.
18// Bezier Curve = A polynomial curve defined by a list of control points. The curve starts at the first control point and ends at the last one. The other control points define the shape of the curve and they are often *NOT* on the curve
19// Control Point = A point that influences the shape of the Bezier curve.
20// Degree = The degree of the polynomial used to make the bezier curve. A bezier curve of degree N will have N+1 control points. Most beziers are cubic (degree 3). The higher the degree, the more the curve can wiggle.
21// Bezier Parameter = A parameter, usually `u` below, that ranges from 0 to 1 to trace out the bezier curve. When `u=0` you get the first control point and when `u=1` you get the last control point. Intermediate points are traced out *non-uniformly*.
22// Bezier Path = A list of bezier control points corresponding to a series of Bezier curves that connect together, end to end. Because they connect, the endpoints are shared between control points and are not repeated, so a degree 3 bezier path representing two bezier curves will have seven entries to represent two sets of four control points. **NOTE:** A "bezier path" is *NOT* a standard path
23// Bezier Patch = A two-dimensional arrangement of Bezier control points that generate a bounded curved Bezier surface. A Bezier patch is a (N+1) by (M+1) grid of control points, which defines surface with four edges (in the non-degenerate case).
24// Bezier Surface = A surface defined by a list of one or more bezier patches.
25// Spline Steps = The number of straight-line segments used to approximate a Bezier curve. The more spline steps, the better the approximation to the curve, but the slower it will be to generate. This plays a role analogous to `$fn` for circles. Usually defaults to 16.
26
27
28// Section: Bezier Curves
29
30// Function: bezier_points()
31// Usage:
32// pt = bezier_points(bezier, u);
33// ptlist = bezier_points(bezier, RANGE);
34// ptlist = bezier_points(bezier, LIST);
35// Topics: Bezier Curves
36// Description:
37// Computes points on a bezier curve with control points specified by `bezier` at parameter values
38// specified by `u`, which can be a scalar or a list. The value `u=0` gives the first endpoint; `u=1` gives the final endpoint,
39// and intermediate values of `u` fill in the curve in a non-uniform fashion. This function uses an optimized method which
40// is best when `u` is a long list and the bezier degree is 10 or less. The degree of the bezier
41// curve is `len(bezier)-1`.
42// Arguments:
43// bezier = The list of endpoints and control points for this bezier curve.
44// u = Parameter values for evaluating the curve, given as a single value, a list or a range.
45// Example(2D): Quadratic (Degree 2) Bezier.
46// bez = [[0,0], [30,30], [80,0]];
47// debug_bezier(bez, N=len(bez)-1);
48// translate(bezier_points(bez, 0.3)) color("red") sphere(1);
49// Example(2D): Cubic (Degree 3) Bezier
50// bez = [[0,0], [5,35], [60,-25], [80,0]];
51// debug_bezier(bez, N=len(bez)-1);
52// translate(bezier_points(bez, 0.4)) color("red") sphere(1);
53// Example(2D): Degree 4 Bezier.
54// bez = [[0,0], [5,15], [40,20], [60,-15], [80,0]];
55// debug_bezier(bez, N=len(bez)-1);
56// translate(bezier_points(bez, 0.8)) color("red") sphere(1);
57// Example(2D): Giving a List of `u`
58// bez = [[0,0], [5,35], [60,-25], [80,0]];
59// debug_bezier(bez, N=len(bez)-1);
60// pts = bezier_points(bez, [0, 0.2, 0.3, 0.7, 0.8, 1]);
61// rainbow(pts) move($item) sphere(1.5, $fn=12);
62// Example(2D): Giving a Range of `u`
63// bez = [[0,0], [5,35], [60,-25], [80,0]];
64// debug_bezier(bez, N=len(bez)-1);
65// pts = bezier_points(bez, [0:0.2:1]);
66// rainbow(pts) move($item) sphere(1.5, $fn=12);
67
68// Ugly but speed optimized code for computing bezier curves using the matrix representation
69// See https://pomax.github.io/bezierinfo/#matrix for explanation.
70//
71// All of the loop unrolling makes and the use of the matrix lookup table make a big difference
72// in the speed of execution. For orders 10 and below this code is 10-20 times faster than
73// the recursive code using the de Casteljau method depending on the bezier order and the
74// number of points evaluated in one call (more points is faster). For orders 11 and above without the
75// lookup table or hard coded powers list the code is about twice as fast as the recursive method.
76// Note that everything I tried to simplify or tidy this code made is slower, sometimes a lot slower.
77function bezier_points(curve, u) =
78 is_num(u) ? bezier_points(curve,[u])[0] :
79 let(
80 N = len(curve)-1,
81 M = _bezier_matrix(N)*curve
82 )
83 N==0 ? [for(uval=u)[1]*M] :
84 N==1 ? [for(uval=u)[1, uval]*M] :
85 N==2 ? [for(uval=u)[1, uval, uval^2]*M] :
86 N==3 ? [for(uval=u)[1, uval, uval^2, uval^3]*M] :
87 N==4 ? [for(uval=u)[1, uval, uval^2, uval^3, uval^4]*M] :
88 N==5 ? [for(uval=u)[1, uval, uval^2, uval^3, uval^4, uval^5]*M] :
89 N==6 ? [for(uval=u)[1, uval, uval^2, uval^3, uval^4, uval^5,uval^6]*M] :
90 N==7 ? [for(uval=u)[1, uval, uval^2, uval^3, uval^4, uval^5,uval^6, uval^7]*M] :
91 N==8 ? [for(uval=u)[1, uval, uval^2, uval^3, uval^4, uval^5,uval^6, uval^7, uval^8]*M] :
92 N==9 ? [for(uval=u)[1, uval, uval^2, uval^3, uval^4, uval^5,uval^6, uval^7, uval^8, uval^9]*M] :
93 N==10? [for(uval=u)[1, uval, uval^2, uval^3, uval^4, uval^5,uval^6, uval^7, uval^8, uval^9, uval^10]*M] :
94 /* N>=11 */ [for(uval=u)[for (i=[0:1:N]) uval^i]*M];
95
96
97// Not public.
98function _signed_pascals_triangle(N,tri=[[-1]]) =
99 len(tri)==N+1 ? tri :
100 let(last=tri[len(tri)-1])
101 _signed_pascals_triangle(N,concat(tri,[[-1, for(i=[0:1:len(tri)-2]) (i%2==1?-1:1)*(abs(last[i])+abs(last[i+1])),len(last)%2==0? -1:1]]));
102
103
104// Not public.
105function _compute_bezier_matrix(N) =
106 let(tri = _signed_pascals_triangle(N))
107 [for(i=[0:N]) concat(tri[N][i]*tri[i], repeat(0,N-i))];
108
109
110// The bezier matrix, which is related to Pascal's triangle, enables nonrecursive computation
111// of bezier points. This method is much faster than the recursive de Casteljau method
112// in OpenScad, but we have to precompute the matrices to reap the full benefit.
113
114// Not public.
115_bezier_matrix_table = [
116 [[1]],
117 [[ 1, 0],
118 [-1, 1]],
119 [[1, 0, 0],
120 [-2, 2, 0],
121 [1, -2, 1]],
122 [[ 1, 0, 0, 0],
123 [-3, 3, 0, 0],
124 [ 3,-6, 3, 0],
125 [-1, 3,-3, 1]],
126 [[ 1, 0, 0, 0, 0],
127 [-4, 4, 0, 0, 0],
128 [ 6,-12, 6, 0, 0],
129 [-4, 12,-12, 4, 0],
130 [ 1, -4, 6,-4, 1]],
131 [[ 1, 0, 0, 0, 0, 0],
132 [ -5, 5, 0, 0, 0, 0],
133 [ 10,-20, 10, 0, 0, 0],
134 [-10, 30,-30, 10, 0, 0],
135 [ 5,-20, 30,-20, 5, 0],
136 [ -1, 5,-10, 10,-5, 1]],
137 [[ 1, 0, 0, 0, 0, 0, 0],
138 [ -6, 6, 0, 0, 0, 0, 0],
139 [ 15,-30, 15, 0, 0, 0, 0],
140 [-20, 60,-60, 20, 0, 0, 0],
141 [ 15,-60, 90,-60, 15, 0, 0],
142 [ -6, 30,-60, 60,-30, 6, 0],
143 [ 1, -6, 15,-20, 15,-6, 1]],
144 [[ 1, 0, 0, 0, 0, 0, 0, 0],
145 [ -7, 7, 0, 0, 0, 0, 0, 0],
146 [ 21, -42, 21, 0, 0, 0, 0, 0],
147 [-35, 105,-105, 35, 0, 0, 0, 0],
148 [ 35,-140, 210,-140, 35, 0, 0, 0],
149 [-21, 105,-210, 210,-105, 21, 0, 0],
150 [ 7, -42, 105,-140, 105,-42, 7, 0],
151 [ -1, 7, -21, 35, -35, 21,-7, 1]],
152 [[ 1, 0, 0, 0, 0, 0, 0, 0, 0],
153 [ -8, 8, 0, 0, 0, 0, 0, 0, 0],
154 [ 28, -56, 28, 0, 0, 0, 0, 0, 0],
155 [-56, 168,-168, 56, 0, 0, 0, 0, 0],
156 [ 70,-280, 420,-280, 70, 0, 0, 0, 0],
157 [-56, 280,-560, 560,-280, 56, 0, 0, 0],
158 [ 28,-168, 420,-560, 420,-168, 28, 0, 0],
159 [ -8, 56,-168, 280,-280, 168,-56, 8, 0],
160 [ 1, -8, 28, -56, 70, -56, 28,-8, 1]],
161 [[1, 0, 0, 0, 0, 0, 0, 0, 0, 0], [-9, 9, 0, 0, 0, 0, 0, 0, 0, 0], [36, -72, 36, 0, 0, 0, 0, 0, 0, 0], [-84, 252, -252, 84, 0, 0, 0, 0, 0, 0],
162 [126, -504, 756, -504, 126, 0, 0, 0, 0, 0], [-126, 630, -1260, 1260, -630, 126, 0, 0, 0, 0], [84, -504, 1260, -1680, 1260, -504, 84, 0, 0, 0],
163 [-36, 252, -756, 1260, -1260, 756, -252, 36, 0, 0], [9, -72, 252, -504, 630, -504, 252, -72, 9, 0], [-1, 9, -36, 84, -126, 126, -84, 36, -9, 1]],
164 [[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [-10, 10, 0, 0, 0, 0, 0, 0, 0, 0, 0], [45, -90, 45, 0, 0, 0, 0, 0, 0, 0, 0], [-120, 360, -360, 120, 0, 0, 0, 0, 0, 0, 0],
165 [210, -840, 1260, -840, 210, 0, 0, 0, 0, 0, 0], [-252, 1260, -2520, 2520, -1260, 252, 0, 0, 0, 0, 0],
166 [210, -1260, 3150, -4200, 3150, -1260, 210, 0, 0, 0, 0], [-120, 840, -2520, 4200, -4200, 2520, -840, 120, 0, 0, 0],
167 [45, -360, 1260, -2520, 3150, -2520, 1260, -360, 45, 0, 0], [-10, 90, -360, 840, -1260, 1260, -840, 360, -90, 10, 0],
168 [1, -10, 45, -120, 210, -252, 210, -120, 45, -10, 1]]
169];
170
171
172// Not public.
173function _bezier_matrix(N) =
174 N>10 ? _compute_bezier_matrix(N) :
175 _bezier_matrix_table[N];
176
177
178// Function: bezier_curve()
179// Usage:
180// path = bezier_curve(bezier, [splinesteps], [endpoint]);
181// Topics: Bezier Curves
182// See Also: bezier_curvature(), bezier_tangent(), bezier_derivative(), bezier_points()
183// Description:
184// Takes a list of bezier control points and generates splinesteps segments (splinesteps+1 points)
185// along the bezier curve they define.
186// Points start at the first control point and are sampled uniformly along the bezier parameter.
187// The endpoints of the output will be *exactly* equal to the first and last bezier control points
188// when endpoint is true. If endpoint is false the sampling stops one step before the final point
189// of the bezier curve, but you still get the same number of (more tightly spaced) points.
190// The distance between the points will *not* be equidistant.
191// The degree of the bezier curve is one less than the number of points in `curve`.
192// Arguments:
193// bezier = The list of control points that define the Bezier curve.
194// splinesteps = The number of segments to create on the bezier curve. Default: 16
195// endpoint = if false then exclude the endpoint. Default: True
196// Example(2D): Quadratic (Degree 2) Bezier.
197// bez = [[0,0], [30,30], [80,0]];
198// move_copies(bezier_curve(bez, 8)) sphere(r=1.5, $fn=12);
199// debug_bezier(bez, N=len(bez)-1);
200// Example(2D): Cubic (Degree 3) Bezier
201// bez = [[0,0], [5,35], [60,-25], [80,0]];
202// move_copies(bezier_curve(bez, 8)) sphere(r=1.5, $fn=12);
203// debug_bezier(bez, N=len(bez)-1);
204// Example(2D): Degree 4 Bezier.
205// bez = [[0,0], [5,15], [40,20], [60,-15], [80,0]];
206// move_copies(bezier_curve(bez, 8)) sphere(r=1.5, $fn=12);
207// debug_bezier(bez, N=len(bez)-1);
208function bezier_curve(bezier,splinesteps=16,endpoint=true) =
209 bezier_points(bezier, lerpn(0,1,splinesteps+1,endpoint));
210
211
212// Function: bezier_derivative()
213// Usage:
214// deriv = bezier_derivative(bezier, u, [order]);
215// derivs = bezier_derivative(bezier, LIST, [order]);
216// derivs = bezier_derivative(bezier, RANGE, [order]);
217// Topics: Bezier Curves
218// See Also: bezier_curvature(), bezier_tangent(), bezier_points()
219// Description:
220// Evaluates the derivative of the bezier curve at the given parameter value or values, `u`. The `order` gives the order of the derivative.
221// The degree of the bezier curve is one less than the number of points in `bezier`.
222// Arguments:
223// bezier = The list of control points that define the Bezier curve.
224// u = Parameter values for evaluating the curve, given as a single value, a list or a range.
225// order = The order of the derivative to return. Default: 1 (for the first derivative)
226function bezier_derivative(bezier, u, order=1) =
227 assert(is_int(order) && order>=0)
228 order==0? bezier_points(bezier, u) : let(
229 N = len(bezier) - 1,
230 dpts = N * deltas(bezier)
231 ) order==1? bezier_points(dpts, u) :
232 bezier_derivative(dpts, u, order-1);
233
234
235
236// Function: bezier_tangent()
237// Usage:
238// tanvec = bezier_tangent(bezier, u);
239// tanvecs = bezier_tangent(bezier, LIST);
240// tanvecs = bezier_tangent(bezier, RANGE);
241// Topics: Bezier Curves
242// See Also: bezier_curvature(), bezier_derivative(), bezier_points()
243// Description:
244// Returns the unit tangent vector at the given parameter values on a bezier curve with control points `bezier`.
245// Arguments:
246// bezier = The list of control points that define the Bezier curve.
247// u = Parameter values for evaluating the curve, given as a single value, a list or a range.
248function bezier_tangent(bezier, u) =
249 let(
250 res = bezier_derivative(bezier, u)
251 ) is_vector(res)? unit(res) :
252 [for (v=res) unit(v)];
253
254
255
256// Function: bezier_curvature()
257// Usage:
258// crv = bezier_curvature(curve, u);
259// crvlist = bezier_curvature(curve, LIST);
260// crvlist = bezier_curvature(curve, RANGE);
261// Topics: Bezier Curves
262// See Also: bezier_tangent(), bezier_derivative(), bezier_points()
263// Description:
264// Returns the curvature value for the given parameters `u` on the bezier curve with control points `bezier`.
265// The curvature is the inverse of the radius of the tangent circle at the given point.
266// Thus, the tighter the curve, the larger the curvature value. Curvature will be 0 for
267// a position with no curvature, since 1/0 is not a number.
268// Arguments:
269// bezier = The list of control points that define the Bezier curve.
270// u = Parameter values for evaluating the curve, given as a single value, a list or a range.
271function bezier_curvature(bezier, u) =
272 is_num(u) ? bezier_curvature(bezier,[u])[0] :
273 let(
274 d1 = bezier_derivative(bezier, u, 1),
275 d2 = bezier_derivative(bezier, u, 2)
276 ) [
277 for(i=idx(d1))
278 sqrt(
279 sqr(norm(d1[i])*norm(d2[i])) -
280 sqr(d1[i]*d2[i])
281 ) / pow(norm(d1[i]),3)
282 ];
283
284
285
286// Function: bezier_closest_point()
287// Usage:
288// u = bezier_closest_point(bezier, pt, [max_err]);
289// Topics: Bezier Curves
290// See Also: bezier_points()
291// Description:
292// Finds the closest part of the given bezier curve to point `pt`.
293// The degree of the curve, N, is one less than the number of points in `curve`.
294// Returns `u` for the closest position on the bezier curve to the given point `pt`.
295// Arguments:
296// bezier = The list of control points that define the Bezier curve.
297// pt = The point to find the closest curve point to.
298// max_err = The maximum allowed error when approximating the closest approach.
299// Example(2D):
300// pt = [40,15];
301// bez = [[0,0], [20,40], [60,-25], [80,0]];
302// u = bezier_closest_point(bez, pt);
303// debug_bezier(bez, N=len(bez)-1);
304// color("red") translate(pt) sphere(r=1);
305// color("blue") translate(bezier_points(bez,u)) sphere(r=1);
306function bezier_closest_point(bezier, pt, max_err=0.01, u=0, end_u=1) =
307 let(
308 steps = len(bezier)*3,
309 uvals = [u, for (i=[0:1:steps]) (end_u-u)*(i/steps)+u, end_u],
310 path = bezier_points(bezier,uvals),
311 minima_ranges = [
312 for (i = [1:1:len(uvals)-2]) let(
313 d1 = norm(path[i-1]-pt),
314 d2 = norm(path[i ]-pt),
315 d3 = norm(path[i+1]-pt)
316 ) if (d2<=d1 && d2<=d3) [uvals[i-1],uvals[i+1]]
317 ]
318 ) len(minima_ranges)>1? (
319 let(
320 min_us = [
321 for (minima = minima_ranges)
322 bezier_closest_point(bezier, pt, max_err=max_err, u=minima.x, end_u=minima.y)
323 ],
324 dists = [for (v=min_us) norm(bezier_points(bezier,v)-pt)],
325 min_i = min_index(dists)
326 ) min_us[min_i]
327 ) : let(
328 minima = minima_ranges[0],
329 pp = bezier_points(bezier, minima),
330 err = norm(pp[1]-pp[0])
331 ) err<max_err? mean(minima) :
332 bezier_closest_point(bezier, pt, max_err=max_err, u=minima[0], end_u=minima[1]);
333
334
335// Function: bezier_length()
336// Usage:
337// pathlen = bezier_length(bezier, [start_u], [end_u], [max_deflect]);
338// Topics: Bezier Curves
339// See Also: bezier_points()
340// Description:
341// Approximates the length of the portion of the bezier curve between start_u and end_u.
342// Arguments:
343// bezier = The list of control points that define the Bezier curve.
344// start_u = The Bezier parameter to start measuring measuring from. Between 0 and 1.
345// end_u = The Bezier parameter to end measuring at. Between 0 and 1. Greater than start_u.
346// max_deflect = The largest amount of deflection from the true curve to allow for approximation.
347// Example:
348// bez = [[0,0], [5,35], [60,-25], [80,0]];
349// echo(bezier_length(bez));
350function bezier_length(bezier, start_u=0, end_u=1, max_deflect=0.01) =
351 let(
352 segs = len(bezier) * 2,
353 uvals = lerpn(start_u, end_u, segs+1),
354 path = bezier_points(bezier,uvals),
355 defl = max([
356 for (i=idx(path,e=-3)) let(
357 mp = (path[i] + path[i+2]) / 2
358 ) norm(path[i+1] - mp)
359 ]),
360 mid_u = lerp(start_u, end_u, 0.5)
361 )
362 defl <= max_deflect? path_length(path) :
363 sum([
364 for (i=[0:1:segs-1]) let(
365 su = lerp(start_u, end_u, i/segs),
366 eu = lerp(start_u, end_u, (i+1)/segs)
367 ) bezier_length(bezier, su, eu, max_deflect)
368 ]);
369
370
371
372// Function: bezier_line_intersection()
373// Usage:
374// u = bezier_line_intersection(bezier, line);
375// Topics: Bezier Curves, Geometry, Intersection
376// See Also: bezier_points(), bezier_length(), bezier_closest_point()
377// Description:
378// Finds the intersection points of the 2D Bezier curve with control points `bezier` and the given line, specified as a pair of points.
379// Returns the intersection as a list of `u` values for the Bezier.
380// Arguments:
381// bezier = The list of control points that define a 2D Bezier curve.
382// line = a list of two distinct 2d points defining a line
383function bezier_line_intersection(bezier, line) =
384 assert(is_path(bezier,2), "The input ´bezier´ must be a 2d bezier")
385 assert(_valid_line(line,2), "The input `line` is not a valid 2d line")
386 let(
387 a = _bezier_matrix(len(bezier)-1)*bezier, // bezier algebraic coeffs.
388 n = [-line[1].y+line[0].y, line[1].x-line[0].x], // line normal
389 q = [for(i=[len(a)-1:-1:1]) a[i]*n, (a[0]-line[0])*n] // bezier SDF to line
390 )
391 [for(u=real_roots(q)) if (u>=0 && u<=1) u];
392
393
394
395
396// Section: Bezier Path Functions
397// To contruct more complicated curves you can connect a sequence of Bezier curves end to end.
398// A Bezier path is a flattened list of control points that, along with the degree, represents such a sequence of bezier curves where all of the curves have the same degree.
399// A Bezier path looks like a regular path, since it is just a list of points, but it is not a regular path. Use {{bezpath_curve()}} to convert a Bezier path to a regular path.
400// We interpret a degree N Bezier path as groups of N+1 control points that
401// share endpoints, so they overlap by one point. So if you have an order 3 bezier path `[p0,p1,p2,p3,p4,p5,p6]` then the first
402// Bezier curve control point set is `[p0,p1,p2,p3]` and the second one is `[p3,p4,p5,p6]`. The endpoint, `p3`, is shared between the control point sets.
403// The Bezier degree, which must be known to interpret the Bezier path, defaults to 3.
404
405
406// Function: bezpath_points()
407// Usage:
408// pt = bezpath_points(bezpath, curveind, u, [N]);
409// ptlist = bezpath_points(bezpath, curveind, LIST, [N]);
410// path = bezpath_points(bezpath, curveind, RANGE, [N]);
411// Topics: Bezier Paths
412// See Also: bezier_points(), bezier_curve()
413// Description:
414// Extracts from the Bezier path `bezpath` the control points for the Bezier curve whose index is `curveind` and
415// computes the point or points on the corresponding Bezier curve specified by `u`. If `curveind` is zero you
416// get the first curve. The number of curves is `(len(bezpath)-1)/N` so the maximum index is that number minus one.
417// Arguments:
418// bezpath = A Bezier path path to approximate.
419// curveind = Curve number along the path.
420// u = Parameter values for evaluating the curve, given as a single value, a list or a range.
421// N = The degree of the Bezier path curves. Default: 3
422function bezpath_points(bezpath, curveind, u, N=3) =
423 bezier_points(select(bezpath,curveind*N,(curveind+1)*N), u);
424
425
426// Function: bezpath_curve()
427// Usage:
428// path = bezpath_curve(bezpath, [splinesteps], [N], [endpoint])
429// Topics: Bezier Paths
430// See Also: bezier_points(), bezier_curve()
431// Description:
432// Takes a bezier path and converts it into a path of points.
433// Arguments:
434// bezpath = A bezier path to approximate.
435// splinesteps = Number of straight lines to split each bezier curve into. default=16
436// N = The degree of the bezier curves. Cubic beziers have N=3. Default: 3
437// endpoint = If true, include the very last point of the bezier path. Default: true
438// Example(2D):
439// bez = [
440// [0,0], [-5,30],
441// [20,60], [50,50], [110,30],
442// [60,25], [70,0], [80,-25],
443// [80,-50], [50,-50]
444// ];
445// debug_bezier(bez, N=3, width=2);
446function bezpath_curve(bezpath, splinesteps=16, N=3, endpoint=true) =
447 assert(is_path(bezpath))
448 assert(is_int(N))
449 assert(is_int(splinesteps) && splinesteps>0)
450 assert(len(bezpath)%N == 1, str("A degree ",N," bezier path should have a multiple of ",N," points in it, plus 1."))
451 let(
452 segs = (len(bezpath)-1) / N,
453 step = 1 / splinesteps
454 ) [
455 for (seg = [0:1:segs-1])
456 each bezier_points(select(bezpath, seg*N, (seg+1)*N), [0:step:1-step/2]),
457 if (endpoint) last(bezpath)
458 ];
459
460
461// Function: bezpath_closest_point()
462// Usage:
463// res = bezpath_closest_point(bezpath, pt, [N], [max_err]);
464// Topics: Bezier Paths
465// See Also: bezier_points(), bezier_curve(), bezier_closest_point()
466// Description:
467// Finds an approximation to the closest part of the given bezier path to point `pt`.
468// Returns [segnum, u] for the closest position on the bezier path to the given point `pt`.
469// Arguments:
470// bezpath = A bezier path to approximate.
471// pt = The point to find the closest curve point to.
472// N = The degree of the bezier curves. Cubic beziers have N=3. Default: 3
473// max_err = The maximum allowed error when approximating the closest approach.
474// Example(2D):
475// pt = [100,0];
476// bez = [[0,0], [20,40], [60,-25], [80,0],
477// [100,25], [140,25], [160,0]];
478// pos = bezpath_closest_point(bez, pt);
479// xy = bezpath_points(bez,pos[0],pos[1]);
480// debug_bezier(bez, N=3);
481// color("red") translate(pt) sphere(r=1);
482// color("blue") translate(xy) sphere(r=1);
483function bezpath_closest_point(bezpath, pt, N=3, max_err=0.01, seg=0, min_seg=undef, min_u=undef, min_dist=undef) =
484 assert(is_vector(pt))
485 assert(is_int(N))
486 assert(is_num(max_err))
487 assert(len(bezpath)%N == 1, str("A degree ",N," bezier path shound have a multiple of ",N," points in it, plus 1."))
488 let(curve = select(bezpath,seg*N,(seg+1)*N))
489 (seg*N+1 >= len(bezpath))? (
490 let(curve = select(bezpath, min_seg*N, (min_seg+1)*N))
491 [min_seg, bezier_closest_point(curve, pt, max_err=max_err)]
492 ) : (
493 let(
494 curve = select(bezpath,seg*N,(seg+1)*N),
495 u = bezier_closest_point(curve, pt, max_err=0.05),
496 dist = norm(bezier_points(curve, u)-pt),
497 mseg = (min_dist==undef || dist<min_dist)? seg : min_seg,
498 mdist = (min_dist==undef || dist<min_dist)? dist : min_dist,
499 mu = (min_dist==undef || dist<min_dist)? u : min_u
500 )
501 bezpath_closest_point(bezpath, pt, N, max_err, seg+1, mseg, mu, mdist)
502 );
503
504
505
506// Function: bezpath_length()
507// Usage:
508// plen = bezpath_length(path, [N], [max_deflect]);
509// Topics: Bezier Paths
510// See Also: bezier_points(), bezier_curve(), bezier_length()
511// Description:
512// Approximates the length of the bezier path.
513// Arguments:
514// path = A bezier path to approximate.
515// N = The degree of the bezier curves. Cubic beziers have N=3. Default: 3
516// max_deflect = The largest amount of deflection from the true curve to allow for approximation.
517function bezpath_length(bezpath, N=3, max_deflect=0.001) =
518 assert(is_int(N))
519 assert(is_num(max_deflect))
520 assert(len(bezpath)%N == 1, str("A degree ",N," bezier path shound have a multiple of ",N," points in it, plus 1."))
521 sum([
522 for (seg=[0:1:(len(bezpath)-1)/N-1]) (
523 bezier_length(
524 select(bezpath, seg*N, (seg+1)*N),
525 max_deflect=max_deflect
526 )
527 )
528 ]);
529
530
531
532// Function: path_to_bezpath()
533// Usage:
534// bezpath = path_to_bezpath(path, [closed], [tangents], [uniform], [size=]|[relsize=]);
535// Topics: Bezier Paths, Rounding
536// See Also: path_tangents()
537// Description:
538// Given a 2d or 3d input path and optional list of tangent vectors, computes a cubic (degree 3) bezier
539// path that passes through every point on the input path and matches the tangent vectors. If you do
540// not supply the tangent it will be computed using `path_tangents()`. If the path is closed specify this
541// by setting `closed=true`. The size or relsize parameter determines how far the curve can deviate from
542// the input path. In the case where the curve has a single hump, the size specifies the exact distance
543// between the specified path and the bezier. If you give relsize then it is relative to the segment
544// length (e.g. 0.05 means 5% of the segment length). In 2d when the bezier curve makes an S-curve
545// the size parameter specifies the sum of the deviations of the two peaks of the curve. In 3-space
546// the bezier curve may have three extrema: two maxima and one minimum. In this case the size specifies
547// the sum of the maxima minus the minimum. If you do not supply the tangents then they are computed
548// using `path_tangents()` with `uniform=false` by default. Tangents computed on non-uniform data tend
549// to display overshoots. See `smooth_path()` for examples.
550// Arguments:
551// path = 2D or 3D point list or 1-region that the curve must pass through
552// closed = true if the curve is closed . Default: false
553// tangents = tangents constraining curve direction at each point
554// uniform = set to true to compute tangents with uniform=true. Default: false
555// ---
556// size = absolute size specification for the curve, a number or vector
557// relsize = relative size specification for the curve, a number or vector. Default: 0.1.
558function path_to_bezpath(path, closed, tangents, uniform=false, size, relsize) =
559 is_1region(path) ? path_to_bezpath(path[0], default(closed,true), tangents, uniform, size, relsize) :
560 let(closed=default(closed,false))
561 assert(is_bool(closed))
562 assert(is_bool(uniform))
563 assert(num_defined([size,relsize])<=1, "Can't define both size and relsize")
564 assert(is_path(path,[2,3]),"Input path is not a valid 2d or 3d path")
565 assert(is_undef(tangents) || is_path(tangents,[2,3]),"Tangents must be a 2d or 3d path")
566 assert(is_undef(tangents) || len(path)==len(tangents), "Input tangents must be the same length as the input path")
567 let(
568 curvesize = first_defined([size,relsize,0.1]),
569 relative = is_undef(size),
570 lastpt = len(path) - (closed?0:1)
571 )
572 assert(is_num(curvesize) || len(curvesize)==lastpt, str("Size or relsize must have length ",lastpt))
573 let(
574 sizevect = is_num(curvesize) ? repeat(curvesize, lastpt) : curvesize,
575 tangents = is_def(tangents) ? [for(t=tangents) let(n=norm(t)) assert(!approx(n,0),"Zero tangent vector") t/n] :
576 path_tangents(path, uniform=uniform, closed=closed)
577 )
578 assert(min(sizevect)>0, "Size and relsize must be greater than zero")
579 [
580 for(i=[0:1:lastpt-1])
581 let(
582 first = path[i],
583 second = select(path,i+1),
584 seglength = norm(second-first),
585 dummy = assert(seglength>0, str("Path segment has zero length from index ",i," to ",i+1)),
586 segdir = (second-first)/seglength,
587 tangent1 = tangents[i],
588 tangent2 = -select(tangents,i+1), // Need this to point backwards, in direction of the curve
589 parallel = abs(tangent1*segdir) + abs(tangent2*segdir), // Total component of tangents parallel to the segment
590 Lmax = seglength/parallel, // May be infinity
591 size = relative ? sizevect[i]*seglength : sizevect[i],
592 normal1 = tangent1-(tangent1*segdir)*segdir, // Components of the tangents orthogonal to the segment
593 normal2 = tangent2-(tangent2*segdir)*segdir,
594 p = [ [-3 ,6,-3 ], // polynomial in power form
595 [ 7,-9, 2 ],
596 [-5, 3, 0 ],
597 [ 1, 0, 0 ] ]*[normal1*normal1, normal1*normal2, normal2*normal2],
598 uextreme = approx(norm(p),0) ? []
599 : [for(root = real_roots(p)) if (root>0 && root<1) root],
600 distlist = [for(d=bezier_points([normal1*0, normal1, normal2, normal2*0], uextreme)) norm(d)],
601 scale = len(distlist)==0 ? 0 :
602 len(distlist)==1 ? distlist[0]
603 : sum(distlist) - 2*min(distlist),
604 Ldesired = size/scale, // This will be infinity when the polynomial is zero
605 L = min(Lmax, Ldesired)
606 )
607 each [
608 first,
609 first + L*tangent1,
610 second + L*tangent2
611 ],
612 select(path,lastpt)
613 ];
614
615
616
617
618// Function: bezpath_close_to_axis()
619// Usage:
620// bezpath = bezpath_close_to_axis(bezpath, [axis], [N]);
621// Topics: Bezier Paths
622// See Also: bezpath_offset()
623// Description:
624// Takes a 2D bezier path and closes it to the specified axis.
625// Arguments:
626// bezpath = The 2D bezier path to close to the axis.
627// axis = The axis to close to, "X", or "Y". Default: "X"
628// N = The degree of the bezier curves. Cubic beziers have N=3. Default: 3
629// Example(2D):
630// bez = [[50,30], [40,10], [10,50], [0,30],
631// [-10, 10], [-30,10], [-50,20]];
632// closed = bezpath_close_to_axis(bez);
633// debug_bezier(closed);
634// Example(2D):
635// bez = [[30,50], [10,40], [50,10], [30,0],
636// [10, -10], [10,-30], [20,-50]];
637// closed = bezpath_close_to_axis(bez, axis="Y");
638// debug_bezier(closed);
639function bezpath_close_to_axis(bezpath, axis="X", N=3) =
640 assert(is_path(bezpath,2), "bezpath_close_to_axis() can only work on 2D bezier paths.")
641 assert(is_int(N))
642 assert(len(bezpath)%N == 1, str("A degree ",N," bezier path shound have a multiple of ",N," points in it, plus 1."))
643 let(
644 sp = bezpath[0],
645 ep = last(bezpath)
646 ) (axis=="X")? concat(
647 lerpn([sp.x,0], sp, N, false),
648 list_head(bezpath),
649 lerpn(ep, [ep.x,0], N, false),
650 lerpn([ep.x,0], [sp.x,0], N+1)
651 ) : (axis=="Y")? concat(
652 lerpn([0,sp.y], sp, N, false),
653 list_head(bezpath),
654 lerpn(ep, [0,ep.y], N, false),
655 lerpn([0,ep.y], [0,sp.y], N+1)
656 ) : (
657 assert(in_list(axis, ["X","Y"]))
658 );
659
660
661// Function: bezpath_offset()
662// Usage:
663// bezpath = bezpath_offset(offset, bezier, [N]);
664// Topics: Bezier Paths
665// See Also: bezpath_close_to_axis()
666// Description:
667// Takes a 2D bezier path and closes it with a matching reversed path that is offset by the given `offset` [X,Y] distance.
668// Arguments:
669// offset = Amount to offset second path by.
670// bezier = The 2D bezier path.
671// N = The degree of the bezier curves. Cubic beziers have N=3. Default: 3
672// Example(2D):
673// bez = [[50,30], [40,10], [10,50], [0,30], [-10, 10], [-30,10], [-50,20]];
674// closed = bezpath_offset([0,-5], bez);
675// debug_bezier(closed);
676// Example(2D):
677// bez = [[30,50], [10,40], [50,10], [30,0], [10, -10], [10,-30], [20,-50]];
678// closed = bezpath_offset([-5,0], bez);
679// debug_bezier(closed);
680function bezpath_offset(offset, bezier, N=3) =
681 assert(is_vector(offset,2))
682 assert(is_path(bezier,2), "bezpath_offset() can only work on 2D bezier paths.")
683 assert(is_int(N))
684 assert(len(bezier)%N == 1, str("A degree ",N," bezier path shound have a multiple of ",N," points in it, plus 1."))
685 let(
686 backbez = reverse([ for (pt = bezier) pt+offset ]),
687 bezend = len(bezier)-1
688 ) concat(
689 list_head(bezier),
690 lerpn(bezier[bezend], backbez[0], N, false),
691 list_head(backbez),
692 lerpn(backbez[bezend], bezier[0], N+1)
693 );
694
695
696
697// Section: Cubic Bezier Path Construction
698
699// Function: bez_begin()
700// Topics: Bezier Paths
701// See Also: bez_tang(), bez_joint(), bez_end()
702// Usage:
703// pts = bez_begin(pt, a, r, [p=]);
704// pts = bez_begin(pt, VECTOR, [r], [p=]);
705// Description:
706// This is used to create the first endpoint and control point of a cubic bezier path.
707// Arguments:
708// pt = The starting endpoint for the bezier path.
709// a = If given a scalar, specifies the theta (XY plane) angle in degrees from X+. If given a vector, specifies the direction and possibly distance of the first control point.
710// r = Specifies the distance of the control point from the endpoint `pt`.
711// ---
712// p = If given, specifies the number of degrees away from the Z+ axis.
713// Example(2D): 2D Bezier Path by Angle
714// bezpath = flatten([
715// bez_begin([-50, 0], 45,20),
716// bez_tang ([ 0, 0],-135,20),
717// bez_joint([ 20,-25], 135, 90, 10, 15),
718// bez_end ([ 50, 0], -90,20),
719// ]);
720// debug_bezier(bezpath);
721// Example(2D): 2D Bezier Path by Vector
722// bezpath = flatten([
723// bez_begin([-50,0],[0,-20]),
724// bez_tang ([-10,0],[0,-20]),
725// bez_joint([ 20,-25], [-10,10], [0,15]),
726// bez_end ([ 50,0],[0, 20]),
727// ]);
728// debug_bezier(bezpath);
729// Example(2D): 2D Bezier Path by Vector and Distance
730// bezpath = flatten([
731// bez_begin([-30,0],FWD, 30),
732// bez_tang ([ 0,0],FWD, 30),
733// bez_joint([ 20,-25], 135, 90, 10, 15),
734// bez_end ([ 30,0],BACK,30),
735// ]);
736// debug_bezier(bezpath);
737// Example(3D,FlatSpin,VPD=200): 3D Bezier Path by Angle
738// bezpath = flatten([
739// bez_begin([-30,0,0],90,20,p=135),
740// bez_tang ([ 0,0,0],-90,20,p=135),
741// bez_joint([20,-25,0], 135, 90, 15, 10, p1=135, p2=45),
742// bez_end ([ 30,0,0],-90,20,p=45),
743// ]);
744// debug_bezier(bezpath);
745// Example(3D,FlatSpin,VPD=225): 3D Bezier Path by Vector
746// bezpath = flatten([
747// bez_begin([-30,0,0],[0,-20, 20]),
748// bez_tang ([ 0,0,0],[0,-20,-20]),
749// bez_joint([20,-25,0],[0,10,-10],[0,15,15]),
750// bez_end ([ 30,0,0],[0,-20,-20]),
751// ]);
752// debug_bezier(bezpath);
753// Example(3D,FlatSpin,VPD=225): 3D Bezier Path by Vector and Distance
754// bezpath = flatten([
755// bez_begin([-30,0,0],FWD, 20),
756// bez_tang ([ 0,0,0],DOWN,20),
757// bez_joint([20,-25,0],LEFT,DOWN,r1=20,r2=15),
758// bez_end ([ 30,0,0],DOWN,20),
759// ]);
760// debug_bezier(bezpath);
761function bez_begin(pt,a,r,p) =
762 assert(is_finite(r) || is_vector(a))
763 assert(len(pt)==3 || is_undef(p))
764 is_vector(a)? [pt, pt+(is_undef(r)? a : r*unit(a))] :
765 is_finite(a)? [pt, pt+spherical_to_xyz(r,a,default(p,90))] :
766 assert(false, "Bad arguments.");
767
768
769// Function: bez_tang()
770// Topics: Bezier Paths
771// See Also: bez_begin(), bez_joint(), bez_end()
772// Usage:
773// pts = bez_tang(pt, a, r1, r2, [p=]);
774// pts = bez_tang(pt, VECTOR, [r1], [r2], [p=]);
775// Description:
776// This creates a smooth joint in a cubic bezier path. It creates three points, being the
777// approaching control point, the fixed bezier control point, and the departing control
778// point. The two control points will be collinear with the fixed point, making for a
779// smooth bezier curve at the fixed point. See {{bez_begin()}} for examples.
780// Arguments:
781// pt = The fixed point for the bezier path.
782// a = If given a scalar, specifies the theta (XY plane) angle in degrees from X+. If given a vector, specifies the direction and possibly distance of the departing control point.
783// r1 = Specifies the distance of the approching control point from the fixed point. Overrides the distance component of the vector if `a` contains a vector.
784// r2 = Specifies the distance of the departing control point from the fixed point. Overrides the distance component of the vector if `a` contains a vector. If `r1` is given and `r2` is not, uses the value of `r1` for `r2`.
785// ---
786// p = If given, specifies the number of degrees away from the Z+ axis.
787function bez_tang(pt,a,r1,r2,p) =
788 assert(is_finite(r1) || is_vector(a))
789 assert(len(pt)==3 || is_undef(p))
790 let(
791 r1 = is_num(r1)? r1 : norm(a),
792 r2 = default(r2,r1),
793 p = default(p, 90)
794 )
795 is_vector(a)? [pt-r1*unit(a), pt, pt+r2*unit(a)] :
796 is_finite(a)? [
797 pt-spherical_to_xyz(r1,a,p),
798 pt,
799 pt+spherical_to_xyz(r2,a,p)
800 ] :
801 assert(false, "Bad arguments.");
802
803
804// Function: bez_joint()
805// Topics: Bezier Paths
806// See Also: bez_begin(), bez_tang(), bez_end()
807// Usage:
808// pts = bez_joint(pt, a1, a2, r1, r2, [p1=], [p2=]);
809// pts = bez_joint(pt, VEC1, VEC2, [r1=], [r2=], [p1=], [p2=]);
810// Description:
811// This creates a disjoint corner joint in a cubic bezier path. It creates three points, being
812// the aproaching control point, the fixed bezier control point, and the departing control point.
813// The two control points can be directed in different arbitrary directions from the fixed bezier
814// point. See {{bez_begin()}} for examples.
815// Arguments:
816// pt = The fixed point for the bezier path.
817// a1 = If given a scalar, specifies the theta (XY plane) angle in degrees from X+. If given a vector, specifies the direction and possibly distance of the approaching control point.
818// a2 = If given a scalar, specifies the theta (XY plane) angle in degrees from X+. If given a vector, specifies the direction and possibly distance of the departing control point.
819// r1 = Specifies the distance of the approching control point from the fixed point. Overrides the distance component of the vector if `a1` contains a vector.
820// r2 = Specifies the distance of the departing control point from the fixed point. Overrides the distance component of the vector if `a2` contains a vector.
821// ---
822// p1 = If given, specifies the number of degrees away from the Z+ axis of the approaching control point.
823// p2 = If given, specifies the number of degrees away from the Z+ axis of the departing control point.
824function bez_joint(pt,a1,a2,r1,r2,p1,p2) =
825 assert(is_finite(r1) || is_vector(a1))
826 assert(is_finite(r2) || is_vector(a2))
827 assert(len(pt)==3 || (is_undef(p1) && is_undef(p2)))
828 let(
829 r1 = is_num(r1)? r1 : norm(a1),
830 r2 = is_num(r2)? r2 : norm(a2),
831 p1 = default(p1, 90),
832 p2 = default(p2, 90)
833 ) [
834 if (is_vector(a1)) (pt+r1*unit(a1))
835 else if (is_finite(a1)) (pt+spherical_to_xyz(r1,a1,p1))
836 else assert(false, "Bad Arguments"),
837 pt,
838 if (is_vector(a2)) (pt+r2*unit(a2))
839 else if (is_finite(a2)) (pt+spherical_to_xyz(r2,a2,p2))
840 else assert(false, "Bad Arguments")
841 ];
842
843
844// Function: bez_end()
845// Topics: Bezier Paths
846// See Also: bez_tang(), bez_joint(), bez_end()
847// Usage:
848// pts = bez_end(pt, a, r, [p=]);
849// pts = bez_end(pt, VECTOR, [r], [p=]);
850// Description:
851// This is used to create the approaching control point, and the endpoint of a cubic bezier path.
852// See {{bez_begin()}} for examples.
853// Arguments:
854// pt = The starting endpoint for the bezier path.
855// a = If given a scalar, specifies the theta (XY plane) angle in degrees from X+. If given a vector, specifies the direction and possibly distance of the first control point.
856// r = Specifies the distance of the control point from the endpoint `pt`.
857// p = If given, specifies the number of degrees away from the Z+ axis.
858function bez_end(pt,a,r,p) =
859 assert(is_finite(r) || is_vector(a))
860 assert(len(pt)==3 || is_undef(p))
861 is_vector(a)? [pt+(is_undef(r)? a : r*unit(a)), pt] :
862 is_finite(a)? [pt+spherical_to_xyz(r,a,default(p,90)), pt] :
863 assert(false, "Bad arguments.");
864
865
866
867// Section: Bezier Surfaces
868
869
870// Function: is_bezier_patch()
871// Usage:
872// bool = is_bezier_patch(x);
873// Topics: Bezier Patches, Type Checking
874// Description:
875// Returns true if the given item is a bezier patch.
876// Arguments:
877// x = The value to check the type of.
878function is_bezier_patch(x) =
879 is_list(x) && is_list(x[0]) && is_vector(x[0][0]) && len(x[0]) == len(x[len(x)-1]);
880
881
882// Function: bezier_patch_flat()
883// Usage:
884// patch = bezier_patch_flat(size, [N=], [spin=], [orient=], [trans=]);
885// Topics: Bezier Patches
886// See Also: bezier_patch_points()
887// Description:
888// Returns a flat rectangular bezier patch of degree `N`, centered on the XY plane.
889// Arguments:
890// size = 2D XY size of the patch.
891// ---
892// N = Degree of the patch to generate. Since this is flat, a degree of 1 should usually be sufficient.
893// orient = The orientation to rotate the edge patch into. Given as an [X,Y,Z] rotation angle list.
894// trans = Amount to translate patch, after rotating to `orient`.
895// Example(3D):
896// patch = bezier_patch_flat(size=[100,100], N=3);
897// debug_bezier_patches([patch], size=1, showcps=true);
898function bezier_patch_flat(size=[100,100], N=4, spin=0, orient=UP, trans=[0,0,0]) =
899 let(
900 patch = [
901 for (x=[0:1:N]) [
902 for (y=[0:1:N])
903 v_mul(point3d(size), [x/N-0.5, 0.5-y/N, 0])
904 ]
905 ],
906 m = move(trans) * rot(a=spin, from=UP, to=orient)
907 ) [for (row=patch) apply(m, row)];
908
909
910
911// Function: bezier_patch_reverse()
912// Usage:
913// rpatch = bezier_patch_reverse(patch);
914// Topics: Bezier Patches
915// See Also: bezier_patch_points(), bezier_patch_flat()
916// Description:
917// Reverses the patch, so that the faces generated from it are flipped back to front.
918// Arguments:
919// patch = The patch to reverse.
920function bezier_patch_reverse(patch) =
921 [for (row=patch) reverse(row)];
922
923
924// Function: bezier_patch_points()
925// Usage:
926// pt = bezier_patch_points(patch, u, v);
927// ptgrid = bezier_patch_points(patch, LIST, LIST);
928// ptgrid = bezier_patch_points(patch, RANGE, RANGE);
929// Topics: Bezier Patches
930// See Also: bezier_patch_normals(), bezier_points(), bezier_curve(), bezpath_curve()
931// Description:
932// Sample a bezier patch on a listed point set. The bezier patch must be a rectangular array of
933// points, and it will be sampled at all the (u,v) pairs that you specify. If you give u and v
934// as single numbers you'll get a single point back. If you give u and v as lists or ranges you'll
935// get a 2d rectangular array of points. If one but not both of u and v is a list or range then you'll
936// get a list of points.
937// Arguments:
938// patch = The 2D array of control points for a Bezier patch.
939// u = The bezier u parameter (inner list of patch). Generally between 0 and 1. Can be a list, range or value.
940// v = The bezier v parameter (outer list of patch). Generally between 0 and 1. Can be a list, range or value.
941// Example(3D):
942// patch = [
943// [[-50, 50, 0], [-16, 50, 20], [ 16, 50, 20], [50, 50, 0]],
944// [[-50, 16, 20], [-16, 16, 40], [ 16, 16, 40], [50, 16, 20]],
945// [[-50,-16, 20], [-16,-16, 40], [ 16,-16, 40], [50,-16, 20]],
946// [[-50,-50, 0], [-16,-50, 20], [ 16,-50, 20], [50,-50, 0]]
947// ];
948// debug_bezier_patches(patches=[patch], size=1, showcps=true);
949// pt = bezier_patch_points(patch, 0.6, 0.75);
950// translate(pt) color("magenta") sphere(d=3, $fn=12);
951// Example(3D): Getting Multiple Points at Once
952// patch = [
953// [[-50, 50, 0], [-16, 50, 20], [ 16, 50, 20], [50, 50, 0]],
954// [[-50, 16, 20], [-16, 16, 40], [ 16, 16, 40], [50, 16, 20]],
955// [[-50,-16, 20], [-16,-16, 40], [ 16,-16, 40], [50,-16, 20]],
956// [[-50,-50, 0], [-16,-50, 20], [ 16,-50, 20], [50,-50, 0]]
957// ];
958// debug_bezier_patches(patches=[patch], size=1, showcps=true);
959// pts = bezier_patch_points(patch, [0:0.2:1], [0:0.2:1]);
960// for (row=pts) move_copies(row) color("magenta") sphere(d=3, $fn=12);
961function bezier_patch_points(patch, u, v) =
962 assert(is_range(u) || is_vector(u) || is_finite(u), "Input u is invalid")
963 assert(is_range(v) || is_vector(v) || is_finite(v), "Input v is invalid")
964 !is_num(u) && !is_num(v) ?
965 let(
966 vbezes = [for (i = idx(patch[0])) bezier_points(column(patch,i), u)]
967 )
968 [for (i = idx(vbezes[0])) bezier_points(column(vbezes,i), v)]
969 : is_num(u) && is_num(v)? bezier_points([for (bez = patch) bezier_points(bez, v)], u)
970 : is_num(u) ? bezier_patch_points(patch,force_list(u),v)[0]
971 : column(bezier_patch_points(patch,u,force_list(v)),0);
972
973
974
975
976function _bezier_rectangle(patch, splinesteps=16, style="default") =
977 let(
978 uvals = lerpn(0,1,splinesteps.x+1),
979 vvals = lerpn(1,0,splinesteps.y+1),
980 pts = bezier_patch_points(patch, uvals, vvals)
981 )
982 vnf_vertex_array(pts, style=style, reverse=false);
983
984
985// Function: bezier_vnf()
986// Usage:
987// vnf = bezier_vnf(patches, [splinesteps], [style]);
988// Topics: Bezier Patches
989// See Also: bezier_patch_points(), bezier_patch_flat()
990// Description:
991// Convert a patch or list of patches into the corresponding Bezier surface, representing the
992// result as a [VNF structure](vnf.scad). The `splinesteps` argument specifies the sampling grid of
993// the surface for each patch by specifying the number of segments on the borders of the surface.
994// It can be a scalar, which gives a uniform grid, or
995// it can be [USTEPS, VSTEPS], which gives difference spacing in the U and V parameters.
996// Note that the surface you produce may be disconnected and is not necessarily a valid manifold in OpenSCAD.
997// Arguments:
998// patches = The bezier patch or list of bezier patches to convert into a vnf.
999// splinesteps = Number of segments on the border of the bezier surface. You can specify [USTEPS,VSTEPS]. Default: 16
1000// style = The style of subdividing the quads into faces. Valid options are "default", "alt", "min_edge", "quincunx", "convex" and "concave". See {{vnf_vertex_array()}}. Default: "default"
1001// Example(3D):
1002// patch = [
1003// // u=0,v=0 u=1,v=0
1004// [[-50,-50, 0], [-16,-50, 20], [ 16,-50, -20], [50,-50, 0]],
1005// [[-50,-16, 20], [-16,-16, 20], [ 16,-16, -20], [50,-16, 20]],
1006// [[-50, 16, 20], [-16, 16, -20], [ 16, 16, 20], [50, 16, 20]],
1007// [[-50, 50, 0], [-16, 50, -20], [ 16, 50, 20], [50, 50, 0]],
1008// // u=0,v=1 u=1,v=1
1009// ];
1010// vnf = bezier_vnf(patch, splinesteps=16);
1011// vnf_polyhedron(vnf);
1012// Example(3D,FlatSpin,VPD=444): Combining multiple patches
1013// patch = [
1014// // u=0,v=0 u=1,v=0
1015// [[0, 0,0], [33, 0, 0], [67, 0, 0], [100, 0,0]],
1016// [[0, 33,0], [33, 33, 33], [67, 33, 33], [100, 33,0]],
1017// [[0, 67,0], [33, 67, 33], [67, 67, 33], [100, 67,0]],
1018// [[0,100,0], [33,100, 0], [67,100, 0], [100,100,0]],
1019// // u=0,v=1 u=1,v=1
1020// ];
1021// tpatch = translate([-50,-50,50], patch);
1022// vnf = bezier_vnf([
1023// tpatch,
1024// xrot(90, tpatch),
1025// xrot(-90, tpatch),
1026// xrot(180, tpatch),
1027// yrot(90, tpatch),
1028// yrot(-90, tpatch)]);
1029// vnf_polyhedron(vnf);
1030// Example(3D):
1031// patch1 = [
1032// [[18,18,0], [33, 0, 0], [ 67, 0, 0], [ 82, 18,0]],
1033// [[ 0,40,0], [ 0, 0,100], [100, 0, 20], [100, 40,0]],
1034// [[ 0,60,0], [ 0,100,100], [100,100, 20], [100, 60,0]],
1035// [[18,82,0], [33,100, 0], [ 67,100, 0], [ 82, 82,0]],
1036// ];
1037// patch2 = [
1038// [[18,82,0], [33,100, 0], [ 67,100, 0], [ 82, 82,0]],
1039// [[ 0,60,0], [ 0,100,-50], [100,100,-50], [100, 60,0]],
1040// [[ 0,40,0], [ 0, 0,-50], [100, 0,-50], [100, 40,0]],
1041// [[18,18,0], [33, 0, 0], [ 67, 0, 0], [ 82, 18,0]],
1042// ];
1043// vnf = bezier_vnf(patches=[patch1, patch2], splinesteps=16);
1044// vnf_polyhedron(vnf);
1045// Example(3D): Connecting Patches with asymmetric splinesteps. Note it is fastest to join all the VNFs at once, which happens in vnf_polyhedron, rather than generating intermediate joined partial surfaces.
1046// steps = 8;
1047// edge_patch = [
1048// // u=0, v=0 u=1,v=0
1049// [[-60, 0,-40], [0, 0,-40], [60, 0,-40]],
1050// [[-60, 0, 0], [0, 0, 0], [60, 0, 0]],
1051// [[-60,40, 0], [0,40, 0], [60,40, 0]],
1052// // u=0, v=1 u=1,v=1
1053// ];
1054// corner_patch = [
1055// // u=0, v=0 u=1,v=0
1056// [[ 0, 40,-40], [ 0, 0,-40], [40, 0,-40]],
1057// [[ 0, 40, 0], [ 0, 0, 0], [40, 0, 0]],
1058// [[40, 40, 0], [40, 40, 0], [40, 40, 0]],
1059// // u=0, v=1 u=1,v=1
1060// ];
1061// face_patch = bezier_patch_flat([120,120],orient=LEFT);
1062// edges = [
1063// for (axrot=[[0,0,0],[0,90,0],[0,0,90]], xang=[-90:90:180])
1064// bezier_vnf(
1065// splinesteps=[steps,1],
1066// rot(a=axrot,
1067// p=rot(a=[xang,0,0],
1068// p=translate(v=[0,-100,100],p=edge_patch)
1069// )
1070// )
1071// )
1072// ];
1073// corners = [
1074// for (xang=[0,180], zang=[-90:90:180])
1075// bezier_vnf(
1076// splinesteps=steps,
1077// rot(a=[xang,0,zang],
1078// p=translate(v=[-100,-100,100],p=corner_patch)
1079// )
1080// )
1081// ];
1082// faces = [
1083// for (axrot=[[0,0,0],[0,90,0],[0,0,90]], zang=[0,180])
1084// bezier_vnf(
1085// splinesteps=1,
1086// rot(a=axrot,
1087// p=zrot(zang,move([-100,0,0], face_patch))
1088// )
1089// )
1090// ];
1091// vnf_polyhedron(concat(edges,corners,faces));
1092function bezier_vnf(patches=[], splinesteps=16, style="default") =
1093 assert(is_num(splinesteps) || is_vector(splinesteps,2))
1094 assert(all_positive(splinesteps))
1095 let(splinesteps = force_list(splinesteps,2))
1096 is_bezier_patch(patches)? _bezier_rectangle(patches, splinesteps=splinesteps,style=style)
1097 : assert(is_list(patches),"Invalid patch list")
1098 vnf_join(
1099 [
1100 for (patch=patches)
1101 is_bezier_patch(patch)? _bezier_rectangle(patch, splinesteps=splinesteps,style=style)
1102 : assert(false,"Invalid patch list")
1103 ]
1104 );
1105
1106
1107
1108// Function: bezier_vnf_degenerate_patch()
1109// Usage:
1110// vnf = bezier_vnf_degenerate_patch(patch, [splinesteps], [reverse]);
1111// vnf_edges = bezier_vnf_degenerate_patch(patch, [splinesteps], [reverse], return_edges=true);
1112// Description:
1113// Returns a VNF for a degenerate rectangular bezier patch where some of the corners of the patch are
1114// equal. If the resulting patch has no faces then returns an empty VNF. Note that due to the degeneracy,
1115// the shape of the surface can be triangular even though the underlying patch is a rectangle.
1116// If you specify return_edges then the return is a list whose first element is the vnf and whose second
1117// element lists the edges in the order [left, right, top, bottom], where each list is a list of the actual
1118// point values, but possibly only a single point if that edge is degenerate.
1119// The method checks for various types of degeneracy and uses a triangular or partly triangular array of sample points.
1120// See examples below for the types of degeneracy detected and how the patch is sampled for those cases.
1121// Note that splinesteps is the same for both directions of the patch, so it cannot be an array.
1122// Arguments:
1123// patch = Patch to process
1124// splinesteps = Number of segments to produce on each side. Default: 16
1125// reverse = reverse direction of faces. Default: false
1126// return_edges = if true return the points on the four edges: [left, right, top, bottom]. Default: false
1127// Example(3D,NoAxes): This quartic patch is degenerate at one corner, where a row of control points are equal. Processing this degenerate patch normally produces excess triangles near the degenerate point.
1128// splinesteps=8;
1129// patch=[
1130// repeat([-12.5, 12.5, 15],5),
1131// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]],
1132// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]],
1133// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]],
1134// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]]
1135// ];
1136// vnf_wireframe((bezier_vnf(patch, splinesteps)),width=0.1);
1137// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1138// Example(3D,NoAxes): With bezier_vnf_degenerate_patch the degenerate point does not have excess triangles. The top half of the patch decreases the number of sampled points by 2 for each row.
1139// splinesteps=8;
1140// patch=[
1141// repeat([-12.5, 12.5, 15],5),
1142// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]],
1143// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]],
1144// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]],
1145// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]]
1146// ];
1147// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1148// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1149// Example(3D,NoAxes): With splinesteps odd you get one "odd" row where the point count decreases by 1 instead of 2. You may prefer even values for splinesteps to avoid this.
1150// splinesteps=7;
1151// patch=[
1152// repeat([-12.5, 12.5, 15],5),
1153// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]],
1154// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]],
1155// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]],
1156// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]]
1157// ];
1158// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1159// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1160// Example(3D,NoAxes): A more extreme degeneracy occurs when the top half of a patch is degenerate to a line. (For odd length patches the middle row must be degenerate to trigger this style.) In this case the number of points in each row decreases by 1 for every row. It doesn't matter of splinesteps is odd or even.
1161// splinesteps=8;
1162// patch = [[[10, 0, 0], [10, -10.4, 0], [10, -20.8, 0], [1.876, -14.30, 0], [-6.24, -7.8, 0]],
1163// [[5, 0, 0], [5, -5.2, 0], [5, -10.4, 0], [0.938, -7.15, 0], [-3.12, -3.9, 0]],
1164// repeat([0,0,0],5),
1165// repeat([0,0,5],5),
1166// repeat([0,0,10],5)
1167// ];
1168// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1169// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1170// Example(3D,NoScales): Here is a degenerate cubic patch.
1171// splinesteps=8;
1172// patch = [ [ [-20,0,0], [-10,0,0],[0,10,0],[0,20,0] ],
1173// [ [-20,0,10], [-10,0,10],[0,10,10],[0,20,10]],
1174// [ [-10,0,20], [-5,0,20], [0,5,20], [0,10,20]],
1175// repeat([0,0,30],4)
1176// ];
1177// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1178// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1179// Example(3D,NoScales): A more extreme degenerate cubic patch, where two rows are equal.
1180// splinesteps=8;
1181// patch = [ [ [-20,0,0], [-10,0,0],[0,10,0],[0,20,0] ],
1182// [ [-20,0,10], [-10,0,10],[0,10,10],[0,20,10] ],
1183// repeat([-10,10,20],4),
1184// repeat([-10,10,30],4)
1185// ];
1186// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1187// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1188// Example(3D,NoScales): Quadratic patch degenerate at the right side:
1189// splinesteps=8;
1190// patch = [[[0, -10, 0],[10, -5, 0],[20, 0, 0]],
1191// [[0, 0, 0], [10, 0, 0], [20, 0, 0]],
1192// [[0, 0, 10], [10, 0, 5], [20, 0, 0]]];
1193// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1194// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1195// Example(3D,NoAxes): Cubic patch degenerate at both ends. In this case the point count changes by 2 at every row.
1196// splinesteps=8;
1197// patch = [
1198// repeat([10,-10,0],4),
1199// [ [-20,0,0], [-1,0,0],[0,10,0],[0,20,0] ],
1200// [ [-20,0,10], [-10,0,10],[0,10,10],[0,20,10] ],
1201// repeat([-10,10,20],4),
1202// ];
1203// vnf_wireframe(bezier_vnf_degenerate_patch(patch, splinesteps),width=0.1);
1204// color("red")move_copies(flatten(patch)) sphere(r=0.3,$fn=9);
1205function bezier_vnf_degenerate_patch(patch, splinesteps=16, reverse=false, return_edges=false) =
1206 !return_edges ? bezier_vnf_degenerate_patch(patch, splinesteps, reverse, true)[0] :
1207 assert(is_bezier_patch(patch), "Input is not a Bezier patch")
1208 assert(is_int(splinesteps) && splinesteps>0, "splinesteps must be a positive integer")
1209 let(
1210 row_degen = [for(row=patch) all_equal(row)],
1211 col_degen = [for(col=transpose(patch)) all_equal(col)],
1212 top_degen = row_degen[0],
1213 bot_degen = last(row_degen),
1214 left_degen = col_degen[0],
1215 right_degen = last(col_degen),
1216 samplepts = lerpn(0,1,splinesteps+1)
1217 )
1218 all(row_degen) && all(col_degen) ? // fully degenerate case
1219 [EMPTY_VNF, repeat([patch[0][0]],4)] :
1220 all(row_degen) ? // degenerate to a line (top to bottom)
1221 let(pts = bezier_points(column(patch,0), samplepts))
1222 [EMPTY_VNF, [pts,pts,[pts[0]],[last(pts)]]] :
1223 all(col_degen) ? // degenerate to a line (left to right)
1224 let(pts = bezier_points(patch[0], samplepts))
1225 [EMPTY_VNF, [[pts[0]], [last(pts)], pts, pts]] :
1226 !top_degen && !bot_degen && !left_degen && !right_degen ? // non-degenerate case
1227 let(pts = bezier_patch_points(patch, samplepts, samplepts))
1228 [
1229 vnf_vertex_array(pts, reverse=!reverse),
1230 [column(pts,0), column(pts,len(pts)-1), pts[0], last(pts)]
1231 ] :
1232 top_degen && bot_degen ?
1233 let(
1234 rowcount = [
1235 each list([3:2:splinesteps]),
1236 if (splinesteps%2==0) splinesteps+1,
1237 each reverse(list([3:2:splinesteps]))
1238 ],
1239 bpatch = [for(i=[0:1:len(patch[0])-1]) bezier_points(column(patch,i), samplepts)],
1240 pts = [
1241 [bpatch[0][0]],
1242 for(j=[0:splinesteps-2]) bezier_points(column(bpatch,j+1), lerpn(0,1,rowcount[j])),
1243 [last(bpatch[0])]
1244 ],
1245 vnf = vnf_tri_array(pts, reverse=!reverse)
1246 ) [
1247 vnf,
1248 [
1249 column(pts,0),
1250 [for(row=pts) last(row)],
1251 pts[0],
1252 last(pts),
1253 ]
1254 ] :
1255 bot_degen ? // only bottom is degenerate
1256 let(
1257 result = bezier_vnf_degenerate_patch(reverse(patch), splinesteps=splinesteps, reverse=!reverse, return_edges=true)
1258 )
1259 [
1260 result[0],
1261 [reverse(result[1][0]), reverse(result[1][1]), (result[1][3]), (result[1][2])]
1262 ] :
1263 top_degen ? // only top is degenerate
1264 let(
1265 full_degen = len(patch)>=4 && all(select(row_degen,1,ceil(len(patch)/2-1))),
1266 rowmax = full_degen ? count(splinesteps+1) :
1267 [for(j=[0:splinesteps]) j<=splinesteps/2 ? 2*j : splinesteps],
1268 bpatch = [for(i=[0:1:len(patch[0])-1]) bezier_points(column(patch,i), samplepts)],
1269 pts = [
1270 [bpatch[0][0]],
1271 for(j=[1:splinesteps]) bezier_points(column(bpatch,j), lerpn(0,1,rowmax[j]+1))
1272 ],
1273 vnf = vnf_tri_array(pts, reverse=!reverse)
1274 ) [
1275 vnf,
1276 [
1277 column(pts,0),
1278 [for(row=pts) last(row)],
1279 pts[0],
1280 last(pts),
1281 ]
1282 ] :
1283 // must have left or right degeneracy, so transpose and recurse
1284 let(
1285 result = bezier_vnf_degenerate_patch(transpose(patch), splinesteps=splinesteps, reverse=!reverse, return_edges=true)
1286 )
1287 [result[0],
1288 select(result[1],[2,3,0,1])
1289 ];
1290
1291
1292// Function: bezier_patch_normals()
1293// Usage:
1294// n = bezier_patch_normals(patch, u, v);
1295// ngrid = bezier_patch_normals(patch, LIST, LIST);
1296// ngrid = bezier_patch_normals(patch, RANGE, RANGE);
1297// Topics: Bezier Patches
1298// See Also: bezier_patch_points(), bezier_points(), bezier_curve(), bezpath_curve()
1299// Description:
1300// Compute the unit normal vector to a bezier patch at the listed point set. The bezier patch must be a rectangular array of
1301// points, and the normal will be computed at all the (u,v) pairs that you specify. If you give u and v
1302// as single numbers you'll get a single point back. If you give u and v as lists or ranges you'll
1303// get a 2d rectangular array of points. If one but not both of u and v is a list or range then you'll
1304// get a list of points.
1305// .
1306// This function works by computing the cross product of the tangents. In some degenerate cases the one of the tangents
1307// can be zero, so the normal vector does not exist. In this case, undef is returned. Another degenerate case
1308// occurs when the tangents are parallel, or nearly parallel. In this case you will get a unit vector returned but it will not
1309// be the correct normal vector. This can happen if you use a degenerate patch, or if you give two of the edges of your patch a smooth "corner"
1310// so that the u and v directions are parallel at the corner.
1311// Arguments:
1312// patch = The 2D array of control points for a Bezier patch.
1313// u = The bezier u parameter (inner list of patch). Generally between 0 and 1. Can be a list, range or value.
1314// v = The bezier v parameter (outer list of patch). Generally between 0 and 1. Can be a list, range or value.
1315// Example(3D,Med,VPR=[71.1,0,155.9],VPD=292.705,VPT=[20.4724,38.7273,22.7683],NoAxes): Normal vectors on a patch
1316// patch = [
1317// // u=0,v=0 u=1,v=0
1318// [[-50,-50, 0], [-16,-50, 20], [ 16,-50, -20], [50,-50, 0]],
1319// [[-50,-16, 40], [-16,-16, 20], [ 16,-16, -20], [50,-16, 70]],
1320// [[-50, 16, 20], [-16, 16, -20], [ 16, 37, 20], [70, 16, 20]],
1321// [[-50, 50, 0], [73, 50, -40], [ 16, 50, 20], [50, 50, 0]],
1322// // u=0,v=1 u=1,v=1
1323// ];
1324// vnf_polyhedron(bezier_vnf(patch,splinesteps=30));
1325// uv = lerpn(0,1,12);
1326// pts = bezier_patch_points(patch, uv, uv);
1327// normals = bezier_patch_normals(patch, uv, uv);
1328// for(i=idx(uv),j=idx(uv)){
1329// stroke([pts[i][j],pts[i][j]-6*normals[i][j]], width=0.5,
1330// endcap1="dot",endcap2="arrow2",color="blue");
1331// }
1332// Example(3D,NoAxes,Med,VPR=[72.5,0,288.9],VPD=192.044,VPT=[51.6089,48.118,5.89088]): This example gives invalid normal vectors at the four corners of the patch where the u and v directions are parallel. You can see how the line of triangulation is approaching parallel at the edge, and the invalid vectors in red point in a completely incorrect direction.
1333// patch = [
1334// [[18,18,0], [33, 0, 0], [ 67, 0, 0], [ 82, 18,0]],
1335// [[ 0,40,0], [ 0, 0,100], [100, 0, 20], [100, 40,0]],
1336// [[ 0,60,0], [ 0,100,100], [100,100, 20], [100, 60,0]],
1337// [[18,82,0], [33,100, 0], [ 67,100, 0], [ 82, 82,0]],
1338// ];
1339// vnf_polyhedron(bezier_vnf(patch,splinesteps=30));
1340// uv = lerpn(0,1,7);
1341// pts = bezier_patch_points(patch, uv, uv);
1342// normals = bezier_patch_normals(patch, uv, uv);
1343// for(i=idx(uv),j=idx(uv)){
1344// color=((uv[i]==0 || uv[i]==1) && (uv[j]==0 || uv[j]==1))
1345// ? "red" : "blue";
1346// stroke([pts[i][j],pts[i][j]-8*normals[i][j]], width=0.5,
1347// endcap1="dot",endcap2="arrow2",color=color);
1348// }
1349// Example(3D,Med,NoAxes,VPR=[56.4,0,71.9],VPD=66.9616,VPT=[10.2954,1.33721,19.4484]): This degenerate patch has normals everywhere, but computational of the normal fails at the point of degeneracy, the top corner.
1350// patch=[
1351// repeat([-12.5, 12.5, 15],5),
1352// [[-6.25, 11.25, 15], [-6.25, 8.75, 15], [-6.25, 6.25, 15], [-8.75, 6.25, 15], [-11.25, 6.25, 15]],
1353// [[0, 10, 15], [0, 5, 15], [0, 0, 15], [-5, 0, 15], [-10, 0, 15]],
1354// [[0, 10, 8.75], [0, 5, 8.75], [0, 0, 8.75], [-5, 0, 8.75], [-10, 0, 8.75]],
1355// [[0, 10, 2.5], [0, 5, 2.5], [0, 0, 2.5], [-5, 0, 2.5], [-10, 0, 2.5]]
1356// ];
1357// vnf_polyhedron(bezier_vnf(patch, 32));
1358// uv = lerpn(0,1,8);
1359// pts = bezier_patch_points(patch, uv, uv);
1360// normals = bezier_patch_normals(patch, uv, uv);
1361// for(i=idx(uv),j=idx(uv)){
1362// if (is_def(normals[i][j]))
1363// stroke([pts[i][j],pts[i][j]-2*normals[i][j]], width=0.1,
1364// endcap1="dot",endcap2="arrow2",color="blue");
1365// }
1366// Example(3D,Med,NoAxes,VPR=[48,0,23.6],VPD=32.0275,VPT=[-0.145727,-0.0532125,1.74224]): This example has a singularities where the tangent lines don't exist, so the normal will be undef at those points.
1367// pts1 = [ [-5,0,0], [5,0,5], [-5,0,5], [5,0,0] ];
1368// pts2 = [ [0,-5,0], [0,5,5], [0,-5,5], [0,5,0] ];
1369// patch = [for(i=[0:3])
1370// [for(j=[0:3]) pts1[i]+pts2[j] ] ];
1371// vnf_polyhedron(bezier_vnf(patch, 163));
1372// uv = [0,.1,.2,.3,.7,.8,.9,1];//lerpn(0,1,8);
1373// pts = bezier_patch_points(patch, uv, uv);
1374// normals = bezier_patch_normals(patch, uv, uv);
1375// for(i=idx(uv),j=idx(uv))
1376// stroke([pts[i][j],pts[i][j]+2*normals[i][j]], width=0.08,
1377// endcap1="dot",endcap2="arrow2",color="blue");
1378
1379function bezier_patch_normals(patch, u, v) =
1380 assert(is_range(u) || is_vector(u) || is_finite(u), "Input u is invalid")
1381 assert(is_range(v) || is_vector(v) || is_finite(v), "Input v is invalid")
1382 !is_num(u) && !is_num(v) ?
1383 let(
1384 vbezes = [for (i = idx(patch[0])) bezier_points(column(patch,i), u)],
1385 dvbezes = [for (i = idx(patch[0])) bezier_derivative(column(patch,i), u)],
1386 v_tangent = [for (i = idx(vbezes[0])) bezier_derivative(column(vbezes,i), v)],
1387 u_tangent = [for (i = idx(vbezes[0])) bezier_points(column(dvbezes,i), v)]
1388 )
1389 [for(i=idx(u_tangent)) [for(j=idx(u_tangent[0])) unit(cross(u_tangent[i][j],v_tangent[i][j]),undef)]]
1390 : is_num(u) && is_num(v)?
1391 let(
1392 du = bezier_derivative([for (bez = patch) bezier_points(bez, v)], u),
1393 dv = bezier_points([for (bez = patch) bezier_derivative(bez, v)], u)
1394 )
1395 unit(cross(du,dv),undef)
1396 : is_num(u) ? bezier_patch_normals(patch,force_list(u),v)[0]
1397 : column(bezier_patch_normals(patch,u,force_list(v)),0);
1398
1399
1400
1401// Section: Debugging Beziers
1402
1403
1404// Module: debug_bezier()
1405// Usage:
1406// debug_bezier(bez, [size], [N=]);
1407// Topics: Bezier Paths, Debugging
1408// See Also: bezpath_curve()
1409// Description:
1410// Renders 2D or 3D bezier paths and their associated control points.
1411// Useful for debugging bezier paths.
1412// Arguments:
1413// bez = the array of points in the bezier.
1414// size = diameter of the lines drawn.
1415// ---
1416// N = Mark the first and every Nth vertex after in a different color and shape.
1417// Example(2D):
1418// bez = [
1419// [-10, 0], [-15, -5],
1420// [ -5, -10], [ 0, -10], [ 5, -10],
1421// [ 14, -5], [ 15, 0], [16, 5],
1422// [ 5, 10], [ 0, 10]
1423// ];
1424// debug_bezier(bez, N=3, width=0.5);
1425module debug_bezier(bezpath, width=1, N=3) {
1426 no_children($children);
1427 check =
1428 assert(is_path(bezpath),"bezpath must be a path")
1429 assert(is_int(N) && N>0, "N must be a positive integer")
1430 assert(len(bezpath)%N == 1, str("A degree ",N," bezier path shound have a multiple of ",N," points in it, plus 1."));
1431 $fn=8;
1432 stroke(bezpath_curve(bezpath, N=N), width=width, color="cyan");
1433 color("green")
1434 if (N!=3)
1435 stroke(bezpath, width=width);
1436 else
1437 for(i=[1:3:len(bezpath)]) stroke(select(bezpath,max(0,i-2), min(len(bezpath)-1,i)), width=width);
1438 twodim = len(bezpath[0])==2;
1439 color("red") move_copies(bezpath)
1440 if ($idx % N !=0)
1441 if (twodim){
1442 rect([width/2, width*3]);
1443 rect([width*3, width/2]);
1444 } else {
1445 zcyl(d=width/2, h=width*3);
1446 xcyl(d=width/2, h=width*3);
1447 ycyl(d=width/2, h=width*3);
1448 }
1449 color("blue") move_copies(bezpath)
1450 if ($idx % N ==0)
1451 if (twodim) circle(d=width*2.25); else sphere(d=width*2.25);
1452 if (twodim) color("red") move_copies(bezpath)
1453 if ($idx % N !=0) circle(d=width/2);
1454}
1455
1456
1457// Module: debug_bezier_patches()
1458// Usage:
1459// debug_bezier_patches(patches, [size=], [splinesteps=], [showcps=], [showdots=], [showpatch=], [convexity=], [style=]);
1460// Topics: Bezier Patches, Debugging
1461// See Also: bezier_patch_points(), bezier_patch_flat(), bezier_vnf()
1462// Description:
1463// Shows the surface, and optionally, control points of a list of bezier patches.
1464// Arguments:
1465// patches = A list of rectangular bezier patches.
1466// ---
1467// splinesteps = Number of segments to divide each bezier curve into. Default: 16
1468// showcps = If true, show the controlpoints as well as the surface. Default: true.
1469// showdots = If true, shows the calculated surface vertices. Default: false.
1470// showpatch = If true, shows the surface faces. Default: true.
1471// size = Size to show control points and lines.
1472// style = The style of subdividing the quads into faces. Valid options are "default", "alt", and "quincunx".
1473// convexity = Max number of times a line could intersect a wall of the shape.
1474// Example:
1475// patch1 = [
1476// [[15,15,0], [33, 0, 0], [ 67, 0, 0], [ 85, 15,0]],
1477// [[ 0,33,0], [33, 33, 50], [ 67, 33, 50], [100, 33,0]],
1478// [[ 0,67,0], [33, 67, 50], [ 67, 67, 50], [100, 67,0]],
1479// [[15,85,0], [33,100, 0], [ 67,100, 0], [ 85, 85,0]],
1480// ];
1481// patch2 = [
1482// [[15,85,0], [33,100, 0], [ 67,100, 0], [ 85, 85,0]],
1483// [[ 0,67,0], [33, 67,-50], [ 67, 67,-50], [100, 67,0]],
1484// [[ 0,33,0], [33, 33,-50], [ 67, 33,-50], [100, 33,0]],
1485// [[15,15,0], [33, 0, 0], [ 67, 0, 0], [ 85, 15,0]],
1486// ];
1487// debug_bezier_patches(patches=[patch1, patch2], splinesteps=8, showcps=true);
1488module debug_bezier_patches(patches=[], size, splinesteps=16, showcps=true, showdots=false, showpatch=true, convexity=10, style="default")
1489{
1490 no_children($children);
1491 assert(is_undef(size)||is_num(size));
1492 assert(is_int(splinesteps) && splinesteps>0);
1493 assert(is_list(patches) && all([for (patch=patches) is_bezier_patch(patch)]));
1494 assert(is_bool(showcps));
1495 assert(is_bool(showdots));
1496 assert(is_bool(showpatch));
1497 assert(is_int(convexity) && convexity>0);
1498 for (patch = patches) {
1499 size = is_num(size)? size :
1500 let( bounds = pointlist_bounds(flatten(patch)) )
1501 max(bounds[1]-bounds[0])*0.01;
1502 if (showcps) {
1503 move_copies(flatten(patch)) color("red") sphere(d=size*2);
1504 color("cyan")
1505 for (i=[0:1:len(patch)-1], j=[0:1:len(patch[i])-1]) {
1506 if (i<len(patch)-1) extrude_from_to(patch[i][j], patch[i+1][j]) circle(d=size);
1507 if (j<len(patch[i])-1) extrude_from_to(patch[i][j], patch[i][j+1]) circle(d=size);
1508 }
1509 }
1510 if (showpatch || showdots){
1511 vnf = bezier_vnf(patch, splinesteps=splinesteps, style=style);
1512 if (showpatch) vnf_polyhedron(vnf, convexity=convexity);
1513 if (showdots) color("blue") move_copies(vnf[0]) sphere(d=size);
1514 }
1515 }
1516}
1517
1518
1519// vim: expandtab tabstop=4 shiftwidth=4 softtabstop=4 nowrap