数字信号处理实验指导书

1、实验一 基于Matlab的数字信号处理基本操作一、试验目的1、熟悉Matlab软件的的运行环境，2、了解Matlab基本命令与操作，熟悉数据输入、输出与编辑方法； 3、用Matlab函数产生离散信号；4、用Matlab进行离散信号的表示与运算；5、运用基本函数绘制常用函数的图形。二、实验环境：1、Matlab 6.52、PC机三、实验方法：Matlab操作界面的默认窗口如图1.1所示，Matlab的通用操作界面包括9个常用的窗口，分别是：命令窗口、历史命令窗口、当前目录浏览器窗口、工作空间浏览器窗口、数组编辑窗口、交互界面分类目录窗口、M文件编辑/调试器窗口、帮助导航/浏览器窗口和程序性能剖析

2、窗口。图1-1进入Matlab 6.5命令窗口（Command Window）,命令窗口会默认出现在Matlab界面的右侧，在命令窗口可键入各种Matlab的命令、函数和表达式，并显示除图形以外的所有运算结果。命令串口也可单独显示，选择菜单“VIEW”“Undock Command Window”命令，可使单独的命令窗口返回Matlab界面。Matlab 6.5版运行时，命令窗口中的每个命令行前会出现提示符“”。与以前的版本不同，命令窗口内显示的字符和数值采用不同的颜色，在默认的情况下，输入的命令、表达式以及计算结果等采用黑色字体，字符串采用红色，“if”、“for”等关键词采用蓝色。命令窗口

3、中命令行的编辑，由于Matlab把命令窗口中输入的所有命令都记录在专门的“历史命令”（Command History）空间中，因此Matlab命令窗口不仅可以对输入的命令进行编辑和运行，而且可以对已输入的命令进行回调、编辑和重运行。命令窗编辑功能：键入和修改程序的方法与通常的文字处理相仿。特殊的功能键为：ESC：恢复命令输入的空白状态：调出下一行命令：调出上一行（历史）命令命令窗编辑功能在程序调试时十分有用。对于已执行过的命令，如要做些修改后重新执行，就可不必重新键入，用键调出原命令做修改即可。主菜单中的编辑项功能：用它可以把屏幕上加深了的文字裁减（Cut）或复制（Copy）下来，然后粘贴（P

4、aste）到任一其他视窗的任何位置上去。主菜单中的视图项功能：用它可以改变屏幕上显示的视窗布局。例如，我们希望只显示命令窗，使它占整个屏幕，依次引出View的下拉菜单，即【View】【Desktop Layout】【Command Window Only】键入“help”，屏幕上将显示系统中已装入的函数库（即子目录）的名称。§键入help子目录名，如help elfun，即得出elfun库中各函数名。§键入help函数名，如help tan2，即得到tan2函数的意义及用法。退出MATLAB有两种方法。一种是键入exit或quit，另一种是用鼠标双击左上角的小方块或单击右上

5、角的×号，后者是非正常退出，该次的输入命令将不记录在“历史命令”中，所以应当尽量避免使用。在命令窗口中，默认情况下数值计算结果的显示格式为：当数值为整数，以整数显示；当数值为实数，以小数后4位的精度情况显示，即以“短（Short）”格式显示；如果数值的有效数字超出了这一范围，则以科学计数法显示结果。注意：数值的显示精度并不代表数值的存储精度。在命令窗口可以直接输入操作命令，对命令窗口进行控制。例如：（1）用于清空命令窗口中的显示内容语法：clc（2）在命令窗口中控制其后每页的显示内容行数。语法：more通常，只要执行了一种绘图命令，就会自动产生图形窗，以后的绘图都在这一个图形窗中进行

6、。如想再建一个或几个图形窗，则可键入figure，MATLAB会新建一个图形窗，并自动给它依次排序。如果要认为规定新图为图3，则可键入figure（3）。如果调看已经存在的图形窗3，也应键入figure（3）。MATLAB程序编制有两种方式。一种称为行命令方式，这就是在命令窗中一行一行地输入程序，计算机每次对一行命令作出反应，像计算器那样。这只能编辑简单的程序，在入门时可以用这种方式。程序稍复杂一些，就应把程序写成一个有多行语句组成的文件，让MATLAB来执行这个文件。编写和修改这种文件程序就要用到文本编辑器。命令窗上方最左边的按钮是用来打开文本编辑器空白页的，左边第二个按钮是用来打开原有程序

7、文件的。在命令窗中键入demo，将出现MATLAB的演示窗，演示窗的左侧是库目录。双击目录中的名称或选中该项目后点击右下角的【Run 】方框，即出现该项目的演示界面。通常，演示画面的右侧是一些功能按钮，左上半部是图形，而左下半部则是相应的MATLAB程序语句。还可以在界面上直接修改这些语句并重新执行。因此演示程序也是一个很好的学习过程。四、实验用MATLAB语言工具函数简介本次试验在熟悉软件的运行环境的基础上，要求完成信号的产生以及图形的输出。MATLAB可以根据给出的数据，用绘图命令在屏幕上画出其图形，通过图形对科学计算进行描述。Plot命令用来绘制xy坐标中的曲线。它是一个功能很强的命令。

8、输入变量不同，可以产生很多不同的结果。1、 plot（y）：输入一个数组的情况如果y是一个数组，函数plot（y）给出线性直角坐标的二维图，以y中元素的下标作为X坐标，y中元素的值作为Y坐标，一一画出在XY坐标平面图上，而且将各点以直线相连。例如，要画出是个随机数的曲线。可列出：y5*(rand(1,10)-.5)由Rand函数产生的随机数的最大值为1，最小数为0，平均值为0.5。所以y的最大值为2.5，最小值为2.5，平均值为0。键入plot（y），MATLAB会产生一个图形窗，自动规定最合适的坐标比例绘图。X方向是下标，从110，Y方向范围则是44，并自动标出刻度。可以用title命令给图

9、加上标题，用xlabel，ylabel命令给坐标轴加上说明，用text或gext命令可在图上任何位置加标注，也可用grid命令在图上打上坐标网格线。MATLAB实现程序和形成的图如下：y=5*(rand(1,10)-.5)y = Columns 1 through 9 2.2506 -1.3443 0.5342 -0.0701 1.9565 1.3105 -0.2177 -2.4075 1.6070 Column 10 -0.2765>> plot(y)>> title('my first plot')>> xlabel('x'

10、;),ylabel('y')>> grid2、 plot（x，y）：输入两个数组的情况如果数组x和y具有相同长度，命令plot（x，y）将绘出以x元素为横坐标，y元素为纵坐标的曲线。3、 用plot（t，y1，y2，）命令该语句中t是向量，y=y1，y2，是矩阵，若t是列（行）向量，则y的列（行）长与t长度相同。Y的行（列）数就是曲线的根树。例如，>> t=0:0.5:4*pi;>> y=exp(-0.1*t).*sin(t);>> y1=exp(-0.1*t).*sin(t+1);>> plot(t,y;y1)4、

11、用hold命令在画完前一张图后用hold命令保持住，再画下一张曲线。如键入Plot（t，y），hold on，plot（t，y1，g）执行此命令时，图形窗产生第一幅图形，同时，命令屏幕显示Current plot held，图形处于保持状态。再执行plot（t，y1，g），就把第二幅图以绿色的曲线迭合在同一张图上。5、 图形控制（1）figure：打开图形窗口。MATLAB中的第一幅图随plot命令自动打开，以后的plot命令都画在同一张图上。如要画在另一张新图上，就要用figure命令打开新的图形窗口。有了顺序为1，2，3，的几个图形窗后，再用plot语句，即键入figure（i），表示打开

12、第i幅图。否则，所有的图都会画在最后显示的那张图上。（2）clf：清除当前图形窗的内容。（3）hold：保持当前图形窗的内容，再键入hold，就解除冻结。这种拉线开关式的控制有时会造成混乱，可以用hold on和hold off命令以得到确定的状态。（4）close：关闭当前图形窗。close all：关闭所有图形窗。（5）subplot（n，m，p）命令：将图形窗口分为nm个子图，在第p个子图处绘制图形。五、实验内容： 1、利用Matlab产生单位脉冲序列，并画图。2、利用Matlab产生单位阶跃序列，并画图。3、利用Matlab产生实指数序列：，并画图。4、产生被噪声污染的正弦信号。5、实

13、现一个由50 H和120H正弦信号构成的和信号。6、周期为N的正弦序列，且。7、矩形序列选出45个题目作为实验题目。六、实验报告1、列出实验程序清单，并附上必要的程序说明。2、记录实验结果，实验后，对结果进行分析。实验二 用Matlab进行信号的一维快速傅立叶变换（FFT）一、试验目的1、了解傅立叶变换的物理意义：序列的离散傅立叶变换就是序列的离散频谱；2、了解FFT，掌握FFT函数的调用格式；3、应用FFT计算信号的频谱。二、实验原理：1、离散傅里叶变换（DFT）及其主要性质DFT表示离散信号的离散频谱，DFT的主要性质中有奇偶对称特性，虚实特性等。通过实验可以加深理解。例如：实序列的DFT

14、具有偶对称的实部和奇对称的虚部。这可以证明如下：由定义 实序列DFT的这个特性，在本实验中可以通过实指数序列及三角序列看出来。对于单一频率的三角序列来说它的DFT谱线也是单一的，这个物理意义我们可以从实验中得到验证，在理论上可以推导如下：设：，其DFT为： 从而 以上这串式中反映了的直流分量，是的一次谐波，又根据虚实特性，而其它分量均为零。当周期减小时显然的谱只应该在及才有分量，实验者可以通过和上述相同的步骤加以理论证明。由于与相位差，所以它的DFT只包括实部而没有虚部，以上这些性质可在本实验中得到验证。2、利用DFT对信号进行频谱分析DFT的重要应用之一是对时域连续信号的频谱进行分析，称为傅

15、里叶分析，时域连续信号离散傅里叶分析的基本步骤如图2.1所示。0 L P F A/D X D F T 图2.1 时域连续信号离散傅里叶分析的处理步骤其中消混叠低通滤波器LPF（预滤波器）的引入，是为了消除或减少时域连续信号转换成序列时可能出现的频谱混叠的影响。实际工作中，时域离散信号的时宽是很长的甚至是无限长的（例如语音或音乐信号）。由于DFT的需要，必须把限制在一定的时间间隔之内，即进行数据截断。数据的截断相当于加窗处理。因此，在计算的DFT之前，用一个时域有限的窗函数加到上是非常必要的。通过A/D变换器转换成采样序列一时域离散序列。其频谱用表示，它是频率的周期函数，即其中或为的频谱。在实际

16、应用中，消混叠低通滤波器的阻带不可能是无限衰减的，故由周期延拓得到的有非零重叠，即出现混叠现象。由于进行DFT的需要，必须对序列进行加窗处理，即加窗对频域的影响，用周期卷积表示。最后是进行DFT运算。加窗后的DFT是，其中假设窗函数长L小于或等于DFT长度N。有限长序列的DFT相当于傅里叶变换的等间隔采样。便是的离散频率函数因为DFT频率间隔为，且模拟频率和数字频率间的关系为，所以离散的频率函数各频率点对应的模拟频率为：显然频率分辨率为：利用DFT计算频谱，只给出频谱或的频率分量，即频率的采样值，而不可能得到连续的频谱函数。如果在两个离散的谱线之间有一个特别大的频谱分量，就无法检测出来了。为了

17、在保持原来频谱形状不变的情况下，使谱线加密，即使频域采样点数增加，从而使原来看不到的频谱分量变得可以看到，可以通过在信号数据的末端补加一些零值点，使DFT计算周期内点数增加，但又不改变原有的记录数据的方法来实现。3、快速离散傅里叶变换（FFT）快速离散傅里叶变换是计算离散傅里叶变换的一种快速算法，为了提高运算速度，FFT将DFT的计算逐次分解成较小点数的DFT。按时间抽取（Decimation-In-Time(DIT)FFT算法把输入序列按其值为偶数或是奇数分解成越来越短的序列。按频域抽取（Decimation-In-Frequency(DIF)FFT算法是把输出序列按其值是偶数或是奇数来分解

18、成越来越短的序列。三、实验用MATLAB语言工具函数简介MATLAB中计算序列的离散傅里叶变换和逆变换是采用快速算法，利用fft和ifft函数实现。函数fft用来求序列的DFT,调用格式为：Xk=fft(xn,N)其中xn为有限长序列，N为序列xn的长度，Xk为序列xn的DFT.函数ifft用来求IDFT,调用格式为：xn=ifft(Xk,N)其中，Xk为有限长序列，N为序列Xk的长度，xn为序列Xk的IDFT.例2.1：已知序列，试绘制x（n）及它的离散傅里叶变换图。MATLAB实现程序：clear allN=100;n=0:N-1;xn=2*sin(0.48*pi*n)+cos(0.52*

19、pi*n);XK=fft(xn,N);magXK=abs(XK);phaXK=angle(XK);subplot(1,2,1)plot(n,xn)xlabel('n');ylabel('x(n)');title('x(n) N=100');subplot(1,2,2)k=0:length(magXK)-1;stem(k,magXK,'.');xlabel('k');ylabel('|X(k)|');title('X(k) N=100');运行结果： 序列x(n)及利用FFT求解的相应

20、频谱四、实验要求与实验内容：1、一个由50Hz和120Hz正弦信号构成的信号，受零均值随机信号的干扰，数据采样率为1000Hz。现通过fft函数来分析其信号频率成分。2、一50Hz受噪声污染的正弦信号，对其进行频谱分析。3、分别对下列序列进行频谱分析：(1) 实指数序列(2) 复指数序列(3) 周期为N的正弦序列，且(4) 周期为N的余弦序列，且(5) 复合函数列(6) 矩形序列实验中你可以随意选择这六种输入序列中的3种作为实验内容进行实验。五、实验报告1、列出实验程序清单，并附上必要的程序说明。2、记录实验结果，实验后，对结果进行分析。实验三 用MATLAB设计IIR数字滤波器一、试验目的熟

21、悉模拟Batterworth滤波器设计和用双线性变换法设计数字IIR滤波器的方法。二、实验原理：利用双线性变换法设计IIR滤波器，首先要设计出满足指标要求的模拟滤波器的传递函数，然后由通过双线性变换可得所要设计的IIR滤波器的系统函数。如果给定的指标为数字滤波器的指标，则首先要转换成模拟域指标。1、低通数字滤波器设计数字域与模拟域关系 （3.1） （3.2）2、高通数字滤波器设计数字域与模拟域关系 （3.3） （3.4）3、带通数字滤波器设计数字域与模拟域关系 （3.5） （3.6） （3.7）4、带阻数字滤波器设计 （3.8） （3.9） （3.10）模拟滤波器设计：巴特沃兹滤波器的振幅平方

22、函数为 （3.11）其传输函数为 （3.12） （3.13）首先确定技术指标：(1)通带中允许的最大衰减和通带截止频率；(2)阻带允许的最小衰减和阻带起始频率。由式（11）可得： （3.14） （3.15）得到 （3.16） （3.17）再利用上面两式得到令， 则 （3.18）已知，可由式(3.18)求出滤波器的阶数N。求出的N可能有小数部分一般取大于等于N的最小整数。关于3dB截止频率，有时在技术指标中给出，如果没有给出可以按照式（3.16）或式（3.17）求出。根据以上所述，巴特沃兹滤波器的设计步骤为：(1)根据要求，由式(3.18)求出阶数N；(2)由式（3.16）或式(5.17)求出3

23、dB截止频率；(3)由式（3.13）求出N个极点；(4)由式（3.12）写出传递函数。实际设计中，第（3）、（4）两步由以下两步代替：(5)由N可查下表，得归一化低通巴特沃兹滤波器：阶数N分母多项式的因式123456(6)去归一化，即将用代替，得到实际。三、实验用MATLAB语言工具函数简介用双线性变换法设计IIR滤波器是IIR滤波器设计的经典方法，首先根据模拟滤波器的指标设计出相应的模拟滤波器，然后将设计好的模拟滤波器转换成满足给定技术指标的数字滤波器。在MATLAB的数字信号处理工具箱中提供了相应的设计函数，常用的有：1Butterworth滤波器阶数选择函数N,Wn=buttord(Wp

24、,Ws,Rp,Rs)输入参数：Wp通带截止频率，Ws阻带截止频率，Rp通带最大衰减，Rs阻带最小衰减；输出参数：N符合要求的滤波器最小阶数，Wn为Butterworth滤波器固有频率（3dB）。2零极点增益模型到传递函数模型的转换 num,den=zp2tf(Z,P,K);输入参数：Z,P,K分别表示零极点增益模型的零点、极点和增益；输出参数：num,den分别为传递函数分子和分母的多项式系数。3从低通向低通的转换b,a=lp2lp(Bap,Aap,Wn);功能：把模拟滤波器原型转换成截止频率为Wn的低通滤波器。4双线性变换函数bz,az=bilinear(b,a,Fs); 功能：把模拟滤波器

25、的零极点模型转换为数字滤波器的零极点模型，其中Fs是采样频率。例31：用双线性变换法设计一个Butterworth低通滤波器，要求其通带截止频率100Hz，阻带截止频率200Hz，通带衰减小于2dB，阻带衰减大于15dB，采样频率Fs=500HZ。MATLAB实现程序：%把数字滤波器参数转换成模拟滤波器参数wp=100*2*pi; %利用 ws=200*2*pi;Rp=2;Rs=15;Fs=500;Ts=1/Fs;%选择滤波器的最小阶数N,Wn=buttord(wp,ws,Rp,Rs,'s');%创建butterworth模拟滤波器Z,P,K=buttap(N);%把滤波器零极

26、点模型转化为传递函数模型Bap,Aap=zp2tf(Z,P,K);%把模拟滤波器原型转换成截至频率为Wn的低通滤波器b,a=lp2lp(Bap,Aap,Wn);%用双线性变换法实现模拟滤波器到数字滤波器的转换bz,az=bilinear(b,a,Fs);%绘制频率响应曲线H,W=freqz(bz,az);plot(W*Fs/(2*pi),abs(H);gridxlabel('频率/Hz')ylabel('幅度')运行结果：图5.1 Butterworth低通滤波器的频率响应四、实验内容及要求1、用双线性变换法设计Batterworth低通IIR数字滤波器，要求通

27、带内频率低于时，容许幅度误差在1dB之内，频率在到之间的阻带衰减大与10dB。2、以为采样间隔，在屏幕上打印出数字滤波器的频率区间上的幅频响应特性由线。五、实验报告1、列出实验程序清单，并附上必要的程序说明。2、记录实验结果，实验后，对结果进行分析。实验四 用窗口法设计FIR滤波器一、试验目的了解一个实际滤波器设计过程，加深掌握用窗口法设计FIR滤波器的原理和窗函数对滤波器性能的影响。二、实验原理：设所希望得到的滤波器的理想频率响应为。那么FIR滤波器的设计就在于寻找一个传递函数去逼近。在这种逼近中最直接的一种方法是从单位取样响应序列着手，使逼近理想的单位取样响应。我们知道可以从理想频率响应通

28、过傅里叶反变换来得到，即： （4.1）但是一般来说，这样得到的单位取样响应往往都是无限长序列；而且是非因果的。我们以一个截止频率为的线性相应位理想低通为例来说明。设低通滤波器的时延为，即： （4.2）则 这是一个以为中心的偶对称的无限长非因果序列。这样一个无限长的序列怎样用一个有限长序列去近似呢？最简单的办法就是直接截取它的一段来代替它。例如把到的一段截取来作为，但是为要保证所得到的是线性相位滤波器。必须满足的对称性，所以时延应该取长度的一半，即这种直接截取的办法可以形象地想象为，好比是通过一个“窗口”所看到的一段。中表达为和一个“窗口函数”的乘积。在这里，窗口函数就是矩形脉冲函数，即 但是一

29、般来说，窗口函数并不一定是矩形函数，可以在矩形以内还对作一定的加权处理，因此，一般可以表示为 这里就是窗口函数。这种对理想单位取样响应加窗的处理对频率响应会产生以下三点影响：（1）使理想特性不连续的边沿加宽，形成一过渡带，过渡带的宽度取决于窗口频谱的主瓣宽度。（2）在过渡带两旁产生肩峰和余振，它们取决于窗口频谱的旁瓣；旁瓣越多，余振也越多；旁瓣相对值越大，肩峰则越强。（3）增加截取长度N，只能缩小窗口频谱的主瓣宽度而不能改变旁瓣的相对值；旁瓣与主瓣的相对关系只决定于窗口函灵敏的形状。因此增加N，只能相对应减小过渡带宽。而不能改变肩峰值。肩峰值的大小直接决定通带内的平稳和阻带的衰减，对滤波器性能有很大关系。例如矩形窗的情况下，肩峰达8.95%，致使阻带最小衰减只有21分贝，这在工程上往往是不够的。怎样才能改善阻带的衰减特性呢？只能从改善窗口函数的形状上找出路，所以希望的窗口频谱中应该减少旁瓣，使能量集中在主瓣，这样可以减少肩峰和余振，提高阻带的衰减。而且要求主瓣宽度尽量窄，以获得较陡的过渡带，然而这两个要求总不能兼得，往往需要用增加主瓣宽度带换取决瓣的抑制，于是提出了海明窗、凯宽窗、切比雪夫窗等窗口函数。三、实验用MATLAB语言工