File indexing completed on 2025-01-22 07:34:14
0001
0002
0003
0004 #include "HLTriggerOffline/Egamma/interface/EmDQM.h"
0005
0006 using namespace ROOT::Math::VectorUtil;
0007
0008
0009
0010
0011 EmDQM::EmDQM(const edm::ParameterSet &pset) {
0012
0013
0014 autoConfMode_ = pset.getUntrackedParameter<bool>("autoConfMode", false);
0015
0016
0017 triggerObject_ = pset.getParameter<edm::InputTag>("triggerobject");
0018 verbosity_ = pset.getUntrackedParameter<unsigned int>("verbosity", 0);
0019 genEtaAcc_ = pset.getParameter<double>("genEtaAcc");
0020 genEtAcc_ = pset.getParameter<double>("genEtAcc");
0021 ptMax_ = pset.getUntrackedParameter<double>("PtMax", 1000.);
0022 ptMin_ = pset.getUntrackedParameter<double>("PtMin", 0.);
0023 etaMax_ = pset.getUntrackedParameter<double>("EtaMax", 2.7);
0024 phiMax_ = pset.getUntrackedParameter<double>("PhiMax", 3.15);
0025 nbins_ = pset.getUntrackedParameter<unsigned int>("Nbins", 40);
0026 eta2DMax_ = pset.getUntrackedParameter<double>("Eta2DMax", 2.8);
0027 phi2DMax_ = pset.getUntrackedParameter<double>("Phi2DMax", 3.2);
0028 nbins2D_ = pset.getUntrackedParameter<unsigned int>("Nbins2D", 16);
0029 minEtForEtaEffPlot_ = pset.getUntrackedParameter<unsigned int>("minEtForEtaEffPlot", 15);
0030 useHumanReadableHistTitles_ = pset.getUntrackedParameter<bool>("useHumanReadableHistTitles", false);
0031 mcMatchedOnly_ = pset.getUntrackedParameter<bool>("mcMatchedOnly", true);
0032 noPhiPlots_ = pset.getUntrackedParameter<bool>("noPhiPlots", true);
0033 noIsolationPlots_ = pset.getUntrackedParameter<bool>("noIsolationPlots", true);
0034
0035 if (!autoConfMode_) {
0036 paramSets.push_back(pset);
0037 isData_ = false;
0038 } else {
0039 isData_ = pset.getParameter<bool>("isData");
0040 }
0041
0042 histoFillerEle = new HistoFiller<reco::ElectronCollection>(this);
0043 histoFillerPho = new HistoFiller<reco::RecoEcalCandidateCollection>(this);
0044 histoFillerClu = new HistoFiller<reco::RecoEcalCandidateCollection>(this);
0045 histoFillerL1Iso = new HistoFiller<l1extra::L1EmParticleCollection>(this);
0046 histoFillerL1NonIso = new HistoFiller<l1extra::L1EmParticleCollection>(this);
0047
0048
0049 genParticles_token = consumes<edm::View<reco::Candidate>>(edm::InputTag("genParticles"));
0050 triggerObject_token = consumes<trigger::TriggerEventWithRefs>(triggerObject_);
0051 hltResults_token = consumes<edm::TriggerResults>(edm::InputTag("TriggerResults", "", triggerObject_.process()));
0052 if (autoConfMode_) {
0053 gencutColl_fidWenu_token = mayConsume<edm::View<reco::Candidate>>(edm::InputTag("fiducialWenu"));
0054 gencutColl_fidZee_token = mayConsume<edm::View<reco::Candidate>>(edm::InputTag("fiducialZee"));
0055 gencutColl_fidTripleEle_token = mayConsume<edm::View<reco::Candidate>>(edm::InputTag("fiducialTripleEle"));
0056 gencutColl_fidGammaJet_token = mayConsume<edm::View<reco::Candidate>>(edm::InputTag("fiducialGammaJet"));
0057 gencutColl_fidDiGamma_token = mayConsume<edm::View<reco::Candidate>>(edm::InputTag("fiducialDiGamma"));
0058 } else {
0059 gencutColl_manualConf_token =
0060 consumes<edm::View<reco::Candidate>>(pset.getParameter<edm::InputTag>("cutcollection"));
0061 }
0062 }
0063
0064 void EmDQM::dqmBeginRun(edm::Run const &iRun, edm::EventSetup const &iSetup) {
0065 bool changed(true);
0066 if (hltConfig_.init(iRun, iSetup, triggerObject_.process(), changed)) {
0067
0068
0069 if (autoConfMode_) {
0070 if (verbosity_ >= OUTPUT_ALL) {
0071
0072 edm::LogPrint("EmDQM") << "inited=" << hltConfig_.inited();
0073 edm::LogPrint("EmDQM") << "changed=" << hltConfig_.changed();
0074 edm::LogPrint("EmDQM") << "processName=" << hltConfig_.processName();
0075 edm::LogPrint("EmDQM") << "tableName=" << hltConfig_.tableName();
0076 edm::LogPrint("EmDQM") << "size=" << hltConfig_.size();
0077 edm::LogInfo("EmDQM") << "The following filter types are not analyzed: \n"
0078 << "\tHLTGlobalSumsMET\n"
0079 << "\tHLTHtMhtFilter\n"
0080 << "\tHLTMhtFilter\n"
0081 << "\tHLTJetTag\n"
0082 << "\tHLT1CaloJet\n"
0083 << "\tHLT1CaloMET\n"
0084 << "\tHLT1CaloBJet\n"
0085 << "\tHLT1Tau\n"
0086 << "\tHLT1PFTau\n"
0087 << "\tPFTauSelector\n"
0088 << "\tHLT1PFJet\n"
0089 << "\tHLTPFJetCollectionsFilter\n"
0090 << "\tHLTPFJetCollectionsVBFFilter\n"
0091 << "\tHLTPFJetTag\n"
0092 << "\tEtMinCaloJetSelector\n"
0093 << "\tEtMinPFJetSelector\n"
0094 << "\tLargestEtCaloJetSelector\n"
0095 << "\tLargestEtPFJetSelector\n"
0096 << "\tHLTEgammaTriggerFilterObjectWrapper\n"
0097 << "\tHLTEgammaDoubleLegCombFilter\n"
0098 << "\tHLT2ElectronTau\n"
0099 << "\tHLT2ElectronMET\n"
0100 << "\tHLT2ElectronPFTau\n"
0101 << "\tHLTPMMassFilter\n"
0102 << "\tHLTHcalTowerFilter\n"
0103 << "\tHLT1Photon\n"
0104 << "\tHLTRFilter\n"
0105 << "\tHLTRHemisphere\n"
0106 << "\tHLTElectronPFMTFilter\n"
0107 << "\tPrimaryVertexObjectFilter\n"
0108 << "\tHLTEgammaAllCombMassFilter\n"
0109 << "\tHLTMuon*\n";
0110 }
0111
0112
0113 std::vector<std::vector<std::string>> egammaPaths = findEgammaPaths();
0114
0115
0116
0117
0118
0119
0120
0121
0122
0123
0124 std::vector<std::string> filterModules;
0125
0126 for (unsigned int j = 0; j < egammaPaths.size(); j++) {
0127 if (verbosity_ >= OUTPUT_ALL) {
0128 switch (j) {
0129 case TYPE_SINGLE_ELE:
0130 edm::LogPrint("EmDQM") << "////////////////////////////////////////"
0131 "/\nSingle electron paths: ";
0132 break;
0133 case TYPE_DOUBLE_ELE:
0134 edm::LogPrint("EmDQM") << "////////////////////////////////////////"
0135 "/\nDouble electron paths: ";
0136 break;
0137 case TYPE_TRIPLE_ELE:
0138 edm::LogPrint("EmDQM") << "////////////////////////////////////////"
0139 "/\nTriple electron paths: ";
0140 break;
0141 case TYPE_SINGLE_PHOTON:
0142 edm::LogPrint("EmDQM") << "////////////////////////////////////////"
0143 "/\nSingle photon paths: ";
0144 break;
0145 case TYPE_DOUBLE_PHOTON:
0146 edm::LogPrint("EmDQM") << "////////////////////////////////////////"
0147 "/\nDouble photon paths: ";
0148 break;
0149 }
0150 }
0151
0152 for (unsigned int i = 0; i < egammaPaths.at(j).size(); i++) {
0153
0154 const std::string pathName = egammaPaths.at(j).at(i);
0155 if (verbosity_ >= OUTPUT_ALL)
0156 edm::LogPrint("EmDQM") << "Path: " << pathName;
0157
0158
0159 filterModules = getFilterModules(pathName);
0160
0161
0162 edm::ParameterSet paramSet;
0163
0164 paramSet.addUntrackedParameter("pathIndex", hltConfig_.triggerIndex(pathName));
0165 paramSet.addParameter("@module_label", hltConfig_.removeVersion(pathName) + "_DQM");
0166
0167
0168
0169
0170 double genEtMin = getPrimaryEtCut(pathName);
0171 if (genEtMin >= 0) {
0172 paramSet.addUntrackedParameter("genEtMin", genEtMin);
0173 } else {
0174 if (verbosity_ >= OUTPUT_WARNINGS)
0175 edm::LogWarning("EmDQM") << "Pathname: '" << pathName
0176 << "': Unable to determine a minimum Et. Will not include "
0177 "this path in the validation.";
0178 continue;
0179 }
0180
0181
0182
0183 double ptMax = ptMax_;
0184 double ptMin = ptMin_;
0185 if (ptMax < (1.2 * genEtMin)) {
0186 paramSet.addUntrackedParameter<double>("PtMax", (10 * ceil(0.12 * genEtMin)));
0187 paramSet.addUntrackedParameter<double>("PtMin", (10 * ceil(0.12 * genEtMin) - ptMax + ptMin));
0188 } else {
0189 paramSet.addUntrackedParameter<double>("PtMax", ptMax);
0190 paramSet.addUntrackedParameter<double>("PtMin", ptMin);
0191 }
0192
0193
0194 switch (j) {
0195 case TYPE_SINGLE_ELE:
0196 paramSet.addParameter<unsigned>("reqNum", 1);
0197 paramSet.addParameter<int>("pdgGen", 11);
0198 paramSet.addParameter<edm::InputTag>("cutcollection", edm::InputTag("fiducialWenu"));
0199 paramSet.addParameter<int>("cutnum", 1);
0200 break;
0201 case TYPE_DOUBLE_ELE:
0202 paramSet.addParameter<unsigned>("reqNum", 2);
0203 paramSet.addParameter<int>("pdgGen", 11);
0204 paramSet.addParameter<edm::InputTag>("cutcollection", edm::InputTag("fiducialZee"));
0205 paramSet.addParameter<int>("cutnum", 2);
0206 break;
0207 case TYPE_TRIPLE_ELE:
0208 paramSet.addParameter<unsigned>("reqNum", 3);
0209 paramSet.addParameter<int>("pdgGen", 11);
0210 paramSet.addParameter<edm::InputTag>("cutcollection", edm::InputTag("fiducialTripleEle"));
0211 paramSet.addParameter<int>("cutnum", 3);
0212 break;
0213 case TYPE_SINGLE_PHOTON:
0214 paramSet.addParameter<unsigned>("reqNum", 1);
0215 paramSet.addParameter<int>("pdgGen", 22);
0216 paramSet.addParameter<edm::InputTag>("cutcollection", edm::InputTag("fiducialGammaJet"));
0217 paramSet.addParameter<int>("cutnum", 1);
0218 break;
0219 case TYPE_DOUBLE_PHOTON:
0220 paramSet.addParameter<unsigned>("reqNum", 2);
0221 paramSet.addParameter<int>("pdgGen", 22);
0222 paramSet.addParameter<edm::InputTag>("cutcollection", edm::InputTag("fiducialDiGamma"));
0223 paramSet.addParameter<int>("cutnum", 2);
0224 }
0225
0226
0227
0228 if (isData_)
0229 paramSet.addParameter<edm::InputTag>("cutcollection", edm::InputTag("gsfElectrons"));
0230
0231 std::vector<edm::ParameterSet> filterVPSet;
0232 edm::ParameterSet filterPSet;
0233 std::string moduleLabel;
0234
0235
0236 for (std::vector<std::string>::iterator filter = filterModules.begin(); filter != filterModules.end();
0237 ++filter) {
0238 std::string moduleType = hltConfig_.modulePSet(*filter).getParameter<std::string>("@module_type");
0239 moduleLabel = hltConfig_.modulePSet(*filter).getParameter<std::string>("@module_label");
0240
0241
0242 if (moduleType == "Pythia6GeneratorFilter" || moduleType == "HLTTriggerTypeFilter" ||
0243 moduleType == "HLTLevel1Activity" || moduleType == "HLTPrescaler" || moduleType == "HLTBool")
0244 continue;
0245
0246
0247 if (moduleType == "HLTLevel1GTSeed") {
0248 filterPSet = makePSetForL1SeedFilter(moduleLabel);
0249 } else if (moduleType == "HLTEgammaL1MatchFilterRegional") {
0250 filterPSet = makePSetForL1SeedToSuperClusterMatchFilter(moduleLabel);
0251 } else if (moduleType == "HLTEgammaEtFilter") {
0252 filterPSet = makePSetForEtFilter(moduleLabel);
0253 } else if (moduleType == "HLTElectronOneOEMinusOneOPFilterRegional") {
0254 filterPSet = makePSetForOneOEMinusOneOPFilter(moduleLabel);
0255 } else if (moduleType == "HLTElectronPixelMatchFilter") {
0256 filterPSet = makePSetForPixelMatchFilter(moduleLabel);
0257 } else if (moduleType == "HLTEgammaGenericFilter") {
0258 filterPSet = makePSetForEgammaGenericFilter(moduleLabel);
0259 } else if (moduleType == "HLTEgammaGenericQuadraticFilter") {
0260 filterPSet = makePSetForEgammaGenericQuadraticFilter(moduleLabel);
0261 } else if (moduleType == "HLTElectronGenericFilter") {
0262 filterPSet = makePSetForElectronGenericFilter(moduleLabel);
0263 } else if (moduleType == "HLTEgammaDoubleEtDeltaPhiFilter") {
0264 filterPSet = makePSetForEgammaDoubleEtDeltaPhiFilter(moduleLabel);
0265 } else if (moduleType == "HLTGlobalSumsMET" || moduleType == "HLTHtMhtFilter" ||
0266 moduleType == "HLTMhtFilter" || moduleType == "HLTJetTag" || moduleType == "HLT1CaloJet" ||
0267 moduleType == "HLT1CaloMET" || moduleType == "HLT1CaloBJet" || moduleType == "HLT1Tau" ||
0268 moduleType == "HLT1PFTau" || moduleType == "PFTauSelector" || moduleType == "HLT1PFJet" ||
0269 moduleType == "HLTPFJetCollectionsFilter" || moduleType == "HLTPFJetCollectionsVBFFilter" ||
0270 moduleType == "HLTPFJetTag" || moduleType == "EtMinCaloJetSelector" ||
0271 moduleType == "EtMinPFJetSelector" || moduleType == "LargestEtCaloJetSelector" ||
0272 moduleType == "LargestEtPFJetSelector" ||
0273 moduleType == "HLTEgammaTriggerFilterObjectWrapper"
0274 || moduleType == "HLTEgammaDoubleLegCombFilter"
0275
0276
0277 || moduleType == "HLT2ElectronMET" || moduleType == "HLT2ElectronTau" ||
0278 moduleType == "HLT2ElectronPFTau" || moduleType == "HLTPMMassFilter" ||
0279 moduleType == "HLTHcalTowerFilter" || moduleType == "HLT1Photon" || moduleType == "HLTRFilter" ||
0280 moduleType == "HLTRHemisphere" || moduleType == "HLTElectronPFMTFilter" ||
0281 moduleType == "PrimaryVertexObjectFilter" || moduleType == "HLTEgammaAllCombMassFilter" ||
0282 moduleType.find("HLTMuon") != std::string::npos)
0283 continue;
0284 else {
0285 if (verbosity_ >= OUTPUT_WARNINGS)
0286 edm::LogWarning("EmDQM") << "No parameter set for filter '" << moduleLabel << "' with filter type '"
0287 << moduleType << "' added. Module will not be analyzed.";
0288 continue;
0289 }
0290
0291
0292 if (!filterPSet.empty()) {
0293 if (!hltConfig_.modulePSet(moduleLabel).exists("saveTags")) {
0294
0295
0296
0297 if (moduleLabel.find("Unseeded") != std::string::npos &&
0298 (j == TYPE_DOUBLE_PHOTON || j == TYPE_SINGLE_PHOTON)) {
0299 filterVPSet.back().addParameter<int>("theHLTOutputTypes", trigger::TriggerPhoton);
0300 }
0301 }
0302
0303
0304
0305 if (filterPSet.getParameter<int>("ncandcut") < 0) {
0306 edm::LogPrint("EmDQM") << "No number of candidates for filter " << moduleLabel << " found. Set to "
0307 << paramSet.getParameter<int>("cutnum") << ", determined from path name.";
0308 filterPSet.addParameter<int>("ncandcut", paramSet.getParameter<int>("cutnum"));
0309 } else if (filterPSet.getParameter<int>("ncandcut") > paramSet.getParameter<int>("cutnum")) {
0310 edm::LogInfo("EmDQM") << "Changed required number of candidates from "
0311 << paramSet.getParameter<int>("cutnum") << " to "
0312 << filterPSet.getParameter<int>("ncandcut") << " for filter " << moduleLabel;
0313 paramSet.addParameter<int>("cutnum", filterPSet.getParameter<int>("ncandcut"));
0314 paramSet.addParameter<unsigned>("reqNum", (unsigned)filterPSet.getParameter<int>("ncandcut"));
0315 }
0316
0317 filterVPSet.push_back(filterPSet);
0318 } else
0319 break;
0320
0321 }
0322
0323
0324 if (!filterPSet.empty()) {
0325 std::string lastModuleName = filterPSet.getParameter<edm::InputTag>("HLTCollectionLabels").label();
0326 if (!hltConfig_.modulePSet(lastModuleName).exists("saveTags")) {
0327
0328
0329
0330 if ((j == TYPE_SINGLE_PHOTON || j == TYPE_DOUBLE_PHOTON) && pathName.rfind("Ele") == std::string::npos) {
0331 filterVPSet.back().addParameter<int>("theHLTOutputTypes", trigger::TriggerPhoton);
0332 }
0333 }
0334 paramSet.addParameter<std::vector<edm::ParameterSet>>("filters", filterVPSet);
0335 } else {
0336 if (verbosity_ >= OUTPUT_ALL)
0337 edm::LogPrint("EmDQM") << "Will not include this path in the validation due to "
0338 "errors while generating the parameter set.";
0339 continue;
0340 }
0341
0342
0343
0344
0345 paramSets.push_back(paramSet);
0346 }
0347
0348 }
0349 }
0350
0351 hltCollectionLabelsFoundPerPath.reserve(paramSets.size());
0352 hltCollectionLabelsMissedPerPath.reserve(paramSets.size());
0353
0354
0355
0356
0357 for (std::vector<edm::ParameterSet>::iterator psetIt = paramSets.begin(); psetIt != paramSets.end(); ++psetIt) {
0358 hltCollectionLabelsFoundPerPath.push_back(hltCollectionLabelsFound);
0359 hltCollectionLabelsMissedPerPath.push_back(hltCollectionLabelsMissed);
0360 }
0361
0362 if (changed) {
0363
0364
0365 }
0366 } else {
0367
0368
0369 if (verbosity_ >= OUTPUT_ERRORS)
0370 edm::LogError("EmDQM") << " HLT config extraction failure with process name '" << triggerObject_.process()
0371 << "'.";
0372
0373 }
0374 }
0375
0376 void EmDQM::bookHistograms(DQMStore::IBooker &iBooker, edm::Run const &iRun, edm::EventSetup const &iSetup) {
0377
0378
0379
0380 for (std::vector<edm::ParameterSet>::iterator psetIt = paramSets.begin(); psetIt != paramSets.end(); ++psetIt) {
0381 SetVarsFromPSet(psetIt);
0382
0383 iBooker.setCurrentFolder(dirname_);
0384
0385
0386
0387
0388
0389
0390
0391 std::vector<MonitorElement *> etahist;
0392 std::vector<MonitorElement *> phihist;
0393 std::vector<MonitorElement *> ethist;
0394 std::vector<MonitorElement *> etahistmatch;
0395 std::vector<MonitorElement *> phihistmatch;
0396 std::vector<MonitorElement *> ethistmatch;
0397 std::vector<MonitorElement *> histEtOfHltObjMatchToGen;
0398 std::vector<MonitorElement *> histEtaOfHltObjMatchToGen;
0399 std::vector<MonitorElement *> histPhiOfHltObjMatchToGen;
0400 std::vector<MonitorElement *> etaphihist;
0401 std::vector<MonitorElement *> etaphihistmatch;
0402 std::vector<MonitorElement *> histEtaPhiOfHltObjMatchToGen;
0403
0404 MonitorElement *total;
0405 MonitorElement *totalmatch;
0406
0407 MonitorElement *etgen;
0408 MonitorElement *etagen;
0409 MonitorElement *phigen;
0410 MonitorElement *etaphigen;
0411
0412 std::string histName = "total_eff";
0413 std::string histTitle = "total events passing";
0414 if (!mcMatchedOnly_) {
0415
0416
0417 total = iBooker.book1D(
0418 histName.c_str(), histTitle.c_str(), numOfHLTCollectionLabels + 2, 0, numOfHLTCollectionLabels + 2);
0419 total->setBinLabel(numOfHLTCollectionLabels + 1, "Total");
0420 total->setBinLabel(numOfHLTCollectionLabels + 2, "Gen");
0421 for (unsigned int u = 0; u < numOfHLTCollectionLabels; u++) {
0422 total->setBinLabel(u + 1, theHLTCollectionLabels[u].label());
0423 }
0424 }
0425
0426 histName = "total_eff_MC_matched";
0427 histTitle = "total events passing (mc matched)";
0428 totalmatch = iBooker.book1D(
0429 histName.c_str(), histTitle.c_str(), numOfHLTCollectionLabels + 2, 0, numOfHLTCollectionLabels + 2);
0430 totalmatch->setBinLabel(numOfHLTCollectionLabels + 1, "Total");
0431 totalmatch->setBinLabel(numOfHLTCollectionLabels + 2, "Gen");
0432 for (unsigned int u = 0; u < numOfHLTCollectionLabels; u++) {
0433 totalmatch->setBinLabel(u + 1, theHLTCollectionLabels[u].label());
0434 }
0435
0436 MonitorElement *tmphisto;
0437
0438
0439
0440
0441
0442 std::string pdgIdString;
0443 switch (pdgGen) {
0444 case 11:
0445 pdgIdString = "Electron";
0446 break;
0447 case 22:
0448 pdgIdString = "Photon";
0449 break;
0450 default:
0451 pdgIdString = "Particle";
0452 }
0453
0454 histName = "gen_et";
0455 histTitle = "E_{T} of " + pdgIdString + "s";
0456 etgen = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, plotPtMin, plotPtMax);
0457 histName = "gen_eta";
0458 histTitle = "#eta of " + pdgIdString + "s ";
0459 etagen = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -etaMax_, etaMax_);
0460 histName = "gen_phi";
0461 histTitle = "#phi of " + pdgIdString + "s ";
0462 if (!noPhiPlots_)
0463 phigen = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -phiMax_, phiMax_);
0464 histName = "gen_etaphi";
0465 histTitle = "#eta-#phi of " + pdgIdString + "s ";
0466 etaphigen = iBooker.book2D(
0467 histName.c_str(), histTitle.c_str(), nbins2D_ - 2, -eta2DMax_, eta2DMax_, nbins2D_, -phi2DMax_, phi2DMax_);
0468
0469
0470
0471
0472 std::vector<std::string> HltHistTitle;
0473 if (theHLTCollectionHumanNames.size() == numOfHLTCollectionLabels && useHumanReadableHistTitles_) {
0474 HltHistTitle = theHLTCollectionHumanNames;
0475 } else {
0476 for (unsigned int i = 0; i < numOfHLTCollectionLabels; i++) {
0477 HltHistTitle.push_back(theHLTCollectionLabels[i].label());
0478 }
0479 }
0480
0481 for (unsigned int i = 0; i < numOfHLTCollectionLabels; i++) {
0482 if (!mcMatchedOnly_) {
0483
0484 histName = theHLTCollectionLabels[i].label() + "et_all";
0485 histTitle = HltHistTitle[i] + " Et (ALL)";
0486 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, plotPtMin, plotPtMax);
0487 ethist.push_back(tmphisto);
0488
0489
0490 histName = theHLTCollectionLabels[i].label() + "eta_all";
0491 histTitle = HltHistTitle[i] + " #eta (ALL)";
0492 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -etaMax_, etaMax_);
0493 etahist.push_back(tmphisto);
0494
0495 if (!noPhiPlots_) {
0496
0497 histName = theHLTCollectionLabels[i].label() + "phi_all";
0498 histTitle = HltHistTitle[i] + " #phi (ALL)";
0499 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -phiMax_, phiMax_);
0500 phihist.push_back(tmphisto);
0501 }
0502
0503 histName = theHLTCollectionLabels[i].label() + "etaphi_all";
0504 histTitle = HltHistTitle[i] + " #eta-#phi (ALL)";
0505 tmphisto = iBooker.book2D(
0506 histName.c_str(), histTitle.c_str(), nbins2D_ - 2, -eta2DMax_, eta2DMax_, nbins2D_, -phi2DMax_, phi2DMax_);
0507 etaphihist.push_back(tmphisto);
0508
0509
0510
0511 histName = theHLTCollectionLabels[i].label() + "et";
0512 histTitle = HltHistTitle[i] + " Et";
0513 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, plotPtMin, plotPtMax);
0514 histEtOfHltObjMatchToGen.push_back(tmphisto);
0515
0516
0517
0518 histName = theHLTCollectionLabels[i].label() + "eta";
0519 histTitle = HltHistTitle[i] + " eta";
0520 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -etaMax_, etaMax_);
0521 histEtaOfHltObjMatchToGen.push_back(tmphisto);
0522
0523 if (!noPhiPlots_) {
0524
0525
0526 histName = theHLTCollectionLabels[i].label() + "phi";
0527 histTitle = HltHistTitle[i] + " phi";
0528 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -phiMax_, phiMax_);
0529 histPhiOfHltObjMatchToGen.push_back(tmphisto);
0530 }
0531
0532 histName = theHLTCollectionLabels[i].label() + "etaphi";
0533 histTitle = HltHistTitle[i] + " eta-phi";
0534 tmphisto = iBooker.book2D(
0535 histName.c_str(), histTitle.c_str(), nbins2D_ - 2, -eta2DMax_, eta2DMax_, nbins2D_, -phi2DMax_, phi2DMax_);
0536 histEtaPhiOfHltObjMatchToGen.push_back(tmphisto);
0537 }
0538
0539
0540 histName = theHLTCollectionLabels[i].label() + "et_MC_matched";
0541 histTitle = HltHistTitle[i] + " Et (MC matched)";
0542 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, plotPtMin, plotPtMax);
0543 ethistmatch.push_back(tmphisto);
0544
0545
0546 histName = theHLTCollectionLabels[i].label() + "eta_MC_matched";
0547 histTitle = HltHistTitle[i] + " #eta (MC matched)";
0548 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -etaMax_, etaMax_);
0549 etahistmatch.push_back(tmphisto);
0550
0551 if (!noPhiPlots_) {
0552
0553 histName = theHLTCollectionLabels[i].label() + "phi_MC_matched";
0554 histTitle = HltHistTitle[i] + " #phi (MC matched)";
0555 tmphisto = iBooker.book1D(histName.c_str(), histTitle.c_str(), nbins_, -phiMax_, phiMax_);
0556 phihistmatch.push_back(tmphisto);
0557 }
0558
0559 histName = theHLTCollectionLabels[i].label() + "etaphi_MC_matched";
0560 histTitle = HltHistTitle[i] + " #eta-#phi (MC matched)";
0561 tmphisto = iBooker.book2D(
0562 histName.c_str(), histTitle.c_str(), nbins2D_ - 2, -eta2DMax_, eta2DMax_, nbins2D_, -phi2DMax_, phi2DMax_);
0563 etaphihistmatch.push_back(tmphisto);
0564 }
0565
0566
0567 etahists.push_back(etahist);
0568 phihists.push_back(phihist);
0569 ethists.push_back(ethist);
0570 etahistmatchs.push_back(etahistmatch);
0571 phihistmatchs.push_back(phihistmatch);
0572 ethistmatchs.push_back(ethistmatch);
0573 histEtOfHltObjMatchToGens.push_back(histEtOfHltObjMatchToGen);
0574 histEtaOfHltObjMatchToGens.push_back(histEtaOfHltObjMatchToGen);
0575 histPhiOfHltObjMatchToGens.push_back(histPhiOfHltObjMatchToGen);
0576 etaphihists.push_back(etaphihist);
0577 etaphihistmatchs.push_back(etaphihistmatch);
0578 histEtaPhiOfHltObjMatchToGens.push_back(histEtaPhiOfHltObjMatchToGen);
0579
0580
0581
0582
0583
0584
0585
0586
0587
0588
0589
0590
0591 totals.push_back(total);
0592 totalmatchs.push_back(totalmatch);
0593 etgens.push_back(etgen);
0594 etagens.push_back(etagen);
0595 phigens.push_back(phigen);
0596 etaphigens.push_back(etaphigen);
0597 }
0598 }
0599
0600
0601
0602
0603 EmDQM::~EmDQM() {}
0604
0605
0606 bool EmDQM::checkGeneratedParticlesRequirement(const edm::Event &event) {
0607
0608
0609
0610
0611
0612 edm::Handle<edm::View<reco::Candidate>> genParticles;
0613 event.getByToken(genParticles_token, genParticles);
0614 if (!genParticles.isValid()) {
0615 if (verbosity_ >= OUTPUT_WARNINGS)
0616 edm::LogWarning("EmDQM") << "genParticles invalid.";
0617 return false;
0618 }
0619
0620 std::vector<reco::LeafCandidate> allSortedGenParticles;
0621
0622 for (edm::View<reco::Candidate>::const_iterator currentGenParticle = genParticles->begin();
0623 currentGenParticle != genParticles->end();
0624 currentGenParticle++) {
0625
0626
0627
0628 if (!(abs((*currentGenParticle).pdgId()) == pdgGen && (*currentGenParticle).status() == 1 &&
0629 (*currentGenParticle).et() > 2.0))
0630 continue;
0631
0632 reco::LeafCandidate tmpcand(*(currentGenParticle));
0633
0634 if (tmpcand.et() < plotEtMin)
0635 continue;
0636
0637 allSortedGenParticles.push_back(tmpcand);
0638 }
0639
0640 std::sort(allSortedGenParticles.begin(), allSortedGenParticles.end(), pTGenComparator_);
0641
0642
0643 if (allSortedGenParticles.size() < gencut_)
0644 return false;
0645
0646
0647
0648
0649
0650
0651
0652
0653 for (unsigned int i = 0; i < gencut_; i++) {
0654 bool inECALgap = fabs(allSortedGenParticles[i].eta()) > 1.4442 && fabs(allSortedGenParticles[i].eta()) < 1.556;
0655 if ((fabs(allSortedGenParticles[i].eta()) > genEtaAcc_) || inECALgap) {
0656
0657
0658
0659
0660 return false;
0661 }
0662 }
0663
0664
0665 return true;
0666 }
0667
0668
0669 bool EmDQM::checkRecoParticlesRequirement(const edm::Event &event) {
0670
0671
0672
0673
0674 edm::Handle<edm::View<reco::Candidate>> referenceParticles;
0675
0676 if (autoConfMode_) {
0677 switch (reqNum) {
0678 case 1:
0679 if (pdgGen == 11)
0680 event.getByToken(gencutColl_fidWenu_token, referenceParticles);
0681 else
0682 event.getByToken(gencutColl_fidGammaJet_token, referenceParticles);
0683 break;
0684 case 2:
0685 if (pdgGen == 11)
0686 event.getByToken(gencutColl_fidZee_token, referenceParticles);
0687 else
0688 event.getByToken(gencutColl_fidDiGamma_token, referenceParticles);
0689 break;
0690 case 3:
0691 event.getByToken(gencutColl_fidTripleEle_token, referenceParticles);
0692 break;
0693 }
0694 } else {
0695 event.getByToken(gencutColl_manualConf_token, referenceParticles);
0696 }
0697 if (!referenceParticles.isValid()) {
0698 if (verbosity_ >= OUTPUT_WARNINGS)
0699 edm::LogWarning("EmDQM") << "referenceParticles invalid.";
0700 return false;
0701 }
0702
0703 std::vector<const reco::Candidate *> allSortedReferenceParticles;
0704
0705 for (edm::View<reco::Candidate>::const_iterator currentReferenceParticle = referenceParticles->begin();
0706 currentReferenceParticle != referenceParticles->end();
0707 currentReferenceParticle++) {
0708 if (currentReferenceParticle->et() <= 2.0)
0709 continue;
0710
0711
0712
0713
0714
0715
0716 if (currentReferenceParticle->et() < plotEtMin)
0717 continue;
0718
0719
0720 allSortedReferenceParticles.push_back(&(*currentReferenceParticle));
0721 }
0722
0723
0724
0725
0726
0727 return allSortedReferenceParticles.size() >= gencut_;
0728 }
0729
0730
0731
0732
0733 void EmDQM::analyze(const edm::Event &event, const edm::EventSetup &setup) {
0734
0735 unsigned int vPos = 0;
0736 for (std::vector<edm::ParameterSet>::iterator psetIt = paramSets.begin(); psetIt != paramSets.end();
0737 ++psetIt, ++vPos) {
0738 SetVarsFromPSet(psetIt);
0739
0740 hltCollectionLabelsFound = hltCollectionLabelsFoundPerPath.at(vPos);
0741 hltCollectionLabelsMissed = hltCollectionLabelsMissedPerPath.at(vPos);
0742
0743
0744
0745
0746
0747
0748 edm::Handle<edm::View<reco::Candidate>> cutCounter;
0749 if (autoConfMode_) {
0750 switch (reqNum) {
0751 case 1:
0752 if (pdgGen == 11)
0753 event.getByToken(gencutColl_fidWenu_token, cutCounter);
0754 else
0755 event.getByToken(gencutColl_fidGammaJet_token, cutCounter);
0756 break;
0757 case 2:
0758 if (pdgGen == 11)
0759 event.getByToken(gencutColl_fidZee_token, cutCounter);
0760 else
0761 event.getByToken(gencutColl_fidDiGamma_token, cutCounter);
0762 break;
0763 case 3:
0764 event.getByToken(gencutColl_fidTripleEle_token, cutCounter);
0765 break;
0766 }
0767 } else {
0768 event.getByToken(gencutColl_manualConf_token, cutCounter);
0769 }
0770 if (cutCounter->size() < (unsigned int)gencut_) {
0771
0772
0773 continue;
0774 }
0775
0776
0777
0778 edm::Handle<trigger::TriggerEventWithRefs> triggerObj;
0779 event.getByToken(triggerObject_token, triggerObj);
0780 if (!triggerObj.isValid()) {
0781 if (verbosity_ >= OUTPUT_WARNINGS)
0782 edm::LogWarning("EmDQM") << "parameter triggerobject (" << triggerObject_
0783 << ") does not corresond to a valid TriggerEventWithRefs product. "
0784 "Please check especially the process name (e.g. when running "
0785 "over reprocessed datasets)";
0786 continue;
0787 }
0788
0789
0790 if (event.isRealData()) {
0791
0792
0793
0794
0795
0796 if (!checkRecoParticlesRequirement(event))
0797 continue;
0798 } else {
0799
0800 if (!checkGeneratedParticlesRequirement(event))
0801 continue;
0802 }
0803
0804
0805
0806
0807
0808
0809
0810 if (!mcMatchedOnly_)
0811 totals.at(vPos)->Fill(numOfHLTCollectionLabels + 0.5);
0812 totalmatchs.at(vPos)->Fill(numOfHLTCollectionLabels + 0.5);
0813
0814
0815
0816
0817
0818
0819
0820 std::vector<reco::Particle> sortedGen;
0821 for (edm::View<reco::Candidate>::const_iterator genpart = cutCounter->begin(); genpart != cutCounter->end();
0822 genpart++) {
0823 reco::Particle tmpcand(genpart->charge(), genpart->p4(), genpart->vertex(), genpart->pdgId(), genpart->status());
0824 if (tmpcand.et() >= plotEtMin) {
0825 sortedGen.push_back(tmpcand);
0826 }
0827 }
0828 std::sort(sortedGen.begin(), sortedGen.end(), pTComparator_);
0829
0830
0831
0832
0833 if (sortedGen.size() < gencut_) {
0834 continue;
0835 }
0836 sortedGen.erase(sortedGen.begin() + gencut_, sortedGen.end());
0837
0838 for (unsigned int i = 0; i < gencut_; i++) {
0839 etgens.at(vPos)->Fill(sortedGen[i].et());
0840
0841 if (sortedGen[i].et() > minEtForEtaEffPlot_) {
0842 etagens.at(vPos)->Fill(sortedGen[i].eta());
0843 if (!noPhiPlots_)
0844 phigens.at(vPos)->Fill(sortedGen[i].phi());
0845 etaphigens.at(vPos)->Fill(sortedGen[i].eta(), sortedGen[i].phi());
0846 }
0847 }
0848 if (gencut_ >= reqNum && !mcMatchedOnly_)
0849 totals.at(vPos)->Fill(numOfHLTCollectionLabels +
0850 1.5);
0851 if (gencut_ >= reqNum)
0852 totalmatchs.at(vPos)->Fill(numOfHLTCollectionLabels +
0853 1.5);
0854
0855 bool accepted = true;
0856 edm::Handle<edm::TriggerResults> hltResults;
0857 event.getByToken(hltResults_token, hltResults);
0858
0859
0860
0861 for (unsigned int n = 0; n < numOfHLTCollectionLabels; n++) {
0862
0863
0864 if (sortedGen.size() < nCandCuts.at(n)) {
0865 if (verbosity_ >= OUTPUT_ERRORS)
0866 edm::LogError("EmDQM") << "There are less generated particles than the module '"
0867 << theHLTCollectionLabels[n].label() << "' requires.";
0868 continue;
0869 }
0870 std::vector<reco::Particle> sortedGenForFilter(sortedGen);
0871 sortedGenForFilter.erase(sortedGenForFilter.begin() + nCandCuts.at(n), sortedGenForFilter.end());
0872
0873
0874 if (pathIndex != 0 &&
0875 hltConfig_.moduleIndex(pathIndex, theHLTCollectionLabels[n].label()) > hltResults->index(pathIndex))
0876 break;
0877
0878
0879 switch (theHLTOutputTypes[n]) {
0880 case trigger::TriggerL1NoIsoEG:
0881 histoFillerL1NonIso->fillHistos(triggerObj, event, vPos, n, sortedGenForFilter, accepted);
0882 break;
0883 case trigger::TriggerL1IsoEG:
0884 histoFillerL1Iso->fillHistos(triggerObj, event, vPos, n, sortedGenForFilter, accepted);
0885 break;
0886 case trigger::TriggerPhoton:
0887 histoFillerPho->fillHistos(triggerObj, event, vPos, n, sortedGenForFilter, accepted);
0888 break;
0889 case trigger::TriggerElectron:
0890 histoFillerEle->fillHistos(triggerObj, event, vPos, n, sortedGenForFilter, accepted);
0891 break;
0892 case trigger::TriggerCluster:
0893 histoFillerClu->fillHistos(triggerObj, event, vPos, n, sortedGenForFilter, accepted);
0894 break;
0895 default:
0896 throw(cms::Exception("Release Validation Error") << "HLT output type not implemented: theHLTOutputTypes[n]");
0897 }
0898 }
0899
0900
0901 hltCollectionLabelsFoundPerPath.at(vPos) = hltCollectionLabelsFound;
0902 hltCollectionLabelsMissedPerPath.at(vPos) = hltCollectionLabelsMissed;
0903 }
0904 }
0905
0906
0907
0908
0909
0910 template <class T>
0911 void HistoFiller<T>::fillHistos(edm::Handle<trigger::TriggerEventWithRefs> &triggerObj,
0912 const edm::Event &iEvent,
0913 unsigned int vPos,
0914 unsigned int n,
0915 std::vector<reco::Particle> &sortedGen,
0916 bool &accepted) {
0917 std::vector<edm::Ref<T>> recoecalcands;
0918 if ((triggerObj->filterIndex(dqm->theHLTCollectionLabels[n]) >= triggerObj->size())) {
0919 dqm->hltCollectionLabelsMissed.insert(dqm->theHLTCollectionLabels[n].encode());
0920 accepted = false;
0921 return;
0922 }
0923
0924 dqm->hltCollectionLabelsFound.insert(dqm->theHLTCollectionLabels[n].encode());
0925
0926
0927
0928
0929 triggerObj->getObjects(
0930 triggerObj->filterIndex(dqm->theHLTCollectionLabels[n]), dqm->theHLTOutputTypes[n], recoecalcands);
0931
0932
0933 if (dqm->theHLTOutputTypes[n] == trigger::TriggerL1NoIsoEG) {
0934 std::vector<edm::Ref<T>> isocands;
0935 triggerObj->getObjects(triggerObj->filterIndex(dqm->theHLTCollectionLabels[n]), trigger::TriggerL1IsoEG, isocands);
0936 if (!isocands.empty()) {
0937 for (unsigned int i = 0; i < isocands.size(); i++)
0938 recoecalcands.push_back(isocands[i]);
0939 }
0940 }
0941
0942 if (recoecalcands.empty()) {
0943 accepted = false;
0944 return;
0945 }
0946
0947
0948 if (recoecalcands.size() >= dqm->nCandCuts.at(n) && !dqm->mcMatchedOnly_)
0949 dqm->totals.at(vPos)->Fill(n + 0.5);
0950
0951
0952
0953
0954
0955 for (unsigned int j = 0; j < recoecalcands.size(); j++) {
0956 if (!(recoecalcands.at(j).isAvailable())) {
0957 if (dqm->verbosity_ >= dqm->OUTPUT_ERRORS)
0958 edm::LogError("EmDQMInvalidRefs") << "Event content inconsistent: TriggerEventWithRefs contains "
0959 "invalid Refs. Invalid refs for: "
0960 << dqm->theHLTCollectionLabels[n].label()
0961 << ". The collection that this module uses may has been dropped in "
0962 "the event.";
0963 return;
0964 }
0965 }
0966
0967 if (!dqm->mcMatchedOnly_) {
0968
0969
0970
0971
0972
0973 for (unsigned int i = 0; i < dqm->nCandCuts.at(n); i++) {
0974 math::XYZVector currentGenParticleMomentum = sortedGen[i].momentum();
0975
0976 float closestDeltaR = 0.5;
0977 int closestEcalCandIndex = -1;
0978 for (unsigned int j = 0; j < recoecalcands.size(); j++) {
0979 float deltaR = DeltaR(recoecalcands[j]->momentum(), currentGenParticleMomentum);
0980
0981 if (deltaR < closestDeltaR) {
0982 closestDeltaR = deltaR;
0983 closestEcalCandIndex = j;
0984 }
0985 }
0986
0987
0988
0989 if (closestEcalCandIndex >= 0) {
0990 dqm->histEtOfHltObjMatchToGens.at(vPos).at(n)->Fill(recoecalcands[closestEcalCandIndex]->et());
0991 dqm->histEtaOfHltObjMatchToGens.at(vPos).at(n)->Fill(recoecalcands[closestEcalCandIndex]->eta());
0992 if (!dqm->noPhiPlots_)
0993 dqm->histPhiOfHltObjMatchToGens.at(vPos).at(n)->Fill(recoecalcands[closestEcalCandIndex]->phi());
0994 dqm->histEtaPhiOfHltObjMatchToGens.at(vPos).at(n)->Fill(recoecalcands[closestEcalCandIndex]->eta(),
0995 recoecalcands[closestEcalCandIndex]->phi());
0996 }
0997 }
0998
0999
1000
1001
1002
1003
1004
1005 for (unsigned int i = 0; i < recoecalcands.size(); i++) {
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025 dqm->ethists.at(vPos).at(n)->Fill(recoecalcands[i]->et());
1026 dqm->etahists.at(vPos).at(n)->Fill(recoecalcands[i]->eta());
1027 if (!dqm->noPhiPlots_)
1028 dqm->phihists.at(vPos).at(n)->Fill(recoecalcands[i]->phi());
1029 dqm->etaphihists.at(vPos).at(n)->Fill(recoecalcands[i]->eta(), recoecalcands[i]->phi());
1030 }
1031 }
1032
1033
1034
1035
1036 unsigned int matchedMcParts = 0;
1037 float mindist = 0.3;
1038 if (n == 0)
1039 mindist = 0.5;
1040 for (unsigned int i = 0; i < dqm->nCandCuts.at(n); ++i) {
1041
1042 bool matchThis = false;
1043 math::XYZVector candDir = sortedGen[i].momentum();
1044
1045 double closestDr = 1000.;
1046 for (unsigned int trigOb = 0; trigOb < recoecalcands.size(); ++trigOb) {
1047 double dr = DeltaR(recoecalcands[trigOb]->momentum(), candDir);
1048 if (dr < closestDr) {
1049 closestDr = dr;
1050
1051 }
1052 if (closestDr > mindist) {
1053
1054 } else {
1055 matchedMcParts++;
1056 matchThis = true;
1057 }
1058 }
1059 if (!matchThis) {
1060 accepted = false;
1061 continue;
1062 }
1063
1064 dqm->ethistmatchs.at(vPos).at(n)->Fill(sortedGen[i].et());
1065 if (sortedGen[i].et() > dqm->minEtForEtaEffPlot_) {
1066 dqm->etahistmatchs.at(vPos).at(n)->Fill(sortedGen[i].eta());
1067 if (!dqm->noPhiPlots_)
1068 dqm->phihistmatchs.at(vPos).at(n)->Fill(sortedGen[i].phi());
1069 dqm->etaphihistmatchs.at(vPos).at(n)->Fill(sortedGen[i].eta(), sortedGen[i].phi());
1070 }
1071 }
1072
1073
1074 if (matchedMcParts >= dqm->nCandCuts.at(n) && accepted == true)
1075 dqm->totalmatchs.at(vPos)->Fill(n + 0.5);
1076 }
1077
1078 void EmDQM::dqmEndRun(edm::Run const &iRun, edm::EventSetup const &iSetup) {
1079
1080 unsigned int vPos = 0;
1081 for (std::vector<edm::ParameterSet>::iterator psetIt = paramSets.begin(); psetIt != paramSets.end();
1082 ++psetIt, ++vPos) {
1083 SetVarsFromPSet(psetIt);
1084
1085
1086
1087
1088
1089 std::vector<std::string> labelsNeverFound;
1090
1091 for (const auto &tag : hltCollectionLabelsMissedPerPath.at(vPos)) {
1092 auto atag = edm::InputTag(tag);
1093 if ((hltCollectionLabelsFoundPerPath.at(vPos)).count(atag.encode()) == 0)
1094
1095 labelsNeverFound.push_back(atag.encode());
1096
1097 }
1098
1099 if (labelsNeverFound.empty())
1100 continue;
1101
1102 std::sort(labelsNeverFound.begin(), labelsNeverFound.end());
1103
1104
1105
1106
1107 if (verbosity_ >= OUTPUT_WARNINGS)
1108 edm::LogWarning("EmDQM") << "There were some HLTCollectionLabels which were never found:";
1109
1110 for (auto const &tag : labelsNeverFound) {
1111 if (verbosity_ >= OUTPUT_ALL)
1112 edm::LogPrint("EmDQM") << " " << tag;
1113 }
1114 }
1115 }
1116
1117
1118 int EmDQM::countSubstring(const std::string &str, const std::string &sub) {
1119 if (sub.empty())
1120 return 0;
1121 int count = 0;
1122 for (size_t offset = str.find(sub); offset != std::string::npos; offset = str.find(sub, offset + sub.length())) {
1123 ++count;
1124 }
1125 return count;
1126 }
1127
1128
1129
1130 std::vector<std::vector<std::string>> EmDQM::findEgammaPaths() {
1131 std::vector<std::vector<std::string>> Paths(5);
1132
1133 for (unsigned int i = 0; i < hltConfig_.size(); i++) {
1134 std::string path = hltConfig_.triggerName(i);
1135
1136
1137 if (int(path.find("HLT_")) == 0) {
1138
1139 int scCount = countSubstring(path, "_SC");
1140 int eleCount = countSubstring(path, "Ele");
1141 int doubleEleCount = countSubstring(path, "DoubleEle");
1142 int doubleSCCount = countSubstring(path, "DiSC");
1143 int tripleEleCount = countSubstring(path, "TripleEle");
1144 int photonCount = countSubstring(path, "hoton");
1145 int doublePhotonCount = countSubstring(path, "DoublePhoton") + countSubstring(path, "Diphoton");
1146
1147 int totEleCount = 2 * tripleEleCount + doubleEleCount + eleCount + scCount + 2 * doubleSCCount;
1148 int totPhotonCount = doublePhotonCount + photonCount;
1149
1150 if (totEleCount + totPhotonCount < 1)
1151 continue;
1152 switch (totEleCount) {
1153 case 1:
1154 Paths[TYPE_SINGLE_ELE].push_back(path);
1155
1156 break;
1157 case 2:
1158 Paths[TYPE_DOUBLE_ELE].push_back(path);
1159
1160 break;
1161 case 3:
1162 Paths[TYPE_TRIPLE_ELE].push_back(path);
1163
1164 break;
1165 }
1166
1167 switch (totPhotonCount) {
1168 case 1:
1169 Paths[TYPE_SINGLE_PHOTON].push_back(path);
1170
1171 break;
1172 case 2:
1173 Paths[TYPE_DOUBLE_PHOTON].push_back(path);
1174
1175 break;
1176 }
1177 }
1178
1179
1180 }
1181 return Paths;
1182 }
1183
1184
1185
1186 std::vector<std::string> EmDQM::getFilterModules(const std::string &path) {
1187 std::vector<std::string> filters;
1188
1189
1190
1191
1192 for (unsigned int i = 0; i < hltConfig_.size(path); i++) {
1193 std::string module = hltConfig_.moduleLabel(path, i);
1194 std::string moduleType = hltConfig_.moduleType(module);
1195 std::string moduleEDMType = hltConfig_.moduleEDMType(module);
1196
1197
1198 if (moduleEDMType == "EDFilter" ||
1199 moduleType.find("Filter") != std::string::npos) {
1200
1201 filters.push_back(module);
1202
1203
1204 }
1205 }
1206 return filters;
1207 }
1208
1209
1210
1211 double EmDQM::getPrimaryEtCut(const std::string &path) {
1212 double minEt = -1;
1213
1214 boost::regex reg("^HLT_.*?(Ele|hoton|EG|SC)([[:digit:]]+).*");
1215
1216 boost::smatch what;
1217 if (boost::regex_match(path, what, reg, boost::match_extra)) {
1218 minEt = std::stod(what[2]);
1219 }
1220
1221 return minEt;
1222 }
1223
1224
1225
1226 edm::ParameterSet EmDQM::makePSetForL1SeedFilter(const std::string &moduleName) {
1227
1228
1229
1230
1231 edm::ParameterSet retPSet;
1232
1233 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1234
1235
1236 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1237 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections",
1238 std::vector<edm::InputTag>(1, std::string("none")));
1239 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerL1NoIsoEG);
1240
1241
1242
1243
1244 int orCount = countSubstring(moduleName, "OR") + countSubstring(moduleName, "Or");
1245 int egCount = countSubstring(moduleName, "EG");
1246 int dEgCount = countSubstring(moduleName, "DoubleEG");
1247 int tEgCount = countSubstring(moduleName, "TripleEG");
1248
1249 int candCount = 2 * tEgCount + dEgCount + egCount;
1250
1251
1252
1253 if (orCount > 0 && candCount > 0) {
1254 if (egCount % (orCount + 1) == 0 && dEgCount % (orCount + 1) == 0 && tEgCount % (orCount + 1) == 0)
1255 candCount /= (orCount + 1);
1256 else if (egCount - dEgCount - tEgCount > 0)
1257 candCount = 1;
1258 else if (dEgCount > 0)
1259 candCount = 2;
1260 else if (tEgCount > 0)
1261 candCount = 3;
1262 else
1263 candCount = -1;
1264 }
1265
1266 switch (candCount) {
1267 case 0:
1268 retPSet.addParameter<int>("ncandcut", 0);
1269 break;
1270 case 1:
1271 retPSet.addParameter<int>("ncandcut", 1);
1272 break;
1273 case 2:
1274 retPSet.addParameter<int>("ncandcut", 2);
1275 break;
1276 case 3:
1277 retPSet.addParameter<int>("ncandcut", 3);
1278 break;
1279 default:
1280 retPSet.addParameter<int>("ncandcut", -1);
1281 }
1282
1283 return retPSet;
1284 }
1285
1286
1287
1288 edm::ParameterSet EmDQM::makePSetForL1SeedToSuperClusterMatchFilter(const std::string &moduleName) {
1289
1290
1291
1292
1293
1294
1295 edm::ParameterSet retPSet;
1296
1297 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1298
1299
1300 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1301 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections",
1302 std::vector<edm::InputTag>(1, std::string("none")));
1303 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerCluster);
1304 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1305
1306 return retPSet;
1307 }
1308
1309
1310
1311 edm::ParameterSet EmDQM::makePSetForEtFilter(const std::string &moduleName) {
1312 edm::ParameterSet retPSet;
1313
1314 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1315
1316
1317 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1318 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections",
1319 std::vector<edm::InputTag>(1, std::string("none")));
1320 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerCluster);
1321 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1322
1323 return retPSet;
1324 }
1325
1326
1327
1328 edm::ParameterSet EmDQM::makePSetForOneOEMinusOneOPFilter(const std::string &moduleName) {
1329 edm::ParameterSet retPSet;
1330
1331 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1332
1333
1334 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1335 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections",
1336 std::vector<edm::InputTag>(1, std::string("none")));
1337 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerElectron);
1338 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1339
1340 return retPSet;
1341 }
1342
1343
1344
1345 edm::ParameterSet EmDQM::makePSetForPixelMatchFilter(const std::string &moduleName) {
1346 edm::ParameterSet retPSet;
1347
1348 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1349
1350
1351 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1352 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections",
1353 std::vector<edm::InputTag>(1, std::string("none")));
1354 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerCluster);
1355 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1356
1357 return retPSet;
1358 }
1359
1360
1361
1362 edm::ParameterSet EmDQM::makePSetForEgammaDoubleEtDeltaPhiFilter(const std::string &moduleName) {
1363 edm::ParameterSet retPSet;
1364
1365 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1366
1367
1368 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1369 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections",
1370 std::vector<edm::InputTag>(1, std::string("none")));
1371 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerCluster);
1372 retPSet.addParameter<int>("ncandcut", 2);
1373
1374 return retPSet;
1375 }
1376
1377
1378
1379 edm::ParameterSet EmDQM::makePSetForEgammaGenericFilter(const std::string &moduleName) {
1380 edm::ParameterSet retPSet;
1381
1382
1383
1384
1385
1386
1387
1388
1389
1390
1391
1392
1393
1394
1395 edm::InputTag varTag = hltConfig_.modulePSet(moduleName).getParameter<edm::InputTag>("varTag");
1396
1397
1398
1399
1400
1401
1402 std::string inputType = hltConfig_.moduleType(varTag.label());
1403
1404
1405
1406
1407
1408
1409
1410
1411
1412
1413
1414
1415
1416
1417
1418
1419
1420
1421
1422
1423 if (hltConfig_.saveTags(moduleName))
1424 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerPhoton);
1425 else
1426 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerCluster);
1427
1428 std::vector<edm::InputTag> isoCollections;
1429 isoCollections.push_back(varTag);
1430
1431
1432
1433
1434
1435
1436
1437 if (inputType == "EgammaHLTR9IDProducer" ||
1438 inputType == "EgammaHLTClusterShapeProducer" ||
1439 inputType == "EgammaHLTGsfTrackVarProducer" ||
1440 inputType == "EgammaHLTBcHcalIsolationProducersRegional" ||
1441 inputType == "EgammaHLTHcalVarProducerFromRecHit" ||
1442 inputType == "EgammaHLTEcalPFClusterIsolationProducer" ||
1443 inputType == "EgammaHLTHcalPFClusterIsolationProducer" ||
1444 inputType == "EgammaHLTElectronTrackIsolationProducers"
1445 ) {
1446 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1447
1448
1449 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1450 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections", isoCollections);
1451 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1452
1453 return retPSet;
1454 }
1455
1456 if (verbosity_ >= OUTPUT_ERRORS)
1457 edm::LogError("EmDQM") << "Can't determine what the HLTEgammaGenericFilter '" << moduleName
1458 << "' should do: uses a collection produced by a module of C++ type '" << inputType << "'.";
1459 return edm::ParameterSet();
1460 }
1461
1462
1463
1464 edm::ParameterSet EmDQM::makePSetForEgammaGenericQuadraticFilter(const std::string &moduleName) {
1465 edm::ParameterSet retPSet;
1466
1467
1468
1469
1470
1471
1472
1473
1474
1475
1476
1477
1478
1479
1480 edm::InputTag varTag = hltConfig_.modulePSet(moduleName).getParameter<edm::InputTag>("varTag");
1481
1482
1483
1484
1485
1486 std::string inputType = hltConfig_.moduleType(varTag.label());
1487
1488
1489
1490
1491
1492
1493
1494
1495
1496
1497
1498
1499
1500
1501
1502
1503
1504
1505
1506
1507 if (hltConfig_.saveTags(moduleName))
1508 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerPhoton);
1509 else
1510 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerCluster);
1511
1512 std::vector<edm::InputTag> isoCollections;
1513 isoCollections.push_back(varTag);
1514
1515
1516
1517
1518
1519
1520
1521 if (inputType == "EgammaHLTR9IDProducer" ||
1522 inputType == "EgammaHLTClusterShapeProducer" ||
1523 inputType == "EgammaHLTBcHcalIsolationProducersRegional" ||
1524 inputType == "EgammaHLTHcalVarProducerFromRecHit" ||
1525 inputType == "EgammaHLTEcalPFClusterIsolationProducer" ||
1526 inputType == "EgammaHLTHcalPFClusterIsolationProducer" ||
1527 inputType == "EgammaHLTPhotonTrackIsolationProducersRegional"
1528 ) {
1529 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1530
1531
1532 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1533 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections", isoCollections);
1534 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1535
1536 return retPSet;
1537 }
1538
1539 if (verbosity_ >= OUTPUT_ERRORS)
1540 edm::LogError("EmDQM") << "Can't determine what the HLTEgammaGenericQuadraticFilter '" << moduleName
1541 << "' should do: uses a collection produced by a module of C++ type '" << inputType << "'.";
1542 return edm::ParameterSet();
1543 }
1544
1545
1546
1547 edm::ParameterSet EmDQM::makePSetForElectronGenericFilter(const std::string &moduleName) {
1548 edm::ParameterSet retPSet;
1549
1550
1551
1552
1553
1554
1555
1556
1557
1558
1559
1560 edm::InputTag varTag = hltConfig_.modulePSet(moduleName).getParameter<edm::InputTag>("varTag");
1561
1562
1563
1564
1565
1566 std::string inputType = hltConfig_.moduleType(varTag.label());
1567
1568
1569
1570
1571
1572
1573
1574
1575
1576
1577
1578
1579
1580
1581
1582
1583
1584
1585
1586
1587
1588 retPSet.addParameter<int>("theHLTOutputTypes", trigger::TriggerElectron);
1589
1590 std::vector<edm::InputTag> isoCollections;
1591 isoCollections.push_back(varTag);
1592
1593
1594
1595
1596
1597
1598
1599
1600
1601 if (inputType == "EgammaHLTElectronTrackIsolationProducers"
1602 ) {
1603 retPSet.addParameter<std::vector<double>>("PlotBounds", std::vector<double>(2, 0.0));
1604
1605
1606 retPSet.addParameter<edm::InputTag>("HLTCollectionLabels", edm::InputTag(moduleName, "", triggerObject_.process()));
1607 retPSet.addParameter<std::vector<edm::InputTag>>("IsoCollections", isoCollections);
1608 retPSet.addParameter<int>("ncandcut", hltConfig_.modulePSet(moduleName).getParameter<int>("ncandcut"));
1609
1610 return retPSet;
1611 }
1612
1613 if (verbosity_ >= OUTPUT_ERRORS)
1614 edm::LogError("EmDQM") << "Can't determine what the HLTElectronGenericFilter '" << moduleName
1615 << "' should do: uses a collection produced by a module of C++ type '" << inputType << "'.";
1616 return edm::ParameterSet();
1617 }
1618
1619
1620
1621 void EmDQM::SetVarsFromPSet(std::vector<edm::ParameterSet>::iterator psetIt) {
1622 dirname_ = "HLT/HLTEgammaValidation/" + psetIt->getParameter<std::string>("@module_label");
1623
1624 pathIndex = psetIt->getUntrackedParameter<unsigned int>("pathIndex", 0);
1625
1626 reqNum = psetIt->getParameter<unsigned int>("reqNum");
1627 pdgGen = psetIt->getParameter<int>("pdgGen");
1628
1629 plotEtMin = psetIt->getUntrackedParameter<double>("genEtMin", 0.);
1630 plotPtMin = psetIt->getUntrackedParameter<double>("PtMin", 0.);
1631 plotPtMax = psetIt->getUntrackedParameter<double>("PtMax", 1000.);
1632
1633
1634 gencutCollection_ = psetIt->getParameter<edm::InputTag>("cutcollection");
1635 gencut_ = psetIt->getParameter<int>("cutnum");
1636
1637
1638
1639
1640
1641 std::vector<edm::ParameterSet> filters = psetIt->getParameter<std::vector<edm::ParameterSet>>("filters");
1642
1643
1644 theHLTCollectionLabels.clear();
1645 theHLTOutputTypes.clear();
1646 theHLTCollectionHumanNames.clear();
1647 plotBounds.clear();
1648 isoNames.clear();
1649 plotiso.clear();
1650 nCandCuts.clear();
1651
1652 int i = 0;
1653 for (std::vector<edm::ParameterSet>::iterator filterconf = filters.begin(); filterconf != filters.end();
1654 filterconf++) {
1655 theHLTCollectionLabels.push_back(filterconf->getParameter<edm::InputTag>("HLTCollectionLabels"));
1656 theHLTOutputTypes.push_back(filterconf->getParameter<int>("theHLTOutputTypes"));
1657
1658
1659 theHLTCollectionHumanNames.push_back(
1660 filterconf->getUntrackedParameter<std::string>("HLTCollectionHumanName", theHLTCollectionLabels[i].label()));
1661
1662 std::vector<double> bounds = filterconf->getParameter<std::vector<double>>("PlotBounds");
1663
1664 assert(bounds.size() == 2);
1665 plotBounds.push_back(std::pair<double, double>(bounds[0], bounds[1]));
1666 isoNames.push_back(filterconf->getParameter<std::vector<edm::InputTag>>("IsoCollections"));
1667
1668
1669
1670
1671
1672
1673
1674
1675
1676
1677
1678
1679
1680
1681
1682
1683
1684
1685
1686
1687
1688
1689
1690
1691
1692 assert(!isoNames.back().empty());
1693 if (isoNames.back().at(0).label() == "none") {
1694 plotiso.push_back(false);
1695 } else {
1696 if (!noIsolationPlots_)
1697 plotiso.push_back(true);
1698 else
1699 plotiso.push_back(false);
1700 }
1701 nCandCuts.push_back(filterconf->getParameter<int>("ncandcut"));
1702 i++;
1703 }
1704
1705
1706 numOfHLTCollectionLabels = theHLTCollectionLabels.size();
1707 }
1708
1709
1710
1711 DEFINE_FWK_MODULE(EmDQM);