table3d.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2020, 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_TABLE3D_H
24 #define O2SCL_TABLE3D_H
25 
26 /** \file table3d.h
27  \brief File defining \ref o2scl::table3d
28 */
29 
30 #include <iostream>
31 #include <fstream>
32 #include <vector>
33 #include <string>
34 #include <cmath>
35 #include <sstream>
36 
37 #include <boost/numeric/ublas/vector.hpp>
38 #include <boost/numeric/ublas/vector_proxy.hpp>
39 #include <boost/numeric/ublas/matrix.hpp>
40 #include <boost/numeric/ublas/matrix_proxy.hpp>
41 
42 #include <o2scl/misc.h>
43 #include <o2scl/err_hnd.h>
44 #include <o2scl/search_vec.h>
45 #include <o2scl/uniform_grid.h>
46 #include <o2scl/interp.h>
47 #include <o2scl/table_units.h>
48 #include <o2scl/contour.h>
49 #include <o2scl/shunting_yard.h>
50 #include <o2scl/vector.h>
51 
52 // Forward definition of the table3d class for HDF I/O
53 namespace o2scl {
54  class table3d;
55 }
56 
57 // Forward definition of HDF I/O to extend friendship
58 namespace o2scl_hdf {
59  class hdf_file;
60  void hdf_input(hdf_file &hf, o2scl::table3d &t, std::string name);
61  void hdf_output(hdf_file &hf, o2scl::table3d &t, std::string name);
62 }
63 
64 #ifndef DOXYGEN_NO_O2NS
65 namespace o2scl {
66 #endif
67 
68  /** \brief A data structure containing one or more slices of
69  two-dimensional data points defined on a grid
70 
71  \future Improve interpolation and derivative caching, possibly
72  through non-const versions of the interpolation functions.
73  \future Should there be a clear_grid() function separate from
74  clear_data() and clear()?
75  \future Allow the user to more clearly probe 'size_set' vs.
76  'xy_set'? (AWS 07/18: This is apparently resolved.)
77  */
78  class table3d {
79 
80  public:
81 
84 
85  // This is used for the interpolation classes
86  typedef boost::numeric::ublas::matrix_row<const ubmatrix> ubmatrix_row;
87  typedef boost::numeric::ublas::matrix_column<const ubmatrix>
88  ubmatrix_column;
89 
90  /** \brief Create a new 3D \table
91  */
92  table3d();
93 
94  virtual ~table3d();
95 
96  /** \brief Create a table3d object from a table, assuming \c scolx
97  and \c scoly store the x- and y-grid data, respectively.
98  */
99  table3d(o2scl::table_units<> &t, std::string colx, std::string coly);
100 
101  /// Copy constructor
102  table3d(const table3d &t);
103 
104  /// Copy constructor
105  table3d &operator=(const table3d &t);
106 
107  /// \name Initialization
108  //@{
109  /** \brief Initialize the x-y grid
110 
111  This function will not allow you to redefine the grid when
112  there is data in the \table if a grid of a different size was
113  already set from a previous call to either set_xy() or
114  set_size(). However, you may freely redefine the grid after a
115  call to clear_data() or clear_table(). You may change
116  individual grid points at any time with set_grid_x() and
117  set_grid_y().
118  */
119  template<class vec_t, class vec2_t>
120  void set_xy(std::string x_name, size_t nx, const vec_t &x,
121  std::string y_name, size_t ny, const vec2_t &y) {
122 
123  if (has_slice && (size_set || xy_set) && (nx!=numx || ny!=numy)) {
124  O2SCL_ERR("Size cannot be reset in table3d::set_xy().",
126  return;
127  }
128  if (xy_set) {
129  xval.clear();
130  yval.clear();
131  }
132  numx=nx;
133  numy=ny;
134  xname=x_name;
135  yname=y_name;
136  xval.resize(nx);
137  yval.resize(ny);
138  for(size_t i=0;i<nx;i++) (xval)[i]=x[i];
139  for(size_t i=0;i<ny;i++) (yval)[i]=y[i];
140  size_set=true;
141  xy_set=true;
142  return;
143  }
144 
145  /** \brief Initialize the x-y grid with \ref uniform_grid objects
146 
147  This function will not allow you to redefine the grid when
148  there is data in the \table if a grid of a different size was
149  already set from a previous call to either set_xy() or
150  set_size(). However, you may freely redefine the grid after a
151  call to clear_data() or clear_table(). You may change
152  individual grid points at any time with set_grid_x() and
153  set_grid_y().
154  */
155  void set_xy(std::string x_name, uniform_grid<double> gx,
156  std::string y_name, uniform_grid<double> gy);
157 
158  /** \brief Initialize \table size
159 
160  This function will not allow you to resize the \table if it
161  already has data or if the size has already been set with the
162  set_xy() function, unless you clear the data with clear_data()
163  or the \table with clear_table() first.
164  */
165  void set_size(size_t nx, size_t ny);
166  //@}
167 
168  // --------------------------------------------------------
169  /// \name On-grid get and set methods
170  //@{
171 
172  /** \brief Set element in slice \c name at location <tt>ix,iy</tt>
173  to value \c val
174  */
175  void set(size_t ix, size_t iy, std::string name, double val);
176 
177  /** \brief Set element in slice of index \c z at location
178  <tt>ix,iy</tt> to value \c val .
179  */
180  void set(size_t ix, size_t iy, size_t z, double val);
181 
182  /** \brief Get element in slice \c name at location
183  <tt>ix,iy</tt>
184  */
185  double &get(size_t ix, size_t iy, std::string name);
186 
187  /** \brief Get element in slice \c name at location
188  <tt>ix,iy</tt> (const version)
189  */
190  const double &get(size_t ix, size_t iy, std::string name) const;
191 
192  /** \brief Get element in slice of index \c z at location
193  <tt>ix,iy</tt>
194  */
195  double &get(size_t ix, size_t iy, size_t z);
196 
197  /** \brief Get element in slice of index \c z at location
198  <tt>ix,iy</tt> (const version)
199  */
200  const double &get(size_t ix, size_t iy, size_t z) const;
201  //@}
202 
203  // --------------------------------------------------------
204  /** \name Off-grid get and set methods
205 
206  These methods return the value of a slice on the grid
207  point nearest to a user-specified location. For
208  interpolation into a point off the grid, use
209  \ref table3d::interp().
210  */
211  //@{
212 
213  /** \brief Set element in slice \c name at the nearest location to
214  <tt>x,y</tt> to value \c val
215  */
216  void set_val(double x, double y, std::string name, double val);
217 
218  /** \brief Set element in slice of index \c z at the nearest
219  location to <tt>x,y</tt> to value \c val
220  */
221  void set_val(double x, double y, size_t z, double val);
222 
223  /** \brief Get element in slice \c name at location closest to
224  <tt>x,y</tt>
225  */
226  double &get_val(double x, double y, std::string name);
227 
228  /** \brief Get element in slice \c name at location closest to
229  <tt>x,y</tt>
230  */
231  const double &get_val(double x, double y, std::string name) const;
232 
233  /** \brief Get element in slice of index \c z at location closest
234  to <tt>x,y</tt>
235  */
236  double &get_val(double x, double y, size_t z);
237 
238  /** \brief Get element in slice of index \c z at location closest
239  to <tt>x,y</tt>
240  */
241  const double &get_val(double x, double y, size_t z) const;
242 
243  /** \brief Set elements in the first <tt>nv</tt> slices at the
244  nearest location to <tt>x,y</tt> to value \c val
245  */
246  template<class vec_t>
247  void set_slices(double x, double y, size_t nv, vec_t &vals) {
248  size_t ix, iy;
249  lookup_x(x,ix);
250  lookup_y(y,iy);
251 
252  for(size_t i=0;i<nv && i<list.size();i++) {
253  list[i](ix,iy)=vals[i];
254  }
255  return;
256  }
257 
258  /** \brief Get the data for every slice at the nearest location to
259  <tt>x,y</tt>
260  */
261  template<class vec_t>
262  void get_slices(double x, double y, size_t nv, vec_t &v) {
263 
264  size_t ix, iy;
265 
266  lookup_x(x,ix);
267  lookup_y(y,iy);
268 
269  for(size_t i=0;i<nv && i<list.size();i++) {
270  v[i]=list[i](ix,iy);
271  }
272 
273  return;
274  }
275  //@}
276 
277  // --------------------------------------------------------
278  /// \name Off-grid get and set methods returning nearest point
279  //@{
280 
281  /** \brief Set element in slice \c name at the nearest location to
282  <tt>x,y</tt> to value \c val
283  */
284  void set_val_ret(double &x, double &y, std::string name, double val);
285 
286  /** \brief Set element in slice of index \c z at the nearest
287  location to <tt>x,y</tt> to value \c val
288  */
289  void set_val_ret(double &x, double &y, size_t z, double val);
290 
291  /** \brief Get element in slice \c name at location closest to
292  <tt>x,y</tt>, and also return the corresponding values of \c x
293  and \c y
294  */
295  double &get_val_ret(double &x, double &y, std::string name);
296 
297  /** \brief Get element in slice \c name at location closest to
298  <tt>x,y</tt>, and also return the corresponding values of \c x
299  and \c y
300  */
301  const double &get_val_ret(double &x, double &y, std::string name) const;
302 
303  /** \brief Get element in slice of index \c z at location closest
304  to <tt>x,y</tt>, and also return the corresponding values of
305  \c x and \c y
306  */
307  double &get_val_ret(double &x, double &y, size_t z);
308 
309  /** \brief Get element in slice of index \c z at location closest
310  to <tt>x,y</tt>, and also return the corresponding values of
311  \c x and \c y
312  */
313  const double &get_val_ret(double &x, double &y, size_t z) const;
314 
315  /** \brief This function adds a slice from a different table3d
316  object, interpolating the results into the current
317  table3d object
318  */
319  void add_slice_from_table(table3d &source, std::string slice,
320  std::string dest_slice="");
321 
322  /** \brief Set elements in the first <tt>nv</tt> slices at the
323  nearest location to <tt>x,y</tt> to values \c vals
324  */
325  template<class vec_t>
326  void set_slices_ret(double &x, double &y, size_t nv, vec_t &vals) {
327  size_t ix, iy;
328  lookup_x(x,ix);
329  lookup_y(y,iy);
330  x=xval[ix];
331  y=yval[iy];
332 
333  for(size_t i=0;i<nv && i<list.size();i++) {
334  list[i](ix,iy)=vals[i];
335  }
336  return;
337  }
338 
339  /** \brief Get elements in the first <tt>nv</tt> slices at the
340  nearest location to <tt>x,y</tt> to value \c val
341  */
342  template<class vec_t>
343  void get_slices_ret(double &x, double &y, size_t nv, vec_t &vals) {
344 
345  size_t ix, iy;
346  lookup_x(x,ix);
347  lookup_y(y,iy);
348  x=xval[ix];
349  y=yval[iy];
350 
351  for(size_t i=0;i<nv && i<list.size();i++) {
352  vals[i]=list[i](ix,iy);
353  }
354  return;
355  }
356 
357  //@}
358 
359  // --------------------------------------------------------
360  /// \name Grid information get and set methods
361  //@{
362 
363  /// Set x grid point at index \c ix
364  void set_grid_x(size_t ix, double val);
365 
366  /// Set y grid point at index \c iy
367  void set_grid_y(size_t iy, double val);
368 
369  /// Get x grid point at index \c ix
370  double get_grid_x(size_t ix) const;
371 
372  /// Get y grid point at index \c iy
373  double get_grid_y(size_t iy) const;
374 
375  /// Get the name of the x grid variable
376  std::string get_x_name() const;
377 
378  /// Get the name of the y grid variable
379  std::string get_y_name() const;
380 
381  /// Set the name of the x grid variable
382  void set_x_name(std::string name);
383 
384  /// Set the name of the y grid variable
385  void set_y_name(std::string name);
386 
387  /// Get a const reference to the full x grid
388  const ubvector &get_x_data() const;
389 
390  /// Get a const reference to the full y grid
391  const ubvector &get_y_data() const;
392  //@}
393 
394  // --------------------------------------------------------
395  /// \name Size get methods
396  //@{
397  /// Get the size of the slices
398  void get_size(size_t &nx, size_t &ny) const;
399 
400  /// Get the x size
401  size_t get_nx() const;
402 
403  /// Get the y size
404  size_t get_ny() const;
405 
406  /// Get the number of slices
407  size_t get_nslices() const;
408 
409  /// True if the size of the table has been set
410  bool is_size_set() const;
411 
412  /// True if the grid has been set
413  bool is_xy_set() const;
414  //@}
415 
416  // --------------------------------------------------------
417  /// \name Slice manipulation
418  //@{
419 
420  /// Create a set of new slices specified in the string \c names
421  void line_of_names(std::string names);
422 
423  /** \brief Returns the name of slice with index \c z
424  */
425  std::string get_slice_name(size_t z) const;
426 
427  /** \brief Add a new slice
428  */
429  void new_slice(std::string name);
430 
431  /** \brief Set all of the values in slice \c name to \c val
432  */
433  void set_slice_all(std::string name, double val);
434 
435  /** \brief Find the index for slice named \c name
436  */
437  size_t lookup_slice(std::string name) const;
438 
439  /** \brief Return true if slice is already present
440  */
441  bool is_slice(std::string name, size_t &ix) const;
442 
443  /** \brief Rename slice named \c olds to \c news
444 
445  This is slow since we have to delete the column and re-insert
446  it. This process in turn mangles all of the iterators in the
447  list.
448  */
449  void rename_slice(std::string olds, std::string news);
450 
451  /** \brief Make a new slice named \c dest which is a copy
452  of the slice with name given in \c src.
453  */
454  void copy_slice(std::string src, std::string dest);
455 
456  /** \brief Initialize all values of slice named \c scol to \c val
457 
458  \note This will call the error handler if the value \c val is
459  not finite (i.e. either <tt>Inf</tt> or <tt>NaN</tt>).
460  */
461  void init_slice(std::string scol, double val);
462 
463  /// Return a constant reference to a slice
464  const ubmatrix &get_slice(std::string scol) const;
465 
466  /// Return a constant reference to a slice
467  const ubmatrix &get_slice(size_t iz) const;
468 
469  /// Return a constant reference to a slice
470  ubmatrix &get_slice(std::string scol);
471 
472  /// Return a constant reference to a slice
473  ubmatrix &get_slice(size_t iz);
474 
475  /** \brief Return a constant reference to all the slice data
476 
477  \comment
478  This isn't designated const, i.e. as
479  const std::vector<ubmatrix> &get_data() const;
480  because it would then have to be
481  const std::vector<const ubmatrix> &get_data() const;
482  \endcomment
483  */
484  const std::vector<ubmatrix> &get_data();
485 
486  /** \brief Copy to a slice from a generic matrix object
487 
488  The type <tt>mat_t</tt> can be any type with an
489  <tt>operator(,)</tt> method.
490  */
491  template<class mat_t>
492  void copy_to_slice(mat_t &m, std::string slice_name) {
493  for(size_t i=0;i<numx;i++) {
494  for(size_t j=0;j<numy;j++) {
495  this->set(i,j,slice_name,m(i,j));
496  }
497  }
498  return;
499  }
500  //@}
501 
502  // --------------------------------------------------------
503  /// \name Lookup and search methods
504  //@{
505  /** \brief Look for a value in the x grid
506  */
507  void lookup_x(double val, size_t &ix) const;
508 
509  /** \brief Look for a value in the y grid
510  */
511  void lookup_y(double val, size_t &iy) const;
512 
513  /** \brief Look for a value in a specified slice
514  */
515  void lookup(double val, std::string slice, size_t &ix,
516  size_t &iy) const;
517  //@}
518 
519  // --------------------------------------------------------
520  /// \name Interpolation, differentiation, and integration
521  //@{
522 
523  /** \brief Specify the interpolation type
524  */
525  void set_interp_type(size_t interp_type);
526 
527  /** \brief Get the interpolation type
528  */
529  size_t get_interp_type() const;
530 
531  /** \brief Interpolate \c x and \c y in slice named \c name
532  */
533  double interp(double x, double y, std::string name) const;
534 
535  /** \brief Interpolate the derivative of the data with respect to
536  the x grid at point \c x and \c y in slice named \c name
537  */
538  double deriv_x(double x, double y, std::string name) const;
539 
540  /** \brief Interpolate the derivative of the data with respect to
541  the y grid at point \c x and \c y in slice named \c name
542  */
543  double deriv_y(double x, double y, std::string name) const;
544 
545  /** \brief Interpolate the mixed second derivative of the data at
546  point \c x and \c y in slice named \c name
547  */
548  double deriv_xy(double x, double y, std::string name) const;
549 
550  /** \brief Interpolate the integral of the data
551  respect to the x grid
552  */
553  double integ_x(double x1, double x2, double y, std::string name) const;
554 
555  /** \brief Interpolate the integral of the data
556  respect to the y grid
557  */
558  double integ_y(double x, double y1, double y2, std::string name) const;
559 
560  /** \brief Fill a vector of interpolated values from each slice at the
561  point <tt>x,y</tt>
562  */
563  template<class vec_t>
564  void interp_slices(double x, double y, size_t nv, vec_t &v) {
565 
566  for (size_t i=0;i<list.size();i++) {
567  std::string name=get_slice_name(i);
568  v[i]=interp(x,y,name);
569  }
570 
571  return;
572  }
573 
574  /** \brief Create a new slice, named \c fpname, containing the
575  derivative of \c fname with respect to the x coordinate
576  */
577  void deriv_x(std::string fname, std::string fpname);
578 
579  /** \brief Create a new slice, named \c fpname, containing the
580  derivative of \c fname with respect to the y coordinate
581  */
582  void deriv_y(std::string fname, std::string fpname);
583  //@}
584 
585  // --------------------------------------------------------
586  /// \name Extract 2-dimensional tables
587  //@{
588  /** \brief Extract a table at a fixed x grid point
589 
590  \note All of the information previously stored in \c t will
591  be lost.
592  */
593  void extract_x(double x, table<> &t);
594 
595  /** \brief Extract a table at a fixed y grid point
596 
597  \note All of the information previously stored in \c t will
598  be lost.
599  */
600  void extract_y(double y, table<> &t);
601  //@}
602 
603  // --------------------------------------------------------
604  /// \name Clear methods
605  //@{
606  /** \brief Zero the data entries but keep the slice names
607  and grid
608  */
609  void zero_table();
610 
611  /** \brief Clear everything
612  */
613  void clear();
614 
615  /** \brief Remove all of the data by setting the number
616  of lines to zero
617 
618  This leaves the column names intact and does not remove
619  the constants.
620  */
621  void clear_data();
622  //@}
623 
624  // --------------------------------------------------------
625  /// \name Summary method
626  //@{
627  /** \brief Output a summary of the information stored
628 
629  Outputs the number of constants, the grid information,
630  and a list of the slice names
631  */
632  void summary(std::ostream *out, int ncol=79) const;
633  //@}
634 
635  // ---------
636  // Allow HDF I/O functions to access table3d data
637 
638  friend void o2scl_hdf::hdf_output(o2scl_hdf::hdf_file &hf,
639  table3d &t,
640  std::string name);
641 
642  friend void o2scl_hdf::hdf_input(o2scl_hdf::hdf_file &hf, table3d &t,
643  std::string name);
644 
645  // --------------------------------------------------------
646  /// \name Contour lines method
647  //@{
648 
649  /** \brief Create contour lines from the slice named \c name
650 
651  This uses \ref contour to compute contour lines (stored in \c
652  clines) from slice \c name given \c nlev contour levels in \c
653  levs .
654  */
655  template<class vec_t>
656  void slice_contours(std::string name, size_t nlev, vec_t &levs,
657  std::vector<contour_line> &clines) {
658 
659  size_t z=lookup_slice(name);
660 
661  contour co;
662  co.set_data(numx,numy,xval,yval,list[z]);
663  co.set_levels(nlev,levs);
664  co.calc_contours(clines);
665 
666  return;
667  }
668  //@}
669 
670  // --------------------------------------------------------
671  /// \name Manipulating constants
672  //@{
673  /** \brief Add a constant, or if the constant already exists, change
674  its value
675  */
676  virtual void add_constant(std::string name, double val);
677 
678  /// Remove a constant
679  virtual void remove_constant(std::string name);
680 
681  /** \brief Set a constant equal to a value, but don't add it if
682  not already present
683 
684  If \c err_on_notfound is <tt>true</tt> (the default), then
685  this function throws an exception if a constant with
686  name \c name is not found. If \c err_on_notfound is
687  <tt>false</tt>, then if a constant with name \c name
688  is not found this function just silently returns
689  \ref o2scl::exc_enotfound.
690  */
691  virtual int set_constant(std::string name, double val,
692  bool err_on_notfound=true);
693 
694  /// Test if \c name is a constant
695  virtual bool is_constant(std::string name) const;
696 
697  /// Get a constant
698  virtual double get_constant(std::string name);
699 
700  /// Get a constant by index
701  virtual void get_constant(size_t ix, std::string &name,
702  double &val) const;
703 
704  /// Get the number of constants
705  virtual size_t get_nconsts() const {
706  return constants.size();
707  }
708  //@}
709 
710  /// \name Miscellaneous methods
711  //@{
712  /** \brief Read a generic table3d object specified as a
713  text file
714 
715  This function reads a set of columns of numerical values,
716  presuming that the first column is the x-grid value, the
717  second column is the y-grid value, and the remaining columns
718  are slices to be added. If the first row appears to be strings
719  rather than numerical values, then the first row is used for
720  the x name, y name, and slice names. Values in the first two
721  columns which differ by less than \c eps are assumed to refer
722  to the same grid point. If not all combinations of x and y are
723  found, then those entries are left unchanged in all slices.
724 
725  \future It would be great to add a function which generates
726  a text file in this format as well.
727  */
728  int read_gen3_list(std::istream &fin, int verbose=0, double eps=1.0e-12);
729 
730  /** \brief Set the current table3d object by reading a
731  \ref o2scl::table
732  */
733  template<class vec_t>
735  std::string xname2="", std::string yname2="",
736  double empty_value=0.0, int verbose=0,
737  bool err_on_fail=true, double eps=1.0e-12) {
738 
739  clear();
740 
741  if (tab.get_ncolumns()<3) {
742  if (err_on_fail) {
743  O2SCL_ERR2("Not enough columns of data in ",
744  "table3d::read_table().",o2scl::exc_einval);
745  } else {
746  return o2scl::exc_einval;
747  }
748  }
749 
750  // Set up xname and yname if unspecified
751  if (xname2.length()==0) {
752  if (yname2==tab.get_column_name(0)) {
753  xname2=tab.get_column_name(1);
754  } else {
755  xname2=tab.get_column_name(0);
756  yname2=tab.get_column_name(1);
757  }
758  }
759  if (yname2.length()==0) {
760  if (xname2==tab.get_column_name(1)) {
761  yname2=tab.get_column_name(0);
762  } else {
763  yname2=tab.get_column_name(1);
764  }
765  }
766 
767  // Setup x and y grid vectors from data
768  const vec_t &xdata=tab[xname2];
769  const vec_t &ydata=tab[yname2];
770 
771  std::vector<double> xgrid, ygrid;
772 
773  // Note it is important that this loop ends at tab.get_nlines()
774  // rather than xdata.size() because the internal vector storage
775  // can be larger than the actual table size
776  for(size_t i=0;i<tab.get_nlines();i++) {
777 
778  // Look for x value in x grid
779  bool found=false;
780  for(size_t j=0;j<xgrid.size();j++) {
781  if (fabs(xdata[i]-xgrid[j])/fabs(xgrid[j])<eps) {
782  found=true;
783  }
784  }
785 
786  // If not found, add it
787  if (found==false) {
788  xgrid.push_back(xdata[i]);
789  }
790 
791  // Now look for y value in y grid
792  found=false;
793  for(size_t j=0;j<ygrid.size();j++) {
794  if (fabs(ydata[i]-ygrid[j])/fabs(ygrid[j])<eps) {
795  found=true;
796  }
797  }
798 
799  // If not found, add it
800  if (found==false) {
801  ygrid.push_back(ydata[i]);
802  }
803  }
804 
805  if (verbose>1) {
806  std::cout << "X grid: " << std::endl;
807  for(size_t k=0;k<xgrid.size();k++) {
808  std::cout << k << " " << xgrid[k] << std::endl;
809  }
810  std::cout << "Y grid: " << std::endl;
811  for(size_t k=0;k<ygrid.size();k++) {
812  std::cout << k << " " << ygrid[k] << std::endl;
813  }
814  }
815 
816  // Sor the grids
817  vector_sort_double(xgrid.size(),xgrid);
818  vector_sort_double(ygrid.size(),ygrid);
819 
820  // Set grid from x and y grid vectors
821  set_xy(xname2,xgrid.size(),xgrid,yname2,ygrid.size(),ygrid);
822 
823  // Create new slices
824  std::vector<std::string> sl_names;
825  for(size_t i=0;i<tab.get_ncolumns();i++) {
826  if (tab.get_column_name(i)!=xname2 && tab.get_column_name(i)!=yname2) {
827  std::string sl=tab.get_column_name(i);
828  if (verbose>0) {
829  std::cout << "New slice: " << sl << std::endl;
830  }
831  sl_names.push_back(sl);
832  new_slice(sl);
833  set_slice_all(sl,empty_value);
834  }
835  }
836 
837  // Set the data
838  for(size_t i=0;i<tab.get_ncolumns();i++) {
839  if (tab.get_column_name(i)!=xname2 && tab.get_column_name(i)!=yname2) {
840  std::string sl=tab.get_column_name(i);
841  for(size_t j=0;j<tab.get_nlines();j++) {
842  set_val(xdata[j],ydata[j],sl,tab.get(i,j));
843  }
844  }
845  }
846 
847  return 0;
848  }
849 
850  /// Return the type, \c "table3d".
851  virtual const char *type() { return "table3d"; }
852  //@}
853 
854  /** \name Parsing mathematical functions specified as strings
855 
856  \comment
857  Note that doxygen doesn't format the documentation propertly
858  if the \name specification covers more than one line
859  \endcomment
860  */
861  //@{
862  /** \brief Fill a matrix from the function specified in \c function
863 
864  \comment
865  This function must return an int rather than void because
866  of the presence of the 'throw_on_err' mechanism
867  \endcomment
868  */
869  template<class resize_mat_t>
870  int function_matrix(std::string function, resize_mat_t &mat,
871  bool throw_on_err=true) {
872 
873  calculator calc;
874  std::map<std::string,double> vars;
875 
876  std::map<std::string,double>::const_iterator mit;
877  for(mit=constants.begin();mit!=constants.end();mit++) {
878  vars[mit->first]=mit->second;
879  }
880 
881  calc.compile(function.c_str(),&vars);
882 
883  if (mat.size1()!=numx || mat.size2()!=numy) {
884  mat.resize(numx,numy);
885  }
886 
887  for(size_t i=0;i<numx;i++) {
888  for(size_t j=0;j<numy;j++) {
889  vars[xname]=xval[i];
890  vars[yname]=yval[j];
891 
892  for(size_t k=0;k<list.size();k++) {
893  vars[get_slice_name(k)]=list[k](i,j);
894  }
895 
896  mat(i,j)=calc.eval(&vars);
897  }
898  }
899 
900 
901  return 0;
902  }
903 
904  /** \brief Make a column from <tt>function</tt> and add it to the table.
905 
906  If the column already exists, the data already present is
907  overwritten with the result.
908  */
909  void function_slice(std::string function, std::string col);
910  //@}
911 
912  /** \brief Copy slice named \c slice to a new \ref o2scl::table3d
913  object with a uniform grid using the current interpolation type
914  */
915  table3d slice_to_uniform_grid(std::string slice, size_t xpts,
916  bool log_x, size_t ypts, bool log_y);
917 
918  /** \brief Copy entire table to a new \ref o2scl::table3d
919  object with a uniform grid using the current interpolation type
920  */
921  table3d table_to_uniform_grid(size_t xpts, bool log_x,
922  size_t ypts, bool log_y);
923 
924 #ifdef O2SCL_NEVER_DEFINED
925 
926  /** \brief Desc
927  */
928  template<class vec_t>
929  void test_uniform_grid_log(std::string slice, bool &log_x, bool &log_y,
930  vec_t &x, vec_t &y, vec_t &s) {
931  vector<double> dev_x;
932  for(size_t i=0;i<numx-1;i++) {
933  dev_x.push_back(xval[i+1]-xval[i]);
934  }
935  double avg=vector_mean(dev_x.size(),dev_x);
936  double stddev=vector_stddev(dev_x.size(),dev_x);
937  bool x_set=false;
938  log_x=false;
939  if (stddev<1.0e-8*fabs(avg)) {
940  x.resize(numx);
941  for(size_t i=0;i<numx;i++) {
942  x[i]=xval[i];
943  }
944  x_set=true;
945  } else {
946  dev_x.clear();
947  for(size_t i=0;i<numx-1;i++) {
948  dev_x.push_back(xval[i+1]/xval[i]);
949  }
950  avg=vector_mean(dev_x.size(),dev_x);
951  stddev=vector_stddev(dev_x.size(),dev_x);
952  if (stddev<1.0e-8*fabs(avg)) {
953  x.resize(numx);
954  for(size_t i=0;i<numx;i++) {
955  x[i]=xval[i];
956  }
957  x_set=true;
958  log_x=true;
959  }
960  }
961  if (x_set==false) {
962  uniform_grid_end
963  }
964  return;
965  }
966 #endif
967 
968  protected:
969 
970 #ifndef DOXYGEN_INTERNAL
971 
972  /// \name Iterator types
973  //@{
974  typedef std::map<std::string,size_t,
975  std::greater<std::string> >::iterator map_iter;
976  typedef std::map<std::string,size_t,
977  std::greater<std::string> >::const_iterator map_const_iter;
978  //@}
979 
980  /// \name Data storage
981  //@{
982  /// The list of constants
983  std::map<std::string,double> constants;
984 
985  /// The size of the x grid
986  size_t numx;
987 
988  /// The size of the y grid
989  size_t numy;
990 
991  /// A tree connecting column names to list indexes
992  std::map<std::string,size_t,std::greater<std::string> > tree;
993 
994  /// The name for the x grid
995  std::string xname;
996 
997  /// The name for the y grid
998  std::string yname;
999 
1000  /// The pointers to the matrices
1001  std::vector<ubmatrix> list;
1002 
1003  /// The x grid
1005 
1006  /// The y grid
1008 
1009  /// True if the grid has been set
1010  bool xy_set;
1011 
1012  /// True if the size of the grid has been set
1013  bool size_set;
1014 
1015  /// True if the table has at least one slice
1017 
1018  /// The interpolation type
1019  size_t itype;
1020  //@}
1021 
1022  /// \name Tree iterator boundaries
1023  //@{
1024  /// Return the beginning of the slice tree
1025  map_iter begin() {return tree.begin();};
1026  /// Return the end of the slice tree
1027  map_iter end() {return tree.end();};
1028  /// Return the beginning of the slice tree
1029  map_const_iter const_begin() const {return tree.begin();};
1030  /// Return the end of the slice tree
1031  map_const_iter const_end() const {return tree.end();};
1032  //@}
1033 
1034 #endif
1035 
1036  };
1037 
1038 #ifndef DOXYGEN_NO_O2NS
1039 }
1040 #endif
1041 
1042 #endif
boost::numeric::ublas::matrix< double >
o2scl::table3d::is_constant
virtual bool is_constant(std::string name) const
Test if name is a constant.
o2scl::table3d::read_table
int read_table(const o2scl::table< vec_t > &tab, std::string xname2="", std::string yname2="", double empty_value=0.0, int verbose=0, bool err_on_fail=true, double eps=1.0e-12)
Set the current table3d object by reading a o2scl::table.
Definition: table3d.h:734
o2scl::contour::set_levels
void set_levels(size_t nlevels, vec_t &ulevels)
Set the contour levels.
Definition: contour.h:337
o2scl::table
Data table table class.
Definition: table.h:49
o2scl::table::get_nlines
size_t get_nlines() const
Return the number of lines.
Definition: table.h:460
o2scl::table3d::end
map_iter end()
Return the end of the slice tree.
Definition: table3d.h:1027
o2scl::table3d::extract_y
void extract_y(double y, table<> &t)
Extract a table at a fixed y grid point.
o2scl::table3d::set
void set(size_t ix, size_t iy, std::string name, double val)
Set element in slice name at location ix,iy to value val.
o2scl::table3d::lookup
void lookup(double val, std::string slice, size_t &ix, size_t &iy) const
Look for a value in a specified slice.
boost::numeric::ublas::vector< double >
o2scl::table3d::extract_x
void extract_x(double x, table<> &t)
Extract a table at a fixed x grid point.
o2scl::table3d::const_end
map_const_iter const_end() const
Return the end of the slice tree.
Definition: table3d.h:1031
o2scl::table3d::line_of_names
void line_of_names(std::string names)
Create a set of new slices specified in the string names.
O2SCL_ERR2
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
Definition: err_hnd.h:281
o2scl::table3d::set_val_ret
void set_val_ret(double &x, double &y, std::string name, double val)
Set element in slice name at the nearest location to x,y to value val.
o2scl::table3d::new_slice
void new_slice(std::string name)
Add a new slice.
o2scl::table3d::get_y_name
std::string get_y_name() const
Get the name of the y grid variable.
o2scl::table3d::get_val
double & get_val(double x, double y, std::string name)
Get element in slice name at location closest to x,y
o2scl::table3d::get_x_name
std::string get_x_name() const
Get the name of the x grid variable.
o2scl
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
Definition: anneal.h:42
o2scl::vector_mean
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
o2scl::table3d::get_constant
virtual double get_constant(std::string name)
Get a constant.
o2scl::table3d::get_nconsts
virtual size_t get_nconsts() const
Get the number of constants.
Definition: table3d.h:705
o2scl::table3d::integ_y
double integ_y(double x, double y1, double y2, std::string name) const
Interpolate the integral of the data respect to the y grid.
o2scl::table3d::deriv_x
double deriv_x(double x, double y, std::string name) const
Interpolate the derivative of the data with respect to the x grid at point x and y in slice named nam...
o2scl::table3d::copy_slice
void copy_slice(std::string src, std::string dest)
Make a new slice named dest which is a copy of the slice with name given in src.
o2scl::contour
Calculate contour lines from a two-dimensional data set.
Definition: contour.h:240
o2scl_hdf::hdf_input
void hdf_input(hdf_file &hf, o2scl::table< vec_t > &t, std::string name)
Input a o2scl::table object from a hdf_file.
Definition: hdf_io.h:148
o2scl::table3d::get_y_data
const ubvector & get_y_data() const
Get a const reference to the full y grid.
o2scl_inte_qng_coeffs::x2
static const double x2[5]
Definition: inte_qng_gsl.h:66
o2scl::table3d::copy_to_slice
void copy_to_slice(mat_t &m, std::string slice_name)
Copy to a slice from a generic matrix object.
Definition: table3d.h:492
o2scl::table3d::integ_x
double integ_x(double x1, double x2, double y, std::string name) const
Interpolate the integral of the data respect to the x grid.
o2scl::table3d::size_set
bool size_set
True if the size of the grid has been set.
Definition: table3d.h:1013
o2scl::calculator::compile
void compile(const char *expr, std::map< std::string, double > *vars=0, bool debug=false, std::map< std::string, int > opPrec=opPrecedence)
Compile expression expr using variables specified in vars and return an integer to indicate success o...
o2scl::table_units
Data table table class with units.
Definition: table_units.h:42
o2scl::table3d::get_slices
void get_slices(double x, double y, size_t nv, vec_t &v)
Get the data for every slice at the nearest location to x,y
Definition: table3d.h:262
o2scl::table3d::table3d
table3d()
Create a new 3D table table .
o2scl::table3d::yname
std::string yname
The name for the y grid.
Definition: table3d.h:998
o2scl::table3d::add_slice_from_table
void add_slice_from_table(table3d &source, std::string slice, std::string dest_slice="")
This function adds a slice from a different table3d object, interpolating the results into the curren...
o2scl::vector_stddev
double vector_stddev(size_t n, const vec_t &data)
Standard deviation with specified mean.
Definition: vec_stats.h:253
o2scl::table3d::numy
size_t numy
The size of the y grid.
Definition: table3d.h:989
o2scl::table3d::xy_set
bool xy_set
True if the grid has been set.
Definition: table3d.h:1010
o2scl_inte_qng_coeffs::x1
static const double x1[5]
Definition: inte_qng_gsl.h:48
o2scl::table3d::zero_table
void zero_table()
Zero the data entries but keep the slice names and grid.
o2scl::table3d::interp_slices
void interp_slices(double x, double y, size_t nv, vec_t &v)
Fill a vector of interpolated values from each slice at the point x,y
Definition: table3d.h:564
o2scl::table3d::init_slice
void init_slice(std::string scol, double val)
Initialize all values of slice named scol to val.
o2scl::table3d::set_size
void set_size(size_t nx, size_t ny)
Initialize table table size.
o2scl::exc_einval
@ exc_einval
invalid argument supplied by user
Definition: err_hnd.h:59
o2scl::table3d::xval
ubvector xval
The x grid.
Definition: table3d.h:1004
o2scl::calculator
Evaluate a mathematical expression in a string.
Definition: shunting_yard.h:119
o2scl::table3d::set_grid_x
void set_grid_x(size_t ix, double val)
Set x grid point at index ix.
o2scl::table3d::constants
std::map< std::string, double > constants
The list of constants.
Definition: table3d.h:983
o2scl::table3d::set_grid_y
void set_grid_y(size_t iy, double val)
Set y grid point at index iy.
o2scl::table3d::deriv_y
double deriv_y(double x, double y, std::string name) const
Interpolate the derivative of the data with respect to the y grid at point x and y in slice named nam...
o2scl::table3d::tree
std::map< std::string, size_t, std::greater< std::string > > tree
A tree connecting column names to list indexes.
Definition: table3d.h:992
o2scl::table3d::get_x_data
const ubvector & get_x_data() const
Get a const reference to the full x grid.
o2scl::table3d::get_nx
size_t get_nx() const
Get the x size.
o2scl_hdf
The O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl namespace ...
Definition: table.h:59
o2scl::contour::calc_contours
void calc_contours(std::vector< contour_line > &clines)
Calculate the contours.
o2scl::table3d::clear_data
void clear_data()
Remove all of the data by setting the number of lines to zero.
o2scl::table3d::get
double & get(size_t ix, size_t iy, std::string name)
Get element in slice name at location ix,iy
o2scl::table3d::get_grid_x
double get_grid_x(size_t ix) const
Get x grid point at index ix.
o2scl::table3d::set_interp_type
void set_interp_type(size_t interp_type)
Specify the interpolation type.
o2scl::table3d::clear
void clear()
Clear everything.
o2scl::table3d::get_slice_name
std::string get_slice_name(size_t z) const
Returns the name of slice with index z.
o2scl::uniform_grid< double >
o2scl::table3d::lookup_x
void lookup_x(double val, size_t &ix) const
Look for a value in the x grid.
o2scl::table3d::begin
map_iter begin()
Return the beginning of the slice tree.
Definition: table3d.h:1025
o2scl::table3d::set_slices
void set_slices(double x, double y, size_t nv, vec_t &vals)
Set elements in the first nv slices at the nearest location to x,y to value val.
Definition: table3d.h:247
o2scl::table3d::set_xy
void set_xy(std::string x_name, size_t nx, const vec_t &x, std::string y_name, size_t ny, const vec2_t &y)
Initialize the x-y grid.
Definition: table3d.h:120
o2scl::table3d::interp
double interp(double x, double y, std::string name) const
Interpolate x and y in slice named name.
o2scl::table3d::deriv_xy
double deriv_xy(double x, double y, std::string name) const
Interpolate the mixed second derivative of the data at point x and y in slice named name.
o2scl::table::get_ncolumns
size_t get_ncolumns() const
Return the number of columns.
Definition: table.h:452
o2scl::table3d::get_nslices
size_t get_nslices() const
Get the number of slices.
o2scl::table3d::slice_contours
void slice_contours(std::string name, size_t nlev, vec_t &levs, std::vector< contour_line > &clines)
Create contour lines from the slice named name.
Definition: table3d.h:656
o2scl::table3d::slice_to_uniform_grid
table3d slice_to_uniform_grid(std::string slice, size_t xpts, bool log_x, size_t ypts, bool log_y)
Copy slice named slice to a new o2scl::table3d object with a uniform grid using the current interpola...
o2scl::table3d::const_begin
map_const_iter const_begin() const
Return the beginning of the slice tree.
Definition: table3d.h:1029
o2scl::table3d::set_y_name
void set_y_name(std::string name)
Set the name of the y grid variable.
o2scl::table3d::is_size_set
bool is_size_set() const
True if the size of the table has been set.
o2scl::table3d::function_matrix
int function_matrix(std::string function, resize_mat_t &mat, bool throw_on_err=true)
Fill a matrix from the function specified in function.
Definition: table3d.h:870
o2scl::table3d::set_constant
virtual int set_constant(std::string name, double val, bool err_on_notfound=true)
Set a constant equal to a value, but don't add it if not already present.
o2scl::table3d::is_slice
bool is_slice(std::string name, size_t &ix) const
Return true if slice is already present.
o2scl::table3d::lookup_slice
size_t lookup_slice(std::string name) const
Find the index for slice named name.
o2scl::table3d::set_val
void set_val(double x, double y, std::string name, double val)
Set element in slice name at the nearest location to x,y to value val.
o2scl::table3d::is_xy_set
bool is_xy_set() const
True if the grid has been set.
o2scl::vector_sort_double
void vector_sort_double(size_t n, vec_t &data)
Sort a vector of doubles (in increasing order)
Definition: vector.h:895
o2scl_hdf::hdf_file
Store data in an O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$sc...
Definition: hdf_file.h:105
o2scl::table3d::get_ny
size_t get_ny() const
Get the y size.
o2scl::table3d::set_slice_all
void set_slice_all(std::string name, double val)
Set all of the values in slice name to val.
O2SCL_ERR
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Definition: err_hnd.h:273
o2scl::table3d::get_size
void get_size(size_t &nx, size_t &ny) const
Get the size of the slices.
o2scl::table3d::table_to_uniform_grid
table3d table_to_uniform_grid(size_t xpts, bool log_x, size_t ypts, bool log_y)
Copy entire table to a new o2scl::table3d object with a uniform grid using the current interpolation ...
o2scl::table3d::summary
void summary(std::ostream *out, int ncol=79) const
Output a summary of the information stored.
o2scl::table3d::get_grid_y
double get_grid_y(size_t iy) const
Get y grid point at index iy.
o2scl::table3d::yval
ubvector yval
The y grid.
Definition: table3d.h:1007
o2scl::calculator::eval
double eval(std::map< std::string, double > *vars=0)
Evalate the previously compiled expression using variables specified in vars.
o2scl::table3d::numx
size_t numx
The size of the x grid.
Definition: table3d.h:986
o2scl::table3d::function_slice
void function_slice(std::string function, std::string col)
Make a column from function and add it to the table.
o2scl::table3d::add_constant
virtual void add_constant(std::string name, double val)
Add a constant, or if the constant already exists, change its value.
o2scl::table3d::get_val_ret
double & get_val_ret(double &x, double &y, std::string name)
Get element in slice name at location closest to x,y, and also return the corresponding values of x a...
o2scl::table3d::get_slice
const ubmatrix & get_slice(std::string scol) const
Return a constant reference to a slice.
o2scl::table3d::xname
std::string xname
The name for the x grid.
Definition: table3d.h:995
o2scl::table::get
double get(std::string scol, size_t row) const
Get value from row row of column named col. .
Definition: table.h:408
o2scl::table3d::has_slice
bool has_slice
True if the table has at least one slice.
Definition: table3d.h:1016
o2scl::table3d::remove_constant
virtual void remove_constant(std::string name)
Remove a constant.
o2scl::table3d::read_gen3_list
int read_gen3_list(std::istream &fin, int verbose=0, double eps=1.0e-12)
Read a generic table3d object specified as a text file.
o2scl::contour::set_data
void set_data(size_t sizex, size_t sizey, const vec_t &x_fun, const vec_t &y_fun, const mat_t &udata)
Set the data.
Definition: contour.h:266
o2scl::table3d::lookup_y
void lookup_y(double val, size_t &iy) const
Look for a value in the y grid.
o2scl::table3d::set_x_name
void set_x_name(std::string name)
Set the name of the x grid variable.
o2scl::table3d::itype
size_t itype
The interpolation type.
Definition: table3d.h:1019
o2scl::table3d::list
std::vector< ubmatrix > list
The pointers to the matrices.
Definition: table3d.h:1001
o2scl::table3d::operator=
table3d & operator=(const table3d &t)
Copy constructor.
o2scl::table3d
A data structure containing one or more slices of two-dimensional data points defined on a grid.
Definition: table3d.h:78
o2scl::table::get_column_name
std::string get_column_name(size_t icol) const
Returns the name of column col .
Definition: table.h:819
o2scl::table3d::set_slices_ret
void set_slices_ret(double &x, double &y, size_t nv, vec_t &vals)
Set elements in the first nv slices at the nearest location to x,y to values vals.
Definition: table3d.h:326
o2scl::table3d::get_interp_type
size_t get_interp_type() const
Get the interpolation type.
o2scl::table3d::rename_slice
void rename_slice(std::string olds, std::string news)
Rename slice named olds to news.
o2scl::table3d::type
virtual const char * type()
Return the type, "table3d".
Definition: table3d.h:851
o2scl::table3d::get_data
const std::vector< ubmatrix > & get_data()
Return a constant reference to all the slice data.
o2scl::table3d::get_slices_ret
void get_slices_ret(double &x, double &y, size_t nv, vec_t &vals)
Get elements in the first nv slices at the nearest location to x,y to value val.
Definition: table3d.h:343

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