DD4hep - The AIDA detector description toolkit for high energy physics experiments
DD4hep  Rev:Unversioneddirectory
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Geant4ExtraParticles.cpp
Go to the documentation of this file.
1 //
2 // Authors: Tomohiko Tanabe <tomohiko@icepp.s.u-tokyo.ac.jp>
3 // Taikan Suehara <suehara@icepp.s.u-tokyo.ac.jp>
4 // Proted from Mokka by A.Sailer (CERN )
5 //
6 // $Id$
7 // $Name: $
8 
22 #include "Geant4ExtraParticles.h"
24 #include "DDG4/Geant4Kernel.h"
25 #include "DDG4/Factories.h"
26 
27 #include "G4ParticleTable.hh"
28 #include "G4ParticleDefinition.hh"
29 #include "G4SystemOfUnits.hh"
30 #include "G4PhysicalConstants.hh"
31 #include "G4Version.hh"
32 
33 #include "DDG4/Factories.h"
34 #include <fstream>
35 #include <sstream>
36 #include <string>
37 
38 
39 #include "G4RunManager.hh"
40 
41 using namespace DD4hep::Simulation;
42 
43 Geant4ExtraParticles::Geant4ExtraParticles(Geant4Context* ctxt, const std::string& nam)
44  : Geant4PhysicsConstructor(ctxt, nam), m_decay(0), m_ionise(0), m_scatter(0)
45 {
46  declareProperty("pdgfile", m_pdgfile);
47 }
48 
53 }
54 
55 // bool Geant4ExtraParticles::FileExists() {
56 // std::ifstream pdgFile( m_pdgfile.c_str(), std::ifstream::in );
57 // return pdgFile.is_open();
58 // }
59 
60 void Geant4ExtraParticles::constructParticle(Constructor& ) {
61  if (m_pdgfile.empty()) return;
62 
63  info("pdgfile: %s",m_pdgfile.c_str());
64 
65  G4ParticleTable *theParticleTable = G4ParticleTable::GetParticleTable();
66  std::ifstream pdgFile( m_pdgfile.c_str(), std::ifstream::in );
67 
68  if (!pdgFile.is_open()) {
69  except("Could not open pdgfile: %s",m_pdgfile.c_str());
70  return;
71  }
72 
73  info("opened pdgfile: %s",m_pdgfile.c_str());
74 
75  while ( !pdgFile.eof() ) {
76  // read line
77  std::string linebuf;
78  getline( pdgFile, linebuf );
79 
80  // ignore comments
81  if (linebuf.substr(0,1) == "#") continue;
82  if (linebuf.substr(0,2) == "//") continue;
83 
84  // ignore empty lines
85  if (linebuf.empty()) continue;
86 
87  // parse line
88  int pdg;
89  std::string nam;
90  double charge;
91  double mass;
92  double width;
93  double lifetime;
94 
95  std::istringstream istr(linebuf);
96 
97  istr >> pdg >> nam >> charge >> mass >> width >> lifetime;
98 
99  // do add particles that don't fly
100  // if (lifetime == 0) continue;
101 
102  if(width<0) width = 0;
103 
104  // normalize to G4 units
105  mass *= GeV;
106 
107  if (charge != 0) {
108  charge /= 3.;
109  }
110 
111  if (lifetime > 0) {
112  lifetime = lifetime*mm/c_light;
113  }
114 
115  if (width == 0 && lifetime > 0) {
116  width = hbar_Planck/lifetime;
117  }
118 
119  // don't add if the particle already exists
120  G4ParticleDefinition* p = theParticleTable->FindParticle(pdg);
121  if ( !p ) {
122 
123  if (abs(pdg)>80 && abs(pdg)<=100) {
124  // don't add generator-specific particles
125  } else {
126  /*
127  if (pdg==5122) {
128  G4cout << "Lambda_b0: " << "PDG=" << pdg << ", name=" << name << ", chrg=" << charge
129  << ", mass=" << mass << ", width=" << width << ", lifetime=" << lifetime << "\n";
130  G4cout << "debug: mass=" << 5.62 << ", width =" << 1.39e-12/6.582e-16 << ", lifetime=" << 1.39e-12 << "\n";
131  }
132  //*/
133  p = new G4ParticleDefinition(
134  nam, // name
135  mass, // mass
136  width, // width
137  charge, // charge
138  0, // 2*spin
139  0, // parity
140  0, // C-conjugation
141  0, // 2*isospin
142  0, // 2*isospin3
143  0, // G-parity
144  "extra", // type
145  0, // lepton number
146  0, // baryon number
147  pdg, // PDG encoding
148  width==0?true:false, // stable
149  lifetime, // lifetime
150  NULL, // decay table
151  false); // short lived
152  }
153  }
154  }
155 
156  G4cout << "Loaded extra particles using file: " << m_pdgfile << G4endl;
157 }
158 
159 void Geant4ExtraParticles::constructProcess(Constructor& ctor) {
160  G4ParticleTable::G4PTblDicIterator* ParticleIterator = ctor.particleIterator();
161 #if G4VERSION_NUMBER < 940
162  if ( 0 == _scatter ) _scatter=new G4hMultipleScattering();
163  if ( 0 == _ionise ) _ionise=new G4hIonisation()
164  if ( 0 == _decay ) _decay=new G4Decay()
165 #endif
166  while((*ParticleIterator)()) {
167  G4ParticleDefinition* pdef = ParticleIterator->value();
168  G4ProcessManager* pmgr = pdef->GetProcessManager();
169  if (pdef->GetParticleType() == "extra") {
170  if (pdef->GetPDGCharge() != 0) {
171 #if G4VERSION_NUMBER < 940
172  pmgr->AddProcess(_scatter, -1, 1, 1); // multiple scattering
173  pmgr->AddProcess(_ionise, -1, 2, 2); // ionisation
174  pmgr->AddProcess(_decay, -1, -1, 2); // decay
175 #else
176  pmgr->AddProcess(new G4hMultipleScattering(), -1, 1, 1); //multiple scattering
177  pmgr->AddProcess(new G4hIonisation(), -1, 2, 2); // ionisation
178  pmgr->AddProcess(new G4Decay(), -1, -1, 2); // decay
179 #endif
180 
181  } else {
182 
183 #if G4VERSION_NUMBER < 940
184  pmgr->AddProcess(_scatter=new G4hMultipleScattering(), -1, 1, 1); // multiple scattering
185  pmgr->AddProcess(_decay=new G4Decay(), -1, -1, 2); // decay
186 #else
187  // pmgr->AddProcess(new G4hMultipleScattering(), -1, 1, 1); // multiple scattering
188  pmgr->AddProcess(new G4Decay(), -1, -1, 2); // decay
189 #endif
190 
191  }
192  }
193  }
194 }
195 
Implementation base of a Geant4 physics constructor.
Geant4ExtraParticles(Geant4Context *ctxt, const std::string &nam)
Standard constructor with initailization parameters.
virtual void constructProcess(Constructor &ctor)
Callback to construct processes (uses the G4 particle table)
void deletePtr(T *&p)
Helper to delete objects from heap and reset the pointer. Saves many many lines of code...
Definition: Primitives.h:234
#define DECLARE_GEANT4ACTION(name)
Plugin defintion to create Geant4Action objects.
Definition: Factories.h:166
void except(const char *fmt,...) const
Support of exceptions: Print fatal message and throw runtime_error.
static const double hbar_Planck
Definition: DD4hepUnits.h:337
static const double GeV
Definition: DD4hepUnits.h:149
static const double c_light
Definition: DD4hepUnits.h:328
Generic context to extend user, run and event information.
void info(const char *fmt,...) const
Support of info messages.
virtual void constructParticle(Constructor &ctor)
Callback to construct particles.
static const double mm
Definition: DD4hepUnits.h:69
Geant4Action & declareProperty(const std::string &nam, T &val)
Declare property.
Definition: Geant4Action.h:358
virtual ~Geant4ExtraParticles()
Default destructor.