数字信号处理第三版上机实验答案

实验一：系统响应及系统稳定性

1.实验目的

（1）掌握求系统响应的方法。

（2）掌握时域离散系统的时域特性。

（3）分析、观察及检验系统的稳定性。

2.实验原理与方法

在时域中，描写系统特性的方法是差分方程和单位脉冲响应，在频域可以用系统函数

描述系统特性。已知输入信号可以由差分方程、单位脉冲响应或系统函数求出系统对于该输

入信号的响应，本实验仅在时域求解在计算机上适合用递推法求差分方程的解，最简单的

方法是采用MATLAB语言的工具箱函数filter函数。也可以用MATLAB语言的工具箱函数

conv函数计算输入信号和系统的单位脉冲响应的线性卷积，求出系统的响应。

系统的时域特性指的是系统的线性时不变性质、因果性和稳定性。重点分析实验系统

的稳定性，包括观察系统的暂态响应和稳定响应。

系统的稳定性是指对任意有界的输入信号，系统都能得到有界的系统响应。或者系统

的单位脉冲响应满足绝对可和的条件。系统的稳定性由其差分方程的系数决定。

实际中检查系统是否稳定，不可能检查系统对所有有界的输入信号，输出是否都是有

界输出，或者检行系统的单位脉冲响应满足绝对可和的系件。可行的方法是在系统的输入端

加入单位阶跃序列，如果系统的输出趋近一个常数（包括零），就可以断定系统是稳定的।⑼。

系统的稳态输出是指当〃-8时，系统的输出。如果系统稳定，信号加入系统后，系统输

出的开始一段称为暂态效应，随n的加大，幅度趋于稳定，达到稳态输出。

注意在以下实验中均假设系统的初始状态为零。

3.实验内容及步骤

（1）编制程序，包括产生输入信号、单位脉冲响应序列的子程序，用filter函数或conv

函数求解系统输出响应的主程序。程序中要有绘制信号波形的功能。

（2）给定一个低通滤波器的差分方程为

），（〃）=005x（〃）+005x（〃-1）+0.9），（〃一1）

输入信号内（〃）=《（〃）

x2（n）=u（n）

a）分别求出系统对七（〃）=凡（〃）和的响应序列，并画出其波形。

b）求出系统的单位冲响应，画出其波形。

（3）给定系统的单位脉冲响应为

/?）（/?）=7?10（77）

/?2（/?）=5（〃）+2.55（〃一1）+2.55（〃-2）+8（n-3）

用线性卷积法分别求系统hi（n）和h2（n）对玉⑺=R式几）的输出响应,并画出波形。

（4）给定一谐振器的差分方程为

y(n)—1.8237y(n—1)—0.9801y(n—2)+box(n)—box(n—2)

令％=1/100.49,谐振器的谐振频率为。

a)用实验方法检查系统是否稳定。输入信号为〃(〃)时，画出系统输出波形。

b)给定输入信号为

x(n)=sin(0.014〃)+sin(0.4/?)

求出系统的输出响应，并画出其波形。

4.思考题

(1)如果输入信号为无限长序列，系统的单位脉冲响应是有限长序列，可否用线性卷积

法求系统的响应？如何求？

(2)如果信号经过低通滤波器，把信号的高频分量滤掉，时域信号会有何变化，用前面

第一个实验结果进行分析说明。

5.实验报告要求

(I)简述在时域求系统响应的方法。

(2)简述通过实验判断系统稳定性的方法。分析上面第三个实验的稳定输出的波形。

(3)对各实验所得结果进行简单分析和解释。

(4)简要回答思考题。

(5)打印程序清单和要求的各信号波形。

实验程序清单

%实验1：系统响应及系统稳定性

xln=[l1111111];

hln=[ones(1,10)zeros(1,10)];

h2n=[11zeros(1,10)];

y21n=conv(hlnzxln);

y22n=conv(h2n,xln);

figure(2)

subplot(2,2,1);y='hl(n)*;stem(hln);

title(•(d)Ipt3p¥i»Ad3^Ii6；hl£-n£©');

boxon

subplot(2,2,2);y=1y21(n),;stem(y21n);

title('(e)hln6eR8(n)pA^i»yy21(n),);

boxon

subplot(2,2,3);y=*y2(n),;stem(h2n);

titled(f)Ipt3p¥I»Ad3^IiO：h2(n)');

boxon

subplot(2,2,4);y=,y22(n)1;stem(y22n);

title('(g)h2(n)6eR8(n)pA^i»yy22(n)');

boxon

实验程序运行结果及分析讨论

程序运行结果如图所示。

实验内容(2)系统的单位冲响应、系统对内(〃)=%(〃)和/2(〃)=〃(〃)的响应序列

分别如图(a)、(b)和(c)所示;

实验内容(3)系统h(n)和h2(n)对$(〃)=R式〃)的输出响应分别如图(e)和(g)所示;

实验内容(4)系统对容〃)和l(力)=sin(0.014〃)+sin(0.4〃)的响应序列分别如图(h)

和⑴所示。由图(h)可见，系统对〃(〃)的响应逐渐衰减到零，所以系统稳定。由图(i)可见，

系统对x(〃)=$山(0.0144)十$加(0.4〃)的稳态响应近似为正弦序列5亩(0.4"),这一结论

验证了该系统的谐振频率是rado

简答思考题

(1)如果输入信号为无限长序列，系统的单位脉冲响应是有限长序列，可否用线性卷积

法求系统的响应。①对输入信号序列分段；②求单位脉冲响应h(n)与各段的卷积；③将各段

卷积结果相加。具体实现方法有第三章介绍的重叠相加法和重叠保留法。

(2)如果信号经过低通滤波器，把信号的高频分量滤掉，时域信号的剧烈变化将被平滑,

由实验内容(1)结果图、(b)和(c)可见,经过系统低通滤波使输入信号6(〃)、$5)=凡(〃)

和%(〃)=履〃)的阶跃变化变得缓慢上升与下降。

实验二时域采样与频域采样

实验指导

1.实验目的

时域采样理论与频域采样理论是数字信号处理中的重要理论。要求掌握模拟信号采样

前后频谱的变化，以及如何选择采样频率才能使采样后的信号不丢失信息；要求掌握撅率域

采样会引起时域周期化的概念，以及频率域采样定理及其对频域采样点数选择的指导作用。

2.实验原理与方法

时域采样定理的要点是：

a)对模拟信号与(7)以间隔T进行时域等间隔理想采样，形成的采样信号的频谱

父(/Q)是原模拟信号频谱X//。)以采样角频率(C,=21/7)为周

期进行周期延拓。公式为：

18

尤(4)=叩W⑺]=ZX”(4-OR)

b)采样频率Q,必须大于等于模拟信号最高频率的两倍以上，才能使采样信号的

频谱不产生频谱混叠。

利用计算机计算上式并不方便，下面我们导出另外一个公式，以便用计算机上进行实验。

理想采样信号白⑺和模拟信号xa(t)之间的关系为：

对上式进行傅立叶变换，得到：

(0)=「凡(/)X^t-nT)]e-^dt

=♦匚演

在上式的积分号内只有当，=〃7时，才有非零值，因此：

£(4)=

上式中，在数值上％(〃7)=x(〃)，再将g=07代入，得到：

用(4)二之](廿5

n="x

上式的右边就是序列的傅立叶变换X("’”)，即

先g)=X©&)|

<y=QT

上式说明理想采样信号的傅立叶变换可用相应的采样序列的傅立叶变换得到，只要将自变量

3用Q7代替即可。

频域采样定理的要点是：

a)对信号x(n)的频谱函数X(eT)在[0,2可上等间隔采样N点，得到

XJk)=Xd)攵=(),1,2,・,，N—1

则N点IDFTLX'(A)]得到的序列就是原序列x(n)以N为周期进行周期延拓后的主值区

序列，公式为：

/(〃)=IDFT[XN(Q]N=lXx(〃+W)]R,V(〃)

b)由上式可知，频域采样点数N必须大于等于时域离散信号的长度M(即N2M),才

能使时域不产生混叠，则N点1DFT[XN(Q]得到的序列/(〃)就是原序列x(n),即

XiV(«)=x(n)o如果N>M,4(〃)比原序列尾部多N-M个零点；如果N<M,z则

xiV(/?)=IDFT[Xt(Z)]发生了时域混叠失真，而且/(〃)的长度N也比x(n)的长度

M短，因此。XN(〃)与x(n)不相同。

在数字信号处理的应用中，只要涉及时域或者频域采样，都必须服从这两个采样理论的

要点。

对比上面叙述的时域采样原理和频域采样原理，得到一个有用的结论，这两个采样理论

具有对偶性：“时域采样频谱周期延拓，频域采样时域信号周期延拓”。因此放在一起进行实

验。

3.实验内容及步骤

(1)时域采样理论的验证。

af

给定模拟信号，xa(t)=Ae~sin(Qor)w(r)

式中4=,<7=5072兀，。0=50后nrad/s,它的幅频特性曲线如图

图的幅频特性曲线

现用DFT(FFT)求该模拟信号的幅频特性，以验证时域采样理论。

安照与”)的幅频特性曲线，选取三种采样频率，即6=1"0，300/0，200也。观

测时间选丁,二50〃乂o

为使用DFT,首先用下面公式产生时域离散信号，对三种采样频率.，采样序列按顺

序用力](〃),工2(力)，尢3(〃)表示。

a,d

x(〃)=xa(nT)=Ae~sin(Q0/?T)M(/?T)

因为采样频率不同，得到的25),%(〃)，刍(〃)的长度不同，长度(点数)用

公式N=，〃xF,》计算。选FFT的变换点数为M=64,序列长度不够64的尾部加零。

X(A)=rTT[M〃)]，Q0,1,2,3,---,M-1

27r

式中A代表的频率为a)k=—k.

M

要求：编写实验程序，计算M(〃)、为(〃)和七(〃)的幅度特性，并绘图显示。观察分

析频谱混叠失真。

(2)频域采样理论的验证。

给定信号如下：

n+10</?<13

x(n)=*27-n14<w<26

、0其它

编写程序分别对频谱函数X(/'”)=FT[M〃)J在区间1。,24]上等间隔采样32

和16点，得到X"幻和七6(公：

X^k)=X(ejO1)2万，攵=0/，2,…31

由五%

Xm(Q=X(")2”，攵=0,1215

-k

16

再分别对招2(幻和伏)进行32点和16点IFFT,得到与2(〃)和司6(〃)：

^2(n)=IFFT[X32(^)]32,/2=0,1,2,,31

X16(W)-IFFT[X16(ZC)]16,n—0,1,2,--,15

分别画出X("")、为32(6和X[6(Z)的幅度谱,并绘图显示Mn)、与2(〃)和七6(〃)的波形,

进行对比和分析，验证总结频域采样理论。

提示：频域采样用以下方法容易变程序实现。

①直接调用MATLAB函数fft计算X,2(Q=FFT[X(//)]32就得至IJX(/)在[0,24]的32

点频率域采样

②抽取X,2(外的偶数点即可得到X(*")在[0,2〃]的16点频率域采样％6伏)，即

Xl6(k)=X32(2k),4=0,1,2,,15。

③当然也可以按照频域采样理论，先将信号Mn)以16为周期进行周期延拓，取其主值

区(16点)，再对其进行16点DFT(FFT),得到的就是X(〃&)在[0,2组的16点频率域采样

XA幻。

4.思考题:

如果序列x(n)的长度为M,希望得到其频谱X("")在[0,2组上的N点等间隔采样，

当N<M时，如何用一次最少点数的DFT得到该频谱采样？

5.实验报告及要求

a)运行程序打印要求显示的图形，。

b)分析比较实验结果，简述由实验得到的主要结论

O简要回答思考题

d)附上程序清单和有关曲线。

实验程序清单

1时域采样理论的验证程序清单

clc;

Tp=64/1000;

Fs=1000;T=l/Fs;

M=Tp*Fs；

n=0:M-1;

A=;alph=pi*50*2A;omega=pi*502A;

xnt=A*exp(-alph*n*T).*sin(omega*n*T);

Xk=T*fft(xntzM);

yn='xa(nT)';subplot(3,2,1);

stem(xnt);

boxon;title(1(a)Fs=1000Hz*);

k=0:M-1;fk=k/Tp;

1

subplot(3,2,2);plot(fk,abs(Xk));title((a)T*FT[xa(nT)]zFS=1000Hz');

xlabel('f(Hz),);ylabel(*-u^E');axis([0,Fs,0,*max(abs(Xk))])

%%

%Fs=300Hz

FS=300;T=l/Fs;

M=Tp*Fs；n=0:M-1;

A=;alph=pi*50*2A;omega=pi*50*2A;

xnt=A*exp(-alph*n*T).*sin(omega*n*T);

Xk=T*fft(xnt,M);

yn=1xa(nT)';subplot(3,2,3);

stem(xnt);

hex(ri)FS=3D0H7.');

k=0:M-1;fk=k/Tp;

subplot(3,2,4);plot(fk,abs(Xk))1(a)7*FT[xa(nT)],FS=300Hz');

xlabel('f(Hz)1);ylabel(**);axis([0,Fs,0,*max(abs(Xk))])

%%

%Fs=200Hz

FS=200;T=l/Fs;

M=Tp*Fs;n=0:M-1;

A=;alph=pi*50*2A;omega=pi*50*2A;

xnt=A*exp<-alph*n*T).*oin(omega*n*T);

Xk=T*fft(xntzM);

yn=1xa(nT)';subplot(3,2,5);

stem(xnt);

boxon;title('(a)Fs=200Hz');

k=0:M-1;fk=k/Tp;

subplot(3,2,6);plot(fk,abs(Xk));title(*(a)7*FT[xa(nT)],FS=200Hz')；

xlabel('f(Hz)');ylabel('-u^IE*);axis([0,Fs,0,*max(abs(Xk))])

2频域采样理论的验证程序清单

M=27;N=32;n=0:M;

xa=O:floor(M/2);

xb=ceil(M/2)-1:-1:0;

xn=[xa,xb];

Xk=fft(xn,1024);

X32k=fft(xnz32);

x32n=ifft(X32k);

XI6k=X32k(1:2:N);

xl6n=ifft(X16k,N/2);

subplot(3,2,2);stem(n,xn);boxon

title('(b)Ey^2"DdADx(n)');xlabel('n1);ylabel(1x(n)');

k=0:1023;wk=2*k/1024;

subplot(3,2,1);plot(wk,abs(Xk));title(*(a)FT[x(n)],);

xlabel('\omega/\pi');ylabel('IX(eAjA\omega)I');

k=0:N/2-l;

subplot(3,2,3);stem(k,abs(X16k));boxon

title('(c)16pa^pOd2ESu,);xlabel('k');ylabel('|X_l_6(k)|;

nl=0:N/2-1;

subplot(3,2,4);stem(nl,xl6n);boxon

title(*(d)16paIDFT：X_l_6(k)]1);xlabel('n1);ylabel(fx_l_6(n)1);

k=0:N-1;

subplot(3,2,5);stem(k,abs(X32k));boxon

title('(e)32paJEp0o2ENu*);xlabel(1k•);ylabel('IX_3_2(k)|1);

nl=0:N-1;

subplot(3,2,6);stem(nl,x32n);boxon

title(1(f)32paiDFT：x_3_2(k)]1;xlabel(*0');ylabel(*x_3_2(n),);

XM«tX««x*I.MM

9--l'WJL口»7<V«!•-

实验程序运行结果

1时域采样理论的验证程序运行结果如图所示。由图可见，采样序列的频谱的确是以

采样频率为周期对模拟信号频谱的周期延拓。当采样频率为1000Hz时频谱混叠很小；当采

样频率为300Hz时，在折叠频率150Hz附近频谱混叠很严重；当采样频率为200Hz时，在

折叠频率110Hz附近频谱混叠更很严重。

(a)Fs=1000Hz(a)rFT|xa(nT)],Fs=1000Hz

f(Hz)

(b)Fs=300Hz

150n

t100

至50i

01」」」——

051015

n

(c)Fs=200Hz

150f--------------------------------------

P100

u

0510

图时域采样理论的验证程序运行结果如图所示。

(a)FT[x(n)](b)三角波序列x(n)

20020

•&10010

x

nn

0.50102030

W/JIn

点

(c)16点频域采样(d)16IDFTIX^k)]

20020

c

■s-

9100910

x7

oo

102030

kn

①32点止邛七⑼

⑻32点频域采样

图该图验证了频域采样理论和频域采样定理。对信号x(n)的频谱函数X(eT)在［o,2ir：上等

间隔采样N-16时，N点IDFT［X,、,(幻］得到的序列正是原序列x(n)以16为周期进行周期

延拓后的主值区序列:

/⑺=IDFT[XN(k)]N=[Jx(〃+W)风(〃)

r=-«o

由于N<M,所以发生了时域混叠失真，因此。XN(〃)与x(n)不相同，如图图和(d)所示。当

N=32时，如图图和(d)所示，由于N>M,频域采样定理，所以不存在时域混叠失真，因此。

乙/〃)与x(n)相同。

简答思考题

先对原序列x(n)以N为周期进行周期延拓后取主值区序列，

xN(n)=[£x(n+iN)]RN(〃)

再计算N点DFT则得到N点领域采样：

jM

XN(k)=DFT[xN(n)]N=X(e)Z=O,1,2,・・・，N-1

7人

实验三：用FFT对信号作频谱分析

实验指导

1.实验目的

学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析

误差及其原因，以便正确应用FFT。

2.实验原理

用FFT对信号作频谱分析是学习数字信号处理的盂要内容。经常需要进行谱分析的信

号是模拟信号和时域离散信号。对信号进行谱分析的重要问题是频谱分辨率D和分析误差。

频谱分辨率直接和FFT的变换区间N有关，因为FFT能够实现的频率分辨率是2%/N,因

此要求2〃/NWO。可以根据此式选择FFT的变换区间N。误差主要来自于用FFT作频谱

分析时，得到的是离散谱，而信号(周期信号除外)是连续谱，只有当N较大时离散谱的

包络才能逼近于连续谱，因此N要适当选择大一些。

周期信号的频谱是离散谱，只有用整数倍周期的长度作FFT,得到的离散谱才能代表周

期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。

对模拟信号进行谱分析时，首先要按照采样定理将其变成时域离散信号。如果是模拟

周期信号，也应该选取整数倍周期的长度，经过采样后形成周期序列，按照周期序列的谱分

析进行。

3.实验步骤及内容

(1)对以下序列进行谱分析。

MS)：&(〃)

〃+1,0<n<3

8-〃，4<??<7

0,其它n

4一〃，0</:<3

七(〃)=<"-3,4<z?<7

0,其它n

选择FFT的变换区间N为8和16两种情况进行频谱分析。分别打印其幅频特性曲线。

并进行对比、分析和讨论，

(2)对以下周期序列进行谱分析。

/、71

X4(72)=COS—H

展(〃)=cos("〃/4)+cos("〃/8)

选择rrr的变换区间N为8和16两种情况分别对以上序列进行频谱分析。分别打印

其幅频特性曲线。并进行对比、分析和讨论。

(3)对模拟周期信号进行谱分析

九6Q)=cos871t+cos164t+cos20R

选择采样频率冗=64Hz,变换区间N=16,32,64三种情况进行谱分析。分别打印其幅频

特性，并进行分析和讨论。

4.思考题

(1)对于周期序列，如果周期不知道，如何用FFT进行谱分析？

(2)如何选择FFT的变换区间？(包括非周期信号和周期信号)

(3)当N=8时，々(〃)和与伽)的幅频特性会相同吗？为什么？N=16呢？

5.实验报告要求

(1)完成各个实验任务和要求。附上程序清单和有关曲线。

(2)简要回答思考题。

实验程序清单

%第10章实验3程序

%用FFT对信号作频谱分析

clearalkcloseall

%实验内容⑴===================================================

xln=|ones(l,4)];%产生序列向量xl(n尸R4(n)

M=8;xa=l:(M/2);xb=(M/2):-l:l;x2n=(xa,xb]；%产生长度为8的三角波序列x2(n)

x3n=(xb,xa];

X1k8=fft(x1n,8);%计算x1n的8点DFT

Xlkl6=fft(xln,16);(^l9xln的16点DFT

X2k8=fft(x2n.8):%计算xln的8点DFT

X2kl6=fft(x2n,i6);%计算xln的16点DFT

X3k8=fft(x3n.8);%计算xln的8点DFT

X3ki6=fft(x3n,16);%计算xln的16点DFT

%以下绘制幅频特性曲线

subplot(2,2,l);mstem(Xlk8);%绘制8点DFT的幅频特性图

titleC(la)8点DFT[x_l(n：i]');xlabelCw/n')；ylabel('幅度')；

axis([0,2,0*max(abs(XIkg))])

subplot(2,2,3)；ms(cm(Xlkl6);%绘制16点DFT的幅频特性图

titleC(lb)l6点DFHx」(r)]');xlabel('3/n');ylabelf幅度，)；

axis([0,2,0,*max(abs(Xlkl6))])

figure⑵

suhplot(2,2,l);mstem(X2k8);%绘制8点DFT的幅频特性图

titleC(2a)8点DFTlx_2(nj]');xlabel('w/n听ylabd('幅度');

axis([0,2,0»max(abs(X2k8))])

subplot(2,2,2);mstem(X2kl6);%绘制16点DFT的幅频特性图

titlc('(2b)16点DrT[x_2(r)r);xlabcl('«/n)ylabcl('幅度，)；

axis([0,2,0,*max(abs(X2k16))])

subplot(2,2,3)；mstem(X3k8);%绘制8点DFT的幅频特性图

titleC(3a)8点DFT[x_3(n)]');xlabel('w/n)；ylabel('幅度);

axis(fO,2,O,*max(abs(X3k8))])

subplot(2,2,4);mstem(X3kl6);%绘制16点DFT的幅频特性图

liUef(3b)16点DFT[x_3(r)『);xlabel('3/兀，)；ylabel('幅度，)；

axis(|O.2.O.*max(abs(X3k16))J)

%实验内容(2)周期序歹ij洲分析==================================

N=8;n=0:N-l;%FFT的变换区间N=8

x4n=cos(pi*n/4);

x5n=cos(pi*n/4)+cos(pi*iV8);

X4k8=fft(x4n);%计算x4n的8点DFT

X5k8=fft(x5n);%计算x5n的8点DFT

N=16;n=0:N-l;%FFT的变换区间N=16

x4n=cos(pi*n/4);

x5n=cos(pi*n/4)+cos(pi*n/8);

X4k!6=fft(x4n);%计算x4n的16点DFT

X5kl6=fft(x5n);%计算x5n的16点DFT

figure⑶

subplot(2.2,l);mstem(X4kS);%绘制8点DFT的幅频特性图

title('(4a)8点DFT[x_4(n：『);xlabel('3/n上ylabcl('幅度')；

axis([0,2,0,*max(abs(X4k8))])

subplot(2,2,3)；mstem(X4kl6);%绘制16点DFT的幅频特性图

litleC(4b)16点DFT[x_4(r)r):xlabd('3/兀')；ylabel('幅度，)；

axis([0.2.0,*max(abs(X4k16))])

subplot(2,2,2);mstem(X5k8);%绘制8点DFT的幅频特性图

titleC(5a)8点DFT[x_5(n)r);xlabelCw/n力ylabelf幅度

axis([0,2,0,*max(abs(X5k8))])

subplot(2,2,4);mstem(X5kl6);%绘制16点DFT的幅频特性图

iitlcC(5b)16点DFnx_5(r))');xlabd('3/jr，)；y]abel('幅度');

axis([0,2,0,*max(abs(X5k16))])

%实验内容⑶模拟周期信号谱分析========================

figure(4)

Fs=64;T=l/Fs;

N=16;n=0:N-l;%FFT的变换区间N=16

x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6(l)16点采样

X6k16=fft(x6nT);%计算x6nT的16点DFT

X6kl6=fftshift(X6kl6);%将零频率移到频谱中心

Tp=N*T;F=l/Tp;%频率分辨率F

k=-N/2:N/2-l;il<=k*F:%产生16点DFT对应的采样点频率(以零频率为中心)

subplot(3,1,1);stcm(fk,absiX6k16)：.');bc)xon%绘制8点DFT的幅频特性图

titic('(6a)16点IDEnx—GSTHDRabcKRIbOXylabclC幅度工

axis([-N*F/2-l,N*F/2-1,0,*max(abs(X6k16))])

N=32;n=0:N-l;%FFT的变换区间N=16

x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6(t)32点采样

X6k32=fft(x6nT);%计算x6nT的32点DFT

X6k32=fftshift(X6k32);%将零频率移到频谱中心

Tp=N叮:F=l/Tp;%频率分辨率F

k=-N/2:N/2-l;fl<=k*F:%产生16点DFT对应的采样点频率(以零频率为中心)

subplot(3,1,2);stem(fk,abs(X6k32),,?);boxon%绘制8点DFT的幅频特性图

title((()b)32A|Ubllx_6(nl)jr):xlabelCt(Hz));ylabel(i^.®,)；

axis([-N*F/2-1,N*F/2-1,0,*max(abs(X6k32))])

N=64:n=0:N-l;%FFT的变换区间N=16

x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x6(t)64点采样

X6k64=fft(x6nT):%计算x6nT的64点DFT

X6k64=fftshift(X6k64);%将零频率移到频谱中心

Tp=N*T;F=l/Tp;%频率分辨率F

k=-N/2:N/2-l;fk=k*F:%产生16点DFT对应的采样点频率(以零频率为中心)

subplot(3,1,3)；stem(fk,abs(X6k64)；.');boxon%绘制8点DFT的幅频特性图

titleC(6a)64点|DFT[x_6(nT)]|'):xlabd('f(Hz)');ylabel('幅度')；

axis(f-N*F/2-1.N*F/2-1.0*max(abs(X6k64))])

实验程序运行结果

实验3程序运行结果如图所示。

(la)8点DFTIXi：n))

(1b)16点DFT[X[(n)]

(2a)8点DFT[x2(n)](2b)16<^DFT[x2(n))

L」.一」

0552

cn/TT

(38)8点DFT[FS)](3b)16点DFTlfS)]

DoLILII

u)/n

°0

°o051152

W/TT3/TT

1odl10--

-30-20-100102030

f(Hz)

(6a)64点|DFTM(nD]l

11I1

■

f(Hz)

图程序运行结果分析讨论；

请读者注意，用DFT(或FFT)分析频谱，绘制频谱图时，最好将X(k)的自变量k换

算成对应的频率，作为横坐标便于观察频谱。

27r

你二犷£A=0,l,2,・・・,N-1

为了便于读取频率值，最好关于兀归一化，即以⑷/万作为横坐标。

1、实验内容(1)

图(1a)和(1b)说明与(〃)=凡(〃)的8点DFT和16点DFT分别是为(〃)的频谱

函数的8点和16点采样；

因为&(〃)=%((〃+3))8%(〃)，所以，工3(〃)与工2(〃)的8点DFT的模相等，

如图(2a)和(3a)。但是，当N=16时，毛(〃)与工2(〃)不满足循环移位关系，所

以图(2b)和(3b)的模不同。

2、实验内容(2),对周期序列谱分析

匕(〃)=COS]〃的周期为8,所以N=8和N=16均是其周期的整数倍，得到正确

的单一频率正弦波的频谱，仅在“处有1根单一谱线。如图(4b)和(4b)所示c

毛(〃)=©0$("〃/4)+(：050〃/8)的周期为16,所以N=8不是其周期的整

数倍，得到的频谱不正确，如图(5a)所示。N=16是其一个周期，得到正确的频谱，

仅在兀和又处有2根单一谱线，如图(5b)所示。

3、实验内容(3),对模拟周期信号谱分析

x6«)=cos87rt+cos167tt+cos20R

%6⑺有3个频率成分,f\=4Hz厅2=8Hz,&=10Hz。所以不⑺的周

期为。采样频率月=64Hz=16/=8力=6.4£。变换区间N=i6时，观察时

间Tp=16T=,不是七。）的整数倍周期，所以所得频谱不正确，如图（6a）所示。变换区

间N=32,64时，观察时间Tp=,1s,是工6。）的整数周期，所以所得频谱正确，如图（6b）

和（6c）所示。图中3根谱线正好位于4Hz,8Hz,l0Hz处。变换区间N=64时频谱幅

度是变换区间N=32时2倍，这种结果正好验证了用DFT对中期序列谱分析的理论。

注意：

（1）用DFT（或FFT）对模拟信号分析频谱时，最好将X（k）的自变量k换算成对应的

模拟频率旅，作为横坐标绘图，便于观察频谱。这样，不管变换区间N取信号周期的几倍,

画出的频谱图中有效离散谐波谱线所在的频率值不变，如图（6b）和（6c）所示。

f.=^k=—k=—k,k=0,l,2,・・.,N-1

卜

NNTTp

（2）本程序直接画出采样序列N点DFT的模值，实际上分析频谱时最好画出归一化幅

度谱，这样就避免了幅度值随变换区间N变化的缺点。本实验程序这样绘图只要是为了验

证了用DFT对中期序列谱分析的理论。

简答思考题

思考题（1）和（2）的答案请读者在教材3.?节找，思考题（3）的答案在程序运行结

果分析讨论已经详细I可答，

实验四IIR数字滤波器设计及软件实现

实验指导

i.实验目的

（1）熟悉用双线性变换法设计HR数字滤波器的原理与方法；

（2）学会调用MATLAB信号处理工具箱中滤波器设计函数（或滤波器设计分析工具

fdatool）设计各种I1R数字滤波器，学会根据滤波需求确定滤波器指标参数。

（3）掌握IIR数字滤波器的MATLAB实现方法。

（3）通过观察滤波器输入输出信号的时域波形及其频谱，建立数字滤波的概念。

2.实验原理

设计HR数字滤波器•般采用间接法（脉冲响应不变法和双线性变换法），应用最广泛

的是双线性变换法。基本设计过程是：①先将给定的数字滤波器的指标转换成过渡模拟滤波

器的指标；②设计过渡模拟滤波器；③将过渡模拟滤波器系统函数转换成数字滤波器的系

统函数。MATLAB信号处理工具箱中的各种IIR数字滤波器设计函数都是采用双线性变换法。

第六章介绍的滤波器设计函数butter、chobyl.choby2和ellip可以分别被调用来直接

设计巴特沃斯、切比雪夫1、切比雪夫2和椭圆模拟和数字滤波器。本实验要求读者调用如

上函数直接设计IIR数字滤波器。

本实验的数字滤波器的MATLAB实现是指调用MATLAB信号处理工具箱函数filter对给

定的输入信号x(n)进行滤波，得到滤波后的输出信号yin)。

3.实验内容及步骤

(1)调用信号产生函数mstg产生由三路抑制载波调幅信号相加构成的复合信号st,

该函数还会自动绘图显示st的时域波形和幅频特性曲线，如图所示。由图可见，三路信号

时域混叠无法在时域分离，但频域是分离的，所以可以通过滤波的方法在频域分禽，这就是

本实验的目的。

(a)s(t)的波形

(b)s(t)的频谱

超：：：：：：：：：

iOS0.5...........+T.......:..........b..........I........."?..........-：...卜...卜...T

IIIIIIIII

口L・—」」一J・・L・L・—一」・■・,一」■一一■•■■一・1...........

0200400600800100012001400160018002030

fZHz

图三路调幅信号st的时域波形和幅频特性曲线

(2)要求将st中三路调幅信号分离，通过观察st的幅频特性曲线，分别确定可以分

离st中三路抑制载波单频调幅信号的三个滤波器(低通滤波器、带通滤波器、高通滤波器)

的通带截止频率和阻带截止频率。要求滤波器的通带最大衰减为，阻带最小衰减为60dB。

提示：抑制载波单频调幅信号的数学表示式为

s(f)=cos(2;ry；j)cos(24。)=g[cos(240)+cos(2i(/.+/)")]

其中，cos(2万〃)称为载波，艮为载波频率，cos(2;r.")称为单频调制信号，f。为调制正

弦波信号频率，且满足。＞人。由上式可见，所谓抑制载波单频调幅信号，就是2个正弦

信号相乘，它有2个频率成分：和频力+4和差频力-玲，这2个频率成分关于我波频率

「对称。所以，1路抑制载波单频调幅信号的频谱图是关于载波频率。对称的2根谱线，其

中没有载频成分，故取名为抑制载波单频调幅信号。容易看出，图中三路调幅信号的载波频

率分别为250Hz、500Hz.lOOOHzo如果调制信号m(t)具有带限连续频谱,无直流成分，则

s(/)=，〃Q)cos(2万/7)就是一般的抑制载波调幅信号。其频谱图是关于载波频率上对称的

2个边带(上下边带)，在专业课通信原理中称为双边带抑制载波(DSB-SC)调幅信号，简称

双边带(DSB)信号。如果调制信号m(t)有直流成分，则sQ)=〃z(f)cos(2/r//)就是一•般的

双边带调幅信号。其频谱可是关于载波频率fc对称的2个边带(上下边带)，并包含载频成

分。

(3)编程序调用MATLAB滤波器设计函数ellipord和ellip分别设计这三个椭圆滤波

器，并绘图显示其幅频响应特性曲线。

(4)调用滤波器实现函数filter,用三个滤波器分别对信号产生函数mstg产生的信

号st进行滤波，分离出st中的三路不同载波频率的调幅信号yKn)、yz(n)和y3(n),并绘

图显示yl(n)、y2(n)和y3(n)的时域波形，观察分离效果。

4.信号产生函数mstg清单

functionst=mstg

生产生信号序列向量St,并显示st的时域波形和频谱

%st=mstg返回三路调幅信号相加形成的混合信号，长度N=1600

N=1600%N为信号st的长度。

Fs=10000;T=l/Fs;Tp=N*T;%采样频率Fs=10kHz,Tp为采样时间

t=0:T:(N-l)*T;k=0:N-l;f=k/Tp;

fcl=Fs/10;%第1路调幅信号的载波频率fcl=1000Hz,

fml=fcl/10;与第1路调幅信号的调制信号频率fml=100Hz

fc2=Fs/20;与第2路调幅信号的载波频率fc2=500Hz

fm2=fc2/10;，第2路调幅信号的调制信号频率fm2=50llz

fc3=Fs/40;与第3路调幅信号的载波频率fc3=250Hz,

fm3=fc3/10;为第3路调幅信号的调制信号频率fm3=25Hz

xtl=cos(2*pi*fml*t).*cos(2*pi*fcl*t);%产生第1路调幅信号

xt2=cos(2*pi*fm2*t).*cos(2*pi*fc2*t);%产生笫2路调幅信号

xt3=cos(2*pi*fm3*t).*cos(2*pi*fc3*t);%产生第3路调幅信号

st=xtl+xt2+x13;%三路调幅信号相加

fxt=fft(st,N);%计算信号st的频谱

$二==以下为绘图部分，绘制st的时域波形和幅频特性曲线===================

subplot(3,1,1)

plot(t,st);grid;xlabelCt/s,)jylabel(,s(t)1);

axis([0,Tp/8,min(st),max(st)]);title()(a)s(t)的波形')

subplot(3,1,2)

stem(f,abs(fxt)/max(abs(fxt))»*.*);grid;title(*(b)s(t)的频谱')

axis([0,Fs/5,(),]);

xlabelCf/Hz*);ylabel(,幅度')

5.实验程序框图如图所示，供读者参考。

图实验4程序框图

6.思考题

（1）请阅读信号产生函数mstg,确定三路调幅信号的载波频率和调制信号频率。

（2）信号产生函数mstg中采样点数N=800,对st进行N点FFT可以得到6根理想谱

线。如果取"10（）0,可否得到6根理想谱线？为什么？N=200（）呢？请改变函数mslg中采样

点数N的值，观察频谱图验证您的判断是否正确。

（3）修改信号产生函数mstg,给每路调幅信号加入载波成分，产生调幅（AM）佶号，

重复本实验，观察AM信号与抑制载波调幅信号的时域波形及其频谱的差别。

提示：AM信号表示式：s（/）=[l+cos（2万/J）]cos（