Fix function names and argument order in catalog guide

Three-tier discovery: pg-orrery-catalog in PATH, sibling dev
checkout, or original build_catalog.py + curl. Indexes use
IF NOT EXISTS for idempotent re-runs.

Makefile

@@ -4,7 +4,8 @@ DATA = sql/pg_orrery--0.1.0.sql sql/pg_orrery--0.2.0.sql sql/pg_orrery--0.1.0--0
        sql/pg_orrery--0.3.0.sql sql/pg_orrery--0.2.0--0.3.0.sql \
        sql/pg_orrery--0.4.0.sql sql/pg_orrery--0.3.0--0.4.0.sql \
        sql/pg_orrery--0.5.0.sql sql/pg_orrery--0.4.0--0.5.0.sql \
-       sql/pg_orrery--0.6.0.sql sql/pg_orrery--0.5.0--0.6.0.sql
+       sql/pg_orrery--0.6.0.sql sql/pg_orrery--0.5.0--0.6.0.sql \
+       sql/pg_orrery--0.7.0.sql sql/pg_orrery--0.6.0--0.7.0.sql
 
 # Our extension C sources
 OBJS = src/pg_orrery.o src/tle_type.o src/eci_type.o src/observer_type.o \
@@ -16,7 +17,8 @@ OBJS = src/pg_orrery.o src/tle_type.o src/eci_type.o src/observer_type.o \
        src/moon_funcs.o src/radio_funcs.o \
        src/lambert.o src/transfer_funcs.o \
        src/de_reader.o src/eph_provider.o src/de_funcs.o \
-       src/od_math.o src/od_iod.o src/od_solver.o src/od_funcs.o
+       src/od_math.o src/od_iod.o src/od_solver.o src/od_funcs.o \
+       src/spgist_tle.o
 
 # Vendored SGP4/SDP4 sources (pure C, from Bill Gray's sat_code, MIT license)
 SGP4_DIR = src/sgp4
@@ -31,7 +33,7 @@ OBJS += $(SGP4_OBJS)
 # Regression tests
 REGRESS = tle_parse sgp4_propagate coord_transforms pass_prediction gist_index convenience \
          star_observe kepler_comet planet_observe moon_observe lambert_transfer \
-         de_ephemeris od_fit vallado_518
+         de_ephemeris od_fit spgist_tle vallado_518
 REGRESS_OPTS = --inputdir=test
 
 # Pure C — no C++ runtime needed. LAPACK for OD solver (dgelss_).

bench/benchmark.sql

@@ -0,0 +1,234 @@
+-- ============================================================
+-- SP-GiST Orbital Trie Benchmark (Phase 3)
+-- CelesTrak active catalog, ~14k satellites
+-- ============================================================
+
+\timing on
+
+-- ============================================================
+-- 1. Catalog distribution analysis
+-- ============================================================
+SELECT
+    CASE
+        WHEN tle_perigee(tle) < 2000 THEN 'LEO (<2000km)'
+        WHEN tle_perigee(tle) < 20000 THEN 'MEO (2000-20000km)'
+        WHEN tle_perigee(tle) < 34000 THEN 'GEO-transfer'
+        ELSE 'GEO/HEO (>34000km)'
+    END AS regime,
+    count(*) AS n,
+    round(100.0 * count(*) / (SELECT count(*) FROM bench_catalog), 1) AS pct
+FROM bench_catalog
+GROUP BY 1
+ORDER BY 2 DESC;
+
+-- ============================================================
+-- 2. Create indexes
+-- ============================================================
+\echo '--- CREATE SP-GiST INDEX ---'
+CREATE INDEX bench_spgist ON bench_catalog USING spgist (tle tle_spgist_ops);
+
+\echo '--- CREATE GiST INDEX ---'
+CREATE INDEX bench_gist ON bench_catalog USING gist (tle tle_ops);
+
+-- Index sizes
+SELECT indexname,
+       pg_size_pretty(pg_relation_size(indexname::regclass)) AS size
+FROM pg_indexes
+WHERE tablename = 'bench_catalog'
+ORDER BY indexname;
+
+-- ============================================================
+-- 3. Benchmark: 2h window, Eagle Idaho (43.7N) — RAAN active
+-- ============================================================
+\echo '--- BENCHMARK 1: 2h window, Eagle Idaho, 10 deg min_el ---'
+
+-- 3a. Sequential scan (baseline)
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT count(*) AS seqscan_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+-- 3b. SP-GiST index scan
+SET enable_seqscan = off;
+SELECT count(*) AS spgist_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 4. Benchmark: 24h window, Eagle Idaho — RAAN bypassed
+-- ============================================================
+\echo '--- BENCHMARK 2: 24h window, Eagle Idaho, 10 deg min_el ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT count(*) AS seqscan_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 00:00:00+00'::timestamptz,
+    '2026-02-18 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+SELECT count(*) AS spgist_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 00:00:00+00'::timestamptz,
+    '2026-02-18 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 5. Benchmark: 2h window, Equatorial observer — all inc pass
+-- ============================================================
+\echo '--- BENCHMARK 3: 2h window, Equator, 10 deg min_el ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT count(*) AS seqscan_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+SELECT count(*) AS spgist_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 6. Benchmark: 45 deg min_el (aggressive altitude filter)
+-- ============================================================
+\echo '--- BENCHMARK 4: 2h window, Eagle Idaho, 45 deg min_el ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT count(*) AS seqscan_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    45.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+SELECT count(*) AS spgist_candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    45.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 7. Consistency check: index and seqscan must agree
+-- ============================================================
+\echo '--- CONSISTENCY CHECK ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+CREATE TEMPORARY TABLE seq_results AS
+SELECT norad_id FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+CREATE TEMPORARY TABLE idx_results AS
+SELECT norad_id FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- These should both be 0
+SELECT count(*) AS in_seq_not_idx FROM seq_results
+WHERE norad_id NOT IN (SELECT norad_id FROM idx_results);
+
+SELECT count(*) AS in_idx_not_seq FROM idx_results
+WHERE norad_id NOT IN (SELECT norad_id FROM seq_results);
+
+DROP TABLE seq_results, idx_results;
+
+-- ============================================================
+-- 8. EXPLAIN ANALYZE for query plan details
+-- ============================================================
+\echo '--- EXPLAIN ANALYZE: SP-GiST scan ---'
+
+SET enable_seqscan = off;
+EXPLAIN ANALYZE
+SELECT count(*)
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+\echo '--- EXPLAIN ANALYZE: Sequential scan ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+EXPLAIN ANALYZE
+SELECT count(*)
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+-- ============================================================
+-- Cleanup
+-- ============================================================
+DROP TABLE bench_catalog;
+
+\timing off

bench/benchmark_results.txt

@@ -0,0 +1,193 @@
+Timing is on.
+       regime       |   n   | pct  
+--------------------+-------+------
+ LEO (<2000km)      | 13587 | 94.5
+ GEO/HEO (>34000km) |   588 |  4.1
+ MEO (2000-20000km) |   111 |  0.8
+ GEO-transfer       |    90 |  0.6
+(4 rows)
+
+Time: 7.139 ms
+--- CREATE SP-GiST INDEX ---
+CREATE INDEX
+Time: 12.351 ms
+ spgist_size 
+-------------
+ 2424 kB
+(1 row)
+
+Time: 1.388 ms
+--- BENCHMARK 1: 2h window, Eagle Idaho, 10 deg min_el ---
+SET
+Time: 0.035 ms
+SET
+Time: 0.021 ms
+Sequential scan:
+                                                          QUERY PLAN                                                          
+------------------------------------------------------------------------------------------------------------------------------
+ Aggregate  (cost=470.74..470.75 rows=1 width=8) (actual time=2.706..2.707 rows=1.00 loops=1)
+   Buffers: shared hit=291
+   ->  Seq Scan on bench_catalog  (cost=0.00..470.70 rows=14 width=0) (actual time=0.015..2.649 rows=2261.00 loops=1)
+         Filter: ('("43.6977N 116.3535W 760m","2026-02-16 19:00:00-07","2026-02-16 21:00:00-07",10)'::observer_window &? tle)
+         Rows Removed by Filter: 12115
+         Buffers: shared hit=291
+ Planning:
+   Buffers: shared hit=8
+ Planning Time: 0.093 ms
+ Execution Time: 2.725 ms
+(10 rows)
+
+Time: 3.621 ms
+RESET
+Time: 0.028 ms
+RESET
+Time: 0.009 ms
+SET
+Time: 0.012 ms
+SP-GiST index scan:
+                                                               QUERY PLAN                                                                
+-----------------------------------------------------------------------------------------------------------------------------------------
+ Aggregate  (cost=1754.89..1754.90 rows=1 width=8) (actual time=3.896..3.897 rows=1.00 loops=1)
+   Buffers: shared hit=1161
+   ->  Bitmap Heap Scan on bench_catalog  (cost=1284.16..1754.86 rows=14 width=0) (actual time=0.957..3.833 rows=2261.00 loops=1)
+         Filter: ('("43.6977N 116.3535W 760m","2026-02-16 19:00:00-07","2026-02-16 21:00:00-07",10)'::observer_window &? tle)
+         Rows Removed by Filter: 12115
+         Heap Blocks: exact=291
+         Buffers: shared hit=1161
+         ->  Bitmap Index Scan on bench_spgist  (cost=0.00..1284.16 rows=14376 width=0) (actual time=0.925..0.925 rows=14376.00 loops=1)
+               Index Searches: 1
+               Buffers: shared hit=870
+ Planning Time: 0.050 ms
+ Execution Time: 3.936 ms
+(12 rows)
+
+Time: 4.150 ms
+RESET
+Time: 0.038 ms
+--- BENCHMARK 2: 24h window, Eagle Idaho, 10 deg min_el ---
+SET
+Time: 0.016 ms
+SET
+Time: 0.009 ms
+Sequential scan:
+ seqscan_24h 
+-------------
+       13562
+(1 row)
+
+Time: 2.867 ms
+RESET
+Time: 0.023 ms
+RESET
+Time: 0.008 ms
+SET
+Time: 0.011 ms
+SP-GiST index scan:
+ spgist_24h 
+------------
+      13562
+(1 row)
+
+Time: 3.832 ms
+RESET
+Time: 0.025 ms
+--- BENCHMARK 3: 2h window, Equator, 10 deg min_el ---
+SET
+Time: 0.011 ms
+SET
+Time: 0.008 ms
+Sequential scan:
+ seqscan_equator 
+-----------------
+            2073
+(1 row)
+
+Time: 2.564 ms
+RESET
+Time: 0.011 ms
+RESET
+Time: 0.007 ms
+SET
+Time: 0.008 ms
+SP-GiST index scan:
+ spgist_equator 
+----------------
+           2073
+(1 row)
+
+Time: 3.566 ms
+RESET
+Time: 0.020 ms
+--- BENCHMARK 4: 2h window, Eagle Idaho, 45 deg min_el ---
+SET
+Time: 0.010 ms
+SET
+Time: 0.006 ms
+Sequential scan:
+ seqscan_45deg 
+---------------
+          1407
+(1 row)
+
+Time: 2.510 ms
+RESET
+Time: 0.021 ms
+RESET
+Time: 0.007 ms
+SET
+Time: 0.008 ms
+SP-GiST index scan:
+ spgist_45deg 
+--------------
+         1407
+(1 row)
+
+Time: 3.591 ms
+RESET
+Time: 0.014 ms
+--- CONSISTENCY CHECK ---
+SET
+Time: 0.017 ms
+SET
+Time: 0.009 ms
+SELECT 2261
+Time: 4.578 ms
+RESET
+Time: 0.023 ms
+RESET
+Time: 0.008 ms
+SET
+Time: 0.010 ms
+SELECT 2261
+Time: 3.583 ms
+RESET
+Time: 0.022 ms
+ in_seq_not_idx 
+----------------
+              0
+(1 row)
+
+Time: 0.586 ms
+ in_idx_not_seq 
+----------------
+              0
+(1 row)
+
+Time: 0.280 ms
+DROP TABLE
+Time: 0.951 ms
+--- PRUNING SUMMARY ---
+     scenario     | catalog_size | candidates | candidate_pct | pruning_pct 
+------------------+--------------+------------+---------------+-------------
+ 2h/Eagle/10deg   |        14376 |       2261 |          15.7 |        84.3
+ 2h/Eagle/45deg   |        14376 |       1407 |           9.8 |        90.2
+ 2h/Equator/10deg |        14376 |       2073 |          14.4 |        85.6
+ 24h/Eagle/10deg  |        14376 |      13562 |          94.3 |         5.7
+(4 rows)
+
+Time: 9.346 ms
+DROP INDEX
+Time: 1.176 ms
+DROP TABLE
+Time: 1.725 ms
+Timing is off.

bench/benchmark_results_30k.txt

@@ -0,0 +1,115 @@
+pg_orrery v0.7.0 SP-GiST Benchmark — 30k Space-Track Catalog
+=============================================================
+Date: 2026-02-17
+Catalog: Space-Track USSPACECOM full catalog (29,784 objects)
+Host: Linux 6.16.5-arch1-1, PostgreSQL 17
+Branch: phase/spgist-orbital-trie
+
+Catalog Composition:
+  LEO (<128 min):          25,641  (86.1%)
+  MEO (128-720 min):        1,801  ( 6.0%)
+  GEO/HEO (720-1500 min):  2,253  ( 7.6%)
+  Super-GEO (>1500 min):      89  ( 0.3%)
+
+Index Build:
+  SP-GiST: 46.7 ms, 4,816 kB
+  GiST:    65.8 ms, 5,856 kB
+  Table:              4,760 kB
+
+==============================================================
+TIMING RESULTS (best of 3 runs, ms)
+==============================================================
+
+Query Pattern                | Seqscan | SP-GiST | Delta
+-----------------------------|---------|---------|-------
+2h window, Eagle ID, 10 deg  |  5.19   |  5.22   | +0.03
+6h window, Eagle ID, 10 deg  |  5.34   |  6.34   | +1.00
+24h window, Eagle ID, 10 deg |  5.42   |  6.68   | +1.26
+2h window, Eagle ID, 30 deg  |  5.54   |  5.64   | +0.10
+2h window, equatorial, 10deg |  5.26   |  5.59   | +0.33
+
+==============================================================
+PRUNING RESULTS
+==============================================================
+
+Query Pattern                | Candidates | % Pass | % Pruned
+-----------------------------|------------|--------|---------
+2h window, Eagle ID, 10 deg  |    7,188   |  24.1% |   75.9%
+6h window, Eagle ID, 10 deg  |   12,373   |  41.5% |   58.5%
+24h window, Eagle ID, 10 deg |   26,971   |  90.6% |    9.4%
+2h window, Eagle ID, 30 deg  |    5,232   |  17.6% |   82.4%
+2h window, equatorial, 10deg |    4,670   |  15.7% |   84.3%
+
+==============================================================
+BUFFER I/O COMPARISON (2h Eagle 10deg)
+==============================================================
+
+Method      | Pages Read | Heap Fetches | Scan Type
+------------|------------|--------------|------------------
+Seqscan     |    595     |     n/a      | Seq Scan
+SP-GiST     |  2,396     |       0      | Index Only Scan
+
+SP-GiST reads 4.0x more pages than seqscan (2,396 vs 595).
+But uses Index Only Scan with zero heap fetches.
+
+==============================================================
+CONSISTENCY CHECK
+==============================================================
+
+False negatives: 0  (index never misses a seqscan result)
+False positives: 0  (index never returns extra results)
+Seqscan count:   7,188
+SP-GiST count:   7,188
+
+==============================================================
+SCALING TREND (2h Eagle 10deg, best-of-3)
+==============================================================
+
+Catalog Size | Seqscan | SP-GiST | Delta  | SP-GiST Pages | Seq Pages
+-------------|---------|---------|--------|---------------|----------
+14,376       |  4.5 ms |  6.1 ms | +1.6ms |      888      |    291
+20,597       |  3.8 ms |  4.7 ms | +0.9ms |    (IOOS)     |   (est)
+29,784       |  5.2 ms |  5.2 ms | +0.0ms |    2,396      |    595
+
+The delta is converging toward zero. At 30k the SP-GiST index is
+essentially tied with seqscan on the 2h/10deg query. For queries
+with fewer survivors (30deg elevation, equatorial observer), the
+index is within 0.1-0.3ms.
+
+==============================================================
+PLANNER BEHAVIOR
+==============================================================
+
+PostgreSQL's query planner CHOOSES the SP-GiST index by default
+at 30k (without any enable_seqscan=off forcing). The planner's
+cost model prefers the Index Only Scan.
+
+EXPLAIN output (default settings):
+  Index Only Scan using bench_spgist on bench_catalog
+    Index Cond: (tle &? ...)
+    Heap Fetches: 0
+    Buffers: shared hit=2396
+    Execution Time: 7.246 ms (with planning)
+
+==============================================================
+NOTES
+==============================================================
+
+1. At 30k objects, the planner voluntarily chooses SP-GiST over
+   seqscan. This is the crossover point where the index becomes
+   the planner's preferred strategy.
+
+2. The Index Only Scan with zero heap fetches means the index
+   contains all information needed — no table access required.
+
+3. The 75.9% pruning rate on the 2h window means only 7,188 of
+   29,784 satellites need SGP4 propagation. This avoids ~22,596
+   unnecessary SGP4 calls in the predict_passes() pipeline.
+
+4. The equatorial observer (84.3% pruned) and high-elevation
+   (82.4% pruned) queries show the strongest filtering because
+   the altitude and RAAN filters are most aggressive there.
+
+5. The 24h window only prunes 9.4% because the RAAN filter
+   self-disables for full Earth rotations, leaving only the
+   altitude and inclination filters active.

bench/benchmark_results_66k.txt

@@ -0,0 +1,114 @@
+pg_orrery v0.7.0 SP-GiST Benchmark — 66k Full Space-Track Catalog
+===================================================================
+Date: 2026-02-17
+Catalog: Space-Track USSPACECOM full catalog including decayed (65,886 objects)
+Host: Linux 6.16.5-arch1-1, PostgreSQL 17
+Branch: phase/spgist-orbital-trie
+Note: After fixing L1 inclination pruning (sma_low -> sma_high)
+
+Catalog Composition:
+  LEO (<128 min):          59,537  (90.4%)
+  MEO (128-720 min):        3,474  ( 5.3%)
+  GEO/HEO (720-1500 min):  2,643  ( 4.0%)
+  Super-GEO (>1500 min):     232  ( 0.4%)
+
+Index Build:
+  SP-GiST: 55.2 ms, 11 MB
+  GiST:   118.2 ms, 13 MB
+  Table:              10 MB
+
+==============================================================
+TIMING RESULTS (best of 2-3 runs, ms)
+==============================================================
+
+Query Pattern                | Seqscan | SP-GiST | Delta
+-----------------------------|---------|---------|-------
+2h window, Eagle ID, 10 deg  |  12.5   |  14.0   | +1.5
+6h window, Eagle ID, 10 deg  |  12.2   |  15.6   | +3.4
+2h window, Tromsø, 10 deg    |  11.3   |  10.9   | -0.4  ★
+24h window, Eagle ID, 10 deg |  12.0   |  16.2   | +4.2
+
+★ Tromsø (69.6°N): SP-GiST beats seqscan. High-latitude observers
+  benefit most from inclination pruning.
+
+==============================================================
+PRUNING RESULTS
+==============================================================
+
+Query Pattern                | Candidates | % Pass | % Pruned
+-----------------------------|------------|--------|---------
+2h window, Eagle ID, 10 deg  |   12,964   |  19.7% |   80.3%
+6h window, Eagle ID, 10 deg  |   24,274   |  36.8% |   63.2%
+24h window, Eagle ID, 10 deg |   60,875   |  92.4% |    7.6%
+2h window, Eagle ID, 30 deg  |    9,680   |  14.7% |   85.3%
+2h window, equatorial, 10deg |    9,699   |  14.7% |   85.3%
+2h window, Tromsø 69.6°N     |    6,529   |   9.9% |   90.1%
+2h window, South Pole 85°S   |    5,248   |   8.0% |   92.0%
+
+==============================================================
+CONSISTENCY CHECKS (all patterns)
+==============================================================
+
+Query Pattern                | False Negatives | False Positives
+-----------------------------|-----------------|----------------
+2h Eagle 10deg               |        0        |       0
+6h Eagle 10deg               |        0        |       0
+24h Eagle 10deg              |        0        |       0
+2h Eagle 30deg               |        0        |       0
+2h Equator 10deg             |        0        |       0
+2h Tromsø 10deg              |        0        |       0
+2h South Pole 10deg          |        0        |       0
+
+==============================================================
+SCALING TREND (2h Eagle 10deg, best-of-N)
+==============================================================
+
+Catalog Size | Seqscan | SP-GiST | Delta  | Notes
+-------------|---------|---------|--------|------
+14,376       |  4.5 ms |  6.1 ms | +1.6ms | Active CelesTrak
+29,784       |  5.2 ms |  5.2 ms | +0.0ms | Active Space-Track (before fix)
+65,886       | 12.5 ms | 14.0 ms | +1.5ms | Full catalog incl decayed (after fix)
+
+The fix (sma_high instead of sma_low for footprint) adds ~1-2ms overhead
+by conservatively keeping more subtrees alive during L1 pruning. This is
+the correct trade-off: zero false negatives is non-negotiable.
+
+==============================================================
+PLANNER BEHAVIOR (66k)
+==============================================================
+
+PostgreSQL still chooses SP-GiST Index Only Scan by default:
+  Index Only Scan using bench_spgist on bench_catalog
+    Index Cond: (tle &? ...)
+    Heap Fetches: 0
+    Buffers: shared hit=4990
+
+Seqscan would read 1,297 pages. Index reads 4,990 pages (3.8x more).
+But Index Only Scan avoids all heap I/O.
+
+==============================================================
+KEY FINDING: HIGH-LATITUDE OBSERVERS
+==============================================================
+
+The SP-GiST index is most valuable for high-latitude observers:
+
+  Tromsø (69.6°N): 90.1% pruned, SP-GiST BEATS seqscan by 0.4ms
+  South Pole (85°S): 92.0% pruned
+
+High-latitude locations eliminate most LEO satellites via the
+inclination filter — only satellites with inc > ~60° can reach
+these latitudes. The SP-GiST trie prunes entire inclination
+subtrees at L1, making the index scan faster than touching
+every page in the table.
+
+==============================================================
+WHAT THE 80-92% PRUNING MEANS IN PRACTICE
+==============================================================
+
+For a 65,886-object catalog with a 2-hour window:
+  - Without &? operator: 65,886 SGP4 predict_passes() calls
+  - With &? operator:    12,964 SGP4 calls (Eagle) or 5,248 (South Pole)
+  - Savings: 52,922-60,638 unnecessary propagation calls avoided
+
+At ~1ms per predict_passes() call (7-day window, 30s resolution),
+that's 53-61 seconds of saved computation per query.

bench/benchmark_results_v2.txt

@@ -0,0 +1,202 @@
+Timing is on.
+       regime       |   n   | pct  
+--------------------+-------+------
+ LEO (<2000km)      | 13587 | 94.5
+ GEO/HEO (>34000km) |   588 |  4.1
+ MEO (2000-20000km) |   111 |  0.8
+ GEO-transfer       |    90 |  0.6
+(4 rows)
+
+Time: 9.226 ms
+--- CREATE SP-GiST INDEX ---
+CREATE INDEX
+Time: 18.724 ms
+--- CREATE GiST INDEX ---
+CREATE INDEX
+Time: 44.994 ms
+     indexname      |  size   
+--------------------+---------
+ bench_catalog_pkey | 336 kB
+ bench_gist         | 2904 kB
+ bench_spgist       | 2344 kB
+(3 rows)
+
+Time: 3.750 ms
+--- BENCHMARK 1: 2h window, Eagle Idaho, 10 deg min_el ---
+SET
+Time: 0.158 ms
+SET
+Time: 0.020 ms
+ seqscan_candidates 
+--------------------
+               2261
+(1 row)
+
+Time: 8.224 ms
+RESET
+Time: 0.102 ms
+RESET
+Time: 0.019 ms
+SET
+Time: 0.023 ms
+ spgist_candidates 
+-------------------
+              2261
+(1 row)
+
+Time: 9.787 ms
+RESET
+Time: 0.142 ms
+--- BENCHMARK 2: 24h window, Eagle Idaho, 10 deg min_el ---
+SET
+Time: 0.044 ms
+SET
+Time: 0.013 ms
+ seqscan_candidates 
+--------------------
+              13562
+(1 row)
+
+Time: 4.272 ms
+RESET
+Time: 0.044 ms
+RESET
+Time: 0.015 ms
+SET
+Time: 0.017 ms
+ spgist_candidates 
+-------------------
+             13562
+(1 row)
+
+Time: 6.832 ms
+RESET
+Time: 0.065 ms
+--- BENCHMARK 3: 2h window, Equator, 10 deg min_el ---
+SET
+Time: 0.025 ms
+SET
+Time: 0.010 ms
+ seqscan_candidates 
+--------------------
+               2073
+(1 row)
+
+Time: 5.868 ms
+RESET
+Time: 1.133 ms
+RESET
+Time: 0.083 ms
+SET
+Time: 0.032 ms
+ spgist_candidates 
+-------------------
+              2073
+(1 row)
+
+Time: 7.401 ms
+RESET
+Time: 0.105 ms
+--- BENCHMARK 4: 2h window, Eagle Idaho, 45 deg min_el ---
+SET
+Time: 0.034 ms
+SET
+Time: 0.010 ms
+ seqscan_candidates 
+--------------------
+               1407
+(1 row)
+
+Time: 5.641 ms
+RESET
+Time: 0.153 ms
+RESET
+Time: 0.018 ms
+SET
+Time: 0.048 ms
+ spgist_candidates 
+-------------------
+              1407
+(1 row)
+
+Time: 6.581 ms
+RESET
+Time: 0.062 ms
+--- CONSISTENCY CHECK ---
+SET
+Time: 0.049 ms
+SET
+Time: 0.012 ms
+SELECT 2261
+Time: 7.979 ms
+RESET
+Time: 0.159 ms
+RESET
+Time: 0.024 ms
+SET
+Time: 0.030 ms
+SELECT 2261
+Time: 7.533 ms
+RESET
+Time: 0.487 ms
+ in_seq_not_idx 
+----------------
+              0
+(1 row)
+
+Time: 1.214 ms
+ in_idx_not_seq 
+----------------
+              0
+(1 row)
+
+Time: 0.864 ms
+DROP TABLE
+Time: 1.814 ms
+--- EXPLAIN ANALYZE: SP-GiST scan ---
+SET
+Time: 0.064 ms
+                                                               QUERY PLAN                                                               
+----------------------------------------------------------------------------------------------------------------------------------------
+ Aggregate  (cost=51.38..51.39 rows=1 width=8) (actual time=7.322..7.325 rows=1.00 loops=1)
+   Buffers: shared hit=1075
+   ->  Bitmap Heap Scan on bench_catalog  (cost=4.38..51.35 rows=14 width=0) (actual time=6.921..7.255 rows=2261.00 loops=1)
+         Recheck Cond: (tle &? '("43.6977N 116.3535W 760m","2026-02-16 19:00:00-07","2026-02-16 21:00:00-07",10)'::observer_window)
+         Heap Blocks: exact=187
+         Buffers: shared hit=1075
+         ->  Bitmap Index Scan on bench_spgist  (cost=0.00..4.38 rows=14 width=0) (actual time=6.887..6.888 rows=2261.00 loops=1)
+               Index Cond: (tle &? '("43.6977N 116.3535W 760m","2026-02-16 19:00:00-07","2026-02-16 21:00:00-07",10)'::observer_window)
+               Index Searches: 1
+               Buffers: shared hit=888
+ Planning Time: 0.143 ms
+ Execution Time: 7.365 ms
+(12 rows)
+
+Time: 7.974 ms
+RESET
+Time: 0.084 ms
+--- EXPLAIN ANALYZE: Sequential scan ---
+SET
+Time: 0.023 ms
+SET
+Time: 0.011 ms
+                                                          QUERY PLAN                                                          
+------------------------------------------------------------------------------------------------------------------------------
+ Aggregate  (cost=470.74..470.75 rows=1 width=8) (actual time=6.037..6.039 rows=1.00 loops=1)
+   Buffers: shared hit=291
+   ->  Seq Scan on bench_catalog  (cost=0.00..470.70 rows=14 width=0) (actual time=0.016..5.952 rows=2261.00 loops=1)
+         Filter: (tle &? '("43.6977N 116.3535W 760m","2026-02-16 19:00:00-07","2026-02-16 21:00:00-07",10)'::observer_window)
+         Rows Removed by Filter: 12115
+         Buffers: shared hit=291
+ Planning Time: 0.130 ms
+ Execution Time: 6.066 ms
+(8 rows)
+
+Time: 6.589 ms
+RESET
+Time: 0.088 ms
+RESET
+Time: 2.471 ms
+DROP TABLE
+Time: 3.314 ms
+Timing is off.

bench/benchmark_spgist_only.sql

@@ -0,0 +1,264 @@
+-- ============================================================
+-- SP-GiST Orbital Trie Benchmark (Phase 3)
+-- CelesTrak active catalog, ~14k satellites
+-- GiST comparison omitted (known crash in gist_tle_picksplit)
+-- ============================================================
+
+\timing on
+
+-- ============================================================
+-- 1. Catalog distribution analysis
+-- ============================================================
+SELECT
+    CASE
+        WHEN tle_perigee(tle) < 2000 THEN 'LEO (<2000km)'
+        WHEN tle_perigee(tle) < 20000 THEN 'MEO (2000-20000km)'
+        WHEN tle_perigee(tle) < 34000 THEN 'GEO-transfer'
+        ELSE 'GEO/HEO (>34000km)'
+    END AS regime,
+    count(*) AS n,
+    round(100.0 * count(*) / (SELECT count(*) FROM bench_catalog), 1) AS pct
+FROM bench_catalog
+GROUP BY 1
+ORDER BY 2 DESC;
+
+-- ============================================================
+-- 2. Create SP-GiST index
+-- ============================================================
+\echo '--- CREATE SP-GiST INDEX ---'
+CREATE INDEX bench_spgist ON bench_catalog USING spgist (tle tle_spgist_ops);
+
+SELECT pg_size_pretty(pg_relation_size('bench_spgist'::regclass)) AS spgist_size;
+
+-- ============================================================
+-- 3. Benchmark: 2h window, Eagle Idaho (43.7N) — RAAN active
+-- ============================================================
+\echo '--- BENCHMARK 1: 2h window, Eagle Idaho, 10 deg min_el ---'
+
+-- 3a. Sequential scan (baseline)
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+\echo 'Sequential scan:'
+EXPLAIN ANALYZE
+SELECT count(*) AS candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+-- 3b. SP-GiST index scan
+SET enable_seqscan = off;
+\echo 'SP-GiST index scan:'
+EXPLAIN ANALYZE
+SELECT count(*) AS candidates
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 4. Benchmark: 24h window, Eagle Idaho — RAAN bypassed
+-- ============================================================
+\echo '--- BENCHMARK 2: 24h window, Eagle Idaho, 10 deg min_el ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+\echo 'Sequential scan:'
+SELECT count(*) AS seqscan_24h
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 00:00:00+00'::timestamptz,
+    '2026-02-18 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+\echo 'SP-GiST index scan:'
+SELECT count(*) AS spgist_24h
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 00:00:00+00'::timestamptz,
+    '2026-02-18 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 5. Benchmark: 2h window, Equatorial observer
+-- ============================================================
+\echo '--- BENCHMARK 3: 2h window, Equator, 10 deg min_el ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+\echo 'Sequential scan:'
+SELECT count(*) AS seqscan_equator
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+\echo 'SP-GiST index scan:'
+SELECT count(*) AS spgist_equator
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 6. Benchmark: High min_el (45 deg)
+-- ============================================================
+\echo '--- BENCHMARK 4: 2h window, Eagle Idaho, 45 deg min_el ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+\echo 'Sequential scan:'
+SELECT count(*) AS seqscan_45deg
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    45.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+\echo 'SP-GiST index scan:'
+SELECT count(*) AS spgist_45deg
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    45.0
+)::observer_window;
+RESET enable_seqscan;
+
+-- ============================================================
+-- 7. Consistency check
+-- ============================================================
+\echo '--- CONSISTENCY CHECK ---'
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+CREATE TEMPORARY TABLE seq_results AS
+SELECT norad_id FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+SET enable_seqscan = off;
+CREATE TEMPORARY TABLE idx_results AS
+SELECT norad_id FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+RESET enable_seqscan;
+
+SELECT count(*) AS in_seq_not_idx FROM seq_results
+WHERE norad_id NOT IN (SELECT norad_id FROM idx_results);
+
+SELECT count(*) AS in_idx_not_seq FROM idx_results
+WHERE norad_id NOT IN (SELECT norad_id FROM seq_results);
+
+DROP TABLE seq_results, idx_results;
+
+-- ============================================================
+-- 8. Pruning summary
+-- ============================================================
+\echo '--- PRUNING SUMMARY ---'
+SELECT
+    '2h/Eagle/10deg' AS scenario,
+    (SELECT count(*) FROM bench_catalog) AS catalog_size,
+    count(*) AS candidates,
+    round(100.0 * count(*) / (SELECT count(*) FROM bench_catalog), 1) AS candidate_pct,
+    round(100.0 * (1.0 - count(*)::numeric / (SELECT count(*) FROM bench_catalog)), 1) AS pruning_pct
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+UNION ALL
+SELECT
+    '24h/Eagle/10deg',
+    (SELECT count(*) FROM bench_catalog),
+    count(*),
+    round(100.0 * count(*) / (SELECT count(*) FROM bench_catalog), 1),
+    round(100.0 * (1.0 - count(*)::numeric / (SELECT count(*) FROM bench_catalog)), 1)
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 00:00:00+00'::timestamptz,
+    '2026-02-18 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+UNION ALL
+SELECT
+    '2h/Equator/10deg',
+    (SELECT count(*) FROM bench_catalog),
+    count(*),
+    round(100.0 * count(*) / (SELECT count(*) FROM bench_catalog), 1),
+    round(100.0 * (1.0 - count(*)::numeric / (SELECT count(*) FROM bench_catalog)), 1)
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+UNION ALL
+SELECT
+    '2h/Eagle/45deg',
+    (SELECT count(*) FROM bench_catalog),
+    count(*),
+    round(100.0 * count(*) / (SELECT count(*) FROM bench_catalog), 1),
+    round(100.0 * (1.0 - count(*)::numeric / (SELECT count(*) FROM bench_catalog)), 1)
+FROM bench_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2026-02-17 02:00:00+00'::timestamptz,
+    '2026-02-17 04:00:00+00'::timestamptz,
+    45.0
+)::observer_window;
+
+-- ============================================================
+-- Cleanup
+-- ============================================================
+DROP INDEX bench_spgist;
+DROP TABLE bench_catalog;
+
+\timing off

bench/load_bench.sh

@@ -0,0 +1,145 @@
+#!/bin/bash
+# Load pg_orrery benchmark catalog into PostgreSQL.
+#
+# Uses pg-orrery-catalog if available, falls back to pre-generated SQL.
+#
+# Usage:
+#   ./bench/load_bench.sh                       # Load from cached SQL or TLE files
+#   ./bench/load_bench.sh --rebuild              # Re-merge from individual source files
+#   ./bench/load_bench.sh --download             # Re-download sources + rebuild + load
+#
+# Environment:
+#   PGPORT       PostgreSQL port (default: 5499)
+#   PGDATABASE   Target database  (default: contrib_regression)
+#   SOCKS_PROXY  SOCKS5 proxy for CelesTrak (default: none)
+#
+set -euo pipefail
+
+BENCH_DIR="$(cd "$(dirname "$0")" && pwd)"
+PGPORT="${PGPORT:-5499}"
+PGDATABASE="${PGDATABASE:-contrib_regression}"
+TABLE="bench_catalog"
+REBUILD=false
+DOWNLOAD=false
+
+for arg in "$@"; do
+    case "$arg" in
+        --rebuild)  REBUILD=true ;;
+        --download) DOWNLOAD=true; REBUILD=true ;;
+        --help|-h)
+            head -14 "$0" | tail -13 | sed 's/^# \?//'
+            exit 0 ;;
+    esac
+done
+
+# ── Check for pg-orrery-catalog ──────────────────────────────
+HAS_CATALOG=false
+if command -v pg-orrery-catalog &>/dev/null; then
+    HAS_CATALOG=true
+elif [ -f "$BENCH_DIR/../pg-orrery-catalog/.venv/bin/pg-orrery-catalog" ]; then
+    # Sibling development checkout
+    export PATH="$BENCH_DIR/../pg-orrery-catalog/.venv/bin:$PATH"
+    HAS_CATALOG=true
+fi
+
+# ── Download sources ─────────────────────────────────────────
+if $DOWNLOAD; then
+    if $HAS_CATALOG; then
+        echo "==> Downloading TLE sources via pg-orrery-catalog..."
+        pg-orrery-catalog download --force
+    else
+        echo "==> pg-orrery-catalog not found, downloading via curl..."
+
+        CURL_PROXY=""
+        [ -n "${SOCKS_PROXY:-}" ] && CURL_PROXY="--socks5-hostname $SOCKS_PROXY"
+
+        # CelesTrak active (no auth needed)
+        CURL_CT="/usr/bin/curl -s $CURL_PROXY --connect-timeout 15 --max-time 120"
+        echo "  CelesTrak active..."
+        $CURL_CT "https://celestrak.org/NORAD/elements/gp.php?GROUP=active&FORMAT=3le" \
+            -o "$BENCH_DIR/celestrak_active.tle" 2>/dev/null || echo "    FAILED"
+
+        # CelesTrak supplemental GP
+        for group in starlink oneweb planet orbcomm; do
+            echo "  CelesTrak SupGP ${group}..."
+            $CURL_CT "https://celestrak.org/NORAD/elements/supplemental/sup-gp.php?FILE=${group}&FORMAT=3le" \
+                -o "$BENCH_DIR/supgp_${group}.tle" 2>/dev/null || true
+        done
+
+        REBUILD=true
+    fi
+fi
+
+# ── Build SQL ────────────────────────────────────────────────
+if $REBUILD; then
+    if $HAS_CATALOG; then
+        echo "==> Building catalog via pg-orrery-catalog..."
+        # Use cached downloads if available, fall back to bench/ TLE files
+        SOURCES=()
+        for f in "$BENCH_DIR"/*.tle; do
+            [ -f "$f" ] && SOURCES+=("$f")
+        done
+        if [ ${#SOURCES[@]} -gt 0 ]; then
+            pg-orrery-catalog build "${SOURCES[@]}" --table "$TABLE" \
+                > "$BENCH_DIR/load_mega_catalog.sql"
+        else
+            pg-orrery-catalog build --table "$TABLE" \
+                > "$BENCH_DIR/load_mega_catalog.sql"
+        fi
+        echo "  Generated load_mega_catalog.sql"
+    else
+        echo "==> Building catalog via build_catalog.py..."
+        SOURCES=()
+        for f in spacetrack_everything.tle celestrak_active.tle satnogs_full.tle \
+                 supgp_starlink.tle supgp_oneweb.tle supgp_planet.tle supgp_orbcomm.tle; do
+            [ -f "$BENCH_DIR/$f" ] && SOURCES+=("$BENCH_DIR/$f")
+        done
+
+        if [ ${#SOURCES[@]} -eq 0 ]; then
+            echo "ERROR: No source TLE files found in $BENCH_DIR" >&2
+            exit 1
+        fi
+
+        python3 "$BENCH_DIR/build_catalog.py" "${SOURCES[@]}" \
+            > "$BENCH_DIR/load_mega_catalog.sql"
+        echo "  Generated load_mega_catalog.sql"
+    fi
+fi
+
+# ── Load into PostgreSQL ─────────────────────────────────────
+if [ ! -f "$BENCH_DIR/load_mega_catalog.sql" ]; then
+    echo "ERROR: $BENCH_DIR/load_mega_catalog.sql not found" >&2
+    echo "  Run with --rebuild or --download first" >&2
+    exit 1
+fi
+
+echo "==> Loading catalog into $PGDATABASE (port $PGPORT)..."
+PGPORT=$PGPORT psql -d "$PGDATABASE" -f "$BENCH_DIR/load_mega_catalog.sql" -q 2>&1 | tail -3
+
+# ── Create indexes ───────────────────────────────────────────
+echo "==> Creating indexes..."
+PGPORT=$PGPORT psql -d "$PGDATABASE" -q << 'SQL'
+\timing on
+CREATE INDEX IF NOT EXISTS bench_spgist_idx ON bench_catalog USING spgist (tle tle_spgist_ops);
+CREATE INDEX IF NOT EXISTS bench_gist_idx ON bench_catalog USING gist (tle);
+\timing off
+SQL
+
+# ── Summary ──────────────────────────────────────────────────
+PGPORT=$PGPORT psql -d "$PGDATABASE" -q << 'SQL'
+SELECT count(*) || ' objects loaded' AS status FROM bench_catalog;
+SELECT
+    CASE
+        WHEN tle_mean_motion(tle) > 11.25 THEN 'LEO'
+        WHEN tle_mean_motion(tle) > 1.8 THEN 'MEO'
+        WHEN tle_mean_motion(tle) > 0.9 THEN 'GEO'
+        ELSE 'HEO'
+    END AS regime,
+    count(*) AS count
+FROM bench_catalog
+GROUP BY 1
+ORDER BY 2 DESC;
+SQL
+
+echo "==> Done.  Run benchmarks with:"
+echo "  PGPORT=$PGPORT psql -d $PGDATABASE -f bench/benchmark.sql"

bench/tle_to_sql.py

@@ -0,0 +1,50 @@
+#!/usr/bin/env python3
+"""Convert 3-line TLE file to SQL COPY format for pg_orrery benchmark."""
+import sys
+
+def main():
+    lines = open(sys.argv[1]).read().strip().split('\n')
+    print("-- CelesTrak active satellite catalog")
+    print("-- Auto-generated for SP-GiST benchmark")
+    print("CREATE TABLE IF NOT EXISTS bench_catalog (")
+    print("    norad_id integer PRIMARY KEY,")
+    print("    name text NOT NULL,")
+    print("    tle tle NOT NULL")
+    print(");")
+    print("TRUNCATE bench_catalog;")
+    print()
+
+    count = 0
+    errors = 0
+    i = 0
+    while i + 2 < len(lines):
+        name = lines[i].strip()
+        line1 = lines[i + 1].strip()
+        line2 = lines[i + 2].strip()
+
+        # Validate TLE format
+        if not line1.startswith('1 ') or not line2.startswith('2 '):
+            i += 1  # skip and try to resync
+            errors += 1
+            continue
+
+        # Extract NORAD ID from line 1 (cols 3-7)
+        try:
+            norad_id = int(line1[2:7].strip())
+        except ValueError:
+            i += 3
+            errors += 1
+            continue
+
+        # Escape single quotes in name
+        name_escaped = name.replace("'", "''")
+        tle_str = line1 + '\n' + line2
+
+        print(f"INSERT INTO bench_catalog VALUES ({norad_id}, '{name_escaped}', E'{tle_str}') ON CONFLICT (norad_id) DO NOTHING;")
+        count += 1
+        i += 3
+
+    print(f"\n-- Loaded {count} satellites ({errors} parse errors skipped)")
+
+if __name__ == '__main__':
+    main()

docs/astro.config.mjs

@@ -69,6 +69,8 @@ export default defineConfig({
             { label: "Conjunction Screening", slug: "guides/conjunction-screening" },
             { label: "JPL DE Ephemeris", slug: "guides/de-ephemeris" },
             { label: "Orbit Determination", slug: "guides/orbit-determination" },
+            { label: "Satellite Pass Prediction", slug: "guides/pass-prediction" },
+            { label: "Building TLE Catalogs", slug: "guides/catalog-management" },
           ],
         },
         {
@@ -95,7 +97,7 @@ export default defineConfig({
             { label: "Functions: Transfers", slug: "reference/functions-transfers" },
             { label: "Functions: DE Ephemeris", slug: "reference/functions-de" },
             { label: "Functions: Orbit Determination", slug: "reference/functions-od" },
-            { label: "Operators & GiST Index", slug: "reference/operators-gist" },
+            { label: "Operators & Indexes", slug: "reference/operators-gist" },
             { label: "Body ID Reference", slug: "reference/body-ids" },
             { label: "Constants & Accuracy", slug: "reference/constants-accuracy" },
           ],

docs/src/content/docs/architecture/design-principles.mdx

@@ -104,7 +104,7 @@ For v0.1.0/v0.2.0 functions, there are no file-scope variables, no static locals
 
 ### Fixed-size types
 
-All seven pg_orrery types use `STORAGE = plain` and fixed `INTERNALLENGTH`. No TOAST, no detoasting, no variable-length headers. The `tle` type is exactly 112 bytes. Direct pointer access via `PG_GETARG_POINTER(n)` --- no copies, no allocations on read.
+All seven pg_orrery base types use `STORAGE = plain` and fixed `INTERNALLENGTH`. No TOAST, no detoasting, no variable-length headers. The `tle` type is exactly 112 bytes. Direct pointer access via `PG_GETARG_POINTER(n)` --- no copies, no allocations on read. The eighth type, `observer_window`, is a SQL composite used only as a query-time parameter --- it is never stored in table columns.
 
 ### Deterministic memory
 

docs/src/content/docs/getting-started/installation.mdx

@@ -97,13 +97,13 @@ import { Tabs, TabItem, Steps, Aside } from "@astrojs/starlight/components";
 
 ## Running the test suite
 
-If building from source, the regression tests verify all functions across 14 test suites:
+If building from source, the regression tests verify all functions across 15 test suites:
 
 ```bash
 make installcheck PG_CONFIG=/usr/bin/pg_config
 ```
 
-This runs the tests listed in the `REGRESS` variable: TLE parsing, SGP4 propagation, coordinate transforms, pass prediction, GiST indexing, convenience functions, star observation, Keplerian propagation, planet observation, moon observation, Lambert transfers, and DE ephemeris.
+This runs the tests listed in the `REGRESS` variable: TLE parsing, SGP4 propagation, coordinate transforms, pass prediction, GiST indexing, convenience functions, star observation, Keplerian propagation, planet observation, moon observation, Lambert transfers, DE ephemeris, orbit determination, SP-GiST visibility index, and the 518 Vallado test vectors.
 
 ## Upgrading
 

docs/src/content/docs/guides/catalog-management.mdx

@@ -0,0 +1,309 @@
+---
+title: Building TLE Catalogs
+sidebar:
+  order: 12
+---
+
+import { Steps, Aside, Tabs, TabItem, Code } from "@astrojs/starlight/components";
+
+Every pg_orrery workflow starts with TLEs in a table. The [Tracking Satellites](/guides/tracking-satellites/) guide shows how to insert a few satellites by hand --- but a real catalog has tens of thousands of objects from multiple sources, each with different freshness and coverage. `pg-orrery-catalog` handles the download, merge, and load pipeline.
+
+## The problem with multiple TLE sources
+
+Three major sources provide TLE data, each with trade-offs:
+
+| Source | Auth | Coverage | Freshness |
+|--------|------|----------|-----------|
+| [Space-Track](https://www.space-track.org) | Login required | Full catalog (~30k+ on-orbit) | Hours to days |
+| [CelesTrak](https://celestrak.org) | None | Active sats + operator supplemental GP | Minutes to hours |
+| [SatNOGS](https://db.satnogs.org) | None | Community-tracked objects | Varies |
+
+The same satellite often appears in all three. CelesTrak's supplemental GP (SupGP) data is particularly valuable --- operators like SpaceX submit Starlink ephemerides that are often hours fresher than Space-Track's own catalog.
+
+The question is which entry to keep. `pg-orrery-catalog` answers with epoch-based deduplication: when the same NORAD ID appears in multiple sources, the entry with the newest epoch wins. This means SupGP data automatically overrides stale Space-Track entries where available.
+
+## Install
+
+```bash
+# Run directly (no install needed)
+uvx pg-orrery-catalog --help
+
+# Or install permanently
+uv pip install pg-orrery-catalog
+
+# For direct database loading (adds psycopg)
+uv pip install "pg-orrery-catalog[pg]"
+```
+
+## Download, build, load
+
+The typical workflow is three steps. Each can run independently.
+
+<Steps>
+1. **Download** TLE data from remote sources into the local cache:
+
+   ```bash
+   pg-orrery-catalog download
+   ```
+
+   This fetches from all configured sources (CelesTrak by default, Space-Track if credentials are set). Files are cached in `~/.cache/pg-orrery-catalog/` and reused unless stale (>24h) or `--force` is passed.
+
+   To download from a specific source:
+
+   ```bash
+   pg-orrery-catalog download --source celestrak
+   pg-orrery-catalog download --source spacetrack --force
+   ```
+
+2. **Build** a merged catalog and output it:
+
+   <Tabs>
+     <TabItem label="Pipe to psql">
+     ```bash
+     pg-orrery-catalog build | psql -d mydb
+     ```
+     </TabItem>
+     <TabItem label="Save SQL file">
+     ```bash
+     pg-orrery-catalog build --table satellites -o catalog.sql
+     ```
+     </TabItem>
+     <TabItem label="Export 3LE">
+     ```bash
+     pg-orrery-catalog build --format 3le -o merged.tle
+     ```
+     </TabItem>
+     <TabItem label="JSON output">
+     ```bash
+     pg-orrery-catalog build --format json -o catalog.json
+     ```
+     </TabItem>
+   </Tabs>
+
+   With no arguments, `build` merges all cached files. You can also pass specific TLE files:
+
+   ```bash
+   pg-orrery-catalog build /path/to/spacetrack.tle /path/to/celestrak.tle
+   ```
+
+   The merge reports what happened:
+
+   ```
+   spacetrack_everything: 33053 objects (33053 new, 0 updated)
+   celestrak_active: 14376 objects (2 new, 0 updated)
+   satnogs_full: 1488 objects (121 new, 5 updated)
+   supgp_starlink: 9703 objects (77 new, 7398 updated)
+   Total: 33253 unique objects
+   Regimes: LEO: 31542, GEO: 1203, MEO: 385, HEO: 123
+   ```
+
+   Notice how SupGP updated 7,398 Starlink entries --- those are fresher epochs from SpaceX overriding stale Space-Track data.
+
+3. **Load** directly into PostgreSQL (requires `[pg]` extra):
+
+   ```bash
+   pg-orrery-catalog load \
+     --database-url postgresql:///mydb \
+     --table satellites \
+     --create-index
+   ```
+
+   The `--create-index` flag creates both GiST and SP-GiST indexes on the `tle` column, ready for spatial queries and KNN ordering.
+</Steps>
+
+## Configuration
+
+Three layers, highest precedence first:
+
+1. **CLI flags** --- `--table`, `--source`, `--database-url`
+2. **Environment variables** --- `SPACETRACK_USER`, `SOCKS_PROXY`, `DATABASE_URL`
+3. **Config file** --- `~/.config/pg-orrery-catalog/config.toml`
+
+### Space-Track credentials
+
+Space-Track requires a free account. Set credentials via environment variables:
+
+```bash
+export SPACETRACK_USER="you@example.com"
+export SPACETRACK_PASSWORD="secret"
+pg-orrery-catalog download --source spacetrack
+```
+
+Or in the config file:
+
+```toml
+[spacetrack]
+user = "you@example.com"
+password = "secret"
+```
+
+### SOCKS proxy
+
+CelesTrak is sometimes unreachable from certain networks. Route through a SOCKS5 proxy:
+
+```bash
+export SOCKS_PROXY="localhost:1080"
+pg-orrery-catalog download
+```
+
+### Full config reference
+
+```toml
+[spacetrack]
+user = "you@example.com"
+password = "secret"
+
+[celestrak]
+proxy = "localhost:1080"
+supgp_groups = ["starlink", "oneweb", "planet", "orbcomm"]
+
+[output]
+table = "satellites"
+
+[database]
+url = "postgresql://localhost/mydb"
+```
+
+## Working with the generated SQL
+
+The SQL output creates a table with three columns:
+
+```sql
+CREATE TABLE satellites (
+    id   serial,
+    name text,
+    tle  tle
+);
+```
+
+Once loaded, the full pg_orrery function set is available:
+
+```sql
+-- Where is every LEO satellite right now?
+SELECT name, observe(tle, '40.0N 105.3W 1655m'::observer, now()) AS topo
+FROM satellites
+WHERE tle_mean_motion(tle) > 11.25;
+
+-- Which satellites are overhead right now?
+SELECT name,
+       round(topo_elevation(
+         observe_safe(tle, '40.0N 105.3W 1655m'::observer, now())
+       )::numeric, 1) AS el
+FROM satellites
+WHERE topo_elevation(
+        observe_safe(tle, '40.0N 105.3W 1655m'::observer, now())
+      ) > 10
+ORDER BY el DESC;
+
+-- Predict ISS passes for the next 24 hours
+SELECT pass_aos_time(p)::timestamp(0) AS rise,
+       round(pass_max_elevation(p)::numeric, 1) AS max_el,
+       pass_los_time(p)::timestamp(0) AS set
+FROM satellites,
+     predict_passes(tle, '40.0N 105.3W 1655m'::observer,
+                    now(), now() + interval '24 hours', 10.0) p
+WHERE tle_norad_id(tle) = 25544;
+```
+
+## NORAD ID encoding
+
+TLE files use a 5-character field for the NORAD catalog number. With more than 100,000 tracked objects, the original 5-digit numeric format ran out of space. The encoding has evolved through four cases:
+
+| Case | Format | Range | Example |
+|------|--------|-------|---------|
+| Traditional | `ddddd` | 0 -- 99,999 | `25544` (ISS) |
+| Alpha-5 | `Ldddd` | 100,000 -- 339,999 | `T0002` = 270,002 |
+| Super-5 case 3 | `xxxxX` | 340,000 -- 906,309,663 | `0000A` = 340,000 |
+| Super-5 case 4 | `xxxXd` | 906,309,664+ | `000A0` = 906,309,664 |
+
+Alpha-5 skips the letters I and O (they look like 1 and 0). Super-5 uses a base-64 alphabet: digits 0--9, uppercase A--Z, lowercase a--z, plus `+` and `-`.
+
+`pg-orrery-catalog` decodes all four cases, matching the `get_norad_number()` implementation in pg_orrery's vendored SGP4 library. This means Alpha-5 objects like Starlink satellites (NORAD IDs above 100,000) load correctly.
+
+<Aside type="note" title="Alpha-5 verification">
+You can verify the decoding independently:
+
+```bash
+python3 -c "
+from pg_orrery_catalog.tle import decode_norad
+print(f'T0002 = {decode_norad(\"T0002\")}')  # 270002
+print(f'A0001 = {decode_norad(\"A0001\")}')  # 100001
+print(f'Z9999 = {decode_norad(\"Z9999\")}')  # 339999
+"
+```
+</Aside>
+
+## Cache management
+
+Downloaded TLE files are stored under `~/.cache/pg-orrery-catalog/`, organized by source:
+
+```
+~/.cache/pg-orrery-catalog/
+  celestrak/
+    celestrak_active.tle
+    supgp_starlink.tle
+    supgp_oneweb.tle
+    ...
+  satnogs/
+    satnogs_full.tle
+  spacetrack/
+    spacetrack_everything.tle
+```
+
+Check what's cached:
+
+```bash
+pg-orrery-catalog info --cache
+```
+
+Files older than 24 hours are considered stale and re-downloaded automatically. Use `--force` to override fresh cache entries.
+
+## Automating catalog updates
+
+For a regularly-updated catalog, a cron job or systemd timer works well:
+
+```bash
+# Update catalog daily at 03:00
+0 3 * * * pg-orrery-catalog download && pg-orrery-catalog build --table satellites | psql -d mydb
+```
+
+Or with the direct load command:
+
+```bash
+0 3 * * * pg-orrery-catalog download && pg-orrery-catalog load --database-url postgresql:///mydb --table satellites --create-index
+```
+
+<Aside type="caution" title="Table replacement">
+The default SQL output includes `DROP TABLE IF EXISTS` before `CREATE TABLE`. This replaces the entire table on each load. If you need to preserve the table and upsert, export as JSON and handle the merge in your application logic.
+</Aside>
+
+## Using as a library
+
+`pg-orrery-catalog` can also be imported as a Python library:
+
+```python
+from pg_orrery_catalog.tle import decode_norad, parse_3le_file
+from pg_orrery_catalog.catalog import merge_sources
+from pg_orrery_catalog.regime import regime_summary
+from pg_orrery_catalog.output.sql import generate_sql
+
+# Parse and merge
+merged, stats = merge_sources(["spacetrack.tle", "celestrak.tle"])
+print(f"{stats.total_unique} unique objects")
+
+# Classify
+regimes = regime_summary(merged)
+print(regimes)  # {'LEO': 31542, 'MEO': 385, 'GEO': 1203, 'HEO': 123}
+
+# Generate SQL
+sql = generate_sql(merged, table="my_catalog")
+```
+
+## What's next
+
+With a catalog loaded, see:
+
+- [Tracking Satellites](/guides/tracking-satellites/) --- observe, predict passes, screen conjunctions
+- [Satellite Pass Prediction](/guides/pass-prediction/) --- detailed pass prediction workflows
+- [Conjunction Screening](/guides/conjunction-screening/) --- find close approaches using GiST indexes
+- [Benchmarks](/performance/benchmarks/) --- performance data with catalogs of 33k--66k objects

docs/src/content/docs/guides/pass-prediction.mdx

@@ -0,0 +1,224 @@
+---
+title: Satellite Pass Prediction
+sidebar:
+  order: 11
+---
+
+import { Steps, Aside, Tabs, TabItem } from "@astrojs/starlight/components";
+
+Satellite pass prediction answers a deceptively simple question: "which satellites will fly over me tonight?" The brute-force approach -- propagating every object in a 30,000-satellite catalog with SGP4 at 10-second intervals over a 24-hour window -- requires billions of floating-point operations. pg_orrery solves this with the `&?` visibility cone operator, which applies three geometric filters (altitude, inclination, RAAN) to eliminate 80-90% of candidates without any SGP4 propagation. An optional SP-GiST index provides tree-level pruning for large catalogs.
+
+## How you do it today
+
+Ground station operators and amateur observers use several approaches:
+
+- **Heavens-Above / N2YO**: Web tools that compute passes for a single observer. Great for casual use. Not queryable from SQL; you scrape or check manually.
+- **Skyfield / PyEphem / predict**: Python libraries that propagate TLEs and compute topocentric coordinates. You write a loop over the catalog and check each object. Works, but scales linearly with catalog size.
+- **GPredict / Orbitron**: Desktop applications for pass prediction. Local install, single observer, no database integration.
+- **Custom scripts**: Propagate everything, compute elevation, filter. For a full catalog this takes minutes per observer per day.
+
+The common thread: every approach propagates every satellite. There is no pre-filtering. A 30,000-object catalog takes the same time whether 2 satellites or 2,000 are visible from your location.
+
+## What changes with pg_orrery
+
+pg_orrery attacks the problem in two stages:
+
+**Stage 1: SP-GiST index eliminates impossible candidates.** The `&?` operator checks whether a satellite _could_ be visible from a ground observer during a time window, using three geometric filters applied without SGP4 propagation:
+
+| Filter | What it checks | What it eliminates |
+|---|---|---|
+| Altitude | Perigee too high for the observer's minimum elevation angle | MEO/GEO satellites for high-elevation queries |
+| Inclination | Inclination + footprint angle must reach the observer's latitude | Equatorial satellites from high-latitude observers |
+| RAAN | Right Ascension of Ascending Node alignment with observer's local sidereal time | Satellites whose orbital plane isn't overhead during the query window |
+
+**Stage 2: SGP4 propagation on survivors.** The handful of candidates that pass the geometric filter are propagated with `predict_passes()` to find exact AOS/LOS times and maximum elevation.
+
+The key type for queries is `observer_window`:
+
+| Field | Type | Meaning |
+|---|---|---|
+| `obs` | `observer` | Ground location (lat, lon, altitude) |
+| `t_start` | `timestamptz` | Start of observation window |
+| `t_end` | `timestamptz` | End of observation window |
+| `min_el` | `float8` | Minimum elevation angle in degrees |
+
+## What pg_orrery does not replace
+
+<Aside type="caution" title="Geometric filter, not a propagator">
+The `&?` operator answers "could this satellite possibly be visible?" -- not "will it definitely pass overhead." It is a conservative superset: it may include satellites that do not actually produce a visible pass (false positives), but it will never exclude one that does (no false negatives). Always follow with `predict_passes()` for ground truth.
+</Aside>
+
+- **Not a pass schedule.** The `&?` operator does not compute AOS, LOS, or maximum elevation. Use `predict_passes()` on the candidates for precise timing.
+- **No optical visibility.** The filter considers only geometric visibility (above the horizon at the required elevation). It does not check whether the satellite is sunlit, whether the observer is in darkness, or whether the pass is bright enough to see. Use `pass_visible()` to check illumination.
+- **J2-only RAAN.** The RAAN filter projects the ascending node using only the J2 zonal harmonic. For short query windows (< 4 hours) the error is small. For long windows (> 12 hours) the RAAN filter automatically disables itself (full Earth rotation makes it meaningless).
+
+## Try it
+
+### Set up a catalog with the SP-GiST index
+
+```sql
+CREATE TABLE catalog (
+    norad_id  integer PRIMARY KEY,
+    name      text NOT NULL,
+    tle       tle NOT NULL
+);
+
+-- Insert your TLEs (from CelesTrak, Space-Track, or any provider)
+-- ISS example:
+INSERT INTO catalog VALUES (25544, 'ISS',
+    '1 25544U 98067A   24001.50000000  .00016717  00000-0  10270-3 0  9025
+2 25544  51.6400 208.9163 0006703  30.1694  61.7520 15.50100486 00001');
+
+-- Create the SP-GiST orbital trie index
+CREATE INDEX catalog_spgist ON catalog USING spgist (tle tle_spgist_ops);
+```
+
+<Aside type="note">
+The SP-GiST index is **not** the default operator class. You must explicitly specify `tle_spgist_ops`. The existing GiST index (`tle_ops`, used by `&&` and `<->`) is unaffected. Both indexes can coexist on the same table.
+</Aside>
+
+### Query: which satellites might be visible tonight?
+
+```sql
+-- Eagle, Idaho: 43.7N 116.4W, 760m elevation
+-- Tonight's 6-hour window, minimum 10 deg elevation
+SELECT name
+FROM catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 01:00:00+00'::timestamptz,
+    '2024-01-01 07:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+```
+
+This query runs the three geometric filters (altitude, inclination, RAAN) on every TLE in the catalog. With the SP-GiST index, PostgreSQL prunes entire subtrees of the index without examining individual TLEs.
+
+### Full pass prediction workflow
+
+<Steps>
+1. **Filter candidates with `&?`:**
+
+   ```sql
+   CREATE TEMPORARY TABLE candidates AS
+   SELECT norad_id, name, tle
+   FROM catalog
+   WHERE tle &? ROW(
+       observer('43.6977N 116.3535W 760m'),
+       '2024-01-01 01:00:00+00'::timestamptz,
+       '2024-01-01 07:00:00+00'::timestamptz,
+       10.0
+   )::observer_window;
+   ```
+
+   For a 30,000-object catalog, this typically returns a few hundred candidates.
+
+2. **Compute actual passes on survivors:**
+
+   ```sql
+   SELECT c.name,
+          (p).aos_time, (p).los_time,
+          round((p).max_elevation::numeric, 1) AS max_el,
+          round((p).aos_azimuth::numeric, 1) AS aos_az
+   FROM candidates c,
+   LATERAL predict_passes(
+       c.tle,
+       observer('43.6977N 116.3535W 760m'),
+       '2024-01-01 01:00:00+00'::timestamptz,
+       '2024-01-01 07:00:00+00'::timestamptz,
+       10.0
+   ) AS p
+   ORDER BY (p).aos_time;
+   ```
+
+   Only the geometric survivors go through full SGP4 propagation.
+
+3. **Filter for optical visibility (optional):**
+
+   ```sql
+   SELECT c.name,
+          (p).aos_time, (p).los_time,
+          round((p).max_elevation::numeric, 1) AS max_el
+   FROM candidates c,
+   LATERAL predict_passes(
+       c.tle,
+       observer('43.6977N 116.3535W 760m'),
+       '2024-01-01 01:00:00+00'::timestamptz,
+       '2024-01-01 07:00:00+00'::timestamptz,
+       10.0
+   ) AS p
+   WHERE pass_visible(c.tle, observer('43.6977N 116.3535W 760m'), (p).aos_time)
+   ORDER BY (p).max_elevation DESC;
+   ```
+
+   This checks whether the satellite is sunlit while the observer is in darkness -- the condition for naked-eye visibility.
+</Steps>
+
+### Comparing query windows
+
+<Tabs>
+<TabItem label="Short window (2h)">
+```sql
+-- Short window: RAAN filter is most aggressive
+-- Only satellites whose orbital plane is currently
+-- aligned with the observer's meridian pass through
+SELECT count(*) AS candidates
+FROM catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 02:00:00+00'::timestamptz,
+    '2024-01-01 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+```
+</TabItem>
+<TabItem label="Full day (24h)">
+```sql
+-- 24-hour window: RAAN filter bypassed (full Earth rotation)
+-- Only altitude and inclination filters apply
+SELECT count(*) AS candidates
+FROM catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+```
+</TabItem>
+<TabItem label="Equatorial observer">
+```sql
+-- Equatorial observer: all inclinations reach latitude 0
+-- Only altitude and RAAN filters apply
+SELECT count(*) AS candidates
+FROM catalog
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2024-01-01 02:00:00+00'::timestamptz,
+    '2024-01-01 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+```
+</TabItem>
+</Tabs>
+
+### Using both GiST and SP-GiST indexes
+
+The two index types serve different purposes and coexist on the same table:
+
+```sql
+-- SP-GiST: "which satellites could I see tonight?"
+CREATE INDEX catalog_spgist ON catalog USING spgist (tle tle_spgist_ops);
+
+-- GiST: "which satellites share an orbital shell?"
+CREATE INDEX catalog_gist ON catalog USING gist (tle tle_ops);
+
+-- Both queries use their respective indexes
+SELECT name FROM catalog WHERE tle &? ROW(...)::observer_window;  -- SP-GiST
+SELECT name FROM catalog WHERE tle && other_tle;                   -- GiST
+```
+
+<Aside type="tip" title="When to use which index">
+Use the **SP-GiST** index (`&?` operator) for observer-to-catalog queries: "what can I see from here?" Use the **GiST** index (`&&` and `<->` operators) for catalog-to-catalog queries: "which satellites are in the same orbital shell?" The SP-GiST index partitions by SMA and inclination with a query-time RAAN filter. The GiST index partitions by altitude band and inclination range for overlap detection.
+</Aside>

docs/src/content/docs/guides/tracking-satellites.mdx

@@ -43,7 +43,7 @@ pg_orrery propagates TLEs and computes look angles. It does not replace the full
 
 - **No real-time GUI.** GPredict and STK provide map displays, polar plots, and Doppler displays. pg_orrery returns numbers. Use any visualization tool to render its output.
 - **No rotator control.** Hamlib drives antenna rotators. pg_orrery computes the azimuth and elevation values Hamlib would consume, but it has no hardware interface.
-- **No TLE fetching.** Bring your own TLEs from Space-Track, CelesTrak, or any provider. pg_orrery parses and propagates them.
+- **TLE fetching via companion tool.** pg_orrery itself doesn't download TLEs, but [`pg-orrery-catalog`](/guides/catalog-management/) handles the full pipeline: download from Space-Track, CelesTrak, and SatNOGS, merge with epoch-based dedup, and load into PostgreSQL.
 - **Orbit determination available.** Since v0.4.0, pg_orrery can fit TLEs from ECI, topocentric, or angles-only observations via differential correction. See the [Orbit Determination guide](/guides/orbit-determination/).
 - **No high-precision propagation.** SGP4/SDP4 accuracy degrades with TLE age. For operational conjunction assessment, use SP ephemerides or owner/operator-provided state vectors. pg_orrery's GiST screening finds candidates; you verify with better data.
 

docs/src/content/docs/performance/benchmarks.mdx

@@ -17,6 +17,7 @@ All benchmarks use PostgreSQL's `EXPLAIN (ANALYZE, BUFFERS)` for timing. The num
 | Operation | Count | Time | Rate | Notes |
 |-----------|-------|------|------|-------|
 | TLE propagation (SGP4) | 12,000 | 17 ms | 706K/sec | Mixed LEO/MEO/GEO |
+| Visibility cone filter (`&?`) | 65,886 | 12.5 ms | 5.3M/sec | 80% pruned, no SGP4 |
 | Planet observation (VSOP87) | 875 | 57 ms | 15.4K/sec | All 7 non-Earth planets, 125 times each |
 | Galilean moon observation | 1,000 | 63 ms | 15.9K/sec | L1.2 + VSOP87 pipeline |
 | Saturn moon observation | 800 | 53 ms | 15.1K/sec | TASS17 + VSOP87 |
@@ -219,12 +220,71 @@ FROM predict_passes(
 
 A 7-day window at 30-second coarse scan resolution requires ~20,160 propagation calls for the coarse scan, plus bisection and ternary search calls for each pass found. Typical ISS result: 25--35 passes found in ~40 ms.
 
+## Visibility cone filtering (`&?` operator)
+
+The `&?` operator answers "could this satellite possibly be visible from this observer?" using three geometric filters (altitude, inclination, RAAN) without any SGP4 propagation. This is the first stage of the pass prediction pipeline, reducing the number of satellites that need full propagation.
+
+```sql
+-- Benchmark: filter a 66,000-object catalog
+-- Eagle, Idaho: 2-hour window, 10 deg minimum elevation
+EXPLAIN (ANALYZE, BUFFERS)
+SELECT count(*)
+FROM satellite_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 02:00:00+00'::timestamptz,
+    '2024-01-01 04:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+```
+
+**65,886 TLEs filtered in 12.5 ms --- 80% pruned, 12,964 candidates survive.**
+
+The operator evaluates three geometric conditions per TLE: perigee altitude vs. maximum visible altitude, inclination + ground footprint vs. observer latitude, and RAAN alignment via J2 secular precession. Each check is a few floating-point operations --- no SGP4 initialization, no Kepler equation, no trigonometric series.
+
+### Pruning rate by query pattern
+
+The pruning rate depends on observer latitude, query window duration, and minimum elevation. Shorter windows and higher latitudes prune more aggressively.
+
+| Query | Candidates | Pruned | Notes |
+|-------|-----------|--------|-------|
+| 2h, Eagle ID (43.7°N), 10° | 12,964 | 80.3% | Typical mid-latitude evening |
+| 2h, Tromsoe (69.6°N), 10° | 6,529 | 90.1% | High latitude: inclination filter strongest |
+| 2h, South Pole (85°S), 10° | 5,248 | 92.0% | Only polar-orbit satellites survive |
+| 2h, Equator (0°N), 10° | 9,699 | 85.3% | All inclinations pass latitude check; RAAN filter dominates |
+| 2h, Eagle ID, 30° | 9,680 | 85.3% | Higher elevation: altitude filter tighter |
+| 6h, Eagle ID, 10° | 24,274 | 63.2% | Wider RAAN window admits more candidates |
+| 24h, Eagle ID, 10° | 60,875 | 7.6% | RAAN filter bypassed (full Earth rotation) |
+
+### SP-GiST index performance
+
+The optional SP-GiST index (`tle_spgist_ops`) builds a 2-level trie partitioned by semi-major axis and inclination. At 66,000 objects, the index adds 1--2 ms overhead compared to a sequential scan for mid-latitude observers, but **beats the sequential scan for high-latitude queries** where inclination pruning eliminates entire subtrees:
+
+| Query | Seqscan | SP-GiST | Difference |
+|-------|---------|---------|------------|
+| 2h, Eagle ID, 10° | 12.5 ms | 14.0 ms | +1.5 ms |
+| 2h, Tromsoe, 10° | 11.3 ms | 10.9 ms | **-0.4 ms** |
+
+The planner chooses the SP-GiST Index Only Scan by default at this catalog size, with zero heap fetches (all data served from index pages).
+
+<Aside type="tip" title="Where the index shines">
+The SP-GiST index is most valuable for high-latitude observers (60°+) and for catalogs larger than 30,000 objects. At typical CelesTrak catalog sizes (12--15,000 active satellites), the `&?` operator's sequential evaluation is fast enough that the index overhead exceeds the pruning benefit.
+</Aside>
+
+### What the pruning means for predict_passes()
+
+For a 65,886-object catalog and a 2-hour window from Eagle, Idaho:
+
+- **Without `&?`:** 65,886 `predict_passes()` calls (each ~1 ms for a 7-day window)
+- **With `&?`:** 12,964 calls --- **52,922 unnecessary propagations avoided**
+- **Time saved:** ~53 seconds per query at typical propagation cost
+
 ## Reproducing these benchmarks
 
 <Tabs>
   <TabItem label="Requirements">
     - PostgreSQL 17 with pg_orrery installed
-    - A satellite catalog table with ~12,000 TLEs (available from CelesTrak)
+    - A satellite catalog table with ~12,000 TLEs (see [Building TLE Catalogs](/guides/catalog-management/) or download directly from CelesTrak)
     - A star catalog table (any subset of Hipparcos or Yale BSC)
     - No concurrent queries during measurement
     - `shared_buffers` and `work_mem` at default or higher
@@ -233,9 +293,10 @@ A 7-day window at 30-second coarse scan resolution requires ~20,160 propagation
     ```sql
     CREATE EXTENSION pg_orrery;
 
-    -- Load a TLE catalog
+    -- Load a TLE catalog (pg-orrery-catalog handles this)
+    -- pg-orrery-catalog build --table satellite_catalog | psql -d mydb
     CREATE TABLE satellite_catalog (tle tle);
-    -- (COPY from CelesTrak bulk TLE file)
+    -- (or COPY from CelesTrak bulk TLE file)
 
     -- Verify catalog size
     SELECT count(*) FROM satellite_catalog;
@@ -265,4 +326,6 @@ A 7-day window at 30-second coarse scan resolution requires ~20,160 propagation
 
 The benchmarks demonstrate that pg_orrery's computation cost is low enough to treat orbital mechanics as a SQL primitive. Propagating an entire satellite catalog takes less time than a typical index scan on a moderately-sized table. Planet observation is fast enough to generate ephemeris tables with `generate_series`. Pork chop plots are feasible as interactive queries rather than batch jobs.
 
-The numbers also show where the bottlenecks are: VSOP87 series evaluation dominates everything except star observation and raw SGP4 propagation. If a future optimization effort targets one component, it should be the VSOP87 evaluation loop.
+The visibility cone filter (`&?`) is the fastest operation per evaluation --- three floating-point comparisons vs. the full SGP4 pipeline --- and its 80--92% pruning rate means the most expensive operation in a pass prediction pipeline (SGP4 propagation) only runs on the small fraction of the catalog that could actually produce a visible pass.
+
+The numbers also show where the bottlenecks are: VSOP87 series evaluation dominates everything except star observation, raw SGP4 propagation, and the visibility cone filter. If a future optimization effort targets one component, it should be the VSOP87 evaluation loop.

docs/src/content/docs/reference/operators-gist.mdx

@@ -1,12 +1,12 @@
 ---
-title: "Operators & GiST Index"
+title: "Operators & Indexes"
 sidebar:
   order: 8
 ---
 
 import { Aside, Tabs, TabItem } from "@astrojs/starlight/components";
 
-pg_orrery defines two operators on the `tle` type and a GiST operator class that enables indexed conjunction screening over large satellite catalogs. The operators work on the orbital altitude band and inclination range extracted from TLE elements, providing a fast necessary-condition filter for proximity analysis.
+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).
 
 ---
 
@@ -188,3 +188,116 @@ TRUNCATE satellite_catalog;
 COPY satellite_catalog FROM '/path/to/catalog.csv' WITH (FORMAT csv);
 REINDEX INDEX idx_tle_gist;
 ```
+
+---
+
+### &? (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).
+
+#### Signature
+
+```sql
+tle &? observer_window → boolean
+```
+
+#### The `observer_window` Type
+
+The right argument is a composite type constructed with `ROW(...)::observer_window`:
+
+| Field | Type | Description |
+|---|---|---|
+| `obs` | `observer` | Ground location (latitude, longitude, altitude) |
+| `t_start` | `timestamptz` | Start of observation window |
+| `t_end` | `timestamptz` | End of observation window |
+| `min_el` | `float8` | Minimum elevation angle in degrees (default 10.0 if NULL) |
+
+#### Description
+
+Applies three geometric filters without SGP4 propagation:
+
+1. **Altitude filter:** Rejects satellites whose perigee altitude exceeds the maximum visible altitude for the given minimum elevation angle
+2. **Inclination filter:** Rejects satellites whose inclination + ground footprint angle cannot reach the observer's latitude
+3. **RAAN filter:** Projects the ascending node to the query midpoint via J2 secular precession and checks alignment with the observer's local sidereal time. Automatically bypassed for query windows spanning a full Earth rotation (>= ~12 hours)
+
+Returns `true` if the satellite passes all three filters. Returns `false` for degenerate TLEs with zero mean motion.
+
+<Aside type="note">
+The `&?` operator is designed as the first stage of a two-stage pipeline. Use it to reduce a 30,000-object catalog to a few hundred candidates, then run `predict_passes()` on the survivors for exact pass times and elevations.
+</Aside>
+
+#### Example
+
+```sql
+-- Which satellites might be visible from Eagle, Idaho tonight?
+SELECT name
+FROM satellite_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 02:00:00+00'::timestamptz,
+    '2024-01-01 08:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+```
+
+---
+
+## SP-GiST Operator Class: tle_spgist_ops
+
+The `tle_spgist_ops` operator class enables SP-GiST indexing on `tle` columns. The index builds a 2-level space-partitioning trie: semi-major axis at level 0 (altitude regime) and inclination at level 1 (latitude coverage). Equal-population splits ensure balanced trees across the dense LEO region.
+
+### Creating the Index
+
+```sql
+CREATE INDEX idx_tle_spgist ON satellite_catalog USING spgist (tle tle_spgist_ops);
+```
+
+<Aside type="note">
+The SP-GiST operator class is **not** the default for the `tle` type. The GiST operator class (`tle_ops`) remains the default. Both can coexist on the same table -- they serve different query patterns.
+</Aside>
+
+### What Gets Indexed
+
+The SP-GiST trie partitions TLEs by two static orbital elements:
+
+- **Level 0: Semi-major axis** (computed from mean motion via Kepler's 3rd law). Separates LEO, MEO, HEO, and GEO objects into altitude bins.
+- **Level 1: Inclination** (in radians). Within each altitude bin, objects are further partitioned by orbital inclination.
+
+Equal-population splits (`floor(sqrt(n))` bins, clamped to [2, 16]) ensure dense orbital regimes like LEO get finer partitioning.
+
+<Aside type="note" title="Eccentricity and highly elliptical orbits">
+The SP-GiST inner nodes store only SMA bin boundaries — they do not carry eccentricity information. This means the index cannot prune HEO (Highly Elliptical Orbit) satellites by altitude at the tree level. A satellite like CLUSTER II (SMA ~70,000 km, eccentricity 0.88) sits in a high-SMA bin alongside GEO satellites, but its actual perigee is only ~2,000 km. The leaf-level filter correctly identifies these using the full TLE (including eccentricity), so **no false negatives occur** — but HEO bins may produce slightly more candidates than a purely circular-orbit catalog would.
+</Aside>
+
+### Index-Accelerated Queries
+
+```sql
+-- The &? operator uses the SP-GiST index when available
+SELECT name
+FROM satellite_catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 02:00:00+00'::timestamptz,
+    '2024-01-01 08:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+```
+
+During the index scan, inner nodes are pruned by altitude band (level 0) and inclination range (level 1). The RAAN filter is applied at the leaf level. This avoids examining individual TLEs in entire subtrees that cannot produce visible passes.
+
+### Performance
+
+The `&?` operator eliminates 80-90% of a satellite catalog without SGP4 propagation --- this is the primary value, regardless of whether a sequential scan or index scan evaluates it. At typical catalog sizes (10-30k objects), the operator evaluates the full catalog in under 10 ms, and PostgreSQL's query planner may choose a sequential scan over the index.
+
+The SP-GiST index becomes advantageous at larger catalog sizes (100k+ objects) where tree-level pruning avoids examining individual TLEs in entire subtrees. The index is most effective for:
+
+- **Short query windows** (1-6 hours): The RAAN filter aggressively eliminates satellites whose orbital plane is not currently aligned with the observer
+- **Mid-latitude observers** (30-60 degrees): The inclination filter eliminates equatorial and low-inclination satellites
+- **High minimum elevation** (> 20 degrees): The altitude filter eliminates distant MEO/GEO objects
+
+For 24-hour query windows, the RAAN filter self-disables (full Earth rotation makes it meaningless), and only the altitude and inclination filters apply.
+
+### Index Maintenance
+
+Like the GiST index, the SP-GiST index is maintained automatically by PostgreSQL on `INSERT`, `UPDATE`, and `DELETE`. No manual `REINDEX` is needed under normal operation.

docs/src/content/docs/reference/types.mdx

@@ -6,7 +6,7 @@ sidebar:
 
 import { Aside, Tabs, TabItem } from "@astrojs/starlight/components";
 
-pg_orrery defines seven composite types that represent the core data structures of orbital mechanics. Each type has a fixed on-disk size, a text I/O format for readability, and accessor functions for extracting individual fields.
+pg_orrery defines seven fixed-size base types and one SQL composite type that represent the core data structures of orbital mechanics. Each base type has a fixed on-disk size, a text I/O format for readability, and accessor functions for extracting individual fields.
 
 ## tle
 
@@ -267,3 +267,40 @@ SELECT body_id,
 FROM generate_series(1, 8) AS body_id,
      planet_heliocentric(body_id, now()) AS h;
 ```
+
+---
+
+## observer_window
+
+**Type:** SQL composite (variable-length)
+
+A query-time parameter that defines a ground observer's visibility window. Unlike the seven base types above, `observer_window` is a SQL composite type --- it is not designed for table storage, but as an argument to the `&?` visibility cone operator.
+
+### Fields
+
+| Field | Type | Description |
+|-------|------|-------------|
+| `obs` | `observer` | Ground location (latitude, longitude, altitude) |
+| `t_start` | `timestamptz` | Start of observation window |
+| `t_end` | `timestamptz` | End of observation window |
+| `min_el` | `float8` | Minimum elevation angle in degrees |
+
+### Construction
+
+Construct with `ROW(...)::observer_window` syntax:
+
+```sql
+-- Eagle, Idaho: 6-hour window, 10 deg minimum elevation
+SELECT tle
+FROM catalog
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 01:00:00+00'::timestamptz,
+    '2024-01-01 07:00:00+00'::timestamptz,
+    10.0
+)::observer_window;
+```
+
+<Aside type="tip">
+See [Satellite Pass Prediction](/guides/pass-prediction/) for the full two-stage workflow using `&?` to filter candidates before SGP4 propagation with `predict_passes()`.
+</Aside>

pg_orrery.control

@@ -1,4 +1,4 @@
 comment = 'A database orrery — celestial mechanics types and functions for PostgreSQL'
-default_version = '0.6.0'
+default_version = '0.7.0'
 module_pathname = '$libdir/pg_orrery'
 relocatable = true

sql/pg_orrery--0.6.0--0.7.0.sql

@@ -0,0 +1,78 @@
+-- pg_orrery 0.6.0 -> 0.7.0 migration
+--
+-- Adds SP-GiST orbital trie index for satellite pass prediction.
+-- 2-level trie: SMA (L0) + inclination (L1) with query-time RAAN filter.
+-- The &? operator answers "might this satellite be visible?"
+
+-- ============================================================
+-- observer_window composite type (query parameter bundle)
+-- ============================================================
+
+CREATE TYPE observer_window AS (
+    obs         observer,
+    t_start     timestamptz,
+    t_end       timestamptz,
+    min_el      float8
+);
+
+COMMENT ON TYPE observer_window IS
+    'Observation query parameters: observer location, time window, and minimum elevation angle (degrees). Used with the &? visibility cone operator.';
+
+-- ============================================================
+-- Visibility cone operator function
+-- ============================================================
+
+CREATE FUNCTION tle_visibility_possible(tle, observer_window) RETURNS boolean
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+
+COMMENT ON FUNCTION tle_visibility_possible(tle, observer_window) IS
+    'Could this satellite be visible from the observer during the time window? Combines altitude, inclination, and RAAN checks. Conservative superset — survivors need SGP4 propagation for ground truth.';
+
+-- ============================================================
+-- &? operator (visibility cone check)
+-- ============================================================
+-- The indexed column (tle) MUST be the left argument so PostgreSQL
+-- can form a ScanKey and pass it to inner_consistent for pruning.
+
+CREATE OPERATOR &? (
+    LEFTARG   = tle,
+    RIGHTARG  = observer_window,
+    FUNCTION  = tle_visibility_possible,
+    RESTRICT  = contsel,
+    JOIN      = contjoinsel
+);
+
+COMMENT ON OPERATOR &? (tle, observer_window) IS
+    'Visibility cone check: could this satellite be visible from the observer during the time window? Index-accelerated via SP-GiST orbital trie.';
+
+-- ============================================================
+-- SP-GiST support functions
+-- ============================================================
+
+CREATE FUNCTION spgist_tle_config(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_choose(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_picksplit(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_inner_consistent(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_leaf_consistent(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+-- ============================================================
+-- SP-GiST operator class (opt-in, not DEFAULT)
+-- ============================================================
+
+CREATE OPERATOR CLASS tle_spgist_ops
+    FOR TYPE tle USING spgist AS
+        OPERATOR  1   &? (tle, observer_window),
+        FUNCTION  1   spgist_tle_config(internal, internal),
+        FUNCTION  2   spgist_tle_choose(internal, internal),
+        FUNCTION  3   spgist_tle_picksplit(internal, internal),
+        FUNCTION  4   spgist_tle_inner_consistent(internal, internal),
+        FUNCTION  5   spgist_tle_leaf_consistent(internal, internal);

sql/pg_orrery--0.7.0.sql

@@ -0,0 +1,963 @@
+-- pg_orrery -- Orbital mechanics types and functions for PostgreSQL
+--
+-- Types: tle, eci_position, geodetic, topocentric, observer, pass_event
+-- Provides SGP4/SDP4 propagation, coordinate transforms, pass prediction,
+-- and GiST indexing on altitude bands for conjunction screening.
+--
+-- All propagation uses WGS-72 constants (matching TLE mean element fitting).
+-- Coordinate output uses WGS-84 (matching modern geodetic standards).
+
+-- ============================================================
+-- Shell types (forward declarations)
+-- ============================================================
+
+CREATE TYPE tle;
+CREATE TYPE eci_position;
+CREATE TYPE geodetic;
+CREATE TYPE topocentric;
+CREATE TYPE observer;
+CREATE TYPE pass_event;
+
+
+-- ============================================================
+-- TLE type: Two-Line Element set
+-- ============================================================
+
+CREATE FUNCTION tle_in(cstring) RETURNS tle
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION tle_out(tle) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION tle_recv(internal) RETURNS tle
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION tle_send(tle) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE tle (
+    INPUT = tle_in,
+    OUTPUT = tle_out,
+    RECEIVE = tle_recv,
+    SEND = tle_send,
+    INTERNALLENGTH = 112,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE tle IS 'Two-Line Element set — parsed mean orbital elements for SGP4/SDP4 propagation';
+
+-- TLE accessor functions
+
+CREATE FUNCTION tle_epoch(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_epoch(tle) IS 'TLE epoch as Julian date (UTC)';
+
+CREATE FUNCTION tle_norad_id(tle) RETURNS int4
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_norad_id(tle) IS 'NORAD catalog number';
+
+CREATE FUNCTION tle_inclination(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_inclination(tle) IS 'Orbital inclination in degrees';
+
+CREATE FUNCTION tle_eccentricity(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_eccentricity(tle) IS 'Orbital eccentricity (dimensionless)';
+
+CREATE FUNCTION tle_raan(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_raan(tle) IS 'Right ascension of ascending node in degrees';
+
+CREATE FUNCTION tle_arg_perigee(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_arg_perigee(tle) IS 'Argument of perigee in degrees';
+
+CREATE FUNCTION tle_mean_anomaly(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_mean_anomaly(tle) IS 'Mean anomaly in degrees';
+
+CREATE FUNCTION tle_mean_motion(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_mean_motion(tle) IS 'Mean motion in revolutions per day';
+
+CREATE FUNCTION tle_bstar(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_bstar(tle) IS 'B* drag coefficient (1/earth-radii)';
+
+CREATE FUNCTION tle_period(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_period(tle) IS 'Orbital period in minutes';
+
+CREATE FUNCTION tle_age(tle, timestamptz) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_age(tle, timestamptz) IS 'TLE age in days (positive = past epoch, negative = before epoch)';
+
+CREATE FUNCTION tle_perigee(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_perigee(tle) IS 'Perigee altitude in km above WGS-72 ellipsoid';
+
+CREATE FUNCTION tle_apogee(tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_apogee(tle) IS 'Apogee altitude in km above WGS-72 ellipsoid';
+
+CREATE FUNCTION tle_intl_desig(tle) RETURNS text
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_intl_desig(tle) IS 'International designator (COSPAR ID)';
+
+CREATE FUNCTION tle_from_lines(text, text) RETURNS tle
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_from_lines(text, text) IS
+    'Construct TLE from separate line1/line2 text columns';
+
+
+-- ============================================================
+-- ECI position type: True Equator Mean Equinox (TEME) frame
+-- ============================================================
+
+CREATE FUNCTION eci_in(cstring) RETURNS eci_position
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_out(eci_position) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_recv(internal) RETURNS eci_position
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_send(eci_position) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE eci_position (
+    INPUT = eci_in,
+    OUTPUT = eci_out,
+    RECEIVE = eci_recv,
+    SEND = eci_send,
+    INTERNALLENGTH = 48,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE eci_position IS 'Earth-Centered Inertial position and velocity in TEME frame (km, km/s)';
+
+-- ECI accessor functions
+
+CREATE FUNCTION eci_x(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_y(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_z(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_vx(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_vy(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_vz(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION eci_speed(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION eci_speed(eci_position) IS 'Velocity magnitude in km/s';
+
+CREATE FUNCTION eci_altitude(eci_position) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION eci_altitude(eci_position) IS 'Approximate geocentric altitude in km (radius - WGS72_AE)';
+
+
+-- ============================================================
+-- Geodetic type: WGS-84 latitude/longitude/altitude
+-- ============================================================
+
+CREATE FUNCTION geodetic_in(cstring) RETURNS geodetic
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION geodetic_out(geodetic) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION geodetic_recv(internal) RETURNS geodetic
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION geodetic_send(geodetic) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE geodetic (
+    INPUT = geodetic_in,
+    OUTPUT = geodetic_out,
+    RECEIVE = geodetic_recv,
+    SEND = geodetic_send,
+    INTERNALLENGTH = 24,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE geodetic IS 'Geodetic coordinates on WGS-84 ellipsoid (lat/lon in degrees, altitude in km)';
+
+CREATE FUNCTION geodetic_lat(geodetic) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION geodetic_lon(geodetic) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION geodetic_alt(geodetic) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+
+-- ============================================================
+-- Topocentric type: observer-relative az/el/range
+-- ============================================================
+
+CREATE FUNCTION topocentric_in(cstring) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION topocentric_out(topocentric) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION topocentric_recv(internal) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION topocentric_send(topocentric) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE topocentric (
+    INPUT = topocentric_in,
+    OUTPUT = topocentric_out,
+    RECEIVE = topocentric_recv,
+    SEND = topocentric_send,
+    INTERNALLENGTH = 32,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE topocentric IS 'Topocentric coordinates relative to observer (azimuth, elevation, range, range rate)';
+
+CREATE FUNCTION topo_azimuth(topocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION topo_azimuth(topocentric) IS 'Azimuth in degrees (0=N, 90=E, 180=S, 270=W)';
+
+CREATE FUNCTION topo_elevation(topocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION topo_elevation(topocentric) IS 'Elevation in degrees (0=horizon, 90=zenith)';
+
+CREATE FUNCTION topo_range(topocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION topo_range(topocentric) IS 'Slant range in km';
+
+CREATE FUNCTION topo_range_rate(topocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION topo_range_rate(topocentric) IS 'Range rate in km/s (positive = receding)';
+
+
+-- ============================================================
+-- Observer type: ground station location
+-- ============================================================
+
+CREATE FUNCTION observer_in(cstring) RETURNS observer
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION observer_out(observer) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION observer_recv(internal) RETURNS observer
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION observer_send(observer) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE observer (
+    INPUT = observer_in,
+    OUTPUT = observer_out,
+    RECEIVE = observer_recv,
+    SEND = observer_send,
+    INTERNALLENGTH = 24,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE observer IS 'Ground observer location (accepts: 40.0N 105.3W 1655m or decimal degrees)';
+
+CREATE FUNCTION observer_lat(observer) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION observer_lat(observer) IS 'Latitude in degrees (positive = North)';
+
+CREATE FUNCTION observer_lon(observer) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION observer_lon(observer) IS 'Longitude in degrees (positive = East)';
+
+CREATE FUNCTION observer_alt(observer) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION observer_alt(observer) IS 'Altitude in meters above WGS-84 ellipsoid';
+
+CREATE FUNCTION observer_from_geodetic(float8, float8, float8 DEFAULT 0.0) RETURNS observer
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION observer_from_geodetic(float8, float8, float8) IS
+    'Construct observer from lat (deg), lon (deg), altitude (meters). Avoids text formatting round-trips.';
+
+
+-- ============================================================
+-- Pass event type: satellite visibility window
+-- ============================================================
+
+CREATE FUNCTION pass_event_in(cstring) RETURNS pass_event
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION pass_event_out(pass_event) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION pass_event_recv(internal) RETURNS pass_event
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION pass_event_send(pass_event) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE pass_event (
+    INPUT = pass_event_in,
+    OUTPUT = pass_event_out,
+    RECEIVE = pass_event_recv,
+    SEND = pass_event_send,
+    INTERNALLENGTH = 48,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE pass_event IS 'Satellite pass event (AOS/MAX/LOS times, max elevation, AOS/LOS azimuths)';
+
+CREATE FUNCTION pass_aos_time(pass_event) RETURNS timestamptz
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_aos_time(pass_event) IS 'Acquisition of signal time';
+
+CREATE FUNCTION pass_max_el_time(pass_event) RETURNS timestamptz
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_max_el_time(pass_event) IS 'Maximum elevation time';
+
+CREATE FUNCTION pass_los_time(pass_event) RETURNS timestamptz
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_los_time(pass_event) IS 'Loss of signal time';
+
+CREATE FUNCTION pass_max_elevation(pass_event) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_max_elevation(pass_event) IS 'Maximum elevation in degrees';
+
+CREATE FUNCTION pass_aos_azimuth(pass_event) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_aos_azimuth(pass_event) IS 'AOS azimuth in degrees (0=N)';
+
+CREATE FUNCTION pass_los_azimuth(pass_event) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_los_azimuth(pass_event) IS 'LOS azimuth in degrees (0=N)';
+
+CREATE FUNCTION pass_duration(pass_event) RETURNS interval
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_duration(pass_event) IS 'Pass duration (LOS - AOS)';
+
+
+-- ============================================================
+-- SGP4/SDP4 propagation functions
+-- ============================================================
+
+CREATE FUNCTION sgp4_propagate(tle, timestamptz) RETURNS eci_position
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION sgp4_propagate(tle, timestamptz) IS
+    'Propagate TLE to a point in time using SGP4 (near-earth) or SDP4 (deep-space). Returns TEME ECI position/velocity.';
+
+CREATE FUNCTION sgp4_propagate_safe(tle, timestamptz) RETURNS eci_position
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE PARALLEL SAFE;
+COMMENT ON FUNCTION sgp4_propagate_safe(tle, timestamptz) IS
+    'Like sgp4_propagate but returns NULL on error instead of raising an exception. For batch queries with potentially invalid TLEs.';
+
+CREATE FUNCTION sgp4_propagate_series(tle, timestamptz, timestamptz, interval)
+    RETURNS TABLE(t timestamptz, x float8, y float8, z float8, vx float8, vy float8, vz float8)
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE
+    ROWS 100;
+COMMENT ON FUNCTION sgp4_propagate_series(tle, timestamptz, timestamptz, interval) IS
+    'Propagate TLE over a time range at regular intervals. Returns time series of TEME positions.';
+
+CREATE FUNCTION tle_distance(tle, tle, timestamptz) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_distance(tle, tle, timestamptz) IS
+    'Euclidean distance in km between two TLEs at a reference time';
+
+
+-- ============================================================
+-- Coordinate transform functions
+-- ============================================================
+
+CREATE FUNCTION eci_to_geodetic(eci_position, timestamptz) RETURNS geodetic
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION eci_to_geodetic(eci_position, timestamptz) IS
+    'Convert TEME ECI position to WGS-84 geodetic coordinates at given time';
+
+CREATE FUNCTION eci_to_topocentric(eci_position, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION eci_to_topocentric(eci_position, observer, timestamptz) IS
+    'Convert TEME ECI position to topocentric (az/el/range) relative to observer';
+
+CREATE FUNCTION subsatellite_point(tle, timestamptz) RETURNS geodetic
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION subsatellite_point(tle, timestamptz) IS
+    'Subsatellite (nadir) point on WGS-84 ellipsoid at given time';
+
+CREATE FUNCTION ground_track(tle, timestamptz, timestamptz, interval)
+    RETURNS TABLE(t timestamptz, lat float8, lon float8, alt float8)
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE
+    ROWS 100;
+COMMENT ON FUNCTION ground_track(tle, timestamptz, timestamptz, interval) IS
+    'Ground track as time series of subsatellite points (lat/lon in degrees, alt in km)';
+
+CREATE FUNCTION observe(tle, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION observe(tle, observer, timestamptz) IS
+    'Propagate TLE and compute observer-relative look angles in one call. Returns topocentric (az/el/range/range_rate).';
+
+CREATE FUNCTION observe_safe(tle, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE PARALLEL SAFE;
+COMMENT ON FUNCTION observe_safe(tle, observer, timestamptz) IS
+    'Like observe() but returns NULL on propagation error. For batch queries with potentially invalid/decayed TLEs.';
+
+
+-- ============================================================
+-- Pass prediction functions
+-- ============================================================
+
+CREATE FUNCTION next_pass(tle, observer, timestamptz) RETURNS pass_event
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION next_pass(tle, observer, timestamptz) IS
+    'Find the next satellite pass over observer (searches up to 7 days ahead)';
+
+CREATE FUNCTION predict_passes(tle, observer, timestamptz, timestamptz, float8 DEFAULT 0.0)
+    RETURNS SETOF pass_event
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE
+    ROWS 10;
+COMMENT ON FUNCTION predict_passes(tle, observer, timestamptz, timestamptz, float8) IS
+    'Predict all satellite passes over observer in time window. Optional min_elevation in degrees.';
+
+CREATE FUNCTION pass_visible(tle, observer, timestamptz, timestamptz) RETURNS boolean
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION pass_visible(tle, observer, timestamptz, timestamptz) IS
+    'True if any pass occurs over observer in the time window';
+
+
+-- ============================================================
+-- GiST operator support functions
+-- ============================================================
+
+-- Overlap operator: do orbital keys overlap in altitude AND inclination?
+CREATE FUNCTION tle_overlap(tle, tle) RETURNS boolean
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+-- Altitude distance operator (altitude-only, for KNN ordering)
+CREATE FUNCTION tle_alt_distance(tle, tle) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE OPERATOR && (
+    LEFTARG = tle,
+    RIGHTARG = tle,
+    FUNCTION = tle_overlap,
+    COMMUTATOR = &&,
+    RESTRICT = areasel,
+    JOIN = areajoinsel
+);
+
+COMMENT ON OPERATOR && (tle, tle) IS 'Orbital key overlap (altitude band AND inclination range) — necessary condition for conjunction';
+
+CREATE OPERATOR <-> (
+    LEFTARG = tle,
+    RIGHTARG = tle,
+    FUNCTION = tle_alt_distance,
+    COMMUTATOR = <->
+);
+
+COMMENT ON OPERATOR <-> (tle, tle) IS 'Minimum altitude-band separation in km (0 if overlapping). Altitude-only — does not account for inclination. Use && for 2-D filtering.';
+
+
+-- ============================================================
+-- GiST operator class for 2-D orbital indexing (altitude + inclination)
+-- ============================================================
+
+-- GiST internal support functions
+CREATE FUNCTION gist_tle_compress(internal) RETURNS internal
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_decompress(internal) RETURNS internal
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_consistent(internal, tle, smallint, oid, internal) RETURNS boolean
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_union(internal, internal) RETURNS internal
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_penalty(internal, internal, internal) RETURNS internal
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_picksplit(internal, internal) RETURNS internal
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_same(internal, internal, internal) RETURNS internal
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION gist_tle_distance(internal, tle, smallint, oid, internal) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE OPERATOR CLASS tle_ops
+    DEFAULT FOR TYPE tle USING gist AS
+        OPERATOR  3   && ,
+        OPERATOR  15  <-> (tle, tle) FOR ORDER BY float_ops,
+        FUNCTION  1   gist_tle_consistent(internal, tle, smallint, oid, internal),
+        FUNCTION  2   gist_tle_union(internal, internal),
+        FUNCTION  3   gist_tle_compress(internal),
+        FUNCTION  4   gist_tle_decompress(internal),
+        FUNCTION  5   gist_tle_penalty(internal, internal, internal),
+        FUNCTION  6   gist_tle_picksplit(internal, internal),
+        FUNCTION  7   gist_tle_same(internal, internal, internal),
+        FUNCTION  8   gist_tle_distance(internal, tle, smallint, oid, internal);
+
+
+-- ============================================================
+-- Heliocentric type: ecliptic J2000 position in AU
+-- ============================================================
+
+CREATE TYPE heliocentric;
+
+CREATE FUNCTION heliocentric_in(cstring) RETURNS heliocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION heliocentric_out(heliocentric) RETURNS cstring
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION heliocentric_recv(internal) RETURNS heliocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION heliocentric_send(heliocentric) RETURNS bytea
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE TYPE heliocentric (
+    INPUT = heliocentric_in,
+    OUTPUT = heliocentric_out,
+    RECEIVE = heliocentric_recv,
+    SEND = heliocentric_send,
+    INTERNALLENGTH = 24,
+    ALIGNMENT = double,
+    STORAGE = plain
+);
+
+COMMENT ON TYPE heliocentric IS 'Heliocentric position in ecliptic J2000 frame (x, y, z in AU)';
+
+CREATE FUNCTION helio_x(heliocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION helio_x(heliocentric) IS 'X component in AU (ecliptic J2000, vernal equinox direction)';
+
+CREATE FUNCTION helio_y(heliocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION helio_y(heliocentric) IS 'Y component in AU (ecliptic J2000)';
+
+CREATE FUNCTION helio_z(heliocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION helio_z(heliocentric) IS 'Z component in AU (ecliptic J2000, north ecliptic pole)';
+
+CREATE FUNCTION helio_distance(heliocentric) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION helio_distance(heliocentric) IS 'Heliocentric distance in AU';
+
+
+-- ============================================================
+-- Star observation functions
+-- ============================================================
+
+CREATE FUNCTION star_observe(float8, float8, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION star_observe(float8, float8, observer, timestamptz) IS
+    'Observe a star from (ra_hours J2000, dec_degrees J2000, observer, time). Returns topocentric az/el. Range is 0 (infinite distance).';
+
+CREATE FUNCTION star_observe_safe(float8, float8, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE PARALLEL SAFE;
+COMMENT ON FUNCTION star_observe_safe(float8, float8, observer, timestamptz) IS
+    'Like star_observe but returns NULL on invalid inputs. For batch queries over star catalogs.';
+
+
+-- ============================================================
+-- Keplerian propagation functions
+-- ============================================================
+
+CREATE FUNCTION kepler_propagate(
+    q_au float8, eccentricity float8,
+    inclination_deg float8, arg_perihelion_deg float8,
+    long_asc_node_deg float8, perihelion_jd float8,
+    t timestamptz
+) RETURNS heliocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION kepler_propagate(float8, float8, float8, float8, float8, float8, timestamptz) IS
+    'Two-body Keplerian propagation from classical orbital elements. Returns heliocentric ecliptic J2000 position in AU. Handles elliptic, parabolic, and hyperbolic orbits.';
+
+
+-- ============================================================
+-- Comet observation
+-- ============================================================
+
+CREATE FUNCTION comet_observe(
+    q_au float8, eccentricity float8,
+    inclination_deg float8, arg_perihelion_deg float8,
+    long_asc_node_deg float8, perihelion_jd float8,
+    earth_x_au float8, earth_y_au float8, earth_z_au float8,
+    obs observer, t timestamptz
+) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION comet_observe(float8, float8, float8, float8, float8, float8, float8, float8, float8, observer, timestamptz) IS
+    'Observe a comet/asteroid from orbital elements. Requires Earth heliocentric ecliptic J2000 position (AU). Returns topocentric az/el with geocentric range in km.';
+
+
+-- ============================================================
+-- VSOP87 planets, ELP82B Moon, Sun observation
+-- ============================================================
+
+CREATE FUNCTION planet_heliocentric(int4, timestamptz) RETURNS heliocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION planet_heliocentric(int4, timestamptz) IS
+    'VSOP87 heliocentric ecliptic J2000 position (AU). Body IDs: 0=Sun, 1=Mercury, ..., 8=Neptune.';
+
+CREATE FUNCTION planet_observe(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION planet_observe(int4, observer, timestamptz) IS
+    'Observe a planet from (body_id 1-8, observer, time). Returns topocentric az/el with geocentric range in km.';
+
+CREATE FUNCTION sun_observe(observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION sun_observe(observer, timestamptz) IS
+    'Observe the Sun from (observer, time). Returns topocentric az/el with geocentric range in km.';
+
+CREATE FUNCTION moon_observe(observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION moon_observe(observer, timestamptz) IS
+    'Observe the Moon via ELP2000-82B from (observer, time). Returns topocentric az/el with geocentric range in km.';
+
+
+-- ============================================================
+-- Planetary moon observation
+-- ============================================================
+
+CREATE FUNCTION galilean_observe(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION galilean_observe(int4, observer, timestamptz) IS
+    'Observe a Galilean moon of Jupiter via L1.2 theory. Body IDs: 0=Io, 1=Europa, 2=Ganymede, 3=Callisto.';
+
+CREATE FUNCTION saturn_moon_observe(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION saturn_moon_observe(int4, observer, timestamptz) IS
+    'Observe a Saturn moon via TASS 1.7. Body IDs: 0=Mimas, 1=Enceladus, 2=Tethys, 3=Dione, 4=Rhea, 5=Titan, 6=Iapetus, 7=Hyperion.';
+
+CREATE FUNCTION uranus_moon_observe(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION uranus_moon_observe(int4, observer, timestamptz) IS
+    'Observe a Uranus moon via GUST86. Body IDs: 0=Miranda, 1=Ariel, 2=Umbriel, 3=Titania, 4=Oberon.';
+
+CREATE FUNCTION mars_moon_observe(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION mars_moon_observe(int4, observer, timestamptz) IS
+    'Observe a Mars moon via MarsSat. Body IDs: 0=Phobos, 1=Deimos.';
+
+
+-- ============================================================
+-- Jupiter decametric radio burst prediction
+-- ============================================================
+
+CREATE FUNCTION io_phase_angle(timestamptz) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION io_phase_angle(timestamptz) IS
+    'Io orbital phase angle in degrees [0,360). 0=superior conjunction (behind Jupiter). Standard Radio JOVE convention.';
+
+CREATE FUNCTION jupiter_cml(observer, timestamptz) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION jupiter_cml(observer, timestamptz) IS
+    'Jupiter Central Meridian Longitude, System III (1965.0), in degrees [0,360). Light-time corrected.';
+
+CREATE FUNCTION jupiter_burst_probability(float8, float8) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION jupiter_burst_probability(float8, float8) IS
+    'Estimated Jupiter decametric burst probability (0-1) from (io_phase_deg, cml_deg). Based on Carr et al. (1983) source regions A, B, C, D.';
+
+
+-- ============================================================
+-- Interplanetary transfer orbits (Lambert solver)
+-- ============================================================
+
+CREATE FUNCTION lambert_transfer(
+    dep_body_id int4, arr_body_id int4,
+    dep_time timestamptz, arr_time timestamptz,
+    OUT c3_departure float8, OUT c3_arrival float8,
+    OUT v_inf_departure float8, OUT v_inf_arrival float8,
+    OUT tof_days float8, OUT transfer_sma float8
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION lambert_transfer(int4, int4, timestamptz, timestamptz) IS
+    'Solve Lambert transfer between two planets. Returns C3 (km^2/s^2), v_infinity (km/s), TOF (days), SMA (AU). Body IDs 1-8.';
+
+CREATE FUNCTION lambert_c3(int4, int4, timestamptz, timestamptz) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION lambert_c3(int4, int4, timestamptz, timestamptz) IS
+    'Departure C3 (km^2/s^2) for a Lambert transfer. Returns NULL if solver fails. For pork chop plots.';
+
+
+-- ============================================================
+-- DE ephemeris functions (optional high-precision)
+-- ============================================================
+
+CREATE FUNCTION planet_heliocentric_de(int4, timestamptz) RETURNS heliocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION planet_heliocentric_de(int4, timestamptz) IS
+    'Heliocentric ecliptic J2000 position via JPL DE (sub-arcsecond). Falls back to VSOP87 if DE unavailable. Body IDs: 0=Sun, 1-8=Mercury-Neptune.';
+
+CREATE FUNCTION planet_observe_de(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION planet_observe_de(int4, observer, timestamptz) IS
+    'Observe planet via JPL DE. Falls back to VSOP87. Body IDs: 1-8 (Mercury-Neptune).';
+
+CREATE FUNCTION sun_observe_de(observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION sun_observe_de(observer, timestamptz) IS
+    'Observe Sun via JPL DE. Falls back to VSOP87.';
+
+CREATE FUNCTION moon_observe_de(observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION moon_observe_de(observer, timestamptz) IS
+    'Observe Moon via JPL DE. Falls back to ELP2000-82B.';
+
+CREATE FUNCTION lambert_transfer_de(
+    dep_body_id int4, arr_body_id int4,
+    dep_time timestamptz, arr_time timestamptz,
+    OUT c3_departure float8, OUT c3_arrival float8,
+    OUT v_inf_departure float8, OUT v_inf_arrival float8,
+    OUT tof_days float8, OUT transfer_sma float8
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION lambert_transfer_de(int4, int4, timestamptz, timestamptz) IS
+    'Lambert transfer via JPL DE positions. Falls back to VSOP87. Body IDs 1-8.';
+
+CREATE FUNCTION lambert_c3_de(int4, int4, timestamptz, timestamptz) RETURNS float8
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION lambert_c3_de(int4, int4, timestamptz, timestamptz) IS
+    'Departure C3 via JPL DE. Falls back to VSOP87. For pork chop plots.';
+
+CREATE FUNCTION galilean_observe_de(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION galilean_observe_de(int4, observer, timestamptz) IS
+    'Observe Galilean moon with JPL DE parent position. L1.2 moon offsets. Body IDs: 0-3 (Io-Callisto).';
+
+CREATE FUNCTION saturn_moon_observe_de(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION saturn_moon_observe_de(int4, observer, timestamptz) IS
+    'Observe Saturn moon with JPL DE parent position. TASS 1.7 moon offsets. Body IDs: 0-7 (Mimas-Hyperion).';
+
+CREATE FUNCTION uranus_moon_observe_de(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION uranus_moon_observe_de(int4, observer, timestamptz) IS
+    'Observe Uranus moon with JPL DE parent position. GUST86 moon offsets. Body IDs: 0-4 (Miranda-Oberon).';
+
+CREATE FUNCTION mars_moon_observe_de(int4, observer, timestamptz) RETURNS topocentric
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION mars_moon_observe_de(int4, observer, timestamptz) IS
+    'Observe Mars moon with JPL DE parent position. MarsSat moon offsets. Body IDs: 0-1 (Phobos-Deimos).';
+
+
+-- Diagnostic function
+
+CREATE FUNCTION pg_orrery_ephemeris_info(
+    OUT provider text, OUT file_path text,
+    OUT start_jd float8, OUT end_jd float8,
+    OUT version int4, OUT au_km float8
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE PARALLEL SAFE;
+COMMENT ON FUNCTION pg_orrery_ephemeris_info() IS
+    'Returns current ephemeris provider status: VSOP87 or JPL_DE with file path, JD range, version, and AU value.';
+
+
+-- ============================================================
+-- Orbit determination (TLE fitting from observations)
+-- ============================================================
+
+-- Fit TLE from ECI position/velocity ephemeris
+-- weights: per-observation weighting (NULL = uniform)
+
+CREATE FUNCTION tle_from_eci(
+    positions eci_position[], times timestamptz[],
+    seed tle DEFAULT NULL, fit_bstar boolean DEFAULT false,
+    max_iter int4 DEFAULT 15,
+    weights float8[] DEFAULT NULL,
+    OUT fitted_tle tle, OUT iterations int4,
+    OUT rms_final float8, OUT rms_initial float8, OUT status text,
+    OUT condition_number float8, OUT covariance float8[], OUT nstate int4
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE PARALLEL SAFE;
+COMMENT ON FUNCTION tle_from_eci(eci_position[], timestamptz[], tle, boolean, int4, float8[]) IS
+    'Fit a TLE from ECI position/velocity observations via differential correction. Optional per-observation weights for heterogeneous sensor fusion. Returns fitted TLE, RMS residuals, convergence status, condition number, and formal covariance matrix.';
+
+-- Fit TLE from topocentric observations (az/el/range) — single observer
+-- fit_range_rate: include range_rate as 4th residual component
+-- weights: per-observation weighting (NULL = uniform)
+
+CREATE FUNCTION tle_from_topocentric(
+    observations topocentric[], times timestamptz[],
+    obs observer,
+    seed tle DEFAULT NULL, fit_bstar boolean DEFAULT false,
+    max_iter int4 DEFAULT 15,
+    fit_range_rate boolean DEFAULT false,
+    weights float8[] DEFAULT NULL,
+    OUT fitted_tle tle, OUT iterations int4,
+    OUT rms_final float8, OUT rms_initial float8, OUT status text,
+    OUT condition_number float8, OUT covariance float8[], OUT nstate int4
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE PARALLEL SAFE;
+COMMENT ON FUNCTION tle_from_topocentric(topocentric[], timestamptz[], observer, tle, boolean, int4, boolean, float8[]) IS
+    'Fit a TLE from topocentric (az/el/range) observations via differential correction. Optional range_rate fitting and per-observation weights. Returns fitted TLE, RMS residuals, convergence status, condition number, and formal covariance matrix.';
+
+-- Fit TLE from topocentric observations — multiple observers
+
+CREATE FUNCTION tle_from_topocentric(
+    observations topocentric[], times timestamptz[],
+    observers observer[], observer_ids int4[],
+    seed tle DEFAULT NULL, fit_bstar boolean DEFAULT false,
+    max_iter int4 DEFAULT 15,
+    fit_range_rate boolean DEFAULT false,
+    weights float8[] DEFAULT NULL,
+    OUT fitted_tle tle, OUT iterations int4,
+    OUT rms_final float8, OUT rms_initial float8, OUT status text,
+    OUT condition_number float8, OUT covariance float8[], OUT nstate int4
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME', 'tle_from_topocentric_multi'
+    LANGUAGE C STABLE PARALLEL SAFE;
+COMMENT ON FUNCTION tle_from_topocentric(topocentric[], timestamptz[], observer[], int4[], tle, boolean, int4, boolean, float8[]) IS
+    'Fit a TLE from topocentric observations collected by multiple ground stations. observer_ids[i] indexes into observers[]. Optional range_rate fitting and per-observation weights.';
+
+-- Per-observation residuals diagnostic
+
+CREATE FUNCTION tle_fit_residuals(
+    fitted tle,
+    positions eci_position[],
+    times timestamptz[]
+) RETURNS TABLE (
+    t timestamptz,
+    dx_km float8,
+    dy_km float8,
+    dz_km float8,
+    pos_err_km float8
+)
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+COMMENT ON FUNCTION tle_fit_residuals(tle, eci_position[], timestamptz[]) IS
+    'Compute per-observation position residuals (km) between a TLE and ECI observations. Useful for fit quality assessment.';
+
+-- Fit TLE from RA/Dec observations — single observer
+-- Uses Gauss method for initial orbit determination when no seed is provided.
+-- RA in hours [0,24), Dec in degrees [-90,90] (matches star_observe convention).
+-- RMS output is in radians for angles-only mode.
+
+CREATE FUNCTION tle_from_angles(
+    ra_hours float8[], dec_degrees float8[],
+    times timestamptz[],
+    obs observer,
+    seed tle DEFAULT NULL, fit_bstar boolean DEFAULT false,
+    max_iter int4 DEFAULT 15,
+    weights float8[] DEFAULT NULL,
+    OUT fitted_tle tle, OUT iterations int4,
+    OUT rms_final float8, OUT rms_initial float8, OUT status text,
+    OUT condition_number float8, OUT covariance float8[], OUT nstate int4
+) RETURNS RECORD AS 'MODULE_PATHNAME' LANGUAGE C STABLE PARALLEL SAFE;
+COMMENT ON FUNCTION tle_from_angles(float8[], float8[], timestamptz[], observer, tle, boolean, int4, float8[]) IS
+    'Fit a TLE from angles-only (RA/Dec) observations via Gauss IOD + differential correction. RA in hours [0,24), Dec in degrees [-90,90]. RMS output in radians. Uses Gauss method for seed-free initial guess.';
+
+-- Fit TLE from RA/Dec observations — multiple observers
+
+CREATE FUNCTION tle_from_angles(
+    ra_hours float8[], dec_degrees float8[],
+    times timestamptz[],
+    observers observer[], observer_ids int4[],
+    seed tle DEFAULT NULL, fit_bstar boolean DEFAULT false,
+    max_iter int4 DEFAULT 15,
+    weights float8[] DEFAULT NULL,
+    OUT fitted_tle tle, OUT iterations int4,
+    OUT rms_final float8, OUT rms_initial float8, OUT status text,
+    OUT condition_number float8, OUT covariance float8[], OUT nstate int4
+) RETURNS RECORD
+    AS 'MODULE_PATHNAME', 'tle_from_angles_multi'
+    LANGUAGE C STABLE PARALLEL SAFE;
+COMMENT ON FUNCTION tle_from_angles(float8[], float8[], timestamptz[], observer[], int4[], tle, boolean, int4, float8[]) IS
+    'Fit a TLE from angles-only (RA/Dec) observations from multiple ground stations via Gauss IOD + differential correction.';
+-- pg_orrery 0.6.0 -> 0.7.0 migration
+--
+-- Adds SP-GiST orbital trie index for satellite pass prediction.
+-- 2-level trie: SMA (L0) + inclination (L1) with query-time RAAN filter.
+-- The &? operator answers "might this satellite be visible?"
+
+-- ============================================================
+-- observer_window composite type (query parameter bundle)
+-- ============================================================
+
+CREATE TYPE observer_window AS (
+    obs         observer,
+    t_start     timestamptz,
+    t_end       timestamptz,
+    min_el      float8
+);
+
+COMMENT ON TYPE observer_window IS
+    'Observation query parameters: observer location, time window, and minimum elevation angle (degrees). Used with the &? visibility cone operator.';
+
+-- ============================================================
+-- Visibility cone operator function
+-- ============================================================
+
+CREATE FUNCTION tle_visibility_possible(tle, observer_window) RETURNS boolean
+    AS 'MODULE_PATHNAME' LANGUAGE C STABLE STRICT PARALLEL SAFE;
+
+COMMENT ON FUNCTION tle_visibility_possible(tle, observer_window) IS
+    'Could this satellite be visible from the observer during the time window? Combines altitude, inclination, and RAAN checks. Conservative superset — survivors need SGP4 propagation for ground truth.';
+
+-- ============================================================
+-- &? operator (visibility cone check)
+-- ============================================================
+-- The indexed column (tle) MUST be the left argument so PostgreSQL
+-- can form a ScanKey and pass it to inner_consistent for pruning.
+
+CREATE OPERATOR &? (
+    LEFTARG   = tle,
+    RIGHTARG  = observer_window,
+    FUNCTION  = tle_visibility_possible,
+    RESTRICT  = contsel,
+    JOIN      = contjoinsel
+);
+
+COMMENT ON OPERATOR &? (tle, observer_window) IS
+    'Visibility cone check: could this satellite be visible from the observer during the time window? Index-accelerated via SP-GiST orbital trie.';
+
+-- ============================================================
+-- SP-GiST support functions
+-- ============================================================
+
+CREATE FUNCTION spgist_tle_config(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_choose(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_picksplit(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_inner_consistent(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+CREATE FUNCTION spgist_tle_leaf_consistent(internal, internal) RETURNS void
+    AS 'MODULE_PATHNAME' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
+
+-- ============================================================
+-- SP-GiST operator class (opt-in, not DEFAULT)
+-- ============================================================
+
+CREATE OPERATOR CLASS tle_spgist_ops
+    FOR TYPE tle USING spgist AS
+        OPERATOR  1   &? (tle, observer_window),
+        FUNCTION  1   spgist_tle_config(internal, internal),
+        FUNCTION  2   spgist_tle_choose(internal, internal),
+        FUNCTION  3   spgist_tle_picksplit(internal, internal),
+        FUNCTION  4   spgist_tle_inner_consistent(internal, internal),
+        FUNCTION  5   spgist_tle_leaf_consistent(internal, internal);

src/gist_tle.c

@@ -327,6 +327,10 @@ gist_tle_consistent(PG_FUNCTION_ARGS)
  * gist_tle_union -- compute 2-D bounding box for a set of entries
  *
  * The union is [min(alt_low), max(alt_high)] x [min(inc_low), max(inc_high)].
+ *
+ * GiST convention: entryvec->vector[] is 1-based (FirstOffsetNumber),
+ * vector[0] is unused.  entryvec->n includes the unused slot, so
+ * valid indices are 1 .. entryvec->n - 1.
  */
 Datum
 gist_tle_union(PG_FUNCTION_ARGS)
@@ -338,10 +342,10 @@ gist_tle_union(PG_FUNCTION_ARGS)
 	tle_orbital_key	*cur;
 
 	result = (tle_orbital_key *) palloc(sizeof(tle_orbital_key));
-	cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[0].key);
+	cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[FirstOffsetNumber].key);
 	*result = *cur;
 
-	for (i = 1; i < entryvec->n; i++)
+	for (i = FirstOffsetNumber + 1; i < entryvec->n; i++)
 	{
 		cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[i].key);
 		key_merge(result, cur);
@@ -413,35 +417,43 @@ picksplit_cmp(const void *a, const void *b)
  * Standard R-tree approach: compute spread in both dimensions, split
  * along whichever dimension has the greater spread.  This prevents
  * the tree from becoming biased toward one dimension.
+ *
+ * GiST convention: entryvec->vector[] is 1-based (FirstOffsetNumber),
+ * vector[0] is unused/uninitialized.  entryvec->n includes the unused
+ * slot, so the actual entry count is (entryvec->n - 1) and valid
+ * indices are FirstOffsetNumber .. entryvec->n - 1.  The OffsetNumbers
+ * placed into spl_left[] and spl_right[] must also be 1-based.
  */
 Datum
 gist_tle_picksplit(PG_FUNCTION_ARGS)
 {
 	GistEntryVector	*entryvec = (GistEntryVector *) PG_GETARG_POINTER(0);
 	GIST_SPLITVEC	*splitvec = (GIST_SPLITVEC *) PG_GETARG_POINTER(1);
-	int				 nentries = entryvec->n;
+	OffsetNumber	 maxoff = entryvec->n - 1;
+	int				 nentries = maxoff - FirstOffsetNumber + 1;
 	picksplit_item	*items;
 	tle_orbital_key	*left_union;
 	tle_orbital_key	*right_union;
 	tle_orbital_key	*cur;
 	int				 split_at;
 	int				 i;
+	OffsetNumber	 off;
 	double			 alt_min, alt_max, inc_min, inc_max;
 	double			 alt_spread, inc_spread;
 	bool			 split_on_alt;
 
 	/* First pass: compute spread in both dimensions */
-	cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[0].key);
+	cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[FirstOffsetNumber].key);
 	alt_min = (cur->alt_low + cur->alt_high) / 2.0;
 	alt_max = alt_min;
 	inc_min = (cur->inc_low + cur->inc_high) / 2.0;
 	inc_max = inc_min;
 
-	for (i = 1; i < nentries; i++)
+	for (off = FirstOffsetNumber + 1; off <= maxoff; off++)
 	{
 		double alt_mid, inc_mid;
 
-		cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[i].key);
+		cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[off].key);
 		alt_mid = (cur->alt_low + cur->alt_high) / 2.0;
 		inc_mid = (cur->inc_low + cur->inc_high) / 2.0;
 
@@ -462,10 +474,10 @@ gist_tle_picksplit(PG_FUNCTION_ARGS)
 
 	/* Second pass: compute sort values in the chosen dimension */
 	items = (picksplit_item *) palloc(sizeof(picksplit_item) * nentries);
-	for (i = 0; i < nentries; i++)
+	for (i = 0, off = FirstOffsetNumber; off <= maxoff; i++, off++)
 	{
-		cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[i].key);
-		items[i].index = i;
+		cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[off].key);
+		items[i].index = off;	/* store 1-based OffsetNumber directly */
 		if (split_on_alt)
 			items[i].sortval = (cur->alt_low + cur->alt_high) / 2.0;
 		else
@@ -476,7 +488,7 @@ gist_tle_picksplit(PG_FUNCTION_ARGS)
 
 	split_at = nentries / 2;
 
-	/* Allocate offset arrays (GiST uses OffsetNumber, 1-based) */
+	/* Allocate offset arrays */
 	splitvec->spl_left  = (OffsetNumber *) palloc(sizeof(OffsetNumber) * nentries);
 	splitvec->spl_right = (OffsetNumber *) palloc(sizeof(OffsetNumber) * nentries);
 	splitvec->spl_nleft  = 0;
@@ -497,21 +509,19 @@ gist_tle_picksplit(PG_FUNCTION_ARGS)
 
 	for (i = 0; i < nentries; i++)
 	{
-		int idx = items[i].index;
+		OffsetNumber idx = items[i].index;	/* already 1-based */
 
 		cur = (tle_orbital_key *) DatumGetPointer(
 				entryvec->vector[idx].key);
 
 		if (i < split_at)
 		{
-			splitvec->spl_left[splitvec->spl_nleft++] =
-				(OffsetNumber)(idx + 1);	/* 1-based */
+			splitvec->spl_left[splitvec->spl_nleft++] = idx;
 			key_merge(left_union, cur);
 		}
 		else
 		{
-			splitvec->spl_right[splitvec->spl_nright++] =
-				(OffsetNumber)(idx + 1);
+			splitvec->spl_right[splitvec->spl_nright++] = idx;
 			key_merge(right_union, cur);
 		}
 	}

src/spgist_tle.c

@@ -0,0 +1,814 @@
+/*
+ * spgist_tle.c -- SP-GiST operator class for orbital trie on TLE
+ *
+ * 2-level space-partitioning trie:
+ *   L0: Semi-major axis (altitude regime, from Kepler's 3rd law)
+ *   L1: Inclination (ground-track latitude coverage)
+ *
+ * Query-time RAAN filter at leaf level projects ascending node to
+ * the query midpoint via J2 secular precession and checks alignment
+ * with the observer's local sidereal position.
+ *
+ * The &? operator answers "could this satellite be visible from this
+ * observer during this time window?" -- a conservative superset of
+ * the actual answer. Survivors go through SGP4 propagation.
+ *
+ * Equal-population splits: picksplit sorts by the level's element
+ * and divides into floor(sqrt(n)) bins, clamped [2,16]. Dense LEO
+ * gets finer SMA bins than sparse MEO/GEO.
+ */
+
+#include "postgres.h"
+#include "fmgr.h"
+#include "access/spgist.h"
+#include "access/htup_details.h"
+#include "catalog/pg_type_d.h"
+#include "utils/builtins.h"
+#include "utils/timestamp.h"
+#include "executor/executor.h"
+#include "types.h"
+#include "astro_math.h"
+#include <math.h>
+#include <float.h>
+
+PG_FUNCTION_INFO_V1(spgist_tle_config);
+PG_FUNCTION_INFO_V1(spgist_tle_choose);
+PG_FUNCTION_INFO_V1(spgist_tle_picksplit);
+PG_FUNCTION_INFO_V1(spgist_tle_inner_consistent);
+PG_FUNCTION_INFO_V1(spgist_tle_leaf_consistent);
+PG_FUNCTION_INFO_V1(tle_visibility_possible);
+
+/* Max trie depth: L0 (SMA) + L1 (inclination) */
+#define SPGIST_TLE_MAX_LEVEL  2
+
+/* Earth angular rotation rate in radians/day */
+#define EARTH_ROT_RAD_PER_DAY  (2.0 * M_PI)
+
+/* Seconds per day */
+#define SECONDS_PER_DAY  86400.0
+
+/* Minutes per day */
+#define MINUTES_PER_DAY  1440.0
+
+
+/* ----------------------------------------------------------------
+ * Helper: semi-major axis in km from mean motion
+ *
+ * Kepler's 3rd law with WGS-72: a = (KE / n)^(2/3) * AE
+ * where n is in radians/minute (TLE internal units).
+ * ----------------------------------------------------------------
+ */
+static inline double
+tle_sma_km(const pg_tle *tle)
+{
+	double n = tle->mean_motion;
+
+	if (n <= 0.0)
+		return 0.0;
+	return pow(WGS72_KE / n, 2.0 / 3.0) * WGS72_AE;
+}
+
+
+/* ----------------------------------------------------------------
+ * Helper: perigee altitude in km above Earth's surface
+ * ----------------------------------------------------------------
+ */
+static inline double
+tle_perigee_alt_km(const pg_tle *tle)
+{
+	double a = tle_sma_km(tle);
+
+	return a * (1.0 - tle->eccentricity) - WGS72_AE;
+}
+
+
+/* ----------------------------------------------------------------
+ * Helper: apogee altitude in km above Earth's surface
+ * ----------------------------------------------------------------
+ */
+static inline double
+tle_apogee_alt_km(const pg_tle *tle)
+{
+	double a = tle_sma_km(tle);
+
+	return a * (1.0 + tle->eccentricity) - WGS72_AE;
+}
+
+
+/* ----------------------------------------------------------------
+ * Helper: maximum satellite altitude visible at a given min elevation
+ *
+ * At min_el degrees elevation, the observer can see a satellite at
+ * most this far above the surface.  Conservative upper bound using
+ * the Earth-center angle geometry:
+ *   h_max = Re * (1/cos(el) - 1) roughly, but for a safe upper
+ *   bound we use the slant range limit.
+ *
+ * For min_el = 10 deg, practical limit is ~2500 km for LEO passes.
+ * For min_el = 0 deg, theoretical limit extends to GEO+.
+ * We use a generous bound: if min_el < 5, return 50000 (no filter).
+ * Otherwise compute from geometry.
+ * ----------------------------------------------------------------
+ */
+static inline double
+max_visible_altitude_km(double min_el_deg)
+{
+	double el_rad, sin_el;
+	double rho, h_max;
+
+	/*
+	 * Below 5 deg elevation the horizon geometry allows visibility to
+	 * GEO+ altitudes.  Disable the altitude filter rather than compute
+	 * an impractically large bound.  50000 km exceeds GEO (35786 km).
+	 */
+	if (min_el_deg < 5.0)
+		return 50000.0;
+
+	el_rad = min_el_deg * DEG_TO_RAD;
+	sin_el = sin(el_rad);
+
+	/*
+	 * Maximum slant range rho where a satellite at altitude h is visible
+	 * at elevation el.  From the geometry:
+	 *   rho = Re * (sqrt((h/Re + 1)^2 - cos^2(el)) - sin(el))
+	 * Invert for h given a max practical rho.  We take rho_max = 5000 km
+	 * (well beyond any LEO pass) and solve for h.
+	 *
+	 * Simpler conservative bound: h_max = rho_max / sin(el) for el > 0.
+	 */
+	rho = 5000.0;	/* max practical slant range; well beyond LEO/MEO */
+	h_max = sqrt(rho * rho + 2.0 * WGS72_AE * rho * sin_el
+				 + WGS72_AE * WGS72_AE) - WGS72_AE;
+
+	return h_max;
+}
+
+
+/* ----------------------------------------------------------------
+ * Helper: angular radius of ground visibility footprint
+ *
+ * For a satellite at altitude h km, the half-angle of the visibility
+ * cone (Earth-center angle) at min_el elevation is:
+ *   lambda = arccos(Re / (Re + h) * cos(el)) - el
+ * Returns degrees.
+ * ----------------------------------------------------------------
+ */
+static inline double
+ground_footprint_deg(double sma_km, double min_el_deg)
+{
+	double h_km = sma_km - WGS72_AE;
+	double el_rad, cos_ratio, lambda;
+
+	if (h_km <= 0.0)
+		return 0.0;
+
+	el_rad = min_el_deg * DEG_TO_RAD;
+	cos_ratio = WGS72_AE / (WGS72_AE + h_km) * cos(el_rad);
+
+	if (cos_ratio >= 1.0)
+		return 0.0;
+
+	lambda = acos(cos_ratio) - el_rad;
+
+	return lambda * RAD_TO_DEG;
+}
+
+
+/* ----------------------------------------------------------------
+ * Helper: J2 secular RAAN precession rate in radians/day
+ *
+ *   dOmega/dt = -1.5 * n * J2 * (Re/a)^2 * cos(i)
+ *
+ * where n is mean motion in rad/s, J2 and Re are WGS-72.
+ * Result in rad/day (multiply rad/s by 86400).
+ *
+ * Uses only the J2 zonal harmonic (no J3/J4 short-period terms).
+ * For LEO this can accumulate ~0.5 deg/day error from J3 short-
+ * period oscillations.  Acceptable because (a) the RAAN window
+ * includes footprint + Earth rotation padding, and (b) any query
+ * window >= ~12 hours bypasses the RAAN filter entirely.
+ * ----------------------------------------------------------------
+ */
+static inline double
+j2_raan_rate(double sma_km, double inc_rad)
+{
+	double a = sma_km;
+	double ratio, n_rad_s;
+
+	if (a <= 0.0)
+		return 0.0;
+
+	ratio = WGS72_AE / a;
+	n_rad_s = sqrt(WGS72_MU / (a * a * a));
+
+	return -1.5 * n_rad_s * WGS72_J2 * ratio * ratio * cos(inc_rad)
+		   * SECONDS_PER_DAY;
+}
+
+
+/* ----------------------------------------------------------------
+ * Traversal state carried down the tree during index scans.
+ * Accumulates SMA and inclination ranges from L0 and L1.
+ * ----------------------------------------------------------------
+ */
+typedef struct OrbitalTraversal
+{
+	double	sma_low;
+	double	sma_high;
+	double	inc_low;
+	double	inc_high;
+} OrbitalTraversal;
+
+
+/* ----------------------------------------------------------------
+ * Query parameter extraction from observer_window composite type
+ *
+ * observer_window is: (obs observer, t_start timestamptz,
+ *                      t_end timestamptz, min_el float8)
+ * ----------------------------------------------------------------
+ */
+typedef struct ObserverWindow
+{
+	pg_observer obs;
+	double		jd_start;
+	double		jd_end;
+	double		jd_mid;
+	double		min_el_deg;
+} ObserverWindow;
+
+
+static void
+extract_observer_window(HeapTupleHeader composite, ObserverWindow *win)
+{
+	bool	isnull;
+	Datum	val;
+	int64	ts_start, ts_end;
+
+	/*
+	 * Field 1: obs (observer -- fixed-size 24B, pass-by-reference).
+	 * GetAttributeByNum returns a Datum that is a direct pointer to
+	 * the pg_observer bytes in the composite tuple.  No varlena header.
+	 */
+	val = GetAttributeByNum(composite, 1, &isnull);
+	if (isnull)
+		ereport(ERROR,
+				(errcode(ERRCODE_NULL_VALUE_NOT_ALLOWED),
+				 errmsg("observer_window.obs must not be NULL")));
+
+	memcpy(&win->obs, DatumGetPointer(val), sizeof(pg_observer));
+
+	/* Field 2: t_start (timestamptz) */
+	val = GetAttributeByNum(composite, 2, &isnull);
+	if (isnull)
+		ereport(ERROR,
+				(errcode(ERRCODE_NULL_VALUE_NOT_ALLOWED),
+				 errmsg("observer_window.t_start must not be NULL")));
+	ts_start = DatumGetTimestampTz(val);
+	win->jd_start = timestamptz_to_jd(ts_start);
+
+	/* Field 3: t_end (timestamptz) */
+	val = GetAttributeByNum(composite, 3, &isnull);
+	if (isnull)
+		ereport(ERROR,
+				(errcode(ERRCODE_NULL_VALUE_NOT_ALLOWED),
+				 errmsg("observer_window.t_end must not be NULL")));
+	ts_end = DatumGetTimestampTz(val);
+	win->jd_end = timestamptz_to_jd(ts_end);
+
+	/* Validate time window ordering */
+	if (win->jd_end < win->jd_start)
+		ereport(ERROR,
+				(errcode(ERRCODE_INVALID_PARAMETER_VALUE),
+				 errmsg("observer_window.t_end must not precede t_start")));
+
+	win->jd_mid = (win->jd_start + win->jd_end) / 2.0;
+
+	/* Field 4: min_el (float8), clamped to [0, 90] */
+	val = GetAttributeByNum(composite, 4, &isnull);
+	if (isnull)
+		win->min_el_deg = 10.0;	/* default */
+	else
+	{
+		win->min_el_deg = DatumGetFloat8(val);
+		if (win->min_el_deg < 0.0)
+			win->min_el_deg = 0.0;
+		if (win->min_el_deg > 90.0)
+			win->min_el_deg = 90.0;
+	}
+}
+
+
+/* ----------------------------------------------------------------
+ * Comparison function for picksplit qsort.
+ * Sorts by extracted key value (SMA or inclination).
+ * ----------------------------------------------------------------
+ */
+typedef struct
+{
+	int		orig_index;
+	double	key;
+} PicksplitEntry;
+
+static int
+picksplit_entry_cmp(const void *a, const void *b)
+{
+	double ka = ((const PicksplitEntry *) a)->key;
+	double kb = ((const PicksplitEntry *) b)->key;
+
+	if (ka < kb)
+		return -1;
+	if (ka > kb)
+		return 1;
+	return 0;
+}
+
+
+/* ================================================================
+ * SP-GiST support functions (5 required)
+ * ================================================================
+ */
+
+
+/*
+ * spgist_tle_config -- declare trie type system
+ *
+ * No prefix data (VOIDOID): bin ranges encoded entirely in node labels.
+ * Labels are float8 (bin boundary values, sorted ascending).
+ * Leaves store full TLE unchanged.
+ */
+Datum
+spgist_tle_config(PG_FUNCTION_ARGS)
+{
+	spgConfigIn  *in  = (spgConfigIn *) PG_GETARG_POINTER(0);
+	spgConfigOut *out = (spgConfigOut *) PG_GETARG_POINTER(1);
+
+	out->prefixType    = VOIDOID;
+	out->labelType     = FLOAT8OID;
+	out->leafType      = in->attType;	/* tle type */
+	out->canReturnData = true;
+	out->longValuesOK  = false;
+
+	PG_RETURN_VOID();
+}
+
+
+/*
+ * spgist_tle_choose -- route a new TLE to the correct child node
+ *
+ * L0: route by SMA.  L1: route by inclination.
+ * restDatum = leafDatum (full TLE unchanged), matching the quad-tree
+ * precedent where the tree terminates by depth, not value exhaustion.
+ *
+ * At level >= 2, all nodes are allTheSame (no further partitioning).
+ */
+Datum
+spgist_tle_choose(PG_FUNCTION_ARGS)
+{
+	spgChooseIn  *in  = (spgChooseIn *)  PG_GETARG_POINTER(0);
+	spgChooseOut *out = (spgChooseOut *) PG_GETARG_POINTER(1);
+
+	pg_tle *tle = (pg_tle *) DatumGetPointer(in->leafDatum);
+	int		level = in->level;
+	double	val;
+
+	/* Extract the routing key for this level */
+	if (level == 0)
+		val = tle_sma_km(tle);
+	else
+		val = tle->inclination;
+
+	if (in->allTheSame || level >= SPGIST_TLE_MAX_LEVEL)
+	{
+		out->resultType = spgMatchNode;
+		out->result.matchNode.nodeN = 0;
+		out->result.matchNode.levelAdd = 1;
+		out->result.matchNode.restDatum = in->leafDatum;
+		PG_RETURN_VOID();
+	}
+
+	/*
+	 * Find the bin: labels are sorted ascending, each label marks the
+	 * lower bound of a bin.  We want the last bin whose label <= val.
+	 */
+	{
+		int		best = 0;
+		int		i;
+
+		for (i = 1; i < in->nNodes; i++)
+		{
+			if (DatumGetFloat8(in->nodeLabels[i]) <= val)
+				best = i;
+			else
+				break;
+		}
+
+		out->resultType = spgMatchNode;
+		out->result.matchNode.nodeN = best;
+		out->result.matchNode.levelAdd = 1;
+		out->result.matchNode.restDatum = in->leafDatum;
+	}
+
+	PG_RETURN_VOID();
+}
+
+
+/*
+ * spgist_tle_picksplit -- split a leaf page using equal-population strategy
+ *
+ * Sort by the current level's element (SMA at L0, inclination at L1),
+ * divide into floor(sqrt(nTuples)) bins clamped to [2, 16].
+ * At level >= 2, emit a single allTheSame node.
+ */
+Datum
+spgist_tle_picksplit(PG_FUNCTION_ARGS)
+{
+	spgPickSplitIn  *in  = (spgPickSplitIn *)  PG_GETARG_POINTER(0);
+	spgPickSplitOut *out = (spgPickSplitOut *) PG_GETARG_POINTER(1);
+
+	int		level    = in->level;
+	int		nTuples  = in->nTuples;
+	int		nBins, perBin, remainder;
+	int		i, bin, pos;
+	PicksplitEntry *entries;
+
+	/*
+	 * At level >= 2 we have no further partitioning dimension.
+	 * Emit a single allTheSame node that accepts everything.
+	 */
+	if (level >= SPGIST_TLE_MAX_LEVEL)
+	{
+		out->nNodes = 1;
+		out->hasPrefix = false;
+		out->prefixDatum = (Datum) 0;
+		out->nodeLabels = (Datum *) palloc(sizeof(Datum));
+		out->nodeLabels[0] = Float8GetDatum(0.0);
+		out->mapTuplesToNodes = (int *) palloc(sizeof(int) * nTuples);
+		out->leafTupleDatums = (Datum *) palloc(sizeof(Datum) * nTuples);
+
+		for (i = 0; i < nTuples; i++)
+		{
+			out->mapTuplesToNodes[i] = 0;
+			out->leafTupleDatums[i] = in->datums[i];
+		}
+
+		PG_RETURN_VOID();
+	}
+
+	/* Extract and sort by the current level's element */
+	entries = (PicksplitEntry *) palloc(sizeof(PicksplitEntry) * nTuples);
+
+	for (i = 0; i < nTuples; i++)
+	{
+		pg_tle *tle = (pg_tle *) DatumGetPointer(in->datums[i]);
+
+		entries[i].orig_index = i;
+		if (level == 0)
+			entries[i].key = tle_sma_km(tle);
+		else
+			entries[i].key = tle->inclination;
+	}
+
+	qsort(entries, nTuples, sizeof(PicksplitEntry), picksplit_entry_cmp);
+
+	/* Equal-population split: floor(sqrt(n)) bins, clamped [2, 16] */
+	nBins = (int) floor(sqrt((double) nTuples));
+	if (nBins < 2)
+		nBins = 2;
+	if (nBins > 16)
+		nBins = 16;
+	/* Prevent over-read: never more bins than tuples */
+	if (nBins > nTuples)
+		nBins = nTuples;
+
+	perBin = nTuples / nBins;
+	remainder = nTuples % nBins;
+
+	out->nNodes = nBins;
+	out->hasPrefix = false;
+	out->prefixDatum = (Datum) 0;
+	out->nodeLabels = (Datum *) palloc(sizeof(Datum) * nBins);
+	out->mapTuplesToNodes = (int *) palloc(sizeof(int) * nTuples);
+	out->leafTupleDatums = (Datum *) palloc(sizeof(Datum) * nTuples);
+
+	pos = 0;
+	for (bin = 0; bin < nBins; bin++)
+	{
+		int count = perBin + (bin < remainder ? 1 : 0);
+
+		/* Node label = key value of the first entry in this bin */
+		out->nodeLabels[bin] = Float8GetDatum(entries[pos].key);
+
+		for (i = 0; i < count; i++)
+		{
+			int orig = entries[pos + i].orig_index;
+
+			out->mapTuplesToNodes[orig] = bin;
+			out->leafTupleDatums[orig] = in->datums[orig];
+		}
+
+		pos += count;
+	}
+
+	pfree(entries);
+
+	PG_RETURN_VOID();
+}
+
+
+/*
+ * spgist_tle_inner_consistent -- prune child nodes during index scan
+ *
+ * L0: skip bins whose perigee altitude exceeds max visible altitude.
+ *     The bin's SMA range is [label[i], label[i+1]).
+ * L1: skip bins whose inclination is too low to reach observer latitude.
+ *     A satellite with inclination i has ground track bounded by [-i, +i].
+ *     Observer at latitude phi needs i + footprint >= |phi|.
+ *
+ * Propagates OrbitalTraversal state to surviving children via
+ * traversalMemoryContext for the RAAN filter at leaf level.
+ */
+Datum
+spgist_tle_inner_consistent(PG_FUNCTION_ARGS)
+{
+	spgInnerConsistentIn  *in  = (spgInnerConsistentIn *)  PG_GETARG_POINTER(0);
+	spgInnerConsistentOut *out = (spgInnerConsistentOut *) PG_GETARG_POINTER(1);
+
+	int		level = in->level;
+	int		nkeys = in->nkeys;
+	int		i;
+	ObserverWindow win;
+	bool	have_query = false;
+
+	/* Extract query from scankeys -- we need the &? operator's arg */
+	for (i = 0; i < nkeys; i++)
+	{
+		if (in->scankeys[i].sk_strategy == 1)
+		{
+			HeapTupleHeader composite;
+
+			composite = DatumGetHeapTupleHeader(in->scankeys[i].sk_argument);
+			extract_observer_window(composite, &win);
+			have_query = true;
+			break;
+		}
+	}
+
+	/* Allocate output arrays */
+	out->nodeNumbers = (int *) palloc(sizeof(int) * in->nNodes);
+	out->levelAdds = (int *) palloc(sizeof(int) * in->nNodes);
+	out->reconstructedValues = NULL;
+	out->traversalValues = (void **) palloc(sizeof(void *) * in->nNodes);
+	out->nNodes = 0;
+
+	for (i = 0; i < in->nNodes; i++)
+	{
+		OrbitalTraversal *parent_trav;
+		OrbitalTraversal *child_trav;
+		double	bin_low, bin_high;
+		bool	dominated = false;
+
+		/* Decode bin range from labels */
+		bin_low = DatumGetFloat8(in->nodeLabels[i]);
+		if (i < in->nNodes - 1)
+			bin_high = DatumGetFloat8(in->nodeLabels[i + 1]);
+		else
+			bin_high = INFINITY;
+
+		/* Inherit parent traversal state or initialize */
+		if (in->traversalValue)
+			parent_trav = (OrbitalTraversal *) in->traversalValue;
+		else
+			parent_trav = NULL;
+
+		/* Pruning logic per level */
+		if (have_query && level == 0)
+		{
+			/*
+			 * L0: SMA range narrowing only — no altitude pruning.
+			 *
+			 * We cannot prune SMA bins by altitude because eccentricity
+			 * is not available at the inner node level.  A satellite
+			 * at SMA 70,000 km with e=0.88 has perigee ~2,000 km —
+			 * well within typical max_alt.  Without knowing e, any SMA
+			 * bin could contain satellites with perigee near Earth's
+			 * surface.
+			 *
+			 * L0 still helps by narrowing the SMA range passed to L1
+			 * for computing a tighter ground footprint.
+			 */
+		}
+		else if (have_query && level == 1)
+		{
+			/*
+			 * L1: Inclination pruning.
+			 * bin_high is the upper bound on inclination in this bin.
+			 * A satellite with inclination i has ground track [-i, +i].
+			 * The observer at latitude phi can see it if:
+			 *   i + footprint >= |phi|
+			 *
+			 * Use the parent SMA range to compute a conservative footprint.
+			 * The largest footprint comes from the HIGHEST altitude (footprint
+			 * grows with altitude: GEO sees 71+ degrees, LEO sees ~7 degrees).
+			 * Use sma_high for conservatism — never prune objects that the
+			 * leaf filter would accept.
+			 */
+			double obs_lat = fabs(win.obs.lat);
+			double sma_for_footprint;
+			double footprint;
+
+			if (parent_trav)
+				sma_for_footprint = parent_trav->sma_high;
+			else
+				sma_for_footprint = 50000.0;	/* above GEO — maximum footprint */
+
+			footprint = ground_footprint_deg(sma_for_footprint,
+											 win.min_el_deg) * DEG_TO_RAD;
+
+			if (bin_high + footprint < obs_lat)
+				dominated = true;
+		}
+
+		if (!dominated)
+		{
+			int idx = out->nNodes;
+
+			/* Build child traversal state */
+			child_trav = (OrbitalTraversal *)
+				MemoryContextAlloc(in->traversalMemoryContext,
+								   sizeof(OrbitalTraversal));
+
+			if (parent_trav)
+				memcpy(child_trav, parent_trav, sizeof(OrbitalTraversal));
+			else
+			{
+				child_trav->sma_low = 0.0;
+				child_trav->sma_high = INFINITY;
+				child_trav->inc_low = 0.0;
+				child_trav->inc_high = M_PI;
+			}
+
+			/* Narrow bounds based on current level */
+			if (level == 0)
+			{
+				child_trav->sma_low = bin_low;
+				child_trav->sma_high = bin_high;
+			}
+			else if (level == 1)
+			{
+				child_trav->inc_low = bin_low;
+				child_trav->inc_high = bin_high;
+			}
+
+			out->nodeNumbers[idx] = i;
+			out->levelAdds[idx] = 1;
+			out->traversalValues[idx] = child_trav;
+			out->nNodes++;
+		}
+	}
+
+	PG_RETURN_VOID();
+}
+
+
+/* ----------------------------------------------------------------
+ * Shared filter: three-stage visibility check on a single TLE.
+ *
+ * 1. Perigee altitude check (with eccentricity)
+ * 2. Inclination + ground footprint vs observer latitude
+ * 3. RAAN query-time filter (J2 precession to query midpoint)
+ *
+ * Called from both leaf_consistent (index scan) and
+ * tle_visibility_possible (sequential scan / standalone operator).
+ * ----------------------------------------------------------------
+ */
+static bool
+tle_passes_visibility_filter(const pg_tle *tle, const ObserverWindow *win)
+{
+	double	sma, perigee_alt, max_alt;
+	double	obs_lat_abs, footprint_rad;
+	double	dt_days, raan_projected, lst;
+	double	earth_rot_rad, raan_window_half, raan_diff;
+
+	/* Reject degenerate TLEs (decay, error data) */
+	if (tle->mean_motion <= 0.0)
+		return false;
+
+	sma = tle_sma_km(tle);
+
+	/* Filter 1: perigee altitude */
+	perigee_alt = sma * (1.0 - tle->eccentricity) - WGS72_AE;
+	max_alt = max_visible_altitude_km(win->min_el_deg);
+	if (perigee_alt > max_alt)
+		return false;
+
+	/* Filter 2: inclination + footprint vs observer latitude */
+	obs_lat_abs = fabs(win->obs.lat);
+	footprint_rad = ground_footprint_deg(sma, win->min_el_deg) * DEG_TO_RAD;
+
+	if (tle->inclination + footprint_rad < obs_lat_abs)
+		return false;
+
+	/* Filter 3: RAAN alignment via J2 secular precession */
+	dt_days = win->jd_mid - tle->epoch;
+	raan_projected = tle->raan
+		+ j2_raan_rate(sma, tle->inclination) * dt_days;
+	raan_projected = fmod(raan_projected, 2.0 * M_PI);
+	if (raan_projected < 0.0)
+		raan_projected += 2.0 * M_PI;
+
+	/* Observer LST at query midpoint */
+	lst = gmst_from_jd(win->jd_mid) + win->obs.lon;
+	lst = fmod(lst, 2.0 * M_PI);
+	if (lst < 0.0)
+		lst += 2.0 * M_PI;
+
+	/* RAAN window: Earth rotation during query + footprint pad */
+	earth_rot_rad = (win->jd_end - win->jd_start) * EARTH_ROT_RAD_PER_DAY;
+	raan_window_half = earth_rot_rad / 2.0
+		+ ground_footprint_deg(sma, win->min_el_deg) * DEG_TO_RAD;
+
+	if (raan_window_half >= M_PI)
+		return true;	/* full rotation -- pass everything */
+
+	raan_diff = fabs(raan_projected - lst);
+	if (raan_diff > M_PI)
+		raan_diff = 2.0 * M_PI - raan_diff;
+
+	return (raan_diff <= raan_window_half);
+}
+
+
+/*
+ * spgist_tle_leaf_consistent -- final check on a leaf TLE
+ *
+ * Delegates to tle_passes_visibility_filter() for the &? operator.
+ * recheck = false: the &? operator IS the superset filter.
+ * The user runs predict_passes() on survivors for SGP4 ground truth.
+ */
+Datum
+spgist_tle_leaf_consistent(PG_FUNCTION_ARGS)
+{
+	spgLeafConsistentIn  *in  = (spgLeafConsistentIn *)  PG_GETARG_POINTER(0);
+	spgLeafConsistentOut *out = (spgLeafConsistentOut *) PG_GETARG_POINTER(1);
+
+	pg_tle *tle;
+	int		i;
+	bool	result = true;
+
+	tle = (pg_tle *) DatumGetPointer(in->leafDatum);
+	out->leafValue = in->leafDatum;
+	out->recheck = false;
+
+	for (i = 0; i < in->nkeys; i++)
+	{
+		if (in->scankeys[i].sk_strategy == 1)
+		{
+			ObserverWindow win;
+			HeapTupleHeader composite;
+
+			composite = DatumGetHeapTupleHeader(in->scankeys[i].sk_argument);
+			extract_observer_window(composite, &win);
+
+			if (!tle_passes_visibility_filter(tle, &win))
+			{
+				result = false;
+				break;
+			}
+		}
+	}
+
+	PG_RETURN_BOOL(result);
+}
+
+
+/* ================================================================
+ * Operator function: &? (visibility cone check)
+ * ================================================================
+ */
+
+
+/*
+ * tle_visibility_possible(tle, observer_window) -> bool
+ *
+ * Standalone operator: can the satellite possibly be visible from
+ * this observer during this time window?  Combines altitude check,
+ * latitude/inclination check, and RAAN filter.
+ *
+ * This is the same logic as leaf_consistent, callable directly
+ * as a SQL operator for sequential scans or WHERE clauses.
+ *
+ * The indexed column (tle) MUST be the left argument so that
+ * PostgreSQL can form a ScanKey and pass it to inner_consistent
+ * for tree-level pruning.  See skey.h:23-26.
+ */
+Datum
+tle_visibility_possible(PG_FUNCTION_ARGS)
+{
+	pg_tle		   *tle = (pg_tle *) PG_GETARG_POINTER(0);
+	HeapTupleHeader composite = PG_GETARG_HEAPTUPLEHEADER(1);
+	ObserverWindow	win;
+
+	extract_observer_window(composite, &win);
+
+	PG_RETURN_BOOL(tle_passes_visibility_filter(tle, &win));
+}

test/expected/spgist_tle.out

@@ -0,0 +1,354 @@
+-- Test SP-GiST orbital trie index and &? visibility cone operator
+SET client_min_messages = warning;
+CREATE EXTENSION IF NOT EXISTS pg_orrery;
+RESET client_min_messages;
+-- ============================================================
+-- Test table with mixed orbital regimes
+-- ============================================================
+CREATE TABLE test_spgist (
+    id serial,
+    name text,
+    tle tle
+);
+-- ISS (LEO, ~400km, 51.64 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('ISS',
+    '1 25544U 98067A   24001.50000000  .00016717  00000-0  10270-3 0  9025
+2 25544  51.6400 208.9163 0006703  30.1694  61.7520 15.50100486 00001');
+-- Hubble (LEO, ~540km, 28.47 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('Hubble',
+    '1 20580U 90037B   24001.50000000  .00000790  00000+0  39573-4 0  9992
+2 20580  28.4705  61.4398 0002797 317.3115  42.7577 15.09395228 00008');
+-- GPS IIR-M (MEO, ~20200km, 55.44 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('GPS-IIR',
+    '1 28874U 05038A   24001.50000000  .00000012  00000+0  00000+0 0  9993
+2 28874  55.4408 300.3467 0117034  51.6543 309.5420  2.00557079 00006');
+-- Equatorial-LEO (same altitude as ISS, 5 deg inclination)
+INSERT INTO test_spgist (name, tle) VALUES ('Equatorial-LEO',
+    '1 99901U 24999A   24001.50000000  .00016717  00000-0  10270-3 0  9990
+2 99901   5.0000 208.9163 0006703  30.1694  61.7520 15.50100486 00001');
+-- SSO-800 (Sun-synchronous, ~800km, 98.7 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('SSO-800',
+    '1 99902U 24999B   24001.50000000  .00000100  00000+0  50000-4 0  9991
+2 99902  98.7000 120.0000 0001000  90.0000 270.0000 14.19553000 00001');
+-- GEO-SAT (Geostationary, ~35786km, 0.04 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('GEO-SAT',
+    '1 99903U 24999C   24001.50000000  .00000000  00000+0  00000+0 0  9992
+2 99903   0.0400 270.0000 0003000   0.0000 180.0000  1.00273791 00001');
+-- ============================================================
+-- Test 1: Operator standalone — ISS from Eagle Idaho (2h window)
+-- Eagle Idaho: 43.6977N 116.3535W, 760m elevation
+-- ISS passes altitude and inclination checks, but RAAN filter
+-- rejects it — the orbital plane isn't overhead during this
+-- specific 2-hour window (correct physics, see Test 5 for 24h).
+-- ============================================================
+SELECT name,
+    tle &? ROW(
+        observer('43.6977N 116.3535W 760m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-01 02:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'ISS';
+ name | visible 
+------+---------
+ ISS  | f
+(1 row)
+
+-- ============================================================
+-- Test 2: Equatorial-LEO NOT visible from Eagle Idaho
+-- 5 deg inc + ~12 deg footprint = 17 deg < 43.7 deg latitude
+-- ============================================================
+SELECT name,
+    tle &? ROW(
+        observer('43.6977N 116.3535W 760m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-01 02:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'Equatorial-LEO';
+      name      | visible 
+----------------+---------
+ Equatorial-LEO | f
+(1 row)
+
+-- ============================================================
+-- Test 3: Create SP-GiST index, verify index scan with positive
+-- results.  Equatorial observer at 0E — SSO-800 RAAN (120 deg)
+-- aligns with LST near 0E at this epoch, so it passes.
+-- ============================================================
+CREATE INDEX test_spgist_idx ON test_spgist USING spgist (tle tle_spgist_ops);
+SET enable_seqscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-01 02:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ SSO-800
+(1 row)
+
+RESET enable_seqscan;
+-- ============================================================
+-- Test 4: Seqscan vs index scan consistency — same query must
+-- return identical results regardless of scan method.
+-- ============================================================
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-01 02:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ SSO-800
+(1 row)
+
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+-- ============================================================
+-- Test 5: 24-hour window — RAAN filter bypassed (full Earth
+-- rotation).  Only ISS and SSO-800 pass inclination from Eagle
+-- Idaho (43.7 deg).  Hubble (28.5+14.8=43.3 deg) barely fails.
+-- GPS-IIR and GEO-SAT filtered by altitude.
+-- ============================================================
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ ISS
+ SSO-800
+(2 rows)
+
+-- ============================================================
+-- Test 6: High min_el (45 deg) changes footprint — wider
+-- footprint lets more inclinations through.  Same 24h window.
+-- ============================================================
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    45.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ ISS
+ SSO-800
+(2 rows)
+
+-- ============================================================
+-- Test 7: GiST coexistence — both index types on same table
+-- ============================================================
+CREATE INDEX test_gist_idx ON test_spgist USING gist (tle);
+-- GiST overlap query still works
+SELECT a.name AS sat_a, b.name AS sat_b, a.tle && b.tle AS overlaps
+FROM test_spgist a, test_spgist b
+WHERE a.name = 'ISS' AND b.name = 'Hubble';
+ sat_a | sat_b  | overlaps 
+-------+--------+----------
+ ISS   | Hubble | f
+(1 row)
+
+-- SP-GiST query still works alongside GiST
+SET enable_seqscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ ISS
+ SSO-800
+(2 rows)
+
+RESET enable_seqscan;
+-- ============================================================
+-- Test 8: NULL TLE handling — NULLs should be excluded
+-- ============================================================
+INSERT INTO test_spgist (name, tle) VALUES ('NULL-SAT', NULL);
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ ISS
+ SSO-800
+(2 rows)
+
+-- ============================================================
+-- Test 9: Degenerate TLE (mean_motion = 0) — rejected by filter
+-- ============================================================
+INSERT INTO test_spgist (name, tle) VALUES ('DECAYED',
+    '1 99904U 24999D   24001.50000000  .00000000  00000+0  00000+0 0  9993
+2 99904   0.0000   0.0000 0000000   0.0000   0.0000  0.00000000 00001');
+SELECT name,
+    tle &? ROW(
+        observer('0.0N 0.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-02 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'DECAYED';
+  name   | visible 
+---------+---------
+ DECAYED | f
+(1 row)
+
+-- ============================================================
+-- Test 10: Polar observer (90N) — only ISS and SSO-800 reach
+-- the pole.  ISS (51.6 + footprint) < 90, so only SSO-800
+-- (retrograde, 98.7 deg inc > 90 deg) passes.  24h window.
+-- ============================================================
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('90.0N 0.0E 0m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ SSO-800
+(1 row)
+
+-- ============================================================
+-- Test 11: Zero-duration window — sees only what is directly
+-- overhead at the instant.  RAAN window = footprint only.
+-- ============================================================
+SELECT name,
+    tle &? ROW(
+        observer('0.0N 0.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-01 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'ISS';
+ name | visible 
+------+---------
+ ISS  | f
+(1 row)
+
+-- ============================================================
+-- Test 12: Index-vs-seqscan consistency on 24h Eagle Idaho
+-- (the primary correctness test, now after all inserts)
+-- ============================================================
+SET enable_seqscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ ISS
+ SSO-800
+(2 rows)
+
+RESET enable_seqscan;
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+  name   
+---------
+ ISS
+ SSO-800
+(2 rows)
+
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+-- ============================================================
+-- Test 13: HEO at high latitude — GTO-class orbit (low inc,
+-- high SMA, high eccentricity) from Tromsø (69.6°N).
+-- The large SMA gives a huge footprint that compensates for the
+-- low inclination.  Must pass the seqscan operator check.
+-- Regression test for the L1 pruning bug (sma_low vs sma_high).
+-- ============================================================
+-- GTO debris: inc 5 deg, perigee ~250 km, apogee ~35786 km
+INSERT INTO test_spgist (name, tle) VALUES ('GTO-DEBRIS',
+    '1 99905U 24999E   24001.50000000  .00000100  00000+0  10000-3 0  9994
+2 99905   5.0000 210.0000 7300000  30.0000  61.0000  2.25600000 00001');
+-- Seqscan: GTO-DEBRIS from Tromsø — must be visible
+-- inc 5 deg + footprint(SMA ~25000) ~65 deg = 70 > 69.6
+SELECT name,
+    tle &? ROW(
+        observer('69.6N 19.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-02 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'GTO-DEBRIS';
+    name    | visible 
+------------+---------
+ GTO-DEBRIS | t
+(1 row)
+
+-- Index scan: same query, must return the same result
+SET enable_seqscan = off;
+SELECT name,
+    tle &? ROW(
+        observer('69.6N 19.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-02 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'GTO-DEBRIS';
+    name    | visible 
+------------+---------
+ GTO-DEBRIS | t
+(1 row)
+
+RESET enable_seqscan;
+-- ============================================================
+-- Cleanup
+-- ============================================================
+DROP TABLE test_spgist;

test/sql/spgist_tle.sql

@@ -0,0 +1,316 @@
+-- Test SP-GiST orbital trie index and &? visibility cone operator
+SET client_min_messages = warning;
+CREATE EXTENSION IF NOT EXISTS pg_orrery;
+RESET client_min_messages;
+
+-- ============================================================
+-- Test table with mixed orbital regimes
+-- ============================================================
+CREATE TABLE test_spgist (
+    id serial,
+    name text,
+    tle tle
+);
+
+-- ISS (LEO, ~400km, 51.64 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('ISS',
+    '1 25544U 98067A   24001.50000000  .00016717  00000-0  10270-3 0  9025
+2 25544  51.6400 208.9163 0006703  30.1694  61.7520 15.50100486 00001');
+
+-- Hubble (LEO, ~540km, 28.47 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('Hubble',
+    '1 20580U 90037B   24001.50000000  .00000790  00000+0  39573-4 0  9992
+2 20580  28.4705  61.4398 0002797 317.3115  42.7577 15.09395228 00008');
+
+-- GPS IIR-M (MEO, ~20200km, 55.44 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('GPS-IIR',
+    '1 28874U 05038A   24001.50000000  .00000012  00000+0  00000+0 0  9993
+2 28874  55.4408 300.3467 0117034  51.6543 309.5420  2.00557079 00006');
+
+-- Equatorial-LEO (same altitude as ISS, 5 deg inclination)
+INSERT INTO test_spgist (name, tle) VALUES ('Equatorial-LEO',
+    '1 99901U 24999A   24001.50000000  .00016717  00000-0  10270-3 0  9990
+2 99901   5.0000 208.9163 0006703  30.1694  61.7520 15.50100486 00001');
+
+-- SSO-800 (Sun-synchronous, ~800km, 98.7 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('SSO-800',
+    '1 99902U 24999B   24001.50000000  .00000100  00000+0  50000-4 0  9991
+2 99902  98.7000 120.0000 0001000  90.0000 270.0000 14.19553000 00001');
+
+-- GEO-SAT (Geostationary, ~35786km, 0.04 deg)
+INSERT INTO test_spgist (name, tle) VALUES ('GEO-SAT',
+    '1 99903U 24999C   24001.50000000  .00000000  00000+0  00000+0 0  9992
+2 99903   0.0400 270.0000 0003000   0.0000 180.0000  1.00273791 00001');
+
+
+-- ============================================================
+-- Test 1: Operator standalone — ISS from Eagle Idaho (2h window)
+-- Eagle Idaho: 43.6977N 116.3535W, 760m elevation
+-- ISS passes altitude and inclination checks, but RAAN filter
+-- rejects it — the orbital plane isn't overhead during this
+-- specific 2-hour window (correct physics, see Test 5 for 24h).
+-- ============================================================
+SELECT name,
+    tle &? ROW(
+        observer('43.6977N 116.3535W 760m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-01 02:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'ISS';
+
+
+-- ============================================================
+-- Test 2: Equatorial-LEO NOT visible from Eagle Idaho
+-- 5 deg inc + ~12 deg footprint = 17 deg < 43.7 deg latitude
+-- ============================================================
+SELECT name,
+    tle &? ROW(
+        observer('43.6977N 116.3535W 760m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-01 02:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'Equatorial-LEO';
+
+
+-- ============================================================
+-- Test 3: Create SP-GiST index, verify index scan with positive
+-- results.  Equatorial observer at 0E — SSO-800 RAAN (120 deg)
+-- aligns with LST near 0E at this epoch, so it passes.
+-- ============================================================
+CREATE INDEX test_spgist_idx ON test_spgist USING spgist (tle tle_spgist_ops);
+
+SET enable_seqscan = off;
+
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-01 02:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+
+RESET enable_seqscan;
+
+
+-- ============================================================
+-- Test 4: Seqscan vs index scan consistency — same query must
+-- return identical results regardless of scan method.
+-- ============================================================
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('0.0N 0.0E 0m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-01 02:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+
+-- ============================================================
+-- Test 5: 24-hour window — RAAN filter bypassed (full Earth
+-- rotation).  Only ISS and SSO-800 pass inclination from Eagle
+-- Idaho (43.7 deg).  Hubble (28.5+14.8=43.3 deg) barely fails.
+-- GPS-IIR and GEO-SAT filtered by altitude.
+-- ============================================================
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+
+
+-- ============================================================
+-- Test 6: High min_el (45 deg) changes footprint — wider
+-- footprint lets more inclinations through.  Same 24h window.
+-- ============================================================
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    45.0
+)::observer_window
+ORDER BY name;
+
+
+-- ============================================================
+-- Test 7: GiST coexistence — both index types on same table
+-- ============================================================
+CREATE INDEX test_gist_idx ON test_spgist USING gist (tle);
+
+-- GiST overlap query still works
+SELECT a.name AS sat_a, b.name AS sat_b, a.tle && b.tle AS overlaps
+FROM test_spgist a, test_spgist b
+WHERE a.name = 'ISS' AND b.name = 'Hubble';
+
+-- SP-GiST query still works alongside GiST
+SET enable_seqscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+RESET enable_seqscan;
+
+
+-- ============================================================
+-- Test 8: NULL TLE handling — NULLs should be excluded
+-- ============================================================
+INSERT INTO test_spgist (name, tle) VALUES ('NULL-SAT', NULL);
+
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+
+
+-- ============================================================
+-- Test 9: Degenerate TLE (mean_motion = 0) — rejected by filter
+-- ============================================================
+INSERT INTO test_spgist (name, tle) VALUES ('DECAYED',
+    '1 99904U 24999D   24001.50000000  .00000000  00000+0  00000+0 0  9993
+2 99904   0.0000   0.0000 0000000   0.0000   0.0000  0.00000000 00001');
+
+SELECT name,
+    tle &? ROW(
+        observer('0.0N 0.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-02 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'DECAYED';
+
+
+-- ============================================================
+-- Test 10: Polar observer (90N) — only ISS and SSO-800 reach
+-- the pole.  ISS (51.6 + footprint) < 90, so only SSO-800
+-- (retrograde, 98.7 deg inc > 90 deg) passes.  24h window.
+-- ============================================================
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('90.0N 0.0E 0m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+
+
+-- ============================================================
+-- Test 11: Zero-duration window — sees only what is directly
+-- overhead at the instant.  RAAN window = footprint only.
+-- ============================================================
+SELECT name,
+    tle &? ROW(
+        observer('0.0N 0.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-01 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'ISS';
+
+
+-- ============================================================
+-- Test 12: Index-vs-seqscan consistency on 24h Eagle Idaho
+-- (the primary correctness test, now after all inserts)
+-- ============================================================
+SET enable_seqscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+RESET enable_seqscan;
+
+SET enable_indexscan = off;
+SET enable_bitmapscan = off;
+SELECT name
+FROM test_spgist
+WHERE tle &? ROW(
+    observer('43.6977N 116.3535W 760m'),
+    '2024-01-01 00:00:00+00'::timestamptz,
+    '2024-01-02 00:00:00+00'::timestamptz,
+    10.0
+)::observer_window
+ORDER BY name;
+RESET enable_indexscan;
+RESET enable_bitmapscan;
+
+
+-- ============================================================
+-- Test 13: HEO at high latitude — GTO-class orbit (low inc,
+-- high SMA, high eccentricity) from Tromsø (69.6°N).
+-- The large SMA gives a huge footprint that compensates for the
+-- low inclination.  Must pass the seqscan operator check.
+-- Regression test for the L1 pruning bug (sma_low vs sma_high).
+-- ============================================================
+
+-- GTO debris: inc 5 deg, perigee ~250 km, apogee ~35786 km
+INSERT INTO test_spgist (name, tle) VALUES ('GTO-DEBRIS',
+    '1 99905U 24999E   24001.50000000  .00000100  00000+0  10000-3 0  9994
+2 99905   5.0000 210.0000 7300000  30.0000  61.0000  2.25600000 00001');
+
+-- Seqscan: GTO-DEBRIS from Tromsø — must be visible
+-- inc 5 deg + footprint(SMA ~25000) ~65 deg = 70 > 69.6
+SELECT name,
+    tle &? ROW(
+        observer('69.6N 19.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-02 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'GTO-DEBRIS';
+
+-- Index scan: same query, must return the same result
+SET enable_seqscan = off;
+SELECT name,
+    tle &? ROW(
+        observer('69.6N 19.0E 0m'),
+        '2024-01-01 00:00:00+00'::timestamptz,
+        '2024-01-02 00:00:00+00'::timestamptz,
+        10.0
+    )::observer_window AS visible
+FROM test_spgist
+WHERE name = 'GTO-DEBRIS';
+RESET enable_seqscan;
+
+
+-- ============================================================
+-- Cleanup
+-- ============================================================
+DROP TABLE test_spgist;