diff --git a/sql/pg_orbit--0.1.0.sql b/sql/pg_orbit--0.1.0.sql index 72e1365..917e663 100644 --- a/sql/pg_orbit--0.1.0.sql +++ b/sql/pg_orbit--0.1.0.sql @@ -444,11 +444,11 @@ COMMENT ON FUNCTION pass_visible(tle, observer, timestamptz, timestamptz) IS -- GiST operator support functions -- ============================================================ --- Overlap operator: do altitude bands intersect? +-- 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 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; @@ -461,7 +461,7 @@ CREATE OPERATOR && ( JOIN = areajoinsel ); -COMMENT ON OPERATOR && (tle, tle) IS 'Altitude band overlap — necessary condition for conjunction'; +COMMENT ON OPERATOR && (tle, tle) IS 'Orbital key overlap (altitude band AND inclination range) — necessary condition for conjunction'; CREATE OPERATOR <-> ( LEFTARG = tle, @@ -470,11 +470,11 @@ CREATE OPERATOR <-> ( COMMUTATOR = <-> ); -COMMENT ON OPERATOR <-> (tle, tle) IS 'Minimum altitude-band separation in km (0 if overlapping)'; +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 altitude-band indexing +-- GiST operator class for 2-D orbital indexing (altitude + inclination) -- ============================================================ -- GiST internal support functions diff --git a/src/gist_tle.c b/src/gist_tle.c index d9ef94a..5e3e12e 100644 --- a/src/gist_tle.c +++ b/src/gist_tle.c @@ -1,14 +1,18 @@ /* - * gist_tle.c -- GiST operator class for altitude-band indexing on TLE + * gist_tle.c -- GiST operator class for 2-D orbital indexing on TLE * - * Every TLE defines an orbit whose perigee and apogee altitudes form a - * 1-D range (the "altitude band"). Two orbits can only be in proximity - * if their altitude bands overlap -- a necessary but not sufficient - * condition for conjunction. + * Every TLE defines an orbit whose perigee/apogee altitudes and + * inclination form a 2-D bounding box in (altitude, inclination) space. + * Two orbits can only be in proximity if their altitude bands overlap + * AND their inclination ranges overlap -- both necessary but not + * sufficient conditions for conjunction. * - * The GiST index stores [perigee_km, apogee_km] ranges as internal - * keys, enabling fast coarse filtering. The && (overlap) operator is - * always rechecked: real conjunction screening requires propagation. + * The GiST index stores [alt_low, alt_high] x [inc_low, inc_high] keys, + * enabling fast coarse filtering. For leaf entries inc_low == inc_high + * (exact value); internal nodes widen to the bounding box of children. + * + * The && (overlap) operator is always rechecked: real conjunction + * screening requires propagation. * * Semi-major axis from Kepler's third law using WGS-72 constants: * a = (KE / n)^(2/3) [earth radii] @@ -37,30 +41,37 @@ PG_FUNCTION_INFO_V1(gist_tle_picksplit); PG_FUNCTION_INFO_V1(gist_tle_same); PG_FUNCTION_INFO_V1(gist_tle_distance); -/* Floating-point comparison tolerance (km) */ -#define ALT_EPSILON 1.0e-9 +/* Floating-point comparison tolerance (km and radians) */ +#define KEY_EPSILON 1.0e-9 /* - * Altitude band extracted from a TLE's mean elements. + * 2-D orbital key extracted from a TLE's mean elements. + * Altitude band (perigee/apogee) plus inclination range. * This is the GiST internal key -- much cheaper to compare * than propagating two full state vectors. + * + * For leaf entries: inc_low == inc_high (exact inclination). + * For internal nodes: bounding box over all children. */ -typedef struct tle_alt_range +typedef struct tle_orbital_key { - double low; /* perigee altitude, km */ - double high; /* apogee altitude, km */ -} tle_alt_range; + double alt_low; /* perigee altitude, km */ + double alt_high; /* apogee altitude, km */ + double inc_low; /* inclination lower bound, radians */ + double inc_high; /* inclination upper bound, radians */ +} tle_orbital_key; /* ---------------------------------------------------------------- - * tle_to_alt_range -- compute [perigee, apogee] from mean elements + * tle_to_orbital_key -- compute [perigee, apogee] x [inc, inc] * * Uses WGS-72 KE and AE (the only constants valid for SGP4 elements). - * Degenerate TLEs with n <= 0 map to a zero-width range at 0 km. + * Degenerate TLEs with n <= 0 map to zero-width ranges at 0. + * Inclination is stored in radians (same as pg_tle). * ---------------------------------------------------------------- */ static void -tle_to_alt_range(const pg_tle *tle, tle_alt_range *range) +tle_to_orbital_key(const pg_tle *tle, tle_orbital_key *key) { double n = tle->mean_motion; /* rad/min */ double e = tle->eccentricity; @@ -68,74 +79,91 @@ tle_to_alt_range(const pg_tle *tle, tle_alt_range *range) if (n <= 0.0) { - range->low = 0.0; - range->high = 0.0; + key->alt_low = 0.0; + key->alt_high = 0.0; + key->inc_low = 0.0; + key->inc_high = 0.0; return; } a_er = pow(WGS72_KE / n, 2.0 / 3.0); - range->low = a_er * (1.0 - e) * WGS72_AE - WGS72_AE; - range->high = a_er * (1.0 + e) * WGS72_AE - WGS72_AE; + key->alt_low = a_er * (1.0 - e) * WGS72_AE - WGS72_AE; + key->alt_high = a_er * (1.0 + e) * WGS72_AE - WGS72_AE; /* Guard against numerical inversion from near-zero eccentricity */ - if (range->low > range->high) + if (key->alt_low > key->alt_high) { - double tmp = range->low; - range->low = range->high; - range->high = tmp; + double tmp = key->alt_low; + key->alt_low = key->alt_high; + key->alt_high = tmp; } + + /* Leaf entry: exact inclination (radians) */ + key->inc_low = tle->inclination; + key->inc_high = tle->inclination; } /* ---------------------------------------------------------------- - * range_overlaps -- do two altitude bands share any interval? + * key_overlaps -- do two orbital keys overlap in BOTH dimensions? + * + * Altitude bands AND inclination ranges must both overlap. * ---------------------------------------------------------------- */ static inline bool -range_overlaps(const tle_alt_range *a, const tle_alt_range *b) +key_overlaps(const tle_orbital_key *a, const tle_orbital_key *b) { - return (a->low <= b->high) && (b->low <= a->high); + return (a->alt_low <= b->alt_high) && (b->alt_low <= a->alt_high) + && (a->inc_low <= b->inc_high) && (b->inc_low <= a->inc_high); } /* ---------------------------------------------------------------- - * range_contains -- does outer fully contain inner? + * key_contains -- does outer fully contain inner in both dimensions? * ---------------------------------------------------------------- */ static inline bool -range_contains(const tle_alt_range *outer, const tle_alt_range *inner) +key_contains(const tle_orbital_key *outer, const tle_orbital_key *inner) { - return (outer->low <= inner->low) && (inner->high <= outer->high); + return (outer->alt_low <= inner->alt_low) + && (inner->alt_high <= outer->alt_high) + && (outer->inc_low <= inner->inc_low) + && (inner->inc_high <= outer->inc_high); } /* ---------------------------------------------------------------- - * range_merge -- expand dst to encompass src + * key_merge -- expand dst bounding box to encompass src * ---------------------------------------------------------------- */ static inline void -range_merge(tle_alt_range *dst, const tle_alt_range *src) +key_merge(tle_orbital_key *dst, const tle_orbital_key *src) { - if (src->low < dst->low) - dst->low = src->low; - if (src->high > dst->high) - dst->high = src->high; + if (src->alt_low < dst->alt_low) + dst->alt_low = src->alt_low; + if (src->alt_high > dst->alt_high) + dst->alt_high = src->alt_high; + if (src->inc_low < dst->inc_low) + dst->inc_low = src->inc_low; + if (src->inc_high > dst->inc_high) + dst->inc_high = src->inc_high; } /* ---------------------------------------------------------------- - * range_separation -- minimum gap between two non-overlapping ranges + * alt_separation -- minimum altitude gap between two keys * - * Returns 0 if the ranges overlap. + * Returns 0 if the altitude bands overlap. + * Used for KNN distance (altitude-dominant ordering). * ---------------------------------------------------------------- */ static inline double -range_separation(const tle_alt_range *a, const tle_alt_range *b) +alt_separation(const tle_orbital_key *a, const tle_orbital_key *b) { - if (a->high < b->low) - return b->low - a->high; - if (b->high < a->low) - return a->low - b->high; + if (a->alt_high < b->alt_low) + return b->alt_low - a->alt_high; + if (b->alt_high < a->alt_low) + return a->alt_low - b->alt_high; return 0.0; } @@ -149,8 +177,8 @@ range_separation(const tle_alt_range *a, const tle_alt_range *b) /* * tle_overlap(tle, tle) -> bool [the && operator] * - * True if the altitude bands of two TLEs share any interval. - * This is a fast pre-filter: overlapping bands are necessary + * True if two TLEs overlap in both altitude AND inclination. + * This is a fast pre-filter: overlapping keys are necessary * but not sufficient for actual conjunction. */ Datum @@ -158,12 +186,12 @@ tle_overlap(PG_FUNCTION_ARGS) { pg_tle *a = (pg_tle *) PG_GETARG_POINTER(0); pg_tle *b = (pg_tle *) PG_GETARG_POINTER(1); - tle_alt_range ra, rb; + tle_orbital_key ka, kb; - tle_to_alt_range(a, &ra); - tle_to_alt_range(b, &rb); + tle_to_orbital_key(a, &ka); + tle_to_orbital_key(b, &kb); - PG_RETURN_BOOL(range_overlaps(&ra, &rb)); + PG_RETURN_BOOL(key_overlaps(&ka, &kb)); } @@ -174,18 +202,21 @@ tle_overlap(PG_FUNCTION_ARGS) * overlap. This is not the physical distance between the objects -- * it is the gap between their orbital shells, useful for ordering * nearest-neighbor queries without propagation. + * + * Altitude-only: inclination weighting adds complexity without + * meaningful benefit for KNN conjunction screening. */ Datum tle_alt_distance(PG_FUNCTION_ARGS) { pg_tle *a = (pg_tle *) PG_GETARG_POINTER(0); pg_tle *b = (pg_tle *) PG_GETARG_POINTER(1); - tle_alt_range ra, rb; + tle_orbital_key ka, kb; - tle_to_alt_range(a, &ra); - tle_to_alt_range(b, &rb); + tle_to_orbital_key(a, &ka); + tle_to_orbital_key(b, &kb); - PG_RETURN_FLOAT8(range_separation(&ra, &rb)); + PG_RETURN_FLOAT8(alt_separation(&ka, &kb)); } @@ -196,10 +227,10 @@ tle_alt_distance(PG_FUNCTION_ARGS) /* - * gist_tle_compress -- extract altitude range from a leaf TLE + * gist_tle_compress -- extract orbital key from a leaf TLE * - * Leaf entries carry the full pg_tle; we compress to tle_alt_range. - * Internal entries are already tle_alt_range from union operations. + * Leaf entries carry the full pg_tle; we compress to tle_orbital_key. + * Internal entries are already tle_orbital_key from union operations. */ Datum gist_tle_compress(PG_FUNCTION_ARGS) @@ -209,18 +240,18 @@ gist_tle_compress(PG_FUNCTION_ARGS) if (entry->leafkey) { - pg_tle *tle = (pg_tle *) DatumGetPointer(entry->key); - tle_alt_range *range = (tle_alt_range *) palloc(sizeof(tle_alt_range)); + pg_tle *tle = (pg_tle *) DatumGetPointer(entry->key); + tle_orbital_key *key = (tle_orbital_key *) palloc(sizeof(tle_orbital_key)); - tle_to_alt_range(tle, range); + tle_to_orbital_key(tle, key); retval = (GISTENTRY *) palloc(sizeof(GISTENTRY)); - gistentryinit(*retval, PointerGetDatum(range), + gistentryinit(*retval, PointerGetDatum(key), entry->rel, entry->page, entry->offset, false); } else { - /* Internal node: already a tle_alt_range */ + /* Internal node: already a tle_orbital_key */ retval = entry; } @@ -241,48 +272,44 @@ gist_tle_decompress(PG_FUNCTION_ARGS) /* * gist_tle_consistent -- can this subtree contain matches for the query? * - * Strategy RTOverlapStrategyNumber (&&): altitude bands must overlap. - * Always sets recheck = true because altitude overlap is only a necessary + * Checks overlap in both altitude AND inclination dimensions. + * Always sets recheck = true because 2-D overlap is only a necessary * condition -- the real conjunction test requires propagation. */ Datum gist_tle_consistent(PG_FUNCTION_ARGS) { - GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); - pg_tle *query = (pg_tle *) PG_GETARG_POINTER(1); - StrategyNumber strategy = (StrategyNumber) PG_GETARG_UINT16(2); + GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); + pg_tle *query = (pg_tle *) PG_GETARG_POINTER(1); + StrategyNumber strategy = (StrategyNumber) PG_GETARG_UINT16(2); /* arg 3 is the query type OID, unused */ - bool *recheck = (bool *) PG_GETARG_POINTER(4); - tle_alt_range *key = (tle_alt_range *) DatumGetPointer(entry->key); - tle_alt_range query_range; - bool result; + bool *recheck = (bool *) PG_GETARG_POINTER(4); + tle_orbital_key *key = (tle_orbital_key *) DatumGetPointer(entry->key); + tle_orbital_key query_key; + bool result; - tle_to_alt_range(query, &query_range); + tle_to_orbital_key(query, &query_key); - /* - * Altitude overlap is necessary, not sufficient. - * The actual operator (if exact) would need propagation, so always recheck. - */ *recheck = true; switch (strategy) { - case RTOverlapStrategyNumber: /* && */ - result = range_overlaps(key, &query_range); + case RTOverlapStrategyNumber: /* && */ + result = key_overlaps(key, &query_key); break; case RTContainedByStrategyNumber: /* <@ */ if (GIST_LEAF(entry)) - result = range_contains(&query_range, key); + result = key_contains(&query_key, key); else - result = range_overlaps(key, &query_range); + result = key_overlaps(key, &query_key); break; case RTContainsStrategyNumber: /* @> */ if (GIST_LEAF(entry)) - result = range_contains(key, &query_range); + result = key_contains(key, &query_key); else - result = range_overlaps(key, &query_range); + result = key_overlaps(key, &query_key); break; default: @@ -297,31 +324,30 @@ gist_tle_consistent(PG_FUNCTION_ARGS) /* - * gist_tle_union -- compute bounding altitude range for a set of entries + * gist_tle_union -- compute 2-D bounding box for a set of entries * - * The union of N altitude ranges is simply [min(low), max(high)]. + * The union is [min(alt_low), max(alt_high)] x [min(inc_low), max(inc_high)]. */ Datum gist_tle_union(PG_FUNCTION_ARGS) { - GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0); - int *sizep = (int *) PG_GETARG_POINTER(1); - int i; - tle_alt_range *result; - tle_alt_range *cur; + GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0); + int *sizep = (int *) PG_GETARG_POINTER(1); + int i; + tle_orbital_key *result; + tle_orbital_key *cur; - result = (tle_alt_range *) palloc(sizeof(tle_alt_range)); - cur = (tle_alt_range *) DatumGetPointer(entryvec->vector[0].key); - result->low = cur->low; - result->high = cur->high; + result = (tle_orbital_key *) palloc(sizeof(tle_orbital_key)); + cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[0].key); + *result = *cur; for (i = 1; i < entryvec->n; i++) { - cur = (tle_alt_range *) DatumGetPointer(entryvec->vector[i].key); - range_merge(result, cur); + cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[i].key); + key_merge(result, cur); } - *sizep = sizeof(tle_alt_range); + *sizep = sizeof(tle_orbital_key); PG_RETURN_POINTER(result); } @@ -330,26 +356,28 @@ gist_tle_union(PG_FUNCTION_ARGS) /* * gist_tle_penalty -- cost of inserting a new entry into an existing subtree * - * Penalty is the increase in altitude span (km) that results from - * expanding the subtree's bounding range to include the new entry. - * A good penalty function keeps the tree balanced and minimizes - * unnecessary range expansion. + * Uses margin (half-perimeter) instead of area. Leaf entries have + * inc_low == inc_high, giving zero area -- area-based penalty would + * degenerate to 0 for every insertion, making the tree random. + * Margin remains non-zero for degenerate (zero-width) bounding boxes. */ Datum gist_tle_penalty(PG_FUNCTION_ARGS) { - GISTENTRY *origentry = (GISTENTRY *) PG_GETARG_POINTER(0); - GISTENTRY *newentry = (GISTENTRY *) PG_GETARG_POINTER(1); - float *penalty = (float *) PG_GETARG_POINTER(2); - tle_alt_range *orig = (tle_alt_range *) DatumGetPointer(origentry->key); - tle_alt_range *add = (tle_alt_range *) DatumGetPointer(newentry->key); - double orig_span; - double merged_span; + GISTENTRY *origentry = (GISTENTRY *) PG_GETARG_POINTER(0); + GISTENTRY *newentry = (GISTENTRY *) PG_GETARG_POINTER(1); + float *penalty = (float *) PG_GETARG_POINTER(2); + tle_orbital_key *orig = (tle_orbital_key *) DatumGetPointer(origentry->key); + tle_orbital_key *add = (tle_orbital_key *) DatumGetPointer(newentry->key); + double orig_margin; + double merged_margin; - orig_span = orig->high - orig->low; - merged_span = fmax(orig->high, add->high) - fmin(orig->low, add->low); + orig_margin = (orig->alt_high - orig->alt_low) + + (orig->inc_high - orig->inc_low); + merged_margin = (fmax(orig->alt_high, add->alt_high) - fmin(orig->alt_low, add->alt_low)) + + (fmax(orig->inc_high, add->inc_high) - fmin(orig->inc_low, add->inc_low)); - *penalty = (float)(merged_span - orig_span); + *penalty = (float)(merged_margin - orig_margin); PG_RETURN_POINTER(penalty); } @@ -357,19 +385,19 @@ gist_tle_penalty(PG_FUNCTION_ARGS) /* * Comparison callback for qsort in picksplit. - * Sorts entries by the midpoint of their altitude range. + * Sorts entries by a sort key chosen at call time. */ typedef struct { int index; /* position in the original entry vector */ - double midpoint; /* (low + high) / 2 */ + double sortval; /* midpoint in the chosen dimension */ } picksplit_item; static int picksplit_cmp(const void *a, const void *b) { - double ma = ((const picksplit_item *) a)->midpoint; - double mb = ((const picksplit_item *) b)->midpoint; + double ma = ((const picksplit_item *) a)->sortval; + double mb = ((const picksplit_item *) b)->sortval; if (ma < mb) return -1; @@ -382,29 +410,66 @@ picksplit_cmp(const void *a, const void *b) /* * gist_tle_picksplit -- split an overfull page into two groups * - * Strategy: sort entries by altitude midpoint, split at the median. - * This keeps nearby altitude bands in the same subtree, which is - * optimal for a 1-D range index. + * 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. */ Datum gist_tle_picksplit(PG_FUNCTION_ARGS) { - GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0); - GIST_SPLITVEC *splitvec = (GIST_SPLITVEC *) PG_GETARG_POINTER(1); + GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0); + GIST_SPLITVEC *splitvec = (GIST_SPLITVEC *) PG_GETARG_POINTER(1); int nentries = entryvec->n; picksplit_item *items; - tle_alt_range *left_union; - tle_alt_range *right_union; - tle_alt_range *cur; + tle_orbital_key *left_union; + tle_orbital_key *right_union; + tle_orbital_key *cur; int split_at; int i; + 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); + 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++) + { + double alt_mid, inc_mid; + + cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[i].key); + alt_mid = (cur->alt_low + cur->alt_high) / 2.0; + inc_mid = (cur->inc_low + cur->inc_high) / 2.0; + + if (alt_mid < alt_min) alt_min = alt_mid; + if (alt_mid > alt_max) alt_max = alt_mid; + if (inc_mid < inc_min) inc_min = inc_mid; + if (inc_mid > inc_max) inc_max = inc_mid; + } + + /* + * Normalize spreads so they're comparable. Altitude is km above + * surface (range ~0-35786), inclination in radians (range 0-pi). + * Divide each by its domain width to get a 0-1 fraction. + */ + alt_spread = (alt_max - alt_min) / 35786.0; /* GEO altitude above surface */ + inc_spread = (inc_max - inc_min) / M_PI; + split_on_alt = (alt_spread >= inc_spread); + + /* Second pass: compute sort values in the chosen dimension */ items = (picksplit_item *) palloc(sizeof(picksplit_item) * nentries); for (i = 0; i < nentries; i++) { - cur = (tle_alt_range *) DatumGetPointer(entryvec->vector[i].key); - items[i].index = i; - items[i].midpoint = (cur->low + cur->high) / 2.0; + cur = (tle_orbital_key *) DatumGetPointer(entryvec->vector[i].key); + items[i].index = i; + if (split_on_alt) + items[i].sortval = (cur->alt_low + cur->alt_high) / 2.0; + else + items[i].sortval = (cur->inc_low + cur->inc_high) / 2.0; } qsort(items, nentries, sizeof(picksplit_item), picksplit_cmp); @@ -417,39 +482,37 @@ gist_tle_picksplit(PG_FUNCTION_ARGS) splitvec->spl_nleft = 0; splitvec->spl_nright = 0; - /* Compute union ranges and assign entries */ - left_union = (tle_alt_range *) palloc(sizeof(tle_alt_range)); - right_union = (tle_alt_range *) palloc(sizeof(tle_alt_range)); + /* Compute union keys and assign entries */ + left_union = (tle_orbital_key *) palloc(sizeof(tle_orbital_key)); + right_union = (tle_orbital_key *) palloc(sizeof(tle_orbital_key)); /* Seed the unions from the first entry in each half */ - cur = (tle_alt_range *) DatumGetPointer( + cur = (tle_orbital_key *) DatumGetPointer( entryvec->vector[items[0].index].key); - left_union->low = cur->low; - left_union->high = cur->high; + *left_union = *cur; - cur = (tle_alt_range *) DatumGetPointer( + cur = (tle_orbital_key *) DatumGetPointer( entryvec->vector[items[split_at].index].key); - right_union->low = cur->low; - right_union->high = cur->high; + *right_union = *cur; for (i = 0; i < nentries; i++) { int idx = items[i].index; - cur = (tle_alt_range *) DatumGetPointer( + cur = (tle_orbital_key *) DatumGetPointer( entryvec->vector[idx].key); if (i < split_at) { splitvec->spl_left[splitvec->spl_nleft++] = (OffsetNumber)(idx + 1); /* 1-based */ - range_merge(left_union, cur); + key_merge(left_union, cur); } else { splitvec->spl_right[splitvec->spl_nright++] = (OffsetNumber)(idx + 1); - range_merge(right_union, cur); + key_merge(right_union, cur); } } @@ -465,19 +528,20 @@ gist_tle_picksplit(PG_FUNCTION_ARGS) /* * gist_tle_same -- equality test on compressed keys * - * Two altitude ranges are "same" if both endpoints match within - * a small tolerance. This lets the GiST machinery detect - * duplicate keys. + * Two orbital keys are "same" if all four bounds match within + * tolerance. This lets the GiST machinery detect duplicate keys. */ Datum gist_tle_same(PG_FUNCTION_ARGS) { - tle_alt_range *a = (tle_alt_range *) PG_GETARG_POINTER(0); - tle_alt_range *b = (tle_alt_range *) PG_GETARG_POINTER(1); - bool *result = (bool *) PG_GETARG_POINTER(2); + tle_orbital_key *a = (tle_orbital_key *) PG_GETARG_POINTER(0); + tle_orbital_key *b = (tle_orbital_key *) PG_GETARG_POINTER(1); + bool *result = (bool *) PG_GETARG_POINTER(2); - *result = (fabs(a->low - b->low) < ALT_EPSILON && - fabs(a->high - b->high) < ALT_EPSILON); + *result = (fabs(a->alt_low - b->alt_low) < KEY_EPSILON + && fabs(a->alt_high - b->alt_high) < KEY_EPSILON + && fabs(a->inc_low - b->inc_low) < KEY_EPSILON + && fabs(a->inc_high - b->inc_high) < KEY_EPSILON); PG_RETURN_POINTER(result); } @@ -489,17 +553,20 @@ gist_tle_same(PG_FUNCTION_ARGS) * Returns the minimum altitude-band separation in km. * For overlapping ranges this is 0, making the entry a candidate. * The planner uses this to drive ORDER BY <-> queries. + * + * Altitude-only: conjunction screening is altitude-dominant. + * Inclination weighting can be added later if needed. */ Datum gist_tle_distance(PG_FUNCTION_ARGS) { - GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); - pg_tle *query = (pg_tle *) PG_GETARG_POINTER(1); + GISTENTRY *entry = (GISTENTRY *) PG_GETARG_POINTER(0); + pg_tle *query = (pg_tle *) PG_GETARG_POINTER(1); /* strategy number at arg 2, unused for single-distance class */ - tle_alt_range *key = (tle_alt_range *) DatumGetPointer(entry->key); - tle_alt_range query_range; + tle_orbital_key *key = (tle_orbital_key *) DatumGetPointer(entry->key); + tle_orbital_key query_key; - tle_to_alt_range(query, &query_range); + tle_to_orbital_key(query, &query_key); - PG_RETURN_FLOAT8(range_separation(key, &query_range)); + PG_RETURN_FLOAT8(alt_separation(key, &query_key)); } diff --git a/test/expected/gist_index.out b/test/expected/gist_index.out index 87d0505..e091787 100644 --- a/test/expected/gist_index.out +++ b/test/expected/gist_index.out @@ -1,4 +1,4 @@ --- Test GiST index and operators +-- Test GiST index and operators (2-D: altitude + inclination) CREATE EXTENSION IF NOT EXISTS pg_orbit; NOTICE: extension "pg_orbit" already exists, skipping -- Create test table with mixed orbit types @@ -7,46 +7,59 @@ CREATE TABLE test_orbits ( name text, tle tle ); --- ISS (LEO, ~400km) +-- ISS (LEO, ~400km, inclination 51.64 deg) INSERT INTO test_orbits (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) +-- Hubble (LEO, ~540km, inclination 28.47 deg) INSERT INTO test_orbits (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) +-- GPS IIR-M (MEO, ~20200km, inclination 55.44 deg) INSERT INTO test_orbits (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 but near-equatorial inclination (5 deg). +-- Same mean motion (15.501 rev/day) and eccentricity as ISS → same altitude band. +-- Different inclination → should NOT overlap with ISS under 2-D key. +INSERT INTO test_orbits (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'); -- Create GiST index CREATE INDEX test_orbits_gist ON test_orbits USING gist (tle); --- Overlap: ISS && Hubble = false (different altitude bands: ~400km vs ~540km) +-- 2-D overlap: ISS && Equatorial-LEO should be false (same altitude, different inclination) SELECT a.name AS sat_a, b.name AS sat_b, a.tle && b.tle AS overlaps FROM test_orbits a, test_orbits b WHERE a.id < b.id ORDER BY a.name, b.name; - sat_a | sat_b | overlaps ---------+---------+---------- - Hubble | GPS-IIR | f - ISS | GPS-IIR | f - ISS | Hubble | f -(3 rows) + sat_a | sat_b | overlaps +---------+----------------+---------- + GPS-IIR | Equatorial-LEO | f + Hubble | Equatorial-LEO | f + Hubble | GPS-IIR | f + ISS | Equatorial-LEO | f + ISS | GPS-IIR | f + ISS | Hubble | f +(6 rows) --- Altitude distance between different orbit regimes +-- Altitude distance: ISS <-> Equatorial-LEO should be ~0 (same altitude shell) SELECT a.name AS sat_a, b.name AS sat_b, round((a.tle <-> b.tle)::numeric, 0) AS alt_dist_km FROM test_orbits a, test_orbits b WHERE a.id < b.id ORDER BY a.name, b.name; - sat_a | sat_b | alt_dist_km ---------+---------+------------- - Hubble | GPS-IIR | 19332 - ISS | GPS-IIR | 19451 - ISS | Hubble | 115 -(3 rows) + sat_a | sat_b | alt_dist_km +---------+----------------+------------- + GPS-IIR | Equatorial-LEO | 19451 + Hubble | Equatorial-LEO | 115 + Hubble | GPS-IIR | 19332 + ISS | Equatorial-LEO | 0 + ISS | GPS-IIR | 19451 + ISS | Hubble | 115 +(6 rows) --- GiST index scan: find all LEO sats (overlap with ISS) +-- GiST index scan: find all sats overlapping ISS (altitude AND inclination) +-- Should return only ISS itself (Equatorial-LEO is at same altitude but wrong inclination) SET enable_seqscan = off; SELECT name FROM test_orbits WHERE tle && (SELECT tle FROM test_orbits WHERE name = 'ISS') @@ -63,21 +76,23 @@ SELECT name, round((tle <-> (SELECT tle FROM test_orbits WHERE name = 'ISS'))::n FROM test_orbits WHERE name != 'ISS' ORDER BY tle <-> (SELECT tle FROM test_orbits WHERE name = 'ISS'); - name | dist ----------+------- - Hubble | 115 - GPS-IIR | 19451 -(2 rows) + name | dist +----------------+------- + Equatorial-LEO | 0 + Hubble | 115 + GPS-IIR | 19451 +(3 rows) RESET enable_seqscan; -- Self-overlap is always true SELECT name, tle && tle AS self_overlap FROM test_orbits ORDER BY name; - name | self_overlap ----------+-------------- - GPS-IIR | t - Hubble | t - ISS | t -(3 rows) + name | self_overlap +----------------+-------------- + Equatorial-LEO | t + GPS-IIR | t + Hubble | t + ISS | t +(4 rows) -- Cleanup DROP TABLE test_orbits; diff --git a/test/sql/gist_index.sql b/test/sql/gist_index.sql index 09f50db..1d1c895 100644 --- a/test/sql/gist_index.sql +++ b/test/sql/gist_index.sql @@ -1,4 +1,4 @@ --- Test GiST index and operators +-- Test GiST index and operators (2-D: altitude + inclination) CREATE EXTENSION IF NOT EXISTS pg_orbit; -- Create test table with mixed orbit types @@ -8,38 +8,46 @@ CREATE TABLE test_orbits ( tle tle ); --- ISS (LEO, ~400km) +-- ISS (LEO, ~400km, inclination 51.64 deg) INSERT INTO test_orbits (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) +-- Hubble (LEO, ~540km, inclination 28.47 deg) INSERT INTO test_orbits (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) +-- GPS IIR-M (MEO, ~20200km, inclination 55.44 deg) INSERT INTO test_orbits (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 but near-equatorial inclination (5 deg). +-- Same mean motion (15.501 rev/day) and eccentricity as ISS → same altitude band. +-- Different inclination → should NOT overlap with ISS under 2-D key. +INSERT INTO test_orbits (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'); + -- Create GiST index CREATE INDEX test_orbits_gist ON test_orbits USING gist (tle); --- Overlap: ISS && Hubble = false (different altitude bands: ~400km vs ~540km) +-- 2-D overlap: ISS && Equatorial-LEO should be false (same altitude, different inclination) SELECT a.name AS sat_a, b.name AS sat_b, a.tle && b.tle AS overlaps FROM test_orbits a, test_orbits b WHERE a.id < b.id ORDER BY a.name, b.name; --- Altitude distance between different orbit regimes +-- Altitude distance: ISS <-> Equatorial-LEO should be ~0 (same altitude shell) SELECT a.name AS sat_a, b.name AS sat_b, round((a.tle <-> b.tle)::numeric, 0) AS alt_dist_km FROM test_orbits a, test_orbits b WHERE a.id < b.id ORDER BY a.name, b.name; --- GiST index scan: find all LEO sats (overlap with ISS) +-- GiST index scan: find all sats overlapping ISS (altitude AND inclination) +-- Should return only ISS itself (Equatorial-LEO is at same altitude but wrong inclination) SET enable_seqscan = off; SELECT name FROM test_orbits WHERE tle && (SELECT tle FROM test_orbits WHERE name = 'ISS')