Kameleon Interpolator Classes

Model-Independent Classes

Interpolator

class ccmc::Interpolator
TODO: Brief description of Interpolator class.

TODO: Full description of Inteprolator class

Public Functions

Interpolator()
Default constructor

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2) = 0
Interpolates the variable specified at position (c0,c1,c2).

TODO: add detailed description

virtual float getConversionFactor(const std::string &variable)
Return
Parameters
variable -
virtual float getConversionFactor(const long &variable_id)
Return
Parameters
variable_id -
virtual ~Interpolator()
Destructor

virtual void convertCoordinates(const std::string &source, const std::string &dest, const long time_et, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Convert from source coordinates to destination coordinates results placed in dc0, dc2, dc2

KameleonInterpolator

class ccmc::KameleonInterpolator
Provides interpolation to derived variables for models accessed via the Kameleon class.

The KameleonInterpolator class inherits from the Interpolator, but computes model-dependent derived variables. Currently these are all hand-coded.

Public Functions

KameleonInterpolator(Model *modelReader)
Parameters
model -
KameleonInterpolator(Model *modelReader, const std::string &preferred_coordinates)
OpenGGCM and BATSRUS: c0,c1,c2 corresponds to x,y,z, respectively. ENLIL and MAS: c0,c1,c2 corresponds to r,phi(latitude), theta(longitude), respectively

Parameters
variable -
c0 -
c1 -
c2 -
long getEphemTime()
Get ephemeris time

void setEphemTime(long time)
Set Ephemeris time in seconds past 12h on 1 January 2000

Parameters
time -
- seconds past 12h January 2000

void setEphemTime()
Use the default times stored in the file

void convertCoordinates(Position *v_in, Position *v_out)
Assume conversion from preferred to model coordinates

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable -
c0 -
c1 -
c2 -
virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
OpenGGCM and BATSRUS: c0,c1,c2 corresponds to x,y,z, respectively. ENLIL and MAS: c0,c1,c2 corresponds to r,phi(latitude), theta(longitude), respectively

Parameters
variable -
c0 -
c1 -
c2 -
dc0 -
dc1 -
dc2 -
variable -
c0 -
c1 -
c2 -
dc0 -
dc1 -
dc2 -
return -
virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2)
OpenGGCM and BATSRUS: c0,c1,c2 corresponds to x,y,z, respectively. ENLIL and MAS: c0,c1,c2 corresponds to r,phi(latitude), theta(longitude), respectively

Return
Parameters
variable_id -
c0 -
c1 -
c2 -
Interpolates using the base model interpolators. Does not perform derived variable calculations.

variable_id -
c0 -
c1 -
c2 -
virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
A variable_id won’t work well, since derived variables can be requested, which do not exist in the data.

OpenGGCM and BATSRUS: c0,c1,c2 maps to x,y,z, respectively. ENLIL and MAS: c0,c1,c2 maps to r,theta(latitude), phi(longitude), respectively

Parameters
variable -
c0 -
c1 -
c2 -
dc0 -
dc1 -
dc2 -

TimeInterpolator

class ccmc::TimeInterpolator
The TimeInterpolator collects multiple Kameleon objects and permits time interpolation between them.

Public Functions

TimeInterpolator()
virtual ~TimeInterpolator()
void addTimestep(const std::string &filename)
void removeTimestep(const Time &time)
std::vector<Time> getTimesteps()
Interpolator *getInterpolator(const Time &time)
Kameleon *getKameleon(const Time &time)
bool manageMemory(const Time &time, const std::vector<std::string> &variable)
Manages the memory and the opened files based on the time and variable specified.

bool manageMemory(const Time &time, const char *variable[], int numVars)
Manages the memory and the opened files based on the time and variable specified.

void closeAndCleanupMemory()
float interpolate(const Time &time, const std::string &variable, const float &c0, const float &c1, const float &c2)
float interpolate(const Time &time, int variable, const float &c0, const float &c1, const float &c2)
void clearAll()
Public Static Functions

Time parseTime(const std::string &timeString)
std::string toString(const Time &time)

Point classes

Point and Point3f are two different implementations of a point class. Point uses templating and leaves the operator implementation to the Vector subclass. Point3f is for floats only, but includes the operators. Point3f is used extensively by the Fieldline class.

template <class T>
class ccmc::Point
TODO: Brief description of Point class.

TODO: Full description of Point class

Public Functions

Point()
Default constructor

Point(T &c0, T &c1, T &c2)
void setComponents(T &c1, T &c2, T &c3)
const T &c0()
const T &c1()
const T &c2()
std::string toString()
virtual ~Point()
class ccmc::Point3f
Public Types

enum Coordinates
Values:

SPHERICAL
CARTESIAN
Public Functions

Point3f(const float &component1, const float &component2, const float &component3)
TODO: finish documentation

Parameters
component1 -
component2 -
component3 -
Point3f(const float &component1, const float &component2, const float &component3, Coordinates c)
TODO: finish documentation

Parameters
component1 -
component2 -
component3 -
Point3f(const Point3f &p)
Point3f()
TODO: finish documentation

Point3f(Coordinates c)
TODO: finish documentation

void normalize()
TODO: finish documentation

std::string toString() const
TODO: finish documentation

float magnitude()
Computes the magnitude of the cartesian vector

float distance(const Point3f &p) const
compute the distance between this point, and Point p TODO: finish documentation

Point3f operator+(const Point3f &p) const
TODO: finish documentation

Point3f operator-(const Point3f &p) const
Minus operator where sender is on the right TODO: test this

Point3f operator*(float value) const
TODO: finish documentation

Point3f operator*(double value) const
TODO: finish documentation

void setCoordinates(Coordinates c)
TODO: finish documentation

Point3f::Coordinates getCoordinates()
TODO: finish documentation

Point3f getCartesian()
virtual ~Point3f()
TODO: finish documentation

Public Members

float component1
float component2
float component3

Vector

template <class T>
class ccmc::Vector
Inherits from Point class

Public Functions

Vector()
Vector(T c0, T c1, T c2)
T length()
Vector<T> operator+(const Vector<T> &v) const
Vector<T> operator-(const Vector<T> &v) const
Vector<T> operator*(const T &value) const
Vector<T> operator/(const T &value) const
void operator+=(const Vector<T> &v)
void operator-=(const Vector<T> &v)
void operator*=(const T &value)
void operator/=(const T &value)
void norm()
virtual ~Vector()
Public Static Functions

T dot(Vector<T> &a, Vector<T> &b)
T dot(Vector<T> *a, Vector<T> *b)
Vector<T> cross(Vector<T> &a, Vector<T> &b)
Vector<T> cross(Vector<T> *a, Vector<T> *b)
void cross(Vector<T> &c, Vector<T> &a, Vector<T> &b)
T triple(Vector<T> &a, Vector<T> &b, Vector<T> &c)
Vector<T> norm(Vector<T> &a)
Vector<T> norm(Vector<T> *a)
T angle(Vector<T> &a, Vector<T> &b)
void angle(T &angle, Vector<T> &a, Vector<T> &b)
Public Static Attributes

Vector<T> dummy

Polyhedron

This class implements Spherical Barycentric Coordinates for interpolation. In the limit that a query point approaches one of the faces, the scheme converges to Mean Value Coordinates.

template <class T>
class ccmc::Polyhedron
Class to perform interpolations on arbitrary convex polyedra with planar faces.

Polyhedron is a class that can represent any arbitary polyhedron where the faces are planar. Each face can have any number of vertices. This class also implements Spherical Barycentric Coordinates, which may be used for linear interpolation for points inside or on the polyhedron’s boundary. Points outside the polyhedron are given zero weight.

Public Functions

Polyhedron()
Constructor initializes interpolation variables.

void setPositions(std::vector<Vector<T>> const &pos)
Sets positions as a set of vectors.

void setLoopsFaces(std::vector<int> const &loop, std::vector<int> const &face)
initializes loops and faces of polyhedron

void setNeighbors(std::vector<int> const &neighbor)
initializes neighbors cell indexes. If a face has no neighbor, set its neighbor value to currentPolyhedron;

Vector<T> centroid()
Calculate this polyhedron’s centroid (average of positions).

Vector<T> faceCentroid(int face)
Calculates face centroid as average of face vertices

Return
face’s centroid.
Parameters
face -
which face to average

void setBarycentricCoordinates(Vector<T> const &point)
Given a starting point, calculates Spherical Barycentric Coordinates, such that every position is given a weight in [0,1]. The weights are stored in lambda.

If the point is outside the polyhedron, the weights will be set to zero. The weights will be linearly precise: The sum of each weight times the corresponding position yields the original point. point = SUM_i(lambda_i x position[i])

See Langer et. al, Spherical Barycentric Coordinates, 2006 Eurographics Symposium on Geometry

Parameters
point -
3D position for which to calculate barycentric weights.

void setFloaterCoordinates(int face, Vector<T> const &point)
If the point lies on the boundary, spherical barycentric coordinates approach Floater coordinates (citation needed), which are implemented here.

void setLinearInterpolants(int v1, int v2, Vector<T> const &point)
If the point lies on an edge, spherical barycentric coordinates and Floater coordinates approach linear interpolation between the vertices that define the edge.

void saveAsDXObject(string filename)
Save as a human-readable Open-DX data object

Vector<T> getFaceVectorError()
Check that face vectors sum to zero according to Stoke’s theorem

void setCurrentPolyhedron(int a)
Sets currentPolyhedron

virtual Polyhedron<T> *getNextPolyhedron()
Points to the next polyhedron. This prototype just returns a pointer to itself.

Vector<T> maxDistanceToCentroid()
Calculates maximum distance from vertices to centroid

Vector<T> minDistanceToCentroid()
Calculates minimum distance from vertices to centroid

int getClosestFace(Vector<T> const &point)
Find closest face by comparing distances to face centroids

Vector<T> testLinearity(Vector<T> const &point)
Test the property that the sum of the vertices times their weights recovers the interpolation point.

bool isPointInside(Vector<T> const &point)
Checks if a new point is inside this polyhedron. Stores the result in isInside.

virtual int getType()
Returns the type of Polyhedron this is. This prototype just returns 0.

void merge(Polyhedron<T> const *p)
Merges this polyhedron with another. Use this to create a group of polyhedra to render as a single object using saveAsDXObject.

Public Members

std::vector<T> lambda
barycentric weights associated with each vertex, determined by setBarycentricCoordinates. Each weight is in [0,1]

std::vector<Vector<T>> positions
3D positions of polyhedron’s vertices

std::vector<int> globalVertex
global index of polyhedron’s vertices

int closestFace
index of face closest to the last interpolation query point

int currentPolyhedron
global index of this polyhedron

std::vector<int> loops
each loop defines a face of polyhedron. loops must index into positions, with face normals determined by right-hand-rule.

std::vector<int> faces
faces must index into the start of each loop in loops, except the last face element will be the index of the last loop element

bool isInside
true if the last interpolation was inside the polyhedron

std::vector<int> neighbors
global index of neighboring polyhedra, matching faces.

Public Static Attributes

T piOver2
used for zeroing vectors

Nanoflann Classes

The Nanoflann header-only library is used for aiding in search queries.

template <typename T>
struct PointCloud
Public Functions

size_t kdtree_get_point_count() const
T kdtree_distance(const T *p1, const size_t idx_p2, size_t size) const
T kdtree_get_pt(const size_t idx, int dim) const
template <class BBOX>
bool kdtree_get_bbox(BBOX &bb) const
Public Members

std::vector<Point> pts
struct Point
Public Members

T x
T y
T z
template <typename DistanceType, typename IndexType = size_t, typename CountType = size_t>
class nanoflann::KNNResultSet
The KNNResultSet is a container for storing the results of a KD-tree query.

Public Functions

KNNResultSet(CountType capacity_)
void init(IndexType *indices_, DistanceType *dists_)
CountType size() const
bool full() const
void addPoint(DistanceType dist, IndexType index)
DistanceType worstDist() const
template <class num_t>
class NanoKdTree
Public Types

typedef KDTreeSingleIndexAdaptor<L2_Simple_Adaptor<num_t, PointCloud<num_t>>, PointCloud<num_t>, 3> my_kd_tree
Public Functions

void build()
void nearest(const num_t *query_point, KNNResultSet<num_t> &resultSet)
~NanoKdTree()
Public Members

PointCloud<num_t> cloud
my_kd_tree *tree

Fieldline

class ccmc::Fieldline
Public Functions

Fieldline()
TODO: finish documentation

Fieldline(int initialSize)
Parameters
initialSize -
~Fieldline()
TODO: finish documentation

void insertPointData(const Point3f &p, const float &d)
TODO: finish documentation

Parameters
p -
d -
void insertVectorData(const Point3f &vector)
void removePoint(int index)
TODO: finish documentation

Parameters
index -
Fieldline reverseOrder()
TODO: finish documentation

Return
void reverseOrderInPlace()
TODO: finish documentation

Return
const std::vector<Point3f> &getPositions()
TODO: finish documentation

Return
const std::vector<float> &getData()
TODO: finish documentation

Return
int size()
TODO: finish documentation

Return
const Point3f &getPosition(int i)
TODO: finish documentation

Return
Parameters
i -
float getData(int i)
TODO: finish documentation

Return
Parameters
i -
Point3f getStartPoint()
TODO: finish documentation

int getStartIndex()
void reserve(int size)
Parameters
size -
void setStartPoint(Point3f p)
TODO: finish documentation

Parameters
startPoint -
void setStartIndex(int index)
void setVariable(std::string variable)
const std::string &getVariable()
const std::vector<float> &getDs()
Calculate the forward difference elements for a field line with ordered positions. Output has length fieldline.size()-1 TODO: Add backward and higher-order differencing

const std::vector<Point3f> &getElements()
const Point3f &getElement(int i)
const std::vector<float> &integrate()
Calculate the integral of ds*values over the length of the field line

const std::vector<float> &integrateVector()
Calculate the integral of (vector) dot dl over the length of the field line

const std::vector<float> &derivative()
Calculate the derivative of d values/ ds over the length of the field line

const std::vector<float> &measure()
Measure the length of the field line up to point i

float getLength(int i)
Get the length up to position i

float getIntegral(int i)
Get the integral up to position i

Fieldline interpolate(int option, int Npoints)
const std::vector<int> &getNearest()
const std::vector<float> &getTlocal()
void minmax()
Public Members

int mincount
int maxcount
std::vector<int> minima
std::vector<int> maxima
int GlobMinIndex
int GlobMaxIndex

Model-specific Interpolators

All Kameleon Interpolators inherit from the Interpolator class

Adapt3D

class ccmc::Adapt3DInterpolator
TODO: brief description of Adapt3DInterpolator class.

TODO: full description of Adapt3DInterpolator class

Public Functions

Adapt3DInterpolator(Model *modelReader)
Parameters
modelReader -
Pointer to the Model object containing the appropriate variable maps. Adapt3DInterpolator should be returned by a Adapt3D::createNewInterpolator() call.

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable -
c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Interpolation method. Note that using the variable ID is significantly faster than using the variable string.

Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable -
The input variable.

c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable_id -
c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Interpolation method. Note that using the variable ID is significantly faster than using the variable string.

Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable_id -
A long representing the variable ID.

c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual ~Adapt3DInterpolator()
Destructor

BATSRUS

class ccmc::BATSRUSInterpolator
TODO: brief description of BATSRUSInterpolator class.

TODO: full description of BATSRUSInterpolator class

Public Functions

BATSRUSInterpolator(Model *modelReader)
Parameters
modelReader -
Pointer to the Model object containing the appropriate variable maps. BATSRUSInterpolator should be returned by a BATSRUS::createNewInterpolator() call.

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable -
c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Interpolation method. Note that using the variable ID is significantly faster than using the variable string.

Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable -
The input variable.

c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable_id -
c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Interpolation method. Note that using the variable ID is significantly faster than using the variable string.

Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable_id -
A long representing the variable ID.

c0 -
X component of the position

c1 -
Y component of the position

c2 -
Z component of the position

dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual ~BATSRUSInterpolator()
Destructor

LFM

class ccmc::LFMInterpolator
LFM interpolator class, combining kd-tree with spherical barycentric coordinates.

The LFMInterpolator makes use of the nanoflann library’s kd-tree implementation NanoKdTree for fast lookup of LFM positions. Once the nearest cell center is found, interpolation is performed using spherical barycentric coordinates via the Polyhedron class.

TODO: Move kd-tree initialization to LFM object instantiation.

Public Functions

LFMInterpolator(Model *modelReader)
Constructor

Parameters
modelReader -
pointer to LFM model object

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable -
c0 -
c1 -
c2 -
virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable -
c0 -
c1 -
c2 -
dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable_id -
c0 -
c1 -
c2 -
virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable_id -
c0 -
c1 -
c2 -
dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual float getConversionFactor(const std::string &variable)
Return
Parameters
variable -
virtual float getConversionFactor(const long &variable_id)
Return
Parameters
variable_id -
virtual ~LFMInterpolator()
Destructor

void setPolyhedralCells()
This function creates and stores the polyhedra used for interpolation. Each polyhedron has vertices set by LFM cell centers. In total, there are (ni-1)*(nj-1)*nk hexahedra (6 faces), 2*(ni-1) axis cells having nk+2 faces, 1 inner boundary cell having (nj-1)*nk + 2 faces. This function also calculates polyhedral cell centers, which are put in a point cloud and used for kd-tree searches

Polyhedron<float> *getCell(Vector<float> point)
Return a pointer to the cell containing the query point.

Polyhedron<float> getInterpolationPolys()
Constructs a single polyhedron out of all the polys used in interpolation

Public Members

PointCloud<float> cell_centers
size ni,nj,nk - nonstatic: only used once for polyhedra

PointCloud<float> helper_nodes
size ni-1,nj-1,nk (+axis later) - nonstatic: copied to lfmtree.cloud

ENLIL

class ccmc::ENLIL
TODO: Brief description of ENLIL class.

TODO: Full description of ENLIL class

Public Functions

virtual long open(const std::string &filename)
Opens a file.

Opens a file and performs any necessary initialization required to work with the data.

Parameters
filename -
virtual Interpolator *createNewInterpolator()
Returns an Interpolator object for the currently opened file.

This returns an Interpolator object that contains all the necessary local variables required to interpolate independent of any other Interpolator object. The pointer must be deleted from the calling program.

Return
A pointer to an Interpolator object.
virtual const std::vector<std::string> getLoadedVariables()
Returns the list of variables that have been loaded into memory, using the loadVariable or loadVectorVariable methods.

MAS

class ccmc::MASInterpolator
TODO: Brief description of MASInterpolator class.

TODO: Full description of MASInteprolator class

Public Functions

MASInterpolator(Model *model)
Parameters
model -
virtual float interpolate(const std::string &variable, const float &r, const float &lat, const float &lon)
Return
Parameters
variable -
r -
lat -
lon -
virtual float interpolate(const std::string &variable, const float &r, const float &lat, const float &lon, float &dr, float &dlat, float &dlon)
Return
Parameters
variable -
r -
lat -
lon -
dr -
dlat -
dlon -
virtual float interpolate(const long &variable_id, const float &r, const float &lat, const float &lon)
Return
Parameters
variableID -
r -
lat -
lon -
virtual float interpolate(const long &variable_id, const float &r, const float &lat, const float &lon, float &dr, float &dlat, float &dlon)
Return
Parameters
variableID -
r -
lat -
lon -
dr -
dlat -
dlon -

OpenGGCM

class ccmc::OpenGGCMInterpolator
TODO: Brief description of OpenGGCMInterpolator class.

TODO: Full description of OpenGGCMInterpolator class

Public Functions

OpenGGCMInterpolator(Model *modelReader)
Constructor

Parameters
modelReader -
virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable -
c0 -
c1 -
c2 -
virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable -
c0 -
c1 -
c2 -
dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2)
Return
Parameters
variable_id -
c0 -
c1 -
c2 -
virtual float interpolate(const long &variable_id, const float &c0, const float &c1, const float &c2, float &dc0, float &dc1, float &dc2)
Return
The interpolated value at position (c0,c1,c2) with deltas (dc0,dc1,dc2)
Parameters
variable_id -
c0 -
c1 -
c2 -
dc0 -
Reference to a variable to store the delta for component 0

dc1 -
Reference to a variable to store the delta for component 1

dc2 -
Reference to a variable to store the delta for component 2

virtual float getConversionFactor(const std::string &variable)
Creates the cell with the necessary values to interpolate the correct value at position (c0,c1,c2).

6_____7
/|    /|
4_|___5 |
| |   | |
| 2 - - 3
|/    |/
0 - - 1
Return
The cell values necessary for inteprolating the value at position (c0,c1,c2).
Return
Parameters
variable -
c0 -
c1 -
c2 -
Parameters
variable -
virtual float getConversionFactor(const long &variable_id)
Return
Parameters
variable_id -
virtual ~OpenGGCMInterpolator()
Destructor

SWMF Ionosphere

class ccmc::SWMFIonoInterpolator
Public Functions

virtual float interpolate(const std::string &variable, const float &c0, const float &c1, const float &c2)
Interpolates the variable specified at position (c0,c1,c2).

TODO: add detailed description

virtual float interpolate(const long &variable_id, const float &r, const float &lat, const float &lon)
Parameters
r -
This parameter is unused, since the grid is 2D.

lat -
The latitude component of the position.

lon -
The longitude component of the position.