00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #ifndef ROOT_Math_GSLInterpolator
00032 #define ROOT_Math_GSLInterpolator
00033
00034 #include <vector>
00035 #include <string>
00036 #include <cassert>
00037
00038 #include "Math/InterpolationTypes.h"
00039
00040 #include "gsl/gsl_interp.h"
00041 #include "gsl/gsl_spline.h"
00042
00043 #include "gsl/gsl_errno.h"
00044 #include "Math/Error.h"
00045
00046 namespace ROOT {
00047 namespace Math {
00048
00049
00050
00051
00052
00053
00054
00055 class GSLInterpolator {
00056
00057 public:
00058
00059 GSLInterpolator(unsigned int ndata, Interpolation::Type type);
00060
00061 GSLInterpolator(const Interpolation::Type type, const std::vector<double> & x, const std::vector<double> & y );
00062 virtual ~GSLInterpolator();
00063
00064 private:
00065
00066 GSLInterpolator(const GSLInterpolator &);
00067 GSLInterpolator & operator = (const GSLInterpolator &);
00068
00069 public:
00070
00071 bool Init(unsigned int ndata, const double *x, const double * y);
00072
00073 double Eval( double x ) const
00074 {
00075 assert(fAccel);
00076 double y = 0;
00077 int ierr = gsl_spline_eval_e(fSpline, x, fAccel, &y );
00078 if (ierr) MATH_WARN_MSG("GSLInterpolator::Eval",gsl_strerror(ierr) )
00079 return y;
00080 }
00081
00082 double Deriv( double x ) const
00083 {
00084 assert(fAccel);
00085 double deriv = 0;
00086 int ierr = gsl_spline_eval_deriv_e(fSpline, x, fAccel, &deriv );
00087 if (ierr) MATH_WARN_MSG("GSLInterpolator::Deriv",gsl_strerror(ierr) )
00088 return deriv;
00089 }
00090
00091 double Deriv2( double x ) const {
00092 assert(fAccel);
00093 double deriv2 = 0;
00094 int ierr = gsl_spline_eval_deriv2_e(fSpline, x, fAccel, &deriv2 );
00095 if (ierr) MATH_WARN_MSG("GSLInterpolator::Deriv2",gsl_strerror(ierr) )
00096 return deriv2;
00097 }
00098
00099 double Integ( double a, double b) const {
00100 if ( a > b) return -Integ(b,a);
00101 assert(fAccel);
00102 double result = 0;
00103 int ierr = gsl_spline_eval_integ_e(fSpline, a, b, fAccel, &result );
00104 if (ierr) MATH_WARN_MSG("GSLInterpolator::Integ",gsl_strerror(ierr) )
00105 return result;
00106 }
00107
00108 std::string Name() {
00109 return fInterpType->name;
00110 }
00111
00112
00113 protected:
00114
00115
00116 private:
00117
00118 gsl_interp_accel * fAccel;
00119 gsl_spline * fSpline;
00120 const gsl_interp_type * fInterpType;
00121
00122 };
00123
00124 }
00125 }
00126
00127
00128 #endif