simulation.h#
Functions
-
std::vector<std::vector<real>> reconstruct_stm(const std::vector<real> &stm)#
Reconstruct the STM matrix from the flattened STM vector.
- Parameters:
stm – [in] The flattened STM vector.
- Returns:
The reconstructed STM matrix.
-
void get_baseBodyFrame(const int &spiceId, const real &tMjdTDB, std::string &baseBodyFrame)#
Get the name of body-fixed frame for a given SPICE ID.
- Parameters:
spiceId – [in] SPICE ID of the body.
tMjdTDB – [in] Time of observation in Modified Julian Date (TDB).
baseBodyFrame – [out] Name of the body-fixed frame.
-
void get_observer_state(const real &tObsMjd, const std::vector<real> &observerInfo, PropSimulation *propSim, const bool tObsInUTC, std::vector<real> &observerState)#
Get the observer state for a given time.
- Parameters:
tObsMjd – [in] Observation time in Modified Julian Date.
observerInfo – [in] Observer information (base body SPICE ID, lon [rad], lat [rad], alt [m])
propSim – [in] PropSimulation object for querying Ephemeris of base body.
tObsInUTC – [in] Flag to indicate if the observation time is in UTC (true) or TDB (false).
observerState – [out] Output observer state (position [DU], velocity [DU/TU]).
-
struct Constants#
- #include <simulation.h>
Structure to hold constants used in a PropSimulation.
- Param du2m:
Distance unit conversion factor (default: AU to meters).
- Param tu2s:
Time unit conversion factor (default: days to seconds).
- Param duptu2mps:
Distance unit per time unit conversion factor (default: AU/day to m/s).
- Param G:
Gravitational constant (default: kg AU^3/day^2).
- Param clight:
Speed of light (default: AU/day).
- Param j2000Jd:
Julian date of J2000 epoch.
- Param JdMinusMjd:
Julian date minus modified Julian date.
-
struct IntegrationParameters#
- #include <simulation.h>
Structure to hold integration parameters used in a PropSimulation.
- Param nInteg:
Number of integrated bodies.
- Param nSpice:
Number of SPICE bodies (queried from PropSimulation Ephemeris object).
- Param nTotal:
Total number of bodies (nInteg + nSpice).
- Param n2Derivs:
Number of second derivatives (usually 3 for each body unless STM is propagated).
- Param t0:
Initial time.
- Param tf:
Final time.
- Param dt0:
Initial timestep.
- Param dtMin:
Minimum timestep.
- Param dtChangeFactor:
Maximum factor by which to change timestep.
- Param adaptiveTimestep:
Flag to use adaptive timestep.
- Param timestepCounter:
Timestep counter.
- Param tolPC:
Tolerance for Gauss-Radau predictor-corrector loop.
- Param tolInteg:
Tolerance for integrator.
-
class Body#
- #include <simulation.h>
Parent class for all bodies in a PropSimulation.
- Param t0:
Initial time.
- Param mass:
Mass of the body.
- Param radius:
Radius of the body.
- Param J2:
J2 coefficient for oblateness.
- Param poleRA:
Right ascension of the pole.
- Param poleDec:
Declination of the pole.
- Param name:
Name of the body.
- Param spiceId:
SPICE ID of the body.
- Param pos:
Position of the body.
- Param vel:
Velocity of the body.
- Param acc:
Acceleration of the body.
- Param isPPN:
Flag to indicate if the body is a PPN body (for EIH relativity loop 1).
- Param isJ2:
Flag to indicate if the body has J2 oblateness.
- Param isNongrav:
Flag to indicate if the body has nongravitational forces.
- Param isMajor:
Flag to indicate if the body is a major body (for EIH relativity loop 2).
- Param caTol:
Distance tolerance for close approach detection (default: 0.1 AU).
Subclassed by IntegBody, SpiceBody
Public Functions
-
void set_J2(real J2, real poleRA, real poleDec)#
Set the J2 for the body.
- Parameters:
J2 – [in] Oblateness coefficient.
poleRA – [in] Right ascension of the pole.
poleDec – [in] Declination of the pole.
-
void set_harmonics(real poleRA, real poleDec, int nZon, int nTes, std::vector<real> J, std::vector<std::vector<real>> C, std::vector<std::vector<real>> S)#
Set the harmonic coefficients for the body.
- Parameters:
poleRA – [in] Right ascension of the pole.
poleDec – [in] Declination of the pole.
nZon – [in] Degree of the zonal coefficients.
nTes – [in] Order of the tesseral coefficients.
J – [in] Zonal coefficient vector.
C – [in] Sectoral coefficient array.
S – [in] Tesseral coefficient array.
Public Members
-
real t0#
-
real mass#
-
real radius#
-
real J2 = 0.0L#
-
real poleRA = 0.0L#
-
real poleDec = 90.0L#
-
size_t nZon = 0#
-
size_t nTes = 0#
-
std::vector<real> J = {}#
-
std::vector<std::vector<real>> C = {}#
-
std::vector<std::vector<real>> S = {}#
-
std::string name#
-
int spiceId#
-
real pos[3]#
-
real vel[3]#
-
real acc[3]#
-
bool isPPN = false#
-
bool isJ2 = false#
-
bool isHarmonic = false#
-
bool isNongrav = false#
-
bool isMajor = false#
-
real caTol = 0.1#
-
class SpiceBody : public Body#
- #include <simulation.h>
Class for SPICE bodies in a PropSimulation.
- Param isSpice:
Flag to indicate if the body is a SPICE body.
Public Functions
Public Members
-
bool isSpice = true#
-
struct NongravParameters#
- #include <simulation.h>
Structure to hold parameters for asteroi/comet nongravitational forces.
- Param a1:
Parameter for radial nongravitational forces.
- Param a2:
Parameter for transverse nongravitational forces.
- Param a3:
Parameter for normal nongravitational forces.
- Param a1Est:
Flag to indicate if a1 is estimated (for orbit determination).
- Param a2Est:
Flag to indicate if a2 is estimated (for orbit determination).
- Param a3Est:
Flag to indicate if a3 is estimated (for orbit determination).
- Param alpha:
Normalizing factor for nongravitational forces.
- Param k:
Exponent for nongravitational forces.
- Param m:
Exponent for nongravitational forces.
- Param n:
Exponent for nongravitational forces.
- Param r0_au:
Normalizing distance for nongravitational forces.
-
class IntegBody : public Body#
- #include <simulation.h>
Class for integrated bodies in a PropSimulation.
- Param spiceId:
SPICE ID of the body.
- Param logCA:
Flag to indicate if close approaches are logged (default: true).
- Param isCometary:
Flag to indicate if the body is initialized with a cometary state (as opposed to Cartesian state).
- Param initState:
Initial state of the body.
- Param isInteg:
Flag to indicate if the body is integrated.
- Param isThrusting:
Flag to indicate if the body is thrusting.
- Param ngParams:
Nongravitational parameters for the body.
- Param n2Derivs:
Number of second derivatives for the body.
- Param propStm:
Flag to indicate if the state transition matrix is propagated.
- Param stm:
State transition matrix.
- Param dCartdState:
Derivatives of initial Cartesian state with respect to initial state.
Public Functions
-
IntegBody(std::string name, real t0, real mass, real radius, std::vector<real> cometaryState, NongravParameters ngParams)#
Construct a new IntegBody object.
- Parameters:
name – [in] Name of the body.
t0 – [in] Initial time.
mass – [in] Mass of the body [kg].
radius – [in] Radius of the body [m].
cometaryState – [in] Initial heliocentric ecliptic cometary state of the body.
ngParams – [in] Nongravitational parameters for the body.
-
IntegBody(std::string name, real t0, real mass, real radius, std::vector<real> pos, std::vector<real> vel, NongravParameters ngParams)#
Construct a new IntegBody object.
- Parameters:
name – [in] Name of the body.
t0 – [in] Initial time.
mass – [in] Mass of the body [kg].
radius – [in] Radius of the body [m].
pos – [in] Initial barycentric position of the body [AU].
vel – [in] Initial barycentric velocity of the body [AU/day].
ngParams – [in] Nongravitational parameters for the body.
-
void prepare_stm()#
Prepare the state transition matrix for propagation.
Public Members
-
int spiceId = -99999#
-
bool logCA = true#
-
bool isCometary = false#
-
std::vector<real> initState#
-
std::vector<real> initCart#
-
bool isInteg = true#
-
bool isThrusting = false#
-
NongravParameters ngParams#
-
size_t n2Derivs = 3#
-
bool propStm = false#
-
std::vector<real> stm#
-
std::vector<std::vector<real>> dCartdState#
-
class Event#
- #include <simulation.h>
Class for events in a PropSimulation.
- Param t:
Time of the event.
- Param bodyName:
Name of the IntegBody for the event
- Param bodyIndex:
Index of the IntegBody in the PropSimulation.
- Param xIntegIndex:
Starting index of the body’s state in the flattened state vector.
- Param deltaV:
Delta-V for the event.
- Param multiplier:
Multiplier for the Delta-V.
- Param dt:
Time duration for the continuous event.
- Param isContinuous:
Flag to indicate if the event is continuous.
- Param isHappening:
Flag to indicate if the continuous event is happening.
Public Functions
-
void apply_impulsive(PropSimulation *propSim, const real &t, std::vector<real> &xInteg)#
Apply the impulse event to the body.
- Parameters:
propSim – [in] PropSimulation object.
t – [in] Time of the event.
xInteg – [inout] State of the body.
Public Members
-
real t#
-
std::string bodyName#
-
bool isContinuous = false#
-
bool eventEst = false#
-
size_t bodyIndex#
-
size_t xIntegIndex#
-
bool hasStarted = false#
-
std::vector<real> deltaV = std::vector<real>(3, std::numeric_limits<real>::quiet_NaN())#
-
real multiplier = std::numeric_limits<real>::quiet_NaN()#
-
bool deltaVEst = false#
-
bool multiplierEst = false#
-
std::vector<real> expAccel0 = std::vector<real>(3, std::numeric_limits<real>::quiet_NaN())#
-
real tau = std::numeric_limits<real>::quiet_NaN()#
-
bool expAccel0Est = false#
-
bool tauEst = false#
-
class EventManager#
- #include <simulation.h>
-
struct BPlaneParameters#
- #include <simulation.h>
Structure to hold parameters for a close approach B-plane.
- Param x:
X-coordinate of the B-plane.
- Param y:
Y-coordinate of the B-plane.
- Param z:
Z-coordinate of the B-plane (usually 0).
- Param dx:
Derivatives of x with respect to the state at the close approach.
- Param dy:
Derivatives of y with respect to the state at the close approach.
-
class CloseApproachParameters#
- #include <simulation.h>
Class for close approach parameters in a PropSimulation.
- Param t:
Time of the close approach or impact (from child class).
- Param xRel:
Relative state at the close approach or impact (from child class).
- Param tCA:
Time of the close approach.
- Param xRelCA:
Relative state at the close approach.
- Param tMap:
Time of the B-plane map.
- Param xRelMap:
Relative state at the B-plane map time.
- Param dist:
Distance at the close approach.
- Param vel:
Relative velocity at the close approach.
- Param vInf:
Hyperbolic excess velocity of the close approach.
- Param flybyBody:
Name of the flyby body.
- Param flybyBodyIdx:
Index of the flyby body in the PropSimulation.
- Param centralBody:
Name of the central body.
- Param centralBodyIdx:
Index of the central body in the PropSimulation.
- Param centralBodySpiceId:
SPICE ID of the central body.
- Param impact:
Flag to indicate if the close approach is an impact.
- Param tPeri:
Time of periapsis.
- Param tLin:
Linearized time of intersection on the B-plane.
- Param bVec:
B-vector of the B-plane.
- Param bMag:
Magnitude of the B-vector.
- Param gravFocusFactor:
Gravitational focus factor to account for nonlinear effects from central body.
- Param kizner:
Kizner B-plane parameters.
- Param opik:
Öpik B-plane parameters.
- Param scaled:
Scaled Kizner B-plane parameters.
- Param mtp:
Modified Target Plane parameters.
- Param dtLin:
Partial derivatives of the (linearized intersection minus map) time with respect to the state at the close approach.
- Param dt:
Partial derivatives of the time of closest approach with respect to the state at the close approach.
Subclassed by ImpactParameters
Public Functions
-
void get_ca_parameters(PropSimulation *propSim, const real &tMap)#
Get the close approach parameters.
- Parameters:
propSim – [in] PropSimulation object.
tMap – [in] Time of the B-plane map.
-
void print_summary(int prec = 8)#
Print a summary of the close approach parameters.
- Parameters:
prec – [in] Precision of the output.
Public Members
-
real t#
-
std::vector<real> xRel#
-
real tCA#
-
std::vector<real> xRelCA#
-
real tMap#
-
std::vector<real> xRelMap#
-
real dist#
-
real vel#
-
real vInf#
-
std::string flybyBody#
-
int flybyBodyIdx#
-
std::string centralBody#
-
int centralBodyIdx#
-
int centralBodySpiceId#
-
bool impact#
-
real tPeri#
-
real tLin#
-
std::vector<real> bVec = std::vector<real>(3, 0.0L)#
-
real bMag#
-
real gravFocusFactor#
-
BPlaneParameters kizner#
-
BPlaneParameters opik#
-
BPlaneParameters scaled#
-
std::vector<real> dtLin = std::vector<real>(6, 0.0L)#
-
std::vector<real> dt = std::vector<real>(6, 0.0L)#
-
class ImpactParameters : public CloseApproachParameters#
- #include <simulation.h>
Class for impact parameters in a PropSimulation.
- Param xRelBodyFixed:
Relative state at the impact time in body-fixed frame.
- Param lon:
Longitude of the impact [rad].
- Param lat:
Latitude of the impact [rad].
- Param alt:
Altitude of the impact [km].
Public Functions
-
void get_impact_parameters(PropSimulation *propSim)#
Get the impact parameters.
- Parameters:
propSim – [in] PropSimulation object.
-
void print_summary(int prec = 8)#
Print a summary of the impact parameters.
- Parameters:
prec – [in] Precision of the output.
-
struct InterpolationParameters#
- #include <simulation.h>
Structure to hold interpolation parameters for a PropSimulation.
- Param interpIdx:
Index of the interpolation parameters.
- Param t0:
Initial time of the interpolation window.
- Param dt0:
Timespan of the interpolation window.
- Param b0:
Interpolation coefficients for the interpolation window.
- Param xInteg0:
States at the initial time of the interpolation window.
- Param accInteg0:
Accelerations at the initial time of the interpolation window.
- Param tStack:
Stack of integrator epochs.
- Param xIntegStack:
Stack of states at integrator epochs.
- Param bStack:
Stack of interpolation coefficients at integrator epochs.
- Param accIntegStack:
Stack of accelerations at integrator epochs.
-
class PropSimulation#
- #include <simulation.h>
Class for a PropSimulation.
- Param name:
Name of the simulation.
- Param DEkernelPath:
Path to DE kernels.
- Param spkEphem:
SPK ephemeris object for the simulation.
- Param pckEphem:
PCK ephemeris object for the simulation.
- Param unsafePersistentMemoryMap:
Flag to indicate if ephemeris memory maps should not be unmapped (unsafe, default: false).
- Param consts:
Constants object for the simulation.
- Param integParams:
Integration parameters for the simulation.
- Param spiceBodies:
Vector of SpiceBody objects in the simulation.
- Param integBodies:
Vector of IntegBody objects in the simulation.
- Param events:
Vector of ImpulseEvent objects in the simulation.
- Param caParams:
Vector of CloseApproachParameters objects in the simulation.
- Param impactParams:
Vector of ImpactParameters objects in the simulation.
- Param isPreprocessed:
Flag to indicate if the simulation is preprocessed.
- Param t:
Current time of the simulation.
- Param xInteg:
Current state of the simulation.
- Param interpParams:
Interpolation parameters for the simulation.
- Param tEvalUTC:
Flag to indicate if the evaluation times are in UTC.
- Param evalApparentState:
Flag to indicate if the apparent state is evaluated.
- Param evalMeasurements:
Flag to indicate if measurements are evaluated.
- Param convergedLightTime:
Flag to indicate if converged light time needs to be computed.
- Param xObserver:
Observer states array.
- Param observerInfo:
Observer information array.
- Param tEvalMargin:
Margin for interpolation outside integration arc.
- Param tEval:
Vector of times at which to evaluate the integrated state.
- Param obsType:
Vector of observation type flags (0=none, 1=delay, 2=doppler, 3=gaia).
- Param lightTimeEval:
Vector of computed light times.
- Param xIntegEval:
Vector of integrated states at evaluation times.
- Param opticalObs:
Vector of optical observations.
- Param opticalPartials:
Vector of optical observation partials.
- Param opticalObsCorr:
Vector of photocenter-barycenter corrections for optical observations.
- Param radarObs:
Vector of radar observations.
- Param radarPartials:
Vector of radar observation partials.
Public Functions
-
PropSimulation(std::string name, real t0, const int defaultSpiceBodies, std::string DEkernelPath)#
Construct a new PropSimulation object.
- Parameters:
name – [in] Name of the simulation.
t0 – [in] Initial time.
defaultSpiceBodies – [in] Default SPICE bodies to load. (0=empty sim with DE440 ephemeris, 430/421=DE430 Sun+planets+Moon+Pluto+big16 asteroids, 440/441=DE440 Sun+planets+Moon+Pluto+big16 asteroids)
DEkernelPath – [in] Path to SPICE metakernel
-
PropSimulation(std::string name, const PropSimulation &simRef)#
Construct a new PropSimulation object from a reference simulation.
- Parameters:
name – [in] Name of the simulation.
simRef – [in] Reference simulation to copy.
-
void map_ephemeris()#
Memory map the ephemeris files.
-
void unmap_ephemeris()#
Unmap the ephemeris files.
-
std::vector<real> get_spiceBody_state(const real t, const std::string &bodyName)#
Get the state of a SpiceBody in the PropSimulation at a given time.
- Parameters:
t – [in] Time at which to get the state.
bodyName – [in] Name of the body.
- Returns:
std::vector<real> State of the body at time t.
-
std::vector<real> interpolate(const real t)#
Get the integrated state from the simulation at a given time.
- Parameters:
t – [in] Time to interpolate the integrated state from the PropSimulation at.
- Returns:
std::vector<real> The interpolated geometric barycentric state at time t.
-
void add_spice_body(SpiceBody body)#
Add a SpiceBody to the simulation.
- Parameters:
body – [in] SpiceBody object to add to the simulation.
-
void add_integ_body(IntegBody body)#
Add an IntegBody to the simulation.
- Parameters:
body – [in] IntegBody object to add to the simulation.
-
void remove_body(std::string name)#
Remove a body from the simulation.
- Parameters:
name – [in] Name of the body to remove.
-
void add_event(Event event)#
Add an event to the simulation.
- Parameters:
event – [in] Event object to add to the simulation.
-
void set_sim_constants(real du2m = 149597870700.0L, real tu2s = 86400.0L, real G = 6.6743e-11L / (149597870700.0L * 149597870700.0L * 149597870700.0L) * 86400.0L * 86400.0L, real clight = 299792458.0L / 149597870700.0L * 86400.0L)#
Set the values of the PropSimulation Constants object.
- Parameters:
du2m – [in] Distance unit conversion factor (default: AU to meters).
tu2s – [in] Time unit conversion factor (default: days to seconds).
G – [in] Gravitational constant (default: kg AU^3/day^2).
clight – [in] Speed of light (default: AU/day).
-
void set_integration_parameters(real tf, std::vector<real> tEval = std::vector<real>(), bool tEvalUTC = false, bool evalApparentState = false, bool convergedLightTime = false, std::vector<std::vector<real>> observerInfo = std::vector<std::vector<real>>(), bool adaptiveTimestep = true, real dt0 = 1.0L, real dtMin = 1.0e-4L, real dtChangeFactor = 0.25L, real tolInteg = 1.0e-11L, real tolPC = 1.0e-16L)#
Set the values of the PropSimulation IntegrationParameters object.
- Parameters:
tf – [in] Final time.
tEval – [in] Vector of times at which to evaluate the integrated state (default: empty).
tEvalUTC – [in] Flag to indicate if the evaluation times are in UTC (default: false).
evalApparentState – [in] Flag to indicate if the apparent state is evaluated (default: false).
convergedLightTime – [in] Flag to indicate if converged light time needs to be computed (default: false).
observerInfo – [in] Observer information array (default: empty).
adaptiveTimestep – [in] Flag to use adaptive timestep (default: true).
dt0 – [in] Initial timestep (default: 0.0).
dtMin – [in] Minimum timestep (default: 0.005).
dtChangeFactor – [in] Maximum factor by which to change timestep (default: 0.25).
tolInteg – [in] Tolerance for integrator (default: 1e-11).
tolPC – [in] Tolerance for Gauss-Radau predictor-corrector loop (default: 1e-16).
-
std::vector<real> get_sim_constants()#
Get the constants used in the simulation.
- Returns:
std::vector<real> Vector of simulation constants.
-
std::vector<real> get_integration_parameters()#
Get the integration parameters used in the simulation.
- Returns:
std::vector<real> Vector of simulation integration parameters.
-
void integrate()#
Integrate the simulation to the final time.
-
void extend(real tf, std::vector<real> tEvalNew = std::vector<real>(), std::vector<std::vector<real>> observerInfoNew = std::vector<std::vector<real>>())#
Extend the simulation to a new final time.
- Parameters:
tf – [in] New final time.
tEvalNew – [in] New vector of times at which to evaluate the integrated state.
observerInfoNew – [in] New vector of observer information.
-
void save(std::string filename)#
Save the simulation to a file.
- Parameters:
filename – [in] Name of the file to save the simulation to.
Public Members
-
std::string name#
-
std::string DEkernelPath#
-
SpkEphemeris spkEphem#
-
PckEphemeris pckEphem#
-
bool unsafePersistentMemoryMap = false#
-
IntegrationParameters integParams#
-
EventManager eventMngr#
-
std::vector<CloseApproachParameters> caParams#
-
std::vector<ImpactParameters> impactParams#
-
real t#
-
std::vector<real> xInteg#
-
InterpolationParameters interpParams#
-
bool tEvalUTC = false#
-
bool evalApparentState = false#
-
bool evalMeasurements = false#
-
bool convergedLightTime = false#
-
std::vector<std::vector<real>> xObserver#
-
std::vector<std::vector<real>> observerInfo#
-
real tEvalMargin = 0.0L#
-
std::vector<real> tEval#
-
std::vector<int> obsType#
-
std::vector<std::vector<real>> lightTimeEval#
-
std::vector<std::vector<real>> xIntegEval#
-
std::vector<std::vector<real>> opticalObs#
-
std::vector<std::vector<real>> opticalObsDot#
-
std::vector<std::vector<real>> opticalPartials#
-
std::vector<std::vector<real>> opticalObsCorr#
-
std::vector<std::vector<real>> radarObs#
-
std::vector<std::vector<real>> radarPartials#
Private Functions
-
void prepare_for_evaluation(std::vector<real> &tEval, std::vector<std::vector<real>> &observerInfo)#
Prepare the simulation for evaluation.
- Parameters:
tEval – [in] Vector of times at which to evaluate the integrated state.
observerInfo – [in] Observer information array.
-
void preprocess()#
Preprocess the simulation before integration.
Private Members
-
bool isPreprocessed = false#