]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HIJING/hipyset1_35/luprep_hijing.F
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fixing the decoding of regional header bits.
[u/mrichter/AliRoot.git] / HIJING / hipyset1_35 / luprep_hijing.F
CommitLineData
e74335a4 1* $Id$
2
3C*********************************************************************
4
5 SUBROUTINE LUPREP_HIJING(IP)
6
7C...Purpose: to rearrange partons along strings, to allow small systems
8C...to collapse into one or two particles and to check flavours.
9 IMPLICIT DOUBLE PRECISION(D)
10#include "lujets_hijing.inc"
11#include "ludat1_hijing.inc"
12#include "ludat2_hijing.inc"
13#include "ludat3_hijing.inc"
14 DIMENSION DPS(5),DPC(5),UE(3)
15
16C...Rearrange parton shower product listing along strings: begin loop.
17 I1=N
18 DO 130 MQGST=1,2
19 DO 120 I=MAX(1,IP),N
20 IF(K(I,1).NE.3) GOTO 120
21 KC=LUCOMP_HIJING(K(I,2))
22 IF(KC.EQ.0) GOTO 120
23 KQ=KCHG(KC,2)
24 IF(KQ.EQ.0.OR.(MQGST.EQ.1.AND.KQ.EQ.2)) GOTO 120
25
26C...Pick up loose string end.
27 KCS=4
28 IF(KQ*ISIGN(1,K(I,2)).LT.0) KCS=5
29 IA=I
30 NSTP=0
31 100 NSTP=NSTP+1
32 IF(NSTP.GT.4*N) THEN
33 CALL LUERRM_HIJING(14
34 $ ,'(LUPREP_HIJING:) caught in infinite loop')
35 RETURN
36 ENDIF
37
38C...Copy undecayed parton.
39 IF(K(IA,1).EQ.3) THEN
40 IF(I1.GE.MSTU(4)-MSTU(32)-5) THEN
41 CALL LUERRM_HIJING(11
42 $ ,'(LUPREP_HIJING:) no more memory left in LUJETS_HIJING'
43 $ )
44 RETURN
45 ENDIF
46 I1=I1+1
47 K(I1,1)=2
48 IF(NSTP.GE.2.AND.IABS(K(IA,2)).NE.21) K(I1,1)=1
49 K(I1,2)=K(IA,2)
50 K(I1,3)=IA
51 K(I1,4)=0
52 K(I1,5)=0
53 DO 110 J=1,5
54 P(I1,J)=P(IA,J)
55 110 V(I1,J)=V(IA,J)
56 K(IA,1)=K(IA,1)+10
57 IF(K(I1,1).EQ.1) GOTO 120
58 ENDIF
59
60C...Go to next parton in colour space.
61 IB=IA
62 IF(MOD(K(IB,KCS)/MSTU(5)**2,2).EQ.0.AND.MOD(K(IB,KCS),MSTU(5)).
63 &NE.0) THEN
64 IA=MOD(K(IB,KCS),MSTU(5))
65 K(IB,KCS)=K(IB,KCS)+MSTU(5)**2
66 MREV=0
67 ELSE
68 IF(K(IB,KCS).GE.2*MSTU(5)**2.OR.MOD(K(IB,KCS)/MSTU(5),MSTU(5)).
69 & EQ.0) KCS=9-KCS
70 IA=MOD(K(IB,KCS)/MSTU(5),MSTU(5))
71 K(IB,KCS)=K(IB,KCS)+2*MSTU(5)**2
72 MREV=1
73 ENDIF
74 IF(IA.LE.0.OR.IA.GT.N) THEN
75 CALL LUERRM_HIJING(12
76 $ ,'(LUPREP_HIJING:) colour rearrangement failed')
77 RETURN
78 ENDIF
79 IF(MOD(K(IA,4)/MSTU(5),MSTU(5)).EQ.IB.OR.MOD(K(IA,5)/MSTU(5),
80 &MSTU(5)).EQ.IB) THEN
81 IF(MREV.EQ.1) KCS=9-KCS
82 IF(MOD(K(IA,KCS)/MSTU(5),MSTU(5)).NE.IB) KCS=9-KCS
83 K(IA,KCS)=K(IA,KCS)+2*MSTU(5)**2
84 ELSE
85 IF(MREV.EQ.0) KCS=9-KCS
86 IF(MOD(K(IA,KCS),MSTU(5)).NE.IB) KCS=9-KCS
87 K(IA,KCS)=K(IA,KCS)+MSTU(5)**2
88 ENDIF
89 IF(IA.NE.I) GOTO 100
90 K(I1,1)=1
91 120 CONTINUE
92 130 CONTINUE
93 N=I1
94
95C...Find lowest-mass colour singlet jet system, OK if above threshold.
96 IF(MSTJ(14).LE.0) GOTO 320
97 NS=N
98 140 NSIN=N-NS
99 PDM=1.+PARJ(32)
100 IC=0
101 DO 190 I=MAX(1,IP),NS
102 IF(K(I,1).NE.1.AND.K(I,1).NE.2) THEN
103 ELSEIF(K(I,1).EQ.2.AND.IC.EQ.0) THEN
104 NSIN=NSIN+1
105 IC=I
106 DO 150 J=1,4
107 150 DPS(J)=P(I,J)
108 MSTJ(93)=1
109 DPS(5)=ULMASS_HIJING(K(I,2))
110 ELSEIF(K(I,1).EQ.2) THEN
111 DO 160 J=1,4
112 160 DPS(J)=DPS(J)+P(I,J)
113 ELSEIF(IC.NE.0.AND.KCHG(LUCOMP_HIJING(K(I,2)),2).NE.0) THEN
114 DO 170 J=1,4
115 170 DPS(J)=DPS(J)+P(I,J)
116 MSTJ(93)=1
117 DPS(5)=DPS(5)+ULMASS_HIJING(K(I,2))
118 PD=SQRT(MAX(0D0,DPS(4)**2-DPS(1)**2-DPS(2)**2-DPS(3)**2))-DPS(5)
119 IF(PD.LT.PDM) THEN
120 PDM=PD
121 DO 180 J=1,5
122 180 DPC(J)=DPS(J)
123 IC1=IC
124 IC2=I
125 ENDIF
126 IC=0
127 ELSE
128 NSIN=NSIN+1
129 ENDIF
130 190 CONTINUE
131 IF(PDM.GE.PARJ(32)) GOTO 320
132
133C...Fill small-mass system as cluster.
134 NSAV=N
135 PECM=SQRT(MAX(0D0,DPC(4)**2-DPC(1)**2-DPC(2)**2-DPC(3)**2))
136 K(N+1,1)=11
137 K(N+1,2)=91
138 K(N+1,3)=IC1
139 K(N+1,4)=N+2
140 K(N+1,5)=N+3
141 P(N+1,1)=DPC(1)
142 P(N+1,2)=DPC(2)
143 P(N+1,3)=DPC(3)
144 P(N+1,4)=DPC(4)
145 P(N+1,5)=PECM
146
147C...Form two particles from flavours of lowest-mass system, if feasible.
148 K(N+2,1)=1
149 K(N+3,1)=1
150 IF(MSTU(16).NE.2) THEN
151 K(N+2,3)=N+1
152 K(N+3,3)=N+1
153 ELSE
154 K(N+2,3)=IC1
155 K(N+3,3)=IC2
156 ENDIF
157 K(N+2,4)=0
158 K(N+3,4)=0
159 K(N+2,5)=0
160 K(N+3,5)=0
161 IF(IABS(K(IC1,2)).NE.21) THEN
162 KC1=LUCOMP_HIJING(K(IC1,2))
163 KC2=LUCOMP_HIJING(K(IC2,2))
164 IF(KC1.EQ.0.OR.KC2.EQ.0) GOTO 320
165 KQ1=KCHG(KC1,2)*ISIGN(1,K(IC1,2))
166 KQ2=KCHG(KC2,2)*ISIGN(1,K(IC2,2))
167 IF(KQ1+KQ2.NE.0) GOTO 320
168 200 CALL LUKFDI_HIJING(K(IC1,2),0,KFLN,K(N+2,2))
169 CALL LUKFDI_HIJING(K(IC2,2),-KFLN,KFLDMP,K(N+3,2))
170 IF(K(N+2,2).EQ.0.OR.K(N+3,2).EQ.0) GOTO 200
171 ELSE
172 IF(IABS(K(IC2,2)).NE.21) GOTO 320
173 210 CALL LUKFDI_HIJING(1+INT((2.+PARJ(2))*RLU_HIJING(0)),0,KFLN
174 $ ,KFDMP)
175 CALL LUKFDI_HIJING(KFLN,0,KFLM,K(N+2,2))
176 CALL LUKFDI_HIJING(-KFLN,-KFLM,KFLDMP,K(N+3,2))
177 IF(K(N+2,2).EQ.0.OR.K(N+3,2).EQ.0) GOTO 210
178 ENDIF
179 P(N+2,5)=ULMASS_HIJING(K(N+2,2))
180 P(N+3,5)=ULMASS_HIJING(K(N+3,2))
181 IF(P(N+2,5)+P(N+3,5)+PARJ(64).GE.PECM.AND.NSIN.EQ.1) GOTO 320
182 IF(P(N+2,5)+P(N+3,5)+PARJ(64).GE.PECM) GOTO 260
183
184C...Perform two-particle decay of jet system, if possible.
185 IF(PECM.GE.0.02*DPC(4)) THEN
186 PA=SQRT((PECM**2-(P(N+2,5)+P(N+3,5))**2)*(PECM**2-
187 & (P(N+2,5)-P(N+3,5))**2))/(2.*PECM)
188 UE(3)=2.*RLU_HIJING(0)-1.
189 PHI=PARU(2)*RLU_HIJING(0)
190 UE(1)=SQRT(1.-UE(3)**2)*COS(PHI)
191 UE(2)=SQRT(1.-UE(3)**2)*SIN(PHI)
192 DO 220 J=1,3
193 P(N+2,J)=PA*UE(J)
194 220 P(N+3,J)=-PA*UE(J)
195 P(N+2,4)=SQRT(PA**2+P(N+2,5)**2)
196 P(N+3,4)=SQRT(PA**2+P(N+3,5)**2)
197 CALL LUDBRB_HIJING(N+2,N+3,0.,0.,DPC(1)/DPC(4),DPC(2)/DPC(4),
198 & DPC(3)/DPC(4))
199 ELSE
200 NP=0
201 DO 230 I=IC1,IC2
202 230 IF(K(I,1).EQ.1.OR.K(I,1).EQ.2) NP=NP+1
203 HA=P(IC1,4)*P(IC2,4)-P(IC1,1)*P(IC2,1)-P(IC1,2)*P(IC2,2)-
204 & P(IC1,3)*P(IC2,3)
205 IF(NP.GE.3.OR.HA.LE.1.25*P(IC1,5)*P(IC2,5)) GOTO 260
206 HD1=0.5*(P(N+2,5)**2-P(IC1,5)**2)
207 HD2=0.5*(P(N+3,5)**2-P(IC2,5)**2)
208 HR=SQRT(MAX(0.,((HA-HD1-HD2)**2-(P(N+2,5)*P(N+3,5))**2)/
209 & (HA**2-(P(IC1,5)*P(IC2,5))**2)))-1.
210 HC=P(IC1,5)**2+2.*HA+P(IC2,5)**2
211 HK1=((P(IC2,5)**2+HA)*HR+HD1-HD2)/HC
212 HK2=((P(IC1,5)**2+HA)*HR+HD2-HD1)/HC
213 DO 240 J=1,4
214 P(N+2,J)=(1.+HK1)*P(IC1,J)-HK2*P(IC2,J)
215 240 P(N+3,J)=(1.+HK2)*P(IC2,J)-HK1*P(IC1,J)
216 ENDIF
217 DO 250 J=1,4
218 V(N+1,J)=V(IC1,J)
219 V(N+2,J)=V(IC1,J)
220 250 V(N+3,J)=V(IC2,J)
221 V(N+1,5)=0.
222 V(N+2,5)=0.
223 V(N+3,5)=0.
224 N=N+3
225 GOTO 300
226
227C...Else form one particle from the flavours available, if possible.
228 260 K(N+1,5)=N+2
229 IF(IABS(K(IC1,2)).GT.100.AND.IABS(K(IC2,2)).GT.100) THEN
230 GOTO 320
231 ELSEIF(IABS(K(IC1,2)).NE.21) THEN
232 CALL LUKFDI_HIJING(K(IC1,2),K(IC2,2),KFLDMP,K(N+2,2))
233 ELSE
234 KFLN=1+INT((2.+PARJ(2))*RLU_HIJING(0))
235 CALL LUKFDI_HIJING(KFLN,-KFLN,KFLDMP,K(N+2,2))
236 ENDIF
237 IF(K(N+2,2).EQ.0) GOTO 260
238 P(N+2,5)=ULMASS_HIJING(K(N+2,2))
239
240C...Find parton/particle which combines to largest extra mass.
241 IR=0
242 HA=0.
243 DO 280 MCOMB=1,3
244 IF(IR.NE.0) GOTO 280
245 DO 270 I=MAX(1,IP),N
246 IF(K(I,1).LE.0.OR.K(I,1).GT.10.OR.(I.GE.IC1.AND.I.LE.IC2.
247 &AND.K(I,1).GE.1.AND.K(I,1).LE.2)) GOTO 270
248 IF(MCOMB.EQ.1) KCI=LUCOMP_HIJING(K(I,2))
249 IF(MCOMB.EQ.1.AND.KCI.EQ.0) GOTO 270
250 IF(MCOMB.EQ.1.AND.KCHG(KCI,2).EQ.0.AND.I.LE.NS) GOTO 270
251 IF(MCOMB.EQ.2.AND.IABS(K(I,2)).GT.10.AND.IABS(K(I,2)).LE.100)
252 &GOTO 270
253 HCR=DPC(4)*P(I,4)-DPC(1)*P(I,1)-DPC(2)*P(I,2)-DPC(3)*P(I,3)
254 IF(HCR.GT.HA) THEN
255 IR=I
256 HA=HCR
257 ENDIF
258 270 CONTINUE
259 280 CONTINUE
260
261C...Shuffle energy and momentum to put new particle on mass shell.
262 HB=PECM**2+HA
263 HC=P(N+2,5)**2+HA
264 HD=P(IR,5)**2+HA
265C******************CHANGES BY HIJING************
266 HK2=0.0
267 IF(HA**2-(PECM*P(IR,5))**2.EQ.0.0.OR.HB+HD.EQ.0.0) GO TO 285
268C******************
269 HK2=0.5*(HB*SQRT(((HB+HC)**2-4.*(HB+HD)*P(N+2,5)**2)/
270 &(HA**2-(PECM*P(IR,5))**2))-(HB+HC))/(HB+HD)
271 285 HK1=(0.5*(P(N+2,5)**2-PECM**2)+HD*HK2)/HB
272 DO 290 J=1,4
273 P(N+2,J)=(1.+HK1)*DPC(J)-HK2*P(IR,J)
274 P(IR,J)=(1.+HK2)*P(IR,J)-HK1*DPC(J)
275 V(N+1,J)=V(IC1,J)
276 290 V(N+2,J)=V(IC1,J)
277 V(N+1,5)=0.
278 V(N+2,5)=0.
279 N=N+2
280
281C...Mark collapsed system and store daughter pointers. Iterate.
282 300 DO 310 I=IC1,IC2
283 IF((K(I,1).EQ.1.OR.K(I,1).EQ.2).AND.KCHG(LUCOMP_HIJING(K(I,2))
284 $ ,2).NE.0)THEN
285 K(I,1)=K(I,1)+10
286 IF(MSTU(16).NE.2) THEN
287 K(I,4)=NSAV+1
288 K(I,5)=NSAV+1
289 ELSE
290 K(I,4)=NSAV+2
291 K(I,5)=N
292 ENDIF
293 ENDIF
294 310 CONTINUE
295 IF(N.LT.MSTU(4)-MSTU(32)-5) GOTO 140
296
297C...Check flavours and invariant masses in parton systems.
298 320 NP=0
299 KFN=0
300 KQS=0
301 DO 330 J=1,5
302 330 DPS(J)=0.
303 DO 360 I=MAX(1,IP),N
304 IF(K(I,1).LE.0.OR.K(I,1).GT.10) GOTO 360
305 KC=LUCOMP_HIJING(K(I,2))
306 IF(KC.EQ.0) GOTO 360
307 KQ=KCHG(KC,2)*ISIGN(1,K(I,2))
308 IF(KQ.EQ.0) GOTO 360
309 NP=NP+1
310 IF(KQ.NE.2) THEN
311 KFN=KFN+1
312 KQS=KQS+KQ
313 MSTJ(93)=1
314 DPS(5)=DPS(5)+ULMASS_HIJING(K(I,2))
315 ENDIF
316 DO 340 J=1,4
317 340 DPS(J)=DPS(J)+P(I,J)
318 IF(K(I,1).EQ.1) THEN
319 IF(NP.NE.1.AND.(KFN.EQ.1.OR.KFN.GE.3.OR.KQS.NE.0)) CALL
320 & LUERRM_HIJING(2
321 $ ,'(LUPREP_HIJING:) unphysical flavour combination')
322 IF(NP.NE.1.AND.DPS(4)**2-DPS(1)**2-DPS(2)**2-DPS(3)**2.LT.
323 & (0.9*PARJ(32)+DPS(5))**2) CALL LUERRM_HIJING(3,
324 & '(LUPREP_HIJING:) too small mass in jet system')
325 NP=0
326 KFN=0
327 KQS=0
328 DO 350 J=1,5
329 350 DPS(J)=0.
330 ENDIF
331 360 CONTINUE
332
333 RETURN
334 END