24 Commits

Author SHA1 Message Date
9158c7c55a Add v0.11.0: make_orbital_elements constructors + moon equatorial functions
6 new SQL functions (114 -> 120):
- make_orbital_elements(): construct from 9 floats, angles in radians
- make_orbital_elements_deg(): same with angles in degrees, matches
  text I/O convention and typical catalog column layouts
- galilean_equatorial(): geocentric RA/Dec for Io/Europa/Ganymede/Callisto
- saturn_moon_equatorial(): geocentric RA/Dec for Mimas through Hyperion
- uranus_moon_equatorial(): geocentric RA/Dec for Miranda through Oberon
- mars_moon_equatorial(): geocentric RA/Dec for Phobos/Deimos

Constructors requested by astrolock-api to replace fragile
format(9 args)::orbital_elements cast pattern. Moon equatorial
functions fill the last NULL RA/Dec gaps in their unified sky query.

All 20 regression suites pass.
2026-02-23 14:07:39 -07:00
b0741c553b Add v0.10.0: aberration, DE apparent, angular separation, stellar parallax
Annual stellar aberration (~20 arcsec) added to all 6 existing _apparent()
functions via classical first-order v/c projection (Ron & Vondrak). Earth
velocity sourced from VSOP87 xyz[3..5] (analytic) or DE numerical
differentiation.

New functions (106 -> 114):
- eq_angular_distance(): Vincenty formula, stable at 0 and 180 deg
- eq_within_cone(): cosine shortcut for fast cone-search predicate
- <-> operator on equatorial type
- 6 DE apparent variants with VSOP87 fallback:
  planet/sun/moon_observe_apparent_de(),
  planet/moon_equatorial_apparent_de(),
  small_body_observe_apparent_de()

Stellar parallax now functional in star_observe_pm() and
star_equatorial_pm() — Green (1985) Eq. 11.3 displacement using
Earth heliocentric position from VSOP87.

All 19 regression suites pass (18 existing + new aberration suite).
2026-02-21 21:47:42 -07:00
b33d63034b Add v0.9.0 apparent position features: equatorial type, refraction, proper motion, light-time
New equatorial type (24 bytes: RA/Dec/distance) captures apparent coordinates
of date — what the observation pipeline computes at precession step 3 but was
discarding before hour angle conversion. Matches telescope GoTo mount conventions.

24 new SQL functions (82 → 106 total):
- equatorial type I/O + 3 accessors (eq_ra, eq_dec, eq_distance)
- Satellite RA/Dec: eci_to_equatorial (topocentric), eci_to_equatorial_geo (geocentric)
- Solar system equatorial: planet/sun/moon/small_body_equatorial
- Atmospheric refraction: Bennett (1982) with domain clamp at -1 deg
- Refracted pass prediction: predict_passes_refracted (horizon at -0.569 deg)
- Stellar proper motion: star_observe_pm, star_equatorial_pm (Hipparcos/Gaia convention)
- Light-time correction: planet/sun/small_body_observe_apparent, *_equatorial_apparent
- DE equatorial variants: planet_equatorial_de, moon_equatorial_de

Also includes v0.8.0 orbital_elements type (MPC parser, small_body_observe),
GiST 0-based indexing fix, llms.txt updates, and doc improvements.

All 18 regression suites pass. Zero build warnings (GCC + Clang).
2026-02-21 15:31:46 -07:00
747b7ae60a Fix L1 inclination pruning for HEO orbits, add 66k benchmark
Bug: inner_consistent used sma_low for footprint calculation, but
ground footprint grows with altitude. High-SMA bins (GTO, HEO)
need sma_high to compute the maximum footprint — using sma_low
caused 453 false negatives at high-latitude observers (Tromsoe).

Fix: use sma_high (not sma_low) in L1 inclination pruning.

Added regression test: GTO-debris (inc 5 deg, e=0.73) at Tromsoe
must return identical results from seqscan and index scan.

Benchmark on 65,886-object catalog (full Space-Track including
decayed): 80-92% pruning, zero false negatives across 7 query
patterns. SP-GiST beats seqscan for high-latitude observers.
2026-02-17 23:05:49 -07:00
e1c22cb873 Fix GiST picksplit crash and SP-GiST operator argument order
GiST: entryvec->vector[] uses 1-based indexing (FirstOffsetNumber),
not 0-based. Reading vector[0] hit uninitialized memory, causing
SIGSEGV on large catalogs (14k+ satellites). Fixed in gist_tle_union
and gist_tle_picksplit.

SP-GiST: PostgreSQL requires the indexed column as the LEFT argument
of the operator to form a ScanKey (skey.h:23-26). Flipped &? from
(observer_window, tle) to (tle, observer_window) so inner_consistent
receives scankeys for tree-level pruning.

Removed L0 altitude pruning from inner_consistent — SMA bins don't
carry eccentricity, so HEO satellites (e.g. CLUSTER II, e=0.88,
SMA ~70000 km, perigee ~2000 km) were falsely pruned. L0 now only
narrows SMA range for L1 footprint computation.

All 15 regression tests pass. Consistency check on 14,376 satellites
confirms 0 false negatives, 0 false positives.
2026-02-17 21:30:28 -07:00
2a7240e739 Harden SP-GiST orbital trie after Hamilton review
Fix M2: clamp picksplit nBins to nTuples to prevent out-of-bounds
read on the entries array when called with a single tuple.

Fix H2: use WGS72_AE as effective bin_low for the first bin in L0
inner_consistent, preventing false negatives when objects with lower
SMA than the first label are inserted after index creation.

Fix H3: reject degenerate TLEs (mean_motion <= 0) early in the
visibility filter rather than propagating nonsensical values.

Fix L1: extract shared tle_passes_visibility_filter() to eliminate
duplicated 3-stage filter logic between leaf_consistent and the
standalone tle_visibility_possible operator.

Add boundary tests: degenerate TLE, polar observer, zero-duration
window, and post-insert index-vs-seqscan consistency check.
2026-02-17 20:36:47 -07:00
6f9be428f9 Add SP-GiST orbital trie index for satellite pass prediction (v0.7.0)
2-level SP-GiST index on TLE data: SMA at L0, inclination at L1, with
query-time RAAN filter via J2 secular precession.  New &? operator
(observer_window &? tle) returns true when a satellite might be visible
from a ground observer during a time window.

Index prunes by altitude band, inclination+footprint vs observer
latitude, and RAAN alignment against local sidereal time.  Operator
class tle_spgist_ops is opt-in (not default), coexists with existing
GiST tle_ops.  Equal-population picksplit with sqrt(n) bins.
2026-02-17 20:27:54 -07:00
adfb6949e1 Add range rate fitting, weighted observations, and Gauss angles-only IOD (v0.6.0)
Range rate: topocentric residuals now include an optional 4th component
(dot(Δr, v_ecef) / |Δr|) with OD_RR_SCALE=10.0 for unit balancing.
Controlled via fit_range_rate parameter on tle_from_topocentric().

Weighted observations: per-observation weights applied as √w scaling
to both residuals and Jacobian rows, producing the weighted normal
equations H'WH without explicit W construction. Weights parameter
added to tle_from_eci, tle_from_topocentric, and tle_from_angles.

Gauss angles-only IOD: Vallado Algorithm 52 implementation for
seed-free orbit recovery from 3+ RA/Dec observations. New RA/Dec
residual function with cos(dec) scaling and wrap-around handling.
New tle_from_angles() and tle_from_angles_multi() SQL functions
accepting RA in hours [0,24), Dec in degrees [-90,90].

New standalone test suite: test_od_gauss (17 assertions).
New regression tests: Tests 18-25 covering range rate, weights,
angles-only with/without seed, and error cases.
2026-02-17 17:48:13 -07:00
bca8b3e7eb Add covariance output and condition number to OD solver (v0.5.0)
Computes formal covariance (H^T·H)^{-1} via LAPACK dpotrf_/dpotri_
after DC convergence. Returns upper-triangle array (21 elements for
6-state, 28 for 7-state with B*), condition number from SVD, and
nstate count. Covariance is computed even for perfect-seed fits.

Bumps extension to v0.5.0 with full install SQL and migration path.
2026-02-17 16:15:44 -07:00
6e57071970 Add Gibbs IOD bootstrap for seed-free orbit determination
Eliminates the seed TLE requirement for topocentric fitting by
computing an initial orbit estimate from 3 well-spaced observations
using the Gibbs method. ECI fitting retains the single-observation
r,v approach (exact for two-body) with Gibbs as fallback.
2026-02-17 16:06:05 -07:00
59fd8ba743 Add multi-observer support for topocentric fitting
Extend od_observation_t with observer_idx so each observation can
reference a different ground station. Config now holds an array of
observers instead of a single pointer. The existing single-observer
tle_from_topocentric() is unchanged (sets observer_idx=0 for all obs).

New overload: tle_from_topocentric(topo[], ts[], observer[], int4[], ...)
accepts parallel observer_ids array indexing into the observers array.
PG function overloading resolves by argument types.

Tests 9-11: two-station fit converges, single-station via multi-observer
API matches, out-of-range observer_id raises error.
2026-02-17 15:59:11 -07:00
9b0634725b Add adaptive step limiting to DC solver
Scale step limits by a trust-region factor that halves on divergence
(RMS increases > 1%) and relaxes toward full step on good convergence
(RMS decreases > 10%). Prevents oscillation with poor initial guesses
without affecting well-seeded fits. Also stores SVD condition number
for diagnostic use in upcoming covariance output.

Existing 8 OD regression tests + 67 standalone math tests unaffected
(adaptive_factor starts at 1.0, round-trip tests never trigger
divergence).
2026-02-17 15:55:20 -07:00
87ab81e7d0 Add observation-to-TLE fitting (orbit determination) for v0.4.0
Batch weighted least-squares differential correction using equinoctial
elements, LAPACK dgelss_() for SVD solve, vendored SGP4/SDP4 as the
propagation engine. Per Vallado & Crawford (2008) AIAA 2008-6770.

New SQL functions:
  - tle_from_eci(): fit TLE from ECI position/velocity ephemeris
  - tle_from_topocentric(): fit TLE from az/el/range observations
  - tle_fit_residuals(): per-observation position residuals diagnostic

Solver features: 6-state (orbital) or 7-state (+ B*) fitting,
equinoctial elements for singularity-free optimization, tiered step
limiting, Brouwer/Kozai Newton-Raphson conversion, auto initial guess
from first ECI observation when no seed TLE provided.

Tested: 8 regression tests (LEO/MEO/near-circular round-trips,
B* recovery, topocentric, seedless, error handling, diagnostics),
67 standalone math unit tests, all 14 suites pass.
2026-02-17 15:44:48 -07:00
b18cded4c2 Add PG version test matrix (14-18)
Shell script drives the Dockerfile builder stage across PG versions,
capturing pass/fail + timing per version. Makefile targets: test-matrix,
test-pg%, test-matrix-clean. Also runs standalone DE reader test in the
builder stage to catch compiler-version regressions.

Fix pork chop grid test: add ORDER BY to CROSS JOIN (optimizer chooses
different join nesting across PG versions, reordering rows).
2026-02-17 14:53:32 -07:00
3915d1784f Rename pg_orbit to pg_orrery
An existing product called PG Orbit (a mobile PostgreSQL client)
creates a naming conflict. pg_orrery — a database orrery built from
Keplerian parameters and SQL instead of brass gears.

Build system: control file, Makefile, Dockerfile, docker init script.
C source: GUC prefix, PG_FUNCTION_INFO_V1 symbol, header guards,
ereport prefixes, comments across ~30 files including vendored SGP4.
SQL: all 5 install/migration scripts, function name pg_orrery_ephemeris_info.
Tests: 9 SQL suites, 8 expected outputs, standalone DE reader test.
Documentation: CLAUDE.md, README.md, DESIGN.md, Starlight site infra,
36 MDX pages, OG renderer, logo SVG, docker-compose, agent threads.

All 13 regression suites pass. Docs site builds (37 pages).
2026-02-17 13:36:22 -07:00
cd338c3c64 Vendor SGP4/SDP4 as pure C, drop g++ dependency
Replace the sat_code git submodule (lib/sat_code/) with vendored
sources in src/sgp4/. The upstream .cpp files are renamed to .c —
the code is valid C99 with zero C++ features. This eliminates the
g++ and -lstdc++ build dependencies.

Adds 518 Vallado test vectors (AIAA 2006-6753-Rev1) as a 13th
regression suite to verify byte-identical numerical output.

Updates all documentation (CLAUDE.md, DESIGN.md, 11 MDX pages,
Dockerfile) to reflect the new layout and pure-C compilation.
2026-02-17 12:45:39 -07:00
bcf460a6ff Add standalone DE reader unit test and aliasing safety docs
Standalone test (test/test_de_reader.c): generates a synthetic 12KB
DE binary with known Chebyshev coefficients, exercises header parsing,
polynomial evaluation at 5 domain points, Earth derivation via
center=99, Moon geocentric, layout validation, range/body error paths,
NAN sentinel, and garbage file rejection. 55 tests, no PG dependency.

Makefile: add `make test-de-reader` target; also includes the earlier
sat_code C++ -> pure C sgp4 migration that was missed from prior commits.

astro_math.h: document that ecliptic_to_equatorial/equatorial_to_ecliptic
are NOT safe for aliased (in-place) calls — equ[1] is written before
equ[2] reads ecl[1]. The vendored sgp4/sdp4.c has separate in-place
versions using a temp variable.
2026-02-17 01:56:44 -07:00
15fa553c0e Add optional JPL DE440/441 ephemeris support (v0.3.0)
Clean-room DE binary reader (~400 lines C) with Chebyshev/Clenshaw
evaluation — no GPL dependency on jpl_eph. Per-backend lazy
initialization preserves PARALLEL SAFE. Existing VSOP87/ELP82B
functions stay IMMUTABLE; new _de() variants are STABLE with
automatic fallback to compiled-in ephemerides on any DE failure.

Implementation:
- de_reader.c: header parse, record seek, Clenshaw recurrence
- eph_provider.c: GUC (pg_orbit.ephemeris_path), lazy init,
  ICRS-to-ecliptic frame rotation, on_proc_exit cleanup
- de_funcs.c: 11 new SQL functions (_de variants + diagnostics)
- Constant chain of custody rules 6-8 (frame rotation,
  same-provider, AU consistency)

Extract observe_from_geocentric() to astro_math.h for shared use
by planet_funcs.c, moon_funcs.c, and de_funcs.c.

57 → 68 functions, 11 → 12 regression test suites, all passing.
2026-02-16 19:54:48 -07:00
70420c3b4f Phase 4: Lambert transfer orbit solver for interplanetary trajectories
Add Universal Variable Lambert solver for computing transfer orbits
between any two planets. Enables pork chop plot generation as SQL:

  SELECT dep_date, arr_date, lambert_c3(3, 4, dep_date, arr_date)
  FROM generate_series(...) dep CROSS JOIN generate_series(...) arr;

New functions:
- lambert_transfer(dep_body, arr_body, dep_time, arr_time) → RECORD
  Returns C3 departure/arrival (km^2/s^2), v_infinity (km/s),
  time of flight (days), and transfer orbit SMA (AU).
- lambert_c3(dep_body, arr_body, dep_time, arr_time) → float8
  Convenience: departure C3 only, NULL on solver failure.

The solver uses Stumpff functions for unified elliptic/parabolic/hyperbolic
handling, with Newton-Raphson iteration and bisection fallback.
Each solve is sub-millisecond; PARALLEL SAFE for batch computation.

All 11 regression tests pass.
2026-02-16 02:00:09 -07:00
ad7209d0db Phase 3: Planetary moons and Jupiter radio burst prediction
Add observation functions for 19 planetary moons across four systems:
- Galilean moons (Io, Europa, Ganymede, Callisto) via clean-room L1.2 theory
- Saturn moons (Mimas through Hyperion) via TASS 1.7
- Uranus moons (Miranda through Oberon) via GUST86
- Mars moons (Phobos, Deimos) via MarsSat

Add Jupiter decametric radio burst prediction for Radio JOVE operators:
- io_phase_angle() — Io orbital phase from superior conjunction
- jupiter_cml() — System III Central Meridian Longitude with light-time correction
- jupiter_burst_probability() — Carr et al. (1983) source regions A, B, C, D

L1.2 Galilean theory is a clean-room MIT implementation from the published
IMCCE FORTRAN coefficients. All other ephemeris libraries are MIT-licensed
extractions from Stellarium with static caching removed for PARALLEL SAFE.

All 10 regression tests pass. Extension .so grows from 2.4MB to 2.5MB.
2026-02-16 01:55:13 -07:00
0544a78276 pg_orbit 0.2.0: Full solar system computation at the SQL layer
Phase 1 — Stars, comets, Keplerian propagation:
- star_observe() / star_observe_safe(): fixed star alt/az via IAU 1976
  precession, equatorial-to-horizontal transform
- kepler_propagate(): two-body Keplerian orbit propagation for
  elliptic, parabolic, and hyperbolic orbits
- comet_observe(): observe comets/asteroids from orbital elements
- heliocentric type: ecliptic J2000 position (x, y, z in AU)

Phase 2 — VSOP87 planets, ELP82B Moon, Sun:
- planet_heliocentric(): VSOP87 heliocentric ecliptic J2000 positions
  for Mercury through Neptune (Bretagnon & Francou, MIT)
- planet_observe(): full observation pipeline for any planet
- sun_observe(): Sun position from negated Earth VSOP87
- moon_observe(): ELP2000-82B lunar position (Chapront-Touzé, MIT)
- Clean-room precession (IAU 2006) and sidereal time (IERS 2010)
- elliptic_to_rectangular utility (Stellarium, MIT)

All Stellarium extractions are MIT-licensed, thread-safe (static
caching removed for PARALLEL SAFE), zero external data files.

All 9 regression tests pass (90ms total).
2026-02-16 01:36:27 -07:00
a792e7e083 Extend GiST index to 2-D: altitude + inclination
The 1-D altitude-band index only pruned ~25% of the 22k satellite
catalog (eliminates MEO/GEO/HEO but 75% is LEO).  Adding inclination
as a second indexed dimension prunes an additional ~40% of remaining
candidates — objects in equatorial or low-inclination orbits that
geometrically cannot pass over the observer's latitude.

Key changes:
- tle_alt_range (16 bytes) → tle_orbital_key (32 bytes) with
  inc_low/inc_high fields
- All 8 GiST support functions updated for 2-D bounding boxes
- Penalty uses margin (half-perimeter) not area to avoid degeneracy
  when leaf entries have zero-width inclination ranges
- Picksplit selects split dimension by normalized spread
- && operator now checks altitude AND inclination overlap
- <-> operator remains altitude-only (conjunction screening is
  altitude-dominant)
- SQL operator comments updated for 2-D semantics
- Test adds Equatorial-LEO satellite at ISS altitude but 5° inclination
  to validate inclination-based pruning
2026-02-15 18:10:19 -07:00
5552bf3280 Add convenience functions for Craft integration
Implements 5 new C functions requested by the Craft (Astrolock) API team:
- tle_from_lines(text, text): two-argument TLE constructor
- observer_from_geodetic(float8, float8, float8): numeric observer constructor
- observe(tle, observer, timestamptz): single-call propagate + topocentric
- sgp4_propagate_safe(tle, timestamptz): returns NULL on propagation error
- observe_safe(tle, observer, timestamptz): returns NULL on propagation error

Refactors do_propagate() into safe/unsafe variants to support NULL returns.
Adds regression test (convenience.sql) covering all new functions including
an equivalence test verifying observe() matches the manual two-step pipeline.
All 6 regression tests pass.
2026-02-15 17:44:41 -07:00
15a830dc40 Initial implementation of pg_orbit PostgreSQL extension
6 custom types (tle, eci_position, geodetic, topocentric, observer,
pass_event), 67 SQL functions, 2 operators (&&, <->), and a GiST
operator class for altitude-band indexing. Wraps Bill Gray's sat_code
for SGP4/SDP4 propagation with WGS-72 constants for propagation and
WGS-84 for coordinate output. All 5 regression tests pass on PG 18.
2026-02-15 17:07:07 -07:00