• 语音信号的产生模型和清浊音特性
  空气从肺部排除行成气流，气流激励升到从嘴唇或鼻孔或同时从嘴唇和鼻孔辐射出来而产生语音。当气流通过声门时，如果声带的张力刚好使声带张弛震荡，产生一股准周期的气流，这一气流激励声道就产生浊音语音。当气流通过声门时，声道的某处发生了收缩而产生一个狭窄的通道，当空气流到达此处时被迫以高速通过紧缩部分而在附近产生气流湍流，这种湍流通过声道便行成清音语音。
  浊音语音信号为周期冲击序列，其周期为。清音语音信号为白噪声序列，满足均值为0，方差为1的高斯分布。
  语音产生的“源-系统”模型（根据发声器官和语音产生过程）：
  
  可实现的数字模型（激励源、声道模型和辐射模型）：
  
• 检测基音周期的基本原理
  基音是指发浊音时声带振动所引起的周期性，基音周期是声带振动频率的倒数。基音检测是语音处理中一个非常重要的问题，找到一个完善的适用于不同的讲话者、不同要求和环境、准确和可靠的基音检测的方法是语音处理领域最具挑战性的任务之一。目前,大多数基音检测方法都采用的是自相关法。
  对于浊音信号的互相关函数，在基音周期的整数倍位置上出现峰值，而清音的自相关函数没有明显的峰值出现；因此检测是否有峰值就可判断是清音或浊音，检测峰值的位置就可以提取基音周期值。
• 程序框图
• 读入语音信号
  已知语音数据文件是经8KHz采样频率采样的信号。

clear;clc;
fid=fopen('female.dat');
x1=fread(fid,16000,'short');
figure;plot(x1);title('读信号');
grid

采用fread函数以二进制数据读入数据。设定信号长度为16000，精度为16位整数型。语音信号图如下：

• 语音信号分帧
  设20ms的语音信号为一帧，计算出帧长为160个样点，并设置帧移为160。用矩形窗作为分帧的窗函数。为方便后续计算，对分帧后的矩阵作转置，计算出帧数。

%分帧
wlen=8000*20*10.^(-3); inc=160;          % 给出帧长和帧移
win=boxcar(wlen);         % 给出矩形窗
N=length(x1);               % 信号长度
time=(0:N-1)/8000;           % 计算出信号的时间刻度
X=enframe(x1,win,inc);     % 分帧
X=X';
fn=size(X,2);              % 求出帧数

• 去均值
  当语音信号在分析窗里有非零均值或有非常低的低频噪声出现时，归一化自相关函数在所要求的所有延迟上均产生高的相关，这不利于计算依赖归一化自相关函数进行清浊分类的安静段语音或低频幅度清音语音段。为此，在计算归一化自相关函数时首先要从分析窗数据中减掉均值。利用mean函数可以计算矩阵特定区域的均值。

%去均值
for p=1:89
    ave=mean(X(:,p));
    X(:,p)=X(:,p)-ave;
    for o=1:160
        XX((p-1)*160+o)=X(o,p);
    end
end
figure;plot(1:14240,XX,'k'); axis([5000 5400 ,-inf,inf]);xlabel('时间(s)');ylabel('时间波形');title('展开信号');
hold on;

• 800Hz低通滤波
  为减少高频共振峰和外来高频噪声的影响，对去均值的语音信号进行800Hz低通滤波，以去掉大部分共振峰的影响，又可以当基音频率为最高500Hz时，仍能保持其一二次谐波。依据双线性变换法设计的800Hz5阶椭圆低通滤波器传递函数为：
  
  默认512点插值，b1和a1为系统传递函数分子和分母的系数向量。构成滤波器并对去均值后的信号矩阵X进行滤波，生成矩阵sf_1。

b1=[0.008233 -0.004879 0.007632 -0.004879 0.008233];
a1=[1 -3.6868 5.8926 -5.0085 2.2518 -0.4271];
[H,w]=freqz(b1,a1,512);
for r=1:89
     sf_1(:,r)=filter(b1,a1,X(:,r));
end
plot(1:14240,sf,'b');xlabel('时间（s）');ylabel('幅值');title('椭圆滤波');axis([5000 5400 ,-inf,inf]) %改变横轴范围
hold on;

• 数值滤波
  经800Hz低通滤波的语音信号，主要去除了第三和第四个高频共振峰及高频噪声的影响，但第一和第二个共振峰有时仍然存在，是的浊音段语音信号的周期性模糊，产生了错误的基音估计。实验发现，若在800Hz低通滤波器之后级联一个宽度N=9的数值滤波器可有效去除这一影响，突出了浊音语音信号的周期性，使基音估计可靠。数值滤波器的传递函数为：
  
  默认512点插值，b2和a2为系统传递函数分子和分母的系数向量。构成滤波器并对去均值后的信号矩阵sf_1进行滤波，生成矩阵sf_2。下图分别为语音信号预处理图：

b2=[1 1 1 1 1 1 1 1 1];
a2=[9 0 0 0 0 0 0 0 0];
[H,w]=freqz(b2,a2,512); 
for r=1:89
    sf_2(:,r)=filter(b2,a2,sf_1(:,r));
end
plot(1:14240,sf,'r');xlabel('时间（s）');ylabel('幅值');title('数值滤波');axis([5000 5400 ,-inf,inf]) %改变横轴范围
legend('去均值','低通滤波','数值滤波');

• 计算语音信号短时平均能量ELP
  语音信号短时平均能量由以下等式求得：
  

由于
是一个极小值，因此在此程序里任取为0.0000000001。在此给出各帧短时平均能量ELP计算代码和结果图：

for i=1 : fn
   for z=2:wlen
     S_LPN(1,1)=sf(1,(1+(i-1)*inc))*sf(1,(1+(i-1)*inc)); 
     S_LPN(1,z)=S_LPN(1,(z-1))+sf(1,(z+(i-1)*inc))*sf(1,(z+(i-1)*inc));
    end
    En(i)=10*log10(S_LPN(1,160)/160+0.0000000001);         
end

• 计算归一化互相关函数NCCF、后处理、清浊判决、基音周期的确定和输出
• 互相关函数NCCF
  令语音信号表示为s(n)，定义波形相关测度为：
  
  这里ρ(τ)称为归一化互相关函数，并作为基音检测的依据。
• 后处理

为了避免选择实际基音的倍数作为估计值，可采用如下后处理措施。

在三个基音区域内计算ρ(τ)，分别得到各自区域的最大值及对应的延迟,记为ρmax1,ρmax2,ρmax3和T1，T2，T3。三个区域分别为：
区域1：80 ~ 147
区域2：40 ~ 79
区域3：20 ~ 39
按如下方法确定最优基音延迟Topt
Topt = T1
ρmax = ρmax1
if ρmax2 ≥ cρmax
ρmax = ρmax2
Topt = T2
if ρmax3 ≥ cρmax
ρmax = ρmax3
Topt = T3
end
end
c为经验因子，经过测试，发现男性取值为0.92，女性取值为0.94时比较合适。

• 清浊判决
  利用前面定义的归一化互相关函数可进行有效的清浊判决。
  第一步，确定低通数值滤波后的语音对数能量ELP（dB）：
  
  第二步，确定周期性水平量Zperiod，设前三个峰的平均值：
  
  则周期性水平量Zperiod定义为：
  
  第三步，取ELP和Zperiod作为清浊判决。当ELP≤Eth时，语音能量太小直接判决为清音。当ELP>Eth但Zperiod≤Zth时，说明周期性不强判决为清音。如果ELP>Eth，Zperiod>Zth且Topt在允许的基频范围内，判为浊音，并输出对应的最优化基音延迟Topt。通过实验测试，取男性Eth=32.62dB，女性Eth=55.5dB，男性Zth=0.79，女性Zth=0.62时能正确分别清浊音。
• 基音周期的确定和输出

for m=2: fn            
%取一帧 
    if En(1,m)>32.62
        for tao=20:39        %对于不同的延迟tao
            rxx3(1,1)=sf(1+(m-1)*inc)*sf(1-tao+(m-1)*inc);
            s_n(1,1)=sf(1+(m-1)*inc)*sf(1+(m-1)*inc);
            s_n_tao3(1,1)=sf(1-tao+(m-1)*inc)*sf(1-tao+(m-1)*inc);
           for j=2:wlen   %分别计算分子和分母
               rxx3(1,j)=rxx3(j-1)+sf(j+(m-1)*inc)*sf(j-tao+(m-1)*inc);  %分子
               s_n(1,j)=s_n(j-1)+sf(j+(m-1)*inc)*sf(j+(m-1)*inc);
               s_n_tao3(1,j)=s_n_tao3(j-1)+sf(j-tao+(m-1)*inc)*sf(j-tao+(m-1)*inc);
           end
           Rxx3(1,tao-19)=rxx3(1,160)/((s_n(1,160)*s_n_tao3(1,160))^0.5);
        end
       [Rxx_tao_max3,index3]=max(Rxx3);  %最大互相关函数峰值=Rxx_tao-max、延迟tao=index
        for tao=40:79        %对于不同的延迟tao
            rxx2(1,1)=sf(1+(m-1)*inc)*sf(1-tao+(m-1)*inc);
           s_n_tao2(1,1)=sf(1-tao+(m-1)*inc)*sf(1-tao+(m-1)*inc);
           for j=2:wlen   %分别计算分子和分母           
               rxx2(j)=rxx2(j-1)+sf(j+(m-1)*inc)*sf(j-tao+(m-1)*inc);  %分子
               s_n_tao2(j)=s_n_tao2(j-1)+sf(j-tao+(m-1)*inc)*sf(j-tao+(m-1)*inc);
           end
           Rxx2(tao-39)=rxx2(160)/((s_n(160)*s_n_tao2(160))^0.5);       
        end
        [Rxx_tao_max2,index2]=max(Rxx2);  %最大互相关函数峰值=Rxx_tao-max、延迟tao=index
        for tao=80:147        %对于不同的延迟tao
            rxx1(1,1)=sf(1+(m-1)*inc)*sf(1-tao+(m-1)*inc);
            s_n_tao1(1,1)=sf(1-tao+(m-1)*inc)*sf(1-tao+(m-1)*inc);
           for j=2:wlen   %分别计算分子和分母
               rxx1(j)=rxx1(j-1)+sf(j+(m-1)*inc)*sf(j-tao+(m-1)*inc);  %分子
               s_n_tao1(j)=s_n_tao1(j-1)+sf(j-tao+(m-1)*inc)*sf(j-tao+(m-1)*inc);
           end
           Rxx1(tao-79)=rxx1(160)/((s_n(160)*s_n_tao1(160))^0.5);       
        end
        [Rxx_tao_max1,index1]=max(Rxx1); %最大互相关函数峰值=Rxx_tao-max、延迟tao=index
        %最优基因延迟Topt
        Rxx_tao_max=Rxx_tao_max1;Topt=index1+79;
        c=0.92;
        if (Rxx_tao_max2>=c*Rxx_tao_max)
            Rxx_tao_max=Rxx_tao_max2;Topt=index2+39;
        end
        if (Rxx_tao_max3>=c*Rxx_tao_max)
            Rxx_tao_max=Rxx_tao_max3;Topt=index3+19;
        end
        P_max(1,m)=Rxx_tao_max;
        %清浊判决
        P_avr=(Rxx_tao_max1+Rxx_tao_max2+Rxx_tao_max3)/3;     %前三个峰的平均值
        Z_period=Rxx_tao_max+P_avr;       %周期性水平量
        if (Z_period>0.79) %Topt在允许的基频范围内
           Topt_tao(1,m)=Topt;
        else
           Topt_tao(1,m)=0;
        end
    else 
        Topt_tao(1,m)=0;
    end
end
figure;plot(1:97,Topt_tao)     % 画出基音周期
title('基音周期');ylabel('基音周期值'); xlabel('帧数');
axis([0 100,0 120]) %改变横轴范围
hold on;

更多精彩内容