基于MATLAB的语音信号滤波处理毕业设计(论文).doc

学毕业设计(论文)题目: 基于matlab的语音信号滤波处理作 者 夏杰 届 别 2015届 院 别 信息与通信工程学院 专 业 电子信息工程 指导教师 罗朝明 职 称 副教授 完成时间 2015年5月17日 摘要随着计算机和信息科学的飞速发展，信号处理逐渐发展成为一门独立的学科，成为信息科学的重要组成部分，在语音处理、雷达、图像处理、通信、生物医学工程等众多领域中得到广泛应应用。本论文研究了滤波器的设计，并通过设计好的滤波器语来分析语音信号。在计算机中录入一段语音信号，通过matlab对语音信号进行时频变换和分析；再通过matlab构造一定频率的噪音信号，再将两个信号进行合成；然后设计不同类型的数字滤波器对其噪音进行滤波处理。关键词：matlab；语音信号；滤波器；仿真；abstractwith the rapid development of computers and information science, signal processing is becoming an independent discipline, and playing an important part in information science. it has been widely applied in speech processing, radar, image processing, communication, bio-medical engineering and many other fields. this paper has studied the filter design, and it analyses the voice signal through the designed filter. after recording a voice signal in the computer, then doing the time-frequency transformation and analyses for the voice signal through matlab. we construct a fixed-frequency signal, and synthesis the two signals, and it deals with filter processing on the speech signal with noise, and analyses the signals in the time domain and frequency domain.keyword: matlab；voice signal；filter；simulation；目 录摘要iabstractii第一章 绪论11.1 语音信号滤波处理的简介11.2 数字滤波器的概念11.3 matlab软件的介绍31.4 论文各章节的安排4第二章 语音信号的频谱分析52.1 快速傅里叶算法的原理52.2语音信号的采集及频谱分析102.3 语音信号的加噪处理及频谱分析12第三章 数字滤波器的设计153.1数字滤波器的设计原理153.2 用窗函数法设计fir数字滤波器183.3 用双线性法设计iir数字滤波器21第四章 语音信号的滤波244.1 语音信号的滤波测试244.2 加噪语音信号的滤波处理27第五章 总结和展望295.1 课题总结295.2 课题展望29参考文献30致谢31ii第1章 绪论1.1 语音信号滤波处理的简介语音信号的基本组成单位是“音素”。音素主要由“浊音”和“清音”这两大类组成。如果将语音只有背景噪声称为“无声”，那么可以将音素分成“清音”、“浊音”、“无声”这三类。一个音节由辅音和元音两部分组成。因为元音都是浊音，所以元音在音节中占据着主要的部分。在普通话中，“辅音元音”构成了每个音节。语音信号的滤波处理是数字信号处理领域目前发展最为迅速的信息科学研究领域的核心技术之一，通过语音传递信息是人类最重要、最有效、最常用和最方便的交换信息形式。语音信号处理是现代信息处理的基本内容，数字信号的处理更是重中之重。数字信号处理的研究内容主要是语音信号和图像信号，而研究语音信号对于现代语音通信有着积极的意义。研究语音信号又分为时域和频域两个方面。（1） 语音信号的时域分析处理：对语音信号进行时域分析的本质就是分析和提取语音信号的时域参数。当进行语音分析时，它的时域波形图是很直观的。其实语音信号的实质就是时域信号，所以使用时域分析是最早也是最广泛的一种分析方法。时域分析一般用于参数的分析和应用，如语音的大分类、分割、预处理等。这种分析方法有以下特点：第一，比较直观的语音信号和明确其物理意义。第二，可以简单的实现并且运算少。第三，获得语音信号的重要参数。最后，它的实现用到了示波器等常规设备，所以应用起来较简单。（2）语音信号的频谱分析处理：语音信号的频域分析就是对语音信号的频域持征进行分析。语音信号的频域分析包括语音信号的功率谱、倒频谱、频谱、频谱包络分析等，而一般常用的频域分析方法包括带傅里叶变换法、通滤波器组法等几种。本文介绍的是语音信号的傅里叶分析法。1.2 数字滤波器的概念数字信号处理是一门发展迅速、应用广泛的前沿性学科，其理论和实践性都很强。为了系统的掌握和理解信号的收集、处理、传输、显示和储存，本课题介绍了用matlab来处理语音信号。在21世纪这个信息技术高速发展的时代，信号处理逐渐趋向于数字化、软件化这方面发展。数字信号处理主要是对语音信号和图像信号的研究，其中，语音信号的研究可以分为从时域和频域这两个方面来进行。时域的分析处理有以下两种方法：一种是对语音信号进行分析，这一分析方法属于线性处理的范畴，其原理主要是通过信号的卷积、倍乘、时移、加减、以及求相关函数等来实现；另外一种是通过处理生成和变换成各种调制信号，这一分析方法属于非线性的范畴，其原理是压缩和扩张信号平均累加器的动态范围，用门限法对噪音进行抑制。通过对信号的频域分析处理，可以对信号的频率特性在频谱中加以分析研究，这是对不确定信号分析的主要方法，主要是扩展了信号分析的范围。通常在世纪应用中，信号的时、频域分析经常同时进行。数字滤波是数字信号处理理论的一部分。具体来说，凡是用数字方式对信号进行滤波、变换、调制、压缩、增强、均衡、解调、估值、识别、产生等加工处理，都可以纳入数字信号处理的领域。 数字信号处理这门学科主要是对数字滤波器设计方法进行研究。对于数字滤波器，早在49年代末期，就有人对其可能性进行过研究讨论，在50年代的时候，也有人研究讨论过数字滤波这个具有跨时代意义的问题。但可惜的是直到60年代中期，人们才是开始慢慢形成关于数字滤波器的一套完整的理论，在这一时期，通过不同形式的数字滤波结构（有的以最小运算误差为特点，有的以最快运算速度为特点，而有的则囊括以上两者）。在数字滤波器领域中，有限冲激响应（fir）和无限冲激响应（iir）的关系认识转化是一个重大的发现。在初期，人们一般认为iir滤波器的运算效率比fir滤波器更高，所以明显青睐于前者，但当有人提出可以用快速傅里叶变换（fft）来实现其卷积运算的理念后，发现高阶fir滤波器统一可以用比较高的运算效率来实现，这就激励人们大量的研究和探索数字滤波器的fir高性能滤波器的设计方法和频域设计方法，从而有了数字滤波器设计中频域方法和时域方法的出世。然而，这些数字滤波器都是一些人们的早期研究。早期的数字滤波器尽管在声纳、语音、地震和医学和的信号处理中曾经起过非常大的作用，但是出于当时计算机主机的价格比较昂贵，在一定程度上阻碍了个别专门用途数字滤波器的发展。在科学技术蓬勃发展的70年代，引进了计算机辅助的设计方法，这使数字滤波器的设计以及实现仅仅是逼近响应模拟滤波器。通常来说，变换模拟滤波器函数来设计数字滤波器，几乎很难逼近任意冲激响应或频率响应，然而采用计算机辅助设计这一方法则可能实现时域或频域的最佳逼近，或时域频域的联合最佳逼近，这样，不但可以使数字滤波器的分析和设计内容变的丰富起来，而且各种新的数字信号处理系统也都可以用专用数字硬件加以实现。1.3 matlab软件的介绍matlab是由美国mathworks公司推出的用于图形处理和数值计算的数学计算环境。在matlab环境下，用户可以集成地进行程序设计、数值计算、输入输出、图形绘制、文件管理等各项操作。所以相比其他同类软件，它卓越的数据可视化和优秀的数字计算能力使其脱颖而出。最初的matlab系统是由cleve moler公司用fortran语言设计的，现在matlab程序是基于mathworks公司在c语言环境下开发的。1984年开始发行它的第一版（dos版本1.0），而后经过10多年的不断完善和改进，版本不断升级 ，其功能越来越丰富，工具也越来越丰富，因此它应用的范围也越来越广发。matlab语言是一种对数据分析和处理功能十分强大的据算计应用软件，它可以将声音文件变换为离散的数据文件，然后利用其强大的矩阵运算能力处理数据，如数字滤波、时域和频域分析、傅立叶变换、声音回放以及各种图的呈现等，信号处理是matlab的最重要的应用领域之一。matlab语言之所以能如此迅速地普及是由于它有着不同于其他语言的特点：1.简洁紧凑的语言使用非常方便，并且库函数非常丰富。matlab程序书写形式自由，利用其丰富的库函数巧妙的避开了子程序编程任务的复杂，压缩了一切不必要的编程工作。2.运算符丰富。由于matlab是用c语言编写的，所以matlab提供了和c语言几乎一样多的运算符，灵活使用matlab的运算符将使程序变得极为简短。3.matlab既具有结构化的控制语句，又有面向对象编程的特性。4.语法限制不严格，程序设计自由度大。5.程序的可移植性很好，基本上不做什么修改就可以运行在各种型号的计算机和操作系统上。6.matlab的图形功能强大。在c和fortran语言里，绘图都不是很容易，但是在matlab里通过数据的可视化则显得非常简单。此外，matlab编辑图形界面的能力非常强。7.matlab的一项特色是拥有功能强大的工具箱。matlab包含两个部分：核心部分和各种可选的工具箱。核心部分有数百个核心内部函数。本设计只使用了上面的部分功能，就是用matlab对含噪的语音信号同时在时域和频域进行滤波处理和分析。1.4 论文各章节的安排第一章 介绍了数字信号处理的发展趋势和所用的的开发软件matlb。第二章 介绍了语音信号的采集和频谱分析的方法和原理。第三章 介绍了数字滤波器的设计方法。第四章 介绍了用不同的数字滤波器对语音信号进行滤波处理。第五章 介绍了本人对这次毕业设计的总结和展望未来数字滤波器的发展。第2章 语音信号的频谱分析2.1 快速傅里叶算法的原理dft是信号分析与处理中的一种重要变换。因为dft直接计算的计算量与变换区间长度（n的平方）成正比，当n较大时，造成计算量太大，因而在快速傅里叶变换（简称fft）出现以前，直接用dft算法对信号进行频谱分析和它的实时处理是不切实际的。直到1965年发现了dft这种快速算法，这种情况才发生了逆转。自1965年图基（j.w.tuky）和库力(t.w.coody)在书名为计算机数学（math.computation,vol,19.1965)杂志上发表了著名的一篇机器计算傅里叶级数的一种算法论文后，桑德（g.sand）-图基等快速傅里叶算法的相继出现，在人们的不断进行改进和不断完善，便很快形成了一套高效运算方法，从而才有了我们现在的快速傅里叶变换，简称为fft（fast fourier transform）。dft的运算效率在这种算法的基础上可以提高12个数量级，这为数字信号处理技术应用于各种信号的实时处理创造了良好的条件，极大的推动了数字信号处理技术在当今社会的发展。fft算法的概念：(1)“级”概念：将n 点dft先分成两个n/2点dft，再是四个n/4点dft直至n/2个两点dft每分一次称为“一”级运算。因为n=2m所以n点dft可分成m级。(2)“组”概念：每一级都有n/2个蝶形单元，例如：n=16,则每级都有8个蝶形单元。每一级的n/2个蝶形单元可以分成若干组，每一组具有相同的结构，相同的 因子分布，第m级的组数为： n/2m+1 n=16=24，分4级。(3)因子的分布: 每由后向前（m由m-1-0级）推进一级，则此系数为后级系数中偶数序号的那一半。(4)按时间抽取法：由于每一步分解都是基于在每级按输入时间序列的次序是属于偶数还是奇数来分解为两个更短的序列，所以称为“按时间抽取法”。如图2.1所示fft算法基本上分为两大类：时域抽取法fft（decimation-in-time fft,简称dit-fft）和频域抽取法fft（decimation-in-frequency fft,简称dif-fft）。下面就介绍dif-fft算法。设序列长度为,且满足 (为正整数）按的奇偶性把分解成为两个点的子序列：2点2点2点2点2点2点2点2点两个2点dft两个2点dft两个2点dft两个2点dft两个4点dtf两个4点dft两个n/2点的dft2.1按时间抽取法 （2.1）则可将dft化为 （2.2）由于故式2.2可表示成 （2.3）与分别是及的点dft： （2.4） （2.5）由此，我们可以看到，一个点dft已分解成两个的dft。这两个点的dft再按照2.4和2.5式组合成一个点dft。 （2.6）这里应该看到,只有个点，即。而x(k)却有个点，即,故用上式计算得到的只是的前一半结果，要用,来表达全部的值，还必须应用系数的周期性，即这样可得到 （2.7）同理可得 （2.8）前两式说明了后半部分值（）所对应的，分别等于前半部分值（）所对应的，。再考虑的以下性质： 这样，把上述各式带入，就可以将表达为前后两部分： k=0,1, （2.9）因此，只要求出0到区间的所有，值，即可求出0到区间内的所有值，显然节省了运算量。可以用蝶形信号流图表示。 一个点dft可以分解为两个点dft，每一个点dft只需次复数乘法以及次复数加法。两个点dft共需次复数乘法和次复数加法。把两个点dft和成为点dft时，其中有个蝶形运算还需要次复数乘法及次复数加法。因而通过第一次分解后，总共还需要次复数乘法以及次复数加法。由因此，通过这样分解后，它的运算工作量至少节省了一半。既然这样分解是有效的，由于，因而仍是偶数，所以可以再一次把每个点子序列再其奇偶部分分解为两个点的子序列。 （2.10） （2.11）且： 式中： （2.12） （2.13）根据上面同样的分析可以得出：利用四个点的dft及两级蝶形组合运算来计算点dft，相比只用一次分解蝶形组合方式的计算量又减少了一大半。它的算法步骤如下：（1）分组，变量置换：dft变换： 先将按的奇偶分为两组，作变量置换:n=偶数时，令n=2r;得到：则其dft可化为两部分：前半部分： 后半部分：代入dft中 （2）代入dft变换式：（3）求出子序列的dft： （2.14）（4）由上面的公式我们可以知道：一个点的dft被分解为两个点dft。，这两个点的dft按照式2.15又合成点dft中的前半部分再应用系数的周期性，求出用，表达的后半部的值。 （2.15）（5）求出后半部的表示式： （2.16）看出：后半部的值所对应的，则完全重复了前半部分的值所对应的，的值。又（6）最后得出：频域中的个点频率成分为：（2.17）结论：只要求出区间内的各个整数值所对结应的，值，即可以求出在整个区间内全部值，这就是为什么fft能大量节省计算的关键所在。由于，因此仍为偶数，依照上面的方法可以进一步把每个点子序列，再按输入的奇偶分解为两个点的子序列，再按这种方法不断划分下去，直到最后剩下的是2点dft，两点dft实际上只是加减运算。上面频域中前/后半部分表示式可以用蝶形信号流图表示。如图2.2所示 2.2 蝶形信号流图2.2语音信号的采集及频谱分析采样频率在进行模拟数字信号的转换过程中，当采样频率fs.max大于信号中最好频率fmax的2倍时，即：fs.max=2fmax，则采样之后的数字信号完整的保留了原始信号中的信息，一般实际应用中保证采样频率为信号最高频率的510倍，采样定理又称为耐奎斯特定理。1924年耐奎斯特（nyquist）就推导出在理想低通信道的最高大码元传输速率的公式：理想低通信道的最高大码元传输速率=2w*log2n（其中w是理想低通信道的带宽，n是电平强度） 采样位数即采样值或取样值，用来衡量声音波动变化的参数，是指声卡在采集和播放声音文件时所使用数字声音信号的二进制位数。采样频率是指录音设备在一秒钟内对声音信号的采集次数，采样频率越高，声音的还原就越真实越自然。采样位数和采样频率对于音频接口来说事最为重要的两个指标，也是选择音频接口的两个重要标准。无论采样频率如何，理论上来说采样的位数决定了音频数据的最大力度范围。每增加一个采样位数，相当于力度范围增加了6db。采样位数越多则捕捉到的信号越精确。对于采样率来说你可以想象它类似于一个照相机，44.1khz意味着音频流进入计算机时计算机每秒会对其拍照达到441000次。显然采样频率越高，计算机摄取的图片越多，对于原始音频的还原也越加精确。 该设计采用的是优美的铃声为分析样本，利用windows下的录音机录制一小段语音，然后在matlab软件平台下，利用wavread函数对语音信号进行采样，采样频率为fs=44100hz，得到的采样频率和采样点数为：y,fs,bit=wavread(c:代码wy.wav); 采样频率：fs=44100hz bits=16 采样点数：n=1329408 对原始语音信号进行采集并播放： x1,fs,bits=wavread(f:yuyinpsb.wav); sound(x1,fs,bits); 信号频谱分析是将信号源发出的信号强度按频率顺序展开，使其成为频率的函数，并考察变化规律，称为频谱分析。频谱分析主要分析信号是由哪些频率的正弦信号叠加得到的，以及这些正弦信号的振幅。信号频谱分析是将信号源发出的信号强度按频率顺序展开，使其成为频率的函数，并考察变化规律，称为频谱分析。频谱分析主要分析信号是由哪些频率的正弦信号叠加得到的，以及这些正弦信号的振幅。首先画出语音信号的时域波形，然后用函数fft对语音信号进行傅立叶变换，得到信号的频谱特性。然后通过sound函数回放原语音信号，我们听到一段优美的英文歌。其程序如下clearclcx1,fs,bits=wavread(c:代码wy.wav);size(x1)%求y的行数和列数y1=x1(:,1);%对信号进行分列处理(单声道处理)，前面表示所有行,后面表示取单列n1=length(y1);%取y的长度，也就是总采样点数t=1/fs;%采样间隔ts=t*1000%化成mst1=(0:n1-1)*ts;%设置波形图横坐标为时间figure(1);%创建绘图窗口plot(t1,y1);%画出音乐信号时域波形图xlabel(时间（ms）);%横坐标命名ylabel(化到-1到1之间的幅值（v）);%纵坐标命名 title(原始音乐信号时域图)%标题n=length(x1) %计算音频信号的长度figure(2);ftemp=-n/2:n/2-1.*fs/n;y1=fft(x1); %做length(x1)点的fftplot(ftemp,fftshift(abs(y1); %画出原始语音信号的频谱图title(原始语音信号的频谱);grid on;sound(x1,fs,bits);如图3.1和3.2所示 图3.1 原信号时域图 图3.2 原信号频谱图2.3 语音信号的加噪处理及频谱分析在matlab中，运用randn函数产生与原信号等长度的随机噪音信号。在通过高通滤波器滤除4000hz以前的噪音信号，在将其频谱与原信号频谱叠加得到最终的噪音信号。其程序如下：n=length(x1) %计算音频信号的长度noise=0.1*randn(n,2); %产生等长度的随机噪声信号(这里的噪声的大小取决于随机函数的幅度倍数）x1_z=noise; %将两个信号叠加成一个新的信号加噪声处理 n=length(x1); %选取变换的点数 x1_zp=fft(x1_z,n); %对n点进行傅里叶变换到频域f=fs*(0:n/2-1)/n; % 对应点的频率figure(2)subplot(2,1,1);plot(x1_z); %加噪语音信号的时域波形图title(加噪语音信号时域波形);xlabel(时间轴)ylabel(幅值 a)subplot(2,1,2);plot(f,abs(x1_zp(1:n/2); %加噪语音信号的频谱图title(加噪语音信号频谱图);xlabel(频率hz);ylabel(频率幅值);%提取4000hz以后的噪音信号%wp1=2*pi*3800/fs;ws1=2*pi*4000/fs;wc1=(wp1+ws1)/2/pi;bt=ws1-wp1;n1=ceil(1.8*pi/bt);n2=n1+mod(n1+1,2);window=boxcar(n2); %长度为n1的矩形窗windowb1=fir1(n2-1,wc1,high,window);x1_high = filter(b1,1, x1_z); %对信号进行低通滤波figure(5);ftemp=-n/2:n/2-1.*fs/n;plot(ftemp,fftshift(abs(fft(x1_high);title(信号经过高通滤波器（频域）);z1=y1+fft(x1_high);figure(6);plot(ftemp,fftshift(abs(z1);title(叠加信号的频谱)grid on;如图3.3和3.4所示 图3.3 噪音信号的频谱图 图3.4 叠加信号的频谱图第3章 数字滤波器的设计3.1数字滤波器的设计原理滤波器分为模拟滤波器和数字滤波器。模拟滤波器由电阻，电容，电感等原器件构成；而数字滤波器则可由数字电路实现或软件实现。按脉冲响应来分类，数字滤波器分为iir和fir，即无限冲激响应滤波器和有限冲激响应滤波器；其中iir网络中有反馈回路，fir网络中没有反馈回路。从小的方面分：a.按所处理的信号分为两种：模拟滤波器和数字滤波器。b.按所通过信号的频段分为低通、高通、带通和带阻滤波器四种。 数字滤波器按单位取样响应时域特性可分为:iir(infinite impulse response) 滤波器(即无限长单位冲激响应滤波器)fir(finite impulse response)滤波器(即有限长单位冲激响应滤波器)。低通、高通、带通和带阻滤波器都称作是经典滤波器，它们每一种又有模拟滤波器和数字滤波器两种形式。 数字滤波器是数字信号处理的基础,可以对信号进行过滤、参数估计和检测等处理，在图像、通信、雷达、语音等许多领域都有着广泛的应用。特别是在数据压缩和图像处理等方面取得了非常大的成就。所以数字滤波器的设计显得非常重要。matlab的出现为数字滤波的研究和应用提供了一个很好的平台。它在以矩阵运算的基础上，把程序设计、计算和可视化融合到了一个交叉式的环境中。数字滤波器用来对信号进行过滤、检测与参数估计等处理,在通信、图像、语音、雷达等许多领域都有着十分广泛的应用。尤其在图像处理、数据压缩等方面取得非常大的成就。这使各个领域的研究人员可以直观方便地进行科学研究与工程应用。当然设计一数字滤波器首先要清楚其原理。fir线性相位滤波器的原理如下:设是长度为的线性相位滤波器的脉冲响应，它的系统函数为：它的频率响应函数为：它具有一个线性相位约束：其中，为一个常数。则有以下面的特性：称之为对称脉冲响应；或者称之为反对称脉冲响应。 于是,根据值的不同以及m的奇偶性,就分别产生了四种类型的线性fir相位滤波器.它们是:1型线性相位fir滤波器:为奇数, 以中心点对称;2型线性相位fir滤波器:为偶数, 以点与点之间的中心对称;3型线性相位fir滤波器:为奇数, 以点为中心反对称;4型线性相位fir滤波器:为偶数, 以点与点之间的中心反对称;根据这四类fir滤波器,可得到其响应的频率响应的特征:其中为振幅响应,与之区别的是它与幅值特性不同,前者的取值可正可负,而后者只能为非负值。iir数字滤波器具有以下特点：(1) iir数字滤波器的系统函数可以写成封闭函数的形式。 (2) iir数字滤波器采用递归型结构结构上带有反馈环路。iir滤波器运算结构一般通常由相加、延时和乘以系数等基本运算组成，即可以构成级联型、直接型、并联型、正准型四种结构形式，它们都是具有反馈回路这一结构。由于在运算中会出现舍入处理，会造成误差不断累积，偶尔会产生微弱的寄生振荡。 (3) iir数字滤波器在设计上可以借鉴成熟的模拟滤波器成果，如椭圆滤波器、巴特沃斯和契比雪夫等，所以在设计模拟器的时候可以查看它们的设计数据以及图表，这样我们的设计工作量相对而言就少了许多，其次，对计算工具的要求也不高。通过根据指标写出模拟滤波器的公式以及通过一定变换，可以将模拟滤波器转换为数字滤波器的公式。 (4) iir数字滤波器是非线性的，所以其相位特性不好控制，当对相位要求较高时，则需要添加相位校准网络。 iir数字滤波器可用一个n阶差分方程表示： （2.2） 或用它的z域系统函数: （2.3） 对照模拟滤波器的传递函数: （2.4） 可以看出,设计数字滤波器的思路与模拟滤波器相仿,设计实质都是寻找一组系数b,a,去逼近所求的频率响应,以便在性能上满足预定的技术要求；不同的是数字滤波器是在z平面寻找合适的h(z), 模拟滤波器则是在s平面用数学逼近法寻找近似的所需特性h(s)。由于模拟滤波器的设计在理论上已十分成熟,因此数字滤波器设计的关键是将h(s)变换至h(z),即利用复值映射将模拟滤波器离散化。3.2 用窗函数法设计fir数字滤波器 窗函数设计法是一种通过截短和计权的方法使无限长非因果序列成为有限长脉冲应响应序列的设计方法。一般在设计滤波器之前，滤波器的技术指标是要根据具体的工程应用确定的。通常在运用中，因为数字滤波器往往被用来实现选频操作，所以指标形式通常为在频域中以分贝值给出的相应响应和幅度响应。窗函数法设计fir滤波器的步骤如下：(1) 根据阻带衰减及过度带宽，选择窗函数的类型并估计其窗口长度n（或阶数m=n-1）。窗函数的类型可以根据最小阻带衰减as独立选择，原因是因为窗口长度n的大小对最小阻带衰减as几乎没有影响。一旦窗函数类型确定以后，可以根据过度带宽小于给定指标确定所拟用的窗函数的窗口长度n，假设待设计的滤波器的过渡带宽为，它与窗口长度n的大小近似成反比，所以只要窗函数类型确定后，那么它的计算公式也就确定了。但是这些公式都是近似的，得出的窗口长度n还需要在计算中不断修正以达到最好的要求。原则是：在保证阻带衰减满足其规定要求的情况下，选择比较小的n，在n和窗函数类型确定后，就可以调用matlab软件中的窗函数求出窗函数w（n）。 (2) 根据待求滤波器的理想频率响应求出理想单位脉冲响应hd（n），如果给出基于窗函数法的fir数字低通滤波器设计。在一般情况下，hd（n）是不能用封闭公式表示的，需要采用数值方法表示；从w=0到w=2采样n点，采用离散傅里叶反变换（idft）即可求出。 (3) 计算滤波器的单位脉冲响应h（n），它是理想单位脉冲响应和窗函数的乘积。 (4) 计算技术指标是否满足要求，为了计算数字滤波器在频域中的特性，可调用matlab中的freqz函数，如果不满足要求的情况下，可根据具体情况适当调整窗函数长度或类型，直到满足要求为止。 通常上述的几点难以同时满足。实际中设计fir数字滤波器往往要求是线性相位的，因此要求)(nw满足线性相位的条件，即要求w(n)满足：w(n)=w(n-1-n) （2） 所以，窗函数不仅有截短的作用，而且能够起到平滑的作用，在很多领域得到了应用。以下为利用窗函数法设计低通、高通和带通滤波器。低通滤波器设计程序：x1,fs,bits=wavread(c:代码wy.wav);wp1=2*pi*3800/fs;ws1=2*pi*4000/fs;wc1=(wp1+ws1)/2/pi;bt=ws1-wp1;n1=ceil(1.8*pi/bt);window=boxcar(n1+1); %长度为n1的矩形窗windowb1=fir1(n1,wc1,window);figure(1);freqz(b1,1,512);title(低通滤波器的频率响应);如3.5图所示图3.5 低通滤波器频率响应高通滤波器设计程序：x1,fs,bits=wavread(c:代码wy.wav);wp1=2*pi*3800/fs;ws1=2*pi*4000/fs;wc1=(wp1+ws1)/2/pi;bt=ws1-wp1;n1=ceil(1.8*pi/bt);n2=n1+mod(n1+1,2);window=boxcar(n2); %长度为n1的矩形窗windowb1=fir1(n2-1,wc1,high,window);figure(3);freqz(b1,1,512);title(高通滤波器的频率响应)如图3.6所示图3.6 高通滤波器频率响应带通滤波器程序：clearclcfs=44100;wp1=2*pi*3205/fs;wp2=2*pi*4410/fs; %带通滤波器的通带范围ws1=2*pi*3000/fs;ws2=2*pi*4600/fs;bt=wp1-ws1;wc=(wp1+ws1)/2/pi,(wp2+ws2)/2/pi;n3=ceil(1.8*pi/bt);window=boxcar(n3); %长度为n1的矩形窗windowb1=fir1(n3-1,wc,window);figure(9);freqz(b1,1,512); %数字滤波器频率响应title(带通滤波器的频率响应);如3.7图所示图3.7 带通滤波器频率响应3.3 用双线性法设计iir数字滤波器采用非线性频率压缩方法，将整个频率轴上的频率范围压缩到-/t/t之间，再用z=est转换到z平面上。也就是说，第一步先将整个s平面压缩映射到s1平面的-/t/t一条横带里；第二步再通过标准变换关系z=es1t将横带转换到z平面上去。这样就使z平面与s平面对应的建立了单值关系，也就是说，消除了多值变换性就代表着消除了频谱混叠现象，映射关系如图3.8所示。图3.8 双线性变换的映射关系为了将s平面的整个虚轴j压缩到s1平面j1轴上的-/t到/t段上，可以通过以下的正切变换实现 (3.6)式中,t仍是采样间隔。当1由-/t经过0变化到/t时，由-经过0变化到+，也即映射了整个j轴。将式（2-6）写成 (3.7)将此关系解析延拓到整个s平面和s1平面，令j=s，j1=s1，则得 (3.8)再将s1平面通过以下标准变换关系映射到z平面 (3.9)从而得到s平面和z平面的单值映射关系为： (3.10) (3.11)式（3.10）与式（3.10）是s平面与z平面之间的单值映射关系，这种变换都是两个线性函数之比，因此称为双线性变换式（3.6）与式（3.10的双线性变换符合映射变换应满足的两点要求。首先,把z=ej，可得 （3.12）即s平面的虚轴映射到z平面的单位圆。其次，将s=+j代入式（2-12），得 （3.13）因此 （3.14）由此看出，当0时，|z|0时，|z|1。也就是说，s平面的左半平面映射到z平面的单位圆内，s平面的右半平面映射到z平面的单位圆外，s平面的虚轴映射到z平面的单位圆上。因此，稳定的模拟滤波器经双线性变换后所得的数字滤波器也一定是稳定的。本文只对巴特沃斯滤波器进行设计。用双线性变换法设计巴特沃斯低通滤波器的程序如下:x1,fs,bits=wavread(c:代码wy.wav);wp=2*3800/fs;ws=2*4000/fs;rp=3;as=10;n,wc=buttord(wp,ws,rp,as,s);b,a=butter(n,wc,s);bz,az=bilinear(b,a,fs);nd,wdc=buttord(wp,ws,rp,as);b,a=butter(nd,wdc);freqz(b,a,1000,fs);图3.9 巴特沃斯低通滤波器的频率响应 巴特沃斯低通滤波器零极点图第4章 语音信号的滤波4.1 语音信号的滤波测试本节是通过将采集的原始语音信号分别通过已设计好的低通、高通和带通滤波器，并观察和分析它们各自的频谱图，并利用函数sound函数播放滤波后的语音号，仔细聆听其中声音的不同。将信号通过4000hz的低通滤波器后，因为只剩下低频段的信号了，我们可以明显发现其声音与原始信号的声音相比声音较低沉。其信号经过低通滤波器程序：x1,fs,bits=wavread(c:代码wy.wav);n=length(x1);wp1=2*pi*3800/fs;ws1=2*pi*4000/fs;wc1=(wp1+ws1)/2/pi;bt=ws1-wp1;n1=ceil(1.8*pi/bt);window=boxcar(n1+1); %长度为n1的矩形窗windowb1=fir1(n1,wc1,window);z1 = filter(b1,1, x1); %对信号进行低通滤波figure(9);ftemp=-n/2:n/2-1.*fs/n;plot(ftemp,fftshift(abs(fft(z1);title(信号经过低通滤波器（频域）);sound(z1,fs);wavwrite(z1,fs,c:代码语音播放低通.wav);%将低通滤波后的音乐信号保存在语音播放中。如4.1图所示。将信号通过4000hz的高通滤波器后，因为只剩下高频段的信号了，我们可以明显发现其声音与原始信号的声音相比声音比较刺耳，让人听着感觉难受。其信号进过高通滤波器程序：x1,fs,bits=wavread(c:代码wy.wav);n=length(x1);图4.1 经过低通滤波器的信号频谱wp1=2*pi*3800/fs;ws1=2*pi*4000/fs;wc1=(wp1+ws1)/2/pi;bt=ws1-wp1;n1=ceil(1.8*pi/bt);n2=n1+mod(n1+1,2);window=boxcar(n2); %长度为n1的矩形窗windowb1=fir1(n2-1,wc1,high,window);z1 = filter(b1,1, x1); %对信号进行高通滤波ftemp=-n/2:n/2-1.*fs/n;plot(ftemp,fftshift(abs(fft(z1);title(信号经过高通滤波器（频域）);sound(z1,fs);wavwrite(z1,fs,c:代码语音播放高通.wav);%如图4.2所示 图4.2 经过高通滤波器的信号频谱将信号通过下面的带通滤波器后，其播放的声音让人听不出有音乐的感觉了。其信号进过带通滤波器程序：x1,fs,bits=wavread(c:代码wy.wav);n=length(x1);wp1=2*pi*3205/fs;wp2=2*pi*4410/fs; %带通滤波器的通带范围ws1=2*pi*3000/fs;ws2=2*pi*4600/fs;bt=wp1-ws1;wc=(wp1+ws1)/2/pi,(wp2+ws2)/2/pi;n3=ceil(1.8*pi/bt);window=boxcar(n3); %长度为n1的矩形窗windowb1=fir1(n3-1,wc,window);z1 = filter(b1,1, x1); %对信号进行高通滤波figure(9);ftemp=-n/2:n/2-1.*fs/n;plot(ftemp,fftshift(abs(fft(z1);title(信号经过带通滤波器（频域）);sound(z1,fs);wavwrite(z1,fs,c:代码语音播放带通.wav);%将带通滤波后的音乐信号保存在语音播放中如图4.3所示 图4.3 经过带通滤波器的信号频谱将信号通过4000hz的巴特沃斯低通滤波器后，因为只剩下低频段的信号了，我们可以明显发现其声音与原始信号的声音相比声音较低沉。其信号经过低通滤波器程序：clcclearx1,fs,bits=wavread(c:代码wy.wav);n=length(x1) %计算音频信号的长度wp=2*3800/fs;ws=2*4000/fs;rp=3;as=10;n,wc=buttord(wp,ws,rp,as,s);b,a=butter(n,wc,s);bz,az=bilinear(b,a,fs);nd,wdc=buttord(wp,ws,rp,as);b,a=butter(nd,wdc);freqz(b,a,1000,fs);z2= filter(b,a,x1); ftemp=-n/2:n/2-1.*fs/n;figure(1);plot(ftemp,fftshift(abs(fft(z2);title(信号经过低通滤波器（频域）);grid on;如图4.4所示：图4.44.2 加噪语音信号的滤波处理在以上的程序中，我们已对原始信号进行了加噪处理(其加杂的是4000hz以后的高频成分)，所以我们现在要想把语音信号进行还原，则必须要通过4000hz的低通滤波器对其进行滤波处理，然后我们就可以得到与原信号无限逼近的语音信号，通过以下程序，我们可以看到其频谱图与原信号频谱图几乎一样，声音频率也很相似。这样，我们的加噪语音信号的滤波处理就实现了。%4000hz的低通滤波器对噪音信号进行滤波%wp1=2*pi*3800/fs;ws1=2*pi*4000/fs;wc1=(wp1+ws1)/2/pi;bt=ws1-wp1;n1=ceil(1.8*pi/bt);window=boxcar(n1+1); %长度为n1的矩形窗windowb1=fir1(n1,wc1,window);z2 = filter(b1,1, ifft(z1); %对信号进行低通滤波figure(9);ftemp=-n/2:n/2-1.*fs/n;plot(ftemp,fftshift(abs(fft(z2);title(信号经过低通滤波器（频域）);grid on;so