一、设计的目的和意义数字滤波器和快速傅立叶变换(FFT)等是语音信号数字处理的理论和技术基础，是20世纪60年代形成的一系列数字信号处理的理论和算法。

在数字信号处理中,滤波器的设计占有极其重要的地位。

而其中，FIR数字滤波器和IIR数字滤波器是重要组成部分。

Matlab具有功能强大、简单易学、编程效率高等特点,深受广大科技工作者的喜爱。

特别是Matlab中还具有信号分析工具箱,所以对于使用者，不需要具备很强的编程能力,就可以方便地进行信号分析、处理和设计。

利用Matlab中的信号处理工具箱，可以快速有效的设计各种数字滤波器。

本论文基于Matlab语音信号处理的设计与实现，综合运用数字信号处理的相关理论知识，对加噪声语音信号进行时域、频域分析并滤波。

而后通过理论推导得出相应结论，再利用Matlab作为编程工具进行计算机实现工作。

本次课程设计的课题为《基于DSP的语音信号滤波去噪》，运用麦克风采集一段语音信号，绘制波形并观察其频谱，给定相应技术指标，用脉冲响应不变法设计的一个满足指标的巴特沃斯IIR滤波器，对该语音信号进行滤波去噪处理，比较滤波前后的波形和频谱并进行分析，根据结果和学过的理论得出合理的结论。

二、设计原理：2.1 巴特沃斯滤波器巴特沃斯滤波器是电子滤波器的一种。

巴特沃斯滤波器的特点是通频带的频率响应曲线最平滑。

巴特沃斯滤波器的特性是通频带内的频率响应曲线最大限度平坦，没有起伏，而在组频带则逐渐下降为零。

在振幅的对数对角频率的波得图上，从某一边界角频率开始，振幅随着角频率的增加而逐步减少，趋向负无穷大。

其振幅平方函数具有如2-1式：（2-1）式中，N为整数，称为滤波器的阶数，N越大，通带和阻带的近似性越好，过渡带也越陡。

如下图2.1所示：图2.1 巴特沃兹filter 振幅平方函数过渡带：通带→阻带间过渡的频率范围，Ωc：截止频率。

理想滤波器的过渡带为Ω，阻带|H(jΩ)|=0，通带内幅度|H(jΩ)|=常数，H(jΩ)线性相位。

通带内，分母Ω/Ωc<1,相应(Ω/Ωc)2N随N的增加而趋于0，A(Ω2)→1，在过渡带和阻带，Ω/Ωc >1，随N的增加，Ωe/Ωc>>1，所以A(Ω2)快速下降。

Ω=Ωc时，，幅度衰减，相当于3bd衰减点。

振幅平方函数的极点可写成如式2-2：Ha(-s).Ha(s)= （2-2）可分解为2N个一次因式令分母为零，→可见，Butterworth 滤波器的振幅平方函数有2N个极点，它们均匀对称地分布在|s|=Ωc的圆周上。

2.2 脉冲响应不变法如果从模拟到数字滤波器我们想要保留脉冲响应的形状，那么就得到一种方法称为脉冲不变响应法的变换方法。

脉冲响应不变法是从滤波器的脉冲响应出发，使数字滤波器的单位脉冲响应序h(n)模仿模拟滤波去的冲击响应ha(t)，使h(n)正好等于ha(t)的采样值，即h(n)=ha(nT) （2-3）T为采样周期。

如以Ha (s)及H(z)分别表示ha(t)的拉式变换及h(n)的z变换，即Ha(s)=L[ha(t)] （2-4）H(z)=Z[h(n)] （2-5）则根据采样序列z变换与模拟信号拉式变换的关系，得：（2-6）上式表明，采样脉冲响应不变法将模拟滤波器变换为数字滤波器时，它所完成的S平面到Z平面的变换，正是以前讨论的拉式变换到Z变换的标准变换，即首先对Ha(s)作周期严拓，然后再经过z=e st的映射关系映射到Z平面上。

应当指出，Z=e st的映射关系表明，S平面上每一条宽为2pi/T 的横带部分，都将重叠地映射到Z平面的整个全部平面上。

每一横带的左半部分映射到Z平面单位圆以内，每一横带的右半部分映射到Z平面单位圆以外，jΩ轴映射在单位圆上，但jΩ轴上的每一段2pi/T都应于绕单位圆一周，如下图2.2所示：图2.2 脉冲响应不变法的映射关系Z=e st的映射关系反映的是Ha(s)的周期严拓与H(z)的关系，而不是Ha(s)b本身与H(z)的关系，因此，使用脉冲响应不变法时，从Ha(s)到H(z)并没有一个由S平面到Z平面的简单代数映射关系，即没有一个s=f(z)的代数关系式。

另外，数字滤波器的频响也不是简单的重现模拟滤波器的频响，而是模拟滤波器频响的周期严拓，周期为ΩS =2π/T=2πfs，即（2-7）三、具体设计步骤(1). 语音信号的采集及分析在我的电脑中通过搜索，找到Windows关机.wav文件，并命名为off.wav，粘贴至C:\MATLAB6p5p1\work目录下。

画出语音信号的时域波形频域幅度谱：>> [y,Fs,bits] = wavread('off.wav');>> sound (y,Fs,bits)>> plot(y);title('时域波形')>> t=(1:16000)/8000;>>plot(t,y);xlabel('t');ylabel('y');('时域波形');图3.1时域波形>> Y=fft(y);>> magY=abs(Y);angY=angle(Y);>> w=(1:16000)/16000*2*pi;>> plot(w/pi,magY);图3.2频域幅度谱(2). 给原始信号加上一个高斯白噪声 >> [y,Fs,bit]=wavread ('off.wav'); %读入语音信号 >>sound(y,Fs,bit); %回放语音信号>>n=length(y); %求出语音信号的长度>>Y=fft(y,n); %对采样得到的语音信号进行DFT 变换>>f=(0:length(y)-1)'*Fs/length(y);%加高斯白噪声 g=awgn(y,20); %给语音信号加上高斯白噪声sound(g,Fs,bit); %回放加噪信号G=fft(g,n); %对加噪后的语音信号进行DFT 变换(3). 设计一个滤波器，滤除高频噪声将数字滤波器的设计指标设为通带截止频率fb=1100HZ,阻带频率fc=1200HZ,通带波纹Ap=1dB ，阻带波纹As=20dB,要求确定H(z)。

设计步骤如下：（1） 选取某一T 。

这是任意的，我们选T=1。

并确定模拟频率如（2） 利用设参数p Ω，s Ω，p A 和s A 设计一个模拟滤波器)(s H a 。

Tpp ω=ΩT s s ω=Ω和 式3-1（3） 利用部分分式展开，将)(s H a 展开为如式3-2：式3-2 （4）现在将模拟极点{p k }变换为数字极点{e p k T }，得到数字滤波器如式3-3：式3-3并作化简得出作为z -1有理函数的H(z)。

四、程序实现及实验结果1．用MATLAB 对原始语音信号进行分析，加入高斯白噪声后播放合成的语音信号，并进行滤波前后，时域波形和频谱分析。

∑=-=N k kk a p s R s H 1)(∑=--=N k T p k z e R z H k 111)(图4.1 滤波前后信号的波形、频谱对比图4.2滤波前后波形相位对比结果分析：由图4.1中滤波前后波形比较可看出，经过滤波后的波形比原波形的振幅有所减小，去除了很多由于噪声所产生的干扰；从滤波前后的频谱比较可以看出经过滤波后除了原本的声音外，中间由于噪声产生的频谱波形已经滤除；由图4.2滤波前后相位比较图可看出由于经过滤波，相位变得稀疏；经过MATLAB仿真，听滤波前后的声音，可以听出有明显的滤波效果。

因此利用脉冲响应不变法设计的巴特沃斯滤波器已经达到了设计的要求。

四、总结这次课程设计，给我留下了很深的印象。

虽然只是短暂的一周，但在这期间，却让我受益匪浅。

通过这次课程设计，使我对语音信号有了全面的认识，对数字信号处理的知识又有了深刻的理解，对时域频域的信号特点有了了解，了解了许多音频频谱方面的基本知识，也复习了数字滤波器的实现原理，加深了印象。

虽然在以前的课程中使用过MATLAB，但只是那些基础的操作，让我在实现这个题目时我从下手，通过DSP与MATLAB的完美结合，我们更熟悉了数字信号处理的基本知识和MATLAB的m语言，让我们把课上的理论知识运用到实际中去，让我们更近一步地巩固了课堂上所学的理论知识，并能很好地理解与掌握数字信号处理中的基本概念、基本原理、基本分析方法。

遗憾的是，课设作品中，采用的是系统自带的WAV文件，并不是自己录制的，去噪效果不如人的声音那样明显，有机会我一定会加以改进，完善这次设计。

随着对课题的深入了解，我也意识到自己对MATLAB的函数知识掌握匮乏，要加以深入学习，才能在以后的学习中更好的运用MATLAB处理数字信号，达到我们预期的目的。

五、参考文献1. 程佩青.数字信号处理:清华大学出版社,2007年2月2. 从玉良.数字信号处理原理及MATLAB实现:电子工业出版社,2005 年7月3. 刘波. MATLAB信号处理:电子工业出版社 ,2006年1月[y,Fs,bit]=wavread ('off.wav'); %读入语音信号%其中向量y为采样值，Fs为采样频率（Hz），bit为采样点数sound(y,Fs,bit); %回放语音信号n=length(y); %求出语音信号的长度Y=fft(y,n); %对采样得到的语音信号进行DFT变换f=(0:length(y)-1)'*Fs/length(y);%加高斯白噪声g=awgn(y,20); %给语音信号加上高斯白噪声sound(g,Fs,bit); %回放加噪信号G=fft(g,n); %对加噪后的语音信号进行DFT变换Fp=1200;%阻带截止频率Fs=1100;%通带截止频率Ft=8000;%采集频率As=20;%通带波纹Ap=1Ap=1;%阻带波纹As=20wp=2*pi*Fp/Ft;ws=2*pi*Fs/Ft;fp=2*Ft*tan(wp/2);fs=2*Ft*tan(ws/2);[n,wn]=buttord(wp,ws,Ap,As,'s'); %求低通滤波器的阶数和截止频率[b,a]=butter(n,wn,'s'); %求S域的频率响应的参数[num,den]=bilinear(b,a,1); %双线性变换实现S域到Z域的变换[h,w]=freqz(num,den); %根据参数求出频率响应figure(3)plot(w*8000*0.5/pi,abs(h));z=filter(num,den,y);sound(z);m=z; %求滤波后的信号figure(1)subplot(2,2,3);plot(abs(m),'r');title('滤波后信号的频谱'); grid;subplot(2,2,4);plot(z,'b');title('滤波后的信号波形');grid;subplot(2,2,2);plot(y,'b');title('滤波前信号的波形');grid;subplot(2,2,1);plot(abs(Y),'r');title('滤波前信号的频谱');grid;figure(2);p=angle(m);q=angle(Y);subplot(2,1,1);plot(q,'b');title('滤波前相位');grid;subplot(2,1,2);plot(p,'b');title('滤波后相位');grid;[y,Fs,bits] = wavread('off.wav');（本资料素材和资料部分来自网络，仅供参考。