COMPUTING THE BOUNDS FOR A VECTOR SHAPE CONSISTING OF LINES, QUADRATIC AND CUBIC BEZIER CURVES

April 2010 - Michiel Kamermans.

NOTE: this page has been superceded by this Primer on Bezier Curves. That said, I'll keep the original article up for historical purposes.
Ever wanted to compute the real bounds of a vector shape? Okay, rephrase: ever needed the real bounds for a vector shape? K, let's do some math. The obvious way to get the genuine bounds is to look at a vector shape as a series of line segments:
vector = (line*, curve*)+
(where * and + have the usual meaning of '0 or more' and '1 or more')
The idea of bounds computation is simply to do this:
keep shape min x, min y, max x and max y unassigned (in case we have no segments)
for each segment in some vector shape do:
	if segment's min x < shape min x, set shape min x to segment's min x
	if segment's min y < shape min y, set shape min y to segment's min y
	if segment's max x > shape max x, set shape max x to segment's max x
	if segment's max y > shape max y, set shape max y to segment's max y
After doing this, either our boundary variables are empty (in which case the shape was just one or more disjoint points), or they now represent the bounding box for our shape. So far so simple. But how do we compute the min x/y and max x/y for segments?

linear bounds are simple

FINDING BOUNDING BOXES FOR LINE SEGMENTS

Lines are pretty simple. They run from x1/y1 to x2/y2, and their bounding box is the almost trivial
{min(x1,x2), min(y1,y2), max(x1,x2), max(y1,y2)}
Curves are a bit trickier. You can't just rely on their coordinates, or the coordinates of the control points... you have to do math! (ohnoes!). Thankfully, the math is simple. So simple in fact, that you can pretty much copy-paste what I am about to tell you.

QUADRATIC (OR "2nd ORDER") BEZIER CURVE SEGMENTS

The quadratic bezier curve is a parametric curve (which means it actually uses different functions using the same parameter, for each dimension) of order 2 (meaning the most complex term of the parameter, let's call it "t", will be something involving t2 but never t3 or higher). I'll be giving you the solutions for 2 dimensions here (namely, x and y), but you are probably smart enough to figure out how to add a third dimension (and indeed, why stop at three?).

the bounding box for a quadratic bezier curve
Let's assume the following: our bezier curve runs from coordinate a/b over control point c/d to end point e/f. This gives us the following Bezier function:
F(t)x = a(1-t)2 + 2c(1-t)t + et2
F(t)y = b(1-t)2 + 2d(1-t)t + ft2
These are simple functions. Things get simpler - to find the bounds, we need to find those points (for each dimension) where the curve changes direction: minima and maxima. This is lower high school math tricks 101 - minima and maxima for a function are those points where the derivative equals zero, so let's derive:
F'(t)x = -2a(1-t) + (2c - 4ct) + 2et
F'(t)y = -2b(1-t) + (2d - 4dt) + 2ft
That wasn't so hard. Lets rewrite them so that they become ridiculously easy to do computation with:
F'(t)x = 2t(a-2c+e) + 2(c-a)
F'(t)y = 2t(b-2d+f) + 2(d-b)
and then solve them for 0, so that we find the value of t that we need to plug into the normal function to get the minimum/maximum
0 = 2t(a-2c+e) + 2(c-a)
-2t(a-2c+e) = 2(c-a)
-t(a-2c+e) = (c-a)
t(a-2c+e) = (a-c)
t = ((a-c) / (a-2c+e))
That was certainly easy. The derivative is zero when t is some value that can be computed with simple numbers and arithmetic! There is one snag, of course: if (a-2c+e) is zero, we're in trouble, because then t is undefined (because division by zero is mathematically entirely possible - you get every possible value in the numerical domain you're working in as possible answer. At the same time. Not "some value in that set", but "every value in that set, at the same time". Welcome to division by zero).
When can this happen? Well, first off, we're not going to compute the curve's bounding box unless the control points extend beyond the start/end coordinates, for the very simple reason that the bounding box is the same as for the line segment from start to end coordinate if they are. Why waste good math on a solution that doesn't require it?
Turns out that solves the entire problem, because the only case where (a-2c+e) is zero is when c = ½(a-e). If you already lost track of the math, that means that the control point has to be exactly in between the start and end x coordinate. Clearly, when the control point is between start and end, the bounds are start and end, and we don't need to do any complex math at all. We're safe!
Of course the same goes for the y direction:
0 = 2t(b-2d+f) + 2(d-b)
-2t(b-2d+f) = 2(d-b)
-t(b-2d+f) = (d-b)
t(b-2d+f) = (b-d)
t = (b-d) / (b-2d+f)
Again, this goes horrendously wrong if the control point is exactly between start and end, but in that case we're not using this function anyway. So what do we do to find the bounds?
if segment is quadratic bezier do:
	for both directions do:
		if control between start and end, compute linear bounding box
		otherwise, compute
			bound = u(1-t)2 + 2v(1-t)t + wt2
				(with t = ((u-v) / (u-2v+w)), with {u = start, v = control, w = end})
			if control precedes start, min = bound, otherwise max = bound
And that's it. Bounding box for 2nd order Bezier curve found.

CUBIC (OR "3rd ORDER") BEZIER CURVE SEGMENTS


the bounding box for a cubic bezier curve
The story is virtually the same for cubic curves, except the math is one order higher, so we see a parametric function with as highest order term t3, and we have four points: start point a/b, 'left' control point c/d, 'right' control point e/f and end point g/h.
F(t)x = a(1-t)3 + 3ct(1-t)2 + 3e(1-t)t2 + gt3
F(t)y = b(1-t)3 + 3dt(1-t)2 + 3f(1-t)t2 + ht3
So... minima maxima computerisation. Derivatives:
F'(t)x = 3(-a(t2-2t+1) + c(3t2-4t+1) + e(2t-3t2) + gt2)
F'(t)y = 3(-b(t2-2t+1) + d(3t2-4t+1) + f(2t-3t2) + ht2)
Alright, snag time. These are quadratic functions, which tend to have two possible roots. In order to get to those, we first have to rewrite the derivatives as a "simplified" quadratic formula:
F'(t)x = 3(-(a-3c-g+3e)t2  + 2(a-2c+e)t - a + c)
F'(t)y = 3(-(b-3d-h+3f)t2  + 2(b-2d+f)t - b + d)
And the zero points are now the incredibly wonderful
t = (-u ± sqrt( u2 - 4vw)) / 2w
where
u = (2a - 4c + 2e),
v = (c-a), and
w = (-a+3c+g-3e)
Lovely. Except of course this gives us two possible answers:
t1 = (-u + sqrt( u2 - 4vw)) / 2w
and
t2 = (-u - sqrt( u2 - 4vw)) / 2w
And only one of these will be relevant. You may have noticed this is a lot of computation.
Let's simplify instead!
Rather than using the derivative of the four coordinate function, we can use a three coordinate approximate function. If we use the *delta* coordinates, so that a/b, c/d, e/f and g/h become c-a/d-b, (c-a)-e/(d-b)-f, and ((c-a)-e)-g/((d-b)-f)-h, things get easier. This looks lengthy, but it's an easy computation: if we start at 0,0, with control points at 6,2 and -1,6, and an end point at 1,1 then the rewritten delta coordinates {dx1/dy1, dx2/dy2, dx3/dy3} are {+6/+2, -7/+4, +2/-5}. Sorted. Using those coordinates instead, we get different "derivative" functions:
f(t)x = 3(dx1(1-t)2 + 2dx2(1-t)t + dx3t2)
f(t)x = 3(dy1(1-t)2 + 2dy2(1-t)t + dy3t2)
This is a lot simpler. It gets even better - if we want to find where these functions are zero, that factor 3 is irrelevant (because 3 times 0 is 0), so we just have to solve:
dx1(1-t)2 + 2dx2(1-t)t + dx3t2 = 0
dy1(1-t)2 + 2dy2(1-t)t + dy3t2 = 0
Which is equivalent to solving:
(dx1-2dx2+dx3)t2 - 2(dx1-dx2)t dx1 = 0
(dy1-2dy2+dy3)t2 - 2(dy1-dy2)t dy1 = 0
Hurray for simpler math! High school abc-rule time:
t0x = ((2dx1-2dx2) ± sqrt((2dx2-2dx1)2 - 4dx1(dx1-2dx2+dx3)))/2(dx1-2dx2+dx3)
t0y = ((2dy1-2dy2) ± sqrt((2dy2-2dy1)2 - 4dy1(dy1-2dy2+dy3)))/2(dy1-2dy2+dy3)
And, of course, we observe that dx1+dx3 may not equal 2 dx2 (and likewise for dy) or we will be dividing by zero. However, if we remember why we're doing this math, namely to compute the bounding box of a bezier curve, we can always "salt" our values ever so slightly to get round this slight mathematical problem.
Before we continue, any value of t that we get from this that is negative, or greater than 1, can immediately be discarded, because the whole idea of bezier curves is that they are defined only for t ∊ [0,1].
With the remaining values of t, we simply test what comes out of the original functions,
F(t)x = a(1-(t0x))3 + 3c(t0x)(1-(t0x))2 + 3e(1-(t0x))(t0x)2 + g(t0x)3
F(t)y = b(1-(t0y))3 + 3d(t0y)(1-(t0y))2 + 3f(1-(t0y))(t0y)2 + h(t0y)3
And do a min/max comparison. Hurray, we have now determined the bounds for an arbitrary bezier curve segment of order three in two dimensions!
Our 3rd order procedure is:
if segment is cubic bezier do:
	for both directions do:
		if control between start and end, compute linear bounding box
		otherwise, if one control is outside the start-end interval,
		and the other is at the start or end point:
			flip, recompute normally, then flip-correct the computed outcome
		otherwise, compute bounds using:
			a(1-(t0))3 + 3c(t0)(1-(t0))2 + 3e(1-(t0))(t0)2 + g(t0)3
To see what this code lets us do, let's have a look at an animated "de Casteljau" algorithm, which subdivides any bezier curve (of any order!) into two subcurves at some value of t. processingjs.nihongoresources.com/decasteljau - this page even comes with links to the source code, in which you'll find the theory you've just read actually implemented in a reasonably understandable programming language (BezierCurve class definition, the calculate_straightened_bbox() and calculate_standard_bbox() functions)