cbremErgLPM.f File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	cbremerglpm (ee, rho, eg)

Function/Subroutine Documentation

cbremerglpm()

subroutine cbremerglpm	(	real*8	ee,
		real*8	rho,
		real*8	eg
	)

Definition at line 20 of file cbremErgLPM.f.

References d0, masele, rndc(), and true.

Referenced by cbrems().

      implicit none
#include "Zmass.h"
      real*8 ee  !input.  electron energy in GeV
      real*8 rho !input.  air density in kg/m3
      real*8 eg  !output.  sampled gamma ray energy. GeV
!
      real*8  u, v, temp, s
      real*8 cmigdg, cmigdpsi, cmigdsforbrem, cmigdgzai

      common /cmigdc/ s1, logs1, vm, logvm, smpa1, smpa2, suma,
     *   a2prob, x0
      real*8 vm, logvm, smpa1, smpa2, suma, a2prob, x0
      real*8 s1  ! /1.1185e-4/  ! (Z^0.333/183)^2 ; Z=7.25
      real*8 logs1   !  /-9.0983/  ! ln(s1)

!      real*8 a1, a2, vm, logvm, suma, a2prob
!      parameter (vm = 1.d-4, logvm = -4*2.3025851, 
!     *           a1 = 0.5, a2 = -4.d0*logvm/3.,
!     *           suma = a1 + a2, a2prob = a2/suma)
!
      do while (.true.)
         call rndc(u)
         if(u .lt. a2prob) then
!             sampling by  1/v
            call rndc(u)
            v = exp(u * logvm)  ! fractional energy
!            rejection
            call rndc(u)
            if(u .lt. 1.-v) then
!                further check; get s
               s = cmigdsforbrem(ee, rho, v)
               if(u .lt. (1.-v) * cmigdgzai(s)*cmigdpsi(s)) then
!                  accepted
                  eg =min( v, 1.d0- masele/ee) * ee
                  goto 100
               endif
            endif
         else
!              sampling by 2vdv
            call rndc(u)
            call rndc(v)
            v = max(u, v)
!            rejection, get s
            s = cmigdsforbrem(ee, rho, v)
            temp = cmigdgzai(s) * (cmigdg(s) + 2* cmigdpsi(s))/3.
            call rndc(u)
            if(u .lt. temp) then
!               accepted
               eg = min(v, 1.d0-masele/ee) * ee
               goto 100
            endif
         endif
      enddo
 100  continue

Here is the call graph for this function:

Here is the caller graph for this function: