winds_engine/
field.rs

1//! The [`WindField`] trait and the concrete field types.
2
3use crate::{Error, Result};
4
5/// World-frame air velocity, m/s: +X is downrange (tailwind), +Y is up, +Z
6/// is the horizontal cross axis to the right of +X (a right-handed frame,
7/// ballistics-engine's McCoy frame). Vertical wind is zero for the horizontal
8/// field types in this crate.
9#[derive(Debug, Clone, Copy, PartialEq)]
10pub struct Wind {
11    /// Downrange component, m/s.
12    pub x_mps: f64,
13    /// Vertical component, m/s.
14    pub y_mps: f64,
15    /// Cross-range component, m/s.
16    pub z_mps: f64,
17}
18
19impl Wind {
20    /// Calm air.
21    pub const ZERO: Self = Self {
22        x_mps: 0.0,
23        y_mps: 0.0,
24        z_mps: 0.0,
25    };
26
27    /// Horizontal wind from a speed and a meteorological direction: the
28    /// direction the air comes FROM, degrees clockwise, where 0° means from
29    /// downrange (+X, so the wind blows toward −X, a headwind) and 90° means
30    /// from +Z. Vertical component is zero.
31    pub fn from_speed_direction(speed_mps: f64, direction_from_deg: f64) -> Self {
32        let toward = (direction_from_deg + 180.0).to_radians();
33        Self {
34            x_mps: speed_mps * toward.cos(),
35            y_mps: 0.0,
36            z_mps: speed_mps * toward.sin(),
37        }
38    }
39
40    /// Horizontal speed.
41    pub fn speed_mps(&self) -> f64 {
42        self.x_mps.hypot(self.z_mps)
43    }
44
45    /// Meteorological direction the air comes FROM, degrees in [0, 360).
46    pub fn direction_from_deg(&self) -> f64 {
47        let toward = self.z_mps.atan2(self.x_mps).to_degrees() + 180.0;
48        toward.rem_euclid(360.0)
49    }
50}
51
52/// A deterministic, stateless wind field: same time and position always
53/// give the same wind, so integrator substeps, snapshots and seeded
54/// reproduction are exact.
55pub trait WindField {
56    /// World-frame air velocity at a time (seconds) and position
57    /// (`[x, y, z]` metres, Y up).
58    fn wind(&self, time_s: f64, position_m: [f64; 3]) -> Wind;
59}
60
61/// Uniform wind everywhere.
62#[derive(Debug, Clone, Copy)]
63pub struct ConstantField {
64    /// The constant air velocity.
65    pub wind: Wind,
66}
67
68impl WindField for ConstantField {
69    fn wind(&self, _time_s: f64, _position_m: [f64; 3]) -> Wind {
70        self.wind
71    }
72}
73
74/// One level of a vertical profile.
75#[derive(Debug, Clone, Copy, PartialEq)]
76pub struct ProfileLevel {
77    /// Height above the profile datum, metres.
78    pub altitude_m: f64,
79    /// Horizontal wind speed, m/s (nonnegative).
80    pub speed_mps: f64,
81    /// Meteorological direction the air comes FROM, degrees (see
82    /// [`Wind::from_speed_direction`]).
83    pub direction_from_deg: f64,
84}
85
86/// Piecewise-linear vertical profile: speed interpolates linearly in
87/// altitude, direction along the shortest circular arc. Below the first
88/// level and above the last, the nearest level holds.
89#[derive(Debug, Clone)]
90pub struct ProfileField {
91    levels: Vec<ProfileLevel>,
92}
93
94impl ProfileField {
95    /// Build from levels, which must be finite, nonnegative in speed and
96    /// strictly ascending in altitude.
97    pub fn new(mut levels: Vec<ProfileLevel>) -> Result<Self> {
98        if levels.is_empty() {
99            return Err(Error::invalid("profile", "needs at least one level"));
100        }
101        levels.sort_by(|a, b| a.altitude_m.total_cmp(&b.altitude_m));
102        for pair in levels.windows(2) {
103            if pair[0].altitude_m == pair[1].altitude_m {
104                return Err(Error::invalid(
105                    "profile",
106                    "altitudes must be strictly ascending",
107                ));
108            }
109        }
110        for level in &levels {
111            if !level.altitude_m.is_finite()
112                || !level.speed_mps.is_finite()
113                || !level.direction_from_deg.is_finite()
114            {
115                return Err(Error::invalid("profile", "levels must be finite"));
116            }
117            if level.speed_mps < 0.0 {
118                return Err(Error::invalid("profile", "speeds must be nonnegative"));
119            }
120        }
121        Ok(Self { levels })
122    }
123
124    /// The sorted levels.
125    pub fn levels(&self) -> &[ProfileLevel] {
126        &self.levels
127    }
128
129    /// Wind at an altitude, with circular direction interpolation.
130    pub fn wind_at(&self, altitude_m: f64) -> Wind {
131        let levels = &self.levels;
132        if altitude_m <= levels[0].altitude_m {
133            let l = levels[0];
134            return Wind::from_speed_direction(l.speed_mps, l.direction_from_deg);
135        }
136        for pair in levels.windows(2) {
137            let (lo, hi) = (pair[0], pair[1]);
138            if altitude_m <= hi.altitude_m {
139                let f = (altitude_m - lo.altitude_m) / (hi.altitude_m - lo.altitude_m);
140                let speed = lo.speed_mps + f * (hi.speed_mps - lo.speed_mps);
141                let mut delta = (hi.direction_from_deg - lo.direction_from_deg) % 360.0;
142                if delta > 180.0 {
143                    delta -= 360.0;
144                } else if delta < -180.0 {
145                    delta += 360.0;
146                }
147                let direction = lo.direction_from_deg + f * delta;
148                return Wind::from_speed_direction(speed, direction);
149            }
150        }
151        let l = levels[levels.len() - 1];
152        Wind::from_speed_direction(l.speed_mps, l.direction_from_deg)
153    }
154}
155
156impl WindField for ProfileField {
157    fn wind(&self, _time_s: f64, position_m: [f64; 3]) -> Wind {
158        self.wind_at(position_m[1])
159    }
160}
winds_engine/field.rs

winds_engine/
field.rs
