This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / fluka / pmprab.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:19:58  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.44  by  S.Giani
11 *-- Author :
12 *$ CREATE PMPRAB.FOR
13 *COPY PMPRAB
14 *
15 *=== pmprab ===========================================================*
16 *
17       SUBROUTINE PMPRAB ( KPROJ, EKIN, PPROJ, TXX, TYY, TZZ, WEE )
18  
19 #include "geant321/dblprc.inc"
20 #include "geant321/dimpar.inc"
21 #include "geant321/iounit.inc"
22 *
23 *----------------------------------------------------------------------*
24 *                                                                      *
25 *     Created on 22 september 1991 by    Alfredo Ferrari & Paola Sala  *
26 *                                                   Infn - Milan       *
27 *                                                                      *
28 *     Last change on 22-sep-91     by    Alfredo Ferrari               *
29 *                                                                      *
30 *                                                                      *
31 *----------------------------------------------------------------------*
32 *
33 #include "geant321/balanc.inc"
34 #include "geant321/finuc.inc"
35 #include "geant321/nucdat.inc"
36 #include "geant321/nucgeo.inc"
37 #include "geant321/parevt.inc"
38 #include "geant321/parnuc.inc"
39 #include "geant321/part.inc"
40 #include "geant321/resnuc.inc"
41 *
42       REAL RNDM(1)
43 *
44       IF ( KPROJ .NE. 14 .OR. EKIN .GT. 2.D+00 * GAMMIN .OR. IBTAR .NE.
45      &     1 .OR. ICHTAR .NE. 1 ) THEN
46          WRITE (LUNOUT,*)' **** Pmprab: kproj,ekin,ibtar,ichtar',
47      &                     KPROJ,EKIN,IBTAR,ICHTAR
48          WRITE (LUNERR,*)' **** Pmprab: kproj,ekin,ibtar,ichtar',
49      &                     KPROJ,EKIN,IBTAR,ICHTAR
50       END IF
51       PXRES  = PXTTOT
52       PYRES  = PYTTOT
53       PZRES  = PZTTOT
54       PTRES  = PTTOT
55       CALL GRNDM(RNDM,1)
56       RNDPAN = RNDM (1)
57       IF ( RNDPAN .GE. 1.D+00 / PNFRAT ) THEN
58          ERES   = EKIN + AM (KPROJ) + EKFERM + AM (1)
59          UMO2   = ( ERES - PTRES ) * ( ERES + PTRES )
60          UMO    = SQRT (UMO2)
61          GAMCM = ERES  / UMO
62          ETAX  = PXRES / UMO
63          ETAY  = PYRES / UMO
64          ETAZ  = PZRES / UMO
65          ECMSNU = 0.5D+00 * ( UMO2 + AMNUSQ (2) ) / UMO
66          PCMS   = UMO - ECMSNU
67          CALL RACO ( PCMSX, PCMSY, PCMSZ )
68          PCMSX = PCMS * PCMSX
69          PCMSY = PCMS * PCMSY
70          PCMSZ = PCMS * PCMSZ
71          ETAPCM = ETAX * PCMSX + ETAY * PCMSY + ETAZ * PCMSZ
72          PHELP  = PCMS + ETAPCM / ( GAMCM + 1.D+00 )
73          ENPHOT = GAMCM * PCMS + ETAPCM
74          PXHELP = PCMSX + ETAX * PHELP
75          PYHELP = PCMSY + ETAY * PHELP
76          PZHELP = PCMSZ + ETAZ * PHELP
77          PXRES = PXRES - PXHELP
78          PYRES = PYRES - PYHELP
79          PZRES = PZRES - PZHELP
80          ERES  = ERES  - ENPHOT
81          NP = NP + 1
82          TKI   (NP) = ENPHOT
83          KPART (NP) = 7
84          PLR   (NP) = ENPHOT
85          CXR   (NP) = PXHELP / PLR (NP)
86          CYR   (NP) = PYHELP / PLR (NP)
87          CZR   (NP) = PZHELP / PLR (NP)
88          WEI   (NP) = WEE
89          IOTHER = IOTHER + 1
90          PXNUCR = PXNUCR + PXHELP
91          PYNUCR = PYNUCR + PYHELP
92          PZNUCR = PZNUCR + PZHELP
93          ENUCR  = ENUCR  + TKI (NP)
94          IBNUCR = IBNUCR + IBAR (KPART(NP))
95          ICNUCR = ICNUCR + ICH  (KPART(NP))
96          ETAPCM = - ETAPCM
97          PHELP  = ECMSNU + ETAPCM / ( GAMCM + 1.D+00 )
98          ENNEU  = GAMCM * ECMSNU + ETAPCM
99          PXHELP = -PCMSX + ETAX * PHELP
100          PYHELP = -PCMSY + ETAY * PHELP
101          PZHELP = -PCMSZ + ETAZ * PHELP
102          NP = NP + 1
103          TKI   (NP) = ENNEU - AM (8)
104          KPART (NP) = 8
105          PLR (NP) = SQRT ( PXHELP**2 + PYHELP**2 + PZHELP**2 )
106          CXR   (NP) = PXHELP / PLR (NP)
107          CYR   (NP) = PYHELP / PLR (NP)
108          CZR   (NP) = PZHELP / PLR (NP)
109          WEI   (NP) = WEE
110          ERES  = ERES  - ENNEU
111          PXRES = PXRES - PXHELP
112          PYRES = PYRES - PYHELP
113          PZRES = PZRES - PZHELP
114          IBRES = 0
115          ICRES = 0
116          PTRES = 0.D+00
117          ANOW  = 0.D+00
118          ZNOW  = 0.D+00
119       ELSE
120          ERES   = EKIN + AM (KPROJ) + EKFERM + AM (1)
121          UMO2   = ( ERES - PTRES ) * ( ERES + PTRES )
122          UMO    = SQRT (UMO2)
123          GAMCM = ERES  / UMO
124          ETAX  = PXRES / UMO
125          ETAY  = PYRES / UMO
126          ETAZ  = PZRES / UMO
127          ECMSNU = 0.5D+00 * ( UMO2 + AM (8)* AM (8) - AM (23) * AM (23)
128      &          ) / UMO
129          ECMSP0 = UMO - ECMSNU
130          PCMS = SQRT ( ( ECMSP0 - AM (23) ) * ( ECMSP0 + AM (23) ) )
131          CALL RACO ( PCMSX, PCMSY, PCMSZ )
132          PCMSX = PCMS * PCMSX
133          PCMSY = PCMS * PCMSY
134          PCMSZ = PCMS * PCMSZ
135          ETAPCM = ETAX * PCMSX + ETAY * PCMSY + ETAZ * PCMSZ
136          PHELP  = ECMSP0 + ETAPCM / ( GAMCM + 1.D+00 )
137          ENPIO0 = GAMCM * ECMSP0 + ETAPCM
138          PXHELP = PCMSX + ETAX * PHELP
139          PYHELP = PCMSY + ETAY * PHELP
140          PZHELP = PCMSZ + ETAZ * PHELP
141          PXRES = PXRES - PXHELP
142          PYRES = PYRES - PYHELP
143          PZRES = PZRES - PZHELP
144          ERES  = ERES  - ENPIO0
145          NP = NP + 1
146          TKI   (NP) = ENPIO0 - AM (23)
147          KPART (NP) = 23
148          PLR (NP) = SQRT ( PXHELP**2 + PYHELP**2 + PZHELP**2 )
149          CXR   (NP) = PXHELP / PLR (NP)
150          CYR   (NP) = PYHELP / PLR (NP)
151          CZR   (NP) = PZHELP / PLR (NP)
152          WEI   (NP) = WEE
153          IOTHER = IOTHER + 1
154          PXNUCR = PXNUCR + PXHELP
155          PYNUCR = PYNUCR + PYHELP
156          PZNUCR = PZNUCR + PZHELP
157          ENUCR  = ENUCR  + TKI (NP)
158          IBNUCR = IBNUCR + IBAR (KPART(NP))
159          ICNUCR = ICNUCR + ICH  (KPART(NP))
160          ETAPCM = - ETAPCM
161          PHELP  = ECMSNU + ETAPCM / ( GAMCM + 1.D+00 )
162          ENNEU  = GAMCM * ECMSNU + ETAPCM
163          PXHELP = -PCMSX + ETAX * PHELP
164          PYHELP = -PCMSY + ETAY * PHELP
165          PZHELP = -PCMSZ + ETAZ * PHELP
166          NP = NP + 1
167          TKI   (NP) = ENNEU - AM (8)
168          KPART (NP) = 8
169          PLR (NP) = SQRT ( PXHELP**2 + PYHELP**2 + PZHELP**2 )
170          CXR   (NP) = PXHELP / PLR (NP)
171          CYR   (NP) = PYHELP / PLR (NP)
172          CZR   (NP) = PZHELP / PLR (NP)
173          WEI   (NP) = WEE
174          ERES  = ERES  - ENNEU
175          PXRES = PXRES - PXHELP
176          PYRES = PYRES - PYHELP
177          PZRES = PZRES - PZHELP
178          IBRES = 0
179          ICRES = 0
180          PTRES = 0.D+00
181          ANOW  = 0.D+00
182          ZNOW  = 0.D+00
183       END IF
184       RETURN
185 *=== End of subroutine PMPRAB =========================================*
186       END