Skip to content

Preserve unbounded hyperbolic eccentric anomalies - #1152

Open
antonsoo wants to merge 1 commit into
skyfielders:masterfrom
antonsoo:fix/unbounded-hyperbolic-anomaly
Open

antonsoo wants to merge 1 commit into
skyfielders:masterfrom
antonsoo:fix/unbounded-hyperbolic-anomaly

Conversation

@antonsoo

Copy link
Copy Markdown

A hyperbolic eccentric anomaly can exceed pi in either direction. eccentric_anomaly() currently wraps it into [-pi, pi], which can move the recovered perihelion date by years. This affects scalar and array positions.

For example, a two-body orbit followed ten years past perihelion:

import pandas as pd
from skyfield.api import load
from skyfield.constants import GM_SUN_DE440_km3_s2 as GM_SUN
from skyfield.data import mpc
from skyfield.elementslib import osculating_elements_of

ts = load.timescale()
row = pd.Series(dict(
    eccentricity=1.2, perihelion_distance_au=0.255,
    inclination_degrees=122.7, longitude_of_ascending_node_degrees=24.6,
    argument_of_perihelion_degrees=241.7, perihelion_year=2017,
    perihelion_month=9, perihelion_day=9.0, designation="hyperbolic orbit"))
orbit = mpc.comet_orbit(row, ts, GM_SUN)
elements = osculating_elements_of(orbit.at(orbit.epoch + 3650), gm_km3_s2=GM_SUN)
print(round(float(elements.eccentric_anomaly.radians), 6))
print(round(float(elements.periapsis_time.tdb - orbit.epoch.tdb), 3))

Before: -1.90112 and 3819.499 days of perihelion-date error. After: 4.382066 and 0.0 days at the displayed precision.

Remove the wrapping in both hyperbolic branches. The new regressions parameterize the Cartesian hyperbola independently and check both the anomaly and recovered perihelion, including inbound positions and the +/-pi boundaries.

I also replayed archived JPL elements for 1I/ʻOumuamua and 2I/Borisov against CSPICE conics and an independent scalar Kepler solve. Twenty incorrect cases out of sixty become correct across scalar and batch calls. Rebuilding the orbit from the recovered perihelion and following it thirty more days agrees within 98 metres; all complete position and apparent-observation controls, and the elliptical-element controls, stay identical. This validates two-body processing, not the physical accuracy of a long-term comet ephemeris.

Validation: full Python 3.8 suite and documentation doctests pass. Python 3.12/NumPy 2.5 has the same 41 deprecation-as-error failures on unchanged upstream and this branch; the two new tests pass. The native workload also passes on NumPy 2.2 and 2.5.

@brandon-rhodes

Copy link
Copy Markdown
Member

Interesting, I've not tried an orbit with a hyperbolic eccentric anomaly before. I'm busy this month but will plan to look at this in November, thanks!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants