Realtime Biquad Dsp Agent Skill

A real-time biquad in pure Kotlin

A biquad is five coefficients and four state values per channel. No library is needed and no

allocation happens per sample:

val y = b0 * x + b1 * x1 + b2 * x2 - a1 * y1 - a2 * y2
x2 = x1; x1 = x
y2 = y1; y1 = y

The coefficients come from the standard audio-EQ-cookbook formulas (Robert Bristow-Johnson).

For a low-pass, with ω = 2π·fc/fs, α = sin(ω)/(2Q):

val v = (1.0 - cos(omega)) / 2.0
b0 = v;  b1 = 1.0 - cos(omega);  b2 = v
val a0 = 1.0 + alpha
a1 = -2.0 * cos(omega);  a2 = 1.0 - alpha
// then divide all five by a0

High-pass is the same shape with v = (1 + cos ω)/2 and b1 = −(1 + cos ω); a0, a1, a2 are

identical. Reach for this when you need a *sweepable* filter — a transition effect, tone control,

crossover — inside an audio callback that must not allocate.

Traps

Every channel needs its own four state values. Sharing x1/x2/y1/y2 between left and right feeds

each channel the other's history: not "slightly wrong" — a mono-ish smear where the stereo image

collapses and moving content sounds phasey. Name the fields so the mistake is visible — x1L, x1R

— and treat "one filter object, N channel state sets" as the shape, not "N filter objects".

Cascading N sections needs N *state sets* as well as N coefficient sets. Two stages sharing delay

lines is not a 4th-order filter — it is some other filter, each stage feeding on the other's

history. Coefficients may legitimately be identical across stages; the state never is:

val mid = b0_1 * input + b1_1 * s1_x1L + b2_1 * s1_x2L - a1_1 * s1_y1L - a2_1 * s1_y2L
s1_x2L = s1_x1L; s1_x1L = input;  s1_y2L = s1_y1L; s1_y1L = mid
// stage 2 is the identical shape, reading `mid` and writing its OWN s2_* state — never s1_*

"Two cascaded Butterworth sections" is not a 4th-order Butterworth. A single Q = 0.707 section

is −3 dB at the cutoff; two of them in series is −6 dB there (reproduced below). A textbook

4th-order Butterworth uses two sections with *different* Q values and stays −3 dB at cutoff.

Cascading identical 0.707 sections is a perfectly good choice — it is the Linkwitz–Riley shape used

for crossovers — but the cutoff you set is no longer the −3 dB point, misleading anyone who reads

"4th-order Butterworth" into hunting a bug that isn't there.

Keep every stage in the chain at neutral gain — do not skip it. In a bank of peaking sections it

is tempting to build only the bands the user moved. Then every adjustment resizes the coefficient

array *and the state array with it*, and a resized state array is a cleared one: the delay lines

vanish mid-stream and the user hears a click on every drag. With A = 10^(gainDb/40):

b = [1 + alpha*A, -2*cos(w0), 1 - alpha*A]      a = [1 + alpha/A, -2*cos(w0), 1 - alpha/A]

At 0 dB, A is exactly 1, so alpha*A and alpha/A are the same number and b equals a term

for term — b0 normalises to exactly 1.0 in IEEE 754, and b1 == a1, b2 == a2 (reproduced

below). So H(z) ≡ 1, and the array sizes depend only on the format. It is not *bit*-identical:

accumulating left to right does not cancel exactly: 200 000 random samples from cleared state give

one neutral stage (1 kHz, 48 kHz, Q 1.41) max |y−x| = 1.0e-14 (190 dB below the 16-bit LSB), a

ten-stage flat bank 1.6e-12. Claim "identity transfer function" in comments, never "bit-identical"

— checkable, and false. The same fix applies when a band's centre exceeds Nyquist (16 kHz on a 22.05

kHz-Nyquist stream): substitute the identity stage rather than a continue, for the identical

fixed-array-size reason.

Clear the history on bypass, on a stream change, and on disable — any discontinuity. Delay lines

shaped by the *previous* setting splice their old shape onto the new one if left in place, even when

every stage is neutral, and even on an ordinary seek. Example: seed a 500 Hz stage (48 kHz, Q 1.41)

at +9 dB with a full-scale sine, switch to 0 dB and feed silence — residue starts at 1.58 full

scale, still 1.5e-2 two hundred samples later (figures are setup-specific). reset() zeroes state

only — never coefficients, or the first buffer after every seek runs at pass-through.

Capture a0 before you overwrite anything with it. The normalisation divides all five coefficients

by the raw a0, and a0 is built from alpha, which several formulas also reuse. Capture it in a

local, then divide: doing it in place, one coefficient at a time, in the wrong order, produces

plausible-looking coefficients and a filter that either goes silent or runs away.

Clamp the cutoff into 20 .. (sampleRate / 2 − 1). Past half the sample rate the design is

meaningless: the recursion goes unstable and overflows into NaN — and once NaN enters y1, every

sample after it is NaN forever, since the state feeds back. Silence that never recovers after one

bad parameter write is the signature.

val clamped = cutoffHz.coerceIn(20f, (sampleRate / 2f) - 1f)

Recompute lazily, on a dirty flag, and only at a buffer boundary. Trigonometry per sample is

wasteful; worse, coefficients written from another thread mid-buffer leave the filter running half

one set, half another. Set @Volatile coefficientsDirty from every parameter setter — cutoff

and filter type, forgetting the type is the classic miss — and consult it once per buffer.

Sweep exponentially, not linearly. Hearing is logarithmic: a linear ramp from 20 kHz to 200 Hz spends

most of its time in a region nobody can hear it move through and then lurches. Interpolate in the log domain:

fun exponentialInterpolate(start: Float, end: Float, t: Float): Float {
    if (start <= 0f || end <= 0f) return end
    return exp(ln(start) + (ln(end) - ln(start)) * t).toFloat()
}

Coefficients are cross-thread, state is not. Mark the coefficient fields @Volatile because a UI

or animation thread writes them; leave the state fields plain, because only the audio thread touches

them and volatile reads in the inner loop are pure cost.

Verifying it

• **Two numeric claims above, reproduced directly: the identity-stage coefficients and the −3 dB vs

−6 dB cascade:**

   python3 -c "
   import math, cmath
   w=2*math.pi*1000/48000; c=math.cos(w); al=math.sin(w)/(2*1.41); A=1.0
   b0,b1,b2 = 1+al*A, -2*c, 1-al*A
   a0,a1,a2 = 1+al/A, -2*c, 1-al/A
   print('identity:', b0/a0, b1/a0, b2/a0, '| b1==a1:', b1==a1, 'b2==a2:', b2==a2)
   al2=math.sin(w)/(2*0.707); b0=(1-c)/2/(1+al2); b1=(1-c)/(1+al2); a1=-2*c/(1+al2); a2=(1-al2)/(1+al2)
   z=cmath.exp(-1j*w); H=(b0+b1*z+b0*z*z)/(1+a1*z+a2*z*z)
   print('cascade dB:', 20*math.log10(abs(H)), 20*math.log10(abs(H*H)))"

Pass condition: identity prints 1.0 -1.8951700997919003 0.9115234478756887 (b0 exactly 1.0)

plus b1==a1: True b2==a2: True; cascade dB prints roughly -3.01 -6.02 — one stage is

textbook, two aren't.

• By hand: sweep a filter's cutoff slowly while audio plays, then separately toggle it off and back

on — a click or zipper noise while sweeping means coefficients aren't gated on a dirty flag or are

recomputed mid-buffer, a thump after re-enabling means the delay lines weren't cleared. Separately, feed

a steady sine at the cutoff through one 0.707-Q stage then the same stage cascaded twice — the second

should sound roughly twice as attenuated, confirming −3/−6 dB by ear.