c     constants

      double precision pi, m_n, h_par, k_b, avog
      double precision he4mol, gam, rgas
      parameter	(pi		= 3.14159,
     +		 m_n            = 0.1045407751,	! meV/(Angstrom/psec)^2
     +		 h_par          = 0.6582122,	! meV*psec
     +		 k_b            = 0.08617,      ! Boltzman const (meV/K)
     +		 avog		= 0.6022,       ! Avodadro's no. (10^24 atom/mol)
     +		 he4mol		= 4.0,		! helium4 (gm/mol)
     +		 gam		= 1.6666666667,	! adiabatic const
     +		 rgas		= 8.314)	! gas constant (J/mol*K)

      double precision gy, mag_b
      parameter(gy     =18.324718,   ! 1/(mGaus*sec)
     +		mag_b  = 0.1)       ! mGaus
      double precision k_0, bb
      parameter	(k_0   = 0.75,                  ! 1/Ang
     +		 bb    = 3.07729087d-15)        ! scattering length (m)
c
c     ****Number of Events*****
      integer nevts
	  parameter (nevts = 1e5)
c
c     *********Helium Target properties**********
      double precision T 
      parameter (T  = 2.10)     ! helium temp. (K)
c     parameter (T  = 4.25)     ! helium temp. (K)
c
      double precision tarx1, tarx2, tarx3, tarx4
	  double precision tary1, tary2, tarz
	  double precision ndenhegas
	  logical hegastar, airtar

c  !!two of the following four should be commented out!!
c  !! for liquid He, both should be false below!!	  
      parameter (hegastar = .false.)  ! false if liquid He target or air
!      parameter (hegastar = .true.)   ! true is gas He target
      parameter (airtar   = .false.)  ! false if liquid He target or He gas
!      parameter (airtar   = .true.)   ! true if air in target

	  parameter (ndenhegas = 0.0000250 )  !! He numb. density (1e24 atm/cc)
      parameter	(tarx1 = -3.0,    ! cm
     +           tarx2 = -0.35,         ! cm
     +           tarx3 = 0.35,          ! cm
     +           tarx4 = 3.0,           ! cm
     +           tary1 = -2.5,          ! cm
     +           tary2 = 2.5,           ! cm
     +           tarz  = 41.6)          ! cm
c     ************************************
c
c     *********Aluminum Window properties**********
      double precision ndenAl, m_al
      parameter (ndenAl  = 0.06026)     ! Al numb. density (1e24 atm/cc)
      parameter (m_al = 2.82)				! Al mass - meV/(Angstrom/psec)^2
c
      double precision talx1, talx2, talx3, talx4, taly1, taly2
	  double precision talz1, talz2, talz3
      parameter	(talx1 = -3.0,         ! cm
     +           talx2 = 3.0,          ! cm
     +           talx3 = -2.8,          ! cm
     +           talx4 = 2.8,           ! cm
     +           taly1 = -3.0,          ! cm
     +           taly2 = 3.0,           ! cm
     +           talz1 = 0.1,           ! cm
     +           talz2 = 0.01)          ! cm
c     ************************************
c
c    *******Air Gaps***********
      double precision m_air, ndenair
	  parameter (ndenair = 0.000054)  !! air numb. density (1e24 atm/cc)
	  parameter (m_air = 1.57)     ! air mass - meV/(Angstrom/psec)^2

      double precision tairx1,tairx2,tairy1,tairy2,tairz1,tairz2
	  double precision tairz3,tairz4,tairz5
      parameter	(tairx1 = -4.0,          ! cm
     +           tairx2 = 4.0,           ! cm
     +           tairy1 = -4.0,          ! cm
     +           tairy2 = 4.0,           ! cm
     +           tairz1 = 3.0,         ! gap after PSM  cm
     +	         tairz2 = 10.0,    ! gap between incoil and target (cm)
     +           tairz3 = 8.0,          ! cm
     +           tairz4 = 13.0, 
     +           tairz5 = 10.0 )         ! cm
c
c	  **BEAM GEOMETRY*******************************
c
c     ******Initial x, y-range and guide length***
	  double precision gl,xin,yin,ng6index
	  parameter	(xin = 2.75,           ! starting x range (E-W) after PSM
     +           yin = 2.25,           ! starting y range (high-low) after PSM
     +           gl = 111.0,           ! Guide length in z-direction (cm)
     +           ng6index = 0.002)    ! NG6 index rad/ang
!     +           ng6index = 0.0000001) ! for no beam divergence index rad/ang
c

c    *******Input Guide********
      double precision gx1, gx2, gy1, gy2
      parameter	(gx1 = -2.20,		! cm
     +           gx2 = 2.20,          ! cm
     +           gy1 = -2.25,         ! cm  
     +           gy2 = 2.25)          ! cm


c    *******Input Coil*********
      double precision ix1, ix2, iy1, iy2, iz0, gin
	  parameter	(gin = 0.0012,      ! Input critical index (rad/Ang)
     +	         ix1 = -2.2,
     +           ix2 = 2.2,
     +           iy1 = -2.25,         ! cm
     +           iy2 = 2.25,          ! cm
     +           iz0 = 89.0)          ! cm
c
c    ****************************
c
c    *******Output Coil**********
      double precision ox1, ox2, ox3, ox4, oy1, oy2, oz, got
      parameter (got = 0.0012,       ! Output critical index (rad/Ang)
     +           ox1 = -3.3255,       ! cm
     +           ox2 = -0.175,       ! cm
     +           ox3 = 0.175,        ! cm
     +           ox4 = 3.3255,       ! cm
     +           oy1 = -2.575,       ! cm
     +           oy2 = 2.575,        ! cm
     +           oz  = 125.)         ! cm
c
c    ******Analyzing Super Mirror (ASM)****
      double precision sx1, sx2, sy1, sy2, sz, gsm
      parameter	(sx1 = -2.85,        ! cm 
     +           sx2 = 2.85,         ! cm
     +           sy1 = -2.25,        ! cm
     +           sy2 = 2.25,         ! cm
     +           sz  = 28.0,         ! cm
     +           gsm = 0.003)        ! (Pyrex glass+ Ti + Ga) rad/Ang