语音信号滤波去噪 - 使用脉冲响应不变法设计的巴特沃斯滤波器

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第1页 共20页

语音信号滤波去噪——使用脉冲响应不变

法设计的巴特沃斯滤波器

学生姓名：王雄杰 指导老师：喻小虎

摘 要 本课程设计主要内容是设计利用脉冲响应不变法设计的巴特沃斯滤波器，对一段含噪语音信号进行滤波去噪处理并根据滤波前后的波形和频谱分析滤波性能。本课程设计仿真平台为MATLAB7.0，开发工具是M语言编程。录制一段语音信号，并人为加入一单频噪声，然后对信号进行频谱分析以确定所加噪声频率，即 fn=3000Hz，并设计滤波器进行滤波去噪处理，最后比较滤波前后的波形和频谱，得出结果为，滤波器后的语音信号与原始信号基本一致，即设计的巴特沃斯滤波器能够去除信号中所加单频噪声，达到了设计目的。

关键词 课程设计；滤波去噪；巴特沃斯滤波器；脉冲响应不变法；MATLAB

1 引 言

本课程设计需要录制一段频率为8000Hz，8位的单声道语音信号，并绘制波形观察其频谱，再对其进行加噪处理，同样要绘制加噪后的频谱图，再用MATLAB利用脉冲响应不变法设计一个巴特沃斯滤波器，将该语音信号进行滤波去噪处理。

1.1 课程设计目的

课程设计是教学的最后一个步骤，课程设计有利于基础知识的理解，我们掌握了基础知识和基本技能，但是要真正接触才能真正理解课程的深入部分；还有利于逻辑思维的锻炼，在许多常规学科的日常教学中，我们不难发现这样一个现象，不少学生的思维常常处于混乱的状态，写起作文来前言不搭后语，解起数学题来步骤混乱，这些都是缺乏思维训练的结果，所以我们可以通过实践来分析问题、解决问题、预测目标等目的；同时也有利于与其他学科的整合，例如我们这次的课程设计就要运用MATLAB软件的

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第2页 共20页

帮助才能实现；最重要的有利于治学态度的培养，在课程设计中，我们可能经常犯很多小错误，可能要通过好几次的反复修改、调试才能成功，但这种现象会随着学校的深入而慢慢改观。这当中就有一个严谨治学、一丝不苟的科学精神的培养，又有一个不怕失败、百折不饶品格的锻炼。

《数字信号处理》课程设计是在学生完成数字信号处理和MATLAB的结合后的基本实验以后开设的。本课程设计的目的是为了让学生综合数字信号处理和MATLAB并实现一个较为完整的小型滤波系统。这一点与验证性的基本实验有本质性的区别。开设课程设计环节的主要目的是通过系统设计、软件仿真、程序安排与调试、写实习报告等步骤，使学生初步掌握工程设计的具体步骤和方法，提高分析问题和解决问题的能力，提高实际应用水平。

1.2 课程设计的要求

（1）滤波器指标必须符合工程设计。

（2）设计完后应检查其频率响应曲线是否满足指标。 （3）处理结果和分析结论应该一致，而且应符合理论。 （4）独立完成课程设计并按要求编写课程设计报告。

1.3 设计平台

MATLAB是美国MathWorks公司生产的一个为科学和工程计算专门设计的交互式大型软件，是一个可以完成各种精确计算和数据处理的、可视化的、强大的计算工具。它集图示和精确计算于一身，在应用数学、物理、化工、机电工程、医学、金融和其他需要进行复杂数值计算的领域得到了广泛应用。它不仅是一个在各类工程设计中便于使用的计算工具，而且也是一个在数学、数值分析和工程计算等课程教学中的优秀的数学工具，在世界各地的高和大型计算机上运行，适用于Windows、UNIX等多种系统平台。

MATLAB作为一种科学计算的高级语言之所以受欢迎，就是因为它有丰富的函数资源和工具箱资源，编程人员可以根据自己的需要选择函数，而无需再去编写大量繁琐的程序代码，从而减轻了编程人员的工作负担，被称为第四代编程语言的MATLAB最大的特点就是简洁开放的程序代码和直观实用的开发环境[1]。

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第3页 共20页

2 设计原理

本课程设计需要录制的一段频率为8000Hz，8位的单声道语音信号，并绘制波形观察其频谱，再对其进行加噪处理，同样要绘制加噪后的频谱图，再用MATLAB利用脉冲响应不变法设计一个巴特沃斯滤波器，将该语音信号进行滤波去噪处理。

2.1 IIR滤波器

IIR滤波器具有无限长脉冲响应，因此能够与模拟滤波器相匹敌；一般来说，所有的模拟滤波器都有无限长脉冲响应。因此，IIR滤波器设计的基本方法是利用复值映射将大家熟知的模拟滤波器变换为数字滤波器。这一方法的优势在于各种模拟滤波器设计（AFD）表格和映射在文献中普遍都能获得。这种IIR滤波器设计的基本方法存在两种途径：

途径1： 设计模拟低通 滤波器

途径2： 设计模拟低 通滤波器

实行频带变换 S → S 实行滤波器 变换S→Z 期望的IIR滤波器

实行滤波器变 换S→Z 实行频带变 换Z→Z 期望的IIR滤波器

我们将研究途径二的设计基本方法，在这种IIR滤波器设计方法将按下列步骤进行： （1）设计模拟低通滤波器。

（2）研究并实行滤波器变换以得到数字低通滤波器。

（3）研究并实行频带变换以便从数字低通滤波器得到其他数字滤波器[2]。

2.2 巴特沃斯滤波器

巴特沃斯滤波器是用上述IIR滤波器的性质表征的：它的幅度响应在通带和阻带都是平坦的。一个N阶低通滤波器的幅度平方响应给出为：

（2-1）

式中，N为整数，是滤波器的阶，?c(rad/s)是截止频率。称为滤波器的阶数，N越大，

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第4页 共20页

通带和阻带的近似性越好，过渡带也越陡。幅度平方响应的图如下所示。

1.21N=100 0.8N=1 0.6N=2 N=200 0.40.200102030405060708090100 图2.1 幅度平方响应图

从这张图可看出下面几个性质：

（1）??0,Ha(j0)2?1=1，对全部N。 （2）???c,Ha(j?c)?2212，对全部N，这意味着在?c有3dB衰减。

（3）Ha(j?)是?的单调下降的函数。

（4）Ha(j?)随N??向一个理想低通滤波器趋近。

（5）Ha(j?)在??0世最大平坦，因为在这里所有阶的导数存在且等于零。 为了确定系统函数Ha(s)现将得到

Ha(s)Ha(?s)?12N22

?s?1???j?c??(j?c)2N(j?c)2N?2N?2Ns?(j?c)2N?(s?pk)k?1（2-2）

由（2-2）式分母多项式的根（或Ha(s)Ha(?s)）的极点）给出为

pk?(?1)2N(j?c)??ce

1j2?N(2k?N?1)??ce?12k?1?j?????22N?,k?0,1,,2N?1 （2-3）

（2-3）式的解释是：

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第5页 共20页

（1）2N个极点均匀分布在半径为Ωc的圆上，相隔 pi/N （2）N 为奇数时， pk= Ωc ej2pik/N （3）N 为偶数时， pk= Ωc ej(pi/2N+kpi/N) （4）极点对虚轴(jΩ)对称分布

（5）极点永远不会落在虚轴上，且仅当N为奇数时才会落在实轴上 （6）选取极点全部落在左半平面可得到因果稳定的Ha(s) 作为一个例子，三阶和四阶巴特沃斯滤波器的极点分布如图所示：

图2.2巴特沃斯滤波器的极点图

通过选取在左半面的极点就能给出一个稳定和因果的Ha(s)，并且能将Ha(s)写成：

N?cHa(s)?

?(s?s)kk?1N（2-4）

2.3 脉冲响应不变法

脉冲响应不变法是从滤波器的脉冲响应出发，使数字滤波器的单位脉冲响应序h(n)模仿模拟滤波去的冲击响应ha(t)，使h(n)正好等于ha(t)的采样值，即

h(n)=ha(nT) （2-5）

T为采样周期。如以Ha(s)及H(z)分别表示ha(t)的拉式变换及h(n)的z变换。

变换前：Ha(s)??k?1NRks?pk （2-6）

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第6页 共20页

变换后：H(z)?Rk ?pkT?11?ezk?1 （2-7）

N由于z=ejw是在单位圆上，s=jj?是在虚轴上，所以有下面从s平面到z平面的变换 z=esT

（2-8）

系统函数H（z）和Ha（s）是经过频域混叠公式联系的：

1?2???H(z)|z?eST??Ha?s?jm?Tm???T? （2-9）

在（2-8）式的映射关系下，复平面的变换如图2.3

图2.3脉冲响应不变法中的复平面映射

映射特征：

（1）记σ = Re(s):σ < 0, 映射到 |z|<1 (单位圆内)，σ = 0,映射到 |z|=1 (单位圆上)。σ >0,映射到 |z|>1 (单位圆外)。

（2）多个 s 映射到一个 z : 多点对一点的映射，每个宽度为2π/T 的全部半无限带都映射到单位圆内。

（3）因果稳定的AF仍然映射为因果稳定的DF。 （4）如果滤波器不是真正带限会产生混叠误差(Aliasing )。

3.设计步骤

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第7页 共20页

3.1设计流程图

语音信号滤波去噪——使用脉冲响应不变法设计的巴特沃斯滤波器，设计流程图如下3.1所示：

开始 1、语音信号信号的采集（wavread函数），画时域图

快速傅里叶变换，并且画频谱图 2、加入3000Hz的噪声 画出加噪后频谱图 3、设定滤波器性能指标，通带截止频率 fp=2600，阻带截止频率fs=2950，通带波纹Rp=3，阻带波纹As=16 脉冲响应不变法法设计巴特沃斯滤波器 验证并进行频谱分析 4、设计好的滤波器进行滤波处 比较滤波前后语音信号的波形及频 回放语音信号 结束 图3.1 脉冲响应不变法巴特沃斯滤波器对语音信号去噪流程图

3.2语言信号的采集

人工录制一段为.wav的语音信号，时间在8-9s左右，要求为8000Hz，8位单声道的

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第8页 共20页

音频格式。然后在Matlab软件平台下，利用函数wavread对语音信号进行采样，函数为[x,fs,bits]=wavread('music.wav'),记住采样频率fs和采样点数。

所下载音乐的属性设置如图3.2：

图3.2 语音信号参数设置

3.3语音信号的频谱分析

首先使用MATLAB分别画出原始语音信号的时域波形，再对原始语音信号进行快速傅里叶变换，画出原始信号的频率——幅度谱，然后再原始语音信号的基础上加入单频噪声后，噪声频率为fn=3000Hz,再所得到的音乐做时间——幅度谱，同样对它做傅里叶变换，得到此信号的频率——频谱特性。可得到图形3.3和图3.4。

其中傅里叶变换表示为：X=abs(fft(x)); Y=abs(fft(y))

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第9页 共20页

图3.3原始信号与加噪后的程序运行图

由上图可知，单频噪声频率fn=3000Hz，此数据是在对原始语音信号的频率——幅度谱得出，此语音信号的主频带在2000Hz内，我们取3000Hz是为了达到更好的干扰效果。因此带上耳塞后可以听到原始语音和加噪后的带有尖锐单频啸叫声的语音。可以从听觉上感知噪声的存在，经过试验当fn取2000Hz以内的频率时，听到的干扰信号无明显的尖锐单频啸叫声。

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第10页 共20页

图3.4原始信号与加噪后的波形图与频谱图

由图我们可以看出，在加入噪声后，原始语音信号和加噪后有不同的时间——幅度谱和频率——幅度谱，证明噪声对原始信号造成了影响。

3.4 滤波器设计

将数字滤波器的设计指标设为通带截止频率fp=2600Hz,阻带频率fc=2950Hz,通带波纹Rp=3dB，阻带波纹As=16dB,要求确定H(z)。设计步骤如下：

（1） 选取T=1。

（2） 利用设参数?p，?s，Rp和As设计一个模拟滤波器Ha(s)。 （3） 利用部分分式将Ha(s)展开。

（4） 现在将模拟极点{pk}变换为数字极点{epkT}，得到数字滤波器。并作化简得出作为z-1有理函数的H(z)。

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第11页 共20页

根据上述步骤，按要求我们用脉冲响应不变法设计好了巴特沃斯滤波器。 根据原理课设计滤波器，可得到如下图3.5：

3.5 利用脉冲响应不变法设计的数字巴特沃斯低通滤波器程序运行图

由图可知，设置fn=3000Hz，即确定了通带截止频率为，fp=2600Hz，阻带截止频率为，fc=2950Hz，过渡带宽为450Hz，时间T=1s，为过渡带宽和时间T的值对于滤波器的设计十分重要。

[cs,ds]=afd_butt(OmegaP,OmegaS,Rp,As); *** Butterworth Filter Order = 15 由此句可得出：滤波器为16阶。 >> wp/pi

ans = 0.6500 >> ws/pi

ans =

0.7375

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第12页 共20页

图3.6 利用脉冲响应不变法设计的数字巴特沃斯滤波器

由滤波器幅度响应图可以看出，设计的通带截止频率wp/pi=0.6500,阻带截止频率ws/pi=0.7375，纵坐标差不多刚好在As=16db处，所以设计的滤波器达到要求。

3.5 信号的滤波

前面已经用脉冲响应不变法设计好了我们要的巴特沃斯低通滤波器，接着就对语音信号进行滤波处理，看自己设计的巴特沃斯滤波器有没有对我们的加噪后的语音信号进行处理。所以就用filter函数进行滤波，即y_fil=filter(b,a,y);。我们将滤波前后的时域波形进行比较，并对其进行快速傅里叶变换，即Y_fil=abs(fft(y_fil));，目的是对比前后的频域频谱，具体分析设计的滤波器是否达到设计要求，再将滤波后的信号与原始信号进行对比，可看出滤波效果。如图3.7是原始语音信号、加干扰后的信号和滤波后信号的时间——幅度谱与频率谱图：

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第13页 共20页

图3.7原始信号、加噪后信号和滤波后信号的波形图和频谱图

观察分析图3.7，由原始信号与加噪后的频率——幅度谱分析可看出在fn=3000Hz的频率上，加了单频噪声，由于加入噪声，对原始信号造成了干扰，在运行sound(y,fs,bits);后，可以明显听出有尖锐的单频啸叫声，对加噪后的信号和滤波后的信号进行频率——幅度谱分析，可得出，设计的巴特沃斯低通滤波器滤除了大部分的噪声，运行sound(y_fil,fs,bits);后，原来的尖锐声微乎其微，由此证明了设计的合理性。最后对原始和滤波后的频率——幅度谱进行分析，滤波后的信号基本与原始信号无异，因此滤波成功。

3.6 结果分析

首先采集语音信号，再按照步骤用脉冲响应不变法设计巴特沃斯滤波器，得到图3.4。并且由图3.4可知，设计的通带截止频率wp/pi=0.6500，阻带截止频率ws/pi=0.7375，纵坐标在As=16db之下，所以设计的滤波器达到要求。我们观察到图3.7滤波前后语音

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第14页 共20页

信号的波形对比图，发现时域波形中的变化不明显，可能是因为我们采集的语音信号噪声不是很大，但是还是有滤去噪声的；但是可以看到在频域波形中，很明显地反应出设计的滤波器滤去了我们采集的语音信号中的噪声。所以，运用脉冲响应不变法设计的巴特沃斯滤波器达到了设计要求。

带上耳塞，运行程序，第一段为原音，第二段为加噪后的音乐，第三段为滤波后的音乐，通过反复听这三段音乐，能够听出第一段音乐和第二段音乐确有不同，能感觉到噪声，第二段与第三段相对比，几乎听不到噪声，第一段与第三段相比，几乎无差距。

图3.8，与相同滤波器（脉冲响应不变法的巴特沃斯低通滤波器）不同音乐（人工语音），通过观察自己的图形与同学的设计结果图，在交流中发现，我们如果设置相同的fn=3000Hz,我们滤波器阶数相同，都为15阶，但是我们的频率——幅度谱的幅度值有差异，但是滤波效果是相同的。

图3.8语音信号的脉冲响应不变法巴特沃斯滤波去噪图

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第15页 共20页

4 出现的问题及解决方法

在这次的课程设计中我们设计的滤波器对语音信号的滤噪处理或多或少会因为我们的疏忽或者实践能力不够还是会遇到问题的，当遇到问题时我们不但运用自己学到的知识还得到同学和老师的帮助将其解决，不断完善和修改。设计中出现的问题及解决方法如下：

（1）设计的滤波器达不到要求，因为在采集语音信号时，录音机的音频格式应该是PCM 8000Hz，8位，单声道。

（2）在运用Matlab设计滤波器时，当编辑完前面两条程序时无法放出声音，后来发现我们应当把采集的语音信号wav文件放到Matlab的work文件夹中，如果要放在其他文件夹在调用此音乐时需标明路径，如： [x,fs,bits]=wavread('d:\\music.wav'); 则标明名字为music的音乐信号时放在d盘。

（3）还要在滤波器性能曲线的wp、ws、Ap处画一根竖线，这样更方便看出结果，其中通带截止频率处线的值由wp/pi，阻带截止频率处的值为ws/pi。

（4）所有的时间波形横坐标都要化为时间，滤波前后频谱的横坐标应是频率，这样在观察通带截止频率和阻带截止频率时更加精确，误差较小。

（5）由于通带截止频率fp=fn-400和阻带截止频率fc=fn-50的值会影响滤波器的过渡带宽，而过渡带宽如果太小会造成滤波器无法达到指标，所以适当的在已知fn的条件下，基本保证过度带宽在0.1*pi附近。

（6）时间T也会影响滤波器的设计，开始时我设定T=1/fs，fs=8000Hz，此时采样频率过大，造成了音乐信号放大，经改正将T设为1s，此时就达到了预知的效果。

（7）wp=fp/fs*2*pi;ws=fc/fs*2*pi; OmegaP=wp/T; OmegaS=ws/T; 在用脉冲响应不变法设计低通滤波器时，这些公式的错误运用也会对滤波器设计造成严重的影响。

（8）对于设计的滤波器我们如何检验其是否可用，我们要对滤波器的频率——幅度谱进行分析，在ws/pi=0.7375与图中曲线的交接处以下的值是否大于As=16db，如果大于，则此滤波器设计完成。

（9）在观察中我们发现，时间——幅度谱对于一些细节的变化并不明显，因此我们可以从频率——幅度谱着手分析。

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第16页 共20页

5 结束语

在这次的课程设计中，我们用一周的时间来做设计，其余时间是自己写与修改报告， 在开始前我们需要复习课本上滤波器设计方面的知识，后来在设计中我们又发现了很多的问题，比如，在课本学习过程中我们以为不太重要的环节在实现滤波器效果时，其实是很重要的，任何一个小小的错误，比如一个参数设置的合理性都会造成无法设计出符合标准的滤波器。并且在这次DSP课程设计，不但巩固了以前学的知识，还学习了新内容，在实践中收获知识。原来滤去语音中的噪声我们运用自己学到的知识能做到，并更加巩固了MATLAB软件的使用。

MATLAB软件的熟练使用是很重要的，整个课程设计过程中，我们都要用到它。当然由于我们的实践能力有限，在过程中，我们需要老师的帮助，在老师的帮助下问题就迎刃而解啦。

这次做的滤波器要滤去语音信号中的噪声，觉得很有意思，做了1次的课程设计了，我发现其实在这样的课程设计中一些细微之处往往是决定成败的关键。其实用MATLAB软件做实验是要细心的，因为很多的语法和常量变量的定义我们都要仔细，一个不小心看错了或者输入不认真是容易出错误。在发现错误时，问了很多同学也问了老师，当发现错误原因时，很多次都是因为自己的不细心造成。但经过思考并仔细检查，最后得以编译成功。在做实验是一定要小心谨慎，也许一个小小的语法错误和常量变量的定义的错误就造成整个程序出现问题，得不到所需的波形，导致实验结果不正确。

实验中，同学间的相互交流也是很重要的，比如，一个同学是跟我用相同音乐的，但是她课题为“语音信号滤波去噪——使用脉冲响应不变法设计的切比雪夫I型滤波器”我们只是滤波器不同，在结束后我们彼此交流了两个不同滤波器的性能，明显能够看出切比雪夫I型滤波器滤波效果更好。在这样的环境下，我们能够收获更多。

这次实验能够顺利完成，有自己的努力，也有同学和老师的帮助，因此，在这里我要衷心地感谢他们的帮助，也感谢学校给我们的机会，让我们能够将自己学到的知识运用到实际中！

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第17页 共20页

参考文献

[1] 张圣勤．MATLAB7.0实用教程〔M〕．北京：机械工程出版社，2006

[2] (美)维纳?K?英格尔，（美）约翰?G?普罗克斯（著）．刘树棠（译）．.数字信号处理（MATLAB版）〔M〕．第2版．西安：西安交通大学出版社，2008 [3] 程佩青．数字信号处理教程〔M〕．北京：清华大学出版社，2002 [4] 刘敏，魏玲．Matlab通信仿真与应用〔M〕．北京：国防工业出版社，2001 [5] Sanjit K.Mitra，孙洪，余翔宇．数字信号处理实验指导书〔M〕．电子工业出版社，2005

[6] 陈怀琛．数字信号处理教程－MATLAB释义与实现〔M〕．电子工业出版社，2008 [7] 丁玉美, 高西全．数字信号处理〔M〕．西安：电子科技大学出版社，2010

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第18页 共20页

附录1：语音信号滤波去噪设计源程序清单

% 程序名称：untitled.m

% 程序功能：采用基于脉冲响应不变法设计法，设计巴特沃斯滤波器对含噪语音信号进行滤波去噪处理。 % 程序作者：王雄杰 /% 最后修改日期：2015-3-27

[x,fs,bits]=wavread('d:music.wav');% 输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。

sound(x,fs,bits); % 按指定的采样率和每样本编码位数回放 N=length(x); % 计算信号x的长度

t=0:1/fs:(N-1)/fs;% 计算时间范围，样本数除以采样频率 deltaf=fs/N; % 计算频谱的谱线间隔 f=0:deltaf:fs/2-deltaf; % 计算频谱频率范围

X=abs(fft(x)); % 对原始信号和加噪信号进行fft变换,取幅度谱 X=X(1:N/2);% 截取前半部分

subplot(2,2,1) ;plot(t,x) %画出原始语音信号的时间幅度谱图 xlabel('时间（单位：s）');ylabel('幅度');title('原始语音信号') subplot(2,2,2) ;plot(f,X) %画出原始信号的频率幅度谱图

xlabel('频率（单位：HZ）');ylabel('幅度谱');title('原始语音信号幅度谱图') fn=3000;% 单频噪声频率

x=x'; y=x+0.02*sin(fn*2*pi*t); %在原始语音信号基础上加单频噪声后的每个样本的值 Y=abs(fft(y));% 对原始信号和加噪信号进行fft变换,取幅度谱 Y=Y(1:N/2); %截取前半部分

subplot(2,2,3) ;plot(t,y) %画出加入单频干扰后的的时间幅度谱图 axis([0 10 -1 1])% 横坐标从0到10 ，纵坐标从-1到1

xlabel('时间（单位：s）');ylabel('幅度');title('加入单频干扰后的语音信号') subplot(2,2,4) ;plot(f,Y)%画出加入单频干扰后的频率幅度谱图

xlabel('频率（单位：HZ）');ylabel('幅度谱');title('加入单频干扰后的语音信号幅度谱图')

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第19页 共20页

sound(y,fs,bits); % 可以明显听出有尖锐的单频啸叫声 % 程序名称：filter.m

fp=fn-400;fc=fn-50; %定义通带和阻带截止频率 T=1; %定义采样间隔

Rp=3;As=16; % 定义通带波纹和阻带衰减

wp=fp/fs*2*pi;ws=fc/fs*2*pi; %计算对应的数字频率 OmegaP=wp/T; OmegaS=ws/T; %截止频率线性变换

[cs,ds]=afd_butt(OmegaP,OmegaS,Rp,As); %计算滤波器阶数和截止频率和滤波器系统函数分子分母系数

[b,a]=imp_invr(cs,ds,T);% 脉冲响应不变法得到数字滤波器系统函数分子分母系数 [db,mag,pha,grd,w]=freqz_m(b,a); %验证滤波器是否达到指定性能 delta=[1,zeros(1,99)];ha=filter(b,a,delta);% 计算脉冲响应 Subplot(221);plot(w/pi,db)

xlabel('w/pi');ylabel('db');title('滤波器幅度响应图') Subplot(222);plot(w/pi,mag)

xlabel('w/pi');ylabel('幅度mag');title('滤波器幅度响应图') Subplot(223);plot(w/pi,pha)

xlabel('w/pi');ylabel('相位pha');title('滤波器相位响应图') Subplot(224);plot(ha)

xlabel('n');ylabel('h(n)');title('滤波器脉冲响应图') % 程序名称：voisenoisefilter.m

[x,fs,bits]=wavread('d:music.wav');% 输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。

sound(x,fs,bits); % 按指定的采样率和每样本编码位数回放 N=length(x); % 计算信号x的长度

t=0:1/fs:(N-1)/fs;% 计算时间范围，样本数除以采样频率 deltaf=fs/N; % 计算频谱的谱线间隔 f=0:deltaf:fs/2-deltaf; % 计算频谱频率范围

X=abs(fft(x)); % 对原始信号和加噪信号进行fft变换,取幅度谱

王雄杰《语音信号滤波去噪--使用脉冲响应不变法设计的巴特沃斯滤波器》第20页 共20页

X=X(1:N/2);% 截取前半部分

subplot(3,2,1) ;plot(t,x) %画出原始语音信号的时间幅度谱图 xlabel('时间t');ylabel('幅度');title('原始语音信号时间x') subplot(3,2,2) ;plot(f,X) %画出原始信号的频率幅度谱图 xlabel('频率f');ylabel('幅度谱');title('原语音信号幅度谱图X') fn=3000;% 单频噪声频率

x=x'; y=x+0.02*sin(fn*2*pi*t); %在原始语音信号基础上加单频噪声后的每个样本的值 Y=abs(fft(y));%对原始信号和加噪信号进行fft变换,取幅度谱 Y=Y(1:N/2); % 截取前半部分

subplot(3,2,3) ;plot(t,y) %画出加入单频干扰后的的时间幅度谱图 xlabel('时间t');ylabel('幅度');title('加入单频干扰后的语音信号时间y') subplot(3,2,4) ;plot(f,Y) %画出加入单频干扰后的频率幅度谱图 axis([0 10 -1 1]) %横坐标从0到10 ，纵坐标从-1到1

xlabel('频率f');ylabel('幅度谱');title('加入单频干扰后的语音信号幅度谱图Y') sound(y,fs,bits); % 可以明显听出有尖锐的单频啸叫声

y_fil=filter(b*T,a*T,y); % IIR滤波器对信号signal进行滤波处理，其中b,a为上面设计好的滤波器参数。

Y_fil=abs(fft(y_fil));Y_fil=Y_fil(1:N/2);% 计算频谱取前一半 subplot(3,2,5) ;plot(t,y_fil) %画出滤波后语音信号时间幅度谱图 axis([0 10 -1 1]) %横坐标从0到10 ，纵坐标从-1到1 xlabel('时间t');ylabel('幅度');title('滤波后语音信号时间y_fil') subplot(3,2,6) ;plot(f,Y_fil)%画出滤波后语音信号频率幅度谱图 xlabel('频率f');ylabel('幅度');title('滤波后语音信号幅度谱Y_fil') sound(y_fil，fs, bits);

本文来源：https://www.bwwdw.com/article/6q2w.html