# IIR biquad real-time filter just output noises

I've implemented a IIR filter but it outputs only noise (like FM radio out of tuning) and I can't see what's wrong.

This is the processing function:

void Filter::process(int size, double *in, double *out)
{
double xb[2] = { 0, 0 };
double yb[2] = { 0, 0 };

for (int n = 0; n < size; n++)
{
out[n] = (coeffs->a0 * in[n]) + (coeffs->a1 * xb[0])
+ (coeffs->a2 * xb[1]) - (coeffs->b1 * yb[0])
- (coeffs->b2 * yb[1]);

xb[1] = xb[0];
xb[0] = in[n];
yb[1] = yb[0];
yb[0] = out[n];
}
}


in and out are the buffers of sampled data at 48kHz; and size is 1024.

The coefficients are for ${Fc \over Fs} = {500 \over 48000 }$ and $Q = 0.7071$

a0:  0.00107054
a1:  0.00214108
a2:  0.00107054
b1:  -1.99572
b2:  0.953753

• you b1 and b2 coefs don't look like they result in stable poles. |b1| needs to be smaller than 1+b2 . |b1+b2| should be less than 1. – robert bristow-johnson Mar 13 '17 at 4:14
• i'm just curious: what does the letter b in your xb[2] and yb[2] states signify? – robert bristow-johnson Mar 13 '17 at 4:18
• @robertbristow-johnson In my crazy head the b stands for buffer – Victor Aurélio Mar 15 '17 at 5:02

I was unable to get a stable filter using your coefficients. You did not state, but am assuming this is a lowpass filter. Using the calculator here, I got:

coeffs = dict(
a0=0.0010232172047183973,
a1=0.0020464344094367946,
a2=0.0010232172047183973,
b1=-1.9075008174364765,
b2=0.91159368625535,
)


Filter Code:

I recast to a Transposed Direct Form 2 (Not necessary, but often recommended for floating point) as:

def filter_process(input_data):
zb = [0., 0.]

output_data = []
for data in input_data:

output_data.append(coeffs['a0'] * data + zb[0])
zb[0] = zb[1] + coeffs['a1'] * data - coeffs['b1'] * output_data[-1]
zb[1] = coeffs['a2'] * data - coeffs['b2'] * output_data[-1]

return output_data


Test Code:

import matplotlib.pyplot as plt
import numpy as np

def plot_f(f):
Fs = 44000
sample = 1000
x = np.arange(sample)
y = np.sin(2 * np.pi * f * x / Fs)
out_y = np.array(filter_process(y))
plt.plot(x, out_y)

for freq in (200, 500, 1000):
plot_f(freq)

plt.show()


Results:

If the input and output buffers are in fact the same buffer, that definitely will be problematic. The transfer function as currently organized, writes a new output to the buffer prior to finishing it use as an input. But that is easily remedied:

void Filter::process(int size, double *in, double *out)
{
double xb[2] = { 0, 0 };
double yb[2] = { 0, 0 };

for (int n = 0; n < size; n++)
{
double in_n = in[n]
out[n] = (coeffs->a0 * in_n) + (coeffs->a1 * xb[0])
+ (coeffs->a2 * xb[1]) - (coeffs->b1 * yb[0])
- (coeffs->b2 * yb[1]);

xb[1] = xb[0];
xb[0] = in_n;
yb[1] = yb[0];
yb[0] = out[n];
}
}

• Yep, low pass filter, I have used formulas on the EQ Cookbook[0] to generate the coefficients [0]: musicdsp.org/files/Audio-EQ-Cookbook.txt – Victor Aurélio Mar 12 '17 at 20:09
• Using coefficients generated from calculator you pointed out and your transfer function and got same result. – Victor Aurélio Mar 12 '17 at 20:36
• I do not have a C compiler handy, but I hand translated your C to python, and after I changed coefficients, your transfer function worked for me. That leaves as a key difference the input/output data. You did not specify what you were using. – Stephen Rauch Mar 12 '17 at 20:42
• The input comes from JackAudio to JackAudio. That's, its processing the system audio (music). – Victor Aurélio Mar 12 '17 at 20:46
• Would suggest some test data until you get this running, or an all-pass filter to verify your data pipeline. – Stephen Rauch Mar 12 '17 at 20:47

It looks like you got your a and b coefficients swapped. The feedback terms: b0,b1,b2, should be tiny (around 0.00x), and a0,a1,a2 can be close to 1 or 2 in magnitude.

• I cant decide whats a or b some people use a as zeros others uses poles. – Victor Aurélio Mar 12 '17 at 20:07
• @VictorAurélio The feedback terms or poles usually need to be small to keep thing from eventually blowing up (becoming NaNs or Inf). – hotpaw2 Mar 12 '17 at 22:52
• What I mean is the variable name, I find some to use b as poles, other people use b as zeros... – Victor Aurélio Mar 13 '17 at 16:43

The states should be declared outside the method and be members of your class so they are actually state-full:

class Filter
{
double xb[2] = { 0, 0 };
double yb[2] = { 0, 0 };

public:
void process(int size, double *in, double *out)
{
for (int n = 0; n < size; n++)
{
out[n] = (coeffs->a0 * in[n]) + (coeffs->a1 * xb[0])
+ (coeffs->a2 * xb[1]) - (coeffs->b1 * yb[0])
- (coeffs->b2 * yb[1]);

xb[1] = xb[0];
xb[0] = in[n];
yb[1] = yb[0];
yb[0] = out[n];
}
}
};

• That didn't help me, same result. – Victor Aurélio Mar 14 '17 at 13:56