]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | *CMZ : 17/07/98 15.44.31 by Federico Carminati |
2 | *-- Author : | |
3 | C********************************************************************* | |
4 | ||
5 | SUBROUTINE LU4ENT(IP,KF1,KF2,KF3,KF4,PECM,X1,X2,X4,X12,X14) | |
6 | ||
7 | C...Purpose: to store four partons or particles in their CM frame, with | |
8 | C...the first along the +z axis, the last in the xz plane with x > 0 | |
9 | C...and the second having y < 0 and y > 0 with equal probability. | |
10 | *KEEP,LUJETS. | |
11 | COMMON /LUJETS/ N,K(200000,5),P(200000,5),V(200000,5) | |
12 | SAVE /LUJETS/ | |
13 | *KEEP,LUDAT1. | |
14 | COMMON /LUDAT1/ MSTU(200),PARU(200),MSTJ(200),PARJ(200) | |
15 | SAVE /LUDAT1/ | |
16 | *KEEP,LUDAT2. | |
17 | COMMON /LUDAT2/ KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4) | |
18 | SAVE /LUDAT2/ | |
19 | *KEND. | |
20 | ||
21 | C...Standard checks. | |
22 | MSTU(28)=0 | |
23 | IF(MSTU(12).GE.1) CALL LULIST(0) | |
24 | IPA=MAX(1,IABS(IP)) | |
25 | IF(IPA.GT.MSTU(4)-3) CALL LUERRM(21, | |
26 | &'(LU4ENT:) writing outside LUJETS momory') | |
27 | KC1=LUCOMP(KF1) | |
28 | KC2=LUCOMP(KF2) | |
29 | KC3=LUCOMP(KF3) | |
30 | KC4=LUCOMP(KF4) | |
31 | IF(KC1.EQ.0.OR.KC2.EQ.0.OR.KC3.EQ.0.OR.KC4.EQ.0) CALL LUERRM(12, | |
32 | &'(LU4ENT:) unknown flavour code') | |
33 | ||
34 | C...Find masses. Reset K, P and V vectors. | |
35 | PM1=0. | |
36 | IF(MSTU(10).EQ.1) PM1=P(IPA,5) | |
37 | IF(MSTU(10).GE.2) PM1=ULMASS(KF1) | |
38 | PM2=0. | |
39 | IF(MSTU(10).EQ.1) PM2=P(IPA+1,5) | |
40 | IF(MSTU(10).GE.2) PM2=ULMASS(KF2) | |
41 | PM3=0. | |
42 | IF(MSTU(10).EQ.1) PM3=P(IPA+2,5) | |
43 | IF(MSTU(10).GE.2) PM3=ULMASS(KF3) | |
44 | PM4=0. | |
45 | IF(MSTU(10).EQ.1) PM4=P(IPA+3,5) | |
46 | IF(MSTU(10).GE.2) PM4=ULMASS(KF4) | |
47 | DO 100 I=IPA,IPA+3 | |
48 | DO 100 J=1,5 | |
49 | K(I,J)=0 | |
50 | P(I,J)=0. | |
51 | 100 V(I,J)=0. | |
52 | ||
53 | C...Check flavours. | |
54 | KQ1=KCHG(KC1,2)*ISIGN(1,KF1) | |
55 | KQ2=KCHG(KC2,2)*ISIGN(1,KF2) | |
56 | KQ3=KCHG(KC3,2)*ISIGN(1,KF3) | |
57 | KQ4=KCHG(KC4,2)*ISIGN(1,KF4) | |
58 | IF(KQ1.EQ.0.AND.KQ2.EQ.0.AND.KQ3.EQ.0.AND.KQ4.EQ.0) THEN | |
59 | ELSEIF(KQ1.NE.0.AND.KQ2.EQ.2.AND.KQ3.EQ.2.AND.(KQ1+KQ4.EQ.0.OR. | |
60 | &KQ1+KQ4.EQ.4)) THEN | |
61 | ELSEIF(KQ1.NE.0.AND.KQ1+KQ2.EQ.0.AND.KQ3.NE.0.AND.KQ3+KQ4.EQ.0.) | |
62 | &THEN | |
63 | ELSE | |
64 | CALL LUERRM(2,'(LU4ENT:) unphysical flavour combination') | |
65 | ENDIF | |
66 | K(IPA,2)=KF1 | |
67 | K(IPA+1,2)=KF2 | |
68 | K(IPA+2,2)=KF3 | |
69 | K(IPA+3,2)=KF4 | |
70 | ||
71 | C...Store partons/particles in K vectors for normal case. | |
72 | IF(IP.GE.0) THEN | |
73 | K(IPA,1)=1 | |
74 | IF(KQ1.NE.0.AND.(KQ2.NE.0.OR.KQ3.NE.0.OR.KQ4.NE.0)) K(IPA,1)=2 | |
75 | K(IPA+1,1)=1 | |
76 | IF(KQ2.NE.0.AND.KQ1+KQ2.NE.0.AND.(KQ3.NE.0.OR.KQ4.NE.0)) | |
77 | & K(IPA+1,1)=2 | |
78 | K(IPA+2,1)=1 | |
79 | IF(KQ3.NE.0.AND.KQ4.NE.0) K(IPA+2,1)=2 | |
80 | K(IPA+3,1)=1 | |
81 | ||
82 | C...Store partons for parton shower evolution from q-g-g-qbar or | |
83 | C...g-g-g-g event. | |
84 | ELSEIF(KQ1+KQ2.NE.0) THEN | |
85 | IF(KQ1.EQ.0.OR.KQ2.EQ.0.OR.KQ3.EQ.0.OR.KQ4.EQ.0) CALL LUERRM(2, | |
86 | & '(LU4ENT:) requested flavours can not develop parton shower') | |
87 | K(IPA,1)=3 | |
88 | K(IPA+1,1)=3 | |
89 | K(IPA+2,1)=3 | |
90 | K(IPA+3,1)=3 | |
91 | KCS=4 | |
92 | IF(KQ1.EQ.-1) KCS=5 | |
93 | K(IPA,KCS)=MSTU(5)*(IPA+1) | |
94 | K(IPA,9-KCS)=MSTU(5)*(IPA+3) | |
95 | K(IPA+1,KCS)=MSTU(5)*(IPA+2) | |
96 | K(IPA+1,9-KCS)=MSTU(5)*IPA | |
97 | K(IPA+2,KCS)=MSTU(5)*(IPA+3) | |
98 | K(IPA+2,9-KCS)=MSTU(5)*(IPA+1) | |
99 | K(IPA+3,KCS)=MSTU(5)*IPA | |
100 | K(IPA+3,9-KCS)=MSTU(5)*(IPA+2) | |
101 | ||
102 | C...Store partons for parton shower evolution from q-qbar-q-qbar event. | |
103 | ELSE | |
104 | IF(KQ1.EQ.0.OR.KQ2.EQ.0.OR.KQ3.EQ.0.OR.KQ4.EQ.0) CALL LUERRM(2, | |
105 | & '(LU4ENT:) requested flavours can not develop parton shower') | |
106 | K(IPA,1)=3 | |
107 | K(IPA+1,1)=3 | |
108 | K(IPA+2,1)=3 | |
109 | K(IPA+3,1)=3 | |
110 | K(IPA,4)=MSTU(5)*(IPA+1) | |
111 | K(IPA,5)=K(IPA,4) | |
112 | K(IPA+1,4)=MSTU(5)*IPA | |
113 | K(IPA+1,5)=K(IPA+1,4) | |
114 | K(IPA+2,4)=MSTU(5)*(IPA+3) | |
115 | K(IPA+2,5)=K(IPA+2,4) | |
116 | K(IPA+3,4)=MSTU(5)*(IPA+2) | |
117 | K(IPA+3,5)=K(IPA+3,4) | |
118 | ENDIF | |
119 | ||
120 | C...Check kinematics. | |
121 | MKERR=0 | |
122 | IF(0.5*X1*PECM.LE.PM1.OR.0.5*X2*PECM.LE.PM2.OR.0.5*(2.-X1-X2-X4)* | |
123 | &PECM.LE.PM3.OR.0.5*X4*PECM.LE.PM4) MKERR=1 | |
124 | PA1=SQRT(MAX(0.,(0.5*X1*PECM)**2-PM1**2)) | |
125 | PA2=SQRT(MAX(0.,(0.5*X2*PECM)**2-PM2**2)) | |
126 | PA3=SQRT(MAX(0.,(0.5*(2.-X1-X2-X4)*PECM)**2-PM3**2)) | |
127 | PA4=SQRT(MAX(0.,(0.5*X4*PECM)**2-PM4**2)) | |
128 | X24=X1+X2+X4-1.-X12-X14+(PM3**2-PM1**2-PM2**2-PM4**2)/PECM**2 | |
129 | CTHE4=(X1*X4-2.*X14)*PECM**2/(4.*PA1*PA4) | |
130 | IF(ABS(CTHE4).GE.1.002) MKERR=1 | |
131 | CTHE4=MAX(-1.,MIN(1.,CTHE4)) | |
132 | STHE4=SQRT(1.-CTHE4**2) | |
133 | CTHE2=(X1*X2-2.*X12)*PECM**2/(4.*PA1*PA2) | |
134 | IF(ABS(CTHE2).GE.1.002) MKERR=1 | |
135 | CTHE2=MAX(-1.,MIN(1.,CTHE2)) | |
136 | STHE2=SQRT(1.-CTHE2**2) | |
137 | CPHI2=((X2*X4-2.*X24)*PECM**2-4.*PA2*CTHE2*PA4*CTHE4)/ | |
138 | &(4.*PA2*STHE2*PA4*STHE4) | |
139 | IF(ABS(CPHI2).GE.1.05) MKERR=1 | |
140 | CPHI2=MAX(-1.,MIN(1.,CPHI2)) | |
141 | IF(MKERR.EQ.1) CALL LUERRM(13, | |
142 | &'(LU4ENT:) unphysical kinematical variable setup') | |
143 | ||
144 | C...Store partons/particles in P vectors. | |
145 | P(IPA,3)=PA1 | |
146 | P(IPA,4)=SQRT(PA1**2+PM1**2) | |
147 | P(IPA,5)=PM1 | |
148 | P(IPA+3,1)=PA4*STHE4 | |
149 | P(IPA+3,3)=PA4*CTHE4 | |
150 | P(IPA+3,4)=SQRT(PA4**2+PM4**2) | |
151 | P(IPA+3,5)=PM4 | |
152 | P(IPA+1,1)=PA2*STHE2*CPHI2 | |
153 | P(IPA+1,2)=PA2*STHE2*SQRT(1.-CPHI2**2)*(-1.)**INT(RLU(0)+0.5) | |
154 | P(IPA+1,3)=PA2*CTHE2 | |
155 | P(IPA+1,4)=SQRT(PA2**2+PM2**2) | |
156 | P(IPA+1,5)=PM2 | |
157 | P(IPA+2,1)=-P(IPA+1,1)-P(IPA+3,1) | |
158 | P(IPA+2,2)=-P(IPA+1,2) | |
159 | P(IPA+2,3)=-P(IPA,3)-P(IPA+1,3)-P(IPA+3,3) | |
160 | P(IPA+2,4)=SQRT(P(IPA+2,1)**2+P(IPA+2,2)**2+P(IPA+2,3)**2+PM3**2) | |
161 | P(IPA+2,5)=PM3 | |
162 | ||
163 | C...Set N. Optionally fragment/decay. | |
164 | N=IPA+3 | |
165 | IF(IP.EQ.0) CALL LUEXEC | |
166 | ||
167 | RETURN | |
168 | END |