Skip to content

Commit

Permalink
Implement 06 CIO locator calculation (#35)
Browse files Browse the repository at this point in the history
  • Loading branch information
@AngusGMorrison
AngusGMorrison authored Dec 16, 2023
1 parent 9d046a5 commit 7f637e3
Show file tree
Hide file tree
Showing 4 changed files with 254 additions and 0 deletions.
1 change: 1 addition & 0 deletions crates/lox_core/src/bodies.rs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@ use crate::time::constants::f64::{SECONDS_PER_DAY, SECONDS_PER_JULIAN_CENTURY};
mod generated;
pub use generated::*;

mod cio;
mod cip;
pub mod fundamental;
pub mod nutation;
Expand Down
12 changes: 12 additions & 0 deletions crates/lox_core/src/bodies/cio.rs
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,12 @@
/*
* Copyright (c) 2023. Helge Eichhorn and the LOX contributors
*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, you can obtain one at https://mozilla.org/MPL/2.0/.
*/

// TODO: Remove this once all module components are actively used.
#![allow(dead_code)]

mod s06;
123 changes: 123 additions & 0 deletions crates/lox_core/src/bodies/cio/s06.rs
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,123 @@
/*
* Copyright (c) 2023. Helge Eichhorn and the LOX contributors
*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, you can obtain one at https://mozilla.org/MPL/2.0/.
*/

mod terms;

use crate::bodies::cip::xy06::XY;
use crate::bodies::fundamental::iers03::{
general_accum_precession_in_longitude_iers03, mean_moon_sun_elongation_iers03,
};
use crate::bodies::{Earth, Moon, Sun, Venus};
use crate::math::arcsec_to_rad;
use crate::time::intervals::TDBJulianCenturiesSinceJ2000;
use crate::types::Radians;

/// l, l', F, D, Ω, LVe, LE and pA.
type FundamentalArgs = [Radians; 8];

/// Computes the Celestial Intermediate Origin (CIO) locator s, in radians, given the (X, Y) coordinates of
/// the Celestial Intermediate Pole (CIP). Based on IAU 2006 precession and IAU 2000A nutation.
pub fn s(t: TDBJulianCenturiesSinceJ2000, xy: XY) -> Radians {
let fundamental_args = fundamental_args(t);
let evaluated_terms = evaluate_terms(&fundamental_args);
let arcsec = fast_polynomial::poly_array(t, &evaluated_terms);
let radians = arcsec_to_rad(arcsec);
radians - xy[0] * xy[1] / 2.0
}

fn fundamental_args(t: TDBJulianCenturiesSinceJ2000) -> FundamentalArgs {
// The output of the CIO calculation is dependent on the ordering of these arguments. DO NOT
// EDIT.
[
Moon.mean_anomaly_iers03(t),
Sun.mean_anomaly_iers03(t),
Moon.mean_longitude_minus_ascending_node_mean_longitude_iers03(t),
mean_moon_sun_elongation_iers03(t),
Moon.ascending_node_mean_longitude_iers03(t),
Venus.mean_longitude_iers03(t),
Earth.mean_longitude_iers03(t),
general_accum_precession_in_longitude_iers03(t),
]
}

fn evaluate_terms(args: &FundamentalArgs) -> [f64; 6] {
[
evaluate_single_order_terms(args, terms::COEFFICIENTS[0], &terms::ZERO_ORDER),
evaluate_single_order_terms(args, terms::COEFFICIENTS[1], &terms::FIRST_ORDER),
evaluate_single_order_terms(args, terms::COEFFICIENTS[2], &terms::SECOND_ORDER),
evaluate_single_order_terms(args, terms::COEFFICIENTS[3], &terms::THIRD_ORDER),
evaluate_single_order_terms(args, terms::COEFFICIENTS[4], &terms::FOURTH_ORDER),
terms::COEFFICIENTS[5],
]
}

fn evaluate_single_order_terms(
args: &FundamentalArgs,
term_coefficient: f64,
terms: &[terms::Term],
) -> f64 {
terms.iter().rev().fold(term_coefficient, |acc, term| {
let a = term
.fundamental_arg_coeffs
.iter()
.zip(args)
.fold(0.0, |acc, (coeff, arg)| acc + coeff * arg);

acc + term.sin_coeff * a.sin() + term.cos_coeff * a.cos()
})
}

#[cfg(test)]
mod tests {
use super::*;
use crate::bodies::cip::xy06::xy;
use float_eq::assert_float_eq;

const TOLERANCE: f64 = 1e-11;

#[test]
fn test_s_jd0() {
let jd0: TDBJulianCenturiesSinceJ2000 = -67.11964407939767;
let xy = xy(jd0);
assert_float_eq!(s(jd0, xy), -0.0723985415686306, rel <= TOLERANCE);
}

#[test]
fn test_s_j2000() {
let j2000: TDBJulianCenturiesSinceJ2000 = 0.0;
let xy = xy(j2000);
assert_float_eq!(s(j2000, xy), -0.00000001013396519178, rel <= TOLERANCE);
}

#[test]
fn test_s_j2100() {
let j2100: TDBJulianCenturiesSinceJ2000 = 1.0;
let xy = xy(j2100);
assert_float_eq!(s(j2100, xy), -0.00000000480511934533, rel <= TOLERANCE);
}

#[test]
fn test_fundamental_args_ordering() {
let j2000: TDBJulianCenturiesSinceJ2000 = 0.0;
let actual = fundamental_args(j2000);
let expected = [
Moon.mean_anomaly_iers03(j2000),
Sun.mean_anomaly_iers03(j2000),
Moon.mean_longitude_minus_ascending_node_mean_longitude_iers03(j2000),
mean_moon_sun_elongation_iers03(j2000),
Moon.ascending_node_mean_longitude_iers03(j2000),
Venus.mean_longitude_iers03(j2000),
Earth.mean_longitude_iers03(j2000),
general_accum_precession_in_longitude_iers03(j2000),
];

expected.iter().enumerate().for_each(|(i, expected)| {
assert_float_eq!(*expected, actual[i], rel <= TOLERANCE);
});
}
}
118 changes: 118 additions & 0 deletions crates/lox_core/src/bodies/cio/s06/terms.rs
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,118 @@
/*
* Copyright (c) 2023. Helge Eichhorn and the LOX contributors
*
* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, you can obtain one at https://mozilla.org/MPL/2.0/.
*/

use crate::types::Radians;

pub(super) const COEFFICIENTS: [f64; 6] = [
94.00e-6,
3808.65e-6,
-122.68e-6,
-72574.11e-6,
27.98e-6,
15.62e-6,
];

/// Coefficients of l, l', F, D, Ω, LVe, LE and pA.
pub(super) type FundamentalArgCoeffs = [Radians; 8];

pub(super) struct Term {
pub fundamental_arg_coeffs: FundamentalArgCoeffs,
pub sin_coeff: f64,
pub cos_coeff: f64,
}

#[rustfmt::skip]
// @formatter:off (sometimes RustRover ignores the rustfmt skip)
pub(super) const ZERO_ORDER: [Term; 33] = [
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -2640.73e-6, cos_coeff: 0.39e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -63.53e-6, cos_coeff: 0.02e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 3.0, 0.0, 0.0, 0.0], sin_coeff: -11.75e-6, cos_coeff: -0.01e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -11.21e-6, cos_coeff: -0.01e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 4.57e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 3.0, 0.0, 0.0, 0.0], sin_coeff: -2.02e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -1.98e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 3.0, 0.0, 0.0, 0.0], sin_coeff: 1.72e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 1.41e-6, cos_coeff: 0.01e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: 1.26e-6, cos_coeff: 0.01e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: 0.63e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 0.63e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 2.0, -2.0, 3.0, 0.0, 0.0, 0.0], sin_coeff: -0.46e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 2.0, -2.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -0.45e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 4.0, -4.0, 4.0, 0.0, 0.0, 0.0], sin_coeff: -0.36e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 1.0, -1.0, 1.0, -8.0, 12.0, 0.0], sin_coeff: 0.24e-6, cos_coeff: 0.12e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -0.32e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -0.28e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 2.0, 0.0, 3.0, 0.0, 0.0, 0.0], sin_coeff: -0.27e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 2.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -0.26e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: 0.21e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, -2.0, 2.0, -3.0, 0.0, 0.0, 0.0], sin_coeff: -0.19e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, -2.0, 2.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: -0.18e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 0.0, 8.0,-13.0, -1.0], sin_coeff: 0.10e-6, cos_coeff: -0.05e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -0.15e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [2.0, 0.0, -2.0, 0.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: 0.14e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 2.0, -2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 0.14e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, -2.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -0.14e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, -2.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: -0.14e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 4.0, -2.0, 4.0, 0.0, 0.0, 0.0], sin_coeff: -0.13e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 4.0, 0.0, 0.0, 0.0], sin_coeff: 0.11e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, -2.0, 0.0, -3.0, 0.0, 0.0, 0.0], sin_coeff: -0.11e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, -2.0, 0.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: -0.11e-6, cos_coeff: 0.00e-6 },
];

#[rustfmt::skip]
// @formatter:off (sometimes RustRover ignores the rustfmt skip)
pub(super) const FIRST_ORDER: [Term; 3] = [
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -0.07e-6, cos_coeff: 3.57e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 1.73e-6, cos_coeff: -0.03e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 3.0, 0.0, 0.0, 0.0], sin_coeff: 0.00e-6, cos_coeff: 0.48e-6 },
];

#[rustfmt::skip]
// @formatter:off (sometimes RustRover ignores the rustfmt skip)
pub(super) const SECOND_ORDER: [Term; 25] = [
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 743.52e-6, cos_coeff: -0.17e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 56.91e-6, cos_coeff: 0.06e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 9.84e-6, cos_coeff: -0.01e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -8.85e-6, cos_coeff: 0.01e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -6.38e-6, cos_coeff: -0.05e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -3.07e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, 2.0, -2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 2.23e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 1.67e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 2.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 1.30e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 1.0, -2.0, 2.0, -2.0, 0.0, 0.0, 0.0], sin_coeff: 0.93e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: 0.68e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -0.55e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, -2.0, 0.0, -2.0, 0.0, 0.0, 0.0], sin_coeff: 0.53e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -0.27e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -0.27e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, -2.0, -2.0, -2.0, 0.0, 0.0, 0.0], sin_coeff: -0.26e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: -0.25e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 2.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 0.22e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [2.0, 0.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -0.21e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [2.0, 0.0, -2.0, 0.0, -1.0, 0.0, 0.0, 0.0], sin_coeff: 0.20e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 0.17e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [2.0, 0.0, 2.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 0.13e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [2.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -0.13e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [1.0, 0.0, 2.0, -2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -0.12e-6, cos_coeff: 0.00e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 0.0, 0.0], sin_coeff: -0.11e-6, cos_coeff: 0.00e-6 },
];

#[rustfmt::skip]
// @formatter:off (sometimes RustRover ignores the rustfmt skip)
pub(super) const THIRD_ORDER: [Term; 4] = [
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: 0.30e-6, cos_coeff: -23.42e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, -2.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -0.03e-6, cos_coeff: -1.46e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 2.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: -0.01e-6, cos_coeff: -0.25e-6 },
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0], sin_coeff: 0.00e-6, cos_coeff: 0.23e-6 },
];

#[rustfmt::skip]
// @formatter:off (sometimes RustRover ignores the rustfmt skip)
pub(super) const FOURTH_ORDER: [Term; 1] = [
Term{ fundamental_arg_coeffs: [0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0], sin_coeff: -0.26e-6, cos_coeff: -0.01e-6 }
];

0 comments on commit 7f637e3

Please sign in to comment.