diff --git a/CLAUDE.md b/CLAUDE.md index c648060..6f0fea3 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -1,9 +1,9 @@ # pg_orrery — A Database Orrery for PostgreSQL ## What This Is -A database orrery — celestial mechanics types and functions for PostgreSQL. Native C extension using PGXS, 114 SQL functions, 9 custom types, covering satellites (SGP4/SDP4), planets (VSOP87 + optional JPL DE441), Moon (ELP2000-82B), 19 planetary moons (L1.2/TASS17/GUST86/MarsSat), stars (with proper motion and annual parallax), comets, asteroids (MPC catalog), Jupiter radio bursts, interplanetary Lambert transfers, equatorial RA/Dec coordinates with angular separation, atmospheric refraction, annual stellar aberration, and light-time correction. +A database orrery — celestial mechanics types and functions for PostgreSQL. Native C extension using PGXS, 132 SQL objects (124 user-visible functions + 8 GiST support), 9 custom types, covering satellites (SGP4/SDP4), planets (VSOP87 + optional JPL DE441), Moon (ELP2000-82B), 19 planetary moons (L1.2/TASS17/GUST86/MarsSat), stars (with proper motion and annual parallax), comets, asteroids (MPC catalog), Jupiter radio bursts, interplanetary Lambert transfers, equatorial RA/Dec coordinates with GiST-indexed angular separation, atmospheric refraction, annual stellar aberration, and light-time correction. -**Current version:** 0.10.0 +**Current version:** 0.12.0 **Repository:** https://git.supported.systems/warehack.ing/pg_orrery **Documentation:** https://pg-orrery.warehack.ing @@ -11,7 +11,7 @@ A database orrery — celestial mechanics types and functions for PostgreSQL. Na ```bash make PG_CONFIG=/usr/bin/pg_config # Compile with PGXS sudo make install PG_CONFIG=/usr/bin/pg_config # Install extension -make installcheck PG_CONFIG=/usr/bin/pg_config # Run 19 regression test suites +make installcheck PG_CONFIG=/usr/bin/pg_config # Run 22 regression test suites ``` Requires: PostgreSQL 17 development headers, GCC, Make. @@ -27,7 +27,7 @@ Image: `git.supported.systems/warehack.ing/pg_orrery:pg17` ## Project Layout ``` -pg_orrery.control # Extension metadata (version 0.9.0) +pg_orrery.control # Extension metadata (version 0.12.0) Makefile # PGXS build + Docker targets sql/ pg_orrery--0.1.0.sql # v0.1.0: satellite types/functions/operators @@ -39,6 +39,9 @@ sql/ pg_orrery--0.7.0.sql # v0.7.0: GiST improvements pg_orrery--0.8.0.sql # v0.8.0: orbital_elements type + MPC parser (82 functions) pg_orrery--0.9.0.sql # v0.9.0: equatorial type, refraction, proper motion, light-time (106 functions) + pg_orrery--0.10.0.sql # v0.10.0: angular separation, cone search, apparent functions (114 functions) + pg_orrery--0.11.0.sql # v0.11.0: orbital_elements constructors, moon equatorial (120 functions) + pg_orrery--0.12.0.sql # v0.12.0: equatorial GiST, DE moon equatorial (132 objects) pg_orrery--0.1.0--0.2.0.sql # Migration: v0.1.0 → v0.2.0 (adds solar system) pg_orrery--0.2.0--0.3.0.sql # Migration: v0.2.0 → v0.3.0 (adds DE ephemeris) pg_orrery--0.3.0--0.4.0.sql # Migration: v0.3.0 → v0.4.0 @@ -47,6 +50,9 @@ sql/ pg_orrery--0.6.0--0.7.0.sql # Migration: v0.6.0 → v0.7.0 pg_orrery--0.7.0--0.8.0.sql # Migration: v0.7.0 → v0.8.0 (orbital_elements type) pg_orrery--0.8.0--0.9.0.sql # Migration: v0.8.0 → v0.9.0 (equatorial, refraction, proper motion, light-time) + pg_orrery--0.9.0--0.10.0.sql # Migration: v0.9.0 → v0.10.0 (angular separation, cone search) + pg_orrery--0.10.0--0.11.0.sql # Migration: v0.10.0 → v0.11.0 (constructors, moon equatorial) + pg_orrery--0.11.0--0.12.0.sql # Migration: v0.11.0 → v0.12.0 (equatorial GiST, DE moon equatorial) src/ pg_orrery.c # PG_MODULE_MAGIC + _PG_init() (GUC registration) types.h # All struct definitions + constants + DE body ID mapping @@ -58,7 +64,8 @@ src/ sgp4_funcs.c # sgp4_propagate(), _safe(), _series(), tle_distance() coord_funcs.c # eci_to_geodetic(), eci_to_topocentric(), ground_track() pass_funcs.c # next_pass(), predict_passes(), predict_passes_refracted(), pass_visible() - gist_tle.c # GiST operator class (&&, <->) + gist_tle.c # GiST operator class for TLE (&&, <->) + gist_equatorial.c # GiST operator class for equatorial (KNN <->) # --- Solar System (v0.2.0) --- vsop87.c / vsop87.h # VSOP87 planetary ephemeris (Bretagnon 1988) elp82b.c / elp82b.h # ELP2000-82B lunar ephemeris (Chapront 1988) @@ -83,7 +90,7 @@ src/ # --- JPL DE Ephemeris (v0.3.0) --- de_reader.h / de_reader.c # Clean-room JPL DE binary reader (Chebyshev/Clenshaw) eph_provider.h / eph_provider.c # Provider dispatch, GUC, lazy init, frame rotation - de_funcs.c # All _de() SQL function implementations + de_funcs.c # All _de() SQL function implementations (incl. moon equatorial DE) sgp4/ # Vendored SGP4/SDP4 (Bill Gray's sat_code, MIT license) sgp4.c # Near-earth propagator (period < 225 min) sdp4.c # Deep-space propagator (period >= 225 min) @@ -96,7 +103,7 @@ src/ PROVENANCE.md # Vendoring decision, modifications, verification LICENSE # MIT license (Bill Gray / Project Pluto) test/ - sql/ # 16 regression test suites + sql/ # 22 regression test suites expected/ # Expected output data/vallado_518.json # 518 Vallado test vectors (AIAA 2006-6753-Rev1) docs/ @@ -123,22 +130,23 @@ All types are fixed-size, `STORAGE = plain`, `ALIGNMENT = double`. No TOAST over | `orbital_elements` | 72 | Classical Keplerian elements for comets/asteroids (epoch, q, e, inc, omega, Omega, tp, H, G) | | `equatorial` | 24 | Apparent RA (hours), Dec (degrees), distance (km) — of date | -## Function Domains (114 total) +## Function Domains (132 SQL objects) | Domain | Theory | Key Functions | Count | |--------|--------|---------------|-------| | Satellite | SGP4/SDP4 (Brouwer 1959) | `observe()`, `predict_passes()`, `eci_to_equatorial()` | 25 | | Planets | VSOP87 (Bretagnon 1988) | `planet_observe()`, `planet_equatorial()`, `planet_observe_apparent()` | 7 | | Sun/Moon | VSOP87 + ELP2000-82B | `sun_observe()`, `moon_observe()`, `sun/moon_equatorial()` | 6 | -| Planetary moons | L1.2, TASS17, GUST86, MarsSat | `galilean_observe()`, `saturn_moon_observe()` | 4 | +| Planetary moons | L1.2, TASS17, GUST86, MarsSat | `galilean_observe()`, `saturn_moon_observe()`, `*_equatorial()` | 12 | | Stars | J2000 + IAU 1976 precession | `star_observe()`, `star_equatorial()`, `star_observe_pm()` | 5 | | Comets/asteroids | Two-body Keplerian + MPC | `small_body_observe()`, `small_body_equatorial()`, `oe_from_mpc()` | 19 | | Refraction | Bennett (1982) | `atmospheric_refraction()`, `predict_passes_refracted()` | 4 | | Equatorial spatial | Vincenty formula | `eq_angular_distance()`, `eq_within_cone()`, `<->` | 2 | | Jupiter radio | Carr et al. (1983) | `jupiter_burst_probability()` | 3 | | Transfers | Lambert (Izzo 2015) | `lambert_transfer()`, `lambert_c3()` | 2 | -| DE ephemeris | JPL DE440/441 (optional) | `planet_observe_de()`, `planet_equatorial_de()`, `*_apparent_de()` | 19 | -| GiST index | Altitude-band approximation | `&&` (overlap), `<->` (distance) | 8 | +| DE ephemeris | JPL DE440/441 (optional) | `planet_observe_de()`, `*_equatorial_de()`, `*_apparent_de()` | 23 | +| GiST index (TLE) | Altitude-band approximation | `&&` (overlap), `<->` (distance) | 8 | +| GiST index (equatorial) | Spherical bounding box | `<->` (KNN ordering) | 8 | | Diagnostics | -- | `pg_orrery_ephemeris_info()` | 1 | All functions are `PARALLEL SAFE`. VSOP87/ELP82B functions are `IMMUTABLE` (compiled-in coefficients). DE functions are `STABLE` (external file dependency). @@ -241,6 +249,10 @@ Every `_de()` function mirrors an existing VSOP87 function: | `mars_moon_observe_de()` | `mars_moon_observe()` | STABLE | | `planet_equatorial_de()` | `planet_equatorial()` | STABLE | | `moon_equatorial_de()` | `moon_equatorial()` | STABLE | +| `galilean_equatorial_de()` | `galilean_equatorial()` | STABLE | +| `saturn_moon_equatorial_de()` | `saturn_moon_equatorial()` | STABLE | +| `uranus_moon_equatorial_de()` | `uranus_moon_equatorial()` | STABLE | +| `mars_moon_equatorial_de()` | `mars_moon_equatorial()` | STABLE | | `planet_observe_apparent_de()` | `planet_observe_apparent()` | STABLE | | `sun_observe_apparent_de()` | `sun_observe_apparent()` | STABLE | | `moon_observe_apparent_de()` | `moon_observe_apparent()` | STABLE | @@ -268,7 +280,7 @@ All numerical logic is byte-identical to upstream. Verified against 518 Vallado ## Testing -19 regression test suites via `make installcheck`: +22 regression test suites via `make installcheck`: | Suite | What it tests | |-------|--------------| @@ -276,7 +288,7 @@ All numerical logic is byte-identical to upstream. Verified against 518 Vallado | sgp4_propagate | SGP4/SDP4, propagation series, tle_distance | | coord_transforms | TEME-to-geodetic, TEME-to-topocentric, ground_track | | pass_prediction | predict_passes, next_pass, pass_visible, min elevation filter | -| gist_index | `&&` overlap, `<->` distance, GiST index scan, KNN ordering | +| gist_index | `&&` overlap, `<->` distance, GiST index scan, KNN ordering (TLE) | | convenience | observe(), observe_safe(), tle_from_lines(), observer_from_geodetic() | | star_observe | Star observation, IAU 1976 precession, heliocentric type I/O | | kepler_comet | Keplerian propagation (elliptic/parabolic/hyperbolic), comet_observe | @@ -291,10 +303,13 @@ All numerical logic is byte-identical to upstream. Verified against 518 Vallado | refraction | Bennett refraction, P/T correction, apparent elevation, refracted pass prediction | | aberration | Annual aberration magnitude, DE apparent fallback, angular distance, cone search, stellar parallax | | vallado_518 | 518 Vallado test vectors (AIAA 2006-6753-Rev1), per-satellite breakdown | +| v011_features | make_orbital_elements constructors, moon equatorial functions | +| gist_equatorial | Equatorial GiST KNN ordering, RA wrapping, cone search, EXPLAIN index scan | +| v012_features | DE moon equatorial fallback to VSOP87, invalid body_id rejection | ### PG Version Matrix -Test all 19 regression suites + DE reader unit test across PostgreSQL 14-18 using Docker: +Test all 22 regression suites + DE reader unit test across PostgreSQL 14-18 using Docker: ```bash make test-matrix # Full matrix (PG 14-18) @@ -318,9 +333,9 @@ Logs saved to `test/matrix-logs/pg${ver}.log`. The script reuses the Dockerfile **Live:** https://pg-orrery.warehack.ing -Starlight docs at `docs/` — 44 MDX pages covering all domains. +Starlight docs at `docs/` — 44+ MDX pages covering all domains. -Sections: Getting Started, Guides (9 domain walkthroughs incl. DE ephemeris), Workflow Translation (Skyfield/Horizons/GMAT/Radio Jupiter Pro comparisons), Reference (all 106 functions incl. DE variants, equatorial, refraction), Architecture (Hamilton's principles, constant custody, observation pipeline), Performance (benchmarks). +Sections: Getting Started, Guides (9 domain walkthroughs incl. DE ephemeris), Workflow Translation (Skyfield/Horizons/GMAT/Radio Jupiter Pro comparisons), Reference (all 132 SQL objects incl. DE variants, equatorial GiST, refraction), Architecture (Hamilton's principles, constant custody, observation pipeline), Performance (benchmarks). ### Local Development ```bash diff --git a/docs/src/content/docs/reference/functions-de.mdx b/docs/src/content/docs/reference/functions-de.mdx index 5fb2d56..eb880a3 100644 --- a/docs/src/content/docs/reference/functions-de.mdx +++ b/docs/src/content/docs/reference/functions-de.mdx @@ -368,6 +368,122 @@ FROM moon_equatorial_de(now()) AS e; --- +## galilean_equatorial_de + +Computes the geocentric equatorial coordinates (RA/Dec) of a Galilean moon using JPL DE ephemeris for Jupiter's position. Falls back to VSOP87 for both Jupiter and Earth when DE is unavailable. Moon offsets always come from Lieske L1.2 theory. + +### Signature + +```sql +galilean_equatorial_de(moon_id int4, t timestamptz) → equatorial +``` + +### Parameters + +| Parameter | Type | Description | +|-----------|------|-------------| +| `moon_id` | `int4` | 0=Io, 1=Europa, 2=Ganymede, 3=Callisto | +| `t` | `timestamptz` | Evaluation time | + +### Returns + +An `equatorial` with RA (hours), Dec (degrees), and distance (km) from Earth's center. + +### Example + +```sql +-- All 4 Galilean moons via DE +SELECT moon_id, + round(eq_ra(galilean_equatorial_de(moon_id, now()))::numeric, 4) AS ra, + round(eq_dec(galilean_equatorial_de(moon_id, now()))::numeric, 4) AS dec +FROM generate_series(0, 3) AS moon_id; +``` + +--- + +## saturn_moon_equatorial_de + +Computes the geocentric equatorial coordinates (RA/Dec) of a Saturn moon using JPL DE ephemeris for Saturn's position. Falls back to VSOP87 when DE is unavailable. Moon offsets come from TASS 1.7 theory. + +### Signature + +```sql +saturn_moon_equatorial_de(moon_id int4, t timestamptz) → equatorial +``` + +### Parameters + +| Parameter | Type | Description | +|-----------|------|-------------| +| `moon_id` | `int4` | 0=Mimas, 1=Enceladus, 2=Tethys, 3=Dione, 4=Rhea, 5=Titan, 6=Iapetus, 7=Hyperion | +| `t` | `timestamptz` | Evaluation time | + +### Example + +```sql +-- Titan's position via DE +SELECT round(eq_ra(saturn_moon_equatorial_de(5, now()))::numeric, 4) AS ra, + round(eq_dec(saturn_moon_equatorial_de(5, now()))::numeric, 4) AS dec; +``` + +--- + +## uranus_moon_equatorial_de + +Computes the geocentric equatorial coordinates (RA/Dec) of a Uranus moon using JPL DE ephemeris for Uranus's position. Falls back to VSOP87 when DE is unavailable. Moon offsets come from GUST86 theory. + +### Signature + +```sql +uranus_moon_equatorial_de(moon_id int4, t timestamptz) → equatorial +``` + +### Parameters + +| Parameter | Type | Description | +|-----------|------|-------------| +| `moon_id` | `int4` | 0=Miranda, 1=Ariel, 2=Umbriel, 3=Titania, 4=Oberon | +| `t` | `timestamptz` | Evaluation time | + +### Example + +```sql +-- Titania's position via DE +SELECT round(eq_ra(uranus_moon_equatorial_de(3, now()))::numeric, 4) AS ra, + round(eq_dec(uranus_moon_equatorial_de(3, now()))::numeric, 4) AS dec; +``` + +--- + +## mars_moon_equatorial_de + +Computes the geocentric equatorial coordinates (RA/Dec) of a Mars moon using JPL DE ephemeris for Mars's position. Falls back to VSOP87 when DE is unavailable. Moon offsets come from MarsSat analytical theory. + +### Signature + +```sql +mars_moon_equatorial_de(moon_id int4, t timestamptz) → equatorial +``` + +### Parameters + +| Parameter | Type | Description | +|-----------|------|-------------| +| `moon_id` | `int4` | 0=Phobos, 1=Deimos | +| `t` | `timestamptz` | Evaluation time | + +### Example + +```sql +-- Both Mars moons via DE +SELECT moon_id, + round(eq_ra(mars_moon_equatorial_de(moon_id, now()))::numeric, 4) AS ra, + round(eq_dec(mars_moon_equatorial_de(moon_id, now()))::numeric, 4) AS dec +FROM generate_series(0, 1) AS moon_id; +``` + +--- + ## pg_orrery_ephemeris_info Returns diagnostic information about the current ephemeris provider. diff --git a/docs/src/content/docs/reference/operators-gist.mdx b/docs/src/content/docs/reference/operators-gist.mdx index 8e8229c..8b26cb6 100644 --- a/docs/src/content/docs/reference/operators-gist.mdx +++ b/docs/src/content/docs/reference/operators-gist.mdx @@ -6,7 +6,7 @@ sidebar: import { Aside, Tabs, TabItem } from "@astrojs/starlight/components"; -pg_orrery defines four operators on the `tle` type and two operator classes (GiST and SP-GiST) that enable indexed satellite queries over large catalogs. The GiST index accelerates conjunction screening (orbit-to-orbit overlap). The SP-GiST index accelerates pass prediction (observer-to-orbit visibility). +pg_orrery defines operators on the `tle` and `equatorial` types with three operator classes (two GiST + one SP-GiST) that enable indexed queries over large catalogs. The TLE GiST index accelerates conjunction screening (orbit-to-orbit overlap). The SP-GiST index accelerates pass prediction (observer-to-orbit visibility). The equatorial GiST index accelerates nearest-neighbor sky queries (angular distance KNN). --- @@ -205,6 +205,81 @@ REINDEX INDEX idx_tle_gist; --- +## GiST Operator Class: eq_gist_ops + +The `eq_gist_ops` operator class enables GiST indexing on `equatorial` columns. With this index, the `<->` operator (angular distance in degrees via the Vincenty formula) supports index-ordered KNN queries — PostgreSQL traverses the tree by increasing angular distance without computing all distances upfront. + + + +### Creating the Index + +```sql +CREATE INDEX idx_sky_eq ON sky_cache USING gist (eq); +``` + +### What Gets Indexed + +The GiST index stores a 24-byte float-precision spherical bounding box for each entry: +- **RA range:** `[ra_low, ra_high]` in radians. When `ra_low > ra_high`, the box wraps across 0h (covers `[ra_low, 2pi) union [0, ra_high]`) +- **Dec range:** `[dec_low, dec_high]` in radians + +Float precision (~0.12 arcsec bounding error at RA = 2pi) is more than sufficient for index pruning. Actual angular distance is computed in double precision via the Vincenty formula during recheck. + +### Index-Accelerated Queries + + + +```sql +-- Find the 10 nearest sky objects to Jupiter +SELECT name, + round((eq <-> planet_equatorial_apparent(5, NOW()))::numeric, 4) AS dist_deg +FROM sky_cache +ORDER BY eq <-> planet_equatorial_apparent(5, NOW()) +LIMIT 10; +``` + + +```sql +-- Everything within 15 degrees of Vega, sorted by distance +SELECT name, + round((eq <-> star_equatorial(18.616, 38.784, NOW()))::numeric, 2) AS dist_deg +FROM sky_cache +WHERE eq_within_cone(eq, star_equatorial(18.616, 38.784, NOW()), 15.0) +ORDER BY eq <-> star_equatorial(18.616, 38.784, NOW()); +``` + + +```sql +-- Confirm the planner uses the GiST index +SET enable_seqscan = off; +EXPLAIN (COSTS OFF) +SELECT name FROM sky_cache +ORDER BY eq <-> '(12.00000000,0.00000000,0.000)'::equatorial +LIMIT 3; +-- Should show: Index Scan using idx_sky_eq on sky_cache +RESET enable_seqscan; +``` + + + +### RA Wrapping + +Objects near 0h (e.g., in Pisces/Aquarius) and objects near 24h are correctly identified as neighbors. The bounding box merge and distance functions handle the RA discontinuity at 0h/24h explicitly. An object at RA = 23.9h and another at RA = 0.1h are approximately 3 degrees apart (at moderate declination), and the KNN traversal finds them as neighbors. + +### Polar Regions + +Near the celestial poles (Dec approaching +/-90 degrees), RA becomes degenerate — a small patch of sky spans a wide RA range. Bounding boxes for polar objects may cover the full RA circle. This does not affect **correctness** (the Vincenty formula handles pole convergence naturally) but can degrade **index selectivity** for dense polar catalogs. For typical sky catalogs (<10,000 objects), the effect is negligible. + +### Design Notes + +- **KNN only** (strategy 15, `<->` ordering). No `&&` overlap operator — meaningless for point types. +- **Distance unit:** degrees, matching `eq_angular_distance()`. +- **Lower-bound contract hardened:** box boundaries widened by float epsilon before distance computation to guarantee KNN correctness under float-to-double promotion. + +--- + ### &? (Visibility Cone) Tests whether a satellite could possibly be visible from a ground observer during a time window. This is a geometric superset filter -- it may include satellites that do not produce an actual pass (false positives), but will never exclude one that does (no false negatives).