]> git.uio.no Git - u/mrichter/AliRoot.git/blame - GEANT321/fluka/evdeex.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Remove AliTRDconst.h
[u/mrichter/AliRoot.git] / GEANT321 / fluka / evdeex.F
CommitLineData
fe4da5cc 1*
2* $Id$
3*
4* $Log$
5* Revision 1.1.1.1 1995/10/24 10:19:55 cernlib
6* Geant
7*
8*
9#include "geant321/pilot.h"
10*CMZ : 3.21/02 29/03/94 15.41.42 by S.Giani
11*-- Author :
12*$ CREATE EVDEEX.FOR
13*COPY EVDEEX
14*
15*=== evdeex ===========================================================*
16*
17 SUBROUTINE EVDEEX ( WEE )
18
19#include "geant321/dblprc.inc"
20#include "geant321/dimpar.inc"
21#include "geant321/iounit.inc"
22*
23*----------------------------------------------------------------------*
24* *
25* Evdeex : created on 5-10-1990 by A. Ferrari & P. Sala, INFN Milan *
26* *
27* Last change on 28-jan-93 by Alfredo Ferrari, INFN-Milan *
28* *
29* This routine provides a simple model for sampling nuclear deexci- *
30* tation gammas following the evaporation step *
31* *
32*----------------------------------------------------------------------*
33*
34#include "geant321/balanc.inc"
35#include "geant321/eva0.inc"
36#include "geant321/finuc.inc"
37#include "geant321/nucdat.inc"
38#include "geant321/parevt.inc"
39#include "geant321/resnuc.inc"
40*
41*----------------------------------------------------------------------*
42* Entering the routine we have: *
43* Ammres is the atomic mass of the residual nucleus *
44* Ibres (Anow) its mass number *
45* Icres (Znow) its atomic number *
46* Eres its total energy *
47* Ptres its momentum *
48* Pxres x-component of the momentum *
49* Pyres y-component of the momentum *
50* Pzres z-component of the momentum *
51* Tvrecl kinetic energy *
52* Tvcms excitation energy *
53*----------------------------------------------------------------------*
54*
55 PARAMETER ( C0M1E1 = 0.306 D+00 )
56 PARAMETER ( C0E2E1 = 7.1 D-01 )
57 PARAMETER ( HNDFE1 = 5.0 D-03 )
58 PARAMETER ( HNDFM1 = 5.0 D-02 )
59 PARAMETER ( HNDFE2 = 1.0 D+01 )
60*
61 DIMENSION COSGAM (3)
62 REAL RNDM(2)
63*
64 IDEEXG = 0
65 IF ( IBRES .LE. 0 .OR. TVCMS .LE. GAMMIN ) THEN
66 TVCMS = 0.D+00
67 RETURN
68 END IF
69 ECHCK = ERES - AMMRES
70 ECHCK0 = ECHCK
71 IBHELP = IBRES / 2
72 IF ( IBHELP * 2 .LT. IBRES ) THEN
73 CALL EEXLVL ( IBRES, ICRES, DELTA, EEX2ND, EEXDUM )
74 IPAR = 1
75 JPAR = 1
76 ELSE
77 ICHELP = ICRES / 2
78 JPAR = 0
79 IF ( ICHELP * 2 .LT. ICRES ) THEN
80 CALL EEXLVL ( IBRES, ICRES, DELTA, EEX1ST, EEXDUM )
81 IPAR = 1
82 ELSE
83 CALL EEXLVL ( IBRES, ICRES, EEX1ST, DELTA, EEXDUM )
84 IPAR = 2
85 END IF
86 END IF
87 DELPAI = MIN ( EEXDUM, DELTA )
88 RNUCL = R0NUCL * RMASS (IBRES)
89 AINERM = 0.24D+00 * AMMRES * RNUCL * RNUCL
90 ROTEN0 = 0.5D+00 * PLABRC * PLABRC / AINERM
91 ENMIN = MAX ( DELTA, TWOTWO * ROTEN0 )
92 RNMASS = AMMRES + TVCMS
93 UMO = RNMASS
94 GAMCM = ERES / RNMASS
95 ETAX = PXRES / RNMASS
96 ETAY = PYRES / RNMASS
97 ETAZ = PZRES / RNMASS
98 IF ( TVCMS .LE. ENMIN .OR. IBRES .LE. 4 ) GO TO 3000
99 AHELP = HNDFE1 * RMASS (IBRES) * RMASS (IBRES)
100 CFM1E1 = C0M1E1 * HNDFM1 / AHELP
101 CFE2E1 = C0E2E1 * HNDFE2 / AHELP
102 1000 CONTINUE
103 UMEV = 1.D+03 * TVCMS
104 ASMALL = GETA ( UMEV , ICRES , IBRES-ICRES, ILVMOD, ISDUM,
105 & AOGMAX, AOGMIN )
106 TEMPSQ = 1.D-03 * ( TVCMS - DELPAI ) / ASMALL
107 TEMPER = SQRT ( TEMPSQ )
108 HHH = TVCMS / TEMPER
109 HHHSQ = HHH * HHH
110 IF(HHH.GT.88) THEN
111 AHELP=0.0
112 ELSE
113 AHELP = EXP ( - HHH )
114 ENDIF
115 BRE1M1 = 6.D+00 - AHELP * ( HHHSQ * HHH + 3.D+00 * HHHSQ
116 & + 6.D+00 * HHH )
117 BRE2 = 20.D+00 * BRE1M1 - AHELP * ( HHHSQ * HHHSQ * HHH
118 & + 5.D+00 * HHHSQ * HHHSQ )
119 BRE1M1 = BRE1M1 * ( 1.D+00 + CFM1E1 )
120 BRE2 = BRE2 * TEMPSQ * CFE2E1
121 CALL GRNDM(RNDM,1)
122 IF ( RNDM (1) .LT. BRE1M1 / ( BRE1M1 + BRE2 ) ) THEN
123 LMULT = 1
124 ELSE
125 LMULT = 2
126 END IF
127 LEXPN = 2 * LMULT + 1
128 DDLEXP = LEXPN
129 XDISMX = MIN ( DDLEXP, HHH )
130 IF(XDISMX.GT.88) THEN
131 DISMX=0.
132 ELSE
133 DISMX = XDISMX ** LEXPN * EXP ( - XDISMX )
134 ENDIF
135 2000 CONTINUE
136 CALL GRNDM(RNDM,2)
137 XTENT = RNDM (1) * HHH
138 IF(XTENT.GT.88) THEN
139 FREJE = 0.
140 ELSE
141 FREJE = XTENT ** LEXPN * EXP ( - XTENT ) / DISMX
142 ENDIF
143 IF ( RNDM (2) .GE. FREJE ) GO TO 2000
144 ENERG0 = XTENT * TEMPER
145 TVCMS = TVCMS - ENERG0
146 RNMASS = AMMRES + TVCMS
147 CALL RACO ( COSGAM (1), COSGAM (2), COSGAM (3) )
148 ERNCM = 0.5D+00 * ( UMO * UMO + RNMASS * RNMASS ) / UMO
149 EEGCM = UMO - ERNCM
150 PCMS = EEGCM
151 PCMSX = PCMS * COSGAM (1)
152 PCMSY = PCMS * COSGAM (2)
153 PCMSZ = PCMS * COSGAM (3)
154 ETAPCM = PCMSX * ETAX + PCMSY * ETAY + PCMSZ * ETAZ
155 ENERG = GAMCM * EEGCM + ETAPCM
156 PHELP = ETAPCM / ( GAMCM + 1.D+00 ) + EEGCM
157 PLBGX = PCMSX + ETAX * PHELP
158 PLBGY = PCMSY + ETAY * PHELP
159 PLBGZ = PCMSZ + ETAZ * PHELP
160 ERES = GAMCM * ERNCM - ETAPCM
161 EKRES = ERES - RNMASS
162 PHELP = - ETAPCM / ( GAMCM + 1.D+00 ) + ERNCM
163 PXRES = - PCMSX + ETAX * PHELP
164 PYRES = - PCMSY + ETAY * PHELP
165 PZRES = - PCMSZ + ETAZ * PHELP
166 ECHCK = ECHCK - ENERG
167 IDEEXG = IDEEXG + 1
168 NP = NP + 1
169 WEI (NP) = WEE
170 KPART (NP) = 7
171 TKI (NP) = ENERG
172 PLR (NP) = ENERG
173 CXR (NP) = PLBGX / ENERG
174 CYR (NP) = PLBGY / ENERG
175 CZR (NP) = PLBGZ / ENERG
176 GAMCM = ERES / RNMASS
177 ETAX = PXRES / RNMASS
178 ETAY = PYRES / RNMASS
179 ETAZ = PZRES / RNMASS
180 UMO = RNMASS
181 IF ( TVCMS .GT. ENMIN ) GO TO 1000
182 IF ( NP .GE. MXP ) THEN
183 WRITE ( LUNOUT, * )' **** Finuc overflow in Evdeex,',
184 & ' program stopped ****'
185 WRITE ( LUNERR, * )' **** Finuc overflow in Evdeex,',
186 & ' program stopped ****'
187 STOP
188 END IF
189 3000 CONTINUE
190 IF ( TVCMS .LE. ( 6 + 2 * JPAR ) * ROTEN0 ) THEN
191 ENERG0 = TVCMS
192 TVCMS = 0.D+00
193 CALL RACO ( COSGAM (1), COSGAM (2), COSGAM (3) )
194 ERNCM = 0.5D+00 * ( UMO * UMO + AMMRES * AMMRES ) / UMO
195 EEGCM = UMO - ERNCM
196 PCMS = EEGCM
197 PCMSX = PCMS * COSGAM (1)
198 PCMSY = PCMS * COSGAM (2)
199 PCMSZ = PCMS * COSGAM (3)
200 ETAPCM = PCMSX * ETAX + PCMSY * ETAY + PCMSZ * ETAZ
201 ENERG = GAMCM * EEGCM + ETAPCM
202 PHELP = ETAPCM / ( GAMCM + 1.D+00 ) + EEGCM
203 PLBGX = PCMSX + ETAX * PHELP
204 PLBGY = PCMSY + ETAY * PHELP
205 PLBGZ = PCMSZ + ETAZ * PHELP
206 ERES = GAMCM * ERNCM - ETAPCM
207 EKRES = ERES - AMMRES
208 PHELP = - ETAPCM / ( GAMCM + 1.D+00 ) + ERNCM
209 PXRES = - PCMSX + ETAX * PHELP
210 PYRES = - PCMSY + ETAY * PHELP
211 PZRES = - PCMSZ + ETAZ * PHELP
212 ECHCK = ECHCK - ENERG
213 IDEEXG = IDEEXG + 1
214 NP = NP + 1
215 WEI (NP) = WEE
216 KPART (NP) = 7
217 TKI (NP) = ENERG
218 PLR (NP) = ENERG
219 CXR (NP) = PLBGX / ENERG
220 CYR (NP) = PLBGY / ENERG
221 CZR (NP) = PLBGZ / ENERG
222 ELSE
223 AIAMOM = SQRT ( 0.25D+00 + TVCMS / ROTEN0
224 & + 0.75D+00 * JPAR ) - 0.25D+00
225 IF ( IPAR .EQ. 1 ) THEN
226 IF ( JPAR .EQ. 1 ) THEN
227 JAMIN = 2 * INT ( AIAMOM ) + 1
228 IAMIN = JAMIN / 2
229 IF ( MOD ( IAMIN, 2 ) .GT. 0 ) THEN
230 EGROUN = 3.75D+00 * ROTEN0
231 ELSE
232 EGROUN = 0.75D+00 * ROTEN0
233 END IF
234 ELSE
235 IAMIN = INT ( AIAMOM )
236 JAMIN = 2 * IAMIN
237 IF ( MOD ( IAMIN, 2 ) .GT. 0 ) THEN
238 EGROUN = 2.D+00 * ROTEN0
239 ELSE
240 EGROUN = 0.D+00
241 END IF
242 END IF
243 ELSE
244 IAMIN = INT ( AIAMOM ) / IPAR
245 IAMIN = IAMIN * IPAR
246 JAMIN = 2 * IAMIN
247 EGROUN = 0.D+00
248 END IF
249 DELTAE = TVCMS + EGROUN - 0.25D+00 * ROTEN0 * JAMIN *
250 & ( JAMIN + 2 )
251 DO 4000 JAMOM = JAMIN, 4, -4
252 ENERG0 = ROTEN0 * ( 2 * JAMOM - 2 )
253 & + DELTAE
254 DELTAE = 0.D+00
255 TVCMS = TVCMS - ENERG0
256 RNMASS = AMMRES + TVCMS
257 CALL RACO ( COSGAM (1), COSGAM (2), COSGAM (3) )
258 ERNCM = 0.5D+00 * ( UMO * UMO + RNMASS * RNMASS ) / UMO
259 EEGCM = UMO - ERNCM
260 PCMS = EEGCM
261 PCMSX = PCMS * COSGAM (1)
262 PCMSY = PCMS * COSGAM (2)
263 PCMSZ = PCMS * COSGAM (3)
264 ETAPCM = PCMSX * ETAX + PCMSY * ETAY + PCMSZ * ETAZ
265 ENERG = GAMCM * EEGCM + ETAPCM
266 PHELP = ETAPCM / ( GAMCM + 1.D+00 ) + EEGCM
267 PLBGX = PCMSX + ETAX * PHELP
268 PLBGY = PCMSY + ETAY * PHELP
269 PLBGZ = PCMSZ + ETAZ * PHELP
270 ERES = GAMCM * ERNCM - ETAPCM
271 EKRES = ERES - RNMASS
272 PHELP = - ETAPCM / ( GAMCM + 1.D+00 ) + ERNCM
273 PXRES = - PCMSX + ETAX * PHELP
274 PYRES = - PCMSY + ETAY * PHELP
275 PZRES = - PCMSZ + ETAZ * PHELP
276 ECHCK = ECHCK - ENERG
277 IDEEXG = IDEEXG + 1
278 NP = NP + 1
279 WEI (NP) = WEE
280 KPART (NP) = 7
281 TKI (NP) = ENERG
282 PLR (NP) = ENERG
283 CXR (NP) = PLBGX / ENERG
284 CYR (NP) = PLBGY / ENERG
285 CZR (NP) = PLBGZ / ENERG
286 GAMCM = ERES / RNMASS
287 ETAX = PXRES / RNMASS
288 ETAY = PYRES / RNMASS
289 ETAZ = PZRES / RNMASS
290 UMO = RNMASS
291 4000 CONTINUE
292 END IF
293 ECHCK = ECHCK - EKRES
294 IF ( ABS ( ECHCK ) .GT. 1.D-7 * ECHCK0 ) THEN
295 WRITE ( LUNERR, * )' **** No energy conservation in Evdeex',
296 & ' ****', ECHCK, IDEEXG
297 END IF
298 TVCMS = 0.D+00
299 IF ( NP .GT. MXP ) THEN
300 WRITE ( LUNOUT, * )' **** Finuc overflow in Evdeex,',
301 & ' program stopped ****'
302 WRITE ( LUNERR, * )' **** Finuc overflow in Evdeex,',
303 & ' program stopped ****'
304 STOP
305 END IF
306 TVRECL = EKRES
307 P2RES = PXRES * PXRES + PYRES * PYRES + PZRES * PZRES
308 PTRES = SQRT (P2RES)
309*=== End of subroutine Evdeex =========================================*
310 RETURN
311 END