Evaluating Discrete Transfer Functions using Z-Transforms

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.

Signal Flow Graph for FIR filter
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.

Signal Flow Graph, IIR filter, Direct Form I
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)$.

$$ H(z) = H_1(z)\,H_2(z)=H_2(z)\,H_1(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.

Signal Flow Graph, IIR filter, Direct Form II

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

$$ Y(z) = X(z)\,H(z) = 1\,H(z) $$

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.

$$ \shaded{ Y(z)=H(z) } $$

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$

$$ z = \left.\mathrm{e}^{sT}\right|_{s=\omega T}=\mathrm{e}^{j\omega T} \label{eq:zunitcircle} $$

Substitute $\eqref{eq:zunitcircle}$ in $\eqref{eq:tf_factors}$

$$ \begin{align} H(\mathrm{e}^{j\omega T}) &= K\,\mathrm{e}^{j(\small N-\small M)\omega T}\frac{(\mathrm{e}^{j\omega T}-q_1)(\mathrm{e}^{j\omega T}-q_2)\dots(\mathrm{e}^{j\omega T}-q_{\small M})}{(\mathrm{e}^{j\omega T}-p_1)(\mathrm{e}^{j\omega T}-p_2)\dots(\mathrm{e}^{j\omega T}-p_{\small N})},&K=\frac{b_M}{a_N} \label{eq:tf_unitcircle} \end{align} $$

this will only converge when the ROC includes the unit circle where $|z|=1$.

Expressing equation $\eqref{eq:tf_unitcircle}$ in polar form, helps us distinguish between the amplitude and phase response

$$ \shaded{ \begin{align} H\left(\mathrm{e}^{j\omega T}\right) &= r\,\mathrm{e}^{j\varphi} \nonumber \\ \text{where}\quad r&=\left|H\left(\mathrm{e}^{j\omega T}\right)\right|&\text{amplitude response} \nonumber \\ \text{and}\quad\varphi&=\angle{H\left(\mathrm{e}^{j\omega T}\right)}&\text{phase response} \nonumber \\ \end{align} } \label{eq:tf_polarform} $$

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}$

$$ \begin{align} \left|H(\mathrm{e}^{j\omega T})\right| =& \left|K\,\mathrm{e}^{j(\small N-\small M)\omega T}\frac{(\mathrm{e}^{j\omega T}-q_1)(\mathrm{e}^{j\omega T}-q_2)\dots(\mathrm{e}^{j\omega T}-q_{\small M})}{(\mathrm{e}^{j\omega T}-p_1)(\mathrm{e}^{j\omega T}-p_2)\dots(\mathrm{e}^{j\omega T}-p_{\small N})}\right|,&K=\frac{b_M}{a_N}\nonumber\\[14mu] =& \left|K\right|\,\left|\mathrm{e}^{j(\small N-\small M)\omega T}\right| \frac{\left|\mathrm{e}^{j\omega T}-q_1\right|\cdot\left|\mathrm{e}^{j\omega T}-q_2\right|\dots\left|\mathrm{e}^{j\omega T}-q_{\small M}\right|} {\left|\mathrm{e}^{j\omega T}-p_1\right|\cdot\left|\mathrm{e}^{j\omega T}-p_2\right|\dots\left|\mathrm{e}^{j\omega T}-p_{\small N}\right|},&K=\frac{b_M}{a_N} \end{align} $$

With $\left|\mathrm{e}^{j(\small N-\small M)\omega T}\right|=1$, according to Euler’s formula, the amplitude response follows as

$$ \begin{align} \shaded{ \left|H(\mathrm{e}^{j\omega T})\right| = |K|\,\frac{\prod_{i=1}^{M}\left|\mathrm{e}^{j\omega T}-q_i\right| } {\prod_{i=1}^{N}\left|\mathrm{e}^{j\omega T}-p_i\right|}}, & & K=\frac{b_M}{a_N} \end{align} $$

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}$

$$ \begin{align} \angle H(\mathrm{e}^{j\omega T}) =\,& \angle\left(K\,\mathrm{e}^{j(\small N-\small M)\omega T}\frac{(\mathrm{e}^{j\omega T}-q_1)(\mathrm{e}^{j\omega T}-q_2)\dots(\mathrm{e}^{j\omega T}-q_{\small M})}{(\mathrm{e}^{j\omega T}-p_1)(\mathrm{e}^{j\omega T}-p_2)\dots(\mathrm{e}^{j\omega T}-p_{\small N})}\right),&K=\frac{b_M}{a_N} \nonumber \\ =\,& \angle K+\angle\mathrm{e}^{j(\small N-\small M)\omega T} \nonumber \\ & +\angle\left(\mathrm{e}^{j\omega T}-q_1\right) + \angle\left(\mathrm{e}^{j\omega T}-q_2\right)+\dots +\angle\left(\mathrm{e}^{j\omega T}-q_{\small M}\right) \nonumber \\ & -\angle\left(\mathrm{e}^{j\omega T}-p_1\right) – \angle\left(\mathrm{e}^{j\omega T}-p_2\right) -\dots -\angle\left(\mathrm{e}^{j\omega T}-p_{\small N}\right),&K=\frac{b_M}{a_N} \nonumber \end{align} $$

With $K$ a real valued scalar and $\angle\mathrm{e}^{j\phi}=\phi$, according to Euler’s formula

$$ \begin{align} \angle H(\mathrm{e}^{j\omega T}) =&(N-M)\omega T\nonumber\\ &+\angle\left(\mathrm{e}^{j\omega T}-q_1\right) + \angle\left(\mathrm{e}^{j\omega T}-q_2\right)+\dots +\angle\left(\mathrm{e}^{j\omega T}-q_{\small M}\right) \nonumber \\ &-\angle\left(\mathrm{e}^{j\omega T}-p_1\right) – \angle\left(\mathrm{e}^{j\omega T}-p_2\right) -\dots -\angle\left(\mathrm{e}^{j\omega T}-p_{\small N}\right) \nonumber \end{align} $$

The phase response $\angle H(\mathrm{e}^{j\omega T})$ follows as

$$ \shaded{\angle H(\mathrm{e}^{j\omega T}) =(N-M)\omega T +\sum_{i=1}^{M}\angle\left(\mathrm{e}^{j\omega T}-q_i\right) -\sum_{i=1}^{N}\angle\left(\mathrm{e}^{j\omega T}-p_i\right) } $$

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)$

$$ Y(z) = F(z)+G(z) $$

if $Y(z)$ was already proper ($M\lt N$), we can skip the long division and set the term $F(z)$ to $0$.

Determine $F(z)$, the FIR part

The $F(z)$ part will be a polynomial in $z^{-1}$ of the order $M-N$.

$$ F(z) = f_0+f_1z^{-1}+f_2z^{-2}+\ldots+f_{\color{purple}{M-N}}z^{-(\color{purple}{M-N})}=\sum_{k=0}^{M-N}f_k\,z^{-k} \label{eq:firpart5} $$

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

$$ \begin{align} \shaded{f[n] = \sum_{k=0}^{M-N}f_k\ \delta[n-k]} & & \text{FIR part} \label{eq:firpart} \end{align} $$

we call this the FIR part, because it does not depend on any value of the output.

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

$$ \begin{align} G(z) &= \dot{K}\,\frac{1+\dot{b}_{1}z^{-1}+\dot{b}_{2}z^{-2}+\ldots+\dot{b}_{\color{red}{N-1}}z^{-(\color{red}{N-1})}} {1+a_{1}z^{-1}+a_{2}z^{-2}+\ldots+a_{\color{red}{N}}z^{-N}},&\dot{K}=\frac{\dot{b}_0}{\dot{a}_0} \nonumber \\[10mu] &= \dot{K}\,\frac{1+\dot{b}_1z^{-1}+\dot{b}_2z^{-2}+\ldots+\dot{b}_{\color{red}{N-1}}z^{-(\color{red}{N-1})}}{(z-r_1)(z-r_2)(z-r_3)\ldots(z-r_N)},&\shaded{\dot{K}=\frac{\dot{b}_0}{\dot{a}_0}} \label{eq:gfactors} \end{align} $$

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}$

$$ G(z) = \frac{c_1}{z-r_1}+\frac{c_2}{z-r_2}+\cdots+\frac{c_N}{z-r_N}=\sum_{k=1}^N\,\frac{c_k}{z-r_k} \label{eq:iir1} $$

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

$$ \shaded{ \begin{align} g[n]&=\left(c_1+c_2r_2+\cdots+c_{\small N}(r_{\small N})^{\,\small N-1}\right)\,\color{grey}{\gamma[n\color{black}{-1}]} \nonumber \\ &= \sum_{k=1}^{N}c_k(r_k)^{k-1}\, \color{grey}{\gamma[n\color{black}{-1}]} \nonumber \end{align} } \label{eq:iir0} $$

note the $-1$ in the step function $\gamma[n\color{black}{-1}]$. For the next choice, will prevent that delay.

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$.

$$ \begin{align} \frac{G(z)}{\color{blue}{z}}&=\dot{K}\,\frac{1+\dot{b}_1z^{-1}+\dot{b}_2z^{-2}+\ldots+\dot{b}_{\color{red}{N-1}}z^{-(\color{red}{N-1})}}{\color{blue}{z}(z-r_1)(z-r_2)(z-r_3)\ldots(z-r_N)} \label{eq:choice2} \end{align} $$

When $\frac{G(z)}{z}$ has only single poles, according to Heaviside, it expands to the summation

$$ \begin{align} \frac{G(z)}{\color{blue}z}&=\color{blue}{\frac{c_0}{z}}+\frac{c_1}{z-r_1}+\frac{c_2}{z-r_2}+\cdots+\frac{c_{\small N}}{z-r_N}=\color{blue}{\frac{c_0}{z}}+\sum_{k=1}^N\,c_k\,\frac{1}{z-r_k}\quad\Rightarrow \nonumber \\[6mu] G(z) &= \color{blue}{c_0}+\frac{c_1z}{z-r_1}+\frac{c_2z}{z-r_2}+\cdots+\frac{c_{\small N}z}{z-r_N}=\color{blue}{c_0}+\sum_{k=1}^N\,c_k\,\frac{\color{blue}{z}}{z-r_k} \end{align} $$

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

$$ \shaded{ \begin{align} g[n] &= c_0\,\color{grey}{\delta[n]}+\left(c_1r_1+c_2(r_2)^2+\cdots+c_{\small N}(r_{\small N})^{\,\small N}\right)\,\color{grey}{\gamma[n]} \nonumber \\ &= c_0\color{grey}{\delta[n]}+\sum_{k=1}^{N}c_k(r_k)^{\small N}\,\color{grey}{\gamma[n]} \nonumber \end{align} } \label{eq:choise2b} $$

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}$

$$ \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}} \nonumber \\[10mu] &= \color{blue}{z^{-M}}\,\frac{b_0z^{\color{blue}M}+b_1z^{\color{blue}{M}-1}+b_2z^{\color{blue}{M}-2}+\cdots+b_M}{a_0+a_1z^{-1}+a_2z^{-2}+\cdots+a_Nz^{-N}} \end{align} $$

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$.

$$ F(z) = f_0+f_1z+f_2z^{2}+\ldots+f_{\color{purple}{M-N}}z^{(\color{purple}{M-N})} = \sum_{k=0}^{M-N}f_k\,z^{k} $$

Determine the $G(z)$ part, by first bringing it back to a polynomial in $z^{−1}$

$$ \begin{align} G(z) &= \dot{K}\,\color{blue}{z^{-M}}\,\frac{\ddot{b}_0+\ddot{b}_1z^{-1}+\ddot{b}_2z^{-2}+\ldots+\ddot{b}_{{N-1}}z^{-({N-1})}} {1+a_1z^{-1}+a_2z^{-2}+\ldots+a_{N}z^{-N}},&\dot{K}=\frac{\ddot{b}_0}{\dot{a}_0} \nonumber \\[10mu] &= \dot{K}\,\color{blue}{z^{-M}}\,\frac{\ddot{b}_0+\ddot{b}_1z^{-1}+\ddot{b}_2z^{-2}+\ldots+\ddot{b}_{{N-1}}z^{-({N-1})}}{(z-r_1)(z-r_2)(z-r_3)\ldots(z-r_N)},&\shaded{\dot{K}=\frac{\ddot{b}_0}{\dot{a}_0}} \end{align} $$

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

$$ \begin{align} G(z)&=\color{blue}{z^{-M}}\left(\frac{c_1}{z-r_1}+\frac{c_2}{z-r_2}+\cdots+\frac{c_N}{z-r_N}\right) \nonumber \\[8mu] &=\color{blue}{z^{-M}}\left(\sum_{k=1}^N\,\frac{c_k}{z-r_k}\right) \end{align} $$

Therefore, $G(z)$ transforms to a parallel combination of delayed scaled step functions in the time-domain

$$ \shaded{ \begin{align} g[n]&=\left(c_1+c_2r_2+\cdots+c_{\small N}(r_{\small N})^{\,\small N-1}\right)\,\color{grey}{\gamma[n\color{black}{-M-1}]} \nonumber \\ &= \sum_{k=1}^{N}c_k(r_k)^{k-1}\, \color{grey}{\gamma[n\color{black}{-M-1}]} \nonumber \end{align} } $$

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}$

$$ \begin{align} F(z) &= f_0+f_1z^{-1}+f_2z^{-2}+\ldots+f_Kz^{-K},\ \ \ K=M-N& \forall_{M\geq N} \end{align} $$

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})$.

$$ \def\lfz#1{\overset{\Large#1}{\,\circ\kern-6mu-\kern-7mu-\kern-7mu-\kern-6mu\bullet\,}} \def\ztransform{\lfz{\mathcal{Z}}} y_{\tiny{IIR}}[n]=\sum_{k=0}^{N}\,c_k\,(r_k)^{n-1}\,\gamma[n-1] \ztransform \sum_{k=0}^{N}c_k\,\frac{1}{z-r_k} $$

where $c_k$ is some finite-order polynomial and $r_k$ is the $k$^th pole of the filter.

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$.

$$ Y(z)=\frac{2z^3+z^2-z+4}{(z-2)^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

$$ \begin{align} Y(z) &= \color{blue}{-\frac{1}{2}}+\color{blue}{\frac{5}{2}}\frac{z}{(z-2)}+\color{blue}{8}\frac{z}{(z-2)^2}+\color{blue}{11}\frac{z}{(z-2)^3} \end{align} $$

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

$$ \begin{align} y[n] &= \color{blue}{-\tfrac{1}{2}}\delta[n]+\color{blue}{\tfrac{5}{2}}2^n\gamma[n]+\tfrac{1}{2}\color{blue}{8}n2^n\gamma[n]+\tfrac{1}{2^2}\color{blue}{11}\tfrac{1}{2}n(n-1)2^n, & n\geq0\nonumber\\[8mu] &= -\tfrac{1}{2}\delta[n]+\tfrac{5}{2}2^n\gamma[n]+4n2^n\gamma[n]+\tfrac{11}{8}n(n-1)2^n, &n\geq0 \end{align} $$

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

$$ H(z) = \frac{\color{teal}{2}+\color{teal}{6}z^{-1}+\color{teal}{6}z^{-2}+\color{teal}{2}z^{-3}}{\color{olive}{1}\color{olive}{-2}z^{-1}+\color{olive}{1}z^{-2}} \label{eq:example2_def} $$

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$)

$$ H(z) = \underbrace{\color{purple}{10}+\color{purple}{2}z^{-1}}_{F(z)}+\underbrace{\frac{\color{green}{-8}+\color{green}{24}z^{-1}}{1-2z^{-1}+z^{-2}}}_{G(z)} $$

Examine the IIR part, $G(z)$

$$ G(z) = \frac{\color{green}{-8}+\color{green}{24}z^{-1}}{1-2z^{-1}+z^{-2}} \label{eq:example2g} $$

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

$$ G(z) = \frac{\color{blue}{-24}}{1-z^{-1}}+\frac{\color{blue}{16}}{\left(1-z^{-1}\right)^{2}} $$

With $H(z) = F(z)+G(z)$, and the impulse response $Y(z) = \Delta(z)\,H(z)=H(z)$

$$ \shaded{ Y(z)=\color{purple}{10}+\color{purple}{2}z^{-1}\color{blue}{-}\frac{\color{blue}{24}}{1-z^{-1}}+\frac{\color{blue}{16}}{\left(1-z^{-1}\right)^{2}} } $$

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

$$ \begin{align} y[n] &= \color{purple}{10}\delta[n]+\color{purple}{2}\gamma[n-1]\color{blue}{-24}\gamma[n]+\color{blue}{16}(n+1)\gamma[n] \\[8mu] &= \{2,10,26,42,\dots\} \end{align} $$

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}$

$$ H(z) = \frac{\color{grey}{2}z^3+\color{grey}{6}z^{2}+\color{grey}{6}z+\color{grey}{2}}{z(\color{purple}{1}z^2\color{purple}{-2}z+\color{purple}{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)

$$ \begin{align} H(z) &= \frac{\color{purple}{2}z+\color{purple}{10}}{z}+\frac{\color{green}{24}z\color{green}{-8}}{z(\color{purple}{1}z^2\color{purple}{-2}z+\color{purple}{1})} \nonumber \\[6mu] &= \underbrace{\color{purple}{2}+\color{purple}{10}z^{-1}}_{\triangleq F(z)}+\underbrace{z^{-1}\frac{\color{green}{24}z\color{green}{-8}}{z^2\color{purple}{-2}z+\color{purple}{1}}}_{\triangleq G(z)} \end{align} $$

Examine the IIR part $G(z)$, by bringing it back to a polynomial in $z^{-1}$

$$ \begin{align} G(z)&=z^{-1}\frac{\color{green}{24}z\color{green}{-8}}{z^2\color{purple}{-2}z+\color{purple}{1}} =z^{-2}\frac{\color{green}{24}\color{green}{-8}z^{-1}}{\color{purple}{1}\color{purple}{-2}z^{-1}+\color{purple}{1}z^{-2}} \end{align} $$

Partial fraction expansion, left as an exercise to the reader

$$ \begin{align} G(z) &= z^{-2}\frac{\color{green}{24}\color{green}{-8}z^{-1}}{\color{purple}{1}\color{purple}{-2}z^{-1}+\color{purple}{1}z^{-2}} = z^{-2}\left(\frac{\color{blue}{8}}{1-z^{-1}}+\frac{\color{blue}{16}}{(1-z^{-1})^{\color{pink}{2}}}\right) \end{align} $$

With $H(z)=F(z)+G(z)$, and the impulse response $Y(z)=\Delta(z)\,H(z)=H(z)$

$$ \shaded{ Y(z)=F(z)+G(z)=\color{purple}{2}+\color{purple}{10}z^{-1}+z^{-2}\left(\frac{\color{blue}{8}}{1-z^{-1}}+\frac{\color{blue}{16}}{(1-z^{-1})^{\color{magenta}{2}}}\right) } $$

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

$$ H(z) = \color{purple}{10}+\color{purple}{2}z^{-1}+\frac{\color{blue}{-24}}{\left(\color{brown}{1}-z^{-1}\right)^{\color{magenta}{1}}}+\frac{\color{blue}{16}}{\left(\color{brown}1-z^{-1}\right)^{\color{magenta}{2}}} $$

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