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, ObservationReliability,
25 RangeReliabilityRow, ReliabilityOptions, ReliabilityReport, ReliabilitySummary,
26};
27
28use crate::astro::frames::transforms::geodetic_from_ecef_proj;
29use crate::dop::{ecef_to_enu_rotation, LineOfSight};
30use crate::frame::Wgs84Geodetic;
31use crate::id::{GnssSatelliteId, GnssSystem};
32use crate::spp::{EphemerisSource, ReceiverSolution};
33
34#[derive(Debug, Clone, Copy, PartialEq)]
36pub struct AraimRow {
37 pub id: GnssSatelliteId,
39 pub line_of_sight: LineOfSight,
41 pub system: GnssSystem,
43 pub elevation_rad: f64,
45}
46
47#[derive(Debug, Clone, PartialEq)]
49pub struct AraimGeometry {
50 pub rows: Vec<AraimRow>,
52 pub receiver: Wgs84Geodetic,
54 pub clock_systems: Vec<GnssSystem>,
56}
57
58impl AraimGeometry {
59 pub fn from_receiver_solution(
64 solution: &ReceiverSolution,
65 eph: &dyn EphemerisSource,
66 t_j2000_s: f64,
67 ) -> Result<Self, AraimError> {
68 if !t_j2000_s.is_finite() {
69 return Err(AraimError::InsufficientGeometry);
70 }
71 let receiver = match solution.geodetic {
72 Some(receiver) => receiver,
73 None => geodetic_from_position(solution.position.as_array())?,
74 };
75 let clock_systems = receiver_solution_clock_systems(solution)?;
76 if solution.used_sats.len() < 3 + clock_systems.len() {
77 return Err(AraimError::InsufficientGeometry);
78 }
79
80 let rx_ecef_m = solution.position.as_array();
81 let enu = ecef_to_enu_rotation(receiver.lat_rad, receiver.lon_rad);
82 let mut rows = Vec::with_capacity(solution.used_sats.len());
83 for &id in &solution.used_sats {
84 let (sat_ecef_m, _) = eph
85 .position_clock_at_j2000_s(id, t_j2000_s)
86 .ok_or(AraimError::InsufficientGeometry)?;
87 let dx = sat_ecef_m[0] - rx_ecef_m[0];
88 let dy = sat_ecef_m[1] - rx_ecef_m[1];
89 let dz = sat_ecef_m[2] - rx_ecef_m[2];
90 let range_m = (dx * dx + dy * dy + dz * dz).sqrt();
91 if !range_m.is_finite() || range_m <= 0.0 {
92 return Err(AraimError::InsufficientGeometry);
93 }
94
95 let line_of_sight = LineOfSight::new(dx / range_m, dy / range_m, dz / range_m);
96 let up = enu[2][0] * line_of_sight.e_x
97 + enu[2][1] * line_of_sight.e_y
98 + enu[2][2] * line_of_sight.e_z;
99 let elevation_rad = up.clamp(-1.0, 1.0).asin();
100 if !elevation_rad.is_finite() {
101 return Err(AraimError::InsufficientGeometry);
102 }
103
104 rows.push(AraimRow {
105 id,
106 line_of_sight,
107 system: id.system,
108 elevation_rad,
109 });
110 }
111
112 Ok(Self {
113 rows,
114 receiver,
115 clock_systems,
116 })
117 }
118}
119
120fn geodetic_from_position(position_m: [f64; 3]) -> Result<Wgs84Geodetic, AraimError> {
121 let [lon_deg, lat_deg, height_m] =
122 geodetic_from_ecef_proj(position_m[0], position_m[1], position_m[2])
123 .map_err(|_| AraimError::InsufficientGeometry)?;
124 Wgs84Geodetic::new(lat_deg.to_radians(), lon_deg.to_radians(), height_m)
125 .map_err(|_| AraimError::InsufficientGeometry)
126}
127
128fn receiver_solution_clock_systems(
129 solution: &ReceiverSolution,
130) -> Result<Vec<GnssSystem>, AraimError> {
131 let clock_systems = if !solution.system_clocks_s.is_empty() {
132 solution
133 .system_clocks_s
134 .iter()
135 .map(|&(system, _)| system)
136 .collect()
137 } else if !solution.metadata.systems.is_empty() {
138 solution.metadata.systems.clone()
139 } else {
140 crate::spp::clock_systems(&solution.used_sats)
141 };
142 if clock_systems.is_empty() {
143 return Err(AraimError::InsufficientGeometry);
144 }
145 for (idx, system) in clock_systems.iter().enumerate() {
146 if clock_systems[..idx].contains(system) {
147 return Err(AraimError::InsufficientGeometry);
148 }
149 }
150 Ok(clock_systems)
151}
152
153#[derive(Debug, Clone, Copy, PartialEq)]
155pub struct IntegrityAllocation {
156 pub phmi_total: f64,
158 pub phmi_vert: f64,
160 pub phmi_hor: f64,
162 pub pfa_vert: f64,
164 pub pfa_hor: f64,
166 pub p_threshold_unmonitored: f64,
168 pub p_emt: f64,
170 pub max_fault_order: usize,
172}
173
174impl IntegrityAllocation {
175 pub const fn lpv_200() -> Self {
177 Self {
178 phmi_total: 1.0e-7,
179 phmi_vert: 9.8e-8,
180 phmi_hor: 2.0e-9,
181 pfa_vert: 3.9e-6,
182 pfa_hor: 9.0e-8,
183 p_threshold_unmonitored: 8.0e-8,
185 p_emt: 1.0e-5,
187 max_fault_order: 2,
188 }
189 }
190}
191
192#[derive(Debug, Clone, Copy, PartialEq, Eq, thiserror::Error)]
194pub enum AraimError {
195 #[error("insufficient ARAIM geometry")]
197 InsufficientGeometry,
198 #[error("unmonitorable ARAIM fault probability exceeds allocation")]
200 UnmonitorableFaultMass,
201 #[error("ARAIM numerical failure")]
203 NumericalFailure,
204 #[error("invalid ARAIM ISM")]
206 InvalidIsm,
207 #[error("invalid ARAIM allocation")]
209 InvalidAllocation,
210}
211
212pub(crate) fn clock_system_for_row(system: GnssSystem) -> GnssSystem {
213 match system {
214 GnssSystem::Sbas => GnssSystem::Gps,
215 other => other,
216 }
217}
218
219pub(crate) fn validate_probability(value: f64, allow_zero: bool) -> bool {
220 value.is_finite()
221 && if allow_zero {
222 (0.0..1.0).contains(&value) || value == 0.0
223 } else {
224 (0.0..1.0).contains(&value)
225 }
226}
227
228pub(crate) fn validate_nonneg_finite(value: f64) -> bool {
229 value.is_finite() && value >= 0.0
230}
231
232pub(crate) fn validate_positive_finite(value: f64) -> bool {
233 value.is_finite() && value > 0.0
234}