Update to pythi8.170
[u/mrichter/AliRoot.git] / PYTHIA8 / pythia8170 / xmldoc / HadronLevelStandalone.xml
1 <chapter name="Hadron-Level Standalone">
2
3 <h2>Hadron-Level Standalone</h2>
4
5 The Les Houches Accord allows external process-level configurations
6 to be fed in, for subsequent parton-level and hadron-level generation
7 to be handled internally by PYTHIA. There is no correspondingly 
8 standardized interface if you have external events that have also 
9 been generated through the parton-level stage, so that only the 
10 hadron-level remains to be handled. A non-standard way to achieve this
11 exists, however, and can be useful both for real applications and
12 for various tests of the hadronization model on its own. 
13
14 <p/>
15 The key trick is to set the flag <code>ProcessLevel:all = off</code>.
16 When <code>pythia.next()</code> is called it then does not try to
17 generate a hard process. Since there are no beams, it is also not 
18 possible to perform the normal <code>PartonLevel</code> step.
19 (It is still possible to generate final-state radiation, but this
20 is not automatic. It would have to be done by hand, using the 
21 <code>pythia.forceTimeShower(...)</code> method, before 
22 <code>pythia.next()</code> is called.) 
23 Thus only the <code>HadronLevel</code> methods are 
24 called, to take the current content of the event record stored in
25 <code>pythia.event</code> as a starting point for any hadronization 
26 and decays that are allowed by the normal parameters of this step. 
27 Often the input would consist solely of partons grouped into colour 
28 singlets, but also (colour-singlet) particles are allowed.
29
30 <p/>
31 To set up all the parameters, a <code>pythia.init()</code> call has
32 to be used, without any arguments. In brief, the structure of the
33 main program therefore should be something like
34 <pre>
35   Pythia pythia;                               // Declare generator.
36   Event& event = pythia.event                  // Convenient shorthand.
37   pythia.readString("ProcessLevel:all = off"); // The trick!
38   pythia.init();                               // Initialization.
39   for (int iEvent = 0; iEvent < nEvent; ++iEvent) {
40     // Insert filling of event here!
41     pythia.next();                             // Do the hadron level.
42   }
43 </pre> 
44 Of course this should be supplemented by analysis of events, error checks,
45 and so on, as for a normal PYTHIA run. The unique aspect is how to fill
46 the <code>event</code> inside the loop, before <code>pythia.next()</code> 
47 is called.
48
49 <h3>Input configuration</h3>
50
51 To set up a new configuration the first step is to throw away the current
52 one, with <code>event.reset()</code>. This routine will also reserve
53 the zeroth entry in the even record to represent the event as a whole. 
54
55 <p/>
56 With the <code>event.append(...)</code> methods a new entry is added at the
57 bottom of the current record, i.e. the first time it is called entry 
58 number 1 is filled, and so on. The <code>append</code>  method basically
59 exists in four variants, either without or with history information,
60 and with four-momentum provided either as a <code>Vec4</code> four-vector 
61 or as four individual components:
62 <pre>
63   append( id, status, col, acol, p, m)
64   append( id, status, col, acol, px, py, pz, e, m)
65   append( id, status, mother1, mother2, daughter1, daughter2, col, acol, p, m)
66   append( id, status, mother1, mother2, daughter1, daughter2, col, acol, px, py, pz, e, m)
67 </pre> 
68 The methods return the index at which the entry has been stored,
69 but normally you would not use this feature. 
70
71 <p/>
72 You can find descriptions of the input variables 
73 <aloc href="ParticleProperties">here</aloc>. 
74 The PDG particle code <code>id</code> and the Les Houches Accord colour 
75 <code>col</code> and anticolour <code>acol</code> tags must be set
76 correctly. The four-momentum and mass have to be provided in units of GeV; 
77 if you omit the mass it defaults to 0. 
78
79 <p/>
80 Outgoing particles that should hadronize should be given status code 23.
81 Often this is the only status code you need. You could e.g. also fill in 
82 incoming partons with -21 and intermediate ones with -22, if you so wish.
83 Usually the choice of status codes is not crucial, so long as you recall 
84 that positive numbers correspond to particles that are still around, while
85 negative numbers denote ones that already hadronized or decayed. However,
86 so as not to run into contradictions with the internal PYTHIA checks
87 (when <code>Check:event = on</code>), or with external formats such as 
88 HepMC, we do recommend the above codes. When <code>pythia.next()</code> 
89 is called the positive-status particles that hadronize/decay get the sign 
90 of the status code flipped to negative but the absolute value is retained. 
91 The new particles are added with normal PYTHIA status codes.
92
93 <p/>
94 For normal hadronization/decays in <code>pythia.next()</code> the
95 history encoded in the mother and daughter indices is not used. 
96 Therefore the first two <code>append</code> methods, which set all these 
97 indices vanishing, should suffice. The subsequent hadronization/decays 
98 will still be properly documented.
99
100 <p/>
101 The exception is when you want to include junctions in your string 
102 topology, i.e. have three string pieces meet. Then you must insert in
103 your event record the (decayed) particle that is the reason for the
104 presence of a junction, e.g. a baryon beam remnant from which several
105 valence quarks have been kicked out, or a neutralino that underwent a
106 baryon-number-violating decay. This particle must have as daughters 
107 the three partons that together carry the baryon number. 
108  
109 <p/>
110 The sample program in <code>main21.cc</code> illustrates how you can work
111 with this facility, both for simple parton configurations and for more
112 complicated ones with junctions.
113
114 <p/>
115 As an alternative to setting up a topology with the methods above, 
116 a <aloc href="LesHouchesAccord">Les Houches Event File</aloc> (LHEF) 
117 can also provide the configurations. Since no beams or processes are 
118 defined, only the <code>&lt;event&gt;....&lt;/event&gt;</code> blocks 
119 need to be present, one for each event, so strictly it is not a true 
120 LHEF. You need to select <code>Beams:frameType = 4</code>, provide 
121 the file name in <code>Beams:LHEF</code> and, as above, set
122 <code>ProcessLevel:all = off</code>. Needless to say, an externally 
123 linked <code>LHAup</code> class works as well as an LHEF,
124 with <code>Beams:frameType = 5</code>.
125
126 <p/>
127 The event information to store in the LHEF, or provide by the 
128 <code>LHAup</code>, is essentially the same as above. The only 
129 difference is in status codes: outgoing particles should have 1 
130 instead of 23, and intermediate resonances 2 instead of -22. 
131 Incoming partons, if any, are -1 instead of -21.
132
133 <h3>Extensions to resonance decays</h3>
134
135 With the above scheme, <code>pythia.next()</code> will generate
136 hadronization, i.e. string fragmentation and subsequent decays of
137 normal unstable particles. It will not decay 
138 <aloc href="ResonanceDecays">resonances</aloc>, i.e.
139 <ei>W, Z</ei>, top, Higgs, SUSY and other massive particles.
140
141 <p/>
142 If the decay products are already provided, of course those 
143 products will be hadronized, but without any of the showers
144 that one would expect in <ei>Z^0 -> q qbar</ei>, say. That is, the
145 presence of a decayed resonance in the event record can be nice
146 for documentation purposes, but otherwise it plays no role. 
147 It is possible to change this behaviour with the following flag. 
148
149 <flag name="Standalone:allowResDec" default="off">
150 If off, then resonances are stable if final-state particles,
151 and irrelevant if intermediate particles. If on, then 
152 resonances are decayed, sequentially if necessary. After each
153 decay step a parton shower is applied, and subsequent decay 
154 kinematics is shifted by the recoil of such branchings. 
155 If the decay chain and kinematics of the resonance is provided
156 at input, this information is used, otherwise it is generated
157 internally according to built-in branching ratios etc.
158 </flag>
159
160 <p/>
161 The input configuration has to follow the rules described above, 
162 with <code>ProcessLevel:all = off</code>. (Terminology not quite 
163 consistent, since resonance decays normally is part of the 
164 process-level step.) It is possible to combine several resonances, 
165 and other coloured or uncoloured particles into the same event.
166  
167 <p/>
168 Final-state radiation is applied to each step of the decay 
169 sequence by default, but can be switched off with 
170 <code>PartonLevel:FSR = off</code> or 
171 <code>PartonLevel:FSRinResonances</code>. 
172
173 <h3>Repeated hadronization or decay</h3>
174
175 An alternative approach is possible with the 
176 <code>pythia.forceHadronLevel()</code> routine. This method does
177 a call to the <code>HadronLevel</code> methods, irrespective of the
178 value of the <code>HadronLevel:all</code> flag. If you hadronize
179 externally generated events it is equivalent to a 
180 <code>pythia.next()</code> call with 
181 <code>ProcessLevel:all = off</code>. 
182
183 <p/>
184 This method truly sticks to the hadron level, and thus cannot handle 
185 resonance decays. You therefore <b>must not</b> mix it with the
186 <code>Standalone:allowResDec = on</code> framework. 
187
188 <p/>
189 The similarity of names indicates that 
190 <code>pythia.forceTimeShower( int iBeg, int iEnd, double pTmax, 
191 int nBranchMax = 0)</code> is intended to belong to the same set of 
192 work-by-hand methods. Here <code>iBeg</code> and <code>iEnd</code>
193 give the range of partons that should be allowed to shower, 
194 <code>pTmax</code> the maximum <ei>pT</ei> scale of emissions,
195 and a nonzero <code>nBranchMax</code> a maximum number of allowed
196 branchings. Additionally, a <code>scale</code> has to be set for each
197 parton that should shower, which requires an additional final argument 
198 to the <code>append</code> methods above. This scale limits the maximum 
199 <ei>pT</ei> allowed for each parton, in addition to the global 
200 <code>pTmax</code>. When not set the scale defaults to 0, meaning no 
201 radiation for that parton.
202
203 <p/>
204 The real application instead is for repeated hadronization of the same
205 PYTHIA process- and parton-level event. This may for some studies
206 help to save time, given that these two first step are more 
207 time-consuming than the hadronization one. 
208
209 <p/>
210 For repeated hadronization you should first generate an event as usual, 
211 but with <code>HadronLevel:all = off</code>. This event you can save
212 in a temporary copy, e.g. <code>Event savedEvent = pythia.event</code>.
213 Inside a loop you copy back with <code>pythia.event = savedEvent</code>, 
214 and call <code>pythia.forceHadronLevel()</code> to obtain a new 
215 hadronization history.
216
217 <p/>
218 A more limited form of repetition is if you want to decay a given kind 
219 of particle repeatedly, without having to generate the rest of the event 
220 anew. This could be the case e.g. in <ei>B</ei> physics applications. 
221 Then you can use the <code>pythia.moreDecays()</code> method, which 
222 decays all particles in the event record that have not been decayed 
223 but should have been done so. The 
224 <code>pythia.particleData.mayDecay( id, false/true)</code> method 
225 may be used to switch off/on the decays of a particle species 
226 <code>id</code>, so that it is not decayed in the 
227 <code>pythia.next()</code> call but only inside a loop over a number of
228 tries. 
229
230 <p/>
231 Between each loop the newly produced decay products must be 
232 removed and the decayed particle status restored to undecayed.
233 The former is simple, since the new products are appended to the
234 end of the event record: <code>event.saveSize()</code> saves the
235 initial size of the event record, and <code>event.restoreSize()</code>
236 can later be used repeatedly to restore this original size, which means
237 that the new particles at the end are thrown away. The latter is more
238 complicated, and requires the user to identify the positions of all
239 particles of the species and restore a positive status code with
240 <code>event[i].statusPos()</code>.
241
242 <p/>
243 The <code>main15.cc</code> program illustrates both these methods,
244 i.e. either repeated hadronization or repeated decay of PYTHIA
245 events.
246
247 </chapter>
248
249 <!-- Copyright (C) 2012 Torbjorn Sjostrand -->