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