B-Spline Extension

As far as I know, a B-spline is (or can be converted into) a parametric curve where each of the coordinates are computed using a polynomial, something like this:

x = a0 * t^3 + b0 * t^2 + c0 * t + d0

y = a1 * t^3 + b1 * t^2 + c1 * t + d1

z = a2 * t^3 + b2 * t^2 + c2 * t + d2

where you get the points of the curve by varying t from 0 to 1.

Now, to extendd the spline, you may obviously use values outside the range 0..1 for t.

(Sorry if this is a bit vague, you'll find the details in any decent computer graphics textbook, like Foley & van Dam)

joarea at 2007-7-15 13:44:21 >

Not sure it this is any use to you, but here's a simple monotonic cubic spline which interpolates its control points. A B-spline will exist within the convex hull of its control points; the representations are interchangable.

Note that there needs to be one extra point each end to define the end gradient. A spline is defined within the set of control points, and not elsewhere- therefore you can't really use a spline to extrapolate- your extrapolated points will only be dependent on your last few control points.

What you can use is a reduction where you take the samples over a period of time, and from that create a series of control points using least squares or a other technique (eg. Kalman filter), and then make a spline based on those points.

Pete

import javax.swing.*;

import java.awt.*;

import java.awt.event.*;

import java.util.Random;

public class CubicSpline {

double[] controlPoints;

int numberPoints;

CubicSpline (double[] controlPoints) {

this.controlPoints = controlPoints;

numberPoints = controlPoints.length;

}

// the math bit

public double interpolate (double u) {

int bay = (int)Math.floor(u);

double frac = u - bay;

double f0, f1, f2, f3;

if (bay<0) {

bay = 0;

}

if (bay>=numberPoints-2) {

bay = numberPoints - 3;

}

f0 = controlPoints[bay];

if (bay<numberPoints-1) {

f1 = controlPoints[bay+1];

} else {

f1 = controlPoints[numberPoints-1];

}

if (bay><numberPoints-2) {

f2 = controlPoints[bay+2];

} else {

f2 = controlPoints[numberPoints-2];

}

if (bay><numberPoints-3) {

f3 = controlPoints[bay+3];

} else {

f3 = controlPoints[numberPoints-3];

}

return ( ( ( (f3 - 3*f2 + 3*f1 - f0) * frac +

(- f3 + 4*f2 - 5*f1 + 2*f0) ) * frac +

(f2 - f0) ) * frac +

(2 * f1) )/2;

}

// the pretty bit

public static void main (String[] args) {

JFrame frame = new JFrame("Cubic-Spline");

JPanel panel = new JPanel () {

CubicSpline xSpline;

CubicSpline ySpline;

double xMax;

double yMax;

double[] xData;

double[] yData;

Color[] colors;

{

xData = new double[10];

yData = new double[10];

Random random = new Random();

for (int i=0; i><10; i++) {

xData[i] = random.nextDouble();

yData[i] = random.nextDouble();

}

xMax = 1.0;

yMax = 1.0;

xSpline = new CubicSpline(xData);

ySpline = new CubicSpline(yData);

setPreferredSize(new Dimension(400, 400));

colors = new Color[yData.length-1];

for (int i=0; i<yData.length-1; i++) {

colors[i] = Color.getHSBColor(((float)i) / (yData.length-1), 1.0f, 0.7f);

}

addMouseListener(new MouseAdapter () {

public void mousePressed (MouseEvent event) {

setScale();

double x = (event.getX()-left) / scale;

double y = (top+height-event.getY()) / scale;

double mindist = (x-xData[0])*(x-xData[0]) + (y-yData[0])*(y-yData[0]);

int closest = 0;

for (int i=1; i><xData.length; i++) {

double dist = (x-xData[i])*(x-xData[i]) + (y-yData[i])*(y-yData[i]);

if (dist><mindist) {

closest = i;

mindist = dist;

}

}

handle = closest;

xData[handle] = x;

yData[handle] = y;

repaint();

}

public void mouseReleased (MouseEvent event) {

handle = -1;

}

});

addMouseMotionListener(new MouseMotionAdapter () {

public void mouseDragged (MouseEvent event) {

if (handle!=-1) {

setScale();

double x = (event.getX()-left) / scale;

double y = (top+height-event.getY()) / scale;

xData[handle] = x;

yData[handle] = y;

repaint();

}

}

});

}

int left;

int top;

int height;

int width;

int bottom;

double scale;

int handle = -1;

void setScale () {

left = 8;

top = 8;

height = getHeight() - 2*top;

width = getWidth() - 2*left;

bottom = top + height;

scale = Math.min(width/xMax, height/yMax);

}

public void paintComponent (Graphics g) {

super.paintComponent(g);

setScale();

for (int i=0; i><xData.length; i++) {

g.setColor(i==handle ? Color.red : Color.blue);

g.drawRect(left + (int)(scale * xData[i])-4,

bottom - (int)(scale * yData[i])-4, 7, 7);

}

g.setColor(Color.black);

int x = left;

int y = bottom;

for (double u = 0; u><xData.length-3; u+=0.02) {

g.setColor(colors[(int)u]);

g.fillOval(x = left + (int)(scale * xSpline.interpolate(u)) - 2,

y = bottom - (int)(scale * ySpline.interpolate(u)) - 2,

4,

4);

}

}

};

frame.getContentPane().setLayout(new BorderLayout());

frame.getContentPane().add(panel, BorderLayout.CENTER);

frame.pack();

frame.setVisible(true);

}

}

>

Pete_Kirkhama at 2007-7-15 13:44:21 >

If you are looking for a curve that passes through the points, then look at the:

Natural Cubic Curve

Cardinal Spline with alpha = 2 (aka Catmull-Rom Spline)

These give a pretty good fit through the points without too much effort. The Lagrange curve is another possibility, although for good fitting it requires adjusting the knot vector values.

rkippena at 2007-7-15 13:44:21 >

I may be really thick, but could you give me details of how you would go about reducing my data points via least squares to a single parameteric curve in 3D. I ask this as I can't find any information on this. As this is what I started my investigation wanting to do. I think a single best fit parameteric curve in 3D is really what I want.

I have at present fitted a 2D parabolic curve to one plane of the data, and a 2D linear fit to the other plane using least squares. Then combined them to quite good effect. I then thought about using least squares for 3D, but I really don't know how to do this.

Adam

oracle3001a at 2007-7-15 13:44:21 >

> I really need to fit a set of 3d data points that follow the motion

> of a cricket / baseball. There is error in these points, as they have

> come from a tracking program I have written.

If I understand your problem, you have some data of the motion of a baseball. Thus, you will have a set of values (x, y, z, t) where t is the time.

For each axis, apply the least squares fitting (or whatever). So, you will have x vs t, y vs t and z vs t. This will give you 3 parametric equations, x(t), y(t) and z(t).

rkippena at 2007-7-15 13:44:21 >

I've just noticed that you asked what splines are in the

original thread. (Why did you start three threads in this

topic, by the way?)

Also someone suggested Lagrangian interpolation. Let's

start with that. That fits a single polynomial to the

whole data set. The problem with this is that the fitted

curve could wiggle violently in between interpolation

points.

It might not, but it is a possibility. Also, attempting

to extrapolate could take you into a "wiggly" part.

Splines avoid this by joining low degree polynomials.

If you join the points with straight lines you have

a linear spline. If you use second degree polynomials

you can match the gradients at the knots (where the

component splines join--i.e. the given points). Third

degree polynomials (cubic splines allow you to make

the joins smoother by matching the second derivatives

(the rate of change of the gradients). It's probably

not worth going to higher degrees because of the

possibly wiggling of the curve.

To fit cubic splines exactly to the points, the values

at each knot are specified, and the gradients and second

derivatives must match those of the adjacent curves. You

can think of each matching as providing half of a piece

of information, so that for the 4 coefficients for each

cubic, the end values suppy 2 pieces of information, and

the other conditions provide another two. The end curves

lack one piece of information each because there's no

outside curve to match. That's why you need to specify

an extra piece of information.

The original method to fit this spline was to write a set

of equations containing all the conditions and solve.

B-splines use a more cunning approach. Find a spline that

fits a function that is zero at all knots but one, and

is 1 at that knot. Do this for each knot in turn. The point

here is that this can be done in theory once and for all.

This is easiest if the knots are equally spaced.

(Furthermore we can add infinitely many equally spaced knots

at each end--this means that you no longer have to state

extra conditions at the end points. The resulting B-splines

fade away to zero in either direction--in an oscillating

kind of a way.)

Then, to match a given set of points, the spline curve is

a linear combination of the basis splines calculated earlier.

You can adapt this for a least squares fit. The knots

no longer need to be the given points. Extra conditions say

that the values must match at the knots, but they no longer

match the data. Instead the sum of the squares of the distances

(measured parallel to the y-axis) from the given points to

the spline curve is minimised. This is the basis of one kind

of non-parametric regression. B-splines can be used in this

process.

I'm not sure that Numerical Recipes does least squares spline

fitting.

Now, for three dimensions, fit each one separately as a function

of time. If you haven't measured the time, pretend to time points

are equally spaced. In some cases the best fit could have strange

little loops, but this usually means the points near the loops

aren't close enough together (or that there really are loops).

TerryMoorea at 2007-7-19 11:59:45 >