Skip to main content

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(crate) use interpolatable::INTERPOLATION_SAMPLES;
28pub use interpolatable::Interpolatable;
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(
49        "Probable bug: Requested epoch {req_epoch}, corresponding to an offset of {req_dur} in a spline of duration {spline_dur}"
50    ))]
51    OutOfSpline {
52        req_epoch: Epoch,
53        req_dur: Duration,
54        spline_dur: Duration,
55    },
56    #[snafu(display("Interpolation failed: {source}"))]
57    Interpolation { source: InterpolationError },
58}
59
60/// Smooth the RIC differences using an in-line median filter.
61/// This avoids allocations and operates directly on the Cartesian components.
62fn smooth_state_diff_in_place(ric_diff: &mut [CartesianState], window_size: usize) {
63    assert!(
64        window_size % 2 == 1,
65        "Window size must be odd for proper median calculation"
66    );
67    let half_window = window_size / 2;
68
69    // Temporary buffer to store sorted values for median calculation
70    let mut temp_buffer = vec![0.0; window_size];
71
72    // Iterate over each state in the array
73    for i in 0..ric_diff.len() {
74        let start = i.saturating_sub(half_window);
75        let end = (i + half_window + 1).min(ric_diff.len());
76
77        // Smooth each component independently
78        for component in 0..6 {
79            // Fill the temporary buffer with values from the current window
80            for (j, idx) in (start..end).enumerate() {
81                temp_buffer[j] = match component {
82                    0 => ric_diff[idx].radius_km.x,
83                    1 => ric_diff[idx].radius_km.y,
84                    2 => ric_diff[idx].radius_km.z,
85                    3 => ric_diff[idx].velocity_km_s.x,
86                    4 => ric_diff[idx].velocity_km_s.y,
87                    5 => ric_diff[idx].velocity_km_s.z,
88                    _ => unreachable!(),
89                };
90            }
91
92            // Sort the buffer to find the median
93            temp_buffer[..end - start].sort_unstable_by(|a, b| a.partial_cmp(b).unwrap());
94
95            // Replace the current value with the median
96            let median = temp_buffer[(end - start) / 2];
97            match component {
98                0 => ric_diff[i].radius_km.x = median,
99                1 => ric_diff[i].radius_km.y = median,
100                2 => ric_diff[i].radius_km.z = median,
101                3 => ric_diff[i].velocity_km_s.x = median,
102                4 => ric_diff[i].velocity_km_s.y = median,
103                5 => ric_diff[i].velocity_km_s.z = median,
104                _ => unreachable!(),
105            }
106        }
107    }
108}