nyx_space/md/trajectory/
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::math::{cartesian::CartesianState, interpolation::InterpolationError};
20use snafu::prelude::*;
21
22mod interpolatable;
23mod sc_traj;
24mod traj;
25mod traj_it;
26
27pub use interpolatable::Interpolatable;
28pub(crate) use interpolatable::INTERPOLATION_SAMPLES;
29pub use traj::Traj;
30
31pub use crate::io::ExportCfg;
32
33use super::StateParameter;
34use crate::time::{Duration, Epoch};
35
36#[derive(Clone, PartialEq, Debug, Snafu)]
37pub enum TrajError {
38    #[snafu(display("Event {event} not found between {start} and {end}"))]
39    EventNotFound {
40        start: Epoch,
41        end: Epoch,
42        event: String,
43    },
44    #[snafu(display("No interpolation data at {epoch}"))]
45    NoInterpolationData { epoch: Epoch },
46    #[snafu(display("Failed to create trajectory: {msg}"))]
47    CreationError { msg: String },
48    #[snafu(display("Probable bug: Requested epoch {req_epoch}, corresponding to an offset of {req_dur} in a spline of duration {spline_dur}"))]
49    OutOfSpline {
50        req_epoch: Epoch,
51        req_dur: Duration,
52        spline_dur: Duration,
53    },
54    #[snafu(display("Interpolation failed: {source}"))]
55    Interpolation { source: InterpolationError },
56}
57
58/// Smooth the RIC differences using an in-line median filter.
59/// This avoids allocations and operates directly on the Cartesian components.
60fn smooth_state_diff_in_place(ric_diff: &mut [CartesianState], window_size: usize) {
61    assert!(
62        window_size % 2 == 1,
63        "Window size must be odd for proper median calculation"
64    );
65    let half_window = window_size / 2;
66
67    // Temporary buffer to store sorted values for median calculation
68    let mut temp_buffer = vec![0.0; window_size];
69
70    // Iterate over each state in the array
71    for i in 0..ric_diff.len() {
72        let start = i.saturating_sub(half_window);
73        let end = (i + half_window + 1).min(ric_diff.len());
74
75        // Smooth each component independently
76        for component in 0..6 {
77            // Fill the temporary buffer with values from the current window
78            for (j, idx) in (start..end).enumerate() {
79                temp_buffer[j] = match component {
80                    0 => ric_diff[idx].radius_km.x,
81                    1 => ric_diff[idx].radius_km.y,
82                    2 => ric_diff[idx].radius_km.z,
83                    3 => ric_diff[idx].velocity_km_s.x,
84                    4 => ric_diff[idx].velocity_km_s.y,
85                    5 => ric_diff[idx].velocity_km_s.z,
86                    _ => unreachable!(),
87                };
88            }
89
90            // Sort the buffer to find the median
91            temp_buffer[..end - start].sort_unstable_by(|a, b| a.partial_cmp(b).unwrap());
92
93            // Replace the current value with the median
94            let median = temp_buffer[(end - start) / 2];
95            match component {
96                0 => ric_diff[i].radius_km.x = median,
97                1 => ric_diff[i].radius_km.y = median,
98                2 => ric_diff[i].radius_km.z = median,
99                3 => ric_diff[i].velocity_km_s.x = median,
100                4 => ric_diff[i].velocity_km_s.y = median,
101                5 => ric_diff[i].velocity_km_s.z = median,
102                _ => unreachable!(),
103            }
104        }
105    }
106}