function pdisc(per) * * This is the discriminant for the periods (=0 at a period) * implicit double precision(a-h,o-z) common/trap/per0,perH,perHe,alpha2,ane,l pi=3.14159265358979 ell=float(l) per0l=per0/sqrt(ell*(ell+1.)) perHl=perH/sqrt(ell*(ell+1.)) x1=(2.*per/perHl-2*per/per0l+ell)*pi/2. x2=(2.*per/perHl-ane-0.5)*pi/2. pdisc=cos(x1)/sin(x1)-alpha2*cos(x2)/sin(x2) return end ************************************************************************