71 lines
2.2 KiB
Plaintext
71 lines
2.2 KiB
Plaintext
[/
|
|
Copyright (c) 2019 Nick Thompson
|
|
Use, modification and distribution are subject to the
|
|
Boost Software License, Version 1.0. (See accompanying file
|
|
LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
|
|
]
|
|
|
|
[section:whittaker_shannon Whittaker-Shannon interpolation]
|
|
|
|
[heading Synopsis]
|
|
``
|
|
#include <boost/math/interpolators/whittaker_shannon.hpp>
|
|
``
|
|
|
|
namespace boost { namespace math { namespace interpolators {
|
|
|
|
template <class RandomAccessContainer>
|
|
class whittaker_shannon
|
|
{
|
|
public:
|
|
|
|
using Real = RandomAccessContainer::value_type;
|
|
|
|
whittaker_shannon(RandomAccessContainer&& v, Real left_endpoint, Real step_size);
|
|
|
|
Real operator()(Real x) const;
|
|
|
|
Real prime(Real x) const;
|
|
};
|
|
|
|
}}} // namespaces
|
|
|
|
|
|
[heading Whittaker-Shannon Interpolation]
|
|
|
|
The Whittaker-Shannon interpolator takes equispaced data and interpolates between them via a sum of sinc functions.
|
|
This interpolation is stable and infinitely smooth, but has linear complexity in the data,
|
|
making it slow relative to compactly-supported b-splines.
|
|
In addition, we cannot pass an infinite amount of data into the class,
|
|
and must truncate the (perhaps) infinite sinc series to a finite number of terms.
|
|
Since the sinc function has slow /1/x/ decay, the truncation of the series can incur large error.
|
|
Hence this interpolator works best when operating on samples of compactly-supported functions.
|
|
Here is an example of interpolating a smooth "bump function":
|
|
|
|
auto bump = [](double x) { if (std::abs(x) >= 1) { return 0.0; } return std::exp(-1.0/(1.0-x*x)); };
|
|
|
|
double t0 = -1;
|
|
size_t n = 2049;
|
|
double h = 2.0/(n-1.0);
|
|
|
|
std::vector<double> v(n);
|
|
for(size_t i = 0; i < n; ++i) {
|
|
double t = t0 + i*h;
|
|
v[i] = bump(t);
|
|
}
|
|
|
|
auto ws = whittaker_shannon(std::move(v), t0, h);
|
|
|
|
double y = ws(0.3);
|
|
|
|
The derivative of the interpolant can also be evaluated, but the accuracy is not as high:
|
|
|
|
double yp = ws.prime(0.3);
|
|
|
|
[heading Complexity and Performance]
|
|
|
|
The call to the constructor requires [bigo](1) operations, simply moving data into the class.
|
|
Each call to the interpolant is [bigo](/n/), where /n/ is the number of points to interpolate.
|
|
|
|
[endsect] [/section:whittaker_shannon]
|