SSJ
V. labo.

umontreal.iro.lecuyer.functionfit
Class BSpline

java.lang.Object
  extended by umontreal.iro.lecuyer.functionfit.BSpline
All Implemented Interfaces:
MathFunction, MathFunctionWithDerivative, MathFunctionWithFirstDerivative, MathFunctionWithIntegral

public class BSpline
extends Object
implements MathFunction, MathFunctionWithIntegral, MathFunctionWithDerivative, MathFunctionWithFirstDerivative

Represents a B-spline with control points at (Xi, Yi). Let Pi = (Xi, Yi), for i = 0,…, n - 1, be a control point and let tj, for j = 0,…, m - 1 be a knot. A B-spline of degree p = m - n - 1 is a parametric curve defined as

P(t) = ∑i=0n-1Ni, p(t)Pi, for tp <= t <= tm-p-1.

Here,
Ni, p(t) = ((t - ti)/(ti+p - ti))Ni, p-1(t) + ((ti+p+1 - t)/(ti+p+1 - ti+1))Ni+1, p-1(t)  
Ni, 0(t) = {$\displaystyle \begin{array}{ll}
 1&\mbox{ for }t_i\le t\le t_{i+1},\\
 0&\mbox{ elsewhere}.
 \end{array}$  

This class provides methods to evaluate P(t) = (X(t), Y(t)) at any value of t, for a B-spline of any degree p >= 1. Note that the evaluate method of this class can be slow, since it uses a root finder to determine the value of t* for which X(t*) = x before it computes Y(t*).


Constructor Summary
BSpline(double[] x, double[] y, double[] knots)
          Constructs a new uniform B-spline with control points at (x[i], y[i]), and knot vector given by the array knots.
BSpline(double[] x, double[] y, int degree)
          Constructs a new uniform B-spline of degree degree with control points at (x[i], y[i]).
 
Method Summary
static BSpline createApproxBSpline(double[] x, double[] y, int degree, int h)
          Returns a B-spline curve of degree degree smoothing (xi, yi), for i = 0,…, n points.
static BSpline createInterpBSpline(double[] x, double[] y, int degree)
          Returns a B-spline curve of degree degree interpolating the (xi, yi) points.
 double derivative(double u)
          Computes (or estimates) the first derivative of the function at point x.
 double derivative(double u, int n)
          .
 BSpline derivativeBSpline()
          Returns the derivative B-spline object of the current variable.
 BSpline derivativeBSpline(int i)
          Returns the ith derivative B-spline object of the current variable; i must be less than the degree of the original B-spline.
 double evaluate(double u)
          Returns the value of the function evaluated at x.
 double evalX(double u)
           
 double evalY(double u)
           
 double[] getKnots()
          Returns an array containing the knot vector (t0, tm-1).
 double getMaxKnot()
          Returns the knot maximal value.
 double getMinKnot()
          Returns the knot minimal value.
 double[] getX()
          Returns the Xi coordinates for this spline.
 double[] getY()
          Returns the Yi coordinates for this spline.
 double integral(double a, double b)
          Computes (or estimates) the integral of the function over the interval [a, b].
 
Methods inherited from class java.lang.Object
equals, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
 

Constructor Detail

BSpline

public BSpline(double[] x,
               double[] y,
               int degree)
Constructs a new uniform B-spline of degree degree with control points at (x[i], y[i]). The knots of the resulting B-spline are set uniformly from x[0] to x[n-1].

Parameters:
x - the values of X.
y - the values of Y.
degree - the degree of the B-spline.

BSpline

public BSpline(double[] x,
               double[] y,
               double[] knots)
Constructs a new uniform B-spline with control points at (x[i], y[i]), and knot vector given by the array knots.

Parameters:
x - the values of X.
y - the values of Y.
knots - the knots of the B-spline.
Method Detail

getX

public double[] getX()
Returns the Xi coordinates for this spline.

Returns:
the Xi coordinates.

getY

public double[] getY()
Returns the Yi coordinates for this spline.

Returns:
the Yi coordinates.

getMaxKnot

public double getMaxKnot()
Returns the knot maximal value.

Returns:
the Yi coordinates.

getMinKnot

public double getMinKnot()
Returns the knot minimal value.

Returns:
the Yi coordinates.

getKnots

public double[] getKnots()
Returns an array containing the knot vector (t0, tm-1).

Returns:
the knot vector.

createInterpBSpline

public static BSpline createInterpBSpline(double[] x,
                                          double[] y,
                                          int degree)
Returns a B-spline curve of degree degree interpolating the (xi, yi) points. This method uses the uniformly spaced method for interpolating points with a B-spline curve, and a uniformed clamped knot vector, as described in http://www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/.

Parameters:
x - the values of X.
y - the values of Y.
degree - the degree of the B-spline.
Returns:
the B-spline curve.

createApproxBSpline

public static BSpline createApproxBSpline(double[] x,
                                          double[] y,
                                          int degree,
                                          int h)
Returns a B-spline curve of degree degree smoothing (xi, yi), for i = 0,…, n points. The precision depends on the parameter h: 1 <= degree <= h < n, which represents the number of control points used by the new B-spline curve, minimizing the quadratic error

L = ∑i=0n($\displaystyle {\frac{{Y_i - S_i(X_i)}}{{W_i}}}$)2.

This method uses the uniformly spaced method for interpolating points with a B-spline curve and a uniformed clamped knot vector, as described in http://www.cs.mtu.edu/~shene/COURSES/cs3621/NOTES/.

Parameters:
x - the values of X.
y - the values of Y.
degree - the degree of the B-spline.
h - the desired number of control points.
Returns:
the B-spline curve.

derivativeBSpline

public BSpline derivativeBSpline()
Returns the derivative B-spline object of the current variable. Using this function and the returned object, instead of the derivative method, is strongly recommended if one wants to compute many derivative values.

Returns:
the derivative B-spline of the current variable.

derivativeBSpline

public BSpline derivativeBSpline(int i)
Returns the ith derivative B-spline object of the current variable; i must be less than the degree of the original B-spline. Using this function and the returned object, instead of the derivative method, is strongly recommended if one wants to compute many derivative values.

Parameters:
i - the degree of the derivative.
Returns:
the ith derivative.

evaluate

public double evaluate(double u)
Description copied from interface: MathFunction
Returns the value of the function evaluated at x.

Specified by:
evaluate in interface MathFunction
Parameters:
u - value at which the function is evaluated
Returns:
function evaluated at x

integral

public double integral(double a,
                       double b)
Description copied from interface: MathFunctionWithIntegral
Computes (or estimates) the integral of the function over the interval [a, b].

Specified by:
integral in interface MathFunctionWithIntegral
Parameters:
a - the starting point of the interval.
b - the ending point of the interval.
Returns:
the value of the integral.

derivative

public double derivative(double u)
Description copied from interface: MathFunctionWithFirstDerivative
Computes (or estimates) the first derivative of the function at point x.

Specified by:
derivative in interface MathFunctionWithFirstDerivative
Parameters:
u - the point to evaluate the derivative to.
Returns:
the value of the derivative.

derivative

public double derivative(double u,
                         int n)
Description copied from interface: MathFunctionWithDerivative
. \begin{tabb}
 Computes (or estimates) the $n$th derivative
 of the function at p...
 ...hod{umontreal.iro.lecuyer.functions}{MathFunction}{evaluate}{double}.
 \end{tabb}
xthe point to evaluate the derivate to. nthe order of the derivative. the resulting derivative. IllegalArgumentExceptionif n is negative or 0.

Specified by:
derivative in interface MathFunctionWithDerivative

evalX

public double evalX(double u)

evalY

public double evalY(double u)

SSJ
V. labo.

To submit a bug or ask questions, send an e-mail to Pierre L'Ecuyer.