LArSoft  v07_13_02
Liquid Argon Software toolkit - http://larsoft.org/
g4b::ConvertMCTruthToG4 Class Reference

#include "ConvertMCTruthToG4.h"

Inheritance diagram for g4b::ConvertMCTruthToG4:

Public Member Functions

 ConvertMCTruthToG4 ()
 Standard constructor and destructor. More...
 
virtual ~ConvertMCTruthToG4 ()
 
void Reset ()
 Get ready to convert a new set of MCTruth objects. More...
 
void Append (art::Ptr< simb::MCTruth > &mct)
 
void Append (const simb::MCTruth *mct)
 
virtual void GeneratePrimaries (G4Event *)
 

Private Attributes

std::vector< const simb::MCTruth * > fConvertList
 List of MCTruth objects to convert for this spill. More...
 
std::map< G4int, G4int > fUnknownPDG
 map of unknown PDG codes to instances More...
 

Static Private Attributes

static G4ParticleTable * fParticleTable = 0
 Geant4's table of particle definitions. More...
 

Detailed Description

Definition at line 36 of file ConvertMCTruthToG4.h.

Constructor & Destructor Documentation

g4b::ConvertMCTruthToG4::ConvertMCTruthToG4 ( )

Standard constructor and destructor.

Definition at line 39 of file ConvertMCTruthToG4.cxx.

40  {
41  }
g4b::ConvertMCTruthToG4::~ConvertMCTruthToG4 ( )
virtual

Definition at line 44 of file ConvertMCTruthToG4.cxx.

References fUnknownPDG.

45  {
46  // Print out a list of "unknown" PDG codes we saw in the input.
47  if ( ! fUnknownPDG.empty() ){
48  std::ostringstream badtxt;
50  for ( ; i != fUnknownPDG.end(); ++i ){
51  int pdg = (*i).first;
52  badtxt << "\n Unknown PDG code = " << pdg
53  << ", seen " << (*i).second << " times.";
54 
55  const int genieLo = 2000000001;
56  const int genieHi = 2000000202;
57  if ( pdg >= genieLo && pdg <= genieHi ) {
58  badtxt << " (GENIE specific)";
59  }
60  }
61 
62  mf::LogWarning("ConvertPrimaryToGeant4")
63  << "The following unknown PDG codes were present in the "
64  << "simb::MCTruth input.\n"
65  << "They were not processed by Geant4."
66  << badtxt.str();
67  }
68 
69  }
intermediate_table::iterator iterator
std::map< G4int, G4int > fUnknownPDG
map of unknown PDG codes to instances
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning

Member Function Documentation

void g4b::ConvertMCTruthToG4::Append ( art::Ptr< simb::MCTruth > &  mct)

Add a new MCTruth object to the list of primary particles to be appended to the Geant4 event.

Definition at line 78 of file ConvertMCTruthToG4.cxx.

References art::Ptr< T >::get().

Referenced by g4b::G4Helper::G4Run().

79  {
80  this->Append( mct.get() );
81  }
T const * get() const
Definition: Ptr.h:321
void Append(art::Ptr< simb::MCTruth > &mct)
void g4b::ConvertMCTruthToG4::Append ( const simb::MCTruth mct)

Add a new MCTruth object to the list of primary particles to be appended to the Geant4 event.

Definition at line 84 of file ConvertMCTruthToG4.cxx.

References fConvertList.

85  {
86  fConvertList.push_back( mct );
87  }
std::vector< const simb::MCTruth * > fConvertList
List of MCTruth objects to convert for this spill.
void g4b::ConvertMCTruthToG4::GeneratePrimaries ( G4Event *  event)
virtual

Required of any class that inherits from G4VUserPrimaryGeneratorAction; append primary particles to a G4Event object. This method is invoked by Geant4, and is not directly called by the user application.

Definition at line 90 of file ConvertMCTruthToG4.cxx.

References fConvertList, fParticleTable, fUnknownPDG, LOG_DEBUG, simb::MCParticle::Momentum(), simb::MCParticle::PdgCode(), simb::MCParticle::Polarization(), g4b::PrimaryParticleInformation::SetMCTruth(), simb::MCParticle::StatusCode(), simb::MCParticle::T(), simb::MCParticle::Vx(), simb::MCParticle::Vy(), simb::MCParticle::Vz(), x, y, Z, and z.

91  {
92  // Take the contents of MCTruth objects and use them to
93  // initialize the G4Event.
94 
95  // A G4Event organizes its particles in terms of "vertices" and
96  // "particles", like HepMC. Unfortunately, ROOT doesn't use
97  // HepMC, so the MCTruth objects aren't organized that way.
98  // For most of the work that we'll ever do, there'll be only one
99  // vertex in the event. However, just in case there are multiple
100  // vertices (e.g., overlays, double vertex studies) I want the
101  // code to function properly.
102 
103  // So create a map of particle positions and associated
104  // G4PrimaryVertex*. Note that the map must use CLHEP's
105  // LorentzVector, and not ROOT's, since ROOT does not define an
106  // operator< for its physics vectors.
107  std::map< CLHEP::HepLorentzVector, G4PrimaryVertex* > vertexMap;
109 
110  // For each MCTruth (probably only one, but you never know):
111  // index keeps track of which MCTruth object you are using
112  size_t index = 0;
113  for( auto const* mct : fConvertList){
114 
115  // For each simb::MCParticle in the MCTruth:
116  for ( int p = 0; p != mct->NParticles(); ++p ){
117  // Implementation note: the following statement copies the
118  // MCParticle from the MCTruth, instead of getting a
119  // const reference.
120  simb::MCParticle const& particle = mct->GetParticle(p);
121 
122  // status code == 1 means "track this particle." Any
123  // other status code should be ignored by the Monte Carlo.
124  if ( particle.StatusCode() != 1 ) continue;
125 
126  // Get the Particle Data Group code for the particle.
127  G4int pdgCode = particle.PdgCode();
128 
129  // Get the vertex. Note that LArSoft/ROOT uses cm, but
130  // Geant4/CLHEP uses mm.
131  G4double x = particle.Vx() * CLHEP::cm;
132  G4double y = particle.Vy() * CLHEP::cm;
133  G4double z = particle.Vz() * CLHEP::cm;
134  G4double t = particle.T() * CLHEP::ns;
135 
136  // Create a CLHEP four-vector from the particle's vertex.
137  CLHEP::HepLorentzVector fourpos(x,y,z,t);
138 
139  // Is this vertex already in our map?
140  G4PrimaryVertex* vertex = 0;
142  if ( result == vertexMap.end() ){
143  // No, it's not, so create a new vertex and add it to the
144  // map.
145  vertex = new G4PrimaryVertex(x, y, z, t);
146  vertexMap[ fourpos ] = vertex;
147 
148  // Add the vertex to the G4Event.
149  event->AddPrimaryVertex( vertex );
150  }
151  else{
152  // Yes, it is, so use the existing vertex.
153  vertex = (*result).second;
154  }
155 
156  // Get additional particle information.
157  TLorentzVector momentum = particle.Momentum(); // (px,py,pz,E)
158  TVector3 polarization = particle.Polarization();
159 
160  // Get the particle table if necessary. (Note: we're
161  // doing this "late" because I'm not sure at what point
162  // the G4 particle table is initialized in the loading process.
163  if ( fParticleTable == 0 ){
164  fParticleTable = G4ParticleTable::GetParticleTable();
165  }
166 
167  // Get Geant4's definition of the particle.
168  G4ParticleDefinition* particleDefinition;
169 
170  if(pdgCode==0){
171  particleDefinition = fParticleTable->FindParticle("opticalphoton");
172  }
173  else
174  particleDefinition = fParticleTable->FindParticle(pdgCode);
175 
176  if ( pdgCode > 1000000000) { // If the particle is a nucleus
177  LOG_DEBUG("ConvertPrimaryToGeant4") << ": %%% Nuclear PDG code = " << pdgCode
178  << " (x,y,z,t)=(" << x
179  << "," << y
180  << "," << z
181  << "," << t << ")"
182  << " P=" << momentum.P()
183  << ", E=" << momentum.E();
184  // If the particle table doesn't have a definition yet, ask the ion
185  // table for one. This will create a new ion definition as needed.
186  if (!particleDefinition) {
187  int Z = (pdgCode % 10000000) / 10000; // atomic number
188  int A = (pdgCode % 10000) / 10; // mass number
189  particleDefinition = fParticleTable->GetIonTable()->GetIon(Z, A, 0.);
190  }
191  }
192 
193  // What if the PDG code is unknown? This has been a known
194  // issue with GENIE.
195  if ( particleDefinition == 0 ){
196  LOG_DEBUG("ConvertPrimaryToGeant4") << ": %%% Code not found = " << pdgCode;
197  fUnknownPDG[ pdgCode ] += 1;
198  continue;
199  }
200 
201  // Create a Geant4 particle to add to the vertex.
202  G4PrimaryParticle* g4particle = new G4PrimaryParticle( particleDefinition,
203  momentum.Px() * CLHEP::GeV,
204  momentum.Py() * CLHEP::GeV,
205  momentum.Pz() * CLHEP::GeV);
206 
207  // Add more particle information the Geant4 particle.
208  G4double charge = particleDefinition->GetPDGCharge();
209  g4particle->SetCharge( charge );
210  g4particle->SetPolarization( polarization.x(),
211  polarization.y(),
212  polarization.z() );
213 
214  // Add the particle to the vertex.
215  vertex->SetPrimary( g4particle );
216 
217  // Create a PrimaryParticleInformation object, and save
218  // the MCTruth pointer in it. This will allow the
219  // ParticleActionList class to access MCTruth
220  // information during Geant4's tracking.
221  PrimaryParticleInformation* primaryParticleInfo = new PrimaryParticleInformation;
222  primaryParticleInfo->SetMCTruth( mct, index, p );
223 
224  // Save the PrimaryParticleInformation in the
225  // G4PrimaryParticle for access during tracking.
226  g4particle->SetUserInformation( primaryParticleInfo );
227 
228  LOG_DEBUG("ConvertPrimaryToGeant4") << ": %%% primary PDG=" << pdgCode
229  << ", (x,y,z,t)=(" << x
230  << "," << y
231  << "," << z
232  << "," << t << ")"
233  << " P=" << momentum.P()
234  << ", E=" << momentum.E();
235 
236  } // for each particle in MCTruth
237  ++index;
238  } // for each MCTruth
239  }// GeneratePrimaries
Float_t x
Definition: compare.C:6
const TVector3 & Polarization() const
Definition: MCParticle.h:218
int PdgCode() const
Definition: MCParticle.h:216
Float_t y
Definition: compare.C:6
Double_t z
Definition: plot.C:279
static G4ParticleTable * fParticleTable
Geant4&#39;s table of particle definitions.
int StatusCode() const
Definition: MCParticle.h:215
Float_t Z
Definition: plot.C:39
intermediate_table::const_iterator const_iterator
double T(const int i=0) const
Definition: MCParticle.h:228
double Vx(const int i=0) const
Definition: MCParticle.h:225
std::map< G4int, G4int > fUnknownPDG
map of unknown PDG codes to instances
const TLorentzVector & Momentum(const int i=0) const
Definition: MCParticle.h:224
double Vz(const int i=0) const
Definition: MCParticle.h:227
#define LOG_DEBUG(id)
std::vector< const simb::MCTruth * > fConvertList
List of MCTruth objects to convert for this spill.
double Vy(const int i=0) const
Definition: MCParticle.h:226
vertex reconstruction
void g4b::ConvertMCTruthToG4::Reset ( )

Get ready to convert a new set of MCTruth objects.

Definition at line 72 of file ConvertMCTruthToG4.cxx.

References fConvertList.

Referenced by g4b::G4Helper::G4Run().

73  {
74  fConvertList.clear();
75  }
std::vector< const simb::MCTruth * > fConvertList
List of MCTruth objects to convert for this spill.

Member Data Documentation

std::vector<const simb::MCTruth*> g4b::ConvertMCTruthToG4::fConvertList
private

List of MCTruth objects to convert for this spill.

Definition at line 62 of file ConvertMCTruthToG4.h.

Referenced by Append(), GeneratePrimaries(), and Reset().

G4ParticleTable * g4b::ConvertMCTruthToG4::fParticleTable = 0
staticprivate

Geant4's table of particle definitions.

Definition at line 61 of file ConvertMCTruthToG4.h.

Referenced by GeneratePrimaries().

std::map<G4int, G4int> g4b::ConvertMCTruthToG4::fUnknownPDG
private

map of unknown PDG codes to instances

Definition at line 63 of file ConvertMCTruthToG4.h.

Referenced by GeneratePrimaries(), and ~ConvertMCTruthToG4().


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