数字信号处理设计报告.doc

1、 目录1信号课程设计内容、目的31.1课程设计的内容、目的1.2课程设计程序语言的选用2课程设计5 2.1课程设计原理2.2课程设计结构图2.3程序运行及分析3调试及结果讨论323.1调试出现的问题3.2讨论分析结果4课程设计心得体会335参考资料336.附录34信号课程设计报告1.信号课程设计内容、目的1.1 课程设计的内容及目的开学初，鉴于上学期的信号处理课程的学习，学院安排了为期两周的信号课程设计。此次课程设计的目的是：1全面复习课程所学理论知识，巩固所学知识重点和难点，将理论与实践很好地结合起来。2提高综合运用所学知识独立分析和解决问题的能力。3熟练使用一种高级语言进行编程实现。 课程

2、设计的内容包括下几个大块：（一）信号分析1编制信号生成程序，产生下述各序列，绘出它们的时域波形，其中包括：单位抽样序列、矩形序列、三角波序列、反三角波序列、Gaussian（高斯）序列、正弦序列、衰减正弦序列。2. 对上述信号完成下列信号分析1）对三角波序列和反三角波序列，作N=8点的FFT，观察比较它们的幅频特性，说明它们有什么异同？绘出两序列及其它们的幅频特性曲线。在和的尾部补零，作N=16点的FFT，观察它们的幅频特性发生了什么变化？分析说明原因。2)观察高斯序列，固定信号中的参数p=8，令q分别等于2，4，8，观察它们的时域和幅频特性，了解当q取不同值时，对信号序列的时域幅频特性的影响

3、；固定q=8，令p分别等于8，13，14，观察参数p变化对信号序列的时域及幅频特性的影响，观察p等于多少时，会发生明显的泄漏现象，混叠是否也随之出现？记录实验中观察到的现象，绘出相应的时域序列和幅频特性曲线。3）对于正弦序列，取数据长度N分别等于8，16，32，分别作N点FFT，观察它们的的时域和幅频特性，说明它们的差别 ，简要说明原因。4)观察衰减正弦序列的时域和幅频特性，绘出幅频特性曲线，改变采样频率fs,使fs=300Hz， 观察此时的频谱的形状和谱峰出现位置？说明产生现象的原因。3设有一连续时间信号s(t),其由20Hz、220Hz和750Hz的正弦信号叠加而成，分析确定采样频率及数据

4、分析长度，计算并绘出信号的频谱，指出各个频率份量。 （二）、数字滤波器的设计与实现设计和实现IIR或FIR数字滤波器：要求：输入数字滤波器的滤波指标,包括通带截止频率,通带最大衰减,阻带截止频率,阻带最小衰减。设计出相应的数字滤波器。显示得到的滤波器的系数和阶数，求出系统的单位脉冲响应h(n)。并画出设计得到的滤波器的幅频特性图(要有坐标标度)和相频特性。 1对于IIR数字滤波器：使用双线性变换法。可选用BW型、CB型和椭圆型等滤波器类型来进行设计2对于FIR数字滤波器：要对多种窗口（三种及以上）、不同阶数设计结果进行分析比较。（三）、综合应用对于给出的任意一段语音信号，进行分析处理。1.读出

5、语音信号，并显示语音信号的时域波形，然后对读出的信号进行FFT变换，得到信号的频谱特性。2.分析其频谱特性，应用设计的滤波器将主要的频率分量分离出来。3.回放此语音信号，感觉滤波前后的声音有什么变化1.2课程设计程序语言的选用本次信号课程设计要求使用一门高级计算机语言来编程，鉴于在许多高等院校，Matlab已经成为许多课程的基本教学工具，我们指导老师就要求我们使用Matlab来实现这次的课程设计。在开始设计之前，我们先要对matlab有一番的认识。Matlab自1984年由MathsWorks公司推向市场以来，历经20多年的发展和竞争，现已风靡世界。在科研与工程应用领域，Matlab已被广泛的

6、应用于科学研究和解决各种具体的实际问题。可靠的数值计算和符号计算功能、强大的绘图功能、简单易学的语言体系以及为数众多的应用工具箱是Matlab区别于其他科技应用软件的显著标志。1. 数值计算和符号计算功能Matlab以矩阵作为数据操作的基本单位，这使得矩阵计算变得非常简洁、方便、高效。Matlab还提供了十分丰富的数值计算函数，而且采用的数值计算算法都是国际公认的最先进、可靠的算法，其程序由世界一流专家编制和高度优化。2. 绘图功能利用Matlab绘图十分方便，它既可以绘制各种图形，包括二维图形和三维图形，还可以对图形进行修饰和控制，以增强图形的表现效果。Matlab提供了两个层次的绘图操作：

7、一种是对图形句柄进行的低层操作，另一种是建立在低层操作之上的高层绘图操作。3. 语言体系Matlab具有程序结构控制、函数调用、数据结构、输入/输出等程序语言特征。因此，对于从事数值计算、计算机辅助设计和系统仿真等领域的人员来说，用Matlab编程的确是一个理想选择。4. Matlab工具箱Matlab工具箱包含两部分内容：基本部分和各种可选的工具箱。Matlab工具箱分为两大类：功能性工具箱和学科性工具箱。功能性工具箱主要用来扩充其符号计算功能、可视建模仿真功能及文字处理功能等。学科性工具箱专业性比较强，如控制系统工具箱、信号处理工具箱、神经网络工具箱、最优化工具箱等，这些工具箱都是由该领域

8、内学术水平很高的专家编写的，用户可以直接利用这些工具箱进行相关领域的科学研究。2.课程设计2.1课程设计原理根据课程设计的目的和要求，指导老师要求我们使用Matlab语言。在课程设计分析，我们可以使用Matlab的工具箱来实现课程设计所要求的内容。在命令窗口中输入Guide，然后出现界面操作。截图如下所示：然后通过switch语句来实现不同功能函数之间的转换，来实现设计要求。在这次的设计中，包括各种序列函数、序列分析函数、滤波器设计函数等，通过一个主菜单选择所要的操作，编写的程序根据你的选择执行对应的操作。序列分析函数可以使用序列的函数表达式来编写程序，序列分析函数多是通过对已知序列进行傅里叶

9、变换来进行幅频分析，滤波器的部分通过各种不同滤波器的函数来实现，例如巴特沃思滤波器可以通过buttord函数来实现巴特沃思数字滤波器的设计，切比雪夫滤波器通过cheb1ord函数等。2.2课程设计结构图选择操作序列单位脉冲序列矩形脉冲序列三角波序列反三角波序列高斯序列正弦序列衰减正弦序列序列分析下拉菜单滤波器设计巴特沃思滤波器切比雪夫滤波器椭圆滤波器语音信号分析2.3程序运行及分析总菜单菜单界面，只要选择想要的操作，就可以之间调用程序中的编码，执行要求。单位抽样序列单位抽样序列的程序如下所示：n=1:50; x=zeros(1,50); x(1)=1;stem(n,x,.);title(x1)

10、;xlabel(n);ylabel(x(n);这个程序编码中，取单位抽样序列的长度为50，n=1：50表示n的取值从1到50，而且之间间隔都是一个单位长度（即为1）。用zeros语句来实现赋零值，zeros函数是产生一个零矩阵。然后只对x=1时赋1值（注意：MATLAB中数组下标从1开始）。矩形序列矩形序列的程序编码如下所示：n=1:50;x=sign(sign(10-n)+1);stem(n,x,.);title(x2);xlabel(n);ylabel(x(n)自行定义序列长度为50个单位，n=1：50表示n的取值从1到50，而且之间间隔都是一个单位长度（即为1）。这个程序中sign函数实

11、现。三角波序列三角波序列的程序编码如下所示：n1=0:3;x1=n1+1;n2=4:7;x2=8-n2;n=n1 n2x=x1 x2stem(n,x,.);title(x3);xlabel(n);ylabel(x(n);text(1,2.5,x=n+1);text(6,2.2,x=8-n);三角波序列是分段序列，程序中n的取值从0到7，其中一段由0到3，另外一段由4到7，两段分别由不同的表达式，这段编码使用stem函数来实现绘图的功能。还有，使用title、xlabel、ylabel、text等函数对图形进行标注。反三角波序列反三角波序列的程序如下所示：n1=0:3;x1=4-n1;n2=4:

12、7;x2=n2-3;n=n1 n2;x=x1 x2;stem(n,x,.);title(x4);xlabel(n);ylabel(x(n);text(2,2.2,x=4-n);text(5,2.2,x=n-3);和三角波序列的程序一样，反三角波序列也是分段序列，同样分成0到3和4到7两段，两段有不同的算术表达式。同样的，程序中也调用stem、title、text等函数，功能和三角波序列的程序中一样。高斯序列高斯序列的程序如下所示：n=0:15; p=8;q=2;x=exp(-1*(n-p).2/q);stem(n,x,.);title(x5);xlabel(n);ylabel(x(n);tex

13、t(11,0.4,x=exp(-1*(n-8).2/2);高斯序列在编程时自行定义序列长度为16，而且确定p=8，q=2，然后输入高斯序列的表达式x=exp(-1*(n-p).2/q)。同样的，该程序也调用了stem等函数。正弦序列正弦序列的程序如下所示：fs=64;N=16;n=0:N-1;x=sin(n*pi/4);stem(n,x,.);title(x6);xlabel(n);ylabel(x(n);text(11,0.8,x=sin(16*pi*n);取fs=64Hz，N=16。那么n的取值是从0到N-1。调用text函数进行图形标注。衰减正弦序列衰减正弦序列程序如下所示：n=0:15

14、; a=0.1;f=0.0625;x=50*exp(-a*n).*sin(2*pi*f*n); stem(n,x,.);title(x7);xlabel(n);ylabel(x(n);text(9,10,x=50*exp(-a*n).*sin(2*pi*f*n);题目中取A=50,采样周期T=1ms,即fs=1000Hz,f0=62.5,a=100。本程序中去序列长度为16，n的取值从0到15.同时调用text函数进行标注。三角波与反三角波的幅频特性这步幅频特性的程序如下所示：n1=0:3;x1=n1+1n2=4:7;x2=8-n2;n=n1 n2；x=x1 x2m1=0:3;y1=4-m1m

15、2=4:7;y2=m2-3;m=m1 m2；y=y1 y2subplot(4,2,3);stem(n,x,.);title(三角波序列x3);xlabel(n);ylabel(x(n);subplot(4,2,5);stem(abs(fft(x,8),.);title(三角波FFT(N=8);xlabel(n);ylabel(幅度值);subplot(4,2,7);stem(abs(fft(x,16),.);title(三角波FFT(N=16)xlabel(n);ylabel(幅度值);subplot(4,2,4);stem(m,y,.);title(反三角波序列x4);xlabel(n);y

16、label(x(n);subplot(4,2,6);stem(abs(fft(y,8),.);title(反三角波FFT(N=8);xlabel(n);ylabel(幅度值);subplot(4,2,8);stem(abs(fft(y,16),.);title(反三角波FFT(N=16);xlabel(n);ylabel(幅度值);如截图所示，程序中实现的功能包含了不同N值（8，16，20）的时域图以及相对应的幅频分析图。利用subplot函数实现图形的分区显示，subplot（m,n,p）函数将当前图形窗口分成m*n个绘图区，即m行，每行n个绘图区，区号按行优先编号，且选定第p个区为当前活动

17、区。在每一个绘图区允许以不同的坐标系单独绘制图形。Matlab提供了对向量或直接对矩阵进行离散傅里叶变换的函数。在这段程序中调用了fft（x,N）函数计算N点离散傅里叶变换。它是限定向量的长度为N，若x的长度小于N，则不足部分补上0；若大于N，则删去超出N的那些元素。观察演示图形，三角波序列和反三角波序列的8点傅里叶变换图形相同，但N=16时，二者就出现了不同。其中的过程大致如下：序列傅里叶变换将x（n）进行周期延拓，得到周期为N的周期序列，对周期序列进行DFS变换，取其主值序列。P=8，q取不同值的时域和幅频特性这段程序如下所示：n=0:15; x1=exp(-1*(n-8).2/2);x2

18、=exp(-1*(n-8).2/4);x3=exp(-1*(n-8).2/8);subplot(4,2,3);stem(n,x1,.);title(x5(p=8,q=2);xlabel(n);ylabel(x(n);subplot(4,2,4);stem(abs(fft(x1),.);title(FFT(p=8,q=2); xlabel(n);ylabel(幅度值);subplot(4,2,5);stem(n,x2,.);title(x5(p=8,q=4);subplot(4,2,6);stem(abs(fft(x2),.); title(FFT(p=8,q=4);xlabel(n);ylab

19、el(幅度值);subplot(4,2,7); stem(n,x3,.);title(x5(p=8,q=8);xlabel(n);ylabel(x(n);subplot(4,2,8);stem(abs(fft(x3),.);title(FFT(p=8,q=8);xlabel(n);ylabel(幅度值)；如截图所示，当q取不同的值时的时域图以及幅频特性图。通过stem和subplot函数进行图形的显示。利用fft函数对高斯序列进行离散傅里叶变换，从而进行幅频分析。通过比较，当固定p=8,而q在取不同值时它们的时域和幅频图形是不一样，q=8，p取不同值时的时域和幅频特性这段程序如下所示：n=0:

20、15;x1=exp(-1*(n-8).2/8);x2=exp(-1*(n-13).2/8);x3=exp(-1*(n-14).2/8);subplot(4,2,3);stem(n,x1,.);title(x5(p=8,q=8);xlabel(n);ylabel(x(n);subplot(4,2,4);stem(abs(fft(x1),.);title(FFT(p=8,q=8);xlabel(n);ylabel(幅度值);subplot(4,2,5);stem(n,x2,.);title(x5(p=13,q=8);subplot(4,2,6);stem(abs(fft(x2),.); title

21、(FFT(p=13,q=8);xlabel(n);ylabel(幅度值);subplot(4,2,7); stem(n,x3,.);title(x5(p=14,q=8);xlabel(n);ylabel(x(n);subplot(4,2,8);stem(abs(fft(x3),.);title(FFT(p=14,q=8);xlabel(n);ylabel(幅度值);如截图所示，当p取不同的值时的时域图以及幅频特性图。同样通过stem和subplot函数进行图形的显示。利用fft函数对高斯序列进行离散傅里叶变换，从而进行幅频分析。通过比较，当p=8，q=8时发生的泄露现象最明显。当p=14，q=

22、8时出现混叠现象。正弦序列不同N的时域和幅频特性这段程序如下：fs=64;n1=0:8-1;x1=sin(n1*pi/16);n2=0:16-1;x2=sin(n2*pi/16);n3=0:20-1;x3=sin(n3*pi/16);subplot(4,2,3);stem(n1,x1,.);title(x6(N=8);xlabel(n);ylabel(x(n);subplot(4,2,4);stem(abs(fft(x1),.);title(FFT(N=8);xlabel(n);ylabel(幅度值);subplot(4,2,5);stem(n2,x2,.);title(x6(N=16);xl

23、abel(n);ylabel(x(n);subplot(4,2,6);stem(abs(fft(x2),.); title(FFT(N=16);xlabel(n);ylabel(幅度值);subplot(4,2,7); stem(n3,x3,.);title(x6(N=20);xlabel(n);ylabel(x(n);subplot(4,2,8);stem(abs(fft(x3),.);title(FFT(N=20);xlabel(n);ylabel(幅度值);这段程序也是通过相同的函数来实现功能。通过观察截图，可以得到不同N值的情况下的幅频特性，一般根据对频谱分析的分辨率确定N，频率分辨率

24、是指频谱分析中能够分辨的两个相邻频率点谱线的最小距离，若图形变换较缓，分辨率可以低一些，相反地，分辨率要高一点。且N越大，分辨率越高。衰减正弦序列不同频率下的时域和幅频特性程序编码如下：n=0:15; a1=0.1;f1=0.0625;a2=1/3;f2=62.5/300;x1=50*exp(-a1*n).*sin(2*pi*f1*n);x2=50*exp(-a2*n).*sin(2*pi*f2*n);subplot(3,2,3);stem(n,x1,.);title(x7(fs=1000);xlabel(n);ylabel(x(n);subplot(3,2,4);stem(abs(fft(x

25、1),.);title(FFT(fs=1000);xlabel(n);ylabel(幅度值);subplot(3,2,5);stem(n,x2,.);title(x7(fs=300);xlabel(n);ylabel(x(n);subplot(3,2,6);stem(abs(fft(x2),.);title(FFT(fs=300);xlabel(n);ylabel(幅度值);分别取不同的频率fs=1000Hz和fs=300Hz，比较图形可得对图形进行采样是对信号进行时间上量化，时域和频域特性图和采样频都有关系。频率比较大时，时域的衰减时间比较长。观察幅频特性，频率比较大时一个完整信号的得长度比

26、较长。连续时间信号信号频谱程序如下：n=0:20; x1=sin(2*pi*20*n/1500);x2=sin(2*pi*220*n/1500);x3=sin(2*pi*750*n/1500);y=x1+x2+x3;subplot(4,1,2);stem(y,.);title(叠加正弦序列);xlabel(x);ylabel(y);subplot(4,1,3);stem(abs(y),.);title(正弦序列幅度谱);xlabel(n);ylabel(幅度值);subplot(4,1,4);plot(angle(y),.);title(正弦序列相位谱);xlabel(n);ylabel(相位

27、);自行定义序列长度为21，取采样频率为1500Hz（三个信号频率分别是20Hz、200Hz以及750Hz，根据采样定理：采样频率至少要大于信号最高频率的2倍，所以可以取采样频率1500Hz。）。观察截图，第一幅图可以看出由20Hz.200Hz.750Hz的信号叠加，幅度发生变化。椭圆滤波器滤波器选择和数据输入界面滤波器的幅频特性系统单位脉冲响应h(n)参数n,wc,a,b这次操作是选择椭圆低通滤波器，程序如下：wp=input(Please enter wp（比ws小）:);ws=input(Please enter ws:);rp=input(Please enter rp:);rs=in

28、put(Please enter rs:);n,wc=ellipord(wp,ws,rp,rs)b,a=ellip(n,rp,rs,wc,low)freqz(b,a);num=b;den=a;figure(1)impulse(num,den);title(系统单位脉冲响应h(n);程序是通过双线性变换法来实现椭圆滤波器的设计，利用ellip函数设计出对应的数字滤波器，利用ellipord函数求出阶数n和wc。设计时还引用impluse函数来求得系统的单位脉冲响应。b,a =ellip(n,rp,rs,wc,x)设计椭圆数字滤波器，并且返回n和wc值。当x取不同的选项时（包含low,high,b

29、andpass,bandstop）设计出不同的滤波器。n,wc=ellipord(wp,ws,rp,rs)返回设计出滤波器的阶数n和wc。滤波器幅频特性截图中的第一个图形是椭圆低通数字滤波器的损耗函数，可以看出所设计的滤波器是低通滤波器，阶数为3，它的截止频率为0.1rad。巴特沃思滤波器滤波器选择和参数输入界面滤波器的幅频特性系统单位脉冲响应h(n)参数n,wc,a,b这次操作是选择椭圆高通滤波器，程序如下：wp=input(Please enter wp（比ws大）:);ws=input(Please enter ws:);rp=input(Please enter rp:);rs=inp

30、ut(Please enter rs:);n,wc=buttord(wp,ws,rp,rs)b,a=butter(n,wc,high)freqz(b,a);num=b;den=a;figure(1)impulse(num,den);title(系统单位脉冲响应h(n);程序通过双线性变换法实现巴特沃思滤波器的设计。程序中调用了buttord和butter函数来实现滤波器的设计和参数返回。其中，b,a=butter(n,wc,x)实现巴特沃思滤波器的设计，改变x成low,high,bandpass和bandstop可以实现低通等各种滤波器的设计。n,wc=buttord(wp,ws,rp,rs)

31、实现巴特沃思滤波器的阶数n和wc的返回。从滤波器幅频特性的截图以及各种参数的截图可以得到所设计的滤波器是高通滤波器，阶数为5，且它的截止频率为0.5524rad。切比雪夫滤波器滤波器选择和参数输入界面滤波器的幅频特性系统单位脉冲响应h(n)参数n,wc,a,b这次操作是选择切比雪夫带通滤波器的设计，程序如下所示：wp=input(Please enter wp（范围比ws小）:);ws=input(Please enter ws:);rp=input(Please enter rp:);rs=input(Please enter rs:);n,wc=cheb1ord(wp,ws,rp,rs)b

32、,a=cheby1(n,rp,wc,bandpass)freqz(b,a);num=b;den=a;figure(1)impulse(num,den);title(系统单位脉冲响应h(n);程序通过双线性变换法实现切比雪夫数字滤波器的设计。程序中调用了cheb1ord和cheby1函数来实现滤波器的设计和参数返回。其中，b,a=cheby1(n,rp,wc,x)实现切比雪夫滤波器的设计，改变x成low,high,bandpass和bandstop可以实现低通等各种滤波器的设计。n,wc=cheb1ord(wp,ws,rp,rs)实现切比雪夫滤波器的阶数n和wc的返回。从滤波器幅频特性的截图以及

33、各种参数的截图可以得到所设计的滤波器是带通滤波器，阶数为3，且它的截止频率分别为0.4rad和0.6rad。切比雪夫滤波器滤波器选择和参数输入界面滤波器的幅频特性系统单位脉冲响应h(n)参数n,wc,a,b这次选择的是切比雪夫带阻滤波器，程序如下所示：wp=input(Please enter wp（范围比ws大）:);ws=input(Please enter ws:);rp=input(Please enter rp:);rs=input(Please enter rs:);n,wc=cheb1ord(wp,ws,rp,rs)b,a=cheby1(n,rp,wc,stop)freqz(b,

34、a);num=b;den=a;figure(1)impulse(num,den);title(系统单位脉冲响应h(n);程序通过双线性变换法实现切比雪夫数字滤波器的设计。程序中调用了cheb1ord和cheby1函数来实现滤波器的设计和参数返回。从滤波器幅频特性的截图以及各种参数的截图可以得到所设计的滤波器是带阻滤波器，阶数为4，且它的截止频率分别为0.2rad和0.8rad。FIR数字滤波器汉宁窗高通数字滤波器幅频特性图参数输入阶数显示FIR数字滤波器是多编写的程序，该程序如下所示：wp=input(Please enter wp:);ws=input(Please enter ws:);D

35、b=wp-ws;N0=ceil(6.2*pi/Db);n=N0+mod(N0+1,2)wc=(wp+ws)/2/pi;hn=fir1(n-1,wc,high,hanning(n);subplot(4,1,2);stem(hn,.);title(系统单位脉冲响应);xlabel(n);ylabel(hn);subplot(4,1,3);stem(abs(fft(hn),.);title(hn的幅频特性);xlabel(n);ylabel(幅度值);subplot(4,1,4);stem(angle(fft(hn),.);title(hn相位特性); xlabel(w);ylabel(相位);设计

36、出的滤波器是汉宁窗高通滤波器，在程序中通过计算先确定滤波器的阶数，然后在选用适合的滤波器设计方式。由于不是要求设计，所以只是设计出了高通的滤波器，其他几种滤波器的设计方法也是类似的。其中还要注意到，该程序中利用hn=fir1(n-1,wc,high,hanning(n)来实现对系统单位脉冲响应的求解。由截图可以看出此次操作的高通滤波器的阶数是33。音频信号的读取和分析该部分程序如下所示：fs=40000; x1=wavread(bird.wav);sound(x1,40000); y1=fft(x1,1024); f=fs*(0:511)/1024;subplot(4,1,2);plot(x1

37、) title(语音信号的时域波形);xlabel(时间);ylabel(幅值);subplot(4,1,3);plot(abs(y1(1:512) title(语音信号FFT频谱)ylabel(|y1|);subplot(4,1,4);plot(f,abs(y1(1:512);title(语音信号频谱特性)xlabel(时间);ylabel(幅值);这段程序通过wavread函数来读取音频文件，然后调用sound函数来播放音频，然后对所读取的音频信号调用fft函数来对信号进行频谱分析。必须注意的是在读取音频时必须要采用适当的采样频率，音频文件的路径也要注意。在这个程序中，我是把音频文件和M文

38、件放在同一个文件夹中，这样就只要写文件名就行，不必麻烦地改文件的读取路径。语音信号过滤低通滤波器y,fs,bits=wavread(bird.wav);fp=3000;fs=6000;rp=0.5;rs=40;fc=40000;wp=2*fp/fc;ws=2*fs/fc;n,wc=buttord(wp,ws,rp,rs);num,den=butter(n,wc);x=filter(num,den,y);sound(x,fc,bits);X=fft(x);subplot(3,2,3);plot(y);title(滤波前语音信号)subplot(3,2,4);plot(x);title(滤波后语音

39、信号)subplot(3,2,5);plot(abs(X);title(滤波后语音信号频谱)h,f=freqz(num,den,256,fc);subplot(3,2,6);plot(f,abs(h),k);title(滤波器频率响应) 程序编写时先读取音频文件，然后对截止频率归一化的转换，然后使用filter函数将信号通过滤波器，对滤波后的信号进行分析，如截图所示。高通滤波器y,fs,bits=wavread(bird.wav);fp=2000;fs=1000;rp=0.5;rs=40;fc=40000;wp=2*fp/fc;ws=2*fs/fc;n,wc=buttord(wp,ws,rp,

40、rs);num,den=butter(n,wc,high);x=filter(num,den,y);sound(x,fc,bits);X=fft(x);subplot(3,2,3);plot(y);title(滤波前语音信号)subplot(3,2,4);plot(x);title(滤波后语音信号)subplot(3,2,5);plot(abs(X);title(滤波后语音信号频谱)h,f=freqz(num,den,256,fc);subplot(3,2,6);plot(f,abs(h),k);title(滤波器频率响应) 这段函数对音频信号的处理相似于上面所说的，只是滤波器是高通滤波器。带

41、通滤波器y,fs,bits=wavread(bird.wav);fp=1500,4500;fs=1000,6000;rp=0.5;rs=40;fc=40000;wp=2*fp/fc;ws=2*fs/fc;n,wc=buttord(wp,ws,rp,rs);num,den=butter(n,wc,bandpass); x=filter(num,den,y);sound(x,fc,bits);X=fft(x);subplot(3,2,3);plot(y);title(滤波前语音信号)subplot(3,2,4);plot(x);title(滤波后语音信号)subplot(3,2,5);plot(a

42、bs(X);title(滤波后语音信号频谱)h,f=freqz(num,den,256,fc);subplot(3,2,6);plot(f,abs(h),k);title(滤波器频率响应)这次语音信号通过的是带通滤波器，只有频率在1500Hz4500Hz的语音信号才能通过滤波器。带阻滤波器y,fs,bits=wavread(bird.wav);fp=1000,6000;fs=1500,4500;rp=0.5;rs=40;fc=40000;wp=2*fp/fc;ws=2*fs/fc;n,wc=buttord(wp,ws,rp,rs);num,den=butter(n,wc,stop);x=fil

43、ter(num,den,y);sound(x,fc,bits);X=fft(x); subplot(3,2,3);plot(y);title(滤波前语音信号)subplot(3,2,4);plot(x);title(滤波后语音信号)subplot(3,2,5);plot(abs(X);title(滤波后语音信号频谱)h,f=freqz(num,den,256,fc);subplot(3,2,6);plot(f,abs(h),k);title(滤波器频率响应)这次设计的滤波器是带阻滤波器，频率在1500Hz4500Hz的信号都被滤掉了，剩下的信号都可以通过滤波器。3.调试及结果讨论3.1调试出现的问题写完程序之后就要对程序进行调试，出现的问题可以总结为一下几点：（1） M