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    /// Returns the polar motion rotation matrix for the given pole coordinates.
13    pub fn polar_motion_matrix(&self, time: Time<Tt>, pole_coords: PoleCoords) -> DMat3 {
14        if pole_coords.is_zero() {
15            return DMat3::IDENTITY;
16        };
17        match self {
18            ReferenceSystem::Iers1996 => pole_coords.polar_motion_matrix(),
19            ReferenceSystem::Iers2003(_) | ReferenceSystem::Iers2010 => {
20                let sp = TioLocator::iau2000(time);
21                pole_coords.polar_motion_matrix_iau2000(sp)
22            }
23        }
24    }
25}
26
27/// Pole coordinates (xp, yp) describing polar motion.
28#[derive(Debug, Clone, Copy, Default, PartialEq)]
29pub struct PoleCoords {
30    /// x-coordinate of the pole.
31    pub xp: Angle,
32    /// y-coordinate of the pole.
33    pub yp: Angle,
34}
35
36impl PoleCoords {
37    /// Computes the classical polar motion matrix W = R_y(-xp) * R_x(-yp).
38    pub fn polar_motion_matrix(&self) -> DMat3 {
39        (-self.xp).rotation_y() * (-self.yp).rotation_x()
40    }
41
42    /// Computes the polar motion matrix including the TIO locator s' (IAU 2000+).
43    pub fn polar_motion_matrix_iau2000(&self, sp: TioLocator) -> DMat3 {
44        self.polar_motion_matrix() * sp.0.rotation_z()
45    }
46
47    /// Returns `true` if both pole coordinates are zero.
48    pub fn is_zero(&self) -> bool {
49        self.xp.is_zero() && self.yp.is_zero()
50    }
51}
52
53#[cfg(test)]
54mod tests {
55    use lox_core::units::AngleUnits;
56    use lox_test_utils::assert_approx_eq;
57
58    use super::*;
59
60    #[test]
61    fn test_polar_motion_matrix_iau2000() {
62        let pole = PoleCoords {
63            xp: 2.55060238e-7.rad(),
64            yp: 1.860359247e-6.rad(),
65        };
66        let sp = TioLocator(-1.367_174_580_728_891_5e-11.rad());
67        let exp = DMat3::from_cols_array(&[
68            0.999_999_999_999_967_5,
69            -1.367_174_580_728_847e-11,
70            2.550_602_379_999_972e-7,
71            1.414_624_947_957_029_7e-11,
72            0.999_999_999_998_269_5,
73            -1.860_359_246_998_866_3e-6,
74            -2.550_602_379_741_215_3e-7,
75            1.860_359_247_002_414e-6,
76            0.999_999_999_998_237,
77        ])
78        .transpose();
79        let act = pole.polar_motion_matrix_iau2000(sp);
80        assert_approx_eq!(act, exp, atol <= 1e-12);
81    }
82}