1、1基于 matlab的语音信号分析与处理摘要:滤波器设计在数字信号处理中占有极其重要的地位,FIR 数字滤波器和IIR 滤波器是滤波器设计的重要组成部分。Matlab 功能强大、编程效率高, 特别是 Matlab 具有信号分析工具箱,不需具备很强的编程能力,就可以很方便地进行信号分析、处理和设计。基于 MATLAB 有噪音语音信号处理的设计与实现,综合运用数字信号处理的理论知识对加噪声语音信号进行时域、频域分析和滤波。使用窗函数法来设计 FIR 数字滤波器,用巴特沃斯、切比雪夫和双线性变法设计 IIR 数字滤波器,并利用 MATLAB 作为辅助工具完成设计中的计算与图形的绘制。关键词:数字滤波
2、器;MATLAB;切比雪夫Abstract: Filter design in digital signal processing plays an extremely important role, FIR digital filters and IIR filter is an important part of filter design. Matlab is powerful, programming efficiency, Matlab also has a particular signal analysis toolbox, it need not have strong pro
3、gramming skills can be easily signal analysis, processing and design. MATLAB based on the noise issue speech signal processing design and implementation of digital signal processing integrated use of the theoretical knowledge of the speech signal plus noise, time domain, frequency domain analysis an
4、d filtering. The corresponding results obtained through theoretical derivation, and then use MATLAB as a programming tool for computer implementation.Implemented in the design process, using the window function method to design FIR digital filters with Butterworth, Chebyshev and bilinear Reform IIR
5、digital filter design and use of MATLAB as a supplementary tool to complete the calculation and graphic design Drawing.Keywords:digital filter; MATLAB; Chebyshev2语音信号处理是研究用数字信号处理技术和语音学知识对语音信号进行处理的新兴的学科,是目前发展最为迅速的信息科学研究领域的核心技术之一。通过语音传递信息是人类最重要、最有效、最常用和最方便的交换信息形式。Matlab 语言是一种数据分析和处理功能十分强大的计算机应用软件,它可以将
6、声音文件变换为离散的数据文件,然后利用其强大的矩阵运算能力处理数据,如数字滤波、傅里叶变换、时域和频域分析、声音回放以及各种图的呈现等,它的信号处理与分析工具箱为语音信号分析提供了十分丰富的功能函数,利用这些功能函数可以快捷而又方便地完成语音信号的处理和分析以及信号的可视化,使人机交互更加便捷。1.语音信号处理的概念及现状语音是语言的声学表现,是人类交流信息最自然、最有效、最方便的手段。随着社会文化的进步和科学技术的发展,人类开始进入了信息化时代,用现代手段研究语音处理技术,使人们能更加有效地产生、传输、存储、和获取语音信息,这对于促进社会的发展具有十分重要的意义,因此,语音信号处理正越来越受
7、到人们的关注和广泛的研究。1.1语音信号处理的概念语音是人类获取信息的重要来源和利用信息的重要手段。通过语言相互传递信息是人类最重要的基本功能之一。语音是语言的声学表现,是相互传递信息的最重要的手段,是人类最重要、最有效、最常用和最方便的交换信息的形式。1.2语音信号处理的现状20 世纪 60 年代中期形成的一系列数字信号处理的理论和算法,如数字滤波器、快速傅立叶变换(FFT)等是语音信号数字处理的理论和技术基础。随着信息科学技术的飞速发展,语音信号处理取得了重大的进展:进入 70 年代之后,提出了用于语音信号的信息压缩和特征提取的线性预测技术(LPC) ,并已成为3语音信号处理最强有力的工具
8、,广泛应用于语音信号的分析、合成及各个应用领域,以及用于输入语音与参考样本之间时间匹配的动态规划方法;80 年代初一种新的基于聚类分析的高效数据压缩技术矢量量化(VQ)应用于语音信号处理中;而用隐马尔可夫模型(HMM)描述语音信号过程的产生是 80 年代语音信号处理技术的重大发展,目前 HMM 已构成了现代语音识别研究的重要基石。近年来人工神经网络(ANN)的研究取得了迅速发展,语音信号处理的各项课题是促进其发展的重要动力之一,同时,它的许多成果也体现在有关语音信号处理的各项技术之中。2.语音信号处理的内容和方法2.1语音信号处理的内容用Matlab对含噪的语音信号同时在时域和频域进行滤波处理
9、和分析,在Matlab应用软件下设计一个简单易用的图形用户界面(GUI) ,来解决一般应用条件下的各种语音信号的处理。主要是通过用带有录音功能的手机或计算机录取一段语音信息,把已录取的语音信息存储为.wav 格式文件,用 matlab 读取语音文件,运用数字信号学基本原理实现语音信号的处理,在 matlab 环境下综合运用信号提取,幅频变换以及傅里叶变换、滤波等技术来进行语音信号处理,能对语音信号进行采集,并对其进行各种处理,达到简单的语音信号处理的目的。2.2语音信号处理的方法在图形用户界面(Graphical User Interface,简称 GUI,又称图形用户接口)是指采用图形方式显
10、示的计算机操作用户界面。与早期计算机使用的命令行界面相比,图形界面对于用户来说在视觉上更易于接受。G UI 的 广 泛 应 用是 当 今 计 算 机 发 展 的 重 大 成 就 之 一 , 他 极 大 地 方 便 了 非 专 业 用 户 的 使 用 人 们从 此 不 再 需 要 死 记 硬 背 大 量 的 命 令 , 取 而 代 之 的 是 可 以 通 过 窗 口 、 菜 单 、 按键 等 方 式 来 方 便 地 进 行 操 作 。4基于 MATLAB 有噪音语音信号处理的设计与实现,综合运用数字信号处理的理论知识对加噪声语音信号进行时域、频域分析和滤波。通过理论推导得出相应结论,再利用 MA
11、TLAB 作为编程工具进行计算机实现。在设计实现的过程中,使用窗函数法来设计 FIR 数字滤波器,用巴特沃斯、切比雪夫和双线性变法设计 IIR 数字滤波器,并利用 MATLAB 作为辅助工具完成设计中的计算与图形的绘制。通过对对所设计滤波器的仿真和频率特性分析,可知利用 MATLAB 信号处理工具箱可以有效快捷地设计 FIR 和 IIR 数字滤波器。3.语音信号处理3.1原始语音信号采集与处理使用带有录音功能的手机或电脑的声卡设备采集一段语音信号,并将其保存在电脑中,语音信息文件为 *.wav 格式。语音信号的处理主要包括信号的提取、信号的调整、信号的变换和滤波等。通过用户图形界面的输出功能,
12、将处理后的信号的语音进行播放,试听处理后的效果。语音信号采集过程如图 3-1所示。麦克风 声卡 滤波 采样 A / D 转换W a v声 音 自 带 录 音 机indowsW图 3-1 语音信号采集过程(1)语音信号的时域分析语音信号是一种非平稳的时变信号,它携带着各种信息。在语音编码、语音合成、语音识别和语音增强等语音处理中无一例外需要提取语音中包含的各种信息。语音信号分析的目的就在与方便有效的提取并表示语音信号所携带的信息。语音信号分析可以分为时域和变换域等处理方法,其中时域分析是最简单的方法,直接对语音信号的时域波形进行分析,提取的特征参数主要有语音的短时能量,短时平均过零率,短时自相关
13、函数等。 提取:通过图形用户界面上的菜单功能按键采集电脑设备上的一段音频信号,完成音频信号的频率,幅度等信息的提取,并得到该语音信号的波形5图。 调整:在设计的用户图形界面下对输入的音频信号进行各种变化,如变化幅度、改变频率等操作,以实现对语音信号的调整。(2)语音信号的频域分析信号的傅立叶表示在信号的分析与处理中起着重要的作用。因为对于线性系统来说,可以很方便地确定其对正弦或复指数和的响应,所以傅立叶分析方法能完善地解决许多信号分析和处理问题。另外,傅立叶表示使信号的某些特性变得更明显,因此,它能更深入地说明信号的各项红物理现象。由于语音信号是随着时间变化的,通常认为,语音是一个受准周期脉冲
14、或随机噪声源激励的线性系统的输出。输出频谱是声道系统频率响应与激励源频谱的乘积。声道系统的频率响应及激励源都是随时间变化的,因此一般标准的傅立叶表示虽然适用于周期及平稳随机信号的表示,但不能直接用于语音信号。由于语音信号可以认为在短时间内,近似不变,因而可以采用短时分析法。 变换:在用户图形界面下对采集的语音信号进行 Fourier 等变换,并画出变换前后的频谱图和变换后的倒谱图。 滤波:滤除语音信号中的噪音部分,可采用低通滤波、高通滤波、带通滤波和帯阻滤波,并比较各种滤波后的效果。(3)语音信号处理流程图语音信号处理的过程包括语音信号的采集、信息提取、信号调整、信号变换、信号滤波。其中信号调
15、整又包括幅度和频率的任意倍数变化,语音信号处理流程图如图 3-2 所示。语音信号采集信息提取 信号调整 信号变换 信号滤波效果显示图 3-2 语音信号处理流程图6信号的滤波采用了四种滤波方式,来观察各种滤波性能的优缺点:如图 3-3 所示:信号滤波切比雪夫 I 型低通滤波切比雪夫 型高通滤波切比雪夫 型带阻滤波椭圆数字带通滤波图 3-3 语音信号滤波的方式在以上两图中,可以看到整个语音信号处理系统的流程大概分为三步,首先要读入待处理的语音信号,然后进行语音信号的处理,包括信息的提取、幅度和频率的变换以及语音信号的傅里叶变换、滤波等;滤波又包括低通滤波、高通滤波、带通滤波和带阻滤波等方式。最后对
16、处理过的语音信号进行处理后的效果显示。3.2 语音的录入与打开在 MATLAB 中,y,fs,bits=wavread(Blip,N1 N2);用于读取语音,采样值放在向量 y 中,fs 表示采样频率(Hz),bits 表示采样位数。N1 N2表示读取从 N1 点到 N2 点的值(若只有一个 N 的点则表示读取前 N 点的采样值)。sound(x,fs,bits); 用于对声音的回放。向量 y 则就代表了一个信号(也即一个复杂的“函数表达式”)也就是说可以像处理一个信号表达式一样处理这个声音信号。3.3时域信号的 FFT分析FFT 即 为 快 速 傅 氏 变 换 , 是 离 散 傅 氏 变 换
17、 的 快 速 算 法 , 它 是 根 据 离 散 傅氏 变 换 的 奇 、 偶 、 虚 、 实 等 特 性 , 对 离 散 傅 立 叶 变 换 的 算 法 进 行 改 进 获 得 的 。在 MATLAB 的信号处理工具箱中函数 FFT 和 IFFT 用于快速傅立叶变换和逆变换。函数 FFT 用于序列快速傅立叶变换,其调用格式为 y=fft(x),其中,x 是序列,y 是序列的 FFT,x 可以为一向量或矩阵,若 x 为一向量,y 是 x 的 FFT 且和 x相同长度;若 x 为一矩阵,则 y 是对矩阵的每一列向量进行 FFT。如果 x 长度是 2 的幂次方,函数 fft 执行高速基2FFT 算
18、法,否则 fft 执行一种混合基的7离散傅立叶变换算法,计算速度较慢。函数 FFT 的另一种调用格式为y=fft(x,N),式中,x,y 意义同前,N 为正整数。函数执行 N 点的 FFT,若 x为向量且长度小于 N,则函数将 x 补零至长度 N;若向量 x 的长度大于 N,则函数截短 x 使之长度为 N;若 x 为矩阵,按相同方法对 x 进行处理。3.4数字滤波器设计原理数字滤波器的作用是利用离散时间系统的特性对输入信号波形(或频谱)进行加工处理,或者说利用数字方法按预定的要求对信号进行变换。数字滤波器可以理解为是一个计算程序或算法,将代表输入信号的数字时间序列转化为代表输出信号的数字时间序
19、列,并在转化过程中,使信号按预定的形式变化。数字滤波器有多种分类,根据数字滤波器冲激响应的时域特征,可将数字滤波器分为两种,即无限长冲激响应(IIR)滤波器和有限长冲激响应(FIR)滤波器。从性能上来说,IIR 滤波器传输函数的极点可位于单位圆内的任何地方,因此可用较低的阶数获得高的选择性,所用的存贮单元少,所以经济而效率高。但是这个高效率是以相位的非线性为代价的。选择性越好,则相位非线性越严重。相反,FIR 滤波器却可以得到严格的线性相位,然而由于 FIR滤波器传输函数的极点固定在原点,所以只能用较高的阶数达到高的选择性;对于同样的滤波器设计指标,FIR 滤波器所要求的阶数可以比 IIR 滤
20、波器高510 倍,结果,成本较高,信号延时也较大;如果按相同的选择性和相同的线性要求来说,则 IIR 滤波器就必须加全通网络进行相位较正,同样要大增加滤波器的节数和复杂性。整体来看,IIR 滤波器达到同样效果阶数少,延迟小,但是有稳定性问题,非线性相位;FIR 滤波器没有稳定性问题,线性相位,但阶数多,延迟大。3.5倒谱的概念定义:倒谱定义为信号短时振幅谱的对数傅里叶反变换。特点:具有可近似地分离并能提取出频谱包络信息和细微结构信息的特点用途:提取声道特征信息:提取频谱包络特征,以此作为描述音韵的特征参数而应用于语音识别。8提取音源信息:提取基音特征,以此作为描述音韵特征的辅助参数而应用于语音
21、识别。求法:D P T L o g I I x D P T D P T峰值检测nx时 间 窗 ABClifter倒 谱 窗 EFA:短时信号;B:短时频谱;C:对数频谱; D:倒谱系数;E:对数频谱包络;F:基本周期4. 语音信号处理实例分析4.1图形用户界面设计在 MATLAB 主窗口中,选择 File 菜单中的 New 菜单项,再选择其中的 GUI命令,就会显示图形用户界面的设计模板。MATLAB 为 GUI 设计一共准备了 4 种模板,分别是 Blank GUI(默认) 、GUI with Uicontrols(带控件对象的 GUI 模板) 、GUI with Axes and Menu
22、(带坐标轴与菜单的 GUI 模板)与 Modal Question Dialog(带模式问话对话框的 GUI 模板)。设计语音信号处理系统的用户图形操作界面(GUI)SoundProcess,其中菜单主要包括 File、Process 和 Output 三大主要部分,其中 File 菜单包括输入(Input) 、保存(Save)和退出(Quit)等功能;Process 菜单主要包括提取(Extract) 、调整(Extract) 、变换(Transform)和滤波(Filter)菜单,其中调整(Extract)包括幅度调整(Range)和频率调整(Frequency) ,滤波(Filter)菜
23、单包含低通滤波(LowpassFilter) 、高通滤波(HighpassFilter) 、带通滤波(BandpassFilter)和帯阻滤波(BandstopFilter)等功能菜单。4.2信号的采集该系统是以一段简短的的语音信号做为分析样本,通过计算机系统将一段“主人,信息收到了”的语音信号保存到到计算机中,并且保存格式为9“*.wav”。4.3 语音信号的处理设计(1)语音信号的提取在 Matlab 中使用 Wavread 函数,可得出信号的采样频率为 22500,并且声音是单声道的。利用 Sound 函数可以清晰的听到“主人,信息收到了”的语音。采集数据并画出波形图。其中声音的采样频率
24、 Fs=22050Hz,y 为采样 数据,NBITS 表示量化阶数。部分程序如下:fn=input( Enter WAV filename:,s); %获取一个*.wav 的文件x,fs,nb=wavread(fn);ms2=floor(fs*0.002);ms10=floor(fs*0.01);ms20=floor(fs*0.02);ms30=floor(fs*0.03);t=(0:length(x)-1)/fs; %计算样本时刻 subplot(2,1,1); %确定显示位置plot(t,x); %画波形图legend(Waveform);xlabel( Time(s);ylabel(Am
25、plitude); 运行后弹出语音信号处理系统的操作界面如图 4-1:图 4-1 语音信号处理系统的操作界面然后点击 File 菜单中的子菜单 Input,回到 Matlab 软件的输入界面如图 4-2:10图 4-2 输入界面输入要处理的语音信号的名称,便可得到语音语音的波形图如图 4-3:图 4-3 语音语音的波形图如图中提取的语音的波形图所示,整段音频数据中得声音高低起伏与录入的声音信号基本一致,并且可以观察到其中包含部分高频噪声。(2)语音信号的调整在语音信号的研究中,经常会对语音信号进行进行多倍频率以及多倍幅度变换调整,日常应用中,这种变换调整也经常要用到。所以在设计中也添加了这种功
26、能,并能够观察调整后的信号的波形图得变化, 而且能通过语音处理界面的输出功能试听处理后的语音信号。语音信号的频率调整11在设计中,可以将语音信号的采样频率提高或降低,来实现语音信号的调整,得到理想的语音信号。例如将采样频率提高一倍,即可得到语音信号频率为原频率 2 倍新的语音信号。运行 ProcessAdjustFrequency,得到如图 4-4 的信号波形图,并试听调整后的效果。图 4-4 频率调整后波形图与原语音信号相比,经过调整后的信号周期变为原来的 1/2,此时的语速明显变快,即实现了信号的 2 倍频功能。语音信号的振幅调整在设计中,可以将语音信号的幅度进行提高或降低操作,来实现语音
27、信号的调整,得到声音音量大小不同的语音信号,例如将原语音信号的幅度提高一倍,得到如下图 4-5 的信号波形图,可以通过 GUI 操作界面的输出功能试听调整后的效果。图 4-5 幅度调整后波形图此时听到的调整后声音声调变高,但不是很明显,可以将幅度的变化值设12置的比较大,那样的话就可以得到效果相当明显的语音信号了。(3)语音信号的傅里叶变换倒谱分析是指信号短时振幅谱的对数进行傅里叶反变换。它具有可近似地分离并提取出频谱包络信息和细微结构信息的特点。对语音信号进行频谱分析,在Matlab中可以利用函数fft对信号行快速傅里叶变换,得到信号的频谱图,并进行倒谱分析,得到倒谱图。傅里叶变换的部分程序
28、如下:x=y(44101:55050,1); %提取原语音信号的一部分t=(0:length(x)-1)/fs; %计算样本时刻subplot(3,1,1); %确定显示位置plot(t,x); %画波形图legend(波形图);xlabel( Time(s);ylabel(Amplitude);Y=fft(x,hamming(length(x); %做加窗傅里叶变换fm=5000*length(Y)/fs; %限定频率范围f=(0:fm)*fs/length(Y); %确定频率刻度subplot(3,1,2);plot(f,20*log10(abs(Y(1:length(f)+eps);le
29、gend(频谱图); %画频谱图ylabel(幅度(db);xlabel(频率(Hz);c=fft(log(abs(x)+eps); %倒频谱计算ms1=fs/1000;ms20=fs/50q=(ms1:ms20)/fs; %确定倒频刻度subplot(3,1,3);plot(q,abs(c(ms1:ms20); %画倒谱图legend(倒谱图);13xlabel(倒频(s));ylabel(倒频谱幅度(Hz));运行Process Transform,对语音信号的一部分进行傅里叶变换,并进行倒谱分析,得到如图4-6:图 4-6 声音样本波形图、频谱图和倒谱图从上面的倒谱图可以看出当读“主人,
30、信息收到了”时,所对应的频率大概在200Hz左右。这与人的语音信号频率集中在200 Hz到4.5 kHz之间是相一致的。而在未发声的时间段内,相对的小高频部分(200500Hz)应该属于背景噪声。(4)语音信号的滤波从图 4-4 中发现,语音信号中包含背景噪声,这些噪声的频率一般较高。所以可以利用 MATLAB 软件中的滤波器进行滤波处理,得到较为理想的语音信号。语音信号的低通滤波系统中设计了一个截止频率为 200Hz 切比雪夫I 型低通滤波器,它的幅频特性如下图 4-7:14图 4-7 低通滤波器的幅频特性低通滤波器性能指标:wp=0.075pi,ws =0.125pi,Rp=0.25;As
31、=50dB;经过低通滤波器处理后,比较处理前后的波形图的变化,如下图 4-8:图 4-8 低通滤波后波形和频谱的变化低通滤波后,声音稍微有些发闷、低沉,原因是高频分量被低通滤波器衰减。但是很接近原来的声音。 语音信号的高通滤波运用切比雪夫型数字高通滤波器,对语音信号进行滤波处理。高通滤波器性能指标:wp=0.375pi,ws=0.425pi,Rp=0.25;As=50dB;然后将其与原信号的比较图如下图 4-9:15图 4-9 高通滤波后波形和频谱的变化高通滤波后,此时只有少许杂音,原因是低频分量被高通滤波器衰减,而人声部分正好是低频部分,所以只剩下杂音,或者发出高频杂音但人的耳朵听不到。语音
32、信号的带通滤波运用椭圆数字带通滤波器函数,对语音信号进行滤波处理后其与原信号的比较图如下图 4-10:图 4-10 带通滤波后波形和频谱的变化16语音信号的带阻滤波运用切比雪夫型数字带阻滤波器,对语音信号进行滤波处理后其与原信号的比较图如下图 4-11:图 4-11 帯阻滤波后波形和频谱的变化从以上各种数字滤波器经过滤波后得出的语音信号相比较,低通滤波后,声音稍微有些发闷,但是很接近原来的声音;高通滤波后听不到人的声音;带通滤波后声音有点像机器人小叮当发出的声音。带阻滤波后,声音比较接近原来的声音。从频谱图中我们可以看出声音的能量主要集中在低频(0.2pi即22045Hz以内)部分。4.4 语
33、音信号的输出可以将处理后的语音信号在 Matlab 软件先播放,体验处理后的语音信号的效果。还可以将处理后的语音信号保存在电脑上。运行 FileSave,保存处理后的语音信号。如果没有语音信号被处理,则系统会出现提示如下图 4-12:17图 4-12 保存提示界面如果有语音信号被处理,运行 FileSave,系统会出现提示如下图 4-13:图 4-13 保存界面保存后,整个操作过程就完成了。18参考文献1 余胜威,吴婷,罗建桥MATLAB GUI 设计入门与实践北京:清华大学出版社,20162 高西全,丁玉美数字信号处理第 3 版西安:西安电子科技大学出版社,20083 罗华飞MATLAB G
34、UI 设计学习手记第 3 版北京:北京航空航天大学出版社,20144 刘泉,阙大顺数字信号处理原理与实现北京:电子工业出版社,20105 张磊,毕靖,郭莲英MATLAB实用教程北京:人民邮电出版社,20086 张威MATLAB基础与编程入门西安:西安电子科技大学出版社,20137 刘帅奇,李会雅,赵杰MATLAB 程序设计基础与应用北京:清华大学出版社,20168 丁伟雄MATLAB R2015a 数字图像处理北京:清华大学出版社,20169 刘维精通 MATLAB 与 C/C+混合程序设计第 4 版北京:北京航空航天大学出版社,201519附录(I) 设计 FIR和 IIR数字滤波器%=II
35、R 低通滤波器=Ft=8000;Fp=1000;Fs=1200;wp=2*pi*Fp/Ft;ws=2*pi*Fs/Ft;fp=2*Ft*tan(wp/2);fs=2*Fs*tan(wp/2);n11,wn11=buttord(wp,ws,1,50,s); b11,a11=butter(n11,wn11,s); num11,den11=bilinear(b11,a11,0.5);h,w=freqz(num11,den11);figure;plot(w*8000*0.5/pi,abs(h);legend(IIR 低通滤波器,Location,NorthWest);grid;程序结果如下图:0 50
36、0 1000 1500 2000 2500 3000 3500 400000.20.40.60.811.21.4 IIR低低低低低20%=IIR 带通滤波器=Fp1=1200;Fp2=3000;Fs1=1000;Fs2=3200;Ft=8000;wp1=tan(pi*Fp1/Ft);wp2=tan(pi*Fp2/Ft);ws1=tan(pi*Fs1/Ft);ws2=tan(pi*Fs2/Ft);w=wp1*wp2/ws2;bw=wp2-wp1; %有效通带频率wp=1;ws=(wp1*wp2-w.2)/(bw*w);n12,wn12=buttord(wp,ws,1,50,s); b12,a12
37、=butter(n12,wn12,s);num2,den2=lp2bp(b12,a12,sqrt(wp1*wp2),bw);num12,den12=bilinear(num2,den2,0.5);h,w=freqz(num12,den12);figure;plot(w*8000*0.5/pi,abs(h);axis(0 4500 0 1.5);legend(IIR 带通滤波器,Location,NorthWest);grid;程序结果如下图:210 500 1000 1500 2000 2500 3000 3500 4000 450000.511.5 IIR低低低低低%=IIR 高通滤波器=F
38、t=8000;Fp=4000;Fs=3500;wp1=tan(pi*Fp/Ft);ws1=tan(pi*Fs/Ft);wp=1;ws=wp1*wp/ws1;n13,wn13=cheb1ord(wp,ws,1,50,s);b13,a13=cheby1(n13,1,wn13,s); num,den=lp2hp(b13,a13,wn13);num13,den13=bilinear(num,den,0.5); h,w=freqz(num13,den13);figure;plot(w*21000*0.5/pi,abs(h);legend(IIR 高通滤波器,Location,NorthWest);axi
39、s(0 11000 0 1.5);grid;22程序结果如下图:0 1000 2000 3000 4000 5000 6000 7000 8000 9000 100001100000.511.5 IIR低低低低低%*FIR 低通滤波*Ft=8000; Fp=1000; Fs=1200; wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.(-rp/20); s=10.(-rs/20);fpts=wp ws;mag=1 0;dev=p s;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);b21=fir1(n21,wn21,kai
40、ser(n21+1,beta);h,w=freqz(b21,1); figure;23plot(w*8000*0.5/pi,abs(h);title(FIR 低通滤波器,fontweight,bold);grid;程序结果如下图:0 500 1000 1500 2000 2500 3000 3500 400000.20.40.60.811.21.4 FIR低低低低低%*FIR 带通滤波器*Fp1=1200; %通带边界频率 Fp2=3000;Fs1=1000; %阻带截止频率Fs2=3200;Ft=8000;wp1=tan(pi*Fp1/Ft); wp2=tan(pi*Fp2/Ft);ws1=
41、tan(pi*Fs1/Ft);ws2=tan(pi*Fs2/Ft); w=wp1*wp2/ws2;bw=wp2-wp1;wp=1;ws=(wp*wp2-w.2)/(bw*w);24n22,wn22=buttord(wp,ws,1,50,s); b22,a22=butter(n22,wn22,s); num2,den2=lp2bp(b22,a22,sqrt(wp1*wp2),bw); num22,den22=bilinear(num2,den2,0.5);h,w=freqz(num22,den22); figure;plot(w*8000*0.5/pi,abs(h);axis(0 4500 0
42、1.5);legend(FIR 带通滤波器,Location,NorthWest);grid;程序结果如下图:0 500 1000 1500 2000 2500 3000 3500 4000 450000.511.5FIR低低低低低%*%FIR 高通滤波器 *Ft=8001;Fp=4000; %通带边界频率 Fs=3500; %阻带截止频率wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;25rs=50;p=1-10.(-rp/20); s=10.(-rs/20);fpts=ws wp;mag=0 1;dev=p s;n23,wn23,beta,ftype=kaiserord(fpts,m
43、ag,dev);b23=fir1(n23,wn23,high,kaiser(n23+1,beta);h,w=freqz(b23,1); figure;plot(w*12000*0.5/pi,abs(h);title(FIR 高通滤波器,fontweight,bold);axis(2500 5500 0 1.2);grid;程序结果如下图:2500 3000 3500 4000 4500 5000 550000.20.40.60.81FIR低低低低低26附录(II)比较滤波前后语音信号的波形及频谱% =双线性变换法=%*低通滤波器*y,fs,nbits=wavread (OriSound); %
44、IIR 低通n = length (y) ; %求出语音信号的长度Noise=0.2*randn(n,2); %随机函数产生噪声s=y+Noise; %语音信号加入噪声S=fft(s); Ft=8000;Fp=1000;Fs=1200;wp=2*pi*Fp/Ft;ws=2*pi*Fs/Ft;n11,wn11=buttord(wp,ws,1,50,s);%求低通滤波器的阶数和截止频率b11,a11=butter(n11,wn11,s); %求 S 域的频率响应的参数 num11,den11=bilinear(b11,a11,0.5); %利用双线性变换实现频率响应 S 域到 Z域的变换 z11=
45、filter(num11,den11,s);sound(z11);m11=fft(z11); %求滤波后的信号figure;subplot(2,2,1);plot(abs(S),g);title(滤波前信号的频谱,fontweight,bold);axis( 0 150000 0 4000);grid;subplot(2,2,2);plot(abs(m11),r);title(滤波后信号的频谱,fontweight,bold);27axis( 0 150000 0 4000);grid;subplot(2,2,3);plot(s);title(滤波前信号的波形,fontweight,bold)
46、;axis(95000 100000 -1 1);grid;subplot(2,2,4);plot(z11);title(滤波后的信号波形,fontweight,bold);axis(95000 100000 -1 1);grid;程序结果如下图:0 5 10 15x 10401000200030004000 低 低 低 低 低 低 低 低0 5 10 15x 10401000200030004000 低 低 低 低 低 低 低 低9.5 9.6 9.7 9.8 9.9 10x 104-1-0.500.51 低 低 低 低 低 低 低 低9.5 9.6 9.7 9.8 9.9 10x 104-
47、1-0.500.51 低 低 低 低 低 低 低 低附 II-1 双线性低通滤波器比较%*带通滤波器*y,fs,nbits=wavread (OriSound); %IIR 带通28n = length (y) ; %求出语音信号的长度Noise=0.2*randn(n,2); %随机函数产生噪声s=y+Noise; %语音信号加入噪声 S=fft(s); %傅里叶变换Ft=8000;Fp=1000;Fs=1200;wp=2*Fp/Ft;ws=2*Fs/Ft;rp=1;rs=50;p=1-10.(-rp/20); %通带阻带波纹q=10.(-rs/20);fpts=wp ws;mag=1 0;dev=p q;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);%由