男女频谱识别

1、如果原序列具有偶对称性，则DFT结果是原时域序列的 N倍。男女频谱识别1. DFT总结（离散傅里叶变换）DFT的定义是针对任意的离散序列 x（nTs）中的有限个离散抽样（°0：N）的，它并不要求 该序列具有周期性。由DFT求出的离散谱H(k) =H k =HNTsk EZ)是离散的周期函数，周期为N 1Nfo =N /T=二 fs、离散间隔为NTs Ts1 fs _丄_怙To 0。离散谱关于变元k的周期为N。NTs N To如果称离散谱经过IDFT所得到的序列为重建信号，x'（nTs）（nZ）,则重建信号是离散的周期1NTs To f函数，周期为0 （对应离散谱的离散间隔的倒

2、数）、离散间隔为Td1Ts =NTs/N 0 :N Nf0 （对应离散谱周期的倒数）。foJ1To NTso经IDFT重建信号的基频就是频域的离散间隔，或时域周期的倒数，为实序列的离散谱关于原点和2 （如果N是偶数）是共轭对称和幅度对称的。因此，真正有用的频谱信息可以从在时域和频域0N0 2 范围获得，从低频到咼频。N范围内的N点分别是各自的主值区间或主值周期。DFT性质线性性：对任意常数amMMDFT 瓦 amXm(n)匕迟 amDFT IXm(n)】,mT_mT奇偶虚实性：DFT的反褶、平移：先把有限长序列周期延拓，再作相应反褶或平移，最后取主值区间的 序列作为最终结果。DFT有如下的奇偶

3、虚实特性：奇二奇；偶二偶；实偶二实偶；实奇二虚奇；实二（实偶）+ j（实奇）；实 （实偶） EXP（实奇）。反褶和共轭性：时域频域反褶反褶共轭共轭+反褶共轭+反褶共轭对偶性：X（n）：=Nx（-k）把离散谱序列当成时域序列进行DFT，结果是原时域序列反褶的N倍;时移性：x(n s) = X(k)WN。序列的时移不影响 DFT离散谱的幅度。频移性：x(n)WNuX(k)时域离散圆卷积定理：x(n) ： y(n)二X(k)丫(k)圆卷积：周期均为N的序列x(n)与y(n)之间的圆卷积为x(n) ： y(n)仍是n的序列，周期为 N。非周期序列之间只可能存在线卷积，不存在圆卷积；周期序列之间存在圆卷

4、积，但不存在线卷积。x(n)y(n)二 1 X(k)®Y(k)频域离散圆卷积定理：N _N J(n) x(i)y*(i n) i=0时域离散圆相关定理：只辭("=X(k)Y (k)R(P)(x(n),y(n)周期为n的序列x(n)和y(n)的圆相关：是n的序列，周期为N。h(n) 1 DFTk H * (k) L*N。其中DFTk1表示按k进行DFT运算。N 4N -1x(n)' X(k)2帕斯瓦尔定理:n=0N k =D2. FFT (离散傅氏变换的快速算法)FFT( Fast Fourier Transformation)，即为快速傅氏变换，是离散傅氏变换的快速

5、算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的FFT算法图(Bufferfly算法) FFT算法图(Bufferfly算法)发现，但是对于在计算机系统或者说数字系统中应用离散傅立叶变换，可以说是进了一大步。设x(n)为N项的复数序列，由 DFT变换，任一 X( m)的计算都需要 N次复数乘法和 N-1 次复数加法，而一次复数乘法等于四次实数乘法和两次实数加法，一次复数加法等于两次实数加法，即使把一次复数乘法和一次复数加法定义成一次“运算”(四次实数乘法和四次实数加法)，那么求出N项复数序列的X(m),即N点DFT变换大约就需

6、要 N2次运算。当 N=1024点甚至更多的时候，需要N2=1048576次运算，在FFT中，禾U用 WN的周期性和对称性，把一个N项序列(设N=2k,k为正整数)，分为两个N/2项的子序列，每个N/2点DFT 变换需要(N/2) 2次运算，再用N次运算把两个 N/2点的DFT变换组合成一个 N点的DFT 变换。这样变换以后，总的运算次数就变成N+2* (N/2)A2=N+ ( NA2) /2。继续上面的例子，N=1024时，总的运算次数就变成了525312次，节省了大约 50%的运算量。而如果我们将这种“一分为二”的思想不断进行下去，直到分成两两一组的DFT运算单元，那么 N点的DFT变换就

7、只需要 Nlog2N次的运算，N在1024点时，运算量仅有10240次，是先前的直接算法的1%，点数越多，运算量的节约就越大，这就是 FFT的优越性。2.源码表示在C环境下的源码源码(1)： 快速傅立叶变换/入口参数：/I: 1=0,傅立叶变换；1=1,逆傅立叶变换/ il: iI=0,不计算傅立叶变换或逆变换模和幅角；il=1,计算模和幅角/ n:输入的点数，为偶数，一般为 32，64, 128，.,1024等/ k:满足n=2Ak(k>0),实质上k是n个采样数据可以分解为偶次幕和奇次幕的次数/ pr: l=0时，存放N点采样数据的实部/ 1=1时，存放傅立叶变换的N个实部/ pi:

8、 l=0时，存放N点采样数据的虚部/ l=1时，存放傅立叶变换的N个虚部/出口参数：/ fr: l=0,返回傅立叶变换的实部/ l=1,返回逆傅立叶变换的实部/ fi: l=0,返回傅立叶变换的虚部/ l=1,返回逆傅立叶变换的虚部/ pr: il=1,i=0 时，返回傅立叶变换的模/ il=1,i=1时，返回逆傅立叶变换的模/ pi: il=1,i=0 时，返回傅立叶变换的辐角/ il=1,i=1时，返回逆傅立叶变换的辐角void fft(double pr, double pi, i nt n, i nt k, double fr, double fi, i nt l, i nt il)

9、int it,m,is,i,j, nv,l0;double p,q,s,vr,vi,poddr,poddi;for(it=0;it<=n-1;m=it+)is=0;for(i=0;i<=k-1;i+)j=m/2;is=2*is+(m-2*j);m=j;frit=pris;fiit=piis;/pr0=1.0;pi0=0.0;p=6.283185306/n;pr1=cos(p);pi1=-si n(p);if (l)pi1=-pi1;for(i=2;i<=n-1;i+)p=pri-1*pr1;q=pii-1*pi1;s=(pri-1+pii-1)*(pr1+pi1);pr=p-

10、q;pi=s-p-q;for(it=0;it<=n-2;it+=2)vr=frit;vi=fiit;frit=vr+frit+1;fiit=vi+fiit+1;frit+1=vr-frit+1;fiit+1=vi-fiit+1;m=n/2;n v=2;for(IO=k-2;IO>=0;IO-)m/=2;nv<<=1;for(it=0;it<=(m-1)* nv;it+=nv)for(j=0;j<=( nv/2)-1;j+)p=prm*j*frit+j+nv/2;q=pim*j*fiit+j+nv/2;s=prm*j+pim*j;s*=(frit+j+nv/2

11、+fiit+j+nv/2);poddr=p_q;poddi=s_p_q;frit+j+nv/2=frit+j-poddr;fiit+j+nv/2=fiit+j-poddi;frit+j+=poddr;fiit+j+=poddi;if(l)for(i=0;i<=n _1;fr/=n, fii+/=n);if(il)for(i=0;i<=n-1;i+)pr=sqrt(fr*fr+fi*fi);if(fabs(fr)<0.000001*fabs(fi)pi=fi*fr>0?90.0-90.0;elsepi=ata n( fi/fr)*360.0/6.283185306;ret

12、urn;源码(2)ps:可以运行的/ The followi ng line must be defi ned before in cludi ng math.h to correctly defi ne M_PI#defi ne _USE_MATH_DEFINES#in clude <math.h>#i nclude <stdio.h>#in clude <stdlib.h>#defi ne PI M_PI /* pi to mach ine precisi on, defi ned in math.h */#defi ne TWOPI (2.0*PI)/

13、*FFT/IFFT routi ne. (see pages 507-508 of Numerical Recipes in C)In puts:data : array of complex* data points of size 2*NFFT+1.data0 is unu sed,* the n 'th complex nu mber x(n), for 0 <= n <= len gth(x)-1, is stored as: data2* n+1 = real(x( n)data2* n+2 = imag(x( n)if len gth(Nx) < NFFT

14、, the rema in der of the array must be padded with zeros nn : FFT order NFFT. This MUST be a power of 2 and >= len gth(x).isig n: if set to 1,computes the forward FFTif set to -1,computes In verse FFT - in this case the output values haveto be man ually no rmalized by multiplyi ng with 1/NFFT.Out

15、puts:data : The FFT or IFFT results are stored in data, overwriti ng the in put.*/void four1(double data, int nn, int isig n)int n, mmax, m, j, istep, i;double wtemp, wr, wpr, wpi, wi, theta;double tempr, tempi;n = nn << 1;j = 1;for (i = 1; i < n; i += 2) if (j > i) tempr = dataj; dataj

16、= datai; datai = tempr;tempr = dataj+1; dataj+1 = datai+1; datai+1 = tempr;m = n >> 1;while (m >= 2 && j > m) j -= m;m >>= 1;j += m;mmax = 2;while (n > mmax) istep = 2*mma x;theta = TWOPI/(isig n*mmax);wtemp = sin( 0.5*theta);wpr = -2.0*wtemp*wtemp;wpi = sin( theta);wr =

17、 1.0;wi = 0.0;for (m = 1; m < mma x; m += 2) for (i = m; i <= n; i += istep) j =i + mma x;tempr = wr*dataj - wi*dataj+1;tempi = wr*dataj+1 + wi*dataj;dataj = datai - tempr;dataj+1 = datai+1 - tempi;datai += tempr;datai+1 += tempi;wr = (wtemp = wr)*wpr - wi*wpi + wr;wi = wi*wpr + wtemp*wpi + wi

18、;mmax = istep;在C+环境下的源码bool FFT(complex<double> * TD, complex<double> * FD, i nt r)一维快速Fourier变换。指向频域数complex<double> * TD指向时域数组的指针 ；complex<double> * FD组的指针；r 2的幕数，即迭代次数LONG count; / Fourier 变换点数int i,j,k; /循环变量int bfsize,p; / 中间变量double angle; / 角度complex<double> *W,*

19、X1,*X2,*X;count = 1 << r; /计算Fourier变换点数为1左移r位W = new complex<double>co unt / 2;X1 = new complex<double>co un t;X2 = new complex<double>cou nt; /分配运算所需存储器/计算加权系数(旋转因子 w的i次幕表)for(i = 0; i < count / 2; i+)an gle = -i * PI * 2 / cou nt;W i = complex<double> (cos(a ngle)

20、, sin(an gle);/将时域点写入X1memcpy(X1, TD, sizeof(complex<double>) * coun t);/采用蝶形算法进行快速Fourier变换for(k = 0; k < r; k+)for(j = 0; j < 1 << k; j+)bfsize = 1 << (r-k);for(i = 0; i < bfsize / 2; i+)p = j * bfsize;X2i + p = X1i + p + X1i + p + bfsize / 2 * Wi * (1<<k);X2i + p

21、+ bfsize / 2 = X1i + p - X1i + p + bfsize / 2 * Wi * (1<<k);X = X1;XI = X2;X2 = X;/重新排序for(j = 0; j < coun t; j+)p = 0;for(i = 0; i < r; i+)if (j&(1<<i)p+=1<<(r-i_1)；FDj=X1p;/释放内存delete W;delete X1;delete X2;return true;在Matlab环境下的源码fun cti on X=myfft(x)%myfft函数用递归实现N=le

22、ngth(x);t=log2(N);t1=floor(t);t2=ceil(t);if t1=t2; %若x的长度N不为2的整数次幕，则补 0至最接近的2的整数次幕x=x zeros(1,2At2-N);N=2At2;endw0=exp(-j*2*pi/N);X=zeros(1,N);if N=2X(1)=x(1)+x(2);X(2)=x(1)-x(2);elsen=1:N/2;xe( n)=x(2* n-1);xo( n)=x(2* n);XE=myfft(xe); % 递归调用XO=myfft(xo);for n=1:N/2X( n)=XE( n)+XO( n)*(wO( n-1);X(n

23、+N/2)=XE( n)-XO( n)*(wO( n-1);endEnd4.MATLAB中FFT的使用方法一 调用方法X=FFT(x);X=FFT(x , N);x=IFFT(X);x=IFFT(X,N)用MATLAB进行谱分析时注意：(1) 函数FFT返回值的数据结构具有对称性。例：N=8;n=0:N-1;xn=4 3 2 6 7 8 9 0;Xk=fft(x n)Xk =39.0000 -10.7782 + 6.2929i 0 - 5.0000i 4.7782 - 7.7071i 5.0000 4.7782 + 7.7071i 0 + 5.0000i -10.7782 - 6.2929iX

24、k与xn的维数相同，共有8个元素。Xk的第一个数对应于直流分量，即频率值为0。(2) 做FFT分析时，幅值大小与 FFT选择的点数有关，但不影响分析结果。在IFFT时已经做了处理。要得到真实的振幅值的大小，只要将得到的变换后结果乘以2除以N即可。二.FFT应用举例例 1: x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t)。采样频率 fs=100Hz，分别绘制 N=128、1024点幅频图。elf;fs=100;N=128; %采样频率和数据点数n=0:N-1;t =n/fs; % 时间序列x=0.5*si n(2*pi*15*t)+2*si n(2*pi*40*t);

25、 %信号y=fft(x,N); %对信号进行快速 Fourier变换 mag=abs(y); %求得Fourier变换后的振幅 f=n*fs/N; %频率序列subplot(2,2,1),plot(f,mag); %绘出随频率变化的振幅xlabel('频率 /Hz');ylabel('振幅');title('N=128');grid on;subplot(2,2,2),plot(f(1:N/2),mag(1:N/2); % 绘出 Nyquist 频率之前随频率变化的振幅 xlabel('频率 /Hz');ylabel('振

26、幅');title('N=128');grid on;%对信号采样数据为1024点的处理fs=100;N=1024; n=0:N-1;t=n/fs;x=0.5*si n(2*pi*15*t)+2*si n(2*pi*40*t); %信号y=fft(x,N); %对信号进行快速 Fourier变换 mag=abs(y); %求取Fourier变换的振幅 f=n *fs/N;subplot(2,2,3),plot(f,mag); %绘出随频率变化的振幅xlabel('频率 /Hz');ylabel('振幅');title('N=102

27、4');grid on;subplot(2,2,4)plot(f(1:N/2),mag(1:N/2); % 绘出Nyquist频率之前随频率变化的振幅 xlabel('频率 /Hz');ylabel('振幅');title('N=1024');grid on;运行结果：fs=100Hz , Nyquist频率为fs/2=50Hz。整个频谱图是以 Nyquist频率为对称轴的。并且可以 明显识别出信号中含有两种频率成分：15Hz和40Hz。由此可以知道FFT变换数据的对称性。因此用FFT对信号做谱分析，只需考察0Nyquist频率范围内的福

28、频特性。若没有给出采样频率和采样间隔，则分析通常对归一化频率01进行。另外，振幅的大小与所用采样点数有关，采用128点和1024点的相同频率的振幅是有不同的表现值，但在同一幅图中，40Hz与15Hz振动幅值之比均为 4： 1，与真实振幅0.5： 2是一致的。为了与真实振幅对应，需要将 变换后结果乘以2除以N。例 3： x=cos(2*pi*0.24*n)+cos(2*pi*0.26*n)(1) 数据点过少，几乎无法看出有关信号频谱的详细信息；(2) 中间的图是将x(n)补90个零，幅度频谱的数据相当密，称为高密度频谱图。但从 图中很难看出信号的频谱成分。(3) 信号的有效数据很长，可以清楚地看

29、出信号的频率成分，一个是0.24Hz，一个是 0.26Hz，称为高分辨率频谱。可见，采样数据过少，运用FFT变换不能分辨出其中的频率成分。添加零后可增加频谱中的数据个数，谱的密度增高了，但仍不能分辨其中的频率成分，即谱的分辨率没有提高。只有数据点数足够多时才能分辨其中的频率成分。3. 男女性语言频谱我们知道男性基音频谱在 60200HZ之间，女性的基音频谱在200450HZ之间。wenku.baiduxonriHJI-T -il研11暹田立iJIBL百JT亠貝中：生音裁暹诵8L百“.忍干H号才氐誓卅上手99科41的囲对声昔信号的处理分桁L " J1”频域上的分析在这之后，我们对信号进

30、行快速傅里叶变换.得到相应的 频谱分析图。就此，我们可以得到相应的(男女芦)频谱图对比如下：4.隐马尔可夫模型(Hidden Markov Model , HMM )是统计模型，它用来描述一个含有隐含未知参数的马尔可夫过程。其难点是从可观察的参数中确定该过程的隐含参数。然后利用这些参数来作进一步的分析，例如模式识别。基本理论隐马尔可夫模型 是马尔可夫链 的一种，它的状态不能直接观察到，但能通过观测向量 序列观察到，每个观测向量都是通过某些 概率密度分布表现为各种状态， 每一个观测向量是 由一个具有相应概率密度分布的状态序列产生。所以，隐马尔可夫模型是一个双重随机过程-具有一定状态数的隐 马尔可

31、夫链和显示随机函数集。自 20世纪80年代以来，HMM被 应用于语音识别，取得重大成功。到了 90年代，HMM还被引入计算机文字识别和移动通 信核心技术 多用户的检测” HMM在生物信息科学、故障诊断等领域也开始得到应用。基本算法针对以下三个问题，人们提出了相应的算法1评估问题前向算法2解码问题Viterbi 算法3学习问题Baum-Welch 算法(向前向后算法)模型表达隐马尔可夫模型(HMM )可以用五个元素来描述，包括2个状态集合和3个概率矩阵：1. 隐含状态(S)这些状态之间满足 马尔可夫性质，是马尔可夫模型 中实际所隐含的状态。这些状态通常 无法通过直接观测而得到。(例如 S1、S2

32、、S3等等)2. 可观测状态(0)在模型中与隐含状态相关联，可通过直接观测而得到。(例如O1、O2、O3等等，可观测状态的数目不一定要和隐含状态的数目一致。3. 初始状态概率矩阵(n)表示隐含状态在初始时刻t=1的概率矩阵，侧如t=1时，P(S1)=p1、P(S2)=P2、P(S3)=p3 ,则初始状态概率矩阵n = p1 p2 p3 .4. 隐含状态转移概率矩阵(A)描述了 HMM模型中各个状态之间的转移概率。其中Aij = P( Sj | Si ),1<表示崔Nt.时刻、状态为 Si的条件下，在 t+1时刻状态是 Sj的概率。5.观测状态转移概率矩阵B (英文名为Confusion

33、Matrix ，直译为混淆矩阵不太易于从字面理解)。令N代表隐含状态数目，M代表可观测状态数目，则：Bij = P( 0i | Sj ), 1 w i < M,1<j < N.表示在t时刻、隐含状态是 Sj条件下，观察状态为Oi的概率。总结：一般的，可以用 入=(A,B, t三元组来简洁的表示一个隐 马尔可夫模型。隐马尔可夫模型实际上是标准马尔可夫模型的扩展，添加了可观测状态集合和这些状态与隐含状态 之间的概率关系。基本问题1. 评估问题。给定观测序列 O=O1O2O3Ot和模型参数 入=(A,B, n)怎样有效计算某一观测序列的概率，进而可对该 HMM做出相关评估。例如，已

34、有一些模型参数各异的HMM，给定观测序列O=O1O2O3Ot，我们想知道哪个 HMM模型最可能生成该观测序列。通常我们利用forward算法分别计算每个 HMM产生给定观测序列 O的概率，然后从中选出最优的 HMM模型。这类评估的问题的一个经典例子是 语音识别。在描述语言识别的隐马尔科夫模型中， 每个单词生成一个对应的 HMM，每个观测序列由一个单词的语音构成，单词的识别是通过 评估进而选出最有可能产生观测序列所代表的读音的HMM而实现的。2. 解码问题给定观测序列 O=O1O2O3Ot 和模型参数入=(A,B, n,)怎样寻找某种意义上最优的隐 状态序列。在这类问题中，我们感兴趣的是马尔科夫

35、模型中隐含状态， 这些状态不能直接观 测但却更具有价值，通常利用 Viterbi算法来寻找。这类问题的一个实际例子是中文分词，即把一个句子如何划分其构成才合适。例如，句子 发展中国家”是划分成 发展-中-国家”还是 发展-中国-家”这个问题可以用隐马尔科夫 模型来解决。句子的分词方法可以看成是隐含状态，而句子则可以看成是给定的可观测状态，从而通过建HMM来寻找出最可能正确的分词方法。5. MFCC(梅尔倒频谱系数)在语音辨识(Speech Recognition)和语者辨识(Speaker Recognition)方面，最常用到的 语音特征就是梅尔倒频谱系数(Mel-scale Frequen

36、cy Cepstral Coefficients，简称 MFCC ),此参数考虑到人耳对不同频率的感受程度，因此特别适合用在语音辨识。下面简单的介绍一下求解MFCC的过程。(1) .预强调(Pre-emphasis ):将语音讯号 s(n)通过一个高通滤波器。 H(z)=1-a* (z-1)系数其中a介于0.9和1.0之间。若以时域的表达式来表示，预强调后的 讯号s2(n)为s2(n) = s(n) - a*s(n-1)。这个目的就是为了消除发声过程中声带和嘴唇的效应，来补偿语音信号受到发音系统所压抑的高频部分。(另一种说法则是要突显在高频的共振峰。)(2) 音框化(Frame blockin

37、g ):先将N个取样点集合成一个观测单位，称为音框(Frame),通常N的值是256或512,涵盖的时间约为2030 ms左右。为了避免相邻两音框的变化过大，所以我们会让两相邻因框之间有一段重迭区域，此重迭区域包含了M个取样点，通常M的值约是 N的一半或1/3。通常语音辨识所用的音讯的取样频率为8 KHz或16KHz，以8 KHz来说，若音框长度为256个取样点，则对应的时间长度是256/8000*1000=32 ms。(3) .汉明窗(Hamming window ):将每一个音框(frame)乘上汉明窗，以增加音框左端和右端的连续性(请见下一个步骤的说明)。假设音框化的讯号为S(n), n

38、 = 0,N-1。N为frame的大小，那么乘上汉明窗后为S'(n) = S(n)*W(n)，此 W(n)形式如下：W(n, a) = (1 - a) - a*cos(2pn/(N-1) , 0W nW N-1 ? ? 不同的a值会产生不同的汉明窗。一般我们都取a =0.46。(4) .快速傅里叶转换(Fast Fourier Transform, or FFT ):由于讯号在时域(Time domain )上的变化通常很难看出讯号的特性，所以通常将它转换成频域(Frequency domain )上的能量分布来观察，不同的能量分布，就能代表不同语音的特性。所以在乘上汉明窗后，每个音框

39、还必需再经过 FFT以得到在频谱上的能量分布.乘上汉明窗的主要目的，是要加强音框左端和右端的连续性，这是因为在进行FFT时，都是假设一个音框内的讯号是代表一个周期性讯号，如果这个周期性不存在， FFT会为了要符合左右端不连续的变化，而产生一些不存 在原讯号的能量分布，造成分析上的误差。当然，如果我们在取音框时，能够使音框中的讯 号就已经包含基本周期的整数倍，这时候的音框左右端就会是连续的，那就可以不需要乘上汉明窗了。但是在实作上，由于基本周期的计算会需要额外的时间，而且也容易算错，因此 我们都用汉明窗来达到类似的效果。(5) .三角带通滤波器(Triangular Bandpass Filte

40、rs ):将能量频谱能量乘以一组20个三角带通滤波器，求得每一个滤波器输出的对数能量(Log Energy ),共20个。必须注意的是：这20个三角带通滤波器在梅尔频率 (Mel Frequency)上是平均分布的，而梅尔频率和一 般频率 f 的关系如下：mel(f)=2595*log10(1+f/700)或是mel(f)=1125*l n( 1+f/700)梅尔频率代表一般人耳对于频率的感受度，由此也可以看出人耳对于频率f的感受是呈对数变化的：在低频部分，人耳感受是比较敏锐。在高频部分，人耳的感受就会越来越粗糙。三角带通滤波器有两个主要目的：对频谱进行平滑化， 并消除谐波的作用，突显原先语音

41、的共振峰。(因此一段语音的音调或音高，是不会呈现在MFCC参数内，换句话说，以MFCC为特征的语音辨识系统，并不会受到输入语音的音调不同而有所影响。) 降低资料量。(6) .离散余弦转换(Discrete cosine transform, or DCT):将上述的 20个对数能量 Ek 带入离散余弦转换， 求出L阶的Mel- scale Cepstrum参数，这里L通常取12。离散余弦 转换公式如下：Cm=Sk=1Ncosm*(k-0.5)*p/N*Ek, m=1,2, ., L 其中 Ek 是由前一个步骤所算出来的三角滤波器和频谱能量的内积值，这里N是三角滤波器的个数。由于之前作了 FFT

42、，所以采用 DCT转换是期望能转回类似Time Domain的情况来看，又称 Quefrency Domain，其实也就是 Cepstrum (倒谱)。又因为之前采用 Mel- Frequency 来 转换至梅尔频率，所以才称之Mel-scale Cepstrum。(7) .对数能量(Log energy): 一个音框的音量(即能量)，也是语音的重要特征，而且非 常容易计算。因此我们通常再加上一个音框的对数能量(定义为一个音框内讯号的平方和，再取以10为底的对数值，再乘以10)，使得每一个音框基本的语音特征就有13维，包含了 1个对数能量和 12个倒频谱参数。(若要加入其它语音特征以测试辨识率

43、， 也可以在此 阶段加入，这些常用的其它语音特征，包含音高、过零率、共振峰等。)(8) .差量倒频谱参数(Delta cepstrum):虽然已经求出 13个特征参数，然而在实际应用 于语音辨识时，我们通常会再加上差量倒频谱参数，以显示倒频谱参数对时间的变化。6. MFCC的特征提取用人类生物特征(声纹、指纹、人脸、虹膜、DNA等)并结合计算机技术进行安全验证是当今的热门课题。声纹识别技术是生物识别技术的一种，与其他生物识别技术相比，声纹识别具有更为简便、准确、经济及非接触性识别等众多优势。声纹识别是一项根据声纹中反映说话人生理和行为特征的语音参数，自动识别说话人身份的技术。声纹识别过程中最主

44、要 的两部分内容是特征提取和模式匹配。特征提取，就是从声音中选取唯一表现说话人身份的有效且稳定可靠的特征；模式匹配就是对训练和鉴别时的特征模式做相似性匹配。基于人耳听觉特性的语音特征 Mel频率倒谱系数(MFCC),可以作为声纹识别的重要特征参数之一。 MFCC是根据人的听觉感知机理提出的能反映人的声音特征的参数，因此它没有任何前提假 设，且具有良好的识别能力和抗噪声能力。由于人耳所听到的声音的高低与声音的频率并不成正比关系，而用 Mel频率尺寸则更符合人耳的听觉特性。 Mel频率与实际频率的具体关 系可用公式表示同:Mel ( f)=2 595 lg( 1+f /700) (1)下图是MFC

45、C参数的提取框图：图1 MFCC的提取过程MFCC参数提取过程如下1:(1)原始语音信号在经过采样量化后需经过预加重处理以提高信号的高频部分的幅度，计算机里用具有6dB/倍频程的提升高频特性的预加重数字滤波器来实现，它一般是一阶的数字滤波器：JH =1z(2)研究发现，当 卩=0.950.97时，说话人识别系统的识别率最高，因此本文采用预加重系 数为0.97。(2)根据语音信号短时平稳的特性,通过分帧操作提取语音短时特性便于建模，一般取帧长30 ms ,此时帧移10 ms;为了平滑信号以减少每帧信号两端的预检测误差，避免频谱出现“破碎”现象采用加窗处理。每帧用窗函数相乘以减少帧起始和结束处的信号不连续性。 在语音信号的时域分析来说，窗函数的形状很重要，虽然矩形窗平滑比较好，但容易使波形细节丢失，并且会产生泄露现象，而Hamming窗则能有效地克服泄露现象，所以，在声纹 识别中加窗的方法大致以Hamming窗最为常见，本系统也采用 Hamming窗(系统采用N=410点):工2二nQ.54 0.46 cos(),0 兰 nEN -1W(n) = $N -10,其他(3) 端点检测是语音识别中的一个重要环节，有效的端点检测方法不仅可以减少数据的存储量和处理时