[u/mrichter/AliRoot.git] / GEANT321 / gheisha / stpair.F

*
* $Id$
*
* $Log$
* Revision 1.1.1.1  1995/10/24 10:20:59  cernlib
* Geant
*
*
#include "geant321/pilot.h"
*CMZ :  3.21/02 29/03/94  15.41.38  by  S.Giani
*-- Author :
      SUBROUTINE STPAIR
C
C *** STRANGE PARTICLE PAIR PRODUCTION ***
C *** NVE 14-MAR-1988 CERN GENEVA ***
C
C ORIGIN : H.FESEFELDT 16-DEC-1987
C
C THE SAME FORMULA FOR <K KB> VS AVAILABLE ENERGY
C                  AND <K Y>  VS AVAILABLE ENERGY
C FOR ALL REACTIONS.
C CHOOSE CHARGE COMBINATIONS K+ K- , K+ K0B, K0 K0B OR K0 K-
C                            K+ Y0, K0 Y+, K0 Y-
C FOR ANTIBARYON INDUCED REACTIONS HALF OF THE CROSS SECTIONS
C KB YB PAIRS ARE PRODUCED
C CHARGE IS NOT CONSERVED , NO EXPERIMENTAL DATA AVAILABLE FOR
C EXCLUSIVE REACTIONS, THEREFORE SOME AVERAGE BEHAVIOUR ASSUMED.
C THE RATIO L/SIGMA IS TAKEN AS 3:1 (FROM EXPERIMENTAL LOW ENERGY)
C
#include "geant321/s_defcom.inc"
C
      REAL KKB,KY
      DIMENSION KKB(9),KY(12),IPAKKB(2,9),IPAKY(2,12),IPAKYB(2,12)
      DIMENSION AVKKB(12),AVKY(12),AVNNB(12),AVRS(12)
      DIMENSION RNDM(1)
      DATA KKB/0.2500,0.3750,0.5000,0.5625,0.6250,0.6875,0.7500,
     *         0.8750,1.000/
      DATA KY /0.200,0.300,0.400,0.550,0.625,0.700,0.800,0.850,
     *         0.900,0.950,0.975,1.000/
      DATA IPAKKB/10,13,10,11,10,12,11,11,11,12,12,11,12,12,
     *            11,13,12,13/
      DATA IPAKY /18,10,18,11,18,12,20,10,20,11,20,12,21,10,
     *            21,11,21,12,22,10,22,11,22,12/
      DATA IPAKYB/19,13,19,12,19,11,23,13,23,12,23,11,24,13,
     *            24,12,24,11,25,13,25,12,25,11/
      DATA AVRS/3.,4.,5.,6.,7.,8.,9.,10.,20.,30.,40.,50./
      DATA AVKKB/0.0015,0.005,0.012,0.0285,0.0525,0.075,0.0975,
     *           0.123,0.28,0.398,0.495,0.573/
      DATA AVKY /0.005,0.03,0.064,0.095,0.115,0.13,0.145,0.155,
     *           0.20,0.205,0.210,0.212/
      DATA AVNNB/0.00001,0.0001,0.0006,0.0025,0.01,0.02,0.04,
     $           0.05,0.12,0.15,0.18,0.20/
C
      IF(IPA(3).LE.0) GO TO 9999
      IER(50)=IER(50)+1
      IPA1=ABS(IPA(1))
      IPA2=ABS(IPA(2))
C --- PROTECTION AGAINST ANNIHILATION PROCESSES ---
      IF ((IPA1 .EQ. 0) .OR. (IPA2 .EQ. 0)) GO TO 9999
      EAB=RS-ABS(RMASS(IPA1))-ABS(RMASS(IPA2))
      IF(EAB.LT.1.) GO TO 9999
C**
C** CHOOSE RANDOM REPLACEMENT OF PRODUCED KAONS (16.12.87)
      DO 111 I=1,60
      IF(IPA(I).EQ.0) GOTO 112
  111 CONTINUE
  112 I=I-3
      CALL GRNDM(RNDM,1)
      I3=3+IFIX(RNDM(1)*I)
  114 CALL GRNDM(RNDM,1)
      I4=3+IFIX(RNDM(1)*I)
      IF(I.EQ.1) I4=4
      IF(I3.EQ.I4) GOTO 114
C
C *** CHOOSE RANDOM REPLACEMENT OF PRODUCED KAONS (16.12.87) ***
C --- GET RS BIN ---
      DO 1 I=2,12
      IF (RS .LE. AVRS(I)) GO TO 2
 1    CONTINUE
      I1=11
      I2=12
      GO TO 3
C
 2    CONTINUE
      I1=I-1
      I2=I
C
C *** USE LINEAR INTERPOLATION OR EXTRAPOLATION BY Y=RC*X+B ***
 3    CONTINUE
      DXNVE=AVRS(I2)-AVRS(I1)
      DYNVE=LOG(AVKKB(I2))-LOG(AVKKB(I1))
      RCNVE=DYNVE/DXNVE
      BNVE=LOG(AVKKB(I1))-RCNVE*AVRS(I1)
      AVK=RCNVE*RS+BNVE
      DYNVE=LOG(AVKY(I2))-LOG(AVKY(I1))
      RCNVE=DYNVE/DXNVE
      BNVE=LOG(AVKY(I1))-RCNVE*AVRS(I1)
      AVY=RCNVE*RS+BNVE
      DYNVE=LOG(AVNNB(I2))-LOG(AVNNB(I1))
      RCNVE=DYNVE/DXNVE
      BNVE =LOG(AVNNB(I1))-RCNVE*AVRS(I1)
      AVN  =RCNVE*RS+BNVE
C
      AVK=EXP(AVK)
      AVY=EXP(AVY)
      AVN=EXP(AVN)
      IF(AVK+AVY+AVN.LE.0.) GOTO 9999
      IF(IPA1.LT.14) AVY=AVY/2.
      IF(IPA2.LT.14) AVY=0.
      AVY=AVY+AVK+AVN
      AVK=    AVK+AVN
      CALL GRNDM(RNDM,1)
      RAN=RNDM(1)
      IF(RAN.LT.AVN) GOTO 5
      IF(RAN.LT.AVK) GOTO 10
      IF(RAN.LT.AVY) GOTO 20
      GO TO 9999
    5 IF((EAB-2.).LT.0.) GO TO 9999
      CALL GRNDM(RNDM,1)
      IF(RNDM(1).LT.0.5) GO TO 6
      IPA(I3)=14
      IPA(I4)=15
      GOTO 30
    6 IPA(I3)=16
      IPA(I4)=17
      GOTO 30
   10 IF((EAB-1.).LT.0.) GO TO 9999
      CALL GRNDM(RNDM,1)
      RAN=RNDM(1)
      DO 11 I=1,9
      IF(RAN.LT.KKB(I)) GOTO 12
   11 CONTINUE
      GO TO 9999
   12 IPA(I3)=IPAKKB(1,I)
      IPA(I4)=IPAKKB(2,I)
      GOTO 30
   20 IF((EAB-1.6).LT.0.) GO TO 9999
      CALL GRNDM(RNDM,1)
      RAN=RNDM(1)
      DO 21 I=1,12
      IF(RAN.LT.KY(I)) GOTO 22
   21 CONTINUE
      GO TO 9999
   22 IF(IPA(1).LT.14) GOTO 23
      CALL GRNDM(RNDM,1)
      IF(RNDM(1).LT.0.5) GOTO 23
      IPA1=ABS(IPA(1))
      IPA(1)=IPAKY(1,I)
      IF(IPA1.EQ.15) GOTO 25
      IF(IPA1.EQ.17) GOTO 25
      IF(IPA1.EQ.19) GOTO 25
      IF(IPA1.GT.22) GOTO 25
      GOTO 24
   25 IPA(1)=IPAKYB(1,I)
      IPA(I3)=IPAKYB(2,I)
      GOTO 30
   23 IPA(2)=IPAKY(1,I)
   24 IPA(I3)=IPAKY(2,I)
C** CHECK THE AVAILABLE ENERGY
   30 EAB=RS
      IJ=0
      DO 31 I=1,60
      IF(IPA(I).EQ.0) GOTO 31
      IPA1=ABS(IPA(I))
      EAB=EAB-ABS(RMASS(IPA1))
      IJ=IJ+1
      IF(EAB.LT.0.) GOTO 35
   31 CONTINUE
      IF (NPRT(4)) WRITE(NEWBCD,1003) (IPA(J),J=1,IJ)
      GO TO 9999
   35 I=I-1
      L=I-1
      IF(L.LE.0) GO TO 9999
      DO 36 J=I,60
   36 IPA(J)=0
      IF (NPRT(4)) WRITE(NEWBCD,1002) (IPA(J),J=1,L)
C
 1002 FORMAT(' *STPAIR* KKB/KY PAIR PRODUCTION NOT ENOUGH ENERGY',
     $ ' REDUCE NUMBER OF PARTICLES ',2X,20I3)
 1003 FORMAT(' *STPAIR* KKB/KY PAIR PRODUCTION ENERGY SUFFICIENT',
     $ ' NUMBER OF PARTICLES ',2X,20I3)
C
 9999 CONTINUE
C
      END
Commit	Line	Data
fe4da5cc	1	*
	2	* $Id$
	3	*
	4	* $Log$
	5	* Revision 1.1.1.1 1995/10/24 10:20:59 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 STPAIR
	13	C
	14	C * STRANGE PARTICLE PAIR PRODUCTION *
	15	C * NVE 14-MAR-1988 CERN GENEVA *
	16	C
	17	C ORIGIN : H.FESEFELDT 16-DEC-1987
	18	C
	19	C THE SAME FORMULA FOR <K KB> VS AVAILABLE ENERGY
	20	C AND <K Y> VS AVAILABLE ENERGY
	21	C FOR ALL REACTIONS.
	22	C CHOOSE CHARGE COMBINATIONS K+ K- , K+ K0B, K0 K0B OR K0 K-
	23	C K+ Y0, K0 Y+, K0 Y-
	24	C FOR ANTIBARYON INDUCED REACTIONS HALF OF THE CROSS SECTIONS
	25	C KB YB PAIRS ARE PRODUCED
	26	C CHARGE IS NOT CONSERVED , NO EXPERIMENTAL DATA AVAILABLE FOR
	27	C EXCLUSIVE REACTIONS, THEREFORE SOME AVERAGE BEHAVIOUR ASSUMED.
	28	C THE RATIO L/SIGMA IS TAKEN AS 3:1 (FROM EXPERIMENTAL LOW ENERGY)
	29	C
	30	#include "geant321/s_defcom.inc"
	31	C
	32	REAL KKB,KY
	33	DIMENSION KKB(9),KY(12),IPAKKB(2,9),IPAKY(2,12),IPAKYB(2,12)
	34	DIMENSION AVKKB(12),AVKY(12),AVNNB(12),AVRS(12)
	35	DIMENSION RNDM(1)
	36	DATA KKB/0.2500,0.3750,0.5000,0.5625,0.6250,0.6875,0.7500,
	37	* 0.8750,1.000/
	38	DATA KY /0.200,0.300,0.400,0.550,0.625,0.700,0.800,0.850,
	39	* 0.900,0.950,0.975,1.000/
	40	DATA IPAKKB/10,13,10,11,10,12,11,11,11,12,12,11,12,12,
	41	* 11,13,12,13/
	42	DATA IPAKY /18,10,18,11,18,12,20,10,20,11,20,12,21,10,
	43	* 21,11,21,12,22,10,22,11,22,12/
	44	DATA IPAKYB/19,13,19,12,19,11,23,13,23,12,23,11,24,13,
	45	* 24,12,24,11,25,13,25,12,25,11/
	46	DATA AVRS/3.,4.,5.,6.,7.,8.,9.,10.,20.,30.,40.,50./
	47	DATA AVKKB/0.0015,0.005,0.012,0.0285,0.0525,0.075,0.0975,
	48	* 0.123,0.28,0.398,0.495,0.573/
	49	DATA AVKY /0.005,0.03,0.064,0.095,0.115,0.13,0.145,0.155,
	50	* 0.20,0.205,0.210,0.212/
	51	DATA AVNNB/0.00001,0.0001,0.0006,0.0025,0.01,0.02,0.04,
	52	$ 0.05,0.12,0.15,0.18,0.20/
	53	C
	54	IF(IPA(3).LE.0) GO TO 9999
	55	IER(50)=IER(50)+1
	56	IPA1=ABS(IPA(1))
	57	IPA2=ABS(IPA(2))
	58	C --- PROTECTION AGAINST ANNIHILATION PROCESSES ---
	59	IF ((IPA1 .EQ. 0) .OR. (IPA2 .EQ. 0)) GO TO 9999
	60	EAB=RS-ABS(RMASS(IPA1))-ABS(RMASS(IPA2))
	61	IF(EAB.LT.1.) GO TO 9999
	62	C**
	63	C** CHOOSE RANDOM REPLACEMENT OF PRODUCED KAONS (16.12.87)
	64	DO 111 I=1,60
65	IF(IPA(I).EQ.0) GOTO 112
66	111 CONTINUE
67	112 I=I-3
68	CALL GRNDM(RNDM,1)
69	I3=3+IFIX(RNDM(1)*I)
70	114 CALL GRNDM(RNDM,1)
71	I4=3+IFIX(RNDM(1)*I)
72	IF(I.EQ.1) I4=4
73	IF(I3.EQ.I4) GOTO 114
74	C
75	C * CHOOSE RANDOM REPLACEMENT OF PRODUCED KAONS (16.12.87) *
76	C --- GET RS BIN ---
77	DO 1 I=2,12
78	IF (RS .LE. AVRS(I)) GO TO 2
79	1 CONTINUE
80	I1=11
81	I2=12
82	GO TO 3
83	C
84	2 CONTINUE
85	I1=I-1
86	I2=I
87	C
88	C *** USE LINEAR INTERPOLATION OR EXTRAPOLATION BY Y=RCX+B **
89	3 CONTINUE
90	DXNVE=AVRS(I2)-AVRS(I1)
91	DYNVE=LOG(AVKKB(I2))-LOG(AVKKB(I1))
92	RCNVE=DYNVE/DXNVE
93	BNVE=LOG(AVKKB(I1))-RCNVE*AVRS(I1)
94	AVK=RCNVE*RS+BNVE
95	DYNVE=LOG(AVKY(I2))-LOG(AVKY(I1))
96	RCNVE=DYNVE/DXNVE
97	BNVE=LOG(AVKY(I1))-RCNVE*AVRS(I1)
98	AVY=RCNVE*RS+BNVE
99	DYNVE=LOG(AVNNB(I2))-LOG(AVNNB(I1))
100	RCNVE=DYNVE/DXNVE
101	BNVE =LOG(AVNNB(I1))-RCNVE*AVRS(I1)
102	AVN =RCNVE*RS+BNVE
103	C
104	AVK=EXP(AVK)
105	AVY=EXP(AVY)
106	AVN=EXP(AVN)
107	IF(AVK+AVY+AVN.LE.0.) GOTO 9999
108	IF(IPA1.LT.14) AVY=AVY/2.
109	IF(IPA2.LT.14) AVY=0.
110	AVY=AVY+AVK+AVN
111	AVK= AVK+AVN
112	CALL GRNDM(RNDM,1)
113	RAN=RNDM(1)
114	IF(RAN.LT.AVN) GOTO 5
115	IF(RAN.LT.AVK) GOTO 10
116	IF(RAN.LT.AVY) GOTO 20
117	GO TO 9999
118	5 IF((EAB-2.).LT.0.) GO TO 9999
119	CALL GRNDM(RNDM,1)
120	IF(RNDM(1).LT.0.5) GO TO 6
121	IPA(I3)=14
122	IPA(I4)=15
123	GOTO 30
124	6 IPA(I3)=16
125	IPA(I4)=17
126	GOTO 30
127	10 IF((EAB-1.).LT.0.) GO TO 9999
128	CALL GRNDM(RNDM,1)
129	RAN=RNDM(1)
130	DO 11 I=1,9
131	IF(RAN.LT.KKB(I)) GOTO 12
132	11 CONTINUE
133	GO TO 9999
134	12 IPA(I3)=IPAKKB(1,I)
135	IPA(I4)=IPAKKB(2,I)
136	GOTO 30
137	20 IF((EAB-1.6).LT.0.) GO TO 9999
138	CALL GRNDM(RNDM,1)
139	RAN=RNDM(1)
140	DO 21 I=1,12
141	IF(RAN.LT.KY(I)) GOTO 22
142	21 CONTINUE
143	GO TO 9999
144	22 IF(IPA(1).LT.14) GOTO 23
145	CALL GRNDM(RNDM,1)
146	IF(RNDM(1).LT.0.5) GOTO 23
147	IPA1=ABS(IPA(1))
148	IPA(1)=IPAKY(1,I)
149	IF(IPA1.EQ.15) GOTO 25
150	IF(IPA1.EQ.17) GOTO 25
151	IF(IPA1.EQ.19) GOTO 25
152	IF(IPA1.GT.22) GOTO 25
153	GOTO 24
154	25 IPA(1)=IPAKYB(1,I)
155	IPA(I3)=IPAKYB(2,I)
156	GOTO 30
157	23 IPA(2)=IPAKY(1,I)
158	24 IPA(I3)=IPAKY(2,I)
159	C** CHECK THE AVAILABLE ENERGY
160	30 EAB=RS
161	IJ=0
162	DO 31 I=1,60
163	IF(IPA(I).EQ.0) GOTO 31
164	IPA1=ABS(IPA(I))
165	EAB=EAB-ABS(RMASS(IPA1))
166	IJ=IJ+1
167	IF(EAB.LT.0.) GOTO 35
168	31 CONTINUE
169	IF (NPRT(4)) WRITE(NEWBCD,1003) (IPA(J),J=1,IJ)
170	GO TO 9999
171	35 I=I-1
172	L=I-1
173	IF(L.LE.0) GO TO 9999
174	DO 36 J=I,60
175	36 IPA(J)=0
176	IF (NPRT(4)) WRITE(NEWBCD,1002) (IPA(J),J=1,L)
177	C
178	1002 FORMAT(' STPAIR KKB/KY PAIR PRODUCTION NOT ENOUGH ENERGY',
179	$ ' REDUCE NUMBER OF PARTICLES ',2X,20I3)
180	1003 FORMAT(' STPAIR KKB/KY PAIR PRODUCTION ENERGY SUFFICIENT',
181	$ ' NUMBER OF PARTICLES ',2X,20I3)
182	C
183	9999 CONTINUE
184	C
185	END