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/he3tube/simulation/SensitiveDetector.h>
12 #include <beast/he3tube/dataobjects/He3tubeSimHit.h>
13 #include <beast/he3tube/dataobjects/HE3G4TrackInfo.h>
14 
15 #include <framework/datastore/StoreArray.h>
16 #include <framework/datastore/RelationArray.h>
17 
18 #include <G4Track.hh>
19 #include <G4Step.hh>
20 #include "G4VProcess.hh"
21 
22 
23 namespace Belle2 {
29  namespace he3tube {
30 
32  Simulation::SensitiveDetectorBase("He3tubeSensitiveDetector", Const::invalidDetector)
33  {
34  m_trackID = 0;
35  //Make sure all collections are registered
36  StoreArray<MCParticle> mcParticles;
37  StoreArray<He3tubeSimHit> simHits;
38  StoreArray<HE3G4TrackInfo> HE3G4TrackInfos;
39  RelationArray relMCSimHit(mcParticles, simHits);
40 
41  //Register all collections we want to modify and require those we want to use
42  mcParticles.registerInDataStore();
43  simHits.registerInDataStore();
44  relMCSimHit.registerInDataStore();
45  HE3G4TrackInfos.registerInDataStore();
46  //Register the Relation so that the TrackIDs get replaced by the actual
47  //MCParticle indices after simulating the events. This is needed as
48  //secondary particles might not be stored so everything relating to those
49  //particles will be attributed to the last saved mother particle
50  registerMCParticleRelation(relMCSimHit);
51  }
52 
53  bool SensitiveDetector::step(G4Step* step, G4TouchableHistory*)
54  {
55  //Get Track information
56  const G4Track& track = *step->GetTrack();
57  const int trackID = track.GetTrackID();
58  const double depEnergy = step->GetTotalEnergyDeposit() * CLHEP::MeV;
59  const double nielEnergy = step->GetNonIonizingEnergyDeposit() * CLHEP::MeV;
60  const G4ThreeVector G4tkPos = step->GetTrack()->GetPosition();
61  float tkPos[3];
62  tkPos[0] = G4tkPos.x() * CLHEP::cm;
63  tkPos[1] = G4tkPos.y() * CLHEP::cm;
64  tkPos[2] = G4tkPos.z() * CLHEP::cm;
65  const G4ThreeVector G4tkMom = step->GetTrack()->GetMomentum();
66  float tkMom[3];
67  tkMom[0] = G4tkMom.x() * CLHEP::MeV;
68  tkMom[1] = G4tkMom.y() * CLHEP::MeV;
69  tkMom[2] = G4tkMom.z() * CLHEP::MeV;
70  const G4ThreeVector G4tkMomDir = step->GetTrack()->GetMomentumDirection();
71  float tkMomDir[3];
72  tkMomDir[0] = G4tkMomDir.x() * CLHEP::MeV;
73  tkMomDir[1] = G4tkMomDir.y() * CLHEP::MeV;
74  tkMomDir[2] = G4tkMomDir.z() * CLHEP::MeV;
75  const int tkPDG = step->GetTrack()->GetDefinition()->GetPDGEncoding();
76  const double tkKEnergy = step->GetTrack()->GetKineticEnergy();
77  const int detNb = step->GetTrack()->GetVolume()->GetCopyNo();
78  const double GlTime = step->GetPreStepPoint()->GetGlobalTime();
79  //Ignore everything below 1eV
80  //if (depEnergy < CLHEP::eV) return false;
81 
82  //Get the datastore arrays
83  StoreArray<MCParticle> mcParticles;
85  RelationArray relMCSimHit(mcParticles, simHits);
86 
87  //find out if the process that created the particle was a neutron process
88  bool neuProc = false;
89  G4String CPName;
90  if (step->GetTrack()->GetCreatorProcess() != 0) {
91  const G4VProcess* creator = step->GetTrack()->GetCreatorProcess();
92  CPName = creator->GetProcessName();
93  if (CPName.contains("Neutron")) neuProc = true;
94  }
95 
96  if (m_trackID != track.GetTrackID()) {
97  //TrackID changed, store track informations
98  m_trackID = track.GetTrackID();
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<He3tubeSimHit> He3tubeHits;
120  He3tubeSimHit* hit = He3tubeHits.appendNew(
121  depEnergy,
122  nielEnergy,
123  tkPDG,
124  tkKEnergy,
125  detNb,
126  GlTime,
127  tkPos,
128  tkMom,
129  tkMomDir,
130  neuProc
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 
140  int SensitiveDetector::saveG4TrackInfo(
141  int trackID,
142  int PDG,
143  float Mass,
144  float Energy,
145  float vtx[3],
146  float mom[3]
147  )
148  {
149  //Get the datastore arrays
150  StoreArray<HE3G4TrackInfo> HE3G4TrackInfos;
151 
152  HE3G4TrackInfos.appendNew(HE3G4TrackInfo(trackID, PDG, Mass, Energy, vtx, mom));
153 
154  int simhitNumber = HE3G4TrackInfos.getEntries() - 1;
155 
156  return (simhitNumber);
157  }//saveG4TrackInfo
158 
159 
160 
161 
162  } //he3tube namespace
164 } //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::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
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::he3tube::SensitiveDetector::SensitiveDetector
SensitiveDetector()
Constructor.
Definition: SensitiveDetector.cc:39
Belle2::fangs::SensitiveDetector::step
bool step(G4Step *step, G4TouchableHistory *) override
Step processing method.
Definition: SensitiveDetector.cc:55
Belle2::StoreArray< MCParticle >
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226
Belle2::HE3G4TrackInfo
Class HE3G4TrackInfo - Geant4 simulated hit for the Microtpc detector.
Definition: HE3G4TrackInfo.h:42