diff --git a/Cargo.toml b/Cargo.toml index b765d36..465bafd 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -15,6 +15,10 @@ path = "src/main.rs" name = "test_net" path = "src/bin/test_net.rs" +[[bin]] +name = "specview" +path = "src/bin/specview.rs" + [dependencies] anyhow = "1.0.102" async-trait = "0.1.89" diff --git a/src/bin/specview.rs b/src/bin/specview.rs new file mode 100644 index 0000000..685661f --- /dev/null +++ b/src/bin/specview.rs @@ -0,0 +1,258 @@ +//! `specview` — a developer CLI for the `dsp` toolkit: conjure test signals, +//! draw spectrograms in the terminal, and measure echo-cancellation quality +//! (ERLE) on WAV files. Not shipped in the GUI; a measurement aid. +//! +//! Usage: +//! specview show [--width N] [--height N] [--max-hz HZ] [--floor DB] [--ascii] +//! specview gen -o out.wav [--secs S] [--amp A] [--freq F] [--f0 F] [--f1 F] [--seed N] +//! kinds: sine | sweep | white | pink | impulse | silence +//! specview erle [--show] +//! +//! Examples: +//! specview gen sweep -o sweep.wav --secs 3 --f0 100 --f1 12000 +//! specview show ~/peerspeak-recordings/2026-06-14/peerspeak-*.wav +//! specview erle echo_only.wav after_aec.wav --show + +use std::collections::HashMap; +use std::path::Path; +use std::process::ExitCode; + +use peerspeak::dsp::{generators, metrics, render, stft, wav}; + +const SAMPLE_RATE: u32 = 48_000; + +fn main() -> ExitCode { + let args: Vec = std::env::args().skip(1).collect(); + let Some(cmd) = args.first() else { + eprintln!("{}", USAGE); + return ExitCode::FAILURE; + }; + + let rest = &args[1..]; + let result = match cmd.as_str() { + "show" => cmd_show(rest), + "gen" => cmd_gen(rest), + "erle" => cmd_erle(rest), + "-h" | "--help" | "help" => { + println!("{}", USAGE); + Ok(()) + } + other => Err(format!("unknown command {other:?}\n\n{USAGE}")), + }; + + match result { + Ok(()) => ExitCode::SUCCESS, + Err(e) => { + eprintln!("error: {e}"); + ExitCode::FAILURE + } + } +} + +const USAGE: &str = "specview — spectrogram + AEC measurement tool\n\n\ + specview show [--width N] [--height N] [--max-hz HZ] [--floor DB] [--ascii]\n\ + specview gen -o out.wav [--secs S] [--amp A] [--freq F] [--f0 F] [--f1 F] [--seed N]\n\ + kinds: sine sweep white pink impulse silence\n\ + specview erle [--show]"; + +/// Renders a WAV as a spectrogram plus a numeric summary. +fn cmd_show(args: &[String]) -> Result<(), String> { + let (positional, flags) = parse_args(args); + let path = positional.first().ok_or("show needs a ")?; + let data = wav::read(Path::new(path))?; + + let opts = render::RenderOpts { + width: flags.usize_or("width", 100), + height: flags.usize_or("height", 30), + max_hz: flags.f32_or("max-hz", 12_000.0), + floor_db: flags.f32_or("floor", -80.0), + ascii: flags.present("ascii"), + }; + let spec = stft::analyze(&data.samples, data.sample_rate, 2048, 512); + print!("{}", render::render(&spec, &opts)); + print_summary(path, &data.samples, &spec); + Ok(()) +} + +/// Generates a test signal and writes it to a WAV file. +fn cmd_gen(args: &[String]) -> Result<(), String> { + let (positional, flags) = parse_args(args); + let kind = positional.first().ok_or("gen needs a ")?; + let out = flags + .get("o") + .or_else(|| flags.get("out")) + .ok_or("gen needs -o ")?; + + let secs = flags.f32_or("secs", 2.0); + let amp = flags.f32_or("amp", 0.8); + let len = (secs * SAMPLE_RATE as f32) as usize; + let seed = flags.u64_or("seed", 1); + + let samples = match kind.as_str() { + "sine" => generators::sine(flags.f32_or("freq", 1000.0), amp, len, SAMPLE_RATE), + "sweep" => generators::log_sweep( + flags.f32_or("f0", 100.0), + flags.f32_or("f1", 12_000.0), + amp, + len, + SAMPLE_RATE, + ), + "white" => generators::white_noise(amp, len, seed), + "pink" => generators::pink_noise(amp, len, seed), + "impulse" => generators::impulse(amp, len), + "silence" => generators::silence(len), + other => return Err(format!("unknown kind {other:?} (sine sweep white pink impulse silence)")), + }; + + wav::write(Path::new(out), &samples, SAMPLE_RATE)?; + println!( + "wrote {} — {} ({:.2}s, {} samples, rms {:.1} dBFS)", + out, + kind, + secs, + samples.len(), + metrics::dbfs(metrics::rms(&samples)), + ); + Ok(()) +} + +/// Measures ERLE between a "before" (echo present) and "after" (post-AEC) file — +/// broadband and per voice band — so you can see how much echo was removed and +/// in which frequency range any residual lives. +fn cmd_erle(args: &[String]) -> Result<(), String> { + let (positional, flags) = parse_args(args); + let before = positional.first().ok_or("erle needs ")?; + let after = positional.get(1).ok_or("erle needs ")?; + + let b = wav::read(Path::new(before))?; + let a = wav::read(Path::new(after))?; + // Compare over the overlapping length. + let n = b.samples.len().min(a.samples.len()); + let bs = &b.samples[..n]; + let as_ = &a.samples[..n]; + + let broadband = metrics::erle(bs, as_); + println!("ERLE (broadband): {broadband:+.1} dB (higher = more echo removed)"); + println!( + " before: {:.1} dBFS rms after: {:.1} dBFS rms", + metrics::dbfs(metrics::rms(bs)), + metrics::dbfs(metrics::rms(as_)), + ); + + // Per-band ERLE via averaged STFT magnitude energy. + let sb = stft::analyze(bs, b.sample_rate, 2048, 512); + let sa = stft::analyze(as_, a.sample_rate, 2048, 512); + let avg_b = average_frame(&sb); + let avg_a = average_frame(&sa); + println!(" per band:"); + for band in metrics::VOICE_BANDS { + let eb = metrics::band_energy(&avg_b, |bin| sb.bin_hz(bin), band.low_hz, band.high_hz); + let ea = metrics::band_energy(&avg_a, |bin| sa.bin_hz(bin), band.low_hz, band.high_hz); + let db = if eb > 1e-12 && ea > 1e-12 { + 10.0 * (eb / ea).log10() + } else if eb > 1e-12 { + 120.0 + } else { + 0.0 + }; + println!(" {:<18} {:+6.1} dB", band.label, db); + } + + if flags.present("show") { + println!("\n--- before ---"); + print!("{}", render::render(&sb, &render::RenderOpts::default())); + println!("\n--- after ---"); + print!("{}", render::render(&sa, &render::RenderOpts::default())); + } + Ok(()) +} + +/// Prints RMS/peak/dBFS plus the per-band energy distribution of a signal. +fn print_summary(path: &str, samples: &[f32], spec: &stft::Spectrogram) { + println!( + "\n{}\n duration {:.2}s rms {:.1} dBFS peak {:.1} dBFS", + path, + samples.len() as f32 / spec.sample_rate as f32, + metrics::dbfs(metrics::rms(samples)), + metrics::dbfs(metrics::peak(samples)), + ); + let avg = average_frame(spec); + let total: f32 = avg.iter().map(|m| m * m).sum::().max(1e-12); + println!(" energy by band:"); + for band in metrics::VOICE_BANDS { + let e = metrics::band_energy(&avg, |bin| spec.bin_hz(bin), band.low_hz, band.high_hz); + let pct = 100.0 * e / total; + let bar = "█".repeat((pct / 2.5).round() as usize); + println!(" {:<18} {:5.1}% {}", band.label, pct, bar); + } +} + +/// Mean magnitude spectrum across all STFT frames (the time-averaged spectrum). +fn average_frame(spec: &stft::Spectrogram) -> Vec { + let bins = spec.bins(); + let mut avg = vec![0.0f32; bins]; + if spec.frames.is_empty() { + return avg; + } + for frame in &spec.frames { + for (bin, a) in avg.iter_mut().enumerate() { + *a += frame[bin]; + } + } + let n = spec.frames.len() as f32; + for a in &mut avg { + *a /= n; + } + avg +} + +/// A trivial `--flag value` / `--flag` (boolean) parser. Positional args are +/// everything not consumed as a flag value. Keeps the binary dependency-free. +struct Flags { + map: HashMap, + bools: Vec, +} + +impl Flags { + fn get(&self, key: &str) -> Option<&String> { + self.map.get(key) + } + fn present(&self, key: &str) -> bool { + self.bools.iter().any(|b| b == key) || self.map.contains_key(key) + } + fn f32_or(&self, key: &str, default: f32) -> f32 { + self.map.get(key).and_then(|v| v.parse().ok()).unwrap_or(default) + } + fn usize_or(&self, key: &str, default: usize) -> usize { + self.map.get(key).and_then(|v| v.parse().ok()).unwrap_or(default) + } + fn u64_or(&self, key: &str, default: u64) -> u64 { + self.map.get(key).and_then(|v| v.parse().ok()).unwrap_or(default) + } +} + +/// Splits args into positionals and flags. A token like `--width`/`-o` followed +/// by a non-flag token captures that token as its value; otherwise it's a bool. +fn parse_args(args: &[String]) -> (Vec, Flags) { + let mut positional = Vec::new(); + let mut map = HashMap::new(); + let mut bools = Vec::new(); + let mut i = 0; + while i < args.len() { + let tok = &args[i]; + if let Some(key) = tok.strip_prefix("--").or_else(|| tok.strip_prefix('-')) { + let next_is_value = args.get(i + 1).is_some_and(|n| !n.starts_with('-')); + if next_is_value { + map.insert(key.to_string(), args[i + 1].clone()); + i += 2; + } else { + bools.push(key.to_string()); + i += 1; + } + } else { + positional.push(tok.clone()); + i += 1; + } + } + (positional, Flags { map, bools }) +} diff --git a/src/dsp/fft.rs b/src/dsp/fft.rs new file mode 100644 index 0000000..7ec15ff --- /dev/null +++ b/src/dsp/fft.rs @@ -0,0 +1,186 @@ +//! A small, dependency-free radix-2 Cooley-Tukey FFT. +//! +//! We hand-roll this rather than pull in `rustfft`/`realfft` because (a) the +//! whole DSP toolkit is a developer/measurement aid, not a hot real-time path, +//! and (b) it keeps the supply-chain surface at zero new crates. The transform +//! is the textbook iterative in-place algorithm (bit-reversal permutation + +//! log2(N) butterfly stages), computed in `f64` for headroom even though the +//! audio it analyses is `f32`. +//! +//! Only power-of-two lengths are supported; the STFT layer always pads frames +//! up to a power of two before calling in. + +use std::f64::consts::PI; + +/// A minimal complex number for the transform. Kept local (rather than pulling a +/// `num-complex` dependency) since the FFT is the only thing that needs it. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct Complex { + pub re: f64, + pub im: f64, +} + +impl Complex { + pub const ZERO: Complex = Complex { re: 0.0, im: 0.0 }; + + pub fn new(re: f64, im: f64) -> Self { + Complex { re, im } + } + + /// Magnitude `sqrt(re^2 + im^2)`. + pub fn magnitude(self) -> f64 { + self.re.hypot(self.im) + } + + fn add(self, o: Complex) -> Complex { + Complex::new(self.re + o.re, self.im + o.im) + } + + fn sub(self, o: Complex) -> Complex { + Complex::new(self.re - o.re, self.im - o.im) + } + + fn mul(self, o: Complex) -> Complex { + Complex::new( + self.re * o.re - self.im * o.im, + self.re * o.im + self.im * o.re, + ) + } +} + +/// In-place forward FFT. `buf.len()` must be a power of two. +/// +/// Uses the standard sign convention `X[k] = sum_n x[n] * exp(-2πi·kn/N)`. +pub fn fft(buf: &mut [Complex]) { + transform(buf, false); +} + +/// In-place inverse FFT (normalized by `1/N`), the exact inverse of [`fft`]. +pub fn ifft(buf: &mut [Complex]) { + transform(buf, true); + let n = buf.len() as f64; + for c in buf.iter_mut() { + c.re /= n; + c.im /= n; + } +} + +fn transform(buf: &mut [Complex], inverse: bool) { + let n = buf.len(); + assert!(n.is_power_of_two(), "FFT length {n} must be a power of two"); + if n <= 1 { + return; + } + + // Bit-reversal permutation: reorder so the iterative butterflies can run + // bottom-up in place. + let mut j = 0usize; + for i in 1..n { + let mut bit = n >> 1; + while j & bit != 0 { + j ^= bit; + bit >>= 1; + } + j ^= bit; + if i < j { + buf.swap(i, j); + } + } + + // Butterfly stages: combine length-`len` DFTs from length-`len/2` halves, + // doubling `len` each pass. + let sign = if inverse { 1.0 } else { -1.0 }; + let mut len = 2; + while len <= n { + let ang = sign * 2.0 * PI / len as f64; + let wlen = Complex::new(ang.cos(), ang.sin()); + let mut i = 0; + while i < n { + let mut w = Complex::new(1.0, 0.0); + for k in 0..len / 2 { + let u = buf[i + k]; + let v = buf[i + k + len / 2].mul(w); + buf[i + k] = u.add(v); + buf[i + k + len / 2] = u.sub(v); + w = w.mul(wlen); + } + i += len; + } + len <<= 1; + } +} + +/// Forward FFT of a real signal, returning the **one-sided** magnitude spectrum: +/// bins `0..=N/2` (DC through Nyquist), where `N` is the next power of two ≥ +/// `samples.len()`. The input is zero-padded up to `N`. +/// +/// Magnitudes are raw (un-normalized) linear amplitudes; callers convert to dB +/// or normalize as needed. +pub fn real_magnitude_spectrum(samples: &[f32]) -> Vec { + let n = samples.len().next_power_of_two().max(2); + let mut buf = vec![Complex::ZERO; n]; + for (i, &s) in samples.iter().enumerate() { + buf[i].re = s as f64; + } + fft(&mut buf); + buf[..=n / 2].iter().map(|c| c.magnitude() as f32).collect() +} + +#[cfg(test)] +mod tests { + use super::*; + + fn approx(a: f64, b: f64, eps: f64) -> bool { + (a - b).abs() <= eps + } + + #[test] + fn impulse_transforms_to_flat_spectrum() { + // FFT of a unit impulse at n=0 is all-ones (flat spectrum). + let mut buf = vec![Complex::ZERO; 8]; + buf[0] = Complex::new(1.0, 0.0); + fft(&mut buf); + for c in &buf { + assert!(approx(c.magnitude(), 1.0, 1e-9), "expected flat 1.0, got {c:?}"); + } + } + + #[test] + fn single_bin_sine_peaks_in_that_bin() { + // A cosine at exactly bin k=2 over N=16 should put all energy in bin 2 + // (and its mirror N-2). Check the one-sided spectrum peaks at bin 2. + let n = 16; + let k = 2; + let samples: Vec = (0..n) + .map(|i| (2.0 * PI * k as f64 * i as f64 / n as f64).cos() as f32) + .collect(); + let mag = real_magnitude_spectrum(&samples); + let peak_bin = mag + .iter() + .enumerate() + .max_by(|a, b| a.1.partial_cmp(b.1).unwrap()) + .unwrap() + .0; + assert_eq!(peak_bin, k, "energy should land in bin {k}, got {peak_bin}"); + } + + #[test] + fn ifft_inverts_fft() { + let original: Vec = (0..32) + .map(|i| Complex::new((i as f64 * 0.3).sin(), (i as f64 * 0.1).cos())) + .collect(); + let mut buf = original.clone(); + fft(&mut buf); + ifft(&mut buf); + for (a, b) in original.iter().zip(&buf) { + assert!(approx(a.re, b.re, 1e-9) && approx(a.im, b.im, 1e-9)); + } + } + + #[test] + #[should_panic(expected = "power of two")] + fn non_power_of_two_panics() { + let mut buf = vec![Complex::ZERO; 6]; + fft(&mut buf); + } +} diff --git a/src/dsp/generators.rs b/src/dsp/generators.rs new file mode 100644 index 0000000..b6fa57b --- /dev/null +++ b/src/dsp/generators.rs @@ -0,0 +1,156 @@ +//! Reproducible test-signal generators — the "conjured noises" the measurement +//! harness drives AEC/analysis with. Everything is deterministic: the noise +//! sources take an explicit seed so a run is byte-for-byte repeatable, which is +//! what makes A/B comparisons (old build vs new, tuning X vs Y) meaningful. +//! +//! All generators emit `f32` samples nominally in `[-amp, amp]` at the given +//! sample rate. Convert to `i16` at the I/O edge. + +use std::f64::consts::PI; + +/// A tiny deterministic PRNG (xorshift64*) so the noise generators need no +/// external crate and produce identical output for a given seed across runs. +struct XorShift64 { + state: u64, +} + +impl XorShift64 { + fn new(seed: u64) -> Self { + // Avoid the all-zero fixed point. + XorShift64 { + state: seed.wrapping_mul(0x9E37_79B9_7F4A_7C15) | 1, + } + } + + fn next_u64(&mut self) -> u64 { + let mut x = self.state; + x ^= x >> 12; + x ^= x << 25; + x ^= x >> 27; + self.state = x; + x.wrapping_mul(0x2545_F491_4F6C_DD1D) + } + + /// Uniform `f32` in `[-1.0, 1.0)`. + fn next_bipolar(&mut self) -> f32 { + // Top 24 bits -> [0,1), then map to [-1,1). + let u = (self.next_u64() >> 40) as f32 / (1u32 << 24) as f32; + u * 2.0 - 1.0 + } +} + +/// `len` samples of silence. +pub fn silence(len: usize) -> Vec { + vec![0.0; len] +} + +/// A pure sine tone at `freq` Hz. +pub fn sine(freq: f32, amp: f32, len: usize, sample_rate: u32) -> Vec { + (0..len) + .map(|i| amp * (2.0 * PI * freq as f64 * i as f64 / sample_rate as f64).sin() as f32) + .collect() +} + +/// A unit impulse: `amp` at sample 0, silence after. Its flat spectrum makes it +/// the natural probe for an echo path's impulse response. +pub fn impulse(amp: f32, len: usize) -> Vec { + let mut v = vec![0.0; len]; + if len > 0 { + v[0] = amp; + } + v +} + +/// A logarithmic (exponential) sine sweep from `f0` to `f1` Hz over the buffer — +/// the standard excitation for measuring a system's frequency response, since it +/// spends equal time per octave. +pub fn log_sweep(f0: f32, f1: f32, amp: f32, len: usize, sample_rate: u32) -> Vec { + if len == 0 { + return Vec::new(); + } + let t_total = len as f64 / sample_rate as f64; + let (f0, f1) = (f0 as f64, f1 as f64); + let k = (f1 / f0).ln(); + (0..len) + .map(|i| { + let t = i as f64 / sample_rate as f64; + // Instantaneous-phase integral of an exponential chirp. + let phase = 2.0 * PI * f0 * t_total / k * ((k * t / t_total).exp() - 1.0); + amp * phase.sin() as f32 + }) + .collect() +} + +/// Uniform white noise in `[-amp, amp]`, seeded for reproducibility. +pub fn white_noise(amp: f32, len: usize, seed: u64) -> Vec { + let mut rng = XorShift64::new(seed); + (0..len).map(|_| amp * rng.next_bipolar()).collect() +} + +/// Pink (1/f) noise via the Voss-McCartney algorithm — perceptually flatter than +/// white noise and a closer stand-in for room/voice spectra. Seeded. +pub fn pink_noise(amp: f32, len: usize, seed: u64) -> Vec { + let mut rng = XorShift64::new(seed); + const ROWS: usize = 16; + let mut rows = [0.0f32; ROWS]; + let mut running = 0.0f32; + let mut out = Vec::with_capacity(len); + for i in 0..len { + // Each row updates half as often as the previous: the row to refresh is + // the index of the lowest set bit of the sample counter. + let n = i + 1; + let row = (n & n.wrapping_neg()).trailing_zeros() as usize; + if row < ROWS { + running -= rows[row]; + rows[row] = rng.next_bipolar(); + running += rows[row]; + } + // Normalize: sum of ROWS unit-bipolar rows spans roughly [-ROWS, ROWS]. + out.push(amp * running / ROWS as f32); + } + out +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::dsp::metrics::rms; + + #[test] + fn sine_has_expected_length_and_range() { + let s = sine(1000.0, 0.5, 480, 48_000); + assert_eq!(s.len(), 480); + assert!(s.iter().all(|&x| x.abs() <= 0.5001)); + } + + #[test] + fn impulse_is_single_spike() { + let s = impulse(1.0, 10); + assert_eq!(s[0], 1.0); + assert!(s[1..].iter().all(|&x| x == 0.0)); + } + + #[test] + fn white_noise_is_deterministic_for_a_seed() { + let a = white_noise(1.0, 256, 42); + let b = white_noise(1.0, 256, 42); + let c = white_noise(1.0, 256, 43); + assert_eq!(a, b, "same seed -> identical output"); + assert_ne!(a, c, "different seed -> different output"); + } + + #[test] + fn noise_energy_is_nonzero_and_bounded() { + let w = white_noise(0.5, 4096, 7); + let r = rms(&w); + assert!(r > 0.0 && r < 0.5, "rms {r} should be in (0, amp)"); + } + + #[test] + fn log_sweep_spans_the_buffer() { + let s = log_sweep(100.0, 8000.0, 1.0, 4800, 48_000); + assert_eq!(s.len(), 4800); + // Non-trivial energy present. + assert!(rms(&s) > 0.1); + } +} diff --git a/src/dsp/metrics.rs b/src/dsp/metrics.rs new file mode 100644 index 0000000..8285566 --- /dev/null +++ b/src/dsp/metrics.rs @@ -0,0 +1,132 @@ +//! Objective signal metrics — the numbers that turn "sounds better" into a +//! measurement. The headline one for echo cancellation is **ERLE** (Echo Return +//! Loss Enhancement): how much echo energy the canceller removed, in dB. + +/// Root-mean-square level of a signal (linear amplitude). `0.0` for empty input. +pub fn rms(samples: &[f32]) -> f32 { + if samples.is_empty() { + return 0.0; + } + let sum_sq: f64 = samples.iter().map(|&s| s as f64 * s as f64).sum(); + (sum_sq / samples.len() as f64).sqrt() as f32 +} + +/// Peak absolute amplitude. `0.0` for empty input. +pub fn peak(samples: &[f32]) -> f32 { + samples.iter().fold(0.0, |m, &s| m.max(s.abs())) +} + +/// Converts a linear amplitude (e.g. an RMS or peak value, relative to full-scale +/// `1.0`) to decibels below full scale. A floor of -120 dBFS is returned for +/// silence so the result is always finite. +pub fn dbfs(linear: f32) -> f32 { + if linear <= 1e-6 { + return -120.0; + } + 20.0 * linear.log10() +} + +/// **Echo Return Loss Enhancement**, in dB: `10·log10(E[echo²] / E[residual²])`. +/// +/// `echo` is the signal *before* cancellation (the echo the mic picked up); +/// `residual` is what's *left after* the canceller ran. A larger number is +/// better — e.g. +30 dB means the canceller removed 99.9% of the echo energy. +/// Returns a +120 dB ceiling if the residual is effectively silent (perfect +/// cancellation) and 0.0 if there was no echo energy to begin with. +pub fn erle(echo: &[f32], residual: &[f32]) -> f32 { + let e_echo = mean_square(echo); + let e_res = mean_square(residual); + if e_echo <= 1e-12 { + return 0.0; + } + if e_res <= 1e-12 { + return 120.0; + } + 10.0 * (e_echo / e_res).log10() as f32 +} + +/// Mean square (average energy per sample) of a signal. `0.0` for empty input. +fn mean_square(samples: &[f32]) -> f64 { + if samples.is_empty() { + return 0.0; + } + samples.iter().map(|&s| s as f64 * s as f64).sum::() / samples.len() as f64 +} + +/// A frequency band for per-band energy analysis, in Hz. +pub struct Band { + pub label: &'static str, + pub low_hz: f32, + pub high_hz: f32, +} + +/// Voice-relevant bands for spotting *where* residual echo or noise lives. +pub const VOICE_BANDS: &[Band] = &[ + Band { label: "low (80-300)", low_hz: 80.0, high_hz: 300.0 }, + Band { label: "low-mid (300-1k)", low_hz: 300.0, high_hz: 1000.0 }, + Band { label: "mid (1k-3k)", low_hz: 1000.0, high_hz: 3000.0 }, + Band { label: "high-mid (3k-6k)", low_hz: 3000.0, high_hz: 6000.0 }, + Band { label: "high (6k-12k)", low_hz: 6000.0, high_hz: 12000.0 }, +]; + +/// Sums the linear magnitude energy within `[low_hz, high_hz)` across a single +/// STFT magnitude frame. `bin_hz` maps a bin index to its centre frequency. +pub fn band_energy(frame: &[f32], bin_hz: impl Fn(usize) -> f32, low_hz: f32, high_hz: f32) -> f32 { + frame + .iter() + .enumerate() + .filter(|&(bin, _)| { + let hz = bin_hz(bin); + hz >= low_hz && hz < high_hz + }) + .map(|(_, &m)| m * m) + .sum() +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn rms_of_constant_is_that_constant() { + assert!((rms(&[0.5; 100]) - 0.5).abs() < 1e-6); + assert_eq!(rms(&[]), 0.0); + } + + #[test] + fn peak_finds_max_magnitude() { + assert_eq!(peak(&[0.1, -0.9, 0.3]), 0.9); + } + + #[test] + fn dbfs_landmarks() { + assert!((dbfs(1.0) - 0.0).abs() < 1e-4, "full scale = 0 dBFS"); + assert!((dbfs(0.5) - -6.0206).abs() < 1e-3, "half = ~-6 dB"); + assert_eq!(dbfs(0.0), -120.0, "silence floors"); + } + + #[test] + fn erle_halving_energy_is_about_3db() { + // residual amplitude = echo/sqrt(2) -> half the energy -> ~3.01 dB. + let echo = vec![1.0f32; 1000]; + let residual = vec![std::f32::consts::FRAC_1_SQRT_2; 1000]; + let e = erle(&echo, &residual); + assert!((e - 3.0103).abs() < 0.01, "expected ~3 dB, got {e}"); + } + + #[test] + fn erle_perfect_cancellation_ceils() { + assert_eq!(erle(&[1.0; 10], &[0.0; 10]), 120.0); + assert_eq!(erle(&[0.0; 10], &[0.0; 10]), 0.0, "no echo -> 0"); + } + + #[test] + fn band_energy_selects_the_right_bins() { + // 5-bin frame, 100 Hz per bin: bins at 0,100,200,300,400 Hz. + let frame = [1.0, 2.0, 3.0, 4.0, 5.0]; + let bin_hz = |b: usize| b as f32 * 100.0; + // [150,350) -> bins 200,300 Hz -> 3^2 + 4^2 = 25. + let e = band_energy(&frame, bin_hz, 150.0, 350.0); + assert!((e - 25.0).abs() < 1e-4, "got {e}"); + } +} diff --git a/src/dsp/mod.rs b/src/dsp/mod.rs new file mode 100644 index 0000000..96f8ae6 --- /dev/null +++ b/src/dsp/mod.rs @@ -0,0 +1,18 @@ +//! Signal-analysis toolkit: FFT, STFT spectrograms, reproducible test-signal +//! generators, objective metrics (RMS/dBFS/ERLE/per-band energy), a minimal WAV +//! reader/writer, and a 24-bit-colour terminal spectrogram renderer. +//! +//! This is a **measurement/developer aid**, not part of the real-time audio +//! path. It exists so audio work — especially echo-cancellation tuning — can be +//! evaluated with reproducible numbers and a visible spectrogram instead of +//! ear alone. The pure layers (`fft`, `window`, `stft`, `generators`, `metrics`, +//! `render`) have no I/O; `wav` is the only edge that touches the filesystem. +//! Driven from the `specview` binary (`src/bin/specview.rs`). + +pub mod fft; +pub mod generators; +pub mod metrics; +pub mod render; +pub mod stft; +pub mod wav; +pub mod window; diff --git a/src/dsp/render.rs b/src/dsp/render.rs new file mode 100644 index 0000000..161c921 --- /dev/null +++ b/src/dsp/render.rs @@ -0,0 +1,260 @@ +//! Renders a [`Spectrogram`] to a string for the terminal. +//! +//! The trick that gives real resolution in a text grid: each character cell is +//! the Unicode upper-half block `▀` with a 24-bit-colour **foreground** (top +//! pixel) and **background** (bottom pixel), so one row of text shows two +//! frequency bins. Magnitudes are mapped to dB against the loudest cell and +//! coloured with a magma-style ramp (dark = quiet, bright = loud). An ASCII +//! fallback ramp is available for non-truecolor / piped output. + +use super::stft::Spectrogram; + +/// Rendering options for [`render`]. +pub struct RenderOpts { + /// Output width in character columns (time axis). + pub width: usize, + /// Output height in character rows; vertical resolution is `2×height` bins. + pub height: usize, + /// Dynamic-range floor in dB below the loudest cell (e.g. `-80.0`). + pub floor_db: f32, + /// Top of the frequency axis in Hz (clamped to Nyquist). + pub max_hz: f32, + /// Use a plain-ASCII brightness ramp instead of 24-bit colour blocks. + pub ascii: bool, +} + +impl Default for RenderOpts { + fn default() -> Self { + RenderOpts { + width: 100, + height: 30, + floor_db: -80.0, + max_hz: 12_000.0, + ascii: false, + } + } +} + +const RESET: &str = "\x1b[0m"; +/// Brightness ramp (dark→bright) for the ASCII fallback. +const ASCII_RAMP: &[u8] = b" .:-=+*#%@"; + +/// Renders `spec` to a multi-line string: a header, the frequency-labelled +/// spectrogram body, a time axis, and a dB colour legend. +pub fn render(spec: &Spectrogram, opts: &RenderOpts) -> String { + if spec.frames.is_empty() { + return "(empty signal — nothing to plot)\n".to_string(); + } + + let nyquist = spec.sample_rate as f32 / 2.0; + let max_hz = opts.max_hz.min(nyquist).max(spec.bin_hz(1)); + // Number of bins from DC up to max_hz (at least 2 so we can aggregate). + let bins_used = (0..spec.bins()) + .take_while(|&b| spec.bin_hz(b) <= max_hz) + .count() + .max(2); + + let n_frames = spec.frames.len(); + let width = opts.width.clamp(1, n_frames); + let rows = opts.height.max(1) * 2; // 2 frequency pixels per character row + + // Stage 1: aggregate frames into `width` time columns (mean magnitude/bin). + let mut columns: Vec> = Vec::with_capacity(width); + for c in 0..width { + let f_lo = c * n_frames / width; + let f_hi = ((c + 1) * n_frames / width).max(f_lo + 1).min(n_frames); + let span = (f_hi - f_lo) as f32; + let mut acc = vec![0.0f32; bins_used]; + for f in f_lo..f_hi { + let frame = &spec.frames[f]; + for (bin, a) in acc.iter_mut().enumerate() { + *a += frame[bin]; + } + } + for a in &mut acc { + *a /= span; + } + columns.push(acc); + } + + // Stage 2: aggregate each column's bins into `rows` vertical pixels. + // grid[col][pixel], pixel 0 = lowest freq. + let mut grid: Vec> = Vec::with_capacity(width); + let mut max_mag = 1e-9f32; + for col in &columns { + let mut pixels = vec![0.0f32; rows]; + for (p, px) in pixels.iter_mut().enumerate() { + let b_lo = p * bins_used / rows; + let b_hi = ((p + 1) * bins_used / rows).max(b_lo + 1).min(bins_used); + let mean = col[b_lo..b_hi].iter().copied().sum::() / (b_hi - b_lo) as f32; + *px = mean; + max_mag = max_mag.max(mean); + } + grid.push(pixels); + } + + // Stage 3: draw. Char row r (0 = top) shows pixels (top=high freq, bottom). + let mut out = String::new(); + out.push_str(&format!( + "Spectrogram {:.0} Hz–{:.0} Hz · {:.2}s · {} frames, fft {}, hop {}\n", + 0.0, + max_hz, + spec.frame_time(n_frames), + n_frames, + spec.fft_size, + spec.hop, + )); + + let label_every = (opts.height.max(1) / 6).max(1); + for r in 0..opts.height { + let top_pixel = rows - 1 - 2 * r; + let bot_pixel = rows.saturating_sub(2 + 2 * r); + // Frequency label for this char row (use the top pixel's centre freq). + let label = if r % label_every == 0 { + let frac = top_pixel as f32 / (rows - 1) as f32; + format!("{:>6.0}", frac * max_hz) + } else { + " ".to_string() + }; + out.push_str(&label); + out.push_str(" \u{2502}"); // " │" + for col in &grid { + let top = norm_db(col[top_pixel], max_mag, opts.floor_db); + let bot = norm_db(col[bot_pixel], max_mag, opts.floor_db); + push_cell(&mut out, top, bot, opts.ascii); + } + if !opts.ascii { + out.push_str(RESET); + } + out.push('\n'); + } + + // Time axis. + out.push_str(" \u{2514}"); // " └" + out.push_str(&"\u{2500}".repeat(width)); + out.push('\n'); + out.push_str(&format!( + " 0.00s{:>width$.2}s\n", + spec.frame_time(n_frames), + width = width.saturating_sub(5), + )); + + out.push_str(&legend(opts)); + out +} + +/// Maps a linear magnitude to `[0,1]` over `[floor_db, 0]` relative to `max_mag`. +fn norm_db(mag: f32, max_mag: f32, floor_db: f32) -> f32 { + if mag <= 1e-9 { + return 0.0; + } + let db = 20.0 * (mag / max_mag).log10(); + ((db - floor_db) / -floor_db).clamp(0.0, 1.0) +} + +/// Appends one character cell (two stacked pixels) to `out`. +fn push_cell(out: &mut String, top: f32, bot: f32, ascii: bool) { + if ascii { + // One pixel per cell in ASCII mode (use the brighter of the pair). + let v = top.max(bot); + let idx = ((v * (ASCII_RAMP.len() - 1) as f32).round() as usize).min(ASCII_RAMP.len() - 1); + out.push(ASCII_RAMP[idx] as char); + } else { + let (tr, tg, tb) = magma(top); + let (br, bg, bb) = magma(bot); + // fg = top pixel, bg = bottom pixel, glyph = upper half block. + out.push_str(&format!( + "\x1b[38;2;{tr};{tg};{tb}m\x1b[48;2;{br};{bg};{bb}m\u{2580}" + )); + } +} + +/// A magma-style colormap: `t` in `[0,1]` → `(r,g,b)`. Piecewise-linear through +/// control colours sampled from matplotlib's `magma`. +fn magma(t: f32) -> (u8, u8, u8) { + const STOPS: &[(f32, (f32, f32, f32))] = &[ + (0.0, (0.0, 0.0, 4.0)), + (0.25, (80.0, 18.0, 123.0)), + (0.5, (182.0, 54.0, 121.0)), + (0.75, (252.0, 136.0, 97.0)), + (1.0, (252.0, 253.0, 191.0)), + ]; + let t = t.clamp(0.0, 1.0); + for win in STOPS.windows(2) { + let (t0, c0) = win[0]; + let (t1, c1) = win[1]; + if t <= t1 { + let f = if t1 > t0 { (t - t0) / (t1 - t0) } else { 0.0 }; + return ( + (c0.0 + (c1.0 - c0.0) * f) as u8, + (c0.1 + (c1.1 - c0.1) * f) as u8, + (c0.2 + (c1.2 - c0.2) * f) as u8, + ); + } + } + (252, 253, 191) +} + +/// A horizontal dB colour legend strip. +fn legend(opts: &RenderOpts) -> String { + let mut s = format!(" dB: {:.0}", opts.floor_db); + let steps = 32; + if opts.ascii { + for i in 0..steps { + let v = i as f32 / (steps - 1) as f32; + let idx = ((v * (ASCII_RAMP.len() - 1) as f32).round() as usize).min(ASCII_RAMP.len() - 1); + s.push(ASCII_RAMP[idx] as char); + } + } else { + s.push(' '); + for i in 0..steps { + let v = i as f32 / (steps - 1) as f32; + let (r, g, b) = magma(v); + s.push_str(&format!("\x1b[48;2;{r};{g};{b}m ")); + } + s.push_str(RESET); + } + s.push_str(" 0\n"); + s +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::dsp::{generators, stft}; + + #[test] + fn norm_db_maps_floor_and_peak() { + // At max magnitude -> 0 dB -> top of range (1.0). + assert!((norm_db(1.0, 1.0, -80.0) - 1.0).abs() < 1e-5); + // Silence -> 0.0. + assert_eq!(norm_db(0.0, 1.0, -80.0), 0.0); + // -80 dB (mag 1e-4 of max) -> ~0.0 bottom of range. + assert!(norm_db(1e-4, 1.0, -80.0) < 0.01); + } + + #[test] + fn magma_endpoints() { + assert_eq!(magma(0.0), (0, 0, 4)); + assert_eq!(magma(1.0), (252, 253, 191)); + } + + #[test] + fn render_produces_grid_of_expected_height() { + let sig = generators::sine(2000.0, 0.8, 48_000, 48_000); + let spec = stft::analyze(&sig, 48_000, 1024, 512); + let opts = RenderOpts { width: 40, height: 10, ..Default::default() }; + let out = render(&spec, &opts); + // Header + 10 body rows + time axis (2) + legend = non-trivial. + let lines = out.lines().count(); + assert!(lines >= 13, "expected >=13 lines, got {lines}"); + assert!(out.contains("Spectrogram")); + } + + #[test] + fn empty_spectrogram_renders_message() { + let spec = stft::analyze(&[], 48_000, 512, 256); + let out = render(&spec, &RenderOpts::default()); + assert!(out.contains("empty")); + } +} diff --git a/src/dsp/stft.rs b/src/dsp/stft.rs new file mode 100644 index 0000000..edf778c --- /dev/null +++ b/src/dsp/stft.rs @@ -0,0 +1,105 @@ +//! Short-Time Fourier Transform: slice a signal into overlapping windowed +//! frames and take each frame's magnitude spectrum. The result is the +//! time×frequency matrix a spectrogram draws. + +use super::fft::real_magnitude_spectrum; +use super::window::{apply, hann}; + +/// One STFT analysis: a sequence of magnitude frames plus the geometry needed to +/// label axes (sample rate, fft size, hop). +pub struct Spectrogram { + /// `frames[t][bin]` = linear magnitude of frequency `bin` at time-step `t`. + /// Each inner vec has `fft_size/2 + 1` bins (DC..Nyquist). + pub frames: Vec>, + pub sample_rate: u32, + pub fft_size: usize, + pub hop: usize, +} + +impl Spectrogram { + /// Number of one-sided frequency bins per frame (`fft_size/2 + 1`). + pub fn bins(&self) -> usize { + self.fft_size / 2 + 1 + } + + /// Centre frequency (Hz) of bin index `bin`. + pub fn bin_hz(&self, bin: usize) -> f32 { + bin as f32 * self.sample_rate as f32 / self.fft_size as f32 + } + + /// Time (seconds) at the start of frame `t`. + pub fn frame_time(&self, t: usize) -> f32 { + (t * self.hop) as f32 / self.sample_rate as f32 + } +} + +/// Computes the STFT of `samples` with the given `fft_size` (rounded up to a +/// power of two) and `hop` (frame advance in samples). A Hann window is applied +/// to each frame. The final partial frame is zero-padded so trailing audio is +/// not dropped. +pub fn analyze(samples: &[f32], sample_rate: u32, fft_size: usize, hop: usize) -> Spectrogram { + let fft_size = fft_size.next_power_of_two().max(2); + let hop = hop.max(1); + let window = hann(fft_size); + + let mut frames = Vec::new(); + if !samples.is_empty() { + let mut start = 0; + while start < samples.len() { + let end = (start + fft_size).min(samples.len()); + let mut frame = vec![0.0f32; fft_size]; + frame[..end - start].copy_from_slice(&samples[start..end]); + let windowed = apply(&frame, &window); + frames.push(real_magnitude_spectrum(&windowed)); + start += hop; + } + } + + Spectrogram { + frames, + sample_rate, + fft_size, + hop, + } +} + +#[cfg(test)] +mod tests { + use super::*; + use std::f64::consts::PI; + + #[test] + fn frame_count_follows_hop() { + // 1000 samples, hop 250 -> frames start at 0,250,500,750 = 4 frames. + let sig = vec![0.0f32; 1000]; + let s = analyze(&sig, 48_000, 512, 250); + assert_eq!(s.frames.len(), 4); + assert_eq!(s.bins(), 512 / 2 + 1); + } + + #[test] + fn tone_lands_in_expected_bin() { + // A 3 kHz tone at 48 kHz with a 1024-pt FFT -> bin ≈ 3000/(48000/1024) = 64. + let sr = 48_000; + let freq = 3000.0; + let sig: Vec = (0..4096) + .map(|i| (2.0 * PI * freq * i as f64 / sr as f64).sin() as f32) + .collect(); + let s = analyze(&sig, sr, 1024, 512); + let mid = &s.frames[s.frames.len() / 2]; + let peak_bin = mid + .iter() + .enumerate() + .max_by(|a, b| a.1.partial_cmp(b.1).unwrap()) + .unwrap() + .0; + let peak_hz = s.bin_hz(peak_bin); + assert!((peak_hz - freq as f32).abs() < 100.0, "peak at {peak_hz} Hz, want {freq}"); + } + + #[test] + fn empty_input_yields_no_frames() { + let s = analyze(&[], 48_000, 512, 256); + assert!(s.frames.is_empty()); + } +} diff --git a/src/dsp/wav.rs b/src/dsp/wav.rs new file mode 100644 index 0000000..37d5975 --- /dev/null +++ b/src/dsp/wav.rs @@ -0,0 +1,139 @@ +//! A minimal canonical-PCM WAV reader (the analysis-side counterpart to +//! `audio::recorder::WavWriter`). Handles 16-bit integer PCM — the only format +//! PeerSpeak writes — in mono or interleaved multi-channel, downmixing to mono +//! `f32` in `[-1, 1]` for analysis. I/O lives here, at the edge; the rest of the +//! `dsp` module is pure. + +use std::fs; +use std::path::Path; + +/// A decoded WAV: mono `f32` samples plus the original sample rate. +pub struct WavData { + pub samples: Vec, + pub sample_rate: u32, +} + +/// Reads a 16-bit PCM WAV file, downmixing any channels to mono `f32`. +/// +/// Returns `Err` with a human-readable reason if the file is missing, truncated, +/// not RIFF/WAVE, or not 16-bit PCM. The parser walks the chunk list rather than +/// assuming a fixed 44-byte header, so files with extra chunks (`LIST`, `fact`, +/// …) still read. +pub fn read(path: &Path) -> Result { + let bytes = fs::read(path).map_err(|e| format!("cannot read {}: {e}", path.display()))?; + if bytes.len() < 12 || &bytes[0..4] != b"RIFF" || &bytes[8..12] != b"WAVE" { + return Err("not a RIFF/WAVE file".to_string()); + } + + let mut channels = 0u16; + let mut sample_rate = 0u32; + let mut bits = 0u16; + let mut data: Option<&[u8]> = None; + + // Walk chunks starting after the 12-byte RIFF/WAVE header. + let mut pos = 12usize; + while pos + 8 <= bytes.len() { + let id = &bytes[pos..pos + 4]; + let size = u32::from_le_bytes([bytes[pos + 4], bytes[pos + 5], bytes[pos + 6], bytes[pos + 7]]) as usize; + let body_start = pos + 8; + let body_end = (body_start + size).min(bytes.len()); + match id { + b"fmt " if size >= 16 => { + let fmt = &bytes[body_start..body_end]; + let audio_format = u16::from_le_bytes([fmt[0], fmt[1]]); + channels = u16::from_le_bytes([fmt[2], fmt[3]]); + sample_rate = u32::from_le_bytes([fmt[4], fmt[5], fmt[6], fmt[7]]); + bits = u16::from_le_bytes([fmt[14], fmt[15]]); + if audio_format != 1 { + return Err(format!("unsupported WAV format tag {audio_format} (need PCM=1)")); + } + } + b"data" => { + data = Some(&bytes[body_start..body_end]); + } + _ => {} + } + // Chunks are word-aligned: an odd size is padded with one byte. + pos = body_start + size + (size & 1); + } + + if bits != 16 { + return Err(format!("only 16-bit PCM supported, got {bits}-bit")); + } + let channels = channels.max(1); + let data = data.ok_or("no data chunk")?; + + // Interleaved S16LE -> per-frame channel average -> mono f32. + let frame_bytes = 2 * channels as usize; + let mut samples = Vec::with_capacity(data.len() / frame_bytes.max(1)); + for frame in data.chunks_exact(frame_bytes) { + let mut acc = 0i32; + for ch in frame.chunks_exact(2) { + acc += i16::from_le_bytes([ch[0], ch[1]]) as i32; + } + let avg = acc as f32 / channels as f32; + samples.push(avg / 32768.0); + } + + Ok(WavData { samples, sample_rate }) +} + +/// Writes mono `f32` samples (clamped to `[-1, 1]`) as a 16-bit PCM WAV. Used by +/// the `specview gen` command to materialize conjured test signals. +pub fn write(path: &Path, samples: &[f32], sample_rate: u32) -> Result<(), String> { + let pcm: Vec = samples + .iter() + .map(|&s| (s.clamp(-1.0, 1.0) * 32767.0).round() as i16) + .collect(); + let data_bytes = (pcm.len() * 2) as u32; + let byte_rate = sample_rate * 2; // mono, 2 bytes/sample + let mut out = Vec::with_capacity(44 + pcm.len() * 2); + out.extend_from_slice(b"RIFF"); + out.extend_from_slice(&(36 + data_bytes).to_le_bytes()); + out.extend_from_slice(b"WAVE"); + out.extend_from_slice(b"fmt "); + out.extend_from_slice(&16u32.to_le_bytes()); + out.extend_from_slice(&1u16.to_le_bytes()); // PCM + out.extend_from_slice(&1u16.to_le_bytes()); // mono + out.extend_from_slice(&sample_rate.to_le_bytes()); + out.extend_from_slice(&byte_rate.to_le_bytes()); + out.extend_from_slice(&2u16.to_le_bytes()); // block align + out.extend_from_slice(&16u16.to_le_bytes()); // bits + out.extend_from_slice(b"data"); + out.extend_from_slice(&data_bytes.to_le_bytes()); + for s in pcm { + out.extend_from_slice(&s.to_le_bytes()); + } + fs::write(path, out).map_err(|e| format!("cannot write {}: {e}", path.display())) +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn write_then_read_round_trips() { + let dir = std::env::temp_dir(); + let path = dir.join(format!("dsp-wav-rt-{}.wav", std::process::id())); + let sig: Vec = (0..1000).map(|i| (i as f32 * 0.01).sin() * 0.5).collect(); + write(&path, &sig, 48_000).unwrap(); + + let back = read(&path).unwrap(); + assert_eq!(back.sample_rate, 48_000); + assert_eq!(back.samples.len(), sig.len()); + // 16-bit quantization tolerance. + for (a, b) in sig.iter().zip(&back.samples) { + assert!((a - b).abs() < 1e-3, "{a} vs {b}"); + } + let _ = fs::remove_file(&path); + } + + #[test] + fn rejects_non_wave() { + let dir = std::env::temp_dir(); + let path = dir.join(format!("dsp-wav-bad-{}.bin", std::process::id())); + fs::write(&path, b"not a wav at all").unwrap(); + assert!(read(&path).is_err()); + let _ = fs::remove_file(&path); + } +} diff --git a/src/dsp/window.rs b/src/dsp/window.rs new file mode 100644 index 0000000..85a3685 --- /dev/null +++ b/src/dsp/window.rs @@ -0,0 +1,55 @@ +//! Analysis windows for the STFT. +//! +//! A raw rectangular frame leaks spectral energy across bins (the abrupt edges +//! look like discontinuities to the FFT). A Hann window tapers each frame to +//! zero at its edges, trading a little main-lobe width for much lower side-lobe +//! leakage — the standard choice for a spectrogram. + +use std::f64::consts::PI; + +/// Returns a length-`n` periodic Hann window, `w[i] = 0.5·(1 - cos(2πi/n))`. +/// +/// The *periodic* form (denominator `n`, not `n-1`) is used because STFT frames +/// tile the signal; it gives perfect overlap-add reconstruction at 50% hop. +pub fn hann(n: usize) -> Vec { + if n <= 1 { + return vec![1.0; n]; + } + (0..n) + .map(|i| (0.5 - 0.5 * (2.0 * PI * i as f64 / n as f64).cos()) as f32) + .collect() +} + +/// Applies `window` to `frame` element-wise into a new buffer. Lengths must match. +pub fn apply(frame: &[f32], window: &[f32]) -> Vec { + debug_assert_eq!(frame.len(), window.len()); + frame.iter().zip(window).map(|(&s, &w)| s * w).collect() +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn hann_endpoints_are_zero_and_center_is_one() { + let w = hann(8); + assert!(w[0].abs() < 1e-6, "first sample tapers to ~0, got {}", w[0]); + // Periodic Hann peaks at the midpoint n/2. + assert!((w[4] - 1.0).abs() < 1e-6, "center peaks at 1, got {}", w[4]); + } + + #[test] + fn hann_is_symmetric_about_center() { + let w = hann(16); + // Periodic window is symmetric across indices 1..n-1. + for i in 1..8 { + assert!((w[i] - w[16 - i]).abs() < 1e-6); + } + } + + #[test] + fn apply_scales_samples() { + let out = apply(&[2.0, 2.0], &[0.5, 0.25]); + assert_eq!(out, vec![1.0, 0.5]); + } +} diff --git a/src/lib.rs b/src/lib.rs index c85a211..77c8fa3 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -1,5 +1,6 @@ pub mod audio; pub mod codec; +pub mod dsp; pub mod network; pub mod core; pub mod app;