This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / gheisha / caspip.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 CASPIP(K,INT,NFL)
13 C
14 C *** CASCADE OF PI+ ***
15 C *** NVE 04-MAY-1988 CERN GENEVA ***
16 C
17 C ORIGIN : H.FESEFELDT (18-SEP-1987)
18 C
19 C PI+  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.2,0.45,0.55,0.65,0.75,0.85,0.90,0.94,0.98/
41       DATA CECH/0.33,0.27,0.29,0.31,0.27,0.18,0.13,0.10,0.09,0.07/
42       DATA B/0.7,0.7/,C/1.25/
43 C
44 C --- INITIALIZATION INDICATED BY KGINIT(18) ---
45       IF (KGINIT(18) .NE. 0) GO TO 10
46       KGINIT(18)=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 P AS TARGET
58 C
59       L=0
60       DO 1 NP1=1,20
61       NP=NP1-1
62       NMM1=NP1-2
63       IF(NMM1.LE.1) NMM1=1
64       DO 1 NM1=NMM1,NP1
65       NM=NM1-1
66       DO 1 NZ1=1,20
67       NZ=NZ1-1
68       L=L+1
69       IF(L.GT.1200) GOTO 1
70       NT=NP+NM+NZ
71       IF(NT.LE.0.OR.NT.GT.60) GOTO 1
72       PMUL(1,L)=PMLTPC(NP,NM,NZ,NT,B(1),C)
73       ANORM(1,NT)=ANORM(1,NT)+PMUL(1,L)
74     1 CONTINUE
75 C** FOR N AS TARGET
76       L=0
77       DO 2 NP1=1,20
78       NP=NP1-1
79       NMM1=NP1-1
80       IF(NMM1.LE.1) NMM1=1
81       NPP1=NP1+1
82       DO 2 NM1=NMM1,NPP1
83       NM=NM1-1
84       DO 2 NZ1=1,20
85       NZ=NZ1-1
86       L=L+1
87       IF(L.GT.1200) GOTO 2
88       NT=NP+NM+NZ
89       IF(NT.LE.0.OR.NT.GT.60) GOTO 2
90       PMUL(2,L)=PMLTPC(NP,NM,NZ,NT,B(2),C)
91       ANORM(2,NT)=ANORM(2,NT)+PMUL(2,L)
92     2 CONTINUE
93       DO 3 I=1,60
94       IF(ANORM(1,I).GT.0.) ANORM(1,I)=1./ANORM(1,I)
95       IF(ANORM(2,I).GT.0.) ANORM(2,I)=1./ANORM(2,I)
96     3 CONTINUE
97       IF(.NOT.NPRT(10)) GOTO 10
98       WRITE(NEWBCD,2001)
99       DO 4 NFL=1,2
100       WRITE(NEWBCD,2002) NFL
101       WRITE(NEWBCD,2003) (ANORM(NFL,I),I=1,60)
102       WRITE(NEWBCD,2003) (PMUL(NFL,I),I=1,1200)
103     4 CONTINUE
104 C**  CHOOSE PROTON OR NEUTRON AS TARGET
105    10 NFL=2
106       CALL GRNDM(RNDM,1)
107       IF(RNDM(1).LT.ZNO2/ATNO2) NFL=1
108       TARMAS=RMASS(14)
109       IF (NFL .EQ. 2) TARMAS=RMASS(16)
110       S=AMASQ+TARMAS**2+2.0*TARMAS*EN
111       RS=SQRT(S)
112       ENP(8)=AMASQ+TARMAS**2+2.0*TARMAS*ENP(6)
113       ENP(9)=SQRT(ENP(8))
114       EAB=RS-TARMAS-RMASS(7)
115 C
116 C**  ELASTIC SCATTERING
117       NP=0
118       NM=0
119       NZ=0
120       N=0.
121       IPA(1)=7
122       IPA(2)=14
123       IF(NFL.EQ.2) IPA(2)=16
124       IF(INT.EQ.2) GOTO 20
125 C**  FOR PI+ N REACTIONS CHANGE SOME OF THE ELASTIC CROSS SECTION
126 C**  TO PI+ N --> PI0 P
127       IF(NFL.EQ.1) GOTO 100
128       IPLAB=IFIX(P   *5.)+1
129       IF(IPLAB.GT.10) IPLAB=10
130       CALL GRNDM(RNDM,1)
131       IF(RNDM(1).GT.CECH(IPLAB)/ATNO2**0.42) GOTO 100
132       IPA(1)=8
133       IPA(2)=14
134       GOTO 100
135 C**  CHECK IF ENERGETICALLY POSSIBLE TO PRODUCE ONE EXTRA PION IN REACT.
136   20  IF (EAB .LE. RMASS(7)) GOTO 55
137 C**  SUPPRESSION OF HIGH MULTIPLICITY EVENTS AT LOW MOMENTUM
138       IEAB=IFIX(EAB*5.)+1
139       IF(IEAB.GT.10) GOTO 22
140       CALL GRNDM(RNDM,1)
141       IF(RNDM(1).LT.SUPP(IEAB)) GOTO 22
142       N=1.
143       GOTO (23,24),NFL
144  23   CONTINUE
145       TEST=-(1+B(1))**2/(2.0*C**2)
146       IF (TEST .LE. EXPXL) TEST=EXPXL
147       IF (TEST .GE. EXPXU) TEST=EXPXU
148       W0=EXP(TEST)
149       WP=EXP(TEST)
150       CALL GRNDM(RNDM,1)
151       RAN=RNDM(1)
152       NP=0
153       NM=0
154       NZ=1
155       IF(RAN.LT.W0/(W0+WP)) GOTO 50
156       NP=1
157       NM=0
158       NZ=0
159       GOTO 50
160  24   CONTINUE
161       TEST=-(1+B(2))**2/(2.0*C**2)
162       IF (TEST .LE. EXPXL) TEST=EXPXL
163       IF (TEST .GE. EXPXU) TEST=EXPXU
164       W0=EXP(TEST)
165       WP=EXP(TEST)
166       TEST=-(-1+B(2))**2/(2.0*C**2)
167       IF (TEST .LE. EXPXL) TEST=EXPXL
168       IF (TEST .GE. EXPXU) TEST=EXPXU
169       WM=EXP(TEST)
170       WT=W0+WP+WM
171       WP=W0+WP
172       CALL GRNDM(RNDM,1)
173       RAN=RNDM(1)
174       NP=0
175       NM=0
176       NZ=1
177       IF(RAN.LT.W0/WT) GOTO 50
178       NP=1
179       NM=0
180       NZ=0
181       IF(RAN.LT.WP/WT) GOTO 50
182       NP=0
183       NM=1
184       NZ=0
185       GOTO 50
186 C
187    22 ALEAB=LOG(EAB)
188 C** NO. OF TOTAL PARTICLES VS SQRT(S)-2*MP
189       N=3.62567+0.665843*ALEAB+0.336514*ALEAB*ALEAB
190      * +0.117712*ALEAB*ALEAB*ALEAB+0.0136912*ALEAB*ALEAB*ALEAB*ALEAB
191       N=N-2.
192 C** NORMALIZATION CONSTANT FOR  KNO-DISTRIBUTION
193       ANPN=0.
194       DO 21 NT=1,60
195       TEST=-(PI/4.0)*(NT/N)**2
196       IF (TEST .LE. EXPXL) TEST=EXPXL
197       IF (TEST .GE. EXPXU) TEST=EXPXU
198       DUM1=PI*NT/(2.0*N*N)
199       DUM2=ABS(DUM1)
200       DUM3=EXP(TEST)
201       ADDNVE=0.0
202       IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
203       IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GE. 1.0E-10)) ADDNVE=DUM1*DUM3
204       ANPN=ANPN+ADDNVE
205    21 CONTINUE
206       ANPN=1./ANPN
207 C** P OR N AS TARGET
208       CALL GRNDM(RNDM,1)
209       RAN=RNDM(1)
210       EXCS=0.
211       GOTO (30,40),NFL
212 C** FOR P AS TARGET
213    30 L=0
214       DO 31 NP1=1,20
215       NP=NP1-1
216       NMM1=NP1-2
217       IF(NMM1.LE.1) NMM1=1
218       DO 31 NM1=NMM1,NP1
219       NM=NM1-1
220       DO 31 NZ1=1,20
221       NZ=NZ1-1
222       L=L+1
223       IF(L.GT.1200) GOTO 31
224       NT=NP+NM+NZ
225       IF(NT.LE.0.OR.NT.GT.60) GOTO 31
226       TEST=-(PI/4.0)*(NT/N)**2
227       IF (TEST .LE. EXPXL) TEST=EXPXL
228       IF (TEST .GE. EXPXU) TEST=EXPXU
229       DUM1=ANPN*PI*NT*PMUL(1,L)*ANORM(1,NT)/(2.0*N*N)
230       DUM2=ABS(DUM1)
231       DUM3=EXP(TEST)
232       ADDNVE=0.0
233       IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
234       IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GE. 1.0E-10)) ADDNVE=DUM1*DUM3
235       EXCS=EXCS+ADDNVE
236       IF(RAN.LT.EXCS) GOTO 50
237    31 CONTINUE
238       GOTO 80
239 C** FOR N AS TARGET
240    40 L=0
241       DO 41 NP1=1,20
242       NP=NP1-1
243       NMM1=NP1-1
244       IF(NMM1.LE.1) NMM1=1
245       NPP1=NP1+1
246       DO 41 NM1=NMM1,NPP1
247       NM=NM1-1
248       DO 41 NZ1=1,20
249       NZ=NZ1-1
250       L=L+1
251       IF(L.GT.1200) GOTO 41
252       NT=NP+NM+NZ
253       IF(NT.LE.0.OR.NT.GT.60) GOTO 41
254       TEST=-(PI/4.0)*(NT/N)**2
255       IF (TEST .LE. EXPXL) TEST=EXPXL
256       IF (TEST .GE. EXPXU) TEST=EXPXU
257       DUM1=ANPN*PI*NT*PMUL(2,L)*ANORM(2,NT)/(2.0*N*N)
258       DUM2=ABS(DUM1)
259       DUM3=EXP(TEST)
260       ADDNVE=0.0
261       IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
262       IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GT. 1.0E-10)) ADDNVE=DUM1*DUM3
263       EXCS=EXCS+ADDNVE
264       IF(RAN.LT.EXCS) GOTO 50
265    41 CONTINUE
266       GOTO 80
267    50 GOTO (60,65),NFL
268    60 IF(NP.EQ.1+NM) GOTO 61
269       IF(NP.EQ.2+NM) GOTO 63
270       IPA(1)=7
271       IPA(2)=14
272       GOTO 100
273    61 CALL GRNDM(RNDM,1)
274       IF(RNDM(1).LT.0.5) GOTO 62
275       IPA(1)=7
276       IPA(2)=16
277       GOTO 100
278    62 IPA(1)=8
279       IPA(2)=14
280       GOTO 100
281    63 IPA(1)=8
282       IPA(2)=16
283       GOTO 100
284    65 IF(NP.EQ.NM) GOTO 66
285       IF(NP.EQ.1+NM) GOTO 68
286       IPA(1)=7
287       IPA(2)=14
288       GOTO 100
289    66 CALL GRNDM(RNDM,1)
290       IF(RNDM(1).LT.0.25) GOTO 67
291       IPA(1)=7
292       IPA(2)=16
293       GOTO 100
294    67 IPA(1)=8
295       IPA(2)=14
296       GOTO 100
297    68 IPA(1)=8
298       IPA(2)=16
299       GOTO 100
300    70 IF(NPRT(4))
301      *WRITE(NEWBCD,1003) EAB,N,NFL,NP,NM,NZ
302       CALL STPAIR
303       IF(INT.EQ.1) CALL TWOB(7,NFL,N)
304       IF(INT.EQ.2) CALL GENXPT(7,NFL,N)
305       GO TO 9999
306    55 IF(NPRT(4))
307      *WRITE(NEWBCD,1001)
308       GOTO 53
309 C** EXCLUSIVE REACTION NOT FOUND
310    80 IF(NPRT(4))
311      *WRITE(NEWBCD,1004) RS,N
312    53 INT=1
313       NP=0
314       NM=0
315       NZ=0
316       N=0.
317       IPA(1)=7
318       IPA(2)=14
319       IF(NFL.EQ.2) IPA(2)=16
320   100 DO 101 I=3,60
321   101 IPA(I)=0
322       IF(INT.LE.0) GOTO 131
323   120 NT=2
324       IF(NP.EQ.0) GOTO 122
325       DO 121 I=1,NP
326       NT=NT+1
327   121 IPA(NT)=7
328   122 IF(NM.EQ.0) GOTO 124
329       DO 123 I=1,NM
330       NT=NT+1
331   123 IPA(NT)=9
332   124 IF(NZ.EQ.0) GOTO 130
333       DO 125 I=1,NZ
334       NT=NT+1
335   125 IPA(NT)=8
336   130 IF(NPRT(4))
337      *WRITE(NEWBCD,2004) NT,(IPA(I),I=1,20)
338       IF(IPA(1).EQ.7) NP=NP+1
339       IF(IPA(1).EQ.8) NZ=NZ+1
340       IF(IPA(1).EQ.9) NM=NM+1
341       GOTO 70
342   131 IF(NPRT(4))
343      *WRITE(NEWBCD,2005)
344 C
345 1001  FORMAT('0*CASPIP* CASCADE ENERGETICALLY NOT POSSIBLE',
346      $ ' CONTINUE WITH QUASI-ELASTIC SCATTERING')
347 1003  FORMAT(' *CASPIP* PION+ -INDUCED CASCADE,',
348      $ ' AVAIL. ENERGY',2X,F8.4,
349      $ 2X,'<NTOT>',2X,F8.4,2X,'FROM',4(2X,I3),2X,'PARTICLES')
350 1004  FORMAT(' *CASPIP* PION+ -INDUCED CASCADE,',
351      $ ' EXCLUSIVE REACTION NOT FOUND',
352      $ ' TRY ELASTIC SCATTERING  AVAIL. ENERGY',2X,F8.4,2X,
353      $ '<NTOT>',2X,F8.4)
354 2001  FORMAT('0*CASPIP* TABLES FOR MULT. DATA PION+  INDUCED REACTION',
355      $ ' FOR DEFINITION OF NUMBERS SEE FORTRAN CODING')
356 2002  FORMAT(' *CASPIP* TARGET PARTICLE FLAG',2X,I5)
357 2003  FORMAT(1H ,10E12.4)
358 2004  FORMAT(' *CASPIP* ',I3,2X,'PARTICLES , MASS INDEX ARRAY',2X,20I4)
359 2005  FORMAT(' *CASPIP* NO PARTICLES PRODUCED')
360 C
361  9999 CONTINUE
362       END