First-order digital lowpass filters

This article shows the frequency filtering properties of the simplest recursive digital filter. The filter has a lowpass characteristic, meaning low frequencies are unaffected and high frequencies are filtered out. We will only consider a first-order filter, i.e. a filter with only one delay element.

Frequency response

The magnitude response versus frequency of the first-order digital lowpass filter is shown in Figure 1. The figure shows multiple lines; one for each chosen filter coefficient ‘c’. Useful values for $c$ range from $0.5$ to almost zero, where closer to zero means a lower cutoff frequency (more frequencies get filtered out).

The figure shows lines for $c=\frac{1}{2}$, $\frac{1}{4}$, $\frac{1}{8}$, $\frac{1}{16}$, $\frac{1}{32}$, $\frac{1}{64}$, $\frac{1}{128}$, $\frac{1}{256}$, $\frac{1}{512}$ and $\frac{1}{1024}$. These values can be implemented efficiently and without the use of floating-point arithmetic if desired.

Every real-world frequency depends on the sample rate used. In order to present information independent of the sample rate, frequencies are given as normalised frequencies. Here, a normalised frequency of $1.0$ means half the sample rate, also called the Nyquist frequency.

Implementation

Figure 2 shows the block diagram of the first-order digital lowpass filter. The block marked $z^{-1}$ is a delay of one sample, usually called the ‘filter state’.

The following code show how to implement the filter in C. It plots the first 10 samples of the impulse response.

#include <stdio.h>

typedef struct
{
    float coeff;
    float state;
} iir1;

void init_iir1(iir1 *filt, float c)
{
    filt->coeff = c;
    filt->state = 0.0f;
}

/* calculate a new output and update the filter */
float update_iir1(iir1 *filt, float input)
{
    filt->state = filt->coeff*input + (1.0f - filt->coeff)*filt->state;
    return filt->state;
}

int main()
{
    const float coeff = 0.5f;
    
    iir1 filter1;
    init_iir1(&filter1, coeff);

    printf("Filter coefficient = %f\n", coeff);
    /* print the first 10 samples of the impulse response */
    float in = 1.0f;
    for(int i=0; i<10; i++)
    {
        float out = update_iir1(&filter1, in);
        printf("[n=%d] in=%f out=%f\n", i, in, out);
        in = 0.0f;
    }
    
    return 0;
};

As can be seen from the source code above, the actual filter update function ‘update_iir1’ is quite simple.

The resulting output:

Filter coefficient = 0.500000
[n=0] in=1.000000 out=0.500000
[n=1] in=0.000000 out=0.250000
[n=2] in=0.000000 out=0.125000
[n=3] in=0.000000 out=0.062500
[n=4] in=0.000000 out=0.031250
[n=5] in=0.000000 out=0.015625
[n=6] in=0.000000 out=0.007812
[n=7] in=0.000000 out=0.003906
[n=8] in=0.000000 out=0.001953
[n=9] in=0.000000 out=0.000977

Cutoff frequency

In this section we present the cutoff frequency determination and it contains a lot of mathematics. If this isn’t your thing, skip to the end where the final answer is given.

The z-domain transfer function of the filter is: $$H(z) = \dfrac{c}{1+ (1-c)\cdot z^{-1}}$$

A good start to get the magnitude response $|H(z)|$ is to multiply $H(z)$ with its complex conjugate $H^{\ast}(z)$, as $ H(z) \cdot H^{\ast}(z)$ equals $|H(z)|^2$. $$|H(z)|^2 = \dfrac{c^2}{1 - 2(1-c)\Re \lbrace z^{-1} \rbrace + (1-c)^2|z^{-1}|^2},$$ where $\Re \lbrace x \rbrace$ takes the real part of $x$.

Checking the result

The filter should pass DC unaltered, so we can check if our formula is correct by pluging in $z=1$. We expect the magnitude $|H(z)|$ to be ‘1’ and therefore $|H(z)|^2$ to be 1, when $c>0$.

$$|H(z=1)|^2 = \dfrac{c^2}{1 - 2(1-c)\Re \lbrace 1 \rbrace + (1-c)^2|1|^2}$$

Which can be written as:

$$|H(z=1)|^2 = \dfrac{c^2}{1 - 2(1-c) + (1-c)^2}$$

Expanding $-2(1-c)$, we get:

$$|H(z=1)|^2 = \dfrac{c^2}{1 - 2 + 2c + (1-c)^2}$$

And because $(1-c)^2 = 1 - 2c + c^2$, we obtain:

$$|H(z=1)|^2 = \dfrac{c^2}{1 - 2 + 2c + 1 -2c + c^2}$$

Which, by elimination of terms in the denominator, simplifies to:

$$|H(z=1)|^2 = \dfrac{c^2}{c^2}$$

So:

$$|H(z=1)|^2 = 1$$

Which at proves that our formula is at least correct for DC.

Determining cutoff when the coefficient ‘c’ is given

The -3 dB cutoff frequency corresponds to $z=e^{j \theta}$ for which $|H(z)| = \frac{1}{\sqrt{2}}$, or equivalently $|H(z)|^2 = \frac{1}{2}$. Plugging $z=e^{j \theta}$ into our formula for $|H(z)|^2$ gives:

$$|H(e^{j \theta})|^2 = \dfrac{c^2}{1 - 2(1-c)\Re \lbrace \exp(-j \theta) \rbrace + (1-c)^2|1|^2}$$

Given Euler’s identity $e^{j \theta} = \cos(\theta) -j\cdot \sin(\theta)$, we can write $\Re \lbrace \exp(-j \theta) \rbrace = \cos(\theta)$, and get:

$$|H(e^{j \theta})|^2 = \dfrac{c^2}{1 - 2(1-c)\cos(\theta) + (1-c)^2}$$

And therefore:

$$|H(e^{j \theta})|^2 = \dfrac{c^2}{2 - 2(1-c)\cos(\theta) - 2c + c^2}$$

To find the cutoff frequency $\theta$ (in radians) given $c$, we must have $|H(e^{j \theta})|^2 = \frac{1}{2}$ as stated earlier. This means that the left part of the denominator $2 - 2(1-c)\cos(\theta) - 2c$ must equal $c^2$ so we get: $$|H(e^{j \theta})|^2 = \dfrac{c^2}{c^2 + c^2}$$

Thus we continue with $2 - 2(1-c)\cos(\theta) - 2c = c^2$ and try so solve for $\theta$.

$$-2(1-c)\cos(\theta) = c^2 + 2c - 2$$ $$\cos(\theta) = -\dfrac{c^2 + 2c - 2}{2(1-c)}$$ $$\theta = \arccos \left( -\dfrac{c^2 + 2c - 2}{2(1-c)} \right) $$

Going from radians to normalised frequency: $$f_{cut} = \dfrac{\theta}{\pi} $$

Now we can finally calculate the cutoff frequencies for our list of coefficients:

Determining the coefficient ‘c’ when the cutoff is given

Given a cutoff frequency $f_{c}$ in Hertz and the sample rate $f_{s}$ in samples per second, the normalised cutoff frequency $f$ is simply: $$f = 2\dfrac{f_{c}}{f_{s}},$$ and $\theta$ is then: $$\theta = \pi f \enspace \Longrightarrow \enspace \theta = 2\pi\dfrac{f_{c}}{f_{s}}$$

Again we use the formula for $|H(z)|^2$ but now we know $\theta$ and want to determine $c$. For the -3 dB frequency, we plug in $z=e^{-j\theta}$ and require that $|H(z=e^{-j\theta})|^2 = \frac{1}{2}$. From earlier, we have: $$|H(e^{j \theta})|^2 = \dfrac{c^2}{1 - 2(1-c)\cos(\theta) + (1-c)^2}$$ We must then solve: $$\dfrac{c^2}{1 - 2(1-c)\cos(\theta) + (1-c)^2} = \frac{1}{2}$$

This seems like a daunting task but remember that $\cos(\theta)$ is already known and is therefore a constant. A good starting step is to expand all the terms involving powers of ‘c’.

$$\dfrac{c^2}{1 - (2-2c)\cos(\theta) + 1 -2 c + c^2} = \frac{1}{2}$$

$$\dfrac{c^2}{2 - 2\cos(\theta) +2c\cos(\theta) -2c + c^2} = \frac{1}{2}$$

Now we can collect all the powers of ‘c’:

$$\dfrac{c^2}{2(1 - \cos(\theta)) -2c(1-\cos(\theta)) + c^2} = \frac{1}{2}$$

Multiplying both sides of the equation by the denominator gives: $$ c^2 = {(1 - \cos(\theta))-c(1-\cos(\theta)) + \frac{1}{2}c^2} $$

$$ 0 = {(1 - \cos(\theta))-c(1-\cos(\theta)) - \frac{1}{2}c^2} $$

$$ {2(1 - \cos(\theta))-2c(1-\cos(\theta)) - c^2} = 0 $$

$$ {c^2 + 2c(1-\cos(\theta)) - 2(1 - \cos(\theta))} = 0 $$

Remember, $\cos(\theta)$ is constant so $1 - \cos(\theta)$ is too. so we can solve it using the ABC formula, with: $$ \hat{a} = 1 $$ $$ \hat{b} = 2(1-\cos(\theta)) $$ $$ \hat{c} = -2(1-\cos(\theta)) $$

Then: $$ c = \dfrac{-\hat{b} \pm \sqrt{\hat{b}^2 -4\hat{a}\hat{c}}}{2\hat{a}} $$

And because $\hat{c} = -\hat{b}$, we actually have: $$ c = \dfrac{-\hat{b} \pm \sqrt{\hat{b}^2 +4\hat{b}}}{2} $$

We could go further unfolding this but there is hardly any simplification to be gained.

Finale

To summerise, we have the following set of formulas to determine the filter coefficient c: $$ \theta = 2\pi\dfrac{f_c}{f_s} $$ To calculate $\theta$ based on the cutoff frequency in Hertz and the sample rate.

Then we calculate our ABC formula coefficient $\hat{b}$: $$ \hat{b} = 2(1 - cos(\theta)) $$

And evaluate the simplified ABC formula to get the filter coefficient ‘c’: $$ c = \dfrac{-\hat{b} + \sqrt{\hat{b}^2 +4\hat{b}}}{2} $$