语音信号的数字滤波处理课设

1、 语音信号的数字滤波处理（4）目录一 绪论.21.1设计目的.21.2 设计内容.21.3 设计题目.3二 设计指标要求.3三 预习题.43.1线性卷积.43.2FFT快速卷积.63.3圆周卷积.83.4采样定理演示.11四 课程设计原理.134.1巴特沃斯低通滤波器 .134.2 用窗函数法设计FIR滤波器.14五 设计的步骤和过程.155.1设计方案.155.2设计步骤.16六 用巴特沃斯设计IIR滤波器.306.1巴特沃斯低通滤波器.196.2巴特沃斯带通滤波器.206.3巴特沃斯高通滤波器.23七 用hamming窗设计FIR滤波器.257.1hamming高通滤波器.257.2ham

2、ming低通滤波器.267.3hamming 带通滤波器.28八 心得体会.30参考文献. 31附录A巴特沃斯低通滤波器程序. .32附录B巴特沃斯高通滤波器程序. 34附录C巴特沃斯带通滤波器程序. 37附录D：hamming高通滤波器程序.40附录E：hamming低通滤波器程序 42附录F：hamming带通滤波器程序 44一 绪论本次数字信号处理实习过程中，要求实际用自己设计的数字滤波器分别对被不同噪声污染的信号进行滤波。数字滤波是数字信号处理的重要内容，数字滤波器可分为IIR和FIR两大类。对于IIR数字滤波器的设计，需要借助模拟原型滤波器，再将模拟滤波器转化为数字滤波器，文中采用的

3、设计方法是脉冲响应不变法、双向性变换法和完全函数设计法；对于FIR数字滤波器的设计，可以根据所给定的频率特性直接设计，文中采用的设计方法是窗函数法。根据IIR滤波器和FIR滤波器的特点，在MATLAB坏境下分别用双线性变换法设计IIR和用窗函数设计FIR数字滤波器，并对采集的语音信号进行分析，最后给出了IIR和FIR对语音滤波的效果。1.1 设计目的这次课程设计的目的是为了让我们综合运用数字信号处理和MATLAB去设计一个语音信号去噪的数字滤波器，使学生能够通过实践来分析问题、解决问题、预测目标等。从而不但能够使学生对本门课程有一个深入的了解与运用而且能够让学生学会工程设计的具体步骤和方法，提

4、高分析问题和解决问题的能力，提高实际应用水平。1.2 设计内容（1）设计卷积运算的演示程序可输入任意两个序列x1(n)、x2(n)，指定x1(n)为自己的学号, x1(n)=2,0,1,15,7,0,5,0,2,1,2例如x2(n)= 1, 2.43, 6.17,12.93,22.17,32.25,40.88, 45.87, 45.87, 40.88, 132.25, 22.17, 12.93, 6.17, 2.43,1.0000。分别动态演示两个序列进行线性卷积x1(n)x2(n)和圆周卷积x1(n)x2 (n)的过程；要求分别动态演示翻转、移位、乘积、求和的过程；圆周卷积默认使用

5、两个序列中的最大长度，但卷积前可以指定卷积长度N用以进行混叠分析；改变圆周卷积长度N，根据实验结果分析2类卷积的关系。在计算机操作系统下选一段声音文件(XP系统在“C:WINDOWSMedia”)，读取文件取10ms的声音数据产生时域序列x1(n),序列内容自定义。利用x2(n)= 1, 2.43, 6.17,12.93,22.17,32.25,40.88, 45.87, 45.87, 40.88, 32.25, 22.17, 12.93, 6.17, 2.43,1。利用FFT实现快速卷积，验证时域卷积定理，并与直接卷积进行效率对比（验证时采用matlab子函数）。（2）编写程序演示采样定理（

6、时域采样、频谱周期延拓），同时演示采样频率小于2fc时，产生的混叠效应对下面连续信号进行采样：，A为幅度因子，a为衰减因子，为模拟角频率，其中n为学号(例如，王墨同学n=23)要求输入采样频率fs（根据程序处理需要指定范围）后，在时域演示信号波形、采样脉冲及采样后信号；在频域演示不同采样频率下对应信号的频谱。1.3设计题目利用Windows下的录音机或其他软件，进行语音信号的采集（*.wav）；语音信号的频谱分析，画出采样后语音信号的时域波形和频谱图；产生噪声信号并加到语音信号中，得到被污染的语音信号，并回放语音信号；污染信号的频谱分析，画出被污染的语音信号时域波形和频谱； 根据有关的频谱特性

7、，采用间接法设计IIR数字滤波器，并画出相应滤波器的幅频、相频图（设计3个IIR滤波器）a. 模拟滤波器类型：巴特沃思滤波器（低通、带通、高通）b总体要求：Matlab原程序+仿真波形+技术指标 根据有关的频谱特性，采用直接法设计FIR数字滤波器，并画出相应滤波器的幅频、相频图（设计3个FIR滤波器）a. 滤波器类型：hamming窗（低通、带通、高通）b总体要求：Matlab原程序+仿真波形+技术指标+窗函数 用自己设计的这些滤波器分别对被不同噪声污染的信号进行滤波用自己设计的这些滤波器分别对被不同噪声污染的信号进行滤波；分析得到信号的频谱，画出滤波后信号的时域波形和频谱，并对滤波前后的信号

8、进行对比，分析信号的变化；回放语音信号。二 设计指标要求1. 掌握数字信号处理的基本概念，基本理论和基本方法。2. 熟悉离散信号和系统的时域特性。3. 掌握序列快速傅里叶变换方法。4. 学会MATLAB的使用，掌握MATLAB的程序设计方法。5. 掌握利用MATLAB对语音信号进行频谱分析。6. 掌握滤波器的网络结构。三 预习题3.1 线性卷积x1(n)=2,0,1,1,5,7,0,5,0,2,1,2x2(n)= 2,4,5,6,2,7,4,8,9,6,4,8,3,5,7,8 ，其线性卷积动态演示Matlab程序如下：clf;n1=1:40;n2=1:40;x1=zeros(1,12),2,0

9、,1,1,5,7,0,5,0,2,1,2,zeros(1,16); x2=zeros(1,12),2,4,5,6,2,7,4,8,9,6,4,8,3,5,7,8,zeros(1,12);x11=zeros(size(x1);subplot(3,1,1)stem(n1,x1)subplot(3,1,2)stem(n2 ,x2)pause(2);for i=1:24 x11(i)=x1(25-i);endfor i=25:40 x11(i)=0;endsubplot(3,1,1)stem(n1,x11);pause(2)x3=zeros(size(x11);for i=1:28x11=zeros(

10、1,1),x11(1:(40-1);x3(i)=sum(x11.*x2);subplot(3,1,1)stem(n1,x11)subplot(3,1,3)x4=zeros(1,12),x3(1:i),zeros(1,28-i);stem(n1,x4)pause(2)end 运行结果如图3.1：图3.1 线性卷积3.2 FFT快速卷积x2(n)= 1, 2.43, 6.17,12.93,22.17,32.25,40.88, 45.87, 45.87, 40.88, 32.25, 22.17, 12.93, 6.17, 2.43,1.0000，其FFT快速卷积Matlab程序如下：clf;x1=w

11、avread('c:WindowsMediading.wav',10);x1=x1'x1=x1(1,:);x2=1 2.43 6.17 12.93 22.17 32.25 40.88 45.87 45.87 40.88 32.25 22.17 12.93 6.17 2.43 1;subplot(7,1,1);stem(x1);ylabel('x1(n)');title('x1(n)');subplot(7,1,2);stem(x2);ylabel('x2(n)');title('x2(n)');y=conv

12、(x1,x2);subplot(7,1,3);stem(y);ylabel('y');title('x1(n)与x2(n)的卷积');N1=length(x1);N2=length(x2);N=N1+N2-1;X1=fft(x1,N);X2=fft(x2,N);subplot(7,1,4);stem(X1);ylabel('X1');title('x1(n)的N点DFT');subplot(7,1,5);stem(X2);ylabel('X2');title('x2(n)的N点DFT');Y1=X1

13、.*X2;subplot(7,1,6);stem(Y1);ylabel('Y1');title('X1与X2相乘的结果');Y2=ifft(Y1);subplot(7,1,7);stem(Y2);ylabel('Y2');title('Y1的IDFT结果');运行结果如图3.2：图3.2 FFT快速卷积3.3圆周卷积x1(n)=2,0,1,1,5,7,0,5,0,2,1,2，x2(n)= 2,4,5,6,2,7,4,8,9,6,4,8,3,5,7,8，其圆周卷积Matlab程序如下：clf;x=2,0,1,1,5,7,0,5,0,

14、2,1,2;y=2,4,5,6,2,7,4,8,9,6,4,8,3,5,7,8;N1=length(x);N2=length(y);N=N1;if N2>N1 N=N2;endM=input('请输入圆周卷积长度n');if M>=N N=M;endif N=N1 x=x;else x=x,zeros(1,N-N1);endX=zeros(N,N);x1=zeros(1,N);for n=1:N if n=1 x1(n)=x(n); else x1(n)=x(N-n+2); endendX(1,1:N)=x1;for i=2:N for n=1:N if n=1 X

15、(i,n)=X(i-1,N); else X(i,n)=X(i-1,n-1); end endendif N2<Ny=y,zeros(1,N-N2);endY=y'W=X*Y;W=W'n=1:N;stem(n,W);当圆周卷积长度为89、30时，结果如图3.3、3.4：图3.3 圆周卷积（N=89）图3.4 圆周卷积（N=30）3.4采样定理演示连续信号：，A为幅度因子，a为衰减因子，为模拟角频率，其中n为16，其采样Matlab程序如下：fs=input('请输入采样频率:');n=input('请输入两位学号:');t=0:1/fs:1

16、;t1=0:0.0005:1;A=n;a=sqrt(2)*n*pi;w=a;x=A*exp(-a*t).*sin(w*t);x1=A*exp(-a*t1).*sin(w*t1);subplot(5,1,1)plot(x1);xn=ones(1,length(t);subplot(512)stem(xn)yn=xn.*x;subplot(5,1,3)stem(yn);y=fft(yn,fs);f=fs*(0:fs/2)/fs;subplot(5,1,4)plot(f,y(1:fs/2+1);yn=ifft(y);subplot(5,1,5)plot(yn)当fs= 50 时运行结果如图3.5：图

17、3.5 采样定理演示（fs=50）当fs= 20 时运行结果如图3.6：图3.6采样定理演示（fs=20）四 课程设计原理4.1 巴特沃斯低通滤波器（1）特性:巴特沃斯滤波器的特点是通频带内的频率响应曲线最大限度平坦，没有起伏，而在阻频带则逐渐下降为零。 在振幅的对数对角频率的波特图上,从某一边界角频率开始,振幅随着角频率的增加而逐步减少,趋向负无穷大。一阶巴特沃斯滤波器的衰减率为每倍频6分贝，每十倍频20分贝。二阶巴特沃斯滤波器的衰减率为每倍频12分贝、三阶巴特沃斯滤波器的衰减率为每倍频18分贝、如此类推。巴特沃斯滤波器的振幅对角频率 单调下降，并且也是唯一的无论阶数，振幅对角频率曲线都保持

18、同样的形状的滤波器。只不过滤波器阶数越高，在阻频带振幅衰减速度越快。其他滤波器高阶的振幅对角频率图和低级数的振幅对角频率有不同的形状。(2)传递函数：巴特沃斯低通滤波器可用如下振幅的平方对频率的公式表示：其中,= 滤波器的阶数= 截止频率 = 振幅下降为 -3分贝时的 频率= 通频带边缘频率在通频带边缘的数值 4.2 用窗函数法设计FIR滤波器根据过渡带宽及阻带衰减要求，选择窗函数的类型并估计窗口长度N（或阶数M=N-1），窗函数类型可根据最小阻带衰减As独立选择，因为窗口长度N对最小阻带衰减As没有影响，在确定窗函数类型以后，可根据过渡带宽小于给定指标确定所拟用的窗函数的窗口长度N

19、，设待求滤波器的过渡带宽为w，它与窗口长度N近似成反比，窗函数类型确定后，其计算公式也确定了，不过这些公式是近似的，得出的窗口长度还要在计算中逐步修正，原则是在保证阻带衰减满足要求的情况下，尽量选择较小的N，在N和窗函数类型确定后，即可调用MATLAB中的窗函数求出窗函数wd（n）。根据待求滤波器的理想频率响应求出理想单位脉冲响应hd(n)，如果给出待求滤波器频率应为Hd，则理想的单位脉冲响应可以用下面的傅里叶反变换式求出：在一般情况下，hd(n)是不能用封闭公式表示的，需要采用数值方法表示；从w=0到w=2采样N点，采用离散傅里叶反变换（IDFT）即可求出。用窗函数wd(n)将hd(n)截断

20、，并进行加权处理，得到如果要求线性相位特性， 则h(n)还必须满足：根据上式中的正、 负号和长度N的奇偶性又将线性相位FIR滤波器分成四类。 要根据所设计的滤波特性正确选择其中一类。 例如， 要设计线性相位低通特性可选择h(n)=h(N-1-n)一类，而不能选h(n)=-h(N-1-n)一类。 验算技术指标是否满足要求，为了计算数字滤波器在频域中的特性，可调用freqz子程序，如果不满足要求，可根据具体情况，调整窗函数类型或长度，直到满足要求为止。五 设计方案及步骤5.1 设计方案1.要求利用windows下的录音机（开始程序附件娱乐录音机，文件属性立即转换8000Hz，8位，单声道）录制一段

21、自己的话音，或者采用Windows自带的声音文件（默认为22050Hz），时间控制在几秒左右。然后在MATLAB软件平台下，利用函数wavread对语音信号进行采样，记住采样频率和采样点数。通过wavread函数的使用，要求理解采样频率、采样位数等概念。wavread函数调用格式：y=wavread(file)，读取file所规定的wav文件，返回采样值放在向量y中。y,fs,nbits=wavread(file)，采样值放在向量y中，fs表示采样频率（Hz），nbits表示采样位数。y=wavread(file,N)，读取前N点的采样值放在向量y中。y=wavread(file,N1,N2)

22、，读取从N1点到N2点的采样值放在向量y中。2语音信号的频谱分析要求首先画出语音信号的时域波形；然后对语音信号进行频谱分析，在MATLAB中，可以利用函数fft对信号进行快速付立叶变换，得到信号的频谱特性；从而加深对频谱特性的理解。3设计数字滤波器和画出频率响应根据语音信号的特点给出有关滤波器的性能指标，例如：1）低通滤波器性能指标，fp=1000Hz，fc=1200 Hz， As=100dB，Ap=1dB；2）高通滤波器性能指标，fc=2800 Hz，fp=3000 Hz，As=100dB，Ap=1dB；3）带通滤波器性能指标，fp1=1200 Hz，fp2=3000 Hz，fc1=1000

23、 Hz，fc2=3200 Hz，As=100dB，Ap=1dB。4.回放语音信号在Matlab中,函数sound可以对声音进行回放。其调用格式：sound(x,fs,bits);可以感觉滤波前后的声音有变化。5.方案设计流程图如图3.5.1所示。（图为巴特沃思滤波器，hamming窗滤波器与巴特沃思类似）图5.1 方案设计流程图5.2设计步骤（1）语音信号的采集Matlab程序如下:运行结果如图3.5.2x,FS,bits=wavread('C:WINDOWSMediading.wav');x=x(:,1);figure(1);subplot(2,1,1);plot(x);so

24、und(x,FS,bits); title('语音信号时域波形图')y=fft(x,3260);f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title('语音信号频谱图');运行结果如图3.5.2。图5.2 语音信号（2）对语音信号加噪Matlab程序如下：t=0:length(x)-1;zs=0.05*cos(2*pi*10000*t/22050);zs0=0.05*cos(2*pi*10000*t/22050000);figure(2);subplot(2,1,1)plot(

25、zs0)title('噪声信号波形');zs1=fft(zs,1200);sound(zs,FS,bits); subplot(2,1,2)plot(f(1:600),abs(zs1(1:600);title('噪声信号频谱');x1=x+zs'sound(x1,FS,bits); y1=fft(x1,1200);figure(3);subplot(2,1,1);plot(x1);title('加入噪声后的信号波形');subplot(2,1,2);plot(f(1:600),abs(y1(1:600);title('加入噪声后的

26、信号频谱');运行结果如图3.5.3(a)(b)图5.3 对语音信号加噪六 用巴特沃斯设计IIR滤波器6.1巴特沃斯低通滤波器Matlab程序见附录A巴特沃斯低通滤波器程序所示，运行程序波形图分析如下：图6.1 6.2. 巴特沃斯带通滤波器Matlab程序见附录B切比雪夫带通滤波器程序所示，运行程序波形图分析如下：图6.26.3 巴特沃斯高通滤波器Matlab程序见附录C比雪夫高通滤波器程序所示，运行程序波形图分析如下：图6.3七 hamming窗设计FIR滤波器7.1 hamming高通滤波器调试程序见附录Dhamming窗高通滤波器程序所示，运行程序波形图分析如下：图7.1 7.2

27、 hamming低通滤波器调试程序见附录hamming窗低通滤波器程序所示，运行程序波形图分析如下：图7.2 7.3 Hamming带通滤波器调试程序见附录hamming窗带通滤波器程序所示，运行程序波形图分析如下：图7.3五 心得体会课程设计是我们专业课程知识综合应用的实践训练，着是我们迈向社会，从事职业工作前一个必不少的过程”千里之行始于足下”，通过这次课程设计，我深深体会到这句千古名言的真正含义我今天认真的进行课程设计，学会脚踏实地迈开这一步，就是为明天能稳健地在社会大潮中奔跑打下坚实的基础在这次的课程设计中，犯了些错误，发现自己还不够认真不够细心以前的知识学习的不是很牢固，因为数字信号

28、处理学的不是很好，很多相关知识的运用不是很熟悉，还有MATLAB 的很多相关指令的运用不熟练，导致课程设计做起来有点吃力，但是只要努力去认真做就会得到你想要的，通过这次DSP课程设计，不但巩固了以前学的知识，还学习了新内容，收获颇丰。原来滤去语音中的噪声我们运用自己学到的知识能做到，并更加巩固了MATLAB软件的使用，在从这次的课程设计中体会到，我觉得什么事都要我们真真正正用心去做了，才会使自己更加成长，没有学习就不可能有实践的机会，没有自己的实践就不会有所突破，希望这次的经历能让我们在以后的学习生活中不断成长，在以后的学习中更能与实际相结合，为自己今后在社会上能立足打好坚实的基础。想起大二的

29、实习和这次的课程设计，至今我仍感慨颇多，从理论到实践，在这段日子里，可以说得是苦多于甜，但是可以学到很多很多的东西，同时不仅可以巩固了以前所学过的知识，而且学到了很多在书本上所没有学到过的知识。通过这次课程设计使我懂得了理论与实际相结合是很重要的，只有理论知识是远远不够的，只有把所学的理论知识与实践相结合起来，从理论中得出结论，才能真正为社会服务，从而提高自己的实际动手能力和独立思考的能力。在设计的过程中遇到问题，可以说得是困难重重，但可喜的是最终都得到了解决。最后在此感谢我们的黄亚飞老师.，老师严谨细致的作风一直是我学习中的榜样。这次课程设计的每个实验细节和每个数据，都离不开老师您的细心指导

30、。而您开朗的个性和宽容的态度，帮助我能够很顺利的完成了这次课程设计。同时感谢对我帮助过的同学们，谢谢你们对我的帮助和支持，让我感受到同学的友谊参考文献（1）数字信号处理（第二版），丁玉美等，西安电子科技大学出版社；（2）数字信号处理试验指导书王创新、文卉编 长沙理工大学印刷（内部使用）（3）数字信号处理及其MATLAB实现，陈怀琛等译，电子工业出版社；（4）MATLAB及在电子信息课程中的应用，陈怀琛等，电子工业出版社。（5）数字信号处理A.V.奥本海姆，R.W.谢弗著，北京：科学出版社（6）数字信号处理理论、算法与实现（第二版）胡广书编著，北京：电子工业出版社附录A：巴特沃斯低通滤波器程序F

31、s=22050;x,FS,bits=wavread('C:WINDOWSMediaWindows XP 关机.wav');x=x(:,1); figure(1); subplot(2,1,1); plot(x); sound(x,FS,bits); %回放语音 title('语音信号时域波形图')y=fft(x,3260); f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title('语音信号频谱图'); %产生噪声信号并加到语音信号 t=0:length(x)-1

32、; zs=0.002*cos(2*pi*10000*t/22050); zs0=0.002*cos(2*pi*10000*t/22050000); figure(2); subplot(2,1,1); plot(zs0) title('噪声信号波形');zs1=fft(zs,1200); sound(zs,FS,bits); %回放噪音 subplot(2,1,2);plot(f(1:600),abs(zs1(1:600);title('噪声信号频谱'); x1=x+zs' sound(x1,FS,bits); %回放加入噪声后的语音 y1=fft(x1

33、,1200); figure(3); subplot(2,1,1);plot(x1); title('加入噪声后的信号波形'); subplot(2,1,2); plot(f(1:600),abs(y1(1:600); title('加入噪声后的信号频谱'); %低通滤波 fp=3000;fs=3500;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'

34、); B,A=butter(N,wc,'s'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯低通滤波器');xlabel('频率HZ');ylabel('耗损dB'); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); title('滤波后信号波形'); ydd=fft(yd,1200);

35、subplot(2,1,2);plot(f(1:600),abs(ydd(1:600); title('滤波后信号频谱'); sound(yd,FS,bits);附录B：巴特沃斯高通滤波器程序Fs=22050; x,FS,bits=wavread('C:WINDOWSMediaWindows XP 关机.wav'); x=x(:,1); figure(1); subplot(2,1,1); plot(x);sound(x,FS,bits); %回放语音 title('语音信号时域波形图') y=fft(x,3260); f=(FS/1630)*1

36、:1630; subplot(2,1,2); plot(f(1:1630),abs(y(1:1630); title('语音信号频谱图'); %产生噪声信号并加到语音信号 t=0:length(x)-1; zs0=0.002*cos(2*pi*100*t/22050); figure(2); subplot(2,1,1); plot(zs0) title('噪声信号波形'); zs1=fft(zs0,1200); sound(zs,FS,bits); %回放噪音 subplot(2,1,2); plot(f(1:600),abs(zs1(1:600); titl

37、e('噪声信号频谱'); x1=x+zs0' sound(x1,FS,bits); %回放加入噪声后的语音 y1=fft(x1,1200); figure(3); subplot(2,1,1);plot(x1); title('加入噪声后的信号波形'); subplot(2,1,2); plot(f(1:600),abs(y1(1:600); title('加入噪声后的信号频谱'); %高通滤波 fp=600;fs=400;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs

38、1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'high','s'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯高通滤波器'); xlabel('频率HZ');ylabel('耗损dB'); grid on; yd=filte

39、r(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); title('滤波后信号波形'); ydd=fft(yd,1200); subplot(2,1,2);plot(f(1:600),abs(ydd(1:600); title('滤波后信号频谱'); sound(yd,FS,bits);附录C：巴特沃斯带通滤波器程序Fs=22050; x,FS,bits=wavread('C:WINDOWSMediaWindows XP 关机.wav'); x=x(:,1); figure(1); subplot(2,

40、1,1); plot(x); sound(x,FS,bits); %回放语音 title('语音信号时域波形图'); y=fft(x,3260); f=(FS/1630)*1:1630; subplot(2,1,2); plot(f(1:1630),abs(y(1:1630); title('语音信号频谱图');%产生噪声信号并加到语音信号 t=0:length(x)-1; zs0=0.002*cos(2*pi*100*t/22050); figure(2); subplot(2,1,1);plot(zs0) title('噪声信号波形'); z

41、s1=fft(zs0,1200); sound(zs,FS,bits); %回放噪音 subplot(2,1,2);plot(f(1:600),abs(zs1(1:600); title('噪声信号频谱'); x1=x+zs0' sound(x1,FS,bits); %回放加入噪声后的语音 y1=fft(x1,1200); figure(3); subplot(2,1,1);plot(x1); title('加入噪声后的信号波形'); subplot(2,1,2); plot(f(1:600),abs(y1(1:600); title('加入噪声

42、后的信号频谱'); %带通滤波 fp=600,6000;fs=400,7000;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'s'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title

43、('巴特沃斯带通滤波器'); xlabel('频率HZ');ylabel('耗损dB'); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); ydd=fft(yd,1200); subplot(2,1,2);plot(f(1:600),abs(ydd(1:600);sound(yd,FS,bits);附录D：hamming高通滤波器程序x,FS,bits=wavread('C:WINDOWSMediading.wav');x=x(:,1);figur

44、e(1);subplot(2,1,1);plot(x);sound(x,FS,bits)title('语音信号时域波形图')y=fft(x,3260);f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title('语音信号频谱图'); t=0:length(x)-1;zs0=0.05*cos(2*pi*100*t/22050);figure(2);subplot(2,1,1)plot(zs0)title('噪声信号波形');zs1=fft(zs0,1200); sub

45、plot(2,1,2)plot(f(1:600),abs(zs1(1:600);title('噪声信号频谱');x1=x+zs0'sound(x1,FS,bits)y1=fft(x1,1200);figure(3);subplot(2,1,1);plot(x1);title('加入噪声后的信号波形');subplot(2,1,2);plot(f(1:600),abs(y1(1:600);title('加入噪声后的信号频谱'); fp=600,fc=400; wp=2*pi*fp/FS;ws=2*pi*fc/FS;Bt=wp-ws; N0=

46、ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wp+ws)/2/pi; hn=fir1(N-1,wc,'high',hamming(N); X=conv(hn,x); sound(X,FS,bits);X1=fft(X,1200);figure(4);subplot(211);plot(X);title('滤波后的信号波形');subplot(212);plot(f(1:600),abs(X1(1:600); title('滤波后的信号频谱')附录E：hamming低通滤波器程序x,FS,bits=wavread(&#

47、39;C:WINDOWSMediading.wav');x=x(:,1);figure(1);subplot(2,1,1);plot(x);sound(x,FS,bits)title('语音信号时域波形图')y=fft(x,3260);f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title('语音信号频谱图'); t=0:length(x)-1;zs=0.05*cos(2*pi*10000*t/22050);zs0=0.05*cos(2*pi*10000*t/22050000);figure(2);subplot(2,1,1)plot(zs0)title(