This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / gphys / gpairm.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:21:28  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.22  by  S.Giani
11 *-- Author :
12       SUBROUTINE GPAIRM
13 C.
14 C.    ******************************************************************
15 C.    *                                                                *
16 C.    *  Simulates direct pair production by muons                     *
17 C.    *                                                                *
18 C.    *    ==>Called by : GTMUON                                       *
19 C.    *       Author    L.Urban  *********                             *
20 C.    *                                                                *
21 C.    ******************************************************************
22 C.
23 #include "geant321/gcbank.inc"
24 #include "geant321/gconsp.inc"
25 #include "geant321/gctrak.inc"
26 #include "geant321/gcmate.inc"
27 #include "geant321/gckine.inc"
28 #include "geant321/gcking.inc"
29 #include "geant321/gcphys.inc"
30 #include "geant321/gccuts.inc"
31       DIMENSION RNDM(2)
32       LOGICAL ROTATE
33 C
34 C         CVM= 3*SQRT(e)*EMMU/4
35 C         EM6= 6*EMMU**2
36       DATA CVM,EM6/0.130652,0.066983/
37       DATA AL10T/9.212/
38 C.
39 C.    ------------------------------------------------------------------
40 C.
41       IF(GEKIN.LE.PPCUTM)GO TO 900
42       EEM1=GETOT
43       KCASE=NAMEC(6)
44 C
45       VMIN=4*EMASS/EEM1
46       VMAX=1.-CVM*Z**0.333333/EEM1
47       IF(VMAX.LE.VMIN)GO TO 900
48       VC  = PPCUTM/EEM1
49       ALE=LOG(EEM1)
50       ALFA=1.+ALE/AL10T
51       V0=0.18*(4.+ALE/AL10T)*ALFA*(ALFA*VMIN)**0.6666667
52       BETA=0.1*(1.+3.*ALE/AL10T)
53       B=0.9/(1.+0.4*ALE+0.022*ALE*ALE)
54       AA=1.+2.*B*LOG(VC/V0)
55       IF(AA.LE.1.) AA=1.05
56       A1=1.-AA
57       CC=EXP(-0.25*A1*A1/B)
58       A1R=1./A1
59       C1=VMAX**A1
60       C2=VC**A1
61 C
62 C     SAMPLE V AND RO
63 C
64   50  CALL GRNDM(RNDM,2)
65       R=RNDM(1)
66       V=(R*C1+(1.-R)*C2)**A1R
67       IF(V.LE.VMIN) GOTO 50
68       IF(V.LT.V0) THEN
69         SCREJ=CC*((V-VMIN)/(V0-VMIN))**BETA*(V0/V)**A1
70       ELSE
71         SCREJ=CC*(V0/V)**(A1+B*LOG(V/V0))
72       ENDIF
73       IF(RNDM(2).GT.SCREJ) GOTO 50
74       R0MAX= SCREJ*(1.-EM6/(EEM1**2*(1.-V)))
75       CALL GRNDM(RNDM,2)
76       R0   = R0MAX*(2.*RNDM(1)-1.)
77 C
78 C           Energies
79 C
80       EPP  = V*EEM1
81       IF(IPAIR.NE.1)THEN
82          NGKINE=0
83          DESTEP=DESTEP+EPP
84          GO TO 60
85       ENDIF
86       EPOS = 0.5*EPP*(1.+R0)
87       EMIN = EPP-EPOS
88 C
89 C           Angles
90 C
91       THETA = AMASS/EEM1
92       SINTH = SIN(THETA)
93       COSTH = COS(THETA)
94       PHI   = TWOPI*RNDM(2)
95       COSPHI= COS(PHI)
96       SINPHI= SIN(PHI)
97 C
98       CALL GFANG(VECT(4),COSAL,SINAL,COSBT,SINBT,ROTATE)
99 C
100 C           Positron
101 C
102       NGKINE = 0
103       TPOS   = EPOS-EMASS
104       IF(TPOS.GT.CUTELE)THEN
105          PPOS  = SQRT((EPOS+EMASS)*TPOS)
106          NGKINE= NGKINE+1
107          GKIN(1,NGKINE)=PPOS*SINTH*COSPHI
108          GKIN(2,NGKINE)=PPOS*SINTH*SINPHI
109          GKIN(3,NGKINE)=PPOS*COSTH
110          GKIN(4,NGKINE)=EPOS
111          GKIN(5,NGKINE)=2.
112          TOFD(NGKINE)=0.
113          GPOS(1,NGKINE) = VECT(1)
114          GPOS(2,NGKINE) = VECT(2)
115          GPOS(3,NGKINE) = VECT(3)
116          IF(ROTATE)
117      +   CALL GDROT(GKIN(1,NGKINE),COSAL,SINAL,COSBT,SINBT)
118       ELSE
119          DESTEP=DESTEP+TPOS
120          CALL GANNI2
121       ENDIF
122 C
123 C           Electron
124 C
125       TMIN=EMIN-EMASS
126       IF(TMIN.GT.CUTELE)THEN
127          PMIN  = SQRT((EMIN+EMASS)*TMIN)
128          NGKINE= NGKINE+1
129          GKIN(1,NGKINE)=-PMIN*SINTH*COSPHI
130          GKIN(2,NGKINE)=-PMIN*SINTH*SINPHI
131          GKIN(3,NGKINE)=PMIN*COSTH
132          GKIN(4,NGKINE)=EMIN
133          GKIN(5,NGKINE)=3.
134          TOFD(NGKINE)=0.
135          GPOS(1,NGKINE) = VECT(1)
136          GPOS(2,NGKINE) = VECT(2)
137          GPOS(3,NGKINE) = VECT(3)
138          IF(ROTATE)
139      +   CALL GDROT(GKIN(1,NGKINE),COSAL,SINAL,COSBT,SINBT)
140       ELSE
141          DESTEP=DESTEP+TMIN
142       ENDIF
143 C
144 C           Correct muon track
145 C
146   60  GEKIN  = GEKIN-EPP
147       GETOT  = GEKIN+AMASS
148       VECT(7)= SQRT((GETOT+AMASS)*GEKIN)
149       CALL GEKBIN
150 C
151 C           Update probabilities
152 C
153  900  CALL GRNDM(RNDM,1)
154       ZINTPA = -LOG(RNDM(1))
155       SLPAIR = SLENG
156       STEPPA = BIG
157 C
158       END