[1] | 1 | <chapter name="The Event Record"> |
---|
| 2 | |
---|
| 3 | <h2>The Event Record</h2> |
---|
| 4 | |
---|
| 5 | A <code>Pythia</code> instance contains two members of the |
---|
| 6 | <code>Event</code> class. The one called <code>process</code> provides |
---|
| 7 | a brief summary of the main steps of the hard process, while the |
---|
| 8 | one called <code>event</code> contains the full history. The |
---|
| 9 | user would normally interact mainly with the second one, so |
---|
| 10 | we will examplify primarily with that one. |
---|
| 11 | |
---|
| 12 | <p/> |
---|
| 13 | The <code>Event</code> class to first approximation is a vector of |
---|
| 14 | <code>Particle</code>s, so that it can expand to fit the current |
---|
| 15 | event size. The index operator is overloaded, so that e.g. |
---|
| 16 | <code>event[i]</code> corresponds to the <ei>i</ei>'th particle |
---|
| 17 | of the object <code>event</code>. Thus <code>event[i].id()</code> |
---|
| 18 | returns the identity of the <ei>i</ei>'th particle, and so on. |
---|
| 19 | Therefore the methods of the |
---|
| 20 | <code><aloc href="ParticleProperties">Particle</aloc></code> class |
---|
| 21 | are at least as essential as those of the <code>Event</code> class |
---|
| 22 | itself. |
---|
| 23 | |
---|
| 24 | <p/> |
---|
| 25 | As used inside PYTHIA, some conventions are imposed on the structure |
---|
| 26 | of the event record. Entry 0 of the <code>vector<Particle></code> |
---|
| 27 | is used to represent the event as a whole, with its total four-momentum |
---|
| 28 | and invariant mass, but does not form part of the event history. |
---|
| 29 | Lines 1 and 2 contains the two incoming beams, and only from here on |
---|
| 30 | history tracing works as could be expected. That way unassigned mother |
---|
| 31 | and daughter indices can be put 0 without ambiguity. Depending on the |
---|
| 32 | task at hand, a loop may therefore start at index 1 rather than 0 |
---|
| 33 | without any loss. Specifically, for translation to other event record |
---|
| 34 | formats such as HepMC <ref>Dob01</ref>, where the first index is 1, the |
---|
| 35 | Pythia entry 0 definitely ought to be skipped in order to minimize the |
---|
| 36 | danger of indexing errors. |
---|
| 37 | |
---|
| 38 | <p/> |
---|
| 39 | In the following we will list the methods available. |
---|
| 40 | Only a few of them have a function to fill in normal user code. |
---|
| 41 | |
---|
| 42 | <h3>Basic output methods</h3> |
---|
| 43 | |
---|
| 44 | Some methods are available to read out information on the |
---|
| 45 | current event record: |
---|
| 46 | |
---|
| 47 | <method name="Particle& Event::operator[](int i)"> |
---|
| 48 | </method> |
---|
| 49 | <methodmore name="const Particle& Event::operator[](int i)"> |
---|
| 50 | </methodmore> |
---|
| 51 | <methodmore name="Particle& Event::at(int i)"> |
---|
| 52 | returns a (<code>const</code>) reference to the <ei>i</ei>'th particle |
---|
| 53 | in the event record, which can be used to get (or set) all the |
---|
| 54 | <aloc href="ParticleProperties">properties</aloc> of this particle. |
---|
| 55 | </methodmore> |
---|
| 56 | |
---|
| 57 | <method name="int Event::size()"> |
---|
| 58 | The event size, i.e. the sie of the <code>vector<Particle></code>. |
---|
| 59 | Thus valid particles, to be accessed by the above indexing operator, |
---|
| 60 | are stored in the range <ei>0 <= i < size()</ei>. See comment |
---|
| 61 | above about the (ir)relevance of entry 0. |
---|
| 62 | </method> |
---|
| 63 | |
---|
| 64 | <method name="void Event::list()"> |
---|
| 65 | </method> |
---|
| 66 | <methodmore name="void Event::list(ostream& os)"> |
---|
| 67 | </methodmore> |
---|
| 68 | <methodmore name="void Event::list(bool showScaleAndVertex, |
---|
| 69 | bool showMothersAndDaughters = false)"> |
---|
| 70 | </methodmore> |
---|
| 71 | <methodmore name="void Event::list(bool showScaleAndVertex, |
---|
| 72 | bool showMothersAndDaughters, ostream& os)"> |
---|
| 73 | Provide a listing of the whole event, i.e. of the |
---|
| 74 | <code>vector<Particle></code>. The methods with fewer arguments |
---|
| 75 | call the final one with the respective default values, and are |
---|
| 76 | non-inlined so they can be used in a debugger. The basic identity |
---|
| 77 | code, status, mother, daughter, colour, four-momentum and mass data |
---|
| 78 | are always given, but the methods can also be called with a few |
---|
| 79 | optional arguments for further information: |
---|
| 80 | <argument name="showScaleAndVertex" default="false"> optionally give a |
---|
| 81 | second line for each particle, with the production scale (in GeV), |
---|
| 82 | the particle polarization (dimensionless), the production vertex |
---|
| 83 | (in mm or mm/c) and the invariant lifetime (also in mm/c). |
---|
| 84 | </argument> |
---|
| 85 | <argument name="showMothersAndDaughters" default="false"> |
---|
| 86 | gives a list of all daughters and mothers of a particle, as defined by |
---|
| 87 | the <code>motherList(i)</code> and <code>daughterList(i)</code> methods |
---|
| 88 | described below. It is mainly intended for debug purposes. |
---|
| 89 | </argument> |
---|
| 90 | <argument name="os" default="cout"> a reference to the <code>ostream</code> |
---|
| 91 | object to which the event listing will be directed. |
---|
| 92 | </argument> |
---|
| 93 | |
---|
| 94 | </method> |
---|
| 95 | |
---|
| 96 | <p/> |
---|
| 97 | Each <code>Particle</code> has two mother and two daughter indices. |
---|
| 98 | These may be used to encode zero, one, two or more mothers/daughters, |
---|
| 99 | depending on the combination of values and status code, according to |
---|
| 100 | well-defined <aloc href="ParticleProperties">rules</aloc>. The |
---|
| 101 | two methods below can do this job easier for you. |
---|
| 102 | |
---|
| 103 | <method name="vector<int> Event::motherList(int i)"> |
---|
| 104 | returns a vector of all the mother indices of the particle at index |
---|
| 105 | <ei>i</ei>. This list is empty for entries 0, 1 and 2, |
---|
| 106 | i.e. the "system" in line 0 is not counted as part of the history. |
---|
| 107 | Normally the list contains one or two mothers, but it can also be more, |
---|
| 108 | e.g. in string fragmentation the whole fragmenting system is counted |
---|
| 109 | as mothers to the primary hadrons. Many particles may have the same |
---|
| 110 | <code>motherList</code>. Mothers are listed in ascending order. |
---|
| 111 | </method> |
---|
| 112 | |
---|
| 113 | <method name="vector<int> Event::daughterList(int i)"> |
---|
| 114 | returns a vector of all the daughter indices of the particle at index |
---|
| 115 | <ei>i</ei>. This list is empty for a particle that did |
---|
| 116 | not decay (or, if the evolution is stopped early enough, a parton |
---|
| 117 | that did not branch), while otherwise it can contain a list of |
---|
| 118 | varying length, from one to many. For the two incoming beam particles, |
---|
| 119 | all shower initiators and beam remnants are counted as daughters, |
---|
| 120 | with the one in slot 0 being the one leading up to the hardest |
---|
| 121 | interaction. The "system" in line 0 does not have any daughters, |
---|
| 122 | i.e. is not counted as part of the history. Many partons may have the |
---|
| 123 | same <code>daughterList</code>. Daughters are listed in ascending order. |
---|
| 124 | </method> |
---|
| 125 | |
---|
| 126 | <method name="int Event::statusHepMC(int i)"> |
---|
| 127 | returns the status code according to the HepMC conventions agreed in |
---|
| 128 | February 2009. This convention does not preserve the full information |
---|
| 129 | provided by the internal PYTHIA status code, as obtained by |
---|
| 130 | <code>Particle::status()</code>, but comes reasonably close. |
---|
| 131 | The allowed output values are: |
---|
| 132 | <ul> |
---|
| 133 | <li>0 : an empty entry, with no meaningful information and therefore |
---|
| 134 | to be skipped unconditionally (should not occur in PYTHIA);</li> |
---|
| 135 | <li>1 : a final-state particle, i.e. a particle that is not decayed |
---|
| 136 | further by the generator (may also include unstable particles that |
---|
| 137 | are to be decayed later, as part of the detector simulation);</li> |
---|
| 138 | <li>2 : a decayed Standard Model hadron or tau or mu lepton, excepting |
---|
| 139 | virtual intermediate states thereof (i.e. the particle must undergo |
---|
| 140 | a normal decay, not e.g. a shower branching);</li> |
---|
| 141 | <li>3 : a documentation entry (not used in PYTHIA);</li> |
---|
| 142 | <li>4 : an incoming beam particle;</li> |
---|
| 143 | <li>11 - 200 : an intermediate (decayed/branched/...) particle that does |
---|
| 144 | not fulfill the criteria of status code 2, with a generator-dependent |
---|
| 145 | classification of its nature; in PYTHIA the absolute value of the normal |
---|
| 146 | status code is used.</li> |
---|
| 147 | </ul> |
---|
| 148 | |
---|
| 149 | </method> |
---|
| 150 | |
---|
| 151 | <h3>Further output methods</h3> |
---|
| 152 | |
---|
| 153 | The above methods are the main ones that a normal user would make |
---|
| 154 | frequent use of. There are some further methods that also could come |
---|
| 155 | in handy, in the exploration of the history of an event, but where |
---|
| 156 | the outcome is not always obvious if one is not familiar with the |
---|
| 157 | detailed structure of an event record. |
---|
| 158 | |
---|
| 159 | <method name="int Event::iTopCopy(int i)"> |
---|
| 160 | </method> |
---|
| 161 | <methodmore name="int Event::iBotCopy(int i)"> |
---|
| 162 | are used to trace carbon copies of the particle at index <ei>i</ei> up |
---|
| 163 | to its top mother or down to its bottom daughter. If there are no such |
---|
| 164 | carbon copies, <ei>i</ei> itself will be returned. A carbon copy is |
---|
| 165 | when the "same" particle appears several times in the event record, but |
---|
| 166 | with changed momentum owing to recoil effects. |
---|
| 167 | </methodmore> |
---|
| 168 | |
---|
| 169 | <method name="int Event::iTopCopyId(int i)"> |
---|
| 170 | </method> |
---|
| 171 | <methodmore name="int Event::iBotCopyId(int i)"> |
---|
| 172 | also trace top mother and bottom daughter, but do not require carbon |
---|
| 173 | copies, only that one can find an unbroken chain, of mothers or daughters, |
---|
| 174 | with the same flavour <code>id</code> code. When it encounters ambiguities, |
---|
| 175 | say a <ei>g -> g g</ei> branching or a <ei>u u -> u u</ei> hard scattering, |
---|
| 176 | it will stop the tracing and return the current position. It can be confused |
---|
| 177 | by nontrivial flavour changes, e.g. a hard process <ei>u d -> d u</ei> |
---|
| 178 | by <ei>W^+-</ei> exchange will give the wrong answer. These methods |
---|
| 179 | therefore are of limited use for common particles, in particular for the |
---|
| 180 | gluon, but should work well for "rare" particles. |
---|
| 181 | </method> |
---|
| 182 | |
---|
| 183 | <method name="vector<int> Event::sisterList(int i)"> |
---|
| 184 | returns a vector of all the sister indices of the particle at index |
---|
| 185 | <ei>i</ei>, i.e. all the daughters of the first mother, except the |
---|
| 186 | particle itself. |
---|
| 187 | </method> |
---|
| 188 | |
---|
| 189 | <method name="vector<int> Event::sisterListTopBot(int i, |
---|
| 190 | bool widenSearch = true)"> |
---|
| 191 | returns a vector of all the sister indices of the particle at index |
---|
| 192 | <ei>i</ei>, tracking up and back down through carbon copies |
---|
| 193 | if required. That is, the particle is first traced up with |
---|
| 194 | <code>iTopCopy()</code> before its mother is found, and then all |
---|
| 195 | the particles in the <code>daughterList()</code> of this mother are |
---|
| 196 | traced down with <code>iBotCopy()</code>, omitting the original |
---|
| 197 | particle itself. Any non-final particles are removed from the list. |
---|
| 198 | Should this make the list empty the search criterion is widened so that |
---|
| 199 | all final daughters are allowed, not only carbon-copy ones. A second |
---|
| 200 | argument <code>false</code> inhibits the second step, and increases |
---|
| 201 | the risk that an empty list is returned. A typical example of this |
---|
| 202 | is for ISR cascades, e.g. <ei>e -> e gamma</ei> where the photon |
---|
| 203 | may not have any obvious sister in the final state if the bottom copy |
---|
| 204 | of the photon is an electron that annihilates and thus is not part of |
---|
| 205 | the final state. |
---|
| 206 | </method> |
---|
| 207 | |
---|
| 208 | <method name="bool Event::isAncestor(int i, int iAncestor)"> |
---|
| 209 | traces the particle <ei>i</ei> upwards through mother, grandmother, |
---|
| 210 | and so on, until either <ei>iAncestor</ei> is found or the top of |
---|
| 211 | the record is reached. Normally one unique mother is required, |
---|
| 212 | as is the case e.g. in decay chains or in parton showers, so that |
---|
| 213 | e.g. the tracing through a hard scattering would not work. For |
---|
| 214 | hadronization, first-rank hadrons are identified with the respective |
---|
| 215 | string endpoint quark, which may be useful e.g. for <ei>b</ei> physics, |
---|
| 216 | while higher-rank hadrons give <code>false</code>. Currently also |
---|
| 217 | ministrings that collapsed to one single hadron and junction topologies |
---|
| 218 | give <code>false</code>. |
---|
| 219 | </method> |
---|
| 220 | |
---|
| 221 | <p/> |
---|
| 222 | One data member in an <code>Event</code> object is used to keep track |
---|
| 223 | of the largest <code>col()</code> or <code>acol()</code> colour tag set |
---|
| 224 | so far, so that new ones do not clash. |
---|
| 225 | |
---|
| 226 | <modeopen name="Event:startColTag" default="100" min="0" max="1000"> |
---|
| 227 | This sets the initial colour tag value used, so that the first one |
---|
| 228 | assigned is <code>startColTag + 1</code>, etc. The Les Houches accord |
---|
| 229 | <ref>Boo01</ref> suggests this number to be 500, but 100 works equally |
---|
| 230 | well. |
---|
| 231 | </modeopen> |
---|
| 232 | |
---|
| 233 | <method name="void Event::initColTag(int colTag = 0)"> |
---|
| 234 | forces the current colour tag value to be the larger of the input |
---|
| 235 | <code>colTag</code> and the above <code>Event:startColTag</code> |
---|
| 236 | values. |
---|
| 237 | </method> |
---|
| 238 | |
---|
| 239 | <method name="int Event::lastColTag()"> |
---|
| 240 | returns the current maximum colour tag. |
---|
| 241 | </method> |
---|
| 242 | |
---|
| 243 | <method name="int Event::nextColTag()"> |
---|
| 244 | increases the current maximum colour tag by one and returns this |
---|
| 245 | new value. This method is used whenever a new colour tag is needed. |
---|
| 246 | </method> |
---|
| 247 | |
---|
| 248 | <p/> |
---|
| 249 | Many event properties are accessible via the <code>Info</code> class, |
---|
| 250 | <aloc href="EventInformation">see here</aloc>. Since they are used |
---|
| 251 | directly in the event generation, a few are stored directly in the |
---|
| 252 | <code>Event</code> class, however. |
---|
| 253 | |
---|
| 254 | <method name="void Event::scale( double scaleIn)"> |
---|
| 255 | </method> |
---|
| 256 | <methodmore name="double Event::scale()"> |
---|
| 257 | set or get the scale (in GeV) of the hardest process in the event. |
---|
| 258 | Matches the function of the <code>scale</code> variable in the |
---|
| 259 | <aloc href="LesHouchesAccord">Les Houches Accord</aloc>. |
---|
| 260 | </methodmore> |
---|
| 261 | |
---|
| 262 | <method name="void Event::scaleSecond( double scaleSecondIn)"> |
---|
| 263 | </method> |
---|
| 264 | <methodmore name="double Event::scaleSecond()"> |
---|
| 265 | set or get the scale (in GeV) of a second hard process in the event, |
---|
| 266 | in those cases where such a one |
---|
| 267 | <aloc href="SecondHardProcess">has been requested</aloc>. |
---|
| 268 | </methodmore> |
---|
| 269 | |
---|
| 270 | <h3>Constructors and modifications of the event record</h3> |
---|
| 271 | |
---|
| 272 | Although you would not normally need to create your own |
---|
| 273 | <code>Event</code> instance, there may be times where that |
---|
| 274 | could be convenient. The typical exampel would be if you want to |
---|
| 275 | create a new event record that is the sum of a few different ones, |
---|
| 276 | e.g. if you want to simulate pileup events. There may also be cases |
---|
| 277 | where you want to add one or a few particles to an existing event |
---|
| 278 | record. |
---|
| 279 | |
---|
| 280 | <method name="Event::Event(int capacity = 100)"> |
---|
| 281 | creates an empty event record, but with a reserved size |
---|
| 282 | <ei>capacity</ei> for the <code>Particle</code> vector. |
---|
| 283 | </method> |
---|
| 284 | |
---|
| 285 | <method name="Event& Event::operator=(const Event& oldEvent)"> |
---|
| 286 | copies the input event record. |
---|
| 287 | </method> |
---|
| 288 | |
---|
| 289 | <method name="Event& Event::operator+=(const Event& addEvent)"> |
---|
| 290 | appends an event to an existing one. For the appended particles |
---|
| 291 | mother, daughter and colour tags are shifted to make a consistent |
---|
| 292 | record. The zeroth particle of the appended event is not copied, |
---|
| 293 | but the zeroth particle of the combined event is updated to the |
---|
| 294 | full energy-momentum content. |
---|
| 295 | </method> |
---|
| 296 | |
---|
| 297 | <method name="void Event::init(string headerIn = "", |
---|
| 298 | ParticleData* particleDataPtrIn = 0, int startColTagIn = 100)"> |
---|
| 299 | initializes colour, the pointer to the particle database, and the |
---|
| 300 | header specification used for the event listing. We remind that a |
---|
| 301 | <code>Pythia</code> object contains two event records |
---|
| 302 | <code>process</code> and <code>event</code>. Thus one may e.g. |
---|
| 303 | call either <code>pythia.process.list()</code> or |
---|
| 304 | <code>pythia.event.list()</code>. To distinguish those two rapidly |
---|
| 305 | at visual inspection, the <code>"Pythia Event Listing"</code> header |
---|
| 306 | is printed out differently, in one case adding |
---|
| 307 | <code>"(hard process)"</code> and in the other |
---|
| 308 | <code>"(complete event)"</code>. When <code>+=</code> is used to |
---|
| 309 | append an event, the modified event is printed with |
---|
| 310 | <code>"(combination of several events)"</code> as a reminder. |
---|
| 311 | </method> |
---|
| 312 | |
---|
| 313 | <method name="void Event::clear()"> |
---|
| 314 | empties event record. Specifically the <code>Particle</code> vector |
---|
| 315 | size is reset to zero. |
---|
| 316 | </method> |
---|
| 317 | |
---|
| 318 | <method name="void Event::reset()"> |
---|
| 319 | empties the event record, as <code>clear()</code> above, but then |
---|
| 320 | fills the zero entry of the <code>Particle</code> vector with the |
---|
| 321 | pseudoparticle used to represent the event as a whole. At this point |
---|
| 322 | the pseudoparticle is not assigned any momentum or mass. |
---|
| 323 | </method> |
---|
| 324 | |
---|
| 325 | <method name="void Event::popBack(int n = 1)"> |
---|
| 326 | removes the last <ei>n</ei> particle entries; must be a positive |
---|
| 327 | number. |
---|
| 328 | </method> |
---|
| 329 | |
---|
| 330 | <method name="int Event::append(Particle entryIn)"> |
---|
| 331 | appends a particle to the bottom of the event record and |
---|
| 332 | returns the index of this position. |
---|
| 333 | </method> |
---|
| 334 | |
---|
| 335 | <method name="int Event::append(int id, int status, int mother1, |
---|
| 336 | int mother2, int daughter1, int daughter2, int col, int acol, |
---|
| 337 | double px, double py, double pz, double e, double m = 0., |
---|
| 338 | double scale = 0., double pol = 9.)"> |
---|
| 339 | appends a particle to the bottom of the event record and |
---|
| 340 | returns the index of this position; |
---|
| 341 | <aloc href="ParticleProperties">see here</aloc> for the meaning |
---|
| 342 | of the various particle properties. |
---|
| 343 | </method> |
---|
| 344 | |
---|
| 345 | <method name="int Event::append(int id, int status, int mother1, |
---|
| 346 | int mother2, int daughter1, int daughter2, int col, int acol, |
---|
| 347 | Vec4 p, double m = 0., double scale = 0., double pol = 9.)"> |
---|
| 348 | appends a particle to the bottom of the event record and |
---|
| 349 | returns the index of this position, as above but with four-momentum |
---|
| 350 | as a <code>Vec4</code>. |
---|
| 351 | </method> |
---|
| 352 | |
---|
| 353 | <method name="int Event::append(int id, int status, int col, int acol, |
---|
| 354 | double px, double py, double pz, double e, double m = 0., |
---|
| 355 | double scale = 0., double pol = 9.)"> |
---|
| 356 | </method> |
---|
| 357 | <methodmore name="int Event::append(int id, int status, int col, |
---|
| 358 | int acol, Vec4 p, double m = 0., double scale = 0., double pol = 9.)"> |
---|
| 359 | appends a particle to the bottom of the event record and |
---|
| 360 | returns the index of this position, as above but with vanishing |
---|
| 361 | (i.e. zero) mother and daughter indices. |
---|
| 362 | </method> |
---|
| 363 | |
---|
| 364 | <method name="int Event::setPDTPtr(int iSet = -1)"> |
---|
| 365 | send in a pointer to the <code>ParticleData</code> database for |
---|
| 366 | particle <code>iSet</code>, by default the most recently appended |
---|
| 367 | particle. Also generates a pointer to the |
---|
| 368 | <code>ParticleDataEntry</code> object of the identity code |
---|
| 369 | of the particle. |
---|
| 370 | </method> |
---|
| 371 | |
---|
| 372 | <method name="int Event::copy(int iCopy, int newStatus = 0)"> |
---|
| 373 | copies the existing particle in entry <code>iCopy</code> to the |
---|
| 374 | bottom of the event record and returns the index of this position. |
---|
| 375 | By default, i.e. with <code>newStatus = 0</code>, everything is |
---|
| 376 | copied precisely as it is, which means that history information |
---|
| 377 | has to be modified further by hand to make sense. With a positive |
---|
| 378 | <code>newStatus</code>, the new copy is set up to be the daughter of |
---|
| 379 | the old, with status code <code>newStatus</code>, while the status |
---|
| 380 | code of <code>iCopy</code> is negated. With a negative |
---|
| 381 | <code>newStatus</code>, the new copy is instead set up to be the |
---|
| 382 | mother of <code>iCopy</code>. An attempt to copy an out-of-range |
---|
| 383 | entry will return -1. |
---|
| 384 | </method> |
---|
| 385 | |
---|
| 386 | <method name="Particle& Event::back()"> |
---|
| 387 | returns a reference to the last particle in the event record. |
---|
| 388 | </method> |
---|
| 389 | |
---|
| 390 | <method name="void Event::restorePtrs()"> |
---|
| 391 | each particle in the event record has a pointer to the whole database |
---|
| 392 | and another to the particle species itself, used to find some particle |
---|
| 393 | properties. The latter pointer is automatically set/changed whenever |
---|
| 394 | the particle identity is set/changed by one of the normal methods. |
---|
| 395 | (It is the "changed" part that prompts the inclusion of a pointer |
---|
| 396 | to the whole database.) Of course the pointer values are specific to |
---|
| 397 | the memory locations of the current run, and so it has no sense to |
---|
| 398 | save them if events are written to file. Should you use some |
---|
| 399 | persistency scheme that bypasses the normal methods when the event is |
---|
| 400 | read back in, you can use <code>restorePtrs()</code> afterwards to set |
---|
| 401 | these pointers appropriately. |
---|
| 402 | </method> |
---|
| 403 | |
---|
| 404 | <p/> |
---|
| 405 | A few methods exist to rotate and boost events. These derive from the |
---|
| 406 | <aloc href="FourVectors">Vec4</aloc> methods, and affect both the |
---|
| 407 | momentum and the vertex (position) components of all particles. |
---|
| 408 | |
---|
| 409 | <method name="void Event::rot(double theta, double phi)"> |
---|
| 410 | rotate all particles in the event by this polar and azimuthal angle |
---|
| 411 | (expressed in radians). |
---|
| 412 | </method> |
---|
| 413 | |
---|
| 414 | <method name="void Event::bst(double betaX, double betaY, double betaZ)"> |
---|
| 415 | </method> |
---|
| 416 | <methodmore name="void Event::bst(double betaX, double betaY, |
---|
| 417 | double betaZ, double gamma)"> |
---|
| 418 | </methodmore> |
---|
| 419 | <methodmore name="void Event::bst(const Vec4& vec)"> |
---|
| 420 | boost all particles in the event by this three-vector. |
---|
| 421 | Optionally you may provide the <ei>gamma</ei> value as a fourth argument, |
---|
| 422 | which may help avoid roundoff errors for big boosts. You may alternatively |
---|
| 423 | supply a <code>Vec4</code> four-vector, in which case the boost vector |
---|
| 424 | becomes <ei>beta = p/E</ei>. |
---|
| 425 | </methodmore> |
---|
| 426 | |
---|
| 427 | <method name="void Event::rotbst(const RotBstMatrix& M)"> |
---|
| 428 | rotate and boost by the combined action encoded in the |
---|
| 429 | <code><aloc href="FourVectors">RotBstMatrix</aloc> M</code>. |
---|
| 430 | </method> |
---|
| 431 | |
---|
| 432 | <h3>The Junction Class</h3> |
---|
| 433 | |
---|
| 434 | The event record also contains a vector of junctions, which often |
---|
| 435 | is empty or else contains only a very few per event. Methods are |
---|
| 436 | available to add further junctions or query the current junction list. |
---|
| 437 | This is only for the expert user, however, and is not discussed |
---|
| 438 | further here, but only the main points. |
---|
| 439 | |
---|
| 440 | <p/> |
---|
| 441 | A junction stores the properites associated with a baryon number that |
---|
| 442 | is fully resolved, i.e. where three different colour indices are |
---|
| 443 | involved. There are two main applications, |
---|
| 444 | <ol> |
---|
| 445 | <li>baryon beams, where at least two valence quarks are kicked out, |
---|
| 446 | and so the motion of the baryon number is notrivial;</li> |
---|
| 447 | <li>baryon-number violating processes, e.g. in SUSY with broken |
---|
| 448 | <ei>R</ei>-parity.</li> |
---|
| 449 | </ol> |
---|
| 450 | Information on junctions is set, partly in the process generation, |
---|
| 451 | partly in the beam remnants machinery, and used by the fragmentation |
---|
| 452 | routines, but the normal user does not have to know the details. |
---|
| 453 | |
---|
| 454 | <p/> |
---|
| 455 | For each junction, information is stored on the kind of junction, and |
---|
| 456 | on the three (anti)colour indices that are involved in the junction. |
---|
| 457 | The possibilities foreseen are: |
---|
| 458 | <ul> |
---|
| 459 | <li><code>kind = 1</code> : incoming colourless particle to three |
---|
| 460 | outgoing colours (e.g. baryon beam remnant or |
---|
| 461 | <ei>neutralino -> q q q</ei>);</li> |
---|
| 462 | <li><code>kind = 2</code> : incoming colourless particle to three |
---|
| 463 | outgoing anticolours;</li> |
---|
| 464 | <li><code>kind = 3</code> : one incoming anticolour (stored first) |
---|
| 465 | and two outgoing colours (e.g. antisquark decaying to two quarks, or |
---|
| 466 | gluino decay to three quarks);</li> |
---|
| 467 | <li><code>kind = 4</code> : one incoming colour (stored first) and two |
---|
| 468 | outgoing anticolours (e.g. squark decaying to two antiquarks, or |
---|
| 469 | gluino decaying to three antiquarks);</li> |
---|
| 470 | <li><code>kind = 5</code> : two incoming anticolours (stored first) |
---|
| 471 | and one outgoing colour (e.g. resonant squark production through RPV);</li> |
---|
| 472 | <li><code>kind = 6</code> : two incoming colours (stored first) |
---|
| 473 | and one outgoing anticolour (e.g. resonant antisquark production |
---|
| 474 | through RPV); |
---|
| 475 | </li> |
---|
| 476 | </ul> |
---|
| 477 | The odd (even) <code>kind</code> codes corresponds to a +1 (-1) change in |
---|
| 478 | baryon number across the junction. |
---|
| 479 | |
---|
| 480 | <p/> |
---|
| 481 | The kind and colour information in the list of junctions can be set |
---|
| 482 | or read with methods of the <code>Event</code> class, but are not of |
---|
| 483 | common interest and so not described here. |
---|
| 484 | |
---|
| 485 | <p/> |
---|
| 486 | A listing of current junctions can be obtained with the |
---|
| 487 | <code>listJunctions()</code> method. |
---|
| 488 | |
---|
| 489 | <h3>Subsystems</h3> |
---|
| 490 | |
---|
| 491 | Separate from the event record as such, but closely tied to it is the |
---|
| 492 | <code><aloc href="AdvancedUsage">PartonSystems</aloc></code> class, |
---|
| 493 | which mainly stores the parton indices of incoming and outgoing partons, |
---|
| 494 | classified by collision subsystem. Such information is needed to |
---|
| 495 | interleave multiparton interactions, initial-state showers and final-state |
---|
| 496 | showers, and append beam remnants. It could also be used in other places. |
---|
| 497 | It is intended to be accessed only by experts, such as implementors of |
---|
| 498 | <aloc href="ImplementNewShowers">new showering models</aloc>. |
---|
| 499 | |
---|
| 500 | </chapter> |
---|
| 501 | |
---|
| 502 | <!-- Copyright (C) 2012 Torbjorn Sjostrand --> |
---|