GEANT321/gphys/gbremm.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1995/10/24 10:21:22  cernlib
   6 * Geant
   7 *
   8 *
   9 #include "geant321/pilot.h"
  10 *CMZ :  3.21/02 29/03/94  15.41.21  by  S.Giani
  11 *-- Author :
  12       SUBROUTINE GBREMM
  13 C.
  14 C.    ******************************************************************
  15 C.    *                                                                *
  16 C.    *  Simulates discrete hard BREMSSTRAHLUNG by muons.              *
  17 C.    *                                                                *
  18 C.    *                                                                *
  19 C.    *    ==>Called by : GTMUON                                       *
  20 C.    *       Author     L.Urban  *********                            *
  21 C.    *                                                                *
  22 C.    ******************************************************************
  23 C.
  24 #include "geant321/gcbank.inc"
  25 #include "geant321/gcphys.inc"
  26 #include "geant321/gconsp.inc"
  27 #include "geant321/gctrak.inc"
  28 #include "geant321/gcjloc.inc"
  29 #include "geant321/gckine.inc"
  30 #include "geant321/gcking.inc"
  31 #include "geant321/gccuts.inc"
  32       DIMENSION RNDM(2)
  33 C.
  34 C.    ------------------------------------------------------------------
  35 C.
  36 C             Ensure cut-off avoids infra-red catastrophe.
  37 C
  38       EEM1   = GETOT
  39       IF (GEKIN.LE.BCUTM) GO TO 900
  40       KCASE  = NAMEC(9)
  41 C
  42       C1   = Q(JPROB+31)
  43       C2   = Q(JPROB+32)
  44       C3   = Q(JPROB+33)
  45       VC   = BCUTM/EEM1
  46       VM   = 1.-C3/EEM1
  47       IF(VM.LE.0.)GO TO 900
  48       CREJ = C2/EEM1
  49 C
  50   50  CALL GRNDM(RNDM,2)
  51       R    = RNDM(1)
  52       V    = VC*(VM/VC)**R
  53 C
  54 C             Compute rejection function
  55 C
  56       V1   = 1.-V
  57       F1   = C1-LOG(1.+CREJ*V/V1)
  58       SCREJ= (V1+0.75*V*V)*F1/C1
  59       IF(RNDM(2).GT.SCREJ)GO TO 50
  60 C
  61 C             Photon energy
  62 C
  63       EGAMMA=EEM1*V
  64 C
  65 C        CUT ON ENERGY THRESHOLD ?
  66 C
  67       IF((IBREM.NE.1).OR.(EGAMMA.LE.CUTGAM)) THEN
  68           NGKINE = 0
  69           DESTEP = DESTEP + EGAMMA
  70           GO TO 800
  71       ENDIF
  72 C
  73 C             Generate emitted photon angles with respect to a Z-axis
  74 C             defined along parent track. PHI is generated isotropically
  75 C             and THETA is assigned a universal angular distribution
  76 C
  77       THETA  = GBTETH(EEM1, AMASS, V)*AMASS/EEM1
  78       SINTH  = SIN(THETA)
  79       COSTH  = COS(THETA)
  80       CALL GRNDM(RNDM,1)
  81       PHI    = TWOPI*RNDM(1)
  82       COSPHI = COS(PHI)
  83       SINPHI = SIN(PHI)
  84 C
  85 C             Polar co-ordinates to momentum components.
  86 C
  87       NGKINE = 1
  88       GKIN(1,1)=EGAMMA*SINTH*COSPHI
  89       GKIN(2,1)=EGAMMA*SINTH*SINPHI
  90       GKIN(3,1)=EGAMMA*COSTH
  91       GKIN(4,1)=EGAMMA
  92       GKIN(5,1)=1.
  93       TOFD(1)  =0.
  94       GPOS(1,1) = VECT(1)
  95       GPOS(2,1) = VECT(2)
  96       GPOS(3,1) = VECT(3)
  97 C
  98 C             Rotate photon into GEANT system
  99 C
 100       CALL GVROT(VECT(4),GKIN)
 101 C
 102 C             Correct track for lost energy
 103 C
 104   800 CONTINUE
 105       GEKIN = GEKIN - EGAMMA
 106       GETOT = GEKIN + AMASS
 107       VECT(7)=SQRT (ABS((GETOT+AMASS)*GEKIN))
 108       CALL GEKBIN
 109 C
 110 C             Update probabilities
 111 C
 112   900 CALL GRNDM(RNDM,1)
 113       ZINTBR=-LOG(RNDM(1))
 114       SLBREM=SLENG
 115       STEPBR=BIG
 116 C
 117       END