********************************************************************** ********************************************************************** ** ** ** June 1991 ** ** A Manual to ** ** ** ** The Lund Monte Carlo for Jet Fragmentation and e+e- Physics ** ** ** ** JETSET version 7.3 ** ** ** ** Torbjorn Sjostrand ** ** ** ** CERN/TH, CH-1211 Geneva 23 ** ** BITNET/EARN address TORSJO@CERNVM ** ** Tel. +22 - 767 28 20 ** ** ** ** LUSHOW is written together with Mats Bengtsson ** ** ** ** Copyright Torbjorn Sjostrand ** ** ** ********************************************************************** ********************************************************************** * * * Table of Contents * * * * 1. Introductory Material * * 1.1. Program Objective * * 1.2. Update History * * 1.3. Major Changes from JETSET 6.3 to 7.1 * * 1.4. Major Changes from JETSET 7.1 to 7.2 * * 1.5. Major Changes from JETSET 7.2 to 7.3 * * 1.6. Installation of Program * * 1.7. Random Numbers * * 1.8. Programming Philosophy * * 1.9. The First Steps * * * * 2. The Fragmentation/Decay Package * * 2.1. Particle Codes * * 2.2. The Event Record * * 2.3. Definition of Initial Configuration or Commonblock Variables* * 2.4. The Physics Routines * * 2.5. Event Study and Data Listing Routines * * 2.6. Event Analysis Routines * * 2.7. The General Switches and Parameters * * 2.8. Further Parameters and Particle Data * * 2.9. Miscellaneous Comments * * 2.10. Examples * * 2.11. Translation to/from Standard Commonblock * * * * 3. The e+e- Routines * * 3.1. e+e- Continuum Event Generation * * 3.2. A Routine for "onium" Decay * * 3.3. The Commonblock Variables * * 3.4. Examples * * * * Acknowledgements * * References * * * ********************************************************************** * * * Legend: * * >= larger than or equal to * * <= smaller than or equal to * * /= not equal to * * -> goes to (rightarrow) * * ^ what follows next is to be read as an upper index * * _ what follows next is to be read as a lower index * * (D=..) default value for commonblock parameter * * (R) commonblock variable which user may read but not change * * (I) commonblock variable for purely internal use * * * ********************************************************************** 1. Introductory Material This first section contains the objective of this program, background information and installation notes. The experienced user, at a place where the program is already installed, can probably go directly to section 2 without too much of a loss. ______________________________________________________________________ 1.1. Program Objective Jet fragmentation, in its broadest sense, occurs anytime a high energy physics event involves the production of several hadrons. Indeed, most of the experiments, present or planned, at high energy physics laboratories like CERN, DESY, Fermilab, SLAC, Cornell, KEK or Serpukhov involve jet studies, either for their own sake or as a prerequisite for the "real" physics studies. The theory for strong interactions, QCD, in principle should give all the properties of jet fragmentation. The problem is only that nobody knows how to solve QCD, not even for questions orders of magnitude easier than the real-life situation of tens or hundreds of particles produced in one single event, each particle with three momentum degrees of freedom, plus additional flavour and spin quantum numbers. The natural way out has been the introduction of phenomenological models for jet fragmentation, implemented in terms of computer programs that generate complete events, which can be directly compared with experimental data. For subprocesses involving high momentum transfers Q^2, perturbation theory can be used to give the leading order behaviour. The generation of an event can therefore be subdivided into two steps. In the first step, a parton configuration is selected, using perturbative standard model (electroweak and QCD) results. In the second step, these partons are then allowed to fragment into hadrons, with unstable hadrons decaying further. All elements specific to the process under study are supposedly contained in the first step, whereas the fragmentation is assumed to occur according to rules independent of the primary process, "jet universality". Because of quantum mechanical effects, neither of the two steps is of a deterministic nature. A correct treatment would utilize amplitudes and thus contain interference terms; since these terms are unknown anyhow, at least for fragmentation, it is natural to choose a probabilistic approach. This approach readily lends itself to an implementation in terms of Monte Carlo computer programs. Ideally the events generated with these programs should not only give the same mean behaviour as experimentally observed events, but also contain the same degree of event-by-event fluctuations. The first part of the present program contains the process-independent fragmentation and decay routines, whereas the second part deals specifically with the production of the initial parton configuration in e+e- annihilation events. Hard interactions are also covered in the PYTHIA program [Ben87]. A new version (5.4) of this program, fully compatible with JETSET 7.3 conventions, is now available. Indeed, some of the new features in JETSET have been included specifically to simplify this compatibility. PYTHIA has traditionally been a program for hadron-hadron collisions, but the new version also contains a number of hard processes in lepton-lepton collisions. In the future, the e+e- annihilation routines in JETSET will be moved to PYTHIA. Other programs have been built around JETSET 6.3. These include: - the leptoproduction program LEPTO [Ing80]; - the hadron-hadron, hadron-nucleus and nucleus-nucleus program FRITIOF [Nil87]; - the dipole radiation program ARIADNE [Pet88] (now also available for JETSET 7); - the higher twist program TWISTER [Ing87]; - the photoproduction program LUCIFR [Ing87a]; and - the leptoproduction heavy flavour program AROMA [Ing88]. The Lund Monte Carlo is therefore a useful tool for the exploration of jet universality questions, in that experimentalists can directly compare e.g. flavour production parameters determined from different fields of high energy physics. A program like this may be used for different purposes. In one extreme, it may offer a guideline for reasonable "conventional" physics when planning a detector design, or help in estimating acceptance corrections for a given detector geometry. In the other extreme, Monte Carlo results may be compared with existing data in order to extract and study interesting physics. The user probably wants to obtain a "best estimate" without too much ado in the former case, whereas he/she may want to play around with different schemes and possibilities in the latter case, often even studying "crackpot" alternatives just to build up the physical intuition. Indeed, the inclusion of a specific option in the program does not in any way imply a belief that this option is physically sound, only that it may be useful to have around. A prime example is independent fragmentation, a concept we consider false and which is now strongly counterindicated by data. In order to find genuine, nontrivial predictions of our own string fragmentation scheme, independent fragmentation still offers a good reference. A word of warning may therefore be in place. The program description is fairly lengthy, and certainly could not be absorbed in one sitting. This is not even necessary, since all switches and parameters are provided with sensible default values, based on our best understanding of jet fragmentation. A new user can therefore disregard all the fancy options, and just run the program as it is, to gain some experience. Later on, the options that might seem useful can be tried out. No single user is ever likely to find need for more than a fraction of the total number of switches and parameters available, yet many of them have been added to meet special user requests. At all times, however, it should be remembered that this program represents a model, not a theory. Some parameter values have been determined from experimental data, others are "informed guesses". Some parts of the model/program seem to be supported by the data, whereas the status of others is more uncertain. A lot of ground has to be covered, and the models used may therefore not be equally well developed for every aspect. There are many areas that have not even been touched upon in the present program, like spin and polarization phenomena. Furthermore, surprises should always be expected when performing a new experiment, why else do it? ______________________________________________________________________ 1.2. Update History The Lund Monte Carlo is by now a fairly old and well established program, but has still been steadily improved on. This is a continuous process, with the official numbered versions little more than snapshots of this process. For the record, we below list the official versions, with some brief notes. no date publ. main new or improved features 1 Nov78 [Sjo78] single quark jets 2 May79 [Sjo79] heavy flavour jets 3.1 Aug79 --- two-jets in e+e-, preliminary three-jets 3.2 Apr80 [Sjo80] three-jets in e+e- with full matrix elements, toponium -> 3g decays 3.3 Aug80 --- softer fragmentation spectrum 4.1 Apr81 --- baryon production and diquark fragmentation, fourth generation quarks, larger jet systems 4.2 Nov81 --- low-p_T physics 4.3 Mar82 [Sjo82] four-jets and QFD structure in e+e-, Jul82 [Sjo83] event analysis routines 5.1 Apr83 --- improved string fragmentation scheme, symmetric fragmentation, full 2^nd order QCD for e+e- 5.2 Nov83 --- momentum conservation schemes for IF, initial state photon radiation in e+e- 5.3 May84 --- "popcorn" model for baryon production 6.1 Jan85 --- commonblocks restructured, parton showers 6.2 Oct85 [Sjo86] error detection 6.3 Oct86 [Sjo87] new parton shower scheme 7.1 Feb89 --- new particle codes and commonblock structure, more mesons, improved decays, vertex information, Abelian gluon model, Bose-Einstein effects 7.2 Nov89 --- interface to new standard commonblock 7.3 May90 [this] Artru-Mennessier-Bowler type fragm.func. All versions preceding 5.1 should be considered completely obsolete by now, and there is little reason to use anything up to 6.2. From a physics point of view, version 6.3 is very similar to 7.1 - 7.3. For physics studies close to completion, there is therefore little to be gained in switching to version 7. Any future developments will be based on the latter, however (i.e. there will never be a version 6.4). ______________________________________________________________________ 1.3. Major Changes from JETSET 6.3 to 7.1 Version 7.1 of JETSET, the Lund Monte Carlo for jet fragmentation and e+e- physics, represented a major break in continuity from a programming point of view, with no backwards compatibility to JETSET 6.3. The major rewritings were prompted by the adoption of the new particle numbering scheme developed under the aegis of the Particle Data Group [PDG88], a scheme which is intended to become standard in all event generation and detector simulation Monte Carlo programs. This affects commonblock structure and content, as well as calling sequences. In connection with these necessary changes, a number of further improvements were made. The physics content of the program, however, was only modestly changed compared to JETSET 6.3. Here is a more extensive list of major changes in JETSET 7.1. Some of the changes affect the default operation of the program, while others appear as new options. - Some subroutine names and arguments have been changed. - Switches and parameters have been regrouped for better consistency. - New KF particle codes have been introduced and are consistently used in the program. - The storing of particle status and origin in commonblock LUJETS has been reorganized and extended. - Information on decay vertices for long-lived particles has been added. - The lowest orbital angular momentum one meson multiplets (one scalar, two pseudovector and one tensor) have been implemented. - All particle data has been updated in accordance with the 1988 Review of Particle Properties [PDG88]. - Tau and charm decays are now mainly given in terms of individual decay channels, with measured branching ratios, where known. - Top and fourth generation decays have been updated to include W propagator effects and a more flexible treatment of different mass hierarchies. - Parton showers may develop from decay product quarks, e.g. in heavy flavour decays. - Pion and eta Dalitz decays are performed with the proper matrix elements. - A larger a value may be used for diquark production in the Lund flavour dependent symmetric fragmentation function. - Diquark fragmentation (e.g. in a baryon target remnant) is performed in accordance with the "popcorn" baryon pair production scheme. - Bose-Einstein effects may be simulated according to a simple algorithm. - Matrix elements for two toy models, a scalar gluon and an Abelian vector gluon, have been included for the testing of QCD. - The Abelian vector gluon model is also available as an option in the shower evolution. - Nonisotropic azimuthal angle in shower evolution due to gluon polarization. - The alpha-strong function used for matrix elements has been made continuous at flavour thresholds by using a Lambda value which depends on the number of flavours assumed. - The LUCLUS cluster search routine can be used also with a mass distance measure. A new routine LUJMAS for heavy and light jet mass reconstruction has been added. - A routine LUTABU has been added to provide various event analysis facilities: particle content, factorial moments, energy-energy correlation, etc. - A random number generator package, based on the new algorithm by Marsaglia, Zaman and Tsang [Mar90,Jam90], comes with the program. - The old J- and I-quark fragmentation scheme for baryon jets has been removed. As a consequence, so has the one-string low-p_T collision model. - Also a number of other obsolete or rarely used features have been removed, such as the simplified Montvay scheme treatment. As may be understood from the list above, the program has increased significantly in size. Offset against this, the last two points represent an effort to remove obsolete material, in particular when this may complicate the introduction of more interesting options. Those who wish to use one of the deleted features are heartily recommended to run JETSET 6.3 instead. The question of diquark fragmentation has not received its final solution yet, what is presented here is just a (reasonable) makeshift arrangement until further studies have been made. ______________________________________________________________________ 1.4. Major Changes from JETSET 7.1 to 7.2 The updates from version 7.1 to 7.2 are all minor, and just about any program that ran with version 7.1 will also work with JETSET 7.2. Some minor modifications have been made to bring documentation in closer agreement with the standard developed by the QCD event generators subgroup of the 1989 LEP physics workshop [QEG89], and also to allow for some more switches where they logically belong. Specifically, these are the changes where compatibility problems might exist with programs written to run with 7.1: - The meaning of KF codes 91 - 94 has been brought in agreement with the agreed standard [QEG89]. To make space, the internal diquark KC codes were shifted from 91 and 92 to 90; this latter should not be visible externally. - The codes KF = 91 - 93 now appear as part of the event history (this can partly be circumvented, see MSTU(16)), i.e. primary hadrons now point back to an entry with KF = 92 for a fragmenting string, KF = 93 for an independent fragmentation jet system, and to KF = 91 for a small parton system that collapsed into one or two primary hadrons. For the designation of the CM frame of a shower, KF = 94 works as KF = 93 used to, except that the K(I,1) code has been changed from 14 to 11. - The (rarely used) MSTJ(111) and MSTJ(112) variables have been moved to MSTJ(115) and MSTJ(116), respectively, to make room for some further matrix element switches. - In qqbarq'qbar' four-jet events generated in the matrix element option, the secondary q' may now also be b, whereas previously only d, u, s and c were allowed. The parametrized qqbarq'qbar' rate has as a consequence been increased by a factor 5/4 (this neglects the presence of a term for q = q' which does not grow with the number of flavours, but this term is negligible anyway). - For matrix elements, it is no longer checked that the y value is sufficiently large that no negative differential three-jet rates may appear. Instead, the rate is assumed zero in such regions, while it is rescaled uniformly in regions of positive cross-section so as to keep the total three-jet rate "correct". (This change of policy was required to make the implementation of optimized matrix elements practicable.) - The Z, W and t default masses have been updated, together with a few other minor changes of particle data. - An error in the second order corrections to the three-jet rate, introduced by mistake when going from 6.3 to 7.1, has been removed. - A few other minor errors have also been corrected; none of these should really make a difference, however. In addition, a number of extensions have been made, where no backwards incompatibilities are involved. - The Zhu parametrization of second order corrections to the three-jet rate has been introduced as an alternative to the old GKS one. - It is possible to pick an "optimized" Q^2 scale for second order matrix elements (and the R value). - The parton shower can also include emission of on-mass-shell photons. If the photons are emitted in the first step of the shower, a correction is made to the matrix element probability [Gro81]. - It is possible to include nonisotropic azimuthal distributions in shower branchings due to soft gluon interference, in addition to the gluon polarization effects made available in 7.1. - The parton shower can be run in a scalar gluon mode, in addition to the already existing QCD and Abelian vector gluon ones. - An interface LUHEPC exists for conversion between the internal event record commonblock LUJETS and the new standard HEPEVT one [QEG89], see section 2.11. - A routine LUJOIN has been added to help users rapidly to set up colour flow information in user-defined jet systems, as required e.g. for subsequent shower evolution. - It is possible to insert separators between sections in LULIST for improved readability. Also, the LULIST event listing gives information on string systems. - A new option LUEDIT(5) has been introduced. - Calls to ALOG have been replaced with LOG, i.e. generic function names are used throughout. ______________________________________________________________________ 1.5. Major Changes from JETSET 7.2 to 7.3 Changes from version 7.2 to 7.3 are essentially all minor, and should not cause any problems. In the following particulars backwards compatibility is not fully preserved. - The parameter PARU(121), tan^2(beta) for the charged Higgs scenario, has been replaced by PARU(141), which is tan(beta); also the default value has been changed from 1 to 5. - The parameters PARJ(35) and PARJ(40) are no longer used for the 'hybrid' Lund symmetric + SLAC fragmentation function option; instead PARJ(33) and PARJ(38) are used, as in the standard Lund. - Some new decay channels have been introduced and a few moved; thus the old decay channel numbers have changed. - As of March 3, 1991, additional decay channels were introduced in the Higgs sector, with forther decay channel changes as a consequence. - Simultaneously the name of particle 36 was changed from H"0 to A0. - A few minor errors have been corrected. The following new features may be found. - The possibility to pick a softer fragmentation function for heavy flavours (charm and beyond) according to the string inspired shape of Bowler (a la Artru-Mennessier, with further extension by Morris). - A routine for a running alpha_em. - Inclusion of a leptoquark LQ as particle 39. - Extended parametrization of non-standard model couplings (for Z', W', H0, H'0, A0 H+, LQ). - Also a lepton may shower photons in LUSHOW. A number of further additions should appear in later versions of the program. Below is given a few examples of what this could include. - Further updates of charm and bottom decays. - Anisotropic decays of polarized tau leptons. - Fragmentation of strings with a junction, as obtained when the two quarks of an effective diquark have a significant relative transverse momentum. - A retuning of fragmentation parameters, with contributions from the recently included meson multiplets accounted for. - Switch to allow top and heavier quarks to decay before they hadronize. No guarantees are given, however, neither that these will be included, nor when it would happen. The present text is intended to contain a complete manual for the JETSET 7.3 user. It is not a description of the physics of the program, however. A new user should therefore have a look in the following literature, at the very least. - [And83] : the standard introduction to the Lund model, unfortunately with a lot of details now out of date. - [Sjo86] : the latest complete description of the physics that goes into the JETSET program, as well as the complete JETSET 6.2 manual. - [Sjo87] : update notes for JETSET 6.3, and in particular a description of the new parton showering routine. - [Sjo88] : a recent summary of fragmentation models (as seen from the Lund horizon). - [QEG89] : the most recent review of QCD event generators, with a fairly extensive description of JETSET 7.1. Some of the new aspects included in JETSET 7 are still not well documented anywhere, unfortunately. ______________________________________________________________________ 1.6. Installation of Program The program is written entirely in standard Fortran 77, and should run on any machine with such a compiler. Unfortunately, experience with IBM, VAX and Norsk Data compilers has been uniform: the options available for obtaining optimized code actually produce erroneous code (e.g., operations inside DO loops are moved out before them, where some of the variables have not yet been properly set). Therefore the general advice is to use a low optimizing level, like OPTIMIZE(1) on IBM, /NOOPT on VAX, OPT OFF on ND, etc. Note that this is often not the default setting. Specifically on the Apollo, there seems to exist a problem with non-compliance to the Fortran 77 standard for the handling of unformatted write and backspace, as used by RLUGET and RLUSET. This does not form a necessary part of the event generator itself, and so RLUGET and RLUSET could well just be commented out, or some alternative solution found. SAVE statements have been included in accordance with the Fortran standard. Since most ordinary machines take SAVE for granted, this part is not particularly well tried out, however. Users on machines without automatic SAVE are therefore warned to be on the lookout for any variables which may have been missed. All default settings and particle data are stored in the BLOCK DATA LUDATA. This subprogram must be linked for a proper functioning of the other routines. On some machines this is not done automatically but must be forced by the user, in particular if JETSET is maintained as a library from which routines are to be loaded only when they are needed. In this connection we note that the library approach does not give any significant space advantages over a loading of the whole JETSET package as a unit, since a normal run will call on most of the routines anyway, directly or indirectly. Since most machines in current use are 32 bit ones, this is the precision normally assumed. A few pieces of code have therefore had to be written in double precision. All double precision variables have as first character D, and this character is never used for single precision variables. Therefore the declaration is done with IMPLICIT DOUBLE PRECISION(D). The only double precision constants that appear are 0D0, 1D0 and 2D0. A user on a 64 bit machine could therefore easily transform the program to single precision throughout. He (she) is then strongly urged to mark this version clearly, and see to it that this version is not handed on to a 32 bit machine user, who could find him(her)self in deep trouble. For applications at very high energies, like SSC, the use of single precision for any real variable starts to become a problem. It might then be necessary to rewrite the program completely, i.e. extend the range of the declaration to IMPLICIT DOUBLE PRECISION(A-H,O-Z), and change all real constants to double precision. Needless to say, the latter is a major undertaking. In some cases, shortcuts are available. On the IBM, e.g., the AUTODBL compiler option for automatic precision doubling works fine, provided only that an extra dummy integer variable (NPAD, say) is inserted directly after N in the LUJETS commonblock, to make the number of integers even. Some pieces of code will then actually run in quadruple precision; this could be corrected as described above. A test program, LUTEST, is included in the JETSET package. It is disguised as a subroutine, so that the user has to run a main program CALL LUTEST(1) END This program will generate six hundred events of different types, under a variety of conditions. If the program has not been properly installed, this program is likely to crash, or at least generate a number of erroneous events. This will then clearly be marked in the output, which otherwise will just contain a few sample event listings and a table of the number of different particles produced. To switch off the output of normal events and final table, use LUTEST(0) instead of LUTEST(1). The final tally of errors detected should read 0. ______________________________________________________________________ 1.7. Random Numbers The construction of a good, portable (pseudo)random generator is not a trivial task. Therefore the Lund Monte Carlo has traditionally stayed away from that area, and just provided the routine RLU as an interface, which the user could modify to call on an existing routine, implemented on the actual machine being used. In recent years, progress has been made in constructing portable generators with large periods and other good properties; see the review [Jam90]. Therefore the current version contains a random number generator based on the algorithm proposed by Marsaglia, Zaman and Tsang [Mar90]. This routine should work on any machine with an at least 24-digit mantissa, i.e. all common 32-bit (or more) computers. Given the same initial state, the sequence will also be identical on different machines. This need not mean that the same sequence of events will be generated on an IBM and a VAX, say, since the different treatments of roundoff errors in numerical operations will lead to slightly different real numbers being tested against these random numbers in IF-statements. Apart from nomenclature issues, and the coding of RLU as a function rather than a subroutine, the only difference between the JETSET code and the code given in [Jam90] is that slightly different algorithms are used to ensure that the random number is not equal to 0 or 1 within machine precision. The generator has a period of roughly 2*10^43, and the possibility to obtain almost 10^9 different and disjoint subsequences, selected by giving an initial integer numbers. The price to be paid for the long period is that the state of the generator at a given moment can not be described by a single integer, but requires about 100 words. Some of these are real numbers, and are thus not correctly represented in decimal form. The normal procedure, with being able to restart the generation from a seed value written among the run output, is therefore not convenient. The CERN library implementation keeps track of the number of random numbers generated since the start. With this value saved, in a subsequent run the random generator can be asked to skip ahead the corresponding number of numbers. The Lund Monte Carlo is a heavy user of random numbers, however: typically 30% of the full run time is spent on random number generation. Of this, half is overhead coming from the function call administration, but the other half is truly related to the speed of the algorithm. Therefore a skipping ahead would take place with 15% the time cost of the original run, i.e. an uncomfortably high figure. Instead a different solution is chosen here. Two special routines are provided for writing and reading the state of the random number (plus some initialization information) on a sequential file, in a machine dependent internal representation. The file used for this purpose has to be specified by the user, and opened for read and write. A state is written as a single record, in free format. It is possible to write an arbitrary number of states on a file, and a record can be overwritten, if so desired. The event generation loop might then look something like: 1) save the state of the generator on file (using flag set in point 3 below), 2) generate an event, 3) study the event for errors or other reasons why to regenerate it later; set flag to overwrite previous generator state if no errors, else set flag to create new record; 4) loop back to point 1. In the end, the file will then contain the state before each of the problematical events. An alternative approach might be to save the state every 100 events or so. It should be emphasized that this option has not been included because errors in JETSET are likely to be frequent, but because a subsequent detector simulation of this event might fail, or for other similar reasons. (The user may then have to save also other sets of seeds, naturally.) In addition to the service routines, the commonblock which contains the state of the generator is available for manipulation by the user, if he/she so desires. In particular, the initial seed value is by default 19780503, i.e. different from the Marsaglia/CERN default 54217137. It is possible to change this value before any random numbers have been generated, or to force reinitialization in mid-run with any desired new seed. Inside JETSET, or other Lund family programs, some initialization may take place in connection with the very first event generated in a run, so sometimes it may be necessary to generate one ordinary event before reading in a saved state to generate an interesting event. It should should be noted that, of course, the appearance of a random number generator package inside JETSET does in no way preclude a user from removing these routines. He/she can easily revert to the old approach, where RLU is nothing but an interface to an arbitrary external random number generator; e.g., to call a routine RNDM all one needs to have is FUNCTION RLU(IDUM) 100 RLU=RNDM(IDUM) IF(RLU.LE.0..OR.RLU.GE.1.) GOTO 100 RETURN END The random generator subpackage consists of the following components. FUNCTION RLU(IDUM) Purpose: to generate a (pseudo)random number uniformly in the range 0 < RLU < 1, i.e. excluding the endpoints. IDUM : dummy input argument. SUBROUTINE RLUGET(LFN,MOVE) Purpose: to dump the current state of the random number generator on a separate file, using internal representation for real and integer numbers. To be precise, the full contents of the LUDATR commonblock are written on the file, with the exception of MRLU(6). LFN : the file number to which the state is dumped. The user must associate this number with a true file (with a machine-dependent name), and see to it that this file is open for write. MOVE : choice of adding a new record to the file or overwriting old record. Normally only options 0 or -1 should be used. = 0 (or > 0) : add a new record to the end of the file. = -1 : overwrite the last record with a new one (i.e. do one BACSPACE before the new write). = -n : back up n records before writing the new record. The records following after the new record are lost, i.e. the last n old records are lost and one new added. SUBROUTINE RLUSET(LFN,MOVE) Purpose: to read in a state for the random number generator, from which the subsequent generation can proceed. The state must previously have been saved by a RLUGET call. Again the full contents of the LUDATR commonblock are read, with the exception of MRLU(6). LFN : the file number from which the state is read. The user must associate this number with a true file previously written with a RLUGET call, and see to it that this file is open for read. MOVE : positioning in file before a record is read. With zero value, records are read one after the other for each new call, while nonzero values may be used to navigate back and forth, and e.g. return to the same initial state several times. = 0 : read the next record. = +n : skip ahead n records before reading the record that sets the state of the random number generator. = -n : back up n records before reading the record that sets the state of the random number generator. COMMON/LUDATR/MRLU(6),RRLU(100) Purpose: to contain the state of the random number generator at any moment (for communication between RLU, RLUGET and RLUSET), and also to provide the user with the possibility to initialize different random number sequences, and to know how many numbers have been generated. MRLU(1) : (D=19780503) the integer number which specifies which of the possible subsequences will be initialized in the next RLU call for which MRLU(2)=0. Allowed values are 0 <= MRLU(1) <= 900 000 000, the original Marsaglia (and CERN library) seed is 54217137. The MRLU(1) value is not changed by any of the JETSET routines. MRLU(2) : (D=0) initialization flag, put to 1 in the first RLU call of run. A reinitialization of the random number generator can be made in mid-run by resetting MRLU(2) to 0 by hand. In addition, anytime the counter MRLU(3) reaches 1000000000, it is reset to 0 and MRLU(2) is increased by 1. MRLU(3) : (D=0) counter for the number of random numbers generated from the beginning of the run. To avoid overflow when very many numbers are generated, MRLU(2) is used as described above. MRLU(4), MRLU(5) : I97 and J97 of the CERN library implementation; part of the state of the generator. MRLU(6) : (D=0) current position, i.e. how many records after beginning, in the file used by RLUGET and RLUSET. RRLU(1) - RRLU(97) : the U array of the CERN library implementation; part of the state of the generator. RRLU(98) - RRLU(100) : C, CD and CM of the CERN library implementation; the first part of the state of the generator, the latter two constants calculated at initialization. ______________________________________________________________________ 1.8. Programming Philosophy The Monte Carlo program is built as a slave system, i.e. the user supplies the main program, and from this the various subroutines are called on to execute specific tasks, after which control is returned to the main program. Some of these tasks may be very trivial, whereas the "high-level" routines by themselves may make a large number of subroutine calls. It should be noted that, while the physics content is obviously at the center of attention, the Lund Monte Carlo also contains a more extensive setup of auxiliary service routines than any other physics event generator. The hope is that this will provide a comfortable working environment, where not only events are generated, but users also linger on to perform a lot of the subsequent studies. (As for the relatively small attention given to physics in this manual, the reason is that the physics is documented separately in a series of papers, but the program pieces only here.) The general rule is that all routines have names six characters long, beginning with LU, except for real-valued functions, which start off with UL instead. There are three exceptions: KLU, PLU and RLU. The former two functions are strongly coupled to the K and P matrices in the LUJETS commonblock, the latter uses R to emphasize the role as a random number generator. Also commonblocks have names starting with LU. For most of the routines no initialization is necessary, except for the one implied by the presence of BLOCK DATA LUDATA; this subprogram must be linked, however, which does not occur automatically with all loaders. The cases where some initialization may indeed be performed (depending on exactly which options are used), and hence events may have to be generated in some coherent fashion, are LUEEVT (and some routines called by it) and LUONIA. In addition, the random number generator is initialized at the first call; see the preceding section. Apart from writing a header, printing error messages if need be, and responding to explicit requests for listings, all tasks of the program are performed silently. All output is directed to unit MSTU(11), by default 6, and it is up to the user to see to it that this unit is open for write. The one exception is LUUPDA where, for obvious reasons, the input/output file is specified at each call. Here the user again has to see to it that proper read/write access is set. The Lund Monte Carlo is extremely versatile, but the price to be paid for this is a large number of adjustable parameters and switches for alternative modes of operation. No single user is ever likely to have need for more than a fraction of the options available. Since all these parameters and switches are assigned sensible default values, there is no reason to worry about them until the need arises. Unless explicitly stated (or obvious from the context) all switches and parameters can be changed independently of each other. One should note, however, that if only a few switches/parameters are changed, this may result in an artificially bad agreement with data. Many disagreements can often be cured by a subsequent retuning of some other parameters of the model, in particular those that have anyway once been determined by a comparison with data in the context of the default scenario. Typically such a retuning could involve one QCD parameter (alpha_S or Lambda), the longitudinal fragmentation function, and the average transverse fragmentation momentum. The program contains a number of checks that flavours specified for jet systems make sense, that the energy is enough to allow hadronization, that the memory space in LUJETS is enough, etc. If anything goes wrong that the program can catch (obviously that may not always be possible) an error message will be printed and the treatment of the corresponding event will be cut short. So long as no error messages appear on the output, it may not be worth the while to look into the rules for error checking, but if but one message appears, it should be enough cause for alarm to receive prompt attention. Also warnings are sometimes printed. These are less serious, and the experienced used might deliberately do operations which goes against the rules, but still can be made to make sense in their context. Only the first few warnings will be printed, thereafter the program will be quite. By default, the program is set to stop execution after ten errors, however, after printing the last erroneous event. It must be emphasized that not all errors will be caught. In particular, one tricky question is what happens if an integer-valued commonblock switch or subroutine/function argument is used with a value not defined. In some subroutine calls, a prompt return will be expedited, but in most instances the subsequent action is entirely unpredictable, and often completely haywire. The same goes for real-valued input assigned values outside the physically sensible range. One example will here suffice: if PARJ(2) is defined as the s/u suppression factor, a value > 1 will not give more profuse s than u production, but actually a spillover into c production. Users, beware! ______________________________________________________________________ 1.9. The First Steps This section is intended as a brief guided tour through some highlights of the Lund Monte Carlo, for persons with no previous experience in using the Lund Monte Carlo. So, Welcome to the Lund World!, or, Please move on!, as the case may be. As a first example, assume that you want to study the production of uubar two-jet systems at energy 20 GeV. To do this, write a main program CALL LU2ENT(0,2,-2,20.) CALL LULIST(1) END and run this program, linked together with JETSET. The routine LU2ENT is specifically intended for storing two entries (jets or particles). The first argument (0) is a command to perform fragmentation and decay directly after the entries have been stored, the second and third that the two entries are u (2) and ubar (-2), and the last that the CM energy of the pair is 20 GeV. When this is run, the resulting event is stored in the LUJETS commonblock. This information can then be read out by the user. No output is produced by LU2ENT itself, except for the title lines The Lund Monte Carlo - JETSET version 7.3 ** Last date of change: 22 May 1990 ** which appear once for every JETSET run. Instead the second command, to LULIST, provides a simple visible summary of the information stored in LUJETS. The argument (1) indicates that the short version should be used, which is suitable for viewing the listing directly on an 80 column terminal screen. It might look something like Event listing (summary) I particle/jet KS KF orig p_x p_y p_z E m 1 (u) A 12 2 0 0.000 0.000 10.000 10.000 0.006 2 (u~) V 11 -2 0 0.000 0.000 -10.000 10.000 0.006 3 (string) 11 92 1 0.000 0.000 0.000 20.000 20.000 4 (rho+) 11 213 3 0.098 -0.154 2.710 2.856 0.885 5 (rho-) 11 -213 3 -0.227 0.145 6.538 6.590 0.781 6 pi+ 1 211 3 0.125 -0.266 0.097 0.339 0.140 7 (Sigma0) 11 3212 3 -0.254 0.034 -1.397 1.855 1.193 8 (K*+) 11 323 3 -0.124 0.709 -2.753 2.968 0.846 9 p~- 1 -2212 3 0.395 -0.614 -3.806 3.988 0.938 10 pi- 1 -211 3 -0.013 0.146 -1.389 1.403 0.140 11 pi+ 1 211 4 0.109 -0.456 2.164 2.218 0.140 12 (pi0) 11 111 4 -0.011 0.301 0.546 0.638 0.135 13 pi- 1 -211 5 0.089 0.343 2.089 2.124 0.140 14 (pi0) 11 111 5 -0.316 -0.197 4.449 4.467 0.135 15 (Lambda0) 11 3122 7 -0.208 0.014 -1.403 1.804 1.116 16 gamma 1 22 7 -0.046 0.020 0.006 0.050 0.000 17 K+ 1 321 8 -0.084 0.299 -2.139 2.217 0.494 18 (pi0) 11 111 8 -0.040 0.410 -0.614 0.751 0.135 19 gamma 1 22 12 0.059 0.146 0.224 0.274 0.000 20 gamma 1 22 12 -0.070 0.155 0.322 0.364 0.000 21 gamma 1 22 14 -0.322 -0.162 4.027 4.043 0.000 22 gamma 1 22 14 0.006 -0.035 0.422 0.423 0.000 23 p+ 1 2212 15 -0.178 0.033 -1.343 1.649 0.938 24 pi- 1 -211 15 -0.030 -0.018 -0.059 0.156 0.140 25 gamma 1 22 18 -0.006 0.384 -0.585 0.699 0.000 26 gamma 1 22 18 -0.034 0.026 -0.029 0.052 0.000 sum: 0.00 0.000 0.000 0.000 20.000 20.000 (a few blanks have been removed between the columns to make it fit into the 72 column format of this file). Look in the particle/jet column and note that the first two lines are the original u and ubar, where 'bar' is actually written '~' to save space in longer names. The parantheses enclosing these names are there as a reminder that these jets actually have been allowed to fragment. They are retained so that event histories can be studied. Also note that the KF (flavour code) column contains 2 in the first line and -2 in the second. These are the actually stored codes denoting the presence of a u and a ubar, cf. the LU2ENT call, while the names written are just conveniences used when producing visible output. The A and V near the end of particle/jet column indicate the beginning and end of a string (or cluster, or independent fragmentation) parton system; any intermediate entries belonging to the same system would have had an I in that column. (This gives a poor man's vertical representation of a doublesided arrow, cf. <----->.) In the orig (origin) column, the zeros indicate that u and ubar are two initial entries. The subsequent line, number 3, denotes the fragmenting u+ubar string system as a whole, and has orig 1, since the first parton of this string system is entry number 1. The particles in lines 4 - 10 have orig 3 to denote that they come directly from the fragmentation of this string. In string fragmentation it is not meaningful to say that a particle comes from only the u quark or only the ubar one. It is the string system as a whole which gives a rho+, a rho-, a pi+, a Sigma0, a K*+, a p~-, and a pi-. Note that some of the particle names are again enclosed by parantheses, indicating that these particles are also not present in the final state, but have decayed further. Thus the pi- in line 13 and the pi0 in line 14 have orig 5, as an indication that they come from the decay of the rho- in line 5. Only the names not enclosed in parantheses are the ones that remain at the end of the fragmentation/decay chain, and thus are experimentally observable. The actual status code used to distinguish between different classes of entries is given in the KS column; codes in the range 1 - 10 correspond to remaining entries. The columns with p_x, p_y, p_z, E and m are rather self-explanatory. All momenta, energies and masses are given in units of GeV, with c = 1. Note that energy and momentum is conserved at each step of the fragmentation/decay process (although there exist options where this is not true). Also note that the z axis plays the role of preferred direction, along which the original partons are put. The final line is intended as a quick check that nothing funny happened. It contains the summed charge, summed momentum, summed energy and invariant mass of the final entries at the end of the fragmentation/decay chain, and the values should agree with the input ones implied by the LU2ENT arguments. (In fact, warnings would appear on the output if anything untoward happened, but that is another story). The example above has illustrated roughly what information is to be had in the event record, but not so much about how it is stored. This is better seen by using a 132 column format for listing events. Try e.g. the following program CALL LU3ENT(0,1,21,-1,30.,0.9,0.7) CALL LULIST(2) CALL LUEDIT(3) CALL LULIST(2) END where a 3-jet dgdbar event is generated in line 1 and listed in line 2. This listing will contain the numbers as directly stored in the commonblock LUJETS COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5) For particle I, K(I,1) thus gives information on whether a jet/particle has fragmented/decayed or not, K(I,2) gives the particle code, K(I,3) the origin, K(I,4) and K(I,5) the position of fragmentation/decay products, and P(I,1) - P(I,5) momentum, energy and mass. The number of lines in current use is given by N, i.e. 1 <= I <= N. The V matrix contains decay vertices; to view those LULIST(3) has to be used. It is important to learn the rules for how information is stored in LUJETS, see sections 2.1 and 2.2. The third line in the program illustrates another important point about JETSET: a number of routines are available for manipulating the event record after an event has been generated. Thus LUEDIT(3) will remove everything except stable charged particles, as shown by the result of the second LULIST call. More advanced possibilities include things like sphericity or clustering routines. Apart from the input arguments of subroutine calls, control on the doings of JETSET may be imposed via the LUDAT1, LUDAT2, LUDAT3 and LUDAT4 commonblocks. Here sensible default values are always provided. A user might want to switch off all particle decays by putting MSTJ(21) = 0 or increase the s/u ratio in fragmentation by putting PARJ(2) = 0.40, to give but two examples. It is by exploring the possibilities here that JETSET can be turned into an extremely versatile tool, even if all the nice physics is already present in the default values. As a final, semirealistic example, assume that the pT spectrum of pi+ particles is to be studied in 93 GeV e+e- annihilation events, where pT is to be defined with respect to the sphericity axis. Using the HBOOK package (version 4, watch out for version- or installation- specific differences) for histogramming, a complete program might look like COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5) COMMON/PAWC/HMEMOR(10000) CALL HLIMIT(10000) CALL HBOOK1(1,'pT spectrum of pi+',100,0.,5.,0.) NEVT=100 DO 110 IEVT=1,NEVT CALL LUEEVT(0,93.) IF(IEVT.EQ.1) CALL LULIST(1) CALL LUSPHE(SPH,APL) CALL LUEDIT(31) DO 100 I=1,N IF(K(I,2).NE.211) GOTO 100 PT=SQRT(P(I,1)**2+P(I,2)**2) CALL HF1(1,PT,1.) 100 CONTINUE 110 CONTINUE CALL HOPERA(1,'+',1,1,20./NEVT,0.) CALL HISTDO END Study this program, try to understand what happens at each step, and run it to check that it works. After that you should be ready to look at the relevant sections of this manual and start writing your own programs. Good luck! ********************************************************************** 2. The Fragmentation/Decay Package The fragmentation and decay routines form the heart of the Lund family of Monte Carlos. The hard kinematics may be specific for e+e- annihilation (section 3), hadron collisions (as in PYTHIA), leptoproduction (as in LEPTO), etc., but in the end they all use the fragmentation routines described here. In particular, the particle code and the organization of the event record is common to all members of the Lund family. (At the time of reading this, the switchover to the new codes may not have taken place in all programs. In the case of PYTHIA, the work is in progress.) ______________________________________________________________________ 2.1. Particle codes The new particle code now adopted by the Particle Data Group [PDG88] is used consistently throughout the program, and is referred to as the KF particle code. This code the user has to be thoroughly familiar with. It is described below. The KF code is not convenient for a direct storing of masses, decay data, or other particle properties, however, since the KF codes are so spread out. Instead a compressed code KC between 1 and 500 is used here, where the most frequently used particles have a separate code, but many heavy flavour hadrons are lumped together in groups. Normally this code is only used at very specific places in the program, not visible to the user. If need be, the correspondence can always be obtained by using the function LUCOMP, KC = LUCOMP(KF). It is therefore not the intention that a user should ever need to know any KC code at all. Below follows a list of KF particle codes. The list is not complete; a more extensive one may be obtained with CALL LULIST(11). Particles are grouped together, and the basic rules are described for each group. Whenever a distinct antiparticle exists, it is given the same KF code with a minus sign (whereas KC codes are always positive). The particle names printed here also corresponds to the ones obtained with the routine LUNAME, which is used extensively, e.g. in LULIST. Greek characters are spelt out in full, with a capital first letter corresponding to a capital greek letter. Generically the name of a particle is made up out of the following pieces: a) The basic root name. This includes a * for most spin 1 (L = 0) mesons and spin 3/2 baryons, and a ' for some spin 1/2 baryons (where there are two states to be distinguished, cf. Lambda - Sigma0). The rules for heavy baryon naming are in accordance with the 1986 Particle Data Group conventions [PDG86]. For mesons with one unit of orbital angular momentum, K (D, B, ..) is used for quark spin 0 and K* (D*, B*, ..) for quark spin 1 mesons; the convention for * may here deviate slightly from the one used by the PDG. b) Any lower indices, separated from the root with a _. For heavy hadrons, this is the additional heavy flavour content not inherent from the root itself. For a diquark, it is the spin. c) The character ~ (alternatively bar, see MSTU(15)) for an antiparticle, wherever the distinction between particle and antiparticle is not inherent in the charge information. d) Charge information: ++, +, 0, -, or --. Charge is not given for quarks or diquarks. Some neutral particles which customarily are given without a 0 also here lack it, like neutrinos, g, gamma, and flavour diagonal mesons other than pi0 and rho0. Note that charge is included both for the proton and the neutron. While nonstandard, it is helpful in avoiding misunderstandings when looking at an event listing. IMPORTANT NOTE: Apart from a trivial printing error for the Omega- baryon in the PDG tables, a few inconsistencies between the KF and the PDG codes are known. These stem from differences of interpretation of the rules agreed on, that form the basis of the PDG tables and (independently) of the JETSET tables. (Of course, my private opinion is that I follow the original agreement, and PDG deviate from it.) Hopefully, this should have few practical consequences, since only rarely produced particles are affected. Anyway, here is a list of the known discrepancies: 1) PDG has not allowed for the existence of an eta_b, which in JETSET is included with code 551. This code is reserved for chi_b0 in PDG, a particle which appears as 10551 in JETSET. (We agree to have eta_c as 441, which illustrates the basic difference: I use the additional recurrence figure to refer to a whole multiplet, whether all particles of that multiplet have been found or not, while PDG does not reserve space for particles which we know should be there but has not yet been discovered, which means that a multiplet need not go together.) 2) PDG has not allowed for the existence of an h_1c, which in JETSET is represented by 10443. Therefore chi_1c is PDG code 10443 but JETSET code 20443. (Comment as for 1.) 3) Different conventions for spin 1/2 baryons with one heavy flavour (charm, bottom, top), one strange flavour, and one light (u or d). Here two states exist, e.g. Xi_c+ and Xi'_c+, both with flavour content csu. By analogy with the Lambda0-Sigma system, JETSET uses decreasing order of flavour content for the heavier state and inversed order of the two lighter flavours for the lighter state, while the PDG tables use the opposite convention. Thus in JETSET Xi_c+ is 4232 and Xi'_c 4322, while in PDG it is the other way around. There are no imminent plans to change the JETSET rules to agree with the PDG ones in either of the cases above. List of KF codes: 1) Quarks and leptons. This group contains the basic building blocks of matter, arranged according to family, with the lower member of weak isodoublets also having the smaller code (thus d precedes u, contrary to the ordering in previous JETSET versions). A fourth generation is included for future reference. The quark codes are used as building blocks for the diquark, meson and baryon codes below. 1 d 11 e- 2 u 12 nu_e 3 s 13 mu- 4 c 14 nu_mu 5 b 15 tau- 6 t 16 nu_tau 7 l 17 chi- 8 h 18 nu_chi 9 19 10 20 2) Gauge bosons and other fundamental bosons. This group includes all the gauge and Higgs bosons of the standard model, as well as some of the bosons appearing in various extensions of the standard model. The latter are not covered by the standard Particle Data group codes. They correspond to one extra U(1) group and one extra SU(2) one, a further Higgs doublet, a (scalar) leptoquark LQ, and a horizontal gauge boson R (coupling between families). 21 gluon 31 22 gamma 32 Z'0 23 Z0 33 Z"0 24 W+ 34 W'+ 25 H0 35 H'0 26 36 A0 27 37 H+ 28 38 29 39 LQ 30 40 R0 3) Free space The positions 41 - 80 are currently unused. In the future, they might come to be used e.g. for supersymmetric partners to the particles above, or for some other kind of new physics. At the moment, they are at the disposal of the user. 4) Various special codes. In a Monte Carlo, it is always necessary to have codes that do not correspond to any specific particle, but are used to lump together groups of similar particles for decay treatment, or to specify generic decay products. These codes, which again are non-standard, are found between numbers 81 and 100. Several are not found in the event record, and therefore properly belong only to the KC group of codes. In the lines below, auxiliary information is given within parentheses. 81 specflav (spectator flavour; used in decay product listings) 82 rndmflav (a random u, d, or s flavour; possible decay product) 83 phasespa (simple isotropic phase space decay) 84 c-hadron (information on decay of generic charm hadron) 85 b-hadron (information on decay of generic bottom hadron) 86 t-hadron (information on decay of generic top hadron) 87 l-hadron (information on decay of generic low hadron) 88 h-hadron (information on decay of generic high hadron) 89 Wvirt (off-mass-shell W in weak decays of t, l, h or chi) 90 diquark (generic code for colour information) 91 cluster (parton system in cluster fragmentation) 92 string (parton system in string fragmentation) 93 indep. (parton system in independent fragmentation) 94 CMshower (four-momentum of timelike showering system) 95 SPHEaxis (event axis found with LUSPHE) 96 THRUaxis (event axis found with LUTHRU) 97 CLUSjet (jet (cluster) found with LUCLUS) 98 CELLjet (jet (cluster) found with LUCELL) 99 table (tabular output from LUTABU) 100 5) Diquark codes. A diquark made up of a quark with code i and another with code j, where i >= j, and with total spin s, is given the code KF = 1000*i + 100*j + 2s+1 i.e. the tens position is left empty (cf. the baryon code below). Some of the most frequently used codes are listed below. All the lowest-lying spin 0 and 1 diquarks are included in the program. 1103 dd_1 2101 ud_0 2103 ud_1 2203 uu_1 3101 sd_0 3103 sd_1 3201 su_0 3203 su_1 3303 ss_1 The corresponding KC code is 90, and is mainly used to store colour charge. 6) Meson codes. A meson made up of a quark with code i and an antiquark with code -j, and with total spin s, is given the code KF = (100*max(i,j) + 10*min(i,j) + 2s+1) * sign(i-j) * (-1)**max(i,j) i.e. if the heaviest quark is a down-type one, an extra - sign enters. This is in accordance with the particle-antiparticle distinction adopted in the 1986 Review of Particle Properties [PDG86]. (Compared to the JETSET 6.3 conventions, however, it corresponds to a sign change for B mesons and fourth generation L mesons.) The flavour-diagonal states are arranged in order of ascending mass. The standard rule of having the last digit being 2s+1 is broken for the K_S0 - K_L0 system, where it is 0, and this convention should carry over to mixed states in the B meson system. For higher multiplets with the same spin, +-10000, +-20000, etc., are added to provide the extra distinction needed. Some of the most frequently used codes are given below. The full lowest-lying pseudoscalar and vector multiplets are included in the program. 211 pi+ 213 rho+ 311 K0 313 K*0 321 K+ 323 K*+ 411 D+ 413 D*+ 421 D0 423 D*0 431 D_s+ 433 D*_s+ 511 B0 513 B*0 521 B+ 523 B*+ 531 B_s0 533 B*_s0 111 pi0 113 rho0 221 eta 223 omega 331 eta' 333 phi 441 eta_c 443 J/psi 551 eta_b 553 Upsilon 661 eta_t 663 Theta 130 K_L0 310 K_S0 Also the lowest lying orbital angular momentum L = 1 mesons are included : one pseudovector multiplet obtained for total quark spin 0 (L = 1, S = 0 -> J = 1) and one scalar, one pseudovector and one tensor multiplet obtained for total quark spin 1 (L = 1, S = 1 -> J = 0, 1 or 2). Any mixing between the two pseudovector multiplets is not taken into account. Please note that some members of these multiplets have still not been found, and are included here only based on guesswork. Even for known ones, particle data (mass, width, decay modes) is highly incomplete. 10213 b_1+ 10211 a_0+ 10313 K_10 10311 K*_00 10323 K_1+ 10321 K*_0+ 10413 D_1+ 10411 D*_0+ 10423 D_10 10421 D*_00 10433 D_1s+ 10431 D*_0s+ 10113 b_10 10111 a_00 10223 h_10 10221 f_00 10333 h'_10 10331 f'_00 10443 h_1c0 10441 chi_0c0 20213 a_1+ 215 a_2+ 20313 K*_10 315 K*_20 20323 K*_1+ 325 K*_2+ 20413 D*_1+ 415 D*_2+ 20423 D*_10 425 D*_20 20433 D*_1s+ 435 D*_2s+ 20113 a_10 115 a_20 20223 f_10 225 f_20 20333 f'_10 335 f'20 20443 chi_1c0 445 chi_2c0 The corresponding meson KC codes, used for organizing mass and decay data, range between 101 and 240. 7) Baryon codes. A baryon made up of quarks i, j and k, with i >= j >= k, and total spin s, is given the code KF = 1000*i + 100*j + 10*k + 2s+1. An exception is provided by spin 1/2 baryons made up of three different type quarks, where the two lightest quarks form a spin 0 diquark (Lambda-like baryons). Here the order of the j and k quarks is changed, so as to provide a simple means of distinction to baryons with the lightest quarks in a spin 1 diquark (Sigma-like baryons). For hadrons with heavy flavours, the root names are Lambda or Sigma for hadrons with two u or d quarks, Xi for those with one and Omega for those without u or d quarks. Some of the most frequently used codes are given below. The full lowest-lying spin 1/2 and 3/2 multiplets are included in the program. 1114 Delta- 2112 n 2114 Delta0 2212 p 2214 Delta+ 2224 Delta++ 3112 Sigma- 3114 Sigma*- 3122 Lambda0 3212 Sigma0 3214 Sigma*0 3222 Sigma+ 3224 Sigma*+ 3312 Xi- 3314 Xi*- 3322 Xi0 3324 Xi*0 3334 Omega- 4112 Sigma_c0 4114 Sigma*_c0 4122 Lambda_c+ 4212 Sigma_c+ 4214 Sigma*_c+ 4222 Sigma_c++ 4224 Sigma*_c++ 4132 Xi_c0 4312 Xi'_c0 4314 Xi*_c0 4232 Xi_c+ 4322 Xi'_c+ 4324 Xi*_c+ 4332 Omega_c0 4334 Omega*_c0 5112 Sigma_b- 5114 Sigma*_b- 5122 Lambda_b0 5212 Sigma_b0 5214 Sigma*_b0 5222 Sigma_b+ 5224 Sigma*_b+ The corresponding KC codes, used for organizing mass and decay data, range between 301 and 400, with some slots still free. 8) Diffractive states. These codes are not standard ones, but have been defined by analogy to be used for denoting diffractive states in the PYTHIA program, as part of the event history. The first two or three digits give flavour content, while the last one is 0, to denote the somewhat unusual character of the code. Only three codes have been introduced. 210 pi_diffr+ 2110 n_diffr 2210 p_diffr+ 9) Free compressed codes The positions 401 - 500 of mass and decay arrays are left open. Here a user may map any new kind of particle from the ordinary KF codes, which probably are above 10000, into a more manageable KC range for mass and decay data information. The mapping must be implemented in the LUCOMP function. ______________________________________________________________________ 2.2. The Event Record Each new event generated is in its entirety stored in the commonblock LUJETS, which thus forms the event record. Here each jet or particle that appears at some stage of the fragmentation or decay chain will occupy one line in the matrices. The different components of this line will tell which jet/particle it is, from where it originates, its present status (fragmented/decayed or not), its momentum, energy and mass, and the space-time position of its production vertex. Note that K(I,3) - K(I,5) and the P and V vectors may take special meaning for some specific applications (e.g. sphericity or cluster analysis), as described in those connections. COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5) Purpose: to contain the event record, i.e. the complete list of all partons and particles in the current event. N : number of lines in the K, P and V matrices occupied by the current event. N is continuously updated as the definition of the original configuration and the treatment of fragmentation and decay proceed. In the following, the individual parton/particle number, running between 1 and N, is called I. K(I,1) : status code KS, which gives the current status of the parton/particle stored in the line. The ground rule is that codes 1 - 10 correspond to currently existing partons/particles, while larger codes contain partons/particles which no longer exist, or other kinds of event information. = 0 : empty line. = 1 : an undecayed particle or an unfragmented jet, the latter either being a single jet or the last one of a jet system. = 2 : an unfragmented jet, which is followed by more jets in the same colour singlet jet system. = 3 : an unfragmented jet with special colour flow information stored in K(I,4) and K(I,5), such that adjacent partons along the string need not follow after each other in the event record. = 4 : a particle which could have decayed, but did not do it within the allowed volume around the original vertex. = 5 : a particle which is is to be forced to decay in the next LUEXEC call, in the vertex position given (this code is only set by user intervention). = 11 : a decayed particle or a fragmented jet, the latter either being a single jet or the last one of a jet system, cf. =1. = 12 : a fragmented jet, which is followed by more jets in the same colour singlet jet system, cf. =2. = 13 : a jet which has been removed when special colour flow information has been used to rearrange a jet system, cf. =3. = 14 : a parton which has branched into further partons, with special colour flow information provided, cf. =3. = 15 : a particle which has been forced to decay (by user intervention), cf. =5. = 21 : documentation lines used to give a compressed story of the event at the beginning of the event record. = 31 : lines with information on sphericity, thrust or cluster search. = 32 : tabular output, as generated by LUTABU. = 41 : junction (currently not fully implemented). < 0 : these codes are never used by the program, and are therefore usually not affected by operations on the record, like LUROBO, LULIST and event analysis routines (the exception is some LUEDIT calls, where lines are moved but not deleted). Such codes may therefore be useful in some connections. K(I,2) : parton/particle KF code, as described in a section 2.1 above. K(I,3) : line number of parent particle or jet, where known, else 0. Note that the assignment of a particle to a given jet of a jet system is unphysical, and what is given there is only related to the way the event was generated. K(I,4) : normally the line number of the first daughter; is 0 for an undecayed particle or unfragmented jet. For K(I,1) = 3, 13 or 14 it instead contains special colour flow information (for internal use only) on the form K(I,4) = 200000000*MCFR + 100000000*MCTO + 10000*ICFR + ICTO, where ICFR and ICTO give the line numbers of the partons from which the colour comes and to where it goes, respectively, and MCFR and MCTO originally are 0 and are set to 1 when the corresponding colour connection has been traced in the LUPREP rearrangement procedure. (The packing may be changed with MSTU(5).) The 'from' colour position may indicate a parton which branched to produce the current parton, or a parton created together with the current parton but with matched anticolour, while the 'to' normally indicates a parton that the current parton branches into. Thus, for setting up an initial colour configuration, it is normally only the 'from' part that is used, while the 'to' part is added by the program in a subsequent call to parton shower evolution (for final state radiation; it is the other way around for initial state radiation). Note: normally most users never have to worry about the exact rules for colour flow storage, since they are used mainly for internal purposes. However, when it is necessary to define this flow, it is recommended to use the LUJOIN routine, since that likely would reduce the chances of making a mistake. K(I,5) : normally the line number of the last daughter; is 0 for an undecayed particle or unfragmented jet. For K(I,1) = 3, 13 or 14 it instead contains special colour flow information (for internal use only) on the form K(I,5) = 200000000*MCFR + 100000000*MCTO + 10000*ICFR + ICTO, where ICFR and ICTO give the line numbers of the partons from which the anticolour comes and to where it goes, respectively, and MCFR and MCTO originally are 0 and are set to 1 when the corresponding colour connection has been traced in the LUPREP rearrangement procedure. (The packing may be changed with MSTU(5).) Comments and note: as for K(I,4). P(I,1) : p_x, momentum in the x direction, in GeV/c. P(I,2) : p_y, momentum in the y direction, in GeV/c. P(I,3) : p_z, momentum in the z direction. in GeV/c. P(I,4) : E, energy, in GeV. P(I,5) : m, mass, in GeV/c^2. In parton showers, with spacelike virtualities, i.e. where Q^2 = - m^2 > 0, P(I,5) = -Q. V(I,1) : x position of production vertex, in mm. V(I,2) : y position of production vertex, in mm. V(I,3) : z position of production vertex, in mm. V(I,4) : time of production, in mm/c (= 3.33 * 10^-12 s). V(I,5) : proper lifetime of particle, in mm/c (= 3.33 * 10^-12 s). If the particle is not expected to decay, V(I,5) = 0. A line with K(I,1) = 4, i.e. a particle that could have decayed, but did not do so within allowed region, has the proper nonzero V(I,5). In the absence of electric or magnetic fields, or other disturbances, the decay vertex V' of an unstable particle may be calculated as V'(j) = V(I,j) + V(I,5)*P(I,j)/P(I,5), j = 1 - 4. ______________________________________________________________________ 2.3. Definition of Initial Configuration or Commonblock Variables With the use of the conventions described for the event record, it is possible to specify any initial jet/particle configuration. This task is simplified for a number of often occuring situations by the existence of the filling routines below. Several calls can be combined in the specification. In case one call is enough, the complete fragmentation/decay chain may be simulated at the same time. At each call, the value of N is updated to the last line used for information in the call, so if several calls are used, they should be made with increasing IP number, or else N should be redefined by hand afterwards. It should be noted that many users do not come in direct contact with these routines, since that is taken care of by higher-level routines for specific processes (LUEEVT, PYTHIA, LEPTO, etc.). As an experiment, the routine LUGIVE contains a facility for setting various comonblock variables in a controlled and documented fashion. SUBROUTINE LU1ENT(IP,KF,PE,THE,PHI) Purpose: to add one entry to the event record, i.e. either a jet or a particle. IP : line number for the jet/particle. If IP = 0 is used, line number 1 is used and LUEXEC is called. If IP < 0, line -IP is used, with status code KS = 2 rather than 1; thus a jet system may be built up by filling all but the last jet of the system with IP < 0. KF : jet/particle flavour code. PE : jet/particle energy. If PE is smaller than the mass, the jet/particle is taken to be at rest. THE, PHI : polar and azimuthal angle for the momentum vector of the jet/particle. SUBROUTINE LU2ENT(IP,KF1,KF2,PECM) Purpose: to add two entries to the event record, i.e. either a two-jet system or two separate particles. IP : line number for the first jet/particle, with second in line IP+1. If IP = 0, lines 1 and 2 are used and LUEXEC is called. If IP < 0, lines -IP and -IP+1 are used, with status code KS = 3, i.e. with special colour connection information, so that a parton shower can be generated by a LUSHOW call, followed by a LUEXEC call, if so desired (only relevant for jets). KF1, KF2 : flavour codes for the two jets/particles. PECM : (= W) the total energy of the system. Remark: the system is given in the CM frame, with the first jet/particle going out in the +z direction. SUBROUTINE LU3ENT(IP,KF1,KF2,KF3,PECM,X1,X3) Purpose: to add three entries to the event record, i.e. either a three-jet system or three separate particles. IP : line number for the first jet/particle, with other two in IP+1 and IP+2. If IP = 0, lines 1, 2 and 3 are used and LUEXEC is called. If IP < 0, lines -IP through -IP+2 are used, with status code KS = 3, i.e. with special colour connection information, so that a parton shower can be generated by a LUSHOW call, followed by a LUEXEC call, if so desired (only relevant for jets). KF1, KF2, KF3: flavour codes for the three jets/particles. PECM : (= W) the total energy of the system. X1, X3 : x_i = 2E_i/W, i.e. twice the energy fraction taken by the i:th jet. Thus x_2 = 2 - x_1 - x_3, and need not be given. Note that not all combinations of x_i are inside the physically allowed region. Remark : the system is given in the CM frame, in the xz-plane, with the first jet going out in the +z direction and the third one having px > 0. SUBROUTINE LU4ENT(IP,KF1,KF2,KF3,KF4,PECM,X1,X2,X4,X12,X14) Purpose: to add four entries to the event record, i.e. either a four-jet system or four separate particles (or, for qqbarq'qbar' events, two two-jet systems). IP : line number for the first jet/particle, with other three in lines IP+1, IP+2 and IP+3. If IP = 0, lines 1, 2, 3 and 4 are used and LUEXEC is called. If IP < 0, lines -IP through -IP+3 are used, with status code KS = 3, i.e. with special colour connection information, so that a parton shower can be generated by a LUSHOW call, followed by a LUEXEC call, if so desired (only relevant for jets). KF1,KF2,KF3,KF4 : flavour codes for the four jets/particles. PECM : (= W) the total energy of the system. X1,X2,X4 : x_i = 2E_i/W, i.e. twice the energy fraction taken by the i:th jet. Thus x_3 = 2 - x_1 - x_2 - x_4, and need not be given. X12,X14 : x_ij = 2p_ip_j/W^2, i.e. twice the four-vector product of the momenta for jets i and j, properly normalized. With the masses known, other x_ij may be constructed from the x_i and x_ij given. Note that not all combinations of x_i and x_ij are inside the physically allowed region. Remark: the system is given in the CM frame, with the first jet going out in the +z direction and the fourth jet lying in the xz-plane with p_x > 0. The second jet will have p_y > 0 and p_y < 0 with equal probability with the third jet balancing this p_y (this corresponds to a random choice between the two possible stereoisomers). SUBROUTINE LUJOIN(NJOIN,IJOIN) Purpose: to connect a number of previously defined partons into a string configuration. Initially the partons must be given with status codes (KS = K(I,1)) 1, 2 or 3. Afterwards the partons all have status code 3, i.e. are given with full colour flow information. Compared to the normal way of defining a parton system, the partons need therefore not appear in the same sequence in the event record as they are assumed to do along the string. It is also possible to call LUSHOW for all or some of the entries making up the string formed by LUJOIN. NJOIN: the number of entries that are to be joined by one string. IJOIN: an one-dimensional array, of size at least NJOIN. The NJOIN first numbers are the positions of the partons that are to be joined, given in the order the partons are assumed to appear along the string. If the system consists entirely of gluons, the string is closed by connecting back the last to the first entry. Remarks: only one string (i.e. one colour singlet) may be defined per call, but one is at liberty to use any number of LUJOIN calls for a given event. The program will check that the parton configuration specified makes sense, and not take any action unless it is. Note, however, that an initially sensible parton configuration may become nonsensical, if only some of the partons are reconnected, while the others are left unchanged. SUBROUTINE LUGIVE(CHIN) Purpose: to set the value of any variable residing in the commmonblocks LUJETS, LUDAT1, LUDAT2, LUDAT3, LUDAT4, or LUDATR. This is done in a more controlled fashion than by directly including the commonblocks in the user program, in that array bounds are checked and the old and new values for the variable changed is written to the output for reference. CHIN : character expression of length at most 100 characters, with requests for variables to be changed, stored in the form variable1=value1;variable2=value2;variable3=value3 ... Note that an arbitrary number of instructions can be stored in one call if separated by semicolons, and that blanks may be included anyplace. The variable_i may be any single variable in the JETSET commonblocks, and the value_i must be of the correct integer, real or character (without extra quotes) type. Array indices and values must be given explicitly, i.e. can not be variables in their own right. The exception is that the first index can be preceded by a C, signifying that the index should be translated from normal KF to compressed KC code with a LUCOMP call; this is allowed for the KCHG, PMAS, MDCY and CHAF arrays. If a value_i is omitted, i.e. with the construction variable=, the current value is written to the output. Remark 1: The routine may also be used to set values for the PYTHIA commonblocks PYSUBS, PYPARS, PYINT1, PYINT2, PYINT3, PYINT4, PYINT5, and PYINT6, according to the same rules as above. Commands to JETSET and to PYTHIA may be freely mixed. Remark 2: The checks on array bounds are hardwired into this routine. Therefore, if some user changes array dimensions and MSTU(3), MSTU(6) and/or MSTU(7), as allowed by other considerations, these changes will not be known to LUGIVE. Normally this should not be a problem, however. ______________________________________________________________________ 2.4. The Physics Routines The physics routines form the major part of the program, but once the initial jet/particle configuration has been specified and default parameter values changed, if so desired, only a LUEXEC call is necessary to simulate the whole fragmentation and decay chain. The physics involved has been described in a number of different publications, see e.g. [Sjo86,Sjo88] and references therein. We will therefore only give a rather brief overview. The routines RLU, RLUGET and RLUSET have been described in section 1.6, but otherwise would also belong to this section. The routine LUSHOW, for timelike shower evolution, is somewhat of a special classification problem. Its main application is to e+e- annihilation events, as administrated from LUEEVT, and it would be natural to place the routine in that group, as has been done in the past. However, in some decays a qqbar pair is formed with such energy that gluon radiation corrections should be included, and thus LUSHOW becomes intimately intertwined with LUDECY. This is particularly acute with the increased lower mass limits for top. A 60 GeV top, say, would typically give a qqbar pair coming from the virtual W with an invariant mass of 30 - 40 GeV, i.e. comparable to those PETRA/PEP energies where gluon emission corrections are known to be quite significant. The placing of this routine here does not preclude its use elsewhere, however. SUBROUTINE LUEXEC Purpose: to administrate the fragmentation and decay chain. LUEXEC may be called several times, but only entries which have not yet been treated (more precisely, have 1 <= KS <= 10) can be affected by further calls. This may apply if more jets/particles have been added by the user, or if particles previously considered stable are now allowed to decay. The actions that will be taken during a LUEXEC call can be tailored extensively via the LUDAT1 - LUDAT3 commonblocks, in particular by setting the MSTJ values suitably. SUBROUTINE LUPREP(IP) Purpose: to rearrange parton shower endproducts (marked with KS = 3) sequentially along strings; also to (optionally) allow small jet systems to collapse into two particles or one only, in the latter case with energy and momentum to be shuffled elsewhere in the event; also to perform checks that e.g. flavours of colour singlet systems make sense. SUBROUTINE LUSTRF(IP) Purpose: to generate the fragmentation of an arbitrary colour singlet jet system according to the Lund string fragmentation model. In many respects, this routine is the very heart and soul of the Lund family of programs. All the technical details are documented in [Sjo84]. SUBROUTINE LUINDF(IP) Purpose: to handle the fragmentation of a jet system according to independent fragmentation models, and implement energy, momentum and flavour conservation, if so desired. Also the fragmentation of a single jet, not belonging to a jet system, is considered here (this is of course physical nonsense, but may sometimes be convenient for specific tasks). SUBROUTINE LUDECY(IP) Purpose: to perform a particle decay, according to known branching ratios or different kinds of models, depending on our level of knowledge. Various matrix elements are included for specific processes. SUBROUTINE LUKFDI(KFL1,KFL2,KFL3,KF) Purpose: to generate a new quark or diquark flavour and to combine it with an existing flavour to give a hadron. KFL1: incoming flavour. KFL2: extra incoming flavour, e.g. for formation of final particle, where the flavours are completely specified. Is normallly 0. KFL3: newly created flavour; is 0 if KFL2 is nonzero. KF: produced hadron. Is 0 if something went wrong (e.g. inconsistent combination of incoming flavours). SUBROUTINE LUPTDI(KFL,PX,PY) Purpose: to give transverse momentum, e.g. for a qqbar pair created in the field, according to independent Gaussian distributions in p_x and p_y. SUBROUTINE LUZDIS(KFL1,KFL3,PR,Z) Purpose: to generate the longitudinal scaling variable z in jet fragmentation, either according to the Lund symmetric fragmentation function, or according to a choice of other shapes. SUBROUTINE LUSHOW(IP1,IP2,QMAX) Purpose: to generate timelike parton showers, conventional or coherent. The performance of the program is regulated by the switches MSTJ(41) - MSTJ(49) and parameters PARJ(82) - PARJ(84). In order to keep track of the colour flow information, the positions K(I,4) and K(I,5) have to be organized properly for showering partons. Inside JETSET 7.1, this is done automatically, but for external use proper care must be taken. IP1 > 0, IP2 = 0 : generate a timelike parton shower for the parton in line IP1 in commonblock LUJETS, with maximum allowed mass QMAX. With only one parton at hand, one can not simultaneously conserve both energy and momentum: we here choose to conserve energy and jet direction, while longitudinal momentum (along the jet axis) is not conserved. IP1 > 0, IP2 > 0 : generate timelike parton showers for the two partons in lines IP1 and IP2 in the commonblock LUJETS, with maximum allowed mass for each parton QMAX. For shower evolution, the two partons are boosted to their CM frame. Energy and momentum is conserved for the pair of partons, although not for each individually. One of the two partons may be replaced by a nonradiating particle, such as a photon or a diquark; the energy and momentum of this particle will then be modified to conserve the total energy and momentum. IP1 > 0, IP2 < 0 : generate timelike parton showers for the -IP2 (at most 3) partons in lines IP1, IP1+1, ... IPI-IP2-1 in the commonblock LUJETS, with maximum allowed mass for each parton QMAX. The actions for IP2 = -1 and IP2 = -2 correspond to what is described above, but additionally IP2 = -3 may be used to generate the evolution starting from three given partons (e.g. in Upsilon -> ggg). Then the three partons are boosted to their CM frame, energy is conserved for each parton individually and momentum for the system as a whole. QMAX : the maximum allowed mass of a radiating parton, i.e. the starting value for the subsequent evolution. (In addition, the mass of a single parton may not exceed its energy, the mass of a parton in a system may not exceed the invariant mass of the system.) SUBROUTINE LUBOEI Purpose: to include Bose-Einstein effects according to a simple parametrization. By default, this routine is not called. If called, this is done after decay of short-lived resonances, but before decay of long-lived ones. See MSTJ(51) - MSTJ(52). FUNCTION ULMASS(KF) Purpose: to give the mass for a parton/particle. SUBROUTINE LUNAME(KF,CHAU) CHARACTER CHAU*16 Purpose: to give the parton/particle name (as a character string). FUNCTION LUCHGE(KF) Purpose: to give three times the charge for a parton/particle. FUNCTION LUCOMP(KF) Purpose: to give the compressed parton/particle code KC for a given KF code, as required to find entry into mass and decay data tables. Also checks whether the given KF code is actually an allowed one (i.e. known by the program), and returns 0 if not. Note that KF may be positive or negative, while the resulting KC code is never negative. SUBROUTINE LUERRM(MERR,MESSAG) Purpose: to keep track of the number of errors and warnings encountered, write out information on them, and abort the program in case of too many errors. FUNCTION ULALPS(Q2) Purpose: to calculate the running strong coupling constant alpha_strong as a function of the momentum transfer Q^2, in first or second order QCD. See MSTU(111) - MSTU(118) and PARU(111) - PARU(118). Is currently not used in parton showers, but only for matrix elements. FUNCTION ULANGL(X,Y) Purpose: to calculate the angle from the x and y coordinates. BLOCK DATA LUDATA Purpose: to give default values for variables in the LUDAT1 - LUDAT4 and LUDATR commonblocks. ______________________________________________________________________ 2.5. Event Study and Data Listing Routines After a LUEXEC call, the event generated is stored in the LUJETS commonblock, and whatever physical variable is desired may be constructed from this record. An event may be rotated, boosted or listed, and particle data may be listed or modified. Via the functions KLU and PLU the values of some frequently appearing variables may be obtained more easily. As described in section 2.6, also more detailed event shape analyses may be performed simply. SUBROUTINE LUROBO(THE,PHI,BEX,BEY,BEZ) Purpose: to perform rotations and Lorentz boosts (in that order, if both in the same call) of jet/particle momenta and vertex position variables. THE, PHI : standard polar coordinates theta, phi, giving the rotated direction of a momentum vector initially along the +z axis. BEX, BEY, BEZ : gives the direction and size beta of a Lorentz boost, such that a particle initially at rest will have p/E = beta afterwards. Remark: all entries 1 through N are affected by the transformation, unless lower and upper bounds are explicitly given by MSTU(1) and MSTU(2), or if status code KS <= 0. ENTRY LUDBRB(IMI,IMA,THE,PHI,DBEX,DBEY,DBEZ) Purpose: to perform rotations and Lorentz boosts (in that order, if both in the same call) of jet/particle momenta and vertex position variables, for a specific range of entries, and with the boost vector given in double precision. Is entry to LUROBO, mainly intended for internal use. IMI, IMA : range of entries affected by transformation, IMI <= I <= IMA. THE, PHI : standard polar coordinates theta, phi, giving the rotated direction of a momentum vector initially along the +z axis. DBEX, DBEY, DBEZ : gives the direction and size beta of a Lorentz boost, such that a particle initially at rest will have p/E = beta afterwards. Is to be given in double precision. Remark: all entries with status codes KS > 0 in the requested range are affected by the transformation. SUBROUTINE LUEDIT(MEDIT) Purpose: to exclude unstable or undetectable jets/particles from the event record. One may also use LUEDIT to store spare copies of events (specifically initial parton configuration) that can be recalled to allow e.g. different fragmentation schemes to be run through with one and the same parton configuration. Finally, an event which has been analysed with LUSPHE, LUTHRU or LUCLUS (see section 2.6) may be rotated to align the event axis with the z direction. MEDIT : tells which action is to be taken. = 0 : empty (KS = 0) and documentation (KS > 20) lines are removed. The jets/particles remaining are compressed in the beginning of the LUJETS commonblock and the N value is updated accordingly. The event history is lost, so that information stored in K(I,3), K(I,4) and K(I,5) is no longer relevant. = 1 : as =0, but in addition all jets/particles that have fragmented/decayed (KS > 10) are removed. = 2 : as =1, but also all neutrinos and unknown particles (i.e. compressed code KC = 0) are removed. = 3 : as =2, but also all uncharged, colour neutral particles are removed, leaving only charged, stable particles (and unfragmented partons, if fragmentation has not been performed). = 5 : as =0, but also all partons which have branched or been rearranged in a parton shower and all particles which have decayed are removed, leaving only the fragmenting parton configuration and the final state particles. = 11 : remove lines with K(I,1) < 0. Update event history information (in K(I,3) - K(I,5)) to refer to remaining entries. = 12 : remove lines with K(I,1) = 0. Update event history information (in K(I,3) - K(I,5)) to refer to remaining entries. = 13 : remove lines with K(I,1) = 11, 12 or 15, except for any line with K(I,2) = 94. Update event history information (in K(I,3) - K(I,5)) to refer to remaining entries. In particular, try to trace origin of daughters, for which the mother is decayed, back to entries not deleted. = 14 : remove lines with K(I,1) = 13 or 14, and also any line with K(I,2) = 94. Update event history information (in K(I,3) - K(I,5)) to refer to remaining entries. In particular, try to trace origin of rearranged jets back through the parton shower history to the shower initiator. = 15 : remove lines with K(I,1) > 20. Update event history information (in K(I,3) - K(I,5)) to refer to remaining entries. = 21 : all partons/particles in current event record are stored (as a spare copy) in bottom of commonblock LUJETS (is e.g. done to save original partons before calling LUEXEC). = 22 : partons/particles stored in bottom of event record with =21 are placed in beginning of record again, overwriting previous information there (so that e.g. a different fragmentation scheme can be used on the same partons). Since the copy at bottom is unaffected, repeated calls with =22 can be made. = 23 : primary partons/particles in the beginning of event record are marked as not fragmented or decayed, and number of entries N is updated accordingly. Is simpe substitute for =21 plus =22 when no fragmentation/decay products precede any of the original partons/particles. = 31 : rotate largest axis, determined by LUSPHE, LUTHRU or LUCLUS, to sit along the z direction, and the second largest axis into the xz plane. For LUCLUS it can be further specified to +z axis and xz plane with x > 0, respectively. Requires that one of these routines has been called before. = 32 : mainly intended for LUSPHE and LUTHRU, this gives a further alignment of the event, in addition to the one implied by =31. The "slim" jet, defined as the side (z > 0 or z < 0) with the smallest summed p_T over square root of number of particles, is rotated into the +z hemisphere. In the opposite hemisphere (now z < 0), the side of x > 0 and x < 0 which has the largest summed p_z absolute is rotated into the z < 0, x > 0 quadrant. Requires that LUSPHE or LUTHRU has been called before. SUBROUTINE LULIST(MLIST) Purpose: to list an event, jet or particle data, or current parameter values. MLIST : determines what is to be listed. = 0 : writes a header with program version number and last date of change; is mostly for internal use. = 1 : gives a simple list of current event record, in an 80 column format suitable for viewing directly on the computer terminal. For each entry, the following information is given: the entry number I, the parton/particle name (see below), the status code KS (K(I,1)), the flavour code KF (K(I,2)), the line number of the mother (K(I,3)), and the three-momentum, energy and mass (P(I,1) - P(I,5)). If MSTU(3) is nonzero, lines immediately after the event record proper are also listed. A final line contains information on total charge, momentum, energy and invariant mass. The particle name is given by a call to the routine LUNAME. For an entry which has decayed/fragmented (KS = 11 - 20), this particle name is given within parantheses. Similarly, a documentation line (KS = 21 - 30) has the name enclosed in expression signs (!...!) and an event/jet axis information line the name within inequality signs (<...>). If the last character of the name is a ?, it is a signal that the complete name has been truncated to fit in, and can therefore not be trusted; this is very rare. For partons which have been arranged along strings (KS = 1, 2, 11 or 12), the end of the parton name column contains information about the colour string arrangement: a A for the first entry of a string, an I for all intermediate ones, and a V for the final one (a poor man's vertical rendering of the doublesided arrow <---->). It is possible to insert lines just consisting of sequences of ====== to separate different sections of the event record, see MSTU(70) - MSTU(80). = 2 : gives a more extensive list of the current event record, in a 132 column format, suitable for printers or workstations. For each entry, the following information is given: the entry number I, the parton/particle name (with padding as described for MLIST = 1), the status code KS (K(I,1)), the flavour code KF (K(I,2)), the line number of the mother (K(I,3)), the decay product/colour flow pointers (K(I,4), K(I,5)), and the three-momentum, energy and mass (P(I,1) - P(I,5)). If MSTU(3) is nonzero, lines immediately after the event record proper are also listed. A final line contains information on total charge, momentum, energy and invariant mass. Lines with only ====== may be inserted as for MLIST(1). = 3 : gives the same basic listing as = 2, but with an additional line for each entry containing information on production vertex position and time (V(I,1) - V(I,4)) and, for unstable particles, invariant lifetime (V(I,5)). = 11 : provides a simple list of all parton/particle codes defined in the program, with KF code and corresponding particle name. The list is grouped by particle kind, and only within each group in ascending order. = 12 : provides a list of all parton/particle and decay data used in the program. Each parton/particle code is represented by one line containing KF flavour code, KC compressed code, particle name, antiparticle name (where appropriate), electrical and colour charge (stored in KCHG), mass, resonance width and maximum broadening, average invariant lifetime (in PMAS) and whether the particle is considered stable or not (in MDCY). Immediately after a particle, each decay channel gets one line, containing decay channel number (IDC read from MDCY), on/off switch for the channel, matrix element type (MDME), branching ratio (BRAT), and decay products (KFDP). The MSTU(14) flag can be used to set the maximum flavour for which particles are listed, with the default (= 0) corresponding to separately defined ones (KC > 100 if KF > 0). In order to keep the size down, decay modes of heavy hadrons collectively defined are never listed; these have KC codes 84 - 88, where the relevant information may be found. = 13 : gives a list of current parameter values for MSTU, PARU, MSTJ and PARJ, and the first 200 entries of PARF. This is useful to keep check of which default values were changed in a given run. SUBROUTINE LUUPDA(MUPDA,LFN) Purpose: to give the user the ability to update particle data, or to keep several versions of modified particle data for special purposes (e.g. charm studies). MUPDA : gives the type of action to be taken. = 1 : write a table of particle data, that the user then can edit at leisure. For ordinary listing of decay data, LULIST(12) should be used, but that listing could not be read back in by the program. For each compressed flavour code KC = 1 - 500, one line is written containing KC (I5), the basic particle name (pieces a and b in section 2.1) (2X,A8) in CHAF, the electric (I3), colour charge (I3) and particle/antiparticle distinction (I3) codes in KCHG, the mass (F12.5), the mass width (F12.5), maximum broadening (F12.5) and average invariant lifetime (2X,F12.5) in PMAS, and the on/off decay switch (I3) in MDCY(KC,1). After a KC line follows one line for each possible decay channel, containing the MDME codes (5X,2I5), the branching ratio (5X,F12.5) in BRAT, and the KFDP codes for the decay products (5I8), with trailing 0:s if the number of decay products is smaller than 5. = 2 : read in particle data, as written with =1 and thereafter edited by the user, and use this data subsequently in the current run. Reading is done with fixed format, which means that the user has to preserve the format codes described for =1 during the editing. A number of checks will be made to see if input looks reasonable, with warnings if not. If some decay channel is said not to conserve charge, it should be taken seriously. Warnings that decay is kinematically unallowed need not be as serious, since that particular decay mode may not be switched on unless the particle mass is increased. = 3 : write current particle data as data lines, which can be edited into BLOCK DATA LUDATA for a permanent replacement of the particle data. This option should never be used by the ordinary user. LFN : the file number which the data should be written to or read from. The user must see to it that this file is properly opened for read or write (since the definition of file names is machine dependent). FUNCTION KLU(I,J) Purpose: to provide various integer-valued event data. Note that many of the options available (in particular I > 0, J >= 14) which refer to event history will not work after a LUEDIT call. I=0, J= : properties referring to the complete event. = 1 : N, total number of lines in event record. = 2 : total number of partons/particles remaining after fragmentation and decay. = 6 : three times the total charge of remaining (stable) partons and particles. I>0, J= : properties referring to the entry in line no. I of the event record. = 1 - 5 : K(I,1) - K(I,5), i.e. parton/particle status KS, flavour code KF and origin/decay product/colour flow information. = 6 : three times parton/particle charge. = 7 : 1 for a remaining entry, 0 for a decayed/fragmented/ documentation entry. = 8 : KF code (K(I,2)) for a remaining entry, 0 for a decayed/fragmented/documentation entry. = 9 : KF code (K(I,2)) for a parton (i.e. not colour neutral entry), 0 for a particle. = 10 : KF code (K(I,2)) for a particle (i.e. colour neutral entry), 0 for a parton. = 11 : compressed flavour code KC. = 12 : colour information code, i.e. 0 for colour neutral, 1 for colour triplet, -1 for antitriplet and 2 for octet. = 13 : flavour of "heaviest" quark or antiquark (i.e. with largest code) in hadron or diquark (including sign for antiquark), 0 else. = 14 : generation number. Beam particles or virtual exchange particles are generation 0, original jets/particles generation 1 and then 1 is added for each step in the fragmentation/decay chain. = 15 : line number of ancestor, i.e. predecessor in first generation (generation 0 entries are disregarded). = 16 : rank of a hadron in the jet it belongs to. Rank denotes the ordering in flavour space, with hadrons containing the original flavour of the jet having rank 1, increasing by 1 for each step away in flavour ordering. All decay products inherit the rank of their parent. Whereas the meaning of a first-rank hadron in a quark jet is always well-defined, the definition of higher ranks is only meaningful for independently fragmenting quark jets. In other cases, rank refers to the ordering in the actual simulation, which may be of little interest. = 17 : generation number after a collapse of a jet system into one particle, with 0 for an entry not coming from a collapse, and -1 for entry with unknown history. A particle formed in a collapse is generation 1, and then one is added in each decay step. = 18 : number of decay/fragmentation products (only defined in a collective sense for fragmentation). = 19 : origin of colour for showering parton, 0 else. = 20 : origin of anticolour for showering parton, 0 else. = 21 : position of colour daughter for showering parton, 0 else. = 22 : position of anticolour daughter for showering parton, 0 else. FUNCTION PLU(I,J) Purpose: to provide various real-valued event data. Note that some of the options available (I > 0, J = 20-25), which are primarily intended for studies of systems in their respective CM frame, requires that a LUEXEC call has been made for the current initial parton/particle configuration, but that the latest LUEXEC call has not been followed by a LUROBO one. I=0, J= : properties referring to the complete event. = 1 - 4 : sum of p_x, p_y, p_z and E, respectively, for the stable remaining entries. = 5 : invariant mass of the stable remaining entries. = 6 : sum of electric charge of the stable remaining entries. I>0, J= : properties referring to the entry in line no. I of the event record. = 1 - 5 : P(I,1) - P(I,5), i.e. normally p_x, p_y, p_z, E and m for jet/particle. = 6 : electric charge q. = 7 : momentum squared p^2 = p_x^2 + p_y^2 + p_z^2. = 8 : momentum p. = 9 : transverse momentum squared p_T^2 = p_x^2 + p_y^2. = 10 : transverse momentum p_T. = 11 : transverse mass squared m_T^2 = m^2 + p_x^2 + p_y^2. = 12 : transverse mass m_T. = 13 - 14 : polar angle theta in radians (between 0 and pi) or degrees, respectively. = 15 - 16 : azimuthal angle phi in radians (between -pi and pi) or degrees, respectively. = 17 : true rapidity y = 0.5*ln((E+p_z)/(E-p_z)). = 18 : rapidity y_pi obtained by assuming that the particle is a pion when calculating the energy E, to be used in the formula above, from the (assumed known) momentum p. = 19 : pseudorapidity eta = 0.5*ln((p+p_z)/(p-p_z)). = 20 : momentum fraction x_p = 2p/W, where W is the total energy of initial jet/particle configuration. = 21 : x_F = 2p_z/W (Feynman-x if system is studied in CM frame). = 22 : x_T = 2p_T/W. = 23 : x_E = 2E/W. = 24 : z_+ = (E+p_z)/W. = 25 : z_- = (E-p_z)/W. ______________________________________________________________________ 2.6. Event Analysis Routines The six routines LUSPHE, LUTHRU, LUCLUS, LUCELL, LUJMAS and LUFOWO give the user the possibility to find some global event shape properties. The routine LUTABU performs a statistical analysis of a number of different quantities like particle content, factorial moments and the energy-energy correlation. Note that, by default, all remaining partons/particles except neutrinos are used in the analysis. The JETSET 6.3 default procedure of also including neutrinos may be obtained with MSTU(41) = 1. Also note that axes determined are stored in LUJETS, but are not proper four-vectors and, as a general rule (with some exceptions), should therefore not be rotated or boosted. SUBROUTINE LUSPHE(SPH,APL) Purpose: to diagonalize the momentum tensor, i.e. find the eigenvalues lambda_1 > lambda_2 > lambda_3, with sum unity, and the corresponding eigenvectors. Momentum power dependence is given by PARU(41); default corresponds to sphericity, PARU(41)=1. gives measures linear in momenta. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value. SPH : 3(lambda_2 + lambda_3)/2, i.e. sphericity (for PARU(41)=2.). = -1. : analysis not performed because event contained less than two particles (or two exactly back-to-back particles, in which case the two transverse directions would be undefined). APL : 3 lambda_3/2, i.e. aplanarity (for PARU(41)=2.). = -1. : as SPH=-1. Remark: the lines N+1 through N+3 (N-2 through N for MSTU(43) = 2) in LUJETS will, after a call, contain the following information: K(N+i,1) = 31; K(N+i,2) = 95; K(N+i,3) : i, the axis number, i=1,2,3; K(N+i,4), K(N+i,5) = 0; P(N+i,1) - P(N+i,3) : the i:th eigenvector, x,y and z components; P(N+i,4) : lambda_i, the i:th eigenvalue; P(N+i,5) = 0; V(N+i,1) - V(N+i,5) = 0. Also, the number of particles used in the analysis is given in MSTU(62). SUBROUTINE LUTHRU(THR,OBL) Purpose: to find the thrust, major and minor axes and corresponding projected momentum quantities, in particular thrust and oblateness. The performance of the program is affected by MSTU(44), MSTU(45), PARU(42) and PARU(48). In particular, PARU(42) gives the momentum dependence, with the default value 1. corresponding to linear dependence. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value. THR : thrust (for PARU(42)=1.). = -1. : analysis not performed because event contained less than two particles. = -2. : remaining space in LUJETS (partly used as working area) not large enough to allow analysis. OBL : oblateness (for PARU(42)=1.). = -1., -2. : as for THR. Remark: the lines N+1 through N+3 (N-2 through N for MSTU(43) = 2) in LUJETS will, after a call, contain the following information: K(N+i,1) = 31; K(N+i,2) = 96; K(N+i,3) : i, the axis number, i=1,2,3; K(N+i,4), K(N+i,5) = 0; P(N+i,1) - P(N+i,3) : the thrust, major and minor axis, respectively, for i = 1, 2 and 3; P(N+i,4) : corresponding thrust, major and minor value; P(N+i,5) = 0; V(N+i,1) - V(N+i,5) = 0. Also, the number of particles used in the analysis is given in MSTU(62). SUBROUTINE LUCLUS(NJET) Purpose: to reconstruct an arbitrary number of jets using a cluster analysis method based on particle momenta. Two different distance measures are available, the traditional one roughly corresponding to relative transverse momentum [Sjo83], and a new one based on the JADE method of Bethke, which roughly corresponds to invariant mass [Bet86] (in both cases with some important modifications). The choice is controlled by MSTU(46). The distance scale d_join, above which two clusters may not be joined, is normally given by PARU(44). In general, d_join may be varied to describe different "jet resolution powers"; the default value, 2.5 GeV, is fairly well suited for e+e- physics at 30 - 40 GeV. With the mass distance, PARU(44) can be used to set the absolute maximum cluster mass, or PARU(45) to set the scaled one, i.e. in y = m^2/W^2, where W^2 is the total invariant mass-squared of the particles being considered. It is possible to continue the cluster search from the configuration already found, with a new higher d_join scale, by selecting MSTU(48) properly. In MSTU(47) one can also require a minimum number of jets to be reconstructed; combined with an artificially large d_join this can be used to reconstruct a predetermined number of jets. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value, whereas assumptions about particle masses is given by MSTU(42). The parameters PARU(43) and PARU(48) regulate more technical details (for events at high energies and large multiplicities, however, the choice of a larger PARU(43) may be necessary to obtain reasonable reconstruction times). NJET : the number of clusters reconstructed. = -1 : analysis not performed because event contained less than MSTU(47) (normally 1) particles, or analysis failed to reconstruct the requested number of jets. = -2 : remaining space in LUJETS (partly used as working area) not large enough to allow analysis. Remark: if the analysis does not fail, further information is found in MSTU(61) - MSTU(63) and PARU(61) - PARU(63). In particular, PARU(61) contains the invariant mass for the system analyzed, i.e. the number used in determining the denominator of y = m^2/W^2. PARU(62) gives the generalized thrust, i.e. the sum of (absolute values of) cluster momenta divided by the sum of particle momenta (roughly the same as multicity). PARU(63) gives the minimum distance d (in p_T or m) between two clusters in the final cluster configuration, 0 in case of only one cluster. Further, the lines N+1 through N+NJET (N-NJET+1 through N for MSTU(43) = 2) in LUJETS will, after a call, contain the following information: K(N+i,1) = 31; K(N+i,2) = 97; K(N+i,3) : i, the jet number, with the jets arranged in falling order of absolute momentum; K(N+i,4) : the number of particles assigned to jet i; K(N+i,5) = 0; P(N+i,1) - P(N+i,5) : momentum, energy and invariant mass of jet i; V(N+i,1) - V(N+i,5) = 0. Also, for a particle which was used in the analysis, K(I,4) = i, where I is the particle number and i the number of the jet it has ben assigned to. Undecayed particles not used then have K(I,4) = 0. An exception is made for lines with K(I,1) = 3 (which anyhow are not normally interesting for cluster search), where the colour flow information stored in K(I,4) is left intact. SUBROUTINE LUCELL(NJET) Purpose: to provide a simpler cluster routine more in line with what is currently used in the study of high-p_T collider events. A detector is assumed to stretch in pseudorapidity between -PARU(51) and +PARU(51) and be segmented in MSTU(51) equally large eta (pseudorapidity) bins and MSTU(52) phi (azimuthal) bins. Transverse energy E_T for undecayed entries are summed up in each bin. For MSTU(53) nonzero, the energy is smeared by calorimetric resolution effects, cell by cell. This is done according to a Gaussian distribution; if MSTU(53) = 1 the standard deviation for the E_T is PARU(55)*sqrt(E_T), if MSTU(53) = 2 the standard deviation for the E is PARU(55)*sqrt(E), E_T and E expressed in GeV. The Gaussian is cut off at 0 and at a factor PARU(56) times the correct E_T or E. All bins with E_T > PARU(52) are taken to be possible initiators of jets, and are tried in falling E_T sequence to check whether the total E_T summed over cells no more distant than PARU(54) in sqrt((Deltaeta)^2 + (Deltaphi)^2) exceeds PARU(53). If so, these cells define one jet, and are removed from further consideration. Contrary to LUCLUS, not all particles need be assigned to jets. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value. NJET : the number of jets reconstructed (may be 0). = -2 : remaining space in LUJETS (partly used as working area) not large enough to allow analysis. Remark: the lines N+1 through N+NJET (N-NJET+1 through N for MSTU(43) = 2) in LUJETS will, after a call, contain the following information: K(N+i,1) = 31; K(N+i,2) = 98; K(N+i,3) : i, the jet number, with the jets arranged in falling order in E_T; K(N+i,4) : the number of particles assigned to jet i; K(N+i,5) = 0; V(N+i,1) - V(N+i,5) = 0; Further, for MSTU(54) = 1 P(N+i,1), P(N+i,2) = position in eta and phi of the center of the jet initiator cell, i.e. geometrical center of jet; P(N+i,3), P(N+i,4) = position in eta and phi of the E_T-weighted center of the jet, i.e. the center of gravity of the jet; P(N+i,5) = sum E_T of the jet; while for MSTU(54) = 2 P(N+i,1) - P(N+i,5) : the jet momentum vector, constructed from the summed E_T and the eta and phi of the E_T-weighted center of the jet as (p_x, p_y, p_z, E, m) = E_T(cos(phi), sin(phi), sinh(eta), cosh(eta), 0); and for MSTU(54) = 3 P(N+i,1) - P(N+i,5) : the jet momentum vector, constructed by adding vectorially the momentum of each cell assigned to the jet, assuming that all the E_T was deposited at the center of the cell, and with the jet mass in P(N+i,5) calculated from the summed E and p as m^2 = E^2 - p_x^2 - p_y^2 - p_z^2. Also, the number of particles used in the analysis is given in MSTU(62), and the number of cells hit in MSTU(63). SUBROUTINE LUJMAS(PMH,PML) Purpose: to reconstruct high and low jet mass of an event, according to the general suggestion of Clavelli and Smilga [Cla79]. Actually, a simplified algorithm is used, wherein a preliminary division of the event into two hemispheres is done transversely to the sphericity axis. Then one particle at a time is reassigned to the other hemisphere if that reduces the sum of the two jet masses-squared (m_H^2 + m_L^2). The procedure is stopped when no further significant change (see PARU(48)) is obtained. Often, the original assignment is retained as it is. Which particles (or partons) used in the analysis is determined by the MSTU(41) value, whereas assumptions about particle masses is given by MSTU(42). PMH : heavy jet mass (in GeV). = -2. : remaining space in LUJETS (partly used as working area) not large enough to allow analysis. PML : light jet mass (in GeV). = -2. : remaining space in LUJETS (partly used as working area) not large enough to allow analysis. Remark: After a successful call, MSTU(62) contains the number of particles used in the analysis, and PARU(61) the invariant mass of the system analyzed. The latter number is helpful in constructing scaled jet masses. SUBROUTINE LUFOWO(H10,H20,H30,H40) Purpose: to do an event analysis in terms of the Fox-Wolfram moments [Fox79]. The moments H_i are normalized to the lowest one, H_0. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value. H10 : H_1/H_0. Is =0 if momentum is balanced. H20 : H_2/H_0. H30 : H_3/H_0. H40 : H_4/H_0. Remark: the number of particles used in the analysis is given in MSTU(62). SUBROUTINE LUTABU(MTABU) Purpose: to provide a number of event analysis options which can be be used on each new event, with accumulated statistics to be written out on request. When errors are quoted, these refer to the uncertainty in the average value for the event sample as a whole, rather than to the spread of the individual events, i.e. errors decrease like one over the square root of the number of events analyzed. For a correct use of LUTABU, it is not permissible to freely mix generation and analysis of different classes of events, since only one set of statistics counters exists. A single run may still contain sequential "subruns", between which statistics is reset. Whenever an event is analyzed, the number of particles/partons used is given in MSTU(62). MTABU : determines which action is to be taken. Generally, a last digit equal to 0 indicates that the statistics counters for this option is to be reset; since the counters are reset (by DATA statements) at the beginning of a run, this is not used normally. Last digit 1 leads to an analysis of current event with respect to the desired properties. Note that the resulting action may depend on how the event generated has been rotated, boosted or edited before this call. The statistics accumulated is output in tabular form with last digit 2, while it is dumped in the LUJETS commonblock for last digit 3. The latter option may be useful for interfacing to graphics output. = 10 : statistics on parton multiplicity is reset. = 11 : the parton content of the current event is analyzed, classified according to the flavour content of the hard interaction and the total number of partons. The flavour content is assumed given in MSTU(161) and MSTU(162); these are automatically set e.g. in LUEEVT calls. = 12 : gives a table on parton multiplicity distribution. = 13 : stores the parton multiplicity distribution of events in /LUJETS/, using the following format: N = total number of different channels found; K(I,1) = 32; K(I,2) = 99; K(I,3), K(I,4) = the two flavours of the flavour content; K(I,5) = total number of events found with flavour content of K(I,3) and K(I,4); P(I,1) - P(I,5) = relative probability to find given flavour content and a total of 1, 2, 3, 4 or 5 partons, respectively; V(I,1) - V(I,5) = relative probability to find given flavour content and a total of 6 - 7, 8 - 10, 11 - 15, 16 - 25 or above 25 partons, respectively. In addition, MSTU(3) = 1 and K(N+1,1) = 32; K(N+1,2) = 99; K(N+1,5) = number of events analyzed. = 20 : statistics on particle content is reset. = 21 : the particle/parton content of the current event is analyzed, also for particles which have subsequently decayed and partons which have fragmented (unless this has been made impossible by a preceding LUEDIT call). Particles are subdivided into primary and secondary ones, the main principle being that primary particles are those produced in the fragmentation of a string, while secondary come from decay of other particles. Since particles (top, say), may decay into partons, the distinction is not always unique. = 22 : gives a table of particle content in events. = 23 : stores particle content in events in /LUJETS/, using the following format: N = number of different particle species found; K(I,1) = 32; K(I,2) = 99; K(I,3) = particle KF code; K(I,5) = total number of particles and antiparticles of this species; P(I,1) = average number of primary particles per event; P(I,2) = average number of secondary particles per event; P(I,3) = average number of primary antiparticles per event; P(I,4) = average number of secondary antiparticles per event; P(I,5) = average total number of particles or antiparticles per event. In addition, MSTU(3) = 1 and K(N+1,1) = 32; K(N+1,2) = 99; K(N+1,5) = number of events analyzed; P(N+1,1) = average primary multiplicity per event; P(N+1,2) = average final multiplicity per event; P(N+1,3) = average charged multiplicity per event. = 30 : statistics on factorial moments is reset. = 31 : analyzes the factorial moments of the multiplicity distribution in different bins of rapidity and azimuth. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value. The selection between usage of true rapidity, pion rapidity or pseudorapidity is regulated by MSTU(42). The z axis is assumed to be event axis; if this is not desirable find an event axis e.g. with LUSPHE or LUTHRU and use LUEDIT(31). Maximum (pion, pseudo) rapidity, which sets the limit for the rapidity plateau or the experimental acceptance, is given by PARU(57). = 32 : prints a table of the first four factorial moments for various bins of pseudorapidity and azimuth. The moments are properly normalized so that they would be unity (up to statistical fluctuations) for uniform and uncorrelated particle production according to Poissonian statistics, but increasing for decreasing bin size in case of "intermittent" behaviour [Bia86]. The error on the average value is based on the actual statistical sample (i.e. does not use any assumptions on the distribution to relate errors to the average values of higher moments), and decrease as 1/sqrt(number-of-events). Note that for small bin sizes, where the average multiplicity is small and the factorial moment therefore only very rarely is nonvanishing, moment values may fluctuate wildly and the errors given may be too low. = 33 : stores the factorial moments in /LUJETS/, using the format: N = 30, with I = 1 - 10 corresponding to results for slicing the rapidity range in 2**(1-I) bins, I = 11 - 20 to slicing the azimuth in 2**(11-I) bins, and I = 21 - 30 to slicing both rapidity and azimuth, each in 2**(21-I) bins; K(I,1) = 32; K(I,2) = 99; K(I,3) = number of bins in rapidity; K(I,4) = number of bins in azimuth; P(I,1) = rapidity bin size; P(I,2) - P(I,5) = - i.e. mean of second, third, fourth and fifth factorial moment; V(I,1) = azimuthal bin size; V(I,2) - V(I,5) = statistical errors on - . In addition, MSTU(3) = 1 and K(31,1) = 32; K(31,2) = 99; K(31,5) = number of events analyzed. = 40 : statistics on energy-energy correlation is reset. = 41 : the energy-energy correlation (EEC) of the current event is analyzed. Which particles (or partons) are used in the analysis is determined by the MSTU(41) value. Events are assumed given in their CM frame. The weight assigned to a pair i and j is 2*E_i*E_j/W^2, where W is the sum of energies of all analyzed particles in the event. Energies are determined from the momenta of particles, with mass determined according to the MSTU(42) value. Statistics is accumulated for the relative angle theta_ij, ranging between 0 and 180 degrees, subdivided into 50 bins. = 42 : prints a table of the energy-energy correlation (EEC) and its asymmetry (EECA), with errors. The definition of errors is not unique. In our approach each event is viewed as one observation, i.e. an EEC and EECA distribution is obtained by summing over all particle pairs of an event, and then the average and spread of this event-distribution is calculated in the standard fashion, with the quoted error containing a factor of 1/sqrt(number-of-events). It could have been possible to view each single particle pair as one observation, which would have given somewhat lower errors, but then one would also be forced to do a complicated correction procedure to account for the pairs in an event not being uncorrelated (two hard jets separated by a given angle typically corresponds to several pairs at about that angle). Note, however, that in our approach the error-squared on an EECA bin is smaller than the sum of error-squares of the corresponding EEC bins (as it should be). Also note that it is not possible to combine the errors of two nearby bins by hand from the information given, since nearby bins are correlated (again a trivial consequence of the presence of jets). = 43 : stores the EEC and EECA in /LUJETS/, using the format: N = 25; K(I,1) = 32; K(I,2) = 99; P(I,1) = EEC for angles between I-1 and I times 3.6 degrees; P(I,2) = EEC for angles between 50-I and 51-I times 3.6 degrees; P(I,3) = EECA for angles between I-1 and I times 3.6 degrees; P(I,4), P(I,5) : lower and upper edge of angular range of bin I, expressed in radians; V(I,1) - V(I,3) : errors on the EEC and EECA values stored in P(I,1) - P(I,3) (see =42 for comments); V(I,4), V(I,5) : lower and upper edge of angular range of bin I, expressed in degrees. In addition, MSTU(3) = 1 and K(26,1) = 32; K(26,2) = 99; K(26,5) = number of events analyzed. = 50 : statistics on complete final states is reset. = 51 : analyzes the particle content of the final state of the current event record. During the course of the run, statistics is thus accumulated on how often different final states appear. Only final states with up to 8 particles are analyzed, and there is only reserved space for up to 200 different final states. Most high energy events have multiplicities far above 8, so the main use for this tool is to study the effective branching ratios obtained with a given decay model for e.g. charm or bottom hadrons. Then LU1ENT may be used to generate one decaying particle at a time, with a subsequent analysis by LUTABU. Depending on at what level this studied is to be carried out, some particle decays may be switched off, like pi0. = 52 : gives a list of the (at most 200) channels with up to 8 particles in the final state, with their relative branching ratio. The ordering is according to multiplicity, and within each multiplicity according to an ascending order of KF codes. The KF codes of the particles belonging to a given channel are given in descending order. = 53 : stores the final states and branching ratios found in /LUJETS/, using the format: N = number of different explicit final states found (at most 200); K(I,1) = 32; K(I,2) = 99; K(I,5) = multiplicity of given final state, a number between 1 and 8; P(I,1) - P(I,5), V(I,1) - V(I,3) : the KF codes of the up to 8 particles of the given final state, converted to real numbers, with trailing zeroes for positions not used;. V(I,5) : effective branching ratio for the given final state. In addition, MSTU(3) = 1 and K(N+1,1) = 32; K(N+1,2) = 99; K(N+1,5) = number of events analyzed; V(N+1,5) = summed branching ratio for finals states not given above, either because they contained more than 8 particles or because all 200 channels have been used up. ______________________________________________________________________ 2.7. The General Switches and Parameters The commonblock LUDAT1 is, next to LUJETS, the one a user is most likely to access. Here he may control in detail what the program is to do, if the default mode of operation is not satisfactory. COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200) Purpose: to give access to a number of status codes and parameters which regulate the performance of the program as a whole. Here MSTU and PARU are related to utility functions, as well as a few parameters of the standard model, while MSTJ and PARJ affect the underlying physics assumptions. Some of the variables in LUDAT1 deal specifically with e+e- physics, and are described in section 3.3. MSTU(1),MSTU(2) : (D=0,0) can be used to replace the ordinary lower and upper limits (normally 1 and N) for the action of LUROBO, and most LUEDIT and LULIST calls. Are reset to 0 in a LUEXEC call. MSTU(3) : (D=0) number of lines with extra information added after line N. Is reset to 0 in a LUEXEC call, or in an LUEDIT call when particles are removed. MSTU(4) : (D=4000) number of lines available in the commonblock LUJETS. Should always be changed if the dimensions of the K and P arrays are changed by the user, but should otherwise never be touched. Maximum allowed value is 10000, unless MSTU(5) is also changed. MSTU(5) : (D=10000) is used in building up the special colour flow information stored in K(I,4) and K(I,5) for K(I,3) = 3, 13 or 14. The generic form for j = 4 or 5 is K(I,j) = 2*MSTU(5)**2*MCFR + MSTU(5)**2*MCTO + MSTU(5)*ICTO + ICFR with notation as in section 2.2. One should always have MSTU(5) >= MSTU(4). On a 32 bit machine, vaules MSTU(5) > 20000 may lead to overflow problems, and should be avoided. MSTU(6) : (D=500) number of KC codes available in the KCHG, PMAS, MDCY, and CHAF arrays; should be changed if these dimensions are changed. MSTU(7) : (D=2000) number of decay channels available in the MDME, BRAT and KFDP arrays; should be changed if these dimensions are changed. MSTU(10) : (D=2) use of parton/particle masses in filling routines (LU1ENT, LU2ENT, LU3ENT, LU4ENT). = 0 : assume the mass to be zero. = 1 : keep the mass value stored in P(I,5), whatever it is. (This may be used e.g. to describe kinematics with off-mass-shell partons). = 2 : find masses according to mass tables as usual. MSTU(11) : (D=6) file number to which all program output is directed. It is the responsibility of the user to see to it that the corresponding file is also opened for output. MSTU(12) : (D=1) writing of header (version number and last date of change) on output file. = 0 : not done. = 1 : header is written at first occasion, at which time MSTU(12) is set =0. MSTU(13) : (D=1) writing of information on variable values changed by LUGIVE calls. = 0 : no information is provided. = 1 : information is written to standard output. MSTU(14) : (D=0) if nonzero, this gives the maximum flavour for which a LULIST(12) call will give particle data on possible hadrons. With MSTU(14)=5 only known hadrons, i.e. up to bottom, are listed. If =0, only separately specified particles are listed (i.e. either KF <= 100 or else both KF > 100 and KC > 100). MSTU(15) : (D=1) selection for characters used in particle names to denote an antiparticle; appear in LULIST listings or other LUNAME applications. = 1 : the tilde character '~'. = 2 : the characters 'bar'. MSTU(16) : (D=1) choice of mother pointers for the particles produced by a fragmenting parton system. = 1 : all primary particles of a system point to a line with KF = 92 or 93, for string or independent fragmentation, respectively, or to a line with KF = 91 if a jet system has so small a mass that it is forced to decay into one or two particles. The two (or more) shower initiators of a showering parton system point to a line with KF = 94. The entries with KF = 91 - 94 in their turn point back to the predecessor partons, so that the KF = 91 - 94 entries form a part of the event history proper. = 2 : although the lines with KF = 91 - 94 are present, and contain the proper mother and daughter pointers, they are not part of the event history proper, in that particles produced in string fragmentation point directly to either of the two endpoint partons of the string (depending on the side they were generated from), particles produced in independent fragmentation point to the respective parton they were generated from, particles in small mass systems point to either endpoint parton, and shower initiators point to the original on-mass-shell counterparts. Also the daugher pointers bypass the KF = 91 - 94 entries. In independent fragmentation, a parton need not produce any particles at all, and then have daughter pointers 0. Note : MSTU(16) should not be changed between the generation of an event and the translation of this event record with a LUHEPC call, since this may give an erroneous translation of the event history. MSTU(17) : (D=0) storage option for MSTU(90) and associated information on z values for heavy flavour production. = 0 : MSTU(90) is reset to zero at each LUEXEC call. This is the appropriate course if LUEXEC is only called once per event, as is normally the case when the user does not himself call LUEXEC. = 1 : user has to reset MSTU(90) to zero himself before each new event. This is the appropriate course if several LUEXEC calls may appear for one event, i.e. if the user calls LUEXEC directly. MSTU(19) : (D=0) advisory warning for unphysical flavour setups in LU2ENT, LU3ENT or LU4ENT calls. = 0 : yes. = 1 : no; MSTU(19) is reset to 0 in such a call. MSTU(21) : (D=2) check on possible errors during program execution. Obviously no guarantee is given that all errors will be caught, but some of the most trivial user-caused errors may be found. = 0 : errors do not cause any immediate action, rather the program will try to cope, which may mean e.g. that it runs into an infinite loop. = 1 : parton/particle configurations are checked for possible errors. In case of problem, an exit is made from the misbehaving subprogram, but the generation of the event is continued from there on. For the first MSTU(22) errors a a message is printed; after that no messages appear. = 2 : parton/particle configurations are checked for possible errors. In case of problem, an exit is made from the misbehaving subprogram, and subsequently from LUEXEC. The user may then choose to correct the error, and continue the execution by another LUEXEC call. For the first MSTU(22) errors a message is printed, after that the last event is printed and execution is stopped. MSTU(22) : (D=10) max number of errors that are printed. MSTU(23) : (I) count of number of errors experienced to date. MSTU(24) : (R) type of latest error experienced; reason that event was not generated in full. Is reset at each LUEXEC call. = 0 : no error experienced. = 1 : have reached end of or are writing outside LUJETS memory. = 2 : unknown flavour code or unphysical combination of codes; may also be caused by erroneous string connection information. = 3 : energy or mass too small or unphysical kinematical variable setup. = 4 : program is caught in an infinite loop. = 5 : momentum, energy or charge was not conserved (even allowing for machine precision errors, see PARU(11)); is evaluated only after event has been generated in full, and does not apply when independent fragmentation without momentum conservation was used. = 6 : error call from outside the fragmentation/decay package (e.g. the e+e- routines). = 7 : inconsistent particle data input in LUUPDA (MUPDA = 2) or other LUUPDA-related problem. = 8 : problems in more peripheral service routines. = 9 : various other problems. MSTU(25) : (D=1) printing of warning messages. = 0 : no warnings are written. = 1 : first MSTU(26) warnings are printed, thereafter no warnings appear. MSTU(26) : (D=10) max number of warnings that are printed. MSTU(27) : (I) count of number of warnings experienced to date. MSTU(28) : (R) type of latest warning given, with codes paralleling those for MSTU(24), but of a less serious nature. MSTU(31) : (I) number of LUEXEC calls in present run. MSTU(32) : (I) number of entries stored with LUEDIT(-1) call. MSTU(33) : (I) if set 1 before a LUDBRB call, the V vectors (in the particle range to be rotated/boosted) are set 0 before the rotation/boost. MSTU(33) is set back to 0 in the LUDBRB call. Is inactive in a LUROBO call. MSTU(41) : (D=2) partons/particles used in the event analysis routines LUSPHE, LUTHRU, LUCLUS, LUCELL, LUJMAS, LUFOWO and LUTABU (LUTABU(11) excepted) (cf. LUEDIT). = 1 : all partons/particles that have not fragmented/decayed. = 2 : ditto, with the exception of neutrinos and unknown particles. = 3 : only charged, stable particles, plus any partons still not fragmented. MSTU(42) : (D=2) assumed particle masses, used in calculating energies E^2 = p^2 + m^2, as subsequently used in LUCLUS, LUJMAS and LUTABU (in the latter also for pseudorapidity/pion rapidity/true rapidity selection). = 0 : all particles are assumed massless. = 1 : all particles, except the photon, are assumed to have the charged pion mass. = 2 : the true masses are used. MSTU(43) : (D=1) storing of event analysis information (mainly jet axes), in LUSPHE, LUTHRU, LUCLUS and LUCELL. = 1 : stored after the event proper, in positions N+1 through N+MSTU(3). If several of the routines are used in succession, all but the latest information is overwritten. = 2 : stored with the event proper, i.e. at the end of the event listing, with N updated accordingly. If several of the routines are used in succession, all the axes determined are available. MSTU(44) : (D=4) is the number of the fastest (i.e. with largest momentum) particles used to construct the (at most) 10 most promising starting configurations for the thrust axis determination. MSTU(45) : (D=2) is the number of different starting configurations above, which have to converge to the same (best) value before this is accepted as the correct thrust axis. MSTU(46) : (D=1) distance measure used for the joining of clusters in LUCLUS. = 1 : (approximately) relative transverse momentum. Anytime two clusters have been joined, particles are reassigned to the cluster they now are closest to. The distance cutoff d_join is stored in PARU(44). = 2 : distance measure as in =1, but particles are never reassigned to new jets. = 3 : (approximately) total invariant mass. Particles may never be reassigned between clusters. The distance cutoff m_min is stored in PARU(44). = 4 : as =3, but a scaled distance y = m^2/W^2 is used instead of m. The distance cutoff y_min is stored in PARU(45). MSTU(47) : (D=1) the minimum number of clusters to be reconstructed by LUCLUS. MSTU(48) : (D=0) mode of operation of the LUCLUS routine. = 0 : the cluster search is started from scratch. = 1 : the clusters obtained in a previous cluster search on the same event (with MSTU(48)=0) are to be taken as the starting point for subsequent cluster joining. For this call to have any effect, the joining scale in PARU(44) or PARU(45) must have been changed. If the event record has been modified after the last LUCLUS call, or if any other cluster search parameter setting has been changed, the subsequent result is unpredictable. MSTU(51) : (D=25) number of pseudorapidity bins that the range between -PARU(51) and +PARU(51) is divided into to define cell size for LUCELL. MSTU(52) : (D=24) number of azimuthal bins, used to define the cell size for LUCELL. MSTU(53) : (D=0) smearing of correct energy, imposed cell-by-cell in LUCELL, to simulate calorimeter resolution effects. = 0 : no smearing. = 1 : the transverse energy in a cell, E_T, is smeared according to a Gaussian distribution with standard deviation PARU(55)*sqrt(E_T), where E_T is given in GeV. The Gaussian is cut off so that 0 < E_T_smeared < PARU(56)*E_T_true. = 2 : as =1, but it is the energy E rather than the transverse energy E_T that is smeared. MSTU(54) : (D=1) form for presentation of information about reconstructed clusters in LUCELL, as stored in LUJETS according to the MSTU(43) value. = 1 : the P vector in each line contains eta and phi for the geometric origin of the jet, eta and phi for the weighted center of the jet, and jet E_T, respectively. = 2 : the P vector in each line contains a massless four-vector giving the direction of the jet, obtained as (p_x,p_y,p_z,E,m) = E_T(cos(phi),sin(phi),sinh(eta),cosh(eta),0), where eta and phi give the weighted center of a jet and E_T its transverse energy. = 3 : the P vector in each line contains a massive four-vector, obtained by adding the massless four-vectors of all cells that form part of the jet, and calculating the jet mass from m^2 = E^2-p_x^2-p_y^2-p_z^2. For each cell, the total E_T is summed up, and then translated into a massless four-vector assuming that all the E_T was deposited in the center of the cell. MSTU(61) : (I) first entry for storage of event analysis information in last event analyzed with LUSPHE, LUTHRU, LUCLUS or LUCELL. MSTU(62) : (R) number of particles/partons used in the last event analysis with LUSPHE, LUTHRU, LUCLUS, LUCELL, LUJMAS, LUFOWO or LUTABU. MSTU(63) : (R) in a LUCLUS call, the number of preclusters constructed in order to speed up analysis (should be equal to MSTU(62) if PARU(43) = 0.). In a LUCELL call, the number of cells hit. MSTU(70) : (D=0) the number of lines consisting only of equal signs (======) that are inserted in the event listing obtained with LULIST(1), LULIST(2) or LULIST(3), so as to distinguish different sections of the event record on output. At most 10 such lines can be inserted; see MSTU(71) - MSTU(80). Is reset at LUEDIT calls with arguments 0 - 5. MSTU(71) - MSTU(80) : line numbers below which lines consisting only of equal signs (======) are inserted in event listings. Only the first MSTU(70) of the 10 allowed positions are enabled. MSTU(90) : number of heavy flavour hadrons (i.e. hadrons containing charm or heavier flavours) produced in current event, for which the positions in the event record are stored in MSTU(91) - MSTU(98) and the z values in the fragmentation in PARU(91) - PARU(98). At most eight values will be stored (normally this is no problem). No z values can be stored for those heavy hadrons produced when a string has so small mass that it collapses to one or two particles, nor for those produced as one of the final two particles in the fragmentation of a string. If MSTU(17) = 1, MSTU(90) should be reset to zero by the user before each new event, else this is done automatically. MSTU(91) - MSTU(98) : the first MSTU(90) positions will be filled with the line numbers of the heavy flavour hadrons produced in the current event. See MSTU(90) for additional comments. Note that the information is corrupted by calls to LUEDIT with options 0 - 5 and 21 - 23; calls with options 11 - 15 work, however. MSTU(101) : (D=1) procedure for alpha_em evaluation in the ULALEM function. = 0 : alpha_em is taken fixed at the value PARU(101). = 1 : alpha_em is running with the Q^2 scale, taking into account corrections from fermion loops (e, mu, tau, d, u, s, c, b). MSTU(111) : (D=1) order of alpha_strong evaluation in the ULALPS function. Is overwritten in LUEEVT, LUONIA or PYINIT (in the PYTHIA program) calls with the value desired for the process under study. = 0 : alpha_strong is fixed at the value PARU(111). = 1 : first order running alpha_strong is used. = 2 : second order running alpha_strong is used. MSTU(112) : (D=5) the nominal number of flavours assumed in the alpha_strong expression, with respect to which Lambda is defined. MSTU(113) : (D=3) minimum number of flavours that may be assumed in alpha_strong expression, see MSTU(112). MSTU(114) : (D=5) maximum number of flavours that may be assumed in alpha_strong expression, see MSTU(112). MSTU(115) : (D=0) treatment of alpha_strong singularity for Q^2 -> 0. = 0 : allow it to diverge like 1/ln(Q^2/Lambda^2). = 1 : soften the divergence to 1/ln(1 + Q^2/Lambda^2). = 2 : freeze Q^2 evolution below PARU(114), i.e. the effective argument is max(Q^2, PARU(114)). MSTU(118) : (I) number of flavours n_f found and used in latest ULALPS call. MSTU(161), MSTU(162) : hard flavours involved in current event, as used in an analysis with LUTABU(11). Either or both may be set 0, to indicate the presence of one or none hard flavours in event. Is normally set by high-level routines, like LUEEVT, but can also set by user. MSTU(181) : (R) JETSET version number. MSTU(182) : (R) JETSET subversion number. MSTU(183) : (R) last year of change for JETSET. MSTU(184) : (R) last month of change for JETSET. MSTU(185) : (R) last day of change for JETSET. PARU(1) : (R) pi = 3.1415927. PARU(2) : (R) 2*pi = 6.2831854. PARU(3) : (D=0.1973) conversion factor for GeV^-1 -> fm or fm^-1 -> GeV. PARU(4) : (D=5.068) conversion factor for fm -> GeV^-1 or GeV -> fm^-1. PARU(5) : (D=0.3894) conversion factor for GeV^-2 -> mb or mb^-1 -> GeV^2. PARU(6) : (D=2.568) conversion factor for mb -> GeV^-2 or GeV^2 -> mb^-1. PARU(11) : (D=0.001) relative error, i.e. nonconservation of momentum and energy divided by total energy, that may be attributable to machine precision problems before a physics error is suspected (see MSTU(24) = 5). PARU(12) : (D=0.09 GeV^2) effective cutoff in mass-square, below which partons may be recombined to simplify (machine precision limited) kinematics of string fragmentation. PARU(13) : (D=0.01) effective angular cutoff in radians for recombination of partons, used in conjunction with PARU(12). PARU(21) : (I) contains the total energy W of all first generation jets/particles after a LUEXEC call; to be used by the PLU function for I>0, J=20-25. PARU(41) : (D=2.) power of momentum-dependence in LUSPHE, default corresponds to sphericity, =1. to linear event measures. PARU(42) : (D=1.) power of momentum-dependence in LUTHRU, default corresponds to thrust. PARU(43) : (D=0.25 GeV) maximum distance d_init allowed in LUCLUS when forming starting clusters used to speed up reconstruction. The meaning of the parameter is in p_T for MSTU(46) <= 2 and in m for MSTU(46) >= 3. If =0., no preclustering is obtained. If chosen too large, more joining may be generated at this stage than is desirable. The main application is at high energies, where some speedup is imperative, and the small details are not so important anyway. PARU(44) : (D=2.5 GeV) maximum distance d_join, below which it is allowed to join two clusters into one in LUCLUS. Is used for MSTU(46) <= 3, i.e. both for p_T and mass distance measure. PARU(45) : (D=0.05) maximum distance y_join = m^2/W^2, below which it is allowed to join two clusters into one in LUCLUS for MSTU(46) = 4. PARU(48) : (D=0.0001) convergence criterion for thrust (in LUTHRU) or generalized thrust (in LUCLUS), or relative change of m_H^2 + m_L^2 (in LUJMAS), i.e. when the value changes by less than this amount between two iterations the process is stopped. PARU(51) : (D=2.5) defines maximum absolute pseudorapidity used for detector assumed in LUCELL. PARU(52) : (D=1.5 GeV) gives minimum E_T for a cell to be considered as a potential jet initiator by LUCELL. PARU(53) : (D=7.0 GeV) gives minimum summed E_T for a collection of cells to be accepted as a jet. PARU(54) : (D=1.) gives the maximum distance in sqrt((Deltaeta)^2 + (Deltaphi)^2) from cell initiator when grouping cells to check whether they qualify as a jet. PARU(55) : (D=0.5) the calorimeter cell resolution assumed when smearing the transverse energy (or energy) in LUCELL (see MSTU(53)) is taken to be PARU(55)*sqrt(E_T) (or PARU(55)*sqrt(E)) for E_T (or E) in GeV. PARU(56) : (D=2.) maximum factor of upward fluctuation in transverse energy or energy in a given cell when calorimeter resolution is included in LUCELL (see MSTU(53)). PARU(57) : (D=3.2) maximum rapidity (or pseudorapidity or pion rapidity, see MSTU(42)) used in the factorial moments analysis in LUTABU. PARU(61) : (I) invariant mass W of a system analyzed with LUCLUS or LUJMAS, with energies calculated according to the MSTU(42) value. PARU(62) : (R) the generalized thrust obtained after a successful LUCLUS call, i.e. ratio of summed cluster momenta and summed particle momenta. PARU(63) : (R) the minimum distance d between two clusters in the final cluster configuration after a successful LUCLUS call; is 0 if only one cluster left. PARU(91) - PARU(98) : the first MSTU(90) positions will be filled with the fragmentation z values used internally in the generation of heavy flavour hadrons - how these are translated into the actual energies and momenta of the observed hadrons is a complicated function of the string configuration. The particle with z value stored in PARU(i) is to be found in line MSTU(i) of the event record. See MSTU(90) and MSTU(91) - MSTU(98) for additional comments. PARU(101) : (D=0.00729735) alpha_em, the electromagnetic fine structure constant at vanishing momentum transfer. PARU(102) : (D=0.229) sin^2(theta_W), the weak mixing angle of the standard electroweak model. PARU(108) : (I) the alpha_em value obtained in the latest call to the ULALEM function. PARU(111) : (D=0.20) fix alpha-strong value assumed in ULALPS when MSTU(111) = 0 (and also in parton showers when alpha_strong is assumed fix there). PARU(112) : (D=0.25 GeV) Lambda used in running alpha_strong expression in ULALPS. Like MSTU(111), this value is overwritten by the calling physics routines, and is therefore purely nominal. PARU(113) : (D=1.) the flavour thresholds, for the effective number of flavours n_f to use in the alpha_strong expression, are assumed to sit at Q^2 = PARU(113)*m_q^2, where m_q is the quark mass. May be overwritten from the calling physics routine. PARU(114) : (D=4 GeV^2) Q^2 value below which the alpha_strong value is assumed constant for MSTU(115) = 2. PARU(115) : (D=10.) maximum alpha_strong value that ULALPS will ever return; is used as a last resort to avoid singularities. PARU(117) : (I) Lambda value (associated with MSTU(118) effective flavours) obtained in latest ULALPS call. PARU(118) : (I) alpha_strong value obtained in latest ULALPS call. PARU(121) - PARU(130) : couplings of a new Z'; for fermion default values are given by the standard model Z values, assuming sin^2theta_W = 0.23. Note that e.g. the Z' width contains couplings-squared, and thus depends quadratically on the values below. PARU(121), PARU(122) : (D=-0.693,-1.) vector and axial couplings of down type quarks to Z'. PARU(123), PARU(124) : (D=0.387,1.) vector and axial couplings of up type quarks to Z'. PARU(125), PARU(126) : (D=-0.08,-1.) vector and axial couplings of leptons to Z'. PARU(127), PARU(128) : (D=1.,1.) vector and axial couplings of neutrinos to Z'. PARU(129) : (D=1.) the coupling Z' -> W+ + W- is taken to be PARU(129) * (the standard model Z -> W+ + W- coupling) * (m_W/m_Z')^2. This gives a Z' -> W+ + W- partial width that increases proportionately to the Z' mass. PARU(130) : (D=0.) in the decay chain Z' -> W+ + W- -> 4 fermions, the angular distribution in the W decays is supposed to be a mixture, with fraction 1-PARU(130) corresponding to the same angular distribution between the four final fermions as in Z -> W W (mixture of transverse and longitudinal W:s), and fraction PARU(130) corresponding to H -> W W the same way (longitudinal W:s). PARU(131) - PARU(136) : couplings of a new W'; for fermions default values are given by the standard model W values (i.e. V-A). Note that e.g. the W' width contains couplings-squared, and thus depends quadratically on the values below. PARU(131), PARU(132) : (D=1.,-1.) vector and axial couplings of a quark-antiquark pair to W'; is further multiplied by the ordinary CKM factors. PARU(133), PARU(134) : (D=1.,-1.) vector and axial couplings of a lepton-neutrino pair to W'. PARU(135) : (D=1.) the coupling W' -> Z0 + W is taken to be PARU(135) * (the standard model W -> Z0 + W coupling) * (m_W/m_W')^2. This gives a W' -> Z0 + W partial width that increases proportionately to the W' mass. PARU(136) : (D=0.) in the decay chain W' -> W + Z -> 4 fermions, the angular distribution in the W/Z decays is supposed to be a mixture, with fraction 1-PARU(130) corresponding to the same angular distribution between the four final fermions as in W -> W Z (mixture of transverse and longitudinal W/Z:s), and fraction PARU(130) corresponding to H+ -> W Z the same way (longitudinal W/Z:s). PARU(141) : (D=5.) tan(beta) parameter of a two Higgs doublet scenario, i.e. the ratio of vacuum expectation values. This affects mass relations and couplings in the Higgs sector. PARU(142) : (D=1.) the Z -> H+ + H- coupling is taken to be PARU(142) * (the MSSM Z -> H+ + H- coupling). PARU(143) : (D=1.) the Z' -> H+ + H- coupling is taken to be PARU(143) * (the MSSM Z -> H+ + H- coupling). PARU(145) : (D=1.) quadratically multiplicative factor in the Z' -> Z + H partial width in left-right-symmetric models, expected to be unity (see [Coc90]). PARU(146) : (D=1.) sin(2*alpha) parameter, enters quadratically as multiplicative factor in the W' -> W + H partial width in left-right-symmetric models (see [Coc90]). PARU(151) : (D=1.) multiplicative factor in the LQ -> q + l Yukawa coupling-squared, and thereby in the LQ partial width and the q + l -> LQ and other cross-sections. Specifically, lambda**2/(4pi) = PARU(151) * alpha_em, i.e. it corresponds to the k factor of [Hew88]. PARU(153) : (D=0.) anomalous magnetic moment of the W; eta = kappa - 1, where eta = 0 (kappa = 1) is the standard model value. PARU(155) : (D=1000. GeV) compositeness scale Lambda. PARU(156) : (D=1.) sign of interference term between standard cross-section and composite term (eta parameter); should be +-1. PARU(157) - PARU(159) : (D=3*1.) strength of SU(2), U(1) and SU(3) couplings, respectively, in an excited fermion scenario; cf. f, f' and f_s of [Bau90]. PARU(161) - PARU(168) : (D=5*1.,3*0.) multiplicative factors that can be used to modify the default couplings of the H particle in PYTHIA. Note that the factors enter quadratically in the partial widths. The default values correspond to the couplings given in the minimal one-Higgs-doublet standard model. PARU(161) : H coupling to down type quarks. PARU(162) : H coupling to up type quarks. PARU(163) : H coupling to leptons. PARU(164) : H coupling to Z0. PARU(165) : H coupling to W+-. PARU(168) : H coupling to H+- in gamma+gamma -> H0 loops, in MSSM sin(beta-alpha)+cos(2*beta)sin(beta+alpha)/(2*cos^2(theta_W)). PARU(171) - PARU(178) : (D=7*1.,0.) multiplicative factors that can be used to modify the default couplings of the H' particle in PYTHIA. Note that the factors enter quadratically in partial widths. The default values for PARU(171) - PARU(175) correspond to the couplings given to H in the minimal one-Higgs-doublet standard model, and are therefore not realistic in a two-Higgs-doublet scenario. The default values should therefore be changed appropriately by the user. Also the last two default values should be changed; for these the expressions of the minimal supersymmetric standard model (MSSM) are given to show parameter normalization. PARU(171) : H' coupling to down type quarks. PARU(172) : H' coupling to up type quarks. PARU(173) : H' coupling to leptons. PARU(174) : H' coupling to Z0. PARU(175) : H' coupling to W+-. PARU(176) : H' coupling to H0 + H0, in MSSM cos(2*alpha) cos(beta+alpha) - 2 sin(2*alpha) sin(beta+alpha). PARU(177) : H' coupling to A + A, in MSSM cos(2*beta) cos(beta+alpha). PARU(178) : H' coupling to H+- in gamma+gamma -> H' loops, in MSSM cos(beta-alpha)-cos(2*beta)cos(beta+alpha)/(2*cos^2(theta_W)). PARU(181) - PARU(190) : (D=3*1.,2*0.,2*1.,3*0.) multiplicative factors that can be used to modify the default couplings of the A particle in PYTHIA. Note that the factors enter quadratically in partial widths. The default values for PARU(181) - PARU(183) correspond to the couplings given to H in the minimal one-Higgs-doublet standard model, and are therefore not realistic in a two-Higgs-doublet scenario. The default values should therefore be changed appropriately by the user. PARU(184) and PARU(185) should be vanishing at the tree level, and are so set; normalization of these couplings agrees with what is used for H and H'. Also the other default values should be changed; for these the expressions of the minimal supersymmetric standard model (MSSM) are given to show parameter normalization. PARU(181) : A coupling to down type quarks. PARU(182) : A coupling to up type quarks. PARU(183) : A coupling to leptons. PARU(184) : A coupling to Z0. PARU(185) : A coupling to W+-. PARU(186) : A coupling to Z + H0 (or Z* to A + H), in MSSM cos(beta-alpha). PARU(187) : A coupling to Z + H' (or Z* to A + H'), in MSSM sin(beta-alpha). PARU(188) : As PARU(186), but coupling to Z' rather than Z. PARU(189) : As PARU(187), but coupling to Z' rather than Z. PARU(190) : A coupling to H+- in gamma+gamma -> A loops, 0 in MSSM. PARU(191) - PARU(195) : (D=4*0.,1.) multiplicative factors that can be used to modify the couplings of the H+ particle in PYTHIA. Currently only PARU(195) is in use. See above for related comments. PARU(195) : H+ coupling to W + H0 (or W* to H+ + H0), in MSSM cos(beta-alpha). MSTJ(1) : (D=1) choice of fragmentation scheme. = 0 : no jet fragmentation at all. = 1 : string fragmentation according to the Lund model. = 2 : independent fragmentation, according to specification in MSTJ(2) and MSTJ(3). MSTJ(2) : (D=3) gluon jet fragmentation scheme in independent fragmentation. = 1 : a gluon is assumed to fragment like a random d, u or s quark or antiquark. = 2 : as =1, but longitudinal (see PARJ(43), PARJ(44) and PARJ(59)) and transverse (see PARJ(22)) momentum properties of quark or antiquark substituting for gluon may be separately specified. = 3 : a gluon is assumed to fragment like a pair of a d, u or s quark and its antiquark, sharing the gluon energy according to the Altarelli-Parisi splitting function. = 4 : as =3, but longitudinal (see PARJ(43), PARJ(44) and PARJ(59)) and transverse (see PARJ(22)) momentum properties of quark and antiquark substituting for gluon may be separately specified. MSTJ(3) : (D=0) energy, momentum and flavour conservation options in independent fragmentation. Whenever momentum conservation is described below, energy and flavour conservation is also implicitly assumed. For physics details, see e.g. [Sjo84a]. = 0 : no explicit conservation of any kind. = 1 : particles share momentum imbalance compensation according to their energy (roughly equivalent to boosting event to CM frame). This is similar to the approach in the Ali et al. program [Ali80]. = 2 : particles share momentum imbalance compensation according to their longitudinal mass with respect to the imbalance direction. = 3 : particles share momentum imbalance compensation equally. = 4 : transverse momenta are compensated separately within each jet, longitudinal momenta are rescaled so that ratio of final jet to initial parton momentum is the same for all the jets of the event. This is similar to the approach in the Hoyer et al. program [Hoy79]. = 5 : only flavour is explicitly conserved. = 6 - 10 : as =1-5, except that above several colour singlet systems that followed immediately after each other in the event listing (e.g. qqbarqqbar) were treated as one single system, whereas here they are treated as separate systems. = -1 : independent fragmentation, where also particles moving backwards with respect to the jet direction are kept, and thus the amount of energy and momentum mismatch may be large. MSTJ(11) : (D=1) choice of longitudinal fragmentation function, i.e. how large a fraction of the energy available a newly-created hadron takes. = 1 : the Lund symmetric fragmentation function, see PARJ(41) - PARJ(45). = 2 : choice of some different forms for each flavour separately including Field-Feynman [Fie78] and SLAC heavy flavour [Pet83] functions, see PARJ(51) - PARJ(59). = 3 : hybrid scheme, where light flavours are treated with symmetric Lund (=1), but charm and heavier can be separately chosen, e.g. according to the SLAC function (=2). = 4 : the Lund symmetric fragmentation function (=1), for heavy endpoint quarks modified according to the Bowler (+Artru- Mennessier, Morris) space-time picture of string evolution, see PARJ(46). = 5 : as =4, but with possibility to interpolate between Bowler and Lund separately for c, b and t; see PARJ(46) - PARJ(48). MSTJ(12) : (D=2) choice of baryon production model. = 0 : no baryon-antibaryon pair production at all; initial diquark treated as a unit. = 1 : diquark-antidiquark pair production allowed; diquark treated as a unit [And82]. = 2 : diquark-antidiquark pair production allowed;, with possibility for diquark to be split according to the "popcorn" scheme [And85]. = 3 : as =2, but additionally the production of first rank baryons may be suppressed by a factor PARJ(19). MSTJ(13) : (D=0) generation of transverse momentum for endpoint quark(s) of single quark jet or qqbar jet system (in multijet events no endpoint transverse momentum is ever allowed for). = 0 : no transverse momentum for endpoint quarks. = 1 : endpoint quarks obtain transverse momenta like ordinary qqbar pairs produced in the field (see PARJ(21)); for two-jet systems the endpoints obtain balancing transverse momenta. MSTJ(14) : (D=1) treatment of colour singlet jet systems with low invariant masses. = 0 : no precautions are taken, meaning that problems may occur in LUSTRF (or LUINDF) later on. = 1 : small jet systems are allowed to collapse into two particles or, failing that, one single particle. Normally all small systems are treated this way, starting with the smallest one, but some systems would require more work and are left untreated; they include diquark-antidiquark pairs below the two-particle threshold. = -1 : special option for LUPREP calls, where no precautions are taken (as for =0), but, in addition, no checks are made on the presence of small-mass systems; i.e. LUPREP only rearranges colour strings. MSTJ(15) : (D=0) production probability for new flavours. = 0 : according to standard Lund parametrization, as given by PARJ(1) - PARJ(20). = 1 : according to probabilities stored by the user in PARF(201) - PARF(1960); note that no default values exist here, i.e. PARF must be set by the user. The MSTJ(12) switch can still be used to set baryon production mode, with the modification that MSTJ(12) = 2 here allows an arbitrary number of mesons to be produced between a baryon and an antibaryon (since the probability for diquark -> meson + new diquark is assumed independent of prehistory). MSTJ(21) : (D=2) form of particle decays. = 0 : all particle decays are inhibited. = 1 : a particle declared unstable in the MDCY vector, and with decay channels defined, may decay within the region given by MSTJ(22). A particle may decay into jets, which then fragment further according to the MSTJ(1) value. = 2 : as =1, except that a qqbar jet system produced in a decay (e.g. of a B meson) is always allowed to fragment according to string fragmentation, rather than according to the MSTJ(1) value (this means that momentum, energy and charge are conserved in the decay). MSTJ(22) : (D=1) cutoff on decay length for a particle that is allowed to decay according to MSTJ(21) and the MDCY value. = 1 : a particle declared unstable is also forced to decay. = 2 : a particle is decayed only if its average invariant lifetime is larger than PARJ(71). = 3 : a particle is decayed only if the decay vertex is within a distance PARJ(72) of the origin. = 4 : a particle is decayed only if the decay vertex is within a cylindrical volume with radius PARJ(73) in the xy-plane and extent to +-PARJ(74) in the z direction. MSTJ(23) : (D=1) possibility of having a shower evolving from a qqbar pair created as decay products. = 0 : never. = 1 : whenever the decay channel matrix element code is MDME(IDC,2) = 22, 23 or 33, the two first decay products (if they are partons) are allowed to shower, like a colour singlet subsystem, with maximum virtuality given by the invariant mass of the pair. MSTJ(24) : (D=2) particle masses. = 0 : discrete mass values are used. = 1 : particles registered as having a mass width in the PMAS vector are given a mass according to a truncated Breit-Wigner shape, linear in m: P(m) dm = 1/((m-m_0)^2 + Gamma^2/4) dm. = 2 : as =1, but gauge bosons (actually all particles with abs(KF) <= 100) are distributed according to a Breit-Wigner quadratic in m, as obtained from propagators. = 3 : as =1, but Breit-Wigner shape is always quadratic in m: P(m) dm^2 = 1/((m^2 - m_0^2)^2 + m_0^2*Gamma^2) dm^2. MSTJ(25) : (D=1) inclusion of the W propagator, in addition to the standard, "infinitely heavy" weak V-A matrix element, in the decay of a t, l or h quark, or chi lepton. = 0 : not included. = 1 : included. MSTJ(26) : (D=0) inclusion of B-B~ mixing in decays. = 0 : no. = 1 : yes, with mixing parameters given by PARJ(76) and PARJ(77). Mixing decays are not specially marked. = 2 : yes, as =1, but a B (B~) that decays as a B~ (B) is marked as K(I,1)=12 rather than the normal K(I,1)=11. MSTJ(41) : (D=1) type of branchings allowed in shower. = 0 : no branchings at all, i.e. shower is switched off. = 1 : QCD type branchings of quarks and gluons. = 2 : also emission of photons off quarks and leptons; the photons are assumed on mass-shell. MSTJ(42) : (D=2) branching mode for timelike showers. = 1 : conventional branching, i.e. without angular ordering. = 2 : coherent branching, i.e. with angular ordering. MSTJ(43) : (D=4) choice of z definition in branching. = 1 : energy fraction in grandmother's rest frame ("local, constrained"). = 2 : energy fraction in grandmother's rest frame assuming massless daughters, with energy and momentum reshuffled for massive ones ("local, unconstrained"). = 3 : energy fraction in CM frame of the showering partons ("global, constrained"). = 4 : energy fraction in CM frame of the showering partons assuming massless daughters, with energy and momentum reshuffled for massive ones ("global, unconstrained"). MSTJ(44) : (D=2) choice of alpha_strong scale for shower. = 0 : fixed at PARU(111) value. = 1 : running with Q^2 = m^2/4, m mass of decaying parton, Lambda as stored in PARJ(81) (natural choice for conventional showers). = 2 : running with Q^2 = z*(1-z)*m^2, i.e. roughly p_T^2 of branching, Lambda as stored in PARJ(81) (natural choice for coherent showers). MSTJ(45) : (D=5) maximum flavour that can be produced in shower by g -> q + qbar; also used to determine the maximum number of active flavours in the alpha_strong factor in parton showers (here with a minimum of 3). MSTJ(46) : (D=0) nonhomogeneous azimuthal distributions in branchings in shower. = 0 : azimuthal angle is chosen uniformly. = 1 : nonhomogeneous azimuthal angle in gluon decays due to a kinematics-dependent effective gluon polarization. Not meaningful for scalar model, i.e. then same as =0. = 2 : nonhomogeneous azimuthal angle in gluon decay due to interference with nearest neighbour (in colour). Not meaningful for Abelian model, i.e. then same as =0. = 3 : nonhomogeneous azimuthal angle in gluon decay due to both polarization (=1) and interference (=2). Not meaningful for Abelian model, i.e. then same as =1. Not meaningful for scalar model, i.e. then same as =2. MSTJ(47) : (D=1) corrections to the lowest order qqbarg three-jet matrix element at the first branching of either initial parton in a shower. = 0 : no corrections. = 1 : included whenever scattered partons are qqbar. = 2 : always included when shower starts from two partons. MSTJ(48) : (D=0) possibility to impose maximum angle of first branching in shower. = 0 : no explicit maximum angle. = 1 : maximum angle given by PARJ(85) for single showering parton, by PARJ(85) and PARJ(86) for pair of showering partons. MSTJ(49) : (D=0) possibility to change the branching probabilities according to some alternative toy models (note that the alpha-strong Q2 evolution may well be different in these models, but that only the MSTJ(44) options are at the disposal of the user). = 0 : standard QCD branchings. = 1 : branchings according to a scalar gluon theory, i.e. the splitting kernels in the Altarelli-Parisi equations are, with a common factor (alpha_strong)^*/(2*pi) omitted, P(q -> qg) = (2/3) * (1-z), P(g -> gg) = PARJ(87), P(g - qqbar) = PARJ(88) (for each separate flavour). The couplings of the gluon have been left as free parameters, since they depend on the colour structure assumed. Note that, since a spin 0 object decays isotropically, the gluon splitting kernels contain no z dependence. = 2 : brachings according to an Abelian vector gluon theory, i.e. the colour factors are changed (compared to QCD) according to C_F = 4/3 -> 1, N_C = 3 -> 0, T_R = 1/2 -> 3. Note that an Abelian model is not expected to contain any coherence effects between gluons, so that one should normally use MSTJ(42)=1 and MSTJ(46)=0 or =1. Also, alpha_strong is expected to increase with increasing Q^2 scale, rather than decrease. No such alpha_strong option is available; the one that comes closest is MSTJ(44)=0, i.e. a fix value. MSTJ(51) : (D=0) inclusion of Bose-Einstein effects. = 0 : no effects included. = 1 : effects included according to an exponential parametrization C_2(Q) = 1 + PARJ(92)*exp(-Q/PARJ(93)), where C_2(Q) represents the ratio of particle production at Q with Bose-Einstein effects to that without, and the relative momentum Q is defined by Q^2(p_1,p_2) = -(p_1 - p_2)^2 = (p_1 + p_2)^2 - 4m^2. Particles with width broader than PARJ(91) are assumed to have time to decay before Bose-Einstein effects are to be considered. = 2 : effects included according to a Gaussian parametrization C_2(Q) = 1 + PARJ(92)*exp(- (Q/PARJ(93))^2 ), with notation and comments as above. MSTJ(52) : (D=3) number of particle species for which Bose-Einstein correlations are to be included, ranged along the chain pi+, pi-, pi0, K+, K-, K0S, K0L, eta and eta'. Default corresponds to including all pions (pi+, pi-, pi0), 7 to including all Kaons as well, and 9 is maximum. MSTJ(91) : (I) flag when generating gluon jet with options MSTJ(2) = 2 or 4 (then =1, else =0). MSTJ(92) : (I) flag that a qqbar or gg pair or a ggg triplet created in LUDECY should be allowed to shower, is 0 if no pair or triplet, is the entry number of the first parton if a pair indeed exists, is the entry number of the first parton, with a - sign, if a triplet indeed exists. MSTJ(93) : (I) switch for ULMASS action. Is reset to 0 in ULMASS call. = 0 : ordinary action. = 1 : light (d, u, s, c, b) quark masses are taken from PARF(101) - PARF(105) rather than PMAS(1,1) - PMAS(5,1). Diquark masses are given as sum of quark masses, without spin splitting term. = 2 : as = 1. Additionally the constant terms PARF(121) and PARF(122) are subtracted from quark and diquark masses, respectively. PARJ(1) : (D=0.10) is P(qq)/P(q), the suppression of diquark-antidiquark pair production in the colour field, compared to quark-antiquark production. PARJ(2) : (D=0.30) is P(s)/P(u), the suppression of s quark pair production in the field compared to u or d pair production. PARJ(3) : (D=0.4) is (P(us)/P(ud))/(P(s)/P(d)), the extra suppression of strange diquark production compared to the normal suppression of strange quarks. PARJ(4) : (D=0.05) is (1/3)P(ud_1)/P(ud_0), the suppression of spin 1 diquarks compared to spin 0 ones (excluding the factor 3 coming from spin counting). PARJ(5) : (D=0.5) parameter determining relative occurence of baryon production by BMBbar and by BBbar configurations in the popcorn baryon production model, roughly P(BMBbar)/(P(BBbar)+P(BMBbar)) = PARJ(5)/(0.5+PARJ(5)). PARJ(6) : (D=0.5) extra suppression for having a ssbar pair shared by the B and Bbar of a BMBbar situation. PARJ(7) : (D=0.5) extra suppression for having a strange meson M in a BMBbar configuration. PARJ(11) : (D=0.5) is the probability that a light meson (containing u and d quarks only) has spin 1 (with 1-PARJ(11) the probability for spin 0) when formed in fragmentation. PARJ(12) : (D=0.6) is the probability that a strange meson has spin 1. PARJ(13) : (D=0.75) is the probability that a charm or heavier meson has spin 1. PARJ(14) : (D=0.) is the probability that a spin = 0 meson is produced with an orbital angular momentum 1, for a total spin = 1. PARJ(15) : (D=0.) is the probability that a spin = 1 meson is produced with an orbital angular momentum 1, for a total spin = 0. PARJ(16) : (D=0.) is the probability that a spin = 1 meson is produced with an orbital angular momentum 1, for a total spin = 1. PARJ(17) : (D=0.) is the probability that a spin = 1 meson is produced with an orbital angular momentum 1, for a total spin = 2. Note on PARJ(11) - PARJ(17) : the end result of the numbers above is that, with i = 11, 12 or 13, depending on flavour content, P(S = 0, L = 0, J = 0) = (1-PARJ(i)) * (1-PARJ(14)), P(S = 0, L = 1, J = 1) = (1-PARJ(i)) * PARJ(14), P(S = 1, L = 0, J = 1) = PARJ(i) * (1-PARJ(15)-PARJ(16)-PARJ(17)), P(S = 1, L = 1, J = 0) = PARJ(i) * PARJ(15), P(S = 1, L = 1, J = 1) = PARJ(i) * PARJ(16), P(S = 1, L = 1, J = 2) = PARJ(i) * PARJ(17), where S is the quark "true" spin and J is the total spin, usually called the spin of the meson. PARJ(18) : (D=1.) is an extra suppression factor multiplying the ordinary SU(6) weight for spin 3/2 baryons, and hence a means to break SU(6) in addition to the dynamic breaking implied by PARJ(2), PARJ(3), PARJ(4), PARJ(6) and PARJ(7). PARJ(19) : (D=1.) extra baryon suppression factor, which multiplies the ordinary diquark-antidiquark production probability for the breakup closest to the endpoint of a string, but leaves other breaks unaffected. Is only used for MSTJ(12) = 3. PARJ(21) : (D=0.35 GeV) corresponds to the width in the Gaussian p_x and p_y transverse mementum distributions for primary hadrons. PARJ(22) : (D=1.) relative increase in transverse momentum in a gluon jet generated with MSTJ(2) = 2 or 4. PARJ(31) : (D=0.1 GeV) gives the remaining W_+ below which the generation of a single jet is stopped (it is chosen smaller than a pion mass, so that no hadrons moving in the forward direction are missed). PARJ(32) : (D=1. GeV) is, with quark masses added, used to define the minimum allowable energy of a colour singlet jet system. PARJ(33) - PARJ(34) : (D=0.8 GeV, 1.5 GeV) are, together with quark masses, used to define the remaining energy below which the fragmentation of a jet system is stopped and two final hadrons formed. PARJ(33) is normally used, except for MSTJ(11) = 2, when PARJ(34) is used. PARJ(36) : (D=2.) represents the dependence on the mass of the final quark pair for defining the stopping point of the fragmentation. Is strongly correlated to the choice of PARJ(33) - PARJ(35). PARJ(37) : (D=0.2) relative width of the smearing of the stopping point energy. PARJ(38) - PARJ(39) : (D=2.5, 0.6) refers to the probability for reverse rapidity ordering of the final two hadrons, with transverse masses m_T_1 and m_T_2, for a total remaining transverse mass W_rem. The probability is parametrized in the form P(reverse ordering) = 0.5*((m_T_1 + m_T_2)/W_rem)^d where d = PARJ(38)*(m_T_1^2 + m_T_2^2)^2 for MSTJ(11) not 2, and d = PARJ(39) for MSTJ(11) = 2. PARJ(41), PARJ(42) : (D=0.5, 0.9 GeV^-2) give the a and b parameters of the symmetric Lund fragmentation function for MSTJ(11) = 1 (and MSTJ(11) = 3 for ordinary hadrons) f(z) = (1/z) * (1-z)^a * exp(-b*m_T^2/z) where m_T is the transverse mass of the hadron. PARJ(43), PARJ(44) : (D=0.5, 0.9 GeV^-2) give the a and b parameters as above for the special case of a gluon jet generated with IF and MSTJ(2) = 2 or 4. PARJ(45) : (D=0.5) the amount by which the effective a parameter in the Lund flavour dependent symmetric fragmentation function is assumed to be larger than the normal a when diquarks are produced. More specifically, with f(z) = (1/z) * z^a1 * ((1-z)/z)^a2 * exp(-b*m_T^2/z) a1 = PARJ(41) when considering the fragmentation of a quark and = PARJ(41) + PARJ(45) for the fragmentation of a diquark, with corresponding expression for a2 depending on whether the newly created object is a quark or diquark (for an independent gluon jet generated with MSTJ(2) = 2 or 4, replace PARJ(41) by PARJ(43)). In the popcorn model, a meson created in between the baryon and antibaryon has a1 = a2 = PARJ(41) + PARJ(45). PARJ(46) - PARJ(48) : (D=3*1.) modification of the Lund symmetric fragmentation for heavy endpoint quarks according to the recipe by Bowler (based on Artru-Mennessier space-time picture of string evolution, also studied by Morris), available when MSTJ(11) = 4 or = 5 is selected. Specifically, the shape is taken to be f(z) = (1/z**(1 + r*b*m_Q^2)) * (1-z)**a * exp(-b*m_T^2/z), where m_Q is the heavy quark mass, and a, b and m_T as described above (see PARJ(41) - PARJ(42)). If MSTJ(11) = 4 then r = PARJ(46) for all flavours, while if MSTJ(11) = 5 then r = PARJ(46) for c, r = PARJ(47) for b and r = PARJ(48) for t and heavier. PARJ(46) - PARJ(48) thus provide a possibility to interpolate between the 'pure' Bowler shape, r = 1., and the normal Lund one, r = 0. The additional modifications made in PARJ(43) - PARJ(45) are automatically taken into account, if necessary. PARJ(51) - PARJ(58) : (D=3*0.77, 5*0.) give four possible ways to parametrize the fragmentation function for MSTJ(11) = 2 (and MSTJ(11) = 3 for charm and heavier). The fragmentation of each flavour KFL may be chosen separately; for a diquark the flavour of the heaviest quark is used. With c = PARJ(50+KFL) the parametrizations are: 0 <= c <= 1 : Field-Feynman, f(z) = 1 - c+ 3*c*(1-z)^2; -1 <= c < 0 : SLAC, f(z) = 1/(z*(1-1/z-(-c)/(1-z))^2); c > 1 : power peaked at z=0, f(z) = (1-z)^(c-1); c < -1 : power peaked at z=1, f(z) = z^(-c-1). PARJ(59) : (D=1.) replaces PARJ(51) - PARJ(53) for gluon jet generated with MSTJ(2) = 2 or 4. PARJ(61) - PARJ(63) : (D=4.5, 0.7, 0.) parametrizes the energy dependence of the primary multiplicity distribution in phase space decays. PARJ(64) : (0.003 GeV) minimum kinetic energy in decays (safety margin for numerical precision errors). PARJ(65) : (D=0.5 GeV) mass which, in addition to the spectator quark or diquark mass, is not assumed to partake in the weak decay of a heavy quark in a hadron. PARJ(66) : (D=0.5) relative probability that colour is rearranged when two singlets are to be formed from decay products. Only applies for MDME(IDC,2) = 11 - 30, i.e. low mass phase space decays. PARJ(71) : (D=10 mm) maximum average invariant lifetime for particles allowed to decay in the MSTJ(22) = 2 option. With the default value, K_S0, Lambda, Sigma-, Sigma+, Xi-, XI0 and Omega- are stable (in addition to those normally taken to be stable), but charm and bottom do still decay. PARJ(72) : (D=1000 mm) maximum distance from the origin at which a decay is allowed to take place in the MSTJ(22) = 3 option. PARJ(73) : (D=100 mm) maximum cylindrical distance rho = sqrt(x^2 + y^2) from the origin at which a decay is allowed to take place in the MSTJ(22) = 4 option. PARJ(74) : (D=1000 mm) maximum z distance from the origin at which a decay is allowed to take place in the MSTJ(22) = 4 option. PARJ(76) : (D=0.7) mixing parameter x_d = (Delta-M)/Gamma in B0-B~0 system. PARJ(77) : (D=10.) mixing parameter x_s = (Delta-M)/Gamma in B_s0-B_s~0 system. PARJ(81) : (D=0.40 GeV) Lambda value used in running alpha_strong for parton showers (see MSTJ(44)). PARJ(82) : (D=1.0 GeV) invariant mass cutoff m_min of parton showers, below which partons are not assumed to radiate. For Q^2 = p_T^2 (MSTJ(44) = 2) PARJ(82)/2 additionally gives the minimum p_T of a branching. To avoid infinite alpha_strong values, one must have PARJ(82) > 2*PARJ(81) for MSTJ(44) >= 1 (this is automatically checked in the program, with 2.2*PARJ(81) as the lowest value attainable). PARJ(83) : (D=1.0 GeV) invariant mass cutoff m_min used for photon emission in parton showers, below which quarks and leptons are not assumed to radiate. The function of PARJ(83) closely parallels that of PARJ(82) for QCD branchings, but there is a priori no requirement that the two be equal. PARJ(85), PARJ(86) : (D=10.,10.) maximum opening angles allowed in the first branching of parton showers; see MSTJ(48). PARJ(87) : (D=0.) coupling of g -> gg in scalar gluon shower, see MSTJ(49) = 1. PARJ(88) : (D=0.) coupling of g -> qqbar in scalar gluon shower (per quark species), see MSTJ(49) = 1. PARJ(91) : (D=0.020 GeV) minimum particle width in PMAS(KC,2), above which particle decays are assumed to take place before the stage where Bose-Einstein effects are introduced. PARJ(92) : (D=1.) nominal strength of Bose-Einstein effects for Q = 0, see MSTJ(51). This parameter, often denoted lambda, expresses the amount of incoherence in particle production. Due to the simplified picture used for the Bose-Einstein effects, in particular for effects from three nearby identical particles, the actual lambda of the simulated events may be larger than the input value. PARJ(93) : (D=0.20 GeV) size of the Bose-Einstein effect region in terms of the Q variable, see MSTJ(51). The more conventional measure, in terms of the radius R of the production volume, is given by R = hslash/PARJ(93) = 0.2 fm GeV/PARJ(93) = PARU(3)/PARJ(93). ______________________________________________________________________ 2.8. Further Parameters and Particle Data The following commonblocks are maybe of a more peripheral interest, with the exception of the MDCY array, which allows a selective inhibiting of particle decays, and PMAS(6,1), the top quark mass. COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4) Purpose: to give access to a number of flavour treatment constants or parameters and particle/parton data. KCHG(KC,1) : three times particle/parton charge for compressed code KC. KCHG(KC,2) : colour information for compressed code KC. = 0 : colour singlet particle. = 1 : quark or antidiquark. = -1 : antiquark or diquark. = 2 : gluon. KCHG(KC,3) : particle/antiparticle distinction for compressed code KC. = 0 : the particle is its own antiparticle. = 1 : a nonidentical antiparticle exists. PMAS(KC,1) : particle/parton mass (in GeV) for compressed code KC. PMAS(KC,2) : the total width Gamma (in GeV) of an assumed symmetric Breit-Wigner mass shape for compressed particle code KC. PMAS(KC,3) : the maximum deviation (in GeV) from the PMAS(KC,1) value at which the Breit-Wigner shape above is truncated. PMAS(KC,4) : the average lifetime tau for compressed particle code KC, with c*tau in mm, i.e. tau in units of 3.33 * 10^-12 s. PARF(1) - PARF(60) : give a parametrization of the uubar-ddbar-ssbar flavour mixing for pseudoscalar, vector and tensor mesons. PARF(61) - PARF(80) : give flavour SU(6) weights for the production of a spin 1/2 or spin 3/2 baryon from a given diquark-quark combination. PARF(101) - PARF(108) : first five contain d, u, s, c and b constituent masses, as to be used in mass formulae, and should not be changed. For t, l and h masses the current values stored in PMAS(6,1) - PMAS(8,1) are copied in. PARF(111), PARF(112) : (D=0.0, 0.11 GeV) constant terms in the mass formulae for heavy mesons and baryons, respectively (with diquark getting 2/3 of baryon). PARF(113), PARF(114) : (D=0.16,0.048 GeV) factors which, together with Clebsch-Gordan coefficients and quark constituent masses, determine the mass splitting due to spin-spin interactions for heavy mesons and baryons, respectively. The latter factor is also used for the splitting between spin 0 and spin 1 diquarks. PARF(115) - PARF(118) : (D=0.50, 0.45, 0.55, 0.60 GeV), constant mass terms, added to the constituent masses, to get the mass of heavy mesons with orbital angular momentum L = 1. The four numbers are for pseudovector mesons with quark spin 0, and for scalar, pseudovector and tensor mesons with quark spin 1, respectively. PARF(121),PARF(122) : (D=0.1, 0.2 GeV) constant terms, which are subtracted for quark and diquark masses, respectively, in defining the allowed phase space in particle decays into partons. PARF(201) - PARF(1960) : (D=1760*0) relative probabilities for flavour production in the MSTJ(15) = 1 option; to be defined by the user before any JETSET calls. The index in PARF is of the compressed form 120 + 80*KTAB1 + 25*KTABS + KTAB3. Here KTAB1 is the old flavour, fixed by preceding fragmentation history, while KTAB3 is the new flavour, to be selected according to the relevant relative probabilities (except for the very last particle, produced when joining two jets, where both KTAB1 and KTAB3 are known). Only the most frequently appearing quarks/diquarks are defined, according to the code 1 = d, 2 = u, 3 = s, 4 = c, 5 = b, 6 = t, 7 = dd_1, 8 = ud_0, 9 = ud_1, 10 = uu_1, 11 = sd_0, 12 = sd_1, 13 = su_0, 14 = su_1, 15 = ss_1, 16 = cd_0, 17 = cd_1, 18 = cu_0, 19 = cu_1, 20 = cs_0, 21 = cs_1, 22 = cc_1. These are thus the only possibilities for the new flavour to be produced; for an occasional old flavour not on this list, the ordinary relative flavour production probabilities will be used. Given the initial and final flavour, the intermediate hadron that is produced is almost fixed. (Initial and final diquark here corresponds to "popcorn" production of mesons intermediate between a baryon and an antibaryon). The additional index KTABS gives the spin type of this hadron, with 0 = pseudoscalar meson or Lambda-like spin 1/2 baryon, 1 = vector meson or Sigma-like spin 1/2 baryon, 2 = tensor meson or spin 3/2 baryon. (Some meson multiplets, not frequently produced, are not accessible by this parametrization.) Note that some combinations of KTAB1, KTAB3 and KTABS do not correspond to a physical particle (a Lambda-like baryon must contain three different quark flavours, a Sigma-like one at least two), and that the user must see to it that the corresponding PARF entries are vanishing. One additional complication exist when KTAB3 and KTAB1 denote the same flavour content (normally KTAB3 = KTAB1, but for diquarks the spin freedom may give KTAB3 = KTAB1 +- 1): then a flavour neutral meson is to be produced, and here ddbar, uubar and ssbar states mix (heavier flavour states do not, and these are therefore no problem). For these cases the ordinary KTAB3 value gives the total probability to produce either of the mesons possible, while KTAB3 = 23 gives the relative probability to produce the lightest meson state (pi0, rho0, a_20), KTAB3 = 24 relative probability for the middle meson (eta, omega, f_20), and KTAB3 = 25 relative probability for the heaviest one (eta', phi, f'_20). Note that, for simplicity, these relative probabilities are assumed the same whether initial and final diquark have the same spin or not; the total probability may well be assumed different, however. As a general comment, the sum of PARF values for a given KTAB1 need not be normalized to unity, but rather the program will find the sum of relevant weights and normalize to that. The same goes for the KTAB3 = 23 - 25 weights. This makes it straightforward to use one common setup of PARF values and still switch between different MSTJ(12) baryon production modes. VCKM(I,J) : squared matrix elements of the Cabibbo-Kobayashi-Maskawa flavour mixing matrix. Are not currently used in JETSET, but appear here e.g. for usage in PYTHIA. I : up type generation index, i.e. 1 = u, 2 = c, 3 = t and 4 = h. J : down type generation index, i.e. 1 = d, 2 = s, 3 = b and 4 = l. COMMON/LUDAT3/MDCY(500,3),MDME(2000,2),BRAT(2000),KFDP(2000,5) Purpose: to give access to particle decay data and parameters. In particular, MDCY(KC,1) may be used to switch on or off the decay of a given particle species, and MDME(IDC,1) to switch on or off an individual decay channel of a particle. For quarks, leptons and gauge bosons, a number of decay channels are included that are not allowed for on-mass-shell particles, see MDME(IDC,2) = 102. These channels are not currently used in JETSET, but instead find applications in PYTHIA. MDCY(KC,1) : switch to tell whether a particle may be allowed to decay or not. For compressed particle code KC, MDCY(KC,1) = 1 signifies that the particle is allowed to decay (if decay information is defined below for the particle), whereas MDCY(KC,1) = 0 implies that the particle is not allowed to decay. The user is reminded that the way to know the KC value is to use the LUCOMP function, e.g. to switch off Lambda decay: MDCY(LUCOMP(3122),1) = 0. MDCY(KC,2) : gives the entry point into the decay channel table for compressed particle code KC. Is 0 if no decay channels have been defined. MDCY(KC,3) : gives the total number of decay channels defined for compressed particle code KC, independently of whether they have been assigned a nonvanishing branching ratio or not. Thus the last decay channel is MDCY(KC,2)+MDCY(KC,3)-1. MDME(IDC,1) : on/off switch for individual decay channel IDC. In addition, a channel may be left selectively open; this has some special applications in PYTHIA which are not currently used in JETSET. Effective branching ratios are automatically recalculated for the decay channels left open. If a particle is allowed to decay by the MDCY(KC,1) value, at least one channel must be left open by the user. A list of decay channels with current IDC numbers may be obtained with LULIST(12). = -1 : this is a non-standard model decay mode, which by default is assumed not to exist. Normally, this option is used for decays involving fourth generation or H+- particles. = 0 : channel is switched off. = 1 : channel is switched on. = 2 : channel is switched on for a particle but off for an antiparticle. It is also on for a particle its own antiparticle, i.e. here it means the same as =1. = 3 : channel is switched on for an antiparticle but off for a particle. It is off for a particle its own antiparticle. = 4 : in the production of a pair of equal or charge conjugate resonances in PYTHIA, say H0 -> W+ W-, either one of the resonances is allowed to decay according to this group of channels, but not both. If the two particles of the pair are different, the channel is on. Within JETSET, this option only means that the channel is switched off. = 5 : as =4, but an independent group of channels, such that in a pair of equal or charge conjugate resonances the decay of either resonance may be specified independently. If the two particles in the pair are different, the channel is off. Within JETSET, this option only means that the channel is switched off. WARNING: the two values -1 and 0 may look similar, but in fact are quite different. In neither case the channel so set is generated, but in the latter case the channel still contributes to the total width of a resonance, and thus affects both simulated line shape and the generated cross-section when the PYTHIA program is run. The value 0 is appropriate to a channel we assume exists, even if we are not currently simulating it, while -1 should be used for channels we believe do not exist. In particular, users are warned unwittingly to set fourth generation channels 0 (rather than -1), since by now the support for a fourth generation is small. Remark: all the options above may be freely mixed. The difference, for those cases where both make sense, between using values 2 and 3 and using 4 and 5 is that the latter automatically include charge conjugate states, e.g. H0 -> W+ W- -> e+ nue d ubar or dbar u e- nuebar, but the former only one of them. In calculations of the joint branching ratio, this makes a factor 2 difference. MDME(IDC,2) : information on special matrix element treatment for decay channel IDC. In addition to the outline below, special rules apply for the order in which decay products should be given, so that matrix elements and colour flow is properly treated. = 0 : no special matrix element treatment; partons and particles are copied directly to the event record, with momentum distributed according to phase space. = 1 : omega and phi decays into three pions. = 2 : pi0 or eta Dalitz decay to gamma e+ e-. = 3 : used for vector meson decays into two pseudoscalars, to signal that, in the decay chain PS_0 -> PS_1 + V -> PS_1 + PS_2 + PS_3 (PS pseudoscalar, V vector), the momentum of PS_2 and PS_3 have a cos^2(theta_02) distribution in the rest frame of V, and that, in the decay chain PS_0 -> gamma + V -> gamma + PS_2 + PS_3, the momentum of PS_2 and PS_3 have a sin^2(theta_02) distribution in the rest frame of V. = 4 : decay of a spin 1 "onium" resonance to three gluons or to a photon and two gluons. The gluons may subsequently develop a shower if MSTJ(23) = 1. = 11 : phase space production of hadrons from the quarks available. = 12 : as =11, but for onia resonances, with the option of modifying the multiplicity distribution separately. = 13 : as =11, but at least three hadrons to be produced (useful when the two-body decays are given explicitly). = 14 : as =11, but at least four hadrons to be produced. = 15 : as =11, but at least five hadrons to be produced. = 22 - 30 : phase space production of hadrons from the quarks available, with the multiplicity fixed to be MDME(IDC,2)-20, i.e. 2 - 10. = 31 : two or more quarks and particles are distributed according to phase space. If three or more products, the last product is a spectator quark, i.e. sitting at rest with respect to the decaying hadron. = 32 : a qqbar or gg pair, distributed according to phase space (in angle), and allowed to develop a shower if MSTJ(23) = 1. = 33 : a triplet qXqbar, where X is either a gluon or a colour singlet particle, the final particle (qbar) is assumed to sit at rest with respect to the decaying hadron, and the two first particles (qX) are allowed to develop a shower if MSTJ(23) = 1. = 41 : weak decay, where particles are distributed according to phase space, multiplied by a factor from the expected shape of the momentum spectrum of the direct product of the weak decay (the tau neutrino in tau decay). = 42 : weak decay matrix element for quarks and leptons. Products may be given either in terms of quarks or hadrons, or leptons for some channels. If the spectator system is given in terms of quarks, it is assumed to collapse into one particle from the onset. If the virtual W decays into quarks, these quarks are converted to particles, according to phase space in the W rest frame, as in =11. Is intended for tau, charm and bottom. = 43 : as =42, but if the W decays into quarks, these will either appear as jets or, for small masses, collapse into a one- or two-body system. = 44 : weak decay matrix element for quarks and leptons, where the spectator system may collapse into one particle for a small invariant mass. If the first two decay products are a qq'bar pair, they may develop a parton shower, if MSTJ(23) = 1. Is intended for top and beyond, but largely superseded by the following option. = 45 : weak decay q -> W + q' or l -> W + l', where the W is registered as a decay product and subsequently treated with MDME = 46. To distinguish from ordinary on-shell W:s, code KF = +-89 is used. The virtual W mass is selected according to the standard weak decay matrix element, times the W propagator (for MSTJ(25) = 1). There may be two or three decay products; if a third this is a spectator taken to sit at rest. The spectator system may collapse into one particle. Is intended for top and beyond. = 46 : W (KF = 89) decay into qq'bar or lnu_l according to relative probabilities given by couplings (as stored in the BRAT vector) times a dynamical phase space factor given by the current W mass. In the decay, the correct V-A angular distribution is generated if the W origin is known (heavy quark or lepton). This is therefore the second step of a decay with MDME = 45. A qq'bar pair may subsequently develop a shower. = 84 - 88 : map the decay of this particle onto the generic c, b, t, l or h decay modes defined for KC = 84 - 88. = 101 : this is not a proper decay channel, but only to be considered as a continuation line for the decay product listing of the immediately preceding channel. Since the KFDP array can contain five decay products per channel, with this code it is possible to define channels with up to ten decay products. It is not allowed to have several continuation lines after each other. = 102 : this is not a proper decay channel for an on-mass-shell (or nearly so) decaying particle, and is therefore assigned branching ratio 0. As an off-mass-shell particle, this decay mode is allowed, however. By including this channel among the others, the switches MDME(IDC,1) may be used to allow or forbid these channels in hard processes, with cross-sections to be calculated separately. As an example, gamma -> u ubar is not possible for a massless photon, but is an allowed channel in e+e- annihilation. BRAT(IDC) : give branching ratios for the different decay channels. In principle, the sum of branching ratios for a given particle should be unity. Since the program anyway has to calculate the sum of branching ratios left open by the MDME(IDC,1) values and normalize to that, the user need not explicitly ensure this normalization, however. (Warnings are printed in LUUPDA(2) calls if the sum is not unity, but this is entirely intended as a help for finding user mistypings.) For decay channels with MDME(IDC,2) > 80 the BRAT values are dummy. KFDP(IDC,J) : contain the decay products in the different channels, with five positions J = 1 through 5 reserved for each channel IDC. The decay products are given following the standard KF code for jets and particles, with 0 for trailing empty positions. Note that the MDME(IDC+1,2) = 101 option allows the user to double the maximum number of decay product in a given channel from 5 to 10, with the five latter products stored KFDP(IDC+1,J). COMMON/LUDAT4/CHAF(500) CHARACTER CHAF*8 Purpose: to give access to character type variables. CHAF : particle names (excluding charge) according to KC code. ______________________________________________________________________ 2.9. Miscellaneous Comments The previous sections have dealt with the subroutine options and variables one at a time. This is certainly important, but for a full use of the capabilities of the Lund Monte Carlo, it is also necessary to understand how to make different pieces work together. This is something that can not be explained fully in a manual, but must also be learnt by trial and error. This section contains some examples of relationships between subroutines, commonblocks and parameters. It also contains comments on issues that did not fit in naturally anywhere else, but still might be useful to have on record. Very often, the output of the program is to be fed into a subsequent detector simulation program. It therefore becomes necessary to set up an interface between the LUJETS commonblock and the detector model. Preferrably this should be done via the HEPEVT standard commonblock, see section 2.11, but sometimes this may not be convenient. If a LUEDIT(2) call is made, the remaining entries exactly correspond to those an ideal detector could see: all non-decayed particles, with the exception of neutrinos. The translation of momenta should be trivial (if need be, a LUROBO call can be made to rotate the "preferred" z direction to whatever is the longitudinal direction of the detector), and so should the translation of particle codes. In particular, if the detector simulation program also uses the standard Particle Data Group codes, no conversion at all is needed. The problem then is to select which particles are allowed to decay, and how decay vertex information should be used. Several switches regulate which particles are allowed to decay. First, the master switch MSTJ(21) can be used to switch on/off all decays (and it also contains a choice of how fragmentation should be interfaced). Second, a particle must have decay modes defined for it, i.e. the corresponding MDCY(KC,2) and MDCY(KC,3) entries must be nonzero for compressed code KC = LUCOMP(KF). This is true for all colour neutral particles except the neutrinos, the photon, the proton and the neutron. (This statement is actually not fully correct, since irrelevant "decay modes" with MDME(IDC,2) = 102 exist in some cases.) Third, the individual switch in MDCY(KC,1) must be on. Of all the particles with decay modes defined, only mu+-, pi+-, K+- and K_L0 are by default considered stable. Finally, if MSTJ(22) does not have its default value 1, checks are also made on the lifetime of a particle before it is allowed to decay. In the simplest alternative, MSTJ(22) = 2, the comparison is based on the average lifetime, or rather c*tau, measured in mm. Thus if the limit PARJ(71) is (the default) 10 mm, then decays of K_S0, Lambda, Sigma-, Sigma+, Xi-, Xi0 and Omega- are all switched off, but charm and bottom still decay. No c*tau values below 1 micron are defined. With the two options MSTJ(22) = 3 or 4, a spherical or cylindrical volume is defined around the origin, and all decays taking place inside this volume are ignored. Whenever a particle is in principle allowed to decay, i.e. MSTJ(21) and MDCY on, an invariant lifetime is selected once and for all and stored in V(I,5). The K(I,1) is then also changed to 4. For MSTJ(22) = 1, such a particle will also decay, but else it could remain in the event record. It is then possible, at a later stage, to expand the volume inside which decays are allowed, and do a new LUEXEC call to have particles fulfilling the new conditions (but not the old) decay. As a further option, the K(I,1) code may be put to 5, signalling that the particle will definitely decay in the next LUEXEC call, at the vertex position given (by the user) in the V vector. This then allows the Lund decay routines to be used inside a detector simulation program, as follows. For a particle which did not decay before entering the detector, its point of decay is still well defined (in the absence of deflections by electric or magnetic fields): V'(j) = V(I,j) + V(I,5)*P(I,j)/P(I,5), j = 1 - 4. If it interacts before that point, the detector simulation program is left to handle things. If not, the V vector is updated according to the formula above, K(I,1) is set to 5, and LUEXEC is called, to give a set of decay products, that can again be tracked. A further possibility is to force particles to decay into specific decay channels; this may be particularly interesting for charm or bottom physics. The choice of channels left open is determined by the values of the switches MDME(IDC,1) for decay channel IDC (use LULIST(12) to obtain the full listing). One or several channels may be left open; in the latter case effective branching ratios are automatically recalculated without the need for user intervention. It is also possible to differentiate between which channels are left open for particles and which for antiparticles. Lifetimes are not affected by the exclusion of some decay channels. Note that, whereas forced decays can enhance the efficiency for several kinds of studies, it can also introduce unexpected biases, in particular when events may contain several particles with forced decays. A nontrivial question is to know which parameter values to use. The default values stored in the program are based on comparisons with e+e- data at around 30 GeV, using a parton shower picture followed by string fragmentation. If fragmentation is indeed an universal phenomenon, as we would like to think, then the same parameters should also apply at other energies and in other processes. Note, however, that the choice of parameters is intertwined with the choice of perturbative QCD description. If instead matrix elements are used, a best fit to 30 GeV data would require the values PARJ(21) = 0.40, PARJ(41) = 1.0 and PARJ(42) = 0.7. With matrix elements one does not expect an energy independence of the parameters, since the effective minimum invariant mass cutoff is then energy dependent, i.e. so is the amount of soft gluon emission effects lumped together with the fragmentation parameters. A mismatch in the perturbative QCD treatment could also lead to small differences between different processes. It is often said that the string fragmentation model contains a wealth of parameters. This is certainly true, but it must be remembered that most of these deal with flavour properties, and to a large extent factorize from the treatment of the general event shape. In a fit to the latter it is therefore usually enough to consider the parameters of the perturbative QCD treatment, like Lambda in alpha_strong and a shower cutoff m_min (or alpha_strong itself and y_min, if matrix elements are used), the a and b parameter of the Lund symmetric fragmentation function (PARJ(41) and PARJ(42)) and the width of the transverse momentum distribution (sigma = PARJ(21)). In addition, the a and b parameters are very strongly correlated by the requirement of having the correct average multiplicity, such that in a typical chi^2 plot, the allowed region corresponds to a very narrow but very long valley, stretched diagonally from small (a,b) pairs to large ones. As to the flavour parameters, these are certainly many more, but most of them are understood qualitiatively within one single framework, that of tunneling pair production of flavours. Since the use of independent fragmentation has fallen in disrespect, it should be pointed out that the default parameters here are not particularly well tuned to the data. This especially applies if one, in addition to asking for independent fragmentation, also asks for another setup of fragmentation functions, i.e. other than the standard Lund symmetric one. In particular, note that most fits to the popular Peterson et al. (SLAC) heavy flavour fragmentation function are based on the actual observed spectrum. In a Monte Carlo simulation, one must then start out with something harder, to compensate for the energy lost by initial state photon radiation and gluon bremsstrahlung. Since independent fragmentation is not collinear safe (i.e, the emission of a collinear gluon changes the properties of the final event), the tuning is strongly dependent on the perturbative QCD treatment chosen. All the parameters needed for a tuning of independent fragmentation are available, however. The masses of most frequently used particles are taken from tables. For some rare charm and bottom hadrons, and for heavier flavour hadrons, this would be unwieldy, and instead mass formulae are used, based on the quark content. For the known quarks d, u, s, c and b, the masses used for this purpose are actually the ones stored in positions 101 - 105 in the PARF vector, rather than the ones found in PMAS. This means that the PMAS masses can be freely changed by the user, to modify the masses that appear in the event record, without courting disaster elsewhere (since mass formulae typically contain 1/m terms from spin-spin splittings, is is necessary to have the nonzero "constituent" masses here). Thus a user should never touch the mass values stored in PARF. For the heavier flavours top, low and high, the current PMAS values are always used. For these flavours, the only individually defined hadrons are the flavour neutral eta, Theta and chi_2 states. A complete change of top mass in the program thus requires changing PMAS(6,1), PMAS(LUCOMP(661),1), PMAS(LUCOMP(663),1) and PMAS(LUCOMP(665),1). Since the latter heavy flavour diagonal states are not normally produced in fragmentation, it would be no disaster to forget changing their masses. The masses of several responances, like rho, K* and Delta, are by default distributed according to simle Breit-Wigner shapes, suitably truncated, with the truncation chosen so that no problems are encountered in decay chains. It should be emphasized that such a simple approach may be a poor approximation, and should never be used for detailed studies of resonance production. It will, however, give a first approximation to the smearing caused by variable resonance masses. Another option (this one by default off) with a similar level of crudity is the simulation of Bose-Einstein effects. Here the detailed physics is not that well understood, see e.g. the review [Lor88]. What is offered is an algorith, more than just a parametrization (since very specific assumptions and choices have been made), and yet less than a true model (since the underlying physics picture is rather fuzzy). The fragmentation is allowed to proceed as usual, and so is the decay of short-lived particles like rho. Then pairs of identical particles, pi+ say, are considered one by one. The Q value is evaluated, and a shifted (smaller) Q' is found such that the (infinite statistics) ratio of shifted to unshifted Q distributions is given by the requested parametrization in MSTJ(51). (In fact, the distribution dips slightly below unity at Q values outside the Bose enhancement region, from conservation of total multiplicity.) This can be translated into an effective shift of the momenta of the two particles, if one uses as extra constraint that the total three-momentum of each pair be conserved in the CM frame of the event. Only after all pairwise momentum shifts have been evaluated with respect to the original momenta are these momenta actually shifted, for each particle by the sum of evaluated shifts. The total energy of the event is slightly reduced in the process, which is compensated by an overall rescaling of all CM frame momentum vectors. Finally, the decay chain is resumed with more long-lived particles like pi0. Two comments can be made on the approach adopted. First, Bose-Einstein effects are here interpreted almost as a classical force acting on the "final state", rather than as a quantum mechanical phenomenon on the production amplitude. This is not a credo, but just an ansatz to make things manageable. In particular, if an event weight were introduced by the product of all weights for individual pairs, at first sight a more correct procedure, events with higher multiplicities would be weighted up by quite significant amounts. Indeed, the effects are so large that they could not be simply compensated by a modest change of the standard fragmentation parameters. Second, since only pairwise interactions are considered, the effects associated with three or more nearby particles tend to get overestimated (For n identical particles with the same momentum and Lambda = 1, the correct weight is n! but the actually simulated one more like 2^(n(n-1)/2).) Thus the input lambda may have to be chosen smaller than what one wants to get out. This option should therefore be used with caution, and only as a first approximation to what Bose-Einstein effects can mean. The commonblock LUJETS has expanded with time, and can now house 4000 entries. This figure may seem ridiculously large, but actually the previous limit of 2000 was often reached in studies of high-p_T processes at the SSC. The reason for this is that the event record contains not only the final particles, but also all intermediate partons and hadrons, which subsequenty showered, fragmented or decayed. Included are also a wealth of photons coming from pi0 decays; the simplest way of reducing the size of the event record is actually to switch off pi0 decays by MDCY(LUCOMP(111),1) = 0. Also note that some routines, like LUCLUS and LUCELL, use memory after the event record proper as a working area. Still, to change the size of the commonblock, upwards or downwards, is easy: just do a global substitute in the commonblock and change the MSTU(4) value to the new number. If more than 10000 lines are to be used, the packing of colour information should also be changed, see MSTU(5). The program contains space so that additional new particles may be introduced. Although not completely trivial, this should not be beyond the ability of an ordinary user. Basically, three steps are involved. First, a mechanism of production has to be introduced. This production may well take place in an external program, like PYTHIA or some user-written correspondence, where matrix elements are used to select the hard process. In this case the new particle already exists in the LUJETS commonblock when JETSET is called. A new particle, meson, baryon or glueball, may also be a part of the fragmentation process, in which case LUKFDI would have to be suitably modified. The particle might also appear as a decay product from some already existing particle, and then the decay data in /LUDAT3/ would have to be expanded; conceivably also LUDECY would be affected. The second step is to teach to program to recognize the new particle. IF a KF code in the range 41 to 80 is used, this is automatically taken care of, and in particular the compressed code KC coincides with KF. If a whole sequence of particles is to be introduced, with KF codes paralleling that of ordinary mesons/baryons (a supersymmetric "meson" multiplet, made of a squark plus an antiquark, say), then LUCOMP must be modified to include a mapping from these KF values to currently unused KC ones, like the range 401 - 500. It is the presence of such a mapping that the program uses to accept a given KF code as bona fide. The third and final step is to define the properties of this new particle. Thus particle charge information must be given in KCHG, mass, width and lifetime in PMAS, particle name in CHAF, and decay data in the MDCY, MDME, BRAT and KFDP arrays. This process is most conveniently carried out by using LUUPDA(1) to produce a table of particle data, which can then be modified by the user and read back in with LUUPDA(2). Note that the particle data is to be introduced for the compressed code KC, not for KF proper. ______________________________________________________________________ 2.10. Examples What every user has to know something about is the commonblock LUJETS COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5) A short summary follows; for details see section 2.2. The complete event record is stored in LUJETS. This normally includes the original jets, the particles these fragment into and the subsequent decay chains. The number N gives the number of lines (1 through N) in the K, P and V matrices that are actually used to store the current event. For line I, K(I,1) gives information about current status. In particular, anything which has 1 <= K(I,1) <= 10 represents a particle or parton that is either considered stable or has not yet been treated, whereas K(I,1) >= 11 corresponds to partons that have fragmented or branched, particles that have decayed and an assortment of special purpose lines. K(I,2) gives a code for what parton or particle we are dealing with, e.g. 1 = d, 2 = u, 11 = e-, 12 = nu_e, 21 = g, 22 = gamma, 2103 = ud_1 diquark, 111 = pi0, 211 = pi+, 213 = rho+, 2212 = p. Antiparticles, where existing, are given by the negative number, -1 = dbar, -211 = pi-. This code is often referred to as KF code, and is described in section 2.1. Codes K(I,3) - K(I,5) contain event history information, line number where mother are stored, line number range of daughters or colour flow information. The P vector contains momentum information, with P(I,1), P(I,2) and P(I,3) giving the three-momentum, P(I,4) the energy and P(I,5) the mass, in GeV/c, GeV and GeV/c^2, respectively. The components V(I,1), V(I,2) and V(I,3) give the production vertex position of the particle, and V(I,4) the production time, all with respect to a primary vertex situated in the origin at time 0. The space components are given in mm, and the time one in mm/c = 3.33*10^-12 s. The final component V(I,5) gives the invariant lifetime of a particle, again in mm/c. Neglecting bending in a magnetic field etc., the decay vertex V' of the particle is then given by V'(j) = V(I,j) + V(I,5)*P(I,j)/P(I,5), j = 1 - 4. For a stable particle V(I,5) = 0. For an unstable, V(I,5) is selected even if the decay is not actually carried out because it takes place outside the allowed region, see section 2.9. A 10 GeV u quark jet going out along the +z axis is generated with CALL LU1ENT(0,2,10.,0.,0.) Note that such a single jet is not required to conserve energy, momentum or flavour. In the generation scheme, particles with negative p_z are produced as well, but these are automatically rejected unless MSTJ(3) = -1. While frequently used in former days, the one-jet generation option is not of much current interest. In e.g. a leptoproduction event a typical situation could be a u quark going out in the +z direction and a ud_0 target remnant essentially at rest. The simplest procedure is probably to treat it in the CM frame and boost it to the lab frame afterwards. Hence, if the CM energy is 20 GeV and the boost beta = 0.996 (corresponding to x_B = 0.045) CALL LU2ENT(0,2,2101,20.) CALL LUROBO(0.,0.,0.,0.,0.996) The jets could of course also be defined and allowed to fragment in the lab frame with CALL LU1ENT(-1,2,223.15,0.,0.) CALL LU1ENT(2,12,0.6837,3.1416,0.) CALL LUEXEC Note here that the target diquark is required to move in the backwards direction with (E-p_z) = m_proton*(1-x_B) to obtain the correct invariant mass for the system. This is, however, only an artefact of using a fixed diquark mass to represent a varying target remnant mass, and is of no importance for the fragmentation. If one wants a nicer-looking event record, it is possible to use the following CALL LU1ENT(-1,2,223.15,0.,0.) MSTU(10)=1 P(2,5)=0.938*(1.-0.045) CALL LU1ENT(2,2101,0.,0.,0.) MSTU(10)=2 CALL LUEXEC A 30 GeV uubarg event with E_u = 8 GeV and E_ubar = 14 GeV is simulated with CALL LU3ENT(0,2,21,-1,30.,2.*8./30.,2.*14./30.) The event will be given in a standard orientation with the u quark along the +z axis and the ubar in the -z, +x quadrant. Note that the flavours of the three partons have to be given in the order they are found along a string, if string fragmentation options are to work. Also note that, for three-jet events, and particularly four-jet ones, not all setups of kinematical variables x lie within the kinematically allowed regions of phase space. All commonblock variables can obviously be changed by including the corresponding commonblock in the user-written main program. Alternatively, the routine LUGIVE can be used to feed in values, with some additional checks on array bounds then performed. A call CALL LUGIVE('MSTJ(21)=3;PMAS(C663,1)=89.;CHAF(401)=funnyino;'// &'PMAS(21,4)=') will thus change the value of MSTJ(21) to 3, the value of PMAS(LUCOMP(663),1) = PMAS(136,1) to 89., the value of CHAF(401) to 'funnyino', and give the current value of PMAS(21,4). Since old and new values of parameters changed are written to output, this may offer a convenient way of documenting non-default values used in a given run. On the other hand, if a variable is changed back and forth frequently, the resulting voluminous output may be undesirable, and a direct usage of the commonblocks is then to be recommended (the output can also be switched off, see MSTU(13)). A general rule of thumb is that none of the physics routines (LUSTRF, LUINDF, LUDECY, etc.) should ever be called directly, but only via LUEXEC. This routine may be called repeatedly for one single event. At each call only those entries that are allowed to fragment or decay, and have not yet done so, are treated. Thus CALL LU2ENT(1,1,-1,20.) ! fill 2 jets without fragmenting MSTJ(1)=0 ! inhibit jet fragmentation MSTJ(21)=0 ! inhibit particle decay MDCY(LUCOMP(111),1)=0 ! inhibit pi0 decay CALL LUEXEC ! will not do anything MSTJ(1)=1 ! CALL LUEXEC ! jets will fragment, but no decays MSTJ(21)=2 ! CALL LUEXEC ! particles decay, except pi0 CALL LUEXEC ! nothing new can happen MDCY(LUCOMP(111),1)=1 ! CALL LUEXEC ! pi0:s decay A partial exception to the rule above is LUSHOW. Its main application is for internal use by LUEEVT and LUDECY, and by other Lund family programs like PYTHIA, but it can also be directly called by the user. Note that a special format for storing colour flow information in K(I,4) and K(I,5) must then be used. For simple cases, the LU2ENT can be made to take care of that automatically, by calling with the first argument negative. CALL LU2ENT(-1,1,-2,40.) ! store dubar with colour flow CALL LUSHOW(1,2,40.) ! shower partons CALL LUEXEC ! subsequent fragmentation/decay It is always good practice to list one or a few events during a run to check that the program is working as intended. With CALL LULIST(1) all particles will be listed and in addition total charge, momentum and energy of stable entries will be given. For string fragmentation these quantities should be conserved exactly (up to machine precision errors), and the same goes when running independent fragmentation with one of the momentum conservation options. LULIST(1) gives a format that comfortably fits on an 80 column screen, at the price of not giving the complete story. With LULIST(2) a more extensive listing is obtained, and LULIST(3) also gives vertex information. Further options are available, like LULIST(12), which gives a list of particle data. An event, as stored in the LUJETS commonblock, will contain the original jets and the whole decay chain, i.e. also particles which subsequently decayed. If parton showers are used, the amount of parton information is also considerable: first the on-shell partons before showers have been considered, then a KS = 22 line with total energy of the showering subsystem, after that the complete shower history treelike structure, starting off with the same initial partons (now off-shell), and finally the endproducts of the shower rearranged along the string directions. This detailed record is useful in many connections, but if one only wants to retain the final particles, superfluous information may be removed with LUEDIT. Thus e.g. CALL LUEDIT(2) will leave you with the final charged and neutral particles, except for neutrinos. The information in LUJETS may be used directly to study an event. Some useful additional quantities derived from these, such as charge and rapidity, may easily be found via the KLU and PLU functions. Thus electric charge = PLU(I,6) (as integer, 3*charge = KLU(I,6)) and true rapidity y with respect to the z axis = PLU(I,17). Event analysis routines for properties of the event as-a-whole include sphericity, thrust, cluster (2 different algorithms), jet masses and Fox-Wolfram moments. These routines may be called directly after the event generation CALL LU2ENT(0,5,-5,40.) CALL LUTHRU(THR,OBL) and then all stable, final particles (except neutrinos) are used in the analysis. This may be changed with the MSTu(41), or by explicitly using LUEDIT, or by hand setting K(I,1) values for unwanted particles before the analysis. A number of utility (MSTU, PARU) and physics (MSTJ, PARJ) switches and parameters are available in commonblock LUDAT1. All of these have sensible default values. Particle data is stored in commonblocks LUDAT2, LUDAT3 and LUDAT4. Note that the data in the arrays KCHG, PMAS, MDCY and CHAF is not stored by KF code, but by the compressed code KC. This code is not to be learnt by heart, but instead accessed via the conversion function LUCOMP, KC = LUCOMP(KF). In the particle tables, the following particles are considered stable: the photon, e+-, mu+-, pi+-, K+-, K_L0, p, pbar, n, nbar and all the neutrinos. It is, however, always possible to inhibit the decay of any given particle by putting the corresponding MDCY value zero or negative, e.g. MDCY(LUCOMP(310),1) = 0 makes K_S0 and MDCY(LUCOMP(3122),1) = 0 Lambda stable. It is also possible to select stability based on the average lifetime (see MSTJ(22)), or based on whether the decay takes place within a given spherical or cylindrical volume around the origin. This is described in more detail in section 2.9. The Field-Feynman jet model [Fie78] is available in the program by changing the following values: MSTJ(1) = 2 (independent fragmentation), MSTJ(3) = -1 (retain particles with p_z < 0; is not mandatory), MSTJ(11) = 2 (choice of longitudinal fragmentation function, with the a parameter stored in PARJ(51) - PARJ(53)), MSTJ(12) = 0 (no baryon production), MSTJ(13) = 1 (give endpoint quarks p_T as quarks created in the field), MSTJ(24) = 0 (no mass broadening of resonances), PARJ(2)=0.5 (s/u ratio for the production of new qqbar pairs), PARJ(11) = PARJ(12) = 0.5 (probability for mesons to have spin 1) and PARJ(21) = 0.35 (width of Gaussian transverse momentum distribution). In addition only d, u and s single quark jets may be generated following the FF recipe. Today the FF "standard jet" concept is probably dead and buried, so the numbers above should more be taken as an example of the flexibility of the program, than as something to apply in practice. A wide range of independent fragmentation options are implemented, to be accessed with the master switch MSTJ(1) = 2. In particular, with MSTJ(2) = 1 a gluon jet is assumed to fragment like a random d, dbar, u, ubar, s or sbar jet, while with MSTJ(2) = 3 the gluon is split into a ddbar, uubar or ssbar pair of jets sharing the energy according to the Altarelli-Parisi splitting function. Whereas energy, momentum and flavour is not explicitly conserved in independent fragmentation, a number of options are available in MSTJ(3) to ensure this "post facto", e.g. MSTJ(3) = 1 will boost the event to ensure momentum conservation and then (in the CM frame) rescale momenta by a common factor to obtain energy conservation, whereas MSTJ(3)=4 rather uses a method of stretching the jets in longitudinal momentum along the respective jet axis to keep angles between jets fixed. ______________________________________________________________________ 2.11. Translation to/from Standard Commonblock Recently a set of commonblocks has been developed and agreed on within the framework of the 1989 LEP physics study, see [QEG89]. This standard defines an event record structure which should make the interfacing of different event generators much more simple. It would be a major work to rewrite JETSET to agree with this standard event record structure. More importantly, the standard only covers quantities which can be defined unambigously, i.e. which are independent of the particular program used. There are thus no provisions for the need for colour flow information in models based on string fragmentation, etc., so the standard commonblocks would anyway have to be supplemented with additional event information. For the moment, the approach adopted is therefore to retain the LUJETS event record, but supply a routine LUHEPC which can convert to or from the standard event record. Due to somewhat different content in the two records, some ambiguities do exist in the translation procedure. LUHEPC has therefore to be used with some judgment. In this section, the new standard event structure is first presented, i.e. the most important points in [QEG89] are recapitulated. Thereafter the conversion routine is described, with particular attention to ambiguities and limitations. The standard event record is stored in two commonblocks. The second of these is specifically intended for spin information. Since JETSET never (explicitly) makes use of spin information, this latter commonblock is never addressed. In order to make the components of the standard more distinguishable in user programs, the three characters HEP (for High Energy Physics) have been chosen to be a part of all names. PARAMETER (NMXHEP=2000) COMMON/HEPEVT/NEVHEP,NHEP,ISTHEP(NMXHEP),IDHEP(NMXHEP), &JMOHEP(2,NMXHEP),JDAHEP(2,NMXHEP),PHEP(5,NMXHEP),VHEP(4,NMXHEP) Purpose: to contain an event record in a Monte Carlo-independent format. NMXHEP: maximum numbers of entries (partons/particles) that can be stored in the commonblock. The default value of 2000 can be changed via the parameter construction. In the translation, it is checked that this value is not exceeded. NEVHEP: is normally the event number, but may have special meaning, according to description below. > 0 : event number, sequentially increased by 1 for each call to the main event generation routine, starting with 1 for the first event generated. = 0 : for a program which does not keep track of event numbers, like JETSET. = -1 : special initialization record; not used by JETSET. = -2 : special final record; not used by JETSET. NHEP: the actual number of entries stored in current event. These are found in the first NHEP positions of the respective arrays below. Index IHEP, 1 <= IHEP <= NHEP, is below used to denote a given entry. ISTHEP(IHEP): status code for entry IHEP, with following meanings: = 0 : null entry. = 1 : an existing entry, which has not decayed or fragmented. This is the main class of entries which represents the "final state" given by the generator. = 2 : an entry which has decayed or fragmented and therefore is not appearing in the final state, but is retained for event history information. = 3 : a documentation line, defined separately from the event history. This could include the two incoming reacting particles, etc. = 4 - 10 : undefined, but reserved for future standards. = 11 - 200 : at the disposal of each model builder for constructs specific to his program, but equivalent to a null line in the context of any other program. = 201 - : at the disposal of users, in particular for event tracking in the detector. IDHEP(IHEP) : particle identity, according to the Particle Data Group standard. The four additional codes 91 - 94 have been introduced to make the event history more legible, see section 2.1 and the MSTU(16) description. JMOHEP(1,IHEP) : pointer to the position where the mother is stored. The value is 0 for initial entries. JMOHEP(2,IHEP) : pointer to position of second mother. Normally only one mother exist, in which case the value 0 is to be used. In JETSET, entries with codes 91 - 94 are the only ones to have two mothers. The flavour contents of these object, as well as details of momentum sharing, have to be found by looking at the mother partons, i.e. the two partons in positions JMOHEP(1,IHEP) and JMOHEP(2,IHEP) for a cluster or a shower system, and the range JMOHEP(1,IHEP) - JMOHEP(2,IHEP) for a string or an independent fragmentation parton system. JDAHEP(1,IHEP) : pointer to the position of the first daughter. If an entry has not decayed, this is 0. JDAHEP(2,IHEP) : pointer to the position of the last daughter. If an entry has not decayed, this is 0. It is assumed that daughters are stored sequentially, so that the whole range JDAHEP(1,IHEP) - JDAHEP(2,IHEP) contains daughters. This should be done also when only one daughter is present, like in K^0 -> K_S^0 decays. Normally daughters are stored after mothers, but in backwards evolution of initial state radiation the opposite may appear, i.e. that mothers are found below the daughters they branch into. Also, the two daughters need then not appear one after the other, but may be separated in the event record. PHEP(1,IHEP) : momentum in the x direction, in GeV/c. PHEP(2,IHEP) : momentum in the y direction, in GeV/c. PHEP(3,IHEP) : momentum in the z direction, in GeV/c. PHEP(4,IHEP) : energy, in GeV. PHEP(5,IHEP) : mass, in GeV/c**2. For spacelike partons, it is allowed to use a negative mass, according to PHEP(5,IHEP) = -sqrt(-m**2). VHEP(1,IHEP) : production vertex x position, in mm. VHEP(2,IHEP) : production vertex y position, in mm. VHEP(3,IHEP) : production vertex z position, in mm. VHEP(4,IHEP) : production time, in mm/c (= 3.33*10**(-12) s). This completes the brief description of the standard. In JETSET, the routine LUHEPC is provided as an interface. SUBROUTINE LUHEPC(MCONV) Purpose: to convert between the LUJETS event record and the HEPEVT event record. MCONV : direction of conversion. = 1 : translate the current LUJETS record into the HEPEVT one, while leaving the original LUJETS one unaffected. = 2 : translate the current HEPEVT record into the LUJETS one, while leaving the original HEPEVT one unaffected. The conversion of momenta is trival: it is just a matter of exchanging the order of the indices. The vertex information is but little more complicated; the extra fifth component present in LUJETS can be easily reconstructed from other information for particles which have decayed. (Some of the advanced features made possible by this component, like the possibility to consider decays within expanding spatial volumes in subsequent LUEXEC calls, can not be used if the record is translated back and forth, however.) Also, the particle codes K(I,2) and IDHEP(I) are identical, since they are both based on the PDG codes. The remaining, nontrivial areas deal with the status codes and the event history. In moving from LUJETS to HEPEVT, information on colour flow is lost. On the other hand, the position of a second mother, if any, has to be found; this only affects lines with K(I,2) = 91 - 94. Also, for lines with K(I,1) = 13 or 14, the daughter pointers have to be found. By and large, however, the translation from LUJETS to HEPEVT should cause little problem, and there should never be any need for user intervention. (We assume that JETSET is run with the default MSTU(16) = 1, else some discrepancies compared to the proposed standard event history description will be present.) In moving from HEPEVT to LUJETS, information on a second mother is lost. Any codes IDHEP(I) not equal to 1, 2 or 3 are translated into K(I,1) = 0, and so all entries with K(I,1) >= 30 are effectively lost in a translation back and forth. All entries with IDHEP(I) = 2 are translated into K(I,1) = 11, and so entries of type K(I,1) = 12, 13 or 14 or 15 are never found. There is thus no colour flow information available for partons which have fragmented. For partons with IDHEP(I) = 1, i.e. which have not fragmented, an attempt is made to subdivide the partonic system into colour singlets, as required for subsequent string fragmentation. To this end, it is assumed that partons are stored sequentially along strings. Normally, a string would then start at a q (qbar) or qbarqbar (qq) entry, cover a number of intermediate gluons, and end at a qbar (q) or qq (qbarqbar) entry. Particles could be interspersed in this list with no adverse effects, i.e. a u-g-gamma- ubar sequence would be interpreted as a u-g-ubar string plus an additional photon. A closed gluon loop would be assumed to be made up of a sequential listing of the gluons, with the string continuing from the last gluon up back to the first one. Contrary to the previous, open string case, the appearance of any particle but a gluon would therefore signal the end of the gluon loop. For example, a g-g-g-g sequence would be interpreted as one single four-gluon loop, while a g-g-gamma-g-g sequence would be seen as composed of two two-gluon systems. If these interpretations, which are not unique, are not to the liking of the user, it is up to him/her to correct them, e.g. by using LUJOIN to tell exactly which partons should be joined, in which sequence, to give a string. Calls to LUJOIN (or the equivalent) are also necessary if LUSHOW is to be used to have some partons develop a shower. For practical applications, one should note that e+e- events, which have been allowed to shower but not to fragment, do have partons arranged in the order assumed above, so that a translation to HEPEVT and back does not destroy the possibility to perform fragmentation by a simple LUEXEC call. Also the hard interactions in PYTHIA fulfil this condition, while problems may appear in the multiple interaction scenario, where several closed gg loops may appear directly following after each other, and thus would be interpreted as a single multigluon loop after translation back and forth. ********************************************************************** 3. The e+e- Routines The routines in this section are more specific than the ones in section 2: aim is taken on applications to e+e- annihilation events, in the continuum or on an onium resonance. ______________________________________________________________________ 3.1. e+e- Continuum Event Generation The only routine a normal user will call is LUEEVT. The other routines listed below, as well as LUSHOW (see section 2.4), are called by LUEEVT. SUBROUTINE LUEEVT(KFL,ECM) Purpose: to generate a complete event e+e- -> gamma/Z0 -> qqbar -> parton shower -> hadrons according to QFD and QCD cross-sections. As an alternative to parton showers, second order matrix elements are available for qqbar + qqbarg + qqbargg + qqbarq'qbar' production. KFL : flavour of events generated. = 0 : mixture of all allowed flavours according to relevant probabilities. = 1 - 8 : primary quarks are only of the specified flavour KFL. ECM : total CM energy of system. Remark: Each call generates one event, which is independent of preceding ones, with one exception, as follows. If radiative corrections are included, the shape of the hard photon spectrum is recalculated only with each LUXTOT call, which normally is done only if KFL, ECM or MSTJ(102) is changed. A change of e.g. the Z0 mass in midrun has to be followed either by a user call to LUXTOT or by an internal call forced e.g. by putting MSTJ(116)=3. SUBROUTINE LUXTOT(KFL,ECM,XTOT) Purpose : to calculate the total hadronic cross-section, including quark thresholds, weak, beam polarization and QCD effects and radiative corrections. In the process, variables necessary for the treatment of hard photon radiation are calculated and stored. KFL, ECM : as for LUEEVT. XTOT : the calculated total cross-section in nb. SUBROUTINE LURADK(ECM,MK,PAK,THEK,PHIK,ALPK) Purpose: to describe initial state hard photon radiation. SUBROUTINE LUXKFL(KFL,ECM,ECMC,KFLC) Purpose: to generate the primary quark flavour in case this is not specified by user. SUBROUTINE LUXJET(ECM,NJET,CUT) Purpose: to determine the number of jets (2, 3 or 4) to be generated within the kinematically allowed region (characterized by CUT = y_cut) in the matrix element approach; to be chosen such that all probabilities are between 0 and 1. SUBROUTINE LUX3JT(NJET,CUT,KFL,ECM,X1,X2) Purpose: to generate the internal momentum variables of a three-jet event, qqbarg, according to first or second order QCD matrix elements. SUBROUTINE LUX4JT(NJET,CUT,KFL,ECM,KFLN,X1,X2,X4,X12,X14) Purpose: to generate the internal momentum variables for a four-jet event, qqbargg or qqbarqqbar, according to second order QCD matrix elements. SUBROUTINE LUXDIF(NC,NJET,KFL,ECM,CHI,THE,PHI) Purpose: to describe the angular orientation of the jets. In first order QCD the complete QED or QFD formulae are used; in second order three-jets are assumed to have the same orientation as in first, and four-jets are approximated by three-jets. ______________________________________________________________________ 3.2. A Routine for "onium" Decay In LUONIA we have implemented the decays of heavy onia resonances into three gluons or two gluons plus a photon, which are the dominant non-backgroundlike decays of Upsilon, and also would have been it for a reasonably light top. With the present mass limits, actually weak decays are expected to dominate, as is already implemented in the ordinary LUDECY treatment of toponium decay. The inclusion of parton showering for Upsilon decay is probably overkill, but is there for completeness. SUBROUTINE LUONIA(KFL,ECM) Purpose: to simulate the process e+e- -> gamma -> 1-- "onium" resonance -> ggg or gggamma -> shower -> hadrons. KFL : the flavour of the quark giving rise to the resonance. = 0 : generate ggg events alone. = 1 - 8 : generate ggg and gggamma events in mixture determined by the charge-squared of flavour KFL. Normally KFL = 5 or 6. ______________________________________________________________________ 3.3. The Commonblock Variables The status codes and parameters relevant for the e+e- routines are found in the commonblock LUDAT1. This commonblock also contains more general status codes and parameters, which were described in section 2.7. COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200) Purpose: to give access to a number of status codes and parameters regulating the performance of the e+e- event generation routines. MSTJ(101) : (D=5) gives the type of QCD corrections used for continuum events. = 0 : only qqbar events are generated. = 1 : qqbar + qqbarg events are generated according to first order QCD. = 2 : qqbar + qqbarg + qqbargg + qqbarqqbar events are generated according to second order QCD. = 3 : qqbar + qqbarg + qqbargg + qqbarqqbar events are generated, but without second order corrections to the three-jet rate. = 5 : a parton shower is allowed to develop from an original qqbar pair, see MSTJ(41) - MSTJ(49) for details. = -1 : only qqbarg events are generated (within same matrix element cuts as for =1). Since the change in flavour composition from mass cuts or radiative corrections is not taken into account, this option is not intended for quantitative studies. = -2 : only qqbargg and qqbarqqbar events are generated (as for =2). The same warning as for =-1 applies. = -3 : only qqbargg events are generated (as for =2). The same warning as for =-1 applies. = -4 : only qqbarqqbar events are generated (as for =2). The same warning as for =-1 applies. Note 1: MSTJ(101) is also used in LUONIA, with <= 4 : ggg + gammagg events are generated according to lowest order matrix elements. >= 5 : a parton shower is allowed to develop from the original ggg or gggamma configuration, see MSTJ(41) - MSTJ(49) for details. Note 2: The default values of fragmentation parameters have been chosen to work well with the default parton shower approach above. If any of the other options are used, or if the parton shower is used in non-default mode, it may be necessary to retune fragmentation parameters. As an example, we note that the second order matrix element approach (MSTJ(101) = 2) at PETRA/PEP energies gives a better description with the a and b parameters of the symmetric fragmentation function set to a = PARJ(41) = 1, b = PARJ(42) = 0.7, and the width of the transverse momentum distribution to sigma = PARJ(21) = 0.40. In principle, one also ought to change the joining parameter to PARJ(33) = PARJ(35) = 1.1 to preserve a flat rapidity plateau, but if this should be forgotten, it does not make too much of a difference. For applications at TRISTAN or LEP, one must expect to have to change the matrix element approach parameters even more, to make up for additional soft gluon effects not covered in this approach. MSTJ(102) : (D=2) inclusion of weak effects (Z0 exchange) for flavour production, angular orientation, cross-sections and initial state photon radiation in continuum events. = 1 : QED, i.e. no weak effects are included. = 2 : QFD, i.e. including weak effects. = 3 : as =2, but at initialization in LUXTOT the Z0 width is calculated from sin^2(theta_W), alpha_em and Z0 and quark masses (including bottom and top threshold factors for MSTJ(103) odd), assuming three full generations, and the result is stored in PARJ(124). MSTJ(103) : (D=7) mass effects in continuum matrix elements, in the form MSTJ(103) = M_1 + 2*M_2 + 4*M_3, where M_i = 0 if no mass effects and M_i = 1 if mass effects should be included. Here; M_1 : threshold factor for new flavour production according to QFD result; M_2 : gluon emission probability (only applies for abs(MSTJ(101)) <= 1, otherwise no mass effects anyhow); M_3 : angular orientation of event (only applies for abs(MSTJ(101)) <= 1 and MSTJ(102)=1, otherwise no mass effects anyhow). MSTJ(104) : (D=5) number of allowed flavours, i.e. flavours that can be produced in a continuum event if the energy is big enough. A change to 6 makes top production allowed above the threshold, etc. Note that in qqbarqqbar events only the five first flavours are allowed in the secondary pair, produced by a gluon breakup. MSTJ(105) : (D=1) fragmentation and decay in LUEEVT and LUONIA calls. = 0 : no LUEXEC calls, i.e. only matrix element and/or parton shower treatment. = 1 : LUEXEC calls are made to generate fragmentation and decay chain. = -1 : no LUEXEC calls and no collapse of small jet systems into one or two particles (in LUPREP). MSTJ(106) : (D=1) angular orientation in LUEEVT and LUONIA. = 0 : standard orientation of events, i.e. q along +z axis and qbar along -z axis or in xz plane with p_x > 0 for continuum events, and g_1g_2g_3 or gammag_2g_3 in xz plane with g_1 or gamma along the +z axis for continuum events. = 1 : random orientation according to matrix elements. MSTJ(107) : (D=0) radiative corrections to continuum events. = 0 : no radiative corrections. = 1 : initial state radiative corrections (including weak effects for MSTJ(102) = 2 or 3). MSTJ(108) : (D=2) calculation of alpha_strong for matrix element alternatives. The MSTU(111) and PARU(112) values are automatically overwritten in LUEEVT or LUONIA calls accordingly. = 0 : fixed alpha_strong value as given in PARU(111). = 1 : first order formula is always used, with Lambda_QCD given by PARJ(121). = 2 : first or second order formula is used, depending on value of MSTJ(101), with Lambda_QCD given by PARJ(121) or PARJ(122). MSTJ(109) : (D=0) gives a possibility to switch from QCD matrix elements to some alternative toy models. Is not relevant for shower evolution, MSTJ(101) = 5, where one can instead use MSTJ(49). = 0 : standard QCD scenario. = 1 : a scalar gluon model. Since no second order corrections are available in this scenario, one can only use this with MSTJ(101) = 1 or -1. Also note that the event-as-a-whole angular distribution is for photon exchange only (i.e. no weak effects), and that no higher order corrections to the total cross-section are included. = 2 : an Abelian vector gluon theory, with the colour factors C_F = 1 (= 4/3 in QCD), N_C = 0 (= 3 in QCD) and T_R = 3 n_f (= n_f/2 in QCD). If one selects alpha_Abelian = (4/3) * alpha_QCD, the three-jet cross-section will agree with the QCD one, and differences are to be found only in four-jets. The MSTJ(109) = 2 option has to be run with the defaults MSTJ(110) = 1 and MSTJ(111) = 0; if need be, the latter variables will be overwritten by the program. Warning: second order corrections give a large negative contribution to the three-jet cross-section, so large that the whole scenario is of doubtful use. In order to make the second order options work at all, the three-jet cross-section is here by hand set exactly equal to zero for MSTJ(101) = 2. It is here probably better to use the option MSTJ(101) = 3, although this is not a consistent procedure either. MSTJ(110) : (D=1) choice of second order contributions to the three-jet rate. = 1 : the GKS second order matrix elements, i.e. the old JETSET standard. = 2 : the Zhu parametrization of the ERT matrix elements, based on the program of Kunszt and Ali, i.e. in historical sequence ERT/Kunszt/Ali/Zhu. The parametrization is available for y = 0.01, 0.02, 0.03, 0.04 and 0.05. Values outside this range are put at the nearest border, while values inside this range are given by linear interpolation between the two nearest points. Since this procedure is rather primitive, one should try to work at one of the values given above. Note that no Abelian CQD parametrization is available for this option. MSTJ(111) : (D=0) use of optimized perturbation theory for second order matrix elements (it can also be used for first order matrix elements, but here it only corresponds to a trivial rescaling of the alpha_strong argument). = 0 : no optimization procedure; i.e. Q^2 = E_CM^2. = 1 : an optimized Q^2 scale is chosen as Q^2 = y' * E_CM^2, where y' = PARJ(128) for the total cross-section R factor. while y' = PARJ(129) for the three- and four-jet rates, This y' value enters via the alpha_strong, and also via a term proportional to alpha_strong**2*ln(y'). Some constraints are imposed; thus the optimized 'three-jet' contribution to R is assumed positive (for PARJ(128)), the total three-jet rate piece of second order) contribution to R is assumed positive, is not allowed to be negative (for PARJ(129)), etc. However, there is no guarantee that the differential three-jet cross-section is not negative (and truncated to 0) somewhere (this can also happen with y' = 1, but is then less frequent). The actually obtained y' values are stored in PARJ(168) and PARJ(169), respectively. If an optimized Q^2 scale is used, then the Lambda_QCD (and alpha_strong) should also be changed. With the value y' = 0.002, it has been shown [Bet89] that a Lambda_QCD = 0.100 GeV gives a reasonable agreement; the parameter to change is PARJ(122) for a second order running alpha_strong. Note that, since the optimized Q^2 scale is sometimes below the charm threshold, the effective number of flavours used in alpha_strong may well be 4 only. If one feels that it is still appropriate to use 5 flavours (one choice might be as good as the other), it is necessary to put MSTU(113) = 5. MSTJ(115) : (D=1) documentation of continuum or onium events, in increasing order of completeness. = 0 : only the parton shower, the fragmenting partons and the generated hadronic system are stored in the LUJETS commonblock. = 1 : also a radiative photon is stored (for continuum events). = 2 : also the original e+e- are stored (with KS=21). = 3 : also the gamma or gamma/Z0 exchanged for continuum events, the onium state for resonance events is stored (with KS=21). MSTJ(116) : (D=1) initialization of total cross-section and radiative photon spectrum in LUEEVT calls. = 0 : never; can not be used together with radiative corrections. = 1 : calculated at first call and then whenever KFL or MSTJ(102) is changed or ECM is changed by more than PARJ(139). = 2 : calculated at each call. = 3 : everything is reinitialized in next call, but MSTJ(116) is afterwards automatically put =1 for use in subsequent calls. MSTJ(119) : (I) check on need to reinitialize LUXTOT. MSTJ(120) : (R) type of continuum event generated with matrix element option (with the shower one, the result is always =1). = 1 : qqbar. = 2 : qqbarg. = 3 : qqbargg from Abelian (QED-like) graphs in matrix element. = 4 : qqbargg from non-Abelian (i.e. containing three-gluon coupling) graphs in matrix element. = 5 : qqbarqqbar. MSTJ(121) : (R) flag set if a negative differential cross-section was encountered in the latest LUX3JT call. Events are still generated, but maybe not quite according to the distribution one would like (the rate is set zero in the regions of negative cross-section, and the differential rate in the regions of positive cross-section is rescaled to give the 'correct' total three-jet rate). PARJ(121) : (D=1.5 GeV) Lambda value used in first order calculation of alpha_strong in the matrix element alternative. PARJ(122) : (D=0.5 GeV) Lambda values used in second order calculation of alpha_strong in the matrix element alternative. PARJ(123) : (D=91.2 GeV) mass of Z0 as used in propagators for QFD case. PARJ(124) : (D=2.4 GeV) width of Z0 as used in propagators for QFD case. Overwritten at initialization if MSTJ(102) = 3. PARJ(125) : (D=0.02) y_cut, minimum scaled invariant mass-squared of any two partons in 3- or 4-jet events; the main user-controlled matrix element cut. PARJ(126) provides an additional constraint. For each new event, it is additionally checked that the total three- plus four-jet fraction does not exceed unity; if so the effective y cut will be dynamically increased. The actual y cut value is stored in PARJ(150), event by event. PARJ(126) : (D=2. GeV) minimum invariant mass of any two partons in 3- or 4-jet events; a cut in addition to the one above, mainly for the case of a radiative photon lowering the hadronic CM energy significantly. PARJ(127) : (D=1. GeV) is used as a safety margin for small colour singlet jet systems, cf. PARJ(32), specifically qqbar masses in qqbarqqbar 4-jet events and gg mass in onium gammagg events. PARJ(128) : (D=0.25) optimized Q^2 scale for the QCD R (total rate) factor for the MSTJ(111) = 1 option is given by Q^2 = y' * E_CM^2, where y' = PARJ(128). For various reasons the actually used y' value may be raised from the nominal one; while PARJ(128) gives the nominal value PARJ(168) gives the actual one for the current event. PARJ(129) : (D=0.002) optimized Q^2 scale for the three- and four-jet rate for the MSTJ(111) = 1 option is given by Q^2 = y' * E_CM^2, where y' = PARJ(129). For various reasons the actually used y' value may be raised from the nominal one; while PARJ(129) gives the nominal value PARJ(169) gives the actual one for the current event. The default value is in agreement with the studies of Bethke [Bet89]. PARJ(131), PARJ(132) : (D=2*0.) longitudinal polarizations P_L+ and P_L- of incoming e+ and e-. PARJ(133) : (D=0.) transverse polarization P_T = (P_T+*P_T-)^*(1/2) with P_T+ and P_T- transverse polarizations of incoming e+ and e-. PARJ(134) : (D=0.) mean of transverse polarization directions of incoming e+ and e-, Deltaphi = (phi_+ + phi_-)/2, with phi azimuthal angle of polarization, leading to a shift in the phi distribution of jets by Deltaphi. PARJ(135) : (D=0.01) minimum photon energy fraction (of beam energy) in initial state radiation; should normally never be changed (if lowered too much, the fraction of events containing a radiative photon will exceed unity, leading to problems). PARJ(136) : (D=0.99) maximum photon energy fraction (of beam energy) in initial state radiation; may be changed to reflect actual trigger conditions of a detector (but must always be larger than PARJ(135)). PARJ(139) : (D=0.2 GeV) maximum deviation of ECM from the corresponding value at last LUXTOT call, above which a new call is made if MSTJ(116) = 1. PARJ(141) : (R) value of R, the ratio of continuum cross-section to the lowest order muon pair production cross-section, as given in massless QED (i.e. three times the sum of active quark charges-squared, possibly modified for polarization). PARJ(142) : (R) value of R including quark mass effects (for MSTJ(102)=1) and/or weak propagator effects (for MSTJ(102)=2). PARJ(143) : (R) value of R as PARJ(142), but including QCD corrections as given by MSTJ(101). PARJ(144) : (R) value of R as PARJ(143), but additionally including corrections from initial state photon radiation (if MSTJ(107)=1). Since the effects of heavy flavour thresholds are not simply integrable, the initial value of PARJ(144) is updated during the course of the run to improve accuracy. PARJ(145) - PARJ(148) : (R) absolute cross-sections in nb as for the cases PARJ(141) - PARJ(144) above. PARJ(150) : (R) current effective matrix element cutoff y_cut, as given by PARJ(125), PARJ(126) and the requirements of having non-negative cross-sections for two-, three- and four-jet events. Not used in parton showers. PARJ(151) : (R) value of CM energy ECM at last LUXTOT call. PARJ(152) : (R) current first-order contribution to the three-jet fraction; modified by mass effects. Not used in parton showers. PARJ(153) : (R) current second-order contribution to the three-jet fraction; modified by mass effects. Not used in parton showers. PARJ(154) : (R) current second-order contribution to the four-jet fraction; modified by mass effects. Not used in parton showers. PARJ(155) : (R) current fraction of four-jet rate attributable to qqbarqqbar events rather than qqbargg ones; modified by mass effects. Not used in parton showers. PARJ(156) : (R) has two functions when using second order QCD. For a three-jet event, it gives the ratio of the second-order to the total three-jet cross-section in the given kinematical point. For a four-jet event, it gives the ratio of the modified four-jet cross-section, obtained when interference terms with not well defined colour flow are neglected, to the full unmodified one, all evaluated in the given kinematical point. Not used in parton showers. PARJ(157) - PARJ(159) : (I) used for cross-section calculations to include mass threshold effects to radiative photon cross-section. What is stored is basic cross-section, number of events generated and number that passed cuts. PARJ(160) : (R) nominal fraction of events that should contain a radiative photon. PARJ(161) - PARJ(164) : (I) give shape of radiative photon spectrum including weak effects. PARJ(168) : (R) actual y' value of current event in optimized perturbation theory for R; see MSTJ(111) and PARJ(128). PARJ(169) : (R) actual y' value of current event in optimized perturbation theory for three- and four-jet rate; see MSTJ(111) and PARJ(129). PARJ(171) : (R) fraction of cross-section corresponding to the axial coupling of quark pair to the intermediate gamma*/Z0 state; needed for the Abelian gluon model three-jet matrix element. ______________________________________________________________________ 3.4. Examples An ordinary e+e- annihilation event in the continuum, at a CM energy of 40 GeV, may be generated with CALL LUEEVT(0,40.) In this case a qqbar event is generated, including weak effects, followed by parton shower evolution and fragmentation/decay treatment. Before a call to LUEEVT, however, a number of default values may be changed, e.g. MSTJ(101) = 2 to use second order QCD matrix elements, giving a mixture of qqbar, qqbarg, qqbargg and qqbarq'qbar' events, MSTJ(102) = 1 to have QED only, MSTJ(104) = 6 to allow ttbar production as well, MSTJ(107) = 1 to include initial state photon radiation (including a treatment of the Z0 pole), PARJ(123) = 92.0 to change the Z0 mass, PARJ(81) = 0.3 to change the parton shower Lambda value, or PARJ(82) = 1.5 to change the parton shower cutoff. If initial state photon radiation is used, some restrictions apply on how one can alternate the generation of events at different energies or with different Z0 mass etc. These restrictions are not there for efficiency reasons (the extra time for recalculating the extra constants every time is small), but because it ties in with the the cross-section calculations (see PARJ(144)). Most parameters can be changed independently of each other. However, if just one or a few parameters/switches are changed, one should not be surprised to find a rather bad agreement with the data, like e.g. a too low or high average hadron multiplicity. It is therefore usually necessary to retune one parameter related to the perturbative QCD description, like alpha_S or Lambda, one of the two parameters a and b of the Lund symmetric fragmentation function (since they are so strongly correlated, it is often not necessary to retune both of them), and the average fragmentation transverse momentum - see Note 2 of the MSTJ(101) description for one example. For very detailed studies it may be necessary to retune even more parameters. The three-gluon or gluon-gluon-photon decay Upsilon may be simulated by a call CALL LUONIA(5,9.46) Unfortunately, with present top mass limits, this routine will not be of much interest for toponium studies (weak decays will dominate). A typical program for analysis of e+e- annihilation events at 100 GeV might look something like COMMON/LUJETS/N,K(4000,5),P(4000,5),V(4000,5) COMMON/LUDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200) COMMON/LUDAT2/KCHG(500,3),PMAS(500,4),PARF(2000),VCKM(4,4) COMMON/LUDAT3/MDCY(500,3),MDME(2000,2),BRAT(2000),KFDP(2000,5) MDCY(LUCOMP(111),1)=0 ! put pi0 stable MSTJ(104)=6 ! allow top-antitop production PMAS(6,1)=45. ! change top quark mass MSTJ(107)=1 ! include initial state radiation PARU(41)=1. ! use linear sphericity ..... ! other desired changes ..... ! initialize analysis statistics DO 100 IEVENT=1,1000 ! loop over events CALL LUEEVT(0,100.) ! generate new event IF(IEVENT.EQ.1) CALL LULIST(2) ! list first event CALL LUTABU(11) ! save particle composition ! statistics CALL LUEDIT(2) ! remove decayed particles CALL LUSPHE(SPH,APL) ! linear sphericity analysis IF(SPH.LT.0.) GOTO 100 ! too few particles in event for ! LUSPHE to work on it (unusual) CALL LUEDIT(31) ! orient event along axes above IF(IEVENT.EQ.1) CALL LULIST(2) ! list first treated event ..... ! fill analysis statistics CALL LUTHRU(THR,OBL) ! now do thrust analysis ..... ! more analysis statistics 100 CONTINUE ! CALL LUTABU(12) ! print particle composition ! statistics ..... ! print analysis statistics END ********************************************************************** Acknowledgements Constructive criticism from a number of users has been very helpful. In the work with JETSET 7.1 and 7.2, comments from Hans-Uno Bengtsson and Bill Gary have been of particular importance. In some cases, users have submitted pieces of code, with suggested extensions. This correspondence has been valuable, but the new code actually used in the program has been written entirely by the author, to ensure a uniform programming style. Also, the (moral) responsibility for any errors and other shortcomings rests solely with the author. ********************************************************************** References Ali80 A. Ali, J. G. Korner, G. Kramer, J. Willrodt, Nucl. Phys. B168 (1980) 409 A. Ali, E. Pietarinen, G. Kramer, J. Willrodt, Phys. Lett. B93 (1980) 155 And82 B. Andersson, G. Gustafson, T. Sjostrand, Nucl. Phys. B197 (1982) 45 And85 B. Andersson, G. Gustafson, T. Sjostrand, Physica Scripta 32 (1985) 574 And83 B. Andersson, G. Gustafson, G. Ingelman, T. Sjostrand, Phys. Rep. 97 (1983) 31 Bau90 U. Baur, M. Spira, P. M. Zerwas, Phys. Rev. D42 (1990) 815 Ben87 H.-U. Bengtsson, T. Sjostrand, Computer Phys. Comm. 46 (1987) 43 Bet86 JADE Collaboration, W. Bartel et al., Z. Physik C33 (1986) 23 S. Bethke, Habilitation thesis, LBL 50-208 (1987) Bet89 S. Bethke, Z. Physik C43 (1989) 331 Bia86 A. Bialas, R. Peschanski, Nucl. Phys. B273 (1986) 703 Cla79 L. Clavelli, Phys. Lett. B85 (1979) 111 A. V. Smilga, Nucl. Phys. B161 (1979) 449 Coc90 D. Cocolicchio, F. Feruglio, G.L. Fogli, J. Terron, CERN-TH.5909/90 F. Feruglio, private communication Fie78 R. D. Field, R. P. Feynman, Nucl. Phys. B136 (1978) 1 Fox79 G. C. Fox, S. Wolfram, Nucl. Phys. B149 (1979) 413 Gro81 T. R. Grose, K. O. Mikaelian, Phys. Rev. D23 (1981) 123 Hoy79 P. Hoyer, P. Osland, H. G. Sander, T. F. Walsh, P. M. Zerwas, Nucl. Phys. B161 (1979) 349 Ing80 G. Ingelman, T. Sjostrand, LUTP 80-12 (1980) G. Ingelman, LEPTO version 4.3, CERN program pool long writeup, program W5046 (1986) Ing87 G. Ingelman, Computer Phys. Comm. 46 (1987) 217 Ing87a G. Ingelman, A. Weigend, Computer Phys. Comm. 46 (1987) 241 Ing88 G. Ingelman, G. A. Schuler, Z. Phys. C40 (1988) 299 Jam90 F. James, Computer Physics Commun. 60 (1990) 329 Hew88 J.L. Hewett, S. Pakvasa, Phys. Rev. D37 (1988) 3165 Kra86 G. Kramer, B. Lampe, Fortschr. Physik 37 (1989) 161 Lor88 B. Lorstad, Int. J. of Mod. Phys. A4 (1989) 2861 Mar90 G. Marsaglia, A. Zaman, W.-W. Tsang, Stat. Prob. Lett. 9 (1990) 35 Nil87 B. Nilsson-Almqvist, E. Stenlund, Computer Phys. Comm. 43 (1987) 387 PDG86 Particle Data Group, M. Aguilar-Benitez et al., Phys. Lett. B170 (1986) 1 PDG88 Particle Data Group, G. P. Yost et al., Phys. Lett. B204 (1988) 1 Pet83 C. Peterson, D. Schlatter, I. Schmitt, P. Zerwas, Phys. Rev. D27 (1983) 105 Pet88 U. Pettersson, LU TP 88-5 (1988) L. Lonnblad, U. Pettersson, LU TP 88-15 (1988) L. Lonnblad, LU TP 89-10 (1989) QEG89 QCD event generators subgroup of the 1989 LEP physics workshop, B. Bambah et al., CERN-TH.5466/89 (1989), to appear in the proceedings of the workshop Sjo78 T. Sjostrand, B. Soderberg, LU TP 78-18 (1978) Sjo79 T. Sjostrand, LU TP 79-8 (1979) Sjo80 T. Sjostrand, LU TP 80-3 (1980) Sjo82 T. Sjostrand, Computer Phys. Comm. 27 (1982) 243 Sjo83 T. Sjostrand, Computer Phys. Comm. 28 (1983) 227 Sjo84 T. Sjostrand, Nucl. Phys. B248 (1984) 469 Sjo84a T. Sjostrand, Z. Physik C26 (1984) 93 Sjo86 T. Sjostrand, Computer Phys. Comm. 39 (1986) 347 Sjo87 T. Sjostrand, M. Bengtsson, Computer Phys. Comm. 43 (1987) 367 Sjo88 T. Sjostrand, Int. J. of Mod. Phys. A3 (1988) 751 **********************************************************************