This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / gheisha / casn.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:21:00  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.38  by  S.Giani
11 *-- Author :
12       SUBROUTINE CASN(K,INT,NFL)
13 C
14 C *** CASCADE OF NEUTRON ***
15 C *** NVE 04-MAY-1988 CERN GENEVA ***
16 C
17 C ORIGIN : H.FESEFELDT (13-SEP-1987)
18 C
19 C N  UNDERGOES INTERACTION WITH NUCLEON WITHIN NUCLEUS.
20 C CHECK IF ENERGETICALLY POSSIBLE TO PRODUCE PIONS/KAONS.
21 C IF NOT ASSUME NUCLEAR EXCITATION OCCURS AND INPUT PARTICLE
22 C IS DEGRADED IN ENERGY.    NO OTHER PARTICLES PRODUCED.
23 C IF REACTION IS POSSIBLE FIND CORRECT NUMBER OF PIONS/PROTONS/
24 C NEUTRONS PRODUCED USING AN INTERPOLATION TO MULTIPLICITY DATA.
25 C REPLACE SOME PIONS OR PROTONS/NEUTRONS BY KAONS OR STRANGE BARYONS
26 C ACCORDING TO AVERAGE MULTIPLICITY PER INELASTIC REACTIONS.
27 C
28 #include "geant321/mxgkgh.inc"
29 #include "geant321/s_consts.inc"
30 #include "geant321/s_curpar.inc"
31 #include "geant321/s_result.inc"
32 #include "geant321/s_prntfl.inc"
33 #include "geant321/limits.inc"
34 #include "geant321/s_kginit.inc"
35 C
36       REAL N
37       DIMENSION PMUL(2,1200),ANORM(2,60),SUPP(10),CECH(10),B(2)
38       DIMENSION RNDM(1)
39       SAVE PMUL,ANORM
40       DATA SUPP/0.,0.4,0.55,0.65,0.75,0.82,0.86,0.90,0.94,0.98/
41       DATA CECH/0.50,0.45,0.40,0.35,0.30,0.25,0.06,0.04,0.005,0./
42       DATA B/0.35,0.0/,C/1.25/
43 C
44 C --- INITIALIZATION INDICATED BY KGINIT(17) ---
45       IF (KGINIT(17) .NE. 0) GO TO 10
46       KGINIT(17)=1
47 C
48 C --- INITIALIZE PMUL AND ANORM ARRAYS ---
49       DO 9000 J=1,1200
50       DO 9001 I=1,2
51       PMUL(I,J)=0.0
52       IF (J .LE. 60) ANORM(I,J)=0.0
53  9001 CONTINUE
54  9000 CONTINUE
55 C
56 C** COMPUTE NORMALIZATION CONSTANTS
57 C** FOR N AS TARGET
58 C
59       L=0
60       DO 1 NP1=1,20
61       NP=NP1-1
62       NPP1=NP1+2
63       DO 1 NM1=NP1,NPP1
64       NM=NM1-1
65       DO 1 NZ1=1,20
66       NZ=NZ1-1
67       L=L+1
68       IF(L.GT.1200) GOTO 1
69       NPROT= -NP+NM
70       NNEUT=2-NPROT
71       NT=NP+NM+NZ
72       IF(NT.LE.0.OR.NT.GT.60) GOTO 1
73       PMUL(1,L)=PMLTPC(NP,NM,NZ,NT,B(2),C)
74       NPROTF=NFAC(NPROT)
75       NNEUTF=NFAC(NNEUT)
76       PMUL(1,L)=PMUL(1,L)/(NPROTF*NNEUTF)
77       ANORM(1,NT)=ANORM(1,NT)+PMUL(1,L)
78     1 CONTINUE
79 C** FOR P AS TARGET
80       L=0
81       DO 2 NP1=1,20
82       NP=NP1-1
83       NMM1=NP1-1
84       IF(NMM1.LE.1) NMM1=1
85       NPP1=NP1+1
86       DO 2 NM1=NMM1,NPP1
87       NM=NM1-1
88       DO 2 NZ1=1,20
89       NZ=NZ1-1
90       L=L+1
91       IF(L.GT.1200) GOTO 2
92       NPROT=1-NP+NM
93       NNEUT=2-NPROT
94       NT=NP+NM+NZ
95       IF(NT.LE.0.OR.NT.GT.60) GOTO 2
96       PMUL(2,L)=PMLTPC(NP,NM,NZ,NT,B(1),C)
97       NPROTF=NFAC(NPROT)
98       NNEUTF=NFAC(NNEUT)
99       PMUL(2,L)=PMUL(2,L)/(NPROTF*NNEUTF)
100       ANORM(2,NT)=ANORM(2,NT)+PMUL(2,L)
101     2 CONTINUE
102       DO 3 I=1,60
103       IF(ANORM(1,I).GT.0.) ANORM(1,I)=1./ANORM(1,I)
104       IF(ANORM(2,I).GT.0.) ANORM(2,I)=1./ANORM(2,I)
105     3 CONTINUE
106       IF(.NOT.NPRT(10)) GOTO 10
107       WRITE(NEWBCD,2001)
108       DO 4 NFL=1,2
109       WRITE(NEWBCD,2002) NFL
110       WRITE(NEWBCD,2003) (ANORM(NFL,I),I=1,60)
111       WRITE(NEWBCD,2003) (PMUL(NFL,I),I=1,1200)
112     4 CONTINUE
113 C**  CHOOSE PROTON OR NEUTRON AS TARGET
114    10 NFL=2
115       CALL GRNDM(RNDM,1)
116       IF(RNDM(1).LT.ZNO2/ATNO2) NFL=1
117       TARMAS=RMASS(14)
118       IF (NFL .EQ. 2) TARMAS=RMASS(16)
119       S=AMASQ+TARMAS**2+2.0*TARMAS*EN
120       RS=SQRT(S)
121       ENP(8)=AMASQ+TARMAS**2+2.0*TARMAS*ENP(6)
122       ENP(9)=SQRT(ENP(8))
123       EAB=RS-TARMAS-RMASS(16)
124 C**  ELASTIC SCATTERING
125       NP=0
126       NM=0
127       NZ=0
128       N=0.
129       NCECH=0
130       IF(INT.EQ.2) GOTO 20
131 C** INTRODUCE CHARGE EXCHANGE REACTION PN --> NP
132       IF(NFL.EQ.2) GOTO 100
133       IPLAB=IFIX(P*2.5)+1
134       IF(IPLAB.GT.10) IPLAB=10
135       CALL GRNDM(RNDM,1)
136       IF(RNDM(1).GT.CECH(IPLAB)/ATNO2**0.42) GOTO 100
137       NCECH=1
138       GOTO 100
139 C**  CHECK IF ENERGETICALLY POSSIBLE TO PRODUCE ONE EXTRA PION IN REACT.
140   20  IF (EAB .LE. RMASS(7)) GOTO 55
141 C**  SUPPRESSION OF HIGH MULTIPLICITY EVENTS AT LOW MOMENTUM
142       IEAB=IFIX(EAB*5.)+1
143       IF(IEAB.GT.10) GOTO 22
144       CALL GRNDM(RNDM,1)
145       IF(RNDM(1).LT.SUPP(IEAB)) GOTO 22
146       N=1.
147       GOTO (24,23),NFL
148  23   CONTINUE
149       TEST=-(1+B(2))**2/(2.0*C**2)
150       IF (TEST .GT. EXPXU) TEST=EXPXU
151       IF (TEST .LT. EXPXL) TEST=EXPXL
152       W0=EXP(TEST)/2.0
153       WM=EXP(TEST)
154       CALL GRNDM(RNDM,1)
155       RAN=RNDM(1)
156       NP=0
157       NM=0
158       NZ=1
159       IF(RAN.LT.W0/(W0+WM)) GOTO 100
160       NP=0
161       NM=1
162       NZ=0
163       GOTO 100
164  24   CONTINUE
165       TEST=-(1+B(1))**2/(2.0*C**2)
166       IF (TEST .GT. EXPXU) TEST=EXPXU
167       IF (TEST .LT. EXPXL) TEST=EXPXL
168       W0=EXP(TEST)
169       WP=EXP(TEST)/2.0
170       TEST=-(-1+B(1))**2/(2.0*C**2)
171       IF (TEST .GT. EXPXU) TEST=EXPXU
172       IF (TEST .LT. EXPXL) TEST=EXPXL
173       WM=EXP(TEST)/2.0
174       WT=W0+WP+WM
175       WP=W0+WP
176       CALL GRNDM(RNDM,1)
177       RAN=RNDM(1)
178       NP=0
179       NM=0
180       NZ=1
181       IF(RAN.LT.W0/WT) GOTO 100
182       NP=1
183       NM=0
184       NZ=0
185       IF(RAN.LT.WP/WT) GOTO 100
186       NP=0
187       NM=1
188       NZ=0
189       GOTO 100
190 C
191    22 ALEAB=LOG(EAB)
192 C** NO. OF TOTAL PARTICLES VS SQRT(S)-2*MP
193       N=3.62567+0.665843*ALEAB+0.336514*ALEAB*ALEAB
194      * +0.117712*ALEAB*ALEAB*ALEAB+0.0136912*ALEAB*ALEAB*ALEAB*ALEAB
195       N=N-2.
196 C** NORMALIZATION CONSTANT FOR  KNO-DISTRIBUTION
197       ANPN=0.
198       DO 21 NT=1,60
199       TEST=-(PI/4.0)*(NT/N)**2
200       IF (TEST .GT. EXPXU) TEST=EXPXU
201       IF (TEST .LT. EXPXL) TEST=EXPXL
202       ANPN=ANPN+PI*NT*EXP(TEST)/(2.0*N*N)
203    21 CONTINUE
204       ANPN=1./ANPN
205 C** P OR N AS TARGET
206       CALL GRNDM(RNDM,1)
207       RAN=RNDM(1)
208       EXCS=0.
209       GOTO (40,30),NFL
210 C** FOR N AS TARGET
211    30 L=0
212       DO 31 NP1=1,20
213       NP=NP1-1
214       NPP1=NP1+2
215       DO 31 NM1=NP1,NPP1
216       NM=NM1-1
217       DO 31 NZ1=1,20
218       NZ=NZ1-1
219       L=L+1
220       IF(L.GT.1200) GOTO 31
221       NT=NP+NM+NZ
222       IF(NT.LE.0.OR.NT.GT.60) GOTO 31
223       TEST=-(PI/4.0)*(NT/N)**2
224       IF (TEST .GT. EXPXU) TEST=EXPXU
225       IF (TEST .LT. EXPXL) TEST=EXPXL
226       DUM1=ANPN*PI*NT*PMUL(1,L)*ANORM(1,NT)/(2.0*N*N)
227       DUM2=ABS(DUM1)
228       DUM3=EXP(TEST)
229       ADDNVE=0.0
230       IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
231       IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GE. 1.0E-10)) ADDNVE=DUM1*DUM3
232       EXCS=EXCS+ADDNVE
233       IF(RAN.LT.EXCS) GOTO 100
234    31 CONTINUE
235       GOTO 80
236 C** FOR P AS TARGET
237    40 L=0
238       DO 41 NP1=1,20
239       NP=NP1-1
240       NMM1=NP1-1
241       IF(NMM1.LE.1) NMM1=1
242       NPP1=NP1+1
243       DO 41 NM1=NMM1,NPP1
244       NM=NM1-1
245       DO 41 NZ1=1,20
246       NZ=NZ1-1
247       L=L+1
248       IF(L.GT.1200) GOTO 41
249       NT=NP+NM+NZ
250       IF(NT.LE.0.OR.NT.GT.60) GOTO 41
251       TEST=-(PI/4.0)*(NT/N)**2
252       IF (TEST .GT. EXPXU) TEST=EXPXU
253       IF (TEST .LT. EXPXL) TEST=EXPXL
254       DUM1=ANPN*PI*NT*PMUL(2,L)*ANORM(2,NT)/(2.0*N*N)
255       DUM2=ABS(DUM1)
256       DUM3=EXP(TEST)
257       ADDNVE=0.0
258       IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
259       IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GE. 1.0E-10)) ADDNVE=DUM1*DUM3
260       EXCS=EXCS+ADDNVE
261       IF(RAN.LT.EXCS) GOTO 100
262    41 CONTINUE
263       GOTO 80
264    50 IF(NPRT(4))
265      *WRITE(NEWBCD,1003) EAB,N,NFL,NP,NM,NZ
266       CALL STPAIR
267       IF(INT.EQ.1) CALL TWOB(16,NFL,N)
268       IF(INT.EQ.2) CALL GENXPT(16,NFL,N)
269       GO TO 9999
270    55 IF(NPRT(4))
271      *WRITE(NEWBCD,1001)
272       GOTO 53
273 C** EXCLUSIVE REACTION NOT FOUND
274    80 IF(NPRT(4))
275      *WRITE(NEWBCD,1004) RS,N
276    53 INT=1
277       NP=0
278       NM=0
279       NZ=0
280   100 DO 101 I=1,60
281   101 IPA(I)=0
282       IF(INT.LE.0) GOTO 131
283       NPROT=1-NP+NM+(1-NFL)
284       NNEUT=2-NPROT
285       GOTO (112,102),NFL
286   102 GOTO (103,104),INT
287   103 IPA(1)=16
288       IPA(2)=16
289       NT=2
290       GOTO 130
291   104 IF(NNEUT.EQ.1) GOTO 105
292       IF(NNEUT.EQ.2) GOTO 106
293       IPA(1)=14
294       IPA(2)=14
295       GOTO 120
296   105 IPA(1)=14
297       IPA(2)=16
298       CALL GRNDM(RNDM,1)
299       IF(RNDM(1).LT.0.5) GOTO 120
300       IPA(1)=16
301       IPA(2)=14
302       GOTO 120
303   106 IPA(1)=16
304       IPA(2)=16
305       GOTO 120
306   112 GOTO (113,114),INT
307   113 IPA(1)=16
308       IPA(2)=14
309       NT=2
310       IF(NCECH.EQ.0) GOTO 130
311       IPA(1)=14
312       IPA(2)=16
313       GOTO 130
314   114 IF(NNEUT.EQ.1) GOTO 115
315       IF(NNEUT.EQ.2) GOTO 116
316       IPA(1)=14
317       IPA(2)=14
318       GOTO 120
319   115 IPA(1)=14
320       IPA(2)=16
321       CALL GRNDM(RNDM,1)
322       IF(RNDM(1).LT.0.33) GOTO 120
323       IPA(1)=16
324       IPA(2)=14
325       GOTO 120
326   116 IPA(1)=16
327       IPA(2)=16
328   120 NT=2
329       IF(NP.EQ.0) GOTO 122
330       DO 121 I=1,NP
331       NT=NT+1
332   121 IPA(NT)=7
333   122 IF(NM.EQ.0) GOTO 124
334       DO 123 I=1,NM
335       NT=NT+1
336   123 IPA(NT)=9
337   124 IF(NZ.EQ.0) GOTO 130
338       DO 125 I=1,NZ
339       NT=NT+1
340   125 IPA(NT)=8
341   130 IF(NPRT(4))
342      *WRITE(NEWBCD,2004) NT,(IPA(I),I=1,20)
343       GOTO 50
344   131 IF(NPRT(4))
345      *WRITE(NEWBCD,2005)
346 C
347 1001  FORMAT('0*CASN* CASCADE ENERGETICALLY NOT POSSIBLE NUCLEAR',
348      * ' EXCITATION',2X,F8.4,2X,'INCIDENT ENERGY LOST')
349 1003  FORMAT(' *CASN* NEUTRON-INDUCED CASCADE,',
350      $ ' AVAIL. ENERGY',2X,F8.4,
351      $ 2X,'<NTOT>',2X,F8.4,2X,'FROM',4(2X,I3),2X,'PARTICLES')
352 1004  FORMAT(' *CASN* NEUTRON-INDUCED CASCADE,',
353      $ ' EXCLUSIVE REACTION NOT FOUND',
354      $ ' TRY ELASTIC SCATTERING  AVAIL. ENERGY',2X,F8.4,2X,
355      $ '<NTOT>',2X,F8.4)
356 2001  FORMAT('0*CASN* TABLES FOR MULT. DATA NEUTRON INDUCED REACTION',
357      $ ' FOR DEFINITION OF NUMBERS SEE FORTRAN CODING')
358 2002  FORMAT(' *CASN* TARGET PARTICLE FLAG',2X,I5)
359 2003  FORMAT(1H ,10E12.4)
360 2004  FORMAT(' *CASN* ',I3,2X,'PARTICLES , MASS INDEX ARRAY',2X,20I4)
361 2005  FORMAT(' *CASN* NO PARTICLES PRODUCED')
362 C
363  9999 CONTINUE
364       END