Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
GSLInterpolator.h
Go to the documentation of this file.
1// @(#)root/mathmore:$Id$
2// Authors: L. Moneta, A. Zsenei 08/2005
3
4 /**********************************************************************
5 * *
6 * Copyright (c) 2004 ROOT Foundation, CERN/PH-SFT *
7 * *
8 * This library is free software; you can redistribute it and/or *
9 * modify it under the terms of the GNU General Public License *
10 * as published by the Free Software Foundation; either version 2 *
11 * of the License, or (at your option) any later version. *
12 * *
13 * This library is distributed in the hope that it will be useful, *
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
16 * General Public License for more details. *
17 * *
18 * You should have received a copy of the GNU General Public License *
19 * along with this library (see file COPYING); if not, write *
20 * to the Free Software Foundation, Inc., 59 Temple Place, Suite *
21 * 330, Boston, MA 02111-1307 USA, or contact the author. *
22 * *
23 **********************************************************************/
24
25// Header file for class GSLInterpolator
26//
27// Created by: moneta at Fri Nov 26 15:31:41 2004
28//
29// Last update: Fri Nov 26 15:31:41 2004
30//
31#ifndef ROOT_Math_GSLInterpolator
32#define ROOT_Math_GSLInterpolator
33
34#include <vector>
35#include <string>
36#include <cassert>
37
39
40#include "gsl/gsl_interp.h"
41#include "gsl/gsl_spline.h"
42
43#include "gsl/gsl_errno.h"
44#include "Math/Error.h"
45
46namespace ROOT {
47namespace Math {
48
49
50 /**
51 Interpolation class based on GSL interpolation functions
52 @ingroup Interpolation
53 */
54
56
57 public:
58
59 GSLInterpolator(unsigned int ndata, Interpolation::Type type);
60
61 GSLInterpolator(const Interpolation::Type type, const std::vector<double> & x, const std::vector<double> & y );
62 virtual ~GSLInterpolator();
63
64 // usually copying is non trivial, so we make delete this
69
70 bool Init(unsigned int ndata, const double *x, const double *y);
71
72 double Eval(double x) const
73 {
74 assert(fAccel);
75 double y = 0;
76 static unsigned int nErrors = 0;
77 int ierr = gsl_spline_eval_e(fSpline, x, fAccel, &y);
78
79 if (fResetNErrors)
80 nErrors = 0, fResetNErrors = false;
81
82 if (ierr) {
83 ++nErrors;
84 if(nErrors <= 4) {
85 MATH_WARN_MSG("GSLInterpolator::Eval", gsl_strerror(ierr));
86 if(nErrors == 4)
87 MATH_WARN_MSG("GSLInterpolator::Eval", "Suppressing additional warnings");
88 }
89 }
90
91 return y;
92 }
93
94 double Deriv(double x) const
95 {
96 assert(fAccel);
97 double deriv = 0;
98 static unsigned int nErrors = 0;
99 int ierr = gsl_spline_eval_deriv_e(fSpline, x, fAccel, &deriv);
100
101 if (fResetNErrors)
102 nErrors = 0, fResetNErrors = false;
103
104 if (ierr) {
105 ++nErrors;
106 if(nErrors <= 4) {
107 MATH_WARN_MSG("GSLInterpolator::Deriv", gsl_strerror(ierr));
108 if(nErrors == 4)
109 MATH_WARN_MSG("GSLInterpolator::Deriv", "Suppressing additional warnings");
110 }
111 }
112
113 return deriv;
114 }
115
116 double Deriv2(double x) const {
117 assert(fAccel);
118 double deriv2 = 0;
119 static unsigned int nErrors = 0;
120 int ierr = gsl_spline_eval_deriv2_e(fSpline, x, fAccel, &deriv2);
121
122 if (fResetNErrors)
123 nErrors = 0, fResetNErrors = false;
124
125 if (ierr) {
126 ++nErrors;
127 if(nErrors <= 4) {
128 MATH_WARN_MSG("GSLInterpolator::Deriv2", gsl_strerror(ierr));
129 if(nErrors == 4)
130 MATH_WARN_MSG("GSLInterpolator::Deriv2", "Suppressing additional warnings");
131 }
132 }
133
134 return deriv2;
135 }
136
137 double Integ(double a, double b) const {
138 if (a > b) // gsl will report an error in this case
139 return -Integ(b, a);
140
141 assert(fAccel);
142 double result = 0;
143 static unsigned int nErrors = 0;
144 int ierr = gsl_spline_eval_integ_e(fSpline, a, b, fAccel, &result);
145
146 if (fResetNErrors)
147 nErrors = 0, fResetNErrors = false;
148
149 if (ierr) {
150 ++nErrors;
151 if(nErrors <= 4) {
152 MATH_WARN_MSG("GSLInterpolator::Integ", gsl_strerror(ierr));
153 if(nErrors == 4)
154 MATH_WARN_MSG("GSLInterpolator::Integ", "Suppressing additional warnings");
155 }
156 }
157
158 return result;
159 }
160
161 std::string Name() {
162 return fInterpType->name;
163 }
164
165 private:
166
167 mutable bool fResetNErrors; // flag to reset counter for error messages
168 gsl_interp_accel * fAccel;
169 gsl_spline * fSpline;
170 const gsl_interp_type * fInterpType;
171
172 };
173
174} // namespace Math
175} // namespace ROOT
176
177#endif /* ROOT_Math_GSLInterpolator */
#define MATH_WARN_MSG(loc, str)
Definition Error.h:80
#define b(i)
Definition RSha256.hxx:100
#define a(i)
Definition RSha256.hxx:99
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
Interpolation class based on GSL interpolation functions.
GSLInterpolator & operator=(const GSLInterpolator &)=delete
GSLInterpolator(GSLInterpolator &&)=delete
double Deriv2(double x) const
double Eval(double x) const
bool Init(unsigned int ndata, const double *x, const double *y)
const gsl_interp_type * fInterpType
double Integ(double a, double b) const
GSLInterpolator(const Interpolation::Type type, const std::vector< double > &x, const std::vector< double > &y)
double Deriv(double x) const
GSLInterpolator(const GSLInterpolator &)=delete
Type
Enumeration defining the types of interpolation methods availables.
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
Namespace for new Math classes and functions.
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...