Skip to main content

harmos_storage/signal/waveform/
summary.rs

1/* ***
2 *
3 * Harmos
4 *
5 * SPDX-License-Identifier: Apache-2.0
6 *
7 * ***
8 */
9
10/* === Modules === */
11
12/* === Re-Exports === */
13
14/* === Imports === */
15
16// Standard Library Imports
17
18// Third Party Imports
19
20// Local Imports
21
22/* === Constants === */
23
24/* === Macros === */
25
26/* === Type Aliases === */
27
28/* === Structs === */
29
30// The streaming, mergeable reduction: count, sum, extremes, mean and spread
31// kept in constant space, so partial summaries of disjoint runs merge into the
32// summary of their union. Non-finite values are skipped, the same values
33// `Series` skips, so both answer the shared statistics identically.
34#[derive(Clone, Copy, Debug, PartialEq)]
35pub struct Summary {
36    count: u64,
37    sum: f64,
38    mean: f64,
39    lowest: f64,
40    highest: f64,
41    deviations: f64,
42}
43
44/* === Enums === */
45
46/* === Traits === */
47
48/* === Implementations === */
49
50impl Default for Summary {
51    fn default() -> Self {
52        Self { count: 0, sum: 0.0, mean: 0.0, lowest: f64::INFINITY, highest: f64::NEG_INFINITY, deviations: 0.0 }
53    }
54}
55
56impl Summary {
57    #[must_use]
58    pub fn new() -> Self {
59        Self::default()
60    }
61
62    #[must_use]
63    pub fn of(values: &[f64]) -> Self {
64        let mut summary = Self::new();
65        for value in values {
66            summary.push(*value);
67        }
68        summary
69    }
70
71    pub fn push(&mut self, value: f64) {
72        if !value.is_finite() {
73            return;
74        }
75        let before = self.mean;
76
77        self.count += 1;
78        self.sum += value;
79        self.lowest = self.lowest.min(value);
80        self.highest = self.highest.max(value);
81
82        // The running mean never forms the sum, so it stays finite wherever
83        // the values do; opposite signs are weighted apart so the step itself
84        // cannot overflow.
85        let n = self.count as f64;
86        self.mean = if self.count == 1 {
87            value
88        } else if before.signum() == value.signum() {
89            before + (value - before) / n
90        } else {
91            before * ((n - 1.0) / n) + value / n
92        };
93        self.deviations += (value - before) * (value - self.mean);
94    }
95
96    #[must_use]
97    pub fn merge(&self, other: &Self) -> Self {
98        if self.count == 0 {
99            return *other;
100        }
101        if other.count == 0 {
102            return *self;
103        }
104
105        let count = self.count + other.count;
106        let total = count as f64;
107        let mean = self.mean * (self.count as f64 / total) + other.mean * (other.count as f64 / total);
108        let apart = other.mean - self.mean;
109        let weight = (self.count as f64 * other.count as f64) / total;
110
111        Self {
112            count,
113            sum: self.sum + other.sum,
114            mean,
115            lowest: self.lowest.min(other.lowest),
116            highest: self.highest.max(other.highest),
117            deviations: self.deviations + other.deviations + apart * apart * weight,
118        }
119    }
120
121    #[must_use]
122    pub fn count(&self) -> u64 {
123        self.count
124    }
125
126    #[must_use]
127    pub fn is_empty(&self) -> bool {
128        self.count == 0
129    }
130
131    #[must_use]
132    pub fn sum(&self) -> f64 {
133        self.sum
134    }
135
136    #[must_use]
137    pub fn min(&self) -> Option<f64> {
138        (self.count > 0).then_some(self.lowest)
139    }
140
141    #[must_use]
142    pub fn max(&self) -> Option<f64> {
143        (self.count > 0).then_some(self.highest)
144    }
145
146    #[must_use]
147    pub fn span(&self) -> Option<f64> {
148        (self.count > 0).then_some(self.highest - self.lowest)
149    }
150
151    #[must_use]
152    pub fn mean(&self) -> Option<f64> {
153        (self.count > 0).then_some(self.mean)
154    }
155
156    // Population variance: the spread of exactly these values, not an
157    // estimate of a wider population's.
158    #[must_use]
159    pub fn variance(&self) -> Option<f64> {
160        (self.count > 0).then(|| (self.deviations / self.count as f64).max(0.0))
161    }
162
163    #[must_use]
164    pub fn std(&self) -> Option<f64> {
165        self.variance().map(f64::sqrt)
166    }
167
168    // The mean square is the squared mean plus the variance; `hypot` takes the
169    // root of that sum without squaring either term into overflow.
170    #[must_use]
171    pub fn rms(&self) -> Option<f64> {
172        Some(self.mean()?.hypot(self.std()?))
173    }
174}
175
176/* === Public Functions === */
177
178/* === Private Functions === */
179
180/* === Tests === */
181
182#[cfg(test)]
183mod tests {
184    use super::Summary;
185    use crate::signal::Series;
186
187    const CLOSE: f64 = 1e-9;
188
189    fn close(found: Option<f64>, expected: f64) {
190        let found = found.expect("an answer");
191        assert!((found - expected).abs() <= CLOSE * expected.abs().max(1.0), "{found} is not {expected}");
192    }
193
194    // Series answers over the whole run at once; Summary is folded from chunks and
195    // merged, so agreement proves the streaming form matches the exact one.
196    #[test]
197    fn chunked_summaries_agree() {
198        let values: Vec<f64> = (0..97).map(|index| ((index * 37) % 23) as f64 * 1.5 - 11.0).chain([f64::NAN, f64::INFINITY]).collect();
199        let series = Series::new(&values);
200        let summary = values.chunks(10).map(Summary::of).fold(Summary::new(), |held, next| held.merge(&next));
201
202        assert_eq!(series.count(), summary.count());
203        close(Some(series.sum()), summary.sum());
204        close(series.mean(), summary.mean().unwrap());
205        assert_eq!(series.min(), summary.min());
206        assert_eq!(series.max(), summary.max());
207        close(series.rms(), summary.rms().unwrap());
208        close(series.std(), summary.std().unwrap());
209    }
210}