Numerical method for

extracting an arc length

parameterization from
parametric curves
Richard J Sharpe and Richard W Thorne*

Let the arc of the quarter-circle be generated by the

A numerical technique is described which allows spacing o f
rational quadratic polynominal:
points as a function o f arc-length for any parametric
curve. The design o f a subroutine package for performing X(t) = [ R + X c + ( X c - R ) t=]/(1 + t =) (4)
this parameterization is outlined. The technique has been
incorporated into a CAD system and has been shown to be Y(t) = [Yc + 2Rt + Yc t 2 ] / ( 1 + ta) (5)
fast and accurate.
computer-aided design, numerical methods, parametric curves where 0 ~ < t ~ 1.
Figure 1 shows points on an arc generated with this
parametric form. An equal number of points separated by
an equal arc-length are also shown. It is evident that arc-
In computer-aided design and drafting applications it is length does not vary linearly with parameter value for this
frequently necessary to find points on a curve at various parametric representation. For many commonly used
arc-lengths from a reference point. This is true whether curve types in CAD it is either extremely difficult or
the system is a 2D drafting system which uses planar curve impossible to generate an arc-length parameterization
primitives or a 3D system based on surface primitives. For analytically lp . Another example which illustrates the
example, in the former case a draftsperson may wish to disparity between uniform parameter spacing and equal
place symbols representing rivets equally spaced along a arc-length distribution is the parametric cubic spline
curve. In the latter case, the creation of a ruled surface illustrated in Figure 2.
between two curves, with the rulings varying linearly with
arc-length, may be the desired result.

NEED FOR THE NUMERICAL METHOD

It is common practice in CAD to use parametric curve and
surface representatives. A parametric curve takes the follow-
ing form:
Fft) = (Xft), Y(t), Z(t)) (1)
A commonly used 2D parametric form for a circle is
given in equations (2) and (3).
x(t) = Xc + Rcos(t) (2)
Yft) = Yc + Rsin (t) (3)
where (Xc, Yc) is the centre of the circle, R is its radius and
t varies between 0.0 and 21r.
One advantage of the parametric form using these
trigonometric functions is that uniformly distributed
values of the parameter t yield directly points spaced
uniformly on the arc of the circle. For other parametric
forms this simple and convenient relationship is not
obtained as is illustrated in the following example of the
arc of a quarter-circle.

University of Connecticut School of Medicine, Farmington, Figure 1. The arc o f a quarter-circle: circles mark five points
CT 06032, USA spaced uniformly along the arc. Five crosses mark points
*Developmentand Engineering Department, Summagraphics located by uniformly distributed values of t in equations
Corporation, 35 Brentwood Avenue, Fairfield, CT 06430, USA (4) and (S).

volume 14 number 2 march 1982 0010-4485/82/020079-03 $03.00 © 1982 Butterworth & Co (Publishers) Ltd 79
 into the original parametric equation of the curve to find
the coordinates of the point.
X = X ('/finaL)
y = y ('tfinal)
Z = Z #final) (11)
A complication in applying the iteration equation is the
evaluation of the integral. This integral cannot, in general,
be solved analytically. A numerical integration technique
must therefore be applied. The authors have tried several
and have had the most success with the Romberg technique 3.
This technique required less time to converge than the
trapezoidal rule or Simpson's rule.

Figure 2. A parametric cubic spline, Eleven circles

corresponding to equally spaced arc-lengths are shown S O F T W A R E DESIGN
along with e/even crosses corresponding to equally spaced The authors have incorporated the algorithm described
parameter values. above into a CAD system. This was done in a manner such
that new parametric primitives could be accommodated by
the system with minimal software changes.
A hierarchical diagram of the subroutine package which
implements the algorithm is shown in Figure 3.
NUMERICAL METHOD
The entry point to the PT_AT_ARC_LENGTH function
The arc-length between two points on a parametric curve contains the following arguments:
is given by:
Input:
L = f.t (X' (0 2 + Y'(t) 2 + Z'(t) 2) ½ dt (6) CURVE_DATA a data structure which contains data
t ref defining the curve primitive and a
where flag indicating the curve type
(eg parametric cubic, B~zier etc)
tr~f = the value of the parameter corresponding to the
point of reference on the curve T_REF the parameter value of the reference
point on the curve
t = the value of the parameter corresponding to some T_ZERO the parameter value of an initial
general point on the curve guess point on the curve. This guess
In generating an arc-length parameterization the problem
then becomes one of finding the value of t for a given arc-
length L. Thus we must solve the non-linear algebraic
PT_AT_ARC_LENGTH
integral equation, shown above, for t.
Compute a point at
n given arc length
By rearranging (6) we obtain: from a reference point
on o curve
M(t) = f t (X'(t) 2 + Y'(t) ~ + Z'(t) 2) ½ dt - L = 0.0 (7)
t ref
I I
A very efficient technique for finding the value of the COMPUTE_ M_OF--T COMPUTE_M_PRIME_OF_T
parameter which satisfies M(t) is to apply the Newton- Compute M (t) as Compute M'(t) as
Raphson technique 3 . This non-linear equation solving tech- given by equation (7) given by equation (9)
nique converges quadratically, provided a good initial point
I (
is used to start the algorithm. ROMBERG-INTEGRATION COMPUTE_G_OF_T
The Newton-Raphson method generates a sequence of Numerically integrate o Compute G(t) as given
successively improved approximations to the solution of function using the Romberg by equation (10)
the equation M(t) = 0 using the relation: method
I
t.+ I = t - M(t.) M' Ct.) (8) t COMPUTE_DERWAT)VES
COMPUTE_G_OF_T Compute the parametric
Thus to solve the above equation we must first compute Compute G(t) as given derivatives of o curve
M'(t). In this case, M'(t) is: by equation (10)
½
M'(t) = (X'(t): + Y'(t) 2 + Z'(t)=) (9) I
COMPUTE_DERIVATIVES
By substituting (7) and (9) into (8), an iteration equation Compute the parametric
is obtained which will allow successively improved derivatives of o curve
approximation to the parameter t to be found.
tn+1: tn - [ f tn G(t) dt - L / /M' (t n) (10) Figure 3. Hierarchy of the subroutine package which
tref implements the point at an arc-length algorithm. To
½ accommodate new parametric curve types only the
where G(tn)= (X'(tn) z + Y'(tn) 2 + Z'(tn) z)
COMPUTE_DE RI VA TI VES subroutine needs to be
Once a final t is found, it can than be substituted back modified.

80 computer-aided design
 point is used to start the N e w t o n -
Raphson solution and to provide a
direction
ARC_LENGTH the arc-length from the reference
point at which the computed point
will be placed
Output:
T_FINAL the parameter value of the computed
point
POINT_FINAL (3) the Cartesian coordinates of the
computed point
Figure 4. Surface generated from parametric cubic spline
boundary curves. An arc-length parameterization was used
to generate the surface.
The calling program must provide an initial guess point
to start the Newton-Raphson iteration. This guess point
also indicates the direction along the curve to place the
computed point. The guess point in our case is computed hierarchy is similarly designed. It too contains the
using equation (12). COMPUTE_DERIVATIVES routine which is the only
subroutine in the leg which knows anything about the
to = tref + (Z~tl~L) • L (12)
different primitive types. Thus, to add a new curve
where primitive to the system only a single subroutine needs to
2n" = 0.1 be modified.
AL = the arc-length from tref to t = tref + 0.1 Performance studies have been carried out on the
L = the arc-length at which to place the point from tref algorithm and indicate that the average time necessary to
cbmpute a point is about 0.05 s on a Digital Equipment
The sign of L determines whether the computed point
Corporation VAX-11/780 for the worst case B~zier
will be placed at a value of t < tref Inegative L) or at a
curves and B-splines.
value of t > tref (positive L).
The basis for the above equation is that although para-
meter value is not exactly proportional to arc-length it is CONCLUSIONS
approximately proportional for most curve types. To
A numerical technique has been developed to extract an
simplify the equation At is set to the somewhat arbitrary arc-length parameterization from parametric curves. This
constant value of 0.1. Thus this linear relationship provides
algorithm makes it possible to compute points at a given
a simple method for computing the initial guess point for
arc-length along a curve when an analytical solution is
the algorithm.
either difficult or impossible to obtain.
Although the Newton--Raphson technique can diverge
The algorithm has been incorporated into a CAD system
in certain situations 4 it has been remarkably stable in this
and is used for placing points or symbols along curves and
particular application. In over six months of use in a CAD
also as the basis for surface generation. An example of the
system not one non-convergence problem has occurred.
latter is the surface shown in Figure 4. The lines shown for
Even when given poor guess points the algorithm has
display purposes were generated using the implicit arc-
converged in every case to a solution. The explanation for
length parameterization described in this paper.
this extremely good behaviour is probably due to the fact
that arc-length is a monotonically increasing function of
parameter value. As such the function does not have the ACKNOWLEDGEMENTS
local maxima, minima or inflection points which can
The authors are grateful to Dr Bertram Herzog for his
potentially cause problems for the Newton-Raphson
comments on the manuscript and to J une Younger for
technique 4 .
patiently typing it. We would also like to thank Peter
The hierarchy chart in Figure 3 provides an overview
Thomas for his technical assistance.
of the subroutine package used to compute a point at an
arc-length. Basically, the high level PT_AT_ARC_LENGTH
routine applies the iteration equation (10) until a REFERENCES
relative error test is satisfied.
The PT_AT_ARC_LENGTH function calls individual 1 Rogers, D F and Adams, J A Mathematical elements for
routines to compute M(t) and Mr(t) as given by equations computer graphics McGraw-Hill, New York, NY, USA
(7) and ~9) respectively. The subroutine which computes 11976)
M(t) must perform a numerical integration. It does this 2 Hartley, P J and Judd, C J 'Parameterization and shape
by calling ROMBERG_INTEGRATION which in turn of B-spline curves for CAD' Comput. Aided Des. Vol 12
calls COMPUTE_G_OF_T to numerically integrate the No 5 [September 1980) pp 235-238
function (equation (10)) using the Romberg technique. A
3 Burden, R L, Faires, J D and Reynolds, A L Numerical
relative error test is also used in the integration routine.
analysis Prindle, Weber and Schmidt, Boston, MA, USA
The only routine which performs any special casing in this
leg of the hierarchy is the routine at the bottom: 11978)
COMPUTE_DERIVATIVES. This routine computes the 4 Conte, S D and de Boor, C Elementary numerical
parametric derivatives (Xt(t), Y'(t), Z'(t)) and returns the analysis: an algorithmic approach McGraw-Hill, New
result to the calling routine. The second leg of the York, NY, USA (1972)

volume 14 number 2 march 1982 81