1pub mod fault_modes;
11pub mod ism;
12mod mhss;
13pub mod protection;
14pub mod reliability;
15
16#[cfg(test)]
17mod tests;
18
19pub use fault_modes::{enumerate_fault_modes, FaultHypothesis};
20pub use ism::{ConstellationIsm, Ism, SatelliteIsm, SatelliteIsmModel};
21pub use mhss::{araim, AraimResult, FaultMode};
22pub use protection::ProtectionModel;
23pub use reliability::{
24 reliability_araim, reliability_design, wtest_noncentrality, wtest_noncentrality_components,
25 ObservationReliability, RangeReliabilityRow, ReliabilityOptions, ReliabilityReport,
26 ReliabilitySummary, WtestNoncentralityComponents,
27};
28
29use crate::astro::frames::transforms::geodetic_from_ecef_proj;
30use crate::dop::{ecef_to_enu_rotation, LineOfSight};
31use crate::frame::Wgs84Geodetic;
32use crate::id::{GnssSatelliteId, GnssSystem};
33use crate::spp::{EphemerisSource, ReceiverSolution};
34
35#[derive(Debug, Clone, Copy, PartialEq)]
37pub struct AraimRow {
38 pub id: GnssSatelliteId,
40 pub line_of_sight: LineOfSight,
42 pub system: GnssSystem,
44 pub elevation_rad: f64,
46}
47
48#[derive(Debug, Clone, PartialEq)]
50pub struct AraimGeometry {
51 pub rows: Vec<AraimRow>,
53 pub receiver: Wgs84Geodetic,
55 pub clock_systems: Vec<GnssSystem>,
57}
58
59impl AraimGeometry {
60 pub fn from_receiver_solution(
65 solution: &ReceiverSolution,
66 eph: &dyn EphemerisSource,
67 t_j2000_s: f64,
68 ) -> Result<Self, AraimError> {
69 if !t_j2000_s.is_finite() {
70 return Err(AraimError::InsufficientGeometry);
71 }
72 let receiver = match solution.geodetic {
73 Some(receiver) => receiver,
74 None => geodetic_from_position(solution.position.as_array())?,
75 };
76 let clock_systems = receiver_solution_clock_systems(solution)?;
77 if solution.used_sats.len() < 3 + clock_systems.len() {
78 return Err(AraimError::InsufficientGeometry);
79 }
80
81 let rx_ecef_m = solution.position.as_array();
82 let enu = ecef_to_enu_rotation(receiver.lat_rad, receiver.lon_rad);
83 let mut rows = Vec::with_capacity(solution.used_sats.len());
84 for &id in &solution.used_sats {
85 let (sat_ecef_m, _) = eph
86 .position_clock_at_j2000_s(id, t_j2000_s)
87 .ok_or(AraimError::InsufficientGeometry)?;
88 let dx = sat_ecef_m[0] - rx_ecef_m[0];
89 let dy = sat_ecef_m[1] - rx_ecef_m[1];
90 let dz = sat_ecef_m[2] - rx_ecef_m[2];
91 let range_m = (dx * dx + dy * dy + dz * dz).sqrt();
92 if !range_m.is_finite() || range_m <= 0.0 {
93 return Err(AraimError::InsufficientGeometry);
94 }
95
96 let line_of_sight = LineOfSight::new(dx / range_m, dy / range_m, dz / range_m);
97 let up = enu[2][0] * line_of_sight.e_x
98 + enu[2][1] * line_of_sight.e_y
99 + enu[2][2] * line_of_sight.e_z;
100 let elevation_rad = up.clamp(-1.0, 1.0).asin();
101 if !elevation_rad.is_finite() {
102 return Err(AraimError::InsufficientGeometry);
103 }
104
105 rows.push(AraimRow {
106 id,
107 line_of_sight,
108 system: id.system,
109 elevation_rad,
110 });
111 }
112
113 Ok(Self {
114 rows,
115 receiver,
116 clock_systems,
117 })
118 }
119}
120
121fn geodetic_from_position(position_m: [f64; 3]) -> Result<Wgs84Geodetic, AraimError> {
122 let [lon_deg, lat_deg, height_m] =
123 geodetic_from_ecef_proj(position_m[0], position_m[1], position_m[2])
124 .map_err(|_| AraimError::InsufficientGeometry)?;
125 Wgs84Geodetic::new(lat_deg.to_radians(), lon_deg.to_radians(), height_m)
126 .map_err(|_| AraimError::InsufficientGeometry)
127}
128
129fn receiver_solution_clock_systems(
130 solution: &ReceiverSolution,
131) -> Result<Vec<GnssSystem>, AraimError> {
132 let clock_systems = if !solution.system_clocks_s.is_empty() {
133 solution
134 .system_clocks_s
135 .iter()
136 .map(|&(system, _)| system)
137 .collect()
138 } else if !solution.metadata.systems.is_empty() {
139 solution.metadata.systems.clone()
140 } else {
141 crate::spp::clock_systems(&solution.used_sats)
142 };
143 if clock_systems.is_empty() {
144 return Err(AraimError::InsufficientGeometry);
145 }
146 for (idx, system) in clock_systems.iter().enumerate() {
147 if clock_systems[..idx].contains(system) {
148 return Err(AraimError::InsufficientGeometry);
149 }
150 }
151 Ok(clock_systems)
152}
153
154#[derive(Debug, Clone, Copy, PartialEq)]
156pub struct IntegrityAllocation {
157 pub phmi_total: f64,
159 pub phmi_vert: f64,
161 pub phmi_hor: f64,
163 pub pfa_vert: f64,
165 pub pfa_hor: f64,
167 pub p_threshold_unmonitored: f64,
169 pub p_emt: f64,
171 pub max_fault_order: usize,
173}
174
175impl IntegrityAllocation {
176 pub const fn lpv_200() -> Self {
178 Self {
179 phmi_total: 1.0e-7,
180 phmi_vert: 9.8e-8,
181 phmi_hor: 2.0e-9,
182 pfa_vert: 3.9e-6,
183 pfa_hor: 9.0e-8,
184 p_threshold_unmonitored: 8.0e-8,
186 p_emt: 1.0e-5,
188 max_fault_order: 2,
189 }
190 }
191}
192
193#[derive(Debug, Clone, Copy, PartialEq, Eq, thiserror::Error)]
195pub enum AraimError {
196 #[error("insufficient ARAIM geometry")]
198 InsufficientGeometry,
199 #[error("unmonitorable ARAIM fault probability exceeds allocation")]
201 UnmonitorableFaultMass,
202 #[error("ARAIM numerical failure")]
204 NumericalFailure,
205 #[error("invalid ARAIM ISM")]
207 InvalidIsm,
208 #[error("invalid ARAIM allocation")]
210 InvalidAllocation,
211}
212
213pub(crate) fn clock_system_for_row(system: GnssSystem) -> GnssSystem {
214 match system {
215 GnssSystem::Sbas => GnssSystem::Gps,
216 other => other,
217 }
218}
219
220pub(crate) fn validate_probability(value: f64, allow_zero: bool) -> bool {
221 value.is_finite()
222 && if allow_zero {
223 (0.0..1.0).contains(&value) || value == 0.0
224 } else {
225 (0.0..1.0).contains(&value)
226 }
227}
228
229pub(crate) fn validate_nonneg_finite(value: f64) -> bool {
230 value.is_finite() && value >= 0.0
231}
232
233pub(crate) fn validate_positive_finite(value: f64) -> bool {
234 value.is_finite() && value > 0.0
235}