Belle II Software  release-05-02-19
SVDPerformanceTTreeModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giulia Casarosa *
7  * *
8  * Module used to create a TTree to study SVD clusters, *
9  * genfit unbiased residuals and many other properties related *
10  * to the track they belong to. *
11  * *
12  * This software is provided "as is" without any warranty. *
13  **************************************************************************/
14 
15 #include <svd/modules/svdPerformance/SVDPerformanceTTreeModule.h>
16 #include <framework/datastore/StoreArray.h>
17 #include <framework/gearbox/Unit.h>
18 #include <framework/core/Environment.h>
19 #include <svd/dataobjects/SVDCluster.h>
20 #include <svd/dataobjects/SVDRecoDigit.h>
21 #include <svd/dataobjects/SVDTrueHit.h>
22 #include <svd/dataobjects/SVDEventInfo.h>
23 #include <vxd/dataobjects/VxdID.h>
24 #include <vxd/geometry/GeoCache.h>
25 #include <vxd/geometry/SensorInfoBase.h>
26 #include <tracking/dataobjects/RecoTrack.h>
27 #include <tracking/dataobjects/RecoHitInformation.h>
28 #include <genfit/TrackPoint.h>
29 #include <TVector3.h>
30 #include <math.h>
31 #include <iostream>
32 #include <mdst/dataobjects/Track.h>
33 #include <mdst/dataobjects/TrackFitResult.h>
34 
35 
36 using namespace Belle2;
37 using namespace std;
38 
39 
40 //-----------------------------------------------------------------
41 // Register the Module
42 //-----------------------------------------------------------------
43 REG_MODULE(SVDPerformanceTTree)
44 
45 //-----------------------------------------------------------------
46 // Implementation
47 //-----------------------------------------------------------------
48 
50 {
51  //Set module properties
52  setDescription("The module is used to create a TTree to study SVD clusters, genfit unbiased residuals and many other properties related to the track they belong to.");
53  //Parameter to take only specific RecoTracks as input
54  addParam("outputFileName", m_rootFileName, "Name of output root file.", std::string("SVDPerformanceTTree.root"));
55  addParam("recoTracksStoreArrayName", m_recoTracksStoreArrayName, "StoreArray name of the input and output RecoTracks.",
56  m_recoTracksStoreArrayName);
57 }
58 
60 {
61  StoreArray<RecoTrack> recoTracks(m_recoTracksStoreArrayName);
62  recoTracks.isOptional();
63 
64  m_rootFilePtr = new TFile(m_rootFileName.c_str(), "RECREATE");
65 
66  //Tree for SVD u overlapping clusters
67  m_t_U = new TTree("t_U", "Tree for SVD u-clusters");
68  m_t_U->Branch("svdClSNR", &m_svdClSNR, "svdClSNR/F");
69  m_t_U->Branch("svdClCharge", &m_svdClCharge, "svdClCharge/F");
70  m_t_U->Branch("svdStripCharge", &m_svdStripCharge);
71  m_t_U->Branch("svdClTime", &m_svdClTime, "svdClTime/F");
72  m_t_U->Branch("svdStripTime", &m_svdStripTime);
73  m_t_U->Branch("svdRes", &m_svdRes, "svdRes/F");
74  m_t_U->Branch("svdClIntStrPos", &m_svdClIntStrPos, "svdClIntStrPos/F");
75  m_t_U->Branch("svdClPos", &m_svdClPos, "svdClPos/F");
76  m_t_U->Branch("svdClPosErr", &m_svdClPosErr, "svdClPosErr/F");
77  m_t_U->Branch("svdTruePos", &m_svdTruePos, "svdTruePos/F");
78  m_t_U->Branch("svdClPhi", &m_svdClPhi, "svdClPhi/F");
79  m_t_U->Branch("svdClZ", &m_svdClZ, "svdClZ/F");
80  m_t_U->Branch("svdTrkd0", &m_svdTrkd0, "svdTrkd0/F");
81  m_t_U->Branch("svdTrkz0", &m_svdTrkz0, "svdTrkz0/F");
82  m_t_U->Branch("svdTrkpT", &m_svdTrkpT, "svdTrkpT/F");
83  m_t_U->Branch("svdTrkpCM", &m_svdTrkpCM, "svdTrkpCM/F");
84  m_t_U->Branch("svdTrkTraversedLength", &m_svdTrkTraversedLength, "svdTrkTraversedLength/F");
85  m_t_U->Branch("svdTrkPXDHits", &m_svdTrkPXDHits, "svdTrkPXDHits/i");
86  m_t_U->Branch("svdTrkSVDHits", &m_svdTrkSVDHits, "svdTrkSVDHits/i");
87  m_t_U->Branch("svdTrkCDCHits", &m_svdTrkCDCHits, "svdTrkCDCHits/i");
88  m_t_U->Branch("svdTrkPos", &m_svdTrkPos, "svdTrkPos/F");
89  m_t_U->Branch("svdTrkPosOS", &m_svdTrkPosOS, "svdTrkPosOS/F");
90  m_t_U->Branch("svdTrkPosErr", &m_svdTrkPosErr, "svdTrkPosErr/F");
91  m_t_U->Branch("svdTrkPosErrOS", &m_svdTrkPosErrOS, "svdTrkPosErrOS/F");
92  m_t_U->Branch("svdTrkQoP", &m_svdTrkQoP, "svdTrkQoP/F");
93  m_t_U->Branch("svdTrkPrime", &m_svdTrkPrime, "svdTrkPrime/F");
94  m_t_U->Branch("svdTrkPrimeOS", &m_svdTrkPrimeOS, "svdTrkPrimeOS/F");
95  m_t_U->Branch("svdTrkPosUnbiased", &m_svdTrkPosUnbiased, "svdTrkPosUnbiased/F");
96  m_t_U->Branch("svdTrkPosErrUnbiased", &m_svdTrkPosErrUnbiased, "svdTrkPosErrUnbiased/F");
97  m_t_U->Branch("svdTrkQoPUnbiased", &m_svdTrkQoPUnbiased, "svdTrkQoPUnbiased/F");
98  m_t_U->Branch("svdTrkPrimeUnbiased", &m_svdTrkPrimeUnbiased, "svdTrkPrimeUnbiased/F");
99  m_t_U->Branch("svdLayer", &m_svdLayer, "svdLayer/i");
100  m_t_U->Branch("svdLadder", &m_svdLadder, "svdLadder/i");
101  m_t_U->Branch("svdSensor", &m_svdSensor, "svdSensor/i");
102  m_t_U->Branch("svdSize", &m_svdSize, "svdSize/i");
103  m_t_U->Branch("svdTB", &m_svdTB, "svdTB/i");
104  //Tree for SVD v overlapping clusters
105  m_t_V = new TTree("t_V", "Tree for M_SVD v-clusters");
106  m_t_V->Branch("svdClSNR", &m_svdClSNR, "svdClSNR/F");
107  m_t_V->Branch("svdClCharge", &m_svdClCharge, "svdClCharge/F");
108  m_t_V->Branch("svdStripCharge", &m_svdStripCharge);
109  m_t_V->Branch("svdClTime", &m_svdClTime, "svdClTime/F");
110  m_t_V->Branch("svdStripTime", &m_svdStripTime);
111  m_t_V->Branch("svdRes", &m_svdRes, "svdRes/F");
112  m_t_V->Branch("svdClIntStrPos", &m_svdClIntStrPos, "svdClIntStrPos/F");
113  m_t_V->Branch("svdClPos", &m_svdClPos, "svdClPos/F");
114  m_t_V->Branch("svdClPosErr", &m_svdClPosErr, "svdClPosErr/F");
115  m_t_V->Branch("svdTruePos", &m_svdTruePos, "svdTruePos/F");
116  m_t_V->Branch("svdClPhi", &m_svdClPhi, "svdClPhi/F");
117  m_t_V->Branch("svdClZ", &m_svdClZ, "svdClZ/F");
118  m_t_V->Branch("svdTrkd0", &m_svdTrkd0, "svdTrkd0/F");
119  m_t_V->Branch("svdTrkz0", &m_svdTrkz0, "svdTrkz0/F");
120  m_t_V->Branch("svdTrkpT", &m_svdTrkpT, "svdTrkpT/F");
121  m_t_V->Branch("svdTrkpCM", &m_svdTrkpCM, "svdTrkpCM/F");
122  m_t_V->Branch("svdTrkTraversedLength", &m_svdTrkTraversedLength, "svdTrkTraversedLength/F");
123  m_t_V->Branch("svdTrkPXDHits", &m_svdTrkPXDHits, "svdTrkPXDHits/i");
124  m_t_V->Branch("svdTrkSVDHits", &m_svdTrkSVDHits, "svdTrkSVDHits/i");
125  m_t_V->Branch("svdTrkCDCHits", &m_svdTrkCDCHits, "svdTrkCDCHits/i");
126  m_t_V->Branch("svdTrkPos", &m_svdTrkPos, "svdTrkPos/F");
127  m_t_V->Branch("svdTrkPosOS", &m_svdTrkPosOS, "svdTrkPosOS/F");
128  m_t_V->Branch("svdTrkPosErr", &m_svdTrkPosErr, "svdTrkPosErr/F");
129  m_t_V->Branch("svdTrkPosErrOS", &m_svdTrkPosErrOS, "svdTrkPosErrOS/F");
130  m_t_V->Branch("svdTrkQoP", &m_svdTrkQoP, "svdTrkQoP/F");
131  m_t_V->Branch("svdTrkPrime", &m_svdTrkPrime, "svdTrkPrime/F");
132  m_t_V->Branch("svdTrkPrimeOS", &m_svdTrkPrimeOS, "svdTrkPrimeOS/F");
133  m_t_V->Branch("svdTrkPosUnbiased", &m_svdTrkPosUnbiased, "svdTrkPosUnbiased/F");
134  m_t_V->Branch("svdTrkPosErrUnbiased", &m_svdTrkPosErrUnbiased, "svdTrkPosErrUnbiased/F");
135  m_t_V->Branch("svdTrkQoPUnbiased", &m_svdTrkQoPUnbiased, "svdTrkQoPUnbiased/F");
136  m_t_V->Branch("svdTrkPrimeUnbiased", &m_svdTrkPrimeUnbiased, "svdTrkPrimeUnbiased/F");
137  m_t_V->Branch("svdLayer", &m_svdLayer, "svdLayer/i");
138  m_t_V->Branch("svdLadder", &m_svdLadder, "svdLadder/i");
139  m_t_V->Branch("svdSensor", &m_svdSensor, "svdSensor/i");
140  m_t_V->Branch("svdSize", &m_svdSize, "svdSize/i");
141  m_t_V->Branch("svdTB", &m_svdTB, "svdTB/i");
142 
143 }
144 
146 {
147 
148  //first check SVDEventInfo name
149  StoreObjPtr<SVDEventInfo> temp_eventinfo("SVDEventInfo");
150  std::string m_svdEventInfoName = "SVDEventInfo";
151  if (!temp_eventinfo.isValid())
152  m_svdEventInfoName = "SVDEventInfoSim";
153  StoreObjPtr<SVDEventInfo> eventinfo(m_svdEventInfoName);
154  if (!eventinfo) B2ERROR("No SVDEventInfo!");
155  m_svdTB = eventinfo->getModeByte().getTriggerBin();
156 
157  bool isMC = Environment::Instance().isMC();
158 
160  StoreArray<RecoTrack> recoTracks(m_recoTracksStoreArrayName);
161  for (const auto& trk : recoTracks) {
162  if (! trk.wasFitSuccessful()) {
163  continue;
164  }
165  int pionCode = 211;
166 
167  RelationVector<Track> theTK = DataStore::getRelationsWithObj<Track>(&trk);
168 
169  if (theTK.size() == 0) {
170  continue;
171  }
172 
173 
174  const TrackFitResult* tfr = theTK[0]->getTrackFitResult(Const::ChargedStable(pionCode));
175  if (tfr) {
176  m_svdTrkd0 = tfr->getD0();
177  m_svdTrkz0 = tfr->getZ0();
178  m_svdTrkpT = tfr->getMomentum().Perp();
179  TLorentzVector pStar = tfr->get4Momentum();
180  pStar.Boost(0, 0, 3. / 11);
181  m_svdTrkpCM = pStar.P();
182  }
183 
184 
185  const vector<SVDCluster* > svdClusters = trk.getSVDHitList();
186  B2DEBUG(40, "FITTED TRACK: NUMBER OF SVD HITS = " << svdClusters.size());
187 
188  m_svdTrkPXDHits = (trk.getPXDHitList()).size();
189  m_svdTrkSVDHits = (trk.getSVDHitList()).size();
190  m_svdTrkCDCHits = (trk.getCDCHitList()).size();
191 
192  for (unsigned int i = 0; i < svdClusters.size(); i++) {
193 
194  const SVDCluster* svd_1 = svdClusters[i];
195 
196  //get true hits, used only if isMC
197  RelationVector<SVDTrueHit> trueHit_1 = DataStore::getRelationsWithObj<SVDTrueHit>(svd_1);
198 
199  const RecoHitInformation* infoSVD_1 = trk.getRecoHitInformation(svd_1);
200  if (!infoSVD_1) {
201  continue;
202  }
203  const auto* hitTrackPoint_1 = trk.getCreatedTrackPoint(infoSVD_1);
204  const auto* fittedResult_1 = hitTrackPoint_1->getFitterInfo();
205  if (!fittedResult_1) {
206  continue;
207  }
208  const VxdID svd_id_1 = svd_1->getSensorID();
209  const unsigned short svd_Layer_1 = svd_id_1.getLayerNumber();
210  const unsigned short svd_Ladder_1 = svd_id_1.getLadderNumber();
211  const unsigned short svd_Sensor_1 = svd_id_1.getSensorNumber();
212 
213  try {
214  const TVectorD resUnBias_1 = fittedResult_1->getResidual(0, false).getState();
215  genfit::MeasuredStateOnPlane state_unbiased = fittedResult_1->getFittedState(false);
216  const TVectorD& svd_predIntersect_unbiased = state_unbiased.getState();
217  const TMatrixDSym& covMatrix_unbiased = state_unbiased.getCov();
218  genfit::MeasuredStateOnPlane state_1 = trk.getMeasuredStateOnPlaneFromRecoHit(infoSVD_1);
219  const TVectorD& svd_predIntersect_1 = state_1.getState();
220  const TMatrixDSym& covMatrix_1 = state_1.getCov();
221 
222  if (svd_1->isUCluster()) {
223 
224  const int strips_1 = svd_1->getSize();
225 
226  const double res_U_1 = resUnBias_1.GetMatrixArray()[0] * Unit::convertValueToUnit(1.0, "um");
227  const TVector3 svdLocal_1(svd_1->getPosition(), svd_predIntersect_1[4], 0.);
228  const VXD::SensorInfoBase& svdSensor_1 = geo.get(svd_id_1);
229  const TVector3& svdGlobal_1 = svdSensor_1.pointToGlobal(svdLocal_1);
230  double svdPhi_1 = atan2(svdGlobal_1(1), svdGlobal_1(0));
231  double svdZ_1 = svdGlobal_1(2);
232 
233  //Fill SVD tree for u-overlaps if required by the user
234  m_svdRes = res_U_1;
235  m_svdClTime = svd_1->getClsTime();
236  m_svdClSNR = svd_1->getSNR();
237  m_svdClCharge = svd_1->getCharge();
238  m_svdClPos = svd_1->getPosition();
239  m_svdClPosErr = svd_1->getPositionSigma();
240  if (isMC && trueHit_1.size() > 0)
241  m_svdTruePos = trueHit_1[0]->getU();
242  else
243  m_svdTruePos = -99;
244  m_svdClPhi = svdPhi_1;
245  m_svdClZ = svdZ_1;
246  m_svdTrkPos = svd_predIntersect_1[3];
247  m_svdTrkPosOS = svd_predIntersect_1[4];
248  m_svdTrkPosErr = sqrt(covMatrix_1[3][3]);
249  m_svdTrkPosErrOS = sqrt(covMatrix_1[4][4]);
250  m_svdTrkQoP = svd_predIntersect_1[0];
251  m_svdTrkPrime = svd_predIntersect_1[1];
252  m_svdTrkPrimeOS = svd_predIntersect_1[2];
253  m_svdTrkTraversedLength = svdSensor_1.getThickness() * sqrt(1 + m_svdTrkPrimeOS * m_svdTrkPrimeOS + m_svdTrkPrime * m_svdTrkPrime);
254  m_svdTrkPosUnbiased = svd_predIntersect_unbiased[3];
255  m_svdTrkPosErrUnbiased = sqrt(covMatrix_unbiased[3][3]);
256  m_svdTrkQoPUnbiased = svd_predIntersect_unbiased[0];
257  m_svdTrkPrimeUnbiased = svd_predIntersect_unbiased[1];
258  m_svdLayer = svd_Layer_1;
259  m_svdLadder = svd_Ladder_1;
260  m_svdSensor = svd_Sensor_1;
261  m_svdSize = strips_1;
262 
263  float pitch = 50e-4;
264  float halfLength = 1.92;
265  if (m_svdLayer > 3) {
266  pitch = 75e-4;
267  halfLength = 2.88;
268  }
269  m_svdClIntStrPos = fmod(m_svdClPos + halfLength, pitch) / pitch;
270 
271  m_svdStripCharge.clear();
272  m_svdStripTime.clear();
273  //retrieve relations and set strip charges and times
274  RelationVector<SVDRecoDigit> theRecoDigits = DataStore::getRelationsWithObj<SVDRecoDigit>(svd_1);
275  if ((theRecoDigits.size() != m_svdSize) && (m_svdSize != 128)) //virtual cluster
276  B2ERROR(" Inconsistency with cluster size! # recoDigits = " << theRecoDigits.size() << " != " << m_svdSize << " cluster size");
277 
278  if (!(m_svdSize == 128 && theRecoDigits.size() == 0))
279  for (unsigned int d = 0; d < m_svdSize; d++) {
280  m_svdStripCharge.push_back(theRecoDigits[d]->getCharge());
281  m_svdStripTime.push_back(theRecoDigits[d]->getTime());
282  }
283 
284  m_t_U->Fill();
285 
286  } else {
287  const int strips_1 = svd_1->getSize();
288  const double res_V_1 = resUnBias_1.GetMatrixArray()[0] * Unit::convertValueToUnit(1.0, "um");
289  const TVector3 svdLocal_1(svd_predIntersect_1[3], svd_1->getPosition(), 0.);
290  const VXD::SensorInfoBase& svdSensor_1 = geo.get(svd_id_1);
291  const TVector3& svdGlobal_1 = svdSensor_1.pointToGlobal(svdLocal_1);
292  double svdPhi_1 = atan2(svdGlobal_1(1), svdGlobal_1(0));
293  double svdZ_1 = svdGlobal_1(2);
294 
295  m_svdRes = res_V_1;
296  m_svdClTime = svd_1->getClsTime();
297  m_svdClSNR = svd_1->getSNR();
298  m_svdClCharge = svd_1->getCharge();
299  m_svdClPos = svd_1->getPosition();
300  m_svdClPosErr = svd_1->getPositionSigma();
301  if (isMC && trueHit_1.size() > 0)
302  m_svdTruePos = trueHit_1[0]->getV();
303  else
304  m_svdTruePos = -99;
305  m_svdClPhi = svdPhi_1;
306  m_svdClZ = svdZ_1;
307  m_svdTrkPos = svd_predIntersect_1[4];
308  m_svdTrkPosOS = svd_predIntersect_1[3];
309  m_svdTrkPosErr = sqrt(covMatrix_1[4][4]);
310  m_svdTrkPosErrOS = sqrt(covMatrix_1[3][3]);
311  m_svdTrkQoP = svd_predIntersect_1[0];
312  m_svdTrkPrime = svd_predIntersect_1[2];
313  m_svdTrkPrimeOS = svd_predIntersect_1[1];
314  m_svdTrkTraversedLength = svdSensor_1.getThickness() * sqrt(1 + m_svdTrkPrimeOS * m_svdTrkPrimeOS + m_svdTrkPrime * m_svdTrkPrime);
315  m_svdTrkPosUnbiased = svd_predIntersect_unbiased[4];
316  m_svdTrkPosErrUnbiased = sqrt(covMatrix_unbiased[4][4]);
317  m_svdTrkQoPUnbiased = svd_predIntersect_unbiased[0];
318  m_svdTrkPrimeUnbiased = svd_predIntersect_unbiased[2];
319  m_svdLayer = svd_Layer_1;
320  m_svdLadder = svd_Ladder_1;
321  m_svdSensor = svd_Sensor_1;
322  m_svdSize = strips_1;
323 
324  float pitch = 160e-4;
325  float halfLength = 6.144;
326  if (m_svdLayer > 3) {
327  pitch = 240e-4;
328  }
329  m_svdClIntStrPos = fmod(m_svdClPos + halfLength, pitch) / pitch;
330 
331  m_svdStripCharge.clear();
332  m_svdStripTime.clear();
333 
334 
335  //retrieve relations and set strip charges and times
336  RelationVector<SVDRecoDigit> theRecoDigits = DataStore::getRelationsWithObj<SVDRecoDigit>(svd_1);
337  if ((theRecoDigits.size() != m_svdSize) && (m_svdSize != 128)) //virtual cluster
338  B2ERROR(" Inconsistency with cluster size! # recoDigits = " << theRecoDigits.size() << " != " << m_svdSize << " cluster size");
339 
340  if (!(m_svdSize == 128 && theRecoDigits.size() == 0))
341  for (unsigned int d = 0; d < m_svdSize; d++) {
342  m_svdStripCharge.push_back(theRecoDigits[d]->getCharge());
343  m_svdStripTime.push_back(theRecoDigits[d]->getTime());
344  }
345 
346  m_t_V->Fill();
347  }
348  } catch (...) {
349  B2INFO("oops...something went wrong in getting the unbiased state, skipping this cluster.");
350  continue;
351  }
352 
353  }
354  }
355 }
356 
357 
359 {
360 
361  if (m_rootFilePtr != NULL) {
362  m_rootFilePtr->cd();
363  m_t_U->Write();
364  m_t_V->Write();
365  m_rootFilePtr->Close();
366  }
367 }
Belle2::RelationVector::size
size_t size() const
Get number of relations.
Definition: RelationVector.h:98
Belle2::SVDCluster::getCharge
float getCharge() const
Get collected charge.
Definition: SVDCluster.h:152
Belle2::TrackFitResult::get4Momentum
TLorentzVector get4Momentum() const
Getter for the 4Momentum at the closest approach of the track in the r/phi projection.
Definition: TrackFitResult.h:125
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
genfit::MeasuredStateOnPlane
#StateOnPlane with additional covariance matrix.
Definition: MeasuredStateOnPlane.h:39
Belle2::TrackFitResult::getMomentum
TVector3 getMomentum() const
Getter for vector of momentum at closest approach of track in r/phi projection.
Definition: TrackFitResult.h:116
Belle2::VXD::GeoCache::get
static const SensorInfoBase & get(Belle2::VxdID id)
Return a reference to the SensorInfo of a given SensorID.
Definition: GeoCache.h:141
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::VxdID::getLadderNumber
baseType getLadderNumber() const
Get the ladder id.
Definition: VxdID.h:108
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::Environment::isMC
bool isMC() const
Do we have generated, not real data?
Definition: Environment.cc:57
Belle2::SVDCluster::getSNR
float getSNR() const
Get cluster SNR.
Definition: SVDCluster.h:167
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::TrackFitResult::getZ0
double getZ0() const
Getter for z0.
Definition: TrackFitResult.h:200
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::SVDCluster::getPositionSigma
float getPositionSigma() const
Get the error of the reconstructed hit coordinate.
Definition: SVDCluster.h:137
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
Belle2::VXD::GeoCache::getInstance
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:215
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VXD::SensorInfoBase::getThickness
double getThickness() const
Return the thickness of the sensor.
Definition: SensorInfoBase.h:121
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: ParticleList.h:33
Belle2::SVDPerformanceTTreeModule::initialize
void initialize() override
Register input and output data.
Definition: SVDPerformanceTTreeModule.cc:59
Belle2::SVDPerformanceTTreeModule
The module is used to create a TTree to study SVD clusters, genfit unbiased residuals and many other ...
Definition: SVDPerformanceTTreeModule.h:35
Belle2::RecoHitInformation
This class stores additional information to every CDC/SVD/PXD hit stored in a RecoTrack.
Definition: RecoHitInformation.h:48
Belle2::Unit::convertValueToUnit
static double convertValueToUnit(double value, const std::string &unitString)
Converts a floating point value from the standard framework unit to the given unit.
Definition: UnitConst.cc:146
Belle2::VXD::SensorInfoBase::pointToGlobal
TVector3 pointToGlobal(const TVector3 &local, bool reco=false) const
Convert a point from local to global coordinates.
Definition: SensorInfoBase.h:373
Belle2::SVDCluster::isUCluster
bool isUCluster() const
Get the direction of strips.
Definition: SVDCluster.h:118
Belle2::VxdID::getSensorNumber
baseType getSensorNumber() const
Get the sensor id.
Definition: VxdID.h:110
Belle2::SVDCluster::getSensorID
VxdID getSensorID() const
Get the sensor ID.
Definition: SVDCluster.h:110
Belle2::SVDCluster::getPosition
float getPosition(double v=0) const
Get the coordinate of reconstructed hit.
Definition: SVDCluster.h:125
Belle2::SVDCluster::getSize
unsigned short getSize() const
Get cluster size.
Definition: SVDCluster.h:162
Belle2::SVDCluster
The SVD Cluster class This class stores all information about reconstructed SVD clusters.
Definition: SVDCluster.h:38
Belle2::SVDPerformanceTTreeModule::event
void event() override
Compute the variables and fill the tree.
Definition: SVDPerformanceTTreeModule.cc:145
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::Const::ChargedStable
Provides a type-safe way to pass members of the chargedStableSet set.
Definition: Const.h:465
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::VxdID::getLayerNumber
baseType getLayerNumber() const
Get the layer id.
Definition: VxdID.h:106
Belle2::Environment::Instance
static Environment & Instance()
Static method to get a reference to the Environment instance.
Definition: Environment.cc:31
Belle2::SVDPerformanceTTreeModule::terminate
void terminate() override
Write the TTrees to the file.
Definition: SVDPerformanceTTreeModule.cc:358
Belle2::SVDCluster::getClsTime
float getClsTime() const
Get average of waveform maximum times of cluster strip signals.
Definition: SVDCluster.h:142
Belle2::TrackFitResult::getD0
double getD0() const
Getter for d0.
Definition: TrackFitResult.h:178
Belle2::StoreObjPtr::isValid
bool isValid() const
Check whether the object was created.
Definition: StoreObjPtr.h:120