lox_orbits/propagators/
sgp4.rs1use glam::DVec3;
6pub use sgp4::Elements;
7use sgp4::{Constants, ElementsError, MinutesSinceEpoch};
8use thiserror::Error;
9
10use lox_bodies::Earth;
11use lox_core::f64::consts::SECONDS_PER_MINUTE;
12use lox_time::Time;
13use lox_time::deltas::TimeDelta;
14use lox_time::time_scales::Tai;
15use lox_time::utc::{Utc, UtcError};
16
17use crate::propagators::Propagator;
18use crate::states::State;
19use crate::trajectories::TrajectoryError;
20use lox_frames::Icrf;
21
22#[derive(Debug, Clone, Error)]
23pub enum Sgp4Error {
24 #[error(transparent)]
25 ElementsError(#[from] ElementsError),
26 #[error(transparent)]
27 TrajectoryError(#[from] TrajectoryError),
28 #[error(transparent)]
29 Sgp4(#[from] sgp4::Error),
30 #[error(transparent)]
31 Utc(#[from] UtcError),
32}
33
34pub struct Sgp4 {
35 constants: Constants,
36 time: Time<Tai>,
37}
38
39impl Sgp4 {
40 pub fn new(initial_state: Elements) -> Result<Self, Sgp4Error> {
41 let epoch = initial_state.epoch();
42 let time = Utc::from_delta(TimeDelta::from_julian_years(epoch))?.to_time();
43 let constants = Constants::from_elements(&initial_state)?;
44 Ok(Self { constants, time })
45 }
46
47 pub fn time(&self) -> Time<Tai> {
48 self.time
49 }
50}
51
52impl Propagator<Tai, Earth, Icrf> for Sgp4 {
53 type Error = Sgp4Error;
54
55 fn propagate(&self, time: Time<Tai>) -> Result<State<Tai, Earth, Icrf>, Self::Error> {
56 let dt = (time - self.time).as_seconds_f64() / SECONDS_PER_MINUTE;
57 let prediction = self.constants.propagate(MinutesSinceEpoch(dt))?;
58 let position = DVec3::from_array(prediction.position);
59 let velocity = DVec3::from_array(prediction.velocity);
60 Ok(State::new(time, position, velocity, Earth, Icrf))
61 }
62}
63
64#[cfg(test)]
65mod tests {
66 use lox_test_utils::assert_approx_eq;
67
68 use super::*;
69
70 #[test]
71 fn test_sgp4() {
72 let tle = Elements::from_tle(
73 Some("ISS (ZARYA)".to_string()),
74 "1 25544U 98067A 24170.37528350 .00016566 00000+0 30244-3 0 9996".as_bytes(),
75 "2 25544 51.6410 309.3890 0010444 339.5369 107.8830 15.49495945458731".as_bytes(),
76 )
77 .unwrap();
78 let sgp4 = Sgp4::new(tle).unwrap();
79 let orbital_period = 92.821;
80 let t1 = sgp4.time() + TimeDelta::from_minutes(orbital_period);
81 let s1 = sgp4.propagate(t1).unwrap();
82 let k1 = s1.to_keplerian();
83 assert_approx_eq!(
84 k1.orbital_period().as_seconds_f64() / SECONDS_PER_MINUTE,
85 orbital_period,
86 rtol <= 1e-4
87 );
88 }
89}