数字信号处理实验三离散时间信号的频域分析Word文档下载推荐.docx

1、 fft函数最通常的应用是计算信号的频谱。考虑一个由100hz和200hz正弦信号构成的信号，受零均值随机信号的干扰，数据采样频率为1000hz。通过fft函数来分析其信号频率成分。t=0:0.001:1;%采样周期为0.001s，即采样频率为1000hzx=sin（2*pi*100*t）+sin（2*pi*200*t）+1.5*rand（1,length（t）;%产生受噪声污染的正弦波信号subplot（2,1,1）;plot（x（1:50）;%画出时域内的信号y=fft（x,512）;%对x进行512点的fftf=1000*（0:256）/512;%设置频率轴（横轴）坐标，1000为采样频

2、率subplot（2,1,2）;plot（f,y（1:257）;%画出频域内的信号实验内容32：频谱泄漏和谱间干扰假设现有含有三种频率成分的信号x（t）=cos（200t）+sin（100t）+cos（50t）用DFT分析x（t）的频谱结构。选择不同的截取长度，观察DFT进行频谱分析十存在的截断效应。试用加窗的方法减少谱间干扰。请分析截取长度对频谱泄漏和频率分辨率的影响，分析不同窗函数对谱间干扰的影响。提示：截断效应使谱分辨率（能分开的两根谱线间的最小间距）降低，并产生谱间干扰；频谱混叠失真使折叠频率（fs/2）附近的频谱产生较大的失真。理论和实践都已证明，加大截取长度可提高频率分辨率；选择合

3、适的窗函数可降低谱间干扰；而频谱混叠失真要通过提高采样频率fs和预滤波来改善。解：取采样频率fs=400Hz，采样信号序列x（n）= x（t）w（n）, n=0,1.N-1;N为采样点数，N=fs*T，T为截取时间长度，w（n）为窗函数。实验取三种长度T1=0.04s，T2=4*0.04s，T3=8*0.04s窗函数分别用矩形窗函数w（n）=RN（n），Hamming 窗。clear;close allfs=400; T=1/fs; %采样频率为400Tp=0.04;N=Tp*fs; %采样点数N1=N,4*N,8*N; %设定三种截取长度供调用st=|X1（jf）|;|X4（jf）|X8（j

4、f）|;%设定三种标注语句供调用%矩形窗截断for m=1:3n=1:N1（m）;xn=cos（200*pi*n*T）+sin（100*pi*n*T）+cos（50*pi*n*T）;%产生采样序列Xk=fft（xn,4096） %4096点DFT，用FFT实现fk=0:4095/4096/T;subplot（3,2,2*m-1）plot（fk,abs（Xk）/max（abs（Xk）;ylabel（st（m,:）;if m=1 title（矩形窗截取）;endend;%加hamming窗改善谱间干扰wn=hamming（N1（m）;%调用函数hamming产生N长hamming窗序列wnxn=（

5、cos（200*pi*n*T）+sin（100*pi*n*T）+cos（50*pi*n*T）.*wnsubplot（3,2,2*m）hamming窗截取2 一维逆快速傅立叶变换y=ifft（x）y=ifft（x,n）实验内容3-3：频域采样定理的验证 （1） 产生一个三角波序列x（n）（2）对M=40，计算x（n）的64点DFT，并图示x（n）和X（k）=DFTx（n）,k=0,1,63。（3）对（2）中所得X（k）在0,2pi上进行32点抽样得X1（k）=X（2k）, k=0,1,31。并图示。（4）求X1（k）的32点IDFT,并图示。即x1（n）=IDFTX1（k）, k=0,1,31（

6、5）采用周期延拓的方法绘出的波形，评述与x（n）的关系，并根据频域采样理论加以解释。M=40;N=64;n=0:M;xa=0:floor（M/2）;xb=ceil（M/2）-1:-1:0;xn=xa,xb; %产生长度为M的三角波序列xnXk=fft（xn,64）; %计算xn的64点FFTX1k=Xk（1:2:N）; %隔点抽取Xk得到X1（k）x1n=ifft（X1k,N/2）; %计算X1k的32点IFFT得到x1（n）nc=0:3*N/2; %取97点进行观察xc=x1n（mod（nc,N/2）+1）; %x1（n）的周期延拓序列subplot（3,2,1）stem（n,xn,.tit

7、le（40点三角波序列x（n）xlabel（nylabel（x（n）subplot（3,2,2）n1=0:N/2-1;stem（n1,x1n,32点IDFTX1（k）x1（n）subplot（3,2,3）k=0:N-1;stem（k,abs（Xk）,64点DFTx（n）kX（k）subplot（3,2,4）k1=0:stem（k1,abs（X1k）,X1（k）subplot（3,2,5）;stem（nc,xc,x1（n）的周期延拓序列x（mod（n.32）由于频域在 0,2上的采样点数N小于x（n）的长度M，所以产生时域混叠现象，不能由X1（k）恢复序列x（n）。只有满足N大于等于M时，可优频

8、域采样X1（k）得到原序列x（n）。这就是频域采样鼎立，请同学们自己编程验证.3 线性调频Z变换离散傅立叶变换（DFT）可以看作信号在Z域上沿单位圆的均匀采样。但在实际应用中，并非整个单位圆上的频谱都有意义。一些情况下，如对于窄带信号，只希望分析信号所在的一段频带等，采样点的轨迹是一条弧线或圆周。这种需求，就导致了线性调频Z变换（Chirp z变换）的出现。Chirp z变换与DFT计算整个频谱的算法不同，它是一种更为灵活的计算频谱的算法，可以用来计算单位圆上任一段曲线的Z变换，作频谱分析时输入的点数和输出的点数可以不相等，从而达到频域“细化”的目的。y=czt（x,m,w,a）y=czt（x

9、）例3-4：利用Chirp z变换计算滤波器h（且h=fir1（30,125/500,boxcar（31）在100Hz200Hz的频率特性，并用图形文件比较CZT函数和FFT函数。h=fir1（30,125/500,boxcar（31）;fs=1000;f1=100;f2=200;m=1024;y=fft（h,1024）;fy=fs*（0:1023）/1024;plot（fy,abs（y）;axis（0,500,0,1.5）;w=exp（-j*2*pi*（f2-f1）/（m*fs）;a=exp（j*2*pi*f1/fs）;z=czt（h,m,w,a）;%产生1024个z值fz=（f2-f1）*

10、（0:1023）/1024+f1;plot（fz,abs（z）4 实验问题回答（1）完成例31的程序并观察信号序列x（t）在时域和频域上分析的特点。（2）在实验内容3-2和3-3中按要求完成程序，并回答3-2和3-3中提出的问题。（2） 完成例34的程序并观察Chirp z变换与fft变换的特点。四实验报告要求1按照实验内容要求完成相关实验程序，并得出相关的实验结果（包括图形结果）。2 回答实验中提出的问题。3 总结本次实验结果，按照实验报告格式要求，书写实验报告。五实验设备PC机，MATLAB软件附录A MATLAB系统的常用概念1、命令窗口在Windows 2000下启动MATLAB系统后

11、，Windows 2000的工作平台上会弹出一个窗口，如下图所示，这个窗口称为MATLAB的命令窗口。MATLAB的命令窗口是用户与MATLAB解释器进行通信的工作环境，提示符表示MATLAB解释器正等待用户输入命令。所有的MATLAB命令、MATLAB函数，以及MATLAB程序都要在这个窗口下运行。在命令窗口，用户可以发出MATLAB命令。例如，为了生成一个3*3的矩阵，可以在提示符下，键入如下的命令：A=1 2 3;4 5 6; 7 8 9方括号命令表示矩阵，空格或逗号将每行的元素分开，而分号将矩阵的各行数值分开。再键入Enter后，MATLAB将回显如下的矩阵：A=1 2 34 5 67

12、 8 9为了求该矩阵的逆矩阵，则只要键入命令？B=inv（A）;MATLAB就将计算出相应的结果。如果不想在命令窗口中显示计算结果，只要如上所示，在该命令后多键入一个分号即可。此时，MATLAB系统只完成该命令所要求的计算任务，其计算结果不回显。这项功能在程序设计中是非常必要的。MATLAB系统也可以说是一种新的语言，该语言十分容易掌握，其结构非常类似于数学式子的书写格式，用户花上很少的时间即可掌握MATLAB的大部分命令。2、 图形窗口MATLAB系统的强大功能之一是其优秀的图形功能。对于任何作图命令，MATLAB将打开另一个窗口来绘制输出图形，这样的窗口在MATLAB系统中被成为图形窗口。

13、在同一个图形窗口中，可以绘制多个图形，也可以生成多个图形窗口，并选择其中的一个图形窗口，在其中绘制图形。生成图形窗口的方法比较多，在没有图形窗口存在时，每个绘图函数都能自动生成一个图形窗口；也可以用figure 命令生成一个新的图形窗口；还可以用命令窗口的File菜单的New子菜单的Figure项来打开一个新的图形窗口。3、搜索路径MATLAB管理着一条搜索路径，它在搜索路径下寻找与命令相关的函数文件。例如，如果在MATLAB提示符下输入example, MATLAB解释器将按照下面的步骤来处理这条字符串：（1） 检查example是不是一个变量；（2） 如果不是，检查example是不是一个

14、内部函数；（3） 如果不是，检查在当前文件夹下是否存在名为example.mex,example.dll或example.m的文件。MEX文件是MATLAB 的执行文件，将优先执行；（4） 如果不存在，检查在MATLAB 的搜索路径的目录下是否存在名为MEX,example.mex,example.dll,或example.m的文件。MEX文件优先执行。用户可以打开路径浏览器（Path Browse）查看MATLAB系统的当前搜索路径，也可以在其中加入自己的路径。4、 文件类型在MATLAB系统中，根据功能可将MATLAB系统所使用的外部文件分为几类，并用不同的扩展名作为其标识，我们用的主要是

15、M文件。M文件以字母m为其扩展名，例如startup.m。一般说来，M文件是以ASCII码文本文件，可以用任何文本编辑器进行编辑。在MATLAB系统中，有两类M文件。一类称为程序M文件，简称M文件；另一类称为函数M文件，或简称为函数，统称为 M文件。M文件的内容是由符合MATLAB语法的语句构成的，函数M文件的第一行必须是以关键字function开始的函数说明语句。两类M文件的共同特征是：在MATLAB命令窗口中的命令提示符下键入文件名，来执行M文件中的所有语句规定的计算任务或完成一定的功能。它们的区别在于以下两方面：第一，程序M文件中创建的变量都是MATLAB工作空间中的变量，工作空间中的其

16、他程序或函数可以共享，而函数M文件中创建的所有变量除了全程变量外，均局限于函数运行空间内的局部变量；第二，函数M文件的调用式中可以有输入参数和输出参数，而程序M文件则没有这种功能。5、语言语法要素MATLAB只管理一种对象矩阵。可以使用下列的任何一种方法在MATLAB环境下创建或输入一个矩阵：1） 显示的输入一个元素序列；2） 用MATLAB的内部函数创建一个矩阵；3） 在M文件中用MATLAB语句创建一个矩阵；4） 从一个外部数据文件中装载并创建一个矩阵。在MATLAB中有两个基本概念：变量和表达式。变量由变量名表示，函数名作为特殊的变量名看待，每个变量名由一个字母后面跟随任意个字母或数字（

17、包括下划线）组成，但MATLAB只能分辨前19个字符。MATLAB能区分组成变量名的大小写字母。MATLAB的语句则是两种形式之一：变量名=表达式或者表达式。在前一种语句形式下，MATLAB将运算的结果赋给“变量名”；而在第二种语句形式下，将运算的结果赋给MATLAB的永久变量ans，每条语句以回车符结束。一般的，运算的结果在命令窗口中显示出来。如果语句的最后一个字符是分号“；”，那末，MATLAB仅执行赋值运算，不再显示运算的结果。与C语言一样，MATLAB将字符串当作数组或矩阵处理。6、MATLAB的基本运算符矩阵运算符矩阵A的转置。A+B,A-B矩阵A与B的和与差。A*B矩阵A与B的乘法

18、。A.*B矩阵A与B的对应元素相乘关系运算符小于小于或等于=大于大于或等于=等于=不等于7、特殊运算符在MATLAB的M文件中，可以加入解释行。解释行的标识符为“%”，该标识符将被作为注解内容。程序执行时，注解被忽略。方括号“ ”用于生成矩阵。特别的，语句A= 生成空矩阵A。行分隔符“；”用于MATLAB语句后时，表示该语句的执行结果不被回显出来，这可避免显示一些不感性趣的结果。冒号“：“最主要的作用是生成向量，从下面的例子中可以看出它的使用方法：j:k 生成向量j,j+1,j+2,ki: 生成向量j,j+i,j+2*i,k，如果jk，则生成空矩阵 A（:,j） 矩阵A的第j列A（I,:） 矩

19、阵A的第I行A（j:k） 向量A（j）,A（j+1）,A（k）,j: 从第 j列到第k列的矩阵子块换行连接符“”，有时一条MATLAB语句会很长，在命令窗口的一行内很可能写不下，此时只要在该语句中加入三连点，再回车即可在下一行接着写该语句。8、MATLAB的常用数学函数三角函数sin正弦函数Cos余弦函数Tan正切函数Asin反正弦函数Atan反正切函数Sinh双曲正弦函数Cosh双曲余弦函数Tanh双曲正切函数Asinh反双曲正弦函数Acosh反双曲余弦函数Atanh反双曲正切函数Acos反余弦函数初等函数Abs实数的绝对值、复数的模、字符串的ASIIC值Angle复数的幅角Sqrt方根函数

20、Real复数的实部Imag复数的虚部Conj复共轭运算Round最邻近整数截断（四舍五入）Ceil不大于自变量的最大整数Rem不小于自变量的最小整数Exp自然指数函数（以e为底）Log自然对数函数（以e为底）log10以10为底的对数函数9、程序流控制与其他的程序设计语言一样，MATLAB语言也提供了条件语句。下面分别予以介绍。1） for循环语句MATLAB也有自己的for循环语句。如果要反复执行的一组语句的循环次数是已知或预定义的，就可以用for循环语句。例如：for I=1:nx（I）=0;这条语句将向量x的前n个元素赋予零值，这里的变量n必须预先给定。每一个for必须与end配对使用。

21、2） while循环语句MATLAB提供有while循环语句，它的作用是在一定的逻辑条件的控制下，不断的循环执行一条或一组语句，直到逻辑条件不满足为止。While 语句的一般形式是While 表达式语句组3） if条件语句和bread语句break 语句用于退出循环体，if条件语句有两种形式，分别是if语句组1else语句组2和 表达式1elseif 表达式2语句组310、MATLAB的在线帮助用户可以随时利用MATLAB的在线帮助查询自己不懂得用法的函数的具体用法，例如：在命令窗口键入help abs 后的显示如下：ABS Absolute value ABS（X） is the absol

22、ute value of the elements of X. When X is complex, ABS（X） is the complex modulus （magnitude） of the elements of X.See also SIGN,ANGLE,UNWRAP.Overloaded methodsHelp sym/abs.m将abs函数的主要用法和用途都列了出来。11 实验上机的具体过程如下：1） 在windows 2000/XP的桌面上找到MATLAB的图标单击进入MATLAB的命令窗口，或在开始菜单里选择程序再找到MATLAB单击也可进入MATLAB的命令窗口。2） 进

23、入命令窗口后，在菜单File中选择open，打开已存在的文件。如果是新文件，在菜单File中选择New，再找到M-file即可打开MATLAB的编辑窗口，在编辑窗口内输入你的源程序后存盘。3） 为了使程序能够在MATLAB中运行，需要在搜索路径中加入你的路径，加入路径的过程如下，在MATLAB的命令窗口中选择File 菜单的Set path选项，则打开一个Path Browser窗口，如下图所示：在该窗口中单击Browse，选择你的路径加入，然后在File中选择Save Path后退出即可。这样，你可以在MATLAB的命令窗口中键入你的源程序文件名，或在编辑窗口中选择Tools中的Run即可编

24、译运行。若编译无错，则可得出结果；若编译有错，则可根据命令窗口的提示进行修改后再编译运行，直至得出正确的结果。附录B 信号处理工具箱函数MATLAB包含了进行信号处理的许多工具箱函数，有关这些工具箱函数的使用可通过Help命令得到。为使用方便，在这里将给出几个常用到的函数的使用说明。函数形式函数功能关于函数参数的说明X=sawtooth（t,width）产生锯齿波或三角波。width用于确定最大值的位置，即从0到2*width函数从-1上升到+1。X=square（t,duty）产生方波Duty用于指定正半周期的比例Y=abs（x）求绝对值当x为复数时，得到的是复数模（幅值），若x为字符串，得

25、到的是各个字符的ASCII码。C=conv（a,b）求卷积求取矢量a和b的卷积,c的长度为a和b的长度和减去1。h,w=freqs（b,a,n）模拟滤波器的频率响应b,a为滤波器的冲击响应s变换的分子和分母多项式的系数，在n个频率点计算频率响应h,f=freqz（b,a,n,Fs）数字滤波器的频率响应。Fs为采样频率，b,a为滤波器的冲击响应的Z变换的分子和分母多项式的系数，该函数的作用是在0Fs/2频率范围内选取n个点（记在f中），并计算相应的频率响应。h,t=impz（b,a,n）数字滤波器的冲击响应b,a为滤波器的冲击响应s变换的分子和分母多项式的系数，计算出冲击响应h，取样点数为n.n

26、,Wn=buttord（Wp,Ws,Rp,Rs,s）Butterworth滤波器阶的选择。Wp和Ws分别为通带和阻带的截止频率，皆大于0小于1。Rp和Rs分别是通带和阻带的波纹系数，s表示模拟域，也可不加s，则为数字域。b,a=butter（n,Wn,ftype,s）Butterworth模拟和数字滤波器设计。设计阶数为n，3dB截止频率为Wn的滤波器，ftype指滤波器的类型，high是高通，stop是带阻，无此参数则是低通，s指模拟域，无则表示数字域，b,a是对应变换的分子分母多项式的系数。n,Wn=cheblord（chebyshevI滤波器阶的选择。b,a=cheby1（n,RpWn,ftype,s）Chebyshev（切比雪夫）I型模拟和数字滤波器设计（通带等波纹）。设计阶数