1 | <html> |
---|
2 | <head> |
---|
3 | <title>Program Flow</title> |
---|
4 | <link rel="stylesheet" type="text/css" href="pythia.css"/> |
---|
5 | <link rel="shortcut icon" href="pythia32.gif"/> |
---|
6 | </head> |
---|
7 | <body> |
---|
8 | |
---|
9 | <h2>Program Flow</h2> |
---|
10 | |
---|
11 | Recall that, to first order, the event generation process can be |
---|
12 | subdivided into three stages: |
---|
13 | <ol> |
---|
14 | <li>Initializaion.</li> |
---|
15 | <li>The event loop.</li> |
---|
16 | <li>Finishing.</li> |
---|
17 | </ol> |
---|
18 | This is reflected in how the top-level <code>Pythia</code> class should |
---|
19 | be used in the user-supplied main program, further outlined in the |
---|
20 | following. Since the nature of the run is defined at the initialization |
---|
21 | stage, this is where most of the PYTHIA user code has to be written. |
---|
22 | So as not to confuse the reader unduly, the description of initialization |
---|
23 | options has been subdivided into what would normally be used and what is |
---|
24 | intended for more special applications. |
---|
25 | |
---|
26 | <p/> |
---|
27 | At the bottom of this webpge is a complete survey of all public |
---|
28 | <code>Pythia</code> methods and data members, in a more formal style |
---|
29 | than the task-oriented descriptions found in the preceding sections. |
---|
30 | This offers complementary information. |
---|
31 | |
---|
32 | <h3>Initialization - normal usage</h3> |
---|
33 | |
---|
34 | <ol> |
---|
35 | |
---|
36 | <li> |
---|
37 | Already at the top of the main program file, you need to include the proper |
---|
38 | header file |
---|
39 | <pre> |
---|
40 | #include "Pythia.h" |
---|
41 | </pre> |
---|
42 | To simplify typing, it also makes sense to declare |
---|
43 | <pre> |
---|
44 | using namespace Pythia8; |
---|
45 | </pre> |
---|
46 | </li> |
---|
47 | |
---|
48 | <p/> |
---|
49 | <li> |
---|
50 | The first step is to create a generator object, |
---|
51 | e.g. with |
---|
52 | <pre> |
---|
53 | Pythia pythia; |
---|
54 | </pre> |
---|
55 | It is this object that we will use from now on. Normally a run |
---|
56 | will only contain one <code>Pythia</code> object. (But you can |
---|
57 | use several <code>Pythia</code> objects, which then will be |
---|
58 | independent of each other.)<br/> |
---|
59 | By default all output from <code>Pythia</code> will be on the |
---|
60 | <code>cout</code> stream, but the <code>list</code> methods below do |
---|
61 | allow output to alternative streams or files. |
---|
62 | </li> |
---|
63 | |
---|
64 | <p/> |
---|
65 | <li> |
---|
66 | You next want to set up the character of the run. |
---|
67 | The pages under the "Setup Run Tasks" heading in the index |
---|
68 | describe all the options available (with some very few exceptions, |
---|
69 | found on the other pages). |
---|
70 | The default values and your modifications are stored in two databases, |
---|
71 | one for <a href="SettingsScheme.html" target="page">generic settings</a> |
---|
72 | and one for <a href="ParticleDataScheme.html" target="page">particle data</a>. |
---|
73 | Both of these are initialized with their default values by the |
---|
74 | <code>Pythia</code> constructor. The default values can then be |
---|
75 | changed, primarily by one of the two ways below, or by a combination |
---|
76 | of them. |
---|
77 | |
---|
78 | <p/> |
---|
79 | a) You can use the |
---|
80 | <pre> |
---|
81 | pythia.readString(string); |
---|
82 | </pre> |
---|
83 | method repeatedly to do a change of a property at a time. |
---|
84 | The information in the string is case-insensitive, but upper- and |
---|
85 | lowercase can be combined for clarity. The rules are that<br/> |
---|
86 | (i) if the first nonblank character of the string is a letter |
---|
87 | it is assumed to contain a setting, and is sent on to |
---|
88 | <code>pythia.settings.readString(string)</code>;<br/> |
---|
89 | (ii) if instead the string begins with a digit it is assumed to |
---|
90 | contain particle data updates, and so sent on to |
---|
91 | <code>pythia.particleData.readString(string)</code>;<br/> |
---|
92 | (iii) if none of the above, the string is assumed to be a comment, |
---|
93 | i.e. nothing will be done.<br/> |
---|
94 | In the former two cases, a warning is issued whenever a string |
---|
95 | cannot be recognized (maybe because of a spelling mistake).<br/> |
---|
96 | Some examples would be |
---|
97 | <pre> |
---|
98 | pythia.readString("TimeShower:pTmin = 1.0"); |
---|
99 | pythia.readString("111:mayDecay = false"); |
---|
100 | </pre> |
---|
101 | The <code>readString(string)</code> method is intended primarily for |
---|
102 | a few changes. It can also be useful if you want to construct a |
---|
103 | parser for input files that contain commands both to PYTHIA and to |
---|
104 | other libraries.<br/> |
---|
105 | |
---|
106 | <p/> |
---|
107 | b) You can read in a file containing a list of those variables |
---|
108 | you want to see changed, with a |
---|
109 | <pre> |
---|
110 | pythia.readFile(fileName); |
---|
111 | </pre> |
---|
112 | Each line in this file with be processes by the |
---|
113 | <code>readString(string)</code> method introduced above. You can thus |
---|
114 | freely mix comment lines and lines handed on to <code>Settings</code> |
---|
115 | or to <code>ParticleData</code>.<br/> |
---|
116 | This approach is better suited for more extensive changes than a direct |
---|
117 | usage of <code>readString(string)</code>, and can also avoid having to |
---|
118 | recompile and relink your main program between runs.<br/> |
---|
119 | It is also possible to read input from an <code>istream</code>, by |
---|
120 | default <code>cin</code>, rather than from a file. This may be convenient |
---|
121 | if information is generated on-the-fly, within the same run. |
---|
122 | |
---|
123 | <p/> |
---|
124 | Changes are made sequentially in the order the commands are encountered |
---|
125 | during execution, meaning that if a parameter is changed several times |
---|
126 | it is the last one that counts. The two special |
---|
127 | <code><a href="Tunes.html" target="page">Tune:ee</a></code> and |
---|
128 | <code><a href="Tunes.html" target="page">Tune:pp</a></code> |
---|
129 | modes are expanded to change several settings in one go, but these obey |
---|
130 | the same ordering rules. |
---|
131 | <br/> |
---|
132 | </li> |
---|
133 | |
---|
134 | <p/> |
---|
135 | <li> |
---|
136 | Next comes the initialization stage, where all |
---|
137 | remaining details of the generation are to be specified. |
---|
138 | There is one standard method to use for this |
---|
139 | |
---|
140 | <p/> |
---|
141 | <code>pythia.init();</code><br/> |
---|
142 | with no arguments will read all relevant information from the |
---|
143 | <code><a href="SettingsScheme.html" target="page">Settings</a></code> |
---|
144 | and <code><a href="ParticleDataScheme.html" target="page">ParticleData</a></code> |
---|
145 | databases. Specifically the setup of incoming beams and energies |
---|
146 | is governed by the the beam parameters from the |
---|
147 | <code><a href="BeamParameters.html" target="page">Beams</a></code> |
---|
148 | group of variables. If you don't change any of those you will |
---|
149 | default to proton-proton collisions at 14 TeV, i.e. the nominal LHC |
---|
150 | values. |
---|
151 | |
---|
152 | <p/> |
---|
153 | A few alternative forms are available, where the arguments of the |
---|
154 | <code>init(...)</code> call can be used to set the beam parameters. |
---|
155 | These alternatives are now deprecated, and will bew removed for |
---|
156 | PYTHIA 8.2. |
---|
157 | |
---|
158 | <p/> |
---|
159 | a) <code>pythia.init( idA, idB, eCM);</code><br/> |
---|
160 | lets you specify the identities and the CM energy of the two incoming |
---|
161 | beam particles, with A (B) assumed moving in the <i>+z (-z)</i> |
---|
162 | direction. |
---|
163 | |
---|
164 | <p/> |
---|
165 | b) <code>pythia.init( idA, idB, eA, eB);</code><br/> |
---|
166 | is similar, but the two beam energies can be different, so the |
---|
167 | collisions do not occur in the CM frame. If one of the beam energies |
---|
168 | is below the particle mass you obtain a fixed-target topology. |
---|
169 | |
---|
170 | <p/> |
---|
171 | c) <code>pythia.init( idA, idB, pxA, pyA, pzA, pxB, pyB, pzB);</code><br/> |
---|
172 | is similar, but here you provide the three-momenta |
---|
173 | <i>(p_x, p_y, p_z)</i> of the two incoming particles, |
---|
174 | to allow for arbitrary beam directions. |
---|
175 | |
---|
176 | <p/> |
---|
177 | d) <code>pythia.init(fileName);</code> <br/> |
---|
178 | assumes a file in the <a href="LesHouchesAccord.html" target="page">Les Houches |
---|
179 | Event File</a> format is provided. |
---|
180 | |
---|
181 | <p/> |
---|
182 | e) <code>pythia.init( LHAup*);</code> <br/> |
---|
183 | assumes <a href="LesHouchesAccord.html" target="page">Les Houches Accord</a> |
---|
184 | initialization and event information is available in an <code>LHAup</code> |
---|
185 | class object, and that a pointer to this object is handed in. |
---|
186 | |
---|
187 | <p/> |
---|
188 | <li> |
---|
189 | If you want to have a list of the generator and particle data used, |
---|
190 | either only what has been changed or everything, you can use |
---|
191 | <pre> |
---|
192 | pythia.settings.listChanged(); |
---|
193 | pythia.settings.listAll(); |
---|
194 | pythia.particleData.listChanged(); |
---|
195 | pythia.particleData.listAll(); |
---|
196 | </pre> |
---|
197 | </li> |
---|
198 | |
---|
199 | </ol> |
---|
200 | |
---|
201 | <h3>The event loop</h3> |
---|
202 | |
---|
203 | <ol> |
---|
204 | |
---|
205 | <li> |
---|
206 | Inside the event generation loop you generate the |
---|
207 | next event using the <code>next()</code> method, |
---|
208 | <pre> |
---|
209 | pythia.next(); |
---|
210 | </pre> |
---|
211 | This method takes no arguments; everything has already been specified. |
---|
212 | It does return a bool value, however, <code>false</code> when the |
---|
213 | generation failed. This can be a "programmed death" when the |
---|
214 | supply of input parton-level configurations on file is exhausted. |
---|
215 | It can alternatively signal a failure of <code>Pythia</code> to |
---|
216 | generate an event, or unphysical features in the event record at the |
---|
217 | end of the generation step. It makes sense to allow a few <code>false</code> |
---|
218 | values before a run is aborted, so long as the related faulty |
---|
219 | events are skipped. |
---|
220 | </li> |
---|
221 | |
---|
222 | <p/> |
---|
223 | <li> |
---|
224 | The generated event is now stored in the <code>event</code> |
---|
225 | object, of type <code><a href="EventRecord.html" target="page">Event</a></code>, |
---|
226 | which is a public member of <code>pythia</code>. You therefore have |
---|
227 | access to all the tools described on the pages under the "Study Output" |
---|
228 | header in the index. For instance, an event can be listed with |
---|
229 | <code>pythia.event.list()</code>, the identity of the <i>i</i>'th |
---|
230 | <a href="ParticleProperties.html" target="page">particle</a> is given by |
---|
231 | <code>pythia.event[i].id()</code>, and so on.<br/> |
---|
232 | The hard process - roughly the information normally stored in the |
---|
233 | Les Houches Accord event record - is available as a second object, |
---|
234 | <code>process</code>, also of type <code>Event</code>.<br/> |
---|
235 | A third useful public object is |
---|
236 | <code><a href="EventInformation.html" target="page">info</a></code>, which offers |
---|
237 | a set of one-of-a kind pieces of information about the most recent |
---|
238 | event. |
---|
239 | </li> |
---|
240 | |
---|
241 | </ol> |
---|
242 | |
---|
243 | <h3>Finishing</h3> |
---|
244 | |
---|
245 | <ol> |
---|
246 | |
---|
247 | <li>At the end of the generation process, you can call |
---|
248 | <pre> |
---|
249 | pythia.stat(); |
---|
250 | </pre> |
---|
251 | to get some run statistics, on cross sections and the number of errors |
---|
252 | and warnings encountered. The alternative |
---|
253 | <code>pythia.statistics(...);</code> is equivalent but deprecated. |
---|
254 | </li> |
---|
255 | |
---|
256 | </ol> |
---|
257 | |
---|
258 | <h3>Advanced usage, mainly for initialization</h3> |
---|
259 | |
---|
260 | A) Necessary data are automatically loaded when you use the |
---|
261 | default PYTHIA installation directory structure and run the main |
---|
262 | programs in the <code>examples</code> subdirectory. However, in the |
---|
263 | general case, you must provide the path of the <code>xmldoc</code> |
---|
264 | directory, where default settings and particle data are found. |
---|
265 | This can be done in two ways. |
---|
266 | |
---|
267 | <ol> |
---|
268 | |
---|
269 | <li> |
---|
270 | You can set the environment variable <code>PYTHIA8DATA</code> to |
---|
271 | contain the location of the <code>xmldoc</code> directory. In the |
---|
272 | <code>csh</code> and <code>tcsh</code> shells this could e.g. be |
---|
273 | <pre> |
---|
274 | setenv PYTHIA8DATA /home/myname/pythia81xx/xmldoc |
---|
275 | </pre> |
---|
276 | while in other shells it could be |
---|
277 | <pre> |
---|
278 | export PYTHIA8DATA=/home/myname/pythia81xx/xmldoc |
---|
279 | </pre> |
---|
280 | where xx is the subversion number.<br/> |
---|
281 | Recall that environment variables set locally are only defined in the |
---|
282 | current instance of the shell. The above lines should go into your |
---|
283 | <code>.cshrc</code> and <code>.bashrc</code> files, respectively, |
---|
284 | if you want a more permanant assignment. |
---|
285 | </li> |
---|
286 | |
---|
287 | <p/> |
---|
288 | <li> |
---|
289 | You can provide the path as argument to the <code>Pythia</code> |
---|
290 | constructor, e.g. |
---|
291 | <pre> |
---|
292 | Pythia pythia("/home/myname/pythia81xx/xmldoc"); |
---|
293 | </pre> |
---|
294 | </li> |
---|
295 | </ol> |
---|
296 | where again xx is the subversion number.<br/> |
---|
297 | When <code>PYTHIA8DATA</code> is set it takes precedence, else |
---|
298 | the path in the constructor is used, else one defaults to the |
---|
299 | <code>../xmldoc</code> directory. |
---|
300 | |
---|
301 | <p/> |
---|
302 | B) You can override the default behaviour of PYTHIA not only by the |
---|
303 | settings and particle data, but also by replacing some of the |
---|
304 | PYTHIA standard routines by ones of your own. Of course, this is only |
---|
305 | possible if your routines fit into the general PYTHIA framework. |
---|
306 | Therefore they must be coded according to the the rules relevant |
---|
307 | in each case, as a derived class of a PYTHIA base class, and a pointer |
---|
308 | to such an object must be handed in by one of the methods below. |
---|
309 | These calls must be made before the <code>pythia.init(...)</code> call. |
---|
310 | |
---|
311 | <ol> |
---|
312 | |
---|
313 | <li> |
---|
314 | If you are not satisfied with the list of parton density functions that |
---|
315 | are implemented internally or available via the LHAPDF interface |
---|
316 | (see the <a href="PDFSelection.html" target="page">PDF Selection</a> page), you |
---|
317 | can suppy your own by a call to the <code>setPDFPtr(...)</code> method |
---|
318 | <pre> |
---|
319 | pythia.setPDFptr( pdfAPtr, pdfBPtr); |
---|
320 | </pre> |
---|
321 | where <code>pdfAPtr</code> and <code>pdfBPtr</code> are pointers to |
---|
322 | two <code>Pythia</code> <a href="PartonDistributions.html" target="page">PDF |
---|
323 | objects</a>. Note that <code>pdfAPtr</code> and <code>pdfBPtr</code> |
---|
324 | cannot point to the same object; even if the PDF set is the same, |
---|
325 | two copies are needed to keep track of two separate sets of <i>x</i> |
---|
326 | and density values.<br/> |
---|
327 | If you further wish to use separate PDF's for the hard process of an |
---|
328 | event than the ones being used for everything else, the extended form |
---|
329 | <pre> |
---|
330 | pythia.setPDFptr( pdfAPtr, pdfBPtr, pdfHardAPtr, pdfHardBPtr); |
---|
331 | </pre> |
---|
332 | allows you to specify those separately, and then the first two sets |
---|
333 | would only be used for the showers and for multiparton interactions. |
---|
334 | </li> |
---|
335 | |
---|
336 | <p/> |
---|
337 | <li> |
---|
338 | If you want to link to an external generator that feeds in events |
---|
339 | in the LHA format, you can call the <code>setLHAupPtr(...)</code> |
---|
340 | method |
---|
341 | <pre> |
---|
342 | pythia.setLHAupPtr( lhaUpPtr); |
---|
343 | </pre> |
---|
344 | where the <code>lhaUpPtr</code> derives from the |
---|
345 | <a href="LesHouchesAccord.html" target="page">LHAup</a> base class. |
---|
346 | </li> |
---|
347 | |
---|
348 | <p/> |
---|
349 | <li> |
---|
350 | If you want to perform some particle decays with an |
---|
351 | external generator, you can call the <code>setDecayPtr(...)</code> |
---|
352 | method |
---|
353 | <pre> |
---|
354 | pythia.setDecayPtr( decayHandlePtr, particles); |
---|
355 | </pre> |
---|
356 | where the <code>decayHandlePtr</code> derives from the |
---|
357 | <code><a href="ExternalDecays.html" target="page">DecayHandler</a></code> base |
---|
358 | class and <code>particles</code> is a vector of particle codes to be |
---|
359 | handled. |
---|
360 | </li> |
---|
361 | |
---|
362 | <p/> |
---|
363 | <li> |
---|
364 | If you want to use an external random number generator, |
---|
365 | you can call the <code>setRndmEnginePtr(...)</code> method |
---|
366 | <pre> |
---|
367 | pythia.setRndmEnginePtr( rndmEnginePtr); |
---|
368 | </pre> |
---|
369 | where <code>rndmEnginePtr</code> derives from the |
---|
370 | <code><a href="RandomNumbers.html" target="page">RndmEngine</a></code> base class. |
---|
371 | The <code>Pythia</code> default random number generator is perfectly |
---|
372 | good, so this is only intended for consistency in bigger frameworks. |
---|
373 | </li> |
---|
374 | |
---|
375 | <p/> |
---|
376 | <li> |
---|
377 | If you want to interrupt the evolution at various stages, |
---|
378 | to interrogate the event and possibly veto it, or you want to |
---|
379 | reweight the cross section, you can use |
---|
380 | <pre> |
---|
381 | pythia.setUserHooksPtr( userHooksPtr); |
---|
382 | </pre> |
---|
383 | where <code>userHooksPtr</code> derives from the |
---|
384 | <code><a href="UserHooks.html" target="page">UserHooks</a></code> base class. |
---|
385 | </li> |
---|
386 | |
---|
387 | <p/> |
---|
388 | <li> |
---|
389 | If you want to use your own merging scale definition for |
---|
390 | matrix element + parton shower merging, you can call |
---|
391 | <pre> |
---|
392 | pythia.setMergingHooksPtr( mergingHooksPtr); |
---|
393 | </pre> |
---|
394 | where <code>mergingHooksPtr</code> derives from the |
---|
395 | <code><a href="MatrixElementMerging.html" target="page">MergingHooks</a></code> base class. |
---|
396 | </li> |
---|
397 | |
---|
398 | <p/> |
---|
399 | <li> |
---|
400 | If you want to use your own parametrization of beam momentum spread and |
---|
401 | interaction vertex, rather than the provided simple Gaussian |
---|
402 | parametrization (off by default), you can call |
---|
403 | <pre> |
---|
404 | pythia.setBeamShapePtr( beamShapePtr); |
---|
405 | </pre> |
---|
406 | where <code>beamShapePtr</code> derives from the |
---|
407 | <code><a href="BeamShape.html" target="page">BeamShape</a></code> base class. |
---|
408 | </li> |
---|
409 | |
---|
410 | <p/> |
---|
411 | <li> |
---|
412 | If you want to implement a cross section of your own, but still make use |
---|
413 | of the built-in phase space selection machinery, you can use |
---|
414 | <pre> |
---|
415 | pythia.setSigmaPtr( sigmaPtr); |
---|
416 | </pre> |
---|
417 | where <code>sigmaPtr</code> of type <code>SigmaProcess*</code> is an |
---|
418 | instance of a class derived from one of the <code>Sigma1Process</code>, |
---|
419 | <code>Sigma2Process</code> and <code>Sigma3Process</code> base classes |
---|
420 | in their turn derived from |
---|
421 | <code><a href="SemiInternalProcesses.html" target="page">SigmaProcess</a></code>. |
---|
422 | This call can be used repeatedly to hand in several different processes. |
---|
423 | </li> |
---|
424 | |
---|
425 | <p/> |
---|
426 | <li> |
---|
427 | If your cross section contains the production of a new resonance |
---|
428 | with known analytical expression for all the relevant partial widths, |
---|
429 | you can make this resonance available to the program with |
---|
430 | <pre> |
---|
431 | pythia.setResonancePtr( resonancePtr); |
---|
432 | </pre> |
---|
433 | where <code>resonancePtr</code> of type <code>ResonanceWidths*</code> |
---|
434 | is an instance of a class derived from the |
---|
435 | <code><a href="SemiInternalResonances.html" target="page">ResonanceWidths</a></code> |
---|
436 | base class. In addition you need to add the particle to the normal |
---|
437 | <a href="ParticleDataScheme.html" target="page">particle and decay database</a>. |
---|
438 | This procedure can be used repeatedly to hand in several different |
---|
439 | resonances. |
---|
440 | </li> |
---|
441 | |
---|
442 | <p/> |
---|
443 | <li> |
---|
444 | If you are a real expert and want to <a href="ImplementNewShowers.html" target="page">replace |
---|
445 | the PYTHIA initial- and final-state showers</a>, you can use |
---|
446 | <pre> |
---|
447 | pythia.setShowerPtr( timesDecPtr, timesPtr, spacePtr); |
---|
448 | </pre> |
---|
449 | where <code>timesDecPtr</code> and <code>timesPtr</code> |
---|
450 | derive from the <code>TimeShower</code> base class, and |
---|
451 | <code>spacePtr</code> from <code>SpaceShower</code>. |
---|
452 | </li> |
---|
453 | |
---|
454 | </ol> |
---|
455 | |
---|
456 | <p/> |
---|
457 | C) Some comments on collecting several tasks in the same run. |
---|
458 | <ol> |
---|
459 | |
---|
460 | <li> |
---|
461 | PYTHIA has not been written for threadsafe execution on multicore |
---|
462 | processors. If you want to use all cores, |
---|
463 | the most efficient way presumably is to start correspondingly many jobs, |
---|
464 | with different random number seeds, and add the statistics at the end. |
---|
465 | However, note that several instances can be set up in the same main |
---|
466 | program, since instances are completely independent of each other, |
---|
467 | so each instance could be run inside a separate thread. |
---|
468 | </li> |
---|
469 | |
---|
470 | <p/> |
---|
471 | <li> |
---|
472 | In some cases it is convenient to use more than one <code>Pythia</code> |
---|
473 | object. The key example would be the simultaneous generation of signal |
---|
474 | and pileup events, see <code>main19.cc</code>. The two objects are then |
---|
475 | set up and initialized separately, and generate events completely |
---|
476 | independently of each other. It is only afterwards that the event records |
---|
477 | are combined into one single super-event per beam crossing. |
---|
478 | </li> |
---|
479 | |
---|
480 | <p/> |
---|
481 | <li> |
---|
482 | When time is not an issue, it may be that you want to perform several |
---|
483 | separate subruns sequentially inside a run, e.g. to combine results for |
---|
484 | several kinematical regions or to compare results for some different |
---|
485 | tunes of the underlying event. One way to go is to create (and destroy) |
---|
486 | one <code>pythia</code> object for each subrun, in which case they are |
---|
487 | completely separate. You can also use the same <code>pythia</code> object, |
---|
488 | only doing a new <code>init(...)</code> call for each subrun. In that |
---|
489 | case, the settings and particle databases remain as they were in the |
---|
490 | previous subrun, only affected by the specific changes you introduced in |
---|
491 | the meantime. You can put those changes in the main program, with |
---|
492 | <code>pythia.readString(string)</code>, using your own logic to decide |
---|
493 | which ones to execute in which subrun. A corresponding possibility |
---|
494 | exists with <code>pythia.readFile(fileName, subrun)</code> (or an |
---|
495 | <code>istream</code> instead of a <code>fileName</code>), which as second |
---|
496 | argument can take a non-negative subrun number. Then only those |
---|
497 | sections of the file before any <code>Main:subrun = ...</code> line |
---|
498 | or with matching <code>subrun</code> number will be read. That is, the |
---|
499 | file could have a structure like |
---|
500 | <pre> |
---|
501 | ( lines always read, i.e. "default values" always (re)set ) |
---|
502 | Main:subrun = 1 |
---|
503 | ( lines only read with readFile(fileName, 1) ) |
---|
504 | Main:subrun = 2 |
---|
505 | ( lines only read with readFile(fileName, 2) ) |
---|
506 | </pre> |
---|
507 | Both of these possibilities are illustrated in <code>main08.cc</code>. |
---|
508 | </li> |
---|
509 | |
---|
510 | <p/> |
---|
511 | <li> |
---|
512 | When working with Les Houches Event Files, it may well be that your |
---|
513 | intended input event sample is spread over several files, that you all |
---|
514 | want to turn into complete events in one and the same run. There is no |
---|
515 | problem with looping over several subruns, where each new subrun |
---|
516 | is initialized with a new file, with name set in <code>Beams:LHEF</code>. |
---|
517 | However, in that case you will do a complete re-initialization each time |
---|
518 | around. If you want to avoid this, note that the flag |
---|
519 | <code>Beams:newLHEFsameInit = true</code> can be set for the second and |
---|
520 | subsequent subruns. Then the new file will be simulated with the same |
---|
521 | initialization data as already set in a previous |
---|
522 | <code>pythia.init()</code> call. The burden rests on you to ensure |
---|
523 | that this is indeed correct, e.g. that the two event samples have not |
---|
524 | been generated for different beam energies. Also note that cross |
---|
525 | sections for processes will be based on the information in the |
---|
526 | first-read file, when the full initialization is performed. |
---|
527 | </li> |
---|
528 | |
---|
529 | </ol> |
---|
530 | |
---|
531 | <h2>The Pythia Class</h2> |
---|
532 | |
---|
533 | Here follows the complete survey of all public <code>Pythia</code> |
---|
534 | methods and data members. |
---|
535 | |
---|
536 | <h3>Constructor and destructor</h3> |
---|
537 | |
---|
538 | <a name="method1"></a> |
---|
539 | <p/><strong>Pythia::Pythia(string xmlDir = "../xmldoc") </strong> <br/> |
---|
540 | creates an instance of the <code>Pythia</code> event generators, |
---|
541 | and sets initial default values, notably for all settings and |
---|
542 | particle data. You may use several <code>Pythia</code> instances |
---|
543 | in the same run; only when you want to access external static |
---|
544 | libraries could this cause problems. (This includes in particular |
---|
545 | Fortran libraries such as <a href="PDFSelection.html" target="page">LHAPDF</a>.) |
---|
546 | <br/><code>argument</code><strong> xmlDir </strong> (<code>default = <strong>../xmldoc</strong></code>) : allows you to choose |
---|
547 | from which directory the default settings and particle data values |
---|
548 | are read in. If the <code>PYTHIA8DATA</code> environment variable |
---|
549 | has been set it takes precedence. Else this optional argument allows |
---|
550 | you to choose another directory location than the default one. Note |
---|
551 | that it is only the directory location you can change, its contents |
---|
552 | must be the ones of the <code>xmldoc</code> directory in the |
---|
553 | standard distribution. |
---|
554 | |
---|
555 | |
---|
556 | |
---|
557 | <a name="method2"></a> |
---|
558 | <p/><strong>Pythia::~Pythia </strong> <br/> |
---|
559 | the destructor deletes the objects created by the constructor. |
---|
560 | |
---|
561 | <h3>Set up run</h3> |
---|
562 | |
---|
563 | <a name="method3"></a> |
---|
564 | <p/><strong>bool Pythia::readString(string line, bool warn = true) </strong> <br/> |
---|
565 | reads in a single string, that is interpreted as an instruction to |
---|
566 | modify the value of a <a href="SettingsScheme.html" target="page">setting</a> or |
---|
567 | <a href="ParticleDataScheme.html" target="page">particle data</a>, as already described |
---|
568 | above. |
---|
569 | <br/><code>argument</code><strong> line </strong> : |
---|
570 | the string to be interpreted as an instruction. |
---|
571 | |
---|
572 | <br/><code>argument</code><strong> warn </strong> (<code>default = <strong>true</strong></code>) : |
---|
573 | write a warning message or not whenever the instruction does not make |
---|
574 | sense, e.g. if the variable does not exist in the databases. |
---|
575 | |
---|
576 | <br/><b>Note:</b> the method returns false if it fails to |
---|
577 | make sense out of the string. |
---|
578 | |
---|
579 | |
---|
580 | <a name="method4"></a> |
---|
581 | <p/><strong>bool Pythia::readFile(string fileName, bool warn = true, int subrun = SUBRUNDEFAULT) </strong> <br/> |
---|
582 | |
---|
583 | <strong>bool Pythia::readFile(string fileName, int subrun = SUBRUNDEFAULT) </strong> <br/> |
---|
584 | |
---|
585 | <strong>bool Pythia::readFile(istream& inStream = cin, bool warn = true, int subrun = SUBRUNDEFAULT) </strong> <br/> |
---|
586 | |
---|
587 | <strong>bool Pythia::readFile(istream& inStream = cin, int subrun = SUBRUNDEFAULT) </strong> <br/> |
---|
588 | reads in a whole file, where each line is interpreted as an instruction |
---|
589 | to modify the value of a <a href="SettingsScheme.html" target="page">setting</a> or |
---|
590 | <a href="ParticleDataScheme.html" target="page">particle data</a>, cf. the above |
---|
591 | <code>readString</code> method. All four forms of the |
---|
592 | <code>readFile</code> command share code for actually reading a file. |
---|
593 | <br/><code>argument</code><strong> fileName </strong> : |
---|
594 | the file from which instructions are read. |
---|
595 | |
---|
596 | <br/><code>argument</code><strong> inStream </strong> : |
---|
597 | an istream from which instructions are read. |
---|
598 | |
---|
599 | <br/><code>argument</code><strong> warn </strong> (<code>default = <strong>true</strong></code>) : |
---|
600 | write a warning message or not whenever the instruction does not make |
---|
601 | sense, e.g. if the variable does not exist in the databases. In the |
---|
602 | command forms where <code>warn</code> is omitted it is true. |
---|
603 | |
---|
604 | <br/><code>argument</code><strong> subrun </strong> : |
---|
605 | allows you have several optional sets of commands within the same file. |
---|
606 | Only those sections of the file before any <code>Main:subrun = ...</code> |
---|
607 | line or following such a line with matching subrun number will be read. |
---|
608 | The subrun number should not be negative; negative codes like |
---|
609 | <code>SUBRUNDEFAULT</code> corresponds to no specific subrun. |
---|
610 | |
---|
611 | <br/><b>Note:</b> the method returns false if it fails to |
---|
612 | make sense out of any one line. |
---|
613 | |
---|
614 | |
---|
615 | <a name="method5"></a> |
---|
616 | <p/><strong>bool Pythia::setPDFPtr( PDF* pdfAPtr, PDF* pdfBPtr, PDF* pdfHardAPtr = 0, PDF* pdfHardBPtr = 0) </strong> <br/> |
---|
617 | offers the possibility to link in external PDF sets for usage inside |
---|
618 | the program. The rules for constructing your own class from |
---|
619 | the <code>PDF</code> base class are described |
---|
620 | <a href="PartonDistributions.html" target="page">here</a>. |
---|
621 | <br/><code>argument</code><strong> pdfAPtr, pdfBPtr </strong> : |
---|
622 | pointers to two <code>PDF</code>-derived objects, one for each of |
---|
623 | the incoming beams. The two objects have to be instantiated by you |
---|
624 | in your program. Even if the two beam particles are the same |
---|
625 | (protons, say) two separate instances are required, since current |
---|
626 | information is cached in the objects. If both arguments are zero |
---|
627 | then any previous linkage to external PDF's is disconnected, |
---|
628 | see further Note 2 below. |
---|
629 | |
---|
630 | <br/><code>argument</code><strong> pdfHardAPtr, pdfHardBPtr </strong> (<code>default = <strong>0</strong></code>) : |
---|
631 | pointers to two further <code>PDF</code>-derived objects, one for each |
---|
632 | of the incoming beams. Normally only the first two arguments above would |
---|
633 | be used, and then the same PDF sets would be invoked everywhere. If you |
---|
634 | provide these two further pointers then two different sets of PDF's are |
---|
635 | used. This second set is then exclusively for the generation of the hard |
---|
636 | process from the process matrix elements library. The first set above |
---|
637 | is for everything else, notably parton showers and multiparton interactions. |
---|
638 | |
---|
639 | <br/><b>Note 1:</b> The method returns false if the input is obviously |
---|
640 | incorrect, e.g. if two (nonzero) pointers agree. |
---|
641 | <br/><b>Note 2:</b> If you want to combine several subruns you can |
---|
642 | call <code>setPDFPtr</code> with new arguments before each |
---|
643 | <code>Pythia::init(...)</code> call. To revert from external PDF's |
---|
644 | to the normal internal PDF selection you must call |
---|
645 | <code>setPDFPtr(0, 0)</code> before <code>Pythia::init(...)</code>. |
---|
646 | |
---|
647 | |
---|
648 | <a name="method6"></a> |
---|
649 | <p/><strong>bool Pythia::setLHAupPtr( LHAup* lhaUpPtrIn) </strong> <br/> |
---|
650 | offers linkage to an external generator that feeds in events |
---|
651 | in the LHA format, see |
---|
652 | <a href="LesHouchesAccord.html" target="page">Les Houches Accord</a>, |
---|
653 | assuming that |
---|
654 | <code><a href="BeamParameters.html" target="page">Beams:frameType = 5</a></code> |
---|
655 | has been set. |
---|
656 | <br/><code>argument</code><strong> lhaUpPtrIn </strong> : |
---|
657 | pointer to a <code>LHAup</code>-derived object. |
---|
658 | |
---|
659 | <br/><b>Note:</b> The method currently always returns true. |
---|
660 | |
---|
661 | |
---|
662 | <a name="method7"></a> |
---|
663 | <p/><strong>bool Pythia::setDecayPtr( DecayHandler* decayHandlePtr, vector<int> handledParticles) </strong> <br/> |
---|
664 | offers the possibility to link to an external program that can do some |
---|
665 | of the particle decays, instad of using the internal decay machinery. |
---|
666 | With particles we here mean the normal hadrons and leptons, not |
---|
667 | top quarks, electroweak bosons or new particles in BSM scenarios. |
---|
668 | The rules for constructing your own class from the |
---|
669 | <code>DecayHandler</code> base class are described |
---|
670 | <a href="ExternalDecays.html" target="page">here</a>. Note that you can only |
---|
671 | provide one external object, but this object in its turn could |
---|
672 | very well hand on different particles to separate decay libraries. |
---|
673 | <br/><code>argument</code><strong> decayHandlePtr </strong> : |
---|
674 | pointer to a <code>DecayHandler</code>-derived object. This object |
---|
675 | must be instantiated by you in your program. |
---|
676 | |
---|
677 | <br/><code>argument</code><strong> handledParticles </strong> : vector with the PDG identity codes |
---|
678 | of the particles that should be handled by the external decay package. |
---|
679 | You should only give the particle (positive) codes; the respective |
---|
680 | antiparticle is always included as well. |
---|
681 | |
---|
682 | <br/><b>Note:</b> The method currently always returns true. |
---|
683 | |
---|
684 | |
---|
685 | <a name="method8"></a> |
---|
686 | <p/><strong>bool Pythia::setRndmEnginePtr( RndmEngine* rndmEnginePtr) </strong> <br/> |
---|
687 | offers the possibility to link to an external random number generator. |
---|
688 | The rules for constructing your own class from the |
---|
689 | <code>RndmEngine</code> base class are described |
---|
690 | <a href="RandomNumbers.html" target="page">here</a>. |
---|
691 | <br/><code>argument</code><strong> rndmEnginePtr </strong> : |
---|
692 | pointer to a <code>RndmEngine</code>-derived object. This object |
---|
693 | must be instantiated by you in your program. |
---|
694 | |
---|
695 | <br/><b>Note:</b> The method returns true if the pointer is different |
---|
696 | from 0. |
---|
697 | |
---|
698 | |
---|
699 | <a name="method9"></a> |
---|
700 | <p/><strong>bool Pythia::setUserHooksPtr( UserHooks* userHooksPtr) </strong> <br/> |
---|
701 | offers the possibility to interact with the generation process at |
---|
702 | a few different specified points, e.g. to reject undesirable events |
---|
703 | at an early stage to save computer time. The rules for constructing |
---|
704 | your own class from the <code>UserHooks</code> base class are described |
---|
705 | <a href="UserHooks.html" target="page">here</a>. You can only hand in one such |
---|
706 | pointer, but this may be to a class that implements several of the |
---|
707 | different allowed possibilities. |
---|
708 | <br/><code>argument</code><strong> userHooksPtr </strong> : |
---|
709 | pointer to a <code>userHooks</code>-derived object. This object |
---|
710 | must be instantiated by you in your program. |
---|
711 | |
---|
712 | <br/><b>Note:</b> The method currently always returns true. |
---|
713 | |
---|
714 | |
---|
715 | <a name="method10"></a> |
---|
716 | <p/><strong>bool Pythia::setBeamShapePtr( BeamShape* beamShapePtr) </strong> <br/> |
---|
717 | offers the possibility to provide your own shape of the momentum and |
---|
718 | space-time spread of the incoming beams. The rules for constructing |
---|
719 | your own class from the <code>BeamShape</code> base class are described |
---|
720 | <a href="BeamShape.html" target="page">here</a>. |
---|
721 | <br/><code>argument</code><strong> BeamShapePtr </strong> : |
---|
722 | pointer to a <code>BeamShape</code>-derived object. This object |
---|
723 | must be instantiated by you in your program. |
---|
724 | |
---|
725 | <br/><b>Note:</b> The method currently always returns true. |
---|
726 | |
---|
727 | |
---|
728 | <a name="method11"></a> |
---|
729 | <p/><strong>bool Pythia::setSigmaPtr( SigmaProcess* sigmaPtr) </strong> <br/> |
---|
730 | offers the possibility to link your own implementation of a process |
---|
731 | and its cross section, to make it a part of the normal process |
---|
732 | generation machinery, without having to recompile the |
---|
733 | <code>Pythia</code> library itself. The rules for constructing your |
---|
734 | own class from the <code>SigmaProcess</code> base class are described |
---|
735 | <a href="SemiInternalProcesses.html" target="page">here</a>. You may call this |
---|
736 | routine repeatedly, to add as many new processes as you wish. |
---|
737 | <br/><code>argument</code><strong> sigmaPtr </strong> : |
---|
738 | pointer to a <code>SigmaProcess</code>-derived object. This object |
---|
739 | must be instantiated by you in your program. |
---|
740 | |
---|
741 | <br/><b>Note:</b> The method currently always returns true. |
---|
742 | |
---|
743 | |
---|
744 | <a name="method12"></a> |
---|
745 | <p/><strong>bool Pythia::setResonancePtr( ResonanceWidths* resonancePtr) </strong> <br/> |
---|
746 | offers the possibility to link your own implementation of the |
---|
747 | calculation of partial resonance widths, to make it a part of the |
---|
748 | normal process generation machinery, without having to recompile the |
---|
749 | <code>Pythia</code> library itself. This allows the decay of new |
---|
750 | resonances to be handled internally, when combined with new particle |
---|
751 | data. Note that the decay of normal hadrons cannot be modelled here; |
---|
752 | this is for New Physics resonances. The rules for constructing your |
---|
753 | own class from the <code>ResonanceWidths</code> base class are described |
---|
754 | <a href="SemiInternalResonances.html" target="page">here</a>. You may call this |
---|
755 | routine repeatedly, to add as many new resonances as you wish. |
---|
756 | <br/><code>argument</code><strong> resonancePtr </strong> : |
---|
757 | pointer to a <code>ResonanceWidths</code>-derived object. This object |
---|
758 | must be instantiated by you in your program. |
---|
759 | |
---|
760 | <br/><b>Note:</b> The method currently always returns true. |
---|
761 | |
---|
762 | |
---|
763 | <a name="method13"></a> |
---|
764 | <p/><strong>bool Pythia::setShowerPtr( TimeShower* timesDecPtr, TimeShower* timesPtr = 0, SpaceShower* spacePtr = 0) </strong> <br/> |
---|
765 | offers the possibility to link your own parton shower routines as |
---|
766 | replacements for the default ones. This is much more complicated |
---|
767 | since the showers are so central and are so interlinked with other |
---|
768 | parts of the program. Therefore it is also possible to do the |
---|
769 | replacement in stages, from the more independent to the more |
---|
770 | intertwined. The rules for constructing your own classes from the |
---|
771 | <code>TimeShower</code> and <code>SpaceShower</code>base classes |
---|
772 | are described <a href="ImplementNewShowers.html" target="page">here</a>. These |
---|
773 | objects must be instantiated by you in your program. |
---|
774 | <br/><code>argument</code><strong> timesDecPtr </strong> : |
---|
775 | pointer to a <code>TimeShower</code>-derived object for doing |
---|
776 | timelike shower evolution in resonance decays, e.g. of a |
---|
777 | <i>Z^0</i>. This is decoupled from beam remnants and parton |
---|
778 | distributions, and is therefore the simplest kind of shower |
---|
779 | to write. If you provide a value 0 then the internal shower |
---|
780 | routine will be used. |
---|
781 | |
---|
782 | <br/><code>argument</code><strong> timesPtr </strong> (<code>default = <strong>0</strong></code>) : |
---|
783 | pointer to a <code>TimeShower</code>-derived object for doing |
---|
784 | all other timelike shower evolution, which is normally interleaved |
---|
785 | with multiparton interactions and spacelike showers, introducing |
---|
786 | both further physics and further technical issues. If you retain |
---|
787 | the default value 0 then the internal shower routine will be used. |
---|
788 | You are allowed to use the same pointer as above for the |
---|
789 | <code>timesDecPtr</code> if the same shower can fulfill both tasks. |
---|
790 | |
---|
791 | <br/><code>argument</code><strong> spacePtr </strong> (<code>default = <strong>0</strong></code>) : |
---|
792 | pointer to a <code>SpaceShower</code>-derived object for doing |
---|
793 | all spacelike shower evolution, which is normally interleaved |
---|
794 | with multiparton interactions and timelike showers. If you retain |
---|
795 | the default value 0 then the internal shower routine will be used. |
---|
796 | |
---|
797 | <br/><b>Note:</b> The method currently always returns true. |
---|
798 | |
---|
799 | |
---|
800 | <h3>Initialize</h3> |
---|
801 | |
---|
802 | At the initialization stage all the information provided above is |
---|
803 | processed, and the stage is set up for the subsequent generation |
---|
804 | of events. Currently several alterative forms of the <code>init</code> |
---|
805 | method are available for this stage, but only the first one is |
---|
806 | recommended. |
---|
807 | |
---|
808 | <a name="method14"></a> |
---|
809 | <p/><strong>bool Pythia::init() </strong> <br/> |
---|
810 | initialize for collisions, in any of the five separate possibilities |
---|
811 | below. In this option the beams are not specified by input arguments, |
---|
812 | but instead by the settings in the |
---|
813 | <a href="BeamParameters.html" target="page">Beam Parameters</a> section. |
---|
814 | This allows the beams to be specified in the same file as other |
---|
815 | run instructions. The default settings give pp collisions at 14 TeV. |
---|
816 | <br/><b>Note:</b> The method returns false if the |
---|
817 | initialization fails. It is then not possible to generate any |
---|
818 | events. |
---|
819 | |
---|
820 | |
---|
821 | <a name="method15"></a> |
---|
822 | <p/><strong>bool Pythia::init( int idA, int idB, double eCM) </strong> <br/> |
---|
823 | initialize for collisions in the center-of-mass frame, with the |
---|
824 | beams moving in the <i>+-z</i> directions. |
---|
825 | <br/><code>argument</code><strong> idA, idB </strong> : |
---|
826 | particle identity code for the two incoming beams. |
---|
827 | |
---|
828 | <br/><code>argument</code><strong> eCM </strong> : |
---|
829 | the CM energy of the collisions. |
---|
830 | |
---|
831 | <br/><b>Notes:</b> Deprecated. The method returns false if the |
---|
832 | initialization fails. It is then not possible to generate any |
---|
833 | events. |
---|
834 | |
---|
835 | |
---|
836 | <a name="method16"></a> |
---|
837 | <p/><strong>bool Pythia::init( int idA, int idB, double eA, double eB) </strong> <br/> |
---|
838 | initialize for collisions with back-to-back beams, |
---|
839 | moving in the <i>+-z</i> directions, but with different energies. |
---|
840 | <br/><code>argument</code><strong> idA, idB </strong> : |
---|
841 | particle identity code for the two incoming beams. |
---|
842 | |
---|
843 | <br/><code>argument</code><strong> eA, eB </strong> : |
---|
844 | the energies of the two beams. If an energy is set to be below |
---|
845 | the mass of the respective beam particle that particle is taken to |
---|
846 | be at rest. This offers a simple possibility to simulate |
---|
847 | fixed-target collisions. |
---|
848 | |
---|
849 | <br/><b>Notes:</b> Deprecated. The method returns false if the |
---|
850 | initialization fails. It is then not possible to generate any |
---|
851 | events. |
---|
852 | |
---|
853 | |
---|
854 | <a name="method17"></a> |
---|
855 | <p/><strong>bool Pythia::init( int idA, int idB, double pxA, double pyA, double pzA, double pxB, double pyB, double pzB) </strong> <br/> |
---|
856 | initialize for collisions with arbitrary beam directions. |
---|
857 | <br/><code>argument</code><strong> idA, idB </strong> : |
---|
858 | particle identity code for the two incoming beams. |
---|
859 | |
---|
860 | <br/><code>argument</code><strong> pxA, pyA, pzA </strong> : |
---|
861 | the three-momntum vector <i>(p_x, p_y, p_z)</i> of the first |
---|
862 | incoming beam. |
---|
863 | |
---|
864 | <br/><code>argument</code><strong> pxB, pyB, pzB </strong> : |
---|
865 | the three-momntum vector <i>(p_x, p_y, p_z)</i> of the second |
---|
866 | incoming beam. |
---|
867 | |
---|
868 | <br/><b>Notes:</b> Deprecated. The method returns false if the |
---|
869 | initialization fails. It is then not possible to generate any |
---|
870 | events. |
---|
871 | |
---|
872 | |
---|
873 | <a name="method18"></a> |
---|
874 | <p/><strong>bool Pythia::init( string LesHouchesEventFile, bool skipInit = false) </strong> <br/> |
---|
875 | initialize for hard-process collisions fed in from an external file |
---|
876 | with events, written according to the |
---|
877 | <a href="LesHouchesAccord.html" target="page">Les Houches Event File</a> |
---|
878 | standard. |
---|
879 | <br/><code>argument</code><strong> LesHouchesEventFile </strong> : |
---|
880 | the file name (including path, where required) where the |
---|
881 | events are stored, including relevant information on beam |
---|
882 | identities and energies. |
---|
883 | |
---|
884 | <br/><code>argument</code><strong> skipInit </strong> (<code>default = <strong>false</strong></code>) : |
---|
885 | By default this method does a complete reinitialization of the |
---|
886 | generation process. If you set this argument to true then |
---|
887 | no reinitialization will occur, only the pointer to the event |
---|
888 | file is updated. This may come in handy if the full event sample |
---|
889 | is split across several files generated under the same conditions |
---|
890 | (except random numbers, of course). You then do the first |
---|
891 | initialization with the default, and all subsequent ones with |
---|
892 | true. Note that things may go wrong if the files are not created |
---|
893 | under the same conditions. |
---|
894 | |
---|
895 | <br/><b>Notes:</b> Deprecated. The method returns false if the |
---|
896 | initialization fails. It is then not possible to generate any |
---|
897 | events. |
---|
898 | |
---|
899 | |
---|
900 | <a name="method19"></a> |
---|
901 | <p/><strong>bool Pythia::init( LHAup* lhaUpPtr) </strong> <br/> |
---|
902 | initialize for hard-process collisions fed in from an external |
---|
903 | source of events, consistent with the Les Houches Accord standard. |
---|
904 | The rules for constructing your own class from the <code>LHAup</code> |
---|
905 | base class are described <a href="LesHouchesAccord.html" target="page">here</a>. |
---|
906 | This class is also required to provide the beam parameters. |
---|
907 | <br/><code>argument</code><strong> lhaUpPtr </strong> : |
---|
908 | pointer to a <code>LHAup</code>-derived object. This object |
---|
909 | must be instantiated by you in your program. |
---|
910 | |
---|
911 | <br/><b>Notes:</b> Deprecated. The method returns false if the |
---|
912 | initialization fails. It is then not possible to generate any |
---|
913 | events. |
---|
914 | |
---|
915 | |
---|
916 | <h3>Generate events</h3> |
---|
917 | |
---|
918 | The <code>next()</code> method is the main one to generate events. |
---|
919 | In this section we also put a few other specialized methods that |
---|
920 | may be useful in some circumstances. |
---|
921 | |
---|
922 | <a name="method20"></a> |
---|
923 | <p/><strong>bool Pythia::next() </strong> <br/> |
---|
924 | generate the next event. No input parameters are required; all |
---|
925 | instructions have already been set up in the initialization stage. |
---|
926 | <br/><b>Note:</b> The method returns false if the event generation |
---|
927 | fails. The event record is then not consistent and should not be |
---|
928 | studied. When reading in hard collisions from a Les Houches Event File |
---|
929 | the problem may be that the end of the file has been reached. This |
---|
930 | can be checked with the |
---|
931 | <code><a href="EventInformation.html" target="page">Info::atEndOfFile()</a></code> |
---|
932 | method. |
---|
933 | |
---|
934 | |
---|
935 | <a name="method21"></a> |
---|
936 | <p/><strong>int Pythia::forceTimeShower( int iBeg, int iEnd, double pTmax, int nBranchMax = 0) </strong> <br/> |
---|
937 | perform a final-state shower evolution on partons in the |
---|
938 | <code>event</code> event record. This could be used for externally |
---|
939 | provided simple events, or even parts of events, for which |
---|
940 | a complete generation is not foreseen. Since the mother source of |
---|
941 | the parton system is not known, one cannot expect as good accuracy |
---|
942 | as in a normal generation. When two different timelike shower |
---|
943 | instances are set up, it is the one used for showering in resonance |
---|
944 | decays that is used here. The <code>forceTimeShower</code> method |
---|
945 | can be used in conjunction with the <code>forceHadronLevel</code> |
---|
946 | one below. Further comments are found |
---|
947 | <a href="HadronLevelStandalone.html" target="page">here</a>. |
---|
948 | <br/><code>argument</code><strong> iBeg, iEnd </strong> : the first and last entry of the event |
---|
949 | record to be affected by the call. |
---|
950 | |
---|
951 | <br/><code>argument</code><strong> pTmax </strong> : the maximum <i>pT</i> scale of emissions. |
---|
952 | Additionally, as always, the <code>scale</code> variable of each parton |
---|
953 | sets the maximum <i>pT</i> scale of branchings of this parton. |
---|
954 | Recall that this scale defaults to 0 if not set, so that no radiation |
---|
955 | can occur. |
---|
956 | |
---|
957 | <br/><code>argument</code><strong> nBranchMax </strong> (<code>default = <strong>0</strong></code>) : when positive, it sets the |
---|
958 | maximum number of branchings that are allowed to occur in the shower, |
---|
959 | i.e. the shower may stop evolving before reaching the lower cutoff. |
---|
960 | The argument has no effect when zero or negative, i.e. then the shower |
---|
961 | will continue to the lower cutoff. |
---|
962 | |
---|
963 | <br/><b>Note:</b> The method returns the number of branchings that |
---|
964 | has been generated. |
---|
965 | |
---|
966 | |
---|
967 | <a name="method22"></a> |
---|
968 | <p/><strong>bool Pythia::forceHadronLevel(bool findJunctions = true) </strong> <br/> |
---|
969 | hadronize the existing event record, i.e. perform string fragmentation |
---|
970 | and particle decays. There are two main applications. Firstly, |
---|
971 | you can use the same parton-level content as a basis for repeated |
---|
972 | hadronization attempts, in schemes intended to save computer time. |
---|
973 | Secondly, you may have an external program that can simulate the full |
---|
974 | partonic level of the event - hard process, parton showers, multiparton |
---|
975 | interactions, beam remnants, colour flow, and so on - but not |
---|
976 | hadronization. Further details are found |
---|
977 | <a href="HadronLevelStandalone.html" target="page">here</a>. |
---|
978 | <br/><code>argument</code><strong> findJunctions </strong> (<code>default = <strong>true</strong></code>) : |
---|
979 | normally this routine will search through the event record and try to |
---|
980 | figure out if any colour junctions are present. If so, the colour |
---|
981 | topology of such junctions must be sorted out. In tricky cases this |
---|
982 | might fail, and then hadronization will not work. A user who is |
---|
983 | aware of this and knows the intended colour flow can set up the |
---|
984 | junction information in the event record, and then call |
---|
985 | <code>forceHadronLevel(false)</code> so as not to have this information |
---|
986 | overwritten. |
---|
987 | <br/><b>Note:</b> The method returns false if the hadronization |
---|
988 | fails. The event record is then not consistent and should not be |
---|
989 | studied. |
---|
990 | |
---|
991 | |
---|
992 | <a name="method23"></a> |
---|
993 | <p/><strong>bool Pythia::moreDecays() </strong> <br/> |
---|
994 | perform decays of all particles in the event record that have not been |
---|
995 | decayed but should have been done so. This can be used e.g. for |
---|
996 | repeated decay attempts, in schemes intended to save computer time. |
---|
997 | Further details are found <a href="HadronLevelStandalone.html" target="page">here</a>. |
---|
998 | <br/><b>Note:</b> The method returns false if the decays fail. The |
---|
999 | event record is then not consistent and should not be studied. |
---|
1000 | |
---|
1001 | |
---|
1002 | <a name="method24"></a> |
---|
1003 | <p/><strong>bool Pythia::forceRHadronDecays() </strong> <br/> |
---|
1004 | perform decays of R-hadrons that were previously considered stable. |
---|
1005 | This could be if an R-hadron is sufficiently long-lived that |
---|
1006 | it may interact in the detector between production and decay, so that |
---|
1007 | its four-momentum is changed. Further details are found |
---|
1008 | <a href="RHadrons.html" target="page">here</a>. |
---|
1009 | <br/><b>Note:</b> The method returns false if the decays fail. The |
---|
1010 | event record is then not consistent and should not be studied. |
---|
1011 | |
---|
1012 | |
---|
1013 | <a name="method25"></a> |
---|
1014 | <p/><strong>void Pythia::LHAeventList(ostream& os = cout) </strong> <br/> |
---|
1015 | list the Les Houches Accord information on the current event, see |
---|
1016 | <code><a href="LesHouchesAccord.html" target="page">LHAup::listEvent(...)</a></code>. |
---|
1017 | (Other listings are available via the class members below, so this |
---|
1018 | listing is a special case that would not fit elsewhere.) |
---|
1019 | <br/><code>argument</code><strong> os </strong> (<code>default = <strong>cout</strong></code>) : |
---|
1020 | output stream where the listing occurs. |
---|
1021 | |
---|
1022 | |
---|
1023 | |
---|
1024 | <a name="method26"></a> |
---|
1025 | <p/><strong>bool Pythia::LHAeventSkip(int nSkip) </strong> <br/> |
---|
1026 | skip ahead a number of events in the Les Houches generation |
---|
1027 | sequence, without doing anything further with them, see |
---|
1028 | <code><a href="LesHouchesAccord.html" target="page">LHAup::skipEvent(nSkip)</a></code>. |
---|
1029 | Mainly intended for debug purposes, e.g. when an event at a known |
---|
1030 | location in a Les Houches Event File is causing problems. |
---|
1031 | <br/><code>argument</code><strong> nSkip </strong> : |
---|
1032 | number of events to skip. |
---|
1033 | |
---|
1034 | <br/><b>Note:</b> The method returns false if the operation fails, |
---|
1035 | specifically if the end of a LHEF has been reached, cf. |
---|
1036 | <code>next()</code> above. |
---|
1037 | |
---|
1038 | |
---|
1039 | <h3>Finalize</h3> |
---|
1040 | |
---|
1041 | There is no required finalization step; you can stop generating events |
---|
1042 | when and how you want. It is still recommended that you make it a |
---|
1043 | routine to call the following method at the end. A second method provides |
---|
1044 | a deprecated alternative. |
---|
1045 | |
---|
1046 | <a name="method27"></a> |
---|
1047 | <p/><strong>void Pythia::stat() </strong> <br/> |
---|
1048 | list statistics on the event generation, specifically total and partial |
---|
1049 | cross sections and the number of different errors. For more details see |
---|
1050 | <a href="EventStatistics.html" target="page">here</a> and for available options |
---|
1051 | <a href="MainProgramSettings.html" target="page">here</a>. |
---|
1052 | |
---|
1053 | |
---|
1054 | <a name="method28"></a> |
---|
1055 | <p/><strong>void Pythia::statistics(bool all = false, bool reset = false) </strong> <br/> |
---|
1056 | list statistics on the event generation, specifically total and partial |
---|
1057 | cross sections and the number of different errors. For more details see |
---|
1058 | <a href="EventStatistics.html" target="page">here</a>. |
---|
1059 | <br/><code>argument</code><strong> all </strong> (<code>default = <strong>false</strong></code>) : |
---|
1060 | if true also statistics on multiparton interactions is shown, by default not. |
---|
1061 | |
---|
1062 | <br/><code>argument</code><strong> reset </strong> (<code>default = <strong>false</strong></code>) : if true then all counters, |
---|
1063 | e.g on events generated and errors experienced, are reset to zero |
---|
1064 | whenever the routine is called. The default instead is that |
---|
1065 | all stored statistics information is unaffected by the call. Counters |
---|
1066 | are automatically reset in each new <code>Pythia::init(...)</code> |
---|
1067 | call, however, so the only time the <code>reset</code> option makes a |
---|
1068 | difference is if <code>statistics(...)</code> is called several times |
---|
1069 | in a (sub)run. |
---|
1070 | |
---|
1071 | <br/><b>Note:</b> Deprecated. |
---|
1072 | |
---|
1073 | |
---|
1074 | <h3>Interrogate settings</h3> |
---|
1075 | |
---|
1076 | Normally settings are used in the setup and initialization stages |
---|
1077 | to determine the character of a run, e.g. read from a file with the |
---|
1078 | above-described <code>Pythia::readFile(...)</code> method. |
---|
1079 | There is no strict need for a user to interact with the |
---|
1080 | <code>Settings</code> database in any other way. However, as an option, |
---|
1081 | some settings variables have been left free for the user to set in |
---|
1082 | such a file, and then use in the main program to directly affect the |
---|
1083 | performance of that program, see |
---|
1084 | <a href="MainProgramSettings.html" target="page">here</a>. A typical example would |
---|
1085 | be the number of events to generate. For such applications the |
---|
1086 | following shortcuts to some <code>Settings</code> methods may be |
---|
1087 | convenient. |
---|
1088 | |
---|
1089 | <a name="method29"></a> |
---|
1090 | <p/><strong>bool Pythia::flag(string key) </strong> <br/> |
---|
1091 | read in a boolean variable from the <code>Settings</code> database. |
---|
1092 | <br/><code>argument</code><strong> key </strong> : |
---|
1093 | the name of the variable to be read. |
---|
1094 | |
---|
1095 | |
---|
1096 | |
---|
1097 | <a name="method30"></a> |
---|
1098 | <p/><strong>int Pythia::mode(string key) </strong> <br/> |
---|
1099 | read in an integer variable from the <code>Settings</code> database. |
---|
1100 | <br/><code>argument</code><strong> key </strong> : |
---|
1101 | the name of the variable to be read. |
---|
1102 | |
---|
1103 | |
---|
1104 | |
---|
1105 | <a name="method31"></a> |
---|
1106 | <p/><strong>double Pythia::parm(string key) </strong> <br/> |
---|
1107 | read in a double-precision variable from the <code>Settings</code> |
---|
1108 | database. |
---|
1109 | <br/><code>argument</code><strong> key </strong> : |
---|
1110 | the name of the variable to be read. |
---|
1111 | |
---|
1112 | |
---|
1113 | |
---|
1114 | <a name="method32"></a> |
---|
1115 | <p/><strong>string Pythia::word(string key) </strong> <br/> |
---|
1116 | read in a string variable from the <code>Settings</code> database. |
---|
1117 | <br/><code>argument</code><strong> key </strong> : |
---|
1118 | the name of the variable to be read. |
---|
1119 | |
---|
1120 | |
---|
1121 | |
---|
1122 | <h3>Data members</h3> |
---|
1123 | |
---|
1124 | The <code>Pythia</code> class contains a few public data members, |
---|
1125 | several of which play a central role. We list them here, with |
---|
1126 | links to the places where they are further described. |
---|
1127 | |
---|
1128 | <a name="method33"></a> |
---|
1129 | <p/><strong>Event Pythia::process </strong> <br/> |
---|
1130 | the hard-process event record, see <a href="EventRecord.html" target="page">here</a> |
---|
1131 | for further details. |
---|
1132 | |
---|
1133 | |
---|
1134 | <a name="method34"></a> |
---|
1135 | <p/><strong>Event Pythia::event </strong> <br/> |
---|
1136 | the complete event record, see <a href="EventRecord.html" target="page">here</a> |
---|
1137 | for further details. |
---|
1138 | |
---|
1139 | |
---|
1140 | <a name="method35"></a> |
---|
1141 | <p/><strong>Info Pythia::info </strong> <br/> |
---|
1142 | further information on the event-generation process, see |
---|
1143 | <a href="EventInformation.html" target="page">here</a> for further details. |
---|
1144 | |
---|
1145 | |
---|
1146 | <a name="method36"></a> |
---|
1147 | <p/><strong>Settings Pythia::settings </strong> <br/> |
---|
1148 | the settings database, see <a href="SettingsScheme.html" target="page">here</a> |
---|
1149 | for further details. |
---|
1150 | |
---|
1151 | |
---|
1152 | <a name="method37"></a> |
---|
1153 | <p/><strong>ParticleData Pythia::particleData </strong> <br/> |
---|
1154 | the particle properties and decay tables database, see |
---|
1155 | <a href="ParticleDataScheme.html" target="page">here</a> for further details. |
---|
1156 | |
---|
1157 | |
---|
1158 | <a name="method38"></a> |
---|
1159 | <p/><strong>Rndm Pythia::rndm </strong> <br/> |
---|
1160 | the random number generator, see <a href="RandomNumberSeed.html" target="page">here</a> |
---|
1161 | and <a href="RandomNumbers.html" target="page">here</a> for further details. |
---|
1162 | |
---|
1163 | |
---|
1164 | <a name="method39"></a> |
---|
1165 | <p/><strong>CoupSM Pythia::coupSM </strong> <br/> |
---|
1166 | Standard Model couplings and mixing matrices, see |
---|
1167 | <a href="StandardModelParameters.html" target="page">here</a> for further details. |
---|
1168 | |
---|
1169 | |
---|
1170 | <a name="method40"></a> |
---|
1171 | <p/><strong>SusyLesHouches Pythia::slha </strong> <br/> |
---|
1172 | parameters and particle data in the context of supersymmetric models, |
---|
1173 | see <a href="SUSYLesHouchesAccord.html" target="page">here</a> for further details. |
---|
1174 | |
---|
1175 | |
---|
1176 | <a name="method41"></a> |
---|
1177 | <p/><strong>PartonSystems Pythia::partonSystems </strong> <br/> |
---|
1178 | a grouping of the partons in the event record by subsystem, |
---|
1179 | see <a href="AdvancedUsage.html" target="page">here</a> for further details. |
---|
1180 | |
---|
1181 | |
---|
1182 | </body> |
---|
1183 | </html> |
---|
1184 | |
---|
1185 | <!-- Copyright (C) 2012 Torbjorn Sjostrand --> |
---|