Merging the VirtualMC branch to the main development branch (HEAD)
[u/mrichter/AliRoot.git] / HIJING / hipyset1_35 / pyklim_hijing.F
CommitLineData
e74335a4 1* $Id$
2
3C***********************************************************************
4
5 SUBROUTINE PYKLIM_HIJING(ILIM)
6
7C...Checks generated variables against pre-set kinematical limits;
8C...also calculates limits on variables used in generation.
9#include "ludat1_hijing.inc"
10#include "ludat2_hijing.inc"
11#include "ludat3_hijing.inc"
12#include "pypars_hijing.inc"
13#include "pysubs_hijing.inc"
14#include "pyint1_hijing.inc"
15#include "pyint2_hijing.inc"
16
17C...Common kinematical expressions.
18 ISUB=MINT(1)
19 IF(ISUB.EQ.96) GOTO 110
20 SQM3=VINT(63)
21 SQM4=VINT(64)
22 IF(ILIM.NE.1) THEN
23 TAU=VINT(21)
24 RM3=SQM3/(TAU*VINT(2))
25 RM4=SQM4/(TAU*VINT(2))
26 BE34=SQRT((1.-RM3-RM4)**2-4.*RM3*RM4)
27 ENDIF
28 PTHMIN=CKIN(3)
29 IF(MIN(SQM3,SQM4).LT.CKIN(6)**2) PTHMIN=MAX(CKIN(3),CKIN(5))
30
31 IF(ILIM.EQ.0) THEN
32C...Check generated values of tau, y*, cos(theta-hat), and tau' against
33C...pre-set kinematical limits.
34 YST=VINT(22)
35 CTH=VINT(23)
36 TAUP=VINT(26)
37 IF(ISET(ISUB).LE.2) THEN
38 X1=SQRT(TAU)*EXP(YST)
39 X2=SQRT(TAU)*EXP(-YST)
40 ELSE
41 X1=SQRT(TAUP)*EXP(YST)
42 X2=SQRT(TAUP)*EXP(-YST)
43 ENDIF
44 XF=X1-X2
45 IF(TAU*VINT(2).LT.CKIN(1)**2) MINT(51)=1
46 IF(CKIN(2).GE.0..AND.TAU*VINT(2).GT.CKIN(2)**2) MINT(51)=1
47 IF(X1.LT.CKIN(21).OR.X1.GT.CKIN(22)) MINT(51)=1
48 IF(X2.LT.CKIN(23).OR.X2.GT.CKIN(24)) MINT(51)=1
49 IF(XF.LT.CKIN(25).OR.XF.GT.CKIN(26)) MINT(51)=1
50 IF(YST.LT.CKIN(7).OR.YST.GT.CKIN(8)) MINT(51)=1
51 IF(ISET(ISUB).EQ.2.OR.ISET(ISUB).EQ.4) THEN
52 PTH=0.5*BE34*SQRT(TAU*VINT(2)*(1.-CTH**2))
53 Y3=YST+0.5*LOG((1.+RM3-RM4+BE34*CTH)/(1.+RM3-RM4-BE34*CTH))
54 Y4=YST+0.5*LOG((1.+RM4-RM3-BE34*CTH)/(1.+RM4-RM3+BE34*CTH))
55 YLARGE=MAX(Y3,Y4)
56 YSMALL=MIN(Y3,Y4)
57 ETALAR=10.
58 ETASMA=-10.
59 STH=SQRT(1.-CTH**2)
60 IF(STH.LT.1.E-6) GOTO 100
61 EXPET3=((1.+RM3-RM4)*SINH(YST)+BE34*COSH(YST)*CTH+
62 & SQRT(((1.+RM3-RM4)*COSH(YST)+BE34*SINH(YST)*CTH)**2-4.*RM3))/
63 & (BE34*STH)
64 EXPET4=((1.-RM3+RM4)*SINH(YST)-BE34*COSH(YST)*CTH+
65 & SQRT(((1.-RM3+RM4)*COSH(YST)-BE34*SINH(YST)*CTH)**2-4.*RM4))/
66 & (BE34*STH)
67 ETA3=LOG(MIN(1.E10,MAX(1.E-10,EXPET3)))
68 ETA4=LOG(MIN(1.E10,MAX(1.E-10,EXPET4)))
69 ETALAR=MAX(ETA3,ETA4)
70 ETASMA=MIN(ETA3,ETA4)
71 100 CTS3=((1.+RM3-RM4)*SINH(YST)+BE34*COSH(YST)*CTH)/
72 & SQRT(((1.+RM3-RM4)*COSH(YST)+BE34*SINH(YST)*CTH)**2-4.*RM3)
73 CTS4=((1.-RM3+RM4)*SINH(YST)-BE34*COSH(YST)*CTH)/
74 & SQRT(((1.-RM3+RM4)*COSH(YST)-BE34*SINH(YST)*CTH)**2-4.*RM4)
75 CTSLAR=MAX(CTS3,CTS4)
76 CTSSMA=MIN(CTS3,CTS4)
77 IF(PTH.LT.PTHMIN) MINT(51)=1
78 IF(CKIN(4).GE.0..AND.PTH.GT.CKIN(4)) MINT(51)=1
79 IF(YLARGE.LT.CKIN(9).OR.YLARGE.GT.CKIN(10)) MINT(51)=1
80 IF(YSMALL.LT.CKIN(11).OR.YSMALL.GT.CKIN(12)) MINT(51)=1
81 IF(ETALAR.LT.CKIN(13).OR.ETALAR.GT.CKIN(14)) MINT(51)=1
82 IF(ETASMA.LT.CKIN(15).OR.ETASMA.GT.CKIN(16)) MINT(51)=1
83 IF(CTSLAR.LT.CKIN(17).OR.CTSLAR.GT.CKIN(18)) MINT(51)=1
84 IF(CTSSMA.LT.CKIN(19).OR.CTSSMA.GT.CKIN(20)) MINT(51)=1
85 IF(CTH.LT.CKIN(27).OR.CTH.GT.CKIN(28)) MINT(51)=1
86 ENDIF
87 IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) THEN
88 IF(TAUP*VINT(2).LT.CKIN(31)**2) MINT(51)=1
89 IF(CKIN(32).GE.0..AND.TAUP*VINT(2).GT.CKIN(32)**2) MINT(51)=1
90 ENDIF
91
92 ELSEIF(ILIM.EQ.1) THEN
93C...Calculate limits on tau
94C...0) due to definition
95 TAUMN0=0.
96 TAUMX0=1.
97C...1) due to limits on subsystem mass
98 TAUMN1=CKIN(1)**2/VINT(2)
99 TAUMX1=1.
100 IF(CKIN(2).GE.0.) TAUMX1=CKIN(2)**2/VINT(2)
101C...2) due to limits on pT-hat (and non-overlapping rapidity intervals)
102 TM3=SQRT(SQM3+PTHMIN**2)
103 TM4=SQRT(SQM4+PTHMIN**2)
104 YDCOSH=1.
105 IF(CKIN(9).GT.CKIN(12)) YDCOSH=COSH(CKIN(9)-CKIN(12))
106 TAUMN2=(TM3**2+2.*TM3*TM4*YDCOSH+TM4**2)/VINT(2)
107 TAUMX2=1.
108C...3) due to limits on pT-hat and cos(theta-hat)
109 CTH2MN=MIN(CKIN(27)**2,CKIN(28)**2)
110 CTH2MX=MAX(CKIN(27)**2,CKIN(28)**2)
111 TAUMN3=0.
112 IF(CKIN(27)*CKIN(28).GT.0.) TAUMN3=
113 & (SQRT(SQM3+PTHMIN**2/(1.-CTH2MN))+
114 & SQRT(SQM4+PTHMIN**2/(1.-CTH2MN)))**2/VINT(2)
115 TAUMX3=1.
116 IF(CKIN(4).GE.0..AND.CTH2MX.LT.1.) TAUMX3=
117 & (SQRT(SQM3+CKIN(4)**2/(1.-CTH2MX))+
118 & SQRT(SQM4+CKIN(4)**2/(1.-CTH2MX)))**2/VINT(2)
119C...4) due to limits on x1 and x2
120 TAUMN4=CKIN(21)*CKIN(23)
121 TAUMX4=CKIN(22)*CKIN(24)
122C...5) due to limits on xF
123 TAUMN5=0.
124 TAUMX5=MAX(1.-CKIN(25),1.+CKIN(26))
125 VINT(11)=MAX(TAUMN0,TAUMN1,TAUMN2,TAUMN3,TAUMN4,TAUMN5)
126 VINT(31)=MIN(TAUMX0,TAUMX1,TAUMX2,TAUMX3,TAUMX4,TAUMX5)
127 IF(MINT(43).EQ.1.AND.(ISET(ISUB).EQ.1.OR.ISET(ISUB).EQ.2)) THEN
128 VINT(11)=0.99999
129 VINT(31)=1.00001
130 ENDIF
131 IF(VINT(31).LE.VINT(11)) MINT(51)=1
132
133 ELSEIF(ILIM.EQ.2) THEN
134C...Calculate limits on y*
135 IF(ISET(ISUB).EQ.3.OR.ISET(ISUB).EQ.4) TAU=VINT(26)
136 TAURT=SQRT(TAU)
137C...0) due to kinematics
138 YSTMN0=LOG(TAURT)
139 YSTMX0=-YSTMN0
140C...1) due to explicit limits
141 YSTMN1=CKIN(7)
142 YSTMX1=CKIN(8)
143C...2) due to limits on x1
144 YSTMN2=LOG(MAX(TAU,CKIN(21))/TAURT)
145 YSTMX2=LOG(MAX(TAU,CKIN(22))/TAURT)
146C...3) due to limits on x2
147 YSTMN3=-LOG(MAX(TAU,CKIN(24))/TAURT)
148 YSTMX3=-LOG(MAX(TAU,CKIN(23))/TAURT)
149C...4) due to limits on xF
150 YEPMN4=0.5*ABS(CKIN(25))/TAURT
151 YSTMN4=SIGN(LOG(SQRT(1.+YEPMN4**2)+YEPMN4),CKIN(25))
152 YEPMX4=0.5*ABS(CKIN(26))/TAURT
153 YSTMX4=SIGN(LOG(SQRT(1.+YEPMX4**2)+YEPMX4),CKIN(26))
154C...5) due to simultaneous limits on y-large and y-small
155 YEPSMN=(RM3-RM4)*SINH(CKIN(9)-CKIN(11))
156 YEPSMX=(RM3-RM4)*SINH(CKIN(10)-CKIN(12))
157 YDIFMN=ABS(LOG(SQRT(1.+YEPSMN**2)-YEPSMN))
158 YDIFMX=ABS(LOG(SQRT(1.+YEPSMX**2)-YEPSMX))
159 YSTMN5=0.5*(CKIN(9)+CKIN(11)-YDIFMN)
160 YSTMX5=0.5*(CKIN(10)+CKIN(12)+YDIFMX)
161C...6) due to simultaneous limits on cos(theta-hat) and y-large or
162C... y-small
163 CTHLIM=SQRT(1.-4.*PTHMIN**2/(BE34*TAU*VINT(2)))
164 RZMN=BE34*MAX(CKIN(27),-CTHLIM)
165 RZMX=BE34*MIN(CKIN(28),CTHLIM)
166 YEX3MX=(1.+RM3-RM4+RZMX)/MAX(1E-10,1.+RM3-RM4-RZMX)
167 YEX4MX=(1.+RM4-RM3-RZMN)/MAX(1E-10,1.+RM4-RM3+RZMN)
168 YEX3MN=MAX(1E-10,1.+RM3-RM4+RZMN)/(1.+RM3-RM4-RZMN)
169 YEX4MN=MAX(1E-10,1.+RM4-RM3-RZMX)/(1.+RM4-RM3+RZMX)
170 YSTMN6=CKIN(9)-0.5*LOG(MAX(YEX3MX,YEX4MX))
171 YSTMX6=CKIN(12)-0.5*LOG(MIN(YEX3MN,YEX4MN))
172 VINT(12)=MAX(YSTMN0,YSTMN1,YSTMN2,YSTMN3,YSTMN4,YSTMN5,YSTMN6)
173 VINT(32)=MIN(YSTMX0,YSTMX1,YSTMX2,YSTMX3,YSTMX4,YSTMX5,YSTMX6)
174 IF(MINT(43).EQ.1) THEN
175 VINT(12)=-0.00001
176 VINT(32)=0.00001
177 ELSEIF(MINT(43).EQ.2) THEN
178 VINT(12)=0.99999*YSTMX0
179 VINT(32)=1.00001*YSTMX0
180 ELSEIF(MINT(43).EQ.3) THEN
181 VINT(12)=-1.00001*YSTMX0
182 VINT(32)=-0.99999*YSTMX0
183 ENDIF
184 IF(VINT(32).LE.VINT(12)) MINT(51)=1
185
186 ELSEIF(ILIM.EQ.3) THEN
187C...Calculate limits on cos(theta-hat)
188 YST=VINT(22)
189C...0) due to definition
190 CTNMN0=-1.
191 CTNMX0=0.
192 CTPMN0=0.
193 CTPMX0=1.
194C...1) due to explicit limits
195 CTNMN1=MIN(0.,CKIN(27))
196 CTNMX1=MIN(0.,CKIN(28))
197 CTPMN1=MAX(0.,CKIN(27))
198 CTPMX1=MAX(0.,CKIN(28))
199C...2) due to limits on pT-hat
200 CTNMN2=-SQRT(1.-4.*PTHMIN**2/(BE34**2*TAU*VINT(2)))
201 CTPMX2=-CTNMN2
202 CTNMX2=0.
203 CTPMN2=0.
204 IF(CKIN(4).GE.0.) THEN
205 CTNMX2=-SQRT(MAX(0.,1.-4.*CKIN(4)**2/(BE34**2*TAU*VINT(2))))
206 CTPMN2=-CTNMX2
207 ENDIF
208C...3) due to limits on y-large and y-small
209 CTNMN3=MIN(0.,MAX((1.+RM3-RM4)/BE34*TANH(CKIN(11)-YST),
210 & -(1.-RM3+RM4)/BE34*TANH(CKIN(10)-YST)))
211 CTNMX3=MIN(0.,(1.+RM3-RM4)/BE34*TANH(CKIN(12)-YST),
212 & -(1.-RM3+RM4)/BE34*TANH(CKIN(9)-YST))
213 CTPMN3=MAX(0.,(1.+RM3-RM4)/BE34*TANH(CKIN(9)-YST),
214 & -(1.-RM3+RM4)/BE34*TANH(CKIN(12)-YST))
215 CTPMX3=MAX(0.,MIN((1.+RM3-RM4)/BE34*TANH(CKIN(10)-YST),
216 & -(1.-RM3+RM4)/BE34*TANH(CKIN(11)-YST)))
217 VINT(13)=MAX(CTNMN0,CTNMN1,CTNMN2,CTNMN3)
218 VINT(33)=MIN(CTNMX0,CTNMX1,CTNMX2,CTNMX3)
219 VINT(14)=MAX(CTPMN0,CTPMN1,CTPMN2,CTPMN3)
220 VINT(34)=MIN(CTPMX0,CTPMX1,CTPMX2,CTPMX3)
221 IF(VINT(33).LE.VINT(13).AND.VINT(34).LE.VINT(14)) MINT(51)=1
222
223 ELSEIF(ILIM.EQ.4) THEN
224C...Calculate limits on tau'
225C...0) due to kinematics
226 TAPMN0=TAU
227 TAPMX0=1.
228C...1) due to explicit limits
229 TAPMN1=CKIN(31)**2/VINT(2)
230 TAPMX1=1.
231 IF(CKIN(32).GE.0.) TAPMX1=CKIN(32)**2/VINT(2)
232 VINT(16)=MAX(TAPMN0,TAPMN1)
233 VINT(36)=MIN(TAPMX0,TAPMX1)
234 IF(MINT(43).EQ.1) THEN
235 VINT(16)=0.99999
236 VINT(36)=1.00001
237 ENDIF
238 IF(VINT(36).LE.VINT(16)) MINT(51)=1
239
240 ENDIF
241 RETURN
242
243C...Special case for low-pT and multiple interactions:
244C...effective kinematical limits for tau, y*, cos(theta-hat).
245 110 IF(ILIM.EQ.0) THEN
246 ELSEIF(ILIM.EQ.1) THEN
247 IF(MSTP(82).LE.1) VINT(11)=4.*PARP(81)**2/VINT(2)
248 IF(MSTP(82).GE.2) VINT(11)=PARP(82)**2/VINT(2)
249 VINT(31)=1.
250 ELSEIF(ILIM.EQ.2) THEN
251 VINT(12)=0.5*LOG(VINT(21))
252 VINT(32)=-VINT(12)
253 ELSEIF(ILIM.EQ.3) THEN
254 IF(MSTP(82).LE.1) ST2EFF=4.*PARP(81)**2/(VINT(21)*VINT(2))
255 IF(MSTP(82).GE.2) ST2EFF=0.01*PARP(82)**2/(VINT(21)*VINT(2))
256 VINT(13)=-SQRT(MAX(0.,1.-ST2EFF))
257 VINT(33)=0.
258 VINT(14)=0.
259 VINT(34)=-VINT(13)
260 ENDIF
261
262 RETURN
263 END