Belle II Software development
SensitiveDetector.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <beast/pindiode/simulation/SensitiveDetector.h>
10#include <beast/pindiode/dataobjects/PindiodeSimHit.h>
11
12#include <framework/datastore/StoreArray.h>
13#include <framework/datastore/RelationArray.h>
14
15#include <G4Track.hh>
16#include <G4Step.hh>
17
18namespace Belle2 {
24 namespace pindiode {
25
27 Simulation::SensitiveDetectorBase("PindiodeSensitiveDetector", Const::invalidDetector)
28 {
29 m_hitNum = 0;
30 m_EvnetNumber = 0;
32 m_trackID = 0;
33 m_startTime = 0;
34 m_endTime = 0;
35 m_WightedTime = 0;
36 m_startEnergy = 0;
38 m_trackLength = 0;
39 iECLCell = 0;
40 TimeIndex = 0;
41 local_pos = 0;
42 T_ave = 0;
43 firstcall = 0;
44 m_phiID = 0;
45 m_thetaID = 0;
46 m_cellID = 0;
47
48 //Make sure all collections are registered
49 StoreArray<MCParticle> mcParticles;
51 RelationArray relMCSimHit(mcParticles, simHits);
52
53 //Register all collections we want to modify and require those we want to use
54 mcParticles.registerInDataStore();
55 simHits.registerInDataStore();
56 relMCSimHit.registerInDataStore();
57
58 //Register the Relation so that the TrackIDs get replaced by the actual
59 //MCParticle indices after simulating the events. This is needed as
60 //secondary particles might not be stored so everything relating to those
61 //particles will be attributed to the last saved mother particle
62 registerMCParticleRelation(relMCSimHit);
63 }
64
66 {
67
68 }
69
70 bool SensitiveDetector::step(G4Step* step, G4TouchableHistory*)
71 {
72 const G4StepPoint& preStep = *step->GetPreStepPoint();
73 const G4StepPoint& postStep = *step->GetPostStepPoint();
74
75 G4Track& track = *step->GetTrack();
76 if (m_trackID != track.GetTrackID()) {
77 //TrackID changed, store track informations
78 m_trackID = track.GetTrackID();
79 //Get momentum
80 m_momentum = preStep.GetMomentum() ;
81 //Get energy
82 m_startEnergy = preStep.GetKineticEnergy() ;
83 //Reset energy deposit;
85 //Reset Wighted Time;
86 m_WightedTime = 0;
87 //Reset m_WightedPos;
88 m_WightedPos.SetXYZ(0, 0, 0);
89
90 }
91 //Update energy deposit
92 m_energyDeposit += step->GetTotalEnergyDeposit() ;
93
94 m_startTime = preStep.GetGlobalTime();
95 m_endTime = postStep.GetGlobalTime();
96 m_WightedTime += (m_startTime + m_endTime) / 2 * (step->GetTotalEnergyDeposit());
97
98 m_startPos = preStep.GetPosition();
99 m_endPos = postStep.GetPosition();
100 ROOT::Math::XYZVector position((m_startPos.getX() + m_endPos.getX()) / 2 / CLHEP::cm,
101 (m_startPos.getY() + m_endPos.getY()) / 2 / CLHEP::cm,
102 (m_startPos.getZ() + m_endPos.getZ()) / 2 / CLHEP::cm);
103 m_WightedPos += position * (step->GetTotalEnergyDeposit());
104
105 //Save Hit if track leaves volume or is killed
106 if (track.GetNextVolume() != track.GetVolume() || track.GetTrackStatus() >= fStopAndKill) {
107 int pdgCode = track.GetDefinition()->GetPDGEncoding();
108
109 //const G4VPhysicalVolume& v = * track.GetVolume();
110 //G4ThreeVector posCell = v.GetTranslation();
111 // Get layer ID
112
113 //if (v.GetName().find("Crystal") != std::string::npos) {
114 //CsiGeometryPar* eclp = CsiGeometryPar::Instance();
115 m_cellID = step->GetTrack()->GetVolume()->GetCopyNo();
116
117 double dTotalEnergy = 1 / m_energyDeposit; //avoid the error no match for 'operator/'
118 if (m_energyDeposit > 0.) {
120 m_energyDeposit, m_momentum, m_WightedPos * dTotalEnergy);
121 //}
122 }
123
124 //Reset TrackID
125 m_trackID = 0;
126 }
127
128 return true;
129 }
130 /*
131 //Get Track information
132 const G4Track& track = *step->GetTrack();
133 const int trackID = track.GetTrackID();
134 const double depEnergy = step->GetTotalEnergyDeposit() * CLHEP::MeV;
135 const double nielEnergy = step->GetNonIonizingEnergyDeposit() * CLHEP::MeV;
136 const G4ThreeVector G4tkPos = step->GetTrack()->GetPosition();
137 float tkPos[3];
138 tkPos[0] = G4tkPos.x() * CLHEP::cm;
139 tkPos[1] = G4tkPos.y() * CLHEP::cm;
140 tkPos[2] = G4tkPos.z() * CLHEP::cm;
141 const G4ThreeVector G4tkMom = step->GetTrack()->GetMomentum();
142 float tkMom[3];
143 tkMom[0] = G4tkMom.x() * CLHEP::MeV;
144 tkMom[1] = G4tkMom.y() * CLHEP::MeV;
145 tkMom[2] = G4tkMom.z() * CLHEP::MeV;
146 const G4ThreeVector G4tkMomDir = step->GetTrack()->GetMomentumDirection();
147 float tkMomDir[3];
148 tkMomDir[0] = G4tkMomDir.x() * CLHEP::MeV;
149 tkMomDir[1] = G4tkMomDir.y() * CLHEP::MeV;
150 tkMomDir[2] = G4tkMomDir.z() * CLHEP::MeV;
151 const int tkPDG = step->GetTrack()->GetDefinition()->GetPDGEncoding();
152 const double tkKEnergy = step->GetTrack()->GetKineticEnergy();
153 const int detNb = step->GetTrack()->GetVolume()->GetCopyNo();
154 const double GlTime = step->GetPreStepPoint()->GetGlobalTime();
155 //Ignore everything below 1eV
156 if (depEnergy < CLHEP::eV) return false;
157
158 //Get the datastore arrays
159 StoreArray<MCParticle> mcParticles;
160 StoreArray<PindiodeSimHit> simHits;
161 RelationArray relMCSimHit(mcParticles, simHits);
162
163 StoreArray<PindiodeSimHit> PindiodeHits;
164 if (!PindiodeHits.isValid()) PindiodeHits.create();
165 PindiodeSimHit* hit = PindiodeHits.appendNew(
166 depEnergy,
167 nielEnergy,
168 tkPDG,
169 tkKEnergy,
170 detNb,
171 GlTime,
172 tkPos,
173 tkMom,
174 tkMomDir
175 );
176
177 //Add Relation between SimHit and MCParticle with a weight of 1. Since
178 //the MCParticle index is not yet defined we use the trackID from Geant4
179 relMCSimHit.add(trackID, hit->getArrayIndex(), 1.0);
180
181 return true;
182 }
183 */
184
186 const G4int cellId,
187 const G4int trackID,
188 const G4int pid,
189 const G4double tof,
190 const G4double edep,
191 G4ThreeVector mom,
192 ROOT::Math::XYZVector posAve)
193 {
194
195 //Get the datastore arraus
196 StoreArray<MCParticle> mcParticles;
198 RelationArray relMPindiodemHit(mcParticles, simHits);
199
200 StoreArray<PindiodeSimHit> PindiodeHits;
201 RelationArray pindiodeSimHitRel(mcParticles, PindiodeHits);
202 ROOT::Math::XYZVector momentum(mom.getX() / CLHEP::GeV, mom.getY() / CLHEP::GeV, mom.getZ() / CLHEP::GeV);
203 PindiodeHits.appendNew(cellId, trackID, pid, tof / CLHEP::ns, edep / CLHEP::GeV, momentum, posAve);
204 int simhitNumber = PindiodeHits.getEntries() - 1;
205 B2DEBUG(150, "HitNumber: " << simhitNumber);
206 pindiodeSimHitRel.add(trackID, simhitNumber);
207 return (simhitNumber);
208 }//saveSimHit
209
210
211 } //pindiode namespace
213} //Belle2 namespace
This class provides a set of constants for the framework.
Definition: Const.h:34
Low-level class to create/modify relations between StoreArrays.
Definition: RelationArray.h:62
void add(index_type from, index_type to, weight_type weight=1.0)
Add a new element to the relation.
static void registerMCParticleRelation(const std::string &name, RelationArray::EConsolidationAction ignoreAction=RelationArray::c_negativeWeight)
Register an relation involving MCParticles.
bool registerInDataStore(DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut)
Register the object/array in the DataStore.
Accessor to arrays stored in the data store.
Definition: StoreArray.h:113
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:246
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:216
int m_thetaID
The current phi ID in an event.
double m_energyDeposit
particle energy at the entrance in volume
int m_EvnetNumber
The current number of created hits in an event.
int saveSimHit(const G4int cellId, const G4int trackID, const G4int pid, const G4double tof, const G4double edep, G4ThreeVector mom, ROOT::Math::XYZVector WightedPos)
Save PindiodeSimHit into datastore.
int TimeIndex
Hit Energy of StoreArray.
int m_trackID
The current number of created hits in an event.
int m_hitNum
members of SensitiveDetector
int firstcall
flight time to diode sensor
double T_ave
position alongthe vector of crystal axis
double local_pos
vector of crystal axis
G4ThreeVector m_momentum
Wighted step Position.
int m_cellID
The current theta ID in an event.
double m_startTime
momentum of track
int iECLCell
length of the track in the volume
double m_trackLength
energy deposited in volume
int m_oldEvnetNumber
The current number of created hits in an event.
ROOT::Math::XYZVector m_WightedPos
Position of poststep.
G4ThreeVector m_endPos
Position of prestep.
bool step(G4Step *step, G4TouchableHistory *) override
Step processing method.
Abstract base class for different kinds of events.