Belle II Software development
SensitiveDetector.h
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#pragma once
10#ifndef VXD_SIMULATION_SENSITIVEDETECTOR_H
11#define VXD_SIMULATION_SENSITIVEDETECTOR_H
12
13#include <vxd/simulation/SensitiveDetectorBase.h>
14#include <framework/datastore/StoreArray.h>
15#include <framework/datastore/RelationArray.h>
16#include <array>
17
18#ifdef VXD_SENSITIVEDETECTOR_DEBUG
19#include <vxd/simulation/SensitiveDetectorDebugHelper.h>
20#endif
21
22namespace Belle2 {
28 namespace VXD {
29
64 template <class SimHitClass, class TrueHitClass> class SensitiveDetector: public SensitiveDetectorBase {
65 public:
68
69 private:
75 int saveTrueHit(const SensorTraversal& traversal) override;
76
84 int saveSimHit(const SensorTraversal& traversal,
85 const SensorTraversal::range& points) override;
86
94 void saveRelations(const SensorTraversal& traversal, int trueHitIndex,
95 std::vector<std::pair<unsigned int, float>> simHitIndices) override;
96
102 std::array<float, 3> vecToFloat(const G4ThreeVector& vec)
103 {
104 return std::array<float, 3> {{(float)vec.x(), (float)vec.y(), (float)vec.z()}};
105 }
106
119 };
120
121 template <class SimHitClass, class TrueHitClass>
123 VXD::SensitiveDetectorBase(sensorInfo)
124 {
125 //Register output collections.
126 //m_mcparticles.isRequired();
127 m_simhits.registerInDataStore();
128 m_truehits.registerInDataStore();
129 m_relMCSimHits.registerInDataStore();
130 m_relMCTrueHits.registerInDataStore();
131 m_relTrueSimHits.registerInDataStore();
134 }
135
136 template <class SimHitClass, class TrueHitClass>
138 {
139 //We have the full sensor traversal information, we only lack the midpoint so lets get it
140 StepInformation midPoint = findMidPoint(traversal);
141
142 //By now everything is done: just convert the position and momenta in float[3] and create a truehit
143 auto posEntry = vecToFloat(traversal.front().position);
144 auto posMidPoint = vecToFloat(midPoint.position);
145 auto posExit = vecToFloat(traversal.back().position);
146 auto momEntry = vecToFloat(traversal.front().momentum);
147 auto momMidPoint = vecToFloat(midPoint.momentum);
148 auto momExit = vecToFloat(traversal.back().momentum);
149 //And we should convert the electrons back in energy ...
150 const double energyDep = traversal.getElectrons() * Const::ehEnergy;
151
152 //create the truehit
153 int trueHitIndex = m_truehits.getEntries();
154 m_truehits.appendNew(getSensorID(), posEntry.data(), posMidPoint.data(), posExit.data(),
155 momEntry.data(), momMidPoint.data(), momExit.data(), energyDep, midPoint.time);
156#ifdef VXD_SENSITIVEDETECTOR_DEBUG
158#endif
159 return trueHitIndex;
160 }
161
162 template <class SimHitClass, class TrueHitClass>
164 const SensorTraversal::range& points)
165 {
166 //Convert position to floats here
167 auto posIn = vecToFloat(points.first->position);
168 auto posOut = vecToFloat(points.second->position);
169 auto electronProfile = simplifyEnergyDeposit(points);
170
171 //Create the simhit
172 int simHitIndex = m_simhits.getEntries();
173 SimHitClass* simhit = m_simhits.appendNew(getSensorID(), traversal.getPDGCode(), points.first->time,
174 posIn.data(), posOut.data());
175 simhit->setEnergyDeposit(electronProfile);
176#ifdef VXD_SENSITIVEDETECTOR_DEBUG
177 int startPoint = std::distance(traversal.begin(), (SensorTraversal::const_iterator)points.first);
178 int endPoint = std::distance(traversal.begin(), (SensorTraversal::const_iterator)points.second);
179 SensitiveDetectorDebugHelper::getInstance().addSimHit(simhit, startPoint, endPoint);
180#endif
181 return simHitIndex;
182 }
183
184 template <class SimHitClass, class TrueHitClass>
186 int trueHitIndex, std::vector<std::pair<unsigned int, float>> simHitIndices)
187 {
188 m_relMCSimHits.add(traversal.getTrackID(), simHitIndices.begin(), simHitIndices.end());
189 //If there is no truehit there are obviously no relations ;)
190 if (trueHitIndex >= 0) {
191 m_relMCTrueHits.add(traversal.getTrackID(), trueHitIndex, traversal.getElectrons());
192 m_relTrueSimHits.add(trueHitIndex, simHitIndices.begin(), simHitIndices.end());
193 }
194 }
195
196 } //VXD Namespace
198} //Belle2 namespace
199
200#endif
static const double ehEnergy
Energy needed to create an electron-hole pair in Si at std.
Definition Const.h:697
Low-level class to create/modify relations between StoreArrays.
@ c_negativeWeight
Flip the sign of the weight to become negative if the original element got re-attributed.
@ c_deleteElement
Delete the whole relation element if the original element got re-attributed.
Class to keep track of the traversal of the sensitive volume for one track.
int getPDGCode() const
get PDG code of the particle
double getElectrons() const
get total number of deposited electrons so far
int getTrackID() const
get Geant4 trackID
std::pair< iterator, iterator > range
Iterator pair for a set of points.
static void registerMCParticleRelation(const std::string &name, RelationArray::EConsolidationAction ignoreAction=RelationArray::c_negativeWeight)
Register an relation involving MCParticles.
Accessor to arrays stored in the data store.
Definition StoreArray.h:113
SensitiveDetectorBase(SensorInfoBase *info)
Constructor.
StepInformation findMidPoint(const SensorTraversal &traversal)
Find the mid-point of the track traversal.
VxdID getSensorID() const
Return the VxdID belonging to this sensitive detector.
std::vector< unsigned int > simplifyEnergyDeposit(const SensorTraversal::range &points)
Simplify the energy deposition profile using Douglas-Peuker-Algorithm We normally force a Geant4 step...
static SensitiveDetectorDebugHelper & getInstance()
Singleton class: get instance.
void addSimHit(const VXDSimHit *simhit, int startPoint, int endPoint)
Write the information normally used when creating SimHits to the entry.
void addTrueHit(const VXDTrueHit *truehit)
Write the information normally used when creating TrueHits to the entry.
int saveSimHit(const SensorTraversal &traversal, const SensorTraversal::range &points) override
Save a SimHit for a track along the given points.
std::array< float, 3 > vecToFloat(const G4ThreeVector &vec)
Convert G4ThreeVector (aka Hep3Vector) to float array to store in TrueHit/SimHit classes.
int saveTrueHit(const SensorTraversal &traversal) override
Save the actual TrueHit for a given sensor traversal information.
void saveRelations(const SensorTraversal &traversal, int trueHitIndex, std::vector< std::pair< unsigned int, float > > simHitIndices) override
Save the relations between MCParticle, TrueHit and SimHits.
SensitiveDetector(VXD::SensorInfoBase *sensorInfo)
Construct instance and take over ownership of the sensorInfo pointer.
Base class to provide Sensor Information for PXD and SVD.
Namespace to provide code needed by both Vertex Detectors, PXD and SVD, and also testbeam telescopes.
Definition GeoCache.h:34
Abstract base class for different kinds of events.
Simple struct to keep information about steps in the sensitive detector.
G4ThreeVector momentum
Step momentum.
G4ThreeVector position
Step position.
double time
timestamp of the step