《语音信号滤波去噪——使用汉宁窗设计的FIR滤波器.docx》由会员分享，可在线阅读，更多相关《语音信号滤波去噪——使用汉宁窗设计的FIR滤波器.docx（24页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页语音信号滤波去噪使用汉宁窗设计的FIR 滤波器学生姓名： 指导教师：摘要 本课程设计主要是对一段语音信号，参加噪声后，用汉宁窗设计出的FIR 滤波器对参加噪声后的语音信号进展滤波去噪处理。在此次课程设计中，系统操作平台为 Windows XP，程序设计的操作软件为 MATLAB 7.0。此课程设计首先是用麦克风采集一段语音信号， 参加噪声，然后承受汉宁窗函数法设计出FIR 滤波器，再用设计出的滤波器对这段加噪后的语音信号进展滤波去噪，最终对前后时域和频域的波形图进展比照分析，从波形可以看出噪声被完全滤除，到达了语音不失真的效果，说
2、明此次设计格外成功。关键词 程序设计；滤波去噪；FIR 滤波器；汉宁窗；MATLAB 7.01 引 言本课程设计主要是对一段语音信号，进展加噪后，用某种函数法设计出的FIR 滤波器对参加噪声后的语音信号进展滤波去噪处理，并且分析比照前后时域和频域波形的程序设计。1.1 课程设计目的在此次课程中主要的要求是用麦克风采集一段语音信号，绘制波形并观看其频谱，给定相应技术指标，用汉宁窗设计一个满足指标的FIR 滤波器，对该语音信号进展滤波去噪处理，比较滤波前后的波形和频谱并进展分析，依据结果和学过的理论得出合理的结论。与不同信源一样滤波方法的同学比较各种信源的特点，与一样信源不同滤波方法的同学比较各种
3、滤波方法性能的优劣。通过此次课程设计，我们能够学会如何综合运用这些学问，并把这些学问运用于实践当中，使所学学问在综合运用力量上以及分析问题、解决问题力量上得到进一步的进展，语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页让自己对这些学问有更深的了解。通过课程设计培育严谨的科学态度，认真的工作作风和团队协作精神。1.2 课程设计的要求（1） 滤波器指标必需符合工程实际。（2） 设计完后应检查其频率响应曲线是否满足指标。（3） 处理结果和分析结论应当全都，而且应符合理论。（4） 独立完成课程设计并按要求编写课程设计报告书。1.3 工作平台简介课程设计的主要设计平台式 MATLAB
4、7.0。如以下图 1-1 所示：MATLAB 的名称源自Matrix Laboratory ，它是美国 MathWorks 公司生产的一个为科学和工程计算特地设计的交互式大型软件，是一个可以完成各种准确计算和数据处理的、可视化的、强大的计算工具。它集图示和准确计算于一身，在应用数学、物理、化工、机电工程、医药、金融和其他需要进展简单数值计算的领域得到广泛应用。它不仅是一个在各类工程设计中便于使用的计算工具，而且也是一个在数学、数值分析和工程计算等课程教学中的优秀的教学工具，在世界各地的高等院校中格外流行，在各类工业应用中更有不俗的表现。MATLAB 可以在几乎全部的 PC 机和大型计算机上运行
5、，适用于Windows、UNIX 等各种系统平台1。总的来说，该软件有三大特点。一是功能强大。具有数值计算和符号计算、计算结果和编程可视化、数学和文字统一处理、离线和在线计算等功能；二是界面友善、语言自然。MATLAB 以复数处理作为计算单元，指令表达与标准教科书的数学表达式相近；三是开放性强。当学好MATLAB 的同时，会更好的帮助自己去就解决一些难题，而且MATLAB 拥有格外好的进展前途，对我们将来的帮助也是不行限量的。语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页图 1-1 MATLAB 7.0 设计平台2 设计原理2.1 数字信号处理简介数字信号处理Digital
6、Signal Processing，简称DSP是一门涉及很多学科而又广泛应用于很多领域的兴学科2。20 世纪 60 年月以来，随着计算机和信息技术的飞速进展， 数字信号处理技术应运而生并得到快速的进展。在过去的二十多年时间里，数字信号处理已经在通信等领域得到极为广泛的应用。数字信号处理是利用计算机或专用处理设备，以数字形式对信号进展采集、变换、滤波、估值、增加、压缩、识别等处理，以得到符合人们需要的信号形式。随着信息技术的迅猛进展，数字信号处理已成为一个极其重要的学科和技术领域。在通信、语音、图像、自动掌握和家用电器等众多领域得到了广泛的应用。数字滤波是数字信号处理的重要环节，它在数字信号处理
7、中占有着重要的地位,它具有牢靠性好、精度高、敏捷性大、体积小、重量轻等优点。随着数字技术的进展，数字滤波器越来越受到人们的重视，广泛地应用于各个领域。数字滤波器的输入输出信号都是数字信号，它是通过肯定的运算过程转变输入信号所含频率成分的相比照例或者滤除某些频率成分来实现滤波的， 这种运算过程是由乘法器、加法器和单位延迟器组成的。数字滤波器是数字信号处理技术语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页的重要内容，其对数字信号进展的最常见处理是保存数字信号中的有用频率成分和去除信号中的无用频率成分。依据时间域的特性，数字滤波器可以分为无限冲激脉冲响应数字滤波器IIR 滤波器和有
8、限冲激脉冲响应数字滤波器FIR 滤波器3。2.2 FIR 滤波器有限长单位脉冲响应数字滤波器Finite Impulse Response Digital Filter,缩写 FIRDF：有限长单位冲激响应滤波器，是数字信号处理系统中最根本的元件，最大优点是可以实现线性相性滤波，它可以在保证任意幅频特性的同时具有严格的线性相频特性，同时其单位抽样响应是有限长的，因而滤波器是稳定的系统。因此，FIR 滤波器在通信、图像处理、模式识别等领域都有着广泛的应用。FIR 滤波器的设计方法主要分为两类：第一类是基于靠近抱负滤波器器特性的方法包括窗函数法、频率采样法、和等波浪最正确靠近法；其次类是最优设计法
9、。设 FIRDF 的单位脉冲响应 h(n) 的长度为 N ，则其频率响应函数为H (e jw ) = N -1 h(n)e- jwn n=0一般将 H (e jw) 表示成如下形式：2-1H (e jw ) = H (w )e jq (w)2-2g式中，H(w ) 是w 的实函数可以去负值。与前面的表示形式，即 H ( e jw ) =H(w ) e jq ( w )gg相比， Hg(w) 与w 不同。q(w) 与 j(w) 不同。为了区分于幅频响应函数 H (e jw) 和相频响应函数j(w) ，称 Hg(w) 为幅频特性函数，称q(w) 为相频特性函数。第一类线性相位 FIRDF 的相位特
10、性函数是w 的严格线性函数：q (w ) = -wt2-3其次类线性相位 FIRDF 的相位特性函数如下：q (w ) = q 0 - wt2-4式中，t 是常数， q 0 是起始相位。q0 = -p / 2 在信号处理中很有有用价值如希伯尔特变换器，这是 FIRDF 除了线性相位滤波外，还具有真正交变换作用。语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页2.3 窗口设计法()(FIR 滤波器的设计方法有很多种，如窗函数设计法、频率采样设计法和最优化设计法等。窗口设计法的根本思想是用 FIRDF 靠近期望的滤波特性。设期望靠近的滤波器的频率响应为 Hde jw，其单位脉冲响应
11、用 hd(n)表示。为了设计简洁便利，通常选择 Hde jw )为具有片段常数特性的抱负滤波器。因此 hd(n)是无限长非因果序列，不能直接作为 FIRDF的单位脉冲响应。窗口设计法就是截取hd(n)为有限长的一段因果序列，并用适宜的窗口函数进展加权作为 FIRDF 的单位脉冲响应h(n)。窗口设计法根本步骤如下：（1） 构造期望靠近的频率响应函数 H(de jw)。以低通线性相位 FIRDF 设计为例，一般选择 Hde jw)为线性抱负低通滤波器，即() e - jwt ， w wHe jw d= 0，cw w p(c2-5（2） 求出hd(n)。对 Hde jw )进展 IFT 得到( )
12、1p( )sinw (n-t)h n = 2p H ejw ejwndw=p(c -t)2-6d-pdn（3） 加窗得到 FIRDF 的单位脉冲响应h(n)，h(n)= h (n)w(n)2-7d式中， w(n)称为窗口函数，其长度为N 。假设要求第一类线性相位FIRDF，则要求h(n)关于(N - 1)/ 2 点偶对称。而hd(N - 1)/ 2 点偶对称。(n)关于n = t 点偶对称，所t = (N - 1)/ 2 ，同时要求w(n)关于常见的窗函数，可以分为以下主要类型：（1） 幂窗-承受时间变量某种幂次的函数，如矩形、三角形、梯形或其它时间t 的高次幂；（2） 三角函数窗-应用三角函
13、数，即正弦或余弦函数等组合成复合函数，例如汉宁窗、海明窗等；语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页（3） 指数窗-承受指数时间函数，例如高斯窗等4。其性能如表 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/M74dBBLACKMANHARRIS16.1/M1
14、09dBBOHMANWIN5.8/M51.5dBCHEBWIN15.2/M113dBNUTTALLWIN15.4/M108dBTUKEYWIN2.4/M22dB2.4 汉宁窗Hanning window汉宁窗函数是余弦平方函数，又称之为升余弦函数，它的时域形式可以表为：w(k ) = 0.5(1- cos(2pk)2-8n +1其中k = 1,2, n 。它的频域幅度特性函数为：W (w) = 0.5W (w) + 0.25W(w -2p) + W(w +2p) e- j ( N -1)wRRN -1RN -122-9其中 WR(w ) 为矩形窗函数的幅度频率特性函数。汉宁窗函数的最大旁瓣值比
15、主瓣值低31dB，但是主瓣宽度比矩形窗函数的主瓣宽度增加了一倍，为 8p 。汉宁窗函数的时域幅N度与频域幅度特性曲线的 MATLAB 实现的曲线图如图 2-1 所示。语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页图 2-1 汉宁窗函数的时域幅度与频域幅度特性曲线3 设计步骤3.1 设计流程图本课程设计主要是用麦克风采集一段语音信号，通过进展自编函数参加噪声，然后承受汉宁窗函数法设计 FIR 滤波器，并且对这段参加噪声的语音信号函数进展滤波去噪，用绘图的程序画出前后时域和频域的波形图进展比照分析。程序的设计流程图如以下图 3-1 所示：语音信号滤波去噪使用汉宁窗设计的FIR 滤
16、波器第页 共 23 页开头在 Windows 下录制语音将语音格式改参加单频噪声对语音信号进展频谱分析，画出时域和频域波形图用汉宁窗设计FIR 滤波器画出其频率响用 FIR 滤波器对语音信号进展滤波画出语音信号滤波前后波形并且进展比较分析完毕图 3-1 程序设计流程图3.2 录制语音信号单击桌面左下角的“开头”命令，找到“附件”中的“消遣”命令，消灭“录音机” 这一项，单击这一项，弹出窗口，如以下图 3-2 所示，接下来开头录制语音信号“大家好， 我是 xxx”，时间为 2-3 秒，如以下图 3-3 所示：图 3-2 未录制信号时录音机的状态图 3-3 录制信号后录音机的状态语音信号滤波去噪使
17、用汉宁窗设计的FIR 滤波器第页 共 23 页3.3 语音加噪处理采集完成后在信号中参加一个单频噪声，设计的任务即为从含噪信号中滤除单频噪声，复原原始信号。所以我们首先对采集来的信号进展加噪处理，调用的程序如下：x,fs,bits=wavread(”d:luo.wav”);% 输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs 是生成该波形文件时的采样率，bits 是波形文件每样本的编码位数。sound(x,fs,bits);% 按指定的采样率和每样本编码位数回放N=length(x);% 计算信号 x 的长度fn=1700;% 单频噪声频率t=0:1/fs:(N-1)/fs
18、;% 计算时间范围，样本数除以采样频率y=x”+0.01*sin(fn*2*pi*t);% 参加一个单频噪声sound(y,fs,bits);% 应当可以明显听出有锋利的单频啸叫声对语音信号参加单频噪声前后的波形进展分析，首先画出语音信号的时域波形；然后对语音号进展快速傅里叶变换，得到信号的频谱特性。调用的程序如下：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;% 计算频谱频率范围语音信号滤波去
19、噪使用汉宁窗设计的FIR 滤波器第页 共 23 页用绘图命令分别画出加噪前后信号的时域和频域波形，如以下图 3-4 所示：图 3-4 语音信号参加单频噪声前后的时域和频域波形图由上图 3-4 我们可以看到语音信号参加单频噪声后的时域波形比未加之前在幅度范围内有了明显的增加，在频谱方面我们可以看到除了在加了噪声后的频谱图上的 1700hz 有个明显的冲激外，其余地方与未加时的频谱一模一样，这现象表现在语音播放时我们可以听到一声锋利的噪声。3.4 滤波器设计在 Matlab 中，可以利用矩形窗、三角窗、汉宁窗、汉明窗、布莱克曼窗、凯塞窗等设计 FIR 滤波器，在本次课程设计中主要应用汉宁窗设计出
20、FIR 滤波器。利用 Matlab 中的函数 freqz 画出各滤波器的频率响应，首先利用数字信号处理里面学过的学问，依据自己选定的参数，用汉宁窗函数法设计FIR 数字滤波器，得到数字滤波器的参数b，a。其中 b 为系统函数的分子系数，a 为系统函数分母系数。再调用 freqz(b,a,512,fs)即可得到该滤波器的频率响应。主程序如下：fpd=1600;fsd=1650;fsu=1750;fpu=1800;% FIR 滤波器的上下截止频率Rp=1;As=37;语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页% 带阻滤波器设计指标fcd=(fpd+fsd)/2;fcu=(fp
21、u+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(6.2*pi/dw)+1;% 计算汉宁窗设计该滤波器时需要的阶数n=0:M-1;% 定义时间范围w_ham=hanning(M);% 产生 M 阶的汉宁窗hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu
22、,M);% 调用自编函数计算抱负带阻滤波器的脉冲响应h_bs=w_ham”.*hd_bs;% 用窗口法计算实际滤波器脉冲响应db,mag,pha,grd,w=freqz_m(h_bs,1);% 调用自编函数计算滤波器的频率特性通过绘图工具我们可以得出滤波器的波形图，如以下图 3-5 所示:语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页图 3-5FIR 滤波器的频率响应上图 3-5 为用汉宁窗函数法设计出的 FIR 滤波器图，我们可以看出，阻带最大衰减为-75dB,FIR 滤波器的主瓣宽度很小,这样可以使过渡带很陡，旁瓣相对于主瓣也比较小。3.5 信号滤波处理用自己设计的各滤波
23、器分别对采集的信号进展滤波，在Matlab 中，FIR 滤波器利用函数 fftfilt 对信号进展滤波，IIR 滤波器利用函数 filter 对信号进展滤波，对语音信号进展滤波后，认真比照滤波前和滤波后的语音信号图，得出结论。主程序如下：y_fil=filter(h_bs,1,y);% 用设计好的滤波器对 y 进展滤波Y_fil=fft(y_fil);Y_fil=Y_fil(1:N/2);% 计算频谱取前一半由绘图工具我们可以得出滤波后的语音信号波形图、原始语音信号的波形图以及参加噪声后的语音信号波形图，图如下 3-6 所示：语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页图
24、3-6 滤波前后的波形图由上图 3-6 我们可以看出，加噪声的后的语音信号经过 FIR 滤波器的滤噪处理，时域和频域图都几乎完全一样，这说明噪声被完全滤掉，同时也说明 FIR 滤波器设计很抱负， 能满足所需要求。3.6 结果分析语音信号经过 FIR 滤波器的滤除噪声的处理，在 Matlab 中，函数 sound 可以对声音进展回放。其调用格式：sound (x，fs，bits)我们可以明显感觉滤波前后的声音有变化。声音中刺耳的声音没有了，几乎恢复成原始的声音，但较原始的声音更平滑一些。这说明用汉宁窗设计FIR 滤波器滤掉了语音中的噪声同时，也把原始语音的很小的一局部也滤掉了，所以回放语音的时候
25、听起来比以前的更加平滑，说明这段程序设计是成功的。语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页4 消灭的问题及解决方法在设计课程设计时，消灭了以下几个问题：1. 在录音时，由于没能掌握好时间，所以消灭了时间上的错误，经屡次修正后，最终成功。2. 调用语音时，在用语句“x,fs,bits=wavread(”d:luo.wav”)”时，消灭错误，后来知道语句的单引号用了中文输入法，最终将单引号的输入改为英文输入法，之后，还觉察在调用路径出了错误，我又将其放到 d 盘内，这才正确。3. 参加单频噪声时，单频噪声为正弦信号，其系数为 1， 得出的图形在幅度变化范围内很大，几乎不能看
26、到振幅的起伏变化，后来经胡教师的教导我觉察了是系数太大，将其变为 0.01 后，得出正确结果。4. 绘制出加噪声后的波形图，所加单频噪声频率值不在所定义的 1700hz 上，我开头以为是单频噪声的值有问题，所以屡次修改后仍没有改正，后来在教师的帮助下觉察是计算频谱的谱线间隔时出错了，就是“ deltaf=fs/2/length(X)”此处多除了一个 2 ，将其改为“deltaf=fs/length(X)”后，谱线间隔增加一倍，“X”和“Y”在时间范围变得全都，所以才能得出正确结果。5. 由于单频噪声频率的屡次改动，所设计的滤波器的性能指标没有随之变化，所以消灭了滤波上的错误，经过认真检查后，觉
27、察是性能指标的错误，没能以单频噪声的频率为中心频率，所以经过修改指标，滤波器得以正确设计。语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页5 完毕语两周的课程设计完毕了，在此次课程设计中，我收获了很多，不仅在学问学习方面， 而且在编程力量方面、团结合作力量等方面都有了肯定的成就。首先，在学问学习方面，虽然我们开头有学过数字信号处理这一门课程，有肯定的根底，但是我觉察我还没能学懂学通透，学问没能很好的与实践结合，特别在设计指标方面， 没能全面考虑，所以经过这两周的课程设计，我觉察学的远远多于书本上的东西，而且印象深刻。其次，在编程力量方面，此次由于有了胡教师的具体讲解和她供给的具
28、体资料，所以， 我能够很快知道接下来应当怎样编程，有什么作用，得出什么结果，这让我们始终有着明确的目标，让我们更加有信念。最终，在团结合作力量方面，虽然每个人有每个人的题目，但是我们还是充分发挥了我们的团结合作力量，找到资料的同学帮没找到的同学找，有什么错误的大家一起争论， 一起找资料解决，特别是题目相像的同学，大家有资料都一起争论，而且还有我们的指导教师胡教师的细心指导，所以我们合作得很开心。在这次设计过程中，既表达出了自己单独设计的力量以及综合运用学问的力量，又体会了学以致用、突出自己劳动成果的喜悦心情，从中也能觉察自己寻常学习的缺乏和薄弱环节，从而加以弥补。同时，也再次体会到了团结合作的
29、欢快。在此我要格外感谢敬重的胡教师，由于你给了我一个这么好的题目，要我学到了更多的东西，还有你的悉心帮助，让我的设计能够成功完成，同时，也要感谢好同伴们给我 的帮助。语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页参考文献1 张圣勤. MATLAB7.0 有用教程M. 北京：机械工程出版社，20232 维纳K 恩格尔，约翰G 普罗克斯.数字信号处理M. 西安：西安交通大学出版社，2023.63 飞思科技产品研发技术中心.MATLAB7 关心信号处理技术与应用M.北京：电子工业出版社，2023(3).4 百科 ROBOT,zivenwong.窗函数.百度百科, :/baike.b
30、aidu /view/2327611.htm?fr=ala0_1 ：2023-12-26语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页附录 1：用汉宁窗设计 FIR 滤波器程序清单%程序名称： FIR.M%程序功能：承受汉宁窗设计 FIR 滤波器的方法，给语音出去噪声。%程序 %最终修改日期： 2023-3-1 x,fs,bits=wavread(”d:luo.wav”);% 输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs 是生成该波形文件时的采样率，bits 是波形文件每样本的编码位数。sound(x,fs,bits);% 按指定的采样率和每样本编码位
31、数回放N=length(x);% 计算信号 x 的长度fn=1700;% 单频噪声频率t=0:1/fs:(N-1)/fs;% 计算时间范围，样本数除以采样频率y=x”+0.01*sin(fn*2*pi*t);% 参加一个单频噪声sound(y,fs,bits);% 可以明显听出有锋利的单频啸叫声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;% 计算频谱频率范围语音信号滤波去噪使用汉宁窗设计的FI
32、R 滤波器第页 共 23 页subplot(2,2,1);plot(t,x);% 布局为 2*2 的四个小图axis(0,2.2,-0.04,0.07);%转变横纵坐标的范围title(”原始语音信号”); xlabel(”时间(单位:s)”); ylabel(”幅度”);%加上标题和横坐标名称grid on;% 加上网格subplot(2,2,2);plot(f,X);axis(0,5000,0,40);title(”原始语音信号频谱”); xlabel(”频率(单位:hz)”); ylabel(”幅度谱”);grid on; subplot(2,2,3);plot(t,y);axis(0,
33、2.2,-0.04,0.08);title(”参加单频干扰后的语音信号”); xlabel(”时间(单位:s)”);ylabel(”幅度”); grid on;subplot(2,2,4);plot(f,Y);axis(0,5000,0,40);title(”参加单频干扰后的语音信号频谱”); xlabel(”频率(单位:hz)”);ylabel(”幅度谱”); grid on; fpd=1600;fsd=1650; fsu=1750;fpu=1800;% FIR 滤波器的上下截止频率语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页Rp=1;As=37;% 带阻滤波器设计指标f
34、cd=(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(6.2*pi/dw)+1;% 计算汉宁窗设计该滤波器时需要的阶数n=0:M-1;% 定义时间范围w_ham=hanning(M);% 产生 M 阶的汉宁窗hd_bs=ideal_lp(wcd,M)+ideal
35、_lp(pi,M)-ideal_lp(wcu,M);% 调用自编函数计算抱负带阻滤波器的脉冲响应h_bs=w_ham”.*hd_bs;% 用窗口法计算实际滤波器脉冲响应db,mag,pha,grd,w=freqz_m(h_bs,1);% 调用自编函数计算滤波器的频率特性figure(2) subplot(2,2,1);plot(w/pi,db);axis(0,0.5,-80,20);title(”以 db 为单位的幅度特性”); xlabel(”w/pi”);ylabel(”db”);grid on; subplot(2,2,2);plot(w/pi,mag); axis(0,0.5,-0.1
36、,1.25);语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页title(”以线性为单位的幅度特性”); xlabel(”w/pi”);ylabel(”mag”);grid on; subplot(2,2,3);plot(w,pha); title(”滤波器相位响应图”); xlabel(”w/pi”);ylabel(”相位 pha”); grid on;subplot(2,2,4);plot(h_bs);axis(0,1400,-0.2,1.5);title(”滤波器脉冲响应图”); xlabel(”n”);ylabel(”h(n)”); grid on; y_fil=fil
37、ter(h_bs,1,y);% 用设计好的滤波器对 y 进展滤波Y_fil=fft(y_fil);Y_fil=Y_fil(1:N/2);% 计算频谱取前一半figure(3) subplot(3,2,1);plot(t,x);axis(0,2.1,-0.05,0.08);title(”原始语音信号时间 x”); xlabel(”时间 t”);ylabel(”幅度”);grid on; subplot(3,2,2);plot(f,X);axis(0,5000,0,40);title(”原始语音信号幅度谱 X”); xlabel(”频率 f”);ylabel(”幅度”); grid on;subp
38、lot(3,2,3);plot(t,y);axis(0,2.1,-0.05,0.1);title(”加干扰语音信号时间 x1”); xlabel(”时间 t”);ylabel(”幅度”);语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页grid on; subplot(3,2,4);plot(f,Y);axis(0,5000,0,40);title(”加干扰语音信号幅度谱 X1”); xlabel(”频率 f”);ylabel(”幅度”);grid on; subplot(3,2,5);plot(t,y_fil);axis(0,2.1,-0.05,0.1);title(”滤波后
39、语音信号时间 y”); xlabel(”时间 t”);ylabel(”幅度”);grid on; subplot(3,2,6);plot(f,Y_fil);axis(0,5000, 0,40);title(”滤波后语音信号幅度谱 Y”); xlabel(”频率 f”);ylabel(”幅度”); grid on;语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页附录 2：function hd = ideal_lp(wc,M);% 抱负低通滤波器计算% % hd = ideal_lp(wc,M)% hd = 0 to M-1 之间的抱负脉冲响应% wc = 截止频率(弧度)% M
40、 = 抱负滤波器的长度alpha = (M-1)/2; n = 0:1:(M-1);m = n - alpha + eps;hd = sin(wc*m) ./ (pi*m);语音信号滤波去噪使用汉宁窗设计的FIR 滤波器第页 共 23 页附录 3：function db,mag,pha,grd,w = freqz_m(b,a);% freqz 子程序的改进版本% % db,mag,pha,grd,w = freqz_m(b,a);% db = 0 到 pi 弧度区间内的相对振幅(db)% mag = 0 到 pi 弧度区间内确实定振幅% pha = 0 到 pi 弧度区间内的相位响应% grd
41、 = 0 到 pi 弧度区间内的群拖延% w =0 到 pi 弧度区间内的 501 个频率样本向量% b = Ha(z)的分子多项式系数(对 FIR b=h)% a = Ha(z)的分母多项式系数(对 FIR: a=1)% H,w = freqz(b,a,1000,”whole”); H = (H(1:1:501)”; w = (w(1:1:501)”;mag = abs(H);db = 20*log10(mag+eps)/max(mag); pha = angle(H);% pha = unwrap(angle(H); grd = grpdelay(b,a,w);% grd = diff(pha);% grd = grd(1) grd;% grd = 0 grd(1:1:500); grd; grd(2:1:501) 0;% grd = median(grd)*500/pi;