function s = psi(t,beta)
     s = 2*sqrt(t + beta^2);