实验一离散系统的时域分析和复频域分析

实验一离散系统的时域分析和复频域分析

1.实验目的

(1)把握在时域求系统响应的方式。

(2)把握时域离散系统的时域特性。

(3)通过实验判定系统稳固性

(4)把握利用Z变换对系统进行复频域分析。

(5)把握系统零、极点的绘制方式。

(6)通过复频域分析系统稳固性、频率特性。

(7)熟悉Z变换的应用

2.实验设备

运算机

MATLABR2021a仿真软件

3.实验原理

(1)离散系统的时域分析

在时域中，描述系统特性的方式是差分方程和单位脉冲响应，在频域能够用

系统函数描述系统特性。已知输入信号能够由差分方程、单位脉冲响应或系统函

数求出系统关于该输入信号的响应，利用filter函数或conv函数计算输入信号

和系统的单位脉冲响应的线性卷积，求出系统的响应。

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

验系统的稳固性，包括观看系统的暂态响应和稳固响应。

系统的稳固性是指对任意有界的输入信号，系统都能取得有界的系统响应,

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

决定。

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

是不是都是有界输出，或检查系统的单位脉冲响应知足绝对可和的条件。可行的

方式是在系统的输入端加入单位阶跃序列，若是系统的输出趋近一个常数(包括

零)，就能够够判定系统是稳固的。系统的稳态输出是指当n-8时系统的输出。

若是系统稳固，信号加入系统后，系统输出的开始一段称为暂态效应，随着n

的加大，幅度趋于稳固，达到稳态输出。

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

(2)离散系统的复频域分析

离散系统的时域方程为

idy{n一幻=Zpg-q

k=()4=0

其变换域分析如下

y(/7)=x(〃)*力(〃)=名x(m)h5-m)o15"")=4(—

m=-»

频域系统频率响应为

〃(「)

Z域

y(/?)=x(〃)*力(力)=Zx(in)h5-in)<=>Yiz)=Z(z)//(z)

系统的转移函数为

加版

Z夕zUu-Jz，

分解因式〃(z)~-=K牛------：一，其中，J和％.称为零、极点。

(1

zLI-2；z)

J=o/=l

在MATLAB中，能够用函数[z,p,K]=tf2zp(num,den)求得有理分式形式的系

统转移函数的零、极点，用函数zplane(z,p)绘出零、极点散布图；也能够用函

数zplane(num,den)直接绘出有理分式形式的系统转移函数的零、极点散布图。

另外，在MATLAB中，能够用函数[r,p,k]=residuez(num,den)完成部份分

式展开计算；能够用函数sos=zp2sos(z,p,K)完成将高阶系统分解为2阶系统的

串联。

判定Z域因果LTI稳固性：当且仅当系统函数的全数极点位于单位圆内时，

因果LTI系统是稳固的。

4.实验内容及步骤

/r-△一〃!!/吨、+x(〃)=〃(〃)一〃(〃—10)曰吨、*鹏士+h(n)=(0.9)“u(〃)

(1)给定矩形脉冲)"7是脉冲响应为

的LTI系统的输入，求输出。利用filter函数对差分方程进行数值求解，

画出输入序列和脉冲口四)一-D+S9ys—2)=”(加

(2)已知差分方程：

①画出在n=-20,…W0的脉冲响应。

②画出在“="2()，',,，10°的单位阶跃响应。

③判定由“加表征的那个系统的稳固性。

(3)一个线性时不变系统由差分方程描述如下：

subplot。,L2);stem(n,h);

ylabel('h(n)');axis([-550min(h)max(h)]);

subplot(3,l,3);stem(n,y);

xlabel('n');ylabel('y(n)');axis([-5,50〃8]);

(2)涉及的功能函数stepseq

function[x,n]=stepseq(n0,nl,n2)

n=[nl:n2];x=[(n-n0)>=0];

(3)内容2:稳固性分析

b=[l];a=[l,-lj；

x=impseq(0/-20,120);n=[-20:120];

h=filter(b,a,x);

subplot(2,l,l);stem(n/h);

title(脉冲响应');xlabel(r');ylabel(h(n)')

x=stepseq(0,-20,120);

s=filter(b,a,x);

subplot(2,lz2);stem(n,s);

title('阶跃响应');xlabel(H);ylabel('u(n)')

(4)涉及的功能函数impseq

function[x,n]=impseq(n0,nlzn2)

n=[nl:n2];x=[(n-n0)==0];

(5)内容3:稳固性分析

b=[l,2,l];a=[lj;

x=impseq(0,0,100);n=[0:100];

h=filter(b,a,x);

subplot(2,lzl);stem(n,h);

title('脉冲响应');xlabel('n');ylabel('h(n)')

xl=5+3*cos*pi*n)+4*sin*pi*n);

y=filter(b,a,xl);

subplot(2,l/2);stem(n/y);

title(在sin和cos作用下的响应');xlabel('n');ylabel('y(n)')

(6)内容4参考程序

b=[l,2.1,3];%H(z)分子系数

a=[l,0,0,0];%H(z)分母系数

zplane(b,a);%画零、极点散布图

(7)内容5参考程序

num=[l,,,,];

den=[l,,

[z,p,k]=tf2zp(num,den);

m=abs(p);

dispC零点);disp⑵;

disp(极点);disp(p);

dispC增益系数);disp(k);

sos=zp2sos(z,p,k);

disp(二阶节');disp(real(sos));

zplane(num,den)

(8)内容6参考程序

b=[l,O];a=[lz];

figure

zplane(b,a);

[H,w]=freqz(b,a,200,,whole,);

magH=abs(H(l:101));phaH=angle(H(l:101));

hw=w(l:101);

figure

subplot(2,l/l);plot(hw/pi,magH);

xlabelC频■位:n')；ylabel('|H(jw)|,)；

title('幅频响应')；

subplot(2,L2);plot(hw/pi,phaH/pi);

xlabel(濒率单位：Tf);ylabel('相位单位：n');

title('相频响应')；

%脉冲响应

[h,T]=impz(b,a,50);

figure

stem(T,h);

xlabel('rT);ylabel('h(n)');

6.结果分析

(1)挪用filter解差分方程的仿真结果如图所示，稳固性分析方面的仿真结果

如图所示

x⑻

图

n

图

由图所示的/脉冲响应图像能够看出，当n>120时，力(加的值就为零了，

这就意味着系统是稳固的。另一种方式是利用MATLAB的roots函数，如：

z=roots(a);

magz=abs(z);

magz=

0,9487

因为每一个根的幅度都小于1,在z平面的单位圆内，因此系统是稳固的，

稳定性分析仿真结果如图

脉冲响应

3

2-

0-

10102030405060708090100

图

由如图所示力(〃)脉冲响应图像能够看出当n>10时力(加的值就为零了,这

以为着系统是稳固的，第二个图反映了输入为sin和cos函数的加权叠加的响应

y(〃)

(2)内容4结果分析如图所示

tEr

d

A

JO

E

U

a

E

E

-

-1

-2-1.5-1-0.500.5

RealPart

图复频域分析的零、极点图

(3)内容5中，输入〃num"和"den”的别离为分子和分母多项式的系

数，计算求得零,极点增益系数和二阶节的系数，在MATLAB的Command

Windows中观看。

“零点

+

极点

增益系数

1

二阶节

系统函数的二阶节形式为

-1

1-0.3885z।550'1+1.1552Z+0.651修

1+0.2885Z'+0«363Qz'1-1.0552z।+0-7679z」

内容5直接型系统函数的零、极点图如所示

8

6

4

0.

0.2

,eE0.

d0.

A0

J-00.

WE

62

E

E0

-4

。6

b

图内容5系统函数零极点图

内容6因果系统的系统函数零，极点图如图所示。

8

00.

6

0.4

0.

e0

e2

d

A0

J0

UE

0

'2

CET

E

-4

6

-0.8

-1

-0.500.5

RealPart

图内容6系统函数零、极点图

内容6系统的幅频特性和相频特性如所示

幅频响应

15

号10

=

I

-5

0

00.10.20.30.40.50.60.70.80.91

频率单位：n

相频响应

0

二-0.1

一

-0.2

-0.3

-0.4

00.10.20.30.40.50.60.70.80.91

频率单位：n

图内容6系统的频率特性图

内容6脉冲响应序列图如图所示。

图内容6脉冲响应序列图

6.实验报告要求

(1)报告中要求给出实验的MATLAB程序，并对每条语句给出注释，说明语

句作用。

(2)粘贴实验结果图

实验二FIR数字滤波器的设计

1.实验目的

(1)把握用窗函数法设计FIR数字滤波器的原理和方式。

(2)了解各类不同窗函数对滤波器性能的阻碍。

(3)把握FIR滤波器的快速卷积实现原理。

(4)学会挪用MATLAB函数设计与实现FIR滤波器。

2.实验设备

运算机

MATLABR2021a仿真软件

3.实睑原理

目前FIR滤波器的设计方式要紧有三种：窗函数法、频率取样法和切比雪夫

等波纹逼近的最优化设计方式。经常使用的是窗函数法和相梆样法。赫脸中

的窗函数法的立简单，可应用毗展窗函数公式，在技术指标要求襁的时候是

H(k)H(k)H(z)

比较灵活方便的。它是从时域动身用一个窗函数截取理想的取得，

以有限长序列近似理想的；若是从撅域动身，用理想的在单

位圆上等距离取样取得，依照取得将逼近理想的〃，,⑵，

这确实是频率取样法。

几种常见的窗函数如下。

(1)矩形窗(RectangleWindows)

挪用格式：w=boxcar(n),依照长度n产生一个矩形窗。。

(2)三角窗(TriangularWindows)

挪用格式：w=triang(n),依照长度n产生一个三角窗。。

(3)汉宁窗(HanningWindows)

挪用格式：w=hanning(n),依照长度n产生一个汉宁窗”。

(4)海明窗(HammingWindows)

挪用格式：w=hamming(n),依照长度n产生一个海明窗①。

(5)布拉克曼窗(BlackmanWindows)

挪用格式：W=blackman(n),依照长度n产生一个矩形窗①。

(6)恺撒窗(KaiserWindows)

挪用格式：w=kaiser(n),依照长度n和阻碍窗函数旁瓣的夕参数产生一个恺撒窗

4.实验内容及步骤

(1)别离用矩形窗、三角窗、汉宁窗、海明窗设计低通数字滤波器。信号采样

频率为1000Hz,数字滤波器的截止频率为100Hz,滤波器的阶数为80。

用海明窗设计一个低通滤波器截止频率。，窗的长度取

(2FIR=nM=121z

回出加窗前后冲激响应及滤波器幅频特性。

5.实验程序清单

(1)内容1参考程序

clearall;closeall;

passrad=*pi;

wl=boxcar(81);

w2=triang(81);

w3=hanning(81);

w4=hamming(81);

n=l:l:81;

hd=sin(passrad*(n-41))./(pi*(n-41));

hd(41)=passrad/pi;

hl=hd.*rot90(wl);

h2=hd.*rot90(w2);

h3=hd.*rot90(w3);

h4=hd.*rot90(w4);

[MAGl/RAD]=freqz(hl);

[MAG2,RAD]=freqz(h2);

[MAG3,RAD]=freqz(h3);

[MAG4,RAD]=freqz(h4);

subplot(2,2,l);plot(RAD20*logl0(abs(MAGl)));title('矩形窗');gridon;

,

subplot(2/2r2);plot(RAD,20*logl0(abs(MAG2)));title(H^W)；gridon;

subplot(223);plot(RAD,20*logl0(abs(MAG3)));title('汉宁窗');gridon;

subplot(2,2,4);plot(RAD,20*logl0(abs(MAG4A);title('海明窗)gridon;

(2)内容2参考程序

M=121;

n=[0:l:M-l];

wc=*pi;

hd=ideal_lp(wc,M);

w_rec=(boxcar(M)),;

w_ham=(hamming(M))';

h_rec=hd.*w_rec;

h_ham=hd.*w_ham;

[dbl/magl,phal,grdl,wl]=freqz_m(h_rec,[l]);

[db2/mag2,pha2,grd2,w2]=freqz_m(h_ham,[l]);

subplot(2,2,l);stem(n,hd,'k');grid

title('加窗前脉冲响应')

axis([O,M-l,,]);ylabel(hd(n),);

subplot(2,2,2);stem(n,hd,'k,);grid

title('加窗后脉冲响应')

axis([O,M-l〃]);ylabel(h(n)');

subplot(2,2,3);plot(wl/pi/dbl/k');grid

title('加海明窗前幅频特性')

axis([0,L-100,10])；ylabel('幅度')；

,

subplot(2,2,4);plot(w2/pi/db2,'k);grid

title(加海明窗后幅频特性')

axis([0,1,-100,10]);ylabel('幅度')；

(3)内容2设计功能函数freqz_m

function[db,mag/pha/grd/w]=freqz_m(b,a);

[H,w]=freqz(b,a,1000,'whole');

H=(H(l:l:501)),;w=(w(l:l:501)),;

mag=abs(H);

db=20*logl0((mag+eps)/max(mag));

pha=angle(H);

grd=grpdelay(b,a,w);

(4)内容2设计功能函数idealjp

functionhd=ideal_lp(wc,M)

alpha=(M-l)/2;

n=[O:l:(M-l)];

m=n-alpha+eps;

hd=sin(wc*m)./(pi*m);

6.实验程序运行结果

内容1和内容2的程序运行结果如图和所示

图多种窗函数设计低通滤波器

加窗前脉冲响应

图加窗前后冲激响应及幅频特性

7.实验报告要求

(1)附程序清单、打印实验内容，要求画图显示的曲线图。

(2)分析总结实验结果。

实验三IIR数字滤波器的设计

1.实验目的

(1)把握双线性变换法及冲激响应不变法设计方式及其原理。

(2)熟悉用双线性变换法及冲激响应不变法设计IIR数字滤波器的MATLAB

(3)观看双线性变换及冲激响应不变法设计的滤波器的频域特性。

(4)熟悉Butterworth滤波器、Chebyshev滤波器和椭圆滤波器的频率特性。

2.实验设备

运算机

MATLABR2021a仿真软件

3.实验原理

设计HR数字滤波器一样采纳间接法(冲激响应不变法和双线性变换法)，应

用最普遍的是双线性变换法。

(1)数字滤波器设计方式如下。

①给定技术指标转换为模拟滤波器设计性能指标。

②求出知足性能指标的模拟相应滤波器性能阶数和截止频率。

利用中、等函数

MATLABbuttordcheblord.cheb2ordxelliporde

③设计模拟滤波器。

MATLAB信号处置工具箱提供了模拟滤波器设计的完全工具函数：

butter,chebyl,cheby2,elip,besself.用户只需一次挪用就可完成低通，高通，带

通，带阻滤波器设计，挪用格式如：[b,a]=butter(N,Q；ftypc's'),其

中，‘代ype'为滤波器类型

本实验的数字滤波器的MATLAB实现是指挪用MATLAB信号处置工具箱函

数filter对指定的输入信号双力进行滤波，取得滤波后的输出信号丫(〃)。

(2)冲激响应不变法设计IIR数字滤波器的步骤如下。

①假设模拟滤波器的传递函数具有一阶极点，且分母的阶数高于分子的阶数，

将展开成部份分式并对求拉氏反变换。

②利用冲激响应不变法求数字滤波器的冲击响应h(n)o

③求h(n)的Z变换。

(3)冲激响应不变法要紧用于设计某些要求在时域上能仿照滤波器功能的数字

滤波器。要紧特点如下。

①频率变换坐标是线性的。

②具有频谱的周期延拓效应，只能用于限带的频率响应特性，即低通和带通数字

滤波器。

(4)双线性变换法是基于对微分方程的积分，利用对积分的数值逼近取得的。

所谓双线性是指变换公式中s与z的关系不管是分子部份仍是分母部份都期性

的。要紧特点如下。

①排除冲激响应不变法所固定的频率混叠现象。

②缺点是模拟频率。和数字频率0之间的非线性关系。

4.实验内容及步骤

(1)采纳冲激响应不变法设计一个低通切比雪夫I型数字滤波器，设取样频率

为1kHz在通带上限临界频率是400Hz处衰减不大于在阻带临界频率500Hz

处衰减不小于60dBo

要求：绘制幅频响应图。

(2)采纳双线性变换法设计一个高通切比雪夫n型数字滤波器，其通带上限临

界频率为1500HZ，阻带临界频率为1000Hz,抽样频率为2000Hz,在通带内

最大衰减为，阻带内的最小衰减为20dBo(此练习涉及的切比雪夫n型数字滤

波器为提高部份)

要求：绘制幅频响应图。

(3)第一设计一个模拟巴特沃斯低通滤波器，取样频率为10kHz,在通带截止

频率1kHz处衰减不大于IdB,在阻带截止频率处衰减不小于15dBo采纳双线

性变换法设计一个同参数数字巴特沃斯低通滤波器，通带幅度归一化，使其在

w=0处幅度为1。

要求：绘制模拟滤波器及数字滤波器的幅频特性曲线进行对照。

5.参考程序

(1)内容1参考程序

%数字滤波器频率特点转换成模拟滤波器的频率特点

Wp=2*pi*400;

Ws=2*pi*600;

Rp=；

Rs=60;

Fs=1000;

%选择滤波器的最小阶数

,

[N,Wn]=cheblord(Wp/Ws/Rp/Rs,'s);

%创建低通切比雪夫滤波器

[Z,BK]=cheblap(N,Rp);

[A,B,C,D]=zp2ss(Z,RK);

[At,Bt,Ct,Dt]=Ip21P(A,B,C,D,Wn);

[numLdenl]=ss2tf(AtzBt,Ct,Dt);

%冲激响应不变法把模拟滤波器转换为数字滤波器

[num2,den2]=impinvar(numl,denl,rs);

%绘制幅频响应图

[HzW]=freqz(num2,den2);

plot(W*Fs/2/pi/abs(H));grid;

xlabel('Hz)ylabelC幅值)；

(2)内容2参考程序

%数字滤波器的频率特点转换为模拟滤波器的频率特点

Wp=2*pi*1400;Wpl=2*pi*Wp;

Ws=2*pi*1000;Wsl=2*pi*Ws;

Rp=；

Rs=15;

fs=20000;

%选择滤波器的最小阶数

,

[N,Wn]=cheb2ord(Wpl,Wsl,Rp/Rs,'s);

%创建低通巴特沃斯滤波器

[Z,BK]=cheb2ap(N,Rs);

[A,B,C,D]=zp2ss(Z/BK);

%低通到高通的转变

[Atl,Btl,Ctl,Dtl]=lp2hp(A,B,C,D,Wn);

%冲激响应不变法把模拟滤波器转换为数字滤波器

[At2,Bt2,Ct2,Dt2]=bilinear(Atl,Btl/Ctl,Dtl,fs);

%绘制幅频响应图

[numfden]=ss2tf(At2,Bt2,Ct2,Dt2);

[H,W]=freqz(num/den);

plot(W*fs/2/pi,abs(H));

grid;

xlabel('Hz)ylabelC幅值)；

(3)内容3参考程序

%数字滤波器指标

wp=*pi;ws=*pi;Rp=l;As=15;

%转换成模拟域指标

T=l;Fs=l/T;

omegap=(2/T)*tan(wp/2);

omegas=(2/T)*tan(ws/2);

ep=sqrt(10A(Rp/10)-l);

Ripple=sqrt(l/(l+ep*ep));

Attn=1/(10A(As/20));

n=ceil((logl0((10A(Rp/10)-l)/(10A(As/10)-l)))/(2*logl0(omegap/omega

s)))；

omegac=omegap/((10A(Rp/10)-l)A(l/(2*n)));

[B,A]=butter^omegac/s');

W=(0:500)*pi/500;

[H]=freqs(B,A,W);

mag=abs(H);

pha=angle(H);

db=20*logl0((mag+eps)/max(mag));

%双线性变换法数字滤波器

[b,a]=bilinear(B,A,T);

[h,w]=freqz(bra,1000,'whole')；

h=(h(l:501)),;

w=(w(l:501)),;

m=abs(h);

p=angle(h);

db=20*logl0((m+eps)/max(m));

%绘制频率特性曲线

subplot(2,2,l);plot(w/pi,mag);xlabel(7pi');ylabel('ll)§^,);

title(模拟滤波器幅度响应');grid;axis([O01]);

subplot(222);plot(w/pi,pha);xlabel('/pi');ylabel('相位')；

title('模拟滤波器相位响应');grid;axis([O-44]);

subplot(2f2,3);plot(w/pifm);xlabel(7pi');ylabeI('it®®');

title('数字滤波器幅度响应');grid;axis([001]);

subplot(2,2,4);plot(w/pi,p);xlabel('/pr);ylabel('相位】

title('数字滤波器相位响应');grid;axis([0-44]);

6.运行结果

内容1、内容2和内容3的运行结果如图、和图所示

1.4

1.2

1

0.8

0.6

0.4

0.2

0

050100150200250300350400450500

Hz

图冲激响应不变法设计低通切比雪夫数字滤波器

图双线性变换法设计高通切比雪夫数字滤波器

模拟滤波器幅度响应模拟滤波器相位响应

/pi/pi

图对照模拟滤波器和数字滤波器的幅频特性

7.实验报告要求

(1)简述实验原理。

(2)打印程序清单。

(3)绘制相应幅频响应图。

实验四DFT/FFT频谱分析及应用

1.实验目的

(1)加深对DFT/FFT的明白得，熟悉MATLAB中的有关函数。

(2)学习用FFT对持续信号和时域离散信号进行谱分析的方式。

(3)了解FFT在系统分析中的应用。

2.实验设备

运算机

MATLABR2021a仿真软件

3.实睑原理

在各类信号序列中，有限长序列占重腹地位。对有限长序列能够利用离散傅

里叶变换（DFT）进行分析。DFT不但能够专门好地反映序列的频谱特性,而且

易于用快速算法（FFT）在运算机上进行分析。有限长序列的DFT是其Z变换

在单位圆上的等距采样，或说是序列傅里叶变换的等距采样，因此能够用于序列

的谱分析。

所谓谱分析确实是计算信号的频谱，包括振幅谱。相位谱和功率谱。在利用

FFT对模拟信号进行诸分析时，应将模拟信号离散化以取得离散时刻信号，同时

考虑谱分析中参数的选择。对信号进行谱分析的重要问题是频谱分析率和分析误

差。在维持分辨率不变的情形下，假设希望增加所分析的信号的最高频率，或在

维持信号最高频率不变的情形下，提高分辨率，唯一的方式是增加在记录长度内

的取样点数N,那么N必需知足条件：

2f

N＞」

F

其中，广。是持续信号最高频率，F是频率分辨率。能够依照此式选择FFT的

变换区间N。误差要紧来自于用FFT做频谱分析时，取得的是离散谱，而信号

（周期信号除外）是持续谱，只有当N较大时，离散谱的包络才能逼近于持续

谱，因此N要适被选择得大一些。

周期信号的频谱是离散谱，只有效整数倍周期的长度做「「T’取得的离散谱

才能代表周期信号的频谱。若是不明白信号周期，能够尽可能选择信号的观看时

刻长一些。

MATLAB信号处置工具箱中提供了4个FFT内部函数用于计算DFT和IDFT,

它们别离是：

fft(x),fft(x,L),ifft(X),ifft(X,L)

4.实验内容及步骤

(1)对模拟信号X")=2sm(4必+58s(8m)进行采样，用DFT

进彳方普分析，如何确信你¥辐部解W4+8s(2袖

(2)已知一持续信号为其中

fl=120Hz,f2=140Hz,用FFT喔葫近但$墩嬲堆跚避舞卷懑情

峰所需要的点数，绘制频谱。

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

选择采样频率fs=64Hz,变换区间N=16,32,64,在三种情形进行谱分析。

别离打印其幅频特性，并进行分析和讨论。

5.参考程序

(1)内容1参考程序

N=8;

fs=16;

n=O:N-l;

w=fs*n/N;

x=2*sin(4*pi*n/fs)+5*cos(8*pi*n/fs);

X=fft(x,N);

figure(l)

,,,

stem(w(l:N/2)/abs(X(l:N/2)),k,'filled);grid;

xlabel('Hz');ylabel('凶f)|');

title('N=8')；

%N取6或32时进行对照讨论

(2)内容2参考程序

N=21;L=256;

fl=120;f2=140;fs=400;

T=l/fs;ws=2*pi*fs;

n=O:N-l;

x=cos(2*pi*fl*n*T)+cos(2*pi*f2*n*T);

X=fftshift(fft(x,L));

w=(-ws/2+(0:L-l)*ws/L)/(2*pi);

figure(2)

subplot(2,1,1);plot(w,abs(X));

xlabel(Hz');ylabel('幅度谱');axis(卜200,200,0,15]);

title('N=21');

N=ll;n=0:N-l;

x=cos(2*pi*fl*n*T)+cos(2*pi*f2*n*T);

X=fftshift(fft(x,L));

subplot(2,l/2);plot(w,abs(X));

xlabel('Hz)ylabel('幅度谱');axis([-200,200,0,15]);

title('N=ir)；

（3）内容3参考程序

figure(3)

fs=64;T=l/fs;

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

xnT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);%对x(t)16点采样

Xkl6=fft(xnT);%计算xnT的16点DFT

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

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

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

频率为中心)

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

title('16点|DFT[x(nT)]|')；xlabelCf(Hz)');ylabel('幅度')；

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

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

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

Xk32=fft(xnT);%计算xnT的32点DFT

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

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

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

频率为中心）

subplot(3,l,2);stem(fk,abs(Xk32)/.');boxon%绘制8点DFT的幅频特性图

title('32点|DFT[x(nT)]|');xlabel('f(Hz))ylabel('幅度')；

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

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