Constructing a cubic Bezier that passes through four points

1 February 2019

Drag any point $\lbrace P_0, P_1, P_2, P_3 \rbrace$ to refit the spline, or drag either tangent $\lbrace T_0, T_1 \rbrace$ to manually change the spline.

$\alpha =$

Bezier curves are extremely useful constructs that allow users to intuitively create smooth curves by adjusting end points and tangents. Mathematically, they are polynomials. In this article, we’ll focus on the most popular variant—the cubic Bezier curve.

Most applications have composite Bezier curves, i.e. curves constructed by chaining together multiple polynomials. The diagram above consists of three curves. For most of this article, we’ll focus on the middle curve one going from $P_1$ to $P_2$. This curve is defined by the four points $\lbrace P_1, T_1, T_2, P_2 \rbrace$, given by the parametric equation:

$\text{bezier}(t) = {P_1} \cdot (1-t)^3 + {T_1} \cdot 3 \cdot (1-t)^2 \cdot t + {T_2} \cdot 3 \cdot (1-t) \cdot t^2 + {P_2} \cdot t^3$

The other two curves are defined by their respective end points, and their tangents which are omitted from the diagram. I won’t go into too much Bezier math theory.
This is thoroughly covered in *A Primer on Bézier Curves* by Pomax.

It is very useful to edit the curve using its tangent points, and most spline tools support this. However, it is also very useful to just drag/create points and make the curve automatically fit them. Very few applications support this use case.

Bezier curves’ formulation doesn’t describe the curve in terms of four points it needs to pass through. They describe the curve in terms of two end points and two tangent points, as we saw above.

However there does exist a class of curves that does support this four-point formulation—Catmull-Rom curves.
Most of the math in this article is based on *Parameterization and Applications of Catmull-Rom Curves* by Cem Yuksel, Scott Schaefer, John Keyser.
There is a sub-class of these called Centripetal Catmull-Rom curves, that guarantees the absence of cusps (i.e. loops), which is exactly what we want.

Both—Bezier curves and Catmull-Rom curves—are cubic polynomials, and hence interconvertible.
*Splines Are Just Obfuscated Beziers* by Joshua Barczak.

**Therefore, by keeping things in Bezier form while also supporting Catmull-Rom-like behavior, we get the best of both worlds—fine-grained tangent control, and intuitive point-fitting.**

Given four points, $\lbrace P_0, P_1, P_2, P_3 \rbrace$, we want to find the tangent points $\lbrace T_1, T_2 \rbrace$. They are given by:

$\begin{aligned} T_1 & = \frac{d_1^2 P_2 - d_2^2 P_0 + (2d_1^2 + 3 d_1 d_2 + d_2^2)P_1} {3d_1(d_1 + d_2)} \cr T_2 & = \frac{d_3^2 P_1 - d_2^2 P_3 + (2d_3^2 + 3 d_3 d_2 + d_2^2)P_2} {3d_3(d_3 + d_2)} \cr \end{aligned}$

where $d_i = |P_i - P_{i - 1}|^\alpha$.

$\alpha$ is a constant, which dictates the “tightness” of the curve. When $\alpha = 0$, the parameterization is said to be uniform. When $\alpha = 1$, it is said to be chordal. When $\alpha = 0.5$, it is said to be centripetal, which is the one we’re interested in the most, since it guarantees lack of cusps.

You may have noticed that the above formulation requires four points $\lbrace P_0, P_1, P_2, P_3 \rbrace$, but only generates the curve between the middle two $\lbrace P_1, P_2 \rbrace$. So how did we generate the other (greyed out) curves in the figure?

We generate phantom points. *E.g.* the curve between $\lbrace P_0, P1 \rbrace$ is given by the phantom point $P_0 + \varepsilon (P_0 - P_1)$ and the points $\lbrace P_0, P_1, P_2 \rbrace$. $\varepsilon$ is any arbitrarily small number. In this example I’m using 0.0001; the main constraint is that the phantom point should be colinear to its adjacent real points, for best-looking results.

Here’s the pertinent Javascript code that powers the interactive figure up top. While the example supports 2D curves, the technique works with any number of dimensions.

```
function vec2_dist(a, b) {
return Math.sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
function refit_bezier(
passthru_0, passthru_1, passthru_2, passthru_3,
out_tangent_1, out_tangent_2
) {
let d1 = Math.pow(vec2_dist(passthru_1, passthru_0), alpha);
let d2 = Math.pow(vec2_dist(passthru_2, passthru_1), alpha);
let d3 = Math.pow(vec2_dist(passthru_3, passthru_2), alpha);
// Modify tangent 1
{
let a = d1 * d1;
let b = d2 * d2;
let c = (2 * d1 * d1) + (3 * d1 * d2) + (d2 * d2);
let d = 3 * d1 * (d1 + d2);
out_tangent_1.x = (a * passthru_2.x - b * passthru_0.x + c * passthru_1.x) / d;
out_tangent_1.y = (a * passthru_2.y - b * passthru_0.y + c * passthru_1.y) / d;
}
// Modify tangent 2
{
let a = d3 * d3;
let b = d2 * d2;
let c = (2 * d3 * d3) + (3 * d3 * d2) + (d2 * d2);
let d = 3 * d3 * (d3 + d2);
out_tangent_2.x = (a * passthru_1.x - b * passthru_3.x + c * passthru_2.x) / d;
out_tangent_2.y = (a * passthru_1.y - b * passthru_3.y + c * passthru_2.y) / d;
}
}
```