aerocontext-core 0.2.1

Provider-neutral aeronautical-context model and the pluggable ContextProvider contract
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80

//! Spherical geodesy on [`GeoPoint`]s.
//!
//! Aviation-scale great-circle math on a spherical Earth (mean radius in
//! nautical miles). Accurate to well under 0.5% versus the ellipsoid —
//! appropriate for corridor construction and spatial association, not for
//! survey-grade work.

use crate::model::GeoPoint;

/// Mean Earth radius in nautical miles.
pub const EARTH_RADIUS_NM: f64 = 3_440.065;

/// Great-circle distance between two points, nautical miles (haversine).
pub fn distance_nm(a: GeoPoint, b: GeoPoint) -> f64 {
    let (lat_a, lat_b) = (a.lat.to_radians(), b.lat.to_radians());
    let half_dlat = (lat_b - lat_a) / 2.0;
    let half_dlon = (b.lon - a.lon).to_radians() / 2.0;
    let h = half_dlat.sin().powi(2) + lat_a.cos() * lat_b.cos() * half_dlon.sin().powi(2);
    2.0 * EARTH_RADIUS_NM * h.sqrt().asin()
}

/// Initial great-circle bearing from `a` to `b`, degrees in `0..360`.
pub fn initial_bearing_deg(a: GeoPoint, b: GeoPoint) -> f64 {
    let (lat_a, lat_b) = (a.lat.to_radians(), b.lat.to_radians());
    let dlon = (b.lon - a.lon).to_radians();
    let y = dlon.sin() * lat_b.cos();
    let x = lat_a.cos() * lat_b.sin() - lat_a.sin() * lat_b.cos() * dlon.cos();
    (y.atan2(x).to_degrees() + 360.0) % 360.0
}

/// The point `distance_nm` out from `from` on the given initial bearing.
pub fn destination(from: GeoPoint, bearing_deg: f64, distance_nm: f64) -> GeoPoint {
    let angular = distance_nm / EARTH_RADIUS_NM;
    let bearing = bearing_deg.to_radians();
    let lat_from = from.lat.to_radians();
    let lat_to =
        (lat_from.sin() * angular.cos() + lat_from.cos() * angular.sin() * bearing.cos()).asin();
    let lon_to = from.lon.to_radians()
        + (bearing.sin() * angular.sin() * lat_from.cos())
            .atan2(angular.cos() - lat_from.sin() * lat_to.sin());
    GeoPoint {
        lat: lat_to.to_degrees(),
        lon: normalize_lon(lon_to.to_degrees()),
    }
}

/// Signed cross-track distance of `point` from the great circle through
/// `leg_start` → `leg_end`, nautical miles. Negative is left of track.
pub fn cross_track_nm(point: GeoPoint, leg_start: GeoPoint, leg_end: GeoPoint) -> f64 {
    let dist_sp = distance_nm(leg_start, point) / EARTH_RADIUS_NM;
    let bearing_sp = initial_bearing_deg(leg_start, point).to_radians();
    let bearing_se = initial_bearing_deg(leg_start, leg_end).to_radians();
    (dist_sp.sin() * (bearing_sp - bearing_se).sin()).asin() * EARTH_RADIUS_NM
}

/// Signed along-track distance of `point` projected onto the leg
/// `leg_start` → `leg_end`, nautical miles from `leg_start`; negative
/// when the projection falls behind the start.
pub fn along_track_nm(point: GeoPoint, leg_start: GeoPoint, leg_end: GeoPoint) -> f64 {
    let dist_sp = distance_nm(leg_start, point) / EARTH_RADIUS_NM;
    let cross = cross_track_nm(point, leg_start, leg_end) / EARTH_RADIUS_NM;
    let magnitude = (dist_sp.cos() / cross.cos()).clamp(-1.0, 1.0).acos() * EARTH_RADIUS_NM;
    let relative = (initial_bearing_deg(leg_start, point)
        - initial_bearing_deg(leg_start, leg_end))
    .to_radians();
    if relative.cos() < 0.0 {
        -magnitude
    } else {
        magnitude
    }
}

/// Wrap a longitude into `-180..=180`.
pub fn normalize_lon(lon: f64) -> f64 {
    let wrapped = (lon + 180.0).rem_euclid(360.0) - 180.0;
    if wrapped == -180.0 { 180.0 } else { wrapped }
}

#[cfg(test)]
mod tests;