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