diff --git a/src/orbital_mechanics.cpp b/src/orbital_mechanics.cpp index 9258ba7..8e5f7e3 100644 --- a/src/orbital_mechanics.cpp +++ b/src/orbital_mechanics.cpp @@ -166,13 +166,25 @@ OrbitalElements cartesian_to_orbital_elements(Vec3 position, Vec3 velocity, doub } else { double cos_nu = r_dot_e / (r_mag * e * mu); cos_nu = fmax(-1.0, fmin(1.0, cos_nu)); - true_anomaly = acos(cos_nu); - if (vec3_dot(r_vec, v_vec) < 0.0) { - true_anomaly = 2.0 * M_PI - true_anomaly; + double sin_nu; + + if (fabs(cos_nu) > 1.0 - 1e-10) { + Vec3 h_cross_e = vec3_cross(h_vec, e_vec); + double denom = r_mag * e * h; + if (denom > 1e-10) { + sin_nu = vec3_dot(r_vec, h_cross_e) / denom; + } else { + sin_nu = 0.0; + } + } else { + Vec3 r_cross_h = vec3_cross(r_vec, h_vec); + double denom = r_mag * e * h * mu; + sin_nu = (denom > 1e-10) ? vec3_dot(r_cross_h, e_vec) / denom : 0.0; } - // Normalize to (-π, π] range - if (true_anomaly > M_PI) { - true_anomaly -= 2.0 * M_PI; + + true_anomaly = atan2(sin_nu, cos_nu); + if (true_anomaly == -M_PI) { + true_anomaly = M_PI; } }