Home > fcn > glotlf.m

glotlf

PURPOSE ^

GLOTLF Liljencrants-Fant glottal model U=(D,T,P)

SYNOPSIS ^

function u=glotlf(d,t,p)

DESCRIPTION ^

GLOTLF   Liljencrants-Fant glottal model U=(D,T,P)
 d is derivative of flow waveform: must be 0, 1 or 2
 t is in fractions of a cycle
 p has one row per output point
    p(:,1)=open phase [0.6]
    p(:,2)=+ve/-ve slope ratio [0.1]
    p(:,3)=closure time constant/closed phase [0.2]
 Note: this signal has not been low-pass filtered
 and will therefore be aliased

 Usage example:    ncyc=5;
            period=80;
            t=0:1/period:ncyc;
            ug=glotlf(0,t);
            plot(t,ug)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function u=glotlf(d,t,p)
0002 %GLOTLF   Liljencrants-Fant glottal model U=(D,T,P)
0003 % d is derivative of flow waveform: must be 0, 1 or 2
0004 % t is in fractions of a cycle
0005 % p has one row per output point
0006 %    p(:,1)=open phase [0.6]
0007 %    p(:,2)=+ve/-ve slope ratio [0.1]
0008 %    p(:,3)=closure time constant/closed phase [0.2]
0009 % Note: this signal has not been low-pass filtered
0010 % and will therefore be aliased
0011 %
0012 % Usage example:    ncyc=5;
0013 %            period=80;
0014 %            t=0:1/period:ncyc;
0015 %            ug=glotlf(0,t);
0016 %            plot(t,ug)
0017 
0018 
0019 %      Copyright (C) Mike Brookes 1998
0020 %
0021 %      Last modified Thu Apr 30 17:22:00 1998
0022 %
0023 %   VOICEBOX home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
0024 %
0025 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0026 %   This program is free software; you can redistribute it and/or modify
0027 %   it under the terms of the GNU General Public License as published by
0028 %   the Free Software Foundation; either version 2 of the License, or
0029 %   (at your option) any later version.
0030 %
0031 %   This program is distributed in the hope that it will be useful,
0032 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0033 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0034 %   GNU General Public License for more details.
0035 %
0036 %   You can obtain a copy of the GNU General Public License from
0037 %   ftp://prep.ai.mit.edu/pub/gnu/COPYING-2.0 or by writing to
0038 %   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
0039 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0040 
0041 if nargin < 2
0042   tt=(0:99)'/100;
0043 else
0044   tt=t-floor(t);
0045 end
0046 u=zeros(size(tt));
0047 de=[0.6 0.1 0.2]';
0048 if nargin < 3
0049   p=de;
0050 elseif length(p)<2
0051   p=[p(:); de(length(p)+1:2)];
0052 end
0053 
0054 te=p(1);
0055 mtc=te-1;
0056 e0=1;
0057 wa=pi/(te*(1-p(3)));
0058 a=-log(-p(2)*sin(wa*te))/te;
0059 inta=e0*((wa/tan(wa*te)-a)/p(2)+wa)/(a^2+wa^2);
0060 
0061 % if inta<0 we should reduce p(2)
0062 % if inta>0.5*p(2)*(1-te) we should increase p(2)
0063 
0064 rb0=p(2)*inta;
0065 rb=rb0;
0066 
0067 % Use Newton to determine closure time constant
0068 % so that flow starts and ends at zero.
0069 
0070 for i=1:4
0071   kk=1-exp(mtc/rb);
0072   err=rb+mtc*(1/kk-1)-rb0;
0073   derr=1-(1-kk)*(mtc/rb/kk)^2;
0074   rb=rb-err/derr;
0075 end
0076 e1=1/(p(2)*(1-exp(mtc/rb)));
0077 
0078 
0079 ta=tt<te;
0080 tb=~ta;
0081 
0082 if d==0
0083   u(ta)=e0*(exp(a*tt(ta)).*(a*sin(wa*tt(ta))-wa*cos(wa*tt(ta)))+wa)/(a^2+wa^2);
0084   u(tb)=e1*(exp(mtc/rb)*(tt(tb)-1-rb)+exp((te-tt(tb))/rb)*rb);
0085 elseif d==1
0086   u(ta)=e0*exp(a*tt(ta)).*sin(wa*tt(ta));
0087   u(tb)=e1*(exp(mtc/rb)-exp((te-tt(tb))/rb));
0088 elseif d==2
0089   u(ta)=e0*exp(a*tt(ta)).*(a*sin(wa*tt(ta))+wa*cos(wa*tt(ta)));
0090   u(tb)=e1*exp((te-tt(tb))/rb)/rb;
0091 else
0092   error('Derivative must be 0,1 or 2');
0093 end

Generated on Tue 27-Mar-2007 12:06:24 by m2html © 2003