Computing the length of a Bézier path

Bill Zaumen

For printing, see the PDF version.

Introduction

Bézier paths, as used in practice, consist of multiple segments, each of which is either a straight-line segment, a quadratic Bézier curve segment, or a cubic Bézier curve segment. These curves can be represented by the following expressions:

• $\mathbf{p}(t) = \mathbf{p}_0B_{0,1}(t) + \mathbf{p}_1B_{1,1}(t)$ for linear segments,
• $\mathbf{p}(t) = \mathbf{p}_0B_{0,2}(t) + \mathbf{p}_1B_{1,2}(t) + \mathbf{p}_2B_{2,2}(t)$ for quadratic segments,
• $\mathbf{p}(t) = \mathbf{p}_0B_{0,3}(t) + \mathbf{p}_1B_{1,3}(t) + \mathbf{p}_2B_{2,3}(t) + \mathbf{p}_3B_{3,3}(t)$ for cubic segments,

If the X and Y components of $\mathbf{p}(t)$ are $x(t)$ and $y(t)$ (with an obvious generalization to more dimensions), then the distance between two points on a segment is $$ d(t_1,t_2) = \int_{t=t1}^{t2}\sqrt{A(t)^2 + B(t)^2}\,dt$$ where $$ \begin{eqnarray} A(t) & = & \frac{dx}{dt} \\ B(t) & = & \frac{dy}{dt} \end{eqnarray} $$ These derivatives can be computed by using the relation, $\frac{d}{dt}B_{i,n}(t) = n(B_{i-1,n-1}(t) - B_{i,n-1}(t))$ with $B_{-1,n-1}(t) = B_{n,n-1}(t) = 0$. This relation can be used to differentiate $x$ and $y$, representing the derivatives as polynomials using a Bézier basis and the derivatives can then be easily transformed to a monomial basis.

Because $A(t)$ and $B(t)$ are squared, the square root will always be positive unless $A(t)$ and $B(t)$ share a common root. In addition, $A(t)$ and $B(t)$ can be at most quadratic polynomials. When $A(t)$ and $B(t)$ share at least one root, we can write the integral as $$ d(t_1,t_2) = \int_{t=t1}^{t2}|P(t)|\sqrt{Q(t)}\,dt$$ where P(t) and Q(t) are polynomials and where the degree of $Q$ is either 2 or 0.

Straight line segments

For linear line segments, $A(t)$ and $B(t)$ are constants and the integral is trivial:the length of a line segment is a constant multiplied by $t_2 - t_1$, with positive values indicating that the path parameter is increasing and negative values indicating that the path parameter is decreasing.

Quadratic line segments

For quadratic line segments, $A(t)$ and $B(t)$ are constants or first degree polynomials. If $A(t)$ and $B(t)$ have a common root, the integrand is the absolute value of a polynomial. Otherwise, the integral can be written as $$ d(t_1, t_2) = \int_{t=t_1}^{t_2} \sqrt{a + bt + ct^2}\,dt$$ That integral can be found in the CRC Standard Mathematical Tables, which states that it is $$ d(t_1,t_2) = \left. \frac{(2ct+b)\sqrt{X}}{4c} \right|_{t=t_1}^{t_2} + \frac{1}{2k}\int_{t=t_1}^{t_2}\frac{1}{ \sqrt{X}}\,dt$$ where $q = 4ac-b^2$, $k = \frac{4c}{q}$ and $X = a + bt + ct^2$ In addition, $\int_{t=t_1}^{t_2}\frac{1}{\sqrt{X}}\,dt$ has different values for

• $c > 0$, in which case $$\int \frac{dt}{\sqrt{X}} = \frac1{\sqrt{c}} \log (\sqrt{X} + t\sqrt{c} + \frac{b}{2/sqrt{c}}) $$
• $c > 0$ and $q > 0$, in which case $$\int \frac{dt}{\sqrt{X}} = \frac1{\sqrt{c}} \sinh^{-1} \left(\frac{2ct+b}{\sqrt{4ac-b^2}}\right) $$
• $c < 0$, in which case $$\int \frac{dt}{\sqrt{X}} = \frac1{\sqrt{-c}} \sin^{-1} \left(\frac{-2ct-b}{\sqrt{b^2 - 4ac}}\right) $$

For integrals containing absolute values of polynomials, it is necessary to break the integral up into segments whose boundaries are the limits of integration and the roots, with no segment including a root internally. Finally, when the integral uses logarithms, it is possible for $a+bt+ct^2$ to be positive for $t\in[t_1,t_2]$ and the argument for the logarithm to be negative. For this case use $\log x = \log (-x) + i\pi$ for $x < 0$. While the integral is a complex number, the imaginary part is a constant and thus drops out of definite integrals.

Cubic line segments

For cubic line segments, $A(t)$ and $B(t)$ could be constants or first degree polynomials, but are typically second degree polynomials If $A(t)$ and $B(t)$ have a common root, the integrand is the square root of a constant or a second degree polynomial $Q(t)$, multiplied by the absolute value of a polynomial $P(t)$. There are three cases:

1. if $A(t)$ and $B(t)$ share one common root, $P(t)$ is a first degree polynomial and $Q(t)$ is a second degree polynomial
2. if $A(t)$ and $B(t)$ share two roots or both have the same double root, then $P(t)$ is a second-degree polynomial and $Q(t)$ is a constant.
3. If there are no roots in common, $Q(t) = A(t)^2 + B(t)^2$ is a fourth degree polynomial with no real roots, and can be factored into two polynomials $Q_1(t)$ and $(Q_2(t)$ that have no real roots and that have positive values for all values of $t$.

The third case, which requires integrating the square root of a quartic polynomial, uses elliptic integrals. If $Q_1(t)$ and $Q_2(t)$ are written as $$\begin{eqnarray} Q_1(t) & = & f_1 + g_1t+h_1t^2 \\ Q_2(t) & = & f_2 + g_2t + h_2t^2 \ , \end{eqnarray}$$ the integral can be found in B.C. Carlson, " A Table of Elliptic Integrals: Two Quadratic Factors," Mathematics of Computation, Volume 59, Number 199, July 1992, Pages 165—180. Equation 2.45 on Page 169, labeled as [1, 1, 1, 1], contains this integral, and makes use of a number of variables Carlson defined before Equation 2.45. One obtains $$\begin{eqnarray} \int_{t=y}^{x}\sqrt{Q_1(t)Q_2(t)}\,dt & = & (\delta^2_{22}/h_2^2 - \delta_{11}^2/h_1^2)[\psi_0H_0 + (\Delta_0 - \delta_{12}^2)R_f]/8 \\ & - & (3\psi_0^2-4h_1h_2\delta_{12}^2)(\Sigma+\delta_{12}^2R_f) / (24h_1^2h_2^2) \\ & + & [\Delta^2R_f-\psi_0A(1,1,1,1)]/(12h_1h_2) + E/(3h_1) \end{eqnarray}$$ where the quantities in this integral are defined below. In the following list of definitions, $i$ has each of two values (1 and 2) and the functions $R_F$, $R_J$, $R_C$, and $R_D$ are the Carlson symmetric forms of elliptic integrals: $$\begin{eqnarray} H_0 & = & \delta_{11}^2\psi_0 [ R_J(M^2,L_{-}^2,L_{+}^2,\Omega_0^2)/3 +R_C(a_0^2,b_0^2)/2]/h_1^2 -X_0R_C(T_0^2,V_0^2) \\ a_0^2 & = & b_0^2 + \Lambda_0(\Lambda_{+} - \Lambda_0) (\Lambda_0 - \Lambda_{-}) \\ b_0^2 & = & (S^2/U^2 + \Lambda_0)\Omega_0^4 \\ V_0^2 & = & \mu_0^2(S^2 + \Lambda_0U^2) \\ T_0 & = & \mu_0S + 2h_1h_2 \\ \mu_0 & = & h_1/(\xi_1\eta_1) \\ X_0 & = & -(\xi_1^\prime\xi_2 + \eta_1^\prime\eta_2)/(x-y) \\ \Omega_0^2 & = & M^2 + \Lambda_0 \\ \Lambda_0 & = & \delta_{11}^2h_2/h_1 \\ \psi_0 & = & g_1h_2 - g_2h_1 \\ A(1,1,1,1) & = & \xi_1\xi_2 - \eta_1\eta_2 \\ S & = & (\xi_1\eta_1\theta_2 + \xi_2\eta_2\theta_1)/(x-y)^2 \\ \Sigma & = & G-\Delta_{+}R_f + B \\ R_f & = & R_F(M^2,L_{-}^2, L_{+}^2) \\ G& = & 2\Delta\Delta_{+}R_D(M^2,L_{-}^2,L_{+}^2)/3 + \Delta/(2U) + (\delta_{12}^2\theta_1 - \delta_{11}^2\theta_2) /(4\xi_1\eta_1U) \\ L_{\pm}^2 & = & M^2 + \Delta_{\pm} \\ \Delta_{\pm} & = & \delta_{12}^2 \pm \Delta \\ \Delta &=& (\delta_{12}^4 -\delta_{11}^2\delta_{22}^2)^{1/2} \\ \delta_{ij} & = & (2f_ih_j + 2f_jh_i - g_ig_j)^{(1/2)} \\ M & = & \zeta_1\zeta_2/(x - y) \\ U & = & (\xi_1\eta_2+\eta_1\xi_2)/(x - y) \\ \zeta_i & = & [(\xi_i+\eta_i)^2 - h_i(x-y)^2]^{(1/2)} \\ \theta_i & = & \xi_i^2 + \eta_i^2 - h_i(x-y)^2 \\ E & = & \xi_1^\prime\xi_1^2\xi_2-\eta_1^\prime\eta_1^2\eta_2 \\ B & = & \xi_1^\prime\xi_2 - \eta_1^\prime\eta_2 \\ \eta_1^\prime & = & (g_1 + 2h_1y)/(2\eta_1) \\ \xi_1^\prime & = & (g_1 + 2h_1x)/(2\xi_1) \\ \eta_i& = & (f_i + g_iy + h_iy^2)^{1/2} \\ \xi_i & = & (f_i + g_ix + h_ix^2)^{1/2} \end{eqnarray}$$ Programs that define these variables should place them in the reverse order: the list uses the convention in which quantities are used before they are defined, the opposite of that used by Carlson.

An algorithm for factoring $Q(t)$ is described in Brookfield, " Factoring Quartic Polynomials: A Lost Art", Mathematics Magazine, Vol. 80, No. 1, February 2007, Pages 67–70. The easiest way to proceed is to factor $Q(t)$ into three terms: a constant, and two quadratic polynomials whose $t^2$ coefficients are 1.0. There are two special cases.

1. the minimum value of at least one of the two factors is small (e.g., less than $1 \times 10^{-4}$). In this case, tests indicate that numerical accuracy is poor, and one should use some other method such as numeric integration.
2. the two factors are identical. In this case, the integrand is the absolute value of a quadratic polynomial. For this case, the integral given in Carlson's paper can fail: a Java implementation returned Double.NaN, most likely because of division by zero, which suggests that Carlson's integral may not be numerically accurate when the two factors are nearly, but not exactly identical.