1 <chapter name="Program Flow">
5 Recall that, to first order, the event generation process can be
6 subdivided into three stages:
8 <li>Initializaion.</li>
9 <li>The event loop.</li>
12 This is reflected in how the top-level <code>Pythia</code> class should
13 be used in the user-supplied main program, further outlined in the
14 following. Since the nature of the run is defined at the initialization
15 stage, this is where most of the PYTHIA user code has to be written.
16 So as not to confuse the reader unduly, the description of initialization
17 options has been subdivided into what would normally be used and what is
18 intended for more special applications.
21 At the bottom of this webpge is a complete survey of all public
22 <code>Pythia</code> methods and data members, in a more formal style
23 than the task-oriented descriptions found in the preceding sections.
24 This offers complementary information.
26 <h3>Initialization - normal usage</h3>
31 Already at the top of the main program file, you need to include the proper
36 To simplify typing, it also makes sense to declare
38 using namespace Pythia8;
44 The first step is to create a generator object,
49 It is this object that we will use from now on. Normally a run
50 will only contain one <code>Pythia</code> object. (But you can
51 use several <code>Pythia</code> objects, which then will be
52 independent of each other.)<br/>
53 By default all output from <code>Pythia</code> will be on the
54 <code>cout</code> stream, but the <code>list</code> methods below do
55 allow output to alternative streams or files.
60 You next want to set up the character of the run.
61 The pages under the "Setup Run Tasks" heading in the index
62 describe all the options available (with some very few exceptions,
63 found on the other pages).
64 The default values and your modifications are stored in two databases,
65 one for <aloc href="SettingsScheme">generic settings</aloc>
66 and one for <aloc href="ParticleDataScheme">particle data</aloc>.
67 Both of these are initialized with their default values by the
68 <code>Pythia</code> constructor. The default values can then be
69 changed, primarily by one of the two ways below, or by a combination
75 pythia.readString(string);
77 method repeatedly to do a change of a property at a time.
78 The information in the string is case-insensitive, but upper- and
79 lowercase can be combined for clarity. The rules are that<br/>
80 (i) if the first nonblank character of the string is a letter
81 it is assumed to contain a setting, and is sent on to
82 <code>pythia.settings.readString(string)</code>;<br/>
83 (ii) if instead the string begins with a digit it is assumed to
84 contain particle data updates, and so sent on to
85 <code>pythia.particleData.readString(string)</code>;<br/>
86 (iii) if none of the above, the string is assumed to be a comment,
87 i.e. nothing will be done.<br/>
88 In the former two cases, a warning is issued whenever a string
89 cannot be recognized (maybe because of a spelling mistake).<br/>
90 Some examples would be
92 pythia.readString("TimeShower:pTmin = 1.0");
93 pythia.readString("111:mayDecay = false");
95 The <code>readString(string)</code> method is intended primarily for
96 a few changes. It can also be useful if you want to construct a
97 parser for input files that contain commands both to PYTHIA and to
101 b) You can read in a file containing a list of those variables
102 you want to see changed, with a
104 pythia.readFile(fileName);
106 Each line in this file with be processes by the
107 <code>readString(string)</code> method introduced above. You can thus
108 freely mix comment lines and lines handed on to <code>Settings</code>
109 or to <code>ParticleData</code>.<br/>
110 This approach is better suited for more extensive changes than a direct
111 usage of <code>readString(string)</code>, and can also avoid having to
112 recompile and relink your main program between runs.<br/>
113 It is also possible to read input from an <code>istream</code>, by
114 default <code>cin</code>, rather than from a file. This may be convenient
115 if information is generated on-the-fly, within the same run.
118 Changes are made sequentially in the order the commands are encountered
119 during execution, meaning that if a parameter is changed several times
120 it is the last one that counts. The two special
121 <code><aloc href="Tunes">Tune:ee</aloc></code> and
122 <code><aloc href="Tunes">Tune:pp</aloc></code>
123 modes are expanded to change several settings in one go, but these obey
124 the same ordering rules.
130 Next comes the initialization stage, where all
131 remaining details of the generation are to be specified. The
132 <code>init(...)</code> method allows a few different input formats,
133 so you can pick the one convenient for you:
136 a) <code>pythia.init( idA, idB, eCM);</code><br/>
137 lets you specify the identities and the CM energy of the two incoming
138 beam particles, with A (B) assumed moving in the <ei>+z (-z)</ei>
142 b) <code>pythia.init( idA, idB, eA, eB);</code><br/>
143 is similar, but the two beam energies can be different, so the
144 collisions do not occur in the CM frame. If one of the beam energies
145 is below the particle mass you obtain a fixed-target topology.
148 c) <code>pythia.init( idA, idB, pxA, pyA, pzA, pxB, pyB, pzB);</code><br/>
149 is similar, but here you provide the three-momenta
150 <ei>(p_x, p_y, p_z)</ei> of the two incoming particles,
151 to allow for arbitrary beam directions.
154 d) <code>pythia.init(fileName);</code> <br/>
155 assumes a file in the <aloc href="LesHouchesAccord">Les Houches
156 Event File</aloc> format is provided.
159 e) <code>pythia.init();</code><br/>
160 with no arguments will read the beam parameters from the
161 <code><aloc href="MainProgramSettings">Main</aloc></code>
162 group of variables, which provides you with the same possibilities as
163 the above options a, b, c and d. If you don't change any of those you will
164 default to proton-proton collisions at 14 TeV, i.e. the nominal LHC
168 f) <code>pythia.init( LHAup*);</code> <br/>
169 assumes <aloc href="LesHouchesAccord">Les Houches Accord</aloc>
170 initialization and event information is available in an <code>LHAup</code>
171 class object, and that a pointer to this object is handed in.
175 If you want to have a list of the generator and particle data used,
176 either only what has been changed or everything, you can use
178 pythia.settings.listChanged();
179 pythia.settings.listAll();
180 pythia.particleData.listChanged();
181 pythia.particleData.listAll();
187 <h3>The event loop</h3>
192 Inside the event generation loop you generate the
193 next event using the <code>next()</code> method,
197 This method takes no arguments; everything has already been specified.
198 It does return a bool value, however, <code>false</code> when the
199 generation failed. This can be a "programmed death" when the
200 supply of input parton-level configurations on file is exhausted.
201 It can alternatively signal a failure of <code>Pythia</code> to
202 generate an event, or unphysical features in the event record at the
203 end of the generation step. It makes sense to allow a few <code>false</code>
204 values before a run is aborted, so long as the related faulty
210 The generated event is now stored in the <code>event</code>
211 object, of type <code><aloc href="EventRecord">Event</aloc></code>,
212 which is a public member of <code>pythia</code>. You therefore have
213 access to all the tools described on the pages under the "Study Output"
214 header in the index. For instance, an event can be listed with
215 <code>pythia.event.list()</code>, the identity of the <ei>i</ei>'th
216 <aloc href="ParticleProperties">particle</aloc> is given by
217 <code>pythia.event[i].id()</code>, and so on.<br/>
218 The hard process - roughly the information normally stored in the
219 Les Houches Accord event record - is available as a second object,
220 <code>process</code>, also of type <code>Event</code>.<br/>
221 A third useful public object is
222 <code><aloc href="EventInformation">info</aloc></code>, which offers
223 a set of one-of-a kind pieces of information about the most recent
233 <li>At the end of the generation process, you can call
237 to get some run statistics, on cross sections and the number of errors
238 and warnings encountered.
243 <h3>Advanced usage, mainly for initialization</h3>
245 A) Necessary data are automatically loaded when you use the
246 default PYTHIA installation directory structure and run the main
247 programs in the <code>examples</code> subdirectory. However, in the
248 general case, you must provide the path of the <code>xmldoc</code>
249 directory, where default settings and particle data are found.
250 This can be done in two ways.
255 You can set the environment variable <code>PYTHIA8DATA</code> to
256 contain the location of the <code>xmldoc</code> directory. In the
257 <code>csh</code> and <code>tcsh</code> shells this could e.g. be
259 setenv PYTHIA8DATA /home/myname/pythia81xx/xmldoc
261 while in other shells it could be
263 export PYTHIA8DATA=/home/myname/pythia81xx/xmldoc
265 where xx is the subversion number.<br/>
266 Recall that environment variables set locally are only defined in the
267 current instance of the shell. The above lines should go into your
268 <code>.cshrc</code> and <code>.bashrc</code> files, respectively,
269 if you want a more permanant assignment.
274 You can provide the path as argument to the <code>Pythia</code>
277 Pythia pythia("/home/myname/pythia81xx/xmldoc");
281 where again xx is the subversion number.<br/>
282 When <code>PYTHIA8DATA</code> is set it takes precedence, else
283 the path in the constructor is used, else one defaults to the
284 <code>../xmldoc</code> directory.
287 B) You can override the default behaviour of PYTHIA not only by the
288 settings and particle data, but also by replacing some of the
289 PYTHIA standard routines by ones of your own. Of course, this is only
290 possible if your routines fit into the general PYTHIA framework.
291 Therefore they must be coded according to the the rules relevant
292 in each case, as a derived class of a PYTHIA base class, and a pointer
293 to such an object must be handed in by one of the methods below.
294 These calls must be made before the <code>pythia.init(...)</code> call.
299 If you are not satisfied with the list of parton density functions that
300 are implemented internally or available via the LHAPDF interface
301 (see the <aloc href="PDFSelection">PDF Selection</aloc> page), you
302 can suppy your own by a call to the <code>setPDFPtr(...)</code> method
304 pythia.setPDFptr( pdfAPtr, pdfBPtr);
306 where <code>pdfAPtr</code> and <code>pdfBPtr</code> are pointers to
307 two <code>Pythia</code> <aloc href="PartonDistributions">PDF
308 objects</aloc>. Note that <code>pdfAPtr</code> and <code>pdfBPtr</code>
309 cannot point to the same object; even if the PDF set is the same,
310 two copies are needed to keep track of two separate sets of <ei>x</ei>
311 and density values.<br/>
312 If you further wish to use separate PDF's for the hard process of an
313 event than the ones being used for everything else, the extended form
315 pythia.setPDFptr( pdfAPtr, pdfBPtr, pdfHardAPtr, pdfHardBPtr);
317 allows you to specify those separately, and then the first two sets
318 would only be used for the showers and for multiple interactions.
323 If you want to perform some particle decays with an
324 external generator, you can call the <code>setDecayPtr(...)</code>
327 pythia.setDecayPtr( decayHandlePtr, particles);
329 where the <code>decayHandlePtr</code> derives from the
330 <code><aloc href="ExternalDecays">DecayHandler</aloc></code> base
331 class and <code>particles</code> is a vector of particle codes to be
337 If you want to use an external random number generator,
338 you can call the <code>setRndmEnginePtr(...)</code> method
340 pythia.setRndmEnginePtr( rndmEnginePtr);
342 where <code>rndmEnginePtr</code> derives from the
343 <code><aloc href="RandomNumbers">RndmEngine</aloc></code> base class.
344 The <code>Pythia</code> default random number generator is perfectly
345 good, so this is only intended for consistency in bigger frameworks.
350 If you want to interrupt the evolution at various stages,
351 to interrogate the event and possibly veto it, or you want to
352 reweight the cross section, you can use
354 pythia.setUserHooksPtr( userHooksPtr);
356 where <code>userHooksPtr</code> derives from the
357 <code><aloc href="UserHooks">UserHooks</aloc></code> base class.
362 If you want to use your own parametrization of beam momentum spread and
363 interaction vertex, rather than the provided simple Gaussian
364 parametrization (off by default), you can call
366 pythia.setBeamShapePtr( beamShapePtr);
368 where <code>beamShapePtr</code> derives from the
369 <code><aloc href="BeamShape">BeamShape</aloc></code> base class.
374 If you want to implement a cross section of your own, but still make use
375 of the built-in phase space selection machinery, you can use
377 pythia.setSigmaPtr( sigmaPtr);
379 where <code>sigmaPtr</code> of type <code>SigmaProcess*</code> is an
380 instance of a class derived from one of the <code>Sigma1Process</code>,
381 <code>Sigma2Process</code> and <code>Sigma3Process</code> base classes
382 in their turn derived from
383 <code><aloc href="SemiInternalProcesses">SigmaProcess</aloc></code>.
384 This call can be used repeatedly to hand in several different processes.
389 If your cross section contains the production of a new resonance
390 with known analytical expression for all the relevant partial widths,
391 you can make this resonance available to the program with
393 pythia.setResonancePtr( resonancePtr);
395 where <code>resonancePtr</code> of type <code>ResonanceWidths*</code>
396 is an instance of a class derived from the
397 <code><aloc href="SemiInternalResonances">ResonanceWidths</aloc></code>
398 base class. In addition you need to add the particle to the normal
399 <aloc href="ParticleDataScheme">particle and decay database</aloc>.
400 This procedure can be used repeatedly to hand in several different
406 If you are a real expert and want to <aloc href="ImplementNewShowers">replace
407 the PYTHIA initial- and final-state showers</aloc>, you can use
409 pythia.setShowerPtr( timesDecPtr, timesPtr, spacePtr);
411 where <code>timesDecPtr</code> and <code>timesPtr</code>
412 derive from the <code>TimeShower</code> base class, and
413 <code>spacePtr</code> from <code>SpaceShower</code>.
419 C) Some comments on collecting several tasks in the same run.
423 PYTHIA has not been written for threadsafe execution on multicore
424 processors. If you want to use all cores,
425 the most efficient way presumably is to start correspondingly many jobs,
426 with different random number seeds, and add the statistics at the end.
431 In some cases it is convenient to use more than one <code>Pythia</code>
432 object. The key example would be the simultaneous generation of signal
433 and pileup events, see <code>main19.cc</code>. The two objects are then
434 set up and initialized separately, and generate events completely
435 independently of each other. It is only afterwards that the event records
436 are combined into one single super-event per beam crossing.
441 When time is not an issue, it may be that you want to perform several
442 separate subruns sequentially inside a run, e.g. to combine results for
443 several kinematical regions or to compare results for some different
444 tunes of the underlying event. One way to go is to create (and destroy)
445 one <code>pythia</code> object for each subrun, in which case they are
446 completely separate. You can also use the same <code>pythia</code> object,
447 only doing a new <code>init(...)</code> call for each subrun. In that
448 case, the settings and particle databases remain as they were in the
449 previous subrun, only affected by the specific changes you introduced in
450 the meantime. You can put those changes in the main program, with
451 <code>pythia.readString(string)</code>, using your own logic to decide
452 which ones to execute in which subrun. A corresponding possibility
453 exists with <code>pythia.readFile(fileName, subrun)</code> (or an
454 <code>istream</code> instead of a <code>fileName</code>), which as second
455 argument can take a non-negative subrun number. Then only those
456 sections of the file before any <code>Main:subrun = ...</code> line
457 or with matching <code>subrun</code> number will be read. That is, the
458 file could have a structure like
460 ( lines always read, i.e. "default values" always (re)set )
462 ( lines only read with readFile(fileName, 1) )
464 ( lines only read with readFile(fileName, 2) )
466 Both of these possibilities are illustrated in <code>main08.cc</code>.
471 When working with Les Houches Event Files, it may well be that your
472 intended input event sample is spread over several files, that you all
473 want to turn into complete events in one and the same run. There is no
474 problem with looping over several subruns, where each new subrun
475 is initialized with a new file. However, in that case you will do
476 a complete re-initialization each time around. If you want to avoid
477 this, note that there is an optional second argument for LHEF
478 initialization: <code>pythia.init(fileName, skipInit)</code>.
479 Alternatively, the tag <code>Main:LHEFskipInit</code> can be put
480 in a file of commands to obtain the same effect.
481 Here <code>skipInit</code> defaults to <code>false</code>,
482 but if set <code>true</code> then the new file will be simulated
483 with the same initialization data as already set in a previous
484 <code>pythia.init(...)</code> call. The burden rests on you to ensure
485 that this is indeed correct, e.g. that the two event samples have not
486 been generated for different beam energies.
491 <h2>The Pythia Class</h2>
493 Here follows the complete survey of all public <code>Pythia</code>
494 methods and data members.
496 <h3>Constructor and destructor</h3>
498 <method name="Pythia::Pythia(string xmlDir = "../xmldoc")">
499 creates an instance of the <code>Pythia</code> event generators,
500 and sets initial default values, notably for all settings and
501 particle data. You may use several <code>Pythia</code> instances
502 in the same run; only when you want to access external static
503 libraries could this cause problems. (This includes in particular
504 Fortran libraries such as <aloc href="PDFSelection">LHAPDF</aloc>.)
505 <argument name="xmlDir" default="../xmldoc">allows you to choose
506 from which directory the default settings and particle data values
507 are read in. If the <code>PYTHIA8DATA</code> environment variable
508 has been set it takes precedence. Else this optional argument allows
509 you to choose another directory location than the default one. Note
510 that it is only the directory location you can change, its contents
511 must be the ones of the <code>xmldoc</code> directory in the
512 standard distribution.
516 <method name="Pythia::~Pythia">
517 the destructor deletes the objects created by the constructor.
521 <method name="bool Pythia::readString(string line, bool warn = true)">
522 reads in a single string, that is interpreted as an instruction to
523 modify the value of a <aloc href="SettingsScheme">setting</aloc> or
524 <aloc href="ParticleDataScheme">particle data</aloc>, as already described
526 <argument name="line">
527 the string to be interpreted as an instruction.
529 <argument name="warn" default="true">
530 write a warning message or not whenever the instruction does not make
531 sense, e.g. if the variable does not exist in the databases.
533 <note>Note:</note> the method returns false if it fails to
534 make sense out of the string.
537 <method name="bool Pythia::readFile(string fileName, bool warn = true,
538 int subrun = SUBRUNDEFAULT)">
540 <methodmore name="bool Pythia::readFile(string fileName,
541 int subrun = SUBRUNDEFAULT)">
543 <methodmore name="bool Pythia::readFile(istream& inStream = cin,
544 bool warn = true, int subrun = SUBRUNDEFAULT)">
546 <methodmore name="bool Pythia::readFile(istream& inStream = cin,
547 int subrun = SUBRUNDEFAULT)">
548 reads in a whole file, where each line is interpreted as an instruction
549 to modify the value of a <aloc href="SettingsScheme">setting</aloc> or
550 <aloc href="ParticleDataScheme">particle data</aloc>, cf. the above
551 <code>readString</code> method. All four forms of the
552 <code>readFile</code> command share code for actually reading a file.
553 <argument name="fileName">
554 the file from which instructions are read.
556 <argument name="inStream">
557 an istream from which instructions are read.
559 <argument name="warn" default="true">
560 write a warning message or not whenever the instruction does not make
561 sense, e.g. if the variable does not exist in the databases. In the
562 command forms where <code>warn</code> is omitted it is true.
564 <argument name="subrun">
565 allows you have several optional sets of commands within the same file.
566 Only those sections of the file before any <code>Main:subrun = ...</code>
567 line or following such a line with matching subrun number will be read.
568 The subrun number should not be negative; negative codes like
569 <code>SUBRUNDEFAULT</code> corresponds to no specific subrun.
571 <note>Note:</note> the method returns false if it fails to
572 make sense out of any one line.
575 <method name="bool Pythia::setPDFPtr( PDF* pdfAPtr, PDF* pdfBPtr,
576 PDF* pdfHardAPtr = 0, PDF* pdfHardBPtr = 0)">
577 offers the possibility to link in external PDF sets for usage inside
578 the program. The rules for constructing your own class from
579 the <code>PDF</code> base class are described
580 <aloc href="PartonDistributions">here</aloc>.
581 <argument name="pdfAPtr, pdfBPtr">
582 pointers to two <code>PDF</code>-derived objects, one for each of
583 the incoming beams. The two objects have to be instantiated by you
584 in your program. Even if the two beam particles are the same
585 (protons, say) two separate instances are required, since current
586 information is cached in the objects. If both arguments are zero
587 then any previous linkage to external PDF's is disconnected,
588 see further Note 2 below.
590 <argument name="pdfHardAPtr, pdfHardBPtr" default="0">
591 pointers to two further <code>PDF</code>-derived objects, one for each
592 of the incoming beams. Normally only the first two arguments above would
593 be used, and then the same PDF sets would be invoked everywhere. If you
594 provide these two further pointers then two different sets of PDF's are
595 used. This second set is then exclusively for the generation of the hard
596 process from the process matrix elements library. The first set above
597 is for everything else, notably parton showers and multiple interactions.
599 <note>Note 1:</note> The method returns false if the input is obviously
600 incorrect, e.g. if two (nonzero) pointers agree.
601 <note>Note 2:</note> If you want to combine several subruns you can
602 call <code>setPDFPtr</code> with new arguments before each
603 <code>Pythia::init(...)</code> call. To revert from external PDF's
604 to the normal internal PDF selection you must call
605 <code>setPDFPtr(0, 0)</code> before <code>Pythia::init(...)</code>.
608 <method name="bool Pythia::setDecayPtr( DecayHandler* decayHandlePtr,
609 vector<int> handledParticles)">
610 offers the possibility to link to an external program that can do some
611 of the particle decays, instad of using the internal decay machinery.
612 With particles we here mean the normal hadrons and leptons, not
613 top quarks, electroweak bosons or new particles in BSM scenarios.
614 The rules for constructing your own class from the
615 <code>DecayHandler</code> base class are described
616 <aloc href="ExternalDecays">here</aloc>. Note that you can only
617 provide one external object, but this object in its turn could
618 very well hand on different particles to separate decay libraries.
619 <argument name="decayHandlePtr">
620 pointer to a <code>DecayHandler</code>-derived object. This object
621 must be instantiated by you in your program.
623 <argument name="handledParticles"> vector with the PDG identity codes
624 of the particles that should be handled by the external decay package.
625 You should only give the particle (positive) codes; the respective
626 antiparticle is always included as well.
628 <note>Note:</note> The method currently always returns true.
631 <method name="bool Pythia::setRndmEnginePtr( RndmEngine* rndmEnginePtr)">
632 offers the possibility to link to an external random number generator.
633 The rules for constructing your own class from the
634 <code>RndmEngine</code> base class are described
635 <aloc href="RandomNumbers">here</aloc>.
636 <argument name="rndmEnginePtr">
637 pointer to a <code>RndmEngine</code>-derived object. This object
638 must be instantiated by you in your program.
640 <note>Note:</note> The method returns true if the pointer is different
644 <method name="bool Pythia::setUserHooksPtr( UserHooks* userHooksPtr)">
645 offers the possibility to interact with the generation process at
646 a few different specified points, e.g. to reject undesirable events
647 at an early stage to save computer time. The rules for constructing
648 your own class from the <code>UserHooks</code> base class are described
649 <aloc href="UserHooks">here</aloc>. You can only hand in one such
650 pointer, but this may be to a class that implements several of the
651 different allowed possibilities.
652 <argument name="userHooksPtr">
653 pointer to a <code>userHooks</code>-derived object. This object
654 must be instantiated by you in your program.
656 <note>Note:</note> The method currently always returns true.
659 <method name="bool Pythia::setBeamShapePtr( BeamShape* beamShapePtr)">
660 offers the possibility to provide your own shape of the momentum and
661 space-time spread of the incoming beams. The rules for constructing
662 your own class from the <code>BeamShape</code> base class are described
663 <aloc href="BeamShape">here</aloc>.
664 <argument name="BeamShapePtr">
665 pointer to a <code>BeamShape</code>-derived object. This object
666 must be instantiated by you in your program.
668 <note>Note:</note> The method currently always returns true.
671 <method name="bool Pythia::setSigmaPtr( SigmaProcess* sigmaPtr)">
672 offers the possibility to link your own implementation of a process
673 and its cross section, to make it a part of the normal process
674 generation machinery, without having to recompile the
675 <code>Pythia</code> library itself. The rules for constructing your
676 own class from the <code>SigmaProcess</code> base class are described
677 <aloc href="SemiInternalProcesses">here</aloc>. You may call this
678 routine repeatedly, to add as many new processes as you wish.
679 <argument name="sigmaPtr">
680 pointer to a <code>SigmaProcess</code>-derived object. This object
681 must be instantiated by you in your program.
683 <note>Note:</note> The method currently always returns true.
686 <method name="bool Pythia::setResonancePtr( ResonanceWidths* resonancePtr)">
687 offers the possibility to link your own implementation of the
688 calculation of partial resonance widths, to make it a part of the
689 normal process generation machinery, without having to recompile the
690 <code>Pythia</code> library itself. This allows the decay of new
691 resonances to be handled internally, when combined with new particle
692 data. Note that the decay of normal hadrons cannot be modelled here;
693 this is for New Physics resonances. The rules for constructing your
694 own class from the <code>ResonanceWidths</code> base class are described
695 <aloc href="SemiInternalResonances">here</aloc>. You may call this
696 routine repeatedly, to add as many new resonances as you wish.
697 <argument name="resonancePtr">
698 pointer to a <code>ResonanceWidths</code>-derived object. This object
699 must be instantiated by you in your program.
701 <note>Note:</note> The method currently always returns true.
704 <method name="bool Pythia::setShowerPtr( TimeShower* timesDecPtr,
705 TimeShower* timesPtr = 0, SpaceShower* spacePtr = 0)">
706 offers the possibility to link your own parton shower routines as
707 replacements for the default ones. This is much more complicated
708 since the showers are so central and are so interlinked with other
709 parts of the program. Therefore it is also possible to do the
710 replacement in stages, from the more independent to the more
711 intertwined. The rules for constructing your own classes from the
712 <code>TimeShower</code> and <code>SpaceShower</code>base classes
713 are described <aloc href="ImplementNewShowers">here</aloc>. These
714 objects must be instantiated by you in your program.
715 <argument name="timesDecPtr">
716 pointer to a <code>TimeShower</code>-derived object for doing
717 timelike shower evolution in resonance decays, e.g. of a
718 <ei>Z^0</ei>. This is decoupled from beam remnants and parton
719 distributions, and is therefore the simplest kind of shower
720 to write. If you provide a value 0 then the internal shower
721 routine will be used.
723 <argument name="timesPtr" default="0">
724 pointer to a <code>TimeShower</code>-derived object for doing
725 all other timelike shower evolution, which is normally interleaved
726 with multiple interactions and spacelike showers, introducing
727 both further physics and further technical issues. If you retain
728 the default value 0 then the internal shower routine will be used.
729 You are allowed to use the same pointer as above for the
730 <code>timesDecPtr</code> if the same shower can fulfill both tasks.
732 <argument name="spacePtr" default="0">
733 pointer to a <code>SpaceShower</code>-derived object for doing
734 all spacelike shower evolution, which is normally interleaved
735 with multiple interactions and timelike showers. If you retain
736 the default value 0 then the internal shower routine will be used.
738 <note>Note:</note> The method currently always returns true.
743 At the initialization stage all the information provided above is
744 processed, and the stage is set up for the subsequent generation
745 of events. Several alterative forms of the <code>init</code> method
746 are available for this stage; pick the one most convenient.
748 <method name="bool Pythia::init( int idA, int idB, double eCM)">
749 initialize for collisions in the center-of-mass frame, with the
750 beams moving in the <ei>+-z</ei> directions.
751 <argument name="idA, idB">
752 particle identity code for the two incoming beams.
754 <argument name="eCM">
755 the CM energy of the collisions.
757 <note>Note:</note> The method returns false if the
758 initialization fails. It is then not possible to generate any
762 <method name="bool Pythia::init( int idA, int idB, double eA, double eB)">
763 initialize for collisions with back-to-back beams,
764 moving in the <ei>+-z</ei> directions, but with different energies.
765 <argument name="idA, idB">
766 particle identity code for the two incoming beams.
768 <argument name="eA, eB">
769 the energies of the two beams. If an energy is set to be below
770 the mass of the respective beam particle that particle is taken to
771 be at rest. This offers a simple possibility to simulate
772 fixed-target collisions.
774 <note>Note:</note> The method returns false if the
775 initialization fails. It is then not possible to generate any
779 <method name="bool Pythia::init( int idA, int idB, double pxA,
780 double pyA, double pzA, double pxB, double pyB, double pzB)">
781 initialize for collisions with arbitrary beam directions.
782 <argument name="idA, idB">
783 particle identity code for the two incoming beams.
785 <argument name="pxA, pyA, pzA">
786 the three-momntum vector <ei>(p_x, p_y, p_z)</ei> of the first
789 <argument name="pxB, pyB, pzB">
790 the three-momntum vector <ei>(p_x, p_y, p_z)</ei> of the second
793 <note>Note:</note> The method returns false if the
794 initialization fails. It is then not possible to generate any
798 <method name="bool Pythia::init( string LesHouchesEventFile,
799 bool skipInit = false)">
800 initialize for hard-process collisions fed in from an external file
801 with events, written according to the
802 <aloc href="LesHouchesAccord">Les Houches Event File</aloc>
804 <argument name="LesHouchesEventFile">
805 the file name (including path, where required) where the
806 events are stored, including relevant information on beam
807 identities and energies.
809 <argument name="skipInit" default="false">
810 By default this method does a complete reinitialization of the
811 generation process. If you set this argument to true then
812 no reinitialization will occur, only the pointer to the event
813 file is updated. This may come in handy if the full event sample
814 is split across several files generated under the same conditions
815 (except random numbers, of course). You then do the first
816 initialization with the default, and all subsequent ones with
817 true. Note that things may go wrong if the files are not created
818 under the same conditions.
820 <note>Note:</note> The method returns false if the
821 initialization fails. It is then not possible to generate any
825 <method name="bool Pythia::init()">
826 initialize for collisions, in any of the four above possibilities.
827 In this option the beams are not specified by input arguments,
828 but instead by the settings in the
829 <aloc href="BeamParameters">Beam Parameters</aloc> section.
830 This allows the beams to be specified in the same file as other
831 run instructions. The default settings give pp collisions at 14 TeV.
832 <note>Note:</note> The method returns false if the
833 initialization fails. It is then not possible to generate any
837 <method name="bool Pythia::init( LHAup* lhaUpPtr)">
838 initialize for hard-process collisions fed in from an external
839 source of events, consistent with the Les Houches Accord standard.
840 The rules for constructing your own class from the <code>LHAup</code>
841 base class are described <aloc href="LesHouchesAccord">here</aloc>.
842 This class is also required to provide the beam parameters.
843 <argument name="lhaUpPtr">
844 pointer to a <code>LHAup</code>-derived object. This object
845 must be instantiated by you in your program.
847 <note>Note:</note> The method returns false if the
848 initialization fails. It is then not possible to generate any
852 <h3>Generate events</h3>
854 The <code>next()</code> method is the main one to generate events.
855 In this section we also put a few other specialized methods that
856 may be useful in some circumstances.
858 <method name="bool Pythia::next()">
859 generate the next event. No input parameters are required; all
860 instructions have already been set up in the initialization stage.
861 <note>Note:</note> The method returns false if the event generation
862 fails. The event record is then not consistent and should not be
863 studied. When reading in hard collisions from a Les Houches Event File
864 the problem may be that the end of the file has been reached. This
865 can be checked with the
866 <code><aloc href="EventInformation">Info::atEndOfFile()</aloc></code>
870 <method name="bool Pythia::forceHadronLevel()">
871 hadronize the existing event record, i.e. perform string fragmentation
872 and particle decays. There are two main applications. Firstly,
873 you can use the same parton-level content as a basis for repeated
874 hadronization attempts, in schemes intended to save computer time.
875 Secondly, you may have an external program that can simulate the full
876 partonic level of the event - hard process, parton showers, multiple
877 interactions, beam remnants, colour flow, and so on - but not
878 hadronization. Further details are found
879 <aloc href="HadronLevelStandalone">here</aloc>.
880 <note>Note:</note> The method returns false if the hadronization
881 fails. The event record is then not consistent and should not be
885 <method name="bool Pythia::moreDecays()">
886 perform decays of all particles in the event record that have not been
887 decayed but should have been done so. This can be used e.g. for
888 repeated decay attempts, in schemes intended to save computer time.
889 Further details are found <aloc href="HadronLevelStandalone">here</aloc>.
890 <note>Note:</note> The method returns false if the decays fail. The
891 event record is then not consistent and should not be studied.
894 <method name="void Pythia::LHAeventList(ostream& os = cout)">
895 list the Les Houches Accord information on the current event, see
896 <code><aloc href="LesHouchesAccord">LHAup::listEvent(...)</aloc></code>.
897 (Other listings are available via the class members below, so this
898 listing is a special case that would not fit elsewhere.)
899 <argument name="os" default = "cout">
900 output stream where the listing occurs.
904 <method name="bool Pythia::LHAeventSkip(int nSkip)">
905 skip ahead a number of events in the Les Houches generation
906 sequence, without doing anything further with them, see
907 <code><aloc href="LesHouchesAccord">LHAup::skipEvent(nSkip)</aloc></code>.
908 Mainly intended for debug purposes, e.g. when an event at a known
909 location in a Les Houches Event File is causing problems.
910 <argument name="nSkip">
911 number of events to skip.
913 <note>Note:</note> The method returns false if the operation fails,
914 specifically if the end of a LHEF has been reached, cf.
915 <code>next()</code> above.
920 There is no required finalization step; you can stop generating events
921 when and how you want. It is still recommended that you make it a
922 routine to call the following method at the end.
924 <method name="void Pythia::statistics(bool all = false, bool reset = false)">
925 list statistics on the event generation, specifically total and partial
926 cross sections and the number of different errors. For more details see
927 <aloc href="EventStatistics">here</aloc>.
928 <argument name="all" default="false">
929 if true also statistics on multiple interactions is shown, by default not.
931 <argument name="reset" default="false"> if true then all counters,
932 e.g on events generated and errors experienced, are reset to zero
933 whenever the routine is called. The default instead is that
934 all stored statistics information is unaffected by the call. Counters
935 are automatically reset in each new <code>Pythia::init(...)</code>
936 call, however, so the only time the <code>reset</code> option makes a
937 difference is if <code>statistics(...)</code> is called several times
942 <h3>Interrogate settings</h3>
944 Normally settings are used in the setup and initialization stages
945 to determine the character of a run, e.g. read from a file with the
946 above-described <code>Pythia::readFile(...)</code> method.
947 There is no strict need for a user to interact with the
948 <code>Settings</code> database in any other way. However, as an option,
949 some settings variables have been left free for the user to set in
950 such a file, and then use in the main program to directly affect the
951 performance of that program, see
952 <aloc href="MainProgramSettings">here</aloc>. A typical example would
953 be the number of events to generate. For such applications the
954 following shortcuts to some <code>Settings</code> methods may be
957 <method name="bool Pythia::flag(string key)">
958 read in a boolean variable from the <code>Settings</code> database.
959 <argument name="key">
960 the name of the variable to be read.
964 <method name="int Pythia::mode(string key)">
965 read in an integer variable from the <code>Settings</code> database.
966 <argument name="key">
967 the name of the variable to be read.
971 <method name="double Pythia::parm(string key)">
972 read in a double-precision variable from the <code>Settings</code>
974 <argument name="key">
975 the name of the variable to be read.
979 <method name="string Pythia::word(string key)">
980 read in a string variable from the <code>Settings</code> database.
981 <argument name="key">
982 the name of the variable to be read.
986 <h3>Data members</h3>
988 The <code>Pythia</code> class contains a few public data members,
989 several of which play a central role. We list them here, with
990 links to the places where they are further described.
992 <method name="Event Pythia::process">
993 the hard-process event record, see <aloc href="EventRecord">here</aloc>
997 <method name="Event Pythia::event">
998 the complete event record, see <aloc href="EventRecord">here</aloc>
1002 <method name="Info Pythia::info">
1003 further information on the event-generation process, see
1004 <aloc href="EventInformation">here</aloc> for further details.
1007 <method name="Settings Pythia::settings">
1008 the settings database, see <aloc href="SettingsScheme">here</aloc>
1009 for further details.
1012 <method name="ParticleData Pythia::particleData">
1013 the particle properties and decay tables database, see
1014 <aloc href="ParticleDataScheme">here</aloc> for further details.
1017 <method name="Rndm Pythia::rndm">
1018 the random number generator, see <aloc href="RandomNumberSeed">here</aloc>
1019 and <aloc href="RandomNumbers">here</aloc> for further details.
1022 <method name="CoupSM Pythia::coupSM">
1023 Standard Model couplings and mixing matrices, see
1024 <aloc href="StandardModelParameters">here</aloc> for further details.
1027 <method name="SusyLesHouches Pythia::slha">
1028 parameters and particle data in the context of supersymmetric models,
1029 see <aloc href="SUSYLesHouchesAccord">here</aloc> for further details.
1032 <method name="PartonSystems Pythia::partonSystems">
1033 a grouping of the partons in the event record by subsystem,
1034 see <aloc href="AdvancedUsage">here</aloc> for further details.
1039 <!-- Copyright (C) 2010 Torbjorn Sjostrand -->