atagen/headroom
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

stage 2

Browse source
This commit is contained in:
atagen 2026-05-19 16:33:09 +10:00
commit ca1910de60
39 changed files with 6328 additions and 0 deletions
Download patch file Download diff file Expand all files Collapse all files

19
crates/headroom-dsp/Cargo.toml Normal file
View file

@ -0,0 +1,19 @@
[package]
name = "headroom-dsp"
description = "DSP kernels for Headroom: true-peak limiter, compressor, AGC envelope helpers."
version.workspace = true
edition.workspace = true
rust-version.workspace = true
license = "MPL-2.0"
homepage.workspace = true
repository.workspace = true
authors.workspace = true
readme = "README.md"
[dependencies]
# Kept intentionally empty. The DSP crate must build clean on any host
# and is the most reusable piece in the workspace. If you find yourself
# wanting to add a dependency here, think twice.
[features]
default = []

19
crates/headroom-dsp/README.md Normal file
View file

@ -0,0 +1,19 @@
# headroom-dsp
DSP kernels for Headroom. Pure Rust, no dependencies.
- `Limiter` — feed-forward true-peak brickwall with configurable
oversampling (1/2/4/8×), lookahead, hold, and release.
- `Compressor` — log-domain feed-forward with peak or RMS detector,
soft knee, attack/release, and optional auto-makeup.
- `AttackRelease` — exponential envelope follower (peak / inverse-gain
modes).
- `DelayLine`, `SlidingMaxBuffer`, `PolyphaseUpsampler`,
`PolyphaseDownsampler` — supporting building blocks.
All processors are allocation-free in their `process_*` methods.
Construction allocates; do not construct in the audio thread.
## License
MPL-2.0.

268
crates/headroom-dsp/src/compressor.rs Normal file
View file

@ -0,0 +1,268 @@
//! Feed-forward dynamics compressor.
//!
//! Log-domain detector → static curve with soft knee → smoothed
//! envelope → linear gain → apply (no internal delay).
use crate::util::{db_to_lin, lin_to_db, time_to_alpha};
/// Detector type used to build the side-chain signal.
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum Detector {
/// Maximum of `|left|, |right|`. Fast, low CPU, slightly more
/// percussive feel on transients.
Peak,
/// One-pole low-passed mean square. Smoother, more relaxed on
/// percussive material.
Rms,
}
/// Compressor parameters.
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct CompressorConfig {
/// Threshold in dBFS. Inputs above this start compressing.
pub threshold_db: f32,
/// Compression ratio (>= 1.0).
pub ratio: f32,
/// Soft-knee width in dB. `0.0` is a hard knee.
pub knee_db: f32,
/// Attack time in ms.
pub attack_ms: f32,
/// Release time in ms.
pub release_ms: f32,
/// Makeup gain in dB. `None` selects an automatic mild boost.
pub makeup_db: Option<f32>,
/// Detector type.
pub detector: Detector,
/// RMS window in ms (only used when `detector == Rms`).
pub rms_window_ms: f32,
}
impl Default for CompressorConfig {
fn default() -> Self {
Self {
threshold_db: -24.0,
ratio: 2.5,
knee_db: 6.0,
attack_ms: 10.0,
release_ms: 100.0,
makeup_db: None,
detector: Detector::Peak,
rms_window_ms: 5.0,
}
}
}
impl CompressorConfig {
/// Clamp ratio to `>= 1.0`, knee/attack/release/window to `>= 0`.
#[must_use]
pub fn sanitized(mut self) -> Self {
if self.ratio < 1.0 {
self.ratio = 1.0;
}
self.knee_db = self.knee_db.max(0.0);
self.attack_ms = self.attack_ms.max(0.0);
self.release_ms = self.release_ms.max(0.0);
self.rms_window_ms = self.rms_window_ms.max(0.1);
self
}
}
/// Feed-forward compressor (stereo-linked).
pub struct Compressor {
cfg: CompressorConfig,
sample_rate: f32,
envelope_db: f32,
attack_alpha: f32,
release_alpha: f32,
rms_state: f32,
rms_alpha: f32,
last_gr_db: f32,
}
impl Compressor {
/// Construct with the given config and input sample rate.
#[must_use]
pub fn new(cfg: CompressorConfig, sample_rate: f32) -> Self {
let cfg = cfg.sanitized();
Self {
cfg,
sample_rate,
envelope_db: -200.0,
attack_alpha: time_to_alpha(cfg.attack_ms, sample_rate),
release_alpha: time_to_alpha(cfg.release_ms, sample_rate),
rms_state: 0.0,
rms_alpha: time_to_alpha(cfg.rms_window_ms, sample_rate),
last_gr_db: 0.0,
}
}
/// Active configuration.
#[must_use]
pub fn config(&self) -> CompressorConfig {
self.cfg
}
/// Most recent gain reduction in dB (negative when compressing).
#[must_use]
pub fn gain_reduction_db(&self) -> f32 {
self.last_gr_db
}
/// Update parameters. Recomputes alphas. Envelope state is kept,
/// so live tweaks don't pop.
pub fn set_config(&mut self, cfg: CompressorConfig) {
let cfg = cfg.sanitized();
self.cfg = cfg;
self.attack_alpha = time_to_alpha(cfg.attack_ms, self.sample_rate);
self.release_alpha = time_to_alpha(cfg.release_ms, self.sample_rate);
self.rms_alpha = time_to_alpha(cfg.rms_window_ms, self.sample_rate);
}
/// Process one stereo frame.
pub fn process_frame(&mut self, left: f32, right: f32) -> (f32, f32) {
let det_lin = match self.cfg.detector {
Detector::Peak => left.abs().max(right.abs()),
Detector::Rms => {
let sq = 0.5 * left.mul_add(left, right * right);
self.rms_state += self.rms_alpha * (sq - self.rms_state);
self.rms_state.max(0.0).sqrt()
}
};
let det_db = lin_to_db(det_lin.max(1e-20));
if det_db > self.envelope_db {
self.envelope_db += self.attack_alpha * (det_db - self.envelope_db);
} else {
self.envelope_db += self.release_alpha * (det_db - self.envelope_db);
}
let gr_db = static_curve_gain_reduction(
self.envelope_db,
self.cfg.threshold_db,
self.cfg.ratio,
self.cfg.knee_db,
);
let makeup_db = self
.cfg
.makeup_db
.unwrap_or_else(|| auto_makeup(self.cfg.threshold_db, self.cfg.ratio));
let lin_gain = db_to_lin(-gr_db + makeup_db);
self.last_gr_db = -gr_db;
(left * lin_gain, right * lin_gain)
}
/// Reset envelopes and detector state.
pub fn reset(&mut self) {
self.envelope_db = -200.0;
self.rms_state = 0.0;
self.last_gr_db = 0.0;
}
}
/// Static compression curve. Returns the positive gain reduction (in
/// dB) that should be subtracted from the input level.
fn static_curve_gain_reduction(
input_db: f32,
threshold_db: f32,
ratio: f32,
knee_db: f32,
) -> f32 {
let over = input_db - threshold_db;
if knee_db > 0.0 && over > -knee_db * 0.5 && over < knee_db * 0.5 {
// Quadratic soft knee.
let x = over + knee_db * 0.5;
let factor = (x * x) / (2.0 * knee_db);
factor * (1.0 - 1.0 / ratio)
} else if over > 0.0 {
over * (1.0 - 1.0 / ratio)
} else {
0.0
}
}
/// Mild auto-makeup: compensate for half the static gain reduction at
/// 0 dBFS. Conservative on purpose — the limiter is downstream and we
/// don't want to push it.
fn auto_makeup(threshold_db: f32, ratio: f32) -> f32 {
let gr_at_zero = (-threshold_db).max(0.0) * (1.0 - 1.0 / ratio);
gr_at_zero * 0.5
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn below_threshold_is_unity_minus_makeup() {
let mut c = Compressor::new(CompressorConfig::default(), 48_000.0);
// Drive a long, low signal and check we land at expected gain.
let mut last = (0.0_f32, 0.0_f32);
for _ in 0..10_000 {
last = c.process_frame(0.01, 0.01);
}
// Below threshold: gain reduction is zero, only makeup applied.
let makeup_db = auto_makeup(-24.0, 2.5);
let expected = 0.01 * db_to_lin(makeup_db);
assert!(
(last.0 - expected).abs() < 1e-3,
"got {} expected {}",
last.0,
expected
);
assert!(c.gain_reduction_db().abs() < 0.1);
}
#[test]
fn above_threshold_reduces_gain() {
let cfg = CompressorConfig {
threshold_db: -20.0,
ratio: 4.0,
knee_db: 0.0, // hard knee for clean math
attack_ms: 0.1,
release_ms: 0.1,
makeup_db: Some(0.0), // no makeup so we test pure reduction
..CompressorConfig::default()
};
let mut c = Compressor::new(cfg, 48_000.0);
// Drive ~-6 dBFS = 0.5 in linear.
let target = 0.5_f32;
let mut last_out = 0.0;
for _ in 0..2_000 {
let (l, _) = c.process_frame(target, target);
last_out = l;
}
// Input is 14 dB above threshold. With ratio 4, GR = 14*(1-0.25) = 10.5 dB.
// Expected output: -6 - 10.5 = -16.5 dB linear = 0.1496.
let expected_db = -6.0 - 14.0 * (1.0 - 0.25);
let expected_lin = db_to_lin(expected_db);
let got_db = lin_to_db(last_out);
assert!(
(got_db - expected_db).abs() < 0.5,
"got {got_db} expected {expected_db}"
);
assert!(c.gain_reduction_db() < -5.0, "gr was {}", c.gain_reduction_db());
let _ = expected_lin;
}
#[test]
fn ratio_below_one_is_clamped() {
let cfg = CompressorConfig {
ratio: 0.5,
..CompressorConfig::default()
}
.sanitized();
assert_eq!(cfg.ratio, 1.0);
}
#[test]
fn static_curve_at_threshold_with_soft_knee() {
// At exactly threshold, soft knee contributes exactly half the
// ratio's compression amount at the upper knee shoulder.
let gr = static_curve_gain_reduction(-24.0, -24.0, 4.0, 6.0);
// At over==0 inside the knee, x = knee/2, factor = knee/8.
// GR = knee/8 * (1 - 1/4) = 6/8 * 0.75 = 0.5625
assert!((gr - 0.5625).abs() < 1e-4, "gr={gr}");
}
}

76
crates/headroom-dsp/src/delay.rs Normal file
View file

@ -0,0 +1,76 @@
//! A simple fixed-length sample delay line.
/// FIFO sample delay of fixed length.
///
/// `push_pop(x)` writes `x` and returns the sample that was written
/// `len` calls ago (initialized to zero).
pub struct DelayLine {
buf: Vec<f32>,
write_idx: usize,
}
impl DelayLine {
/// Construct a new delay line of `samples` length.
///
/// Lengths of 0 are clamped to 1 so the type always behaves like a
/// one-sample identity at minimum.
#[must_use]
pub fn new(samples: usize) -> Self {
Self {
buf: vec![0.0; samples.max(1)],
write_idx: 0,
}
}
/// Effective delay in samples.
#[must_use]
pub fn len(&self) -> usize {
self.buf.len()
}
/// Always false.
#[must_use]
pub fn is_empty(&self) -> bool {
false
}
/// Write `x` and return the sample written `len` calls ago.
pub fn push_pop(&mut self, x: f32) -> f32 {
let out = self.buf[self.write_idx];
self.buf[self.write_idx] = x;
self.write_idx = (self.write_idx + 1) % self.buf.len();
out
}
/// Clear the delay line to silence.
pub fn reset(&mut self) {
for v in &mut self.buf {
*v = 0.0;
}
self.write_idx = 0;
}
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn delays_exactly_n_samples() {
let mut d = DelayLine::new(4);
let expected = [0.0, 0.0, 0.0, 0.0, 1.0, 2.0, 3.0, 4.0];
let inputs = [1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0];
for (i, &x) in inputs.iter().enumerate() {
let y = d.push_pop(x);
assert!((y - expected[i]).abs() < 1e-9, "i={i} y={y}");
}
}
#[test]
fn zero_length_clamps_to_one() {
let mut d = DelayLine::new(0);
assert_eq!(d.len(), 1);
assert_eq!(d.push_pop(1.0), 0.0);
assert_eq!(d.push_pop(2.0), 1.0);
}
}

98
crates/headroom-dsp/src/envelope.rs Normal file
View file

@ -0,0 +1,98 @@
//! Exponential attack/release envelope follower.
use crate::util::time_to_alpha;
/// One-pole smoother with separate attack and release coefficients.
pub struct AttackRelease {
attack_alpha: f32,
release_alpha: f32,
state: f32,
}
impl AttackRelease {
/// Construct from times in milliseconds and a sample rate (Hz).
#[must_use]
pub fn new(attack_ms: f32, release_ms: f32, sample_rate: f32) -> Self {
Self {
attack_alpha: time_to_alpha(attack_ms, sample_rate),
release_alpha: time_to_alpha(release_ms, sample_rate),
state: 0.0,
}
}
/// Update the coefficients (e.g. after a sample-rate change).
pub fn set_times(&mut self, attack_ms: f32, release_ms: f32, sample_rate: f32) {
self.attack_alpha = time_to_alpha(attack_ms, sample_rate);
self.release_alpha = time_to_alpha(release_ms, sample_rate);
}
/// Peak-detector mode: attack on rising input, release on falling.
/// Typical use: envelope detector for compressors.
pub fn process_peak(&mut self, target: f32) -> f32 {
if target > self.state {
self.state += self.attack_alpha * (target - self.state);
} else {
self.state += self.release_alpha * (target - self.state);
}
self.state
}
/// Gain-follower mode: attack on falling input (gain dropping),
/// release on rising input (gain recovering toward unity). The
/// inverse direction from [`process_peak`](Self::process_peak).
pub fn process_gain(&mut self, target: f32) -> f32 {
if target < self.state {
self.state += self.attack_alpha * (target - self.state);
} else {
self.state += self.release_alpha * (target - self.state);
}
self.state
}
/// Current state.
#[must_use]
pub fn state(&self) -> f32 {
self.state
}
/// Force the state to a given value.
pub fn reset(&mut self, value: f32) {
self.state = value;
}
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn peak_mode_attacks_fast_releases_slow() {
let fs = 48_000.0;
let mut env = AttackRelease::new(0.1, 100.0, fs);
// Drive to 1.0 and let it settle.
for _ in 0..100 {
env.process_peak(1.0);
}
assert!(env.state() > 0.99);
// Drop input to 0.0 and verify slow decay.
env.process_peak(0.0);
assert!(env.state() > 0.999);
for _ in 0..10 {
env.process_peak(0.0);
}
// Still well above zero on the release time scale.
assert!(env.state() > 0.8);
}
#[test]
fn gain_mode_attacks_on_drop() {
let fs = 48_000.0;
let mut env = AttackRelease::new(0.1, 100.0, fs);
env.reset(1.0);
// Demand a gain drop. Should snap down quickly.
for _ in 0..100 {
env.process_gain(0.5);
}
assert!(env.state() < 0.51);
}
}

25
crates/headroom-dsp/src/lib.rs Normal file
View file

@ -0,0 +1,25 @@
//! DSP kernels for Headroom.
//!
//! The contract: every `process_*` method on the public types is
//! allocation-free and bounded-time. Construction (`new`) allocates and
//! is not realtime-safe — do it ahead of time.
//!
//! See `PLAN.md` §3 for the role each kernel plays in the daemon.
#![forbid(unsafe_code)]
#![warn(missing_docs)]
mod compressor;
mod delay;
mod envelope;
mod limiter;
mod oversample;
mod sliding_max;
pub mod util;
pub use compressor::{Compressor, CompressorConfig, Detector};
pub use delay::DelayLine;
pub use envelope::AttackRelease;
pub use limiter::{Limiter, LimiterConfig, SoftTierConfig};
pub use oversample::{design_lowpass_blackman, PolyphaseDownsampler, PolyphaseUpsampler};
pub use sliding_max::SlidingMaxBuffer;

850
crates/headroom-dsp/src/limiter.rs Normal file
View file

@ -0,0 +1,850 @@
//! Two-tier true-peak limiter.
//!
//! Architecture (per channel, stereo-linked gain):
//!
//! ```text
//! input ─► upsample ─► delay ─► × gain ─► clamp ─► downsample ─► clamp ─► out
//! │ ▲
//! └─► peak ─────┤
//! │ │
//! ▼ │
//! ┌────────┐ │
//! │ soft │── × │ smooth attack/release
//! │ ceil │ │ target = program_lufs + max_psr_db
//! └────────┘ │
//! ┌────────┐ │
//! │ hard │── × ┘ instant attack + hold + release
//! │ ceil │ target = ceiling_dbtp (e.g. −0.1)
//! └────────┘
//! ```
//!
//! The **hard tier** enforces the absolute output contract: a
//! configurable ceiling (default `−0.1 dBTP`) with full inter-sample
//! peak handling. Its gain has instant attack, a brief hold, and an
//! exponential release. Two defensive `clamp` stages downstream
//! guarantee the contract numerically — the envelope can misbehave and
//! the contract still holds.
//!
//! The **soft tier** sits in parallel. Its target ceiling is *dynamic*:
//! `program_loudness_lufs + soft.max_psr_db`. It uses a smooth
//! attack/release envelope (musical, not slappy) and pulls transients
//! down to a comfortable peak-to-loudness ratio *before* they ever
//! threaten the hard ceiling. Listeners hear "loud-but-not-shocking"
//! transients instead of bricks landing exactly at the ceiling.
//!
//! The two tiers share the upsampler, downsampler, delay line, and
//! sliding peak buffer. The cost of the soft tier is one extra
//! envelope evaluation per oversampled sample — no additional latency,
//! no additional FIR work.
//!
//! When no program loudness has been provided (typical at startup, or
//! before the AGC's first `ebur128` window completes), the soft tier
//! falls back to a static ceiling. When the soft tier is disabled
//! entirely (`LimiterConfig::soft = None`), the limiter behaves as a
//! pure brickwall — see the `transparent` profile.
use crate::delay::DelayLine;
use crate::envelope::AttackRelease;
use crate::oversample::{design_lowpass_blackman, PolyphaseDownsampler, PolyphaseUpsampler};
use crate::sliding_max::SlidingMaxBuffer;
use crate::util::{db_to_lin, lin_to_db, time_to_alpha};
/// Soft-tier configuration.
///
/// The soft tier targets a *dynamic* ceiling computed as
/// `program_loudness_lufs + max_psr_db`. It is responsible for the
/// listening experience — keeping the peak-to-loudness ratio bounded
/// — without acting as a safety contract (that's the hard tier's job).
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct SoftTierConfig {
/// Maximum allowed peak-to-shortterm-loudness ratio in dB.
/// Effective ceiling becomes `program_lufs + max_psr_db`.
pub max_psr_db: f32,
/// Fallback ceiling in dBTP used when no program loudness has been
/// supplied (e.g. during startup before the AGC has measured the
/// first short-term window).
pub static_ceiling_dbtp: f32,
/// Attack time in ms (smooth, not instant).
pub attack_ms: f32,
/// Release time in ms.
pub release_ms: f32,
}
impl Default for SoftTierConfig {
fn default() -> Self {
Self {
max_psr_db: 14.0,
static_ceiling_dbtp: -6.0,
attack_ms: 5.0,
release_ms: 200.0,
}
}
}
impl SoftTierConfig {
/// Sanitize: positive ceilings clamp to 0, non-finite or negative
/// times clamp to small positives.
#[must_use]
pub fn sanitized(mut self) -> Self {
if self.static_ceiling_dbtp > 0.0 {
self.static_ceiling_dbtp = 0.0;
}
if !self.max_psr_db.is_finite() || self.max_psr_db < 0.0 {
self.max_psr_db = 0.0;
}
if self.attack_ms < 0.0 || !self.attack_ms.is_finite() {
self.attack_ms = 0.0;
}
if self.release_ms < 0.0 || !self.release_ms.is_finite() {
self.release_ms = 0.0;
}
self
}
}
/// Configurable limiter parameters.
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct LimiterConfig {
/// Hard-tier output ceiling in dBTP. Must be `<= 0.0`.
pub ceiling_dbtp: f32,
/// Lookahead time in milliseconds. Sets the delay-line length and
/// the size of the peak-detector sliding window. Shared by both
/// tiers.
pub lookahead_ms: f32,
/// Hard-tier exponential release (toward unity gain) in ms.
pub release_ms: f32,
/// Hard-tier hold time after a gain reduction before release
/// begins, in ms.
pub hold_ms: f32,
/// Oversampling factor. 1 disables ISP detection; 4 is the
/// BS.1770-4 reference.
pub oversample: usize,
/// Number of FIR taps used by the oversampling filter (odd).
pub fir_taps: usize,
/// Soft-tier configuration. `None` disables the soft tier and the
/// limiter behaves as a pure brickwall.
pub soft: Option<SoftTierConfig>,
}
impl Default for LimiterConfig {
fn default() -> Self {
Self {
ceiling_dbtp: -0.1,
lookahead_ms: 2.0,
release_ms: 80.0,
hold_ms: 5.0,
oversample: 4,
fir_taps: 31,
soft: Some(SoftTierConfig::default()),
}
}
}
impl LimiterConfig {
/// Sanitize a user-supplied configuration: clamp ceiling,
/// oversample factor, ensure odd FIR length, sanitize the soft
/// tier if present.
#[must_use]
pub fn sanitized(mut self) -> Self {
if self.ceiling_dbtp > 0.0 {
self.ceiling_dbtp = 0.0;
}
self.oversample = self.oversample.clamp(1, 8);
if self.fir_taps < 5 {
self.fir_taps = 5;
}
if self.fir_taps % 2 == 0 {
self.fir_taps += 1;
}
if let Some(soft) = self.soft {
self.soft = Some(soft.sanitized());
}
self
}
/// Convenience: brickwall only (no soft tier).
#[must_use]
pub fn brickwall_only() -> Self {
Self {
soft: None,
..Self::default()
}
}
}
const MAX_OVERSAMPLE: usize = 8;
/// Two-tier feed-forward true-peak limiter.
pub struct Limiter {
cfg: LimiterConfig,
ceiling_lin: f32,
os: usize,
// Per-channel oversampler / downsampler / delay-line paths.
up_l: PolyphaseUpsampler,
up_r: PolyphaseUpsampler,
down_l: PolyphaseDownsampler,
down_r: PolyphaseDownsampler,
delay_l: DelayLine,
delay_r: DelayLine,
/// Sliding-window peak (oversampled domain), shared across
/// channels and tiers.
peak_buf: SlidingMaxBuffer,
// ---- Hard tier state (instant attack + hold + release) ----
hard_gain: f32,
hold_remaining: u32,
hold_samples_os: u32,
hard_release_alpha: f32,
// ---- Soft tier state (smooth envelope) ----
soft_envelope: Option<AttackRelease>,
soft_max_psr_db: f32,
soft_static_ceiling_lin: f32,
program_loudness_lufs: Option<f32>,
/// Effective soft ceiling in linear gain, recomputed whenever the
/// program loudness changes (or once at startup).
soft_ceiling_lin: f32,
// Scratch buffers (sized for the maximum supported oversample
// factor).
up_buf_l: [f32; MAX_OVERSAMPLE],
up_buf_r: [f32; MAX_OVERSAMPLE],
gained_buf_l: [f32; MAX_OVERSAMPLE],
gained_buf_r: [f32; MAX_OVERSAMPLE],
// Telemetry (sampled per frame).
last_peak_lin: f32,
last_gr_db: f32,
last_soft_gr_db: f32,
last_hard_gr_db: f32,
}
impl Limiter {
/// Construct from a configuration and input sample rate.
///
/// Allocates the FIR coefficients, polyphase tables, and delay
/// buffers. Not realtime-safe.
#[must_use]
pub fn new(cfg: LimiterConfig, sample_rate: f32) -> Self {
let cfg = cfg.sanitized();
let os = cfg.oversample;
let lowpass = if os > 1 {
design_lowpass_blackman(cfg.fir_taps, 0.45 / os as f32)
} else {
vec![1.0]
};
let os_rate = sample_rate * os as f32;
let lookahead_samples_os = (cfg.lookahead_ms * 1e-3 * os_rate).round() as usize;
let lookahead_samples_os = lookahead_samples_os.max(1);
let hold_samples_os = (cfg.hold_ms * 1e-3 * os_rate).round() as u32;
let hard_release_alpha = time_to_alpha(cfg.release_ms, os_rate);
let ceiling_lin = db_to_lin(cfg.ceiling_dbtp);
let (soft_envelope, soft_max_psr_db, soft_static_ceiling_lin, soft_ceiling_lin) =
if let Some(soft) = cfg.soft {
let env = AttackRelease::new(soft.attack_ms, soft.release_ms, os_rate);
let static_ceiling_lin = db_to_lin(soft.static_ceiling_dbtp);
(
Some(env),
soft.max_psr_db,
static_ceiling_lin,
static_ceiling_lin,
)
} else {
(None, 0.0, 1.0, 1.0)
};
let mut me = Self {
cfg,
ceiling_lin,
os,
up_l: PolyphaseUpsampler::new(os, &lowpass),
up_r: PolyphaseUpsampler::new(os, &lowpass),
down_l: PolyphaseDownsampler::new(os, &lowpass),
down_r: PolyphaseDownsampler::new(os, &lowpass),
delay_l: DelayLine::new(lookahead_samples_os),
delay_r: DelayLine::new(lookahead_samples_os),
peak_buf: SlidingMaxBuffer::new(lookahead_samples_os),
hard_gain: 1.0,
hold_remaining: 0,
hold_samples_os,
hard_release_alpha,
soft_envelope,
soft_max_psr_db,
soft_static_ceiling_lin,
program_loudness_lufs: None,
soft_ceiling_lin,
up_buf_l: [0.0; MAX_OVERSAMPLE],
up_buf_r: [0.0; MAX_OVERSAMPLE],
gained_buf_l: [0.0; MAX_OVERSAMPLE],
gained_buf_r: [0.0; MAX_OVERSAMPLE],
last_peak_lin: 0.0,
last_gr_db: 0.0,
last_soft_gr_db: 0.0,
last_hard_gr_db: 0.0,
};
// Seed soft envelope to unity so we don't start with phantom
// gain reduction during the first frames.
if let Some(env) = &mut me.soft_envelope {
env.reset(1.0);
}
me
}
/// Active configuration.
#[must_use]
pub fn config(&self) -> LimiterConfig {
self.cfg
}
/// Hard-tier output ceiling in dBTP.
#[must_use]
pub fn ceiling_dbtp(&self) -> f32 {
self.cfg.ceiling_dbtp
}
/// Most recent total gain reduction in dB (negative when limiting).
/// This is the *applied* reduction: `min(soft_gain, hard_gain)`.
#[must_use]
pub fn gain_reduction_db(&self) -> f32 {
self.last_gr_db
}
/// Most recent soft-tier gain reduction in dB.
#[must_use]
pub fn soft_gain_reduction_db(&self) -> f32 {
self.last_soft_gr_db
}
/// Most recent hard-tier gain reduction in dB.
///
/// A non-zero value here indicates the soft tier did not keep the
/// signal under the absolute ceiling and the brickwall engaged.
/// Routinely non-zero values for benign material suggest the soft
/// tier is under-configured (too high `max_psr_db`, too slow
/// attack, or the lookahead is too short for the chosen attack).
#[must_use]
pub fn hard_gain_reduction_db(&self) -> f32 {
self.last_hard_gr_db
}
/// Most recent observed true-peak in dBTP.
#[must_use]
pub fn true_peak_dbtp(&self) -> f32 {
lin_to_db(self.last_peak_lin.max(1e-20))
}
/// Effective soft ceiling currently in use, in dBTP.
///
/// Equals `program_loudness_lufs + soft.max_psr_db` when both are
/// known, otherwise the configured `soft.static_ceiling_dbtp`.
/// Returns `None` if the soft tier is disabled.
#[must_use]
pub fn effective_soft_ceiling_dbtp(&self) -> Option<f32> {
self.cfg.soft.map(|_| lin_to_db(self.soft_ceiling_lin))
}
/// Update the program loudness used to compute the dynamic soft
/// ceiling. Typically called by the AGC at its tick rate with the
/// short-term BS.1770 loudness; non-finite values are ignored.
pub fn set_program_loudness_lufs(&mut self, lufs: f32) {
if !lufs.is_finite() {
return;
}
self.program_loudness_lufs = Some(lufs);
self.recompute_soft_ceiling();
}
/// Forget the program loudness; soft tier falls back to its static
/// ceiling. Useful when the AGC stalls or is reset.
pub fn clear_program_loudness(&mut self) {
self.program_loudness_lufs = None;
self.recompute_soft_ceiling();
}
fn recompute_soft_ceiling(&mut self) {
self.soft_ceiling_lin = match (self.cfg.soft, self.program_loudness_lufs) {
(Some(_), Some(lufs)) => {
let dynamic_dbtp = (lufs + self.soft_max_psr_db).min(0.0);
db_to_lin(dynamic_dbtp)
}
(Some(_), None) => self.soft_static_ceiling_lin,
(None, _) => 1.0,
};
}
/// Process one stereo frame.
///
/// Allocation-free. Returns `(left, right)` guaranteed to lie
/// within `±ceiling_dbtp` (the hard contract).
pub fn process_frame(&mut self, left: f32, right: f32) -> (f32, f32) {
// Sanitize NaN / Inf to zero defensively; never propagate
// garbage into the limiter state.
let left = if left.is_finite() { left } else { 0.0 };
let right = if right.is_finite() { right } else { 0.0 };
self.up_l.process(left, &mut self.up_buf_l[..self.os]);
self.up_r.process(right, &mut self.up_buf_r[..self.os]);
let mut frame_peak = 0.0_f32;
let mut min_soft_gain = 1.0_f32;
let mut min_total_gain = 1.0_f32;
for k in 0..self.os {
let s_l = self.up_buf_l[k];
let s_r = self.up_buf_r[k];
let peak = s_l.abs().max(s_r.abs());
frame_peak = frame_peak.max(peak);
let window_peak = self.peak_buf.push_and_max(peak);
// ---- Soft tier --------------------------------------
let soft_gain = if let Some(env) = &mut self.soft_envelope {
let target = if window_peak > self.soft_ceiling_lin && window_peak > 1e-20 {
self.soft_ceiling_lin / window_peak
} else {
1.0
};
env.process_gain(target)
} else {
1.0
};
if soft_gain < min_soft_gain {
min_soft_gain = soft_gain;
}
// ---- Hard tier --------------------------------------
// The hard tier defends the ceiling, but it shouldn't do
// redundant work when the soft tier already handles the
// peak. So compute the predicted peak *after* the soft
// tier acts, then size the hard gain against that.
//
// `predicted_post_soft` takes the max of:
// - `immediate`: the peak with the *current* soft gain
// applied (safe if soft hasn't ramped yet)
// - `asymptotic`: the peak after the soft tier converges
// to its target (the steady-state)
// The max is the more conservative (larger) prediction.
let predicted_post_soft = if self.soft_envelope.is_some() {
let asymptotic = window_peak.min(self.soft_ceiling_lin);
let immediate = window_peak * soft_gain;
asymptotic.max(immediate)
} else {
window_peak
};
let hard_target =
if predicted_post_soft > self.ceiling_lin && predicted_post_soft > 1e-20 {
self.ceiling_lin / predicted_post_soft
} else {
1.0
};
if hard_target < self.hard_gain {
self.hard_gain = hard_target;
self.hold_remaining = self.hold_samples_os;
} else if self.hold_remaining > 0 {
self.hold_remaining -= 1;
} else {
self.hard_gain += self.hard_release_alpha * (hard_target - self.hard_gain);
if self.hard_gain > hard_target {
self.hard_gain = hard_target;
}
}
// ---- Combine ----------------------------------------
let total_gain = soft_gain.min(self.hard_gain);
if total_gain < min_total_gain {
min_total_gain = total_gain;
}
let d_l = self.delay_l.push_pop(s_l);
let d_r = self.delay_r.push_pop(s_r);
let mut out_l = d_l * total_gain;
let mut out_r = d_r * total_gain;
// Defense-in-depth #1: brickwall clip in the oversampled
// domain. Prevents extreme overshoots from passing into
// the downsampler.
out_l = out_l.clamp(-self.ceiling_lin, self.ceiling_lin);
out_r = out_r.clamp(-self.ceiling_lin, self.ceiling_lin);
self.gained_buf_l[k] = out_l;
self.gained_buf_r[k] = out_r;
}
let mut out_l = self.down_l.process(&self.gained_buf_l[..self.os]);
let mut out_r = self.down_r.process(&self.gained_buf_r[..self.os]);
// Defense-in-depth #2: brickwall clip at the input sample
// rate, after downsampling. Guards against FIR-induced ringing
// nudging the output above the ceiling in the downsampled
// domain.
out_l = out_l.clamp(-self.ceiling_lin, self.ceiling_lin);
out_r = out_r.clamp(-self.ceiling_lin, self.ceiling_lin);
self.last_peak_lin = frame_peak;
self.last_soft_gr_db = lin_to_db(min_soft_gain.max(1e-12));
self.last_hard_gr_db = lin_to_db(self.hard_gain.max(1e-12));
self.last_gr_db = lin_to_db(min_total_gain.max(1e-12));
(out_l, out_r)
}
/// Process an interleaved stereo buffer in place.
pub fn process_interleaved_stereo(&mut self, buf: &mut [f32]) {
debug_assert!(buf.len() % 2 == 0);
for frame in buf.chunks_exact_mut(2) {
let (l, r) = self.process_frame(frame[0], frame[1]);
frame[0] = l;
frame[1] = r;
}
}
/// Reset all internal state. Program loudness is also cleared.
pub fn reset(&mut self) {
self.up_l.reset();
self.up_r.reset();
self.down_l.reset();
self.down_r.reset();
self.delay_l.reset();
self.delay_r.reset();
self.peak_buf.reset();
self.hard_gain = 1.0;
self.hold_remaining = 0;
if let Some(env) = &mut self.soft_envelope {
env.reset(1.0);
}
self.program_loudness_lufs = None;
self.recompute_soft_ceiling();
self.last_peak_lin = 0.0;
self.last_gr_db = 0.0;
self.last_soft_gr_db = 0.0;
self.last_hard_gr_db = 0.0;
}
}
#[cfg(test)]
mod tests {
use super::*;
use std::f32::consts::PI;
fn run_sine(
limiter: &mut Limiter,
freq: f32,
amp_db: f32,
samples: usize,
sr: f32,
) -> Vec<f32> {
let amp = db_to_lin(amp_db);
let mut out = Vec::with_capacity(samples * 2);
for n in 0..samples {
let t = n as f32 / sr;
let s = amp * (2.0 * PI * freq * t).sin();
let (l, r) = limiter.process_frame(s, s);
out.push(l);
out.push(r);
}
out
}
// ----------------------------------------------------------------
// Hard-tier contract: holds with or without the soft tier present.
// ----------------------------------------------------------------
#[test]
fn passes_signal_below_both_ceilings_unchanged() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
// -18 dBFS is below the default static soft ceiling of -6 dBTP
// and the hard ceiling. Neither tier should engage.
let out = run_sine(&mut l, 440.0, -18.0, 4_800, sr);
let max_abs = out.iter().skip(1_000).fold(0.0_f32, |a, &b| a.max(b.abs()));
let max_db = lin_to_db(max_abs);
assert!(
(max_db - (-18.0)).abs() < 0.5,
"expected ~-18 dB, got {max_db}"
);
assert!(
l.gain_reduction_db().abs() < 0.5,
"expected ~0 GR, got {}",
l.gain_reduction_db()
);
}
#[test]
fn enforces_hard_ceiling_on_hot_signal_with_soft_tier() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
let out = run_sine(&mut l, 440.0, 6.0, 9_600, sr);
let ceiling_lin = db_to_lin(-0.1);
let max_abs = out
.iter()
.skip(2_000)
.fold(0.0_f32, |a, &b| a.max(b.abs()));
assert!(
max_abs <= ceiling_lin + 1e-6,
"above hard ceiling: max_abs={max_abs}, ceiling_lin={ceiling_lin}"
);
}
#[test]
fn enforces_hard_ceiling_on_intersample_peak_with_soft_tier() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
let ceiling_lin = db_to_lin(-0.1);
let mut max_abs = 0.0_f32;
let mut sign = 1.0_f32;
let amp = 0.95_f32;
for n in 0..9_600 {
let s = sign * amp;
sign = -sign;
let (lo, ro) = l.process_frame(s, s);
if n > 1_500 {
max_abs = max_abs.max(lo.abs()).max(ro.abs());
}
}
assert!(
max_abs <= ceiling_lin + 1e-6,
"ISP: above hard ceiling: max_abs={max_abs}, ceiling_lin={ceiling_lin}"
);
}
#[test]
fn enforces_hard_ceiling_on_transient_impulse_with_soft_tier() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
let ceiling_lin = db_to_lin(-0.1);
let mut max_abs = 0.0_f32;
for n in 0..4_800_usize {
let s = if n == 1_000 { 4.0 } else { 0.0 };
let (lo, ro) = l.process_frame(s, s);
max_abs = max_abs.max(lo.abs()).max(ro.abs());
}
assert!(
max_abs <= ceiling_lin + 1e-6,
"impulse: above hard ceiling: max_abs={max_abs}, ceiling_lin={ceiling_lin}"
);
}
#[test]
fn brickwall_only_skips_soft_tier_entirely() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::brickwall_only(), sr);
assert!(l.effective_soft_ceiling_dbtp().is_none());
// Drive a hot signal; brickwall must still hold.
let out = run_sine(&mut l, 440.0, 6.0, 4_800, sr);
let ceiling_lin = db_to_lin(-0.1);
let max_abs = out.iter().skip(800).fold(0.0_f32, |a, &b| a.max(b.abs()));
assert!(max_abs <= ceiling_lin + 1e-6);
// No soft gain reduction should ever have been recorded.
assert!(l.soft_gain_reduction_db().abs() < 1e-6);
}
// ----------------------------------------------------------------
// Soft tier: static fallback ceiling
// ----------------------------------------------------------------
#[test]
fn soft_tier_static_ceiling_engages_before_hard() {
let sr = 48_000.0;
// Static soft ceiling at -6 dBTP, attack short enough to
// settle inside the lookahead.
let cfg = LimiterConfig {
lookahead_ms: 5.0,
soft: Some(SoftTierConfig {
static_ceiling_dbtp: -6.0,
attack_ms: 1.0,
release_ms: 100.0,
..SoftTierConfig::default()
}),
..LimiterConfig::default()
};
let mut l = Limiter::new(cfg, sr);
// Drive a +6 dB sine — well above the soft ceiling.
let out = run_sine(&mut l, 440.0, 6.0, 9_600, sr);
// Output should sit near the soft ceiling, well below hard.
let soft_ceiling_lin = db_to_lin(-6.0);
let max_abs = out
.iter()
.skip(2_000)
.fold(0.0_f32, |a, &b| a.max(b.abs()));
// Allow small overshoot during soft attack (gain hasn't fully
// settled when the peak arrives), but it must be well under
// the hard ceiling.
assert!(
max_abs <= soft_ceiling_lin * 1.1,
"output above soft ceiling: max_abs={max_abs}, soft_lin={soft_ceiling_lin}"
);
// Soft tier should report meaningful GR; hard tier ideally
// does very little once the soft tier has settled.
assert!(
l.soft_gain_reduction_db() < -3.0,
"soft GR too small: {}",
l.soft_gain_reduction_db()
);
}
// ----------------------------------------------------------------
// Soft tier: dynamic ceiling from program loudness
// ----------------------------------------------------------------
#[test]
fn dynamic_ceiling_tracks_program_loudness() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
// Default max_psr_db = 14.
l.set_program_loudness_lufs(-18.0);
let dyn_ceiling = l.effective_soft_ceiling_dbtp().expect("soft tier active");
assert!(
(dyn_ceiling - (-4.0)).abs() < 1e-3,
"expected -4 dBTP, got {dyn_ceiling}"
);
// Move the program louder; ceiling rises (and clamps at 0).
l.set_program_loudness_lufs(-2.0);
let dyn_ceiling = l.effective_soft_ceiling_dbtp().unwrap();
assert!(
(-0.1..=0.0).contains(&dyn_ceiling),
"expected clamp near 0 dBTP, got {dyn_ceiling}"
);
// Clear it; falls back to static.
l.clear_program_loudness();
let fallback = l.effective_soft_ceiling_dbtp().unwrap();
assert!(
(fallback - (-6.0)).abs() < 1e-3,
"expected static -6 dBTP, got {fallback}"
);
}
#[test]
fn dynamic_ceiling_bounds_psr_on_hot_transient() {
let sr = 48_000.0;
// Long lookahead and fast soft attack so the soft tier
// demonstrably catches the transient before the hard tier
// needs to.
let cfg = LimiterConfig {
lookahead_ms: 5.0,
soft: Some(SoftTierConfig {
max_psr_db: 14.0,
static_ceiling_dbtp: -6.0,
attack_ms: 1.0,
release_ms: 100.0,
}),
..LimiterConfig::default()
};
let mut l = Limiter::new(cfg, sr);
l.set_program_loudness_lufs(-18.0);
// Expected dynamic ceiling: -18 + 14 = -4 dBTP ≈ 0.631 lin.
let dyn_ceil_lin = db_to_lin(-4.0);
// Slam a +6 dBFS impulse.
let mut max_after = 0.0_f32;
for n in 0..4_800_usize {
let s = if n == 800 { db_to_lin(6.0) } else { 0.0 };
let (lo, _) = l.process_frame(s, s);
if n > 700 {
max_after = max_after.max(lo.abs());
}
}
// Output should be at or below the dynamic soft ceiling with
// a small ringing margin. Critically, the hard tier should
// *not* be the thing that catches it — its GR should be small.
assert!(
max_after <= dyn_ceil_lin * 1.15,
"soft tier didn't bound the transient: max={max_after}, dyn_ceil={dyn_ceil_lin}"
);
// The hard tier may snap briefly at peak entry (soft envelope
// hasn't ramped yet), then take its release time to recover.
// We don't require zero hard engagement here — only that it
// isn't doing the majority of the work.
assert!(
l.hard_gain_reduction_db().abs() < 4.0,
"hard tier engaged unreasonably: {}",
l.hard_gain_reduction_db()
);
}
// ----------------------------------------------------------------
// Misc
// ----------------------------------------------------------------
#[test]
fn nan_inputs_do_not_propagate_with_soft_tier() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
for _ in 0..1_000 {
let (lo, ro) = l.process_frame(f32::NAN, f32::INFINITY);
assert!(lo.is_finite() && ro.is_finite());
}
}
#[test]
fn ceiling_clamps_positive_config_to_zero() {
let cfg = LimiterConfig {
ceiling_dbtp: 3.0,
..LimiterConfig::default()
}
.sanitized();
assert_eq!(cfg.ceiling_dbtp, 0.0);
}
#[test]
fn set_program_loudness_ignores_non_finite() {
let sr = 48_000.0;
let mut l = Limiter::new(LimiterConfig::default(), sr);
// Establish a baseline.
l.set_program_loudness_lufs(-20.0);
let baseline = l.effective_soft_ceiling_dbtp().unwrap();
// NaN / Inf should be ignored.
l.set_program_loudness_lufs(f32::NAN);
assert_eq!(l.effective_soft_ceiling_dbtp().unwrap(), baseline);
l.set_program_loudness_lufs(f32::INFINITY);
assert_eq!(l.effective_soft_ceiling_dbtp().unwrap(), baseline);
}
#[test]
fn soft_tier_reduces_perceived_peak_to_loudness_ratio() {
// The whole point of the soft tier: a transient on top of a
// quieter program should NOT come out near the hard ceiling.
let sr = 48_000.0;
let cfg = LimiterConfig {
lookahead_ms: 5.0,
soft: Some(SoftTierConfig {
max_psr_db: 12.0,
static_ceiling_dbtp: -8.0,
attack_ms: 1.0,
release_ms: 100.0,
}),
..LimiterConfig::default()
};
let mut brickwall = Limiter::new(LimiterConfig::brickwall_only(), sr);
let mut two_tier = Limiter::new(cfg, sr);
two_tier.set_program_loudness_lufs(-20.0);
let mut bw_peak = 0.0_f32;
let mut tt_peak = 0.0_f32;
for n in 0..4_800_usize {
// Quiet program with a single big spike.
let s = if n == 1_200 { db_to_lin(3.0) } else { 0.01 };
let (lo_bw, _) = brickwall.process_frame(s, s);
let (lo_tt, _) = two_tier.process_frame(s, s);
if n > 1_000 {
bw_peak = bw_peak.max(lo_bw.abs());
tt_peak = tt_peak.max(lo_tt.abs());
}
}
// Brickwall lets the spike through near the hard ceiling.
// Two-tier holds it much lower.
assert!(
tt_peak < bw_peak * 0.6,
"soft tier did not meaningfully reduce peak: bw={bw_peak}, tt={tt_peak}"
);
}
}

230
crates/headroom-dsp/src/oversample.rs Normal file
View file

@ -0,0 +1,230 @@
//! Polyphase FIR up/downsamplers.
//!
//! Used by the true-peak limiter to detect inter-sample peaks via
//! oversampled-domain peak detection (per ITU-R BS.1770-4).
use std::f32::consts::PI;
/// Design a Blackman-windowed sinc lowpass FIR.
///
/// * `taps` — filter length. Odd values give a linear-phase filter
/// with an exact group delay of `(taps - 1) / 2` samples.
/// * `fc` — normalized cutoff in `0.0..0.5` (fraction of sample rate).
///
/// Coefficients are normalized for unity DC gain. Suitable as the
/// prototype lowpass for `M`-times oversampling at `fc = 0.5 / M`
/// (slightly below Nyquist of the lower rate).
#[must_use]
pub fn design_lowpass_blackman(taps: usize, fc: f32) -> Vec<f32> {
let taps = taps.max(1);
let m = (taps as f32 - 1.0).max(1.0);
let mut h = vec![0.0_f32; taps];
let mut sum = 0.0_f32;
for (n, h_n) in h.iter_mut().enumerate() {
let x = n as f32 - m / 2.0;
let sinc = if x.abs() < 1e-9 {
2.0 * fc
} else {
(2.0 * PI * fc * x).sin() / (PI * x)
};
let w = 0.42 - 0.5 * (2.0 * PI * n as f32 / m).cos() + 0.08 * (4.0 * PI * n as f32 / m).cos();
*h_n = sinc * w;
sum += *h_n;
}
if sum.abs() > 1e-12 {
for v in &mut h {
*v /= sum;
}
}
h
}
/// Polyphase FIR upsampler.
///
/// One input sample produces `factor` output samples. Coefficients are
/// pre-scaled by `factor` so the output's DC gain equals the input.
pub struct PolyphaseUpsampler {
factor: usize,
taps_per_phase: usize,
/// `phases[j * taps_per_phase + p] = h[p * factor + j] * factor`.
phases: Vec<f32>,
history: Vec<f32>,
write_idx: usize,
}
impl PolyphaseUpsampler {
/// Construct from a prototype lowpass and an upsample `factor`.
#[must_use]
pub fn new(factor: usize, fir_taps: &[f32]) -> Self {
let factor = factor.max(1);
let taps_per_phase = fir_taps.len().div_ceil(factor);
let mut phases = vec![0.0_f32; factor * taps_per_phase];
for (n, &h) in fir_taps.iter().enumerate() {
let j = n % factor;
let p = n / factor;
phases[j * taps_per_phase + p] = h * factor as f32;
}
Self {
factor,
taps_per_phase,
phases,
history: vec![0.0_f32; taps_per_phase.max(1)],
write_idx: 0,
}
}
/// Upsample factor.
#[must_use]
pub fn factor(&self) -> usize {
self.factor
}
/// Push one input sample and emit `factor` output samples into
/// `out[..factor]`. `out` must have length `>= factor`.
pub fn process(&mut self, x: f32, out: &mut [f32]) {
debug_assert!(out.len() >= self.factor);
let len = self.history.len();
self.history[self.write_idx] = x;
let just_written = self.write_idx;
self.write_idx = (self.write_idx + 1) % len;
for (j, slot) in out.iter_mut().take(self.factor).enumerate() {
let phase = &self.phases[j * self.taps_per_phase..(j + 1) * self.taps_per_phase];
let mut acc = 0.0_f32;
for (p, &h) in phase.iter().enumerate() {
let idx = (just_written + len - p) % len;
acc += h * self.history[idx];
}
*slot = acc;
}
}
/// Clear filter history.
pub fn reset(&mut self) {
for v in &mut self.history {
*v = 0.0;
}
self.write_idx = 0;
}
}
/// FIR downsampler. Takes `factor` input samples and emits one output.
///
/// Uses the same prototype lowpass coefficients as the upsampler.
/// Implementation is straightforward (no polyphase split) — for our
/// filter length the savings are modest and code clarity wins.
pub struct PolyphaseDownsampler {
factor: usize,
taps: Vec<f32>,
history: Vec<f32>,
write_idx: usize,
}
impl PolyphaseDownsampler {
/// Construct from a prototype lowpass and a downsample `factor`.
#[must_use]
pub fn new(factor: usize, fir_taps: &[f32]) -> Self {
let factor = factor.max(1);
Self {
factor,
taps: fir_taps.to_vec(),
history: vec![0.0_f32; fir_taps.len().max(1)],
write_idx: 0,
}
}
/// Downsample factor.
#[must_use]
pub fn factor(&self) -> usize {
self.factor
}
/// Push `factor` input samples and return one filtered output.
pub fn process(&mut self, ins: &[f32]) -> f32 {
debug_assert_eq!(ins.len(), self.factor);
let len = self.history.len();
for &x in ins {
self.history[self.write_idx] = x;
self.write_idx = (self.write_idx + 1) % len;
}
let mut acc = 0.0_f32;
for (k, &h) in self.taps.iter().enumerate() {
let idx = (self.write_idx + len - 1 - k) % len;
acc += h * self.history[idx];
}
acc
}
/// Clear filter history.
pub fn reset(&mut self) {
for v in &mut self.history {
*v = 0.0;
}
self.write_idx = 0;
}
}
#[cfg(test)]
mod tests {
use super::*;
fn fir(taps: usize, factor: usize) -> Vec<f32> {
design_lowpass_blackman(taps, 0.45 / factor as f32)
}
#[test]
fn upsampler_dc_gain_preserved() {
let h = fir(31, 4);
let mut up = PolyphaseUpsampler::new(4, &h);
// Drive DC and let the filter settle, then check unity gain.
let mut buf = [0.0_f32; 8];
let mut last_avg = 0.0;
for _ in 0..200 {
up.process(1.0, &mut buf);
last_avg = buf[..4].iter().sum::<f32>() / 4.0;
}
assert!((last_avg - 1.0).abs() < 1e-3, "got {last_avg}");
}
#[test]
fn down_then_up_roundtrip_is_bounded() {
// Stuff zero-padded input through up then down; output amplitude
// should approximately equal input on smooth signals.
let h = fir(31, 4);
let mut up = PolyphaseUpsampler::new(4, &h);
let mut down = PolyphaseDownsampler::new(4, &h);
let mut max_err = 0.0_f32;
let mut up_buf = [0.0_f32; 8];
// Drive a slow sine well below Nyquist.
for n in 0..2_000 {
let t = n as f32 / 48_000.0;
let x = (2.0 * std::f32::consts::PI * 1_000.0 * t).sin() * 0.5;
up.process(x, &mut up_buf);
let y = down.process(&up_buf[..4]);
// After group-delay warm-up, the error should be small.
if n > 80 {
max_err = max_err.max((x - y).abs());
}
}
// The filter is symmetric, so up/down with the same kernel
// introduces ~6 dB attenuation by design (each pass contributes
// half the gain). What we care about here is finite, bounded
// output and no runaway.
assert!(max_err < 1.0, "max_err {max_err}");
}
#[test]
fn upsampler_handles_impulse() {
let h = fir(15, 4);
let mut up = PolyphaseUpsampler::new(4, &h);
let mut buf = [0.0_f32; 8];
up.process(1.0, &mut buf);
// Some non-zero output expected on first impulse already.
assert!(buf[..4].iter().any(|&v| v.abs() > 1e-6));
// Drive zeros; output decays to zero.
for _ in 0..200 {
up.process(0.0, &mut buf);
}
assert!(buf[..4].iter().all(|&v| v.abs() < 1e-6));
}
}

108
crates/headroom-dsp/src/sliding_max.rs Normal file
View file

@ -0,0 +1,108 @@
//! Amortized-O(1) sliding-window maximum.
//!
//! Uses the standard monotonic-decreasing-deque trick. The deque's
//! capacity is bounded by `window`, so it never reallocates after
//! construction.
use std::collections::VecDeque;
/// Streaming max over a fixed-size sliding window.
pub struct SlidingMaxBuffer {
window: usize,
counter: u64,
/// `(index, value)`, monotonically decreasing in value from front.
deque: VecDeque<(u64, f32)>,
}
impl SlidingMaxBuffer {
/// Construct with the given window size. Lengths of 0 are clamped
/// to 1.
#[must_use]
pub fn new(window: usize) -> Self {
let window = window.max(1);
Self {
window,
counter: 0,
deque: VecDeque::with_capacity(window),
}
}
/// Window length in samples.
#[must_use]
pub fn window(&self) -> usize {
self.window
}
/// Push `value` and return the maximum over the most recent
/// `window` samples (inclusive of the value just pushed).
pub fn push_and_max(&mut self, value: f32) -> f32 {
// Drop entries that have aged out of the window.
let cutoff = self.counter.saturating_sub(self.window as u64 - 1);
while let Some(&(idx, _)) = self.deque.front() {
if idx < cutoff {
self.deque.pop_front();
} else {
break;
}
}
// Drop entries from the back smaller than the new value — they
// can never become the maximum.
while let Some(&(_, v)) = self.deque.back() {
if v <= value {
self.deque.pop_back();
} else {
break;
}
}
self.deque.push_back((self.counter, value));
self.counter += 1;
// SAFETY-ish: deque is non-empty (we just pushed).
self.deque.front().map_or(0.0, |&(_, v)| v)
}
/// Reset to empty state.
pub fn reset(&mut self) {
self.counter = 0;
self.deque.clear();
}
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn tracks_window_max() {
let mut s = SlidingMaxBuffer::new(3);
assert_eq!(s.push_and_max(1.0), 1.0);
assert_eq!(s.push_and_max(3.0), 3.0);
assert_eq!(s.push_and_max(2.0), 3.0);
assert_eq!(s.push_and_max(2.0), 3.0); // 3.0 aged out... actually still in (window=3, last 3 are [3,2,2])
assert_eq!(s.push_and_max(0.5), 2.0); // window is now [2,2,0.5]
assert_eq!(s.push_and_max(0.5), 2.0); // [2,0.5,0.5]
assert_eq!(s.push_and_max(0.5), 0.5); // [0.5,0.5,0.5]
}
#[test]
fn monotonically_decreasing_input() {
let mut s = SlidingMaxBuffer::new(4);
for (i, &v) in [5.0_f32, 4.0, 3.0, 2.0, 1.0, 0.5].iter().enumerate() {
let m = s.push_and_max(v);
// After window is filled, max is the value `window-1` back.
let expected = match i {
0..=3 => 5.0,
4 => 4.0,
_ => 3.0,
};
assert_eq!(m, expected);
}
}
#[test]
fn window_one_is_identity() {
let mut s = SlidingMaxBuffer::new(1);
for v in [1.0, 2.0, 0.5, 9.0_f32, -3.0] {
assert_eq!(s.push_and_max(v), v);
}
}
}

60
crates/headroom-dsp/src/util.rs Normal file
View file

@ -0,0 +1,60 @@
//! Common helpers: dB <-> linear conversions, time constants.
/// Lower bound used to avoid `log10(0)`.
pub const PEAK_FLOOR: f32 = 1e-20;
/// Convert linear amplitude to decibels. Inputs at or below
/// [`PEAK_FLOOR`] clamp to `-200 dB`.
#[must_use]
pub fn lin_to_db(x: f32) -> f32 {
if x <= PEAK_FLOOR {
-200.0
} else {
20.0 * x.log10()
}
}
/// Convert decibels to linear amplitude.
#[must_use]
pub fn db_to_lin(db: f32) -> f32 {
10.0_f32.powf(db / 20.0)
}
/// Convert a time constant in milliseconds to a one-pole smoother
/// coefficient at the given sample rate.
///
/// `y[n] = y[n-1] + alpha * (x[n] - y[n-1])`. The returned alpha is
/// `1 - exp(-1 / (tau * fs))` where `tau` is `time_ms / 1000`. A
/// `time_ms` of 0 or below returns `1.0` (instantaneous).
#[must_use]
pub fn time_to_alpha(time_ms: f32, sample_rate: f32) -> f32 {
if time_ms <= 0.0 || sample_rate <= 0.0 {
1.0
} else {
let tau_samples = (time_ms * 1e-3) * sample_rate;
1.0 - (-1.0 / tau_samples).exp()
}
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn db_round_trips() {
for db in [-60.0, -20.0, -6.0, -0.1, 0.0, 3.0, 6.0_f32] {
let lin = db_to_lin(db);
let back = lin_to_db(lin);
assert!((back - db).abs() < 1e-3, "db={db} back={back}");
}
}
#[test]
fn time_to_alpha_endpoints() {
assert_eq!(time_to_alpha(0.0, 48_000.0), 1.0);
assert!(time_to_alpha(1000.0, 48_000.0) < 0.01);
// Very fast attack: alpha approaches 1.
let a_fast = time_to_alpha(0.01, 48_000.0);
assert!(a_fast > 0.05, "fast alpha was {a_fast}");
}
}