DDG4/src/Geant4Data.cpp

0001 //==========================================================================
0002 //  AIDA Detector description implementation
0003 //--------------------------------------------------------------------------
0004 // Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
0005 // All rights reserved.
0006 //
0007 // For the licensing terms see $DD4hepINSTALL/LICENSE.
0008 // For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
0009 //
0010 // Author     : M.Frank
0011 //
0012 //==========================================================================
0013
0014 // Framework include files
0015 #include <DD4hep/Printout.h>
0016 #include <DD4hep/InstanceCount.h>
0017
0018 #include <DDG4/Geant4Data.h>
0019 #include <DDG4/Geant4StepHandler.h>
0020 #include <DDG4/Geant4FastSimHandler.h>
0021
0022 // Geant4 include files
0023 #include <G4Step.hh>
0024 #include <G4Allocator.hh>
0025 #include <G4OpticalPhoton.hh>
0026
0027 using namespace dd4hep::sim;
0028
0029 /// Default constructor
0030 SimpleRun::SimpleRun()  {
0031   InstanceCount::increment(this);
0032 }
0033
0034 /// Default destructor
0035 SimpleRun::~SimpleRun() {
0036   InstanceCount::decrement(this);
0037 }
0038
0039 /// Default constructor
0040 SimpleEvent::SimpleEvent()  {
0041   InstanceCount::increment(this);
0042 }
0043
0044 /// Default destructor
0045 SimpleEvent::~SimpleEvent() {
0046   InstanceCount::decrement(this);
0047 }
0048
0049 /// Default destructor
0050 DataExtension::~DataExtension() {
0051 }
0052
0053 /// Default constructor
0054 Geant4HitData::Geant4HitData()   {
0055   InstanceCount::increment(this);
0056 }
0057
0058 /// Default destructor
0059 Geant4HitData::~Geant4HitData() {
0060   InstanceCount::decrement(this);
0061 }
0062
0063 /// Extract the MC contribution for a given hit from the step information
0064 Geant4HitData::Contribution Geant4HitData::extractContribution(const G4Step* step) {
0065   Geant4StepHandler h(step);
0066   double deposit =
0067     (h.trackDef() == G4OpticalPhoton::OpticalPhotonDefinition()) ? h.trkEnergy() : h.totalEnergy();
0068   const G4ThreeVector& pre   = h.prePosG4();
0069   const G4ThreeVector& post  = h.postPosG4();
0070   G4ThreeVector        mom   = h.track->GetMomentum();
0071   double               len   = (post-pre).mag() ;
0072   double               position[] = { (pre.x()+post.x())/2.0,(pre.y()+post.y())/2.0,(pre.z()+post.z())/2.0 };
0073   double               momentum[] = { mom.x(), mom.y(), mom.z() };
0074   return Contribution(h.trkID(), h.trkPdgID(), deposit, h.trkTime(), len, position, momentum);
0075 }
0076
0077 /// Extract the MC contribution for a given hit from the step information with BirksLaw effect option
0078 Geant4HitData::Contribution Geant4HitData::extractContribution(const G4Step* step, bool ApplyBirksLaw) {
0079   Geant4StepHandler h(step);
0080   if ( ApplyBirksLaw == true ) h.doApplyBirksLaw();
0081   double deposit =
0082     (h.trackDef() == G4OpticalPhoton::OpticalPhotonDefinition()) ? h.trkEnergy() : h.totalEnergy();
0083   const G4ThreeVector& pre  = h.prePosG4();
0084   const G4ThreeVector& post = h.postPosG4();
0085   G4ThreeVector        mom  = h.track->GetMomentum();
0086   double length = (post-pre).mag() ;
0087   double momentum[] = { mom.x(), mom.y(), mom.z() };
0088   double position[] = { (pre.x()+post.x())/2.0,(pre.y()+post.y())/2.0,(pre.z()+post.z())/2.0 };
0089   return Contribution(h.trkID(), h.trkPdgID(), deposit, h.trkTime(), length, position, momentum);
0090 }
0091
0092 /// Extract the MC contribution for a given hit from the fast simulation spot information
0093 Geant4HitData::Contribution Geant4HitData::extractContribution(const Geant4FastSimSpot* spot) {
0094   Geant4FastSimHandler h(spot);
0095   const G4Track*       t = spot->primary;
0096   G4ThreeVector        mom = t->GetMomentum();
0097   G4ThreeVector        pos = h.avgPositionG4();
0098   double               position[] = { pos.x(), pos.y(), pos.z() };
0099   double               momentum[] = { mom.x(), mom.y(), mom.z() };
0100   return Contribution( h.trkID(), h.trkPdgID(), h.energy(), h. trkTime(), 0e0, position, momentum);
0101 }
0102
0103 /// Default constructor
0104 Geant4Tracker::Hit::Hit()  {
0105   InstanceCount::increment(this);
0106 }
0107
0108 /// Standard initializing constructor
0109 Geant4Tracker::Hit::Hit(int track_id, int pdg_id, double deposit, double time_stamp,
0110             double len, const Position& pos, const Direction& mom)
0111   : Geant4HitData(), position(pos), momentum(mom), length(len), energyDeposit(deposit),
0112     truth(track_id, pdg_id, deposit, time_stamp, len, pos, mom)
0113 {
0114   g4ID = track_id;
0115   InstanceCount::increment(this);
0116 }
0117
0118 /// Optimized constructor for sensitive detectors
0119 Geant4Tracker::Hit::Hit(const Geant4HitData::Contribution& contrib, const Direction& mom, double depo)
0120   : Geant4HitData(), position(contrib.x, contrib.y, contrib.z),
0121     momentum(mom), length(contrib.length), energyDeposit(depo), truth(contrib)
0122 {
0123   g4ID = truth.trackID;
0124   InstanceCount::increment(this);
0125 }
0126
0127 /// Default destructor
0128 Geant4Tracker::Hit::~Hit() {
0129   InstanceCount::decrement(this);
0130 }
0131
0132 /// Explicit assignment operation
0133 void Geant4Tracker::Hit::copyFrom(const Hit& c) {
0134   if ( &c != this )  {
0135     energyDeposit = c.energyDeposit;
0136     position = c.position;
0137     momentum = c.momentum;
0138     length = c.length;
0139     truth = c.truth;
0140   }
0141 }
0142
0143 /// Clear hit content
0144 Geant4Tracker::Hit& Geant4Tracker::Hit::clear() {
0145   position.SetXYZ(0, 0, 0);
0146   momentum.SetXYZ(0, 0, 0);
0147   energyDeposit = 0.0;
0148   length = 0.0;
0149   truth.clear();
0150   return *this;
0151 }
0152
0153 /// Store Geant4 point and step information into tracker hit structure.
0154 Geant4Tracker::Hit& Geant4Tracker::Hit::storePoint(const G4Step* step, const G4StepPoint* pnt) {
0155   G4Track*      trk = step->GetTrack();
0156   G4ThreeVector trm = trk->GetMomentum();
0157   G4ThreeVector pos = pnt->GetPosition();
0158   G4ThreeVector mom = pnt->GetMomentum();
0159   double        dep = step->GetTotalEnergyDeposit();
0160
0161   truth.deposit = dep;
0162   truth.trackID = trk->GetTrackID();
0163   truth.pdgID   = trk->GetDefinition()->GetPDGEncoding();
0164   truth.time    = trk->GetGlobalTime();
0165   truth.setPosition(pos.x(), pos.y(), pos.z());
0166   truth.setMomentum(trm.x(), trm.y(), trm.z());
0167
0168   energyDeposit = dep;
0169   position.SetXYZ(pos.x(), pos.y(), pos.z());
0170   momentum.SetXYZ(mom.x(), mom.y(), mom.z());
0171   length = 0;
0172   return *this;
0173 }
0174
0175 /// Store Geant4 spot information into tracker hit structure.
0176 Geant4Tracker::Hit& Geant4Tracker::Hit::storePoint(const Geant4FastSimSpot* spot)   {
0177   const G4Track* trk  = spot->primary;
0178   double         dep  = spot->hit->GetEnergy();
0179   G4ThreeVector  trm  = trk->GetMomentum();
0180   G4ThreeVector  pos  = spot->hitPosition();
0181   G4ThreeVector  mom  = trk->GetMomentum().unit() * dep;
0182
0183   this->truth.deposit = dep;
0184   this->truth.trackID = trk->GetTrackID();
0185   this->truth.time    = trk->GetGlobalTime();
0186   this->truth.pdgID   = trk->GetDefinition()->GetPDGEncoding();
0187   this->truth.setPosition(pos.x(), pos.y(), pos.z());
0188   this->truth.setMomentum(trm.x(), trm.y(), trm.z());
0189
0190   this->energyDeposit = dep;
0191   this->position.SetXYZ(pos.x(), pos.y(), pos.z());
0192   this->momentum.SetXYZ(mom.x(), mom.y(), mom.z());
0193   this->length = 0;
0194   return *this;
0195 }
0196
0197 /// Default constructor (for ROOT)
0198 Geant4Calorimeter::Hit::Hit()  {
0199   InstanceCount::increment(this);
0200 }
0201
0202 /// Standard constructor
0203 Geant4Calorimeter::Hit::Hit(const Position& pos) : position(pos)  {
0204   InstanceCount::increment(this);
0205 }
0206
0207 /// Default destructor
0208 Geant4Calorimeter::Hit::~Hit() {
0209   InstanceCount::decrement(this);
0210 }