
function psi=phikeyps(zet)
if zet<0
   c=[1 -15*zet 0 0 -1];
   q=roots(c);
   psi=q(4);
else
	%psi=-4.7*zet;
	c=min(50,.35*zet);
	psi=-(zet+.667*(zet-14.28)/exp(c)+9.525);

end;