PROGRAM angle
C======================================================================
C     THIS PROGRAM CALCULATES epsilon, and from there delta and psi  
C======================================================================
      REAL wmn,wmx,x,epsinf,wp,gp
     * ,teta,pi,sn1,cs1,sn2,cs2,rho,psi,delta
      COMPLEX eps,n1,n2,rp,rs,rps
      INTEGER i,nw
      read(*,*) epsinf,wp,gp
      READ(*,*) wmn,wmx,nw
      read(*,*) teta
      pi=4.*atan(1.)
      teta=teta*pi/180
      sn1=sin(teta)
      sn2=sn1**2
      cs2=1.-sn2
      cs1=sqrt(cs2)
      do 10 i=1,nw
       x=wmn+real(I-1)*(wmx-wmn)/real(nw)
       w=abs(x)
       eps=epsinf-wp**2/(w*(w+cmplx(0.,gp)))
       n1=eps*cs1
       n2=csqrt(eps-sn2)
       rp=(n1-n2)/(n1+n2)
       n1=cs1
       rs=(n1-n2)/(n1+n2)
       rps=rp/rs
       rho=real(cabs(rps))
       psi=atan(rho)
       delta=aimag(clog(rps))
       if (x.gt.0) then
         write(*,*) x,psi
       else
         write(*,*) x,delta
       endif
10    continue
      END