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

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

Fixes #466. GSL_error is now thrown on GSL errors. Last thing in this ticket.

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