svf.hpp

Go to the documentation of this file.

// sbl/widgets/proc/svf.hpp — State Variable Filter (Audio Stack — Widgets)
//
// Zero-delay feedback (ZDF) SVF based on Andrew Simper's topology, as
// documented in Vadim Zavalishin's "The Art of VA Filter Design."
// Reference implementation: Mutable Instruments stmlib/dsp/filter.h (MIT).
//
// Unconditionally stable at all cutoff and resonance settings — no clamping
// needed, no oversampling required. This replaces the previous Chamberlin
// SVF which was unstable at high cutoff + resonance near Nyquist.
//
// Uses float state variables — STM32H7 (Cortex-M7) has hardware FPU with
// single-cycle float ops. Fixed-point variant for M0+ is future work.
//
// Resonance compensation: optional input gain reduction as resonance
// increases, preventing output clipping during filter sweeps. This is
// the "volume drop" characteristic of classic analog filters (SH-101,
// MS-20). Disabled by default (k=0) for backward compatibility.
//
// process_modulated() provides per-sample exponential cutoff modulation
// from a normalized signal [-1, 1]. This eliminates the manual per-sample
// coefficient loop that every LFO→filter instrument must write. Uses
// sbl::signal::exp_mod() internally for perceptually correct scaling.
// See FDP-031 for architecture rationale, RPT-014 for the triggering bug.
//
// Usage:
//   sbl::widgets::proc::Svf filter;
//   filter.set_sample_rate(48000);
//   filter.set_cutoff(1000.0f);   // 1 kHz
//   filter.set_resonance(0.75f);  // moderate resonance (0.0 = flat, 1.0 = max)
//   filter.set_resonance_compensation(3.0f);  // analog-style volume drop
//   filter.set_mode(FilterMode::LowPass);
//   filter.process(buf, frames);
//
//   // Per-sample modulated filtering (LFO → cutoff):
//   lfo.process(lfo_buf, frames);
//   filter.process_modulated(buf, frames, lfo_buf, 1000.0f, 2.0f);  // ±2 oct
 
#pragma once
 
#include <cstdint>
 
#include <sbl/dsp/fast_math.hpp>
#include <sbl/signal/exp_mod.hpp>
 
namespace sbl::widgets::proc {
 
enum class FilterMode : uint8_t {
    LowPass,
    HighPass,
    BandPass,
    Notch,
};
 
class Svf {
public:
    /// @note All public methods are ISR-safe — bounded computation, no I/O.
 
    void set_sample_rate(uint32_t sr) { sample_rate_ = sr; }
 
    /**
     * @brief Set cutoff frequency in Hz
     * @param freq_hz Cutoff in Hz (e.g., 1000.0f)
     */
    void set_cutoff(float freq_hz) {
        float f = freq_hz / static_cast<float>(sample_rate_);
        if (f > 0.497f) f = 0.497f;
        g_ = dsp::fast_tan_pif(f);
        h_ = 1.0f / (1.0f + r_ * g_ + g_ * g_);
    }
 
    /**
     * @brief Set resonance
     * @param q 0.0 = no resonance (flat), 1.0 = maximum resonance
     */
    void set_resonance(float q) {
        // Map float to damping: 0.0 → 2.0 (flat), 1.0 → 0.0 (self-oscillation)
        q_ = q;
        r_ = (1.0f - q) * 2.0f;
        h_ = 1.0f / (1.0f + r_ * g_ + g_ * g_);
        update_comp_gain();
    }
 
    /**
     * @brief Set resonance compensation strength
     *
     * Controls how much the input signal is attenuated as resonance
     * increases: gain = 1 / (1 + q * k). This prevents the resonant
     * peak from clipping during filter sweeps.
     *
     *   k = 0   No compensation (default, backward compatible)
     *   k = 2   Gentle drop (~33% volume at full resonance)
     *   k = 4   Strong drop (~20% volume at full resonance)
     *
     * @param k Compensation amount (0.0 = off)
     */
    void set_resonance_compensation(float k) {
        comp_k_ = k;
        update_comp_gain();
    }
 
    /** @brief Set filter output mode */
    void set_mode(FilterMode mode) { mode_ = mode; }
 
    /**
     * @brief Process a block of audio samples in-place
     *
     * ZDF SVF — unconditionally stable, no oversampling needed.
     *
     * @param buf Float audio buffer (modified in-place, [-1.0, 1.0])
     * @param frames Number of samples
     */
    void process(float* buf, uint16_t frames) {
        float g = g_, r = r_, h = h_;
        float comp = comp_gain_;
        float s1 = state_1_, s2 = state_2_;
 
        for (uint16_t i = 0; i < frames; ++i) {
            float x = buf[i] * comp;
            float hp = (x - r * s1 - g * s1 - s2) * h;
            float bp = g * hp + s1;
            s1 = g * hp + bp;
            float lp = g * bp + s2;
            s2 = g * bp + lp;
 
            switch (mode_) {
                case FilterMode::LowPass:  buf[i] = lp; break;
                case FilterMode::HighPass:  buf[i] = hp; break;
                case FilterMode::BandPass:  buf[i] = bp; break;
                case FilterMode::Notch:     buf[i] = lp + hp; break;
                default:                    buf[i] = lp; break;
            }
        }
        state_1_ = s1;
        state_2_ = s2;
    }
 
    /**
     * @brief Process with per-sample exponential cutoff modulation
     *
     * Applies exponential modulation to the cutoff frequency per sample,
     * recomputing ZDF coefficients each sample. The modulation signal is
     * normalized [-1, 1] and depth is in octaves — perceptually uniform
     * regardless of the center frequency.
     *
     * This replaces the manual per-sample loop pattern:
     *   for (i = 0; i < n; ++i) {
     *       filter.set_cutoff(cutoff + lfo[i]);  // WRONG: linear, broken
     *       filter.process(&buf[i], 1);
     *   }
     *
     * With:
     *   filter.process_modulated(buf, n, lfo, cutoff, 2.0f);  // ±2 octaves
     *
     * @param buf     Float audio buffer (modified in-place)
     * @param frames  Number of samples
     * @param mod     Modulation signal buffer (normalized [-1.0, 1.0])
     * @param center  Center cutoff frequency in Hz
     * @param depth   Modulation depth in octaves (e.g., 2.0 = ±2 octaves)
     * @param lo      Minimum cutoff in Hz (default 20)
     * @param hi      Maximum cutoff in Hz (default 20000)
     */
    void process_modulated(float* buf, uint16_t frames,
                           const float* mod, float center, float depth,
                           float lo = 20.0f, float hi = 20000.0f) {
        float r = r_;
        float comp = comp_gain_;
        float s1 = state_1_, s2 = state_2_;
        float sr = static_cast<float>(sample_rate_);
 
        for (uint16_t i = 0; i < frames; ++i) {
            // Exponential modulation: center * 2^(mod * depth)
            float freq = signal::exp_mod(center, mod[i], depth, lo, hi);
 
            // Recompute ZDF coefficients for this sample
            float f = freq / sr;
            if (f > 0.497f) f = 0.497f;
            float g = dsp::fast_tan_pif(f);
            float h = 1.0f / (1.0f + r * g + g * g);
 
            // ZDF SVF tick
            float x = buf[i] * comp;
            float hp = (x - r * s1 - g * s1 - s2) * h;
            float bp = g * hp + s1;
            s1 = g * hp + bp;
            float lp = g * bp + s2;
            s2 = g * bp + lp;
 
            switch (mode_) {
                case FilterMode::LowPass:  buf[i] = lp; break;
                case FilterMode::HighPass:  buf[i] = hp; break;
                case FilterMode::BandPass:  buf[i] = bp; break;
                case FilterMode::Notch:     buf[i] = lp + hp; break;
                default:                    buf[i] = lp; break;
            }
        }
        state_1_ = s1;
        state_2_ = s2;
 
        // Update cached coefficients to reflect last sample's frequency
        // so that a subsequent process() call uses a reasonable state
        float f_last = signal::exp_mod(center, mod[frames - 1], depth, lo, hi) / sr;
        if (f_last > 0.497f) f_last = 0.497f;
        g_ = dsp::fast_tan_pif(f_last);
        h_ = 1.0f / (1.0f + r_ * g_ + g_ * g_);
    }
 
    /** @brief Reset filter state */
    void reset() {
        state_1_ = 0.0f;
        state_2_ = 0.0f;
    }
 
private:
    void update_comp_gain() {
        comp_gain_ = (comp_k_ > 0.0f)
            ? 1.0f / (1.0f + q_ * comp_k_)
            : 1.0f;
    }
 
    float g_ = 0.0f;         // tan(π·f) — frequency coefficient
    float r_ = 1.0f;         // Damping (2.0 = flat, 0.0 = self-oscillation)
    float h_ = 0.0f;         // 1 / (1 + r·g + g²)
    float q_ = 0.0f;         // Resonance (0.0–1.0), cached for compensation
    float comp_k_ = 0.0f;    // Compensation strength (0 = off)
    float comp_gain_ = 1.0f; // Precomputed: 1 / (1 + q * k)
    float state_1_ = 0.0f;
    float state_2_ = 0.0f;
    FilterMode mode_ = FilterMode::LowPass;
    uint32_t sample_rate_ = 48000;
};
 
} // namespace sbl::widgets::proc