在Matlab平台上实现对语音信号的去噪研究和仿真
摘要
语音信号在数字信号处理中占有极其重要的地位，因此选择通过对语音信号的研究来巩固和掌握数字信号处理的基本能力十分具有代表性。对数字信号处理离不开滤波器，因此滤波器的设计在信号处理中占有极其重要的地位。而MATLAB软件工具箱提供了对各种数字滤波器的设计。本论文“在MATLAB平台上实现对语音信号的去噪研究与仿真”综合运用了数字信号处理的各种基本知识，进而对不带噪语音信号进行谱分析以及带噪语音信号进行谱分析和滤波处理。通过理论推导得出相应的结论，再通过利用MATLAB作为编程工具来进行计算机实现比价已验证推导出来的结论。在设计过程中，通过设计FIR数字滤波器和IIR数字滤波器来完成滤波处理。在设计过程中，运用了MATLAB对整个设计中的图形的绘制和一些数据的计算以及仿真。
关键字 滤波器；MATLAB；仿真；滤波
目录
第1章 前言
1.1 研究的意义
语音是语言的声学表现，是人类交流信息最自然、最有效、最方便的手段。随着社会文化的进步和科学技术的发展，人类开始进入了信息化时代，用现代手段研究语音处理技术，使人们能更加有效地产生、传输、存储、和获取语音信息，这对于促进社会的发展具有十分重要的意义，因此，语音信号处理正越来越受到人们的关注和广泛的研究。
语音信号是信息技术处理中最重要的一门科学，是人类社会几步的标志。那么什么是语音？语音是人类特有的功能，也是人类获取外界信息的重要工具，也是人与人交流必不可少的重要手段。那么什么又是信号？那信号是什么呢？信号是传递信息的函数。离散时间信号——序列——可以用图形来表示。
语音信号处理是一门用研究数字信号处理研究信号的科学。它是一新兴的信
息科学，同时又是综合多个学科领域的一门交叉科学。语音在我们的日常生
活中随时可见，也随处可见，语音很大程度上可以影响我们的生活。所以研究语音信号无论是在科学领域上还是日常生活中都有其广泛而重要的意义。
1.2 国内外研究现状
20世纪60年代中期形成的一系列数字信号处理的理论和算法，如数字滤波器、快速傅立叶变换（FFT）等是语音信号数字处理的理论和技术基础。随着信息科学技术的飞速发展，语音信号处理取得了重大的进展：进入70年代之后，提出了用于语音信号的信息压缩和特征提取的线性预测技术（LPC），并已成为语音信号处理最强有力的工具，广泛应用于语音信号的分析、合成及各个应用领域，以及用于输入语音与参考样本之间，时间匹配的动态规划方法；80年代初一种新的基于聚类分析的高效数据压缩技术—矢量量化（VQ）应用于语音信号处理中；而用隐马尔可夫模型（HMM）描述语音信号过程的产生是80年代语音信号处理技术的重大发展，目前HMM已构成了现代语音识别研究的重要基石。近年来人工神经网络(ANN)的研究取得了迅速发展，语音信号处理的各项课题是促进其发展的重要动力之一，他的各项成果也体现在语音信号处理的各项技术之中。
1.3 研究的内容
本论文主要介绍的是的语音信号的简单处理。本论文针对以上问题，运用数字信号学基本原理实现语音信号的处理，在matlab7.0环境下综合运用信号提取，幅频变换以及傅里叶变换、滤波等技术来进行语音信号处理。我所做的工作就是在matlab7.0软件上编写一个处理语音信号的程序，能对语音信号进行采集，并对其进行各种处理，达到简单语音信号处理的目的。
对语音信号的研究，本论文采用了设计两种滤波器的基本研究方法来达到研究语音信号去噪的目的，最终结合图像以及对语音信号的回放，通过对比，得出结论。
本课题的研究基本步骤如下：
- 语音信号的录制。
- 在MATLAB平台上读入语音信号。
- 绘制频谱图并回放原始语音信号。
- 利用MATLAB编程加入一段正弦波噪音，设计滤波器去噪。
- 利用MATLAB编程加入一段随机噪音信号，设计FIR和IIR滤波器去噪，并分别绘制频谱图、回放语音信号。
6 通过仿真后的图像以及对语音信号的回放，对比两种去噪方式的优缺点。
其大概流程框图可如下表示：（图1-1）
图1-1 论文设计的流程
第2章 语音信号去噪方法的研究
2.2 去噪的原理
2.2.1 采样定理
在进行模拟/数字信号的转换过程中，当采样频率fs.max大于信号中，最高频率fmax的2倍时，即：fs.max>=2fmax,则采样之后的数字信号完整地保留了原始信号中的信息，一般实际应用中保证采样频率为信号最高频率的5～10倍；采样定理又称奈奎斯特定理。 1924年奈奎斯特(Nyquist)就推导出在理想低通信道的最高大码元传输速率的公式: 理想低通信道的最高大码元传输速率=2W*log2 N (其中W是理想低通信道的带宽,N是电平强度)为什么把采样频率设为8kHz?在数字通信中，根据采样定理, 最小采样频率为语音信号最高频率的2倍
频带为F的连续信号 f(t)可用一系列离散的采样值f(t1),f(t1±Δt)，f(t1±2Δt)，...来表示,只要这些采样点的时间间隔Δt≤1/2F，便可根据各采样值完全恢复原来的信号f(t)。 这是时域采样定理的一种表述方式。
时域采样定理的另一种表述方式是：当时间信号函数f(t)的最高频率分量为fM时,f(t)的值可由一系列采样间隔小于或等于1/2fM的采样值来确定,即采样点的重复频率f≥2fM。图为模拟信号和采样样本的示意图。
时域采样定理是采样误差理论、随机变量采样理论和多变量采样理论的基础。对于时间上受限制的连续信号f(t)（即当│t│>T时,f(t)=0,这里T=T2-T1是信号的持续时间），若其频谱为F（ω）,则可在频域上用一系列离散的采样值
（2-1）
采样值来表示,只要这些采样点的频率间隔
（2-2）
。
2.2.2 采样频率
采样频率，也称为采样速度或者采样率，定义了每秒从连续信号中提取并组成离散信号的采样个数，它用赫兹（Hz）来表示。采样频率的倒数是采样周期或者叫作采样时间，它是采样之间的时间间隔。通俗的讲采样频率是指计算机每秒钟采集多少个声音样本，是描述声音文件的音质、音调，衡量声卡、声音文件的质量标准。
采样频率只能用于周期性采样的采样器，对于非周期性采样的采样器没有规则限制。 采样频率的常用的表示符号是 fs。 通俗的讲采样频率是指计算机每秒钟采集多少个声音样本，是描述声音文件的音质、音调，衡量声卡、声音文件的质量标准。采样频率越高，即采样的间隔时间越短，则在单位时间内计算机得到的声音样本数据就越多，对声音波形的表示也越精确。采样频率与声音频率之间有一定的关系，根据采样定理，只有采样频率高于声音信号最高频率的两倍时，才能把数字信号表示的声音还原成为原来的声音。这就是说采样频率是衡量声卡采集、记录和还原声音文件的质量标准。
采样位数和采样率对于音频接口来说是最为重要的两个指标，也是选择音频接口的两个重要标准。无论采样频率如何，理论上来说采样的位数决定了音频数据最大的力度范围。每增加一个采样位数相当于力度范围增加了6dB。采样位数越多则捕捉到的信号越精确。对于采样率来说你可以想象它类似于一个照相机，44.1kHz意味着音频流进入计算机时计算机每秒会对其拍照达441000次。显然采样率越高，计算机摄取的图片越多，对于原始音频的还原也越加精确
2.2 去噪的方法
数字信号处理技术经过几十年的发展，在国内外已经取得了很大的成绩。到目前为止，已经比较成熟的去噪方法比较典型的有切比雪夫去噪法、双线性变换去噪法、窗函数去噪法、谱相减去噪法、巴特沃兹去噪法等有名的去噪方法。下面分别对上述去噪方法中比较有代表性的谱相减去噪法和FIR以及IIR滤波法做一个简单的介绍。
（1）谱相减算法去噪法
在我们的日常交流和语音通信系统中，加性宽带噪声严重影响了语音质量和可懂度。从带噪语音中提取原始语音信号的方法很多，在单信道条件下，谱相减算法以其运算量小、原理简单、易于实现并且有不错的增强效果而得到了广泛的应用。谱相减语音增强算法的核心是噪声检测和谱减规则。在分析了语音增强算法理论的基础上，本文首先研究了语音激活检测算法。对基于短时能量和短时过零率双门限法语音激活检测的噪声估计算法做了研究及仿真，同时还研究了一种基于最小子带能量的噪声估计方法。然后，通过分析经典谱减法的原理及其一般改进形式，研究了一种基于噪声残差的谱相减改进算法和一种可以不以噪声是零均值的高斯分布为前提的谱减法改进算法。最后通过大量的仿真实验，验证了所研究的几种改进算法都能有效地提高增强效果。
由于谱相减算法去噪比较复杂，我们在这里就不多做讨论，并且本论文也不采用此种方法。下面着重讲解本论文采用的两种方法。
（2）FIR滤波法
FIR工作原理：在进入FIR滤波器前，首先要将信号通过A/D器件进行模数转换，使之成为8bit的数字信号，一般可用速度较高的逐次逼进式A/D转换器，不论采用乘累加方法还是分布式算法设计FIR滤波器，滤波器输出的数据都是一串序列，要使它能直观地反应出来，还需经过数模转换，因此由FPGA构成的FIR滤波器的输出须外接D/A模块。FPGA有着规整的内部逻辑阵列和丰富的连线资源，特别适合于数字信号处理任务，相对于串行运算为主导的通用DSP芯片来说，其并行性和可扩展性更好，利用FPGA乘累加的快速算法，可以设计出高速的FIR数字滤波器。
FIR的特点：
有限长单位冲激响应（FIR）滤波器有以下特点：
A位冲激响应h (n)在有限个n值处不为零
B系统函数0处收敛，极点全部在z = 0处
C结构上主要是非递归结构，没有输出到输入的反馈，但有些结构中（例如频率抽样结构）也包含有反馈的递归部分。
设FIR滤波器的单位冲激响应h (n)为一个N点序列，0 ≤ n ≤N —1，则滤波器的系统函数为
H(z)=∑h(n)*z^-n （2-3）
就是说，它有（N—1）阶极点在z = 0处，有（N—1）个零点位于有限z平面的任何位置。
FIR滤波器基本结构：
横截型
式的系统的差分方程表达式为
y(n)=∑h(m)x(n-m) ( 2-4）
很明显，这就是线性移不变系统的卷积和公式，也是x (n)的延时链的横向结构，称为横截型结构或卷积型结构，也可称为直接型结构。将转置定理用于可得到的转置直接型结构。
级联型
其中[N/2]表示取N/2的整数部分。若N为偶数，则N—1为奇数，故系数B2K中有一个为零，这是因为，这时有奇数个根，其中复数根成共轭对必为偶数，必然有奇数个实根。画出N为奇数时，FIR滤波器的级联结构，其中每一个二阶因子用图4-11的横型结构。
这种结构的每一节控制一对零点，因而再需要控制传输零点时，可以采用它。但是这种结构所需要的系数B2k（I = 0，1，2，k，= 1，2，．．．，[N/2]）比卷积型的系数h (n)要多，因而所需的乘法次数也比卷积型的要多。
IIRnfinite Impulse Response)数字滤波器，又名“无限脉冲响应数字滤波器”，或“递归滤波器”。递归滤波器，也就是IIR数字滤波器，顾名思义，具有反馈，一般认为具有无限的脉冲响应。
IIR数字滤波器的设计
利用 MATALAB工具箱分析工具(FDATool)可以很方便地设计出符合应用要求的未经量化的IIR数字滤波器。需要将MATLAB设计出的IIR数字滤波器进一步分解和量化，从而获得可用FPGA实现的滤波器系数。
IIR数字滤波器的设计步骤
由于采用了级联结构，因此如何将滤波器的每一个极点和零点相组合，从而使得数字滤波器输出所含的噪声最小是个十分关键的问题。为了产生最优的量化后的IIR数字滤波器，采用如下步骤进行设计。
首先计算整体传递函数的零极点；
选取具有最大幅度的极点以及距离它最近的零点，使用它们组成一个二阶基本节的传递函数；
对于剩下的极点和零点采用与相类似的步骤，直至形成所有的二阶基本节。
通过上面三步法进行的设计可以保证IIR数字滤波器中N位乘法器产生的量化舍入误差最小。
获得最优IIR数字滤波器系数：
为了设计出可用FPGA实现的数字滤波器，需要对上一步分解获得的二阶基本节的滤波器系数进行量化，即用一个固定的字长加以表示。量化过程中由于存在不同程度的量化误差，由此会导致滤波器的频率响应出现偏差，严重时会使IIR滤波器的极点移到单位圆之外，系统因而失去稳定性。为了获得最优的滤波器系数，采用以下步骤进行量化。
A计算每个系数的 绝对值；
B查找出每个系数绝对值中的最大值；
C计算比此绝对值大的最小整数；
D对3.2.3的结果取反获得负整数；
E计算需要表示此整数的最小位数；
F计算用于表示系数值分数部分的余下位数。
除了系数存在量化误差，数字滤波器运算过程中有限字长效应也会造成误差，因此对滤波器中乘法器、加法器及寄存器的数据宽度要也进行合理的设计，以防止产生极限环现象和溢出振荡。
与FIR数字滤波器的设计不同，IIR滤波器设计时的阶数不是由设计者指定，而是根据设计者输入的各个滤波器参数（截止频率、通带滤纹、阻带衰减等），由软件设计出满足这些参数的最低滤波器阶数。在MATLAB下设计不同类型IIR滤波器均有与之对应的函数用于阶数的选择。
第3章 滤波器的设计及实现
3.1数字滤波器设计的基本原理
滤波器是一种对信号有处理作用的器件或电路。
滤波器分为有源滤波和无源滤波，它的主要作用是：让有用信号尽可能无衰减的通过，对无用信号尽可能大的衰减。
滤波器是由电感器和电容器构成的网路，可使混合的交直流电流分开。整流器中，即借助此网路滤净脉动直流中的涟波，而获得比较纯净的直流输出。最基本的滤波器，是由一个电容器和一个电感器构成，称为L型滤波。所有各型的滤波器，都是集合L型单节滤波器而成。基本单节式滤波器由一个串联臂及一个并联臂所组成，串联臂为电感器，并联臂为电容器。在电源及声频电路中之滤波器，最通用者为L型及π型两种。就L型单节滤波器而言，其电感抗XL与电容抗XC，对任一频率为一常数，其关系为
XL·XC=K2 (3-1)
故L型滤波器又称为K常数滤波器。倘若一滤波器的构成部分，较K常数型具有较尖锐的截止频率（即对范围选择性强），而同时对此截止频率以外的其他频率只有较小的衰减率者，称为m常数滤波器。所谓截止频率，亦即与滤波器有尖锐谐振的频率。通带与带阻滤波器都是m常数滤波器，m为截止频率与被衰减的其他频率之衰减比的函数。每一m常数滤波器的阻抗与K常数滤波器之间的关系，均由m常数决定，此常数介于0～1之间。当m接近零值时，截止频率的尖锐度增高，但对于截止频的倍频之衰减率将随着而减小。最合于实用的m值为0.6。至于那一频率需被截止，可调节共振臂以决定之。m常数滤波器对截止频率的衰减度，决定于共振臂的有效Q值之大小。若达K常数及m常数滤波器组成级联电路，可获得尖锐的滤波作用及良好的频率衰减。
简而言之，滤波，本质上是从被噪声畸变和污染了的信号中提取原始信号所携带的信息的过程。
滤波器的分类：
巴特沃斯响应（最平坦响应），巴特沃斯响应能够最大化滤波器的通带平坦度。该响应非常平坦，非常接近DC信号，然后慢慢衰减至截止频率点为-3dB，最终逼近-20ndB/decade的衰减率，其中n为滤波器的阶数。巴特沃斯滤波器特别适用于低频应用，其对于维护增益的平坦性来说非常重要。
贝塞尔响应，除了会改变依赖于频率的输入信号的幅度外，滤波器还会为其引入了一个延迟。延迟使得基于频率的相移产生非正弦信号失真。就像巴特沃斯响应利用通带最大化了幅度的平坦度一样，贝塞尔响应最小化了通带的相位非线性。
切贝雪夫响应，在一些应用当中，最为重要的因素是滤波器截断不必要信号的速度。如果你可以接受通带具有一些纹波，就可以得到比巴特沃斯滤波器更快速的衰减。
本课题主要采样FIR和IIR数字滤波器对加噪语音信号进行滤波。下面就分别论述一下FIR和IIR数字滤波器的设计原理和基本步骤。
3.2 FIR数字滤波器的设计及实现
FIR：有限脉冲响应滤波器。有限说明其脉冲响应是有限的。与IIR相比，它具有线性相位、容易设计的优点。这也就说明，IIR滤波器具有相位不线性，不容易设计的缺点。而另一方面，IIR却拥有FIR所不具有的缺点，那就是设计同样参数的滤波器，FIR比IIR需要更多的参数。这也就说明，要增加DSP的计算量。DSP需要更多的计算时间，对DSP的实时性有影响。FIR滤波器的设计比较简单，就是要设计一个数字滤波器去逼近一个理想的低通滤波器。通常这个理想的低通滤波器在频域上是一个矩形窗。根据傅里叶变换我们可以知道，此函数在时域上是一个采样函数。通常此函数的表达式为：
sa（n）＝sin（n∩）/n（3-2）
但是这个采样序列是无限的，计算机是无法对它进行计算的。故我们需要对此采样函数进行截断处理。也就是加一个窗函数。就是传说中的加窗。也就是把这个时域采样序列去乘一个窗函数，就把这个无限的时域采样序列截成了有限个序列值。但是加窗后对此采样序列的频域也产生了影响：此时的频域便不再是一个理想的矩形窗，而是成了一个有过渡带，阻带有波动的低通滤波器。通常根据所加的窗函数的不同，对采样信号加窗后，在频域所得的低通滤波器的阻带衰减也不同。通常我们就是根据此阻带衰减去选择一个合适的窗函数。如矩形窗、汉宁窗、汉明窗、BLACKMAN窗、凯撒窗等。
面我们着重介绍用窗函数法设计FIR滤波器的步骤。如下：
（1）根据对阻带衰减及过渡带的指标要求，选择串窗数类型（矩形窗、三角窗、汉宁窗、哈明窗、凯塞窗等），并估计窗口长度N。先按照阻带衰减选择窗函数类型。原则是在保证阻带衰减满足要求的情况下，尽量选择主瓣的窗函数。
（2）构造希望逼近的频率响应函数。
（3）计算h(n).。
（4）加窗得到设计结果。
接下来，我们根据语音信号的特点给出有关滤波器的技术指标：
低通滤波器的性能指标：
fp=1000Hz，fc=1200Hz，As=50db ,Ap=1dB
高通滤波器的性能指标：
fp=3500Hz，fc=4000Hz，As=50dB，Ap=1dB
在Matlab中，可以利用函数fir1设计FIR滤波器，利用Matlab中的函数freqz画出各步步器的频率响应。
MATLAB信号处理工具箱函数cheblap,cheblord和cheeby1是切比雪夫I型滤波器设计函数。我们用到的是cheeby1函数，其调用格式如下：
[B,A]=cheby1(N,Rp,wpo,’ftypr’)
[B,A]=cheby1(N,Rp,wpo,’ftypr’,’s’)
下面我们将给出设计FIR数字滤波器的主要程序和图像。
FIR低通滤波器程序见附录1
FIR低通滤波器图像：（图3—1）
图3—1 FIR低通滤波器
FIR高通滤波程序见附录2
FIR高通滤波图像：（图3-2）
图3-2 FIR高通滤波器
3.3 IIR数字滤波器的设计及实现
对于数字高通、带通滤波器的设计，通用方法为双线性变换法。可以借助于模拟滤波器的频率转换设计一个所需类型的过渡模拟滤波器，再经过双线性变换将其转换策划那个所需的数字滤波器。具体设计步骤如下：
（1）确定所需类型数字滤波器的技术指标。
（2）将所需类型数字滤波器的边界频率转换成相应的模拟滤波器的边界频率，转换公式为Ω=2/T tan(0.5ω) （3-1）
(3)将相应类型的模拟滤波器技术指标转换成模拟低通滤波器技术指标。
（4）设计模拟低通滤波器。
（5）通过频率变换将模拟低通转换成相应类型的过渡模拟滤波器。
（6）采用双线性变换法将相应类型的过渡模拟滤波器转换成所需类型的数字滤波器。
我们知道，脉冲响应不变法的主要缺点是会产生频谱混叠现象，使数字滤波器的频响偏离模拟滤波器的频响特性。为了克服之一缺点，可以采用双线性变换法。
下面我们总结一下利用模拟滤波器设计IIR数字低通滤波器的步骤：
（1）确定数字低通滤波器的技术指标：通带边界频率、通带最大衰减，阻带截止频率、阻带最小衰减。
（2）将数字低通滤波器的技术指标转换成相应的模拟低通滤波器的技术指标。
（3）按照模拟低通滤波器的技术指标设计及过渡模拟低通滤波器。
（4）用双线性变换法，模拟滤波器系统函数转换成数字低通滤波器系统函数。
MATLAB信号处理工具箱函数cheblap,cheblord和cheeby1是切比雪夫I型滤波器设计函数。我们用到的是cheeby1函数，其调用格式如下：
[B,A]=cheby1(N,Rp,wpo,’ftypr’)
[B,A]=cheby1(N,Rp,wpo,’ftypr’,’s’)
函数butter,cheby1和ellip设计IIR滤波器时都是默认的双线性变换法，所以在设计滤波器时只需要代入相应的实现函数即可。下面我们将给出IIR数字滤波器的主要程序。
IIR低通滤波器程序见附录3；
IIR低通滤波器图像：（图3—3）
图3-3 IIR低通滤波器
IIR滤波器高通程序
IIR滤波器高通图像：（图3-4）
图3-4 IIR高通滤波器
第四章 去噪及仿真的研究
4.1 语音文件在MATLAB平台上的录入与打开
单击自己的电脑开始程序，选择所有程序，接着选择附件，再选择娱乐，最后选择录音。自己录入“毕业设计”语音信号，然后保存在MATLAB文件夹里面，命名为“wangqingtian.wav”。
利用MATLAB中的wavread命令来读入（采集）语音信号，将它赋值给某一向量。
[y,fs,bits]=wavread(' [N1 N2]);用于读取语音，采样值放在向量y中，fs表示采样频率(Hz)，bits表示采样位数。[N1 N2]表示读取从N1点到N2点的值（若只有一个N的点则表示读取前N点的采样值）。
4.2 原始语音信号频谱分析及仿真
利用MATLAB中的wavread命令来读入（采集）语音信号，将它赋值给某一向量。再将该向量看作一个普通的信号，对其进行FFT变换实现频谱分析，再依据实际情况对它进行滤波。对于波形图与频谱图（包括滤波前后的对比图）都可以用 MATLAB画出。我们还可以通过sound命令来对语音信号进行回放，以便在听觉上来感受声音的变化。选择设计此方案，是对数字信号处理的一次实践。在数字信号处理的课程学习过程中，我们过多的是理论学习，几乎没有进行实践方面的运用。这个课题正好是对数字语音处理的一次有利实践，而且语音处理也可以说是信号处理在实际应用中很大众化的一方面。这个方案用到的软件也是在数字信号处理中非常通用的一个软件——MATLAB软件。所以这个课题的设计过程也是一次数字信号处理在MATLAB中应用的学习过程。课题用到了较多的MATLAB语句，而由于课题研究范围所限，真正与数字信号有关的命令函数却并不多。sound(x,fs,bits); 用于对声音的回放。向量y则就代表了一个信号（也即一个复杂的“函数表达式”）也就是说可以像处理一个信号表达式一样处理这个声音信号。
FFT的MATLAB实现：
在MATLAB的信号处理工具箱中函数FFT和IFFT用于快速傅立叶变换和逆变换。下面介绍这些函数。
函数FFT用于序列快速傅立叶变换。
函数的一种调用格式为 y=fft(x)
其中，x是序列，y是序列的FFT，x可以为一向量或矩阵，若x为一向量，y是x的FFT。且和x相同长度。若x为一矩阵，则y是对矩阵的每一列向量进行FFT。
如果x长度是2的幂次方，函数fft执行高速基－2FFT算法；否则fft执行一种混合基的离散傅立叶变换算法，计算速度较慢。
函数FFT的另一种调用格式为y=fft(x,N)
式中，x，y意义同前，N为正整数。函数执行N点的FFT。若x为向量且长度小于N，则函数将x补零至长度N。若向量x的长度大于N，则函数截短x使之长度为N。若x 为矩阵，按相同方法对x进行处理。
经函数fft求得的序列y一般是复序列，通常要求其幅值和相位。MATLAB提供求复数的幅值和相位函数：abs，angle，这些函数一般和 FFT同时使用。
函数abs(x)用于计算复向量x的幅值，函数angle(x)用于计算复向量的相角，介于 和 之间，以弧度表示。函数unwrap(p)用于展开弧度相位角p ,当相位角绝对变化超过 时，函数把它扩展至 。
用MATLAB工具箱函数fft进行频谱分析时需注意：
（1）函数fft返回值y的数据结构对称性
一般而言，对于N点的x(n)序列的FFT是N点的复数序列，其点n=N/2+1对应Nyquist频率，作频谱分析时仅取序列X(k)的前一半，即前N/2点即可。X(k)的后一半序列和前一半序列时对称的。
（2）频率计算
若N点序列x(n)(n=0,1,…,N-1)是在采样频率 下获得的。它的FFT也是N点序列，即X(k)(k=0,1,2,…,N-1)，则第k点所对应实际频率值为f=k*f /N.
（3）作FFT分析时，幅值大小与FFT选择点数有关，但不影响分析结果。
下面的一段程序是语音信号在MATLAB中的最简单表现，它实现了语音的读入打开，以及绘出了语音信号的波形频谱图。
[x,fs,bits]=wavread('wangqingtian.wav’);
sound(x,fs,bits);
X=fft(x,4096);
magX=abs(X);
angX=angle(X);
subplot(221);plot(x);title('原始信号波形');
subplot(222);plot(X); title('原始语音信号采样后的频谱图‘)
subplot(223);plot(magX);title('原始信号幅值');
subplot(224);plot(angX);title('原始信号相位');
程序运行可以听到声音，得到的图形为：（图4-1、图4-2、图4-3、图4-4）
图4-1 原始信号波形
图4-3 原始信号幅值
图4-4 原始信号相位
4.3 加噪语音信号频谱分析及仿真
(1)正弦波信号加入原始语音信号
前面已经介绍了MATLAB软件相关知识，那么我们怎么在NATLAB平台上实现对一段原始语音信号加入一个正弦波信号呢？
下面一段程序实现了在原始语音信号加入正弦波信号。程序见附录5
分析此段程序可知，此程序是先对原始语音信号做时域波形分析和频谱分析，然后再对加噪的语音信号做时域波形分析和频谱分析。
首先通过MATLAB中调用和回放语音信号命令来实现对原始语音信号的调用和回放，程序如下：
[y,fs,bits]=wavread('wangqingtian.wav');
sound(y,fs)
由于在MATLAB中，如要实现两个信号的相加减，那么两个信号的长度和维度都要一样才能相加减。程序中：n=length(y) 用于计算信号的长度和选取变换的点数。然后用傅里叶变换到频域：y_p=fft(y,n);通过函数 f=fs*(0:n/2-1)/n;计
图4—5 原始语音信号采样后时域波形和频谱图
上段程序中，函数noise是频率为3000Hz的正弦波信号噪音，语句y_z=y+noise实现了两个信号的相加，然后绘制加噪后的语音信号时域波形和频谱图并回放加噪后的语音信号。加噪后的时域波形和频谱图如下：（图4—6）
图4-6 加噪语音信号时域波形和频谱图
如上所示，通过对加噪前和加噪后语音信号的图像的对比和对语音信号回放的人耳感知可以知道，加入正弦波信号后频谱图和时域波形并没有什么明显的变化，而人耳听到的声音也几乎没有什么变化。
（2）随机噪音信号加入原始语音信号
matlab函数randn：产生正态分布的随机数或矩阵的函数
产生均值为0，方差 σ^2 = 1，标准差σ = 1的正态分布的随机数或矩阵的函数。
用法：
Y = randn(n)
返回一个n*n的随机项的矩阵。如果n不是个数量，将返回错误信息。
Y = randn(m,n) 或 Y = randn([m n])
返回一个m*n的随机项矩阵。
Y = randn(m,n,p,...) 或 Y = randn([m n p...])
产生随机数组。
Y = randn(size(A))
返回一个和A有同样维数大小的随机数组。
randn
返回一个每次都变化的数量。
下面一段程序实现了利用randn函数把一段随机噪音信号加入原始语音信号的信号处理过程：见附录6
图4-7 加噪语音信号时域波形和频谱图
通过对两张图片的对比，很明显可以看加噪后的语音信号时域波形比原始语音信号浑浊了许多，在时间轴上可以明显看出0—0.5S的幅值增大了；通过对原始语音信号的频谱图与加噪后的语音信号频谱图的对比，也可以看出在频率5000Hz以后的频率幅值发生了明显的增加。
再通过对原始语音信号的回放效果与加噪后的语音信号回放的效果的对比，人耳可以明显辨别出两种语音信号不一样了，加噪后的语音信号在听觉上比原始语音信号要浑浊很多，而且还有吱吱嘎嘎的混杂音。
4.4 去噪及仿真
（1）FIR滤波器法去噪
通过对上一节中加噪语音信号和原始语音信号频谱图对比可以知道，噪音大部分是Hz大于5000的部分，故设计低通滤波器进行滤波处理。接下来我们要用设计的FIR低通滤波器对上一节中加噪语音信号进行滤波处理。
图4-8 FIR滤波前和滤波后波形及频谱
分析：从以上四图可以很明显和直观的看出原始语音信号和加噪语音信号时域波形和频谱图的区别。加噪后的语音信号的时域波形比原始语音信号要模糊得多，频谱图则是在频率5000Hz以后出现了明显的变化。
再通过滤波前的信号波形和频谱图的对比，可以明显看出滤波后的波形开始变得清晰了，有点接近原始信号的波形图了。滤波后信号的频谱图也在5000Hz以后开始逐渐接近原始语音信号的频谱图。
图4-9 IIR滤波前和滤波后波形及频谱
4.5 结合去噪后的频谱图对比两种方式滤波的优缺点
IIR数字滤波器采用递归型结构，即结构上带有反馈环路。IIR滤波器运算结构通常由延时、乘以系数和相加等基本运算组成，可以组合成直接型、正准型、级联型、并联型四种结构形式，都具有反馈回路。由于运算中的舍入处理，使误差不断累积，有时会产生微弱的寄生振荡。
（1）IIR数字滤波器的相位特性不好控制，对相位要求较高时，需加相位校准网络。FIR滤波器则要求较低。
（2）IIR滤波器运算误差大，有可能出现极限环振荡，FIR相比之下运算误差较小，不会出现极限环振荡。
（3）IIR幅频特性精度很高，不是线性相位的，可以应用于对相位信息不敏感的音频信号上；
（4）与FIR滤波器的设计不同，IIR滤波器设计时的阶数不是由设计者指定，而是根据设计者输入的各个滤波器参数（截止频率、通带滤纹、阻带衰减等），由软件设计出满足这些参数的最低滤波器阶数。在MATLAB下设计不同类型IIR滤波器均有与之对应的函数用于阶数的选择。
（5）IIR单位响应为无限脉冲序列FIR单位响应为有限的
（6）FIR幅频特性精度较之于iir低，但是线性相位，就是不同频率分量的信号经过FIR滤波器后他们的时间差不变。这是很好的性质。
（7）IIR滤波器有噪声反馈，而且噪声较大，FIR滤波器噪声较小。
FIR幅频特性精度较之于iir低，但是线性相位，就是不同频率分量的信号经过FIR滤波器后他们的时间差不变。这是很好的性质。
参考文献
[1] Boll S F.Suppression of Acoustic Noise in Speech Using Spectral Subtraction[J].IEEE Trans.on Acoustics,Speech,and Signal Processing,1979,27(2):113-120.
[2] Berouti M,Schwartz R,Makhoul J.Enhancement of SpeechCorrupted by Acoustic Noise[J].IEEE Trans.on Acoustics,Speech,and Signal Processing,1979,4:208-211.
[3] 胡航，语音信号处理，哈尔滨工业大学出版社，2000 年 5 月
[4]ThomsonDJ.SpectrumEstimationand armonicAnalysis[J].Proc.IEEE,1982,70(9):1 055-1 096.
[5]皇甫堪,陈建文,楼生强.现代数字信号处理[M].北京:电子工业出版社,2003.
[6]Hu Yi,Loizou P C.Speech Enhancement Based on WaveletThresholding the Multitaper Spectrum[J]. IEEE Trans.onSpeech and Audio Processing,2004,12(1):59-67.
[7]吴红卫,吴镇扬,赵力.基于多窗谱的心理声学语音增强[J].声学学 报,2007,32(3):275-281.
[8]潘欣裕，童兴法，赵鹤鸣，基于谱能比例加权的谱减法语音增强研究[J]. 中 国电子学会第十五届信息论学术年会暨第一届全国网络编码学术年会论文集
[9] 程正,赵鹤鸣. 基于多频带谱减法的语音增强算法的研究[J]. 002-8331(2007)36-0040-03.
[10] 曹瑜镠,方元,吕勇.基于最小统计及谱减法的语音增强[J]. 语音技术. 002-8684(2006)12-0043-04
[11］白文雅，黄健群，陈智伶.基于维纳滤波语音增强算法的改进实现[J].电声 技术，2007,31（1） ：44-46.
[12］蔡斌.一种改进型 MMSE 语音增强方法[J].信号处理，2004,20（1） ：70-74.
[13］陈俊，孙洪，董航.基于 MMSE 先验信噪比估计的语音增强[J].武汉大学学 理学版，2005，51（5） ：638-642.
[14］樊昌信，曹丽娜.通信原理[M].北京：国防工业出版社，2008.
[15］姚天任.数字语音处理[M].武汉：华中科技大学出版社，2005.
2万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
