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

Last change on this file since 1655 was 1655, checked in by Jari Häkkinen, 13 years ago

Addresses #466. Removed bsearch from the GSL wrapper.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.4 KB
Line 
1#ifndef _theplu_yat_regression_gslinterpolation_
2#define _theplu_yat_regression_gslinterpolation_
3
4// $Id: GSLInterpolation.h 1655 2008-12-16 23:29:55Z jari $
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 Check the error status
66
67       Some underlying GSL functions return GSL error codes. The error
68       code is stored internally and should be checked for by the
69       caller. Refer to the gsl_errno.h for the error code listing.
70
71       \return The current error status. A zero return values
72       indicates no errors where encountered.
73    */
74    int error_status(void) const;
75
76    /**
77       \brief Calculate the interpolated value for \a x.
78
79       GSL error status is set if evaluation is requested outside the
80       range defined by the interpolation algorithm. Use function
81       error_status to check for errors.
82
83       \return The interpolated value of \f$ y \f$ for a given point
84       \a x.
85
86       \see error_status(void)
87    */
88    double evaluate(double x);
89
90    /**
91       \brief Calculate the derivative of the interpolated function at
92       \a x.
93
94       GSL error status is set if evaluation is requested outside the
95       range defined by the interpolation algorithm. Use function
96       error_status to check for errors.
97
98       \return The derivative.
99
100       \see error_status(void)
101    */
102    double evaluate_derivative(double x);
103
104    /**
105       \brief Calculate the 2nd derivative of the interpolated
106       function at \a x.
107
108       GSL error status is set if evaluation is requested outside the
109       range defined by the interpolation algorithm. Use function
110       error_status to check for errors.
111
112       \return The 2nd derivative.
113
114       \see error_status(void)
115    */
116    double evaluate_derivative2(double x);
117
118    /**
119       \brief Calculate the numerical integral of the interpolated
120       function over the range \f$ [a,b] \f$.
121
122       GSL error status is set if evaluation is requested outside the
123       range defined by the interpolation algorithm. Use function
124       error_status to check for errors.
125
126       \return The integral.
127
128       \see error_status(void)
129    */
130    double evaluate_integral(double a, double b);
131
132    /**
133       \brief The result of the latest evaluaion function call is
134       stored and can be retrieved with this function.
135
136       \return The latest evaluated value.
137     */
138    double evaluation(void) const;
139
140    /**
141       \brief This function returns the minimum number of points
142       required by the interpolation type.
143
144       For example, Akima spline interpolation requires a minimum of 5
145       points.
146
147       \return The minimum number of points required.
148    */
149    unsigned int min_size(void) const;
150
151  protected:
152    /**
153       \brief The default constructor
154
155       Initialization of the interpolation object for the data \f$ (x,
156       y) \f$ where \a x and \a y are vector like objects of the same
157       size. The content of \a x and \a y are copied for internal
158       storage. \a x is always assumed to be strictly ordered, with
159       increasing \a x values; the behavior for other arrangements is
160       not defined.
161
162       Some underlying GSL functions return GSL error codes. The error
163       code is stored internally and should be checked for by the
164       caller. Refer to the gsl_errno.h for the error code listing.
165
166       \see error_status(void)
167    */
168    GSLInterpolation(const gsl_interp_type*, const utility::VectorBase& x,
169                     const utility::VectorBase& y);
170
171    /**
172       \brief The destructor
173    */
174    virtual ~GSLInterpolation(void);
175
176  private:
177    /**
178       \brief Copy Constructor. (not implemented)
179    */
180    GSLInterpolation(const GSLInterpolation&);
181
182    gsl_interp_accel* accelerator_;
183    double evaluation_;
184    gsl_interp* interpolator_;
185    int gsl_status_;
186    double* x_;
187    double* y_;
188  };
189
190}}} // of namespaces regression, yat, and theplu
191
192#endif
Note: See TracBrowser for help on using the repository browser.