The article on Z-Transforms introduced a difference equation for discrete stable causal Linear Time Invariant (LTI) systems, that from here on we will refer to as a LTI system, or system for short \(\)
$$ \begin{align} \sum_{k=0}^N a_k\,y[n-k]&= \sum_{k=0}^M b_k\,x[n-k]\quad\Rightarrow \nonumber \\[10mu] a_0y[n]+a_1y[n-1]+a_2y[n-2]+\ldots &= b_0x[n]+b_1x[n-1]+b_2x[n-2]+\ldots \label{eq:diffequation} \end{align}\nonumber $$
The article Discrete Transfer Functions showed us the discrete transfer function \(H(z)\) for causal LTI systems. $$ \begin{align} H(z) &= \frac{b_0+b_1z^{-1}+b_2z^{-2}+\cdots+b_Mz^{-M}}{a_0+a_1z^{-1}+a_2z^{-2}+\cdots+a_Nz^{-N}} \label{eq:tf_polynomial} \\[10mu] &= K\,z^{\small N-\small M}\frac{(z-q_1)(z-q_2)\cdots(z-q_{\small M})}{(z-p_1)(z-p_2)\cdots(z-p_{\small N})},&K=\frac{b_M}{a_N} \label{eq:tf_factors} \end{align} $$
Here we will evaluate the response of discrete transfer functions to sinusoidal inputs, introduce stability criteria and give methods to transfer a response back to the time-domain.
We will following the notation used in our piece on Z Transforms, where \( \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \ztransform\) denotes a unilateral Z-transform, equivalent to the more common notation \(\mathfrak{Z}\left\{\,f[n]\,\right\}\), and \(f[n]\) is defined as the sample taken at time \(nT\) or \(f(nT)\). The terms filter and system will be used interchangeably.
Filter Types
In general, the term electronic filters refers to circuits that perform signal processing to remove unwanted frequency components from a signal and/or to enhance wanted components. Examples include enhancing X-ray images at airports and extracting radio signals from far away space probes. An electronic control system typically refers to a circuit that processes one signal into another to give the desired system response. E.g. modern thermostats learn the characteristics of your house, and chemical plants measure fluid levels to control flow pumps. In the context of this writing we refer to filters and systems interchangeably.
In equation \(\eqref{eq:diffequation}\), the \(b_i\) coefficients are called feedforward coefficients, and the \(a_i\) coefficients are called feedback coefficients.
We classify filters based on whether or not they use any previous value of the output, in what case we say they have feedback. Based on this we classify filters into two groups: finite impulse response and infinite impulse response filters.
Finite Impulse Response
As the name implies, Finite Impulse Response (FIR) filters have a finite response to an input. If the filter order is \(M\), then the maximum delay for the input to the output will be \(M\) samples. In other words, given an impulse input, the output will return to \(0\) after \(M\) samples.
In FIR filters, the filter output does not depend on any previous value of the output and the coefficients \(a_i=0\) for all \(i\gt 0\) in equation \(\eqref{eq:diffequation}\).
The figure below shows the signal flow of a FIR filter, where the \(z^{-1}\) block represents an one sample delay.
Output \(y[n]\)depends on | \(x[n-M]\cdots x[n]\) |
Impulse response | Has a finite duration |
Coefficients | \(b_i\) |
Infinite Impulse Response
A filter is said to be recursive when \(a_i\neq 0\) for some \(i\gt 0\). Recursive filters are also called Infinite Impulse Response (IIR) filters.
Output \(y[n]\) depends on | \(x[n-M]\ldots x[n],\ y[n-N]\ldots y[n]\) |
Impulse response | Duration depends on feedback |
Coefficients | \(a_i, b_i\) |
We may view an IIR filter \(H(z)\) as a series combination of two subsystems \(H_1(z)\) and \(H_2(z)\).
The commutative property of multiplication, allows the order of the subsystems to be reversed. If we draw such a circuit, it becomes apparent that each delay element \(z^{-1}\) is next to another delay element with the same input. We may then replace each such set of delay elements by one delay element with the same input. The resulting signal flow graph is called “Direct Form II” as depicted below.
Impulse response
The response of a system to any input can be calculated by the time convolution or the frequency product of the impulse response of the system and the input signal.
We derived the Z-Transform for the impulse function \(\delta[n]\) as
$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \delta[n] \ztransform 1 \triangleq \Delta(z)\nonumber \label{eq:impulse} $$
When applying this input \(X(z)\) to the filter with impulse response \(H(z)\), the output \(Y(z)\) is
In other words, the response to an impulse input, is simply the transfer function \(H(z)\) itself. That is not very surprising considering that \(h[n]\) was defined as the response to an impulse input function.
Since the transfer function \(H(z)\) equals the impulse response \(Y(z)\) of the transfer function, these terms are often used interchangeably.
Frequency response
Discrete Transfer Functions introduced the concept of poles and zeros and their effect on the transfer function. Here we will take it a step further by evaluating specific \(z\) values.
The frequency response of a linear time invariant system is defined as the steady state response to a sinusoidal input. In other words the output after all transients have died out.
As part of the Z-transform, we defined:
$$ z\triangleq\mathrm{e}^{sT} \nonumber $$ where \(s=\sigma+j\omega\)
To find the frequency response, we follow the same methodology as we did for the Continuous Frequency Response and evaluate the expression \(F(z)\) along \(s=j\omega\)
Substitute \(\eqref{eq:zunitcircle}\) in \(\eqref{eq:tf_factors}\)
Expressing equation \(\eqref{eq:tf_unitcircle}\) in polar form, helps us distinguish between the amplitude and phase response
Remember, that in the \(z\)-plane, angular frequency are shown in normalized form, where the normalized angular frequency \(\omega T\) is the angle with the positive horizontal axis.
Amplitude response
We find the amplitude response as the magnitude \(|H(\mathrm{e}^{j\omega T})|\) when substituting equation \(\eqref{eq:tf_unitcircle}\) in \(\eqref{eq:tf_polarform}\)
With \(\left|\mathrm{e}^{j(\small N-\small M)\omega T}\right|=1\), according to Euler’s formula, the amplitude response follows as
The amplitude response can be visualized with the length of vectors from the poles and zeros to point \(z\) on the unit circle that corresponding to the natural frequency for which the function is evaluated.
The product of the vector lengths from the zeros divided by that of the poles times \(K\) represents the amplitude response for that normalized angular frequency.
Phase response
Along the same line of thought, the phase response is the angle \(\angle\left(H(\mathrm{e}^{j\omega T})\right)\) when substituting equation \(\eqref{eq:tf_unitcircle}\) in \(\eqref{eq:tf_polarform}\)
With \(K\) a real valued scalar and \(\angle\mathrm{e}^{j\phi}=\phi\), according to Euler’s formula
The phase response \(\angle H(\mathrm{e}^{j\omega T})\) follows as
The phase response can be visualized using the angle of vectors from the poles and zeros to point \(z\) compared to a horizontal line. Where \(z\) is a point on the unit circle (\(|z|=1\)) for which the function is evaluated. The phase response follows as is the sum of the angles from the zeroes minus that of the poles plus \((N-M)\omega T\).
Inverse Z-transform
Eventually there comes a time to return to the time-domain using an inverse Z-transform. The article Z-Transforms metioned some of the techniques for the inverse Z-Transform:
- using the binomial theorem
- using the convolution theorem
- performing long division
- using the initial-value theorem
- partial fractions expansion
We already used the binomial theorem to proof the binomial scaled pair. Here we use long division to reduce the order of the numerator, and use partial fraction expansion to split up the remaining fraction.
Partial Fraction Expansion (and long division)
As we have seen in equation \(\eqref{eq:tf_polynomial}\), the impulse response \(Y(z)=H(z)\) is a rational fraction with \(N\) poles and \(M\) zeroes
$$ \begin{align} Y(z) &= \frac{b_0+b_1z^{-1}+b_2z^{-2}+\ldots+b_Mz^{-M}}{a_0+a_1z^{-1}+a_2z^{-2}+\ldots+a_Nz^{-N}},&a_0=1 \nonumber \end{align} \nonumber $$
This rational fraction is proper when the degree of the numerator polynomial is less than the degree of the denominator polynomial. To make the function proper, we use long division of the denominator/nominator until the order of the numerator is less than that of the denominator.
Let’s call the quotient from the long division \(F(z)\), and the ratio of the remainder/denominator \(G(z)\)
Determine \(F(z)\), the FIR part
The \(F(z)\) part will be a polynomial in \(z^{-1}\) of the order \(M-N\).
Recall the delay from the Z-transform pairs
$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \delta[n-a]\,\ztransform\, z^{-a}\nonumber $$
Using this transform, \(F(z)\) transforms to a parallel combination of delayed impulses in the time-domain
Determine \(G(z)\), the IIR part
The proper rational function \(G(z)\) is the remainder/denominator of the long division. Note that the numerator uses \(\dot{b}_i\) coefficients and the constant \(\dot{K}\) brings the numerator and denominator in unity form
This proper fraction can be split up into a sum of simpler fractions as introduced by Oliver Heaviside and described in Partial Fraction Expansion (PFE).
At this point in time we need to decide on the format of the time-domain function \(g[n]\). We will consider \(G(z)\) with singular poles and show you three forms. The same can be done with multiple poles, but it is a bit more involved.
While all formats lead to the same output sequence, some may be more intuitive than others. For example, if you expect a exponentially decaying response, you may want to work towards that and see how well it matches.
Choice 1: The obvious
When \(G(z)\) only has single poles, PFE gives a summation of partial fractions in the form \(\frac{c}{z-a}\)
Thanks to the linearity property of the Z-transform, these the simpler fraction can each be transformed to and summed up in the time-domain. The Z-transform for the scaled delay pair is then found in the table of Z-transform pairs as
$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \begin{align} a^{n-1}\gamma[n-1] \ztransform \dfrac{1}{z-a}, & & |z|\gt|a| \nonumber \end{align} \nonumber $$
Therefore, \(G(z)\) transforms to a parallel combination of delayed scaled step functions in the time-domain
Choice 2: Work towards \(\frac{z}{z-a}\) partial fractions
This time, we decide to work towards partial fractions in the form \(\frac{z}{z-a}\) that transform to \(a^n\,\gamma[n]\) in the time-domain.
Once more, start with equation \(\eqref{eq:gfactors}\), but this time preserve a power of \(z\) by divide both sides by \(z\).
When \(\frac{G(z)}{z}\) has only single poles, according to Heaviside, it expands to the summation
The Z-transform for the constant and the the scaled pair are found in the table of Z-transform pairs as
$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \begin{align} \small{a\delta[n]\triangleq\begin{cases}a,&n=0\\0,&n\neq0\end{cases}} \ztransform a \nonumber \\ a^{n}\gamma[n] \ztransform \dfrac{z}{z-a},&&|z|\gt|a| \nonumber \end{align} \nonumber $$
Thus \(G(z)\) transforms to a parallel combination of impulse and scaled step functions in the time-domain
Choice 3: delayed form
In this variation, the IRR part begins after the FIR part has finished. This can be more accurate in signal modeling applications, as the IIR part may be delayed so that its impulse response begins where that of the FIR part died out.
Start at the beginning with equation \(\eqref{eq:gfactors}\), and multiply the numerator with the highest power of \(z^{-1}\): \(z^{\small M}\)
Similar to before, split the expression in \(F(z)\) and \(G(z)\) parts using long division. Call the quotient from the long division \(F(z)\), and the ratio of the remainder/denominator \(G(z)\)
The \(F(z)\) part will be a polynomial in \(z\) of the order \(M-N\).
Determine the \(G(z)\) part, by first bringing it back to a polynomial in \(z^{−1}\)
Split up the proper fraction into a sum of simpler fractions using PFE, using the “obvious” choice and assuming \(G(z)\) only has single poles
Therefore, \(G(z)\) transforms to a parallel combination of delayed scaled step functions in the time-domain
Examples for these forms are given in the appendix.
Stability
As we have seen in Z-transforms, the convergence of a transfer function depends on its magnitude, while its phase has no effect. A system is stable if the magnitude of its impulse response \(h[n]\) decays to \(0\) as \(t\to\infty\).
As we have seen, every finite-order LTI filter can be expressed as FIR and IIR parts. We will now examine how these parts contribute to the stability of the system.
FIR part
The FIR part \(F(z)\) from equation \(\eqref{eq:firpart5}\) is a finite-order polynomial in \(z^{-1}\)
In the time-domain this transforms to \(\eqref{eq:firpart}\)
$$ f[n]=\sum_{k=0}^{M-N}f_k\ \delta[n-k]\nonumber $$This is always stable because there are finite terms.
IIR part
As shown in Choice 1 above, the IIR part \(\eqref{eq:choice2}\) can be expressed as a summation of \(\frac{1}{z-a}\) terms \(\eqref{eq:choise2b}\). These terms transform to the time-domain as \((\ref{eq:iir1},\ref{eq:iir0})\).
If all the poles are inside the unit circle in the \(z\)-plane, then the IIR part is stable and consequently the transfer function \(H(z)\) is stable. More formally:
An irreducible transfer function \(H(z)\) is stable if and only if all its poles have a magnitude less than one.
Appendix
Example 1
This example shows an Inverse Z-Transform for a rational function where the numerator and denominator have the same degree \(N=M=3\).
This conveniently matches Example 2 in the Partial Fraction Expansion article.
$$ \frac{\color{green}{2}x^3+\color{green}{x}^2\color{green}{-1}x+\color{green}{4}}{(x-2)^3}=\color{blue}{-\frac{1}{2}}+\color{blue}{11}\frac{x}{(x-2)^3}+\color{blue}{8}\frac{x}{(x-2)^2}+\color{blue}{\frac{5}{2}}\frac{x}{x-2}\nonumber $$
This implies that \(Y(z)\) can be expressed in partial fractions as
These terms are readily found in the Z-transform pairs table
$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \begin{align} \delta[n] &\,\ztransform\, 1 \nonumber \\ a^n\,\gamma[n] &\,\ztransform \frac{z}{z-a},&|z|\gt |a| \nonumber\\[6mu] n\,a^n\,\gamma[n] &\,\ztransform \frac{az}{(z-a)^2},&|z|\gt |a| \nonumber \\ \tfrac{1}{2}{n(n-1)}\,a^n\,\gamma[n] &\,\ztransform \frac{a^2z}{(z-a)},&|z|\gt |a| \nonumber \end{align} \nonumber $$
Using these transform pairs, the time-domain response is
Example 2
This example shows an Inverse Z-Transform for a rational function where the degree of the numerator is 1 more as that of the denominator \(M=3, N=2\). [CCRMA]
Examine the impulse response of a filter with transfer function
Solve it by using long division to bring the order of the numerator down to \(N-1\), so we can use partial fraction expansion on the remaining IIR part.
To help with notation, define \(d\triangleq z^{-1}\) Do the long division
This bought the order of the numerator down to one less than that of the denominator (\(M=1, N=2\))
Examine the IIR part, \(G(z)\)
This matches Example 1 in the partial fraction expansion article.
$$ \frac{\color{green}{-8}+\color{green}{24}x}{1-2x+x^2}=\frac{\color{blue}{-24}}{1-x}+\frac{\color{blue}{16}}{\left(1-x\right)^{2}}\nonumber $$
Analogue to that example, fraction \(G(z)\) can be expressed in partial fractions as
With \(H(z) = F(z)+G(z)\), and the impulse response \(Y(z) = \Delta(z)\,H(z)=H(z)\)
These terms are readily found in the Z-transform pairs table
$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} \begin{align} \delta[n] &\ztransform 1 \nonumber \\ \gamma[n-a]\,\color{grey}{\gamma[n]} &\ztransform \small{\begin{cases}z^{-a},&a\geq0\\0,&a\lt0\end{cases}} &z\neq0 \nonumber \\ a^n\,\color{grey}{\gamma[n]} &\ztransform {\frac{1}{1-az^{-1}}},&|z|\gt|a| \nonumber \\ \small{\left(\begin{array}{c}n+m-1\\m-1\end{array}\right)}\,a^n\,\color{grey}{\gamma[n]} &\ztransform\dfrac{1}{(1-az^{-1})^m}&|z|\gt|a| \nonumber \\ n\,\color{grey}{\gamma[n]} &\ztransform \frac{z^{-1}}{\left(1-z^{-1}\right)^2},&|z|\gt1 \nonumber \end{align}\nonumber $$
The impulse response \(y[n]\) follows as
Delayed version
The delayed version is found by first multiplying numerator and denominator of equation \(\eqref{eq:example2_def}\) with \(z^3\), to make them a polynomial in \(z\) instead of \(z^{-1}\)
Do a long division to reduce the order of the numerator
This bought the order of the numerator down to less than that of the denominator (M=1,N=3)
Examine the IIR part \(G(z)\), by bringing it back to a polynomial in \(z^{-1}\)
Partial fraction expansion, left as an exercise to the reader
With \(H(z)=F(z)+G(z)\), and the impulse response \(Y(z)=\Delta(z)\,H(z)=H(z)\)
Using GNU/Octave
GNU/Octave residuez
function returns the FIR part as \(\color{purple}{f}\), the filter-pole residues as \(\color{blue}{r}\), the filter poles as \(\color{brown}{p}\) and the pole multiplicities as \(\color{magenta}{m}\)
pkg load signal
B=[2 6 6 2]; A=[1 -2 1];
[r,p,f,m] = residuez(B,A)
r =
-24
16
p =
1
1
f =
10 2
m =
1
2
In other words
We can also use GNU/Octave to determine the delayed form of the IIR
For the same example, the residued
function returns
[r,p,f,m] = residued(B,A)
r =
8
16
p =
1
1
f =
2 10
m =
1
2
In other words