GEANT321/gphys/grayl.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.2  1996/11/06 12:53:21  ravndal
   6 * Correct the rotation into the GEANT system
   7 *
   8 * Revision 1.1.1.1  1995/10/24 10:21:33  cernlib
   9 * Geant
  10 *
  11 *
  12 #include "geant321/pilot.h"
  13 *CMZ :  3.21/02 29/03/94  15.41.23  by  S.Giani
  14 *-- Author :
  15       SUBROUTINE GRAYL
  16 C.
  17 C.    ******************************************************************
  18 C.    *                                                                *
  19 C.    *  Generates COHERENT SCATTERING mechanism                       *
  20 C.    *                                                                *
  21 C.    *    ==>Called by : GTGAMA                                       *
  22 C.    *       Author    G.Tromba, P.Bregant                             *
  23 C.    *                                                                *
  24 C.    ******************************************************************
  25 C.
  26 #include "geant321/gcbank.inc"
  27 #include "geant321/gctrak.inc"
  28 #include "geant321/gcphys.inc"
  29 #include "geant321/gconsp.inc"
  30 #include "geant321/gcking.inc"
  31 #include "geant321/gcjloc.inc"
  32 #include "geant321/gcmulo.inc"
  33       DIMENSION VROT(3)
  34       DIMENSION RNDM(2)
  35       LOGICAL ROTATE
  36 C.
  37 C.    ------------------------------------------------------------------
  38 C.
  39       KCASE = NAMEC(25)
  40       NGKINE=0
  41 *
  42       IMAX=NEKBIN
  43       DO 2 I=2,NEK1
  44          IF (Q(JRAYL+I).GT.0.9*BIG) THEN
  45             IMAX = I-1
  46             GO TO 2
  47          ENDIF
  48    2  CONTINUE
  49 C
  50 C             Extract a value for the momentum transfer according
  51 C             with the FORM FACTORS distribution.
  52 C
  53       IF (IEKBIN.LT.IMAX) THEN
  54         GEKRT1=1.-GEKRAT
  55         A2MAX=GEKRT1*Q(JRAYL+NEK1+IEKBIN)+GEKRAT*Q(JRAYL+NEK1+IEKBIN+1)
  56       ELSE
  57         A2MAX=Q(JRAYL+NEK1+IEKBIN)
  58       ENDIF
  59       GEKIN1=1./GEKIN
  60    5  CALL GRNDM(RNDM,2)
  61       A2=RNDM(1)*A2MAX
  62       DO 10 I=1,(IMAX-1)
  63          IF (A2.LE.Q(JRAYL+I+NEK1)) THEN
  64             IE=I
  65             GO TO 20
  66          ENDIF
  67 10    CONTINUE
  68       IE=IMAX
  69 20    IF (IE.GT.1) THEN
  70         GA2RAT =(A2-Q(JRAYL+NEK1+IE-1))/(Q(JRAYL+NEK1+IE)
  71      +          -Q(JRAYL+NEK1+IE-1))
  72         QNEW   = ELOW(IE-1)*(1.-GA2RAT)+ELOW(IE)*GA2RAT
  73       ELSE
  74         GA2RAT  = A2/Q(JRAYL+NEK1+IE)
  75         QNEW    = ELOW(IE)*GA2RAT
  76       ENDIF
  77 C
  78 C            Find out cos(theta)
  79       COSTH=1.-0.5*(QNEW*GEKIN1)**2
  80       IF(RNDM(2).GT.0.5*(1.+COSTH**2)) GO TO 5
  81       SINTH=SQRT((1.-COSTH)*(1.+COSTH))
  82 C
  83 C            Extract PH  (PHI unif. between 0 and 2PI)
  84       CALL GRNDM(RNDM,2)
  85       PH=RNDM(1)*TWOPI
  86       COSPH=COS(PH)
  87       SINPH=SIN(PH)
  88 C
  89 C            Set new VECT(I), i=4,5,6
  90       VROT(1)=GEKIN*SINTH*COSPH
  91       VROT(2)=GEKIN*SINTH*SINPH
  92       VROT(3)=GEKIN*COSTH
  93 C
  94 C             Rotate scattered photon into GEANT system
  95 C
  96       CALL GFANG(VECT(4),COSTH,SINTH,COSPH,SINPH,ROTATE)
  97       IF (ROTATE) CALL GDROT(VROT,COSTH,SINTH,COSPH,SINPH)
  98
  99       VMROT=1./SQRT(VROT(1)*VROT(1)+VROT(2)*VROT(2)+VROT(3)*VROT(3))
 100       VECT(4)=VROT(1)*VMROT
 101       VECT(5)=VROT(2)*VMROT
 102       VECT(6)=VROT(3)*VMROT
 103 C
 104 C            Update probabilities
 105 C
 106 *
 107       ZINTRA=-LOG(RNDM(2))
 108       SLRAYL=SLENG
 109       STEPRA=BIG
 110 *
 111 C
 112       END
 113