WIP

crates/lox-earth/src/lib.rs

@@ -6,13 +6,10 @@
  * 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 cio;
-mod cip;
-mod coordinate_transformations;
-mod nutation;
-mod rotation_angle;
-mod tides;
-mod tio;
+pub mod cio;
+pub mod cip;
+pub mod coordinate_transformations;
+pub mod nutation;
+pub mod rotation_angle;
+pub mod tides;
+pub mod tio;

crates/lox-orbits/Cargo.toml

@@ -10,6 +10,7 @@ repository.workspace = true
 
 [dependencies]
 lox-bodies.workspace = true
+lox-earth.workspace = true
 lox-ephem.workspace = true
 lox-time.workspace = true
 lox-math.workspace = true

crates/lox-orbits/src/elements.rs

@@ -11,21 +11,28 @@ use std::f64::consts::TAU;
 use float_eq::float_eq;
 use glam::{DMat3, DVec3};
 
-use lox_bodies::{DynOrigin, Origin, PointMass};
+use lox_bodies::{DynOrigin, MaybePointMass, Origin, PointMass};
 use lox_time::deltas::TimeDelta;
 use lox_time::TimeLike;
 
-use crate::frames::{CoordinateSystem, Icrf};
-use crate::states::{State, ToCartesian};
+use crate::frames::{CoordinateSystem, DynFrame, Icrf, ReferenceFrame};
+use crate::states::State;
 
-pub trait ToKeplerian<T: TimeLike, O: PointMass> {
-    fn to_keplerian(&self) -> Keplerian<T, O>;
+#[derive(Debug, Clone, PartialEq)]
+pub(crate) struct KeplerianElements {
+    pub semi_major_axis: f64,
+    pub eccentricity: f64,
+    pub inclination: f64,
+    pub longitude_of_ascending_node: f64,
+    pub argument_of_periapsis: f64,
+    pub true_anomaly: f64,
 }
 
 #[derive(Debug, Clone, PartialEq)]
-pub struct Keplerian<T: TimeLike, O: Origin> {
+pub struct Keplerian<T: TimeLike, O: MaybePointMass, R: ReferenceFrame> {
     time: T,
     origin: O,
+    frame: R,
     semi_major_axis: f64,
     eccentricity: f64,
     inclination: f64,
@@ -34,12 +41,12 @@ pub struct Keplerian<T: TimeLike, O: Origin> {
     true_anomaly: f64,
 }
 
-pub type DynKeplerian<T> = Keplerian<T, DynOrigin>;
+pub type DynKeplerian<T> = Keplerian<T, DynOrigin, DynFrame>;
 
-impl<T, O> Keplerian<T, O>
+impl<T, O> Keplerian<T, O, Icrf>
 where
     T: TimeLike,
-    O: Origin,
+    O: PointMass,
 {
     #[allow(clippy::too_many_arguments)]
     pub fn new(
@@ -55,6 +62,7 @@ where
         Self {
             time,
             origin,
+            frame: Icrf,
             semi_major_axis,
             eccentricity,
             inclination,
@@ -63,7 +71,46 @@ where
             true_anomaly,
         }
     }
+}
 
+impl<T> DynKeplerian<T>
+where
+    T: TimeLike,
+{
+    #[allow(clippy::too_many_arguments)]
+    pub fn with_dynamic(
+        time: T,
+        origin: DynOrigin,
+        semi_major_axis: f64,
+        eccentricity: f64,
+        inclination: f64,
+        longitude_of_ascending_node: f64,
+        argument_of_periapsis: f64,
+        true_anomaly: f64,
+    ) -> Result<Self, &'static str> {
+        if origin.maybe_gravitational_parameter().is_none() {
+            return Err("undefined gravitational parameter");
+        }
+        Ok(Self {
+            time,
+            origin,
+            frame: DynFrame::Icrf,
+            semi_major_axis,
+            eccentricity,
+            inclination,
+            longitude_of_ascending_node,
+            argument_of_periapsis,
+            true_anomaly,
+        })
+    }
+}
+
+impl<T, O, R> Keplerian<T, O, R>
+where
+    T: TimeLike,
+    O: MaybePointMass,
+    R: ReferenceFrame,
+{
     pub fn origin(&self) -> O
     where
         O: Clone,
@@ -78,6 +125,12 @@ where
         self.time.clone()
     }
 
+    pub fn gravitational_parameter(&self) -> f64 {
+        self.origin
+            .maybe_gravitational_parameter()
+            .expect("gravitational parameter should be available")
+    }
+
     pub fn semi_major_axis(&self) -> f64 {
         self.semi_major_axis
     }
@@ -109,15 +162,8 @@ where
             self.semi_major_axis * (1.0 - self.eccentricity.powi(2))
         }
     }
-}
-
-impl<T, O> Keplerian<T, O>
-where
-    T: TimeLike,
-    O: PointMass,
-{
     pub fn to_perifocal(&self) -> (DVec3, DVec3) {
-        let grav_param = self.origin.gravitational_parameter();
+        let grav_param = self.gravitational_parameter();
         let semiparameter = self.semiparameter();
         let (sin_nu, cos_nu) = self.true_anomaly.sin_cos();
         let sqrt_mu_p = (grav_param / semiparameter).sqrt();
@@ -130,29 +176,38 @@ where
     }
 
     pub fn orbital_period(&self) -> TimeDelta {
-        let mu = self.origin.gravitational_parameter();
+        let mu = self.gravitational_parameter();
         let a = self.semi_major_axis();
         TimeDelta::from_decimal_seconds(TAU * (a.powi(3) / mu).sqrt()).unwrap()
     }
 }
 
-impl<T: TimeLike, O: PointMass> CoordinateSystem<Icrf> for Keplerian<T, O> {
-    fn reference_frame(&self) -> Icrf {
-        Icrf
+impl<T: TimeLike, O: MaybePointMass, R: ReferenceFrame + Clone> CoordinateSystem<R>
+    for Keplerian<T, O, R>
+{
+    fn reference_frame(&self) -> R {
+        self.frame.clone()
     }
 }
 
-impl<T, O> ToCartesian<T, O, Icrf> for Keplerian<T, O>
+impl<T, O, R> Keplerian<T, O, R>
 where
     T: TimeLike + Clone,
-    O: PointMass + Clone,
+    O: MaybePointMass + Clone,
+    R: ReferenceFrame + Clone,
 {
-    fn to_cartesian(&self) -> State<T, O, Icrf> {
+    pub(crate) fn to_cartesian(&self) -> State<T, O, R> {
         let (pos, vel) = self.to_perifocal();
         let rot = DMat3::from_rotation_z(self.longitude_of_ascending_node)
             * DMat3::from_rotation_x(self.inclination)
             * DMat3::from_rotation_z(self.argument_of_periapsis);
-        State::new(self.time(), rot * pos, rot * vel, self.origin(), Icrf)
+        State::new(
+            self.time(),
+            rot * pos,
+            rot * vel,
+            self.origin(),
+            self.reference_frame(),
+        )
     }
 }
 
@@ -164,41 +219,6 @@ pub fn is_circular(eccentricity: f64) -> bool {
     float_eq!(eccentricity, 0.0, abs <= 1e-8)
 }
 
-enum OrbitShape {
-    SemiMajorAndEccentricity {
-        semi_major_axis: f64,
-        eccentricity: f64,
-    },
-    Radii {
-        periapsis_radius: f64,
-        apoapsis_radius: f64,
-    },
-    Altitudes {
-        periapsis_altitude: f64,
-        apoapsis_altitude: f64,
-    },
-}
-
-enum Anomaly {
-    True(f64),
-    Eccentric(f64),
-    Mean(f64),
-}
-
-impl Default for Anomaly {
-    fn default() -> Self {
-        Anomaly::True(0.0)
-    }
-}
-
-pub struct KeplerianBuilder {
-    shape: OrbitShape,
-    inclination: Option<f64>,
-    longitude_of_ascending_node: Option<f64>,
-    argument_of_periapisis: Option<f64>,
-    anomaly: Option<Anomaly>,
-}
-
 #[cfg(test)]
 mod tests {
     use super::*;

crates/lox-orbits/src/frames.rs

@@ -6,10 +6,8 @@
  * file, you can obtain one at https://mozilla.org/MPL/2.0/.
  */
 
-use std::f64::consts::FRAC_PI_2;
 use std::{convert::Infallible, str::FromStr};
 
-use crate::ground::GroundLocation;
 use glam::{DMat3, DVec3};
 use lox_bodies::{DynOrigin, MaybeRotationalElements, Origin, RotationalElements, Spheroid};
 use lox_math::types::units::Seconds;
@@ -18,6 +16,9 @@ use thiserror::Error;
 
 use crate::rotations::Rotation;
 
+mod iau;
+mod iers;
+
 pub trait ReferenceFrame {
     fn name(&self) -> String;
     fn abbreviation(&self) -> String;
@@ -149,26 +150,6 @@ impl<T: RotationalElements> ReferenceFrame for BodyFixed<T> {
     }
 }
 
-// #[derive(Clone, Debug)]
-// pub struct Topocentric<B: Spheroid>(GroundLocation<B>);
-
-// impl<B: Spheroid> Topocentric<B> {
-//     pub fn new(location: GroundLocation<B>) -> Self {
-//         Topocentric(location)
-//     }
-//
-//     pub fn from_coords(longitude: f64, latitude: f64, altitude: f64, body: B) -> Self {
-//         let location = GroundLocation::new(longitude, latitude, altitude, body);
-//         Topocentric(location)
-//     }
-//
-//     pub fn rotation_from_body_fixed(&self) -> DMat3 {
-//         let r1 = DMat3::from_rotation_z(self.0.longitude()).transpose();
-//         let r2 = DMat3::from_rotation_y(FRAC_PI_2 - self.0.latitude()).transpose();
-//         r2 * r1
-//     }
-// }
-
 #[derive(Clone, Copy, Debug, PartialEq, Eq, PartialOrd, Ord)]
 pub enum DynFrame {
     Icrf,
@@ -256,34 +237,43 @@ impl FromStr for DynFrame {
     }
 }
 
+pub trait RotateTo<T: ReferenceFrame> {
+    type Time;
+
+    fn rotation(&self, frame: &T, time: Self::Time) -> Rotation;
+}
+
+pub trait TryRotateTo<T: ReferenceFrame> {
+    type Time;
+    type Error;
+
+    fn try_rotation(&self, frame: &T) -> Result<Rotation, Self::Error>;
+}
+
+impl<T: ReferenceFrame, U: RotateTo<T>> TryRotateTo<T> for U {
+    type Error = Infallible;
+
+    fn try_rotation(&self, frame: &T) -> Result<Rotation, Self::Error> {
+        Ok(self.rotation(frame))
+    }
+}
+
+impl TryRotateTo<DynFrame> for DynFrame {
+    type Error = ();
+
+    fn try_rotation(&self, frame: &DynFrame) -> Result<Rotation, Self::Error> {
+        match self {
+            _ => todo!(),
+        }
+    }
+}
+
 #[cfg(test)]
 mod tests {
     use super::*;
-    use lox_bodies::Earth;
-    use lox_math::assert_close;
     use lox_math::is_close::IsClose;
     use rstest::rstest;
 
-    // #[test]
-    // fn test_topocentric() {
-    //     let topo = Topocentric::from_coords(8.0, 50.0, 0.0, Earth);
-    //     let r = topo.rotation_from_body_fixed();
-    //     let x_axis = DVec3::new(
-    //         0.038175550084451906,
-    //         -0.9893582466233818,
-    //         -0.14040258976976597,
-    //     );
-    //     let y_axis = DVec3::new(
-    //         -0.25958272521858694,
-    //         -0.14550003380861354,
-    //         0.9546970980000851,
-    //     );
-    //     let z_axis = DVec3::new(-0.9649660284921128, 0.0, -0.2623748537039304);
-    //     assert_close!(r.x_axis, x_axis);
-    //     assert_close!(r.y_axis, y_axis);
-    //     assert_close!(r.z_axis, z_axis);
-    // }
-
     #[rstest]
     #[case("IAU_EARTH", Some(DynFrame::BodyFixed(DynOrigin::Earth)))]
     #[case("FOO_EARTH", None)]

crates/lox-orbits/src/frames/iers.rs

@@ -0,0 +1,15 @@
+/*
+ * Copyright (c) 2024. 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::frames::{Cirf, Icrf, RotateTo};
+use crate::rotations::Rotation;
+
+impl RotateTo<Cirf> for Icrf {
+    fn rotation(&self, _: &Cirf) -> Rotation {
+        todo!()
+    }
+}