source: trunk/yat/regression/GSLInterpolation.h

Last change on this file was 4207, checked in by Peter, 6 weeks ago

update copyright statements

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