一、离散时间系统的响应

离散时间LTI系统可用线性常系数差分方程来描述，即：

$\sum^{N}_{i=0}{a_{i}y(n-1)} = \sum_{j=0}^{M}{b_{j}x(n-j)}$

其中，$a_i(i=0,1,2,...,N)$和$b_i(i=0,1,2,...,M)$为实常数。

MATLAB中函数filter可对差分方程在指定时间范围内的输入序列所产生的响应进行求解。函数filter的语句格式为：y=filter(b,a,x)，其中，x为输入的离散序列：y为输出的离散序列；y的长度与x的长度一样；b与a分别为差分方程左端与右端的系数向量。下面的示例，使用MATLAB命令绘制出当激励信号为$x(n)=(1/2)^nu(n)$时，系统$3y(n)-4y(n-1)+2y(n-2)=x(n)+2x(n-1)$的零状态响应。

1
2
3
4
5
6
7
8
9
10
11
12
13
14

% 离散系统系统零状态响应

a = [3 -4 2];
b = [1 2];

n = 0:30;
x = (1/2) .^n;
y = filter(b, a, x);

% 绘图
stem(n, y, "fill");
grid on;
title("系统响应y(n)");
xlabel("n");

二、离散时间系统的单位取样响应

系统的单位取样响应定义为系统在$\delta{n}$激励下系统的零状态响应，用$h(n)$表示。MATLAB求解单位取样响应可利用函数filter，并将激励设置为单位抽样序列。下面的示例演示了如何求解上述系统的单位取样响应：

1
2
3
4
5
6
7
8
9
10
11

% 离散系统的单位取样响应
a = [3 -4 2];
b = [1 2];

n = 0:30;
x = (n == 0);
h = filter(b, a, x);
stem(n, h, "fill");
grid on;
xlabel("n");
title("系统单位取样响应h(n)");

MATLAB另一种求单位取样响应的方法是利用控制系统工具箱提供的函数impz来实现。impz函数的常用语句格式为：impz(b,a,N)，其中，参数N通常为正整数，代表计算单位采样响应的样值个数。下面的示例演示了如何使用impz函数绘制出系统$3y(n)-4y(n-1)+2y(n-2)=x(n)+2x(n-1)$的单位取样响应。

1
2
3
4
5
6
7

% impz函数求解离散系统单位采样响应

a = [3 -4 2];
b = [1 2];
impz(b, a,30);
grid on;
title("系统单位取样响应h(n)")

三、离散时间信号的卷积和运算

由于系统的零状态响应是激励与系统的单位取样响应的卷积，因此卷积运算在离散时间信号处理领域被广泛应用。离散时间信号的卷积定义为 ：

$y(n)=x(n)^{*}h(n)=\sum_{m=-\infty}^{\infty}{x(m)h(n-m)}$

可见，离散时间信号的卷积运算是求和运算，因而常被成为卷积和。

MATLAB求离散时间信号卷积和的命令为conv，其语句格式为：$y=conv(x,h)$，其中，x和h表示离散时间信号值的向量，y为卷积结果。用MATLAB进行卷积和运算时，无法实现无限的累加，只能计算时有限信号的卷积。下面的示例演示了如何利用MATLAB的conv命令求两个长为4的矩形序列的卷积和，即$g(n)=[u(n)-u(n-4)]^{*}[u(n)-u(n-4)]$，其结果应该是长为7（4+4-1=7）的三角序列。用向量$[1 1 1 1]$表示矩形序列：

1
2
3
4
5
6
7
8

% 离散时间信号的卷积
x1 = [1 1 1 1];
x2 = [1 1 1 1];
g = conv(x1, x2);

% 绘制
n = 0:6;
stem(n, g, "fill");

对于给定函数的卷积和，我们应计算卷积结果的起始点及其长度。两个时限序列的卷积和长度等于两个序列长度的和减 1。

例如：已知某系统的单位取样响应为$h(t)=0.8^{n}[u(n)-u(n-8)]$，试用MATLAB求当激励信号为$x(n)=u(n)-u(n-4)$时，系统的零状态响应。

解：MATLAB 中可通过卷积求解零状态响应,即$x(n)*h(n)$。由题意可知，描述$h(n)$向量的长度至少为8，描述x(n)向量的长度至少为4，因此为了美观，将$h(n)$向量和$x(n)$向量加上一些附加的零值。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31

nx = -1:5; % x(n)向量显示范围（添加了附加的零值）
nh = -2:10; % h(n)向量显示范围（添加了附加的零值）
x = uDT(nx) - uDT(nx - 4); 
h = 0.8 .^nh.*(uDT(nh)-uDT(nh-8));
y = conv(x, h);

ny1 = nx(1) + nh(1); % 卷积结果起始点
ny2 = nx(end) + nh(end); % 卷积结果重点

ny = ny1:ny2;

subplot(3, 1, 1)
stem(nx, x, "filled");
grid on;
xlabel("n");
title("x(n)");
axis([-4 16 0 3]);

subplot(3, 1, 2);
stem(nh, h, "filled");
grid on;
xlabel("n");
title("h(n)");
axis([-4 16 0 3]);

subplot(3, 1, 3);
stem(ny, y, "filled");
grid on;
xlabel("n");
title("x(n)*h(n)");
axis([-4 16 0 3]);