diff --git a/CMakeLists.txt b/CMakeLists.txt index 1c09552425d74a6197d9339bf36d700fb551ef98..77c48a3f34bb659ebf97d085bf9063bd164a5164 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -33,10 +33,12 @@ set(ALE_INSTALL_INCLUDE_DIR "include/ale") add_library(ale SHARED ${CMAKE_CURRENT_SOURCE_DIR}/src/ale.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/Rotation.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/src/States.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/Isd.cpp ${CMAKE_CURRENT_SOURCE_DIR}/src/Util.cpp) set(ALE_HEADERS "${ALE_BUILD_INCLUDE_DIR}/ale.h" "${ALE_BUILD_INCLUDE_DIR}/Rotation.h" + "${ALE_BUILD_INCLUDE_DIR}/States.h" "${ALE_BUILD_INCLUDE_DIR}/Isd.h" "${ALE_BUILD_INCLUDE_DIR}/Util.h") set_target_properties(ale PROPERTIES diff --git a/include/States.h b/include/States.h new file mode 100644 index 0000000000000000000000000000000000000000..c87cc7f059caa4dea35c0a1a244d9379c096afb3 --- /dev/null +++ b/include/States.h @@ -0,0 +1,133 @@ +#ifndef ALE_STATES_H +#define ALE_STATES_H + +#include<vector> +#include <stdexcept> +#include <gsl/gsl_interp.h> +#include <ale.h> + +namespace ale { + +/** A 3D cartesian vector */ +struct Vec3d { + double x; + double y; + double z; + + // Accepts an {x,y,z} vector + Vec3d(const std::vector<double>& vec) { + if (vec.size() != 3) { + throw std::invalid_argument("Input vector must have 3 entries."); + } + x = vec[0]; + y = vec[1]; + z = vec[2]; + }; + + Vec3d(double x, double y, double z) : x(x), y(y), z(z) {}; + Vec3d() : x(0.0), y(0.0), z(0.0) {}; +}; + +/** A state vector with position and velocity*/ +struct State { + Vec3d position; + Vec3d velocity; + + // Accepts a {x, y, z, vx, vy, vz} vector + State(const std::vector<double>& vec) { + if (vec.size() != 6) { + throw std::invalid_argument("Input vector must have 6 entries."); + } + position = {vec[0], vec[1], vec[2]}; + velocity = {vec[3], vec[4], vec[5]}; + }; + + State(ale::Vec3d position, ale::Vec3d velocity) : position(position), velocity(velocity) {}; + State() {}; +}; + +class States { + public: + // Constructors + States(); + States(const std::vector<double>& ephemTimes, const std::vector<ale::Vec3d>& positions, + int refFrame=1); + + States(const std::vector<double>& ephemTimes, const std::vector<ale::Vec3d>& positions, + const std::vector<ale::Vec3d>& velocities, int refFrame=1); + + States(const std::vector<double>& ephemTimes, const std::vector<ale::State>& states, + int refFrame=1); + + ~States(); + + // Getters + std::vector<ale::State> getStates() const; //! Returns state vectors (6-element positions&velocities) + std::vector<ale::Vec3d> getPositions() const; //! Returns the current positions + std::vector<ale::Vec3d> getVelocities() const; //! Returns the current velocities + std::vector<double> getTimes() const; //! Returns the current times + int getReferenceFrame() const; //! Returns reference frame as NAIF ID + bool hasVelocity() const; //! Returns true if any velocities have been provided + + /** + * Returns a single state by interpolating state. + * If the Cache has been minimized, a cubic hermite is used to interpolate the + * position and velocity over the reduced cache. + * If not, a standard gsl interpolation will be done. + * + * @param time Time to get a value at + * @param interp Interpolation type to use. Will be ignored if cache is minimized. + * + * @return ale::State + */ + ale::State getState(double time, ale::interpolation interp=LINEAR) const; + + /** Gets a position at a single time. Operates the same way as getState() **/ + ale::Vec3d getPosition(double time, ale::interpolation interp=LINEAR) const; + + /** Gets a velocity at a single time. Operates the same way as getState() **/ + ale::Vec3d getVelocity(double time, ale::interpolation interp=LINEAR) const; + + /** Returns the first ephemeris time **/ + double getStartTime(); + + /** Returns the last ephemeris time **/ + double getStopTime(); + + /** + * Perform a cache reduction. After running this, getStates(), getPositions(), + * and getVelocities() will return vectors of reduced size, and getState(), + * getPosition(), and getVelocity() will + * returns values interpolated over the reduced cache using a cubic hermite spline + * + * Adapted from Isis::SpicePosition::reduceCache(). + * + * @param tolerance Maximum error between hermite approximation and original value. + */ + States minimizeCache(double tolerance=0.01); + + private: + + /** + * Calculates the points (indicies) which need to be kept for the hermite spline to + * interpolate between to mantain a maximum error of tolerance. + * + * Adapted from Isis::SpicePosition::HermiteIndices. + * + * @param tolerance Maximum error between hermite approximation and original value. + * @param indexList The list of indicies that need to be kept. + * @param baseTime Scaled base time for fit + * @param timeScale Time scale for fit. + * + * @return std::vector<int> + */ + std::vector<int> hermiteIndices(double tolerance, std::vector <int> indexList, + double baseTime, double timeScale); + std::vector<ale::State> m_states; //! Represent at states internally to keep pos, vel together + std::vector<double> m_ephemTimes; //! Time in seconds + int m_refFrame; //! Naif IDs for reference frames + }; +} + +#endif + diff --git a/include/ale.h b/include/ale.h index ee268126c0aac3a4a20ef3e9cb995a3cfb006649..64c5ab41b4bbe336da2144d1dfea7605b7023e3b 100644 --- a/include/ale.h +++ b/include/ale.h @@ -5,7 +5,7 @@ #include <vector> #include <gsl/gsl_interp.h> - + #include <nlohmann/json.hpp> using json = nlohmann::json; @@ -14,11 +14,30 @@ namespace ale { /// Interpolation enum for defining different methods of interpolation enum interpolation { /// Interpolate using linear interpolation - LINEAR, + LINEAR = 0, /// Interpolate using spline interpolation - SPLINE + SPLINE = 1, }; + + // Temporarily moved interpolation and associated helper functions over from States. Will be + // moved to interpUtils in the future. + + /** The following helper functions are used to calculate the reduced states cache and cubic hermite + to interpolate over it. They were migrated, with minor modifications, from + Isis::NumericalApproximation **/ + + /** Determines the lower index for the interpolation interval. */ + int interpolationIndex(const std::vector<double> ×, double interpTime); + + /** Evaluates a cubic hermite at time, interpTime, between the appropriate two points in x. **/ + double evaluateCubicHermite(const double interpTime, const std::vector<double>& derivs, + const std::vector<double>& x, const std::vector<double>& y); + + /** Evaluate velocities using a Cubic Hermite Spline at a time a, within some interval in x, **/ + double evaluateCubicHermiteFirstDeriv(const double interpTime, const std::vector<double>& deriv, + const std::vector<double>& times, const std::vector<double>& y); + /** *@brief Get the position of the spacecraft at a given time based on a set of coordinates, and their associated times *@param coords A vector of double vectors of coordinates diff --git a/src/States.cpp b/src/States.cpp new file mode 100644 index 0000000000000000000000000000000000000000..9cf78300c9889773f8f85b9036485efcb1af3e8b --- /dev/null +++ b/src/States.cpp @@ -0,0 +1,289 @@ +#include "States.h" + +#include <iostream> +#include <algorithm> +#include <gsl/gsl_interp.h> +#include <gsl/gsl_spline.h> +#include <gsl/gsl_poly.h> +#include <cmath> +#include <float.h> + +namespace ale { + + // Empty constructor + States::States() : m_refFrame(0) { + m_states = {}; + m_ephemTimes = {}; + } + + + States::States(const std::vector<double>& ephemTimes, const std::vector<ale::Vec3d>& positions, + int refFrame) : + m_ephemTimes(ephemTimes), m_refFrame(refFrame) { + // Construct State vector from position and velocity vectors + + if (positions.size() != ephemTimes.size()) { + throw std::invalid_argument("Length of times must match number of positions"); + } + + ale::Vec3d velocities = {0.0, 0.0, 0.0}; + for (ale::Vec3d position : positions) { + m_states.push_back(ale::State(position, velocities)); + } + } + + + States::States(const std::vector<double>& ephemTimes, const std::vector<ale::Vec3d>& positions, + const std::vector<ale::Vec3d>& velocities, int refFrame) : + m_ephemTimes(ephemTimes), m_refFrame(refFrame) { + + if ((positions.size() != ephemTimes.size())||(ephemTimes.size() != velocities.size())) { + throw std::invalid_argument("Length of times must match number of positions and velocities."); + } + + for (int i=0; i < positions.size() ;i++) { + m_states.push_back(ale::State(positions[i], velocities[i])); + } + } + + + States::States(const std::vector<double>& ephemTimes, const std::vector<ale::State>& states, + int refFrame) : + m_ephemTimes(ephemTimes), m_states(states), m_refFrame(refFrame) { + if (states.size() != ephemTimes.size()) { + throw std::invalid_argument("Length of times must match number of states."); + } + } + + // Default Destructor + States::~States() {} + + // Getters + std::vector<ale::State> States::getStates() const { + return m_states; + } + + std::vector<ale::Vec3d> States::getPositions() const { + // extract positions from state vector + std::vector<ale::Vec3d> positions; + + for(ale::State state : m_states) { + positions.push_back(state.position); + } + return positions; + } + + + std::vector<ale::Vec3d> States::getVelocities() const { + // extract velocities from state vector + std::vector<ale::Vec3d> velocities; + + for(ale::State state : m_states) { + velocities.push_back(state.velocity); + } + return velocities; + } + + + std::vector<double> States::getTimes() const { + return m_ephemTimes; + } + + + int States::getReferenceFrame() const { + return m_refFrame; + } + + + bool States::hasVelocity() const { + std::vector<ale::Vec3d> velocities = getVelocities(); + bool allZero = std::all_of(velocities.begin(), velocities.end(), [](ale::Vec3d vec) + { return vec.x==0.0 && vec.y==0.0 && vec.z==0.0; }); + return !allZero; + } + + + ale::State States::getState(double time, ale::interpolation interp) const { + int lowerBound = ale::interpolationIndex(m_ephemTimes, time); + int interpStart; + int interpStop; + + if (m_ephemTimes.size() <= 3) { + interpStart = 0; + interpStop = m_ephemTimes.size() - 1; + } + else if (lowerBound == 0) { + interpStart = lowerBound; + interpStop = lowerBound + 3; + } + else if (lowerBound == m_ephemTimes.size() - 1) { + interpStart = lowerBound - 3; + interpStop = lowerBound; + } + else if (lowerBound == m_ephemTimes.size() - 2) { + interpStart = lowerBound - 2; + interpStop = lowerBound + 1; + } + else { + interpStart = lowerBound - 1; + interpStop = lowerBound + 2; + } + + ale::State state; + std::vector<double> xs, ys, zs, vxs, vys, vzs, interpTimes; + + for (int i = interpStart; i <= interpStop; i++) { + state = m_states[i]; + interpTimes.push_back(m_ephemTimes[i]); + xs.push_back(state.position.x); + ys.push_back(state.position.y); + zs.push_back(state.position.z); + vxs.push_back(state.velocity.x); + vys.push_back(state.velocity.y); + vzs.push_back(state.velocity.z); + } + + ale::Vec3d position, velocity; + + if ( interp == LINEAR || (interp == SPLINE && !hasVelocity())) { + position = {interpolate(xs, interpTimes, time, interp, 0), + interpolate(ys, interpTimes, time, interp, 0), + interpolate(zs, interpTimes, time, interp, 0)}; + + velocity = {interpolate(xs, interpTimes, time, interp, 1), + interpolate(ys, interpTimes, time, interp, 1), + interpolate(zs, interpTimes, time, interp, 1)}; + } + else if (interp == SPLINE && hasVelocity()){ + // Do hermite spline if velocities are available + double baseTime = (interpTimes.front() + interpTimes.back()) / 2; + + std::vector<double> scaledEphemTimes; + for(unsigned int i = 0; i < interpTimes.size(); i++) { + scaledEphemTimes.push_back(interpTimes[i] - baseTime); + } + double sTime = time - baseTime; + position.x = ale::evaluateCubicHermite(sTime, vxs, scaledEphemTimes, xs); + position.y = ale::evaluateCubicHermite(sTime, vys, scaledEphemTimes, ys); + position.z = ale::evaluateCubicHermite(sTime, vzs, scaledEphemTimes, zs); + + velocity.x = ale::evaluateCubicHermiteFirstDeriv(sTime, vxs, scaledEphemTimes, xs); + velocity.y = ale::evaluateCubicHermiteFirstDeriv(sTime, vys, scaledEphemTimes, ys); + velocity.z = ale::evaluateCubicHermiteFirstDeriv(sTime, vzs, scaledEphemTimes, zs); + } + return ale::State(position, velocity); + } + + + ale::Vec3d States::getPosition(double time, ale::interpolation interp) const { + ale::State interpState = getState(time, interp); + return interpState.position; + } + + + ale::Vec3d States::getVelocity(double time, ale::interpolation interp) const { + ale::State interpState = getState(time, interp); + return interpState.velocity; + } + + + double States::getStartTime() { + return m_ephemTimes[0]; + } + + + double States::getStopTime() { + int len = m_ephemTimes.size(); + return m_ephemTimes.back(); + } + + + States States::minimizeCache(double tolerance) { + if (m_states.size() <= 2) { + throw std::invalid_argument("Cache size is 2, cannot minimize."); + } + if (!hasVelocity()) { + throw std::invalid_argument("The cache can only be minimized if velocity is provided."); + } + + // Compute scaled time to use for fitting. + double baseTime = (m_ephemTimes.at(0) + m_ephemTimes.back())/ 2.0; + double timeScale = 1.0; + + // Find current size of m_states + int currentSize = m_ephemTimes.size() - 1; + + // Create 3 starting values for the new size-minimized cache + std::vector <int> inputIndices; + inputIndices.push_back(0); + inputIndices.push_back(currentSize / 2); + inputIndices.push_back(currentSize); + + // find all indices needed to make a hermite table within the appropriate tolerance + std::vector <int> indexList = hermiteIndices(tolerance, inputIndices, baseTime, timeScale); + + // Update m_states and m_ephemTimes to only save the necessary indicies in the index list + std::vector<ale::State> tempStates; + std::vector<double> tempTimes; + + for(int i : indexList) { + tempStates.push_back(m_states[i]); + tempTimes.push_back(m_ephemTimes[i]); + } + return States(tempTimes, tempStates, m_refFrame); + } + + std::vector<int> States::hermiteIndices(double tolerance, std::vector<int> indexList, + double baseTime, double timeScale) { + unsigned int n = indexList.size(); + double sTime; + + std::vector<double> x, y, z, vx, vy, vz, scaledEphemTimes; + for(unsigned int i = 0; i < indexList.size(); i++) { + scaledEphemTimes.push_back((m_ephemTimes[indexList[i]] - baseTime) / timeScale); + x.push_back(m_states[indexList[i]].position.x); + y.push_back(m_states[indexList[i]].position.y); + z.push_back(m_states[indexList[i]].position.z); + vx.push_back(m_states[i].velocity.x); + vy.push_back(m_states[i].velocity.y); + vz.push_back(m_states[i].velocity.z); + } + + // loop through the saved indices from the end + for(unsigned int i = indexList.size() - 1; i > 0; i--) { + double xerror = 0; + double yerror = 0; + double zerror = 0; + + // check every value of the original kernel values within interval + for(int line = indexList[i-1] + 1; line < indexList[i]; line++) { + sTime = (m_ephemTimes[line] - baseTime) / timeScale; + + // find the errors at each value + xerror = fabs(ale::evaluateCubicHermite(sTime, vx, scaledEphemTimes, x) - m_states[line].position.x); + yerror = fabs(ale::evaluateCubicHermite(sTime, vy, scaledEphemTimes, y) - m_states[line].position.y); + zerror = fabs(ale::evaluateCubicHermite(sTime, vz, scaledEphemTimes, z) - m_states[line].position.z); + + if(xerror > tolerance || yerror > tolerance || zerror > tolerance) { + // if any error is greater than tolerance, no need to continue looking, break + break; + } + } + + if(xerror < tolerance && yerror < tolerance && zerror < tolerance) { + // if errors are less than tolerance after looping interval, no new point is necessary + continue; + } + else { + // if any error is greater than tolerance, add midpoint of interval to indexList vector + indexList.push_back((indexList[i] + indexList[i-1]) / 2); + } + } + + if(indexList.size() > n) { + std::sort(indexList.begin(), indexList.end()); + indexList = hermiteIndices(tolerance, indexList, baseTime, timeScale); + } + return indexList; + } +} diff --git a/src/ale.cpp b/src/ale.cpp index 5646c3ee3d9e61a37160614e6b87ae46986ba187..3bcdc31c2f9989562a6509f10b4580233124bab4 100644 --- a/src/ale.cpp +++ b/src/ale.cpp @@ -21,6 +21,93 @@ using namespace std; namespace ale { + + // Temporarily moved over from States.cpp. Will be moved into interpUtils in the future. + + /** The following helper functions are used to calculate the reduced states cache and cubic hermite + to interpolate over it. They were migrated, with minor modifications, from + Isis::NumericalApproximation **/ + + /** Determines the lower index for the interpolation interval. */ + int interpolationIndex(const std::vector<double> ×, double interpTime) { + if (times.size() < 2){ + throw std::invalid_argument("There must be at least two times."); + } + auto nextTimeIt = std::upper_bound(times.begin(), times.end(), interpTime); + if (nextTimeIt == times.end()) { + --nextTimeIt; + } + if (nextTimeIt != times.begin()) { + --nextTimeIt; + } + return std::distance(times.begin(), nextTimeIt); + } + + + /** Evaluates a cubic hermite at time, interpTime, between the appropriate two points in x. **/ + double evaluateCubicHermite(const double interpTime, const std::vector<double>& derivs, + const std::vector<double>& x, const std::vector<double>& y) { + if( (derivs.size() != x.size()) || (derivs.size() != y.size()) ) { + throw std::invalid_argument("EvaluateCubicHermite - The size of the first derivative vector does not match the number of (x,y) data points."); + } + + // Find the interval in which "a" exists + int lowerIndex = ale::interpolationIndex(x, interpTime); + + double x0, x1, y0, y1, m0, m1; + // interpTime is contained within the interval (x0,x1) + x0 = x[lowerIndex]; + x1 = x[lowerIndex+1]; + // the corresponding known y-values for x0 and x1 + y0 = y[lowerIndex]; + y1 = y[lowerIndex+1]; + // the corresponding known tangents (slopes) at (x0,y0) and (x1,y1) + m0 = derivs[lowerIndex]; + m1 = derivs[lowerIndex+1]; + + double h, t; + h = x1 - x0; + t = (interpTime - x0) / h; + return (2 * t * t * t - 3 * t * t + 1) * y0 + (t * t * t - 2 * t * t + t) * h * m0 + (-2 * t * t * t + 3 * t * t) * y1 + (t * t * t - t * t) * h * m1; + } + + /** Evaluate velocities using a Cubic Hermite Spline at a time a, within some interval in x, **/ + double evaluateCubicHermiteFirstDeriv(const double interpTime, const std::vector<double>& deriv, + const std::vector<double>& times, const std::vector<double>& y) { + if(deriv.size() != times.size()) { + throw std::invalid_argument("EvaluateCubicHermiteFirstDeriv - The size of the first derivative vector does not match the number of (x,y) data points."); + } + + // find the interval in which "interpTime" exists + int lowerIndex = ale::interpolationIndex(times, interpTime); + + double x0, x1, y0, y1, m0, m1; + + // interpTime is contained within the interval (x0,x1) + x0 = times[lowerIndex]; + x1 = times[lowerIndex+1]; + + // the corresponding known y-values for x0 and x1 + y0 = y[lowerIndex]; + y1 = y[lowerIndex+1]; + + // the corresponding known tangents (slopes) at (x0,y0) and (x1,y1) + m0 = deriv[lowerIndex]; + m1 = deriv[lowerIndex+1]; + + double h, t; + h = x1 - x0; + t = (interpTime - x0) / h; + if(h != 0.0) { + return ((6 * t * t - 6 * t) * y0 + (3 * t * t - 4 * t + 1) * h * m0 + (-6 * t * t + 6 * t) * y1 + (3 * t * t - 2 * t) * h * m1) / h; + + } + else { + throw std::invalid_argument("Error in evaluating cubic hermite velocities, values at" + "lower and upper indicies are exactly equal."); + } + } + // Position Data Functions vector<double> getPosition(vector<vector<double>> coords, vector<double> times, double time, interpolation interp) { diff --git a/tests/ctests/AleTest.cpp b/tests/ctests/AleTest.cpp index 4e690a3392e8edd8eeaf5ac13ce7bf909663ae3b..9ee7d1345f9d3e65ad09891fc713c1ec811a9794 100644 --- a/tests/ctests/AleTest.cpp +++ b/tests/ctests/AleTest.cpp @@ -347,8 +347,14 @@ TEST(VelocityInterpTest, InvalidTime) { EXPECT_THROW(ale::getVelocity(coords, times, 3, ale::LINEAR), invalid_argument); } -TEST(Interpolation, Derivative2) -{ +TEST(Interpolation, Derivative1) { + vector<double> points = {0, 2, 4}; + vector<double> times = {0, 1, 2}; + + EXPECT_DOUBLE_EQ(ale::interpolate(points, times, 1, ale::LINEAR, 1), 2); +} + +TEST(Interpolation, Derivative2) { //only checks that case 2 of the switch block is hit vector<double> points = {0, 0, 0}; vector<double> times = {0, 1, 2}; @@ -362,3 +368,53 @@ TEST(Interpolation, InvalidDerivative) { EXPECT_THROW(ale::interpolate(points, times, 1, ale::LINEAR, 3), invalid_argument); } + +TEST(InterpIndex, InvalidTimes) { + std::vector<double> times = {}; + + EXPECT_THROW(ale::interpolationIndex(times, 0), invalid_argument); +} + +TEST(EvaluateCubicHermite, SimplyPolynomial) { + // Cubic function is y = x^3 - 2x^2 + 1 + // derivative is dy/dx = 3x^2 - 4x + std::vector<double> derivs = {7.0, -1.0}; + std::vector<double> x = {-1.0, 1.0}; + std::vector<double> y = {-2.0, 0.0}; + + EXPECT_DOUBLE_EQ(ale::evaluateCubicHermite(0.0, derivs, x, y), 1.0); +} + +TEST(EvaluateCubicHermite, InvalidDervisXY) { + std::vector<double> derivs = {}; + std::vector<double> x = {1.0}; + std::vector<double> y = {1.0}; + + EXPECT_THROW(ale::evaluateCubicHermite(0.0, derivs, x, y), invalid_argument); +} + +TEST(EvaluateCubicHermiteFirstDeriv, SimplyPolynomial) { + // Cubic function is y = x^3 - 2x^2 + 1 + // derivative is dy/dx = 3x^2 - 4x + std::vector<double> derivs = {7.0, -1.0}; + std::vector<double> x = {-1.0, 1.0}; + std::vector<double> y = {-2.0, 0.0}; + + EXPECT_DOUBLE_EQ(ale::evaluateCubicHermiteFirstDeriv(0.5, derivs, x, y), -1.25); +} + +TEST(EvaluateCubicHermiteFirstDeriv, InvalidDervisTimes) { + std::vector<double> derivs = {}; + std::vector<double> times = {1.0}; + std::vector<double> y = {1.0}; + + EXPECT_THROW(ale::evaluateCubicHermiteFirstDeriv(0.0, derivs, times, y), invalid_argument); +} + +TEST(EvaluateCubicHermiteFirstDeriv, InvalidVelocities) { + std::vector<double> derivs = {5.0, 6.0}; + std::vector<double> times = {1.0, 1.0}; + std::vector<double> y = {1.0}; + + EXPECT_THROW(ale::evaluateCubicHermiteFirstDeriv(0.0, derivs, times, y), invalid_argument); +} diff --git a/tests/ctests/StatesTest.cpp b/tests/ctests/StatesTest.cpp new file mode 100644 index 0000000000000000000000000000000000000000..b96e6650b910eb9a5376c6f0cd2deab4d4cc146f --- /dev/null +++ b/tests/ctests/StatesTest.cpp @@ -0,0 +1,467 @@ +#include "gtest/gtest.h" + +#include "States.h" + +#include <cmath> +#include <exception> + +using namespace std; +using namespace ale; + +TEST(StatesTest, DefaultConstructor) { + States defaultState; + vector<ale::State> states = defaultState.getStates(); + EXPECT_EQ(states.size(), 0); +} + +TEST(StatesTest, ConstructorPositionNoVelocity) { + std::vector<double> ephemTimes = {0.0, 1.0, 2.0, 3.0}; + std::vector<ale::Vec3d> positions = { + ale::Vec3d(4.0, 1.0, 4.0), + ale::Vec3d (5.0, 2.0, 3.0), + ale::Vec3d (6.0, 3.0, 2.0), + ale::Vec3d (7.0, 4.0, 1.0)}; + + States noVelocityState(ephemTimes, positions); + vector<ale::State> states = noVelocityState.getStates(); + EXPECT_EQ(states.size(), 4); + EXPECT_NEAR(states[0].position.x, 4.0, 1e-10); + EXPECT_NEAR(states[0].position.y, 1.0, 1e-10); + EXPECT_NEAR(states[0].position.z, 4.0, 1e-10); + EXPECT_NEAR(states[0].velocity.x, 0.0, 1e-10); + EXPECT_NEAR(states[0].velocity.y, 0.0, 1e-10); + EXPECT_NEAR(states[0].velocity.z, 0.0, 1e-10); + EXPECT_NEAR(states[3].position.x, 7.0, 1e-10); + EXPECT_NEAR(states[3].position.y, 4.0, 1e-10); + EXPECT_NEAR(states[3].position.z, 1.0, 1e-10); + EXPECT_NEAR(states[3].velocity.x, 0.0, 1e-10); + EXPECT_NEAR(states[3].velocity.y, 0.0, 1e-10); + EXPECT_NEAR(states[3].velocity.z, 0.0, 1e-10); +} + +TEST(StatesTest, ConstructorPositionAndVelocity) { + std::vector<double> ephemTimes = {0.0, 1.0, 2.0, 3.0}; + std::vector<ale::Vec3d> positions = { + ale::Vec3d(4.0, 1.0, 4.0), + ale::Vec3d (5.0, 2.0, 3.0), + ale::Vec3d (6.0, 3.0, 2.0), + ale::Vec3d (7.0, 4.0, 1.0)}; + + std::vector<ale::Vec3d> velocities = { + ale::Vec3d(-4.0, -1.0, -4.0), + ale::Vec3d(-5.0, -2.0, -3.0), + ale::Vec3d(-6.0, -3.0, -2.0), + ale::Vec3d(-7.0, -4.0, -1.0)}; + + States positionVelocityState(ephemTimes, positions, velocities); + vector<ale::State> states = positionVelocityState.getStates(); + EXPECT_EQ(states.size(), 4); + EXPECT_NEAR(states[0].position.x, 4.0, 1e-10); + EXPECT_NEAR(states[0].position.y, 1.0, 1e-10); + EXPECT_NEAR(states[0].position.z, 4.0, 1e-10); + EXPECT_NEAR(states[0].velocity.x, -4.0, 1e-10); + EXPECT_NEAR(states[0].velocity.y, -1.0, 1e-10); + EXPECT_NEAR(states[0].velocity.z, -4.0, 1e-10); + EXPECT_NEAR(states[3].position.x, 7.0, 1e-10); + EXPECT_NEAR(states[3].position.y, 4.0, 1e-10); + EXPECT_NEAR(states[3].position.z, 1.0, 1e-10); + EXPECT_NEAR(states[3].velocity.x, -7.0, 1e-10); + EXPECT_NEAR(states[3].velocity.y, -4.0, 1e-10); + EXPECT_NEAR(states[3].velocity.z, -1.0, 1e-10); +} + +TEST(StatesTest, ConstructorStates) { + std::vector<double> ephemTimes = {0.0, 1.0, 2.0, 3.0}; + + std::vector<ale::Vec3d> positions = { + ale::Vec3d(4.0, 1.0, 4.0), + ale::Vec3d (5.0, 2.0, 3.0), + ale::Vec3d (6.0, 3.0, 2.0), + ale::Vec3d (7.0, 4.0, 1.0)}; + + std::vector<ale::Vec3d> velocities = { + ale::Vec3d(-4.0, -1.0, -4.0), + ale::Vec3d(-5.0, -2.0, -3.0), + ale::Vec3d(-6.0, -3.0, -2.0), + ale::Vec3d(-7.0, -4.0, -1.0)}; + + std::vector<ale::State> stateVector; + for (int i=0; i < positions.size(); i++) { + stateVector.push_back(ale::State(positions[i], velocities[i])); + } + + States statesState(ephemTimes, stateVector, 1); + vector<ale::State> states = statesState.getStates(); + EXPECT_EQ(states.size(), 4); + EXPECT_NEAR(states[0].position.x, 4.0, 1e-10); + EXPECT_NEAR(states[0].position.y, 1.0, 1e-10); + EXPECT_NEAR(states[0].position.z, 4.0, 1e-10); + EXPECT_NEAR(states[0].velocity.x, -4.0, 1e-10); + EXPECT_NEAR(states[0].velocity.y, -1.0, 1e-10); + EXPECT_NEAR(states[0].velocity.z, -4.0, 1e-10); + EXPECT_NEAR(states[3].position.x, 7.0, 1e-10); + EXPECT_NEAR(states[3].position.y, 4.0, 1e-10); + EXPECT_NEAR(states[3].position.z, 1.0, 1e-10); + EXPECT_NEAR(states[3].velocity.x, -7.0, 1e-10); + EXPECT_NEAR(states[3].velocity.y, -4.0, 1e-10); + EXPECT_NEAR(states[3].velocity.z, -1.0, 1e-10); +} + +class TestState : public ::testing::Test { +protected: + States *states; + States *statesNoVelocity; + + // fixtures.h had a lot of overrides before the steup and tearDqwn functions + void SetUp() override { + states = NULL; + statesNoVelocity = NULL; + + //define test data values + std::vector<double> ephemTimes = {0.0, 1.0, 2.0, 3.0}; + + /** + X func: x = time + 4 + Y func: y = 4 - (time - 2)^2 + Z func: z = (time / 2.5)^3 + + All postiion values were obtained using the above functions and all velocity + values were obtained using the derivative of the above functions at the defined + ephemeris times in ephemTimes. + **/ + std::vector<ale::Vec3d> positions = { + ale::Vec3d(4, 0, 0), + ale::Vec3d(5, 3, 0.064), + ale::Vec3d(6, 4, 0.512), + ale::Vec3d(7, 3, 1.728)}; + + std::vector<ale::Vec3d> velocities = { + ale::Vec3d(1, 4, 0), + ale::Vec3d(1, 2, 0.48), + ale::Vec3d(1, 0, 1.92), + ale::Vec3d(1, -2, 4.32)}; + + std::vector<ale::State> stateVector; + for (int i=0; i < positions.size(); i++) { + stateVector.push_back(ale::State(positions[i], velocities[i])); + } + + states = new ale::States(ephemTimes, stateVector); + statesNoVelocity = new ale::States(ephemTimes, positions); + } + + void TearDown() override { + if (states) { + delete states; + states = NULL; + } + if (statesNoVelocity) { + delete statesNoVelocity; + statesNoVelocity = NULL; + } + } +}; + + +// Tests GSL spline and linear interp - force to use GSL spline by omitting velocity +TEST_F(TestState, getPosition) { + double time = 1.5; + + ale::Vec3d linear_position = states->getPosition(time, LINEAR); + ale::Vec3d spline_position = states->getPosition(time, SPLINE); + ale::Vec3d linear_no_vel_position = statesNoVelocity->getPosition(time, LINEAR); + ale::Vec3d spline_no_vel_position = statesNoVelocity->getPosition(time, SPLINE); + + double linear_x = 5.5; + double linear_y = 3.5; + double linear_z = 0.288; + + EXPECT_NEAR(linear_position.x, linear_x, 1e-10); + EXPECT_NEAR(linear_position.y, linear_y, 1e-10); + EXPECT_NEAR(linear_position.z, linear_z, 1e-10); + EXPECT_NEAR(spline_position.x, 5.5, 1e-10); + EXPECT_NEAR(spline_position.y, 3.75, 1e-10); + EXPECT_NEAR(spline_position.z, 0.108, 1e-10); + EXPECT_NEAR(linear_no_vel_position.x, linear_x, 1e-10); + EXPECT_NEAR(linear_no_vel_position.y, linear_y, 1e-10); + EXPECT_NEAR(linear_no_vel_position.z, linear_z, 1e-10); + EXPECT_NEAR(spline_no_vel_position.x, 5.5, 1e-10); + EXPECT_NEAR(spline_no_vel_position.y, 3.8, 1e-10); + EXPECT_NEAR(spline_no_vel_position.z, 0.2016, 1e-10); +} + + +TEST_F(TestState, getVelocity) { + double time = 1.5; + ale::Vec3d linear_velocity = states->getVelocity(time, LINEAR); + ale::Vec3d spline_velocity = states->getVelocity(time, SPLINE); + ale::Vec3d linear_no_vel_velocity = statesNoVelocity->getVelocity(time, LINEAR); + ale::Vec3d spline_no_vel_velocity = statesNoVelocity->getVelocity(time, SPLINE); + + EXPECT_NEAR(linear_velocity.x, 1, 1e-10); + EXPECT_NEAR(linear_velocity.y, 1, 1e-10); + EXPECT_NEAR(linear_velocity.z, 0.448, 1e-10); + EXPECT_NEAR(spline_velocity.x, 1, 1e-10); + EXPECT_NEAR(spline_velocity.y, 1, 1e-10); + EXPECT_NEAR(spline_velocity.z, 0.072, 1e-10); + EXPECT_NEAR(linear_no_vel_velocity.x, 1, 1e-10); + EXPECT_NEAR(linear_no_vel_velocity.y, 1, 1e-10); + EXPECT_NEAR(linear_no_vel_velocity.z, 0.448, 1e-10); + EXPECT_NEAR(spline_no_vel_velocity.x, 1, 1e-10); + EXPECT_NEAR(spline_no_vel_velocity.y, 1, 1e-10); + EXPECT_NEAR(spline_no_vel_velocity.z, 0.416, 1e-10); +} + +// getState() and interpolateState() are tested when testing getPosition and +// getVelocity, because they are derived from those methods + +TEST_F(TestState, getStartTime) { + double time = states->getStartTime(); + EXPECT_NEAR(time, 0, 1e-10); +} + +TEST_F(TestState, getStopTime) { + double time = states->getStopTime(); + EXPECT_NEAR(time, 3.0, 1e-10); +} + +TEST_F(TestState, getTimes) { + std::vector<double> time = states->getTimes(); + EXPECT_NEAR(time[0], 0.0, 1e-10); + EXPECT_NEAR(time[3], 3.0, 1e-10); +} + +TEST_F(TestState, getReferenceFrame) { + int frame = states->getReferenceFrame(); + EXPECT_EQ(frame, 1); +} + +TEST_F(TestState, getPositions) { + std::vector<ale::Vec3d> positions = states->getPositions(); + EXPECT_EQ(positions.size(), 4); + EXPECT_NEAR(positions[0].x, 4.0, 1e-10); + EXPECT_NEAR(positions[0].y, 0.0, 1e-10); + EXPECT_NEAR(positions[0].z, 0.0, 1e-10); + EXPECT_NEAR(positions[3].x, 7.0, 1e-10); + EXPECT_NEAR(positions[3].y, 3.0, 1e-10); + EXPECT_NEAR(positions[3].z, 1.728, 1e-10); +} + +// minimizeCache +// one test where the cache can be minimized +// one test where the cache cannot be minimized +TEST(StatesTest, minimizeCache_true) { + //create information that can be reduced + std::vector<double> ephemTimes = {0.0, 1.0, 2.0, 3.0, 4.0}; + std::vector<ale::Vec3d> positions = { + ale::Vec3d(5.0, 3.0, 0.0), + ale::Vec3d(5.6, 5.1875, 0.3), + ale::Vec3d(7.0, 6.0, 1.0), + ale::Vec3d(10.025, 5.5625, 2.5125), + ale::Vec3d(15.0, 4.0, 5.0)}; + + States withVelocityState(ephemTimes, positions, positions, 1); + + vector<ale::State> states = withVelocityState.getStates(); + EXPECT_EQ(states.size(), 5); + States minimizedStates = withVelocityState.minimizeCache(1); + vector<ale::State> states_min = minimizedStates.getStates(); + EXPECT_EQ(states_min.size(), 4); +} + +// Creates a setup in which the cache cannot be minimized +TEST(StatesTest, minimizeCache_false) { + //create information that cannot be reduced + std::vector<double> ephemTimes = {0.0, 1.0, 2.0, 3.0, 4.0}; + std::vector<ale::Vec3d> positions = { + ale::Vec3d(5.0, 3.0, 0.0), + ale::Vec3d(6.5, 4.0, 1.0), + ale::Vec3d(7.0, 6.0, 1.0), + ale::Vec3d(8.7, 3.0, 3.0), + ale::Vec3d(15.0, 4.0, 5.0)}; + + //create states with no velocity information + States noVelocityState(ephemTimes, positions); + vector<ale::State> states = noVelocityState.getStates(); + EXPECT_EQ(states.size(), 5); + + ASSERT_THROW(noVelocityState.minimizeCache(), std::invalid_argument); + + std::vector<ale::Vec3d> velocities = { + ale::Vec3d(5.0, 3.0, 0.0), + ale::Vec3d(6.5, 4.0, 1.0), + ale::Vec3d(7.0, 6.0, 1.0), + ale::Vec3d(8.7, 3.0, 3.0), + ale::Vec3d(15.0, 4.0, 5.0)}; + + States withVelocityState(ephemTimes, positions, velocities, 1); + vector<ale::State> states_min = withVelocityState.getStates(); + EXPECT_EQ(states_min.size(), 5); +} + + +// This test checks to see if the minimized cache looks identical to ISIS's minimized cache for +// a Dawn IR image VIR_IR_1A_1_362681634_1 (located in dawnvir2isis's IR app test.) +// Values were obtained by adding strategic couts to SpicePosition.cpp, running spiceinit, and +// pasting the results in here. +TEST(StatesTest, minimizeCache_matchesISIS) { + // 362681869.7384 + std::vector<double> ephemTimes = {362681633.613412,362681649.355078,362681665.096744,362681680.83841,362681696.580076,362681712.321741,362681728.063407,362681743.805073,362681759.546739,362681775.288405,362681791.030071,362681806.771737,362681822.513403,362681838.255068,362681853.996734,362681869.7384,362681885.480066,362681901.221732,362681916.963398,362681932.705064,362681948.44673,362681964.188395,362681979.930061,362681995.671727,362682011.413393,362682027.155059,362682042.896725,362682058.638391,362682074.380057,362682090.121722,362682105.863388,362682121.605054,362682137.34672,362682153.088386,362682168.830052,362682184.571718,362682200.313384,362682216.055049,362682231.796715,362682247.538381,362682263.280047,362682279.021713,362682294.763379,362682310.505045,362682326.246711,362682341.988376,362682357.730042,362682373.471708,362682389.213374,362682404.95504,362682420.696706,362682436.438372,362682452.180037,362682467.921703,362682483.663369,362682499.405035,362682515.146701,362682530.888367,362682546.630033,362682562.371699,362682578.113365}; + + std::vector<ale::Vec3d> positions = { +ale::Vec3d(-17382.7468835417, 96577.8576989543, 16364.0677257831), +ale::Vec3d(-17382.9888346502, 96576.3340464433, 16363.6197071606), +ale::Vec3d(-17383.2307847157, 96574.8103941707, 16363.1716886571), +ale::Vec3d(-17383.472735884, 96573.2867410342, 16362.7236700492), +ale::Vec3d(-17383.7146860391, 96571.7630881361, 16362.2756515455), +ale::Vec3d(-17383.9566372671, 96570.2394343143, 16361.8276329373), +ale::Vec3d(-17384.1985874818, 96568.715780731, 16361.3796144184), +ale::Vec3d(-17384.4405387696, 96567.1921262539, 16360.9315958249), +ale::Vec3d(-17384.6824890441, 96565.668472045, 16360.4835773057), +ale::Vec3d(-17384.9244393485, 96564.1448175084, 16360.0355588162), +ale::Vec3d(-17385.1663907556, 96562.6211619886, 16359.5875402223), +ale::Vec3d(-17385.4083410899, 96561.0975068265, 16359.1395216879), +ale::Vec3d(-17385.650292527, 96559.5738506513, 16358.6915030937), +ale::Vec3d(-17385.8922429508, 96558.050194804, 16358.2434845738), +ale::Vec3d(-17386.1341944775, 96556.526538033, 16357.7954659943), +ale::Vec3d(-17386.3761449312, 96555.0028814707, 16357.3474474742),// position to test +ale::Vec3d(-17386.6180964878, 96553.4792240146, 16356.8994288944), +ale::Vec3d(-17386.8600470013, 96551.9555667671, 16356.4514103591), +ale::Vec3d(-17387.1019987069, 96550.4319087152, 16356.003391794), +ale::Vec3d(-17387.3439492802, 96548.9082508421, 16355.5553732733), +ale::Vec3d(-17387.5859009562, 96547.3845920454, 16355.1073546781), +ale::Vec3d(-17387.8278515892, 96545.8609335467, 16354.659336187), +ale::Vec3d(-17388.0698022818, 96544.3372746906, 16354.2113176809), +ale::Vec3d(-17388.3117540474, 96542.8136149108, 16353.7632991002), +ale::Vec3d(-17388.5537048296, 96541.2899553993, 16353.3152806087), +ale::Vec3d(-17388.7956566549, 96539.766294994, 16352.8672620277), +ale::Vec3d(-17389.0376074372, 96538.2426347973, 16352.4192435359), +ale::Vec3d(-17389.2795593521, 96536.7189737366, 16351.9712249249), +ale::Vec3d(-17389.5215102239, 96535.1953128845, 16351.5232064478), +ale::Vec3d(-17389.7634621985, 96533.6716511088, 16351.0751878513), +ale::Vec3d(-17390.0054131003, 96532.1479896013, 16350.627169374), +ale::Vec3d(-17390.2473651048, 96530.6243272596, 16350.1791507773), +ale::Vec3d(-17390.4893161259, 96529.100665067, 16349.7311322847), +ale::Vec3d(-17390.7312681604, 96527.5770019805, 16349.2831137176), +ale::Vec3d(-17390.9732192412, 96526.0533391921, 16348.8350952099), +ale::Vec3d(-17391.2151713652, 96524.5296754503, 16348.3870766424), +ale::Vec3d(-17391.4571224462, 96523.0060120065, 16347.9390581494), +ale::Vec3d(-17391.6990736166, 96521.4823481456, 16347.4910396562), +ale::Vec3d(-17391.9410258302, 96519.9586833909, 16347.0430210884), +ale::Vec3d(-17392.1829770008, 96518.435018964, 16346.59500258), +ale::Vec3d(-17392.4249293039, 96516.9113535837, 16346.1469840417), +ale::Vec3d(-17392.666880564, 96515.3876884717, 16345.6989655182), +ale::Vec3d(-17392.9088328673, 96513.8640223764, 16345.2509469647), +ale::Vec3d(-17393.1507842169, 96512.340356609, 16344.8029284856), +ale::Vec3d(-17393.3927366097, 96510.8166898881, 16344.354909917), +ale::Vec3d(-17393.6346879893, 96509.2930234355, 16343.9068914228), +ale::Vec3d(-17393.8766404716, 96507.7693560891, 16343.4588728688), +ale::Vec3d(-17394.1185918513, 96506.2456889215, 16343.0108543743), +ale::Vec3d(-17394.3605443934, 96504.7220208899, 16342.5628358051), +ale::Vec3d(-17394.6024958924, 96503.1983530967, 16342.1148173252), +ale::Vec3d(-17394.8444484644, 96501.6746844097, 16341.6667987558), +ale::Vec3d(-17395.0863999934, 96500.1510159612, 16341.2187802906), +ale::Vec3d(-17395.3283526549, 96498.6273466784, 16340.7707617209), +ale::Vec3d(-17395.5703042734, 96497.1036775447, 16340.3227432554), +ale::Vec3d(-17395.8122558622, 96495.5800080833, 16339.8747247748), +ale::Vec3d(-17396.0542086133, 96494.0563377578, 16339.4267061899), +ale::Vec3d(-17396.2961603214, 96492.5326676411, 16338.9786877239), +ale::Vec3d(-17396.5381131323, 96491.0089966602, 16338.5306691834), +ale::Vec3d(-17396.7800648703, 96489.4853258881, 16338.0826506874), +ale::Vec3d(-17397.0220177409, 96487.9616541922, 16337.6346321317), +ale::Vec3d(-17397.2639695387, 96486.4379827647, 16337.1866136653) + }; + +std::vector<ale::Vec3d> velocities = { +ale::Vec3d(-0.0153700718765843, -0.0967910279657348, -0.0284606828741588), +ale::Vec3d(-0.0153700737931066, -0.096791048945897, -0.0284606828190555), +ale::Vec3d(-0.0153700757092848, -0.0967910699268904, -0.0284606827640637), +ale::Vec3d(-0.0153700776250131, -0.0967910909088181, -0.0284606827092379), +ale::Vec3d(-0.0153700795403742, -0.0967911118915676, -0.0284606826545261), +ale::Vec3d(-0.0153700814553067, -0.0967911328752586, -0.0284606825999784), +ale::Vec3d(-0.015370083369856, -0.0967911538597701, -0.0284606825455435), +ale::Vec3d(-0.0153700852839879, -0.0967911748452176, -0.0284606824912741), +ale::Vec3d(-0.0153700871977456, -0.0967911958314943, -0.0284606824371165), +ale::Vec3d(-0.0153700891110573, -0.0967912168186004, -0.0284606823830493), +ale::Vec3d(-0.0153700910239836, -0.0967912378066924, -0.0284606823291967), +ale::Vec3d(-0.0153700929365416, -0.0967912587956045, -0.0284606822754565), +ale::Vec3d(-0.0153700948486603, -0.0967912797854579, -0.0284606822218836), +ale::Vec3d(-0.0153700967604118, -0.0967913007761374, -0.0284606821684221), +ale::Vec3d(-0.0153700986717228, -0.0967913217677527, -0.0284606821151287),// velocity to test +ale::Vec3d(-0.0153701005826786, -0.0967913427601938, -0.0284606820619469), +ale::Vec3d(-0.0153701024931907, -0.0967913637535689, -0.0284606820089321), +ale::Vec3d(-0.0153701044033458, -0.0967913847477737, -0.0284606819560303), +ale::Vec3d(-0.0153701063130629, -0.0967914057429141, -0.0284606819032941), +ale::Vec3d(-0.0153701082224151, -0.0967914267388804, -0.0284606818506722), +ale::Vec3d(-0.0153701101313204, -0.0967914477357808, -0.0284606817982155), +ale::Vec3d(-0.0153701120398714, -0.0967914687335107, -0.0284606817458717), +ale::Vec3d(-0.0153701139480175, -0.096791489732125, -0.028460681693667), +ale::Vec3d(-0.0153701158557277, -0.0967915107316699, -0.0284606816416301), +ale::Vec3d(-0.0153701177630685, -0.0967915317320424, -0.0284606815897053), +ale::Vec3d(-0.0153701196699842, -0.0967915527333597, -0.0284606815379469), +ale::Vec3d(-0.0153701215765303, -0.0967915737354921, -0.028460681486301), +ale::Vec3d(-0.0153701234826341, -0.0967915947385677, -0.0284606814348225), +ale::Vec3d(-0.0153701253883785, -0.0967916157424725, -0.0284606813834553), +ale::Vec3d(-0.0153701272936794, -0.0967916367473043, -0.0284606813322556), +ale::Vec3d(-0.015370129198612, -0.0967916577529709, -0.0284606812811673), +ale::Vec3d(-0.0153701311031148, -0.0967916787595734, -0.0284606812302468), +ale::Vec3d(-0.0153701330072486, -0.0967916997669981, -0.0284606811794379), +ale::Vec3d(-0.0153701349109495, -0.0967917207753623, -0.0284606811287951), +ale::Vec3d(-0.0153701368142852, -0.0967917417845576, -0.0284606810782657), +ale::Vec3d(-0.0153701387171744, -0.0967917627946783, -0.0284606810279025), +ale::Vec3d(-0.0153701406197154, -0.0967917838056372, -0.0284606809776528), +ale::Vec3d(-0.0153701425218418, -0.0967918048174734, -0.0284606809275409), +ale::Vec3d(-0.0153701444235427, -0.0967918258302474, -0.0284606808775962), +ale::Vec3d(-0.0153701463248632, -0.0967918468438453, -0.0284606808277636), +ale::Vec3d(-0.0153701482257517, -0.0967918678583828, -0.028460680778098), +ale::Vec3d(-0.0153701501262747, -0.0967918888737514, -0.0284606807285448), +ale::Vec3d(-0.0153701520263642, -0.0967919098900505, -0.0284606806791581), +ale::Vec3d(-0.0153701539260803, -0.0967919309071788, -0.0284606806298831), +ale::Vec3d(-0.0153701558253628, -0.0967919519252454, -0.0284606805807748), +ale::Vec3d(-0.0153701577242883, -0.0967919729441372, -0.028460680531779), +ale::Vec3d(-0.0153701596227715, -0.0967919939639688, -0.0284606804829498), +ale::Vec3d(-0.0153701615208887, -0.0967920149846261, -0.0284606804342323), +ale::Vec3d(-0.015370163418563, -0.0967920360062177, -0.028460680385682), +ale::Vec3d(-0.0153701653158791, -0.0967920570286386, -0.0284606803372438), +ale::Vec3d(-0.0153701672127516, -0.0967920780520043, -0.0284606802889809), +ale::Vec3d(-0.0153701691092572, -0.0967920990761938, -0.0284606802408221), +ale::Vec3d(-0.0153701710053287, -0.0967921201013125, -0.0284606801928288), +ale::Vec3d(-0.0153701729010412, -0.096792141127264, -0.0284606801449487), +ale::Vec3d(-0.0153701747963447, -0.0967921621540997, -0.0284606800972071), +ale::Vec3d(-0.0153701766912129, -0.096792183181868, -0.0284606800496334), +ale::Vec3d(-0.0153701785857075, -0.0967922042104672, -0.0284606800021717), +ale::Vec3d(-0.0153701804797747, -0.0967922252399974, -0.0284606799548765), +ale::Vec3d(-0.0153701823734752, -0.0967922462703656, -0.0284606799076943), +ale::Vec3d(-0.0153701842667371, -0.0967922673016664, -0.0284606798606782), +ale::Vec3d(-0.015370186159624, -0.0967922883337961, -0.0284606798137744) +}; + +States testState(ephemTimes, positions, velocities); +States minimizedState = testState.minimizeCache(); + +// Test the ability to recover the original coordinates and velocity within the tolerance +// from the reduced cache. (Aribtrarily selected the 15th index.) +ale::State result = minimizedState.getState(362681869.7384); +EXPECT_NEAR(result.position.x, -17386.3761449312, 1e-5); +EXPECT_NEAR(result.position.y, 96555.0028814707, 1e-4); +EXPECT_NEAR(result.position.z, 16357.3474474742, 1e-5); +EXPECT_NEAR(result.velocity.x, -0.0153700986717228, 1e-6); +EXPECT_NEAR(result.velocity.y, -0.0967913217677527, 1e-6); +EXPECT_NEAR(result.velocity.z, -0.0284606821151287, 1e-6); + +// Get all the states to check that they match exactly with ISIS's reduced cache: +std::vector<ale::State> results = minimizedState.getStates(); +EXPECT_NEAR(results[0].position.x, -17382.7468835417, 1e-8); +EXPECT_NEAR(results[0].position.y, 96577.8576989543, 1e-8); +EXPECT_NEAR(results[0].position.z, 16364.0677257831, 1e-8); +EXPECT_NEAR(results[1].position.x, -17390.0054131003, 1e-8); +EXPECT_NEAR(results[1].position.y, 96532.1479896013, 1e-8); +EXPECT_NEAR(results[1].position.z, 16350.627169374, 1e-8); +EXPECT_NEAR(results[2].position.x, -17397.2639695387, 1e-8); +EXPECT_NEAR(results[2].position.y, 96486.4379827647, 1e-8); +EXPECT_NEAR(results[2].position.z, 16337.1866136653, 1e-8); + +// The cache size should be reduced from 61 to 3 to match ISIS +EXPECT_EQ(results.size(), 3); + +}