4.离散时间系统的Matlab实现

– 调用格式：y=filter(b , a, x)，对于由矢量b, a决 定的数字系统（b和a分别表示系统函数H（z） 对应的分子项和分母项系数构成的数组，而且 分母系数要归一化处理。）当输入信号为x时， 对x中的数据进行滤波，结果存于y中，长度取 max(na , nb).
MATLAB数字信号处理
MATLAB数字信号处理
0.3 0.25 0.2 0.15 0.1 0.05 0 -0.05 -0.1 -0.15 -0.2
0
20
40
60
80
100
120
MATLAB数字信号处理
filtic子函数
• 功能：为filter子函数选择初始条件。 • 调用格式：
– zi=filtic(b ,a ,y ,x)；求给定输入x和y时的初始状 态。 – zi=filtic(b , a, y)；求x=0，给定输入y时的初始状 态。 – 其中，x和y分别是表示过去的输入和输出。
离散时间系统的Matlab实现
MATLAB数字信号处理
Impz函数
• 功能：求解系统的单位冲击响应 • 调用方式：
– [h,t]=impz(b,a)：b、a分别为系统传递函数的分 子和分母的系数向量。返回系统(b,a)的冲激响 应h和相应的时间轴向量t。 – [h,t]=impz(b,a,n)：返回n点冲激响应。 – impz(b,a);在当前窗口用stem(t,h )函数出图。
MATLAB数字信号处理
• 【例5】已知一个因果系统的差分方程为 6y(n)+2y(n-2)=x(n)+3x(n-1)+3x(n-2)+x(n-3)满 足初始条件y(-1)=0,x(-1)=0.求系统的单位冲 激响应和阶跃响应。用filter()函数的MATLAB 程序（取N=32点作图）
MATLAB数字信号处理
MATLAB数字信号处理
freqz
• 功能：用于求离散时间系统的频率响应函 数。 • 调用格式： • 1）[h, w]=freqz(b ,a ,n)。可以得到数字滤波 器的n点复频响应值，这n个点均匀地分布 之[0, pi ]上，并将这n个频点的频率记录在w 中，相应的频响值记录在h中，n缺省时取 512点。
MATLAB数字信号处理
MATLAB数字信号处理
• • • • • • • •
subplot(1,2,1),stem(n,hn,’filled’); title(‘系统的单位冲击响应’)； ylabel(‘h(n)’);xlabel(‘n’); axis([0,N-1,-1.1*min(hn),1.1*max(hn)]); subplot(1,2,2), stem(n,gn,’k’); title(‘系统的单位阶跃响应’)； ylabel(‘g(n)’);xlabel(‘n’); axis([0,N-1,-1.1*min(gn),1.1*max(gn)]);
MATLAB数字信号处理
• 【例2】已知一个因果系统的差分方程为 6y(n)+2y(n-2)=x(n)+3x(n-1)+3x(n-2)+x(n-3)满 足初始条件y(-1)=0,x(-1)=0.求系统的单位冲 激响应和阶跃响应。
– 解：将上述方程对y(n)项系数进行归一化，得 到其系统函数分子和分母系数wk.baidu.com
MATLAB数字信号处理
• • • • • •
x2=[n>=0]; %单位阶跃信号 gn=filter(b,a,x2,xi); subplot(1,2,2), stem(n,gn,’k’); title(‘系统的单位阶跃响应’)； ylabel(‘g(n)’);xlabel(‘n’); axis([0,N-1,-1.1*min(gn),1.1*max(gn)]);
MATLAB数字信号处理
• 2）[h f]=freqz(b, a, n, Fs)；用于对H(jw) 在[0, Fs/2]上等间隔采样n点，采样点频率及相应 频响值分别记录在f和h中，由用户指定Fs（ 以Hz为单位）的值。 • 3）h=freqz(b, a, w)；用于对H(jw) 在[0, 2*pi ] 上进行采样，采样频率点由矢量w指定。 • 4）h=freqz(b, a, f, Fs)；用于对H(jw)在[0, Fs] 上采样，采样频率点由矢量f指定。
0.2
0.1
0
-0.1
-0.2
0
5
10
15
20 25 30 n (samples)
35
40
45
MATLAB数字信号处理
dstep
• 功能：求解数字系统的阶跃响应。 • 调用格式：
– [h,t]=dstep(b,a);求解数字系统的阶跃响应h，取 样点数为缺省值。 – [h,t]=dstep(b,a,n);求解数字系统的阶跃响应h， 取样点数为n值。 – dstep(b,a)；在当前窗口用stairs(t,h)函数出图。
• •
a0=1, a1=0, a2=1/3, a3=0 b0=1/6, b1=1/2, b2=1/2, b3=1/6
MATLAB数字信号处理
用impz()函数的MATLAB程序（取 N=32点作图）
• • • • • • a=[1,0,1/3, 0]; b=[1/6, 1/2, 1/2, 1/6]; N=32; n=0:N-1; hn=impz(b,a,n); gn=dstep(b,a,n);
MATLAB数字信号处理
• subplot(2,2, 2), plot(n/pi, angle(h)); grid %作系统的相位 频响图 • axis([0,1,1.1*min(angle(h)),1.1*max(angle(h))]); • ylabel('相位');xlabel('以pi为单位的频率'); • title('相频响应'); • db=20*log10(abs(h)); • subplot(2, 2, 3), plot(n/pi, db); grid %作系统的相对幅度 频响图 • title('幅频响应（db）'); • subplot(2, 2, 4), zplane(b, a); grid %作零极点分布图 • title('零极点分布');
MATLAB数字信号处理
• 【例6】已知离散时间系统的系统函数
0.2 0.1z 0.3z 0.1z 0.2 z H ( z) 1 1.1z 1 1.5 z 2 0.7 z 3 0.3z 4
1
2
3
4
• 求系统在0-pi 频率范围内，归一化的绝对幅 度频率响应，相对幅度频率响应，相位频 率响应和零极点分部图。
MATLAB数字信号处理
幅 频 响 应 （ V） 1 0.8
幅度 相位
相频响应
2 0 -2 0 0.5 幅 频 响 应 （ db） 1 0 0.5 以 pi为 单 位 的 频 率 零极点分布 1
0.6 0.4 0.2
0 -20 -40 -60 -80
Imaginary Part
1 0.5 0 -0.5 -1 0 0.5 1 -1 -0.5 0 0.5 Real Part 1
MATLAB数字信号处理
系统的单位冲击响应 0.5 0.4 0.3 0.2 0 5 10 15 20 n 系统的单位阶跃响应 25 30
h(n) g(n)
1
0.5
0 0 5 10 15 n 20 25 30
MATLAB数字信号处理
zplane
• zplane(z,p)：绘制系统零极点图，“o”表示 零点，“x”表示极点。z，p分别为零点和极 点向量。 • zplane(b,a)：b、a分别为系统传递函数的分 子和分母系数向量。
MATLAB数字信号处理
状态空间的转换
• tf2ss函数：由传递函数到状态空间的转换 • [A,B,C,D] = TF2SS(NUM,DEN) H(s)=NUM(s)/DEN(s)到
– x = Ax + Bu – y = Cx + Du – 举例：>> b=[0.3 0.2]; >> a=[1 -0.4 -0.7]; >> [A,B,C,D]=tf2ss(b,a);
零点 极点
Imaginary Part
0.2 0 -0.2 -0.4 -0.6 -0.8 -1 -1 -0.5 0 Real Part 0.5 1
MATLAB数字信号处理
filter
• 功能：对数字系统的输入信号进行滤波处 理。因为一个离散系统可以看作是一个滤 波器，系统的输出就是输入经过滤波器滤 波的结果。
MATLAB数字信号处理
• b=[0.2, 0.1, 0.3, 0.1, 0.2]; • a=[1, -1.1, 1.5, -0.7, 0.3]; • n=(0: 500)*pi/500; %在[0,pi]的范围内取 501个采样点 • [h, w]=freqz(b, a, n); %求系统的频率响应 • subplot(2, 2, 1), plot(n/pi, abs(h)); grid %作 系统的绝对幅度频响图 • axis([0,1,1.1*min(abs(h)), 1.1*max(abs(h))]); • ylabel('幅度'); • title('幅频响应（V）');
• • • • • • • • • • •
a=[1,0,1/3, 0]; b=[1/6, 1/2, 1/2, 1/6]; xi=filtic(b,a,0,0); N=32; n=0:N-1; x1=[n==0]; %单位冲激信号 hn=filter(b,a,x1,xi); subplot(1,2,1),stem(n,hn,’filled’); title(‘系统的单位冲激响应’)； ylabel(‘h(n)’);xlabel(‘n’); axis([0,N-1,-1.1*min(hn),1.1*max(hn)]);
• [y , zf]=filter(b , a, x);除得到结果矢量y外， 还得到x 的最终状态矢量zf。 • y=filter(b , a, x, zi)；可在zi中指定x的初始 状态。
MATLAB数字信号处理
• • • • •
【例4】计算低通滤波器的冲激响应。 x=[1 zeros(1,100)]; [b,a] = cheby1(11,1,.4); y=filter(b,a,x); stem(y);
MATLAB数字信号处理
• • • • •
【例3】计算线性系统(b,a)的零点和极点。 b=[0.2 0.1 0.3 0.1 0.15]; a=[1 -1.1 1.45 -0.6 0.3]; zplane(b,a); legend('零点','极点');
MATLAB数字信号处理
1 0.8 0.6 0.4
MATLAB数字信号处理
• • • •
【例1】计算线性系统(b,a)的冲激响应： b=[0.2 0.1 0.3 0.1 0.15]; a=[1 -1.1 1.45 -0.6 0.3]; impz(b,a,50);
MATLAB数字信号处理
Impulse Response 0.5
0.4
0.3
Amplitude