From eeb725f996c1dd6e3e316acd039d08bc8a4a2131 Mon Sep 17 00:00:00 2001 From: Mollusk Date: Sun, 14 Jun 2026 02:39:01 -0400 Subject: [PATCH] feat(dsp): in-process NLMS echo canceller + headless AEC test loop MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Stage A+B of the in-process AEC: a Normalized LMS adaptive-filter echo canceller and a synthetic echo-path simulator, driven end-to-end by a new `specview aec` command so cancellation can be measured fully headless on conjured signals (no speakers/mic needed). - src/dsp/echo_path.rs: EchoPath — synthetic acoustic echo (bulk delay + exponentially-decaying diffuse RIR, energy-normalized to a target attenuation) convolved over a far-end signal. Gives ground-truth echo. - src/dsp/aec.rs: Nlms — sample-at-a-time NLMS adaptive FIR (ring-buffered reference history, energy-normalized update, freezable for double-talk). Cleaned output = mic minus the learned echo estimate. - specview aec: far -> sim echo -> (+ optional near-end) -> cancel -> measure. ERLE via oracle residual (cleaned - near), broadband + early/late (convergence) + per voice band, optional before/after spectrograms. - 6 new tests (path delay/attenuation, NLMS convergence >20 dB on a known path, frozen-filter no-op, near-end passthrough). 33 dsp tests total. Verified: single-talk pink-noise echo cancels +14.5 -> +32.0 dB ERLE as the filter converges; double-talk (no DTD yet) drives ERLE negative as the filter diverges onto the near-end tone — the motivating result for Stage C (DTD). Co-Authored-By: Claude Opus 4.8 --- src/bin/specview.rs | 99 ++++++++++++++++++++++++- src/dsp/aec.rs | 171 +++++++++++++++++++++++++++++++++++++++++++ src/dsp/echo_path.rs | 129 ++++++++++++++++++++++++++++++++ src/dsp/mod.rs | 2 + 4 files changed, 399 insertions(+), 2 deletions(-) create mode 100644 src/dsp/aec.rs create mode 100644 src/dsp/echo_path.rs diff --git a/src/bin/specview.rs b/src/bin/specview.rs index 685661f..bad4007 100644 --- a/src/bin/specview.rs +++ b/src/bin/specview.rs @@ -17,7 +17,7 @@ use std::collections::HashMap; use std::path::Path; use std::process::ExitCode; -use peerspeak::dsp::{generators, metrics, render, stft, wav}; +use peerspeak::dsp::{aec, echo_path, generators, metrics, render, stft, wav}; const SAMPLE_RATE: u32 = 48_000; @@ -33,6 +33,7 @@ fn main() -> ExitCode { "show" => cmd_show(rest), "gen" => cmd_gen(rest), "erle" => cmd_erle(rest), + "aec" => cmd_aec(rest), "-h" | "--help" | "help" => { println!("{}", USAGE); Ok(()) @@ -53,7 +54,9 @@ 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]"; + specview erle [--show]\n\ + specview aec --far far.wav [--near near.wav] [--delay N] [--tail N] [--atten A]\n\ + [--taps N] [--mu M] [--out cleaned.wav] [--show]"; /// Renders a WAV as a spectrogram plus a numeric summary. fn cmd_show(args: &[String]) -> Result<(), String> { @@ -167,6 +170,98 @@ fn cmd_erle(args: &[String]) -> Result<(), String> { Ok(()) } +/// The full headless AEC loop: take a far-end signal, run it through a synthetic +/// echo path to make the echo a mic would hear, optionally add near-end speech, +/// then cancel with the NLMS filter and measure how much echo was removed. +/// +/// Because the simulator hands us the true echo and near-end separately, ERLE is +/// measured against ground truth (residual echo = cleaned − near) — exact even +/// when near-end talk is present. +fn cmd_aec(args: &[String]) -> Result<(), String> { + let (_positional, flags) = parse_args(args); + let far_path = flags.get("far").ok_or("aec needs --far ")?; + let far = wav::read(Path::new(far_path))?; + let sr = far.sample_rate; + + // Optional near-end (the local voice the canceller must preserve). + let near = match flags.get("near") { + Some(p) => { + let d = wav::read(Path::new(p))?; + let mut v = d.samples; + v.resize(far.samples.len(), 0.0); // align length to far-end + v + } + None => vec![0.0f32; far.samples.len()], + }; + + // Build the echo and the mic signal (echo + near-end). + let delay = flags.usize_or("delay", 480); // ~10 ms at 48 kHz + let tail = flags.usize_or("tail", 960); // ~20 ms reverb tail + let atten = flags.f32_or("atten", 0.5); // ~ -6 dB echo + let path = echo_path::EchoPath::synthetic(delay, tail, atten, flags.u64_or("seed", 1)); + let echo = path.apply(&far.samples); + let mic: Vec = echo.iter().zip(&near).map(|(&e, &n)| e + n).collect(); + + // Cancel. Default taps comfortably cover the synthetic path length. + let taps = flags.usize_or("taps", (delay + tail).next_power_of_two()); + let mu = flags.f32_or("mu", 0.5); + let mut canceller = aec::Nlms::new(taps, mu, 1e-6); + let cleaned = canceller.process(&far.samples, &mic); + + // Residual echo (oracle): whatever's left after removing the known near-end. + let residual: Vec = cleaned.iter().zip(&near).map(|(&c, &n)| c - n).collect(); + + let broadband = metrics::erle(&echo, &residual); + let q = echo.len() / 4; + let early = metrics::erle(&echo[..q], &residual[..q]); + let late = metrics::erle(&echo[3 * q..], &residual[3 * q..]); + let has_near = flags.present("near"); + + println!("AEC simulation"); + println!( + " echo path: delay {delay} ({:.0} ms), tail {tail} ({:.0} ms), atten {atten} ({:.1} dB)", + 1000.0 * delay as f32 / sr as f32, + 1000.0 * tail as f32 / sr as f32, + metrics::dbfs(atten), + ); + println!(" filter: {taps} taps, mu {mu}{}", if has_near { " (with near-end / double-talk)" } else { "" }); + println!(" mic before: {:.1} dBFS rms", metrics::dbfs(metrics::rms(&mic))); + println!(" residual echo after: {:.1} dBFS rms", metrics::dbfs(metrics::rms(&residual))); + println!(" ERLE broadband: {broadband:+.1} dB"); + println!(" ERLE early/late: {early:+.1} -> {late:+.1} dB (rise = filter converging)"); + + // Per-band residual echo, to see where any leak lives. + let se = stft::analyze(&echo, sr, 2048, 512); + let sr_res = stft::analyze(&residual, sr, 2048, 512); + let avg_e = average_frame(&se); + let avg_r = average_frame(&sr_res); + println!(" ERLE by band:"); + for band in metrics::VOICE_BANDS { + let ee = metrics::band_energy(&avg_e, |b| se.bin_hz(b), band.low_hz, band.high_hz); + let er = metrics::band_energy(&avg_r, |b| sr_res.bin_hz(b), band.low_hz, band.high_hz); + let db = if ee > 1e-12 && er > 1e-12 { + 10.0 * (ee / er).log10() + } else if ee > 1e-12 { + 120.0 + } else { + 0.0 + }; + println!(" {:<18} {:+6.1} dB", band.label, db); + } + + if let Some(out) = flags.get("out") { + wav::write(Path::new(out), &cleaned, sr)?; + println!(" wrote cleaned output -> {out}"); + } + if flags.present("show") { + println!("\n--- mic (echo present) ---"); + print!("{}", render::render(&stft::analyze(&mic, sr, 2048, 512), &render::RenderOpts::default())); + println!("\n--- cleaned (post-AEC) ---"); + print!("{}", render::render(&stft::analyze(&cleaned, sr, 2048, 512), &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!( diff --git a/src/dsp/aec.rs b/src/dsp/aec.rs new file mode 100644 index 0000000..a54454a --- /dev/null +++ b/src/dsp/aec.rs @@ -0,0 +1,171 @@ +//! An in-process acoustic echo canceller built on a **Normalized Least Mean +//! Squares (NLMS)** adaptive FIR filter. +//! +//! The idea: we know the far-end signal exactly (it's what we're about to play +//! out), and the mic picks up an echo of it shaped by the unknown room. An +//! adaptive filter continuously estimates that room response from the data: +//! each sample it predicts the echo as a weighted sum of recent far-end samples, +//! subtracts that prediction from the mic, and nudges its weights to shrink the +//! leftover. With only far-end present, the leftover (the "error" signal) drives +//! the weights toward the true echo path; once converged, the error *is* the +//! cleaned near-end audio. +//! +//! NLMS over plain LMS divides the update by the reference's instantaneous energy, +//! which makes the convergence speed independent of how loud the far-end is — the +//! reason it's the workhorse of practical echo cancellation. +//! +//! This is the pure, sample-at-a-time core. Double-talk handling and a residual +//! suppressor build on top of it; integration into the live audio path wraps it. + +/// A single-channel NLMS echo canceller. +pub struct Nlms { + /// Adaptive filter coefficients (the estimated echo-path impulse response). + /// `weights[k]` multiplies the far-end sample `k` steps in the past. + weights: Vec, + /// Ring buffer of recent far-end (reference) samples, same length as weights. + history: Vec, + /// Index in `history` where the newest sample currently lives. + pos: usize, + /// Step size (0 < mu < 2). Larger = faster adaptation but more misadjustment. + mu: f32, + /// Regularization added to the energy normalizer to avoid divide-by-zero and + /// runaway updates when the far-end is near silent. + eps: f32, + /// When false, weights are held fixed (used during double-talk). + adapting: bool, +} + +impl Nlms { + /// Creates an NLMS canceller with `taps` filter coefficients (should cover the + /// echo path's length: delay + reverb tail), step size `mu` in `(0, 2)`, and + /// regularization `eps`. Sensible defaults: `mu = 0.5`, `eps = 1e-6`. + pub fn new(taps: usize, mu: f32, eps: f32) -> Self { + let taps = taps.max(1); + Nlms { + weights: vec![0.0; taps], + history: vec![0.0; taps], + pos: 0, + mu, + eps, + adapting: true, + } + } + + /// Number of filter taps. + pub fn taps(&self) -> usize { + self.weights.len() + } + + /// Enables or disables weight adaptation. Disable during double-talk so the + /// near-end voice doesn't corrupt the learned echo path. + pub fn set_adapting(&mut self, on: bool) { + self.adapting = on; + } + + /// A copy of the current estimated impulse response (the filter weights), + /// newest-tap-first. Useful for tests/inspection. + pub fn weights(&self) -> &[f32] { + &self.weights + } + + /// Processes one sample pair and returns the echo-cancelled output. + /// + /// `reference` is the far-end sample about to be (or being) played out; + /// `mic` is the simultaneously-captured near-end mic sample (near-end voice + + /// echo of the far-end + noise). The returned value is `mic` minus the + /// filter's echo estimate — i.e. the cleaned near-end. + pub fn process_sample(&mut self, reference: f32, mic: f32) -> f32 { + let l = self.weights.len(); + // Store newest reference at the current ring position. + self.history[self.pos] = reference; + + // Predict the echo: weights[k] pairs with the sample k steps back, which + // sits at ring index (pos - k). Accumulate the reference energy too. + let mut echo_est = 0.0f32; + let mut energy = 0.0f32; + for k in 0..l { + let idx = (self.pos + l - k) % l; + let s = self.history[idx]; + echo_est += self.weights[k] * s; + energy += s * s; + } + + let error = mic - echo_est; + + // NLMS weight update: w += mu * error * x / (||x||^2 + eps). + if self.adapting { + let g = self.mu * error / (energy + self.eps); + for k in 0..l { + let idx = (self.pos + l - k) % l; + self.weights[k] += g * self.history[idx]; + } + } + + // Advance the ring so the next sample overwrites the oldest entry. + self.pos = (self.pos + 1) % l; + error + } + + /// Convenience: cancels a whole buffer, returning the cleaned output. `mic` + /// and `reference` are processed pairwise over their shared length. + pub fn process(&mut self, reference: &[f32], mic: &[f32]) -> Vec { + let n = reference.len().min(mic.len()); + let mut out = Vec::with_capacity(n); + for i in 0..n { + out.push(self.process_sample(reference[i], mic[i])); + } + out + } +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::dsp::echo_path::EchoPath; + use crate::dsp::generators::white_noise; + use crate::dsp::metrics::erle; + + /// With only far-end present (no near-end talk), the filter should identify + /// the echo path and drive ERLE high — and the late part of the call should + /// be far better cancelled than the early part (it learns over time). + #[test] + fn converges_on_a_known_echo_path() { + let far = white_noise(0.5, 48_000, 1); + let path = EchoPath::synthetic(64, 192, 0.5, 7); // 256-tap path + let echo = path.apply(&far); + // Mic is echo-only (single-talk): the ideal case for adaptation. + let mut aec = Nlms::new(512, 0.5, 1e-6); + let cleaned = aec.process(&far, &echo); + + // Compare ERLE on the first vs last quarter. + let q = echo.len() / 4; + let early = erle(&echo[..q], &cleaned[..q]); + let late = erle(&echo[3 * q..], &cleaned[3 * q..]); + assert!(late > early + 10.0, "should improve markedly: early {early:.1} late {late:.1}"); + assert!(late > 20.0, "converged ERLE should exceed 20 dB, got {late:.1}"); + } + + #[test] + fn frozen_filter_does_not_adapt() { + let far = white_noise(0.5, 4000, 2); + let echo = EchoPath::synthetic(16, 48, 0.5, 3).apply(&far); + let mut aec = Nlms::new(128, 0.5, 1e-6); + aec.set_adapting(false); + aec.process(&far, &echo); + // Weights must remain exactly zero with adaptation off. + assert!(aec.weights().iter().all(|&w| w == 0.0)); + } + + #[test] + fn passes_near_end_through_when_no_echo() { + // No far-end (reference silent) -> filter predicts ~0 -> near-end passes + // through essentially unchanged. + let near = white_noise(0.3, 2000, 9); + let silent_ref = vec![0.0f32; near.len()]; + let mut aec = Nlms::new(128, 0.5, 1e-6); + let out = aec.process(&silent_ref, &near); + for (a, b) in near.iter().zip(&out) { + assert!((a - b).abs() < 1e-6, "near-end should pass through: {a} vs {b}"); + } + } +} diff --git a/src/dsp/echo_path.rs b/src/dsp/echo_path.rs new file mode 100644 index 0000000..e8ca149 --- /dev/null +++ b/src/dsp/echo_path.rs @@ -0,0 +1,129 @@ +//! A synthetic acoustic echo path — the model the AEC test harness uses in place +//! of real speakers + a room + a microphone. +//! +//! When you're on a speakerphone call, the far-end voice comes out of your +//! speaker, bounces around the room, and re-enters your mic delayed, coloured, +//! and attenuated. That transformation is a linear filter: an **impulse +//! response** (IR). Convolving the far-end signal with an IR reproduces the echo +//! a mic would capture — fully deterministically, so the canceller can be +//! measured headless against ground truth. +//! +//! The IR here is `delay` samples of silence (bulk acoustic propagation / +//! buffering latency) followed by an exponentially-decaying diffuse tail (the +//! room reverb), scaled so the whole echo sits `~attenuation` below the far-end. + +use super::generators; + +/// A linear echo path expressed as an impulse response. +pub struct EchoPath { + /// The impulse response; `echo[n] = Σ_k ir[k]·far[n-k]`. + pub ir: Vec, +} + +impl EchoPath { + /// Builds the path directly from a caller-supplied impulse response. + pub fn from_ir(ir: Vec) -> Self { + EchoPath { ir } + } + + /// A synthetic room echo: `delay` samples of pure delay, then a `tail`-sample + /// exponentially-decaying diffuse reflection cluster. `attenuation` is the + /// peak linear gain of the path (e.g. `0.5` ≈ -6 dB echo). `seed` makes the + /// diffuse tail reproducible. + /// + /// A strong early reflection is placed at the delay (the dominant first bounce) + /// followed by a decaying noisy tail, then the whole IR is normalized so its + /// energy gain matches `attenuation` — a realistic, non-trivial path for the + /// adaptive filter to identify. + pub fn synthetic(delay: usize, tail: usize, attenuation: f32, seed: u64) -> Self { + let tail = tail.max(1); + let mut ir = vec![0.0f32; delay + tail]; + let noise = generators::white_noise(1.0, tail, seed); + // Time constant: decay to ~e^-4 (≈1.8%) by the end of the tail. + let tau = tail as f32 / 4.0; + for i in 0..tail { + let env = (-(i as f32) / tau).exp(); + // Dominant direct reflection at i=0, diffuse decaying noise after. + let direct = if i == 0 { 1.0 } else { 0.4 * noise[i] }; + ir[delay + i] = env * direct; + } + // Normalize so the path's RMS gain equals `attenuation` (energy-based, so + // ERLE targets are predictable regardless of delay/tail choices). + let energy: f32 = ir.iter().map(|&x| x * x).sum::().sqrt(); + if energy > 1e-9 { + let scale = attenuation / energy; + for x in &mut ir { + *x *= scale; + } + } + EchoPath { ir } + } + + /// Number of taps in the impulse response (its total length). + pub fn len(&self) -> usize { + self.ir.len() + } + + pub fn is_empty(&self) -> bool { + self.ir.is_empty() + } + + /// Convolves `far` with the impulse response, returning the echo signal the + /// mic would capture. Output length matches `far` (the convolution tail past + /// the input end is truncated — the far-end keeps going in a real call). + pub fn apply(&self, far: &[f32]) -> Vec { + let n = far.len(); + let m = self.ir.len(); + let mut echo = vec![0.0f32; n]; + for (k, &h) in self.ir.iter().enumerate() { + if h == 0.0 { + continue; + } + // echo[i] += h * far[i-k] for all valid i. + for i in k..n { + echo[i] += h * far[i - k]; + } + } + let _ = m; + echo + } +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::dsp::metrics::rms; + + #[test] + fn pure_delay_shifts_the_signal() { + // IR = single unit tap at index 3 -> echo is far delayed by 3 samples. + let path = EchoPath::from_ir({ + let mut ir = vec![0.0; 4]; + ir[3] = 1.0; + ir + }); + let far = vec![1.0, 2.0, 3.0, 4.0, 5.0]; + let echo = path.apply(&far); + assert_eq!(echo, vec![0.0, 0.0, 0.0, 1.0, 2.0]); + } + + #[test] + fn synthetic_path_attenuates_to_target() { + // A white-noise far-end through a path normalized to 0.5 RMS gain should + // produce an echo roughly 0.5x the far-end RMS (≈ -6 dB). + let far = generators::white_noise(0.5, 48_000, 1); + let path = EchoPath::synthetic(480, 480, 0.5, 99); + let echo = path.apply(&far); + let ratio = rms(&echo) / rms(&far); + assert!((0.3..0.7).contains(&ratio), "echo/far rms ratio {ratio} off target"); + } + + #[test] + fn synthetic_path_has_leading_delay() { + let path = EchoPath::synthetic(100, 200, 0.5, 5); + // The first `delay` IR taps are silent (pure propagation delay). + assert!(path.ir[..100].iter().all(|&x| x == 0.0)); + assert!(path.ir[100..].iter().any(|&x| x != 0.0)); + assert_eq!(path.len(), 300); + } +} diff --git a/src/dsp/mod.rs b/src/dsp/mod.rs index 96f8ae6..0fd7967 100644 --- a/src/dsp/mod.rs +++ b/src/dsp/mod.rs @@ -9,6 +9,8 @@ //! `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 aec; +pub mod echo_path; pub mod fft; pub mod generators; pub mod metrics;