Belle II Software development
SVDDQMExpressRecoModule.cc
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#include "svd/modules/svdDQM/SVDDQMExpressRecoModule.h"
10
11#include <hlt/softwaretrigger/core/FinalTriggerDecisionCalculator.h>
12#include <framework/datastore/StoreObjPtr.h>
13#include <framework/datastore/StoreArray.h>
14#include <framework/dataobjects/EventMetaData.h>
15
16#include <svd/dataobjects/SVDShaperDigit.h>
17#include <svd/dataobjects/SVDRecoDigit.h>
18#include <svd/dataobjects/SVDCluster.h>
19
20#include <vxd/geometry/SensorInfoBase.h>
21#include <vxd/geometry/GeoTools.h>
22
23#include <boost/format.hpp>
24
25#include "TDirectory.h"
26
27using namespace std;
28using boost::format;
29using namespace Belle2;
30using namespace SoftwareTrigger;
31
32//-----------------------------------------------------------------
33// Register the Module
34//-----------------------------------------------------------------
35REG_MODULE(SVDDQMExpressReco);
36
37
38//-----------------------------------------------------------------
39// Implementation
40//-----------------------------------------------------------------
41
43{
44 //Set module properties
45 setDescription("Original SVD DQM module for ExpressReco.");
46
47 setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
48 addParam("offlineZSShaperDigits", m_storeSVDShaperDigitsName, "ShaperDigits StoreArray name - usually ZS5 strips.",
49 std::string("SVDShaperDigitsZS5"));
50 addParam("ShaperDigits", m_storeNoZSSVDShaperDigitsName, "not zero-suppressed ShaperDigits StoreArray name.",
51 std::string("SVDShaperDigits"));
52 addParam("Clusters", m_storeSVDClustersName, "Cluster StoreArray name.",
53 std::string("SVDClusters"));
54 addParam("ShowAllHistos", m_ShowAllHistos, "Flag to show all histos in DQM, default = 0.", int(0));
55 addParam("desynchronizeSVDTime", m_desynchSVDTime,
56 "if True, svd time back in SVD time reference.", bool(false));
57 addParam("CutSVDCharge", m_CutSVDCharge,
58 "minimum charge (ADC) to fill the strip-hitmap histogram.", float(0));
59 addParam("CutSVDClusterCharge", m_CutSVDClusterCharge,
60 "minimum charge (in e-) to fill the cluster-hitmap histogram.", float(0));
61 addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed.",
62 std::string("SVDExpReco"));
63 addParam("additionalPlots", m_additionalPlots, "Flag to produce additional plots",
64 bool(false));
65 addParam("useParamFromDB", m_useParamFromDB, "use SVDDQMPlotsConfiguration from DB", bool(true));
66 addParam("skipHLTRejectedEvents", m_skipRejectedEvents, "If True, skip events rejected by HLT.", bool(false));
67 addParam("samples3", m_3Samples, "if True 3 samples histograms analysis is performed", bool(false));
68
69 m_histoList = new TList();
70}
71
72
73SVDDQMExpressRecoModule::~SVDDQMExpressRecoModule()
74{
75}
76
77//------------------------------------------------------------------
78// Function to define histograms
79//-----------------------------------------------------------------
80
82{
83 if (m_useParamFromDB) {
84 if (!m_svdPlotsConfig.isValid())
85 B2FATAL("no valid configuration found for SVD reconstruction");
86 else {
87 B2DEBUG(20, "SVDRecoConfiguration: from now on we are using " << m_svdPlotsConfig->get_uniqueID());
88 m_3Samples = m_svdPlotsConfig->isPlotsFor3SampleMonitoring();
89 m_skipRejectedEvents = m_svdPlotsConfig->isSkipHLTRejectedEvents();
90 }
91 }
92
94 if (gTools->getNumberOfLayers() == 0) {
95 B2FATAL("Missing geometry for VXD, check steering file.");
96 }
97 if (gTools->getNumberOfSVDLayers() == 0) {
98 B2WARNING("Missing geometry for SVD, SVD-DQM is skipped.");
99 return;
100 }
101
102 // Create a separate histogram directories and cd into it.
103 TDirectory* oldDir = gDirectory;
104 if (m_histogramDirectoryName != "") {
105 oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
106 oldDir->cd(m_histogramDirectoryName.c_str());
107 }
108
109 // basic constants presets:
110 int nSVDSensors = gTools->getNumberOfSVDSensors();
111 int nSVDChips = gTools->getTotalSVDChips();
112
113 // number of events counter
114 m_nEvents = new TH1F("SVDDQM_nEvents", "SVD Number of Events", 1, -0.5, 0.5);
115 m_nEvents->GetYaxis()->SetTitle("N events");
117
118 // Create basic histograms:
119 // basic counters per sensor:
120 m_hitMapCountsU = new TH1F("SVDDQM_StripCountsU", "SVD Integrated Number of ZS5 Fired U-Strips per sensor",
121 nSVDSensors, 0, nSVDSensors);
122 m_hitMapCountsU->GetXaxis()->SetTitle("Sensor ID");
123 m_hitMapCountsU->GetYaxis()->SetTitle("counts");
125 m_hitMapCountsV = new TH1F("SVDDQM_StripCountsV", "SVD Integrated Number of ZS5 Fired V-Strips per sensor",
126 nSVDSensors, 0, nSVDSensors);
127 m_hitMapCountsV->GetXaxis()->SetTitle("Sensor ID");
128 m_hitMapCountsV->GetYaxis()->SetTitle("counts");
130 m_hitMapClCountsU = new TH1F("SVDDQM_ClusterCountsU", "SVD Integrated Number of U-Clusters per sensor",
131 nSVDSensors, 0, nSVDSensors);
132 m_hitMapClCountsU->GetXaxis()->SetTitle("Sensor ID");
133 m_hitMapClCountsU->GetYaxis()->SetTitle("counts");
135 m_hitMapClCountsV = new TH1F("SVDDQM_ClusterCountsV", "SVD Integrated Number of V-Clusters per sensor",
136 nSVDSensors, 0, nSVDSensors);
137 m_hitMapClCountsV->GetXaxis()->SetTitle("Sensor ID");
138 m_hitMapClCountsV->GetYaxis()->SetTitle("counts");
140 for (int i = 0; i < nSVDSensors; i++) {
141 VxdID id = gTools->getSensorIDFromSVDIndex(i);
142 int iLayer = id.getLayerNumber();
143 int iLadder = id.getLadderNumber();
144 int iSensor = id.getSensorNumber();
145 TString AxisTicks = Form("%i_%i_%i", iLayer, iLadder, iSensor);
146 m_hitMapCountsU->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
147 m_hitMapCountsV->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
148 m_hitMapClCountsU->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
149 m_hitMapClCountsV->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
150 }
151
152 // basic counters per chip:
153 m_hitMapCountsChip = new TH1F("SVDDQM_StripCountsChip", "SVD Integrated Number of ZS5 Fired Strips per chip",
154 nSVDChips, 0, nSVDChips);
155 m_hitMapCountsChip->GetXaxis()->SetTitle("Chip ID");
156 m_hitMapCountsChip->GetYaxis()->SetTitle("counts");
158 m_hitMapClCountsChip = new TH1F("SVDDQM_ClusterCountsChip", "SVD Integrated Number of Clusters per chip",
159 nSVDChips, 0, nSVDChips);
160 m_hitMapClCountsChip->GetXaxis()->SetTitle("Chip ID");
161 m_hitMapClCountsChip->GetYaxis()->SetTitle("counts");
163
164 if (m_additionalPlots) {
165 m_firedU = new TH1F*[nSVDSensors];
166 m_firedV = new TH1F*[nSVDSensors];
167 m_clustersU = new TH1F*[nSVDSensors];
168 m_clustersV = new TH1F*[nSVDSensors];
169 m_stripSignalU = new TH1F*[nSVDSensors];
170 m_stripSignalV = new TH1F*[nSVDSensors];
171 }
172
173 m_clusterChargeU = new TH1F*[nSVDSensors];
174 m_clusterChargeV = new TH1F*[nSVDSensors];
175 m_clusterSNRU = new TH1F*[nSVDSensors];
176 m_clusterSNRV = new TH1F*[nSVDSensors];
177
178 m_stripCountU = new TH1F*[nSVDSensors];
179 m_stripCountV = new TH1F*[nSVDSensors];
180 m_strip3SampleCountU = new TH1F*[nSVDSensors];
181 m_strip3SampleCountV = new TH1F*[nSVDSensors];
182 m_strip6SampleCountU = new TH1F*[nSVDSensors];
183 m_strip6SampleCountV = new TH1F*[nSVDSensors];
184
185 m_stripCountGroupId0U = new TH1F*[nSVDSensors];
186 m_stripCountGroupId0V = new TH1F*[nSVDSensors];
187
188 m_onlineZSstripCountU = new TH1F*[nSVDSensors];
189 m_onlineZSstripCountV = new TH1F*[nSVDSensors];
190
191 if (m_3Samples) {
192 m_onlineZSstrip3SampleCountU = new TH1F*[nSVDSensors];
193 m_onlineZSstrip3SampleCountV = new TH1F*[nSVDSensors];
194 m_onlineZSstrip6sampleCountU = new TH1F*[nSVDSensors];
195 m_onlineZSstrip6sampleCountV = new TH1F*[nSVDSensors];
196 }
197
198 m_clusterSizeU = new TH1F*[nSVDSensors];
199 m_clusterSizeV = new TH1F*[nSVDSensors];
200 m_clusterTimeU = new TH1F*[nSVDSensors];
201 m_clusterTimeV = new TH1F*[nSVDSensors];
202
203 int ChargeBins = 80;
204 float ChargeMax = 80;
205 int SNRBins = 50;
206 float SNRMax = 100;
207 int TimeBins = 300;
208 float TimeMin = -150;
209 float TimeMax = 150;
210
211 int GroupIdBins = 21;
212 float GroupIdMin = -1.5;
213 float GroupIdMax = 19.5;
214
215 int MaxBinBins = 6;
216 int MaxBinMax = 6;
217
218 TString refFrame = "in FTSW reference";
220 refFrame = "in SVD reference";
221
222
223 //----------------------------------------------------------------
224 // Charge of clusters for all sensors
225 //----------------------------------------------------------------
226 string name = str(format("SVDDQM_ClusterChargeUAll"));
227 string title = str(format("SVD U-Cluster Charge for all sensors"));
228 m_clusterChargeUAll = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
229 m_clusterChargeUAll->GetXaxis()->SetTitle("cluster charge [ke-]");
230 m_clusterChargeUAll->GetYaxis()->SetTitle("count");
232 name = str(format("SVDDQM_ClusterChargeVAll"));
233 title = str(format("SVD V-Cluster Charge for all sensors"));
234 m_clusterChargeVAll = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
235 m_clusterChargeVAll->GetXaxis()->SetTitle("cluster charge [ke-]");
236 m_clusterChargeVAll->GetYaxis()->SetTitle("count");
238 //----------------------------------------------------------------
239 // Charge of clusters for L3/L456 sensors
240 //----------------------------------------------------------------
241 name = str(format("SVDDQM_ClusterChargeU3"));
242 title = str(format("SVD U-Cluster Charge for layer 3 sensors"));
243 m_clusterChargeU3 = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
244 m_clusterChargeU3->GetXaxis()->SetTitle("cluster charge [ke-]");
245 m_clusterChargeU3->GetYaxis()->SetTitle("count");
247 name = str(format("SVDDQM_ClusterChargeV3"));
248 title = str(format("SVD V-Cluster Charge for layer 3 sensors"));
249 m_clusterChargeV3 = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
250 m_clusterChargeV3->GetXaxis()->SetTitle("cluster charge [ke-]");
251 m_clusterChargeV3->GetYaxis()->SetTitle("count");
253
254 name = str(format("SVDDQM_ClusterChargeU456"));
255 title = str(format("SVD U-Cluster Charge for layers 4,5,6 sensors"));
256 m_clusterChargeU456 = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
257 m_clusterChargeU456->GetXaxis()->SetTitle("cluster charge [ke-]");
258 m_clusterChargeU456->GetYaxis()->SetTitle("count");
260
261 name = str(format("SVDDQM_ClusterChargeV456"));
262 title = str(format("SVD V-Cluster Charge for layers 4,5,6 sensors"));
263 m_clusterChargeV456 = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
264 m_clusterChargeV456->GetXaxis()->SetTitle("cluster charge [ke-]");
265 m_clusterChargeV456->GetYaxis()->SetTitle("count");
267
268 //----------------------------------------------------------------
269 // SNR of clusters for all sensors
270 //----------------------------------------------------------------
271 name = str(format("SVDDQM_ClusterSNRUAll"));
272 title = str(format("SVD U-Cluster SNR for all sensors"));
273 m_clusterSNRUAll = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax); // max = ~ 60
274 m_clusterSNRUAll->GetXaxis()->SetTitle("cluster SNR");
275 m_clusterSNRUAll->GetYaxis()->SetTitle("count");
277 name = str(format("SVDDQM_ClusterSNRVAll"));
278 title = str(format("SVD V-Cluster SNR for all sensors"));
279 m_clusterSNRVAll = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
280 m_clusterSNRVAll->GetXaxis()->SetTitle("cluster SNR");
281 m_clusterSNRVAll->GetYaxis()->SetTitle("count");
283 //----------------------------------------------------------------
284 // SNR of clusters for L3/L456 sensors
285 //----------------------------------------------------------------
286 name = str(format("SVDDQM_ClusterSNRU3"));
287 title = str(format("SVD U-Cluster SNR for layer 3 sensors"));
288 m_clusterSNRU3 = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
289 m_clusterSNRU3->GetXaxis()->SetTitle("cluster SNR");
290 m_clusterSNRU3->GetYaxis()->SetTitle("count");
292 name = str(format("SVDDQM_ClusterSNRV3"));
293 title = str(format("SVD V-Cluster SNR for layer 3 sensors"));
294 m_clusterSNRV3 = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
295 m_clusterSNRV3->GetXaxis()->SetTitle("cluster SNR");
296 m_clusterSNRV3->GetYaxis()->SetTitle("count");
298
299 name = str(format("SVDDQM_ClusterSNRU456"));
300 title = str(format("SVD U-Cluster SNR for layers 4,5,6 sensors"));
301 m_clusterSNRU456 = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
302 m_clusterSNRU456->GetXaxis()->SetTitle("cluster SNR");
303 m_clusterSNRU456->GetYaxis()->SetTitle("count");
305 name = str(format("SVDDQM_ClusterSNRV456"));
306 title = str(format("SVD V-Cluster SNR for layers 4,5,6 sensors"));
307 m_clusterSNRV456 = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
308 m_clusterSNRV456->GetXaxis()->SetTitle("cluster SNR");
309 m_clusterSNRV456->GetYaxis()->SetTitle("count");
311 //----------------------------------------------------------------
312 // Cluster time distribution for all sensors
313 //----------------------------------------------------------------
314 TString Name = "SVDDQM_ClusterTimeUAll";
315 TString Title = Form("SVD U-Cluster Time %s for all sensors", refFrame.Data());
316 m_clusterTimeUAll = new TH1F(Name.Data(), Title.Data(), TimeBins, TimeMin, TimeMax);
317 m_clusterTimeUAll->GetXaxis()->SetTitle("cluster time (ns)");
318 m_clusterTimeUAll->GetYaxis()->SetTitle("count");
320 Name = "SVDDQM_ClusterTimeVAll";
321 Title = Form("SVD V-Cluster Time %s for all sensors", refFrame.Data());
322 m_clusterTimeVAll = new TH1F(Name.Data(), Title.Data(), TimeBins, TimeMin, TimeMax);
323 m_clusterTimeVAll->GetXaxis()->SetTitle("cluster time (ns)");
324 m_clusterTimeVAll->GetYaxis()->SetTitle("count");
326 //----------------------------------------------------------------
327 // Time of clusters for L3/L456 sensors
328 //----------------------------------------------------------------
329 Name = "SVDDQM_ClusterTimeU3";
330 Title = Form("SVD U-Cluster Time %s for layer 3 sensors", refFrame.Data());
331 m_clusterTimeU3 = new TH1F(Name.Data(), Title.Data(), TimeBins, TimeMin, TimeMax);
332 m_clusterTimeU3->GetXaxis()->SetTitle("cluster time (ns)");
333 m_clusterTimeU3->GetYaxis()->SetTitle("count");
335 name = str(format("SVDDQM_ClusterTimeV3"));
336 Title = Form("SVD V-Cluster Time %s for layer 3 sensors", refFrame.Data());
337 m_clusterTimeV3 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
338 m_clusterTimeV3->GetXaxis()->SetTitle("cluster time (ns)");
339 m_clusterTimeV3->GetYaxis()->SetTitle("count");
341
342 name = str(format("SVDDQM_ClusterTimeU456"));
343 Title = Form("SVD U-Cluster Time %s for layers 4,5,6 sensors", refFrame.Data());
344 m_clusterTimeU456 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
345 m_clusterTimeU456->GetXaxis()->SetTitle("cluster time (ns)");
346 m_clusterTimeU456->GetYaxis()->SetTitle("count");
348 name = str(format("SVDDQM_ClusterTimeV456"));
349 Title = Form("SVD V-Cluster Time %s for layers 4,5,6 sensors", refFrame.Data());
350 m_clusterTimeV456 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
351 m_clusterTimeV456->GetXaxis()->SetTitle("cluster time (ns)");
352 m_clusterTimeV456->GetYaxis()->SetTitle("count");
354
355 //----------------------------------------------------------------
356 // Time of clusters for L3/L456 sensors for 3 samples
357 //----------------------------------------------------------------
358 if (m_3Samples) {
359 Name = "SVDDQM_Cluster3TimeU3";
360 Title = Form("SVD U-Cluster Time %s for layer 3 sensors for 3 samples", refFrame.Data());
361 m_cluster3SampleTimeU3 = new TH1F(Name.Data(), Title.Data(), TimeBins, TimeMin, TimeMax);
362 m_cluster3SampleTimeU3->GetXaxis()->SetTitle("cluster time (ns)");
363 m_cluster3SampleTimeU3->GetYaxis()->SetTitle("count");
365 name = str(format("SVDDQM_Cluster3TimeV3"));
366 Title = Form("SVD V-Cluster Time %s for layer 3 sensors for 3 samples", refFrame.Data());
367 m_cluster3SampleTimeV3 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
368 m_cluster3SampleTimeV3->GetXaxis()->SetTitle("cluster time (ns)");
369 m_cluster3SampleTimeV3->GetYaxis()->SetTitle("count");
371 name = str(format("SVDDQM_Cluster3TimeU456"));
372 Title = Form("SVD U-Cluster Time %s for layers 4,5,6 sensors for 3 samples", refFrame.Data());
373 m_cluster3SampleTimeU456 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
374 m_cluster3SampleTimeU456->GetXaxis()->SetTitle("cluster time (ns)");
375 m_cluster3SampleTimeU456->GetYaxis()->SetTitle("count");
377 name = str(format("SVDDQM_Cluster3TimeV456"));
378 Title = Form("SVD V-Cluster Time %s for layers 4,5,6 sensors for 3 samples", refFrame.Data());
379 m_cluster3SampleTimeV456 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
380 m_cluster3SampleTimeV456->GetXaxis()->SetTitle("cluster time (ns)");
381 m_cluster3SampleTimeV456->GetYaxis()->SetTitle("count");
383
384 //----------------------------------------------------------------
385 // Time of clusters for L3/L456 sensors for 6 samples
386 //----------------------------------------------------------------
387 Name = "SVDDQM_Cluster6TimeU3";
388 Title = Form("SVD U-Cluster Time %s for layer 3 sensors for 6 samples", refFrame.Data());
389 m_cluster6SampleTimeU3 = new TH1F(Name.Data(), Title.Data(), TimeBins, TimeMin, TimeMax);
390 m_cluster6SampleTimeU3->GetXaxis()->SetTitle("cluster time (ns)");
391 m_cluster6SampleTimeU3->GetYaxis()->SetTitle("count");
393 name = str(format("SVDDQM_Cluster6TimeV3"));
394 Title = Form("SVD V-Cluster Time %s for layer 3 sensors for 6 samples", refFrame.Data());
395 m_cluster6SampleTimeV3 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
396 m_cluster6SampleTimeV3->GetXaxis()->SetTitle("cluster time (ns)");
397 m_cluster6SampleTimeV3->GetYaxis()->SetTitle("count");
399
400 name = str(format("SVDDQM_Cluster6TimeU456"));
401 Title = Form("SVD U-Cluster Time %s for layers 4,5,6 sensors for 6 samples", refFrame.Data());
402 m_cluster6SampleTimeU456 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
403 m_cluster6SampleTimeU456->GetXaxis()->SetTitle("cluster time (ns)");
404 m_cluster6SampleTimeU456->GetYaxis()->SetTitle("count");
406 name = str(format("SVDDQM_Cluster6TimeV456"));
407 Title = Form("SVD V-Cluster Time %s for layers 4,5,6 sensors for 6 samples", refFrame.Data());
408 m_cluster6SampleTimeV456 = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
409 m_cluster6SampleTimeV456->GetXaxis()->SetTitle("cluster time (ns)");
410 m_cluster6SampleTimeV456->GetYaxis()->SetTitle("count");
412 }
413
414 //----------------------------------------------------------------
415 // Cluster time group Id vs cluster time for U/V sensors
416 //----------------------------------------------------------------
417 Name = "SVDDQM_ClusterTimeGroupIdU";
418 Title = Form("SVD cluster Time Group Id %s vs cluster time for U/P Side", refFrame.Data());
419 m_clusterTimeGroupIdU = new TH2F(Name.Data(), Title.Data(), TimeBins / 2, TimeMin, TimeMax, GroupIdBins, GroupIdMin, GroupIdMax);
420 m_clusterTimeGroupIdU->GetXaxis()->SetTitle("cluster time (ns)");
421 m_clusterTimeGroupIdU->GetYaxis()->SetTitle("cluster group id");
423 Name = "SVDDQM_ClusterTimeGroupIdV";
424 Title = Form("SVD cluster Time Group Id %s vs cluster time for V/N Side", refFrame.Data());
425 m_clusterTimeGroupIdV = new TH2F(Name.Data(), Title.Data(), TimeBins / 2, TimeMin, TimeMax, GroupIdBins, GroupIdMin, GroupIdMax);
426 m_clusterTimeGroupIdV->GetXaxis()->SetTitle("cluster time (ns)");
427 m_clusterTimeGroupIdV->GetYaxis()->SetTitle("cluster group id");
429
430 //----------------------------------------------------------------
431 // Cluster time group Id vs cluster time for U/V sensors for coarse and fine trigger
432 //----------------------------------------------------------------
433 Name = "SVDDQM_cluster6TimeGroupIdU";
434 Title = Form("SVD cluster Time Group Id %s vs cluster time for U/P Side for coarse trigger", refFrame.Data());
435 m_clusterTimeCoarseGroupIdU = new TH2F(Name.Data(), Title.Data(), TimeBins / 2, TimeMin, TimeMax, GroupIdBins, GroupIdMin,
436 GroupIdMax);
437 m_clusterTimeCoarseGroupIdU->GetXaxis()->SetTitle("cluster time (ns)");
438 m_clusterTimeCoarseGroupIdU->GetYaxis()->SetTitle("cluster group id");
440 Name = "SVDDQM_cluster6TimeGroupIdV";
441 Title = Form("SVD cluster Time Group Id %s vs cluster time for V/N Side for coarse trigger", refFrame.Data());
442 m_clusterTimeCoarseGroupIdV = new TH2F(Name.Data(), Title.Data(), TimeBins / 2, TimeMin, TimeMax, GroupIdBins, GroupIdMin,
443 GroupIdMax);
444 m_clusterTimeCoarseGroupIdV->GetXaxis()->SetTitle("cluster time (ns)");
445 m_clusterTimeCoarseGroupIdV->GetYaxis()->SetTitle("cluster group id");
447
448 Name = "SVDDQM_cluster3TimeGroupIdU";
449 Title = Form("SVD cluster Time Group Id %s vs cluster time for U/P Side for fine trigger", refFrame.Data());
450 m_clusterTimeFineGroupIdU = new TH2F(Name.Data(), Title.Data(), TimeBins / 2, TimeMin, TimeMax, GroupIdBins, GroupIdMin,
451 GroupIdMax);
452 m_clusterTimeFineGroupIdU->GetXaxis()->SetTitle("cluster time (ns)");
453 m_clusterTimeFineGroupIdU->GetYaxis()->SetTitle("cluster group id");
455 Name = "SVDDQM_cluster3TimeGroupIdV";
456 Title = Form("SVD cluster Time Group Id %s vs cluster time for V/N Side for fine trigger", refFrame.Data());
457 m_clusterTimeFineGroupIdV = new TH2F(Name.Data(), Title.Data(), TimeBins / 2, TimeMin, TimeMax, GroupIdBins, GroupIdMin,
458 GroupIdMax);
459 m_clusterTimeFineGroupIdV->GetXaxis()->SetTitle("cluster time (ns)");
460 m_clusterTimeFineGroupIdV->GetYaxis()->SetTitle("cluster group id");
462
463 //----------------------------------------------------------------
464 // MaxBin of strips for all sensors (offline ZS)
465 //----------------------------------------------------------------
466 name = str(format("SVDDQM_StripMaxBinUAll"));
467 title = str(format("SVD U-Strip MaxBin for all sensors"));
468 m_stripMaxBinUAll = new TH1F(name.c_str(), title.c_str(), MaxBinBins, 0, MaxBinMax);
469 m_stripMaxBinUAll->GetXaxis()->SetTitle("max bin");
470 m_stripMaxBinUAll->GetYaxis()->SetTitle("count");
472 name = str(format("SVDDQM_StripMaxBinVAll"));
473 title = str(format("SVD V-Strip MaxBin for all sensors"));
474 m_stripMaxBinVAll = new TH1F(name.c_str(), title.c_str(), MaxBinBins, 0, MaxBinMax);
475 m_stripMaxBinVAll->GetXaxis()->SetTitle("max bin");
476 m_stripMaxBinVAll->GetYaxis()->SetTitle("count");
478
479 name = str(format("SVDDQM_StripMaxBinU3"));
480 title = str(format("SVD U-Strip MaxBin for layer 3 sensors"));
481 m_stripMaxBinU3 = new TH1F(name.c_str(), title.c_str(), MaxBinBins, 0, MaxBinMax);
482 m_stripMaxBinU3->GetXaxis()->SetTitle("max bin");
483 m_stripMaxBinU3->GetYaxis()->SetTitle("count");
485 name = str(format("SVDDQM_StripMaxBinV3"));
486 title = str(format("SVD V-Strip MaxBin for layer 3 sensors"));
487 m_stripMaxBinV3 = new TH1F(name.c_str(), title.c_str(), MaxBinBins, 0, MaxBinMax);
488 m_stripMaxBinV3->GetXaxis()->SetTitle("max bin");
489 m_stripMaxBinV3->GetYaxis()->SetTitle("count");
491
492 name = str(format("SVDDQM_StripMaxBinU6"));
493 title = str(format("SVD U-Strip MaxBin for layer 6 sensors"));
494 m_stripMaxBinU6 = new TH1F(name.c_str(), title.c_str(), MaxBinBins, 0, MaxBinMax);
495 m_stripMaxBinU6->GetXaxis()->SetTitle("max bin");
496 m_stripMaxBinU6->GetYaxis()->SetTitle("count");
498 name = str(format("SVDDQM_StripMaxBinV6"));
499 title = str(format("SVD V-Strip MaxBin for layer 6 sensors"));
500 m_stripMaxBinV6 = new TH1F(name.c_str(), title.c_str(), MaxBinBins, 0, MaxBinMax);
501 m_stripMaxBinV6->GetXaxis()->SetTitle("max bin");
502 m_stripMaxBinV6->GetYaxis()->SetTitle("count");
504
505 for (int i = 0; i < nSVDSensors; i++) {
506 VxdID id = gTools->getSensorIDFromSVDIndex(i);
507 int iLayer = id.getLayerNumber();
508 int iLadder = id.getLadderNumber();
509 int iSensor = id.getSensorNumber();
510 VxdID sensorID(iLayer, iLadder, iSensor);
511 SVD::SensorInfo SensorInfo = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensorID));
512 string sensorDescr = str(format("%1%_%2%_%3%") % iLayer % iLadder % iSensor);
513
514 if (m_additionalPlots) {
515 //----------------------------------------------------------------
516 // Number of fired strips per sensor
517 //----------------------------------------------------------------
518 name = str(format("SVDDQM_%1%_FiredU") % sensorDescr);
519 title = str(format("SVD Sensor %1% Number of Fired U-Strips") % sensorDescr);
520 m_firedU[i] = new TH1F(name.c_str(), title.c_str(), 50, 0, 50);
521 m_firedU[i]->GetXaxis()->SetTitle("# fired strips");
522 m_firedU[i]->GetYaxis()->SetTitle("count");
523 m_histoList->Add(m_firedU[i]);
524 name = str(format("SVDDQM_%1%_FiredV") % sensorDescr);
525 title = str(format("SVD Sensor %1% Number of Fired V-Strips") % sensorDescr);
526 m_firedV[i] = new TH1F(name.c_str(), title.c_str(), 50, 0, 50);
527 m_firedV[i]->GetXaxis()->SetTitle("# fired strips");
528 m_firedV[i]->GetYaxis()->SetTitle("count");
529 m_histoList->Add(m_firedV[i]);
530 //----------------------------------------------------------------
531 // Number of clusters per sensor
532 //----------------------------------------------------------------
533 name = str(format("SVDDQM_%1%_ClustersU") % sensorDescr);
534 title = str(format("SVD Sensor %1% Number of U-Clusters") % sensorDescr);
535 m_clustersU[i] = new TH1F(name.c_str(), title.c_str(), 20, 0, 20);
536 m_clustersU[i]->GetXaxis()->SetTitle("# clusters");
537 m_clustersU[i]->GetYaxis()->SetTitle("count");
538 m_histoList->Add(m_clustersU[i]);
539 name = str(format("SVDDQM_%1%_ClustersV") % sensorDescr);
540 title = str(format("SVD Sensor %1% Number of V-Clusters") % sensorDescr);
541 m_clustersV[i] = new TH1F(name.c_str(), title.c_str(), 20, 0, 20);
542 m_clustersV[i]->GetXaxis()->SetTitle("# clusters");
543 m_clustersV[i]->GetYaxis()->SetTitle("count");
544 m_histoList->Add(m_clustersV[i]);
545 //----------------------------------------------------------------
546 // Charge of strips
547 //----------------------------------------------------------------
548 name = str(format("SVDDQM_%1%_ADCStripU") % sensorDescr);
549 title = str(format("SVD Sensor %1% U-Strip signal in ADC Counts, all 6 APV samples") % sensorDescr);
550 m_stripSignalU[i] = new TH1F(name.c_str(), title.c_str(), 256, -0.5, 255.5);
551 m_stripSignalU[i]->GetXaxis()->SetTitle("signal ADC");
552 m_stripSignalU[i]->GetYaxis()->SetTitle("count");
554 name = str(format("SVDDQM_%1%_ADCStripV") % sensorDescr);
555 title = str(format("SVD Sensor %1% V-Strip signal in ADC Counts, all 6 APV samples") % sensorDescr);
556 m_stripSignalV[i] = new TH1F(name.c_str(), title.c_str(), 256, -0.5, 255.5);
557 m_stripSignalV[i]->GetXaxis()->SetTitle("signal ADC");
558 m_stripSignalV[i]->GetYaxis()->SetTitle("count");
560 }
561
562 //----------------------------------------------------------------
563 // Charge of clusters
564 //----------------------------------------------------------------
565 name = str(format("SVDDQM_%1%_ClusterChargeU") % sensorDescr);
566 title = str(format("SVD Sensor %1% U-Cluster Charge") % sensorDescr);
567 m_clusterChargeU[i] = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
568 m_clusterChargeU[i]->GetXaxis()->SetTitle("cluster charge [ke-]");
569 m_clusterChargeU[i]->GetYaxis()->SetTitle("count");
571 name = str(format("SVDDQM_%1%_ClusterChargeV") % sensorDescr);
572 title = str(format("SVD Sensor %1% V-Cluster Charge") % sensorDescr);
573 m_clusterChargeV[i] = new TH1F(name.c_str(), title.c_str(), ChargeBins, 0, ChargeMax);
574 m_clusterChargeV[i]->GetXaxis()->SetTitle("cluster charge [ke-]");
575 m_clusterChargeV[i]->GetYaxis()->SetTitle("count");
577 //----------------------------------------------------------------
578 // SNR of clusters
579 //----------------------------------------------------------------
580 name = str(format("SVDDQM_%1%_ClusterSNRU") % sensorDescr);
581 title = str(format("SVD Sensor %1% U-Cluster SNR") % sensorDescr);
582 m_clusterSNRU[i] = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
583 m_clusterSNRU[i]->GetXaxis()->SetTitle("cluster SNR");
584 m_clusterSNRU[i]->GetYaxis()->SetTitle("count");
585 m_histoList->Add(m_clusterSNRU[i]);
586 name = str(format("SVDDQM_%1%_ClusterSNRV") % sensorDescr);
587 title = str(format("SVD Sensor %1% V-Cluster SNR") % sensorDescr);
588 m_clusterSNRV[i] = new TH1F(name.c_str(), title.c_str(), SNRBins, 0, SNRMax);
589 m_clusterSNRV[i]->GetXaxis()->SetTitle("cluster SNR");
590 m_clusterSNRV[i]->GetYaxis()->SetTitle("count");
591 m_histoList->Add(m_clusterSNRV[i]);
592
593 //----------------------------------------------------------------
594 // Strips Counts
595 //----------------------------------------------------------------
596 name = str(format("SVDDQM_%1%_StripCountU") % sensorDescr);
597 title = str(format("SVD Sensor %1% Integrated Number of ZS5 Fired U-Strip vs Strip Number") % sensorDescr);
598 m_stripCountU[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
599 m_stripCountU[i]->GetXaxis()->SetTitle("cellID");
600 m_stripCountU[i]->GetYaxis()->SetTitle("count");
601 m_histoList->Add(m_stripCountU[i]);
602 name = str(format("SVDDQM_%1%_StripCountV") % sensorDescr);
603 title = str(format("SVD Sensor %1% Integrated Number of ZS5 Fired V-Strip vs Strip Number") % sensorDescr);
604 m_stripCountV[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
605 m_stripCountV[i]->GetXaxis()->SetTitle("cellID");
606 m_stripCountV[i]->GetYaxis()->SetTitle("count");
607 m_histoList->Add(m_stripCountV[i]);
608 //----------------------------------------------------------------
609 // Strips Counts with online ZS
610 //----------------------------------------------------------------
611 name = str(format("SVDDQM_%1%_OnlineZSStripCountU") % sensorDescr);
612 title = str(format("SVD Sensor %1% Integrated Number of online-ZS Fired U-Strip vs Strip Number") % sensorDescr);
613 m_onlineZSstripCountU[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
614 m_onlineZSstripCountU[i]->GetXaxis()->SetTitle("cellID");
615 m_onlineZSstripCountU[i]->GetYaxis()->SetTitle("count");
617 name = str(format("SVDDQM_%1%_OnlineZSStripCountV") % sensorDescr);
618 title = str(format("SVD Sensor %1% Integrated Number of online-ZS Fired V-Strip vs Strip Number") % sensorDescr);
619 m_onlineZSstripCountV[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
620 m_onlineZSstripCountV[i]->GetXaxis()->SetTitle("cellID");
621 m_onlineZSstripCountV[i]->GetYaxis()->SetTitle("count");
623
624 //----------------------------------------------------------------
625 // Strips Counts for 3 samples
626 //----------------------------------------------------------------
627 if (m_3Samples) {
628 name = str(format("SVDDQM_%1%_Strip3CountU") % sensorDescr);
629 title = str(format("SVD Sensor %1% Integrated Number of ZS5 Fired U-Strip vs Strip Number for 3 samples") % sensorDescr);
630 m_strip3SampleCountU[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
631 m_strip3SampleCountU[i]->GetXaxis()->SetTitle("cellID");
632 m_strip3SampleCountU[i]->GetYaxis()->SetTitle("count");
634 name = str(format("SVDDQM_%1%_Strip3CountV") % sensorDescr);
635 title = str(format("SVD Sensor %1% Integrated Number of ZS5 Fired V-Strip vs Strip Number for 3 samples") % sensorDescr);
636 m_strip3SampleCountV[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
637 m_strip3SampleCountV[i]->GetXaxis()->SetTitle("cellID");
638 m_strip3SampleCountV[i]->GetYaxis()->SetTitle("count");
640
641 //----------------------------------------------------------------
642 // Strips Counts with online ZS for 3 samples
643 //----------------------------------------------------------------
644 name = str(format("SVDDQM_%1%_OnlineZSStrip3CountU") % sensorDescr);
645 title = str(format("SVD Sensor %1% Integrated Number of online-ZS Fired U-Strip vs Strip Number for 3 samples") % sensorDescr);
646 m_onlineZSstrip3SampleCountU[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
647 m_onlineZSstrip3SampleCountU[i]->GetXaxis()->SetTitle("cellID");
648 m_onlineZSstrip3SampleCountU[i]->GetYaxis()->SetTitle("count");
650 name = str(format("SVDDQM_%1%_OnlineZSStrip3CountV") % sensorDescr);
651 title = str(format("SVD Sensor %1% Integrated Number of online-ZS Fired V-Strip vs Strip Number for 3 samples") % sensorDescr);
652 m_onlineZSstrip3SampleCountV[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
653 m_onlineZSstrip3SampleCountV[i]->GetXaxis()->SetTitle("cellID");
654 m_onlineZSstrip3SampleCountV[i]->GetYaxis()->SetTitle("count");
656
657 //----------------------------------------------------------------
658 // Strips Counts for 6 samples
659 //----------------------------------------------------------------
660 name = str(format("SVDDQM_%1%_Strip6CountU") % sensorDescr);
661 title = str(format("SVD Sensor %1% Integrated Number of ZS5 Fired U-Strip vs Strip Number for 6 samples") % sensorDescr);
662 m_strip6SampleCountU[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
663 m_strip6SampleCountU[i]->GetXaxis()->SetTitle("cellID");
664 m_strip6SampleCountU[i]->GetYaxis()->SetTitle("count");
666 name = str(format("SVDDQM_%1%_strip6CountV") % sensorDescr);
667 title = str(format("SVD Sensor %1% Integrated Number of ZS5 Fired V-Strip vs Strip Number for 6 samples") % sensorDescr);
668 m_strip6SampleCountV[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
669 m_strip6SampleCountV[i]->GetXaxis()->SetTitle("cellID");
670 m_strip6SampleCountV[i]->GetYaxis()->SetTitle("count");
672 //----------------------------------------------------------------
673 // Strips Counts with online ZS for 6 samples
674 //----------------------------------------------------------------
675 name = str(format("SVDDQM_%1%_OnlineZSStrip6CountU") % sensorDescr);
676 title = str(format("SVD Sensor %1% Integrated Number of online-ZS Fired U-Strip vs Strip Number for 6 samples") % sensorDescr);
677 m_onlineZSstrip6sampleCountU[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
678 m_onlineZSstrip6sampleCountU[i]->GetXaxis()->SetTitle("cellID");
679 m_onlineZSstrip6sampleCountU[i]->GetYaxis()->SetTitle("count");
681 name = str(format("SVDDQM_%1%_OnlineZSStrip6CountV") % sensorDescr);
682 title = str(format("SVD Sensor %1% Integrated Number of online-ZS Fired V-Strip vs Strip Number for 6 samples") % sensorDescr);
683 m_onlineZSstrip6sampleCountV[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
684 m_onlineZSstrip6sampleCountV[i]->GetXaxis()->SetTitle("cellID");
685 m_onlineZSstrip6sampleCountV[i]->GetYaxis()->SetTitle("count");
687 }
688
689 //----------------------------------------------------------------
690 // Strips Counts for cluster time group id = 0
691 //----------------------------------------------------------------
692 name = str(format("SVDDQM_%1%_StripCountGroupId0U") % sensorDescr);
693 title = str(format("SVD Sensor %1% Integrated NumberFired U-Strip for group Id = 0 vs Strip Number") % sensorDescr);
694 m_stripCountGroupId0U[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
695 m_stripCountGroupId0U[i]->GetXaxis()->SetTitle("cellID");
696 m_stripCountGroupId0U[i]->GetYaxis()->SetTitle("count");
698 name = str(format("SVDDQM_%1%_StripCountGroupId0V") % sensorDescr);
699 title = str(format("SVD Sensor %1% Integrated Number of Fired V-Strip for group Id = 0 vs Strip Number") % sensorDescr);
700 m_stripCountGroupId0V[i] = new TH1F(name.c_str(), title.c_str(), 768, -0.5, 767.5);
701 m_stripCountGroupId0V[i]->GetXaxis()->SetTitle("cellID");
702 m_stripCountGroupId0V[i]->GetYaxis()->SetTitle("count");
704
705 //----------------------------------------------------------------
706 // Cluster size distribution
707 //----------------------------------------------------------------
708 name = str(format("SVDDQM_%1%_ClusterSizeU") % sensorDescr);
709 title = str(format("SVD Sensor %1% U-Cluster Size") % sensorDescr);
710 m_clusterSizeU[i] = new TH1F(name.c_str(), title.c_str(), 9, 1, 10);
711 m_clusterSizeU[i]->GetXaxis()->SetTitle("cluster size");
712 m_clusterSizeU[i]->GetYaxis()->SetTitle("count");
714 name = str(format("SVDDQM_%1%_ClusterSizeV") % sensorDescr);
715 title = str(format("SVD Sensor %1% V-Cluster Size") % sensorDescr);
716 m_clusterSizeV[i] = new TH1F(name.c_str(), title.c_str(), 9, 1, 10);
717 m_clusterSizeV[i]->GetXaxis()->SetTitle("cluster size");
718 m_clusterSizeV[i]->GetYaxis()->SetTitle("count");
720 //----------------------------------------------------------------
721 // Cluster time distribution
722 //----------------------------------------------------------------
723 name = str(format("SVDDQM_%1%_ClusterTimeU") % sensorDescr);
724 Title = Form("SVD Sensor %s U-Cluster Time %s", sensorDescr.c_str(), refFrame.Data());
725 m_clusterTimeU[i] = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
726 m_clusterTimeU[i]->GetXaxis()->SetTitle("cluster time (ns)");
727 m_clusterTimeU[i]->GetYaxis()->SetTitle("count");
729 name = str(format("SVDDQM_%1%_ClusterTimeV") % sensorDescr);
730 Title = Form("SVD Sensor %s V-Cluster Time %s", sensorDescr.c_str(), refFrame.Data());
731 m_clusterTimeV[i] = new TH1F(name.c_str(), Title.Data(), TimeBins, TimeMin, TimeMax);
732 m_clusterTimeV[i]->GetXaxis()->SetTitle("cluster time (ns)");
733 m_clusterTimeV[i]->GetYaxis()->SetTitle("count");
735 }
736
737 for (int i = 0; i < nSVDChips; i++) {
738 VxdID id = gTools->getChipIDFromSVDIndex(i);
739 int iLayer = id.getLayerNumber();
740 int iLadder = id.getLadderNumber();
741 int iSensor = id.getSensorNumber();
742 int iChip = gTools->getSVDChipNumber(id);
743 int IsU = gTools->isSVDSideU(id);
744 TString AxisTicks = Form("%i_%i_%i_u%i", iLayer, iLadder, iSensor, iChip);
745 if (!IsU)
746 AxisTicks = Form("%i_%i_%i_v%i", iLayer, iLadder, iSensor, iChip);
747 m_hitMapCountsChip->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
748 m_hitMapClCountsChip->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
749 }
750
751
752
753 //----------------------------------------------------------------
754 // Additional histograms for out of ExpressReco
755 //----------------------------------------------------------------
756
757 if (m_ShowAllHistos == 1) {
758 TDirectory* dirShowAll = nullptr;
759 dirShowAll = oldDir->mkdir("SVDDQMAll");
760 dirShowAll->cd();
761
762 m_hitMapU = new TH2F*[nSVDSensors];
763 m_hitMapV = new TH2F*[nSVDSensors];
764 m_hitMapUCl = new TH1F*[nSVDSensors];
765 m_hitMapVCl = new TH1F*[nSVDSensors];
766 for (int i = 0; i < nSVDSensors; i++) {
767 VxdID id = gTools->getSensorIDFromSVDIndex(i);
768 int iLayer = id.getLayerNumber();
769 int iLadder = id.getLadderNumber();
770 int iSensor = id.getSensorNumber();
771 VxdID sensorID(iLayer, iLadder, iSensor);
772 SVD::SensorInfo SensorInfo = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensorID));
773 string sensorDescr = str(format("%1%_%2%_%3%") % iLayer % iLadder % iSensor);
774 //----------------------------------------------------------------
775 // Hitmaps: Number of strips by coordinate
776 //----------------------------------------------------------------
777 name = str(format("SVD_%1%_StripHitmapU") % sensorDescr);
778 title = str(format("SVD Sensor %1% Strip Hitmap in U") % sensorDescr);
779 int nStrips = SensorInfo.getUCells();
780 m_hitMapU[i] = new TH2F(name.c_str(), title.c_str(), nStrips, 0, nStrips, SVDShaperDigit::c_nAPVSamples, 0,
782 m_hitMapU[i]->GetXaxis()->SetTitle("u position [pitch units]");
783 m_hitMapU[i]->GetYaxis()->SetTitle("timebin [time units]");
784 m_hitMapU[i]->GetZaxis()->SetTitle("hits");
785 m_histoList->Add(m_hitMapU[i]);
786 name = str(format("SVD_%1%_StripHitmapV") % sensorDescr);
787 title = str(format("SVD Sensor %1% Strip Hitmap in V") % sensorDescr);
788 nStrips = SensorInfo.getVCells();
789 m_hitMapV[i] = new TH2F(name.c_str(), title.c_str(), nStrips, 0, nStrips, SVDShaperDigit::c_nAPVSamples, 0,
791 m_hitMapV[i]->GetXaxis()->SetTitle("v position [pitch units]");
792 m_hitMapV[i]->GetYaxis()->SetTitle("timebin [time units]");
793 m_hitMapV[i]->GetZaxis()->SetTitle("hits");
794 m_histoList->Add(m_hitMapV[i]);
795 //----------------------------------------------------------------
796 // Hitmaps: Number of clusters by coordinate
797 //----------------------------------------------------------------
798 name = str(format("SVD_%1%_HitmapClstU") % sensorDescr);
799 title = str(format("SVD Sensor %1% Hitmap Clusters in U") % sensorDescr);
800 nStrips = SensorInfo.getUCells();
801 m_hitMapUCl[i] = new TH1F(name.c_str(), title.c_str(), nStrips, 0, nStrips);
802 m_hitMapUCl[i]->GetXaxis()->SetTitle("u position [pitch units]");
803 m_hitMapUCl[i]->GetYaxis()->SetTitle("hits");
804 m_histoList->Add(m_hitMapUCl[i]);
805 name = str(format("SVD_%1%_HitmapClstV") % sensorDescr);
806 title = str(format("SVD Sensor %1% Hitmap Clusters in V") % sensorDescr);
807 nStrips = SensorInfo.getVCells();
808 m_hitMapVCl[i] = new TH1F(name.c_str(), title.c_str(), nStrips, 0, nStrips);
809 m_hitMapVCl[i]->GetXaxis()->SetTitle("v position [pitch units]");
810 m_hitMapVCl[i]->GetYaxis()->SetTitle("hits");
811 m_histoList->Add(m_hitMapVCl[i]);
812 }
813 }
814
815 oldDir->cd();
816}
817
818
820{
821 // Register histograms (calls back defineHisto)
822 REG_HISTOGRAM
823
824 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
825 if (gTools->getNumberOfSVDLayers() != 0) {
826 //Register collections
830
831 storeSVDClusters.isOptional();
832 storeSVDShaperDigits.isOptional();
833 m_svdEventInfo.isOptional();
834 storeNoZSSVDShaperDigits.isOptional();
835
836 //Store names to speed up creation later
837 m_storeSVDShaperDigitsName = storeSVDShaperDigits.getName();
838 }
839
840 m_objTrgSummary.isOptional();
841}
842
844{
845 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
846 if (gTools->getNumberOfSVDLayers() == 0) return;
847
848
849 StoreObjPtr<EventMetaData> evtMetaData;
850 m_expNumber = evtMetaData->getExperiment();
851 m_runNumber = evtMetaData->getRun();
852
853 // Add experiment and run number to the title of selected histograms (CR shifter plots)
854 TString runID = TString::Format(" ~ Exp%d Run%d", m_expNumber, m_runNumber);
855 TObject* obj;
856 TIter nextH(m_histoList);
857 while ((obj = nextH()))
858 if (obj->InheritsFrom("TH1")) {
859
860 TString tmp = (TString)obj->GetTitle();
861 Int_t pos = tmp.Last('~');
862 if (pos == -1) pos = tmp.Length() + 2;
863
864 TString title = tmp(0, pos - 2);
865 ((TH1F*)obj)->SetTitle(title + runID);
866 ((TH1F*)obj)->Reset();
867 }
868}
869
871{
872 //check HLT decision and increase number of events only if the event has been accepted
875 if (!eventAccepted) return;
876 }
877 m_nEvents->Fill(0);
878
879 int nSamples = 0;
880 if (m_svdEventInfo.isValid())
881 nSamples = m_svdEventInfo->getNSamples();
882 else
883 return;
884
885 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
886 if (gTools->getNumberOfSVDLayers() == 0) return;
887
888
889 const StoreArray<SVDShaperDigit> storeNoZSSVDShaperDigits(m_storeNoZSSVDShaperDigitsName);
890 const StoreArray<SVDShaperDigit> storeSVDShaperDigits(m_storeSVDShaperDigitsName);
891 const StoreArray<SVDCluster> storeSVDClusters(m_storeSVDClustersName);
892
893 if (!storeSVDShaperDigits.isValid() || !storeSVDShaperDigits.getEntries()) {
894 return;
895 }
896
897 int firstSVDLayer = gTools->getFirstSVDLayer();
898 int lastSVDLayer = gTools->getLastSVDLayer();
899 int nSVDSensors = gTools->getNumberOfSVDSensors();
900
901 // Fired strips offline ZS
902 vector< set<int> > uStrips(nSVDSensors); // sets to eliminate multiple samples per strip
903 vector< set<int> > vStrips(nSVDSensors);
904 for (const SVDShaperDigit& digitIn : storeSVDShaperDigits) {
905 int iLayer = digitIn.getSensorID().getLayerNumber();
906 if ((iLayer < firstSVDLayer) || (iLayer > lastSVDLayer)) continue;
907 int iLadder = digitIn.getSensorID().getLadderNumber();
908 int iSensor = digitIn.getSensorID().getSensorNumber();
909 VxdID sensorID(iLayer, iLadder, iSensor);
910 int index = gTools->getSVDSensorIndex(sensorID);
911 SVD::SensorInfo SensorInfo = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensorID));
912 if (digitIn.isUStrip()) {
913
914 //fill strip count first
915 if (m_stripCountU[index] != nullptr) m_stripCountU[index]->Fill(digitIn.getCellID());
916
917 if (m_3Samples) {
918 if (nSamples == 3) {
919 if (m_strip3SampleCountU[index] != nullptr) m_strip3SampleCountU[index]->Fill(digitIn.getCellID());
920 } else {
921 if (m_strip6SampleCountU[index] != nullptr) m_strip6SampleCountU[index]->Fill(digitIn.getCellID());
922 }
923 }
924 //fill max bin
925 if (m_stripMaxBinUAll != nullptr) m_stripMaxBinUAll->Fill(digitIn.getMaxTimeBin());
926 if (iLayer == 3)
927 if (m_stripMaxBinU3 != nullptr) m_stripMaxBinU3->Fill(digitIn.getMaxTimeBin());
928 if (iLayer == 6)
929 if (m_stripMaxBinU6 != nullptr) m_stripMaxBinU6->Fill(digitIn.getMaxTimeBin());
930
931 uStrips.at(index).insert(digitIn.getCellID());
932 int Chip = (int)(digitIn.getCellID() / gTools->getSVDChannelsPerChip()) + 1;
933 int indexChip = gTools->getSVDChipIndex(sensorID, kTRUE, Chip);
934 // 6-to-1 relation weights are equal to digit signals, modulo rounding error
935 SVDShaperDigit::APVFloatSamples samples = digitIn.getSamples();
936 int isSample = 0;
937 for (size_t i = 0; i < SVDShaperDigit::c_nAPVSamples; ++i) {
939 if (m_stripSignalU[index] != nullptr) m_stripSignalU[index]->Fill(samples[i]);
940 if (samples[i] > m_CutSVDCharge) {
941 isSample = 1;
942 if (m_ShowAllHistos == 1) {
943 if (m_hitMapU[index] != nullptr) m_hitMapU[index]->Fill(digitIn.getCellID(), i);
944 }
945 }
946 }
947 if (isSample) {
948 if (m_hitMapCountsU != nullptr) m_hitMapCountsU->Fill(index);
949 if (m_hitMapCountsChip != nullptr) m_hitMapCountsChip->Fill(indexChip);
950 }
951 } else {
952 //fill strip count first
953 if (m_stripCountV[index] != nullptr) m_stripCountV[index]->Fill(digitIn.getCellID());
954
955 if (m_3Samples) {
956 if (nSamples == 3) {
957 if (m_strip3SampleCountV[index] != nullptr) m_strip3SampleCountV[index]->Fill(digitIn.getCellID());
958 } else {
959 if (m_strip6SampleCountV[index] != nullptr) m_strip6SampleCountV[index]->Fill(digitIn.getCellID());
960 }
961 }
962
963 //fill max bin
964 if (m_stripMaxBinVAll != nullptr) m_stripMaxBinVAll->Fill(digitIn.getMaxTimeBin());
965
966 if (iLayer == 3)
967 if (m_stripMaxBinV3 != nullptr) m_stripMaxBinV3->Fill(digitIn.getMaxTimeBin());
968 if (iLayer == 6)
969 if (m_stripMaxBinV6 != nullptr) m_stripMaxBinV6->Fill(digitIn.getMaxTimeBin());
970
971 vStrips.at(index).insert(digitIn.getCellID());
972 int Chip = (int)(digitIn.getCellID() / gTools->getSVDChannelsPerChip()) + 1;
973 int indexChip = gTools->getSVDChipIndex(sensorID, kFALSE, Chip);
974 // 6-to-1 relation weights are equal to digit signals, modulo rounding error
975 SVDShaperDigit::APVFloatSamples samples = digitIn.getSamples();
976 int isSample = 0;
977 for (size_t i = 0; i < SVDShaperDigit::c_nAPVSamples; ++i) {
979 if (m_stripSignalV[index] != nullptr) m_stripSignalV[index]->Fill(samples[i]);
980 if (samples[i] > m_CutSVDCharge) {
981 isSample = 1;
982 if (m_ShowAllHistos == 1) {
983 if (m_hitMapV[index] != nullptr) m_hitMapV[index]->Fill(digitIn.getCellID(), i);
984 }
985 }
986 }
987 if (isSample) {
988 if (m_hitMapCountsV != nullptr) m_hitMapCountsV->Fill(index);
989 if (m_hitMapCountsChip != nullptr) m_hitMapCountsChip->Fill(indexChip);
990 }
991 }
992 }
993 if (m_additionalPlots) {
994 for (int i = 0; i < nSVDSensors; i++) {
995 if ((m_firedU[i] != nullptr) && (uStrips[i].size() > 0))
996 m_firedU[i]->Fill(uStrips[i].size());
997 if ((m_firedV[i] != nullptr) && (vStrips[i].size() > 0))
998 m_firedV[i]->Fill(vStrips[i].size());
999 }
1000 }
1001
1002 // Fired strips ONLINE ZS
1003 if (storeNoZSSVDShaperDigits.isValid())
1004 for (const SVDShaperDigit& digitIn : storeNoZSSVDShaperDigits) {
1005 int iLayer = digitIn.getSensorID().getLayerNumber();
1006 if ((iLayer < firstSVDLayer) || (iLayer > lastSVDLayer)) continue;
1007 int iLadder = digitIn.getSensorID().getLadderNumber();
1008 int iSensor = digitIn.getSensorID().getSensorNumber();
1009 VxdID sensorID(iLayer, iLadder, iSensor);
1010 int index = gTools->getSVDSensorIndex(sensorID);
1011 SVD::SensorInfo SensorInfo = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensorID));
1012 if (digitIn.isUStrip()) {
1013 if (m_onlineZSstripCountU[index] != nullptr) m_onlineZSstripCountU[index]->Fill(digitIn.getCellID());
1014 if (m_3Samples) {
1015 if (nSamples == 3) {
1016 if (m_onlineZSstrip3SampleCountU[index] != nullptr) m_onlineZSstrip3SampleCountU[index]->Fill(digitIn.getCellID());
1017 } else {
1018 if (m_onlineZSstrip6sampleCountU[index] != nullptr) m_onlineZSstrip6sampleCountU[index]->Fill(digitIn.getCellID());
1019 }
1020 }
1021 } else {
1022 if (m_onlineZSstripCountV[index] != nullptr) m_onlineZSstripCountV[index]->Fill(digitIn.getCellID());
1023 if (m_3Samples) {
1024 if (nSamples == 3) {
1025 if (m_onlineZSstrip3SampleCountV[index] != nullptr) m_onlineZSstrip3SampleCountV[index]->Fill(digitIn.getCellID());
1026 } else {
1027 if (m_onlineZSstrip6sampleCountV[index] != nullptr) m_onlineZSstrip6sampleCountV[index]->Fill(digitIn.getCellID());
1028 }
1029 }
1030 }
1031 }
1032
1033 vector< set<int> > countsU(nSVDSensors); // sets to eliminate multiple samples per strip
1034 vector< set<int> > countsV(nSVDSensors);
1035 // Hitmaps, Charge, Seed, Size, Time, ...
1036 for (const SVDCluster& cluster : storeSVDClusters) {
1037 if (cluster.getCharge() < m_CutSVDClusterCharge) continue;
1038 int iLayer = cluster.getSensorID().getLayerNumber();
1039 if ((iLayer < firstSVDLayer) || (iLayer > lastSVDLayer)) continue;
1040 int iLadder = cluster.getSensorID().getLadderNumber();
1041 int iSensor = cluster.getSensorID().getSensorNumber();
1042 VxdID sensorID(iLayer, iLadder, iSensor);
1043 int index = gTools->getSVDSensorIndex(sensorID);
1044 SVD::SensorInfo SensorInfo = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensorID));
1045
1046 float time = cluster.getClsTime();
1047 if (m_desynchSVDTime && m_svdEventInfo.isValid())
1048 time = time - m_svdEventInfo->getSVD2FTSWTimeShift(cluster.getFirstFrame());
1049
1050 vector<int> vec = cluster.getTimeGroupId();
1051 auto minElement = min_element(vec.begin(), vec.end());
1052 int groupId = -1;
1053 if (vec.size() > 0) {
1054 groupId = *minElement;
1055
1056 if (cluster.isUCluster()) {
1057 if (m_clusterTimeGroupIdU != nullptr) m_clusterTimeGroupIdU->Fill(time, groupId);
1058 if (m_objTrgSummary.isValid()) {
1059 int trgQuality = m_objTrgSummary->getTimQuality();
1060 if (trgQuality == 1)
1061 if (m_clusterTimeCoarseGroupIdU != nullptr) m_clusterTimeCoarseGroupIdU->Fill(time, groupId);
1062 if (trgQuality == 2)
1063 if (m_clusterTimeFineGroupIdU != nullptr) m_clusterTimeFineGroupIdU->Fill(time, groupId);
1064 }
1065
1066 } else {
1067 if (m_clusterTimeGroupIdV != nullptr) m_clusterTimeGroupIdV->Fill(time, groupId);
1068 if (m_objTrgSummary.isValid()) {
1069 int trgQuality = m_objTrgSummary->getTimQuality();
1070 if (trgQuality == 1)
1071 if (m_clusterTimeCoarseGroupIdV != nullptr) m_clusterTimeCoarseGroupIdV->Fill(time, groupId);
1072 if (trgQuality == 2)
1073 if (m_clusterTimeFineGroupIdV != nullptr) m_clusterTimeFineGroupIdV->Fill(time, groupId);
1074 }
1075 }
1076 }
1077
1078 if (cluster.isUCluster()) {
1079 countsU.at(index).insert(SensorInfo.getUCellID(cluster.getPosition()));
1080 int indexChip = gTools->getSVDChipIndex(sensorID, kTRUE,
1081 (int)(SensorInfo.getUCellID(cluster.getPosition()) / gTools->getSVDChannelsPerChip()) + 1);
1082 if (m_hitMapClCountsU != nullptr) m_hitMapClCountsU->Fill(index);
1083 if (m_hitMapClCountsChip != nullptr) m_hitMapClCountsChip->Fill(indexChip);
1084 if (m_clusterChargeU[index] != nullptr) m_clusterChargeU[index]->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1085 if (m_clusterSNRU[index] != nullptr) m_clusterSNRU[index]->Fill(cluster.getSNR());
1086 if (m_clusterChargeUAll != nullptr) m_clusterChargeUAll->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1087 if (m_clusterSNRUAll != nullptr) m_clusterSNRUAll->Fill(cluster.getSNR());
1088 if (m_clusterSizeU[index] != nullptr) m_clusterSizeU[index]->Fill(cluster.getSize());
1089 if (m_clusterTimeU[index] != nullptr) m_clusterTimeU[index]->Fill(time);
1090 if (m_clusterTimeUAll != nullptr) m_clusterTimeUAll->Fill(time);
1091 if (iLayer == 3) {
1092 if (m_clusterChargeU3 != nullptr) m_clusterChargeU3->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1093 if (m_clusterSNRU3 != nullptr) m_clusterSNRU3->Fill(cluster.getSNR());
1094 if (m_clusterTimeU3 != nullptr) m_clusterTimeU3->Fill(time);
1095 if (m_3Samples) {
1096 if (nSamples == 3) {
1097 if (m_cluster3SampleTimeU3 != nullptr) m_cluster3SampleTimeU3->Fill(time);
1098 } else {
1099 if (m_cluster6SampleTimeU3 != nullptr) m_cluster6SampleTimeU3->Fill(time);
1100 }
1101 }
1102 } else {
1103 if (m_clusterChargeU456 != nullptr) m_clusterChargeU456->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1104 if (m_clusterSNRU456 != nullptr) m_clusterSNRU456->Fill(cluster.getSNR());
1105 if (m_clusterTimeU456 != nullptr) m_clusterTimeU456->Fill(time);
1106 if (m_3Samples) {
1107 if (nSamples == 3) {
1108 if (m_cluster3SampleTimeU456 != nullptr) m_cluster3SampleTimeU456->Fill(time);
1109 } else {
1110 if (m_cluster6SampleTimeU456 != nullptr) m_cluster6SampleTimeU456->Fill(time);
1111 }
1112 }
1113 }
1114
1115 if (m_ShowAllHistos == 1)
1116 if (m_hitMapUCl[index] != nullptr) m_hitMapUCl[index]->Fill(SensorInfo.getUCellID(cluster.getPosition()));
1117
1118 // groupId for U side
1119 if (groupId == 0) {
1120 for (const SVDShaperDigit& digitIn : cluster.getRelationsTo<SVDShaperDigit>(m_storeSVDShaperDigitsName)) {
1121 if (m_stripCountGroupId0U != nullptr) m_stripCountGroupId0U[index]->Fill(digitIn.getCellID());
1122 }
1123 }
1124 } else {
1125 countsV.at(index).insert(SensorInfo.getVCellID(cluster.getPosition()));
1126 int indexChip = gTools->getSVDChipIndex(sensorID, kFALSE,
1127 (int)(SensorInfo.getVCellID(cluster.getPosition()) / gTools->getSVDChannelsPerChip()) + 1);
1128 if (m_hitMapClCountsV != nullptr) m_hitMapClCountsV->Fill(index);
1129 if (m_hitMapClCountsChip != nullptr) m_hitMapClCountsChip->Fill(indexChip);
1130 if (m_clusterChargeV[index] != nullptr) m_clusterChargeV[index]->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1131 if (m_clusterSNRV[index] != nullptr) m_clusterSNRV[index]->Fill(cluster.getSNR());
1132 if (m_clusterChargeVAll != nullptr) m_clusterChargeVAll->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1133 if (m_clusterSNRVAll != nullptr) m_clusterSNRVAll->Fill(cluster.getSNR());
1134 if (m_clusterSizeV[index] != nullptr) m_clusterSizeV[index]->Fill(cluster.getSize());
1135 if (m_clusterTimeV[index] != nullptr) m_clusterTimeV[index]->Fill(time);
1136 if (m_clusterTimeVAll != nullptr) m_clusterTimeVAll->Fill(time);
1137 if (iLayer == 3) {
1138 if (m_clusterChargeV3 != nullptr) m_clusterChargeV3->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1139 if (m_clusterSNRV3 != nullptr) m_clusterSNRV3->Fill(cluster.getSNR());
1140 if (m_clusterTimeV3 != nullptr) m_clusterTimeV3->Fill(time);
1141 if (m_3Samples) {
1142 if (nSamples == 3) {
1143 if (m_cluster3SampleTimeV3 != nullptr) m_cluster3SampleTimeV3->Fill(time);
1144 } else {
1145 if (m_cluster6SampleTimeV3 != nullptr) m_cluster6SampleTimeV3->Fill(time);
1146 }
1147 }
1148 } else {
1149 if (m_clusterChargeV456 != nullptr) m_clusterChargeV456->Fill(cluster.getCharge() / 1000.0); // in kelectrons
1150 if (m_clusterSNRV456 != nullptr) m_clusterSNRV456->Fill(cluster.getSNR());
1151 if (m_clusterTimeV456 != nullptr) m_clusterTimeV456->Fill(time);
1152 if (m_3Samples) {
1153 if (nSamples == 3) {
1154 if (m_cluster3SampleTimeV456 != nullptr) m_cluster3SampleTimeV456->Fill(time);
1155 } else {
1156 if (m_cluster6SampleTimeV456 != nullptr) m_cluster6SampleTimeV456->Fill(time);
1157 }
1158 }
1159 }
1160 if (m_ShowAllHistos == 1)
1161 if (m_hitMapVCl[index] != nullptr) m_hitMapVCl[index]->Fill(SensorInfo.getVCellID(cluster.getPosition()));
1162
1163 // groupId for V side
1164 if (groupId == 0) {
1165 for (const SVDShaperDigit& digitIn : cluster.getRelationsTo<SVDShaperDigit>(m_storeSVDShaperDigitsName)) {
1166 if (m_stripCountGroupId0V != nullptr) m_stripCountGroupId0V[index]->Fill(digitIn.getCellID());
1167 }
1168 }
1169 }
1170 }
1171 if (m_additionalPlots) {
1172 for (int i = 0; i < nSVDSensors; i++) {
1173 if ((m_clustersU[i] != nullptr) && (countsU[i].size() > 0))
1174 m_clustersU[i]->Fill(countsU[i].size());
1175 if ((m_clustersV[i] != nullptr) && (countsV[i].size() > 0))
1176 m_clustersV[i]->Fill(countsV[i].size());
1177 }
1178 }
1179}
1180
1181
1183{
1184 // m_histoList->Delete();
1185 delete m_histoList;
1186
1187}
HistoModule()
Constructor.
Definition HistoModule.h:32
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
The SVD Cluster class This class stores all information about reconstructed SVD clusters.
Definition SVDCluster.h:29
TH1F ** m_hitMapUCl
Hitmaps clusters for u.
TH1F ** m_clustersV
number of v clusters per event
TH1F * m_clusterSNRVAll
v SNR of clusters for all sensors
TH1F ** m_strip3SampleCountV
v strip count for 3 samples
TH1F * m_clusterTimeV456
v Time of clusters for layer 4,5,6 sensors
bool m_useParamFromDB
if true read back from DB configuration parameters
TH1F * m_clusterSNRUAll
u SNR of clusters for all sensors
TH1F * m_clusterChargeU3
u charge of clusters for layer 3 sensors
void initialize() override final
Module function initialize.
TH2F * m_clusterTimeGroupIdV
time group id for V side
TH1F * m_clusterSNRV3
v SNR of clusters for layer 3 sensors
TH1F * m_hitMapCountsV
Hitmaps v of Digits.
TH2F ** m_hitMapU
Hitmaps pixels for u.
TH1F ** m_onlineZSstrip6sampleCountU
u strip count (online Zero Suppression) for 6 samples
std::string m_storeNoZSSVDShaperDigitsName
not zero-suppressed SVDShaperDigits StoreArray name
bool m_skipRejectedEvents
if true skip events rejected by HLT
TH1F * m_stripMaxBinU6
u MaxBin of strips for layer 6 sensors (offline Zero Suppression)
TH1F * m_clusterTimeUAll
u time of clusters for all sensors
TH1F ** m_onlineZSstrip3SampleCountV
v strip count (online Zero Suppression for 3 samples
TH1F * m_hitMapClCountsU
Hitmaps u of Clusters.
TH2F * m_clusterTimeFineGroupIdV
time group id for V side for fine trigger
float m_CutSVDCharge
cut for accepting strips to hitmap histogram default = 0 ADU
TH1F * m_cluster6SampleTimeU3
u Time of clusters for layer 3 sensors for 6 samples
StoreObjPtr< SVDEventInfo > m_svdEventInfo
SVDEventInfo data object.
TH1F * m_hitMapCountsU
Hitmaps u of Digits.
TH1F ** m_clusterSNRV
v SNR of clusters per sensor
TH1F ** m_clusterChargeV
v charge of clusters
TH1F * m_cluster3SampleTimeU3
u Time of clusters for layer 3 sensors for 3 samples
TH1F ** m_stripSignalU
u charge of strips
TH1F * m_clusterChargeUAll
u charge of clusters for all sensors
TH1F * m_cluster3SampleTimeU456
u Time of clusters for layer 4,5,6 sensors for 3 samples
TH1F * m_clusterChargeU456
u charge of clusters for layer 4,5,6 sensors
void defineHisto() override final
Histogram definitions such as TH1(), TH2(), TNtuple(), TTree()....
TH1F * m_cluster3SampleTimeV456
v Time of clusters for layer 4,5,6 sensors for 3 samples
TH1F * m_clusterTimeV3
v Time of clusters for layer 3 sensors
DBObjPtr< SVDDQMPlotsConfiguration > m_svdPlotsConfig
SVD DQM plots configuration.
TH1F ** m_clusterSNRU
u SNR of clusters per sensor
TH1F ** m_hitMapVCl
Hitmaps clusters for v.
TH1F ** m_strip3SampleCountU
u strip count for 3 samples
void terminate() override final
Module function terminate.
TH1F * m_stripMaxBinV3
v MaxBin of strips for layer 3 sensors (offline Zero Suppression)
TH1F * m_clusterTimeVAll
v time of clusters for all sensors
TH2F * m_clusterTimeCoarseGroupIdU
time group id for U side for coarse trigger
void event() override final
Module function event.
TH1F ** m_onlineZSstripCountV
v strip count (online Zero Suppression
TH1F ** m_stripSignalV
v charge of strips
StoreObjPtr< TRGSummary > m_objTrgSummary
Trigger Summary data object.
std::string m_storeSVDShaperDigitsName
SVDShaperDigits StoreArray name.
std::string m_histogramDirectoryName
Name of the histogram directory in ROOT file.
TH1F * m_clusterChargeVAll
v charge of clusters for all sensors
TH1F * m_clusterSNRU3
u SNR of clusters for layer 3 sensors
TH1F * m_clusterSNRV456
v SNR of clusters for layer 4,5,6 sensors
TH1F * m_stripMaxBinUAll
u MaxBin of strips for all sensors (offline Zero Suppression)
bool m_3Samples
if true enable 3 samples histograms analysis
TList * m_histoList
list of cumulative histograms
TH1F ** m_clusterChargeU
u charge of clusters
TH1F * m_hitMapCountsChip
Hitmaps of digits on chips.
TH1F * m_cluster3SampleTimeV3
v Time of clusters for layer 3 sensors for 3 samples
TH1F * m_clusterChargeV3
v charge of clusters for layer 3 sensors
TH2F ** m_hitMapV
Hitmaps pixels for v.
TH1F * m_stripMaxBinV6
v MaxBin of strips for layer 6 sensors (offline Zero Suppression)
TH1F ** m_stripCountGroupId0V
V strip count for cluster time group Id = 0.
TH1F ** m_strip6SampleCountV
v strip count for 3 samples
TH1F * m_cluster6SampleTimeV456
v Time of clusters for layer 4,5,6 sensors for 6 samples
TH1F ** m_firedU
Fired u strips per event.
TH1F * m_clusterTimeU3
u Time of clusters for layer 3 sensors
TH1F * m_cluster6SampleTimeV3
v Time of clusters for layer 3 sensors for 6 samples
TH1F * m_clusterSNRU456
u SNR of clusters for layer 4,5,6 sensors
void beginRun() override final
Module function beginRun.
TH1F ** m_onlineZSstripCountU
u strip count (online Zero Suppression)
int m_ShowAllHistos
Flag to show all histos in DQM, default = 0 (do not show)
bool m_additionalPlots
additional plots flag
TH1F * m_stripMaxBinU3
u MaxBin of strips for layer 3 sensors (offline Zero Suppression)
bool m_desynchSVDTime
if TRUE: svdTime back in SVD time reference
TH1F * m_stripMaxBinVAll
v MaxBin of strips for all sensors (offline Zero Suppression)
TH1F * m_hitMapClCountsChip
Hitmaps of clusters on chips.
TH1F ** m_strip6SampleCountU
u strip count for 6 samples
TH1F ** m_onlineZSstrip3SampleCountU
u strip count (online Zero Suppression) for 3 samples
TH1F ** m_onlineZSstrip6sampleCountV
v strip count (online Zero Suppression for 6 samples
StoreObjPtr< SoftwareTriggerResult > m_resultStoreObjectPointer
Store Object for reading the trigger decision.
float m_CutSVDClusterCharge
cut for accepting clusters to hitmap histogram, default = 0 ke-
TH2F * m_clusterTimeGroupIdU
time group id for U side
std::string m_storeSVDClustersName
SVDClusters StoreArray name.
TH1F * m_clusterTimeU456
u Time of clusters for layer 4,5,6 sensors
TH2F * m_clusterTimeCoarseGroupIdV
time group id for V side for coarse trigger
TH1F ** m_firedV
Fired v strips per event.
TH1F * m_clusterChargeV456
v charge of clusters for layer 4,5,6 sensors
TH1F ** m_stripCountGroupId0U
U strip count for cluster time group Id = 0.
TH1F * m_hitMapClCountsV
Hitmaps v of Clusters.
TH1F ** m_clustersU
number of u clusters per event
TH2F * m_clusterTimeFineGroupIdU
time group id for U side for fine trigger
TH1F * m_cluster6SampleTimeU456
u Time of clusters for layer 4,5,6 sensors for 6 samples
The SVD ShaperDigit class.
static const std::size_t c_nAPVSamples
Number of APV samples stored.
std::array< APVFloatSampleType, c_nAPVSamples > APVFloatSamples
array of APVFloatSampleType objects
Specific implementation of SensorInfo for SVD Sensors which provides additional sensor specific infor...
Definition SensorInfo.h:25
static bool getFinalTriggerDecision(const SoftwareTriggerResult &result, bool forgetTotalResult=false)
Calculate the final cut decision using all "total_results" of all sub triggers in the software trigge...
const std::string & getName() const
Return name under which the object is saved in the DataStore.
bool isOptional(const std::string &name="")
Tell the DataStore about an optional input.
Accessor to arrays stored in the data store.
Definition StoreArray.h:113
bool isValid() const
Check whether the array was registered.
Definition StoreArray.h:288
int getEntries() const
Get the number of objects in the array.
Definition StoreArray.h:216
Type-safe access to single objects in the data store.
Definition StoreObjPtr.h:96
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a reference to the SensorInfo of a given SensorID.
Definition GeoCache.cc:67
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition GeoCache.cc:214
const GeoTools * getGeoTools()
Return a raw pointer to a GeoTools object.
Definition GeoCache.h:142
unsigned short getNumberOfSVDSensors() const
Get number of SVD sensors.
Definition GeoTools.h:138
int getVCells() const
Return number of pixel/strips in v direction.
int getUCells() const
Return number of pixel/strips in u direction.
int getVCellID(double v, bool clamp=false) const
Return the corresponding pixel/strip ID of a given v coordinate.
int getUCellID(double u, double v=0, bool clamp=false) const
Return the corresponding pixel/strip ID of a given u coordinate.
Class to uniquely identify a any structure of the PXD and SVD.
Definition VxdID.h:33
baseType getLayerNumber() const
Get the layer id.
Definition VxdID.h:96
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Abstract base class for different kinds of events.
STL namespace.