Class interp2_planar (o2scl)

O2scl : Class List

template<class vec_t>
class interp2_planar

Interpolate among two independent variables with planes.

This class performs planar interpolation when the data points are not arranged in a specified order (i.e. not on a grid). For a set of data \( {x_i,y_i,f_i} \), the value of \( f \) is predicted given a new value of x and y. This interpolation is performed by finding the plane that goes through three closest points in the data set. Distances are determined with

\[ d_{ij} = \sqrt{\left(\frac{x_i-x_j}{\Delta x}\right)^2 + \left(\frac{y_i-y_j}{\Delta y}\right)^2} \]
The values \( \Delta_x \) and \( \Delta_y \) are specified in x_scale and y_scale, respectively. If these values are negative (the default) then they are computed with \( \Delta x = x_{\mathrm{max}}-x_{\mathrm{min}} \) and \( \Delta y = y_{\mathrm{max}}-y_{\mathrm{min}} \) .

If the x- and y-values of the entire data set lie on a line, then the interpolation will fail and the error handler will be called. Colinearity is defined by a threshold, thresh which defaults to \( 10^{-12} \). If the denominator,

\[ \sum_{k=1}^{3}\varepsilon_{ijk} x_i y_j < \mathrm{thresh} \]
where \( \varepsilon_{ijk} \) is an anti-symmetric Levi-Cevita tensor, then the points are colinear. The value of thresh can be zero, but if it is negative then it will be reset to the default value for the next interpolation.

This class stores pointers to the data, not a copy. The data can be changed between interpolations without an additional call to set_data(), but the scales may need to be recomputed with compute_scale().

The vector type can be any type with a suitably defined operator[].

The interpolation requires at least three points and set_data() will call the error handler if the first argument is less than three.

Idea for Future:

Make a parent class for this and o2scl::interp2_neigh.

Note

This class is experimental.

Note

This class operates by performing a \( {\cal O}(N) \) brute-force search to find the three closest points. If the three closest points are colinear, then the data are sorted by distance [ \( {\cal O}(N \log N) \) ], and the closest triplets are enumerated until a non-colinear triplet is found.

Note

I believe this interpolation is a bit unstable because it doesn’t ensure that the user-specified objective point is inside the region determined by the three closest data points, and this can lead to some strong extrapolation.

Public Types

typedef boost::numeric::ublas::vector<double> ubvector
typedef boost::numeric::ublas::vector<size_t> ubvector_size_t

Public Functions

inline interp2_planar()
inline void set_data(size_t n_points, vec_t &x, vec_t &y, vec_t &f)

Initialize the data for the planar interpolation and compute the scaling factors.

inline void compute_scale()

Find scaling.

inline double eval(double x, double y) const

Perform the planar interpolation.

inline double operator()(double x, double y) const

Perform the planar interpolation.

template<class vec2_t>
inline double operator()(vec2_t &v) const

Perform the planar interpolation using the first two elements of v as input.

inline void eval_points(double x, double y, double &f, size_t &i1, double &x1, double &y1, size_t &i2, double &x2, double &y2, size_t &i3, double &x3, double &y3) const

Planar interpolation returning the closest points.

This function interpolates x and y into the data returning f. It also returns the three closest x- and y-values used for computing the plane.

Public Members

double thresh

Threshold for colinearity (default \( 10^{-12} \))

double x_scale

The user-specified x scale (default -1)

double y_scale

The user-specified y scale (default -1)

Protected Functions

inline int swap(size_t &index_1, double &dist_1, size_t &index_2, double &dist_2) const

Swap points 1 and 2.

Protected Attributes

double dx

The scale in the x direction.

double dy

The scale in the y direction.

size_t np

The number of points.

vec_t *ux

The x-values.

vec_t *uy

The y-values.

vec_t *uf

The f-values.

bool data_set

True if the data has been specified.

Private Functions

interp2_planar(const interp2_planar<vec_t>&)
interp2_planar<vec_t> &operator=(const interp2_planar<vec_t>&)