interpm_idw.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2018, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_INTERPM_IDW_H
24 #define O2SCL_INTERPM_IDW_H
25 
26 /** \file interpm_idw.h
27  \brief File defining \ref o2scl::interpm_idw
28 */
29 
30 #include <iostream>
31 #include <string>
32 #include <cmath>
33 
34 #include <boost/numeric/ublas/matrix.hpp>
35 
36 #include <gsl/gsl_combination.h>
37 
38 #include <o2scl/err_hnd.h>
39 #include <o2scl/vector.h>
40 #include <o2scl/vec_stats.h>
41 #include <o2scl/linear_solver.h>
42 #include <o2scl/columnify.h>
43 
44 #ifndef DOXYGEN_NO_O2NS
45 namespace o2scl {
46 #endif
47 
48  /** \brief Multi-dimensional interpolation by inverse distance
49  weighting
50 
51  This class performs interpolation on a multi-dimensional data
52  set specified as a series of scattered points using the inverse
53  distance-weighted average of nearby points. The function \ref
54  set_data() takes as input: the number of input dimensions, the
55  number of output functions, the number of points which specify
56  the data, and a "vector of vectors" which contains the data for
57  all the points. The vector of vectors must be of a type which
58  allows std::swap on individual elements (which are of type
59  <tt>vec_t</tt>).
60 
61  The "order" of the interpolation, i.e. the number of nearby
62  points which are averaged, defaults to 3 and can be changed
63  by \ref set_order(). To obtain interpolation uncertainties,
64  this class finds the nearest <tt>order+1</tt> points and
65  returns the standard deviation of the interpolated value
66  over all of the subsets of <tt>order</tt> points.
67  The value <tt>order=1</tt> corresponds to nearest-neighbor
68  interpolation.
69 
70  This class requires a distance metric to weight the
71  interpolation, and a Euclidean distance is used. By default, the
72  length scales in each direction are automatically determined by
73  extent of the data (absolute value of max minus min in each
74  direction), but the user can specify the length scales manually
75  in \ref set_scales() .
76 
77  First derivatives can be obtained using \ref derivs_err() , but
78  these derivatives are not obtained from the same approximation
79  used in the interpolation. That is, the derivatives returned are
80  not equal to exact derivatives from the interpolated function
81  (as is the case in, e.g., cubic spline interpolation in one
82  dimension). This will typically only be particularly noticable
83  near discontinuities.
84 
85  If the user specifies an array of pointers, the data can be
86  changed between calls to the interpolation, but data points
87  cannot be added (as set data separately stores the total number
88  of data points) without a new call to \ref set_data(). Also, the
89  automatically-determined length scales may need to be recomputed
90  by calling \ref auto_scale().
91  */
92  template<class vec_t> class interpm_idw {
93 
94  protected:
95 
96  public:
97 
101 
102  interpm_idw() {
103  data_set=false;
104  scales.resize(1);
105  scales[0]=1.0;
106  order=3;
107  verbose=0;
108  }
109 
110  /** \brief Verbosity parameter (default 0)
111  */
112  int verbose;
113 
114  /** \brief Set the number of closest points to use
115  for each interpolation (default 3)
116  */
117  void set_order(size_t n) {
118  if (n==0) {
119  O2SCL_ERR("Order cannot be zero in interpm_idw.",
121  }
122  order=n;
123  return;
124  }
125 
126  /** \brief Set the scales for the distance metric
127 
128  All the scales must be positive and non-zero. The size of the
129  vector \c (specified in \c n) must be larger than zero.
130  */
131  template<class vec2_t> void set_scales(size_t n, vec2_t &v) {
132  if (n==0) {
133  O2SCL_ERR("Scale vector size cannot be zero in interpm_idw.",
135  }
136  for(size_t i=0;i<n;i++) {
137  if (v[i]<=0.0) {
138  O2SCL_ERR("Scale must be positive and non-zero in interpm_idw.",
140  }
141  }
142  scales.resize(n);
144  return;
145  }
146 
147  /** \brief Initialize the data for the interpolation
148 
149  The object \c vecs should be a vector (of size <tt>n_in+n_out</tt>)
150  of vectors (all of size <tt>n_points</tt>). It may have be
151  any time which allows the use of <tt>std::swap</tt> for
152  each vector in the list.
153  */
154  template<class vec_vec_t>
155  void set_data(size_t n_in, size_t n_out, size_t n_points,
156  vec_vec_t &vecs, bool auto_scale_flag=true) {
157 
158  if (n_points<3) {
159  O2SCL_ERR2("Must provide at least three points in ",
160  "interpm_idw::set_data()",exc_efailed);
161  }
162  if (n_in<1) {
163  O2SCL_ERR2("Must provide at least one input column in ",
164  "interpm_idw::set_data()",exc_efailed);
165  }
166  if (n_out<1) {
167  O2SCL_ERR2("Must provide at least one output column in ",
168  "interpm_idw::set_data()",exc_efailed);
169  }
170  np=n_points;
171  nd_in=n_in;
172  nd_out=n_out;
173  ptrs.resize(n_in+n_out);
174  for(size_t i=0;i<n_in+n_out;i++) {
175  std::swap(ptrs[i],vecs[i]);
176  }
177  data_set=true;
178 
179  if (auto_scale_flag) {
180  auto_scale();
181  }
182 
183  return;
184  }
185 
186  /** \brief Get the data used for interpolation
187  */
188  template<class vec_vec_t>
189  void get_data(size_t &n_in, size_t &n_out, size_t &n_points,
190  vec_vec_t &vecs) {
191  n_points=np;
192  n_in=nd_in;
193  n_out=nd_out;
194  vecs.resize(n_in+n_out);
195  for(size_t i=0;i<n_in+n_out;i++) {
196  std::swap(ptrs[i],vecs[i]);
197  }
198  data_set=false;
199  n_points=0;
200  n_in=0;
201  n_out=0;
202  ptrs.clear();
203  return;
204  }
205 
206  /** \brief Automatically determine the length scales from the
207  data
208  */
209  void auto_scale() {
210  scales.resize(nd_in);
211  for(size_t i=0;i<nd_in;i++) {
212  scales[i]=fabs(o2scl::vector_max_value<vec_t,double>(np,ptrs[i])-
213  o2scl::vector_min_value<vec_t,double>(np,ptrs[i]));
214  }
215  return;
216  }
217 
218  /** \brief Initialize the data for the interpolation
219  for only one output function
220 
221  The object \c vecs should be a vector (of size <tt>n_in+1</tt>)
222  of vectors (all of size <tt>n_points</tt>). It may be
223  any type which allows the use of <tt>std::swap</tt> for
224  each vector in the list.
225  */
226  template<class vec_vec_t>
227  void set_data(size_t n_in, size_t n_points,
228  vec_vec_t &vecs, bool auto_scale=true) {
229  set_data(n_in,1,n_points,vecs,auto_scale);
230  return;
231  }
232 
233  /** \brief Perform the interpolation over the first function
234  */
235  template<class vec2_t> double operator()(const vec2_t &x) const {
236  return eval(x);
237  }
238 
239  /** \brief Perform the interpolation over the first function
240  */
241  template<class vec2_t> double eval(const vec2_t &x) const {
242 
243  if (data_set==false) {
244  O2SCL_ERR("Data not set in interpm_idw::eval().",
245  exc_einval);
246  }
247 
248  // Compute distances
249  std::vector<double> dists(np);
250  for(size_t i=0;i<np;i++) {
251  dists[i]=dist(i,x);
252  }
253 
254  // Find closest points
255  std::vector<size_t> index;
256  o2scl::vector_smallest_index<std::vector<double>,double,
257  std::vector<size_t> >(dists,order,index);
258 
259  // Check if the closest distance is zero
260  if (dists[index[0]]<=0.0) {
261  return ptrs[nd_in][index[0]];
262  }
263 
264  // Compute normalization
265  double norm=0.0;
266  for(size_t i=0;i<order;i++) {
267  norm+=1.0/dists[index[i]];
268  }
269 
270  // Compute the inverse-distance weighted average
271  double ret=0.0;
272  for(size_t i=0;i<order;i++) {
273  ret+=ptrs[nd_in][index[i]]/dists[index[i]];
274  }
275  ret/=norm;
276 
277  // Return the average
278  return ret;
279  }
280 
281  /** \brief Perform the interpolation over the first function
282  with uncertainty
283  */
284  template<class vec2_t> void eval_err(const vec2_t &x, double &val,
285  double &err) const {
286 
287  if (data_set==false) {
288  O2SCL_ERR("Data not set in interpm_idw::eval_err().",
289  exc_einval);
290  }
291 
292  // Compute distances
293  std::vector<double> dists(np);
294  for(size_t i=0;i<np;i++) {
295  dists[i]=dist(i,x);
296  }
297 
298  // Find closest points
299  std::vector<size_t> index;
300  o2scl::vector_smallest_index<std::vector<double>,double,
301  std::vector<size_t> >(dists,order+1,index);
302 
303  if (dists[index[0]]<=0.0) {
304 
305  // If the closest distance is zero, just set the value
306  val=ptrs[nd_in][index[0]];
307  err=0.0;
308  return;
309 
310  } else {
311 
312  std::vector<double> vals(order+1);
313 
314  for(size_t j=0;j<order+1;j++) {
315 
316  // Compute normalization
317  double norm=0.0;
318  for(size_t i=0;i<order+1;i++) {
319  if (i!=j) norm+=1.0/dists[index[i]];
320  }
321 
322  // Compute the inverse-distance weighted average
323  vals[j]=0.0;
324  for(size_t i=0;i<order+1;i++) {
325  if (i!=j) {
326  vals[j]+=ptrs[nd_in][index[i]]/dists[index[i]];
327  }
328  }
329  vals[j]/=norm;
330 
331  }
332 
333  val=vals[order];
334  err=o2scl::vector_stddev(vals);
335 
336  }
337 
338  return;
339  }
340 
341  /** \brief Perform the interpolation over all the functions,
342  storing the result in \c y
343  */
344  template<class vec2_t, class vec3_t>
345  void eval(vec2_t &x, vec3_t &y) const {
346 
347  if (data_set==false) {
348  O2SCL_ERR("Data not set in interpm_idw::eval().",
349  exc_einval);
350  }
351 
352  if (verbose>0) {
353  std::cout << "interpm_idw: input: ";
354  for(size_t k=0;k<nd_in;k++) {
355  std::cout << x[k] << " ";
356  }
357  std::cout << std::endl;
358  }
359 
360  // Compute distances
361  std::vector<double> dists(np);
362  for(size_t i=0;i<np;i++) {
363  dists[i]=dist(i,x);
364  }
365 
366  // Find closest points
367  std::vector<size_t> index;
368  o2scl::vector_smallest_index<std::vector<double>,double,
369  std::vector<size_t> >(dists,order,index);
370  if (verbose>0) {
371  for(size_t i=0;i<order;i++) {
372  std::cout << "interpm_idw: closest point: ";
373  for(size_t k=0;k<nd_in;k++) {
374  std::cout << ptrs[k][index[i]] << " ";
375  }
376  std::cout << std::endl;
377  }
378  }
379 
380  // Check if the closest distance is zero, if so, just
381  // return the value
382  if (dists[index[0]]<=0.0) {
383  for(size_t i=0;i<nd_out;i++) {
384  y[i]=ptrs[nd_in+i][index[0]];
385  }
386  if (verbose>0) {
387  std::cout << "interpm_idw: distance zero. "
388  << "Returning values at index: " << index[0] << std::endl;
389  std::cout << "\t";
390  o2scl::vector_out(std::cout,nd_out,y,true);
391  }
392  return;
393  }
394 
395  // Compute normalization
396  double norm=0.0;
397  for(size_t i=0;i<order;i++) {
398  norm+=1.0/dists[index[i]];
399  }
400  if (verbose>0) {
401  std::cout << "interpm_idw: norm is " << norm << std::endl;
402  }
403 
404  // Compute the inverse-distance weighted averages
405  for(size_t j=0;j<nd_out;j++) {
406  y[j]=0.0;
407  for(size_t i=0;i<order;i++) {
408  if (j==0 && verbose>0) {
409  std::cout << "interpm_idw: Point: ";
410  for(size_t k=0;k<nd_in;k++) {
411  std::cout << ptrs[k][index[i]] << " ";
412  }
413  std::cout << std::endl;
414  }
415  y[j]+=ptrs[nd_in+j][index[i]]/dists[index[i]];
416  if (verbose>0) {
417  std::cout << "interpm_idw: j,order,value,1/dist: "
418  << j << " " << i << " "
419  << ptrs[nd_in+j][index[i]] << " "
420  << 1.0/dists[index[i]] << std::endl;
421  }
422  }
423  y[j]/=norm;
424  if (verbose>0) {
425  std::cout << "interpm_idw: y[" << j << "]: " << y[j]
426  << std::endl;
427  }
428  }
429 
430 
431 
432  return;
433  }
434 
435  /** \brief Perform the interpolation over all the functions
436  with uncertainties
437  */
438  template<class vec2_t, class vec3_t, class vec4_t>
439  void eval_err_index(const vec2_t &x, vec3_t &val, vec4_t &err,
440  std::vector<size_t> &index) const {
441 
442  if (data_set==false) {
443  O2SCL_ERR("Data not set in interpm_idw::eval_err().",
444  exc_einval);
445  }
446 
447  // Compute distances
448  std::vector<double> dists(np);
449  for(size_t i=0;i<np;i++) {
450  dists[i]=dist(i,x);
451  }
452 
453  // Find closest points
454  o2scl::vector_smallest_index<std::vector<double>,double,
455  std::vector<size_t> >(dists,order+1,index);
456 
457  if (dists[index[0]]<=0.0) {
458 
459  // If the closest distance is zero, just set the values and
460  // errors
461  for(size_t k=0;k<nd_out;k++) {
462  val[k]=ptrs[nd_in+k][index[0]];
463  err[k]=0.0;
464  }
465  return;
466 
467  } else {
468 
469  for(size_t k=0;k<nd_out;k++) {
470 
471  std::vector<double> vals(order+1);
472 
473  for(size_t j=0;j<order+1;j++) {
474 
475  // Compute normalization
476  double norm=0.0;
477  for(size_t i=0;i<order+1;i++) {
478  if (i!=j) norm+=1.0/dists[index[i]];
479  }
480 
481  // Compute the inverse-distance weighted average
482  vals[j]=0.0;
483  for(size_t i=0;i<order+1;i++) {
484  if (i!=j) {
485  vals[j]+=ptrs[nd_in+k][index[i]]/dists[index[i]];
486  }
487  }
488  vals[j]/=norm;
489 
490  }
491 
492  // Instead of using the average, we report the
493  // value as the last element in the array, which
494  // is the interpolated value from the closest points
495  val[k]=vals[order];
496 
497  err[k]=o2scl::vector_stddev(vals);
498 
499  }
500 
501  }
502 
503  return;
504  }
505 
506  /** \brief Perform the interpolation over all the functions
507  with uncertainties
508  */
509  template<class vec2_t, class vec3_t, class vec4_t>
510  void eval_err(const vec2_t &x, vec3_t &val, vec4_t &err) const {
511  std::vector<size_t> index;
512  return eval_err_index(x,val,err,index);
513  }
514 
515  /** \brief For one of the functions, compute the partial
516  derivatives (and uncertainties) with respect to all of the
517  inputs at one data point
518 
519  \note This function ignores the order chosen by \ref
520  set_order() and always chooses to average derivative
521  calculations determined from \c n_in+1 combinations of \c n_in
522  points .
523 
524  \future This function requires an extra copy from
525  "ders" to "ders2" which could be removed.
526  */
527  template<class vec3_t>
528  void derivs_err(size_t func_index, size_t point_index,
529  vec3_t &derivs, vec3_t &errs) const {
530 
531  if (data_set==false) {
532  O2SCL_ERR("Data not set in interpm_idw::derivs_err().",
533  exc_einval);
534  }
535 
536  // Set x equal to the specified point
537  ubvector x(nd_in);
538  for(size_t i=0;i<nd_in;i++) {
539  x[i]=ptrs[i][point_index];
540  }
541  // Set f equal to the value of the function at the specified point
542  double f=ptrs[nd_in+func_index][point_index];
543 
544  // The linear solver
546 
547  // Compute distances
548  std::vector<double> dists(np);
549  for(size_t i=0;i<np;i++) {
550  dists[i]=dist(i,x);
551  }
552 
553  // Find closest (but not identical) points
554 
555  std::vector<size_t> index;
556  size_t max_smallest=(nd_in+2)*2;
557  if (max_smallest>np) max_smallest=np;
558  if (max_smallest<nd_in+1) {
559  O2SCL_ERR("Couldn't find enough nearby points.",o2scl::exc_einval);
560  }
561 
562  if (verbose>0) {
563  std::cout << "max_smallest: " << max_smallest << std::endl;
564  }
565 
566  o2scl::vector_smallest_index<std::vector<double>,double,
567  std::vector<size_t> >(dists,max_smallest,index);
568 
569  if (verbose>0) {
570  for(size_t i=0;i<index.size();i++) {
571  std::cout << "index[" << i << "] = " << index[i] << " "
572  << dists[index[i]] << std::endl;
573  }
574  }
575 
576  std::vector<size_t> index2;
577  for(size_t i=0;i<max_smallest;i++) {
578  if (dists[index[i]]>0.0) {
579  index2.push_back(index[i]);
580  if (index2.size()==nd_in+1) i=max_smallest;
581  }
582  }
583  if (index2.size()<nd_in+1) {
584  O2SCL_ERR("Couldn't find enough nearby points (2).",
586  }
587 
588  if (verbose>0) {
589  for(size_t i=0;i<index2.size();i++) {
590  std::cout << "index2[" << i << "] = " << index2[i] << " "
591  << dists[index2[i]] << std::endl;
592  }
593  }
594 
595  // Unit vector storage
596  std::vector<ubvector> units(nd_in+1);
597  // Difference vector norms
598  std::vector<double> diff_norms(nd_in+1);
599  // Storage for the derivative estimates
600  std::vector<ubvector> ders(nd_in+1);
601  // Matrix of unit vectors
602  ubmatrix m(nd_in,nd_in);
603  // Vector of function value differences
604  ubvector v(nd_in);
605  // Rearranged derivative object
606  std::vector<ubvector> ders2(nd_in);
607 
608  for(size_t i=0;i<nd_in+1;i++) {
609 
610  // Assign unit vector elements
611  units[i].resize(nd_in);
612  for(size_t j=0;j<nd_in;j++) {
613  units[i][j]=ptrs[j][index2[i]]-x[j];
614  }
615 
616  // Normalize the unit vectors
617  diff_norms[i]=o2scl::vector_norm<ubvector,double>(units[i]);
618  for(size_t j=0;j<nd_in;j++) {
619  units[i][j]/=diff_norms[i];
620  }
621 
622  }
623 
624  // Verbose output of the closest points and their norms
625  if (verbose>0) {
626  std::cout << "Point: ";
627  for(size_t i=0;i<nd_in;i++) {
628  std::cout << x[i] << " ";
629  }
630  std::cout << f << std::endl;
631  for(size_t j=0;j<nd_in+1;j++) {
632  std::cout << "Closest: " << j << " " << index2[j] << " ";
633  for(size_t i=0;i<nd_in;i++) {
634  std::cout << ptrs[i][index2[j]] << " ";
635  }
636  std::cout << ptrs[func_index+nd_in][index2[j]] << " "
637  << diff_norms[j] << std::endl;
638  }
639  for(size_t j=0;j<nd_in+1;j++) {
640  std::cout << "diff_norm: " << j << " " << diff_norms[j]
641  << std::endl;
642  }
643  // End of verbose output
644  }
645 
646  // Go through each set of points
647  for(size_t i=0;i<nd_in+1;i++) {
648 
649  ders[i].resize(nd_in);
650 
651  // Construct the matrix and vector for the solver
652  size_t jj=0;
653  for(size_t j=0;j<nd_in+1;j++) {
654  if (j!=i) {
655  for(size_t k=0;k<nd_in;k++) {
656  m(jj,k)=units[j][k];
657  }
658  v[jj]=(ptrs[func_index+nd_in][index2[j]]-f)/diff_norms[j];
659  jj++;
660  }
661  }
662 
663  // Solve to compute the derivatives
664  if (verbose>0) {
665  std::cout << "m:" << std::endl;
666  o2scl::matrix_out(std::cout,nd_in,nd_in,m);
667  std::cout << "v:" << std::endl;
668  o2scl::vector_out(std::cout,nd_in,v,true);
669  }
670  lshh.solve(nd_in,m,v,ders[i]);
671  if (verbose>0) {
672  std::cout << "Derivs: " << i << " ";
673  std::cout.setf(std::ios::showpos);
674  for(size_t j=0;j<nd_in;j++) {
675  std::cout << ders[i][j] << " ";
676  }
677  std::cout.unsetf(std::ios::showpos);
678  std::cout << std::endl;
679  }
680 
681  // Go to next derivative estimate
682  }
683 
684  for(size_t i=0;i<nd_in;i++) {
685 
686  // Rearrange derivatives
687  ders2[i].resize(nd_in+1);
688  for(size_t j=0;j<nd_in+1;j++) {
689  ders2[i][j]=ders[j][i];
690  }
691 
692  // Compute mean and standard deviation
693  derivs[i]=o2scl::vector_mean(ders2[i]);
694  errs[i]=o2scl::vector_stddev(ders2[i]);
695  }
696 
697  return;
698  }
699 
700 #ifndef DOXYGEN_INTERNAL
701 
702  protected:
703 
704  /// Distance scales for each coordinate
705  ubvector scales;
706 
707  /// The number of points
708  size_t np;
709  /// The number of dimensions of the inputs
710  size_t nd_in;
711  /// The number of dimensions of the outputs
712  size_t nd_out;
713  /// A vector of pointers holding the data
714  std::vector<vec_t> ptrs;
715  /// True if the data has been specified
716  bool data_set;
717  /// Number of points to include in each interpolation (default 3)
718  size_t order;
719 
720  /// Compute the distance between \c x and the point at index \c index
721  template<class vec2_t> double dist(size_t index,
722  const vec2_t &x) const {
723  double ret=0.0;
724  size_t nscales=scales.size();
725  for(size_t i=0;i<nd_in;i++) {
726  ret+=pow((x[i]-ptrs[i][index])/scales[i%nscales],2.0);
727  }
728  return sqrt(ret);
729  }
730 
731 #endif
732 
733  };
734 
735 #ifndef DOXYGEN_NO_O2NS
736 }
737 #endif
738 
739 #endif
740 
741 
742 
void derivs_err(size_t func_index, size_t point_index, vec3_t &derivs, vec3_t &errs) const
For one of the functions, compute the partial derivatives (and uncertainties) with respect to all of ...
Definition: interpm_idw.h:528
double vector_mean(size_t n, const vec_t &data)
Compute the mean of the first n elements of a vector.
Definition: vec_stats.h:55
Multi-dimensional interpolation by inverse distance weighting.
Definition: interpm_idw.h:92
The main O<span style=&#39;position: relative; top: 0.3em; font-size: 0.8em&#39;>2</span>scl O$_2$scl names...
Definition: anneal.h:42
invalid argument supplied by user
Definition: err_hnd.h:59
void get_data(size_t &n_in, size_t &n_out, size_t &n_points, vec_vec_t &vecs)
Get the data used for interpolation.
Definition: interpm_idw.h:189
virtual void solve(size_t n, mat_t &A, vec_t &b, vec_t &x)
Solve square linear system of size n.
Generic Householder linear solver.
int matrix_out(std::ostream &os, size_t nrows, size_t ncols, mat_t &A)
A operator for simple matrix output using operator()
Definition: columnify.h:254
generic failure
Definition: err_hnd.h:61
size_t nd_out
The number of dimensions of the outputs.
Definition: interpm_idw.h:712
void set_order(size_t n)
Set the number of closest points to use for each interpolation (default 3)
Definition: interpm_idw.h:117
double operator()(const vec2_t &x) const
Perform the interpolation over the first function.
Definition: interpm_idw.h:235
void eval_err_index(const vec2_t &x, vec3_t &val, vec4_t &err, std::vector< size_t > &index) const
Perform the interpolation over all the functions with uncertainties.
Definition: interpm_idw.h:439
double vector_stddev(size_t n, const vec_t &data)
Standard deviation with specified mean.
Definition: vec_stats.h:253
void auto_scale()
Automatically determine the length scales from the data.
Definition: interpm_idw.h:209
ubvector scales
Distance scales for each coordinate.
Definition: interpm_idw.h:705
void vector_copy(const vec_t &src, vec2_t &dest)
Simple vector copy.
Definition: vector.h:127
std::vector< vec_t > ptrs
A vector of pointers holding the data.
Definition: interpm_idw.h:714
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
Definition: err_hnd.h:281
double dist(size_t index, const vec2_t &x) const
Compute the distance between x and the point at index index.
Definition: interpm_idw.h:721
size_t order
Number of points to include in each interpolation (default 3)
Definition: interpm_idw.h:718
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Definition: err_hnd.h:273
void eval_err(const vec2_t &x, vec3_t &val, vec4_t &err) const
Perform the interpolation over all the functions with uncertainties.
Definition: interpm_idw.h:510
size_t np
The number of points.
Definition: interpm_idw.h:708
bool data_set
True if the data has been specified.
Definition: interpm_idw.h:716
void set_data(size_t n_in, size_t n_points, vec_vec_t &vecs, bool auto_scale=true)
Initialize the data for the interpolation for only one output function.
Definition: interpm_idw.h:227
void vector_out(std::ostream &os, size_t n, const vec_t &v, bool endline=false)
Output the first n elements of a vector to a stream, os.
Definition: vector.h:2528
size_t nd_in
The number of dimensions of the inputs.
Definition: interpm_idw.h:710
int verbose
Verbosity parameter (default 0)
Definition: interpm_idw.h:112
void set_scales(size_t n, vec2_t &v)
Set the scales for the distance metric.
Definition: interpm_idw.h:131
void set_data(size_t n_in, size_t n_out, size_t n_points, vec_vec_t &vecs, bool auto_scale_flag=true)
Initialize the data for the interpolation.
Definition: interpm_idw.h:155
double eval(const vec2_t &x) const
Perform the interpolation over the first function.
Definition: interpm_idw.h:241
void eval_err(const vec2_t &x, double &val, double &err) const
Perform the interpolation over the first function with uncertainty.
Definition: interpm_idw.h:284
void eval(vec2_t &x, vec3_t &y) const
Perform the interpolation over all the functions, storing the result in y.
Definition: interpm_idw.h:345

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).