#!/usr/bin/env python3
"""
Precalculate expected values for test_inclined_orbits.cpp.

Usage:
    python3 scripts/precalc_inclined_orbits.py

Outputs C++-style comments with precalculated values for embedding in the test.
Uses scripts/sim_engine.py for the physics engine.
"""

import sys, math
sys.path.insert(0, 'scripts')
from sim_engine import orbital_to_cartesian, vmag, OrbitalElements, G

# =============================================================================
# Molniya orbit
# =============================================================================
a = 26540000.0
e = 0.74
inc = 1.107
omega = 4.71
Omega = 0.0
mu = G * 5.972e24

r_peri = a * (1.0 - e)
r_apo = a * (1.0 + e)
r_90 = a * (1.0 - e*e) / (1.0 + e * math.cos(math.pi/2.0))
r_270 = a * (1.0 - e*e) / (1.0 + e * math.cos(3.0*math.pi/2.0))

T = 2 * math.pi * math.sqrt(a**3 / mu)
T_half = T / 2

print("# Molniya radii:")
print(f"#   r_peri = {r_peri:.6f}")
print(f"#   r_90   = {r_90:.6f}")
print(f"#   r_apo  = {r_apo:.6f}")
print(f"#   r_270  = {r_270:.6f}")
print(f"#")
print(f"# Period: {T:.6f} s = {T/3600:.6f} hours")
print(f"# Half period: {T_half:.6f} s = {T_half/3600:.6f} hours")

# =============================================================================
# Generic inclined orbit
# =============================================================================
a2 = 10000000.0
e2 = 0.5
inc2 = math.radians(45)
omega2 = math.pi / 2

elements2 = OrbitalElements(a=a2, e=e2, nu=0.0, inc=inc2, Omega=0.0, omega=omega2)
pos2, vel2 = orbital_to_cartesian(elements2, 5.972e24)
r2 = vmag(pos2)
z2 = pos2[2]

print(f"\n# Generic inclined (a={a2}, e={e2}, i=45deg, omega=90deg):")
print(f"#   r = {r2:.6f} m")
print(f"#   z = {z2:.6f} m")