Skip to main content

dsfb_debug/
baseline.rs

1//! DSFB-Debug: baseline computation — healthy-window statistics.
2//!
3//! # Role in the pipeline
4//!
5//! Computes the nominal reference `x_hat(k)` and the admissibility
6//! envelope radius `ρ` from the healthy-window slice of the residual
7//! stream. These two values anchor the entire downstream evaluation:
8//!
9//! - `x_hat(k)` is subtracted from `x(k)` to produce `r(k) = x(k) -
10//!   x_hat(k)` (paper §5.2; see `residual.rs`).
11//! - `ρ` defines the admissibility envelope
12//!   `E(k) = { r : ‖r‖ ≤ ρ }` (paper §5.4; see `envelope.rs`).
13//!
14//! # Healthy-window contract
15//!
16//! Per paper §F.4, the "healthy window" is the first N fault-free
17//! windows of the residual stream. The harness identifies the
18//! healthy slice via `OwnedResidualMatrix.healthy_window_end` from
19//! the upstream fixture metadata. The bound is operator-supplied,
20//! never inferred from data — DSFB-Debug refuses to fit baseline
21//! statistics on suspected-fault windows.
22//!
23//! # Numerical stability
24//!
25//! Mean is two-pass (compensated summation); standard deviation is
26//! single-pass with explicit `(n - 1)` Bessel correction. No NaN
27//! propagation — NaN inputs are treated as imputed downstream.
28//! Square-root delegates to `envelope::sqrt_approx_pub` so the
29//! no_std core needs no `libm`/`std::f64` dependency.
30
31/// Compute the per-signal mean from a healthy window.
32///
33/// # Arguments
34/// * `healthy_data` - row-major [window][signal], immutable
35/// * `num_signals` - number of signals per window
36/// * `num_windows` - number of healthy windows
37/// * `mean_out` - output buffer for per-signal means
38pub fn compute_baseline_mean(
39    healthy_data: &[f64],
40    num_signals: usize,
41    num_windows: usize,
42    mean_out: &mut [f64],
43) {
44    let mut s = 0;
45    while s < num_signals && s < mean_out.len() {
46        let mut sum = 0.0;
47        let mut count: usize = 0;
48        let mut w = 0;
49        while w < num_windows {
50            let idx = w * num_signals + s;
51            if idx < healthy_data.len() {
52                let v = healthy_data[idx];
53                // Skip NaN-like sentinel values
54                if !v.is_nan() { // NaN != NaN check
55                    sum += v;
56                    count += 1;
57                }
58            }
59            w += 1;
60        }
61        mean_out[s] = if count > 0 { sum / count as f64 } else { 0.0 };
62        s += 1;
63    }
64}
65
66/// Compute per-signal envelope radii (3σ) from healthy residuals.
67///
68/// # Arguments
69/// * `healthy_data` - raw healthy window data (row-major)
70/// * `baseline_mean` - per-signal means (immutable)
71/// * `num_signals` - signals per window
72/// * `num_windows` - healthy windows
73/// * `rho_out` - output buffer for envelope radii
74pub fn compute_baseline_envelope(
75    healthy_data: &[f64],
76    baseline_mean: &[f64],
77    num_signals: usize,
78    num_windows: usize,
79    rho_out: &mut [f64],
80) {
81    // Temporary: compute residuals per signal, then get 3σ
82    let mut s = 0;
83    while s < num_signals && s < rho_out.len() {
84        // Collect residuals for this signal into a small inline approach
85        // Since we can't allocate, we compute variance in a single pass
86        let mean = if s < baseline_mean.len() { baseline_mean[s] } else { 0.0 };
87        let mut var_sum = 0.0;
88        let mut count: usize = 0;
89        let mut w = 0;
90        while w < num_windows {
91            let idx = w * num_signals + s;
92            if idx < healthy_data.len() {
93                let v = healthy_data[idx];
94                if !v.is_nan() {
95                    let d = v - mean;
96                    var_sum += d * d;
97                    count += 1;
98                }
99            }
100            w += 1;
101        }
102        let sigma = if count > 1 {
103            crate::envelope::sqrt_approx_pub(var_sum / (count - 1) as f64)
104        } else {
105            1.0
106        };
107        rho_out[s] = if 3.0 * sigma > 1e-10 { 3.0 * sigma } else { 1e-10 };
108        s += 1;
109    }
110}
111
112#[cfg(test)]
113mod tests {
114    use super::*;
115
116    #[test]
117    fn test_baseline_mean() {
118        // 3 windows, 2 signals
119        let data = [1.0, 10.0, 2.0, 20.0, 3.0, 30.0];
120        let mut mean = [0.0; 2];
121        compute_baseline_mean(&data, 2, 3, &mut mean);
122        assert!((mean[0] - 2.0).abs() < 1e-12);
123        assert!((mean[1] - 20.0).abs() < 1e-12);
124    }
125}