Spline.h

Go to the documentation of this file.

#if !defined(LSST_AFW_MATH_DETAIL_SPLINE)
#define LSST_AFW_MATH_DETAIL_SPLINE 1
#include <cmath>
#include <vector>

namespace lsst { namespace afw { namespace math { namespace detail {

/*****************************************************************************/
/*
 * Splines
 */
class Spline {
public:
    virtual ~Spline() {}

    void interpolate(std::vector<double> const& x, 
                     std::vector<double> &y        
                    ) const;
    void derivative(std::vector<double> const& x, 
                    std::vector<double> &dydx     
                   ) const;

    std::vector<double> roots(double const value, 
                              double const x0,    
                              double const x1     
                             ) const;
    
protected:
    Spline() {}
    void _allocateSpline(int const nknot);

    std::vector<double> _knots;                // positions of knots
    std::vector<std::vector<double> > _coeffs; // and associated coefficients
};

class TautSpline : public Spline {
public:
    enum Symmetry { Unknown, Odd, Even };

    TautSpline(std::vector<double> const& x,
               std::vector<double> const& y,
               double const gamma=0,
               Symmetry type=Unknown
          );
private:
    void calculateTautSpline(std::vector<double> const& x,
                             std::vector<double> const& y,
                             double const gamma0
                            );
    void calculateTautSplineEvenOdd(std::vector<double> const& x,
                                    std::vector<double> const& y,
                                    double const gamma0,
                                    bool even
                                );
};

class SmoothedSpline : public Spline {
public:
    SmoothedSpline(std::vector<double> const& x,  
                   std::vector<double> const& y,  
                   std::vector<double> const& dy, 
                   double s,                      
                   double *chisq=NULL,            
                   std::vector<double> *errs=NULL 
          );
};
}}}}
#endif