Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-29 23:10:48

0001 /**\class HeavyFlavorValidation HeavyFlavorValidation.cc HLTriggerOfflineHeavyFlavor/src/HeavyFlavorValidation.cc
0002 
0003  Description: Analyzer to fill Monitoring Elements for muon, dimuon and trigger path efficiency studies (HLT/RECO, RECO/GEN)
0004 
0005  Implementation:
0006      matching is based on closest in delta R, no duplicates allowed. Generated to Global based on momentum at IP; L1, L2, L2v to Global based on position in muon system, L3 to Global based on momentum at IP.
0007 */
0008 // Original Author:  Zoltan Gecse
0009 
0010 #include <memory>
0011 #include <initializer_list>
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/MakerMacros.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "FWCore/Utilities/interface/InputTag.h"
0018 #include "FWCore/Common/interface/TriggerNames.h"
0019 
0020 #include "DataFormats/Common/interface/RefToBase.h"
0021 #include "DataFormats/TrackReco/interface/Track.h"
0022 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0023 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0024 #include "DataFormats/L1Trigger/interface/L1MuonParticle.h"
0025 #include "DataFormats/L1Trigger/interface/L1MuonParticleFwd.h"
0026 #include "DataFormats/Candidate/interface/Candidate.h"
0027 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0028 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0029 #include "DataFormats/MuonReco/interface/Muon.h"
0030 #include "DataFormats/Math/interface/deltaR.h"
0031 #include "DataFormats/Common/interface/TriggerResults.h"
0032 #include "DataFormats/Candidate/interface/Particle.h"
0033 #include "DataFormats/HLTReco/interface/TriggerEventWithRefs.h"
0034 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0035 
0036 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0037 
0038 #include "DQMServices/Core/interface/DQMStore.h"
0039 #include <DQMServices/Core/interface/DQMEDAnalyzer.h>
0040 
0041 #include "CommonTools/Utils/interface/PtComparator.h"
0042 
0043 #include "TLorentzVector.h"
0044 
0045 using namespace std;
0046 using namespace edm;
0047 using namespace reco;
0048 using namespace l1extra;
0049 using namespace trigger;
0050 
0051 class HeavyFlavorValidation : public DQMEDAnalyzer {
0052 public:
0053   explicit HeavyFlavorValidation(const edm::ParameterSet &);
0054   ~HeavyFlavorValidation() override;
0055 
0056 protected:
0057   void dqmBeginRun(const edm::Run &, const edm::EventSetup &) override;
0058   void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override;
0059   void analyze(const edm::Event &, const edm::EventSetup &) override;
0060 
0061 private:
0062   int getMotherId(const Candidate *p);
0063   void match(MonitorElement *me,
0064              vector<LeafCandidate> &from,
0065              vector<LeafCandidate> &to,
0066              double deltaRMatchingCut,
0067              vector<int> &map);
0068   void myBook2D(DQMStore::IBooker &ibooker,
0069                 TString name,
0070                 vector<double> &xBins,
0071                 TString xLabel,
0072                 vector<double> &yBins,
0073                 TString yLabel,
0074                 TString title);
0075   void myBook2D(DQMStore::IBooker &ibooker,
0076                 TString name,
0077                 vector<double> &xBins,
0078                 TString xLabel,
0079                 vector<double> &yBins,
0080                 TString yLabel) {
0081     myBook2D(ibooker, name, xBins, xLabel, yBins, yLabel, name);
0082   }
0083   void myBookProfile2D(DQMStore::IBooker &ibooker,
0084                        TString name,
0085                        vector<double> &xBins,
0086                        TString xLabel,
0087                        vector<double> &yBins,
0088                        TString yLabel,
0089                        TString title);
0090   void myBookProfile2D(DQMStore::IBooker &ibooker,
0091                        TString name,
0092                        vector<double> &xBins,
0093                        TString xLabel,
0094                        vector<double> &yBins,
0095                        TString yLabel) {
0096     myBookProfile2D(ibooker, name, xBins, xLabel, yBins, yLabel, name);
0097   }
0098   void myBook1D(DQMStore::IBooker &ibooker, TString name, vector<double> &xBins, TString label, TString title);
0099   void myBook1D(DQMStore::IBooker &ibooker, TString name, vector<double> &xBins, TString label) {
0100     myBook1D(ibooker, name, xBins, label, name);
0101   }
0102 
0103   /**
0104      * Get the filter "level" (as it is defined for the use of this module and its corresponding
0105      * harvesting module).
0106      *
0107      * level 1 - 3 -> more or less synonymously to the the trigger levels
0108      * level 4 and 5 -> vertex, dz, track, etc.. filters
0109      *
0110      * See the comments in the definition for some more details.
0111      */
0112   int getFilterLevel(const std::string &moduleName, const HLTConfigProvider &hltConfig);
0113 
0114   string dqmFolder;
0115   string triggerProcessName;
0116   string triggerPathName;
0117 
0118   EDGetTokenT<TriggerEventWithRefs> triggerSummaryRAWTag;
0119   EDGetTokenT<TriggerEvent> triggerSummaryAODTag;
0120   InputTag triggerResultsTag;
0121   EDGetTokenT<TriggerResults> triggerResultsToken;
0122   InputTag recoMuonsTag;
0123   EDGetTokenT<MuonCollection> recoMuonsToken;
0124   InputTag genParticlesTag;
0125   EDGetTokenT<GenParticleCollection> genParticlesToken;
0126 
0127   vector<int> motherIDs;
0128   double genGlobDeltaRMatchingCut;
0129   double globL1DeltaRMatchingCut;
0130   double globL2DeltaRMatchingCut;
0131   double globL3DeltaRMatchingCut;
0132   vector<double> deltaEtaBins;
0133   vector<double> deltaPhiBins;
0134   vector<double> muonPtBins;
0135   vector<double> muonEtaBins;
0136   vector<double> muonPhiBins;
0137   vector<double> dimuonPtBins;
0138   vector<double> dimuonEtaBins;
0139   vector<double> dimuonDRBins;
0140   map<TString, MonitorElement *> ME;
0141   vector<pair<string, int> > filterNamesLevels;
0142   const double muonMass;
0143 };
0144 
0145 HeavyFlavorValidation::HeavyFlavorValidation(const ParameterSet &pset)
0146     :  //get parameters
0147       dqmFolder(pset.getUntrackedParameter<string>("DQMFolder")),
0148       triggerProcessName(pset.getUntrackedParameter<string>("TriggerProcessName")),
0149       triggerPathName(pset.getUntrackedParameter<string>("TriggerPathName")),
0150       motherIDs(pset.getUntrackedParameter<vector<int> >("MotherIDs")),
0151       genGlobDeltaRMatchingCut(pset.getUntrackedParameter<double>("GenGlobDeltaRMatchingCut")),
0152       globL1DeltaRMatchingCut(pset.getUntrackedParameter<double>("GlobL1DeltaRMatchingCut")),
0153       globL2DeltaRMatchingCut(pset.getUntrackedParameter<double>("GlobL2DeltaRMatchingCut")),
0154       globL3DeltaRMatchingCut(pset.getUntrackedParameter<double>("GlobL3DeltaRMatchingCut")),
0155       deltaEtaBins(pset.getUntrackedParameter<vector<double> >("DeltaEtaBins")),
0156       deltaPhiBins(pset.getUntrackedParameter<vector<double> >("DeltaPhiBins")),
0157       muonPtBins(pset.getUntrackedParameter<vector<double> >("MuonPtBins")),
0158       muonEtaBins(pset.getUntrackedParameter<vector<double> >("MuonEtaBins")),
0159       muonPhiBins(pset.getUntrackedParameter<vector<double> >("MuonPhiBins")),
0160       dimuonPtBins(pset.getUntrackedParameter<vector<double> >("DimuonPtBins")),
0161       dimuonEtaBins(pset.getUntrackedParameter<vector<double> >("DimuonEtaBins")),
0162       dimuonDRBins(pset.getUntrackedParameter<vector<double> >("DimuonDRBins")),
0163       muonMass(0.106) {
0164   triggerSummaryRAWTag = consumes<TriggerEventWithRefs>(
0165       InputTag(pset.getUntrackedParameter<string>("TriggerSummaryRAW"), "", triggerProcessName));
0166   triggerSummaryAODTag =
0167       consumes<TriggerEvent>(InputTag(pset.getUntrackedParameter<string>("TriggerSummaryAOD"), "", triggerProcessName));
0168   triggerResultsTag = InputTag(pset.getUntrackedParameter<string>("TriggerResults"), "", triggerProcessName);
0169   triggerResultsToken = consumes<TriggerResults>(triggerResultsTag);
0170   recoMuonsTag = pset.getParameter<InputTag>("RecoMuons");
0171   recoMuonsToken = consumes<MuonCollection>(recoMuonsTag);
0172   genParticlesTag = pset.getParameter<InputTag>("GenParticles");
0173   genParticlesToken = consumes<GenParticleCollection>(genParticlesTag);
0174 }
0175 
0176 void HeavyFlavorValidation::dqmBeginRun(const edm::Run &iRun, const edm::EventSetup &iSetup) {
0177   //discover HLT configuration
0178   HLTConfigProvider hltConfig;
0179   bool isChanged;
0180   if (hltConfig.init(iRun, iSetup, triggerProcessName, isChanged)) {
0181     LogDebug("HLTriggerOfflineHeavyFlavor")
0182         << "Successfully initialized HLTConfigProvider with process name: " << triggerProcessName << endl;
0183   } else {
0184     LogWarning("HLTriggerOfflineHeavyFlavor")
0185         << "Could not initialize HLTConfigProvider with process name: " << triggerProcessName << endl;
0186     return;
0187   }
0188   vector<string> triggerNames = hltConfig.triggerNames();
0189   for (const auto &trigName : triggerNames) {
0190     // TString triggerName = trigName;
0191     if (trigName.find(triggerPathName) != std::string::npos) {
0192       vector<string> moduleNames = hltConfig.moduleLabels(trigName);
0193       for (const auto &moduleName : moduleNames) {
0194         const int level = getFilterLevel(moduleName, hltConfig);
0195         if (level > 0) {
0196           filterNamesLevels.push_back({moduleName, level});
0197         }
0198       }
0199       break;
0200     }
0201   }
0202 
0203   if (filterNamesLevels.empty()) {
0204     LogDebug("HLTriggerOfflineHeavyFlavor") << "Bad Trigger Path: " << triggerPathName << endl;
0205     return;
0206   } else {
0207     std::string str;
0208     str.reserve(
0209         512);  // avoid too many realloctions in the following loop (allows for filter names with roughly 100 chars each)
0210     for (const auto &filters : filterNamesLevels)
0211       str = str + " " + filters.first;
0212     LogDebug("HLTriggerOfflineHeavyFlavor") << "Trigger Path: " << triggerPathName << " has filters:" << str;
0213   }
0214 }
0215 
0216 void HeavyFlavorValidation::bookHistograms(DQMStore::IBooker &ibooker,
0217                                            edm::Run const &iRun,
0218                                            edm::EventSetup const &iSetup) {
0219   ibooker.cd();
0220   ibooker.setCurrentFolder((dqmFolder + "/") + triggerProcessName + "/" + triggerPathName);
0221 
0222   // create Monitor Elements
0223   // Eta Pt Single
0224   myBook2D(ibooker, "genMuon_genEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
0225   myBook2D(ibooker, "globMuon_genEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
0226   myBook2D(ibooker, "globMuon_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
0227 
0228   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0229     myBook2D(ibooker,
0230              TString::Format("filt%dMuon_recoEtaPt", int(i + 1)),
0231              muonEtaBins,
0232              "#mu eta",
0233              muonPtBins,
0234              " #mu pT (GeV)",
0235              filterNamesLevels[i].first);
0236   }
0237   myBook2D(ibooker, "pathMuon_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)", triggerPathName);
0238   myBook2D(ibooker, "resultMuon_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
0239   // Eta Pt Single Resolution
0240   myBookProfile2D(ibooker, "resGlobGen_genEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
0241   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0242     myBookProfile2D(ibooker,
0243                     TString::Format("resFilt%dGlob_recoEtaPt", int(i + 1)),
0244                     muonEtaBins,
0245                     "#mu eta",
0246                     muonPtBins,
0247                     " #mu pT (GeV)",
0248                     filterNamesLevels[i].first);
0249   }
0250   myBookProfile2D(
0251       ibooker, "resPathGlob_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)", triggerPathName);
0252   // Eta Pt Double
0253   myBook2D(ibooker, "genDimuon_genEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
0254   myBook2D(ibooker, "globDimuon_genEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
0255   myBook2D(ibooker, "globDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
0256   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0257     myBook2D(ibooker,
0258              TString::Format("filt%dDimuon_recoEtaPt", int(i + 1)),
0259              dimuonEtaBins,
0260              "#mu#mu eta",
0261              dimuonPtBins,
0262              " #mu#mu pT (GeV)",
0263              filterNamesLevels[i].first);
0264   }
0265   myBook2D(
0266       ibooker, "pathDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)", triggerPathName);
0267   myBook2D(ibooker, "resultDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
0268   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0269     myBook2D(ibooker,
0270              TString::Format("diFilt%dDimuon_recoEtaPt", int(i + 1)),
0271              dimuonEtaBins,
0272              "#mu#mu eta",
0273              dimuonPtBins,
0274              " #mu#mu pT (GeV)",
0275              filterNamesLevels[i].first);
0276   }
0277   myBook2D(
0278       ibooker, "diPathDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)", triggerPathName);
0279   // Eta Phi Single
0280   myBook2D(ibooker, "genMuon_genEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
0281   myBook2D(ibooker, "globMuon_genEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
0282   myBook2D(ibooker, "globMuon_recoEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
0283   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0284     myBook2D(ibooker,
0285              TString::Format("filt%dMuon_recoEtaPhi", int(i + 1)),
0286              muonEtaBins,
0287              "#mu eta",
0288              muonPhiBins,
0289              "#mu phi",
0290              filterNamesLevels[i].first);
0291   }
0292   myBook2D(ibooker, "pathMuon_recoEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi", triggerPathName);
0293   myBook2D(ibooker, "resultMuon_recoEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
0294   // Rap Pt Double
0295   myBook2D(ibooker, "genDimuon_genRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
0296   myBook2D(ibooker, "globDimuon_genRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
0297   myBook2D(ibooker, "globDimuon_recoRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
0298   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0299     myBook2D(ibooker,
0300              TString::Format("filt%dDimuon_recoRapPt", int(i + 1)),
0301              dimuonEtaBins,
0302              "#mu#mu rapidity",
0303              dimuonPtBins,
0304              " #mu#mu pT (GeV)",
0305              filterNamesLevels[i].first);
0306   }
0307   myBook2D(ibooker,
0308            "pathDimuon_recoRapPt",
0309            dimuonEtaBins,
0310            "#mu#mu rapidity",
0311            dimuonPtBins,
0312            " #mu#mu pT (GeV)",
0313            triggerPathName);
0314   myBook2D(ibooker, "resultDimuon_recoRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
0315   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0316     myBook2D(ibooker,
0317              TString::Format("diFilt%dDimuon_recoRapPt", int(i + 1)),
0318              dimuonEtaBins,
0319              "#mu#mu rapidity",
0320              dimuonPtBins,
0321              " #mu#mu pT (GeV)",
0322              filterNamesLevels[i].first);
0323   }
0324   myBook2D(ibooker,
0325            "diPathDimuon_recoRapPt",
0326            dimuonEtaBins,
0327            "#mu#mu rapidity",
0328            dimuonPtBins,
0329            " #mu#mu pT (GeV)",
0330            triggerPathName);
0331   // Pt DR Double
0332   myBook2D(ibooker, "genDimuon_genPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
0333   myBook2D(ibooker, "globDimuon_genPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
0334   myBook2D(ibooker, "globDimuon_recoPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
0335   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0336     myBook2D(ibooker,
0337              TString::Format("filt%dDimuon_recoPtDR", int(i + 1)),
0338              dimuonPtBins,
0339              " #mu#mu pT (GeV)",
0340              dimuonDRBins,
0341              "#mu#mu #Delta R at IP",
0342              filterNamesLevels[i].first);
0343   }
0344   myBook2D(ibooker,
0345            "pathDimuon_recoPtDR",
0346            dimuonPtBins,
0347            " #mu#mu pT (GeV)",
0348            dimuonDRBins,
0349            "#mu#mu #Delta R at IP",
0350            triggerPathName);
0351   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0352     myBook2D(ibooker,
0353              TString::Format("diFilt%dDimuon_recoPtDR", int(i + 1)),
0354              dimuonPtBins,
0355              " #mu#mu pT (GeV)",
0356              dimuonDRBins,
0357              "#mu#mu #Delta R at IP",
0358              filterNamesLevels[i].first);
0359   }
0360   myBook2D(ibooker,
0361            "diPathDimuon_recoPtDR",
0362            dimuonPtBins,
0363            " #mu#mu pT (GeV)",
0364            dimuonDRBins,
0365            "#mu#mu #Delta R at IP",
0366            triggerPathName);
0367   myBook2D(ibooker, "resultDimuon_recoPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
0368   // Pt DRpos Double
0369   myBook2D(ibooker, "globDimuon_recoPtDRpos", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS");
0370   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0371     myBook2D(ibooker,
0372              TString::Format("filt%dDimuon_recoPtDRpos", int(i + 1)),
0373              dimuonPtBins,
0374              " #mu#mu pT (GeV)",
0375              dimuonDRBins,
0376              "#mu#mu #Delta R in MS",
0377              filterNamesLevels[i].first);
0378   }
0379   myBook2D(ibooker,
0380            "pathDimuon_recoPtDRpos",
0381            dimuonPtBins,
0382            " #mu#mu pT (GeV)",
0383            dimuonDRBins,
0384            "#mu#mu #Delta R in MS",
0385            triggerPathName);
0386   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0387     myBook2D(ibooker,
0388              TString::Format("diFilt%dDimuon_recoPtDRpos", int(i + 1)),
0389              dimuonPtBins,
0390              " #mu#mu pT (GeV)",
0391              dimuonDRBins,
0392              "#mu#mu #Delta R in MS",
0393              filterNamesLevels[i].first);
0394   }
0395   myBook2D(ibooker,
0396            "diPathDimuon_recoPtDRpos",
0397            dimuonPtBins,
0398            " #mu#mu pT (GeV)",
0399            dimuonDRBins,
0400            "#mu#mu #Delta R in MS",
0401            triggerPathName);
0402   myBook2D(
0403       ibooker, "resultDimuon_recoPtDRpos", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS");
0404 
0405   // Matching
0406   myBook2D(ibooker, "globGen_deltaEtaDeltaPhi", deltaEtaBins, "#Delta eta", deltaPhiBins, "#Delta phi");
0407   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0408     myBook2D(ibooker,
0409              TString::Format("filt%dGlob_deltaEtaDeltaPhi", int(i + 1)),
0410              deltaEtaBins,
0411              "#Delta eta",
0412              deltaPhiBins,
0413              "#Delta phi",
0414              filterNamesLevels[i].first);
0415   }
0416   myBook2D(
0417       ibooker, "pathGlob_deltaEtaDeltaPhi", deltaEtaBins, "#Delta eta", deltaPhiBins, "#Delta phi", triggerPathName);
0418   // Size of containers
0419   vector<double> sizeBins;
0420   sizeBins.push_back(10);
0421   sizeBins.push_back(0);
0422   sizeBins.push_back(10);
0423   myBook1D(ibooker, "genMuon_size", sizeBins, "container size");
0424   myBook1D(ibooker, "globMuon_size", sizeBins, "container size");
0425   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0426     myBook1D(
0427         ibooker, TString::Format("filt%dMuon_size", int(i + 1)), sizeBins, "container size", filterNamesLevels[i].first);
0428   }
0429   myBook1D(ibooker, "pathMuon_size", sizeBins, "container size", triggerPathName);
0430 }
0431 
0432 void HeavyFlavorValidation::analyze(const Event &iEvent, const EventSetup &iSetup) {
0433   if (filterNamesLevels.empty()) {
0434     return;
0435   }
0436   //access the containers and create LeafCandidate copies
0437   vector<LeafCandidate> genMuons;
0438   Handle<GenParticleCollection> genParticles;
0439   iEvent.getByToken(genParticlesToken, genParticles);
0440   if (genParticles.isValid()) {
0441     for (GenParticleCollection::const_iterator p = genParticles->begin(); p != genParticles->end(); ++p) {
0442       if (p->status() == 1 && std::abs(p->pdgId()) == 13 &&
0443           (find(motherIDs.begin(), motherIDs.end(), -1) != motherIDs.end() ||
0444            find(motherIDs.begin(), motherIDs.end(), getMotherId(&(*p))) != motherIDs.end())) {
0445         genMuons.push_back(*p);
0446       }
0447     }
0448   } else {
0449     LogDebug("HLTriggerOfflineHeavyFlavor") << "Could not access GenParticleCollection" << endl;
0450   }
0451   sort(genMuons.begin(), genMuons.end(), GreaterByPt<LeafCandidate>());
0452   ME["genMuon_size"]->Fill(genMuons.size());
0453   LogDebug("HLTriggerOfflineHeavyFlavor")
0454       << "GenParticleCollection from " << genParticlesTag << " has size: " << genMuons.size() << endl;
0455 
0456   vector<LeafCandidate> globMuons;
0457   vector<LeafCandidate> globMuons_position;
0458   Handle<MuonCollection> recoMuonsHandle;
0459   iEvent.getByToken(recoMuonsToken, recoMuonsHandle);
0460   if (recoMuonsHandle.isValid()) {
0461     for (MuonCollection::const_iterator p = recoMuonsHandle->begin(); p != recoMuonsHandle->end(); ++p) {
0462       if (p->isGlobalMuon()) {
0463         globMuons.push_back(*p);
0464         globMuons_position.push_back(LeafCandidate(p->charge(),
0465                                                    math::XYZTLorentzVector(p->outerTrack()->innerPosition().x(),
0466                                                                            p->outerTrack()->innerPosition().y(),
0467                                                                            p->outerTrack()->innerPosition().z(),
0468                                                                            0.)));
0469       }
0470     }
0471   } else {
0472     LogDebug("HLTriggerOfflineHeavyFlavor") << "Could not access reco Muons" << endl;
0473   }
0474   ME["globMuon_size"]->Fill(globMuons.size());
0475   LogDebug("HLTriggerOfflineHeavyFlavor")
0476       << "Global Muons from " << recoMuonsTag << " has size: " << globMuons.size() << endl;
0477 
0478   // access RAW trigger event
0479   vector<vector<LeafCandidate> > muonsAtFilter;
0480   vector<vector<LeafCandidate> > muonPositionsAtFilter;
0481   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0482     muonsAtFilter.push_back(vector<LeafCandidate>());
0483     muonPositionsAtFilter.push_back(vector<LeafCandidate>());
0484   }
0485   Handle<TriggerEventWithRefs> rawTriggerEvent;
0486   iEvent.getByToken(triggerSummaryRAWTag, rawTriggerEvent);
0487   if (rawTriggerEvent.isValid()) {
0488     for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0489       size_t index = rawTriggerEvent->filterIndex(InputTag(filterNamesLevels[i].first, "", triggerProcessName));
0490       if (index < rawTriggerEvent->size()) {
0491         if (filterNamesLevels[i].second == 1) {
0492           vector<L1MuonParticleRef> l1Cands;
0493           rawTriggerEvent->getObjects(index, TriggerL1Mu, l1Cands);
0494           for (size_t j = 0; j < l1Cands.size(); j++) {
0495             muonsAtFilter[i].push_back(*l1Cands[j]);
0496           }
0497         } else {
0498           vector<RecoChargedCandidateRef> hltCands;
0499           rawTriggerEvent->getObjects(index, TriggerMuon, hltCands);
0500           for (size_t j = 0; j < hltCands.size(); j++) {
0501             muonsAtFilter[i].push_back(*hltCands[j]);
0502             if (filterNamesLevels[i].second == 2) {
0503               muonPositionsAtFilter[i].push_back(
0504                   LeafCandidate(hltCands[j]->charge(),
0505                                 math::XYZTLorentzVector(hltCands[j]->track()->innerPosition().x(),
0506                                                         hltCands[j]->track()->innerPosition().y(),
0507                                                         hltCands[j]->track()->innerPosition().z(),
0508                                                         0.)));
0509             }
0510           }
0511         }
0512       }
0513       ME[TString::Format("filt%dMuon_size", int(i + 1))]->Fill(muonsAtFilter[i].size());
0514       LogDebug("HLTriggerOfflineHeavyFlavor")
0515           << "Filter \"" << filterNamesLevels[i].first << "\" has " << muonsAtFilter[i].size() << " muons" << endl;
0516     }
0517   } else {
0518     LogDebug("HLTriggerOfflineHeavyFlavor") << "Could not access RAWTriggerEvent" << endl;
0519   }
0520 
0521   // access AOD trigger event
0522   vector<LeafCandidate> pathMuons;
0523   Handle<TriggerEvent> aodTriggerEvent;
0524   iEvent.getByToken(triggerSummaryAODTag, aodTriggerEvent);
0525   if (aodTriggerEvent.isValid()) {
0526     TriggerObjectCollection allObjects = aodTriggerEvent->getObjects();
0527     for (int i = 0; i < aodTriggerEvent->sizeFilters(); i++) {
0528       if (aodTriggerEvent->filterTag(i) == InputTag((filterNamesLevels.end() - 1)->first, "", triggerProcessName)) {
0529         Keys keys = aodTriggerEvent->filterKeys(i);
0530         for (size_t j = 0; j < keys.size(); j++) {
0531           pathMuons.push_back(LeafCandidate(
0532               allObjects[keys[j]].id() > 0 ? 1 : -1,
0533               math::PtEtaPhiMLorentzVector(
0534                   allObjects[keys[j]].pt(), allObjects[keys[j]].eta(), allObjects[keys[j]].phi(), muonMass)));
0535         }
0536       }
0537     }
0538     ME["pathMuon_size"]->Fill(pathMuons.size());
0539     LogDebug("HLTriggerOfflineHeavyFlavor")
0540         << "Path \"" << triggerPathName << "\" has " << pathMuons.size() << " muons at last filter \""
0541         << (filterNamesLevels.end() - 1)->first << "\"" << endl;
0542   } else {
0543     LogDebug("HLTriggerOfflineHeavyFlavor") << "Could not access AODTriggerEvent" << endl;
0544   }
0545 
0546   // access Trigger Results
0547   bool triggerFired = false;
0548   Handle<TriggerResults> triggerResults;
0549   iEvent.getByToken(triggerResultsToken, triggerResults);
0550   if (triggerResults.isValid()) {
0551     LogDebug("HLTriggerOfflineHeavyFlavor") << "Successfully initialized " << triggerResultsTag << endl;
0552     const edm::TriggerNames &triggerNames = iEvent.triggerNames(*triggerResults);
0553     bool hlt_exists = false;
0554     for (unsigned int i = 0; i != triggerNames.size(); i++) {
0555       TString hlt_name = triggerNames.triggerName(i);
0556       if (hlt_name.Contains(triggerPathName)) {
0557         triggerFired = triggerResults->accept(i);
0558         hlt_exists = true;
0559         break;
0560       }
0561     }
0562     if (!hlt_exists) {
0563       LogDebug("HLTriggerOfflineHeavyFlavor") << triggerResultsTag << " has no trigger: " << triggerPathName << endl;
0564     }
0565   } else {
0566     LogDebug("HLTriggerOfflineHeavyFlavor") << "Could not initialize " << triggerResultsTag << endl;
0567   }
0568 
0569   //create matching maps
0570   vector<int> glob_gen(genMuons.size(), -1);
0571   match(ME["globGen_deltaEtaDeltaPhi"], genMuons, globMuons, genGlobDeltaRMatchingCut, glob_gen);
0572   vector<vector<int> > filt_glob;
0573   for (size_t i = 0; i < filterNamesLevels.size(); i++) {
0574     filt_glob.push_back(vector<int>(globMuons.size(), -1));
0575     if (filterNamesLevels[i].second == 1) {
0576       match(ME[TString::Format("filt%dGlob_deltaEtaDeltaPhi", int(i + 1))],
0577             globMuons_position,
0578             muonsAtFilter[i],
0579             globL1DeltaRMatchingCut,
0580             filt_glob[i]);
0581     } else if (filterNamesLevels[i].second == 2) {
0582       match(ME[TString::Format("filt%dGlob_deltaEtaDeltaPhi", int(i + 1))],
0583             globMuons_position,
0584             muonPositionsAtFilter[i],
0585             globL2DeltaRMatchingCut,
0586             filt_glob[i]);
0587     } else if (filterNamesLevels[i].second > 2) {
0588       match(ME[TString::Format("filt%dGlob_deltaEtaDeltaPhi", int(i + 1))],
0589             globMuons,
0590             muonsAtFilter[i],
0591             globL3DeltaRMatchingCut,
0592             filt_glob[i]);
0593     }
0594   }
0595   vector<int> path_glob(globMuons.size(), -1);
0596   if ((filterNamesLevels.end() - 1)->second == 1) {
0597     match(ME["pathGlob_deltaEtaDeltaPhi"], globMuons_position, pathMuons, globL1DeltaRMatchingCut, path_glob);
0598   } else if ((filterNamesLevels.end() - 1)->second == 2) {
0599     match(ME["pathGlob_deltaEtaDeltaPhi"], globMuons, pathMuons, globL2DeltaRMatchingCut, path_glob);
0600   } else if ((filterNamesLevels.end() - 1)->second > 2) {
0601     match(ME["pathGlob_deltaEtaDeltaPhi"], globMuons, pathMuons, globL3DeltaRMatchingCut, path_glob);
0602   }
0603 
0604   //fill histos
0605   bool first = true;
0606   for (size_t i = 0; i < genMuons.size(); i++) {
0607     ME["genMuon_genEtaPt"]->Fill(genMuons[i].eta(), genMuons[i].pt());
0608     ME["genMuon_genEtaPhi"]->Fill(genMuons[i].eta(), genMuons[i].phi());
0609     if (glob_gen[i] != -1) {
0610       ME["resGlobGen_genEtaPt"]->Fill(
0611           genMuons[i].eta(), genMuons[i].pt(), (globMuons[glob_gen[i]].pt() - genMuons[i].pt()) / genMuons[i].pt());
0612       ME["globMuon_genEtaPt"]->Fill(genMuons[i].eta(), genMuons[i].pt());
0613       ME["globMuon_genEtaPhi"]->Fill(genMuons[i].eta(), genMuons[i].phi());
0614       ME["globMuon_recoEtaPt"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
0615       ME["globMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
0616       for (size_t f = 0; f < filterNamesLevels.size(); f++) {
0617         if (filt_glob[f][glob_gen[i]] != -1) {
0618           ME[TString::Format("resFilt%dGlob_recoEtaPt", int(f + 1))]->Fill(
0619               globMuons[glob_gen[i]].eta(),
0620               globMuons[glob_gen[i]].pt(),
0621               (muonsAtFilter[f][filt_glob[f][glob_gen[i]]].pt() - globMuons[glob_gen[i]].pt()) /
0622                   globMuons[glob_gen[i]].pt());
0623           ME[TString::Format("filt%dMuon_recoEtaPt", int(f + 1))]->Fill(globMuons[glob_gen[i]].eta(),
0624                                                                         globMuons[glob_gen[i]].pt());
0625           ME[TString::Format("filt%dMuon_recoEtaPhi", int(f + 1))]->Fill(globMuons[glob_gen[i]].eta(),
0626                                                                          globMuons[glob_gen[i]].phi());
0627         } else {
0628           break;
0629         }
0630       }
0631       if (path_glob[glob_gen[i]] != -1) {
0632         ME["resPathGlob_recoEtaPt"]->Fill(
0633             globMuons[glob_gen[i]].eta(),
0634             globMuons[glob_gen[i]].pt(),
0635             (pathMuons[path_glob[glob_gen[i]]].pt() - globMuons[glob_gen[i]].pt()) / globMuons[glob_gen[i]].pt());
0636         ME["pathMuon_recoEtaPt"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
0637         ME["pathMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
0638       }
0639       //highest pt muon
0640       if (first) {
0641         first = false;
0642         if (triggerFired) {
0643           ME["resultMuon_recoEtaPt"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
0644           ME["resultMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
0645         }
0646       }
0647     }
0648   }
0649 
0650   //fill dimuon histograms (highest pT, opposite charge)
0651   int secondMuon = 0;
0652   for (size_t j = 1; j < genMuons.size(); j++) {
0653     if (genMuons[0].charge() * genMuons[j].charge() == -1) {
0654       secondMuon = j;
0655       break;
0656     }
0657   }
0658   if (secondMuon > 0) {
0659     //two generated
0660     double genDimuonPt = (genMuons[0].p4() + genMuons[secondMuon].p4()).pt();
0661     double genDimuonEta = (genMuons[0].p4() + genMuons[secondMuon].p4()).eta();
0662     double genDimuonRap = (genMuons[0].p4() + genMuons[secondMuon].p4()).Rapidity();
0663     double genDimuonDR = deltaR<LeafCandidate, LeafCandidate>(genMuons[0], genMuons[secondMuon]);
0664     bool highPt = genMuons[0].pt() > 7. && genMuons[secondMuon].pt() > 7;
0665     ME["genDimuon_genEtaPt"]->Fill(genDimuonEta, genDimuonPt);
0666     ME["genDimuon_genRapPt"]->Fill(genDimuonRap, genDimuonPt);
0667     if (highPt)
0668       ME["genDimuon_genPtDR"]->Fill(genDimuonPt, genDimuonDR);
0669     //two global
0670     if (glob_gen[0] != -1 && glob_gen[secondMuon] != -1) {
0671       ME["globDimuon_genEtaPt"]->Fill(genDimuonEta, genDimuonPt);
0672       ME["globDimuon_genRapPt"]->Fill(genDimuonRap, genDimuonPt);
0673       if (highPt)
0674         ME["globDimuon_genPtDR"]->Fill(genDimuonPt, genDimuonDR);
0675       double globDimuonPt = (globMuons[glob_gen[0]].p4() + globMuons[glob_gen[secondMuon]].p4()).pt();
0676       double globDimuonEta = (globMuons[glob_gen[0]].p4() + globMuons[glob_gen[secondMuon]].p4()).eta();
0677       double globDimuonRap = (globMuons[glob_gen[0]].p4() + globMuons[glob_gen[secondMuon]].p4()).Rapidity();
0678       double globDimuonDR =
0679           deltaR<LeafCandidate, LeafCandidate>(globMuons[glob_gen[0]], globMuons[glob_gen[secondMuon]]);
0680       double globDimuonDRpos = deltaR<LeafCandidate, LeafCandidate>(globMuons_position[glob_gen[0]],
0681                                                                     globMuons_position[glob_gen[secondMuon]]);
0682       ME["globDimuon_recoEtaPt"]->Fill(globDimuonEta, globDimuonPt);
0683       ME["globDimuon_recoRapPt"]->Fill(globDimuonRap, globDimuonPt);
0684       if (highPt)
0685         ME["globDimuon_recoPtDR"]->Fill(globDimuonPt, globDimuonDR);
0686       if (highPt)
0687         ME["globDimuon_recoPtDRpos"]->Fill(globDimuonPt, globDimuonDRpos);
0688       //two filter objects
0689       for (size_t f = 0; f < filterNamesLevels.size(); f++) {
0690         if (filt_glob[f][glob_gen[0]] != -1 && filt_glob[f][glob_gen[secondMuon]] != -1) {
0691           ME[TString::Format("diFilt%dDimuon_recoEtaPt", int(f + 1))]->Fill(globDimuonEta, globDimuonPt);
0692           ME[TString::Format("diFilt%dDimuon_recoRapPt", int(f + 1))]->Fill(globDimuonRap, globDimuonPt);
0693           if (highPt)
0694             ME[TString::Format("diFilt%dDimuon_recoPtDR", int(f + 1))]->Fill(globDimuonPt, globDimuonDR);
0695           if (highPt)
0696             ME[TString::Format("diFilt%dDimuon_recoPtDRpos", int(f + 1))]->Fill(globDimuonPt, globDimuonDRpos);
0697         } else {
0698           break;
0699         }
0700       }
0701       //one filter object
0702       for (size_t f = 0; f < filterNamesLevels.size(); f++) {
0703         if (filt_glob[f][glob_gen[0]] != -1 || filt_glob[f][glob_gen[secondMuon]] != -1) {
0704           ME[TString::Format("filt%dDimuon_recoEtaPt", int(f + 1))]->Fill(globDimuonEta, globDimuonPt);
0705           ME[TString::Format("filt%dDimuon_recoRapPt", int(f + 1))]->Fill(globDimuonRap, globDimuonPt);
0706           if (highPt)
0707             ME[TString::Format("filt%dDimuon_recoPtDR", int(f + 1))]->Fill(globDimuonPt, globDimuonDR);
0708           if (highPt)
0709             ME[TString::Format("filt%dDimuon_recoPtDRpos", int(f + 1))]->Fill(globDimuonPt, globDimuonDRpos);
0710         } else {
0711           break;
0712         }
0713       }
0714       //two path objects
0715       if (path_glob[glob_gen[0]] != -1 && path_glob[glob_gen[secondMuon]] != -1) {
0716         ME["diPathDimuon_recoEtaPt"]->Fill(globDimuonEta, globDimuonPt);
0717         ME["diPathDimuon_recoRapPt"]->Fill(globDimuonRap, globDimuonPt);
0718         if (highPt)
0719           ME["diPathDimuon_recoPtDR"]->Fill(globDimuonPt, globDimuonDR);
0720         if (highPt)
0721           ME["diPathDimuon_recoPtDRpos"]->Fill(globDimuonPt, globDimuonDRpos);
0722       }
0723       //one path object
0724       if (path_glob[glob_gen[0]] != -1 || path_glob[glob_gen[secondMuon]] != -1) {
0725         ME["pathDimuon_recoEtaPt"]->Fill(globDimuonEta, globDimuonPt);
0726         ME["pathDimuon_recoRapPt"]->Fill(globDimuonRap, globDimuonPt);
0727         if (highPt)
0728           ME["pathDimuon_recoPtDR"]->Fill(globDimuonPt, globDimuonDR);
0729         if (highPt)
0730           ME["pathDimuon_recoPtDRpos"]->Fill(globDimuonPt, globDimuonDRpos);
0731       }
0732       //trigger result
0733       if (triggerFired) {
0734         ME["resultDimuon_recoEtaPt"]->Fill(globDimuonEta, globDimuonPt);
0735         ME["resultDimuon_recoRapPt"]->Fill(globDimuonRap, globDimuonPt);
0736         if (highPt)
0737           ME["resultDimuon_recoPtDR"]->Fill(globDimuonPt, globDimuonDR);
0738         if (highPt)
0739           ME["resultDimuon_recoPtDRpos"]->Fill(globDimuonPt, globDimuonDRpos);
0740       }
0741     }
0742   }
0743 }
0744 
0745 int HeavyFlavorValidation::getMotherId(const Candidate *p) {
0746   const Candidate *mother = p->mother();
0747   if (mother) {
0748     if (mother->pdgId() == p->pdgId()) {
0749       return getMotherId(mother);
0750     } else {
0751       return mother->pdgId();
0752     }
0753   } else {
0754     return 0;
0755   }
0756 }
0757 
0758 void HeavyFlavorValidation::match(MonitorElement *me,
0759                                   vector<LeafCandidate> &from,
0760                                   vector<LeafCandidate> &to,
0761                                   double dRMatchingCut,
0762                                   vector<int> &map) {
0763   vector<double> dR(from.size());
0764   for (size_t i = 0; i < from.size(); i++) {
0765     map[i] = -1;
0766     dR[i] = 10.;
0767     //find closest
0768     for (size_t j = 0; j < to.size(); j++) {
0769       double dRtmp = deltaR<double>(from[i].eta(), from[i].phi(), to[j].eta(), to[j].phi());
0770       if (dRtmp < dR[i]) {
0771         dR[i] = dRtmp;
0772         map[i] = j;
0773       }
0774     }
0775     //fill matching histo
0776     if (map[i] != -1) {
0777       me->Fill(to[map[i]].eta() - from[i].eta(), deltaPhi<double>(to[map[i]].phi(), from[i].phi()));
0778     }
0779     //apply matching cut
0780     if (dR[i] > dRMatchingCut) {
0781       map[i] = -1;
0782     }
0783     //remove duplication
0784     if (map[i] != -1) {
0785       for (size_t k = 0; k < i; k++) {
0786         if (map[k] != -1 && map[i] == map[k]) {
0787           if (dR[i] < dR[k]) {
0788             map[k] = -1;
0789           } else {
0790             map[i] = -1;
0791           }
0792           break;
0793         }
0794       }
0795     }
0796   }
0797 }
0798 
0799 void HeavyFlavorValidation::myBook2D(DQMStore::IBooker &ibooker,
0800                                      TString name,
0801                                      vector<double> &ptBins,
0802                                      TString ptLabel,
0803                                      vector<double> &etaBins,
0804                                      TString etaLabel,
0805                                      TString title) {
0806   //   dqmStore->setCurrentFolder(dqmFolder+"/"+folder);
0807   int ptN = ptBins.size() == 3 ? (int)ptBins[0] + 1 : ptBins.size();
0808   Double_t *pt = new Double_t[ptN];
0809   for (int i = 0; i < ptN; i++) {
0810     pt[i] = ptBins.size() == 3 ? ptBins[1] + i * (ptBins[2] - ptBins[1]) / ptBins[0] : ptBins[i];
0811   }
0812   int etaN = etaBins.size() == 3 ? (int)etaBins[0] + 1 : etaBins.size();
0813   Double_t *eta = new Double_t[etaN];
0814   for (int i = 0; i < etaN; i++) {
0815     eta[i] = etaBins.size() == 3 ? etaBins[1] + i * (etaBins[2] - etaBins[1]) / etaBins[0] : etaBins[i];
0816   }
0817   TH2F *h = new TH2F(name, name, ptN - 1, pt, etaN - 1, eta);
0818   h->SetXTitle(ptLabel);
0819   h->SetYTitle(etaLabel);
0820   h->SetTitle(title);
0821   ME[name] = ibooker.book2D(name.Data(), h);
0822   delete h;
0823 }
0824 
0825 void HeavyFlavorValidation::myBookProfile2D(DQMStore::IBooker &ibooker,
0826                                             TString name,
0827                                             vector<double> &ptBins,
0828                                             TString ptLabel,
0829                                             vector<double> &etaBins,
0830                                             TString etaLabel,
0831                                             TString title) {
0832   //   dqmStore->setCurrentFolder(dqmFolder+"/"+folder);
0833   int ptN = ptBins.size() == 3 ? (int)ptBins[0] + 1 : ptBins.size();
0834   Double_t *pt = new Double_t[ptN];
0835   for (int i = 0; i < ptN; i++) {
0836     pt[i] = ptBins.size() == 3 ? ptBins[1] + i * (ptBins[2] - ptBins[1]) / ptBins[0] : ptBins[i];
0837   }
0838   int etaN = etaBins.size() == 3 ? (int)etaBins[0] + 1 : etaBins.size();
0839   Double_t *eta = new Double_t[etaN];
0840   for (int i = 0; i < etaN; i++) {
0841     eta[i] = etaBins.size() == 3 ? etaBins[1] + i * (etaBins[2] - etaBins[1]) / etaBins[0] : etaBins[i];
0842   }
0843   TProfile2D *h = new TProfile2D(name, name, ptN - 1, pt, etaN - 1, eta);
0844   h->SetXTitle(ptLabel);
0845   h->SetYTitle(etaLabel);
0846   h->SetTitle(title);
0847   ME[name] = ibooker.bookProfile2D(name.Data(), h);
0848   delete h;
0849 }
0850 
0851 void HeavyFlavorValidation::myBook1D(
0852     DQMStore::IBooker &ibooker, TString name, vector<double> &bins, TString label, TString title) {
0853   //   dqmStore->setCurrentFolder(dqmFolder+"/"+folder);
0854   int binsN = bins.size() == 3 ? (int)bins[0] + 1 : bins.size();
0855   Double_t *myBins = new Double_t[binsN];
0856   for (int i = 0; i < binsN; i++) {
0857     myBins[i] = bins.size() == 3 ? bins[1] + i * (bins[2] - bins[1]) / bins[0] : bins[i];
0858   }
0859   TH1F *h = new TH1F(name, name, binsN - 1, myBins);
0860   h->SetXTitle(label);
0861   h->SetTitle(title);
0862   ME[name] = ibooker.book1D(name.Data(), h);
0863   delete h;
0864 }
0865 
0866 int HeavyFlavorValidation::getFilterLevel(const std::string &moduleName, const HLTConfigProvider &hltConfig) {
0867   // helper lambda to check if a string contains a substring
0868   const auto contains = [](const std::string &s, const std::string &sub) -> bool {
0869     return s.find(sub) != std::string::npos;
0870   };
0871 
0872   // helper lambda to check if a string contains any of a list of substrings
0873   const auto containsAny = [](const std::string &s, const std::vector<std::string> &subs) -> bool {
0874     for (const auto &sub : subs) {
0875       if (s.find(sub) != std::string::npos)
0876         return true;
0877     }
0878     return false;
0879   };
0880 
0881   // helper lambda to check if string s is any of the strings in vector ms
0882   const auto isAnyOf = [](const std::string &s, const std::vector<std::string> &ms) -> bool {
0883     for (const auto &m : ms) {
0884       if (s == m)
0885         return true;
0886     }
0887     return false;
0888   };
0889 
0890   // tmadlener, 20.08.2017:
0891   // define the valid module names for the different "levels", to add a little bit more stability
0892   // to the checking compared to just doing some name matching.
0893   // Note, that the name matching is not completely remved, since at level 4 and 5 some of the
0894   // valid modules are the same, so that the name matching is still needed.
0895   // With the current definition this yields the exact same levels as before, but weeds out some
0896   // of the "false" positives at level 3 (naming matches also to some HLTMuonL1TFilter modules due to
0897   // the 'forIterL3' in the name)
0898   const std::string l1Filter = "HLTMuonL1TFilter";
0899   const std::string l2Filter = "HLTMuonL2FromL1TPreFilter";
0900   const std::vector<std::string> l3Filters = {"HLTMuonDimuonL3Filter", "HLTMuonL3PreFilter"};
0901   const std::vector<std::string> l4Filters = {
0902       "HLTDisplacedmumuFilter", "HLTDiMuonGlbTrkFilter", "HLTMuonTrackMassFilter"};
0903   const std::vector<std::string> l5Filters = {"HLTmumutkFilter", "HLT2MuonMuonDZ", "HLTDisplacedmumuFilter"};
0904 
0905   if (contains(moduleName, "Filter") && hltConfig.moduleEDMType(moduleName) == "EDFilter") {
0906     if (contains(moduleName, "L1") && !contains(moduleName, "ForIterL3") &&
0907         hltConfig.moduleType(moduleName) == l1Filter) {
0908       return 1;
0909     }
0910     if (contains(moduleName, "L2") && hltConfig.moduleType(moduleName) == l2Filter) {
0911       return 2;
0912     }
0913     if (contains(moduleName, "L3") && isAnyOf(hltConfig.moduleType(moduleName), l3Filters)) {
0914       return 3;
0915     }
0916     if (containsAny(moduleName, {"DisplacedmumuFilter", "DiMuon", "MuonL3Filtered", "TrackMassFiltered"}) &&
0917         isAnyOf(hltConfig.moduleType(moduleName), l4Filters)) {
0918       return 4;
0919     }
0920     if (containsAny(moduleName, {"Vertex", "Dz"}) && isAnyOf(hltConfig.moduleType(moduleName), l5Filters)) {
0921       return 5;
0922     }
0923   }
0924 
0925   return -1;
0926 }
0927 
0928 HeavyFlavorValidation::~HeavyFlavorValidation() {}
0929 
0930 //define this as a plug-in
0931 DEFINE_FWK_MODULE(HeavyFlavorValidation);