gistools/space/
sat.rs

1use crate::{
2    space::{
3        propagation::{SGP4ErrorOutput, SGP4Output, sgp4, sgp4init},
4        util::{
5            constants::MINUTES_PER_DAY,
6            time::{TimeStamp, days2mdhms, jday},
7        },
8    },
9    util::Date,
10};
11use alloc::{format, string::String, vec, vec::Vec};
12use core::f64::consts::PI;
13use regex::Regex;
14use serde::{Deserialize, Serialize};
15
16/// Classification of TLE
17/// - U: unclassified
18/// - C: classified
19/// - S: secret
20#[derive(Debug, Default, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
21#[repr(u8)]
22pub enum Classification {
23    /// Unclassified
24    #[default]
25    U,
26    /// Classified
27    C,
28    /// Secret
29    S,
30}
31impl From<&str> for Classification {
32    fn from(s: &str) -> Self {
33        match s {
34            "U" | "u" => Classification::U,
35            "C" | "c" => Classification::C,
36            "S" | "s" => Classification::S,
37            _ => Classification::U,
38        }
39    }
40}
41
42/// Mode of operation AFSPC or Improved
43/// - a: afspc
44/// - i: improved
45#[derive(Debug, Default, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
46#[repr(u8)]
47pub enum OperationMode {
48    /// AFSPC
49    A,
50    /// Improved
51    #[default]
52    I,
53}
54/// Method of orbit determination
55/// - d: deep space
56/// - n: near earth
57#[derive(Debug, Default, Clone, Copy, PartialEq, Eq, Serialize, Deserialize)]
58#[repr(u8)]
59pub enum Method {
60    /// Deep Space
61    D,
62    /// Near Earth
63    #[default]
64    N,
65}
66
67/// TLE Data Interface
68#[derive(Debug, Default, Clone, PartialEq, Serialize, Deserialize)]
69pub struct TLEData {
70    /// Name
71    pub name: String,
72    /// Number
73    pub number: f64,
74    /// Classification
75    pub class: Classification,
76    /// Identifier
77    pub id: String,
78    /// Date
79    pub date: Date,
80    /// Epoch Days
81    pub epochdays: f64,
82    /// fd mm
83    pub fdmm: f64,
84    /// sd mm
85    pub sdmm: f64,
86    /// Drag coefficient
87    pub drag: f64,
88    /// Mean motion
89    pub ephemeris: f64,
90    /// Eccentricity
91    pub esn: f64,
92    /// Inclination
93    pub inclination: f64,
94    /// Right ascension
95    pub ascension: f64,
96    /// Eccentricity
97    pub eccentricity: f64,
98    /// Perigee
99    pub perigee: f64,
100    /// Mean anomaly
101    pub anomaly: f64,
102    /// Mean motion
103    pub motion: f64,
104    /// Revolution
105    pub revolution: f64,
106    /// Bstar
107    pub rms: Option<f64>,
108}
109impl From<&str> for TLEData {
110    /// Parse TLE string
111    /// https://en.wikipedia.org/wiki/Two-line_element_set
112    fn from(value: &str) -> Self {
113        let mut lines: Vec<&str> = trim(value).lines().collect();
114        let mut tle = TLEData::default();
115
116        // Line 0: optional name
117        if lines.len() >= 3 {
118            let mut name = trim(lines.remove(0));
119            if name.starts_with("0 ") {
120                name = &name[2..];
121            }
122            tle.name = name.into();
123        }
124
125        // Line 1
126        let line = lines.remove(0);
127        let checksum = check(line);
128        if checksum != line[68..69].parse::<u32>().unwrap() {
129            panic!("Line 1 checksum mismatch: {} != {}: {}", checksum, &line[68..69], line);
130        }
131
132        tle.number = parse_float(&alpha5_converter(&line[2..7]));
133        tle.class = trim(&line[7..9]).into();
134        tle.id = trim(&line[9..18]).into();
135        (tle.date, tle.epochdays) = parse_epoch(&line[18..33]);
136        tle.fdmm = parse_float(&line[33..44]);
137        tle.sdmm = parse_float(&line[44..53]);
138        tle.drag = parse_drag(&line[53..62]);
139        tle.ephemeris = parse_float(&line[62..64]);
140        tle.esn = parse_float(&line[64..68]);
141
142        // Line 2
143        let line = lines.remove(0);
144        let checksum = check(line);
145        if checksum != line[68..69].parse::<u32>().unwrap() {
146            panic!("Line 2 checksum mismatch: {} != {}: {}", checksum, &line[68..69], line);
147        }
148
149        tle.inclination = parse_float(&line[8..17]);
150        tle.ascension = parse_float(&line[17..26]);
151        tle.eccentricity = parse_float(&format!("0.{}", &line[26..34]));
152        tle.perigee = parse_float(&line[34..43]);
153        tle.anomaly = parse_float(&line[43..52]);
154        tle.motion = parse_float(&line[52..63]);
155        tle.revolution = parse_float(&line[63..68]);
156
157        tle
158    }
159}
160
161/// Celestrak TLE Data Interface
162#[derive(Debug, Default, Clone, PartialEq, Serialize, Deserialize)]
163pub struct TLEDataCelestrak {
164    /// Object name
165    #[serde(rename = "OBJECT_NAME")]
166    pub object_name: String,
167    /// Object ID
168    #[serde(rename = "OBJECT_ID")]
169    pub object_id: String,
170    /// Epoch
171    #[serde(rename = "EPOCH")]
172    pub epoch: String,
173    /// Mean Motion
174    #[serde(rename = "MEAN_MOTION")]
175    pub mean_motion: f64,
176    /// Eccentricity
177    #[serde(rename = "ECCENTRICITY")]
178    pub eccentricity: f64,
179    /// Inclination
180    #[serde(rename = "INCLINATION")]
181    pub inclination: f64,
182    /// Right Ascension
183    #[serde(rename = "RA_OF_ASC_NODE")]
184    pub ra_of_asc_node: f64,
185    /// Argument of Peri-center
186    #[serde(rename = "ARG_OF_PERICENTER")]
187    pub arg_of_pericenter: f64,
188    /// Mean Anomaly
189    #[serde(rename = "MEAN_ANOMALY")]
190    pub mean_anomaly: f64,
191    /// Ephemeris Type
192    #[serde(rename = "EPHEMERIS_TYPE")]
193    pub ephemeris_type: f64,
194    /// Classification Type
195    #[serde(rename = "CLASSIFICATION_TYPE")]
196    pub classification_type: String,
197    /// Norad Cat ID
198    #[serde(rename = "NORAD_CAT_ID")]
199    pub norad_cat_id: f64,
200    /// Element Set Number
201    #[serde(rename = "ELEMENT_SET_NO")]
202    pub element_set_no: f64,
203    /// Rev at Epoch
204    #[serde(rename = "REV_AT_EPOCH")]
205    pub rev_at_epoch: f64,
206    /// Bstar
207    #[serde(rename = "BSTAR")]
208    pub bstar: f64,
209    /// Mean Motion Dot
210    #[serde(rename = "MEAN_MOTION_DOT")]
211    pub mean_motion_dot: f64,
212    /// Mean Motion Ddot
213    #[serde(rename = "MEAN_MOTION_DDOT")]
214    pub mean_motion_ddot: f64,
215    /// RMS
216    #[serde(rename = "RMS")]
217    pub rms: String,
218    /// Data Source
219    #[serde(rename = "DATA_SOURCE")]
220    pub data_source: String,
221}
222/// Convert Celestrak TLE data to a standard TLE data object
223/// [JSON example](https://celestrak.org/NORAD/elements/supplemental/index.php?FORMAT=json)
224impl From<&TLEDataCelestrak> for TLEData {
225    fn from(data: &TLEDataCelestrak) -> Self {
226        // convert date to UTC to avoid javascripts localization issues
227        let date: Date = (&*data.epoch).into();
228        let start = Date::new(date.year, 0, 0);
229        TLEData {
230            name: data.object_name.clone(),
231            number: data.norad_cat_id,
232            class: (&*data.classification_type).into(),
233            id: data.object_id.clone(),
234            date,
235            epochdays: jday(&date) - jday(&start),
236            fdmm: data.mean_motion_dot,
237            sdmm: data.mean_motion_ddot,
238            drag: data.bstar,
239            ephemeris: data.ephemeris_type,
240            esn: data.element_set_no,
241            inclination: data.inclination,
242            ascension: data.ra_of_asc_node,
243            eccentricity: data.eccentricity,
244            perigee: data.arg_of_pericenter,
245            anomaly: data.mean_anomaly,
246            motion: data.mean_motion,
247            revolution: data.rev_at_epoch,
248            rms: data.rms.parse().ok(),
249        }
250    }
251}
252
253/// # Satellite Orbit Class
254///
255/// ## Description
256/// A class representing a satellite orbit.
257///
258/// ## Examples
259///
260/// ### Input TLE example
261/// ```txt
262/// STARLINK-1007
263/// 1 44713C 19074A   23048.53451389 -.00009219  00000+0 -61811-3 0   482
264/// 2 44713  53.0512 157.2379 0001140  81.3827  74.7980 15.06382459    15
265/// ```
266///
267/// ### Run example
268/// ```ts
269/// import { Satellite } from 'gis-tools-ts';
270///
271/// const sat = new Satellite(tleString);
272/// // get propagation at time
273/// const { position, velocity } = sat.propagate(new Date());
274/// ```
275///
276/// ## Links
277/// - https://en.wikipedia.org/wiki/Two-line_element_set
278/// - https://celestrak.org/NORAD/documentation/tle-fmt.php
279/// - https://www.space-track.org/documentation#tle-basic
280#[derive(Debug, Default, Clone, PartialEq, Serialize, Deserialize)]
281pub struct Satellite {
282    /// If the satellite is initialized
283    pub init: bool, // = false;
284    // Line 0
285    /// Name of the satellite
286    pub name: String, // = 'default',
287    // Line 1
288    /// (satnum) Satellite catalog number or NORAD (North American Aerospace Defense) Catalog Number
289    pub number: f64,
290    /// Classification (U: unclassified, C: classified, S: secret)
291    pub class: Classification,
292    /// International Designator
293    pub id: String, // = 'null';
294    /// (epochyr + epochdays)
295    pub date: Date, // = new Date();
296    /// Epoch year
297    pub epochyr: f64,
298    /// Epoch days
299    pub epochdays: f64,
300    /// full Sat epoch
301    pub jdsatepoch: f64,
302    /// (ndot) First derivative of mean motion; the ballistic coefficient
303    pub fdmm: f64,
304    /// (nddot) Second derivative of mean motion (decimal point assumed)
305    pub sdmm: f64,
306    /// (bstar) B*, the drag term, or radiation pressure coefficient (decimal point assumed)
307    pub drag: f64,
308    /// Ephemeris type (always zero; only used in undistributed TLE data)
309    pub ephemeris: f64,
310    /// Element set number. Incremented when a new TLE is generated for this object.
311    pub esn: f64,
312    // Line 2
313    /// Inclination (degrees)
314    pub inclination: f64,
315    /// Right ascension of the ascending node (degrees)
316    pub ascension: f64,
317    /// Eccentricity (decimal point assumed)
318    pub eccentricity: f64,
319    /// Argument of perigee (degrees)
320    pub perigee: f64,
321    /// Mean anomaly (degrees)    
322    pub anomaly: f64,
323    /// Mean motion (revolutions per day)
324    pub motion: f64,
325    /// Revolution number at epoch (revolutions)
326    pub revolution: f64,
327    // extra
328    /// Operation Mode
329    pub opsmode: OperationMode,
330    /// RMS
331    pub rms: Option<f64>,
332    // ------------ all near earth variables ------------
333    /// Is IMP
334    pub isimp: f64,
335    /// method
336    pub method: Method,
337    /// ay coefficient
338    pub aycof: f64,
339    /// con 41
340    pub con41: f64,
341    /// cc1
342    pub cc1: f64,
343    /// cc4
344    pub cc4: f64,
345    /// cc5
346    pub cc5: f64,
347    /// d2
348    pub d2: f64,
349    /// d3
350    pub d3: f64,
351    /// d4
352    pub d4: f64,
353    /// delmo
354    pub delmo: f64,
355    /// eta
356    pub eta: f64,
357    /// argpdot
358    pub argpdot: f64,
359    /// omgcof
360    pub omgcof: f64,
361    /// sinmao
362    pub sinmao: f64,
363    /// t2cof
364    pub t2cof: f64,
365    /// t3cof
366    pub t3cof: f64,
367    /// t4cof
368    pub t4cof: f64,
369    /// t5cof
370    pub t5cof: f64,
371    /// x1mth2
372    pub x1mth2: f64,
373    /// x7thm1
374    pub x7thm1: f64,
375    /// mdot
376    pub mdot: f64,
377    /// nodedot
378    pub nodedot: f64,
379    /// xlcof
380    pub xlcof: f64,
381    /// xmcof
382    pub xmcof: f64,
383    /// nodecf
384    pub nodecf: f64,
385    // ------------ all deep space variables ------------
386    /// irez
387    pub irez: f64,
388    /// d2201
389    pub d2201: f64,
390    /// d2211
391    pub d2211: f64,
392    /// d3210
393    pub d3210: f64,
394    /// d3222
395    pub d3222: f64,
396    /// d4410
397    pub d4410: f64,
398    /// d4422
399    pub d4422: f64,
400    /// d5220
401    pub d5220: f64,
402    /// d5232
403    pub d5232: f64,
404    /// d5421
405    pub d5421: f64,
406    /// d5433
407    pub d5433: f64,
408    /// dedt
409    pub dedt: f64,
410    /// del1
411    pub del1: f64,
412    /// del2
413    pub del2: f64,
414    /// del3
415    pub del3: f64,
416    /// didt
417    pub didt: f64,
418    /// dmdt
419    pub dmdt: f64,
420    /// dnodt
421    pub dnodt: f64,
422    /// domdt
423    pub domdt: f64,
424    /// e3
425    pub e3: f64,
426    /// ee2
427    pub ee2: f64,
428    /// peo
429    pub peo: f64,
430    /// pgho
431    pub pgho: f64,
432    /// pho
433    pub pho: f64,
434    /// pinco
435    pub pinco: f64,
436    /// plo
437    pub plo: f64,
438    /// se2
439    pub se2: f64,
440    /// se3
441    pub se3: f64,
442    /// sgh2
443    pub sgh2: f64,
444    /// sgh3
445    pub sgh3: f64,
446    /// sgh4
447    pub sgh4: f64,
448    /// sh2
449    pub sh2: f64,
450    /// sh3
451    pub sh3: f64,
452    /// si2
453    pub si2: f64,
454    /// si3
455    pub si3: f64,
456    /// sl2
457    pub sl2: f64,
458    /// sl3
459    pub sl3: f64,
460    /// sl4
461    pub sl4: f64,
462    /// gsto
463    pub gsto: f64,
464    /// xfact
465    pub xfact: f64,
466    /// xgh2
467    pub xgh2: f64,
468    /// xgh3
469    pub xgh3: f64,
470    /// xgh4
471    pub xgh4: f64,
472    /// xh2
473    pub xh2: f64,
474    /// xh3
475    pub xh3: f64,
476    /// xi2
477    pub xi2: f64,
478    /// xi3
479    pub xi3: f64,
480    /// xl2
481    pub xl2: f64,
482    /// xl3
483    pub xl3: f64,
484    /// xl4
485    pub xl4: f64,
486    /// xlamo
487    pub xlamo: f64,
488    /// zmol
489    pub zmol: f64,
490    /// zmos
491    pub zmos: f64,
492    /// atime
493    pub atime: f64,
494    /// xli
495    pub xli: f64,
496    /// xni
497    pub xni: f64,
498}
499impl Satellite {
500    /**
501     * Constructor
502     * @param data - TLE data or TLE string
503     * @param initialize - initialize the object on creation
504     */
505    pub fn new(data: &TLEData, initialize: Option<bool>) -> Self {
506        let mut this = Self::default();
507        let initialize = initialize.unwrap_or(true);
508        this.rms = data.rms;
509        this.name = data.name.clone();
510        this.number = data.number;
511        this.class = data.class;
512        this.id = data.id.clone();
513        this.date = data.date;
514        this.fdmm = data.fdmm;
515        this.sdmm = data.sdmm;
516        this.drag = data.drag;
517        this.ephemeris = data.ephemeris;
518        this.esn = data.esn;
519        this.inclination = data.inclination.to_radians();
520        this.ascension = data.ascension.to_radians();
521        this.eccentricity = data.eccentricity;
522        this.perigee = data.perigee.to_radians();
523        this.anomaly = data.anomaly.to_radians();
524        // convert revolution from deg/day to rad/minute
525        this.motion = data.motion;
526        this.revolution = data.revolution;
527        this.epochdays = data.epochdays;
528
529        this.epochyr = (this.date.year % 100) as f64;
530        // convert revolution from deg/day to rad/minute
531        this.motion /= 1440. / (2. * PI); // rad/min (229.1831180523293)
532        // find sgp4epoch time of element set
533        // remember that sgp4 uses units of days from 0 jan 1950 (sgp4epoch)
534        // and minutes from the epoch (time)
535        let year = if this.epochyr < 57. { this.epochyr + 2000. } else { this.epochyr + 1900. };
536        let mdhms_result = days2mdhms(year as u16, this.epochdays);
537
538        let TimeStamp { mon, day, hr, min, sec } = mdhms_result;
539        this.jdsatepoch = jday(&Date::new_full(
540            year as u16,
541            mon as u8,
542            day as u8,
543            hr as u8,
544            min as u8,
545            sec as u8,
546        ));
547
548        if initialize {
549            sgp4init(&mut this);
550        }
551
552        this
553    }
554
555    /// Converts satellite state to an array that is readable by the GPU
556    ///
557    /// ## Returns
558    /// Satellite state in an array
559    pub fn gpu(&self) -> Vec<f64> {
560        vec![
561            self.anomaly,
562            self.motion,
563            self.eccentricity,
564            self.inclination,
565            if self.method == Method::D { 0. } else { 1. }, // 0 -> 'd', 1 -> 'n'
566            if self.opsmode == OperationMode::A { 0. } else { 1. }, // 0 -> 'a'; 1 -> 'i'
567            self.drag,
568            self.mdot,
569            self.perigee,
570            self.argpdot,
571            self.ascension,
572            self.nodedot,
573            self.nodecf,
574            self.cc1,
575            self.cc4,
576            self.cc5,
577            self.t2cof,
578            self.isimp,
579            self.omgcof,
580            self.eta,
581            self.xmcof,
582            self.delmo,
583            self.d2,
584            self.d3,
585            self.d4,
586            self.sinmao,
587            self.t3cof,
588            self.t4cof,
589            self.t5cof,
590            self.irez,
591            self.d2201,
592            self.d2211,
593            self.d3210,
594            self.d3222,
595            self.d4410,
596            self.d4422,
597            self.d5220,
598            self.d5232,
599            self.d5421,
600            self.d5433,
601            self.dedt,
602            self.del1,
603            self.del2,
604            self.del3,
605            self.didt,
606            self.dmdt,
607            self.dnodt,
608            self.domdt,
609            self.gsto,
610            self.xfact,
611            self.xlamo,
612            self.atime,
613            self.xli,
614            self.xni,
615            self.aycof,
616            self.xlcof,
617            self.con41,
618            self.x1mth2,
619            self.x7thm1,
620            self.zmos,
621            self.zmol,
622            self.se2,
623            self.se3,
624            self.si2,
625            self.si3,
626            self.sl2,
627            self.sl3,
628            self.sl4,
629            self.sgh2,
630            self.sgh3,
631            self.sgh4,
632            self.sh2,
633            self.sh3,
634            self.ee2,
635            self.e3,
636            self.xi2,
637            self.xi3,
638            self.xl2,
639            self.xl3,
640            self.xl4,
641            self.xgh2,
642            self.xgh3,
643            self.xgh4,
644            self.xh2,
645            self.xh3,
646            self.peo,
647            self.pinco,
648            self.plo,
649            self.pgho,
650            self.pho,
651        ]
652    }
653
654    /// propagate the satellite's position and velocity given a Date input
655    ///
656    /// ## Parameters
657    /// - `time`: Date object
658    ///
659    /// ## Returns
660    /// satellite state at that time given
661    pub fn propagate(&self, time: &Date) -> Result<SGP4Output, SGP4ErrorOutput> {
662        let j = jday(time);
663        sgp4(self, (j - self.jdsatepoch) * MINUTES_PER_DAY)
664    }
665
666    /// time in minutes since epoch
667    ///
668    /// ## Parameters
669    /// - `time`: time in minutes from satellite epoch
670    ///
671    /// ## Returns
672    /// satellite state at that time given
673    pub fn sgp4(&self, time: f64) -> Result<SGP4Output, SGP4ErrorOutput> {
674        sgp4(self, time)
675    }
676}
677
678/// Parse a float number from a string (TLE-style scientific notation)
679fn parse_float(value: &str) -> f64 {
680    let re = Regex::new(r"([-])?([.\d]+)([+-]\d+)?").unwrap();
681
682    if let Some(caps) = re.captures(value) {
683        let sign = if caps.get(1).map_or("", |m| m.as_str()) == "-" { -1.0 } else { 1.0 };
684        let base = caps.get(2).map_or("0", |m| m.as_str());
685        let power = caps.get(3).map_or_else(|| "e0".into(), |m| format!("e{}", m.as_str()));
686        let combined = format!("{}{}", base, power);
687        return sign * combined.parse::<f64>().unwrap();
688    }
689
690    0.0
691}
692
693/// Parse a drag coefficient from TLE-style string
694fn parse_drag(value: &str) -> f64 {
695    let re = Regex::new(r"([-])?([.\d]+)([+-]\d+)?").unwrap();
696    if let Some(caps) = re.captures(value) {
697        let sign = if caps.get(1).map_or("", |m| m.as_str()) == "-" { -1.0 } else { 1.0 };
698        let base = caps.get(2).map_or("0", |m| m.as_str());
699        let base = if base.contains('.') { base.into() } else { format!("0.{}", base) };
700        let power = caps.get(3).map_or_else(|| "e0".into(), |m| format!("e{}", m.as_str()));
701        return sign * format!("{}{}", base, power).parse::<f64>().unwrap();
702    }
703
704    0.0
705}
706
707/// Parse TLE epoch string (e.g., "22345.6789") into your Date struct
708fn parse_epoch(value: &str) -> (Date, f64) {
709    // Trim whitespace
710    let re = Regex::new(r"^\s+|\s+$").unwrap();
711    let value: String = re.replace_all(value, "").into();
712
713    // Parse epoch year and day-of-year
714    let epoch = value[0..2].parse::<u16>().unwrap();
715    let days = value[2..].parse::<f64>().unwrap();
716
717    // Compute full year
718    let now_year = 2025;
719    let current_epoch = (now_year % 100) as u16;
720    let century = now_year - current_epoch as i32;
721    let year = if epoch > current_epoch + 1 {
722        (century - 100 + epoch as i32) as u16
723    } else {
724        (century + epoch as i32) as u16
725    };
726
727    // Convert day-of-year to month/day/hour/minute/second
728    let ts = days2mdhms(year, days);
729
730    // Return Date
731    (
732        Date::new_full(year, ts.mon as u8, ts.day as u8, ts.hr as u8, ts.min as u8, ts.sec as u8),
733        days,
734    )
735}
736
737/// Check TLE checksum
738fn check(line: &str) -> u32 {
739    let mut sum = 0;
740
741    for c in line.chars().take(68) {
742        if c.is_ascii_digit() {
743            sum += c.to_digit(10).unwrap();
744        } else if c == '-' {
745            sum += 1;
746        }
747    }
748
749    sum % 10
750}
751
752/// Converts alpha5 to alpha2
753/// NOTE: Alpha5 skips I and O
754fn alpha5_converter(s: &str) -> String {
755    if let Some(first_char) = s.chars().next() {
756        // if first char is numeric, return unchanged
757        if first_char.is_ascii_digit() {
758            return s.into();
759        }
760
761        let alphabet = "ABCDEFGHJKLMNPQRSTUVWXYZ";
762        if let Some(idx) = alphabet.find(first_char) {
763            // prepend converted index+10, keep rest
764            let rest: String = s.chars().skip(1).collect();
765            return format!("{}{}", idx + 10, rest);
766        }
767    }
768
769    // fallback: return unchanged
770    s.into()
771}
772
773/// Trim leading and trailing whitespace, BOM, and non-breaking spaces
774fn trim(s: &str) -> &str {
775    fn is_trim_char(c: char) -> bool {
776        c.is_whitespace() || c == '\u{FEFF}' || c == '\u{00A0}'
777    }
778    s.trim_matches(is_trim_char)
779}