sound_to_patch/dsp/
frames.rs

1//! Active frames and the robust statistics over them.
2
3use super::StftParams;
4use super::stft::padded_sample;
5use crate::cast;
6use crate::units::Decibels;
7
8/// The gate of the project's convention: a frame is active when its RMS
9/// lies within 40 dB of the loudest frame.
10pub const ACTIVE_GATE_DB: Decibels = Decibels::from_const(40.0);
11
12/// The fewest active frames a statistic over frames needs; with fewer it
13/// gives no estimate.
14///
15/// The first and last two frames of a clip reach into the zero padding, so
16/// a median over few frames follows the edges. On 48 kHz noise, sawtooth
17/// and sine clips with 4096-point frames, the worst error of the centroid,
18/// spread, rolloff and flatness medians against a 10 s run was 4 % at 10
19/// frames, 3 % at 12, 2 % at 16 and up to 47 times the value at 8 or fewer.
20pub const MIN_ACTIVE_FRAMES: usize = 16;
21
22/// Returns the RMS of each centred frame of an STFT layout, over the raw
23/// samples (no window) with the layout's padding. Frame `m` lines up with
24/// STFT frame `m`.
25#[must_use]
26pub fn frame_rms(samples: &[f64], params: StftParams) -> Vec<f64> {
27    let len = cast::usize_to_f64(params.frame_len());
28    (0..params.frame_count(samples.len()))
29        .map(|frame| {
30            let start = params.frame_start(frame);
31            let energy: f64 = (0..params.frame_len())
32                .map(|offset| {
33                    let index = start.saturating_add(isize::try_from(offset).unwrap_or(isize::MAX));
34                    padded_sample(samples, index, params.padding()).powi(2)
35                })
36                .sum();
37            (energy / len).sqrt()
38        })
39        .collect()
40}
41
42/// The frames that carry signal: frame RMS within a gate of the loudest
43/// frame. Digital silence has no active frame.
44#[derive(Clone, Debug, PartialEq, Eq)]
45pub struct ActiveFrames {
46    active: Vec<bool>,
47}
48
49impl ActiveFrames {
50    /// Marks each frame whose RMS lies within `gate` decibels of the largest
51    /// frame RMS. A frame of zero RMS is never active.
52    #[must_use]
53    pub fn from_rms(rms: &[f64], gate: Decibels) -> Self {
54        let peak = rms.iter().copied().fold(0.0, f64::max);
55        let threshold = peak
56            * Decibels::clamped(-gate.get().abs())
57                .to_amplitude_gain()
58                .unwrap_or(0.0);
59        Self {
60            active: rms
61                .iter()
62                .map(|&value| value > 0.0 && value >= threshold)
63                .collect(),
64        }
65    }
66
67    /// The project's convention: frame RMS of an STFT layout, gated at
68    /// [`ACTIVE_GATE_DB`].
69    #[must_use]
70    pub fn of(samples: &[f64], params: StftParams) -> Self {
71        Self::from_rms(&frame_rms(samples, params), ACTIVE_GATE_DB)
72    }
73
74    /// Returns whether frame `m` is active; a frame past the end is not.
75    #[must_use]
76    pub fn is_active(&self, frame: usize) -> bool {
77        self.active.get(frame).copied().unwrap_or(false)
78    }
79
80    /// Returns the indices of the active frames in order.
81    pub fn indices(&self) -> impl Iterator<Item = usize> + '_ {
82        self.active
83            .iter()
84            .enumerate()
85            .filter_map(|(frame, &active)| active.then_some(frame))
86    }
87
88    /// Returns the number of active frames.
89    #[must_use]
90    pub fn count(&self) -> usize {
91        self.active.iter().filter(|&&active| active).count()
92    }
93
94    /// Returns the number of frames, active or not.
95    #[must_use]
96    pub fn len(&self) -> usize {
97        self.active.len()
98    }
99
100    /// Returns whether no frame was given.
101    #[must_use]
102    pub fn is_empty(&self) -> bool {
103        self.active.is_empty()
104    }
105}
106
107/// The median and the interquartile range of a set of values.
108#[derive(Clone, Copy, Debug, PartialEq)]
109pub struct Summary {
110    median: f64,
111    iqr: f64,
112}
113
114impl Summary {
115    /// Summarises finite values. Quantiles use linear interpolation between
116    /// the sorted values, at position `q (n - 1)`, the convention of the
117    /// common numeric libraries. Returns `None` for no value or a NaN or
118    /// infinite one.
119    #[must_use]
120    pub fn of(values: &[f64]) -> Option<Self> {
121        if values.is_empty() || values.iter().any(|v| !v.is_finite()) {
122            return None;
123        }
124        let mut sorted = values.to_vec();
125        sorted.sort_by(f64::total_cmp);
126        Some(Self {
127            median: quantile(&sorted, 0.5)?,
128            iqr: quantile(&sorted, 0.75)? - quantile(&sorted, 0.25)?,
129        })
130    }
131
132    /// Returns the median.
133    #[must_use]
134    pub fn median(self) -> f64 {
135        self.median
136    }
137
138    /// Returns the interquartile range, the 75th minus the 25th percentile.
139    #[must_use]
140    pub fn iqr(self) -> f64 {
141        self.iqr
142    }
143}
144
145impl Summary {
146    /// Carries the median and the IQR into a unit type.
147    #[must_use]
148    pub fn to_stats<T>(self, unit: impl Fn(f64) -> T) -> Stats<T> {
149        Stats {
150            median: unit(self.median),
151            iqr: unit(self.iqr),
152        }
153    }
154}
155
156/// The median and the interquartile range of a value over active frames, in
157/// its unit type.
158#[derive(Clone, Copy, Debug, PartialEq)]
159pub struct Stats<T> {
160    median: T,
161    iqr: T,
162}
163
164impl<T: Copy> Stats<T> {
165    /// Returns the median.
166    #[must_use]
167    pub fn median(&self) -> T {
168        self.median
169    }
170
171    /// Returns the interquartile range.
172    #[must_use]
173    pub fn iqr(&self) -> T {
174        self.iqr
175    }
176}
177
178/// The linear quantile of sorted values, the default of `numpy.quantile`.
179/// Returns `None` for no values.
180#[must_use]
181pub fn quantile(sorted: &[f64], q: f64) -> Option<f64> {
182    let last = sorted.len().checked_sub(1)?;
183    let position = q * cast::usize_to_f64(last);
184    let below = cast::floor_to_usize(position)?.min(last);
185    let above = (below + 1).min(last);
186    let fraction = position - cast::usize_to_f64(below);
187    Some(sorted[below] + (sorted[above] - sorted[below]) * fraction)
188}
189
190#[cfg(test)]
191mod tests {
192    use super::*;
193
194    #[test]
195    fn the_gate_constant_passes_the_unit_check() {
196        assert_eq!(ACTIVE_GATE_DB, Decibels::new(40.0).unwrap());
197    }
198
199    #[test]
200    fn quantiles_interpolate_linearly() {
201        let summary = Summary::of(&[4.0, 1.0, 3.0, 2.0]).unwrap();
202        // Sorted 1 2 3 4: median at 1.5, quartiles at 0.75 and 2.25.
203        assert_eq!(summary.median(), 2.5);
204        assert_eq!(summary.iqr(), 3.25 - 1.75);
205        let single = Summary::of(&[7.0]).unwrap();
206        assert_eq!((single.median(), single.iqr()), (7.0, 0.0));
207        assert!(Summary::of(&[]).is_none());
208        assert!(Summary::of(&[1.0, f64::NAN]).is_none());
209    }
210
211    #[test]
212    fn frames_within_the_gate_are_active() {
213        let rms = [1.0, 0.011, 0.009, 0.0, 0.5];
214        let active = ActiveFrames::from_rms(&rms, Decibels::new(40.0).unwrap());
215        assert_eq!(active.indices().collect::<Vec<_>>(), vec![0, 1, 4]);
216        assert_eq!(active.count(), 3);
217        assert_eq!(active.len(), 5);
218        assert!(!active.is_active(9));
219    }
220
221    #[test]
222    fn silence_has_no_active_frame() {
223        let params = StftParams::new(64).unwrap();
224        let active = ActiveFrames::of(&[0.0; 1000], params);
225        assert_eq!(active.count(), 0);
226        assert_eq!(active.len(), params.frame_count(1000));
227    }
228
229    #[test]
230    fn frame_rms_of_a_constant_counts_the_padding() {
231        let params = StftParams::new(8).unwrap();
232        let rms = frame_rms(&[1.0; 16], params);
233        // Frame 0 covers 4 padded zeros and 4 ones.
234        assert!((rms[0] - 0.5_f64.sqrt()).abs() < 1e-15);
235        assert!((rms[4] - 1.0).abs() < 1e-15);
236    }
237}
sound_to_patch/dsp/frames.rs

sound_to_patch/dsp/
frames.rs
