BSplineBase< T > Class Template Reference

The base class for a spline object containing the nodes for a given domain, cutoff wavelength, and boundary condition. More...

#include <BSpline.h>

Inheritance diagram for BSplineBase< T >:

List of all members.

Public Types
typedef T	datum_type
enum	BoundaryConditionTypes { BC_ZERO_ENDPOINTS = 0, BC_ZERO_FIRST = 1, BC_ZERO_SECOND = 2 }
	Boundary condition types. More...
Public Member Functions
	BSplineBase (const T *x, int nx, double wl, int bc_type=BC_ZERO_SECOND, int num_nodes=0)
	Construct a spline domain for the given set of x values, cutoff wavelength, and boundary condition type.
	BSplineBase (const BSplineBase &)
	Copy constructor.
bool	setDomain (const T *x, int nx, double wl, int bc_type=BC_ZERO_SECOND, int num_nodes=0)
	Change the domain of this base.
BSpline< T > *	apply (const T *y)
	Create a BSpline smoothed curve for the given set of NX y values.
const T *	nodes (int *nnodes)
	Return array of the node coordinates.
int	nNodes ()
	Return the number of nodes (one more than the number of intervals).
int	nX ()
	Number of original x values.
T	Xmin ()
	Minimum x value found.
T	Xmax ()
	Maximum x value found.
double	Alpha (double wavelength)
	Return the Alpha value for a given wavelength.
double	Alpha ()
	Return alpha currently in use by this domain.
bool	ok ()
	Return the current state of the object, either ok or not ok.
Static Public Member Functions
static const char *	ImplVersion ()
	Return a string describing the implementation version.
static const char *	IfaceVersion ()
	Return a string describing the interface version.
static bool	Debug (int on=-1)
	Call this class method with a value greater than zero to enable debug messages, or with zero to disable messages.
Protected Types
typedef BSplineBaseP< T >	Base
Protected Member Functions
bool	Setup (int num_nodes=0)
void	calculateQ ()
double	qDelta (int m1, int m2)
double	Beta (int m)
void	addP ()
bool	factor ()
double	Basis (int m, T x)
double	DBasis (int m, T x)
double	Ratiod (int &, double &, double &)
Protected Attributes
double	waveLength
int	NX
int	K
int	BC
T	xmax
T	xmin
int	M
double	DX
double	alpha
bool	OK
Base *	base
Static Protected Attributes
static const double	BoundaryConditions [3][4]
static const double	PI = 3.1415927

---

Detailed Description

template<class T>
class BSplineBase< T >

The base class for a spline object containing the nodes for a given domain, cutoff wavelength, and boundary condition.

To smooth a single curve, the BSpline interface contains a constructor which both sets up the domain and solves for the spline. Subsequent curves over the same domain can be created by apply()ing them to the BSpline object, where apply() is a BSplineBase method. [See apply().] New curves can also be smoothed within the same BSpline object by calling solve() with the new set of y values. [See BSpline.] A BSplineBase can be created on its own, in which case all of the computations dependent on the x values, boundary conditions, and cutoff wavelength have already been completed.

The solution of the cubic b-spline is divided into two parts. The first is the setup of the domain given the x values, boundary conditions, and wavelength. The second is the solution of the spline for a set of y values corresponding to the x values in the domain. The first part is done in the creation of the BSplineBase object (or when calling the setDomain method). The second part is done when creating a BSpline object (or calling solve() on a BSpline object).

A BSpline object can be created with either one of its constructors, or by calling apply() on an existing BSplineBase object. Once a spline has been solved, it can be evaluated at any x value. The following example creates a spline curve and evaluates it over the domain:

    vector<float> x;
    vector<float> y;
    { ... }
    int bc = BSplineBase<float>::BC_ZERO_SECOND;
    BSpline<float>::Debug = true;
    BSpline<float> spline (x.begin(), x.size(), y.begin(), wl, bc);
    if (spline.ok())
    {
        ostream_iterator<float> of(cout, "\t ");
    	float xi = spline.Xmin();
	float xs = (spline.Xmax() - xi) / 2000.0;
	for (; xi <= spline.Xmax(); xi += xs)
	{
	    *of++ = spline.evaluate (xi);
	}
    }
 

In the usual usage, the BSplineBase can compute a reasonable number of nodes for the spline, balancing between a few desirable factors. There needs to be at least 2 nodes per cutoff wavelength (preferably 4 or more) for the derivative constraint to reliably approximate a lo-pass filter. There should be at least 1 and preferably about 2 data points per node (measured just by their number and not by any check of the density of points across the domain). Lastly, of course, the fewer the nodes then the faster the computation of the spline. The computation of the number of nodes happens in the Setup() method during BSplineBase construction and when setDomain() is called. If the setup fails to find a desirable number of nodes, then the BSplineBase object will return false from ok().

The ok() method returns false when a BSplineBase or BSpline could not complete any operation successfully. In particular, as mentioned above, ok() will return false if some problem was detected with the domain values or if no reasonable number of nodes could be found for the given cutoff wavelength. Also, ok() on a BSpline object will return false if the matrix equation could not be solved, such as after BSpline construction or after a call to apply().

If letting Setup() determine the number of nodes is not acceptable, the constructors and setDomain() accept the parameter num_nodes. By default, num_nodes is passed as zero, forcing Setup() to calculate the number of nodes. However, if num_nodes is passed as 2 or greater, then Setup() will bypass its own algorithm and accept the given number of nodes instead. Obviously, it's up to the programmer to understand the affects of the number of nodes on the representation of the data and on the solution (or non-solution) of the spline. Remember to check the ok() method to detect when the spline solution has failed.

The interface for the BSplineBase and BSpline templates is defined in the header file BSpline.h. The implementation is defined in BSpline.cpp. Source files which will instantiate the template should include the implementation file and not the interface. If the implementation for a specific type will be linked from elsewhere, such as a static library or Windows DLL, source files should only include the interface file. On Windows, applications should link with the import library BSpline.lib and make sure BSpline.dll is on the path. The DLL contains an implementation for BSpline<float> and BSpline<double>. For debugging, an application can include the implementation to get its own instantiation.

The algorithm is based on the cubic spline described by Katsuyuki Ooyama in Montly Weather Review, Vol 115, October 1987. This implementation has benefited from comparisons with a previous FORTRAN implementation by James L. Franklin, NOAA/Hurricane Research Division. In particular, the algorithm in the Setup() method is based mostly on his implementation (VICSETUP). The Setup() method finds a suitable default for the number of nodes given a domain and cutoff frequency. This implementation adopts most of the same constraints, including a constraint that the cutoff wavelength not be greater than the span of the domain values: wl < max(x) - min(x). If this is not an acceptable constraint, then use the num_nodes parameter to specify the number of nodes explicitly.

The cubic b-spline is formulated as the sum of some multiple of the basis function centered at each node in the domain. The number of nodes is determined by the desired cutoff wavelength and a desirable number of x values per node. The basis function is continuous and differentiable up to the second degree. A derivative constraint is included in the solution to achieve the effect of a low-pass frequency filter with the given cutoff wavelength. The derivative constraint can be disabled by specifying a wavelength value of zero, which reduces the analysis to a least squares fit to a cubic b-spline. The domain nodes, boundary constraints, and wavelength determine a linear system of equations, Qa=b, where a is the vector of basis function coefficients at each node. The coefficient vector is solved by first LU factoring along the diagonally banded matrix Q in BSplineBase. The BSpline object then computes the B vector for a set of y values and solves for the coefficient vector with the LU matrix. Only the diagonal bands are stored in memory and calculated during LU factoring and back substitution, and the basis function is evaluated as few times as possible in computing the diagonal matrix and B vector.

Author:

Gary Granger (http://www.atd.ucar.edu/~granger)

Copyright (c) 1998-2003
University Corporation for Atmospheric Research, UCAR

---

Member Enumeration Documentation

template<class T> enum BSplineBase::BoundaryConditionTypes

Boundary condition types. Enumeration values: BC_ZERO_ENDPOINTS  Set the endpoints of the spline to zero. BC_ZERO_FIRST  Set the first derivative of the spline to zero at the endpoints. BC_ZERO_SECOND  Set the second derivative to zero.

---

Constructor & Destructor Documentation

template<class T> BSplineBase< T >::BSplineBase (  const T *  x, int  nx, double  wl, int  bc_type = BC_ZERO_SECOND, int  num_nodes = 0 )

Construct a spline domain for the given set of x values, cutoff wavelength, and boundary condition type. The parameters are the same as for setDomain(). Call ok() to check whether domain setup succeeded after construction.

---

Member Function Documentation

template<class T> double BSplineBase< T >::Alpha (  double  wavelength  )

Return the Alpha value for a given wavelength. Note that this depends on the current node interval length (DX).

template<class T> BSpline< T > * BSplineBase< T >::apply (  const T *  y  )

Create a BSpline smoothed curve for the given set of NX y values. The returned object will need to be deleted by the caller. Parameters: y  The array of y values corresponding to each of the nX() x values in the domain. See also: ok()

template<class T> bool BSplineBase< T >::Debug (  int  on = -1  )  [inline, static]

Call this class method with a value greater than zero to enable debug messages, or with zero to disable messages. Calling with no arguments returns true if debugging enabled, else false.

template<class T> const T * BSplineBase< T >::nodes (  int *  nnodes  )

Return array of the node coordinates. Returns 0 if not ok(). The array of nodes returned by nodes() belongs to the object and should not be deleted; it will also be invalid if the object is destroyed.

template<class T> bool BSplineBase< T >::ok (   )  [inline]

Return the current state of the object, either ok or not ok. Use this method to test for valid state after construction or after a call to setDomain(). ok() will return false if either fail, such as when an appropriate number of nodes and node interval cannot be found for a given wavelength, or when the linear equation for the coefficients cannot be solved.

template<class T> bool BSplineBase< T >::setDomain (  const T *  x, int  nx, double  wl, int  bc_type = BC_ZERO_SECOND, int  num_nodes = 0 )

Change the domain of this base. [If this is part of a BSpline object, this method {will not} change the existing curve or re-apply the smoothing to any set of y values.] The x values can be in any order, but they must be of sufficient density to support the requested cutoff wavelength. The setup of the domain may fail because of either inconsistency between the x density and the cutoff wavelength, or because the resulting matrix could not be factored. If setup fails, the method returns false. Parameters: x  The array of x values in the domain. nx  The number of values in the x array. wl  The cutoff wavelength, in the same units as the x values. A wavelength of zero disables the derivative constraint. bc_type  The enumerated boundary condition type. If omitted it defaults to BC_ZERO_SECOND. num_nodes  The number of nodes to use for the cubic b-spline. If less than 2 a reasonable number will be calculated automatically, if possible, taking into account the given cutoff wavelength. See also: ok().

---

Member Data Documentation

template<class T> const double BSplineBase< T >::BoundaryConditions [static, protected]

Initial value: { { -4, -1, -1, -4 }, { 0, 1, 1, 0 }, { 2, -1, -1, 2 } }

---

The documentation for this class was generated from the following files:

• /code/granger/bspline/bspline-1_3/BSpline/BSpline.h
• BSpline.cpp

---