语音信号滤波去噪――使用双线性变换法设计的巴特沃斯滤波器

月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第１页共17页拉弦乐器（二胡）语音信号滤波去噪——使用双线性变换法设计的巴特沃斯滤波器学生姓名：巢攀指导老师：胡双红摘要本课程设计是采用双线性变换法设计的巴特沃斯滤波器对拉弦音乐滤波去噪。在网上下载一段拉弦音乐，在MATLAB集成环境下，首先用wavread函数求出拉弦音乐的相关参数，对拉弦音乐进行读取和加噪；然后再给定相应技术指标，设计一个满足指标的巴普沃斯滤波器，对该拉弦音乐进行滤波去噪处理，并绘制对比图，比较滤波前后的波形和频谱并进行分析；最后通过回放拉弦音乐，对比滤波前后的信号变换。本课程设计成功的对拉弦音乐进行滤波去噪，初步完成了设计指标。关键词MATLAB；滤波去噪；IIR滤波器；巴特沃斯；双线性变换法1引言本课程设计采用双线性变换法设计的巴特沃斯滤波器对拉弦乐器语音信号滤波去噪。通过课程设计了解IIR滤波器的原理及使用方法，了解使用MATLAB语言实际IIR滤波器的设计方法。通过观察滤波前后的波形图及回放滤波前后的时域和频域图形，加深对IIR滤波器作用的理解。通过对比滤波前后的波形图，可以看出滤波器对有用信号的无失真放大有重要作用。1.1课程设计目的《数字信号处理》课程设计是在学生完成数字信号处理和MATLAB的结合后的基本实验以后开设的。本课程设计的目的是为了让学生综合数字信号处理和MATLAB并实现一个较为完整的小型滤波系统。这一点与验证性的基本实验有本质性的区别。开设课程设计环节的主要目的是通过系统设计、软件仿真、程序安排与调试、写实习报告等步骤，使月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第２页共17页学生初步掌握工程设计的具体步骤和方法，提高分析问题和解决问题的能力，提高实际应用水平。在本次课程设计的目的是如何设计一个双线性变换法的巴普沃斯滤波器滤波器，从而达到对拉弦乐器语音信号滤波去噪的效果。1.2课程设计的要求（1）滤波器指标必须符合工程实际。（2）设计完后应检查其频率响应曲线是否满足指标。（3）处理结果和分析结论应该一致，而且应符合理论。（4）独立完成课程设计并按要求编写课程设计报告书。（5）学会MATLAB的使用，掌握MATLAB的程序设计方法；1.3设计平台MATLAB是矩阵实验室（MatrixLaboratory）的简称，是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分[1]。MATLAB是由美国mathwork公司发布的主要面对科学计算、可视化以及交互式程序设计的高科技计算环境。它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中，为科学研究、工程设计以及必须进行有效数值计算的众多科学领域提供了一种全面的解决方案，并在很大程度上摆脱了传统非交互式程序设计语言（如C、Fortran）的编辑模式，代表了当今国际科学计算软件的先进水平。MATLAB和Mathematica、Maple并称为三大数学软件。它在数学类科技应用软件中在数值计算方面首屈一指。MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连接其他编程语言的程序等，主要应用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域[2]。月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第３页共17页2设计原理2.1IIR滤波器从离散时间来看，若系统的单位抽样(冲激)响应延伸到无穷长，称之为“无限长单位冲激响应系统”，简称为IIR系统。无限长单位冲激响应（IIR）滤波器有以下几个特点:(1)系统的单位冲激响应h(n)是无限长；(2)系统函数H(z)在有限z平面（0z∞）；(3)结构上存在着输出到输入的反馈，也就是结构上是递归型的。IIR滤波器采用递归型结构，即结构上带有反馈环路。同一种系统函数H（z）可以有多种不同的结构，基本网络结构有直接Ⅰ型、直接Ⅱ型、级联型、并联型四种，都具有反馈回路。同时，IIR数字滤波器在设计上可以借助成熟的模拟滤波器的成果，巴特沃斯(Butterworth)滤波器、切比雪夫(Chebyshev)滤波器、椭圆(Cauer)滤波器、贝塞尔(Bessel)滤波器等，这些典型的滤波器各有特点。有现成的设计数据或图表可查，在设计一个IIR数字滤波器时，我们根据指标先写出模拟滤波器的公式，然后通过一定的变换，将模拟滤波器的公式转换成数字滤波器的公式[3]。2.2巴特沃斯滤波器MATLAB信号处理工具箱函数buttpbuttorbutter是巴特沃斯滤波器设计函数，其有5种调用格式，本课程设计中用到的是[N,wc]=butter(N,wc,Rp,As,’s’),该格式用于计算巴特沃斯模拟滤波器的阶数N和3dB截止频率wc。函数butter,cheby1和ellip设计IIR滤波器时都是默认的双线性变换法，所以在设计滤波器时只需要代入相应的实现函数即可。下面我们设计的IIR数字滤波器的主要程序：wp=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');月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第４页共17页[b,a]=butter(n,wn,'s');[num,den]=bilinear(b,a,1);[h,w]=freqz(num,den);2.3双线性双线性变换法是使数字信号滤波器的频率响应与模拟滤波器的频率响应相似的一种变换方法。为了客服多值映射这一缺点，我们首先把整个s平面压缩变换到某一中介的s1平面的一条横带里（宽度为T2,即从T到T),其次再通过上面讨论过的标准变换关系Tsez1将此横带变换到这个z平面上去，这样就使s平面与z平面式一一对应的关系，消除了多值变换性，也就消除了频谱混叠现象。将s平面整个j轴压缩变换到s1平面j轴上的T到T一段，可以采用以下变换关系：)2tan(1T（2-1）这样，变换到T1，0变到01可将（1）式写成22221111TjTjTjTjeeeej（2-2）解析延拓到整个s平面和s1平面，令sj，11sj，则得TsTsTjTjTjTjeeTstheeeej11111111212222（2-3）再将s1平面通过以下标准变换关系映射到z平面：Tsez1（2-4）从而得到s平面和z平面的单值映射的关系为1111zzs（2-5）ssz11（2-6）一般来说，为了使模拟滤波器的某一频率与数字滤波器的任一频率有对应的关系，可月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第５页共17页以引入待定常数c，使（1）式和（2）式变换成)2tan(1Tc（2-7）TsTseecTscths111121（2-8）仍将Tsez1代入（8）式，可得1111zzcs（2-9）scscz（2-10）（9）式和（10）式是s平面与z平面之间的单值映射关系，这种变换就称为双线性变换[4]。3设计步骤3.1设计流程图拉弦乐器单频语音频信号滤波去噪——使用双线性变换法设计的巴普沃斯滤波器的设计流程图如下图3.1所示：开始拉弦乐器单频语音信号的的采集在MATLAB集成环境下用wavread函数对语音信号采样，快速傅里叶变换，并对语音信号加入噪声，画出语音信号加噪声前后的时域和频域波形图并进行比较。月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第６页共17页图3-1脉冲响应不变法巴特沃斯滤波器对语音信号去噪流程图3.2弹拨乐器（吉他）单音频语音信号的采集在网上下载一段格式为.wav的吉他单音频语音信号，时间在3s左右。然后保存在目录d:\下，命名为xyz.wav，点击windows系统桌面的“开始”按钮，点击开始菜单栏里的“附件”，选择“录音机”选项，点击录音机“文件”选项，进入“声音选定”设置，把属性一栏设置成“8000Hz，8位，单声道，7KB/秒”(见图3.2)。点击确定。设定滤波器性能指标，通带截止频率fb=1900，阻带截止频率fc=2050，通带波纹Ap=3，阻带波纹As=20脉冲响应不变法法设计巴特沃斯滤波器验证并进行频谱分析设计好的滤波器进行滤波处理比较滤波前后语音信号的波形及频谱回放语音信号结束月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第７页共17页图3-2图3-3采集声音3.3原语音信号的频谱分析在MATLAB中编辑m函数，使用wavread函数读取采集的声音文件（.wav)将它赋值给某一向量，再对其进行采样,然后使用plot语句画出相关的频谱图形。（1）Wavread函数调用格式：[y,fs,nbits]=wavread(file)功能说明：采样值放在向量x中，fs表示采样频率（Hz），bits表示采样位数。（2）快速傅里叶变换算法FFT计算DFT的函数fft,其调用格式如下：Xk=fft(x,n)参数x为被变换的时域序列向量，N是DFT变换区间长度，当n大于x的长度时，fft函数自动在x后面补零。，当n小于xn的长度时，fft函数计算x的前n个元素，忽略月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第８页共17页其后面的元素[5]。在本次课程设计中，我们利用fft函数对拉弦音乐进行快速傅里叶变换，就可以得到信号的频谱特性。（3）声音采样文件读取的程序（文件名：xyz.wav)；拉弦音乐的提取：[x,fs,bits]=wavread('d:\xyz.wav');%输入参数为文件的全路径和文件名，输出的第一个参数是每个样本的值，fs是生成该波形文件时的采样率，bits是波形文件每样本的编码位数。sound(x,fs,bits);%按指定的采样率和每样本编码位数回放N=length(x);%计算信号x的长度fn=2200;%单频噪声频率，此参数可改t=0:1/fs:(N-1)/fs;%计算时间范围，样本数除以采样频率x=x';y=x+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;%计算频谱频率范围得到的原始拉弦音乐和加上噪音后的拉弦音乐的时域波形和频谱图如图3.4所示。得到如下图形月《拉弦音乐滤波去噪--使用双线性变换法设计的巴特沃斯滤波器》第９页共17页图3-4加噪前后信号的时域和频域波形3.4滤波器设计设计指标：通带截止频率为1900Hz，阻带截止频率为2050Hz，通带波纹为3dB，阻带波纹为20dB，用双线性变换法设计的一个满足上述指标的巴特沃斯滤波器。双线性变换法设计的巴特沃斯滤波器fp=fn-200;fc=fn-50;%定义通带和阻带截止频率Rp=3;As=20;%定义通带波纹和阻带衰减wp=fp/fs*2*pi;ws=fc/fs*2*pi;%计算对应的数字频率T=2;Fs=1/T;%定义采样间隔Omegap=2/T*tan(wp/2);Omegas=2/T*tan(ws/2);%截止频率预畸变[N,wn]=buttord(Omegap,Omegas,Rp,As,'s')%计算滤波器阶数和截止频率[cs,ds]=butter(N,wn,'s')%计算滤波器系统函数分子分母系数[b,a]=bilinear(cs,ds,Fs)%双线性变换得到数字滤波器系统函数分子分母