Go to the documentation of this file.
23 #ifndef O2SCL_INTERP2_PLANAR_H
24 #define O2SCL_INTERP2_PLANAR_H
34 #include <gsl/gsl_combination.h>
36 #include <o2scl/err_hnd.h>
37 #include <o2scl/vector.h>
39 #ifndef DOXYGEN_NO_O2NS
136 void set_data(
size_t n_points, vec_t &x, vec_t &y, vec_t &f) {
138 O2SCL_ERR2(
"Must provide at least three points in ",
155 double minx=(*ux)[0], maxx=(*ux)[0];
156 for(
size_t i=1;i<
np;i++) {
157 if ((*
ux)[i]<minx) minx=(*ux)[i];
158 if ((*
ux)[i]>maxx) maxx=(*ux)[i];
165 double miny=(*uy)[0], maxy=(*uy)[0];
166 for(
size_t i=1;i<
np;i++) {
167 if ((*
uy)[i]<miny) miny=(*uy)[i];
168 if ((*
uy)[i]>maxy) maxy=(*uy)[i];
175 if (
dx<=0.0 ||
dy<=0.0) {
184 double eval(
double x,
double y)
const {
185 double x1,
x2,
x3, y1, y2, y3, f;
187 eval_points(x,y,f,i1,
x1,y1,i2,
x2,y2,i3,
x3,y3);
201 return eval(v[0],v[1]);
211 size_t &i1,
double &
x1,
double &y1,
212 size_t &i2,
double &
x2,
double &y2,
213 size_t &i3,
double &
x3,
double &y3)
const {
216 O2SCL_ERR(
"Data not set in interp2_planar::eval_points().",
225 double c1=sqrt(pow((x-(*
ux)[0])/
dx,2.0)+pow((y-(*
uy)[0])/
dy,2.0));
226 double c2=sqrt(pow((x-(*
ux)[1])/
dx,2.0)+pow((y-(*
uy)[1])/
dy,2.0));
227 double c3=sqrt(pow((x-(*
ux)[2])/
dx,2.0)+pow((y-(*
uy)[2])/
dy,2.0));
255 for(
size_t j=3;j<
np;j++) {
257 double c4=sqrt(pow((x-(*
ux)[i4])/
dx,2.0)+pow((y-(*
uy)[i4])/
dy,2.0));
271 double a, b, c, den, z1, z2, z3;
290 for(
size_t i=0;i<
np;i++) {
291 dist[i]=sqrt(pow((x-(*
ux)[i])/
dx,2.0)+pow((y-(*
uy)[i])/
dy,2.0));
296 gsl_combination *cc=gsl_combination_alloc(
np,3);
298 gsl_combination_init_first(cc);
301 while (done==
false) {
302 int status=gsl_combination_next(cc);
305 i1=order[gsl_combination_get(cc,0)];
306 i2=order[gsl_combination_get(cc,1)];
307 i3=order[gsl_combination_get(cc,2)];
330 gsl_combination_free(cc);
341 a=-(-y2*z1+y3*z1+y1*z2-y3*z2-y1*z3+y2*z3)/den;
343 c=-(
x3*y2*z1-
x2*y3*z1-
x3*y1*z2+
x1*y3*z2+
x2*y1*z3-
x1*y2*z3)/den;
350 #ifndef DOXYGEN_INTERNAL
366 int swap(
size_t &index_1,
double &dist_1,
size_t &index_2,
367 double &dist_2)
const {
372 index_temp=index_1; dist_temp=dist_1;
373 index_1=index_2; dist_1=dist_2;
374 index_2=index_temp; dist_2=dist_temp;
388 #ifndef DOXYGEN_NO_O2NS
double x_scale
The user-specified x scale (default -1)
size_t np
The number of points.
Interpolate among two independent variables with planes.
@ exc_efailed
generic failure
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
double y_scale
The user-specified y scale (default -1)
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.
double operator()(double x, double y) const
Perform the planar interpolation.
static const double x2[5]
double dx
The scale in the x direction.
void vector_sort_index(size_t n, const vec_t &data, vec_size_t &order)
Create a permutation which sorts the first n elements of a vector (in increasing order)
static const double x1[5]
@ exc_einval
invalid argument supplied by user
double operator()(vec2_t &v) const
Perform the planar interpolation using the first two elements of v as input.
double eval(double x, double y) const
Perform the planar interpolation.
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.
double thresh
Threshold for colinearity (default )
void compute_scale()
Find scaling.
bool data_set
True if the data has been specified.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
double dy
The scale in the y direction.
int swap(size_t &index_1, double &dist_1, size_t &index_2, double &dist_2) const
Swap points 1 and 2.
static const double x3[11]
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).