]> git.uio.no Git - u/mrichter/AliRoot.git/blame - GEANT321/fluka/samcst.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Default compile option changed to -g (Alpha)
[u/mrichter/AliRoot.git] / GEANT321 / fluka / samcst.F
CommitLineData
fe4da5cc 1*
2* $Id$
3*
4* $Log$
5* Revision 1.1.1.1 1995/10/24 10:20:04 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 SAMCST.FOR
13*COPY SAMCST
14*
15*=== samcst ===========================================================*
16*
17C
18C*******************************************************************
19C
20 SUBROUTINE SAMCST(KPROJ,EKIN,CST)
21
22#include "geant321/dblprc.inc"
23#include "geant321/dimpar.inc"
24#include "geant321/iounit.inc"
25C
26C HJM 24/10/88
27C
28C SAMPLING OF COS(THETA)
29C FOR NUCLEON-PROTON ELASTIC SCATTERING
30C ACCORDING TO HETKFA2/BERTINI PARAMETRIZATION
31C
32C-------------------------------------------------------------------
33C
34 DIMENSION DCLIN(195),DCHN(143),DCHNA(36),DCHNB(60)
35 DIMENSION PDCI(60),PDCH(55)
36 REAL RNDM(4)
37 SAVE DCLIN, DCHN, DCHNA, DCHNB, PDCI, PDCH
38C
39 DATA (DCLIN(I),I=1,80) /
40C*** DCLN ARRAY
41 * 5.000D-01, 1.000D+00, 0.000D+00, 1.000D+00, 0.000D+00,
42 * 4.993D-01, 9.881D-01, 5.963D-02, 9.851D-01, 5.945D-02,
43 * 4.936D-01, 8.955D-01, 5.224D-01, 8.727D-01, 5.091D-01,
44 * 4.889D-01, 8.228D-01, 8.859D-01, 7.871D-01, 8.518D-01,
45 * 4.874D-01, 7.580D-01, 1.210D+00, 7.207D-01, 1.117D+00,
46 * 4.912D-01, 6.969D-01, 1.516D+00, 6.728D-01, 1.309D+00,
47 * 5.075D-01, 6.471D-01, 1.765D+00, 6.667D-01, 1.333D+00,
48 * 5.383D-01, 6.054D-01, 1.973D+00, 7.059D-01, 1.176D+00,
49 * 5.397D-01, 5.990D-01, 2.005D+00, 7.023D-01, 1.191D+00,
50 * 5.336D-01, 6.083D-01, 1.958D+00, 6.959D-01, 1.216D+00,
51 * 5.317D-01, 6.075D-01, 1.962D+00, 6.897D-01, 1.241D+00,
52 * 5.300D-01, 6.016D-01, 1.992D+00, 6.786D-01, 1.286D+00,
53 * 5.281D-01, 6.063D-01, 1.969D+00, 6.786D-01, 1.286D+00,
54 * 5.280D-01, 5.960D-01, 2.020D+00, 6.667D-01, 1.333D+00,
55 * 5.273D-01, 5.920D-01, 2.040D+00, 6.604D-01, 1.358D+00,
56 * 5.273D-01, 5.862D-01, 2.069D+00, 6.538D-01, 1.385D+00/
57 DATA (DCLIN(I),I=81,160) /
58C*** DCIN ARRAY
59 * 5.223D-01, 5.980D-01, 2.814D+00, 6.538D-01, 1.385D+00,
60 * 5.202D-01, 5.969D-01, 2.822D+00, 6.471D-01, 1.412D+00,
61 * 5.183D-01, 5.881D-01, 2.883D+00, 6.327D-01, 1.469D+00,
62 * 5.159D-01, 5.866D-01, 2.894D+00, 6.250D-01, 1.500D+00,
63 * 5.133D-01, 5.850D-01, 2.905D+00, 6.170D-01, 1.532D+00,
64 * 5.106D-01, 5.833D-01, 2.917D+00, 6.087D-01, 1.565D+00,
65 * 5.084D-01, 5.801D-01, 2.939D+00, 6.000D-01, 1.600D+00,
66 * 5.063D-01, 5.763D-01, 2.966D+00, 5.909D-01, 1.636D+00,
67 * 5.036D-01, 5.730D-01, 2.989D+00, 5.814D-01, 1.674D+00,
68 * 5.014D-01, 5.683D-01, 3.022D+00, 5.714D-01, 1.714D+00,
69 * 4.986D-01, 5.641D-01, 3.051D+00, 5.610D-01, 1.756D+00,
70 * 4.964D-01, 5.580D-01, 3.094D+00, 5.500D-01, 1.800D+00,
71 * 4.936D-01, 5.573D-01, 3.099D+00, 5.431D-01, 1.827D+00,
72 * 4.909D-01, 5.509D-01, 3.144D+00, 5.313D-01, 1.875D+00,
73 * 4.885D-01, 5.512D-01, 3.142D+00, 5.263D-01, 1.895D+00,
74 * 4.857D-01, 5.437D-01, 3.194D+00, 5.135D-01, 1.946D+00/
75 DATA (DCLIN(I),I=161,195) /
76 * 4.830D-01, 5.353D-01, 3.253D+00, 5.000D-01, 2.000D+00,
77 * 4.801D-01, 5.323D-01, 3.274D+00, 4.915D-01, 2.034D+00,
78 * 4.770D-01, 5.228D-01, 3.341D+00, 4.767D-01, 2.093D+00,
79 * 4.738D-01, 5.156D-01, 3.391D+00, 4.643D-01, 2.143D+00,
80 * 4.701D-01, 5.010D-01, 3.493D+00, 4.444D-01, 2.222D+00,
81 * 4.672D-01, 4.990D-01, 3.507D+00, 4.375D-01, 2.250D+00,
82 * 4.634D-01, 4.856D-01, 3.601D+00, 4.194D-01, 2.323D+00/
83C
84 DATA PDCI /
85 * 4.400D+02, 1.896D-01, 1.931D-01, 1.982D-01, 1.015D-01,
86 * 1.029D-01, 4.180D-02, 4.228D-02, 4.282D-02, 4.350D-02,
87 * 2.204D-02, 2.236D-02, 5.900D+02, 1.433D-01, 1.555D-01,
88 * 1.774D-01, 1.000D-01, 1.128D-01, 5.132D-02, 5.600D-02,
89 * 6.158D-02, 6.796D-02, 3.660D-02, 3.820D-02, 6.500D+02,
90 * 1.192D-01, 1.334D-01, 1.620D-01, 9.527D-02, 1.141D-01,
91 * 5.283D-02, 5.952D-02, 6.765D-02, 7.878D-02, 4.796D-02,
92 * 6.957D-02, 8.000D+02, 4.872D-02, 6.694D-02, 1.152D-01,
93 * 9.348D-02, 1.368D-01, 6.912D-02, 7.953D-02, 9.577D-02,
94 * 1.222D-01, 7.755D-02, 9.525D-02, 1.000D+03, 3.997D-02,
95 * 5.456D-02, 9.804D-02, 8.084D-02, 1.208D-01, 6.520D-02,
96 * 8.233D-02, 1.084D-01, 1.474D-01, 9.328D-02, 1.093D-01/
97C
98 DATA PDCH /
99 * 1.000D+03, 9.453D-02, 9.804D-02, 8.084D-02, 1.208D-01,
100 * 6.520D-02, 8.233D-02, 1.084D-01, 1.474D-01, 9.328D-02,
101 * 1.093D-01, 1.400D+03, 1.072D-01, 7.450D-02, 6.645D-02,
102 * 1.136D-01, 6.750D-02, 8.580D-02, 1.110D-01, 1.530D-01,
103 * 1.010D-01, 1.350D-01, 2.170D+03, 4.004D-02, 3.013D-02,
104 * 2.664D-02, 5.511D-02, 4.240D-02, 7.660D-02, 1.364D-01,
105 * 2.300D-01, 1.670D-01, 2.010D-01, 2.900D+03, 1.870D-02,
106 * 1.804D-02, 1.320D-02, 2.970D-02, 2.860D-02, 5.160D-02,
107 * 1.020D-01, 2.400D-01, 2.250D-01, 3.370D-01, 4.400D+03,
108 * 1.196D-03, 8.784D-03, 1.517D-02, 2.874D-02, 2.488D-02,
109 * 4.464D-02, 8.330D-02, 2.008D-01, 2.360D-01, 3.567D-01/
110C
111 DATA (DCHN(I),I=1,90) /
112 * 4.770D-01, 4.750D-01, 4.715D-01, 4.685D-01, 4.650D-01,
113 * 4.610D-01, 4.570D-01, 4.550D-01, 4.500D-01, 4.450D-01,
114 * 4.405D-01, 4.350D-01, 4.300D-01, 4.250D-01, 4.200D-01,
115 * 4.130D-01, 4.060D-01, 4.000D-01, 3.915D-01, 3.840D-01,
116 * 3.760D-01, 3.675D-01, 3.580D-01, 3.500D-01, 3.400D-01,
117 * 3.300D-01, 3.200D-01, 3.100D-01, 3.000D-01, 2.900D-01,
118 * 2.800D-01, 2.700D-01, 2.600D-01, 2.500D-01, 2.400D-01,
119 * 2.315D-01, 2.240D-01, 2.150D-01, 2.060D-01, 2.000D-01,
120 * 1.915D-01, 1.850D-01, 1.780D-01, 1.720D-01, 1.660D-01,
121 * 1.600D-01, 1.550D-01, 1.500D-01, 1.450D-01, 1.400D-01,
122 * 1.360D-01, 1.320D-01, 1.280D-01, 1.250D-01, 1.210D-01,
123 * 1.180D-01, 1.150D-01, 1.120D-01, 1.100D-01, 1.070D-01,
124 * 1.050D-01, 1.030D-01, 1.010D-01, 9.900D-02, 9.700D-02,
125 * 9.550D-02, 9.480D-02, 9.400D-02, 9.200D-02, 9.150D-02,
126 * 9.100D-02, 9.000D-02, 8.990D-02, 8.900D-02, 8.850D-02,
127 * 8.750D-02, 8.700D-02, 8.650D-02, 8.550D-02, 8.500D-02,
128 * 8.499D-02, 8.450D-02, 8.350D-02, 8.300D-02, 8.250D-02,
129 * 8.150D-02, 8.100D-02, 8.030D-02, 8.000D-02, 7.990D-02/
130 DATA (DCHN(I),I=91,143) /
131 * 7.980D-02, 7.950D-02, 7.900D-02, 7.860D-02, 7.800D-02,
132 * 7.750D-02, 7.650D-02, 7.620D-02, 7.600D-02, 7.550D-02,
133 * 7.530D-02, 7.500D-02, 7.499D-02, 7.498D-02, 7.480D-02,
134 * 7.450D-02, 7.400D-02, 7.350D-02, 7.300D-02, 7.250D-02,
135 * 7.230D-02, 7.200D-02, 7.100D-02, 7.050D-02, 7.020D-02,
136 * 7.000D-02, 6.999D-02, 6.995D-02, 6.993D-02, 6.991D-02,
137 * 6.990D-02, 6.870D-02, 6.850D-02, 6.800D-02, 6.780D-02,
138 * 6.750D-02, 6.700D-02, 6.650D-02, 6.630D-02, 6.600D-02,
139 * 6.550D-02, 6.525D-02, 6.510D-02, 6.500D-02, 6.499D-02,
140 * 6.498D-02, 6.496D-02, 6.494D-02, 6.493D-02, 6.490D-02,
141 * 6.488D-02, 6.485D-02, 6.480D-02/
142C
143 DATA DCHNA /
144 * 6.300D+02, 7.810D-02, 1.421D-01, 1.979D-01, 2.479D-01,
145 * 3.360D-01, 5.400D-01, 7.236D-01, 1.000D+00, 1.540D+03,
146 * 2.225D-01, 3.950D-01, 5.279D-01, 6.298D-01, 7.718D-01,
147 * 9.405D-01, 9.835D-01, 1.000D+00, 2.560D+03, 2.625D-01,
148 * 4.550D-01, 5.963D-01, 7.020D-01, 8.380D-01, 9.603D-01,
149 * 9.903D-01, 1.000D+00, 3.520D+03, 4.250D-01, 6.875D-01,
150 * 8.363D-01, 9.163D-01, 9.828D-01, 1.000D+00, 1.000D+00,
151 * 1.000D+00/
152C
153 DATA DCHNB /
154 * 6.300D+02, 3.800D-02, 7.164D-02, 1.275D-01, 2.171D-01,
155 * 3.227D-01, 4.091D-01, 5.051D-01, 6.061D-01, 7.074D-01,
156 * 8.434D-01, 1.000D+00, 2.040D+03, 1.200D-01, 2.115D-01,
157 * 3.395D-01, 5.295D-01, 7.251D-01, 8.511D-01, 9.487D-01,
158 * 9.987D-01, 1.000D+00, 1.000D+00, 1.000D+00, 2.200D+03,
159 * 1.344D-01, 2.324D-01, 3.754D-01, 5.674D-01, 7.624D-01,
160 * 8.896D-01, 9.808D-01, 1.000D+00, 1.000D+00, 1.000D+00,
161 * 1.000D+00, 2.850D+03, 2.330D-01, 4.130D-01, 6.610D-01,
162 * 9.010D-01, 9.970D-01, 1.000D+00, 1.000D+00, 1.000D+00,
163 * 1.000D+00, 1.000D+00, 1.000D+00, 3.500D+03, 3.300D-01,
164 * 5.450D-01, 7.950D-01, 1.000D+00, 1.000D+00, 1.000D+00,
165 * 1.000D+00, 1.000D+00, 1.000D+00, 1.000D+00, 1.000D+00/
166C
167C---------------------------------------------------------------
168 CST=1.0D+0
169C
170C* IS THE KINETIC ENERGY GREATER THAN LIMIT ?
171C
172 IF (EKIN .GT. 3.5D0) RETURN
173C
174 IF(KPROJ.EQ.8) GOTO 101
175 IF(KPROJ.EQ.1) GOTO 102
176C* INVALID REACTION
177 RETURN
178C-------------------------------- NP ELASTIC SCATTERING----------
179101 CONTINUE
180 IF (EKIN.GT.0.740D+0)GOTO 1000
181 IF (EKIN.LT.0.300D+0)THEN
182C EKIN .LT. 300 MEV
183 IDAT=1
184 ELSE
185C 300 MEV < EKIN < 740 MEV
186 IDAT=6
187 END IF
188C
189 ENER=EKIN
190 IE=ABS(ENER/0.020D+0)
191 UNIV=(ENER-IE*0.020D+0)/0.020D+0
192C FORWARD/BACKWARD DECISION
193 K=IDAT+5*IE
194 BWFW=(DCLIN(K+5)-DCLIN(K))*UNIV + DCLIN(K)
195 CALL GRNDM(RNDM,1)
196 RND=RNDM(1)
197 IF (RND.LT.BWFW)THEN
198 VALUE2=-1.0D+0
199 K=K+1
200 ELSE
201 VALUE2=1.0D+0
202 K=K+3
203 END IF
204C
205 COEF=(DCLIN(K+5)-DCLIN(K))*UNIV + DCLIN(K)
206 CALL GRNDM(RNDM,1)
207 RND=RNDM(1)
208C
209 IF(RND.LT.COEF)THEN
210 CALL GRNDM(RNDM,1)
211 CST=RNDM(1)
212 CST=CST*VALUE2
213 ELSE
214 CALL GRNDM(RNDM,4)
215 R1=RNDM(1)
216 R2=RNDM(2)
217 R3=RNDM(3)
218 R4=RNDM(4)
219C
220 IF(VALUE2.GT.0.D0)THEN
221 CST=MAX(R1,R2,R3,R4)
222 GOTO 1500
223 ELSE
224 CALL GRNDM(RNDM,1)
225 R5=RNDM(1)
226C
227 IF (IDAT.EQ.1)THEN
228 CST=-MAX(R1,R2,R3,R4,R5)
229 ELSE
230 CALL GRNDM(RNDM,2)
231 R6=RNDM(1)
232 R7=RNDM(2)
233 CST=-MAX(R1,R2,R3,R4,R5,R6,R7)
234 END IF
235C
236 END IF
237C
238 END IF
239C
240 GOTO 1500
241C
242C******** EKIN .GT. 0.74 GEV
243C
2441000 CONTINUE
245 ENER=EKIN - 0.66D0
246C IE=ABS(ENER/0.02D0)
247 IE=ENER/0.02D0
248 EMEV=EKIN*1.D+03
249C
250 UNIV=(ENER-IE*0.020D+0)/0.020D+0
251 K=IE
252 BWFW=(DCHN(K+1)-DCHN(K))*UNIV + DCHN(K)
253 CALL GRNDM(RNDM,1)
254 RND=RNDM(1)
255C FORWARD NEUTRON
256 IF (RND.GE.BWFW)THEN
257 DO 1200 K=10,36,9
258 IF (DCHNA(K).GT.EMEV) THEN
259 UNIVE=(EMEV-DCHNA(K-9))/(DCHNA(K)-DCHNA(K-9))
260 CALL GRNDM(RNDM,1)
261 UNIV=RNDM(1)
262 DO 1100 I=1,8
263 II=K+I
264 P=(DCHNA(II)-DCHNA(II-9))*UNIVE + DCHNA(II-9)
265C
266 IF (P.GT.UNIV)THEN
267 CALL GRNDM(RNDM,1)
268 UNIV=RNDM(1)
269 FLTI=I-UNIV
270 GOTO(290,290,290,290,330,340,350,360) I
271 END IF
272 1100 CONTINUE
273 END IF
274 1200 CONTINUE
275C
276 ELSE
277C BACKWARD NEUTRON
278 DO 1400 K=13,60,12
279 IF (DCHNB(K).GT.EMEV) THEN
280 UNIVE=(EMEV-DCHNB(K-12))/(DCHNB(K)-DCHNB(K-12))
281 CALL GRNDM(RNDM,1)
282 UNIV=RNDM(1)
283 DO 1300 I=1,11
284 II=K+I
285 P=(DCHNB(II)-DCHNB(II-12))*UNIVE + DCHNB(II-12)
286C
287 IF (P.GT.UNIV)THEN
288 CALL GRNDM(RNDM,1)
289 UNIV=RNDM(1)
290 FLTI=I-UNIV
291 GOTO(120,120,140,150,160,160,180,190,200,210,220) I
292 END IF
293 1300 CONTINUE
294 END IF
295 1400 CONTINUE
296 END IF
297C
298120 CST=1.0D-2*FLTI-1.0D+0
299 GOTO 1500
300140 CST=2.0D-2*UNIV-0.98D+0
301 GOTO 1500
302150 CST=4.0D-2*UNIV-0.96D+0
303 GOTO 1500
304160 CST=6.0D-2*FLTI-1.16D+0
305 GOTO 1500
306180 CST=8.0D-2*UNIV-0.80D+0
307 GOTO 1500
308190 CST=1.0D-1*UNIV-0.72D+0
309 GOTO 1500
310200 CST=1.2D-1*UNIV-0.62D+0
311 GOTO 1500
312210 CST=2.0D-1*UNIV-0.50D+0
313 GOTO 1500
314220 CST=3.0D-1*(UNIV-1.0D+0)
315 GOTO 1500
316C
317290 CST=1.0D0 - 2.5D-2*FLTI
318 GOTO 1500
319330 CST=0.85D0 + 0.5D-1*UNIV
320 GOTO 1500
321340 CST=0.70D0 + 1.5D-1*UNIV
322 GOTO 1500
323350 CST=0.50D0 + 2.0D-1*UNIV
324 GOTO 1500
325360 CST=0.50D0*UNIV
326C
3271500 RETURN
328C
329C----------------------------------- PP ELASTIC SCATTERING -------
330C
331 102 CONTINUE
332 EMEV=EKIN*1.D+03
333C
334 IF (EKIN.LE.0.500D0) THEN
335 CALL GRNDM(RNDM,1)
336 RND=RNDM(1)
337 CST=2.0D0*RND-1.0D0
338 RETURN
339C
340 ELSE IF (EKIN.LT.1.D0) THEN
341 DO 2200 K=13,60,12
342 IF (PDCI(K).GT.EMEV) THEN
343 UNIVE=(EMEV-PDCI(K-12))/(PDCI(K)-PDCI(K-12))
344 CALL GRNDM(RNDM,1)
345 UNIV=RNDM(1)
346 SUM=0.0D0
347 DO 2100 I=1,11
348 II=K+I
349 SUM=SUM + (PDCI(II)-PDCI(II-12))*UNIVE + PDCI(II-12)
350C
351 IF (UNIV.LT.SUM)THEN
352 CALL GRNDM(RNDM,1)
353 UNIV=RNDM(1)
354 FLTI=I-UNIV
355 GOTO(55,55,55,60,60,65,65,65,65,70,70) I
356 END IF
357 2100 CONTINUE
358 END IF
359 2200 CONTINUE
360 ELSE
361 DO 2400 K=12,55,11
362 IF (PDCH(K).GT.EMEV) THEN
363 UNIVE=(EMEV-PDCH(K-11))/(PDCH(K)-PDCH(K-11))
364 CALL GRNDM(RNDM,1)
365 UNIV=RNDM(1)
366 SUM=0.D0
367 DO 2300 I=1,10
368 II=K+I
369 SUM=SUM + (PDCH(II)-PDCH(II-11))*UNIVE + PDCH(II-11)
370C
371 IF (UNIV.LT.SUM)THEN
372 CALL GRNDM(RNDM,1)
373 UNIV=RNDM(1)
374 FLTI=UNIV+I
375 GOTO(50,55,60,60,65,65,65,65,70,70) I
376 END IF
377 2300 CONTINUE
378 END IF
379 2400 CONTINUE
380 END IF
381C
38250 CST=0.4D0*UNIV
383 GOTO 2500
38455 CST=0.2D0*FLTI
385 GOTO 2500
38660 CST=0.3D0 + 0.1D0*FLTI
387 GOTO 2500
38865 CST=0.6D0 + 0.04D0*FLTI
389 GOTO 2500
39070 CST=0.78D0 + 0.02D0*FLTI
391C
3922500 CONTINUE
393 CALL GRNDM(RNDM,1)
394 RND=RNDM(1)
395 IF (RND.GT.0.5D+00)CST=-CST
396C
397 RETURN
398 END