Belle II Software  release-05-01-25
SVDDQMClustersOnTrackModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2017 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giulia Casarosa *
7  * *
8  * Prepared for Belle II geometry *
9  * *
10  * This software is provided "as is" without any warranty. *
11  **************************************************************************/
12 
13 #include "svd/modules/svdDQM/SVDDQMClustersOnTrackModule.h"
14 
15 #include <hlt/softwaretrigger/core/FinalTriggerDecisionCalculator.h>
16 #include <framework/datastore/DataStore.h>
17 #include <framework/datastore/StoreObjPtr.h>
18 #include <framework/datastore/StoreArray.h>
19 #include <framework/dataobjects/EventMetaData.h>
20 #include <svd/dataobjects/SVDShaperDigit.h>
21 #include <svd/dataobjects/SVDRecoDigit.h>
22 #include <svd/dataobjects/SVDCluster.h>
23 #include <tracking/dataobjects/RecoTrack.h>
24 
25 #include <vxd/geometry/GeoTools.h>
26 
27 #include <boost/format.hpp>
28 #include <boost/foreach.hpp>
29 
30 #include "TDirectory.h"
31 
32 using namespace std;
33 using boost::format;
34 using namespace Belle2;
35 using namespace SoftwareTrigger;
36 
37 //-----------------------------------------------------------------
38 // Register the Module
39 //-----------------------------------------------------------------
40 REG_MODULE(SVDDQMClustersOnTrack)
41 
42 
43 //-----------------------------------------------------------------
44 // Implementation
45 //-----------------------------------------------------------------
46 
47 SVDDQMClustersOnTrackModule::SVDDQMClustersOnTrackModule() : HistoModule()
48 {
49  //Set module properties
50  setDescription("SVD DQM module for Clusters related to Tracks.");
51 
52  setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
53  addParam("skipHLTRejectedEvents", m_skipRejectedEvents, "If TRUE skip events rejected by HLT", bool(true));
54  addParam("TriggerBin", m_tb, "select events for a specific trigger bin, if -1 then no seleciton is applied (default)", int(-1));
55  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
56  std::string("SVDClsTrk"));
57  addParam("desynchronizeSVDTime", m_desynchSVDTime,
58  "if TRUE (default is FALSE): svdTime back in SVD time reference", bool(false));
59 
60  m_histoList = new TList();
61 }
62 
63 
64 SVDDQMClustersOnTrackModule::~SVDDQMClustersOnTrackModule()
65 {
66 }
67 
68 //------------------------------------------------------------------
69 // Function to define histograms
70 //-----------------------------------------------------------------
71 
72 void SVDDQMClustersOnTrackModule::defineHisto()
73 {
74  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
75  if (gTools->getNumberOfLayers() == 0) {
76  B2FATAL("Missing geometry for VXD, check steering file.");
77  }
78  if (gTools->getNumberOfSVDLayers() == 0) {
79  B2WARNING("Missing geometry for SVD, SVD-DQM is skiped.");
80  return;
81  }
82 
83  // Create a separate histogram directories and cd into it.
84  TDirectory* oldDir = gDirectory;
85  if (m_histogramDirectoryName != "") {
86  oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
87  oldDir->cd(m_histogramDirectoryName.c_str());
88  }
89 
90  int ChargeBins = 80;
91  float ChargeMax = 160;
92  int SNRBins = 50;
93  float SNRMax = 100;
94  int TimeBins = 300;
95  float TimeMin = -150;
96  float TimeMax = 150;
97 
98  int MaxBinBins = 6;
99  int MaxBinMax = 6;
100 
101  TString refFrame = "in FTSW reference";
102  if (m_desynchSVDTime)
103  refFrame = "in SVD reference";
104 
105 
106  //----------------------------------------------------------------
107  // Charge of clusters for L3/L456 sensors
108  //----------------------------------------------------------------
109  TString name = "SVDTRK_ClusterChargeU3";
110  TString title = "SVD U-Cluster-on-Track Charge for layer 3 sensors";
111  m_clsTrkChargeU3 = new TH1F(name.Data(), title.Data(), ChargeBins, 0, ChargeMax);
112  m_clsTrkChargeU3->GetXaxis()->SetTitle("cluster charge [ke-]");
113  m_clsTrkChargeU3->GetYaxis()->SetTitle("count");
114  m_histoList->Add(m_clsTrkChargeU3);
115  name = "SVDTRK_ClusterChargeV3";
116  title = "SVD V-Cluster-on-Track Charge for layer 3 sensors";
117  m_clsTrkChargeV3 = new TH1F(name.Data(), title.Data(), ChargeBins, 0, ChargeMax);
118  m_clsTrkChargeV3->GetXaxis()->SetTitle("cluster charge [ke-]");
119  m_clsTrkChargeV3->GetYaxis()->SetTitle("count");
120  m_histoList->Add(m_clsTrkChargeV3);
121 
122  name = "SVDTRK_ClusterChargeU456";
123  title = "SVD U-Cluster-on-Track Charge for layers 4,5,6 sensors";
124  m_clsTrkChargeU456 = new TH1F(name.Data(), title.Data(), ChargeBins, 0, ChargeMax);
125  m_clsTrkChargeU456->GetXaxis()->SetTitle("cluster charge [ke-]");
126  m_clsTrkChargeU456->GetYaxis()->SetTitle("count");
127  m_histoList->Add(m_clsTrkChargeU456);
128 
129  name = "SVDTRK_ClusterChargeV456";
130  title = "SVD V-Cluster-on-Track Charge for layers 4,5,6 sensors";
131  m_clsTrkChargeV456 = new TH1F(name.Data(), title.Data(), ChargeBins, 0, ChargeMax);
132  m_clsTrkChargeV456->GetXaxis()->SetTitle("cluster charge [ke-]");
133  m_clsTrkChargeV456->GetYaxis()->SetTitle("count");
134  m_histoList->Add(m_clsTrkChargeV456);
135 
136  //----------------------------------------------------------------
137  // SNR of clusters for L3/L456 sensors
138  //----------------------------------------------------------------
139  name = "SVDTRK_ClusterSNRU3";
140  title = "SVD U-Cluster-on-Track SNR for layer 3 sensors";
141  m_clsTrkSNRU3 = new TH1F(name.Data(), title.Data(), SNRBins, 0, SNRMax);
142  m_clsTrkSNRU3->GetXaxis()->SetTitle("cluster SNR");
143  m_clsTrkSNRU3->GetYaxis()->SetTitle("count");
144  m_histoList->Add(m_clsTrkSNRU3);
145  name = "SVDTRK_ClusterSNRV3";
146  title = "SVD V-Cluster-on-Track SNR for layer 3 sensors";
147  m_clsTrkSNRV3 = new TH1F(name.Data(), title.Data(), SNRBins, 0, SNRMax);
148  m_clsTrkSNRV3->GetXaxis()->SetTitle("cluster SNR");
149  m_clsTrkSNRV3->GetYaxis()->SetTitle("count");
150  m_histoList->Add(m_clsTrkSNRV3);
151 
152  name = "SVDTRK_ClusterSNRU456";
153  title = "SVD U-Cluster-on-Track SNR for layers 4,5,6 sensors";
154  m_clsTrkSNRU456 = new TH1F(name.Data(), title.Data(), SNRBins, 0, SNRMax);
155  m_clsTrkSNRU456->GetXaxis()->SetTitle("cluster SNR");
156  m_clsTrkSNRU456->GetYaxis()->SetTitle("count");
157  m_histoList->Add(m_clsTrkSNRU456);
158  name = "SVDTRK_ClusterSNRV456";
159  title = "SVD V-Cluster-on-Track SNR for layers 4,5,6 sensors";
160  m_clsTrkSNRV456 = new TH1F(name.Data(), title.Data(), SNRBins, 0, SNRMax);
161  m_clsTrkSNRV456->GetXaxis()->SetTitle("cluster SNR");
162  m_clsTrkSNRV456->GetYaxis()->SetTitle("count");
163  m_histoList->Add(m_clsTrkSNRV456);
164 
165  //----------------------------------------------------------------
166  // Time of clusters for L3/L456 sensors
167  //----------------------------------------------------------------
168  name = "SVDTRK_ClusterTimeU3";
169  title = Form("SVD U-Cluster-on-Track Time %s for layer 3 sensors", refFrame.Data());
170  m_clsTrkTimeU3 = new TH1F(name.Data(), title.Data(), TimeBins, TimeMin, TimeMax);
171  m_clsTrkTimeU3->GetXaxis()->SetTitle("cluster time (ns)");
172  m_clsTrkTimeU3->GetYaxis()->SetTitle("count");
173  m_histoList->Add(m_clsTrkTimeU3);
174  name = "SVDTRK_ClusterTimeV3";
175  title = Form("SVD V-Cluster-on-Track Time %s for layer 3 sensors", refFrame.Data());
176  m_clsTrkTimeV3 = new TH1F(name.Data(), title.Data(), TimeBins, TimeMin, TimeMax);
177  m_clsTrkTimeV3->GetXaxis()->SetTitle("cluster time (ns)");
178  m_clsTrkTimeV3->GetYaxis()->SetTitle("count");
179  m_histoList->Add(m_clsTrkTimeV3);
180 
181  name = "SVDTRK_ClusterTimeU456";
182  title = Form("SVD U-Cluster-on-Track Time %s for layers 4,5,6 sensors", refFrame.Data());
183  m_clsTrkTimeU456 = new TH1F(name.Data(), title.Data(), TimeBins, TimeMin, TimeMax);
184  m_clsTrkTimeU456->GetXaxis()->SetTitle("cluster time (ns)");
185  m_clsTrkTimeU456->GetYaxis()->SetTitle("count");
186  m_histoList->Add(m_clsTrkTimeU456);
187  name = "SVDTRK_ClusterTimeV456";
188  title = Form("SVD V-Cluster-on-Track Time %s for layers 4,5,6 sensors", refFrame.Data());
189  m_clsTrkTimeV456 = new TH1F(name.Data(), title.Data(), TimeBins, TimeMin, TimeMax);
190  m_clsTrkTimeV456->GetXaxis()->SetTitle("cluster time (ns)");
191  m_clsTrkTimeV456->GetYaxis()->SetTitle("count");
192  m_histoList->Add(m_clsTrkTimeV456);
193 
194  //----------------------------------------------------------------
195  // EventT0 vs Time of clusters for U and V sides
196  //----------------------------------------------------------------
197  name = "SVDTRK_ClusterTimeUvsEventT0";
198  title = Form("SVD U-Cluster-on-Track Time vs EventT0 %s for layer 3 sensors", refFrame.Data());
199  m_clsTrkTimeUEvtT0 = new TH2F(name.Data(), title.Data(), TimeBins, TimeMin, TimeMax, 100, -50, 50);
200  m_clsTrkTimeUEvtT0->GetXaxis()->SetTitle("clusters time (ns)");
201  m_clsTrkTimeUEvtT0->GetYaxis()->SetTitle("EventT0 (ns)");
202  m_histoList->Add(m_clsTrkTimeUEvtT0);
203  name = "SVDTRK_ClusterTimeVvsEventT0";
204  title = Form("SVD V-Cluster-on-Track Time vs EventT0 %s for layer 3 sensors", refFrame.Data());
205  m_clsTrkTimeVEvtT0 = new TH2F(name.Data(), title.Data(), TimeBins, TimeMin, TimeMax, 100, -50, 50);
206  m_clsTrkTimeVEvtT0->GetXaxis()->SetTitle("cluster time (ns)");
207  m_clsTrkTimeVEvtT0->GetYaxis()->SetTitle("EventT0 (ns)");
208  m_histoList->Add(m_clsTrkTimeVEvtT0);
209 
210  //----------------------------------------------------------------
211  // MaxBin of strips for all sensors (offline ZS)
212  //----------------------------------------------------------------
213  name = "SVDTRK_StripMaxBinUAll";
214  title = "SVD U-Strip-on-Track MaxBin for all sensors";
215  m_stripMaxBinUAll = new TH1F(name.Data(), title.Data(), MaxBinBins, 0, MaxBinMax);
216  m_stripMaxBinUAll->GetXaxis()->SetTitle("max bin");
217  m_stripMaxBinUAll->GetYaxis()->SetTitle("count");
218  m_histoList->Add(m_stripMaxBinUAll);
219  name = "SVDTRK_StripMaxBinVAll";
220  title = "SVD V-Strip-on-Track MaxBin for all sensors";
221  m_stripMaxBinVAll = new TH1F(name.Data(), title.Data(), MaxBinBins, 0, MaxBinMax);
222  m_stripMaxBinVAll->GetXaxis()->SetTitle("max bin");
223  m_stripMaxBinVAll->GetYaxis()->SetTitle("count");
224  m_histoList->Add(m_stripMaxBinVAll);
225 
226 
227  oldDir->cd();
228 }
229 
230 
231 void SVDDQMClustersOnTrackModule::initialize()
232 {
233  // Register histograms (calls back defineHisto)
234  REG_HISTOGRAM
235 
236  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
237  if (gTools->getNumberOfSVDLayers() != 0) {
238 
239  m_svdEventInfo.isOptional();
240  m_eventT0.isOptional();
241  m_storeTracks.isOptional();
242  m_resultStoreObjectPointer.isOptional();
243 
244  }
245 }
246 
247 void SVDDQMClustersOnTrackModule::beginRun()
248 {
249  StoreObjPtr<EventMetaData> evtMetaData;
250  m_expNumber = evtMetaData->getExperiment();
251  m_runNumber = evtMetaData->getRun();
252 
253  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
254  if (gTools->getNumberOfSVDLayers() == 0) return;
255 
256  // Add experiment and run number to the title of selected histograms (CR shifter plots)
257  TString runID = TString::Format(" ~ Exp%d Run%d", m_expNumber, m_runNumber);
258  TObject* obj;
259  TIter nextH(m_histoList);
260  while ((obj = nextH()))
261  if (obj->InheritsFrom("TH1")) {
262  ((TH1F*)obj)->SetTitle(obj->GetTitle() + runID);
263  if (obj != NULL)((TH1F*)obj)->Reset();
264  }
265 
266 }
267 
268 void SVDDQMClustersOnTrackModule::event()
269 {
270 
271  if (!m_storeTracks.isValid()) {
272  B2WARNING("Missing Tracks StoreArray. Skipping SVDDQMClustersOnTrack");
273  return;
274  }
275 
276  if (!m_svdEventInfo.isValid())
277  m_tb = -1;
278  else {
279  if (m_tb != -1)
280  if (m_svdEventInfo->getModeByte().getTriggerBin() != m_tb)
281  return;
282  }
283 
284  // get EventT0 if present and valid
285  double eventT0 = -1000;
286  if (m_eventT0.isOptional())
287  if (m_eventT0.isValid())
288  if (m_eventT0->hasEventT0())
289  eventT0 = m_eventT0->getEventT0();
290 
291  // if svd time in SVD time reference is shown, eventT0 is also synchronized with SVD reference frame, firstFrame = 0
292  if (m_desynchSVDTime && m_svdEventInfo.isValid())
293  eventT0 = eventT0 - m_svdEventInfo->getSVD2FTSWTimeShift(0);
294 
295  //check HLT decision and increase number of events only if the event has been accepted
296 
297  if (m_skipRejectedEvents && (m_resultStoreObjectPointer.isValid())) {
298  const bool eventAccepted = FinalTriggerDecisionCalculator::getFinalTriggerDecision(*m_resultStoreObjectPointer);
299  if (!eventAccepted) return;
300  }
301 
302  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
303  if (gTools->getNumberOfSVDLayers() == 0) return;
304 
305 
306  BOOST_FOREACH(Track & track, m_storeTracks) {
307 
308  const TrackFitResult* tfr = track.getTrackFitResultWithClosestMass(Const::pion);
309  if (!tfr) continue;
310 
311  RelationVector<RecoTrack> theRC = DataStore::getRelationsWithObj<RecoTrack>(&track);
312  RelationVector<SVDCluster> svdClustersTrack = DataStore::getRelationsWithObj<SVDCluster>(theRC[0]);
313 
314  for (int cl = 0 ; cl < (int)svdClustersTrack.size(); cl++) {
315 
316  int iLayer = svdClustersTrack[cl]->getSensorID().getLayerNumber();
317 
318  float time = svdClustersTrack[cl]->getClsTime();
319  if (m_desynchSVDTime && m_svdEventInfo.isValid())
320  time = time - m_svdEventInfo->getSVD2FTSWTimeShift(svdClustersTrack[cl]->getFirstFrame());
321 
322  if (svdClustersTrack[cl]->isUCluster()) {
323 
324  m_clsTrkTimeUEvtT0->Fill(time, eventT0);
325 
326  if (iLayer == 3) {
327  if (m_clsTrkChargeU3 != NULL) m_clsTrkChargeU3->Fill(svdClustersTrack[cl]->getCharge() / 1000.0); // in kelectrons
328  if (m_clsTrkSNRU3 != NULL) m_clsTrkSNRU3->Fill(svdClustersTrack[cl]->getSNR());
329  if (m_clsTrkTimeU3 != NULL) m_clsTrkTimeU3->Fill(time);
330  } else {
331  if (m_clsTrkChargeU456 != NULL) m_clsTrkChargeU456->Fill(svdClustersTrack[cl]->getCharge() / 1000.0); // in kelectrons
332  if (m_clsTrkSNRU456 != NULL) m_clsTrkSNRU456->Fill(svdClustersTrack[cl]->getSNR());
333  if (m_clsTrkTimeU456 != NULL) m_clsTrkTimeU456->Fill(time);
334  }
335 
336  RelationVector<SVDRecoDigit> recoDigits = svdClustersTrack[cl]->getRelationsTo<SVDRecoDigit>();
337  if (recoDigits.size() == 0) continue;
338 
339  for (UInt_t re = 0; re < recoDigits.size(); re++) {
340 
341  RelationVector<SVDShaperDigit> shaper = recoDigits[re]->getRelationsTo<SVDShaperDigit>();
342  if (m_stripMaxBinUAll != NULL) m_stripMaxBinUAll->Fill(shaper[0]->getMaxTimeBin());
343  }
344 
345 
346  } else {
347 
348  m_clsTrkTimeVEvtT0->Fill(time, eventT0);
349 
350  if (iLayer == 3) {
351  if (m_clsTrkChargeV3 != NULL) m_clsTrkChargeV3->Fill(svdClustersTrack[cl]->getCharge() / 1000.0); // in kelectrons
352  if (m_clsTrkSNRV3 != NULL) m_clsTrkSNRV3->Fill(svdClustersTrack[cl]->getSNR());
353  if (m_clsTrkTimeV3 != NULL) m_clsTrkTimeV3->Fill(time);
354  } else {
355  if (m_clsTrkChargeV456 != NULL) m_clsTrkChargeV456->Fill(svdClustersTrack[cl]->getCharge() / 1000.0); // in kelectrons
356  if (m_clsTrkSNRV456 != NULL) m_clsTrkSNRV456->Fill(svdClustersTrack[cl]->getSNR());
357  if (m_clsTrkTimeV456 != NULL) m_clsTrkTimeV456->Fill(time);
358  }
359 
360  RelationVector<SVDRecoDigit> recoDigits = svdClustersTrack[cl]->getRelationsTo<SVDRecoDigit>();
361  if (recoDigits.size() == 0) continue;
362 
363  for (UInt_t re = 0; re < recoDigits.size(); re++) {
364 
365  RelationVector<SVDShaperDigit> shaper = recoDigits[re]->getRelationsTo<SVDShaperDigit>();
366  if (m_stripMaxBinVAll != NULL) m_stripMaxBinVAll->Fill(shaper[0]->getMaxTimeBin());
367 
368  }
369 
370  }
371 
372  }
373  }
374 }
375 
376 
377 void SVDDQMClustersOnTrackModule::terminate()
378 {
379 
380  delete m_histoList;
381 
382 }
Belle2::RelationVector::size
size_t size() const
Get number of relations.
Definition: RelationVector.h:98
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::SVDShaperDigit
The SVD ShaperDigit class.
Definition: SVDShaperDigit.h:46
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: ParticleList.h:33
Belle2::SVDDQMClustersOnTrackModule
SVD DQM Module for Clusters related to Tracks.
Definition: SVDDQMClustersOnTrackModule.h:46
Belle2::Track
Class that bundles various TrackFitResults.
Definition: Track.h:35
Belle2::SVDRecoDigit
The SVD RecoDigit class.
Definition: SVDRecoDigit.h:54
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29