SUBROUTINE CASPIP(K,INT,NFL)
C
C *** CASCADE OF PI+ ***
C *** NVE 04-MAY-1988 CERN GENEVA ***
C
C ORIGIN : H.FESEFELDT (18-SEP-1987)
C
C PI+ UNDERGOES INTER ACTION WITH NUCLEON WITHIN NUCLEUS.
C CHECK IF ENERGETICALLY POSSIBLE TO PRODUCE PIONS/KAONS.
C IF NOT ASSUME NUCLEAR EXCITATION OCCURS AND INPUT PARTICLE
C IS DEGRADED IN ENERGY. NO OTHER PARTICLES PRODUCED.
C IF REACTION IS POSSIBLE FIND CORRECT NUMBER OF PIONS/PROTONS/
C NEUTRONS PRODUCED USING AN INTERPOLATION TO MULTIPLICITY DATA.
C REPLACE SOME PIONS OR PROTONS/NEUTRONS BY KAONS OR STRANGE BARYONS
C ACCORDING TO AVERAGE MULTIPLICITY PER INELASTIC REACTIONS.
C
#include "ghcdes/mxgkgh.inc"
#include "ghcdes/consts.inc"
#include "ghcdes/curpar.inc"
#include "ghcdes/result.inc"
#include "ghcdes/prntfl.inc"
#include "ghcdes/limits.inc"
#include "ghcdes/kginit.inc"
C
REAL N
DIMENSION PMUL(2,1200),ANORM(2,60),SUPP(10),CECH(10),B(2)
DIMENSION RNDM(1)
SAVE PMUL,ANORM
DATA SUPP/0.,0.2,0.45,0.55,0.65,0.75,0.85,0.90,0.94,0.98/
DATA CECH/0.33,0.27,0.29,0.31,0.27,0.18,0.13,0.10,0.09,0.07/
DATA B/0.7,0.7/,C/1.25/
C
INTEGER *4 NFLFORCE/0/
CHARACTER *(*) COMMAND
INTEGER *4 IDUMM
C
C --- INITIALIZATION INDICATED BY KGINIT(18) ---
IF (KGINIT(18) .NE. 0) GO TO 10
KGINIT(18)=1
C
C --- INITIALIZE PMUL AND ANORM ARRAYS ---
DO 9000 J=1,1200
DO 9001 I=1,2
PMUL(I,J)=0.0
IF (J .LE. 60) ANORM(I,J)=0.0
9001 CONTINUE
9000 CONTINUE
C
C** COMPUTE NORMALIZATION CONSTANTS
C** FOR P AS TARGET
C
L=0
DO 1 NP1=1,20
NP=NP1-1
NMM1=NP1-2
IF(NMM1.LE.1) NMM1=1
DO 1 NM1=NMM1,NP1
NM=NM1-1
DO 1 NZ1=1,20
NZ=NZ1-1
L=L+1
IF(L.GT.1200) GOTO 1
NT=NP+NM+NZ
IF(NT.LE.0.OR.NT.GT.60) GOTO 1
PMUL(1,L)=PMLTPC(NP,NM,NZ,NT,B(1),C)
ANORM(1,NT)=ANORM(1,NT)+PMUL(1,L)
1 CONTINUE
C** FOR N AS TARGET
L=0
DO 2 NP1=1,20
NP=NP1-1
NMM1=NP1-1
IF(NMM1.LE.1) NMM1=1
NPP1=NP1+1
DO 2 NM1=NMM1,NPP1
NM=NM1-1
DO 2 NZ1=1,20
NZ=NZ1-1
L=L+1
IF(L.GT.1200) GOTO 2
NT=NP+NM+NZ
IF(NT.LE.0.OR.NT.GT.60) GOTO 2
PMUL(2,L)=PMLTPC(NP,NM,NZ,NT,B(2),C)
ANORM(2,NT)=ANORM(2,NT)+PMUL(2,L)
2 CONTINUE
DO 3 I=1,60
IF(ANORM(1,I).GT.0.) ANORM(1,I)=1./ANORM(1,I)
IF(ANORM(2,I).GT.0.) ANORM(2,I)=1./ANORM(2,I)
3 CONTINUE
IF(.NOT.NPRT(10)) GOTO 10
WRITE(NEWBCD,2001)
DO 4 NFL=1,2
WRITE(NEWBCD,2002) NFL
WRITE(NEWBCD,2003) (ANORM(NFL,I),I=1,60)
WRITE(NEWBCD,2003) (PMUL(NFL,I),I=1,1200)
4 CONTINUE
C** CHOOSE PROTON OR NEUTRON AS TARGET
10 CONTINUE
IF(NFLFORCE.EQ.0) THEN
NFL=2
CALL GRNDM(RNDM,1)
IF(RNDM(1).LT.ZNO2/ATNO2) NFL=1
ELSE
NFL=NFLFORCE
ENDIF !NFLFORCE.EQ.0
TARMAS=RMASS(14)
IF (NFL .EQ. 2) TARMAS=RMASS(16)
S=AMASQ+TARMAS**2+2.0*TARMAS*EN
RS=SQRT(S)
ENP(8)=AMASQ+TARMAS**2+2.0*TARMAS*ENP(6)
ENP(9)=SQRT(ENP(8))
EAB=RS-TARMAS-RMASS(7)
C
C** ELASTIC SCATTERING
NP=0
NM=0
NZ=0
N=0.
IPA(1)=7
IPA(2)=14
IF(NFL.EQ.2) IPA(2)=16
IF(INT.EQ.2) GOTO 20
C** FOR PI+ N REACTIONS CHANGE SOME OF THE ELASTIC CROSS SECTION
C** TO PI+ N --> PI0 P
IF(NFL.EQ.1) GOTO 100
IPLAB=IFIX(P *5.)+1
IF(IPLAB.GT.10) IPLAB=10
CALL GRNDM(RNDM,1)
IF(RNDM(1).GT.CECH(IPLAB)/ATNO2**0.42) GOTO 100
IPA(1)=8
IPA(2)=14
GOTO 100
C** CHECK IF ENERGETICALLY POSSIBLE TO PRODUCE ONE EXTRA PION IN REACT.
20 IF (EAB .LE. RMASS(7)) GOTO 55
C** SUPPRESSION OF HIGH MULTIPLICITY EVENTS AT LOW MOMENTUM
IEAB=IFIX(EAB*5.)+1
IF(IEAB.GT.10) GOTO 22
CALL GRNDM(RNDM,1)
IF(RNDM(1).LT.SUPP(IEAB)) GOTO 22
N=1.
GOTO (23,24),NFL
23 CONTINUE
TEST=-(1+B(1))**2/(2.0*C**2)
IF (TEST .LE. EXPXL) TEST=EXPXL
IF (TEST .GE. EXPXU) TEST=EXPXU
W0=EXP(TEST)
WP=EXP(TEST)
CALL GRNDM(RNDM,1)
RAN=RNDM(1)
NP=0
NM=0
NZ=1
IF(RAN.LT.W0/(W0+WP)) GOTO 50
NP=1
NM=0
NZ=0
GOTO 50
24 CONTINUE
TEST=-(1+B(2))**2/(2.0*C**2)
IF (TEST .LE. EXPXL) TEST=EXPXL
IF (TEST .GE. EXPXU) TEST=EXPXU
W0=EXP(TEST)
WP=EXP(TEST)
TEST=-(-1+B(2))**2/(2.0*C**2)
IF (TEST .LE. EXPXL) TEST=EXPXL
IF (TEST .GE. EXPXU) TEST=EXPXU
WM=EXP(TEST)
WT=W0+WP+WM
WP=W0+WP
CALL GRNDM(RNDM,1)
RAN=RNDM(1)
NP=0
NM=0
NZ=1
IF(RAN.LT.W0/WT) GOTO 50
NP=1
NM=0
NZ=0
IF(RAN.LT.WP/WT) GOTO 50
NP=0
NM=1
NZ=0
GOTO 50
C
22 ALEAB=LOG(EAB)
C** NO. OF TOTAL PARTICLES VS SQRT(S)-2*MP
N=3.62567+0.665843*ALEAB+0.336514*ALEAB*ALEAB
* +0.117712*ALEAB*ALEAB*ALEAB+0.0136912*ALEAB*ALEAB*ALEAB*ALEAB
N=N-2.
C** NORMALIZATION CONSTANT FOR KNO-DISTRIBUTION
ANPN=0.
DO 21 NT=1,60
TEST=-(PI/4.0)*(NT/N)**2
IF (TEST .LE. EXPXL) TEST=EXPXL
IF (TEST .GE. EXPXU) TEST=EXPXU
DUM1=PI*NT/(2.0*N*N)
DUM2=ABS(DUM1)
DUM3=EXP(TEST)
ADDNVE=0.0
IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GE. 1.0E-10)) ADDNVE=DUM1*DUM3
ANPN=ANPN+ADDNVE
21 CONTINUE
ANPN=1./ANPN
C** P OR N AS TARGET
CALL GRNDM(RNDM,1)
RAN=RNDM(1)
EXCS=0.
GOTO (30,40),NFL
C** FOR P AS TARGET
30 L=0
DO 31 NP1=1,20
NP=NP1-1
NMM1=NP1-2
IF(NMM1.LE.1) NMM1=1
DO 31 NM1=NMM1,NP1
NM=NM1-1
DO 31 NZ1=1,20
NZ=NZ1-1
L=L+1
IF(L.GT.1200) GOTO 31
NT=NP+NM+NZ
IF(NT.LE.0.OR.NT.GT.60) GOTO 31
TEST=-(PI/4.0)*(NT/N)**2
IF (TEST .LE. EXPXL) TEST=EXPXL
IF (TEST .GE. EXPXU) TEST=EXPXU
DUM1=ANPN*PI*NT*PMUL(1,L)*ANORM(1,NT)/(2.0*N*N)
DUM2=ABS(DUM1)
DUM3=EXP(TEST)
ADDNVE=0.0
IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GE. 1.0E-10)) ADDNVE=DUM1*DUM3
EXCS=EXCS+ADDNVE
IF(RAN.LT.EXCS) GOTO 50
31 CONTINUE
GOTO 80
C** FOR N AS TARGET
40 L=0
DO 41 NP1=1,20
NP=NP1-1
NMM1=NP1-1
IF(NMM1.LE.1) NMM1=1
NPP1=NP1+1
DO 41 NM1=NMM1,NPP1
NM=NM1-1
DO 41 NZ1=1,20
NZ=NZ1-1
L=L+1
IF(L.GT.1200) GOTO 41
NT=NP+NM+NZ
IF(NT.LE.0.OR.NT.GT.60) GOTO 41
TEST=-(PI/4.0)*(NT/N)**2
IF (TEST .LE. EXPXL) TEST=EXPXL
IF (TEST .GE. EXPXU) TEST=EXPXU
DUM1=ANPN*PI*NT*PMUL(2,L)*ANORM(2,NT)/(2.0*N*N)
DUM2=ABS(DUM1)
DUM3=EXP(TEST)
ADDNVE=0.0
IF (DUM2 .GE. 1.0) ADDNVE=DUM1*DUM3
IF ((DUM2 .LT. 1.0) .AND. (DUM3 .GT. 1.0E-10)) ADDNVE=DUM1*DUM3
EXCS=EXCS+ADDNVE
IF(RAN.LT.EXCS) GOTO 50
41 CONTINUE
GOTO 80
50 GOTO (60,65),NFL
60 IF(NP.EQ.1+NM) GOTO 61
IF(NP.EQ.2+NM) GOTO 63
IPA(1)=7
IPA(2)=14
GOTO 100
61 CALL GRNDM(RNDM,1)
IF(RNDM(1).LT.0.5) GOTO 62
IPA(1)=7
IPA(2)=16
GOTO 100
62 IPA(1)=8
IPA(2)=14
GOTO 100
63 IPA(1)=8
IPA(2)=16
GOTO 100
65 IF(NP.EQ.NM) GOTO 66
IF(NP.EQ.1+NM) GOTO 68
IPA(1)=7
IPA(2)=14
GOTO 100
66 CALL GRNDM(RNDM,1)
IF(RNDM(1).LT.0.25) GOTO 67
IPA(1)=7
IPA(2)=16
GOTO 100
67 IPA(1)=8
IPA(2)=14
GOTO 100
68 IPA(1)=8
IPA(2)=16
GOTO 100
70 IF(NPRT(4))
*WRITE(NEWBCD,1003) EAB,N,NFL,NP,NM,NZ
CALL STPAIR
IF(INT.EQ.1) CALL TWOB(7,NFL,N)
IF(INT.EQ.2) CALL GENXPT(7,NFL,N)
GO TO 9999
55 IF(NPRT(4))
*WRITE(NEWBCD,1001)
GOTO 53
C** EXCLUSIVE REACTION NOT FOUND
80 IF(NPRT(4))
*WRITE(NEWBCD,1004) RS,N
53 INT=1
NP=0
NM=0
NZ=0
N=0.
IPA(1)=7
IPA(2)=14
IF(NFL.EQ.2) IPA(2)=16
100 DO 101 I=3,60
101 IPA(I)=0
IF(INT.LE.0) GOTO 131
120 NT=2
IF(NP.EQ.0) GOTO 122
DO 121 I=1,NP
NT=NT+1
121 IPA(NT)=7
122 IF(NM.EQ.0) GOTO 124
DO 123 I=1,NM
NT=NT+1
123 IPA(NT)=9
124 IF(NZ.EQ.0) GOTO 130
DO 125 I=1,NZ
NT=NT+1
125 IPA(NT)=8
130 IF(NPRT(4))
*WRITE(NEWBCD,2004) NT,(IPA(I),I=1,20)
IF(IPA(1).EQ.7) NP=NP+1
IF(IPA(1).EQ.8) NZ=NZ+1
IF(IPA(1).EQ.9) NM=NM+1
GOTO 70
131 IF(NPRT(4))
*WRITE(NEWBCD,2005)
C
1001 FORMAT('0*CASPIP* CASCADE ENERGETICALLY NOT POSSIBLE',
$ ' CONTINUE WITH QUASI-ELASTIC SCATTERING')
1003 FORMAT(' *CASPIP* PION+ -INDUCED CASCADE,',
$ ' AVAIL. ENERGY',2X,F8.4,
$ 2X,'<NTOT>',2X,F8.4,2X,'FROM',4(2X,I3),2X,'PARTICLES')
1004 FORMAT(' *CASPIP* PION+ -INDUCED CASCADE,',
$ ' EXCLUSIVE REACTION NOT FOUND',
$ ' TRY ELASTIC SCATTERING AVAIL. ENERGY',2X,F8.4,2X,
$ '<NTOT>',2X,F8.4)
2001 FORMAT('0*CASPIP* TABLES FOR MULT. DATA PION+ INDUCED REACTION',
$ ' FOR DEFINITION OF NUMBERS SEE FORTRAN CODING')
2002 FORMAT(' *CASPIP* TARGET PARTICLE FLAG',2X,I5)
2003 FORMAT(1H ,10E12.4)
2004 FORMAT(' *CASPIP* ',I3,2X,'PARTICLES , MASS INDEX ARRAY',2X,20I4)
2005 FORMAT(' *CASPIP* NO PARTICLES PRODUCED')
C
9999 CONTINUE
RETURN
C
ENTRY CASPIPSET(COMMAND,IDUMM)
IF(COMMAND.EQ.'SET:NFLFORCE') THEN
NFLFORCE=IDUMM
ENDIF !COMMAND.EQ.'SET:NFLFORCE'
RETURN
END