Browse Source

fix: Correct true anomaly normalization to handle negative values

- Replace acos() + quadrant handling with atan2(sin_nu, cos_nu)
- Calculate sin_nu using cross products for proper quadrant determination
- atan2 inherently returns values in (-π, π] range
- Add edge case handling for cos_nu near ±1
- Fix true_anomaly == -π normalization to +π
- Test improvement: 78/82 → 80/82 tests passing
main
cinnaboot 5 months ago
parent
commit
7137b51926
  1. 24
      src/orbital_mechanics.cpp

24
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;
}
}

Loading…
Cancel
Save