GHepRecord.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2019, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5  or see $GENIE/LICENSE
6 
7  Author: Costas Andreopoulos <costas.andreopoulos \at stfc.ac.uk>
8  University of Liverpool & STFC Rutherford Appleton Lab
9 
10  For the class documentation see the corresponding header file.
11 
12  Important revisions after version 2.0.0 :
13  @ Mar 08, 2008 - CA
14  Modified the algorithm for compactifying the daughter lists. The new
15  algorithm is more robust and avoids failure modes of the old algorithm
16  which appeared once simulation of nuclear de-excitations was enabled.
17  @ Jun 20, 2008 - CA
18  Fixed memory leak in Print()
19  @ Jan 28, 2009 - CA
20  When checking for energy / momentum conservation in Print(), use the new
21  kIStFinalStateNuclearRemnant status code for nuclear remnants (previously
22  marked as kIStStableFinalState).
23  @ Sep 15, 2009 - CA
24  IsNucleus() is no longer available in GHepParticle. Use pdg::IsIon().
25  Print-out the particle 'rescattering code' (if it was set).
26  @ May 05, 2010 - CR
27  Adding special ctor for ROOT I/O purposes so as to avoid memory leak due to
28  memory allocated in the default ctor when objects of this class are read by
29  the ROOT Streamer.
30  @ Sep 26, 2011 - CA
31  Demote a few messages from `warning' to `notice'.
32  @ Nov 17, 2011 - CA
33  Added `GEvGenMode_t EventGenerationMode(void) const'
34  @ Nov 28, 2011 - CA
35  HitNucleon() can return a hit nucleon cluster too, as needed for MEC.
36  @ Jan 29, 2013 - CA
37  Demote a few messages from `notice' to `info'.
38  @ Jan 31, 2013 - CA
39  Added static SetPrintLevel(int print_level) and corresponding static
40  data member. $GHEPPRINTLEVEL env var is no longer used.
41  @ Feb 01, 2013 - CA
42  The GUNPHYSMASK env. var is no longer used. Added a private data member to
43  store the bit-field mask in the GHEP record. Added GHepRecord::Accept() and
44  GHepRecord::SetUnphysEventMask(). Tweaks in Print().
45  @ Feb 06, 2013 - CA
46  Added KinePhaseSpace_t fDiffXSecPhSp prov data members to specify which
47  differential cross-section value is stored in fDiffXSec. Added method to
48  set it and tweaked Print() accordingly.
49  @ May 02, 2013 - CA
50  Added `KinePhaseSpace_t DiffXSecVars(void) const' to return fDiffXSecPhSp.
51 
52 */
53 //____________________________________________________________________________
54 
55 #include <cassert>
56 #include <algorithm>
57 #include <iomanip>
58 
59 #include <TLorentzVector.h>
60 #include <TVector3.h>
61 #include <TSystem.h>
62 #include <TRootIOCtor.h>
63 
75 
76 using std::endl;
77 using std::setw;
78 using std::setprecision;
79 using std::setfill;
80 using std::ios;
81 
82 using namespace genie;
83 
85 
86 int GHepRecord::fPrintLevel = 3;
87 
88 //___________________________________________________________________________
90  ostream & operator << (ostream & stream, const GHepRecord & rec)
91  {
92  rec.Print(stream);
93  return stream;
94  }
95 }
96 //___________________________________________________________________________
98 TClonesArray("genie::GHepParticle")
99 {
100  this->InitRecord();
101 }
102 //___________________________________________________________________________
104 TClonesArray("genie::GHepParticle", size)
105 {
106  this->InitRecord();
107 }
108 //___________________________________________________________________________
110 TClonesArray("genie::GHepParticle", record.GetEntries())
111 {
112  this->InitRecord();
113  this->Copy(record);
114 }
115 //___________________________________________________________________________
116 GHepRecord::GHepRecord(TRootIOCtor*) :
117 TClonesArray("genie::GHepParticle"),
118 fInteraction(0),
119 fVtx(0),
120 fEventFlags(0),
121 fEventMask(0),
122 fWeight(0.),
123 fProb(0.),
124 fXSec(0.),
125 fDiffXSec(0.)
126 {
127 
128 }
129 //___________________________________________________________________________
131 {
132  this->CleanRecord();
133 }
134 //___________________________________________________________________________
136 {
137  if(!fInteraction) {
138  LOG("GHEP", pWARN) << "Returning NULL interaction";
139  }
140  return fInteraction;
141 }
142 //___________________________________________________________________________
144 {
145  fInteraction = interaction;
146 }
147 //___________________________________________________________________________
149 {
150 // Returns the GHepParticle from the specified position of the event record.
151 
152  if( position >=0 && position < this->GetEntries() ) {
153  GHepParticle * particle = (GHepParticle *) (*this)[position];
154  if(particle) return particle;
155  }
156  LOG("GHEP", pINFO)
157  << "No particle found with: (pos = " << position << ")";
158 
159  return 0;
160 }
161 //___________________________________________________________________________
163  int pdg, GHepStatus_t status, int start) const
164 {
165 // Returns the first GHepParticle with the input pdg-code and status
166 // starting from the specified position of the event record.
167 
168  int nentries = this->GetEntries();
169  for(int i = start; i < nentries; i++) {
170  GHepParticle * p = (GHepParticle *) (*this)[i];
171  if(p->Status() == status && p->Pdg() == pdg) return p;
172  }
173 
174  LOG("GHEP", pINFO)
175  << "No particle found with: (pos >= " << start
176  << ", pdg = " << pdg << ", ist = " << status << ")";
177 
178  return 0;
179 }
180 //___________________________________________________________________________
182  int pdg, GHepStatus_t status, int start) const
183 {
184 // Returns the position of the first GHepParticle with the input pdg-code
185 // and status starting from the specified position of the event record.
186 
187  int nentries = this->GetEntries();
188  for(int i = start; i < nentries; i++) {
189  GHepParticle * p = (GHepParticle *) (*this)[i];
190  if(p->Status() == status && p->Pdg() == pdg) return i;
191  }
192 
193  LOG("GHEP", pINFO)
194  << "No particle found with: (pos >= " << start
195  << ", pdg = " << pdg << ", ist = " << status << ")";
196 
197  return -1;
198 }
199 //___________________________________________________________________________
201 {
202 // Returns the position of the first match with the specified GHepParticle
203 // starting from the specified position of the event record.
204 
205  int nentries = this->GetEntries();
206  for(int i = start; i < nentries; i++) {
207  GHepParticle * p = (GHepParticle *) (*this)[i];
208  if( p->Compare(particle) ) return i;
209  }
210 
211  LOG("GHEP", pINFO)
212  << "No particle found with pos >= " << start
213  << " matching particle: " << *particle;
214 
215  return -1;
216 }
217 //___________________________________________________________________________
219 {
220 // Returns a list of all stable descendants of the GHEP entry in the input
221 // slot. The user adopts the output vector.
222 
223  // return null if particle index is out of range
224  int nentries = this->GetEntries();
225  if(position<0 || position>=nentries) return 0;
226 
227  vector<int> * descendants = new vector<int>;
228 
229  // return itself if it is a stable final state particle
230  if(this->Particle(position)->Status() == kIStStableFinalState) {
231  descendants->push_back(position);
232  return descendants;
233  }
234 
235  for(int i = 0; i < nentries; i++) {
236  if(i==position) continue;
237  GHepParticle * p = (GHepParticle *) (*this)[i];
238  if(p->Status() != kIStStableFinalState) continue;
239  bool is_descendant=false;
240  int mom = p->FirstMother();
241  while(mom>-1) {
242  if(mom==position) is_descendant=true;
243  if(is_descendant) {
244  descendants->push_back(i);
245  break;
246  }
247  mom = this->Particle(mom)->FirstMother();
248  }
249  }
250  return descendants;
251 }
252 //___________________________________________________________________________
254 {
255  GHepParticle * p0 = this->Particle(0);
256  if(!p0) return kGMdUnknown;
257  GHepParticle * p1 = this->Particle(1);
258  if(!p1) return kGMdUnknown;
259 
260  int p0pdg = p0->Pdg();
261  GHepStatus_t p0st = p0->Status();
262  int p1pdg = p1->Pdg();
263  GHepStatus_t p1st = p1->Status();
264 
265  // In lepton+nucleon/nucleus mode, the 1st entry in the event record
266  // is a charged or neutral lepton with status code = kIStInitialState
267  if( pdg::IsLepton(p0pdg) && p0st == kIStInitialState )
268  {
269  return kGMdLeptonNucleus;
270  }
271 
272  // In dark matter mode, the 1st entry in the event record is a dark
273  // matter particle
274  if( pdg::IsDarkMatter(p0pdg) && p0st == kIStInitialState )
275  {
276  return kGMdDarkMatterNucleus;
277  }
278 
279  // In hadron+nucleon/nucleus mode, the 1st entry in the event record
280  // is a hadron with status code = kIStInitialState and the 2nd entry
281  // is a nucleon or nucleus with status code = kIStInitialState
282  if( pdg::IsHadron(p0pdg) && p0st == kIStInitialState )
283  {
284  if( (pdg::IsIon(p1pdg) || pdg::IsNucleon(p1pdg)) && p1st == kIStInitialState)
285  {
286  return kGMdHadronNucleus;
287  }
288  }
289 
290  // As above, with a photon as a probe
291  if( p0pdg == kPdgGamma && p0st == kIStInitialState )
292  {
293  if( (pdg::IsIon(p1pdg) || pdg::IsNucleon(p1pdg)) && p1st == kIStInitialState)
294  {
295  return kGMdPhotonNucleus;
296  }
297  }
298 
299  // In nucleon decay mode,
300  // - [if the decayed nucleon was a bound one] the 1st entry in the event
301  // record is a nucleus with status code = kIStInitialState and the
302  // 2nd entry is a nucleon with code = kIStDecayedState
303  // - [if the decayed nucleon was a free one] the first entry in the event
304  // record is a nucleon with status code = kIStInitialState and it has a
305  // single daughter which is a nucleon with status code = kIStDecayedState.
306 
307  if( pdg::IsIon(p0pdg) && p0st == kIStInitialState &&
308  pdg::IsNucleon(p1pdg) && p1st == kIStDecayedState)
309  {
310  return kGMdNucleonDecay;
311  }
312  if( pdg::IsNucleon(p0pdg) && p0st == kIStInitialState &&
313  pdg::IsNucleon(p1pdg) && p1st == kIStDecayedState)
314  {
315  return kGMdNucleonDecay;
316  }
317 
318  return kGMdUnknown;
319 }
320 //___________________________________________________________________________
322 {
323 // Returns the GHepParticle representing the probe (neutrino, e,...).
324 
325  int ipos = this->ProbePosition();
326  if(ipos>-1) return this->Particle(ipos);
327  return 0;
328 }
329 //___________________________________________________________________________
331 {
332 // Returns the GHepParticle representing the target / initial state nucleus,
333 // or 0 if it does not exist.
334 
335  int ipos = this->TargetNucleusPosition();
336  if(ipos>-1) return this->Particle(ipos);
337  return 0;
338 }
339 //___________________________________________________________________________
341 {
342 // Returns the GHepParticle representing the remnant nucleus,
343 // or 0 if it does not exist.
344 
345  int ipos = this->RemnantNucleusPosition();
346  if(ipos>-1) return this->Particle(ipos);
347  return 0;
348 }
349 //___________________________________________________________________________
351 {
352 // Returns the GHepParticle representing the struck nucleon, or 0 if it does
353 // not exist.
354 
355  int ipos = this->HitNucleonPosition();
356  if(ipos>-1) return this->Particle(ipos);
357  return 0;
358 }
359 //___________________________________________________________________________
361 {
362 // Returns the GHepParticle representing the struck electron, or 0 if it does
363 // not exist.
364 
365  int ipos = this->HitElectronPosition();
366  if(ipos>-1) return this->Particle(ipos);
367  return 0;
368 }
369 //___________________________________________________________________________
371 {
372 // Returns the GHepParticle representing the final state primary lepton.
373 
374  int ipos = this->FinalStatePrimaryLeptonPosition();
375  if(ipos>-1) return this->Particle(ipos);
376  return 0;
377 }
378 //___________________________________________________________________________
380 {
381 // Returns the GHepParticle representing the sum of the DIS pre-fragm f/s
382 // hadronic system, or 0 if it does not exist.
383 
384  int ipos = this->FinalStateHadronicSystemPosition();
385  if(ipos>-1) return this->Particle(ipos);
386  return 0;
387 }
388 //___________________________________________________________________________
390 {
391 // Returns the GHEP position of the GHepParticle representing the probe
392 // (neutrino, e,...).
393 
394  // The probe is *always* at slot 0.
396  if(mode == kGMdLeptonNucleus ||
397  mode == kGMdDarkMatterNucleus ||
398  mode == kGMdHadronNucleus ||
399  mode == kGMdPhotonNucleus)
400  {
401  return 0;
402  }
403  return -1;
404 }
405 //___________________________________________________________________________
407 {
408 // Returns the GHEP position of the GHepParticle representing the target
409 // nucleus - or -1 if the interaction takes place at a free nucleon.
410 
412 
413  if(mode == kGMdLeptonNucleus ||
414  mode == kGMdDarkMatterNucleus ||
415  mode == kGMdHadronNucleus ||
416  mode == kGMdPhotonNucleus)
417  {
418  GHepParticle * p = this->Particle(1); // If exists, it will be at slot 1
419  if(!p) return -1;
420  int pdgc = p->Pdg();
421  if(pdg::IsIon(pdgc) && p->Status()==kIStInitialState) return 1;
422  }
423  if(mode == kGMdNucleonDecay) {
424  GHepParticle * p = this->Particle(0); // If exists, it will be at slot 0
425  if(!p) return -1;
426  int pdgc = p->Pdg();
427  if(pdg::IsIon(pdgc) && p->Status()==kIStInitialState) return 0;
428  }
429 
430  return -1;
431 }
432 //___________________________________________________________________________
434 {
435 // Returns the GHEP position of the GHepParticle representing the remnant
436 // nucleus - or -1 if the interaction takes place at a free nucleon.
437 
438  GHepParticle * p = this->TargetNucleus();
439  if(!p) return -1;
440 
441  int dau1 = p->FirstDaughter();
442  int dau2 = p->LastDaughter();
443 
444  if(dau1==-1 && dau2==-1) return -1;
445 
446  for(int i=dau1; i<=dau2; i++) {
447  GHepParticle * dp = this->Particle(i);
448  int dpdgc = dp->Pdg();
449  if(pdg::IsIon(dpdgc) && dp->Status()==kIStStableFinalState) return i;
450  }
451  return -1;
452 }
453 //___________________________________________________________________________
455 {
456 // Returns the GHEP position of the GHepParticle representing the hit nucleon.
457 // If a struck nucleon is set it will be at slot 2 (for scattering off nuclear
458 // targets) or at slot 1 (for free nucleon scattering).
459 // If the struck nucleon is not set (eg coherent scattering, ve- scattering)
460 // it returns 0.
461 
462  GHepParticle * nucleus = this->TargetNucleus();
463 
464  int ipos = (nucleus) ? 2 : 1;
465  GHepStatus_t ist = (nucleus) ? kIStNucleonTarget : kIStInitialState;
466 
467  GHepParticle * p = this->Particle(ipos);
468  if(!p) return -1;
469 
470 // bool isN = pdg::IsNeutronOrProton(p->Pdg());
471  bool isN = pdg::IsNucleon(p->Pdg()) || pdg::Is2NucleonCluster(p->Pdg());
472  if(isN && p->Status()==ist) return ipos;
473 
474  return -1;
475 }
476 //___________________________________________________________________________
478 {
479 // Returns the GHEP position of the GHepParticle representing a hit electron.
480 // Same as above..
481 
482  GHepParticle * nucleus = this->TargetNucleus();
483 
484  int ipos = (nucleus) ? 2 : 1;
485 
486  GHepParticle * p = this->Particle(ipos);
487  if(!p) return -1;
488 
489  bool ise = pdg::IsElectron(p->Pdg());
490  if(ise && p->Status()==kIStInitialState) return ipos;
491 
492  return -1;
493 }
494 //___________________________________________________________________________
496 {
497 // Returns the GHEP position GHepParticle representing the final state
498 // primary lepton.
499 
500  GHepParticle * probe = this->Probe();
501  if(!probe) return -1;
502 
503  int ifsl = probe->FirstDaughter();
504  return ifsl;
505 }
506 //___________________________________________________________________________
508 {
509  return this->ParticlePosition(
511 }
512 //___________________________________________________________________________
513 unsigned int GHepRecord::NEntries(int pdg, GHepStatus_t ist, int start) const
514 {
515  unsigned int nentries = 0;
516 
517  for(int i = start; i < this->GetEntries(); i++) {
518  GHepParticle * p = (GHepParticle *) (*this)[i];
519  if(p->Pdg()==pdg && p->Status()==ist) nentries++;
520  }
521  return nentries;
522 }
523 //___________________________________________________________________________
524 unsigned int GHepRecord::NEntries(int pdg, int start) const
525 {
526  unsigned int nentries = 0;
527 
528  for(int i = start; i < this->GetEntries(); i++) {
529  GHepParticle * p = (GHepParticle *) (*this)[i];
530  if(p->Pdg()==pdg) nentries++;
531  }
532  return nentries;
533 }
534 //___________________________________________________________________________
536 {
537 // Provides a simplified method for inserting entries in the TClonesArray
538 
539  unsigned int pos = this->GetEntries();
540 
541 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
542  LOG("GHEP", pINFO)
543  << "Adding particle with pdgc = " << p.Pdg() << " at slot = " << pos;
544 #endif
545  new ((*this)[pos]) GHepParticle(p);
546 
547  // Update the mother's daughter list. If the newly inserted particle broke
548  // compactification, then run CompactifyDaughterLists()
549  this->UpdateDaughterLists();
550 }
551 //___________________________________________________________________________
553  int pdg, GHepStatus_t status, int mom1, int mom2, int dau1, int dau2,
554  const TLorentzVector & p, const TLorentzVector & v)
555 {
556 // Provides a 'simplified' method for inserting entries in the TClonesArray
557 
558  unsigned int pos = this->GetEntries();
559 
560 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
561  LOG("GHEP", pINFO)
562  << "Adding particle with pdgc = " << pdg << " at slot = " << pos;
563 #endif
564  new ((*this)[pos]) GHepParticle(pdg,status, mom1,mom2,dau1,dau2, p, v);
565 
566  // Update the mother's daughter list. If the newly inserted particle broke
567  // compactification, then run CompactifyDaughterLists()
568  this->UpdateDaughterLists();
569 }
570 //___________________________________________________________________________
572  int pdg, GHepStatus_t status, int mom1, int mom2, int dau1, int dau2,
573  double px, double py, double pz, double E,
574  double x, double y, double z, double t)
575 {
576 // Provides a 'simplified' method for inserting entries in the TClonesArray
577 
578  unsigned int pos = this->GetEntries();
579 
580 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
581  LOG("GHEP", pINFO)
582  << "Adding particle with pdgc = " << pdg << " at slot = " << pos;
583 #endif
584  new ( (*this)[pos] ) GHepParticle (
585  pdg, status, mom1, mom2, dau1, dau2, px, py, pz, E, x, y, z, t);
586 
587  // Update the mother's daughter list. If the newly inserted particle broke
588  // compactification, then run CompactifyDaughterLists()
589  this->UpdateDaughterLists();
590 }
591 //___________________________________________________________________________
593 {
594  int pos = this->GetEntries() - 1; // position of last entry
595 
596 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
597  LOG("GHEP", pINFO)
598  << "Updating the daughter-list for the mother of particle at: " << pos;
599 #endif
600 
601  GHepParticle * p = this->Particle(pos);
602  assert(p);
603 
604  int mom_pos = p->FirstMother();
605 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
606  LOG("GHEP", pINFO) << "Mother particle is at slot: " << mom_pos;
607 #endif
608  if(mom_pos==-1) return; // may not have mom (eg init state)
609  GHepParticle * mom = this->Particle(mom_pos);
610  if(!mom) return; // may not have mom (eg init state)
611 
612  int dau1 = mom->FirstDaughter();
613  int dau2 = mom->LastDaughter();
614 
615  // handles the case where the daughter list was initially empty
616  if(dau1 == -1) {
617  mom->SetFirstDaughter(pos);
618  mom->SetLastDaughter(pos);
619  #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
620  LOG("GHEP", pINFO)
621  << "Done! Daughter-list is compact: [" << pos << ", " << pos << "]";
622  #endif
623  return;
624  }
625  // handles the case where the new daughter is added at the slot just before
626  // an already compact daughter list
627  if(pos == dau1-1) {
628  mom->SetFirstDaughter(pos);
629  #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
630  LOG("GHEP", pINFO)
631  << "Done! Daughter-list is compact: [" << pos << ", " << dau2 << "]";
632  #endif
633  return;
634  }
635  // handles the case where the new daughter is added at the slot just after
636  // an already compact daughter list
637  if(pos == dau2+1) {
638  mom->SetLastDaughter(pos);
639  #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
640  LOG("GHEP", pINFO)
641  << "Done! Daughter-list is compact: [" << dau1 << ", " << pos << "]";
642  #endif
643  return;
644  }
645 
646  // If you are here, then the last particle insertion broke the daughter
647  // list compactification - Run the compactifier
648  LOG("GHEP", pNOTICE)
649  << "Daughter-list is not compact - Running compactifier";
650  this->CompactifyDaughterLists();
651 }
652 //___________________________________________________________________________
654 {
655  LOG("GHEP", pNOTICE) << "Removing all intermediate particles from GHEP";
656  this->Compress();
657 
658  int i=0;
659  GHepParticle * p = 0;
660 
661  TIter iter(this);
662  while( (p = (GHepParticle *)iter.Next()) ) {
663 
664  if(!p) continue;
665  GHepStatus_t ist = p->Status();
666 
667  bool keep = (ist==kIStInitialState) ||
668  (ist==kIStStableFinalState) || (ist==kIStNucleonTarget);
669  if(keep) {
670  p->SetFirstDaughter(-1);
671  p->SetLastDaughter(-1);
672  p->SetFirstMother(-1);
673  p->SetLastMother(-1);
674  } else {
675  LOG("GHEP", pNOTICE)
676  << "Removing: " << p->Name() << " from slot: " << i;
677  this->RemoveAt(i);
678  }
679  i++;
680  }
681 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
682  LOG("GHEP", pDEBUG) << "Compressing GHEP record to remove empty slots";
683 #endif
684  this->Compress();
685 }
686 //___________________________________________________________________________
688 {
689  int n = this->GetEntries();
690  if(n<1) return;
691 
692  int i = this->Particle(n-1)->FirstMother();
693  if(i<0) return;
694 
695  // for(int i=0; i<n; i++) {
696  bool compact = this->HasCompactDaughterList(i);
697 
698 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
699  LOG("GHEP", pNOTICE)
700  << "Particle's " << i << " daughter list is "
701  << ((compact) ? "compact" : "__not__ compact");
702 #endif
703 
704  if(!compact) {
705  GHepParticle * p = this->Particle(i);
706 
707  int dau1 = p->FirstDaughter();
708  int dau2 = p->LastDaughter();
709  int ndau = dau2-dau1+1;
710  int ndp = dau2+1;
711  if(dau1==-1) {ndau=0;}
712 
713  int curr_pos = n-1;
714  while (curr_pos > ndp) {
715  this->SwapParticles(curr_pos,curr_pos-1);
716  curr_pos--;
717  }
718  if(ndau>0) {
719  this->Particle(i)->SetFirstDaughter(dau1);
720  this->Particle(i)->SetLastDaughter(dau2+1);
721  } else {
722  this->Particle(i)->SetFirstDaughter(-1);
723  this->Particle(i)->SetLastDaughter(-1);
724  }
725  this->FinalizeDaughterLists();
726 
727  } //!compact
728 
729 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
730  LOG("GHEP", pINFO)
731  << "Done ompactifying daughter-list for particle at slot: " << i;
732 #endif
733  // }
734 }
735 //___________________________________________________________________________
737 {
738 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
739  LOG("GHEP", pDEBUG) << "Examining daughter-list of particle at: " << pos;
740 #endif
741  vector<int> daughters;
742  GHepParticle * p = 0;
743  TIter iter(this);
744  int i=0;
745  while( (p = (GHepParticle *)iter.Next()) ) {
746  if(p->FirstMother() == pos) {
747 
748 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
749  LOG("GHEP", pDEBUG) << "Particle at: " << i << " is a daughter";
750 #endif
751  daughters.push_back(i);
752  }
753  i++;
754  }
755 
756  bool is_compact = true;
757  if(daughters.size()>1) {
758  sort(daughters.begin(), daughters.end());
759  vector<int>::iterator diter = daughters.begin();
760  int prev = *diter;
761  for(; diter != daughters.end(); ++diter) {
762  int curr = *diter;
763  is_compact = is_compact && (TMath::Abs(prev-curr)<=1);
764  prev = curr;
765  }
766  }
767 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
768  LOG("GHEP", pINFO)
769  << "Daughter-list of particle at: " << pos << " is "
770  << (is_compact ? "" : "not") << " compact";
771 #endif
772  return is_compact;
773 }
774 //___________________________________________________________________________
776 {
777  GHepParticle * p = 0;
778  TIter iter(this);
779  int pos = 0;
780  while( (p = (GHepParticle *)iter.Next()) ) {
781  int ist = p->Status();
782  if(ist != kIStInitialState && ist != kIStNucleonTarget) return pos;
783  pos++;
784  }
785  return pos;
786 }
787 //___________________________________________________________________________
789 {
790 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
791  LOG("GHEP", pINFO) << "Swapping GHepParticles : " << i << " <--> " << j;
792 #endif
793  int n = this->GetEntries();
794  assert(i>=0 && j>=0 && i<n && j<n);
795 
796  if(i==j) return;
797 
798  GHepParticle * pi = this->Particle(i);
799  GHepParticle * pj = this->Particle(j);
800  GHepParticle * tmp = new GHepParticle(*pi);
801 
802  pi->Copy(*pj);
803  pj->Copy(*tmp);
804 
805  delete tmp;
806 
807  // tell their daughters
808  if(pi->HasDaughters()) {
809 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
810  LOG("GHEP", pINFO)
811  << pi->Name() << "(previously at pos: " << j
812  << ") is now at pos: " << i << " -> Notify daughters";
813 #endif
814  for(int k=0; k<n; k++) {
815  if(this->Particle(k)->FirstMother()==j) {
816  this->Particle(k)->SetFirstMother(i);
817  }
818  }
819  }
820 
821  if(pj->HasDaughters()) {
822 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
823  LOG("GHEP", pINFO)
824  << pj->Name() << "(previously at pos: " << i
825  << ") is now at pos: " << j << " -> Notify daughters";
826 #endif
827  for(int k=0; k<n; k++) {
828  if(this->Particle(k)->FirstMother()==i) {
829  this->Particle(k)->SetFirstMother(j);
830  }
831  }
832  }
833 }
834 //___________________________________________________________________________
836 {
837 // Update all daughter-lists based on particle 'first mother' field.
838 // To work correctly, the daughter-lists must have been compactified first.
839 
840  GHepParticle * p1 = 0;
841  TIter iter1(this);
842  int i1=0;
843  while( (p1 = (GHepParticle *)iter1.Next()) ) {
844  int dau1 = -1;
845  int dau2 = -1;
846  GHepParticle * p2 = 0;
847  TIter iter2(this);
848  int i2=0;
849  while( (p2 = (GHepParticle *)iter2.Next()) ) {
850 
851  if(p2->FirstMother() == i1) {
852  dau1 = (dau1<0) ? i2 : TMath::Min(dau1,i2);
853  dau2 = (dau2<0) ? i2 : TMath::Max(dau2,i2);
854  }
855  i2++;
856  }
857  i1++;
858  p1 -> SetFirstDaughter (dau1);
859  p1 -> SetLastDaughter (dau2);
860  }
861 }
862 //___________________________________________________________________________
863 void GHepRecord::SetVertex(double x, double y, double z, double t)
864 {
865  fVtx->SetXYZT(x,y,z,t);
866 }
867 //___________________________________________________________________________
868 void GHepRecord::SetVertex(const TLorentzVector & vtx)
869 {
870  fVtx->SetXYZT(vtx.X(),vtx.Y(),vtx.Z(),vtx.T());
871 }
872 //___________________________________________________________________________
874 {
875 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
876  LOG("GHEP", pDEBUG) << "Initializing GHepRecord";
877 #endif
878  fInteraction = 0;
879  fWeight = 1.;
880  fProb = 1.;
881  fXSec = 0.;
882  fDiffXSec = 0.;
884  fVtx = new TLorentzVector(0,0,0,0);
885 
886  fEventFlags = new TBits(GHepFlags::NFlags());
887  fEventFlags -> ResetAllBits(false);
888 
889  fEventMask = new TBits(GHepFlags::NFlags());
890 //fEventMask -> ResetAllBits(true);
891  for(unsigned int i = 0; i < GHepFlags::NFlags(); i++) {
892  fEventMask->SetBitNumber(i, true);
893  }
894 
895  LOG("GHEP", pINFO)
896  << "Initialised unphysical event mask (bits: " << GHepFlags::NFlags() - 1
897  << " -> 0) : " << *fEventMask;
898 
899  this->SetOwner(true);
900 }
901 //___________________________________________________________________________
903 {
904 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
905  LOG("GHEP", pDEBUG) << "Cleaning up GHepRecord";
906 #endif
907  this->Clear("C");
908 }
909 //___________________________________________________________________________
911 {
912 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
913  LOG("GHEP", pDEBUG) << "Reseting GHepRecord";
914 #endif
915  this->CleanRecord();
916  this->InitRecord();
917 }
918 //___________________________________________________________________________
919 void GHepRecord::Clear(Option_t * opt)
920 {
921  if (fInteraction) delete fInteraction;
922  fInteraction=0;
923 
924  if (fVtx) delete fVtx;
925  fVtx=0;
926 
927  if(fEventFlags) delete fEventFlags;
928  fEventFlags=0;
929 
930  if(fEventMask) delete fEventMask;
931  fEventMask=0;
932 
933  TClonesArray::Clear(opt);
934 
935 // if (fInteraction) delete fInteraction;
936 // delete fVtx;
937 //
938 // delete fEventFlags;
939 //
940 // TClonesArray::Clear(opt);
941 }
942 //___________________________________________________________________________
943 void GHepRecord::Copy(const GHepRecord & record)
944 {
945  // clean up
946  this->ResetRecord();
947 
948  // copy event record entries
949  unsigned int ientry = 0;
950  GHepParticle * p = 0;
951  TIter ghepiter(&record);
952  while ( (p = (GHepParticle *) ghepiter.Next()) )
953  new ( (*this)[ientry++] ) GHepParticle(*p);
954 
955  // copy summary
956  fInteraction = new Interaction( *record.fInteraction );
957 
958  // copy flags & mask
959  *fEventFlags = *(record.EventFlags());
960  *fEventMask = *(record.EventMask());
961 
962  // copy vtx position
963  TLorentzVector * v = record.Vertex();
964  fVtx->SetXYZT(v->X(),v->Y(),v->Z(),v->T());
965 
966  // copy weights & xsecs
967  fWeight = record.fWeight;
968  fProb = record.fProb;
969  fXSec = record.fXSec;
970  fDiffXSec = record.fDiffXSec;
971  fDiffXSecPhSp = record.fDiffXSecPhSp;
972 }
973 //___________________________________________________________________________
975 {
976  *fEventMask = mask;
977 
978  LOG("GHEP", pINFO)
979  << "Setting unphysical event mask (bits: " << GHepFlags::NFlags() - 1
980  << " -> 0) : " << *fEventMask;
981 }
982 //___________________________________________________________________________
983 bool GHepRecord::Accept(void) const
984 {
985  TBits flags = *fEventFlags;
986  TBits mask = *fEventMask;
987  TBits bitwiseand = flags & mask;
988  bool accept = (bitwiseand.CountBits() == 0);
989  return accept;
990 }
991 //___________________________________________________________________________
992 void GHepRecord::SetPrintLevel(int print_level)
993 {
994  fPrintLevel = print_level;
995 }
997 {
998  return fPrintLevel;
999 }
1000 //___________________________________________________________________________
1001 void GHepRecord::Print(ostream & stream) const
1002 {
1003  // Print levels:
1004  // 0 -> prints particle list
1005  // 1 -> prints particle list + event flags
1006  // 2 -> prints particle list + event flags + wght/xsec
1007  // 3 -> prints particle list + event flags + wght/xsec + summary
1008  // 10 -> as in level 0 but showing particle positions too
1009  // 11 -> as in level 1 but showing particle positions too
1010  // 12 -> as in level 2 but showing particle positions too
1011  // 13 -> as in level 3 but showing particle positions too
1012 
1013  bool accept_input_print_level =
1014  (fPrintLevel >= 0 && fPrintLevel <= 3) ||
1015  (fPrintLevel >=10 && fPrintLevel <=13);
1016 
1017  int printlevel = (accept_input_print_level) ? fPrintLevel : 3;
1018  int printlevel_orig = printlevel;
1019 
1020  bool showpos = false;
1021  if(printlevel >= 10) {
1022  printlevel-=10;
1023  showpos=true;
1024  }
1025 
1026  stream << "\n\n|";
1027  stream << setfill('-') << setw(115) << "|";
1028 
1029  stream << "\n|GENIE GHEP Event Record [print level: "
1030  << setfill(' ') << setw(3) << printlevel_orig << "]"
1031  << setfill(' ') << setw(73) << "|";
1032 
1033  stream << "\n|";
1034  stream << setfill('-') << setw(115) << "|";
1035 
1036  stream << "\n| ";
1037  stream << setfill(' ') << setw(6) << "Idx | "
1038  << setfill(' ') << setw(16) << "Name | "
1039  << setfill(' ') << setw(6) << "Ist | "
1040  << setfill(' ') << setw(13) << "PDG | "
1041  << setfill(' ') << setw(12) << "Mother | "
1042  << setfill(' ') << setw(12) << "Daughter | "
1043  << setfill(' ') << setw(10) << ((showpos) ? "Px(x) |" : "Px | ")
1044  << setfill(' ') << setw(10) << ((showpos) ? "Py(y) |" : "Py | ")
1045  << setfill(' ') << setw(10) << ((showpos) ? "Pz(z) |" : "Pz | ")
1046  << setfill(' ') << setw(10) << ((showpos) ? "E(t) |" : "E | ")
1047  << setfill(' ') << setw(10) << "m | ";
1048 
1049  stream << "\n|";
1050  stream << setfill('-') << setw(115) << "|";
1051 
1052  GHepParticle * p = 0;
1053  TObjArrayIter piter(this);
1054  TVector3 polarization(0,0,0);
1055 
1056  unsigned int idx = 0;
1057 
1058  double sum_E = 0;
1059  double sum_px = 0;
1060  double sum_py = 0;
1061  double sum_pz = 0;
1062 
1063  while( (p = (GHepParticle *) piter.Next()) ) {
1064 
1065  stream << "\n| ";
1066  stream << setfill(' ') << setw(3) << idx++ << " | ";
1067  stream << setfill(' ') << setw(13) << p->Name() << " | ";
1068  stream << setfill(' ') << setw(3) << p->Status() << " | ";
1069  stream << setfill(' ') << setw(10) << p->Pdg() << " | ";
1070  stream << setfill(' ') << setw(3) << p->FirstMother() << " | ";
1071  stream << setfill(' ') << setw(3) << p->LastMother() << " | ";
1072  stream << setfill(' ') << setw(3) << p->FirstDaughter() << " | ";
1073  stream << setfill(' ') << setw(3) << p->LastDaughter() << " | ";
1074  stream << std::fixed << setprecision(3);
1075  stream << setfill(' ') << setw(7) << p->Px() << " | ";
1076  stream << setfill(' ') << setw(7) << p->Py() << " | ";
1077  stream << setfill(' ') << setw(7) << p->Pz() << " | ";
1078  stream << setfill(' ') << setw(7) << p->E() << " | ";
1079 
1080  if (p->IsOnMassShell())
1081  stream << setfill(' ') << setw(7) << p->Mass() << " | ";
1082  else
1083  stream << setfill('*') << setw(7) << p->Mass() << " | M = "
1084  << p->P4()->M() << " ";
1085 
1086  if (p->PolzIsSet()) {
1087  p->GetPolarization(polarization);
1088  stream << "P = (" << polarization.x() << "," << polarization.y()
1089  << "," << polarization.z() << ")";
1090  }
1091 
1092  if (p->RescatterCode() != -1) {
1093  stream << "FSI = " << p->RescatterCode();
1094  }
1095 
1096  // plot particle position if requested
1097  if(showpos) {
1098  stream << "\n| ";
1099  stream << setfill(' ') << setw(6) << " | ";
1100  stream << setfill(' ') << setw(16) << " | ";
1101  stream << setfill(' ') << setw(6) << " | ";
1102  stream << setfill(' ') << setw(13) << " | ";
1103  stream << setfill(' ') << setw(6) << " | ";
1104  stream << setfill(' ') << setw(6) << " | ";
1105  stream << setfill(' ') << setw(6) << " | ";
1106  stream << setfill(' ') << setw(6) << " | ";
1107  stream << std::fixed << setprecision(3);
1108  stream << setfill(' ') << setw(7) << p->Vx() << " | ";
1109  stream << setfill(' ') << setw(7) << p->Vy() << " | ";
1110  stream << setfill(' ') << setw(7) << p->Vz() << " | ";
1111  stream << setfill(' ') << setw(7) << p->Vt() << " | ";
1112  stream << setfill(' ') << setw(10) << " | ";
1113  }
1114 
1115  // compute P4_{final} - P4_{nitial}
1116 
1117  if(p->Status() == kIStStableFinalState ||
1119 
1120  sum_E += p->E();
1121  sum_px += p->Px();
1122  sum_py += p->Py();
1123  sum_pz += p->Pz();
1124  } else
1125  if(p->Status() == kIStInitialState) {
1126  /*
1127  if(p->Status() == kIStInitialState || p->Status() == kIStNucleonTarget) {
1128  */
1129  sum_E -= p->E();
1130  sum_px -= p->Px();
1131  sum_py -= p->Py();
1132  sum_pz -= p->Pz();
1133  }
1134 
1135  } // loop over particles
1136 
1137  stream << "\n|";
1138  stream << setfill('-') << setw(115) << "|";
1139 
1140  // Print SUMS
1141  stream << "\n| ";
1142  stream << setfill(' ') << setw(17) << "Fin-Init: "
1143  << setfill(' ') << setw(6) << " "
1144  << setfill(' ') << setw(18) << " "
1145  << setfill(' ') << setw(12) << " "
1146  << setfill(' ') << setw(12) << " | ";
1147  stream << std::fixed << setprecision(3);
1148  stream << setfill(' ') << setw(7) << sum_px << " | ";
1149  stream << setfill(' ') << setw(7) << sum_py << " | ";
1150  stream << setfill(' ') << setw(7) << sum_pz << " | ";
1151  stream << setfill(' ') << setw(7) << sum_E << " | ";
1152  stream << setfill(' ') << setw(10) << " | ";
1153 
1154  stream << "\n|";
1155  stream << setfill('-') << setw(115) << "|";
1156 
1157  // Print vertex
1158 
1159  GHepParticle * probe = this->Probe();
1160  if(probe){
1161  stream << "\n| ";
1162  stream << setfill(' ') << setw(17) << "Vertex: ";
1163  stream << setfill(' ') << setw(11)
1164  << ((probe) ? probe->Name() : "unknown probe") << " @ (";
1165 
1166  stream << std::fixed << setprecision(5);
1167  stream << "x = " << setfill(' ') << setw(11) << this->Vertex()->X() << " m, ";
1168  stream << "y = " << setfill(' ') << setw(11) << this->Vertex()->Y() << " m, ";
1169  stream << "z = " << setfill(' ') << setw(11) << this->Vertex()->Z() << " m, ";
1170  stream << std::scientific << setprecision(6);
1171  stream << "t = " << setfill(' ') << setw(15) << this->Vertex()->T() << " s) ";
1172  stream << std::fixed << setprecision(3);
1173  stream << setfill(' ') << setw(2) << "|";
1174 
1175  stream << "\n|";
1176  stream << setfill('-') << setw(115) << "|";
1177  }
1178 
1179  // Print FLAGS
1180 
1181  if(printlevel>=1) {
1182  stream << "\n| ";
1183  stream << "Err flag [bits:" << fEventFlags->GetNbits()-1 << "->0] : "
1184  << *fEventFlags << " | "
1185  << "1st set: " << setfill(' ') << setw(56)
1186  << ( this->IsUnphysical() ?
1187  GHepFlags::Describe(GHepFlag_t(fEventFlags->FirstSetBit())) :
1188  "none") << " | ";
1189  stream << "\n| ";
1190  stream << "Err mask [bits:" << fEventMask->GetNbits()-1 << "->0] : "
1191  << *fEventMask << " | "
1192  << "Is unphysical: " << setfill(' ') << setw(5)
1193  << utils::print::BoolAsYNString(this->IsUnphysical()) << " | "
1194  << "Accepted: " << setfill(' ') << setw(5)
1196  << " |";
1197  stream << "\n|";
1198  stream << setfill('-') << setw(115) << "|";
1199  }
1200 
1201  if(printlevel>=2) {
1202  stream << "\n| ";
1203  stream << std::scientific << setprecision(5);
1204 
1205  stream << "sig(Ev) = "
1206  << setfill(' ') << setw(17) << fXSec/units::cm2
1207  << " cm^2 |";
1208 
1209  switch(fDiffXSecPhSp) {
1210  case ( kPSyfE ) :
1211  stream << " dsig(y;E)/dy = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2 |";
1212  break;
1213  case ( kPSxyfE ) :
1214  stream << " d2sig(x,y;E)/dxdy = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2 |";
1215  break;
1216  case ( kPSxytfE ) :
1217  stream << " d3sig(x,y,t;E)/dxdydt = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2/GeV^2 |";
1218  break;
1219  case ( kPSQ2fE ) :
1220  stream << " dsig(Q2;E)/dQ2 = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2/GeV^2 |";
1221  break;
1222  case ( kPSQ2vfE ) :
1223  stream << " dsig(Q2,v;E)/dQ2dv = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2/GeV^3 |";
1224  break;
1225  case ( kPSWQ2fE ) :
1226  stream << " d2sig(W,Q2;E)/dWdQ2 = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2/GeV^3 |";
1227  break;
1228  default :
1229  stream << " dsig(Ev;{K_s})/dK = " << setfill(' ') << setw(13) << fDiffXSec/units::cm2 << " cm^2/{K} |";
1230  }
1231  stream << " Weight = "
1232  << setfill(' ') << setw(16)
1233  << std::fixed << setprecision(5)
1234  << fWeight
1235  << " |";
1236 
1237  stream << "\n|";
1238  stream << setfill('-') << setw(115) << "|";
1239  }
1240 
1241  stream << "\n";
1242  stream << setfill(' ');
1243 
1244  if(printlevel==3) {
1245  if(fInteraction) stream << *fInteraction;
1246  else stream << "NULL Interaction!" << endl;
1247  }
1248  stream << "\n";
1249 }
1250 //___________________________________________________________________________
static void SetPrintLevel(int print_level)
Definition: GHepRecord.cxx:992
void SetFirstMother(int m)
Definition: GHepParticle.h:133
int RescatterCode(void) const
Definition: GHepParticle.h:66
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int ParticlePosition(int pdg, GHepStatus_t i, int start=0) const
Definition: GHepRecord.cxx:181
virtual void Copy(const GHepRecord &record)
Definition: GHepRecord.cxx:943
virtual void UpdateDaughterLists(void)
Definition: GHepRecord.cxx:592
double fProb
event probability (for given flux neutrino and density-weighted path-length for target element) ...
Definition: GHepRecord.h:180
void Print(ostream &stream) const
int status
Definition: fabricate.py:1613
bool IsOnMassShell(void) const
virtual void SwapParticles(int i, int j)
Definition: GHepRecord.cxx:788
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
double E(void) const
Get energy.
Definition: GHepParticle.h:92
::xsd::cxx::tree::flags flags
Definition: Database.h:214
virtual GHepParticle * HitElectron(void) const
Definition: GHepRecord.cxx:360
virtual Interaction * Summary(void) const
Definition: GHepRecord.cxx:135
virtual GHepParticle * FindParticle(int pdg, GHepStatus_t ist, int start) const
Definition: GHepRecord.cxx:162
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
int FirstDaughter(void) const
Definition: GHepParticle.h:69
bool IsNucleon(int pdgc)
Definition: PDGUtils.cxx:309
virtual int RemnantNucleusPosition(void) const
Definition: GHepRecord.cxx:433
enum genie::EGHepStatus GHepStatus_t
string BoolAsYNString(bool b)
Definition: PrintUtils.cxx:115
const char * p
Definition: xmltok.h:285
enum genie::EGHepFlag GHepFlag_t
virtual void RemoveIntermediateParticles(void)
Definition: GHepRecord.cxx:653
void InitRecord(void)
Definition: GHepRecord.cxx:873
virtual void FinalizeDaughterLists(void)
Definition: GHepRecord.cxx:835
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
virtual int HitNucleonPosition(void) const
Definition: GHepRecord.cxx:454
bool IsDarkMatter(int pdgc)
Definition: PDGUtils.cxx:125
caf::StandardRecord * rec
Definition: tutCAFMacro.C:20
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
virtual unsigned int NEntries(int pdg, GHepStatus_t ist, int start=0) const
Definition: GHepRecord.cxx:513
Float_t tmp
Definition: plot.C:36
double Mass(void) const
Mass that corresponds to the PDG code.
static const double cm2
Definition: Units.h:77
cout<< t1-> GetEntries()
Definition: plottest35.C:29
double Pz(void) const
Get Pz.
Definition: GHepParticle.h:91
static const char * Describe(GHepFlag_t flag)
Definition: GHepFlags.h:43
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
virtual void AttachSummary(Interaction *interaction)
Definition: GHepRecord.cxx:143
virtual ~GHepRecord()
Definition: GHepRecord.cxx:130
virtual TLorentzVector * Vertex(void) const
Definition: GHepRecord.h:141
double Px(void) const
Get Px.
Definition: GHepParticle.h:89
bool Compare(const GHepParticle *p) const
virtual GHepParticle * Probe(void) const
Definition: GHepRecord.cxx:321
int LastMother(void) const
Definition: GHepParticle.h:68
double Vt(void) const
Get production time.
Definition: GHepParticle.h:98
int Pdg(void) const
Definition: GHepParticle.h:64
int FirstMother(void) const
Definition: GHepParticle.h:67
string Name(void) const
Name that corresponds to the PDG code.
bool IsHadron(int pdgc)
Definition: PDGUtils.cxx:355
Enumeration of GENIE event generation modes.
double fXSec
cross section for selected event
Definition: GHepRecord.h:181
Summary information for an interaction.
Definition: Interaction.h:56
Interaction * fInteraction
attached summary information
Definition: GHepRecord.h:169
int LastDaughter(void) const
Definition: GHepParticle.h:70
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
static unsigned int NFlags(void)
Definition: GHepFlags.h:77
virtual bool Accept(void) const
Definition: GHepRecord.cxx:983
void prev()
Definition: show_event.C:91
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
void SetLastDaughter(int d)
Definition: GHepParticle.h:136
Float_t E
Definition: plot.C:20
Long64_t nentries
bool HasDaughters(void) const
Definition: GHepParticle.h:71
virtual bool HasCompactDaughterList(int pos)
Definition: GHepRecord.cxx:736
void Copy(const GHepParticle &particle)
const int kPdgGamma
Definition: PDGCodes.h:166
virtual int FirstNonInitStateEntry(void)
Definition: GHepRecord.cxx:775
virtual GHepParticle * FinalStatePrimaryLepton(void) const
Definition: GHepRecord.cxx:370
GEvGenMode_t EventGenerationMode(void) const
Definition: GHepRecord.cxx:253
ClassImp(EventRecord) namespace genie
Definition: EventRecord.cxx:25
static int GetPrintLevel()
Definition: GHepRecord.cxx:996
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
virtual vector< int > * GetStableDescendants(int position) const
Definition: GHepRecord.cxx:218
const double j
Definition: BetheBloch.cxx:29
#define pINFO
Definition: Messenger.h:63
virtual void ResetRecord(void)
Definition: GHepRecord.cxx:910
virtual int HitElectronPosition(void) const
Definition: GHepRecord.cxx:477
virtual GHepParticle * FinalStateHadronicSystem(void) const
Definition: GHepRecord.cxx:379
z
Definition: test.py:28
bool PolzIsSet(void) const
#define pWARN
Definition: Messenger.h:61
void GetPolarization(TVector3 &polz)
bool Is2NucleonCluster(int pdgc)
Definition: PDGUtils.cxx:365
virtual void CompactifyDaughterLists(void)
Definition: GHepRecord.cxx:687
virtual TBits * EventFlags(void) const
Definition: GHepRecord.h:118
void SetLastMother(int m)
Definition: GHepParticle.h:134
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
virtual bool IsUnphysical(void) const
Definition: GHepRecord.h:120
double Vz(void) const
Get production z.
Definition: GHepParticle.h:97
ostream & operator<<(ostream &stream, const AlgConfigPool &config_pool)
virtual void Clear(Option_t *opt="")
Definition: GHepRecord.cxx:919
virtual GHepParticle * HitNucleon(void) const
Definition: GHepRecord.cxx:350
virtual void SetVertex(double x, double y, double z, double t)
Definition: GHepRecord.cxx:863
KinePhaseSpace_t fDiffXSecPhSp
specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)
Definition: GHepRecord.h:183
virtual GHepParticle * RemnantNucleus(void) const
Definition: GHepRecord.cxx:340
assert(nhit_max >=nhit_nbins)
virtual void AddParticle(const GHepParticle &p)
Definition: GHepRecord.cxx:535
double Vy(void) const
Get production y.
Definition: GHepParticle.h:96
#define pNOTICE
Definition: Messenger.h:62
double fWeight
event weight
Definition: GHepRecord.h:179
virtual int ProbePosition(void) const
Definition: GHepRecord.cxx:389
void SetFirstDaughter(int d)
Definition: GHepParticle.h:135
const int kPdgHadronicSyst
Definition: PDGCodes.h:187
void SetUnphysEventMask(const TBits &mask)
Definition: GHepRecord.cxx:974
GENIE&#39;s GHEP MC event record.
Definition: GHepRecord.h:46
Most commonly used PDG codes. A set of utility functions to handle PDG codes is provided in PDGUtils...
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
virtual int FinalStateHadronicSystemPosition(void) const
Definition: GHepRecord.cxx:507
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
bool IsLepton(int pdgc)
Definition: PDGUtils.cxx:84
virtual int FinalStatePrimaryLeptonPosition(void) const
Definition: GHepRecord.cxx:495
virtual int TargetNucleusPosition(void) const
Definition: GHepRecord.cxx:406
static int fPrintLevel
Definition: GHepRecord.h:197
bool IsElectron(int pdgc)
Definition: PDGUtils.cxx:180
virtual TBits * EventMask(void) const
Definition: GHepRecord.h:119
double Vx(void) const
Get production x.
Definition: GHepParticle.h:95
#define pDEBUG
Definition: Messenger.h:64
void CleanRecord(void)
Definition: GHepRecord.cxx:902
double Py(void) const
Get Py.
Definition: GHepParticle.h:90