【信号处理教程】基于倒谱图判断浊音的基音周期附MATLAB代码

clcclear allclose allwaveFile = 'qinghua.wav ';[y, fs] =audioread(waveFile);time=(1:length(y))/fs;frameSize=floor(40*fs/1000);     %帧长startIndex=round(15000);         %起始序号endIndex=startIndex+frameSize-1; %结束序号frame = y(startIndex:endIndex);  %取出该帧frameSize=length(frame);frame2=frame.*hamming(length(frame));  % 加hamming窗rwy = rceps(frame2);                   % 求倒谱ylen=length(rwy);cepstrum=rwy(1:ylen/2); %基音检测LF=floor(fs/500);     %设置基音搜索的范围HF=floor(fs/70);cn=cepstrum(LF:HF);[mx_cep ind]=max(cn); %设置门限，找到峰值位置if mx_cep > 0.08 & ind >LF    a= fs/(LF+ind);else    a=0;endfigure(1);plot(time, y); title(waveFile); axis tightylim=get(gca, 'ylim');line([time(startIndex), time(startIndex)], ylim, 'color', 'r');line([time(endIndex), time(endIndex)], ylim, 'color', 'r');title('语音波形');figure(2);subplot(2,1,1);plot(frame);title('取出帧的波形');subplot(2,1,2);plot(cepstrum);title('倒谱图');[x,sr]=audioread('qinghua.wav');  meen=mean(x);x= x - meen;updRate=floor(20*sr/1000);          %每20ms更新fRate=floor(40*sr/1000);            %40ms一帧n_samples=length(x);nFrames=floor(n_samples/updRate)-1; %帧数k=1;pitch=zeros(1,nFrames);f0=zeros(1,nFrames);LF=floor(sr/500);HF=floor(sr/70);m=1;avgF0=0;for t=1:nFrames        yin=x(k:k+fRate-1);        cn1=rceps(yin);        cn=cn1(LF:HF);        [mx_cep ind]=max(cn);        if mx_cep > 0.08 & ind >LF           a= sr/(LF+ind);        else            a=0;        end         f0(t)=a;        if t>2 & nFrames>3   %中值滤波对基音轨迹图进行平滑           z=f0(t-2:t);           md=median(z);           pitch(t-2)=md;           if md > 0             avgF0=avgF0+md;             m=m+1;           end        else            if nFrames<=3            pitch(t)=a;            avgF0=avgF0+a;            m=m+1;           end         end     k=k+updRate;endfigure(3)subplot(211);plot((1:length(x))/sr, x);ylabel('幅度');xlabel('时间');subplot(212);xt=1:nFrames;xt=20*xt;plot(xt,pitch)xlim([0,3]);axis([xt(1) xt(nFrames) 0 max(pitch)+50]);ylabel('基音频率/HZ');xlabel('时间');