Belle II Software  release-08-01-10
SVDLocalCalibrationsMonitorModule.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/svdCalibration/SVDLocalCalibrationsMonitorModule.h>
10 #include <vxd/geometry/GeoCache.h>
11 #include <svd/geometry/SensorInfo.h>
12 #include <framework/datastore/StoreObjPtr.h>
13 #include <framework/dataobjects/EventMetaData.h>
14 
15 using namespace Belle2;
16 
17 //-----------------------------------------------------------------
18 // Register the Module
19 //-----------------------------------------------------------------
20 REG_MODULE(SVDLocalCalibrationsMonitor);
21 
22 //-----------------------------------------------------------------
23 // Implementation
24 //-----------------------------------------------------------------
25 
27 {
28  // Set module properties
29  setDescription("Module to produce a list of histograms showing the uploaded calibration constants");
30 
31  // Parameter definitions
32  addParam("outputFileName", m_rootFileName, "Name of output root file.", std::string("SVDLocalCalibrationMonitor_output.root"));
33 }
34 
36 {
37 
38  // create new root file
39  m_rootFilePtr = new TFile(m_rootFileName.c_str(), "RECREATE");
40 
41  //tree initialization
42  m_tree = new TTree("calibLocal", "RECREATE");
43  b_exp = m_tree->Branch("exp", &m_exp, "exp/i");
44  b_run = m_tree->Branch("run", &m_run, "run/i");
45  b_date = m_tree->Branch("date", m_date, "date/C");
46  b_hv = m_tree->Branch("hv", &m_hv, "hv/F");
47  b_layer = m_tree->Branch("layer", &m_layer, "layer/i");
48  b_ladder = m_tree->Branch("ladder", &m_ladder, "ladder/i");
49  b_sensor = m_tree->Branch("sensor", &m_sensor, "sensor/i");
50  b_side = m_tree->Branch("side", &m_side, "side/i");
51  b_maskAVE = m_tree->Branch("maskAVE", &m_maskAVE, "maskAVE/F");
52  b_hotstripsAVE = m_tree->Branch("hotstripsAVE", &m_hotstripsAVE, "hotstripsAVE/F");
53  b_pedestalAVE = m_tree->Branch("pedestalAVE", &m_pedestalAVE, "pedestalAVE/F");
54  b_pedestalRMS = m_tree->Branch("pedestalRMS", &m_pedestalRMS, "pedestalRMS/F");
55  b_noiseAVE = m_tree->Branch("noiseAVE", &m_noiseAVE, "noiseAVE/F");
56  b_noiseRMS = m_tree->Branch("noiseRMS", &m_noiseRMS, "noiseRMS/F");
57  b_noiseElAVE = m_tree->Branch("noiseElAVE", &m_noiseElAVE, "noiseElAVE/F");
58  b_noiseElRMS = m_tree->Branch("noiseElRMS", &m_noiseElRMS, "noiseElRMS/F");
59  b_occupancyAVE = m_tree->Branch("occupancyAVE", &m_occupancyAVE, "occupancyAVE/F");
60  b_occupancyRMS = m_tree->Branch("occupancyRMS", &m_occupancyRMS, "occupancyRMS/F");
61  b_gainAVE = m_tree->Branch("gainAVE", &m_gainAVE, "gainAVE/F");
62  b_gainRMS = m_tree->Branch("gainRMS", &m_gainRMS, "gainRMS/F");
63  b_calPeakADCAVE = m_tree->Branch("calPeakADCAVE", &m_calPeakADCAVE, "calPeakADCAVE/F");
64  b_calPeakADCRMS = m_tree->Branch("calPeakADCRMS", &m_calPeakADCRMS, "calPeakADCRMS/F");
65  b_calPeakTimeAVE = m_tree->Branch("calPeakTimeAVE", &m_calPeakTimeAVE, "calPeakTimeAVE/F");
66  b_calPeakTimeRMS = m_tree->Branch("calPeakTimeRMS", &m_calPeakTimeRMS, "calPeakTimeRMS/F");
67  b_pulseWidthAVE = m_tree->Branch("pulseWidthAVE", &m_pulseWidthAVE, "pulseWidthAVE/F");
68  b_pulseWidthRMS = m_tree->Branch("pulseWidthRMS", &m_pulseWidthRMS, "pulseWidthRMS/F");
69 
70  m_treeDetailed = new TTree("calibLocalDetailed", "RECREATE");
71  b_exp = m_treeDetailed->Branch("exp", &m_exp, "exp/i");
72  b_run = m_treeDetailed->Branch("run", &m_run, "run/i");
73  b_date = m_treeDetailed->Branch("date", m_date, "date/C");
74  b_hv = m_treeDetailed->Branch("hv", &m_hv, "hv/F");
75  b_layer = m_treeDetailed->Branch("layer", &m_layer, "layer/i");
76  b_ladder = m_treeDetailed->Branch("ladder", &m_ladder, "ladder/i");
77  b_sensor = m_treeDetailed->Branch("sensor", &m_sensor, "sensor/i");
78  b_side = m_treeDetailed->Branch("side", &m_side, "side/i");
79  b_strip = m_treeDetailed->Branch("strip", &m_strip, "strip/i");
80  b_mask = m_treeDetailed->Branch("mask", &m_mask, "mask/F");
81  b_hotstrips = m_treeDetailed->Branch("hotstrips", &m_hotstrips, "hotstrips/F");
82  b_noise = m_treeDetailed->Branch("noise", &m_noise, "noise/F");
83  b_occupancy = m_treeDetailed->Branch("occupancy", &m_occupancy, "occupancy/F");
84  b_noiseEl = m_treeDetailed->Branch("noiseEl", &m_noiseEl, "noiseEl/F");
85  b_gain = m_treeDetailed->Branch("gain", &m_gain, "gain/F");
86  b_pedestal = m_treeDetailed->Branch("pedestal", &m_pedestal, "pedestal/F");
87  b_calPeakTime = m_treeDetailed->Branch("calPeakTime", &m_calPeakTime, "calPeakTime/F");
88  b_calPeakADC = m_treeDetailed->Branch("calPeakADC", &m_calPeakADC, "calPeakADC/F");
89  b_pulseWidth = m_treeDetailed->Branch("pulseWidth", &m_pulseWidth, "pulseWidth/F");
90 
91 
92  if (!m_MaskedStr.isValid())
93  B2WARNING("No valid SVDFADCMaskedStrip for the requested IoV");
94  if (!m_NoiseCal.isValid())
95  B2WARNING("No valid SVDNoiseCalibration for the requested IoV");
96  if (!m_svdLocalConfig.isValid())
97  B2FATAL("No valid SVDLocalConfigParameters for the requested IoV");
98  if (!m_svdGlobalConfig.isValid())
99  B2FATAL("No valid SVDGlobalConfigParameters for the requested IoV");
100  if (!m_PedestalCal.isValid())
101  B2WARNING("No valid SVDPedestalCalibration for the requested IoV");
102  if (! m_PulseShapeCal.isValid())
103  B2WARNING("No valid SVDPulseShapeCalibrations for the requested IoV");
104  /* if (!m_OccupancyCal.isValid())
105  B2WARNING("No valid SVDOccupancyCalibrations for the requested IoV");
106  if (!m_HotStripsCal.isValid())
107  B2WARNING("No valid SVDHotStripsCalibrations for the requested IoV");
108  */
109 
111  TH1F hOccupancy("occupancy_L@layerL@ladderS@sensor@view",
112  "occupancy in hits/evt in @layer.@ladder.@sensor @view/@side",
113  1500, 0.0, 0.006);
114  hOccupancy.GetXaxis()->SetTitle("strip occupancy ()");
115  m_hOccupancy = new SVDHistograms<TH1F>(hOccupancy);
116 
117  TH2F h2Occupancy_512("occupancy2D_512_L@layerL@ladderS@sensor@view",
118  "occupancy in HITS/EVT in @layer.@ladder.@sensor @view/@side VS cellID",
119  128 * 4, -0.5, 128 * 4 - 0.5, 1500, 0.0, 0.006);
120  h2Occupancy_512.GetYaxis()->SetTitle("strip occupancy (HITS/EVT)");
121  h2Occupancy_512.GetXaxis()->SetTitle("cellID");
122 
123  TH2F h2Occupancy_768("occupancy2D_768_L@layerL@ladderS@sensor@view",
124  "occupancy in HITS/EVT in @layer.@ladder.@sensor @view/@side VS cellID",
125  128 * 6, -0.5, 128 * 6 - 0.5, 1500, 0.0, 0.006);
126  h2Occupancy_768.GetYaxis()->SetTitle("strip occupancy (HITS/EVT)");
127  h2Occupancy_768.GetXaxis()->SetTitle("cellID");
128 
129  m_h2Occupancy = new SVDHistograms<TH2F>(h2Occupancy_768, h2Occupancy_768, h2Occupancy_768, h2Occupancy_512);
130 
132  TH1F hHotstrips("hotstrips_L@layerL@ladderS@sensor@view",
133  "hot strips in @layer.@ladder.@sensor @view/@side",
134  2, -0.5, 1.5);
135  hHotstrips.GetXaxis()->SetTitle("isHotStrips");
136  m_hHotstrips = new SVDHistograms<TH1F>(hHotstrips);
137 
138  //imported from SVDHSfinder module
139  TH1F hHotStrips768("HotStrips768_L@layerL@ladderS@sensor@view", "Hot Strips of @layer.@ladder.@sensor @view/@side side", 768, 0,
140  768);
141  hHotStrips768.GetXaxis()->SetTitle("cellID");
142  TH1F hHotStrips512("HotStrips512_L@layerL@ladderS@sensor@view", "Hot Strips of @layer.@ladder.@sensor @view/@side side", 512, 0,
143  512);
144  hHotStrips512.GetXaxis()->SetTitle("cellID");
145  hm_hot_strips = new SVDHistograms<TH1F>(hHotStrips768, hHotStrips768, hHotStrips768, hHotStrips512);
146 
147  TH2F h2Hotstrips_512("hotstrips2D_512_L@layerL@ladderS@sensor@view",
148  "hot strips in @layer.@ladder.@sensor @view/@side VS cellID",
149  128 * 4, -0.5, 128 * 4 - 0.5, 2, -0.5, 1.5);
150  h2Hotstrips_512.GetYaxis()->SetTitle("isHotStrips");
151  h2Hotstrips_512.GetXaxis()->SetTitle("cellID");
152 
153  TH2F h2Hotstrips_768("hotstrips2D_768_L@layerL@ladderS@sensor@view",
154  "hot strips in @layer.@ladder.@sensor @view/@side VS cellID",
155  128 * 6, -0.5, 128 * 6 - 0.5, 2, -0.5, 1.5);
156  h2Hotstrips_768.GetYaxis()->SetTitle("isHotStrips");
157  h2Hotstrips_768.GetXaxis()->SetTitle("cellID");
158 
159  m_h2Hotstrips = new SVDHistograms<TH2F>(h2Hotstrips_768, h2Hotstrips_768, h2Hotstrips_768, h2Hotstrips_512);
160 
161 
162  //summary plot of the hot strips per sensor
163  m_hHotStripsSummary = new SVDSummaryPlots("hotStripsSummary@view", "Number of HotStrips on @view/@side Side");
164 
166  TH1F hMask("masked_L@layerL@ladderS@sensor@view",
167  "masked strip in @layer.@ladder.@sensor @view/@side",
168  2, -0.5, 1.5);
169  hMask.GetXaxis()->SetTitle("isMasked");
170  m_hMask = new SVDHistograms<TH1F>(hMask);
171 
172  TH2F h2Mask_512("masked2D_512_L@layerL@ladderS@sensor@view",
173  "masked strip in @layer.@ladder.@sensor @view/@side VS cellID",
174  128 * 4, -0.5, 128 * 4 - 0.5, 2, -0.5, 1.5);
175  h2Mask_512.GetYaxis()->SetTitle("isMasked");
176  h2Mask_512.GetXaxis()->SetTitle("cellID");
177 
178  TH2F h2Mask_768("masked2D_768_L@layerL@ladderS@sensor@view",
179  "masked strip in @layer.@ladder.@sensor @view/@side VS cellID",
180  128 * 6, -0.5, 128 * 6 - 0.5, 2, -0.5, 1.5);
181  h2Mask_768.GetYaxis()->SetTitle("isMasked");
182  h2Mask_768.GetXaxis()->SetTitle("cellID");
183 
184  m_h2Mask = new SVDHistograms<TH2F>(h2Mask_768, h2Mask_768, h2Mask_768, h2Mask_512);
185 
187  TH1F hNoise("noiseADC_L@layerL@ladderS@sensor@view",
188  "noise in ADC in @layer.@ladder.@sensor @view/@side",
189  160, -0.5, 19.5);
190  hNoise.GetXaxis()->SetTitle("strip noise (ADC)");
191  m_hNoise = new SVDHistograms<TH1F>(hNoise);
192 
193  TH2F h2Noise_512("noise2D_512_L@layerL@ladderS@sensor@view",
194  "noise in ADC in @layer.@ladder.@sensor @view/@side VS cellID",
195  128 * 4, -0.5, 128 * 4 - 0.5, 80, -0.5, 9.5);
196  h2Noise_512.GetYaxis()->SetTitle("strip noise (ADC)");
197  h2Noise_512.GetXaxis()->SetTitle("cellID");
198 
199  TH2F h2Noise_768("noise2D_768_L@layerL@ladderS@sensor@view",
200  "noise in ADC in @layer.@ladder.@sensor @view/@side VS cellID",
201  128 * 6, -0.5, 128 * 6 - 0.5, 80, -0.5, 9.5);
202  h2Noise_768.GetYaxis()->SetTitle("strip noise (ADC)");
203  h2Noise_768.GetXaxis()->SetTitle("cellID");
204 
205  m_h2Noise = new SVDHistograms<TH2F>(h2Noise_768, h2Noise_768, h2Noise_768, h2Noise_512);
206 
207 
209  TH1F hNoiseEl("noiseEl_L@layerL@ladderS@sensor@view",
210  "noise in e- in @layer.@ladder.@sensor @view/@side",
211  600, -199.5, 1499.5);
212  hNoiseEl.GetXaxis()->SetTitle("strip noise (e-)");
213  m_hNoiseEl = new SVDHistograms<TH1F>(hNoiseEl);
214 
215  TH2F h2NoiseEl_512("noiseEl2D_512_L@layerL@ladderS@sensor@view",
216  "noise in e- in @layer.@ladder.@sensor @view/@side VS cellID",
217  128 * 4, -0.5, 128 * 4 - 0.5, 600, -199.5, 1499.5);
218  h2NoiseEl_512.GetYaxis()->SetTitle("strip noise (e-)");
219  h2NoiseEl_512.GetXaxis()->SetTitle("cellID");
220 
221  TH2F h2NoiseEl_768("noiseEl2D_768_L@layerL@ladderS@sensor@view",
222  "noise in e- in @layer.@ladder.@sensor @view/@side VS cellID",
223  128 * 6, -0.5, 128 * 6 - 0.5, 600, -199.5, 1499.5);
224  h2NoiseEl_768.GetYaxis()->SetTitle("strip noise (e-)");
225  h2NoiseEl_768.GetXaxis()->SetTitle("cellID");
226 
227  m_h2NoiseEl = new SVDHistograms<TH2F>(h2NoiseEl_768, h2NoiseEl_768, h2NoiseEl_768, h2NoiseEl_512);
228 
229 
231  TH1F hPedestal("pedestalADC_L@layerL@ladderS@sensor@view",
232  "pedestal in ADC in @layer.@ladder.@sensor @view/@side",
233  200, -199.5, 599.5);
234  hPedestal.GetXaxis()->SetTitle("strip pedestal (ADC)");
235  m_hPedestal = new SVDHistograms<TH1F>(hPedestal);
236 
237  TH2F h2Pedestal_512("pedestal2D_512_L@layerL@ladderS@sensor@view",
238  "pedestal in ADC in @layer.@ladder.@sensor @view/@side VS cellID",
239  128 * 4, -0.5, 128 * 4 - 0.5, 200, -199.5, 599.5);
240  h2Pedestal_512.GetYaxis()->SetTitle("strip pedestal (ADC)");
241  h2Pedestal_512.GetXaxis()->SetTitle("cellID");
242 
243  TH2F h2Pedestal_768("pedestal2D_768_L@layerL@ladderS@sensor@view",
244  "pedestal in ADC in @layer.@ladder.@sensor @view/@side VS cellID",
245  128 * 6, -0.5, 128 * 6 - 0.5, 200, -199.5, 599.5);
246  h2Pedestal_768.GetYaxis()->SetTitle("strip pedestal (ADC)");
247  h2Pedestal_768.GetXaxis()->SetTitle("cellID");
248 
249  m_h2Pedestal = new SVDHistograms<TH2F>(h2Pedestal_768, h2Pedestal_768, h2Pedestal_768, h2Pedestal_512);
250 
252  TH1F hGain("gainADC_L@layerL@ladderS@sensor@view",
253  "1/gain in @layer.@ladder.@sensor @view/@side",
254  300, -0.5, 499.5);
255  hGain.GetXaxis()->SetTitle("strip 1/gain (e-/ADC)");
256  m_hGain = new SVDHistograms<TH1F>(hGain);
257 
258  TH2F h2Gain_512("gain2D_512_L@layerL@ladderS@sensor@view",
259  "1/gain in @layer.@ladder.@sensor @view/@side VS cellID",
260  128 * 4, -0.5, 128 * 4 - 0.5, 300, -0.5, 499.5);
261  h2Gain_512.GetYaxis()->SetTitle("strip 1/gain (e-/ADC)");
262  h2Gain_512.GetXaxis()->SetTitle("cellID");
263 
264  TH2F h2Gain_768("gain2D_768_L@layerL@ladderS@sensor@view",
265  "1/gain in @layer.@ladder.@sensor @view/@side VS cellID",
266  128 * 6, -0.5, 128 * 6 - 0.5, 300, -0.5, 499.5);
267  h2Gain_768.GetYaxis()->SetTitle("strip 1/gain (e-/ADC)");
268  h2Gain_768.GetXaxis()->SetTitle("cellID");
269 
270  m_h2Gain = new SVDHistograms<TH2F>(h2Gain_768, h2Gain_768, h2Gain_768, h2Gain_512);
271 
272  // PEAKTIME (ns)
273  TH1F hCalPeakTime("calPeakTime_L@layerL@ladderS@sensor@view",
274  "calPeakTime in @layer.@ladder.@sensor @view/@side",
275  255, -0.5, 254.5);
276  hCalPeakTime.GetXaxis()->SetTitle("strip calPeakTime (ns)");
277  m_hCalPeakTime = new SVDHistograms<TH1F>(hCalPeakTime);
278 
279  TH2F h2CalPeakTime_512("calPeakTime2D_512_L@layerL@ladderS@sensor@view",
280  "calPeakTime in @layer.@ladder.@sensor @view/@side VS cellID",
281  128 * 4, -0.5, 128 * 4 - 0.5, 255, -0.5, 254.5);
282  h2CalPeakTime_512.GetYaxis()->SetTitle("strip calPeakTime (ns)");
283  h2CalPeakTime_512.GetXaxis()->SetTitle("cellID");
284 
285  TH2F h2CalPeakTime_768("calPeakTime2D_768_L@layerL@ladderS@sensor@view",
286  "calPeakTime in @layer.@ladder.@sensor @view/@side VS cellID",
287  128 * 6, -0.5, 128 * 6 - 0.5, 255, -0.5, 254.5);
288  h2CalPeakTime_768.GetYaxis()->SetTitle("strip calPeakTime (ns)");
289  h2CalPeakTime_768.GetXaxis()->SetTitle("cellID");
290 
291  m_h2CalPeakTime = new SVDHistograms<TH2F>(h2CalPeakTime_768, h2CalPeakTime_768, h2CalPeakTime_768, h2CalPeakTime_512);
292 
293  //CALPEAK ADC
294  TH1F hCalPeakADC("calPeakADC_L@layerL@ladderS@sensor@view",
295  "calPeakADC in @layer.@ladder.@sensor @view/@side",
296  80, 44.5, 124.5);
297  hCalPeakADC.GetXaxis()->SetTitle("strip calPeakADC (ADC)");
298  m_hCalPeakADC = new SVDHistograms<TH1F>(hCalPeakADC);
299 
300  TH2F h2CalPeakADC_512("calPeakADC2D_512_L@layerL@ladderS@sensor@view",
301  "calPeakADC in @layer.@ladder.@sensor @view/@side VS cellID",
302  128 * 4, -0.5, 128 * 4 - 0.5, 80, 44.5, 124.5);
303  h2CalPeakADC_512.GetYaxis()->SetTitle("strip calPeakADC (ADC)");
304  h2CalPeakADC_512.GetXaxis()->SetTitle("cellID");
305 
306  TH2F h2CalPeakADC_768("calPeakADC2D_768_L@layerL@ladderS@sensor@view",
307  "calPeakADC in @layer.@ladder.@sensor @view/@side VS cellID",
308  128 * 6, -0.5, 128 * 6 - 0.5, 80, 44.5, 124.5);
309  h2CalPeakADC_768.GetYaxis()->SetTitle("strip calPeakADC (ADC)");
310  h2CalPeakADC_768.GetXaxis()->SetTitle("cellID");
311 
312  m_h2CalPeakADC = new SVDHistograms<TH2F>(h2CalPeakADC_768, h2CalPeakADC_768, h2CalPeakADC_768, h2CalPeakADC_512);
313 
314  // PULSE WIDTH (ns)
315  TH1F hPulseWidth("pulseWidth_L@layerL@ladderS@sensor@view",
316  "pulseWidth in @layer.@ladder.@sensor @view/@side",
317  255, -0.5, 254.5);
318  hPulseWidth.GetXaxis()->SetTitle("strip pulseWidth (ns)");
319  m_hPulseWidth = new SVDHistograms<TH1F>(hPulseWidth);
320 
321  TH2F h2PulseWidth_512("pulseWidth2D_512_L@layerL@ladderS@sensor@view",
322  "pulseWidth in @layer.@ladder.@sensor @view/@side VS cellID",
323  128 * 4, -0.5, 128 * 4 - 0.5, 255, -0.5, 254.5);
324  h2PulseWidth_512.GetYaxis()->SetTitle("strip pulseWidth (ns)");
325  h2PulseWidth_512.GetXaxis()->SetTitle("cellID");
326 
327  TH2F h2PulseWidth_768("pulseWidth2D_768_L@layerL@ladderS@sensor@view",
328  "pulseWidth in @layer.@ladder.@sensor @view/@side VS cellID",
329  128 * 6, -0.5, 128 * 6 - 0.5, 255, -0.5, 254.5);
330  h2PulseWidth_768.GetYaxis()->SetTitle("strip pulseWidth (ns)");
331  h2PulseWidth_768.GetXaxis()->SetTitle("cellID");
332 
333  m_h2PulseWidth = new SVDHistograms<TH2F>(h2PulseWidth_768, h2PulseWidth_768, h2PulseWidth_768, h2PulseWidth_512);
334 
335 }
336 
338 {
339 
341  m_exp = meta->getExperiment();
342  m_run = meta->getRun();
343 
344  m_hv = m_svdGlobalConfig->getHV();
345  m_svdLocalConfig->getCalibDate().copy(m_date, 10);
346  m_date[10] = '\0';
347 
348  //call for a geometry instance
350  std::set<Belle2::VxdID> svdLayers = aGeometry.getLayers(VXD::SensorInfoBase::SVD);
351  std::set<Belle2::VxdID>::iterator itSvdLayers = svdLayers.begin();
352 
353  while ((itSvdLayers != svdLayers.end()) && (itSvdLayers->getLayerNumber() != 7)) { //loop on Layers
354 
355  std::set<Belle2::VxdID> svdLadders = aGeometry.getLadders(*itSvdLayers);
356  std::set<Belle2::VxdID>::iterator itSvdLadders = svdLadders.begin();
357 
358  while (itSvdLadders != svdLadders.end()) { //loop on Ladders
359 
360  std::set<Belle2::VxdID> svdSensors = aGeometry.getSensors(*itSvdLadders);
361  std::set<Belle2::VxdID>::iterator itSvdSensors = svdSensors.begin();
362  B2DEBUG(1, " svd sensor info " << * (svdSensors.begin()));
363 
364  while (itSvdSensors != svdSensors.end()) { //loop on sensors
365  B2DEBUG(1, " svd sensor info " << *itSvdSensors);
366 
367  int layer = itSvdSensors->getLayerNumber();
368  int ladder = itSvdSensors->getLadderNumber();
369  int sensor = itSvdSensors->getSensorNumber();
370  Belle2::VxdID theVxdID(layer, ladder, sensor);
371  const SVD::SensorInfo* currentSensorInfo = dynamic_cast<const SVD::SensorInfo*>(&VXD::GeoCache::get(theVxdID));
372  m_layer = layer;
373  m_ladder = ladder;
374  m_sensor = sensor;
375 
376  for (m_side = 0; m_side < 2; m_side++) {
377 
378  int Ncells = currentSensorInfo->getUCells();
379  if (m_side == 0)
380  Ncells = currentSensorInfo->getVCells();
381 
382  for (m_strip = 0; m_strip < Ncells; m_strip++) {
383  m_occupancy = -1;
384  /* if (m_OccupancyCal.isValid()) {
385  m_occupancy = m_OccupancyCal.getOccupancy(theVxdID, m_side, m_strip);
386  }*/
387  m_hOccupancy->fill(theVxdID, m_side, m_occupancy);
389 
390 
391  m_hotstrips = -1;
392  /* if (m_HotStripsCal.isValid())
393  m_hotstrips = m_HotStripsCal.isHot(theVxdID, m_side, m_strip);*/
394 
395  //aux histo for hotStripSummary table
396  hm_hot_strips->getHistogram(*itSvdSensors, m_side)->SetBinContent(m_strip + 1, m_hotstrips);
397  m_hHotstrips->fill(theVxdID, m_side, m_hotstrips);
399 
400  m_mask = -1;
401  if (m_MaskedStr.isValid())
402  m_mask = m_MaskedStr.isMasked(theVxdID, m_side, m_strip);
403  m_hMask->fill(theVxdID, m_side, m_mask);
404  m_h2Mask->fill(theVxdID, m_side, m_strip, m_mask);
405 
406  m_noise = -1;
407  m_noiseEl = -1;
408  if (m_NoiseCal.isValid()) {
409  m_noise = m_NoiseCal.getNoise(theVxdID, m_side, m_strip);
411  }
412  m_hNoise->fill(theVxdID, m_side, m_noise);
413  m_h2Noise->fill(theVxdID, m_side, m_strip, m_noise);
414  m_hNoiseEl->fill(theVxdID, m_side, m_noiseEl);
415  m_h2NoiseEl->fill(theVxdID, m_side, m_strip, m_noiseEl);
416 
417  m_pedestal = -1;
418  if (m_PedestalCal.isValid())
420  m_hPedestal->fill(theVxdID, m_side, m_pedestal);
421  m_h2Pedestal->fill(theVxdID, m_side, m_strip, m_pedestal);
422 
423  m_gain = -1;
424  if (m_PulseShapeCal.isValid()) {
425  m_gain = m_PulseShapeCal.getChargeFromADC(theVxdID, m_side, m_strip, 1/*ADC*/);
426  m_calPeakADC = 22500. / m_PulseShapeCal.getChargeFromADC(theVxdID, m_side, m_strip, 1/*ADC*/);
429  }
430  m_hGain->fill(theVxdID, m_side, m_gain);
431  m_h2Gain->fill(theVxdID, m_side, m_strip, m_gain);
434  m_hCalPeakADC->fill(theVxdID, m_side, m_calPeakADC);
436  m_hPulseWidth->fill(theVxdID, m_side, m_pulseWidth);
438 
439  m_treeDetailed->Fill();
440 
441  }
442  }
443  ++itSvdSensors;
444  }
445  ++itSvdLadders;
446  }
447  ++itSvdLayers;
448  }
449 
450  B2INFO("now computing Mean and RMS of local calibration constants");
451 
452  //compute averages and RMS
453 
454  itSvdLayers = svdLayers.begin();
455 
456  while ((itSvdLayers != svdLayers.end()) && (itSvdLayers->getLayerNumber() != 7)) { //loop on Layers
457 
458  std::set<Belle2::VxdID> svdLadders = aGeometry.getLadders(*itSvdLayers);
459  std::set<Belle2::VxdID>::iterator itSvdLadders = svdLadders.begin();
460 
461  while (itSvdLadders != svdLadders.end()) { //loop on Ladders
462 
463  std::set<Belle2::VxdID> svdSensors = aGeometry.getSensors(*itSvdLadders);
464  std::set<Belle2::VxdID>::iterator itSvdSensors = svdSensors.begin();
465  B2DEBUG(1, " svd sensor info " << * (svdSensors.begin()));
466 
467  while (itSvdSensors != svdSensors.end()) { //loop on sensors
468  B2DEBUG(1, " svd sensor info " << *itSvdSensors);
469 
470  m_layer = itSvdSensors->getLayerNumber();
471  m_ladder = itSvdSensors->getLadderNumber();
472  m_sensor = itSvdSensors->getSensorNumber();
474 
475 
476  for (m_side = 0; m_side < 2; m_side++) {
477  m_maskAVE = (m_hMask->getHistogram(theVxdID, m_side))->GetMean();
478  m_hotstripsAVE = (m_hHotstrips->getHistogram(theVxdID, m_side))->GetMean();
479  m_pedestalAVE = (m_hPedestal->getHistogram(theVxdID, m_side))->GetMean();
480  m_pedestalRMS = (m_hPedestal->getHistogram(theVxdID, m_side))->GetRMS();
481  m_noiseAVE = (m_hNoise->getHistogram(theVxdID, m_side))->GetMean();
482  m_noiseRMS = (m_hNoise->getHistogram(theVxdID, m_side))->GetRMS();
483  m_noiseElAVE = (m_hNoiseEl->getHistogram(theVxdID, m_side))->GetMean();
484  m_noiseElRMS = (m_hNoiseEl->getHistogram(theVxdID, m_side))->GetRMS();
485  m_occupancyAVE = (m_hOccupancy->getHistogram(theVxdID, m_side))->GetMean();
486  m_occupancyRMS = (m_hOccupancy->getHistogram(theVxdID, m_side))->GetRMS();
487  m_gainAVE = (m_hGain->getHistogram(theVxdID, m_side))->GetMean();
488  m_gainRMS = (m_hGain->getHistogram(theVxdID, m_side))->GetRMS();
489  m_calPeakTimeAVE = (m_hCalPeakTime->getHistogram(theVxdID, m_side))->GetMean();
490  m_calPeakTimeRMS = (m_hCalPeakTime->getHistogram(theVxdID, m_side))->GetRMS();
491  m_calPeakADCAVE = (m_hCalPeakADC->getHistogram(theVxdID, m_side))->GetMean();
492  m_calPeakADCRMS = (m_hCalPeakADC->getHistogram(theVxdID, m_side))->GetRMS();
493  m_pulseWidthAVE = (m_hPulseWidth->getHistogram(theVxdID, m_side))->GetMean();
494  m_pulseWidthRMS = (m_hPulseWidth->getHistogram(theVxdID, m_side))->GetRMS();
495 
496 
497  // for (int s = 0; s < hm_hot_strips->getHistogram(*itSvdSensors, m_side)->GetEntries(); s++)
498  // m_hHotStripsSummary->fill(*itSvdSensors, m_side, 1);
499 
500  m_tree->Fill();
501 
502  }
503  ++itSvdSensors;
504  }
505  ++itSvdLadders;
506  }
507  ++itSvdLayers;
508  }
509 
510 
511 }
512 
514 {
515  B2RESULT("******************************************");
516  B2RESULT("** UNIQUE IDs of calibration DB objects **");
517  B2RESULT("");
518 
519  /* if (m_OccupancyCal.isValid())
520  B2RESULT(" - SVDOccupancyCalibrations:" << m_OccupancyCal.getUniqueID());
521  else
522  B2WARNING("No valid SVDOccupancyCalibrations for the requested IoV");
523 
524  if (m_HotStripsCal.isValid())
525  B2RESULT(" - SVDHotStripsCalibrations:" << m_HotStripsCal.getUniqueID());
526  else
527  B2WARNING("No valid SVDHotStripsCalibrations for the requested IoV");*/
528 
529 
530  if (m_MaskedStr.isValid())
531  B2RESULT(" - SVDFADCMaskedStrips:" << m_MaskedStr.getUniqueID());
532  else
533  B2WARNING("No valid SVDFADCMaskedStrips for the requested IoV");
534 
535  if (m_NoiseCal.isValid())
536  B2RESULT(" - SVDNoiseCalibrations:" << m_NoiseCal.getUniqueID());
537  else
538  B2WARNING("No valid SVDNoiseCalibrations for the requested IoV");
539 
540  if (m_PedestalCal.isValid())
541  B2RESULT(" - SVDPedestalCalibrations:" << m_PedestalCal.getUniqueID());
542  else
543  B2WARNING("No valid SVDPedestalCalibrations for the requested IoV");
544 
545  if (m_PulseShapeCal.isValid())
546  B2RESULT(" - SVDPulseShapeCalibrations:" << m_PulseShapeCal.getUniqueID());
547  else
548  B2WARNING("No valid SVDPulseShapeCalibrations for the requested IoV");
549  //}
550 
551  //void SVDLocalCalibrationsMonitorModule::terminate()
552  //{
553 
554  if (m_rootFilePtr != nullptr) {
555 
556  m_rootFilePtr->cd();
557 
558  //write the tree
559  m_treeDetailed->Write();
560  m_tree->Write();
561 
562  m_rootFilePtr->mkdir("hotstrips");
563  m_rootFilePtr->mkdir("masked_strips");
564  m_rootFilePtr->mkdir("pedestal_ADCunits");
565  m_rootFilePtr->mkdir("noise_ADCunits");
566  m_rootFilePtr->mkdir("occupancy");
567  m_rootFilePtr->mkdir("noise_electronsCharge");
568  m_rootFilePtr->mkdir("gain_electronsCharge");
569  m_rootFilePtr->mkdir("calPeakTime");
570  m_rootFilePtr->mkdir("calPeakADC");
571  m_rootFilePtr->mkdir("pulseWidth");
572 
573 
575 
576  for (auto layer : geoCache.getLayers(VXD::SensorInfoBase::SVD))
577  for (auto ladder : geoCache.getLadders(layer))
578  for (Belle2::VxdID sensor : geoCache.getSensors(ladder))
579  for (int view = SVDHistograms<TH1F>::VIndex ; view < SVDHistograms<TH1F>::UIndex + 1; view++) {
580 
581  //writing the histogram list for the noises in ADC units
582 
583  m_rootFilePtr->cd("occupancy");
584  (m_hOccupancy->getHistogram(sensor, view))->Write();
585  (m_h2Occupancy->getHistogram(sensor, view))->Write();
586 
587  //writing the histogram list for the hotstrips
588  m_rootFilePtr->cd("hotstrips");
589  //------imported from SVDHSfinder module
590  hm_hot_strips->getHistogram(sensor, view)->SetLineColor(kBlack);
591  hm_hot_strips->getHistogram(sensor, view)->SetMarkerColor(kBlack);
592  hm_hot_strips->getHistogram(sensor, view)->SetFillStyle(3001);
593  hm_hot_strips->getHistogram(sensor, view)->SetFillColor(kBlack);
594  hm_hot_strips->getHistogram(sensor, view)->Write();
595 
596  //--------------------
597  (m_hHotstrips->getHistogram(sensor, view))->Write();
598  (m_h2Hotstrips->getHistogram(sensor, view))->Write();
599 
600 
601 
602  //writing the histogram list for the masks in ADC units
603  m_rootFilePtr->cd("masked_strips");
604  (m_hMask->getHistogram(sensor, view))->Write();
605  (m_h2Mask->getHistogram(sensor, view))->Write();
606 
607  //writing the histogram list for the pedestals in ADC units
608  m_rootFilePtr->cd("pedestal_ADCunits");
609  (m_hPedestal->getHistogram(sensor, view))->Write();
610  (m_h2Pedestal->getHistogram(sensor, view))->Write();
611 
612  //writing the histogram list for the noises in ADC units
613  m_rootFilePtr->cd("noise_ADCunits");
614  (m_hNoise->getHistogram(sensor, view))->Write();
615  (m_h2Noise->getHistogram(sensor, view))->Write();
616 
617  //writing the histogram list for the noises in electron charge
618  m_rootFilePtr->cd("noise_electronsCharge");
619  (m_hNoiseEl->getHistogram(sensor, view))->Write();
620  (m_h2NoiseEl->getHistogram(sensor, view))->Write();
621 
622  //writing the histogram list for the gains in electron charge
623  m_rootFilePtr->cd("gain_electronsCharge");
624  (m_hGain->getHistogram(sensor, view))->Write();
625  (m_h2Gain->getHistogram(sensor, view))->Write();
626 
627  //writing the histogram list for the peak times in ns
628  m_rootFilePtr->cd("calPeakTime");
629  (m_hCalPeakTime->getHistogram(sensor, view))->Write();
630  (m_h2CalPeakTime->getHistogram(sensor, view))->Write();
631 
632  //writing the histogram list for the peak in ADC
633  m_rootFilePtr->cd("calPeakADC");
634  (m_hCalPeakADC->getHistogram(sensor, view))->Write();
635  (m_h2CalPeakADC->getHistogram(sensor, view))->Write();
636 
637  //writing the histogram list for the pulse widths in ns
638  m_rootFilePtr->cd("pulseWidth");
639  (m_hPulseWidth->getHistogram(sensor, view))->Write();
640  (m_h2PulseWidth->getHistogram(sensor, view))->Write();
641 
642  }
643  m_rootFilePtr->mkdir("expert");
644 
645  m_rootFilePtr->cd("expert");
646  m_h2Noise->Write("h2Noise");
647  m_h2Occupancy->Write("h2Occupancy");
648  m_h2PulseWidth->Write("h2PulseShape");
649  m_h2Pedestal->Write("h2Pedestal");
650  m_h2Gain->Write("h2Gain");
651  m_h2CalPeakADC->Write("h2CalPeakADC");
652  m_h2CalPeakTime->Write("h2CalPeakTime");
653 
654  m_rootFilePtr->Close();
655  B2RESULT("The rootfile containing the list of histograms has been filled and closed [Local].");
656 
657 
658  }
659 }
Base class for Modules.
Definition: Module.h:72
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:214
float isMasked(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
This is the method for getting the comprehensive list of masked strips at FADC level.
TString getUniqueID()
returns the unique ID of the payload
bool isValid()
returns true if the m_aDBObtPtr is valid in the requested IoV
void fill(const VxdID &vxdID, int view, Types ... args)
fill the histogram for
Definition: SVDHistograms.h:77
H * getHistogram(const VxdID &vxdID, int view)
get a reference to the histogram for
Definition: SVDHistograms.h:56
SVDHistograms< TH2F > * m_h2Pedestal
pedestal (ADC) VS strip 2D histo
SVDHistograms< TH2F > * m_h2CalPeakADC
calPeakADC (ns) VS strip 2D histo
SVDFADCMaskedStrips m_MaskedStr
FADC masked strip payload.
SVDHistograms< TH2F > * m_h2Noise
noise (ADC) VS strip 2D histo
SVDLocalCalibrationsMonitorModule()
Constructor: Sets the description, the properties and the parameters of the module.
TTree * m_tree
pointer at tree containing the mean and RMS of calibration constants
TBranch * b_date
date of the noise local run in yyyy-mm-dd format
SVDHistograms< TH2F > * m_h2CalPeakTime
calPeakTime (ns) VS strip 2D histo
TTree * m_treeDetailed
pointer at tree containing the calibration constants of each strip
SVDHistograms< TH1F > * hm_hot_strips
hot strips per sensor
virtual void event() override
fill trees and histograms
SVDHistograms< TH1F > * m_hGain
gain (e-/ADC) histo
virtual void endRun() override
print the payloads uniqueID and write trees and histograms to the rootfile
SVDSummaryPlots * m_hHotStripsSummary
hot strip summary histo
SVDHistograms< TH1F > * m_hCalPeakTime
calPeakTime (ns) histo
SVDHistograms< TH2F > * m_h2Mask
mask VS strip 2D histo
SVDPulseShapeCalibrations m_PulseShapeCal
pulse shape payload
SVDHistograms< TH1F > * m_hNoiseEl
noise in e- histo
virtual void beginRun() override
initialize the TTrees and check validities of payloads
SVDHistograms< TH2F > * m_h2NoiseEl
noise in e- VS strip 2D histo
DBObjPtr< SVDLocalConfigParameters > m_svdLocalConfig
SVD Local Configuration payload.
char m_date[11]
date of the noise local run in yyyy-mm-dd format
SVDHistograms< TH1F > * m_hHotstrips
hot strips histo
SVDHistograms< TH1F > * m_hNoise
noise (ADC) histo
DBObjPtr< SVDGlobalConfigParameters > m_svdGlobalConfig
SVD Global Configuration payload.
SVDHistograms< TH1F > * m_hPulseWidth
calPeakTime (ns) histo
SVDHistograms< TH1F > * m_hCalPeakADC
calPeakADC (ns) histo
SVDHistograms< TH2F > * m_h2PulseWidth
calPeakTime (ns) VS strip 2D histo
SVDHistograms< TH1F > * m_hPedestal
pedestal (ADC) histo
SVDHistograms< TH2F > * m_h2Occupancy
occupancy (hits/evt) VS strip 2D histo
TBranch * b_occupancyAVE
sensor occupancy average (ADC)
TFile * m_rootFilePtr
pointer at root file used for storing histograms
SVDPedestalCalibrations m_PedestalCal
pedestal payload
SVDHistograms< TH2F > * m_h2Gain
gain (e-/ADC) VS strip 2D histo
SVDHistograms< TH2F > * m_h2Hotstrips
hotstrips VS strip 2D histo
SVDHistograms< TH1F > * m_hOccupancy
occupancy (hits/evt) histo
float getNoiseInElectrons(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
This method provides the correct noise conversion into electrons, taking into account that the noise ...
float getNoise(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
This is the method for getting the noise.
TString getUniqueID()
returns the unique ID of the payload
bool isValid()
returns true if the m_aDBObtPtr is valid in the requested IoV
float getPedestal(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
This is the method for getting the pedestal.
TString getUniqueID()
returns the unique ID of the payload
bool isValid()
returns true if the m_aDBObtPtr is valid in the requested IoV
double getChargeFromADC(const Belle2::VxdID &sensorID, const bool &isU, const unsigned short &strip, const double &pulseADC) const
Return the charge (number of electrons/holes) collected on a specific strip, given the number of ADC ...
float getPeakTime(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
Return the peaking time of the strip.
TString getUniqueID()
returns the unique ID of the payload
float getWidth(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
Return the width of the pulse shape for a given strip.
bool isValid()
returns true if the m_aDBObtPtr is valid in the requested IoV
class to summarize SVD quantities per sensor and side
Specific implementation of SensorInfo for SVD Sensors which provides additional sensor specific infor...
Definition: SensorInfo.h:25
Type-safe access to single objects in the data store.
Definition: StoreObjPtr.h:96
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:39
const std::set< Belle2::VxdID > getLayers(SensorInfoBase::SensorType sensortype=SensorInfoBase::VXD)
Return a set of all known Layers.
Definition: GeoCache.cc:176
const std::set< Belle2::VxdID > & getSensors(Belle2::VxdID ladder) const
Return a set of all sensor IDs belonging to a given ladder.
Definition: GeoCache.cc:204
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:214
static const SensorInfoBase & get(Belle2::VxdID id)
Return a reference to the SensorInfo of a given SensorID.
Definition: GeoCache.h:139
const std::set< Belle2::VxdID > & getLadders(Belle2::VxdID layer) const
Return a set of all ladder IDs belonging to a given layer.
Definition: GeoCache.cc:193
int getVCells() const
Return number of pixel/strips in v direction.
int getUCells() const
Return number of pixel/strips in u direction.
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:33
REG_MODULE(arichBtest)
Register the Module.
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:560
Abstract base class for different kinds of events.