35 #include <blitz/array.h> 46 #ifdef ATBSPLINE_PROFILE 67 template<
typename ControlPo
intT>
285 double knot(
size_t index)
const;
294 std::vector<double>
const &
knots()
const;
309 void setKnot(
size_t index,
double value);
380 size_t segment,
size_t index,
int derivative = 0)
const;
485 ControlPointT
operator()(
double u,
double lBound,
double uBound)
const;
529 double curveIntegral(
double uFrom,
double uTo,
double uStep = 0.001)
const;
578 double u,
size_t derivative = 1,
double lBound = 0.0,
579 double uBound = 1.0)
const;
604 double uFrom,
double uTo,
double lBound = 0.0,
double uBound = 1.0,
605 double uStep = 0.001)
const;
625 void save(std::string
const &fileName, std::string
const &groupName,
626 bool throwExceptions =
false)
const;
644 void load(std::string
const &fileName, std::string
const &groupName,
645 bool throwExceptions =
false);
685 void _setNKnots(
size_t nKnots);
692 void _clearCache()
const;
709 void _updateCache()
const;
712 std::vector<ControlPointT> _controlPoints;
713 mutable std::vector<double> _knots;
716 std::vector< std::vector< std::vector<Polynomial<double>* > > > >
717 _basePolynomialDerivativeCache;
719 std::vector< std::vector< std::vector<Polynomial<double>* > > > >
720 _basePolynomialIntegralCache;
729 BSpline< blitz::TinyVector<double,Dim> > &spline,
730 BSpline< blitz::TinyVector<double,Dim> >
const &reference);
734 std::vector<double>
const &points);
738 BSpline< blitz::TinyVector<double,Dim> > &spline,
739 std::vector<double>
const &
u,
740 std::vector< blitz::TinyVector<double,Dim> >
const &points);
746 BSpline< blitz::TinyVector<double,Dim> >
const &spline,
747 blitz::TinyVector<double,Dim>
const &x,
double &
u);
751 BSpline< blitz::TinyVector<double,Dim> >
const &spline,
752 blitz::TinyVector<double,Dim>
const &x,
double &
u,
753 double lBound,
double uBound);
768 template<
typename ControlPo
intT>
797 BSpline< blitz::TinyVector<double,Dim> > &spline,
798 BSpline< blitz::TinyVector<double,Dim> >
const &reference);
820 std::vector<double>
const &points);
842 BSpline< blitz::TinyVector<double,Dim> > &spline,
843 std::vector<double>
const &
u,
844 std::vector< blitz::TinyVector<double,Dim> >
const &points);
876 BSpline< blitz::TinyVector<double,Dim> >
const &spline,
877 blitz::TinyVector<double,Dim>
const &x,
double &
u);
899 double lBound = 0.0,
double uBound = 1.0);
921 BSpline< blitz::TinyVector<double,Dim> >
const &spline,
922 blitz::TinyVector<double,Dim>
const &x,
double &
u,
923 double lBound = 0.0,
double uBound = 1.0);
927 #include "ATBSpline.icc" ATBPolynomial.hh provides the ATB::Polynomial class to store and do arithmetics within the Polynomial...
BSpline< ControlPointT > & operator=(BSpline< ControlPointT > const &spline)
Copy assignment operator.
The BSpline class provides functions for fitting B-Splines to point clouds and evaluating them at arb...
Exception specialization for error handling within libArrayToolbox.
std::ostream & operator<<(std::ostream &os, BasicTreeNode< KeyT, ContentT > const &n)
std::vector< double > const & knots() const
Read accessor to the knot vector of the B-Spline.
ControlPointT derivative(double u, size_t derivative=1) const
Computation of the derivative of the spline function at position u.
friend 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...
void setDegree(int degree)
Set the degree of the B-Spline.
friend 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...
void setNControlPoints(size_t nControlPoints)
Set the number of control points of the B-Spline.
The BaseCacheIndex struct provides a sortable quadrupel for uniquely identifying the evaluation of th...
double extendedCurveLengthToU(double origin, double distance) const
Get the curve parameter u at the specified distance from the specified origin.
size_t curvePositionToSegment(double u) const
Compute the index of the spline segment, that corresponds to the curve position u.
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 t...
Lightweight alternative to libBlitzHDF5 providing its basic functionality.
void setKnot(size_t index, double value)
Set the value of the knot at the specified index.
size_t nControlPoints() const
Returns the number of control points of the B-Spline.
std::vector< ControlPointT > const & controlPoints() const
Read accessor to the control point vector of the B-Spline.
BSpline(int degree=3, size_t nControlPoints=0)
Default constructor.
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 pro...
size_t nKnots() const
Returns the number of knots of the B-Spline.
BaseCacheIndex(double u, size_t index, int degree, int derivative)
double knot(size_t index) const
Read accessor to the k'th knot of the B-Spline.
ATBTiming.hh provides the MyDateTime and Timer classes for high accuracy profiling.
The BasePolynomialCacheIndex struct provides a sortable quadrupel to uniquely identify the polynomial...
friend void fitSplineToSpline(BSpline< double > &spline, BSpline< double > const &reference)
Fitting of a spline curve to another spline curve.
ControlPointT operator()(double u) const
Evaluation of the spline curve at position u.
friend void fitSplineToPointCloud(BSpline< double > &spline, const std::vector< double > &u, std::vector< double > const &points)
Fitting of a spline curve to scalar data.
bool operator<(BaseCacheIndex const &i) const
void setControlPoints(std::vector< ControlPointT > const &controlPoints)
Replace the current control points by the provided ones.
int degree() const
Returns the B-Spline degree.
ControlPointT const & controlPoint(size_t index) const
Read the jth control point.
void setKnots(std::vector< double > const &knots)
With this function the knot Array can be set manually.
std::string printBasePolynomial(BasePolynomialCacheIndex idx) const
Get a string representation of the (ell'th derivative/integral) of the j'th basis function polynomial...
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 p...
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.
double curveLengthToU(double origin, double distance) const
Get the curve parameter u at the specified distance from the specified origin.
void setOpenUniformKnots(double uMin=0.0, double uMax=1.0)
Generation of an open uniform B-Spline knot configuration.
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.
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.
double curveIntegral(double uFrom, double uTo, double uStep=0.001) const
Estimate the curve length from curve position uFrom to curve position uTo.
void setControlPoint(size_t index, ControlPointT const &point)
Set the ControlPoint at the specified index.