function f = function::interpolate_linear("filename",function(_1,_2),_3);

vector x(10), y(10);
// fill x and y here
function lin = function::interpolate_linear(x,y);

function sp = function::spline("filename",_2, _5+_6);

delaunay_interpolator d;
d.add_point(x,y,z); // many times
cerr<<"Value at (1.2, 3.4) = "<<d(1.2,3.4)<<endl;
// create a blop-function from this interpolator:
function f(d);
cerr<<"Value at (1.2, 3.4) = "<<f(1.2,3.4)<<endl;
// read from file:
function fff = function::interpolate_delaunay("filename",function(_1,_2),_3);

Interpolation in blop is intended to be as generic as possible. There are template interpolator classes being able to interpolate in 1-dimension and also in multi-dimensions, for any types for which addition/subtraction and multiplication/division-by-a-number is valid. Currently linear, spline, and 'shape-preserving piecewise cubic hermite polynomial' (sppchip) are implemented in 1D, linear and delaunay interpolations are implemented in 2D, and linear interpolation is implemented in multi-dimensions.

You can use the low-level interpolation classes directly, or use interpolated blop-functions which use the same base classes.

• Interpolating functions
• The interpolator class
• The 1D interpolator base class (common member functions of 1D interpolators)
• Polynomial interpolators in 1D
  (linear, spline, etc)

Use the following functions (static member functions of the function class) to create blop-functions that interpolate discrete data provided either in a stream (file, pipe, etc) or in arrays. XXX should be replaced by linear, spline or sppchip

function::interpolate_XXX(const var &filename, const function &x=_1, const function &y=_2);
function::interpolate_XXX(const array &x, const array &y);
function::interpolate_XXX(const vector &x, const vector &y);
function::interpolate_XXX(const double x[], const double y[], int n);

An example for their usage:

function f = function::interpolate_linear("filename.dat",function(_1,_2),_3);

This creates a 2D linear interpolator, (x,y) values (independent variables) are the first two columns of the file, the interpolated value is the 3rd column. For spline and sppchip interpolation only 1D is possible (i.e. the independent variable must be a single value, unlike in this example).

Delaunay interpolation can be used to interpolate z-values on a set of irregular (x,y) points. This is therefore a 2D interpolation. The (x,y) plane is Delaunay-triangulated, and a triangulated surface is created having nodes at the z-(x,y) points.

function::interpolate_delaunay(const var &filename);   // (x,y)=(1st,2nd)column, z=3rd column
function::interpolate_delaunay(const var &filename, const function &xy, const function &z);
     // exapmle for the second: function f = function::interpolate_delaunay(filename,function(_2,_3),_1);

If you don't want to have the overhead from the usage of blop-functions, or if you want to programmatically fill up your interpolator with the discrete values, or you want to have access to these discrete values later, you can also use the interpolator classes directly.

---

content below this line should be checked

---

The base class interpolator provides some common interface functions to be used by blop classes (i.e. no interest for the user), and a name() member function returning a name of the interpolation method (for example LINEAR, SPLINE, etc).

This is the base class for 1D interpolation. This class does nothing more than storing the (x,y) point pairs used for the interpolation. Derived classes will implement the different interpolation methods.

Member functions:

points(const vector &x, const vector &y)

points(T x, T y, int n) [T is an indexeable template like double*]

provide the discrete points for the interpolation

add_point(double x, double y)

Append an (x,y) point pair to the end of the stored points. It will be improved to make an ordered addition, but currently it is just appended to the end, the user should make sure the increasing order in x is kept

double operator()(double x)

Evaluate the interpolated value at x

int n() const

Return the number of points

double x(int i) const
double y(int i) const

Return the ith x/y value

interpolator_1d &x(int i, double xval)
interpolator_1d &y(int i, double yval)

Set the ith x or y value (if i is a non-existent index, the vectors storing the datapoints are automatically expanded, missing elements are initialized to 0. This can screw up the orderedness in the x array, it is the user's responsibility that the final array (upon completion of setting up all the points) is ordered

This derived class, polynomial_interpolator_1d interpolates within each interval [x[i];x[i+1]] using a polynomial of some order, the coefficients of which are stored in this function. Methods include linear, spline, shape-preserving piecewise cubic hermite, etc. interpolations. To make an interpolated blop-function, use these functions (where XXX can be linear, spline or sppchip (shape-preserving piecewise cubic hermite interpolating polynomial, see http://www.mathworks.com/moler/interp.pdf for details) More details about the interpolator classes are given below.

static function function::XXX(const var &filename, const function &f1=unset,... const function &f4=unset);	
    // read discrete values from the given file, possibly transform
    // the values, for example first column is x, 2nd+3rd column is y
    // function linint = function::linear("filename",_1,_2+_3);
static function function::XXX(const array &x, const array &y);
static function function::XXX(const vector &x, const vector &y);
static function function::XXX(double x[], double y[], int n);

The coefficients of the interpolation are calculated by calling one of the functions linear(), spline(), etc. Note that you need to call one of these functions every time you change the discrete points, or add new ones.

Member functions (besides those inherited from interpolator_1d, see above)

int order() const

Return the order of the polynomial

bool linear()

Calculate interpolation coefficients to make linear interpolation. Returns false if too few (i.e. <2) points are given

bool spline()

Calculate interpolation coefficients to make spline interpolation. Returns false if too few (<3) points are given

bool sppchip()

Calculate interpolation coefficients to make a shape-preserving piecewise cubic hermite interpolating polynomial (huh...), see http://www.mathworks.com/moler/interp.pdf

Example:

vector<double> x(10), y(10);
x[0] = 1;
y[0] = 1.2;
.... etc, fill the vectors x,y
polynomial_interpolator_1d interp(x,y);
interp.spline();
cerr<<"Value at 2.3 = "<<interp(2.3)<<endl;

---

 Source files:
   interpolate.h
   interpolate.cc

---