layup.iod
Pluggable initial orbit determination (IOD) layer for layup.
An IOD method is a callable that proposes one or more seed orbits for the Marquardt fitter to refine. The expected signature is
iod(observations, seq) -> list[FitResult] | None
where observations is the time-ordered list of layup Observation`s, `seq is the per-segment index list (seq[0] is the primary segment used for the IOD), and the return is either a list of candidate seed orbits (each a FitResult with at least state and epoch populated) or None / empty if no candidate could be produced.
Multiple IOD candidates are returned when the underlying method is multi-valued (Gauss’s polynomial in r₂ has up to eight real roots); do_fit runs LM from each candidate and picks the best converged fit (smallest χ² subject to a sanity bound on heliocentric distance).
Methods register themselves at import time via the module-level registry; use register_iod(name, callable) to add new methods, get_iod(name) to look one up, and iod_methods() to list all available names.
Why a registry instead of subclassing: IOD methods are stateless strategies whose entire interface fits on one line. A function-pointer registry is the smallest abstraction that supports drop-in replacements (e.g. a Lambert-based method, a motion-rate prior, a prelim from BK’s tangent-plane linear fit) without forcing each implementation through a class hierarchy.
Attributes
Functions
|
Register an IOD method under name. Overwrites an existing entry. |
|
Look up an IOD method by name. Raises ValueError if unknown. |
|
Return the sorted list of registered IOD method names. |
|
Gauss's method on the first/middle/last observation of seq[0]. |
|
Cheap algebraic feasibility check on an IOD candidate state. |
|
Propagate state to obs.epoch via full ASSIST and return the |
|
Smallest |candidate - observer| over the observation arc, treating |
|
Drop IOD candidates whose predicted positions miss the observations |
Module Contents
- _REGISTRY: dict[str, IODCallable][source]
- register_iod(name: str, func: IODCallable) None[source]
Register an IOD method under name. Overwrites an existing entry.
- get_iod(name: str) IODCallable[source]
Look up an IOD method by name. Raises ValueError if unknown.
- gauss_iod(observations, seq)[source]
Gauss’s method on the first/middle/last observation of seq[0].
The C++ gauss binding returns up to eight candidate seed orbits (corresponding to the real roots of the 8th-degree polynomial in r₂); we pass them all upstream so the picker can pick the right one rather than committing to solns[0] blindly.
- _passes_physical_bounds(candidate, min_r_au: float = _MIN_R_AU, max_r_au: float = _MAX_R_AU) bool[source]
Cheap algebraic feasibility check on an IOD candidate state.
Rejects candidates with non-positive r² or |r| outside [min, max] AU. Deliberately does not reject hyperbolic-looking velocities: Gauss’s velocity can be wildly wrong even for the correct geometric root, and LM walks those to convergence routinely.
- _predict_rho_hat(ephem, state, state_epoch, obs)[source]
Propagate state to obs.epoch via full ASSIST and return the predicted apparent unit direction (no light-time correction; coarse filter only).
- _inertial_min_geocentric_AU(state, state_epoch, observations) float[source]
Smallest |candidate - observer| over the observation arc, treating candidate motion as inertial (position + velocity·Δt).
Used to detect candidates whose trajectory passes close to Earth (or the ground observer); full ASSIST integration would then spend most of its time resolving the close encounter. Inertial-extrapolation is OK for the detection (we just need an order of magnitude); the actual close approach with gravity could be different.
- filter_candidates_by_residual(candidates, observations, ephem, threshold_sigma: float = 1000.0, residual_percentile: float = 80.0, min_obs_for_filter: int = 4, close_earth_AU: float = _CLOSE_EARTH_AU)[source]
Drop IOD candidates whose predicted positions miss the observations by more than threshold_sigma times the per-axis astrometric σ.
The right Gauss root predicts the bulk of the observations within a few σ; phantom roots are typically off by 10⁵-10⁶ σ on essentially every observation. A loose threshold (1000σ default) keeps the right root in every realistic case while throwing out the obviously-wrong ones before LM ever runs on them.
The per-candidate metric is the residual_percentile-th percentile of the per-observation residuals (in σ), not the worst single point. A max-residual criterion is brittle: one contaminating observation (or a rough multi-point seed propagated across a long arc) can throw a single large residual that exceeds the threshold even for the correct root. Using a high percentile (80th by default) tolerates a minority of bad points — which the downstream robust LM fit cleans up — while still rejecting candidates that miss the majority of observations. Keeping it a percentile rather than the median means a candidate must still fit most of the arc, so a seed that only matches one of two mis-linked tracklet groups is not waved through.
Candidates whose inertial trajectory passes within close_earth_AU of the observer at any obs time are passed through unfiltered. Full ASSIST integration gets stuck on close Earth encounters (tens of seconds per propagation), and replacing it with a 2-body approximation would silently mishandle the real physics of NEO close passes — those are valid science targets that need a different solution. Until that solution exists, we just skip the filter for such candidates and let LM handle them (slowly, on the same close encounters, but that’s a separate known issue).
- Parameters:
candidates (list[FitResult]) – Output of an IOD method (states + epoch filled in).
observations (sequence[Observation]) – Full observation list; we evaluate the candidate against every one of them.
ephem (assist.Ephem) – The Python ASSIST ephemeris handle (e.g. assist.Ephem(planets_path, sb_path)). Not the C struct from layup.routines.get_ephem.
threshold_sigma (float) – Reject candidates whose residual_percentile-th-percentile angular residual exceeds this multiple of the per-observation σ.
residual_percentile (float) – Percentile (0-100) of the per-observation residuals used as the per-candidate rejection metric. 80.0 by default; 50.0 gives the median (tolerant of up to half the points being outliers), 100.0 recovers the legacy worst-point behavior.
min_obs_for_filter (int) – Bypass the filter (return all candidates that pass the physical bounds) when there are fewer than this many observations.
close_earth_AU (float) – Pass-through threshold for the close-Earth-approach check described above.
- Returns:
Filtered list. If every candidate fails the residual test, the list of physical-bound-passing candidates is returned instead — we’d rather hand a bad seed to LM than no seed at all.
- Return type:
list[FitResult]