MATLAB实验教案_matlab实验教案

MATLAB实验教案由刀豆文库小编整理，希望给你工作、学习、生活带来方便，猜你可能喜欢“matlab实验教案”。

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

1.实验目的（1）掌握在时域求系统响应的方法。（2）掌握时域离散系统的时域特性。（3）通过实验判断系统稳定性

（4）掌握利用Z变换对系统进行复频域分析。（5）掌握系统零、极点的绘制方法。

（6）通过复频域分析系统稳定性、频率特性。（7）熟悉Z变换的应用

2.实验设备

 计算机

 MATLAB R2012a仿真软件

3.实验原理

（1）离散系统的时域分析

在时域中，描写系统特性的方法是差分方程和单位脉冲响应，在频域可以用系统函数描述系统特性。已知输入信号可以由差分方程、单位脉冲响应或系统函数求出系统对于该输入信号的响应，利用filter函数或conv函数计算输入信号和系统的单位脉冲响应的线性卷积，求出系统的响应。

系统的时域特性是指系统的线性移不变性质、因果性和稳定性。重点分析实验系统的稳定性，包括观察系统的暂态响应和稳定响应。

系统的稳定性是指对任意有界的输入信号，系统都能得到有界的系统响应，或者系统的单位脉冲响应满足绝对可和的条件。系统的稳定性由其差分方程的系数决定。

实际中检查系统是否稳定，不可能检查系统对所有有界的输入信号、输出是否都是有界输出，或者检查系统的单位脉冲响应满足绝对可和的条件。可行的方法是在系统的输入端加入单位阶跃序列，如果系统的输出趋近一个常数（包括零），就可以断定系统是稳定的。系统的稳态输出是指当n→∞时系统的输出。如果系统稳定，信号加入系统后，系统输出的开始一段称为暂态效应，随着n的加大，幅度趋于稳定，达到稳态输出。

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

（2）离散系统的复频域分析

离散系统的时域方程为 dk0Nky(nk)x(nk)pkK0M其变换域分析如下

y(n)x(n)*h(n)频域系统频率响应为

mx(m)h(nm)Y(e)X(e)H(e)

jjjH(e)Z域 jY(e)X(e)jj

y(n)x(n)*h(n)系统的转移函数为

mx(m)h(nm)Y(z)Y(z)X(Z)i X(z)H（z）H(z)0分解因式H(z)iNpkzdkzi0M=Ki(1(1i1i1NMiz1)),其中，和

izi1i称为零、极点。

在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.实验内容及步骤

h(n)(0.9)u(n)x(n)u(n)u(n10)（1）给定矩形脉冲是脉冲响应为

(n)filter函数对差分方程进行数值求解，画出输入的LTI系统的输入，求输出。y利用序列和脉冲响应。（2）已知差分方程：y(n)y(n1)0.9y(n2)x(n)n20,...,①画出在的脉冲响应。100n20,...,100②画出在的单位阶跃响应。h(n)③判断由表征的这个系统的稳定性。

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

y(n)0.5y(n1)0.25y(n2)x(n)2x(n1)x(n3)n

0n100①画出内的系统的脉冲响应，并确定稳定性。

0②如果这个系统输入为,求在范围内的响应。0.6n)]y()x(n)n[5n3200cos(0.2n)4sin(u(n)（4）求下列序列的Z变换并用“zplane”画出零、极点分布图。，该序列的变换为： }x(n)Z{1,2,1,332 z2zz3123X(z)12zz3z3z

（5）求下列直接型系统函数的零、极点，并将它转换成二阶节形式

10.1z10.3z20.3z30.2z4H(z) 10.1z10.2z20.2z30.5z4y(n（6）给定因果系统要求:)0.9y(n1)x(n)①求系统函数H（z）并画出零、极点示意图。②画出系统的幅频特性和相频特性。()|H(ej)| ③求脉冲响应并画出序列图。h(n)

5.程序清单如下

(1)内容1:调用filter解差分方程 b=[1];a=[1,-0.9];n=-5:50;x=stepseq(0,-5,50)-stepseq(10,-5,50);h=0.9.^n;h(1:5)=0;y=filter(b,a,x);subplot(3,1,1);stem(n,x);ylabel('x(n)');axis([-5 50 min(x)max(x)]);subplot(3,1,2);stem(n,h);ylabel('h(n)');axis([-5 50 min(h)max(h)]);subplot(3,1,3);stem(n,y);xlabel('n');ylabel('y(n)');axis([-5,50,-0.5,8]);（2）涉及的功能函数stepseq function[x,n]=stepseq(n0,n1,n2)n=[n1:n2];x=[(n-n0)>=0];(3)内容2：稳定性分析 b=[1];a=[1,-1,0.9];x=impseq(0,-20,120);n=[-20:120];h=filter(b,a,x);subplot(2,1,1);stem(n,h);title('脉冲响应');xlabel('n');ylabel('h(n)')x=stepseq(0,-20,120);s=filter(b,a,x);subplot(2,1,2);stem(n,s);title('阶跃响应');xlabel('n');ylabel('u(n)')（4）涉及的功能函数impseq function[x,n]=impseq(n0,n1,n2)n=[n1:n2];x=[(n-n0)==0];(5)内容3：稳定性分析 b=[1,2,1];a=[1,-0.5,0.25];x=impseq(0,0,100);n=[0:100];h=filter(b,a,x);subplot(2,1,1);stem(n,h);title('脉冲响应');xlabel('n');ylabel('h(n)')x1=5+3*cos(0.2*pi*n)+4*sin(0.6*pi*n);y=filter(b,a,x1);subplot(2,1,2);stem(n,y);title('在sin和cos作用下的响应');xlabel('n');ylabel('y(n)')（6）内容4参考程序

b=[1, 2, 1, 3];%H(z)分子系数 a=[1, 0, 0, 0];%H(z)分母系数 zplane(b,a);%画零、极点分布图（7）内容5参考程序

num=[1,-0.1,-0.3,-0.3,-0.2];den=[1, 0.1, 0.2 ,0.2 ,0.5];[z,p,k]=tf2zp(num,den);m=abs(p);disp('零点');disp(z);disp('极点');disp(p);disp('增益系数');disp(k);sos=zp2sos(z,p,k);disp('二阶节');disp(real(sos));zplane(num,den)（8）内容6参考程序

b=[1,0];a=[1,-0.9];figure zplane(b,a);[H,w]=freqz(b,a,200,'whole');magH=abs(H(1:101));phaH=angle(H(1:101));hw=w(1:101);figure subplot(2,1,1);plot(hw/pi,magH);xlabel('频率单位：π');ylabel('|H(jw)|');title('幅频响应');subplot(2,1,2);plot(hw/pi,phaH/pi);xlabel('频率单位：π');ylabel('相位单位：π');title('相频响应');%脉冲响应

[h,T]=impz(b,a,50);figure stem(T,h);xlabel('n');ylabel('h(n)');

6.结果分析

（1）调用filter解差分方程的仿真结果如图1.1所示，稳定性分析方面的仿真结

果如图1.2所示

10.50x(n)01h(n)10203040500.508642001020n***50y(n)

图1.1

脉冲响应10.5h(n)0-0.5-1-***20n阶跃响应32u(n)10-2002040n6080100120图1.2

h(n)由图1.2所示的脉冲响应图像可以看出，当n>120时，h(n)的值就为零了，这就意味着系统是稳定的。另一种方法是利用MATLAB的roots函数，如：

z=roots(a);

magz=abs(z);

magz=

0.9487

0,9487 因为每个根的幅度都小于1，在z平面的单位圆内，所以系统是稳定的，稳定性分析仿真结果如图1.3

脉冲响应32h(n)10-10102030506070n在sin和cos作用下的响应4080901006040y(n)20001020304050n60708090100图1.3 由如图1.3所示h(n)脉冲响应图像可以看出当n>10时h(n)的值就为零了，这以为着系统是稳定的，第二个图反映了输入为sin和cos函数的加权叠加的响应y(n)

（2）内容4结果分析如图1.4所示

10.5Imaginary Part03-0.5-1-2-1.5-1-0.5Real Part00.51图1.4复频域分析的零、极点图

（3）内容5中，输入“num”和“den”的分别为分子和分母多项式的系数，计算求得零，极点增益系数和二阶节的系数，在MATLAB的Command Windows中观察。

>>零点

0.9615

-0.5730

-0.1443 + 0.5850i

-0.14430.6997i

-0.5776 + 0.5635i

-0.5776-0.5635i

增益系数

二阶节

1.0000

-0.3885

-0.5509

1.0000

1.1552

0.6511

1.0000

0.2885

0.3630

1.0000

-1.0552

0.7679

系统函数的二阶节形式为

1-0.3885z10.5509z211.1552z10.6511z210.2885z10.3630z2.11.0552z10.7679z2

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

10.80.60.4rtaP0.2 yran0igam-0.2I-0.4-0.6-0.8-1-1-0.500.51Real Part1.5 图

内容5系统函数零极点图

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

10.80.60.4Imaginary Part0.20-0.2-0.4-0.6-0.8-1-1-0.50Real Part0.51图1.6 内容6系统函数零、极点图

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

幅频响应151050

|H(jw)|00.10.20.30.40.50.6频率单位：π相频响应0.70.80.910相位单位：π-0.1-0.2-0.3-0.400.10.20.30.40.50.6频率单位：π0.70.80.91图1.7 内容6系统的频率特性图

内容6脉冲响应序列图如图1.8所示。10.90.80.70.6h(n)0.50.40.30.20.100510152025n3035404550图1.8 内容6脉冲响应序列图

6.实验报告要求

（1）报告中要求给出实验的MATLAB程序，并对每条语句给出注释，说明语句作用。

（2）粘贴实验结果图

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

1.实验目的（1）掌握用窗函数法设计FIR数字滤波器的原理和方法。（2）了解各种不同窗函数对滤波器性能的影响。（3）掌握FIR滤波器的快速卷积实现原理。

（4）学会调用MATLAB函数设计与实现FIR滤波器。

2.实验设备

 计算机

 MATLAB R2012a仿真软件

3.实验原理 目前FIR滤波器的设计方法主要有三种：窗函数法、频率取样法和切比雪夫等波纹逼近的最优化设计方法。常用的是窗函数法和频率取样法。本实验中的窗函数法比较简单，可应用现成的窗函数公式，在技术指标要求不高的时候是比较灵活方便的。它是从时域出发，用一个窗函数截取理想的得到，h(n)hd(n)以有限长序列近

hd(n)似理想的；如果从频域出发，用理想的在单位圆上等间隔取样得到，h(n)Hd根据得到将(ej)逼近理想的HdH(k)H(k)(z)，这就是频率取样法。

H(z)几种常见的窗函数如下。

（1）矩形窗（Rectangle Windows）

调用格式：w=boxcar(n),根据长度n产生一个矩形窗。（2）三角窗（Triangular Windows）

调用格式：w=triang(n),根据长度n产生一个三角窗。（3）汉宁窗（Hanning Windows）

调用格式：w=hanning(n),根据长度n产生一个汉宁窗。（4）海明窗（Hamming Windows）

调用格式：w=hamming(n),根据长度n产生一个海明窗。（5）布拉克曼窗（Blackman Windows）

调用格式：w=blackman(n),根据长度n产生一个矩形窗。（6）恺撒窗（KaiserWindows）

调用格式：w=kaiser(n),根据长度n和影响窗函数旁瓣的参数产生一个恺撒窗。

4.实验内容及步骤

（1）分别用矩形窗、三角窗、汉宁窗、海明窗设计低通数字滤波器。信号采样频率为1000Hz，数字滤波器的截止频率为100Hz，滤波器的阶数为80。（2）用海明窗设计一个FIR低通滤波器，截止频率画出加窗前后冲激响应及滤波器幅频特性。

c=0.45，窗的长度取M=121，5．实验程序清单

（1）内容1 参考程序

clear all;close all;parad=0.2*pi;w1=boxcar(81);w2=triang(81);w3=hanning(81);w4=hamming(81);n=1:1:81;hd=sin(parad*(n-41))./(pi*(n-41));hd(41)=parad/pi;h1=hd.*rot90(w1);h2=hd.*rot90(w2);h3=hd.*rot90(w3);h4=hd.*rot90(w4);[MAG1,RAD]=freqz(h1);[MAG2,RAD]=freqz(h2);[MAG3,RAD]=freqz(h3);[MAG4,RAD]=freqz(h4);subplot(2,2,1);plot(RAD,20*log10(abs(MAG1)));title('矩形窗');grid on;subplot(2,2,2);plot(RAD,20*log10(abs(MAG2)));title('三角窗');grid on;subplot(2,2,3);plot(RAD,20*log10(abs(MAG3)));title('汉宁窗');grid on;subplot(2,2,4);plot(RAD,20*log10(abs(MAG4)));title('海明窗');grid on;（2）内容2 参考程序 M=121;n=[0:1:M-1];wc=0.45*pi;hd=ideal_lp(wc,M);w_rec=(boxcar(M))';w_ham=(hamming(M))';h_rec=hd.*w_rec;h_ham=hd.*w_ham;[db1,mag1,pha1,grd1,w1]=freqz_m(h_rec,[1]);[db2,mag2,pha2,grd2,w2]=freqz_m(h_ham,[1]);subplot(2,2,1);stem(n,hd,'k');grid title('加窗前脉冲响应')axis([0,M-1,-0.2,0.5]);ylabel('hd(n)');subplot(2,2,2);stem(n,hd,'k');grid title('加窗后脉冲响应')axis([0,M-1,-0.2,0.5]);ylabel('h(n)');subplot(2,2,3);plot(w1/pi,db1,'k');grid title('加海明窗前幅频特性')axis([0,1,-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(1:1:501))';w=(w(1:1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);grd=grpdelay(b,a,w);（4）内容2设计功能函数 ideal_lp function hd=ideal_lp(wc,M)alpha=(M-1)/2;n=[0:1:(M-1)];m=n-alpha+eps;hd=sin(wc*m)./(pi*m);

6.实验程序运行结果

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

矩形窗500-50-100-150012汉宁窗1000500-50-100-10001234-1500123434-400-20三角窗-60012海明窗34-200

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

加窗前脉冲响应0.40.40.20-0.2加窗后脉冲响应hd(n)0-0.2050100h(n)0.2050100加海明窗前幅频特性00加海明窗后幅频特性幅度-50幅度00.51-50-100-10000.51图4.2 加窗前后冲激响应及幅频特性

7.实验报告要求（1）附程序清单、打印实验内容，要求绘图显示的曲线图。（2）分析总结实验结果。

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

1.实验目的（1）掌握双线性变换法及冲激响应不变法设计方法及其原理。（2）熟悉用双线性变换法及冲激响应不变法设计IIR数字滤波器的MATLAB编程。（3）观察双线性变换及冲激响应不变法设计的滤波器的频域特性。

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

2.实验设备

 计算机

 MATLAB R2012a仿真软件

3.实验原理

设计IIR数字滤波器一般采用间接法（冲激响应不变法和双线性变换法），应用最广泛的是双线性变换法。（1）数字滤波器设计方法如下。

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

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

利用MATLAB中buttord、cheb1ord、cheb2ord、ellipord等函数。③设计模拟滤波器。

MATLAB信号处理工具箱提供了模拟滤波器设计的完全工具函数：butter,cheby1,cheby2,elip,beelf.用户只需一次调用就可完成低通，高通，带通，带阻滤波器设计。调用格式如：[b,a]=butter(N,c,’ftypc’,’s’),其中，’ftype’为滤波器类型

本实验的数字滤波器的MATLAB实现是指调用MATLAB信号处理工具箱函数filter对指定的输入信号x(t)进行滤波，得到滤波后的输出信号y(n)。

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

①假设模拟滤波器的传递函数具有一阶极点，且分母的阶数高于分子的阶数，将展开成部分分式并对求拉氏反变换。

②使用冲激响应不变法求数字滤波器的冲击响应h(n)。③求h(n)的Z变换。

（3）冲激响应不变法主要用于设计某些要求在时域上能模仿滤波器功能的数字滤波器。主要特点如下。①频率变换坐标是线性的。②具有频谱的周期延拓效应，只能用于限带的频率响应特性，即低通和带通数字滤波器。

（4）双线性变换法是基于对微分方程的积分，利用对积分的数值逼近得到的。所谓双线性是指变换公式中s与z的关系无论是分子部分还是分母部分都是线性的。主要特点如下。

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

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

4.实验内容及步骤

（1）采用冲激响应不变法设计一个低通切比雪夫Ⅰ型数字滤波器，设取样频率为1kHz，在通带上限临界频率是400Hz处衰减不大于0.3dB，在阻带临界频率500Hz处衰减不小于60dB。要求：绘制幅频响应图。

（2）采用双线性变换法设计一个高通切比雪夫Ⅱ型数字滤波器，其通带上限临界频率为1500Hz，阻带临界频率为1000 Hz，抽样频率为2000Hz，在通带内最大衰减为0.3dB，阻带内的最小衰减为20dB。（此练习涉及的切比雪夫Ⅱ型数字滤波器为提高部分）要求：绘制幅频响应图。

（3）首先设计一个模拟巴特沃斯低通滤波器，取样频率为10kHz，在通带截止频率1kHz处衰减不大于1dB，在阻带截止频率1.5kHz处衰减不小于15dB。采用双线性变换法设计一个同参数数字巴特沃斯低通滤波器，通带幅度归一化，使其在w=0处幅度为1。

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

5.参考程序

（1）内容1参考程序

%数字滤波器频率特征转换成模拟滤波器的频率特征 Wp=2*pi*400;Ws=2*pi*600;Rp=0.3;Rs=60;Fs=1000;%选择滤波器的最小阶数

[N,Wn]=cheb1ord(Wp,Ws,Rp,Rs,'s');%创建低通切比雪夫滤波器 [Z,P,K]=cheb1ap(N,Rp);[A,B,C,D]=zp2(Z,P,K);[At,Bt,Ct,Dt]=lp2lp(A,B,C,D,Wn);[num1,den1]=2tf(At,Bt,Ct,Dt);%冲激响应不变法把模拟滤波器转换为数字滤波器 [num2,den2]=impinvar(num1,den1,Fs);%绘制幅频响应图

[H,W]=freqz(num2,den2);plot(W*Fs/2/pi,abs(H));grid;xlabel('Hz');ylabel('幅值');（2）内容2参考程序

%数字滤波器的频率特征转换为模拟滤波器的频率特征 Wp=2*pi*1400;Wp1=2*pi*Wp;Ws=2*pi*1000;Ws1=2*pi*Ws;Rp=0.3;Rs=15;fs=20000;%选择滤波器的最小阶数

[N,Wn]=cheb2ord(Wp1,Ws1,Rp,Rs,'s');%创建低通巴特沃斯滤波器 [Z,P,K]=cheb2ap(N,Rs);[A,B,C,D]=zp2(Z,P,K);%低通到高通的转变

[At1,Bt1,Ct1,Dt1]=lp2hp(A,B,C,D,Wn);%冲激响应不变法把模拟滤波器转换为数字滤波器 [At2,Bt2,Ct2,Dt2]=bilinear(At1,Bt1,Ct1,Dt1,fs);%绘制幅频响应图

[num,den]=2tf(At2,Bt2,Ct2,Dt2);[H,W]=freqz(num,den);plot(W*fs/2/pi,abs(H));grid;xlabel('Hz');ylabel('幅值');（3）内容3参考程序 %数字滤波器指标

wp=0.2*pi;ws=0.3*pi;Rp=1;As=15;%转换成模拟域指标 T=1;Fs=1/T;omegap=(2/T)*tan(wp/2);omegas=(2/T)*tan(ws/2);ep=sqrt(10^(Rp/10)-1);Ripple=sqrt(1/(1+ep*ep));Attn=1/(10^(As/20));n=ceil((log10((10^(Rp/10)-1)/(10^(As/10)-1)))/(2*log10(omegap/omegas)));omegac=omegap/((10^(Rp/10)-1)^(1/(2*n)));[B,A]=butter(n,omegac,'s');W=(0:500)*pi/500;[H]=freqs(B,A,W);mag=abs(H);pha=angle(H);db=20*log10((mag+eps)/max(mag));%双线性变换法数字滤波器

[b,a]=bilinear(B,A,T);[h,w]=freqz(b,a,1000,'whole');h=(h(1:501))';w=(w(1:501))';m=abs(h);p=angle(h);db=20*log10((m+eps)/max(m));%绘制频率特性曲线

subplot(2,2,1);plot(w/pi,mag);xlabel('/pi');ylabel('幅度');title('模拟滤波器幅度响应');grid;axis([0 0.8 0 1]);subplot(2,2,2);plot(w/pi,pha);xlabel('/pi');ylabel('相位');title('模拟滤波器相位响应');grid;axis([0 0.8-4 4]);subplot(2,2,3);plot(w/pi,m);xlabel('/pi');ylabel('幅度');title('数字滤波器幅度响应');grid;axis([0 0.8 0 1]);subplot(2,2,4);plot(w/pi,p);xlabel('/pi');ylabel('相位');title('数字滤波器相位响应');grid;axis([0 0.8-4 4]);

6.运行结果

内容

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

1.41.210.8幅值0.60.40.***50Hz***图4.1 冲激响应不变法设计低通切比雪夫数字滤波器 10.90.80.70.6幅值0.50.40.30.20.***40005000Hz***010000图4.2 双线性变换法设计高通切比雪夫数字滤波器

模拟滤波器幅度响应142模拟滤波器相位响应幅度0.5相位00.40.6/pi数字滤波器幅度响应0.20.80-20-400.40.6/pi数字滤波器相位响应0.20.8142幅度0.5相位00.20.4/pi0.60.80-20-400.20.4/pi0.60.8图4.3 对比模拟滤波器和数字滤波器的幅频特性

7.实验报告要求（1）简述实验原理。（2）打印程序清单。

（3）绘制相应幅频响应图。

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

1.实验目的（1）加深对DFT/FFT的理解，熟悉MATLAB中的有关函数。

（2）学习用FFT对连续信号和时域离散信号进行谱分析的方法。（3）了解FFT在系统分析中的应用。

2.实验设备

 计算机

 MATLAB R2012a仿真软件

3.实验原理

在各种信号序列中，有限长序列占重要地位。对有限长序列可以利用离散傅里叶变换（DFT）进行分析。DFT不但可以很好地反映序列的频谱特性，而且易于用快速算法（FFT）在计算机上进行分析。有限长序列的DFT是其Z变换在单位圆上的等距采样，或者说是序列傅里叶变换的等距采样，因此可以用于序列的谱分析。

所谓谱分析就是计算信号的频谱，包括振幅谱。相位谱和功率谱。在利用FFT对模拟信号进行谱分析时，应将模拟信号离散化以得到离散时间信号，同时考虑谱分析中参数的选择。对信号进行谱分析的重要问题是频谱分析率和分析误差。在保持分辨率不变的情况下，若希望增加所分析的信号的最高频率，或在保持信号最高频率不变的情况下，提高分辨率，唯一的方法是增加在记录长度内的取样点数N，那么N必须满足条件：

N其中，2f0F

f0是连续信号最高频率，F是频率分辨率。可以根据此式选择FFT的变换区间N。误差主要来自于用FFT做频谱分析时，得到的是离散谱，而信号（周期信号除外）是连续谱，只有当N较大时，离散谱的包络才能逼近于连续谱，因此N要适当选择得大一些。

周期信号的频谱是离散谱，只有用整数倍周期的长度做FFT，得到的离散谱才能代表周期信号的频谱。如果不知道信号周期，可以尽量选择信号的观察时间长一些。

MATLAB信号处理工具箱中提供了4个FFT内部函数用于计算DFT和IDFT，它们分别是：

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

4.实验内容及步骤

x(t)2sin(4πt)5cos(8πt)DFT采样率和采样点（1）对模拟信号进行采样，用DFT进行谱分析，如何确定数？

x(t)f1=120Hzcos(2f1t)cos(2πf（2）已知一连续信号为其中，πf2=140Hz，用FFT进行分析时，选取取2t)样频率及能够分辨两个谱峰所需要的点数，绘制频谱。（3）对模拟周期信号进行谱分析

x(t)cos8πtcos16πtcos20πt

选择采样频率fs=64Hz，变换区间N=16，32，64，在三种情况进行谱分析。分别打印其幅频特性，并进行分析和讨论。

5.参考程序

（1）内容1参考程序 N=8;fs=16;n=0:N-1;w=fs*n/N;x=2*sin(4*pi*n/fs)+5*cos(8*pi*n/fs);X=fft(x,N);figure(1)stem(w(1:N/2),abs(X(1:N/2)),'k','filled');grid;xlabel('Hz');ylabel('|X(f)|');title('N=8');% N取6或32时进行对比讨论

（2）内容2参考程序 N=21;L=256;f1=120;f2=140;fs=400;T=1/fs;ws=2*pi*fs;n=0:N-1;x=cos(2*pi*f1*n*T)+cos(2*pi*f2*n*T);X=fftshift(fft(x,L));w=(-ws/2+(0:L-1)*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=11;n=0:N-1;x=cos(2*pi*f1*n*T)+cos(2*pi*f2*n*T);X=fftshift(fft(x,L));subplot(2,1,2);plot(w,abs(X));xlabel('Hz');ylabel('幅度谱');axis([-200,200,0,15]);title('N=11');（3）内容3参考程序 figure(3)fs=64;T=1/fs;N=16;n=0:N-1;

%FFT的变换区间N=16 xnT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);

%对x(t)16点采样 Xk16=fft(xnT);

%计算xnT的16点DFT Xk16=fftshift(Xk16);

%将零频率移到频谱中心 Tp=N*T;F=1/Tp;

%频率分辨率F k=-N/2:N/2-1;fk=k*F;

%产生16点DFT对应的采样点频率（以零频率为中心）

subplot(3,1,1);stem(fk,abs(Xk16),'.');box on %绘制8点DFT的幅频特性图 title('16点|DFT[x(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(Xk16))])N=32;n=0:N-1;

%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=1/Tp;

%频率分辨率F k=-N/2:N/2-1;fk=k*F;

%产生16点DFT对应的采样点频率（以零频率为中心）

subplot(3,1,2);stem(fk,abs(Xk32),'.');box on %绘制8点DFT的幅频特性图 title('32点|DFT[x(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(Xk32))])N=64;n=0:N-1;

%FFT的变换区间N=16 xnT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T);

%对x(t)64点采样 Xk64=fft(xnT);

%计算xnT的64点DFT Xk64=fftshift(Xk64);

%将零频率移到频谱中心 Tp=N*T;F=1/Tp;

%频率分辨率F k=-N/2:N/2-1;fk=k*F;

%产生16点DFT对应的采样点频率（以零频率为中心）

subplot(3,1,3);stem(fk,abs(Xk64),'.');box on%绘制8点DFT的幅频特性图 title(' 64点|DFT[x(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(Xk64))])

6.实验程序运行结果

内容

1、内容2以及内容3的实验程序的运行结果分别如图2.1、图2.2和图2.3，图2.4及图2.5所示

N=82018161412|X(f)|10864200123Hz456

图

2.1 实验内容1采样点数N=8时的频谱

N=612108|X(f)|64200123Hz456图

2.2 实验内容1采样点数N=6时的频谱

N=3280706050|X(f)|40302010001234Hz5678图

2.3 实验内容1采样点数N=32时的频谱

N=21151050-200

幅度谱-150-100-500HzN=1***0151050-200幅度谱-150-100-500Hz50100150200图2.4 实验内容2频谱图

7.分析讨论

（1）实验内容1，根据给定模拟f1=2Hz，f=4Hz,所以满足采样定理

2fs2fh，采取样频率fs=16Hz,F=4-2=2Hz,根据N`2

f0/F,所以N8，当N=6时，不满足采样点数，所以频谱在时域混叠，而当N=32时，满足采样点数,频谱中能够看到有两个频率成分，在f1=2Hz,f2=4Hz处有谱线，如图2.1，2.2，2.3所示。

(2)实验内容2,f1=120Hz,f2=140Hz,取采样频率

fs=400Hz,F=140-120=20Hz,所以N=20，显然，当N=21时，能把f1和f2两个谱峰分辨出来，如图2.4中第一个图所示，当N=11时，两个谱峰分辨不出来，如图2.4第二幅图所示。

16点|DFT[x(nT)]|10幅度50-30-20-100f(Hz)32点|DFT[x(nT)]|102030幅度100-30-20-100f(Hz)64点|DFT[x(nT)]|102030幅度200-30-20-100f(Hz)102030图2.5 实验内容3频谱图

（3）实验内容3，对莫逆周期信号进行谱分析

x(t)cos8tcos16tcos20t x(t)有3个频率成分，f1样频率

=4Hz,f2=8Hz，f3=10Hz,所以x(t)的周期为0.5s，采fs=64Hz。变换区间N=16时，观察时间Tp=16T=0.25s,不是x(t)的整数倍周期，所以所得频谱不正确，如图2.5的第一幅图所示，变换区间N=32，64时，观察时间Tp=0.5s,1s,是x(t)的整数周期，所以所得频谱正确，如图2.5的后两幅图所示，图中3根谱线正好位于4Hz,8Hz,10Hz处。变换区间N=64时，频谱幅度是变换区间N=32时的2倍，这种结果正好验证了用DFT对中期序列谱分析的理论。

8.实验报告要求

（1）完成各个实验任务和要求，附上程序清单和有关曲线。（2）简要回答思考问题