ocnodes(sys),

function [Tp,ThtP,rmax] = ocnodes(sys),
%OCNODES	Computes performance period and time offset
%
%	 [Tp,ThtP] = OCNODES(sys)
%
%		Tp is the performance period
%		ThtP is the time offset 
%       sys is the structure as given by OCstruct

%Extract initial conditions
n=1;xo=sys.xo(n);vo=sys.vo(n);
%Compute relations
wo = (n.*pi).^2*sys.Alpha;
R = -sys.Gamma*wo;
u = (wo.^2-R.^2).^.5;
Tp = pi./u;%performance period
if xo~=0,
    V=vo/xo/wo;
else,
    V=inf;
end
%root 1
r1=mod(atan(u*vo/(R^2*xo+xo*u^2-R*vo))/u,Tp);
if r1==0,r1=Tp;end
u1=1;u2=1;
%root 2
r2=mod(-atan(u*(R^2*u2*xo-2*R*u2*vo-xo*u1+xo*u2*u^2)/(u2*xo*R^3-R^2*u2*vo+...
    R*xo*u2*u^2+R*xo*u1-u1*vo+u2*vo*u^2))/u,Tp);
if r2==0,r2=Tp;end
%Find min and max root
if wo>=(u1/u2)^.5,
   ThtP=r1;rmax=r2;
else
   ThtP=r2;rmax=r1;
end