From 8f6c0f3d10553388afea2649fa4d09c10d4e37e1 Mon Sep 17 00:00:00 2001 From: cinnaboot Date: Sat, 14 May 2022 10:13:55 -0400 Subject: [PATCH] update orbitGetTimeofFlight for orbits that pass through periapsis --- src/orbits.cpp | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/src/orbits.cpp b/src/orbits.cpp index 01baf43..dbf5fcc 100644 --- a/src/orbits.cpp +++ b/src/orbits.cpp @@ -273,21 +273,21 @@ orbitGetTimeOfFlight(system_2body sys, double theta_begin, double theta_end) { double e = sys.ep.e; - // get mean motion (n) double n = getMeanMotion(sys.body.mu, sys.ep.a); - // get mean anomaly for theta_0 (M0) - double ecc_1 = getEccAnomFromTrueAnom(sys.ep.e, theta_0); + double ecc_begin = getEccAnomFromTrueAnom(sys.ep.e, theta_begin); + double ecc_end = getEccAnomFromTrueAnom(sys.ep.e, theta_end); - // get mean anomaly for theta_1 (M1) - double ecc_2 = getEccAnomFromTrueAnom(sys.ep.e, theta_1); + // NOTE: test if flight passes through perisapsis + if (ecc_begin > ecc_end) + ecc_end += 2 * M_PI; - // test if orbit passes through periapsis, and if so, add 2_PI - // TODO: could also check for multiple revolutions, but that would be - // pretty unlikely if we're updating the simulation every 33ms - // use Kepler's equation for time of flight - double tof_1 = 1 / n * (ecc_1 - e * sin(ecc_1)); - double tof_2 = 1 / n * (ecc_2 - e * sin(ecc_2)); + double M1 = getMeanAnomFromEccAnom(ecc_begin, e); + double M2 = getMeanAnomFromEccAnom(ecc_end, e); - return tof_2 - tof_1; + // NOTE: Kepler's equation for time of flight + double tof_begin = 1 / n * M1; + double tof_end = 1 / n * M2; + + return tof_end - tof_begin; }