Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | Friends | List of all members
genie::GHepRecord Class Reference

GENIE's GHEP MC event record. More...

#include "/cvmfs/nova.opensciencegrid.org/externals/genie/v3_00_06_p01/Linux64bit+2.6-2.12-e17-debug/GENIE-Generator/src/Framework/GHEP/GHepRecord.h"

Inheritance diagram for genie::GHepRecord:
genie::EventRecord

Public Member Functions

 GHepRecord ()
 
 GHepRecord (int size)
 
 GHepRecord (const GHepRecord &record)
 
 GHepRecord (TRootIOCtor *)
 
virtual ~GHepRecord ()
 
virtual InteractionSummary (void) const
 
virtual void AttachSummary (Interaction *interaction)
 
virtual void AddParticle (const GHepParticle &p)
 
virtual void AddParticle (int pdg, GHepStatus_t ist, int mom1, int mom2, int dau1, int dau2, const TLorentzVector &p, const TLorentzVector &v)
 
virtual void AddParticle (int pdg, GHepStatus_t ist, int mom1, int mom2, int dau1, int dau2, double px, double py, double pz, double E, double x, double y, double z, double t)
 
virtual GHepParticleParticle (int position) const
 
virtual GHepParticleFindParticle (int pdg, GHepStatus_t ist, int start) const
 
virtual int ParticlePosition (int pdg, GHepStatus_t i, int start=0) const
 
virtual int ParticlePosition (GHepParticle *particle, int start=0) const
 
virtual vector< int > * GetStableDescendants (int position) const
 
GEvGenMode_t EventGenerationMode (void) const
 
virtual GHepParticleProbe (void) const
 
virtual GHepParticleTargetNucleus (void) const
 
virtual GHepParticleRemnantNucleus (void) const
 
virtual GHepParticleHitNucleon (void) const
 
virtual GHepParticleHitElectron (void) const
 
virtual GHepParticleFinalStatePrimaryLepton (void) const
 
virtual GHepParticleFinalStateHadronicSystem (void) const
 
virtual int ProbePosition (void) const
 
virtual int TargetNucleusPosition (void) const
 
virtual int RemnantNucleusPosition (void) const
 
virtual int HitNucleonPosition (void) const
 
virtual int HitElectronPosition (void) const
 
virtual int FinalStatePrimaryLeptonPosition (void) const
 
virtual int FinalStateHadronicSystemPosition (void) const
 
virtual unsigned int NEntries (int pdg, GHepStatus_t ist, int start=0) const
 
virtual unsigned int NEntries (int pdg, int start=0) const
 
virtual TBits * EventFlags (void) const
 
virtual TBits * EventMask (void) const
 
virtual bool IsUnphysical (void) const
 
virtual bool Accept (void) const
 
virtual double Weight (void) const
 
virtual double Probability (void) const
 
virtual double XSec (void) const
 
virtual double DiffXSec (void) const
 
virtual KinePhaseSpace_t DiffXSecVars (void) const
 
virtual void SetWeight (double wght)
 
virtual void SetProbability (double prob)
 
virtual void SetXSec (double xsec)
 
virtual void SetDiffXSec (double xsec, KinePhaseSpace_t ps)
 
virtual TLorentzVector * Vertex (void) const
 
virtual void SetVertex (double x, double y, double z, double t)
 
virtual void SetVertex (const TLorentzVector &vtx)
 
virtual void Copy (const GHepRecord &record)
 
virtual void Clear (Option_t *opt="")
 
virtual void ResetRecord (void)
 
virtual void CompactifyDaughterLists (void)
 
virtual void RemoveIntermediateParticles (void)
 
void SetUnphysEventMask (const TBits &mask)
 
void Print (ostream &stream) const
 

Static Public Member Functions

static void SetPrintLevel (int print_level)
 
static int GetPrintLevel ()
 

Protected Member Functions

void InitRecord (void)
 
void CleanRecord (void)
 
virtual void UpdateDaughterLists (void)
 
virtual bool HasCompactDaughterList (int pos)
 
virtual void SwapParticles (int i, int j)
 
virtual void FinalizeDaughterLists (void)
 
virtual int FirstNonInitStateEntry (void)
 

Protected Attributes

InteractionfInteraction
 attached summary information More...
 
TLorentzVector * fVtx
 vertex in the detector coordinate system More...
 
TBits * fEventFlags
 event flags indicating various pathologies or an unphysical event More...
 
TBits * fEventMask
 an input bit-field mask allowing one to ignore bits set in fEventFlags More...
 
double fWeight
 event weight More...
 
double fProb
 event probability (for given flux neutrino and density-weighted path-length for target element) More...
 
double fXSec
 cross section for selected event More...
 
double fDiffXSec
 differential cross section for selected event kinematics More...
 
KinePhaseSpace_t fDiffXSecPhSp
 specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...) More...
 

Static Protected Attributes

static int fPrintLevel
 

Friends

ostream & operator<< (ostream &stream, const GHepRecord &event)
 

Detailed Description

GENIE's GHEP MC event record.

Author
Costas Andreopoulos <costas.andreopoulos stfc.ac.uk> University of Liverpool & STFC Rutherford Appleton Lab

October 1, 2004

Copyright (c) 2003-2019, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org or see $GENIE/LICENSE

Definition at line 46 of file GHepRecord.h.

Constructor & Destructor Documentation

GHepRecord::GHepRecord ( )

Definition at line 97 of file GHepRecord.cxx.

References InitRecord().

97  :
98 TClonesArray("genie::GHepParticle")
99 {
100  this->InitRecord();
101 }
void InitRecord(void)
Definition: GHepRecord.cxx:873
GHepRecord::GHepRecord ( int  size)

Definition at line 103 of file GHepRecord.cxx.

References InitRecord().

103  :
104 TClonesArray("genie::GHepParticle", size)
105 {
106  this->InitRecord();
107 }
void InitRecord(void)
Definition: GHepRecord.cxx:873
GHepRecord::GHepRecord ( const GHepRecord record)

Definition at line 109 of file GHepRecord.cxx.

References Copy(), and InitRecord().

109  :
110 TClonesArray("genie::GHepParticle", record.GetEntries())
111 {
112  this->InitRecord();
113  this->Copy(record);
114 }
virtual void Copy(const GHepRecord &record)
Definition: GHepRecord.cxx:943
void InitRecord(void)
Definition: GHepRecord.cxx:873
GHepRecord::GHepRecord ( TRootIOCtor *  )

Definition at line 116 of file GHepRecord.cxx.

116  :
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 }
double fProb
event probability (for given flux neutrino and density-weighted path-length for target element) ...
Definition: GHepRecord.h:180
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
double fXSec
cross section for selected event
Definition: GHepRecord.h:181
Interaction * fInteraction
attached summary information
Definition: GHepRecord.h:169
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
double fWeight
event weight
Definition: GHepRecord.h:179
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
GHepRecord::~GHepRecord ( )
virtual

Definition at line 130 of file GHepRecord.cxx.

References CleanRecord().

131 {
132  this->CleanRecord();
133 }
void CleanRecord(void)
Definition: GHepRecord.cxx:902

Member Function Documentation

bool GHepRecord::Accept ( void  ) const
virtual

Definition at line 983 of file GHepRecord.cxx.

References fEventFlags, fEventMask, and lem_server::mask.

Referenced by IsUnphysical(), Print(), and genie::EventGenerator::ProcessEventRecord().

984 {
985  TBits flags = *fEventFlags;
986  TBits mask = *fEventMask;
987  TBits bitwiseand = flags & mask;
988  bool accept = (bitwiseand.CountBits() == 0);
989  return accept;
990 }
::xsd::cxx::tree::flags flags
Definition: Database.h:214
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
void GHepRecord::AddParticle ( const GHepParticle p)
virtual

Definition at line 535 of file GHepRecord.cxx.

References GetEntries(), LOG, genie::GHepParticle::Pdg(), pINFO, elec2geo::pos, and UpdateDaughterLists().

Referenced by genie::HNIntranuke2018::AbsorbHN(), genie::NuETargetRemnantGenerator::AddElectronNeutrino(), genie::HadronicSystemGenerator::AddFinalHadronicSyst(), genie::AMNuGammaGenerator::AddFinalStateNeutrino(), genie::DISHadronicSystemGenerator::AddFragmentationProducts(), genie::InitialStateAppender::AddNeutrino(), genie::InitialStateAppender::AddNucleus(), genie::AMNuGammaGenerator::AddPhoton(), genie::NucDeExcitationSim::AddPhoton(), genie::QELHadronicSystemGenerator::AddRecoilBaryon(), genie::IBDHadronicSystemGenerator::AddRecoilBaryon(), genie::AMNuGammaGenerator::AddRecoilNucleon(), genie::RSPPResonanceSelector::AddResonance(), genie::RESHadronicSystemGenerator::AddResonance(), genie::RSPPHadronicSystemGenerator::AddResonanceDecayProducts(), genie::InitialStateAppender::AddStruckParticle(), genie::NuETargetRemnantGenerator::AddTargetNucleusRemnant(), genie::HadronicSystemGenerator::AddTargetNucleusRemnant(), genie::FermiMover::AddTargetNucleusRemnant(), genie::QELEventGenerator::AddTargetNucleusRemnant(), genie::QELEventGeneratorSM::AddTargetNucleusRemnant(), genie::AMNuGammaGenerator::AddTargetRemnant(), genie::PrimaryLeptonGenerator::AddToEventRecord(), genie::OutgoingDarkGenerator::AddToEventRecord(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_AlvarezRuso(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_BergerSehgalFM(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_ReinSehgal(), ConvertToGHepMock(), genie::HAIntranuke::ElasHA(), genie::HAIntranuke2018::ElasHA(), genie::HNIntranuke2018::ElasHN(), genie::FermiMover::Emit2ndNucleonFromSRC(), genie::HNIntranuke2018::GammaInelasticHN(), genie::HNIntranuke2018::HandleCompoundNucleus(), genie::HAIntranuke::Inelastic(), genie::HAIntranuke2018::Inelastic(), genie::HAIntranuke::InelasticHA(), genie::HAIntranuke2018::InelasticHA(), genie::HNIntranuke2018::InelasticHN(), InitializeEvent(), genie::utils::intranuke::PhaseSpaceDecay(), genie::utils::intranuke2018::PhaseSpaceDecay(), genie::utils::intranuke::PionProduction(), genie::utils::intranuke2018::PionProduction(), genie::COHElHadronicSystemGenerator::ProcessEventRecord(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::NucBindEnergyAggregator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), genie::HAIntranuke::SimulateHadronicFinalState(), genie::HAIntranuke2018::SimulateHadronicFinalState(), genie::HNIntranuke2018::SimulateHadronicFinalState(), genie::Intranuke::TransportHadrons(), genie::Intranuke2018::TransportHadrons(), and genie::HadronTransporter::TransportInTransparentNuc().

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 }
virtual void UpdateDaughterLists(void)
Definition: GHepRecord.cxx:592
cout<< t1-> GetEntries()
Definition: plottest35.C:29
int Pdg(void) const
Definition: GHepParticle.h:64
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::AddParticle ( int  pdg,
GHepStatus_t  ist,
int  mom1,
int  mom2,
int  dau1,
int  dau2,
const TLorentzVector &  p,
const TLorentzVector &  v 
)
virtual

Definition at line 552 of file GHepRecord.cxx.

References GetEntries(), LOG, pINFO, elec2geo::pos, and UpdateDaughterLists().

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 }
virtual void UpdateDaughterLists(void)
Definition: GHepRecord.cxx:592
int status
Definition: fabricate.py:1613
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::AddParticle ( int  pdg,
GHepStatus_t  ist,
int  mom1,
int  mom2,
int  dau1,
int  dau2,
double  px,
double  py,
double  pz,
double  E,
double  x,
double  y,
double  z,
double  t 
)
virtual

Definition at line 571 of file GHepRecord.cxx.

References GetEntries(), LOG, pINFO, elec2geo::pos, and UpdateDaughterLists().

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 }
virtual void UpdateDaughterLists(void)
Definition: GHepRecord.cxx:592
int status
Definition: fabricate.py:1613
cout<< t1-> GetEntries()
Definition: plottest35.C:29
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
Float_t E
Definition: plot.C:20
#define pINFO
Definition: Messenger.h:63
z
Definition: test.py:28
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::AttachSummary ( Interaction interaction)
virtual
void GHepRecord::CleanRecord ( void  )
protected

Definition at line 902 of file GHepRecord.cxx.

References Clear(), LOG, and pDEBUG.

Referenced by ResetRecord(), and ~GHepRecord().

903 {
904 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
905  LOG("GHEP", pDEBUG) << "Cleaning up GHepRecord";
906 #endif
907  this->Clear("C");
908 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
virtual void Clear(Option_t *opt="")
Definition: GHepRecord.cxx:919
#define pDEBUG
Definition: Messenger.h:64
void GHepRecord::Clear ( Option_t *  opt = "")
virtual

Definition at line 919 of file GHepRecord.cxx.

References fEventFlags, fEventMask, fInteraction, and fVtx.

Referenced by CleanRecord(), and Vertex().

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 }
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
Interaction * fInteraction
attached summary information
Definition: GHepRecord.h:169
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
void GHepRecord::CompactifyDaughterLists ( void  )
virtual

compact

Definition at line 687 of file GHepRecord.cxx.

References FinalizeDaughterLists(), genie::GHepParticle::FirstDaughter(), genie::GHepParticle::FirstMother(), GetEntries(), HasCompactDaughterList(), MECModelEnuComparisons::i, genie::GHepParticle::LastDaughter(), LOG, getGoodRuns4SAM::n, Particle(), pINFO, pNOTICE, genie::GHepParticle::SetFirstDaughter(), genie::GHepParticle::SetLastDaughter(), and SwapParticles().

Referenced by UpdateDaughterLists(), and Vertex().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual void SwapParticles(int i, int j)
Definition: GHepRecord.cxx:788
int FirstDaughter(void) const
Definition: GHepParticle.h:69
const char * p
Definition: xmltok.h:285
virtual void FinalizeDaughterLists(void)
Definition: GHepRecord.cxx:835
cout<< t1-> GetEntries()
Definition: plottest35.C:29
int FirstMother(void) const
Definition: GHepParticle.h:67
int LastDaughter(void) const
Definition: GHepParticle.h:70
#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
virtual bool HasCompactDaughterList(int pos)
Definition: GHepRecord.cxx:736
#define pINFO
Definition: Messenger.h:63
#define pNOTICE
Definition: Messenger.h:62
void SetFirstDaughter(int d)
Definition: GHepParticle.h:135
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::Copy ( const GHepRecord record)
virtual

Definition at line 943 of file GHepRecord.cxx.

References EventFlags(), EventMask(), fDiffXSec, fDiffXSecPhSp, fEventFlags, fEventMask, fInteraction, fProb, fVtx, fWeight, fXSec, ResetRecord(), registry_explorer::v, and Vertex().

Referenced by genie::EventRecord::Copy(), GHepRecord(), genie::EventGenerator::ProcessEventRecord(), and Vertex().

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 }
double fProb
event probability (for given flux neutrino and density-weighted path-length for target element) ...
Definition: GHepRecord.h:180
const char * p
Definition: xmltok.h:285
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
virtual TLorentzVector * Vertex(void) const
Definition: GHepRecord.h:141
double fXSec
cross section for selected event
Definition: GHepRecord.h:181
Interaction * fInteraction
attached summary information
Definition: GHepRecord.h:169
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
virtual void ResetRecord(void)
Definition: GHepRecord.cxx:910
virtual TBits * EventFlags(void) const
Definition: GHepRecord.h:118
KinePhaseSpace_t fDiffXSecPhSp
specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)
Definition: GHepRecord.h:183
double fWeight
event weight
Definition: GHepRecord.h:179
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
virtual TBits * EventMask(void) const
Definition: GHepRecord.h:119
virtual double genie::GHepRecord::DiffXSec ( void  ) const
inlinevirtual

Definition at line 128 of file GHepRecord.h.

References fDiffXSec.

128 { return fDiffXSec; }
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
virtual KinePhaseSpace_t genie::GHepRecord::DiffXSecVars ( void  ) const
inlinevirtual

Definition at line 129 of file GHepRecord.h.

References fDiffXSecPhSp.

129 { return fDiffXSecPhSp; }
KinePhaseSpace_t fDiffXSecPhSp
specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)
Definition: GHepRecord.h:183
virtual TBits* genie::GHepRecord::EventFlags ( void  ) const
inlinevirtual

Definition at line 118 of file GHepRecord.h.

References fEventFlags.

Referenced by genie::DISHadronicSystemGenerator::AddFragmentationProducts(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), genie::SKKinematicsGenerator::CalculateKin_AtharSingleKaon(), Copy(), genie::FermiMover::KickHitNucleon(), genie::QELEventGeneratorSM::MaxDiffv(), genie::KineGeneratorWithCache::MaxXSec(), genie::QELEventGeneratorSM::MaxXSec2(), genie::DISPrimaryLeptonGenerator::ProcessEventRecord(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::DFRPrimaryLeptonGenerator::ProcessEventRecord(), genie::NuEKinematicsGenerator::ProcessEventRecord(), genie::DFRKinematicsGenerator::ProcessEventRecord(), genie::SKPrimaryLeptonGenerator::ProcessEventRecord(), genie::COHElKinematicsGenerator::ProcessEventRecord(), genie::DMDISOutgoingDarkGenerator::ProcessEventRecord(), genie::QELKinematicsGenerator::ProcessEventRecord(), genie::RESKinematicsGenerator::ProcessEventRecord(), genie::DMELKinematicsGenerator::ProcessEventRecord(), genie::IBDKinematicsGenerator::ProcessEventRecord(), genie::DISKinematicsGenerator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::DMDISKinematicsGenerator::ProcessEventRecord(), genie::PauliBlocker::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), genie::QELKinematicsGenerator::SpectralFuncExperimentalCode(), genie::DMELKinematicsGenerator::SpectralFuncExperimentalCode(), and genie::COHKinematicsGenerator::throwOnTooManyIterations().

118 { return fEventFlags; }
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
GEvGenMode_t GHepRecord::EventGenerationMode ( void  ) const

Definition at line 253 of file GHepRecord.cxx.

References genie::pdg::IsDarkMatter(), genie::pdg::IsHadron(), genie::pdg::IsIon(), genie::pdg::IsLepton(), genie::pdg::IsNucleon(), genie::kGMdDarkMatterNucleus, genie::kGMdHadronNucleus, genie::kGMdLeptonNucleus, genie::kGMdNucleonDecay, genie::kGMdPhotonNucleus, genie::kGMdUnknown, genie::kIStDecayedState, genie::kIStInitialState, genie::kPdgGamma, PandAna.Demos.pi0_spectra::p0, plot_validation_datamc::p1, Particle(), genie::GHepParticle::Pdg(), and genie::GHepParticle::Status().

Referenced by ProbePosition(), genie::Intranuke2018::ProcessEventRecord(), genie::Intranuke::ProcessEventRecord(), and TargetNucleusPosition().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
bool IsNucleon(int pdgc)
Definition: PDGUtils.cxx:309
enum genie::EGHepStatus GHepStatus_t
bool IsDarkMatter(int pdgc)
Definition: PDGUtils.cxx:125
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
bool IsHadron(int pdgc)
Definition: PDGUtils.cxx:355
const int kPdgGamma
Definition: PDGCodes.h:166
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
bool IsLepton(int pdgc)
Definition: PDGUtils.cxx:84
virtual TBits* genie::GHepRecord::EventMask ( void  ) const
inlinevirtual

Definition at line 119 of file GHepRecord.h.

References fEventMask.

Referenced by Copy().

119 { return fEventMask; }
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
void GHepRecord::FinalizeDaughterLists ( void  )
protectedvirtual

Definition at line 835 of file GHepRecord.cxx.

References genie::GHepParticle::FirstMother(), plot_validation_datamc::p1, and plot_validation_datamc::p2.

Referenced by CompactifyDaughterLists().

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 }
int FirstMother(void) const
Definition: GHepParticle.h:67
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
GHepParticle * GHepRecord::FinalStateHadronicSystem ( void  ) const
virtual

Definition at line 379 of file GHepRecord.cxx.

References FinalStateHadronicSystemPosition(), and Particle().

Referenced by genie::DISHadronicSystemGenerator::SimulateFormationZone().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int FinalStateHadronicSystemPosition(void) const
Definition: GHepRecord.cxx:507
int GHepRecord::FinalStateHadronicSystemPosition ( void  ) const
virtual

Definition at line 507 of file GHepRecord.cxx.

References genie::kIStDISPreFragmHadronicState, genie::kPdgHadronicSyst, and ParticlePosition().

Referenced by genie::DISHadronicSystemGenerator::AddFragmentationProducts(), and FinalStateHadronicSystem().

508 {
509  return this->ParticlePosition(
511 }
virtual int ParticlePosition(int pdg, GHepStatus_t i, int start=0) const
Definition: GHepRecord.cxx:181
const int kPdgHadronicSyst
Definition: PDGCodes.h:187
GHepParticle * GHepRecord::FinalStatePrimaryLepton ( void  ) const
virtual

Definition at line 370 of file GHepRecord.cxx.

References FinalStatePrimaryLeptonPosition(), and Particle().

Referenced by genie::NuETargetRemnantGenerator::AddElectronNeutrino(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_AlvarezRuso(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_BergerSehgalFM(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_ReinSehgal(), novarwgt::ConvertGenieEvent(), genie::HadronicSystemGenerator::Hadronic4pLAB(), supernova::SnovaGen::MakeMCTruth(), genie::HadronicSystemGenerator::MomentumTransferLAB(), genie::DFRPrimaryLeptonGenerator::ProcessEventRecord(), genie::DISPrimaryLeptonGenerator::ProcessEventRecord(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::SKPrimaryLeptonGenerator::ProcessEventRecord(), genie::DMDISOutgoingDarkGenerator::ProcessEventRecord(), genie::PrimaryLeptonGenerator::SetPolarization(), and genie::OutgoingDarkGenerator::SetPolarization().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int FinalStatePrimaryLeptonPosition(void) const
Definition: GHepRecord.cxx:495
int GHepRecord::FinalStatePrimaryLeptonPosition ( void  ) const
virtual

Definition at line 495 of file GHepRecord.cxx.

References genie::GHepParticle::FirstDaughter(), and Probe().

Referenced by FinalStatePrimaryLepton().

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 }
int FirstDaughter(void) const
Definition: GHepParticle.h:69
virtual GHepParticle * Probe(void) const
Definition: GHepRecord.cxx:321
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
GHepParticle * GHepRecord::FindParticle ( int  pdg,
GHepStatus_t  ist,
int  start 
) const
virtual

Definition at line 162 of file GHepRecord.cxx.

References GetEntries(), MECModelEnuComparisons::i, LOG, nentries, genie::GHepParticle::Pdg(), make_root_from_grid_output::pdg, pINFO, and genie::GHepParticle::Status().

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 }
int status
Definition: fabricate.py:1613
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
Long64_t nentries
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
int GHepRecord::FirstNonInitStateEntry ( void  )
protectedvirtual

Definition at line 775 of file GHepRecord.cxx.

References genie::kIStInitialState, genie::kIStNucleonTarget, elec2geo::pos, and genie::GHepParticle::Status().

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 }
const char * p
Definition: xmltok.h:285
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
int GHepRecord::GetPrintLevel ( )
static

Definition at line 996 of file GHepRecord.cxx.

References fPrintLevel.

Referenced by Vertex().

997 {
998  return fPrintLevel;
999 }
static int fPrintLevel
Definition: GHepRecord.h:197
vector< int > * GHepRecord::GetStableDescendants ( int  position) const
virtual

Definition at line 218 of file GHepRecord.cxx.

References genie::GHepParticle::FirstMother(), GetEntries(), MECModelEnuComparisons::i, genie::kIStStableFinalState, nentries, Particle(), and genie::GHepParticle::Status().

Referenced by genie::utils::intranuke::PreEquilibrium(), and genie::utils::intranuke2018::PreEquilibrium().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int FirstMother(void) const
Definition: GHepParticle.h:67
Long64_t nentries
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
bool GHepRecord::HasCompactDaughterList ( int  pos)
protectedvirtual

Definition at line 736 of file GHepRecord.cxx.

References febshutoff_auto::curr, genie::GHepParticle::FirstMother(), MECModelEnuComparisons::i, LOG, pDEBUG, pINFO, elec2geo::pos, and prev().

Referenced by CompactifyDaughterLists().

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 }
const char * p
Definition: xmltok.h:285
int FirstMother(void) const
Definition: GHepParticle.h:67
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
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
#define pDEBUG
Definition: Messenger.h:64
GHepParticle * GHepRecord::HitElectron ( void  ) const
virtual

Definition at line 360 of file GHepRecord.cxx.

References HitElectronPosition(), and Particle().

Referenced by genie::NuETargetRemnantGenerator::AddElectronNeutrino().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int HitElectronPosition(void) const
Definition: GHepRecord.cxx:477
int GHepRecord::HitElectronPosition ( void  ) const
virtual

Definition at line 477 of file GHepRecord.cxx.

References genie::pdg::IsElectron(), genie::kIStInitialState, Particle(), genie::GHepParticle::Pdg(), genie::GHepParticle::Status(), and TargetNucleus().

Referenced by genie::NuETargetRemnantGenerator::AddElectronNeutrino(), and HitElectron().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
const char * p
Definition: xmltok.h:285
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
bool IsElectron(int pdgc)
Definition: PDGUtils.cxx:180
GHepParticle * GHepRecord::HitNucleon ( void  ) const
virtual

Definition at line 350 of file GHepRecord.cxx.

References HitNucleonPosition(), and Particle().

Referenced by genie::AMNuGammaGenerator::AddPhoton(), genie::AMNuGammaGenerator::AddRecoilNucleon(), genie::HadronicSystemGenerator::AddTargetNucleusRemnant(), genie::AMNuGammaGenerator::AddTargetRemnant(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), genie::FermiMover::Emit2ndNucleonFromSRC(), genie::HadronicSystemGenerator::Hadronic4pLAB(), genie::FermiMover::KickHitNucleon(), supernova::SnovaGen::MakeMCTruth(), genie::NucDeExcitationSim::OxygenTargetSim(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::QELKinematicsGenerator::ProcessEventRecord(), genie::FermiMover::ProcessEventRecord(), genie::DMELKinematicsGenerator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::PauliBlocker::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), cafrwgt::CAFReweight::RetrieveGHEP(), genie::QELKinematicsGenerator::SpectralFuncExperimentalCode(), and genie::DMELKinematicsGenerator::SpectralFuncExperimentalCode().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int HitNucleonPosition(void) const
Definition: GHepRecord.cxx:454
int GHepRecord::HitNucleonPosition ( void  ) const
virtual

Definition at line 454 of file GHepRecord.cxx.

References genie::pdg::Is2NucleonCluster(), genie::pdg::IsNucleon(), genie::kIStInitialState, genie::kIStNucleonTarget, Particle(), genie::GHepParticle::Pdg(), genie::GHepParticle::Status(), and TargetNucleus().

Referenced by genie::HadronicSystemGenerator::AddFinalHadronicSyst(), genie::QELHadronicSystemGenerator::AddRecoilBaryon(), genie::IBDHadronicSystemGenerator::AddRecoilBaryon(), genie::AMNuGammaGenerator::AddRecoilNucleon(), genie::RSPPResonanceSelector::AddResonance(), genie::RESHadronicSystemGenerator::AddResonance(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), HitNucleon(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), and cafrwgt::CAFReweight::RetrieveGHEP().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
bool IsNucleon(int pdgc)
Definition: PDGUtils.cxx:309
enum genie::EGHepStatus GHepStatus_t
const char * p
Definition: xmltok.h:285
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
bool Is2NucleonCluster(int pdgc)
Definition: PDGUtils.cxx:365
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::InitRecord ( void  )
protected

Definition at line 873 of file GHepRecord.cxx.

References fDiffXSec, fDiffXSecPhSp, fEventFlags, fEventMask, fInteraction, fProb, fVtx, fWeight, fXSec, MECModelEnuComparisons::i, genie::kPSNull, LOG, genie::GHepFlags::NFlags(), pDEBUG, and pINFO.

Referenced by GHepRecord(), and ResetRecord().

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 }
double fProb
event probability (for given flux neutrino and density-weighted path-length for target element) ...
Definition: GHepRecord.h:180
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
double fXSec
cross section for selected event
Definition: GHepRecord.h:181
Interaction * fInteraction
attached summary information
Definition: GHepRecord.h:169
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
static unsigned int NFlags(void)
Definition: GHepFlags.h:77
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pINFO
Definition: Messenger.h:63
KinePhaseSpace_t fDiffXSecPhSp
specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)
Definition: GHepRecord.h:183
double fWeight
event weight
Definition: GHepRecord.h:179
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
#define pDEBUG
Definition: Messenger.h:64
virtual bool genie::GHepRecord::IsUnphysical ( void  ) const
inlinevirtual

Definition at line 120 of file GHepRecord.h.

References Accept(), and fEventFlags.

Referenced by GenerateEvent(), and Print().

120 { return (fEventFlags->CountBits()>0); }
TBits * fEventFlags
event flags indicating various pathologies or an unphysical event
Definition: GHepRecord.h:175
unsigned int GHepRecord::NEntries ( int  pdg,
GHepStatus_t  ist,
int  start = 0 
) const
virtual

Definition at line 513 of file GHepRecord.cxx.

References GetEntries(), MECModelEnuComparisons::i, nentries, genie::GHepParticle::Pdg(), and genie::GHepParticle::Status().

Referenced by novarwgt::ConvertGenieEvent().

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 }
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
Long64_t nentries
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
unsigned int GHepRecord::NEntries ( int  pdg,
int  start = 0 
) const
virtual

Definition at line 524 of file GHepRecord.cxx.

References GetEntries(), MECModelEnuComparisons::i, nentries, genie::GHepParticle::Pdg(), and make_root_from_grid_output::pdg.

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 }
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
int Pdg(void) const
Definition: GHepParticle.h:64
Long64_t nentries
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
GHepParticle * GHepRecord::Particle ( int  position) const
virtual

Definition at line 148 of file GHepRecord.cxx.

References GetEntries(), LOG, pINFO, and makeBrightnessMap::position.

Referenced by AcceptEvent(), genie::AMNuGammaGenerator::AddFinalStateNeutrino(), genie::NucDeExcitationSim::AddPhoton(), genie::QELHadronicSystemGenerator::AddRecoilBaryon(), genie::IBDHadronicSystemGenerator::AddRecoilBaryon(), genie::RSPPHadronicSystemGenerator::AddResonanceDecayProducts(), genie::FermiMover::AddTargetNucleusRemnant(), genie::QELEventGenerator::AddTargetNucleusRemnant(), genie::QELEventGeneratorSM::AddTargetNucleusRemnant(), genie::AMNuGammaGenerator::AddTargetRemnant(), CompactifyDaughterLists(), EventGenerationMode(), FinalStateHadronicSystem(), FinalStatePrimaryLepton(), GetStableDescendants(), HitElectron(), HitElectronPosition(), HitNucleon(), HitNucleonPosition(), main(), genie::utils::ghep::NeutReactionCode(), genie::utils::intranuke::PreEquilibrium(), genie::utils::intranuke2018::PreEquilibrium(), Probe(), genie::COHElHadronicSystemGenerator::ProcessEventRecord(), genie::PauliBlocker::ProcessEventRecord(), RemnantNucleus(), RemnantNucleusPosition(), genie::HAIntranuke::SimulateHadronicFinalState(), genie::HAIntranuke2018::SimulateHadronicFinalState(), genie::HNIntranuke2018::SimulateHadronicFinalState(), genie::HAIntranuke::SimulateHadronicFinalStateKinematics(), genie::HAIntranuke2018::SimulateHadronicFinalStateKinematics(), SwapParticles(), TargetNucleus(), TargetNucleusPosition(), genie::Intranuke::TransportHadrons(), genie::Intranuke2018::TransportHadrons(), and UpdateDaughterLists().

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 }
cout<< t1-> GetEntries()
Definition: plottest35.C:29
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
int GHepRecord::ParticlePosition ( int  pdg,
GHepStatus_t  i,
int  start = 0 
) const
virtual

Definition at line 181 of file GHepRecord.cxx.

References GetEntries(), MECModelEnuComparisons::i, LOG, nentries, genie::GHepParticle::Pdg(), make_root_from_grid_output::pdg, pINFO, and genie::GHepParticle::Status().

Referenced by FinalStateHadronicSystemPosition(), genie::utils::intranuke::PhaseSpaceDecay(), and genie::utils::intranuke2018::PhaseSpaceDecay().

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 }
int status
Definition: fabricate.py:1613
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
Long64_t nentries
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
int GHepRecord::ParticlePosition ( GHepParticle particle,
int  start = 0 
) const
virtual

Definition at line 200 of file GHepRecord.cxx.

References genie::GHepParticle::Compare(), GetEntries(), MECModelEnuComparisons::i, LOG, nentries, and pINFO.

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 }
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
bool Compare(const GHepParticle *p) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
Long64_t nentries
#define pINFO
Definition: Messenger.h:63
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::Print ( ostream &  stream) const

Definition at line 1001 of file GHepRecord.cxx.

References Accept(), genie::utils::print::BoolAsYNString(), genie::units::cm2, genie::GHepFlags::Describe(), genie::GHepParticle::E(), allTimeWatchdog::endl, fDiffXSec, fDiffXSecPhSp, fEventFlags, fEventMask, fInteraction, genie::GHepParticle::FirstDaughter(), genie::GHepParticle::FirstMother(), fPrintLevel, fWeight, fXSec, genie::GHepParticle::GetPolarization(), compare_h5_caf::idx, genie::GHepParticle::IsOnMassShell(), IsUnphysical(), genie::kIStFinalStateNuclearRemnant, genie::kIStInitialState, genie::kIStStableFinalState, genie::kPSQ2fE, genie::kPSQ2vfE, genie::kPSWQ2fE, genie::kPSxyfE, genie::kPSxytfE, genie::kPSyfE, genie::GHepParticle::LastDaughter(), genie::GHepParticle::LastMother(), genie::GHepParticle::Mass(), genie::GHepParticle::Name(), genie::GHepParticle::P4(), genie::GHepParticle::Pdg(), genie::GHepParticle::PolzIsSet(), Probe(), genie::GHepParticle::Px(), genie::GHepParticle::Py(), genie::GHepParticle::Pz(), genie::GHepParticle::RescatterCode(), genie::GHepParticle::Status(), Vertex(), genie::GHepParticle::Vt(), genie::GHepParticle::Vx(), genie::GHepParticle::Vy(), and genie::GHepParticle::Vz().

Referenced by genie::operator<<(), genie::EventRecord::Print(), and Vertex().

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) {
1246  else stream << "NULL Interaction!" << endl;
1247  }
1248  stream << "\n";
1249 }
int RescatterCode(void) const
Definition: GHepParticle.h:66
bool IsOnMassShell(void) const
double E(void) const
Get energy.
Definition: GHepParticle.h:92
const TLorentzVector * P4(void) const
Definition: GHepParticle.h:79
int FirstDaughter(void) const
Definition: GHepParticle.h:69
string BoolAsYNString(bool b)
Definition: PrintUtils.cxx:115
const char * p
Definition: xmltok.h:285
enum genie::EGHepFlag GHepFlag_t
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
double Mass(void) const
Mass that corresponds to the PDG code.
static const double cm2
Definition: Units.h:77
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 TLorentzVector * Vertex(void) const
Definition: GHepRecord.h:141
double Px(void) const
Get Px.
Definition: GHepParticle.h:89
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.
double fXSec
cross section for selected event
Definition: GHepRecord.h:181
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
virtual bool Accept(void) const
Definition: GHepRecord.cxx:983
bool PolzIsSet(void) const
void GetPolarization(TVector3 &polz)
virtual bool IsUnphysical(void) const
Definition: GHepRecord.h:120
double Vz(void) const
Get production z.
Definition: GHepParticle.h:97
KinePhaseSpace_t fDiffXSecPhSp
specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)
Definition: GHepRecord.h:183
double Vy(void) const
Get production y.
Definition: GHepParticle.h:96
double fWeight
event weight
Definition: GHepRecord.h:179
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
static int fPrintLevel
Definition: GHepRecord.h:197
double Vx(void) const
Get production x.
Definition: GHepParticle.h:95
double Py(void) const
Get Py.
Definition: GHepParticle.h:90
virtual double genie::GHepRecord::Probability ( void  ) const
inlinevirtual

Definition at line 126 of file GHepRecord.h.

References fProb.

126 { return fProb; }
double fProb
event probability (for given flux neutrino and density-weighted path-length for target element) ...
Definition: GHepRecord.h:180
GHepParticle * GHepRecord::Probe ( void  ) const
virtual

Definition at line 321 of file GHepRecord.cxx.

References Particle(), and ProbePosition().

Referenced by genie::NuETargetRemnantGenerator::AddElectronNeutrino(), genie::MECGenerator::AddFinalStateLepton(), genie::AMNuGammaGenerator::AddFinalStateNeutrino(), genie::DISHadronicSystemGenerator::AddFragmentationProducts(), genie::AMNuGammaGenerator::AddPhoton(), genie::QELHadronicSystemGenerator::AddRecoilBaryon(), genie::IBDHadronicSystemGenerator::AddRecoilBaryon(), genie::AMNuGammaGenerator::AddRecoilNucleon(), genie::RESHadronicSystemGenerator::AddResonance(), genie::PrimaryLeptonGenerator::AddToEventRecord(), genie::OutgoingDarkGenerator::AddToEventRecord(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_AlvarezRuso(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_BergerSehgalFM(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_ReinSehgal(), novarwgt::ConvertGenieEvent(), FinalStatePrimaryLeptonPosition(), FindhAFate(), genie::Intranuke::GenerateVertex(), genie::Intranuke2018::GenerateVertex(), genie::HadronicSystemGenerator::Hadronic4pLAB(), genie::HAIntranuke::InelasticHA(), genie::HAIntranuke2018::InelasticHA(), supernova::SnovaGen::MakeMCTruth(), genie::HadronicSystemGenerator::MomentumTransferLAB(), Print(), genie::NuEPrimaryLeptonGenerator::ProcessEventRecord(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::PrimaryLeptonGenerator::ProcessEventRecord(), genie::OutgoingDarkGenerator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), cafrwgt::CAFReweight::RetrieveGHEP(), and genie::MECGenerator::SelectNSVLeptonKinematics().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int ProbePosition(void) const
Definition: GHepRecord.cxx:389
int GHepRecord::ProbePosition ( void  ) const
virtual

Definition at line 389 of file GHepRecord.cxx.

References EventGenerationMode(), genie::kGMdDarkMatterNucleus, genie::kGMdHadronNucleus, genie::kGMdLeptonNucleus, genie::kGMdPhotonNucleus, and submit_nova_art::mode.

Referenced by genie::PrimaryLeptonGenerator::AddToEventRecord(), genie::OutgoingDarkGenerator::AddToEventRecord(), Probe(), genie::QELEventGenerator::ProcessEventRecord(), and genie::QELEventGeneratorSM::ProcessEventRecord().

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 }
Enumeration of GENIE event generation modes.
GEvGenMode_t EventGenerationMode(void) const
Definition: GHepRecord.cxx:253
GHepParticle * GHepRecord::RemnantNucleus ( void  ) const
virtual

Definition at line 340 of file GHepRecord.cxx.

References Particle(), and RemnantNucleusPosition().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int RemnantNucleusPosition(void) const
Definition: GHepRecord.cxx:433
int GHepRecord::RemnantNucleusPosition ( void  ) const
virtual

Definition at line 433 of file GHepRecord.cxx.

References genie::GHepParticle::FirstDaughter(), MECModelEnuComparisons::i, genie::pdg::IsIon(), genie::kIStStableFinalState, genie::GHepParticle::LastDaughter(), Particle(), genie::GHepParticle::Pdg(), genie::GHepParticle::Status(), and TargetNucleus().

Referenced by RemnantNucleus(), genie::Intranuke::TransportHadrons(), and genie::Intranuke2018::TransportHadrons().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
int FirstDaughter(void) const
Definition: GHepParticle.h:69
const char * p
Definition: xmltok.h:285
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
int LastDaughter(void) const
Definition: GHepParticle.h:70
virtual GHepParticle * TargetNucleus(void) const
Definition: GHepRecord.cxx:330
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::RemoveIntermediateParticles ( void  )
virtual

Definition at line 653 of file GHepRecord.cxx.

References MECModelEnuComparisons::i, sanity_check_grl::keep, genie::kIStInitialState, genie::kIStNucleonTarget, genie::kIStStableFinalState, LOG, genie::GHepParticle::Name(), pDEBUG, pNOTICE, genie::GHepParticle::SetFirstDaughter(), genie::GHepParticle::SetFirstMother(), genie::GHepParticle::SetLastDaughter(), genie::GHepParticle::SetLastMother(), and genie::GHepParticle::Status().

Referenced by Vertex().

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 }
void SetFirstMother(int m)
Definition: GHepParticle.h:133
enum genie::EGHepStatus GHepStatus_t
const char * p
Definition: xmltok.h:285
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
string Name(void) const
Name that corresponds to the PDG code.
#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
void SetLastMother(int m)
Definition: GHepParticle.h:134
#define pNOTICE
Definition: Messenger.h:62
void SetFirstDaughter(int d)
Definition: GHepParticle.h:135
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
#define pDEBUG
Definition: Messenger.h:64
void GHepRecord::ResetRecord ( void  )
virtual

Definition at line 910 of file GHepRecord.cxx.

References CleanRecord(), InitRecord(), LOG, and pDEBUG.

Referenced by Copy(), genie::EventGenerator::ProcessEventRecord(), and Vertex().

911 {
912 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
913  LOG("GHEP", pDEBUG) << "Reseting GHepRecord";
914 #endif
915  this->CleanRecord();
916  this->InitRecord();
917 }
void InitRecord(void)
Definition: GHepRecord.cxx:873
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pDEBUG
Definition: Messenger.h:64
void CleanRecord(void)
Definition: GHepRecord.cxx:902
virtual void genie::GHepRecord::SetDiffXSec ( double  xsec,
KinePhaseSpace_t  ps 
)
inlinevirtual

Definition at line 134 of file GHepRecord.h.

References fDiffXSec, fDiffXSecPhSp, and nd_projection_maker::ps.

Referenced by genie::COHKinematicsGenerator::CalculateKin_AlvarezRuso(), genie::SKKinematicsGenerator::CalculateKin_AtharSingleKaon(), genie::COHKinematicsGenerator::CalculateKin_BergerSehgal(), genie::COHKinematicsGenerator::CalculateKin_BergerSehgalFM(), genie::COHKinematicsGenerator::CalculateKin_ReinSehgal(), genie::QELEventGeneratorSM::MaxDiffv(), genie::KineGeneratorWithCache::MaxXSec(), genie::QELEventGeneratorSM::MaxXSec2(), genie::DFRKinematicsGenerator::ProcessEventRecord(), genie::NuEKinematicsGenerator::ProcessEventRecord(), genie::QELKinematicsGenerator::ProcessEventRecord(), genie::COHElKinematicsGenerator::ProcessEventRecord(), genie::RESKinematicsGenerator::ProcessEventRecord(), genie::IBDKinematicsGenerator::ProcessEventRecord(), genie::DMELKinematicsGenerator::ProcessEventRecord(), genie::DISKinematicsGenerator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::DMDISKinematicsGenerator::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), cafrwgt::CAFReweight::RetrieveGHEP(), genie::QELKinematicsGenerator::SpectralFuncExperimentalCode(), and genie::DMELKinematicsGenerator::SpectralFuncExperimentalCode().

135  { fDiffXSecPhSp = ps;
136  fDiffXSec = (xsec>0) ? xsec : 0.;
137  }
Double_t xsec[nknots]
Definition: testXsec.C:47
KinePhaseSpace_t fDiffXSecPhSp
specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)
Definition: GHepRecord.h:183
double fDiffXSec
differential cross section for selected event kinematics
Definition: GHepRecord.h:182
void GHepRecord::SetPrintLevel ( int  print_level)
static

Definition at line 992 of file GHepRecord.cxx.

References fPrintLevel.

Referenced by Initialize(), main(), and Vertex().

993 {
994  fPrintLevel = print_level;
995 }
static int fPrintLevel
Definition: GHepRecord.h:197
virtual void genie::GHepRecord::SetProbability ( double  prob)
inlinevirtual

Definition at line 132 of file GHepRecord.h.

References fProb.

Referenced by cafrwgt::CAFReweight::RetrieveGHEP().

132 { fProb = (prob>0) ? prob : 0.; }
double fProb
event probability (for given flux neutrino and density-weighted path-length for target element) ...
Definition: GHepRecord.h:180
void GHepRecord::SetUnphysEventMask ( const TBits &  mask)

Definition at line 974 of file GHepRecord.cxx.

References fEventMask, LOG, lem_server::mask, genie::GHepFlags::NFlags(), and pINFO.

Referenced by Vertex().

975 {
976  *fEventMask = mask;
977 
978  LOG("GHEP", pINFO)
979  << "Setting unphysical event mask (bits: " << GHepFlags::NFlags() - 1
980  << " -> 0) : " << *fEventMask;
981 }
TBits * fEventMask
an input bit-field mask allowing one to ignore bits set in fEventFlags
Definition: GHepRecord.h:176
static unsigned int NFlags(void)
Definition: GHepFlags.h:77
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pINFO
Definition: Messenger.h:63
void GHepRecord::SetVertex ( double  x,
double  y,
double  z,
double  t 
)
virtual

Definition at line 863 of file GHepRecord.cxx.

References fVtx.

Referenced by cafrwgt::CAFReweight::RetrieveGHEP(), evgen::GENIENeutronOscGen::setRandomEventVertexPosition(), and Vertex().

864 {
865  fVtx->SetXYZT(x,y,z,t);
866 }
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
z
Definition: test.py:28
void GHepRecord::SetVertex ( const TLorentzVector &  vtx)
virtual

Definition at line 868 of file GHepRecord.cxx.

References fVtx.

869 {
870  fVtx->SetXYZT(vtx.X(),vtx.Y(),vtx.Z(),vtx.T());
871 }
TLorentzVector * fVtx
vertex in the detector coordinate system
Definition: GHepRecord.h:172
virtual void genie::GHepRecord::SetWeight ( double  wght)
inlinevirtual
virtual void genie::GHepRecord::SetXSec ( double  xsec)
inlinevirtual

Definition at line 133 of file GHepRecord.h.

References fXSec.

Referenced by cafrwgt::CAFReweight::RetrieveGHEP(), and genie::PhysInteractionSelector::SelectInteraction().

133 { fXSec = (xsec>0) ? xsec : 0.; }
double fXSec
cross section for selected event
Definition: GHepRecord.h:181
Double_t xsec[nknots]
Definition: testXsec.C:47
Interaction * GHepRecord::Summary ( void  ) const
virtual

Definition at line 135 of file GHepRecord.cxx.

References fInteraction, LOG, and pWARN.

Referenced by genie::NuETargetRemnantGenerator::AddElectronNeutrino(), genie::HadronicSystemGenerator::AddFinalHadronicSyst(), genie::DISHadronicSystemGenerator::AddFragmentationProducts(), genie::InitialStateAppender::AddNeutrino(), genie::InitialStateAppender::AddNucleus(), genie::QELHadronicSystemGenerator::AddRecoilBaryon(), genie::IBDHadronicSystemGenerator::AddRecoilBaryon(), genie::RSPPResonanceSelector::AddResonance(), genie::RSPPHadronicSystemGenerator::AddResonanceDecayProducts(), genie::InitialStateAppender::AddStruckParticle(), genie::NuETargetRemnantGenerator::AddTargetNucleusRemnant(), genie::PrimaryLeptonGenerator::AddToEventRecord(), genie::OutgoingDarkGenerator::AddToEventRecord(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_AlvarezRuso(), genie::SKHadronicSystemGenerator::CalculateHadronicSystem_AtharSingleKaon(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_BergerSehgalFM(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_ReinSehgal(), genie::COHKinematicsGenerator::CalculateKin_AlvarezRuso(), genie::SKKinematicsGenerator::CalculateKin_AtharSingleKaon(), genie::COHKinematicsGenerator::CalculateKin_BergerSehgal(), genie::COHKinematicsGenerator::CalculateKin_BergerSehgalFM(), genie::COHKinematicsGenerator::CalculateKin_ReinSehgal(), genie::COHPrimaryLeptonGenerator::CalculatePrimaryLepton_AlvarezRuso(), novarwgt::ConvertGenieEvent(), genie::RESHadronicSystemGenerator::GetResonancePdgCode(), genie::HadronicSystemGenerator::HadronShowerCharge(), main(), supernova::SnovaGen::MakeMCTruth(), genie::QELEventGeneratorSM::MaxDiffv(), genie::KineGeneratorWithCache::MaxXSec(), genie::QELEventGeneratorSM::MaxXSec2(), genie::PrimaryLeptonGenerator::NucRestFrame2Lab(), genie::OutgoingDarkGenerator::NucRestFrame2Lab(), genie::NuEPrimaryLeptonGenerator::ProcessEventRecord(), genie::DFRHadronicSystemGenerator::ProcessEventRecord(), genie::DFRKinematicsGenerator::ProcessEventRecord(), genie::VertexGenerator::ProcessEventRecord(), genie::NuEKinematicsGenerator::ProcessEventRecord(), genie::PrimaryLeptonGenerator::ProcessEventRecord(), genie::COHElKinematicsGenerator::ProcessEventRecord(), genie::QELKinematicsGenerator::ProcessEventRecord(), genie::OutgoingDarkGenerator::ProcessEventRecord(), genie::FermiMover::ProcessEventRecord(), genie::RESKinematicsGenerator::ProcessEventRecord(), genie::DMELKinematicsGenerator::ProcessEventRecord(), genie::IBDKinematicsGenerator::ProcessEventRecord(), genie::RSPPResonanceSelector::ProcessEventRecord(), genie::DISKinematicsGenerator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::DMDISKinematicsGenerator::ProcessEventRecord(), genie::PauliBlocker::ProcessEventRecord(), genie::QELEventGeneratorSM::ProcessEventRecord(), genie::RSPPResonanceSelector::SelectResonance(), genie::QELKinematicsGenerator::SpectralFuncExperimentalCode(), and genie::DMELKinematicsGenerator::SpectralFuncExperimentalCode().

136 {
137  if(!fInteraction) {
138  LOG("GHEP", pWARN) << "Returning NULL interaction";
139  }
140  return fInteraction;
141 }
Interaction * fInteraction
attached summary information
Definition: GHepRecord.h:169
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
#define pWARN
Definition: Messenger.h:61
void GHepRecord::SwapParticles ( int  i,
int  j 
)
protectedvirtual

Definition at line 788 of file GHepRecord.cxx.

References ana::assert(), genie::GHepParticle::Copy(), GetEntries(), genie::GHepParticle::HasDaughters(), calib::j, LOG, getGoodRuns4SAM::n, genie::GHepParticle::Name(), Particle(), python.hepunit::pi, pINFO, genie::GHepParticle::SetFirstMother(), and tmp.

Referenced by CompactifyDaughterLists().

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 }
void SetFirstMother(int m)
Definition: GHepParticle.h:133
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
Float_t tmp
Definition: plot.C:36
cout<< t1-> GetEntries()
Definition: plottest35.C:29
string Name(void) const
Name that corresponds to the PDG code.
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool HasDaughters(void) const
Definition: GHepParticle.h:71
void Copy(const GHepParticle &particle)
const double j
Definition: BetheBloch.cxx:29
#define pINFO
Definition: Messenger.h:63
assert(nhit_max >=nhit_nbins)
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
GHepParticle * GHepRecord::TargetNucleus ( void  ) const
virtual

Definition at line 330 of file GHepRecord.cxx.

References Particle(), and TargetNucleusPosition().

Referenced by genie::HNIntranuke2018::AbsorbHN(), genie::AMNuGammaGenerator::AddRecoilNucleon(), genie::HadronicSystemGenerator::AddTargetNucleusRemnant(), genie::FermiMover::AddTargetNucleusRemnant(), genie::QELEventGenerator::AddTargetNucleusRemnant(), genie::QELEventGeneratorSM::AddTargetNucleusRemnant(), genie::AMNuGammaGenerator::AddTargetRemnant(), genie::PrimaryLeptonGenerator::AddToEventRecord(), genie::OutgoingDarkGenerator::AddToEventRecord(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_AlvarezRuso(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_BergerSehgalFM(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_ReinSehgal(), genie::HAIntranuke::ElasHA(), genie::HAIntranuke2018::ElasHA(), genie::HNIntranuke2018::ElasHN(), genie::utils::intranuke::Equilibrium(), genie::HNIntranuke2018::GammaInelasticHN(), genie::Intranuke::GenerateVertex(), genie::Intranuke2018::GenerateVertex(), HitElectronPosition(), HitNucleonPosition(), genie::HAIntranuke::Inelastic(), genie::HAIntranuke2018::Inelastic(), genie::HAIntranuke::InelasticHA(), genie::HAIntranuke2018::InelasticHA(), genie::FermiMover::KickHitNucleon(), genie::utils::intranuke::PreEquilibrium(), genie::VertexGenerator::ProcessEventRecord(), genie::HadronTransporter::ProcessEventRecord(), genie::NucDeExcitationSim::ProcessEventRecord(), genie::NucBindEnergyAggregator::ProcessEventRecord(), genie::QELEventGenerator::ProcessEventRecord(), genie::PauliBlocker::ProcessEventRecord(), genie::Intranuke2018::ProcessEventRecord(), genie::Intranuke::ProcessEventRecord(), RemnantNucleusPosition(), cafrwgt::CAFReweight::RetrieveGHEP(), genie::DISHadronicSystemGenerator::SimulateFormationZone(), genie::utils::intranuke::ThreeBodyKinematics(), genie::utils::intranuke2018::ThreeBodyKinematics(), genie::utils::intranuke::TwoBodyCollision(), and genie::utils::intranuke2018::TwoBodyCollision().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
virtual int TargetNucleusPosition(void) const
Definition: GHepRecord.cxx:406
int GHepRecord::TargetNucleusPosition ( void  ) const
virtual

Definition at line 406 of file GHepRecord.cxx.

References EventGenerationMode(), genie::pdg::IsIon(), genie::kGMdDarkMatterNucleus, genie::kGMdHadronNucleus, genie::kGMdLeptonNucleus, genie::kGMdNucleonDecay, genie::kGMdPhotonNucleus, genie::kIStInitialState, submit_nova_art::mode, Particle(), genie::GHepParticle::Pdg(), and genie::GHepParticle::Status().

Referenced by genie::NuETargetRemnantGenerator::AddTargetNucleusRemnant(), genie::HadronicSystemGenerator::AddTargetNucleusRemnant(), genie::FermiMover::AddTargetNucleusRemnant(), genie::QELEventGenerator::AddTargetNucleusRemnant(), genie::QELEventGeneratorSM::AddTargetNucleusRemnant(), genie::AMNuGammaGenerator::AddTargetRemnant(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_AlvarezRuso(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_BergerSehgalFM(), genie::COHHadronicSystemGenerator::CalculateHadronicSystem_ReinSehgal(), genie::FermiMover::Emit2ndNucleonFromSRC(), genie::COHElHadronicSystemGenerator::ProcessEventRecord(), cafrwgt::CAFReweight::RetrieveGHEP(), TargetNucleus(), genie::Intranuke::TransportHadrons(), and genie::Intranuke2018::TransportHadrons().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
const char * p
Definition: xmltok.h:285
GHepStatus_t Status(void) const
Definition: GHepParticle.h:65
int Pdg(void) const
Definition: GHepParticle.h:64
Enumeration of GENIE event generation modes.
GEvGenMode_t EventGenerationMode(void) const
Definition: GHepRecord.cxx:253
bool IsIon(int pdgc)
Definition: PDGUtils.cxx:40
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
void GHepRecord::UpdateDaughterLists ( void  )
protectedvirtual

Definition at line 592 of file GHepRecord.cxx.

References ana::assert(), CompactifyDaughterLists(), genie::GHepParticle::FirstDaughter(), genie::GHepParticle::FirstMother(), GetEntries(), genie::GHepParticle::LastDaughter(), LOG, Particle(), pINFO, pNOTICE, elec2geo::pos, genie::GHepParticle::SetFirstDaughter(), and genie::GHepParticle::SetLastDaughter().

Referenced by AddParticle().

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 }
virtual GHepParticle * Particle(int position) const
Definition: GHepRecord.cxx:148
int FirstDaughter(void) const
Definition: GHepParticle.h:69
const char * p
Definition: xmltok.h:285
cout<< t1-> GetEntries()
Definition: plottest35.C:29
int FirstMother(void) const
Definition: GHepParticle.h:67
int LastDaughter(void) const
Definition: GHepParticle.h:70
#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
#define pINFO
Definition: Messenger.h:63
virtual void CompactifyDaughterLists(void)
Definition: GHepRecord.cxx:687
assert(nhit_max >=nhit_nbins)
#define pNOTICE
Definition: Messenger.h:62
void SetFirstDaughter(int d)
Definition: GHepParticle.h:135
STDHEP-like event record entry that can fit a particle or a nucleus.
Definition: GHepParticle.h:40
virtual TLorentzVector* genie::GHepRecord::Vertex ( void  ) const
inlinevirtual
virtual double genie::GHepRecord::Weight ( void  ) const
inlinevirtual
virtual double genie::GHepRecord::XSec ( void  ) const
inlinevirtual

Friends And Related Function Documentation

ostream& operator<< ( ostream &  stream,
const GHepRecord event 
)
friend

Definition at line 90 of file GHepRecord.cxx.

Referenced by Vertex().

91  {
92  rec.Print(stream);
93  return stream;
94  }
caf::StandardRecord * rec
Definition: tutCAFMacro.C:20

Member Data Documentation

double genie::GHepRecord::fDiffXSec
protected

differential cross section for selected event kinematics

Definition at line 182 of file GHepRecord.h.

Referenced by Copy(), DiffXSec(), InitRecord(), Print(), and SetDiffXSec().

KinePhaseSpace_t genie::GHepRecord::fDiffXSecPhSp
protected

specifies which differential cross-section (dsig/dQ2, dsig/dQ2dW, dsig/dxdy,...)

Definition at line 183 of file GHepRecord.h.

Referenced by Copy(), DiffXSecVars(), InitRecord(), Print(), and SetDiffXSec().

TBits* genie::GHepRecord::fEventFlags
protected

event flags indicating various pathologies or an unphysical event

Definition at line 175 of file GHepRecord.h.

Referenced by Accept(), Clear(), Copy(), EventFlags(), InitRecord(), IsUnphysical(), and Print().

TBits* genie::GHepRecord::fEventMask
protected

an input bit-field mask allowing one to ignore bits set in fEventFlags

Definition at line 176 of file GHepRecord.h.

Referenced by Accept(), Clear(), Copy(), EventMask(), InitRecord(), Print(), and SetUnphysEventMask().

Interaction* genie::GHepRecord::fInteraction
protected

attached summary information

Definition at line 169 of file GHepRecord.h.

Referenced by AttachSummary(), Clear(), Copy(), InitRecord(), Print(), and Summary().

int genie::GHepRecord::fPrintLevel
staticprotected

Definition at line 197 of file GHepRecord.h.

Referenced by GetPrintLevel(), Print(), and SetPrintLevel().

double genie::GHepRecord::fProb
protected

event probability (for given flux neutrino and density-weighted path-length for target element)

Definition at line 180 of file GHepRecord.h.

Referenced by Copy(), InitRecord(), Probability(), and SetProbability().

TLorentzVector* genie::GHepRecord::fVtx
protected

vertex in the detector coordinate system

Definition at line 172 of file GHepRecord.h.

Referenced by Clear(), Copy(), InitRecord(), SetVertex(), and Vertex().

double genie::GHepRecord::fWeight
protected

event weight

Definition at line 179 of file GHepRecord.h.

Referenced by Copy(), InitRecord(), Print(), SetWeight(), and Weight().

double genie::GHepRecord::fXSec
protected

cross section for selected event

Definition at line 181 of file GHepRecord.h.

Referenced by Copy(), InitRecord(), Print(), SetXSec(), and XSec().


The documentation for this class was generated from the following files: