LPCRF2AR Convert reflection coefs to autoregressive coefs [AR,ARP,ARU,G]=(RF) Input: RF(:,p+1) gives reflection coefficients of one or more p-section lossless tubes Ouputs: G is the gain of the all-pole AR filter AR/G is the transfer function from the output flow to the glottal input wave. ARP*K is the transfor function from the output flow to the pressure at the glottis where K = rho*c/Alips: rho = air density 1.23 kg/m^3, c=sound speed 340 m/s, Alips = effective area of free space beyond the lips. ARU is the transfer function from the output flow to the total input flow Energy into the vcal tract is equal to K*filter(ARP,1,Ulips).*filter(ARU,1,Ulips) reverse glottal flows divided by 1-r0 where r0 is the glottal reflection coefficient. The scale factor is included to avoid a zero answer when the glottis is closed giving r0=1. The transfer functions have ar(:,1)=art(:,1)=1 They should both be multiplied by z^(p/2)/prod(1+rf) to correct the absolute gain and to compensate for the delay of p/2 samples along the length of the tube. The energy into the vocal tract is given by ars(speech) * are(speech)
0001 function [ar,arp,aru,g]=lpcrf2ar(rf) 0002 %LPCRF2AR Convert reflection coefs to autoregressive coefs [AR,ARP,ARU,G]=(RF) 0003 % 0004 % Input: RF(:,p+1) gives reflection coefficients of one or more p-section lossless tubes 0005 % Ouputs: G is the gain of the all-pole AR filter 0006 % AR/G is the transfer function from the output flow to the glottal input wave. 0007 % ARP*K is the transfor function from the output flow to the pressure at the glottis 0008 % where K = rho*c/Alips: rho = air density 1.23 kg/m^3, c=sound speed 340 m/s, 0009 % Alips = effective area of free space beyond the lips. 0010 % ARU is the transfer function from the output flow to the total input flow 0011 % 0012 % Energy into the vcal tract is equal to K*filter(ARP,1,Ulips).*filter(ARU,1,Ulips) 0013 % reverse glottal flows divided by 1-r0 where r0 is the glottal reflection coefficient. 0014 % The scale factor is included to avoid a zero answer when the glottis is closed giving r0=1. 0015 % 0016 % The transfer functions have ar(:,1)=art(:,1)=1 0017 % They should both be multiplied by z^(p/2)/prod(1+rf) to correct the absolute 0018 % gain and to compensate for the delay of p/2 samples along the length of the tube. 0019 % 0020 % The energy into the vocal tract is given by ars(speech) * are(speech) 0021 0022 0023 % Copyright (C) Mike Brookes 1997 0024 % 0025 % Last modified Tue May 12 16:19:13 1998 0026 % 0027 % VOICEBOX home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html 0028 % 0029 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0030 % This program is free software; you can redistribute it and/or modify 0031 % it under the terms of the GNU General Public License as published by 0032 % the Free Software Foundation; either version 2 of the License, or 0033 % (at your option) any later version. 0034 % 0035 % This program is distributed in the hope that it will be useful, 0036 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0037 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0038 % GNU General Public License for more details. 0039 % 0040 % You can obtain a copy of the GNU General Public License from 0041 % ftp://prep.ai.mit.edu/pub/gnu/COPYING-2.0 or by writing to 0042 % Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA. 0043 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0044 0045 [nf,p1]=size(rf); 0046 p2=p1+1; 0047 p=p1-1; 0048 pm=p-1; 0049 arf=[ones(nf,1) zeros(nf,p)]; 0050 arr=[zeros(nf,p) rf(:,p1)]; 0051 cr=zeros(nf,p); 0052 for k=1:p-1 0053 rk=rf(:,(p1-k)*ones(1,k)); 0054 cr(:,1:k)=arr(:,p2-k:p1); 0055 arr(:,p1-k:p)=arr(:,p1-k:p)+rk.*arf(:,1:k); 0056 arf(:,2:k+1)=arf(:,2:k+1)+rk.*cr(:,1:k); 0057 end 0058 r1=rf(:,1); 0059 ar=arf+r1(:,ones(1,p1)).*arr; 0060 if nargout>1 0061 kp=prod(1-rf(:,2:p1),2); 0062 arp=(arf-arr)./kp(:,ones(1,p1)); 0063 if nargout>2 0064 g=prod(1+rf(:,2:p1),2); 0065 aru=(arf+arr)./g(:,ones(1,p1)); 0066 if nargout>3 0067 g=g.*(1+rf(:,1)); 0068 end 0069 end 0070 end 0071 0072