HLTObjectsMonitor.cc

Go to the documentation of this file.

// system include files
#include <memory>
#include <sys/time.h>
#include <cstdlib>
#include <unordered_map>

// user include files
#include "DataFormats/Common/interface/TriggerResults.h"
#include "DataFormats/HLTReco/interface/TriggerObject.h"
#include "DataFormats/HLTReco/interface/TriggerEvent.h"
#include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"

#include "FWCore/Utilities/interface/EDGetToken.h"
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include <DQMServices/Core/interface/DQMEDAnalyzer.h>

#include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"

//for collections
#include "HLTrigger/JetMET/interface/AlphaT.h"

#include "DataFormats/Math/interface/deltaR.h"
#include "DataFormats/BeamSpot/interface/BeamSpot.h"
#include "DataFormats/BTauReco/interface/JetTag.h"
#include "DataFormats/METReco/interface/MET.h"
#include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
#include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
#include "DataFormats/JetReco/interface/PFJet.h"
#include "DataFormats/JetReco/interface/CaloJet.h"

#include "DQMServices/Core/interface/DQMStore.h"
#include "DQMServices/Core/interface/DQMEDAnalyzer.h"

#include "TLorentzVector.h"

struct hltPlot {
  typedef dqm::reco::MonitorElement MonitorElement;

  std::pair<MonitorElement*, bool> nME;
  std::pair<MonitorElement*, bool> etaME;
  std::pair<MonitorElement*, bool> phiME;
  std::pair<MonitorElement*, bool> ptME;
  std::pair<MonitorElement*, bool> massME;
  std::pair<MonitorElement*, bool> energyME;
  std::pair<MonitorElement*, bool> csvME;
  std::pair<MonitorElement*, bool> etaVSphiME;
  std::pair<MonitorElement*, bool> ptMEhep17;
  std::pair<MonitorElement*, bool> ptMEhem17;  // in harvesting step ratio
  std::pair<MonitorElement*, bool> mrME;
  std::pair<MonitorElement*, bool> rsqME;
  std::pair<MonitorElement*, bool> dxyME;
  std::pair<MonitorElement*, bool> dzME;
  std::pair<MonitorElement*, bool> dimassME;
  std::pair<MonitorElement*, bool> dRME;
  std::pair<MonitorElement*, bool> dRetaVSphiME;
  std::pair<MonitorElement*, bool> q1q2ME;

  std::string label;
  std::string pathNAME;
  int pathIDX;
  std::string moduleNAME;

  std::string xTITLE;
  std::vector<double> etaBINNING;
  std::vector<double> ptBINNING;
  std::vector<double> phiBINNING;
  std::vector<double> massBINNING;
  std::vector<double> dxyBINNING;
  std::vector<double> dzBINNING;
  std::vector<double> dimassBINNING;

  bool doPlot2D;
  bool doPlotETA;
  bool doPlotMASS;
  bool doPlotENERGY;
  bool doPlotHEP17;
  bool doPlotCSV;
  bool doCALO;
  bool doPF;
  bool doPlotRazor;
  bool doPlotDXY;
  bool doPlotDZ;
  bool doPlotDiMass;
};
//
// class declaration
//

//using namespace edm;
using std::unordered_map;

class HLTObjectsMonitor : public DQMEDAnalyzer {
public:
  explicit HLTObjectsMonitor(const edm::ParameterSet&);
  ~HLTObjectsMonitor() override = default;

  //      static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);

private:
  void analyze(const edm::Event&, const edm::EventSetup&) override;
  void bookHistograms(DQMStore::IBooker& i, edm::Run const&, edm::EventSetup const&) override;
  void dqmBeginRun(edm::Run const&, edm::EventSetup const&) override;

  static hltPlot getPlotPSet(edm::ParameterSet pset);
  void getPSet();
  bool isHEP17(double eta, double phi);
  bool isHEM17(double eta, double phi);

  double dxyFinder(
      double, double, edm::Handle<reco::RecoChargedCandidateCollection>, edm::Handle<reco::BeamSpot>, double);
  double dzFinder(double, double, double, double, edm::Handle<reco::RecoChargedCandidateCollection>, double);
  // ----------member data ---------------------------

  std::string TopFolder_;
  std::string label_;
  std::string processName_;
  std::vector<edm::ParameterSet> plotPSETS_;

  HLTConfigProvider hltConfig_;

  MonitorElement* eventsPlot_;
  std::vector<hltPlot> hltPlots_;
  std::unordered_map<std::string, int> plotNamesToBins_;

  bool debug_;

  std::string mainFolder_;
  std::string backupFolder_;

  edm::EDGetTokenT<edm::TriggerResults> triggerResultsToken_;
  edm::EDGetTokenT<trigger::TriggerEvent> triggerEventToken_;

  edm::InputTag beamSpot_;
  edm::EDGetTokenT<reco::BeamSpot> beamSpotToken_;

  edm::EDGetTokenT<reco::JetTagCollection> caloJetBTagsToken_;
  edm::EDGetTokenT<reco::JetTagCollection> pfJetBTagsToken_;

  edm::InputTag muCandidates_;
  edm::EDGetTokenT<std::vector<reco::RecoChargedCandidate>> muCandidatesToken_;
  edm::InputTag eleCandidates_;
  edm::EDGetTokenT<std::vector<reco::RecoChargedCandidate>> eleCandidatesToken_;

  const double MASS_MU = .105658;

  struct MEbinning {
    int nbins;
    double xmin;
    double xmax;
  };

  double MAX_PHI = 3.2;
  int N_PHI = 64;
  const MEbinning phi_binning_{N_PHI, -MAX_PHI, MAX_PHI};

  double MAX_CSV = 1.;
  int N_CSV = 20;
  const MEbinning csv_binning_{N_CSV, -MAX_CSV, MAX_CSV};

  std::vector<double> phi_variable_binning_;

  /*
  HEP17 covers 
  - phi between 310° and 330° (-50° to -30°, or -0.87 t.52 rad)
  - eta between +1.3 and +3.0 (positive side only)
*/
  double MAX_PHI_HEP17 = -0.52;
  double MIN_PHI_HEP17 = -0.87;
  int N_PHI_HEP17 = 7;
  const MEbinning phi_binning_hep17_{N_PHI_HEP17, MIN_PHI_HEP17, MAX_PHI_HEP17};
  double MAX_ETA_HEP17 = 3.0;
  double MIN_ETA_HEP17 = 1.3;
  int N_ETA_HEP17 = 6;
  const MEbinning eta_binning_hep17_{N_ETA_HEP17, MIN_ETA_HEP17, MAX_ETA_HEP17};

  const MEbinning eta_binning_hem17_{N_ETA_HEP17, -MAX_ETA_HEP17, MIN_ETA_HEP17};
};

//
// constants, enums and typedefs
//

//
// static data member definitions
//
hltPlot HLTObjectsMonitor::getPlotPSet(edm::ParameterSet pset) {
  return hltPlot{std::make_pair<MonitorElement*, bool>(nullptr, false),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_eta")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_phi")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_pt")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_mass")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_energy")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_csv")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_etaVSphi")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_pt_HEP17")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_pt_HEM17")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_MR")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_RSQ")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_dxy")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_dz")),
                 std::make_pair<MonitorElement*, bool>(nullptr, pset.getParameter<bool>("displayInPrimary_dimass")),
                 std::make_pair<MonitorElement*, bool>(nullptr, false),
                 std::make_pair<MonitorElement*, bool>(nullptr, false),
                 std::make_pair<MonitorElement*, bool>(nullptr, false),
                 pset.getParameter<std::string>("label"),
                 pset.getParameter<std::string>("pathNAME"),
                 -1,
                 pset.getParameter<std::string>("moduleNAME"),
                 pset.getParameter<std::string>("xTITLE"),
                 pset.getParameter<std::vector<double>>("etaBINNING"),
                 pset.getParameter<std::vector<double>>("ptBINNING"),
                 pset.getParameter<std::vector<double>>("phiBINNING"),
                 pset.getParameter<std::vector<double>>("massBINNING"),
                 pset.getParameter<std::vector<double>>("dxyBINNING"),
                 pset.getParameter<std::vector<double>>("dzBINNING"),
                 pset.getParameter<std::vector<double>>("dimassBINNING"),
                 pset.getUntrackedParameter<bool>("doPlot2D", false),
                 pset.getUntrackedParameter<bool>("doPlotETA", true),
                 pset.getUntrackedParameter<bool>("doPlotMASS", false),
                 pset.getUntrackedParameter<bool>("doPlotENERGY", false),
                 pset.getUntrackedParameter<bool>("doPlotHEP17", true),
                 pset.getUntrackedParameter<bool>("doPlotCSV", false),
                 pset.getUntrackedParameter<bool>("doCALO", false),
                 pset.getUntrackedParameter<bool>("doPF", false),
                 pset.getUntrackedParameter<bool>("doPlotRazor", false),
                 pset.getUntrackedParameter<bool>("doPlotDXY", false),
                 pset.getUntrackedParameter<bool>("doPlotDZ", false),
                 pset.getUntrackedParameter<bool>("doPlotDiMass", false)};
}

void HLTObjectsMonitor::getPSet() {
  for (const auto& pset : plotPSETS_)
    hltPlots_.push_back(getPlotPSet(pset));
}

bool HLTObjectsMonitor::isHEP17(double eta, double phi) {
  if ((eta >= MIN_ETA_HEP17 && eta <= MAX_ETA_HEP17) && (phi >= MIN_PHI_HEP17 && phi <= MAX_PHI_HEP17))
    return true;
  else
    return false;
}
bool HLTObjectsMonitor::isHEM17(double eta, double phi) {
  if ((eta >= -MAX_ETA_HEP17 && eta <= -MIN_ETA_HEP17) && (phi >= MIN_PHI_HEP17 && phi <= MAX_PHI_HEP17))
    return true;
  else
    return false;
}
//
// constructors and destructor
//
HLTObjectsMonitor::HLTObjectsMonitor(const edm::ParameterSet& iConfig)
    : TopFolder_(iConfig.getParameter<std::string>("TopFolder")),
      label_(iConfig.getParameter<std::string>("label")),
      processName_(iConfig.getParameter<std::string>("processName")),
      plotPSETS_(iConfig.getParameter<std::vector<edm::ParameterSet>>("plots")),
      debug_(iConfig.getUntrackedParameter<bool>("debug", false)),
      triggerResultsToken_(consumes<edm::TriggerResults>(iConfig.getParameter<edm::InputTag>("TriggerResults"))),
      triggerEventToken_(consumes<trigger::TriggerEvent>(iConfig.getParameter<edm::InputTag>("TriggerSummary"))),
      beamSpot_(iConfig.getParameter<edm::InputTag>("beamspot")),
      beamSpotToken_(consumes<reco::BeamSpot>(beamSpot_)),
      caloJetBTagsToken_(consumes<reco::JetTagCollection>(iConfig.getParameter<edm::InputTag>("caloJetBTags"))),
      pfJetBTagsToken_(consumes<reco::JetTagCollection>(iConfig.getParameter<edm::InputTag>("pfJetBTags"))),
      muCandidates_(iConfig.getParameter<edm::InputTag>("muCandidates")),
      muCandidatesToken_(consumes<std::vector<reco::RecoChargedCandidate>>(muCandidates_)),
      eleCandidates_(iConfig.getParameter<edm::InputTag>("eleCandidates")),
      eleCandidatesToken_(consumes<std::vector<reco::RecoChargedCandidate>>(eleCandidates_)) {
  getPSet();

  //now do what ever initialization is needed
  mainFolder_ = TopFolder_ + "/MainShifter";
  backupFolder_ = TopFolder_ + "/Backup";

  //set Token(s)

  double step = 2 * MAX_PHI / double(N_PHI);
  for (int i = 0; i <= N_PHI; i++)
    phi_variable_binning_.push_back(-MAX_PHI + step * i);
}

//
// member functions
//

// ------------ method called for each event  ------------
void HLTObjectsMonitor::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
  //  if ( debug_ )
  //    std::cout << "[HLTObjectsMonitor::analyze]" << std::endl;

  // access trigger results
  edm::Handle<edm::TriggerResults> triggerResults;
  iEvent.getByToken(triggerResultsToken_, triggerResults);
  if (!triggerResults.isValid())
    return;

  edm::Handle<trigger::TriggerEvent> triggerEvent;
  iEvent.getByToken(triggerEventToken_, triggerEvent);
  if (!triggerEvent.isValid())
    return;

  edm::Handle<reco::JetTagCollection> caloJetBTags;
  iEvent.getByToken(caloJetBTagsToken_, caloJetBTags);

  edm::Handle<reco::JetTagCollection> pfJetBTags;
  iEvent.getByToken(pfJetBTagsToken_, pfJetBTags);

  edm::Handle<std::vector<reco::RecoChargedCandidate>> muCandidates;
  iEvent.getByToken(muCandidatesToken_, muCandidates);

  edm::Handle<std::vector<reco::RecoChargedCandidate>> eleCandidates;
  iEvent.getByToken(eleCandidatesToken_, eleCandidates);

  edm::Handle<reco::BeamSpot> beamspot;
  iEvent.getByToken(beamSpotToken_, beamspot);

  // loop over path
  int ibin = -1;
  std::vector<bool> plottedPathIndices(plotNamesToBins_.size(), false);
  for (auto& plot : hltPlots_) {
    ibin++;
    if (plot.pathIDX <= 0)
      continue;

    if (triggerResults->accept(plot.pathIDX)) {
      //We only want to fill this once per pathNAME per event
      auto index = plotNamesToBins_[plot.pathNAME];
      if (not plottedPathIndices[index]) {
        plottedPathIndices[index] = true;
        if (debug_)
          std::cout << plot.pathNAME << " --> bin: " << ibin << std::endl;
        eventsPlot_->Fill(index);
      }
      const trigger::TriggerObjectCollection objects = triggerEvent->getObjects();
      if (hltConfig_.saveTags(plot.moduleNAME)) {
        if (debug_)
          std::cout << "objects: " << objects.size() << std::endl;

        bool moduleFOUND = false;
        std::vector<std::string> moduleNames = hltConfig_.moduleLabels(plot.pathIDX);
        for (const auto& module : moduleNames) {
          if (module == plot.moduleNAME)
            moduleFOUND = true;
        }
        if (debug_)
          std::cout << plot.moduleNAME << (moduleFOUND ? "" : "NOT") << " found in the list of modules" << std::endl;

        if (debug_)
          for (const auto& module : moduleNames) {
            unsigned int idx = triggerEvent->filterIndex(edm::InputTag(module, "", processName_));
            std::cout << "module: " << module;
            if (idx < triggerEvent->sizeFilters())
              std::cout << " --> " << idx;
            std::cout << std::endl;
          }
        //
        // trigger accepted and collection w/ objects is available
        edm::InputTag moduleName = edm::InputTag(plot.moduleNAME, "", processName_);
        unsigned int moduleIDX = triggerEvent->filterIndex(moduleName);
        if (debug_)
          std::cout << "moduleNAME: " << plot.moduleNAME << " --> " << moduleIDX << std::endl;

        if (moduleIDX >= triggerEvent->sizeFilters()) {
          LogDebug("HLTObjectsMonitor")
              << plot.pathNAME << " " << plot.moduleNAME
              << " is not available ! please, fix update DQM/HLTEvF/python/HLTObjectsMonitor_cfi.py";
          return;
        }

        const trigger::Keys& keys = triggerEvent->filterKeys(moduleIDX);
        if (debug_)
          std::cout << "keys: " << keys.size() << std::endl;

        plot.nME.first->Fill(keys.size());

        double MR = 0.;
        double RSQ = 0.;
        for (const auto& key : keys) {
          double pt = objects[key].pt();
          double eta = objects[key].eta();
          double phi = objects[key].phi();
          double mass = objects[key].mass();
          double energy = objects[key].energy();
          int id = objects[key].id();
          if (debug_)
            std::cout << "object ID " << id << " mass: " << mass << std::endl;

          // single-object plots
          plot.ptME.first->Fill(pt);
          if (plot.doPlotETA)
            plot.etaME.first->Fill(eta);
          plot.phiME.first->Fill(phi);

          if (plot.doPlotCSV) {
            if (plot.doCALO) {
              if (!caloJetBTags.isValid())
                plot.csvME.first->Fill(-10.);
              else {
                for (auto it = caloJetBTags->begin(); it != caloJetBTags->end(); ++it) {
                  double dR = deltaR(eta, phi, it->first->eta(), it->first->phi());
                  if (debug_)
                    std::cout << "[HLTObjectsMonitor::analyze] deltaR: " << dR << " matched ? "
                              << (dR <= 0.4 ? "YEAP" : "NOPE") << std::endl;
                  plot.csvME.first->Fill(it->second);
                }
              }

            } else if (plot.doPF) {
              if (!pfJetBTags.isValid())
                plot.csvME.first->Fill(-10.);
              else {
                for (auto it = pfJetBTags->begin(); it != pfJetBTags->end(); ++it) {
                  double dR = deltaR(eta, phi, it->first->eta(), it->first->phi());
                  if (debug_)
                    std::cout << "[HLTObjectsMonitor::analyze] deltaR: " << dR << " matched ? "
                              << (dR <= 0.4 ? "YEAP" : "NOPE") << std::endl;
                  plot.csvME.first->Fill(it->second);
                }
              }
            }
          }
          if (plot.doPlotMASS)
            plot.massME.first->Fill(mass);
          if (plot.doPlotENERGY)
            plot.energyME.first->Fill(energy);
          if (plot.doPlot2D)
            plot.etaVSphiME.first->Fill(eta, phi);
          if (plot.doPlotHEP17) {
            if (isHEP17(eta, phi))
              plot.ptMEhep17.first->Fill(pt);
            if (isHEM17(eta, phi))
              plot.ptMEhem17.first->Fill(pt);
          }

          if (id == 0) {             //the MET object containing MR and Rsq will show up with ID = 0
            MR = objects[key].px();  //razor variables stored in dummy reco::MET objects
            RSQ = objects[key].py();
          }

          if (plot.doPlotDXY) {
            double dxy = -99.;
            if (abs(id) == 13)
              dxy = dxyFinder(eta, phi, muCandidates, beamspot, 0.1);  // dRcut = 0.1
            else
              dxy = dxyFinder(eta, phi, eleCandidates, beamspot, 0.1);  // dRcut = 0.1
            plot.dxyME.first->Fill(dxy);
          }
        }  // end loop on keys
        if (plot.doPlotRazor) {
          plot.mrME.first->Fill(MR);
          plot.rsqME.first->Fill(RSQ);
        }

        if (keys.size() < 2) {
          if (plot.doPlotDiMass || plot.doPlotDZ)
            LogDebug("HLTObjectsMonitor") << plot.pathNAME << " " << plot.moduleNAME << " # objects is (" << keys.size()
                                          << ") less than 2 ! you probably want to either change the moduleNAME or "
                                             "switch off di-object system plots (doPlotDZ: "
                                          << plot.doPlotDZ << " doPlotDiMass: " << plot.doPlotDiMass
                                          << ") in DQM/HLTEvF/python/HLTObjectsMonitor_cfi.py)";
        } else {
          for (const auto& key : keys) {
            double pt = objects[key].pt();
            double eta = objects[key].eta();
            double phi = objects[key].phi();
            int id = objects[key].id();

            unsigned int kCnt0 = 0;

            TLorentzVector v1;
            if (abs(id) == 13)  // check if it is a muon
              v1.SetPtEtaPhiM(pt, eta, phi, MASS_MU);
            else
              v1.SetPtEtaPhiM(pt, eta, phi, 0);

            unsigned int kCnt1 = 0;
            for (const auto& key1 : keys) {
              if (key != key1 &&
                  kCnt1 > kCnt0) {  // avoid filling hists with same objs && avoid double counting separate objs

                double pt2 = objects[key1].pt();
                double eta2 = objects[key1].eta();
                double phi2 = objects[key1].phi();
                int id2 = objects[key1].id();

                double dR = deltaR(eta, phi, eta2, phi2);
                plot.dRME.first->Fill(dR);
                plot.dRetaVSphiME.first->Fill(eta, phi, dR);

                int q1 = (id == 0 ? 0 : id / abs(id));
                int q2 = (id2 == 0 ? 0 : id2 / abs(id2));
                int q1q2 = q1 * q2;
                plot.q1q2ME.first->Fill(q1q2);

                if ((id + id2) == 0) {  // check di-object system charge and flavor

                  TLorentzVector v2;
                  if (abs(id2) == 13)  // check if it is a muon
                    v2.SetPtEtaPhiM(pt2, eta2, phi2, MASS_MU);
                  else
                    v2.SetPtEtaPhiM(pt2, eta2, phi2, 0);

                  if (plot.doPlotDiMass) {
                    TLorentzVector v = v1 + v2;
                    plot.dimassME.first->Fill(v.M());
                  }

                  if (plot.doPlotDZ) {
                    double dz = -99.;
                    if (abs(id) == 13)
                      dz = dzFinder(eta, phi, eta2, phi2, muCandidates, 0.1);  // dRcut = 0.1
                    else
                      dz = dzFinder(eta, phi, eta2, phi2, eleCandidates, 0.1);  // dRcut = 0.1
                    plot.dzME.first->Fill(dz);
                  }
                }
              }
              kCnt1++;
            }
            kCnt0++;
          }
        }
      }
    }
  }
}

// ------------ method called when starting to processes a run  ------------
void HLTObjectsMonitor::dqmBeginRun(edm::Run const& iRun, edm::EventSetup const& iSetup) {
  bool changed = true;
  if (hltConfig_.init(iRun, iSetup, processName_, changed))
    if (debug_)
      std::cout << "[HLTObjectsMonitor::dqmBeginRun] extracting HLTconfig" << std::endl;

  //get path indicies from menu
  std::string pathName_noVersion;
  std::vector<std::string> triggerPaths = hltConfig_.triggerNames();

  if (debug_)
    std::cout << "[HLTObjectsMonitor::dqmBeginRun] triggerPaths: " << triggerPaths.size() << " <--> "
              << hltPlots_.size() << std::endl;

  for (const auto& pathName : triggerPaths) {
    if (debug_)
      std::cout << "[HLTObjectsMonitor::dqmBeginRun] " << pathName << std::endl;

    pathName_noVersion = hltConfig_.removeVersion(pathName);
    //    std::cout << "pathName: " << pathName << " --> " << pathName_noVersion << std::endl;
    for (auto& plot : hltPlots_) {
      if (plot.pathNAME == pathName_noVersion) {
        plot.pathIDX = hltConfig_.triggerIndex(pathName);
        // check that the index makes sense, otherwise force pathIDX = -1
        if (plot.pathIDX <= 0 || plot.pathIDX == int(triggerPaths.size()))
          plot.pathIDX = -1;
      }
    }
  }

  if (debug_) {
    for (const auto& plot : hltPlots_)
      std::cout << "plot: " << plot.pathNAME << " --> pathIDX: " << plot.pathIDX << std::endl;
    std::cout << "[HLTObjectsMonitor::dqmBeginRun] DONE" << std::endl;
  }
}

void HLTObjectsMonitor::bookHistograms(DQMStore::IBooker& ibooker,
                                       edm::Run const& iRun,
                                       edm::EventSetup const& iSetup) {
  if (debug_)
    std::cout << "[HLTObjectsMonitor::bookHistograms]" << std::endl;

  ibooker.setCurrentFolder(TopFolder_);

  std::string name = "eventsPerPath_" + label_;
  std::string title = " events per path";

  //We must avoid repeating the same pathNAME
  {
    std::unordered_map<std::string, int> uniqueNames;
    for (auto const& p : hltPlots_) {
      plotNamesToBins_[p.pathNAME] = -1;
    }
    int nbins = plotNamesToBins_.size();
    eventsPlot_ = ibooker.book1D(name, title, nbins, -0.5, double(nbins) - 0.5);
    eventsPlot_->setAxisTitle("HLT path");
    int i = 0;
    //keep the bin order the same as hltPlots_
    for (auto const& p : hltPlots_) {
      //only add a bin if this is the first time we've seen the name
      if (-1 == plotNamesToBins_[p.pathNAME]) {
        plotNamesToBins_[p.pathNAME] = ++i;
        eventsPlot_->setBinLabel(i, p.pathNAME);
        if (debug_)
          std::cout << p.pathNAME << " --> bin: " << i << std::endl;
      }
    }
  }

  for (auto& plot : hltPlots_) {
    if (debug_)
      std::cout << "booking plots for " << plot.label << std::endl;

    if (plot.pathIDX <= 0) {
      LogDebug("HLTObjectsMonitor") << plot.pathNAME
                                    << " is not available in the HLT menu ! no plots are going to be booked for it "
                                       "(update DQM/HLTEvF/python/HLTObjectsMonitor_cfi.py)";
      continue;
    }
    if (debug_)
      std::cout << "booking histograms for " << plot.pathNAME << std::endl;

    {
      if (plot.nME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_nobjects";
      title = plot.pathNAME + " # objects";
      plot.nME.first = ibooker.book1D(name, title, 20, -0.5, 19.5);
      plot.nME.first->setAxisTitle(plot.xTITLE + " # objects");
    }

    if (plot.ptME.second)
      ibooker.setCurrentFolder(mainFolder_);
    else
      ibooker.setCurrentFolder(backupFolder_);

    name = plot.label + "_pt";
    title = plot.pathNAME + " p_T";
    int nbins = (plot.ptBINNING).size() - 1;
    std::vector<float> fbinning((plot.ptBINNING).begin(), (plot.ptBINNING).end());
    float* arr = &fbinning[0];
    plot.ptME.first = ibooker.book1D(name, title, nbins, arr);
    plot.ptME.first->setAxisTitle(plot.xTITLE + " p_{T} [GeV]");

    {
      if (plot.phiME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_phi";
      title = plot.pathNAME + " #phi";
      int nbins = (plot.phiBINNING).size() - 1;
      std::vector<float> fbinning((plot.phiBINNING).begin(), (plot.phiBINNING).end());
      float* arr = &fbinning[0];
      plot.phiME.first = ibooker.book1D(name, title, nbins, arr);
      plot.phiME.first->setAxisTitle(plot.xTITLE + " #phi [rad]");
    }

    if (plot.doPlotETA) {
      if (plot.etaME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_eta";
      title = plot.pathNAME + " #eta";
      int nbins = (plot.etaBINNING).size() - 1;
      std::vector<float> fbinning((plot.etaBINNING).begin(), (plot.etaBINNING).end());
      float* arr = &fbinning[0];

      plot.etaME.first = ibooker.book1D(name, title, nbins, arr);
      plot.etaME.first->setAxisTitle(plot.xTITLE + " #eta");
    }

    if (plot.doPlotMASS) {
      if (plot.massME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_mass";
      title = plot.pathNAME + " mass";
      int nbins = (plot.massBINNING).size() - 1;
      std::vector<float> fbinning((plot.massBINNING).begin(), (plot.massBINNING).end());
      float* arr = &fbinning[0];

      plot.massME.first = ibooker.book1D(name, title, nbins, arr);
      plot.massME.first->setAxisTitle(plot.xTITLE + " mass [GeV]");
    }

    if (plot.doPlotENERGY) {
      if (plot.energyME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_energy";
      title = plot.pathNAME + " energy";
      int nbins = (plot.ptBINNING).size() - 1;
      std::vector<float> fbinning((plot.ptBINNING).begin(), (plot.ptBINNING).end());
      float* arr = &fbinning[0];

      plot.energyME.first = ibooker.book1D(name, title, nbins, arr);
      plot.energyME.first->setAxisTitle(plot.xTITLE + " energy [GeV]");
    }

    if (plot.doPlotCSV) {
      if (plot.csvME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_csv";
      title = plot.pathNAME + " CSV";

      plot.csvME.first = ibooker.book1D(name, title, csv_binning_.nbins, csv_binning_.xmin, csv_binning_.xmax);
      plot.csvME.first->setAxisTitle(plot.xTITLE + " CSV discriminator");
    }

    if (plot.doPlot2D) {
      if (plot.etaVSphiME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_etaVSphi";
      title = plot.pathNAME + " #eta vs #phi";
      int nbinsX = (plot.etaBINNING).size() - 1;
      std::vector<float> fbinningX((plot.etaBINNING).begin(), (plot.etaBINNING).end());
      float* arrX = &fbinningX[0];
      int nbinsY = (plot.phiBINNING).size() - 1;
      ;
      std::vector<float> fbinningY((plot.phiBINNING).begin(), (plot.phiBINNING).end());
      float* arrY = &fbinningY[0];
      plot.etaVSphiME.first = ibooker.book2D(name, title, nbinsX, arrX, nbinsY, arrY);
      plot.etaVSphiME.first->setAxisTitle(plot.xTITLE + " #eta", 1);
      plot.etaVSphiME.first->setAxisTitle(plot.xTITLE + " #phi", 2);
    }

    if (plot.doPlotHEP17) {
      if (plot.ptMEhep17.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      int nbins = (plot.ptBINNING).size() - 1;
      std::vector<float> fbinning((plot.ptBINNING).begin(), (plot.ptBINNING).end());
      float* arr = &fbinning[0];

      name = plot.label + "_pt_HEP17";
      title = plot.pathNAME + " p_{T} HEP17";
      plot.ptMEhep17.first = ibooker.book1D(name, title, nbins, arr);
      plot.ptMEhep17.first->setAxisTitle(plot.xTITLE + " p_{T} [GeV]", 1);

      if (plot.ptMEhem17.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_pt_HEM17";
      title = plot.pathNAME + " p_{T} HEM17";
      plot.ptMEhem17.first = ibooker.book1D(name, title, nbins, arr);
      plot.ptMEhem17.first->setAxisTitle(plot.xTITLE + " p_{T} [GeV]", 1);
    }

    if (plot.doPlotRazor) {
      if (plot.mrME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_mr";
      title = plot.pathNAME + " M_{R}";
      plot.mrME.first = ibooker.book1D(name, title, 100, 0., 100.);
      plot.mrME.first->setAxisTitle(plot.xTITLE + " M_{R} [GeV]", 1);

      if (plot.rsqME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_rsq";
      title = plot.pathNAME + " RSQ";
      plot.rsqME.first = ibooker.book1D(name, title, 100, 0., 100.);
      plot.rsqME.first->setAxisTitle(plot.xTITLE + " RSQ [GeV]", 1);
    }

    if (plot.doPlotDXY) {
      if (plot.dxyME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_dxy";
      title = plot.pathNAME + " d_{xy}";
      int nbins = (plot.dxyBINNING).size() - 1;
      std::vector<float> fbinning((plot.dxyBINNING).begin(), (plot.dxyBINNING).end());
      float* arr = &fbinning[0];
      plot.dxyME.first = ibooker.book1D(name, title, nbins, arr);
      plot.dxyME.first->setAxisTitle(plot.xTITLE + " d_{xy} [cm]");
    }

    if (plot.doPlotDZ) {
      if (plot.dzME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_dz";
      title = plot.pathNAME + " d_{z}";
      int nbins = (plot.dzBINNING).size() - 1;
      std::vector<float> fbinning((plot.dzBINNING).begin(), (plot.dzBINNING).end());
      float* arr = &fbinning[0];
      plot.dzME.first = ibooker.book1D(name, title, nbins, arr);
      plot.dzME.first->setAxisTitle(plot.xTITLE + " d_{z} [cm]");
    }

    if (plot.dRME.second)
      ibooker.setCurrentFolder(mainFolder_);
    else
      ibooker.setCurrentFolder(backupFolder_);

    name = plot.label + "_dR";
    title = plot.pathNAME + " di-object dR";
    plot.dRME.first = ibooker.book1D(name, title, 50, 0., 5.);
    plot.dRME.first->setAxisTitle(plot.xTITLE + " dR_{obj,obj}");

    if (plot.dRetaVSphiME.second)
      ibooker.setCurrentFolder(mainFolder_);
    else
      ibooker.setCurrentFolder(backupFolder_);

    name = plot.label + "_dR_etaVSphi";
    title = plot.pathNAME + " di-object dR in the #eta-#phi plane (of 1st obj)";
    plot.dRetaVSphiME.first = ibooker.bookProfile2D(name, title, 60, -3., 3., 64, -3.2, 3.2, 0., 5.);
    plot.dRetaVSphiME.first->setAxisTitle(plot.xTITLE + " #eta", 1);
    plot.dRetaVSphiME.first->setAxisTitle(plot.xTITLE + " #phi", 2);
    plot.dRetaVSphiME.first->setAxisTitle(plot.xTITLE + " dR_{obj,obj}", 3);

    if (plot.q1q2ME.second)
      ibooker.setCurrentFolder(mainFolder_);
    else
      ibooker.setCurrentFolder(backupFolder_);

    name = plot.label + "_q1q2";
    title = plot.pathNAME + " di-object q1xq2";
    plot.q1q2ME.first = ibooker.book1D(name, title, 3, -1., 1.);
    plot.q1q2ME.first->setAxisTitle(plot.xTITLE + " q_{obj1} x q_{obj2}");

    if (plot.doPlotDiMass) {
      if (plot.dimassME.second)
        ibooker.setCurrentFolder(mainFolder_);
      else
        ibooker.setCurrentFolder(backupFolder_);

      name = plot.label + "_dimass";
      title = plot.pathNAME + " di-object mass";
      int nbins = (plot.dimassBINNING).size() - 1;
      std::vector<float> fbinning((plot.dimassBINNING).begin(), (plot.dimassBINNING).end());
      float* arr = &fbinning[0];
      plot.dimassME.first = ibooker.book1D(name, title, nbins, arr);
      plot.dimassME.first->setAxisTitle(plot.xTITLE + " m_{obj,obj} [GeV]");
    }
  }

  if (debug_)
    std::cout << "[HLTObjectsMonitor::bookHistograms] DONE" << std::endl;
}

double HLTObjectsMonitor::dxyFinder(double eta,
                                    double phi,
                                    edm::Handle<reco::RecoChargedCandidateCollection> candidates,
                                    edm::Handle<reco::BeamSpot> beamspot,
                                    double dRcut = 0.1) {
  double dxy = -99.;
  if (!candidates.isValid()) {
    LogDebug("HLTObjectsMonitor") << "either " << muCandidates_ << " or " << eleCandidates_
                                  << " is not valid ! please, update DQM/HLTEvF/python/HLTObjectsMonitor_cfi.py"
                                  << " by switching OFF doPlotDXY or updating the InputTag collection";
    return dxy;
  }
  if (!beamspot.isValid()) {
    LogDebug("HLTObjectsMonitor") << beamSpot_
                                  << " is not valid ! please, update DQM/HLTEvF/python/HLTObjectsMonitor_cfi.py"
                                  << " by switching OFF doPlotDXY or updating the InputTag collection";
    return dxy;
  }

  bool matched = false;
  for (reco::RecoChargedCandidateCollection::const_iterator candidate = candidates->begin();
       candidate != candidates->end();
       ++candidate) {
    if (deltaR(eta, phi, candidate->eta(), candidate->phi()) < dRcut) {
      matched = true;
      dxy = (-(candidate->vx() - beamspot->x0()) * candidate->py() +
             (candidate->vy() - beamspot->y0()) * candidate->px()) /
            candidate->pt();
      break;
    }
  }
  if (!matched)
    edm::LogWarning("HLTObjectsMonitor") << "trigger object does not match ( dR > " << dRcut
                                         << ") to any of the candidates in either " << muCandidates_ << " or "
                                         << eleCandidates_;

  return dxy;
}

double HLTObjectsMonitor::dzFinder(double eta1,
                                   double phi1,
                                   double eta2,
                                   double phi2,
                                   edm::Handle<reco::RecoChargedCandidateCollection> candidates,
                                   double dRcut = 0.1) {
  double dz = -99.;
  if (!candidates.isValid()) {
    LogDebug("HLTObjectsMonitor") << "either " << muCandidates_ << " or " << eleCandidates_
                                  << " is not valid ! please, update DQM/HLTEvF/python/HLTObjectsMonitor_cfi.py"
                                  << " by switching OFF doPlotDZ or updating the InputTag collection";
    return dz;
  }

  const reco::RecoChargedCandidate* cand1;
  const reco::RecoChargedCandidate* cand2;
  bool matched1 = false;
  bool matched2 = false;
  for (reco::RecoChargedCandidateCollection::const_iterator candidate = candidates->begin();
       candidate != candidates->end();
       ++candidate) {
    if (deltaR(eta1, phi1, candidate->eta(), candidate->phi()) < dRcut) {
      matched1 = true;
      cand1 = &*candidate;
      if (debug_)
        std::cout << "cand1: " << cand1->pt() << " " << cand1->eta() << " " << cand1->phi() << std::endl;
      break;
    }
  }
  if (!matched1) {
    LogDebug("HLTObjectsMonitor") << "trigger object1 does not match ( dR > " << dRcut
                                  << ") to any of the candidates in either " << muCandidates_ << " or "
                                  << eleCandidates_;
    return dz;
  }

  for (reco::RecoChargedCandidateCollection::const_iterator candidate = candidates->begin();
       candidate != candidates->end();
       ++candidate) {
    if (debug_) {
      std::cout << "candidate: " << candidate->pt() << " cand1: " << cand1->pt() << std::endl;
      std::cout << "candidate: " << candidate->eta() << " cand1: " << cand1->eta() << std::endl;
      std::cout << "candidate: " << candidate->phi() << " cand1: " << cand1->phi() << std::endl;
    }
    if (&*candidate == cand1)
      continue;

    if (deltaR(eta2, phi2, candidate->eta(), candidate->phi()) < dRcut) {
      matched2 = true;
      cand2 = &*candidate;
      if (debug_)
        std::cout << "cand2: " << cand2->pt() << " " << cand2->eta() << " " << cand2->phi() << std::endl;
      break;
    }
  }
  if (!matched2) {
    LogDebug("HLTObjectsMonitor") << "trigger object2 does not match ( dR > " << dRcut
                                  << ") to any of the candidates in either " << muCandidates_ << " or "
                                  << eleCandidates_;
    return dz;
  }

  dz = cand1->vz() - cand2->vz();
  return dz;
}

//define this as a plug-in
DEFINE_FWK_MODULE(HLTObjectsMonitor);