1 <chapter name="Event Analysis">
3 <h2>Event Analysis</h2>
5 <h3>Introduction</h3>
7 The routines in this section are intended to be used to analyze
8 event properties. As such they are not part of the main event
9 generation chain, but can be used in comparisons between Monte
10 Carlo events and real data. They are rather free-standing, but
11 assume that input is provided in the PYTHIA 8
12 <code>Event</code> format, and use a few basic facilities such
13 as four-vectors.
15 <h3>Sphericity</h3>
17 The standard sphericity tensor is
18 <eq>
19     S^{ab} = (sum_i p_i^a p_i^b) / (sum_i p_i^2)
20 </eq>
21 where the <ei>sum i</ei> runs over the particles in the event,
22 <ei>a, b = x, y, z,</ei> and <ei>p</ei> without such an index is
23 the absolute size of the three-momentum . This tensor can be
24 diagonalized to find eigenvalues and eigenvectors.
26 <p/>
27 The above tensor can be generalized by introducing a power
28 <ei>r</ei>, such that
29 <eq>
30     S^{ab} = (sum_i p_i^a p_i^b p_i^{r-2}) / (sum_i p_i^r)
31 </eq>
32 In particular, <ei>r = 1</ei> gives a linear dependence on momenta
33 and thus a collinear safe definition, unlike sphericity.
35 <p/>
36 A sphericity analysis object is declared by
37 <class name="Sphericity sph( power, select)">
38 where
39 <argument name="power" default="2.">
40 is the power <ei>r</ei> defined above, i.e.
41 <argoption value="2.">gives Spericity, and</argoption>
42 <argoption value="1.">gives the linear form.</argoption>
43 </argument>
44 <argument name="select" default="2">
45 tells which particles are analyzed,
46 <argoption value="1">all final-state particles,</argoption>
47 <argoption value="2">all observable final-state particles,
48 i.e. excluding neutrinos and other particles without strong or
49 electromagnetic interactions (the <code>isVisible()</code>
50 particle method), and
51 </argoption>
52 <argoption value="3">only charged final-state particles.</argoption>
53 </argument>
54 </class>
56 <p/>
57 The analysis is performed by a call to the method
58 <method name="analyze( event)">
59 where
60 <argument name="event">is an object of the <code>Event</code> class,
61 most likely the <code>pythia.event</code> one.
62 </argument>
63 <br/>If the routine returns <code>false</code> the analysis failed,
64 e.g. if too few particles are present to analyze.
65 </method>
67 <p/>
68 After the analysis has been performed, a few <code>Sphericity</code>
69 class methods are available to return the result of the analysis:
70 <method name="sphericity()">
71 gives the sphericity (or equivalent if <ei>r</ei> is not 2),
72 </method>
73 <method name="aplanarity()">
74 gives the aplanarity (with the same comment),
75 </method>
76 <method name="eigenValue(i)">
77 gives one of the three eigenvalues for <ei>i</ei> = 1, 2 or 3, in
78 descending order,
79 </method>
80 <method name="EventAxis(i)">
81 gives the matching eigenvector, as a <code>Vec4</code> with vanishing
82 time/energy component.
83 </method>
84 <method name="list()">
85 provides a listing of the above information.
86 </method>
87 <method name="nError()">
88 tells the number of times <code>analyze</code> failed to analyze events.
89 </method>
91 <h3>Thrust</h3>
93 Thrust is obtained by varying the thrust axis so that the longitudinal
94 momentum component projected onto it is maximized, and thrust itself is
95 then defined as the sum of absolute longitudinal momenta divided by
96 the sum of absolute momenta. The major axis is found correspondingly
97 in the plane transverse to thrust, and the minor one is then defined
98 to be transverse to both. Oblateness is the difference between the major
99 and the minor values.
101 <p/>
102 The calculation of thrust is more computer-time-intensive than e.g.
103 linear sphericity, introduced above, and has no specific advantages except
104 historical precedent. In the PYTHIA 6 implementation the search was
105 speeded up at the price of then not being guaranteed to hit the absolute
106 maximum. The current implementation studies all possibilities, but at
107 the price of being slower, with time consumption for an event with
108 <ei>n</ei> particles growing like <ei>n^3</ei>.
110 <p/>
111 A thrust analysis object is declared by
112 <class name="Thrust thr( select)">
113 where
114 <argument name="select" default="2">
115 tells which particles are analyzed,
116 <argoption value="1">all final-state particles,</argoption>
117 <argoption value="2">all observable final-state particles,
118 i.e. excluding neutrinos and other particles without strong or
119 electromagnetic interactions (the <code>isVisible()</code>
120 particle method), and
121 </argoption>
122 <argoption value="3">only charged final-state particles.</argoption>
123 </argument>
124 </class>
126 <p/>
127 The analysis is performed by a call to the method
128 <method name="analyze( event)">
129 where
130 <argument name="event">is an object of the <code>Event</code> class,
131 most likely the <code>pythia.event</code> one.
132 </argument>
133 <br/>If the routine returns <code>false</code> the analysis failed,
134 e.g. if too few particles are present to analyze.
135 </method>
137 <p/>
138 After the analysis has been performed, a few <code>Thrust</code>
139 class methods are available to return the result of the analysis:
140 <method name="thrust(), tMajor(), tMinor(), oblateness()">
141 gives the thrust, major, minor and oblateness values, respectively,
142 <method name="EventAxis(i)">
143 gives the matching event-axis vectors, for <ei>i</ei> = 1, 2 or 3
144 corresponding to thrust, major or minor, as a <code>Vec4</code> with
145 vanishing time/energy component.
146 </method>
147 <method name="list()">
148 provides a listing of the above information.
149 </method>
150 <method name="nError()">
151 tells the number of times <code>analyze</code> failed to analyze events.
152 </method>
154 <h3>ClusterJet</h3>
156 <code>ClusterJet</code> (a.k.a. <code>LUCLUS</code> and
157 <code>PYCLUS</code>) is a clustering algorithm of the type used for
158 analyses of <ei>e^+e^-</ei> events, see the PYTHIA 6 manual. A few
159 options are available for some well-known distance measures. Cutoff
160 distances can either be given in terms of a scaled quadratic quantity
161 like <ei>y = pT^2/E^2</ei> or an unscaled linear one like <ei>pT</ei>.
163 <p/>
164 A cluster-jet analysis object is declared by
165 <class name="ClusterJet clusterJet( measure, select, massSet,
166 precluster, reassign)">
167 where
168 <argument name="measure" default="Lund">distance measure, to be provided
169 as a character string (actually, only the first character is necessary)
170 <argoption value="Lund">the Lund <ei>pT</ei> distance,</argoption>
171 <argoption value="JADE">the JADE mass distance, and</argoption>
172 <argoption value="Durham">the Durham <ei>kT</ei> measure.</argoption>
173 </argument>
174 <argument name="select" default="2">
175 tells which particles are analyzed,
176 <argoption value="1">all final-state particles,</argoption>
177 <argoption value="2">all observable final-state particles,
178 i.e. excluding neutrinos and other particles without strong or
179 electromagnetic interactions (the <code>isVisible()</code> particle
180 method), and
181 </argoption>
182 <argoption value="3">only charged final-state particles.</argoption>
183 </argument>
184 <argument name="massSet" default="2">masses assumed for the particles
185 used in the analysis
186 <argoption value="0">all massless,</argoption>
187 <argoption value="1">photons are massless while all others are
188 assigned the <ei>pi+-</ei> mass, and
189 </argoption>
190 <argoption value="2">all given their correct masses.</argoption>
191 </argument>
192 <argument name="precluster" default="false">
193 perform or not a early preclustering step, where nearby particles
194 are lumped together so as to speed up the subsequent normal clustering.
195 </argument>
196 <argument name="reassign" default="false">
197 reassign all particles to the nearest jet each time after two jets
198 have been joined.
199 </argument>
200 </class>
202 <p/>
203 The analysis is performed by a
204 <method name="analyze( event, yScale, pTscale, nJetMin, nJetMax)">
205 where
206 <argument name="event">is an object of the <code>Event</code> class,
207 most likely the <code>pythia.event</code> one.
208 </argument>
209 <argument name="yScale">
210 is the cutoff joining scale, below which jets are joined. Is given
211 in quadratic dimensionless quantities. Either <code>yScale</code>
212 or <code>pTscale</code> must be set nonvanishing, and the larger of
213 the two dictates the actual value.
214 </argument>
215 <argument name="pTscale">
216 is the cutoff joining scale, below which jets are joined. Is given
217 in linear quantities, such as <ei>pT</ei> or <ei>m</ei> depending on
218 the measure used, but always in units of GeV. Either <code>yScale</code>
219 or <code>pTscale</code> must be set nonvanishing, and the larger of
220 the two dictates the actual value.
221 </argument>
222 <argument name="nJetMin" default="1">
223 the minimum number of jets to be reconstructed. If used, it can override
224 the <code>yScale</code> and <code>pTscale</code> values.
225 </argument>
226 <argument name="nJetMax" default="0">
227 the maximum number of jets to be reconstructed. Is not used if below
228 <code>nJetMin</code>. If used, it can override the <code>yScale</code>
229 and <code>pTscale</code> values. Thus e.g.
230 <code>nJetMin = nJetMax = 3</code> can be used to reconstruct exactly
231 3 jets.
232 </argument>
233 <br/>If the routine returns <code>false</code> the analysis failed,
234 e.g. because the number of particles was smaller than the minimum number
235 of jets requested.
236 </method>
238 <p/>
239 After the analysis has been performed, a few <code>ClusterJet</code>
240 class methods are available to return the result of the analysis:
241 <method name="size()">
242 gives the number of jets found, with jets numbered 0 through
243 <code>size() - 1</code>,
244 </method>
245 <method name="p(i)">
246 gives a <code>Vec4</code> corresponding to the four-momentum defined by
247 the sum of all the contributing particles to the <ei>i</ei>'th jet,
248 </method>
249 <method name="jetAssignment(i)">
250 gives the index of the jet that the particle <ei>i</ei> of the event
251 record belongs to,
252 </method>
253 <method name="list()">
254 provides a listing of the reconstructed jets.
255 </method>
256 <method name="nError()">
257 tells the number of times <code>analyze</code> failed to analyze events.
258 </method>
260 <h3>CellJet</h3>
262 <code>CellJet</code> (a.k.a. <code>PYCELL</code>) is a simple cone jet
263 finder in the UA1 spirit, see the PYTHIA 6 manual. It works in an
264 <ei>(eta, phi, eT)</ei> space, where <ei>eta</ei> is pseudorapidity,
265 <ei>phi</ei> azimuthal angle and <ei>eT</ei> transverse energy.
266 It will draw cones in <ei>R = sqrt(Delta-eta^2 + Delta-phi^2)</ei>
267 around seed cells. If the total <ei>eT</ei> inside the cone exceeds
268 the threshold, a jet is formed, and the cells are removed from further
269 analysis. There are no split or merge procedures, so later-found jets
270 may be missing some of the edge regions already used up by previous
271 ones.
273 <p/>
274 A cell-jet analysis object is declared by
275 <class name="CellJet cellJet( etaMax, nEta, nPhi, select, smear,
276 resolution, upperCut, threshold)">
277 where
278 <argument name="etaMax" default="5.">
279 the maximum +-pseudorapidity that the detector is assumed to cover.
280 </argument>
281 <argument name="nEta" default="50">
282 the number of equal-sized bins that the <ei>+-etaMax</ei> range
283 is assumed to be divided into.
284 </argument>
285 <argument name="nPhi" default="32">
286 the number of equal-sized bins that the <ei>phi</ei> range
287 <ei>+-pi</ei> is assumed to be divided into.
288 </argument>
289 <argument name="select" default="2">
290 tells which particles are analyzed,
291 <argoption value="1">all final-state particles,</argoption>
292 <argoption value="2">all observable final-state particles,
293 i.e. excluding neutrinos and other particles without strong or
294 electromagnetic interactions (the <code>isVisible()</code> particle
295 method),
296 and</argoption>
297 <argoption value="3">only charged final-state particles.</argoption>
298 </argument>
299 <argument name="smear" default="0">
300 strategy to smear the actual <ei>eT</ei> bin by bin,
301 <argoption value="0">no smearing,</argoption>
302 <argoption value="1">smear the <ei>eT</ei> according to a Gaussian
303 with width <ei>resolution * sqrt(eT)</ei>, with the Gaussian truncated
304 at 0 and <ei>upperCut * eT</ei>,</argoption>
305 <argoption value="2">smear the <ei>e = eT * cosh(eta)</ei> according
306 to a Gaussian with width <ei>resolution * sqrt(e)</ei>, with the
307 Gaussian truncated at 0 and <ei>upperCut * e</ei>.</argoption>
308 </argument>
309 <argument name="resolution" default="0.5">
310 see above
311 </argument>
312 <argument name="upperCut" default="2.">
313 see above
314 </argument>
315 <argument name="threshold" default="0 GeV">
316 completely neglect all bins with an <ei>eT &lt; threshold</ei>.
317 </argument>
318 </class>
320 <p/>
321 The analysis is performed by a
322 <method name="analyze( event, eTjetMin, coneRadius, eTseed)">
323 where
324 <argument name="event">is an object of the <code>Event</code> class,
325 most likely the <code>pythia.event</code> one.
326 </argument>
327 <argument name="eTjetMin" default="20. GeV">
328 is the minimum transverse energy inside a cone for this to be
329 accepted as a jet.
330 </argument>
331 <argument name="coneRadius" default="0.7">
332  is the size of the cone in <ei>(eta, phi)</ei> space drawn around
333 the geometric center of the jet.
334 </argument>
335 <argument name="eTseed" default="1.5 GeV">
336 the mimimum <ei>eT</ei> in a cell for this to be acceptable as
337 the trial center of a jet.
338 </argument>
339 <br/>If the routine returns <code>false</code> the analysis failed,
340 but currently this is not foreseen ever to happen.
341 </method>
343 <p/>
344 After the analysis has been performed, a few <code>CellJet</code>
345 class methods are available to return the result of the analysis:
346 <method name="size()">
347 gives the number of jets found, with jets numbered 0 through
348 <code>size() - 1</code>,
349 </method>
350 <method name="eT(i)">
351 gives the <ei>eT</ei> of the <ei>i</ei>'th jet, where jets have been
352 ordered with decreasing <ei>eT</ei> values,
353 </method>
354 <method name="etaCenter(i), phiCenter(i)">
355 gives the <ei>eta</ei> and <ei>phi</ei> coordinates of the geometrical
356 center of the <ei>i</ei>'th jet,
357 </method>
358 <method name="etaWeighted(i), phiWeighted(i)">
359 gives the <ei>eta</ei> and <ei>phi</ei> coordinates of the
360 <ei>eT</ei>-weighted center of the <ei>i</ei>'th jet,
361 </method>
362 <method name="multiplicity(i)">
363 gives the number of particles clustered into the <ei>i</ei>'th jet,
364 </method>
365 <method name="pMassless(i)">
366 gives a Vec4 corresponding to the four-momentum defined by the
367 <ei>eT</ei> and the weighted center of the <ei>i</ei>'th jet,
368 </method>
369 <method name="pMassive(i)">
370 gives a <code>Vec4</code> corresponding to the four-momentum defined by
371 the sum of all the contributing cells to the <ei>i</ei>'th jet, where
372 each cell contributes a four-momentum as if all the <ei>eT</ei> is
373 deposited in the center of the cell,
374 </method>
375 <method name="m(i)">
376 gives the invariant mass of the <ei>i</ei>'th jet, defined by the
377 <code>pMassive</code> above,
378 </method>
379 <method name="list()">
380 provides a listing of the above information (except <code>pMassless</code>,
381 for reasons of space).
382 </method>
383 <method name="nError()">
384 tells the number of times <code>analyze</code> failed to analyze events.
385 </method>
387 </chapter>
389 <!-- Copyright (C) 2008 Torbjorn Sjostrand -->