1 | SUBROUTINE HDPM
|
---|
2 |
|
---|
3 | C-----------------------------------------------------------------------
|
---|
4 | C H(ADRONIC) D(UAL) P(ARTON) M(ODEL)
|
---|
5 | C
|
---|
6 | C GENERATOR OF HADRONIC COLLISION INSPIRED BY DUAL PARTON MODEL
|
---|
7 | C THIS SUBROUTINE IS CALLED FROM SDPM
|
---|
8 | C-----------------------------------------------------------------------
|
---|
9 |
|
---|
10 | IMPLICIT DOUBLE PRECISION (A-H,O-Z)
|
---|
11 | *KEEP,CONST.
|
---|
12 | COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
|
---|
13 | DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
|
---|
14 | *KEEP,DPMFLG.
|
---|
15 | COMMON /DPMFLG/ NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
|
---|
16 | INTEGER NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
|
---|
17 | *KEEP,ELADPM.
|
---|
18 | COMMON /ELADPM/ ELMEAN,ELMEAA,IELDPM,IELDPA
|
---|
19 | DOUBLE PRECISION ELMEAN(37),ELMEAA(37)
|
---|
20 | INTEGER IELDPM(37,13),IELDPA(37,13)
|
---|
21 | *KEEP,ELASTY.
|
---|
22 | COMMON /ELASTY/ ELAST,IELIS,IELHM,IELNU,IELPI
|
---|
23 | DOUBLE PRECISION ELAST
|
---|
24 | INTEGER IELIS(20),IELHM(20),IELNU(20),IELPI(20)
|
---|
25 | *KEEP,INDICE.
|
---|
26 | COMMON /INDICE/ NNUCN,NKA0,NHYPN,NETA,NETAS,NPIZER,
|
---|
27 | * NNC,NKC,NHC,NPC,NCH,NNN,NKN,NHN,NET,NPN
|
---|
28 | INTEGER NNUCN(2:3),NKA0(2:3),NHYPN(2:3),NETA(2:3,1:4),
|
---|
29 | * NETAS(2:3),NPIZER(2:3),
|
---|
30 | * NNC,NKC,NHC,NPC,NCH,NNN,NKN,NHN,NET,NPN
|
---|
31 | *KEEP,INTER.
|
---|
32 | COMMON /INTER/ AVCH,AVCH3,DC0,DLOG,DMLOG,ECMDIF,ECMDPM,ELAB,
|
---|
33 | * FNEUT,FNEUT2,GNU,PLAB,POSC2,POSC3,POSN2,POSN3,
|
---|
34 | * RC3TO2,S,SEUGF,SEUGP,SLOG,SLOGSQ,SMLOG,
|
---|
35 | * WIDC2,WIDC3,WIDN2,WIDN3,YCM,YY0,ZN,
|
---|
36 | * IDIF,ITAR
|
---|
37 | DOUBLE PRECISION AVCH,AVCH3,DC0,DLOG,DMLOG,ECMDIF,ECMDPM,ELAB,
|
---|
38 | * FNEUT,FNEUT2,GNU,PLAB,POSC2,POSC3,POSN2,POSN3,
|
---|
39 | * RC3TO2,S,SEUGF,SEUGP,SLOG,SLOGSQ,SMLOG,
|
---|
40 | * WIDC2,WIDC3,WIDN2,WIDN3,YCM,YY0,ZN
|
---|
41 | INTEGER IDIF,ITAR
|
---|
42 | *KEEP,ISTA.
|
---|
43 | COMMON /ISTA/ IFINET,IFINNU,IFINKA,IFINPI,IFINHY
|
---|
44 | INTEGER IFINET,IFINNU,IFINKA,IFINPI,IFINHY
|
---|
45 | *KEEP,LEPAR.
|
---|
46 | COMMON /LEPAR/ LEPAR1,LEPAR2,LASTPI,NRESPC,NRESPN,NCPLUS
|
---|
47 | INTEGER LEPAR1,LEPAR2,LASTPI,NRESPC,NRESPN,NCPLUS
|
---|
48 | *KEEP,MULT.
|
---|
49 | COMMON /MULT/ EKINL,MSMM,MULTMA,MULTOT
|
---|
50 | DOUBLE PRECISION EKINL
|
---|
51 | INTEGER MSMM,MULTMA(37,13),MULTOT(37,13)
|
---|
52 | *KEEP,NEWPAR.
|
---|
53 | COMMON /NEWPAR/ EA,PT2,PX,PY,TMAS,YR,ITYP,
|
---|
54 | * IA1,IA2,IB1,IB2,IC1,IC2,ID1,ID2,IE1,IE2,IF1,IF2,
|
---|
55 | * IG1,IG2,IH1,IH2,II1,II2,IJ1,NTOT
|
---|
56 | DOUBLE PRECISION EA(3000),PT2(3000),PX(3000),PY(3000),TMAS(3000),
|
---|
57 | * YR(3000)
|
---|
58 | INTEGER ITYP(3000),
|
---|
59 | * IA1,IA2,IB1,IB2,IC1,IC2,ID1,ID2,IE1,IE2,IF1,IF2,
|
---|
60 | * IG1,IG2,IH1,IH2,II1,II2,IJ1,NTOT
|
---|
61 | *KEEP,PAM.
|
---|
62 | COMMON /PAM/ PAMA,SIGNUM
|
---|
63 | DOUBLE PRECISION PAMA(6000),SIGNUM(6000)
|
---|
64 | *KEEP,PARPAR.
|
---|
65 | COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
|
---|
66 | * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
|
---|
67 | DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
|
---|
68 | * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
|
---|
69 | INTEGER ITYPE,LEVL
|
---|
70 | *KEEP,PARPAE.
|
---|
71 | DOUBLE PRECISION GAMMA,COSTHE,PHI,H,T,X,Y,CHI,BETA,GCM,ECM
|
---|
72 | EQUIVALENCE (CURPAR(2),GAMMA), (CURPAR(3),COSTHE),
|
---|
73 | * (CURPAR(4), PHI ), (CURPAR(5), H ),
|
---|
74 | * (CURPAR(6), T ), (CURPAR(7), X ),
|
---|
75 | * (CURPAR(8), Y ), (CURPAR(9), CHI ),
|
---|
76 | * (CURPAR(10),BETA), (CURPAR(11),GCM ),
|
---|
77 | * (CURPAR(12),ECM )
|
---|
78 | *KEEP,RANDPA.
|
---|
79 | COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
|
---|
80 | DOUBLE PRECISION FAC,U1,U2
|
---|
81 | REAL RD(3000)
|
---|
82 | INTEGER ISEED(103,10),NSEQ
|
---|
83 | LOGICAL KNOR
|
---|
84 | *KEEP,RATIOS.
|
---|
85 | COMMON /RATIOS/ RPI0R,RPIER,RPEKR,RPEKNR,PPICH,PPINCH,PPNKCH,
|
---|
86 | * ISEL,NEUTOT,NTOTEM
|
---|
87 | DOUBLE PRECISION RPI0R,RPIER,RPEKR,RPEKNR,PPICH,PPINCH,PPNKCH
|
---|
88 | INTEGER ISEL,NEUTOT,NTOTEM
|
---|
89 | *KEEP,RESON.
|
---|
90 | COMMON /RESON/ RDRES,RESRAN,IRESPAR
|
---|
91 | REAL RDRES(2),RESRAN(1000)
|
---|
92 | INTEGER IRESPAR
|
---|
93 |
|
---|
94 | *KEEP,REST.
|
---|
95 | COMMON /REST/ CONTNE,TAR,LT
|
---|
96 | DOUBLE PRECISION CONTNE(3),TAR
|
---|
97 | INTEGER LT
|
---|
98 | *KEEP,RUNPAR.
|
---|
99 | COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
|
---|
100 | * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
|
---|
101 | * MONIOU,MDEBUG,NUCNUC,
|
---|
102 | * CETAPE,
|
---|
103 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
---|
104 | * N1STTR,MDBASE,
|
---|
105 | * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
---|
106 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
---|
107 | * ,GHEISH,GHESIG
|
---|
108 | COMMON /RUNPAC/ DSN,HOST,USER
|
---|
109 | DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
|
---|
110 | REAL STEPFC
|
---|
111 | INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
|
---|
112 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
---|
113 | * N1STTR,MDBASE
|
---|
114 | INTEGER CETAPE
|
---|
115 | CHARACTER*79 DSN
|
---|
116 | CHARACTER*20 HOST,USER
|
---|
117 |
|
---|
118 | LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
---|
119 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
---|
120 | * ,GHEISH,GHESIG
|
---|
121 | *KEND.
|
---|
122 |
|
---|
123 | C-----------------------------------------------------------------------
|
---|
124 |
|
---|
125 | IF ( DEBUG ) WRITE(MDEBUG,444) (CURPAR(I),I=1,9)
|
---|
126 | 444 FORMAT(' HDPM : CURPAR=',1P,9E10.3)
|
---|
127 |
|
---|
128 | C SET ANTI-LEADER TO PROTON OR NEUTRON; TARGET IS ALWAYS NUCLEON
|
---|
129 | CALL RMMAR( RD,1,1 )
|
---|
130 | IF ( RD(1) .LE. CONTNE(LT) ) THEN
|
---|
131 | ITAR = 13
|
---|
132 | ELSE
|
---|
133 | ITAR = 14
|
---|
134 | ENDIF
|
---|
135 |
|
---|
136 |
|
---|
137 | C CALCULATE LAB AND CM ENERGY
|
---|
138 | IF ( ITYPE .NE. 1 ) THEN
|
---|
139 | ELAB = PAMA(ITYPE) * GAMMA
|
---|
140 | PLAB = ELAB * BETA
|
---|
141 | S = PAMA(ITYPE)**2 + PAMA(ITAR)**2 + 2.D0*PAMA(ITAR)*ELAB
|
---|
142 | ELSE
|
---|
143 | C FOR GAMMA-INDUCED REACTION TAKE PI(0) AS LEADING PARTICLE
|
---|
144 | ITYPE = 7
|
---|
145 | ELAB = GAMMA
|
---|
146 | PLAB = ELAB
|
---|
147 | S = PAMA(ITAR)**2 + 2.D0*PAMA(ITAR)*ELAB
|
---|
148 | ENDIF
|
---|
149 |
|
---|
150 | ECMDPM = SQRT(S)
|
---|
151 | IF ( DEBUG ) WRITE(MDEBUG,*) 'HDPM : ITYPE,ELAB,PLAB,S,ECMDPM=',
|
---|
152 | * ITYPE,SNGL(ELAB),SNGL(PLAB),SNGL(S),SNGL(ECMDPM)
|
---|
153 |
|
---|
154 | C LN(S), LN(S)**2 AND RAPIDITY OF C. M. SYSTEM IN LAB
|
---|
155 | SLOG = LOG(S)
|
---|
156 | SLOGSQ = SLOG**2
|
---|
157 | SMLOG = LOG( 2.D0 * PAMA(ITAR) * ELAB )
|
---|
158 | ELABLG = LOG(ELAB)
|
---|
159 | EPLUSP = ELAB + PLAB
|
---|
160 | * YCM = 0.5D0 * LOG( (ELAB+PAMA(ITAR)+PLAB)/(ELAB+PAMA(ITAR)-PLAB) )
|
---|
161 | YCM = 0.5D0 * LOG( (EPLUSP**2 +PAMA(ITAR)*EPLUSP)/
|
---|
162 | * (PAMA(ITYPE)**2+PAMA(ITAR)*EPLUSP) )
|
---|
163 | IF ( DEBUG ) WRITE(MDEBUG,*)'HDPM : SLOG,SLOGSQ,YCM=',
|
---|
164 | * SNGL(SLOG),SNGL(SLOGSQ),SNGL(YCM)
|
---|
165 |
|
---|
166 | C-----------------------------------------------------------------------
|
---|
167 | C RETURN POINT IF CALCULATION OF PARTICLES GOES WRONG
|
---|
168 | 1 CONTINUE
|
---|
169 |
|
---|
170 | IF ( ITYPE .NE. 7 ) THEN
|
---|
171 | C CHOOSE NUMBER OF INTERACTIONS IN TARGET
|
---|
172 | CALL TARINT
|
---|
173 | ELSE
|
---|
174 | C FOR GAMMA-INDUCED REACTIONS TAKE ALWAYS ONE COLLISION
|
---|
175 | GNU = 1.D0
|
---|
176 | ENDIF
|
---|
177 |
|
---|
178 | C-----------------------------------------------------------------------
|
---|
179 | C NO DIFFRACTION IF
|
---|
180 | C OR THE NUMBER OF INTERACTIONS IN TARGET IS CHOSEN RANDOMLY
|
---|
181 | C AND MORE THAN ONE INTERACTION TAKES PLACE
|
---|
182 | C OR PRIMARY PARTICLE IS GAMMA (PI0)
|
---|
183 | C NOW NFLDIF DECIDES WHETHER DIFFRACTIVE PROCESS POSSIBLE OR NOT
|
---|
184 | IF ( ( NFLAIN.EQ.0 .AND. GNU.GT.1.D0 .AND. NFLDIF.EQ.0 )
|
---|
185 | * .OR. ( ITYPE .EQ. 7 ) ) THEN
|
---|
186 | IDIF = 0
|
---|
187 | ELSE
|
---|
188 | C SET DIFFRACTION FLAG IF RANDOM NUMBER < PROBABILITY
|
---|
189 | CALL RMMAR( RD,1,1 )
|
---|
190 | C IDIF IS 0 : NO DIFFRACTION ; IDIF IS 1 : DIFFRACTION
|
---|
191 | C DIFFRACTION RISES WITH ENERGY AND SATURATES AT 10000 GEV
|
---|
192 | C ### DAS TUT ES ABER NICHT: ES IST KONSTANT 0.15 (SIEHE DPFUNC) !!!!
|
---|
193 | IF ( RD(1) .GT. DPFUNC(ECMDPM) ) THEN
|
---|
194 | IDIF = 0
|
---|
195 | ELSE
|
---|
196 | IDIF = 1
|
---|
197 | ENDIF
|
---|
198 | ENDIF
|
---|
199 |
|
---|
200 |
|
---|
201 | C PRINTOUT FOR DEBUG
|
---|
202 | IF ( DEBUG ) THEN
|
---|
203 | WRITE(MDEBUG,*) ' DIFFRACTIVE INTERACTION (0/1) = ',IDIF
|
---|
204 | ENDIF
|
---|
205 |
|
---|
206 | C SET COUNTER FOR REPEAT TO 0
|
---|
207 | NREPRD = 0
|
---|
208 |
|
---|
209 | C GENERATION OF INTERACTION
|
---|
210 | 1919 CONTINUE
|
---|
211 |
|
---|
212 | C FLAG TO CHECK NUMBER OF SECONDARIES;
|
---|
213 | C IS CHANGED TO 1 IF SECONDARY MULTIPLICITY IS LOW
|
---|
214 | ISEL = 0
|
---|
215 | C SET LEADING PARTICLE TO INCOMING PARTICLE AND ANTI-LEADER TO NUCLEON
|
---|
216 | C (AS IT COMES FROM TARGET NUCLEUS) BOTH MAY BE CHANGED BY LEPACX
|
---|
217 | LEPAR1 = ITYPE
|
---|
218 | LEPAR2 = ITAR
|
---|
219 |
|
---|
220 | IF ( IDIF .EQ. 0 ) THEN
|
---|
221 | C-----------------------------------------------------------------------
|
---|
222 | C NON SINGLE DIFFRACTIVE PROCESS STARTS HERE
|
---|
223 |
|
---|
224 | CALL NSD
|
---|
225 | C CHARGE EXCHANGE ENABLED? EXCHANGE LEADER AND ANTI-LEADER
|
---|
226 | LASTPI = 0
|
---|
227 | NRESPC = 0
|
---|
228 | NRESPN = 0
|
---|
229 | NCPLUS = 0
|
---|
230 | IF ( NFLCHE .EQ. 0 ) THEN
|
---|
231 | CALL LEPACX( ECMDPM,ELABLG,LEPAR1,1 )
|
---|
232 | CALL LEPACX( ECMDPM,ELABLG,LEPAR2,2 )
|
---|
233 | ENDIF
|
---|
234 | 1921 CONTINUE
|
---|
235 | CALL RNEGBI( NCH,AVCH,ECMDPM )
|
---|
236 | C NCH IS # OF ALL CHARGED PARTICLES INCLUDING EXCESS FROM TARGET
|
---|
237 | IF ( NCH .LT. 1 ) THEN
|
---|
238 | IF ( LEPAR1 .LT. 50 .OR. LEPAR2 .LT. 50 ) THEN
|
---|
239 | NREPRD = NREPRD + 1
|
---|
240 | IF ( NREPRD .GT. 10 ) GOTO 1
|
---|
241 | GOTO 1921
|
---|
242 | ELSE
|
---|
243 | C INTERACTION IS ONLY RESONANCE PRODUCTION
|
---|
244 | ISEL = 1
|
---|
245 | ENDIF
|
---|
246 | ENDIF
|
---|
247 | C WIDTH PLATEAU FOR CLUSTERS AND FOR CALCULATION OF CENTR.RAP.DENSITY
|
---|
248 | DELRAP = 0.6722D0 * (2.95D0 + 0.0302D0 * SLOG)
|
---|
249 | C SET RSLOG FOR CALCULATION OF PARTICLE RATIOS
|
---|
250 | RSLOG = SLOG
|
---|
251 | C AVERAGE TRANSVERSE MOMENTUM
|
---|
252 | CALL AVEPT( ECMDPM,SLOG )
|
---|
253 |
|
---|
254 | ELSE
|
---|
255 | C-----------------------------------------------------------------------
|
---|
256 | C SINGLE DIFFRACTIVE PROCESS STARTS HERE
|
---|
257 |
|
---|
258 | 1920 CONTINUE
|
---|
259 | CALL DIFRAC( NRETDF )
|
---|
260 | IF ( NRETDF .EQ. 1 ) GOTO 1
|
---|
261 | C CHARGE EXCHANGE ENABLED? EXCHANGE CHARGE OF DIFFRACTING PARTICLE
|
---|
262 | LASTPI = 0
|
---|
263 | NRESPC = 0
|
---|
264 | NRESPN = 0
|
---|
265 | NCPLUS = 0
|
---|
266 | IF ( NFLCHE .EQ. 0 ) THEN
|
---|
267 | IF ( YY0 .GT. 0.D0 ) THEN
|
---|
268 | C PROJECTILE DIFFRACTION
|
---|
269 | CALL LEPACX( ECMDIF,DMLOG,LEPAR1,1 )
|
---|
270 | ELSE
|
---|
271 | C TARGET DIFFRACTION
|
---|
272 | CALL LEPACX( ECMDIF,DMLOG,LEPAR2,2 )
|
---|
273 | ENDIF
|
---|
274 | ENDIF
|
---|
275 | C FLUCTUATION OF MULTIPLICITY ACCORDING TO NEG.BIN. DISTRIBUTION
|
---|
276 | CALL RNEGBI( NCH,AVCH,ECMDIF )
|
---|
277 | C REPEAT CALCULATION AS SOMETHING WENT WRONG
|
---|
278 | IF ( NCH .LT. 1 ) THEN
|
---|
279 | IF ( (YY0 .GT. 0.D0 .AND. LEPAR1 .LT. 50) .OR.
|
---|
280 | * (YY0 .LT. 0.D0 .AND. LEPAR2 .LT. 50) ) THEN
|
---|
281 | NREPRD = NREPRD + 1
|
---|
282 | IF ( NREPRD .GT. 10 ) GOTO 1
|
---|
283 | GOTO 1920
|
---|
284 | ELSE
|
---|
285 | C DIFFRACTIVE INTERACTION IS ONLY RESONANCE PRODUCTION
|
---|
286 | ISEL = 1
|
---|
287 | ENDIF
|
---|
288 | ENDIF
|
---|
289 | C SET RSLOG FOR CALCULATION OF PARTICLE RATIOS
|
---|
290 | RSLOG = DLOG
|
---|
291 | C HERE WE USE ECMDPM, BECAUSE THE MOMENTUM TRANSFER IS DEPENDENT
|
---|
292 | C ON THE ENERGY OF THE TOTAL SYSTEM AND NOT ON THE DIFFRACTING MASS
|
---|
293 | CALL AVEPT( ECMDPM,SLOG )
|
---|
294 |
|
---|
295 | ENDIF
|
---|
296 |
|
---|
297 | C-----------------------------------------------------------------------
|
---|
298 | C NOW FOR NSD AND DIFFRACTIVE PROCESSES
|
---|
299 |
|
---|
300 | C IN CASE OF LOW MULTIPLICITY SET FLAG ISEL
|
---|
301 | IF ( NCH .LE. 2 ) ISEL=1
|
---|
302 | C FNCH IS FLUCTUATING TOT.NUMBER OF CHARGED PARTICLES FOR ALL 3 STRINGS
|
---|
303 | FNCH = DBLE(NCH)
|
---|
304 | C RATIO ALL CHARGED PARTICLES WITH FLUCTUATION/WITHOUT FLUCTUATION
|
---|
305 | XZ = FNCH / AVCH
|
---|
306 | C FNCH3 IS FLUCTUATING NUMBER OF CHARGED PARTICLES FOR 3RD STRING
|
---|
307 | FNCH3 = XZ * AVCH3
|
---|
308 | C FNCH2 IS FLUCTUATING NUMBER OF CHARGED PARTICLES 1ST AND 2ND STRING
|
---|
309 | FNCH2 = FNCH - FNCH3
|
---|
310 | C RC3TO2 IS RATIO (CHARGED 3RD STRING)/(CHARGED 1ST AND 2ND STRING)
|
---|
311 | IF ( FNCH2 .NE. 0.D0 ) THEN
|
---|
312 | RC3TO2 = FNCH3 / FNCH2
|
---|
313 | ELSE
|
---|
314 | RC3TO2 = 0.D0
|
---|
315 | ENDIF
|
---|
316 | IF ( DEBUG ) WRITE(MDEBUG,*) ' FNCH,FNCH2,FNCH3,RC3TO2=',
|
---|
317 | * SNGL(FNCH),SNGL(FNCH2),SNGL(FNCH3),SNGL(RC3TO2)
|
---|
318 |
|
---|
319 | C IS NUMBER OF NEUTRALS FLUCTUATING AS NUMBER OF CHARGED ?
|
---|
320 | IF ( NFLPIF .EQ. 0 .OR. IDIF .EQ. 1 .OR. ECMDPM .LT. 60.D0 ) THEN
|
---|
321 | C SET NUMBER OF GAMMAS ACCORDING TO NEG. BIN. VARIABLE XZ
|
---|
322 | C AS NUMBER OF NEUTRALS FLUCTUATES AS CHARGED.
|
---|
323 | SEUGF = SEUGP * XZ
|
---|
324 | ZG = XZ
|
---|
325 | ELSE
|
---|
326 | C NFLPIF IS 1 MEANS: # OF PI(0) FLUCTUATES AS MEASURED AT COLLIDER
|
---|
327 | IF ( ECMDPM .LT. 200.D0 ) THEN
|
---|
328 | SEUGF = SEUGP * XZ
|
---|
329 | * SEUGF = (0.0786D0*SLOG-0.010D0)*FNCH2 + (0.391D0*SLOG+0.305D0)
|
---|
330 | ELSE
|
---|
331 | C DETERMINE NEW NUMBER OF GAMMAS WITH FLUCTUATION AROUND SEUGP*XZ
|
---|
332 | AGR = EXP(-XZ)
|
---|
333 | DGR = SEUGP * XZ * (0.9823D0 - 0.3756D0 * AGR)
|
---|
334 | SGS = DGR * (0.14718D0 + 2.53213D0 * AGR)
|
---|
335 | 723 CONTINUE
|
---|
336 | SEUGF = 0.88D0 * RANNOR(DGR,SGS)
|
---|
337 | IF ( SEUGF .LT. 1.D0 ) GOTO 723
|
---|
338 | ENDIF
|
---|
339 | C SET NEGATIVE BINOMIAL VARIABLE ZG FOR GAMMAS
|
---|
340 | ZG = SEUGF / SEUGP
|
---|
341 | ENDIF
|
---|
342 | SEUGF = MAX( 1.D0, SEUGF )
|
---|
343 | IF ( DEBUG ) WRITE(MDEBUG,*) 'HDPM :XZ,ZG,SEUGF=',
|
---|
344 | * SNGL(XZ),SNGL(ZG),SNGL(SEUGF)
|
---|
345 |
|
---|
346 | C-----------------------------------------------------------------------
|
---|
347 | C RATIO ALL-NUCLEON/ALL-CHARGED
|
---|
348 | C PARAMETRISATION FROM UA5, NUCL. PHYS. B291 (1987) 445 EQ.(2.4)
|
---|
349 | RNUCCH = MAX( 0.D0, -0.008D0 + 0.00865D0 * RSLOG )
|
---|
350 | C NUMBER FOR DIRECT NEUTRON/ANTINEUTRON PRODUCTION 1ST AND 2ND STRING
|
---|
351 | C MULTIPLY BY 0.5 BECAUSE RATIO RNUCCH GIVES (ALL-NUCL)/(ALL-CHARGED)
|
---|
352 | C AND HERE ONLY THE NEUTRON-ANTINEUTRONS ARE COUNTED
|
---|
353 | FNUCN = 0.5D0 * RNUCCH * FNCH2
|
---|
354 | C RATIO (ALL CHARGED SIGMAS)/(ALL CHARGED) IS 1/3 OF ALL STRANGE BARYON
|
---|
355 | C PARAMETRISATION FORM UA5, NUCL. PHYS. B291 (1987) 445 EQ.(2.5)
|
---|
356 | RHYPCH = MAX( 0.D0, (-0.007D0 + 0.0028D0 * RSLOG) * OB3 )
|
---|
357 | C NEUTRAL STRANGE BARYONS ARE DOUBLE OF CHARGED STRANGE BARYONS
|
---|
358 | FHYPN = 2.D0 * RHYPCH * FNCH2
|
---|
359 | C CORRECT NUMBER OF GAMMAS FROM NEUTRAL HYPERON DECAY S0-->L+GAMMA
|
---|
360 | SEUGFC = MAX( 0.D0, SEUGF - 0.5D0 * FHYPN )
|
---|
361 | C RATIO CHARGED-KAON/CHARGED PIONS
|
---|
362 | C PARAMETRISATION FROM UA5, NUCL. PHYS. B291 (1987) 445 EQ.(2.7)
|
---|
363 | RKPI = MAX (0.D0, 0.024D0 + 0.0062D0 * RSLOG )
|
---|
364 | C RKCH IS RATIO (CHARGED-KAON)/(ALL-CHARGED) DERIVED FROM RKPI;
|
---|
365 | C THE FACTOR 0.5 IN FRONT OF RNUCCH IS BECAUSE ONLY HALF OF NUCLEONS
|
---|
366 | C ARE P/PBAR. THE 1.17 IS AN APROXIMATE CONVERSION FACTOR FROM
|
---|
367 | C (ALL-NUCL)/(ALL-CHARGED) TO (ALL-NUCL)/(CHARGED-PI), WHICH IS A BIT
|
---|
368 | C ENERGY DEPENDENT (1.14 ...1.21) SEE GEICH-GIMBEL TABLE 7.1
|
---|
369 | RKCH = RKPI / (1.D0 + RKPI + (0.5D0*RNUCCH+RHYPCH) * 1.17D0)
|
---|
370 | C K0/K0-BAR FOR 1ST AND 2ND STRING
|
---|
371 | C NEUTRAL KAONS ARE PRODUCED WITH THE SAME RATE AS CHARGED KAONS
|
---|
372 | FKA0 = RKCH * FNCH2
|
---|
373 | C RATIO ETA/PI(0) IS ASSUMED TO BE INDEPENDENT OF ENERGY = 0.19
|
---|
374 | C SEE: ANSORGE ET AL. (UA5-COLLABORATION) Z.PHYS.C43(1989)75
|
---|
375 | * RETPI0 = 0.19D0
|
---|
376 | C RATIO ETA/PI(0) IS ASSUMED TO BE DEPENDENT ON ENERGY
|
---|
377 | C SEE: GEICH-GIMBEL,INT.J.MOD.PHYS.A4(1989)1527 TAB.7.1
|
---|
378 | C HECK'S FIT: RETPI0 IS 0.06 + 0.006*RSLOG + 0.0011*RSLOG**2
|
---|
379 | RETPI0 = 0.06D0 + 0.006D0 * RSLOG + 0.0011D0 * RSLOG**2
|
---|
380 | C AUXIL1 IS FRACTION OF PI(0)/(PI(0)+ETA)
|
---|
381 | AUXIL1 = 1.D0 / (1.D0 + RETPI0)
|
---|
382 | C NUMBER OF GAMMAS FROM PI(0) IS 2, FROM ETA IS 3.216 IN AVERAGE;
|
---|
383 | C AUXIL2 IS NUMBER OF GAMMA-PRODUCING-PARTICLES: PI(0) AND ETA
|
---|
384 | AUXIL2 = SEUGFC / ( AUXIL1 * 2.D0 + (1.D0 - AUXIL1) * 3.216D0 )
|
---|
385 | FETA = (1.D0 - AUXIL1) * AUXIL2
|
---|
386 | FPI0 = AUXIL1 * AUXIL2
|
---|
387 | C CORRECT FPI0 BY DECAYS OF STRANGE BARYONS; NEUTRAL: FHYPN*0.357
|
---|
388 | C CHARGED: 0.5*FNCH2*RHYPCH*0.5157; IT YIELDS FHYPN*(0.357+0.12893)
|
---|
389 | FPI0 = MAX( 0.D0, FPI0 - FHYPN * 0.486D0 )
|
---|
390 | C SUMMED NEUTRAL PARTICLES FOR 1ST AND 2ND STRING
|
---|
391 | FNEUT2 = FNUCN + FKA0 + FHYPN + FETA + FPI0
|
---|
392 | C NEUTRAL PARTICLES FROM 3RD STRING
|
---|
393 | FNEUT3 = RC3TO2 * FNEUT2
|
---|
394 | C TOTAL NUMBER OF NEUTRALS
|
---|
395 | FNEUT = FNEUT2 + FNEUT3
|
---|
396 | NEUTOT = NINT( FNEUT )
|
---|
397 | C CALCULATE TOTAL NUMBER OF PARTICLES TO BE CREATED
|
---|
398 | NTOTEM = NCH + NEUTOT
|
---|
399 | IF ( DEBUG ) WRITE(MDEBUG,*)
|
---|
400 | * ' FNUCN,FKA0,FHYPN,FETA,FPI0,FNEUT2,FNEUT3,NTOTEM=',
|
---|
401 | * SNGL(FNUCN),SNGL(FKA0),SNGL(FHYPN),SNGL(FETA),SNGL(FPI0),
|
---|
402 | * SNGL(FNEUT2),SNGL(FNEUT3),NTOTEM
|
---|
403 | C LIMIT OF SECONDARIES PRODUCED (GIVEN BY SIZE OF ARRAY : 3000)
|
---|
404 | C LIMIT IS ARRAY SIZE - SIZE OF LARGEST TARGET NUCLEUS(40)
|
---|
405 | IF ( NTOTEM .GE. 2956 ) THEN
|
---|
406 | WRITE(MONIOU,*) 'HDPM : REJECT EVENT WITH ',NTOTEM,
|
---|
407 | * ' SECONDARIES'
|
---|
408 | GOTO 1
|
---|
409 | ENDIF
|
---|
410 | C SPECIAL TREATMENT IF MULTIPLICITY IS TOO LOW
|
---|
411 | IF ( NTOTEM .LE. 3 ) ISEL = 1
|
---|
412 |
|
---|
413 | C FRACTION OF THE VARIOUS NEUTRAL PARTICLES (NN, K(0), L+S0 AS PAIRS)
|
---|
414 | C NORMALIZE WITH THE SUM OF ALL NEUTRAL PARTICLES
|
---|
415 | FNORML = 1.D0 / ( 0.5D0 * (FNUCN+FKA0+FHYPN) + FETA + FPI0 )
|
---|
416 | RNUCNR = FNUCN * FNORML * 0.5D0
|
---|
417 | RKA0R = FKA0 * FNORML * 0.5D0
|
---|
418 | RHYPNR = FHYPN * FNORML * 0.5D0
|
---|
419 | RETAR = FETA * FNORML
|
---|
420 | RPI0R = FPI0 * FNORML
|
---|
421 | C CUMULATED RATIOS (NN, K(0), LAMBDA+SIGMA0 AS PAIRS)
|
---|
422 | RPIER = RPI0R + RETAR
|
---|
423 | RPEKR = RPIER + RKA0R
|
---|
424 | RPEKNR = RPEKR + RNUCNR
|
---|
425 | C THEN THE REMAINDER (1-RPEKNR) MUST BE NEUTRAL HYPERON PAIRS
|
---|
426 | IF ( DEBUG ) WRITE(MDEBUG,*)
|
---|
427 | * ' RPI0R,RETAR,RKA0R,RNUCNR,RHYPNR,FNORML=',
|
---|
428 | * SNGL(RPI0R),SNGL(RETAR),SNGL(RKA0R),SNGL(RNUCNR),SNGL(RHYPNR),
|
---|
429 | * SNGL(FNORML)
|
---|
430 |
|
---|
431 | C PROBABILITY TO PRODUCE CHARGED PIONS IS PROBABILITY NOT TO PRODUCE
|
---|
432 | C CHARGED KAONS OR PROTONS OR CHARGED HYPERONS, WHERE PROTON/ANTIPROTON
|
---|
433 | C IS HALF OF (ALL-NUCL)/(ALL-CHARGED)
|
---|
434 | AUXIL = RKCH + 0.5D0 * RNUCCH + RHYPCH
|
---|
435 | AUXIL3 = 1.D0 - AUXIL
|
---|
436 | C RENORMALIZATION AS P/P_BAR, K+-, AND HYPERONS ARE PRODUCED IN PAIRS
|
---|
437 | C AUXIL2 IS INVERSE OF NORMALISATION
|
---|
438 | AUXIL2 = 1.D0 / (1.D0 - 0.5D0 * AUXIL)
|
---|
439 | C CUMULATED PROBABILITIES (PP, K+-, SIGMA+- AS PAIRS)
|
---|
440 | PPICH = AUXIL3 * AUXIL2
|
---|
441 | PPINCH = PPICH + 0.25D0 * RNUCCH * AUXIL2
|
---|
442 | PPNKCH = PPINCH + 0.5D0 * RKCH * AUXIL2
|
---|
443 | C THEN THE REMAINDER (1-PPNKCH) MUST BE CHARGED HYPERON PAIRS
|
---|
444 | IF ( DEBUG ) WRITE(MDEBUG,*) ' PPICH,PPINCH,PPNKCH=',
|
---|
445 | * SNGL(PPICH),SNGL(PPINCH),SNGL(PPNKCH)
|
---|
446 |
|
---|
447 | C NOW SELECT HOW MANY PARTICLES OF EACH TYPE ARE PRODUCED
|
---|
448 | CALL PARNUM( INUMFL )
|
---|
449 | IF ( INUMFL .NE. 0 ) GOTO 1919
|
---|
450 |
|
---|
451 | C DEFINE PARTICLE NUMBERS WHERE SPECIAL RAPIDITY IS CALCULATED
|
---|
452 | C FOR PARTICLES FROM TARGET (THIRD STRING)
|
---|
453 | PPP = RC3TO2 / (1.D0+RC3TO2)
|
---|
454 | C NUMBER OF PARTICLES IN PROTON ANTIPROTON PAIRS FROM TARGET
|
---|
455 | ITA = NINT(PPP * 2.D0 * NNC)
|
---|
456 | C NUMBER OF PARTICLES IN K+ K- PAIRS FROM TARGET
|
---|
457 | ITB = NINT(PPP * 2.D0 * NKC)
|
---|
458 | C NUMBER OF PARTICLES IN SIGMA+ SIGMA- PAIRS FROM TARGET
|
---|
459 | ITC = NINT(PPP * 2.D0 * NHC)
|
---|
460 | C NUMBER OF PI+ PI- FROM TARGET
|
---|
461 | ITD = NINT(PPP * NPC )
|
---|
462 | C CALCULATE BOUNDARIES
|
---|
463 | IA1 = 2
|
---|
464 | IA2 = IA1 + ITA
|
---|
465 | IB1 = IA1 + 2 * NNC
|
---|
466 | IB2 = IB1 + ITB
|
---|
467 | IC1 = IB1 + 2 * NKC
|
---|
468 | IC2 = IC1 + ITC
|
---|
469 | ID1 = IC1 + 2 * NHC
|
---|
470 | ID2 = ID1 + ITD
|
---|
471 | IE1 = ID1 + NPC
|
---|
472 | C NUMBER OF PARTICLES IN NEUTRON ANTINEUTRON PAIRS FROM TARGET
|
---|
473 | IE2 = IE1 + 2 * NNUCN(3)
|
---|
474 | IF1 = IE1 + 2 * NNN
|
---|
475 | C NUMBER OF PARTICLES IN K0S K0L PAIRS FROM TARGET
|
---|
476 | IF2 = IF1 + 2 * NKA0(3)
|
---|
477 | IG1 = IF1 + 2 * NKN
|
---|
478 | C NUMBER OF PARTICLES IN NEUTRAL HYPERON PAIRS FROM TARGET
|
---|
479 | IG2 = IG1 + 2 * NHYPN(3)
|
---|
480 | IH1 = IG1 + 2 * NHN
|
---|
481 | C NUMBER OF ETA FROM TARGET
|
---|
482 | IH2 = IH1 + NETAS(3)
|
---|
483 | II1 = IH1 + NET
|
---|
484 | C NUMBER OF PI(0) FROM TARGET
|
---|
485 | II2 = II1 + NPIZER(3)
|
---|
486 | IJ1 = II1 + NPN
|
---|
487 | IF ( DEBUG ) THEN
|
---|
488 | WRITE(MDEBUG,*) ' CHARGED FROM TARGET:',ITA,ITB,ITC,ITD
|
---|
489 | WRITE(MDEBUG,*) ' NEUTRAL FROM TARGET:',
|
---|
490 | * 2*NNUCN(3),2*NKA0(3),2*NHYPN(3),NETAS(3),NPIZER(3)
|
---|
491 | WRITE(MDEBUG,*) ' NTOTEM,IJ1=',NTOTEM,IJ1
|
---|
492 | ENDIF
|
---|
493 | C REDEFINE TOTAL NUMBER OF SECONDARY PARTICLES : NTOTEM
|
---|
494 | C BY CHARGE EXCHANGE AND RESONANCE FORMATION THIS NUMBER MAY BE ALTERED
|
---|
495 | NTOTEM = IJ1 - 2
|
---|
496 |
|
---|
497 | C-----------------------------------------------------------------------
|
---|
498 | C RATIO OF RAPIDITY DENSITY TO MEAN PSEUDORAPIDITY IN CENTER
|
---|
499 | C PARAMETRISATION SEE CAPDEVIELLE, J.PHYS.G:NUCL.PHYS.15(1989)909,EQ.6
|
---|
500 | IF ( XZ .LT. 1.5D0 ) THEN
|
---|
501 | RDS = (0.24396D0 + 0.70150424D0 * XZ)**2
|
---|
502 | ELSE
|
---|
503 | RDS = (0.55685D0 + 0.48664753D0 * XZ)**2
|
---|
504 | ENDIF
|
---|
505 | C CALCULATE NOW: DN/DY AT Y = 0; DC0 IS AVERAGE PSEUDORAPIDITY DENSITY
|
---|
506 | C TRAP IS RATIO (RAPID.DENS.)/(PSEUDORAP.DENS.) IN CENTER OF RAPIDITY
|
---|
507 | TRAP = 1.25D0
|
---|
508 | IF ( IDIF .EQ. 0 .AND. ECMDPM .GT. 19.4D0 )
|
---|
509 | * TRAP = MAX( 1.D0, 1.28852D0 - 0.0065D0 * SMLOG )
|
---|
510 | DCN2 = DC0 * RDS * TRAP
|
---|
511 | IF ( DEBUG ) WRITE(MDEBUG,*) ' RDS,TRAP,DCN2=',
|
---|
512 | * SNGL(RDS),SNGL(TRAP),SNGL(DCN2)
|
---|
513 | C AMPLITUDE OF GAUSSIAN 1ST AND 2ND STRING
|
---|
514 | ATG2 = FNCH2 / (5.0132566D0 * WIDC2)
|
---|
515 | C NEW DEFINITION OF POSITION BASED ON SEMI INCLUSIVE DATA
|
---|
516 | SQ2 = 2.D0 * ATG2 / DCN2
|
---|
517 | C FINAL POSITION OF GAUSSIAN; WIDTH WIDC2 IS UNCHANGED
|
---|
518 | IF ( SQ2 .GT. 1.D0 ) POSC2 = WIDC2 * SQRT( 2.D0*LOG(SQ2) )
|
---|
519 | C DENSITY OF CHARGED IN EXCESS FROM TARGET IN CENTER OF RAPIDITY
|
---|
520 | DCN3 = 0.5D0 * (GNU - 1.D0) * DCN2
|
---|
521 | IF (DEBUG) WRITE(MDEBUG,*) ' SQ2,POSC2,DCN3=',
|
---|
522 | * SNGL(SQ2),SNGL(POSC2),SNGL(DCN3)
|
---|
523 | IF ( DCN3 .GT. 0.D0 ) THEN
|
---|
524 | C AMPLITUDE 3RD GAUSSIAN
|
---|
525 | ATG3 = FNCH3 / (5.0132566D0 * WIDC3)
|
---|
526 | C AMPLITUDE IS DIVIDED BY DENSITY FOR GETTING CENTER OF 3RD GAUSSIAN
|
---|
527 | SQ3 = 2.D0 * ATG3 / DCN3
|
---|
528 | C CHECK IF ADDITIVE MULTIPLICITY IS TOO LOW
|
---|
529 | IF ( SQ3 .GT. 1.D0 ) POSC3 = WIDC3 * SQRT( 2.D0*LOG(SQ3) )
|
---|
530 | IF (DEBUG) WRITE(MDEBUG,*)' SQ3,POSC3=',SNGL(SQ3),SNGL(POSC3)
|
---|
531 | ENDIF
|
---|
532 |
|
---|
533 | C NFLPI0 .EQ. 0 MEANS TREAT PI(0) RAPIDITY ACCORDING TO COLLIDER DATA
|
---|
534 | IF ( NFLPI0 .EQ. 0 ) THEN
|
---|
535 | C RATIO OF RAPIDITY DENSITY TO MEAN PSEUDORAPIDITY AT CENTER WITH Z<1.5
|
---|
536 | IF ( ZG .LT. 1.5D0 ) THEN
|
---|
537 | RDG = (0.24396D0 + 0.70150424D0 * ZG)**2
|
---|
538 | ELSE
|
---|
539 | RDG = (0.55685D0 + 0.48664753D0 * ZG)**2
|
---|
540 | ENDIF
|
---|
541 | C GAMMAS USE RATIO TRAG TO CALCULATE RATIO OF RAPIDITY TO
|
---|
542 | C PSEUDO RAPIDITY DENSITY IN CENTER (TRAG = 1.1 * 0.5 ).
|
---|
543 | C FACTOR 0.5 COMES FROM RATIO NEUTRAL/CHARGED, AS WE USE DC0, WHICH
|
---|
544 | C IS AVERAGE PSEUDORAPIDITY DENSITY FOR CHARGED PIONS
|
---|
545 | TRAG = 0.55D0
|
---|
546 | IF ( IDIF .EQ. 0 ) THEN
|
---|
547 | IF ( ECMDPM .GT. 19.4D0 )
|
---|
548 | * TRAG = MAX( 0.4D0, 0.6658D0 - 0.01954D0 * SMLOG )
|
---|
549 | IF ( ECMDPM .LE. 50.D0 ) THEN
|
---|
550 | DCG = DC0 * RDG * TRAG
|
---|
551 | ELSEIF ( ECMDPM .LE. 200.D0 ) THEN
|
---|
552 | DCG = DC0 * RDG * TRAG * (1.D0 + 0.18D0 * LOG(ECMDPM/50.D0))
|
---|
553 | ELSE
|
---|
554 | DCG = DC0 * RDG * TRAG * 1.25D0
|
---|
555 | ENDIF
|
---|
556 | ELSE
|
---|
557 | DCG = DC0 * RDG * TRAG
|
---|
558 | ENDIF
|
---|
559 | C DEFINE WIDTH OF STRINGS FOR NEUTRAL PIONS AND ETAS
|
---|
560 | WIDN2 = WIDC2 * MIN( 1.D0, 1.12275D0 - 0.0208D0 * RSLOG )
|
---|
561 | C NEW DEFINITION OF CENTER OF GAUSSIAN BASED ON SEMI INCLUSIVE DATA
|
---|
562 | C USING AMPLITUDE OF THE GAUSSIAN FOR NEUTRALS
|
---|
563 | AUXIL = 2.D0 / (5.0132566D0 * WIDN2 * DCG)
|
---|
564 | C TOTAL MULTIPLICITY USED FOR 1ST AND 2ND STRING OF PI(0) AND ETA
|
---|
565 | C IS GIVEN BY THEIR NUMBERS. ANALOGOUS FOR 3RD STRING
|
---|
566 | SP2 = DBLE ( NPIZER(2)+NETAS(2)) * AUXIL
|
---|
567 | C FINAL CENTER OF GAUSSIANS FOR PI(0) AND ETA (WIDC2 IS UNCHANGED)
|
---|
568 | IF ( SP2 .GT. 1.D0 ) THEN
|
---|
569 | POSN2 = WIDN2 * SQRT( 2.D0 * LOG(SP2) )
|
---|
570 | ELSE
|
---|
571 | POSN2 = POSC2
|
---|
572 | ENDIF
|
---|
573 | WIDN3 = WIDN2
|
---|
574 | SP3 = DBLE(NPIZER(3)+NETAS(3)) * AUXIL
|
---|
575 | IF ( SP3 .GT. 1.D0 ) THEN
|
---|
576 | POSN3 = WIDN3 * SQRT( 2.D0 * LOG(SP3) )
|
---|
577 | ELSE
|
---|
578 | POSN3 = POSC3
|
---|
579 | ENDIF
|
---|
580 | ELSE
|
---|
581 | C NFLPI0 .EQ. 1 MEANS RAPIDITY OF PI(0) AND ETA SAME AS THAT OF CHARGED
|
---|
582 | POSN2 = POSC2
|
---|
583 | WIDN2 = WIDC2
|
---|
584 | POSN3 = POSC3
|
---|
585 | WIDN3 = WIDC3
|
---|
586 | ENDIF
|
---|
587 | IF ( DEBUG ) WRITE(MDEBUG,*)
|
---|
588 | * ' ZG,RDG,DCG,SP2,SP3,POSN2,POSN3,WIDN2 =',
|
---|
589 | * SNGL(ZG),SNGL(RDG),SNGL(DCG),SNGL(SP2),SNGL(SP3),SNGL(POSN2),
|
---|
590 | * SNGL(POSN3),SNGL(WIDN2)
|
---|
591 |
|
---|
592 | C-----------------------------------------------------------------------
|
---|
593 | NREPR1 = 0
|
---|
594 | C RETURN POINT. NUMBERS OF PARTICLES REMAIN UNCHANGED FOR NEXT TRY,
|
---|
595 | C BUT INDIVIDUAL RAPIDITIES GET NEW VALUES.
|
---|
596 | C START FROM BEGINNING IF NO MATCH AFTER 20 TRIES
|
---|
597 | 30 CONTINUE
|
---|
598 | NREPR1 = NREPR1 + 1
|
---|
599 | IF ( NREPR1 .GT. 20 ) THEN
|
---|
600 | IF ( IDIF .EQ. 1 .AND. NREPRD .LE. 10 ) GOTO 1919
|
---|
601 | GOTO 1
|
---|
602 | ENDIF
|
---|
603 |
|
---|
604 | C FOR TOTAL NUMBER OF PARTICLES ADD 2 FOR LEADER AND ANTILEADER
|
---|
605 | NTOT = NTOTEM + 2
|
---|
606 |
|
---|
607 | C PRODUCTION OF INDIVIDUAL RAPIDITIES FOR ALL SECONDARY PARTICLES
|
---|
608 | CALL PARRAP
|
---|
609 | CC IF ( DEBUG ) THEN
|
---|
610 | CC WRITE (MDEBUG,*) ' RAPIDITIES:'
|
---|
611 | CC WRITE (MDEBUG,134) (I,YR(I), I=3,NTOT)
|
---|
612 | C134 FORMAT(' ',1P, (1X, I4, 5X, G13.6 ))
|
---|
613 | CC ENDIF
|
---|
614 |
|
---|
615 |
|
---|
616 | C CALCULATION OF CENTRAL RAPIDITY WITHOUT (ANTI)LEADER
|
---|
617 | C MULTIPLICITY IN CENTER OF RAPIDITY DISTRIBUTION
|
---|
618 | IZN = 0.D0
|
---|
619 | IF ( IDIF .EQ. 0 ) THEN
|
---|
620 | DO 111 I = 3,NTOT
|
---|
621 | IF ( ABS(YR(I)) .LT. DELRAP ) IZN = IZN + 1
|
---|
622 | 111 CONTINUE
|
---|
623 | IF ( IZN .LT. 1 ) THEN
|
---|
624 | IF ( ISEL .EQ. 0 ) GOTO 30
|
---|
625 | C IN CASE OF FEW PARTICLES, SET PARTICLE NUMBER IN PLATEAU TO 1
|
---|
626 | IZN = 1
|
---|
627 | ENDIF
|
---|
628 | C CENTRAL RAPIDITY DENSITY FOR CHARGED PARTICLES
|
---|
629 | IF ( NTOTEM .GE. 1 ) THEN
|
---|
630 | ZNC = MAX( 1.1D0, DBLE(NCH)*IZN/(DBLE(NTOTEM)*2.D0*DELRAP) )
|
---|
631 | ELSE
|
---|
632 | ZNC = 1.1D0
|
---|
633 | ENDIF
|
---|
634 | ELSE
|
---|
635 | C DIFFRACTION: SHIFT RAPIDITIES + TAKE CENT.RAP.DENS. FROM PARAMETRISAT
|
---|
636 | DO 112 I = 3,NTOT
|
---|
637 | YR(I) = YR(I) + YY0
|
---|
638 | 112 CONTINUE
|
---|
639 | ZNC = MAX( 1.1D0, DCN2 )
|
---|
640 | ENDIF
|
---|
641 |
|
---|
642 | C ZN ACCOUNTS FOR THE RISE OF PT WITH CENTRAL RAP.DENSITY. THE FORMULA
|
---|
643 | C IS A FIT TO UA1 VALUES OF ARNISON ET AL, PHYS.LETT.B118(1982)167
|
---|
644 | C REGARD, THAT OUR ZN IS DEFINED DIFFERENT FROM LITERATURE N BY 1
|
---|
645 | C - - - - - -
|
---|
646 | C MODIFICATION AFTER J.N. CAPDEVIELLE, (DEC.96)
|
---|
647 | IF ( ECMDPM .LE. 500.D0 ) THEN
|
---|
648 | ZN = MAX( 1.00001D0, 3.669D0 / ZNC**0.435D0 + 6.4D0 )
|
---|
649 | ELSE
|
---|
650 | C TAKE INTO ACCOUNT THE RESULTS OF UA1/MIMI EXPERIMENT
|
---|
651 | IF ( ZNC .GE. 3.D0 ) THEN
|
---|
652 | ZN = 1.D0 /(ZNC*0.004111D0 + 0.145D0)+3.D0
|
---|
653 | ELSE
|
---|
654 | C FOR ROCH < 3.00 (MIMI) (TO BE USED IN PTRAM)
|
---|
655 | ZM = 0.0033D0 * (ZNC-1.56D0)**2 + 0.406D0
|
---|
656 | ZN = 2.64D0/ZM + 3.D0
|
---|
657 | ENDIF
|
---|
658 | ENDIF
|
---|
659 | C - - - - - -
|
---|
660 | C NOW SET PARTICLE TYPE AND TRANSV. MOMENTA FOR NEW PARTICLES IN PPARAM
|
---|
661 | C SET ALSO TRANSVERSE MASS FOR ALL PARTICLES (INCL. LEADER+ANTILEADER)
|
---|
662 | CALL PPARAM
|
---|
663 |
|
---|
664 | IF ( IDIF .EQ. 0 ) THEN
|
---|
665 | C NOW SET THE RAPIDITY OF THE ANTILEADER ACCORDING TO THE DISTRIBUTION
|
---|
666 | C IN THE FEYNMAN X VARIABLE; SET THE RAPIDITY OF LEADER TO CONSUME
|
---|
667 | C THE REMAINDER OF ENERGY
|
---|
668 | CALL LEDENY( LEDEFL )
|
---|
669 | IF ( LEDEFL .NE. 0 ) THEN
|
---|
670 | IF ( DEBUG ) WRITE(MDEBUG,*) ' LEDEFL=',LEDEFL
|
---|
671 | GOTO 30
|
---|
672 | ENDIF
|
---|
673 |
|
---|
674 | C CALCULATE FOR SINGLE COLLISION SYSTEM C.M. ENERGY + RAPIDITY SHIFT
|
---|
675 | IF ( GNU .LE. 1.D0 ) THEN
|
---|
676 | JGNU = 0.D0
|
---|
677 | DYGNU = 0.D0
|
---|
678 | ECMJAD = ECMDPM
|
---|
679 | ELSE
|
---|
680 | C MULTIPLE COLLISION IN TARGET
|
---|
681 | JGNU = NINT(GNU-1.D0)
|
---|
682 | C ADD ADDITIONALLY INTERACTING
|
---|
683 | C TARGET NUCLEONS TO GET CORRECT JADACH FILTERING
|
---|
684 | C CHOSE RANDOMLY WHETHER PROTON OR NEUTRON
|
---|
685 | CALL RMMAR( RD,JGNU,1 )
|
---|
686 | IPR = 0
|
---|
687 | INE = 0
|
---|
688 | TARMAS = PAMA(ITYP(2))
|
---|
689 | DO 114 I = 1,JGNU
|
---|
690 | NTOT = NTOT + 1
|
---|
691 | IF ( RD(I) .LE. .5D0 ) THEN
|
---|
692 | ITYP(NTOT) = 13
|
---|
693 | INE = INE + 1
|
---|
694 | ELSE
|
---|
695 | ITYP(NTOT) = 14
|
---|
696 | IPR = IPR + 1
|
---|
697 | ENDIF
|
---|
698 | TMAS(NTOT) = PAMA(ITYP(NTOT))
|
---|
699 | TARMAS = TARMAS + TMAS(NTOT)
|
---|
700 | EA(NTOT) = TMAS(NTOT)
|
---|
701 | PX(NTOT) = 0.D0
|
---|
702 | PY(NTOT) = 0.D0
|
---|
703 | PT2(NTOT) = 0.D0
|
---|
704 | 114 CONTINUE
|
---|
705 |
|
---|
706 | C CALCULATE C.M. ENERGY + RAPIDITY SHIFT
|
---|
707 | * YCMGNU = 0.5D0 * LOG( (ELAB+TARMAS+PLAB)/(ELAB+TARMAS-PLAB) )
|
---|
708 | YCMGNU = 0.5D0 * LOG( (EPLUSP**2 +TARMAS*EPLUSP)/
|
---|
709 | * (PAMA(ITYPE)**2+TARMAS*EPLUSP) )
|
---|
710 | DYGNU = YCM - YCMGNU
|
---|
711 |
|
---|
712 | C CALCULATE RAPIDITIES OF ADDITIONALLY INTERACTING TARGET NUCLEONS
|
---|
713 | C IN THE CM SYSTEM OF NUCLEON-NUCLEON SYSTEM
|
---|
714 | DO 115 I = NTOT-JGNU+1,NTOT
|
---|
715 | YR(I) = - YCM
|
---|
716 | 115 CONTINUE
|
---|
717 | C SHIFT RAPIDITIES INTO CM SYSTEM OF GNU+1 MASSES
|
---|
718 | DO 113 I = 1,NTOT
|
---|
719 | YR(I) = YR(I) + DYGNU
|
---|
720 | 113 CONTINUE
|
---|
721 |
|
---|
722 | C CENTER OF MASS ENERGY OF 1 PROJECTILE AND GNU TARGET NUCLEONS TO
|
---|
723 | C BE USED IN THE JADACH FILTER.
|
---|
724 | ECMJAD = SQRT( PAMA(ITYPE)**2 + TARMAS**2 + 2.D0*TARMAS*ELAB )
|
---|
725 |
|
---|
726 | ENDIF
|
---|
727 |
|
---|
728 | ELSE
|
---|
729 | C IN CASE OF DIFFRACTION SET THE RAPIDITY OF LEADER AND ANTILEADER
|
---|
730 | C IN SUBROUTINE LEADDF
|
---|
731 | DYGNU = 0.D0
|
---|
732 | ECMJAD = ECMDPM
|
---|
733 | CALL LEADDF( IFLGLD )
|
---|
734 | IF ( IFLGLD .NE. 0 ) THEN
|
---|
735 | IF ( DEBUG ) WRITE(MDEBUG,*) ' IFLGLD=',IFLGLD
|
---|
736 | GOTO 30
|
---|
737 | ENDIF
|
---|
738 | ENDIF
|
---|
739 |
|
---|
740 | C CORRECT THE RAPIDITIES TO CONSERVE LONGITUDINAL MOMENTA AND ENERGY
|
---|
741 | C USING THE ALGORITHM OF JADACH (SIMPLIFIED VERSION BY R. ATTALLAH)
|
---|
742 | CALL JADACH( ECMJAD,JADFLG )
|
---|
743 | IF ( JADFLG .NE. 0 ) THEN
|
---|
744 | IF ( DEBUG ) WRITE(MDEBUG,*) ' JADFLG=', JADFLG
|
---|
745 | IF ( JADFLG .GT. 0 ) GOTO 30
|
---|
746 | IF ( JADFLG .LT. 0 ) THEN
|
---|
747 | NREPRD = NREPRD + 1
|
---|
748 | IF ( NREPRD .GT. 10 ) GOTO 1
|
---|
749 | GOTO 1919
|
---|
750 | ENDIF
|
---|
751 | ENDIF
|
---|
752 |
|
---|
753 |
|
---|
754 | C CALCULATE LAB ENERGIES OF SECONDARY PARTICLES FROM THE RAPIDITIES
|
---|
755 | C INCLUDING THE ADDITIONAL TARGET NUCLEONS
|
---|
756 | ETOT = 0.D0
|
---|
757 | DO 135 I = 1,NTOT
|
---|
758 | YR(I) = YR(I) + YCM - DYGNU
|
---|
759 | EA(I) = TMAS(I) * COSH( YR(I) )
|
---|
760 | ETOT = ETOT + EA(I)
|
---|
761 | 135 CONTINUE
|
---|
762 |
|
---|
763 | IF ( DEBUG ) WRITE(MDEBUG,136)
|
---|
764 | * (I,ITYP(I),PX(I),PY(I),YR(I),EA(I),I=1,NTOT)
|
---|
765 | 136 FORMAT(' NO ITYP PX PY YR EA'/
|
---|
766 | * (' ',I4,I3,1X,1P,4G13.6) )
|
---|
767 |
|
---|
768 | C-----------------------------------------------------------------------
|
---|
769 | C LOOP OVER ALL SECONDARY PARTICLES AND THE LEADING PARTICLE
|
---|
770 | C PUT THEM ON THE STACK
|
---|
771 | DO 139 LK = 5,8
|
---|
772 | SECPAR(LK) = CURPAR(LK)
|
---|
773 | 139 CONTINUE
|
---|
774 |
|
---|
775 | C PROCESS LOOP
|
---|
776 | DO 140 J = 1,NTOT
|
---|
777 | C REJECTION OF BACKWARD GOING PARTICLES
|
---|
778 | IF ( YR(J) .LE. 0.D0 ) THEN
|
---|
779 | IF ( DEBUG ) WRITE(MDEBUG,*)'HDPM : YR REJECT PARTICLE ',J
|
---|
780 | GOTO 140
|
---|
781 | ENDIF
|
---|
782 | C CALCULATE THE PROPERTIES OF ALL SECONDARIES
|
---|
783 | C PARTICLE TYPE
|
---|
784 | SECPAR(1) = ITYP(J)
|
---|
785 | C CALCULATE GAMMA FACTOR
|
---|
786 | SECPAR(2) = EA(J) / PAMA(ITYP(J))
|
---|
787 | C TOTAL MOMENTUM SQUARED
|
---|
788 | PTM = EA(J)**2 - PAMA(ITYP(J))**2
|
---|
789 | IF ( PT2(J) .GT. PTM ) THEN
|
---|
790 | IF ( DEBUG ) WRITE(MDEBUG,*)'HDPM : PT REJECT PARTICLE ',J
|
---|
791 | GOTO 140
|
---|
792 | ENDIF
|
---|
793 | C EMISSION ZENITH ANGLE AGAINST TRAJECTORY OF PROJECTILE
|
---|
794 | IF ( PTM .EQ. 0.D0 ) THEN
|
---|
795 | COSTET = 1.D0
|
---|
796 | ELSE
|
---|
797 | COSTET = SQRT( 1.D0 - PT2(J) / PTM )
|
---|
798 | ENDIF
|
---|
799 | C EMISSION AZIMUTH ANGLE
|
---|
800 | IF ( PX(J) .NE. 0.D0 .OR. PY(J) .NE. 0.D0 ) THEN
|
---|
801 | PHIJ = ATAN2( PY(J), PX(J) )
|
---|
802 | ELSE
|
---|
803 | PHIJ = 0.D0
|
---|
804 | ENDIF
|
---|
805 | CALL ADDANG( COSTHE,PHI, COSTET,PHIJ, SECPAR(3),SECPAR(4) )
|
---|
806 | IF ( SECPAR(3) .LT. C(29) ) THEN
|
---|
807 | C OMIT UPWARD GOING PARTICLES
|
---|
808 | IF (DEBUG) WRITE(MDEBUG,*)'HDPM : ANGLE REJECT PARTICLE ',J
|
---|
809 | GOTO 140
|
---|
810 | ENDIF
|
---|
811 | C PUT SECONDARY PARTICLES ON STACK, IF NOT GOING UPWARDS
|
---|
812 | IF ( J .GT. 2 ) THEN
|
---|
813 | CALL TSTACK
|
---|
814 | ELSE
|
---|
815 | C PUT LEADER OR ANTI-LEADER ON STACK, IF NOT GOING UPWARDS
|
---|
816 | IF ( ITYP(J) .GT. 50 ) THEN
|
---|
817 | C LEADER OR ANTI LEADER ARE RESONANCES AND DECAY
|
---|
818 | IRESPAR = IRESPAR + 1
|
---|
819 | IF ( IRESPAR .GE. 1000 ) THEN
|
---|
820 | WRITE(MONIOU,*) 'STACK OF RESDEC RANDOM NUMBERS FULL'
|
---|
821 | IRESPAR = 999
|
---|
822 | ENDIF
|
---|
823 | RESRAN(IRESPAR) = RDRES(J)
|
---|
824 | C COUNTER FOR ENERGY-MULTIPLICITY MATRIX
|
---|
825 | MSMM = MSMM + 1
|
---|
826 | ENDIF
|
---|
827 | CALL TSTACK
|
---|
828 |
|
---|
829 | C CALCULATE ELASTICITY FROM ENERGY OF LEADER (REST OF RESONANCE DECAY)
|
---|
830 | IF ( J. EQ. 1 ) ELASTI = SECPAR(2)*PAMA(NINT(SECPAR(1)))/ELAB
|
---|
831 | ENDIF
|
---|
832 | C COUNTERS FOR FIRST INTERACTION
|
---|
833 | IF ( FIRSTI ) THEN
|
---|
834 | IF ( SECPAR(1) .EQ. 7.D0 .OR. SECPAR(1) .EQ. 8.D0
|
---|
835 | * .OR. SECPAR(1) .EQ. 9.D0 ) THEN
|
---|
836 | IFINPI = IFINPI + 1
|
---|
837 | ELSEIF ( SECPAR(1) .EQ. 13.D0 .OR. SECPAR(1) .EQ. 14.D0
|
---|
838 | * .OR. SECPAR(1) .EQ. 15.D0 .OR. SECPAR(1) .EQ. 25.D0 ) THEN
|
---|
839 | IFINNU = IFINNU + 1
|
---|
840 | ELSEIF ( SECPAR(1) .EQ. 10.D0 .OR. SECPAR(1) .EQ. 11.D0
|
---|
841 | * .OR. SECPAR(1) .EQ. 12.D0 .OR. SECPAR(1) .EQ. 16.D0 ) THEN
|
---|
842 | IFINKA = IFINKA + 1
|
---|
843 | ELSEIF ( SECPAR(1) .GE. 71.D0 .AND. SECPAR(1) .LE. 74.D0) THEN
|
---|
844 | IFINET = IFINET + 1
|
---|
845 | ELSEIF ((SECPAR(1) .GE. 18.D0 .AND. SECPAR(1) .LE. 24.D0)
|
---|
846 | * .OR. (SECPAR(1) .GE. 26.D0 .AND. SECPAR(1) .LE. 32.D0))THEN
|
---|
847 | IFINHY = IFINHY + 1
|
---|
848 | ENDIF
|
---|
849 | ENDIF
|
---|
850 |
|
---|
851 | 140 CONTINUE
|
---|
852 |
|
---|
853 | C COUNTER FOR ENERGY-MULTIPLICITY MATRIX
|
---|
854 | MSMM = MSMM + NTOT - 2
|
---|
855 |
|
---|
856 | C FILL ELASTICITY IN MATRICES
|
---|
857 | MEL = MIN ( 1.D0+10.D0* MAX( 0.D0, ELASTI ) , 11.D0 )
|
---|
858 | MEN = MIN ( 4.D0+ 3.D0*LOG10(MAX( .1D0, EKINL )), 37.D0 )
|
---|
859 | IELDPM(MEN,MEL) = IELDPM(MEN,MEL) + 1
|
---|
860 | IELDPA(MEN,MEL) = IELDPA(MEN,MEL) + 1
|
---|
861 | IF ( ELASTI .LT. 1.D0 ) THEN
|
---|
862 | ELMEAN(MEN) = ELMEAN(MEN) + ELASTI
|
---|
863 | ELMEAA(MEN) = ELMEAA(MEN) + ELASTI
|
---|
864 | ENDIF
|
---|
865 |
|
---|
866 | IF ( FIRSTI ) THEN
|
---|
867 | ELAST = ELASTI
|
---|
868 | FIRSTI = .FALSE.
|
---|
869 | ENDIF
|
---|
870 |
|
---|
871 | IF ( DEBUG ) WRITE(MDEBUG,*) 'HDPM : ELAST=',SNGL(ELASTI),
|
---|
872 | * SNGL(ETOT),SNGL(ELAB)
|
---|
873 |
|
---|
874 | RETURN
|
---|
875 | END
|
---|