语音信号降噪分析

中北大学

课程设计任务书

课程设计题目： 信息处理实践： 语音信号降噪分析

下达任务书日期: 2010 年6 月 6 日

课 程 设 计 任 务 书

1．设计目的： 采集一段语音信号，分析该信号的时域和频域特征，了解噪声在语音信号中的表现特征以及对语音信号的影响，并运用有关滤波器对该信号中的噪声进行去除。这培养学生根据实际信号，采用合适的方法对信号进行处理的能力。 2．设计内容和要求（包括原始数据、技术参数、条件、设计要求等）： 1、用MATLAB对原始语音信号进行分析，画出它的时域波形和频谱； 2、掌握语音信号与噪声的特点； 3、对滤除噪声后的语音信号进行时域和频谱分析； 4、比较滤波前后的波形。 3．设计工作任务及工作量的要求〔包括课程设计计算说明书(论文)、图纸、实物样品等〕： 课程设计说明书一份； 仿真结果及分析。 课 程 设 计 任 务 书

4．主要参考文献： [1] 张平.matlab基础与应用. 第2版. 北京：北京航空航天大学出版社，2007 [2] 谢德芳. 数字信号处理. 第1版.北京：科学出版社，2005 [3] 陈后金，薛健，胡健. 数字信号处理. 第1版. 北京：高等教育出版社，2004 [4] 维纳K英格尔，约翰G普罗克斯. 刘树棠. 数字信号处理（MATLAB版）. 第2版. 西安：西安交通大学出版社，2008 [5] 张小虹. 信号系统与数字信号处理. 第1版. 西安：西安电子科技出版社，2002 [6] 蔡启仲. 控制系统计算机辅助设计（MATLAB版）. 第1版. 四川：重庆大学出版社， 5．设计成果形式及要求： 课程设计说明书 仿真结果 6．工作计划及进度： 2010年 6月6日 ~ 6月10日：查资料； 6月11日 ~ 6月16 日：在指导教师指导下设计方案； 6月17日 ~ 6月24日：在指导教师辅导下完成实验；撰写课程设计说明书； 6月25日：答辩 系主任审查意见： 签字： 年 月 日

目 录

封面…………………………………………………….. 设计任务书说明……………………………………….. 一．引言……………………………………………..… （1）课程设计目的…………………………………… （2）课程设计要求…………………………………… （3）设计平台………………………………………… 二．设计原理………………………………………….. （1）FIR滤波器………………………………………. （2）窗口设计法………………………………………. （3）PARZENWIN窗………………………………… 三、设计步骤…………………………………………… （1）设计流程图………………………………………

（2）录制语音信号…………………………………… （3）滤波器设计……………………………………… （4）信号滤波处理…………………………………… （5）结果分析……………………………………………

四、出现问题及解决办法……………………………

五、实验体会……………………………………… 参考文献………………………………………………

1 引 言

本课程设计是采用parzenwin窗设计的FIR滤波器对语音信号滤波去噪。通过课程设计了解FIR滤波器的原理及使用方法，了解使用Matlab语言设计FIR滤波器的方法，了解DSP对FIR滤波器的设计及编程方法。通过观察滤波前后的时域图形，加深对FIR滤波器作用的理解。通过对比滤波前后的波形图及回放滤波前后的语音信号，可以看出滤波器对有用信号的无失真放大具有重要作用。

1.1 课程设计目的

利用Matlab环境下的M文件，用parzenwin窗设计的FIR滤波器来实现对音乐信号去噪，并绘制出滤波前后的时域和频域波形及播放声音的变化，根据运行结果和波形来分析设计过程的正确性。通过这次课程设计，加深对parzenwin窗设计的FIR滤波器的理解，掌握Matlab软件在滤波器设计中的应用，锻炼逻辑思维能力、动手能力以及独立解决问题的能力，对以后更深入地学习和应用数字信号处理及相关知识作准备。

1.2 课程设计的要求

（1）滤波器指标必须符合工程实际。

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

1.3 设计平台

MATLAB是美国Math works公司开发的新一代科学计算软件，是一套高性能的数值计算和可视化软件，功能强大，编程简单，开放性强，广泛应用于计算机辅助分析设计、仿真、数据处理等领域，是当今国际上公认的在科技领域方面最为优秀的应用软件和开发环境。在欧美各高等院校，已经成为应用线性代数、自动控制理论、数据统计、数字信号处理、时间序列分析、动态系统仿真、图像处理等高级课程的基本教学工具。MATLAB 是一种既可交互使用又能解释执行的计算机编程语言，使用接近数学表达式的自然化语言，简单易学，具有可靠的数值、符号运算能力。和强大的图形和可视化功能；此外，MATLAB 内部包括许多专业性较强的工具包并与其他高级语言有接口。MATLAB 的功能和特点使它具备了对应用学科(特别是边缘学科和交叉学科)的极强适应力，并很快成为应用学科计算机辅助分析、设计、仿真、数学乃至科技文字处理不可缺少的基础软件[6]。

2 设计原理

2.1 FIR滤波器

根据冲激响应的时域特性，数字滤波器可分为无限长冲激响应（IIR）和有限

长冲激响应滤波器（FIR），FIR的突出优点是：系统总是稳定的、易于实现线性相位、允许设计多通带（或多阻带）滤波器，但与IIR相比，在满足同样阻带衰减的情况下需要的阶数较高，滤波器的阶数越高，占用的运算时间越多，因此在

满足指标要求的情况下应尽量减少滤波器的阶数。

FIR滤波器的基本结构可以理解为一个分节的延时线，把每一节的输出加权累加，可得到滤波器的输出，FIR滤波器的冲激响应h（n）是有限长的，数学上M阶FIR滤波器可以表示为：

y(n)=

h(i)x(ni);H(z)h(n)zi0n0N1N1n

FIR数字滤波器设计的基本步骤如下：

（1）确定指标

在设计一个滤波器之前，必须首先根据工程实际的需要确定滤波器的技术指标。在很多实际应用中，数字滤波器常常被用来实现选频操作。因此，指标的形式一般在频域中给出幅度和相位响应。幅度指标主要以两种方式给出。第一种是绝对指标。它提供对幅度响应函数的要求，一般应用于FIR滤波器的设计。第二种指标是相对指标。它以分贝值的形式给出要求。在工程实际中，这种指标最受欢迎。对于相位响应指标形式，通常希望系统在通频带中人有线性相位。运用线性相位响应指标进行滤波器设计具有如下优点：①只包含实数算法，不涉及复数运算；②不存在延迟失真，只有固定数量的延迟；③长度为N的滤波器（阶数为N-1），计算量为N/2数量级。因此，本文中滤波器的设计就以线性相位FIR滤波器的设计为例。

（2）逼近

确定了技术指标后，就可以建立一个目标的数字滤波器模型。通常采用理想的数字滤波器模型。之后，利用数字滤波器的设计方法，设计出一个实际滤波器模型来逼近给定的目标。

（3）性能分析和计算机仿真

上两步的结果是得到以差分或系统函数或冲激响应描述的滤波器。根据这个描述就可以分析其频率特性和相位特性，以验证设计结果是否满足指标要求；或者利用计算机仿真实现设计的滤波器，再分析滤波结果来判断FIR滤波器的设计问题实质上是确定能满足所要求的转移序列或脉冲响应的常数的问题，设计方法主要有窗函数，频率采样法和等波纹最佳逼近法等。

2.2 窗口设计法

窗函数设计法是一种通过截短和计权的方法使无限长非因果序列成为有限长脉冲响应序列的设计方法。通常在设计滤波器之前，应该先根据具体的工程应用确定滤波器的技术指标。在大多数实际应用中，数字滤波器常常被用来实现选频操作，所以指标的形式一般为在频域中以分贝值给出的相对幅度响应和相位响应。

用窗函数法设计FIR滤波器的步骤如下：

（1）根据过渡带宽及阻带衰减要求，选择窗函数的类型并估计窗口长度N（或阶数M=N-1）。窗函数类型可根据最小阻带衰减AS独立选择，因为窗口长度N对最小阻带衰减AS没有影响。在确定窗函数类型以后，可根据过渡带宽小于给定指标确定所拟用的窗函数的窗口长度N。设待求滤波器的过渡带宽为△ω，它与窗口长度N近似成反比。窗函数类型确定后，其计算公式也确定了，不过这些公式是近似的，得出的窗口长度还要在计算中逐步修正。原则是在保证阻带衰减满足要求的情况下，尽量选择较小的N。在N和窗函数类型确定后，即可调用MATLAB中的窗函数求出窗函数wd(n)。

（2）根据待求滤波器的理想频率响应求出理想单位脉冲响应hd(n)。如果给出待求滤波器的频率响应为Hd(ejω)，则理想的单位脉冲响应可以用下面的傅里叶反变换式求出：

在一般情况下，hd(n)是不能用封闭公式表示的，需要采用数值方法表示。从ω=0到ω=2π采样N点，采用离散傅里叶反变换(IDFT)即可求出。

（3）计算滤波器的单位脉冲响应h(n)。它是理想单位脉冲响应和窗函数的乘积，即h(n)=hd(n)·wd(n)，在MATLAB中用点乘命令表示为h=hd·wd。 （4） 验算技术指标是否满足要求。为了计算数字滤波器在频域中的特性，可调用freqz子程序，如果不满足要求，可根据具体情况，调整窗函数类型或长度，直到满足要求为止

使用窗函数法设计时要满足以下两个条件： （1）窗谱主瓣尽可能地窄，以获得较陡的过渡带；

（2）尽量减少窗谱的最大旁瓣的相对幅度，也就是使能量尽量集中于主瓣，减小峰肩和纹波，进而增加阻带的衰减。在实际工程中常用的窗函数有五种，即矩形窗、三角窗、汉宁窗、海明窗和凯泽窗。 实际应用的窗函数，可分为以下主要类型：

（1）幂窗--采用时间变量某种幂次的函数，如矩形、三角形、梯形或其它时间（t）的高次幂；

（2）三角函数窗--应用三角函数，即正弦或余弦函数等组合成复合函数，例如汉宁窗、海明窗等；

（3）指数窗--采用指数时间函数，如e-st形式，例如高斯窗等。

下面介绍几种常用窗函数的性质和特点： （1）矩形窗

矩形窗属于时间变量的零次幂窗。矩形窗使用最多，习惯上不加窗就是使信号通过了矩形窗。这种窗的优点是主瓣比较集中，缺点是旁瓣较高，并有负旁瓣，导致变换中带进了高频干扰和泄漏，甚至出现负谱现象。 （2）三角窗

三角窗亦称费杰（Fejer）窗，是幂窗的一次方形式。与矩形窗比较，主瓣宽约等于矩形窗的两倍，但旁瓣小，而且无负旁瓣。 （3）汉宁（Hanning）窗

汉宁窗又称升余弦窗，汉宁窗可以看作是3个矩形时间窗的频谱之和，或者说是 3个 sine（t）型函数之和，而括号中的两项相对于第一个谱窗向左、右各移动了 π/T，从而使旁瓣互相抵消，消去高频干扰和漏能。可以看出，汉宁窗主瓣加宽并降低，旁瓣则显著减小，从减小泄漏观点出发，汉宁窗优于矩形窗．但汉宁窗主瓣加宽，相当于分析带宽加宽，频率分辨力下 （4）海明（Hamming）窗

海明窗与汉宁窗都是余弦窗，只是加权系数不同。海明窗加权的系数能使旁瓣达到更小。分析表明，海明窗的第一旁瓣衰减为一42dB．海明窗的频谱也是由3个矩形时窗的频谱合成，但其旁瓣衰减速度为20dB／（10oct），这比汉宁窗衰减速度慢。海明窗与汉宁窗都是很有用的窗函数。

（5）高斯窗

高斯窗是一种指数窗。高斯窗谱无负的旁瓣，第一旁瓣衰减达一55dB。高斯富谱的主瓣较宽，故而频率分辨力低．高斯窗函数常被用来截断一些非周期信号，如指数衰减信号等。

不同的窗函数对信号频谱的影响是不一样的，这主要是因为不同的窗函数，产生泄漏的大小不一样，频率分辨能力也不一样。信号的截断产生了能量泄漏，而用FFT算法计算频谱又产生了栅栏效应，从原理上讲这两种误差都是不能消除的，但是我们可以通过选择不同的窗函数对它们的影响进行抑制。

对于窗函数的选择，应考虑被分析信号的性质与处理要求。如果仅要求精确读出主瓣频率，而不考虑幅值精度，则可选用主瓣宽度比较窄而便于分辨的矩形窗，例如测量物体的自振频率等；如果分析窄带信号，且有较强的干扰噪声，则应选用旁瓣幅度小的窗函数，如汉宁窗、三角窗等；对于随时间按指数衰减函

数，可采用指数窗来提高信噪比

2.3 PARZENWIN窗

w = parzenwin(n) 返回N -点Parzen（日拉瓦莱-普桑）窗口中列向量。Parzenwin窗分段立方米近似高斯窗口。窗口旁瓣脱落的 1w4 。 Parzenwin窗的定义为：

nwn1.06N22nN1.0 0nN24nNNwn21.0 nN24233设计步骤

3.1 设计流程图

开始 录制一段语音信号：“大家好，我是阮正杰”, 命名为ruan.wav，绘制出其时域波形和频谱图 根据其频谱图计算出通带截止频率wp和阻带截止频率ws,设定指标As,Rp 自己计算parzenwin窗的过度带宽？pi/M,设计滤波器 将语音信号通过自己所设计的FIR滤波器，进行滤波去噪 比较滤波前后语音信号的时域波形和频谱图，并回放，验证是否达到去噪效果 是 结束

3.2 录制语音信号

用windows附件中的录音机录一段自己说的话（语音信号），设定其采样率fs为8000HZ，比特数8bits，单声道。语音为：“大家好，我是阮正杰”，时间为2.0s。将语音信号的文件名命名为ruan.wav，将语音文件保存后，在MATLAB软件平台下，首先调用wavread函数可采集到录制的音乐信号。具体调用如下： >> [y,fs,bits]=wavread('C:\\MATLAB7\\work\\ruan.wav')%读取语音信号

运行后得出fs=8000，bits=16。输出的第一个参数x是每个样本的值，fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。 >> sound(y,fs,bits); %播放语音信号 >> N=length(y); %计算语音信号的长度 >> t=[0:N-1]/fs;

>> plot(t,y); xlabel('t in s.');ylabel('y(t)');title('每个样本的值');grid; %画原始语音信号的时域图

用plot函数画出原始信号的时域图，在原始信号的时域图中应该能看出原始语音信号的大致内容。在对应的位置上应该有相应的比较明显的凸起，而且重音的凸起较大，那么经过用plot函数画出相应的图形后，我们得到的图形是运行后的波形如图3.1：

图3.2 原始语音信号时域波形图

那么由原始信号的时域波形图可以看出语音信号的大致内容。再对原始语音信号做快速傅立叶变换。实现程序如下： >> ft=8000/N; >> F=ft:ft:8000;

>> Y=fft(y); %对原始语音信号做快速傅立叶变化 >> magY=abs(Y);phaY=angle(Y);

>> subplot(2,1,1);plot(F,magY);title('幅度响应');grid; axis([0,8000,0,100]) %绘制原始语音信号的幅度响应

>> subplot(2,1,2);plot(F,phaY);title('相位响应');grid ; axis([0,8000,-4,4]) %绘制原始语音信号的相位响应运行后我们得到语音信号的幅度响应和相位响应如图3.3所示 ：

幅度响应10080604020001000200030004000相位响应420-2-40100020003000400050006000700080005000600070008000图3.3 原始语音信号的频谱图

3.3 滤波器设计

滤波器设计就是要找到一组能满足特定滤波要求的系数向量a和b,而它主要是通过设计指标来实现的。滤波器设计的要求或指标一般是在频域上给出的，常用的滤波器频域指标有：通带截止频率fp，阻带截止频率fs1，通带波纹Rp，阻带衰减As。要达到最佳的滤波效果，则需要对fp，fs1和As进行适当的调整。

因为采样频率f为8000Hz，由采样定理，fs1>=2fc,由图3.3的相位响应可知，取截止频率为2000Hz。所以语音信号可以选择fp= 1800Hz;fs1=2000Hz。设定指

标为：Rp=1,As=55

在MATLAB中，通常采用1/2采样频率进行归一化处理，如果将频率转化为角频率，则需将归一化频率乘以pi。设计的程序如下： >> fp=1800;fs1=2000;

>> wp=2*pi*fp/fs ;ws=2*pi*fs1/fs; %将模拟频率转换为数字角频率 >> tr_width=ws-wp; %过渡带宽

>> M=ceil(13.0725*pi/tr_width)+1%滤波器所用窗函数的最小长度 M =263

>> wc=(ws+wp)/2; >> hd=ideal_lp(wc,M);

>> w_parzen=(parzenwin(M))'; %采用parzenwin窗 >> h=hd.*w_parzen; %在时间域乘积对应于频率域的卷积 >> [db,mag,pha,grd,w]=freqz_m(h,1);delta_w=2*pi/1000; >> Rp=-(min(db(1:1:wp/delta_w+1))) )%计算通带波纹 Rp =0.0136

>> As=-round(max(db(ws/delta_w+1:1:501))) %计算阻带衰减 As = 56 >> b=h;

>> freqz(b,1,512,fs); %绘制滤波器的频率响应

通过运行上面一段程序，MATLAB设计出来的滤波器的真正参数为Rp =0.0136，As =56，符合指标。所设计的滤波器的频率响应如图3.4：

0Magnitude (dB)-50-100-1500x 1045001000150020002500Frequency (Hz)3000350040000Phase (degrees)-1-2-305001000150020002500Frequency (Hz)300035004000

图3.4 滤波器的频率响应

3.4 信号滤波处理

滤波器设计完成后，在MATLAB平台上用函数filter实现滤波，那么程序如下：

>> x=fftfilt(b,y); % FIR滤波器对语音信号进行滤波处理 >> sound(x,fs,bits) %回放滤波之后的语音信号

可以听出滤波之后的信号比原始语音信号更清晰一些，清除了环境噪声。 >> subplot(2,1,1);plot(t,y);title('原信号的时域波形') >> subplot(2,1,2);plot(t,x);title('滤波后的时域波形')

得到的滤波前后语音信号的时域波形对比图如图：

原信号的时域波形0.40.20-0.2-0.400.20.40.60.811.21.41.61.82滤波后的时域波形0.40.20-0.2-0.400.20.40.60.811.21.41.61.82图3.5 滤波前后的时域波形图

对滤波以后的语音信号作快速傅立叶变化，画出其相应的频谱图，实现的程序如下

>> X=fft(x); %对滤波后的语音信号做快速傅立叶变化 >> magX=abs(X);phaX=angle(X);

>> subplot(2,2,1);plot(F,magY);title('原信号的幅度响应') >> subplot(2,2,3);plot(F,magX);title('滤波后的幅度响应') >> subplot(2,2,2);plot(F,phaY);title('原信号的相位响应') >> subplot(2,2,4);plot(F,phaX);title('滤波后的相位响应')

得到滤波前后的语音信号的频谱图如图3.6：

原信号的幅度响应15042050-2-4原信号的相位响应10000200040006000800002000400060008000滤波后的幅度响应15010042050-202000400060008000-40滤波后的相位响应02000400060008000

图3.6 滤波前后的语音信号的频谱图

3.5 结果分析

在MATLAB中，经过sound(x,fs,bits)函数，对经过parzenwin窗设计的FIR滤波器之后的语音信号进行回放，可以听出滤波之后的信号比原始语音信号更清晰一些，清除了环境噪声，通过以下语句来进行语音信号回放比较： >> sound(y,fs,bits); >> sound(x,fs,bits);

所得结果，证明了用parzenwin窗设计的FIR滤波器和语音信号去噪设计是成功的。

4 出现问题及解决方法

在本次课程设计中我遇到了以下问题：

1、不知用什么软件去录制声音，上网搜索也没找的合适的软件。 2、在进行语音信号的提取时没设采样频率，结果频率没达到要求。 3、对提取的信号做频域图和时域图时效果不理想。

4、对parzenwin窗设计法很生疏，不知如何计算M的值。

5、在采用parzenwin窗设计的FIR滤波器时得不到理想的滤波器，把通带截止频率与阻带起始频率之间的差值设置的太小或者太大。

针对以上问题，相应的解决方案如下：

1、通过老师的指导，从windows自带的附件中找到录音软件。 2、在录音机软件中设定好采样频率的属性，使时域图达到要求。

3、通过多次录音的尝试，最终得到了理想的语音信号时域图和频域图的幅度响应。

4、通过同学的帮助，利用matlab编写关于parzenwin的M文件得到text5函数，再通过4个不同参数代入text5函数，画出图形并计算出过渡带宽为13.0725*pi/M，

可求出M。

5、通过适当的选择参数fp和fs1，绘制出来的图形效果比较明显，基本符合设计指标。最后得到As为56，这样通过MATLAB运算出来的滤波器的阻带波纹达到了滤波器设计的要求，得到了比较理想的滤波器。

5 实验体会

本次课程设计，我的任务是利用parzenwin窗设计的FIR滤波器对语音信号滤波去噪。开始我对parzenwin窗了解特别少，通过问懂得的同学，我掌握了用parzenwin窗设计FIR滤波器的方法，了解了窗函数的基本设计流程。经过几天忙碌的课程设计我体会到了很多。

首先通过这次课程设计使我明白了自己原来知识还比较欠缺。自己要学习的东西还太多，以前老是觉得自己什么东西都会，什么东西都懂，有点眼高手低。通过这次课程设计，我才明白学习是一个长期积累的过程，在以后的生活中都应该不断的学习，努力提高自己知识和综合素质。

其次让我明白谦虚的重要性：课程过程中遇到难题是在所难免的事情，当我们无法解决的时候我们应该问问其他的同学或老师，这样可以节省很多的时间。

最后我认识到理论运用到实践的重要性，正所谓“纸上得来终觉浅，绝知此事要躬行”。学习任何知识，仅从理论上去求知，而不去实践、探索是不够的。所以在学完数字信号处理之际，紧接着来一次数字信号处理的课程设计事很有必要的。这样不仅加深我们对数字信号处理的认识，而且还及时真正做到了学以致

用。

在此要感谢我们的指导老师胡老师对我们悉心的指导，感谢胡老师给我们的帮助。虽然这个设计做的也不太好，但是在设计过程中所学到的东西是这次课程设计的最大收获和财富，使我终身受益。

参考文献

[1] 张平.matlab基础与应用. 第2版. 北京：北京航空航天大学出版社，2007 [2] 谢德芳. 数字信号处理. 第1版.北京：科学出版社，2005

[3] 陈后金，薛健，胡健. 数字信号处理. 第1版. 北京：高等教育出版社，2004 [4] 维纳K英格尔，约翰G普罗克斯. 刘树棠. 数字信号处理（MATLAB版）. 第2版. 西安：西安交通大学出版社，2008

[5] 张小虹. 信号系统与数字信号处理. 第1版. 西安：西安电子科技出版社，2002

[6] 蔡启仲. 控制系统计算机辅助设计（MATLAB版）. 第1版. 四川：重庆大学出版社，2003