Mantid
|
The ChebfunBase class provides a base for function approximation with Chebyshev polynomials. More...
#include <ChebfunBase.h>
Public Member Functions | |
template<class FunctionType > | |
ChebfunBase_sptr | bestFitTempl (double start, double end, FunctionType f, std::vector< double > &p, std::vector< double > &a, double maxA, double tolerance, size_t maxSize) |
Fit a function until full convergence. More... | |
std::vector< double > | calcA (const std::vector< double > &p) const |
Calculate expansion coefficients. More... | |
std::vector< double > | calcP (const std::vector< double > &a) const |
Calculate function values. More... | |
ChebfunBase (ChebfunBase &&) noexcept=default | |
Move constructor. More... | |
ChebfunBase (const ChebfunBase &)=default | |
Copy constructor. More... | |
ChebfunBase (size_t n, double start, double end, double tolerance=0.0) | |
Constructor. More... | |
void | derivative (const std::vector< double > &a, std::vector< double > &aout) const |
Calculate the derivative. More... | |
double | endX () const |
End of the interval. More... | |
double | eval (double x, const std::vector< double > &p) const |
Evaluate a function. More... | |
std::vector< double > | evalVector (const std::vector< double > &x, const std::vector< double > &p) const |
Evaluate a function. More... | |
void | evalVector (const std::vector< double > &x, const std::vector< double > &p, std::vector< double > &res) const |
Evaluate a function. More... | |
std::vector< double > | fit (ChebfunFunctionType f) const |
Calculate function values at chebfun x-points. More... | |
std::vector< double > | fit (const API::IFunction &f) const |
Calculate function values at chebfun x-points. More... | |
std::shared_ptr< ChebfunBase > | integral (const std::vector< double > &a, std::vector< double > &aout) const |
Calculate the integral. More... | |
double | integrate (const std::vector< double > &p) const |
Calculate an integral. More... | |
const std::vector< double > & | integrationWeights () const |
Get a reference to the integration weights. More... | |
std::vector< double > | linspace (size_t n) const |
Create a vector of x values linearly spaced on the approximation interval. More... | |
size_t | order () const |
Get the polynomial order of this base. More... | |
std::vector< double > | roots (const std::vector< double > &a) const |
Find all roots of a function on this interval. More... | |
size_t | size () const |
Get the size of the base which is the number of x-points. More... | |
std::vector< double > | smooth (const std::vector< double > &xvalues, const std::vector< double > &yvalues) const |
Smooth some data. More... | |
double | startX () const |
Start of the interval. More... | |
double | tolerance () |
Tolerance for comparing doubles. More... | |
double | width () const |
Get the width of the interval. More... | |
const std::vector< double > & | xPoints () const |
Get a reference to the x-points. More... | |
Static Public Member Functions | |
static std::shared_ptr< ChebfunBase > | bestFit (double start, double end, ChebfunFunctionType, std::vector< double > &p, std::vector< double > &a, double maxA=0.0, double tolerance=0.0, size_t maxSize=0) |
Fit a function until full convergence. More... | |
static std::shared_ptr< ChebfunBase > | bestFit (double start, double end, const API::IFunction &, std::vector< double > &p, std::vector< double > &a, double maxA=0.0, double tolerance=0.0, size_t maxSize=0) |
Fit a function until full convergence. More... | |
template<class FunctionType > | |
static std::shared_ptr< ChebfunBase > | bestFitAnyTolerance (double start, double end, FunctionType f, std::vector< double > &p, std::vector< double > &a, double maxA=0.0, double tolerance=0.0, size_t maxSize=0) |
Find best fit with highest possible tolerance (to be used with noisy data). More... | |
Private Member Functions | |
void | calcIntegrationWeights () const |
Calculate the integration weights. More... | |
void | calcX () |
Calculate the x-values based on the (start,end) interval. More... | |
std::vector< double > | fitOdd (const API::IFunction &f, std::vector< double > &pEven) const |
Calculate function values at odd-valued indices of the base x-points. More... | |
std::vector< double > | fitOdd (const ChebfunFunctionType &f, std::vector< double > &p) const |
Calculate function values at odd-valued indices of the base x-points. More... | |
ChebfunBase & | operator= (ChebfunBase &&other)=delete |
ChebfunBase & | operator= (const ChebfunBase &other)=delete |
Private assingment operator to stress the immutability of ChebfunBase. More... | |
Static Private Member Functions | |
template<class FunctionType > | |
static std::shared_ptr< ChebfunBase > | bestFitTempl (double start, double end, FunctionType f, std::vector< double > &p, std::vector< double > &a, double maxA, double tolerance, size_t maxSize) |
Templated implementation of bestFit method. More... | |
static bool | hasConverged (const std::vector< double > &a, double maxA, double tolerance, size_t shift=0) |
Test an array of Chebyshev coefficients for convergence. More... | |
Private Attributes | |
std::vector< double > | m_bw |
The barycentric weights. More... | |
double | m_end |
End of the interval. More... | |
std::vector< double > | m_integrationWeights |
The integration weights. More... | |
size_t | m_n |
Polynomial order. More... | |
double | m_start |
Start of the interval. More... | |
const double | m_tolerance |
Actual tolerance in comparing doubles. More... | |
std::vector< double > | m_x |
The x-points. More... | |
Static Private Attributes | |
static const size_t | g_maxNumberPoints = 1026 |
Maximum number of (x) points in a base. More... | |
static const double | g_tolerance = 1e-15 |
Maximum tolerance in comparing doubles. More... | |
The ChebfunBase class provides a base for function approximation with Chebyshev polynomials.
A smooth function on a finite interval [a,b] can be approximated by a Chebyshev expansion of order n. Finding an approximation is very easy: the function needs to be evaluated at n+1 specific x- points. These n+1 values can be used to interpolate the function at any x-point in interval [a,b]. This is done by calling the fit(...) method.
Different functions require different polynomial orders to reach the same accuracy of approximation. Static method bestFit(...) tries to find the smallest value of n that provides the required accuracy. If it fails to find an n smaller than some maximum number it returns an empty shared pointer.
Knowing the vector of the function values (P) at the n+1 base x-points and the related vector of the Chebyshev expansion coefficients (A) (claculated by calcA(...) method) allows one to perform various manipulations on the approximation:
This calss doesn't represent a function approximation itself but keeps proerties that can be shared by multiple approximations.
This class is based on the ideas from the Chebfun matlab package (http://www.chebfun.org/).
Definition at line 66 of file ChebfunBase.h.
Mantid::CurveFitting::Functions::ChebfunBase::ChebfunBase | ( | size_t | n, |
double | start, | ||
double | end, | ||
double | tolerance = 0.0 |
||
) |
Constructor.
n | :: Polynomial order == number of points - 1. |
start | :: The start (lower bound) of an interval on the x-axis. |
end | :: The end (upper bound) of an interval on the x-axis. |
tolerance | :: Tolerance in comparing the expansion coefficients. Setting this tolerance is a way to specify the accuracy of the approximation. |
Definition at line 51 of file ChebfunBase.cpp.
|
default |
Copy constructor.
|
defaultnoexcept |
Move constructor.
|
static |
Fit a function until full convergence.
Template specialization for a generic function type.
Definition at line 445 of file ChebfunBase.cpp.
References bestFitTempl(), and tolerance().
Referenced by bestFitAnyTolerance().
|
static |
Fit a function until full convergence.
Template specialization for IFunction.
Definition at line 451 of file ChebfunBase.cpp.
References tolerance().
|
static |
Find best fit with highest possible tolerance (to be used with noisy data).
Definition at line 178 of file ChebfunBase.h.
References bestFit(), g_tolerance, and tolerance().
Referenced by Mantid::CurveFitting::Algorithms::ChiSlice::makeApprox().
|
staticprivate |
Templated implementation of bestFit method.
Referenced by bestFit().
ChebfunBase_sptr Mantid::CurveFitting::Functions::ChebfunBase::bestFitTempl | ( | double | start, |
double | end, | ||
FunctionType | f, | ||
std::vector< double > & | p, | ||
std::vector< double > & | a, | ||
double | maxA, | ||
double | tolerance, | ||
size_t | maxSize | ||
) |
Fit a function until full convergence.
Increases size of the base until full conversion or a size limit is reached. If size limit is reached returns an empty pointer. In this case the calling method can divide the interval and fit each part separately.
start | :: Lower limit of the x-interval. |
end | :: Upper limit of the x-interval. |
f | :: Function to fit. |
p | :: Function values at the found x-points. |
a | :: Chebyshev expansion coefficients. |
maxA | :: A maximum value of the coefficients to compare against. |
tolerance | :: A tolerance for comparing values. Defines the accuracy of the approximation. |
maxSize | :: A maximum size of the expansion. If a large base is required to reach the desired accuracy returns an empty pointer. |
Definition at line 358 of file ChebfunBase.cpp.
References calcA(), fabs, fit(), fitOdd(), g_maxNumberPoints, g_tolerance, hasConverged(), Mantid::Geometry::m, n, tmp, and tolerance().
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::calcA | ( | const std::vector< double > & | p | ) | const |
Calculate expansion coefficients.
Calculate the chebyshev expansion coefficients given function values at the x-points.
p | :: Function values at chebyshev points. |
Definition at line 478 of file ChebfunBase.cpp.
References m_n, Mantid::CurveFitting::HalfComplex::real(), tmp, wavetable, and workspace.
Referenced by bestFitTempl(), and smooth().
|
private |
Calculate the integration weights.
Definition at line 124 of file ChebfunBase.cpp.
References m_end, m_integrationWeights, m_n, m_start, and n.
Referenced by integrate(), and integrationWeights().
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::calcP | ( | const std::vector< double > & | a | ) | const |
Calculate function values.
Calculate function values at chebyshev points given chebyshev expansion coefficiens (inverse of calcA()).
a | :: Chebyshev expansion coefficients. |
Definition at line 538 of file ChebfunBase.cpp.
References Mantid::Geometry::d, m_n, Mantid::CurveFitting::HalfComplex::set(), tmp, wavetable, and workspace.
Referenced by smooth().
|
private |
Calculate the x-values based on the (start,end) interval.
Definition at line 105 of file ChebfunBase.cpp.
References m_end, m_n, m_start, and m_x.
Referenced by ChebfunBase().
void Mantid::CurveFitting::Functions::ChebfunBase::derivative | ( | const std::vector< double > & | a, |
std::vector< double > & | aout | ||
) | const |
|
inline |
double Mantid::CurveFitting::Functions::ChebfunBase::eval | ( | double | x, |
const std::vector< double > & | p | ||
) | const |
Evaluate a function.
x | :: Point of evaluation. |
p | :: The function y-points |
Definition at line 192 of file ChebfunBase.cpp.
References m_bw, m_end, m_x, and Mantid::Geometry::x.
Referenced by roots().
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::evalVector | ( | const std::vector< double > & | x, |
const std::vector< double > & | p | ||
) | const |
Evaluate a function.
Evaluate function on a vector of x-values.
x | :: A vector of x-values. |
p | :: The y-points of a function. |
Definition at line 276 of file ChebfunBase.cpp.
References evalVector(), and Mantid::Geometry::x.
void Mantid::CurveFitting::Functions::ChebfunBase::evalVector | ( | const std::vector< double > & | x, |
const std::vector< double > & | p, | ||
std::vector< double > & | res | ||
) | const |
Evaluate a function.
Evaluate function on a vector of x-values.
x | :: A vector of x-values. |
p | :: The y-points of a function. |
res | :: Output result. res.size() == x.size() |
Definition at line 222 of file ChebfunBase.cpp.
References m_bw, m_start, m_x, value, and Mantid::Geometry::x.
Referenced by evalVector().
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::fit | ( | ChebfunFunctionType | f | ) | const |
Calculate function values at chebfun x-points.
Approximate a function using this base.
f | :: A function pointer. |
Definition at line 576 of file ChebfunBase.cpp.
Referenced by bestFitTempl().
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::fit | ( | const API::IFunction & | f | ) | const |
Calculate function values at chebfun x-points.
Approximate a function using this base.
f | :: A reference to an IFunction |
Definition at line 587 of file ChebfunBase.cpp.
References Mantid::API::IFunction::function(), m_x, Mantid::Geometry::x, and Mantid::Geometry::y.
|
private |
Calculate function values at odd-valued indices of the base x-points.
This method is used by bestFit to minimize the number of calls to the approximated function.
f | :: Function to calculate. |
pEven | :: Values of function f at the even-valued indices of m_x. |
Definition at line 625 of file ChebfunBase.cpp.
References Mantid::API::IFunction::function(), m_x, size(), Mantid::Geometry::x, and Mantid::Geometry::y.
|
private |
Calculate function values at odd-valued indices of the base x-points.
This method is used by bestFit to minimize the number of calls to the approximated function.
f | :: Function to calculate. |
p | :: Values of function f at the even-valued indices of m_x. |
Definition at line 601 of file ChebfunBase.cpp.
References size(), Mantid::Geometry::x, and xPoints().
Referenced by bestFitTempl().
|
staticprivate |
Test an array of Chebyshev coefficients for convergence.
Test if an array of Chebyshev expansion coefficients converged to the specified tolerance.
a | :: A vector of coefficients. |
maxA | :: A maximum value of the coefficients to compare against. |
tolerance | :: Convergence tolerance. |
shift | :: Displacement from the back of the coeffs vector (a) from which the coefficients are considered for testing. |
Definition at line 161 of file ChebfunBase.cpp.
References fabs, tmp, and tolerance().
Referenced by bestFitTempl().
ChebfunBase_sptr Mantid::CurveFitting::Functions::ChebfunBase::integral | ( | const std::vector< double > & | a, |
std::vector< double > & | aout | ||
) | const |
Calculate the integral.
Calculate the first integral of a function.
a | :: Chebyshev coefficients of the integrated function. |
aout | :: Output coeffs of the integral. |
Definition at line 320 of file ChebfunBase.cpp.
References Mantid::Geometry::d, m_end, m_n, m_start, and m_x.
double Mantid::CurveFitting::Functions::ChebfunBase::integrate | ( | const std::vector< double > & | p | ) | const |
Calculate an integral.
Calculate an integral of a function given its values at the base x-points.
p | :: Function values at the x-points. |
Definition at line 84 of file ChebfunBase.cpp.
References calcIntegrationWeights(), m_integrationWeights, m_x, and tmp.
const std::vector< double > & Mantid::CurveFitting::Functions::ChebfunBase::integrationWeights | ( | ) | const |
Get a reference to the integration weights.
Return the integration weights that can be used in function manipulations involving integration.
Definition at line 72 of file ChebfunBase.cpp.
References calcIntegrationWeights(), m_integrationWeights, and m_x.
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::linspace | ( | size_t | n | ) | const |
Create a vector of x values linearly spaced on the approximation interval.
Return a vector of linearly spaced values in the domain interval m_start <= x <= m_end.
n | :: Number of points in the output vector. |
Definition at line 461 of file ChebfunBase.cpp.
References m_end, m_start, n, width(), and Mantid::Geometry::x.
|
privatedelete |
|
privatedelete |
Private assingment operator to stress the immutability of ChebfunBase.
|
inline |
Get the polynomial order of this base.
Definition at line 74 of file ChebfunBase.h.
Referenced by roots().
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::roots | ( | const std::vector< double > & | a | ) | const |
Find all roots of a function on this interval.
Find all roots of this chebfun.
a | :: A vector with the Chebyshev expansion coefficients. |
Definition at line 656 of file ChebfunBase.cpp.
References endX(), eval(), fabs, Mantid::CurveFitting::getGSLMatrixView(), Mantid::CurveFitting::EigenMatrix::mutator(), order(), Mantid::CurveFitting::EigenMatrix::set(), startX(), tmp, Mantid::CurveFitting::EigenMatrix::tr(), workspace, Mantid::Geometry::x, and Mantid::CurveFitting::EigenMatrix::zero().
|
inline |
Get the size of the base which is the number of x-points.
Definition at line 76 of file ChebfunBase.h.
std::vector< double > Mantid::CurveFitting::Functions::ChebfunBase::smooth | ( | const std::vector< double > & | xvalues, |
const std::vector< double > & | yvalues | ||
) | const |
Smooth some data.
xvalues | :: X-values of the data to smooth. |
yvalues | :: Y-values of the data to smooth. xvalues.size() == yvalues.size() |
Definition at line 730 of file ChebfunBase.cpp.
References calcA(), calcP(), g_tolerance, m_x, n, sigma, size(), tmp, Mantid::Geometry::x, and Mantid::Geometry::y.
|
inline |
|
inline |
Tolerance for comparing doubles.
Definition at line 120 of file ChebfunBase.h.
Referenced by bestFit(), bestFitAnyTolerance(), bestFitTempl(), and hasConverged().
|
inline |
Get the width of the interval.
Definition at line 82 of file ChebfunBase.h.
Referenced by linspace().
|
inline |
Get a reference to the x-points.
Definition at line 84 of file ChebfunBase.h.
Referenced by fitOdd().
|
staticprivate |
Maximum number of (x) points in a base.
Definition at line 171 of file ChebfunBase.h.
Referenced by bestFitTempl().
|
staticprivate |
Maximum tolerance in comparing doubles.
Definition at line 169 of file ChebfunBase.h.
Referenced by bestFitAnyTolerance(), bestFitTempl(), and smooth().
|
private |
The barycentric weights.
Definition at line 165 of file ChebfunBase.h.
Referenced by ChebfunBase(), eval(), and evalVector().
|
private |
End of the interval.
Definition at line 161 of file ChebfunBase.h.
Referenced by calcIntegrationWeights(), calcX(), derivative(), eval(), integral(), and linspace().
|
mutableprivate |
The integration weights.
Definition at line 167 of file ChebfunBase.h.
Referenced by calcIntegrationWeights(), integrate(), and integrationWeights().
|
private |
Polynomial order.
Definition at line 157 of file ChebfunBase.h.
Referenced by calcA(), calcIntegrationWeights(), calcP(), calcX(), derivative(), and integral().
|
private |
Start of the interval.
Definition at line 159 of file ChebfunBase.h.
Referenced by calcIntegrationWeights(), calcX(), derivative(), evalVector(), integral(), and linspace().
|
private |
Actual tolerance in comparing doubles.
Definition at line 155 of file ChebfunBase.h.
|
private |
The x-points.
Definition at line 163 of file ChebfunBase.h.
Referenced by calcX(), ChebfunBase(), derivative(), eval(), evalVector(), fit(), fitOdd(), integral(), integrate(), integrationWeights(), and smooth().