function [L, xp] = PitchOL (x, LSubframe, PitchOLpar)% Find the open loop pitch lag% $Id: PitchOL.m 1.1 2003/11/21 G.723.1-v2r1a $% Form the extended vectorxe = [PitchOLpar.xp; x];% Form open-loop estimatesNSubframe = length (LSubframe);LMem = length (PitchOLpar.xp);j = LMem;for (i = 1:NSubframe)  Nx = LSubframe(i);  L(i) = PitchEst (xe, j, Nx, PitchOLpar);  j = j + Nx;end% Memory for next framexp = xe(end-LMem+1:end);return%---------------------function Lo = PitchEst (x, iref, N, PitchOLpar)PMin = PitchOLpar.PMin;PMax = PitchOLpar.PMax;PMultThr = PitchOLpar.PMultThr;ELo = 1;CLo = 0;Lo = PMin;L = PMin - 1;i1 = iref - L;i2 = i1 + N - 1;EL = x(i1+1:i2+1)' * x(i1+1:i2+1);	% Energy for pitch lag PMin-1for (L = PMin:PMax)  i1 = i1 - 1;  i2 = i2 - 1;  % Recursive calculation of the energy term  EL = EL + (x(i1+1)^2 - x(i2+2)^2);  % Cross correlation term  CL = x(i1+1:i2+1)' * x(iref+1:iref+N);  % The error for a pitch predictor is  %   e(n) = x(n) - p x(n-L)  % The squared error for a frame is  %   ES(L) = E(0) - 2p C(L) + p^2 E(L),  % where E(L) = SUM x(n-L)^2 and C(L) = SUM x(n) x(n-L).  % For a given L, the optimal coefficient p is  %   popt = C(L) / E(L).  % For this choice of p, the squared error is  %  ES(L) = E(0) - C(L)^2 / E(L).  % In the following code we find L which maximizes C(L)^2 / E(L). For  % each L, we compare the largest values found so far (at Lo) with a  % value at L and choose the new value of L if  %  C(L)^2   C(Lo)^2  %  ------ > -------  or C(L)^2 E(Lo) > C(Lo)^2 E(L)  %   E(L)     E(Lo)  % The search is conducted from small values of L up. Only positive values  % of C(L) are candidates for a pitch lag. If the better value of squared  % error is within Pmin of the Lo, it is kept. If it further away, the new  % value must be at least 1/A times better.  %  A C(L)^2 E(Lo) > C(Lo)^2 E(L)  % This second check is used to avoid locking onto pitch multiples.  if (EL > 0 & CL > 0)    %  ( [ ------------------------- & --------------- ] | ...    if ( ( (CL^2 * ELo > EL * CLo^2) & ((L-Lo) < PMin) ) | ...           (PMultThr * CL^2 * ELo > EL * CLo^2) )    %      ------------------------------------ )      Lo = L;      ELo = EL;      CLo = CL;    end  endendreturn