ndesmic

Posted on

# Splines from Scratch: Bézier Curves

The Catmull-rom spline I found one of the easiest to interpret and it has a lot of uses for creating smooth paths. However, the Bézier Curve is perhaps one of the most used in actual practice.

# Intuition

Unlike Catmull Rom splines which have 2 control points per segment, Bézier curves can have a number of control points, from two to infinity. In most cases you'll hear something like "Cubic Bézier" which notes how many control points we have. Usually we don't have too many because that gets a bit unwieldy to deal with. If you were to draw a line from point to point you'd wind up with a polygon (a "hull" if you will) and the algorithm is basically trying to fit the curve inside that polygon. It's also possible to have non-convex polygons and these tend to yield interesting results like loops but are harder to wrap your head around. The cubic bézier is probably the most common because it has 4 points and thus the forms a quadrilateral. This is "cubic" because the exponent of the function is going to be equal to the number of points minus 1.

# The confusing way

If you look around for the equation for a Bézier curve you're likely to find a lot of math jargon.

The big sigma `Σ` always causes my brain to shut-down even if I know exactly what it means, and I have to really focus to parse it. It's just a summation, from the bottom value to the top so you can think of it like a for loop:

``````let result = 0;
for(let i = 0; i < n; i++){
result += /*rest of equation*/
}
``````

The next piece is the combinatorial notation "N chooses K". This is the number of possible combinations you could get from reaching into a bag of N things and picking K of them. It's usually written like:

With a vertical set of numbers in parens (the middle representation). And the answer is `n! / k!(n-k)!`. Knowing this we can convert it into something more usable:

``````function factorial(n){
if(n < 0) throw `Cannot get factorial of negative number \${n}`;
if(n === 0) return 1;
let result = n;
for(let i = 1; i < n; i++){
result *= n - i;
}
return result;
}

function combinations(n, k){
return factorial(n) / (factorial(k) * factorial(n - k));
}
``````

Two notes here: this will get really large really quick. Not a problem since we won't be going into 20 control points but if you try it you might have issues. Second is that this is a very slow implementation (only a smidge faster than a recursive factorial) meant to be easy to read and understand. There's a lot you can do with memoization or look-up tables to speed this up.

The rest of the equation is relatively straight-forward if a bit like reverse engineering minified code `((1 - t) ** (n - i)) * (t ** i) * P[i]`. `P` in this case is a vector (or in javascript terms an array) with the component indexed by `i`. But actually there's a little more. Because we're doing things in 2D, each component of `P` is a point with x,y and potentially more dimensions so we actually do this whole thing once per dimension.

Putting this into code:

``````function bezierSpline(points, t){
const dimensions = points[0].length;
const n = points.length - 1;
const result = new Array(dimensions);
for(let d = 0; d < dimensions; d++){
let dimensionResult = 0;
for(let i = 0; i < points.length; i++){
dimensionResult += combinations(n,i) * ((1 - t) ** (n - i)) * (t ** i) * points[i][d];
}
result[d] = dimensionResult;
}
return result;
}
``````

It took a while to get the loop right due to off by one. We need to iterate over 4 points starting at 0, but the `n` term is the "order" of a polynomial which is 1 less than the number of points.

If we plug in value 0 and 1 for `t` you can see that the other terms fall out and we're just left with the starting and ending point which is what we expect.

The result:

If the points aren't convex, that is if say point 2 is further right than point 3 so the edges cross each other you can get some stranger curves:

# The more intuitive way

Even if I can translate the equations, it's really hard to make sense of what all of that means. Why are we using combinatorics? Well, luckily there's a second way to go about this that doesn't involve exponents, factorials or any of that and it also much easier to visualize. It's called De Casteljau's algorithm. The basic algorithm is this:

1) Take the control points and form lines (you will get points.length - 1) lines.
2) Linearly interpolate each line by `t` to get a value
3) If you have more than 1 point, do #1 again
4) The final point is the point along the curve at `t`

It turns out to be easier to implement:

``````function lerp(pointA, pointB, normalValue) {
const result = [];
for(let d = 0; d < pointA.length; d++){
result.push(pointA[d] + (pointB[d] - pointA[d]) * normalValue)
}
return result;
}

function bezierSpline2(points, t){
const segs = pointsToSegments(points);
const values = segs.map(seg => lerp(seg[0], seg[1], t));
if(values.length === 1) return values[0];
return bezierSpline2(values, t);
}
``````

It will produce the exact same result and I think it's faster too because it doesn't have crazy factorials to compute.

# Rational Bézier curves

So before we throw the first algorithm out, we can still use it to build a new, more powerful curve drawing function. I'm not sure if there's an easy way to implement it with De Casteljau's algorithm so perhaps there is value yet. This isn't too hard to implement either. What we will do is add an extra term for the point's weight. I think it makes most sense to encode that as the final dimension in the point. Then we just multiply it:

``````const dimensions = points[0].length - 1;
//...
const result = combinations(n, i) * ((1 - t) ** (n - i)) * (t ** i) * points[i][dimensions] * points[i][d];
``````

So now every point is contributing a weight but we actually need to divide the entire thing by the weighted value without the point to get a correct result. That is:

``````const weight = combinations(n, i) * ((1 - t) ** (n - i)) * (t ** i) * points[i][dimensions];
``````

We aggregate both the results and the weights and then divide. The final code looks like this:

``````function rationalBezierSpline(points, t) {
const dimensions = points[0].length - 1;
const n = points.length - 1;
const result = new Array(dimensions);
for (let d = 0; d < dimensions; d++) {
let dimensionResult = 0;
let dimensionWeight = 0;
for (let i = 0; i < points.length; i++) {
const weight = combinations(n, i) * ((1 - t) ** (n - i)) * (t ** i) * points[i][dimensions];
dimensionResult += weight * points[i][d];
dimensionWeight += weight;
}
result[d] = dimensionResult / dimensionWeight;
}
return result;
}
``````

Since we're reserving the last value as the weight the dimensions are now the tuple length minus 1. But in practice if you are only ever drawing 2d, you can simply assume for length 2 tuples the weight is always 1 and you can assume dimensions is always 2.

So what does this do? It can give use the ability to add more pull to a point. Here's the grumpy face spline from earlier:

And here it is with the 2nd point's weight changed to `3`:

And with the weight changed to `10`:

This gives us even more ability to describe the curve. Since it's the ratio of weights, we call it a "rational bézier."

# A slightly improved drag and drop

The previous article on Catmull-Rom splines explained how to setup a component with nice drag-and-drop for control points to make it interactive. Since then I discovered a new, handy API `setPointerCapture`.

``````#selectedIndex;
#lastPointer;
#currentOffset;
#selectedPoint;

onPointerDown(e) {
const rect = this.dom.canvas.getBoundingClientRect();
this.#lastPointer = [
e.offsetX,
rect.height - e.offsetY
];
this.#selectedIndex = this.#points.findIndex(p => inRadius(this.#lastPointer[0], this.#lastPointer[1], 4, p[0], p[1]));

if (this.#selectedIndex > -1) {
this.#selectedPoint = [
this.#points[this.#selectedIndex][0],
this.#points[this.#selectedIndex][1]
];
this.dom.canvas.setPointerCapture(e.pointerId);
}
}
onPointerMove(e) {
const rect = this.dom.canvas.getBoundingClientRect();
const currentPointer = [
e.offsetX,
rect.height - e.offsetY
];
this.#currentOffset = [
currentPointer[0] - this.#lastPointer[0],
this.#lastPointer[1] - currentPointer[1]
];

this.#points[this.#selectedIndex] = [
this.#selectedPoint[0] + this.#currentOffset[0],
this.#selectedPoint[1] - this.#currentOffset[1]
];

this.renderSpline();
}
onPointerUp(e) {
this.#points[this.#selectedIndex] = [
this.#selectedPoint[0] + this.#currentOffset[0],
this.#selectedPoint[1] - this.#currentOffset[1]
];
this.renderSpline();
this.dom.canvas.removeEventListener("pointermove", this.onPointerMove);
this.dom.canvas.removeEventListener("pointerup", this.onPointerUp);
this.dom.canvas.releasePointerCapture(e.pointerId);
}
``````

This works the same as last time. On `pointerdown` we get the the closest point (if there is one) and set a `pointermove` and `pointerup` event. I've also converted the points into tuples rather than X/Y variables because it was cleaner. What you will notice is the the `this.dom.canvas.setPointerCapture` call. This lets the canvas element continue receiving pointer events even when the pointer is not over the canvas. This let's us have a nice exclusive access to the mouse during manipulation. In `pointermove` we're taking the offset of the new pointer position since we started, adding that to the original point value, and then redrawing. On `pointerup` we're committing the change, unregistering the move and up events and releasing the pointer lock.

Now even if you start dragging beyond the canvas everything should still work as expected and we don't need to set event listeners on the document to do it.

Resources: