Skip to main content

lox_frames/iers/
polar_motion.rs

1// SPDX-FileCopyrightText: 2025 Helge Eichhorn <git@helgeeichhorn.de>
2//
3// SPDX-License-Identifier: MPL-2.0
4
5use glam::DMat3;
6use lox_core::units::Angle;
7use lox_time::{Time, time_scales::Tt};
8
9use crate::iers::{ReferenceSystem, tio::TioLocator};
10
11impl ReferenceSystem {
12    pub fn polar_motion_matrix(&self, time: Time<Tt>, pole_coords: PoleCoords) -> DMat3 {
13        if pole_coords.is_zero() {
14            return DMat3::IDENTITY;
15        };
16        match self {
17            ReferenceSystem::Iers1996 => pole_coords.polar_motion_matrix(),
18            ReferenceSystem::Iers2003(_) | ReferenceSystem::Iers2010 => {
19                let sp = TioLocator::iau2000(time);
20                pole_coords.polar_motion_matrix_iau2000(sp)
21            }
22        }
23    }
24}
25
26#[derive(Debug, Clone, Copy, Default, PartialEq)]
27pub struct PoleCoords {
28    pub xp: Angle,
29    pub yp: Angle,
30}
31
32impl PoleCoords {
33    pub fn polar_motion_matrix(&self) -> DMat3 {
34        (-self.xp).rotation_y() * (-self.yp).rotation_x()
35    }
36
37    pub fn polar_motion_matrix_iau2000(&self, sp: TioLocator) -> DMat3 {
38        self.polar_motion_matrix() * sp.0.rotation_z()
39    }
40
41    pub fn is_zero(&self) -> bool {
42        self.xp.is_zero() && self.yp.is_zero()
43    }
44}
45
46#[cfg(test)]
47mod tests {
48    use lox_core::units::AngleUnits;
49    use lox_test_utils::assert_approx_eq;
50
51    use super::*;
52
53    #[test]
54    fn test_polar_motion_matrix_iau2000() {
55        let pole = PoleCoords {
56            xp: 2.55060238e-7.rad(),
57            yp: 1.860359247e-6.rad(),
58        };
59        let sp = TioLocator(-1.367_174_580_728_891_5e-11.rad());
60        let exp = DMat3::from_cols_array(&[
61            0.999_999_999_999_967_5,
62            -1.367_174_580_728_847e-11,
63            2.550_602_379_999_972e-7,
64            1.414_624_947_957_029_7e-11,
65            0.999_999_999_998_269_5,
66            -1.860_359_246_998_866_3e-6,
67            -2.550_602_379_741_215_3e-7,
68            1.860_359_247_002_414e-6,
69            0.999_999_999_998_237,
70        ])
71        .transpose();
72        let act = pole.polar_motion_matrix_iau2000(sp);
73        assert_approx_eq!(act, exp, atol <= 1e-12);
74    }
75}