Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:09:55

0001 // -*- C++ -*-
0002 //
0003 // Package:     HigPhotonJetHLTOfflineSource
0004 // Class:       HigPhotonJetHLTOfflineSource
0005 //
0006 
0007 //
0008 // Author: Xin Shi <Xin.Shi@cern.ch>
0009 // Created: 2014.07.22
0010 //
0011 
0012 // system include files
0013 #include <memory>
0014 #include <iostream>
0015 
0016 // user include files
0017 #include "DQMServices/Core/interface/DQMOneEDAnalyzer.h"
0018 
0019 #include "FWCore/Framework/interface/Frameworkfwd.h"
0020 #include "FWCore/Framework/interface/ConsumesCollector.h"
0021 #include "FWCore/Framework/interface/MakerMacros.h"
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0024 #include "FWCore/ServiceRegistry/interface/Service.h"
0025 #include "FWCore/Common/interface/TriggerNames.h"
0026 
0027 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0028 
0029 #include "DataFormats/Common/interface/TriggerResults.h"
0030 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0031 #include "DataFormats/VertexReco/interface/Vertex.h"
0032 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0033 #include "DataFormats/EgammaCandidates/interface/Photon.h"
0034 #include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
0035 #include "DataFormats/METReco/interface/PFMETCollection.h"
0036 #include "DataFormats/METReco/interface/PFMET.h"
0037 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0038 #include "DataFormats/JetReco/interface/PFJet.h"
0039 #include "DataFormats/Math/interface/deltaPhi.h"
0040 
0041 #include <TLorentzVector.h>
0042 #include <TH2F.h>
0043 
0044 //  Define the interface
0045 class HigPhotonJetHLTOfflineSource : public DQMOneEDAnalyzer<> {
0046 public:
0047   explicit HigPhotonJetHLTOfflineSource(const edm::ParameterSet&);
0048 
0049 private:
0050   // Analyzer Methods
0051   void dqmBeginRun(const edm::Run&, const edm::EventSetup&) override;
0052   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0053   void analyze(const edm::Event&, const edm::EventSetup&) override;
0054   void dqmEndRun(const edm::Run&, const edm::EventSetup&) override;
0055   bool isMonitoredTriggerAccepted(const edm::TriggerNames&, const edm::Handle<edm::TriggerResults>&);
0056 
0057   // Input from Configuration File
0058   edm::ParameterSet pset_;
0059   std::string hltProcessName_;
0060   std::vector<std::string> hltPathsToCheck_;
0061   std::string dirname_;
0062   bool verbose_;
0063   bool triggerAccept_;
0064   bool perLSsaving_;  //to avoid nanoDQMIO crashing, driven by  DQMServices/Core/python/DQMStore_cfi.py
0065 
0066   edm::EDGetTokenT<edm::TriggerResults> triggerResultsToken_;
0067   edm::EDGetTokenT<reco::VertexCollection> pvToken_;
0068   edm::EDGetTokenT<reco::PhotonCollection> photonsToken_;
0069   edm::EDGetTokenT<reco::PFMETCollection> pfMetToken_;
0070   edm::EDGetTokenT<reco::PFJetCollection> pfJetsToken_;
0071 
0072   double pfjetMinPt_;
0073   double photonMinPt_;
0074 
0075   // Member Variables
0076 
0077   MonitorElement* nvertices_reco_;
0078   MonitorElement* nvertices_;
0079   MonitorElement* nphotons_reco_;
0080   MonitorElement* nphotons_;
0081   MonitorElement* photonpt_reco_;
0082   MonitorElement* photonpt_;
0083   MonitorElement* photonrapidity_reco_;
0084   MonitorElement* photonrapidity_;
0085   MonitorElement* pfmet_reco_;
0086   MonitorElement* pfmet_;
0087   MonitorElement* pfmetphi_reco_;
0088   MonitorElement* pfmetphi_;
0089   MonitorElement* npfjets_reco_;
0090   MonitorElement* npfjets_;
0091   MonitorElement* delphiphomet_reco_;
0092   MonitorElement* delphiphomet_;
0093   MonitorElement* delphijetmet_reco_;
0094   MonitorElement* delphijetmet_;
0095   MonitorElement* invmassjj_reco_;
0096   MonitorElement* invmassjj_;
0097   MonitorElement* deletajj_reco_;
0098   MonitorElement* deletajj_;
0099   MonitorElement* triggers_reco_;
0100   MonitorElement* triggers_;
0101   MonitorElement* trigvsnvtx_reco_;
0102   MonitorElement* trigvsnvtx_;
0103 
0104   double evtsrun_;
0105 };
0106 
0107 // Class Methods
0108 
0109 HigPhotonJetHLTOfflineSource::HigPhotonJetHLTOfflineSource(const edm::ParameterSet& pset) : pset_(pset) {
0110   hltProcessName_ = pset.getParameter<std::string>("hltProcessName");
0111   hltPathsToCheck_ = pset.getParameter<std::vector<std::string>>("hltPathsToCheck");
0112   verbose_ = pset.getUntrackedParameter<bool>("verbose", false);
0113   perLSsaving_ = pset.getUntrackedParameter<bool>("perLSsaving", false);
0114   triggerAccept_ = pset.getUntrackedParameter<bool>("triggerAccept", true);
0115   triggerResultsToken_ = consumes<edm::TriggerResults>(pset.getParameter<edm::InputTag>("triggerResultsToken"));
0116   dirname_ = pset.getUntrackedParameter<std::string>("dirname", std::string("HLT/Higgs/PhotonJet/"));
0117   pvToken_ = consumes<reco::VertexCollection>(pset.getParameter<edm::InputTag>("pvToken"));
0118   photonsToken_ = consumes<reco::PhotonCollection>(pset.getParameter<edm::InputTag>("photonsToken"));
0119   pfMetToken_ = consumes<reco::PFMETCollection>(pset.getParameter<edm::InputTag>("pfMetToken"));
0120   pfJetsToken_ = consumes<reco::PFJetCollection>(pset.getParameter<edm::InputTag>("pfJetsToken"));
0121   pfjetMinPt_ = pset.getUntrackedParameter<double>("pfjetMinPt", 0.0);
0122   photonMinPt_ = pset.getUntrackedParameter<double>("photonMinPt", 0.0);
0123 }
0124 
0125 void HigPhotonJetHLTOfflineSource::dqmBeginRun(const edm::Run& iRun,
0126                                                const edm::EventSetup& iSetup) {  // Initialize hltConfig
0127   HLTConfigProvider hltConfig;
0128   bool changedConfig;
0129   if (!hltConfig.init(iRun, iSetup, hltProcessName_, changedConfig)) {
0130     edm::LogError("HLTPhotonJetVal") << "Initialization of HLTConfigProvider failed!!";
0131     return;
0132   }
0133 
0134   evtsrun_ = 0;
0135 }
0136 
0137 void HigPhotonJetHLTOfflineSource::bookHistograms(DQMStore::IBooker& iBooker,
0138                                                   edm::Run const& iRun,
0139                                                   edm::EventSetup const& iSetup) {
0140   iBooker.setCurrentFolder(dirname_);
0141   nvertices_reco_ = iBooker.book1D("nvertices_reco", "Reco: Number of vertices", 100, 0, 100);
0142   nvertices_ = iBooker.book1D("nvertices", "Number of vertices", 100, 0, 100);
0143   nphotons_reco_ = iBooker.book1D("nphotons_reco", "Reco: Number of photons", 100, 0, 10);
0144   nphotons_ = iBooker.book1D("nphotons", "Number of photons", 100, 0, 10);
0145   photonpt_reco_ = iBooker.book1D("photonpt_reco", "Reco: Photons pT", 100, 0, 500);
0146   photonpt_ = iBooker.book1D("photonpt", "Photons pT", 100, 0, 500);
0147   photonrapidity_reco_ = iBooker.book1D("photonrapidity_reco", "Reco: Photons rapidity;y_{#gamma}", 100, -2.5, 2.5);
0148   photonrapidity_ = iBooker.book1D("photonrapidity", "Photons rapidity;y_{#gamma}", 100, -2.5, 2.5);
0149   pfmet_reco_ = iBooker.book1D("pfmet_reco", "Reco: PF MET", 100, 0, 250);
0150   pfmet_ = iBooker.book1D("pfmet", "PF MET", 100, 0, 250);
0151   pfmetphi_reco_ = iBooker.book1D("pfmetphi_reco", "Reco: PF MET phi;#phi_{PFMET}", 100, -4, 4);
0152   pfmetphi_ = iBooker.book1D("pfmetphi", "PF MET phi;#phi_{PFMET}", 100, -4, 4);
0153   delphiphomet_reco_ =
0154       iBooker.book1D("delphiphomet_reco", "Reco: #Delta#phi(photon, MET);#Delta#phi(#gamma,MET)", 100, 0, 4);
0155   delphiphomet_ = iBooker.book1D("delphiphomet", "#Delta#phi(photon, MET);#Delta#phi(#gamma,MET)", 100, 0, 4);
0156   npfjets_reco_ = iBooker.book1D("npfjets_reco", "Reco: Number of PF Jets", 100, 0, 20);
0157   npfjets_ = iBooker.book1D("npfjets", "Number of PF Jets", 100, 0, 20);
0158   delphijetmet_reco_ =
0159       iBooker.book1D("delphijetmet_reco", "Reco: #Delta#phi(PFJet, MET);#Delta#phi(Jet,MET)", 100, 0, 4);
0160   delphijetmet_ = iBooker.book1D("delphijetmet", "#Delta#phi(PFJet, MET);#Delta#phi(Jet,MET)", 100, 0, 4);
0161   invmassjj_reco_ = iBooker.book1D("invmassjj_reco", "Reco: Inv mass two leading jets;M_{jj}[GeV]", 100, 0, 2000);
0162   invmassjj_ = iBooker.book1D("invmassjj", "Inv mass two leading jets;M_{jj}[GeV]", 100, 0, 2000);
0163   deletajj_reco_ = iBooker.book1D("deletajj_reco", "Reco: #Delta#eta(jj);|#Delta#eta_{jj}|", 100, 0, 6);
0164   deletajj_ = iBooker.book1D("deletajj", "#Delta#eta(jj);|#Delta#eta_{jj}|", 100, 0, 6);
0165   triggers_reco_ =
0166       iBooker.book1D("triggers_reco", "Reco: Triggers", hltPathsToCheck_.size(), 0, hltPathsToCheck_.size());
0167   triggers_ = iBooker.book1D("triggers", "Triggers", hltPathsToCheck_.size(), 0, hltPathsToCheck_.size());
0168   trigvsnvtx_reco_ = iBooker.book2D("trigvsnvtx_reco",
0169                                     "Reco: Trigger vs. # vertices;N_{vertices};Trigger",
0170                                     100,
0171                                     0,
0172                                     100,
0173                                     hltPathsToCheck_.size(),
0174                                     0,
0175                                     hltPathsToCheck_.size());
0176   trigvsnvtx_ = iBooker.book2D("trigvsnvtx",
0177                                "Trigger vs. # vertices;N_{vertices};Trigger",
0178                                100,
0179                                0,
0180                                100,
0181                                hltPathsToCheck_.size(),
0182                                0,
0183                                hltPathsToCheck_.size());
0184 }
0185 
0186 void HigPhotonJetHLTOfflineSource::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0187   // Count total number of events in one run
0188   evtsrun_++;
0189 
0190   edm::Handle<edm::TriggerResults> triggerResults;
0191   iEvent.getByToken(triggerResultsToken_, triggerResults);
0192   if (!triggerResults.isValid()) {
0193     edm::LogError("HigPhotonJetHLT") << "Missing triggerResults collection" << std::endl;
0194     return;
0195   }
0196 
0197   // Check whether contains monitored trigger and accepted
0198   const edm::TriggerNames& triggerNames = iEvent.triggerNames(*triggerResults);
0199   bool triggered = isMonitoredTriggerAccepted(triggerNames, triggerResults);
0200 
0201   // if (!triggered) return;
0202 
0203   // Test scale
0204   // if (evtsrun_ > 10) return;
0205 
0206   // N Vertices
0207   edm::Handle<reco::VertexCollection> vertices;
0208   iEvent.getByToken(pvToken_, vertices);
0209   if (!vertices.isValid())
0210     return;
0211   if (verbose_)
0212     std::cout << "xshi:: N vertices : " << vertices->size() << std::endl;
0213 
0214   // Set trigger name labels
0215   for (size_t i = 0; i < hltPathsToCheck_.size(); i++) {
0216     triggers_->setBinLabel(i + 1, hltPathsToCheck_[i]);
0217   }
0218 
0219   // Fill trigger info
0220   for (unsigned int itrig = 0; itrig < triggerResults->size(); itrig++) {
0221     const std::string& triggername = triggerNames.triggerName(itrig);
0222     for (size_t i = 0; i < hltPathsToCheck_.size(); i++) {
0223       if (triggername.find(hltPathsToCheck_[i]) != std::string::npos) {
0224         triggers_reco_->Fill(i);
0225         trigvsnvtx_reco_->Fill(vertices->size(), i);
0226         if (triggered)
0227           triggers_->Fill(i);
0228         if (triggered)
0229           trigvsnvtx_->Fill(vertices->size(), i);
0230       }
0231     }
0232   }
0233 
0234   nvertices_reco_->Fill(vertices->size());
0235   if (triggered)
0236     nvertices_->Fill(vertices->size());
0237 
0238   // PF MET
0239   edm::Handle<reco::PFMETCollection> pfmets;
0240   iEvent.getByToken(pfMetToken_, pfmets);
0241   if (!pfmets.isValid())
0242     return;
0243   const reco::PFMET pfmet = pfmets->front();
0244   pfmet_reco_->Fill(pfmet.et());
0245   if (triggered)
0246     pfmet_->Fill(pfmet.et());
0247   if (verbose_)
0248     std::cout << "xshi:: number of pfmets: " << pfmets->size() << std::endl;
0249 
0250   pfmetphi_reco_->Fill(pfmet.phi());
0251   if (triggered)
0252     pfmetphi_->Fill(pfmet.phi());
0253 
0254   // Photons
0255   edm::Handle<reco::PhotonCollection> photons;
0256   iEvent.getByToken(photonsToken_, photons);
0257   if (!photons.isValid())
0258     return;
0259   int nphotons = 0;
0260   for (auto const& phoIter : *photons) {
0261     if (phoIter.pt() < photonMinPt_)
0262       continue;
0263     nphotons++;
0264     photonpt_reco_->Fill(phoIter.pt());
0265     photonrapidity_reco_->Fill(phoIter.rapidity());
0266     if (triggered)
0267       photonpt_->Fill(phoIter.pt());
0268     if (triggered)
0269       photonrapidity_->Fill(phoIter.rapidity());
0270     double tmp_delphiphomet = fabs(deltaPhi(phoIter.phi(), pfmet.phi()));
0271     delphiphomet_reco_->Fill(tmp_delphiphomet);
0272     if (triggered)
0273       delphiphomet_->Fill(tmp_delphiphomet);
0274   }
0275   nphotons_reco_->Fill(nphotons);
0276   if (triggered)
0277     nphotons_->Fill(nphotons);
0278 
0279   // PF Jet
0280   edm::Handle<reco::PFJetCollection> pfjets;
0281   iEvent.getByToken(pfJetsToken_, pfjets);
0282   if (!pfjets.isValid())
0283     return;
0284   if (verbose_)
0285     std::cout << "xshi:: N pfjets : " << pfjets->size() << std::endl;
0286 
0287   double min_delphijetmet = 6.0;
0288   TLorentzVector p4jet1, p4jet2, p4jj;
0289   // Two leading jets eta
0290   double etajet1(0), etajet2(0);
0291   int njet = 0;
0292   for (auto const& jetIter : *pfjets) {
0293     if (jetIter.pt() < pfjetMinPt_)
0294       continue;
0295     njet++;
0296 
0297     double tmp_delphijetmet = fabs(deltaPhi(jetIter.phi(), pfmet.phi()));
0298     if (tmp_delphijetmet < min_delphijetmet)
0299       min_delphijetmet = tmp_delphijetmet;
0300 
0301     if (njet == 1) {
0302       p4jet1.SetXYZM(jetIter.px(), jetIter.py(), jetIter.pz(), jetIter.mass());
0303       etajet1 = jetIter.eta();
0304     }
0305     if (njet == 2) {
0306       p4jet2.SetXYZM(jetIter.px(), jetIter.py(), jetIter.pz(), jetIter.mass());
0307       etajet2 = jetIter.eta();
0308     }
0309   }
0310   npfjets_reco_->Fill(njet);
0311   if (triggered)
0312     npfjets_->Fill(njet);
0313 
0314   delphijetmet_reco_->Fill(min_delphijetmet);
0315   if (triggered)
0316     delphijetmet_->Fill(min_delphijetmet);
0317   p4jj = p4jet1 + p4jet2;
0318   double deletajj = etajet1 - etajet2;
0319   if (verbose_)
0320     std::cout << "xshi:: invmass jj " << p4jj.M() << std::endl;
0321 
0322   invmassjj_reco_->Fill(p4jj.M());
0323   deletajj_reco_->Fill(deletajj);
0324   if (triggered)
0325     invmassjj_->Fill(p4jj.M());
0326   if (triggered)
0327     deletajj_->Fill(deletajj);
0328 }
0329 
0330 void HigPhotonJetHLTOfflineSource::dqmEndRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0331   // Normalize to the total number of events in the run
0332   if (!perLSsaving_) {
0333     TH2F* h = trigvsnvtx_->getTH2F();
0334     double integral = h->Integral();
0335     double norm = (integral > 0.) ? evtsrun_ * hltPathsToCheck_.size() / integral : 1.;
0336     h->Scale(norm);
0337     if (verbose_) {
0338       std::cout << "xshi:: endRun total number of events: " << evtsrun_ << ", integral = " << h->Integral()
0339                 << ", norm = " << norm << std::endl;
0340     }
0341   }
0342 }
0343 
0344 bool HigPhotonJetHLTOfflineSource::isMonitoredTriggerAccepted(const edm::TriggerNames& triggerNames,
0345                                                               const edm::Handle<edm::TriggerResults>& triggerResults) {
0346   for (unsigned int itrig = 0; itrig < triggerResults->size(); itrig++) {
0347     // Only consider the triggered case.
0348     if (triggerAccept_ && ((*triggerResults)[itrig].accept() != 1))
0349       continue;
0350     const std::string& triggername = triggerNames.triggerName(itrig);
0351     for (auto const& i : hltPathsToCheck_) {
0352       if (triggername.find(i) != std::string::npos) {
0353         return true;
0354       }
0355     }
0356   }
0357 
0358   return false;
0359 }
0360 
0361 //define this as a plug-in
0362 DEFINE_FWK_MODULE(HigPhotonJetHLTOfflineSource);