Double check if SM is running added. Some redundant output removed from SM
[u/mrichter/AliRoot.git] / DPMJET / pythia6115dpm3.f
CommitLineData
9aaba0d6 1C*********************************************************************
2C*********************************************************************
3C* **
4C* March 1997 **
5C* **
6C* The Lund Monte Carlo for Hadronic Processes **
7C* **
8C* PYTHIA version 6.1 **
9C* **
10C* Torbjorn Sjostrand **
11C* Department of Theoretical Physics 2 **
12C* Lund University **
13C* Solvegatan 14A, S-223 62 Lund, Sweden **
14C* phone +46 - 46 - 222 48 16 **
15C* E-mail torbjorn@thep.lu.se **
16C* **
17C* SUSY parts by **
18C* Stephen Mrenna **
19C* Argonne National Laboratory **
20C* 9700 South Cass Avenue, Argonne, IL 60439, USA **
21C* phone + 1 - 630 - 252 - 7615 **
22C* E-mail mrenna@hep.anl.gov **
23C* **
24C* Several parts are written by Hans-Uno Bengtsson **
25C* PYSHOW is written together with Mats Bengtsson **
26C* advanced popcorn baryon production written by Patrik Eden **
27C* CTEQ 3 parton distributions are by the CTEQ collaboration **
28C* GRV 94 parton distributions are by Glueck, Reya and Vogt **
29C* SaS photon parton distributions together with Gerhard Schuler **
30C* g + g and q + qbar -> t + tbar + H code by Zoltan Kunszt **
31C* MSSM Higgs mass calculation code by M. Carena, **
32C* J.R. Espinosa, M. Quiros and C.E.M. Wagner **
33C* PYGAUS adapted from CERN library (K.S. Kolbig) **
34C* **
35C* The latest program version and documentation is found on WWW **
36C* http://www.thep.lu.se/tf2/staff/torbjorn/Pythia.html **
37C* **
38C* Copyright Torbjorn Sjostrand, Lund 1997 **
39C* **
40C*********************************************************************
41C*********************************************************************
42C *
43C List of subprograms in order of appearance, with main purpose *
44C (S = subroutine, F = function, B = block data) *
45C *
46C B PYDATA to contain all default values *
47C S PYTEST to test the proper functioning of the package *
48C S PYHEPC to convert between /PYJETS/ and /HEPEVT/ records *
49C *
50C S PYINIT to administer the initialization procedure *
51C S PYEVNT to administer the generation of an event *
52C S PYSTAT to print cross-section and other information *
53C S PYINRE to initialize treatment of resonances *
54C S PYINBM to read in beam, target and frame choices *
55C S PYINKI to initialize kinematics of incoming particles *
56C S PYINPR to set up the selection of included processes *
57C S PYXTOT to give total, elastic and diffractive cross-sect. *
58C S PYMAXI to find differential cross-section maxima *
59C S PYPILE to select multiplicity of pileup events *
60C S PYSAVE to save alternatives for gamma-p and gamma-gamma *
61C S PYRAND to select subprocess and kinematics for event *
62C S PYSCAT to set up kinematics and colour flow of event *
63C S PYSSPA to simulate initial state spacelike showers *
64C S PYRESD to perform resonance decays *
65C S PYMULT to generate multiple interactions *
66C S PYREMN to add on target remnants *
67C S PYDIFF to set up kinematics for diffractive events *
68C S PYDOCU to compute cross-sections and handle documentation *
69C S PYFRAM to perform boosts between different frames *
70C S PYWIDT to calculate full and partial widths of resonances *
71C S PYOFSH to calculate partial width into off-shell channels *
72C S PYRECO to handle colour reconnection in W+W- events *
73C S PYKLIM to calculate borders of allowed kinematical region *
74C S PYKMAP to construct value of kinematical variable *
75C S PYSIGH to calculate differential cross-sections *
76C S PYPDFU to evaluate parton distributions *
77C S PYPDFL to evaluate parton distributions at low x and Q^2 *
78C S PYPDEL to evaluate electron parton distributions *
79C S PYPDGA to evaluate photon parton distributions (generic) *
80C S PYGGAM to evaluate photon parton distributions (SaS sets) *
81C S PYGVMD to evaluate VMD part of photon parton distributions *
82C S PYGANO to evaluate anomalous part of photon pdf's *
83C S PYGBEH to evaluate Bethe-Heitler part of photon pdf's *
84C S PYGDIR to evaluate direct contribution to photon pdf's *
85C S PYPDPI to evaluate pion parton distributions *
86C S PYPDPR to evaluate proton parton distributions *
87C F PYCTEQ to evaluate the CTEQ 3 proton parton distributions *
88C S PYGRVL to evaluate the GRV 94L pronton parton distributions *
89C S PYGRVM to evaluate the GRV 94M pronton parton distributions *
90C S PYGRVD to evaluate the GRV 94D pronton parton distributions *
91C F PYGRVV auxiliary to the PYGRV* routines *
92C F PYGRVW auxiliary to the PYGRV* routines *
93C F PYGRVS auxiliary to the PYGRV* routines *
94C F PYHFTH to evaluate threshold factor for heavy flavour *
95C S PYSPLI to find flavours left in hadron when one removed *
96C F PYGAMM to evaluate ordinary Gamma function Gamma(x) *
97C S PYWAUX to evaluate auxiliary functions W1(s) and W2(s) *
98C S PYI3AU to evaluate auxiliary function I3(s,t,u,v) *
99C F PYSPEN to evaluate Spence (dilogarithm) function Sp(x) *
100C S PYQQBH to evaluate matrix element for g + g -> Q + Qbar + H *
101C *
102C S PYMSIN to initialize the supersymmetry simulation *
103C S PYAPPS to determine MSSM parameters from SUGRA input *
104C F PYRNMQ to determine running quark masses *
105C F PYRNMT to determine running top mass *
106C S PYTHRG to calculate sfermion third-gen. mass eigenstates *
107C S PYINOM to calculate neutralino/chargino mass eigenstates *
108C F PYRNM3 to determine running M3, gluino mass *
109C S PYEIG4 to calculate eigenvalues and -vectors in 4*4 matrix *
110C S PYHGGM to determine Higgs mass spectrum *
111C S PYSUBH to determine Higgs masses in the MSSM *
112C S PYPOLE to determine Higgs masses in the MSSM *
113C S PYVACU to determine Higgs masses in the MSSM *
114C S PYRGHM auxiliary to PYVACU *
115C S PYGFXX auxiliary to PYRGHM *
116C F PYFINT auxiliary to PYVACU *
117C F PYFISB auxiliary to PYFINT *
118C S PYSFDC to calculate sfermion decay partial widths *
119C S PYGLUI to calculate gluino decay partial widths *
120C S PYTBBN to calculate 3-body decay of gluino to neutralino *
121C S PYTBBC to calculate 3-body decay of gluino to chargino *
122C S PYNJDC to calculate neutralino decay partial widths *
123C S PYCJDC to calculate chargino decay partial widths *
124C F PYXXZ5 auxiliary for neutralino 3-body decay *
125C F PYXXW5 auxiliary for ino charge change 3-body decay *
126C F PYXXGA auxiliary for ino -> ino + gamma decay *
127C F PYX2XG auxiliary for ino -> ino + gauge boson decay *
128C F PYX2XH auxiliary for ino -> ino + Higgs decay *
129C F PYXXZ2 auxiliary for chargino 3-body decay *
130C S PYHEXT to calculate non-SM Higgs decay partial widths *
131C F PYH2XX auxiliary for H -> ino + ino decay *
132C F PYGAUS to perform Gaussian integration *
133C F PYSIMP to perform Simpson integration *
134C F PYLAMF to evaluate the lambda kinematics function *
135C S PYTBDY to perform 3-body decay of gauginos *
136C *
137C S PY1ENT to fill one entry (= parton or particle) *
138C S PY2ENT to fill two entries *
139C S PY3ENT to fill three entries *
140C S PY4ENT to fill four entries *
141C S PYJOIN to connect entries with colour flow information *
142C S PYGIVE to fill (or query) commonblock variables *
143C S PYEXEC to administrate fragmentation and decay chain *
144C S PYPREP to rearrange showered partons along strings *
145C S PYSTRF to do string fragmentation of jet system *
146C S PYINDF to do independent fragmentation of one or many jets *
147C S PYDECY to do the decay of a particle *
148C S PYDCYK to select parton and hadron flavours in decays *
149C S PYKFDI to select parton and hadron flavours in fragm *
150C S PYNMES to select number of popcorn mesons *
151C S PYKFIN to calculate falvour prod. ratios from input params. *
152C S PYPTDI to select transverse momenta in fragm *
153C S PYZDIS to select longitudinal scaling variable in fragm *
154C S PYSHOW to do timelike parton shower evolution *
155C S PYBOEI to include Bose-Einstein effects (crudely) *
156C F PYMASS to give the mass of a particle or parton *
157C S PYNAME to give the name of a particle or parton *
158C F PYCHGE to give three times the electric charge *
159C F PYCOMP to compress standard KF flavour code to internal KC *
160C S PYERRM to write error messages and abort faulty run *
161C F PYALEM to give the alpha_electromagnetic value *
162C F PYALPS to give the alpha_strong value *
163C F PYANGL to give the angle from known x and y components *
164C F PYR to provide a random number generator *
165C S PYRGET to save the state of the random number generator *
166C S PYRSET to set the state of the random number generator *
167C S PYROBO to rotate and/or boost an event *
168C S PYEDIT to remove unwanted entries from record *
169C S PYLIST to list event record or particle data *
170C S PYLOGO to write a logo *
171C S PYUPDA to update particle data *
172C F PYK to provide integer-valued event information *
173C F PYP to provide real-valued event information *
174C S PYSPHE to perform sphericity analysis *
175C S PYTHRU to perform thrust analysis *
176C S PYCLUS to perform three-dimensional cluster analysis *
177C S PYCELL to perform cluster analysis in (eta, phi, E_T) *
178C S PYJMAS to give high and low jet mass of event *
179C S PYFOWO to give Fox-Wolfram moments *
180C S PYTABU to analyze events, with tabular output *
181C *
182C S PYEEVT to administrate the generation of an e+e- event *
183C S PYXTEE to give the total cross-section at given CM energy *
184C S PYRADK to generate initial state photon radiation *
185C S PYXKFL to select flavour of primary qqbar pair *
186C S PYXJET to select (matrix element) jet multiplicity *
187C S PYX3JT to select kinematics of three-jet event *
188C S PYX4JT to select kinematics of four-jet event *
189C S PYXDIF to select angular orientation of event *
190C S PYONIA to perform generation of onium decay to gluons *
191C *
192C S PYBOOK to book a histogram *
193C S PYFILL to fill an entry in a histogram *
194C S PYFACT to multiply histogram contents by a factor *
195C S PYOPER to perform operations between histograms *
196C S PYHIST to print and reset all histograms *
197C S PYPLOT to print a single histogram *
198C S PYNULL to reset contents of a single histogram *
199C S PYDUMP to dump histogram contents onto a file *
200C *
201C S PYKCUT dummy routine for user kinematical cuts *
202C S PYEVWT dummy routine for weighting events *
203C S PYUPIN dummy routine to initialize a user process *
204C S PYUPEV dummy routine to generate a user process event *
205C S PDFSET dummy routine to be removed when using PDFLIB *
206C S STRUCTM dummy routine to be removed when using PDFLIB *
207C S PYTAUD dummy routine for interface to tau decay libraries *
208C S PYTIME dummy routine for giving date and time *
209C *
210C*********************************************************************
211
212*$ CREATE PYDATA.FOR
213*COPY PYDATA
214C...PYDATA
215C...Default values for switches and parameters,
216C...and particle, decay and process data.
217
218 BLOCK DATA PYDATA
219
220C...Double precision and integer declarations.
221 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
222 INTEGER PYK,PYCHGE,PYCOMP
223C...Commonblocks.
224 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
225 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
226 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
227 COMMON/PYDAT4/CHAF(500,2)
228 CHARACTER CHAF*16
229 COMMON/PYDATR/MRPY(6),RRPY(100)
230 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
231 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
232 COMMON/PYINT1/MINT(400),VINT(400)
233 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
234 COMMON/PYINT3/XSFX(2,-40:40),ISIG(1000,3),SIGH(1000)
235 COMMON/PYINT4/MWID(500),WIDS(500,5)
236 COMMON/PYINT5/NGENPD,NGEN(0:500,3),XSEC(0:500,3)
237 COMMON/PYINT6/PROC(0:500)
238 CHARACTER PROC*28
239 COMMON/PYINT7/SIGT(0:6,0:6,0:5)
240 COMMON/PYMSSM/IMSS(0:99),RMSS(0:99)
241 COMMON/PYSSMT/ZMIX(4,4),UMIX(2,2),VMIX(2,2),SMZ(4),SMW(2),
242 &SFMIX(16,4)
243 COMMON/PYBINS/IHIST(4),INDX(1000),BIN(20000)
244 SAVE /PYDAT1/,/PYDAT2/,/PYDAT3/,/PYDAT4/,/PYDATR/,/PYSUBS/,
245 &/PYPARS/,/PYINT1/,/PYINT2/,/PYINT3/,/PYINT4/,/PYINT5/,
246 &/PYINT6/,/PYINT7/,/PYMSSM/,/PYSSMT/,/PYBINS/
247
248C...PYDAT1, containing status codes and most parameters.
249 DATA MSTU/
250 & 0, 0, 0, 4000,10000, 500, 4000, 0, 0, 2,
251 1 6, 1, 1, 0, 1, 1, 0, 0, 0, 0,
252 2 2, 10, 0, 0, 1, 10, 0, 0, 0, 0,
253 3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
254 4 2, 2, 1, 4, 2, 1, 1, 0, 0, 0,
255 5 25, 24, 0, 1, 0, 0, 0, 0, 0, 0,
256 6 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
257 7 30*0,
258 1 1, 0, 0, 0, 0, 0, 0, 0, 0, 0,
259 2 1, 5, 3, 5, 0, 0, 0, 0, 0, 0,
260 & 80*0/
261 DATA PARU/
262 & 3.141592653589793D0, 6.283185307179586D0,
263 & 0.197327D0, 5.06773D0, 0.389380D0, 2.56819D0, 4*0D0,
264 1 0.001D0, 0.09D0, 0.01D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
265 2 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
266 3 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
267 4 2.0D0, 1.0D0, 0.25D0, 2.5D0, 0.05D0,
268 4 0D0, 0D0, 0.0001D0, 0D0, 0D0,
269 5 2.5D0,1.5D0,7.0D0,1.0D0,0.5D0,2.0D0,3.2D0, 0D0, 0D0, 0D0,
270 6 40*0D0,
271 & 0.00729735D0, 0.232D0, 0.007764D0, 1.0D0, 1.16639D-5,
272 & 0D0, 0D0, 0D0, 0D0, 0D0,
273 1 0.20D0, 0.25D0, 1.0D0, 4.0D0, 10D0, 0D0, 0D0, 0D0, 0D0, 0D0,
274 2 -0.693D0, -1.0D0, 0.387D0, 1.0D0, -0.08D0,
275 2 -1.0D0, 1.0D0, 1.0D0, 1.0D0, 0D0,
276 3 1.0D0,-1.0D0, 1.0D0,-1.0D0, 1.0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
277 4 5.0D0, 1.0D0, 1.0D0, 0D0, 1.0D0, 1.0D0, 0D0, 0D0, 0D0, 0D0,
278 5 1.0D0, 0D0, 0D0, 0D0, 1000D0, 1.0D0, 1.0D0, 1.0D0, 1.0D0,0D0,
279 6 1.0D0, 1.0D0, 1.0D0, 1.0D0, 1.0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
280 7 1.0D0, 1.0D0, 1.0D0, 1.0D0, 1.0D0, 1.0D0, 1.0D0, 0D0,0D0,0D0,
281 8 1.0D0, 1.0D0, 1.0D0, 0.0D0, 0.0D0, 1.0D0, 1.0D0, 0D0,0D0,0D0,
282 9 0D0, 0D0, 0D0, 0D0, 1.0D0, 0D0, 0D0, 0D0, 0D0, 0D0/
283 DATA MSTJ/
284 & 1, 3, 0, 0, 0, 0, 0, 0, 0, 0,
285 1 4, 2, 0, 1, 0, 0, 0, 0, 0, 0,
286 2 2, 1, 1, 2, 1, 2, 2, 0, 0, 0,
287 3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
288 4 2, 2, 4, 2, 5, 3, 3, 0, 0, 3,
289 5 0, 3, 0, 0, 0, 0, 0, 0, 0, 0,
290 6 40*0,
291 & 5, 2, 7, 5, 1, 1, 0, 2, 0, 2,
292 1 0, 0, 0, 0, 1, 1, 0, 0, 0, 0,
293 2 80*0/
294 DATA PARJ/
295 & 0.10D0, 0.30D0, 0.40D0, 0.05D0, 0.50D0,
296 & 0.50D0, 0.50D0, 0.6D0, 1.2D0, 0.6D0,
297 1 0.50D0,0.60D0,0.75D0, 0D0, 0D0, 0D0, 0D0, 1.0D0, 1.0D0, 0D0,
298 2 0.36D0, 1.0D0,0.01D0, 2.0D0,1.0D0,0.4D0, 0D0, 0D0, 0D0, 0D0,
299 3 0.10D0, 1.0D0, 0.8D0, 1.5D0,0D0,2.0D0,0.2D0,2.5D0,0.6D0,0D0,
300 4 0.3D0, 0.58D0, 0.5D0, 0.9D0,0.5D0,1.0D0,1.0D0,1.0D0,0D0,0D0,
301 5 0.77D0, 0.77D0, 0.77D0, -0.05D0, -0.005D0,
302 5 -0.00001D0, -0.00001D0, -0.00001D0, 1.0D0, 0D0,
303 6 4.5D0, 0.7D0, 0D0,0.003D0, 0.5D0, 0.5D0, 0D0, 0D0, 0D0, 0D0,
304 7 10D0, 1000D0, 100D0, 1000D0, 0D0, 0.7D0,10D0, 0D0, 0D0, 0D0,
305 8 0.29D0, 1.0D0, 1.0D0, 0D0, 10D0, 10D0, 0D0, 0D0, 0D0, 0D0,
306 9 0.02D0, 1.0D0, 0.2D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
307 & 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
308 1 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
309 2 1.0D0, 0.25D0,91.187D0,2.489D0, 0.01D0,
310 2 2.0D0, 1.0D0, 0.25D0,0.002D0, 0D0,
311 3 0D0, 0D0, 0D0, 0D0, 0.01D0, 0.99D0, 0D0, 0D0, 0.2D0, 0D0,
312 4 60*0D0/
313
314C...PYDAT2, with particle data and flavour treatment parameters.
315 DATA (KCHG(I,1),I= 1, 500)/-1,2,-1,2,-1,2,-1,2,2*0,-3,0,-3,0,
316 &-3,0,-3,6*0,3,9*0,3,2*0,3,0,-1,12*0,3,2*0,3,28*0,2,-1,20*0,4*3,
317 &8*0,3*3,4*0,3*3,3*0,3*3,7*0,3*3,3*0,3*3,3*0,-2,-3,2*1,3*0,4,3*3,
318 &6,2*-2,2*-3,0,2*1,2*0,2*3,-2,2*-3,2*0,-3,2*1,2*0,3,0,2*4,2*3,2*6,
319 &3,2*1,2*0,2*3,2*0,4,2*3,2*6,2*3,6,2*-2,2*-3,0,-3,0,2*1,2*0,2*3,0,
320 &3,2*-2,2*-3,2*0,2*-3,0,2*1,2*0,2*3,2*0,2*3,-2,2*-3,2*0,2*-3,2*0,
321 &-3,2*0,2*3,4*0,2*3,2*0,2*3,2*0,2*3,4*0,2*3,2*0,2*3,3*0,3,2*0,3,0,
322 &3,0,3,2*0,3,0,3,3*0,-1,2,-1,2,-1,2,-3,0,-3,0,-3,4*0,3,2*0,3,0,-1,
323 &2,-1,2,-1,2,-3,0,-3,0,-3,0,-1,2,-3,164*0/
324 DATA (KCHG(I,2),I= 1, 500)/8*1,12*0,2,16*0,2,1,113*0,-1,0,2*-1,
325 &3*0,-1,4*0,2*-1,3*0,2*-1,4*0,-1,5*0,2*-1,4*0,2*-1,5*0,2*-1,6*0,
326 &-1,7*0,2*-1,5*0,2*-1,6*0,2*-1,7*0,2*-1,8*0,-1,56*0,6*1,6*0,2,7*0,
327 &6*1,6*0,2*1,165*0/
328 DATA (KCHG(I,3),I= 1, 500)/8*1,2*0,8*1,5*0,1,9*0,1,2*0,1,0,2*1,
329 &11*0,1,2*0,1,26*0,1,0,2*1,20*0,4*1,5*0,6*1,4*0,9*1,4*0,12*1,3*0,
330 &102*1,2*0,2*1,2*0,4*1,2*0,6*1,2*0,8*1,3*0,1,0,2*1,0,3*1,0,4*1,
331 &3*0,12*1,3*0,1,2*0,1,0,16*1,163*0/
332 DATA (KCHG(I,4),I= 1, 293)/1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,
333 &16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,
334 &37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,
335 &58,59,60,61,62,63,64,65,66,67,68,69,70,71,72,73,74,75,76,77,78,
336 &79,80,81,82,83,84,85,86,87,88,89,90,91,92,93,94,95,96,97,98,99,
337 &100,110,111,113,115,130,210,211,213,215,220,221,223,225,310,311,
338 &313,315,321,323,325,330,331,333,335,411,413,415,421,423,425,431,
339 &433,435,440,441,443,445,511,513,515,521,523,525,531,533,535,541,
340 &543,545,551,553,555,1103,1114,2101,2103,2110,2112,2114,2203,2210,
341 &2212,2214,2224,3101,3103,3112,3114,3122,3201,3203,3212,3214,3222,
342 &3224,3303,3312,3314,3322,3324,3334,4101,4103,4112,4114,4122,4132,
343 &4201,4203,4212,4214,4222,4224,4232,4301,4303,4312,4314,4322,4324,
344 &4332,4334,4403,4412,4414,4422,4424,4432,4434,4444,5101,5103,5112,
345 &5114,5122,5132,5142,5201,5203,5212,5214,5222,5224,5232,5242,5301,
346 &5303,5312,5314,5322,5324,5332,5334,5342,5401,5403,5412,5414,5422,
347 &5424,5432,5434,5442,5444,5503,5512,5514,5522,5524,5532,5534,5542,
348 &5544,5554,10111,10113,10211,10213,10221,10223,10311,10313,10321,
349 &10323,10331,10333,10411,10413,10421,10423,10431,10433,10441,
350 &10443,10511,10513,10521,10523,10531,10533,10541,10543,10551,
351 &10553,20113,20213,20223,20313,20323,20333,20413,20423,20433/
352 DATA (KCHG(I,4),I= 294, 500)/20443,20513,20523,20533,20543,20553,
353 &100443,100553,1000001,1000002,1000003,1000004,1000005,1000006,
354 &1000011,1000012,1000013,1000014,1000015,1000016,1000021,1000022,
355 &1000023,1000024,1000025,1000035,1000037,1000039,2000001,2000002,
356 &2000003,2000004,2000005,2000006,2000011,2000012,2000013,2000014,
357 &2000015,2000016,4000001,4000002,4000011,4000012,163*0/
358 DATA (PMAS(I,1),I= 1, 214)/0.0099D0,0.0056D0,0.199D0,1.35D0,
359 &5D0,175D0,2*400D0,2*0D0,0.00051D0,0D0,0.10566D0,0D0,1.777D0,0D0,
360 &400D0,5*0D0,91.187D0,80.33D0,80D0,6*0D0,500D0,900D0,500D0,
361 &3*300D0,350D0,200D0,5000D0,10*0D0,3*100D0,3*200D0,26*0D0,1D0,2D0,
362 &5D0,16*0D0,0.13498D0,0.7685D0,1.318D0,0.49767D0,0D0,0.13957D0,
363 &0.7669D0,1.318D0,0D0,0.54745D0,0.78194D0,1.275D0,2*0.49767D0,
364 &0.8961D0,1.432D0,0.4936D0,0.8916D0,1.425D0,0D0,0.95777D0,
365 &1.0194D0,1.525D0,1.8693D0,2.01D0,2.46D0,1.8645D0,2.0067D0,2.46D0,
366 &1.9685D0,2.1124D0,2.5735D0,0D0,2.9798D0,3.09688D0,3.5562D0,
367 &5.2792D0,5.3248D0,5.83D0,5.2789D0,5.3248D0,5.83D0,5.3693D0,
368 &5.4163D0,6.07D0,6.594D0,6.602D0,7.35D0,9.4D0,9.4603D0,9.9132D0,
369 &0.77133D0,1.234D0,0.57933D0,0.77133D0,0D0,0.93957D0,1.233D0,
370 &0.77133D0,0D0,0.93827D0,1.232D0,1.231D0,0.80473D0,0.92953D0,
371 &1.19744D0,1.3872D0,1.11568D0,0.80473D0,0.92953D0,1.19255D0,
372 &1.3837D0,1.18937D0,1.3828D0,1.09361D0,1.3213D0,1.535D0,1.3149D0,
373 &1.5318D0,1.67245D0,1.96908D0,2.00808D0,2.4521D0,2.5D0,2.2849D0,
374 &2.4703D0,1.96908D0,2.00808D0,2.4535D0,2.5D0,2.4529D0,2.5D0,
375 &2.4656D0,2.15432D0,2.17967D0,2.55D0,2.63D0,2.55D0,2.63D0,2.704D0,
376 &2.8D0,3.27531D0,3.59798D0,3.65648D0,3.59798D0,3.65648D0,
377 &3.78663D0,3.82466D0,4.91594D0,5.38897D0,5.40145D0,5.8D0,5.81D0/
378 DATA (PMAS(I,1),I= 215, 500)/5.641D0,5.84D0,7.00575D0,5.38897D0,
379 &5.40145D0,5.8D0,5.81D0,5.8D0,5.81D0,5.84D0,7.00575D0,5.56725D0,
380 &5.57536D0,5.96D0,5.97D0,5.96D0,5.97D0,6.12D0,6.13D0,7.19099D0,
381 &6.67143D0,6.67397D0,7.03724D0,7.0485D0,7.03724D0,7.0485D0,
382 &7.21101D0,7.219D0,8.30945D0,8.31325D0,10.07354D0,10.42272D0,
383 &10.44144D0,10.42272D0,10.44144D0,10.60209D0,10.61426D0,
384 &11.70767D0,11.71147D0,15.11061D0,0.9835D0,1.231D0,0.9835D0,
385 &1.231D0,1D0,1.17D0,1.429D0,1.29D0,1.429D0,1.29D0,2*1.4D0,2.272D0,
386 &2.424D0,2.272D0,2.424D0,2.5D0,2.536D0,3.4151D0,3.46D0,5.68D0,
387 &5.73D0,5.68D0,5.73D0,5.92D0,5.97D0,7.25D0,7.3D0,9.8598D0,9.875D0,
388 &2*1.23D0,1.282D0,2*1.402D0,1.427D0,2*2.372D0,2.56D0,3.5106D0,
389 &2*5.78D0,6.02D0,7.3D0,9.8919D0,3.686D0,10.0233D0,32*500D0,
390 &4*400D0,163*0D0/
391 DATA (PMAS(I,2),I= 1, 500)/5*0D0,1.4D0,16*0D0,2.47833D0,
392 &2.069D0,0.00295D0,6*0D0,14.67788D0,0D0,16.79392D0,8.45231D0,
393 &4.93534D0,5.80468D0,19.1898D0,0.39162D0,417.35283D0,62*0D0,
394 &0.151D0,0.107D0,3*0D0,0.149D0,0.107D0,2*0D0,0.00843D0,0.185D0,
395 &2*0D0,0.0505D0,0.109D0,0D0,0.0498D0,0.098D0,0D0,0.0002D0,
396 &0.00443D0,0.076D0,2*0D0,0.023D0,2*0D0,0.023D0,2*0D0,0.015D0,0D0,
397 &0.0013D0,0D0,0.002D0,2*0D0,0.02D0,2*0D0,0.02D0,2*0D0,0.02D0,
398 &2*0D0,0.02D0,4*0D0,0.12D0,4*0D0,0.12D0,3*0D0,2*0.12D0,3*0D0,
399 &0.0394D0,4*0D0,0.036D0,0D0,0.0358D0,2*0D0,0.0099D0,0D0,0.0091D0,
400 &74*0D0,0.06D0,0.142D0,0.06D0,0.142D0,0D0,0.36D0,0.287D0,0.09D0,
401 &0.287D0,0.09D0,0.25D0,0.08D0,0.05D0,0.02D0,0.05D0,0.02D0,0.05D0,
402 &0D0,0.014D0,0.01D0,8*0.05D0,0D0,0.01D0,2*0.4D0,0.025D0,2*0.174D0,
403 &0.053D0,3*0.05D0,0.0009D0,4*0.05D0,3*0D0,19*1D0,0D0,7*1D0,0D0,
404 &1D0,0D0,1D0,0D0,2.60511D0,2.60839D0,0.42904D0,0.41921D0,163*0D0/
405 DATA (PMAS(I,3),I= 1, 500)/5*0D0,14D0,16*0D0,24.78326D0,
406 &20.69D0,0.02954D0,6*0D0,146.77876D0,0D0,167.93924D0,84.52308D0,
407 &49.35344D0,58.04675D0,191.89803D0,3.91624D0,4173.5283D0,62*0D0,
408 &0.4D0,0.25D0,3*0D0,0.4D0,0.25D0,2*0D0,0.1D0,0.17D0,2*0D0,0.2D0,
409 &0.12D0,0D0,0.2D0,0.12D0,0D0,0.002D0,0.015D0,0.2D0,2*0D0,0.12D0,
410 &2*0D0,0.12D0,2*0D0,0.05D0,0D0,0.005D0,0D0,0.01D0,2*0D0,0.05D0,
411 &2*0D0,0.05D0,2*0D0,0.05D0,2*0D0,0.05D0,4*0D0,0.14D0,4*0D0,0.14D0,
412 &3*0D0,2*0.14D0,3*0D0,0.04D0,4*0D0,0.035D0,0D0,0.035D0,2*0D0,
413 &0.05D0,0D0,0.05D0,74*0D0,0.05D0,0.25D0,0.05D0,0.25D0,0D0,0.2D0,
414 &0.4D0,0.005D0,0.4D0,0.01D0,0.35D0,0.001D0,0.1D0,0.08D0,0.1D0,
415 &0.08D0,0.1D0,0D0,0.05D0,0.02D0,6*0.1D0,0.05D0,0.1D0,0D0,0.02D0,
416 &2*0.3D0,0.05D0,2*0.3D0,0.02D0,2*0.1D0,0.03D0,0.001D0,4*0.1D0,
417 &3*0D0,19*10D0,0.00001D0,7*10D0,0.00001D0,10D0,0.00001D0,10D0,
418 &0.00001D0,26.05109D0,26.08388D0,4.29043D0,4.19206D0,163*0D0/
419 DATA (PMAS(I,4),I= 1, 500)/12*0D0,658654D0,0D0,0.0872D0,68*0D0,
420 &0.1D0,0.387D0,16*0D0,0.00003D0,2*0D0,15500D0,0D0,7804.5D0,6*0D0,
421 &26.762D0,3*0D0,3709D0,6*0D0,0.317D0,2*0D0,0.1244D0,2*0D0,0.14D0,
422 &6*0D0,0.468D0,2*0D0,0.462D0,2*0D0,0.483D0,2*0D0,0.15D0,19*0D0,
423 &44.34D0,0D0,78.88D0,4*0D0,23.96D0,2*0D0,49.1D0,0D0,87.1D0,0D0,
424 &24.6D0,4*0D0,0.0618D0,0.029D0,6*0D0,0.106D0,6*0D0,0.019D0,2*0D0,
425 &7*0.1D0,4*0D0,0.342D0,2*0.387D0,6*0D0,2*0.387D0,6*0D0,0.387D0,
426 &0D0,0.387D0,2*0D0,8*0.387D0,0D0,9*0.387D0,83*0D0,163*0D0/
427 DATA PARF/
428 & 0.5D0,0.25D0, 0.5D0,0.25D0, 1D0, 0.5D0, 0D0, 0D0, 0D0, 0D0,
429 1 0.5D0, 0D0, 0.5D0, 0D0, 1D0, 1D0, 0D0, 0D0, 0D0, 0D0,
430 2 0.5D0, 0D0, 0.5D0, 0D0, 1D0, 1D0, 0D0, 0D0, 0D0, 0D0,
431 3 0.5D0, 0D0, 0.5D0, 0D0, 1D0, 1D0, 0D0, 0D0, 0D0, 0D0,
432 4 0.5D0, 0D0, 0.5D0, 0D0, 1D0, 1D0, 0D0, 0D0, 0D0, 0D0,
433 5 0.5D0, 0D0, 0.5D0, 0D0, 1D0, 1D0, 0D0, 0D0, 0D0, 0D0,
434 6 0.75D0, 0.5D0, 0D0,0.1667D0,0.0833D0,0.1667D0,0D0,0D0,0D0, 0D0,
435 7 0D0, 0D0, 1D0,0.3333D0,0.6667D0,0.3333D0,0D0,0D0,0D0, 0D0,
436 8 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
437 9 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
438 & 0.325D0,0.325D0,0.5D0,1.6D0, 5.0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
439 1 0D0,0.11D0,0.16D0,0.048D0,0.50D0,0.45D0,0.55D0,0.60D0,0D0,0D0,
440 2 0.2D0, 0.1D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0, 0D0,
441 3 60*0D0,
442 4 0.2D0, 0.5D0, 8*0D0,
443 5 1800*0D0/
444 DATA ((VCKM(I,J),J=1,4),I=1,4)/
445 & 0.95113D0, 0.04884D0, 0.00003D0, 0.00000D0,
446 & 0.04884D0, 0.94940D0, 0.00176D0, 0.00000D0,
447 & 0.00003D0, 0.00176D0, 0.99821D0, 0.00000D0,
448 & 0.00000D0, 0.00000D0, 0.00000D0, 1.00000D0/
449
450C...PYDAT3, with particle decay parameters and data.
451 DATA (MDCY(I,1),I= 1, 500)/5*0,3*1,6*0,1,0,1,5*0,3*1,6*0,1,0,
452 &7*1,10*0,2*1,0,3*1,26*0,3*1,16*0,3*1,3*0,2*1,0,7*1,0,2*1,0,12*1,
453 &0,18*1,0,1,4*0,1,3*0,2*1,2*0,3*1,2*0,4*1,0,5*1,2*0,4*1,2*0,5*1,
454 &2*0,6*1,0,7*1,2*0,5*1,2*0,6*1,2*0,7*1,2*0,8*1,0,75*1,0,7*1,0,1,0,
455 &1,0,4*1,163*0/
456 DATA (MDCY(I,2),I= 1, 500)/1,9,17,25,33,41,54,64,2*0,74,78,80,
457 &85,87,141,143,148,2*0,151,160,172,188,208,6*0,287,0,309,332,414,
458 &494,521,524,525,10*0,534,539,0,544,564,588,26*0,606,607,611,16*0,
459 &620,622,627,636,0,645,647,649,0,656,664,670,679,681,683,686,696,
460 &702,705,0,716,722,733,739,802,805,813,874,876,884,917,919,0,923,
461 &924,927,929,965,966,974,1010,1011,1019,1058,1059,1063,1094,1095,
462 &1099,1100,1109,0,1111,4*0,1112,3*0,1115,1118,2*0,1119,1121,1124,
463 &2*0,1128,1129,1132,1135,0,1138,1143,1145,1148,1150,2*0,1154,1155,
464 &1156,1232,2*0,1236,1237,1238,1239,1240,2*0,1244,1245,1247,1248,
465 &1250,1254,0,1255,1259,1263,1267,1271,1275,1279,2*0,1283,1284,
466 &1285,1302,1311,2*0,1320,1321,1322,1323,1324,1333,2*0,1342,1343,
467 &1344,1345,1346,1355,1356,2*0,1365,1374,1383,1392,1401,1410,1419,
468 &1428,0,1437,1446,1455,1464,1473,1482,1491,1500,1509,1518,1519,
469 &1520,1521,1522,1527,1530,1532,1537,1539,1544,1551,1555,1557,1559,
470 &1561,1563,1565,1567,1569,1570,1572,1574,1576,1578,1580,1582,1584,
471 &1586,1588,1589,1591,1593,1607,1609,1611,1615,1617,1619,1621,1623,
472 &1625,1627,1629,1631,1633,1644,1658,1670,1682,1694,1706,1718,1731,
473 &1742,1753,1764,1775,1786,1797,1858,1863,1965,2021,2139,2273,0,
474 &2344,2360,2376,2392,2408,2424,2440,0,2455,0,2470,0,2485,2489,
475 &2493,2496,163*0/
476 DATA (MDCY(I,3),I= 1, 500)/5*8,13,2*10,2*0,4,2,5,2,54,2,5,3,
477 &2*0,9,12,16,20,79,6*0,22,0,23,82,80,27,3,1,9,10*0,2*5,0,20,24,18,
478 &26*0,1,4,9,16*0,2,5,2*9,0,2*2,7,0,8,6,9,2*2,3,10,6,3,11,0,6,11,6,
479 &63,3,8,61,2,8,33,2,4,0,1,3,2,36,1,8,36,1,8,39,1,4,31,1,4,1,9,2,0,
480 &1,4*0,3,3*0,3,1,2*0,2,3,4,2*0,1,3*3,0,5,2,3,2,4,2*0,2*1,76,4,2*0,
481 &4*1,4,2*0,1,2,1,2,4,1,0,7*4,2*0,2*1,17,2*9,2*0,4*1,2*9,2*0,4*1,9,
482 &1,9,2*0,8*9,0,9*9,4*1,5,3,2,5,2,5,7,4,7*2,1,9*2,1,2*2,14,2*2,4,
483 &9*2,11,14,5*12,13,6*11,61,5,102,56,118,134,71,0,6*16,15,0,15,0,
484 &15,0,2*4,3,2,163*0/
485 DATA (MDME(I,1),I= 1,4000)/6*1,-1,7*1,-1,7*1,-1,7*1,-1,7*1,-1,
486 &7*1,-1,1,-1,12*1,2*-1,8*1,2*-1,73*1,-1,2*1,-1,6*1,2*-1,7*1,2*-1,
487 &3*1,-1,6*1,2*-1,6*1,2*-1,3*1,-1,3*1,-1,3*1,5*-1,3*1,-1,85*1,2*-1,
488 &6*1,8*-1,3*1,-1,3*1,-1,3*1,5*-1,3*1,4*-1,197*1,2*-1,2*1,-1,20*1,
489 &2*-1,6*1,2*-1,7*1,-1,3*1,-1,3*1,5*-1,3*1,-1,1,-1,6*1,2*-1,6*1,
490 &2*-1,1892*1,1503*0/
491 DATA (MDME(I,2),I= 1,4000)/43*102,4*0,102,0,4*53,3*102,4*0,102,
492 &2*0,3*102,4*0,102,2*0,6*102,42,6*102,2*42,2*0,8*41,2*0,36*41,
493 &8*102,0,102,0,102,2*0,21*102,8*32,8*0,16*32,21*0,62*53,8*32,14*0,
494 &16*32,27*0,62*53,18*0,62*53,9*0,18*53,3*32,0,6*32,3*0,2*32,3*0,
495 &2*32,7*0,8*32,12*0,16*32,6*0,8*32,8*0,12,2*42,2*11,9*42,0,2,3,
496 &15*0,4*42,5*0,3,12*0,2,3*0,1,0,3,16*0,2*3,15*0,2*42,2*3,18*0,2*3,
497 &3*0,1,11*0,22*42,41*0,2*3,9*0,16*42,45*0,3,10*0,10*42,20*0,2*13,
498 &6*0,12,2*0,12,0,12,14*42,16*0,48,3*13,2*42,9*0,14*42,16*0,48,
499 &3*13,2*42,9*0,14*42,19*0,48,3*13,2*42,6*0,2*11,28*42,5*0,32,3*0,
500 &4*32,2*4,0,32,45*0,14*42,52*0,10*13,2*42,2*11,4*0,2*42,2*11,6*0,
501 &2*42,2*11,0,2*42,2*11,2*42,2*11,2*42,2*11,2*42,2*11,2*42,2*11,
502 &2*42,2*11,2*42,2*11,2*0,3*42,8*0,48,3*13,20*42,4*0,18*42,4*0,
503 &9*42,0,162*42,50*0,2*12,17*0,2*32,33*0,12,9*0,32,2*0,12,11*0,
504 &4*32,2*4,5*0,828*53,1515*0/
505 DATA (BRAT(I) ,I= 1, 418)/43*0D0,0.00003D0,0.00177D0,0.9982D0,
506 &33*0D0,1D0,6*0D0,0.1783D0,0.1735D0,0.1131D0,0.2494D0,0.003D0,
507 &0.09D0,0.0027D0,0.01D0,0.0014D0,0.0012D0,2*0.00025D0,0.0071D0,
508 &0.012D0,0.0004D0,0.00075D0,0.00006D0,2*0.00078D0,0.0034D0,0.08D0,
509 &0.011D0,0.0191D0,0.00006D0,0.005D0,0.0133D0,0.0067D0,0.0005D0,
510 &0.0035D0,0.0006D0,0.0015D0,0.00021D0,0.0002D0,0.00075D0,0.0001D0,
511 &0.0002D0,0.0011D0,3*0.0002D0,0.00022D0,0.0004D0,0.0001D0,
512 &2*0.00205D0,2*0.00069D0,0.00025D0,0.00051D0,0.00025D0,35*0D0,
513 &0.15403D0,0.11945D0,0.15402D0,0.11931D0,0.15215D0,3*0D0,
514 &0.03357D0,0.0668D0,0.03357D0,0.0668D0,0.0335D0,0.0668D0,2*0D0,
515 &0.32139D0,0.0165D0,2*0D0,0.0165D0,0.32067D0,2*0D0,0.00001D0,
516 &0.00059D0,6*0D0,2*0.10814D0,0.10806D0,3*0D0,0.00031D0,0.04438D0,
517 &0.88031D0,4*0D0,0.0002D0,0.05531D0,0D0,0.01838D0,0.00071D0,0D0,
518 &0.00009D0,0.00032D0,62*0D0,0.14449D0,0.11223D0,0.14449D0,
519 &0.11223D0,0.14443D0,0.05782D0,2*0D0,0.03172D0,0.06305D0,
520 &0.03172D0,0.06305D0,0.03172D0,0.06305D0,8*0D0,0.24928D0,0.0128D0,
521 &0.00001D0,0D0,0.0128D0,0.24882D0,0.00039D0,0D0,0.00001D0,
522 &0.00046D0,0.22153D0,5*0D0,2*0.08464D0,0.08463D0,7*0D0,0.00005D0,
523 &0.00097D0,5*0D0,0.00007D0,0D0,0.00049D0,0.00001D0,0.00006D0,
524 &0.30591D0,0.68863D0,0D0,0.0038D0,66*0D0,0.00008D0,0.00167D0/
525 DATA (BRAT(I) ,I= 419, 722)/5*0D0,0.00013D0,0D0,0.00294D0,
526 &0.00001D0,3*0D0,0.99517D0,63*0D0,0.00002D0,0.07231D0,2*0D0,
527 &0.00001D0,0.00269D0,0D0,0.92497D0,18*0D0,0.0024D0,0.99483D0,
528 &0.00278D0,1D0,3*0.21511D0,0.21478D0,2*0D0,2*0.06995D0,2*0D0,1D0,
529 &3*0D0,0.95D0,0.05D0,3*0D0,4*0.25D0,16*0D0,4*0.25D0,20*0D0,1D0,
530 &17*0D0,1D0,2*0.08D0,0.76D0,0.08D0,2*0.105D0,0.04D0,0.5D0,0.08D0,
531 &0.14D0,0.01D0,0.015D0,0.005D0,0.988D0,0.012D0,0.998739D0,
532 &0.00079D0,0.00038D0,0.000046D0,0.000045D0,2*0.34725D0,0.144D0,
533 &0.104D0,0.0245D0,2*0.01225D0,0.0028D0,0.0057D0,0.2112D0,0.1256D0,
534 &2*0.1939D0,2*0.1359D0,0.002D0,0.001D0,0.0006D0,0.999877D0,
535 &0.000123D0,0.99955D0,0.00045D0,2*0.34725D0,0.144D0,0.104D0,
536 &0.049D0,0.0028D0,0.0057D0,0.3923D0,0.321D0,0.2317D0,0.0478D0,
537 &0.0049D0,0.0013D0,0.0003D0,0.0007D0,0.89D0,0.08693D0,0.0221D0,
538 &0.00083D0,2*0.00007D0,0.564D0,0.282D0,0.072D0,0.028D0,0.023D0,
539 &2*0.0115D0,0.005D0,0.003D0,0.6861D0,0.3139D0,2*0.5D0,0.665D0,
540 &0.333D0,0.002D0,0.333D0,0.166D0,0.168D0,0.084D0,0.087D0,0.043D0,
541 &0.059D0,2*0.029D0,0.002D0,0.6352D0,0.2116D0,0.0559D0,0.0173D0,
542 &0.0482D0,0.0318D0,0.666D0,0.333D0,0.001D0,0.332D0,0.166D0,
543 &0.168D0,0.084D0,0.086D0,0.043D0,0.059D0,2*0.029D0,2*0.002D0,
544 &0.437D0,0.208D0,0.302D0,0.0302D0,0.0212D0,0.0016D0,0.48947D0/
545 DATA (BRAT(I) ,I= 723, 897)/0.34D0,3*0.043D0,0.027D0,0.0126D0,
546 &0.0013D0,0.0003D0,0.00025D0,0.00008D0,0.444D0,2*0.222D0,0.104D0,
547 &2*0.004D0,0.07D0,0.065D0,2*0.005D0,2*0.011D0,5*0.001D0,0.07D0,
548 &0.065D0,2*0.005D0,2*0.011D0,5*0.001D0,0.026D0,0.019D0,0.066D0,
549 &0.041D0,0.045D0,0.076D0,0.0073D0,2*0.0047D0,0.026D0,0.001D0,
550 &0.0006D0,0.0066D0,0.005D0,2*0.003D0,2*0.0006D0,2*0.001D0,0.006D0,
551 &0.005D0,0.012D0,0.0057D0,0.067D0,0.008D0,0.0022D0,0.027D0,
552 &0.004D0,0.019D0,0.012D0,0.002D0,0.009D0,0.0218D0,0.001D0,0.022D0,
553 &0.087D0,0.001D0,0.0019D0,0.0015D0,0.0028D0,0.683D0,0.306D0,
554 &0.011D0,0.3D0,0.15D0,0.16D0,0.08D0,0.13D0,0.06D0,0.08D0,0.04D0,
555 &0.034D0,0.027D0,2*0.002D0,2*0.004D0,2*0.002D0,0.034D0,0.027D0,
556 &2*0.002D0,2*0.004D0,2*0.002D0,0.0365D0,0.045D0,0.073D0,0.062D0,
557 &3*0.021D0,0.0061D0,0.015D0,0.025D0,0.0088D0,0.074D0,0.0109D0,
558 &0.0041D0,0.002D0,0.0035D0,0.0011D0,0.001D0,0.0027D0,2*0.0016D0,
559 &0.0018D0,0.011D0,0.0063D0,0.0052D0,0.018D0,0.016D0,0.0034D0,
560 &0.0036D0,0.0009D0,0.0006D0,0.015D0,0.0923D0,0.018D0,0.022D0,
561 &0.0077D0,0.009D0,0.0075D0,0.024D0,0.0085D0,0.067D0,0.0511D0,
562 &0.017D0,0.0004D0,0.0028D0,0.619D0,0.381D0,0.3D0,0.15D0,0.16D0,
563 &0.08D0,0.13D0,0.06D0,0.08D0,0.04D0,0.01D0,2*0.02D0,0.03D0,
564 &2*0.005D0,2*0.02D0,0.03D0,2*0.005D0,0.015D0,0.037D0,0.028D0/
565 DATA (BRAT(I) ,I= 898,1063)/0.079D0,0.095D0,0.052D0,0.0078D0,
566 &4*0.001D0,0.028D0,0.033D0,0.026D0,0.05D0,0.01D0,4*0.005D0,0.25D0,
567 &0.0952D0,0.94D0,0.06D0,2*0.4D0,2*0.1D0,1D0,0.0602D0,0.0601D0,
568 &0.8797D0,0.135D0,0.865D0,0.02D0,0.055D0,2*0.005D0,0.008D0,
569 &0.012D0,0.02D0,0.055D0,2*0.005D0,0.008D0,0.012D0,0.01D0,0.03D0,
570 &0.0035D0,0.011D0,0.0055D0,0.0042D0,0.009D0,0.018D0,0.015D0,
571 &0.0185D0,0.0135D0,0.025D0,0.0004D0,0.0007D0,0.0008D0,0.0014D0,
572 &0.0019D0,0.0025D0,0.4291D0,0.08D0,0.07D0,0.02D0,0.015D0,0.005D0,
573 &1D0,0.3D0,0.15D0,0.16D0,0.08D0,0.13D0,0.06D0,0.08D0,0.04D0,
574 &0.02D0,0.055D0,2*0.005D0,0.008D0,0.012D0,0.02D0,0.055D0,
575 &2*0.005D0,0.008D0,0.012D0,0.01D0,0.03D0,0.0035D0,0.011D0,
576 &0.0055D0,0.0042D0,0.009D0,0.018D0,0.015D0,0.0185D0,0.0135D0,
577 &0.025D0,0.0004D0,0.0007D0,0.0008D0,0.0014D0,0.0019D0,0.0025D0,
578 &0.4291D0,0.08D0,0.07D0,0.02D0,0.015D0,0.005D0,1D0,0.3D0,0.15D0,
579 &0.16D0,0.08D0,0.13D0,0.06D0,0.08D0,0.04D0,0.02D0,0.055D0,
580 &2*0.005D0,0.008D0,0.012D0,0.02D0,0.055D0,2*0.005D0,0.008D0,
581 &0.012D0,0.01D0,0.03D0,0.0035D0,0.011D0,0.0055D0,0.0042D0,0.009D0,
582 &0.018D0,0.015D0,0.0185D0,0.0135D0,0.025D0,2*0.0002D0,0.0007D0,
583 &2*0.0004D0,0.0014D0,0.001D0,0.0009D0,0.0025D0,0.4291D0,0.08D0,
584 &0.07D0,0.02D0,0.015D0,0.005D0,1D0,2*0.3D0,2*0.2D0,0.047D0/
585 DATA (BRAT(I) ,I=1064,1254)/0.122D0,0.006D0,0.012D0,0.035D0,
586 &0.012D0,0.035D0,0.003D0,0.007D0,0.15D0,0.037D0,0.008D0,0.002D0,
587 &0.05D0,0.015D0,0.003D0,0.001D0,0.014D0,0.042D0,0.014D0,0.042D0,
588 &0.24D0,0.065D0,0.012D0,0.003D0,0.001D0,0.002D0,0.001D0,0.002D0,
589 &0.014D0,0.003D0,1D0,2*0.3D0,2*0.2D0,1D0,0.0252D0,0.0248D0,
590 &0.0267D0,0.015D0,0.045D0,0.015D0,0.045D0,0.7743D0,0.029D0,0.22D0,
591 &0.78D0,1D0,0.331D0,0.663D0,0.006D0,0.663D0,0.331D0,0.006D0,1D0,
592 &0.999D0,0.001D0,0.88D0,2*0.06D0,0.639D0,0.358D0,0.002D0,0.001D0,
593 &1D0,0.88D0,2*0.06D0,0.516D0,0.483D0,0.001D0,0.88D0,2*0.06D0,
594 &0.9988D0,0.0001D0,0.0006D0,0.0004D0,0.0001D0,0.667D0,0.333D0,
595 &0.9954D0,0.0011D0,0.0035D0,0.333D0,0.667D0,0.676D0,0.234D0,
596 &0.085D0,0.005D0,2*1D0,0.018D0,2*0.005D0,0.003D0,0.002D0,
597 &2*0.006D0,0.018D0,2*0.005D0,0.003D0,0.002D0,2*0.006D0,0.0066D0,
598 &0.025D0,0.016D0,0.0088D0,2*0.005D0,0.0058D0,0.005D0,0.0055D0,
599 &4*0.004D0,2*0.002D0,2*0.004D0,0.003D0,0.002D0,2*0.003D0,
600 &3*0.002D0,2*0.001D0,0.002D0,2*0.001D0,2*0.002D0,0.0013D0,
601 &0.0018D0,5*0.001D0,4*0.003D0,2*0.005D0,2*0.002D0,2*0.001D0,
602 &2*0.002D0,2*0.001D0,0.2432D0,0.057D0,2*0.035D0,0.15D0,2*0.075D0,
603 &0.03D0,2*0.015D0,2*0.08D0,0.76D0,0.08D0,4*1D0,2*0.08D0,0.76D0,
604 &0.08D0,1D0,2*0.5D0,1D0,2*0.5D0,2*0.08D0,0.76D0,0.08D0,1D0/
605 DATA (BRAT(I) ,I=1255,1447)/2*0.08D0,0.76D0,3*0.08D0,0.76D0,
606 &3*0.08D0,0.76D0,3*0.08D0,0.76D0,3*0.08D0,0.76D0,3*0.08D0,0.76D0,
607 &3*0.08D0,0.76D0,0.08D0,2*1D0,2*0.105D0,0.04D0,0.0077D0,0.02D0,
608 &0.0235D0,0.0285D0,0.0435D0,0.0011D0,0.0022D0,0.0044D0,0.4291D0,
609 &0.08D0,0.07D0,0.02D0,0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,
610 &0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,
611 &0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,4*1D0,2*0.105D0,0.04D0,
612 &0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,0.04D0,
613 &0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,4*1D0,2*0.105D0,
614 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,1D0,2*0.105D0,
615 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
616 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
617 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
618 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
619 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
620 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
621 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
622 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
623 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0,
624 &0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,0.015D0,0.005D0,2*0.105D0/
625 DATA (BRAT(I) ,I=1448,1648)/0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
626 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
627 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
628 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
629 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
630 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
631 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
632 &0.015D0,0.005D0,2*0.105D0,0.04D0,0.5D0,0.08D0,0.14D0,0.01D0,
633 &0.015D0,0.005D0,4*1D0,0.52D0,0.26D0,0.11D0,2*0.055D0,0.333D0,
634 &0.334D0,0.333D0,0.667D0,0.333D0,0.28D0,0.14D0,0.313D0,0.157D0,
635 &0.11D0,0.667D0,0.333D0,0.28D0,0.14D0,0.313D0,0.157D0,0.11D0,
636 &0.36D0,0.18D0,0.03D0,2*0.015D0,2*0.2D0,4*0.25D0,0.667D0,0.333D0,
637 &0.667D0,0.333D0,0.667D0,0.333D0,0.667D0,0.333D0,4*0.5D0,0.007D0,
638 &0.993D0,1D0,0.667D0,0.333D0,0.667D0,0.333D0,0.667D0,0.333D0,
639 &0.667D0,0.333D0,8*0.5D0,0.02D0,0.98D0,1D0,4*0.5D0,3*0.146D0,
640 &3*0.05D0,0.15D0,2*0.05D0,4*0.024D0,0.066D0,0.667D0,0.333D0,
641 &0.667D0,0.333D0,4*0.25D0,0.667D0,0.333D0,0.667D0,0.333D0,2*0.5D0,
642 &0.273D0,0.727D0,0.667D0,0.333D0,0.667D0,0.333D0,4*0.5D0,0.35D0,
643 &0.65D0,2*0.0083D0,0.1866D0,0.324D0,0.184D0,0.027D0,0.001D0,
644 &0.093D0,0.087D0,0.078D0,0.0028D0,3*0.014D0,0.008D0,0.024D0/
645 DATA (BRAT(I) ,I=1649,4000)/0.008D0,0.024D0,0.425D0,0.02D0,
646 &0.185D0,0.088D0,0.043D0,0.067D0,0.066D0,827*0D0,0.8516D0,
647 &0.00539D0,0.04483D0,0.09819D0,0.85053D0,0.02152D0,0.02989D0,
648 &0.09806D0,0.29439D0,0.10943D0,0.59618D0,0.38983D0,0.61017D0,
649 &1503*0D0/
650 DATA (KFDP(I,1),I= 1, 375)/21,22,23,4*-24,25,21,22,23,4*24,25,
651 &21,22,23,4*-24,25,21,22,23,4*24,25,21,22,23,4*-24,25,21,22,23,
652 &4*24,25,37,1000022,1000023,1000025,1000035,21,22,23,4*-24,25,
653 &2*-37,21,22,23,4*24,25,2*37,22,23,-24,25,23,24,-12,22,23,-24,25,
654 &23,24,-12,-14,48*16,22,23,-24,25,23,24,22,23,-24,25,-37,23,24,37,
655 &1,2,3,4,5,6,7,8,21,1,2,3,4,5,6,7,8,11,13,15,17,1,2,3,4,5,6,7,8,
656 &11,12,13,14,15,16,17,18,4*-1,4*-3,4*-5,4*-7,-11,-13,-15,-17,1,2,
657 &3,4,5,6,7,8,11,13,15,17,21,2*22,23,24,1000022,2*1000023,
658 &3*1000025,4*1000035,2*1000024,2*1000037,1000001,2000001,1000001,
659 &-1000001,1000002,2000002,1000002,-1000002,1000003,2000003,
660 &1000003,-1000003,1000004,2000004,1000004,-1000004,1000005,
661 &2000005,1000005,-1000005,1000006,2000006,1000006,-1000006,
662 &1000011,2000011,1000011,-1000011,1000012,2000012,1000012,
663 &-1000012,1000013,2000013,1000013,-1000013,1000014,2000014,
664 &1000014,-1000014,1000015,2000015,1000015,-1000015,1000016,
665 &2000016,1000016,-1000016,1,2,3,4,5,6,7,8,11,12,13,14,15,16,17,18,
666 &24,37,2*23,25,35,4*-1,4*-3,4*-5,4*-7,-11,-13,-15,-17,3*24,1,2,3,
667 &4,5,6,7,8,11,13,15,17,21,2*22,23,24,23,25,36,1000022,2*1000023,
668 &3*1000025,4*1000035,2*1000024,2*1000037,1000001,2000001,1000001,
669 &-1000001,1000002,2000002,1000002,-1000002,1000003,2000003/
670 DATA (KFDP(I,1),I= 376, 606)/1000003,-1000003,1000004,2000004,
671 &1000004,-1000004,1000005,2000005,1000005,-1000005,1000006,
672 &2000006,1000006,-1000006,1000011,2000011,1000011,-1000011,
673 &1000012,2000012,1000012,-1000012,1000013,2000013,1000013,
674 &-1000013,1000014,2000014,1000014,-1000014,1000015,2000015,
675 &1000015,-1000015,1000016,2000016,1000016,-1000016,1,2,3,4,5,6,7,
676 &8,11,13,15,17,21,2*22,23,24,23,1000022,2*1000023,3*1000025,
677 &4*1000035,2*1000024,2*1000037,1000001,2000001,1000001,-1000001,
678 &1000002,2000002,1000002,-1000002,1000003,2000003,1000003,
679 &-1000003,1000004,2000004,1000004,-1000004,1000005,2000005,
680 &1000005,-1000005,1000006,2000006,1000006,-1000006,1000011,
681 &2000011,1000011,-1000011,1000012,2000012,1000012,-1000012,
682 &1000013,2000013,1000013,-1000013,1000014,2000014,1000014,
683 &-1000014,1000015,2000015,1000015,-1000015,1000016,2000016,
684 &1000016,-1000016,-1,-3,-5,-7,-11,-13,-15,-17,24,2*1000022,
685 &2*1000023,2*1000025,2*1000035,1000006,2000006,1000006,2000006,
686 &-1000001,-1000003,-1000011,-1000013,-1000015,-2000015,5,6,21,2,1,
687 &2,3,4,5,6,11,13,15,4,5,11,13,15,2*4,-11,-13,-15,2*24,2*52,1,2,3,
688 &4,5,6,7,8,11,12,13,14,15,16,17,18,2*24,2*52,4*-1,4*-3,4*-5,4*-7,
689 &-11,-13,-15,-17,22,23,1,2,3,4,5,6,7,8,11,12,13,14,15,16,17,18,82/
690 DATA (KFDP(I,1),I= 607,1001)/-11,-13,2*2,-12,-14,-16,2*-2,2*-4,
691 &-2,-4,2*22,211,111,221,13,11,213,-213,221,223,321,130,310,111,
692 &331,111,211,-12,12,-14,14,211,111,22,-13,-11,2*211,213,113,221,
693 &223,321,211,331,22,111,211,2*22,211,22,111,211,22,211,221,111,11,
694 &211,111,2*211,321,130,310,221,111,211,111,130,310,321,2*311,321,
695 &311,323,313,323,313,321,3*311,-13,3*211,12,14,311,2*321,311,321,
696 &313,323,313,323,311,4*321,211,111,3*22,111,321,130,-213,113,213,
697 &211,22,111,11,13,211,321,130,310,221,211,111,11*-11,11*-13,-311,
698 &-313,-311,-313,-20313,2*-311,-313,-311,-313,2*111,2*221,2*331,
699 &2*113,2*223,2*333,-311,-313,2*-321,211,-311,-321,333,-311,-313,
700 &-321,211,2*-321,2*-311,-321,211,113,421,2*411,421,411,423,413,
701 &423,413,421,411,8*-11,8*-13,-321,-323,-321,-323,-311,2*-313,-311,
702 &-313,2*-311,-321,-10323,-321,-323,-321,-311,2*-313,211,111,333,
703 &3*-321,-311,-313,-321,-313,310,333,211,2*-321,-311,-313,-311,211,
704 &-321,3*-311,211,113,321,2*421,411,421,413,423,413,423,411,421,
705 &-15,5*-11,5*-13,221,331,333,221,331,333,10221,211,213,211,213,
706 &321,323,321,323,2212,221,331,333,221,2*2,2*431,421,411,423,413,
707 &82,11,13,82,443,82,6*12,6*14,2*16,3*-411,3*-413,2*-411,2*-413,
708 &2*441,2*443,2*20443,2*2,2*4,2,4,511,521,511,523,513,523,513,521,
709 &511,6*12,6*14,2*16,3*-421,3*-423,2*-421,2*-423,2*441,2*443/
710 DATA (KFDP(I,1),I=1002,1428)/2*20443,2*2,2*4,2,4,521,511,521,513,
711 &523,513,523,511,521,6*12,6*14,2*16,3*-431,3*-433,2*-431,2*-433,
712 &3*441,3*443,3*20443,2*2,2*4,2,4,531,521,511,523,513,16,2*4,2*12,
713 &2*14,2*16,4*2,4*4,2*-11,2*-13,2*-1,2*-3,2*-11,2*-13,2*-1,541,511,
714 &521,513,523,21,11,13,15,1,2,3,4,21,22,553,21,2112,2212,2*2112,
715 &2212,2112,2*2212,2112,-12,3122,3212,3112,2212,2*2112,-12,2*3122,
716 &3222,3112,2212,2112,2212,3122,3222,3212,3122,3112,-12,-14,-12,
717 &3322,3312,2*3122,3212,3322,3312,3122,3322,3312,-12,2*4122,7*-11,
718 &7*-13,2*2224,2*2212,2*2214,2*3122,2*3212,2*3214,5*3222,4*3224,
719 &2*3322,3324,2*2224,7*2212,5*2214,2*2112,2*2114,2*3122,2*3212,
720 &2*3214,2*3222,2*3224,4*2,3,2*2,1,2*2,-11,-13,2*2,4*4122,-11,-13,
721 &2*2,3*4132,3*4232,-11,-13,2*2,4332,-11,-13,2*2,-11,-13,2*2,-11,
722 &-13,2*2,-11,-13,2*2,-11,-13,2*2,-11,-13,2*2,-11,-13,2*2,2*5122,
723 &-12,-14,-16,5*4122,441,443,20443,2*-2,2*-4,-2,-4,-12,-14,-16,
724 &2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,4*5122,-12,-14,-16,
725 &2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,2*5132,2*5232,-12,
726 &-14,-16,2*-2,2*-4,-2,-4,5332,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,
727 &-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,
728 &2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,
729 &-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12/
730 DATA (KFDP(I,1),I=1429,1710)/-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,
731 &2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,
732 &2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,
733 &-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,-12,
734 &-14,-16,2*-2,2*-4,-2,-4,-12,-14,-16,2*-2,2*-4,-2,-4,221,223,221,
735 &223,211,111,321,130,310,213,113,-213,321,311,321,311,323,313,
736 &2*311,321,311,321,313,323,321,211,111,321,130,310,2*211,313,-313,
737 &323,-323,421,411,423,413,411,421,413,423,411,421,423,413,443,
738 &2*82,521,511,523,513,511,521,513,523,521,511,523,513,511,521,513,
739 &523,553,2*21,213,-213,113,213,10211,10111,-10211,2*221,213,2*113,
740 &-213,2*321,2*311,113,323,2*313,323,313,-313,323,-323,423,2*413,
741 &2*423,413,443,82,523,2*513,2*523,2*513,523,553,21,11,13,82,4*443,
742 &10441,20443,445,441,11,13,15,1,2,3,4,21,22,2*553,10551,20553,555,
743 &1000039,-1000024,-1000037,1000022,1000023,1000025,1000035,
744 &1000002,2000002,1000002,2000002,1000021,1000039,1000024,1000037,
745 &1000022,1000023,1000025,1000035,1000001,2000001,1000001,2000001,
746 &1000021,1000039,-1000024,-1000037,1000022,1000023,1000025,
747 &1000035,1000004,2000004,1000004,2000004,1000021,1000039,1000024,
748 &1000037,1000022,1000023,1000025,1000035,1000003,2000003,1000003,
749 &2000003,1000021,1000039,-1000024,-1000037,1000022,1000023/
750 DATA (KFDP(I,1),I=1711,1900)/1000025,1000035,1000006,2000006,
751 &1000006,2000006,1000021,1000039,1000024,1000037,1000022,1000023,
752 &1000025,1000035,1000005,2000005,1000005,2000005,1000021,1000022,
753 &1000039,-1000024,-1000037,1000022,1000023,1000025,1000035,
754 &1000012,2000012,1000012,2000012,1000039,1000024,1000037,1000022,
755 &1000023,1000025,1000035,1000011,2000011,1000011,2000011,1000039,
756 &-1000024,-1000037,1000022,1000023,1000025,1000035,1000014,
757 &2000014,1000014,2000014,1000039,1000024,1000037,1000022,1000023,
758 &1000025,1000035,1000013,2000013,1000013,2000013,1000039,-1000024,
759 &-1000037,1000022,1000023,1000025,1000035,1000016,2000016,1000016,
760 &2000016,1000039,1000024,1000037,1000022,1000023,1000025,1000035,
761 &1000015,2000015,1000015,2000015,1000039,1000001,-1000001,2000001,
762 &-2000001,1000002,-1000002,2000002,-2000002,1000003,-1000003,
763 &2000003,-2000003,1000004,-1000004,2000004,-2000004,1000005,
764 &-1000005,2000005,-2000005,1000006,-1000006,2000006,-2000006,
765 &6*1000022,6*1000023,6*1000025,6*1000035,1000024,-1000024,1000024,
766 &-1000024,1000024,-1000024,1000037,-1000037,1000037,-1000037,
767 &1000037,-1000037,10*1000039,16*1000022,1000024,-1000024,1000024,
768 &-1000024,1000024,-1000024,1000024,-1000024,1000024,-1000024,
769 &1000024,-1000024,1000037,-1000037,1000037,-1000037,1000037/
770 DATA (KFDP(I,1),I=1901,2095)/-1000037,1000037,-1000037,1000037,
771 &-1000037,1000037,-1000037,1000024,-1000024,1000037,-1000037,
772 &1000001,-1000001,2000001,-2000001,1000002,-1000002,2000002,
773 &-2000002,1000003,-1000003,2000003,-2000003,1000004,-1000004,
774 &2000004,-2000004,1000005,-1000005,2000005,-2000005,1000006,
775 &-1000006,2000006,-2000006,1000011,-1000011,2000011,-2000011,
776 &1000012,-1000012,2000012,-2000012,1000013,-1000013,2000013,
777 &-2000013,1000014,-1000014,2000014,-2000014,1000015,-1000015,
778 &2000015,-2000015,1000016,-1000016,2000016,-2000016,5*1000021,
779 &2*1000039,6*1000022,6*1000023,6*1000025,6*1000035,1000022,
780 &1000023,1000025,1000035,1000002,2000002,-1000001,-2000001,
781 &1000004,2000004,-1000003,-2000003,1000006,2000006,-1000005,
782 &-2000005,1000012,2000012,-1000011,-2000011,1000014,2000014,
783 &-1000013,-2000013,1000016,2000016,-1000015,-2000015,2*1000021,
784 &5*1000039,16*1000022,16*1000023,1000024,-1000024,1000024,
785 &-1000024,1000024,-1000024,1000024,-1000024,1000024,-1000024,
786 &1000024,-1000024,1000037,-1000037,1000037,-1000037,1000037,
787 &-1000037,1000037,-1000037,1000037,-1000037,1000037,-1000037,
788 &1000024,-1000024,1000037,-1000037,1000001,-1000001,2000001,
789 &-2000001,1000002,-1000002,2000002,-2000002,1000003,-1000003/
790 DATA (KFDP(I,1),I=2096,2323)/2000003,-2000003,1000004,-1000004,
791 &2000004,-2000004,1000005,-1000005,2000005,-2000005,1000006,
792 &-1000006,2000006,-2000006,1000011,-1000011,2000011,-2000011,
793 &1000012,-1000012,2000012,-2000012,1000013,-1000013,2000013,
794 &-2000013,1000014,-1000014,2000014,-2000014,1000015,-1000015,
795 &2000015,-2000015,1000016,-1000016,2000016,-2000016,5*1000021,
796 &5*1000039,16*1000022,16*1000023,16*1000025,1000024,-1000024,
797 &1000024,-1000024,1000024,-1000024,1000024,-1000024,1000024,
798 &-1000024,1000024,-1000024,1000037,-1000037,1000037,-1000037,
799 &1000037,-1000037,1000037,-1000037,1000037,-1000037,1000037,
800 &-1000037,1000024,-1000024,1000037,-1000037,1000001,-1000001,
801 &2000001,-2000001,1000002,-1000002,2000002,-2000002,1000003,
802 &-1000003,2000003,-2000003,1000004,-1000004,2000004,-2000004,
803 &1000005,-1000005,2000005,-2000005,1000006,-1000006,2000006,
804 &-2000006,1000011,-1000011,2000011,-2000011,1000012,-1000012,
805 &2000012,-2000012,1000013,-1000013,2000013,-2000013,1000014,
806 &-1000014,2000014,-2000014,1000015,-1000015,2000015,-2000015,
807 &1000016,-1000016,2000016,-2000016,5*1000021,2*1000039,15*1000024,
808 &6*1000022,6*1000023,6*1000025,6*1000035,1000022,1000023,1000025,
809 &1000035,1000002,2000002,-1000001,-2000001,1000004,2000004/
810 DATA (KFDP(I,1),I=2324,4000)/-1000003,-2000003,1000006,2000006,
811 &-1000005,-2000005,1000012,2000012,-1000011,-2000011,1000014,
812 &2000014,-1000013,-2000013,1000016,2000016,-1000015,-2000015,
813 &2*1000021,1000039,-1000024,-1000037,1000022,1000023,1000025,
814 &1000035,4*1000001,1000002,2000002,1000002,2000002,1000021,
815 &1000039,1000024,1000037,1000022,1000023,1000025,1000035,
816 &4*1000002,1000001,2000001,1000001,2000001,1000021,1000039,
817 &-1000024,-1000037,1000022,1000023,1000025,1000035,4*1000003,
818 &1000004,2000004,1000004,2000004,1000021,1000039,1000024,1000037,
819 &1000022,1000023,1000025,1000035,4*1000004,1000003,2000003,
820 &1000003,2000003,1000021,1000039,-1000024,-1000037,1000022,
821 &1000023,1000025,1000035,4*1000005,1000006,2000006,1000006,
822 &2000006,1000021,1000039,1000024,1000037,1000022,1000023,1000025,
823 &1000035,4*1000006,1000005,2000005,1000005,2000005,1000021,
824 &1000039,-1000024,-1000037,1000022,1000023,1000025,1000035,
825 &4*1000011,1000012,2000012,1000012,2000012,1000039,-1000024,
826 &-1000037,1000022,1000023,1000025,1000035,4*1000013,1000014,
827 &2000014,1000014,2000014,1000039,-1000024,-1000037,1000022,
828 &1000023,1000025,1000035,4*1000015,1000016,2000016,1000016,
829 &2000016,21,22,23,-24,21,22,23,24,22,23,-24,23,24,1503*0/
830 DATA (KFDP(I,2),I= 1, 337)/3*1,2,4,6,8,1,3*2,1,3,5,7,2,3*3,2,4,
831 &6,8,3,3*4,1,3,5,7,4,3*5,2,4,6,8,5,3*6,1,3,5,7,6,5,4*1000006,3*7,
832 &2,4,6,8,7,4,6,3*8,1,3,5,7,8,5,7,2*11,12,11,12,2*11,2*13,14,13,14,
833 &13,11,13,-211,-213,-211,-213,-211,-213,-211,-213,2*-211,-321,
834 &-323,-321,2*-323,3*-321,4*-211,-213,-211,-213,-211,-213,-211,
835 &-213,-211,-213,3*-211,-213,4*-211,-323,-321,2*-211,2*-321,3*-211,
836 &2*15,16,15,16,15,2*17,18,17,2*18,2*17,-1,-2,-3,-4,-5,-6,-7,-8,21,
837 &-1,-2,-3,-4,-5,-6,-7,-8,-11,-13,-15,-17,-1,-2,-3,-4,-5,-6,-7,-8,
838 &-11,-12,-13,-14,-15,-16,-17,-18,2,4,6,8,2,4,6,8,2,4,6,8,2,4,6,8,
839 &12,14,16,18,-1,-2,-3,-4,-5,-6,-7,-8,-11,-13,-15,-17,21,22,2*23,
840 &-24,2*1000022,1000023,1000022,1000023,1000025,1000022,1000023,
841 &1000025,1000035,-1000024,-1000037,-1000024,-1000037,-1000001,
842 &2*-2000001,2000001,-1000002,2*-2000002,2000002,-1000003,
843 &2*-2000003,2000003,-1000004,2*-2000004,2000004,-1000005,
844 &2*-2000005,2000005,-1000006,2*-2000006,2000006,-1000011,
845 &2*-2000011,2000011,-1000012,2*-2000012,2000012,-1000013,
846 &2*-2000013,2000013,-1000014,2*-2000014,2000014,-1000015,
847 &2*-2000015,2000015,-1000016,2*-2000016,2000016,-1,-2,-3,-4,-5,-6,
848 &-7,-8,-11,-12,-13,-14,-15,-16,-17,-18,-24,-37,22,25,2*36,2,4,6,8,
849 &2,4,6,8,2,4,6,8,2,4,6,8,12,14,16,18,23,22,25,-1,-2,-3,-4,-5,-6/
850 DATA (KFDP(I,2),I= 338, 524)/-7,-8,-11,-13,-15,-17,21,22,2*23,
851 &-24,2*25,36,2*1000022,1000023,1000022,1000023,1000025,1000022,
852 &1000023,1000025,1000035,-1000024,-1000037,-1000024,-1000037,
853 &-1000001,2*-2000001,2000001,-1000002,2*-2000002,2000002,-1000003,
854 &2*-2000003,2000003,-1000004,2*-2000004,2000004,-1000005,
855 &2*-2000005,2000005,-1000006,2*-2000006,2000006,-1000011,
856 &2*-2000011,2000011,-1000012,2*-2000012,2000012,-1000013,
857 &2*-2000013,2000013,-1000014,2*-2000014,2000014,-1000015,
858 &2*-2000015,2000015,-1000016,2*-2000016,2000016,-1,-2,-3,-4,-5,-6,
859 &-7,-8,-11,-13,-15,-17,21,22,2*23,-24,25,2*1000022,1000023,
860 &1000022,1000023,1000025,1000022,1000023,1000025,1000035,-1000024,
861 &-1000037,-1000024,-1000037,-1000001,2*-2000001,2000001,-1000002,
862 &2*-2000002,2000002,-1000003,2*-2000003,2000003,-1000004,
863 &2*-2000004,2000004,-1000005,2*-2000005,2000005,-1000006,
864 &2*-2000006,2000006,-1000011,2*-2000011,2000011,-1000012,
865 &2*-2000012,2000012,-1000013,2*-2000013,2000013,-1000014,
866 &2*-2000014,2000014,-1000015,2*-2000015,2000015,-1000016,
867 &2*-2000016,2000016,2,4,6,8,12,14,16,18,25,1000024,1000037,
868 &1000024,1000037,1000024,1000037,1000024,1000037,2*-1000005,
869 &2*-2000005,1000002,1000004,1000012,1000014,2*1000016,-5,-6,21,11/
870 DATA (KFDP(I,2),I= 525, 940)/-3,-4,-5,-6,-7,-8,-13,-15,-17,-4,-5,
871 &-11,-13,-15,-5,-3,12,14,16,-24,-52,-24,-52,-1,-2,-3,-4,-5,-6,-7,
872 &-8,-11,-12,-13,-14,-15,-16,-17,-18,23,51,23,51,2,4,6,8,2,4,6,8,2,
873 &4,6,8,2,4,6,8,12,14,16,18,2*51,-1,-2,-3,-4,-5,-6,-7,-8,-11,-12,
874 &-13,-14,-15,-16,-17,-18,-82,12,14,-1,-3,11,13,15,1,4,3,4,1,3,22,
875 &11,-211,2*22,-13,-11,-211,211,111,211,-321,130,310,22,2*111,-211,
876 &11,-11,13,-13,-211,111,22,14,12,111,22,111,3*211,-311,22,211,22,
877 &111,-211,211,11,-211,13,22,-211,111,-211,22,111,-11,-211,111,
878 &2*-211,-321,130,310,221,111,-211,111,2*0,-211,111,22,-211,111,
879 &-211,111,-211,211,-213,113,223,221,14,111,211,111,-11,-13,211,
880 &111,22,211,111,211,111,2*211,213,113,223,221,22,-211,111,113,223,
881 &22,111,-321,310,211,111,2*-211,221,22,-11,-13,-211,-321,130,310,
882 &221,-211,111,11*12,11*14,2*211,2*213,211,20213,2*321,2*323,211,
883 &213,211,213,211,213,211,213,211,213,211,213,3*211,213,211,2*321,
884 &8*211,2*113,3*211,111,22,211,111,211,111,4*211,8*12,8*14,2*211,
885 &2*213,2*111,221,2*113,223,333,20213,211,2*321,323,2*311,313,-211,
886 &111,113,2*211,321,2*211,311,321,310,211,-211,4*211,321,4*211,113,
887 &2*211,-321,111,22,-211,111,-211,111,-211,211,-211,211,16,5*12,
888 &5*14,3*211,3*213,211,2*111,2*113,2*-311,2*-313,-2112,3*321,323,
889 &2*-1,22,111,321,311,321,311,-82,-11,-13,-82,22,-82,6*-11,6*-13/
890 DATA (KFDP(I,2),I= 941,1318)/2*-15,211,213,20213,211,213,20213,
891 &431,433,431,433,311,313,311,313,311,313,-1,-4,-3,-4,-1,-3,22,
892 &-211,111,-211,111,-211,211,-211,211,6*-11,6*-13,2*-15,211,213,
893 &20213,211,213,20213,431,433,431,433,321,323,321,323,321,323,-1,
894 &-4,-3,-4,-1,-3,22,211,111,211,111,4*211,6*-11,6*-13,2*-15,211,
895 &213,20213,211,213,20213,431,433,431,433,221,331,333,221,331,333,
896 &221,331,333,-1,-4,-3,-4,-1,-3,22,-321,-311,-321,-311,-15,-3,-1,
897 &2*-11,2*-13,2*-15,-1,-4,-3,-4,-3,-4,-1,-4,2*12,2*14,2,3,2,3,2*12,
898 &2*14,2,1,22,411,421,411,421,21,-11,-13,-15,-1,-2,-3,-4,2*21,22,
899 &21,2*-211,111,22,111,211,22,211,-211,11,2*-211,111,-211,111,22,
900 &11,22,111,-211,211,111,211,22,211,111,211,-211,22,11,13,11,-211,
901 &2*111,2*22,111,211,-321,-211,111,11,2*-211,7*12,7*14,-321,-323,
902 &-311,-313,-311,-313,211,213,211,213,211,213,111,221,331,113,223,
903 &111,221,113,223,321,323,321,-211,-213,111,221,331,113,223,333,
904 &10221,111,221,331,113,223,211,213,211,213,321,323,321,323,321,
905 &323,311,313,311,313,2*-1,-3,-1,2203,3201,3203,2203,2101,2103,12,
906 &14,-1,-3,2*111,2*211,12,14,-1,-3,22,111,2*22,111,22,12,14,-1,-3,
907 &22,12,14,-1,-3,12,14,-1,-3,12,14,-1,-3,12,14,-1,-3,12,14,-1,-3,
908 &12,14,-1,-3,12,14,-1,-3,2*-211,11,13,15,-211,-213,-20213,-431,
909 &-433,3*3122,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1/
910 DATA (KFDP(I,2),I=1319,1774)/3,2*111,2*211,11,13,15,1,4,3,4,1,3,
911 &11,13,15,1,4,3,4,1,3,4*22,11,13,15,1,4,3,4,1,3,22,11,13,15,1,4,3,
912 &4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,
913 &1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,
914 &3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,
915 &11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,
916 &11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,
917 &11,13,15,1,4,3,4,1,3,11,13,15,1,4,3,4,1,3,2*111,2*211,-211,111,
918 &-321,130,310,-211,111,211,-211,111,-213,113,-211,111,223,211,111,
919 &213,113,211,111,223,-211,111,-321,130,310,2*-211,-311,311,-321,
920 &321,211,111,211,111,-211,111,-211,111,311,2*321,311,22,2*-82,
921 &-211,111,-211,111,211,111,211,111,-321,-311,-321,-311,411,421,
922 &411,421,22,2*21,-211,2*211,111,-211,111,2*211,111,-211,211,111,
923 &211,-321,2*-311,-321,22,-211,111,211,111,-311,311,-321,321,211,
924 &111,-211,111,321,311,22,-82,-211,111,211,111,-321,-311,411,421,
925 &22,21,-11,-13,-82,211,111,221,111,4*22,-11,-13,-15,-1,-2,-3,-4,
926 &2*21,211,111,3*22,1,2*2,4*1,2*-24,2*-37,1,2,2*1,4*2,2*24,2*37,2,
927 &3,2*4,4*3,2*-24,2*-37,3,4,2*3,4*4,2*24,2*37,4,5,2*6,4*5,2*-24,
928 &2*-37,5,6,2*5,4*6,2*24,2*37,6,4,11,2*12,4*11,2*-24,2*-37,12,2*11,
929 &4*12,2*24,2*37,13,2*14,4*13,2*-24,2*-37,14,2*13,4*14,2*24,2*37/
930 DATA (KFDP(I,2),I=1775,2218)/15,2*16,4*15,2*-24,2*-37,16,2*15,
931 &4*16,2*24,2*37,21,-1,1,-1,1,-2,2,-2,2,-3,3,-3,3,-4,4,-4,4,-5,5,
932 &-5,5,-6,6,-6,6,1,3,5,2,4,6,1,3,5,2,4,6,1,3,5,2,4,6,1,3,5,2,4,6,1,
933 &-1,3,-3,5,-5,1,-1,3,-3,5,-5,22,23,25,35,36,22,23,25,35,36,22,23,
934 &11,13,15,12,14,16,1,3,5,2,4,25,35,36,-24,24,11,-11,13,-13,15,-15,
935 &1,-1,3,-3,-24,24,11,-11,13,-13,15,-15,1,-1,3,-3,-37,37,-37,37,-1,
936 &1,-1,1,-2,2,-2,2,-3,3,-3,3,-4,4,-4,4,-5,5,-5,5,-6,6,-6,6,-11,11,
937 &-11,11,-12,12,-12,12,-13,13,-13,13,-14,14,-14,14,-15,15,-15,15,
938 &-16,16,-16,16,1,3,5,2,4,24,37,24,-11,-13,-15,-1,-3,24,-11,-13,
939 &-15,-1,-3,24,-11,-13,-15,-1,-3,24,-11,-13,-15,-1,-3,4*37,2*-1,
940 &2*2,2*-3,2*4,2*-5,2*6,2*-11,2*12,2*-13,2*14,2*-15,2*16,-1,-3,22,
941 &23,25,35,36,22,23,11,13,15,12,14,16,1,3,5,2,4,25,35,36,22,23,11,
942 &13,15,12,14,16,1,3,5,2,4,25,35,36,-24,24,11,-11,13,-13,15,-15,1,
943 &-1,3,-3,-24,24,11,-11,13,-13,15,-15,1,-1,3,-3,-37,37,-37,37,-1,1,
944 &-1,1,-2,2,-2,2,-3,3,-3,3,-4,4,-4,4,-5,5,-5,5,-6,6,-6,6,-11,11,
945 &-11,11,-12,12,-12,12,-13,13,-13,13,-14,14,-14,14,-15,15,-15,15,
946 &-16,16,-16,16,1,3,5,2,4,22,23,25,35,36,22,23,11,13,15,12,14,16,1,
947 &3,5,2,4,25,35,36,22,23,11,13,15,12,14,16,1,3,5,2,4,25,35,36,22,
948 &23,11,13,15,12,14,16,1,3,5,2,4,25,35,36,-24,24,11,-11,13,-13,15,
949 &-15,1,-1,3,-3,-24,24,11,-11,13,-13,15,-15,1,-1,3,-3,-37,37,-37/
950 DATA (KFDP(I,2),I=2219,4000)/37,-1,1,-1,1,-2,2,-2,2,-3,3,-3,3,-4,
951 &4,-4,4,-5,5,-5,5,-6,6,-6,6,-11,11,-11,11,-12,12,-12,12,-13,13,
952 &-13,13,-14,14,-14,14,-15,15,-15,15,-16,16,-16,16,1,3,5,2,4,24,37,
953 &23,11,13,15,12,14,16,1,3,5,2,4,25,35,36,24,-11,-13,-15,-1,-3,24,
954 &-11,-13,-15,-1,-3,24,-11,-13,-15,-1,-3,24,-11,-13,-15,-1,-3,4*37,
955 &2*-1,2*2,2*-3,2*4,2*-5,2*6,2*-11,2*12,2*-13,2*14,2*-15,2*16,-1,
956 &-3,1,2*2,4*1,23,25,35,36,2*-24,2*-37,1,2,2*1,4*2,23,25,35,36,
957 &2*24,2*37,2,3,2*4,4*3,23,25,35,36,2*-24,2*-37,3,4,2*3,4*4,23,25,
958 &35,36,2*24,2*37,4,5,2*6,4*5,23,25,35,36,2*-24,2*-37,5,6,2*5,4*6,
959 &23,25,35,36,2*24,2*37,6,11,2*12,4*11,23,25,35,36,2*-24,2*-37,13,
960 &2*14,4*13,23,25,35,36,2*-24,2*-37,15,2*16,4*15,23,25,35,36,2*-24,
961 &2*-37,3*1,4*2,1,2*11,2*12,11,1503*0/
962 DATA (KFDP(I,3),I= 1,1087)/79*0,14,6*0,2*16,2*0,6*111,310,130,
963 &2*0,3*111,310,130,321,113,211,223,221,2*113,2*211,2*223,2*221,
964 &2*113,221,2*113,2*213,-213,113,2*111,310,130,310,130,2*310,130,
965 &470*0,4*3,4*4,1,4,3,2*2,0,-11,8*0,-211,5*0,2*111,211,-211,211,
966 &-211,10*0,111,4*0,2*111,-211,-11,11,-13,22,111,3*0,22,3*0,111,
967 &211,4*0,111,11*0,111,-211,6*0,-211,3*111,7*0,111,-211,5*0,2*221,
968 &3*0,111,5*0,111,11*0,-311,-313,-311,-321,-313,-323,111,221,331,
969 &113,223,-311,-313,-311,-321,-313,-323,111,221,331,113,223,22*0,
970 &111,113,2*211,-211,-311,211,111,3*211,-211,7*211,7*0,111,-211,
971 &111,-211,-321,-323,-311,-321,-313,-323,-211,-213,-321,-323,-311,
972 &-321,-313,-323,-211,-213,22*0,111,113,-311,2*-211,211,-211,310,
973 &-211,2*111,211,2*-211,-321,-211,2*211,-211,111,-211,2*211,6*0,
974 &111,-211,111,-211,0,221,331,333,321,311,221,331,333,321,311,20*0,
975 &3,13*0,-411,-413,-10413,-10411,-20413,-415,-411,-413,-10413,
976 &-10411,-20413,-415,-411,-413,16*0,-4,-1,-4,-3,2*-2,5*0,111,-211,
977 &111,-211,-421,-423,-10423,-10421,-20423,-425,-421,-423,-10423,
978 &-10421,-20423,-425,-421,-423,16*0,-4,-1,-4,-3,2*-2,5*0,111,-211,
979 &111,-211,-431,-433,-10433,-10431,-20433,-435,-431,-433,-10433,
980 &-10431,-20433,-435,-431,-433,19*0,-4,-1,-4,-3,2*-2,8*0,441,443,
981 &441,443,441,443,-4,-1,-4,-3,-4,-3,-4,-1,531,533,531,533,3,2,3,2/
982 DATA (KFDP(I,3),I=1088,2186)/511,513,511,513,1,2,13*0,2*21,11*0,
983 &2112,6*0,2212,12*0,2*3122,3212,10*0,3322,2*0,3122,3212,3214,2112,
984 &2114,2212,2112,3122,3212,3214,2112,2114,2212,2112,52*0,3*3,1,6*0,
985 &4*3,4*0,4*3,6*0,4*3,0,28*3,2*0,3*4122,8*0,4,1,4,3,2*2,4*4,1,4,3,
986 &2*2,4*4,1,4,3,2*2,4*0,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*0,4*4,1,4,3,
987 &2*2,0,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,
988 &4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,
989 &3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,
990 &4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,3,2*2,4*4,1,4,
991 &3,2*2,31*0,211,111,45*0,-211,2*111,-211,3*111,-211,111,211,30*0,
992 &-211,111,13*0,2*21,-211,111,167*0,-1,-3,-5,-2,-4,-6,-1,-3,-5,-2,
993 &-4,-6,-1,-3,-5,-2,-4,-6,-1,-3,-5,-2,-4,-6,-2,2,-4,4,-6,6,-2,2,-4,
994 &4,-6,6,12*0,-11,-13,-15,-12,-14,-16,-1,-3,-5,-2,-4,5*0,-12,12,
995 &-14,14,-16,16,-2,2,-4,4,2*0,-12,12,-14,14,-16,16,-2,2,-4,4,52*0,
996 &-1,-3,-5,-2,-4,3*0,12,14,16,2,4,0,12,14,16,2,4,0,12,14,16,2,4,0,
997 &12,14,16,2,4,28*0,2,4,7*0,-11,-13,-15,-12,-14,-16,-1,-3,-5,-2,-4,
998 &5*0,-11,-13,-15,-12,-14,-16,-1,-3,-5,-2,-4,5*0,-12,12,-14,14,-16,
999 &16,-2,2,-4,4,2*0,-12,12,-14,14,-16,16,-2,2,-4,4,52*0,-1,-3,-5,-2,
1000 &-4,7*0,-11,-13,-15,-12,-14,-16,-1,-3,-5,-2,-4,5*0,-11,-13,-15,
1001 &-12,-14,-16,-1,-3,-5,-2,-4,5*0,-11,-13,-15,-12,-14,-16,-1,-3,-5/
1002 DATA (KFDP(I,3),I=2187,4000)/-2,-4,5*0,-12,12,-14,14,-16,16,-2,2,
1003 &-4,4,2*0,-12,12,-14,14,-16,16,-2,2,-4,4,52*0,-1,-3,-5,-2,-4,3*0,
1004 &-11,-13,-15,-12,-14,-16,-1,-3,-5,-2,-4,4*0,12,14,16,2,4,0,12,14,
1005 &16,2,4,0,12,14,16,2,4,0,12,14,16,2,4,28*0,2,4,1657*0/
1006 DATA (KFDP(I,4),I= 1,4000)/92*0,4*111,6*0,111,2*0,-211,0,-211,
1007 &3*0,111,2*-211,0,111,0,2*111,113,221,2*111,-213,-211,211,113,
1008 &6*111,310,2*130,470*0,13*81,41*0,-11,10*0,111,-211,4*0,111,62*0,
1009 &111,211,111,211,7*0,111,211,111,211,35*0,2*-211,2*111,211,111,
1010 &-211,2*211,2*-211,13*0,-211,111,-211,111,4*0,-211,111,-211,111,
1011 &34*0,111,-211,3*111,3*-211,2*111,3*-211,14*0,-321,-311,3*0,-321,
1012 &-311,20*0,-3,43*0,6*1,39*0,6*2,42*0,6*3,14*0,8*4,4*0,4*-5,4*0,
1013 &2*-5,67*0,-211,111,5*0,-211,111,52*0,2101,2103,2*2101,6*0,4*81,
1014 &4*0,4*81,6*0,4*81,0,28*81,13*0,6*2101,18*81,4*0,18*81,4*0,9*81,0,
1015 &162*81,31*0,-211,111,2450*0/
1016 DATA (KFDP(I,5),I= 1,4000)/94*0,2*111,17*0,111,7*0,2*111,0,
1017 &3*111,0,111,665*0,-211,2*111,-211,111,-211,111,65*0,111,-211,
1018 &3*111,-211,111,3127*0/
1019
1020C...PYDAT4, with particle names (character strings).
1021 DATA (CHAF(I,1),I= 1, 190)/'d','u','s','c','b','t','b''','t''',
1022 &2*' ','e-','nu_e','mu-','nu_mu','tau-','nu_tau','tau''-',
1023 &'nu''_tau',2*' ','g','gamma','Z0','W+','h0',2*' ','reggeon',
1024 &'pomeron',2*' ','Z''0','Z"0','W''+','H0','A0','H+','eta_tech0',
1025 &'LQ_ue','R0',10*' ','pi_tech0','pi_tech+','pi''_tech0',
1026 &'rho_tech0','rho_tech+','omega_tech',24*' ','specflav',
1027 &'rndmflav','phasespa','c-hadron','b-hadron',5*' ','cluster',
1028 &'string','indep.','CMshower','SPHEaxis','THRUaxis','CLUSjet',
1029 &'CELLjet','table',' ','rho_diff0','pi0','rho0','a_20','K_L0',
1030 &'pi_diffr+','pi+','rho+','a_2+','omega_di','eta','omega','f_2',
1031 &'K_S0','K0','K*0','K*_20','K+','K*+','K*_2+','phi_diff','eta''',
1032 &'phi','f''_2','D+','D*+','D*_2+','D0','D*0','D*_20','D_s+',
1033 &'D*_s+','D*_2s+','J/psi_di','eta_c','J/psi','chi_2c','B0','B*0',
1034 &'B*_20','B+','B*+','B*_2+','B_s0','B*_s0','B*_2s0','B_c+',
1035 &'B*_c+','B*_2c+','eta_b','Upsilon','chi_2b','dd_1','Delta-',
1036 &'ud_0','ud_1','n_diffr0','n0','Delta0','uu_1','p_diffr+','p+',
1037 &'Delta+','Delta++','sd_0','sd_1','Sigma-','Sigma*-','Lambda0',
1038 &'su_0','su_1','Sigma0','Sigma*0','Sigma+','Sigma*+','ss_1','Xi-',
1039 &'Xi*-','Xi0','Xi*0','Omega-','cd_0','cd_1','Sigma_c0',
1040 &'Sigma*_c0','Lambda_c+','Xi_c0','cu_0','cu_1','Sigma_c+'/
1041 DATA (CHAF(I,1),I= 191, 317)/'Sigma*_c+','Sigma_c++',
1042 &'Sigma*_c++','Xi_c+','cs_0','cs_1','Xi''_c0','Xi*_c0','Xi''_c+',
1043 &'Xi*_c+','Omega_c0','Omega*_c0','cc_1','Xi_cc+','Xi*_cc+',
1044 &'Xi_cc++','Xi*_cc++','Omega_cc+','Omega*_cc+','Omega*_ccc++',
1045 &'bd_0','bd_1','Sigma_b-','Sigma*_b-','Lambda_b0','Xi_b-',
1046 &'Xi_bc0','bu_0','bu_1','Sigma_b0','Sigma*_b0','Sigma_b+',
1047 &'Sigma*_b+','Xi_b0','Xi_bc+','bs_0','bs_1','Xi''_b-','Xi*_b-',
1048 &'Xi''_b0','Xi*_b0','Omega_b-','Omega*_b-','Omega_bc0','bc_0',
1049 &'bc_1','Xi''_bc0','Xi*_bc0','Xi''_bc+','Xi*_bc+','Omega''_bc0',
1050 &'Omega*_bc0','Omega_bcc+','Omega*_bcc+','bb_1','Xi_bb-',
1051 &'Xi*_bb-','Xi_bb0','Xi*_bb0','Omega_bb-','Omega*_bb-',
1052 &'Omega_bbc0','Omega*_bbc0','Omega*_bbb-','a_00','b_10','a_0+',
1053 &'b_1+','f_0','h_1','K*_00','K_10','K*_0+','K_1+','f''_0','h''_1',
1054 &'D*_0+','D_1+','D*_00','D_10','D*_0s+','D_1s+','chi_0c','h_1c',
1055 &'B*_00','B_10','B*_0+','B_1+','B*_0s0','B_1s0','B*_0c+','B_1c+',
1056 &'chi_0b','h_1b','a_10','a_1+','f_1','K*_10','K*_1+','f''_1',
1057 &'D*_1+','D*_10','D*_1s+','chi_1c','B*_10','B*_1+','B*_1s0',
1058 &'B*_1c+','chi_1b','psi''','Upsilon''','~d_L','~u_L','~s_L',
1059 &'~c_L','~b_1','~t_1','~e_L-','~nu_eL','~mu_L-','~nu_muL',
1060 &'~tau_1-','~nu_tauL','~g','~chi_10','~chi_20','~chi_1+'/
1061 DATA (CHAF(I,1),I= 318, 500)/'~chi_30','~chi_40','~chi_2+',
1062 &'~gravitino','~d_R','~u_R','~s_R','~c_R','~b_2','~t_2','~e_R-',
1063 &'~nu_eR','~mu_R-','~nu_muR','~tau_2-','~nu_tauR','d*','u*','e*-',
1064 &'nu*_e0',163*' '/
1065 DATA (CHAF(I,2),I= 1, 206)/'dbar','ubar','sbar','cbar','bbar',
1066 &'tbar','b''bar','t''bar',2*' ','e+','nu_ebar','mu+','nu_mubar',
1067 &'tau+','nu_taubar','tau''+','nu''_taubar',5*' ','W-',9*' ',
1068 &'W''-',2*' ','H-',' ','LQ_uebar','Rbar0',11*' ','pi_tech-',2*' ',
1069 &'rho_tech-',26*' ','rndmflavbar',' ','c-hadronbar','b-hadronbar',
1070 &20*' ','pi_diffr-','pi-','rho-','a_2-',5*' ','Kbar0','K*bar0',
1071 &'K*_2bar0','K-','K*-','K*_2-',4*' ','D-','D*-','D*_2-','Dbar0',
1072 &'D*bar0','D*_2bar0','D_s-','D*_s-','D*_2s-',4*' ','Bbar0',
1073 &'B*bar0','B*_2bar0','B-','B*-','B*_2-','B_sbar0','B*_sbar0',
1074 &'B*_2sbar0','B_c-','B*_c-','B*_2c-',3*' ','dd_1bar','Deltabar+',
1075 &'ud_0bar','ud_1bar','n_diffrbar0','nbar0','Deltabar0','uu_1bar',
1076 &'p_diffrbar-','pbar-','Deltabar-','Deltabar--','sd_0bar',
1077 &'sd_1bar','Sigmabar+','Sigma*bar+','Lambdabar0','su_0bar',
1078 &'su_1bar','Sigmabar0','Sigma*bar0','Sigmabar-','Sigma*bar-',
1079 &'ss_1bar','Xibar+','Xi*bar+','Xibar0','Xi*bar0','Omegabar+',
1080 &'cd_0bar','cd_1bar','Sigma_cbar0','Sigma*_cbar0','Lambda_cbar-',
1081 &'Xi_cbar0','cu_0bar','cu_1bar','Sigma_cbar-','Sigma*_cbar-',
1082 &'Sigma_cbar--','Sigma*_cbar--','Xi_cbar-','cs_0bar','cs_1bar',
1083 &'Xi''_cbar0','Xi*_cbar0','Xi''_cbar-','Xi*_cbar-','Omega_cbar0',
1084 &'Omega*_cbar0','cc_1bar','Xi_ccbar-','Xi*_ccbar-','Xi_ccbar--'/
1085 DATA (CHAF(I,2),I= 207, 324)/'Xi*_ccbar--','Omega_ccbar-',
1086 &'Omega*_ccbar-','Omega*_cccbar-','bd_0bar','bd_1bar',
1087 &'Sigma_bbar+','Sigma*_bbar+','Lambda_bbar0','Xi_bbar+',
1088 &'Xi_bcbar0','bu_0bar','bu_1bar','Sigma_bbar0','Sigma*_bbar0',
1089 &'Sigma_bbar-','Sigma*_bbar-','Xi_bbar0','Xi_bcbar-','bs_0bar',
1090 &'bs_1bar','Xi''_bbar+','Xi*_bbar+','Xi''_bbar0','Xi*_bbar0',
1091 &'Omega_bbar+','Omega*_bbar+','Omega_bcbar0','bc_0bar','bc_1bar',
1092 &'Xi''_bcbar0','Xi*_bcbar0','Xi''_bcbar-','Xi*_bcbar-',
1093 &'Omega''_bcba','Omega*_bcbar0','Omega_bccbar-','Omega*_bccbar-',
1094 &'bb_1bar','Xi_bbbar+','Xi*_bbbar+','Xi_bbbar0','Xi*_bbbar0',
1095 &'Omega_bbbar+','Omega*_bbbar+','Omega_bbcbar0','Omega*_bbcbar0',
1096 &'Omega*_bbbbar+',2*' ','a_0-','b_1-',2*' ','K*_0bar0','K_1bar0',
1097 &'K*_0-','K_1-',2*' ','D*_0-','D_1-','D*_0bar0','D_1bar0',
1098 &'D*_0s-','D_1s-',2*' ','B*_0bar0','B_1bar0','B*_0-','B_1-',
1099 &'B*_0sbar0','B_1sbar0','B*_0c-','B_1c-',3*' ','a_1-',' ',
1100 &'K*_1bar0','K*_1-',' ','D*_1-','D*_1bar0','D*_1s-',' ',
1101 &'B*_1bar0','B*_1-','B*_1sbar0','B*_1c-',3*' ','~d_Lbar',
1102 &'~u_Lbar','~s_Lbar','~c_Lbar','~b_1bar','~t_1bar','~e_L+',
1103 &'~nu_eLbar','~mu_L+','~nu_muLbar','~tau_1+','~nu_tauLbar',3*' ',
1104 &'~chi_1-',2*' ','~chi_2-',' ','~d_Rbar','~u_Rbar','~s_Rbar'/
1105 DATA (CHAF(I,2),I= 325, 500)/'~c_Rbar','~b_2bar','~t_2bar',
1106 &'~e_R+','~nu_eRbar','~mu_R+','~nu_muRbar','~tau_2+',
1107 &'~nu_tauRbar','d*bar','u*bar','e*bar+','nu*_ebar0',163*' '/
1108
1109C...PYDATR, with initial values for the random number generator.
1110 DATA MRPY/19780503,0,0,97,33,0/
1111
1112C...Default values for allowed processes and kinematics constraints.
1113 DATA MSEL/1/
1114 DATA MSUB/500*0/
1115 DATA ((KFIN(I,J),J=-40,40),I=1,2)/16*0,4*1,4*0,6*1,5*0,5*1,0,
1116 &5*1,5*0,6*1,4*0,4*1,16*0,16*0,4*1,4*0,6*1,5*0,5*1,0,5*1,5*0,
1117 &6*1,4*0,4*1,16*0/
1118 DATA CKIN/
1119 & 2.0D0, -1.0D0, 0.0D0, -1.0D0, 1.0D0,
1120 & 1.0D0, -10D0, 10D0, -10D0, 10D0,
1121 1 -10D0, 10D0, -10D0, 10D0, -10D0,
1122 1 10D0, -1.0D0, 1.0D0, -1.0D0, 1.0D0,
1123 2 0.0D0, 1.0D0, 0.0D0, 1.0D0, -1.0D0,
1124 2 1.0D0, -1.0D0, 1.0D0, 0D0, 0D0,
1125 3 2.0D0, -1.0D0, 0D0, 0D0, 0.0D0,
1126 3 -1.0D0, 0.0D0, -1.0D0, 4.0D0, -1.0D0,
1127 4 12.0D0, -1.0D0, 12.0D0, -1.0D0, 12.0D0,
1128 4 -1.0D0, 12.0D0, -1.0D0, 0D0, 0D0,
1129 5 0.0D0, -1.0D0, 0.0D0, -1.0D0, 0.0D0,
1130 5 -1.0D0, 0D0, 0D0, 0D0, 0D0,
1131 6 140*0D0/
1132
1133C...Default values for main switches and parameters. Reset information.
1134 DATA (MSTP(I),I=1,100)/
1135 & 3, 1, 2, 0, 0, 0, 0, 0, 0, 0,
1136 1 1, 0, 1, 0, 5, 0, 0, 0, 0, 0,
1137 2 1, 0, 1, 0, 0, 0, 0, 0, 0, 1,
1138 3 1, 2, 0, 1, 0, 2, 1, 5, 2, 0,
1139 4 1, 1, 3, 7, 3, 1, 1, 0, 1, 0,
1140 5 4, 1, 3, 1, 5, 1, 1, 6, 1, 7,
1141 6 1, 3, 2, 2, 1, 1, 2, 0, 0, 0,
1142 7 1, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1143 8 1, 1, 100, 0, 0, 0, 0, 0, 0, 0,
1144 9 1, 4, 1, 2, 0, 0, 0, 0, 0, 0/
1145 DATA (MSTP(I),I=101,200)/
1146 & 3, 1, 0, 0, 0, 0, 0, 0, 0, 0,
1147 1 1, 1, 1, 0, 0, 0, 0, 0, 0, 0,
1148 2 0, 1, 2, 1, 1, 50, 0, 0, 10, 0,
1149 3 0, 4, 0, 1, 0, 0, 0, 0, 0, 0,
1150 4 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1151 5 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1152 6 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1153 7 0, 2, 0, 0, 0, 0, 0, 0, 0, 0,
1154 8 6, 115, 1998, 01, 27, 0, 0, 0, 0, 0,
1155 9 0, 0, 0, 0, 0, 0, 0, 0, 0, 0/
1156 DATA (PARP(I),I=1,100)/
1157 & 0.25D0, 10D0, 8*0D0,
1158 1 0D0, 0D0, 1.0D0, 0.01D0, 0.6D0, 1.0D0, 1.0D0, 3*0D0,
1159 2 10*0D0,
1160 3 1.5D0,2.0D0,0.075D0,1.0D0,0.2D0,0D0,2.0D0,0.70D0,0.006D0,0D0,
1161 4 0.02D0,2.0D0,0.10D0,1000D0,2054D0, 123D0, 246D0, 50D0, 2*0D0,
1162 5 1.0D0, 9*0D0,
1163 6 0.25D0, 1.0D0,0.25D0, 1.0D0, 2.0D0,1D-3, 4.0D0,1D-3,2*0D0,
1164 7 4.0D0, 0.25D0, 8*0D0,
1165 8 1.40D0,1.55D0,0.5D0, 0.2D0,0.33D0,0.66D0, 0.7D0, 0.5D0,2*0D0,
1166 9 0.44D0,0.20D0,2.0D0,1.0D0,0D0,3.0D0,1.0D0,0.75D0,0.44D0,2.0D0/
1167 DATA (PARP(I),I=101,200)/
1168 & 0.5D0, 0.28D0, 1.0D0, 0.8D0, 6*0D0,
1169 1 2.0D0, 3*0D0, 1.5D0, 0.5D0, 0.6D0, 2.5D0, 2.0D0, 1.0D0,
1170 2 1.0D0, 0.4D0, 8*0D0,
1171 3 0.01D0, 9*0D0,
1172 4 0.33333D0, 82D0, 1D0, 4D0, 200D0, 5*0D0,
1173 5 0D0, 0D0, 0D0, 0D0, 6*0D0,
1174 6 2.20D0, 23.6D0, 18.4D0, 11.5D0, 6*0D0,
1175 7 0D0, 0D0, 0D0, 1.0D0, 6*0D0,
1176 8 20*0D0/
1177 DATA MSTI/200*0/
1178 DATA PARI/200*0D0/
1179 DATA MINT/400*0/
1180 DATA VINT/400*0D0/
1181
1182C...Constants for the generation of the various processes.
1183 DATA (ISET(I),I=1,100)/
1184 & 1, 1, 1, -1, 3, -1, -1, 3, -2, 2,
1185 1 2, 2, 2, 2, 2, 2, -1, 2, 2, 2,
1186 2 -1, 2, 2, 2, 2, 2, -1, 2, 2, 2,
1187 3 2, -1, 2, 2, 2, 2, -1, -1, -1, -1,
1188 4 -1, -1, -1, -1, -1, -1, -1, -1, -1, -1,
1189 5 -1, -1, 2, 2, -1, -1, -1, 2, -1, -1,
1190 6 -1, -1, -1, -1, -1, -1, -1, 2, 2, 2,
1191 7 4, 4, 4, -1, -1, 4, 4, -1, -1, 2,
1192 8 2, 2, 2, 2, 2, 2, 2, 2, 2, -2,
1193 9 0, 0, 0, 0, 0, 9, -2, -2, -2, -2/
1194 DATA (ISET(I),I=101,200)/
1195 & -1, 1, 1, -2, -2, 2, 2, 2, -2, 2,
1196 1 2, 2, 2, 2, 2, -1, -1, -1, -2, -2,
1197 2 5, 5, 5, 5, -2, -2, -2, -2, -2, -2,
1198 3 -1, -2, -2, -2, -2, -2, -2, -2, -2, -2,
1199 4 1, 1, 1, 1, 1, -2, 1, 1, 1, -2,
1200 5 1, 1, 1, -2, -2, 1, 1, 1, -2, -2,
1201 6 2, 2, 2, 2, 2, 2, 2, 2, -2, -2,
1202 7 2, 2, 5, 5, -2, 2, 2, 5, 5, -2,
1203 8 5, 5, -2, -2, -2, 5, 5, -2, -2, -2,
1204 9 1, 1, 1, 2, -2, -2, -2, -2, -2, -2/
1205 DATA (ISET(I),I=201,300)/
1206 & 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
1207 1 2, 2, 2, 2, -2, 2, 2, 2, 2, 2,
1208 2 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
1209 3 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
1210 4 2, 2, 2, 2, -1, 2, 2, 2, 2, 2,
1211 5 2, 2, 2, 2, -1, 2, -1, 2, 2, -2,
1212 6 2, 2, 2, 2, 2, -1, -1, -1, -1, -1,
1213 7 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
1214 8 -2, -2, -2, -2, -2, -2, -2, -2, -2, -2,
1215 9 -2, -2, -2, -2, -2, -2, -2, -2, -2, -2/
1216 DATA (ISET(I),I=301,500)/200*-2/
1217 DATA ((KFPR(I,J),J=1,2),I=1,50)/
1218 & 23, 0, 24, 0, 25, 0, 24, 0, 25, 0,
1219 & 24, 0, 23, 0, 25, 0, 0, 0, 0, 0,
1220 1 0, 0, 0, 0, 21, 21, 21, 22, 21, 23,
1221 1 21, 24, 21, 25, 22, 22, 22, 23, 22, 24,
1222 2 22, 25, 23, 23, 23, 24, 23, 25, 24, 24,
1223 2 24, 25, 25, 25, 0, 21, 0, 22, 0, 23,
1224 3 0, 24, 0, 25, 0, 21, 0, 22, 0, 23,
1225 3 0, 24, 0, 25, 0, 21, 0, 22, 0, 23,
1226 4 0, 24, 0, 25, 0, 21, 0, 22, 0, 23,
1227 4 0, 24, 0, 25, 0, 21, 0, 22, 0, 23/
1228 DATA ((KFPR(I,J),J=1,2),I=51,100)/
1229 5 0, 24, 0, 25, 0, 0, 0, 0, 0, 0,
1230 5 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1231 6 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1232 6 0, 0, 0, 0, 21, 21, 24, 24, 23, 24,
1233 7 23, 23, 24, 24, 23, 24, 23, 25, 22, 22,
1234 7 23, 23, 24, 24, 24, 25, 25, 25, 0, 211,
1235 8 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1236 8 443, 21,10441, 21,20443, 21, 445, 21, 0, 0,
1237 9 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1238 9 0, 0, 0, 0, 0, 0, 0, 0, 0, 0/
1239 DATA ((KFPR(I,J),J=1,2),I=101,150)/
1240 & 23, 0, 25, 0, 25, 0, 0, 0, 0, 0,
1241 & 443, 22, 443, 21, 443, 22, 0, 0, 22, 25,
1242 1 21, 25, 0, 25, 21, 25, 22, 22, 21, 22,
1243 1 22, 23, 23, 23, 24, 24, 0, 0, 0, 0,
1244 2 25, 6, 25, 6, 25, 0, 25, 0, 0, 0,
1245 2 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1246 3 23, 5, 0, 0, 0, 0, 0, 0, 0, 0,
1247 3 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1248 4 32, 0, 34, 0, 37, 0, 40, 0, 39, 0,
1249 4 0, 0, 4000001, 0, 4000002, 0, 38, 0, 0, 0/
1250 DATA ((KFPR(I,J),J=1,2),I=151,200)/
1251 5 35, 0, 35, 0, 35, 0, 0, 0, 0, 0,
1252 5 36, 0, 36, 0, 36, 0, 0, 0, 0, 0,
1253 6 6, 37, 39, 0, 39, 39, 39, 39, 11, 0,
1254 6 11, 0, 0, 4000001, 0, 4000002, 0, 0, 0, 0,
1255 7 23, 35, 24, 35, 35, 0, 35, 0, 0, 0,
1256 7 23, 36, 24, 36, 36, 0, 36, 0, 0, 0,
1257 8 35, 6, 35, 6, 0, 0, 0, 0, 0, 0,
1258 8 36, 6, 36, 6, 0, 0, 0, 0, 0, 0,
1259 9 54, 0, 55, 0, 56, 0, 11, 0, 0, 0,
1260 9 0, 0, 0, 0, 0, 0, 0, 0, 0, 0/
1261 DATA ((KFPR(I,J),J=1,2),I=201,240)/
1262 & 1000011, 1000011, 2000011, 2000011, 1000011,
1263 & 2000011, 1000013, 1000013, 2000013, 2000013,
1264 & 1000013, 2000013, 1000015, 1000015, 2000015,
1265 & 2000015, 1000015, 2000015, 1000011, 1000012,
1266 1 1000015, 1000016, 2000015, 1000016, 1000012,
1267 1 1000012, 1000016, 1000016, 0, 0,
1268 1 1000022, 1000022, 1000023, 1000023, 1000025,
1269 1 1000025, 1000035, 1000035, 1000022, 1000023,
1270 2 1000022, 1000025, 1000022, 1000035, 1000023,
1271 2 1000025, 1000023, 1000035, 1000025, 1000035,
1272 2 1000024, 1000024, 1000037, 1000037, 1000024,
1273 2 1000037, 1000022, 1000024, 1000023, 1000024,
1274 3 1000025, 1000024, 1000035, 1000024, 1000022,
1275 3 1000037, 1000023, 1000037, 1000025, 1000037,
1276 3 1000035, 1000037, 1000021, 1000022, 1000021,
1277 3 1000023, 1000021, 1000025, 1000021, 1000035/
1278 DATA ((KFPR(I,J),J=1,2),I=241,280)/
1279 4 1000021, 1000024, 1000021, 1000037, 1000021,
1280 4 1000021, 1000021, 1000021, 0, 0,
1281 4 1000002, 1000022, 2000002, 1000022, 1000002,
1282 4 1000023, 2000002, 1000023, 1000002, 1000025,
1283 5 2000002, 1000025, 1000002, 1000035, 2000002,
1284 5 1000035, 1000001, 1000024, 2000005, 1000024,
1285 5 1000001, 1000037, 2000005, 1000037, 1000002,
1286 5 1000021, 2000002, 1000021, 0, 0,
1287 6 1000006, 1000006, 2000006, 2000006, 1000006,
1288 6 2000006, 1000006, 1000006, 2000006, 2000006,
1289 6 0, 0, 0, 0, 0,
1290 6 0, 0, 0, 0, 0,
1291 7 1000002, 1000002, 2000002, 2000002, 1000002,
1292 7 2000002, 1000002, 1000002, 2000002, 2000002,
1293 7 1000002, 2000002, 1000002, 1000002, 2000002,
1294 7 2000002, 1000002, 1000002, 2000002, 2000002/
1295 DATA ((KFPR(I,J),J=1,2),I=281,500)/440*0/
1296 DATA COEF/10000*0D0/
1297 DATA (((ICOL(I,J,K),K=1,2),J=1,4),I=1,40)/
1298 &4,0,3,0,2,0,1,0,3,0,4,0,1,0,2,0,2,0,0,1,4,0,0,3,3,0,0,4,1,0,0,2,
1299 &3,0,0,4,1,4,3,2,4,0,0,3,4,2,1,3,2,0,4,1,4,0,2,3,4,0,3,4,2,0,1,2,
1300 &3,2,1,0,1,4,3,0,4,3,3,0,2,1,1,0,3,2,1,4,1,0,0,2,2,4,3,1,2,0,0,1,
1301 &3,2,1,4,1,4,3,2,4,2,1,3,4,2,1,3,3,4,4,3,1,2,2,1,2,0,3,1,2,0,0,0,
1302 &4,2,1,0,0,0,1,0,3,0,0,3,1,2,0,0,4,0,0,4,0,0,1,2,2,0,0,1,4,4,3,3,
1303 &2,2,1,1,4,4,3,3,3,3,4,4,1,1,2,2,3,2,1,3,1,2,0,0,4,2,1,4,0,0,1,2,
1304 &4,0,0,0,4,0,1,3,0,0,3,0,2,4,3,0,3,4,0,0,1,0,0,1,0,0,3,4,2,0,0,2,
1305 &3,0,0,0,1,0,0,0,0,0,3,0,2,0,0,0,2,0,3,1,2,0,0,0,3,2,1,0,1,0,0,0,
1306 &4,4,3,3,2,2,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
1307 &0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0/
1308
1309C...Treatment of resonances.
1310 DATA (MWID(I) ,I= 1, 500)/5*0,3*1,8*0,1,5*0,3*1,6*0,1,0,7*1,
1311 &10*0,2*1,0,3*1,245*0,19*2,0,7*2,0,2,0,2,0,4*1,163*0/
1312
1313C...Character constants: name of processes.
1314 DATA PROC(0)/ 'All included subprocesses '/
1315 DATA (PROC(I),I=1,20)/
1316 &'f + fbar -> gamma*/Z0 ', 'f + fbar'' -> W+/- ',
1317 &'f + fbar -> h0 ', 'gamma + W+/- -> W+/- ',
1318 &'Z0 + Z0 -> h0 ', 'Z0 + W+/- -> W+/- ',
1319 &' ', 'W+ + W- -> h0 ',
1320 &' ', 'f + f'' -> f + f'' (QFD) ',
1321 1'f + f'' -> f + f'' (QCD) ','f + fbar -> f'' + fbar'' ',
1322 1'f + fbar -> g + g ', 'f + fbar -> g + gamma ',
1323 1'f + fbar -> g + Z0 ', 'f + fbar'' -> g + W+/- ',
1324 1'f + fbar -> g + h0 ', 'f + fbar -> gamma + gamma ',
1325 1'f + fbar -> gamma + Z0 ', 'f + fbar'' -> gamma + W+/- '/
1326 DATA (PROC(I),I=21,40)/
1327 2'f + fbar -> gamma + h0 ', 'f + fbar -> Z0 + Z0 ',
1328 2'f + fbar'' -> Z0 + W+/- ', 'f + fbar -> Z0 + h0 ',
1329 2'f + fbar -> W+ + W- ', 'f + fbar'' -> W+/- + h0 ',
1330 2'f + fbar -> h0 + h0 ', 'f + g -> f + g ',
1331 2'f + g -> f + gamma ', 'f + g -> f + Z0 ',
1332 3'f + g -> f'' + W+/- ', 'f + g -> f + h0 ',
1333 3'f + gamma -> f + g ', 'f + gamma -> f + gamma ',
1334 3'f + gamma -> f + Z0 ', 'f + gamma -> f'' + W+/- ',
1335 3'f + gamma -> f + h0 ', 'f + Z0 -> f + g ',
1336 3'f + Z0 -> f + gamma ', 'f + Z0 -> f + Z0 '/
1337 DATA (PROC(I),I=41,60)/
1338 4'f + Z0 -> f'' + W+/- ', 'f + Z0 -> f + h0 ',
1339 4'f + W+/- -> f'' + g ', 'f + W+/- -> f'' + gamma ',
1340 4'f + W+/- -> f'' + Z0 ', 'f + W+/- -> f'' + W+/- ',
1341 4'f + W+/- -> f'' + h0 ', 'f + h0 -> f + g ',
1342 4'f + h0 -> f + gamma ', 'f + h0 -> f + Z0 ',
1343 5'f + h0 -> f'' + W+/- ', 'f + h0 -> f + h0 ',
1344 5'g + g -> f + fbar ', 'g + gamma -> f + fbar ',
1345 5'g + Z0 -> f + fbar ', 'g + W+/- -> f + fbar'' ',
1346 5'g + h0 -> f + fbar ', 'gamma + gamma -> f + fbar ',
1347 5'gamma + Z0 -> f + fbar ', 'gamma + W+/- -> f + fbar'' '/
1348 DATA (PROC(I),I=61,80)/
1349 6'gamma + h0 -> f + fbar ', 'Z0 + Z0 -> f + fbar ',
1350 6'Z0 + W+/- -> f + fbar'' ', 'Z0 + h0 -> f + fbar ',
1351 6'W+ + W- -> f + fbar ', 'W+/- + h0 -> f + fbar'' ',
1352 6'h0 + h0 -> f + fbar ', 'g + g -> g + g ',
1353 6'gamma + gamma -> W+ + W- ', 'gamma + W+/- -> Z0 + W+/- ',
1354 7'Z0 + Z0 -> Z0 + Z0 ', 'Z0 + Z0 -> W+ + W- ',
1355 7'Z0 + W+/- -> Z0 + W+/- ', 'Z0 + Z0 -> Z0 + h0 ',
1356 7'W+ + W- -> gamma + gamma ', 'W+ + W- -> Z0 + Z0 ',
1357 7'W+/- + W+/- -> W+/- + W+/- ', 'W+/- + h0 -> W+/- + h0 ',
1358 7'h0 + h0 -> h0 + h0 ', 'q + gamma -> q'' + pi+/- '/
1359 DATA (PROC(I),I=81,100)/
1360 8'q + qbar -> Q + Qbar, mass ', 'g + g -> Q + Qbar, massive ',
1361 8'f + q -> f'' + Q, massive ', 'g + gamma -> Q + Qbar, mass ',
1362 8'gamma + gamma -> F + Fbar, m', 'g + g -> J/Psi + g ',
1363 8'g + g -> chi_0c + g ', 'g + g -> chi_1c + g ',
1364 8'g + g -> chi_2c + g ', ' ',
1365 9'Elastic scattering ', 'Single diffractive (XB) ',
1366 9'Single diffractive (AX) ', 'Double diffractive ',
1367 9'Low-pT scattering ', 'Semihard QCD 2 -> 2 ',
1368 9' ', ' ',
1369 9' ', ' '/
1370 DATA (PROC(I),I=101,120)/
1371 &'g + g -> gamma*/Z0 ', 'g + g -> h0 ',
1372 &'gamma + gamma -> h0 ', ' ',
1373 &' ', 'g + g -> J/Psi + gamma ',
1374 &'gamma + g -> J/Psi + g ', 'gamma+gamma -> J/Psi + gamma',
1375 &' ', 'f + fbar -> gamma + h0 ',
1376 1'f + fbar -> g + h0 ', 'q + g -> q + h0 ',
1377 1'g + g -> g + h0 ', 'g + g -> gamma + gamma ',
1378 1'g + g -> g + gamma ', 'g + g -> gamma + Z0 ',
1379 1'g + g -> Z0 + Z0 ', 'g + g -> W+ + W- ',
1380 1' ', ' '/
1381 DATA (PROC(I),I=121,140)/
1382 2'g + g -> Q + Qbar + h0 ', 'q + qbar -> Q + Qbar + h0 ',
1383 2'f + f'' -> f + f'' + h0 ',
1384 2'f + f'' -> f" + f"'' + h0 ',
1385 2' ', ' ',
1386 2' ', ' ',
1387 2' ', ' ',
1388 3'g + g -> Z0 + q + qbar ', ' ',
1389 3' ', ' ',
1390 3' ', ' ',
1391 3' ', ' ',
1392 3' ', ' '/
1393 DATA (PROC(I),I=141,160)/
1394 4'f + fbar -> gamma*/Z0/Z''0 ', 'f + fbar'' -> W''+/- ',
1395 4'f + fbar'' -> H+/- ', 'f + fbar'' -> R ',
1396 4'q + l -> LQ ', ' ',
1397 4'd + g -> d* ', 'u + g -> u* ',
1398 4'g + g -> eta_techni ', ' ',
1399 5'f + fbar -> H0 ', 'g + g -> H0 ',
1400 5'gamma + gamma -> H0 ', ' ',
1401 5' ', 'f + fbar -> A0 ',
1402 5'g + g -> A0 ', 'gamma + gamma -> A0 ',
1403 5' ', ' '/
1404 DATA (PROC(I),I=161,180)/
1405 6'f + g -> f'' + H+/- ', 'q + g -> LQ + lbar ',
1406 6'g + g -> LQ + LQbar ', 'q + qbar -> LQ + LQbar ',
1407 6'f + fbar -> f'' + fbar'' (g/Z)',
1408 6'f +fbar'' -> f" + fbar"'' (W) ',
1409 6'q + q'' -> q" + d* ', 'q + q'' -> q" + u* ',
1410 6' ', ' ',
1411 7'f + fbar -> Z0 + H0 ', 'f + fbar'' -> W+/- + H0 ',
1412 7'f + f'' -> f + f'' + H0 ',
1413 7'f + f'' -> f" + f"'' + H0 ',
1414 7' ', 'f + fbar -> Z0 + A0 ',
1415 7'f + fbar'' -> W+/- + A0 ',
1416 7'f + f'' -> f + f'' + A0 ',
1417 7'f + f'' -> f" + f"'' + A0 ',
1418 7' '/
1419 DATA (PROC(I),I=181,200)/
1420 8'g + g -> Q + Qbar + H0 ', 'q + qbar -> Q + Qbar + H0 ',
1421 8' ', ' ',
1422 8' ', 'g + g -> Q + Qbar + A0 ',
1423 8'q + qbar -> Q + Qbar + A0 ', ' ',
1424 8' ', ' ',
1425 9'f + fbar -> rho_tech0 ', 'f + f'' -> rho_tech+/- ',
1426 9'f + fbar -> omega_tech0 ', 'f+fbar -> f''+fbar'' (technic)',
1427 9' ', ' ',
1428 9' ', ' ',
1429 9' ', ' '/
1430 DATA (PROC(I),I=201,220)/
1431 &'f + fbar -> ~e_L + ~e_Lbar ', 'f + fbar -> ~e_R + ~e_Rbar ',
1432 &'f + fbar -> ~e_R + ~e_Lbar ', 'f + fbar -> ~mu_L + ~mu_Lbar',
1433 &'f + fbar -> ~mu_R + ~mu_Rbar', 'f + fbar -> ~mu_L + ~mu_Rbar',
1434 &'f+fbar -> ~tau_1 + ~tau_1bar', 'f+fbar -> ~tau_2 + ~tau_2bar',
1435 &'f+fbar -> ~tau_1 + ~tau_2bar', 'q + qbar'' -> ~l_L + ~nulbar ',
1436 1'q+qbar''-> ~tau_1 + ~nutaubar', 'q+qbar''-> ~tau_2 + ~nutaubar',
1437 1'f + fbar -> ~nul + ~nulbar ', 'f+fbar -> ~nutau + ~nutaubar',
1438 1' ', 'f + fbar -> ~chi1 + ~chi1 ',
1439 1'f + fbar -> ~chi2 + ~chi2 ', 'f + fbar -> ~chi3 + ~chi3 ',
1440 1'f + fbar -> ~chi4 + ~chi4 ', 'f + fbar -> ~chi1 + ~chi2 '/
1441 DATA (PROC(I),I=221,240)/
1442 2'f + fbar -> ~chi1 + ~chi3 ', 'f + fbar -> ~chi1 + ~chi4 ',
1443 2'f + fbar -> ~chi2 + ~chi3 ', 'f + fbar -> ~chi2 + ~chi4 ',
1444 2'f + fbar -> ~chi3 + ~chi4 ', 'f+fbar -> ~chi+-1 + ~chi-+1 ',
1445 2'f+fbar -> ~chi+-2 + ~chi-+2 ', 'f+fbar -> ~chi+-1 + ~chi-+2 ',
1446 2'q + qbar'' -> ~chi1 + ~chi+-1', 'q + qbar'' -> ~chi2 + ~chi+-1',
1447 3'q + qbar'' -> ~chi3 + ~chi+-1', 'q + qbar'' -> ~chi4 + ~chi+-1',
1448 3'q + qbar'' -> ~chi1 + ~chi+-2', 'q + qbar'' -> ~chi2 + ~chi+-2',
1449 3'q + qbar'' -> ~chi3 + ~chi+-2', 'q + qbar'' -> ~chi4 + ~chi+-2',
1450 3'q + qbar -> ~chi1 + ~g ', 'q + qbar -> ~chi2 + ~g ',
1451 3'q + qbar -> ~chi3 + ~g ', 'q + qbar -> ~chi4 + ~g '/
1452 DATA (PROC(I),I=241,260)/
1453 4'q + qbar'' -> ~chi+-1 + ~g ', 'q + qbar'' -> ~chi+-2 + ~g ',
1454 4'q + qbar -> ~g + ~g ', 'g + g -> ~g + ~g ',
1455 4' ', 'qj + g -> ~qj_L + ~chi1 ',
1456 4'qj + g -> ~qj_R + ~chi1 ', 'qj + g -> ~qj_L + ~chi2 ',
1457 4'qj + g -> ~qj_R + ~chi2 ', 'qj + g -> ~qj_L + ~chi3 ',
1458 5'qj + g -> ~qj_R + ~chi3 ', 'qj + g -> ~qj_L + ~chi4 ',
1459 5'qj + g -> ~qj_R + ~chi4 ', 'qj + g -> ~qk_L + ~chi+-1 ',
1460 5'qj + g -> ~qk_R + ~chi+-1 ', 'qj + g -> ~qk_L + ~chi+-2 ',
1461 5'qj + g -> ~qk_R + ~chi+-2 ', 'qj + g -> ~qj_L + ~g ',
1462 5'qj + g -> ~qj_R + ~g ', ' '/
1463 DATA (PROC(I),I=261,280)/
1464 6'f + fbar -> ~t_1 + ~t_1bar ', 'f + fbar -> ~t_2 + ~t_2bar ',
1465 6'f + fbar -> ~t_1 + ~t_2bar ', 'g + g -> ~t_1 + ~t_1bar ',
1466 6'g + g -> ~t_2 + ~t_2bar ', ' ',
1467 6' ', ' ',
1468 6' ', ' ',
1469 7'qi + qj -> ~qi_L + ~qj_L ', 'qi + qj -> ~qi_R + ~qj_R ',
1470 7'qi + qj -> ~qi_L + ~qj_R ', 'qi+qjbar -> ~qi_L + ~qj_Lbar',
1471 7'qi+qjbar -> ~qi_R + ~qj_Rbar', 'qi+qjbar -> ~qi_L + ~qj_Rbar',
1472 7'f + fbar -> ~qi_L + ~qi_Lbar', 'f + fbar -> ~qi_R + ~qi_Rbar',
1473 7'g + g -> ~qi_L + ~qi_Lbar ', 'g + g -> ~qi_R + ~qi_Rbar '/
1474 DATA (PROC(I),I=281,500)/220*' '/
1475
1476C...Cross sections and slope offsets.
1477 DATA SIGT/294*0D0/
1478
1479C...Supersymmetry switches and parameters.
1480 DATA IMSS/0,
1481 & 0, 0, 0, 1, 0, 0, 0, 1, 0, 0,
1482 1 89*0/
1483 DATA RMSS/0D0,
1484 & 80D0,160D0,500D0,800D0,2D0,250D0,200D0,800D0,700D0,800D0,
1485 1 700D0,500D0,250D0,200D0,800D0,400D0,0D0,0.1D0,850D0,0.041D0,
1486 2 1D0,800D0,1D4,1D4,1D4,0D0,0D0,24D17,2*0D0,
1487 3 69*0D0/
1488
1489C...Data for histogramming routines.
1490 DATA IHIST/1000,20000,55,1/
1491 DATA INDX/1000*0/
1492
1493 END
1494
1495C*********************************************************************
1496
1497*$ CREATE PYTEST.FOR
1498*COPY PYTEST
1499C...PYTEST
1500C...A simple program (disguised as subroutine) to run at installation
1501C...as a check that the program works as intended.
1502
1503 SUBROUTINE PYTEST(MTEST)
1504
1505C...Double precision and integer declarations.
1506 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
1507 INTEGER PYK,PYCHGE,PYCOMP
1508C...Commonblocks.
1509 COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5)
1510 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
1511 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
1512 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
1513 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
1514 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
1515 SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYDAT3/,/PYSUBS/,/PYPARS/
1516C...Local arrays.
1517 DIMENSION PSUM(5),PINI(6),PFIN(6)
1518
1519C...Save defaults for values that are changed.
1520 MSTJ1=MSTJ(1)
1521 MSTJ3=MSTJ(3)
1522 MSTJ11=MSTJ(11)
1523 MSTJ42=MSTJ(42)
1524 MSTJ43=MSTJ(43)
1525 MSTJ44=MSTJ(44)
1526 PARJ17=PARJ(17)
1527 PARJ22=PARJ(22)
1528 PARJ43=PARJ(43)
1529 PARJ54=PARJ(54)
1530 MST101=MSTJ(101)
1531 MST104=MSTJ(104)
1532 MST105=MSTJ(105)
1533 MST107=MSTJ(107)
1534 MST116=MSTJ(116)
1535
1536C...First part: loop over simple events to be generated.
1537 IF(MTEST.GE.1) CALL PYTABU(20)
1538 NERR=0
1539 DO 180 IEV=1,500
1540
1541C...Reset parameter values. Switch on some nonstandard features.
1542 MSTJ(1)=1
1543 MSTJ(3)=0
1544 MSTJ(11)=1
1545 MSTJ(42)=2
1546 MSTJ(43)=4
1547 MSTJ(44)=2
1548 PARJ(17)=0.1D0
1549 PARJ(22)=1.5D0
1550 PARJ(43)=1D0
1551 PARJ(54)=-0.05D0
1552 MSTJ(101)=5
1553 MSTJ(104)=5
1554 MSTJ(105)=0
1555 MSTJ(107)=1
1556 IF(IEV.EQ.301.OR.IEV.EQ.351.OR.IEV.EQ.401) MSTJ(116)=3
1557
1558C...Ten events each for some single jets configurations.
1559 IF(IEV.LE.50) THEN
1560 ITY=(IEV+9)/10
1561 MSTJ(3)=-1
1562 IF(ITY.EQ.3.OR.ITY.EQ.4) MSTJ(11)=2
1563 IF(ITY.EQ.1) CALL PY1ENT(1,1,15D0,0D0,0D0)
1564 IF(ITY.EQ.2) CALL PY1ENT(1,3101,15D0,0D0,0D0)
1565 IF(ITY.EQ.3) CALL PY1ENT(1,-2203,15D0,0D0,0D0)
1566 IF(ITY.EQ.4) CALL PY1ENT(1,-4,30D0,0D0,0D0)
1567 IF(ITY.EQ.5) CALL PY1ENT(1,21,15D0,0D0,0D0)
1568
1569C...Ten events each for some simple jet systems; string fragmentation.
1570 ELSEIF(IEV.LE.130) THEN
1571 ITY=(IEV-41)/10
1572 IF(ITY.EQ.1) CALL PY2ENT(1,1,-1,40D0)
1573 IF(ITY.EQ.2) CALL PY2ENT(1,4,-4,30D0)
1574 IF(ITY.EQ.3) CALL PY2ENT(1,2,2103,100D0)
1575 IF(ITY.EQ.4) CALL PY2ENT(1,21,21,40D0)
1576 IF(ITY.EQ.5) CALL PY3ENT(1,2101,21,-3203,30D0,0.6D0,0.8D0)
1577 IF(ITY.EQ.6) CALL PY3ENT(1,5,21,-5,40D0,0.9D0,0.8D0)
1578 IF(ITY.EQ.7) CALL PY3ENT(1,21,21,21,60D0,0.7D0,0.5D0)
1579 IF(ITY.EQ.8) CALL PY4ENT(1,2,21,21,-2,40D0,
1580 & 0.4D0,0.64D0,0.6D0,0.12D0,0.2D0)
1581
1582C...Seventy events with independent fragmentation and momentum cons.
1583 ELSEIF(IEV.LE.200) THEN
1584 ITY=1+(IEV-131)/16
1585 MSTJ(2)=1+MOD(IEV-131,4)
1586 MSTJ(3)=1+MOD((IEV-131)/4,4)
1587 IF(ITY.EQ.1) CALL PY2ENT(1,4,-5,40D0)
1588 IF(ITY.EQ.2) CALL PY3ENT(1,3,21,-3,40D0,0.9D0,0.4D0)
1589 IF(ITY.EQ.3) CALL PY4ENT(1,2,21,21,-2,40D0,
1590 & 0.4D0,0.64D0,0.6D0,0.12D0,0.2D0)
1591 IF(ITY.GE.4) CALL PY4ENT(1,2,-3,3,-2,40D0,
1592 & 0.4D0,0.64D0,0.6D0,0.12D0,0.2D0)
1593
1594C...A hundred events with random jets (check invariant mass).
1595 ELSEIF(IEV.LE.300) THEN
1596 100 DO 110 J=1,5
1597 PSUM(J)=0D0
1598 110 CONTINUE
1599 NJET=2D0+6D0*PYR(0)
1600 DO 130 I=1,NJET
1601 KFL=21
1602 IF(I.EQ.1) KFL=INT(1D0+4D0*PYR(0))
1603 IF(I.EQ.NJET) KFL=-INT(1D0+4D0*PYR(0))
1604 EJET=5D0+20D0*PYR(0)
1605 THETA=ACOS(2D0*PYR(0)-1D0)
1606 PHI=6.2832D0*PYR(0)
1607 IF(I.LT.NJET) CALL PY1ENT(-I,KFL,EJET,THETA,PHI)
1608 IF(I.EQ.NJET) CALL PY1ENT(I,KFL,EJET,THETA,PHI)
1609 IF(I.EQ.1.OR.I.EQ.NJET) MSTJ(93)=1
1610 IF(I.EQ.1.OR.I.EQ.NJET) PSUM(5)=PSUM(5)+PYMASS(KFL)
1611 DO 120 J=1,4
1612 PSUM(J)=PSUM(J)+P(I,J)
1613 120 CONTINUE
1614 130 CONTINUE
1615 IF(PSUM(4)**2-PSUM(1)**2-PSUM(2)**2-PSUM(3)**2.LT.
1616 & (PSUM(5)+PARJ(32))**2) GOTO 100
1617
1618C...Fifty e+e- continuum events with matrix elements.
1619 ELSEIF(IEV.LE.350) THEN
1620 MSTJ(101)=2
1621 CALL PYEEVT(0,40D0)
1622
1623C...Fifty e+e- continuum event with varying shower options.
1624 ELSEIF(IEV.LE.400) THEN
1625 MSTJ(42)=1+MOD(IEV,2)
1626 MSTJ(43)=1+MOD(IEV/2,4)
1627 MSTJ(44)=MOD(IEV/8,3)
1628 CALL PYEEVT(0,90D0)
1629
1630C...Fifty e+e- continuum events with coherent shower.
1631 ELSEIF(IEV.LE.450) THEN
1632 CALL PYEEVT(0,500D0)
1633
1634C...Fifty Upsilon decays to ggg or gammagg with coherent shower.
1635 ELSE
1636 CALL PYONIA(5,9.46D0)
1637 ENDIF
1638
1639C...Generate event. Find total momentum, energy and charge.
1640 DO 140 J=1,4
1641 PINI(J)=PYP(0,J)
1642 140 CONTINUE
1643 PINI(6)=PYP(0,6)
1644 CALL PYEXEC
1645 DO 150 J=1,4
1646 PFIN(J)=PYP(0,J)
1647 150 CONTINUE
1648 PFIN(6)=PYP(0,6)
1649
1650C...Check conservation of energy, momentum and charge;
1651C...usually exact, but only approximate for single jets.
1652 MERR=0
1653 IF(IEV.LE.50) THEN
1654 IF((PFIN(1)-PINI(1))**2+(PFIN(2)-PINI(2))**2.GE.4D0)
1655 & MERR=MERR+1
1656 EPZREM=PINI(4)+PINI(3)-PFIN(4)-PFIN(3)
1657 IF(EPZREM.LT.0D0.OR.EPZREM.GT.2D0*PARJ(31)) MERR=MERR+1
1658 IF(ABS(PFIN(6)-PINI(6)).GT.2.1D0) MERR=MERR+1
1659 ELSE
1660 DO 160 J=1,4
1661 IF(ABS(PFIN(J)-PINI(J)).GT.0.0001D0*PINI(4)) MERR=MERR+1
1662 160 CONTINUE
1663 IF(ABS(PFIN(6)-PINI(6)).GT.0.1D0) MERR=MERR+1
1664 ENDIF
1665 IF(MERR.NE.0) WRITE(MSTU(11),5000) (PINI(J),J=1,4),PINI(6),
1666 & (PFIN(J),J=1,4),PFIN(6)
1667
1668C...Check that all KF codes are known ones, and that partons/particles
1669C...satisfy energy-momentum-mass relation. Store particle statistics.
1670 DO 170 I=1,N
1671 IF(K(I,1).GT.20) GOTO 170
1672 IF(PYCOMP(K(I,2)).EQ.0) THEN
1673 WRITE(MSTU(11),5100) I
1674 MERR=MERR+1
1675 ENDIF
1676 PD=P(I,4)**2-P(I,1)**2-P(I,2)**2-P(I,3)**2-P(I,5)**2
1677 IF(ABS(PD).GT.MAX(0.1D0,0.001D0*P(I,4)**2).OR.P(I,4).LT.0D0)
1678 & THEN
1679 WRITE(MSTU(11),5200) I
1680 MERR=MERR+1
1681 ENDIF
1682 170 CONTINUE
1683 IF(MTEST.GE.1) CALL PYTABU(21)
1684
1685C...List all erroneous events and some normal ones.
1686 IF(MERR.NE.0.OR.MSTU(24).NE.0.OR.MSTU(28).NE.0) THEN
1687 IF(MERR.GE.1) WRITE(MSTU(11),6400)
1688 CALL PYLIST(2)
1689 ELSEIF(MTEST.GE.1.AND.MOD(IEV-5,100).EQ.0) THEN
1690 CALL PYLIST(1)
1691 ENDIF
1692
1693C...Stop execution if too many errors.
1694 IF(MERR.NE.0) NERR=NERR+1
1695 IF(NERR.GE.10) THEN
1696 WRITE(MSTU(11),6300)
1697 CALL PYLIST(1)
1698 STOP
1699 ENDIF
1700 180 CONTINUE
1701
1702C...Summarize result of run.
1703 IF(MTEST.GE.1) CALL PYTABU(22)
1704
1705C...Reset commonblock variables changed during run.
1706 MSTJ(1)=MSTJ1
1707 MSTJ(3)=MSTJ3
1708 MSTJ(11)=MSTJ11
1709 MSTJ(42)=MSTJ42
1710 MSTJ(43)=MSTJ43
1711 MSTJ(44)=MSTJ44
1712 PARJ(17)=PARJ17
1713 PARJ(22)=PARJ22
1714 PARJ(43)=PARJ43
1715 PARJ(54)=PARJ54
1716 MSTJ(101)=MST101
1717 MSTJ(104)=MST104
1718 MSTJ(105)=MST105
1719 MSTJ(107)=MST107
1720 MSTJ(116)=MST116
1721
1722C...Second part: complete events of various kinds.
1723C...Common initial values. Loop over initiating conditions.
1724 MSTP(122)=MAX(0,MIN(2,MTEST))
1725 MDCY(PYCOMP(111),1)=0
1726 DO 230 IPROC=1,8
1727
1728C...Reset process type, kinematics cuts, and the flags used.
1729 MSEL=0
1730 DO 190 ISUB=1,500
1731 MSUB(ISUB)=0
1732 190 CONTINUE
1733 CKIN(1)=2D0
1734 CKIN(3)=0D0
1735 MSTP(2)=1
1736 MSTP(11)=0
1737 MSTP(33)=0
1738 MSTP(81)=1
1739 MSTP(82)=1
1740 MSTP(111)=1
1741 MSTP(131)=0
1742 MSTP(133)=0
1743 PARP(131)=0.01D0
1744
1745C...Prompt photon production at fixed target.
1746 IF(IPROC.EQ.1) THEN
1747 PZSUM=300D0
1748 PESUM=SQRT(PZSUM**2+PYMASS(211)**2)+PYMASS(2212)
1749 PQSUM=2D0
1750 MSEL=10
1751 CKIN(3)=5D0
1752 CALL PYINIT('FIXT','pi+','p',PZSUM)
1753
1754C...QCD processes at ISR energies.
1755 ELSEIF(IPROC.EQ.2) THEN
1756 PESUM=63D0
1757 PZSUM=0D0
1758 PQSUM=2D0
1759 MSEL=1
1760 CKIN(3)=5D0
1761 CALL PYINIT('CMS','p','p',PESUM)
1762
1763C...W production + multiple interactions at CERN Collider.
1764 ELSEIF(IPROC.EQ.3) THEN
1765 PESUM=630D0
1766 PZSUM=0D0
1767 PQSUM=0D0
1768 MSEL=12
1769 CKIN(1)=20D0
1770 MSTP(82)=4
1771 MSTP(2)=2
1772 MSTP(33)=3
1773 CALL PYINIT('CMS','p','pbar',PESUM)
1774
1775C...W/Z gauge boson pairs + pileup events at the Tevatron.
1776 ELSEIF(IPROC.EQ.4) THEN
1777 PESUM=1800D0
1778 PZSUM=0D0
1779 PQSUM=0D0
1780 MSUB(22)=1
1781 MSUB(23)=1
1782 MSUB(25)=1
1783 CKIN(1)=200D0
1784 MSTP(111)=0
1785 MSTP(131)=1
1786 MSTP(133)=2
1787 PARP(131)=0.04D0
1788 CALL PYINIT('CMS','p','pbar',PESUM)
1789
1790C...Higgs production at LHC.
1791 ELSEIF(IPROC.EQ.5) THEN
1792 PESUM=15400D0
1793 PZSUM=0D0
1794 PQSUM=2D0
1795 MSUB(3)=1
1796 MSUB(102)=1
1797 MSUB(123)=1
1798 MSUB(124)=1
1799 PMAS(25,1)=300D0
1800 CKIN(1)=200D0
1801 MSTP(81)=0
1802 MSTP(111)=0
1803 CALL PYINIT('CMS','p','p',PESUM)
1804
1805C...Z' production at SSC.
1806 ELSEIF(IPROC.EQ.6) THEN
1807 PESUM=40000D0
1808 PZSUM=0D0
1809 PQSUM=2D0
1810 MSEL=21
1811 PMAS(32,1)=600D0
1812 CKIN(1)=400D0
1813 MSTP(81)=0
1814 MSTP(111)=0
1815 CALL PYINIT('CMS','p','p',PESUM)
1816
1817C...W pair production at 1 TeV e+e- collider.
1818 ELSEIF(IPROC.EQ.7) THEN
1819 PESUM=1000D0
1820 PZSUM=0D0
1821 PQSUM=0D0
1822 MSUB(25)=1
1823 MSUB(69)=1
1824 MSTP(11)=1
1825 CALL PYINIT('CMS','e+','e-',PESUM)
1826
1827C...Deep inelastic scattering at a LEP+LHC ep collider.
1828 ELSEIF(IPROC.EQ.8) THEN
1829 P(1,1)=0D0
1830 P(1,2)=0D0
1831 P(1,3)=8000D0
1832 P(2,1)=0D0
1833 P(2,2)=0D0
1834 P(2,3)=-80D0
1835 PESUM=8080D0
1836 PZSUM=7920D0
1837 PQSUM=0D0
1838 MSUB(10)=1
1839 CKIN(3)=50D0
1840 MSTP(111)=0
1841 CALL PYINIT('USER','p','e-',PESUM)
1842 ENDIF
1843
1844C...Generate 20 events of each required type.
1845 DO 220 IEV=1,20
1846 CALL PYEVNT
1847 PESUMM=PESUM
1848 IF(IPROC.EQ.4) PESUMM=MSTI(41)*PESUM
1849
1850C...Check conservation of energy/momentum/flavour.
1851 PINI(1)=0D0
1852 PINI(2)=0D0
1853 PINI(3)=PZSUM
1854 PINI(4)=PESUMM
1855 PINI(6)=PQSUM
1856 DO 200 J=1,4
1857 PFIN(J)=PYP(0,J)
1858 200 CONTINUE
1859 PFIN(6)=PYP(0,6)
1860 MERR=0
1861 DEVE=ABS(PFIN(4)-PINI(4))+ABS(PFIN(3)-PINI(3))
1862 DEVT=ABS(PFIN(1)-PINI(1))+ABS(PFIN(2)-PINI(2))
1863 DEVQ=ABS(PFIN(6)-PINI(6))
1864 IF(DEVE.GT.2D-3*PESUM.OR.DEVT.GT.MAX(0.01D0,1D-4*PESUM).OR.
1865 & DEVQ.GT.0.1D0) MERR=1
1866 IF(MERR.NE.0) WRITE(MSTU(11),5000) (PINI(J),J=1,4),PINI(6),
1867 & (PFIN(J),J=1,4),PFIN(6)
1868
1869C...Check that all KF codes are known ones, and that partons/particles
1870C...satisfy energy-momentum-mass relation.
1871 DO 210 I=1,N
1872 IF(K(I,1).GT.20) GOTO 210
1873 IF(PYCOMP(K(I,2)).EQ.0) THEN
1874 WRITE(MSTU(11),5100) I
1875 MERR=MERR+1
1876 ENDIF
1877 PD=P(I,4)**2-P(I,1)**2-P(I,2)**2-P(I,3)**2-P(I,5)**2*
1878 & SIGN(1D0,P(I,5))
1879 IF(ABS(PD).GT.MAX(0.1D0,0.002D0*P(I,4)**2,0.002D0*P(I,5)**2)
1880 & .OR.(P(I,5).GE.0D0.AND.P(I,4).LT.0D0)) THEN
1881 WRITE(MSTU(11),5200) I
1882 MERR=MERR+1
1883 ENDIF
1884 210 CONTINUE
1885
1886C...Listing of erroneous events, and first event of each type.
1887 IF(MERR.GE.1) NERR=NERR+1
1888 IF(NERR.GE.10) THEN
1889 WRITE(MSTU(11),6300)
1890 CALL PYLIST(1)
1891 STOP
1892 ENDIF
1893 IF(MTEST.GE.1.AND.(MERR.GE.1.OR.IEV.EQ.1)) THEN
1894 IF(MERR.GE.1) WRITE(MSTU(11),6400)
1895 CALL PYLIST(1)
1896 ENDIF
1897 220 CONTINUE
1898
1899C...List statistics for each process type.
1900 IF(MTEST.GE.1) CALL PYSTAT(1)
1901 230 CONTINUE
1902
1903C...Summarize result of run.
1904 IF(NERR.EQ.0) WRITE(MSTU(11),6500)
1905 IF(NERR.GT.0) WRITE(MSTU(11),6600) NERR
1906
1907C...Format statements for output.
1908 5000 FORMAT(/' Momentum, energy and/or charge were not conserved ',
1909 &'in following event'/' sum of',9X,'px',11X,'py',11X,'pz',11X,
1910 &'E',8X,'charge'/' before',2X,4(1X,F12.5),1X,F8.2/' after',3X,
1911 &4(1X,F12.5),1X,F8.2)
1912 5100 FORMAT(/5X,'Entry no.',I4,' in following event not known code')
1913 5200 FORMAT(/5X,'Entry no.',I4,' in following event has faulty ',
1914 &'kinematics')
1915 6300 FORMAT(/5X,'This is the tenth error experienced! Something is ',
1916 &'wrong.'/5X,'Execution will be stopped after listing of event.')
1917 6400 FORMAT(5X,'Faulty event follows:')
1918 6500 FORMAT(//5X,'End result of PYTEST: no errors detected.')
1919 6600 FORMAT(//5X,'End result of PYTEST:',I2,' errors detected.'/
1920 &5X,'This should not have happened!')
1921
1922 RETURN
1923 END
1924
1925C*********************************************************************
1926
1927*$ CREATE PYHEPC.FOR
1928*COPY PYHEPC
1929C...PYHEPC
1930C...Converts PYTHIA event record contents to or from
1931C...the standard event record commonblock.
1932
1933 SUBROUTINE PYHEPC(MCONV)
1934
1935C...Double precision and integer declarations.
1936 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
1937 INTEGER PYK,PYCHGE,PYCOMP
1938C...Commonblocks.
1939 COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5)
1940 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
1941 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
1942 SAVE /PYJETS/,/PYDAT1/,/PYDAT2/
1943C...HEPEVT commonblock.
1944 PARAMETER (NMXHEP=4000)
1945 COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP),
1946 &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP)
1947 DOUBLE PRECISION PHEP,VHEP
1948 SAVE /HEPEVT/
1949
1950C...Conversion from PYTHIA to standard, the easy part.
1951 IF(MCONV.EQ.1) THEN
1952 NEVHEP=0
1953 IF(N.GT.NMXHEP) CALL PYERRM(8,
1954 & '(PYHEPC:) no more space in /HEPEVT/')
1955 NHEP=MIN(N,NMXHEP)
1956 DO 140 I=1,NHEP
1957 ISTHEP(I)=0
1958 IF(K(I,1).GE.1.AND.K(I,1).LE.10) ISTHEP(I)=1
1959 IF(K(I,1).GE.11.AND.K(I,1).LE.20) ISTHEP(I)=2
1960 IF(K(I,1).GE.21.AND.K(I,1).LE.30) ISTHEP(I)=3
1961 IF(K(I,1).GE.31.AND.K(I,1).LE.100) ISTHEP(I)=K(I,1)
1962 IDHEP(I)=K(I,2)
1963 JMOHEP(1,I)=K(I,3)
1964 JMOHEP(2,I)=0
1965 IF(K(I,1).NE.3.AND.K(I,1).NE.13.AND.K(I,1).NE.14) THEN
1966 JDAHEP(1,I)=K(I,4)
1967 JDAHEP(2,I)=K(I,5)
1968 ELSE
1969 JDAHEP(1,I)=0
1970 JDAHEP(2,I)=0
1971 ENDIF
1972 DO 100 J=1,5
1973 PHEP(J,I)=P(I,J)
1974 100 CONTINUE
1975 DO 110 J=1,4
1976 VHEP(J,I)=V(I,J)
1977 110 CONTINUE
1978
1979C...Check if new event (from pileup).
1980 IF(I.EQ.1) THEN
1981 INEW=1
1982 ELSE
1983 IF(K(I,1).EQ.21.AND.K(I-1,1).NE.21) INEW=I
1984 ENDIF
1985
1986C...Fill in missing mother information.
1987 IF(I.GE.INEW+2.AND.K(I,1).EQ.21.AND.K(I,3).EQ.0) THEN
1988 IMO1=I-2
1989 IF(I.GE.INEW+3.AND.K(I-1,1).EQ.21.AND.K(I-1,3).EQ.0)
1990 & IMO1=IMO1-1
1991 JMOHEP(1,I)=IMO1
1992 JMOHEP(2,I)=IMO1+1
1993 ELSEIF(K(I,2).GE.91.AND.K(I,2).LE.93) THEN
1994 I1=K(I,3)-1
1995 120 I1=I1+1
1996 IF(I1.GE.I) CALL PYERRM(8,
1997 & '(PYHEPC:) translation of inconsistent event history')
1998 IF(I1.LT.I.AND.K(I1,1).NE.1.AND.K(I1,1).NE.11) GOTO 120
1999 KC=PYCOMP(K(I1,2))
2000 IF(I1.LT.I.AND.KC.EQ.0) GOTO 120
2001 IF(I1.LT.I.AND.KCHG(KC,2).EQ.0) GOTO 120
2002 JMOHEP(2,I)=I1
2003 ELSEIF(K(I,2).EQ.94) THEN
2004 NJET=2
2005 IF(NHEP.GE.I+3.AND.K(I+3,3).LE.I) NJET=3
2006 IF(NHEP.GE.I+4.AND.K(I+4,3).LE.I) NJET=4
2007 JMOHEP(2,I)=MOD(K(I+NJET,4)/MSTU(5),MSTU(5))
2008 IF(JMOHEP(2,I).EQ.JMOHEP(1,I)) JMOHEP(2,I)=
2009 & MOD(K(I+1,4)/MSTU(5),MSTU(5))
2010 ENDIF
2011
2012C...Fill in missing daughter information.
2013 IF(K(I,2).EQ.94.AND.MSTU(16).NE.2) THEN
2014 DO 130 I1=JDAHEP(1,I),JDAHEP(2,I)
2015 I2=MOD(K(I1,4)/MSTU(5),MSTU(5))
2016 JDAHEP(1,I2)=I
2017 130 CONTINUE
2018 ENDIF
2019 IF(K(I,2).GE.91.AND.K(I,2).LE.94) GOTO 140
2020 I1=JMOHEP(1,I)
2021 IF(I1.LE.0.OR.I1.GT.NHEP) GOTO 140
2022 IF(K(I1,1).NE.13.AND.K(I1,1).NE.14) GOTO 140
2023 IF(JDAHEP(1,I1).EQ.0) THEN
2024 JDAHEP(1,I1)=I
2025 ELSE
2026 JDAHEP(2,I1)=I
2027 ENDIF
2028 140 CONTINUE
2029 DO 150 I=1,NHEP
2030 IF(K(I,1).NE.13.AND.K(I,1).NE.14) GOTO 150
2031 IF(JDAHEP(2,I).EQ.0) JDAHEP(2,I)=JDAHEP(1,I)
2032 150 CONTINUE
2033
2034C...Conversion from standard to PYTHIA, the easy part.
2035 ELSE
2036 IF(NHEP.GT.MSTU(4)) CALL PYERRM(8,
2037 & '(PYHEPC:) no more space in /PYJETS/')
2038 N=MIN(NHEP,MSTU(4))
2039 NKQ=0
2040 KQSUM=0
2041 DO 180 I=1,N
2042 K(I,1)=0
2043 IF(ISTHEP(I).EQ.1) K(I,1)=1
2044 IF(ISTHEP(I).EQ.2) K(I,1)=11
2045 IF(ISTHEP(I).EQ.3) K(I,1)=21
2046 K(I,2)=IDHEP(I)
2047 K(I,3)=JMOHEP(1,I)
2048 K(I,4)=JDAHEP(1,I)
2049 K(I,5)=JDAHEP(2,I)
2050 DO 160 J=1,5
2051 P(I,J)=PHEP(J,I)
2052 160 CONTINUE
2053 DO 170 J=1,4
2054 V(I,J)=VHEP(J,I)
2055 170 CONTINUE
2056 V(I,5)=0D0
2057 IF(ISTHEP(I).EQ.2.AND.PHEP(4,I).GT.PHEP(5,I)) THEN
2058 I1=JDAHEP(1,I)
2059 IF(I1.GT.0.AND.I1.LE.NHEP) V(I,5)=(VHEP(4,I1)-VHEP(4,I))*
2060 & PHEP(5,I)/PHEP(4,I)
2061 ENDIF
2062
2063C...Fill in missing information on colour connection in jet systems.
2064 IF(ISTHEP(I).EQ.1) THEN
2065 KC=PYCOMP(K(I,2))
2066 KQ=0
2067 IF(KC.NE.0) KQ=KCHG(KC,2)*ISIGN(1,K(I,2))
2068 IF(KQ.NE.0) NKQ=NKQ+1
2069 IF(KQ.NE.2) KQSUM=KQSUM+KQ
2070 IF(KQ.NE.0.AND.KQSUM.NE.0) THEN
2071 K(I,1)=2
2072 ELSEIF(KQ.EQ.2.AND.I.LT.N) THEN
2073 IF(K(I+1,2).EQ.21) K(I,1)=2
2074 ENDIF
2075 ENDIF
2076 180 CONTINUE
2077 IF(NKQ.EQ.1.OR.KQSUM.NE.0) CALL PYERRM(8,
2078 & '(PYHEPC:) input parton configuration not colour singlet')
2079 ENDIF
2080
2081 END
2082
2083C*********************************************************************
2084
2085*$ CREATE PYINIT.FOR
2086*COPY PYINIT
2087C...PYINIT
2088C...Initializes the generation procedure; finds maxima of the
2089C...differential cross-sections to be used for weighting.
2090
2091 SUBROUTINE PYINIT(FRAME,BEAM,TARGET,WIN)
2092
2093C...Double precision and integer declarations.
2094 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
2095 INTEGER PYK,PYCHGE,PYCOMP
2096C...Commonblocks.
2097 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
2098 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
2099 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
2100 COMMON/PYDAT4/CHAF(500,2)
2101 CHARACTER CHAF*16
2102 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
2103 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
2104 COMMON/PYINT1/MINT(400),VINT(400)
2105 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
2106 COMMON/PYINT5/NGENPD,NGEN(0:500,3),XSEC(0:500,3)
2107 SAVE /PYDAT1/,/PYDAT2/,/PYDAT3/,/PYDAT4/,/PYSUBS/,/PYPARS/,
2108 &/PYINT1/,/PYINT2/,/PYINT5/
2109C...Local arrays and character variables.
2110 DIMENSION ALAMIN(20),NFIN(20)
2111 CHARACTER*(*) FRAME,BEAM,TARGET
2112 CHARACTER CHFRAM*8,CHBEAM*8,CHTARG*8,CHLH(2)*6
2113
2114C...Interface to PDFLIB.
2115 COMMON/W50512/QCDL4,QCDL5
2116 SAVE /W50512/
2117 DOUBLE PRECISION VALUE(20),QCDL4,QCDL5
2118 CHARACTER*20 PARM(20)
2119 DATA VALUE/20*0D0/,PARM/20*' '/
2120
2121C...Data:Lambda and n_f values for parton distributions; months.
2122 DATA ALAMIN/0.177D0,0.239D0,0.247D0,0.2322D0,0.248D0,0.248D0,
2123 &14*0.2D0/,NFIN/20*4/
2124 DATA CHLH/'lepton','hadron'/
2125
2126C...Reset MINT and VINT arrays. Write headers.
2127 DO 100 J=1,400
2128 MINT(J)=0
2129 VINT(J)=0D0
2130 100 CONTINUE
2131 IF(MSTU(12).GE.1) CALL PYLIST(0)
2132 IF(MSTP(122).GE.1) WRITE(MSTU(11),5100)
2133
2134C...Maximum 4 generations; set maximum number of allowed flavours.
2135 MSTP(1)=MIN(4,MSTP(1))
2136 MSTU(114)=MIN(MSTU(114),2*MSTP(1))
2137 MSTP(58)=MIN(MSTP(58),2*MSTP(1))
2138
2139C...Sum up Cabibbo-Kobayashi-Maskawa factors for each quark/lepton.
2140 DO 120 I=-20,20
2141 VINT(180+I)=0D0
2142 IA=IABS(I)
2143 IF(IA.GE.1.AND.IA.LE.2*MSTP(1)) THEN
2144 DO 110 J=1,MSTP(1)
2145 IB=2*J-1+MOD(IA,2)
2146 IF(IB.GE.6.AND.MSTP(9).EQ.0) GOTO 110
2147 IPM=(5-ISIGN(1,I))/2
2148 IDC=J+MDCY(IA,2)+2
2149 IF(MDME(IDC,1).EQ.1.OR.MDME(IDC,1).EQ.IPM) VINT(180+I)=
2150 & VINT(180+I)+VCKM((IA+1)/2,(IB+1)/2)
2151 110 CONTINUE
2152 ELSEIF(IA.GE.11.AND.IA.LE.10+2*MSTP(1)) THEN
2153 VINT(180+I)=1D0
2154 ENDIF
2155 120 CONTINUE
2156
2157C...Initialize parton distributions: PDFLIB.
2158 IF(MSTP(52).EQ.2) THEN
2159 PARM(1)='NPTYPE'
2160 VALUE(1)=1
2161 PARM(2)='NGROUP'
2162 VALUE(2)=MSTP(51)/1000
2163 PARM(3)='NSET'
2164 VALUE(3)=MOD(MSTP(51),1000)
2165 PARM(4)='TMAS'
2166 VALUE(4)=PMAS(6,1)
2167 CALL PDFSET(PARM,VALUE)
2168 MINT(93)=1000000+MSTP(51)
2169 ENDIF
2170
2171C...Choose Lambda value to use in alpha-strong.
2172 MSTU(111)=MSTP(2)
2173 IF(MSTP(3).GE.2) THEN
2174 ALAM=0.2D0
2175 NF=4
2176 IF(MSTP(52).EQ.1.AND.MSTP(51).GE.1.AND.MSTP(51).LE.10) THEN
2177 ALAM=ALAMIN(MSTP(51))
2178 NF=NFIN(MSTP(51))
2179 ELSEIF(MSTP(52).EQ.2) THEN
2180 ALAM=QCDL4
2181 NF=4
2182 ENDIF
2183 PARP(1)=ALAM
2184 PARP(61)=ALAM
2185 PARP(72)=ALAM
2186 PARU(112)=ALAM
2187 MSTU(112)=NF
2188 IF(MSTP(3).EQ.3) PARJ(81)=ALAM
2189 ENDIF
2190
2191C...Initialize the SUSY generation: couplings, masses,
2192C...decay modes, branching ratios, and so on.
2193 CALL PYMSIN
2194
2195C...Initialize widths and partial widths for resonances.
2196 CALL PYINRE
2197C...Set Z0 mass and width for e+e- routines.
2198 PARJ(123)=PMAS(23,1)
2199 PARJ(124)=PMAS(23,2)
2200
2201C...Identify beam and target particles and frame of process.
2202 CHFRAM=FRAME//' '
2203 CHBEAM=BEAM//' '
2204 CHTARG=TARGET//' '
2205 CALL PYINBM(CHFRAM,CHBEAM,CHTARG,WIN)
2206 IF(MINT(65).EQ.1) GOTO 170
2207
2208C...For gamma-p or gamma-gamma allow many (3 or 6) alternatives.
2209C...For e-gamma allow 2 alternatives.
2210 MINT(121)=1
2211 MINT(123)=MSTP(14)
2212 IF(MSTP(14).EQ.10.AND.(MSEL.EQ.1.OR.MSEL.EQ.2)) THEN
2213 IF((MINT(11).EQ.22.OR.MINT(12).EQ.22).AND.
2214 & (IABS(MINT(11)).GE.28.OR.IABS(MINT(12)).GE.28)) MINT(121)=3
2215 IF(MINT(11).EQ.22.AND.MINT(12).EQ.22) MINT(121)=6
2216 IF((MINT(11).EQ.22.OR.MINT(12).EQ.22).AND.
2217 & (IABS(MINT(11)).EQ.11.OR.IABS(MINT(12)).EQ.11)) MINT(121)=2
2218 ENDIF
2219
2220C...Set up kinematics of process.
2221 CALL PYINKI(0)
2222
2223C...Precalculate flavour selection weights
2224 CALL PYKFIN
2225
2226C...Loop over gamma-p or gamma-gamma alternatives.
2227 DO 160 IGA=1,MINT(121)
2228 MINT(122)=IGA
2229
2230C...Select partonic subprocesses to be included in the simulation.
2231 CALL PYINPR
2232
2233C...Count number of subprocesses on.
2234 MINT(48)=0
2235 DO 130 ISUB=1,500
2236 IF(MINT(50).EQ.0.AND.ISUB.GE.91.AND.ISUB.LE.96.AND.
2237 & MSUB(ISUB).EQ.1) THEN
2238 WRITE(MSTU(11),5200) ISUB,CHLH(MINT(41)),CHLH(MINT(42))
2239 STOP
2240 ELSEIF(MSUB(ISUB).EQ.1.AND.ISET(ISUB).EQ.-1) THEN
2241 WRITE(MSTU(11),5300) ISUB
2242 STOP
2243 ELSEIF(MSUB(ISUB).EQ.1.AND.ISET(ISUB).LE.-2) THEN
2244 WRITE(MSTU(11),5400) ISUB
2245 STOP
2246 ELSEIF(MSUB(ISUB).EQ.1) THEN
2247 MINT(48)=MINT(48)+1
2248 ENDIF
2249 130 CONTINUE
2250 IF(MINT(48).EQ.0) THEN
2251 WRITE(MSTU(11),5500)
2252 STOP
2253 ENDIF
2254 MINT(49)=MINT(48)-MSUB(91)-MSUB(92)-MSUB(93)-MSUB(94)
2255
2256C...Reset variables for cross-section calculation.
2257 DO 150 I=0,500
2258 DO 140 J=1,3
2259 NGEN(I,J)=0
2260 XSEC(I,J)=0D0
2261 140 CONTINUE
2262 150 CONTINUE
2263
2264C...Find parametrized total cross-sections.
2265 CALL PYXTOT
2266
2267C...Maxima of differential cross-sections.
2268 IF(MSTP(121).LE.1) CALL PYMAXI
2269
2270C...Initialize possibility of pileup events.
2271 IF(MINT(121).GT.1) MSTP(131)=0
2272 IF(MSTP(131).NE.0) CALL PYPILE(1)
2273
2274C...Initialize multiple interactions with variable impact parameter.
2275 IF(MINT(50).EQ.1.AND.(MINT(49).NE.0.OR.MSTP(131).NE.0).AND.
2276 & MSTP(82).GE.2) CALL PYMULT(1)
2277
2278C...Save results for gamma-p and gamma-gamma alternatives.
2279 IF(MINT(121).GT.1) CALL PYSAVE(1,IGA)
2280 160 CONTINUE
2281
2282C...Initialization finished.
2283 170 IF(MSTP(122).GE.1) WRITE(MSTU(11),5600)
2284
2285C...Formats for initialization information.
2286 5100 FORMAT('1',18('*'),1X,'PYINIT: initialization of PYTHIA ',
2287 &'routines',1X,17('*'))
2288 5200 FORMAT(1X,'Error: process number ',I3,' not meaningful for ',A6,
2289 &'-',A6,' interactions.'/1X,'Execution stopped!')
2290 5300 FORMAT(1X,'Error: requested subprocess',I4,' not implemented.'/
2291 &1X,'Execution stopped!')
2292 5400 FORMAT(1X,'Error: requested subprocess',I4,' not existing.'/
2293 &1X,'Execution stopped!')
2294 5500 FORMAT(1X,'Error: no subprocess switched on.'/
2295 &1X,'Execution stopped.')
2296 5600 FORMAT(/1X,22('*'),1X,'PYINIT: initialization completed',1X,
2297 &22('*'))
2298
2299 RETURN
2300 END
2301
2302C*********************************************************************
2303
2304*$ CREATE PYEVNT.FOR
2305*COPY PYEVNT
2306C...PYEVNT
2307C...Administers the generation of a high-pT event via calls to
2308C...a number of subroutines.
2309
2310 SUBROUTINE PYEVNT
2311
2312C...Double precision and integer declarations.
2313 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
2314 INTEGER PYK,PYCHGE,PYCOMP
2315C...Commonblocks.
2316 COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5)
2317 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
2318 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
2319 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
2320 COMMON/PYINT1/MINT(400),VINT(400)
2321 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
2322 COMMON/PYINT4/MWID(500),WIDS(500,5)
2323 COMMON/PYINT5/NGENPD,NGEN(0:500,3),XSEC(0:500,3)
2324 COMMON/PYUPPR/NUP,KUP(20,7),NFUP,IFUP(10,2),PUP(20,5),Q2UP(0:10)
2325 SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYPARS/,/PYINT1/,/PYINT2/,
2326 &/PYINT4/,/PYINT5/,/PYUPPR/
2327C...Local array.
2328 DIMENSION VTX(4)
2329
2330C...Initial values for some counters.
2331 N=0
2332 MINT(5)=MINT(5)+1
2333 MINT(7)=0
2334 MINT(8)=0
2335 MINT(83)=0
2336 MINT(84)=MSTP(126)
2337 MSTU(24)=0
2338 MSTU70=0
2339 MSTJ14=MSTJ(14)
2340
2341C...If variable energies: redo incoming kinematics and cross-section.
2342 MSTI(61)=0
2343 IF(MSTP(171).EQ.1) THEN
2344 CALL PYINKI(1)
2345 IF(MSTI(61).EQ.1) THEN
2346 MINT(5)=MINT(5)-1
2347 RETURN
2348 ENDIF
2349 IF(MINT(121).GT.1) CALL PYSAVE(3,1)
2350 CALL PYXTOT
2351 ENDIF
2352
2353C...Loop over number of pileup events; check space left.
2354 IF(MSTP(131).LE.0) THEN
2355 NPILE=1
2356 ELSE
2357 CALL PYPILE(2)
2358 NPILE=MINT(81)
2359 ENDIF
2360 DO 260 IPILE=1,NPILE
2361 IF(MINT(84)+100.GE.MSTU(4)) THEN
2362 CALL PYERRM(11,
2363 & '(PYEVNT:) no more space in PYJETS for pileup events')
2364 IF(MSTU(21).GE.1) GOTO 270
2365 ENDIF
2366 MINT(82)=IPILE
2367
2368C...Generate variables of hard scattering.
2369 MINT(51)=0
2370 MSTI(52)=0
2371 100 CONTINUE
2372 IF(MINT(51).NE.0.OR.MSTU(24).NE.0) MSTI(52)=MSTI(52)+1
2373 MINT(31)=0
2374 MINT(51)=0
2375 MINT(57)=0
2376 CALL PYRAND
2377 IF(MSTI(61).EQ.1) THEN
2378 MINT(5)=MINT(5)-1
2379 RETURN
2380 ENDIF
2381 IF(MINT(51).EQ.2) RETURN
2382 ISUB=MINT(1)
2383 IF(MSTP(111).EQ.-1) GOTO 250
2384
2385 IF(ISUB.LE.90.OR.ISUB.GE.95) THEN
2386C...Hard scattering (including low-pT):
2387C...reconstruct kinematics and colour flow of hard scattering.
2388 110 MINT(51)=0
2389 CALL PYSCAT
2390 IF(MINT(51).EQ.1) GOTO 100
2391 IPU1=MINT(84)+1
2392 IPU2=MINT(84)+2
2393 IF(ISUB.EQ.95) GOTO 130
2394
2395C...Showering of initial state partons (optional).
2396 ALAMSV=PARJ(81)
2397 PARJ(81)=PARP(72)
2398 IF(MSTP(61).GE.1.AND.MINT(47).GE.2) CALL PYSSPA(IPU1,IPU2)
2399 PARJ(81)=ALAMSV
2400 IF(MINT(51).EQ.1) GOTO 100
2401
2402C...Showering of final state partons (optional).
2403 ALAMSV=PARJ(81)
2404 PARJ(81)=PARP(72)
2405 IF(MSTP(71).GE.1.AND.ISET(ISUB).GE.2.AND.ISET(ISUB).LE.10)
2406 & THEN
2407 IPU3=MINT(84)+3
2408 IPU4=MINT(84)+4
2409 IF(ISET(ISUB).EQ.5) IPU4=-3
2410 QMAX=VINT(55)
2411 IF(ISET(ISUB).EQ.2) QMAX=SQRT(PARP(71))*VINT(55)
2412 CALL PYSHOW(IPU3,IPU4,QMAX)
2413 ELSEIF(MSTP(71).GE.1.AND.ISET(ISUB).EQ.11.AND.NFUP.GE.1) THEN
2414 DO 120 IUP=1,NFUP
2415 IPU3=IFUP(IUP,1)+MINT(84)
2416 IPU4=IFUP(IUP,2)+MINT(84)
2417 QMAX=SQRT(MAX(0D0,Q2UP(IUP)))
2418 CALL PYSHOW(IPU3,IPU4,QMAX)
2419 120 CONTINUE
2420 ENDIF
2421 PARJ(81)=ALAMSV
2422
2423C...Decay of final state resonances.
2424 MINT(32)=0
2425 IF(MSTP(41).GE.1.AND.ISET(ISUB).LE.10) CALL PYRESD(0)
2426 IF(MINT(51).EQ.1) GOTO 100
2427 MINT(52)=N
2428
2429C...Multiple interactions.
2430 IF(MSTP(81).GE.1.AND.MINT(50).EQ.1) CALL PYMULT(6)
2431 MINT(53)=N
2432
2433C...Hadron remnants and primordial kT.
2434 130 CALL PYREMN(IPU1,IPU2)
2435 IF(MINT(51).EQ.1.AND.MINT(57).GE.1.AND.MINT(57).LE.5) GOTO 110
2436 IF(MINT(51).EQ.1) GOTO 100
2437
2438 ELSE
2439C...Diffractive and elastic scattering.
2440 CALL PYDIFF
2441 ENDIF
2442
2443C...Check that no odd resonance left undecayed.
2444 IF(MSTP(111).GE.1) THEN
2445 NFIX=N
2446 DO 140 I=MINT(84)+1,NFIX
2447 IF(K(I,1).GE.1.AND.K(I,1).LE.10.AND.K(I,2).NE.21.AND.
2448 & K(I,2).NE.22) THEN
2449 IF(MWID(PYCOMP(K(I,2))).NE.0) THEN
2450 CALL PYRESD(I)
2451 IF(MINT(51).EQ.1) GOTO 100
2452 ENDIF
2453 ENDIF
2454 140 CONTINUE
2455 ENDIF
2456
2457C...Recalculate energies from momenta and masses (if desired).
2458 IF(MSTP(113).GE.1) THEN
2459 DO 150 I=MINT(83)+1,N
2460 IF(K(I,1).GT.0.AND.K(I,1).LE.10) P(I,4)=SQRT(P(I,1)**2+
2461 & P(I,2)**2+P(I,3)**2+P(I,5)**2)
2462 150 CONTINUE
2463 NRECAL=N
2464 ENDIF
2465
2466C...Rearrange partons along strings, check invariant mass cuts.
2467 MSTU(28)=0
2468 IF(MSTP(111).LE.0) MSTJ(14)=-1
2469 CALL PYPREP(MINT(84)+1)
2470 MSTJ(14)=MSTJ14
2471 IF(MSTP(112).EQ.1.AND.MSTU(28).EQ.3) GOTO 100
2472 IF(MSTP(125).EQ.0.OR.MSTP(125).EQ.1) THEN
2473 DO 180 I=MINT(84)+1,N
2474 IF(K(I,2).EQ.94) THEN
2475 DO 170 I1=I+1,MIN(N,I+3)
2476 IF(K(I1,3).EQ.I) THEN
2477 K(I1,3)=MOD(K(I1,4)/MSTU(5),MSTU(5))
2478 IF(K(I1,3).EQ.0) THEN
2479 DO 160 II=MINT(84)+1,I-1
2480 IF(K(II,2).EQ.K(I1,2)) THEN
2481 IF(MOD(K(II,4),MSTU(5)).EQ.I1.OR.
2482 & MOD(K(II,5),MSTU(5)).EQ.I1) K(I1,3)=II
2483 ENDIF
2484 160 CONTINUE
2485 IF(K(I+1,3).EQ.0) K(I+1,3)=K(I,3)
2486 ENDIF
2487 ENDIF
2488 170 CONTINUE
2489 ENDIF
2490 180 CONTINUE
2491 CALL PYEDIT(12)
2492 CALL PYEDIT(14)
2493 IF(MSTP(125).EQ.0) CALL PYEDIT(15)
2494 IF(MSTP(125).EQ.0) MINT(4)=0
2495 DO 200 I=MINT(83)+1,N
2496 IF(K(I,1).EQ.11.AND.K(I,4).EQ.0.AND.K(I,5).EQ.0) THEN
2497 DO 190 I1=I+1,N
2498 IF(K(I1,3).EQ.I.AND.K(I,4).EQ.0) K(I,4)=I1
2499 IF(K(I1,3).EQ.I) K(I,5)=I1
2500 190 CONTINUE
2501 ENDIF
2502 200 CONTINUE
2503 ENDIF
2504
2505C...Introduce separators between sections in PYLIST event listing.
2506 IF(IPILE.EQ.1.AND.MSTP(125).LE.0) THEN
2507 MSTU70=1
2508 MSTU(71)=N
2509 ELSEIF(IPILE.EQ.1) THEN
2510 MSTU70=3
2511 MSTU(71)=2
2512 MSTU(72)=MINT(4)
2513 MSTU(73)=N
2514 ENDIF
2515
2516C...Go back to lab frame (needed for vertices, also in fragmentation).
2517 CALL PYFRAM(1)
2518
2519C...Set nonvanishing production vertex (optional).
2520 IF(MSTP(151).EQ.1) THEN
2521 DO 210 J=1,4
2522 VTX(J)=PARP(150+J)*SQRT(-2D0*LOG(MAX(1D-10,PYR(0))))*
2523 & SIN(PARU(2)*PYR(0))
2524 210 CONTINUE
2525 DO 230 I=MINT(83)+1,N
2526 DO 220 J=1,4
2527 V(I,J)=V(I,J)+VTX(J)
2528 220 CONTINUE
2529 230 CONTINUE
2530 ENDIF
2531
2532C...Perform hadronization (if desired).
2533 IF(MSTP(111).GE.1) THEN
2534 CALL PYEXEC
2535 IF(MSTU(24).NE.0) GOTO 100
2536 ENDIF
2537 IF(MSTP(113).GE.1) THEN
2538 DO 240 I=NRECAL,N
2539 IF(P(I,5).GT.0D0) P(I,4)=SQRT(P(I,1)**2+
2540 & P(I,2)**2+P(I,3)**2+P(I,5)**2)
2541 240 CONTINUE
2542 ENDIF
2543 IF(MSTP(125).EQ.0.OR.MSTP(125).EQ.1) CALL PYEDIT(14)
2544
2545C...Store event information and calculate Monte Carlo estimates of
2546C...subprocess cross-sections.
2547 250 IF(IPILE.EQ.1) CALL PYDOCU
2548
2549C...Set counters for current pileup event and loop to next one.
2550 MSTI(41)=IPILE
2551 IF(IPILE.GE.2.AND.IPILE.LE.10) MSTI(40+IPILE)=ISUB
2552 IF(MSTU70.LT.10) THEN
2553 MSTU70=MSTU70+1
2554 MSTU(70+MSTU70)=N
2555 ENDIF
2556 MINT(83)=N
2557 MINT(84)=N+MSTP(126)
2558 IF(IPILE.LT.NPILE) CALL PYFRAM(2)
2559 260 CONTINUE
2560
2561C...Generic information on pileup events. Reconstruct missing history.
2562 IF(MSTP(131).EQ.1.AND.MSTP(133).GE.1) THEN
2563 PARI(91)=VINT(132)
2564 PARI(92)=VINT(133)
2565 PARI(93)=VINT(134)
2566 IF(MSTP(133).GE.2) PARI(93)=PARI(93)*XSEC(0,3)/VINT(131)
2567 ENDIF
2568 CALL PYEDIT(16)
2569
2570C...Transform to the desired coordinate frame.
2571 270 CALL PYFRAM(MSTP(124))
2572 MSTU(70)=MSTU70
2573 PARU(21)=VINT(1)
2574
2575 RETURN
2576 END
2577
2578C***********************************************************************
2579
2580*$ CREATE PYSTAT.FOR
2581*COPY PYSTAT
2582C...PYSTAT
2583C...Prints out information about cross-sections, decay widths, branching
2584C...ratios, kinematical limits, status codes and parameter values.
2585
2586 SUBROUTINE PYSTAT(MSTAT)
2587
2588C...Double precision and integer declarations.
2589 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
2590 INTEGER PYK,PYCHGE,PYCOMP
2591C...Parameter statement to help give large particle numbers.
2592 PARAMETER (KSUSY1=1000000,KSUSY2=2000000,KEXCIT=4000000)
2593C...Commonblocks.
2594 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
2595 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
2596 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
2597 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
2598 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
2599 COMMON/PYINT1/MINT(400),VINT(400)
2600 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
2601 COMMON/PYINT4/MWID(500),WIDS(500,5)
2602 COMMON/PYINT5/NGENPD,NGEN(0:500,3),XSEC(0:500,3)
2603 COMMON/PYINT6/PROC(0:500)
2604 CHARACTER PROC*28
2605 COMMON/PYMSSM/IMSS(0:99),RMSS(0:99)
2606 SAVE /PYDAT1/,/PYDAT2/,/PYDAT3/,/PYSUBS/,/PYPARS/,/PYINT1/,
2607 &/PYINT2/,/PYINT4/,/PYINT5/,/PYINT6/,/PYMSSM/
2608C...Local arrays, character variables and data.
2609 DIMENSION WDTP(0:200),WDTE(0:200,0:5)
2610 CHARACTER PROGA(6)*28,CHAU*16,CHKF*16,CHD1*16,CHD2*16,CHD3*16,
2611 &CHIN(2)*12,STATE(-1:5)*4,CHKIN(21)*18,DISGA(2)*28
2612 DATA PROGA/
2613 &'VMD/hadron * VMD ','VMD/hadron * direct ',
2614 &'VMD/hadron * anomalous ','direct * direct ',
2615 &'direct * anomalous ','anomalous * anomalous '/
2616 DATA DISGA/'e * VMD','e * anomalous'/
2617 DATA STATE/'----','off ','on ','on/+','on/-','on/1','on/2'/,
2618 &CHKIN/' m_hard (GeV/c^2) ',' p_T_hard (GeV/c) ',
2619 &'m_finite (GeV/c^2)',' y*_subsystem ',' y*_large ',
2620 &' y*_small ',' eta*_large ',' eta*_small ',
2621 &'cos(theta*)_large ','cos(theta*)_small ',' x_1 ',
2622 &' x_2 ',' x_F ',' cos(theta_hard) ',
2623 &'m''_hard (GeV/c^2) ',' tau ',' y* ',
2624 &'cos(theta_hard^-) ','cos(theta_hard^+) ',' x_T^2 ',
2625 &' tau'' '/
2626
2627C...Cross-sections.
2628 IF(MSTAT.LE.1) THEN
2629 IF(MINT(121).GT.1) CALL PYSAVE(5,0)
2630 WRITE(MSTU(11),5000)
2631 WRITE(MSTU(11),5100)
2632 WRITE(MSTU(11),5200) 0,PROC(0),NGEN(0,3),NGEN(0,1),XSEC(0,3)
2633 DO 100 I=1,500
2634 IF(MSUB(I).NE.1) GOTO 100
2635 WRITE(MSTU(11),5200) I,PROC(I),NGEN(I,3),NGEN(I,1),XSEC(I,3)
2636 100 CONTINUE
2637 IF(MINT(121).GT.1) THEN
2638 WRITE(MSTU(11),5300)
2639 DO 110 IGA=1,MINT(121)
2640 CALL PYSAVE(3,IGA)
2641 IF(MINT(121).EQ.2) THEN
2642 WRITE(MSTU(11),5200) IGA,DISGA(IGA),NGEN(0,3),NGEN(0,1),
2643 & XSEC(0,3)
2644 ELSE
2645 WRITE(MSTU(11),5200) IGA,PROGA(IGA),NGEN(0,3),NGEN(0,1),
2646 & XSEC(0,3)
2647 ENDIF
2648 110 CONTINUE
2649 CALL PYSAVE(5,0)
2650 ENDIF
2651 WRITE(MSTU(11),5400) 1D0-DBLE(NGEN(0,3))/
2652 & MAX(1D0,DBLE(NGEN(0,2)))
2653
2654C...Decay widths and branching ratios.
2655 ELSEIF(MSTAT.EQ.2) THEN
2656 WRITE(MSTU(11),5500)
2657 WRITE(MSTU(11),5600)
2658 DO 140 KC=1,500
2659 KF=KCHG(KC,4)
2660 CALL PYNAME(KF,CHKF)
2661 IOFF=0
2662 IF(KC.LE.22) THEN
2663 IF(KC.GT.2*MSTP(1).AND.KC.LE.10) GOTO 140
2664 IF(KC.GT.10+2*MSTP(1).AND.KC.LE.20) GOTO 140
2665 IF(KC.LE.5.OR.(KC.GE.11.AND.KC.LE.16)) IOFF=1
2666 IF(KC.EQ.18.AND.PMAS(18,1).LT.1D0) IOFF=1
2667 IF(KC.EQ.21.OR.KC.EQ.22) IOFF=1
2668 ELSE
2669 IF(MWID(KC).LE.0) GOTO 140
2670 IF(IMSS(1).LE.0.AND.(KF/KSUSY1.EQ.1.OR.
2671 & KF/KSUSY1.EQ.2)) GOTO 140
2672 ENDIF
2673C...Off-shell branchings.
2674 IF(IOFF.EQ.1) THEN
2675 NGP=0
2676 IF(KC.LE.20) NGP=(MOD(KC,10)+1)/2
2677 IF(NGP.LE.MSTP(1)) WRITE(MSTU(11),5700) KF,CHKF(1:10),
2678 & PMAS(KC,1),0D0,0D0,STATE(MDCY(KC,1)),0D0
2679 DO 120 J=1,MDCY(KC,3)
2680 IDC=J+MDCY(KC,2)-1
2681 NGP1=0
2682 IF(IABS(KFDP(IDC,1)).LE.20) NGP1=
2683 & (MOD(IABS(KFDP(IDC,1)),10)+1)/2
2684 NGP2=0
2685 IF(IABS(KFDP(IDC,2)).LE.20) NGP2=
2686 & (MOD(IABS(KFDP(IDC,2)),10)+1)/2
2687 CALL PYNAME(KFDP(IDC,1),CHD1)
2688 CALL PYNAME(KFDP(IDC,2),CHD2)
2689 IF(KFDP(IDC,3).EQ.0) THEN
2690 IF(MDME(IDC,2).EQ.102.AND.NGP1.LE.MSTP(1).AND.
2691 & NGP2.LE.MSTP(1)) WRITE(MSTU(11),5800) IDC,CHD1(1:10),
2692 & CHD2(1:10),0D0,0D0,STATE(MDME(IDC,1)),0D0
2693 ELSE
2694 CALL PYNAME(KFDP(IDC,3),CHD3)
2695 IF(MDME(IDC,2).EQ.102.AND.NGP1.LE.MSTP(1).AND.
2696 & NGP2.LE.MSTP(1)) WRITE(MSTU(11),5900) IDC,CHD1(1:10),
2697 & CHD2(1:10),CHD3(1:10),0D0,0D0,STATE(MDME(IDC,1)),0D0
2698 ENDIF
2699 120 CONTINUE
2700C...On-shell decays.
2701 ELSE
2702 CALL PYWIDT(KF,PMAS(KC,1)**2,WDTP,WDTE)
2703 BRFIN=1D0
2704 IF(WDTE(0,0).LE.0D0) BRFIN=0D0
2705 WRITE(MSTU(11),5700) KF,CHKF(1:10),PMAS(KC,1),WDTP(0),1D0,
2706 & STATE(MDCY(KC,1)),BRFIN
2707 DO 130 J=1,MDCY(KC,3)
2708 IDC=J+MDCY(KC,2)-1
2709 NGP1=0
2710 IF(IABS(KFDP(IDC,1)).LE.20) NGP1=
2711 & (MOD(IABS(KFDP(IDC,1)),10)+1)/2
2712 NGP2=0
2713 IF(IABS(KFDP(IDC,2)).LE.20) NGP2=
2714 & (MOD(IABS(KFDP(IDC,2)),10)+1)/2
2715 BRFIN=0D0
2716 IF(WDTE(0,0).GT.0D0) BRFIN=WDTE(J,0)/WDTE(0,0)
2717 CALL PYNAME(KFDP(IDC,1),CHD1)
2718 CALL PYNAME(KFDP(IDC,2),CHD2)
2719 IF(KFDP(IDC,3).EQ.0) THEN
2720 IF(NGP1.LE.MSTP(1).AND.NGP2.LE.MSTP(1))
2721 & WRITE(MSTU(11),5800) IDC,CHD1(1:10),
2722 & CHD2(1:10),WDTP(J),WDTP(J)/WDTP(0),
2723 & STATE(MDME(IDC,1)),BRFIN
2724 ELSE
2725 CALL PYNAME(KFDP(IDC,3),CHD3)
2726 IF(NGP1.LE.MSTP(1).AND.NGP2.LE.MSTP(1))
2727 & WRITE(MSTU(11),5900) IDC,CHD1(1:10),
2728 & CHD2(1:10),CHD3(1:10),WDTP(J),WDTP(J)/WDTP(0),
2729 & STATE(MDME(IDC,1)),BRFIN
2730 ENDIF
2731 130 CONTINUE
2732 ENDIF
2733 140 CONTINUE
2734 WRITE(MSTU(11),6000)
2735
2736C...Allowed incoming partons/particles at hard interaction.
2737 ELSEIF(MSTAT.EQ.3) THEN
2738 WRITE(MSTU(11),6100)
2739 CALL PYNAME(MINT(11),CHAU)
2740 CHIN(1)=CHAU(1:12)
2741 CALL PYNAME(MINT(12),CHAU)
2742 CHIN(2)=CHAU(1:12)
2743 WRITE(MSTU(11),6200) CHIN(1),CHIN(2)
2744 DO 150 I=-20,22
2745 IF(I.EQ.0) GOTO 150
2746 IA=IABS(I)
2747 IF(IA.GT.MSTP(58).AND.IA.LE.10) GOTO 150
2748 IF(IA.GT.10+2*MSTP(1).AND.IA.LE.20) GOTO 150
2749 CALL PYNAME(I,CHAU)
2750 WRITE(MSTU(11),6300) CHAU,STATE(KFIN(1,I)),CHAU,
2751 & STATE(KFIN(2,I))
2752 150 CONTINUE
2753 WRITE(MSTU(11),6400)
2754
2755C...User-defined limits on kinematical variables.
2756 ELSEIF(MSTAT.EQ.4) THEN
2757 WRITE(MSTU(11),6500)
2758 WRITE(MSTU(11),6600)
2759 SHRMAX=CKIN(2)
2760 IF(SHRMAX.LT.0D0) SHRMAX=VINT(1)
2761 WRITE(MSTU(11),6700) CKIN(1),CHKIN(1),SHRMAX
2762 PTHMIN=MAX(CKIN(3),CKIN(5))
2763 PTHMAX=CKIN(4)
2764 IF(PTHMAX.LT.0D0) PTHMAX=0.5D0*SHRMAX
2765 WRITE(MSTU(11),6800) CKIN(3),PTHMIN,CHKIN(2),PTHMAX
2766 WRITE(MSTU(11),6900) CHKIN(3),CKIN(6)
2767 DO 160 I=4,14
2768 WRITE(MSTU(11),6700) CKIN(2*I-1),CHKIN(I),CKIN(2*I)
2769 160 CONTINUE
2770 SPRMAX=CKIN(32)
2771 IF(SPRMAX.LT.0D0) SPRMAX=VINT(1)
2772 WRITE(MSTU(11),6700) CKIN(31),CHKIN(15),SPRMAX
2773 WRITE(MSTU(11),7000)
2774
2775C...Status codes and parameter values.
2776 ELSEIF(MSTAT.EQ.5) THEN
2777 WRITE(MSTU(11),7100)
2778 WRITE(MSTU(11),7200)
2779 DO 170 I=1,100
2780 WRITE(MSTU(11),7300) I,MSTP(I),PARP(I),100+I,MSTP(100+I),
2781 & PARP(100+I)
2782 170 CONTINUE
2783
2784C...List of all processes implemented in the program.
2785 ELSEIF(MSTAT.EQ.6) THEN
2786 WRITE(MSTU(11),7400)
2787 WRITE(MSTU(11),7500)
2788 DO 180 I=1,500
2789 IF(ISET(I).LT.0) GOTO 180
2790 WRITE(MSTU(11),7600) I,PROC(I),ISET(I),KFPR(I,1),KFPR(I,2)
2791 180 CONTINUE
2792 WRITE(MSTU(11),7700)
2793 ENDIF
2794
2795C...Formats for printouts.
2796 5000 FORMAT('1',9('*'),1X,'PYSTAT: Statistics on Number of ',
2797 &'Events and Cross-sections',1X,9('*'))
2798 5100 FORMAT(/1X,78('=')/1X,'I',34X,'I',28X,'I',12X,'I'/1X,'I',12X,
2799 &'Subprocess',12X,'I',6X,'Number of points',6X,'I',4X,'Sigma',3X,
2800 &'I'/1X,'I',34X,'I',28X,'I',12X,'I'/1X,'I',34('-'),'I',28('-'),
2801 &'I',4X,'(mb)',4X,'I'/1X,'I',34X,'I',28X,'I',12X,'I'/1X,'I',1X,
2802 &'N:o',1X,'Type',25X,'I',4X,'Generated',9X,'Tried',1X,'I',12X,
2803 &'I'/1X,'I',34X,'I',28X,'I',12X,'I'/1X,78('=')/1X,'I',34X,'I',28X,
2804 &'I',12X,'I')
2805 5200 FORMAT(1X,'I',1X,I3,1X,A28,1X,'I',1X,I12,1X,I13,1X,'I',1X,1P,
2806 &D10.3,1X,'I')
2807 5300 FORMAT(1X,'I',34X,'I',28X,'I',12X,'I'/1X,78('=')/
2808 &1X,'I',34X,'I',28X,'I',12X,'I')
2809 5400 FORMAT(1X,'I',34X,'I',28X,'I',12X,'I'/1X,78('=')//
2810 &1X,'********* Fraction of events that fail fragmentation ',
2811 &'cuts =',1X,F8.5,' *********'/)
2812 5500 FORMAT('1',27('*'),1X,'PYSTAT: Decay Widths and Branching ',
2813 &'Ratios',1X,27('*'))
2814 5600 FORMAT(/1X,98('=')/1X,'I',49X,'I',13X,'I',12X,'I',6X,'I',12X,'I'/
2815 &1X,'I',5X,'Mother --> Branching/Decay Channel',8X,'I',1X,
2816 &'Width (GeV)',1X,'I',7X,'B.R.',1X,'I',1X,'Stat',1X,'I',2X,
2817 &'Eff. B.R.',1X,'I'/1X,'I',49X,'I',13X,'I',12X,'I',6X,'I',12X,'I'/
2818 &1X,98('='))
2819 5700 FORMAT(1X,'I',49X,'I',13X,'I',12X,'I',6X,'I',12X,'I'/1X,'I',1X,
2820 &I8,2X,A10,3X,'(m =',F10.3,')',2X,'-->',5X,'I',2X,1P,D10.3,0P,1X,
2821 &'I',1X,1P,D10.3,0P,1X,'I',1X,A4,1X,'I',1X,1P,D10.3,0P,1X,'I')
2822 5800 FORMAT(1X,'I',1X,I8,2X,A10,1X,'+',1X,A10,15X,'I',2X,
2823 &1P,D10.3,0P,1X,'I',1X,1P,D10.3,0P,1X,'I',1X,A4,1X,'I',1X,
2824 &1P,D10.3,0P,1X,'I')
2825 5900 FORMAT(1X,'I',1X,I8,2X,A10,1X,'+',1X,A10,1X,'+',1X,A10,2X,'I',2X,
2826 &1P,D10.3,0P,1X,'I',1X,1P,D10.3,0P,1X,'I',1X,A4,1X,'I',1X,
2827 &1P,D10.3,0P,1X,'I')
2828 6000 FORMAT(1X,'I',49X,'I',13X,'I',12X,'I',6X,'I',12X,'I'/1X,98('='))
2829 6100 FORMAT('1',7('*'),1X,'PYSTAT: Allowed Incoming Partons/',
2830 &'Particles at Hard Interaction',1X,7('*'))
2831 6200 FORMAT(/1X,78('=')/1X,'I',38X,'I',37X,'I'/1X,'I',1X,
2832 &'Beam particle:',1X,A12,10X,'I',1X,'Target particle:',1X,A12,7X,
2833 &'I'/1X,'I',38X,'I',37X,'I'/1X,'I',1X,'Content',6X,'State',19X,
2834 &'I',1X,'Content',6X,'State',18X,'I'/1X,'I',38X,'I',37X,'I'/1X,
2835 &78('=')/1X,'I',38X,'I',37X,'I')
2836 6300 FORMAT(1X,'I',1X,A9,5X,A4,19X,'I',1X,A9,5X,A4,18X,'I')
2837 6400 FORMAT(1X,'I',38X,'I',37X,'I'/1X,78('='))
2838 6500 FORMAT('1',12('*'),1X,'PYSTAT: User-Defined Limits on ',
2839 &'Kinematical Variables',1X,12('*'))
2840 6600 FORMAT(/1X,78('=')/1X,'I',76X,'I')
2841 6700 FORMAT(1X,'I',16X,1P,D10.3,0P,1X,'<',1X,A,1X,'<',1X,1P,D10.3,0P,
2842 &16X,'I')
2843 6800 FORMAT(1X,'I',3X,1P,D10.3,0P,1X,'(',1P,D10.3,0P,')',1X,'<',1X,A,
2844 &1X,'<',1X,1P,D10.3,0P,16X,'I')
2845 6900 FORMAT(1X,'I',29X,A,1X,'=',1X,1P,D10.3,0P,16X,'I')
2846 7000 FORMAT(1X,'I',76X,'I'/1X,78('='))
2847 7100 FORMAT('1',12('*'),1X,'PYSTAT: Summary of Status Codes and ',
2848 &'Parameter Values',1X,12('*'))
2849 7200 FORMAT(/3X,'I',4X,'MSTP(I)',9X,'PARP(I)',20X,'I',4X,'MSTP(I)',9X,
2850 &'PARP(I)'/)
2851 7300 FORMAT(1X,I3,5X,I6,6X,1P,D10.3,0P,18X,I3,5X,I6,6X,1P,D10.3)
2852 7400 FORMAT('1',13('*'),1X,'PYSTAT: List of implemented processes',
2853 &1X,13('*'))
2854 7500 FORMAT(/1X,65('=')/1X,'I',34X,'I',28X,'I'/1X,'I',12X,
2855 &'Subprocess',12X,'I',1X,'ISET',2X,'KFPR(I,1)',2X,'KFPR(I,2)',1X,
2856 &'I'/1X,'I',34X,'I',28X,'I'/1X,65('=')/1X,'I',34X,'I',28X,'I')
2857 7600 FORMAT(1X,'I',1X,I3,1X,A28,1X,'I',1X,I4,1X,I10,1X,I10,1X,'I')
2858 7700 FORMAT(1X,'I',34X,'I',28X,'I'/1X,65('='))
2859
2860 RETURN
2861 END
2862
2863C*********************************************************************
2864
2865*$ CREATE PYINRE.FOR
2866*COPY PYINRE
2867C...PYINRE
2868C...Calculates full and effective widths of gauge bosons, stores
2869C...masses and widths, rescales coefficients to be used for
2870C...resonance production generation.
2871
2872 SUBROUTINE PYINRE
2873
2874C...Double precision and integer declarations.
2875 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
2876 INTEGER PYK,PYCHGE,PYCOMP
2877C...Parameter statement to help give large particle numbers.
2878 PARAMETER (KSUSY1=1000000,KSUSY2=2000000,KEXCIT=4000000)
2879C...Commonblocks.
2880 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
2881 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
2882 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
2883 COMMON/PYDAT4/CHAF(500,2)
2884 CHARACTER CHAF*16
2885 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
2886 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
2887 COMMON/PYINT1/MINT(400),VINT(400)
2888 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
2889 COMMON/PYINT4/MWID(500),WIDS(500,5)
2890 COMMON/PYINT6/PROC(0:500)
2891 CHARACTER PROC*28
2892 COMMON/PYMSSM/IMSS(0:99),RMSS(0:99)
2893 SAVE /PYDAT1/,/PYDAT2/,/PYDAT3/,/PYDAT4/,/PYSUBS/,/PYPARS/,
2894 &/PYINT1/,/PYINT2/,/PYINT4/,/PYINT6/,/PYMSSM/
2895C...Local arrays and data.
2896 DIMENSION WDTP(0:200),WDTE(0:200,0:5),WDTPM(0:200),
2897 &WDTEM(0:200,0:5),KCORD(500),PMORD(500)
2898
2899C...Born level couplings in MSSM Higgs doublet sector.
2900 XW=PARU(102)
2901 XWV=XW
2902 IF(MSTP(8).GE.2) XW=1D0-(PMAS(24,1)/PMAS(23,1))**2
2903 XW1=1D0-XW
2904 IF(MSTP(4).EQ.2) THEN
2905 TANBE=PARU(141)
2906 RATBE=((1D0-TANBE**2)/(1D0+TANBE**2))**2
2907 SQMZ=PMAS(23,1)**2
2908 SQMW=PMAS(24,1)**2
2909 SQMH=PMAS(25,1)**2
2910 SQMA=SQMH*(SQMZ-SQMH)/(SQMZ*RATBE-SQMH)
2911 SQMHP=0.5D0*(SQMA+SQMZ+SQRT((SQMA+SQMZ)**2-4D0*SQMA*SQMZ*RATBE))
2912 SQMHC=SQMA+SQMW
2913 IF(SQMH.GE.SQMZ.OR.MIN(SQMA,SQMHP,SQMHC).LE.0D0) THEN
2914 WRITE(MSTU(11),5000)
2915 STOP
2916 ENDIF
2917 PMAS(35,1)=SQRT(SQMHP)
2918 PMAS(36,1)=SQRT(SQMA)
2919 PMAS(37,1)=SQRT(SQMHC)
2920 ALSU=0.5D0*ATAN(2D0*TANBE*(SQMA+SQMZ)/((1D0-TANBE**2)*
2921 & (SQMA-SQMZ)))
2922 BESU=ATAN(TANBE)
2923 PARU(142)=1D0
2924 PARU(143)=1D0
2925 PARU(161)=-SIN(ALSU)/COS(BESU)
2926 PARU(162)=COS(ALSU)/SIN(BESU)
2927 PARU(163)=PARU(161)
2928 PARU(164)=SIN(BESU-ALSU)
2929 PARU(165)=PARU(164)
2930 PARU(168)=SIN(BESU-ALSU)+0.5D0*COS(2D0*BESU)*SIN(BESU+ALSU)/XW
2931 PARU(171)=COS(ALSU)/COS(BESU)
2932 PARU(172)=SIN(ALSU)/SIN(BESU)
2933 PARU(173)=PARU(171)
2934 PARU(174)=COS(BESU-ALSU)
2935 PARU(175)=PARU(174)
2936 PARU(176)=COS(2D0*ALSU)*COS(BESU+ALSU)-2D0*SIN(2D0*ALSU)*
2937 & SIN(BESU+ALSU)
2938 PARU(177)=COS(2D0*BESU)*COS(BESU+ALSU)
2939 PARU(178)=COS(BESU-ALSU)-0.5D0*COS(2D0*BESU)*COS(BESU+ALSU)/XW
2940 PARU(181)=TANBE
2941 PARU(182)=1D0/TANBE
2942 PARU(183)=PARU(181)
2943 PARU(184)=0D0
2944 PARU(185)=PARU(184)
2945 PARU(186)=COS(BESU-ALSU)
2946 PARU(187)=SIN(BESU-ALSU)
2947 PARU(188)=PARU(186)
2948 PARU(189)=PARU(187)
2949 PARU(190)=0D0
2950 PARU(195)=COS(BESU-ALSU)
2951 ENDIF
2952
2953C...Reset effective widths of gauge bosons.
2954 DO 110 I=1,500
2955 DO 100 J=1,5
2956 WIDS(I,J)=1D0
2957 100 CONTINUE
2958 110 CONTINUE
2959
2960C...Order resonances by increasing mass (except Z0 and W+/-).
2961 NRES=0
2962 DO 140 KC=1,500
2963 KF=KCHG(KC,4)
2964 IF(KF.EQ.0) GOTO 140
2965 IF(MWID(KC).EQ.0) GOTO 140
2966 IF(KC.EQ.7.OR.KC.EQ.8.OR.KC.EQ.17.OR.KC.EQ.18) THEN
2967 IF(MSTP(1).LE.3) GOTO 140
2968 ENDIF
2969 IF(KF/KSUSY1.EQ.1.OR.KF/KSUSY1.EQ.2) THEN
2970 IF(IMSS(1).LE.0) GOTO 140
2971 ENDIF
2972 NRES=NRES+1
2973 PMRES=PMAS(KC,1)
2974 IF(KC.EQ.23.OR.KC.EQ.24) PMRES=0D0
2975 DO 120 I1=NRES-1,1,-1
2976 IF(PMRES.GE.PMORD(I1)) GOTO 130
2977 KCORD(I1+1)=KCORD(I1)
2978 PMORD(I1+1)=PMORD(I1)
2979 120 CONTINUE
2980 130 KCORD(I1+1)=KC
2981 PMORD(I1+1)=PMRES
2982 140 CONTINUE
2983
2984C...Loop over possible resonances.
2985 DO 180 I=1,NRES
2986 KC=KCORD(I)
2987 KF=KCHG(KC,4)
2988
2989C...Check that no fourth generation channels on by mistake.
2990 IF(MSTP(1).LE.3) THEN
2991 DO 150 J=1,MDCY(KC,3)
2992 IDC=J+MDCY(KC,2)-1
2993 KFA1=IABS(KFDP(IDC,1))
2994 KFA2=IABS(KFDP(IDC,2))
2995 IF(KFA1.EQ.7.OR.KFA1.EQ.8.OR.KFA1.EQ.17.OR.KFA1.EQ.18.OR.
2996 & KFA2.EQ.7.OR.KFA2.EQ.8.OR.KFA2.EQ.17.OR.KFA2.EQ.18)
2997 & MDME(IDC,1)=-1
2998 150 CONTINUE
2999 ENDIF
3000
3001C...Check that no supersymmetric channels on by mistake.
3002 IF(IMSS(1).LE.0) THEN
3003 DO 160 J=1,MDCY(KC,3)
3004 IDC=J+MDCY(KC,2)-1
3005 KFA1S=IABS(KFDP(IDC,1))/KSUSY1
3006 KFA2S=IABS(KFDP(IDC,2))/KSUSY1
3007 IF(KFA1S.EQ.1.OR.KFA1S.EQ.2.OR.KFA2S.EQ.1.OR.KFA2S.EQ.2)
3008 & MDME(IDC,1)=-1
3009 160 CONTINUE
3010 ENDIF
3011
3012C...Find mass and evaluate width.
3013 PMR=PMAS(KC,1)
3014 IF(KF.EQ.25.OR.KF.EQ.35.OR.KF.EQ.36) MINT(62)=1
3015 IF(MWID(KC).EQ.3) MINT(63)=1
3016 CALL PYWIDT(KF,PMR**2,WDTP,WDTE)
3017 MINT(51)=0
3018
3019C...Evaluate suppression factors due to non-simulated channels.
3020 IF(KCHG(KC,3).EQ.0) THEN
3021 WIDS(KC,1)=((WDTE(0,1)+WDTE(0,2))**2+
3022 & 2D0*(WDTE(0,1)+WDTE(0,2))*(WDTE(0,4)+WDTE(0,5))+
3023 & 2D0*WDTE(0,4)*WDTE(0,5))/WDTP(0)**2
3024 WIDS(KC,2)=(WDTE(0,1)+WDTE(0,2)+WDTE(0,4))/WDTP(0)
3025 WIDS(KC,3)=0D0
3026 WIDS(KC,4)=0D0
3027 WIDS(KC,5)=0D0
3028 ELSE
3029 IF(MWID(KC).EQ.3) MINT(63)=1
3030 CALL PYWIDT(-KF,PMR**2,WDTPM,WDTEM)
3031 MINT(51)=0
3032 WIDS(KC,1)=((WDTE(0,1)+WDTE(0,2))*(WDTEM(0,1)+WDTEM(0,3))+
3033 & (WDTE(0,1)+WDTE(0,2))*(WDTEM(0,4)+WDTEM(0,5))+
3034 & (WDTE(0,4)+WDTE(0,5))*(WDTEM(0,1)+WDTEM(0,3))+
3035 & WDTE(0,4)*WDTEM(0,5)+WDTE(0,5)*WDTEM(0,4))/WDTP(0)**2
3036 WIDS(KC,2)=(WDTE(0,1)+WDTE(0,2)+WDTE(0,4))/WDTP(0)
3037 WIDS(KC,3)=(WDTEM(0,1)+WDTEM(0,3)+WDTEM(0,4))/WDTP(0)
3038 WIDS(KC,4)=((WDTE(0,1)+WDTE(0,2))**2+
3039 & 2D0*(WDTE(0,1)+WDTE(0,2))*(WDTE(0,4)+WDTE(0,5))+
3040 & 2D0*WDTE(0,4)*WDTE(0,5))/WDTP(0)**2
3041 WIDS(KC,5)=((WDTEM(0,1)+WDTEM(0,3))**2+
3042 & 2D0*(WDTEM(0,1)+WDTEM(0,3))*(WDTEM(0,4)+WDTEM(0,5))+
3043 & 2D0*WDTEM(0,4)*WDTEM(0,5))/WDTP(0)**2
3044 ENDIF
3045
3046C...Set resonance widths and branching ratios;
3047C...also on/off switch for decays.
3048 IF(MWID(KC).EQ.1.OR.MWID(KC).EQ.3) THEN
3049 PMAS(KC,2)=WDTP(0)
3050 PMAS(KC,3)=MIN(0.9D0*PMAS(KC,1),10D0*PMAS(KC,2))
3051 MDCY(KC,1)=MSTP(41)
3052 DO 170 J=1,MDCY(KC,3)
3053 IDC=J+MDCY(KC,2)-1
3054 BRAT(IDC)=0D0
3055 IF(WDTP(0).GT.0D0) BRAT(IDC)=WDTP(J)/WDTP(0)
3056 170 CONTINUE
3057 ENDIF
3058 180 CONTINUE
3059
3060C...Flavours of leptoquark: redefine charge and name.
3061 KFLQQ=KFDP(MDCY(39,2),1)
3062 KFLQL=KFDP(MDCY(39,2),2)
3063 KCHG(39,1)=KCHG(PYCOMP(KFLQQ),1)*ISIGN(1,KFLQQ)+
3064 &KCHG(PYCOMP(KFLQL),1)*ISIGN(1,KFLQL)
3065 LL=1
3066 IF(IABS(KFLQL).EQ.13) LL=2
3067 IF(IABS(KFLQL).EQ.15) LL=3
3068 CHAF(39,1)='LQ_'//CHAF(IABS(KFLQQ),1)(1:1)//
3069 &CHAF(IABS(KFLQL),1)(1:LL)//' '
3070 CHAF(39,2)=CHAF(39,2)(1:4+LL)//'bar '
3071
3072C...Special cases in treatment of gamma*/Z0: redefine process name.
3073 IF(MSTP(43).EQ.1) THEN
3074 PROC(1)='f + fbar -> gamma*'
3075 PROC(15)='f + fbar -> g + gamma*'
3076 PROC(19)='f + fbar -> gamma + gamma*'
3077 PROC(30)='f + g -> f + gamma*'
3078 PROC(35)='f + gamma -> f + gamma*'
3079 ELSEIF(MSTP(43).EQ.2) THEN
3080 PROC(1)='f + fbar -> Z0'
3081 PROC(15)='f + fbar -> g + Z0'
3082 PROC(19)='f + fbar -> gamma + Z0'
3083 PROC(30)='f + g -> f + Z0'
3084 PROC(35)='f + gamma -> f + Z0'
3085 ELSEIF(MSTP(43).EQ.3) THEN
3086 PROC(1)='f + fbar -> gamma*/Z0'
3087 PROC(15)='f + fbar -> g + gamma*/Z0'
3088 PROC(19)='f + fbar -> gamma + gamma*/Z0'
3089 PROC(30)='f + g -> f + gamma*/Z0'
3090 PROC(35)='f + gamma -> f + gamma*/Z0'
3091 ENDIF
3092
3093C...Special cases in treatment of gamma*/Z0/Z'0: redefine process name.
3094 IF(MSTP(44).EQ.1) THEN
3095 PROC(141)='f + fbar -> gamma*'
3096 ELSEIF(MSTP(44).EQ.2) THEN
3097 PROC(141)='f + fbar -> Z0'
3098 ELSEIF(MSTP(44).EQ.3) THEN
3099 PROC(141)='f + fbar -> Z''0'
3100 ELSEIF(MSTP(44).EQ.4) THEN
3101 PROC(141)='f + fbar -> gamma*/Z0'
3102 ELSEIF(MSTP(44).EQ.5) THEN
3103 PROC(141)='f + fbar -> gamma*/Z''0'
3104 ELSEIF(MSTP(44).EQ.6) THEN
3105 PROC(141)='f + fbar -> Z0/Z''0'
3106 ELSEIF(MSTP(44).EQ.7) THEN
3107 PROC(141)='f + fbar -> gamma*/Z0/Z''0'
3108 ENDIF
3109
3110C...Special cases in treatment of WW -> WW: redefine process name.
3111 IF(MSTP(45).EQ.1) THEN
3112 PROC(77)='W+ + W+ -> W+ + W+'
3113 ELSEIF(MSTP(45).EQ.2) THEN
3114 PROC(77)='W+ + W- -> W+ + W-'
3115 ELSEIF(MSTP(45).EQ.3) THEN
3116 PROC(77)='W+/- + W+/- -> W+/- + W+/-'
3117 ENDIF
3118
3119C...Format for error information.
3120 5000 FORMAT(1X,'Error: unphysical input tan^2(beta) and m_H ',
3121 &'combination'/1X,'Execution stopped!')
3122
3123 RETURN
3124 END
3125
3126C*********************************************************************
3127
3128*$ CREATE PYINBM.FOR
3129*COPY PYINBM
3130C...PYINBM
3131C...Identifies the two incoming particles and the choice of frame.
3132
3133 SUBROUTINE PYINBM(CHFRAM,CHBEAM,CHTARG,WIN)
3134
3135C...Double precision and integer declarations.
3136 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
3137 INTEGER PYK,PYCHGE,PYCOMP
3138C...Commonblocks.
3139 COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5)
3140 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
3141 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
3142 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
3143 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
3144 COMMON/PYINT1/MINT(400),VINT(400)
3145 SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYSUBS/,/PYPARS/,/PYINT1/
3146C...Local arrays, character variables and data.
3147 CHARACTER CHFRAM*8,CHBEAM*8,CHTARG*8,CHCOM(3)*8,CHALP(2)*26,
3148 &CHIDNT(3)*8,CHTEMP*8,CHCDE(29)*8,CHINIT*76
3149 DIMENSION LEN(3),KCDE(29),PM(2)
3150 DATA CHALP/'abcdefghijklmnopqrstuvwxyz',
3151 &'ABCDEFGHIJKLMNOPQRSTUVWXYZ'/
3152 DATA CHCDE/'e- ','e+ ','nu_e ','nu_ebar ',
3153 &'mu- ','mu+ ','nu_mu ','nu_mubar','tau- ',
3154 &'tau+ ','nu_tau ','nu_tauba','pi+ ','pi- ',
3155 &'n0 ','nbar0 ','p+ ','pbar- ','gamma ',
3156 &'lambda0 ','sigma- ','sigma0 ','sigma+ ','xi- ',
3157 &'xi0 ','omega- ','pi0 ','reggeon ','pomeron '/
3158 DATA KCDE/11,-11,12,-12,13,-13,14,-14,15,-15,16,-16,
3159 &211,-211,2112,-2112,2212,-2212,22,3122,3112,3212,3222,
3160 &3312,3322,3334,111,28,29/
3161
3162C...Store initial energy. Default frame.
3163 VINT(290)=WIN
3164 MINT(111)=0
3165
3166C...Convert character variables to lowercase and find their length.
3167 CHCOM(1)=CHFRAM
3168 CHCOM(2)=CHBEAM
3169 CHCOM(3)=CHTARG
3170 DO 130 I=1,3
3171 LEN(I)=8
3172 DO 110 LL=8,1,-1
3173 IF(LEN(I).EQ.LL.AND.CHCOM(I)(LL:LL).EQ.' ') LEN(I)=LL-1
3174 DO 100 LA=1,26
3175 IF(CHCOM(I)(LL:LL).EQ.CHALP(2)(LA:LA)) CHCOM(I)(LL:LL)=
3176 & CHALP(1)(LA:LA)
3177 100 CONTINUE
3178 110 CONTINUE
3179 CHIDNT(I)=CHCOM(I)
3180
3181C...Fix up bar, underscore and charge in particle name (if needed).
3182 DO 120 LL=1,6
3183 IF(CHIDNT(I)(LL:LL).EQ.'~') THEN
3184 CHTEMP=CHIDNT(I)
3185 CHIDNT(I)=CHTEMP(1:LL-1)//'bar'//CHTEMP(LL+1:6)//' '
3186 ENDIF
3187 120 CONTINUE
3188 IF(CHIDNT(I)(7:7).EQ.'~') CHIDNT(I)(7:8)='ba'
3189 IF(CHIDNT(I)(1:2).EQ.'nu'.AND.CHIDNT(I)(3:3).NE.'_') THEN
3190 CHTEMP=CHIDNT(I)
3191 CHIDNT(I)='nu_'//CHTEMP(3:7)
3192 ELSEIF(CHIDNT(I)(1:2).EQ.'n ') THEN
3193 CHIDNT(I)(1:3)='n0 '
3194 ELSEIF(CHIDNT(I)(1:4).EQ.'nbar') THEN
3195 CHIDNT(I)(1:5)='nbar0'
3196 ELSEIF(CHIDNT(I)(1:2).EQ.'p ') THEN
3197 CHIDNT(I)(1:3)='p+ '
3198 ELSEIF(CHIDNT(I)(1:4).EQ.'pbar'.OR.
3199 & CHIDNT(I)(1:2).EQ.'p-') THEN
3200 CHIDNT(I)(1:5)='pbar-'
3201 ELSEIF(CHIDNT(I)(1:6).EQ.'lambda') THEN
3202 CHIDNT(I)(7:7)='0'
3203 ELSEIF(CHIDNT(I)(1:3).EQ.'reg') THEN
3204 CHIDNT(I)(1:7)='reggeon'
3205 ELSEIF(CHIDNT(I)(1:3).EQ.'pom') THEN
3206 CHIDNT(I)(1:7)='pomeron'
3207 ENDIF
3208 130 CONTINUE
3209
3210C...Identify free initialization.
3211 IF(CHCOM(1)(1:2).EQ.'no') THEN
3212 MINT(65)=1
3213 RETURN
3214 ENDIF
3215
3216C...Identify incoming beam and target particles.
3217 DO 150 I=1,2
3218 DO 140 J=1,29
3219 IF(CHIDNT(I+1).EQ.CHCDE(J)) MINT(10+I)=KCDE(J)
3220 140 CONTINUE
3221 PM(I)=PYMASS(MINT(10+I))
3222 VINT(2+I)=PM(I)
3223 150 CONTINUE
3224 IF(MINT(11).EQ.0) WRITE(MSTU(11),5000) CHBEAM(1:LEN(2))
3225 IF(MINT(12).EQ.0) WRITE(MSTU(11),5100) CHTARG(1:LEN(3))
3226 IF(MINT(11).EQ.0.OR.MINT(12).EQ.0) STOP
3227
3228C...Identify choice of frame and input energies.
3229 CHINIT=' '
3230
3231C...Events defined in the CM frame.
3232 IF(CHCOM(1)(1:2).EQ.'cm') THEN
3233 MINT(111)=1
3234 S=WIN**2
3235 IF(MSTP(122).GE.1) THEN
3236 IF(CHCOM(2)(1:1).NE.'e') THEN
3237 LOFFS=(31-(LEN(2)+LEN(3)))/2
3238 CHINIT(LOFFS+1:76)='PYTHIA will be initialized for a '//
3239 & CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
3240 & ' collider'//' '
3241 ELSE
3242 LOFFS=(30-(LEN(2)+LEN(3)))/2
3243 CHINIT(LOFFS+1:76)='PYTHIA will be initialized for an '//
3244 & CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
3245 & ' collider'//' '
3246 ENDIF
3247 WRITE(MSTU(11),5200) CHINIT
3248 WRITE(MSTU(11),5300) WIN
3249 ENDIF
3250
3251C...Events defined in fixed target frame.
3252 ELSEIF(CHCOM(1)(1:3).EQ.'fix') THEN
3253 MINT(111)=2
3254 S=PM(1)**2+PM(2)**2+2D0*PM(2)*SQRT(PM(1)**2+WIN**2)
3255 IF(MSTP(122).GE.1) THEN
3256 LOFFS=(29-(LEN(2)+LEN(3)))/2
3257 CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//
3258 & CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
3259 & ' fixed target'//' '
3260 WRITE(MSTU(11),5200) CHINIT
3261 WRITE(MSTU(11),5400) WIN
3262 WRITE(MSTU(11),5500) SQRT(S)
3263 ENDIF
3264
3265C...Frame defined by user three-vectors.
3266 ELSEIF(CHCOM(1)(1:3).EQ.'use') THEN
3267 MINT(111)=3
3268 P(1,5)=PM(1)
3269 P(2,5)=PM(2)
3270 P(1,4)=SQRT(P(1,1)**2+P(1,2)**2+P(1,3)**2+P(1,5)**2)
3271 P(2,4)=SQRT(P(2,1)**2+P(2,2)**2+P(2,3)**2+P(2,5)**2)
3272 S=(P(1,4)+P(2,4))**2-(P(1,1)+P(2,1))**2-(P(1,2)+P(2,2))**2-
3273 & (P(1,3)+P(2,3))**2
3274 IF(MSTP(122).GE.1) THEN
3275 LOFFS=(12-(LEN(2)+LEN(3)))/2
3276 CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//
3277 & CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
3278 & ' user-specified configuration'//' '
3279 WRITE(MSTU(11),5200) CHINIT
3280 WRITE(MSTU(11),5600)
3281 WRITE(MSTU(11),5700) CHCOM(2),P(1,1),P(1,2),P(1,3),P(1,4)
3282 WRITE(MSTU(11),5700) CHCOM(3),P(2,1),P(2,2),P(2,3),P(2,4)
3283 WRITE(MSTU(11),5500) SQRT(MAX(0D0,S))
3284 ENDIF
3285
3286C...Frame defined by user four-vectors.
3287 ELSEIF(CHCOM(1)(1:4).EQ.'four') THEN
3288 MINT(111)=4
3289 PMS1=P(1,4)**2-P(1,1)**2-P(1,2)**2-P(1,3)**2
3290 P(1,5)=SIGN(SQRT(ABS(PMS1)),PMS1)
3291 PMS2=P(2,4)**2-P(2,1)**2-P(2,2)**2-P(2,3)**2
3292 P(2,5)=SIGN(SQRT(ABS(PMS2)),PMS2)
3293 S=(P(1,4)+P(2,4))**2-(P(1,1)+P(2,1))**2-(P(1,2)+P(2,2))**2-
3294 & (P(1,3)+P(2,3))**2
3295 IF(MSTP(122).GE.1) THEN
3296 LOFFS=(12-(LEN(2)+LEN(3)))/2
3297 CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//
3298 & CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
3299 & ' user-specified configuration'//' '
3300 WRITE(MSTU(11),5200) CHINIT
3301 WRITE(MSTU(11),5600)
3302 WRITE(MSTU(11),5700) CHCOM(2),P(1,1),P(1,2),P(1,3),P(1,4)
3303 WRITE(MSTU(11),5700) CHCOM(3),P(2,1),P(2,2),P(2,3),P(2,4)
3304 WRITE(MSTU(11),5500) SQRT(MAX(0D0,S))
3305 ENDIF
3306
3307C...Frame defined by user five-vectors.
3308 ELSEIF(CHCOM(1)(1:4).EQ.'five') THEN
3309 MINT(111)=5
3310 S=(P(1,4)+P(2,4))**2-(P(1,1)+P(2,1))**2-(P(1,2)+P(2,2))**2-
3311 & (P(1,3)+P(2,3))**2
3312 IF(MSTP(122).GE.1) THEN
3313 LOFFS=(12-(LEN(2)+LEN(3)))/2
3314 CHINIT(LOFFS+1:76)='PYTHIA will be initialized for '//
3315 & CHCOM(2)(1:LEN(2))//' on '//CHCOM(3)(1:LEN(3))//
3316 & ' user-specified configuration'//' '
3317 WRITE(MSTU(11),5200) CHINIT
3318 WRITE(MSTU(11),5600)
3319 WRITE(MSTU(11),5700) CHCOM(2),P(1,1),P(1,2),P(1,3),P(1,4)
3320 WRITE(MSTU(11),5700) CHCOM(3),P(2,1),P(2,2),P(2,3),P(2,4)
3321 WRITE(MSTU(11),5500) SQRT(MAX(0D0,S))
3322 ENDIF
3323
3324C...Unknown frame. Error for too low CM energy.
3325 ELSE
3326 WRITE(MSTU(11),5800) CHFRAM(1:LEN(1))
3327 STOP
3328 ENDIF
3329 IF(S.LT.PARP(2)**2) THEN
3330 WRITE(MSTU(11),5900) SQRT(S)
3331 STOP
3332 ENDIF
3333
3334C...Formats for initialization and error information.
3335 5000 FORMAT(1X,'Error: unrecognized beam particle ''',A,'''D0'/
3336 &1X,'Execution stopped!')
3337 5100 FORMAT(1X,'Error: unrecognized target particle ''',A,'''D0'/
3338 &1X,'Execution stopped!')
3339 5200 FORMAT(/1X,78('=')/1X,'I',76X,'I'/1X,'I',A76,'I')
3340 5300 FORMAT(1X,'I',18X,'at',1X,F10.3,1X,'GeV center-of-mass energy',
3341 &19X,'I'/1X,'I',76X,'I'/1X,78('='))
3342 5400 FORMAT(1X,'I',22X,'at',1X,F10.3,1X,'GeV/c lab-momentum',22X,'I')
3343 5500 FORMAT(1X,'I',76X,'I'/1X,'I',11X,'corresponding to',1X,F10.3,1X,
3344 &'GeV center-of-mass energy',12X,'I'/1X,'I',76X,'I'/1X,78('='))
3345 5600 FORMAT(1X,'I',76X,'I'/1X,'I',18X,'px (GeV/c)',3X,'py (GeV/c)',3X,
3346 &'pz (GeV/c)',6X,'E (GeV)',9X,'I')
3347 5700 FORMAT(1X,'I',8X,A8,4(2X,F10.3,1X),8X,'I')
3348 5800 FORMAT(1X,'Error: unrecognized coordinate frame ''',A,'''D0'/
3349 &1X,'Execution stopped!')
3350 5900 FORMAT(1X,'Error: too low CM energy,',F8.3,' GeV for event ',
3351 &'generation.'/1X,'Execution stopped!')
3352
3353 RETURN
3354 END
3355
3356C*********************************************************************
3357
3358*$ CREATE PYINKI.FOR
3359*COPY PYINKI
3360C...PYINKI
3361C...Sets up kinematics, including rotations and boosts to/from CM frame.
3362
3363 SUBROUTINE PYINKI(MODKI)
3364
3365C...Double precision and integer declarations.
3366 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
3367 INTEGER PYK,PYCHGE,PYCOMP
3368C...Commonblocks.
3369 COMMON/PYJETS/N,NPAD,K(4000,5),P(4000,5),V(4000,5)
3370 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
3371 COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
3372 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
3373 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
3374 COMMON/PYINT1/MINT(400),VINT(400)
3375 SAVE /PYJETS/,/PYDAT1/,/PYDAT2/,/PYSUBS/,/PYPARS/,/PYINT1/
3376
3377C...Set initial flavour state.
3378 N=2
3379 DO 100 I=1,2
3380 K(I,1)=1
3381 K(I,2)=MINT(10+I)
3382 100 CONTINUE
3383
3384C...Reset boost. Do kinematics for various cases.
3385 DO 110 J=6,10
3386 VINT(J)=0D0
3387 110 CONTINUE
3388
3389C...Set up kinematics for events defined in CM frame.
3390 IF(MINT(111).EQ.1) THEN
3391 WIN=VINT(290)
3392 IF(MODKI.EQ.1) WIN=PARP(171)*VINT(290)
3393 S=WIN**2
3394 P(1,5)=VINT(3)
3395 P(2,5)=VINT(4)
3396 P(1,1)=0D0
3397 P(1,2)=0D0
3398 P(2,1)=0D0
3399 P(2,2)=0D0
3400 P(1,3)=SQRT(((S-P(1,5)**2-P(2,5)**2)**2-(2D0*P(1,5)*P(2,5))**2)/
3401 & (4D0*S))
3402 P(2,3)=-P(1,3)
3403 P(1,4)=SQRT(P(1,3)**2+P(1,5)**2)
3404 P(2,4)=SQRT(P(2,3)**2+P(2,5)**2)
3405
3406C...Set up kinematics for fixed target events.
3407 ELSEIF(MINT(111).EQ.2) THEN
3408 WIN=VINT(290)
3409 IF(MODKI.EQ.1) WIN=PARP(171)*VINT(290)
3410 P(1,5)=VINT(3)
3411 P(2,5)=VINT(4)
3412 P(1,1)=0D0
3413 P(1,2)=0D0
3414 P(2,1)=0D0
3415 P(2,2)=0D0
3416 P(1,3)=WIN
3417 P(1,4)=SQRT(P(1,3)**2+P(1,5)**2)
3418 P(2,3)=0D0
3419 P(2,4)=P(2,5)
3420 S=P(1,5)**2+P(2,5)**2+2D0*P(2,4)*P(1,4)
3421 VINT(10)=P(1,3)/(P(1,4)+P(2,4))
3422 CALL PYROBO(0,0,0D0,0D0,0D0,0D0,-VINT(10))
3423
3424C...Set up kinematics for events in user-defined frame.
3425 ELSEIF(MINT(111).EQ.3) THEN
3426 P(1,5)=VINT(3)
3427 P(2,5)=VINT(4)
3428 P(1,4)=SQRT(P(1,1)**2+P(1,2)**2+P(1,3)**2+P(1,5)**2)
3429 P(2,4)=SQRT(P(2,1)**2+P(2,2)**2+P(2,3)**2+P(2,5)**2)
3430 DO 120 J=1,3
3431 VINT(7+J)=(P(1,J)+P(2,J))/(P(1,4)+P(2,4))
3432 120 CONTINUE
3433 CALL PYROBO(0,0,0D0,0D0,-VINT(8),-VINT(9),-VINT(10))
3434 VINT(7)=PYANGL(P(1,1),P(1,2))
3435 CALL PYROBO(0,0,0D0,-VINT(7),0D0,0D0,0D0)
3436 VINT(6)=PYANGL(P(1,3),P(1,1))
3437 CALL PYROBO(0,0,-VINT(6),0D0,0D0,0D0,0D0)
3438 S=P(1,5)**2+P(2,5)**2+2D0*(P(1,4)*P(2,4)-P(1,3)*P(2,3))
3439
3440C...Set up kinematics for events with user-defined four-vectors.
3441 ELSEIF(MINT(111).EQ.4) THEN
3442 PMS1=P(1,4)**2-P(1,1)**2-P(1,2)**2-P(1,3)**2
3443 P(1,5)=SIGN(SQRT(ABS(PMS1)),PMS1)
3444 PMS2=P(2,4)**2-P(2,1)**2-P(2,2)**2-P(2,3)**2
3445 P(2,5)=SIGN(SQRT(ABS(PMS2)),PMS2)
3446 DO 130 J=1,3
3447 VINT(7+J)=(P(1,J)+P(2,J))/(P(1,4)+P(2,4))
3448 130 CONTINUE
3449 CALL PYROBO(0,0,0D0,0D0,-VINT(8),-VINT(9),-VINT(10))
3450 VINT(7)=PYANGL(P(1,1),P(1,2))
3451 CALL PYROBO(0,0,0D0,-VINT(7),0D0,0D0,0D0)
3452 VINT(6)=PYANGL(P(1,3),P(1,1))
3453 CALL PYROBO(0,0,-VINT(6),0D0,0D0,0D0,0D0)
3454 S=(P(1,4)+P(2,4))**2
3455
3456C...Set up kinematics for events with user-defined five-vectors.
3457 ELSEIF(MINT(111).EQ.5) THEN
3458 DO 140 J=1,3
3459 VINT(7+J)=(P(1,J)+P(2,J))/(P(1,4)+P(2,4))
3460 140 CONTINUE
3461 CALL PYROBO(0,0,0D0,0D0,-VINT(8),-VINT(9),-VINT(10))
3462 VINT(7)=PYANGL(P(1,1),P(1,2))
3463 CALL PYROBO(0,0,0D0,-VINT(7),0D0,0D0,0D0)
3464 VINT(6)=PYANGL(P(1,3),P(1,1))
3465 CALL PYROBO(0,0,-VINT(6),0D0,0D0,0D0,0D0)
3466 S=(P(1,4)+P(2,4))**2
3467 ENDIF
3468
3469C...Return or error for too low CM energy.
3470 IF(MODKI.EQ.1.AND.S.LT.PARP(2)**2) THEN
3471 IF(MSTP(172).LE.1) THEN
3472 CALL PYERRM(23,
3473 & '(PYINKI:) too low invariant mass in this event')
3474 ELSE
3475 MSTI(61)=1
3476 RETURN
3477 ENDIF
3478 ENDIF
3479
3480C...Save information on incoming particles.
3481 VINT(1)=SQRT(S)
3482 VINT(2)=S
3483 IF(MINT(111).GE.4) VINT(3)=P(1,5)
3484 IF(MINT(111).GE.4) VINT(4)=P(2,5)
3485 VINT(5)=P(1,3)
3486 IF(MODKI.EQ.0) VINT(289)=S
3487 DO 150 J=1,5
3488 V(1,J)=0D0
3489 V(2,J)=0D0
3490 VINT(290+J)=P(1,J)
3491 VINT(295+J)=P(2,J)
3492 150 CONTINUE
3493
3494C...Store pT cut-off and related constants to be used in generation.
3495 IF(MODKI.EQ.0) VINT(285)=CKIN(3)
3496 IF(MSTP(82).LE.1) THEN
3497 IF(MINT(121).GT.1) PARP(81)=1.30D0+0.15D0*LOG(VINT(1)/200D0)/
3498 & LOG(900D0/200D0)
3499 PTMN=PARP(81)
3500 ELSE
3501 IF(MINT(121).GT.1) PARP(82)=1.25D0+0.15D0*LOG(VINT(1)/200D0)/
3502 & LOG(900D0/200D0)
3503 PTMN=PARP(82)
3504 ENDIF
3505 VINT(149)=4D0*PTMN**2/S
3506
3507 RETURN
3508 END
3509
3510C*********************************************************************
3511
3512*$ CREATE PYINPR.FOR
3513*COPY PYINPR
3514C...PYINPR
3515C...Selects partonic subprocesses to be included in the simulation.
3516
3517 SUBROUTINE PYINPR
3518
3519C...Double precision and integer declarations.
3520 IMPLICIT DOUBLE PRECISION(A-H, O-Z)
3521 INTEGER PYK,PYCHGE,PYCOMP
3522C...Commonblocks.
3523 COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
3524 COMMON/PYDAT3/MDCY(500,3),MDME(4000,2),BRAT(4000),KFDP(4000,5)
3525 COMMON/PYSUBS/MSEL,MSELPD,MSUB(500),KFIN(2,-40:40),CKIN(200)
3526 COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
3527 COMMON/PYINT1/MINT(400),VINT(400)
3528 COMMON/PYINT2/ISET(500),KFPR(500,2),COEF(500,20),ICOL(40,4,2)
3529 SAVE /PYDAT1/,/PYDAT3/,/PYSUBS/,/PYPARS/,/PYINT1/,/PYINT2/
3530
3531C...Reset processes to be included.
3532 IF(MSEL.NE.0) THEN
3533 DO 100 I=1,500
3534 MSUB(I)=0
3535 100 CONTINUE
3536 ENDIF
3537
3538C...For e-gamma witn MSTP(14)=10 allow mixture of VMD and anomalous.
3539 IF(MINT(121).EQ.2) THEN
3540 MSUB(10)=1
3541 MINT(123)=MINT(122)+1
3542
3543C...For gamma-p or gamma-gamma with MSTP(14)=10 allow mixture.
3544C...Here also set a few parameters otherwise normally not touched.
3545 ELSEIF(MINT(121).GT.1) THEN
3546
3547C...Parton distributions dampened at small Q2; go to low energies,
3548C...alpha_s <1; no minimum pT cut-off a priori.
3549 MSTP(57)=3
3550 MSTP(85)=0
3551 PARP(2)=2D0
3552 PARU(115)=1D0
3553 CKIN(5)=0.2D0
3554 CKIN(6)=0.2D0
3555
3556C...Define pT cut-off parameters and whether run involves low-pT.
3557 IF(MSTP(82).LE.1) THEN
3558 PTMVMD=1.30D0+0.15D0*LOG(VINT(1)/200D0)/LOG(900D0/200D0)
3559 ELSE
3560 PTMVMD=1.25D0+0.15D0*LOG(VINT(1)/200D0)/LOG(900D0/200D0)
3561 ENDIF
3562 PTMDIR=PARP(15)
3563 PTMANO=PTMVMD
3564 IF(MSTP(15).EQ.5) PTMANO=0.60D0+
3565 & 0.125D0*LOG(1D0+0.10D0*VINT(1))**2
3566 IPTL=1
3567 IF(VINT(285).GT.MAX(PTMVMD,PTMDIR,PTMANO)) IPTL=0
3568 IF(MSEL.EQ.2) IPTL=1
3569
3570C...Set up for p/VMD * VMD.
3571 IF(MINT(122).EQ.1) THEN
3572 MINT(123)=2
3573 MSUB(11)=1
3574 MSUB(12)=1
3575 MSUB(13)=1
3576 MSUB(28)=1
3577 MSUB(53)=1
3578 MSUB(68)=1
3579 IF(IPTL.EQ.1) MSUB(95)=1
3580 IF(MSEL.EQ.2) THEN
3581 MSUB(91)=1
3582 MSUB(92)=1
3583 MSUB(93)=1
3584 MSUB(94)=1
3585 ENDIF
3586 PARP(81)=PTMVMD
3587 PARP(82)=PTMVMD
3588 IF(IPTL.EQ.1) CKIN(3)=0D0
3589
3590C...Set up for p/VMD * direct gamma.
3591 ELSEIF(MINT(122).EQ.2) THEN
3592 MINT(123)=0
3593 IF(MINT(121).EQ.6) MINT(123)=5
3594 MSUB(33)=1
3595 MSUB(54)=1
3596 IF(IPTL.EQ.1) CKIN(3)=PTMDIR
3597
3598C...Set up for p/VMD * anomalous gamma.
3599 ELSEIF(MINT(122).EQ.3) THEN
3600 MINT(123)=3
3601 IF(MINT(121).EQ.6) MINT(123)=7
3602 MSUB(11)=1
3603 MSUB(12)=1
3604 MSUB(13)=1
3605 MSUB(28)=1
3606 MSUB(53)=1
3607 MSUB(68)=1
3608 IF(MSTP(82).GE.2) MSTP(85)=1
3609 IF(IPTL.EQ.1) CKIN(3)=PTMANO
3610
3611C...Set up for direct * direct gamma (switch off leptons).
3612 ELSEIF(MINT(122).EQ.4) THEN
3613 MINT(123)=0
3614 MSUB(58)=1
3615 DO 110 II=MDCY(22,2),MDCY(22,2)+MDCY(22,3)-1
3616 IF(IABS(KFDP(II,1)).GE.10) MDME(II,1)=MIN(0,MDME(II,1))
3617 110 CONTINUE
3618 IF(IPTL.EQ.1) CKIN(3)=PTMDIR
3619
3620C...Set up for direct * anomalous gamma.
3621 ELSEIF(MINT(122).EQ.5) THEN
3622 MINT(123)=6
3623 MSUB(33)=1
3624 MSUB(54)=1
3625 IF(IPTL.EQ.1) CKIN(3)=PTMANO
3626
3627C...Set up for anomalous * anomalous gamma.
3628 ELSEIF(MINT(122).EQ.6) THEN
3629 MINT(123)=3
3630 MSUB(11)=1
3631 MSUB(12)=1
3632 MSUB(13)=1
3633 MSUB(28)=1
3634 MSUB(53)=1
3635 MSUB(68)=1
3636 IF(MSTP(82).GE.2) MSTP(85)=1
3637 IF(IPTL.EQ.1) CKIN(3)=PTMANO
3638 ENDIF
3639
3640C...End of special set up for gamma-p and gamma-gamma.
3641 CKIN