Belle II Software  release-05-01-25
SensitiveDetector.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - 2014 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Andreas Moll, Zbynek Drasal, Peter Kvasnicka *
7  * Martin Ritter *
8  * *
9  * This software is provided "as is" without any warranty. *
10  **************************************************************************/
11 
12 #pragma once
13 #ifndef VXD_SIMULATION_SENSITIVEDETECTOR_H
14 #define VXD_SIMULATION_SENSITIVEDETECTOR_H
15 
16 #include <vxd/simulation/SensitiveDetectorBase.h>
17 #include <framework/datastore/StoreArray.h>
18 #include <framework/datastore/RelationArray.h>
19 #include <array>
20 
21 #ifdef VXD_SENSITIVEDETECTOR_DEBUG
22 #include <vxd/simulation/SensitiveDetectorDebugHelper.h>
23 #endif
24 
25 namespace Belle2 {
31  namespace VXD {
32 
67  template <class SimHitClass, class TrueHitClass> class SensitiveDetector: public SensitiveDetectorBase {
68  public:
70  explicit SensitiveDetector(VXD::SensorInfoBase* sensorInfo);
71 
72  private:
78  int saveTrueHit(const SensorTraversal& traversal) override;
79 
87  int saveSimHit(const SensorTraversal& traversal,
88  const SensorTraversal::range& points) override;
89 
97  void saveRelations(const SensorTraversal& traversal, int trueHitIndex,
98  std::vector<std::pair<unsigned int, float>> simHitIndices) override;
99 
105  std::array<float, 3> vecToFloat(const G4ThreeVector& vec)
106  {
107  return std::array<float, 3> {{(float)vec.x(), (float)vec.y(), (float)vec.z()}};
108  }
109 
122  };
123 
124  template <class SimHitClass, class TrueHitClass>
126  VXD::SensitiveDetectorBase(sensorInfo)
127  {
128  //Register output collections.
129  //m_mcparticles.isRequired();
130  m_simhits.registerInDataStore();
131  m_truehits.registerInDataStore();
137  }
138 
139  template <class SimHitClass, class TrueHitClass>
141  {
142  //We have the full sensor traversal information, we only lack the midpoint so lets get it
143  StepInformation midPoint = findMidPoint(traversal);
144 
145  //By now everything is done: just convert the position and momenta in float[3] and create a truehit
146  auto posEntry = vecToFloat(traversal.front().position);
147  auto posMidPoint = vecToFloat(midPoint.position);
148  auto posExit = vecToFloat(traversal.back().position);
149  auto momEntry = vecToFloat(traversal.front().momentum);
150  auto momMidPoint = vecToFloat(midPoint.momentum);
151  auto momExit = vecToFloat(traversal.back().momentum);
152  //And we should convert the electrons back in energy ...
153  const double energyDep = traversal.getElectrons() * Const::ehEnergy;
154 
155  //create the truehit
156  int trueHitIndex = m_truehits.getEntries();
157  m_truehits.appendNew(getSensorID(), posEntry.data(), posMidPoint.data(), posExit.data(),
158  momEntry.data(), momMidPoint.data(), momExit.data(), energyDep, midPoint.time);
159 #ifdef VXD_SENSITIVEDETECTOR_DEBUG
160  SensitiveDetectorDebugHelper::getInstance().addTrueHit(m_truehits[trueHitIndex]);
161 #endif
162  return trueHitIndex;
163  }
164 
165  template <class SimHitClass, class TrueHitClass>
167  const SensorTraversal::range& points)
168  {
169  //Convert position to floats here
170  auto posIn = vecToFloat(points.first->position);
171  auto posOut = vecToFloat(points.second->position);
172  auto electronProfile = simplifyEnergyDeposit(points);
173 
174  //Create the simhit
175  int simHitIndex = m_simhits.getEntries();
176  SimHitClass* simhit = m_simhits.appendNew(getSensorID(), traversal.getPDGCode(), points.first->time,
177  posIn.data(), posOut.data());
178  simhit->setEnergyDeposit(electronProfile);
179 #ifdef VXD_SENSITIVEDETECTOR_DEBUG
180  int startPoint = std::distance(traversal.begin(), (SensorTraversal::const_iterator)points.first);
181  int endPoint = std::distance(traversal.begin(), (SensorTraversal::const_iterator)points.second);
182  SensitiveDetectorDebugHelper::getInstance().addSimHit(simhit, startPoint, endPoint);
183 #endif
184  return simHitIndex;
185  }
186 
187  template <class SimHitClass, class TrueHitClass>
189  int trueHitIndex, std::vector<std::pair<unsigned int, float>> simHitIndices)
190  {
191  m_relMCSimHits.add(traversal.getTrackID(), simHitIndices.begin(), simHitIndices.end());
192  //If there is no truehit there are obviously no relations ;)
193  if (trueHitIndex >= 0) {
194  m_relMCTrueHits.add(traversal.getTrackID(), trueHitIndex, traversal.getElectrons());
195  m_relTrueSimHits.add(trueHitIndex, simHitIndices.begin(), simHitIndices.end());
196  }
197  }
198 
199  } //VXD Namespace
201 } //Belle2 namespace
202 
203 #endif
Belle2::RelationArray
Low-level class to create/modify relations between StoreArrays.
Definition: RelationArray.h:72
Belle2::StepInformation::time
double time
timestamp of the step
Definition: SensorTraversal.h:44
Belle2::SensorTraversal
Class to keep track of the traversal of the sensitive volume for one track.
Definition: SensorTraversal.h:54
Belle2::RelationArray::c_negativeWeight
@ c_negativeWeight
Flip the sign of the weight to become negative if the original element got re-attributed.
Definition: RelationArray.h:89
Belle2::VXD::SensitiveDetector::m_relMCSimHits
RelationArray m_relMCSimHits
Relation between MCParticle and SimHits.
Definition: SensitiveDetector.h:117
Belle2::SensorTraversal::getElectrons
double getElectrons() const
get total number of deposited electrons so far
Definition: SensorTraversal.h:79
Belle2::VXD::SensitiveDetectorBase
Base class for Sensitive Detector implementation of PXD and SVD.
Definition: SensitiveDetectorBase.h:47
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::StoreAccessorBase::registerInDataStore
bool registerInDataStore(DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut)
Register the object/array in the DataStore.
Definition: StoreAccessorBase.h:54
Belle2::StepInformation::position
G4ThreeVector position
Step position.
Definition: SensorTraversal.h:38
Belle2::VXD::SensitiveDetectorDebugHelper::addSimHit
void addSimHit(const VXDSimHit *simhit, int startPoint, int endPoint)
Write the information normally used when creating SimHits to the entry.
Definition: SensitiveDetectorDebugHelper.cc:113
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::VXD::SensitiveDetector::m_relMCTrueHits
RelationArray m_relMCTrueHits
Relation between MCParticle and TrueHits.
Definition: SensitiveDetector.h:119
Belle2::VXD::SensitiveDetector::m_simhits
StoreArray< SimHitClass > m_simhits
StoreArray for the SimHits.
Definition: SensitiveDetector.h:111
Belle2::SensorTraversal::getPDGCode
int getPDGCode() const
get PDG code of the particle
Definition: SensorTraversal.h:77
Belle2::VXD::SensitiveDetector::saveRelations
void saveRelations(const SensorTraversal &traversal, int trueHitIndex, std::vector< std::pair< unsigned int, float >> simHitIndices) override
Save the relations between MCParticle, TrueHit and SimHits.
Definition: SensitiveDetector.h:188
Belle2::StepInformation
Simple struct to keep information about steps in the sensitive detector.
Definition: SensorTraversal.h:27
Belle2::SensorTraversal::getTrackID
int getTrackID() const
get Geant4 trackID
Definition: SensorTraversal.h:75
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::RelationArray::c_deleteElement
@ c_deleteElement
Delete the whole relation element if the original element got re-attributed.
Definition: RelationArray.h:91
Belle2::SensorTraversal::range
std::pair< iterator, iterator > range
Iterator pair for a set of points.
Definition: SensorTraversal.h:57
Belle2::VXD::SensitiveDetector::saveTrueHit
int saveTrueHit(const SensorTraversal &traversal) override
Save the actual TrueHit for a given sensor traversal information.
Definition: SensitiveDetector.h:140
Belle2::VXD::SensitiveDetectorDebugHelper::addTrueHit
void addTrueHit(const VXDTrueHit *truehit)
Write the information normally used when creating TrueHits to the entry.
Definition: SensitiveDetectorDebugHelper.cc:88
Belle2::VXD::SensitiveDetector::SensitiveDetector
SensitiveDetector(VXD::SensorInfoBase *sensorInfo)
Construct instance and take over ownership of the sensorInfo pointer.
Definition: SensitiveDetector.h:125
Belle2::VXD::SensitiveDetectorDebugHelper::getInstance
static SensitiveDetectorDebugHelper & getInstance()
Singleton class: get instance.
Definition: SensitiveDetectorDebugHelper.cc:30
Belle2::VXD::SensitiveDetector::vecToFloat
std::array< float, 3 > vecToFloat(const G4ThreeVector &vec)
Convert G4ThreeVector (aka Hep3Vector) to float array to store in TrueHit/SimHit classes.
Definition: SensitiveDetector.h:105
Belle2::VXD::SensitiveDetector::m_truehits
StoreArray< TrueHitClass > m_truehits
StoreArray for the TrueHits.
Definition: SensitiveDetector.h:113
Belle2::StepInformation::momentum
G4ThreeVector momentum
Step momentum.
Definition: SensorTraversal.h:40
Belle2::VXD::SensitiveDetector::saveSimHit
int saveSimHit(const SensorTraversal &traversal, const SensorTraversal::range &points) override
Save a SimHit for a track along the given points.
Definition: SensitiveDetector.h:166
Belle2::VXD::SensitiveDetector::m_mcparticles
StoreArray< MCParticle > m_mcparticles
StoreArray for the MCParticles, needed by relations.
Definition: SensitiveDetector.h:115
Belle2::VXD::SensitiveDetector
Sensitive Detector implementation of PXD and SVD.
Definition: SensitiveDetector.h:67
Belle2::StoreArray< SimHitClass >
Belle2::Const::ehEnergy
static const double ehEnergy
Energy needed to create an electron-hole pair in Si at std.
Definition: Const.h:570
Belle2::VXD::SensitiveDetector::m_relTrueSimHits
RelationArray m_relTrueSimHits
Relation between TrueHits and SimHits.
Definition: SensitiveDetector.h:121