mathr / blog / #

One representation of a biquad filter's transfer function is:

\[ H(z) = \frac{ b_0 + b_1 z^{-1} + b_2 z^{-2} } { a_0 + a_1 z^{-1} + a_2 z^{-2} } \]

Without loss of generality, it can be assumed that \(a_0 = 1\) (otherwise divide both top and bottom by the original \(a_0\) to give a new set of coefficients).

Another representation of the transfer function is:

\[ H(z) = g \frac{ (z - Z) (z - Z^*) } { (z - P) (z - P^*) } \]

where \({}^*\) denotes complex conjugation. Here \(Z\) is a zero and \(P\) is a pole, and both are complex, and \(g\) is gain factor, which is real. It will make things easier to assume \(g = 1\), if it isn't then just multiply the input \(x\) by \(g\) before-hand. Multiplying out the top and bottom, and dividing both by \(z^2\), recovers the first representation:

\[ H(z) = \frac{ 1 - 2 \operatorname{Re}{(Z)} z^{-1} + |Z|^2 z^{-2} } { 1 - 2 \operatorname{Re}{(P)} z^{-1} + |P|^2 z^{-2} } \]

Or more explicitly:

\[ \begin{aligned} a_0 &= 1 & a_1 &= -2 \operatorname{Re}{(P)} & a_2 &= |P|^2 \\ b_0 &= 1 & b_1 &= -2 \operatorname{Re}{(Z)} & b_2 &= |Z|^2 \end{aligned} \]

The utility of the second pole-zero representation is in filter design, the former is useful in filter implementation. The direct form 2 implementation is a pair of recurrence relations, where \(x\) is the input, \(y\) is the output, and \(w\) is only used internally:

\[ \begin{aligned} w_n &= \phantom{a_0} x_n - a_1 w_{n-1} - a_2 w_{n-2} \\ y_n &= b_0 w_n + b_1 w_{n-1} + b_2 w_{n-2} \end{aligned} \]

The software Pure-data implements its biquad~ object using direct form 2, except with \(a_1\) and \(a_2\) negated. So an instantiation biquad~ A B C D E has the correspondence:

\[ \begin{aligned} A &= -a_1 & B &= -a_2 & C &= b_0 & D &= b_1 & E &= b_2 \end{aligned} \]

The final puzzle is the factorization from direct form 2 into pole-zero representation. Assuming \(a_0 = 1 = b_0\) (otherwise factor out the gain coefficient and normalize), we get:

\[ \begin{aligned} \operatorname{Re}{(Z)} &= -\frac{1}{2} b_1 & \operatorname{Im}{(Z)} &= \sqrt{ b_2 - \frac{1}{4} b_1^2 } \\ \operatorname{Re}{(P)} &= -\frac{1}{2} a_1 & \operatorname{Im}{(P)} &= \sqrt{ a_2 - \frac{1}{4} a_1^2 } \end{aligned} \]

The motivation for this post was mkfilter giving code output which is inherently unstable for higher order filters, reducing it to smaller stages (like biquads) in serial is much less likely to blow up. I did some experiments converting the pole-zero output from mkfilter into Pd patches using cpole~ and czero~, but biquads are more efficient as the input and ouput are both real (a biquad needs 5 real multiplications and 4 real additions, whereas a single complex multiplication needs 4 real multiplications and 2 real additions). The next step is to write a program into which I can pipe mkfilter output, and get out a Pd patch using some chained biquads.

In the process I fixed some minor issues in the mkfilter source code (retrieved 2015-03-26) and wrote a new build system to get it to compile in the modern age.