]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | |
2 | C*********************************************************************** | |
3 | ||
4 | SUBROUTINE PYOFSH(MOFSH,KFMO,KFD1,KFD2,PMMO,RET1,RET2) | |
5 | ||
6 | C...Calculates partial width and differential cross-section maxima | |
7 | C...of channels/processes not allowed on mass-shell, and selects | |
8 | C...masses in such channels/processes. | |
9 | COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200) | |
10 | COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4) | |
11 | COMMON/LUDAT3/MDCY(500,3),MDME(2000,2),BRAT(2000),KFDP(2000,5) | |
12 | COMMON/PYSUBS/MSEL,MSUB(200),KFIN(2,-40:40),CKIN(200) | |
13 | COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200) | |
14 | COMMON/PYINT1/MINT(400),VINT(400) | |
15 | COMMON/PYINT2/ISET(200),KFPR(200,2),COEF(200,20),ICOL(40,4,2) | |
16 | COMMON/PYINT5/NGEN(0:200,3),XSEC(0:200,3) | |
17 | SAVE /LUDAT1/,/LUDAT2/,/LUDAT3/ | |
18 | SAVE /PYSUBS/,/PYPARS/,/PYINT1/,/PYINT2/,/PYINT5/ | |
19 | DIMENSION KFD(2),MBW(2),PMD(2),PGD(2),PMG(2),PML(2),PMU(2), | |
20 | &PMH(2),ATL(2),ATU(2),ATH(2),RMG(2),INX1(100),XPT1(100), | |
21 | &FPT1(100),INX2(100),XPT2(100),FPT2(100),WDTP(0:40), | |
22 | &WDTE(0:40,0:5) | |
23 | ||
24 | C...Find if particles equal, maximum mass, matrix elements, etc. | |
25 | MINT(51)=0 | |
26 | ISUB=MINT(1) | |
27 | KFD(1)=IABS(KFD1) | |
28 | KFD(2)=IABS(KFD2) | |
29 | MEQL=0 | |
30 | IF(KFD(1).EQ.KFD(2)) MEQL=1 | |
31 | MLM=0 | |
32 | IF(MOFSH.GE.2.AND.MEQL.EQ.1) MLM=INT(1.5+RLU(0)) | |
33 | IF(MOFSH.LE.2.OR.MOFSH.EQ.7) THEN | |
34 | NOFF=44 | |
35 | PMMX=PMMO | |
36 | ELSE | |
37 | NOFF=40 | |
38 | PMMX=VINT(1) | |
39 | IF(CKIN(2).GT.CKIN(1)) PMMX=MIN(CKIN(2),VINT(1)) | |
40 | ENDIF | |
41 | MMED=0 | |
42 | IF((KFMO.EQ.25.OR.KFMO.EQ.35.OR.KFMO.EQ.36).AND.MEQL.EQ.1.AND. | |
43 | &(KFD(1).EQ.23.OR.KFD(1).EQ.24)) MMED=1 | |
44 | IF((KFMO.EQ.32.OR.IABS(KFMO).EQ.34).AND.(KFD(1).EQ.23.OR. | |
45 | &KFD(1).EQ.24).AND.(KFD(2).EQ.23.OR.KFD(2).EQ.24)) MMED=2 | |
46 | IF((KFMO.EQ.32.OR.IABS(KFMO).EQ.34).AND.(KFD(2).EQ.25.OR. | |
47 | &KFD(2).EQ.35.OR.KFD(2).EQ.36)) MMED=3 | |
48 | LOOP=1 | |
49 | ||
50 | C...Find where Breit-Wigners are required, else select discrete masses. | |
51 | 100 DO 110 I=1,2 | |
52 | KFCA=KFD(I) | |
53 | IF(KFCA.GT.100) KFCA=LUCOMP(KFCA) | |
54 | IF(KFCA.GT.0) THEN | |
55 | PMD(I)=PMAS(KFCA,1) | |
56 | PGD(I)=PMAS(KFCA,2) | |
57 | ELSE | |
58 | PMD(I)=0. | |
59 | PGD(I)=0. | |
60 | ENDIF | |
61 | IF(MSTP(42).LE.0.OR.PGD(I).LT.PARP(41)) THEN | |
62 | MBW(I)=0 | |
63 | PMG(I)=PMD(I) | |
64 | RMG(I)=(PMG(I)/PMMX)**2 | |
65 | ELSE | |
66 | MBW(I)=1 | |
67 | ENDIF | |
68 | 110 CONTINUE | |
69 | ||
70 | C...Find allowed mass range and Breit-Wigner parameters. | |
71 | DO 120 I=1,2 | |
72 | IF(MOFSH.EQ.1.AND.LOOP.EQ.1.AND.MBW(I).EQ.1) THEN | |
73 | PML(I)=PARP(42) | |
74 | PMU(I)=PMMX-PARP(42) | |
75 | IF(MBW(3-I).EQ.0) PMU(I)=MIN(PMU(I),PMMX-PMD(3-I)) | |
76 | IF(PMU(I).LT.PML(I)+PARJ(64)) MBW(I)=-1 | |
77 | ELSEIF((MBW(I).EQ.1.OR.MOFSH.GE.5).AND.MOFSH.NE.7) THEN | |
78 | ILM=I | |
79 | IF(MLM.EQ.2) ILM=3-I | |
80 | PML(I)=MAX(CKIN(NOFF+2*ILM-1),PARP(42)) | |
81 | IF(MOFSH.GE.5.AND.I.EQ.2) PML(I)=MAX(PML(I),2.*PMAS(KFD2,1)) | |
82 | PMU(I)=PMMX-MAX(CKIN(NOFF+5-2*ILM),PARP(42)) | |
83 | IF(MOFSH.GE.5.AND.I.EQ.1) PMU(I)=MIN(PMU(I),PMMX-2.* | |
84 | & PMAS(KFD2,1)) | |
85 | IF(CKIN(NOFF+2*ILM).GT.CKIN(NOFF+2*ILM-1)) PMU(I)=MIN(PMU(I), | |
86 | & CKIN(NOFF+2*ILM)) | |
87 | IF(MBW(3-I).EQ.0) PMU(I)=MIN(PMU(I),PMMX-PMD(3-I)) | |
88 | IF(I.EQ.MLM) PMU(I)=MIN(PMU(I),0.5*PMMX) | |
89 | IF(MEQL.EQ.0) PMH(I)=MIN(PMU(I),0.5*PMMX) | |
90 | IF(PMU(I).LT.PML(I)+PARJ(64)) MBW(I)=-1 | |
91 | IF(MBW(I).EQ.1) THEN | |
92 | ATL(I)=ATAN((PML(I)**2-PMD(I)**2)/(PMD(I)*PGD(I))) | |
93 | ATU(I)=ATAN((PMU(I)**2-PMD(I)**2)/(PMD(I)*PGD(I))) | |
94 | IF(MEQL.EQ.0) ATH(I)=ATAN((PMH(I)**2-PMD(I)**2)/(PMD(I)* | |
95 | & PGD(I))) | |
96 | ENDIF | |
97 | ELSEIF(MBW(I).EQ.1.AND.MOFSH.EQ.7) THEN | |
98 | ILM=I | |
99 | IF(MLM.EQ.2) ILM=3-I | |
100 | PML(I)=PARP(42) | |
101 | PMU(I)=PMMX-PARP(42) | |
102 | IF(MBW(3-I).EQ.0) PMU(I)=MIN(PMU(I),PMMX-PMD(3-I)) | |
103 | IF(I.EQ.MLM) PMU(I)=MIN(PMU(I),0.5*PMMX) | |
104 | IF(MEQL.EQ.0) PMH(I)=MIN(PMU(I),0.5*PMMX) | |
105 | IF(PMU(I).LT.PML(I)+PARJ(64)) MBW(I)=-1 | |
106 | IF(MBW(I).EQ.1) THEN | |
107 | ATL(I)=ATAN((PML(I)**2-PMD(I)**2)/(PMD(I)*PGD(I))) | |
108 | ATU(I)=ATAN((PMU(I)**2-PMD(I)**2)/(PMD(I)*PGD(I))) | |
109 | IF(MEQL.EQ.0) ATH(I)=ATAN((PMH(I)**2-PMD(I)**2)/(PMD(I)* | |
110 | & PGD(I))) | |
111 | ENDIF | |
112 | ENDIF | |
113 | 120 CONTINUE | |
114 | IF(MBW(1).LT.0.OR.MBW(2).LT.0.OR.(MBW(1).EQ.0.AND.MBW(2).EQ.0)) | |
115 | &THEN | |
116 | CALL LUERRM(13,'(PYOFSH:) no allowed decay product masses') | |
117 | MINT(51)=1 | |
118 | RETURN | |
119 | ENDIF | |
120 | ||
121 | C...Calculation of partial width of resonance. | |
122 | IF(MOFSH.EQ.1) THEN | |
123 | ||
124 | C..If only one integration, pick that to be the inner. | |
125 | IF(MBW(1).EQ.0) THEN | |
126 | PM2=PMD(1) | |
127 | PMD(1)=PMD(2) | |
128 | PGD(1)=PGD(2) | |
129 | PML(1)=PML(2) | |
130 | PMU(1)=PMU(2) | |
131 | ELSEIF(MBW(2).EQ.0) THEN | |
132 | PM2=PMD(2) | |
133 | ENDIF | |
134 | ||
135 | C...Start outer loop of integration. | |
136 | IF(MBW(1).EQ.1.AND.MBW(2).EQ.1) THEN | |
137 | ATL2=ATAN((PML(2)**2-PMD(2)**2)/(PMD(2)*PGD(2))) | |
138 | ATU2=ATAN((PMU(2)**2-PMD(2)**2)/(PMD(2)*PGD(2))) | |
139 | NPT2=1 | |
140 | XPT2(1)=1. | |
141 | INX2(1)=0 | |
142 | FMAX2=0. | |
143 | ENDIF | |
144 | 130 IF(MBW(1).EQ.1.AND.MBW(2).EQ.1) THEN | |
145 | PM2S=PMD(2)**2+PMD(2)*PGD(2)*TAN(ATL2+XPT2(NPT2)*(ATU2-ATL2)) | |
146 | PM2=MIN(PMU(2),MAX(PML(2),SQRT(MAX(0.,PM2S)))) | |
147 | ENDIF | |
148 | RM2=(PM2/PMMX)**2 | |
149 | ||
150 | C...Start inner loop of integration. | |
151 | PML1=PML(1) | |
152 | PMU1=MIN(PMU(1),PMMX-PM2) | |
153 | IF(MEQL.EQ.1) PMU1=MIN(PMU1,PM2) | |
154 | ATL1=ATAN((PML1**2-PMD(1)**2)/(PMD(1)*PGD(1))) | |
155 | ATU1=ATAN((PMU1**2-PMD(1)**2)/(PMD(1)*PGD(1))) | |
156 | IF(PML1+PARJ(64).GE.PMU1.OR.ATL1+1E-7.GE.ATU1) THEN | |
157 | FUNC2=0. | |
158 | GOTO 180 | |
159 | ENDIF | |
160 | NPT1=1 | |
161 | XPT1(1)=1. | |
162 | INX1(1)=0 | |
163 | FMAX1=0. | |
164 | 140 PM1S=PMD(1)**2+PMD(1)*PGD(1)*TAN(ATL1+XPT1(NPT1)*(ATU1-ATL1)) | |
165 | PM1=MIN(PMU1,MAX(PML1,SQRT(MAX(0.,PM1S)))) | |
166 | RM1=(PM1/PMMX)**2 | |
167 | ||
168 | C...Evaluate function value - inner loop. | |
169 | FUNC1=SQRT(MAX(0.,(1.-RM1-RM2)**2-4.*RM1*RM2)) | |
170 | IF(MMED.EQ.1) FUNC1=FUNC1*((1.-RM1-RM2)**2+8.*RM1*RM2) | |
171 | IF(MMED.EQ.2) FUNC1=FUNC1**3*(1.+10.*RM1+10.*RM2+RM1**2+ | |
172 | & RM2**2+10.*RM1*RM2) | |
173 | IF(FUNC1.GT.FMAX1) FMAX1=FUNC1 | |
174 | FPT1(NPT1)=FUNC1 | |
175 | ||
176 | C...Go to next position in inner loop. | |
177 | IF(NPT1.EQ.1) THEN | |
178 | NPT1=NPT1+1 | |
179 | XPT1(NPT1)=0. | |
180 | INX1(NPT1)=1 | |
181 | GOTO 140 | |
182 | ELSEIF(NPT1.LE.8) THEN | |
183 | NPT1=NPT1+1 | |
184 | IF(NPT1.LE.4.OR.NPT1.EQ.6) ISH1=1 | |
185 | ISH1=ISH1+1 | |
186 | XPT1(NPT1)=0.5*(XPT1(ISH1)+XPT1(INX1(ISH1))) | |
187 | INX1(NPT1)=INX1(ISH1) | |
188 | INX1(ISH1)=NPT1 | |
189 | GOTO 140 | |
190 | ELSEIF(NPT1.LT.100) THEN | |
191 | ISN1=ISH1 | |
192 | 150 ISH1=ISH1+1 | |
193 | IF(ISH1.GT.NPT1) ISH1=2 | |
194 | IF(ISH1.EQ.ISN1) GOTO 160 | |
195 | DFPT1=ABS(FPT1(ISH1)-FPT1(INX1(ISH1))) | |
196 | IF(DFPT1.LT.PARP(43)*FMAX1) GOTO 150 | |
197 | NPT1=NPT1+1 | |
198 | XPT1(NPT1)=0.5*(XPT1(ISH1)+XPT1(INX1(ISH1))) | |
199 | INX1(NPT1)=INX1(ISH1) | |
200 | INX1(ISH1)=NPT1 | |
201 | GOTO 140 | |
202 | ENDIF | |
203 | ||
204 | C...Calculate integral over inner loop. | |
205 | 160 FSUM1=0. | |
206 | DO 170 IPT1=2,NPT1 | |
207 | FSUM1=FSUM1+0.5*(FPT1(IPT1)+FPT1(INX1(IPT1)))* | |
208 | & (XPT1(INX1(IPT1))-XPT1(IPT1)) | |
209 | 170 CONTINUE | |
210 | FUNC2=FSUM1*(ATU1-ATL1)/PARU(1) | |
211 | 180 IF(MBW(1).EQ.1.AND.MBW(2).EQ.1) THEN | |
212 | IF(FUNC2.GT.FMAX2) FMAX2=FUNC2 | |
213 | FPT2(NPT2)=FUNC2 | |
214 | ||
215 | C...Go to next position in outer loop. | |
216 | IF(NPT2.EQ.1) THEN | |
217 | NPT2=NPT2+1 | |
218 | XPT2(NPT2)=0. | |
219 | INX2(NPT2)=1 | |
220 | GOTO 130 | |
221 | ELSEIF(NPT2.LE.8) THEN | |
222 | NPT2=NPT2+1 | |
223 | IF(NPT2.LE.4.OR.NPT2.EQ.6) ISH2=1 | |
224 | ISH2=ISH2+1 | |
225 | XPT2(NPT2)=0.5*(XPT2(ISH2)+XPT2(INX2(ISH2))) | |
226 | INX2(NPT2)=INX2(ISH2) | |
227 | INX2(ISH2)=NPT2 | |
228 | GOTO 130 | |
229 | ELSEIF(NPT2.LT.100) THEN | |
230 | ISN2=ISH2 | |
231 | 190 ISH2=ISH2+1 | |
232 | IF(ISH2.GT.NPT2) ISH2=2 | |
233 | IF(ISH2.EQ.ISN2) GOTO 200 | |
234 | DFPT2=ABS(FPT2(ISH2)-FPT2(INX2(ISH2))) | |
235 | IF(DFPT2.LT.PARP(43)*FMAX2) GOTO 190 | |
236 | NPT2=NPT2+1 | |
237 | XPT2(NPT2)=0.5*(XPT2(ISH2)+XPT2(INX2(ISH2))) | |
238 | INX2(NPT2)=INX2(ISH2) | |
239 | INX2(ISH2)=NPT2 | |
240 | GOTO 130 | |
241 | ENDIF | |
242 | ||
243 | C...Calculate integral over outer loop. | |
244 | 200 FSUM2=0. | |
245 | DO 210 IPT2=2,NPT2 | |
246 | FSUM2=FSUM2+0.5*(FPT2(IPT2)+FPT2(INX2(IPT2)))* | |
247 | & (XPT2(INX2(IPT2))-XPT2(IPT2)) | |
248 | 210 CONTINUE | |
249 | FSUM2=FSUM2*(ATU2-ATL2)/PARU(1) | |
250 | IF(MEQL.EQ.1) FSUM2=2.*FSUM2 | |
251 | ELSE | |
252 | FSUM2=FUNC2 | |
253 | ENDIF | |
254 | ||
255 | C...Save result; second integration for user-selected mass range. | |
256 | IF(LOOP.EQ.1) WIDW=FSUM2 | |
257 | WID2=FSUM2 | |
258 | IF(LOOP.EQ.1.AND.(CKIN(46).GE.CKIN(45).OR.CKIN(48).GE.CKIN(47) | |
259 | & .OR.MAX(CKIN(45),CKIN(47)).GE.1.01*PARP(42))) THEN | |
260 | LOOP=2 | |
261 | GOTO 100 | |
262 | ENDIF | |
263 | RET1=WIDW | |
264 | RET2=WID2/WIDW | |
265 | ||
266 | C...Select two decay product masses of a resonance. | |
267 | ELSEIF(MOFSH.EQ.2.OR.MOFSH.EQ.7) THEN | |
268 | 220 DO 230 I=1,2 | |
269 | IF(MBW(I).EQ.0) GOTO 230 | |
270 | PMBW=PMD(I)**2+PMD(I)*PGD(I)*TAN(ATL(I)+RLU(0)*(ATU(I)-ATL(I))) | |
271 | PMG(I)=MIN(PMU(I),MAX(PML(I),SQRT(MAX(0.,PMBW)))) | |
272 | RMG(I)=(PMG(I)/PMMX)**2 | |
273 | 230 CONTINUE | |
274 | IF((MEQL.EQ.1.AND.PMG(MAX(1,MLM)).GT.PMG(MIN(2,3-MLM))).OR. | |
275 | & PMG(1)+PMG(2)+PARJ(64).GT.PMMX) GOTO 220 | |
276 | ||
277 | C...Weight with matrix element (if none known, use beta factor). | |
278 | FLAM=SQRT(MAX(0.,(1.-RMG(1)-RMG(2))**2-4.*RMG(1)*RMG(2))) | |
279 | IF(MMED.EQ.1) THEN | |
280 | WTBE=FLAM*((1.-RMG(1)-RMG(2))**2+8.*RMG(1)*RMG(2)) | |
281 | ELSEIF(MMED.EQ.2) THEN | |
282 | WTBE=FLAM**3*(1.+10.*RMG(1)+10.*RMG(2)+RMG(1)**2+ | |
283 | & RMG(2)**2+10.*RMG(1)*RMG(2)) | |
284 | ELSEIF(MMED.EQ.3) THEN | |
285 | WTBE=FLAM*(RMG(1)+FLAM**2/12.) | |
286 | ELSE | |
287 | WTBE=FLAM | |
288 | ENDIF | |
289 | IF(WTBE.LT.RLU(0)) GOTO 220 | |
290 | RET1=PMG(1) | |
291 | RET2=PMG(2) | |
292 | ||
293 | C...Find suitable set of masses for initialization of 2 -> 2 processes. | |
294 | ELSEIF(MOFSH.EQ.3) THEN | |
295 | IF(MBW(1).NE.0.AND.MBW(2).EQ.0) THEN | |
296 | PMG(1)=MIN(PMD(1),0.5*(PML(1)+PMU(1))) | |
297 | PMG(2)=PMD(2) | |
298 | ELSEIF(MBW(2).NE.0.AND.MBW(1).EQ.0) THEN | |
299 | PMG(1)=PMD(1) | |
300 | PMG(2)=MIN(PMD(2),0.5*(PML(2)+PMU(2))) | |
301 | ELSE | |
302 | IDIV=-1 | |
303 | 240 IDIV=IDIV+1 | |
304 | PMG(1)=MIN(PMD(1),0.1*(IDIV*PML(1)+(10-IDIV)*PMU(1))) | |
305 | PMG(2)=MIN(PMD(2),0.1*(IDIV*PML(2)+(10-IDIV)*PMU(2))) | |
306 | IF(IDIV.LE.9.AND.PMG(1)+PMG(2).GT.0.9*PMMX) GOTO 240 | |
307 | ENDIF | |
308 | RET1=PMG(1) | |
309 | RET2=PMG(2) | |
310 | ||
311 | C...Evaluate importance of excluded tails of Breit-Wigners. | |
312 | IF(MEQL.EQ.0.AND.MBW(1).EQ.1.AND.MBW(2).EQ.1.AND.PMD(1)+PMD(2). | |
313 | & GT.PMMX.AND.PMH(1).GT.PML(1).AND.PMH(2).GT.PML(2)) MEQL=2 | |
314 | IF(MEQL.LE.1) THEN | |
315 | VINT(80)=1. | |
316 | DO 250 I=1,2 | |
317 | IF(MBW(I).NE.0) VINT(80)=VINT(80)*1.25*(ATU(I)-ATL(I))/PARU(1) | |
318 | 250 CONTINUE | |
319 | ELSE | |
320 | VINT(80)=(1.25/PARU(1))**2*MAX((ATU(1)-ATL(1))* | |
321 | & (ATH(2)-ATL(2)),(ATH(1)-ATL(1))*(ATU(2)-ATL(2))) | |
322 | ENDIF | |
323 | IF((ISUB.EQ.15.OR.ISUB.EQ.19.OR.ISUB.EQ.30.OR.ISUB.EQ.35).AND. | |
324 | & MSTP(43).NE.2) VINT(80)=2.*VINT(80) | |
325 | IF(ISUB.EQ.22.AND.MSTP(43).NE.2) VINT(80)=4.*VINT(80) | |
326 | IF(MEQL.GE.1) VINT(80)=2.*VINT(80) | |
327 | ||
328 | C...Pick one particle to be the lighter (if improves efficiency). | |
329 | ELSEIF(MOFSH.EQ.4) THEN | |
330 | IF(MEQL.EQ.0.AND.MBW(1).EQ.1.AND.MBW(2).EQ.1.AND.PMD(1)+PMD(2). | |
331 | & GT.PMMX.AND.PMH(1).GT.PML(1).AND.PMH(2).GT.PML(2)) MEQL=2 | |
332 | 260 IF(MEQL.EQ.2) MLM=INT(1.5+RLU(0)) | |
333 | ||
334 | C...Select two masses according to Breit-Wigner + flat in s + 1/s. | |
335 | DO 270 I=1,2 | |
336 | IF(MBW(I).EQ.0) GOTO 270 | |
337 | PMV=PMU(I) | |
338 | IF(MEQL.EQ.2.AND.I.EQ.MLM) PMV=PMH(I) | |
339 | ATV=ATU(I) | |
340 | IF(MEQL.EQ.2.AND.I.EQ.MLM) ATV=ATH(I) | |
341 | RBR=RLU(0) | |
342 | IF((ISUB.EQ.15.OR.ISUB.EQ.19.OR.ISUB.EQ.22.OR.ISUB.EQ.30.OR. | |
343 | & ISUB.EQ.35).AND.MSTP(43).NE.2) RBR=2.*RBR | |
344 | IF(RBR.LT.0.8) THEN | |
345 | PMSR=PMD(I)**2+PMD(I)*PGD(I)*TAN(ATL(I)+RLU(0)*(ATV-ATL(I))) | |
346 | PMG(I)=MIN(PMV,MAX(PML(I),SQRT(MAX(0.,PMSR)))) | |
347 | ELSEIF(RBR.LT.0.9) THEN | |
348 | PMG(I)=SQRT(MAX(0.,PML(I)**2+RLU(0)*(PMV**2-PML(I)**2))) | |
349 | ELSEIF(RBR.LT.1.5) THEN | |
350 | PMG(I)=PML(I)*(PMV/PML(I))**RLU(0) | |
351 | ELSE | |
352 | PMG(I)=SQRT(MAX(0.,PML(I)**2*PMV**2/(PML(I)**2+RLU(0)* | |
353 | & (PMV**2-PML(I)**2)))) | |
354 | ENDIF | |
355 | 270 CONTINUE | |
356 | IF((MEQL.GE.1.AND.PMG(MAX(1,MLM)).GT.PMG(MIN(2,3-MLM))).OR. | |
357 | & PMG(1)+PMG(2)+PARJ(64).GT.PMMX) THEN | |
358 | IF(MINT(48).EQ.1) THEN | |
359 | NGEN(0,1)=NGEN(0,1)+1 | |
360 | NGEN(MINT(1),1)=NGEN(MINT(1),1)+1 | |
361 | GOTO 260 | |
362 | ELSE | |
363 | MINT(51)=1 | |
364 | RETURN | |
365 | ENDIF | |
366 | ENDIF | |
367 | RET1=PMG(1) | |
368 | RET2=PMG(2) | |
369 | ||
370 | C...Give weight for selected mass distribution. | |
371 | VINT(80)=1. | |
372 | DO 280 I=1,2 | |
373 | IF(MBW(I).EQ.0) GOTO 280 | |
374 | PMV=PMU(I) | |
375 | IF(MEQL.EQ.2.AND.I.EQ.MLM) PMV=PMH(I) | |
376 | ATV=ATU(I) | |
377 | IF(MEQL.EQ.2.AND.I.EQ.MLM) ATV=ATH(I) | |
378 | F0=PMD(I)*PGD(I)/((PMG(I)**2-PMD(I)**2)**2+ | |
379 | & (PMD(I)*PGD(I))**2)/PARU(1) | |
380 | F1=1. | |
381 | F2=1./PMG(I)**2 | |
382 | F3=1./PMG(I)**4 | |
383 | FI0=(ATV-ATL(I))/PARU(1) | |
384 | FI1=PMV**2-PML(I)**2 | |
385 | FI2=2.*LOG(PMV/PML(I)) | |
386 | FI3=1./PML(I)**2-1./PMV**2 | |
387 | IF((ISUB.EQ.15.OR.ISUB.EQ.19.OR.ISUB.EQ.22.OR.ISUB.EQ.30.OR. | |
388 | & ISUB.EQ.35).AND.MSTP(43).NE.2) THEN | |
389 | VINT(80)=VINT(80)*20./(8.+(FI0/F0)*(F1/FI1+6.*F2/FI2+ | |
390 | & 5.*F3/FI3)) | |
391 | ELSE | |
392 | VINT(80)=VINT(80)*10./(8.+(FI0/F0)*(F1/FI1+F2/FI2)) | |
393 | ENDIF | |
394 | VINT(80)=VINT(80)*FI0 | |
395 | 280 CONTINUE | |
396 | IF(MEQL.GE.1) VINT(80)=2.*VINT(80) | |
397 | ||
398 | ELSEIF(MOFSH.EQ.5) THEN | |
399 | C...Find suitable set of masses for initialization of 2 -> 3 process. | |
400 | IDIV=6 | |
401 | 290 IDIV=IDIV-1 | |
402 | IF(MBW(1).EQ.0) THEN | |
403 | PMG(1)=PMD(1) | |
404 | ELSE | |
405 | PMSR=PMD(1)**2+PMD(1)*PGD(1)*TAN(ATL(1)+0.1*IDIV*(ATU(1)- | |
406 | & ATL(1))) | |
407 | PMG(1)=MIN(PMU(1),MAX(PML(1),SQRT(MAX(0.,PMSR)))) | |
408 | ENDIF | |
409 | PMG(2)=PML(2)*(PMU(2)/PML(2))**(0.1*IDIV) | |
410 | IF(IDIV.GE.1.AND.PMG(1)+PMG(2).GT.0.9*PMMX) GOTO 290 | |
411 | RET1=PMG(1) | |
412 | RET2=PMG(2) | |
413 | ||
414 | C...Evaluate size of selected phase space volume. | |
415 | VINT(80)=2.*LOG(PMU(2)/PML(2)) | |
416 | IF(MBW(1).NE.0) VINT(80)=VINT(80)*1.25*(ATU(1)-ATL(1))/PARU(1) | |
417 | ||
418 | C...Pick decay angles. | |
419 | VINT(81)=0. | |
420 | VINT(82)=0.5*PARU(1) | |
421 | VINT(83)=1. | |
422 | VINT(84)=0. | |
423 | ||
424 | C...Select flavour of resonance decays. | |
425 | KFA=KFPR(ISUB,1) | |
426 | CALL PYWIDT(KFA,PMG(1)**2,WDTP,WDTE) | |
427 | IF(KCHG(KFA,3).EQ.0) THEN | |
428 | IPM=2 | |
429 | ELSE | |
430 | IPM=(5-ISIGN(1,KFA))/2 | |
431 | ENDIF | |
432 | WDTE0S=WDTE(0,1)+WDTE(0,IPM)+WDTE(0,4) | |
433 | IF(WDTE0S.LE.0.) THEN | |
434 | CALL LUERRM(12,'(PYOFSH:) no allowed resonace decay channel') | |
435 | MINT(51)=1 | |
436 | RETURN | |
437 | ENDIF | |
438 | WDTEC=0. | |
439 | DO 300 IDL=1,MDCY(KFA,3) | |
440 | WDTEK=WDTE(IDL,1)+WDTE(IDL,IPM)+WDTE(IDL,4) | |
441 | IF(WDTEK.GT.WDTEC) THEN | |
442 | IDC=IDL+MDCY(KFA,2)-1 | |
443 | WDTEC=WDTEK | |
444 | ENDIF | |
445 | 300 CONTINUE | |
446 | MINT(35)=IDC | |
447 | ||
448 | C...Compensating factor for all flavours. | |
449 | KFL=IABS(KFDP(IDC,1)) | |
450 | QFL=KCHG(KFL,1)/3. | |
451 | AFL=SIGN(1.,QFL+0.1) | |
452 | VFL=AFL-4.*PARU(102)*QFL | |
453 | WDTEK=VFL**2+AFL**2 | |
454 | VINT(80)=VINT(80)*WDTE0S/WDTEK | |
455 | ||
456 | ELSEIF(MOFSH.EQ.6) THEN | |
457 | C...Select two masses, one basically Breit-Wigner, other dm^2/m^2. | |
458 | IF(MBW(1).NE.0) THEN | |
459 | RBR=RLU(0) | |
460 | IF(RBR.LT.0.8) THEN | |
461 | PMSR=PMD(1)**2+PMD(1)*PGD(1)*TAN(ATL(1)+RLU(0)* | |
462 | & (ATU(1)-ATL(1))) | |
463 | PMG(1)=MIN(PMU(1),MAX(PML(1),SQRT(MAX(0.,PMSR)))) | |
464 | ELSEIF(RBR.LT.0.9) THEN | |
465 | PMG(1)=SQRT(MAX(0.,PML(1)**2+RLU(0)*(PMU(1)**2-PML(1)**2))) | |
466 | ELSE | |
467 | PMG(1)=PML(1)*(PMU(1)/PML(1))**RLU(0) | |
468 | ENDIF | |
469 | ENDIF | |
470 | PMG(2)=PML(2)*(PMU(2)/PML(2))**RLU(0) | |
471 | IF(SQRT(MAX(0.,1.-(PML(2)/PMG(2))**2)).LT.RLU(0).OR. | |
472 | & PMG(1)+PMG(2)+PARJ(64).GT.PMMX) THEN | |
473 | MINT(51)=1 | |
474 | RETURN | |
475 | ENDIF | |
476 | RET1=PMG(1) | |
477 | RET2=PMG(2) | |
478 | ||
479 | C...Give weight for selected mass distribution. | |
480 | VINT(80)=2.*LOG(PMU(2)/PML(2)) | |
481 | IF(MBW(1).NE.0) THEN | |
482 | F0=PMD(1)*PGD(1)/((PMG(1)**2-PMD(1)**2)**2+ | |
483 | & (PMD(1)*PGD(1))**2)/PARU(1) | |
484 | F1=1. | |
485 | F2=1./PMG(1)**2 | |
486 | FI0=(ATU(1)-ATL(1))/PARU(1) | |
487 | FI1=PMU(1)**2-PML(1)**2 | |
488 | FI2=2.*LOG(PMU(1)/PML(1)) | |
489 | VINT(80)=VINT(80)*10.*FI0/(8.+(FI0/F0)*(F1/FI1+F2/FI2)) | |
490 | ENDIF | |
491 | ||
492 | C...Select decay angles. | |
493 | VINT(81)=2.*RLU(0)-1. | |
494 | VINT(82)=PARU(2)*RLU(0) | |
495 | VINT(83)=2.*RLU(0)-1. | |
496 | VINT(84)=PARU(2)*RLU(0) | |
497 | ||
498 | C...Select flavour of resonance decays. | |
499 | KFA=KFPR(ISUB,1) | |
500 | CALL PYWIDT(KFA,PMG(1)**2,WDTP,WDTE) | |
501 | IF(KCHG(KFA,3).EQ.0) THEN | |
502 | IPM=2 | |
503 | ELSE | |
504 | IPM=(5-ISIGN(1,KFA))/2 | |
505 | ENDIF | |
506 | WDTE0S=WDTE(0,1)+WDTE(0,IPM)+WDTE(0,4) | |
507 | IF(WDTE0S.LE.0.) THEN | |
508 | CALL LUERRM(12,'(PYOFSH:) no allowed resonace decay channel') | |
509 | MINT(51)=1 | |
510 | RETURN | |
511 | ENDIF | |
512 | RKFL=WDTE0S*RLU(0) | |
513 | IDL=0 | |
514 | 310 IDL=IDL+1 | |
515 | IDC=IDL+MDCY(KFA,2)-1 | |
516 | RKFL=RKFL-(WDTE(IDL,1)+WDTE(IDL,IPM)+WDTE(IDL,4)) | |
517 | IF(IDL.LT.MDCY(KFA,3).AND.RKFL.GT.0.) GOTO 310 | |
518 | MINT(35)=IDC | |
519 | ||
520 | C...Compensating factor for all flavours. | |
521 | KFL=IABS(KFDP(IDC,1)) | |
522 | QFL=KCHG(KFL,1)/3. | |
523 | AFL=SIGN(1.,QFL+0.1) | |
524 | VFL=AFL-4.*PARU(102)*QFL | |
525 | WDTEK=VFL**2+AFL**2 | |
526 | VINT(80)=VINT(80)*WDTE0S/WDTEK | |
527 | ENDIF | |
528 | ||
529 | RETURN | |
530 | END |