nyx_space/od/ground_station/
mod.rs

1/*
2    Nyx, blazing fast astrodynamics
3    Copyright (C) 2018-onwards Christopher Rabotin <christopher.rabotin@gmail.com>
4
5    This program is free software: you can redistribute it and/or modify
6    it under the terms of the GNU Affero General Public License as published
7    by the Free Software Foundation, either version 3 of the License, or
8    (at your option) any later version.
9
10    This program is distributed in the hope that it will be useful,
11    but WITHOUT ANY WARRANTY; without even the implied warranty of
12    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13    GNU Affero General Public License for more details.
14
15    You should have received a copy of the GNU Affero General Public License
16    along with this program.  If not, see <https://www.gnu.org/licenses/>.
17*/
18
19use anise::astro::{Aberration, AzElRange, PhysicsResult};
20use anise::constants::frames::EARTH_J2000;
21use anise::errors::AlmanacResult;
22use anise::prelude::{Almanac, Frame, Orbit};
23use indexmap::{IndexMap, IndexSet};
24use snafu::ensure;
25
26use super::msr::MeasurementType;
27use super::noise::{GaussMarkov, StochasticNoise};
28use super::{ODAlmanacSnafu, ODError, ODTrajSnafu, TrackingDevice};
29use crate::io::ConfigRepr;
30use crate::od::NoiseNotConfiguredSnafu;
31use crate::time::Epoch;
32use hifitime::Duration;
33use rand_pcg::Pcg64Mcg;
34use serde_derive::{Deserialize, Serialize};
35use std::fmt;
36
37pub mod builtin;
38pub mod event;
39pub mod trk_device;
40
41/// GroundStation defines a two-way ranging and doppler station.
42#[derive(Debug, Clone, Serialize, Deserialize, PartialEq)]
43pub struct GroundStation {
44    pub name: String,
45    /// in degrees
46    pub elevation_mask_deg: f64,
47    /// in degrees
48    pub latitude_deg: f64,
49    /// in degrees
50    pub longitude_deg: f64,
51    /// in km
52    pub height_km: f64,
53    pub frame: Frame,
54    pub measurement_types: IndexSet<MeasurementType>,
55    /// Duration needed to generate a measurement (if unset, it is assumed to be instantaneous)
56    pub integration_time: Option<Duration>,
57    /// Whether to correct for light travel time
58    pub light_time_correction: bool,
59    /// Noise on the timestamp of the measurement
60    pub timestamp_noise_s: Option<StochasticNoise>,
61    pub stochastic_noises: Option<IndexMap<MeasurementType, StochasticNoise>>,
62}
63
64impl GroundStation {
65    /// Initializes a point on the surface of a celestial object.
66    /// This is meant for analysis, not for spacecraft navigation.
67    pub fn from_point(
68        name: String,
69        latitude_deg: f64,
70        longitude_deg: f64,
71        height_km: f64,
72        frame: Frame,
73    ) -> Self {
74        Self {
75            name,
76            elevation_mask_deg: 0.0,
77            latitude_deg,
78            longitude_deg,
79            height_km,
80            frame,
81            measurement_types: IndexSet::new(),
82            integration_time: None,
83            light_time_correction: false,
84            timestamp_noise_s: None,
85            stochastic_noises: None,
86        }
87    }
88
89    /// Returns a copy of this ground station with the new measurement type added (or replaced)
90    pub fn with_msr_type(mut self, msr_type: MeasurementType, noise: StochasticNoise) -> Self {
91        if self.stochastic_noises.is_none() {
92            self.stochastic_noises = Some(IndexMap::new());
93        }
94
95        self.stochastic_noises
96            .as_mut()
97            .unwrap()
98            .insert(msr_type, noise);
99
100        self.measurement_types.insert(msr_type);
101
102        self
103    }
104
105    /// Returns a copy of this ground station without the provided measurement type (if defined, else no error)
106    pub fn without_msr_type(mut self, msr_type: MeasurementType) -> Self {
107        if let Some(noises) = self.stochastic_noises.as_mut() {
108            noises.swap_remove(&msr_type);
109        }
110
111        self.measurement_types.swap_remove(&msr_type);
112
113        self
114    }
115
116    pub fn with_integration_time(mut self, integration_time: Option<Duration>) -> Self {
117        self.integration_time = integration_time;
118
119        self
120    }
121
122    /// Returns a copy of this ground station with the measurement type noises' constant bias set to the provided value.
123    pub fn with_msr_bias_constant(
124        mut self,
125        msr_type: MeasurementType,
126        bias_constant: f64,
127    ) -> Result<Self, ODError> {
128        if self.stochastic_noises.is_none() {
129            self.stochastic_noises = Some(IndexMap::new());
130        }
131
132        let stochastics = self.stochastic_noises.as_mut().unwrap();
133
134        let this_noise = stochastics
135            .get_mut(&msr_type)
136            .ok_or(ODError::NoiseNotConfigured {
137                kind: format!("{msr_type:?}"),
138            })
139            .unwrap();
140
141        if this_noise.bias.is_none() {
142            this_noise.bias = Some(GaussMarkov::ZERO);
143        }
144
145        this_noise.bias.unwrap().constant = Some(bias_constant);
146
147        Ok(self)
148    }
149
150    /// Computes the azimuth and elevation of the provided object seen from this ground station, both in degrees.
151    /// This is a shortcut to almanac.azimuth_elevation_range_sez.
152    pub fn azimuth_elevation_of(
153        &self,
154        rx: Orbit,
155        obstructing_body: Option<Frame>,
156        almanac: &Almanac,
157    ) -> AlmanacResult<AzElRange> {
158        let ab_corr = if self.light_time_correction {
159            Aberration::LT
160        } else {
161            Aberration::NONE
162        };
163        almanac.azimuth_elevation_range_sez(
164            rx,
165            self.to_orbit(rx.epoch, almanac).unwrap(),
166            obstructing_body,
167            ab_corr,
168        )
169    }
170
171    /// Return this ground station as an orbit in its current frame
172    pub fn to_orbit(&self, epoch: Epoch, almanac: &Almanac) -> PhysicsResult<Orbit> {
173        Orbit::try_latlongalt(
174            self.latitude_deg,
175            self.longitude_deg,
176            self.height_km,
177            epoch,
178            almanac.frame_info(self.frame).unwrap(),
179        )
180    }
181
182    /// Returns the noises for all measurement types configured for this ground station at the provided epoch, timestamp noise is the first entry.
183    fn noises(&mut self, epoch: Epoch, rng: Option<&mut Pcg64Mcg>) -> Result<Vec<f64>, ODError> {
184        let mut noises = vec![0.0; self.measurement_types.len() + 1];
185
186        if let Some(rng) = rng {
187            ensure!(
188                self.stochastic_noises.is_some(),
189                NoiseNotConfiguredSnafu {
190                    kind: "ground station stochastics".to_string(),
191                }
192            );
193            // Add the timestamp noise first
194
195            if let Some(mut timestamp_noise) = self.timestamp_noise_s {
196                noises[0] = timestamp_noise.sample(epoch, rng);
197            }
198
199            let stochastics = self.stochastic_noises.as_mut().unwrap();
200
201            for (ii, msr_type) in self.measurement_types.iter().enumerate() {
202                noises[ii + 1] = stochastics
203                    .get_mut(msr_type)
204                    .ok_or(ODError::NoiseNotConfigured {
205                        kind: format!("{msr_type:?}"),
206                    })?
207                    .sample(epoch, rng);
208            }
209        }
210
211        Ok(noises)
212    }
213}
214
215impl Default for GroundStation {
216    fn default() -> Self {
217        let mut measurement_types = IndexSet::new();
218        measurement_types.insert(MeasurementType::Range);
219        measurement_types.insert(MeasurementType::Doppler);
220        Self {
221            name: "UNDEFINED".to_string(),
222            measurement_types,
223            elevation_mask_deg: 0.0,
224            latitude_deg: 0.0,
225            longitude_deg: 0.0,
226            height_km: 0.0,
227            frame: EARTH_J2000,
228            integration_time: None,
229            light_time_correction: false,
230            timestamp_noise_s: None,
231            stochastic_noises: None,
232        }
233    }
234}
235
236impl ConfigRepr for GroundStation {}
237
238impl fmt::Display for GroundStation {
239    // Prints the Keplerian orbital elements with units
240    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
241        write!(
242            f,
243            "{} (lat.: {:.4} deg    long.: {:.4} deg    alt.: {:.3} m) [{}]",
244            self.name,
245            self.latitude_deg,
246            self.longitude_deg,
247            self.height_km * 1e3,
248            self.frame,
249        )
250    }
251}
252
253#[cfg(test)]
254mod gs_ut {
255
256    use anise::constants::frames::IAU_EARTH_FRAME;
257    use indexmap::{IndexMap, IndexSet};
258
259    use crate::io::ConfigRepr;
260    use crate::od::prelude::*;
261
262    #[test]
263    fn test_load_single() {
264        use std::env;
265        use std::path::PathBuf;
266
267        use hifitime::TimeUnits;
268
269        let test_data: PathBuf = [
270            env::var("CARGO_MANIFEST_DIR").unwrap(),
271            "data".to_string(),
272            "03_tests".to_string(),
273            "config".to_string(),
274            "one_ground_station.yaml".to_string(),
275        ]
276        .iter()
277        .collect();
278
279        assert!(test_data.exists(), "Could not find the test data");
280
281        let gs = GroundStation::load(test_data).unwrap();
282
283        dbg!(&gs);
284
285        let mut measurement_types = IndexSet::new();
286        measurement_types.insert(MeasurementType::Range);
287        measurement_types.insert(MeasurementType::Doppler);
288
289        let mut stochastics = IndexMap::new();
290        stochastics.insert(
291            MeasurementType::Range,
292            StochasticNoise {
293                bias: Some(GaussMarkov::new(1.days(), 5e-3).unwrap()),
294                ..Default::default()
295            },
296        );
297        stochastics.insert(
298            MeasurementType::Doppler,
299            StochasticNoise {
300                bias: Some(GaussMarkov::new(1.days(), 5e-5).unwrap()),
301                ..Default::default()
302            },
303        );
304
305        let expected_gs = GroundStation {
306            name: "Demo ground station".to_string(),
307            frame: IAU_EARTH_FRAME,
308            measurement_types,
309            elevation_mask_deg: 5.0,
310            stochastic_noises: Some(stochastics),
311            latitude_deg: 2.3522,
312            longitude_deg: 48.8566,
313            height_km: 0.4,
314            light_time_correction: false,
315            timestamp_noise_s: None,
316            integration_time: Some(60 * Unit::Second),
317        };
318
319        println!("{}", serde_yml::to_string(&expected_gs).unwrap());
320
321        assert_eq!(expected_gs, gs);
322    }
323
324    #[test]
325    fn test_load_many() {
326        use hifitime::TimeUnits;
327        use std::env;
328        use std::path::PathBuf;
329
330        let test_file: PathBuf = [
331            env::var("CARGO_MANIFEST_DIR").unwrap(),
332            "data".to_string(),
333            "03_tests".to_string(),
334            "config".to_string(),
335            "many_ground_stations.yaml".to_string(),
336        ]
337        .iter()
338        .collect();
339
340        let stations = GroundStation::load_many(test_file).unwrap();
341
342        dbg!(&stations);
343
344        let mut measurement_types = IndexSet::new();
345        measurement_types.insert(MeasurementType::Range);
346        measurement_types.insert(MeasurementType::Doppler);
347
348        let mut stochastics = IndexMap::new();
349        stochastics.insert(
350            MeasurementType::Range,
351            StochasticNoise {
352                bias: Some(GaussMarkov::new(1.days(), 5e-3).unwrap()),
353                ..Default::default()
354            },
355        );
356        stochastics.insert(
357            MeasurementType::Doppler,
358            StochasticNoise {
359                bias: Some(GaussMarkov::new(1.days(), 5e-5).unwrap()),
360                ..Default::default()
361            },
362        );
363
364        let expected = vec![
365            GroundStation {
366                name: "Demo ground station".to_string(),
367                frame: IAU_EARTH_FRAME.with_mu_km3_s2(398600.435436096),
368                measurement_types: measurement_types.clone(),
369                elevation_mask_deg: 5.0,
370                stochastic_noises: Some(stochastics.clone()),
371                latitude_deg: 2.3522,
372                longitude_deg: 48.8566,
373                height_km: 0.4,
374                light_time_correction: false,
375                timestamp_noise_s: None,
376                integration_time: None,
377            },
378            GroundStation {
379                name: "Canberra".to_string(),
380                frame: IAU_EARTH_FRAME.with_mu_km3_s2(398600.435436096),
381                measurement_types,
382                elevation_mask_deg: 5.0,
383                stochastic_noises: Some(stochastics),
384                latitude_deg: -35.398333,
385                longitude_deg: 148.981944,
386                height_km: 0.691750,
387                light_time_correction: false,
388                timestamp_noise_s: None,
389                integration_time: None,
390            },
391        ];
392
393        assert_eq!(expected, stations);
394
395        // Serialize back
396        let reser = serde_yml::to_string(&expected).unwrap();
397        dbg!(reser);
398    }
399}