语音中提取基音频率matlab程序.zip

%%corr.m % correlation for pitch estimation % flag =1 from left to right, otherwise, from right to left % r=corr(s,L,flag) % % %flag==1时 % % % %|--i----| %s:..................................................... % . % . % . % . % . % . %s:..................................................... % |_______d________|__________________L-d_____________| %返回值是一个序列 %注意到在做自相关时，错位值从20开始 %考虑到在用这个自相关序列时可能是用于求周期，这个20时有利的 %语音的基音周期一般在20到140，因此小于20没有周期 %而错开20以上可以避免出现不是所求的峰值避免出错 %求到结果后放到rr序列中时，注意到下标为i时实际上已经是错开 %i+20-1,因此pitch.m在调用这个时需要有一句p=19+ind_pitch; function r=corr(s,d,flag) epsilon=1.0e-10; % to avoid divided by zero L=length(s); NM=120; if flag==1 % From left to right for i=20:min(L-d,NM) temp=0.0; temp1=0.0; temp2=0.0; for j=1:d temp=temp+s(j)*s(j+i); temp1=temp1+s(j+i)*s(j+i); temp2=temp2+s(j)*s(j); end pw1=sqrt(temp1)+epsilon; pw2=sqrt(temp2)+epsilon; r(i-19)=temp/(pw1*pw2); end else % From right to left for i=20:min(L-d,NM) temp=0.0; temp1=0.0; temp2=0.0; for j=L-d+1:L temp=temp+s(j-i)*s(j); temp1=temp1+s(j-i)*s(j-i); temp2=temp2+s(j)*s(j); end pw1=sqrt(temp1)+epsilon; pw2=sqrt(temp2)+epsilon; r(i-19)=temp/(pw1*pw2); end end %%maxx.m % Find all the local maxima of a function which are greater than % t times the global maximum % max_index=maxx(s,t) %输入一个序列和一个标量t %首先求出序列s的最大值max_value %然后如果序列中的一个点比相邻两个点都大 %而且还大于max_value的t倍 %就把这个点的位置存放在序列mi中 %在序列mi最后再补上序列的最后一个点 %考虑到第一个点和最后一个点有可能满足要求 %但这两个点都不能跟旁边的比较 %所以还是把这两个点保留了，放在mi序列的第一个和最后一个 %实际上就是以最大值点的t倍为界找出这个界以上的极值点的位置 function mi=maxx(s,t) max_value=max(s); oldstate=0; L=length(s); j=1; mi=[1]; for i=1:L-1 if s(i)>s(i+1) newstate=-1; end if s(i)t*max_value j=j+1; mi(j)=i; end end oldstate=newstate; end ll=length(mi); mi(ll+1)=length(s); %%pitch.m % This is a function to find the pitch period of a long enough speech % interval. In order to find correct pitch period, it is in favor of % selecting short period. i.e., if p1p2*Thr, we select p1 % instead of p2. % % usage [p,mx]=pitch(s) %输入一个序列s,首先求出其自相关序列rr %要求出第一个周期的极大值点 %首先求出序列的最大值maxim所在的位置，如果这个位置在THR_largepitch以内 %就认为这就是要求的点。如果这个位置大于THR_largepitch，就考虑小于 %THR_largepitch的范围内有没有点满足以下条件：大于maxim*THR_MAXX的极值点， %大于max(THR_corr,THR_pitch*maxim). %如果有满足条件的，第一个点就是所求点 %至于为什么可以这样做，还没搞懂 %总的说来，这个函数返回一个序列的第一个周期 function [p,mx]=pitch(s) THR_pitch=0.75; THR_maxx=0.7; THR_corr=0.33; THR_largePitch=60; epsilon=1.0e-10; % used to avoid dividing by zero L=length(s); d=100; rr=corr(s,d,0);%最后一个参数是0，求自相关时从右到左 [maxim,mxi]=max(rr); %当rr是一个行向量时，[Y,I]=max(s);Y是s中的最大值，I是标量 %存放了最大值的位置 ind_pitch=mxi; if mxi>THR_largePitch mi=maxx(rr,THR_maxx);%实际上就是以最大值点的t倍为界找出这个界以上的左右极值点的位置 %这个序列的第一个值是1，最后一个是序列rr的长度 %一开始看到这里的时候会觉得比较多余了，直接一次maxim*THR_pitch分界比较就可以了 %事实上，下面的比较没有只针对极值点处理。 %maxx的作用是筛选出大于maxim*THR_maxx的极值点出来，所以这个还是必须的 for i=1:length(mi) if rr(mi(i))>max(THR_corr,maxim*THR_pitch); ind_pitch=mi(i); break;%break说明只取第一个满足条件的 end end end p=19+ind_pitch;%这个19参考corr.m中的解释，其实就是rr(i)对应的是i+19 mx=rr(ind_pitch); %%zerocros.m % This is a function to calculate the zero-cross of a signal % used to determine a signal is voiced or not. % % usage zc=zerocros(s) %输入一个序列，如果相邻两个点异号，统计数字加1 %这个结果可以作为判断是否是噪音的一个依据 %当count足够大时就可以认为是噪音 function zc=zerocros(s) L=length(s); count=0; for i=1:L-1 if s(i)*s(i+1)<0 count=count+1; end end zc=count;