3 C*********************************************************************
5 SUBROUTINE LUKFDI_HIJING(KFL1,KFL2,KFL3,KF)
7 C...Purpose: to generate a new flavour pair and combine off a hadron.
8 #include "ludat1_hijing.inc"
9 #include "ludat2_hijing.inc"
11 C...Default flavour values. Input consistency checks.
18 IF(KF1A.LE.10.AND.KF2A.LE.10.AND.KFL1*KFL2.GT.0) RETURN
19 IF(KF1A.GT.10.AND.KF2A.GT.10) RETURN
20 IF((KF1A.GT.10.OR.KF2A.GT.10).AND.KFL1*KFL2.LT.0) RETURN
23 C...Check if tabulated flavour probabilities are to be used.
24 IF(MSTJ(15).EQ.1) THEN
26 IF(KF1A.GE.1.AND.KF1A.LE.6) KTAB1=KF1A
27 KFL1A=MOD(KF1A/1000,10)
28 KFL1B=MOD(KF1A/100,10)
30 IF(KFL1A.GE.1.AND.KFL1A.LE.4.AND.KFL1B.GE.1.AND.KFL1B.LE.4)
31 & KTAB1=6+KFL1A*(KFL1A-2)+2*KFL1B+(KFL1S-1)/2
32 IF(KFL1A.GE.1.AND.KFL1A.LE.4.AND.KFL1A.EQ.KFL1B) KTAB1=KTAB1-1
33 IF(KF1A.GE.1.AND.KF1A.LE.6) KFL1A=KF1A
37 IF(KF2A.GE.1.AND.KF2A.LE.6) KTAB2=KF2A
38 KFL2A=MOD(KF2A/1000,10)
39 KFL2B=MOD(KF2A/100,10)
41 IF(KFL2A.GE.1.AND.KFL2A.LE.4.AND.KFL2B.GE.1.AND.KFL2B.LE.4)
42 & KTAB2=6+KFL2A*(KFL2A-2)+2*KFL2B+(KFL2S-1)/2
43 IF(KFL2A.GE.1.AND.KFL2A.LE.4.AND.KFL2A.EQ.KFL2B) KTAB2=KTAB2-1
45 IF(KTAB1.GE.0.AND.KTAB2.GE.0) GOTO 140
48 C...Parameters and breaking diquark parameter combinations.
52 IF(MSTJ(12).GE.2) THEN
54 PAR4M=1./(3.*SQRT(PARJ(4)))
55 PARDM=PARJ(7)/(PARJ(7)+PAR3M*PARJ(6))
56 PARS0=PARJ(5)*(2.+(1.+PAR2*PAR3M*PARJ(7))*(1.+PAR4M))
57 PARS1=PARJ(7)*PARS0/(2.*PAR3M)+PARJ(5)*(PARJ(6)*(1.+PAR4M)+
58 & PAR2*PAR3M*PARJ(6)*PARJ(7))
59 PARS2=PARJ(5)*2.*PARJ(6)*PARJ(7)*(PAR2*PARJ(7)+(1.+PAR4M)/PAR3M)
60 PARSM=MAX(PARS0,PARS1,PARS2)
61 PAR4=PAR4*(1.+PARSM)/(1.+PARSM/(3.*PAR4M))
64 C...Choice of whether to generate meson or baryon.
68 IF(KF2A.EQ.0.AND.MSTJ(12).GE.1.AND.(1.+PARJ(1))*RLU_HIJING(0)
70 IF(KF2A.GT.10) MBARY=2
71 IF(KF2A.GT.10.AND.KF2A.LE.10000) KFDA=KF2A
74 IF(KF1A.LE.10000) KFDA=KF1A
77 C...Possibility of process diquark -> meson + new diquark.
78 IF(KFDA.NE.0.AND.MSTJ(12).GE.2) THEN
79 KFLDA=MOD(KFDA/1000,10)
80 KFLDB=MOD(KFDA/100,10)
83 IF(MAX(KFLDA,KFLDB).EQ.3) WTDQ=PARS1
84 IF(MIN(KFLDA,KFLDB).EQ.3) WTDQ=PARS2
85 IF(KFLDS.EQ.1) WTDQ=WTDQ/(3.*PAR4M)
86 IF((1.+WTDQ)*RLU_HIJING(0).GT.1.) MBARY=-1
87 IF(MBARY.EQ.-1.AND.KF2A.NE.0) RETURN
90 C...Flavour for meson, possibly with new flavour.
94 IF(KF2A.EQ.0) KFL3=ISIGN(1+INT((2.+PAR2)*RLU_HIJING(0)),-KFL1)
95 KFLA=MAX(KF1A,KF2A+IABS(KFL3))
96 KFLB=MIN(KF1A,KF2A+IABS(KFL3))
97 IF(KFLA.NE.KF1A) KFS=-KFS
99 C...Splitting of diquark into meson plus new diquark.
101 KFL1A=MOD(KF1A/1000,10)
102 KFL1B=MOD(KF1A/100,10)
103 110 KFL1D=KFL1A+INT(RLU_HIJING(0)+0.5)*(KFL1B-KFL1A)
104 KFL1E=KFL1A+KFL1B-KFL1D
105 IF((KFL1D.EQ.3.AND.RLU_HIJING(0).GT.PARDM).OR.(KFL1E.EQ.3.AND.
106 & RLU_HIJING(0).LT.PARDM)) THEN
107 KFL1D=KFL1A+KFL1B-KFL1D
108 KFL1E=KFL1A+KFL1B-KFL1E
110 KFL3A=1+INT((2.+PAR2*PAR3M*PARJ(7))*RLU_HIJING(0))
111 IF((KFL1E.NE.KFL3A.AND.RLU_HIJING(0).GT.(1.+PAR4M)/MAX(2.,1.
112 $ +PAR4M)).OR.(KFL1E.EQ.KFL3A.AND.RLU_HIJING(0).GT.2./MAX(2
113 $ .,1.+PAR4M)))GOTO 110
115 IF(KFL1E.NE.KFL3A) KFLDS=2*INT(RLU_HIJING(0)+1./(1.+PAR4M))+1
116 KFL3=ISIGN(10000+1000*MAX(KFL1E,KFL3A)+100*MIN(KFL1E,KFL3A)+
118 KFLA=MAX(KFL1D,KFL3A)
119 KFLB=MIN(KFL1D,KFL3A)
120 IF(KFLA.NE.KFL1D) KFS=-KFS
123 C...Form meson, with spin and flavour mixing for diagonal states.
124 IF(KFLA.LE.2) KMUL=INT(PARJ(11)+RLU_HIJING(0))
125 IF(KFLA.EQ.3) KMUL=INT(PARJ(12)+RLU_HIJING(0))
126 IF(KFLA.GE.4) KMUL=INT(PARJ(13)+RLU_HIJING(0))
127 IF(KMUL.EQ.0.AND.PARJ(14).GT.0.) THEN
128 IF(RLU_HIJING(0).LT.PARJ(14)) KMUL=2
129 ELSEIF(KMUL.EQ.1.AND.PARJ(15)+PARJ(16)+PARJ(17).GT.0.) THEN
131 IF(RMUL.LT.PARJ(15)) KMUL=3
132 IF(KMUL.EQ.1.AND.RMUL.LT.PARJ(15)+PARJ(16)) KMUL=4
133 IF(KMUL.EQ.1.AND.RMUL.LT.PARJ(15)+PARJ(16)+PARJ(17)) KMUL=5
136 IF(KMUL.EQ.0.OR.KMUL.EQ.3) KFLS=1
138 IF(KFLA.NE.KFLB) THEN
139 KF=(100*KFLA+10*KFLB+KFLS)*KFS*(-1)**KFLA
143 IF(KFLA.LE.3) KF=110*(1+INT(RMIX+PARF(IMIX-1))+
144 & INT(RMIX+PARF(IMIX)))+KFLS
145 IF(KFLA.GE.4) KF=110*KFLA+KFLS
147 IF(KMUL.EQ.2.OR.KMUL.EQ.3) KF=KF+ISIGN(10000,KF)
148 IF(KMUL.EQ.4) KF=KF+ISIGN(20000,KF)
150 C...Generate diquark flavour.
152 120 IF(KF1A.LE.10.AND.KF2A.EQ.0) THEN
154 130 KFLB=1+INT((2.+PAR2*PAR3)*RLU_HIJING(0))
155 KFLC=1+INT((2.+PAR2*PAR3)*RLU_HIJING(0))
157 IF(KFLB.GE.KFLC) KFLDS=3
158 IF(KFLDS.EQ.1.AND.PAR4*RLU_HIJING(0).GT.1.) GOTO 130
159 IF(KFLDS.EQ.3.AND.PAR4.LT.RLU_HIJING(0)) GOTO 130
160 KFL3=ISIGN(1000*MAX(KFLB,KFLC)+100*MIN(KFLB,KFLC)+KFLDS,KFL1)
162 C...Take diquark flavour from input.
163 ELSEIF(KF1A.LE.10) THEN
165 KFLB=MOD(KF2A/1000,10)
166 KFLC=MOD(KF2A/100,10)
169 C...Generate (or take from input) quark to go with diquark.
171 IF(KF2A.EQ.0) KFL3=ISIGN(1+INT((2.+PAR2)*RLU_HIJING(0)),KFL1)
173 KFLB=MOD(KF1A/1000,10)
174 KFLC=MOD(KF1A/100,10)
178 C...SU(6) factors for formation of baryon. Try again if fails.
180 IF(KFLDS.EQ.3.AND.KFLB.NE.KFLC) KBARY=5
181 IF(KFLA.NE.KFLB.AND.KFLA.NE.KFLC) KBARY=KBARY+1
182 WT=PARF(60+KBARY)+PARJ(18)*PARF(70+KBARY)
183 IF(MBARY.EQ.1.AND.MSTJ(12).GE.2) THEN
185 IF(MAX(KFLB,KFLC).EQ.3) WTDQ=PARS1
186 IF(MIN(KFLB,KFLC).EQ.3) WTDQ=PARS2
187 IF(KFLDS.EQ.1) WTDQ=WTDQ/(3.*PAR4M)
188 IF(KFLDS.EQ.1) WT=WT*(1.+WTDQ)/(1.+PARSM/(3.*PAR4M))
189 IF(KFLDS.EQ.3) WT=WT*(1.+WTDQ)/(1.+PARSM)
191 IF(KF2A.EQ.0.AND.WT.LT.RLU_HIJING(0)) GOTO 120
193 C...Form baryon. Distinguish Lambda- and Sigmalike baryons.
194 KFLD=MAX(KFLA,KFLB,KFLC)
195 KFLF=MIN(KFLA,KFLB,KFLC)
196 KFLE=KFLA+KFLB+KFLC-KFLD-KFLF
198 IF((PARF(60+KBARY)+PARJ(18)*PARF(70+KBARY))*RLU_HIJING(0).GT.
199 & PARF(60+KBARY)) KFLS=4
201 IF(KFLS.EQ.2.AND.KFLD.GT.KFLE.AND.KFLE.GT.KFLF) THEN
202 IF(KFLDS.EQ.1.AND.KFLA.EQ.KFLD) KFLL=1
203 IF(KFLDS.EQ.1.AND.KFLA.NE.KFLD) KFLL=INT(0.25+RLU_HIJING(0))
204 IF(KFLDS.EQ.3.AND.KFLA.NE.KFLD) KFLL=INT(0.75+RLU_HIJING(0))
206 IF(KFLL.EQ.0) KF=ISIGN(1000*KFLD+100*KFLE+10*KFLF+KFLS,KFL1)
207 IF(KFLL.EQ.1) KF=ISIGN(1000*KFLD+100*KFLF+10*KFLE+KFLS,KFL1)
211 C...Use tabulated probabilities to select new flavour and hadron.
212 140 IF(KTAB2.EQ.0.AND.MSTJ(12).LE.0) THEN
215 ELSEIF(KTAB2.EQ.0.AND.KTAB1.GE.7.AND.MSTJ(12).LE.1) THEN
218 ELSEIF(KTAB2.EQ.0) THEN
228 RFL=RFL+PARF(120+80*KTAB1+25*KTS+KT3)
230 RFL=RLU_HIJING(0)*RFL
235 RFL=RFL-PARF(120+80*KTAB1+25*KTS+KT3)
236 160 IF(RFL.LE.0.) GOTO 170
239 C...Reconstruct flavour of produced quark/diquark.
243 KFL3=ISIGN(KFL3A,KFL1*(2*KTAB1-13))
246 IF(KTAB3.GE.8) KFL3A=2
247 IF(KTAB3.GE.11) KFL3A=3
248 IF(KTAB3.GE.16) KFL3A=4
249 KFL3B=(KTAB3-6-KFL3A*(KFL3A-2))/2
250 KFL3=1000*KFL3A+100*KFL3B+1
251 IF(KFL3A.EQ.KFL3B.OR.KTAB3.NE.6+KFL3A*(KFL3A-2)+2*KFL3B) KFL3=
253 KFL3=ISIGN(KFL3,KFL1*(13-2*KTAB1))
256 C...Reconstruct meson code.
257 IF(KFL3A.EQ.KFL1A.AND.KFL3B.EQ.KFL1B.AND.(KFL3A.LE.3.OR.
259 RFL=RLU_HIJING(0)*(PARF(143+80*KTAB1+25*KTABS)+PARF(144+80
260 $ *KTAB1+25*KTABS)+PARF(145+80*KTAB1+25*KTABS))
262 IF(RFL.GT.PARF(143+80*KTAB1+25*KTABS)) KF=220+2*KTABS+1
263 IF(RFL.GT.PARF(143+80*KTAB1+25*KTABS)+PARF(144+80*KTAB1+
264 & 25*KTABS)) KF=330+2*KTABS+1
265 ELSEIF(KTAB1.LE.6.AND.KTAB3.LE.6) THEN
266 KFLA=MAX(KTAB1,KTAB3)
267 KFLB=MIN(KTAB1,KTAB3)
269 IF(KFLA.NE.KF1A) KFS=-KFS
270 KF=(100*KFLA+10*KFLB+2*KTABS+1)*KFS*(-1)**KFLA
271 ELSEIF(KTAB1.GE.7.AND.KTAB3.GE.7) THEN
273 IF(KFL1A.EQ.KFL3A) THEN
274 KFLA=MAX(KFL1B,KFL3B)
275 KFLB=MIN(KFL1B,KFL3B)
276 IF(KFLA.NE.KFL1B) KFS=-KFS
277 ELSEIF(KFL1A.EQ.KFL3B) THEN
281 ELSEIF(KFL1B.EQ.KFL3A) THEN
284 ELSEIF(KFL1B.EQ.KFL3B) THEN
285 KFLA=MAX(KFL1A,KFL3A)
286 KFLB=MIN(KFL1A,KFL3A)
287 IF(KFLA.NE.KFL1A) KFS=-KFS
290 $ ,'(LUKFDI_HIJING:) no matching flavours for qq -> qq')
293 KF=(100*KFLA+10*KFLB+2*KTABS+1)*KFS*(-1)**KFLA
295 C...Reconstruct baryon code.
306 KFLD=MAX(KFLA,KFLB,KFLC)
307 KFLF=MIN(KFLA,KFLB,KFLC)
308 KFLE=KFLA+KFLB+KFLC-KFLD-KFLF
309 IF(KTABS.EQ.0) KF=ISIGN(1000*KFLD+100*KFLF+10*KFLE+2,KFL1)
310 IF(KTABS.GE.1) KF=ISIGN(1000*KFLD+100*KFLE+10*KFLF+2*KTABS,KFL1)
313 C...Check that constructed flavour code is an allowed one.
318 $ ,'(LUKFDI_HIJING:) user-defined flavour probabilities'