astrodyn_time 0.1.0

Time scales (TAI/UTC/UT1/TDB/TT/GMST) and converters for the astrodyn orbital-dynamics pipeline
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148

//! TAI ↔ TDB time-scale conversion.
//!
//! Ports
//! [`models/environment/time/src/time_converter_tai_tdb.cc`](https://github.com/nasa/jeod/blob/jeod_v5.4.0/models/environment/time/src/time_converter_tai_tdb.cc)
//! from JEOD v5.4.0. Uses the simplified Astronomical-Almanac formula
//! for the TDB − TT offset (accurate to ±0.1 µs), then chains through
//! TT ↔ TAI in [`crate::time_converter_tai_tt`].

use crate::epoch::{J2000_TAI_TJT, TAI_TT_OFFSET};
use crate::time_converter_tai_tt::tai_to_tt;
use crate::{SecondsSince, TAI, TDB};
use std::f64::consts::PI;

/// Compute the TDB-TT offset in seconds at a given TAI truncated Julian time.
///
/// Uses the simplified formula from the Explanatory Supplement to the
/// Astronomical Almanac (2006), accurate to ±0.1 microseconds.
///
/// Ported from JEOD `time_converter_tai_tdb.cc`:
/// ```text
/// g = (π/180) * (357.53 + 0.9856003 * (tai_tjt - tai_tjt_at_J2000))
/// offset = 0.001658 * sin(g) + 0.000014 * sin(2g)
/// ```
pub fn tdb_tt_offset(tai_tjt: f64) -> f64 {
    let dt_days = tai_tjt - J2000_TAI_TJT;
    let g = (PI / 180.0) * (357.53 + 0.9856003 * dt_days);
    0.001658 * g.sin() + 0.000014 * (2.0 * g).sin()
}

/// Convert TAI seconds-since-epoch to TDB seconds-since-epoch.
///
/// TDB = TT + periodic_offset = TAI + 32.184 + periodic_offset
pub fn tai_to_tdb(tai_seconds: f64, tai_tjt: f64) -> f64 {
    let tt = tai_to_tt(tai_seconds);
    let offset = tdb_tt_offset(tai_tjt);
    tt + offset
}

/// Convert TDB seconds-since-epoch to TAI seconds-since-epoch.
///
/// Iterative solution (matches JEOD's 5-iteration convergence).
pub fn tdb_to_tai(tdb_seconds: f64, tai_tjt_initial: f64) -> f64 {
    let mut tai = tdb_seconds - TAI_TT_OFFSET;
    let mut tai_tjt = tai_tjt_initial;

    for _ in 0..5 {
        let offset = tdb_tt_offset(tai_tjt);
        let new_tai = tdb_seconds - TAI_TT_OFFSET - offset;
        let dtai = new_tai - tai;
        tai = new_tai;
        tai_tjt = tai_tjt_initial + dtai / 86400.0;
        if dtai.abs() < 1e-15 || (tai.abs() > 0.0 && (dtai / tai).abs() < 1.0e-15) {
            break;
        }
    }
    tai
}

// --- Typed variants (Phase 1 #103) ------------------------------------------
//
// Additive typed siblings. They delegate to the f64 free functions so
// behavior is bit-identical. `tai_tjt` is a scale-neutral truncated Julian
// time anchor carried through untouched.

/// Typed TAI → TDB converter (delegates to [`tai_to_tdb`]).
///
/// `tai_tjt` anchors the calendar date for the periodic TDB-TT offset; it is
/// the same scalar the f64 free function takes (truncated Julian time in TAI).
pub fn tai_to_tdb_typed(t: SecondsSince<TAI>, tai_tjt: f64) -> SecondsSince<TDB> {
    SecondsSince::from_seconds(tai_to_tdb(t.as_seconds(), tai_tjt))
}

/// Typed TDB → TAI converter (delegates to [`tdb_to_tai`]).
pub fn tdb_to_tai_typed(t: SecondsSince<TDB>, tai_tjt_initial: f64) -> SecondsSince<TAI> {
    SecondsSince::from_seconds(tdb_to_tai(t.as_seconds(), tai_tjt_initial))
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn tdb_offset_at_j2000() {
        // At J2000, g = 357.53 degrees, offset should be small
        let offset = tdb_tt_offset(J2000_TAI_TJT);
        // sin(357.53°) ≈ -0.0431, so offset ≈ -7.1e-5 s
        assert!(
            offset.abs() < 0.002,
            "TDB-TT offset at J2000: {} (should be < 2ms)",
            offset
        );
    }

    #[test]
    fn tai_tdb_round_trip() {
        let tai = 1_000_000.0;
        let tai_tjt = J2000_TAI_TJT + tai / 86400.0;
        let tdb = tai_to_tdb(tai, tai_tjt);
        let back = tdb_to_tai(tdb, tai_tjt);
        assert!(
            (back - tai).abs() < 1e-10,
            "TAI-TDB round trip: {} -> {} -> {}, err={}",
            tai,
            tdb,
            back,
            (back - tai).abs()
        );
    }

    #[test]
    fn tai_tdb_round_trip_near_zero() {
        // C4: convergence must work when tai ≈ 0 (near simulation start)
        let tai = 0.0;
        let tai_tjt = J2000_TAI_TJT;
        let tdb = tai_to_tdb(tai, tai_tjt);
        let back = tdb_to_tai(tdb, tai_tjt);
        assert!(
            (back - tai).abs() < 1e-10,
            "TAI-TDB round trip near zero: {} -> {} -> {}, err={}",
            tai,
            tdb,
            back,
            (back - tai).abs()
        );
    }

    #[test]
    fn tai_tdb_typed_matches_f64() {
        let tai = 1_000_000.0_f64;
        let tai_tjt = J2000_TAI_TJT + tai / 86400.0;
        let tdb_typed = tai_to_tdb_typed(SecondsSince::<TAI>::from_seconds(tai), tai_tjt);
        assert_eq!(tdb_typed.as_seconds(), tai_to_tdb(tai, tai_tjt));
    }

    #[test]
    fn tai_tdb_typed_round_trip() {
        // Non-trivial input, typed round-trip.
        let tai_s = 1_000_000.0_f64;
        let tai_tjt = J2000_TAI_TJT + tai_s / 86400.0;
        let tai = SecondsSince::<TAI>::from_seconds(tai_s);
        let tdb = tai_to_tdb_typed(tai, tai_tjt);
        let back = tdb_to_tai_typed(tdb, tai_tjt);
        let err = (back.as_seconds() - tai.as_seconds()).abs();
        // Same 1e-10 bound the iterative f64 round-trip achieves; well below
        // the 1e-14 TAI↔TT spec target (TDB adds ~ms-scale periodic noise).
        assert!(err < 1e-10, "TAI->TDB->TAI typed round-trip err={}", err);
    }
}