I’m still optimizing memory and speed so nothing exciting to report under the lighting topic. Maybe later this week or next week.

So let’s look at a bit of code from the book, shall we?

- float
- IvBezier::FindParameterByDistance( float t1, float s )
- {
- // ensure that we remain within valid parameter space
- if ( s > ArcLength(t1, mTimes[mCount-1]) )
- return mTimes[mCount-1];
- // make first guess
- float p = t1 + s*(mTimes[mCount-1]-mTimes[0])/mTotalLength;
- for ( u_int i = 0; i < 32; ++i )
- {
- // compute function value and test against zero
- float func = ArcLength(t1, p) - s;
- if ( ::IvAbs(func) < 1.0e-03f )
- {
- return p;
- }
- // perform Newton-Raphson iteration step
- float speed = Velocity(p).Length();
- ASSERT( !::IsZero( speed ) );
- p -= func/speed;
- }
- // done iterating, return failure case
- return FLT_MAX;
- } // End of IvBezier::FindParameterByDistance()

For those who aren’t intimately familiar with curves, this code is used to help generate a constant velocity along a curve. Under normal circumstances, if you were to take the parameter t and vary it at a constant rate, the resulting points along the curve are separated by a non-constant distance. So, for example, if you were to plug in 0.5 as your parameter, you wouldn’t necessarily end up with the point which is halfway along the curve in terms of length. This means that by using the parameter directly, you are moving at non-constant speed.

However, in most cases what you want to is to travel along the curve at a constant speed. What the function above does is to take the distance along the curve (starting at a parameter t1, which is usually set to 0) and generate the parameter for the point that lies at that distance. So to move at constant speed, you start with value s = 0 and increase that at a constant rate up to the maximum length. At each step you can plug s into this function to get the curve parameter, which you can use to generate the actual point.

The analytical way to do this is to solve a length integral. However for all practical purposes the particular integral for this case is impossible to solve. So instead a numerical technique is used — the Newton-Raphson method — which finds the roots, or zeros, of functions. In this case we’re trying to solve for the parameter p that gives us ArcLength(t1, p) = s, so the function we’re trying to find the zero for is ArcLength(t1, p) – s. The code above computes the current function value, determines if we’re close enough to zero yet, and if not, performs a Newton-Raphson step.

Which finally brings us to the point of this blog entry. In most cases, this works pretty well, but there is a flaw. In particular, the Newton-Raphson step is problematic:

- float speed = Velocity(p).Length();
- ASSERT( !::IsZero( speed ) );
- p -= func/speed;

If the speed is zero, or close enough so that floating point error bites us in the backside, then we will end up subtracting NaN from p, which gives us a garbage result. In my initial tests this didn’t happen, but when I used it in production code the designers managed to create some camera path curves that triggered this case. So… back to the drawing board.

The key to solving this is to note that there is another method to finding roots other than Newton-Raphson: the bisection method. This is basically a binary search — cutting the problem interval in half and searching again in the more fruitful half. The problem with bisection is that it converges much slower than Newton-Raphson, so sometimes we want to use Newton-Raphson, and sometimes bisection.

It’s getting late, so I’ll wrap this up later: hopefully tomorrow, by Friday at the latest.