lox_earth/
providers.rs

1// SPDX-FileCopyrightText: 2025 Helge Eichhorn <git@helgeeichhorn.de>
2//
3// SPDX-License-Identifier: MPL-2.0
4
5use crate::eop::{EopProvider, EopProviderError};
6use lox_frames::iers::polar_motion::PoleCoords;
7use lox_frames::iers::{Corrections, ReferenceSystem};
8use lox_frames::rotations::RotationProvider;
9use lox_time::Time;
10use lox_time::deltas::TimeDelta;
11use lox_time::offsets::{OffsetProvider, TryOffset};
12use lox_time::time_scales::{Tai, TimeScale};
13use lox_time::utc::Utc;
14use lox_time::utc::leap_seconds::{DefaultLeapSecondsProvider, LeapSecondsProvider};
15use lox_time::utc::transformations::ToUtc;
16
17impl LeapSecondsProvider for EopProvider {
18    fn delta_tai_utc(&self, tai: Time<Tai>) -> TimeDelta {
19        self.get_lsk().map_or_else(
20            || DefaultLeapSecondsProvider.delta_tai_utc(tai),
21            |lsk| lsk.delta_tai_utc(tai),
22        )
23    }
24
25    fn delta_utc_tai(&self, utc: Utc) -> TimeDelta {
26        self.get_lsk().map_or_else(
27            || DefaultLeapSecondsProvider.delta_utc_tai(utc),
28            |lsk| lsk.delta_utc_tai(utc),
29        )
30    }
31
32    fn is_leap_second_date(&self, date: lox_time::calendar_dates::Date) -> bool {
33        self.get_lsk().map_or_else(
34            || DefaultLeapSecondsProvider.is_leap_second_date(date),
35            |lsk| lsk.is_leap_second_date(date),
36        )
37    }
38
39    fn is_leap_second(&self, tai: Time<Tai>) -> bool {
40        self.get_lsk().map_or_else(
41            || DefaultLeapSecondsProvider.is_leap_second(tai),
42            |lsk| lsk.is_leap_second(tai),
43        )
44    }
45}
46
47impl OffsetProvider for EopProvider {
48    type Error = EopProviderError;
49
50    fn tai_to_ut1(&self, delta: TimeDelta) -> Result<TimeDelta, Self::Error> {
51        self.delta_ut1_tai(delta)
52    }
53
54    fn ut1_to_tai(&self, delta: TimeDelta) -> Result<TimeDelta, Self::Error> {
55        self.delta_tai_ut1(delta)
56    }
57}
58
59impl<T> RotationProvider<T> for EopProvider
60where
61    T: TimeScale + Copy,
62    Self: TryOffset<T, Tai>,
63{
64    type EopError = EopProviderError;
65
66    fn corrections(
67        &self,
68        time: Time<T>,
69        sys: ReferenceSystem,
70    ) -> Result<Corrections, EopProviderError> {
71        let utc = time
72            .try_to_scale(Tai, self)
73            .map_err(|err| EopProviderError::Offset(err.to_string()))?
74            .to_utc();
75        match sys {
76            ReferenceSystem::Iers1996 => self.nutation_precession_iau1980(utc),
77            ReferenceSystem::Iers2003(_) | ReferenceSystem::Iers2010 => {
78                self.nutation_precession_iau2000(utc)
79            }
80        }
81    }
82
83    fn pole_coords(&self, time: Time<T>) -> Result<PoleCoords, EopProviderError> {
84        let utc = time
85            .try_to_scale(Tai, self)
86            .map_err(|err| EopProviderError::Offset(err.to_string()))?
87            .to_utc();
88        self.polar_motion(utc)
89    }
90}
91
92#[cfg(test)]
93mod tests {
94    use super::*;
95
96    use std::sync::OnceLock;
97
98    use crate::eop::EopParser;
99
100    use lox_io::spice::lsk::LeapSecondsKernel;
101    use lox_test_utils::assert_approx_eq;
102    use lox_test_utils::data_file;
103    use lox_time::Time;
104    use lox_time::deltas::ToDelta;
105    use lox_time::offsets::TryOffset;
106    use lox_time::subsecond::Subsecond;
107    use lox_time::time;
108    use lox_time::time_scales::DynTimeScale;
109    use lox_time::time_scales::{Tai, Ut1};
110    use lox_time::utc::transformations::ToUtc;
111    use lox_time::{DynTime, calendar_dates::Date, time_of_day::TimeOfDay};
112    use rstest::{fixture, rstest};
113
114    #[rstest]
115    #[case(536414400, -36.40775963091942)]
116    #[case(536415400, -36.40776991477994)]
117    #[case(536416400, -36.407780212136835)]
118    #[case(536417400, -36.40779052210695)]
119    #[case(536418400, -36.407800844763194)]
120    #[case(536419400, -36.407811180178484)]
121    #[case(536420400, -36.40782152842571)]
122    #[case(536421400, -36.40783188957778)]
123    #[case(536422400, -36.4078422637076)]
124    #[case(536423400, -36.407852650888074)]
125    #[case(536424400, -36.4078630511921)]
126    #[case(536425400, -36.4078734646926)]
127    #[case(536426400, -36.40788389146246)]
128    #[case(536427400, -36.407894331574596)]
129    #[case(536428400, -36.407904785101906)]
130    #[case(536429400, -36.40791525211729)]
131    #[case(536430400, -36.40792573269367)]
132    #[case(536431400, -36.407936226903935)]
133    #[case(536432400, -36.40794673482099)]
134    #[case(536433400, -36.40795725651775)]
135    #[case(536434400, -36.407967792067105)]
136    #[case(536435400, -36.40797834154198)]
137    #[case(536436400, -36.40798890501525)]
138    #[case(536437400, -36.407999482559845)]
139    #[case(536438400, -36.40801007424866)]
140    #[case(536439400, -36.4080206801546)]
141    #[case(536440400, -36.40803130035057)]
142    #[case(536441400, -36.40804193490947)]
143    #[case(536442400, -36.408052583904215)]
144    #[case(536443400, -36.408063247407696)]
145    #[case(536444400, -36.40807392549283)]
146    #[case(536445400, -36.408084618232515)]
147    #[case(536446400, -36.40809532569965)]
148    #[case(536447400, -36.40810604796715)]
149    #[case(536448400, -36.408116785107914)]
150    #[case(536449400, -36.408127537194844)]
151    #[case(536450400, -36.408138304300856)]
152    #[case(536451400, -36.40814908649884)]
153    #[case(536452400, -36.40815988386171)]
154    #[case(536453400, -36.40817069646236)]
155    #[case(536454400, -36.40818152437371)]
156    #[case(536455400, -36.408192367668654)]
157    #[case(536456400, -36.4082032264201)]
158    #[case(536457400, -36.408214100700945)]
159    #[case(536458400, -36.4082249905841)]
160    #[case(536459400, -36.40823589614247)]
161    #[case(536460400, -36.40824681744896)]
162    #[case(536461400, -36.408257754576475)]
163    #[case(536462400, -36.40826870759791)]
164    #[case(536463400, -36.40827967658618)]
165    #[case(536464400, -36.408290661614195)]
166    #[case(536465400, -36.40830166275484)]
167    #[case(536466400, -36.40831268008103)]
168    #[case(536467400, -36.40832371366567)]
169    #[case(536468400, -36.40833476358167)]
170    #[case(536469400, -36.408345829901926)]
171    #[case(536470400, -36.40835691269934)]
172    #[case(536471400, -36.40836801204682)]
173    #[case(536472400, -36.40837912801728)]
174    #[case(536473400, -36.40839026068361)]
175    #[case(536474400, -36.40840141011872)]
176    #[case(536475400, -36.40841257639551)]
177    #[case(536476400, -36.4084237595869)]
178    #[case(536477400, -36.40843495976578)]
179    #[case(536478400, -36.408446177005054)]
180    #[case(536479400, -36.40845741137764)]
181    #[case(536480400, -36.40846866295642)]
182    #[case(536481400, -36.40847993181432)]
183    #[case(536482400, -36.40849121802424)]
184    #[case(536483400, -36.408502521659074)]
185    #[case(536484400, -36.40851384279173)]
186    #[case(536485400, -36.40852518149512)]
187    #[case(536486400, -36.40853653784214)]
188    #[case(536487400, -36.4085479119057)]
189    #[case(536488400, -36.40855930375871)]
190    #[case(536489400, -36.408570713474056)]
191    #[case(536490400, -36.40858214112466)]
192    #[case(536491400, -36.40859358678342)]
193    #[case(536492400, -36.408605050523235)]
194    #[case(536493400, -36.408616532417014)]
195    #[case(536494400, -36.40862803253767)]
196    #[case(536495400, -36.40863955095809)]
197    #[case(536496400, -36.408651087751196)]
198    #[case(536497400, -36.40866264298989)]
199    #[case(536498400, -36.40867421674706)]
200    #[case(536499400, -36.40868580909562)]
201    #[case(536500400, -36.40869742010849)]
202    fn test_delta_ut1_tai_orekit(
203        #[case] seconds: i64,
204        #[case] expected: f64,
205        provider: &EopProvider,
206    ) {
207        let tai = Time::new(Tai, seconds, Subsecond::default());
208        let ut1 = Time::new(Ut1, seconds, Subsecond::default());
209        let actual = provider
210            .delta_ut1_tai(tai.to_delta())
211            .unwrap()
212            .to_seconds()
213            .to_f64();
214        assert_approx_eq!(actual, expected, rtol <= 1e-6);
215        let actual = provider
216            .delta_tai_ut1(ut1.to_delta())
217            .unwrap()
218            .to_seconds()
219            .to_f64();
220        assert_approx_eq!(actual, -expected, rtol <= 1e-6);
221    }
222
223    #[rstest]
224    #[case(time!(Tai, 1973, 1, 1).unwrap())]
225    #[case(time!(Tai, 2100, 1, 1).unwrap())]
226    fn test_delta_ut1_tai_extrapolation(#[case] time: Time<Tai>, provider: &EopProvider) {
227        let act = provider.delta_ut1_tai(time.to_delta()).unwrap_err();
228        assert!(matches!(act, EopProviderError::ExtrapolatedValue(_)));
229        let ut1 = time.with_scale(Ut1);
230        let act = provider.delta_tai_ut1(ut1.to_delta()).unwrap_err();
231        assert!(matches!(act, EopProviderError::ExtrapolatedValue(_)));
232    }
233
234    const UT1_TOL: f64 = 1e-2;
235
236    // Reference values from Orekit
237    //
238    // Since we use a different algorithm for TCB UT1 we need to
239    // adjust the tolerance.
240    //
241    #[rstest]
242    #[case::tai_ut1("TAI", "UT1", -36.949521832072996)]
243    #[case::tcb_ut1("TCB", "UT1", -92.61803559995818)]
244    #[case::tcg_ut1("TCG", "UT1", -70.1891114139689)]
245    #[case::tdb_ut1("TDB", "UT1", -69.13340440689674)]
246    #[case::tt_ut1("TT", "UT1", -69.13352209269237)]
247    #[case::ut1_tai("UT1", "TAI", 36.949521532869305)]
248    #[case::ut1_tcb("UT1", "TCB", 92.61803631703046)]
249    #[case::ut1_tcg("UT1", "TCG", 70.18911089451464)]
250    #[case::ut1_tdb("UT1", "TDB", 69.13340387022173)]
251    #[case::ut1_tt("UT1", "TT", 69.13352153286931)]
252    #[case::ut1_ut1("UT1", "UT1", 0.0)]
253    fn test_dyn_time_scale_ut1(
254        #[case] scale1: &str,
255        #[case] scale2: &str,
256        #[case] exp: f64,
257        provider: &EopProvider,
258    ) {
259        let scale1: DynTimeScale = scale1.parse().unwrap();
260        let scale2: DynTimeScale = scale2.parse().unwrap();
261        let date = Date::new(2024, 12, 30).unwrap();
262        let time = TimeOfDay::from_hms(10, 27, 13.145).unwrap();
263        let dt = DynTime::from_date_and_time(scale1, date, time)
264            .unwrap()
265            .to_delta();
266        let act = provider
267            .try_offset(scale1, scale2, dt)
268            .unwrap()
269            .to_seconds()
270            .to_f64();
271        assert_approx_eq!(act, exp, rtol <= UT1_TOL);
272    }
273
274    #[test]
275    fn test_ut1_to_utc() {
276        let tai = time!(Tai, 2024, 5, 17, 12, 13, 14.0).unwrap();
277        let exp = tai.to_utc();
278        let ut1 = tai.try_to_scale(Ut1, provider()).unwrap();
279        let act = ut1.try_to_scale(Tai, provider()).unwrap().to_utc();
280        assert_eq!(act, exp);
281    }
282
283    #[fixture]
284    fn provider() -> &'static EopProvider {
285        static PROVIDER: OnceLock<EopProvider> = OnceLock::new();
286        PROVIDER.get_or_init(|| {
287            EopParser::new()
288                .from_paths(
289                    data_file("iers/finals.all.csv"),
290                    data_file("iers/finals2000A.all.csv"),
291                )
292                .with_leap_seconds_kernel(
293                    LeapSecondsKernel::from_file(data_file("spice/naif0012.tls")).unwrap(),
294                )
295                .parse()
296                .unwrap()
297        })
298    }
299}
lox_earth/providers.rs

lox_earth/
providers.rs
