sidereon_core/estimation/recipe.rs
1//! Named operation-order recipes for the GNSS estimation substrate.
2//!
3//! Phase-2 collapses the three thick estimator stacks (`spp`, `rtk`/`rtk_filter`,
4//! `precise_positioning`) onto one shared substrate plus thin, runtime-selectable
5//! strategies. The single hard constraint is that each external reference's
6//! bit-exactness (Skyfield for SPP, RTKLIB for RTK, the PPP oracle for PPP) must
7//! be preserved to 0 ULP. Different references need different floating-point
8//! operation orders for the *same* physical quantity, so the substrate never
9//! "simplifies" a parity-sensitive formula into one shared form. Instead every
10//! such choice is a NAMED variant: a strategy selects the op-order it needs by
11//! enum value rather than by owning a copy of the helper.
12//!
13//! This module *names* the recipes; the substrate and strategies route every
14//! caller through them. Each reference-faithful strategy resolves to the single
15//! op-order it was already using, so threading the recipe through the shared
16//! spine reproduces the prior code path bit-for-bit and leaves every existing
17//! 0-ULP golden unchanged.
18//!
19//! The `Canonical*` variants belong to the single consistent IERS-rigorous
20//! model (the bounded-tolerance canonical strategy, P6). They are NOT used by
21//! any reference-faithful strategy; canonical is an additional selectable
22//! strategy that changes nothing about the references. The SPP canonical range
23//! and frame variants ([`RangeRecipe::CanonicalLightTimeClosedFormSagnac`],
24//! [`FrameRecipe::CanonicalWgs84`]) are implemented and driven by
25//! [`EstimationRecipe::canonical_spp`]; the RTK and PPP canonical square-root
26//! solve ([`NormalRecipe::CanonicalSquareRoot`] on
27//! [`SolverRecipe::OwnedDeterministicCholesky`]) by
28//! [`EstimationRecipe::canonical_rtk`] and [`EstimationRecipe::canonical_ppp`].
29//! Canonical SPP, RTK, and PPP are all wired.
30
31/// Estimation technique: which physical observation model and parameter set a
32/// strategy estimates.
33#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
34pub enum Technique {
35 /// Single-point positioning: undifferenced pseudorange PVT.
36 #[default]
37 Spp,
38 /// Real-time kinematic: double-differenced code/phase baseline.
39 Rtk,
40 /// Precise point positioning: undifferenced ionosphere-free code/phase.
41 Ppp,
42}
43
44/// The reference a reference-faithful strategy is bit-exact against. The
45/// external oracles (Skyfield, RTKLIB, the PPP oracle) are CI validation targets
46/// whose goldens stay 0-ULP unchanged through P0-P5; [`Self::OwnedDeterministic`]
47/// is instead pinned to the owned solver's own frozen-bits golden (P5).
48#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
49pub enum ReferenceTarget {
50 /// Skyfield (the SPP geometry/clock/Sagnac reference).
51 #[default]
52 Skyfield,
53 /// RTKLIB (the RTK double-difference baseline reference).
54 Rtklib,
55 /// scipy least-squares host solve (the SPP solver-agreement reference).
56 /// Named for the `trust-region-least-squares` host-LAPACK fingerprint study;
57 /// not a runtime
58 /// estimation strategy (it is not wired into the SPP solve path), so it is
59 /// not a valid [`StrategyId`] target.
60 Scipy,
61 /// The PPP float/fixed oracle arc.
62 PppOracle,
63 /// The SPP owned deterministic trust-region solver
64 /// ([`SolverRecipe::OwnedDeterministicTrf`]): a fixed-reduction-order dense
65 /// subproblem factorization with no nalgebra LU and no black-box BLAS in that
66 /// solve. It is pinned to its own frozen-bits golden rather than to an
67 /// external library, and is selectable only for [`Technique::Spp`]. The owned
68 /// kernel owns only the subproblem factorization (the surrounding
69 /// normal-matrix / gradient / norm reductions stay on nalgebra), so its
70 /// cross-platform bit guarantee is scoped to the factorization; see
71 /// [`SolverRecipe::OwnedDeterministicTrf`] for the precise scope.
72 OwnedDeterministic,
73}
74
75/// Runtime-selectable strategy identity. `Reference` strategies are 0-ULP
76/// bit-exact to an external reference and remain the validation oracles;
77/// `Canonical` is the single bounded-tolerance "best" model (P6). Canonical SPP,
78/// RTK, and PPP are all wired.
79#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
80pub enum StrategyId {
81 /// A reference-faithful strategy: 0-ULP to `target` for `technique`.
82 Reference {
83 technique: Technique,
84 target: ReferenceTarget,
85 },
86 /// The canonical strategy for `technique` (bounded-tolerance, truth-gated).
87 Canonical { technique: Technique },
88}
89
90impl Default for StrategyId {
91 fn default() -> Self {
92 Self::Reference {
93 technique: Technique::Spp,
94 target: ReferenceTarget::Skyfield,
95 }
96 }
97}
98
99impl StrategyId {
100 /// SPP, bit-exact to Skyfield (`spp::solve`).
101 pub const fn spp_reference() -> Self {
102 Self::Reference {
103 technique: Technique::Spp,
104 target: ReferenceTarget::Skyfield,
105 }
106 }
107
108 /// RTK, bit-exact to RTKLIB (`rtk` / `rtk_filter`).
109 pub const fn rtk_reference() -> Self {
110 Self::Reference {
111 technique: Technique::Rtk,
112 target: ReferenceTarget::Rtklib,
113 }
114 }
115
116 /// PPP, bit-exact to the PPP oracle arc (`precise_positioning`).
117 pub const fn ppp_reference() -> Self {
118 Self::Reference {
119 technique: Technique::Ppp,
120 target: ReferenceTarget::PppOracle,
121 }
122 }
123
124 /// SPP via the owned deterministic trust-region solver
125 /// ([`SolverRecipe::OwnedDeterministicTrf`]): the owned dense subproblem
126 /// factorization, pinned to its own frozen-bits golden (its cross-platform
127 /// bit guarantee is scoped to the factorization). Selecting this through
128 /// [`crate::estimation::strategies::estimate`] drives the owned solver
129 /// rather than the legacy nalgebra LU path.
130 pub const fn spp_owned_deterministic() -> Self {
131 Self::Reference {
132 technique: Technique::Spp,
133 target: ReferenceTarget::OwnedDeterministic,
134 }
135 }
136}
137
138/// Geometric range / light-time / transmit-time operation order. Each variant
139/// names an existing range model; the substrate selects the op-order rather
140/// than copying the helper.
141#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
142pub enum RangeRecipe {
143 /// SPP closed-form light-time with a fixed transmit-time iteration count and
144 /// a measured-pseudorange seed (`spp/mod.rs` `sat_model`).
145 #[default]
146 SppMeasuredPseudorangeFixedIter,
147 /// `observables::predict` rounded-microsecond transmit time with a fixed
148 /// light-time iteration count (PPP / forward-prediction model).
149 ObservableRoundedMicrosecondFixedIter,
150 /// RTK provided-transmit-position range with the RTKLIB first-order Sagnac
151 /// scalar (`rtk_filter::model` line-of-sight / geometric range).
152 RtkProvidedTxFirstOrderSagnac,
153 /// Canonical: full iterative light-time (iterated to convergence, not a
154 /// fixed truncation) with the closed-form Sagnac Z-rotation, never a
155 /// first-order scalar Sagnac. Driven by [`EstimationRecipe::canonical_spp`]
156 /// in the SPP measurement model; not used by any reference strategy.
157 CanonicalLightTimeClosedFormSagnac,
158}
159
160/// Earth-rotation (Sagnac) correction operation order.
161#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
162pub enum SagnacRecipe {
163 /// Closed-form z-axis rotation of the satellite ECEF position by
164 /// `OMEGA_E_DOT * tau` (SPP / observables).
165 #[default]
166 ClosedFormZRotation,
167 /// RTKLIB first-order scalar Sagnac term added to the range
168 /// (`rtk_filter::model`).
169 RtklibFirstOrderScalar,
170 /// No Sagnac correction (synthetic / ECI-consistent inputs).
171 Off,
172}
173
174/// Local-frame / ENU / az-el basis construction operation order.
175#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
176pub enum FrameRecipe {
177 /// SPP Skyfield-parity ECEF->geodetic with the three-iteration AU-scaled
178 /// latitude solve (`spp` geodetic conversion).
179 #[default]
180 SppSkyfieldAuThreeIter,
181 /// Geocentric-up local frame used by the RTK elevation reference
182 /// (`rtk_filter` elevation mask / antenna projection).
183 GeocentricUpRtkReference,
184 /// Geodetic NEU basis built from the cross-product convention
185 /// (`precise_positioning::model` troposphere geometry).
186 GeodeticNeuCrossProduct,
187 /// DOP ENU rotation basis (`dop`).
188 DopEnuRotation,
189 /// Canonical: one consistent meters-native WGS84/ITRF geodetic basis under
190 /// IERS conventions (the core PROJ-pinned closed-form solve), never a
191 /// reference-specific AU-scaled path. Driven by
192 /// [`EstimationRecipe::canonical_spp`]; not used by any reference strategy.
193 CanonicalWgs84,
194}
195
196/// Normal-equation assembly tie-breaking / fold order. The tie order is the
197/// pivot/elimination convention that fixes the bit pattern of the reduced
198/// system.
199#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
200pub enum NormalRecipe {
201 /// SPP weighted-residual rows with a finite-difference design matrix
202 /// (`spp` least-squares).
203 #[default]
204 SppWeightedResidualFiniteDifference,
205 /// RTK double-difference block assembly with the first-tie covariance fold
206 /// (`rtk_filter::normal` first-tie block).
207 RtkDoubleDifferenceBlockFirstTie,
208 /// PPP weighted normal equations with epoch-local receiver clocks eliminated
209 /// and the reduced static system solved last-tie.
210 PppDenseLastTie,
211 /// Canonical square-root-information solve, shared by canonical RTK and
212 /// canonical PPP: the SPD normal system is solved by the owned deterministic
213 /// Cholesky factorization `Λ = L Lᵀ` plus forward/back substitution, where
214 /// `L` is the information-matrix square root. For RTK this is the
215 /// double-difference information system `Λ x = η` assembled by the same shared
216 /// block fold the RTK reference uses; for PPP it is the weighted normal
217 /// system assembled from the same undifferenced rows the PPP reference uses,
218 /// with epoch-local receiver clocks eliminated before factorization. This is
219 /// the numerically rigorous op-order for an SPD normal
220 /// matrix (no pivoting; exploits symmetry), distinct from the reference RTK
221 /// general first-tie Gaussian elimination
222 /// ([`Self::RtkDoubleDifferenceBlockFirstTie`]) and the reference PPP last-tie
223 /// Gaussian elimination ([`Self::PppDenseLastTie`]). Driven by
224 /// [`EstimationRecipe::canonical_rtk`] and [`EstimationRecipe::canonical_ppp`]
225 /// on the owned [`SolverRecipe::OwnedDeterministicCholesky`] kernel; not used
226 /// by any reference strategy.
227 CanonicalSquareRoot,
228}
229
230/// Linear-solve / factorization operation order. Determinism note: the legacy
231/// SPP path is nalgebra LU (not bit-portable end-to-end), preserved as a named
232/// variant; the owned deterministic kernel (P5) owns the dense subproblem
233/// factorization with its own goldens. Its determinism scope is the
234/// factorization, not the surrounding nalgebra reductions that build the
235/// subproblem -- see [`Self::OwnedDeterministicTrf`].
236#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
237pub enum SolverRecipe {
238 /// nalgebra trust-region least squares, the current SPP solver
239 /// (`spp` / `crate::astro::math::least_squares`). Existing SPP goldens use
240 /// this; kept unchanged.
241 #[default]
242 NalgebraTrfLegacy,
243 /// Flat first-tie Gaussian elimination (RTK baseline/filter solve).
244 FlatGaussianFirstTie,
245 /// Dense last-tie Gaussian elimination (PPP solve,
246 /// `crate::astro::math::linear::solve_linear_last_tie`).
247 DenseGaussianLastTie,
248 /// scipy host LAPACK reference solve (machine-dependent; only as a
249 /// fingerprinted CI reference, never canonical).
250 ScipyHostLapackReference,
251 /// Owned deterministic Cholesky (square-root) linear solve, the canonical RTK
252 /// (P6 increment 2) and canonical PPP (P6 increment 3) solver: the SPD normal
253 /// system is factored `Λ = L Lᵀ` and solved by forward/back substitution
254 /// through the owned
255 /// [`crate::astro::math::linear::solve_flat_normal_square_root_into`] kernel,
256 /// with no nalgebra LU and no black-box BLAS. Paired with
257 /// [`NormalRecipe::CanonicalSquareRoot`]. Both the RTK and PPP canonical paths
258 /// are owned scalar arithmetic and f64 sqrt is IEEE-754 correctly rounded, so
259 /// unlike [`Self::OwnedDeterministicTrf`] (whose surrounding reductions ride
260 /// nalgebra) its bit guarantee covers the full solve and is portable across
261 /// platforms.
262 OwnedDeterministicCholesky,
263 /// Owned deterministic trust-region subproblem solve added in P5: a
264 /// fixed-reduction-order dense Gaussian elimination (the
265 /// `OwnedGaussianFirstTie` kernel) with no nalgebra LU and no black-box BLAS
266 /// in the factorization, pinned to its OWN frozen-bits goldens. Scope: it
267 /// owns ONLY the subproblem factorization; the normal-matrix / gradient /
268 /// norm reductions that build the subproblem still flow through nalgebra's
269 /// CPU-dispatched dense algebra, so the cross-platform bit guarantee is
270 /// scoped to the factorization, not the full solve.
271 OwnedDeterministicTrf,
272}
273
274/// The full operation-order recipe a strategy composes: one variant per stage.
275/// `Default` and the named constructors reproduce the CURRENT behavior of each
276/// existing strategy, so selecting a recipe never changes a reference golden
277/// (PPP goldens were re-frozen once when static PPP moved to clock-eliminated
278/// reduced normals; see `estimation::strategies`).
279#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash, Default)]
280pub struct EstimationRecipe {
281 pub range: RangeRecipe,
282 pub sagnac: SagnacRecipe,
283 pub frame: FrameRecipe,
284 pub normal: NormalRecipe,
285 pub solver: SolverRecipe,
286}
287
288impl EstimationRecipe {
289 /// The current SPP reference recipe (`spp::solve`, Skyfield-parity).
290 pub const fn spp() -> Self {
291 Self {
292 range: RangeRecipe::SppMeasuredPseudorangeFixedIter,
293 sagnac: SagnacRecipe::ClosedFormZRotation,
294 frame: FrameRecipe::SppSkyfieldAuThreeIter,
295 normal: NormalRecipe::SppWeightedResidualFiniteDifference,
296 solver: SolverRecipe::NalgebraTrfLegacy,
297 }
298 }
299
300 /// The current RTK reference recipe (`rtk` / `rtk_filter`, RTKLIB-parity).
301 pub const fn rtk() -> Self {
302 Self {
303 range: RangeRecipe::RtkProvidedTxFirstOrderSagnac,
304 sagnac: SagnacRecipe::RtklibFirstOrderScalar,
305 frame: FrameRecipe::GeocentricUpRtkReference,
306 normal: NormalRecipe::RtkDoubleDifferenceBlockFirstTie,
307 solver: SolverRecipe::FlatGaussianFirstTie,
308 }
309 }
310
311 /// The current PPP reference recipe (`precise_positioning`, oracle-parity).
312 pub const fn ppp() -> Self {
313 Self {
314 range: RangeRecipe::ObservableRoundedMicrosecondFixedIter,
315 sagnac: SagnacRecipe::ClosedFormZRotation,
316 frame: FrameRecipe::GeodeticNeuCrossProduct,
317 normal: NormalRecipe::PppDenseLastTie,
318 solver: SolverRecipe::DenseGaussianLastTie,
319 }
320 }
321
322 /// The SPP recipe driving the owned deterministic trust-region solver: the
323 /// SPP reference model with [`SolverRecipe::OwnedDeterministicTrf`] swapped
324 /// in for the legacy nalgebra LU linear-solve stage. Every other stage is the
325 /// SPP reference op-order, so only the factorization changes.
326 pub const fn spp_owned_deterministic() -> Self {
327 let mut recipe = Self::spp();
328 recipe.solver = SolverRecipe::OwnedDeterministicTrf;
329 recipe
330 }
331
332 /// The canonical SPP recipe: the single consistent IERS-rigorous SPP
333 /// measurement model. It diverges from [`Self::spp`] (the Skyfield-faithful
334 /// reference) only where the physics says to:
335 /// - range: [`RangeRecipe::CanonicalLightTimeClosedFormSagnac`] iterates the
336 /// light-time loop to convergence (vs the reference's fixed
337 /// transmit-time truncation), with the closed-form Sagnac Z-rotation (never
338 /// a first-order scalar Sagnac).
339 /// - frame: [`FrameRecipe::CanonicalWgs84`] solves ECEF->geodetic directly in
340 /// meters on the WGS84 ellipsoid (vs the reference's Skyfield AU-scaled
341 /// three-iteration latitude loop).
342 /// - solver: [`SolverRecipe::OwnedDeterministicTrf`] owns the trust-region
343 /// subproblem factorization so canonical is deterministic run-to-run on a
344 /// pinned build (its cross-platform bit guarantee is scoped to the
345 /// factorization; the surrounding reductions ride nalgebra).
346 ///
347 /// The Sagnac stage is the closed-form Z-rotation the SPP reference already
348 /// uses (the rigorous form), and the normal stage is the SPP
349 /// weighted-residual finite-difference assembly the trust-region solver
350 /// consumes; neither needs a separate canonical variant for SPP.
351 pub const fn canonical_spp() -> Self {
352 Self {
353 range: RangeRecipe::CanonicalLightTimeClosedFormSagnac,
354 sagnac: SagnacRecipe::ClosedFormZRotation,
355 frame: FrameRecipe::CanonicalWgs84,
356 normal: NormalRecipe::SppWeightedResidualFiniteDifference,
357 solver: SolverRecipe::OwnedDeterministicTrf,
358 }
359 }
360
361 /// The canonical RTK recipe: the double-difference baseline under the
362 /// numerically rigorous square-root-information solve. It keeps the RTK
363 /// reference's double-difference measurement physics (the provided-transmit
364 /// range with the RTKLIB first-order Sagnac scalar, the geocentric-up
365 /// elevation frame), because the canonical RTK divergence the physics calls
366 /// for is in the linear algebra, not the observation model: the same SPD
367 /// information system the reference assembles is solved by the owned
368 /// deterministic Cholesky square-root factorization
369 /// ([`NormalRecipe::CanonicalSquareRoot`] on
370 /// [`SolverRecipe::OwnedDeterministicCholesky`]) instead of the reference's
371 /// general first-tie Gaussian elimination. The square-root solve needs no
372 /// pivoting, exploits the symmetry of the SPD normal matrix, and is entirely
373 /// owned scalar arithmetic (no nalgebra, no BLAS), so canonical RTK is
374 /// well-conditioned and bit-reproducible across platforms.
375 pub const fn canonical_rtk() -> Self {
376 Self {
377 range: RangeRecipe::RtkProvidedTxFirstOrderSagnac,
378 sagnac: SagnacRecipe::RtklibFirstOrderScalar,
379 frame: FrameRecipe::GeocentricUpRtkReference,
380 normal: NormalRecipe::CanonicalSquareRoot,
381 solver: SolverRecipe::OwnedDeterministicCholesky,
382 }
383 }
384
385 /// The canonical PPP recipe: the undifferenced ionosphere-free PPP arc under
386 /// the numerically rigorous square-root-information solve. Like
387 /// [`Self::canonical_rtk`] it keeps the PPP reference's measurement physics
388 /// (the rounded-microsecond fixed-iteration light-time with the rigorous
389 /// closed-form Sagnac Z-rotation, and the geodetic NEU antenna frame), because
390 /// the canonical PPP divergence the physics calls for is in the linear
391 /// algebra, not the observation model: the same weighted normal equations
392 /// the reference assembles from the undifferenced rows are reduced by
393 /// eliminating epoch-local receiver clocks, then solved by the owned
394 /// deterministic Cholesky square-root factorization
395 /// ([`NormalRecipe::CanonicalSquareRoot`] on
396 /// [`SolverRecipe::OwnedDeterministicCholesky`]) instead of the reference's
397 /// dense last-tie Gaussian elimination ([`SolverRecipe::DenseGaussianLastTie`]).
398 /// The square-root solve needs no pivoting, exploits the symmetry of the SPD
399 /// normal matrix, and is entirely owned scalar arithmetic (no nalgebra, no
400 /// BLAS), so it is well-conditioned and the solve itself is bit-portable.
401 /// Determinism scope (calibrated, not overstated): unlike canonical RTK, the PPP
402 /// measurement model that builds the rows evaluates troposphere mapping,
403 /// antenna, and geodetic-frame transcendentals through the platform math
404 /// library, so canonical PPP's overall output is bit-reproducible run-to-run on
405 /// a pinned build but is not claimed bit-portable across platforms; only the
406 /// owned Cholesky solve carries the cross-platform guarantee.
407 pub const fn canonical_ppp() -> Self {
408 Self {
409 range: RangeRecipe::ObservableRoundedMicrosecondFixedIter,
410 sagnac: SagnacRecipe::ClosedFormZRotation,
411 frame: FrameRecipe::GeodeticNeuCrossProduct,
412 normal: NormalRecipe::CanonicalSquareRoot,
413 solver: SolverRecipe::OwnedDeterministicCholesky,
414 }
415 }
416
417 /// The canonical recipe for a `technique`. Canonical SPP (P6 increment 1),
418 /// canonical RTK (P6 increment 2), and canonical PPP (P6 increment 3) are all
419 /// wired, so every technique has a canonical strategy. Returns `Option` to keep
420 /// the resolver's "not yet implemented" surface stable.
421 pub const fn for_canonical(technique: Technique) -> Option<Self> {
422 match technique {
423 Technique::Spp => Some(Self::canonical_spp()),
424 Technique::Rtk => Some(Self::canonical_rtk()),
425 Technique::Ppp => Some(Self::canonical_ppp()),
426 }
427 }
428
429 /// The reference recipe for an explicit `(technique, target)` pair, or `None`
430 /// if the pair is not a supported reference strategy. This is the single
431 /// source of truth for which targets each technique can run: only the wired
432 /// reference oracles (Skyfield for SPP, RTKLIB for RTK, the PPP oracle for
433 /// PPP) and the SPP owned deterministic solver are valid. Every other pair
434 /// (a cross-technique oracle, or the unwired scipy host-LAPACK reference) is
435 /// rejected so an impossible strategy can never silently run a mismatched
436 /// recipe.
437 pub const fn for_reference(technique: Technique, target: ReferenceTarget) -> Option<Self> {
438 match (technique, target) {
439 (Technique::Spp, ReferenceTarget::Skyfield) => Some(Self::spp()),
440 (Technique::Spp, ReferenceTarget::OwnedDeterministic) => {
441 Some(Self::spp_owned_deterministic())
442 }
443 (Technique::Rtk, ReferenceTarget::Rtklib) => Some(Self::rtk()),
444 (Technique::Ppp, ReferenceTarget::PppOracle) => Some(Self::ppp()),
445 _ => None,
446 }
447 }
448}
449
450/// How a strategy forms its integer-ambiguity identifiers, and against what they
451/// are referenced. Naming this lets the RTK and PPP fixed solvers share one
452/// LAMBDA resolution kernel
453/// ([`crate::estimation::substrate::ambiguity::resolve_integer_lattice`]) and
454/// differ only in DATA rather than in separate algorithm trees.
455#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
456pub enum DifferencingMode {
457 /// Double-differenced ambiguities, one reference satellite per constellation
458 /// (the RTK baseline / sequential-filter convention: each non-reference
459 /// satellite is differenced against its own system's reference).
460 DoubleDifferencePerSystemReference,
461 /// Undifferenced ambiguities, one per satellite per receiver (the PPP
462 /// convention: no reference satellite, all satellites carry their own
463 /// ionosphere-free ambiguity).
464 Undifferenced,
465}
466
467/// Whether partial ambiguity resolution is attempted when the full-set integer
468/// fix fails its ratio test, and with what floor on the retained subset size.
469#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
470pub enum PartialResolution {
471 /// Full-set only: a failed ratio test means "not fixed" (PPP, and the RTK
472 /// sequential filter, both take the full set or nothing).
473 Disabled,
474 /// Confidence-ranked then exhaustive subset fallback down to
475 /// `min_ambiguities` retained (the RTK static fixed solver,
476 /// `rtk_filter::search::search_partial_fixed_ambiguities`).
477 Exhaustive { min_ambiguities: usize },
478}
479
480/// The integer-ambiguity identity/eligibility policy a strategy resolves under:
481/// the strategy DATA that replaces the RTK-vs-PPP algorithm-tree split. The
482/// LAMBDA resolution kernel is common; only these fields differ between the
483/// reference strategies. Named in P3; consumed by the runtime selector in P4.
484#[derive(Debug, Clone, Copy, PartialEq)]
485pub struct AmbiguityIdPolicy {
486 pub differencing: DifferencingMode,
487 /// Exclude float-only constellations from the integer search set.
488 pub float_only_gating: bool,
489 pub partial: PartialResolution,
490 /// Ratio-test acceptance threshold passed to the LAMBDA kernel.
491 pub ratio_threshold: f64,
492}
493
494impl AmbiguityIdPolicy {
495 /// The static RTK fixed-baseline policy (`rtk_filter::fixed`): per-system
496 /// double differences, float-only constellations excluded from the search,
497 /// partial resolution down to `partial_min_ambiguities`.
498 pub const fn rtk_static(ratio_threshold: f64, partial_min_ambiguities: usize) -> Self {
499 Self {
500 differencing: DifferencingMode::DoubleDifferencePerSystemReference,
501 float_only_gating: true,
502 partial: PartialResolution::Exhaustive {
503 min_ambiguities: partial_min_ambiguities,
504 },
505 ratio_threshold,
506 }
507 }
508
509 /// The sequential RTK filter policy (`rtk_filter::update`): per-system double
510 /// differences, float-only constellations excluded, full set or nothing.
511 pub const fn rtk_sequential(ratio_threshold: f64) -> Self {
512 Self {
513 differencing: DifferencingMode::DoubleDifferencePerSystemReference,
514 float_only_gating: true,
515 partial: PartialResolution::Disabled,
516 ratio_threshold,
517 }
518 }
519
520 /// The static PPP fixed policy (`precise_positioning::fixed`): undifferenced
521 /// per-satellite ambiguities, no constellation gating, full set or nothing.
522 pub const fn ppp(ratio_threshold: f64) -> Self {
523 Self {
524 differencing: DifferencingMode::Undifferenced,
525 float_only_gating: false,
526 partial: PartialResolution::Disabled,
527 ratio_threshold,
528 }
529 }
530}
531
532/// The operation order used to normalize one residual against its weight before
533/// the sigma comparison in a per-residual screen. Naming the order keeps each
534/// screen bit-identical while the formula lives in exactly one place
535/// ([`crate::estimation::substrate::qc::normalized_residual`]).
536#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
537pub enum ResidualNormRecipe {
538 /// `value · weight` where `weight` is an inverse double-difference *variance*
539 /// (`1/(sigma_sat^2 + sigma_ref^2)`), so the normalized innovation is
540 /// `value / sigma^2`. The RTK sequential information filter, whose DD rows
541 /// weight by inverse variance, screens its predicted innovations this way.
542 RtkInverseVarianceInnovation,
543 /// `value · weight` where `weight` is an inverse *sigma*
544 /// (`1/sqrt(sigma_sat^2 + sigma_ref^2)`), so the normalized residual is
545 /// `value / sigma`. The RTK static float/fixed least-squares baselines, whose
546 /// DD rows weight by inverse sigma, screen their post-fit residuals this way.
547 RtkInverseSigmaResidual,
548 /// `|value| · sqrt(weight)` where `weight` is an inverse *sigma* (`1/sigma`):
549 /// the residual magnitude times the square root of the inverse-sigma weight.
550 /// The PPP float leave-one-out screen (PPP rows weight by inverse sigma, as
551 /// `MeasurementWeights` documents).
552 PppInverseSigmaMagnitude,
553}
554
555/// The residual-screen family a strategy applies after (or, for the filter,
556/// before) a solve. Strategy DATA for the P4 selector; the chi-square variant is
557/// the SPP RAIM aggregate test, the rest are per-residual sigma screens that
558/// share [`crate::estimation::substrate::qc::normalized_residual`].
559#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
560pub enum ScreenKind {
561 /// SPP RAIM: aggregate chi-square on the weighted residual sum, then FDE
562 /// leave-one-out exclusion (`quality::raim` / `quality::fde`).
563 RaimChiSquare,
564 /// RTK static fixed: worst information-weighted residual vs a sigma gate,
565 /// excluding the worst satellite within a budget
566 /// (`rtk_filter::fixed::solve_fixed_baseline_validated`).
567 RtkFixedResidualValidation,
568 /// RTK sequential filter: information-weighted innovation screen on predicted
569 /// DD rows, masking rejected rows and coasting (`rtk_filter::update`).
570 RtkSequentialInnovation,
571 /// PPP float: worst studentized residual vs a sigma gate, leave-one-out prune
572 /// and re-solve while WRMS improves (`precise_positioning::float`).
573 PppFloatLeaveOneOut,
574}
575
576impl ScreenKind {
577 /// The per-residual normalization op-order this screen uses, or `None` for
578 /// the aggregate chi-square RAIM screen (which scores the weighted residual
579 /// sum, not individual residuals).
580 pub const fn residual_norm(self) -> Option<ResidualNormRecipe> {
581 match self {
582 Self::RaimChiSquare => None,
583 Self::RtkFixedResidualValidation => Some(ResidualNormRecipe::RtkInverseSigmaResidual),
584 Self::RtkSequentialInnovation => Some(ResidualNormRecipe::RtkInverseVarianceInnovation),
585 Self::PppFloatLeaveOneOut => Some(ResidualNormRecipe::PppInverseSigmaMagnitude),
586 }
587 }
588}
589
590#[cfg(test)]
591mod tests {
592 use super::*;
593
594 #[test]
595 fn defaults_name_current_spp_behavior() {
596 // The per-stage defaults are the SPP reference op-orders, so an
597 // unspecified recipe reproduces the current SPP path.
598 assert_eq!(EstimationRecipe::default(), EstimationRecipe::spp());
599 assert_eq!(
600 RangeRecipe::default(),
601 RangeRecipe::SppMeasuredPseudorangeFixedIter
602 );
603 assert_eq!(SagnacRecipe::default(), SagnacRecipe::ClosedFormZRotation);
604 assert_eq!(FrameRecipe::default(), FrameRecipe::SppSkyfieldAuThreeIter);
605 assert_eq!(
606 NormalRecipe::default(),
607 NormalRecipe::SppWeightedResidualFiniteDifference
608 );
609 assert_eq!(SolverRecipe::default(), SolverRecipe::NalgebraTrfLegacy);
610 assert_eq!(StrategyId::default(), StrategyId::spp_reference());
611 }
612
613 #[test]
614 fn strategy_constructors_match_reference_targets() {
615 assert_eq!(
616 StrategyId::spp_reference(),
617 StrategyId::Reference {
618 technique: Technique::Spp,
619 target: ReferenceTarget::Skyfield,
620 }
621 );
622 assert_eq!(
623 StrategyId::rtk_reference(),
624 StrategyId::Reference {
625 technique: Technique::Rtk,
626 target: ReferenceTarget::Rtklib,
627 }
628 );
629 assert_eq!(
630 StrategyId::ppp_reference(),
631 StrategyId::Reference {
632 technique: Technique::Ppp,
633 target: ReferenceTarget::PppOracle,
634 }
635 );
636 }
637
638 #[test]
639 fn for_reference_selects_each_supported_pairs_recipe() {
640 assert_eq!(
641 EstimationRecipe::for_reference(Technique::Spp, ReferenceTarget::Skyfield),
642 Some(EstimationRecipe::spp())
643 );
644 assert_eq!(
645 EstimationRecipe::for_reference(Technique::Rtk, ReferenceTarget::Rtklib),
646 Some(EstimationRecipe::rtk())
647 );
648 assert_eq!(
649 EstimationRecipe::for_reference(Technique::Ppp, ReferenceTarget::PppOracle),
650 Some(EstimationRecipe::ppp())
651 );
652 }
653
654 #[test]
655 fn owned_deterministic_recipe_swaps_only_the_solver() {
656 let owned = EstimationRecipe::spp_owned_deterministic();
657 assert_eq!(owned.solver, SolverRecipe::OwnedDeterministicTrf);
658 // Every non-solver stage is the SPP reference op-order.
659 assert_eq!(
660 EstimationRecipe {
661 solver: SolverRecipe::NalgebraTrfLegacy,
662 ..owned
663 },
664 EstimationRecipe::spp()
665 );
666 assert_eq!(
667 EstimationRecipe::for_reference(Technique::Spp, ReferenceTarget::OwnedDeterministic),
668 Some(owned)
669 );
670 }
671
672 #[test]
673 fn canonical_spp_recipe_uses_the_rigorous_op_orders() {
674 let canonical = EstimationRecipe::canonical_spp();
675 // Range: full iterative light-time with closed-form Sagnac, not the SPP
676 // reference's fixed-iteration measured-pseudorange recipe.
677 assert_eq!(
678 canonical.range,
679 RangeRecipe::CanonicalLightTimeClosedFormSagnac
680 );
681 assert_ne!(canonical.range, EstimationRecipe::spp().range);
682 // Frame: one consistent meters-native WGS84 basis, not the Skyfield AU
683 // path.
684 assert_eq!(canonical.frame, FrameRecipe::CanonicalWgs84);
685 assert_ne!(canonical.frame, EstimationRecipe::spp().frame);
686 // Sagnac stays the closed-form Z-rotation (the rigorous form the SPP
687 // reference already uses); the canonical divergence is never a
688 // first-order scalar Sagnac.
689 assert_eq!(canonical.sagnac, SagnacRecipe::ClosedFormZRotation);
690 assert_ne!(canonical.sagnac, SagnacRecipe::RtklibFirstOrderScalar);
691 // Solver: the owned deterministic factorization, for run-to-run
692 // determinism on a pinned build.
693 assert_eq!(canonical.solver, SolverRecipe::OwnedDeterministicTrf);
694 assert_eq!(
695 EstimationRecipe::for_canonical(Technique::Spp),
696 Some(canonical)
697 );
698 }
699
700 #[test]
701 fn canonical_rtk_recipe_uses_the_square_root_solve() {
702 let canonical = EstimationRecipe::canonical_rtk();
703 // Normal + solver: the owned Cholesky square-root information solve, not
704 // the reference RTK first-tie Gaussian elimination.
705 assert_eq!(canonical.normal, NormalRecipe::CanonicalSquareRoot);
706 assert_eq!(canonical.solver, SolverRecipe::OwnedDeterministicCholesky);
707 assert_ne!(canonical.normal, EstimationRecipe::rtk().normal);
708 assert_ne!(canonical.solver, EstimationRecipe::rtk().solver);
709 // Measurement physics stays the RTK reference double-difference model: the
710 // canonical RTK divergence is in the linear algebra, not the observation
711 // model, so range/sagnac/frame match the reference.
712 assert_eq!(canonical.range, EstimationRecipe::rtk().range);
713 assert_eq!(canonical.sagnac, EstimationRecipe::rtk().sagnac);
714 assert_eq!(canonical.frame, EstimationRecipe::rtk().frame);
715 assert_eq!(
716 EstimationRecipe::for_canonical(Technique::Rtk),
717 Some(canonical)
718 );
719 }
720
721 #[test]
722 fn canonical_ppp_recipe_uses_the_square_root_solve() {
723 let canonical = EstimationRecipe::canonical_ppp();
724 // Normal + solver: the owned Cholesky square-root information solve, not
725 // the reference PPP dense last-tie Gaussian elimination.
726 assert_eq!(canonical.normal, NormalRecipe::CanonicalSquareRoot);
727 assert_eq!(canonical.solver, SolverRecipe::OwnedDeterministicCholesky);
728 assert_ne!(canonical.normal, EstimationRecipe::ppp().normal);
729 assert_ne!(canonical.solver, EstimationRecipe::ppp().solver);
730 // Measurement physics stays the PPP reference undifferenced model: the
731 // canonical PPP divergence is in the linear algebra, not the observation
732 // model, so range/sagnac/frame match the reference.
733 assert_eq!(canonical.range, EstimationRecipe::ppp().range);
734 assert_eq!(canonical.sagnac, EstimationRecipe::ppp().sagnac);
735 assert_eq!(canonical.frame, EstimationRecipe::ppp().frame);
736 // Canonical RTK and PPP share the square-root normal + owned Cholesky
737 // solver (the same numerically rigorous SPD op-order).
738 assert_eq!(canonical.normal, EstimationRecipe::canonical_rtk().normal);
739 assert_eq!(canonical.solver, EstimationRecipe::canonical_rtk().solver);
740 assert_eq!(
741 EstimationRecipe::for_canonical(Technique::Ppp),
742 Some(canonical)
743 );
744 }
745
746 #[test]
747 fn for_canonical_wires_all_three_techniques() {
748 assert_eq!(
749 EstimationRecipe::for_canonical(Technique::Spp),
750 Some(EstimationRecipe::canonical_spp())
751 );
752 assert_eq!(
753 EstimationRecipe::for_canonical(Technique::Rtk),
754 Some(EstimationRecipe::canonical_rtk())
755 );
756 assert_eq!(
757 EstimationRecipe::for_canonical(Technique::Ppp),
758 Some(EstimationRecipe::canonical_ppp())
759 );
760 }
761
762 #[test]
763 fn for_reference_rejects_impossible_pairs() {
764 // Cross-technique oracles and the unwired scipy reference are not
765 // supported reference strategies.
766 for (technique, target) in [
767 (Technique::Spp, ReferenceTarget::Rtklib),
768 (Technique::Spp, ReferenceTarget::PppOracle),
769 (Technique::Spp, ReferenceTarget::Scipy),
770 (Technique::Rtk, ReferenceTarget::Skyfield),
771 (Technique::Rtk, ReferenceTarget::OwnedDeterministic),
772 (Technique::Rtk, ReferenceTarget::PppOracle),
773 (Technique::Ppp, ReferenceTarget::Skyfield),
774 (Technique::Ppp, ReferenceTarget::OwnedDeterministic),
775 ] {
776 assert_eq!(
777 EstimationRecipe::for_reference(technique, target),
778 None,
779 "{technique:?} + {target:?} must be rejected"
780 );
781 }
782 }
783
784 #[test]
785 fn reference_ambiguity_policies_name_current_behavior() {
786 let rtk_static = AmbiguityIdPolicy::rtk_static(3.0, 4);
787 assert_eq!(
788 rtk_static.differencing,
789 DifferencingMode::DoubleDifferencePerSystemReference
790 );
791 assert!(rtk_static.float_only_gating);
792 assert_eq!(
793 rtk_static.partial,
794 PartialResolution::Exhaustive { min_ambiguities: 4 }
795 );
796
797 let rtk_seq = AmbiguityIdPolicy::rtk_sequential(3.0);
798 assert_eq!(
799 rtk_seq.differencing,
800 DifferencingMode::DoubleDifferencePerSystemReference
801 );
802 assert!(rtk_seq.float_only_gating);
803 assert_eq!(rtk_seq.partial, PartialResolution::Disabled);
804
805 let ppp = AmbiguityIdPolicy::ppp(2.5);
806 assert_eq!(ppp.differencing, DifferencingMode::Undifferenced);
807 assert!(!ppp.float_only_gating);
808 assert_eq!(ppp.partial, PartialResolution::Disabled);
809 }
810
811 #[test]
812 fn rtk_and_ppp_id_policies_differ_only_in_data() {
813 // Same LAMBDA kernel, different identity/eligibility data: the two stacks
814 // are no longer separate algorithm trees, only different policy values.
815 let rtk = AmbiguityIdPolicy::rtk_static(3.0, 1);
816 let ppp = AmbiguityIdPolicy::ppp(3.0);
817 assert_ne!(rtk.differencing, ppp.differencing);
818 assert_ne!(rtk.float_only_gating, ppp.float_only_gating);
819 assert_ne!(rtk.partial, ppp.partial);
820 }
821
822 #[test]
823 fn screen_kinds_select_their_normalization_order() {
824 assert_eq!(ScreenKind::RaimChiSquare.residual_norm(), None);
825 assert_eq!(
826 ScreenKind::RtkFixedResidualValidation.residual_norm(),
827 Some(ResidualNormRecipe::RtkInverseSigmaResidual)
828 );
829 assert_eq!(
830 ScreenKind::RtkSequentialInnovation.residual_norm(),
831 Some(ResidualNormRecipe::RtkInverseVarianceInnovation)
832 );
833 assert_eq!(
834 ScreenKind::PppFloatLeaveOneOut.residual_norm(),
835 Some(ResidualNormRecipe::PppInverseSigmaMagnitude)
836 );
837 }
838
839 #[test]
840 fn each_strategy_selects_a_distinct_solver_order() {
841 // The three reference strategies must not collapse onto one solver
842 // op-order; that distinction is what preserves their separate goldens.
843 assert_ne!(
844 EstimationRecipe::spp().solver,
845 EstimationRecipe::rtk().solver
846 );
847 assert_ne!(
848 EstimationRecipe::rtk().solver,
849 EstimationRecipe::ppp().solver
850 );
851 assert_ne!(
852 EstimationRecipe::spp().solver,
853 EstimationRecipe::ppp().solver
854 );
855 }
856}