Belle II Software  release-05-01-25
SensitiveDetector.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Ritter, Igal Jaegle *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <beast/microtpc/simulation/SensitiveDetector.h>
12 #include <beast/microtpc/dataobjects/MicrotpcSimHit.h>
13 #include <beast/microtpc/dataobjects/TPCG4TrackInfo.h>
14 
15 #include <framework/datastore/StoreArray.h>
16 #include <framework/datastore/RelationArray.h>
17 
18 #include <G4Track.hh>
19 #include <G4Step.hh>
20 
21 namespace Belle2 {
27  namespace microtpc {
28 
30  Simulation::SensitiveDetectorBase("MicrotpcSensitiveDetector", Const::invalidDetector)
31  {
32  m_trackID = 0;
33  //Make sure all collections are registered
34  StoreArray<MCParticle> mcParticles;
35  StoreArray<MicrotpcSimHit> simHits;
36  StoreArray<TPCG4TrackInfo> TPCG4TrackInfos;
37  RelationArray relMCSimHit(mcParticles, simHits);
38 
39  //Register all collections we want to modify and require those we want to use
40  mcParticles.registerInDataStore();
41  simHits.registerInDataStore();
42  relMCSimHit.registerInDataStore();
43  TPCG4TrackInfos.registerInDataStore();
44 
45  //Register the Relation so that the TrackIDs get replaced by the actual
46  //MCParticle indices after simulating the events. This is needed as
47  //secondary particles might not be stored so everything relating to those
48  //particles will be attributed to the last saved mother particle
49  registerMCParticleRelation(relMCSimHit);
50  }
51 
52  bool SensitiveDetector::step(G4Step* step, G4TouchableHistory*)
53  {
54  //Get Track information
55  const G4Track& track = *step->GetTrack();
56  const int trackID = track.GetTrackID();
57  const double depEnergy = step->GetTotalEnergyDeposit() * CLHEP::MeV;
58  const double nielEnergy = step->GetNonIonizingEnergyDeposit() * CLHEP::MeV;
59  const G4ThreeVector G4tkPos = step->GetTrack()->GetPosition();
60  float tkPos[3];
61  tkPos[0] = G4tkPos.x() * CLHEP::cm;
62  tkPos[1] = G4tkPos.y() * CLHEP::cm;
63  tkPos[2] = G4tkPos.z() * CLHEP::cm;
64  const G4ThreeVector G4tkMom = step->GetTrack()->GetMomentum();
65  float tkMom[3];
66  tkMom[0] = G4tkMom.x() * CLHEP::MeV;
67  tkMom[1] = G4tkMom.y() * CLHEP::MeV;
68  tkMom[2] = G4tkMom.z() * CLHEP::MeV;
69  const G4ThreeVector G4tkMomDir = step->GetTrack()->GetMomentumDirection();
70  float tkMomDir[3];
71  tkMomDir[0] = G4tkMomDir.x() * CLHEP::MeV;
72  tkMomDir[1] = G4tkMomDir.y() * CLHEP::MeV;
73  tkMomDir[2] = G4tkMomDir.z() * CLHEP::MeV;
74  const int tkPDG = step->GetTrack()->GetDefinition()->GetPDGEncoding();
75  const double tkKEnergy = step->GetTrack()->GetKineticEnergy();
76  const int detNb = step->GetTrack()->GetVolume()->GetCopyNo();
77  const double GlTime = step->GetPreStepPoint()->GetGlobalTime();
78  //Ignore everything below 1eV
79  //if (depEnergy < CLHEP::eV) return false;
80 
81  //Get the datastore arrays
82  StoreArray<MCParticle> mcParticles;
84  RelationArray relMCSimHit(mcParticles, simHits);
85  /*
86  //find out if the process that created the particle was a neutron process
87  bool neuProc = false;
88  G4String CPName;
89  if (step->GetTrack()->GetCreatorProcess() != 0) {
90  const G4VProcess* creator = step->GetTrack()->GetCreatorProcess();
91  CPName = creator->GetProcessName();
92  if (CPName.contains("Neutron")) neuProc = true;
93  }
94  */
95  if (m_trackID != track.GetTrackID()) {
96  //TrackID changed, store track informations
97  m_trackID = track.GetTrackID();
98  }
99 
100  //Save Hit if track leaves volume or is killed
101  if (track.GetNextVolume() != track.GetVolume() || track.GetTrackStatus() >= fStopAndKill) {
102  auto& mcparticle = Simulation::TrackInfo::getInfo(track);
103  int PDG = mcparticle.getPDG();
104  float Mass = mcparticle.getMass();
105  float Energy = mcparticle.getEnergy();
106  float vtx[3];
107  vtx[0] = mcparticle.getProductionVertex().X();
108  vtx[1] = mcparticle.getProductionVertex().Y();
109  vtx[2] = mcparticle.getProductionVertex().Z();
110  float mom[3];
111  mom[0] = mcparticle.getMomentum().X();
112  mom[1] = mcparticle.getMomentum().Y();
113  mom[2] = mcparticle.getMomentum().Z();
114  saveG4TrackInfo(m_trackID, PDG, Mass, Energy, vtx, mom);
115  //Reset TrackID
116  m_trackID = 0;
117  }
118 
119  StoreArray<MicrotpcSimHit> MicrotpcHits;
120  MicrotpcSimHit* hit = MicrotpcHits.appendNew(
121  trackID,
122  depEnergy,
123  nielEnergy,
124  tkPDG,
125  tkKEnergy,
126  detNb,
127  GlTime,
128  tkPos,
129  tkMom,
130  tkMomDir
131  );
132 
133  //Add Relation between SimHit and MCParticle with a weight of 1. Since
134  //the MCParticle index is not yet defined we use the trackID from Geant4
135  relMCSimHit.add(trackID, hit->getArrayIndex(), 1.0);
136 
137  return true;
138  }
139 
141  (
142  int trackID,
143  int PDG,
144  float Mass,
145  float Energy,
146  float vtx[3],
147  float mom[3]
148  )
149  {
150  //Get the datastore arrays
151  StoreArray<TPCG4TrackInfo> TPCG4TrackInfos;
152 
153  TPCG4TrackInfos.appendNew(TPCG4TrackInfo(trackID, PDG, Mass, Energy, vtx, mom));
154 
155  int simhitNumber = TPCG4TrackInfos.getEntries() - 1;
156 
157  return (simhitNumber);
158  }//saveG4TrackInfo
159 
160  } //microtpc namespace
162 } //Belle2 namespace
Belle2::StoreArray::appendNew
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:256
Belle2::RelationArray
Low-level class to create/modify relations between StoreArrays.
Definition: RelationArray.h:72
Belle2::Simulation::UserInfo::getInfo
static Payload getInfo(Carrier &obj)
Static function to just return UserInformation attached to the obj of type Carrier.
Definition: UserInfo.h:110
Belle2::he3tube::SensitiveDetector::saveG4TrackInfo
int saveG4TrackInfo(int trackID, int PDG, float Mass, float Energy, float vtx[3], float mom[3])
Save saveG4TrackInfo into datastore.
Definition: SensitiveDetector.cc:148
Belle2::MicrotpcSimHit
ClassMicrotpcSimHit - Geant4 simulated hit for the Microtpc detector.
Definition: MicrotpcSimHit.h:40
Belle2::he3tube::SensitiveDetector::step
bool step(G4Step *step, G4TouchableHistory *) override
Step processing method.
Definition: SensitiveDetector.cc:61
Belle2::Simulation::SensitiveDetectorBase::registerMCParticleRelation
static void registerMCParticleRelation(const std::string &name, RelationArray::EConsolidationAction ignoreAction=RelationArray::c_negativeWeight)
Register an relation involving MCParticles.
Definition: SensitiveDetectorBase.cc:24
Belle2::microtpc::SensitiveDetector::SensitiveDetector
SensitiveDetector()
Constructor.
Definition: SensitiveDetector.cc:37
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TPCG4TrackInfo
Class TPCG4TrackInfo - Geant4 simulated hit for the Microtpc detector.
Definition: TPCG4TrackInfo.h:42
Belle2::StoreArray< MCParticle >
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226