kshana 0.25.0

Open, reproducible PNT-resilience simulator with quantum-sensor performance models
Documentation
// SPDX-License-Identifier: AGPL-3.0-only
//! D0.7 — Deep-space radiometric **light-time** reference test, gating kshana's
//! [`kshana::radiometric::light_time_solution`] against a committed **ANISE/SPICE DE440** oracle.
//!
//! ## What this gates, and how it stays ANISE-free
//!
//! The independent oracle is [ANISE](https://github.com/nyx-space/anise)'s converged-Newtonian
//! aberration light time (`Almanac::transform(.., Aberration::CN)`, a rewrite of NAIF SPICE's
//! `spkapo`) over JPL DE440. ANISE is MPL-2.0 / edition-2024 and cannot enter this crate's
//! dependency graph (it would break the `cargo deny` license gate and the MSRV-1.75 build), so it
//! lives in the workspace-EXCLUDED `xval/anise-mars-od` crate. That crate compiled+ran ANISE over
//! DE440 (`de440s.bsp`) and **pinned** its per-leg light times into the committed fixture
//! [`FIXTURE`] — together with, per leg, the exact DE440 SSB-relative geometry kshana needs:
//!
//! * the frozen observer (Earth) SSB position at the reception epoch (`obs_r`), and
//! * the target's SSB-relative quadratic Taylor coefficients about that epoch
//!   (`tgt_r`, `tgt_v`, `tgt_a`).
//!
//! This test reads only that text fixture — **no ANISE at test time** — reconstructs the target's
//! curved SSB motion `r(t) = tgt_r + tgt_v·Δ + ½·tgt_a·Δ²` (Δ = t − t_rx, seconds) via a local
//! [`EphemerisProvider`], re-runs kshana's *own* `light_time_solution` against it, and asserts the
//! resulting τ matches the committed ANISE τ to **≤ 1e-6 s** (sub-microsecond) on every leg.
//!
//! Because the `xval` crate proved (before writing the fixture) that this quadratic geometry
//! reproduces the full-DE440 kshana solve to < 1.1e-10 s — four orders below this gate — the gate
//! measures the genuine kshana-vs-ANISE light-time residual, not a degraded-fixture artefact.
//!
//! The fixture is regenerated by `cd xval/anise-mars-od && KSHANA_ANISE_DE440S=kernels/de440s.bsp
//! cargo run --release --bin lighttime-xval`.

use kshana::body::Body;
use kshana::ephem_provider::EphemerisProvider;
use kshana::radiometric::light_time_solution;
use kshana::timegeo::C_M_PER_S;
use kshana::timescales::{TwoPartJd, SECONDS_PER_DAY};

type Vec3 = [f64; 3];

/// The committed ANISE/DE440 light-time fixture (real `cargo run` output of the `xval` crate).
const FIXTURE: &str =
    include_str!("fixtures/deep_space_mars_radiometric/anise_lighttime_de440.txt");

/// The sub-microsecond gate: kshana's retarded τ must match the pinned ANISE converged-Newtonian τ
/// to within this bound on every leg.
const TAU_GATE_S: f64 = 1e-6;

/// One parsed fixture row: the pinned ANISE oracle plus the DE440 geometry kshana re-solves against.
#[derive(Debug, Clone)]
struct Leg {
    target: String,
    jd_tdb: f64,
    anise_tau_s: f64,
    /// Frozen observer (Earth) SSB position at the reception epoch (m).
    obs_r: Vec3,
    /// Target SSB Taylor coefficients about t_rx: position (m), velocity (m/s), acceleration (m/s²).
    tgt_r: Vec3,
    tgt_v: Vec3,
    tgt_a: Vec3,
}

/// An [`EphemerisProvider`] that returns the target's SSB-relative position from the committed
/// quadratic Taylor model about the reception epoch — the exact geometry ANISE's CN branch
/// differences, reconstructed with no ANISE present. Dispatch is on the target name; `center` is
/// ignored (positions are SSB-relative by construction, matching the fixture).
#[derive(Debug)]
struct TaylorProvider {
    target: String,
    t_rx_jd: f64,
    r0: Vec3,
    v: Vec3,
    a: Vec3,
}

impl EphemerisProvider for TaylorProvider {
    fn relative_position(&self, target: &Body, _center: &Body, jd_tdb: f64) -> Option<Vec3> {
        if target.name != self.target {
            return None;
        }
        let dt_s = (jd_tdb - self.t_rx_jd) * SECONDS_PER_DAY;
        Some([
            self.r0[0] + self.v[0] * dt_s + 0.5 * self.a[0] * dt_s * dt_s,
            self.r0[1] + self.v[1] * dt_s + 0.5 * self.a[1] * dt_s * dt_s,
            self.r0[2] + self.v[2] * dt_s + 0.5 * self.a[2] * dt_s * dt_s,
        ])
    }
}

/// Map a fixture target name to the kshana [`Body`] the provider dispatches on.
fn target_body(name: &str) -> Body {
    match name {
        "Mars" => Body::mars(),
        "Sun" => Body::sun(),
        "Moon" => Body::moon(),
        other => panic!("fixture has unexpected target {other}"),
    }
}

/// Parse the whitespace-delimited fixture rows (comment lines start with `#`). The column order is
/// documented in the fixture header and mirrored here exactly.
fn parse_legs() -> Vec<Leg> {
    let mut legs = Vec::new();
    for line in FIXTURE.lines() {
        let line = line.trim();
        if line.is_empty() || line.starts_with('#') {
            continue;
        }
        let f: Vec<&str> = line.split_whitespace().collect();
        assert_eq!(
            f.len(),
            20,
            "fixture row must have 20 columns, got {}: {line}",
            f.len()
        );
        let p = |i: usize| -> f64 {
            f[i].parse::<f64>()
                .unwrap_or_else(|_| panic!("non-numeric fixture field {i}: {:?}", f[i]))
        };
        legs.push(Leg {
            target: f[0].to_string(),
            jd_tdb: p(1),
            // f[2] kshana_tau (informational), f[3] anise_tau (the oracle we re-assert).
            anise_tau_s: p(3),
            obs_r: [p(8), p(9), p(10)],
            tgt_r: [p(11), p(12), p(13)],
            tgt_v: [p(14), p(15), p(16)],
            tgt_a: [p(17), p(18), p(19)], // f[19] = tgt_az, the last column.
        });
    }
    legs
}

/// The core gate: for every committed leg, re-run kshana's `light_time_solution` against the DE440
/// Taylor geometry and assert τ matches the pinned ANISE converged-Newtonian τ to ≤ 1e-6 s.
#[test]
fn kshana_light_time_matches_anise_de440_to_sub_microsecond() {
    let legs = parse_legs();
    assert!(
        legs.len() >= 8,
        "expected at least 8 reference legs, found {}",
        legs.len()
    );

    // The kshana `center` placeholder; the SSB-relative provider ignores it.
    let center = Body::earth();

    let mut worst_d_tau = 0.0_f64;
    let mut worst_leg = String::new();
    let mut n_mars = 0usize;

    for leg in &legs {
        let body = target_body(&leg.target);
        let provider = TaylorProvider {
            target: leg.target.clone(),
            t_rx_jd: leg.jd_tdb,
            r0: leg.tgt_r,
            v: leg.tgt_v,
            a: leg.tgt_a,
        };
        let t_rx = TwoPartJd::from_f64(leg.jd_tdb);

        let lt =
            light_time_solution(leg.obs_r, t_rx, &body, &center, &provider).unwrap_or_else(|| {
                panic!(
                    "kshana light_time_solution returned None for {} @ JD {}",
                    leg.target, leg.jd_tdb
                )
            });

        let d_tau = (lt.tau_s - leg.anise_tau_s).abs();
        if d_tau > worst_d_tau {
            worst_d_tau = d_tau;
            worst_leg = format!("{} @ JD {}", leg.target, leg.jd_tdb);
        }
        if leg.target == "Mars" {
            n_mars += 1;
        }

        assert!(
            d_tau <= TAU_GATE_S,
            "{} @ JD {}: kshana τ = {} s vs ANISE τ = {} s (Δ = {:.3e} s) exceeds the {:.0e} s gate",
            leg.target,
            leg.jd_tdb,
            lt.tau_s,
            leg.anise_tau_s,
            d_tau,
            TAU_GATE_S
        );

        // Self-consistency of the retarded solution: |obs − retarded target|/c == τ.
        let range_to_retarded = ((leg.obs_r[0] - lt.tx_pos[0]).powi(2)
            + (leg.obs_r[1] - lt.tx_pos[1]).powi(2)
            + (leg.obs_r[2] - lt.tx_pos[2]).powi(2))
        .sqrt()
            / C_M_PER_S;
        assert!(
            (range_to_retarded - lt.tau_s).abs() < 1e-6,
            "{} @ JD {}: retarded geometry inconsistent (|obs−tx|/c = {} s, τ = {} s)",
            leg.target,
            leg.jd_tdb,
            range_to_retarded,
            lt.tau_s
        );
    }

    // The Earth–Mars deep-space leg is the headline case the gate exists for.
    assert!(
        n_mars >= 4,
        "expected several Earth–Mars deep-space legs, found {n_mars}"
    );

    println!(
        "kshana vs ANISE CN/DE440 light time: {} legs, worst |Δτ| = {:.3e} s ({}), gate {:.0e} s — PASS",
        legs.len(),
        worst_d_tau,
        worst_leg,
        TAU_GATE_S
    );
}

/// Sanity bands on the pinned ANISE light times, independent of kshana — so a corrupted fixture is
/// caught regardless of the solver: Earth–Mars τ in the ~5–22 min range (≈ 0.6–2.5 AU), Earth–Sun
/// ≈ 1 AU/c ≈ 490–491 s, Earth–Moon ≈ 1.2–1.4 s.
#[test]
fn pinned_anise_light_times_are_physically_banded() {
    for leg in parse_legs() {
        match leg.target.as_str() {
            "Mars" => assert!(
                (300.0..=1300.0).contains(&leg.anise_tau_s),
                "Earth–Mars τ {} s out of the ~0.6–2.5 AU band",
                leg.anise_tau_s
            ),
            "Sun" => assert!(
                (488.0..=493.0).contains(&leg.anise_tau_s),
                "Earth–Sun τ {} s not near 1 AU/c",
                leg.anise_tau_s
            ),
            "Moon" => assert!(
                (1.15..=1.40).contains(&leg.anise_tau_s),
                "Earth–Moon τ {} s out of the perigee–apogee band",
                leg.anise_tau_s
            ),
            other => panic!("unexpected target {other}"),
        }
    }
}