- Email: [email protected]
Circle approximation for CADCAM using orthogonal C2 cubic B-splines
International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229 www.elsevier.com/locate/ijmactool
Circle approximation for CADCAM using orthogonal C2 cubic B-splines Robert J. Cripps*, Peter S. Lockyer Geometric Modelling Group, School of Engineering, Mechanical and Manufacturing, University of Birmingham, Edgbaston, Birmingham B15 2TT, UK Received 18 October 2004; accepted 3 November 2004 Available online 9 March 2005
Abstract A method to approximate a circular arc of any subtended angle with orthogonal C2 cubic B-splines is presented. The approach is to specify end points and end tangent directions and to give extra points lying on the span to ensure the required level of accuracy. The critical elements are to use a parameterisation that is optimum for circular data and to determine the tangent magnitudes at the end points. The data can be sampled in any way, although evenly spaced data requires less points to achieve a given accuracy. The method is general and does not depend on any iterative schemes to determine the parameterisation, knot values or tangent magnitudes. q 2004 Elsevier Ltd. All rights reserved.
1. Introduction Circular arcs play an important role in the design and construction of computer-based components. Most current CADCAM systems support both analytic and parametric elements. Thus it is a straightforward task to join two curves together. A typical scenario may be that a designer has to construct a boundary curve that is constrained at each end by specified circular arcs, generally of different radii. An example of such a construction comes from BAE SYSTEMS. Air intakes and ducts are generally constructed by first defining a 2D profile which is then swept around a 2D guide to form a skin. The 2D profile is typically constructed using a collection of ellipses, splines and circular arcs blended with G1 continuity. A typical intake profile is shown in Fig. 1(a). The upper (outside) part is a spline section, the lower (inside) part is an ellipse. Fig. 1(b) shows a close up of the leading edge circular arc and the tangent continuous connections. Fig. l(c) is a 3D view of the section and the guide (composed of three filleted lines) around which the section will be swept. The resulting intake surface is shown in Fig. 2. It is possible to construct the B-spline blend curve between two circular arcs but difficulties arise when * Corresponding author. Tel.: C44 121 414 4223; fax: C44 121 414 3515. E-mail address: [email protected] (R.J. Cripps). 0890-6955/$ - see front matter q 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijmachtools.2004.11.010
the CADCAM system has to perform operations on the resulting boundary curve. Two approaches are common. First, the arcs and B-spline are stored as separate entities and any operations on the curve will use separate algorithms for circles and splines. The second approach is to approximate the circular arcs with B-spline curves so that only one entity needs be supported. Of course, the designer does not care which is adopted so long as the computer-based model corresponds to their design. From a software point of view it is far more efficient to minimise the number of entities that need supporting. It is therefore more common to approximate the circular arcs with B-splines. This gives rise to a problem of how well we can approximate the required circular arc. If the CADCAM system supports true NURBS entities then the circle can be represented exactly using the weights of the rational form. However, most CADCAM systems and most users of CADCAM systems use non-rational NUBS in preference to NURBS. This introduces a tolerance issue and eventually leads to a compromise between the amount of data, the number of segments and the degree of the B-splines used to achieve a specified tolerance. Clearly data proliferation is undesirable in terms of processing, data management and data manipulation all of which can have adverse effects on follow-on activities. Thus the goal is to approximate a circular arc to within a user specified tolerance with the minimum of data points or curve segments. The issue of specifying a suitable tolerance
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
1223
cubic B-spline curves, which are the lowest degree splines that match C2 continuity. There have been many attempts to provide an acceptable solution to the problem of approximating a circular arc for use in CADCAM systems. Specifically the problem is to find a polynomial-based representation of the 2D planar function: uðuÞ Z ðcosðuÞ; sinðuÞÞ
Fig. 1. Aircraft engine intake profile definition.
is by no means straightforward and is highly dependent on the application. For example in the automotive industries very high tolerances are required to achieve A-class surfaces for exterior car panels whereas a toolmaker for the casting industry will probably use more relaxed tolerances. Experience with companies operating across the spectrum of tolerance requirements [1] indicates that matching curvatures to within 5% is sufficient to satisfy the exacting tolerances required by aeronautical and automotive industries. Further, these industries desire C2 continuity between elements. Thus our aim is to approximate circular arcs such that the curvature profiles are within 5% using non-rational
0% u% 2p
As reported by Piegl and Tiller [9] previous methods have used cubic splines, Hermite polynomials and quadratic, cubic and quintic Bezier curves. Fang [3] gives a comparison of several of these methods. Piegl and Tiller [9] proposed a new method to approximate circular arcs for evenly spaced data by interpolating end points and end derivatives using nonrational B-splines. The approach is to specify the order of the B-spline and the required level, k, of continuity by specifying the first k end derivatives. Their approach is algorithmic rather than a closed form solution aimed at achieving a high level of precision of approximation in terms of distance from the circular arc. A key problem with all methods is the specification of a suitable parameterisation. The common parameterisation schemes are distance-based, i.e. the parameterisation is based on proportional arc lengths. Some, like the approach of Piegl and Tiller [9], use some optimisation process to improve the proportional arc lengths by adjusting the knot vector or the end tangent magnitudes. However, rather than resort to iterative schemes, an alternative approach has been proposed by Cripps and Ball [2]. It is a purely mathematical construction to determine a closed form of parameterisation that is optimal for circular arcs. The basic idea is to construct a parameterisation that tries to make the first and second derivatives orthogonal at the knots of the resulting nonrational cubic B-spline curve. The rest of the article is organised as follows. We give the basics of NUBS interpolation and briefly review distance-based parameterisation. Then we introduce the orthogonal parameterisation scheme and illustrate its effectiveness in approximating circular arcs of varying spans. We conclude with some suggestions for extending the method for approximating non-circular data.
2. NUBS basics A NUBS curve is defined by: CðuÞ Z
n X Ni;3 ðuÞPi ; ðuk % u% unC1 Þ iZ0
Fig. 2. Aircraft engine intake surface.
where Pi are the control points, u0%u%unC1 is the parameter and Ni,k(u) are the kth normalised B-spline
1224
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
basis functions defined by the knot vector: UZ½u0 Z/Zuk Z0;ukC1 ;.;un ;unC1 Z/ZunCkC1 Z1 (1) Given a string of points pi; iZ0,1,.,n along with two end tangent directions, t0 and tn, find the control points Pi such that piZC(ui), i.e. the cubic B-spline curve interpolates the data points at the corresponding parameter values. In other words, assign the knot vector to be the parameter values. Specifying U and the end tangent magnitudes results in a system of linear equations which can be readily solved [8] giving the unique clamped cubic B-spline interpolation. Given that the basis functions are determined by the knots, which in turn are defined by the parameterisation, it is clear that the parameterisation is key in controlling the resulting NUBS curve. In the following the quality of the approximation is measured by the percentage discrepancies in the curvature profiles, k(u) which is a more demanding measuring than position discrepancy. The number of points required to achieve a tolerance of about 5% error in the resulting k(u) is also considered since our aim is to avoid data proliferation.
_ CðuÞ is a unit vector, so distance-based parameterisation schemes try to keep s 0 z constant, i.e. the segments of the B-spline curve have approximately equal proportional arc lengths. For evenly spaced data the span lengths will be nearly equal. For unevenly spaced data the span lengths between knots will clearly vary according to the spacing of the data. In some distance-based schemes, for example [9] an optimisation process is employed to improve the proportional arc lengths by adjusting the knot vector or the end tangent magnitudes. We propose an alternative approach.
4. Orthogonal parameterisation We wish to determine a parameterisation using a purely mathematical construction that is optimal for circular arcs. Since, for circular arcs: U 0 ðuÞ Z ðKsinðuÞ; cosðuÞÞ U 00 ðuÞ Z ðKcosðuÞ; K sinðuÞÞ thus
3. Parameterisation The most popular approach to constructing a parameterisation is to use a distance-based method that approximates the ‘ideal’ arc length parameterisation. Examples of distance based parameterisation schemes can be found in [4–8]. Of all the distance-based methods, the most commonly implemented is chord length parameterisation [9]. Let the approximate arc length be given by: dZ
n X
jjpi K piK1 jj
iZ1
Then the chord length parameterisation assigns: ui Z uiK1 C
jjpi K piK1 jj ; i Z 1; .; n K 1 d
with u0Z0 and unZ1. The end tangent magnitudes are given by the approximate arc length scaled by the appropriate parameter interval at each end, i.e. a0 Z
d ukC1 K uk
an Z
d : unC1 K un
U 00 0ðuÞ$U 0 ðuÞ Z ðKsinðuÞ; cosðuÞÞ$ðKcosðuÞ; K sinðuÞÞ Z sinðuÞ cosðuÞ K cosðuÞ sinðuÞ Z 0
c u 2½0; 2p;
i.e. the first and second parametric derivatives are orthogonal along the entire arc. Thus the idea is to construct a parameterisation that tries to maintain this orthogonality condition at the knots of the resulting cubic B-spline curve. Given pi (iZ0,.,1), let ai Z pi K piK1 ; biZaiC1 and ci Z ai C bi Z ai C aiC1 , as illustrated in Fig. 3, then our parameterisation [2] is given by: ½ðbi $ci Þ2 C gðai $ci Þ uiC1 Z ui C ðui K uiK1 Þ ½ðai $ci Þ2 C dðbi $ci Þ where 1 g Z jjbi jj2 jjci jj2 2 1 d Z jjai jj2 jjci jj2 2 are scalars determined to ensure C 00 (u)$C 0 (u)Z0 at the knots.
Now, distance based schemes try to keep the proportional distances between the resulting curve segments the same. This is evident since: 0 _ C 0 ðuÞ Z CðuÞ$s
where C 0 (u) is the derivative with respect to the parameter u _ and CðuÞ is the derivative with respect to distance, s. Now,
Fig. 3. Orthogonal construction.
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
The end tangent magnitudes are evaluated from: C 0 ðuj Þ Z
aj t ujC1 K uj j
where ½3ða1 $t0 Þjja1 jj2 ½3ðan $tn Þjjan jj2 ¼ a0 ¼ and a n ½2ða1 $t0 Þ þ jja1 jj2 ½2ðan $tn Þ þ jjan jj2 are scalars that ensure orthogonality at uZui; jZ0,n. Clearly when the data is evenly spaced, ui which is approximating proportional distance along the arc, will be evenly spaced and C 00 (u)$C 0 (u) is precisely zero at the knots. As the data becomes more unevenly spaced, the exact orthogonality condition at the knots will be lost.
5. Unit circular arc approximation—evenly spaced data Following Piegl and Tiller [9] the effect of arc span on the quality of the approximation is considered. However, we measure the quality of the approximation by the curvature error. The maximum curvature error for the B-spline approximation occurs between the knot values and a simple search along each spline segment determines the maximum error over the entire curve. This is a more demanding requirement than radial error. For example for an arc of 308 using 2 points the maximum curvature error is 24%. The corresponding maximum radial error between the circular arc and the spline is 0.00021 which matches the results for the unit circle quoted by [9]. Table 1 shows that as the arc span increases the number of internal points increases in order to match the curvature profile k(u) to within the specified tolerance. This increase is not linear with arc span. For a span of (p/3) requires 3 points, a span of (p/2) nZ4 points and a span of (2p/3), nZ6 points. However, as can be seen in Table 1 the curvature error for the (2p/3) span is about half that for the (p/3) span. From the results in Table 1 a heuristic rule for the number of points to use to approximate a unit circle of a given span q to within 5% curvature error is given by: 6q n Z 1 C int : p
1225
included in Table 1 is the curvature error for each span when the number of points is reduced by 1. This shows how sensitive the B-spline approximation is to the number of points and span length. To illustrate the quality of the approximation, the curvature profiles of the resulting B-spline curves are given. Fig. 4 gives arc spans of (p/12), (p/6), (p/4), nZ5, 6, 7 evenly spaced data points, respectively. Fig. 5 gives the corresponding profiles for spans (p/3), (p/2), (3p/4), p again for nZ5, 6, 7 evenly spaced data points. For spans of less than (p/4) the approximations are well within tolerance and little can be seen in the profiles. For arc spans greater than (p/4) the approximations are having to work much harder and whilst still within acceptable tolerance, the oscillations of curvature are becoming more marked. Note the nature of the curvature profiles and the oscillations about the target curvature k(u)Z1. Even when the curvature errors are within the desired tolerance, i.e. 5% and hence can be ignored in practical applications, Fig. 6 illustrates the oscillatory nature of using polynomial based interpolation and approximation schemes and clearly demonstrates the need to control both position and curvature accuracy to achieve high quality curve approximations. The orthogonality condition along the approximation spline can also be measured by considering the discrepancy between the angle between the first and second derivatives and 908: eangle Z cosK1 ðjÞ K
p 2
Encouragingly only 10 points are required on an arc to match the 5% curvature error for a full 2p span. Also Table 1 Number of points to achieve 5% curvature error for given arc spans Span (deg) n % Curvature error (n) % Curvature error (nK1)
30 2 2.4
60 3 2.4
90 4 2.4
120 6 1.5
180 6 3.5
270 8 4.1
360 10 4.4
n/a
11.1
5.7
2.4
5.7
5.7
5.7
Fig. 4. Curvature profile for arc spans (p/12), (p/6), (p/4), nZ5, 6, 7 evenly spaced data points.
1226
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
Fig. 5. Curvature profile for arc spans (p/3), (p/2), (p/4), p, nZ5, 6, 7: evenly spaced data points.
where j is given by jZ
½C 00 ðuÞ$C 0 ðuÞ sC 00 ðuÞs$sC 0 ðuÞs
For evenly spaced data, the angular error will be zero at the knots, i.e. the data points, and oscillate between the knots. For unevenly spaced data the error profiles still oscillate about zero, but the position of the roots is more difficult to predict. For extreme unevenly spaced data some of the roots may not occur in the knot interval [uj, ujC1). A simple search determines the maximum angular error along the entire spline. In general C 00 (u)$C 0 (u) is a cubic (B-spline) polynomial in the parameter u that has at most three roots. The degenerate cases where C(u) is quadratic or linear will affect the number of roots in C 00 (u)$C 0 (u). When C(u) is quadratic there is at most 1 root and when C(u) is linear there are an infinity of roots. These degenerate cases are unlikely to occur when sampling data from circular arcs of reasonable span lengths. Given C 00 (u)$C 0 (u) is a cubic, each segment of the Bspline curve on the knot interval [uj, ujC1) will have at most three roots, i.e. the angular error eangle is zero at three parameter values in each segment. For evenly spaced data, these are at uj, ujC1 and by symmetry u * Z ð1=2Þðuj C ujC1 Þ Further, there must be turning points in the intervals (uj, u*) and (u*, uj). Again, symmetry tells us that these turning points must be at u*Gt(ujC1Kuj) for 0%t%1. Numerical evidence suggests that the turning points occur when tZ (1/4) with the minimum angular discrepancy at (1/4) (3ujC ujC1) and the maximum at (1/4) (ujC3ujC1). This phenomenon can be seen in Fig. 7. Geometrically Fig. 7 suggests that for each segment of the B-spline the approximation curve starts with the exact orthogonality condition. Moving away from the knot, the angle between C 00 (u) and C 0 (u) is less than 908, and reaches a maximum at a quarter of the knot interval, passes through
Fig. 6. Effect of span on curvature.
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
1227
Fig. 7. % Angular error for range of unit arc spans.
zero at (1/2) (ujCujC1), then increases, i.e. becomes greater than 908 before returning to zero at (ujC1,). Clearly as the number of points increase, the length of each segment reduces and the oscillations about 908 reduce.
6. Unit circular arc approximation—unevenly spaced data The orthogonal parameterisation is optimum when the data is evenly spaced since this guarantees that the orthogonality conditions is maintained at the knots. Since the method is general, it can be used when the data is unevenly spaced. In this case the orthogonality is lost at the knots, i.e. some or all of the roots of C 00 (u)$C 0 (u) may lie outside [uj, ujC1). The effect by approximating a unit circular arc of span 908 is illustrated by choosing four points separated by 21, 28 and 418, respectively. These were chosen to avoid any bias in the span lengths. The resulting B-spline approximation has three segments of proportional arc lengths 0.233, 0.311 and 0.455 respectively. The resulting curvature and angular
errors at the knots are given in Table 2. The maximum curvature error along the B-spline is 5.1% which is just outside our acceptable tolerance. Compare this with a maximum curvature error of just 2.4% for four evenly spaced data points (left hand of Table 2). Table 2 indicates that the greater the length between adjacent points the higher the curvature discrepancy. This is intuitively as expected since the spline segment is having to do more work over the longer span with the same amount of information. The corresponding curvature plots in Fig. 8 highlight this phenomenon as do the angular discrepancies in Figs. 9 and 10. Note also how the unevenly spaced data has affected the symmetry of the orthogonality error. As the data becomes more extreme in spacing, there may be segments on the B-spline curve that never achieve orthogonality between C 00 (u) and C 0 (u). The parameterisation has been optimised for evenly spaced data. Clearly unevenly spaced data can be used but to
Table 2 Data spacing: Unit circle span (p/2)4 points ui
0.000 0 333 0.666 1.000
Evenly
ui
% Curvature error
% Angular error
2.39 2.39 2.39 2.39
0.0 0.0 0.0 0.0
0.000 0.233 0.544 1.000
Unevenly % Curvature error
% Angular error
1.04 1.26 3.63 5.07
0.04 0.13 0.13 0.18
Fig. 8. Curvature profile for unit span (p/2) (a) evenly (b) unevenly spaced data.
1228
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
Fig. 11. Radius of curvature profiles r(u): spiral, evenly spaced data. Orthogonal (nZ4, 5, 6).
Fig. 9. % angular error for unit span (p/2): evenly spaced data.
Fig. 10. % angular error for unit span (p/2): unevenly spaced data.
match such a span to a given tolerance will require more points than for evenly spaced data.
7. Development The table of errors in curvature profiles show that the orthogonal interpolation is able to approximate the unit circle of specified span to within a curvature error of G0.05 using at most 10 evenly spaced data points. The circle is used in CADCAM and the design of components because of its simple curvature profile, k(s)Zc. However, this is quite limiting and a curve with a more general curvature profile is of more use. Applying the orthogonal parameterisation to approximate a spiral gives acceptable results at the cost of using extra data points [2]. Fig. 11 illustrates the resulting radius of curvature of the orthogonal cubic B-spline curve using nZ4, 5 and 6 evenly spaced points from a 1808 spiral
arc from radius rZ1 to rZ2 of span lengthZ(p/2). However, as Fig. 11 clearly shows, the errors grow as we move along the curve and suggests that evenly spaced data is not appropriate for the spiral. Spacing the data based on the turning angle may help to distribute the errors more evenly. The orthogonal parameterisation was optimised for circular segments and therefore could not be expected to work as well for spiral data. There is therefore scope to improve the approximation, especially when approximating spiral data. One approach is to modify the resulting B-spline curve. It is possible to raise the orthogonal cubic B-spline to a quintic C2 B-spline. This generates multiple knots in U, but this is not a problem for CADCAM systems. The C2 quintic also gives two extra vector degrees of freedom. These can be manipulated independently without affecting the continuity across spans and thus gives flexibility to improve the curvature profiles for the spline. Thus, the optimum choice of internal control points with respect to spiral approximation can be determined. Whilst this approach will work, it relies on optimisation or iterative techniques.
8. Conclusions A parameterisation scheme has been presented that is optimal with respect to approximating circular arcs. Our aim was to approximate circular arcs with constrained piecewise cubic C2 interpolating curves. By constructing the parameterisation to approximate the feature of the circle, i.e. orthogonality of first and second derivatives at the knots we are able to closely approximate circular arcs by measuring the % error between curvature profiles. Our approximations are within 5% for the case studies shown. Encouragingly, the approach is suitable for both evenly and unevenly spaced data. The approach performs better than nonoptimised chord length parameterisation and compares with the method of Piegl and Tiller [9] for equivalent amounts of data. The proposed method is therefore worthy of further attention since it is stable and direct, requiring no iterative methods or optimisation. As expected it performs less well on spiral segments, but is good enough to suggest that an acceptable approximation to spiral segments can be achieved by using the orthogonal parameterised spline,
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
raising it to quintic and adjusting the degrees of freedom to suit.
Acknowledgements We are please to be able to acknowledge the support of EPSRC (Grant reference GR/R35735/01) and Delcam, UK. The aircraft cowling definition was kindly supplied by BAE SYSTEMS.
References [1] R.J. Cripps, Algorithms to support point-based CADCAM, Int. J. Machine Tool Manf. 43 (2003) 425–432.
1229
[2] R.J. Cripps, A.A. Ball. Orthogonal Cubic C2 B-splines. Proceedings of the 10th SIAM International Conference on Geometric Design. Seattle November; 2003. [3] L. Fang, Circular arc approximation by quintic polynomial curves, Comput. Aided Geom. Des. 15 (1998) 843–861. [4] T.A. Foley, G.M. Nielson, Knot selection for parametric spline interpolation in: T. Lyche, L.L. Schumaker (Eds.), Mathematical Methods in Computer Aided Geometric Design (1989), pp. 261–272. [5] E.T.Y. Lee, Choosing nodes in parametric curve interpolation, Comput. Aided Geom. Des. 21 (1989) 363–370. [6] C-G. Lim, A universal parametrization in B-spline curve and surface interpolation, Comput. Aided Geom. Des. 16 (1999) 407–422. [7] H. Park, Choosing nodes and knots in closed B-spline curve interpolation to point data, Comput. Aided Des. 33 (2001) 967–974. [8] L. Piegl, W. Tiller, The NURBS Book, Monographs in Visual Communication, Springer, New York, 1995. [9] L. Piegl, W. Tiller, Circular approximation using integral B-splines, Comput. Aided Geom. Des. 35 (2003) 601–607.
Circle approximation for CADCAM using orthogonal C2 cubic B-splines Robert J. Cripps*, Peter S. Lockyer Geometric Modelling Group, School of Engineering, Mechanical and Manufacturing, University of Birmingham, Edgbaston, Birmingham B15 2TT, UK Received 18 October 2004; accepted 3 November 2004 Available online 9 March 2005
Abstract A method to approximate a circular arc of any subtended angle with orthogonal C2 cubic B-splines is presented. The approach is to specify end points and end tangent directions and to give extra points lying on the span to ensure the required level of accuracy. The critical elements are to use a parameterisation that is optimum for circular data and to determine the tangent magnitudes at the end points. The data can be sampled in any way, although evenly spaced data requires less points to achieve a given accuracy. The method is general and does not depend on any iterative schemes to determine the parameterisation, knot values or tangent magnitudes. q 2004 Elsevier Ltd. All rights reserved.
1. Introduction Circular arcs play an important role in the design and construction of computer-based components. Most current CADCAM systems support both analytic and parametric elements. Thus it is a straightforward task to join two curves together. A typical scenario may be that a designer has to construct a boundary curve that is constrained at each end by specified circular arcs, generally of different radii. An example of such a construction comes from BAE SYSTEMS. Air intakes and ducts are generally constructed by first defining a 2D profile which is then swept around a 2D guide to form a skin. The 2D profile is typically constructed using a collection of ellipses, splines and circular arcs blended with G1 continuity. A typical intake profile is shown in Fig. 1(a). The upper (outside) part is a spline section, the lower (inside) part is an ellipse. Fig. 1(b) shows a close up of the leading edge circular arc and the tangent continuous connections. Fig. l(c) is a 3D view of the section and the guide (composed of three filleted lines) around which the section will be swept. The resulting intake surface is shown in Fig. 2. It is possible to construct the B-spline blend curve between two circular arcs but difficulties arise when * Corresponding author. Tel.: C44 121 414 4223; fax: C44 121 414 3515. E-mail address: [email protected] (R.J. Cripps). 0890-6955/$ - see front matter q 2004 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijmachtools.2004.11.010
the CADCAM system has to perform operations on the resulting boundary curve. Two approaches are common. First, the arcs and B-spline are stored as separate entities and any operations on the curve will use separate algorithms for circles and splines. The second approach is to approximate the circular arcs with B-spline curves so that only one entity needs be supported. Of course, the designer does not care which is adopted so long as the computer-based model corresponds to their design. From a software point of view it is far more efficient to minimise the number of entities that need supporting. It is therefore more common to approximate the circular arcs with B-splines. This gives rise to a problem of how well we can approximate the required circular arc. If the CADCAM system supports true NURBS entities then the circle can be represented exactly using the weights of the rational form. However, most CADCAM systems and most users of CADCAM systems use non-rational NUBS in preference to NURBS. This introduces a tolerance issue and eventually leads to a compromise between the amount of data, the number of segments and the degree of the B-splines used to achieve a specified tolerance. Clearly data proliferation is undesirable in terms of processing, data management and data manipulation all of which can have adverse effects on follow-on activities. Thus the goal is to approximate a circular arc to within a user specified tolerance with the minimum of data points or curve segments. The issue of specifying a suitable tolerance
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
1223
cubic B-spline curves, which are the lowest degree splines that match C2 continuity. There have been many attempts to provide an acceptable solution to the problem of approximating a circular arc for use in CADCAM systems. Specifically the problem is to find a polynomial-based representation of the 2D planar function: uðuÞ Z ðcosðuÞ; sinðuÞÞ
Fig. 1. Aircraft engine intake profile definition.
is by no means straightforward and is highly dependent on the application. For example in the automotive industries very high tolerances are required to achieve A-class surfaces for exterior car panels whereas a toolmaker for the casting industry will probably use more relaxed tolerances. Experience with companies operating across the spectrum of tolerance requirements [1] indicates that matching curvatures to within 5% is sufficient to satisfy the exacting tolerances required by aeronautical and automotive industries. Further, these industries desire C2 continuity between elements. Thus our aim is to approximate circular arcs such that the curvature profiles are within 5% using non-rational
0% u% 2p
As reported by Piegl and Tiller [9] previous methods have used cubic splines, Hermite polynomials and quadratic, cubic and quintic Bezier curves. Fang [3] gives a comparison of several of these methods. Piegl and Tiller [9] proposed a new method to approximate circular arcs for evenly spaced data by interpolating end points and end derivatives using nonrational B-splines. The approach is to specify the order of the B-spline and the required level, k, of continuity by specifying the first k end derivatives. Their approach is algorithmic rather than a closed form solution aimed at achieving a high level of precision of approximation in terms of distance from the circular arc. A key problem with all methods is the specification of a suitable parameterisation. The common parameterisation schemes are distance-based, i.e. the parameterisation is based on proportional arc lengths. Some, like the approach of Piegl and Tiller [9], use some optimisation process to improve the proportional arc lengths by adjusting the knot vector or the end tangent magnitudes. However, rather than resort to iterative schemes, an alternative approach has been proposed by Cripps and Ball [2]. It is a purely mathematical construction to determine a closed form of parameterisation that is optimal for circular arcs. The basic idea is to construct a parameterisation that tries to make the first and second derivatives orthogonal at the knots of the resulting nonrational cubic B-spline curve. The rest of the article is organised as follows. We give the basics of NUBS interpolation and briefly review distance-based parameterisation. Then we introduce the orthogonal parameterisation scheme and illustrate its effectiveness in approximating circular arcs of varying spans. We conclude with some suggestions for extending the method for approximating non-circular data.
2. NUBS basics A NUBS curve is defined by: CðuÞ Z
n X Ni;3 ðuÞPi ; ðuk % u% unC1 Þ iZ0
Fig. 2. Aircraft engine intake surface.
where Pi are the control points, u0%u%unC1 is the parameter and Ni,k(u) are the kth normalised B-spline
1224
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
basis functions defined by the knot vector: UZ½u0 Z/Zuk Z0;ukC1 ;.;un ;unC1 Z/ZunCkC1 Z1 (1) Given a string of points pi; iZ0,1,.,n along with two end tangent directions, t0 and tn, find the control points Pi such that piZC(ui), i.e. the cubic B-spline curve interpolates the data points at the corresponding parameter values. In other words, assign the knot vector to be the parameter values. Specifying U and the end tangent magnitudes results in a system of linear equations which can be readily solved [8] giving the unique clamped cubic B-spline interpolation. Given that the basis functions are determined by the knots, which in turn are defined by the parameterisation, it is clear that the parameterisation is key in controlling the resulting NUBS curve. In the following the quality of the approximation is measured by the percentage discrepancies in the curvature profiles, k(u) which is a more demanding measuring than position discrepancy. The number of points required to achieve a tolerance of about 5% error in the resulting k(u) is also considered since our aim is to avoid data proliferation.
_ CðuÞ is a unit vector, so distance-based parameterisation schemes try to keep s 0 z constant, i.e. the segments of the B-spline curve have approximately equal proportional arc lengths. For evenly spaced data the span lengths will be nearly equal. For unevenly spaced data the span lengths between knots will clearly vary according to the spacing of the data. In some distance-based schemes, for example [9] an optimisation process is employed to improve the proportional arc lengths by adjusting the knot vector or the end tangent magnitudes. We propose an alternative approach.
4. Orthogonal parameterisation We wish to determine a parameterisation using a purely mathematical construction that is optimal for circular arcs. Since, for circular arcs: U 0 ðuÞ Z ðKsinðuÞ; cosðuÞÞ U 00 ðuÞ Z ðKcosðuÞ; K sinðuÞÞ thus
3. Parameterisation The most popular approach to constructing a parameterisation is to use a distance-based method that approximates the ‘ideal’ arc length parameterisation. Examples of distance based parameterisation schemes can be found in [4–8]. Of all the distance-based methods, the most commonly implemented is chord length parameterisation [9]. Let the approximate arc length be given by: dZ
n X
jjpi K piK1 jj
iZ1
Then the chord length parameterisation assigns: ui Z uiK1 C
jjpi K piK1 jj ; i Z 1; .; n K 1 d
with u0Z0 and unZ1. The end tangent magnitudes are given by the approximate arc length scaled by the appropriate parameter interval at each end, i.e. a0 Z
d ukC1 K uk
an Z
d : unC1 K un
U 00 0ðuÞ$U 0 ðuÞ Z ðKsinðuÞ; cosðuÞÞ$ðKcosðuÞ; K sinðuÞÞ Z sinðuÞ cosðuÞ K cosðuÞ sinðuÞ Z 0
c u 2½0; 2p;
i.e. the first and second parametric derivatives are orthogonal along the entire arc. Thus the idea is to construct a parameterisation that tries to maintain this orthogonality condition at the knots of the resulting cubic B-spline curve. Given pi (iZ0,.,1), let ai Z pi K piK1 ; biZaiC1 and ci Z ai C bi Z ai C aiC1 , as illustrated in Fig. 3, then our parameterisation [2] is given by: ½ðbi $ci Þ2 C gðai $ci Þ uiC1 Z ui C ðui K uiK1 Þ ½ðai $ci Þ2 C dðbi $ci Þ where 1 g Z jjbi jj2 jjci jj2 2 1 d Z jjai jj2 jjci jj2 2 are scalars determined to ensure C 00 (u)$C 0 (u)Z0 at the knots.
Now, distance based schemes try to keep the proportional distances between the resulting curve segments the same. This is evident since: 0 _ C 0 ðuÞ Z CðuÞ$s
where C 0 (u) is the derivative with respect to the parameter u _ and CðuÞ is the derivative with respect to distance, s. Now,
Fig. 3. Orthogonal construction.
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
The end tangent magnitudes are evaluated from: C 0 ðuj Þ Z
aj t ujC1 K uj j
where ½3ða1 $t0 Þjja1 jj2 ½3ðan $tn Þjjan jj2 ¼ a0 ¼ and a n ½2ða1 $t0 Þ þ jja1 jj2 ½2ðan $tn Þ þ jjan jj2 are scalars that ensure orthogonality at uZui; jZ0,n. Clearly when the data is evenly spaced, ui which is approximating proportional distance along the arc, will be evenly spaced and C 00 (u)$C 0 (u) is precisely zero at the knots. As the data becomes more unevenly spaced, the exact orthogonality condition at the knots will be lost.
5. Unit circular arc approximation—evenly spaced data Following Piegl and Tiller [9] the effect of arc span on the quality of the approximation is considered. However, we measure the quality of the approximation by the curvature error. The maximum curvature error for the B-spline approximation occurs between the knot values and a simple search along each spline segment determines the maximum error over the entire curve. This is a more demanding requirement than radial error. For example for an arc of 308 using 2 points the maximum curvature error is 24%. The corresponding maximum radial error between the circular arc and the spline is 0.00021 which matches the results for the unit circle quoted by [9]. Table 1 shows that as the arc span increases the number of internal points increases in order to match the curvature profile k(u) to within the specified tolerance. This increase is not linear with arc span. For a span of (p/3) requires 3 points, a span of (p/2) nZ4 points and a span of (2p/3), nZ6 points. However, as can be seen in Table 1 the curvature error for the (2p/3) span is about half that for the (p/3) span. From the results in Table 1 a heuristic rule for the number of points to use to approximate a unit circle of a given span q to within 5% curvature error is given by: 6q n Z 1 C int : p
1225
included in Table 1 is the curvature error for each span when the number of points is reduced by 1. This shows how sensitive the B-spline approximation is to the number of points and span length. To illustrate the quality of the approximation, the curvature profiles of the resulting B-spline curves are given. Fig. 4 gives arc spans of (p/12), (p/6), (p/4), nZ5, 6, 7 evenly spaced data points, respectively. Fig. 5 gives the corresponding profiles for spans (p/3), (p/2), (3p/4), p again for nZ5, 6, 7 evenly spaced data points. For spans of less than (p/4) the approximations are well within tolerance and little can be seen in the profiles. For arc spans greater than (p/4) the approximations are having to work much harder and whilst still within acceptable tolerance, the oscillations of curvature are becoming more marked. Note the nature of the curvature profiles and the oscillations about the target curvature k(u)Z1. Even when the curvature errors are within the desired tolerance, i.e. 5% and hence can be ignored in practical applications, Fig. 6 illustrates the oscillatory nature of using polynomial based interpolation and approximation schemes and clearly demonstrates the need to control both position and curvature accuracy to achieve high quality curve approximations. The orthogonality condition along the approximation spline can also be measured by considering the discrepancy between the angle between the first and second derivatives and 908: eangle Z cosK1 ðjÞ K
p 2
Encouragingly only 10 points are required on an arc to match the 5% curvature error for a full 2p span. Also Table 1 Number of points to achieve 5% curvature error for given arc spans Span (deg) n % Curvature error (n) % Curvature error (nK1)
30 2 2.4
60 3 2.4
90 4 2.4
120 6 1.5
180 6 3.5
270 8 4.1
360 10 4.4
n/a
11.1
5.7
2.4
5.7
5.7
5.7
Fig. 4. Curvature profile for arc spans (p/12), (p/6), (p/4), nZ5, 6, 7 evenly spaced data points.
1226
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
Fig. 5. Curvature profile for arc spans (p/3), (p/2), (p/4), p, nZ5, 6, 7: evenly spaced data points.
where j is given by jZ
½C 00 ðuÞ$C 0 ðuÞ sC 00 ðuÞs$sC 0 ðuÞs
For evenly spaced data, the angular error will be zero at the knots, i.e. the data points, and oscillate between the knots. For unevenly spaced data the error profiles still oscillate about zero, but the position of the roots is more difficult to predict. For extreme unevenly spaced data some of the roots may not occur in the knot interval [uj, ujC1). A simple search determines the maximum angular error along the entire spline. In general C 00 (u)$C 0 (u) is a cubic (B-spline) polynomial in the parameter u that has at most three roots. The degenerate cases where C(u) is quadratic or linear will affect the number of roots in C 00 (u)$C 0 (u). When C(u) is quadratic there is at most 1 root and when C(u) is linear there are an infinity of roots. These degenerate cases are unlikely to occur when sampling data from circular arcs of reasonable span lengths. Given C 00 (u)$C 0 (u) is a cubic, each segment of the Bspline curve on the knot interval [uj, ujC1) will have at most three roots, i.e. the angular error eangle is zero at three parameter values in each segment. For evenly spaced data, these are at uj, ujC1 and by symmetry u * Z ð1=2Þðuj C ujC1 Þ Further, there must be turning points in the intervals (uj, u*) and (u*, uj). Again, symmetry tells us that these turning points must be at u*Gt(ujC1Kuj) for 0%t%1. Numerical evidence suggests that the turning points occur when tZ (1/4) with the minimum angular discrepancy at (1/4) (3ujC ujC1) and the maximum at (1/4) (ujC3ujC1). This phenomenon can be seen in Fig. 7. Geometrically Fig. 7 suggests that for each segment of the B-spline the approximation curve starts with the exact orthogonality condition. Moving away from the knot, the angle between C 00 (u) and C 0 (u) is less than 908, and reaches a maximum at a quarter of the knot interval, passes through
Fig. 6. Effect of span on curvature.
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
1227
Fig. 7. % Angular error for range of unit arc spans.
zero at (1/2) (ujCujC1), then increases, i.e. becomes greater than 908 before returning to zero at (ujC1,). Clearly as the number of points increase, the length of each segment reduces and the oscillations about 908 reduce.
6. Unit circular arc approximation—unevenly spaced data The orthogonal parameterisation is optimum when the data is evenly spaced since this guarantees that the orthogonality conditions is maintained at the knots. Since the method is general, it can be used when the data is unevenly spaced. In this case the orthogonality is lost at the knots, i.e. some or all of the roots of C 00 (u)$C 0 (u) may lie outside [uj, ujC1). The effect by approximating a unit circular arc of span 908 is illustrated by choosing four points separated by 21, 28 and 418, respectively. These were chosen to avoid any bias in the span lengths. The resulting B-spline approximation has three segments of proportional arc lengths 0.233, 0.311 and 0.455 respectively. The resulting curvature and angular
errors at the knots are given in Table 2. The maximum curvature error along the B-spline is 5.1% which is just outside our acceptable tolerance. Compare this with a maximum curvature error of just 2.4% for four evenly spaced data points (left hand of Table 2). Table 2 indicates that the greater the length between adjacent points the higher the curvature discrepancy. This is intuitively as expected since the spline segment is having to do more work over the longer span with the same amount of information. The corresponding curvature plots in Fig. 8 highlight this phenomenon as do the angular discrepancies in Figs. 9 and 10. Note also how the unevenly spaced data has affected the symmetry of the orthogonality error. As the data becomes more extreme in spacing, there may be segments on the B-spline curve that never achieve orthogonality between C 00 (u) and C 0 (u). The parameterisation has been optimised for evenly spaced data. Clearly unevenly spaced data can be used but to
Table 2 Data spacing: Unit circle span (p/2)4 points ui
0.000 0 333 0.666 1.000
Evenly
ui
% Curvature error
% Angular error
2.39 2.39 2.39 2.39
0.0 0.0 0.0 0.0
0.000 0.233 0.544 1.000
Unevenly % Curvature error
% Angular error
1.04 1.26 3.63 5.07
0.04 0.13 0.13 0.18
Fig. 8. Curvature profile for unit span (p/2) (a) evenly (b) unevenly spaced data.
1228
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
Fig. 11. Radius of curvature profiles r(u): spiral, evenly spaced data. Orthogonal (nZ4, 5, 6).
Fig. 9. % angular error for unit span (p/2): evenly spaced data.
Fig. 10. % angular error for unit span (p/2): unevenly spaced data.
match such a span to a given tolerance will require more points than for evenly spaced data.
7. Development The table of errors in curvature profiles show that the orthogonal interpolation is able to approximate the unit circle of specified span to within a curvature error of G0.05 using at most 10 evenly spaced data points. The circle is used in CADCAM and the design of components because of its simple curvature profile, k(s)Zc. However, this is quite limiting and a curve with a more general curvature profile is of more use. Applying the orthogonal parameterisation to approximate a spiral gives acceptable results at the cost of using extra data points [2]. Fig. 11 illustrates the resulting radius of curvature of the orthogonal cubic B-spline curve using nZ4, 5 and 6 evenly spaced points from a 1808 spiral
arc from radius rZ1 to rZ2 of span lengthZ(p/2). However, as Fig. 11 clearly shows, the errors grow as we move along the curve and suggests that evenly spaced data is not appropriate for the spiral. Spacing the data based on the turning angle may help to distribute the errors more evenly. The orthogonal parameterisation was optimised for circular segments and therefore could not be expected to work as well for spiral data. There is therefore scope to improve the approximation, especially when approximating spiral data. One approach is to modify the resulting B-spline curve. It is possible to raise the orthogonal cubic B-spline to a quintic C2 B-spline. This generates multiple knots in U, but this is not a problem for CADCAM systems. The C2 quintic also gives two extra vector degrees of freedom. These can be manipulated independently without affecting the continuity across spans and thus gives flexibility to improve the curvature profiles for the spline. Thus, the optimum choice of internal control points with respect to spiral approximation can be determined. Whilst this approach will work, it relies on optimisation or iterative techniques.
8. Conclusions A parameterisation scheme has been presented that is optimal with respect to approximating circular arcs. Our aim was to approximate circular arcs with constrained piecewise cubic C2 interpolating curves. By constructing the parameterisation to approximate the feature of the circle, i.e. orthogonality of first and second derivatives at the knots we are able to closely approximate circular arcs by measuring the % error between curvature profiles. Our approximations are within 5% for the case studies shown. Encouragingly, the approach is suitable for both evenly and unevenly spaced data. The approach performs better than nonoptimised chord length parameterisation and compares with the method of Piegl and Tiller [9] for equivalent amounts of data. The proposed method is therefore worthy of further attention since it is stable and direct, requiring no iterative methods or optimisation. As expected it performs less well on spiral segments, but is good enough to suggest that an acceptable approximation to spiral segments can be achieved by using the orthogonal parameterised spline,
R.J. Cripps, P.S. Lockyer / International Journal of Machine Tools & Manufacture 45 (2005) 1222–1229
raising it to quintic and adjusting the degrees of freedom to suit.
Acknowledgements We are please to be able to acknowledge the support of EPSRC (Grant reference GR/R35735/01) and Delcam, UK. The aircraft cowling definition was kindly supplied by BAE SYSTEMS.
References [1] R.J. Cripps, Algorithms to support point-based CADCAM, Int. J. Machine Tool Manf. 43 (2003) 425–432.
1229
[2] R.J. Cripps, A.A. Ball. Orthogonal Cubic C2 B-splines. Proceedings of the 10th SIAM International Conference on Geometric Design. Seattle November; 2003. [3] L. Fang, Circular arc approximation by quintic polynomial curves, Comput. Aided Geom. Des. 15 (1998) 843–861. [4] T.A. Foley, G.M. Nielson, Knot selection for parametric spline interpolation in: T. Lyche, L.L. Schumaker (Eds.), Mathematical Methods in Computer Aided Geometric Design (1989), pp. 261–272. [5] E.T.Y. Lee, Choosing nodes in parametric curve interpolation, Comput. Aided Geom. Des. 21 (1989) 363–370. [6] C-G. Lim, A universal parametrization in B-spline curve and surface interpolation, Comput. Aided Geom. Des. 16 (1999) 407–422. [7] H. Park, Choosing nodes and knots in closed B-spline curve interpolation to point data, Comput. Aided Des. 33 (2001) 967–974. [8] L. Piegl, W. Tiller, The NURBS Book, Monographs in Visual Communication, Springer, New York, 1995. [9] L. Piegl, W. Tiller, Circular approximation using integral B-splines, Comput. Aided Geom. Des. 35 (2003) 601–607.