Approximating the helix with rational cubic Bézier curves

Approximating the helix with rational cubic Bézier curves

computer-Aided Design, vol. 27, NO. 9, pp. 597-593, 1995 w 0 1995 Elswier Science Ltd EINEMANN printed in Great Britain.AHtights resewed 0010~4485(...

683KB Sizes 0 Downloads 155 Views

computer-Aided Design, vol. 27, NO. 9, pp. 597-593, 1995 w 0 1995 Elswier Science Ltd

EINEMANN

printed in Great Britain.AHtights resewed

0010~4485(p4)ooo12-3

0010-4485/95

510.00+0.00

Approxi mating the helix with rational cubic B6zier curves lmre Juh&z

The paper is on the approximation of a cylindrical helix by cubic rational Btzier curves. Two methods are given, one which guarantees exact slopes at the end-points of the approximating curve, and another which ensures 2nd-order continuity at the junction of two approximating curves. In both cases, there are error bounds which enable the helix to be approximated within any prescribed tolerance. Keywords: rational polynomials, B&ier curves, helices

this paper, we propose two rather elementary methods for the solution of the problem above, and we include error calculations. For the sake of simplicity, we consider a helix with a parameter p which is on a circular cylinder of radius r and axis z. This is not a restriction, because one can always achieve such a situation by applying an appropriate coordinate system transformation. The parametric form of the helix under consideration is h(v) = (r cos(u), r sin(v), pu)

Because of their flexibility, the most currently widespread curve and surface descriptions are rational curves and surfaces, especially NURBS.Increasing numbers of CAD systems and graphics libraries offer this method of curve and surface definition. That is why it is meaningful to deal with the approximation of known curves and surfaces by rational ones, even if the former can be described by exact mathematical formulae. The objective of this paper is to approximate a cylindrical helix by cubic rational BCzier curves. The helix is a unique curve, since its curvature and torsion are constant; the helix is the only spatial curve with this property. There are only two plane curves which share this property; the circle and the straight line. A practical consequence of this property is that the helix can be moved onto itself, i.e. by applying the appropriate helical motion to any arc of the curve, the moving arc remains on the original helix at any of its positions. This property can be utilized in the approximation of a helix, since, if we can approximate any arc of the helix, this approximating arc can be used repeatedly. Another important characteristic of our approximation is the symmetry of the approximating curves, which is also due to the constant curvature and torsion of the helix. Although the helix is a special curve with advantageous properties, it is an irrational curve. Note that a straight line can have an infinite number of points of intersection with a helix. Because of this, there is no exact rational representation of the helix, i.e. it is meaningful to investigate rational approximations. In an earlier paper’, a helix was approximated by cubic and quartic rational Bezier curves, following a kinematic approach, but without error calculations. In Department of Descriptive Geometry, Miskolc-Egyetemvlros, Hungary H3515

University

RATIONAL

BkZIER

form of rational

Bezier

B,“(t)w,b,

i =

(1)

CURVE

We will use the parametric curves, which is

r(t)

uE R

‘=,o

t E LO, 11

C B:(t)wi i=O

where the bi are control points, 0 < wi E [ware weights, and = (;)ti(l

B/(t)

- t)n-i

are the Bernstein polynomials of degree n, 0 I i I n. We are going to examine the cubic case, i.e. n = 3, where w,b,B;(t) r(t)

=

+ w,b,B:W

+ w,b,B;(t) +w,b&(t)

w,B;(t)

+ w,B;(t)

+ w,B;(t)

+ w,B;(t)

The first and second derivatives at the end-points will also be needed. These are

i(o) = 2(b,

- b,)

i(l) = ?(b,

-b,)

(2)

Computer-Aided Design Volume 27 Number 8 August 1995

587

of Miskolc,

Paper received: 22 March 1994. Revised: 23 September 1994

the helix with rational cubic B6zier curves: I Juhkz

Approximating

6(w,w, P(O) =

(b, -b,)

4 6(3w;

f(1)

- 3~;)

- w2w3)

=

t 2(b2

(b,-b,)--(b,-

2

w3

6w,

w3

- b,)

b ,>

\ (3) Figure 1 Control points of cubic ratlonal describes circular arc under consideration

Btzier

cuxve which

BASIC IDEA OF APPROXIMATION w,, =

We approximate an arc of the helix. The approximation process consists of the following steps.

Mi,=;(l+2cosu)

(1) Describe a circular arc in the (x, y> coordinate plane degree 2.

as a rational

_ ;(1+2cos kVl=

BCzier curve of

(3)

Using the notation w = wr = w2, we find that

Lift out the control points of the cubic rational BCzier curve from the (x,y> coordinate plane along a line parallel to the z axis.

h, = aI) b, =

The process described above results in a cubic rational BCzier curve on the cylinder of the helix to be approximated. A circular arc on the (x, y) coordinate plane with radius r, central angle ~CX(0 < (Y< n/2) and centre at the origin can be described as a rational BCzier curve2X3 of degree 2. Its control points are given by a,=(rcos

a,--rsin

a1

a,=(&-p) a,=(rcos

(r,rsin

(Y)

u’3 = 1

(2) Elevate the degree of the previous rational BCzier curve by 1.

1

b, =

a, + 2 cos cfa, 3w a2 + 2 cos ffa,

3w

b, = a2 The result is shown in Figure I. The next step is to lift out the control points from the (x, y) coordinate plane along a line parallel to the z axis. The first and last control points, i.e. b, and b,, are fitted onto the helix, that is

bo,=

(~1

-"P

b,, = 'YP

and its weights are given by

without loss of generality. Owing to this symmetry, b,, and b,, differ only in their sign. Thus we apply the notation b,; = -b

By elevating the degree4 of this rational BCzier curve, we produce the control points and weights of that cubic rational Bezier curve which describes the previous circular arc. Denoting the new control points by bj and the new weights by wi, we get

Now b is the only free parameter for the unique determination of the approximating twisted curve. The coordinates of the control points are finally

i

i

Wi=~i_,-+~~ 3

b2z = b

l-i

3I

b, = (r cos a, -r sin (Y, --up) r

b,= i

&,-I bi =

b,=

w,

1 + 2 cos (Y’ 2 + cos CY

-r

sin ff 1 + 2 cos (Y’ sin (Y

b r ,r 1 + 2 cos (Y’ 1 1 + 2 cos a!

b, = (r cos a, r sin (Y, ap)

where 588 Computer-Aided

2 + cos (Y

Design Volume 27 Number 8 August 1995

-b

Approximating the helix with rational cubic BQier curves: I Juhk

The coordinate functions rational BCzier curve are

of the approximating

cubic

r (t) = r (2t2 - 2t + 1) cos (Y+ 2t(l - t) * 2t(l -t> cos a + 2t2 - 2t + 1 2t-

r,(t) = r sin a ap(P

1

t> cos (Y+ 2t2 - 2t + 1

2t(l-

- (1 - tj3) + b(1 + 2 cos a)t(l - t) (2t - 1)

r,(t) =

2t(l-

Figure 2

t)cos (Y+ 2t2 - 2t + 1

Error surface in case of exact slopes

This results in the equalities

(4)

h,(u) = rJu)

EXACT SLOPES AT END-POINTS

P

b3z ((b3x

-

u2(1 - cos a) + (Y2(1 + cos a)

This yields a 2nd-degree polynomial in u. From its two roots we choose

The free parameter b can be determined in several ways. One option is to determine b in so as to ensure exact slopes at the end-points, i.e. to enforce the equality of the angle between the tangent line of the rational Btzier curve and the (x, y> coordinate plane, and the angle between the tangent line of the helix and the same coordinate plane at the first and last points of the arc. This means that the helix and the approximating curve have common tangent lines, but not necessarily common first derivatives, at the end-points. Using the previously introduced notation, this requirement implies the equality

-= r

2au

r sin v = r sin a

u=(

lasinra)

(6)

(7;;““)

by the preliminary conditions u E [-a, [-a, aI. The error function is given by 6(a,u)

=r,(u)

i

b2J2 + (b3y - b2y)2)L/2

and u E

-pu

u3(sin ff - cr cos a> +u(u~(~Lx+(Y cos a-sin

=P

b,,

a]

u2(a

-

a cos a) +

(Y)

a3(1 + cos (Y)

-V

1

cIE(o,T/2);UE[--(Y,LY]

(7)

from which

b=b,z=p

(

a-

2 sin a 1+2cos

a

1

since the error depends not only on u but on (Yas well (u has to be substituted for using Equation 6). One can see that the error is proportional to p, but independent of the radius. Moreover, a((~, - (Y)= 0, a((~, 0) = 0, S(cu, (Y)= 0 and

(5)

This choice of b results in a cubic rational BCzier curve with the following properties: . l

l

lili08(cy,u)

All of its points are on the cylinder of the helix. The points at the parameter values t = 0, t = 0.5 and t = 1 are on the helix. The tangent lines at the end-points coincide with the tangent lines of the helix.

=0

This error surface is shown in Figure 2. Because of symmetry, it is sufficient to examine the error on the domain (YE (0, n/2), v E [O,CX.].At any fixed value of a(cl! = G), we can calculate the u that gives the maximum error. In order to find this u, the equation

The deviation of the approximating curve from the helix is also of interest. A difference between the approximating rational Bezier curve and the helix can only exist in terms of their z coordinates. The arc of the helix is h(u), u E [-a, a] (see Equation 1) and the approximating curve is r(t), t E [0, l] (b is substituted in Equations 4 with Equation 5). We apply a domain transformation to the approximating curve by changing [O,11to [-a, a], i.e. by the function t(u) = (u + CX)/~CX, u E [ - a, a]. Then we have to find the matching u and u values, i.e. those which correspond to the same generator of the cylinder.

as(&,u) du

= 0

has to be solved, which results in cos ucos2 (Y+ 2cos2 cr + (Ycos (Ysin (Y- cos u - 2 + (Ysin (Y+ (Ycos u sin (Y= 0 2sin (Y-~(l+cos cos u =

0~)

LY- sin (Y

Computer-Aided Design Volume 27 Number 8 August 1995

589

Approximating the helix with rational cubic Bkier

curves: I Juh&sz

Now we show that 6( (Y,u) is monotone increasing in a on the specified domain. In order to prove this, it is sufficient to verify that the inequality

(1 - cos u)(cos (Y- cos u)((Y(l + 2cos (r)

-sin =P

(1 + cos a)(cos

20

(YE(O,rr/2);

(~(2 + cos (u))

(Y- 1)’ sin u

uE[O,al

Figure 4

holds. Since the denominator and the first factor of the numerator are positive, and the second factor of the numerator is nonpositive, it is enough to justify the inequality (~(1+2cos

cr)lsin

(Y(~+COS a)

(8)

By Maclaurin’s theorem, the inequalities

sin ff 2 ff - 5 6 ff2 cos (Y2 1 - 2 (Y2

(Y4

2

24

cosall--+-

are valid. Decreasing the right-hand side and increasing the left-hand side of Expression 8, we get

cr(lf2cos

a)53a-u”+~5sin

(w(2+cos

a)

i.e. Expression 8 is justified. A consequence of these properties of 6(a, u) is that the maximum error is monotone decreasing to zero as LXdecreases to 0. This property of the maximum error enables us to choose an cy which guarantees that the error does not exceed a prescribed tolerance. Assume that we are given any tolerance E > 0 E R. To find the maximum value of (Y which satisfies our requirements, we have to solve the equation 6(a, u) = E. This can be rearranged in the form 6( LX,u( (Y)) = E/P. In other words, the optimal value of (Y can be expressed as a function of E/P. Figure 3 shows a graph of this function.

Control

‘YP _=b

my b2,

from which b=L

5

i/--E/P

Figure 3

590

Optimal

w

0.01

0.02

0.03

values of a for exact slopes

Computer-Aided Design Volume 27 Number 8 August 1995

curves

In practice, it is quite common that the approximating curve has to be composed of several arcs. Considering our initial conditions, this is the case in which the helical arc to be approximated is greater than half of the pitch. In these cases, it is advisable to ensure a certain degree of continuity at junctions. In the method described above, there can only be lst-order (C’) continuity at junctions, because even the continuity of osculating planes is not guaranteed, which is a necessary condition for C2 continuity. Our first aim is to determine b in order to ensure common osculating planes at junctions. Let us consider the cubic rational BCzier curve r(t), r E [O,11, determined by the control points b,, b,, b,, b, and the weights w0 = 1, wi = W, w2 = w, w = 1, and the connecting cubic rational BCzier curve kt), t E [l, 21, b, = b,, b,, b2, b3, Go = w,,, 6, = wl, 6,, = w2, G3 = wj (see Figure 4). The curve f(t) is obtained from r(t) by a helical motion. This helical motion is determined by CX, p, and its axis is the t axis, i.e. the axis of the helix to be approximated. These two curves have a common osculating plane at the junc$on_point b, = b, if the control points b,, b,, b,, b,, b,, b, are coplanar. Owing to the symmetry, fhese control points are coplanar if the lines b,, b, and b,, b, are intersecting, and the distance from their point of intersection m to the (x, y) coordinate plane is equal to the distance from the control point b, = b, to the same coordinate plane. Using the notation of Figure 5,

a

0 , 0

Btzier

C 2 CONTINUITY AT JUNCTIONS

b2,

‘YP mY

*

points of two connecting

Figure 5

Approximatingthe helix with rational cubic B6zier curves: I Juhbz

By the equality

curvature of r(t) at t = 1 by c, 2r+r

my=bzx

tana=

cos ff sin LY

2w,w, l(b3 -b,) K=w

1+ 2 cos (Y cos ff

Ire,- b,13

we get the expression b = ap

_

K=T

cm a

2w,w,

I& - a,) X (h2 - 8,)l

b1- &I3

2 + cos (Y

The resulting rational BCzier curves have a common osculating plane at the junction, i.e. at the parameter value t = 1. Let us examine whether there is a C2 continuity. The condition for C’ continuity is

~(b,_b2)=3w, ‘b _i wo(l 0)

K = 2 since w. = wg = 1 and w2 = wr = w. This means that these curves not only have common osculating planes but also equal curvatures at the junction. Summarizing our results, we can see that the approximating cubic rational BCzier curve has the following properties: l l

(see Equations 2). This equality holds_because w. = w3 = 1, w2 = w1 = w, and b, - b, = b, - b,. The further condition for C2 continuity is the equality k2 - 2 cos c&

-a,)

= b, + 2 cos db,

- b,)

which can be derived from Equations 3. Owing to symmetry (see Figure 6), this equality holds if the equalities 6,, - 2 cos LX(6,, - &%,J= 6,, b,, + 2 cos db,,

- b2J = b,,

are valid. This results in cos a = 1, i.e. there can be C2 continuity only in the case in which (Y= 0. Neither can C2 continuity be gained by the application of an affine parameter transformation, since, for C1 continuity, the equality of Equations 2 requires equal-length parameter ranges under the previous circumstances. Denoting the curvature of r(t) at t = 1 by K, and the

x (b, - b,)l

.

All of its points are on the cylinder of the helix. The points at the parameter values t = 0, t = 0.5 and t = 1 are on the helix. There is a 2nd-order (G2> geometric continuity at the junction.

In order to achieve C2 continuity, we carry out a projective parameter transformation3p5 which, as is already known, results in the transformation of weights. Thereby, the rational BCzier curve determined by the control points bi and weights w,(i = 0, 1,. . . , n> and the curve determined by the control points bi and weights c’w,(i = 0, 1,. . . , n), 0 < c E R, are two different representations of the same curve. Different points belong to the same parameter value in different representations. We apply a projective parameter transformation for the cubic rational BCzier curves of Figure 4. Consider the cubic rational BCzier curve r(t), t E [O,11, determined by the control points b,, b,, b,, b, and the weights w. = 1/c3, w1 = w/c2, w2 = w/c, w3 = 1, and the curve Xt), t E [l, 21, of the same type determined by the control points cb,, hr, cb2,b3 and the weights Go = 1, 6i = Ew, G2 = Z2w, w3 = E3. We examine whether there are positive constants c and E which guarantee C2 continuity of the curves I(t) and if;(t) at the parameter value t = 1. For Co continuity, b, = ho must hold. For C’ continuity,

This results in l/c = E, since b, - b, = b1 - ho. Making use of this, we get the weights w. = c3, w1 = c2w, w =cw, w3=1, and i+),=l, $i=cw, i?,,=c2w, G3= C?

For C2 continuity, 3Wz - W2W3

4

(b, -b,)

- $(b,

-b,)

Figure 6 Control points of two cubic rational B6zier curves with common osculating plane at junction

Computer-Aided Design Volume 27 Number 8 August 1995

591

Approximating the helix with rational cubic B6zier curves: I Juh6sz

After some simplification, we get 1

@JZ - h,l 2b,-b,(

-=3wC

Considering the similar triangles of Figure 7, we find that 1

cc-

cos a By this choice of c, a C* continuity is ensured at the junction. The error calculations are analogous to that of the previous section. We apply the same domain transformation, and we get the error function txa,u)=p

Au” + Bu @$D --u

where A = 1 - cos* LY B=a2(3+2cos a+cos2

Figure 8

Surface dNcr, u)/du,

where A, B, C and D are as defined above, and E = a sin (~/cl - cos a). The root should be in the range (cos (Y,l), since v E (0, (Y) for the nontrivial solution. In Figure 8, there is a plot of the surface &?(cx, v)/av, (YE (0,7r/2), v E [O,a]. The surface is illustrated by the perpendicular projection of its lines (Y= constant on the coordinate plane (u,(&Y(a,v)/dv)). The error function increases monotonically in (Y on the specified domain, since the inequality

a)

asca,v1

c = 2 - cos a - cos* (Y D=(r*(2+3cos CX+COS~(u) and u is determined by Equation 6. This error function has the same properties as the function in Equation 7. Owing to the symmetry, it is enough to examine the error on the domain (YE (0,7r/2), LJE [O, ff 1. At any fixed value of (Y(CX = &), we can compute the u that gives the maximum error. For this, we have to find a roof of &%&, v)/a v = 0. This means finding a root of

da

(1 -cos

u)

= (cos* IY- l)(cos (Y+ 2)* sin v

( cy(5 cos2 a + 2 cos u cos* (Y+ 2 cos (Ycos v -cosv+8cos -sin tcos”

a+2)

(~(~0s acosv+2cosv+10cos

ff

(Yf4))20

can be proved. To show this, we must justify the inequality. a(5cos2

cos3 v(E4C2 - 2E*CD + D2)

(Y+2COsvCOs2

-cosv+8cos

+ cos2 v( -E4C2 - 2E*CD + 3D2 - E’AC

5 -sin

+E3(3AD - BC) - EBD) + cos v( -E4C2 + 2E*CD + 3D2 + 2E5AC - 2EBD) + E4C2 + 2E*CD + D* - E5AC -E3(3AD - BC) - EBD) = 0

CIE (0, n-/2), v E [O,a]

a+2cos

ctcosv

a+2)

cw(cos crcosv+2cosu+10cos

a

+cos3 a!+4) This can be proved in a similar way to Expression 8, i.e. one can use the inequalities derived from Maclaurin’s theorem (but with more terms) and the inequality cos v r cos CY.A consequence of these properties of 6( (Y,v) is that the maximum error decreases monotonically to zero as CYdecreases to 0. Making use of this property of the maximum error, the optimal (Y can be determined for any prescribed error limit, as in the previous section. However, in this case, the evaluation of a(&, v), i.e. the computation of the maximum value of 6 for a fixed (Y= &‘, is more complicated. The optimal value of (Ycan be considered as a function of E/P, i.e. Na, v(a)) = E/P. Figure 9 shows a graph of this function.

CONCLUDING REMARKS

Figure 7

592

In the previous section, we have shown that two approximating segments can be joined to form a C* curve. Unfortunately, it is impossible to join more than

Computer-Aided Design Volume 27 Number 8 August 1995

Approximating the helix with rational cubic Bdzier curves: I Juh&z

0

e/P W

, 0

Figure 9

0.01

0.02

0.03

Optimal values of a for C2 continuity

Q lr

z *-

cc--

___-___--

degree of continuity at junctions decreases the accuracy of the approximation, i.e. for the same accuracy, more approximating arcs have to be applied. Since the error of the approximation is proportional to p, the increase of p decreases the optimal (Y (see Figures 3 and 9), i.e. for a very large p, (Y is too small to be applicable in practice. However, in practical applications such as screw threads, p is small. By means of the method described above, helicoids can be approximated by rational BCzier patches. If the helicoid is generated by a helical motion with parameter p of a rational BCzier curve, then the vertical error of this approximation is equal to the error of the approximation of a helix of the same parameter. Figure II shows an axonometric view of a helicoid approximated by rational BCzier patches.

I

ie/P ”

0

Fire

ACKNOWLEDGMENTS

I

0:01

0:02

-

exact slopes

-__

Cz or Gs continuity

0:03

-

10 Optimal values of a

The author would like to thank Z S Bancsik for discussions with him, and the referees of this paper for their helpful comments and suggestions.

REFERENCES Mick, S and Riischel, 0 ‘Interpolation of helical patches by kinematic rational Bezier patches’ Compui. & Graph. Vol 14 No 2 (1990) pp 275-280 B&m, W, Farin, G and Kahmann, J ‘A survey of curve and surface methods in CAGD’ Comput. Aided Geom. Des. Vol 1 No 1 (1984) pp l-60 Farin, G Curves and Surfaces for Computer Aided Geometric Design Academic Press (1988) Farin, G ‘Algorithms for rational B&zier curves’ Compur.-Aided Des. Vol 15 No 2 (1983) pp 73-77 Patterson, R R ‘Projective transformation of the parameter of a Bernstein-BCzier curve’ ACM Trans. Graph. Vol 4 No 4 (1985) pp 276-290

Figure 11

Axonometric view of helicoid approximated with rational Bezier patches

two segements to form a C2 curve, but any number of segments (with the same (Y)can be joined to form a G2 curve. It is obvious from Figure 10 that the increase in the

Imre Juhdsz is a lecturer at the Department of Descriptioe Geometry at the Uniuersity of Miskolc in Hungary. He graduated from Kossuth Lajos University of Sciences and Arts in Debrecen, from where he also obtained his doctorate. His research interests are in geomem’c moa’eiling, especially the description of curves and surfaces, algorithms for computer graphics, and graphics standards.

Computer-Aided Design Volume 27 Number 8 August 1995

593