倒谱计算与分析

1、视频语音处理技术倒谱计算与分析学院名称： 计算机与信息工程学院专业名称： 计算机科学与技术年级班级：姓 名：学 号：计算机与信息技术学院综合性、设计性实验报告专业：计算机科学技术年级/班级：2011级 2012 2013学年第一学期课程名称视频语音处理技术指导教师张新明本组成员 学号姓名1108114153 王彦秋实验地点计科楼324实验时间项目名称倒谱计算与分析实验类型设计性一、实验目的：对语音信号进行同态分析可得到语音信号的倒谱参数。语音的倒谱是将语音的短 时谱取对数后再进行IDFT得到的，所以浊音信号的激励反映在倒谱上是同样周期的冲 激，借此，可从倒谱波形中估计出基音周期。对倒谱进行低时

2、窗选，通过语音倒谱分析 的最后一级，进行DFT后的输出即为平滑后的对数模函数，这个平滑的对数谱显示了特 定输入语音段的谐振结构，即谱的峰值基本上对应于共振峰频率，对于平滑过的对数谱 中的峰值进行定位，即可估计共振峰。对于倒谱计算与分析的设计实验可作如下训练：1、复倒谱的几种计算方法：2、最小相位信号法和递归法；3、基音检测；4、共振峰检测。二、实验仪器或设备：windowsXP下的Matlab编程环境三、总体设计(设计原理、设计方案及流程等)1. 复倒谱的几种计算方法：在复倒谱分析中，z变换后得到的是复数，所以取对数时要进行复对数运算。这时 存在相位的多值性问题，称为“相位卷绕”。设信号为x(

3、n)为(n) X2(n)则其傅里叶变换为X(ej ) X,ej ) X2(ej )对上式取复对数为ln X(ej ) ln Xj(ej ) ln X2(ej )则其幅度和相位分别为：ln X(ej ) ln Xj(ej ) ln X2 )()1( )2()()1( )2()上式中，虽然 1( ) ,2()的范围均在 ， 内，但 ()的值可能超过,范围。计算机处理时总相位值只能用其主值()表示，然后把这个相位主值“展开”，得到连续相位。所以存在下面的情况：()()2k (K为整数)此时即产生了相位卷绕。下面介绍几种避免相位卷绕求复倒谱的方法。最小相位信号法这是解决相位卷绕的一种较好的方法。但它有

4、一个限制条件：被处理的信号想x(n)必须是最小相位信号。实际上许多信号就是最小相位信号，或可以看作是最小相位信号。语音信号的模型就是极点都在 z平面单位圆内的全极点模型，或者极零点都在z平面单位圆内的极零点模型。设信号x(n)的z变换为X(z)=N (z)/D(z)，则有ln X(z)lnN(z)D(z)根据z变换的微分特性有nX(n)z爭z lndzN(z)D(z)D(z)N (z) N(z)D (z) zN(z) D(z)若x( n)是最小相位信号，则X(必然是稳定的因果序列由Hilbert变换的性质可知,任一因果复倒谱序列都可分解为偶对称分量和奇对称分量之和：X(n) Xe(n) X)(

5、n)其中Xe(n)X(n) X( n) /2X(n)X(n) ? n) /2这两个分量的傅里叶变换分别为？(的傅里叶变换的实部和虚部Xej )x(n)e jn:n*R(ej ) j)?I(ej )从而可得0n 0x(n)鬼(n)n 02xe( n)n 0此即复倒谱的性质3,也就是说一个因果序列可由其偶对称分量来恢复。如果引入一个X(n) g(n) Xe(n)辅助因子g(n)，上式可写作其中0n0g( n)1n02n0最小相位信号法求复倒谱原理框图如下递归法这种方法仅限于是最小相位信号的情况。根据z变换的微分特性得zX(z)：Xz?X(z)dzdz对上式求逆z变换，根据z变换的微分特性，有n ?

6、(n) x(n) n x(n)所以x(n)kk?(k)x(nnk)n 0设x( n)是最小相位序列，而最小相位信号序列一定为因果序列，所以有n kn 1k小cx(n)( )X(k)x(n k)5?(k)x( n k) 5?(n) x(0)k 0 nk 0n由于x(k) 0(k0)及5?(n k) 0( kn)可得递推公式x(n)x(n) n1 k x(k)x(n k)x(0) k 0 nx(0)递归运算后由复倒谱定义x(n)z 1 l n z x(n)z 1 lnx(n)z nn1可知x(0) z In z x(0) In x(0) (n)ln x(0)如果x(n)是最大相位序列，则变为0n

7、0g(n)1n 02n 0x(n) x(：)0 k()?(k)x(n k) n 0x(0)k n 1 nx(0)其中5?(0) In x(0)2、基音检测；语音的倒谱是将语音的短时谱取对数后再进行IDFT得到的，所以浊音信号的周期性激励反映在倒谱上是同样周期的冲激。借此，可从倒谱波形中估计出基音周期。一般 把倒谱波形中第二个冲激，认为是对应激励源的基频。下面给出一种倒谱法求基音周期 的框图及流程图如下w( n)s(n)HDFT冷丽刑In |S(k)卜|lDFT卜W值检测基音检测图5.8一种倒谱法求基音周期的实现框图图5.9 一种倒谱法求基音周期的流程图3. 共振峰检测对倒谱进行低时窗选，通过语

8、音倒谱倒谱将基音谐波和声道的频谱包络分离开来分析系统的最后一级，进行 DFT后的输出即为平滑后的对数模函数，这个平滑的对数 谱显示了特定输入语音段的谐振结构， 即谱的峰值基本上对应于共振蜂频率，对平滑过 的对数谱中的峰值进行定位，即可估计共振峰。原理框图及流程图如下。|fft|*log*IFFT加窗*|fft|分帧加窗log平滑处理一共振峰图5.12共振峰检测框图图5.13共振峰检测流程图s,fs ,n bit=wavread(beiji ng.wav); % b=s;%x=b(5000:5399);%N=le ngth(x);%S=fft(x);%Sa=log(abs(S);%logsa=i

9、fft(Sa);%对Sa进行傅里叶逆变换四、实验步骤（包括主要步骤、代码分析等）1.倒谱MATLAB实现代码段clear all; %倒谱读入一段语音将s转置取400点语音读入语音的长度对x进行傅里叶变换为以e为底的对数yle n=len gth(sa);for i=1:yle n/2;sal(i)=sa(yle n/2+1-i);endfor i=(yle n/2+1):yle n; sal(i)=sa(i+1-yle n/2);end淤图figure(1);subplot(2,1,1);plot(x);%axis(0,400,-0.5,0.5) title(截取的语音段);xlabel(样

10、点数);ylabel(幅度);subplot(2,1,2);time2=-199:1:-1,0:1:200; plot(time2,sa1);%axis(-200,200,-0.5,0.5) title(截取语音的倒谱);xlabel(样点数);ylabel(幅度);运行结果如图：皆5050r帧长起始序号结束序号取出该帧2.倒谱法求浊音、清音的基音周期fun cti on s=p5_2pitchdetect waveFile=beiji ng.wav;y,fs, nbits=wavread(waveFile);time1=1:le ngth(y);time=(1:le ngth(y)/fs;

11、frameSize=floor(50*fs/1000);%startl ndex=rou nd(5000);%en di ndex=startl ndex+frameSize-1;% frame=y(start In dex:e ndln dex); frameSize=le ngth(frame);frame2=frame.*ham min g(le ngth(frame); rwy=rceps(frame2);%yle n=len gth(rwy);cepstrum=rwy(1:yle n/2);for i=1:yle n/2;cepstrum1(i)=rwy(yle n/2+1-i);e

12、ndfor i=(yle n/2+1):yle n;cepstrum1(i)=rwy(i+1-yle n/2);end%基因检测LF=floor(fs/500);%HF=floor(fs/70);cn=cepstrum(LF:HF);mx_cep in d=max(c n);if mx_cep0.08&in dLF;a=fs/(LF+i nd);elsea=0;end加汉明窗求倒谱基因周期的范围是70500Hzpitch=a%画图figure(1);subplot(3,1,1);plot(time1,y);title(语音波形);%axis tightylim=get(gca,ylim);li

13、n e(time1(startI ndex),time1(startl ndex),ylim,color,r); lin e(time1(e ndl ndex),time1(e ndl ndex),ylim,color,r); xlabel(样点数);ylabel(幅度);subplot(3,1,2); plot(frame);%axis(0,400,-0.5,0.5) title(帧语音);xlabel(样点数); ylabel(幅度) subplot(3,1,3);time2=-199:1:-1,0:1:200;plot(time2,cepstrum1);%axis(-200,200,-0

14、.5,0.5) ti tle(帧语音的倒谱);xlabel（样点数）;ylabel（幅度）;运行结果如下图：倒谱法求浊音的基音周期OO05TD.511.522.533.5样点数恢语音-帧语音的倒i普2D0 -Iffl -100石 0050100150200拝点数C.4r 2 工2 -O.dat 11.B2? t 33巧禅翩XID*一帧诲音匚 5 iaiiii小 J0 10019Q2003D03504C0样点数嵋蒔音的闻清fenJU.5Q 5-?CC J 和 -10D-5Dn样点皱inn !和 知口清音的倒谱3.共振峰检测程序waveFile=qi nghua.wav;y,fs, nbits=w

15、avread(waveFile);time=(1:le ngth(y)/fs;frameSize=floor(40*fs/1000);%帧长startl ndex=rou nd(15000);%起始序号结束序号取出该帧加汉明窗倒谱求en di ndex=startl ndex+frameSize-1;%frame=y(start In dex:e ndln dex);%frameSize=le ngth(frame);frame2=frame.*ham min g(le ngth(frame); % rwy=rceps(frame2);% yle n=len gth(rwy);cepstrum

16、=rwy(1:yle n/2);%基音检测LF=floor(fs/500);HF=floor(fs/70);cn=cepstrum(LF:HF);mx_cepi nd=max(c n);%找到最大的突起的位置%共振峰检测核心代码NN=i nd+LF;ham=hammi ng(NN);cep=cepstrum(1:NN);ceps=cep.*ham;%汉明窗forma nt1=20*log(abs(fft(ceps);forma nt(1:2)=forma nt1(1:2);for t=3:NN% do some median filteringz=forma nt1(t-2:t);md=med

17、ia n(z);forma nt2(t)=md;endfor t=1:NN-1if t=2forma nt(t)=forma nt1(t);elseforma nt(t)=forma nt2(t-1)*0.25+forma nt2(t)*0.5+forma nt2(t+1)*0.25;endendsubplot(3,1,1);plot(cepstrum);title( 倒谱);xlabel(样点数);ylabel(幅度)%axis(0,220,-0.5,0.5)spectral=20*log(abs(fft(frame2);subplot(3,1,2);xj=(1:le ngth(spectral)/2)*fs/le ngth(spectral);plot(xj,spectral(1:le ngth(spectral)/2);title(频谱);xlabel(频率 /Hz);ylabel(幅度 /dB) %axis(0,500,-100,50) subplot(3,1,3);xi=(1:NN /2)*fs/NN; plot(xi,forma nt(1:NN/2); title( 平滑对数幅度谱); xlabel(频率 /Hz); ylabel(幅度 /dB);%axis(0,5500,-80,0)运行结果如图所示：倒谱-0.5r1II3Lr jJ|2 -J1 (rrrrJ