基于语音信号去噪处理的FIR低通滤波器设计

1、摘要本次课程设计分析了FIR数字滤波器的基本原理，在MATLAB环境下利用窗函数设计FIR低通滤波器，实现了FIR低通滤波器的设计仿真。本文根据滤波后的时域图和原始语音信号时域图的比较，以及滤波后信号的频谱图和原始语音信号频谱图的比较，最后回放滤波后语音信号，滤波后的语音信号与原始语音信号一样清晰，仿真结果表明，设计的FIR滤波器的各项性能指标均达到了指定要求，设计过程简便易行。该方法为快速、高效地设计FIR滤波器提供了一个可靠而有效的途径。关键词：DSP ；FIR；低通滤波器；语音信号；MATLAB 目录第一章 引言11.1 设计目的及意义11.2 设计任务及要求21.3 课程设计平台2第二

2、章 基本原理32.1 FIR滤波器的基本概念32.2 FIR滤波器的特点32.3 FIR滤波器的种类4第三章 FIR数字低通滤波器的设计53.1 FIR低通滤波器设计原理53.2 FIR低通滤波器的设计方法53.2.1 频率采样法53.2.2 最优化设计63.2.3 窗函数法63.3 窗函数法设计步骤8第四章 详细设计94.1 语音信号的采集94.2 语音信号的读入与打开104.3 语音信号的FFT变换114.4 含噪信号的合成124.5 利用FIR滤波器滤波134.6 结果分析14总结15参考文献16附录17致谢210第一章 引言随着信息科学和计算机技术的不断发展，数字信号处理(DSP，Di

3、gital Signal Processing)的理论和技术也得到了飞速的发展，并逐渐成为一门重要的学科，它的重要性在日常通信、图像处理、遥感、声纳、生物医学、地震、消费电子、国防军事、医疗方面等显得尤为突出。在我们面临的信息革命中，数字信号处理几乎涉及了所有的工程技术领域。 数字信号处理是一种将信号以数字形式进行处理的一种理论和技术，它的目的是将真实世界中的一些信号进行分析并滤波，最后得出其中的有用的信号。数字滤波器是数字信号处理的一种，一般根据单位脉冲响应h(n)分为无限脉冲响应(IIR)和有限脉冲响应(FIR)系统。IIR数字滤波器的设计方法简单，特别是采用

4、双线性变换法来设计的数字滤波器不存在频域混叠的现象，但是IIR滤波器存在一个较为明显的缺憾，就是它的相位响应一般都是非线性的，而在传输频带内的相位响应如果不是线性的，就会造成有用信号的传输失真，而FIR数字滤波器不仅可以设计成任意的幅度响应，而且可以设计成在通频带内具有良好的线性相位响应。FIR数字滤波器的单位脉冲响应h(n)有限长，所以FIR数字滤波器是稳定的，不存在稳定性的问题，且可以通过快速傅里叶变换(FFT)的算法来实现信号滤波，大大的提高的运算效率。1因此，FIR数字滤波器日益引起了人们的关注。本课程设计是采用kaiser窗设计的FIR滤波器对语音信号进行滤波去噪。通过课程设计了解F

5、IR滤波器设计的原理和步骤，掌握用Matlab语言设计滤波器的方法，了解DSP对FIR滤波器的设计及编程方法。通过观察语音信号滤波前后的时域波形的比较，加深对滤波器作用的理解。通过对比滤波前后波形图的比较和放滤波前后语音信号的对比 ，可以看出滤波器对有用信号无失真放大具有重大意义。1.1 设计目的及意义信号处理课程设计是现代信号处理技术课程的有效补充部分，通过课程设计，使得学生在设计过程中了解完整的现代信号处理技术的工程实现方法和流程，从而对现代信号处理技术的理论有更深入的认识。本课程设计的目的是通过学生使用MATLAB等工具，采用窗函数法设计符合一定参数要求的FIR滤波器，并用所设计的滤波器

6、对加噪语音信号进行滤波去噪处理。1.2 设计任务及要求利用MATLAB语言及其工具箱来完成如下的工作： 1、学会MATLAB的使用，掌握MATLAB的程序设计方法； 2、录制一段语音信号，完成对信号的采样，画出信号的时域波形和频谱图，确定信号的频谱范围； 3、给信号叠加噪声（噪声类型分为如下几种：a白噪声；b单频噪色（正弦干扰）；c多频噪声（多正弦干扰）；d其它干扰。），画出受噪声干扰的信号时域波形和频谱图； 4、采用窗函数法设计FIR低通滤波器，低通滤波器性能指标 fb1000 Hz，fc1200 Hz，as100 dB，ap1 dB。画出滤波器的频响特性图； 5、用所设计的滤波器对受噪声影

7、响的信号进行滤波，画出滤波后语音信号的时域波形图和频谱图； 6、对滤波前后的信号进行对比，分析信号的变化；回放语音信号，并与原始语音信号对比；7、掌握窗函数法设计FIR数字滤波器的方法。1.3 课程设计平台20世纪70年代后期，时任美国新墨西哥大学计算机科学系主任的Cleve Moler教授出于减轻学生编程负担的动机，为学生设计了一组调用LINPACK和EISPACK库程序的“通俗易用”的接口，此即用FORTRAN编写的萌芽状态的MATLAB。MATLAB软件包括五大通用功能：数值计算功能（Nemeric）；符号运算功能（Symbolic）；数据可视化功能（Graphic）；数据图形文字统一处

8、理功能（Notebook）和建模仿真可视化功能（Simulink）。其中，符号运算功能的实现是通过请求MAPLE内核计算并将结果返回到MATLAB命令窗口。该软件有三大特点：一是功能强大；二是界面友善、语言自然；三是开放性强。目前，Mathworks公司已推出30多个应用工具箱。MATLAB在线性代数、矩阵分析、数理统计和随机信号分析、电路与系统、系统动力学、信号和图像处理、建模和仿真、通信系统、以及财政金融等众多领域的理论研究和工程设计中得到了广泛应用。第二章 基本原理2.1 FIR滤波器的基本概念数字滤波器从实现的网络结构或者从单位脉冲响应可分为无限脉冲响应（IIR）滤波器和有限脉冲响应（

9、FIR）滤波器。FIR滤波器是有限长单位冲激响应滤波器，是数字信号处理系统中最基本的元件，它可以在保证任意幅频特性的同时具有严格的线性相频特性，同时其单位抽样响应是有限长的，因而滤波器是稳定的系统。2因此，FIR滤波器在通信、图像处理、模式识别等领域都有着广泛的应用。从性能上说，IIR滤波器以非线性相位为代价以较低的阶数获得较高的选择性。而FIR滤波器想要获得相同的选择性阶数是IIR滤波器的5-10倍，结果成本较高、信号时延也较大：从结构上说，IIR采用递归结构，FIR采用非递归结构；从设计工具上说；IIR可以借助于模拟滤波器的成果，FIR滤波器一般采用没有封闭形式的设计公式；从使用场合上来看

10、，在对相位要求不敏感的场合，如语音通讯等，选用IIR较为合适，可以充分发挥经济高效的特点。对图像处理、数据传输等以波形携带信息的系统，使用FIR较好。2.2 FIR滤波器的特点有限长单位冲激响应（FIR）滤波器有以下特点： （1) 系统的单位冲激响应在有限个n值处不为零； （2) 系统函数H（z）在|z|>0处收敛，极点全部在z = 0处（因果系统）； （3） 结构上主要是非递归结构，没有输出到输入的反馈，但有些结构中（例如频率抽样结构）也包含有反馈的递归部分。 设FIR滤波器的单位冲激响应为一个N点序列，则滤波器的系统函数为 式（2-1）就是说，它有（N1）阶极点在z = 0处，有（N

11、1）个零点位于有限z平面的任何位置。优点 ： （1）很容易获得严格的线性相位，避免被处理的信号产生相位失真，这一特点在宽频带信号处理、阵列信号处理、数据传输等系统中非常重要； （2）可得到多带幅频特性； （3）极点全部在原点（永远稳定），无稳定性问题； （4）任何一个非因果的有限长序列，总可以通过一定的延时，转变为因果序列，所以因果性总是满足； （5）无反馈运算，运算误差小。缺点： （6）因为无极点，要获得好的过渡带特性，需以较高的阶数为代价； （7）无法利用模拟滤波器的设计结果，一般无解析设计公式，要借助计算机辅助设计程序完成。32.3 FIR滤波器的种类（1）数字集成电路FIR滤波器4一种

12、是使用单片通用数字滤波器集成电路，这种电路使用简单，但是由于字长和阶数的规格较少，不易完全满足实际需要。虽然可采用多片扩展来满足要求，但会增加体积和功耗，因而在实际应用中受到限制。设计数字滤波器的任务就是寻求一个因果稳定的线性时不变系统，使其系统函数H(z)具有指定的频率特性。（2）DSP芯片FIR滤波器另一种是使用DSP芯片。DSP芯片有专用的数字信号处理函数可调用，实现FIR滤波器相对简单，但是由于程序顺序执行，速度受到限制。而且，就是同一公司的不同系统的DSP芯片，其编程指令也会有所不同，开发周期较长。（3）可编程FIR滤波器还有一种是使用可编程逻辑器件，FPGACPLD。FPGA有着规

13、整的内部逻辑块整列和丰富的连线资源，特别适合用于细粒度和高并行度结构的FIR滤波器的实现，相对于串行运算主导的通用DSP芯片来说，并行性和可扩展性都更好。第三章 FIR数字低通滤波器的设计3.1 FIR低通滤波器设计原理由于FIR 滤波器冲击响应hn是有限长序列，因此这种结构可用非递归结构来实现1。FIR 数字滤波器系统函数一般形式为如式（3-1）所示： 式（3-1）FIR 滤波器数学表达式可用差分方程（3-2）来表示： 式（3-2）式中:y（n）输出序列；h（k）滤波器系数；n滤波器阶数；x（k）输入序列。应用Matlab 设计FIR滤波器的主要任务就是根据给定的性能指标，设计一个H（z）

14、,使其逼近这一指标，进而计算并确定滤波器的系数b(n)，再将所设计滤波器的幅频响应、相频响应曲线作为输出,与设计要求进行比较，对设计的滤波器进行优化。3.2 FIR低通滤波器的设计方法3.2.1 频率采样法频率采样法是从频域出发,根据频域采样定理，对给定的理想滤波器的频率响应加以等间隔的抽样 ，得到： ，k=0,1,N-1 式（3-3）再利用可求得FIR滤波器的系统函数及频率响应。 而在各采样点间的频率响应则是其加权内插函数延伸叠加的结果。但对于一个无限长的序列，用频率采样法必然有一定的逼近误差,误差的大小取决于理想频响曲线的形状,理想频响特性变换越平缓,则内插函数值越接近理想值,误差越小。5

15、为了提高逼近的质量，可以通过在频率相应的过渡带内插入比较连续的采样点，扩展过渡带使其比较连续，从而使得通带和阻带之间变换比较缓慢，以达到减少逼近误差的目的。 选取w0,2内N个采样点的约束条件为如式（34）： 0kN-1 式（3-4） 增大阻带衰减三种方法：1）加宽过渡带宽，以牺牲过渡带换取阻带衰减的增加。2）过渡带的优化设计利用线性最优化的方法确定过渡带采样点的值，得到要求的滤波器的最佳逼近（而不是盲目地设定一个过渡带值）。3）增大N。如果要进一步增加阻带衰减，但又不增加过渡带宽，可增加采样点数N。代价是滤波器阶数增加，运算量增加。直接从频域进行设计，物理概念清楚，直观方便；适合于窄带滤波器

16、设计，这时频率响应只有少数几个非零值，但是截止频率难以控制。3.2.2 最优化设计最优化设计方法是指采用最优化准则来设计的方法。在 FIR DF的最优化设计中 ,最优化准则有均方误差最小化准则和等波纹切比雪夫逼近(也称最大误差最小化)准则两种。但是对于本次的课程设计来说存在一定的局限性和难度。3.2.3 窗函数法设计FIR数字滤波器的最简单的方法是窗函数法，通常也称之为傅立叶级数法。FIR数字滤波器的设计首先给出要求的理想滤波器的频率响应,设计一个FIR数字滤波器频率响应，去逼近理想的滤波响应。窗函数法设计FIR数字滤波器是在时域进行的，因而必须由理想的频率响应推导出对应的单位取样响应，再设计

17、一个FIR数字滤波器的单位取样响应去逼近。设计过程如图3.1： 图3.1加窗的作用是通过把理想滤波器的无限长脉冲响应乘以窗函数来产生一个被截断的脉冲响应，即并且对频率响应进行平滑。MATLAB工具箱提供的窗函数有：矩形窗(Rectangularwindow)、三角窗(Triangular window)、布拉克曼窗(Blackman window)、汉宁窗(Hanningwindow)、海明窗(Hamming window)、凯塞窗(Kaiser window)、切比雪夫窗(Chebyshev window)。窗函数主要用来减少序列因截断而产生的Gibbs效应。但当这个窗函数为矩形时，得到的F

18、IR滤波器幅频响应会有明显的Gibbs效应，并且任意增加窗函数的长度（即FIR滤波器的抽头数）Gibbs效应也不能得到改善。为了克服这种现象,窗函数应该使设计的滤波器：（1） 频率特性的主瓣宽度应尽量窄,且尽可能将能量集中在主瓣内；（2） 窗函数频率特性的旁瓣趋于 的过程中，其能量迅速减小为零。窗函数法设计FIR滤波器由于其运算简单、精度高，已成为工程中应用最广泛的设计方法。因此，本次设计采用窗函数法设计滤波器。其基本思想是让待设计的实际滤波器逼近理想特性。理想低通滤波器的频率特性应为：振幅特性在通带内为1，阻带内为0；在通带内的相位特性与w成线性关系，如式（35）： 式（35）它对应的理想单

19、位冲激相应为式（36）: 式（36)可见，理想低通滤波器的单位冲激响应hd (n)是无限长的非因果序列。而我们要设计的是FIR滤波器，其h( n)必然是有限长的，为了构造物理上可以实现的长度为N 的因果线性相位滤波器，必须将hd( n) 截取长度为N 的一段，或者说用一个有限长度的窗口函数序列w ( n) 来截取hd ( n)，即h( n) = w ( n) hd ( n)，因而窗函数序列的形状及长度的选择就是关键。工程实际中常用的窗函数有6种，即矩形窗、三角形窗、汉宁(Hanning)窗、哈明(Hamming)窗、布莱克曼(Blackman) 窗和凯泽( Kaiser) 窗。它们之间的性能比

20、较如表2所示。实际的FIR数字滤波器的频率响应H() 逼近理想滤波器频率响应Hd()的好坏，取决于窗函数的频率特性W( w)。6各种窗函数的特性如表1所示： 表3-1 6种窗函数的特性窗函数旁瓣峰值/dB近似过渡带宽精确过渡带宽阻带最小衰减/dB矩形窗-134/N1.8/N21三角形窗-258/N6.1/N25汉宁窗-318/N6.2/N44哈明窗-418/N6.6/N53布莱克曼窗-5712/N11/N74凯泽窗（）-5710/N803.3 窗函数法设计步骤用窗函数设计FIR滤波器的步骤为: (1)选择窗函数类型能够和长度，写出窗函数w（n）表达式。根据阻带最小衰减选择窗函数w(n)的类型，

21、再根据过渡带宽确定所选窗函数的长度N。(2)构造希望逼近的频率响应函数。选择理想低通滤波器，应当注意，理想滤波器的截止频率wc点近似为最终设计的FIRDFD的-6dB频率。(3)根据式（37)计算： 式（37)（4）加窗得到设计结果式（38)： 式（38)21第四章 详细设计本次课程设计先完成语音信号的采集，并对所采集的语音信号加入不同的干扰噪声，对加入噪声的信号进行频谱分析，针对受干扰语音信号的特点设计不同的滤波器，然后利用窗函数法设计低通，高通，带通等滤波器对采集到的语音信号进行滤波处理，分析语音信号各频率段的特性。对加噪信号进行滤波，恢复原信号。把原始语音信号、加噪语音信号和滤波后的信号

22、进行时域变换和频域变换，画出它们的时域波形和频域波形图，从视觉角度比较分析滤波的效果。实现框图如图4-1所示：信号的采集 语音信号分析 结果显示及分析 含噪信号合成滤 波 数字滤波器的设计 图4-1 整体设计流程图4.1 语音信号的采集 要求录制一段语音信号，并对录制的信号进行采样。在MATLAB中可以使用wavrecord（）函数完成，例如：y=wavrecord(13*fs,fs,);可以用来录制13秒的语音，wavplay（y,fs）; 用来播放录制的语音，wavwrite(y,fs,'111.wav')；wavwrite(y,'111.wav');两句则

23、是将录制的语音保存下来。74.2 语音信号的读入与打开 在MATLAB中，y,fs,bits=wavread('111.wav');用于读取语音，采样值放在向量y中，fs表示采样频率(Hz)，bits表示采样位数。wavplay（y,fs）；用于语音的播放，也可以使用sound(y)； 用于对声音的回放。向量y则就代表了一个信号，也即一个复杂的“函数表达式”，也可以说像处理一个信号的表达式一样处理这个声音信号。在MATLAB中，实现语音的读入，并绘出了语音信号的原始信号图可采用fs=11025; 语音信号的采样频率为11025，x1=wavread('aa.wav

24、9;); 读取语音信号的数据赋给变量，sound(x1,11025); 用于播放语音信号。语音信号是一种非平稳的时变信号，它携带着各种信息。在语音编码、语音合成、语音识别和语音增强等语音处理中无一例外需要提取语音中包含的各种信息。语音信号分析的目的就在与方便有效的提取并表示语音信号所携带的信息。语音信号分析可以分为时域和变换域等处理方法，其中时域分析是最简单的方法，直接对语音信号的时域波形进行分析，提取的特征参数主要有语音的短时能量，短时平均过零率，短时自相关函数等。原始语音信号图如图4-2所示： 图4-2 原始语音信号图如图中提取的语音的波形图所示，整段音频数据中得声音高低起伏与录入的声音信

25、号基本一致，从图中可以看出语音信号的截止频率为4000Hz。4.3 语音信号的FFT变换FFT即为快速傅氏变换，是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。在MATLAB的信号处理工具箱中函数FFT和IFFT用于快速傅立叶变换和逆变换。函数FFT用于序列快速傅立叶变换，其调用格式为y=fft(x)，其中，x是序列，y是序列的FFT，x可以为一向量或矩阵，若x为一向量，y是x的FFT且和x相同长度；若x为一矩阵，则y是对矩阵的每一列向量进行FFT。8如果x长度是2的幂次方，函数fft执行高速基2FFT算法，否则fft执行一种混合基的

26、离散傅立叶变换算法，计算速度较慢。函数FFT的另一种调用格式为y=fft(x,N)，式中，x，y意义同前，N为正整数。函数执行N点的FFT，若x为向量且长度小于N，则函数将x补零至长度N；若向量x的长度大于N，则函数截短x使之长度为N；若x 为矩阵，按相同方法对x进行处理。对语音信号进行FFT变换就是为了得到它的频域的图形，便于从图中观察出信号的幅度等特性。如图4-3所示即为原始信号时域波形和频谱图： 图4-3 为原始信号时域波形和频谱图 信号的傅立叶表示在信号的分析与处理中起着重要的作用。因为对于线性系统来说，可以很方便地确定其对正弦或复指数和的响应，所以傅立叶分析方法能完善地解决许多信号分

27、析和处理问题。语音信号经过FFT变换之后，就可以得到信号的频域图形，另外，傅立叶表示使信号的某些特性变得更明显，因此，它能更深入地说明信号的各项宏物理现象。4.4 含噪信号的合成在MATLAB软件平台下，给原始的语音信号叠加上噪声，噪声类型分为如下几种：白噪声、单频噪声（正弦干扰）、多频噪声（多正弦干扰）、其他干扰，可设置为低频、高频、带限噪声、Chirp干扰、冲激干扰。绘制出叠加噪声之后的语音信号时域图形及频域图形，在视觉上与原始信号图形对比，也可以通过Windows播放软件或者是利用MATLAB中的sound()函数从听觉上进行对比，分析并体会含噪语音信号频谱和时域波形的改变。加噪后语音信

28、号时域及频域图如图4-4所示： 图4-4 加噪后信号的时域及频域图与原始信号对比，加入噪声后的波形与原始信号相比出现了一些区别：先原始信号没加噪音之前0到2000有幅值，在4000到6000之间没有幅值，但是在加了噪音之后4000到6000之间出现超过最大幅值15，超出了正常值。4.5 利用FIR滤波器滤波利用凯瑟窗设计低通FIR滤波器，利用凯瑟窗设计高通FIR滤波器对合成的含噪信号进行滤波，绘制出FIR滤波器的频率响应，绘出滤波后的时域波形和频谱图，并对滤波前后的信号进行对比，分析信号的变化；在相同的性能指标下比较个方法的滤波效果，并从理论上进行分析。语音信号中包含背景噪声，这些噪声的频率一

29、般较高。所以可以利用MATLAB软件中设计的FIR滤波器进行滤波处理，得到较为理想的语音信号。从以上各种数字滤波器经过滤波后得出的语音信号相比较，使用低通滤波器滤波以后，声音稍微有些发闷，但是很接近原来录制的声音；使用高通滤波器滤波后则听不到人的声音。（1） 用凯瑟窗设计FIR低通滤波器利用凯瑟窗设计的FIR低通滤波器，滤波器的性能指标为： fb1000 Hz，fc1200 Hz，as100 dB，ap1 dB。凯瑟窗是全面地反映主瓣与旁瓣衰减之间的交换关系，可以在它们两者之间自由地选择主瓣与旁瓣衰减的比重。图示的4-5及4-6分别为利用凯瑟窗设计的低通滤波器和滤波前后的语音信号图： 图4-5

30、 凯瑟窗设计低通滤波器的频谱图 图4-6 滤波前后的信号的比较加噪声后的语音信号经过FIR滤波器滤波后，时域和频域图与原始信号几乎一样，说明噪声几乎都被滤掉了，也说明FIR滤波器设计很理想，能满足所需要求。使用凯瑟窗设计的FIR低通滤波器滤波以后，听到的声音稍微有些发闷、低沉，原因是高频分量被FIR低通滤波器衰减。但是得到的已经是很接近原来的声音了。4.6 结果分析通过对比分析可知，滤波后的输出波形和原始语音加噪声信号的图形发生了一些变化：滤波后的输出波形明显在对应时间幅度比原语音加噪声信号的要小，而且滤波的效果也与滤波器的选择有关，可以看出滤波器的性能差异以及参数的调节会对滤波器产生一定的影

31、响，并且通过回放可以发现滤波前后的声音有变化.低通滤波后,已很接近原来的声音。从频谱图中我们还可以看出声音的能量信号主要集中在低频部分，说明高频语音信号被滤出，滤波器达到既定要求。在采样时要注意音频要为单声道，否则会提示出错，无法利用freqz绘制原始语音信号的频率响应图。整个过程的程序代码如附录所示：总结语音信号是一种非平稳的时变信号，它携带着各种信息。语音信号分析的目的就在与方便有效的提取并表示语音信号所携带的信息。语音信号分析可以分为时域和变换域等处理方法，其中时域分析是最简单的方法，直接对语音信号的时域波形进行分析，提取的特征参数。语音信号经过FFT变换之后，就可以得到信号的频域图形，

32、另外，频域图形使信号的某些特性变得更明显，信号的傅立叶表示在信号的分析与处理中起着重要的作用。因为对于线性系统来说，可以很方便地确定其对正弦或复指数和的响应，所以傅立叶分析方法能完善地解决许多信号分析和处理问题。本课程设计的主要目的是在MATLAB环境下，使用窗函数设计法设计一个FIR低通滤波器，并对叠加了噪声后的语音信号进行滤波去噪。应用MATLAB 根据给定的技术指标进行FIR 数字滤波器设计时，可根据设计要求随时改变参数，以使滤波器达到最优化。 使用凯瑟窗设计低通滤波器，并对叠加噪声后的语音信号进行滤波。比较滤波前后的时域波形和频谱并进行分析，低通滤波后，声音稍微有些发闷、低沉，原因是高

33、频分量被低通滤波器衰减。但是很接近原来的声音。总之，加噪声的后的语音信号经过FIR滤波器的滤波处理，时域和频域图与原始信号都几乎一样，这说明噪声几乎全被滤掉了，同时也说明FIR滤波器设计还算理想，能满足所需要求。也可看出滤波之后仍有一小段噪声没有滤除，设计还有些许的不足。 参考文献1丁美玉.数字信号处理M. 西安电子科技大学出版社 151:2222朱冰莲.数字信号处理M.电子工业出版社 87:933刘树棠译.数字信号处理使用MATLABM.西安交通大学出版社 42:494李勇、徐震.MATLAB辅助现代工程数字信号处理M.西安电子科技大学出版社 48:565约翰··普罗克斯

34、数字信号处理M.西安交通大学出版社，20096张小虹.信号系统与数字信号处理M第1版西安：西安电子科技出版社，20027罗军辉等编著.MATLAB7.0在数字信号处理中的应用M.北京：机械工业出版社.28:348张志涌.精通MATLAB 6.5版M北京：北京航空航天大学出版社，2003附录%声音采集及绘制原始语音信号的时域波形和FFT频谱图x,fs=wavread('G:111.wav')x=x(:,1);FS=length(x);f=0;fs/FS;(FS-1)*fs/FS;sound(x); %播放语音信号X=fft(x,4096);magX=abs(X);angX=ang

35、le(X);y1=fft(x,1024); %对信号做1024点FFT变换f=fs*(0:511)/1024;figure(1)subplot(211);plot(magX);title('原始信号幅值');grid on;subplot(212);plot(angX);title('原始信号相位');grid on;figure(2)subplot(211);plot(x); %绘制原始语音信号的时域波形图title('原始语音信号时域波形图');xlabel('time n');ylabel('fuzhi n'

36、); grid on; Subplot(212); %绘制原始语音信号的频率响应图 plot(f,abs(y1(1:512);title('原始语音信号频谱图')xlabel('Hz');ylabel('fudu'); grid on;%添加随机噪声及添加噪声后的时域图和频谱图noise_mu=0;noise_var=0.05; x0=randn(size(x).*sqrt(noise_var)+noise_mu; x1=x+x0; ts=1/fs; % 绘制在原始信号上加随机噪声的信号图 ta=(length(x)-1)/fs; t=0:ts:ta;figure(3);subplot(211); plot(t,x1); title('加随机噪声后语音信号时域图'); xlabel('t');ylabel('x1');grid on;y2=fft(x1,1024); %对信号做1024点FFT变换f=fs*(