搜索
华南理工大学《语音信号处理》实验报告实验名称：DTW算法实现及语音模板匹配姓名：学号：班级：11级电信7班日期：2014年5月一、实验目的利用DTW(DynamicTimeWarping，动态时间规整)算法，进行说话者的语音识别。二、实验原理1、语音识别系统概述一个完整特定人语音识别系统的方案框图如图1所示。输入的模拟语音信号首先要进行预处理，包括预滤波、采样和量化、加窗、端点检测、预加重等,然后是参数特征量的提取。提取的特征参数满足如下要求：(1)特征参数能有效地代表语音特征,具有很好的区分性；(2)参数间有良好的独立性；(3)特征参数要计算方便,要考虑到语音识别的实时实现。图1语音识别系统方案框图语音识别的过程可以被看作模式匹配的过程，模式匹配是指根据一定的准则，使未知模式与模型库中的某一个模型获得最佳匹配的过程。模式匹配中需要用到的参考模板通过模板训练获得。在训练阶段，将特征参数进行一定的处理后，为每个词条建立一个模型，保存为模板库。在识别阶段，语音信号经过相同的通道得到语音特征参数，生成测试模板，与参考模板进行匹配，将匹配分数最高的参考模板作为识别结果。2、语音信号的处理1、语音识别的DTW算法本设计中，采用DTW算法，该算法基于动态规划(DP)的思想解决了发音长短不一的模板匹配问题，在训练和建立模板以及识别阶段，都先采用端点检测算法确定语音的起点和终点。在本设计当中，我们建立的参考模板，m为训练语音帧的时序标号，M为该模板所包含的语音帧总数，R(m)为第m帧的语音特征矢量。所要识别的输入词条语音称为测试模板，n为测试语音帧的时序标号，N为该模板所包含的语音帧总数，T(n)为第n帧的语音特征矢量。参考模板和测试模板一般都采用相同类型的特征矢量(如LPCC系数)、相同的帧长、相同的窗函数和相同的帧移。通常，规整函数被限制在一个平行四边形的网格内，如图2所示。它的一条边斜率为2，另一条边斜率为1/2。规整函数的起点是（1,1），终点为（N，M）。DTW算法的目的是在此平行四边形内由起点到终点寻找一个规整函数，使其具有最小的代价函数，保证了测试模板与参考模板之间具有最大的声学相似特性。图2匹配路径约束示意图由于在模板匹配过程中限定了弯折的斜率，因此平行四边形之外的格点对应的帧匹配距离是不需要计算的。另外，因为每一列各格点上的匹配计算只用到了前一列的3个网格，所以没有必要保存所有的帧匹配距离矩阵和累积距离矩阵。充分利用这两个特点可以减少计算量和存储空间的需求，形成一种高效的DTW算法。图2中，把实际的动态弯折分为三段，（1，xa），（xa+1，xb），（xb+1，N），其中：xa=(2M-N)/3，xb=2(2N-M)/3xa和xb都取最相近的整数，由此可得出对M和N长度的限制条件：2M-N≥3，2N-M≥2当不满足以上条件时，认为两者差别太大，则无法进行动态弯折匹配。在x轴上的每一帧不再需要与y轴上的每一帧进行比较，而只是与y轴上[ymin，ymax]间的帧进行比较，ymin和ymax的计算公式为：ymin=x/2，0≤x≤xb，2x+(M-2N)，xbx≤Nymax=2x，0≤x≤xa，x/2+(M-N/2)，xax≤N如果出现xaxb的情况，则弯折匹配的三段为（1，xb），（xb+1，xa），（xa+1，N）。对于x轴上每前进一帧，虽然所要比较的y轴上的帧数不同，但弯折特性是一样的，累积距离的更新都是用下式实现的：D(x，y)=d(x，y)+min[D(x-1，y)，D(x-1，y-1)，D(x-1，y-2)]号的短时能量或短时平均幅度就能够把语音段和噪声背景区分开。这是仅基于短时能量的端点检测方法。信号{x(n)}的短时能量定义为:语音信号的短时平均幅度定义为:其中w(n)为窗函数。2、短时平均过零率短时过零表示一帧语音信号波形穿过横轴(零电平)的次数。过零分析是语音时域分析中最简单的一种。信号{x(n)}的短时平均过零率定义为:式中，sgn为符号函数，即:过零率有两类重要的应用:第一，用于粗略地描述信号的频谱特性;第二，用于判别清音和浊音、有话和无话解决这个问题的办法，一个是做高通滤波器或带通滤波，减小随机噪声的影响；另一个有效方法是对上述定义做一点修改，设一个门限T，将过零率的含义修改为跨过正负门限。于是，有定义:3、检测方法此时整个端点检测可分为四段：静音段、过渡段、语音段、结束。实验时使用一个变量表示当前状态。静音段，如果能量或过零率超过低门限，就开始标记起始点，进入过渡段。过渡段当两个参数值都回落到低门限以下，就将当前状态恢复到静音状态。而如果过渡段中两个参数中的任一个超过高门限，即被认为进入语音段。处于语音段时，如果两参数降低到门限以下，而且总的计时长度小于最短时间门限，则认为是一段噪音，继续扫描以后的语音数据，否则标一记结束端点。三、实验数据及平台本实验所采用的数据是录制的语音文件，平台是MATLAB。四、实验过程（步骤）1、语音信号预处理语音信号的预处理包括预滤波、采样和量化、加窗、预加重、端点检测等过程。由于语音信号在帧长为10ms~30ms之内是相对平稳的，同时为了便于计算FFT，本系统选取帧长N为256个语音点，帧移M为128点。本文采用汉明窗对语音信号进行分帧处理，如下式：ω(n)=0.54-0.46cos(2πn/(N-1))，0≤n≤N-1预加重用具有6dB/倍频程的提升高频特性的一阶数字滤波器实现：H(z)=1-0.9375/z端点检测采用基于短时能量和短时平均过零率法，利用已知为“静态”的最初十帧信号为短时能量设置2个门限ampl和amph，以及过零率阈值zcr。语音起始点从第11帧开始检测，其流程图如图4。语音结束点的检测方法与检测起点相似，但此时从后向前搜索。图3语音起点检测流程图2、特征参数提取及语音识别本文选取能够反映人对语音的感知特性的Mel频率倒谱系数(MFCC)作为特征参数，阶数为12。经过MFCC特征参数提取后，各帧语音信号就形成了一个个特征矢量。识别时，将待测语音与模板库中的每一个模板进行模式匹配，找到距离最小的模板作为输出结果。五、实验结果及讨论如图4为语音信号“8”的处理结果，其他语音信号处理结果图不在此一一给出：图4语音信号“8”的处理结果图以下为得到的最终处理结果，10个数字识别正确。经测试，程序等到了较好的语音识别效果。六、实验代码主程序：yuyinshibie.mdisp('正在计算参考模板的参数...')fori=1:10fname=sprintf('%da.wav',i-1);x=wavread(fname);[x1x2]=vad(x);m=mfcc(x);m=m(x1-2:x2-4,:);ref(i).mfcc=m;enddisp('正在分析语音信号...')fori=1:10fname=sprintf('%da.wav',i-1);[x,fs,bit]=wavread(fname,[2000,2512]);%采样%%sound(x,fs);%播放语音信号figure(i);subplot(3,3,1);plot(x(1:256));%原始语音信号的时域图形%title('原始信号')subplot(3,3,2)[h,w]=freqz(x,fs)%原始语音信号的频率响应图hr=abs(h);%求系统幅频响应plot(w,hr);title('幅频图');xlabel('Frequencyinradmple')ylabel('MagnitudeindB')subplot(3,3,3)hphase=angle(h);hphase=unwrap(hphase);%求系统相频响应plot(w,hphase);title('相频图');xlabel('Frequencyinradmple')ylabel('Phaseindegrees')y=fft(x,512);%傅立叶变换%mag=abs(y);mag1=10*log10(mag);f=fs*(0:255)/512;subplot(3,3,4)plot(f,mag(1:256));%FFT频谱图%title('fft变换后信号')iff=ifft(y,512);%反傅立叶变换%ifm=abs(iff);subplot(3,3,5)plot(f,ifm(1:256))title('ifft后信号')%短时傅里叶变换Ts=1/fs;%N=T/Ts;N=512;Nw=20;%窗函数长L=Nw/2;%窗函数每次移动的样点数Tn=(N-Nw)/L+1;%计算把数据x共分成多少段nfft=32;%FFT的长度TF=zeros(Tn,nfft);%将存放三维谱图，先清零fori=1:Tnxw=x((i-1)*10+1:i*10+10);%取一段数据temp=fft(xw,nfft);%FFT变换temp=fftshift(temp);%频谱以0频为中心forj=1:nfft;TF(i,j)=temp(j);%把谱图存放在TF中endendsubplot(3,3,6)fnew=((1:nfft)-nfft/2)*fs/nfft;tnew=(1:Tn)*L*Ts;[F,T]=meshgrid(fnew,tnew);mesh(F,T,abs(TF))title('短时傅立叶变换时频图')subplot(3,3,7)contour(F,T,abs(TF))title('等高线表示')enddisp('正在计算测试模板的参数...')fori=1:10fname=sprintf('%db.wav',i-1);x=wavread(fname);[x1x2]=vad(x);m=mfcc(x);m=m(x1-2:x2-4,:);test(i).mfcc=m;enddisp('正在进行模板匹配...')dist=zeros(10,10);fori=1:10forj=1:10dist(i,j)=dtw(test(i).mfcc,ref(j).mfcc);endenddisp('正在计算匹配结果...')fori=1:10[d,j]=min(dist(i,:));fprintf('测试信号%d的识别结果为:%d\n',i-1,j-1);end各子程序模块：dtw.mfunctiondist=dtw(t,r)n=size(t,1);m=size(r,1);%帧匹配距离矩阵d=zeros(n,m);fori=1:nforj=1:md(i,j)=sum((t(i,:)-r(j,:)).^2);endend%累积距离矩阵D=ones(n,m)*realmax;D(1,1)=d(1,1);%动态规划fori=2:nforj=1:mD1=D(i-1,j);ifj1D2=D(i-1,j-1);elseD2=realmax;endifj2D3=D(i-1,j-2);elseD3=realmax;endD(i,j)=d(i,j)+min([D1,D2,D3]);endenddist=D(n,m);enframe.mfunctionf=enframe(x,win,inc)nx=length(x(:));nwin=length(win);if(nwin==1)len=win;elselen=nwin;endif(nargin3)inc=len;endnf=fix((nx-len+inc)/inc);f=zeros(nf,len);indf=inc*(0:(nf-1)).';inds=(1:len);f(:)=x(indf(:,ones(1,len))+inds(ones(nf,1),:));if(nwin1)w=win(:)';f=f.*w(ones(nf,1),:);endmelbankm.mfunction[x,mn,mx]=melbankm(p,n,fs,fl,fh,w)ifnargin6w='tz';ifnargin5fh=0.5;ifnargin4fl=0;endendendf0=700/fs;fn2=floor(n/