sva_samples/measure/
bands.rs

1// Concern: traces one envelope per ERB band, beside that band's own impulse floor | Non-concern: the discrete spectrum (spectrum.rs), what a band means to a listener | IO: (&[f64], sample rate) -> bands
2
3use sva_formula::filter::Shape;
4
5use crate::biquad::{Coeffs, State, design};
6
7/// Glasberg and Moore's auditory filter width in Hz.
8pub fn erb_hz(hz: f64) -> f64 {
9    24.7 * (4.37 * hz / 1000.0 + 1.0)
10}
11
12/// One Cam is one filter width along the cochlea.
13pub fn cam(hz: f64) -> f64 {
14    21.4 * (1.0 + 4.37 * hz / 1000.0).log10()
15}
16
17pub fn hz_at_cam(cam: f64) -> f64 {
18    (10f64.powf(cam / 21.4) - 1.0) * 1000.0 / 4.37
19}
20
21/// 1 to 40 Cam is 51 Hz to 16.4 kHz.
22pub const FIRST_CAM: f64 = 1.0;
23pub const BAND_COUNT: usize = 40;
24
25/// Under the 2-3 ms gap-detection limit, a hundredth of the samples it came from.
26pub const DECIMATED_HZ: f64 = 2000.0;
27
28/// Sixteen time constants is past every band's peak.
29fn floor_secs(erb: f64) -> f64 {
30    (16.0 / erb).clamp(0.02, 2.0)
31}
32
33/// A filter cannot report an onset faster than its own, so a caller subtracts this.
34#[derive(Clone, Copy, Debug, PartialEq)]
35pub struct BandFloor {
36    pub peak: f64,
37    pub time_to_peak_secs: f64,
38    pub rise_10_90_secs: Option<f64>,
39}
40
41#[derive(Clone, Debug, PartialEq)]
42pub struct BandTrack {
43    pub centre_hz: f64,
44    pub cam: f64,
45    pub erb_hz: f64,
46    pub q: f64,
47    pub peak: f64,
48    pub time_to_peak_secs: Option<f64>,
49    pub rise_10_90_secs: Option<f64>,
50    pub floor: BandFloor,
51    pub rms: Vec<f64>,
52}
53
54#[derive(Clone, Debug, PartialEq)]
55pub struct Bands {
56    pub rate_hz: f64,
57    pub start_secs: f64,
58    pub bands: Vec<BandTrack>,
59}
60
61/// `hz: None` is the band's own ERB, capped at a quarter of the decimated rate.
62#[derive(Clone, Copy, Debug, PartialEq)]
63pub struct Smoothing {
64    pub hz: Option<f64>,
65    pub poles: usize,
66}
67
68const OWN_WIDTH: Smoothing = Smoothing { hz: None, poles: 1 };
69
70/// Squaring rectifies; the root at decimation makes the track an RMS.
71struct Chain {
72    pass: Coeffs,
73    smooth: Coeffs,
74    a: State,
75    b: State,
76    energy: Vec<State>,
77}
78
79impl Chain {
80    fn new(centre_hz: f64, sr: f64, smoothing: Smoothing) -> Chain {
81        let erb = erb_hz(centre_hz);
82        let hz = smoothing.hz.unwrap_or(erb.min(DECIMATED_HZ / 4.0));
83        Chain {
84            pass: design(Shape::Bandpass, centre_hz, centre_hz / erb, 0.0, sr),
85            smooth: design(Shape::OnePole, hz, 0.0, 0.0, sr),
86            a: State::default(),
87            b: State::default(),
88            energy: vec![State::default(); smoothing.poles],
89        }
90    }
91
92    #[inline]
93    fn step(&mut self, x: f64) -> f64 {
94        let banded = self.b.step(&self.pass, self.a.step(&self.pass, x));
95        self.energy
96            .iter_mut()
97            .fold(banded * banded, |e, pole| pole.step(&self.smooth, e))
98    }
99}
100
101fn trace(
102    centre_hz: f64,
103    sr: f64,
104    smoothing: Smoothing,
105    hop: usize,
106    n: usize,
107    at: impl Fn(usize) -> f64,
108) -> Vec<f64> {
109    let mut chain = Chain::new(centre_hz, sr, smoothing);
110    let mut out = Vec::with_capacity(n / hop + 1);
111    for i in 0..n {
112        let energy = chain.step(at(i));
113        if i % hop == 0 {
114            out.push(energy.max(0.0).sqrt());
115        }
116    }
117    out
118}
119
120/// Interpolated between the frames straddling it.
121fn crossing(track: &[f64], until: usize, level: f64, rate: f64) -> Option<f64> {
122    let hit = track[..=until].iter().position(|&v| v >= level)?;
123    if hit == 0 {
124        return Some(0.0);
125    }
126    let (lo, hi) = (track[hit - 1], track[hit]);
127    let frac = match hi > lo {
128        true => (level - lo) / (hi - lo),
129        false => 0.0,
130    };
131    Some((hit as f64 - 1.0 + frac) / rate)
132}
133
134fn stats(track: &[f64], rate: f64) -> (f64, Option<f64>, Option<f64>) {
135    let peak = track.iter().copied().fold(0.0, f64::max);
136    if peak <= 0.0 {
137        return (peak, None, None);
138    }
139    let at = track
140        .iter()
141        .position(|&v| v == peak)
142        .expect("the peak is one of the frames");
143    let t10 = crossing(track, at, 0.1 * peak, rate);
144    let t90 = crossing(track, at, 0.9 * peak, rate);
145    let rise = t10.zip(t90).map(|(a, b)| (b - a).max(0.0));
146    (peak, Some(at as f64 / rate), rise)
147}
148
149fn floor_of(centre_hz: f64, sr: f64, smoothing: Smoothing, hop: usize, rate: f64) -> BandFloor {
150    let n = (floor_secs(erb_hz(centre_hz)) * sr).round().max(4.0) as usize;
151    let track = trace(centre_hz, sr, smoothing, hop, n, |i| {
152        f64::from(u8::from(i == 0))
153    });
154    let (peak, at, rise) = stats(&track, rate);
155    BandFloor {
156        peak,
157        time_to_peak_secs: at.unwrap_or(0.0),
158        rise_10_90_secs: rise,
159    }
160}
161
162/// A centre the design cannot place is left out, never clamped onto the last it could.
163fn centres(sr: f64) -> Vec<f64> {
164    (0..BAND_COUNT)
165        .map(|k| hz_at_cam(FIRST_CAM + k as f64))
166        .filter(|hz| *hz < 0.45 * sr)
167        .collect()
168}
169
170pub fn analyze(samples: &[f64], sr: f64, start_secs: f64) -> Bands {
171    smoothed(samples, sr, start_secs, OWN_WIDTH)
172}
173
174pub fn smoothed(samples: &[f64], sr: f64, start_secs: f64, smoothing: Smoothing) -> Bands {
175    let hop = (sr / DECIMATED_HZ).round().max(1.0) as usize;
176    let rate = sr / hop as f64;
177    Bands {
178        rate_hz: rate,
179        start_secs,
180        bands: centres(sr)
181            .into_iter()
182            .map(|centre_hz| {
183                let track = trace(centre_hz, sr, smoothing, hop, samples.len(), |i| samples[i]);
184                let (peak, at, rise) = stats(&track, rate);
185                let erb = erb_hz(centre_hz);
186                BandTrack {
187                    centre_hz,
188                    cam: cam(centre_hz),
189                    erb_hz: erb,
190                    q: centre_hz / erb,
191                    peak,
192                    time_to_peak_secs: at,
193                    rise_10_90_secs: rise,
194                    floor: floor_of(centre_hz, sr, smoothing, hop, rate),
195                    rms: track,
196                }
197            })
198            .collect(),
199    }
200}
sva_samples/measure/bands.rs

sva_samples/measure/
bands.rs
