Skip to main content

thrust/data/eurocontrol/
database.rs

1//! EUROCONTROL airway database handling.
2//!
3//! This module provides functionality to load and query an airway database
4
5use std::hash::Hash;
6use std::{collections::HashMap, path};
7
8use geodesy::prelude::*;
9use once_cell::sync::Lazy;
10use serde::Serialize;
11
12use crate::data::field15::{Connector, Field15Element, Modifier, Point};
13use crate::data::{
14    eurocontrol::aixm::{
15        airport_heliport::{parse_airport_heliport_zip_file, AirportHeliport},
16        arrival_leg::{parse_arrival_leg_zip_file, ArrivalLeg},
17        departure_leg::{parse_departure_leg_zip_file, DepartureLeg},
18        designated_point::{parse_designated_point_zip_file, DesignatedPoint},
19        navaid::{parse_navaid_zip_file, Navaid},
20        route::{parse_route_zip_file, Route},
21        route_segment::{parse_route_segment_zip_file, PointReference, RouteSegment},
22        standard_instrument_arrival::{parse_standard_instrument_arrival_zip_file, StandardInstrumentArrival},
23        standard_instrument_departure::{parse_standard_instrument_departure_zip_file, StandardInstrumentDeparture},
24    },
25    field15::{Altitude, Speed},
26};
27
28/**
29 * An airway database containing navaids, designated points, route segments, and routes.
30 */
31pub struct AirwayDatabase {
32    airports: HashMap<String, AirportHeliport>,
33    navaids: HashMap<String, Navaid>,
34    designated_points: HashMap<String, DesignatedPoint>,
35    route_segments: HashMap<String, RouteSegment>,
36    routes: HashMap<String, Route>,
37    arrival_legs: HashMap<String, ArrivalLeg>,
38    departure_legs: HashMap<String, DepartureLeg>,
39    standard_instrument_arrivals: HashMap<String, StandardInstrumentArrival>,
40    standard_instrument_departures: HashMap<String, StandardInstrumentDeparture>,
41}
42
43impl AirwayDatabase {
44    /// Load the airway database from the specified directory path.
45    pub fn new(path: &path::Path) -> Result<Self, Box<dyn std::error::Error>> {
46        Ok(AirwayDatabase {
47            airports: parse_airport_heliport_zip_file(path.join("AirportHeliport.BASELINE.zip"))?,
48            navaids: parse_navaid_zip_file(path.join("Navaid.BASELINE.zip"))?,
49            designated_points: parse_designated_point_zip_file(path.join("DesignatedPoint.BASELINE.zip"))?,
50            route_segments: parse_route_segment_zip_file(path.join("RouteSegment.BASELINE.zip"))?,
51            routes: parse_route_zip_file(path.join("Route.BASELINE.zip"))?,
52            arrival_legs: if path.join("ArrivalLeg.BASELINE.zip").exists() {
53                parse_arrival_leg_zip_file(path.join("ArrivalLeg.BASELINE.zip"))?
54            } else {
55                HashMap::new()
56            },
57            departure_legs: if path.join("DepartureLeg.BASELINE.zip").exists() {
58                parse_departure_leg_zip_file(path.join("DepartureLeg.BASELINE.zip"))?
59            } else {
60                HashMap::new()
61            },
62            standard_instrument_arrivals: if path.join("StandardInstrumentArrival.BASELINE.zip").exists() {
63                parse_standard_instrument_arrival_zip_file(path.join("StandardInstrumentArrival.BASELINE.zip"))?
64            } else {
65                HashMap::new()
66            },
67            standard_instrument_departures: if path.join("StandardInstrumentDeparture.BASELINE.zip").exists() {
68                parse_standard_instrument_departure_zip_file(path.join("StandardInstrumentDeparture.BASELINE.zip"))?
69            } else {
70                HashMap::new()
71            },
72        })
73    }
74
75    /// Resolve SID connecting points by procedure designator.
76    pub fn resolve_sid_points(&self, name: &str) -> Vec<ResolvedPoint> {
77        let sid_ids = self
78            .standard_instrument_departures
79            .values()
80            .filter(|sid| sid.designator.trim().eq_ignore_ascii_case(name))
81            .map(|sid| sid.identifier.clone())
82            .collect::<std::collections::HashSet<_>>();
83
84        let fallback_points = self
85            .standard_instrument_departures
86            .values()
87            .filter(|sid| sid.designator.trim().eq_ignore_ascii_case(name))
88            .flat_map(|sid| sid.connecting_points.iter().cloned())
89            .collect::<Vec<_>>();
90
91        let exit_points = self.procedure_exit_points(&sid_ids, true);
92        let points_to_resolve = if exit_points.is_empty() {
93            fallback_points
94        } else {
95            exit_points
96        };
97
98        let mut points = points_to_resolve
99            .iter()
100            .map(|p| ResolvedPoint::from_db(p, self))
101            .filter(|p| !matches!(p, ResolvedPoint::None))
102            .collect::<Vec<_>>();
103        points.sort_by_key(|p| format!("{p:?}"));
104        points.dedup();
105        points
106    }
107
108    /// Resolve STAR connecting points by procedure designator.
109    pub fn resolve_star_points(&self, name: &str) -> Vec<ResolvedPoint> {
110        let star_ids = self
111            .standard_instrument_arrivals
112            .values()
113            .filter(|star| star.designator.trim().eq_ignore_ascii_case(name))
114            .map(|star| star.identifier.clone())
115            .collect::<std::collections::HashSet<_>>();
116
117        let fallback_points = self
118            .standard_instrument_arrivals
119            .values()
120            .filter(|star| star.designator.trim().eq_ignore_ascii_case(name))
121            .flat_map(|star| star.connecting_points.iter().cloned())
122            .collect::<Vec<_>>();
123
124        let entry_points = self.procedure_entry_points(&star_ids, false);
125        let points_to_resolve = if entry_points.is_empty() {
126            fallback_points
127        } else {
128            entry_points
129        };
130
131        let mut points = points_to_resolve
132            .iter()
133            .map(|p| ResolvedPoint::from_db(p, self))
134            .filter(|p| !matches!(p, ResolvedPoint::None))
135            .collect::<Vec<_>>();
136        points.sort_by_key(|p| format!("{p:?}"));
137        points.dedup();
138        points
139    }
140
141    /// Resolve SID procedures by designator as route-like segments.
142    pub fn resolve_sid_routes(&self, name: &str) -> Vec<ResolvedRoute> {
143        self.standard_instrument_departures
144            .values()
145            .filter(|sid| sid.designator.trim().eq_ignore_ascii_case(name))
146            .map(|sid| {
147                let segments = self
148                    .departure_legs
149                    .values()
150                    .filter(|leg| leg.departure.as_ref().is_some_and(|id| id == &sid.identifier))
151                    .filter_map(|leg| {
152                        let start = ResolvedPoint::from_db(&leg.start, self);
153                        let end = ResolvedPoint::from_db(&leg.end, self);
154                        if matches!(start, ResolvedPoint::None) || matches!(end, ResolvedPoint::None) {
155                            None
156                        } else {
157                            Some(ResolvedRouteSegment {
158                                start,
159                                end,
160                                name: Some(sid.designator.clone()),
161                                altitude: None,
162                                speed: None,
163                            })
164                        }
165                    })
166                    .collect::<Vec<_>>();
167                ResolvedRoute {
168                    segments: order_route_segments(segments),
169                    name: sid.designator.clone(),
170                }
171            })
172            .collect()
173    }
174
175    /// Resolve STAR procedures by designator as route-like segments.
176    pub fn resolve_star_routes(&self, name: &str) -> Vec<ResolvedRoute> {
177        self.standard_instrument_arrivals
178            .values()
179            .filter(|star| star.designator.trim().eq_ignore_ascii_case(name))
180            .map(|star| {
181                let segments = self
182                    .arrival_legs
183                    .values()
184                    .filter(|leg| leg.arrival.as_ref().is_some_and(|id| id == &star.identifier))
185                    .filter_map(|leg| {
186                        let start = ResolvedPoint::from_db(&leg.start, self);
187                        let end = ResolvedPoint::from_db(&leg.end, self);
188                        if matches!(start, ResolvedPoint::None) || matches!(end, ResolvedPoint::None) {
189                            None
190                        } else {
191                            Some(ResolvedRouteSegment {
192                                start,
193                                end,
194                                name: Some(star.designator.clone()),
195                                altitude: None,
196                                speed: None,
197                            })
198                        }
199                    })
200                    .collect::<Vec<_>>();
201                ResolvedRoute {
202                    segments: order_route_segments(segments),
203                    name: star.designator.clone(),
204                }
205            })
206            .collect()
207    }
208
209    fn procedure_exit_points(
210        &self,
211        procedure_ids: &std::collections::HashSet<String>,
212        is_departure: bool,
213    ) -> Vec<PointReference> {
214        let legs: Vec<(PointReference, PointReference)> = if is_departure {
215            self.departure_legs
216                .values()
217                .filter(|leg| leg.departure.as_ref().is_some_and(|id| procedure_ids.contains(id)))
218                .map(|leg| (leg.start.clone(), leg.end.clone()))
219                .collect()
220        } else {
221            self.arrival_legs
222                .values()
223                .filter(|leg| leg.arrival.as_ref().is_some_and(|id| procedure_ids.contains(id)))
224                .map(|leg| (leg.start.clone(), leg.end.clone()))
225                .collect()
226        };
227
228        let mut indegree: HashMap<String, usize> = HashMap::new();
229        let mut outdegree: HashMap<String, usize> = HashMap::new();
230        let mut refs: HashMap<String, PointReference> = HashMap::new();
231
232        for (start, end) in &legs {
233            let s = start.name();
234            let e = end.name();
235
236            if !s.is_empty() {
237                refs.insert(s.clone(), start.clone());
238                outdegree.entry(s).and_modify(|v| *v += 1).or_insert(1);
239            }
240            if !e.is_empty() {
241                refs.insert(e.clone(), end.clone());
242                indegree.entry(e).and_modify(|v| *v += 1).or_insert(1);
243            }
244        }
245
246        refs.into_iter()
247            .filter_map(|(name, point_ref)| {
248                if point_ref.is_airport_heliport() {
249                    return None;
250                }
251                let in_d = *indegree.get(&name).unwrap_or(&0);
252                let out_d = *outdegree.get(&name).unwrap_or(&0);
253                if out_d == 0 && in_d > 0 {
254                    Some(point_ref)
255                } else {
256                    None
257                }
258            })
259            .collect()
260    }
261
262    fn procedure_entry_points(
263        &self,
264        procedure_ids: &std::collections::HashSet<String>,
265        is_departure: bool,
266    ) -> Vec<PointReference> {
267        let legs: Vec<(PointReference, PointReference)> = if is_departure {
268            self.departure_legs
269                .values()
270                .filter(|leg| leg.departure.as_ref().is_some_and(|id| procedure_ids.contains(id)))
271                .map(|leg| (leg.start.clone(), leg.end.clone()))
272                .collect()
273        } else {
274            self.arrival_legs
275                .values()
276                .filter(|leg| leg.arrival.as_ref().is_some_and(|id| procedure_ids.contains(id)))
277                .map(|leg| (leg.start.clone(), leg.end.clone()))
278                .collect()
279        };
280
281        let mut indegree: HashMap<String, usize> = HashMap::new();
282        let mut outdegree: HashMap<String, usize> = HashMap::new();
283        let mut refs: HashMap<String, PointReference> = HashMap::new();
284
285        for (start, end) in &legs {
286            let s = start.name();
287            let e = end.name();
288
289            if !s.is_empty() {
290                refs.insert(s.clone(), start.clone());
291                outdegree.entry(s).and_modify(|v| *v += 1).or_insert(1);
292            }
293            if !e.is_empty() {
294                refs.insert(e.clone(), end.clone());
295                indegree.entry(e).and_modify(|v| *v += 1).or_insert(1);
296            }
297        }
298
299        refs.into_iter()
300            .filter_map(|(name, point_ref)| {
301                if point_ref.is_airport_heliport() {
302                    return None;
303                }
304                let in_d = *indegree.get(&name).unwrap_or(&0);
305                let out_d = *outdegree.get(&name).unwrap_or(&0);
306                if in_d == 0 && out_d > 0 {
307                    Some(point_ref)
308                } else {
309                    None
310                }
311            })
312            .collect()
313    }
314}
315
316fn order_route_segments(segments: Vec<ResolvedRouteSegment>) -> Vec<ResolvedRouteSegment> {
317    let mut out_map: HashMap<ResolvedPoint, Vec<ResolvedRouteSegment>> = HashMap::new();
318    let mut indegree: HashMap<ResolvedPoint, usize> = HashMap::new();
319
320    for segment in segments {
321        out_map.entry(segment.start.clone()).or_default().push(segment.clone());
322        indegree.entry(segment.start.clone()).or_insert(0);
323        indegree.entry(segment.end.clone()).and_modify(|v| *v += 1).or_insert(1);
324    }
325
326    for edges in out_map.values_mut() {
327        edges.sort_by_key(|e| format!("{}", e.end));
328    }
329
330    let mut starts = indegree
331        .iter()
332        .filter_map(|(point, deg)| if *deg == 0 { Some(point.clone()) } else { None })
333        .collect::<Vec<_>>();
334    starts.sort_by_key(|p| format!("{p}"));
335
336    let mut ordered = Vec::new();
337
338    let walk_from = |start: ResolvedPoint,
339                     out_map: &mut HashMap<ResolvedPoint, Vec<ResolvedRouteSegment>>,
340                     ordered: &mut Vec<ResolvedRouteSegment>| {
341        let mut current = start;
342        loop {
343            let next = out_map
344                .get_mut(&current)
345                .and_then(|edges| if edges.is_empty() { None } else { Some(edges.remove(0)) });
346            if let Some(segment) = next {
347                current = segment.end.clone();
348                ordered.push(segment);
349            } else {
350                break;
351            }
352        }
353    };
354
355    for start in starts {
356        walk_from(start, &mut out_map, &mut ordered);
357    }
358
359    loop {
360        let mut keys = out_map
361            .iter()
362            .filter_map(|(k, v)| if v.is_empty() { None } else { Some(k.clone()) })
363            .collect::<Vec<_>>();
364        if keys.is_empty() {
365            break;
366        }
367        keys.sort_by_key(|p| format!("{p}"));
368        walk_from(keys.remove(0), &mut out_map, &mut ordered);
369    }
370
371    ordered
372}
373
374const VALID_ROUTE_PREFIXES: [&str; 32] = [
375    "UN", "UM", "UL", "UT", "UZ", "UY", "UP", "UA", "UB", "UG", "UH", "UJ", "UQ", "UR", "UV", "UW", "L", "A", "B", "G",
376    "H", "J", "Q", "R", "T", "V", "W", "Y", "Z", "M", "N", "P",
377];
378
379/// The WGS84 ellipsoid.
380static WGS84: Lazy<Ellipsoid> = Lazy::new(|| Ellipsoid::named("WGS84").unwrap());
381
382/**
383 * A resolved route consisting of candidate route segments, based on their name.
384 */
385#[derive(Debug, Clone, Serialize)]
386pub struct ResolvedRoute {
387    pub segments: Vec<ResolvedRouteSegment>,
388    pub name: String,
389}
390
391/**
392 * A resolved route segment consisting of start and end points.
393 * Optionally, altitude and speed constraints can be included, propagated from modifiers.
394 */
395#[derive(Debug, Clone, Serialize)]
396pub struct ResolvedRouteSegment {
397    pub start: ResolvedPoint,
398    pub end: ResolvedPoint,
399    #[serde(skip_serializing_if = "Option::is_none")]
400    pub name: Option<String>,
401    #[serde(skip_serializing_if = "Option::is_none")]
402    pub altitude: Option<Altitude>,
403    #[serde(skip_serializing_if = "Option::is_none")]
404    pub speed: Option<Speed>,
405}
406
407/**
408 * A resolved point (based on their name), which can be a navaid, designated point, coordinates, or None.
409 */
410#[derive(Debug, Clone, Serialize)]
411#[serde(untagged)]
412pub enum ResolvedPoint {
413    AirportHeliport(AirportHeliport),
414    Navaid(Navaid),
415    DesignatedPoint(DesignatedPoint),
416    Coordinates { latitude: f64, longitude: f64 },
417    None,
418}
419
420impl std::fmt::Display for ResolvedPoint {
421    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
422        match self {
423            ResolvedPoint::AirportHeliport(airport) => {
424                write!(
425                    f,
426                    "AirportHeliport({}: {:.3},{:.3})",
427                    airport.icao, airport.latitude, airport.longitude
428                )
429            }
430            ResolvedPoint::Navaid(navaid) => write!(
431                f,
432                "Navaid({}: {:.3},{:.3})",
433                navaid.name.as_ref().unwrap(),
434                navaid.latitude,
435                navaid.longitude
436            ),
437            ResolvedPoint::DesignatedPoint(dp) => write!(
438                f,
439                "DesignatedPoint({}: {:.3}, {:.3})",
440                dp.designator, dp.latitude, dp.longitude
441            ),
442            ResolvedPoint::Coordinates { latitude, longitude } => {
443                write!(f, "Coordinates(lat: {}, lon: {})", latitude, longitude)
444            }
445            ResolvedPoint::None => write!(f, "None"),
446        }
447    }
448}
449
450impl From<&ResolvedPoint> for Coor2D {
451    fn from(val: &ResolvedPoint) -> Self {
452        match val {
453            ResolvedPoint::AirportHeliport(airport) => Coor2D::geo(airport.latitude, airport.longitude),
454            ResolvedPoint::Navaid(navaid) => Coor2D::geo(navaid.latitude, navaid.longitude),
455            ResolvedPoint::DesignatedPoint(dp) => Coor2D::geo(dp.latitude, dp.longitude),
456            ResolvedPoint::Coordinates { latitude, longitude } => Coor2D::geo(*latitude, *longitude),
457            ResolvedPoint::None => Coor2D::default(),
458        }
459    }
460}
461
462impl PartialEq for ResolvedPoint {
463    fn eq(&self, other: &Self) -> bool {
464        match (self, other) {
465            (ResolvedPoint::AirportHeliport(a), ResolvedPoint::AirportHeliport(b)) => a.identifier == b.identifier,
466            (ResolvedPoint::Navaid(a), ResolvedPoint::Navaid(b)) => a.identifier == b.identifier,
467            (ResolvedPoint::DesignatedPoint(a), ResolvedPoint::DesignatedPoint(b)) => a.identifier == b.identifier,
468            (
469                ResolvedPoint::Coordinates {
470                    latitude: a_lat,
471                    longitude: a_lon,
472                },
473                ResolvedPoint::Coordinates {
474                    latitude: b_lat,
475                    longitude: b_lon,
476                },
477            ) => (a_lat - b_lat).abs() < f64::EPSILON && (a_lon - b_lon).abs() < f64::EPSILON,
478            (ResolvedPoint::None, ResolvedPoint::None) => true,
479            _ => false,
480        }
481    }
482}
483
484impl Eq for ResolvedPoint {}
485
486// The hash trait implementation is needed for the DFS algorithm to reconstruct paths
487impl Hash for ResolvedPoint {
488    fn hash<H: std::hash::Hasher>(&self, state: &mut H) {
489        match self {
490            ResolvedPoint::AirportHeliport(airport) => {
491                airport.identifier.hash(state);
492            }
493            ResolvedPoint::Navaid(navaid) => {
494                navaid.identifier.hash(state);
495            }
496            ResolvedPoint::DesignatedPoint(dp) => {
497                dp.identifier.hash(state);
498            }
499            ResolvedPoint::Coordinates { latitude, longitude } => {
500                latitude.to_bits().hash(state);
501                longitude.to_bits().hash(state);
502            }
503            ResolvedPoint::None => {
504                0.hash(state);
505            }
506        }
507    }
508}
509
510impl ResolvedPoint {
511    /// Resolve a point from the database.
512    pub fn from_db(point: &PointReference, db: &AirwayDatabase) -> Self {
513        match point {
514            PointReference::AirportHeliport(id) => {
515                if let Some(airport) = db.airports.get(id) {
516                    ResolvedPoint::AirportHeliport(airport.clone())
517                } else {
518                    ResolvedPoint::None
519                }
520            }
521            PointReference::Navaid(id) => {
522                if let Some(navaid) = db.navaids.get(id) {
523                    ResolvedPoint::Navaid(navaid.clone())
524                } else {
525                    ResolvedPoint::None
526                }
527            }
528            PointReference::DesignatedPoint(id) => {
529                if let Some(dp) = db.designated_points.get(id) {
530                    ResolvedPoint::DesignatedPoint(dp.clone())
531                } else {
532                    ResolvedPoint::None
533                }
534            }
535            PointReference::None => ResolvedPoint::None,
536        }
537    }
538    /// Resolve a point by its name from the database.
539    pub fn lookup(name: &str, db: &AirwayDatabase) -> Vec<Self> {
540        let candidates = db
541            .navaids
542            .values()
543            .filter(|n| {
544                n.name
545                    .as_deref()
546                    .is_some_and(|n_name| n_name.trim().eq_ignore_ascii_case(name))
547            })
548            .collect::<Vec<_>>();
549        if !candidates.is_empty() {
550            return candidates.iter().map(|n| ResolvedPoint::Navaid((*n).clone())).collect();
551        }
552        let candidates = db
553            .designated_points
554            .values()
555            .filter(|dp| dp.designator.trim().eq_ignore_ascii_case(name))
556            .collect::<Vec<_>>();
557        if !candidates.is_empty() {
558            return candidates
559                .iter()
560                .map(|dp| ResolvedPoint::DesignatedPoint((*dp).clone()))
561                .collect();
562        }
563        vec![]
564    }
565}
566
567impl ResolvedRouteSegment {
568    /// Resolve a route segment from the database.
569    pub fn from_db(segment: &RouteSegment, db: &AirwayDatabase) -> Self {
570        ResolvedRouteSegment {
571            start: ResolvedPoint::from_db(&segment.start, db),
572            end: ResolvedPoint::from_db(&segment.end, db),
573            name: None,
574            altitude: None,
575            speed: None,
576        }
577    }
578}
579
580impl ResolvedRoute {
581    /// Resolve a route from the database.
582    pub fn from_db(route: &Route, db: &AirwayDatabase) -> Self {
583        let segments = db
584            .route_segments
585            .values()
586            .filter(|segment| segment.route_formed.as_deref() == Some(&route.identifier))
587            .map(|segment| ResolvedRouteSegment::from_db(segment, db))
588            .collect::<Vec<_>>();
589        ResolvedRoute {
590            segments,
591            name: format!(
592                "{}{}{}",
593                route.prefix.as_deref().unwrap_or(""),
594                route.second_letter.as_deref().unwrap_or(""),
595                route.number.as_deref().unwrap_or("")
596            ),
597        }
598    }
599
600    /// Lookup routes by their name from the database.
601    pub fn lookup(name: &str, db: &AirwayDatabase) -> Vec<Self> {
602        if VALID_ROUTE_PREFIXES.iter().any(|prefix| name.starts_with(prefix)) {
603            // First decompose the name into its components
604            // Another approach would be to make a single string match,
605            // but this serves as sanity check as well.
606            let last = name.chars().last().unwrap();
607            let (name, multiple) = if last.is_alphabetic() {
608                (&name[..name.len() - 1], Some(last.to_string()))
609            } else {
610                (name, None)
611            };
612            let (prefix, second_letter, number) = if name.starts_with('U') && name.len() >= 3 {
613                (
614                    Some("U".to_string()),
615                    name.chars().nth(1).unwrap().to_string(),
616                    name[2..].to_string(),
617                )
618            } else if name.len() >= 2 {
619                (None, name.chars().next().unwrap().to_string(), name[1..].to_string())
620            } else {
621                (None, "".to_string(), "".to_string())
622            };
623            let candidates = db
624                .routes
625                .values()
626                .filter(|route| {
627                    route.prefix.as_deref() == prefix.as_deref()
628                        && route.second_letter.as_deref() == Some(&second_letter)
629                        && route.number.as_deref() == Some(&number)
630                        && route.multiple_identifier.as_deref() == multiple.as_deref()
631                })
632                .collect::<Vec<_>>();
633            return candidates
634                .iter()
635                .map(|route| ResolvedRoute::from_db(route, db))
636                .collect();
637        }
638        vec![]
639    }
640
641    /// Check if the route contains the specified point.
642    pub fn contains(&self, point: &ResolvedPoint) -> bool {
643        self.segments
644            .iter()
645            .any(|segment| &segment.start == point || &segment.end == point)
646    }
647
648    /// Find a sub-route between two points, if it exists.
649    /// The implementation uses a depth-first search (DFS) algorithm to find a path
650    /// between the start and end points within the route segments.
651    pub fn between(&self, start: &ResolvedPoint, end: &ResolvedPoint) -> Option<ResolvedRoute> {
652        // Build adjacency map: point -> list of (next_point, segment_index, is_forward)
653        let mut graph: HashMap<&ResolvedPoint, Vec<(&ResolvedPoint, usize, bool)>> = HashMap::new();
654
655        for (i, segment) in self.segments.iter().enumerate() {
656            // Forward direction: start -> end
657            graph.entry(&segment.start).or_default().push((&segment.end, i, true));
658
659            // Backward direction: end -> start
660            graph.entry(&segment.end).or_default().push((&segment.start, i, false));
661        }
662
663        // Try DFS from start to end
664        if let Some(path) = Self::dfs_helper(
665            &graph,
666            start,
667            end,
668            &mut Vec::new(),
669            &mut std::collections::HashSet::new(),
670        ) {
671            return Some(self.build_route_from_path(path));
672        }
673
674        None
675    }
676
677    fn dfs_helper<'a>(
678        graph: &HashMap<&'a ResolvedPoint, Vec<(&'a ResolvedPoint, usize, bool)>>,
679        current: &'a ResolvedPoint,
680        target: &'a ResolvedPoint,
681        path: &mut Vec<(usize, bool)>,
682        visited: &mut std::collections::HashSet<usize>,
683    ) -> Option<Vec<(usize, bool)>> {
684        if current == target {
685            return Some(path.clone());
686        }
687
688        if let Some(neighbors) = graph.get(current) {
689            for (next_point, segment_idx, is_forward) in neighbors {
690                if !visited.contains(segment_idx) {
691                    visited.insert(*segment_idx);
692                    path.push((*segment_idx, *is_forward));
693
694                    if let Some(result) = Self::dfs_helper(graph, next_point, target, path, visited) {
695                        return Some(result);
696                    }
697
698                    path.pop();
699                    visited.remove(segment_idx);
700                }
701            }
702        }
703
704        None
705    }
706
707    fn build_route_from_path(&self, path: Vec<(usize, bool)>) -> ResolvedRoute {
708        let mut segments = Vec::new();
709
710        for (segment_idx, is_forward) in path {
711            let segment = &self.segments[segment_idx];
712            if is_forward {
713                segments.push(segment.clone());
714            } else {
715                // Reverse the segment
716                segments.push(ResolvedRouteSegment {
717                    start: segment.end.clone(),
718                    end: segment.start.clone(),
719                    name: Some(self.name.clone()),
720                    altitude: segment.altitude.clone(),
721                    speed: segment.speed.clone(),
722                });
723            }
724        }
725
726        ResolvedRoute {
727            segments,
728            name: self.name.clone(),
729        }
730    }
731}
732
733#[derive(Debug)]
734enum EnrichedCandidates {
735    Point((Vec<ResolvedPoint>, Option<Altitude>, Option<Speed>)),
736    PointCoords((ResolvedPoint, Option<Altitude>, Option<Speed>)),
737    Airway((Vec<ResolvedRoute>, String, Option<Altitude>, Option<Speed>)),
738    Direct(),
739}
740
741impl AirwayDatabase {
742    /// Enrich a sequence of Field15Elements into resolved route segments.
743    /// A resolved route segment consists of start and end points,
744    /// along with optional altitude and speed constraints.
745    /// All points and airways are resolved against the database.
746    pub fn enrich_route(&self, elements: Vec<Field15Element>) -> Vec<ResolvedRouteSegment> {
747        let mut altitude = None;
748        let mut speed = None;
749
750        // First, resolve all candidates
751        let mut resolved: Vec<EnrichedCandidates> = Vec::new();
752        for element in &elements {
753            match element {
754                Field15Element::Modifier(m) => {
755                    let Modifier {
756                        speed: s, altitude: a, ..
757                    } = m;
758                    altitude = a.clone();
759                    speed = s.clone();
760                }
761                Field15Element::Point(Point::Waypoint(name)) => {
762                    let lookup = ResolvedPoint::lookup(name, self);
763                    if lookup.is_empty() {
764                        tracing::warn!("No point found for identifier '{}'", name);
765                    }
766                    resolved.push(EnrichedCandidates::Point((
767                        ResolvedPoint::lookup(name, self),
768                        altitude.clone(),
769                        speed.clone(),
770                    )));
771                }
772                Field15Element::Point(Point::Coordinates((lat, lon))) => {
773                    resolved.push(EnrichedCandidates::PointCoords((
774                        ResolvedPoint::Coordinates {
775                            latitude: *lat,
776                            longitude: *lon,
777                        },
778                        altitude.clone(),
779                        speed.clone(),
780                    )));
781                }
782                Field15Element::Connector(Connector::Airway(name)) => {
783                    let lookup = ResolvedRoute::lookup(name, self);
784                    if lookup.is_empty() {
785                        tracing::warn!("No airway found for identifier '{}'", name);
786                        resolved.push(EnrichedCandidates::Direct());
787                    } else {
788                        resolved.push(EnrichedCandidates::Airway((
789                            lookup,
790                            name.to_string(),
791                            altitude.clone(),
792                            speed.clone(),
793                        )));
794                    }
795                }
796                Field15Element::Connector(Connector::Sid(name)) => {
797                    let lookup = self.resolve_sid_routes(name);
798                    if lookup.is_empty() {
799                        tracing::warn!("No SID found for identifier '{}'", name);
800                        resolved.push(EnrichedCandidates::Direct());
801                    } else {
802                        resolved.push(EnrichedCandidates::Airway((
803                            lookup,
804                            name.to_string(),
805                            altitude.clone(),
806                            speed.clone(),
807                        )));
808                    }
809                }
810                Field15Element::Connector(Connector::Star(name)) => {
811                    let lookup = self.resolve_star_routes(name);
812                    if lookup.is_empty() {
813                        tracing::warn!("No STAR found for identifier '{}'", name);
814                        resolved.push(EnrichedCandidates::Direct());
815                    } else {
816                        resolved.push(EnrichedCandidates::Airway((
817                            lookup,
818                            name.to_string(),
819                            altitude.clone(),
820                            speed.clone(),
821                        )));
822                    }
823                }
824                Field15Element::Connector(Connector::Direct) => {
825                    resolved.push(EnrichedCandidates::Direct());
826                }
827                Field15Element::Connector(Connector::Nat(_)) | Field15Element::Connector(Connector::Pts(_)) => {
828                    // NAT and PTS are not handled yet
829                    resolved.push(EnrichedCandidates::Direct());
830                }
831                _ => {}
832            }
833        }
834
835        // 1. For each candidate airway, retain only those that contain both the previous and next point.
836        for i in 1..resolved.len() - 1 {
837            let (before_i, i_and_after) = &mut resolved.split_at_mut(i);
838            if let (EnrichedCandidates::Airway((routes, _, _, _)), after_i) = i_and_after.split_first_mut().unwrap() {
839                tracing::debug!("Filtering airway candidates: {:?}", routes);
840                if let Some(EnrichedCandidates::Point((points, _, _))) = before_i.last() {
841                    routes.retain(|r| points.iter().any(|p| r.contains(p)));
842                    tracing::debug!("Filtering airway candidates with point {:?}: {:?}", points, routes);
843                }
844                if let Some(EnrichedCandidates::Point((points, _, _))) = after_i.first() {
845                    routes.retain(|r| points.iter().any(|p| r.contains(p)));
846                    tracing::debug!("Filtering airway candidates with point {:?}: {:?}", points, routes);
847                }
848            }
849        }
850
851        // 2. Transform now empty airway candidates to Direct
852        for candidate in resolved.iter_mut() {
853            if let EnrichedCandidates::Airway((routes, name, _, _)) = candidate {
854                if routes.is_empty() {
855                    tracing::warn!("No valid airway remaining for '{}'", name);
856                    *candidate = EnrichedCandidates::Direct();
857                }
858            }
859        }
860
861        // 3. For each point, retain only those that are present in the adjacent airway segments.
862        for i in 0..resolved.len() {
863            let (before_i, i_and_after) = &mut resolved.split_at_mut(i);
864            if let (EnrichedCandidates::Point((points, _, _)), after_i) = i_and_after.split_first_mut().unwrap() {
865                tracing::debug!("Filtering point candidates: {:?}", points);
866                if let Some(EnrichedCandidates::Airway((routes, _, _, _))) = before_i.last() {
867                    points.retain(|p| routes.iter().any(|r| r.contains(p)));
868                    tracing::debug!("Filtering point candidates with airway {:?}: {:?}", routes, points);
869                }
870                if let Some(EnrichedCandidates::Airway((routes, _, _, _))) = after_i.first() {
871                    points.retain(|p| routes.iter().any(|r| r.contains(p)));
872                    tracing::debug!("Filtering point candidates with airway {:?}: {:?}", routes, points);
873                }
874            }
875        }
876
877        // 4. Trim airways to the segments between the before and after points.
878        for i in 1..resolved.len() - 1 {
879            let (before_i, i_and_after) = &mut resolved.split_at_mut(i);
880            if let (EnrichedCandidates::Airway((routes, _, _, _)), after_i) = i_and_after.split_first_mut().unwrap() {
881                if let Some(EnrichedCandidates::Point((before, _, _))) = before_i.last() {
882                    if let Some(EnrichedCandidates::Point((after, _, _))) = after_i.first() {
883                        if let Some(before) = before.first() {
884                            if let Some(after) = after.first() {
885                                for route in routes.iter_mut() {
886                                    if let Some(trimmed) = route.between(before, after) {
887                                        *route = trimmed;
888                                        tracing::debug!(
889                                            "Trimmed airway '{}' between points {} and {}: {:?}",
890                                            route.name,
891                                            before,
892                                            after,
893                                            *route
894                                        );
895                                    }
896                                }
897                            }
898                        }
899                    }
900                }
901            }
902        }
903
904        // 5. Replace empty routes with Direct
905        for candidate in resolved.iter_mut() {
906            if let EnrichedCandidates::Airway((routes, name, _, _)) = candidate {
907                if routes.iter().all(|r| r.segments.is_empty()) {
908                    tracing::warn!("No valid segments remaining for airway '{}'", name);
909                    *candidate = EnrichedCandidates::Direct();
910                }
911            }
912        }
913
914        // 6. Break the tie for remaining multiple candidate points
915        let mut last_known: Option<ResolvedPoint> = None;
916
917        for i in 0..resolved.len() {
918            if let EnrichedCandidates::Point((points, _, _)) = &resolved[i] {
919                if points.len() > 1 {
920                    // Find the next definitive point ahead
921                    let mut next_definitive: Option<&ResolvedPoint> = None;
922                    for candidate in resolved[i..].iter() {
923                        match candidate {
924                            EnrichedCandidates::Point((pts, _, _)) if pts.len() == 1 => {
925                                next_definitive = pts.first();
926                                break;
927                            }
928                            EnrichedCandidates::PointCoords((pt, _, _)) => {
929                                next_definitive = Some(pt);
930                                break;
931                            }
932                            _ => {}
933                        }
934                    }
935
936                    match (&last_known, next_definitive) {
937                        (None, None) => {
938                            tracing::warn!("Cannot disambiguate point {:?}: no reference points available", points);
939                        }
940                        (None, Some(_)) => {
941                            tracing::info!("Disambiguating point {:?} using only next definitive point", points);
942                        }
943                        (Some(a), None) => {
944                            tracing::info!("Disambiguating point {:?} using only last known point", points);
945
946                            // Only last known point is available, pick the closest candidate
947                            let mut best_idx = 0;
948                            let mut best_distance = f64::INFINITY;
949
950                            for (idx, candidate) in points.iter().enumerate() {
951                                let distance =
952                                    WGS84.distance(&Into::<Coor2D>::into(a), &Into::<Coor2D>::into(candidate));
953                                if distance < best_distance {
954                                    best_distance = distance;
955                                    best_idx = idx;
956                                }
957                            }
958
959                            // Keep only the best candidate
960                            if let EnrichedCandidates::Point((points, _, _)) = &mut resolved[i] {
961                                let best = points[best_idx].clone();
962                                points.clear();
963                                points.push(best);
964                            }
965                        }
966                        (Some(a), Some(b)) => {
967                            tracing::info!("Disambiguating point {:?} using both reference points", points);
968
969                            let mut best_idx = 0;
970                            let mut best_score = f64::INFINITY;
971
972                            for (idx, candidate) in points.iter().enumerate() {
973                                tracing::info!("Scoring candidate {}: {} ({}-{})", idx, candidate, a, b);
974                                let score = score_hybrid(&a.into(), &b.into(), &candidate.into());
975                                if score < best_score {
976                                    best_score = score;
977                                    best_idx = idx;
978                                }
979                            }
980
981                            // Keep only the best candidate
982                            if let EnrichedCandidates::Point((points, _, _)) = &mut resolved[i] {
983                                let best = points[best_idx].clone();
984                                points.clear();
985                                points.push(best);
986                            }
987                        }
988                    }
989                }
990
991                // Update last_known point
992                if let EnrichedCandidates::Point((pts, _, _)) = &resolved[i] {
993                    if let Some(pt) = pts.first() {
994                        last_known = Some(pt.clone());
995                    }
996                }
997            } else if let EnrichedCandidates::PointCoords((pt, _, _)) = &resolved[i] {
998                last_known = Some(pt.clone());
999            }
1000        }
1001
1002        // 7. Build the final sequence of resolved route segments.
1003        let mut segments = Vec::new();
1004        let mut previous_point: Option<ResolvedPoint> = None;
1005
1006        for enriched in resolved {
1007            match enriched {
1008                EnrichedCandidates::Point((points, alt, spd)) => {
1009                    if let Some(point) = points.first() {
1010                        if let Some(prev) = &previous_point {
1011                            if prev == point {
1012                                continue;
1013                            }
1014                            segments.push(ResolvedRouteSegment {
1015                                start: prev.clone(),
1016                                end: point.clone(),
1017                                name: None,
1018                                altitude: alt,
1019                                speed: spd,
1020                            });
1021                        }
1022                        previous_point = Some(point.clone());
1023                    }
1024                }
1025                EnrichedCandidates::PointCoords((point, alt, spd)) => {
1026                    if let Some(prev) = previous_point {
1027                        segments.push(ResolvedRouteSegment {
1028                            start: prev,
1029                            end: point.clone(),
1030                            name: None,
1031                            altitude: alt,
1032                            speed: spd,
1033                        });
1034                    }
1035                    previous_point = Some(point.clone());
1036                }
1037                EnrichedCandidates::Airway((routes, name, alt, spd)) => {
1038                    if let Some(route) = routes.first() {
1039                        for segment in &route.segments {
1040                            segments.push(ResolvedRouteSegment {
1041                                start: segment.start.clone(),
1042                                end: segment.end.clone(),
1043                                name: Some(name.clone()),
1044                                altitude: alt.clone(),
1045                                speed: spd.clone(),
1046                            });
1047                        }
1048                        previous_point = Some(route.segments.last().unwrap().end.clone());
1049                    }
1050                }
1051                EnrichedCandidates::Direct() => {
1052                    // Direct segments are handled by just carrying forward the previous point
1053                }
1054            }
1055        }
1056        segments
1057    }
1058}
1059
1060fn score_hybrid(a: &Coor2D, b: &Coor2D, x: &Coor2D) -> f64 {
1061    // Ideally gap_ration is close to 1.0 and the bearing difference close to 0.0
1062    let ab = WGS84.geodesic_inv(a, b).to_degrees();
1063    let ax = WGS84.geodesic_inv(a, x).to_degrees();
1064    let xb = WGS84.geodesic_inv(x, b).to_degrees();
1065
1066    // Think about triangular inequality, we want x to be "between" a and b
1067    let gap_ratio = (ax[2] + xb[2]) / ab[2].max(1e-9);
1068
1069    let delta_a = (ax[0] - ab[0]).abs().min(360.0 - (ax[0] - ab[0]).abs());
1070    let delta_b = (xb[0] - ab[0]).abs().min(360.0 - (xb[0] - ab[0]).abs());
1071    let bearing_diff = (delta_a + delta_b) / 2.0; // Normalize to [0,1]
1072
1073    tracing::info!(
1074        "Scoring point: {} = {} + {}; bearing_diff = {:.3}, gap_ratio = {:.3}",
1075        ab[2],
1076        ax[2],
1077        xb[2],
1078        bearing_diff,
1079        gap_ratio
1080    );
1081    // Combine the two metrics into a score
1082    bearing_diff / 180. + (gap_ratio - 1.0).max(0.)
1083}