13.472J/1.128J/2.158J/16.940J COMPUTATIONAL GEOMETRY Lecture 6 N. M. Patrikalakis Massachusetts Institute of Technology Cambridge, MA 02139-4307, USA Copyright c 2003 Massachusetts Institute of Technology Contents 6 B-splines (Uniform and Non-uniform) 2 6.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 6.2 De nition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 6.2.1 Knot vector . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 6.2.2 Properties and de nition of basis function Ni;k(u) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 6.2.3 Example: 2nd order B-spline basis function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 6.2.4 Example: 3rd order B-spline basis function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 6.2.5 Example: 4th order basis function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 6.2.6 Special case n = k 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 6.3 Derivatives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 6.4 Approaches to design with B-spline curves . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 6.5 Interpolation of data points with cubic B-splines . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 6.6 Evaluation and subdivision of B-splines . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 6.6.1 De Boor algorithm for B-spline curve evaluation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 6.6.2 Knot insertion: Boehm’s algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 6.7 Tensor product piecewise polynomial surface patches . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 6.7.1 Example: B ezier patch . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 6.7.2 Example: B-Spline patch . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 6.7.3 Example: Composite B ezier patches . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 6.8 Generalization of B-splines to NURBS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 6.8.1 Curves and Surface Patches . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 6.8.2 Example: Representation of a quarter circle as a rational polynomial . . . . . . . . . . . . . . . . . . . . . . 20 6.8.3 Trimmed patches . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21 6.9 Comparison of free-form curve/surface representation methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 Bibliography 23 Reading in the Textbook Chapter 1, pp. 6 - 33 1 Lecture 6 B-splines (Uniform and Non-uniform) 6.1 Introduction The formulation of uniform B-splines can be generalized to accomplish certain objectives. These include Non-uniform parameterization. Greater general exibility. Change of one polygon vertex in a B ezier curve or of one data point in a cardinal (or interpolatory) spline curve changes entire curve (global schemes). Remove necessity to increase degree of B ezier curves or construct composite B ezier curves using previous schemes in order to increase degrees of freedom. Obtain a \local" approximation scheme. The development extends the B ezier curve formulation to a piecewise polynomial curve with easy continuity control. 2 6.2 De nition A parametric non-uniform B-spline curve is de ned by P(u) = nX i=0 PiNi;k(u) (6.1) where, Pi are n + 1 control points; Ni;k(u) are piecewise polynomial B-spline basis functions of order k (or degree k 1) with n k 1. Therefore n; k are independent, unlike B ezier curves. The parameter u obeys the inequality to u tn+k (6.2) 6.2.1 Knot vector For open (non-periodic) curves, it is usual to de ne a set T of non-decreasing real numbers which is called the knot vector, as follows: T = fto = t1 = = tk 1| {z } k equal values < tk tk+1 tn| {z } n k+1 internal knots < tn+1 = = tn+k| {z } k equal values g (6.3) At each knot value, the curve P(u) has some degree of discontinuity in its derivatives above a certain order as we will see later. The total number of knots is n + k + 1 which equals the number of control points in the curve plus the curve’s order. 6.2.2 Properties and de nition of basis function Ni;k(u) 1. Pni=0 Ni;k(u) = 1 (partition of unity). 2. Ni;k(u) 0 (positivity). 3. Ni;k(u) = 0 if u 62 [ti; ti+k] (local support). 4. Ni;k(u) is (k 2) times continuously di erentiable at simple knots. If a knot tj is of multiplicity ( k), ie. if tj = tj+1 = = tj+ 1 (6.4) then Ni;k(u) is (k 1) times continuously di erentiable, ie. it is of class Ck 1. 3 5. Recursive de nition (Cox-de Boor algorithm) Ni;1(u) = ( 1 u 2 [t i; ti+1) 0 u 62 [ti; ti+1) (6.5) Ni;k(u) = u tit i+k 1 ti Ni;k 1(u) + ti+k ut i+k ti+1 Ni+1;k 1(u) (6.6) (set 00 = 0 above when it occurs) Properties 1-4 or property 5 by itself (given a known knot vector T) de ne the B-spline basis. 6.2.3 Example: 2nd order B-spline basis function (Piecewise linear case{ see Figure 6.1) k = 2; Ck 2 = C2 2 = C0 ti 1 ti+1 ti+2 ti+3ti Ni;2(u) Figure 6.1: Plot of 2nd order B-spline basis functions. Ni;2(u) consists of two piecewise linear polynomials: Ni;2(u) = ( u ti ti+1 ti ti u ti+1t i+2 u ti+2 ti+1 ti+1 u ti+2 4 6.2.4 Example: 3rd order B-spline basis function (Piecewise quadratic case{ See Figure 6.2) k = 3; Ck 2 = C3 2 = C1 Ni;3(ti) consists of three piecewise quadratic polynomials y1(u), y2(u) and y3(u). We want to nd y1(u), y2(u) and y3(u) such that Figure 6.2: Plot of 3rd order B-spline basis function. Ni;3(ti) = N0i;3(ti) = 0: Ni;3(ti+3) = N0i;3(ti+3) = 0: (6.7) Position continuity: y1(u) = A u t i ti+1 ti !2 (6.8) y3(u) = B t i+3 u ti+3 ti+2 !2 (6.9) y2(u) = A(1 s)2 + Bs2 + C2s(1 s) (6.10) s = u ti+1t i+2 ti+1 (6.11) Note that y1(ti+1) = y2(ti+1) = A y3(ti+2) = y2(ti+2) = B 5 at s = 0; s = 1(u = ti+1; u = ti+2). First derivative continuity: 2(C A) 1t i+2 ti+1 = 2A 1t i+1 ti (6.12) 2(B C) 1t i+2 ti+1 = 2B 1t i+3 ti+2 (6.13) hence, A " 1 ti+1 ti + 1 ti+2 ti+1 # = C 1t i+2 ti+1 (6.14) B " 1 ti+2 ti+1 + 1 ti+3 ti+2 # = C 1t i+2 ti+1 (6.15) where A; B = functions(C) Need one more condition (normalization): Z ti+k or +1 ti or 1 Ni;k(u)du = 1k(ti+k ti) We obtain A(ti+2 ti) + B(ti+3 ti+1) + C(ti+2 ti+1) = ti+3 ti (6.16) From Equation 6.14 to 6.16, we obtain A = ti+1 tit i+2 ti ; B = ti+3 ti+2t i+3 ti+1 ; C = 1 (6.17) Finally, Ni;3(u) = y1(u)Ni;1(u) + y2(u)Ni+1;1(u) + y3(u)Ni+2;1(u) = u tit i+2 ti Ni;2(u) + ti+3 ut i+3 ti+1 Ni+1;2(u) (6.18) 6 6.2.5 Example: 4th order basis function (Cubic B-spline case{ K = 4 see Figure 6.3) n = 6 ! 7 control points n + k + 1 = 11 knots t 1 o t1 t2 t3 t4 t7 t8 t9 t10 t6t5 N0,4 N1,4 N2,4 N3,4 N4,4 N 5,4 N6,4 N0,4 is C-1 N0,4 is C2 N1,4 is C2 N5,4 is C2 N2,4 is C2 N6,4 is C0N 1,4 is C 2 N4,4 is C2 N1,4 is C0 N2,4 is C1 N3,4 is C2 N3,4 is C2 N4,4 is C1 N5,4 is C0 N6,4 is C-1 Figure 6.3: Cubic B-spline functions. From property 1: _N0;4(t0) + _N1;4(t0) = 0 (6.19) Therefore, _P(t0) = (P1 P0) _N1;4(t0) (6.20) Similarly, _P(t10) = (P6 P5) _N6;4(t10) (6.21) 7 u=t0=t1=t2=t3 t4 t5 t6 P0 P1 P2 P3 P4 P5 u=t7=t8=t9=t10 P6 span 1 affected by P0,P1,P2,P3 span 2 affected by span 3 affected by span 4 affected by P1,P2,P3,P4 P2,P3,P4,P5 P3,P4,P5,P6 Figure 6.4: Example of local convex hull property of B-spline curve 8 Properties: Local support (eg. P6 a ects span 4 only), see Figure 6.4 i.e. Pi a ects [ti; ti+k] (k spans) Convex hull (stronger that B ezier) let u 2 [ti; ti+1], then Nj;k(u) 6= 0 for j 2 i k + 1; ; i (k values) iX j=i k+1 Nj;k(u) = 1; Nj;k(u) 0 Each span is in the convex hull of the k vertices contributing to its de nition Consequence: k consecutive vertices are collinear ! span is a straight line segment Variation diminishing property as for B ezier curves Exploit knot multiplicity to make complex curves 6.2.6 Special case n = k 1 The B-spline curve is also a B ezier curve in this case. T = fto = t1 = = tk 1| {z } k equal knots < tk = tk+1 = = t2k 1| {z } k equal knots g 6.3 Derivatives P(u) = nX i=0 PiNi;k(u) (6.22) P0(u) = nX i=1 diNi;k 1(u) (6.23) di = (k 1)Pi Pi 1t i+k 1 ti (6.24) 9 6.4 Approaches to design with B-spline curves Design procedure A 1. Designer chooses knot vector and control points. 2. Designer displays curve and tweaks control points to improve curve. Design procedure B 1. Designer starts with data points on or near curve. 2. Construct an interpolating/approximating B-spline curve. 3. Display curve and tweak control points to improve curve. 10 6.5 Interpolation of data points with cubic B-splines Given N data points: Pi; i = 1; 2; ; N, and no other derivative data at the boundaries. The problem is to construct a cubic B-spline curve which interpolates (precisely matches) these data points. Construction of knot vector Let ^u1 = 0 ^ui+1 = ^ui + di+1 di+1 = jPi+1 Pij i = 1; 2; ; N 1 d = NX i=2 di ! ui = ^uid i = 1; 2; ; N 1 Remove two knot values u2 and uN 1 from knot vector to obtain proper number of degrees of freedom (instead of having to prescribe boundary conditions, which could be done if needed) T = fu1 = u1 = u1 = u1| {z } 4 times < u3 u4 uN 2| {z } internal knots=N 4 < uN = uN = uN = uN| {z } 4 times g knots = n + 4 + 1 = (N 4) + 4 + 4 = N + 4 ! n = N 1 R(u) = nX i=0 RiBi;4(u) 0 u 1 (6.25) n = N 1 Require that R(uj) = N 1X i=0 RiBi;4(uj) = Pj j = 1; 2; ; N (6.26) Solve for Ri (system is banded). There are also other more sophisticated ways to choose knot vector and parameterization attempting to make u proportional to arc length. 11 6.6 Evaluation and subdivision of B-splines 6.6.1 De Boor algorithm for B-spline curve evaluation P12 P22 P32 t4 t t4 t1 P11 P0 P1 P2 P3 P13 t t1 t4 t1 t5 t t5 t3 t t3 t5 t3 t4 t t4 t3 t t3 t4 t3 t t2 t4 t2 t5 t t5 t2 t t2 t5 t2 t6 t t6 t3 t t3 t6 t3 t4 t t4 t2 P33 Figure 6.5: Diagram of de Boor subdivision algorithm over a cubic spline segment, where t 2 [t3; t4]. Evaluation of B-spline curves can be performed as follows (de Boor’s algorithm): P(u) = nX i=0 PiNi;k(u) t0 u tn+k (6.27) Let u 2 [tl; tl+1) be a particular span. Ni;k(u) 6= 0 for i = l k + 1; ; l Let P0i = Pi (6.28) De Boor’s recursive formula: Pji = (1 ji)Pj 1i 1 + jiPj 1i ; i l k + 2 (6.29) 12 where, ji = u tit i+k j ti (6.30) ) Pk 1k 1 = P(t) This algorithm is shown graphically in Figure 6.5. An example of the evaluation of a point on an cubic (k = 4, l = 3) B-spline curve is shown in Figure 6.6. P11 P22 P3 3 P32 P31 P30 P20 P10 P2 1 P00 t0=t1=t2=t3 t4 t 5 t6u Figure 6.6: Evaluation of a point on a B-spline curve with the de Boor algorithm. The algorithm shown in Figure 6.5 also permits the splitting of the segment into two smaller segments of the same order: left polygon: P 00 P11 P22 P33 right polygon: P 33 P23 P13 P03 13 6.6.2 Knot insertion: Boehm’s algorithm nX i=0 PiNi;k(t) | {z } over [ti]=[t0; ;tl;tl+1; ] = n+1X i=0 ^Pi ^Ni;k(t) | {z } over [ti]=[t0; ;tl;^t;tl+1; ] (6.31) ^Pi = aiPi + (1 ai)Pi 1 (6.32) ai = 8 >< >: 1 i l k + 1 0 l + 1 i t ti ti+k 1 ti l k + 2 i l P0 = ^P0 P3 = ^P4 ^P3 P2 ^P1 t2 t 3 ^t < t4 P1 ^P2 t4^t t3 t 5 Figure 6.7: Control polygon of B-spline curve before and after knot insertion in interval [t3; t4]. By repeatedly applying knot insertion algorithm, multiple knots within the knot vector can be created. P1 P2 P3 P0 P1 P2 P3 t ?t4 t ?t4 1 t ?t1 t ?t4 1 t ?t5 t ?t5 2 t ?t2 t ?t5 2 t ?t6 t ?t6 3 t ?t t ?t6 3 3 Figure 6.8: Boehm’s algorithm diagram with added knot t in interval [t3; t4]. 14 6.7 Tensor product piecewise polynomial surface patches Let R(u) = nX i=0 RiFi(u) A u B (6.33) be 3-D (or 2-D) curve expressed as linear combinations of basis functions Fi(u). Let this curve sweep a surface by moving and possibly deforming. This can be described by letting each Ri trace a curve Ri(v) = mX k=0 aikGk(v) C v D (6.34) The resulting surface is a tensor product surface (see Figure 6.9). R(u; v) = nX i=0 mX k=0 aikFi(u)Gk(v) (u; v) 2 [A; B] [C; D] (6.35) A B u v C D (u,v) v=D u=A u=BR(u,v) Patch v=C Figure 6.9: A tensor product patch. 15 Matrix (tensor) form: R(u; v) = [F0(u) Fn(u)][aik][G0(v) Gm(v)]T (6.36) The basis functions can be: Monomials ! Ferguson patch Hermite ! Hermite-Coons patch Bernstein ! B ezier patch Lagrange ! Lagrange patch Uniform B-splines ! Uniform B-spline patch Non-uniform B-splines ! Non-uniform B-spline patch. 16 6.7.1 Example: B ezier patch R(u; v) = nX i=0 mX j=0 RijBi;n(u)Bj;m(v) 0 u; v 1 (6.37) where Rij are the control points creating a control polyhedron (net), see Figure 6.10. 0 1 v u R00 1 R10 R30 R33 Figure 6.10: A B ezier patch. Properties: Lines of v= v = constant (isoparameter lines) are B ezier curves of degree n with control points Qi = mX j=0 RijBj;m( v): (6.38) The boundary isoparameter lines have the same control points as the corresponding polyhedron points. The relation between the patch and B ezier net is a nely invariant (translation, rotation, scaling). Convex hull. No known variation diminishing property. The procedure to create piecewise continuous surfaces with B ezier patches is complex. 17 6.7.2 Example: B-Spline patch R(u; v) = nX i=0 mX j=0 RijNi;k(u)Qj;l(v) n k 1 and m l 1 U = fu0 = u1 = = uk 1| {z } k equal knots < uk uk+1 un < un+1 = = un+k| {z } k equal knots g V = fv0 = v1 = = vl 1| {z } l equal knots < vl vl+1 vm < vm+1 = = vm+l| {z } l equal knots g Properties: Obeys same properties as a B ezier patch to which it reduces for n = k 1 and m = l 1 (isoparameter lines, boundaries, a ne invariance, convex hull). Easy construction of complex piecewise continuous geometries. Local control: 1. Rij a ects [ui; ui+k] [vj; vj+l]; 2. Subpatch [ui; ui+1] [vj; vj+1] a ected by Rp;q where (p; q) 2 [i k + 1; i] [j l + 1; j]. Strong convex hull { each subpatch lies in the convex hull of the vertices contributing to its de nition. 18 6.7.3 Example: Composite B ezier patches u v R(1)(u,v) (u,v)R(2) 0 ≤ u,v ≤ 1 Figure 6.11: Composite B ezier patches. 1. Positional continuity R(1)(1; v) = R(2)(0; v) (6.39) R(1)3i = R(2)0i i = 0; 1; 2; 3 (6.40) 2. Tangent plane (or normal) continuity R(2)u (0; v) R(2)v (0; v) = (v)R(1)u (1; v) R(1)v (1; v) (6.41) (direction of surface normal continuous for (v) > 0) Since R(2)v (0; v) = R(1)v (1; v); (6.42) use R(2)u (0; v)| {z } cubic in v = (v)|{z} constant R(1)u (1; v)| {z } cubic in v (6.43) to show R(2)1i R(2)0i = (R(1)3i R(1)2i ): (6.44) Therefore, collinearity of above polyhedron edges is required 19 6.8 Generalization of B-splines to NURBS The acronym NURBS stands for non-uniform rational B-splines. These functions have the Same properties as B-splines, and Are capable of representing a wider class of geometries. 6.8.1 Curves and Surface Patches NURBS curves are de ned by R(u) = Pn i=0 wiRiNi;k(u)P n i=0 wiNi;k(u) (6.45) where weights wi > 0; if all wi = 1, the integral piecewise polynomial spline case is recovered. This formulation permits exact representation of conics, eq. circle, ellipse, hyperbola. NURBS surface patches are de ned by R(u; v) = Pn i=0 Pm j=0 wijRijNi;k(u)Qj;l(v)P n i=0 Pm j=0 wijNi;k(u)Qj;l(v) (6.46) where weights wij 0. This formulation allows for exact representation of quadrics, tori, surfaces of revolution and very general free-form surfaces. If all wij = 1, the integral piecewise polynomial case is recovered. 6.8.2 Example: Representation of a quarter circle as a rational polynomial y x R Figure 6.12: Quarter circle. Consider a quarter circle (see Figure 6.12) described in terms of trigonometric functions by x = R cos( ) y = R sin( ) ) for 0 2: 20 Setting t = tan( 2) and using basic identities from trigonometry, we can express x and y as functions of t: x(t) = R1 t21+t2 y(t) = R 2t1+t2 9> = >; for 0 t 1: (6.47) For the conversion of Equation 6.47 to the B ezier representation, apply [t2 t 1] 2 64 c2c 1 c0 3 75 = [t2 t 1] 2 64 1 2 1 2 2 0 1 0 0 3 75 2 64 b0b 1 b2 3 75 separately to numerators and denominators to obtain the B ezier form. 6.8.3 Trimmed patches R(u; v) is an untrimmed patch in the parametric domain (u; v) 2 [A; B] [C; D]. Describe external loop as a set of edges (ie. curves in parameter space ri(t) = [ui(ti)vi(ti)] { eg. external loop the Figure 6.13 if made up of fr1;r2;r3;r4;r5g, while the internal loop is made up of curve fr6g. r6 Trimming lines r2 r4 r3 r5 r1 C A B D Figure 6.13: Trimmed surface patch. 21 6.9 Comparison of free-form curve/surface representa- tion methods Single span/patch Composite Ferguson (monomial or power basis) B ezier Hermite Cardinal or interpolatory spline B ezier B-spline Lagrange Table 6.1: Classi cation of free-form curve/surface representation. Property Ferguson Hermite- B ezier Lagrange Composite Cardinal Coons B ezier Spline B-Splines NURBS Easy geometric representation low med high med high Medium high high Convex hull no no yes no yes no yes yes Variation diminishing no no yes no yes no yes yes Easiness for creation low med med inappr. med high high high Local yes but control no no no no complex no yes yes Approximation ease med med high low high medium high high Interpolation high but ease med med med inappr. med high high high Generality med med med med med med med high Popularity low low med low med med high very high Table 6.2: Comparison of curve/surface representation methods. Variation diminishing property does not apply to surfaces. Popularity in industry and STEP/PDES standards. 22 Bibliography [1] G. Farin. Curves and Surfaces for Computer Aided Geometric Design: A Practical Guide. Academic Press, Boston, MA, 3rd edition, 1993. [2] I. D. Faux and M. J. Pratt. Computational Geometry for Design and Manufacture. Ellis Horwood, Chichester, England, 1981. [3] J. Hoschek and D. Lasser. Fundamentals of Computer Aided Geometric Design. A. K. Pe- ters, Wellesley, MA, 1993. Translated by L. L. Schumaker. [4] J. Owen. STEP: An Introduction. Information Geometers, Winchester, UK, 1993. [5] L. A. Piegl and W. Tiller. The NURBS Book. Springer, New York, 1995. [6] F. Yamaguchi. Curves and Surfaces in Computer Aided Geometric Design. Springer-Verlag, NY, 1988. 23