[1] | 1 | <html> |
---|
| 2 | <head> |
---|
| 3 | <title>Event Analysis</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 | <script language=javascript type=text/javascript> |
---|
| 10 | function stopRKey(evt) { |
---|
| 11 | var evt = (evt) ? evt : ((event) ? event : null); |
---|
| 12 | var node = (evt.target) ? evt.target :((evt.srcElement) ? evt.srcElement : null); |
---|
| 13 | if ((evt.keyCode == 13) && (node.type=="text")) |
---|
| 14 | {return false;} |
---|
| 15 | } |
---|
| 16 | |
---|
| 17 | document.onkeypress = stopRKey; |
---|
| 18 | </script> |
---|
| 19 | <?php |
---|
| 20 | if($_POST['saved'] == 1) { |
---|
| 21 | if($_POST['filepath'] != "files/") { |
---|
| 22 | echo "<font color='red'>SETTINGS SAVED TO FILE</font><br/><br/>"; } |
---|
| 23 | else { |
---|
| 24 | echo "<font color='red'>NO FILE SELECTED YET.. PLEASE DO SO </font><a href='SaveSettings.php'>HERE</a><br/><br/>"; } |
---|
| 25 | } |
---|
| 26 | ?> |
---|
| 27 | |
---|
| 28 | <form method='post' action='EventAnalysis.php'> |
---|
| 29 | |
---|
| 30 | <h2>Event Analysis</h2> |
---|
| 31 | |
---|
| 32 | <h3>Introduction</h3> |
---|
| 33 | |
---|
| 34 | The routines in this section are intended to be used to analyze |
---|
| 35 | event properties. As such they are not part of the main event |
---|
| 36 | generation chain, but can be used in comparisons between Monte |
---|
| 37 | Carlo events and real data. They are rather free-standing, but |
---|
| 38 | assume that input is provided in the PYTHIA 8 |
---|
| 39 | <code>Event</code> format, and use a few basic facilities such |
---|
| 40 | as four-vectors. Their ordering is mainly by history; for current |
---|
| 41 | LHC applications the final one, <code>SlowJet</code>, is of |
---|
| 42 | special interest. |
---|
| 43 | |
---|
| 44 | <p/> |
---|
| 45 | In addition to the methods presented here, there is also the |
---|
| 46 | possibility to make use of <?php $filepath = $_GET["filepath"]; |
---|
| 47 | echo "<a href='JetFinders.php?filepath=".$filepath."' target='page'>";?>external |
---|
| 48 | jet finders </a>. |
---|
| 49 | |
---|
| 50 | <h3>Sphericity</h3> |
---|
| 51 | |
---|
| 52 | The standard sphericity tensor is |
---|
| 53 | <br/><i> |
---|
| 54 | S^{ab} = (sum_i p_i^a p_i^b) / (sum_i p_i^2) |
---|
| 55 | </i><br/> |
---|
| 56 | where the <i>sum i</i> runs over the particles in the event, |
---|
| 57 | <i>a, b = x, y, z,</i> and <i>p</i> without such an index is |
---|
| 58 | the absolute size of the three-momentum . This tensor can be |
---|
| 59 | diagonalized to find eigenvalues and eigenvectors. |
---|
| 60 | |
---|
| 61 | <p/> |
---|
| 62 | The above tensor can be generalized by introducing a power |
---|
| 63 | <i>r</i>, such that |
---|
| 64 | <br/><i> |
---|
| 65 | S^{ab} = (sum_i p_i^a p_i^b p_i^{r-2}) / (sum_i p_i^r) |
---|
| 66 | </i><br/> |
---|
| 67 | In particular, <i>r = 1</i> gives a linear dependence on momenta |
---|
| 68 | and thus a collinear safe definition, unlike sphericity. |
---|
| 69 | |
---|
| 70 | <p/> |
---|
| 71 | To do sphericity analyses you have to set up a <code>Sphericity</code> |
---|
| 72 | instance, and then feed in events to it, one at a time. The results |
---|
| 73 | for the latest event are available as output from a few methods. |
---|
| 74 | |
---|
| 75 | <a name="method1"></a> |
---|
| 76 | <p/><strong>Sphericity::Sphericity(double power = 2., int select = 2) </strong> <br/> |
---|
| 77 | create a sphericity analysis object, where |
---|
| 78 | <br/><code>argument</code><strong> power </strong> (<code>default = <strong>2.</strong></code>) : |
---|
| 79 | is the power <i>r</i> defined above, i.e. |
---|
| 80 | <br/><code>argumentoption </code><strong> 2.</strong> : gives Spericity, and |
---|
| 81 | <br/><code>argumentoption </code><strong> 1.</strong> : gives the linear form. |
---|
| 82 | |
---|
| 83 | <br/><code>argument</code><strong> select </strong> (<code>default = <strong>2</strong></code>) : |
---|
| 84 | tells which particles are analyzed, |
---|
| 85 | <br/><code>argumentoption </code><strong> 1</strong> : all final-state particles, |
---|
| 86 | <br/><code>argumentoption </code><strong> 2</strong> : all observable final-state particles, |
---|
| 87 | i.e. excluding neutrinos and other particles without strong or |
---|
| 88 | electromagnetic interactions (the <code>isVisible()</code> |
---|
| 89 | particle method), and |
---|
| 90 | |
---|
| 91 | <br/><code>argumentoption </code><strong> 3</strong> : only charged final-state particles. |
---|
| 92 | |
---|
| 93 | |
---|
| 94 | |
---|
| 95 | <a name="method2"></a> |
---|
| 96 | <p/><strong>bool Sphericity::analyze( const Event& event, ostream& os = cout) </strong> <br/> |
---|
| 97 | perform a sphericity analysis, where |
---|
| 98 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 99 | most likely the <code>pythia.event</code> one. |
---|
| 100 | |
---|
| 101 | <br/><code>argument</code><strong> os </strong> (<code>default = <strong>cout</strong></code>) : is the output stream for |
---|
| 102 | error messages. (The method does not rely on the <code>Info</code> |
---|
| 103 | mchinery for error messages.) |
---|
| 104 | |
---|
| 105 | <br/>If the routine returns <code>false</code> the |
---|
| 106 | analysis failed, e.g. if too few particles are present to analyze. |
---|
| 107 | |
---|
| 108 | |
---|
| 109 | <p/> |
---|
| 110 | After the analysis has been performed, a few methods are available |
---|
| 111 | to return the result of the analysis of the latest event: |
---|
| 112 | |
---|
| 113 | <a name="method3"></a> |
---|
| 114 | <p/><strong>double Sphericity::sphericity() </strong> <br/> |
---|
| 115 | gives the sphericity (or equivalent if <i>r</i> is not 2), |
---|
| 116 | |
---|
| 117 | |
---|
| 118 | <a name="method4"></a> |
---|
| 119 | <p/><strong>double Sphericity::aplanarity() </strong> <br/> |
---|
| 120 | gives the aplanarity (with the same comment), |
---|
| 121 | |
---|
| 122 | |
---|
| 123 | <a name="method5"></a> |
---|
| 124 | <p/><strong>double Sphericity::eigenValue(int i) </strong> <br/> |
---|
| 125 | gives one of the three eigenvalues for <i>i</i> = 1, 2 or 3, in |
---|
| 126 | descending order, |
---|
| 127 | |
---|
| 128 | |
---|
| 129 | <a name="method6"></a> |
---|
| 130 | <p/><strong>Vec4 Sphericity::eventAxis(i) </strong> <br/> |
---|
| 131 | gives the matching normalized eigenvector, as a <code>Vec4</code> |
---|
| 132 | with vanishing time/energy component. |
---|
| 133 | |
---|
| 134 | |
---|
| 135 | <a name="method7"></a> |
---|
| 136 | <p/><strong>void Sphericity::list(ostream& os = cout) </strong> <br/> |
---|
| 137 | provides a listing of the above information. |
---|
| 138 | |
---|
| 139 | |
---|
| 140 | <p/> |
---|
| 141 | There is also one method that returns information accumulated for all |
---|
| 142 | the events analyzed so far. |
---|
| 143 | |
---|
| 144 | <a name="method8"></a> |
---|
| 145 | <p/><strong>int Sphericity::nError() </strong> <br/> |
---|
| 146 | tells the number of times <code>analyze(...)</code> failed to analyze |
---|
| 147 | events, i.e. returned <code>false</code>. |
---|
| 148 | |
---|
| 149 | |
---|
| 150 | <h3>Thrust</h3> |
---|
| 151 | |
---|
| 152 | Thrust is obtained by varying the thrust axis so that the longitudinal |
---|
| 153 | momentum component projected onto it is maximized, and thrust itself is |
---|
| 154 | then defined as the sum of absolute longitudinal momenta divided by |
---|
| 155 | the sum of absolute momenta. The major axis is found correspondingly |
---|
| 156 | in the plane transverse to thrust, and the minor one is then defined |
---|
| 157 | to be transverse to both. Oblateness is the difference between the major |
---|
| 158 | and the minor values. |
---|
| 159 | |
---|
| 160 | <p/> |
---|
| 161 | The calculation of thrust is more computer-time-intensive than e.g. |
---|
| 162 | linear sphericity, introduced above, and has no specific advantages except |
---|
| 163 | historical precedent. In the PYTHIA 6 implementation the search was |
---|
| 164 | speeded up at the price of then not being guaranteed to hit the absolute |
---|
| 165 | maximum. The current implementation studies all possibilities, but at |
---|
| 166 | the price of being slower, with time consumption for an event with |
---|
| 167 | <i>n</i> particles growing like <i>n^3</i>. |
---|
| 168 | |
---|
| 169 | <p/> |
---|
| 170 | To do thrust analyses you have to set up a <code>Thrust</code> |
---|
| 171 | instance, and then feed in events to it, one at a time. The results |
---|
| 172 | for the latest event are available as output from a few methods. |
---|
| 173 | |
---|
| 174 | <a name="method9"></a> |
---|
| 175 | <p/><strong>Thrust::Thrust(int select = 2) </strong> <br/> |
---|
| 176 | create a thrust analysis object, where |
---|
| 177 | <br/><code>argument</code><strong> select </strong> (<code>default = <strong>2</strong></code>) : |
---|
| 178 | tells which particles are analyzed, |
---|
| 179 | <br/><code>argumentoption </code><strong> 1</strong> : all final-state particles, |
---|
| 180 | <br/><code>argumentoption </code><strong> 2</strong> : all observable final-state particles, |
---|
| 181 | i.e. excluding neutrinos and other particles without strong or |
---|
| 182 | electromagnetic interactions (the <code>isVisible()</code> |
---|
| 183 | particle method), and |
---|
| 184 | |
---|
| 185 | <br/><code>argumentoption </code><strong> 3</strong> : only charged final-state particles. |
---|
| 186 | |
---|
| 187 | |
---|
| 188 | |
---|
| 189 | <a name="method10"></a> |
---|
| 190 | <p/><strong>bool Thrust::analyze( const Event& event, ostream& os = cout) </strong> <br/> |
---|
| 191 | perform a thrust analysis, where |
---|
| 192 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 193 | most likely the <code>pythia.event</code> one. |
---|
| 194 | |
---|
| 195 | <br/><code>argument</code><strong> os </strong> (<code>default = <strong>cout</strong></code>) : is the output stream for |
---|
| 196 | error messages. (The method does not rely on the <code>Info</code> |
---|
| 197 | mchinery for error messages.) |
---|
| 198 | |
---|
| 199 | <br/>If the routine returns <code>false</code> the |
---|
| 200 | analysis failed, e.g. if too few particles are present to analyze. |
---|
| 201 | |
---|
| 202 | |
---|
| 203 | <p/> |
---|
| 204 | After the analysis has been performed, a few methods are available |
---|
| 205 | to return the result of the analysis of the latest event: |
---|
| 206 | |
---|
| 207 | <a name="method11"></a> |
---|
| 208 | <p/><strong>double Thrust::thrust() </strong> <br/> |
---|
| 209 | |
---|
| 210 | <strong>double Thrust::tMajor() </strong> <br/> |
---|
| 211 | |
---|
| 212 | <strong>double Thrust::tMinor() </strong> <br/> |
---|
| 213 | |
---|
| 214 | <strong>double Thrust::oblateness() </strong> <br/> |
---|
| 215 | gives the thrust, major, minor and oblateness values, respectively, |
---|
| 216 | |
---|
| 217 | |
---|
| 218 | <a name="method12"></a> |
---|
| 219 | <p/><strong>Vec4 Thrust::eventAxis(int i) </strong> <br/> |
---|
| 220 | gives the matching normalized event-axis vectors, for <i>i</i> = 1, 2 or 3 |
---|
| 221 | corresponding to thrust, major or minor, as a <code>Vec4</code> with |
---|
| 222 | vanishing time/energy component. |
---|
| 223 | |
---|
| 224 | |
---|
| 225 | <a name="method13"></a> |
---|
| 226 | <p/><strong>void Thrust::list(ostream& os = cout) </strong> <br/> |
---|
| 227 | provides a listing of the above information. |
---|
| 228 | |
---|
| 229 | |
---|
| 230 | <p/> |
---|
| 231 | There is also one method that returns information accumulated for all |
---|
| 232 | the events analyzed so far. |
---|
| 233 | |
---|
| 234 | <a name="method14"></a> |
---|
| 235 | <p/><strong>int Thrust::nError() </strong> <br/> |
---|
| 236 | tells the number of times <code>analyze(...)</code> failed to analyze |
---|
| 237 | events, i.e. returned <code>false</code>. |
---|
| 238 | |
---|
| 239 | |
---|
| 240 | <h3>ClusterJet</h3> |
---|
| 241 | |
---|
| 242 | <code>ClusterJet</code> (a.k.a. <code>LUCLUS</code> and |
---|
| 243 | <code>PYCLUS</code>) is a clustering algorithm of the type used for |
---|
| 244 | analyses of <i>e^+e^-</i> events, see the PYTHIA 6 manual. All |
---|
| 245 | visible particles in the events are clustered into jets. A few options |
---|
| 246 | are available for some well-known distance measures. Cutoff |
---|
| 247 | distances can either be given in terms of a scaled quadratic quantity |
---|
| 248 | like <i>y = pT^2/E^2</i> or an unscaled linear one like <i>pT</i>. |
---|
| 249 | |
---|
| 250 | <p/> |
---|
| 251 | Note that we have deliberately chosen not to include the <i>e^+e^-</i> |
---|
| 252 | equivalents of the Cambridge/Aachen and anti-<i>kRT</i> algorithms. |
---|
| 253 | These tend to be good at clustering the densely populated (in angle) |
---|
| 254 | cores of jets, but less successful for the sparsely populated transverse |
---|
| 255 | regions, where many jets may come to consist of a single low-momentum |
---|
| 256 | particle. In hadron collisions such jets could easily be disregarded, |
---|
| 257 | while in <i>e^+e^-</i> annihilation all particles derive back from the |
---|
| 258 | hard process. |
---|
| 259 | |
---|
| 260 | <p/> |
---|
| 261 | To do jet finding analyses you have to set up a <code>ClusterJet</code> |
---|
| 262 | instance, and then feed in events to it, one at a time. The results |
---|
| 263 | for the latest event are available as output from a few methods. |
---|
| 264 | |
---|
| 265 | <a name="method15"></a> |
---|
| 266 | <p/><strong>ClusterJet::ClusterJet(string measure = "Lund", int select = 2, int massSet = 2, bool precluster = false, bool reassign = false) </strong> <br/> |
---|
| 267 | create a <code>ClusterJet</code> instance, where |
---|
| 268 | <br/><code>argument</code><strong> measure </strong> (<code>default = <strong>"Lund"</strong></code>) : distance measure, |
---|
| 269 | to be provided as a character string (actually, only the first character |
---|
| 270 | is necessary) |
---|
| 271 | <br/><code>argumentoption </code><strong> "Lund"</strong> : the Lund <i>pT</i> distance, |
---|
| 272 | |
---|
| 273 | <br/><code>argumentoption </code><strong> "JADE"</strong> : the JADE mass distance, and |
---|
| 274 | |
---|
| 275 | <br/><code>argumentoption </code><strong> "Durham"</strong> : the Durham <i>kT</i> measure. |
---|
| 276 | |
---|
| 277 | |
---|
| 278 | <br/><code>argument</code><strong> select </strong> (<code>default = <strong>2</strong></code>) : |
---|
| 279 | tells which particles are analyzed, |
---|
| 280 | <br/><code>argumentoption </code><strong> 1</strong> : all final-state particles, |
---|
| 281 | <br/><code>argumentoption </code><strong> 2</strong> : all observable final-state particles, |
---|
| 282 | i.e. excluding neutrinos and other particles without strong or |
---|
| 283 | electromagnetic interactions (the <code>isVisible()</code> particle |
---|
| 284 | method), and |
---|
| 285 | |
---|
| 286 | <br/><code>argumentoption </code><strong> 3</strong> : only charged final-state particles. |
---|
| 287 | |
---|
| 288 | <br/><code>argument</code><strong> massSet </strong> (<code>default = <strong>2</strong></code>) : masses assumed for the particles |
---|
| 289 | used in the analysis |
---|
| 290 | <br/><code>argumentoption </code><strong> 0</strong> : all massless, |
---|
| 291 | <br/><code>argumentoption </code><strong> 1</strong> : photons are massless while all others are |
---|
| 292 | assigned the <i>pi+-</i> mass, and |
---|
| 293 | |
---|
| 294 | <br/><code>argumentoption </code><strong> 2</strong> : all given their correct masses. |
---|
| 295 | |
---|
| 296 | <br/><code>argument</code><strong> precluster </strong> (<code>default = <strong>false</strong></code>) : |
---|
| 297 | perform or not a early preclustering step, where nearby particles |
---|
| 298 | are lumped together so as to speed up the subsequent normal clustering. |
---|
| 299 | |
---|
| 300 | <br/><code>argument</code><strong> reassign </strong> (<code>default = <strong>false</strong></code>) : |
---|
| 301 | reassign all particles to the nearest jet each time after two jets |
---|
| 302 | have been joined. |
---|
| 303 | |
---|
| 304 | |
---|
| 305 | |
---|
| 306 | <a name="method16"></a> |
---|
| 307 | <p/><strong>ClusterJet::analyze( const Event& event, double yScale, double pTscale, int nJetMin = 1, int nJetMax = 0, ostream& os = cout) </strong> <br/> |
---|
| 308 | performs a jet finding analysis, where |
---|
| 309 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 310 | most likely the <code>pythia.event</code> one. |
---|
| 311 | |
---|
| 312 | <br/><code>argument</code><strong> yScale </strong> : |
---|
| 313 | is the cutoff joining scale, below which jets are joined. Is given |
---|
| 314 | in quadratic dimensionless quantities. Either <code>yScale</code> |
---|
| 315 | or <code>pTscale</code> must be set nonvanishing, and the larger of |
---|
| 316 | the two dictates the actual value. |
---|
| 317 | |
---|
| 318 | <br/><code>argument</code><strong> pTscale </strong> : |
---|
| 319 | is the cutoff joining scale, below which jets are joined. Is given |
---|
| 320 | in linear quantities, such as <i>pT</i> or <i>m</i> depending on |
---|
| 321 | the measure used, but always in units of GeV. Either <code>yScale</code> |
---|
| 322 | or <code>pTscale</code> must be set nonvanishing, and the larger of |
---|
| 323 | the two dictates the actual value. |
---|
| 324 | |
---|
| 325 | <br/><code>argument</code><strong> nJetMin </strong> (<code>default = <strong>1</strong></code>) : |
---|
| 326 | the minimum number of jets to be reconstructed. If used, it can override |
---|
| 327 | the <code>yScale</code> and <code>pTscale</code> values. |
---|
| 328 | |
---|
| 329 | <br/><code>argument</code><strong> nJetMax </strong> (<code>default = <strong>0</strong></code>) : |
---|
| 330 | the maximum number of jets to be reconstructed. Is not used if below |
---|
| 331 | <code>nJetMin</code>. If used, it can override the <code>yScale</code> |
---|
| 332 | and <code>pTscale</code> values. Thus e.g. |
---|
| 333 | <code>nJetMin = nJetMax = 3</code> can be used to reconstruct exactly |
---|
| 334 | 3 jets. |
---|
| 335 | |
---|
| 336 | <br/><code>argument</code><strong> os </strong> (<code>default = <strong>cout</strong></code>) : is the output stream for |
---|
| 337 | error messages. (The method does not rely on the <code>Info</code> |
---|
| 338 | mchinery for error messages.) |
---|
| 339 | |
---|
| 340 | <br/>If the routine returns <code>false</code> the analysis failed, |
---|
| 341 | e.g. because the number of particles was smaller than the minimum number |
---|
| 342 | of jets requested. |
---|
| 343 | |
---|
| 344 | |
---|
| 345 | <p/> |
---|
| 346 | After the analysis has been performed, a few <code>ClusterJet</code> |
---|
| 347 | class methods are available to return the result of the analysis: |
---|
| 348 | |
---|
| 349 | <a name="method17"></a> |
---|
| 350 | <p/><strong>int ClusterJet::size() </strong> <br/> |
---|
| 351 | gives the number of jets found, with jets numbered 0 through |
---|
| 352 | <code>size() - 1</code>. |
---|
| 353 | |
---|
| 354 | |
---|
| 355 | <a name="method18"></a> |
---|
| 356 | <p/><strong>Vec4 ClusterJet::p(int i) </strong> <br/> |
---|
| 357 | gives a <code>Vec4</code> corresponding to the four-momentum defined by |
---|
| 358 | the sum of all the contributing particles to the <i>i</i>'th jet. |
---|
| 359 | |
---|
| 360 | |
---|
| 361 | <a name="method19"></a> |
---|
| 362 | <p/><strong>int ClusterJet::mult(int i) </strong> <br/> |
---|
| 363 | the number of particles that have been clustered into the <i>i</i>'th jet. |
---|
| 364 | |
---|
| 365 | |
---|
| 366 | <a name="method20"></a> |
---|
| 367 | <p/><strong>int ClusterJet::jetAssignment(int i) </strong> <br/> |
---|
| 368 | gives the index of the jet that the particle <i>i</i> of the event |
---|
| 369 | record belongs to, |
---|
| 370 | |
---|
| 371 | |
---|
| 372 | <a name="method21"></a> |
---|
| 373 | <p/><strong>void ClusterJet::list(ostream& os = cout) </strong> <br/> |
---|
| 374 | provides a listing of the reconstructed jets. |
---|
| 375 | |
---|
| 376 | |
---|
| 377 | <a name="method22"></a> |
---|
| 378 | <p/><strong>int ClusterJet::distanceSize() </strong> <br/> |
---|
| 379 | the number of most recent clustering scales that have been stored |
---|
| 380 | for readout with the next method. Normally this would be five, |
---|
| 381 | but less if fewer clustering steps occured. |
---|
| 382 | |
---|
| 383 | |
---|
| 384 | <a name="method23"></a> |
---|
| 385 | <p/><strong>double ClusterJet::distance(int i) </strong> <br/> |
---|
| 386 | clustering scales, with <code>distance(0)</code> being the most |
---|
| 387 | recent one, i.e. normally the highest, up to <code>distance(4)</code> |
---|
| 388 | being the fifth most recent. That is, with <i>n</i> being the final |
---|
| 389 | number of jets, <code>ClusterJet::size()</code>, the scales at which |
---|
| 390 | <i>n+1</i> jets become <i>n</i>, <i>n+2</i> become <i>n+1</i>, |
---|
| 391 | and so on till <i>n+5</i> become <i>n+4</i>. Nonexisting clustering |
---|
| 392 | scales are returned as zero. The physical interpretation of a scale is |
---|
| 393 | as provided by the respective distance measure (Lund, JADE, Durham). |
---|
| 394 | |
---|
| 395 | |
---|
| 396 | <p/> |
---|
| 397 | There is also one method that returns information accumulated for all |
---|
| 398 | the events analyzed so far. |
---|
| 399 | |
---|
| 400 | <a name="method24"></a> |
---|
| 401 | <p/><strong>int ClusterJet::nError() </strong> <br/> |
---|
| 402 | tells the number of times <code>analyze(...)</code> failed to analyze |
---|
| 403 | events, i.e. returned <code>false</code>. |
---|
| 404 | |
---|
| 405 | |
---|
| 406 | <h3>CellJet</h3> |
---|
| 407 | |
---|
| 408 | <code>CellJet</code> (a.k.a. <code>PYCELL</code>) is a simple cone jet |
---|
| 409 | finder in the UA1 spirit, see the PYTHIA 6 manual. It works in an |
---|
| 410 | <i>(eta, phi, eT)</i> space, where <i>eta</i> is pseudorapidity, |
---|
| 411 | <i>phi</i> azimuthal angle and <i>eT</i> transverse energy. |
---|
| 412 | It will draw cones in <i>R = sqrt(Delta-eta^2 + Delta-phi^2)</i> |
---|
| 413 | around seed cells. If the total <i>eT</i> inside the cone exceeds |
---|
| 414 | the threshold, a jet is formed, and the cells are removed from further |
---|
| 415 | analysis. There are no split or merge procedures, so later-found jets |
---|
| 416 | may be missing some of the edge regions already used up by previous |
---|
| 417 | ones. Not all particles in the event are assigned to jets; leftovers |
---|
| 418 | may be viewed as belonging to beam remnants or the underlying event. |
---|
| 419 | It is not used by any experimental collaboration, but is closely |
---|
| 420 | related to the more recent and better theoretically motivated |
---|
| 421 | anti-<i>kT</i> algorithm [<a href="Bibliography.php" target="page">Cac08</a>]. |
---|
| 422 | |
---|
| 423 | <p/> |
---|
| 424 | To do jet finding analyses you have to set up a <code>CellJet</code> |
---|
| 425 | instance, and then feed in events to it, one at a time. Especially note |
---|
| 426 | that, if you want to use the options where energies are smeared in |
---|
| 427 | order so emulate detector imperfections, you must hand in an external |
---|
| 428 | random number generator, preferably the one residing in the |
---|
| 429 | <code>Pythia</code> class. The results for the latest event are |
---|
| 430 | available as output from a few methods. |
---|
| 431 | |
---|
| 432 | <a name="method25"></a> |
---|
| 433 | <p/><strong>CellJet::CellJet(double etaMax = 5., int nEta = 50, int nPhi = 32, int select = 2, int smear = 0, double resolution = 0.5, double upperCut = 2., double threshold = 0., Rndm* rndmPtr = 0) </strong> <br/> |
---|
| 434 | create a <code>CellJet</code> instance, where |
---|
| 435 | <br/><code>argument</code><strong> etaMax </strong> (<code>default = <strong>5.</strong></code>) : |
---|
| 436 | the maximum +-pseudorapidity that the detector is assumed to cover. |
---|
| 437 | |
---|
| 438 | <br/><code>argument</code><strong> nEta </strong> (<code>default = <strong>50</strong></code>) : |
---|
| 439 | the number of equal-sized bins that the <i>+-etaMax</i> range |
---|
| 440 | is assumed to be divided into. |
---|
| 441 | |
---|
| 442 | <br/><code>argument</code><strong> nPhi </strong> (<code>default = <strong>32</strong></code>) : |
---|
| 443 | the number of equal-sized bins that the <i>phi</i> range |
---|
| 444 | <i>+-pi</i> is assumed to be divided into. |
---|
| 445 | |
---|
| 446 | <br/><code>argument</code><strong> select </strong> (<code>default = <strong>2</strong></code>) : |
---|
| 447 | tells which particles are analyzed, |
---|
| 448 | <br/><code>argumentoption </code><strong> 1</strong> : all final-state particles, |
---|
| 449 | <br/><code>argumentoption </code><strong> 2</strong> : all observable final-state particles, |
---|
| 450 | i.e. excluding neutrinos and other particles without strong or |
---|
| 451 | electromagnetic interactions (the <code>isVisible()</code> particle |
---|
| 452 | method), |
---|
| 453 | and |
---|
| 454 | <br/><code>argumentoption </code><strong> 3</strong> : only charged final-state particles. |
---|
| 455 | |
---|
| 456 | <br/><code>argument</code><strong> smear </strong> (<code>default = <strong>0</strong></code>) : |
---|
| 457 | strategy to smear the actual <i>eT</i> bin by bin, |
---|
| 458 | <br/><code>argumentoption </code><strong> 0</strong> : no smearing, |
---|
| 459 | <br/><code>argumentoption </code><strong> 1</strong> : smear the <i>eT</i> according to a Gaussian |
---|
| 460 | with width <i>resolution * sqrt(eT)</i>, with the Gaussian truncated |
---|
| 461 | at 0 and <i>upperCut * eT</i>, |
---|
| 462 | <br/><code>argumentoption </code><strong> 2</strong> : smear the <i>e = eT * cosh(eta)</i> according |
---|
| 463 | to a Gaussian with width <i>resolution * sqrt(e)</i>, with the |
---|
| 464 | Gaussian truncated at 0 and <i>upperCut * e</i>. |
---|
| 465 | |
---|
| 466 | <br/><code>argument</code><strong> resolution </strong> (<code>default = <strong>0.5</strong></code>) : |
---|
| 467 | see above. |
---|
| 468 | |
---|
| 469 | <br/><code>argument</code><strong> upperCut </strong> (<code>default = <strong>2.</strong></code>) : |
---|
| 470 | see above. |
---|
| 471 | |
---|
| 472 | <br/><code>argument</code><strong> threshold </strong> (<code>default = <strong>0 GeV</strong></code>) : |
---|
| 473 | completely neglect all bins with an <i>eT < threshold</i>. |
---|
| 474 | |
---|
| 475 | <br/><code>argument</code><strong> rndmPtr </strong> (<code>default = <strong>0</strong></code>) : |
---|
| 476 | the random-number generator used to select the smearing described |
---|
| 477 | above. Must be handed in for smearing to be possible. If your |
---|
| 478 | <code>Pythia</code> class instance is named <code>pythia</code>, |
---|
| 479 | then <code>&pythia.rndm</code> would be the logical choice. |
---|
| 480 | |
---|
| 481 | |
---|
| 482 | |
---|
| 483 | <a name="method26"></a> |
---|
| 484 | <p/><strong>bool CellJet::analyze( const Event& event, double eTjetMin = 20., double coneRadius = 0.7, double eTseed = 1.5, ostream& os = cout) </strong> <br/> |
---|
| 485 | performs a jet finding analysis, where |
---|
| 486 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 487 | most likely the <code>pythia.event</code> one. |
---|
| 488 | |
---|
| 489 | <br/><code>argument</code><strong> eTjetMin </strong> (<code>default = <strong>20. GeV</strong></code>) : |
---|
| 490 | is the minimum transverse energy inside a cone for this to be |
---|
| 491 | accepted as a jet. |
---|
| 492 | |
---|
| 493 | <br/><code>argument</code><strong> coneRadius </strong> (<code>default = <strong>0.7</strong></code>) : |
---|
| 494 | is the size of the cone in <i>(eta, phi)</i> space drawn around |
---|
| 495 | the geometric center of the jet. |
---|
| 496 | |
---|
| 497 | <br/><code>argument</code><strong> eTseed </strong> (<code>default = <strong>1.5 GeV</strong></code>) : |
---|
| 498 | the mimimum <i>eT</i> in a cell for this to be acceptable as |
---|
| 499 | the trial center of a jet. |
---|
| 500 | |
---|
| 501 | <br/><code>argument</code><strong> os </strong> (<code>default = <strong>cout</strong></code>) : is the output stream for |
---|
| 502 | error messages. (The method does not rely on the <code>Info</code> |
---|
| 503 | mchinery for error messages.) |
---|
| 504 | |
---|
| 505 | <br/>If the routine returns <code>false</code> the analysis failed, |
---|
| 506 | but currently this is not foreseen ever to happen. |
---|
| 507 | |
---|
| 508 | |
---|
| 509 | <p/> |
---|
| 510 | After the analysis has been performed, a few <code>CellJet</code> |
---|
| 511 | class methods are available to return the result of the analysis: |
---|
| 512 | |
---|
| 513 | <a name="method27"></a> |
---|
| 514 | <p/><strong>int CellJet::size() </strong> <br/> |
---|
| 515 | gives the number of jets found, with jets numbered 0 through |
---|
| 516 | <code>size() - 1</code>, |
---|
| 517 | |
---|
| 518 | |
---|
| 519 | <a name="method28"></a> |
---|
| 520 | <p/><strong>double CellJet::eT(i) </strong> <br/> |
---|
| 521 | gives the <i>eT</i> of the <i>i</i>'th jet, where jets have been |
---|
| 522 | ordered with decreasing <i>eT</i> values, |
---|
| 523 | |
---|
| 524 | |
---|
| 525 | <a name="method29"></a> |
---|
| 526 | <p/><strong>double CellJet::etaCenter(int i) </strong> <br/> |
---|
| 527 | |
---|
| 528 | <strong>double CellJet::phiCenter(int i) </strong> <br/> |
---|
| 529 | gives the <i>eta</i> and <i>phi</i> coordinates of the geometrical |
---|
| 530 | center of the <i>i</i>'th jet, |
---|
| 531 | |
---|
| 532 | |
---|
| 533 | <a name="method30"></a> |
---|
| 534 | <p/><strong>double CellJet::etaWeighted(int i) </strong> <br/> |
---|
| 535 | |
---|
| 536 | <strong>double CellJet::phiWeighted(int i) </strong> <br/> |
---|
| 537 | gives the <i>eta</i> and <i>phi</i> coordinates of the |
---|
| 538 | <i>eT</i>-weighted center of the <i>i</i>'th jet, |
---|
| 539 | |
---|
| 540 | |
---|
| 541 | <a name="method31"></a> |
---|
| 542 | <p/><strong>int CellJet::multiplicity(int i) </strong> <br/> |
---|
| 543 | gives the number of particles clustered into the <i>i</i>'th jet, |
---|
| 544 | |
---|
| 545 | |
---|
| 546 | <a name="method32"></a> |
---|
| 547 | <p/><strong>Vec4 CellJet::pMassless(int i) </strong> <br/> |
---|
| 548 | gives a <code>Vec4</code> corresponding to the four-momentum defined |
---|
| 549 | by the <i>eT</i> and the weighted center of the <i>i</i>'th jet, |
---|
| 550 | |
---|
| 551 | |
---|
| 552 | <a name="method33"></a> |
---|
| 553 | <p/><strong>Vec4 CellJet::pMassive(int i) </strong> <br/> |
---|
| 554 | gives a <code>Vec4</code> corresponding to the four-momentum defined by |
---|
| 555 | the sum of all the contributing cells to the <i>i</i>'th jet, where |
---|
| 556 | each cell contributes a four-momentum as if all the <i>eT</i> is |
---|
| 557 | deposited in the center of the cell, |
---|
| 558 | |
---|
| 559 | |
---|
| 560 | <a name="method34"></a> |
---|
| 561 | <p/><strong>double CellJet::m(int i) </strong> <br/> |
---|
| 562 | gives the invariant mass of the <i>i</i>'th jet, defined by the |
---|
| 563 | <code>pMassive</code> above, |
---|
| 564 | |
---|
| 565 | |
---|
| 566 | <a name="method35"></a> |
---|
| 567 | <p/><strong>void CellJet::list() </strong> <br/> |
---|
| 568 | provides a listing of the above information (except <code>pMassless</code>, |
---|
| 569 | for reasons of space). |
---|
| 570 | |
---|
| 571 | |
---|
| 572 | <p/> |
---|
| 573 | There is also one method that returns information accumulated for all |
---|
| 574 | the events analyzed so far. |
---|
| 575 | <a name="method36"></a> |
---|
| 576 | <p/><strong>int CellJet::nError() </strong> <br/> |
---|
| 577 | tells the number of times <code>analyze(...)</code> failed to analyze |
---|
| 578 | events, i.e. returned <code>false</code>. |
---|
| 579 | |
---|
| 580 | |
---|
| 581 | <h3>SlowJet</h3> |
---|
| 582 | |
---|
| 583 | <code>SlowJet</code> is a simple program for doing jet finding according |
---|
| 584 | to either of the <i>kT</i>, anti-<i>kT</i>, and Cambridge/Aachen |
---|
| 585 | algorithms, in a cylindrical coordinate frame. The name is obviously |
---|
| 586 | an homage to the <code>FastJet</code> program [<a href="Bibliography.php" target="page">Cac06</a>]. |
---|
| 587 | That package contains many more algorithms, with many more options, |
---|
| 588 | and, above all, is <i>much</i> faster. Therefore <code>SlowJet</code> |
---|
| 589 | is not so much intended for massive processing of data or Monte Carlo |
---|
| 590 | files as for simple first tests. Nevertheless, within the advertised |
---|
| 591 | capabilities of <code>SlowJet</code>, it has been checked to find |
---|
| 592 | identically the same jets as <code>FastJet</code>. The time consumption |
---|
| 593 | typically is around or below that to generate an LHC <i>pp</i> event |
---|
| 594 | in the first place, so is not prohibitive. But the time rises rapidly |
---|
| 595 | for large multiplicities, so obviously <code>SlowJet</code> can not |
---|
| 596 | be used for tricks like distributing a dense grid of pseudoparticles |
---|
| 597 | to be able to define jet areas, like <code>FastJet</code> can, and also |
---|
| 598 | not for events with much pileup or other noise. |
---|
| 599 | |
---|
| 600 | <p/> |
---|
| 601 | The first step is to decide which particles should be included in the |
---|
| 602 | analysis, and with what four-momenta. The <code>SlowJet</code> constructor |
---|
| 603 | allows to pick a maximum pseudorapidity defined by the extent of the |
---|
| 604 | assumed detector, to pick among some standard options of which particles |
---|
| 605 | to analyze, and to allow for some standard mass assumptions, like that |
---|
| 606 | all charged particles have the pion mass. Obviously this is only a |
---|
| 607 | restricted set of possibilities. |
---|
| 608 | |
---|
| 609 | <p/> |
---|
| 610 | Full flexibility can be obtained by deriving from the base class |
---|
| 611 | <code>SlowJetHook</code> to write your own <code>include</code> method. |
---|
| 612 | This will be presented with one final-state particle at a time, and |
---|
| 613 | should return <code>true</code> for those particles that should be |
---|
| 614 | analyzed. It is also possible to return modified four-momenta and masses, |
---|
| 615 | to take into account detector smearing effects or particle identity |
---|
| 616 | misassignments, but you must respect <i>E^2 - p^2 = m^2</i>. |
---|
| 617 | |
---|
| 618 | <p/> |
---|
| 619 | Alternatively you can modify the event record itself, or a copy of it |
---|
| 620 | (if you want to keep the original intact). For instance, only final |
---|
| 621 | particles are considered in the analysis, i.e. particles with positive |
---|
| 622 | status code, so negative status code should then be assigned to those |
---|
| 623 | particles that you do not want to see analyzed. Again four-momenta and |
---|
| 624 | masses can be modified, subject to <i>E^2 - p^2 = m^2</i>. |
---|
| 625 | |
---|
| 626 | <p/> |
---|
| 627 | The jet reconstructions is then based on sequential recombination with |
---|
| 628 | progressive removal, using the <i>E</i> recombination scheme. To be |
---|
| 629 | more specific, the algorithm works as follows. |
---|
| 630 | <ol> |
---|
| 631 | <li>Each particle to be analyzed defines an original cluster. It has a |
---|
| 632 | well-defined four-momentum and mass at input. From this information the |
---|
| 633 | triplet <i>(pT, y, phi)</i> is calculated, i.e. the transverse momentum, |
---|
| 634 | rapidity and azimuthal angle of the cluster.</li> |
---|
| 635 | <li>Define distance measures of all clusters <i>i</i> to the beam |
---|
| 636 | <br/><i>d_iB = pT_i^2p</i><br/> |
---|
| 637 | and of all pairs <i>(i,j)</i> relative to each other |
---|
| 638 | <br/><i>d_ij = min( pT_i^2p, pT_j^2p) DeltaR_ij^2 / R^2 </i><br/> |
---|
| 639 | where |
---|
| 640 | <br/><i>DeltaR_ij^2 = (y_i - y_j)^2 + (phi_i - phi_j)^2.</i><br/> |
---|
| 641 | The jet algorithm is determined by the user-selected <i>p</i> value, |
---|
| 642 | where <i>p = -1</i> corresponds to the anti-<i>kT</i> one, |
---|
| 643 | <i>p = 0</i> to the Cambridge/Aachen one and <i>p = 1</i> to the |
---|
| 644 | <i>kT</i> one. Also <i>R</i> is chosen by the user, to give an |
---|
| 645 | approximate measure of the size of jets. However, note that jets need |
---|
| 646 | not have a circular shape in <i>(y, phi)</i> space, so <i>R</i> |
---|
| 647 | can not straight off be interpreted as a jet radius.</li> |
---|
| 648 | <li>Find the smallest of all <i>d_iB</i> and <i>d_ij</i>.</li> |
---|
| 649 | <li>If this is a <i>d_iB</i> then cluster <i>i</i> is removed from |
---|
| 650 | the clusters list and instead added to the jets list. |
---|
| 651 | Optionally, a <i>pTjetMin</i> requirement is imposed, where only |
---|
| 652 | clusters with <i>pT > pTjetMin</i> are added to the jets list. |
---|
| 653 | If so, some of the analyzed particles will not be part of any final |
---|
| 654 | jet.</li> |
---|
| 655 | <li>If instead te smallest measure is a <i>d_ij</i> then the |
---|
| 656 | four-momenta of the <i>i</i> and <i>j</i> clusters are added |
---|
| 657 | to define a single new cluster. Convert this four-momentum to a new |
---|
| 658 | <i>(pT, y, phi)</i> triplet and update the list of <i>d_iB</i> |
---|
| 659 | and <i>d_ij</i>.</li> |
---|
| 660 | <li>Return to step 3 until no clusters remain.</li> |
---|
| 661 | </ol> |
---|
| 662 | |
---|
| 663 | <p/> |
---|
| 664 | To do jet finding analyses you first have to set up a <code>SlowJet</code> |
---|
| 665 | instance, where the arguments of the constructor specifies the details |
---|
| 666 | of the subsequent analyses. Thereafter you can feed in events to it, |
---|
| 667 | one at a time, and have them analyzed by the <code>analyze</code> method. |
---|
| 668 | Information on the resulting jets can be extracted by a few different methods. |
---|
| 669 | The minimal procedure only requires one call per event to do the analysis. |
---|
| 670 | We will begin by presenting it, and only afterwards some extensions. |
---|
| 671 | |
---|
| 672 | <a name="method37"></a> |
---|
| 673 | <p/><strong>SlowJet::SlowJet(double power, double R, double pTjetMin = 0.,double etaMax = 25., int select = 2, int massSet = 2, SlowJetHook* sjHookPtr = 0) </strong> <br/> |
---|
| 674 | create a <code>SlowJet</code> instance, where |
---|
| 675 | <br/><code>argument</code><strong> power </strong> : |
---|
| 676 | tells (half) the power of the transverse-momentum dependence of the |
---|
| 677 | distance measure, |
---|
| 678 | <br/><code>argumentoption </code><strong> -1</strong> : the anti-<i>kT</i> algorithm, |
---|
| 679 | <br/><code>argumentoption </code><strong> 0</strong> : the Cambridge/Aachen algorithm, and |
---|
| 680 | <br/><code>argumentoption </code><strong> 1</strong> : the <i>kT</i> algorithm. |
---|
| 681 | |
---|
| 682 | <br/><code>argument</code><strong> R </strong> : |
---|
| 683 | the <i>R</i> size parameter, which is crudely related to the radius of |
---|
| 684 | the jet cone in <i>(y, phi)</i> space around the center of the jet. |
---|
| 685 | |
---|
| 686 | <br/><code>argument</code><strong> pTjetMin </strong> (<code>default = <strong>0.0 GeV</strong></code>) : |
---|
| 687 | the minimum transverse momentum required for a cluster |
---|
| 688 | to become a jet. By default all clusters become jets, and therefore |
---|
| 689 | all analyzed particles are assigned to a jet. |
---|
| 690 | For comparisons with perturbative QCD, however, it is only meaningful |
---|
| 691 | to consider jets with a significant <i>pT</i>. |
---|
| 692 | |
---|
| 693 | <br/><code>argument</code><strong> etaMax </strong> (<code>default = <strong>25.</strong></code>) : |
---|
| 694 | the maximum +-pseudorapidity that the detector is assumed to cover. |
---|
| 695 | If you pick a value above 20 there is assumed to be full coverage |
---|
| 696 | (obviously only meaningful for theoretical studies). |
---|
| 697 | |
---|
| 698 | <br/><code>argument</code><strong> select </strong> (<code>default = <strong>2</strong></code>) : |
---|
| 699 | tells which particles are analyzed, |
---|
| 700 | <br/><code>argumentoption </code><strong> 1</strong> : all final-state particles, |
---|
| 701 | <br/><code>argumentoption </code><strong> 2</strong> : all observable final-state particles, |
---|
| 702 | i.e. excluding neutrinos and other particles without strong or |
---|
| 703 | electromagnetic interactions (the <code>isVisible()</code> particle |
---|
| 704 | method), |
---|
| 705 | and |
---|
| 706 | <br/><code>argumentoption </code><strong> 3</strong> : only charged final-state particles. |
---|
| 707 | |
---|
| 708 | <br/><code>argument</code><strong> massSet </strong> (<code>default = <strong>2</strong></code>) : masses assumed for the particles |
---|
| 709 | used in the analysis |
---|
| 710 | <br/><code>argumentoption </code><strong> 0</strong> : all massless, |
---|
| 711 | <br/><code>argumentoption </code><strong> 1</strong> : photons are massless while all others are |
---|
| 712 | assigned the <i>pi+-</i> mass, and |
---|
| 713 | |
---|
| 714 | <br/><code>argumentoption </code><strong> 2</strong> : all given their correct masses. |
---|
| 715 | |
---|
| 716 | <br/><code>argument</code><strong> sjHookPtr </strong> (<code>default = <strong>0</strong></code>) : |
---|
| 717 | gives the possibility to send in your own selection routine for which |
---|
| 718 | particles should be part of the analysis; see further below on the |
---|
| 719 | <code>SlowJetHook</code> class. If this pointer is sent in nonzero, |
---|
| 720 | <code>etaMax</code> and <code>massSet</code> are disregarded, |
---|
| 721 | and <code>select</code> only gives the basic selection, to which |
---|
| 722 | the user can add further requirements. |
---|
| 723 | |
---|
| 724 | |
---|
| 725 | |
---|
| 726 | <a name="method38"></a> |
---|
| 727 | <p/><strong>bool SlowJet::analyze( const Event& event) </strong> <br/> |
---|
| 728 | performs a jet finding analysis, where |
---|
| 729 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 730 | most likely the <code>pythia.event</code> one. |
---|
| 731 | |
---|
| 732 | <br/>If the routine returns <code>false</code> the analysis failed, |
---|
| 733 | but currently this is not foreseen ever to happen. |
---|
| 734 | |
---|
| 735 | |
---|
| 736 | <p/> |
---|
| 737 | After the analysis has been performed, a few <code>SlowJet</code> |
---|
| 738 | class methods are available to return the result of the analysis: |
---|
| 739 | |
---|
| 740 | <a name="method39"></a> |
---|
| 741 | <p/><strong>int SlowJet::sizeOrig() </strong> <br/> |
---|
| 742 | gives the original number of particles (and thus clusters) that the |
---|
| 743 | analysis starts with. |
---|
| 744 | |
---|
| 745 | |
---|
| 746 | <a name="method40"></a> |
---|
| 747 | <p/><strong>int SlowJet::sizeJet() </strong> <br/> |
---|
| 748 | gives the number of jets found, with jets numbered 0 through |
---|
| 749 | <code>sizeJet() - 1</code>, and ordered in terms of decreasing |
---|
| 750 | transverse momentum values w.r.t. the beam axis, |
---|
| 751 | |
---|
| 752 | |
---|
| 753 | <a name="method41"></a> |
---|
| 754 | <p/><strong>double SlowJet::pT(i) </strong> <br/> |
---|
| 755 | gives the transverse momentum <i>pT</i> of the <i>i</i>'th jet, |
---|
| 756 | |
---|
| 757 | |
---|
| 758 | <a name="method42"></a> |
---|
| 759 | <p/><strong>double SlowJet::y(int i) </strong> <br/> |
---|
| 760 | |
---|
| 761 | <strong>double SlowJet::phi(int i) </strong> <br/> |
---|
| 762 | gives the rapidity <i>y</i> and azimuthal angle <i>phi</i> |
---|
| 763 | of the center of the <i>i</i>'th jet (defined by the vector sum |
---|
| 764 | of constituent four-momenta), |
---|
| 765 | |
---|
| 766 | |
---|
| 767 | <a name="method43"></a> |
---|
| 768 | <p/><strong>Vec4 SlowJet::p(int i) </strong> <br/> |
---|
| 769 | |
---|
| 770 | <strong>double SlowJet::m(int i) </strong> <br/> |
---|
| 771 | gives a <code>Vec4</code> corresponding to the four-momentum |
---|
| 772 | sum of the particles assigned to the <i>i</i>'th jet, and |
---|
| 773 | the invariant mass of this four-vector, |
---|
| 774 | |
---|
| 775 | |
---|
| 776 | <a name="method44"></a> |
---|
| 777 | <p/><strong>int SlowJet::multiplicity(int i) </strong> <br/> |
---|
| 778 | gives the number of particles clustered into the <i>i</i>'th jet, |
---|
| 779 | |
---|
| 780 | |
---|
| 781 | <a name="method45"></a> |
---|
| 782 | <p/><strong>int SlowJet::jetAssignment(int i) </strong> <br/> |
---|
| 783 | gives the index of the jet that the particle <i>i</i> of the event |
---|
| 784 | record belongs to, or -1 if there is no jet containing particle |
---|
| 785 | <i>i</i>, |
---|
| 786 | |
---|
| 787 | |
---|
| 788 | <a name="method46"></a> |
---|
| 789 | <p/><strong>void SlowJet::removeJet(int i) </strong> <br/> |
---|
| 790 | removes the <i>i</i>'th jet, |
---|
| 791 | |
---|
| 792 | |
---|
| 793 | <a name="method47"></a> |
---|
| 794 | <p/><strong>void SlowJet::list() </strong> <br/> |
---|
| 795 | provides a listing of the above information. |
---|
| 796 | |
---|
| 797 | |
---|
| 798 | <p/> |
---|
| 799 | These are the basic methods. For more sophisticated usage |
---|
| 800 | it is possible to trace the clustering, one step at a time. If so, the |
---|
| 801 | <code>setup</code> method should be used to read in the event and find |
---|
| 802 | the initial smallest distance. Each subsequent <code>doStep</code> |
---|
| 803 | will then do one cluster joining and find the new smallest distance. |
---|
| 804 | You can therefore interrogate which clusters will be joined next |
---|
| 805 | before the joining actually is performed. Alternatively you can take |
---|
| 806 | several steps in one go, or take steps down to a predetermined number |
---|
| 807 | of jets plus remaining clusters. |
---|
| 808 | |
---|
| 809 | <a name="method48"></a> |
---|
| 810 | <p/><strong>bool SlowJet::setup( const Event& event) </strong> <br/> |
---|
| 811 | selects the particles to be analyzed, calculates initial distances, |
---|
| 812 | and finds the initial smallest distance. |
---|
| 813 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 814 | most likely the <code>pythia.event</code> one. |
---|
| 815 | |
---|
| 816 | <br/>If the routine returns <code>false</code> the setup failed, |
---|
| 817 | but currently this is not foreseen ever to happen. |
---|
| 818 | |
---|
| 819 | |
---|
| 820 | <a name="method49"></a> |
---|
| 821 | <p/><strong>bool SlowJet::doStep() </strong> <br/> |
---|
| 822 | do the next step of the clustering. This can either be that two |
---|
| 823 | clusters are joined to one, or that a cluster is promoted to a jet |
---|
| 824 | (which is discarded if its <i>pT</i> value is below |
---|
| 825 | <code>pTjetMin</code>). |
---|
| 826 | <br/>The routine will only return <code>false</code> if there are no |
---|
| 827 | clusters left. |
---|
| 828 | |
---|
| 829 | |
---|
| 830 | <a name="method50"></a> |
---|
| 831 | <p/><strong>bool SlowJet::doNSteps(int nStep) </strong> <br/> |
---|
| 832 | calls the <code>doStep()</code> method <code>nStep</code> times, |
---|
| 833 | if possible. Will return <code>false</code> if the list of clusters |
---|
| 834 | is emptied before then. The stored jet information is still perfectly |
---|
| 835 | fine; it is only the number of steps that is wrong. |
---|
| 836 | |
---|
| 837 | |
---|
| 838 | <a name="method51"></a> |
---|
| 839 | <p/><strong>bool SlowJet::stopAtN(int nStop) </strong> <br/> |
---|
| 840 | calls the <code>doStep()</code> method until a total of <code>nStop</code> |
---|
| 841 | jet and cluster objects remain. Will return <code>false</code> if this |
---|
| 842 | is not possible, specifically if the number of objects already is smaller |
---|
| 843 | than <code>nStop</code> when the method is called. The stored jet and |
---|
| 844 | cluster information is still perfectly fine; it only does not have the |
---|
| 845 | expected multiplicity. |
---|
| 846 | |
---|
| 847 | |
---|
| 848 | <a name="method52"></a> |
---|
| 849 | <p/><strong>int SlowJet::sizeAll() </strong> <br/> |
---|
| 850 | gives the total current number of jets and clusters. The jets are |
---|
| 851 | numbered 0 through <code>sizeJet() - 1</code>, while the clusters |
---|
| 852 | are numbered <code>sizeJet()</code> through <code>sizeAll() - 1</code>. |
---|
| 853 | (Internally jets and clusters are represented by two separate arrays, |
---|
| 854 | but are here presented in one flat range.) Note that the jets are ordered |
---|
| 855 | in decreasing <i>pT</i>, while the clusters are not ordered. |
---|
| 856 | |
---|
| 857 | |
---|
| 858 | <p/> |
---|
| 859 | With this extension, the methods <code>double pT(int i)</code>, |
---|
| 860 | <code>double y(int i)</code>, <code>double phi(int i)</code>, |
---|
| 861 | <code>Vec4 p(int i)</code>, <code>double m(int i)</code> and |
---|
| 862 | <code>int multiplicity(int i)</code> can be used as before. |
---|
| 863 | Furthermore, <code>list()</code> generalizes |
---|
| 864 | |
---|
| 865 | <a name="method53"></a> |
---|
| 866 | <p/><strong>void SlowJet::list(bool listAll = false, ostream& os = cout) </strong> <br/> |
---|
| 867 | provides a listing of the above information. |
---|
| 868 | <br/><code>argument</code><strong> listAll </strong> : lists both jets and clusters if <code>true</code>, |
---|
| 869 | else only jets. |
---|
| 870 | |
---|
| 871 | |
---|
| 872 | |
---|
| 873 | <p/> |
---|
| 874 | Three further methods can be used to check what will happen next. |
---|
| 875 | |
---|
| 876 | <a name="method54"></a> |
---|
| 877 | <p/><strong>int SlowJet::iNext() </strong> <br/> |
---|
| 878 | |
---|
| 879 | <strong>int SlowJet::jNext() </strong> <br/> |
---|
| 880 | |
---|
| 881 | <strong>double SlowJet::dNext() </strong> <br/> |
---|
| 882 | if the next step is to join two clusters, then the methods give |
---|
| 883 | the <i>(i,j, d_ij)</i> values, if instead to promote |
---|
| 884 | a cluster to a jet then <i>(i, -1, d_iB)</i>. |
---|
| 885 | If no clusters remain then <i>(-1, -1, 0.)</i>. Note that |
---|
| 886 | the cluster numbers are offset as described above, i.e. they begin at |
---|
| 887 | <code>sizeJet()</code>, which of course easily could be subtracted off. |
---|
| 888 | Also note that the jet and cluster lists are (moderately) reshuffled |
---|
| 889 | in each new step. |
---|
| 890 | |
---|
| 891 | |
---|
| 892 | <p/> |
---|
| 893 | Finally, and separately, the <code>SlowJetHook</code> class can be used |
---|
| 894 | for a more smart selection of which particles to include in the analysis. |
---|
| 895 | For instance, isolated and/or high-<i>pT</i> muons, electrons and |
---|
| 896 | photons should presumably be identified separately at an early stage, |
---|
| 897 | and then not clustered to jets. |
---|
| 898 | |
---|
| 899 | <p/> |
---|
| 900 | Technically, it works like with <?php $filepath = $_GET["filepath"]; |
---|
| 901 | echo "<a href='UserHooks.php?filepath=".$filepath."' target='page'>";?>User Hooks</a>. |
---|
| 902 | That is, PYTHIA contains the base class. You write a derived class. |
---|
| 903 | In the main program you create an instance of this class, and hand it |
---|
| 904 | in to <code>SlowJet</code>; in this case already as part of the |
---|
| 905 | constructor. |
---|
| 906 | |
---|
| 907 | <p/> |
---|
| 908 | The following methods should be defined in your derived class. |
---|
| 909 | |
---|
| 910 | <a name="method55"></a> |
---|
| 911 | <p/><strong>SlowJetHook::SlowJetHook() </strong> <br/> |
---|
| 912 | |
---|
| 913 | <strong>virtual SlowJetHook::~SlowJetHook() </strong> <br/> |
---|
| 914 | the constructor and destructor need not do anything, and if so you |
---|
| 915 | need not write your own destructor. |
---|
| 916 | |
---|
| 917 | |
---|
| 918 | <a name="method56"></a> |
---|
| 919 | <p/><strong>virtual bool SlowJetHook::include(int iSel, const Event& event, Vec4& pSel, double& mSel) </strong> <br/> |
---|
| 920 | is the main method that you will need to write. It will be called |
---|
| 921 | once for each final-state particle in an event, subject to the |
---|
| 922 | value of the <code>select</code> switch in the <code>SlowJet</code> |
---|
| 923 | constructor. The value <code>select = 2</code> may be convenient |
---|
| 924 | since then you do not need to remove e.g. neutrinos yourself, but |
---|
| 925 | use <code>select = 1</code> for full control. The method should then |
---|
| 926 | return <code>true</code> if you want to see particle included in the |
---|
| 927 | jet clustering, and <code>false</code> if not. |
---|
| 928 | <br/><code>argument</code><strong> iSel </strong> : is the index in the event record of the |
---|
| 929 | currently studied particle. |
---|
| 930 | |
---|
| 931 | <br/><code>argument</code><strong> event </strong> : is an object of the <code>Event</code> class, |
---|
| 932 | most likely the <code>pythia.event</code> one, where the currently |
---|
| 933 | studied particle is found. |
---|
| 934 | |
---|
| 935 | <br/><code>argument</code><strong> pSel </strong> : is at input the four-momentum of the |
---|
| 936 | currently studied particle. You can change the values, e.g. to take |
---|
| 937 | into account energy smearing in the detector, to define the initial |
---|
| 938 | cluster value, without corrupting the event record itself. |
---|
| 939 | |
---|
| 940 | <br/><code>argument</code><strong> mSel </strong> : is at input the mass of the currently studied |
---|
| 941 | particle. You can change the value, e.g. to take into account |
---|
| 942 | particle misidentification, to define the initial cluster value, |
---|
| 943 | without corrupting the event record itself. Note that the changes of |
---|
| 944 | <code>pSel</code> and <code>mSel</code> must be coordinated such that |
---|
| 945 | <i>E^2 - p^2 = m^2</i> holds. |
---|
| 946 | |
---|
| 947 | |
---|
| 948 | |
---|
| 949 | <p/> |
---|
| 950 | It is also possible to define further methods of your own. |
---|
| 951 | One such could e.g. be called directly in the main program before the |
---|
| 952 | <code>analyze</code> method is called, to identify and bookkeep |
---|
| 953 | some event properties you may not want to reanalyze for each |
---|
| 954 | individual particle. |
---|
| 955 | |
---|
| 956 | </body> |
---|
| 957 | </html> |
---|
| 958 | |
---|
| 959 | <!-- Copyright (C) 2012 Torbjorn Sjostrand --> |
---|