6.1. Interpolators

For many applications, an interpolation routine will be an important part of your pre- and/or post-processing. Tudat comes with a number of interpolation routines, each of which is handled through the same class architecture.

As is the case for most other model that you create in Tudat, you provide settings to a create function to set up an interpolator. Before going into the details of this process, below is a list of the interpolation routines in Tudat:

  • One-dimensional:

    • Linear interpolation

    • Piece-wise constant interpolator

    • Cubic spline interpolation

    • Lagrange interpolation

    • Hermite interpolation

  • Multi-dimensional:

    • Multi-linear interpolation

Clearly, the first five options are interpolation routines with a single independent variable (i.e., they are one-dimensional), and the final one can be used for interpolation with an arbitrary number of independent variables. In this chapter, you will first introduced to how the one-dimensional interpolators are used, and then to how the multi-dimensional ones are implemented. These will be followed by a description of the various methods that can be used in case the requested independent variable for interpolation is outside the defined domain.

6.1.1. One-dimensional Interpolator Creation

Creating interpolators of a single independent variable is most easily done in Tudat by creating an InterpolatorSettings object, and passing this to the createOneDimensionalInterpolator function.

class InterpolatorSettings

Class containing settings to setup an interpolator.

// Load data.
std::map< double, Eigen::Vector6d > stateMap;
stateMap = ....

// Create interpolator
std::shared_ptr< interpolators::InterpolatorSettings > interpolatorSettings =
   std::make_shared< interpolators::InterpolatorSettings >( linear_interpolator );
std::shared_ptr< OneDimensionalInterpolator< double, Eigen::Vector6d > > interpolator =
     interpolators::createOneDimensionalInterpolator(
             stateMap, interpolatorSettings );

// Interpolate
Eigen::Vector6d interpolatedResult = interpolator->interpolate( 0.0 );

Note that if the data you interpolate is a different type, this should be reflected in how the interpolator is created. For instance, for std::map< double, Eigen::VectorXd > input, use

// Load data.
std::map< double, Eigen::VectorXd > stateMap;
stateMap = ....

// Create interpolator
std::shared_ptr< interpolators::InterpolatorSettings > interpolatorSettings =
   std::make_shared< interpolators::InterpolatorSettings >( linear_interpolator );
std::shared_ptr< OneDimensionalInterpolator< double, Eigen::VectorXd > > interpolator =
     interpolators::createOneDimensionalInterpolator(
             stateMap, interpolatorSettings );

// Interpolate
Eigen::VectorXd interpolatedResult = interpolator->interpolate( 0.0 );

In this example, the stateMap contains the data that is interpolated, using the double key (time) as independent variable and the Eigen::Vector6d value (state) as dependent variable. The interpolation type is linear, and the interpolatorSettings object is created by passing only the argument linear_interpolator.

The interpolator itself is of the type std::shared_ptr< OneDimensionalInterpolator< double, Eigen::Vector6d > >, where the double and Eigen::Vector6d denote the (in)dependent variables of the interpolation. For you application they may be different (for instance double and double, or long double and Time). The interpolator is created by calling the createOneDimensionalInterpolator function with the map containing the data and the interpolator settings as input.

The different interpolator types are handled in a similar manner:

  • Linear: Requires no additional information to the createOneDimensionalInterpolator and is defined by using the InterpolatorSettings base class (with linear_interpolator as argument).

  • Piece-wise constant: Requires no additional information to the createOneDimensionalInterpolator and is defined by using the InterpolatorSettings base class (with piecewise_constant_interpolator as argument).

  • Cubic spline: Requires no additional information to the createOneDimensionalInterpolator and is defined by using the InterpolatorSettings base class (with cubic_spline_interpolator as argument).

    Note

    Note that for the cubic spline implementation, natural boundary conditions are imposed (2nd derivatives at the boundaries equal to zero) and the first derivatives are continuous throughout the curve.

  • Hermite spline: This interpolation requires the definition of the values, as well as the derivatives of the curve at each of the nodes. The interpolator is defined by using the InterpolatorSettings base class (with hermite_spline_interpolator) as argument.

    // Load data for values at node points
std::map< double, Eigen::Vector6d > stateMap;
stateMap = ....

// Load data for first derivatives at node points
std::map< double, Eigen::Vector6d > stateDerivativeMap;
stateDerivativeMap = ....

// Load pair of default values for extrapolation (empty by default)
std::pair< DependentVariableType, DependentVariableType > defaultExtrapolationValues;
defaultExtrapolationValues = ...

// Create interpolator
std::shared_ptr< interpolators::InterpolatorSettings > interpolatorSettings =
 std::make_shared< interpolators::InterpolatorSettings >( hermite_spline_interpolator )
std::shared_ptr< OneDimensionalInterpolator< double, Eigen::Vector6d > > interpolator =
         interpolators::createOneDimensionalInterpolator(
                 stateMap, interpolatorSettings, defaultExtrapolationValues, stateDerivativeMap );

// Interpolate
Eigen::Vector6d interpolatedResult = interpolator->interpolate( 0.0 );

6.1.1.1. Lagrange Interpolator

This interpolation routine uses an \(n\)-th degree polynomial to approximate a function from \(n+1\) data points. In our implementation, you can use a large data set of \(m\) data points (with \(m > n\)) to generate a set of interpolating polynomials. When interpolating a data point, the interpolation routine will automatically select the polynomial where the requested data point lies between the two middle points, to prevent wild oscillations (which occur at the edge of the polynomial). At the boundaries of the full interval, a cubic spline interpolator is used.

To create a Lagrange interpolator, the number of data points used for each interpolating polynomial should be defined, using the dedicated derived class LagrangeInterpolatorSettings. The input argument for this class is the amount of points per polynomial. An example, when using 8 data points per polynomial, is described below:

class LagrangeInterpolatorSettings

Derived class used for the settings of a lagrange interpolator.

// Load data for values at node points
std::map< double, Eigen::Vector6d > stateMap;
stateMap = ....

// Create interpolator
std::shared_ptr< interpolators::InterpolatorSettings > interpolatorSettings =
        std::make_shared< interpolators::LagrangeInterpolatorSettings >( 8 )
std::shared_ptr< OneDimensionalInterpolator< double, Eigen::Vector6d > > interpolator =
                interpolators::createOneDimensionalInterpolator(
                        stateMap, interpolatorSettings );

// Interpolate
Eigen::Vector6d interpolatedResult = interpolator->interpolate( 0.0 );

6.1.2. Multi-dimensional Interpolator Creation

Currently, only one multi-dimensional interpolator is supported in Tudat. This is the multi-linear interpolator, which is created by also using the InterpolatorSettings class, and by specifying multi_linear_interpolator as first input. To create the interpolator object, a different function needs to be used. This is the createMultiDimensionalInterpolator function, whose inputs are as follows:

  • independentValues

    This is a vector of vectors, where each external vector corresponds to one independent variable, and the internal vector simply lists the independent variable points.

  • dependentData

    This variable is of type boost::multi_array, which is very similar to a multi-dimensional matrix in MATLAB and is generally created by reading data from a file by using the MultiArrayFileReader class.

  • interpolatorSettings

    This is simply a pointer to the InterpolatorSettings class described above.

6.1.3. Interpolation When Independent Variables are Out-of-range

In some instances, it can happen that the independent variable at which interpolation has to take place, falls outside the list of independent variables defined in the interpolator. In this case, the default behavior is extrapolation, but as you can imagine, this may lead to very wrong and strange behaviors, especially if a method like cubic spline is selected. For this reason, the user has the option of changing this so called, boundary handling method, or boundaryHandling, by specifying an extra input to the InterpolatorSettings object. The allowed methods for independent variable out-of-range are the following:

  • throw_exception_at_boundary: throw a runtime error

  • use_boundary_value: use the closest defined value for these conditions

  • extrapolate_at_boundary: apply extrapolation (default behavior)

  • use_default_value: use the value provided by defaultExtrapolationValue (defined as an input to the function createOneDimensionalInterpolator or createMultiDimensionalInterpolator)

Tip

Each of the methods above, with the exception of the first one (i.e., the one that gives the runtime error), comes with the option of also giving a warning in case the variable goes out-of-range. You just need to add _with_warning at the end of the enumeration. For example, you can use extrapolate_at_boundary_with_warning to extrapolate and give a warning.

For a one-dimensional interpolator, the input should be a single DependentVariableType value, whereas for a multi-dimensional one, you can also specify one method per each independent variable (i.e., for each dimension). Unfortunately, it is not yet possible to define a boundaryHandling method for each extremity of the independent variable range (e.g., give an error in case an independent variable is requested below lowest value and extrapolate if beyond the highest value).

In case you have specified use_default_value: or use_default_value_with_warning: as the boundary handling method, you should also specify the default value that the extrapolator has to output. As previously mentioned, this is done by adding an extra input to the createOneDimensionalInterpolator or createMultiDimensionalInterpolator functions. The format of this input, which is named defaultExtrapolationValue, depends on the dimensionality, and is as follows:

std::pair< DependentVariableType, DependentVariableType >

For one-dimenisonal interpolators, this variable is defined as a pair. The first value of the pair (accessed with the method .first) denotes the default extrapolation value in case the independent variable is requested below the lower domain boundary. The second value (.second), on the other hand, is output in case the independent variable is asked above the higher domain boundary.

std::vector< std::pair< DependentVariableType, DependentVariableType > >

For multi-dimensional interpolators, it is defined as a vector of pairs. In this case, the pairs are used in the same manner as above, but a pair is defined for each independent variable.