polym1.m

这是一个用于语音信号处理的工具箱

%FUNCTION: model and generate the differential glottal pulse for 6th order 
%          polynomial model.
%        [cof,dgpuls]=polym1(dgf) returns the extimated differential glottal
%        waveform and coefficient vector that fits the polynomial model 
%        with amplitude normalization. This method emphasizes the constraint: 
%        f(0)=f(1);
%
% INPUT: dgf = the differential glottal flow , i.e. the integral of 
%              residual signal found by inverse filtering.
% OUTPUT:  cof = coefficients for fitting a 5-th order polynomial
%       dgpuls = an estimated differential glottal waveform
%
% SEE ALSO: gen_dgf1.

function [cof,dgpuls]=polym1(dgf);

%*********************************************
% 1. set the value of the starting point as zero
%    and normalize the amplitude to 100
%*********************************************
xlens=length(dgf)-1;

% set starting point as zero -- step(1)
 startp=dgf(1);
 dgf=dgf-startp;
 
% reset the value of final point -- step(2)
 xx=0:1/xlens:1;
 lastp=dgf(xlens+1);
 dgf=dgf-lastp*xx;

%normalize the amplitude -- step(3)
 amp1=max(dgf);
 dgf=dgf*100/amp1;
 
%*************************************************
% 2. weighting function which emphasizes f(0)==f(1) -- step(4) 
%*************************************************
 
 for i=1:xlens+1
     x=xx(i);
     if x<0.1
        weit(i)=200*x*x-40*x+3;
     elseif x<0.8
        weit(i)=1;
     else
        weit(i)=25*x*x-40*x+17;
     end
 end
 dgf1=dgf.*weit;

%******************************************************
% 3. smooth the curve using a 6th order polynomial fit -- step(5)
%******************************************************
  
 cof=polyfit( xx,  dgf1, 6);

if nargout<2
   return;
end

%********************************************
% 3. modify C0 such that the sum will be zero -- step(6)
%********************************************
%cof(7)=0;
%for i=1:6
%    cof(7)=cof(7)-cof(i)/(8-i);
%end
 
%*****************************************************
% 4. create the glottal pulse with this polynomial model

% recovering for step (5)
 dgpuls=polyval(cof,xx);

% recovering for step (4)
 dgpuls=dgpuls./weit;

% recovering from step (3)
 dgpuls=dgpuls/100*amp1;

% recovering from step (2)
 dgpuls=dgpuls+lastp*xx;

%recovering from step (1)
 dgpuls=dgpuls+startp;

% Produce a high-energy pulse at glottal closure, follow Dr.Hu's approach
maxg=max(dgpuls(1:xlens/5));
dgpuls(3)=(maxg+1)*.42-1;
dgpuls(4)=.5*dgpuls(3)+.5*dgpuls(4);

% Remove the spectral tilt of the excitation pulse by inverse filtering
ss=dgpuls(2:xlens+1);
ss=ss(:);
energy=ss'*ss;
rc1=ss(1:xlens-1)'*ss(2:xlens)/energy; % first order reflection coefficient
%dgpuls=filter([1 -.7*rc1],1,dgpuls);  % first-order inverse filtering

%retreat the startint point
% dgpulse(1)=dgpuls(xlens+1);