A tower crane carries a tilt sensor on its jib, read fs times per second. The jib really does lean slowly as loads are lifted, but the reading also shakes with the cable and the wind. The safety computer smooths it with a biquad: a recursive (IIR) filter that uses its own last two outputs as well as the last inputs, so a sharp low-pass needs only five weights.
From the cutoff fc (hertz) and the quality factor q (how sharply the response turns the corner; q = 0.7071... gives the flattest passband), the coefficients are:
w0 = 2 * pi * fc / fs
alpha = sin(w0) / (2 * q)
b0 = (1 - cos(w0)) / 2 b1 = 1 - cos(w0) b2 = (1 - cos(w0)) / 2
a0 = 1 + alpha a1 = -2 * cos(w0) a2 = 1 - alpha
Divide all five of b0, b1, b2, a1, a2 by a0. Then each new reading x[n] gives the output
y[n] = b0 * x[n] + b1 * x[n-1] + b2 * x[n-2] - a1 * y[n-1] - a2 * y[n-2]
where every input and output before the first sample is 0. Compute the outputs in order: y[n] needs y[n-1] and y[n-2].
Write biquad_lowpass(x, fc, fs, q) that returns the list of outputs, one per reading. Press Run with plot(x, label="raw") and plot(biquad_lowpass(x, 2, 100, 0.7071), label="smoothed") to compare.
Examples
Input: x = [1, 0, 0, 0], fc = 25, fs = 100, q = 0.7071067811865476
Output: [0.2928932188134524, 0.5857864376269049, 0.24264068711928521, -0.10050506338833858]
Explanation: w0 = pi / 2, so cos(w0) = 0 and alpha = 0.7071. After dividing by a0 = 1.7071,
b0 = 0.2929, b1 = 0.5858, b2 = 0.2929, a1 = 0, a2 = 0.1716. The outputs are the filter's
impulse response: y[2] = 0.2929 - 0.1716 * 0.2929, y[3] = -0.1716 * 0.5858.
Input: x = [2, 2, 2], fc = 25, fs = 100, q = 0.7071067811865476
Output: [0.5857864376269049, 1.7573593128807146, 2.2426406871192857]
Explanation: a steady reading of 2 is passed (after a short overshoot) at its full size,
because b0 + b1 + b2 = 1 + a1 + a2: the gain at 0 Hz is exactly 1.
Constraints
0 <= len(x) <= 10**5,0 < fc < fs / 2,0.3 <= q <= 10; answers are compared with a tolerance of1e-6
Goals
- Compute the five coefficients of a biquad low-pass from its cutoff and Q
- Run a recursive filter sample by sample, feeding back its own past outputs
- See that a filter with feedback needs only five multiplications per sample