TMDlib is hosted by Hepforge, IPPP Durham
Interpolation2D Class Reference

#include <TMDlib_InterpolationKS.h>

Inheritance diagram for Interpolation2D:

Public Member Functions

 Interpolation2D (const vector< double > &x1vv, const vector< double > &x2vv, boost::multi_array< double, 2 > ym, gsl_interp_type int_type=*gsl_interp_linear)
 ~Interpolation2D ()
 destructor More...
vector< pair< double, double > > grid_limits () const
double interp (const vector< double > &x)
 return value from 2d spline interpolation at (x1,x2) More...
- Public Member Functions inherited from BaseInterpolation
 BaseInterpolation ()
 default constructor More...
virtual ~BaseInterpolation ()
 default destructor More...
unsigned int eval_status ()

Detailed Description

2-dim interpolation - linear or cubic spline

The class uses 1d gsl interpolation function and depends on boost multi_array so both of those libraries need to be linked.

The evalaution is done with the function interp(const vector<double> &x) and the status of the last evaluation can be checked with eval_status().

Example of usage:

boost::multi_array<double, 2> matrix_A(boost::extents[n][m]);

for (unsigned i=0;i<n;i++) {         
  for (unsigned j=0;j<m;j++) {      
      vec_x1[i] = arg1;
      vec_x2[j] = arg2;
      matrix_A[i][j] = val;

Interpolation2D int2d(vec_x1, vec_x2,  matrix_A, *gsl_interp_cspline);

Constructor & Destructor Documentation

Interpolation2D::Interpolation2D ( const vector< double > &  x1vv,
const vector< double > &  x2vv,
boost::multi_array< double, 2 >  ym,
gsl_interp_type  int_type = *gsl_interp_linear 

constructor - read the input and interpolate in the innermost dimension x2

x1vv,x2vvinput vectors with argument data; entries in both vectors are required to grow monotonically
yminput array with values data
int_typetype of interpolation, it takes the values values *gsl_interp_linear or *gsl_interp_cspline
Interpolation2D::~Interpolation2D ( )


Member Function Documentation

vector< pair< double, double > > Interpolation2D::grid_limits ( ) const

return grid limits in form of 2d vector of pairs: [[x1min,x1max], [x2min, x2max]]

Min and max values in each dimension can be obtained by using first or second pair data members.

Reimplemented from BaseInterpolation.

double Interpolation2D::interp ( const vector< double > &  x)

return value from 2d spline interpolation at (x1,x2)

Reimplemented from BaseInterpolation.

The documentation for this class was generated from the following files: