De Casteljau's algorithm explained

In the mathematical field of numerical analysis, De Casteljau's algorithm is a recursive method to evaluate polynomials in Bernstein form or Bézier curves, named after its inventor Paul de Casteljau. De Casteljau's algorithm can also be used to split a single Bézier curve into two Bézier curves at an arbitrary parameter value.

Although the algorithm is slower for most architectures when compared with the direct approach, it is more numerically stable.

Definition

A Bézier curve

B

(of degree

n

, with control points

\beta0,\ldots,\betan

) can be written in Bernstein form as followsB(t) = \sum_^\beta_b_(t),where

b

is a Bernstein basis polynomialb_(t) = (1-t)^t^i.The curve at point

t0

can be evaluated with the recurrence relation\begin\beta_i^ &:= \beta_i, && i=0,\ldots,n \\\beta_i^ &:= \beta_i^ (1-t_0) + \beta_^ t_0, && i = 0,\ldots,n-j,\ \ j= 1,\ldots,n\end

Then, the evaluation of

B

at point

t0

can be evaluated in \binom operations. The result

B(t0)

is given byB(t_0) = \beta_0^.

Moreover, the Bézier curve

B

can be split at point

t0

into two curves with respective control points:\begin&\beta_0^,\beta_0^,\ldots,\beta_0^ \\[1ex]&\beta_0^,\beta_1^,\ldots,\beta_n^\end

Geometric interpretation

The geometric interpretation of De Casteljau's algorithm is straightforward.

P0,...,Pn

. Connecting the consecutive points we create the control polygon of the curve.

t:(1-t)

and connect the points you get. This way you arrive at the new polygon having one fewer segment.

t

.The following picture shows this process for a cubic Bézier curve:

Note that the intermediate points that were constructed are in fact the control points for two new Bézier curves, both exactly coincident with the old one. This algorithm not only evaluates the curve at

t

, but splits the curve into two pieces at

t

, and provides the equations of the two sub-curves in Bézier form.

The interpretation given above is valid for a nonrational Bézier curve. To evaluate a rational Bézier curve in

Rn

, we may project the point into

Rn+1

; for example, a curve in three dimensions may have its control points

\{(xi,yi,zi)\}

and weights

\{wi\}

projected to the weighted control points

\{(wixi,wiyi,wizi,wi)\}

. The algorithm then proceeds as usual, interpolating in

R4

. The resulting four-dimensional points may be projected back into three-space with a perspective divide.

In general, operations on a rational curve (or surface) are equivalent to operations on a nonrational curve in a projective space. This representation as the "weighted control points" and weights is often convenient when evaluating rational curves.

Notation

When doing the calculation by hand it is useful to write down the coefficients in a triangle scheme as\begin\beta_0 & = \beta_0^ & & & \\ & & \beta_0^ & & \\\beta_1 & = \beta_1^ & & & \\ & & & \ddots & \\\vdots & & \vdots & & \beta_0^ \\ & & & & \\\beta_ & = \beta_^ & & & \\ & & \beta_^ & & \\\beta_n & = \beta_n^ & & & \\\endWhen choosing a point t0 to evaluate a Bernstein polynomial we can use the two diagonals of the triangle scheme to construct a division of the polynomialB(t) = \sum_^n \beta_i^ b_(t), \quad t \in [0,1]intoB_1(t) = \sum_^n \beta_0^ b_\left(\frac\right)\!, \quad t \in [0,t_0]andB_2(t) = \sum_^n \beta_i^ b_\left(\frac\right)\!, \quad t \in [t_0,1].

Bézier curve

When evaluating a Bézier curve of degree n in 3-dimensional space with n + 1 control points Pi\mathbf(t) = \sum_^ \mathbf_i b_(t),\ t \in [0,1]with\mathbf_i := \begin x_i \\ y_i \\ z_i \end,we split the Bézier curve into three separate equations\beginB_1(t) &= \sum_^ x_i b_(t), & t \in [0,1] \\[1ex]B_2(t) &= \sum_^ y_i b_(t), & t \in [0,1] \\[1ex]B_3(t) &= \sum_^ z_i b_(t), & t \in [0,1]\endwhich we evaluate individually using De Casteljau's algorithm.

Example

We want to evaluate the Bernstein polynomial of degree 2 with the Bernstein coefficients\begin\beta_0^ &= \beta_0 \\[1ex]\beta_1^ &= \beta_1 \\[1ex]\beta_2^ &= \beta_2\endat the point t0.

We start the recursion with\begin\beta_0^ &&=&& \beta_0^ (1-t_0) + \beta_1^t_0 &&=&& \beta_0(1-t_0) + \beta_1 t_0 \\[1ex]\beta_1^ &&=&& \beta_1^ (1-t_0) + \beta_2^t_0 &&=&& \beta_1(1-t_0) + \beta_2 t_0\endand with the second iteration the recursion stops with \begin\beta_0^ & = \beta_0^ (1-t_0) + \beta_1^ t_0 \\\ & = \beta_0(1-t_0) (1-t_0) + \beta_1 t_0 (1-t_0) + \beta_1(1-t_0)t_0 + \beta_2 t_0 t_0 \\\ & = \beta_0 (1-t_0)^2 + \beta_1 2t_0(1-t_0) + \beta_2 t_0^2\endwhich is the expected Bernstein polynomial of degree 2.

Implementations

Here are example implementations of De Casteljau's algorithm in various programming languages.

deCasteljau :: Double -> [(Double, Double)] -> (Double, Double)deCasteljau t [b] = bdeCasteljau t coefs = deCasteljau t reduced where reduced = zipWith (lerpP t) coefs (tail coefs) lerpP t (x0, y0) (x1, y1) = (lerp t x0 x1, lerp t y0 y1) lerp t a b = t * b + (1 - t) * a

def de_casteljau(t: float, coefs: list[float]) -> float: """De Casteljau's algorithm.""" beta = coefs.copy # values in this list are overridden n = len(beta) for j in range(1, n): for k in range(n - j): beta[k] = beta[k] * (1 - t) + beta[k + 1] * t return beta[0]

public double deCasteljau(double t, double[] coefficients)

The following function applies De Casteljau's algorithm to an array of, resolving the final midpoint with the additional properties and (for the midpoint's "in" and "out" tangents, respectively).

function deCasteljau(points, position = 0.5)

The following example calls this function with the points below, exactly halfway along the curve. The resulting coordinates should equal

(192,32)

, or the position of the centremost point.

See also

References

External links