3
$\begingroup$

I've implemented second order biquad filters using the famous RBJ Cookbook. I also want to implement first order lowpass and highpass filters. I was hoping to use the same code, so I need to define these first order filters in Direct Form I.

I found examples of first order filters in Direct Form II, but nothing for Direct Form I. I tried fudging the Direct Form II coefficients into Direct Form I, but I can't get it right.

What would be the correct way to calculating the coefficients? Thank you!

EDIT: I tried converting the Direct Form II filters to Direct Form I, but this makes them hard to automate because a0 is hardcoded to 1:

const K = Math.tan(Math.PI * freq / sampleRate);
const norm = 1 / (1 / K + 1);
b0 = norm;
b1 = norm;
b2 = 0;
a0 = 1;
a1 = (1 - 1 / K) * norm;
a2 = 0;

Is there a way to convert this filter (listed as first order lowpass in Biquad Calculator v3) to a standard Direct Form I filter?

$\endgroup$
4
  • $\begingroup$ You can use first-order filters only for LPF (and Low Shelf), HPF (and High Shelf), and APF. They're not hard to design, there are only 3 coefficients and two degrees of freedom to worry about (excluding overall gain) or one for APF. Sorry I didn't bother to include these designs. I also didn't include the Orfanidis variation of the Peaking EQ either. $\endgroup$ Mar 31 at 4:16
  • 1
    $\begingroup$ First get the analog prototypes with normalized $s$ for these first-order buggers. Then use the same Bilinear Transform: $$ \text{(normalized) } \qquad s \leftarrow \frac{1}{\tan \frac{\omega_0}{2}} \frac{1-z^{-1}}{1+z^{-1}}$$ $\endgroup$ Mar 31 at 4:20
  • $\begingroup$ Thank you for this answer! I've tried using the bilinear transform, but it's a bit too hard for me to understand as I don't come from an engineering background. I found out how to translate the Direct Form II filters, but since a0 is hardcoded to 1 they can't be automated cleanly like the original filters $\endgroup$ Mar 31 at 6:07
  • 1
    $\begingroup$ LPF1: H(s) = 1 / (1 + s/w0) --> num = 1, den = [1 w0]. HPF1: H(s) = s/w0 / (1 + s/w0) --> num = [w0 0] and den = [1 w0]. w0 = 2*pi*f . $\endgroup$
    – Juha P
    Mar 31 at 17:02

2 Answers 2

3
$\begingroup$

Okay, I'm gonna try to toss out here my first guess of general 1st-order filters done in the style of the Audio EQ Cookbook. The cookbook assumes the same EE definition of resonant frequency and normalizes that parameter out of $s$ in the analog transfer functions and then uses this version of the Bilinear Transform with prewarping:

$$\begin{align} \text{(normalized)} \qquad s &\leftarrow \frac{1}{\tan \frac{\omega_0}{2}} \times \frac{1-z^{-1}}{1+z^{-1}} \\ \\ &= \frac{1+\cos(\omega_0)}{\sin(\omega_0)} \times \frac{1-z^{-1}}{1+z^{-1}} \\ \end{align}$$

where

$$ \omega_0 \triangleq 2 \pi \frac{f_0}{F_\mathrm{s}} $$

and $f_0$ is the "significant frequency", which is the -3 dB frequency for LPF and HPF and the shelf-midpoint frequency for the shelf filters (which is how the Cookbook defines it for its 2nd-order shelf filters). $F_\mathrm{s}$ is the sample rate.

For the shelf filters

$$ A \triangleq \sqrt{10^{dB_\mathrm{gain}/20}} = 10^{dB_\mathrm{gain}/40} $$

The $z$-plane transfer function is:

$$\begin{align} H(z) &= \frac{b_0 + b_1 z^{-1}}{a_0 + a_1 z^{-1}} \\ \\ &= \frac{\frac{b_0}{a_0} + \frac{b_1}{a_0} z^{-1}}{1 + \frac{a_1}{a_0} z^{-1}} \\ \end{align}$$

The Direct I Form implementation is:

$$ y[n] = \tfrac{b_0}{a_0} x[n] + \tfrac{b_1}{a_0} x[n-1] - \tfrac{a_1}{a_0} y[n-1]$$


LPF

$$ H(s) = \frac{1}{1+s} $$

$$\begin{align} a_0 &= \sin(\omega_0) + 1 + \cos(\omega_0) \\ a_1 &= \sin(\omega_0) - 1 - \cos(\omega_0)\\ b_0 &= \sin(\omega_0) \\ b_1 &= \sin(\omega_0) \\ \end{align}$$

HPF

$$ H(s) = \frac{s}{1+s} $$

$$\begin{align} a_0 &= \sin(\omega_0) + 1 + \cos(\omega_0) \\ a_1 &= \sin(\omega_0) - 1 - \cos(\omega_0)\\ b_0 &= 1 + \cos(\omega_0) \\ b_1 &= -1 - \cos(\omega_0) \\ \end{align}$$

APF

$$ H(s) = \frac{1-s}{1+s} $$

$$\begin{align} a_0 &= \sin(\omega_0) + 1 + \cos(\omega_0) \\ a_1 &= \sin(\omega_0) - 1 - \cos(\omega_0) \\ b_0 &= \sin(\omega_0) - 1 - \cos(\omega_0) \\ b_1 &= \sin(\omega_0) + 1 + \cos(\omega_0) \\ \end{align}$$

LowShelf

$$ H(s) = \frac{A+s}{\frac{1}{A}+s} $$

$$\begin{align} a_0 &= \tfrac{1}{A}\sin(\omega_0) + 1 + \cos(\omega_0) \\ a_1 &= \tfrac{1}{A}\sin(\omega_0) - 1 - \cos(\omega_0) \\ b_0 &= A\sin(\omega_0) + 1 + \cos(\omega_0) \\ b_1 &= A\sin(\omega_0) - 1 - \cos(\omega_0) \\ \end{align}$$

HighShelf

$$ H(s) = \frac{1+As}{1+\frac{s}{A}} $$

$$\begin{align} a_0 &= \sin(\omega_0) + A + A\cos(\omega_0) \\ a_1 &= \sin(\omega_0) - A - A\cos(\omega_0) \\ b_0 &= \sin(\omega_0) + \tfrac{1}{A} + \tfrac{1}{A}\cos(\omega_0) \\ b_1 &= \sin(\omega_0) - \tfrac{1}{A} - \tfrac{1}{A}\cos(\omega_0) \\ \end{align}$$

This maybe should get added to the Cookbook. I guess I never really thought of doing that. Ain't that weird?

$\endgroup$
3
  • $\begingroup$ Comments have been moved to chat; please do not continue the discussion here. Before posting a comment below this one, please review the purposes of comments. Comments that do not request clarification or suggest improvements usually belong as an answer, on Signal Processing Meta, or in Signal Processing Chat. Comments continuing discussion may be removed. $\endgroup$
    – Peter K.
    Apr 3 at 18:38
  • $\begingroup$ I think you might be a little bit trigger happy about moving comments to chat, Peter. You probably knew that. I wanted to finish the discussion with @DanBoschen about what happens with phase for the APF. $\endgroup$ Apr 3 at 19:55
  • $\begingroup$ SE.SP is not a discussion site. I've taken a zero tolerance approach to comments: when the system flags more than 20 comments on post, the comments get moved to chat. Zero tolerance. $\endgroup$
    – Peter K.
    Apr 6 at 22:58
3
$\begingroup$

A solution with the DC or Nyquist bin zeroed for a high pass or low pass respectively with real coefficients is:

LOW PASS

$b_0 = \frac{1+a_1}{2}$
$b_1 = b_0$
$b_2 = 0$
$a_0 = 1$
$a_1 = -\frac{\cos(\omega_c)}{\sin(\omega_c)+1}$
$a_2 = 0$

HIGH PASS

$b_0 = \frac{1-a_1}{2}$
$b_1 = -b_0$
$b_2 = 0$
$a_0 = 1$
$a_1 = -\frac{\cos(\omega_c)}{\sin(\omega_c)+1}$
$a_2 = 0$

Where $\omega_c$ is the cutoff frequency as the normalized radian frequency in units of radians/sample. For a sampling frequency in Hz as $f_s$ and cutoff frequency in Hz as $f$, the cutoff frequency in normalied radian units is given as:

$$\omega_c = \frac{2\pi f}{f_s}$$

For the above block diagrams, $a_0$ not shown explicitly is normalized such that $a_0=1$.

Validating this, I tested with $\omega_c=\pi/4$ and $\omega_c=3\pi/4$ using the functions in Octave as (and plotted the magnitude responses by passing num and den into freqz):

function [num, den] = lowpass(omega)
  a0 = 1;
  a1 = -cos(omega)/(sin(omega)+1);
  b0= (1+a1)/2;
  b1 = b0;
  num = [b0 b1];
  den = [a0 a1];
endfunction

function [num, den] = highpass(omega)
  a0 = 1;
  a1 = -cos(omega)/(sin(omega)+1);
  b0= (1-a1)/2;
  b1 = -b0;
  num = [b0 b1];
  den = [a0 a1];
endfunction

magnitude response

The coefficients don't change between Direct Form I and Direct Form II: Direct Form II is simply Direct Form I with the two sections reordered:

Direct Form I

Direct Form I reordered

Direct Form II

DETAILS / DERIVATION

The OP had originally posted the question using $a1=0$ and $a2=0$. This would have been a simple first order FIR filters given as:

$$H(z) = K(1 - b_1z^{-1})$$

Where $K$ is a scaling constant for unity gain. If we assume real coefficients only, then for $b_1>0$ we get a low pass response, and for $b_1<0$ we get a high pass response with a single real zero at $z=b_1$ and a trivial pole at the origin $z=0$. This structure is of limited use since the dynamic range between the frequency at $\omega=0$ (DC) and $\omega = \pi$ ($f_s/2$ where $f_s$ is the sampling rate) will diminish rapidly as the zero moves away from $z=1$ or $z=-1$.

A more useful but still first order filter, as an IIR filter would be given as:

$$H(z) = K\frac{1+b_1z^{-1}}{1+a_1z^{-1}}$$

For this case, the filter has a single zero at $z=-b_1$ and pole at $z=-a_1$. For the low-pass case, set the zero to be $z=-1$ placing a frequency response null at Nyquist resulting in $b_1=1$. For the high-pass case, set the zero to be $z=1$ placing a frequency response null at DC and resulting in $b_1=-1$. We can set the bandwidth as a function of the pole placement in each case as follows:

LOW PASS FILTER

$$H(z) = K\frac{1+z^{-1}}{1+a_1z^{-1}}= K\frac{z+1}{z+a_1}$$

Set $K$ to normalize gain = 1 at DC ($z=1$):

$$H(z) = \bigg(\frac{1+a_1}{2}\bigg)\frac{z+1}{z+a_1}\label{1}\tag{1}$$

The squared magnitude response is given as:

$$|H(e^{j\omega})|^2 = H(e^{j\omega})H(e^{-j\omega}) $$

$$= \bigg(\frac{1+a_1}{2}\bigg)^2\bigg(\frac{e^{j\omega}+1}{e^{j\omega}+a_1}\bigg)\bigg(\frac{e^{-j\omega}+1}{e^{-j\omega}+a_1}\bigg) =\bigg(\frac{1+a_1}{2}\bigg)^2\frac{2+2\cos(\omega)}{1+2a_1\cos(\omega)+a_1^2}\label{2}\tag{2}$$

The -3 dB cutoff frequency $\omega_c$ is where the magnitude squared goes to 1/2 power:

$$|H(e^{j\omega_c})|^2 = \bigg(\frac{1+a_1}{2}\bigg)^2\frac{2+2\cos(\omega_c)}{1+2a_1\cos(\omega_c)+a_1^2} = 1/2 \label{3}\tag{3}$$

Solving \ref{3} for $\omega_c$ as a function of real $a_1$, we get:

$$\omega_c = \cos^{-1}\bigg(\frac{ -2a_1}{a_1^2+1} \bigg)\label{4}\tag{4}$$

Rearranging \ref{4} to solve for $a_1$ as a function of $\omega_c$:

$$a_1 = -\frac{\cos(\omega)}{\sqrt{\sin^2(\omega)}+1}\label{5}\tag{5}$$

Since $0<\omega<\pi$, ref{5} reduces further to be:

$$a_1 = -\frac{\cos(\omega)}{\sin(\omega)+1}\label{6}\tag{6}$$

HIGH PASS FILTER

Similarly solving for the case of the High Pass Filter, where:

$$H(z) = \bigg(\frac{1-a_1}{2}\bigg)\frac{z-1}{z+a_1}\label{7}\tag{7}$$

$$|H(e^{j\omega_c})|^2 = \bigg(\frac{1-a_1}{2}\bigg)^2\frac{2-2\cos(\omega_c)}{1+2a_1\cos(\omega_c)+a_1^2} = 1/2 \label{8}\tag{8}$$

Which results in the same relationship for $a_1$ and $\omega_c$ as the Low Pass case.


$\endgroup$
6
  • $\begingroup$ Thank you! I noticed something seems different with these coefficients, I think a1 should be negative. If I set a1 = -a1 or set a0 = -1 it gives the same coefficients as my original code, which seems to produce the correct frequency response. $\endgroup$ Mar 31 at 9:40
  • $\begingroup$ @MysteryPancake, Ah good catch, thank you! I am finding other mistakes from my late night attempts at math...fixing now and will verify. $\endgroup$ Mar 31 at 13:17
  • 1
    $\begingroup$ I got other fish to fry at the moment (Easter), but you should do the first-order recursive APF and then low-shelf and high-shelf. $\endgroup$ Mar 31 at 14:37
  • 2
    $\begingroup$ Happy Easter @robertbristow-johnson !! I will be frying my own fish very soon. $\endgroup$ Mar 31 at 14:48
  • 1
    $\begingroup$ Fish? It's chocolate in my house. Rabbits for preference. :-) $\endgroup$
    – Peter K.
    Mar 31 at 17:20

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.