CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
HcalDigisValidation.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: HcalDigisValidation
4 // Class: HcalDigisValidation
5 //
13 //
14 // Original Author: Ali Fahim,22 R-013,+41227672649,
15 // Created: Wed Mar 23 11:42:34 CET 2011
16 //
17 //
18 
22 
24  subdet_ = iConfig.getUntrackedParameter<std::string>("subdetector", "all");
25  outputFile_ = iConfig.getUntrackedParameter<std::string>("outputFile", "");
26  // inputLabel_ = iConfig.getParameter<std::string > ("digiLabel");
27  inputTag_ = iConfig.getParameter<edm::InputTag>("digiTag");
28  QIE10inputTag_ = iConfig.getParameter<edm::InputTag>("QIE10digiTag");
29  QIE11inputTag_ = iConfig.getParameter<edm::InputTag>("QIE11digiTag");
30  emulTPsTag_ = iConfig.getParameter<edm::InputTag>("emulTPs");
31  dataTPsTag_ = iConfig.getParameter<edm::InputTag>("dataTPs");
32  mc_ = iConfig.getUntrackedParameter<std::string>("mc", "no");
33  mode_ = iConfig.getUntrackedParameter<std::string>("mode", "multi");
34  dirName_ = iConfig.getUntrackedParameter<std::string>("dirName", "HcalDigisV/HcalDigiTask");
35  testNumber_ = iConfig.getParameter<bool>("TestNumber");
36  hep17_ = iConfig.getParameter<bool>("hep17");
37  HEPhase1_ = iConfig.getParameter<bool>("HEPhase1");
38  HBPhase1_ = iConfig.getParameter<bool>("HBPhase1");
39  Plot_TP_ver_ = iConfig.getParameter<bool>("Plot_TP_ver");
40 
41  // register for data access
42  if (iConfig.exists("simHits")) {
43  tok_mc_ = consumes<edm::PCaloHitContainer>(iConfig.getUntrackedParameter<edm::InputTag>("simHits"));
44  }
45  tok_hbhe_ = consumes<HBHEDigiCollection>(inputTag_);
46  tok_ho_ = consumes<HODigiCollection>(inputTag_);
47  tok_hf_ = consumes<HFDigiCollection>(inputTag_);
48  tok_emulTPs_ = consumes<HcalTrigPrimDigiCollection>(emulTPsTag_);
49  if (dataTPsTag_ == edm::InputTag(""))
50  skipDataTPs = true;
51  else {
52  skipDataTPs = false;
53  tok_dataTPs_ = consumes<HcalTrigPrimDigiCollection>(dataTPsTag_);
54  }
55 
56  tok_qie10_hf_ = consumes<QIE10DigiCollection>(QIE10inputTag_);
57  tok_qie11_hbhe_ = consumes<QIE11DigiCollection>(QIE11inputTag_);
58 
59  tok_HRNDC_ = esConsumes<HcalDDDRecConstants, HcalRecNumberingRecord, edm::Transition::BeginRun>();
60  tok_Geom_ = esConsumes<CaloGeometry, CaloGeometryRecord, edm::Transition::BeginRun>();
61  tok_Decoder_ = esConsumes<CaloTPGTranscoder, CaloTPGRecord>();
62  tok_TPGeom_ = esConsumes<HcalTrigTowerGeometry, CaloGeometryRecord>();
63  tok_Topo_ = esConsumes<HcalTopology, HcalRecNumberingRecord>();
64  tok_Cond_ = esConsumes<HcalDbService, HcalDbRecord>();
65 
66  nevent1 = 0;
67  nevent2 = 0;
68  nevent3 = 0;
69  nevent4 = 0;
70  nevtot = 0;
71 
72  msm_ = new std::map<std::string, MonitorElement*>();
73 
74  if (!outputFile_.empty())
75  edm::LogVerbatim("OutputInfo") << " Hcal Digi Task histograms will be saved to '" << outputFile_.c_str() << "'";
76  else
77  edm::LogVerbatim("OutputInfo") << " Hcal Digi Task histograms will NOT be saved";
78 }
79 
81 
83  hcons_ = &es.getData(tok_HRNDC_);
84 
85  maxDepth_[1] = hcons_->getMaxDepth(0); // HB
86  maxDepth_[2] = hcons_->getMaxDepth(1); // HE
87  maxDepth_[3] = hcons_->getMaxDepth(3); // HO
88  maxDepth_[4] = hcons_->getMaxDepth(2); // HF
89  maxDepth_[0] = (maxDepth_[1] > maxDepth_[2] ? maxDepth_[1] : maxDepth_[2]);
90  maxDepth_[0] = (maxDepth_[0] > maxDepth_[3] ? maxDepth_[0] : maxDepth_[3]);
91  maxDepth_[0] = (maxDepth_[0] > maxDepth_[4] ? maxDepth_[0] : maxDepth_[4]); // any of HB/HE/HO/HF
92 
93  const CaloGeometry* geo = &es.getData(tok_Geom_);
94  const HcalGeometry* gHB = static_cast<const HcalGeometry*>(geo->getSubdetectorGeometry(DetId::Hcal, HcalBarrel));
95  const HcalGeometry* gHE = static_cast<const HcalGeometry*>(geo->getSubdetectorGeometry(DetId::Hcal, HcalEndcap));
96  const HcalGeometry* gHO = static_cast<const HcalGeometry*>(geo->getSubdetectorGeometry(DetId::Hcal, HcalOuter));
97  const HcalGeometry* gHF = static_cast<const HcalGeometry*>(geo->getSubdetectorGeometry(DetId::Hcal, HcalForward));
98 
99  nChannels_[1] = gHB->getHxSize(1);
100  nChannels_[2] = gHE->getHxSize(2);
101  nChannels_[3] = gHO->getHxSize(3);
102  nChannels_[4] = gHF->getHxSize(4);
103 
104  nChannels_[0] = nChannels_[1] + nChannels_[2] + nChannels_[3] + nChannels_[4];
105 }
106 
109 
110  // book
111  book1D(ib, "nevtot", 1, 0, 1);
112  int bnoise = 0;
113  int bmc = 0;
114  if (subdet_ == "noise")
115  bnoise = 1;
116  if (mc_ == "yes")
117  bmc = 1;
118  if (subdet_ == "noise" || subdet_ == "all") {
119  booking(ib, "HB", bnoise, bmc);
120  booking(ib, "HO", bnoise, bmc);
121  booking(ib, "HF", bnoise, bmc);
122  booking(ib, "HE", bnoise, bmc);
123  } else {
124  booking(ib, subdet_, 0, bmc);
125  }
126 
127  if (skipDataTPs)
128  return;
129 
130  HistLim tp_hl_et(260, -10, 250);
131  HistLim tp_hl_ntp(640, -20, 3180);
132  HistLim tp_hl_ntp_sub(404, -20, 2000);
133  HistLim tp_hl_ieta(85, -42.5, 42.5);
134  HistLim tp_hl_iphi(74, -0.5, 73.5);
135 
136  book1D(ib, "HcalDigiTask_tp_et", tp_hl_et);
137  book1D(ib, "HcalDigiTask_tp_et_HB", tp_hl_et);
138  book1D(ib, "HcalDigiTask_tp_et_HE", tp_hl_et);
139  book1D(ib, "HcalDigiTask_tp_et_HF", tp_hl_et);
140  book1D(ib, "HcalDigiTask_tp_ntp", tp_hl_ntp);
141  book1D(ib, "HcalDigiTask_tp_ntp_HB", tp_hl_ntp_sub);
142  book1D(ib, "HcalDigiTask_tp_ntp_HE", tp_hl_ntp_sub);
143  book1D(ib, "HcalDigiTask_tp_ntp_HF", tp_hl_ntp_sub);
144  book1D(ib, "HcalDigiTask_tp_ntp_ieta", tp_hl_ieta);
145  book1D(ib, "HcalDigiTask_tp_ntp_iphi", tp_hl_iphi);
146  book1D(ib, "HcalDigiTask_tp_ntp_10_ieta", tp_hl_ieta);
147  book2D(ib, "HcalDigiTask_tp_et_ieta", tp_hl_ieta, tp_hl_et);
148  book2D(ib, "HcalDigiTask_tp_ieta_iphi", tp_hl_ieta, tp_hl_iphi);
149  bookPf(ib, "HcalDigiTask_tp_ave_et_ieta", tp_hl_ieta, tp_hl_et, " ");
150  if (Plot_TP_ver_) {
151  book1D(ib, "HcalDigiTask_tp_et_v0", tp_hl_et);
152  book1D(ib, "HcalDigiTask_tp_et_v1", tp_hl_et);
153  book1D(ib, "HcalDigiTask_tp_et_HF_v0", tp_hl_et);
154  book1D(ib, "HcalDigiTask_tp_et_HF_v1", tp_hl_et);
155  book1D(ib, "HcalDigiTask_tp_ntp_v0", tp_hl_ntp);
156  book1D(ib, "HcalDigiTask_tp_ntp_v1", tp_hl_ntp);
157  book1D(ib, "HcalDigiTask_tp_ntp_HF_v0", tp_hl_ntp_sub);
158  book1D(ib, "HcalDigiTask_tp_ntp_HF_v1", tp_hl_ntp_sub);
159  book1D(ib, "HcalDigiTask_tp_ntp_ieta_v0", tp_hl_ieta);
160  book1D(ib, "HcalDigiTask_tp_ntp_ieta_v1", tp_hl_ieta);
161  book1D(ib, "HcalDigiTask_tp_ntp_iphi_v0", tp_hl_iphi);
162  book1D(ib, "HcalDigiTask_tp_ntp_iphi_v1", tp_hl_iphi);
163  book1D(ib, "HcalDigiTask_tp_ntp_10_ieta_v0", tp_hl_ieta);
164  book1D(ib, "HcalDigiTask_tp_ntp_10_ieta_v1", tp_hl_ieta);
165  book2D(ib, "HcalDigiTask_tp_et_ieta_v0", tp_hl_ieta, tp_hl_et);
166  book2D(ib, "HcalDigiTask_tp_et_ieta_v1", tp_hl_ieta, tp_hl_et);
167  book2D(ib, "HcalDigiTask_tp_ieta_iphi_v0", tp_hl_ieta, tp_hl_iphi);
168  book2D(ib, "HcalDigiTask_tp_ieta_iphi_v1", tp_hl_ieta, tp_hl_iphi);
169  bookPf(ib, "HcalDigiTask_tp_ave_et_ieta_v0", tp_hl_ieta, tp_hl_et, " ");
170  bookPf(ib, "HcalDigiTask_tp_ave_et_ieta_v1", tp_hl_ieta, tp_hl_et, " ");
171  }
172 }
173 
174 void HcalDigisValidation::booking(DQMStore::IBooker& ib, const std::string bsubdet, int bnoise, int bmc) {
175  // Adjust/Optimize binning (JR Dittmann, 16-JUL-2015)
176 
177  HistLim Ndigis(2600, 0., 2600.);
178  HistLim sime(200, 0., 1.0);
179 
180  HistLim digiAmp(360, -100., 7100.);
181  HistLim digiAmpWide(2410, -3000., 720000.); //300 fC binning
182  HistLim ratio(2000, -100., 3900.);
183  HistLim sumAmp(100, -500., 1500.);
184 
185  HistLim nbin(11, -0.5, 10.5);
186 
187  HistLim pedestal(80, -1.0, 15.);
188  HistLim pedestalfC(400, -10., 30.);
189 
190  HistLim frac(80, -0.20, 1.40);
191 
192  HistLim pedLim(80, 0., 8.);
193  HistLim pedWidthLim(100, 0., 2.);
194 
195  HistLim gainLim(120, 0., 0.6);
196  HistLim gainWidthLim(160, 0., 0.32);
197 
198  HistLim ietaLim(85, -42.5, 42.5);
199  HistLim iphiLim(74, -0.5, 73.5);
200 
201  HistLim depthLim(15, -0.5, 14.5);
202 
203  //...TDC
204  HistLim tdcLim(250, 0., 250.);
205  HistLim adcLim(256, 0., 256.);
206 
207  if (bsubdet == "HB") {
208  Ndigis = HistLim(((int)(nChannels_[1] / 100) + 1) * 100, 0., (float)((int)(nChannels_[1] / 100) + 1) * 100);
209  } else if (bsubdet == "HE") {
210  sime = HistLim(200, 0., 1.0);
211  Ndigis = HistLim(((int)(nChannels_[2] / 100) + 1) * 100, 0., (float)((int)(nChannels_[2] / 100) + 1) * 100);
212  } else if (bsubdet == "HF") {
213  sime = HistLim(100, 0., 100.);
214  pedLim = HistLim(100, 0., 20.);
215  pedWidthLim = HistLim(100, 0., 5.);
216  frac = HistLim(400, -4.00, 4.00);
217  Ndigis = HistLim(((int)(nChannels_[4] / 100) + 1) * 100, 0., (float)((int)(nChannels_[4] / 100) + 1) * 100);
218  } else if (bsubdet == "HO") {
219  sime = HistLim(200, 0., 1.0);
220  gainLim = HistLim(160, 0., 1.6);
221  Ndigis = HistLim(((int)(nChannels_[3] / 100) + 1) * 100, 0., (float)((int)(nChannels_[3] / 100) + 1) * 100);
222  }
223 
224  int isubdet = 0;
225  if (bsubdet == "HB")
226  isubdet = 1;
227  else if (bsubdet == "HE")
228  isubdet = 2;
229  else if (bsubdet == "HO")
230  isubdet = 3;
231  else if (bsubdet == "HF")
232  isubdet = 4;
233  else
234  edm::LogWarning("HcalDigisValidation") << "HcalDigisValidation Warning: not HB/HE/HF/HO " << bsubdet << std::endl;
235 
236  Char_t histo[100];
237  const char* sub = bsubdet.c_str();
238  if (bnoise == 0) {
239  // number of digis in each subdetector
240  sprintf(histo, "HcalDigiTask_Ndigis_%s", sub);
241  book1D(ib, histo, Ndigis);
242 
243  // maps of occupancies
244  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
245  sprintf(histo, "HcalDigiTask_ieta_iphi_occupancy_map_depth%d_%s", depth, sub);
246  book2D(ib, histo, ietaLim, iphiLim);
247  }
248 
249  //Depths
250  sprintf(histo, "HcalDigiTask_depths_%s", sub);
251  book1D(ib, histo, depthLim);
252 
253  // occupancies vs ieta
254  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
255  sprintf(histo, "HcalDigiTask_occupancy_vs_ieta_depth%d_%s", depth, sub);
256  book1D(ib, histo, ietaLim);
257  }
258 
259  // just 1D of all cells' amplitudes
260  sprintf(histo, "HcalDigiTask_sum_all_amplitudes_%s", sub);
261  if ((HBPhase1_ && bsubdet == "HB") || (HEPhase1_ && bsubdet == "HE"))
262  book1D(ib, histo, digiAmpWide);
263  else
264  book1D(ib, histo, digiAmp);
265 
266  sprintf(histo, "HcalDigiTask_number_of_amplitudes_above_10fC_%s", sub);
267  book1D(ib, histo, Ndigis);
268 
269  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
270  sprintf(histo, "HcalDigiTask_ADC0_adc_depth%d_%s", depth, sub);
271  book1D(ib, histo, pedestal);
272  }
273 
274  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
275  sprintf(histo, "HcalDigiTask_ADC0_fC_depth%d_%s", depth, sub);
276  book1D(ib, histo, pedestalfC);
277  }
278 
279  sprintf(histo, "HcalDigiTask_signal_amplitude_%s", sub);
280  if ((HBPhase1_ && bsubdet == "HB") || (HEPhase1_ && bsubdet == "HE"))
281  book1D(ib, histo, digiAmpWide);
282  else
283  book1D(ib, histo, digiAmp);
284 
285  if (hep17_ && bsubdet == "HE") {
286  sprintf(histo, "HcalDigiTask_signal_amplitude_HEP17");
287  book1D(ib, histo, digiAmpWide);
288  }
289  //
290  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
291  sprintf(histo, "HcalDigiTask_signal_amplitude_depth%d_%s", depth, sub);
292  if ((HBPhase1_ && bsubdet == "HB") || (HEPhase1_ && bsubdet == "HE"))
293  book1D(ib, histo, digiAmpWide);
294  else
295  book1D(ib, histo, digiAmp);
296  if (hep17_ && bsubdet == "HE") {
297  sprintf(histo, "HcalDigiTask_signal_amplitude_depth%d_HEP17", depth);
298  book1D(ib, histo, digiAmpWide);
299  }
300  }
301 
302  sprintf(histo, "HcalDigiTask_signal_amplitude_vs_bin_all_depths_%s", sub);
303  if ((HBPhase1_ && bsubdet == "HB") || (HEPhase1_ && bsubdet == "HE"))
304  book2D(ib, histo, nbin, digiAmpWide);
305  else
306  book2D(ib, histo, nbin, digiAmp);
307  if (hep17_ && bsubdet == "HE") {
308  sprintf(histo, "HcalDigiTask_signal_amplitude_vs_bin_all_depths_HEP17");
309  book2D(ib, histo, nbin, digiAmpWide);
310  }
311 
312  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
313  sprintf(histo, "HcalDigiTask_all_amplitudes_vs_bin_1D_depth%d_%s", depth, sub);
314  book1D(ib, histo, nbin);
315  if (hep17_ && bsubdet == "HE") {
316  sprintf(histo, "HcalDigiTask_all_amplitudes_vs_bin_1D_depth%d_HEP17", depth);
317  book1D(ib, histo, nbin);
318  }
319  }
320 
321  sprintf(histo, "HcalDigiTask_SOI_frac_%s", sub);
322  book1D(ib, histo, frac);
323  sprintf(histo, "HcalDigiTask_postSOI_frac_%s", sub);
324  book1D(ib, histo, frac);
325 
326  if (bmc == 1) {
327  sprintf(histo, "HcalDigiTask_amplitude_vs_simhits_%s", sub);
328  book2D(ib, histo, sime, digiAmp);
329  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
330  sprintf(histo, "HcalDigiTask_amplitude_vs_simhits_depth%d_%s", depth, sub);
331  book2D(ib, histo, sime, digiAmp);
332  }
333 
334  sprintf(histo, "HcalDigiTask_amplitude_vs_simhits_profile_%s", sub);
335  bookPf(ib, histo, sime, digiAmp);
336  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
337  sprintf(histo, "HcalDigiTask_amplitude_vs_simhits_profile_depth%d_%s", depth, sub);
338  bookPf(ib, histo, sime, digiAmp);
339  }
340 
341  sprintf(histo, "HcalDigiTask_ratio_amplitude_vs_simhits_%s", sub);
342  book1D(ib, histo, ratio);
343  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
344  sprintf(histo, "HcalDigiTask_ratio_amplitude_vs_simhits_depth%d_%s", depth, sub);
345  book1D(ib, histo, ratio);
346  }
347 
348  //...TDC
349  if (bsubdet == "HB" || bsubdet == "HE") {
350  sprintf(histo, "HcalDigiTask_TDCtime_%s", sub);
351  book1D(ib, histo, tdcLim);
352 
353  sprintf(histo, "HcalDigiTask_TDCtime_vs_ADC_%s", sub);
354  book2D(ib, histo, adcLim, tdcLim);
355  }
356 
357  } //mc only
358 
359  } else { // noise only
360 
361  // EVENT "1" distributions of all cells properties
362 
363  //KH
364  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
365  sprintf(histo, "HcalDigiTask_gain_capId0_Depth%d_%s", depth, sub);
366  book1D(ib, histo, gainLim);
367  sprintf(histo, "HcalDigiTask_gain_capId1_Depth%d_%s", depth, sub);
368  book1D(ib, histo, gainLim);
369  sprintf(histo, "HcalDigiTask_gain_capId2_Depth%d_%s", depth, sub);
370  book1D(ib, histo, gainLim);
371  sprintf(histo, "HcalDigiTask_gain_capId3_Depth%d_%s", depth, sub);
372  book1D(ib, histo, gainLim);
373 
374  sprintf(histo, "HcalDigiTask_gainWidth_capId0_Depth%d_%s", depth, sub);
375  book1D(ib, histo, gainWidthLim);
376  sprintf(histo, "HcalDigiTask_gainWidth_capId1_Depth%d_%s", depth, sub);
377  book1D(ib, histo, gainWidthLim);
378  sprintf(histo, "HcalDigiTask_gainWidth_capId2_Depth%d_%s", depth, sub);
379  book1D(ib, histo, gainWidthLim);
380  sprintf(histo, "HcalDigiTask_gainWidth_capId3_Depth%d_%s", depth, sub);
381  book1D(ib, histo, gainWidthLim);
382 
383  sprintf(histo, "HcalDigiTask_pedestal_capId0_Depth%d_%s", depth, sub);
384  book1D(ib, histo, pedLim);
385  sprintf(histo, "HcalDigiTask_pedestal_capId1_Depth%d_%s", depth, sub);
386  book1D(ib, histo, pedLim);
387  sprintf(histo, "HcalDigiTask_pedestal_capId2_Depth%d_%s", depth, sub);
388  book1D(ib, histo, pedLim);
389  sprintf(histo, "HcalDigiTask_pedestal_capId3_Depth%d_%s", depth, sub);
390  book1D(ib, histo, pedLim);
391 
392  sprintf(histo, "HcalDigiTask_pedestal_width_capId0_Depth%d_%s", depth, sub);
393  book1D(ib, histo, pedWidthLim);
394  sprintf(histo, "HcalDigiTask_pedestal_width_capId1_Depth%d_%s", depth, sub);
395  book1D(ib, histo, pedWidthLim);
396  sprintf(histo, "HcalDigiTask_pedestal_width_capId2_Depth%d_%s", depth, sub);
397  book1D(ib, histo, pedWidthLim);
398  sprintf(histo, "HcalDigiTask_pedestal_width_capId3_Depth%d_%s", depth, sub);
399  book1D(ib, histo, pedWidthLim);
400  }
401 
402  //KH
403  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
404  sprintf(histo, "HcalDigiTask_gainMap_Depth%d_%s", depth, sub);
405  book2D(ib, histo, ietaLim, iphiLim);
406  sprintf(histo, "HcalDigiTask_pwidthMap_Depth%d_%s", depth, sub);
407  book2D(ib, histo, ietaLim, iphiLim);
408  }
409 
410  } //end of noise-only
411 } //book
412 
414  conditions_ = &iSetup.getData(tok_Cond_);
415 
416  //TP Code
417  const auto& decoder = &iSetup.getData(tok_Decoder_);
418  const auto& tp_geometry = &iSetup.getData(tok_TPGeom_);
419  htopo_ = &iSetup.getData(tok_Topo_);
420 
421  //Get all handles
423  iEvent.getByToken(tok_emulTPs_, emulTPs);
424 
426  if (!skipDataTPs)
427  iEvent.getByToken(tok_dataTPs_, dataTPs);
428  //iEvent.getByLabel("hcalDigis", dataTPs);
429 
430  //~TP Code
431 
432  if (subdet_ != "all") {
433  noise_ = 0;
434  if (subdet_ == "HB") {
435  reco<HBHEDataFrame>(iEvent, iSetup, tok_hbhe_);
436  reco<QIE11DataFrame>(iEvent, iSetup, tok_qie11_hbhe_);
437  }
438  if (subdet_ == "HE") {
439  reco<HBHEDataFrame>(iEvent, iSetup, tok_hbhe_);
440  reco<QIE11DataFrame>(iEvent, iSetup, tok_qie11_hbhe_);
441  }
442  if (subdet_ == "HO")
443  reco<HODataFrame>(iEvent, iSetup, tok_ho_);
444  if (subdet_ == "HF") {
445  reco<HFDataFrame>(iEvent, iSetup, tok_hf_);
446  reco<QIE10DataFrame>(iEvent, iSetup, tok_qie10_hf_);
447  }
448 
449  if (subdet_ == "noise") {
450  noise_ = 1;
451  subdet_ = "HB";
452  reco<HBHEDataFrame>(iEvent, iSetup, tok_hbhe_);
453  reco<QIE11DataFrame>(iEvent, iSetup, tok_qie11_hbhe_);
454  subdet_ = "HE";
455  reco<HBHEDataFrame>(iEvent, iSetup, tok_hbhe_);
456  reco<QIE11DataFrame>(iEvent, iSetup, tok_qie11_hbhe_);
457  subdet_ = "HO";
458  reco<HODataFrame>(iEvent, iSetup, tok_ho_);
459  subdet_ = "HF";
460  reco<HFDataFrame>(iEvent, iSetup, tok_hf_);
461  reco<QIE10DataFrame>(iEvent, iSetup, tok_qie10_hf_);
462  subdet_ = "noise";
463  }
464  } // all subdetectors
465  else {
466  noise_ = 0;
467 
468  subdet_ = "HB";
469  reco<HBHEDataFrame>(iEvent, iSetup, tok_hbhe_);
470  reco<QIE11DataFrame>(iEvent, iSetup, tok_qie11_hbhe_);
471  subdet_ = "HE";
472  reco<HBHEDataFrame>(iEvent, iSetup, tok_hbhe_);
473  reco<QIE11DataFrame>(iEvent, iSetup, tok_qie11_hbhe_);
474  subdet_ = "HO";
475  reco<HODataFrame>(iEvent, iSetup, tok_ho_);
476  subdet_ = "HF";
477  reco<HFDataFrame>(iEvent, iSetup, tok_hf_);
478  reco<QIE10DataFrame>(iEvent, iSetup, tok_qie10_hf_);
479  subdet_ = "all";
480  }
481 
482  fill1D("nevtot", 0);
483  nevtot++;
484 
485  //TP Code
486  //Counters
487  int c = 0, chb = 0, che = 0, chf = 0, cv0 = 0, cv1 = 0, chfv0 = 0, chfv1 = 0;
488 
489  if (skipDataTPs)
490  return;
491 
492  for (HcalTrigPrimDigiCollection::const_iterator itr = dataTPs->begin(); itr != dataTPs->end(); ++itr) {
493  int ieta = itr->id().ieta();
494  int iphi = itr->id().iphi();
495 
496  HcalSubdetector subdet = (HcalSubdetector)0;
497 
498  if (abs(ieta) <= 16)
500  else if (abs(ieta) < tp_geometry->firstHFTower(itr->id().version()))
502  else if (abs(ieta) <= 42)
504 
505  //Right now, the only case where version matters is in HF
506  //If the subdetector is not HF, set version to -1
507  int tpVersion = (subdet == HcalSubdetector::HcalForward ? itr->id().version() : -1);
508 
509  float en = decoder->hcaletValue(itr->id(), itr->t0());
510 
511  if (en < 0.00001)
512  continue;
513 
514  //Plot the variables
515  //Retain classic behavior (include all tps)
516  //Additional plots that only include HF 3x2 or HF 1x1
517 
518  //Classics
519  fill1D("HcalDigiTask_tp_et", en);
520  fill2D("HcalDigiTask_tp_et_ieta", ieta, en);
521  fill2D("HcalDigiTask_tp_ieta_iphi", ieta, iphi);
522  fillPf("HcalDigiTask_tp_ave_et_ieta", ieta, en);
523  fill1D("HcalDigiTask_tp_ntp_ieta", ieta);
524  fill1D("HcalDigiTask_tp_ntp_iphi", iphi);
525  if (en > 10.)
526  fill1D("HcalDigiTask_tp_ntp_10_ieta", ieta);
527 
528  //3x2 Trig Primitives (tpVersion == 0)
529  if ((subdet != HcalSubdetector::HcalForward || tpVersion == 0) && Plot_TP_ver_) {
530  fill1D("HcalDigiTask_tp_et_v0", en);
531  fill2D("HcalDigiTask_tp_et_ieta_v0", ieta, en);
532  fill2D("HcalDigiTask_tp_ieta_iphi_v0", ieta, iphi);
533  fillPf("HcalDigiTask_tp_ave_et_ieta_v0", ieta, en);
534  fill1D("HcalDigiTask_tp_ntp_ieta_v0", ieta);
535  fill1D("HcalDigiTask_tp_ntp_iphi_v0", iphi);
536  if (en > 10.)
537  fill1D("HcalDigiTask_tp_ntp_10_ieta_v0", ieta);
538  }
539 
540  //1x1 Trig Primitives (tpVersion == 1)
541  if ((subdet != HcalSubdetector::HcalForward || tpVersion == 1) && Plot_TP_ver_) {
542  fill1D("HcalDigiTask_tp_et_v1", en);
543  fill2D("HcalDigiTask_tp_et_ieta_v1", ieta, en);
544  fill2D("HcalDigiTask_tp_ieta_iphi_v1", ieta, iphi);
545  fillPf("HcalDigiTask_tp_ave_et_ieta_v1", ieta, en);
546  fill1D("HcalDigiTask_tp_ntp_ieta_v1", ieta);
547  fill1D("HcalDigiTask_tp_ntp_iphi_v1", iphi);
548  if (en > 10.)
549  fill1D("HcalDigiTask_tp_ntp_10_ieta_v1", ieta);
550  }
551 
552  ++c;
553  if (subdet == HcalSubdetector::HcalBarrel) {
554  fill1D("HcalDigiTask_tp_et_HB", en);
555  ++chb;
556  if (Plot_TP_ver_) {
557  ++cv0;
558  ++cv1;
559  }
560  }
561  if (subdet == HcalSubdetector::HcalEndcap) {
562  fill1D("HcalDigiTask_tp_et_HE", en);
563  ++che;
564  if (Plot_TP_ver_) {
565  ++cv0;
566  ++cv1;
567  }
568  }
569  if (subdet == HcalSubdetector::HcalForward) {
570  fill1D("HcalDigiTask_tp_et_HF", en);
571  ++chf;
572 
573  if (tpVersion == 0 && Plot_TP_ver_) {
574  fill1D("HcalDigiTask_tp_et_HF_v0", en);
575  ++chfv0;
576  ++cv0;
577  }
578 
579  if (tpVersion == 1 && Plot_TP_ver_) {
580  fill1D("HcalDigiTask_tp_et_HF_v1", en);
581  ++chfv1;
582  ++cv1;
583  }
584  }
585 
586  } //end data TP collection
587 
588  fill1D("HcalDigiTask_tp_ntp", c);
589  fill1D("HcalDigiTask_tp_ntp_HB", chb);
590  fill1D("HcalDigiTask_tp_ntp_HE", che);
591  fill1D("HcalDigiTask_tp_ntp_HF", chf);
592  if (Plot_TP_ver_) {
593  fill1D("HcalDigiTask_tp_ntp_v0", cv0);
594  fill1D("HcalDigiTask_tp_ntp_v1", cv1);
595  fill1D("HcalDigiTask_tp_ntp_HF_v0", chfv0);
596  fill1D("HcalDigiTask_tp_ntp_HF_v1", chfv1);
597  }
598 
599  //~TP Code
600 }
601 
602 template <class Digi>
604  const edm::EventSetup& iSetup,
606  // HistLim =============================================================
607 
608  std::string strtmp;
609 
610  // ======================================================================
613 
614  // ADC2fC
615  CaloSamples tool;
616  iEvent.getByToken(tok, digiCollection);
617  if (!digiCollection.isValid())
618  return;
619  int isubdet = 0;
620  if (subdet_ == "HB")
621  isubdet = 1;
622  if (subdet_ == "HE")
623  isubdet = 2;
624  if (subdet_ == "HO")
625  isubdet = 3;
626  if (subdet_ == "HF")
627  isubdet = 4;
628 
629  if (isubdet == 1)
630  nevent1++;
631  if (isubdet == 2)
632  nevent2++;
633  if (isubdet == 3)
634  nevent3++;
635  if (isubdet == 4)
636  nevent4++;
637 
638  int indigis = 0;
639  // amplitude for signal cell at diff. depths
640  std::vector<double> v_ampl_c(maxDepth_[isubdet] + 1, 0);
641 
642  // is set to 1 if "seed" SimHit is found
643  int seedSimHit = 0;
644 
645  int ieta_Sim = 9999;
646  int iphi_Sim = 9999;
647  double emax_Sim = -9999.;
648 
649  // SimHits MC only
650  if (mc_ == "yes") {
652  iEvent.getByToken(tok_mc_, hcalHits);
653  const edm::PCaloHitContainer* simhitResult = hcalHits.product();
654 
655  if (isubdet != 0 && noise_ == 0) { // signal only SimHits
656 
657  for (std::vector<PCaloHit>::const_iterator simhits = simhitResult->begin(); simhits != simhitResult->end();
658  ++simhits) {
659  unsigned int id = simhits->id();
660  int sub, ieta, iphi;
661  HcalDetId hid;
662  if (testNumber_)
664  else
665  hid = HcalDetId(id);
666  sub = hid.subdet();
667  ieta = hid.ieta();
668  iphi = hid.iphi();
669 
670  double en = simhits->energy();
671 
672  if (en > emax_Sim && sub == isubdet) {
673  emax_Sim = en;
674  ieta_Sim = ieta;
675  iphi_Sim = iphi;
676  // to limit "seed" SimHit energy in case of "multi" event
677  if (mode_ == "multi" && ((sub == 4 && en < 100. && en > 1.) || ((sub != 4) && en < 1. && en > 0.02))) {
678  seedSimHit = 1;
679  break;
680  }
681  }
682 
683  } // end of SimHits cycle
684 
685  // found highest-energy SimHit for single-particle
686  if (mode_ != "multi" && emax_Sim > 0.)
687  seedSimHit = 1;
688  } // end of SimHits
689  } // end of mc_ == "yes"
690 
691  // CYCLE OVER CELLS ========================================================
692  int Ndig = 0;
693 
694  for (digiItr = digiCollection->begin(); digiItr != digiCollection->end(); digiItr++) {
695  HcalDetId cell(digiItr->id());
696  int depth = cell.depth();
697  int iphi = cell.iphi();
698  int ieta = cell.ieta();
699  int sub = cell.subdet();
700 
701  if (depth > maxDepth_[isubdet] && sub == isubdet) {
702  edm::LogWarning("HcalDetId") << "HcalDetID presents conflicting information. Depth: " << depth
703  << ", iphi: " << iphi << ", ieta: " << ieta
704  << ". Max depth from geometry is: " << maxDepth_[isubdet]
705  << ". TestNumber = " << testNumber_;
706  continue;
707  }
708 
709  // amplitude for signal cell at diff. depths
710  std::vector<double> v_ampl(maxDepth_[isubdet] + 1, 0);
711 
712  // Gains, pedestals (once !) and only for "noise" case
713  if (((nevent1 == 1 && isubdet == 1) || (nevent2 == 1 && isubdet == 2) || (nevent3 == 1 && isubdet == 3) ||
714  (nevent4 == 1 && isubdet == 4)) &&
715  noise_ == 1 && sub == isubdet) {
716  HcalGenericDetId hcalGenDetId(digiItr->id());
717  const HcalPedestal* pedestal = conditions_->getPedestal(hcalGenDetId);
718  const HcalGain* gain = conditions_->getGain(hcalGenDetId);
719  const HcalGainWidth* gainWidth = conditions_->getGainWidth(hcalGenDetId);
720  const HcalPedestalWidth* pedWidth = conditions_->getPedestalWidth(hcalGenDetId);
721 
722  for (int i = 0; i < 4; i++) {
723  fill1D("HcalDigiTask_gain_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_, gain->getValue(i));
724  fill1D("HcalDigiTask_gainWidth_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_, gainWidth->getValue(i));
725  fill1D("HcalDigiTask_pedestal_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_, pedestal->getValue(i));
726  fill1D("HcalDigiTask_pedestal_width_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_,
727  pedWidth->getWidth(i));
728  }
729 
730  fill2D("HcalDigiTask_gainMap_Depth" + str(depth) + "_" + subdet_, double(ieta), double(iphi), gain->getValue(0));
731  fill2D("HcalDigiTask_pwidthMap_Depth" + str(depth) + "_" + subdet_,
732  double(ieta),
733  double(iphi),
734  pedWidth->getWidth(0));
735 
736  } // end of event #1
737 
738  if (sub == isubdet)
739  Ndig++; // subdet number of digi
740 
741  // No-noise case, only single subdet selected ===========================
742 
743  if (sub == isubdet && noise_ == 0) {
744  HcalCalibrations calibrations = conditions_->getHcalCalibrations(cell);
745 
746  const HcalQIECoder* channelCoder = conditions_->getHcalCoder(cell);
747  const HcalQIEShape* shape = conditions_->getHcalShape(channelCoder);
748  HcalCoderDb coder(*channelCoder, *shape);
749  coder.adc2fC(*digiItr, tool);
750 
751  // for dynamic digi time sample analysis
752  int soi = tool.presamples();
753  int lastbin = tool.size() - 1;
754 
755  double noiseADC = (*digiItr)[0].adc();
756  double noisefC = tool[0];
757  // noise evaluations from "pre-samples"
758  fill1D("HcalDigiTask_ADC0_adc_depth" + str(depth) + "_" + subdet_, noiseADC);
759  fill1D("HcalDigiTask_ADC0_fC_depth" + str(depth) + "_" + subdet_, noisefC);
760 
761  // OCCUPANCY maps fill
762  fill2D("HcalDigiTask_ieta_iphi_occupancy_map_depth" + str(depth) + "_" + subdet_, double(ieta), double(iphi));
763 
764  fill1D("HcalDigiTask_depths_" + subdet_, double(depth));
765 
766  // Cycle on time slices
767  // - for each Digi
768  // - for one Digi with max SimHits E in subdet
769 
770  int closen = 0; // =1 if 1) seedSimHit = 1 and 2) the cell is the same
771  if (ieta == ieta_Sim && iphi == iphi_Sim)
772  closen = seedSimHit;
773 
774  for (int ii = 0; ii < tool.size(); ii++) {
775  int capid = (*digiItr)[ii].capid();
776  // single ts amplitude
777  double val = (tool[ii] - calibrations.pedestal(capid));
778 
779  if (val > 100.) {
780  fill1D("HcalDigiTask_ADC0_adc_depth" + str(depth) + "_" + subdet_, noiseADC);
781  strtmp = "HcalDigiTask_all_amplitudes_vs_bin_1D_depth" + str(depth) + "_" + subdet_;
782  fill1D(strtmp, double(ii), val);
783  }
784 
785  if (closen == 1) {
786  strtmp = "HcalDigiTask_signal_amplitude_vs_bin_all_depths_" + subdet_;
787  fill2D(strtmp, double(ii), val);
788  }
789 
790  // all detectors
791  if (ii >= soi && ii <= lastbin) {
792  v_ampl[0] += val;
793  v_ampl[depth] += val;
794 
795  if (closen == 1) {
796  v_ampl_c[0] += val;
797  v_ampl_c[depth] += val;
798  }
799  }
800  }
801  // end of time bucket sample
802 
803  // maps of sum of amplitudes (sum lin.digis(4,5,6,7) - ped) all depths
804  // just 1D of all cells' amplitudes
805  strtmp = "HcalDigiTask_sum_all_amplitudes_" + subdet_;
806  fill1D(strtmp, v_ampl[0]);
807 
808  std::vector<int> v_ampl_sub(v_ampl.begin() + 1, v_ampl.end()); // remove element 0, which is the sum of any depth
809  double ampl_max = *std::max_element(v_ampl_sub.begin(), v_ampl_sub.end());
810  if (ampl_max > 10.)
811  indigis++;
812  //KH if (ampl1 > 10. || ampl2 > 10. || ampl3 > 10. || ampl4 > 10.) indigis++;
813 
814  // fraction 5,6 bins if ampl. is big.
815  if (v_ampl[depth] > 30.) {
816  double fbinSOI = tool[soi] - calibrations.pedestal((*digiItr)[soi].capid());
817  double fbinPS = 0;
818 
819  for (int j = soi + 1; j <= lastbin; j++)
820  fbinPS += tool[j] - calibrations.pedestal((*digiItr)[j].capid());
821 
822  fbinSOI /= v_ampl[depth];
823  fbinPS /= v_ampl[depth];
824  strtmp = "HcalDigiTask_SOI_frac_" + subdet_;
825  fill1D(strtmp, fbinSOI);
826  strtmp = "HcalDigiTask_postSOI_frac_" + subdet_;
827  fill1D(strtmp, fbinPS);
828  }
829 
830  strtmp = "HcalDigiTask_signal_amplitude_" + subdet_;
831  fill1D(strtmp, v_ampl[0]);
832  strtmp = "HcalDigiTask_signal_amplitude_depth" + str(depth) + "_" + subdet_;
833  fill1D(strtmp, v_ampl[depth]);
834  }
835  } // End of CYCLE OVER CELLS =============================================
836 
837  if (isubdet != 0 && noise_ == 0) { // signal only, once per event
838  strtmp = "HcalDigiTask_number_of_amplitudes_above_10fC_" + subdet_;
839  fill1D(strtmp, indigis);
840 
841  // SimHits once again !!!
842  double eps = 1.e-3;
843  std::vector<double> v_ehits(maxDepth_[isubdet] + 1, 0);
844 
845  if (mc_ == "yes") {
847  iEvent.getByToken(tok_mc_, hcalHits);
848  const edm::PCaloHitContainer* simhitResult = hcalHits.product();
849  for (std::vector<PCaloHit>::const_iterator simhits = simhitResult->begin(); simhits != simhitResult->end();
850  ++simhits) {
851  unsigned int id = simhits->id();
852  int sub, depth, ieta, iphi;
853  HcalDetId hid;
854  if (testNumber_)
856  else
857  hid = HcalDetId(id);
858  sub = hid.subdet();
859  depth = hid.depth();
860  ieta = hid.ieta();
861  iphi = hid.iphi();
862 
863  if (depth > maxDepth_[isubdet] && sub == isubdet) {
864  edm::LogWarning("HcalDetId") << "HcalDetID(SimHit) presents conflicting information. Depth: " << depth
865  << ", iphi: " << iphi << ", ieta: " << ieta
866  << ". Max depth from geometry is: " << maxDepth_[isubdet]
867  << ". TestNumber = " << testNumber_;
868  continue;
869  }
870 
871  // take cell already found to be max energy in a particular subdet
872  if (sub == isubdet && ieta == ieta_Sim && iphi == iphi_Sim) {
873  double en = simhits->energy();
874 
875  v_ehits[0] += en;
876  v_ehits[depth] += en;
877  }
878  } // simhit loop
879 
880  strtmp = "HcalDigiTask_amplitude_vs_simhits_" + subdet_;
881  if (v_ehits[0] > eps)
882  fill2D(strtmp, v_ehits[0], v_ampl_c[0]);
883  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
884  strtmp = "HcalDigiTask_amplitude_vs_simhits_depth" + str(depth) + "_" + subdet_;
885  if (v_ehits[depth] > eps)
886  fill2D(strtmp, v_ehits[depth], v_ampl_c[depth]);
887  }
888 
889  strtmp = "HcalDigiTask_amplitude_vs_simhits_profile_" + subdet_;
890  if (v_ehits[0] > eps)
891  fillPf(strtmp, v_ehits[0], v_ampl_c[0]);
892  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
893  strtmp = "HcalDigiTask_amplitude_vs_simhits_profile_depth" + str(depth) + "_" + subdet_;
894  if (v_ehits[depth] > eps)
895  fillPf(strtmp, v_ehits[depth], v_ampl_c[depth]);
896  }
897 
898  strtmp = "HcalDigiTask_ratio_amplitude_vs_simhits_" + subdet_;
899  if (v_ehits[0] > eps)
900  fill1D(strtmp, v_ampl_c[0] / v_ehits[0]);
901  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
902  strtmp = "HcalDigiTask_amplitude_vs_simhits_profile_depth" + str(depth) + "_" + subdet_;
903  if (v_ehits[depth] > eps)
904  fillPf(strtmp, v_ehits[depth], v_ampl_c[depth]);
905  strtmp = "HcalDigiTask_ratio_amplitude_vs_simhits_depth" + str(depth) + "_" + subdet_;
906  if (v_ehits[depth] > eps)
907  fill1D(strtmp, v_ampl_c[depth] / v_ehits[depth]);
908  }
909 
910  } // end of if(mc_ == "yes")
911 
912  strtmp = "HcalDigiTask_Ndigis_" + subdet_;
913  fill1D(strtmp, double(Ndig));
914 
915  } // end of if( subdet != 0 && noise_ == 0) { // signal only
916 }
917 template <class dataFrameType>
919  const edm::EventSetup& iSetup,
921  // HistLim =============================================================
922 
923  std::string strtmp;
924 
925  // ======================================================================
927  //typename HcalDataFrameContainer<dataFrameType>::const_iterator digiItr;
928 
929  // ADC2fC
930  CaloSamples tool;
931  iEvent.getByToken(tok, digiHandle);
932  if (!digiHandle.isValid())
933  return;
935  int isubdet = 0;
936  if (subdet_ == "HB")
937  isubdet = 1;
938  if (subdet_ == "HE")
939  isubdet = 2;
940  if (subdet_ == "HO")
941  isubdet = 3;
942  if (subdet_ == "HF")
943  isubdet = 4;
944 
945  if (isubdet == 1)
946  nevent1++;
947  if (isubdet == 2)
948  nevent2++;
949  if (isubdet == 3)
950  nevent3++;
951  if (isubdet == 4)
952  nevent4++;
953 
954  int indigis = 0;
955  // amplitude for signal cell at diff. depths
956  std::vector<double> v_ampl_c(maxDepth_[isubdet] + 1, 0);
957 
958  // is set to 1 if "seed" SimHit is found
959  int seedSimHit = 0;
960 
961  int ieta_Sim = 9999;
962  int iphi_Sim = 9999;
963  double emax_Sim = -9999.;
964 
965  // SimHits MC only
966  if (mc_ == "yes") {
968  iEvent.getByToken(tok_mc_, hcalHits);
969  const edm::PCaloHitContainer* simhitResult = hcalHits.product();
970 
971  if (isubdet != 0 && noise_ == 0) { // signal only SimHits
972 
973  for (std::vector<PCaloHit>::const_iterator simhits = simhitResult->begin(); simhits != simhitResult->end();
974  ++simhits) {
975  unsigned int id = simhits->id();
976  int sub, ieta, iphi;
977  HcalDetId hid;
978  if (testNumber_)
980  else
981  hid = HcalDetId(id);
982  sub = hid.subdet();
983  ieta = hid.ieta();
984  iphi = hid.iphi();
985 
986  double en = simhits->energy();
987 
988  if (en > emax_Sim && sub == isubdet) {
989  emax_Sim = en;
990  ieta_Sim = ieta;
991  iphi_Sim = iphi;
992  // to limit "seed" SimHit energy in case of "multi" event
993  if (mode_ == "multi" && ((sub == 4 && en < 100. && en > 1.) || ((sub != 4) && en < 1. && en > 0.02))) {
994  seedSimHit = 1;
995  break;
996  }
997  }
998 
999  } // end of SimHits cycle
1000 
1001  // found highest-energy SimHit for single-particle
1002  if (mode_ != "multi" && emax_Sim > 0.)
1003  seedSimHit = 1;
1004  } // end of SimHits
1005  } // end of mc_ == "yes"
1006 
1007  // CYCLE OVER CELLS ========================================================
1008  int Ndig = 0;
1009 
1010  for (typename HcalDataFrameContainer<dataFrameType>::const_iterator digiItr = digiCollection->begin();
1011  digiItr != digiCollection->end();
1012  digiItr++) {
1013  dataFrameType dataFrame = *digiItr;
1014 
1015  HcalDetId cell(digiItr->id());
1016  int depth = cell.depth();
1017  int iphi = cell.iphi();
1018  int ieta = cell.ieta();
1019  int sub = cell.subdet();
1020 
1021  //Is this in HEP17
1022  bool isHEP17 = (iphi >= 63) && (iphi <= 66) && (ieta > 0) && (sub == 2);
1023 
1024  if (depth > maxDepth_[isubdet] && sub == isubdet) {
1025  edm::LogWarning("HcalDetId") << "HcalDetID presents conflicting information. Depth: " << depth
1026  << ", iphi: " << iphi << ", ieta: " << ieta
1027  << ". Max depth from geometry is: " << maxDepth_[isubdet]
1028  << ". TestNumber = " << testNumber_;
1029  continue;
1030  }
1031 
1032  // amplitude for signal cell at diff. depths
1033  std::vector<double> v_ampl(maxDepth_[isubdet] + 1, 0);
1034 
1035  // Gains, pedestals (once !) and only for "noise" case
1036  if (((nevent1 == 1 && isubdet == 1) || (nevent2 == 1 && isubdet == 2) || (nevent3 == 1 && isubdet == 3) ||
1037  (nevent4 == 1 && isubdet == 4)) &&
1038  noise_ == 1 && sub == isubdet) {
1039  HcalGenericDetId hcalGenDetId(digiItr->id());
1040  const HcalPedestal* pedestal = conditions_->getPedestal(hcalGenDetId);
1041  const HcalGain* gain = conditions_->getGain(hcalGenDetId);
1042  const HcalGainWidth* gainWidth = conditions_->getGainWidth(hcalGenDetId);
1043  const HcalPedestalWidth* pedWidth = conditions_->getPedestalWidth(hcalGenDetId);
1044 
1045  for (int i = 0; i < 4; i++) {
1046  fill1D("HcalDigiTask_gain_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_, gain->getValue(i));
1047  fill1D("HcalDigiTask_gainWidth_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_, gainWidth->getValue(i));
1048  fill1D("HcalDigiTask_pedestal_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_, pedestal->getValue(i));
1049  fill1D("HcalDigiTask_pedestal_width_capId" + str(i) + "_Depth" + str(depth) + "_" + subdet_,
1050  pedWidth->getWidth(i));
1051  }
1052 
1053  fill2D("HcalDigiTask_gainMap_Depth" + str(depth) + "_" + subdet_, double(ieta), double(iphi), gain->getValue(0));
1054  fill2D("HcalDigiTask_pwidthMap_Depth" + str(depth) + "_" + subdet_,
1055  double(ieta),
1056  double(iphi),
1057  pedWidth->getWidth(0));
1058 
1059  } // end of event #1
1060  //edm::LogVerbatim("OutputInfo") << "==== End of event noise block in cell cycle";
1061 
1062  if (sub == isubdet)
1063  Ndig++; // subdet number of digi
1064 
1065  // No-noise case, only single subdet selected ===========================
1066 
1067  if (sub == isubdet && noise_ == 0) {
1068  HcalCalibrations calibrations = conditions_->getHcalCalibrations(cell);
1069 
1070  const HcalQIECoder* channelCoder = conditions_->getHcalCoder(cell);
1071  const HcalQIEShape* shape = conditions_->getHcalShape(channelCoder);
1072  HcalCoderDb coder(*channelCoder, *shape);
1073  coder.adc2fC(dataFrame, tool);
1074 
1075  // for dynamic digi time sample analysis
1076  int soi = tool.presamples();
1077  int lastbin = tool.size() - 1;
1078 
1079  double noiseADC = (dataFrame)[0].adc();
1080  double noisefC = tool[0];
1081  // noise evaluations from "pre-samples"
1082  fill1D("HcalDigiTask_ADC0_adc_depth" + str(depth) + "_" + subdet_, noiseADC);
1083  fill1D("HcalDigiTask_ADC0_fC_depth" + str(depth) + "_" + subdet_, noisefC);
1084 
1085  // OCCUPANCY maps fill
1086  fill2D("HcalDigiTask_ieta_iphi_occupancy_map_depth" + str(depth) + "_" + subdet_, double(ieta), double(iphi));
1087 
1088  fill1D("HcalDigiTask_depths_" + subdet_, double(depth));
1089 
1090  // Cycle on time slices
1091  // - for each Digi
1092  // - for one Digi with max SimHits E in subdet
1093 
1094  int closen = 0; // =1 if 1) seedSimHit = 1 and 2) the cell is the same
1095  if (ieta == ieta_Sim && iphi == iphi_Sim)
1096  closen = seedSimHit;
1097 
1098  for (int ii = 0; ii < tool.size(); ii++) {
1099  int capid = (dataFrame)[ii].capid();
1100  // single ts amplitude
1101  double val = (tool[ii] - calibrations.pedestal(capid));
1102 
1103  if (val > 100.) {
1104  fill1D("HcalDigiTask_ADC0_adc_depth" + str(depth) + "_" + subdet_, noiseADC);
1105  if (hep17_) {
1106  if (!isHEP17) {
1107  strtmp = "HcalDigiTask_all_amplitudes_vs_bin_1D_depth" + str(depth) + "_" + subdet_;
1108  fill1D(strtmp, double(ii), val);
1109  } else {
1110  strtmp = "HcalDigiTask_all_amplitudes_vs_bin_1D_depth" + str(depth) + "_HEP17";
1111  fill1D(strtmp, double(ii), val);
1112  }
1113  } else {
1114  strtmp = "HcalDigiTask_all_amplitudes_vs_bin_1D_depth" + str(depth) + "_" + subdet_;
1115  fill1D(strtmp, double(ii), val);
1116  }
1117  }
1118 
1119  if (closen == 1) {
1120  if (hep17_) {
1121  if (!isHEP17) {
1122  strtmp = "HcalDigiTask_signal_amplitude_vs_bin_all_depths_" + subdet_;
1123  fill2D(strtmp, double(ii), val);
1124  } else {
1125  strtmp = "HcalDigiTask_signal_amplitude_vs_bin_all_depths_HEP17";
1126  fill2D(strtmp, double(ii), val);
1127  }
1128  } else {
1129  strtmp = "HcalDigiTask_signal_amplitude_vs_bin_all_depths_" + subdet_;
1130  fill2D(strtmp, double(ii), val);
1131  }
1132  }
1133 
1134  // all detectors
1135  if (ii >= soi && ii <= lastbin) {
1136  v_ampl[0] += val;
1137  v_ampl[depth] += val;
1138 
1139  if (closen == 1) {
1140  v_ampl_c[0] += val;
1141  v_ampl_c[depth] += val;
1142  }
1143  }
1144 
1145  //...TDC
1146 
1147  if ((HBPhase1_ && sub == 1) || (HEPhase1_ && sub == 2)) {
1148  double digiADC = (dataFrame)[ii].adc();
1149  const QIE11DataFrame dataFrameHBHE = static_cast<const QIE11DataFrame>(*digiItr);
1150  double digiTDC = (dataFrameHBHE)[ii].tdc();
1151  if (digiTDC < 50) {
1152  double time = ii * 25. + (digiTDC * 0.5);
1153  strtmp = "HcalDigiTask_TDCtime_" + subdet_;
1154  fill1D(strtmp, time);
1155 
1156  strtmp = "HcalDigiTask_TDCtime_vs_ADC_" + subdet_;
1157  fill2D(strtmp, digiADC, time);
1158  }
1159  }
1160  }
1161  // end of time bucket sample
1162 
1163  // just 1D of all cells' amplitudes
1164  strtmp = "HcalDigiTask_sum_all_amplitudes_" + subdet_;
1165  fill1D(strtmp, v_ampl[0]);
1166 
1167  std::vector<int> v_ampl_sub(v_ampl.begin() + 1, v_ampl.end()); // remove element 0, which is the sum of any depth
1168  double ampl_max = *std::max_element(v_ampl_sub.begin(), v_ampl_sub.end());
1169  if (ampl_max > 10.)
1170  indigis++;
1171  //KH if (ampl1 > 10. || ampl2 > 10. || ampl3 > 10. || ampl4 > 10.) indigis++;
1172 
1173  // fraction 5,6 bins if ampl. is big.
1174  //histogram names have not been changed, but it should be understood that bin_5 is soi, and bin_6_7 is latter TS'
1175  if ((v_ampl[depth] > 30. && (subdet_ != "HE" || subdet_ != "HB")) ||
1176  (v_ampl[depth] > 300.)) { //300 fC cut for QIE-11 HB & HE
1177  double fbinSOI = tool[soi] - calibrations.pedestal((dataFrame)[soi].capid());
1178  double fbinPS = 0;
1179 
1180  for (int j = soi + 1; j <= lastbin; j++)
1181  fbinPS += tool[j] - calibrations.pedestal((dataFrame)[j].capid());
1182 
1183  fbinSOI /= v_ampl[depth];
1184  fbinPS /= v_ampl[depth];
1185  strtmp = "HcalDigiTask_SOI_frac_" + subdet_;
1186  fill1D(strtmp, fbinSOI);
1187  strtmp = "HcalDigiTask_postSOI_frac_" + subdet_;
1188  fill1D(strtmp, fbinPS);
1189  }
1190 
1191  if (hep17_) {
1192  if (!isHEP17) {
1193  strtmp = "HcalDigiTask_signal_amplitude_" + subdet_;
1194  fill1D(strtmp, v_ampl[0]);
1195  strtmp = "HcalDigiTask_signal_amplitude_depth" + str(depth) + "_" + subdet_;
1196  fill1D(strtmp, v_ampl[depth]);
1197  } else {
1198  strtmp = "HcalDigiTask_signal_amplitude_HEP17";
1199  fill1D(strtmp, v_ampl[0]);
1200  strtmp = "HcalDigiTask_signal_amplitude_depth" + str(depth) + "_HEP17";
1201  fill1D(strtmp, v_ampl[depth]);
1202  }
1203  } else {
1204  strtmp = "HcalDigiTask_signal_amplitude_" + subdet_;
1205  fill1D(strtmp, v_ampl[0]);
1206  strtmp = "HcalDigiTask_signal_amplitude_depth" + str(depth) + "_" + subdet_;
1207  fill1D(strtmp, v_ampl[depth]);
1208  }
1209  }
1210  } // End of CYCLE OVER CELLS =============================================
1211 
1212  if (isubdet != 0 && noise_ == 0) { // signal only, once per event
1213  strtmp = "HcalDigiTask_number_of_amplitudes_above_10fC_" + subdet_;
1214  fill1D(strtmp, indigis);
1215 
1216  // SimHits once again !!!
1217  double eps = 1.e-3;
1218  std::vector<double> v_ehits(maxDepth_[isubdet] + 1, 0);
1219 
1220  if (mc_ == "yes") {
1222  iEvent.getByToken(tok_mc_, hcalHits);
1223  const edm::PCaloHitContainer* simhitResult = hcalHits.product();
1224  for (std::vector<PCaloHit>::const_iterator simhits = simhitResult->begin(); simhits != simhitResult->end();
1225  ++simhits) {
1226  unsigned int id = simhits->id();
1227  int sub, depth, ieta, iphi;
1228  HcalDetId hid;
1229  if (testNumber_)
1231  else
1232  hid = HcalDetId(id);
1233  sub = hid.subdet();
1234  depth = hid.depth();
1235  ieta = hid.ieta();
1236  iphi = hid.iphi();
1237 
1238  if (depth > maxDepth_[isubdet] && sub == isubdet) {
1239  edm::LogWarning("HcalDetId") << "HcalDetID(SimHit) presents conflicting information. Depth: " << depth
1240  << ", iphi: " << iphi << ", ieta: " << ieta
1241  << ". Max depth from geometry is: " << maxDepth_[isubdet]
1242  << ". TestNumber = " << testNumber_;
1243  continue;
1244  }
1245 
1246  // take cell already found to be max energy in a particular subdet
1247  if (sub == isubdet && ieta == ieta_Sim && iphi == iphi_Sim) {
1248  double en = simhits->energy();
1249 
1250  v_ehits[0] += en;
1251  v_ehits[depth] += en;
1252  }
1253  } // simhit loop
1254 
1255  strtmp = "HcalDigiTask_amplitude_vs_simhits_" + subdet_;
1256  if (v_ehits[0] > eps)
1257  fill2D(strtmp, v_ehits[0], v_ampl_c[0]);
1258  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
1259  strtmp = "HcalDigiTask_amplitude_vs_simhits_depth" + str(depth) + "_" + subdet_;
1260  if (v_ehits[depth] > eps)
1261  fill2D(strtmp, v_ehits[depth], v_ampl_c[depth]);
1262  }
1263 
1264  strtmp = "HcalDigiTask_amplitude_vs_simhits_profile_" + subdet_;
1265  if (v_ehits[0] > eps)
1266  fillPf(strtmp, v_ehits[0], v_ampl_c[0]);
1267  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
1268  strtmp = "HcalDigiTask_amplitude_vs_simhits_profile_depth" + str(depth) + "_" + subdet_;
1269  if (v_ehits[depth] > eps)
1270  fillPf(strtmp, v_ehits[depth], v_ampl_c[depth]);
1271  }
1272 
1273  strtmp = "HcalDigiTask_ratio_amplitude_vs_simhits_" + subdet_;
1274  if (v_ehits[0] > eps)
1275  fill1D(strtmp, v_ampl_c[0] / v_ehits[0]);
1276  for (int depth = 1; depth <= maxDepth_[isubdet]; depth++) {
1277  strtmp = "HcalDigiTask_amplitude_vs_simhits_profile_depth" + str(depth) + "_" + subdet_;
1278  if (v_ehits[depth] > eps)
1279  fillPf(strtmp, v_ehits[depth], v_ampl_c[depth]);
1280  strtmp = "HcalDigiTask_ratio_amplitude_vs_simhits_depth" + str(depth) + "_" + subdet_;
1281  if (v_ehits[depth] > eps)
1282  fill1D(strtmp, v_ampl_c[depth] / v_ehits[depth]);
1283  }
1284 
1285  } // end of if(mc_ == "yes")
1286 
1287  strtmp = "HcalDigiTask_Ndigis_" + subdet_;
1288  fill1D(strtmp, double(Ndig));
1289 
1290  } // end of if( subdet != 0 && noise_ == 0) { // signal only
1291 
1292 } //HcalDataFrameContainer
1293 
1295  if (!msm_->count(name))
1296  (*msm_)[name] = ib.book1D(name.c_str(), name.c_str(), n, min, max);
1297 }
1298 
1300  if (!msm_->count(name))
1301  (*msm_)[name] = ib.book1D(name.c_str(), name.c_str(), limX.n, limX.min, limX.max);
1302 }
1303 
1305  msm_->find(name)->second->Fill(X, weight);
1306 }
1307 
1309  if (!msm_->count(name))
1310  (*msm_)[name] = ib.book2D(name.c_str(), name.c_str(), limX.n, limX.min, limX.max, limY.n, limY.min, limY.max);
1311 }
1312 
1313 void HcalDigisValidation::fill2D(std::string name, double X, double Y, double weight) {
1314  msm_->find(name)->second->Fill(X, Y, weight);
1315 }
1316 
1318  if (!msm_->count(name))
1319  (*msm_)[name] = ib.bookProfile(name.c_str(), name.c_str(), limX.n, limX.min, limX.max, limY.n, limY.min, limY.max);
1320 }
1321 
1323  DQMStore::IBooker& ib, std::string name, const HistLim& limX, const HistLim& limY, const char* option) {
1324  if (!msm_->count(name))
1325  (*msm_)[name] =
1326  ib.bookProfile(name.c_str(), name.c_str(), limX.n, limX.min, limX.max, limY.n, limY.min, limY.max, option);
1327 }
1328 
1329 void HcalDigisValidation::fillPf(std::string name, double X, double Y) { msm_->find(name)->second->Fill(X, Y); }
1330 
1332  if (!msm_->count(name))
1333  return nullptr;
1334  else
1335  return msm_->find(name)->second;
1336 }
1337 
1339  std::stringstream out;
1340  out << x;
1341  return out.str();
1342 }
1343 
1344 //define this as a plug-in
Log< level::Info, true > LogVerbatim
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:34
const HcalGainWidth * getGainWidth(const HcalGenericDetId &fId) const
T getUntrackedParameter(std::string const &, T const &) const
edm::EDGetTokenT< HODigiCollection > tok_ho_
std::vector< PCaloHit > PCaloHitContainer
const edm::EventSetup & c
int ib
Definition: cuy.py:661
edm::ESGetToken< CaloTPGTranscoder, CaloTPGRecord > tok_Decoder_
const HcalDDDRecConstants * hcons_
WorkSpace int float eps
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:32
const HcalDbService * conditions_
unsigned int getHxSize(const int type) const
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
void dqmBeginRun(const edm::Run &run, const edm::EventSetup &c) override
std::vector< T >::const_iterator const_iterator
bool exists(std::string const &parameterName) const
checks if a parameter exists
#define X(str)
Definition: MuonsGrabber.cc:38
int presamples() const
access presample information
Definition: CaloSamples.h:36
edm::ESGetToken< HcalDDDRecConstants, HcalRecNumberingRecord > tok_HRNDC_
void reco(const edm::Event &iEvent, const edm::EventSetup &iSetup, const edm::EDGetTokenT< edm::SortedCollection< Digi > > &tok)
const_iterator begin() const
The iterator returned can not safely be used across threads.
edm::EDGetTokenT< QIE11DigiCollection > tok_qie11_hbhe_
float getValue(int fCapId) const
get value for capId = 0..3
Definition: HcalGain.h:21
int ii
Definition: cuy.py:589
const HcalTopology * htopo_
float getValue(int fCapId) const
get value for capId = 0..3
Definition: HcalGainWidth.h:20
void analyze(const edm::Event &, const edm::EventSetup &) override
bool getData(T &iHolder) const
Definition: EventSetup.h:128
const HcalPedestalWidth * getPedestalWidth(const HcalGenericDetId &fId) const
constexpr HcalSubdetector subdet() const
get the subdetector
Definition: HcalDetId.h:138
edm::ESGetToken< HcalDbService, HcalDbRecord > tok_Cond_
void adc2fC(const HBHEDataFrame &df, CaloSamples &lf) const override
Definition: HcalCoderDb.cc:73
int iEvent
Definition: GenABIO.cc:224
constexpr int iphi() const
get the cell iphi
Definition: HcalDetId.h:157
MonitorElement * bookProfile(TString const &name, TString const &title, int nchX, double lowX, double highX, int, double lowY, double highY, char const *option="s", FUNC onbooking=NOOP())
Definition: DQMStore.h:322
void fillPf(std::string name, double X, double Y)
void fill2D(std::string name, double X, double Y, double weight=1)
HcalSubdetector
Definition: HcalAssistant.h:31
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
uint16_t const *__restrict__ x
Definition: gpuClustering.h:39
constexpr int ieta() const
get the cell ieta
Definition: HcalDetId.h:155
HcalDigisValidation(const edm::ParameterSet &)
void bookPf(DQMStore::IBooker &ib, std::string name, const HistLim &limX, const HistLim &limY)
T min(T a, T b)
Definition: MathUtil.h:58
edm::ESGetToken< CaloGeometry, CaloGeometryRecord > tok_Geom_
bool isValid() const
Definition: HandleBase.h:70
edm::EDGetTokenT< HBHEDigiCollection > tok_hbhe_
edm::InputTag QIE10inputTag_
void book2D(DQMStore::IBooker &ib, std::string name, const HistLim &limX, const HistLim &limY)
std::map< std::string, MonitorElement * > * msm_
constexpr double pedestal(int fCapId) const
get pedestal for capid=0..3
T const * product() const
Definition: Handle.h:70
void fill1D(std::string name, double X, double weight=1)
int size() const
get the size
Definition: CaloSamples.h:24
edm::EDGetTokenT< edm::PCaloHitContainer > tok_mc_
MonitorElement * monitor(std::string name)
int getMaxDepth(const int &type) const
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:177
void booking(DQMStore::IBooker &ib, std::string subdetopt, int bnoise, int bmc)
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
const HcalGain * getGain(const HcalGenericDetId &fId) const
float getWidth(int fCapId) const
get width (sqrt(sigma_i_i)) for capId = 0..3
const HcalQIECoder * getHcalCoder(const HcalGenericDetId &fId) const
edm::EDGetTokenT< HcalTrigPrimDigiCollection > tok_emulTPs_
edm::InputTag QIE11inputTag_
boost::transform_iterator< IterHelp, boost::counting_iterator< int > > const_iterator
const HcalQIEShape * getHcalShape(const HcalGenericDetId &fId) const
std::string str(int x)
const_iterator end() const
edm::EDGetTokenT< HcalTrigPrimDigiCollection > tok_dataTPs_
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
edm::ESGetToken< HcalTrigTowerGeometry, CaloGeometryRecord > tok_TPGeom_
void book1D(DQMStore::IBooker &ib, std::string name, int n, double min, double max)
DetId relabel(const uint32_t testId) const
constexpr int depth() const
get the tower depth
Definition: HcalDetId.h:164
Log< level::Warning, false > LogWarning
int weight
Definition: histoStyle.py:51
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
const HcalCalibrations & getHcalCalibrations(const HcalGenericDetId &fId) const
edm::EDGetTokenT< QIE10DigiCollection > tok_qie10_hf_
const HcalPedestal * getPedestal(const HcalGenericDetId &fId) const
Definition: Run.h:45
uint16_t *__restrict__ uint16_t const *__restrict__ adc
edm::ESGetToken< HcalTopology, HcalRecNumberingRecord > tok_Topo_
edm::EDGetTokenT< HFDigiCollection > tok_hf_