Belle II Software  release-05-02-19
SensitiveDetectorDebugHelper.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 *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <vxd/simulation/SensitiveDetectorDebugHelper.h>
12 #include <framework/gearbox/Unit.h>
13 #include <cassert>
14 
15 namespace Belle2 {
20  namespace VXD {
21 
22  SensitiveDetectorDebugHelper& SensitiveDetectorDebugHelper::getInstance()
23  {
24  static SensitiveDetectorDebugHelper instance;
25  return instance;
26  }
27 
29  {
30  m_file = new TFile("debug-vxdsensitivedetector.root", "RECREATE");
31  m_tree = new TTree("sensortraversal", "Debug info for VXD Sensitive Detector implementation");
32 //Due to laziness, define some macros to create branches with less typing
33 #define ADDBRANCH__(name,var,type) m_tree->Branch(#name,&m_info.var,#type)
34 #define ADDBARRAY__(name,var,n,type) ADDBRANCH__(name,var,name[n]/type)
35 #define ADDBRANCH(x,type) ADDBRANCH__(x,x,x/type)
36 #define ADDBARRAY(x,n,type) ADDBARRAY__(x,x,n,type)
37  //and create all branches
38  ADDBRANCH(sensorID, I);
39  ADDBRANCH(pdg, I);
40  ADDBRANCH(contained, I);
41  ADDBRANCH(primary, I);
42  ADDBRANCH(length, D);
43  ADDBRANCH(stepN, I);
44  ADDBARRAY(stepInfo, stepN, D);
45  ADDBRANCH(trueN, I);
46  ADDBARRAY(trueInfo, trueN, D);
47  ADDBRANCH(edepN, I);
48  ADDBARRAY(edepInfo, edepN, D);
49  ADDBRANCH(simhitN, I);
50  ADDBARRAY(simhitInfo, simhitN, D);
51  }
53  {
54  m_info.sensorID = sensorID;
55  m_info.pdg = traversal.getPDGCode();
56  m_info.contained = traversal.isContained();
57  m_info.primary = traversal.isPrimary();
58  m_info.length = traversal.getLength();
59  m_info.stepN = 0;
61  m_info.edepN = 0;
62  m_info.simhitN = 0;
63 
64  //add all steps if there is enough space
65  assert(traversal.size() <= MAX_STEPS);
66  for (const StepInformation& step : traversal) {
67  m_info.stepInfo[m_info.stepN + 0] = step.position.x();
68  m_info.stepInfo[m_info.stepN + 1] = step.position.y();
69  m_info.stepInfo[m_info.stepN + 2] = step.position.z();
70  m_info.stepInfo[m_info.stepN + 3] = step.momentum.x();
71  m_info.stepInfo[m_info.stepN + 4] = step.momentum.y();
72  m_info.stepInfo[m_info.stepN + 5] = step.momentum.z();
73  m_info.stepInfo[m_info.stepN + 6] = step.electrons;
74  m_info.stepInfo[m_info.stepN + 7] = step.time;
75  m_info.stepInfo[m_info.stepN + 8] = step.length;
77  }
78  }
79 
81  {
82  m_info.trueInfo[0] = truehit->getEntryU();
83  m_info.trueInfo[1] = truehit->getEntryV();
84  m_info.trueInfo[2] = truehit->getEntryW();
85  m_info.trueInfo[3] = truehit->getU();
86  m_info.trueInfo[4] = truehit->getV();
87  m_info.trueInfo[5] = truehit->getW();
88  m_info.trueInfo[6] = truehit->getExitU();
89  m_info.trueInfo[7] = truehit->getExitV();
90  m_info.trueInfo[8] = truehit->getExitW();
91  m_info.trueInfo[9] = truehit->getEntryMomentum()[0];
92  m_info.trueInfo[10] = truehit->getEntryMomentum()[1];
93  m_info.trueInfo[11] = truehit->getEntryMomentum()[2];
94  m_info.trueInfo[12] = truehit->getMomentum()[0];
95  m_info.trueInfo[13] = truehit->getMomentum()[1];
96  m_info.trueInfo[14] = truehit->getMomentum()[2];
97  m_info.trueInfo[15] = truehit->getExitMomentum()[0];
98  m_info.trueInfo[16] = truehit->getExitMomentum()[1];
99  m_info.trueInfo[17] = truehit->getExitMomentum()[2];
100  m_info.trueInfo[18] = truehit->getEnergyDep();
101  m_info.trueInfo[19] = truehit->getGlobalTime();
103  }
104 
105  void SensitiveDetectorDebugHelper::addSimHit(const VXDSimHit* simhit, int startPoint, int endPoint)
106  {
107  assert(m_info.simhitN + SIZE_HITS < MAX_HITS * SIZE_HITS);
108  m_info.simhitInfo[m_info.simhitN + 0] = simhit->getPosIn()[0];
109  m_info.simhitInfo[m_info.simhitN + 1] = simhit->getPosIn()[1];
110  m_info.simhitInfo[m_info.simhitN + 2] = simhit->getPosIn()[2];
111  m_info.simhitInfo[m_info.simhitN + 3] = simhit->getPosOut()[0];
112  m_info.simhitInfo[m_info.simhitN + 4] = simhit->getPosOut()[1];
113  m_info.simhitInfo[m_info.simhitN + 5] = simhit->getPosOut()[2];
114  m_info.simhitInfo[m_info.simhitN + 6] = simhit->getGlobalTime();
115  m_info.simhitInfo[m_info.simhitN + 7] = startPoint;
116  m_info.simhitInfo[m_info.simhitN + 8] = endPoint;
117  //Add optimized energy deposition profile
119  for (auto step : simhit->getElectronProfile()) {
120  m_info.edepInfo[m_info.edepN] = step.first;
121  ++m_info.edepN;
122  m_info.edepInfo[m_info.edepN] = step.second;
123  ++m_info.edepN;
124  }
126  //And also add electron profile sampled by distance
127  for (auto step : simhit->getElectronsConstantDistance(50 * Unit::um)) {
128  m_info.edepInfo[m_info.edepN] = step.first;
129  ++m_info.edepN;
130  m_info.edepInfo[m_info.edepN] = step.second;
131  ++m_info.edepN;
132  }
134  //And sampled by electrons
135  for (auto step : simhit->getElectronsConstantNumber(1000)) {
136  m_info.edepInfo[m_info.edepN] = step.first;
137  ++m_info.edepN;
138  m_info.edepInfo[m_info.edepN] = step.second;
139  ++m_info.edepN;
140  }
143  }
144  }
146 } //Belle2 namespace
Belle2::VXD::SensitiveDetectorDebugHelper::MAX_STEPS
@ MAX_STEPS
assume maximum number of steps per sensor traversal
Definition: SensitiveDetectorDebugHelper.h:77
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::VXD::SensitiveDetectorDebugHelper::MAX_HITS
@ MAX_HITS
assume maximum number of simhits per sensor traversal
Definition: SensitiveDetectorDebugHelper.h:79
Belle2::VXDTrueHit::getExitU
float getExitU() const
Return local u coordinate of hit at the endpoint of the track.
Definition: VXDTrueHit.h:88
Belle2::SensorTraversal
Class to keep track of the traversal of the sensitive volume for one track.
Definition: SensorTraversal.h:54
Belle2::VXD::SensitiveDetectorDebugHelper::info::trueN
int trueN
number of values used in trueInfo
Definition: SensitiveDetectorDebugHelper.h:112
Belle2::VXD::SensitiveDetectorDebugHelper::SensitiveDetectorDebugHelper
SensitiveDetectorDebugHelper()
Singleton class: private constructor.
Definition: SensitiveDetectorDebugHelper.cc:36
Belle2::VXDTrueHit::getEntryV
float getEntryV() const
Return local v coordinate of the start point of the track.
Definition: VXDTrueHit.h:84
Belle2::VXD::SensitiveDetectorDebugHelper::m_file
TFile * m_file
ROOT File to write information to.
Definition: SensitiveDetectorDebugHelper.h:71
Belle2::VXDTrueHit
Class VXDTrueHit - Records of tracks that either enter or leave the sensitive volume.
Definition: VXDTrueHit.h:38
Belle2::SensorTraversal::getLength
double getLength() const
get flight length so far
Definition: SensorTraversal.h:81
Belle2::VXDTrueHit::getEntryW
float getEntryW() const
Return local w coordinate of the start point of the track.
Definition: VXDTrueHit.h:86
Belle2::VXD::SensitiveDetectorDebugHelper::info::length
double length
track length
Definition: SensitiveDetectorDebugHelper.h:104
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::VXD::SensitiveDetectorDebugHelper::info::stepInfo
double stepInfo[SIZE_STEP *MAX_STEPS]
values for the step points
Definition: SensitiveDetectorDebugHelper.h:109
Belle2::VXD::SensitiveDetectorDebugHelper::SIZE_HITS
@ SIZE_HITS
number of values per simhit
Definition: SensitiveDetectorDebugHelper.h:88
Belle2::VXD::SensitiveDetectorDebugHelper::SIZE_TRUE
@ SIZE_TRUE
number of values per truehit
Definition: SensitiveDetectorDebugHelper.h:86
Belle2::VXDSimHit::getElectronsConstantNumber
std::vector< std::pair< float, float > > getElectronsConstantNumber(double electronsPerStep) const
Get the electron deposition with constant number of electrons between sampling points.
Definition: VXDSimhit.cc:68
Belle2::SensorTraversal::getPDGCode
int getPDGCode() const
get PDG code of the particle
Definition: SensorTraversal.h:77
Belle2::VXD::SensitiveDetectorDebugHelper::info::sensorID
int sensorID
id of the sensor for the current track
Definition: SensitiveDetectorDebugHelper.h:96
Belle2::VXD::SensitiveDetectorDebugHelper::m_tree
TTree * m_tree
ROOT Tree to write information to.
Definition: SensitiveDetectorDebugHelper.h:73
Belle2::StepInformation
Simple struct to keep information about steps in the sensitive detector.
Definition: SensorTraversal.h:27
Belle2::VXD::SensitiveDetectorDebugHelper::info::edepN
int edepN
number of values used in edepInfo
Definition: SensitiveDetectorDebugHelper.h:117
Belle2::VXDSimHit::getPosOut
TVector3 getPosOut() const
Return the end point of the electron deposition in local coordinates.
Definition: VXDSimHit.h:74
Belle2::VXD::SensitiveDetectorDebugHelper::SIZE_STEP
@ SIZE_STEP
number of values per step
Definition: SensitiveDetectorDebugHelper.h:84
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VXDTrueHit::getV
float getV() const
Return local v coordinate of hit.
Definition: VXDTrueHit.h:78
Belle2::VXDTrueHit::getW
float getW() const
Return local w coordinate of hit.
Definition: VXDTrueHit.h:80
Belle2::VXD::SensitiveDetectorDebugHelper::info::edepInfo
double edepInfo[SIZE_EDEP *MAX_EDEP]
values for the energy deposition
Definition: SensitiveDetectorDebugHelper.h:119
Belle2::VXDSimHit::getPosIn
TVector3 getPosIn() const
Return the start point of the electron deposition in local coordinates.
Definition: VXDSimHit.h:72
Belle2::VXDTrueHit::getExitV
float getExitV() const
Return local v coordinate of hit at the endpoint of the track.
Definition: VXDTrueHit.h:90
Belle2::VXDSimHit::getElectronProfile
std::vector< std::pair< float, float > > getElectronProfile() const
Get the decoded electron profile.
Definition: VXDSimhit.cc:23
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::SensitiveDetectorDebugHelper::startTraversal
void startTraversal(VxdID sensorID, const SensorTraversal &traversal)
Start writing a new sensor traversal: Geant4 steps will be added to the ROOT entry.
Definition: SensitiveDetectorDebugHelper.cc:60
Belle2::VXD::SensitiveDetectorDebugHelper::info::trueInfo
double trueInfo[SIZE_TRUE]
values for the truehit
Definition: SensitiveDetectorDebugHelper.h:114
Belle2::SensorTraversal::isPrimary
bool isPrimary() const
return whether the track belongs to a primary particle
Definition: SensorTraversal.h:85
Belle2::VXDSimHit::getGlobalTime
float getGlobalTime() const override
Return the time of the electron deposition.
Definition: VXDSimHit.h:80
Belle2::VXD::SensitiveDetectorDebugHelper::getInstance
static SensitiveDetectorDebugHelper & getInstance()
Singleton class: get instance.
Definition: SensitiveDetectorDebugHelper.cc:30
Belle2::VXD::SensitiveDetectorDebugHelper::m_info
struct Belle2::VXD::SensitiveDetectorDebugHelper::info m_info
object to store all variables
Belle2::VXD::SensitiveDetectorDebugHelper::info::stepN
int stepN
number of values used in stepInfo
Definition: SensitiveDetectorDebugHelper.h:107
Belle2::VXDTrueHit::getEntryMomentum
TVector3 getEntryMomentum() const
Return momentum at the start point of the track.
Definition: VXDTrueHit.h:100
Belle2::VXDTrueHit::getEntryU
float getEntryU() const
Return local u coordinate of hit when entering silicon.
Definition: VXDTrueHit.h:82
Belle2::VXDTrueHit::getEnergyDep
float getEnergyDep() const
Return energy deposited during traversal of sensor.
Definition: VXDTrueHit.h:94
Belle2::VXDTrueHit::getExitMomentum
TVector3 getExitMomentum() const
Return momentum at the endpoint of the track.
Definition: VXDTrueHit.h:102
Belle2::VXD::SensitiveDetectorDebugHelper::info::pdg
int pdg
pdg code for the current track
Definition: SensitiveDetectorDebugHelper.h:98
Belle2::VXDTrueHit::getU
float getU() const
Return local u coordinate of hit.
Definition: VXDTrueHit.h:76
Belle2::VXDSimHit::getElectronsConstantDistance
std::vector< std::pair< float, float > > getElectronsConstantDistance(double length) const
Get the electron deposition along constant stepsize.
Definition: VXDSimhit.cc:34
Belle2::VXDTrueHit::getExitW
float getExitW() const
Return local w coordinate of hit at the endpoint of the track.
Definition: VXDTrueHit.h:92
Belle2::VXD::SensitiveDetectorDebugHelper::info::primary
int primary
whether or not the track is from a primary particle
Definition: SensitiveDetectorDebugHelper.h:102
Belle2::VXDTrueHit::getGlobalTime
float getGlobalTime() const
Return the time when the track reached its midpoint.
Definition: VXDTrueHit.h:96
Belle2::VXDTrueHit::getMomentum
TVector3 getMomentum() const
Return momentum at the midpoint of the track.
Definition: VXDTrueHit.h:98
Belle2::VXD::SensitiveDetectorDebugHelper::info::simhitN
int simhitN
number of values used in simhitInfo
Definition: SensitiveDetectorDebugHelper.h:122
Belle2::SensorTraversal::isContained
bool isContained() const
return whether the track was contained in the volume so far
Definition: SensorTraversal.h:83
Belle2::VXD::SensitiveDetectorDebugHelper::info::simhitInfo
double simhitInfo[SIZE_HITS *MAX_HITS]
values for the simhits
Definition: SensitiveDetectorDebugHelper.h:124
Belle2::VXD::SensitiveDetectorDebugHelper::info::contained
int contained
whether or not the track was contained in the senitive volume
Definition: SensitiveDetectorDebugHelper.h:100
Belle2::Unit::um
static const double um
[micrometers]
Definition: Unit.h:81
Belle2::VXDSimHit
Class VXDSimHit - Geant4 simulated hit for the VXD.
Definition: VXDSimHit.h:32