shuzixinhaochuli

1、 徐刚 语音信号滤波去噪使用TUKEYWIN窗设计的FIR滤波器 第20页 共19页 语音信号滤波去噪使用TUKEYWIN窗设计的FIR滤波器学生姓名：徐刚 指导老师：胡双红摘 要 本课程设计是使用TUKEYWIN窗设计的FIR滤波器对语音信号进行滤波去噪。主要使用MATLAB7.0进行实现。先用windows工具中的录音机采集一段语音信号；再用wavread函数求出语音信号的三个参数，读取录取的信号，绘制其时域和频谱图；再对信号进行傅里叶变化，绘制出时域和频域图，并进行对比；最后回放语音信号，对比其滤波前后的信号变化。此次课程设计成功的对语音信号进行了滤波去噪，初步完成设计指标。关键词 滤波

2、去噪；TUKEYWIN窗设计法；FIR滤波器；MATLAB7.01 引 言 滤波去噪是信号处理中的一种重要的技术。利用滤波去噪可以从复杂的信号中提取需要的信号。通过数字滤波器我们可将输入信号的某些频率成分或频带进行压缩放大，从而改变输入信号的频谱结构，也就是一个频率选择器。通过噪音与所需信号的频段的不同进行选择达到滤波去噪的效果。本次课程设计是采用TUKEYWIN窗设计FIR滤波器对语言信号进行滤波去噪。根据FIR滤波器的特点设计滤波器后，通过对比滤波前后的波形图及回放滤波前后的语音信号，可以看出滤波器对有用信号的无失真放大具有重要作用。1.1课程设计目的（1） 了解FIR滤波器的特性和特点；

3、（2） 理解用TUKEYWIN窗设计法设计FIR滤波器，并掌握滤波去噪的方法；（3） 熟悉MATLAB软件下有关函数的调用；（4） 学会对设计指标的分析。1.2课程设计的要求（1）滤波器指标必须符合工程实际。（2）设计完后应检查其频率响应曲线是否满足要求。（3）处理结果和分析结论应该一致，而且应符合理论。（4）独立完成课程设计并按要求编写课程设计报告书。1.3设计平台MATLAB是矩阵实验室（Matrix Laboratory）的简称，是美国MathWorks主要包括公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Si

4、mulink两大部分。MATLAB的基本数据单位是矩阵，它的指令表达式与数学、工程中常用的形式十分相似，故用MATLAB来解算问题要比用C，FORTRAN等语言完成相同的事情简捷得多。MATLAB由一系列工具组成。这些工具方便用户使用MATLAB的函数和文件，其中许多工具采用的是图形用户界面。包括MATLAB桌面和命令窗口、历史命令窗口、编辑器和调试器、路径搜索和用于用户浏览帮助、工作空间、文件的浏览器。随着MATLAB的商业化以及软件本身的不断升级，MATLAB的用户界面也越来越精致，更加接近Windows的标准界面，人机交互性更强，操作更简单。而且新版本的MATLAB提供了完整的联机查询、

5、帮助系统，极大的方便了用户的使用。简单的编程环境提供了比较完备的调试系统，程序不必经过编译就可以直接运行，而且能够及时地报告出现的错误及进行出错原因分析。另外Matlab一个高级的矩阵/阵列语言，它包含控制语句、函数、数据结构、输入和输出和面向对象编程特点。用户可以在命令窗口中将输入语句与执行命令同步，也可以先编写好一个较大的复杂的应用程序（M文件）后再一起运行。新版本的MATLAB语言是基于最为流行的C语言基础上的，因此语法特征与C语言极为相似，而且更加简单，更加符合科技人员对数学表达式的书写格式。使之更利于非计算机专业的科技人员使用。而且这种语言可移植性好、可拓展性极强，这也是MATLAB

6、能够深入到科学研究及工程计算各个领域的重要原因【1】。2设计原理有限长单位脉冲响应数字滤波器（FIRDF）的最大优点是可以实现线性相位滤波。而我们知道，无失真传输与滤波处理的条件是，在信号的有效频谱范围内系统幅频响应为常数，相频相应为频率的线性函数。另外，FIRDF是全零点滤波器，硬件和软件实现结构简单，不用考虑稳定性问题。FIR设计方法主要分为两类：第一类是基于逼近理想滤波器特性的方法，包括窗函数法，频率采样法和等波纹最佳逼近法，第二类是最优设计法。窗函数设计法的基本思想是用FIR逼近希望的滤波特性。设希望逼近的滤波器的频率响应函数为，其单位脉冲响应用表示。为了设计简单方便，通常选择为具有片

7、段常数特性的理想滤波器。因此是无限长非因果序列，不能直接诶作为FIR的单位脉冲响应。窗函数设计法就是截取为有限长的一段因果序列，并用合适的窗函数进行加权作为FIR的单位脉冲响应。截取的长度和加权窗函数的类型都直接影响逼近精度（滤波器指标）【2】。2.1FIR滤波器数字滤波器的设计师数字信号处理中的一个重要内容。数字滤波器设计包括FIR（有限单位脉冲响应）滤波器和IIR（无限单位脉冲响应）滤波器两种。与IIR滤波器相比，FIR滤波器在保证幅度特性满足技术要求的同时，很容易做到严格线性相位特性。设FIR滤波器单位脉冲响应长度为N，其系统函数为：是z1的N1次多项式，在z平面上有N1个零点。远点z=

8、0是N1阶重极点，因此是永远稳定的。稳定和线性相位特性是FIR滤波器突出的优点。FIR滤波器的设计任务是选择优先长度的好。使传输函数满足技术要求。FIR滤波器的设计方法有很多种，如窗函数法、频率采样法及其他各种优化设计方法。本实验介绍TUKEYWIN窗的设计方法。TUKEYWIN窗：w=tukeywin(n),产生一个n点的TUKEYWIN窗函数3。数字滤波器的设计步骤：（1） 确定目标。在设计一个滤波器之前，必须首先确定一些技术指标。这些技术指标需要来制定。在很多实际应用中，例如语音或者音频信号处理中，数字滤波器常用来实现选频操作。因此，指标的形式一般确定为频域中幅度和单位响应。（2） 逼近

9、。确定了礼拜五的技术指标后，就可以利用数学和DSP的基本原理提出一个滤波器模型来逼近给定的目标。这一步是滤波器设计的所要研究的主演问题，是整个设计过程的重中之重。（3） 硬件或软件实现。上两步的结果得到以差分或者系统函数或者冲激响应描述的滤波器,根据这个描述就可以用硬件或者在计算机上涌软件实现4。FIR滤波器的优点（1） 相位响应可以真正线性；（2） 系统绝对稳定,设计相对容易；（3） 高效实现；（4） 可用DFT实现。2.2 窗口设计法窗函数设计法的基本原理是用有限长单位脉冲响应序列逼近。由于往往是无限长序列，而且是非因果的，所以用窗函数将截断，并进行加权处理，得到：就作为实际设计的FIR数

10、字滤波器的单位脉冲响应序列，其频率响应函数为我们知道，用窗函数法设计的滤波器性能取决于窗函数的类型及窗口长度N的取值。设计过程中，要根据对阻带最小衰减和过渡带宽度的要求选择合适的窗函数类型和窗口长度N 。各种类型的窗函数可达到的阻带最小衰减和过渡带宽度见表2-1。表2-1 常见窗函数性能表名称滤波器过渡带宽最小阻带衰减名称滤波器过渡带宽最小阻带衰减矩形1.8/M21dBPARZENWIN6.6/M56db巴特利特6.1/M25dBFLATTOPWIN19.6/M108db汉宁6.2/M44dBGAUSSWIN5.8/M60db汉明6.6/M51dBBARTHANNWIN3.6/M40db布莱克

11、曼11/M74dBBLACKMANHARRIS16.1/M109dbBOHMANWIN5.8/M51.5dbCHEBWIN15.2/M113dbNUTTALLWIN15.4/M108dbTUKEYWIN2.4/M22db在基于逼近理想滤波器特性的方法中，运用窗函数设计FIR滤波器可以利用MATLAB平台上的窗函数，根据窗函数的性能表设置滤波器的基本指标结构化程序设计中关键是如何将问题域中的行为（即操作）抽取出来。而在面向对象程序设计中关键是如何将问题域中的实体（即日常所见的概念）抽取出来，作为C+程序中的类，而属性与行为作为类的两类要素通常是必不可少的，甚至还应考虑类必须满足的约束。2.3TU

12、KEYWIN窗函数的设计方法窗函数的理论依据：窗函数设计滤波器的基本思想，就是依据给定的滤波器的技术指标，选择滤波器的阶数N和合适的窗函数。即用一个有限长度的窗口函数序列来截取一个无限长序列获得一个有限长序列，即，并且要满足以下两个条件：（1）窗谱主瓣尽可能的窄，以获得较陡的过渡带；尽量减少窗谱的最大旁瓣的相对幅度，也就是能量尽量集中于主瓣，使峰肩和纹波减小，就可以增多阻带的衰减。这就给窗函数序列的形状和长度选择提出了严格的要求。3设计步骤3.1 设计流程图设计本课题的流程为：先利用Windows下的录音机，录制语音信号“天气很好”，时间在2-3 s左右。然后在Matlab软件平台下，利用函数

13、wavread对语音信号进行采样，记住采样频率和采样点数。通过wavread函数的使用，求出音乐信号的三个参数，分别为：每个样本的值，生成该音乐波形文件时的采样频率，波形文件样本的码数；再对信号做傅里叶变化，绘制出时域和频域的波形；最后通过滤波绘制滤波前后时域必行对比图和幅频特性对比图，并回放滤波前后的音乐信号来验证是否达到去噪的目的。流程图如图3-1所示：开始录制一段语音信号设计一个TUKEYWIN窗的FIR滤波器As，Ap是否达标否是信号滤波处理滤波前后各图的比较结束图3-1 设计流程图 3.2 录制语音信号用工具中的录音机录制一段语音信号，格式为.wav。将音乐信号的文件命名为xg.wa

14、v,然后在Matlab软件平台下，利用函数wavread对语音信号进行采样，记住采样频率fs和采样点数bits。采集完成后在信号中加入一个单频噪声，设计的任务即为从含噪信号中滤除单频噪声，还原原始信号。具体如下：x,fs,bits=wavread(xg.wav)文件的全路径和文件名，输出的第一个参数是每个样本的值,fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。sound(x,fs,bits); % 按指定的采样率和每样本编码位数回放N=length(x); % 计算信号x的长度fn=6000; % 单频噪声频率，此参数可改t=0:1/fs:(N-1)/fs; % 计算时间

15、范围，样本数除以采样频率x=x; y=x+0.1*sin(fn*2*pi*t); sound(y,fs,bits); % 应该可以明显听出有尖锐的单频啸叫声subplot(2,2,3);plot(y);title(加入单频干扰后的语音信号);xlabel(时间（单位：s);ylabel(幅度);X=abs(fft(x); Y=abs(fft(y); % 对原始信号和加噪信号进行fft变换,取幅度谱X=X(1:N/2); Y=Y(1:N/2); % 截取前半部分deltaf=fs/N; % 计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf; % 计算频谱频率范围 用绘图命令分别画出

16、加噪前后信号的时域和频域波形，如图3-2所示：图3-2 加噪前后信号的时域和频域波形图3.3 滤波器设计滤波器设计就是要找到一组能满足特定滤波要求的系数向量a和b，而它主要是通过设计指标来实现的。滤波器设计的要求或指标一般是在频域上给出的，常用的滤波器频域指标有：通带截止频率fp，阻带截止频率fc，通带波纹Rp，阻带衰减As。要达到最佳的滤波效果，则需要对fp，fc和As进行适当的调整。由图3-2可以看出，音乐信号可以选择fpd=2000；fsd=2050；fsu=2150；fpu=2200；Rp=1；As=50的滤波器。在MATLAB中，通常采用1/2采样频率进行归一化处理，如果将频率转化为

17、角频率，则需将归一化频率乘以pi。设计程序：fpd=2000;fsd=2050;fsu=2150;fpu=2200;Rp=1;As=50; % 带阻滤波器设计指标fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=min(fsd-fpd),(fpu-fsu); % 计算上下边带中心频率，和频率间隔wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi; % 将Hz为单位的模拟频率换算为rad为单位的数字频率wsd=fsd/fs*2*pi;wsu=fsu/fs*2*pi;M=ceil(2.4*pi/dw)+1; % 计算TUKEYWIN窗设计

18、该滤波器时需要的阶数n=0:M-1; % 定义时间范围w_tuk=tukeywin（M）; %产生M阶的TUKEYWIN窗hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); % 调用自编函数计算理想带阻滤波器的脉冲响应h_bs=w_ham.*hd_bs; % 用窗口法计算实际滤波器脉冲响应db,mag,pha,grd,w=freqz_m(h_bs,1); % 调用自编函数计算滤波器的频率特性通过运行上面一段程序，MATLAB设计出来的滤波器的真正参数为Rp=1，As=22，符合FIR滤波器的要求，那么设计的数字滤波器基本实现设计的要求。所

19、设计的滤波器的幅度响应和脉冲响应如图3-3所示：图3-2 滤波器的幅度响应图和脉冲响应图3.4 信号滤波处理滤波器设计完成后，在MATLAB平台上涌函数fftfilt实现滤波，那么程序如下：y_fil=filter(h_bs,1,y); % 用设计好的滤波器对y进行滤波Y_fil=fft(y_fil);Y_fil=Y_fil(1:N/2); % 计算频谱取前一半sound (y_fil,fs,bits); %音乐信号回放通过运行上面一段程序，便得到原始信号x，加噪信号y，滤波去噪信号y_fil的时域波形和频谱，如图3-3所示：图3-3 原始信号x，加噪信号y，滤波去噪信号y_fil的时域波形和

20、频谱图3.5 结果分析在MATLAB中，经过sound（y，fs，bits）函数，可听到经过FIR滤波器之后的语音信号尖锐的单频啸叫声已消失。所得结果，证明了FIR滤波器和音乐信号去噪设计是成功的。经过滤波器的滤波去噪，成功的消去噪音，还原了语音信号“天气很好”。这已完成了设计指标。而对声音开始时的噼啪声，经分析应是延时产生的。4出现的问题及解决方法4.1出现的问题在这次设计中，我遇到了两个问题：（1） 滤波器的性能不如人意，不能完全滤掉高频噪声；（2） 音频开始时有噼啪声。4.2解决方法通过询问老师，我得知改变通频带就可以改善性能。我尝试去改变相应的数据，最终我找到了合适的数据。得到合适的通

21、频带使得滤波器性能得到了改善，消除了噪音，最后完成了设计要求。对开始时一下的噼啪声，一开始我没能寻找到问题产生的原因。所幸我在询问同学是否有相同的状况时，颜知生同学告诉我这是延时的结果，是无法避免的。因而第二个问题也由此得到了解决。5 结束语经过这两周的课程设计，我大体上完成此次的目标。这次的课程设计使我对数字信号处理这门课程有了更深入的认识。让我得到了一次知识结合实践的机会，这是非常宝贵的经验。在课程设计中，我锻炼了自己的能力，加强了自身分析问题、解决问题的训练。在这次设计中自己也遇到一些困难，但很庆幸有老师和同学的大力帮助得以顺利的完成这次设计。总的来说这次设计使我受益匪浅，锻炼了自身的能

22、力，并加强了对知识的回顾与运用。虽说还有不足的地方，设计的也非完美，但我相信以后会做的更好。参考文献1孙祥，徐流美，吴清.MATLAB7.0基础教程.北京：机械工程出版社，20052高西全，丁玉美，阔永红.数字信号处理原理、实现及应用.北京：电子工业出版社，20063王正林，刘明.精通MATLAB7.北京：电子工业出版社，20064Richard G.Lyons.Understanding Digital Signal Processing.北京：科学出版社，2003.3附录1：语音信号采集源程序%程序名称：语音信号的采集%程序功能：实现语音信号的采集，并且画出音频信号的时域和频域图%程序作者

23、：徐刚%最后修改日期：2010-3-11x,fs,bits=wavread(xg.wav); % 输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值,fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。sound(x,fs,bits); % 按指定的采样率和每样本编码位数回放N=length(x); % 计算信号x的长度fn=6000; % 单频噪声频率，此参数可改t=0:1/fs:(N-1)/fs; % 计算时间范围，样本数除以采样频率x=x; y=x+0.1*sin(fn*2*pi*t); sound(y,fs,bits); % 应该可以明显听出有尖锐的单频啸叫

24、声subplot(2,2,3);plot(y);title(加入单频干扰后的语音信号);xlabel(时间（单位：s);ylabel(幅度);X=abs(fft(x); Y=abs(fft(y); % 对原始信号和加噪信号进行fft变换,取幅度谱X=X(1:N/2); Y=Y(1:N/2); % 截取前半部分deltaf=fs/N; % 计算频谱的谱线间隔f=0:deltaf:fs/2-deltaf; % 计算频谱频率范围subplot(2,2,1);plot(t,x);title(原始语音信号);xlabel(时间（单位：s);ylabel(幅度);subplot(2,2,2);plot(f

25、,X);title(语音信号幅度谱图);xlabel(频率（单位：Hz);ylabel(幅度谱);subplot(2,2,3);plot(t,y);title(加入单频干扰后的语音信号);xlabel(时间（单位：s);ylabel(幅度);subplot(2,2,4);plot(f,Y);title(加入单频干扰后的语音信号幅度谱图);xlabel(频率（单位：Hz);ylabel(幅度谱);附录2：滤波器设计源程序%程序名称：滤波器的设计%程序功能：用FIR滤波器原型，再用TUKEYWIN窗函数设计法将其离散化为数字滤波器。%程序作者：徐刚%最后修改日期：2010-3-11fpd=2000

26、;fsd=2010;fsu=2190;fpu=2200;Rp=1;As=22; % 带阻滤波器设计指标fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=min(fsd-fpd),(fpu-fsu); % 计算上下边带中心频率，和频率间隔wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi; % 将Hz为单位的模拟频率换算为rad为单位的数字频率wsd=fsd/fs*2*pi;wsu=fsu/fs*2*pi;M=ceil(2.4*pi/dw)+1; % 计算TUKEYWIN窗设计该滤波器时需要的阶数n=0:M-1; % 定义时间范围w_tuk=tukeywin（M）；产生M阶的TUKEYWIN窗hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); % 调用自编函数计算理想带阻滤波器的脉冲响应h_bs=w_ham.*hd_bs; % 用窗口法计算实际滤波器脉冲响应db,mag,pha,grd,w=freqz_m(h_bs,1); % 调用自编函数计算滤波器的频率特性subplot(2,2,1);plot(w/pi,db);grid on;title(滤波器幅度响应图);xlabe