Skip to main content

sound_to_patch/notes/
quiet.rs

1//! The quiet-note guard: a note that sounds only where the stem is far
2//! below its loudest passage is bleed or noise, not part of the stem.
3//!
4//! The model scales each 2 s window by its own level, so a near-silent
5//! window with bleed fills its full range and gives notes. The guard
6//! drops every note whose span holds no active frame: no frame within the
7//! gate of the loudest frame of the stem.
8
9use std::ops::Range;
10
11use crate::audio::SampleRate;
12use crate::cast;
13use crate::dsp::{ActiveFrames, StftParams, Window, frame_rms};
14use crate::units::{Decibels, Hz, MidiPitch, Seconds};
15
16use super::NoteList;
17
18/// The length of a guard frame: 25 ms, short enough to follow a note.
19pub const QUIET_FRAME: Seconds = Seconds::from_const(0.025);
20
21/// The gate of the guard: a frame is active within this many decibels of
22/// the loudest frame. It is the project's active-frame gate.
23pub const QUIET_GATE: Decibels = crate::dsp::ACTIVE_GATE_DB;
24
25/// The stem's level in [`QUIET_FRAME`] frames at its own rate, with a
26/// quarter-frame hop, frame `m` centred at `m` hops. The quiet-note guard,
27/// the attack split and the split join read the same frames.
28#[derive(Clone, Debug)]
29pub struct StemLevels<'a> {
30    mono: &'a [f64],
31    rate: SampleRate,
32    params: StftParams,
33    rms: Vec<f64>,
34    window: Vec<f64>,
35}
36
37impl<'a> StemLevels<'a> {
38    /// Measures the RMS of every frame of `mono` at `rate`.
39    #[must_use]
40    pub fn new(mono: &'a [f64], rate: SampleRate) -> Self {
41        let frame = cast::round_to_usize(QUIET_FRAME.get() * f64::from(rate.get()))
42            .unwrap_or(4)
43            .max(4);
44        let params = StftParams::fixed(frame);
45        Self {
46            mono,
47            rate,
48            params,
49            rms: frame_rms(mono, params),
50            window: Window::Hann.samples(params.frame_len()),
51        }
52    }
53
54    /// Returns the RMS of each frame.
55    #[must_use]
56    pub fn rms(&self) -> &[f64] {
57        &self.rms
58    }
59
60    /// Returns the RMS of frame `frame` about its own mean: the level with
61    /// any DC offset removed.
62    pub(crate) fn ac_rms(&self, frame: usize) -> f64 {
63        self.ac_level(frame * self.params.hop(), self.params.frame_len())
64    }
65
66    /// Returns the time at which frame `frame` is centred.
67    pub(crate) fn frame_time(&self, frame: usize) -> Seconds {
68        Seconds::clamped(cast::usize_to_f64(frame * self.params.hop()) / f64::from(self.rate.get()))
69    }
70
71    fn hops(&self, time: Seconds) -> f64 {
72        time.get() * f64::from(self.rate.get()) / cast::usize_to_f64(self.params.hop())
73    }
74
75    /// Returns sample `offset` of the `len` samples centred at sample
76    /// `centre`. Samples outside the stem count as silence.
77    fn sample(&self, centre: usize, offset: usize, len: usize) -> f64 {
78        (centre + offset)
79            .checked_sub(len / 2)
80            .and_then(|index| self.mono.get(index))
81            .copied()
82            .unwrap_or(0.0)
83    }
84
85    /// Returns the RMS about its own mean of the `len` samples centred at
86    /// sample `centre`.
87    fn ac_level(&self, centre: usize, len: usize) -> f64 {
88        let (sum, energy) = (0..len)
89            .map(|offset| self.sample(centre, offset, len))
90            .fold((0.0, 0.0), |(sum, energy), x| (sum + x, energy + x * x));
91        let count = cast::usize_to_f64(len);
92        let mean = sum / count;
93        (energy / count - mean * mean).max(0.0).sqrt()
94    }
95
96    /// Returns the RMS about its own mean of the `window` of samples
97    /// centred at `time`: the envelope with any DC offset removed. Samples
98    /// outside the stem count as silence.
99    #[must_use]
100    pub fn envelope(&self, time: Seconds, window: Seconds) -> f64 {
101        let rate = f64::from(self.rate.get());
102        let len = cast::round_to_usize(window.get() * rate)
103            .unwrap_or(1)
104            .max(1);
105        self.ac_level(cast::round_to_usize(time.get() * rate).unwrap_or(0), len)
106    }
107
108    /// Returns the frames centred from `from` to `to`, both included.
109    #[must_use]
110    pub fn frames_between(&self, from: Seconds, to: Seconds) -> Range<usize> {
111        let first = cast::ceil_to_usize(self.hops(from)).unwrap_or(0);
112        let end = cast::floor_to_usize(self.hops(to))
113            .map_or(0, |last| last + 1)
114            .min(self.rms.len());
115        first.min(end)..end
116    }
117
118    /// Returns the level of frame `frame` at the first [`HARMONICS`]
119    /// harmonics of `pitch` below the Nyquist frequency: the root of the
120    /// summed power of each, from a Hann window over the frame with the
121    /// frame's mean removed. Samples outside the stem count as silence.
122    ///
123    /// The window's main lobe spans 160 Hz, so below about A2 (110 Hz) it
124    /// cannot tell a note's fundamental from a tone a fifth or a tritone
125    /// away.
126    #[must_use]
127    pub fn harmonic_level(&self, frame: usize, pitch: MidiPitch) -> f64 {
128        let len = self.params.frame_len();
129        let rate = f64::from(self.rate.get());
130        let centre = frame * self.params.hop();
131        let raw: Vec<f64> = (0..len)
132            .map(|offset| self.sample(centre, offset, len))
133            .collect();
134        let mean = raw.iter().sum::<f64>() / cast::usize_to_f64(len);
135        let windowed: Vec<f64> = raw
136            .iter()
137            .zip(&self.window)
138            .map(|(x, weight)| (x - mean) * weight)
139            .collect();
140        let f0 = Hz::from(pitch).get();
141        (1..=HARMONICS)
142            .map(|h| f0 * f64::from(h))
143            .take_while(|&frequency| frequency < rate / 2.0)
144            .map(|frequency| bin_power(&windowed, frequency, rate))
145            .sum::<f64>()
146            .sqrt()
147    }
148
149    /// Returns the level at the fundamental of `pitch` of the one period of
150    /// samples centred at `time`, with their mean removed: one DFT bin over
151    /// the period, unnormalised, so only ratios of two readings mean
152    /// anything. The bin's zeros fall on every harmonic of the note when
153    /// the period is a whole number of samples. It is rounded to one, so at
154    /// 22,050 Hz and at the worst phase a harmonic leaks about -70 dB at
155    /// A1, -46 to -50 dB from A2 to A5, -23 dB at A6 and -17 dB at A7. The
156    /// bin rejects only the note's harmonics: a key 1 to 3 semitones up
157    /// reads within 1 dB of the note, and a fifth up 2 to 6 dB under it.
158    /// Samples outside the stem count as silence.
159    #[must_use]
160    pub fn fundamental_level(&self, time: Seconds, pitch: MidiPitch) -> f64 {
161        let rate = f64::from(self.rate.get());
162        let f0 = Hz::from(pitch).get();
163        let len = cast::round_to_usize(rate / f0).unwrap_or(1).max(1);
164        let centre = cast::round_to_usize(time.get() * rate).unwrap_or(0);
165        let raw: Vec<f64> = (0..len)
166            .map(|offset| self.sample(centre, offset, len))
167            .collect();
168        let mean = raw.iter().sum::<f64>() / cast::usize_to_f64(len);
169        let centred: Vec<f64> = raw.iter().map(|x| x - mean).collect();
170        bin_power(&centred, f0, rate).sqrt()
171    }
172}
173
174/// The power of one DFT bin of `samples` at `frequency`, at `rate`, by the
175/// Goertzel recurrence, which reads any frequency.
176fn bin_power(samples: &[f64], frequency: f64, rate: f64) -> f64 {
177    let coefficient = 2.0 * (std::f64::consts::TAU * frequency / rate).cos();
178    let (mut s1, mut s2) = (0.0, 0.0);
179    for x in samples {
180        let s0 = x + coefficient * s1 - s2;
181        s2 = s1;
182        s1 = s0;
183    }
184    (s1 * s1 + s2 * s2 - coefficient * s1 * s2).max(0.0)
185}
186
187/// The harmonics of a note that [`StemLevels::harmonic_level`] sums.
188pub const HARMONICS: u8 = 8;
189
190/// Drops every note whose span holds no active frame of the stem, and
191/// returns the kept notes with the count of dropped ones.
192///
193/// A note's frames run from the frame nearest its onset, which can be
194/// centred up to half a hop before it, to the last frame centred before
195/// its offset; a note shorter than a hop takes the frame nearest its
196/// onset.
197#[must_use]
198pub fn drop_quiet(
199    notes: &NoteList,
200    mono: &[f64],
201    rate: SampleRate,
202    gate: Decibels,
203) -> (NoteList, usize) {
204    drop_quiet_with(notes, &StemLevels::new(mono, rate), gate)
205}
206
207/// Drops the notes as [`drop_quiet`] does, over levels already measured.
208#[must_use]
209pub fn drop_quiet_with(
210    notes: &NoteList,
211    levels: &StemLevels<'_>,
212    gate: Decibels,
213) -> (NoteList, usize) {
214    let active = ActiveFrames::from_rms(levels.rms(), gate);
215    let before = notes.notes().len();
216    let kept: Vec<_> = notes
217        .notes()
218        .iter()
219        .filter(|note| {
220            let first = cast::round_to_usize(levels.hops(note.onset())).unwrap_or(0);
221            let last = cast::ceil_to_usize(levels.hops(note.offset()))
222                .unwrap_or(first)
223                .max(first + 1);
224            (first..last).any(|m| active.is_active(m))
225        })
226        .cloned()
227        .collect();
228    let dropped = before - kept.len();
229    (NoteList::new(kept), dropped)
230}
231
232#[cfg(test)]
233mod tests {
234    use super::*;
235    use crate::notes::Note;
236    use crate::units::MidiPitch;
237
238    fn note(onset: f64, offset: f64) -> Note {
239        Note::new(
240            Seconds::new(onset).unwrap(),
241            Seconds::new(offset).unwrap(),
242            MidiPitch::new(60.0).unwrap(),
243        )
244        .unwrap()
245    }
246
247    #[test]
248    fn a_note_in_a_passage_50_db_down_is_dropped_and_one_at_30_db_kept() {
249        let rate = SampleRate::new(8000).unwrap();
250        // 1 s at full scale, 1 s at -50 dB, 1 s at -30 dB.
251        let mono: Vec<f64> = (0..24_000)
252            .map(|n| {
253                let level = match n / 8000 {
254                    0 => 1.0,
255                    1 => 10_f64.powf(-50.0 / 20.0),
256                    _ => 10_f64.powf(-30.0 / 20.0),
257                };
258                let t = f64::from(n) / 8000.0;
259                level * (std::f64::consts::TAU * 220.0 * t).sin()
260            })
261            .collect();
262        let notes = NoteList::new(vec![
263            note(0.2, 0.8),
264            note(1.2, 1.8),
265            note(2.2, 2.8),
266            note(0.9, 1.5),
267        ]);
268        let (kept, dropped) = drop_quiet(&notes, &mono, rate, QUIET_GATE);
269        assert_eq!(dropped, 1);
270        let onsets: Vec<f64> = kept.notes().iter().map(|n| n.onset().get()).collect();
271        assert_eq!(onsets, vec![0.2, 0.9, 2.2]);
272    }
273
274    /// 1 s of silence, then 1 s of a 220 Hz sine at full scale, at 8 kHz:
275    /// frames of 200 samples, a hop of 6.25 ms.
276    fn silence_then_tone() -> Vec<f64> {
277        (0..16_000)
278            .map(|n| {
279                let t = f64::from(n) / 8000.0;
280                if n < 8000 {
281                    0.0
282                } else {
283                    (std::f64::consts::TAU * 220.0 * t).sin()
284                }
285            })
286            .collect()
287    }
288
289    #[test]
290    fn a_note_that_reaches_into_a_loud_passage_is_kept() {
291        let rate = SampleRate::new(8000).unwrap();
292        let notes = NoteList::new(vec![note(0.5, 1.03), note(0.5, 0.9)]);
293        let (kept, dropped) = drop_quiet(&notes, &silence_then_tone(), rate, QUIET_GATE);
294        assert_eq!(dropped, 1);
295        assert!((kept.notes()[0].offset().get() - 1.03).abs() < 1e-9);
296    }
297
298    #[test]
299    fn a_note_in_a_10_ms_gap_is_kept_by_the_25_ms_frames_over_it() {
300        let rate = SampleRate::new(8000).unwrap();
301        // The tone, 10 ms of silence from 1 s, then the tone again.
302        let tone = silence_then_tone().split_off(8000);
303        let mut mono = tone.clone();
304        mono.extend(std::iter::repeat_n(0.0, 80));
305        mono.extend(&tone);
306        let notes = NoteList::new(vec![note(1.001, 1.008)]);
307        let (kept, dropped) = drop_quiet(&notes, &mono, rate, QUIET_GATE);
308        assert_eq!((kept.notes().len(), dropped), (1, 0));
309    }
310
311    #[test]
312    fn a_note_shorter_than_a_hop_takes_the_frame_at_its_onset() {
313        let rate = SampleRate::new(8000).unwrap();
314        // From 0.6 to 0.9 of hop 240, inside the tone: the onset rounds up
315        // to frame 241, and the offset reaches no further.
316        let hop = 0.006_25;
317        let notes = NoteList::new(vec![note(240.6 * hop, 240.9 * hop)]);
318        let (kept, dropped) = drop_quiet(&notes, &silence_then_tone(), rate, QUIET_GATE);
319        assert_eq!((kept.notes().len(), dropped), (1, 0));
320    }
321
322    /// A sine of amplitude 1 at `frequency` from sample `start` of 1 s at
323    /// `rate`.
324    fn sine_from(rate: u32, frequency: f64, start: usize) -> Vec<f64> {
325        (0..rate)
326            .map(|n| {
327                let t = f64::from(n) / f64::from(rate);
328                if usize::try_from(n).unwrap() < start {
329                    0.0
330                } else {
331                    (std::f64::consts::TAU * frequency * t).sin()
332                }
333            })
334            .collect()
335    }
336
337    #[test]
338    fn frames_between_holds_both_ends() {
339        let rate = SampleRate::new(8000).unwrap();
340        let mono = vec![0.0; 8000];
341        let levels = StemLevels::new(&mono, rate);
342        let s = |value| Seconds::new(value).unwrap();
343        assert_eq!(levels.frames_between(s(0.0), s(0.0125)), 0..3);
344        assert_eq!(levels.frames_between(s(0.001), s(0.0125)), 1..3);
345        assert_eq!(
346            levels.frames_between(s(0.99), s(2.0)).end,
347            levels.rms().len()
348        );
349    }
350
351    #[test]
352    fn the_harmonic_level_reads_each_harmonic_in_its_own_frame() {
353        // At 8 kHz a frame holds 200 samples: a full sine of amplitude 1
354        // reads half the Hann window's sum, 50.
355        let rate = SampleRate::new(8000).unwrap();
356        let a4 = MidiPitch::new(69.0).unwrap();
357        for frequency in [440.0, 1320.0] {
358            let mono = sine_from(8000, frequency, 4000);
359            let levels = StemLevels::new(&mono, rate);
360            // Frame 80 is centred on sample 4000, where the sine starts.
361            assert!(levels.harmonic_level(78, a4).abs() < 1e-12);
362            let level = levels.harmonic_level(82, a4);
363            assert!((level - 50.0).abs() < 0.5, "{frequency} Hz reads {level}");
364        }
365    }
366
367    #[test]
368    fn the_envelope_reads_the_rms_of_its_window_about_its_mean() {
369        // At 8 kHz 10 ms hold 80 samples, 4.4 periods of 440 Hz: close to
370        // the RMS of a sine, 0.707, with or without a DC offset.
371        let rate = SampleRate::new(8000).unwrap();
372        let s = |value| Seconds::new(value).unwrap();
373        let mono = sine_from(8000, 440.0, 4000);
374        let offset: Vec<f64> = mono.iter().map(|x| x + 0.5).collect();
375        for stem in [&mono, &offset] {
376            let levels = StemLevels::new(stem, rate);
377            let full = levels.envelope(s(0.75), s(0.01));
378            assert!(
379                (full - std::f64::consts::FRAC_1_SQRT_2).abs() < 0.02,
380                "{full}"
381            );
382        }
383        // Centred on the start of the sine, half the window is silence.
384        let half = StemLevels::new(&mono, rate).envelope(s(0.5), s(0.01));
385        assert!((half - 0.5).abs() < 0.03, "{half}");
386    }
387
388    #[test]
389    fn a_dc_offset_leaves_the_harmonic_level_alone() {
390        let rate = SampleRate::new(8000).unwrap();
391        let a1 = MidiPitch::new(33.0).unwrap();
392        let mono = sine_from(8000, 55.0, 0);
393        let offset: Vec<f64> = mono.iter().map(|x| x + 0.5).collect();
394        let clean = StemLevels::new(&mono, rate).harmonic_level(80, a1);
395        let shifted = StemLevels::new(&offset, rate).harmonic_level(80, a1);
396        assert!(
397            (shifted / clean - 1.0).abs() < 0.05,
398            "{shifted} against {clean}"
399        );
400    }
401
402    #[test]
403    fn below_a2_the_harmonic_level_hears_a_fifth_as_the_note() {
404        // At A1 a tone a fifth up, at 82.4 Hz, reads within 3 dB of A1
405        // itself: the window cannot resolve them.
406        let rate = SampleRate::new(8000).unwrap();
407        let a1 = MidiPitch::new(33.0).unwrap();
408        let level = |frequency| {
409            let mono = sine_from(8000, frequency, 0);
410            StemLevels::new(&mono, rate).harmonic_level(80, a1)
411        };
412        let ratio = level(82.407) / level(55.0);
413        assert!(
414            ratio > 10_f64.powf(-3.0 / 20.0),
415            "the fifth reads {ratio} of A1"
416        );
417    }
418
419    #[test]
420    fn the_harmonic_level_stops_below_the_nyquist_frequency() {
421        // At 7040 Hz the eighth harmonic of A4 lies on the Nyquist frequency.
422        let nyquist: Vec<f64> = (0..7040)
423            .map(|n| if n % 2 == 0 { 1.0 } else { -1.0 })
424            .collect();
425        let levels = StemLevels::new(&nyquist, SampleRate::new(7040).unwrap());
426        assert!(levels.harmonic_level(80, MidiPitch::new(69.0).unwrap()) < 1.0);
427        // At 8 kHz the fifth harmonic of 1250 Hz lies above it, and the
428        // sixth would alias to 1750 Hz.
429        let pitch = MidiPitch::new(69.0 + 12.0 * (1250.0_f64 / 440.0).log2()).unwrap();
430        let mono = sine_from(8000, 1750.0, 0);
431        let levels = StemLevels::new(&mono, SampleRate::new(8000).unwrap());
432        assert!(levels.harmonic_level(80, pitch) < 1.0);
433    }
434
435    #[test]
436    fn silence_drops_every_note() {
437        let rate = SampleRate::new(8000).unwrap();
438        let (kept, dropped) = drop_quiet(
439            &NoteList::new(vec![note(0.1, 0.2)]),
440            &vec![0.0; 8000],
441            rate,
442            QUIET_GATE,
443        );
444        assert!(kept.is_empty());
445        assert_eq!(dropped, 1);
446    }
447}