JAPAn
Just Another Parity Analyzer
Loading...
Searching...
No Matches
QwInterpolator< value_t, value_n > Class Template Reference

A multi-dimensional grid of values with interpolation methods. More...

#include <QwInterpolator.h>

+ Collaboration diagram for QwInterpolator< value_t, value_n >:

Public Member Functions

 QwInterpolator (const unsigned int ndim=1)
 Constructor with number of dimensions.
 
 QwInterpolator (const std::vector< coord_t > &min, const std::vector< coord_t > &max, const std::vector< coord_t > &step)
 Constructor with minimum, maximum, and step size.
 
 QwInterpolator (const std::string &filename)
 Constructor with file name.
 
virtual ~QwInterpolator ()
 Destructor.
 
void SetDimensions (const unsigned int ndim)
 Set the number of coordinate dimensions and resize vectors.
 
void SetMinimumMaximumStep (const coord_t min, const coord_t max, const coord_t step)
 Set minimum, maximum, and step size to single values.
 
void SetMinimumMaximumStep (const std::vector< coord_t > &min, const std::vector< coord_t > &max, const std::vector< coord_t > &step)
 Set minimum, maximum, and step size to different values.
 
coord_t GetMinimum (const unsigned int dim) const
 Get minimum in dimension.
 
coord_t GetMaximum (const unsigned int dim) const
 Get maximum in dimension.
 
coord_t GetStepSize (const unsigned int dim) const
 Get minimum in dimension.
 
unsigned int GetMaximumEntries () const
 Get the maximum number of entries.
 
unsigned int GetCurrentEntries () const
 Get the current number of entries.
 
unsigned int GetWrapCoordinate (const unsigned int dim) const
 Get wrapping coordinate.
 
void SetWrapCoordinate (const unsigned int dim, const size_t wrap=1)
 Set wrapping coordinate.
 
void SetWrapCoordinate (const std::vector< size_t > &wrap)
 
int GetDataReductionFactor (const unsigned int dim) const
 Get data reduction factor.
 
void SetDataReductionFactor (const unsigned int dim, const unsigned int redux)
 Set data reduction factor.
 
void SetDataReductionFactor (const unsigned int redux)
 
void SetDataReductionFactor (const std::vector< unsigned int > &redux)
 
void SetInterpolationMethod (const EQwInterpolationMethod method)
 Set the interpolation method.
 
EQwInterpolationMethod GetInterpolationMethod () const
 Get the interpolation method.
 
void PrintCoverage (const unsigned int dim)
 Print coverage map for all bins in one dimension.
 
bool InBounds (const coord_t *coord) const
 Return true if the coordinate is in bounds.
 
Functions to write grid values
bool Set (const coord_t &coord, const value_t &value)
 Set a single value at a coordinate (false if not possible)
 
bool Set (const coord_t *coord, const value_t &value)
 Set a single value at a coordinate (false if not possible)
 
bool Set (const coord_t *coord, const value_t *value)
 Set a set of values at a coordinate (false if not possible)
 
bool Set (const unsigned int linear_index, const value_t &value)
 Set a single value at a linearized index (false if not possible)
 
bool Set (const unsigned int linear_index, const value_t *value)
 Set a set of values at a linearized index (false if not possible)
 
bool Set (const unsigned int *cell_index, const value_t &value)
 Set a single value at a grid point (false if out of bounds)
 
bool Set (const unsigned int *cell_index, const value_t *value)
 Set a set of values at a grid point (false if out of bounds)
 
Functions to retrieve interpolated values
value_t GetValue (const coord_t &coord) const
 Get the interpolated value at coordinate (zero if out of bounds)
 
value_t GetValue (const coord_t *coord) const
 Get the interpolated value at coordinate (zero if out of bounds)
 
bool GetValue (const coord_t *coord, double &value) const
 Get the interpolated value at coordinate (zero if out of bounds)
 
bool GetValue (const coord_t *coord, double *value) const
 Get the interpolated value at coordinate (zero if out of bounds)
 
File reading and writing
bool WriteText (std::ostream &stream) const
 Write the grid as text.
 
bool WriteTextFile (const std::string &filename) const
 Write the grid to text file.
 
bool WriteTextScreen () const
 Write the grid to screen.
 
bool ReadText (std::istream &stream)
 Read the grid from text.
 
bool ReadTextFile (const std::string &filename)
 Read the grid from text file.
 
bool WriteBinaryFile (const std::string &filename) const
 Write the grid values to binary file.
 
bool ReadBinaryFile (const std::string &filename)
 Read the grid values from binary file.
 

Private Attributes

unsigned int fNDim
 Number of dimensions in coordinates.
 
std::vector< coord_tfMin
 Minimum in each dimension.
 
std::vector< coord_tfMax
 Maximum in each dimension.
 
std::vector< coord_tfStep
 Step size in each dimension.
 
std::vector< size_t > fSize
 Number of points in each dimension.
 
std::vector< size_t > fWrap
 Wrap around this coordinate.
 
std::vector< size_t > fRedux
 Data reduction factor.
 
std::vector< size_t > fExtent
 
std::vector< value_t > fValues [value_n]
 Table with pointers to arrays of values.
 
unsigned int fCurrentEntries
 Number of values read in.
 
unsigned int fMaximumEntries
 Maximum number of values.
 
unsigned int * f_cell_index
 Pre-allocated cell index.
 
double * f_cell_local
 Pre-allocated local coordinates.
 
EQwInterpolationMethod fInterpolationMethod
 Interpolation method.
 

Indexing functions (publicly available and unchecked)

unsigned int Index (const coord_t *coord) const
 Return the linearized index based on the point coordinates (unchecked)
 
unsigned int Index (const unsigned int *cell_index) const
 Return the linearized index based on the cell indices (unchecked)
 
unsigned int Index (const unsigned int *cell_index, const unsigned int offset) const
 Return the linearized index based on the cell indices and offset (unchecked)
 
void Cell (const coord_t coord, unsigned int &cell_index, double &cell_local, const unsigned int dim) const
 Return the cell index and local coordinates in one dimension (unchecked)
 
void Cell (const coord_t *coord, unsigned int *cell_index, double *cell_local) const
 Return the cell index and local coordinates (unchecked)
 
void Cell (const coord_t *coord, unsigned int *cell_index) const
 Return the cell index (unchecked)
 
void Cell (unsigned int linear_index, unsigned int *cell_index) const
 Return the cell index based on the linearized index (unchecked)
 
void Coord (const unsigned int *cell_index, coord_t *coord) const
 Return the coordinates based on the cell index (unchecked)
 
void Coord (const unsigned int linear_index, coord_t *coord) const
 Return the coordinates based on the linearized index (unchecked)
 
void Nearest (const coord_t *coord, unsigned int *cell_index) const
 Return the cell index closest to the coordinate (could be above) (unchecked)
 
bool Linear (const coord_t *coord, value_t *value) const
 Linear interpolation (unchecked)
 
bool NearestNeighbor (const coord_t *coord, value_t *value) const
 Nearest-neighbor (unchecked)
 
bool Check (const coord_t *coord) const
 Check for boundaries with coordinate argument.
 
bool Check (const unsigned int *cell_index) const
 Check for boundaries with cell index argument.
 
bool Check (const unsigned int linear_index) const
 Check for boundaries with linearized index argument.
 
value_t Get (const unsigned int *cell_index) const
 Get a single value by cell index (unchecked)
 
value_t Get (const unsigned int index) const
 Get a single value by linearized index (unchecked)
 
bool Get (const unsigned int index, value_t *value) const
 Get a vector value by linearized index (unchecked)
 

Detailed Description

template<class value_t = float, unsigned int value_n = 1>
class QwInterpolator< value_t, value_n >

A multi-dimensional grid of values with interpolation methods.

This class provides various interpolation methods on a multi-dimensional grid of multi-dimensional values. Linear interpolation and nearest-neighbor are implemented for all dimensions.

The template arguments are the internal storage data type (defaults to float) and the number of dimensions of the stored data at each grid point (defaults to one). The dimension of the grid itself is set through the constructor. To describe a double vector field with 5 components on a 3-dimensional grid, you would write

QwInterpolate<float,5> *field = new QwInterpolator<float,5>(3);
QwInterpolator(const unsigned int ndim=1)
Constructor with number of dimensions.

The minimum, maximum, and step size of the grid have to be known before the values are filled.

Definition at line 54 of file QwInterpolator.h.

Constructor & Destructor Documentation

◆ QwInterpolator() [1/3]

template<class value_t = float, unsigned int value_n = 1>
QwInterpolator< value_t, value_n >::QwInterpolator ( const unsigned int ndim = 1)
inline

Constructor with number of dimensions.

Definition at line 59 of file QwInterpolator.h.

59 {
62 f_cell_index = new unsigned int[fNDim];
63 f_cell_local = new double[fNDim];
64 };
A multi-dimensional grid of values with interpolation methods.
unsigned int * f_cell_index
Pre-allocated cell index.
void SetDimensions(const unsigned int ndim)
Set the number of coordinate dimensions and resize vectors.
unsigned int fNDim
Number of dimensions in coordinates.
void SetInterpolationMethod(const EQwInterpolationMethod method)
Set the interpolation method.
double * f_cell_local
Pre-allocated local coordinates.

References f_cell_index, f_cell_local, fNDim, kMultiLinear, SetDimensions(), and SetInterpolationMethod().

+ Here is the call graph for this function:

◆ QwInterpolator() [2/3]

template<class value_t = float, unsigned int value_n = 1>
QwInterpolator< value_t, value_n >::QwInterpolator ( const std::vector< coord_t > & min,
const std::vector< coord_t > & max,
const std::vector< coord_t > & step )
inline

Constructor with minimum, maximum, and step size.

Definition at line 66 of file QwInterpolator.h.

68 {
69 SetDimensions(min.size());
72 f_cell_index = new unsigned int[fNDim];
73 f_cell_local = new double[fNDim];
74 };
void SetMinimumMaximumStep(const coord_t min, const coord_t max, const coord_t step)
Set minimum, maximum, and step size to single values.

References f_cell_index, f_cell_local, fNDim, kMultiLinear, SetDimensions(), SetInterpolationMethod(), and SetMinimumMaximumStep().

+ Here is the call graph for this function:

◆ QwInterpolator() [3/3]

template<class value_t = float, unsigned int value_n = 1>
QwInterpolator< value_t, value_n >::QwInterpolator ( const std::string & filename)
inline

Constructor with file name.

Definition at line 76 of file QwInterpolator.h.

76 {
79 f_cell_index = new unsigned int[fNDim];
80 f_cell_local = new double[fNDim];
81 };
bool ReadBinaryFile(const std::string &filename)
Read the grid values from binary file.

References f_cell_index, f_cell_local, fNDim, kMultiLinear, ReadBinaryFile(), and SetInterpolationMethod().

+ Here is the call graph for this function:

◆ ~QwInterpolator()

template<class value_t = float, unsigned int value_n = 1>
virtual QwInterpolator< value_t, value_n >::~QwInterpolator ( )
inlinevirtual

Destructor.

Definition at line 83 of file QwInterpolator.h.

83 {
84 delete[] f_cell_index;
85 delete[] f_cell_local;
86 };

References f_cell_index, and f_cell_local.

Member Function Documentation

◆ Cell() [1/4]

template<class value_t, unsigned int value_n>
void QwInterpolator< value_t, value_n >::Cell ( const coord_t * coord,
unsigned int * cell_index ) const
inline

Return the cell index (unchecked)

Return the cell index (unchecked)

Parameters
coordPoint coordinates
cell_indexCell index of the point (reference)

Definition at line 652 of file QwInterpolator.h.

655{
656 // Get cell index and ignore local coordinates
658}
void Cell(const coord_t coord, unsigned int &cell_index, double &cell_local, const unsigned int dim) const
Return the cell index and local coordinates in one dimension (unchecked)

References Cell(), and f_cell_local.

+ Here is the call graph for this function:

◆ Cell() [2/4]

template<class value_t, unsigned int value_n>
void QwInterpolator< value_t, value_n >::Cell ( const coord_t * coord,
unsigned int * cell_index,
double * cell_local ) const
inline

Return the cell index and local coordinates (unchecked)

Return the cell index and local coordinates (unchecked)

Parameters
coordPoint coordinates
cell_indexCell index of the point (reference)
cell_localLocal coordinates in cell (reference)

Definition at line 636 of file QwInterpolator.h.

640{
641 // Get cell index and local coordinates in each dimension
642 for (unsigned int dim = 0; dim < fNDim; dim++)
644}

References Cell(), and fNDim.

+ Here is the call graph for this function:

◆ Cell() [3/4]

template<class value_t, unsigned int value_n>
void QwInterpolator< value_t, value_n >::Cell ( const coord_t coord,
unsigned int & cell_index,
double & cell_local,
const unsigned int dim ) const
inline

Return the cell index and local coordinates in one dimension (unchecked)

Return the cell index and local coordinates in one dimension (unchecked)

Parameters
coordPoint coordinate in one dimension
cell_indexCell index of the point (reference)
cell_localLocal coordinates in cell (reference)
dimDimension

Definition at line 612 of file QwInterpolator.h.

617{
618 // Normalize coordinate (positive, with integers on grid points)
619 double norm_coord = (coord - fMin[dim]) / fStep[dim];
620 // Split in fractional and integer part
621 double frac_part, int_part;
624 cell_index = static_cast<int>(int_part); // cast to integer
625 // Wrap index
626 if (fWrap[dim] > 0) cell_index %= (fSize[dim] - fWrap[dim]);
627}
std::vector< size_t > fWrap
Wrap around this coordinate.
std::vector< coord_t > fMin
Minimum in each dimension.
std::vector< coord_t > fStep
Step size in each dimension.
std::vector< size_t > fSize
Number of points in each dimension.

References fMin, fSize, fStep, and fWrap.

Referenced by Cell(), Cell(), Coord(), Index(), Linear(), Nearest(), and PrintCoverage().

+ Here is the caller graph for this function:

◆ Cell() [4/4]

template<class value_t, unsigned int value_n>
void QwInterpolator< value_t, value_n >::Cell ( unsigned int linear_index,
unsigned int * cell_index ) const
inline

Return the cell index based on the linearized index (unchecked)

Return the cell index based on the linearized index (unchecked)

Parameters
linear_indexLinearized index
cell_indexCell index (reference)

Definition at line 666 of file QwInterpolator.h.

669{
670 // This does not work with unsigned int, because that is always >= 0 and wraps
671 for (int dim = fNDim-1; dim >= 0; dim--) {
674 }
675}
std::vector< size_t > fExtent

References fExtent, and fNDim.

◆ Check() [1/3]

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::Check ( const coord_t * coord) const
inlineprivate

Check for boundaries with coordinate argument.

Check whether the point is inside the valid region

Parameters
coordPoint coordinates
Returns
True if the coordinates are in the valid region

Definition at line 528 of file QwInterpolator.h.

529{
530 for (unsigned int dim = 0; dim < fNDim; dim++)
531 if (fWrap[dim] == 0 && (coord[dim] < fMin[dim] || coord[dim] > fMax[dim]))
532 return false;
533 // Otherwise
534 return true;
535}
std::vector< coord_t > fMax
Maximum in each dimension.

References fMax, fMin, fNDim, and fWrap.

Referenced by GetValue(), InBounds(), Set(), Set(), and Set().

+ Here is the caller graph for this function:

◆ Check() [2/3]

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::Check ( const unsigned int * cell_index) const
inlineprivate

Check for boundaries with cell index argument.

Check whether the cell index is inside the valid region

Parameters
cell_indexCell index
Returns
True if the cell index is in the valid region

Definition at line 543 of file QwInterpolator.h.

544{
545 for (unsigned int dim = 0; dim < fNDim; dim++)
546 if (cell_index[dim] >= fSize[dim])
547 return false;
548 // Otherwise
549 return true;
550}

References fNDim, and fSize.

◆ Check() [3/3]

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::Check ( const unsigned int linear_index) const
inlineprivate

Check for boundaries with linearized index argument.

Check whether the linearized index is inside the valid region

Parameters
linear_indexLinearized index
Returns
True if the cell index is in the valid region

Definition at line 558 of file QwInterpolator.h.

559{
561 return false;
562 // Otherwise
563 return true;
564}
unsigned int fMaximumEntries
Maximum number of values.

References fMaximumEntries.

◆ Coord() [1/2]

template<class value_t, unsigned int value_n>
void QwInterpolator< value_t, value_n >::Coord ( const unsigned int * cell_index,
coord_t * coord ) const
inline

Return the coordinates based on the cell index (unchecked)

Return the coordinates based on the cell index (unchecked)

Parameters
cell_indexCell index
coordPoint coordinates (reference)

Definition at line 683 of file QwInterpolator.h.

686{
687 for (unsigned int dim = 0; dim < fNDim; dim++)
689}

References fMin, fNDim, and fStep.

Referenced by Coord(), and PrintCoverage().

+ Here is the caller graph for this function:

◆ Coord() [2/2]

template<class value_t, unsigned int value_n>
void QwInterpolator< value_t, value_n >::Coord ( const unsigned int linear_index,
coord_t * coord ) const
inline

Return the coordinates based on the linearized index (unchecked)

Return the coordinates based on the linearized index (unchecked)

Parameters
linear_indexLinearized index
coordPoint coordinates (reference)

Definition at line 697 of file QwInterpolator.h.

700{
703}
void Coord(const unsigned int *cell_index, coord_t *coord) const
Return the coordinates based on the cell index (unchecked)

References Cell(), Coord(), and f_cell_index.

+ Here is the call graph for this function:

◆ Get() [1/3]

template<class value_t = float, unsigned int value_n = 1>
value_t QwInterpolator< value_t, value_n >::Get ( const unsigned int * cell_index) const
inlineprivate

Get a single value by cell index (unchecked)

Definition at line 457 of file QwInterpolator.h.

457 {
458 if (value_n != 1) return 0; // only for one-dimensional values
459 return fValues[0][Index(cell_index)];
460 };
std::vector< value_t > fValues[value_n]
Table with pointers to arrays of values.
unsigned int Index(const coord_t *coord) const
Return the linearized index based on the point coordinates (unchecked)

References fValues, and Index().

Referenced by Linear(), NearestNeighbor(), and PrintCoverage().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ Get() [2/3]

template<class value_t = float, unsigned int value_n = 1>
value_t QwInterpolator< value_t, value_n >::Get ( const unsigned int index) const
inlineprivate

Get a single value by linearized index (unchecked)

Definition at line 463 of file QwInterpolator.h.

463 {
464 return fValues[0][index];
465 };

References fValues.

◆ Get() [3/3]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Get ( const unsigned int index,
value_t * value ) const
inlineprivate

Get a vector value by linearized index (unchecked)

Definition at line 467 of file QwInterpolator.h.

467 {
468 for (unsigned int i = 0; i < value_n; i++)
469 value[i] = fValues[i][index];
470 return true;
471 };

References fValues.

◆ GetCurrentEntries()

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::GetCurrentEntries ( ) const
inline

Get the current number of entries.

Definition at line 181 of file QwInterpolator.h.

181{ return fCurrentEntries; };
unsigned int fCurrentEntries
Number of values read in.

References fCurrentEntries.

◆ GetDataReductionFactor()

template<class value_t = float, unsigned int value_n = 1>
int QwInterpolator< value_t, value_n >::GetDataReductionFactor ( const unsigned int dim) const
inline

Get data reduction factor.

Definition at line 195 of file QwInterpolator.h.

196 { return fRedux.at(dim); }
std::vector< size_t > fRedux
Data reduction factor.

References fRedux.

◆ GetInterpolationMethod()

template<class value_t = float, unsigned int value_n = 1>
EQwInterpolationMethod QwInterpolator< value_t, value_n >::GetInterpolationMethod ( ) const
inline

Get the interpolation method.

Definition at line 213 of file QwInterpolator.h.

214 { return fInterpolationMethod; };
EQwInterpolationMethod fInterpolationMethod
Interpolation method.

References fInterpolationMethod.

◆ GetMaximum()

template<class value_t = float, unsigned int value_n = 1>
coord_t QwInterpolator< value_t, value_n >::GetMaximum ( const unsigned int dim) const
inline

Get maximum in dimension.

Definition at line 175 of file QwInterpolator.h.

175{ return fMax.at(dim); };

References fMax.

◆ GetMaximumEntries()

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::GetMaximumEntries ( ) const
inline

Get the maximum number of entries.

Definition at line 179 of file QwInterpolator.h.

179{ return fMaximumEntries; };

References fMaximumEntries.

◆ GetMinimum()

template<class value_t = float, unsigned int value_n = 1>
coord_t QwInterpolator< value_t, value_n >::GetMinimum ( const unsigned int dim) const
inline

Get minimum in dimension.

Definition at line 173 of file QwInterpolator.h.

173{ return fMin.at(dim); };

References fMin.

◆ GetStepSize()

template<class value_t = float, unsigned int value_n = 1>
coord_t QwInterpolator< value_t, value_n >::GetStepSize ( const unsigned int dim) const
inline

Get minimum in dimension.

Definition at line 177 of file QwInterpolator.h.

177{ return fStep.at(dim); };

References fStep.

◆ GetValue() [1/4]

template<class value_t = float, unsigned int value_n = 1>
value_t QwInterpolator< value_t, value_n >::GetValue ( const coord_t & coord) const
inline

Get the interpolated value at coordinate (zero if out of bounds)

Definition at line 332 of file QwInterpolator.h.

332 {
333 if (value_n != 1) return false; // only for one-dimensional values
334 if (fNDim != 1) return false; // only for one-dimensional grids
336 if (GetValue(&coord, &value)) return value;
337 else return 0; // interpolation failed
338 };
value_t GetValue(const coord_t &coord) const
Get the interpolated value at coordinate (zero if out of bounds)

References fNDim, and GetValue().

Referenced by GetValue(), GetValue(), and GetValue().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ GetValue() [2/4]

template<class value_t = float, unsigned int value_n = 1>
value_t QwInterpolator< value_t, value_n >::GetValue ( const coord_t * coord) const
inline

Get the interpolated value at coordinate (zero if out of bounds)

Definition at line 340 of file QwInterpolator.h.

340 {
341 if (value_n != 1) return 0; // only for one-dimensional values
343 if (GetValue(coord, &value)) return value;
344 else return 0; // interpolation failed
345 };

References GetValue().

+ Here is the call graph for this function:

◆ GetValue() [3/4]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::GetValue ( const coord_t * coord,
double & value ) const
inline

Get the interpolated value at coordinate (zero if out of bounds)

Definition at line 347 of file QwInterpolator.h.

347 {
348 if (value_n != 1) return false; // only for one-dimensional values
349 return GetValue(coord, &value);
350 };

References GetValue().

+ Here is the call graph for this function:

◆ GetValue() [4/4]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::GetValue ( const coord_t * coord,
double * value ) const
inline

Get the interpolated value at coordinate (zero if out of bounds)

Definition at line 352 of file QwInterpolator.h.

352 {
353 for (unsigned int i = 0; i < value_n; i++) value[i] = 0.0; // zero
354 if (! Check(coord)) return false; // out of bounds
355 value_t result[value_n]; // we need a local copy of type value_t
356 switch (fInterpolationMethod) {
357 // return interpolation status
358 case kMultiLinear:
359 if (! Linear(coord, result)) return false;
360 break;
361 case kNearestNeighbor:
362 if (! NearestNeighbor(coord, result)) return false;
363 break;
364 default: return false;
365 }
366 for (unsigned int i = 0; i < value_n; i++) value[i] = result[i]; // cast
367 return true;
368 };
bool Check(const coord_t *coord) const
Check for boundaries with coordinate argument.
bool Linear(const coord_t *coord, value_t *value) const
Linear interpolation (unchecked)
bool NearestNeighbor(const coord_t *coord, value_t *value) const
Nearest-neighbor (unchecked)

References Check(), fInterpolationMethod, kMultiLinear, kNearestNeighbor, Linear(), and NearestNeighbor().

+ Here is the call graph for this function:

◆ GetWrapCoordinate()

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::GetWrapCoordinate ( const unsigned int dim) const
inline

Get wrapping coordinate.

Definition at line 184 of file QwInterpolator.h.

185 { return fWrap.at(dim); }

References fWrap.

◆ InBounds()

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::InBounds ( const coord_t * coord) const
inline

Return true if the coordinate is in bounds.

Definition at line 248 of file QwInterpolator.h.

248 {
249 return Check(coord);
250 };

References Check().

+ Here is the call graph for this function:

◆ Index() [1/3]

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::Index ( const coord_t * coord) const
inline

Return the linearized index based on the point coordinates (unchecked)

Definition at line 409 of file QwInterpolator.h.

409 {
411 return Index(f_cell_index);
412 };

References Cell(), f_cell_index, and Index().

Referenced by Get(), Index(), Linear(), NearestNeighbor(), Set(), and Set().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ Index() [2/3]

template<class value_t, unsigned int value_n>
unsigned int QwInterpolator< value_t, value_n >::Index ( const unsigned int * cell_index) const
inline

Return the linearized index based on the cell indices (unchecked)

Return the linearized index based on the cell indices (unchecked)

Parameters
cell_indexIndex in each dimension
Returns
Linearized index

Definition at line 572 of file QwInterpolator.h.

574{
575 unsigned int linear_index = 0;
576 // Loop over dimensions
577 for (unsigned int dim = 0; dim < fNDim; dim++)
578 // ... and use the stored extents for every dimensions
580 return linear_index;
581}

References fExtent, and fNDim.

◆ Index() [3/3]

template<class value_t, unsigned int value_n>
unsigned int QwInterpolator< value_t, value_n >::Index ( const unsigned int * cell_index,
const unsigned int pattern ) const
inline

Return the linearized index based on the cell indices and offset (unchecked)

Return the linearized index based on the cell indices and offset (unchecked)

Parameters
cell_indexIndex in each dimension
patternBit pattern with offsets in each dimension
Returns
Linearized index

Definition at line 590 of file QwInterpolator.h.

593{
594 unsigned int linear_index = 0;
595 // Loop over dimensions
596 for (unsigned int dim = 0; dim < fNDim; dim++) {
597 // Add offset of one based on the bit at position dim
598 unsigned int offset = (pattern >> dim) & 0x1;
600 }
601 return linear_index;
602}

References fExtent, and fNDim.

◆ Linear()

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::Linear ( const coord_t * coord,
value_t * value ) const
inlineprivate

Linear interpolation (unchecked)

Perform the multi-dimensional linear interpolation (unchecked)

Parameters
coordPoint coordinates
valueInterpolated value (reference)

Definition at line 483 of file QwInterpolator.h.

486{
487 // Get cell and local normalized coordinates
489 // Initialize to zero
490 for (unsigned int i = 0; i < value_n; i++) value[i] = 0.0;
491 // Calculate the interpolated value
492 // by summing the 2^fNDim = 1 << fNDim neighbor points (1U is unsigned one)
493 for (unsigned int offset = 0; offset < (1U << fNDim); offset++) {
495 if (! Get(Index(f_cell_index,offset), neighbor)) return false;
496 // ... with appropriate weighting factors (1 - cell_local or cell_local)
497 double fac = 1.0;
498 for (unsigned int dim = 0; dim < fNDim; dim++)
499 fac *= ((offset >> dim) & 0x1)? f_cell_local[dim] : (1 - f_cell_local[dim]);
500 // ... for all vector components
501 for (unsigned int i = 0; i < value_n; i++)
502 value[i] += fac * neighbor[i];
503 }
504 return true;
505}
value_t Get(const unsigned int *cell_index) const
Get a single value by cell index (unchecked)

References Cell(), f_cell_index, f_cell_local, fNDim, Get(), and Index().

Referenced by GetValue().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ Nearest()

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::Nearest ( const coord_t * coord,
unsigned int * cell_index ) const
inlineprivate

Return the cell index closest to the coordinate (could be above) (unchecked)

Definition at line 437 of file QwInterpolator.h.

437 {
439 // Loop over all dimensions and add one if larger than 0.5
440 for (unsigned int dim = 0; dim < fNDim; dim++)
441 if (f_cell_local[dim] > 0.5) cell_index[dim]++;
442 };

References Cell(), f_cell_local, and fNDim.

Referenced by NearestNeighbor(), and Set().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ NearestNeighbor()

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::NearestNeighbor ( const coord_t * coord,
value_t * value ) const
inlineprivate

Nearest-neighbor (unchecked)

Perform the nearest-neighbor interpolation (unchecked)

Parameters
coordPoint coordinates
valueInterpolated value (reference)

Definition at line 513 of file QwInterpolator.h.

516{
517 // Get nearest cell
519 return Get(Index(f_cell_index), value);
520}
void Nearest(const coord_t *coord, unsigned int *cell_index) const
Return the cell index closest to the coordinate (could be above) (unchecked)

References f_cell_index, Get(), Index(), and Nearest().

Referenced by GetValue().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ PrintCoverage()

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::PrintCoverage ( const unsigned int dim)
inline

Print coverage map for all bins in one dimension.

Definition at line 218 of file QwInterpolator.h.

218 {
219 // Initialize coverage
220 unsigned int* cover = new unsigned int[fSize[dim]];
221 unsigned int* total = new unsigned int[fSize[dim]];
222 for (unsigned int i = 0; i < fSize[dim]; i++) {
223 cover[i] = 0;
224 total[i] = 0;
225 }
226 // Loop over all entries
228 for (unsigned int linear_index = 0; linear_index < fMaximumEntries; linear_index++) {
231 if (Get(linear_index,value) && value[0] != 0)
232 cover[f_cell_index[dim]]++; // non-zero field
233 }
234 // Print coverage
235 coord_t* coord = new coord_t[fNDim];
236 for (size_t i = 0; i < fSize[dim]; i++) {
237 f_cell_index[dim] = i;
239 mycout << "bin " << i << ", coord " << coord[dim] << ": "
240 << double(cover[i]) / double(total[i]) * 100 << "%"<< myendl;
241 }
242 delete[] coord;
243 delete[] cover;
244 delete[] total;
245 }

References Cell(), Coord(), f_cell_index, fMaximumEntries, fNDim, fSize, Get(), mycout, and myendl.

+ Here is the call graph for this function:

◆ ReadBinaryFile()

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::ReadBinaryFile ( const std::string & filename)
inline

Read the grid values from binary file.

Read the grid values from binary file (should be 64-bit safe, untested)

Parameters
filenameFile name
Returns
True if read successfully

Definition at line 840 of file QwInterpolator.h.

841{
843 if (! file.is_open()) return false;
844 // Informational message
845 mycout << "Reading binary file: ";
846 // Go to end and store length (could also use std::ios::ate)
847 file.seekg(0, std::ios::end);
848 // Go to begin and start reading template parameters
849 file.seekg(0, std::ios::beg);
850 uint32_t n, size;
851 file.read(reinterpret_cast<char*>(&n), sizeof(n));
852 file.read(reinterpret_cast<char*>(&size), sizeof(size));
853 if (n != value_n) return false; // not same dimensionality
854 if (size != sizeof(value_t)) return false; // not same type
855 // Read grid parameters
857 file.read(reinterpret_cast<char*>(&ndim), sizeof(ndim));
859 // Read sizes
860 for (unsigned int dim = 0; dim < fNDim; dim++) {
861 file.read(reinterpret_cast<char*>(&fMin[dim]),sizeof(fMin[dim]));
862 file.read(reinterpret_cast<char*>(&fMax[dim]),sizeof(fMax[dim]));
863 file.read(reinterpret_cast<char*>(&fStep[dim]),sizeof(fStep[dim]));
864 }
865 SetMinimumMaximumStep(fMin, fMax, fStep); // calculates fMaximumEntries
866 // Read values
868 file.read(reinterpret_cast<char*>(&maximum_entries),sizeof(maximum_entries));
869 if (maximum_entries != fMaximumEntries) return false; // not expected number of entries
870 int value_size = sizeof(value_t);
871 unsigned int entries = fValues[0].size();
872 for (unsigned int index = 0; index < entries; index++) {
873 // Read values
874 for (unsigned int i = 0; i < value_n; i++) {
875 file.read((char*)(&fValues[i][index]),value_size);
876 }
877 // Progress bar
878 if (index % (entries / 10) == 0)
879 mycout << index / (entries / 100) << "%" << QwLog::flush;
880 if (index % (entries / 10) != 0
881 && index % (entries / 40) == 0)
882 mycout << "." << QwLog::flush;
883 }
884 mycout << myendl;
885 file.close();
886 return true;
887}

References QwLog::flush(), fMax, fMaximumEntries, fMin, fNDim, fStep, fValues, mycout, myendl, SetDimensions(), and SetMinimumMaximumStep().

Referenced by QwInterpolator().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ ReadText()

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::ReadText ( std::istream & stream)
inline

Read the grid from text.

Read the grid values from a text stream

Parameters
streamInput stream
Returns
True if successfully read all values

Definition at line 746 of file QwInterpolator.h.

747{
748 // Informational message
749 mycout << "Reading text stream: ";
750 // Read the dimensions
751 unsigned int n;
752 stream >> fNDim >> n;
753 if (n != value_n) return false; // not same dimensionality
755 // Read the grid parameters
756 for (unsigned int dim = 0; dim < fNDim; dim++)
757 stream >> dim >> fMin[dim] >> fMax[dim] >> fStep[dim];
759 // Read the grid values
760 unsigned int entries;
761 stream >> entries;
762 for (unsigned int index = 0; index < entries; index++) {
763 // Read values
764 for (unsigned int i = 0; i < value_n; i++) {
765 stream >> fValues[i][index];
766 }
767 // Progress bar
768 if (index % (entries / 10) == 0)
769 mycout << index / (entries / 100) << "%" << QwLog::flush;
770 if (index % (entries / 10) != 0
771 && index % (entries / 40) == 0)
772 mycout << "." << QwLog::flush;
773 }
774 // Check for end of file
776 stream >> end;
777 // Informational message
778 mycout << myendl;
779 if (end == "end") return true;
780 else return false;
781}

References QwLog::flush(), fMax, fMin, fNDim, fStep, fValues, mycout, myendl, SetDimensions(), and SetMinimumMaximumStep().

Referenced by ReadTextFile().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ ReadTextFile()

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::ReadTextFile ( const std::string & filename)
inline

Read the grid from text file.

Definition at line 392 of file QwInterpolator.h.

392 {
393 std::ifstream file(filename.c_str());
394 if (! file.is_open()) return false;
395 if (! ReadText(file)) return false;
396 file.close();
397 return true;
398 };
bool ReadText(std::istream &stream)
Read the grid from text.

References ReadText().

+ Here is the call graph for this function:

◆ Set() [1/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const coord_t & coord,
const value_t & value )
inline

Set a single value at a coordinate (false if not possible)

Definition at line 256 of file QwInterpolator.h.

256 {
257 if (value_n != 1) return false; // only for one-dimensional values
258 if (fNDim != 1) return false; // only for one-dimensional grids
259 return Set(&coord, &value);
260 };
bool Set(const coord_t &coord, const value_t &value)
Set a single value at a coordinate (false if not possible)

References fNDim, and Set().

Referenced by Set(), Set(), Set(), Set(), Set(), and Set().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ Set() [2/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const coord_t * coord,
const value_t & value )
inline

Set a single value at a coordinate (false if not possible)

Definition at line 262 of file QwInterpolator.h.

262 {
263 if (value_n != 1) return false; // only for one-dimensional values
264 return Set(coord, &value);
265 };

References Set().

+ Here is the call graph for this function:

◆ Set() [3/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const coord_t * coord,
const value_t * value )
inline

Set a set of values at a coordinate (false if not possible)

Definition at line 267 of file QwInterpolator.h.

267 {
268 Nearest(coord, f_cell_index); // nearest cell
269 if (! Check(f_cell_index)) return false; // out of bounds
270 bool status = true;
271 bool written = false;
272 unsigned int linear_index;
273 for (unsigned int dim = 0; dim < fNDim; dim++) {
274 // skip dimensions that are not wrapped around
275 if (fWrap[dim] == 0) continue;
276 // FIXME only one wrapping coordinate correctly supported in Set()
277 if ((f_cell_index[dim] < fWrap[dim]) ||
278 (fSize[dim] - f_cell_index[dim] - 1 < fWrap[dim])) {
279 // there are equivalent grid points
280 for (size_t wrap = 0; wrap < fWrap[dim]; wrap++) {
281 // at the minimum
285 // at the maximum
286 f_cell_index[dim] = fSize[dim] - wrap - 1;
289 // set flag
290 written = true;
291 }
292 }
293 }
294 if (not written) {
295 // this is an unambiguous grid point
298 }
299 return status;
300 };

References Check(), f_cell_index, fNDim, fSize, fWrap, Index(), Nearest(), and Set().

+ Here is the call graph for this function:

◆ Set() [4/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const unsigned int * cell_index,
const value_t & value )
inline

Set a single value at a grid point (false if out of bounds)

Definition at line 317 of file QwInterpolator.h.

317 {
318 if (value_n != 1) return false; // only for one-dimensional values
319 return Set(cell_index, &value);
320 };

References Set().

+ Here is the call graph for this function:

◆ Set() [5/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const unsigned int * cell_index,
const value_t * value )
inline

Set a set of values at a grid point (false if out of bounds)

Definition at line 322 of file QwInterpolator.h.

322 {
323 if (! Check(cell_index)) return false; // out of bounds
324 return Set(Index(cell_index),value);
325 };

References Check(), Index(), and Set().

+ Here is the call graph for this function:

◆ Set() [6/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const unsigned int linear_index,
const value_t & value )
inline

Set a single value at a linearized index (false if not possible)

Definition at line 303 of file QwInterpolator.h.

303 {
304 if (value_n != 1) return false; // only for one-dimensional values
305 return Set(linear_index, &value);
306 };

References Set().

+ Here is the call graph for this function:

◆ Set() [7/7]

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::Set ( const unsigned int linear_index,
const value_t * value )
inline

Set a set of values at a linearized index (false if not possible)

Definition at line 308 of file QwInterpolator.h.

308 {
309 if (! Check(linear_index)) return false; // out of bounds
310 for (unsigned int i = 0; i < value_n; i++)
313 return true;
314 };

References Check(), fCurrentEntries, and fValues.

+ Here is the call graph for this function:

◆ SetDataReductionFactor() [1/3]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetDataReductionFactor ( const std::vector< unsigned int > & redux)
inline

Definition at line 204 of file QwInterpolator.h.

204 {
205 if (redux.size() != fNDim) return;
206 fRedux = redux;
207 }

References fNDim, and fRedux.

◆ SetDataReductionFactor() [2/3]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetDataReductionFactor ( const unsigned int dim,
const unsigned int redux )
inline

Set data reduction factor.

Definition at line 198 of file QwInterpolator.h.

199 { fRedux.at(dim) = redux; }

References fRedux.

Referenced by SetDataReductionFactor().

+ Here is the caller graph for this function:

◆ SetDataReductionFactor() [3/3]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetDataReductionFactor ( const unsigned int redux)
inline

Definition at line 200 of file QwInterpolator.h.

200 {
201 for (unsigned int dim = 0; dim < fNDim; dim++)
203 }
void SetDataReductionFactor(const unsigned int dim, const unsigned int redux)
Set data reduction factor.

References fNDim, and SetDataReductionFactor().

+ Here is the call graph for this function:

◆ SetDimensions()

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetDimensions ( const unsigned int ndim)
inline

Set the number of coordinate dimensions and resize vectors.

Definition at line 131 of file QwInterpolator.h.

131 {
132 if (ndim == 0) {
133 mycerr << "Dimensions of the grid should be strictly positive!" << myendl;
134 return;
135 }
136 fNDim = ndim;
137 fMin.resize(fNDim); fMax.resize(fNDim); fStep.resize(fNDim); fWrap.resize(fNDim);
138 fSize.resize(fNDim); fExtent.resize(fNDim+1);
139 };

References fExtent, fMax, fMin, fNDim, fSize, fStep, fWrap, mycerr, and myendl.

Referenced by QwInterpolator(), QwInterpolator(), ReadBinaryFile(), ReadText(), and SetMinimumMaximumStep().

+ Here is the caller graph for this function:

◆ SetInterpolationMethod()

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetInterpolationMethod ( const EQwInterpolationMethod method)
inline

Set the interpolation method.

Definition at line 210 of file QwInterpolator.h.

References fInterpolationMethod.

Referenced by QwInterpolator(), QwInterpolator(), and QwInterpolator().

+ Here is the caller graph for this function:

◆ SetMinimumMaximumStep() [1/2]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetMinimumMaximumStep ( const coord_t min,
const coord_t max,
const coord_t step )
inline

Set minimum, maximum, and step size to single values.

Definition at line 141 of file QwInterpolator.h.

References fNDim, and SetMinimumMaximumStep().

Referenced by QwInterpolator(), ReadBinaryFile(), ReadText(), and SetMinimumMaximumStep().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ SetMinimumMaximumStep() [2/2]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetMinimumMaximumStep ( const std::vector< coord_t > & min,
const std::vector< coord_t > & max,
const std::vector< coord_t > & step )
inline

Set minimum, maximum, and step size to different values.

Definition at line 147 of file QwInterpolator.h.

149 {
150 // Set the dimensionality
151 if (min.size() != fNDim) SetDimensions(min.size());
152 // Check the dimensionality and assign boundaries and step size vectors
153 if (min.size() != fNDim || min.size() != fNDim || step.size() != fNDim) return;
154 fMin = min; fMax = max; fStep = step;
155 fExtent[0] = 1;
156 for (size_t i = 0; i < fMin.size(); i++) {
157 coord_t int_part; // safer to use modf than a direct static cast
158 coord_t frac_part = modf((fMax[i] - fMin[i]) / fStep[i] + 1.0, &int_part);
159 fSize[i] = static_cast<unsigned int>(int_part) + (frac_part > 0.5 ? 1 : 0);
160 fExtent[i+1] = fExtent[i] * fSize[i];
161 }
162 // Try resizing to allocate memory and initialize with zeroes
163 fMaximumEntries = fExtent[fNDim]; fCurrentEntries = 0; // no entries read yet
164 for (unsigned int i = 0; i < value_n; i++) {
165 try {
166 fValues[i].resize(fMaximumEntries,0);
167 } catch (std::bad_alloc& e) {
168 mycerr << "QwInterpolator could not allocate memory for values!" << myendl;
169 }
170 }
171 };

References fCurrentEntries, fExtent, fMax, fMaximumEntries, fMin, fNDim, fSize, fStep, fValues, mycerr, myendl, and SetDimensions().

+ Here is the call graph for this function:

◆ SetWrapCoordinate() [1/2]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetWrapCoordinate ( const std::vector< size_t > & wrap)
inline

Definition at line 189 of file QwInterpolator.h.

189 {
190 if (wrap.size() != fNDim) return;
191 fWrap = wrap;
192 }

References fNDim, and fWrap.

◆ SetWrapCoordinate() [2/2]

template<class value_t = float, unsigned int value_n = 1>
void QwInterpolator< value_t, value_n >::SetWrapCoordinate ( const unsigned int dim,
const size_t wrap = 1 )
inline

Set wrapping coordinate.

Definition at line 187 of file QwInterpolator.h.

188 { fWrap.at(dim) = wrap; }

References fWrap.

◆ WriteBinaryFile()

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::WriteBinaryFile ( const std::string & filename) const
inline

Write the grid values to binary file.

Write the grid values to binary file (should be 64-bit safe, untested)

Integer data types can be stored differently on 32-bit and 64-bit systems. Fixed-length types uint32_t and u_int32_t are provided in stdint.h and sys/types.h, respectively. The floating point types float and double will always have a length of 32 and 64 bit, per the IEEE convention.

Parameters
filenameFile name
Returns
True if written successfully

Definition at line 795 of file QwInterpolator.h.

796{
798 if (! file.is_open()) return false;
799 // Informational message
800 mycout << "Writing binary file: ";
801 // Write template parameters
802 uint32_t n = value_n; // uint32_t has length of 32 bits on any system
803 uint32_t size = sizeof(value_t);
804 file.write(reinterpret_cast<const char*>(&n), sizeof(n));
805 file.write(reinterpret_cast<const char*>(&size), sizeof(size));
806 // Write grid parameters
808 file.write(reinterpret_cast<const char*>(&ndim), sizeof(ndim));
809 for (unsigned int dim = 0; dim < fNDim; dim++) {
810 file.write(reinterpret_cast<const char*>(&fMin[dim]),sizeof(fMin[dim]));
811 file.write(reinterpret_cast<const char*>(&fMax[dim]),sizeof(fMax[dim]));
812 file.write(reinterpret_cast<const char*>(&fStep[dim]),sizeof(fStep[dim]));
813 }
814 uint32_t entries = fValues[0].size();
815 file.write(reinterpret_cast<const char*>(&entries),sizeof(entries));
816 for (unsigned int index = 0; index < entries; index++) {
817 // Write values
818 for (unsigned int i = 0; i < value_n; i++) {
820 file.write(reinterpret_cast<const char*>(&value),sizeof(value));
821 }
822 // Progress bar
823 if (index % (entries / 10) == 0)
824 mycout << index / (entries / 100) << "%" << QwLog::flush;
825 if (index % (entries / 10) != 0
826 && index % (entries / 40) == 0)
827 mycout << "." << QwLog::flush;
828 }
829 mycout << myendl;
830 file.close();
831 return true;
832}

References QwLog::flush(), fMax, fMin, fNDim, fStep, fValues, mycout, and myendl.

+ Here is the call graph for this function:

◆ WriteText()

template<class value_t, unsigned int value_n>
bool QwInterpolator< value_t, value_n >::WriteText ( std::ostream & stream) const
inline

Write the grid as text.

Write the grid values to a text stream

Parameters
streamOutput stream

Definition at line 710 of file QwInterpolator.h.

711{
712 // Informational message
713 mycout << "Writing text stream: ";
714 // Write the dimensions
715 stream << fNDim << "\t" << value_n << std::endl;
716 // Write the grid parameters
717 for (unsigned int dim = 0; dim < fNDim; dim++)
718 stream << dim << "\t" << fMin[dim] << "\t" << fMax[dim] << "\t" << fStep[dim] << std::endl;
719 // Write the values
720 stream << fValues[0].size() << std::endl;
721 unsigned int entries = fValues[0].size();
722 for (unsigned int index = 0; index < entries; index++) {
723 // Write values
724 for (unsigned int i = 0; i < value_n; i++) {
725 stream << fValues[i][index] << "\t";
726 }
727 stream << std::endl;
728 // Progress bar
729 if (index % (entries / 10) == 0)
730 mycout << index / (entries / 100) << "%" << QwLog::flush;
731 if (index % (entries / 10) != 0
732 && index % (entries / 40) == 0)
733 mycout << "." << QwLog::flush;
734 }
735 stream << "end" << std::endl;
736 mycout << myendl;
737 return true;
738}

References QwLog::flush(), fMax, fMin, fNDim, fStep, fValues, mycout, and myendl.

Referenced by WriteTextFile(), and WriteTextScreen().

+ Here is the call graph for this function:
+ Here is the caller graph for this function:

◆ WriteTextFile()

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::WriteTextFile ( const std::string & filename) const
inline

Write the grid to text file.

Definition at line 377 of file QwInterpolator.h.

377 {
378 std::ofstream file(filename.c_str());
379 if (! file.is_open()) return false;
381 file.close();
382 return true;
383 };
bool WriteText(std::ostream &stream) const
Write the grid as text.

References WriteText().

+ Here is the call graph for this function:

◆ WriteTextScreen()

template<class value_t = float, unsigned int value_n = 1>
bool QwInterpolator< value_t, value_n >::WriteTextScreen ( ) const
inline

Write the grid to screen.

Definition at line 385 of file QwInterpolator.h.

385 {
387 return true;
388 };

References WriteText().

+ Here is the call graph for this function:

Field Documentation

◆ f_cell_index

template<class value_t = float, unsigned int value_n = 1>
unsigned int* QwInterpolator< value_t, value_n >::f_cell_index
private

◆ f_cell_local

template<class value_t = float, unsigned int value_n = 1>
double* QwInterpolator< value_t, value_n >::f_cell_local
private

Pre-allocated local coordinates.

Definition at line 122 of file QwInterpolator.h.

Referenced by Cell(), Linear(), Nearest(), QwInterpolator(), QwInterpolator(), QwInterpolator(), and ~QwInterpolator().

◆ fCurrentEntries

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::fCurrentEntries
private

Number of values read in.

Definition at line 115 of file QwInterpolator.h.

Referenced by GetCurrentEntries(), Set(), and SetMinimumMaximumStep().

◆ fExtent

template<class value_t = float, unsigned int value_n = 1>
std::vector<size_t> QwInterpolator< value_t, value_n >::fExtent
private

Linear extent between neighbor points in each dimension (e.g. for the least significant index this will be 1, for the next index the number of points in the first index, etc...)

Definition at line 109 of file QwInterpolator.h.

Referenced by Cell(), Index(), Index(), SetDimensions(), and SetMinimumMaximumStep().

◆ fInterpolationMethod

template<class value_t = float, unsigned int value_n = 1>
EQwInterpolationMethod QwInterpolator< value_t, value_n >::fInterpolationMethod
private

Interpolation method.

Definition at line 126 of file QwInterpolator.h.

Referenced by GetInterpolationMethod(), GetValue(), and SetInterpolationMethod().

◆ fMax

template<class value_t = float, unsigned int value_n = 1>
std::vector<coord_t> QwInterpolator< value_t, value_n >::fMax
private

Maximum in each dimension.

Definition at line 96 of file QwInterpolator.h.

Referenced by Check(), GetMaximum(), ReadBinaryFile(), ReadText(), SetDimensions(), SetMinimumMaximumStep(), WriteBinaryFile(), and WriteText().

◆ fMaximumEntries

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::fMaximumEntries
private

Maximum number of values.

Definition at line 117 of file QwInterpolator.h.

Referenced by Check(), GetMaximumEntries(), PrintCoverage(), ReadBinaryFile(), and SetMinimumMaximumStep().

◆ fMin

template<class value_t = float, unsigned int value_n = 1>
std::vector<coord_t> QwInterpolator< value_t, value_n >::fMin
private

Minimum in each dimension.

Definition at line 94 of file QwInterpolator.h.

Referenced by Cell(), Check(), Coord(), GetMinimum(), ReadBinaryFile(), ReadText(), SetDimensions(), SetMinimumMaximumStep(), WriteBinaryFile(), and WriteText().

◆ fNDim

template<class value_t = float, unsigned int value_n = 1>
unsigned int QwInterpolator< value_t, value_n >::fNDim
private

◆ fRedux

template<class value_t = float, unsigned int value_n = 1>
std::vector<size_t> QwInterpolator< value_t, value_n >::fRedux
private

Data reduction factor.

Definition at line 104 of file QwInterpolator.h.

Referenced by GetDataReductionFactor(), SetDataReductionFactor(), and SetDataReductionFactor().

◆ fSize

template<class value_t = float, unsigned int value_n = 1>
std::vector<size_t> QwInterpolator< value_t, value_n >::fSize
private

Number of points in each dimension.

Definition at line 100 of file QwInterpolator.h.

Referenced by Cell(), Check(), PrintCoverage(), Set(), SetDimensions(), and SetMinimumMaximumStep().

◆ fStep

template<class value_t = float, unsigned int value_n = 1>
std::vector<coord_t> QwInterpolator< value_t, value_n >::fStep
private

Step size in each dimension.

Definition at line 98 of file QwInterpolator.h.

Referenced by Cell(), Coord(), GetStepSize(), ReadBinaryFile(), ReadText(), SetDimensions(), SetMinimumMaximumStep(), WriteBinaryFile(), and WriteText().

◆ fValues

template<class value_t = float, unsigned int value_n = 1>
std::vector<value_t> QwInterpolator< value_t, value_n >::fValues[value_n]
private

Table with pointers to arrays of values.

Definition at line 112 of file QwInterpolator.h.

Referenced by Get(), Get(), Get(), ReadBinaryFile(), ReadText(), Set(), SetMinimumMaximumStep(), WriteBinaryFile(), and WriteText().

◆ fWrap

template<class value_t = float, unsigned int value_n = 1>
std::vector<size_t> QwInterpolator< value_t, value_n >::fWrap
private

Wrap around this coordinate.

Definition at line 102 of file QwInterpolator.h.

Referenced by Cell(), Check(), GetWrapCoordinate(), Set(), SetDimensions(), SetWrapCoordinate(), and SetWrapCoordinate().


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