Chapter 16: Discrete-Time and Sampled-Data System Dynamics

Lesson 3: z-Transform and Discrete-Time Transfer Functions

This lesson develops the z-transform as the primary analysis tool for linear discrete-time systems. We derive key properties, define the region of convergence (ROC), obtain discrete-time transfer functions from difference equations, and connect poles/zeros to stability and frequency response. The lesson includes implementations in Python, C++, Java, MATLAB/Simulink, and Wolfram Mathematica.

1. Conceptual Overview

In continuous-time analysis, the Laplace transform converts differentiation into algebraic multiplication by \( s \). In discrete-time analysis, the z-transform converts time shifts into multiplication by powers of \( z^{-1} \). This makes linear constant-coefficient difference equations (LCCDEs) algebraically solvable in the transform domain.

The object \( H(z) \) plays the same structural role as \( G(s) \): it is the transfer function under zero initial conditions and exposes poles, zeros, and stability properties.

flowchart TD
  A["Difference equation model"] --> B["Apply z-transform"]
  B --> C["Use shift property for delays"]
  C --> D["Solve for H(z)=Y(z)/X(z)"]
  D --> E["Find poles, zeros, ROC"]
  E --> F["Check stability"]
  F --> G["Compute H(e^(j w)) and time response"]
        

2. Definitions and Region of Convergence

Bilaterial z-transform:

\[ X(z) = \sum_{n=-\infty}^{\infty} x[n] z^{-n}. \]

Unilateral z-transform (causal/IVP-oriented):

\[ X^{+}(z) = \sum_{n=0}^{\infty} x[n] z^{-n}. \]

The region of convergence (ROC) is the set of \( z \) for which the series converges absolutely. Suppose \( |x[n]| \le C \alpha^n \) for \( n \ge 0 \). Then

\[ \sum_{n=0}^{\infty} |x[n] z^{-n}| \le \sum_{n=0}^{\infty} C \left(\frac{\alpha}{|z|}\right)^n, \]

so the unilateral transform converges whenever \( |z| > \alpha \) (comparison with a geometric series).

Example: for \( x[n]=a^n u[n] \),

\[ X(z)=\sum_{n=0}^{\infty} (a z^{-1})^n = \frac{1}{1-a z^{-1}}, \qquad \text{ROC: } |z| > |a|. \]

Multiplying numerator and denominator by \( z \) gives \( X(z)=\frac{z}{z-a} \), making the pole explicit.

3. Core Properties and Short Proofs

Linearity follows directly from sum distributivity:

\[ \mathcal{Z}\{\alpha x[n]+\beta y[n]\}=\alpha X(z)+\beta Y(z). \]

Shift property (bilateral form):

\[ \mathcal{Z}\{x[n-k]\}=z^{-k}X(z). \]

Proof: with \( m=n-k \),

\[ \sum_{n=-\infty}^{\infty}x[n-k]z^{-n} =\sum_{m=-\infty}^{\infty}x[m]z^{-(m+k)} =z^{-k}\sum_{m=-\infty}^{\infty}x[m]z^{-m}=z^{-k}X(z). \]

Convolution property:

\[ (x*h)[n]=\sum_{m=-\infty}^{\infty}x[m]h[n-m] \quad \Longrightarrow \quad \mathcal{Z}\{x*h\}=X(z)H(z). \]

The proof is obtained by substituting the convolution sum and interchanging the order of summation.

Index multiplication property:

\[ \mathcal{Z}\{n\,x[n]\}=-z\frac{dX(z)}{dz}. \]

Initial value theorem (causal):

\[ x[0]=\lim_{z \to \infty} X^{+}(z). \]

Final value theorem (causal and pole-condition satisfied):

\[ \lim_{n \to \infty}x[n]=\lim_{z \to 1}(1-z^{-1})X^{+}(z). \]

Always verify pole locations before using the final value theorem.

4. Discrete-Time Transfer Function from Difference Equations

Consider the LCCDE

\[ y[n] + a_1 y[n-1] + \cdots + a_N y[n-N] = b_0 x[n] + b_1 x[n-1] + \cdots + b_M x[n-M]. \]

Under zero initial conditions, applying the z-transform gives

\[ Y(z)\left(1+a_1 z^{-1}+\cdots+a_N z^{-N}\right) = X(z)\left(b_0+b_1 z^{-1}+\cdots+b_M z^{-M}\right). \]

Hence the transfer function is

\[ H(z)\equiv \frac{Y(z)}{X(z)} = \frac{b_0+b_1 z^{-1}+\cdots+b_M z^{-M}} {1+a_1 z^{-1}+\cdots+a_N z^{-N}}. \]

Multiplying numerator and denominator by \( z^N \) yields a rational function in positive powers of \( z \), which is convenient for explicit pole-zero computation.

flowchart TD
  A["LCCDE coefficients a_i, b_i"] --> B["Transform delays to z^(-k)"]
  B --> C["Collect Y(z) and X(z) terms"]
  C --> D["Form H(z)"]
  D --> E["Factor denominator and numerator"]
  E --> F["Poles / Zeros / Stability"]
        

5. Poles, Zeros, ROC, Stability, and Frequency Response

For a causal rational system,

\[ H(z)=K\frac{\prod_{i=1}^{M}(1-z_i z^{-1})}{\prod_{j=1}^{N}(1-p_j z^{-1})}. \]

ROC for a causal system:

\[ \text{ROC}: \quad |z| > \max_j |p_j|. \]

BIBO stability criterion (causal rational LTI):

\[ |p_j| < 1 \quad \forall j. \]

Equivalently, the unit circle must lie inside the ROC so that the impulse response is absolutely summable.

Frequency response on the unit circle:

\[ H(e^{j\omega}) = H(z)\big|_{z=e^{j\omega}}, \quad \omega \in [-\pi,\pi]. \]

This is the discrete-time analog of sinusoidal steady-state response in continuous-time frequency analysis.

Continuous-to-discrete pole mapping insight:

\[ z=e^{sT_s}. \]

If \( \Re(s) < 0 \), then \( |z|=e^{\Re(s)T_s} < 1 \); hence stable continuous-time poles map inside the unit circle.

6. Worked Example with Pole-Zero Interpretation

Consider

\[ y[n]-1.5y[n-1]+0.56y[n-2]=0.2x[n]+0.1x[n-1]. \]

The transfer function is

\[ H(z)=\frac{0.2+0.1z^{-1}}{1-1.5z^{-1}+0.56z^{-2}} = \frac{0.2z^2+0.1z}{z^2-1.5z+0.56}. \]

Denominator factorization:

\[ z^2-1.5z+0.56=(z-0.8)(z-0.7). \]

Therefore the poles are \( 0.8 \) and \( 0.7 \), both inside the unit circle, so the system is BIBO stable. The numerator yields a zero at \( z=-0.5 \) (and an origin zero appears in the positive-power representation due to degree alignment).

Inverse transform by partial fractions:

\[ H(z)=\frac{A}{1-0.8z^{-1}}+\frac{B}{1-0.7z^{-1}} \;\Longrightarrow\; h[n]=A(0.8)^n u[n]+B(0.7)^n u[n]. \]

7. Python Implementation

File: Chapter16_Lesson3.py


# Chapter16_Lesson3.py
import numpy as np

def impulse_response(b, a, N=20):
    x = np.zeros(N); x[0] = 1.0
    y = np.zeros(N)
    for n in range(N):
        y[n] += sum(b[k] * x[n-k] for k in range(len(b)) if n-k >= 0)
        y[n] -= sum(a[k] * y[n-k] for k in range(1, len(a)) if n-k >= 0)
    return y / a[0]

def H_ejw(b, a, w):
    zinv = np.exp(-1j * w)
    num = sum(b[k] * zinv**k for k in range(len(b)))
    den = sum(a[k] * zinv**k for k in range(len(a)))
    return num / den

if __name__ == "__main__":
    # H(z) = (0.2 + 0.1 z^-1) / (1 - 1.5 z^-1 + 0.56 z^-2)
    b = np.array([0.2, 0.1])
    a = np.array([1.0, -1.5, 0.56])

    h = impulse_response(b, a, 20)
    print("h[n] (first 10 samples):", np.round(h[:10], 6))

    ws = np.linspace(0, np.pi, 5)
    for w in ws:
        H = H_ejw(b, a, w)
        print(f"w={w:.3f}, |H|={abs(H):.6f}, phase={np.angle(H):.6f}")

    poles = np.roots(a)
    zeros = np.roots(np.array([0.2, 0.1]))
    print("poles:", poles, "stable:", np.all(np.abs(poles) < 1))
    print("zeros:", zeros)
      

8. C++ Implementation

File: Chapter16_Lesson3.cpp


// Chapter16_Lesson3.cpp
#include <iostream>
#include <vector>
#include <complex>
#include <cmath>

using cd = std::complex<double>;

std::vector<double> impulseResponse(const std::vector<double>& b, const std::vector<double>& a, int N){
    std::vector<double> x(N,0.0), y(N,0.0); x[0]=1.0;
    for(int n=0;n<N;++n){
        for(int k=0;k<(int)b.size();++k) if(n-k>=0) y[n] += b[k]*x[n-k];
        for(int k=1;k<(int)a.size();++k) if(n-k>=0) y[n] -= a[k]*y[n-k];
        y[n] /= a[0];
    }
    return y;
}

cd Hejw(const std::vector<double>& b, const std::vector<double>& a, double w){
    cd zinv = std::exp(cd(0.0,-w)), num(0,0), den(0,0);
    for(int k=0;k<(int)b.size();++k) num += b[k]*std::pow(zinv,k);
    for(int k=0;k<(int)a.size();++k) den += a[k]*std::pow(zinv,k);
    return num/den;
}

int main(){
    std::vector<double> b{0.2,0.1}, a{1.0,-1.5,0.56};
    auto h = impulseResponse(b,a,20);
    for(int i=0;i<10;++i) std::cout << "h[" << i << "]=" << h[i] << "\n";
    for(int i=0;i<=4;++i){
        double w = M_PI*i/4.0;
        cd H = Hejw(b,a,w);
        std::cout << "w=" << w << " |H|=" << std::abs(H) << " phase=" << std::arg(H) << "\n";
    }
    return 0;
}
      

9. Java Implementation

File: Chapter16_Lesson3.java


// Chapter16_Lesson3.java
public class Chapter16_Lesson3 {
    static class C {
        double re, im;
        C(double re, double im){ this.re = re; this.im = im; }
        C add(C o){ return new C(re + o.re, im + o.im); }
        C mul(C o){ return new C(re*o.re - im*o.im, re*o.im + im*o.re); }
        C mul(double s){ return new C(re*s, im*s); }
        C div(C o){
            double d = o.re*o.re + o.im*o.im;
            return new C((re*o.re + im*o.im)/d, (im*o.re - re*o.im)/d);
        }
        double abs(){ return Math.hypot(re, im); }
        double ang(){ return Math.atan2(im, re); }
    }

    static double[] impulse(double[] b, double[] a, int N){
        double[] x = new double[N], y = new double[N];
        x[0] = 1.0;
        for(int n=0; n<N; n++){
            for(int k=0; k<b.length; k++) if(n-k >= 0) y[n] += b[k]*x[n-k];
            for(int k=1; k<a.length; k++) if(n-k >= 0) y[n] -= a[k]*y[n-k];
            y[n] /= a[0];
        }
        return y;
    }

    static C Hejw(double[] b, double[] a, double w){
        C zinv = new C(Math.cos(w), -Math.sin(w));
        C num = new C(0,0), den = new C(0,0), p = new C(1,0);
        for (double bi : b){ num = num.add(p.mul(bi)); p = p.mul(zinv); }
        p = new C(1,0);
        for (double ai : a){ den = den.add(p.mul(ai)); p = p.mul(zinv); }
        return num.div(den);
    }

    public static void main(String[] args){
        double[] b = {0.2, 0.1}, a = {1.0, -1.5, 0.56};
        double[] h = impulse(b, a, 20);
        for(int i=0; i<10; i++) System.out.printf("h[%d]=%.6f%n", i, h[i]);
        for(int i=0; i<=4; i++){
            double w = Math.PI * i / 4.0;
            C H = Hejw(b, a, w);
            System.out.printf("w=%.3f |H|=%.6f phase=%.6f%n", w, H.abs(), H.ang());
        }
    }
}
      

10. MATLAB / Simulink Implementation

File: Chapter16_Lesson3.m


% Chapter16_Lesson3.m
% y[n] - 1.5 y[n-1] + 0.56 y[n-2] = 0.2 x[n] + 0.1 x[n-1]
b = [0.2 0.1];
a = [1 -1.5 0.56];

N = 20;
x = [1 zeros(1,N-1)];
h = filter(b, a, x);
disp('h[n] first 10 samples:'); disp(h(1:10).');

[H, w] = freqz(b, a, 5);
disp(table(w, abs(H), angle(H), 'VariableNames', {'w','Mag','Phase'}));

Ts = 1;
Gz = tf(b, a, Ts, 'Variable', 'z^-1');
disp(Gz);

p = pole(Gz); z = zero(Gz);
disp('Poles:'), disp(p)
disp('Zeros:'), disp(z)
disp(['Stable = ', num2str(all(abs(p) < 1))])
      

Simulink note: Use a Discrete Transfer Fcn block with numerator \( [0.2\ \ 0.1] \), denominator \( [1\ \ -1.5\ \ 0.56] \), and sample time \( T_s \). This directly implements \( H(z) \).

11. Wolfram Mathematica Implementation

File: Chapter16_Lesson3.nb


Notebook[{
  Cell["Chapter16_Lesson3.nb - z-Transform and Discrete-Time Transfer Functions", "Title"],
  Cell[BoxData@ToBoxes[
    eq = y[n] - 1.5 y[n - 1] + 0.56 y[n - 2] == 0.2 KroneckerDelta[n] + 0.1 KroneckerDelta[n - 1];
  ], "Input"],
  Cell[BoxData@ToBoxes[
    sol = RSolveValue[{eq, y[-1] == 0, y[-2] == 0}, y[n], n];
  ], "Input"],
  Cell[BoxData@ToBoxes[
    h = Table[N[sol /. n -> k], {k, 0, 19}]
  ], "Input"],
  Cell[BoxData@ToBoxes[
    Hz = FullSimplify[ZTransform[sol, n, z]]
  ], "Input"],
  Cell[BoxData@ToBoxes[
    w = Table[k Pi/4, {k, 0, 4}];
    vals = N[(Hz /. z -> Exp[I #]) & /@ w];
    Transpose[{w, Abs[vals], Arg[vals]}]
  ], "Input"]
}]
      

12. Problems and Solutions

Problem 1 (basic transform and ROC): Find the z-transform and ROC of \( x[n]=3(0.4)^n u[n] \).

Solution:

\[ X(z)=\sum_{n=0}^{\infty}3(0.4)^n z^{-n} =3\sum_{n=0}^{\infty}(0.4 z^{-1})^n =\frac{3}{1-0.4 z^{-1}}, \quad |z| > 0.4. \]

In positive-power form: \( X(z)=\frac{3z}{z-0.4} \).

Problem 2 (transfer function derivation): For \( y[n]+0.3y[n-1]-0.1y[n-2]=x[n]-2x[n-1] \), derive \( H(z) \).

Solution:

\[ Y(z)\left(1+0.3z^{-1}-0.1z^{-2}\right)=X(z)\left(1-2z^{-1}\right). \]

\[ H(z)=\frac{1-2z^{-1}}{1+0.3z^{-1}-0.1z^{-2}} =\frac{z^2-2z}{z^2+0.3z-0.1}. \]

Problem 3 (final value theorem): Let \( H(z)=\frac{1}{1-0.5z^{-1}} \) and input \( x[n]=u[n] \). Find the steady-state output.

Solution: Since \( X(z)=\frac{1}{1-z^{-1}} \),

\[ Y(z)=\frac{1}{(1-z^{-1})(1-0.5z^{-1})}. \]

\[ \lim_{n \to \infty}y[n] =\lim_{z \to 1}(1-z^{-1})Y(z) =\frac{1}{1-0.5}=2. \]

Problem 4 (inverse z-transform): Invert \( X(z)=\frac{1}{(1-0.8z^{-1})(1-0.2z^{-1})} \) for a causal sequence.

Solution:

\[ X(z)=\frac{A}{1-0.8z^{-1}}+\frac{B}{1-0.2z^{-1}}. \]

\[ 1=(A+B)+(-0.2A-0.8B)z^{-1}. \]

Hence \( A+B=1 \) and \( 0.2A+0.8B=0 \), so \( A=\frac{4}{3} \), \( B=-\frac{1}{3} \). Therefore

\[ x[n]=\frac{4}{3}(0.8)^n u[n]-\frac{1}{3}(0.2)^n u[n]. \]

Problem 5 (s-plane to z-plane pole mapping): A continuous-time pole is \( s=-3+4j \) and \( T_s=0.1 \). Find the corresponding \( z \)-plane pole and decide if it lies inside the unit circle.

Solution:

\[ z=e^{sT_s}=e^{(-3+4j)0.1}=e^{-0.3}e^{j0.4}. \]

Its magnitude is \( |z|=e^{-0.3}\approx 0.7408 \), hence it lies inside the unit circle.

13. Summary

We introduced unilateral and bilateral z-transforms, proved the key properties used in system analysis, derived discrete-time transfer functions from LCCDEs, and connected ROC, poles, zeros, and unit-circle evaluation to stability and frequency response. These tools directly support the next lesson on exact discretization and continuous-discrete conversion methods.

14. References

  1. Ragazzini, J.R., & Zadeh, L.A. (1952). The analysis of sampled-data systems. Transactions of the American Institute of Electrical Engineers, Part II: Applications and Industry, 71(1), 225–234.
  2. Kalman, R.E., & Bertram, J.E. (1959). A unified approach to the theory of sampling systems. Journal of the Franklin Institute, 267(5), 405–436.
  3. Boxer, R. (1957). A note on numerical transform calculus. Proceedings of the IRE, 45, 1401–1406.
  4. Lago, G.V. (1955). Additions to z-transform theory for sampled-data systems. AIEE Transactions, Part II, 73, 403–408.
  5. Linvill, W.K. (1951). Sampled-data control systems studied through comparison of sampling with amplitude modulation. AIEE Transactions, Part II, 70, 1779–1788.
  6. Thompson, W.E. (1956). On two-sided z-transforms. IRE Transactions on Circuit Theory, 3, 156.
  7. Lindorff, D.P. (1959). Application of pole-zero concepts to design of sampled-data systems. IRE Transactions on Automatic Control, AC-4(31), 173–184.
  8. Blachuta, M.J. (1999). On zeros of pulse transfer functions. IEEE Transactions on Automatic Control, 44(6), 1229–1234.
  9. Bistritz, Y. (1986). Zero location of polynomials with respect to the unit circle by the Schur-Cohn test. Systems & Control Letters, 7(4), 267–272.
  10. Razavi, B. (2020). The z-transform for analog designers [The Analog Mind]. IEEE Solid-State Circuits Magazine, 12(3), 8–14.