CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
EmDQMReco.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <string>
3 #include <memory>
4 
5 #include <fmt/printf.h>
6 
8 // Root include files //
10 #include <TDirectory.h>
11 #include <TFile.h>
12 #include <TH1F.h>
13 #include <Math/VectorUtil.h>
14 
16 // Header file for this //
19 
21 // Collaborating Class Header //
42 
43 using namespace ROOT::Math::VectorUtil;
44 
45 //----------------------------------------------------------------------
46 // class EmDQMReco::FourVectorMonitorElements
47 //----------------------------------------------------------------------
49  DQMStore::IBooker &iBooker,
50  const std::string &histogramNameTemplate,
51  const std::string &histogramTitleTemplate)
52  : parent(_parent) {
53  // introducing variables for better code readability later on
54  std::string histName;
55  std::string histTitle;
56 
57  // et
58  histName = fmt::sprintf(histogramNameTemplate, "et");
59  histTitle = fmt::sprintf(histogramTitleTemplate, "E_{T}");
61  iBooker.book1D(histName.c_str(), histTitle.c_str(), parent->plotBins, parent->plotPtMin, parent->plotPtMax);
62 
63  // eta
64  histName = fmt::sprintf(histogramNameTemplate, "eta");
65  histTitle = fmt::sprintf(histogramTitleTemplate, "#eta");
67  iBooker.book1D(histName.c_str(), histTitle.c_str(), parent->plotBins, -parent->plotEtaMax, parent->plotEtaMax);
68 
69  // phi
70  histName = fmt::sprintf(histogramNameTemplate, "phi");
71  histTitle = fmt::sprintf(histogramTitleTemplate, "#phi");
73  iBooker.book1D(histName.c_str(), histTitle.c_str(), parent->plotBins, -parent->plotPhiMax, parent->plotPhiMax);
74 }
75 
76 //----------------------------------------------------------------------
77 
79  etMonitorElement->Fill(momentum.Et());
80  etaMonitorElement->Fill(momentum.eta());
81  phiMonitorElement->Fill(momentum.phi());
82 }
83 
84 //----------------------------------------------------------------------
85 
87 // Constructor //
91  // Read from configuration file //
93  dirname_ = "HLT/HLTEgammaValidation/" + pset.getParameter<std::string>("@module_label");
94 
95  // parameters for generator study
96  reqNum = pset.getParameter<unsigned int>("reqNum");
97  pdgGen = pset.getParameter<int>("pdgGen");
98  recoEtaAcc = pset.getParameter<double>("genEtaAcc");
99  recoEtAcc = pset.getParameter<double>("genEtAcc");
100  // plotting parameters (untracked because they don't affect the physics)
101  plotPtMin = pset.getUntrackedParameter<double>("PtMin", 0.);
102  plotPtMax = pset.getUntrackedParameter<double>("PtMax", 1000.);
103  plotEtaMax = pset.getUntrackedParameter<double>("EtaMax", 2.7);
104  plotPhiMax = pset.getUntrackedParameter<double>("PhiMax", 3.15);
105  plotBins = pset.getUntrackedParameter<unsigned int>("Nbins", 50);
106  useHumanReadableHistTitles = pset.getUntrackedParameter<bool>("useHumanReadableHistTitles", false);
107 
108  triggerNameRecoMonPath = pset.getUntrackedParameter<std::string>("triggerNameRecoMonPath", "HLT_MinBias");
109  processNameRecoMonPath = pset.getUntrackedParameter<std::string>("processNameRecoMonPath", "HLT");
110 
111  recoElectronsInput = consumes<reco::GsfElectronCollection>(
112  pset.getUntrackedParameter<edm::InputTag>("recoElectrons", edm::InputTag("gsfElectrons")));
113  recoObjectsEBT = consumes<std::vector<reco::SuperCluster>>(edm::InputTag("correctedHybridSuperClusters"));
115  consumes<std::vector<reco::SuperCluster>>(edm::InputTag("correctedMulti5x5SuperClustersWithPreshower"));
116  hltResultsT = consumes<edm::TriggerResults>(edm::InputTag("TriggerResults", "", processNameRecoMonPath));
117  triggerObjT = consumes<trigger::TriggerEventWithRefs>(edm::InputTag("hltTriggerSummaryRAW"));
118 
119  // preselction cuts
120  // recocutCollection_= pset.getParameter<edm::InputTag>("cutcollection");
121  recocut_ = pset.getParameter<int>("cutnum");
122 
123  // prescale = 10;
124  eventnum = 0;
125 
126  // just init
127  isHltConfigInitialized_ = false;
128 
130  // Read in the Vector of Parameter Sets. //
131  // Information for each filter-step //
133  std::vector<edm::ParameterSet> filters = pset.getParameter<std::vector<edm::ParameterSet>>("filters");
134 
135  int i = 0;
136  for (std::vector<edm::ParameterSet>::iterator filterconf = filters.begin(); filterconf != filters.end();
137  filterconf++) {
138  theHLTCollectionLabels.push_back(filterconf->getParameter<edm::InputTag>("HLTCollectionLabels"));
139  theHLTOutputTypes.push_back(filterconf->getParameter<int>("theHLTOutputTypes"));
140  // Grab the human-readable name, if it is not specified, use the Collection
141  // Label
142  theHLTCollectionHumanNames.push_back(
143  filterconf->getUntrackedParameter<std::string>("HLTCollectionHumanName", theHLTCollectionLabels[i].label()));
144 
145  std::vector<double> bounds = filterconf->getParameter<std::vector<double>>("PlotBounds");
146  // If the size of plot "bounds" vector != 2, abort
147  assert(bounds.size() == 2);
148  plotBounds.push_back(std::pair<double, double>(bounds[0], bounds[1]));
149  isoNames.push_back(filterconf->getParameter<std::vector<edm::InputTag>>("IsoCollections"));
150 
151  for (unsigned int i = 0; i < isoNames.back().size(); i++) {
152  switch (theHLTOutputTypes.back()) {
156  isoNames.back()[i]));
157  break;
158  case trigger::TriggerL1IsoEG: // Isolated Level 1
159  histoFillerL1Iso->isoNameTokens_.push_back(
161  isoNames.back()[i]));
162  break;
163  case trigger::TriggerPhoton: // Photon
164  histoFillerPho->isoNameTokens_.push_back(
166  isoNames.back()[i]));
167  break;
168  case trigger::TriggerElectron: // Electron
169  histoFillerEle->isoNameTokens_.push_back(
171  break;
172  case trigger::TriggerCluster: // TriggerCluster
173  histoFillerClu->isoNameTokens_.push_back(
175  isoNames.back()[i]));
176  break;
177  default:
178  throw(cms::Exception("Release Validation Error") << "HLT output type not implemented: theHLTOutputTypes[n]");
179  }
180  }
181 
182  // If the size of the isoNames vector is not greater than zero, abort
183  assert(!isoNames.back().empty());
184  if (isoNames.back().at(0).label() == "none") {
185  plotiso.push_back(false);
186  } else {
187  plotiso.push_back(true);
188  }
189  i++;
190  } // END of loop over parameter sets
191 
192  // Record number of HLTCollectionLabels
194 }
195 
199 void EmDQMReco::dqmBeginRun(const edm::Run &iRun, const edm::EventSetup &iSetup) {
200  bool isHltConfigChanged = false; // change of cfg at run boundaries?
201  isHltConfigInitialized_ = hltConfig_.init(iRun, iSetup, "HLT", isHltConfigChanged);
202 }
203 
205 // book DQM histograms //
207 void EmDQMReco::bookHistograms(DQMStore::IBooker &iBooker, edm::Run const &iRun, edm::EventSetup const &iSetup) {
208  // edm::Service<TFileService> fs;
209  iBooker.setCurrentFolder(dirname_);
210 
212  // Set up Histogram of Effiency vs Step. //
213  // theHLTCollectionLabels is a vector of InputTags //
214  // from the configuration file. //
216 
217  std::string histName = "total_eff";
218  std::string histTitle = "total events passing";
219  // This plot will have bins equal to 2+(number of
220  // HLTCollectionLabels in the config file)
221  totalreco = iBooker.book1D(
222  histName.c_str(), histTitle.c_str(), numOfHLTCollectionLabels + 2, 0, numOfHLTCollectionLabels + 2);
225  for (unsigned int u = 0; u < numOfHLTCollectionLabels; u++) {
227  }
228 
229  histName = "total_eff_RECO_matched";
230  histTitle = "total events passing (Reco matched)";
231  totalmatchreco = iBooker.book1D(
232  histName.c_str(), histTitle.c_str(), numOfHLTCollectionLabels + 2, 0, numOfHLTCollectionLabels + 2);
233  totalmatchreco->setBinLabel(numOfHLTCollectionLabels + 1, "Total");
234  totalmatchreco->setBinLabel(numOfHLTCollectionLabels + 2, "Reco");
235  for (unsigned int u = 0; u < numOfHLTCollectionLabels; u++) {
237  }
238 
239  // MonitorElement* tmphisto;
240  MonitorElement *tmpiso;
241 
243  // Set up generator-level histograms //
245  std::string pdgIdString;
246  switch (pdgGen) {
247  case 11:
248  pdgIdString = "Electron";
249  break;
250  case 22:
251  pdgIdString = "Photon";
252  break;
253  default:
254  pdgIdString = "Particle";
255  }
256 
257  //--------------------
258 
259  // reco
260  // (note that reset(..) must be used to set the value of the scoped_ptr...)
261  histReco = std::make_unique<FourVectorMonitorElements>(this,
262  iBooker,
263  "reco_%s", // pattern for histogram name
264  "%s of " + pdgIdString + "s");
265 
266  //--------------------
267 
268  // monpath
269  histRecoMonpath = std::make_unique<FourVectorMonitorElements>(this,
270  iBooker,
271  "reco_%s_monpath", // pattern for histogram name
272  "%s of " + pdgIdString + "s monpath");
273 
274  //--------------------
275 
276  // TODO: WHAT ARE THESE HISTOGRAMS FOR ? THEY SEEM NEVER REFERENCED ANYWHERE
277  // IN THIS FILE... final X monpath
278  histMonpath = std::make_unique<FourVectorMonitorElements>(this,
279  iBooker,
280  "final_%s_monpath", // pattern for histogram name
281  "Final %s Monpath");
282 
283  //--------------------
284 
286  // Set up histograms of HLT objects //
288 
289  // Determine what strings to use for histogram titles
290  std::vector<std::string> HltHistTitle;
291  if (theHLTCollectionHumanNames.size() == numOfHLTCollectionLabels && useHumanReadableHistTitles) {
292  HltHistTitle = theHLTCollectionHumanNames;
293  } else {
294  for (unsigned int i = 0; i < numOfHLTCollectionLabels; i++) {
295  HltHistTitle.push_back(theHLTCollectionLabels[i].label());
296  }
297  }
298 
299  for (unsigned int i = 0; i < numOfHLTCollectionLabels; i++) {
300  //--------------------
301  // distributions of HLT objects passing filter i
302  //--------------------
303 
304  // // Et
305  // histName = theHLTCollectionLabels[i].label()+"et_all";
306  // histTitle = HltHistTitle[i]+" Et (ALL)";
307  // tmphisto =
308  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,plotPtMin,plotPtMax);
309  // ethist.push_back(tmphisto);
310  //
311  // // Eta
312  // histName = theHLTCollectionLabels[i].label()+"eta_all";
313  // histTitle = HltHistTitle[i]+" #eta (ALL)";
314  // tmphisto =
315  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotEtaMax,plotEtaMax);
316  // etahist.push_back(tmphisto);
317  //
318  // // phi
319  // histName = theHLTCollectionLabels[i].label()+"phi_all";
320  // histTitle = HltHistTitle[i]+" #phi (ALL)";
321  // tmphisto =
322  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotPhiMax,plotPhiMax);
323  // phiHist.push_back(tmphisto);
324 
325  standardHist.push_back(std::make_unique<FourVectorMonitorElements>(
326  this,
327  iBooker,
328  theHLTCollectionLabels[i].label() + "%s_all", // histogram name
329  HltHistTitle[i] + " %s (ALL)" // histogram title
330  ));
331 
332  //--------------------
333  // distributions of reco object matching HLT object passing filter i
334  //--------------------
335 
336  // Et
337  // histName = theHLTCollectionLabels[i].label()+"et_RECO_matched";
338  // histTitle = HltHistTitle[i]+" Et (RECO matched)";
339  // tmphisto =
340  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,plotPtMin,plotPtMax);
341  // ethistmatchreco.push_back(tmphisto);
342 
343  // // Eta
344  // histName = theHLTCollectionLabels[i].label()+"eta_RECO_matched";
345  // histTitle = HltHistTitle[i]+" #eta (RECO matched)";
346  // tmphisto =
347  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotEtaMax,plotEtaMax);
348  // etahistmatchreco.push_back(tmphisto);
349  //
350  // // phi
351  // histName = theHLTCollectionLabels[i].label()+"phi_RECO_matched";
352  // histTitle = HltHistTitle[i]+" #phi (RECO matched)";
353  // tmphisto =
354  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotPhiMax,plotPhiMax);
355  // phiHistMatchReco.push_back(tmphisto);
356  histMatchReco.push_back(std::make_unique<FourVectorMonitorElements>(
357  this,
358  iBooker,
359  theHLTCollectionLabels[i].label() + "%s_RECO_matched", // histogram name
360  HltHistTitle[i] + " %s (RECO matched)" // histogram title
361  ));
362 
363  //--------------------
364  // distributions of reco object matching HLT object passing filter i
365  //--------------------
366 
367  // // Et
368  // histName =
369  // theHLTCollectionLabels[i].label()+"et_RECO_matched_monpath"; histTitle
370  // = HltHistTitle[i]+" Et (RECO matched, monpath)"; tmphisto =
371  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,plotPtMin,plotPtMax);
372  // ethistmatchrecomonpath.push_back(tmphisto);
373  //
374  // // Eta
375  // histName =
376  // theHLTCollectionLabels[i].label()+"eta_RECO_matched_monpath";
377  // histTitle = HltHistTitle[i]+" #eta (RECO matched, monpath)";
378  // tmphisto =
379  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotEtaMax,plotEtaMax);
380  // etahistmatchrecomonpath.push_back(tmphisto);
381  //
382  // // phi
383  // histName =
384  // theHLTCollectionLabels[i].label()+"phi_RECO_matched_monpath";
385  // histTitle = HltHistTitle[i]+" #phi (RECO matched, monpath)";
386  // tmphisto =
387  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotPhiMax,plotPhiMax);
388  // phiHistMatchRecoMonPath.push_back(tmphisto);
389 
390  histMatchRecoMonPath.push_back(std::make_unique<FourVectorMonitorElements>(
391  this,
392  iBooker,
393  theHLTCollectionLabels[i].label() + "%s_RECO_matched_monpath", // histogram name
394  HltHistTitle[i] + " %s (RECO matched, monpath)" // histogram title
395  ));
396  //--------------------
397  // distributions of HLT object that is closest delta-R match to sorted reco
398  // particle(s)
399  //--------------------
400 
401  // Et
402  // histName = theHLTCollectionLabels[i].label()+"et_reco";
403  // histTitle = HltHistTitle[i]+" Et (reco)";
404  // tmphisto =
405  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,plotPtMin,plotPtMax);
406  // histEtOfHltObjMatchToReco.push_back(tmphisto);
407  //
408  // // eta
409  // histName = theHLTCollectionLabels[i].label()+"eta_reco";
410  // histTitle = HltHistTitle[i]+" eta (reco)";
411  // tmphisto =
412  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotEtaMax,plotEtaMax);
413  // histEtaOfHltObjMatchToReco.push_back(tmphisto);
414  //
415  // // phi
416  // histName = theHLTCollectionLabels[i].label()+"phi_reco";
417  // histTitle = HltHistTitle[i]+" phi (reco)";
418  // tmphisto =
419  // iBooker.book1D(histName.c_str(),histTitle.c_str(),plotBins,-plotPhiMax,plotPhiMax);
420  // histPhiOfHltObjMatchToReco.push_back(tmphisto);
421 
422  histHltObjMatchToReco.push_back(std::make_unique<FourVectorMonitorElements>(
423  this,
424  iBooker,
425  theHLTCollectionLabels[i].label() + "%s_reco", // histogram name
426  HltHistTitle[i] + " %s (reco)" // histogram title
427  ));
428 
429  //--------------------
430 
431  if (!plotiso[i]) {
432  tmpiso = nullptr;
433  etahistiso.push_back(tmpiso);
434  ethistiso.push_back(tmpiso);
435  phiHistIso.push_back(tmpiso);
436 
437  etahistisomatchreco.push_back(tmpiso);
438  ethistisomatchreco.push_back(tmpiso);
439  phiHistIsoMatchReco.push_back(tmpiso);
440 
441  histEtaIsoOfHltObjMatchToReco.push_back(tmpiso);
442  histEtIsoOfHltObjMatchToReco.push_back(tmpiso);
443  histPhiIsoOfHltObjMatchToReco.push_back(tmpiso);
444 
445  } else {
446  //--------------------
447  // 2D plot: Isolation values vs X for all objects
448  //--------------------
449 
450  // X = eta
451  histName = theHLTCollectionLabels[i].label() + "eta_isolation_all";
452  histTitle = HltHistTitle[i] + " isolation vs #eta (all)";
453  tmpiso = iBooker.book2D(histName.c_str(),
454  histTitle.c_str(),
455  plotBins,
456  -plotEtaMax,
457  plotEtaMax,
458  plotBins,
459  plotBounds[i].first,
460  plotBounds[i].second);
461  etahistiso.push_back(tmpiso);
462 
463  // X = et
464  histName = theHLTCollectionLabels[i].label() + "et_isolation_all";
465  histTitle = HltHistTitle[i] + " isolation vs Et (all)";
466  tmpiso = iBooker.book2D(histName.c_str(),
467  histTitle.c_str(),
468  plotBins,
469  plotPtMin,
470  plotPtMax,
471  plotBins,
472  plotBounds[i].first,
473  plotBounds[i].second);
474  ethistiso.push_back(tmpiso);
475 
476  // X = phi
477  histName = theHLTCollectionLabels[i].label() + "phi_isolation_all";
478  histTitle = HltHistTitle[i] + " isolation vs #phi (all)";
479  tmpiso = iBooker.book2D(histName.c_str(),
480  histTitle.c_str(),
481  plotBins,
482  -plotPhiMax,
483  plotPhiMax,
484  plotBins,
485  plotBounds[i].first,
486  plotBounds[i].second);
487  phiHistIso.push_back(tmpiso);
488 
489  //--------------------
490  // 2D plot: Isolation values vs X for reco matched objects
491  //--------------------
492 
493  // X = eta
494  histName = theHLTCollectionLabels[i].label() + "eta_isolation_RECO_matched";
495  histTitle = HltHistTitle[i] + " isolation vs #eta (reco matched)";
496  tmpiso = iBooker.book2D(histName.c_str(),
497  histTitle.c_str(),
498  plotBins,
499  -plotEtaMax,
500  plotEtaMax,
501  plotBins,
502  plotBounds[i].first,
503  plotBounds[i].second);
504  etahistisomatchreco.push_back(tmpiso);
505 
506  // X = et
507  histName = theHLTCollectionLabels[i].label() + "et_isolation_RECO_matched";
508  histTitle = HltHistTitle[i] + " isolation vs Et (reco matched)";
509  tmpiso = iBooker.book2D(histName.c_str(),
510  histTitle.c_str(),
511  plotBins,
512  plotPtMin,
513  plotPtMax,
514  plotBins,
515  plotBounds[i].first,
516  plotBounds[i].second);
517  ethistisomatchreco.push_back(tmpiso);
518 
519  // X = eta
520  histName = theHLTCollectionLabels[i].label() + "phi_isolation_RECO_matched";
521  histTitle = HltHistTitle[i] + " isolation vs #phi (reco matched)";
522  tmpiso = iBooker.book2D(histName.c_str(),
523  histTitle.c_str(),
524  plotBins,
525  -plotPhiMax,
526  plotPhiMax,
527  plotBins,
528  plotBounds[i].first,
529  plotBounds[i].second);
530  phiHistIsoMatchReco.push_back(tmpiso);
531 
532  //--------------------
533  // 2D plot: Isolation values vs X for HLT object that
534  // is closest delta-R match to sorted reco particle(s)
535  //--------------------
536 
537  // X = eta
538  histName = theHLTCollectionLabels[i].label() + "eta_isolation_reco";
539  histTitle = HltHistTitle[i] + " isolation vs #eta (reco)";
540  tmpiso = iBooker.book2D(histName.c_str(),
541  histTitle.c_str(),
542  plotBins,
543  -plotEtaMax,
544  plotEtaMax,
545  plotBins,
546  plotBounds[i].first,
547  plotBounds[i].second);
548  histEtaIsoOfHltObjMatchToReco.push_back(tmpiso);
549 
550  // X = et
551  histName = theHLTCollectionLabels[i].label() + "et_isolation_reco";
552  histTitle = HltHistTitle[i] + " isolation vs Et (reco)";
553  tmpiso = iBooker.book2D(histName.c_str(),
554  histTitle.c_str(),
555  plotBins,
556  plotPtMin,
557  plotPtMax,
558  plotBins,
559  plotBounds[i].first,
560  plotBounds[i].second);
561  histEtIsoOfHltObjMatchToReco.push_back(tmpiso);
562 
563  // X = phi
564  histName = theHLTCollectionLabels[i].label() + "phi_isolation_reco";
565  histTitle = HltHistTitle[i] + " isolation vs #phi (reco)";
566  tmpiso = iBooker.book2D(histName.c_str(),
567  histTitle.c_str(),
568  plotBins,
569  -plotPhiMax,
570  plotPhiMax,
571  plotBins,
572  plotBounds[i].first,
573  plotBounds[i].second);
574  histPhiIsoOfHltObjMatchToReco.push_back(tmpiso);
575  //--------------------
576 
577  } // END of HLT histograms
578  }
579 }
580 
582 // Destructor //
585 
587 // method called to for each event //
590  // protect from hlt config failure
592  return;
593 
594  eventnum++;
595  bool plotMonpath = false;
596  bool plotReco = true;
597 
601 
602  if (pdgGen == 11) {
603  event.getByToken(recoElectronsInput, recoObjects);
604 
605  if (recoObjects->size() < (unsigned int)recocut_) {
606  // edm::LogWarning("EmDQMReco") << "Less than "<< recocut_ <<" Reco
607  // particles with pdgId=" << pdgGen << ". Only " <<
608  // cutRecoCounter->size() << " particles.";
609  return;
610  }
611  } else if (pdgGen == 22) {
612  event.getByToken(recoObjectsEBT, recoObjectsEB);
613  event.getByToken(recoObjectsEET, recoObjectsEE);
614 
615  if (recoObjectsEB->size() + recoObjectsEE->size() < (unsigned int)recocut_) {
616  // edm::LogWarning("EmDQMReco") << "Less than "<< recocut_ <<" Reco
617  // particles with pdgId=" << pdgGen << ". Only " << cutRecoCounter.size()
618  // << " particles.";
619  return;
620  }
621  }
622 
624  event.getByToken(hltResultsT, HLTR);
625 
630 
631  /* if (theHLTCollectionHumanNames[0] == "hltL1sRelaxedSingleEgammaEt8"){
632  triggerIndex = hltConfig.triggerIndex("HLT_L1SingleEG8");
633  } else if (theHLTCollectionHumanNames[0] == "hltL1sRelaxedSingleEgammaEt5") {
634  triggerIndex = hltConfig.triggerIndex("HLT_L1SingleEG5");
635  } else if (theHLTCollectionHumanNames[0] == "hltL1sRelaxedDoubleEgammaEt5") {
636  triggerIndex = hltConfig.triggerIndex("HLT_L1DoubleEG5");
637  } else {
638  triggerIndex = hltConfig.triggerIndex("");
639  } */
640 
641  unsigned int triggerIndex;
643 
644  // triggerIndex must be less than the size of HLTR or you get a CMSException
645  bool isFired = false;
646  if (triggerIndex < HLTR->size()) {
647  isFired = HLTR->accept(triggerIndex);
648  }
649 
650  // fill L1 and HLT info
651  // get objects possed by each filter
653  event.getByToken(triggerObjT, triggerObj);
654 
655  if (!triggerObj.isValid()) {
656  edm::LogWarning("EmDQMReco") << "RAW-type HLT results not found, skipping event";
657  return;
658  }
659 
661  // Fill the bin labeled "Total" //
662  // This will be the number of events looked at. //
666 
668  // Fill the bin labeled "Total" //
669  // This will be the number of events looked at. //
671  // total->Fill(numOfHLTCollectionLabels+0.5);
672  // totalmatch->Fill(numOfHLTCollectionLabels+0.5);
673 
675  // Fill reconstruction info //
677  // the recocut_ highest Et generator objects of the preselected type are our
678  // matches
679 
680  std::vector<reco::Particle> sortedReco;
681  if (plotReco == true) {
682  if (pdgGen == 11) {
683  for (edm::View<reco::Candidate>::const_iterator recopart = recoObjects->begin(); recopart != recoObjects->end();
684  recopart++) {
685  reco::Particle tmpcand(
686  recopart->charge(), recopart->p4(), recopart->vertex(), recopart->pdgId(), recopart->status());
687  sortedReco.push_back(tmpcand);
688  }
689  } else if (pdgGen == 22) {
690  for (std::vector<reco::SuperCluster>::const_iterator recopart2 = recoObjectsEB->begin();
691  recopart2 != recoObjectsEB->end();
692  recopart2++) {
693  float en = recopart2->energy();
694  float er = sqrt(pow(recopart2->x(), 2) + pow(recopart2->y(), 2) + pow(recopart2->z(), 2));
695  float px = recopart2->energy() * recopart2->x() / er;
696  float py = recopart2->energy() * recopart2->y() / er;
697  float pz = recopart2->energy() * recopart2->z() / er;
698  reco::Candidate::LorentzVector thisLV(px, py, pz, en);
699  reco::Particle tmpcand(0, thisLV, math::XYZPoint(0., 0., 0.), 22, 1);
700  sortedReco.push_back(tmpcand);
701  }
702  for (std::vector<reco::SuperCluster>::const_iterator recopart2 = recoObjectsEE->begin();
703  recopart2 != recoObjectsEE->end();
704  recopart2++) {
705  float en = recopart2->energy();
706  float er = sqrt(pow(recopart2->x(), 2) + pow(recopart2->y(), 2) + pow(recopart2->z(), 2));
707  float px = recopart2->energy() * recopart2->x() / er;
708  float py = recopart2->energy() * recopart2->y() / er;
709  float pz = recopart2->energy() * recopart2->z() / er;
710  reco::Candidate::LorentzVector thisLV(px, py, pz, en);
711  reco::Particle tmpcand(0, thisLV, math::XYZPoint(0., 0., 0.), 22, 1);
712  sortedReco.push_back(tmpcand);
713  }
714  }
715 
716  std::sort(sortedReco.begin(), sortedReco.end(), pTComparator_);
717 
718  // Now the collection of gen particles is sorted by pt.
719  // So, remove all particles from the collection so that we
720  // only have the top "1 thru recocut_" particles in it
721 
722  sortedReco.erase(sortedReco.begin() + recocut_, sortedReco.end());
723 
724  for (unsigned int i = 0; i < recocut_; i++) {
725  // validity has been implicitily checked by the cut on recocut_ above
726  histReco->fill(sortedReco[i].p4());
727 
728  // etreco ->Fill( sortedReco[i].et() );
729  // etareco->Fill( sortedReco[i].eta() );
730  // phiReco->Fill( sortedReco[i].phi() );
731 
732  if (isFired) {
733  histRecoMonpath->fill(sortedReco[i].p4());
734  plotMonpath = true;
735  }
736 
737  } // END of loop over Reconstructed particles
738 
739  if (recocut_ >= reqNum)
740  totalreco->Fill(numOfHLTCollectionLabels + 1.5); // this isn't really needed anymore keep for backward comp.
741  if (recocut_ >= reqNum)
742  totalmatchreco->Fill(numOfHLTCollectionLabels + 1.5); // this isn't really needed anymore keep for backward comp.
743  }
744 
746  // Loop over filter modules //
748  for (unsigned int n = 0; n < numOfHLTCollectionLabels; n++) {
749  // These numbers are from the Parameter Set, such as:
750  // theHLTOutputTypes = cms.uint32(100)
751  switch (theHLTOutputTypes[n]) {
752  case trigger::TriggerL1NoIsoEG: // Non-isolated Level 1
753  histoFillerL1NonIso->fillHistos(triggerObj, event, n, sortedReco, plotReco, plotMonpath);
754  break;
755  case trigger::TriggerL1IsoEG: // Isolated Level 1
756  histoFillerL1Iso->fillHistos(triggerObj, event, n, sortedReco, plotReco, plotMonpath);
757  break;
758  case trigger::TriggerPhoton: // Photon
759  histoFillerPho->fillHistos(triggerObj, event, n, sortedReco, plotReco, plotMonpath);
760  break;
761  case trigger::TriggerElectron: // Electron
762  histoFillerEle->fillHistos(triggerObj, event, n, sortedReco, plotReco, plotMonpath);
763  break;
764  case trigger::TriggerCluster: // TriggerCluster
765  histoFillerClu->fillHistos(triggerObj, event, n, sortedReco, plotReco, plotMonpath);
766  break;
767  default:
768  throw(cms::Exception("Release Validation Error") << "HLT output type not implemented: theHLTOutputTypes[n]");
769  }
770  } // END of loop over filter modules
771 }
772 
774 // fillHistos //
775 // Called by analyze method. //
777 template <class T>
779  const edm::Event &iEvent,
780  unsigned int n,
781  std::vector<reco::Particle> &sortedReco,
782  bool plotReco,
783  bool plotMonpath) {
784  std::vector<edm::Ref<T>> recoecalcands;
785  if ((triggerObj->filterIndex(dqm->theHLTCollectionLabels[n]) >= triggerObj->size())) { // only process if available
786  return;
787  }
788 
790  // Retrieve saved filter objects //
792  triggerObj->getObjects(
793  triggerObj->filterIndex(dqm->theHLTCollectionLabels[n]), dqm->theHLTOutputTypes[n], recoecalcands);
794  // Danger: special case, L1 non-isolated
795  // needs to be merged with L1 iso
797  std::vector<edm::Ref<T>> isocands;
798  triggerObj->getObjects(triggerObj->filterIndex(dqm->theHLTCollectionLabels[n]), trigger::TriggerL1IsoEG, isocands);
799  if (!isocands.empty()) {
800  for (unsigned int i = 0; i < isocands.size(); i++)
801  recoecalcands.push_back(isocands[i]);
802  }
803  } // END of if theHLTOutputTypes == 82
804 
805  if (recoecalcands.empty()) { // stop if no object passed the previous filter
806  return;
807  }
808 
809  if (recoecalcands.size() >= dqm->reqNum)
810  dqm->totalreco->Fill(n + 0.5);
811 
813  // check for validity //
814  // prevents crash in CMSSW_3_1_0_pre6 //
816  for (unsigned int j = 0; j < recoecalcands.size(); j++) {
817  if (!(recoecalcands.at(j).isAvailable())) {
818  edm::LogError("EmDQMReco") << "Event content inconsistent: TriggerEventWithRefs contains "
819  "invalid Refs"
820  << std::endl
821  << "invalid refs for: " << dqm->theHLTCollectionLabels[n].label();
822  return;
823  }
824  }
825 
827  // Loop over all HLT objects in this filter step, and //
828  // fill histograms. //
830  // bool foundAllMatches = false;
831  // unsigned int numOfHLTobjectsMatched = 0;
832  for (unsigned int i = 0; i < recoecalcands.size(); i++) {
833  dqm->standardHist[n]->fill(recoecalcands[i]->p4());
834 
836  // Plot isolation variables (show the not-yet-cut //
837  // isolation, i.e. associated to next filter) //
839  if (n + 1 < dqm->numOfHLTCollectionLabels) { // can't plot beyond last
840  if (dqm->plotiso[n + 1]) {
841  for (unsigned int j = 0; j < isoNameTokens_.size(); j++) {
843  iEvent.getByToken(isoNameTokens_.at(j), depMap);
844  if (depMap.isValid()) { // Map may not exist if only one candidate
845  // passes a double filter
846  typename edm::AssociationMap<edm::OneToValue<T, float>>::const_iterator mapi =
847  depMap->find(recoecalcands[i]);
848  if (mapi != depMap->end()) { // found candidate in isolation map!
849  dqm->etahistiso[n + 1]->Fill(recoecalcands[i]->eta(), mapi->val);
850  dqm->ethistiso[n + 1]->Fill(recoecalcands[i]->et(), mapi->val);
851  dqm->phiHistIso[n + 1]->Fill(recoecalcands[i]->phi(), mapi->val);
852  }
853  }
854  }
855  }
856  } // END of if n+1 < then the number of hlt collections
857  }
858 
860  // Loop over the Reconstructed Particles, and find the //
861  // closest HLT object match. //
863  if (plotReco == true) {
864  for (unsigned int i = 0; i < dqm->recocut_; i++) {
865  math::XYZVector currentRecoParticleMomentum = sortedReco[i].momentum();
866 
867  // float closestRecoDeltaR = 0.5;
868  float closestRecoDeltaR = 1000.;
869  int closestRecoEcalCandIndex = -1;
870  for (unsigned int j = 0; j < recoecalcands.size(); j++) {
871  float deltaR = DeltaR(recoecalcands[j]->momentum(), currentRecoParticleMomentum);
872 
873  if (deltaR < closestRecoDeltaR) {
874  closestRecoDeltaR = deltaR;
875  closestRecoEcalCandIndex = j;
876  }
877  }
878 
879  // If an HLT object was found within some delta-R
880  // of this reco particle, store it in a histogram
881  if (closestRecoEcalCandIndex >= 0) {
882  // histEtOfHltObjMatchToReco[n] ->Fill(
883  // recoecalcands[closestRecoEcalCandIndex]->et() );
884  // histEtaOfHltObjMatchToReco[n]->Fill(
885  // recoecalcands[closestRecoEcalCandIndex]->eta() );
886  // histPhiOfHltObjMatchToReco[n]->Fill(
887  // recoecalcands[closestRecoEcalCandIndex]->phi() );
888 
889  dqm->histHltObjMatchToReco[n]->fill(recoecalcands[closestRecoEcalCandIndex]->p4());
890 
891  // Also store isolation info
892  if (n + 1 < dqm->numOfHLTCollectionLabels) { // can't plot beyond last
893  if (dqm->plotiso[n + 1]) { // only plot if requested in config
894  for (unsigned int j = 0; j < isoNameTokens_.size(); j++) {
896  iEvent.getByToken(isoNameTokens_.at(j), depMap);
897  if (depMap.isValid()) { // Map may not exist if only one candidate
898  // passes a double filter
899  typename edm::AssociationMap<edm::OneToValue<T, float>>::const_iterator mapi =
900  depMap->find(recoecalcands[closestRecoEcalCandIndex]);
901  if (mapi != depMap->end()) { // found candidate in isolation map!
902  dqm->histEtaIsoOfHltObjMatchToReco[n + 1]->Fill(recoecalcands[closestRecoEcalCandIndex]->eta(),
903  mapi->val);
904  dqm->histEtIsoOfHltObjMatchToReco[n + 1]->Fill(recoecalcands[closestRecoEcalCandIndex]->et(),
905  mapi->val);
906  dqm->histPhiIsoOfHltObjMatchToReco[n + 1]->Fill(recoecalcands[closestRecoEcalCandIndex]->phi(),
907  mapi->val);
908  }
909  }
910  }
911  }
912  }
913  } // END of if closestEcalCandIndex >= 0
914  }
915 
917  // Fill reco matched objects into histograms //
919  unsigned int mtachedRecoParts = 0;
920  float minrecodist = 0.3;
921  if (n == 0)
922  minrecodist = 0.5; // low L1-resolution => allow wider matching
923  for (unsigned int i = 0; i < dqm->recocut_; i++) {
924  // match generator candidate
925  bool matchThis = false;
926  math::XYZVector candDir = sortedReco[i].momentum();
927  unsigned int closest = 0;
928  double closestDr = 1000.;
929  for (unsigned int trigOb = 0; trigOb < recoecalcands.size(); trigOb++) {
930  double dr = DeltaR(recoecalcands[trigOb]->momentum(), candDir);
931  if (dr < closestDr) {
932  closestDr = dr;
933  closest = trigOb;
934  }
935  if (closestDr > minrecodist) { // it's not really a "match" if it's that far away
936  closest = -1;
937  } else {
938  mtachedRecoParts++;
939  matchThis = true;
940  }
941  }
942  if (!matchThis)
943  continue; // only plot matched candidates
944  // fill coordinates of mc particle matching trigger object
945  // ethistmatchreco[n] ->Fill( sortedReco[i].et() );
946  // etahistmatchreco[n]->Fill( sortedReco[i].eta() );
947  // phiHistMatchReco[n]->Fill( sortedReco[i].phi() );
948  dqm->histMatchReco[n]->fill(sortedReco[i].p4());
949 
950  if (plotMonpath) {
951  // ethistmatchrecomonpath[n]->Fill( sortedReco[i].et() );
952  // etahistmatchrecomonpath[n]->Fill( sortedReco[i].eta() );
953  // phiHistMatchRecoMonPath[n]->Fill( sortedReco[i].phi() );
954  dqm->histMatchRecoMonPath[n]->fill(sortedReco[i].p4());
955  }
957  // Plot isolation variables (show the not-yet-cut //
958  // isolation, i.e. associated to next filter) //
960  if (n + 1 < dqm->numOfHLTCollectionLabels) { // can't plot beyond last
961  if (dqm->plotiso[n + 1]) { // only plot if requested in config
962  for (unsigned int j = 0; j < isoNameTokens_.size(); j++) {
964  iEvent.getByToken(isoNameTokens_.at(j), depMapReco);
965  if (depMapReco.isValid()) { // Map may not exist if only one
966  // candidate passes a double filter
967  typename edm::AssociationMap<edm::OneToValue<T, float>>::const_iterator mapi =
968  depMapReco->find(recoecalcands[closest]);
969  if (mapi != depMapReco->end()) { // found candidate in isolation map!
970  dqm->etahistisomatchreco[n + 1]->Fill(sortedReco[i].eta(), mapi->val);
971  dqm->ethistisomatchreco[n + 1]->Fill(sortedReco[i].et(), mapi->val);
972  dqm->phiHistIsoMatchReco[n + 1]->Fill(sortedReco[i].eta(), mapi->val);
973  }
974  }
975  }
976  }
977  } // END of if n+1 < then the number of hlt collections
978  }
979  // fill total reco matched efficiency
980  if (mtachedRecoParts >= dqm->reqNum)
981  dqm->totalmatchreco->Fill(n + 0.5);
982  }
983 }
984 
std::vector< MonitorElement * > histEtIsoOfHltObjMatchToReco
Definition: EmDQMReco.h:192
int pdgGen
Definition: EmDQMReco.h:120
void fill(const math::XYZTLorentzVector &momentum)
Definition: EmDQMReco.cc:78
std::vector< std::unique_ptr< FourVectorMonitorElements > > histMatchReco
Definition: EmDQMReco.h:168
T getUntrackedParameter(std::string const &, T const &) const
std::vector< std::vector< edm::InputTag > > isoNames
Definition: EmDQMReco.h:109
std::vector< std::unique_ptr< FourVectorMonitorElements > > standardHist
Definition: EmDQMReco.h:163
double energy() const
energy
Definition: Particle.h:99
std::unique_ptr< FourVectorMonitorElements > histReco
Definition: EmDQMReco.h:207
std::vector< MonitorElement * > etahistisomatchreco
Definition: EmDQMReco.h:188
HistoFillerReco< reco::RecoEcalCandidateCollection > * histoFillerPho
Definition: EmDQMReco.h:229
HistoFillerReco< reco::ElectronCollection > * histoFillerEle
Definition: EmDQMReco.h:226
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:32
std::string dirname_
Definition: EmDQMReco.h:224
void dqmBeginRun(const edm::Run &, const edm::EventSetup &) override
Definition: EmDQMReco.cc:199
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
int eventnum
Definition: EmDQMReco.h:220
std::vector< MonitorElement * > etahistiso
Definition: EmDQMReco.h:184
EmDQMReco(const edm::ParameterSet &pset)
Constructor.
Definition: EmDQMReco.cc:89
std::vector< int > theHLTOutputTypes
Definition: EmDQMReco.h:107
void fillHistos(edm::Handle< trigger::TriggerEventWithRefs > &triggerObj, const edm::Event &iEvent, unsigned int n, std::vector< reco::Particle > &sortedReco, bool plotReco, bool plotMonpath)
Definition: EmDQMReco.cc:778
std::vector< edm::EDGetTokenT< edm::AssociationMap< edm::OneToValue< T, float > > > > isoNameTokens_
Definition: EmDQMReco.h:40
Log< level::Error, false > LogError
edm::EDGetTokenT< std::vector< reco::SuperCluster > > recoObjectsEET
Definition: EmDQMReco.h:150
assert(be >=bs)
unsigned int plotBins
Definition: EmDQMReco.h:130
std::vector< MonitorElement * > phiHistIsoMatchReco
Definition: EmDQMReco.h:190
std::vector< TPRegexp > filters
Definition: eve_filter.cc:22
bool isHltConfigInitialized_
Definition: EmDQMReco.h:113
double recoEtaAcc
Definition: EmDQMReco.h:121
HLTConfigProvider hltConfig_
Definition: EmDQMReco.h:112
FourVectorMonitorElements(EmDQMReco *_parent, DQMStore::IBooker &iBooker, const std::string &histogramNameTemplate, const std::string &histogramTitleTemplate)
Definition: EmDQMReco.cc:48
unsigned int triggerIndex(const std::string &triggerName) const
slot position of trigger path in trigger table (0 to size-1)
void Fill(long long x)
std::unique_ptr< FourVectorMonitorElements > histRecoMonpath
Definition: EmDQMReco.h:212
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
std::vector< std::unique_ptr< FourVectorMonitorElements > > histMatchRecoMonPath
Definition: EmDQMReco.h:173
char const * label
double recoEtAcc
Definition: EmDQMReco.h:122
std::vector< MonitorElement * > ethistiso
Definition: EmDQMReco.h:185
GreaterByPt< reco::Particle > pTComparator_
Definition: EmDQMReco.h:235
int iEvent
Definition: GenABIO.cc:224
std::vector< std::string > theHLTCollectionHumanNames
Definition: EmDQMReco.h:105
unsigned int numOfHLTCollectionLabels
Definition: EmDQMReco.h:102
HistoFillerReco< l1extra::L1EmParticleCollection > * histoFillerL1Iso
Definition: EmDQMReco.h:230
T sqrt(T t)
Definition: SSEVec.h:19
MonitorElement * totalreco
Definition: EmDQMReco.h:199
std::vector< std::unique_ptr< FourVectorMonitorElements > > histHltObjMatchToReco
Definition: EmDQMReco.h:178
EmDQMReco * dqm
Definition: EmDQMReco.h:43
MonitorElement * etaMonitorElement
Definition: EmDQMReco.h:75
edm::EDGetTokenT< std::vector< reco::SuperCluster > > recoObjectsEBT
Definition: EmDQMReco.h:149
void analyze(const edm::Event &event, const edm::EventSetup &) override
Definition: EmDQMReco.cc:589
MonitorElement * totalmatchreco
Definition: EmDQMReco.h:200
std::vector< MonitorElement * > histPhiIsoOfHltObjMatchToReco
Definition: EmDQMReco.h:194
std::vector< MonitorElement * > phiHistIso
Definition: EmDQMReco.h:186
unsigned int recocut_
Definition: EmDQMReco.h:134
bool isValid() const
Definition: HandleBase.h:70
virtual void setBinLabel(int bin, const std::string &label, int axis=1)
set bin label for x, y or z axis (axis=1, 2, 3 respectively)
std::string triggerNameRecoMonPath
Definition: EmDQMReco.h:138
std::vector< bool > plotiso
Definition: EmDQMReco.h:108
edm::EDGetTokenT< reco::GsfElectronCollection > recoElectronsInput
Definition: EmDQMReco.h:148
bool useHumanReadableHistTitles
Definition: EmDQMReco.h:104
double plotPhiMax
Definition: EmDQMReco.h:127
edm::EDGetTokenT< trigger::TriggerEventWithRefs > triggerObjT
Definition: EmDQMReco.h:152
HistoFillerReco< l1extra::L1EmParticleCollection > * histoFillerL1NonIso
Definition: EmDQMReco.h:228
MonitorElement * phiMonitorElement
Definition: EmDQMReco.h:76
unsigned int reqNum
Definition: EmDQMReco.h:119
std::vector< MonitorElement * > ethistisomatchreco
Definition: EmDQMReco.h:189
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
bool init(const edm::Run &iRun, const edm::EventSetup &iSetup, const std::string &processName, bool &changed)
d&#39;tor
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
std::vector< edm::InputTag > theHLTCollectionLabels
Definition: EmDQMReco.h:100
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:36
std::vector< std::pair< double, double > > plotBounds
Definition: EmDQMReco.h:110
edm::EDGetTokenT< edm::TriggerResults > hltResultsT
Definition: EmDQMReco.h:151
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:86
std::vector< MonitorElement * > histEtaIsoOfHltObjMatchToReco
Definition: EmDQMReco.h:193
double plotPtMin
Definition: EmDQMReco.h:125
double plotPtMax
Definition: EmDQMReco.h:126
~EmDQMReco() override
Destructor.
Definition: EmDQMReco.cc:584
Log< level::Warning, false > LogWarning
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
double plotEtaMax
Definition: EmDQMReco.h:124
std::string processNameRecoMonPath
Definition: EmDQMReco.h:143
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
Definition: EmDQMReco.cc:207
tuple size
Write out results.
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
std::unique_ptr< FourVectorMonitorElements > histMonpath
Definition: EmDQMReco.h:217
Definition: Run.h:45
HistoFillerReco< reco::RecoEcalCandidateCollection > * histoFillerClu
Definition: EmDQMReco.h:227