Implementation of a filters

The discrete-time transfer function H(z) is often obtained from its counterpart, the continuous-time transfer function H(s) via discretization. A discrete-time transfer function has the following form:

(1)   \begin{equation*}H(z) = \frac{Y(z)}{X(z)} = \frac{\sum_{i = 0}^{N} b_i z^i}{\sum_{j = 0}^{M} a_j z^j}\end{equation*}

Herein, X(z) is the input and Y(z) is the output of the system, N and M are the degree of the numerator and denominator, respectively. Where N \leq M, meaning we are dealing with a proper transfer function. While (1) is valid for any order it is not recommended to directly use transfer functions of high order. These can namely introduce numerical problems very quickly. Rather factorize the numerator and denominator of (1) into a cascade of first and second order polynomials.

    \begin{equation*}H(z) = K \frac{\displaystyle \prod_{i = 0}^V \left( z + b_{0i} \right)}{\displaystyle \prod_{i = 0} ^W\left( z + a_{0i} \right)} \frac{\displaystyle \prod _{i = 0}^N \left( z^2 + b_{1i} z + b_{2i} \right)}{\displaystyle\prod_{i = 0}^M \left( z^2 + a_{1i} z + a_{2i} \right)}\end{equation*}

Now let us look at the simple discrete-time transfer function of order two:

    \begin{equation*}H(z) = \frac{b_0 z^2 + b_1 z + b_2}{z^2 + a_1 z + a_2}\end{equation*}

This function is non-causal, because it depends on future inputs. Therefor, both the numerator and denominator are multiplied by reciprocal of the highest order of z occuring in the denominator, in this case z^{-2}, to make the system causal. Hence, we obtain:

    \begin{equation*}H(z) = \frac{b_0 + b_1 z^{-1} + b_2 z^{-2}}{1 + a_1 z^{-1} + a_2 z^{-2}}\end{equation*}

Followingly, using the linearity and time-shifting properties of the \mathcal{Z}-transform, i.e., \mathcal{Z}(a_1 x_1[n] + a_2 x_2[n]) = a_1 X_1(z) + a_2 X_2(z) and \mathcal{Z}(q^{-k}x[n]) = z^{-k}X(z), we obtain the difference equation. Remark that the shift operator q is defined as q^kx[n] = x[n+k], the forward shift operation and q^{-k}x[n] = x[n-k], the backward shift (delay) operator. As a result we obtain,

(2)   \begin{equation*}y[n] = \frac{b_0 + b_1 q^{-1} + b_2 q^{-2}}{1 + a_1 q^{-1} + a_2 q^{-2}} x[n]\end{equation*}

Rewriting (2) gives us:

    \begin{gather*}\left( 1 + a_1 q^{-1} + a_2 q^{-2} \right) y[n] = \left( b_0 + b_1 q^{-1} + b_2 q^{-2} \right) x[n] \\\Leftrightarrow y[n] + a_1 y[n-1] + a_2 y[n-2] = b_0 x[n] + b_1 x[n-1] + b_2 x[n-2] \\\Leftrightarrow y[n] = b_0 x[n] + b_1 x[n-1] + b_2 x[n-2] - a_1 y[n-1] - a_2 y[n-2]\end{gather*}

This last equation is the difference equation which we can easily implement on our digital platform. Numerous methods exists on how to implement a filter. Four of these methods are closely related to each other. These are:

  • Direct form I
  • Direct form II
  • Transposed direct form I
  • Transposed direct form II

Direct form I

The direct form I is an FIR filter followed by an IIR filter. That is to say, it implements Y(z) followed by \frac{1}{X(z)}

(3)   \begin{equation*}y[n] = b_0 x[n] + b_1 x[n-1] + b_2 x[n - 2] - a_1 y[n-1] - a_2 y[n-2]\end{equation*}

In an algorithm you can implement it as:

yk = b0 * xn + b1 * x1 + b2 * x2 - a1 * y1 - a2 * y2
x2 = x1
x1 = xn
y2 = y1
y1 = yn

Where, x1, y1, x2 and y2 are the four state variables.

Direct form II

The direct form I is an IIR filter followed by an FIR filter. Which implements \frac{1}{X(z)} followed by Y(z).

Signal flow diagram of direct form II. Where s[n] a state variable

(4)   \begin{align*}s[n] &= x[n] - a_1 s[n-1] - a_2 s[n-2] \\y[n] &= b_0 s[n] + b_1 s[n-1] + b_2 s[n-2]\end{align*}

In an algorithm you can implement it as:

s0 =      xn - a1 * s1 - a2 * s2
yn = b0 * s0 + b1 * s1 + b2 * s2
s2 = s1
s1 = s0

Where, s0, s1, s2 are the three state variables.

Direct form I transposed

Both direct forms can be converted to an equivalent transposed form via:

  • Reverse direction of each interconnection
  • Reverse direction of each multiplier
  • Change junctions to adders and vice-versa
  • Interchange the input and output signals
Signal flow diagram of direct form I transposed.

(5)   \begin{align*}y[n] &= b_0 v[n] + b_1 v[n-1] + b_2 v[n-2] \\v[n] &= x[n] - a_1 v[n-1] - a_2 v[n-2] \\\end{align*}

In an algorithm you can implement it as:

vn = xn + s2
yn = s4 + b0 * vn
s4 = s3 + b1 * vn
s3 = b2 * vn
s2 = s1 - a1 * vn
s1 = - a2 * vn

Note that this form is very inefficient. You can easily shift the delays to the center, obtaining the direct form II.

Direct form II transposed

Signal flow diagram of direct form II transposed.

(6)   \begin{align*}s_1[n] &= b_2 x[n-1] - a_2 y[n-1] \\s_2[n] &= b_1 x[n-1]  - a_1 y[n-1] + s_1[n-1] \\y[n] &= b_0 x[n] + s_2[n]\end{align*}

In an algorithm you can implement it as:

yn = s2 + b0 * xn
s2 = s1 + b1 * xn - a1 * yn
s1 = b2 * xn - a2 * yn

Symbol to represent Emergency Response Team (ERT) Bedrijfshulpverlening (BHV)

Emergency response officers (ERT), or in Dutch Bedrijfshulpverlening (BHV), within companies are trained in providing both first aid and limiting and combating (small) fires. All to limit the consequences of accidents. Various icons are used to resemble these organisations. One symbol which I have seen quiet often is the combination of first aid and a fire symbol.

However, since there was not yet a vector format available, I created one. Using the established first aid and fire sign from ISO 7010.

Emergency Response Team (ERT), or in Dutch Bedrijfshulpverlening (BHV), icon.

On the left the fire symbol with the red color. On the right half of the first aid cross with the green color. The red and green color form a gradient in the middle.

You can download the vector file here.

Filters in Control Systems

In this blog I want to elaborate on the various types of filters that are being used in Control Systems. These filters are often used in combination with a PI, PD or PID controller to obtain a robust controller. The filters which we will discuss are:

  • Low-Pass filter
  • High-Pass filter
  • Lead-Lag filter
  • Notch filter

Using these four filters we can create other filter types, such as a Band-Stop or Band-Pass filter. The behavior of each filter can be captured by a transfer function in the continuous-time using the Laplace domain s of either first and/or second order. We use the transfer functions to describes the filter H(s) = Y(s) / U(s) and as such the relation between the input U(s) and the output Y(s). Throughout this blog we will write the equations in the form of their angular frequency \omega in [rad/s], also known as radian frequency. However, we will specify filters using their frequency f in [Hz]. Using the radial frequency notation results in a more visual compact formula. Additionally, it is also possible to specify the filters in form of their time constant \tau in [s]. The following relation holds between the angular frequency, frequency and time constant.

    \begin{equation*} \omega \triangleq 2\pi f \triangleq \frac{1}{\tau} \end{equation*}

Many forms are used within literature, one book will use angular frequencies, the other will use time-constants. Finally, for second order filters, the only filter with possibly complex poles or zeros, can be written in various ways. We will specify second order filters in terms of the damping \beta of the corresponding frequency. It is also possible to describe these formulas using the Quality factor Q. Whereas \beta describes how oscillations decay in a system after a disturbance, Q describes how underdamped the system is. The following relation holds between \beta and Q

    \begin{equation*} \beta \triangleq \frac{1}{2Q} \end{equation*}

Low-pass filter

A low-pass filter is used to pass signals with a frequency lower than a certain cut-off frequency f_\text{lp}. Below the formulas for both the first- and second-order low-pass filter is given.

    \begin{equation*} H(s) = K\cdot\frac{\omega_{\text{lp}}}{s + \omega_{\text{lp}}} \quad\quad H(s) = K\cdot\frac{\omega_{\text{lp}}^2}{s^2 + 2\beta\omega_{\text{lp}}s + \omega_{\text{lp}}^2} \end{equation*}

Herein, K denotes the gain, f_\text{lp} denotes the low-pass cut-off frequency and \beta_\text{lp} denotes the damping. Whereas the first-order supresses with 20 [dB/dec], the second-order supresses with 40 [dB/dec]. Low frequent the filter gain is K.

High-pass filter

The complement of a low-pass filter is a high-pass filter. This filter is used to pass signals with a frequency higher than a certain cut-off frequency f_\text{hp}. Below the formulas for both the first- and second-order high-pass filter is given.

    \begin{equation*} H(s) = K\cdot\frac{s}{s + \omega_{\text{hp}}} \quad\quad H(s) = K\cdot\frac{s^2}{s^2 + 2\beta\omega_{\text{hp}}s + \omega_{\text{hp}}^2} \end{equation*}

Herein, K denotes the gain, f_\text{hp} denotes the high-pass cut-off frequency and \beta_\text{hp} denotes the damping. Likewise as the low-pass filter the first-order supresses frequencies with 20 [dB/dec] and the second-order with 40 [dB/dec].

Lead-lag filter

A lead-lag filter, also known as a lead-lag compensator, is often mainly used for phase compensation rather then magnitude. Below the formula for a lead or lag filter is shown.

    \begin{equation*} H(s) = K\cdot\frac{\omega_\text{p}}{\omega_\text{z}}\cdot\frac{s + \omega_\text{z}}{s + \omega_\text{p}} \end{equation*}

Herein, f_\text{p} and f_\text{z} denote the frequency of the pole and zero, respectively. The filter functions as a lead filter if f_\text{p} > f_\text{z}  and otherwise as a lag filter. The filter has its maximum or minimum phase at \sqrt{f_\text{p}f_\text{z}}. Finally, at f = \infty the filter has a gain of Kf_\text{p} / f_\text{z} or Kf_\text{z} / f_\text{p} in case of a lead or lag filter, respectively. Naturally, the filter can be cascaded with itself by which a the filter can be a lead and lag filter simultaneously.

Notch filter

A notch filter is often used to filter undesired resonance peaks. Below the formula for a notch filter is shown.

    \begin{equation*} H(s) = K\cdot\frac{\omega_\text{p}}{\omega_\text{z}}\cdot\frac{s^2 + 2\beta_\text{z}\omega_\text{z}s + \omega_\text{z}^2}{s^2 + 2\beta_\text{p}\omega_\text{p}s + \omega_\text{p}^2} \end{equation*}

Herein, f_\text{p} and f_\text{z} denote the frequency of the pole and zero, respectively. Likewise, \beta_\text{p} and \beta_\text{z} denote the damping of the pole and zero and K is as usual the gain. When f_\text{z} = f_\text{z} the notch will target one specific frequency. The gain at that frequency is given by \beta_\text{z} / \beta_\text{p}. When f_1 \neq f_2 the notch filter is also referred as a skewed notch and the difference between gain at low and high frequencies is given by (f_\text{p} / f_\text{z})^2.