Exercise for the Reader

November 1, 2009

The B-spline Series: Definitions (Abridged)

Filed under: 3D Graphics — Seth Porter @ 4:24 pm

This is the second part of a so-far three part series. See part one for some context.

For a typical background you can look at MathWorld or Wikipedia. These are pretty representative of the descriptions you’ll find on the internet in casual searching: technically correct, with nicely formatted equations, but not a whole lot of actual help if you’re trying to understand the things. In particular, the recursive definition of the basis functions are very clever but incredibly unintuitive as far as what they mean. (I also notice that the links share my own indecisiveness about how to title-case “B-spline”.)

A much more useful writeup can be found in the “B-spline Curves” section in these excellent class notes. This is a nice treatment with lots of pretty pictures, although I think the author makes things somewhat more complex by staying with the raw 0 to 1 parameterization throughout (so the elements of the knot vector end up as fractions rather than integers). Nothing wrong with it, just a different treatment (and I suspect one motivated by uniform representation across many kinds of curves, rather than my focus on uniform representation for all segments of the B-spline to ease high-volume computation).

(I should say at this point that I’m torn about how to present things. It’s very tempting to start with what I know now, tell a nice clean story, and maybe then offer an amusing anecdote about how even I didn’t quite get it at first… but that’s exactly what the textbook I was working from did (well, except for the amusing anecdote), and it didn’t serve me very well. So instead I’ll let the game tell the story, and just lay it out as it happened to me.)

First I have to explain why I was motivated to take a shortcut, which is a mixture of raw intellectual laziness and premature optimization. Starting with the intellectual laziness, here’s the basic definition of a B-spline. First, some definitions and terminology (Note that wherever I say something condescending like “the casual reader could be confused”, I really mean “I was confused for about a day”):

The parametric input variable. In general usage, u definitionally ranges between 0 and 1 — this defines the extent of any parametric curve (one where x=f(u) and y=g(u)). However, spline literature, or at least the text I’m working from, is pretty casual about redefining u to work over various different ranges as convenient. To avoid confusion, I’ll use u in the traditional 0..1 sense, and v and w to track these other ranges.
The list of control points. This is probably the most important input to the B-spline curve. Note that indexing is zero-based, and there are n+1 control points. Note that P is a list of vectors, and their dimensionality wholly determines the dimensionality of the output. In other words, if P is a list of (x,y) points, then the position of a B-spline defined on P at a given u will also be an (x,y) point. Likewise if P is defined in 3-space, that is (x,y,z) then the resulting curve is also defined in 3-space.
One less than the number of control points.
An input to the B-spline curve, this determines the degree of the polynomials used for each curve segment, as well as various other things. Most of this post considers the case of k=3.
Generally: the highest exponent appearing in a polynomial: the degree of u2+3 u + 1 is 2. Specifically, for the B-spline segments, degree = k-1 (so my working example with k=3 produces a quadratic on u).
A variable in the range [ 0 , n-k+2 ]. This range is mapped to the 0 to 1 range of u in the obvious way: v = u (n-k+2)
T0..n-k+2 T0..n+k+1
[Edit 11/03/09 – there are n+k+2 total knots (numbered 0 through n+k+1), even though the highest number appearing in the knot vector is n-k+2.] The knot vector. Note that the knot vector’s indexing is zero based. In the case we’re considering (uniform end-interpolating B-splines) it consists of a list of integers between 0 and n-k+2. The first k entries are 0, then the numbers 1 through n-k+1, then k copies of n-k+2. For example, the knot vector for a 6 control point (n=5), second degree (k=3) B-spline would be T = { 0, 0, 0, 1, 2, 3, 4, 4, 4 }. Note that these are not quite indices into the control point array (for one thing, they wouldn’t reference all n+1 control points). These values actually are in the same parameter space as v, above, and basically they serve to knock out terms in the basis functions (when the same value appears repeatedly in the knot vector).
A B-spline basis function, parameterized on i and k. Defined in two parts, asN_i,1(v) := (if t_i <= v <= t_{i+1} then 1 else 0)


N_i,k(v) := (u- T[i])* N(i, k-1, u) / (T[i+k-1] - T[i]) + (T[i+k] - u) * N(i+1, k-1, u) / (T[i+k] - T[i+1])

where 0/0 is defined to be 0. This was the intimidating expression that inspired me to take any shortcut presented. In retrospect, I should have noticed that the denominator would be 1 in both terms unless you’re near an end of the curve, in which case one and then both terms will go to 0 (thanks to the slightly odd definition of dividing by zero). The rest of the time, the Ni,1 terms will act as switches, picking a set of terms from each basis function which will be non-zero in the current section of the curve. Also, note that in much of the literature, this function is written as a function of u, but we’ve already made the jump to the v range as described above, so I replaced the u’s with v’s.

The actual B-spline curve, given a list of control points p and a degree-controlling k. This function is sometimes referred to as p(u), but this is easily confused with the list of control points. Defined asBSpline(n, k, u) := sum(P[i] * N(i, k, u), i, 0, n);

The key thing you’ll notice in this section is that B-splines come out of the mathematical world, not the computer world. The core clever trick is using the lowest-order basis functions as “switches” to select which terms are non-zero. (To their credit, though, at least the lists all use zero-based indices, which will save some tearing-out-of-hair later. (Well, it would if Maxima, the symbolic math program I’m using, didn’t insist on one-based indices…)) This allows a very pretty representation of the function as a whole (see B(v) above), at the expense of hiding all the conditionals way down deep inside these (recursively defined) basis functions. This makes for pretty math, but a hideously slow and ugly direct implementation. As far as I can tell, the intent was to get the top-level B(v) function to look a lot like the definition of Bezier and Hermite curves. Mission accomplished, but with a twist — the definitions of those functions are quite friendly to naive implementation, while this one is pure pain. See http://www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/spline/B-spline/bspline-curve-coef.html for a discussion of this and a different evaluation approach than I will examine here.


I looked at all that, winced hard, and skipped ahead to where they’d solved the uniform version, giving the nice-ish equation for a single segment of

$$B_{uniform,k=3,i}(u)=\frac{1}{2}\,\left( -{p}_{i}\,\left( 2\,{u}^{2}-2\,u-1\right) +{p}_{i+1}\,{u}^{2}+{p}_{i-1}\,{\left( u-1\right) }^{2}\right)$$

or in matrix form,

$$B_{uniform,k=3,i}(u)=\frac{1}{2}\,\left( \begin{pmatrix}{u}^{2} & u & 1\end{pmatrix}.\begin{pmatrix}1 & -2 & 1\cr -2 & 2 & 0\cr 1 & 1 & 0\end{pmatrix}.\begin{pmatrix}{P}_{i-1}\cr {P}_{i}\cr {P}_{i+1}\end{pmatrix}\right)$$

This was what I wanted anyway, an invariant equation where I could just feed it different parameter values and it would produce an x and y. To prevent completely following the path of my own confusion, I’ll point out here that this function is actually defined in terms of

Yet another parameterization. In this case, w ranges from 0 to 1 within each curve segment. This is a reinterpretation of the v parameter as the sum of an integer part (the segment number) and a fractional part (w).
The curve segment number, not to be confused with the i which goes from 0 to n in the definition of B above. To tie these various interpretations of parameter together, v=u (n-k+2)=j+w-1. j ranges from 1 to n-k+2. Note that j is one-based, NOT zero-based (hence the -1 offset when mapped from the 0 to n-k+2 range). This is because at the end of the curve, j=n-k+2 and w=1 — effectively, w is adding one to the range of j (so we need a -1 correction to get the ends to line up right). At the beginning of the curve, w=0, so j must be 1 or we’d again go out of range. A more intuitive interpretation of this is that j is the count of intervals, so this is a pretty standard fencepost correction (http://en.wikipedia.org/wiki/Off-by-one_error).

In terms of these renamed variables, the uniform equation would be

$$B_{uniform,k=3,j}(w)=\frac{1}{2}\,\left( -{p}_{j}\,\left( 2\,{w}^{2}-2\,w-1\right) +{p}_{j+1}\,{w}^{2}+{p}_{j-1}\,{\left( w-1\right) }^{2}\right)$$

So hurray, there was an easy and uniform equation after all; I could either pass in a different window of three control points for each segment, or I could pass in the segment number in addition to the segment-relative parameter (w) and look up the right control points.

The first problem I encountered was that the curve created this way doesn’t actually interpolate the end points. That is to say, when v=0 the curve isn’t at p0. Here’s where I really took a shortcut. Instead of backing off and trying to figure out what was going on, I noticed two things:

  • The knot vector sort of looks like index indirection lookup table. I think this is an especially Computer Science-ish assumption — adding a layer of indirection is a pretty standard way to make non-uniform processes look uniform (see below for a more detailed discussion of what I mean by “uniform”), and a table of integers looks a lot like a list of indices. Particularly because the knot vector repeats its first and last values.
  • If you feed the Buniform function the points [p0 , p0 , p1] (pretending that the control point list had multiple copies of the end-points), and evaluate it at v=0, it does actually evaluate to p0.

I was actually working with k=4 for these initial experiments; I only dropped back to k=3 in an attempt to simplify things to figure out what was going on. The effect I saw is magnified at higher values of k — in the cubic k=4 case, each Buniform is mixing together 4 control points instead of 3. To get it to interpolate the end-point required two additional segments to get from the last “normal” segment (on [ p0 , p1 , p2 , p3 ]) to actually end up at p0. The last normal segment, evaluated at u=0, gives us 1/6 (P2+4 P1 + P0). When we evaluate on [ p0 , p0 , p1 , p2 ] as the first “weird” term, at u=0 we get 1/6 (P1+5 P0). It’s not until we’re looking at [ p0 , p0 , p0 , p1 ] that we actually get p0 at u=0.

In short, I got a little seduced by algebraic numerology — of course if you have a function which mixes 4 points, and you replace enough of those four points with copies of the same point, you’ll end up with a function which returns that point no matter what you give it. (This is the final reductive case, where all four points are actually copies of p0 .)

The actual problem here is that we’ve wandered pretty far from the original definition of a B-spline. In particular, we have too many curve segments — v is ranging from something like -1 to n-k rather than 1 to n-k+2. (The key difference from the “real” B-spline is that the proper B-spline doesn’t actually have a “normal” segment on [ p0 , p1 , p2 , p3 ] — instead, it’s already started going into special cases.)

The practical upshot of this is that we’re spending far too much of the parameter space near the end-points of the curve, leaving only a small portion for the whole interior of the curve. This could be disastrous for rendering, where it’s very nice to sample at evenly-spaced intervals of u and assume that this will lead to generally evenly-spaced vertices in 3-space: with this scheme we’d end up with a great many sampling points at the beginning and end of the curve, and a paucity in the interior.

Even this difficulty could be overcome, by introducing yet another reparameterization to redistribute the incoming u values more evenly over the length of the curve. However, there was no getting around the problem that the actual 3-space shape of the curve was different than my reference implementation (a 3D modeling tool). I’m always willing to entertain the possibility that the reference is wrong and I’m right (perhaps too willing), but at least I needed to understand what was going on here.


1 Comment »

  1. […] previous three parts of this series, we’d sorted out why B-splines are interesting curves, attempted to sort out B-splines terminology, and decided that it would be really nice to be able to evaluation them with uniform functions. […]

    Pingback by B-Splines: Better Math « Exercise for the Reader — December 5, 2009 @ 9:14 pm

RSS feed for comments on this post. TrackBack URI

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google+ photo

You are commenting using your Google+ account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )


Connecting to %s

Blog at WordPress.com.

%d bloggers like this: