lox_frames/iers/cio/
iau2006.rs

1// SPDX-FileCopyrightText: 2023 Angus Morrison <github@angus-morrison.com>
2// SPDX-FileCopyrightText: 2024 Helge Eichhorn <git@helgeeichhorn.de>
3//
4// SPDX-License-Identifier: MPL-2.0
5
6//! Module s06 exposes a function for calculating the Celestial Intermediate Origin (CIO) locator,
7//! s, using IAU 2006 precession and IAU 2000A nutation.
8
9use fast_polynomial::poly_array;
10use lox_core::types::units::JulianCenturies;
11use lox_core::units::{Angle, AngleUnits};
12use lox_time::Time;
13use lox_time::julian_dates::JulianDate;
14use lox_time::time_scales::Tdb;
15
16use crate::iers::cio::CioLocator;
17use crate::iers::cip::CipCoords;
18use crate::iers::fundamental::iers03::{
19    d_iers03, earth_l_iers03, f_iers03, l_iers03, lp_iers03, omega_iers03, pa_iers03,
20    venus_l_iers03,
21};
22
23mod terms;
24
25impl CioLocator {
26    /// Computes the CIO locator using IAU 2006 precession and IAU 2000A nutation.
27    pub fn iau2006(time: Time<Tdb>, CipCoords { x, y }: CipCoords) -> Self {
28        let t = time.centuries_since_j2000();
29        let fundamental_args = fundamental_args(t);
30        let evaluated_terms = evaluate_terms(&fundamental_args);
31        let s = poly_array(t, &evaluated_terms).arcsec();
32        CioLocator(s - y.to_radians() / 2.0 * x)
33    }
34}
35
36/// l, l', F, D, Ω, LVe, LE and pA.
37type FundamentalArgs = [Angle; 8];
38
39/// Computes the Celestial Intermediate Origin (CIO) locator s, in radians, given the (X, Y)
40/// coordinates of the Celestial Intermediate Pole (CIP).
41pub fn cio_locator(
42    centuries_since_j2000_tdb: JulianCenturies,
43    CipCoords { x, y }: CipCoords,
44) -> Angle {
45    let fundamental_args = fundamental_args(centuries_since_j2000_tdb);
46    let evaluated_terms = evaluate_terms(&fundamental_args);
47    let s = poly_array(centuries_since_j2000_tdb, &evaluated_terms).arcsec();
48    Angle::radians(s.to_radians() - x.to_radians() * y.to_radians() / 2.0)
49}
50
51fn fundamental_args(centuries_since_j2000_tdb: JulianCenturies) -> FundamentalArgs {
52    // The output of the CIO calculation is dependent on the ordering of these arguments. DO NOT
53    // EDIT.
54    [
55        l_iers03(centuries_since_j2000_tdb),
56        lp_iers03(centuries_since_j2000_tdb),
57        f_iers03(centuries_since_j2000_tdb),
58        d_iers03(centuries_since_j2000_tdb),
59        omega_iers03(centuries_since_j2000_tdb),
60        venus_l_iers03(centuries_since_j2000_tdb),
61        earth_l_iers03(centuries_since_j2000_tdb),
62        pa_iers03(centuries_since_j2000_tdb),
63    ]
64}
65
66fn evaluate_terms(args: &FundamentalArgs) -> [f64; 6] {
67    [
68        evaluate_single_order_terms(args, terms::COEFFICIENTS[0], &terms::ZERO_ORDER),
69        evaluate_single_order_terms(args, terms::COEFFICIENTS[1], &terms::FIRST_ORDER),
70        evaluate_single_order_terms(args, terms::COEFFICIENTS[2], &terms::SECOND_ORDER),
71        evaluate_single_order_terms(args, terms::COEFFICIENTS[3], &terms::THIRD_ORDER),
72        evaluate_single_order_terms(args, terms::COEFFICIENTS[4], &terms::FOURTH_ORDER),
73        terms::COEFFICIENTS[5],
74    ]
75}
76
77fn evaluate_single_order_terms(
78    args: &FundamentalArgs,
79    term_coefficient: f64,
80    terms: &[terms::Term],
81) -> f64 {
82    terms.iter().rev().fold(term_coefficient, |acc, term| {
83        let a = term
84            .fundamental_arg_coeffs
85            .iter()
86            .zip(args)
87            .fold(0.0.rad(), |acc, (&coeff, &arg)| acc + coeff * arg);
88
89        acc + term.sin_coeff * a.sin() + term.cos_coeff * a.cos()
90    })
91}
92
93#[cfg(test)]
94mod tests {
95    use lox_test_utils::assert_approx_eq;
96
97    use super::*;
98
99    #[test]
100    fn test_cio_locator_iau2006() {
101        let cip_coords = CipCoords {
102            x: 5.791_308_486_706_011e-4.rad(),
103            y: 4.020_579_816_732_961e-5.rad(),
104        };
105        let time = Time::from_two_part_julian_date(Tdb, 2400000.5, 53736.0);
106        let exp = CioLocator(-1.220_032_213_076_463e-8.rad());
107        let act = CioLocator::iau2006(time, cip_coords);
108        assert_approx_eq!(act, exp, atol <= 1e-18);
109    }
110
111    #[test]
112    fn test_fundamental_args_ordering() {
113        let j2000: JulianCenturies = 0.0;
114        let act = fundamental_args(j2000);
115        let exp = [
116            l_iers03(j2000),
117            lp_iers03(j2000),
118            f_iers03(j2000),
119            d_iers03(j2000),
120            omega_iers03(j2000),
121            venus_l_iers03(j2000),
122            earth_l_iers03(j2000),
123            pa_iers03(j2000),
124        ];
125
126        assert_approx_eq!(act, exp, rtol <= 1e-12)
127    }
128}
lox_frames/iers/cio/iau2006.rs

lox_frames/iers/cio/
iau2006.rs
