iRoCS Toolbox  1.1.0
Data Structures | Public Member Functions | Friends
atb::BSpline< ControlPointT > Class Template Reference

The BSpline class provides functions for fitting B-Splines to point clouds and evaluating them at arbitrary curve position u. More...

#include <ATBSpline.hh>

Inheritance diagram for atb::BSpline< ControlPointT >:
Collaboration diagram for atb::BSpline< ControlPointT >:

Data Structures

struct  BaseCacheIndex
 The BaseCacheIndex struct provides a sortable quadrupel for uniquely identifying the evaluation of the lth indefinite integral of jth B-spline basis function of degree p at position u. More...
 
struct  BasePolynomialCacheIndex
 The BasePolynomialCacheIndex struct provides a sortable quadrupel to uniquely identify the polynomial of the lth indefinite integral of the jth B-spline basis function of degree p at segment m. More...
 

Public Member Functions

 BSpline (int degree=3, size_t nControlPoints=0)
 Default constructor. More...
 
 BSpline (int degree, std::vector< ControlPointT > const &controlPoints)
 This constructor generates a B-Spline object with given degree $p$ and the specified control point vector of length $m$. More...
 
 BSpline (BSpline< ControlPointT > const &spline)
 Copy constructor. More...
 
 ~BSpline ()
 Destructor. More...
 
double base (double u, size_t index, int derivative=0) const
 Evalation of the B-Spline basis function index at position u along the curve. More...
 
double base (BaseCacheIndex const &i) const
 Evaluation of the (derivative/integral) value of a B-Spline basis function. More...
 
Polynomial< double > const & basePolynomial (size_t segment, size_t index, int derivative=0) const
 Get the polynomial representation of the (ell'th derivative/integral) of the j'th basis function on the m'th B-Spline segment. More...
 
Polynomial< double > const & basePolynomial (BasePolynomialCacheIndex const &idx) const
 Get the polynomial representation of the (ell'th derivative/integral) of the j'th basis function on the m'th B-Spline segment. More...
 
ControlPointT const & controlPoint (size_t index) const
 Read the jth control point. More...
 
std::vector< ControlPointT > const & controlPoints () const
 Read accessor to the control point vector of the B-Spline. More...
 
double curveIntegral (double uFrom, double uTo, double uStep=0.001) const
 Estimate the curve length from curve position uFrom to curve position uTo. More...
 
double curveLengthToU (double origin, double distance) const
 Get the curve parameter u at the specified distance from the specified origin. More...
 
size_t curvePositionToSegment (double u) const
 Compute the index of the spline segment, that corresponds to the curve position u. More...
 
size_t curvePositionToSegment (double u, double tolerance) const
 Compute the index of the spline segment, that corresponds to the curve position u. More...
 
int degree () const
 Returns the B-Spline degree. More...
 
ControlPointT derivative (double u, size_t derivative=1) const
 Computation of the derivative of the spline function at position u. More...
 
BSpline< ControlPointT > derivative (size_t derivative=1) const
 Computation of the derivative of the spline function. More...
 
double extendedCurveIntegral (double uFrom, double uTo, double lBound=0.0, double uBound=1.0, double uStep=0.001) const
 Estimate the curve length from curve position uFrom to curve position uTo. More...
 
double extendedCurveLengthToU (double origin, double distance) const
 Get the curve parameter u at the specified distance from the specified origin. More...
 
ControlPointT extendedDerivative (double u, size_t derivative=1, double lBound=0.0, double uBound=1.0) const
 Computation of the derivative of the spline function at position u. More...
 
double knot (size_t index) const
 Read accessor to the k'th knot of the B-Spline. More...
 
std::vector< double > const & knots () const
 Read accessor to the knot vector of the B-Spline. More...
 
void load (std::string const &fileName, std::string const &groupName, bool throwExceptions=false)
 Load the spline with corresponding meta-data from an hdf5 file with given name, using the groupName provided. More...
 
void load (BlitzH5File const &inFile, std::string const &groupName)
 Load the spline with corresponding meta-data from an hdf5 file with given name, using the groupName provided. More...
 
size_t nControlPoints () const
 Returns the number of control points of the B-Spline. More...
 
size_t nKnots () const
 Returns the number of knots of the B-Spline. More...
 
ControlPointT operator() (double u) const
 Evaluation of the spline curve at position u. More...
 
ControlPointT operator() (double u, double lBound, double uBound) const
 Evaluation of the spline curve at position u. More...
 
BSpline< ControlPointT > & operator= (BSpline< ControlPointT > const &spline)
 Copy assignment operator. More...
 
std::string printBasePolynomial (BasePolynomialCacheIndex idx) const
 Get a string representation of the (ell'th derivative/integral) of the j'th basis function polynomial on the m'th B-Spline segment. More...
 
void save (std::string const &fileName, std::string const &groupName, bool throwExceptions=false) const
 Save the spline with corresponding meta-data to an hdf5 file with given name, using the groupName provided. More...
 
void save (BlitzH5File &outFile, std::string const &groupName) const
 Save the spline with corresponding meta-data to an hdf5 file with given name, using the groupName provided. More...
 
void setControlPoint (size_t index, ControlPointT const &point)
 Set the ControlPoint at the specified index. More...
 
void setControlPoints (std::vector< ControlPointT > const &controlPoints)
 Replace the current control points by the provided ones. More...
 
void setDegree (int degree)
 Set the degree of the B-Spline. More...
 
void setKnot (size_t index, double value)
 Set the value of the knot at the specified index. More...
 
void setKnots (std::vector< double > const &knots)
 With this function the knot Array can be set manually. More...
 
void setNControlPoints (size_t nControlPoints)
 Set the number of control points of the B-Spline. More...
 
void setOpenUniformKnots (double uMin=0.0, double uMax=1.0)
 Generation of an open uniform B-Spline knot configuration. More...
 

Friends

double distance (BSpline< double > const &spline, double x, double &u)
 Compute the distance and the corresponding curve position of the given point to the spline. More...
 
template<int Dim>
double distance (BSpline< blitz::TinyVector< double, Dim > > const &spline, blitz::TinyVector< double, Dim > const &x, double &u)
 Compute the distance and the corresponding curve position of the given point to the spline. More...
 
template<int Dim>
double extendedDistance (BSpline< blitz::TinyVector< double, Dim > > const &spline, blitz::TinyVector< double, Dim > const &x, double &u, double lBound, double uBound)
 Compute the distance and the corresponding curve position of the given point to the spline. More...
 
void fitSplineToPointCloud (BSpline< double > &spline, const std::vector< double > &u, std::vector< double > const &points)
 Fitting of a spline curve to scalar data. More...
 
template<int Dim>
void fitSplineToPointCloud (BSpline< blitz::TinyVector< double, Dim > > &spline, std::vector< double > const &u, std::vector< blitz::TinyVector< double, Dim > > const &points)
 Fitting of a spline curve to vectorial data. More...
 
void fitSplineToSpline (BSpline< double > &spline, BSpline< double > const &reference)
 Fitting of a spline curve to another spline curve. More...
 
template<int Dim>
void fitSplineToSpline (BSpline< blitz::TinyVector< double, Dim > > &spline, BSpline< blitz::TinyVector< double, Dim > > const &reference)
 Fitting of a spline curve to another spline curve. More...
 

Detailed Description

template<typename ControlPointT>
class atb::BSpline< ControlPointT >

The BSpline class provides functions for fitting B-Splines to point clouds and evaluating them at arbitrary curve position u.

It also supports computation of derivatives and integrals of the base functions as well as derivatives and line integrals along the full spline.

Within the B-Spline class a Cache for all basis function evaluations is maintained to speed up the computations. Changing the number of control points or the spline degree invalidates the cache and will immediately clear all cache entries.

Definition at line 68 of file ATBSpline.hh.

Constructor & Destructor Documentation

◆ BSpline() [1/3]

template<typename ControlPointT>
atb::BSpline< ControlPointT >::BSpline ( int  degree = 3,
size_t  nControlPoints = 0 
)

Default constructor.

This constructor generates a B-Spline object with given degree $p$ and the specified number of control points $m$. The knot vector will be initialized according to those parameters with a length of $m+p+1$ knots.

Parameters
degreeThe spline degree $p$
nControlPointsThe number of control points $m$

◆ BSpline() [2/3]

template<typename ControlPointT>
atb::BSpline< ControlPointT >::BSpline ( int  degree,
std::vector< ControlPointT > const &  controlPoints 
)

This constructor generates a B-Spline object with given degree $p$ and the specified control point vector of length $m$.

The knot vector will be initialized according to those parameters with a length of $m+p+1$ knots. To really do evaluations the knot-vector needs to be filled with useful values before.

Parameters
degreeThe spline degree $p$
controlPointsThe control point vector of length $m$

◆ BSpline() [3/3]

template<typename ControlPointT>
atb::BSpline< ControlPointT >::BSpline ( BSpline< ControlPointT > const &  spline)

Copy constructor.

This constructor creates a deep copy of a spline including the cache.

Parameters
splineThe spline object to copy

◆ ~BSpline()

template<typename ControlPointT>
atb::BSpline< ControlPointT >::~BSpline ( )

Destructor.

Member Function Documentation

◆ operator=()

template<typename ControlPointT>
BSpline<ControlPointT>& atb::BSpline< ControlPointT >::operator= ( BSpline< ControlPointT > const &  spline)

Copy assignment operator.

This operator deep copies the given spline including the cache.

Parameters
splineThe spline object to copy

◆ degree()

template<typename ControlPointT>
int atb::BSpline< ControlPointT >::degree ( ) const

Returns the B-Spline degree.

Returns
degree of the B-Spline

◆ setDegree()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setDegree ( int  degree)

Set the degree of the B-Spline.

This implies reinitialization of the knot vector and clearing the basis function Cache.

Parameters
degreeThe new B-Spline degree

◆ nControlPoints()

template<typename ControlPointT>
size_t atb::BSpline< ControlPointT >::nControlPoints ( ) const

Returns the number of control points of the B-Spline.

Returns
number of B-Spline control points

◆ setNControlPoints()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setNControlPoints ( size_t  nControlPoints)

Set the number of control points of the B-Spline.

This implies resizing of the knot vector and clearing the basis function cache.

Parameters
nControlPointsThe new number of control points

◆ setControlPoints()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setControlPoints ( std::vector< ControlPointT > const &  controlPoints)

Replace the current control points by the provided ones.

If the number of control points changes that way, the basis function cache will be resetted, leading to a reset of the knot vector to open uniform behaviour.

Parameters
controlPointsThe control points to define the spline shape

◆ controlPoints()

template<typename ControlPointT>
std::vector<ControlPointT> const& atb::BSpline< ControlPointT >::controlPoints ( ) const

Read accessor to the control point vector of the B-Spline.

Returns
The B-Spline's control point vector

◆ controlPoint()

template<typename ControlPointT>
ControlPointT const& atb::BSpline< ControlPointT >::controlPoint ( size_t  index) const

Read the jth control point.

Parameters
indexThe index of the control point to read.
Returns
A reference to the control point (in the scalar case only a double)

◆ setControlPoint()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setControlPoint ( size_t  index,
ControlPointT const &  point 
)

Set the ControlPoint at the specified index.

Parameters
indexThe index of the control point to set.
pointThe control point position

◆ nKnots()

template<typename ControlPointT>
size_t atb::BSpline< ControlPointT >::nKnots ( ) const

Returns the number of knots of the B-Spline.

The number of knots is computed from the number of control points and the degree of the spline, and is thus read-only, and only a convenience-function.

Returns
number of B-Spline knots

◆ setKnots()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setKnots ( std::vector< double > const &  knots)

With this function the knot Array can be set manually.

This has to be done after setting the degree and number of control points. A change in either the degree or the number of control points resets the knot vector to behave like an open uniform B-Spline of the corresponding degree. So make sure, that this function is the last step in B-Spline preparation!

Parameters
knotsThe knot vector to be set
Exceptions
ArrayToolsErrorIf the knot vector has the wrong length or is not monotoneous this error is thrown

◆ knot()

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::knot ( size_t  index) const

Read accessor to the k'th knot of the B-Spline.

Parameters
indexThe knot index
Returns
The k'th knot value

◆ knots()

template<typename ControlPointT>
std::vector<double> const& atb::BSpline< ControlPointT >::knots ( ) const

Read accessor to the knot vector of the B-Spline.

Returns
The B-Spline's knot vector

◆ setKnot()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setKnot ( size_t  index,
double  value 
)

Set the value of the knot at the specified index.

If the value would break the monotonicity of the knot vector the knot vector will not be changed and an ArrayToolsError is thrown.

Parameters
indexThe knot index
valueThe knot value
Exceptions
ArrayToolsErrorIf the monotonicity invariant would be broken after the update this error is thrown

◆ setOpenUniformKnots()

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::setOpenUniformKnots ( double  uMin = 0.0,
double  uMax = 1.0 
)

Generation of an open uniform B-Spline knot configuration.

This implies clearing the basis function cache.

Parameters
uMinThe minimum spline parameter value
uMaxThe maximum spline parameter value

◆ base() [1/2]

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::base ( double  u,
size_t  index,
int  derivative = 0 
) const

Evalation of the B-Spline basis function index at position u along the curve.

You can also evaluate derivatives and integrals of the spline by selecting the degree of the indefinte integral as third parameter.

Parameters
uThe position along the curve to evaluate the B-Spline basis function at.
indexThe index of the basis function to evaluate.
derivativeThe degree of the indefinite integral to evaluate. Zero gives the original spline basis function while positive values give the corresponding integral and negative values the corresponding derivative.
Returns
The B-Spline basis function evaluation at position u

◆ base() [2/2]

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::base ( BaseCacheIndex const &  i) const

Evaluation of the (derivative/integral) value of a B-Spline basis function.

Index, position, spline degree and derivative degree are passed via the BaseCacheIndex struct i.

Parameters
iThe B-Spline basis evaluation parameters.
See also
BaseCacheIndex
Returns
The B-Spline basis function evaluation

◆ basePolynomial() [1/2]

template<typename ControlPointT>
Polynomial<double> const& atb::BSpline< ControlPointT >::basePolynomial ( size_t  segment,
size_t  index,
int  derivative = 0 
) const

Get the polynomial representation of the (ell'th derivative/integral) of the j'th basis function on the m'th B-Spline segment.

Warning: The integral is not the integral of the basis function, it is the integral of the polynomial, that makes up the basis function in segment m. To obtain the integral of the basis function add the contributions of the polynomial pieces for this basis function over all segments.

Parameters
indexThe index j of the base function to get the polynomial for
segmentThe segment m
derivativeThe degree of the indefinite integral to compute. Zero gives the polynomial representation of the original spline basis function while positive values give the corresponding integral and negative values the corresponding derivative.
Returns
The polynomial representation of the j'th basis function on the m'th segment.

◆ basePolynomial() [2/2]

template<typename ControlPointT>
Polynomial<double> const& atb::BSpline< ControlPointT >::basePolynomial ( BasePolynomialCacheIndex const &  idx) const

Get the polynomial representation of the (ell'th derivative/integral) of the j'th basis function on the m'th B-Spline segment.

Parameters
idxThe B-Spline basis polynomial parameters.
See also
BasePolynomialCacheIndex
Returns
The polynomial representation of the j'th basis function on the m'th segment.

◆ printBasePolynomial()

template<typename ControlPointT>
std::string atb::BSpline< ControlPointT >::printBasePolynomial ( BasePolynomialCacheIndex  idx) const

Get a string representation of the (ell'th derivative/integral) of the j'th basis function polynomial on the m'th B-Spline segment.

Parameters
idxThe B-Spline basis polynomial parameters.
See also
BasePolynomialCacheIndex
Returns
The string representation of the j'th basis function polynomial on the m'th segment.

◆ curvePositionToSegment() [1/2]

template<typename ControlPointT>
size_t atb::BSpline< ControlPointT >::curvePositionToSegment ( double  u) const

Compute the index of the spline segment, that corresponds to the curve position u.

If the position given is not within the knot range the number of knots is returned. Make sure to check this value before using it as index for accessing the corresponding knots, otherwise a segmentation fault is very likely to happen or even worse no segmentation fault will occur and you overwrite foreign memory.

Parameters
uThe position along the curve
Returns
The spline segment index, or the number of knots of the spline (one behind the last segment) if u is not in the knot range

◆ curvePositionToSegment() [2/2]

template<typename ControlPointT>
size_t atb::BSpline< ControlPointT >::curvePositionToSegment ( double  u,
double  tolerance 
) const

Compute the index of the spline segment, that corresponds to the curve position u.

If the position given is not within the knot range the number of knots is returned. Make sure to check this value before using it as index for accessing the corresponding knots, otherwise a segmentation fault is very likely to happen or even worse no segmentation fault will occur and you overwrite foreign memory. To allow small round-off errors at the spline boundaries, you can set a tolerance, indicating the allowed deviation from the valid spline range to be accepted as beginning and end of the spline.

Parameters
uThe position along the curve
tolerancePoints within the tolerance interval before or after the valid range will still return the first or last segment
Returns
The spline segment index, or the number of knots of the spline (one behind the last segment) if u is not in the knot range

◆ operator()() [1/2]

template<typename ControlPointT>
ControlPointT atb::BSpline< ControlPointT >::operator() ( double  u) const

Evaluation of the spline curve at position u.

Parameters
uThe position to evaluate the spline at
Returns
The spline value (position on the curve)

◆ operator()() [2/2]

template<typename ControlPointT>
ControlPointT atb::BSpline< ControlPointT >::operator() ( double  u,
double  lBound,
double  uBound 
) const

Evaluation of the spline curve at position u.

This version also takes a lower bound and an upper bound. If u exceeds these bounds, the spline value is extrapolated linearly using the spline tangent at the according bound

\[ f'(u) = \left\{ \begin{array}{l@{\qquad}l} f(\mathrm{lBound}) + \frac{\mathrm{d}}{\mathrm{d}u} f(\mathrm{lBound}) \cdot (u - lBound) & u < \mathrm{lBound} \\ f(u) & \mathrm{lBound} \leq u \leq \mathrm{uBound} \\ f(\mathrm{uBound}) + \frac{\mathrm{d}}{\mathrm{d}u} f(\mathrm{uBound}) \cdot (u - rBound) & \mathrm{uBound} < u \end{array} \right. \]

Parameters
uThe position to evaluate the spline at
lBoundThe spline curve parameter up to which the spline linearly extrapolates
uBoundThe spline curve parameter from which the spline linearly extrapolates
Returns
The spline value (position on the curve)

◆ derivative() [1/2]

template<typename ControlPointT>
ControlPointT atb::BSpline< ControlPointT >::derivative ( double  u,
size_t  derivative = 1 
) const

Computation of the derivative of the spline function at position u.

Parameters
uThe position to get the spline curve derivative at
derivativeThe derivative degree
Returns
The derivative of the spline curve at position u

◆ derivative() [2/2]

template<typename ControlPointT>
BSpline<ControlPointT> atb::BSpline< ControlPointT >::derivative ( size_t  derivative = 1) const

Computation of the derivative of the spline function.

Parameters
derivativeThe degree of the derivative
Returns
The lower degree spline resulting from derivation of the spline by the given derivative

◆ curveIntegral()

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::curveIntegral ( double  uFrom,
double  uTo,
double  uStep = 0.001 
) const

Estimate the curve length from curve position uFrom to curve position uTo.

In fact the lengths of the linear segments of uStep length are added up. If uFrom is greater than uTo the distance will be returned with negative sign.

Parameters
uFromThe position on the curve to start measuring
uToThe position on the curve to stop measuring
uStepThe resolution of the approximation. Instead of the real curve length, the curve is split up into a poly-line with nodes of distance uStep along the curve parameterization. The curve length is then the sum of the segment lengths.
Returns
The curve length between curve positions uFrom and uTo

◆ curveLengthToU()

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::curveLengthToU ( double  origin,
double  distance 
) const

Get the curve parameter u at the specified distance from the specified origin.

Parameters
originThe Reference position on the spline to start measuring the distance
distanceThe distance to travel along the spline
Returns
The spline parameter u at distance from origin

◆ extendedCurveLengthToU()

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::extendedCurveLengthToU ( double  origin,
double  distance 
) const

Get the curve parameter u at the specified distance from the specified origin.

The Spline is extended linearly beyond its end points.

Parameters
originThe Reference position on the spline to start measuring the distance
distanceThe distance to travel along the spline
Returns
The spline parameter u at distance from origin

◆ extendedDerivative()

template<typename ControlPointT>
ControlPointT atb::BSpline< ControlPointT >::extendedDerivative ( double  u,
size_t  derivative = 1,
double  lBound = 0.0,
double  uBound = 1.0 
) const

Computation of the derivative of the spline function at position u.

This function also extrapolates. It returns the usual spline derivative within the [lBound, uBound] interval. Ourside the interval the derivative at the corresponding bound is returned, which resembles linear extrapolation.

Parameters
uThe position to get the spline curve derivative at
derivativeThe derivative degree
lBoundThe lower bound up to which linear extrapolation is performed
uBoundThe upper bound. Starting from this curve parameter the spline is linearly extended
Returns
The derivative of the spline curve at position u

◆ extendedCurveIntegral()

template<typename ControlPointT>
double atb::BSpline< ControlPointT >::extendedCurveIntegral ( double  uFrom,
double  uTo,
double  lBound = 0.0,
double  uBound = 1.0,
double  uStep = 0.001 
) const

Estimate the curve length from curve position uFrom to curve position uTo.

In fact the lengths of the linear segments of uStep length are added up. If uFrom is greater than uTo the distance will be returned with negative sign. This function allows to give out of bounds values and extrapolates linearly.

Parameters
uFromThe position on the curve to start measuring
uToThe position on the curve to stop measuring
lBoundThe lower bound up to which the spline is linearly extrapolated
uBoundThe upper bound from which linear interpolation starts
uStepThe resolution of the approximation. Instead of the real curve length, the curve is split up into a poly-line with nodes of distance uStep along the curve parameterization. The curve length is then the sum of the segment lengths.
Returns
The curve length between curve positions uFrom and uTo

◆ save() [1/2]

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::save ( std::string const &  fileName,
std::string const &  groupName,
bool  throwExceptions = false 
) const

Save the spline with corresponding meta-data to an hdf5 file with given name, using the groupName provided.

Parameters
fileNameThe HDF5 file name. If a file of that name already exists, the spline structures will be added, otherwise a new file will be created
groupNameThe HDF5 group to write the spline data to
throwExceptionsIf this flag is given the method throws exceptions if something goes wrong, otherwise it will report a warning via standard out and return to the caller
Exceptions
BlitzH5ErrorIf the throwExceptions flag is set this error is thrown

◆ load() [1/2]

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::load ( std::string const &  fileName,
std::string const &  groupName,
bool  throwExceptions = false 
)

Load the spline with corresponding meta-data from an hdf5 file with given name, using the groupName provided.

Parameters
fileNameThe HDF5 file name.
groupNameThe HDF5 group to read the spline data from
throwExceptionsIf this flag is given the method throws exceptions if something goes wrong, otherwise it will report a warning via standard out and return to the caller
Exceptions
BlitzH5ErrorIf the throwExceptions flag is set this error is thrown

◆ save() [2/2]

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::save ( BlitzH5File outFile,
std::string const &  groupName 
) const

Save the spline with corresponding meta-data to an hdf5 file with given name, using the groupName provided.

Parameters
outFileThe HDF5 file object the spline structures will be written to.
groupNameThe HDF5 group to write the spline data to
Exceptions
BlitzH5ErrorIf the save fails, this error is thrown

◆ load() [2/2]

template<typename ControlPointT>
void atb::BSpline< ControlPointT >::load ( BlitzH5File const &  inFile,
std::string const &  groupName 
)

Load the spline with corresponding meta-data from an hdf5 file with given name, using the groupName provided.

Parameters
inFileThe HDF5 file object.
groupNameThe HDF5 group to read the spline data from
Exceptions
BlitzH5ErrorOn load error this error is thrown

Friends And Related Function Documentation

◆ fitSplineToSpline [1/2]

template<typename ControlPointT>
void fitSplineToSpline ( BSpline< double > &  spline,
BSpline< double > const &  reference 
)
friend

Fitting of a spline curve to another spline curve.

The function alters the control point positions to make the resulting spline optimally fit the given spline regarding the mean squared error.

Parameters
splineThe output spline
referenceThe spline to be fitted

◆ fitSplineToSpline [2/2]

template<typename ControlPointT>
template<int Dim>
void fitSplineToSpline ( BSpline< blitz::TinyVector< double, Dim > > &  spline,
BSpline< blitz::TinyVector< double, Dim > > const &  reference 
)
friend

Fitting of a spline curve to another spline curve.

The function alters the control point positions to make the resulting spline optimally fit the given spline regarding the mean squared error.

Parameters
splineThe output spline
referenceThe spline to be fitted

◆ fitSplineToPointCloud [1/2]

template<typename ControlPointT>
void fitSplineToPointCloud ( BSpline< double > &  spline,
const std::vector< double > &  u,
std::vector< double > const &  points 
)
friend

Fitting of a spline curve to scalar data.

There must be a functional relation between u and points, thus all entries in u need to be unique otherwise the resulting system of linear equations will be singular. Another necessary condition for the fitting to work is, that each spline segment must at least have one data point to fit it to, thus the number of control points may not exceed the number of data points and they should be well distributed in the spline domain. The function alters the control point positions to make the resulting spline optimally fit the data regarding the mean squared error.

Parameters
splineThe spline to fit to the data values provided
uThe curve position of the corresponding point in the points vector.
pointsThe values to fit the spline to

◆ fitSplineToPointCloud [2/2]

template<typename ControlPointT>
template<int Dim>
void fitSplineToPointCloud ( BSpline< blitz::TinyVector< double, Dim > > &  spline,
std::vector< double > const &  u,
std::vector< blitz::TinyVector< double, Dim > > const &  points 
)
friend

Fitting of a spline curve to vectorial data.

There must be a functional relation between u and points, thus all entries in u need to be unique otherwise the resulting system of linear equations will be singular. Another necessary condition for the fitting to work is, that each spline segment must at least have one data point to fit it to, thus the number of control points may not exceed the number of data points and they should be well distributed in the spline domain. The function alters the control point positions to make the resulting spline optimally fit the data regarding the mean squared error.

Parameters
splineThe spline to fit to the data points provided
uThe curve position of the corresponding point in the points vector.
pointsThe points to fit the spline to

◆ distance [1/2]

template<typename ControlPointT>
double distance ( BSpline< double > const &  spline,
double  x,
double &  u 
)
friend

Compute the distance and the corresponding curve position of the given point to the spline.

Parameters
splineA spline
xA point
uThe curve position with minimum distance to x is returned in this reference
Returns
The distance between the spline and the point

◆ distance [2/2]

template<typename ControlPointT>
template<int Dim>
double distance ( BSpline< blitz::TinyVector< double, Dim > > const &  spline,
blitz::TinyVector< double, Dim > const &  x,
double &  u 
)
friend

Compute the distance and the corresponding curve position of the given point to the spline.

Parameters
splineA spline
xA point
uThe curve position with minimum distance to x is returned in this reference
Returns
The distance between the spline and the point

◆ extendedDistance

template<typename ControlPointT>
template<int Dim>
double extendedDistance ( BSpline< blitz::TinyVector< double, Dim > > const &  spline,
blitz::TinyVector< double, Dim > const &  x,
double &  u,
double  lBound,
double  uBound 
)
friend

Compute the distance and the corresponding curve position of the given point to the spline.

This function behaves as if the spline continues as a straight line beyond the given bounds.

Parameters
splineA spline
xA point
uThe curve position with minimum distance to x is returned in this reference
lBoundThe lower bound up to which the spline behaves as a straight line
uBoundThe lower bound starting from which the spline behaves as a straight line
Returns
The distance between the spline and the point

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