3 C*********************************************************************
5 SUBROUTINE PYREMN_HIJING(IPU1,IPU2)
7 C...Adds on target remnants (one or two from each side) and
8 C...includes primordial kT.
10 #include "histrng.inc"
11 C...COMMON BLOCK FROM HIJING
12 #include "lujets_hijing.inc"
13 #include "ludat1_hijing.inc"
14 #include "ludat2_hijing.inc"
15 #include "pypars_hijing.inc"
16 #include "pyint1_hijing.inc"
17 DIMENSION KFLCH(2),KFLSP(2),CHI(2),PMS(6),IS(2),ROBO(5)
19 C...Special case for lepton-lepton interaction.
20 IF(MINT(43).EQ.1) THEN
30 C...Find event type, set pointers.
31 IF(IPU1.EQ.0.AND.IPU2.EQ.0) RETURN
36 IF(ISUB.EQ.95) ILEP=-1
37 IF(ILEP.EQ.1) IQ=MINT(84)+1
38 IF(ILEP.EQ.2) IQ=MINT(84)+2
43 C...Define initial partons, including primordial kT.
53 ELSEIF(MINT(40+JT).EQ.1.AND.IPU.NE.0) THEN
62 C...No primordial kT or chosen according to truncated Gaussian or
69 SS_W2=(PP(IHNT2(11),4)+PT(IHNT2(12),4))**2
70 & -(PP(IHNT2(11),1)+PT(IHNT2(12),1))**2
71 & -(PP(IHNT2(11),2)+PT(IHNT2(12),2))**2
72 & -(PP(IHNT2(11),3)+PT(IHNT2(12),3))**2
74 C********this is s of the current NN collision
75 IF(SS_W2.LE.4.0*PARP(93)**2) GOTO 1211
77 IF(IHPR2(5).LE.0) THEN
78 120 IF(MSTP(91).LE.0) THEN
80 ELSEIF(MSTP(91).EQ.1) THEN
81 PTL=PARP(91)*SQRT(-LOG(RLU_HIJING(0)))
85 PTL=-PARP(92)*LOG(RPT1*RPT2)
87 IF(PTL.GT.PARP(93)) GOTO 120
88 PHI=PARU(2)*RLU_HIJING(0)
91 ELSE IF(IHPR2(5).EQ.1) THEN
92 IF(JT.EQ.1) JPT=NFP(IHNT2(11),11)
93 IF(JT.EQ.2) JPT=NFT(IHNT2(12),11)
94 1205 PTGS=PARP(91)*SQRT(-LOG(RLU_HIJING(0)))
95 IF(PTGS.GT.PARP(93)) GO TO 1205
96 PHI=2.0*HIPR1(40)*RLU_HIJING(0)
100 PKCSQ=PARP(91)*SQRT(-LOG(RLU_HIJING(0)))
101 PHI=2.0*HIPR1(40)*RLU_HIJING(0)
102 RPT1=RPT1+PKCSQ*COS(PHI)
103 RPT2=RPT2+PKCSQ*SIN(PHI)
105 IF(RPT1**2+RPT2**2.GE.SS_W2/4.0) GO TO 1205
108 C ********When initial interaction among soft partons is
109 C assumed the primordial pt comes from the sum of
110 C pt of JPT-1 number of initial interaction, JPT
111 C is the number of interaction including present
112 C one that nucleon hassuffered
115 PMS(JT)=P(I,5)**2+P(I,1)**2+P(I,2)**2
121 SHS=(1.-VINT(43-JT))*Q2/VINT(43-JT)+VINT(5-JT)**2
125 C...Kinematics construction for initial partons.
128 IF(ILEP.EQ.0) SHS=VINT(141)*VINT(142)*VINT(2)+
129 &(P(I1,1)+P(I2,1))**2+(P(I1,2)+P(I2,2))**2
130 SHR=SQRT(MAX(0.,SHS))
132 IF((SHS-PMS(1)-PMS(2))**2-4.*PMS(1)*PMS(2).LE.0.) GOTO 110
133 P(I1,4)=0.5*(SHR+(PMS(1)-PMS(2))/SHR)
134 P(I1,3)=SQRT(MAX(0.,P(I1,4)**2-PMS(1)))
137 ELSEIF(ILEP.EQ.1) THEN
142 ELSEIF(ILEP.EQ.2) THEN
148 IF(MINT(43).EQ.1) RETURN
150 C...Transform partons to overall CM-frame (not for leptoproduction).
152 ROBO(3)=(P(I1,1)+P(I2,1))/SHR
153 ROBO(4)=(P(I1,2)+P(I2,2))/SHR
154 CALL LUDBRB_HIJING(I1,I2,0.,0.,-DBLE(ROBO(3)),-DBLE(ROBO(4)),0D0
156 ROBO(2)=ULANGL_HIJING(P(I1,1),P(I1,2))
157 CALL LUDBRB_HIJING(I1,I2,0.,-ROBO(2),0D0,0D0,0D0)
158 ROBO(1)=ULANGL_HIJING(P(I1,3),P(I1,1))
159 CALL LUDBRB_HIJING(I1,I2,-ROBO(1),0.,0D0,0D0,0D0)
160 NMAX=MAX(MINT(52),IPU1,IPU2)
161 CALL LUDBRB_HIJING(I1,NMAX,ROBO(1),ROBO(2),DBLE(ROBO(3))
162 $ ,DBLE(ROBO(4)),0D0)
163 ROBO(5)=MAX(-0.999999,MIN(0.999999,(VINT(141)-VINT(142))/
164 & (VINT(141)+VINT(142))))
165 CALL LUDBRB_HIJING(I1,NMAX,0.,0.,0D0,0D0,DBLE(ROBO(5)))
168 C...Check invariant mass of remnant system:
169 C...hadronic events or leptoproduction.
171 IF(MSTP(81).LE.0.OR.MSTP(82).LE.0.OR.ISUB.EQ.95) THEN
175 PEH=P(I1,4)+P(I2,4)+0.5*VINT(1)*(VINT(151)+VINT(152))
176 PZH=P(I1,3)+P(I2,3)+0.5*VINT(1)*(VINT(151)-VINT(152))
177 SHH=(VINT(1)-PEH)**2-(P(I1,1)+P(I2,1))**2-(P(I1,2)+P(I2,2))**2-
179 PMMIN=P(MINT(83)+1,5)+P(MINT(83)+2,5)+ULMASS_HIJING(K(I1,2))+
180 & ULMASS_HIJING(K(I2,2))
181 IF(SHR.GE.VINT(1).OR.SHH.LE.(PMMIN+PARP(111))**2) THEN
185 SHR=SQRT(SHH+(P(I1,1)+P(I2,1))**2+(P(I1,2)+P(I2,2))**2)
189 PMS(ILEP)=MAX(0.,PEI**2-PZI**2)
190 PMMIN=P(ILEPR-2,5)+ULMASS_HIJING(K(ILEPR,2))+SQRT(PMS(ILEP))
191 IF(SHR.LE.PMMIN+PARP(111)) THEN
197 C...Subdivide remnant if necessary, store first parton.
200 IF(JT.EQ.ILEP) GOTO 190
203 CALL PYSPLI_HIJING(MINT(10+JT),MINT(12+JT),KFLCH(JT),KFLSP(JT))
213 P(I,5)=ULMASS_HIJING(K(I,2))
215 C...First parton colour connections and transverse mass.
216 KFLS=(3-KCHG(LUCOMP_HIJING(KFLSP(JT)),2)*ISIGN(1,KFLSP(JT)))/2
218 K(IPU,6-KFLS)=MOD(K(IPU,6-KFLS),MSTU(5))+MSTU(5)*I
219 IF(KFLCH(JT).EQ.0) THEN
220 P(I,1)=-P(MINT(83)+JT+2,1)
221 P(I,2)=-P(MINT(83)+JT+2,2)
222 PMS(JT)=P(I,5)**2+P(I,1)**2+P(I,2)**2
224 C...When extra remnant parton or hadron: find relative pT, store.
226 CALL LUPTDI_HIJING(1,P(I,1),P(I,2))
227 PMS(JT+2)=P(I,5)**2+P(I,1)**2+P(I,2)**2
236 P(I,5)=ULMASS_HIJING(K(I,2))
237 P(I,1)=-P(MINT(83)+JT+2,1)-P(I-1,1)
238 P(I,2)=-P(MINT(83)+JT+2,2)-P(I-1,2)
239 PMS(JT+4)=P(I,5)**2+P(I,1)**2+P(I,2)**2
240 C...Relative distribution of energy for particle into two jets.
242 IF(MOD(MINT(10+JT)/1000,10).NE.0) IMB=2
243 IF(IABS(KFLCH(JT)).LE.10.OR.KFLCH(JT).EQ.21) THEN
245 IF(MSTP(92).LE.1) THEN
246 IF(IMB.EQ.1) CHI(JT)=RLU_HIJING(0)
247 IF(IMB.EQ.2) CHI(JT)=1.-SQRT(RLU_HIJING(0))
248 ELSEIF(MSTP(92).EQ.2) THEN
249 CHI(JT)=1.-RLU_HIJING(0)**(1./(1.+CHIK))
250 ELSEIF(MSTP(92).EQ.3) THEN
252 170 CHI(JT)=RLU_HIJING(0)**2
253 IF((CHI(JT)**2/(CHI(JT)**2+CUT**2))**0.25*(1.-CHI(JT))**CHIK
254 & .LT.RLU_HIJING(0)) GOTO 170
257 CUTR=(1.+SQRT(1.+CUT**2))/CUT
258 180 CHIR=CUT*CUTR**RLU_HIJING(0)
259 CHI(JT)=(CHIR**2-CUT**2)/(2.*CHIR)
260 IF((1.-CHI(JT))**CHIK.LT.RLU_HIJING(0)) GOTO 180
262 C...Relative distribution of energy for particle into jet plus particle.
265 IF(MSTP(92).LE.1) THEN
266 IF(IMB.EQ.1) CHI(JT)=RLU_HIJING(0)
267 IF(IMB.EQ.2) CHI(JT)=1.-SQRT(RLU_HIJING(0))
269 CHI(JT)=1.-RLU_HIJING(0)**(1./(1.+PARP(93+2*IMB)))
271 IF(MOD(KFLCH(JT)/1000,10).NE.0) CHI(JT)=1.-CHI(JT)
272 IF (CHI(JT) .EQ. 1. .OR. CHI(JT) .EQ. 0.) GOTO 111
274 PMS(JT)=PMS(JT+4)/CHI(JT)+PMS(JT+2)/(1.-CHI(JT))
275 KFLS=KCHG(LUCOMP_HIJING(KFLCH(JT)),2)*ISIGN(1,KFLCH(JT))
280 K(IPU,6-KFLS)=MOD(K(IPU,6-KFLS),MSTU(5))+MSTU(5)*I
284 IF(SHR.LE.SQRT(PMS(1))+SQRT(PMS(2))) GOTO 140
287 C...Reconstruct kinematics of remnants.
289 IF(JT.EQ.ILEP) GOTO 200
290 PE=0.5*(SHR+(PMS(JT)-PMS(3-JT))/SHR)
291 PZ=SQRT(PE**2-PMS(JT))
292 IF(KFLCH(JT).EQ.0) THEN
294 P(IS(JT),3)=PZ*(-1)**(JT-1)
297 P(IS(JT)+1,4)=0.5*(PW1+PMS(JT+4)/PW1)
298 P(IS(JT)+1,3)=0.5*(PW1-PMS(JT+4)/PW1)*(-1)**(JT-1)
299 P(IS(JT),4)=PE-P(IS(JT)+1,4)
300 P(IS(JT),3)=PZ*(-1)**(JT-1)-P(IS(JT)+1,3)
304 C...Hadronic events: boost remnants to correct longitudinal frame.
306 CALL LUDBRB_HIJING(NS+1,N,0.,0.,0D0,0D0,-DBLE(PZH/(VINT(1)-PEH)
308 C...Leptoproduction events: boost colliding subsystem.
310 NMAX=MAX(IP,MINT(52))
312 PZF=PZ*(-1)**(ILEP-1)
313 PT2=P(ILEPR,1)**2+P(ILEPR,2)**2
314 PHIPT=ULANGL_HIJING(P(ILEPR,1),P(ILEPR,2))
315 CALL LUDBRB_HIJING(MINT(84)+1,NMAX,0.,-PHIPT,0D0,0D0,0D0)
316 RQP=P(IQ,3)*(PT2+PEI**2)-P(IQ,4)*PEI*PZI
317 SINTH=P(IQ,4)*SQRT(PT2*(PT2+PEI**2)/(RQP**2+PT2*
318 & P(IQ,4)**2*PZI**2))*SIGN(1.,-RQP)
319 CALL LUDBRB_HIJING(MINT(84)+1,NMAX,ASIN(SINTH),0.,0D0,0D0,0D0)
320 BETAX=(-PEI*PZI*SINTH+SQRT(PT2*(PT2+PEI**2-(PZI*SINTH)**2)))/
322 CALL LUDBRB_HIJING(MINT(84)+1,NMAX,0.,0.,DBLE(BETAX),0D0,0D0)
323 CALL LUDBRB_HIJING(MINT(84)+1,NMAX,0.,PHIPT,0D0,0D0,0D0)
326 BETAZ=(-PEM*PZM+PZF*SQRT(PZF**2+PEM**2-PZM**2))/(PZF**2+PEM**2)
327 CALL LUDBRB_HIJING(MINT(84)+1,NMAX,0.,0.,0D0,0D0,DBLE(BETAZ))
328 CALL LUDBRB_HIJING(I1,I2,ASIN(SINTH),0.,DBLE(BETAX),0D0,0D0)
329 CALL LUDBRB_HIJING(I1,I2,0.,PHIPT,0D0,0D0,DBLE(BETAZ))