Orbit Fitting

Orbit fitting and observation ingestion.

Fitting tools:

  • Initial Orbit Determination (IOD) (initial_orbit_determination()) – find candidate orbits from a small number of observations using statistical ranging. Returns scored candidates sorted best-first.

  • Orbit Fitting (fit_orbit()) – refine an orbit guess to best match the data, with automatic outlier rejection. Produces a best-fit state and Gaussian uncertainty estimate (covariance).

  • MCMC Uncertainty Estimation (fit_orbit_mcmc()) – for short arcs where the Gaussian approximation is unreliable, sample the full range of plausible orbits consistent with the data.

Observation ingestion:

class kete.orbit_fitting.MPCObservation(desig, prov_desig, discovery, note1, note2, jd, ra, dec, mag_band, catalog_code, obs_code, sun2sc)

Representation of an 80-character format MPC observation.

import kete
import gzip

# Comet Observations
# url = "https://www.minorplanetcenter.net/iau/ECS/MPCAT-OBS/CmtObs.txt.gz"

# Download the database of unnumbered observations from the MPC
url = "https://www.minorplanetcenter.net/iau/ECS/MPCAT-OBS/UnnObs.txt.gz"
url = "https://www.minorplanetcenter.net/iau/ECS/MPCAT-OBS/NumObs.txt.gz"
path = kete.data.download_file(url)

# Fetch all lines from the file which contain C51 (WISE) observatory code.
obs_code = "C51".encode()
with gzip.open(path) as f:
    lines = [line.decode() for line in f if obs_code == line[77:80]]

# Parse lines into a list of MPCObservations
observations = kete.observations.MPCObservation.from_lines(lines)
Parameters:

  • desig (str)

  • prov_desig (str)

  • discovery (bool)

  • note1 (str)

  • note2 (str)

  • jd (float)

  • ra (float)

  • dec (float)

  • mag_band (str)

  • catalog_code (str)

  • obs_code (str)

  • sun2sc (list[float])

classmethod from_lines(lines, load_sc_pos=True)

Create a list of MPCObservations from a list of single 80 char lines.

property sc2obj
class kete.orbit_fitting.Observation

Astronomical observation for orbit determination.

Observations can be optical (RA/Dec), radar range, or radar range-rate. Each carries the observer state, measured value(s), and 1-sigma uncertainties.

Use the static methods to construct instances:

band

Photometric filter name (optical only, empty string for radar).

dec

Declination in degrees (optical only, None otherwise).

epoch

The observation epoch. For radar this is the receive epoch t_rx.

is_occultation

True if this optical measurement was derived from a stellar occultation (optical only, None for radar).

mag

Apparent magnitude (NaN when unavailable).

observer

The observer state (Sun-centered, Ecliptic).

For radar observations this returns the receiver state at the receive epoch (computed via PCK lookup on demand).

static optical(observer, ra, dec, sigma_ra, sigma_dec, band=Ellipsis, mag=Ellipsis, time_sigma=0.1, sigma_corr=0.0, is_occultation=False)

Create an optical (RA/Dec) observation.

The observer state is automatically re-centered to the solar system barycenter if needed (the fitting engine works internally in SSB-centered coordinates).

Parameters:

  • observer (State) – Observer state (any center / frame - will be converted to SSB-centered Equatorial internally). The observation epoch is taken from the observer’s epoch.

  • ra (float) – Right ascension in degrees.

  • dec (float) – Declination in degrees.

  • sigma_ra (float) – 1-sigma sky-plane RA uncertainty in arcseconds (sigma_ra * cos(dec)). This matches the convention used by all standard astrometric data formats (MPC ADES rmsRA, Gaia DR3 ra_error_*, etc.). The cos(dec) factor is removed internally to convert to the RA coordinate direction used during chi-squared accumulation.

  • sigma_dec (float) – 1-sigma Dec uncertainty in arcseconds.

  • band (str) – Photometric filter name (default "V").

  • mag (float) – Apparent magnitude (default NaN).

  • time_sigma (float) – 1-sigma timing uncertainty in seconds (default 0.5). When non-zero the along-track positional uncertainty is inflated by time_sigma * apparent_speed, producing a non-diagonal weight matrix.

  • sigma_corr (float) – Correlation coefficient between RA and Dec uncertainties (default 0.0, i.e. axis-aligned). Must lie in (-1, 1). Used when ADES-format astrometry reports rmscorr or Gaia DR3 reports a correlation between RA and Dec errors – the position covariance ellipse is tilted accordingly.

  • is_occultation (bool) – Set to True when the measurement was derived from a stellar occultation rather than direct astrometry. Default False. The fitting math is identical for both, but ingestion code uses this flag to skip the per-observatory residual table, EFCC18 debiasing, and over-observation reweighting – which do not apply to occultations.

ra

Right ascension in degrees (optical only, None otherwise).

static radar_range(xmit_lat, xmit_lon, xmit_height, rcvr_lat, rcvr_lon, rcvr_height, epoch, range, sigma_range)

Create a radar range observation.

The transmitter and receiver are stored as WGS84 geodetic coordinates. Their inertial states are computed inside the residual via PCK kernel lookups at the receive epoch (and at the iteratively-refined transmit epoch for the xmit station). Use the same site coordinates for monostatic observations.

Parameters:

  • xmit_lat (float) – Transmitter geodetic latitude in degrees.

  • xmit_lon (float) – Transmitter geodetic longitude in degrees.

  • xmit_height (float) – Transmitter height above the WGS84 ellipsoid in km.

  • rcvr_lat (float) – Receiver geodetic latitude in degrees.

  • rcvr_lon (float) – Receiver geodetic longitude in degrees.

  • rcvr_height (float) – Receiver height above the WGS84 ellipsoid in km.

  • epoch (Time) – Receive epoch t_rx (TDB).

  • range (float) – Measured range in AU.

  • sigma_range (float) – 1-sigma range uncertainty in AU.

static radar_rate(xmit_lat, xmit_lon, xmit_height, rcvr_lat, rcvr_lon, rcvr_height, epoch, range_rate, sigma_range_rate)

Create a radar range-rate (Doppler) observation.

See Observation.radar_range() for the station-coordinate convention.

Parameters:

  • xmit_lat (float) – Transmitter geodetic latitude in degrees.

  • xmit_lon (float) – Transmitter geodetic longitude in degrees.

  • xmit_height (float) – Transmitter height above the WGS84 ellipsoid in km.

  • rcvr_lat (float) – Receiver geodetic latitude in degrees.

  • rcvr_lon (float) – Receiver geodetic longitude in degrees.

  • rcvr_height (float) – Receiver height above the WGS84 ellipsoid in km.

  • epoch (Time) – Receive epoch t_rx (TDB).

  • range_rate (float) – Measured range-rate in AU/day (positive = receding).

  • sigma_range_rate (float) – 1-sigma range-rate uncertainty in AU/day.

range

Measured range in AU (radar range only, None otherwise).

range_rate

Measured range-rate in AU/day (radar rate only, None otherwise).

residual(obj_state)

Compute the residual (observed - predicted) and the predicted measurement for this observation against a given object state.

The state is expected at this observation’s epoch (use propagate_n_body() to bring an arbitrary-epoch state to self.epoch first). Light-time correction is applied internally; for radar the iterative t_tx refinement and relativistic two-way Doppler are handled by the same code path used during fitting.

Parameters:

obj_state (State) – Object state at the observation epoch (any center / frame).

Returns:

(residual, predicted).

  • Optical: 2-element lists. residual is [dRA_sky, dDec] in arcseconds, sky-plane convention (dRA includes cos(dec), matching Observation.optical()’s sigma_ra units). predicted is [ra_deg, dec_deg].

  • Radar range: 1-element lists in AU.

  • Radar range-rate: 1-element lists in AU/day.

Return type:

tuple[list[float], list[float]]

sigma_corr

RA/Dec correlation coefficient (optical only, None otherwise).

sigma_dec

1-sigma Dec uncertainty in arcseconds (optical only, None otherwise).

sigma_ra

1-sigma sky-plane RA uncertainty in arcseconds (sigma_ra * cos(dec)); optical only, None otherwise. Matches the input convention of Observation.optical().

sigma_range

1-sigma range uncertainty in AU (radar range only, None otherwise).

sigma_range_rate

1-sigma range-rate uncertainty in AU/day (radar rate only, None otherwise).

time_sigma

1-sigma timing uncertainty in seconds (optical only, None otherwise).

class kete.orbit_fitting.OrbitFit

Result of fitting an orbit to observations.

Returned by fit_orbit(). Contains the best-fit orbital state, its uncertainty (covariance), and diagnostic information about the fit.

all_observations

All input observations (time-sorted), including rejected outliers.

converged

Whether the solver achieved strict convergence.

When False the fit is the best found within the iteration limit but the correction norm did not drop below tol.

included

Per-observation inclusion mask (time-sorted).

True means the observation was used in the final fit; False means it was rejected as an outlier.

non_grav

Fitted non-gravitational model, or None if not fitted.

Convenience shortcut for self.uncertain_state.non_grav.

observations

Observations included in the final fit (time-sorted).

Rejected outliers are not present in this list. To see all observations (including rejected ones), use all_observations.

Band and magnitude metadata from the original Observation inputs are preserved.

residuals

Post-fit residuals as a list of lists (included observations only, time-sorted order).

For optical observations the two elements are (DeltaRA, DeltaDec) in arcseconds. Radar residuals remain in AU or AU/day.

rms

Reduced weighted RMS of post-fit residuals.

state

Best-fit state at the reference epoch (Sun-centered, Ecliptic).

Convenience shortcut for self.uncertain_state.state.

uncertain_state

The uncertain orbit state (state + covariance + non-grav model).

class kete.orbit_fitting.OrbitSamples

Collection of plausible orbits from MCMC uncertainty estimation.

Returned by fit_orbit_mcmc(). Each orbit (draw) is statistically consistent with the observations; the spread of the collection represents the uncertainty in the orbit.

desig

Designator of the fitted object.

divergent

Per-draw divergence flag.

A divergent sample indicates the sampler had difficulty exploring that region of orbit space. A small fraction of divergences is normal; many divergences suggest the model or data are problematic.

draws

Sampled orbits as a list of State objects.

Each state is Sun-centered Ecliptic at the reference epoch.

Non-gravitational parameters (if fitted) are available via raw_draws.

epoch

Common reference epoch (JD, TDB).

log_posterior

Log-posterior density at each draw (nats, relative only).

Values are only meaningful relative to each other within a single run. Useful for weighting draws or diagnosing chain quality.

raw_draws

Raw orbit samples as a list of lists.

Each inner list is [x, y, z, vx, vy, vz, ng_params...] in the Equatorial frame at the reference epoch. Use np.array(samples.raw_draws) to convert to a 2-D array.

seed_id

Seed index (0-based) that generated each draw.

Each seed produces an independent NUTS chain. This field records which seed produced each draw so per-seed diagnostics can be computed.

class kete.orbit_fitting.RangingSamples

Orbit samples from the weighted ranging grid.

Returned by fit_orbit_ranging().

convergence_warning

Warning message if ESS < 50, else None.

draws

Sampled orbits as State objects (Sun-centered Ecliptic).

effective_sample_size

Effective sample size of the grid before drawing.

epoch

Reference epoch (JD TDB).

log_posterior

Normalized log-posterior weight per draw.

raw_draws

Raw draws as [x, y, z, vx, vy, vz] in AU/AU*day, SSB Equatorial.

class kete.orbit_fitting.UncertainState

Uncertain orbit state: a best-fit Cartesian state together with a covariance matrix that may span the 6 position/velocity components and any fitted non-gravitational parameters.

This is the canonical representation of orbit uncertainty in kete, providing correct handling of non-gravitational model templates (preserving fixed physical parameters such as alpha, r_0, etc.).

Construction

  • from_state() – from a State with isotropic uncertainties.

  • from_cometary() – from cometary orbital elements and an element-space covariance (e.g. from JPL Horizons).

  • Returned as part of OrbitFit from orbit fitting.

cov_matrix

Covariance matrix as a list of lists (use np.array() to convert).

desig

Object designator (shortcut for self.state.desig).

epoch

Reference epoch as a Time (shortcut for self.state.epoch).

static from_cometary(elements, cov_matrix, non_grav=None)

Build an UncertainState from cometary orbital elements and a covariance expressed in element space.

The element-space covariance is transformed to a Cartesian covariance via a numerically evaluated Jacobian.

Parameters:

  • elements (CometElements) – Cometary orbital elements (with desig and epoch).

  • cov_matrix (list[list[float]]) – Covariance matrix in element space, (6+Np)x(6+Np). Element order: [e, q, tp, node, w, i, <nongrav...>].

  • non_grav (NonGravModel, optional) – Non-gravitational model template.

static from_state(state, pos_sigma=0.01, vel_sigma=0.0001, non_grav=None)

Build an UncertainState from a state with isotropic diagonal uncertainties.

The covariance is initialized to a diagonal matrix with the given pos_sigma (AU) and vel_sigma (AU/day) on the diagonal. Useful for seeding MCMC from an IOD candidate.

The input state is automatically re-centered to the solar system barycenter if needed.

Parameters:

  • state (State) – Object state (any center / frame – will be converted to SSB-centered Equatorial internally).

  • pos_sigma (float) – 1-sigma position uncertainty in AU (default 0.01).

  • vel_sigma (float) – 1-sigma velocity uncertainty in AU/day (default 0.0001).

  • non_grav (NonGravModel, optional) – Non-gravitational model template. If provided, the covariance is extended to (6+Np)x(6+Np) with tiny diagonal entries for the non-grav parameters.

non_grav

Non-gravitational model template, or None.

param_names

Names of all parameters in the covariance matrix, in row/column order.

Always starts with ["x", "y", "z", "vx", "vy", "vz"], followed by any non-gravitational parameter names.

sample(n_samples, seed=None)

Draw random samples from the covariance distribution.

Returns a tuple (states, non_gravs) where states is a list of State objects and non_gravs is a list of NonGravModel or None.

Parameters:

  • n_samples (int) – Number of samples to draw.

  • seed (int) – Random seed for reproducibility (optional).

state

Best-fit state at the reference epoch (Sun-centered, Ecliptic).

kete.orbit_fitting.fetch_gaia_observations(desig, update_cache=False)

Fetch Gaia DR3 solar system object observations and convert to fitting Observations.

Queries the gaiadr3.sso_observation table via the Gaia TAP service for the given object and returns one Observation per accepted astrometric transit. Only transits with astrometric_outcome_transit == 1 (good positions) are returned.

The Gaia spacecraft position and velocity are taken directly from the table (barycentric ICRS, SSB-centered) so no SPICE lookup is needed.

Observation epoch is the Gaia-centric TCB epoch stored in the table, converted to TDB via Time with scaling='tcb'.

Positional uncertainties are the quadrature sum of the random and systematic components from the table (in mas, already multiplied by cos(dec) for the RA component).

Results are cached via query_tap(); pass update_cache=True to force a fresh query.

Parameters:

  • desig (str) – Object identifier as recognized by Gaia DR3. If the string parses as an integer it is matched against the number_mp column (recommended for numbered minor planets); otherwise it is matched against the denomination column (e.g. "Apophis", "1999 RQ36").

  • update_cache (bool) – If True, discard any cached result and re-query the TAP service.

Returns:

One Observation.optical per accepted transit observation.

Return type:

list[Observation]

Examples

import kete

observations = kete.observations.fetch_gaia_observations("Apophis")
fit = kete.fitting.fit_orbit(initial_state, observations)
kete.orbit_fitting.fetch_mpc_observations(desig, use_observatory_residuals=True, debias=True, apply_over_obs_reweight=True, update_cache=False)

Fetch observations from the MPC API and convert to fitting Observations.

Queries https://data.minorplanetcenter.net/api/get-obs for the given object designation and returns optical observations ready for orbit fitting. Only optical records with valid RA/Dec are included; radar and other types are silently skipped.

Results are cached under ~/.kete/observations/ so that repeated queries for the same designation do not hit the network.

Uncertainties are first taken from a pre-computed table of residual error by observatory code, if available. If not, the MPC-provided rmsra and rmsdec fields are used, defaulting to 1 arcsecond if not provided.

When debias is True, the EFCC18 star-catalog bias correction is applied to each observation’s (RA, Dec) using the ADES astCat field. Observations with unknown or missing catalog codes are passed through unchanged.

When apply_over_obs_reweight is True, sigma is inflated by sqrt(n/4) for groups of more than 4 observations from the same observatory on the same night, following Veres et al. 2017.

Parameters:

  • desig (str) – Object designation recognized by the MPC (e.g. "Apophis", "101955", "1999 RQ36").

  • use_observatory_residuals (bool) – If True, apply per-observatory sigmas from the pre-computed residual table when available; otherwise, use the MPC-provided rmsra and rmsdec fields, defaulting to 1 arcsecond.

  • debias (bool) – If True, apply the EFCC18 star-catalog debiasing correction. Requires the JPL debias_2018.tgz archive, downloaded on first use.

  • apply_over_obs_reweight (bool) – When True (default), inflate sigma for over-observed nights following Veres et al. 2017.

  • update_cache (bool) – If True, discard any cached result and re-query the MPC.

Returns:

One Observation.optical per valid optical record.

Return type:

list[Observation]

Raises:

RuntimeError – If the MPC API request fails.

Examples

import kete

observations = kete.observations.fetch_mpc_observations("Apophis")
kete.orbit_fitting.fetch_observations(desig, update_cache=False)

Fetch all available observations for an object as fitting Observations.

This is a convenience wrapper around the various specific fetchers that tries to get “everything” for a given object. Currently this means MPC ADES optical, JPL radar, and Gaia DR3 measurements.

MPC observations have the EFCC18 debiasing correction applied and per-observatory sigmas when available.

Parameters:

  • desig (str) – Object designation (e.g. "Apophis", "101955", "1999 RQ36").

  • update_cache (bool) – If True, refresh any cached data before filtering.

Returns:

All available observations for the object.

Return type:

list[Observation]

kete.orbit_fitting.fetch_radar_observations(desig, update_cache=False)

Fetch JPL radar observations for an object as fitting Observations.

Wraps fetch_radar_table() and converts each record into either a Observation.radar_range() or Observation.radar_rate() instance ready for orbit fitting.

Each record is converted using the station as the observer state. For bistatic measurements (transmitter != receiver) the observer is placed at the midpoint of the two stations; this approximates the bistatic round-trip path to second order in (baseline / range) and is accurate to a few meters for Earth-baseline NEO radar. Records with missing coordinates or unrecognized units are skipped.

Conversions assumed:

  • Delay (us): round-trip light time. One-way range in AU is c * (delay_us * 1e-6) / 2 / AU_M.

  • Doppler (Hz): two-way Doppler shift at carrier frequency freq (MHz). Range-rate in AU/day is -c * doppler / (2 * freq_Hz) converted from m/s. JPL convention: positive Doppler means approaching; kete convention: positive range-rate means receding, so the sign is flipped.

Parameters:

  • desig (str) – Object designation as it appears in the JPL des field (e.g. "99942", "1998 KY26").

  • update_cache (bool) – If True, refresh the cached master parquet before filtering.

Returns:

One Observation.radar_range or Observation.radar_rate per valid record.

Return type:

list[Observation]

kete.orbit_fitting.fetch_radar_table(desig=None, update_cache=False)

Fetch JPL Small-Body radar astrometry data.

The complete JPL Small-Body Radar Astrometry data set is downloaded once and cached locally as a parquet file. Subsequent calls reuse the cached file unless update_cache is set. If desig is supplied, only the records whose des field matches that designation are returned.

The returned table includes the default API fields, the observer list for each measurement, and geodetic station coordinates joined onto each record for both the receiver (rcvr_*) and transmitter (xmit_*) stations.

See https://ssd-api.jpl.nasa.gov/doc/sb_radar.html for field definitions.

Parameters:

  • desig (str | None) – Object designation as it appears in the JPL des field (e.g. "99942", "1998 KY26", "45P"). If None, every record in the data set is returned.

  • update_cache (bool) – Force download a new file from JPL, refreshing the cached parquet.

kete.orbit_fitting.fit_orbit(initial_state, observations, non_grav=None, include_asteroids=False, chi2_threshold=4.5, max_reject_passes=10)

Fit an orbit to observations using iterative least squares.

Given an initial guess and a set of observations, this function refines the orbital state until it best matches the data, and estimates the uncertainty of the result via a covariance matrix. It can also automatically identify and reject outlier observations.

This is the standard approach for orbit determination when you have a reasonable initial guess (e.g. from initial_orbit_determination()). It works well for arcs of any length, but the Gaussian uncertainty estimate is most reliable for long, well-sampled arcs. For short arcs where the uncertainty is non-Gaussian, consider fit_orbit_mcmc() instead.

For arcs longer than 180 days, progressively wider time windows are fitted around the reference epoch so that each stage bootstraps from the previous converged solution. The final pass fits the full arc and re-evaluates all observations for outlier rejection (if enabled).

The input state is automatically re-centered to the solar system barycenter internally.

Parameters:

  • initial_state (State) – Initial guess for the object state (any center / frame).

  • observations (list) – List of Observation to fit.

  • non_grav (NonGravModel, optional) – Non-gravitational force model.

  • include_asteroids (bool) – If True, include asteroid masses in the force model (slower but more accurate for near-Earth objects). Default is False.

  • chi2_threshold (float) – Per-observation z-score threshold for outlier rejection following Carpino, Milani & Chesley (2003). An observation is rejected when its weighted chi-squared exceeds this multiple of the leverage- corrected expected value, and recovered when it falls below 6/7 of this threshold. For a 2-component optical observation with calibrated weights, a 3-sigma outlier in one axis gives z = 9/2 = 4.5, so values of 4-5 match per-component 3-sigma rejection. Default is 4.5.

  • max_reject_passes (int) – Maximum number of outlier-rejection passes. Set to 0 to disable outlier rejection entirely. Default is 10.

Returns:

The fitted orbit, including the best-fit state, covariance, residuals, and convergence diagnostics.

Return type:

OrbitFit

kete.orbit_fitting.fit_orbit_mcmc(seeds, observations, include_asteroids=False, num_draws=1000, num_tune=500, non_grav=None, maxdepth=10, target_accept=0.8)

Estimate orbit uncertainty from observations using Markov Chain Monte Carlo.

Given one or more candidate orbital states (seeds) and a set of observations, this function produces a collection of plausible orbits that are statistically consistent with the data. The spread of returned orbits represents the uncertainty in the orbit determination – wider spread means less certainty about the true orbit.

This is most useful for short-arc observations (a few nights) where the usual least-squares approach (fit_orbit()) underestimates the true uncertainty. For well-observed objects with long arcs, fit_orbit() alone is usually sufficient and far cheaper.

Under the hood this uses the No-U-Turn Sampler (NUTS), an adaptive variant of Hamiltonian Monte Carlo that efficiently explores the space of possible orbits. Each draw requires a full numerical propagation, so this is orders of magnitude more expensive than least squares – but the result captures non-Gaussian and multi-modal uncertainty that least squares cannot represent.

Seeds are raw State objects (typically from initial_orbit_determination()). No prior orbit fit is required – the sampler builds its own internal mass matrix from a linearization at each seed.

Sampling is parallelized automatically across available CPU cores. When there are fewer seeds than cores, each seed spawns multiple independent sub-chains. The seed_id in the returned OrbitSamples identifies which seed (orbital mode) produced each draw.

num_draws is the total number of orbit samples returned across all seeds. Each seed receives roughly num_draws / len(seeds) draws.

Parameters:

  • seeds (list) – List of State candidate states (e.g. from initial_orbit_determination()), one per orbital mode. Seeds at different epochs are automatically propagated to the first seed’s epoch. The input states are re-centered to SSB Equatorial internally.

  • observations (list) – List of Observation to evaluate against.

  • include_asteroids (bool) – If True, include asteroid masses in the force model. Default is False.

  • num_draws (int) – Total orbit samples across all seeds (after tuning). Default is 1000.

  • num_tune (int) – Number of warmup steps per sub-chain used to adapt internal sampling parameters. These draws are discarded. Default is 500.

  • non_grav (NonGravModel, optional) – Shared non-gravitational force model applied to all chains.

  • maxdepth (int) – Maximum tree depth for the sampler. Higher values allow more thorough exploration at greater computational cost. Default is 10.

  • target_accept (float) – Target acceptance probability for step-size adaptation during warmup. Lowering this (e.g. 0.6) makes the sampler take larger leapfrog steps, which helps in poorly constrained situations. Default is 0.8.

Returns:

Collection of plausible orbits sampled from the posterior.

Return type:

OrbitSamples

Raises:

ValueError – If seeds is empty or two-body epoch propagation fails.

kete.orbit_fitting.fit_orbit_ranging(observations, num_draws=1000, temperature=10.0, seed=0)

Generate orbit samples covering the full admissible region from sparse observations.

Scans a grid over topocentric distances at a selected pair of observations, scores each cell via Gaussian chi^2, adaptively refines in high-probability regions, and draws samples with within-cell Gaussian perturbation.

This is the appropriate tool when the arc is short and MCMC cannot explore the ridge, or when multiple orbital families may be consistent with the data.

Parameters:

  • observations (list) – At least 2 Observation objects.

  • num_draws (int) – Number of orbit samples to return. Default 1000.

  • temperature (float) – Likelihood temperature. 1.0 gives the true Bayesian posterior. Higher values produce a softer distribution that is easier to sample but less statistically rigorous. Default is 10.0, producing results similar to JPL Scout.

  • seed (int) – RNG seed for reproducibility. Default 0.

Return type:

RangingSamples

kete.orbit_fitting.initial_orbit_determination(observations)

Compute an initial orbit from observations.

Returns candidate orbital states sorted by ascending residual score (best first). The score is the trimmed-mean angular residual (radians squared) plus a soft eccentricity penalty.

Parameters:

observations (list[Observation]) – At least 2 optical observations.

Returns:

(score, state) pairs sorted by ascending score. States are at the midpoint of the selected observational arc.

Return type:

list[tuple[float, State]]

kete.orbit_fitting.lambert(r1, r2, dt, prograde=True, max_revs=0)

Solve Lambert’s problem for a Keplerian transfer.

Given two heliocentric position vectors and a transfer time, compute the velocity vectors at departure and arrival that connect them via two-body (Keplerian) motion.

When max_revs is 0 only the zero-revolution (direct) solution is returned. For max_revs > 0 the solver also finds multi-revolution solutions up to that limit; the result list contains one pair per solution, ordered by ascending number of revolutions.

Parameters:

  • r1 (Vector) – Heliocentric position at departure (AU).

  • r2 (Vector) – Heliocentric position at arrival (AU).

  • dt (float) – Transfer time in days. Must be positive.

  • prograde (bool) – If True (default), selects the short-way transfer (transfer angle less than 180 degrees for prograde orbits). If False, selects the long-way transfer.

  • max_revs (int) – Maximum number of complete revolutions to include. Default is 0 (single direct transfer only).

Returns:

[(v1, v2), ...] – velocity pairs at r1 and r2 for each solution found.

Return type:

list[tuple[Vector, Vector]]

Raises:

  • ValueError – If dt <= 0, positions are zero-length, or positions are nearly collinear (transfer angle near 0 or 180 degrees).

  • RuntimeError – If the iterative solver fails to converge.

kete.orbit_fitting.mpc_obs_to_observations(mpc_obs, apply_over_obs_reweight=True, debias=True)

Convert a list of MPCObservation objects to fitting Observations.

Only optical (RA/Dec) observations are supported. Each MPCObservation is converted to an Observation.optical with the observer state computed from the MPC observatory code (ground-based) or from the stored spacecraft position.

Per-observatory uncertainties are applied when available from the pre-computed residual table. When no table entry exists for an observatory code, an epoch- and observation-type-based fallback is used (see _time_sigma_for_obs()).

When debias is True, the EFCC18 star-catalog bias correction is applied using the catalog code stored on each observation (column 72 of the 80-char format). Observations with an unknown or blank catalog code are passed through unchanged.

Parameters:

  • mpc_obs (list) – List of MPCObservation objects (see kete.observations).

  • apply_over_obs_reweight (bool) – When True (default), inflate sigma by sqrt(n/4) for groups of more than 4 observations from the same observatory on the same night, following Veres et al. 2017. Spacecraft observations are exempt.

  • debias (bool) – When True (default), apply the EFCC18 star-catalog bias correction. Requires the JPL debias_2018.tgz archive, downloaded on first use.

Returns:

One Observation.optical per input observation.

Return type:

list[Observation]

Examples

import kete

lines = [...]  # 80-char MPC observation lines
mpc_obs = kete.observations.MPCObservation.from_lines(lines)
observations = kete.observations.mpc_obs_to_observations(mpc_obs)
fit = kete.fitting.fit_orbit(initial_state, observations)