program impsamp
      implicit none
      integer nsamp,isamp,k
      real a,sa,s2pia,exa,xm,x2,chi,x,e,rand,dpi,a0
      parameter (nsamp=10000)

      write (*,*) 'input a0'
      read (*,*) a0

      do k=1,10
      a=a0+(k-1)*0.1
      sa = sqrt(a)
      dpi = 8.*atan(1.)
      s2pia=sqrt(dpi*a)
      exa=(1./a-1.)/2.
      
      xm=0.
      x2=0.
      do isamp=1,nsamp
       chi=sqrt(-2.*log(1.-rand()))*cos(dpi*rand())
       x=sa*chi
       e=s2pia*exp(x*x*exa)/(x*x+1.)
!       write (*,*) isamp,e
       xm=xm+e
       x2=x2+e**2
      enddo
      write (*,*) 'a=',a,' mean = ',xm/nsamp,' +/- ',
     .  sqrt( (x2/nsamp-(xm/nsamp)**2)/(nsamp-1) ) 
      write (10,*) a,xm/nsamp,sqrt( (x2/nsamp-(xm/nsamp)**2)/(nsamp-1) )
      enddo
      stop 
      end