source: trunk/yat/regression/GSLInterpolation.h @ 1676

Last change on this file since 1676 was 1676, checked in by Peter, 12 years ago

fixes #473 - declaring a private assignement operator in interpolation classes

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.5 KB
Line 
1#ifndef _theplu_yat_regression_gslinterpolation_
2#define _theplu_yat_regression_gslinterpolation_
3
4// $Id: GSLInterpolation.h 1676 2008-12-24 01:43:04Z peter $
5
6/*
7  Copyright (C) 2008 Jari Häkkinen
8
9  This file is part of the yat library, http://dev.thep.lu.se/yat
10
11  The yat library is free software; you can redistribute it and/or
12  modify it under the terms of the GNU General Public License as
13  published by the Free Software Foundation; either version 3 of the
14  License, or (at your option) any later version.
15
16  The yat library is distributed in the hope that it will be useful,
17  but WITHOUT ANY WARRANTY; without even the implied warranty of
18  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19  General Public License for more details.
20
21  You should have received a copy of the GNU General Public License
22  along with yat. If not, see <http://www.gnu.org/licenses/>.
23*/
24
25#include <gsl/gsl_interp.h>
26
27namespace theplu {
28namespace yat {
29namespace utility {
30  class VectorBase;
31}
32namespace regression {
33
34  /**
35     \brief Base class for interfacing GSL interpolation.
36
37     The GSL interpolation is described in the
38     <A HREF=
39     "http://www.gnu.org/software/gsl/manual/html_node/Interpolation.html">
40     GSL Manual.</A>
41     The GSL library provides a variety of interpolation methods,
42     including Cubic splines and Akima splines. Interpolations can be
43     defined for both normal and periodic boundary
44     conditions. Additional functions are available for computing
45     derivatives and integrals of interpolating functions.
46
47     Given a set of data points \f$ (x_1, y_1) \dots (x_n, y_n) \f$
48     the sub classes compute a continuous interpolating function \f$
49     y(x) \f$ such that \f$ y(x_i) = y_i \f$. The interpolation is
50     piecewise smooth, and its behavior at the end-points is
51     determined by the type of interpolation used.
52
53     Some underlying GSL functions return GSL error codes. The error
54     code is stored internally and should be checked for by the caller
55     using GSLInterpolation::error_status(void). Refer to the
56     gsl_errno.h for the error code listing.
57
58     \since New in yat 0.5
59  */
60  class GSLInterpolation
61  {
62
63  public:
64    /**
65       \brief The destructor
66    */
67    virtual ~GSLInterpolation(void)=0;
68
69    /**
70       \brief Check the error status
71
72       Some underlying GSL functions return GSL error codes. The error
73       code is stored internally and should be checked for by the
74       caller. Refer to the gsl_errno.h for the error code listing.
75
76       \return The current error status. A zero return values
77       indicates no errors where encountered.
78    */
79    int error_status(void) const;
80
81    /**
82       \brief Calculate the interpolated value for \a x.
83
84       GSL error status is set if evaluation is requested outside the
85       range defined by the interpolation algorithm. Use function
86       error_status to check for errors.
87
88       \return The interpolated value of \f$ y \f$ for a given point
89       \a x.
90
91       \see error_status(void)
92    */
93    double evaluate(double x);
94
95    /**
96       \brief Calculate the derivative of the interpolated function at
97       \a x.
98
99       GSL error status is set if evaluation is requested outside the
100       range defined by the interpolation algorithm. Use function
101       error_status to check for errors.
102
103       \return The derivative.
104
105       \see error_status(void)
106    */
107    double evaluate_derivative(double x);
108
109    /**
110       \brief Calculate the 2nd derivative of the interpolated
111       function at \a x.
112
113       GSL error status is set if evaluation is requested outside the
114       range defined by the interpolation algorithm. Use function
115       error_status to check for errors.
116
117       \return The 2nd derivative.
118
119       \see error_status(void)
120    */
121    double evaluate_derivative2(double x);
122
123    /**
124       \brief Calculate the numerical integral of the interpolated
125       function over the range \f$ [a,b] \f$.
126
127       GSL error status is set if evaluation is requested outside the
128       range defined by the interpolation algorithm. Use function
129       error_status to check for errors.
130
131       \return The integral.
132
133       \see error_status(void)
134    */
135    double evaluate_integral(double a, double b);
136
137    /**
138       \brief The result of the latest evaluaion function call is
139       stored and can be retrieved with this function.
140
141       \return The latest evaluated value.
142     */
143    double evaluation(void) const;
144
145    /**
146       \brief This function returns the minimum number of points
147       required by the interpolation type.
148
149       For example, Akima spline interpolation requires a minimum of 5
150       points.
151
152       \return The minimum number of points required.
153    */
154    unsigned int min_size(void) const;
155
156  protected:
157    /**
158       \brief The default constructor
159
160       Initialization of the interpolation object for the data \f$ (x,
161       y) \f$ where \a x and \a y are vector like objects of the same
162       size. The content of \a x and \a y are copied for internal
163       storage. \a x is always assumed to be strictly ordered, with
164       increasing \a x values; the behavior for other arrangements is
165       not defined.
166
167       Some underlying GSL functions return GSL error codes. The error
168       code is stored internally and should be checked for by the
169       caller. Refer to the gsl_errno.h for the error code listing.
170
171       \see error_status(void)
172    */
173    GSLInterpolation(const gsl_interp_type*, const utility::VectorBase& x,
174                     const utility::VectorBase& y);
175
176  private:
177    /**
178       \brief Copy Constructor. (not implemented)
179    */
180    GSLInterpolation(const GSLInterpolation&);
181    // no assignment
182    GSLInterpolation& operator=(const GSLInterpolation&);
183
184    gsl_interp_accel* accelerator_;
185    double evaluation_;
186    gsl_interp* interpolator_;
187    int gsl_status_;
188    double* x_;
189    double* y_;
190  };
191
192}}} // of namespaces regression, yat, and theplu
193
194#endif
Note: See TracBrowser for help on using the repository browser.