Skip to main content

winds_engine/
validate.rs

1//! Compare a wind profile against a measured sounding.
2//!
3//! Metrics are sampled at the sounding's own levels: speed RMSE, direction
4//! RMSE (circular), and mean absolute speed error per altitude layer
5//! (1, 3, 6 and 10 km boundaries, then the rest).
6
7use crate::field::{ProfileField, WindField};
8use crate::sounding::Sounding;
9
10/// Altitude-layer boundaries, metres: 0–1 km, 1–3 km, 3–6 km, 6–10 km, rest.
11const LAYER_TOPS_M: [f64; 4] = [1_000.0, 3_000.0, 6_000.0, 10_000.0];
12
13/// Mean absolute speed error within one altitude layer.
14#[derive(Debug, Clone, Copy, PartialEq)]
15pub struct LayerError {
16    /// Top of the layer, metres (`f64::INFINITY` for the open layer).
17    pub top_m: f64,
18    /// Mean absolute speed error, m/s.
19    pub mean_abs_speed_error_mps: f64,
20    /// Sounding levels in the layer.
21    pub levels: usize,
22}
23
24/// The comparison result.
25#[derive(Debug, Clone, PartialEq)]
26pub struct ValidationReport {
27    /// Sounding levels compared.
28    pub levels: usize,
29    /// RMS of (profile speed − measured speed), m/s.
30    pub speed_rmse_mps: f64,
31    /// RMS of the circular direction difference, degrees.
32    pub direction_rmse_deg: f64,
33    /// Per-layer mean absolute speed errors.
34    pub layers: Vec<LayerError>,
35}
36
37/// Sample `profile` at every sounding level and compute the metrics.
38pub fn validate(profile: &ProfileField, sounding: &Sounding) -> ValidationReport {
39    let mut speed_sq = 0.0;
40    let mut direction_sq = 0.0;
41    let mut layer_abs = vec![0.0; LAYER_TOPS_M.len() + 1];
42    let mut layer_n = vec![0usize; LAYER_TOPS_M.len() + 1];
43    for level in &sounding.levels {
44        let wind = profile.wind(0.0, [0.0, level.altitude_m, 0.0]);
45        let speed_err = wind.speed_mps() - level.speed_mps;
46        speed_sq += speed_err * speed_err;
47        let mut dir_err = (wind.direction_from_deg() - level.direction_from_deg) % 360.0;
48        if dir_err > 180.0 {
49            dir_err -= 360.0;
50        } else if dir_err < -180.0 {
51            dir_err += 360.0;
52        }
53        direction_sq += dir_err * dir_err;
54        let layer = LAYER_TOPS_M
55            .iter()
56            .position(|top| level.altitude_m < *top)
57            .unwrap_or(LAYER_TOPS_M.len());
58        layer_abs[layer] += speed_err.abs();
59        layer_n[layer] += 1;
60    }
61    let n = sounding.levels.len().max(1) as f64;
62    ValidationReport {
63        levels: sounding.levels.len(),
64        speed_rmse_mps: (speed_sq / n).sqrt(),
65        direction_rmse_deg: (direction_sq / n).sqrt(),
66        layers: layer_abs
67            .into_iter()
68            .zip(layer_n)
69            .enumerate()
70            .filter(|(_, (_, count))| *count > 0)
71            .map(|(i, (abs, count))| LayerError {
72                top_m: LAYER_TOPS_M.get(i).copied().unwrap_or(f64::INFINITY),
73                mean_abs_speed_error_mps: abs / count as f64,
74                levels: count,
75            })
76            .collect(),
77    }
78}