Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-09-18 22:46:52

0001 /*
0002   HLTInclusiveVBFSource
0003   Phat Srimanobhas
0004   To monitor VBF DataParking
0005 */
0006 
0007 #include "DQMOffline/Trigger/interface/HLTInclusiveVBFSource.h"
0008 
0009 #include "FWCore/Common/interface/TriggerNames.h"
0010 #include "FWCore/Framework/interface/Run.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 
0015 #include "DQMServices/Core/interface/DQMStore.h"
0016 
0017 #include "DataFormats/Common/interface/Handle.h"
0018 #include "DataFormats/Common/interface/TriggerResults.h"
0019 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0020 #include "DataFormats/HLTReco/interface/TriggerObject.h"
0021 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0022 #include "DataFormats/Math/interface/deltaR.h"
0023 #include "DataFormats/Math/interface/deltaPhi.h"
0024 #include "DataFormats/VertexReco/interface/Vertex.h"
0025 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0026 
0027 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0028 
0029 #include "JetMETCorrections/Objects/interface/JetCorrector.h"
0030 
0031 #include <cmath>
0032 #include "TH1F.h"
0033 #include "TProfile.h"
0034 #include "TH2F.h"
0035 #include "TPRegexp.h"
0036 #include "TMath.h"
0037 
0038 using namespace edm;
0039 using namespace reco;
0040 using namespace std;
0041 
0042 HLTInclusiveVBFSource::HLTInclusiveVBFSource(const edm::ParameterSet& iConfig) {
0043   LogDebug("HLTInclusiveVBFSource") << "constructor....";
0044   nCount_ = 0;
0045 
0046   dirname_ = iConfig.getUntrackedParameter("dirname", std::string("HLT/InclusiveVBF"));
0047   processname_ = iConfig.getParameter<std::string>("processname");
0048   triggerSummaryLabel_ = iConfig.getParameter<edm::InputTag>("triggerSummaryLabel");
0049   triggerResultsLabel_ = iConfig.getParameter<edm::InputTag>("triggerResultsLabel");
0050   triggerSummaryToken = consumes<trigger::TriggerEvent>(triggerSummaryLabel_);
0051   triggerResultsToken = consumes<edm::TriggerResults>(triggerResultsLabel_);
0052   triggerSummaryFUToken = consumes<trigger::TriggerEvent>(
0053       edm::InputTag(triggerSummaryLabel_.label(), triggerSummaryLabel_.instance(), std::string("FU")));
0054   triggerResultsFUToken = consumes<edm::TriggerResults>(
0055       edm::InputTag(triggerResultsLabel_.label(), triggerResultsLabel_.instance(), std::string("FU")));
0056 
0057   //path_                = iConfig.getUntrackedParameter<std::vector<std::string> >("paths");
0058   //l1path_              = iConfig.getUntrackedParameter<std::vector<std::string> >("l1paths");
0059   debug_ = iConfig.getUntrackedParameter<bool>("debug", false);
0060 
0061   caloJetsToken = consumes<reco::CaloJetCollection>(iConfig.getParameter<edm::InputTag>("CaloJetCollectionLabel"));
0062   caloMetToken = consumes<reco::CaloMETCollection>(iConfig.getParameter<edm::InputTag>("CaloMETCollectionLabel"));
0063   pfJetsToken = consumes<edm::View<reco::PFJet> >(iConfig.getParameter<edm::InputTag>("PFJetCollectionLabel"));
0064   pfMetToken = consumes<edm::View<reco::PFMET> >(iConfig.getParameter<edm::InputTag>("PFMETCollectionLabel"));
0065   //jetID                = new reco::helper::JetIDHelper(iConfig.getParameter<ParameterSet>("JetIDParams"));
0066 
0067   minPtHigh_ = iConfig.getUntrackedParameter<double>("minPtHigh", 40.);
0068   minPtLow_ = iConfig.getUntrackedParameter<double>("minPtLow", 40.);
0069   minDeltaEta_ = iConfig.getUntrackedParameter<double>("minDeltaEta", 3.5);
0070   deltaRMatch_ = iConfig.getUntrackedParameter<double>("deltaRMatch", 0.1);
0071   minInvMass_ = iConfig.getUntrackedParameter<double>("minInvMass", 1000.0);
0072   etaOpposite_ = iConfig.getUntrackedParameter<bool>("etaOpposite", true);
0073 
0074   check_mjj650_Pt35_DEta3p5 = false;
0075   check_mjj700_Pt35_DEta3p5 = false;
0076   check_mjj750_Pt35_DEta3p5 = false;
0077   check_mjj800_Pt35_DEta3p5 = false;
0078   check_mjj650_Pt40_DEta3p5 = false;
0079   check_mjj700_Pt40_DEta3p5 = false;
0080   check_mjj750_Pt40_DEta3p5 = false;
0081   check_mjj800_Pt40_DEta3p5 = false;
0082 }
0083 
0084 HLTInclusiveVBFSource::~HLTInclusiveVBFSource() {
0085   //
0086   // do anything here that needs to be done at desctruction time
0087   // (e.g. close files, deallocate resources etc.)
0088 }
0089 
0090 void HLTInclusiveVBFSource::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0091   using namespace std;
0092   using namespace edm;
0093   using namespace trigger;
0094   using namespace reco;
0095 
0096   if (debug_)
0097     cout << "DEBUG-0: Start to analyze" << endl;
0098 
0099   //****************************************************
0100   // Get trigger information.
0101   //****************************************************
0102   //
0103   //---------- triggerResults ----------
0104   iEvent.getByToken(triggerResultsToken, triggerResults_);
0105   if (!triggerResults_.isValid()) {
0106     iEvent.getByToken(triggerResultsFUToken, triggerResults_);
0107     if (!triggerResults_.isValid()) {
0108       edm::LogInfo("HLTInclusiveVBFSource") << "TriggerResults not found, "
0109                                                "skipping event";
0110       return;
0111     }
0112   }
0113 
0114   // Check how many HLT triggers are in triggerResults
0115   triggerNames_ = iEvent.triggerNames(*triggerResults_);
0116 
0117   //---------- triggerSummary ----------
0118   iEvent.getByToken(triggerSummaryToken, triggerObj_);
0119   if (!triggerObj_.isValid()) {
0120     iEvent.getByToken(triggerSummaryFUToken, triggerObj_);
0121     if (!triggerObj_.isValid()) {
0122       edm::LogInfo("HLTInclusiveVBFSource") << "TriggerEvent not found, "
0123                                                "skipping event";
0124       return;
0125     }
0126   }
0127 
0128   if (debug_)
0129     cout << "DEBUG-1: Trigger information" << endl;
0130 
0131   //****************************************************
0132   // Get AOD information
0133   //****************************************************
0134   //
0135   edm::Handle<edm::View<reco::PFMET> > metSrc;
0136   bool ValidPFMET_ = iEvent.getByToken(pfMetToken, metSrc);
0137   if (!ValidPFMET_)
0138     return;
0139 
0140   edm::Handle<edm::View<reco::PFJet> > jetSrc;
0141   bool ValidPFJet_ = iEvent.getByToken(pfJetsToken, jetSrc);
0142   if (!ValidPFJet_)
0143     return;
0144 
0145   if (!metSrc.isValid())
0146     return;
0147   if (!jetSrc.isValid())
0148     return;
0149   const edm::View<reco::PFMET>& mets = *metSrc;
0150   const edm::View<reco::PFJet>& jets = *jetSrc;
0151   if (jets.empty())
0152     return;
0153   if (mets.empty())
0154     return;
0155 
0156   if (debug_)
0157     cout << "DEBUG-2: AOD Information" << endl;
0158 
0159   //****************************************************
0160   // Variable setting
0161   //****************************************************
0162   //
0163   pathname = "dummy";
0164   filtername = "dummy";
0165 
0166   //
0167   reco_ejet1 = 0.;
0168   //reco_etjet1               = 0.;
0169   reco_pxjet1 = 0.;
0170   reco_pyjet1 = 0.;
0171   reco_pzjet1 = 0.;
0172   reco_ptjet1 = 0.;
0173   reco_etajet1 = 0.;
0174   reco_phijet1 = 0.;
0175 
0176   //
0177   reco_ejet2 = 0.;
0178   //reco_etjet2               = 0.;
0179   reco_pxjet2 = 0.;
0180   reco_pyjet2 = 0.;
0181   reco_pzjet2 = 0.;
0182   reco_ptjet2 = 0.;
0183   reco_etajet2 = 0.;
0184   reco_phijet2 = 0.;
0185 
0186   //
0187   hlt_ejet1 = 0.;
0188   //hlt_etjet1                = 0.;
0189   hlt_pxjet1 = 0.;
0190   hlt_pyjet1 = 0.;
0191   hlt_pzjet1 = 0.;
0192   hlt_ptjet1 = 0.;
0193   hlt_etajet1 = 0.;
0194   hlt_phijet1 = 0.;
0195 
0196   //
0197   hlt_ejet2 = 0.;
0198   //hlt_etjet2                = 0.;
0199   hlt_pxjet2 = 0.;
0200   hlt_pyjet2 = 0.;
0201   hlt_pzjet2 = 0.;
0202   hlt_ptjet2 = 0.;
0203   hlt_etajet2 = 0.;
0204   hlt_phijet2 = 0.;
0205 
0206   //
0207   checkOffline = false;
0208   checkHLT = false;
0209   checkHLTIndex = false;
0210 
0211   //
0212   dR_HLT_RECO_11 = 0.;
0213   dR_HLT_RECO_22 = 0.;
0214   dR_HLT_RECO_12 = 0.;
0215   dR_HLT_RECO_21 = 0.;
0216 
0217   //
0218   checkdR_sameOrder = false;
0219   checkdR_crossOrder = false;
0220 
0221   //
0222   reco_deltaetajet = 0.;
0223   reco_deltaphijet = 0.;
0224   reco_invmassjet = 0.;
0225   hlt_deltaetajet = 0.;
0226   hlt_deltaphijet = 0.;
0227   hlt_invmassjet = 0.;
0228 
0229   //****************************************************
0230   // Offline analysis
0231   //****************************************************
0232   //
0233   checkOffline = false;
0234   for (unsigned int ijet1 = 0; ijet1 < jets.size(); ijet1++) {
0235     if (jets[ijet1].neutralHadronEnergyFraction() > 0.99)
0236       continue;
0237     if (jets[ijet1].neutralEmEnergyFraction() > 0.99)
0238       continue;
0239     for (unsigned int ijet2 = ijet1 + 1; ijet2 < jets.size(); ijet2++) {
0240       if (jets[ijet2].neutralHadronEnergyFraction() > 0.99)
0241         continue;
0242       if (jets[ijet2].neutralEmEnergyFraction() > 0.99)
0243         continue;
0244       //
0245       reco_ejet1 = jets[ijet1].energy();
0246       //reco_etjet1  = jets[ijet1].et();
0247       reco_pxjet1 = jets[ijet1].momentum().X();
0248       reco_pyjet1 = jets[ijet1].momentum().Y();
0249       reco_pzjet1 = jets[ijet1].momentum().Z();
0250       reco_ptjet1 = jets[ijet1].pt();
0251       reco_etajet1 = jets[ijet1].eta();
0252       reco_phijet1 = jets[ijet1].phi();
0253       //
0254       reco_ejet2 = jets[ijet2].energy();
0255       //reco_etjet2  = jets[ijet2].et();
0256       reco_pxjet2 = jets[ijet2].momentum().X();
0257       reco_pyjet2 = jets[ijet2].momentum().Y();
0258       reco_pzjet2 = jets[ijet2].momentum().Z();
0259       reco_ptjet2 = jets[ijet2].pt();
0260       reco_etajet2 = jets[ijet2].eta();
0261       reco_phijet2 = jets[ijet2].phi();
0262       //
0263       reco_deltaetajet = reco_etajet1 - reco_etajet2;
0264       reco_deltaphijet = reco::deltaPhi(reco_phijet1, reco_phijet2);
0265       reco_invmassjet = sqrt((reco_ejet1 + reco_ejet2) * (reco_ejet1 + reco_ejet2) -
0266                              (reco_pxjet1 + reco_pxjet2) * (reco_pxjet1 + reco_pxjet2) -
0267                              (reco_pyjet1 + reco_pyjet2) * (reco_pyjet1 + reco_pyjet2) -
0268                              (reco_pzjet1 + reco_pzjet2) * (reco_pzjet1 + reco_pzjet2));
0269 
0270       //
0271       if (reco_ptjet1 < minPtHigh_)
0272         continue;
0273       if (reco_ptjet2 < minPtLow_)
0274         continue;
0275       if (etaOpposite_ == true && reco_etajet1 * reco_etajet2 > 0)
0276         continue;
0277       if (std::abs(reco_deltaetajet) < minDeltaEta_)
0278         continue;
0279       if (std::abs(reco_invmassjet) < minInvMass_)
0280         continue;
0281 
0282       //
0283       if (debug_)
0284         cout << "DEBUG-3" << endl;
0285       checkOffline = true;
0286       break;
0287     }
0288     if (checkOffline == true)
0289       break;
0290   }
0291   if (checkOffline == false)
0292     return;
0293 
0294   //****************************************************
0295   // Trigger efficiency: Loop for all VBF paths
0296   //****************************************************
0297   //const unsigned int numberOfPaths(hltConfig_.size());
0298   const trigger::TriggerObjectCollection& toc(triggerObj_->getObjects());
0299   for (auto& v : hltPathsAll_) {
0300     checkHLT = false;
0301     checkHLTIndex = false;
0302 
0303     //
0304     v.getMEhisto_RECO_deltaEta_DiJet()->Fill(reco_deltaetajet);
0305     v.getMEhisto_RECO_deltaPhi_DiJet()->Fill(reco_deltaphijet);
0306     v.getMEhisto_RECO_invMass_DiJet()->Fill(reco_invmassjet);
0307 
0308     //
0309     if (debug_)
0310       cout << "DEBUG-4-0: Path loops" << endl;
0311 
0312     //
0313     if (isHLTPathAccepted(v.getPath()) == false)
0314       continue;
0315     checkHLT = true;
0316 
0317     //
0318     if (debug_)
0319       cout << "DEBUG-4-1: Path is accepted. Now we are looking for " << v.getLabel() << " module." << endl;
0320 
0321     //
0322     edm::InputTag hltTag(v.getLabel(), "", processname_);
0323     const int hltIndex = triggerObj_->filterIndex(hltTag);
0324     if (hltIndex >= triggerObj_->sizeFilters())
0325       continue;
0326     checkHLT = true;
0327     if (debug_)
0328       cout << "DEBUG-4-2: HLT module " << v.getLabel() << " exists" << endl;
0329     const trigger::Keys& khlt = triggerObj_->filterKeys(hltIndex);
0330     auto kj = khlt.begin();
0331     for (; kj != khlt.end(); kj += 2) {
0332       if (debug_)
0333         cout << "DEBUG-5" << endl;
0334       checkdR_sameOrder = false;
0335       checkdR_crossOrder = false;  //
0336       hlt_ejet1 = toc[*kj].energy();
0337       //hlt_etjet1  = toc[*kj].et();
0338       hlt_pxjet1 = toc[*kj].px();
0339       hlt_pyjet1 = toc[*kj].py();
0340       hlt_pzjet1 = toc[*kj].pz();
0341       hlt_ptjet1 = toc[*kj].pt();
0342       hlt_etajet1 = toc[*kj].eta();
0343       hlt_phijet1 = toc[*kj].phi();
0344       //
0345       hlt_ejet2 = toc[*(kj + 1)].energy();
0346       //hlt_etjet2  = toc[*(kj+1)].et();
0347       hlt_pxjet2 = toc[*(kj + 1)].px();
0348       hlt_pyjet2 = toc[*(kj + 1)].py();
0349       hlt_pzjet2 = toc[*(kj + 1)].pz();
0350       hlt_ptjet2 = toc[*(kj + 1)].pt();
0351       hlt_etajet2 = toc[*(kj + 1)].eta();
0352       hlt_phijet2 = toc[*(kj + 1)].phi();
0353       //
0354       dR_HLT_RECO_11 = reco::deltaR(hlt_etajet1, hlt_phijet1, reco_etajet1, reco_phijet1);
0355       dR_HLT_RECO_22 = reco::deltaR(hlt_etajet2, hlt_phijet2, reco_etajet2, reco_phijet2);
0356       dR_HLT_RECO_12 = reco::deltaR(hlt_etajet1, hlt_phijet1, reco_etajet2, reco_phijet2);
0357       dR_HLT_RECO_21 = reco::deltaR(hlt_etajet2, hlt_phijet2, reco_etajet1, reco_phijet1);
0358       if (dR_HLT_RECO_11 < deltaRMatch_ && dR_HLT_RECO_22 < deltaRMatch_)
0359         checkdR_sameOrder = true;
0360       if (dR_HLT_RECO_12 < deltaRMatch_ && dR_HLT_RECO_21 < deltaRMatch_)
0361         checkdR_crossOrder = true;
0362       if (checkdR_sameOrder == false && checkdR_crossOrder == false)
0363         continue;
0364       checkHLTIndex = true;
0365       //
0366       if (debug_)
0367         cout << "DEBUG-6: Match" << endl;
0368       hlt_deltaetajet = hlt_etajet1 - hlt_etajet2;
0369       hlt_deltaphijet = reco::deltaPhi(hlt_phijet1, hlt_phijet2);
0370       if (checkdR_crossOrder) {
0371         hlt_deltaetajet = (-1) * hlt_deltaetajet;
0372         hlt_deltaphijet = reco::deltaPhi(hlt_phijet2, hlt_phijet1);
0373       }
0374       hlt_invmassjet = sqrt((hlt_ejet1 + hlt_ejet2) * (hlt_ejet1 + hlt_ejet2) -
0375                             (hlt_pxjet1 + hlt_pxjet2) * (hlt_pxjet1 + hlt_pxjet2) -
0376                             (hlt_pyjet1 + hlt_pyjet2) * (hlt_pyjet1 + hlt_pyjet2) -
0377                             (hlt_pzjet1 + hlt_pzjet2) * (hlt_pzjet1 + hlt_pzjet2));
0378       v.getMEhisto_HLT_deltaEta_DiJet()->Fill(hlt_deltaetajet);
0379       v.getMEhisto_HLT_deltaPhi_DiJet()->Fill(hlt_deltaphijet);
0380       v.getMEhisto_HLT_invMass_DiJet()->Fill(hlt_invmassjet);
0381       //
0382       v.getMEhisto_RECO_deltaEta_DiJet_Match()->Fill(reco_deltaetajet);
0383       v.getMEhisto_RECO_deltaPhi_DiJet_Match()->Fill(reco_deltaphijet);
0384       v.getMEhisto_RECO_invMass_DiJet_Match()->Fill(reco_invmassjet);
0385       //
0386       v.getMEhisto_RECOHLT_deltaEta()->Fill(reco_deltaetajet, hlt_deltaetajet);
0387       v.getMEhisto_RECOHLT_deltaPhi()->Fill(reco_deltaphijet, hlt_deltaphijet);
0388       v.getMEhisto_RECOHLT_invMass()->Fill(reco_invmassjet, hlt_invmassjet);
0389       //
0390       if (checkHLTIndex == true)
0391         break;
0392     }
0393 
0394     //****************************************************
0395     // Match information
0396     //****************************************************
0397     if (checkHLT == true && checkHLTIndex == true) {
0398       if (debug_)
0399         cout << "DEBUG-7: Match" << endl;
0400       v.getMEhisto_NumberOfMatches()->Fill(1);
0401     } else {
0402       if (debug_)
0403         cout << "DEBUG-8: Not match" << endl;
0404       v.getMEhisto_NumberOfMatches()->Fill(0);
0405     }
0406   }
0407 
0408   //****************************************************
0409   //
0410   //****************************************************
0411   for (auto& v : hltPathsAll_) {
0412     if (isHLTPathAccepted(v.getPath()) == false)
0413       continue;
0414     if (debug_)
0415       cout << "DEBUG-9: Loop for rate approximation: " << v.getPath() << endl;
0416     check_mjj650_Pt35_DEta3p5 = false;
0417     check_mjj700_Pt35_DEta3p5 = false;
0418     check_mjj750_Pt35_DEta3p5 = false;
0419     check_mjj800_Pt35_DEta3p5 = false;
0420     check_mjj650_Pt40_DEta3p5 = false;
0421     check_mjj700_Pt40_DEta3p5 = false;
0422     check_mjj750_Pt40_DEta3p5 = false;
0423     check_mjj800_Pt40_DEta3p5 = false;
0424     edm::InputTag hltTag(v.getLabel(), "", processname_);
0425     const int hltIndex = triggerObj_->filterIndex(hltTag);
0426     if (hltIndex >= triggerObj_->sizeFilters())
0427       continue;
0428     const trigger::Keys& khlt = triggerObj_->filterKeys(hltIndex);
0429     auto kj = khlt.begin();
0430     for (; kj != khlt.end(); kj += 2) {
0431       checkdR_sameOrder = false;
0432       checkdR_crossOrder = false;
0433       //
0434       hlt_ejet1 = toc[*kj].energy();
0435       //hlt_etjet1  = toc[*kj].et();
0436       hlt_pxjet1 = toc[*kj].px();
0437       hlt_pyjet1 = toc[*kj].py();
0438       hlt_pzjet1 = toc[*kj].pz();
0439       hlt_ptjet1 = toc[*kj].pt();
0440       hlt_etajet1 = toc[*kj].eta();
0441       hlt_phijet1 = toc[*kj].phi();
0442       //
0443       hlt_ejet2 = toc[*(kj + 1)].energy();
0444       //hlt_etjet2  = toc[*(kj+1)].et();
0445       hlt_pxjet2 = toc[*(kj + 1)].px();
0446       hlt_pyjet2 = toc[*(kj + 1)].py();
0447       hlt_pzjet2 = toc[*(kj + 1)].pz();
0448       hlt_ptjet2 = toc[*(kj + 1)].pt();
0449       hlt_etajet2 = toc[*(kj + 1)].eta();
0450       hlt_phijet2 = toc[*(kj + 1)].phi();
0451       //
0452       hlt_deltaetajet = hlt_etajet1 - hlt_etajet2;
0453       hlt_deltaphijet = reco::deltaPhi(hlt_phijet1, hlt_phijet2);
0454       hlt_invmassjet = sqrt((hlt_ejet1 + hlt_ejet2) * (hlt_ejet1 + hlt_ejet2) -
0455                             (hlt_pxjet1 + hlt_pxjet2) * (hlt_pxjet1 + hlt_pxjet2) -
0456                             (hlt_pyjet1 + hlt_pyjet2) * (hlt_pyjet1 + hlt_pyjet2) -
0457                             (hlt_pzjet1 + hlt_pzjet2) * (hlt_pzjet1 + hlt_pzjet2));
0458       //
0459       if (check_mjj650_Pt35_DEta3p5 == false && hlt_ptjet1 > 35. && hlt_ptjet2 >= 35. && hlt_invmassjet > 650 &&
0460           std::abs(hlt_deltaetajet) > 3.5) {
0461         check_mjj650_Pt35_DEta3p5 = true;
0462       }
0463       if (check_mjj700_Pt35_DEta3p5 == false && hlt_ptjet1 > 35. && hlt_ptjet2 >= 35. && hlt_invmassjet > 700 &&
0464           std::abs(hlt_deltaetajet) > 3.5) {
0465         check_mjj700_Pt35_DEta3p5 = true;
0466       }
0467       if (check_mjj750_Pt35_DEta3p5 == false && hlt_ptjet1 > 35. && hlt_ptjet2 >= 35. && hlt_invmassjet > 750 &&
0468           std::abs(hlt_deltaetajet) > 3.5) {
0469         check_mjj750_Pt35_DEta3p5 = true;
0470       }
0471       if (check_mjj800_Pt35_DEta3p5 == false && hlt_ptjet1 > 35. && hlt_ptjet2 >= 35. && hlt_invmassjet > 800 &&
0472           std::abs(hlt_deltaetajet) > 3.5) {
0473         check_mjj800_Pt35_DEta3p5 = true;
0474       }
0475       if (check_mjj650_Pt40_DEta3p5 == false && hlt_ptjet1 > 40. && hlt_ptjet2 >= 40. && hlt_invmassjet > 650 &&
0476           std::abs(hlt_deltaetajet) > 3.5) {
0477         check_mjj650_Pt40_DEta3p5 = true;
0478       }
0479       if (check_mjj700_Pt40_DEta3p5 == false && hlt_ptjet1 > 40. && hlt_ptjet2 >= 40. && hlt_invmassjet > 700 &&
0480           std::abs(hlt_deltaetajet) > 3.5) {
0481         check_mjj700_Pt40_DEta3p5 = true;
0482       }
0483       if (check_mjj750_Pt40_DEta3p5 == false && hlt_ptjet1 > 40. && hlt_ptjet2 >= 40. && hlt_invmassjet > 750 &&
0484           std::abs(hlt_deltaetajet) > 3.5) {
0485         check_mjj750_Pt40_DEta3p5 = true;
0486       }
0487       if (check_mjj800_Pt40_DEta3p5 == false && hlt_ptjet1 > 40. && hlt_ptjet2 >= 40. && hlt_invmassjet > 800 &&
0488           std::abs(hlt_deltaetajet) > 3.5) {
0489         check_mjj800_Pt40_DEta3p5 = true;
0490       }
0491     }
0492     if (check_mjj650_Pt35_DEta3p5 == true)
0493       v.getMEhisto_NumberOfEvents()->Fill(0);
0494     if (check_mjj700_Pt35_DEta3p5 == true)
0495       v.getMEhisto_NumberOfEvents()->Fill(1);
0496     if (check_mjj750_Pt35_DEta3p5 == true)
0497       v.getMEhisto_NumberOfEvents()->Fill(2);
0498     if (check_mjj800_Pt35_DEta3p5 == true)
0499       v.getMEhisto_NumberOfEvents()->Fill(3);
0500     if (check_mjj650_Pt40_DEta3p5 == true)
0501       v.getMEhisto_NumberOfEvents()->Fill(4);
0502     if (check_mjj700_Pt40_DEta3p5 == true)
0503       v.getMEhisto_NumberOfEvents()->Fill(5);
0504     if (check_mjj750_Pt40_DEta3p5 == true)
0505       v.getMEhisto_NumberOfEvents()->Fill(6);
0506     if (check_mjj800_Pt40_DEta3p5 == true)
0507       v.getMEhisto_NumberOfEvents()->Fill(7);
0508   }
0509 }
0510 
0511 // BeginRun
0512 void HLTInclusiveVBFSource::bookHistograms(DQMStore::IBooker& iBooker, edm::Run const& run, edm::EventSetup const& c) {
0513   iBooker.setCurrentFolder(dirname_);
0514 
0515   //--- htlConfig_
0516   bool changed(true);
0517   if (!hltConfig_.init(run, c, processname_, changed)) {
0518     LogDebug("HLTInclusiveVBFSource") << "HLTConfigProvider failed to initialize.";
0519   }
0520 
0521   const unsigned int numberOfPaths(hltConfig_.size());
0522   for (unsigned int i = 0; i != numberOfPaths; ++i) {
0523     bool numFound = false;
0524     pathname = hltConfig_.triggerName(i);
0525     filtername = "dummy";
0526     unsigned int usedPrescale = 1;
0527     unsigned int objectType = 0;
0528     std::string triggerType = "";
0529 
0530     if (pathname.find("HLT_Di") == std::string::npos)
0531       continue;
0532     if (pathname.find("Jet") == std::string::npos)
0533       continue;
0534     if (pathname.find("MJJ") == std::string::npos)
0535       continue;
0536     if (pathname.find("VBF_v") == std::string::npos)
0537       continue;
0538 
0539     if (debug_) {
0540       cout << " - Startup:Path = " << pathname << endl;
0541       //cout<<" - Startup:PS = "<<hltConfig_.prescaleSize()<<endl;
0542     }
0543 
0544     triggerType = "DiJet_Trigger";
0545     objectType = trigger::TriggerJet;
0546 
0547     // Checking if the trigger exist in HLT table or not
0548     for (unsigned int i = 0; i != numberOfPaths; ++i) {
0549       std::string HLTname = hltConfig_.triggerName(i);
0550       if (HLTname == pathname)
0551         numFound = true;
0552     }
0553 
0554     if (numFound == false)
0555       continue;
0556     std::vector<std::string> numpathmodules = hltConfig_.moduleLabels(pathname);
0557     auto numpathmodule = numpathmodules.begin();
0558     for (; numpathmodule != numpathmodules.end(); ++numpathmodule) {
0559       edm::InputTag testTag(*numpathmodule, "", processname_);
0560       if (hltConfig_.moduleType(*numpathmodule) == "HLTCaloJetVBFFilter" ||
0561           hltConfig_.moduleType(*numpathmodule) == "HLTPFJetVBFFilter") {
0562         filtername = *numpathmodule;
0563         if (debug_)
0564           cout << " - Startup:Module = " << hltConfig_.moduleType(*numpathmodule) << ", FilterName = " << filtername
0565                << endl;
0566       }
0567     }
0568     if (debug_)
0569       cout << " - Startup:Final filter = " << filtername << endl;
0570 
0571     if (objectType == 0 || numFound == false)
0572       continue;
0573     //if(debug_){
0574     //cout<<"Pathname = "<<pathname
0575     //    <<", Filtername = "<<filtername
0576     //    <<", ObjectType = "<<objectType<<endl;
0577     //}
0578     hltPathsAll_.push_back(PathInfo(usedPrescale, pathname, filtername, processname_, objectType, triggerType));
0579   }  //Loop over paths
0580 
0581   //if(debug_) cout<<"== end hltPathsEff_.push_back ======" << endl;
0582 
0583   std::string dirName = dirname_ + "/MonitorInclusiveVBFTrigger/";
0584   for (auto& v : hltPathsAll_) {
0585     if (debug_)
0586       cout << "Storing: " << v.getPath() << ", Prescale = " << v.getprescaleUsed() << endl;
0587     //if(v->getprescaleUsed()!=1) continue;
0588 
0589     std::string subdirName = dirName + v.getPath();
0590     std::string trigPath = "(" + v.getPath() + ")";
0591     iBooker.setCurrentFolder(subdirName);
0592 
0593     MonitorElement* RECO_deltaEta_DiJet;
0594     MonitorElement* RECO_deltaPhi_DiJet;
0595     MonitorElement* RECO_invMass_DiJet;
0596     MonitorElement* HLT_deltaEta_DiJet;
0597     MonitorElement* HLT_deltaPhi_DiJet;
0598     MonitorElement* HLT_invMass_DiJet;
0599     MonitorElement* RECO_deltaEta_DiJet_Match;
0600     MonitorElement* RECO_deltaPhi_DiJet_Match;
0601     MonitorElement* RECO_invMass_DiJet_Match;
0602     MonitorElement* RECOHLT_deltaEta;
0603     MonitorElement* RECOHLT_deltaPhi;
0604     MonitorElement* RECOHLT_invMass;
0605     MonitorElement* NumberOfMatches;
0606     MonitorElement* NumberOfEvents;
0607 
0608     //dummy                     = iBooker.bookFloat("dummy");
0609     RECO_deltaEta_DiJet = iBooker.bookFloat("RECO_deltaEta_DiJet");
0610     RECO_deltaPhi_DiJet = iBooker.bookFloat("RECO_deltaPhi_DiJet");
0611     RECO_invMass_DiJet = iBooker.bookFloat("RECO_invMass_DiJet");
0612     HLT_deltaEta_DiJet = iBooker.bookFloat("HLT_deltaEta_DiJet");
0613     HLT_deltaPhi_DiJet = iBooker.bookFloat("HLT_deltaPhi_DiJet ");
0614     HLT_invMass_DiJet = iBooker.bookFloat("HLT_invMass_DiJet");
0615     RECO_deltaEta_DiJet_Match = iBooker.bookFloat("RECO_deltaEta_DiJet_Match");
0616     RECO_deltaPhi_DiJet_Match = iBooker.bookFloat("RECO_deltaPhi_DiJet_Match");
0617     RECO_invMass_DiJet_Match = iBooker.bookFloat("RECO_invMass_DiJet_Match");
0618     RECOHLT_deltaEta = iBooker.bookFloat("RECOHLT_deltaEta");
0619     RECOHLT_deltaPhi = iBooker.bookFloat("RECOHLT_deltaPhi ");
0620     RECOHLT_invMass = iBooker.bookFloat("RECOHLT_invMass");
0621     NumberOfMatches = iBooker.bookFloat("NumberOfMatches");
0622     NumberOfEvents = iBooker.bookFloat("NumberOfEvents");
0623 
0624     std::string labelname("ME");
0625     std::string histoname(labelname + "");
0626     std::string title(labelname + "");
0627 
0628     //RECO_deltaEta_DiJet
0629     histoname = labelname + "_RECO_deltaEta_DiJet";
0630     title = labelname + "_RECO_deltaEta_DiJet " + trigPath;
0631     RECO_deltaEta_DiJet = iBooker.book1D(histoname.c_str(), title.c_str(), 50, -10., 10.);
0632     RECO_deltaEta_DiJet->getTH1F();
0633 
0634     //RECO_deltaPhi_DiJet
0635     histoname = labelname + "_RECO_deltaPhi_DiJet";
0636     title = labelname + "_RECO_deltaPhi_DiJet " + trigPath;
0637     RECO_deltaPhi_DiJet = iBooker.book1D(histoname.c_str(), title.c_str(), 35, -3.5, 3.5);
0638     RECO_deltaPhi_DiJet->getTH1F();
0639 
0640     //RECO_invMass_DiJet
0641     histoname = labelname + "_RECO_invMass_DiJet";
0642     title = labelname + "_RECO_invMass_DiJet " + trigPath;
0643     RECO_invMass_DiJet = iBooker.book1D(histoname.c_str(), title.c_str(), 100, 500., 2000.);
0644     RECO_invMass_DiJet->getTH1F();
0645 
0646     //HLT_deltaEta_DiJet
0647     histoname = labelname + "_HLT_deltaEta_DiJet";
0648     title = labelname + "_HLT_deltaEta_DiJet " + trigPath;
0649     HLT_deltaEta_DiJet = iBooker.book1D(histoname.c_str(), title.c_str(), 50, -10., 10.);
0650     HLT_deltaEta_DiJet->getTH1F();
0651 
0652     //HLT_deltaPhi_DiJet
0653     histoname = labelname + "_HLT_deltaPhi_DiJet";
0654     title = labelname + "_HLT_deltaPhi_DiJet " + trigPath;
0655     HLT_deltaPhi_DiJet = iBooker.book1D(histoname.c_str(), title.c_str(), 35, -3.5, 3.5);
0656     HLT_deltaPhi_DiJet->getTH1F();
0657 
0658     //HLT_invMass_DiJet
0659     histoname = labelname + "_HLT_invMass_DiJet";
0660     title = labelname + "_HLT_invMass_DiJet " + trigPath;
0661     HLT_invMass_DiJet = iBooker.book1D(histoname.c_str(), title.c_str(), 100, 500., 2000.);
0662     HLT_invMass_DiJet->getTH1F();
0663 
0664     //RECO_deltaEta_DiJet_Match
0665     histoname = labelname + "_RECO_deltaEta_DiJet_Match";
0666     title = labelname + "_RECO_deltaEta_DiJet_Match " + trigPath;
0667     RECO_deltaEta_DiJet_Match = iBooker.book1D(histoname.c_str(), title.c_str(), 50, -10., 10.);
0668     RECO_deltaEta_DiJet_Match->getTH1F();
0669 
0670     //RECO_deltaPhi_DiJet_Match
0671     histoname = labelname + "_RECO_deltaPhi_DiJet_Match";
0672     title = labelname + "_RECO_deltaPhi_DiJet_Match " + trigPath;
0673     RECO_deltaPhi_DiJet_Match = iBooker.book1D(histoname.c_str(), title.c_str(), 35, -3.5, 3.5);
0674     RECO_deltaPhi_DiJet_Match->getTH1F();
0675 
0676     //RECO_invMass_DiJet_Match
0677     histoname = labelname + "_RECO_invMass_DiJet_Match";
0678     title = labelname + "_RECO_invMass_DiJet_Match " + trigPath;
0679     RECO_invMass_DiJet_Match = iBooker.book1D(histoname.c_str(), title.c_str(), 100, 500., 2000.);
0680     RECO_invMass_DiJet_Match->getTH1F();
0681 
0682     //RECOHLT_deltaEta
0683     histoname = labelname + "_RECOHLT_deltaEta";
0684     title = labelname + "_RECOHLT_deltaEta " + trigPath;
0685     RECOHLT_deltaEta = iBooker.book2D(histoname.c_str(), title.c_str(), 50, -10., 10., 50, -10., 10.);
0686     RECOHLT_deltaEta->getTH2F();
0687 
0688     //RECOHLT_deltaPhi
0689     histoname = labelname + "_RECOHLT_deltaPhi";
0690     title = labelname + "_RECOHLT_deltaPhi " + trigPath;
0691     RECOHLT_deltaPhi = iBooker.book2D(histoname.c_str(), title.c_str(), 35, -3.5, 3.5, 35, -3.5, 3.5);
0692     RECOHLT_deltaPhi->getTH2F();
0693 
0694     //RECOHLT_invMass
0695     histoname = labelname + "_RECOHLT_invMass";
0696     title = labelname + "_RECOHLT_invMass " + trigPath;
0697     RECOHLT_invMass = iBooker.book2D(histoname.c_str(), title.c_str(), 100, 500., 2000., 100, 500., 2000.);
0698     RECOHLT_invMass->getTH2F();
0699 
0700     //NumberOfMatches
0701     histoname = labelname + "_NumberOfMatches ";
0702     title = labelname + "_NumberOfMatches  " + trigPath;
0703     NumberOfMatches = iBooker.book1D(histoname.c_str(), title.c_str(), 2, 0., 2.);
0704     NumberOfMatches->getTH1F();
0705 
0706     //NumberOfEvents
0707     histoname = labelname + "_NumberOfEvents";
0708     title = labelname + "_NumberOfEvents " + trigPath;
0709     NumberOfEvents = iBooker.book1D(histoname.c_str(), title.c_str(), 10, 0., 10.);
0710     NumberOfEvents->getTH1F();
0711 
0712     //}
0713     v.setHistos(RECO_deltaEta_DiJet,
0714                 RECO_deltaPhi_DiJet,
0715                 RECO_invMass_DiJet,
0716                 HLT_deltaEta_DiJet,
0717                 HLT_deltaPhi_DiJet,
0718                 HLT_invMass_DiJet,
0719                 RECO_deltaEta_DiJet_Match,
0720                 RECO_deltaPhi_DiJet_Match,
0721                 RECO_invMass_DiJet_Match,
0722                 RECOHLT_deltaEta,
0723                 RECOHLT_deltaPhi,
0724                 RECOHLT_invMass,
0725                 NumberOfMatches,
0726                 NumberOfEvents);
0727     //break;//We need only the first unprescale paths
0728   }
0729 }
0730 
0731 bool HLTInclusiveVBFSource::isBarrel(double eta) {
0732   bool output = false;
0733   if (fabs(eta) <= 1.3)
0734     output = true;
0735   return output;
0736 }
0737 
0738 bool HLTInclusiveVBFSource::isEndCap(double eta) {
0739   bool output = false;
0740   if (fabs(eta) <= 3.0 && fabs(eta) > 1.3)
0741     output = true;
0742   return output;
0743 }
0744 
0745 bool HLTInclusiveVBFSource::isForward(double eta) {
0746   bool output = false;
0747   if (fabs(eta) > 3.0)
0748     output = true;
0749   return output;
0750 }
0751 
0752 bool HLTInclusiveVBFSource::validPathHLT(std::string pathname) {
0753   // hltConfig_ has to be defined first before calling this method
0754   bool output = false;
0755   for (unsigned int j = 0; j != hltConfig_.size(); ++j) {
0756     if (hltConfig_.triggerName(j) == pathname)
0757       output = true;
0758   }
0759   return output;
0760 }
0761 
0762 bool HLTInclusiveVBFSource::isHLTPathAccepted(std::string pathName) {
0763   // triggerResults_, triggerNames_ has to be defined first before calling this method
0764   bool output = false;
0765   if (triggerResults_.isValid()) {
0766     unsigned index = triggerNames_.triggerIndex(pathName);
0767     if (index < triggerNames_.size() && triggerResults_->accept(index))
0768       output = true;
0769   }
0770   return output;
0771 }
0772 
0773 bool HLTInclusiveVBFSource::isTriggerObjectFound(std::string objectName) {
0774   // processname_, triggerObj_ has to be defined before calling this method
0775   bool output = false;
0776   edm::InputTag testTag(objectName, "", processname_);
0777   const int index = triggerObj_->filterIndex(testTag);
0778   if (index >= triggerObj_->sizeFilters()) {
0779     edm::LogInfo("HLTInclusiveVBFSource") << "no index " << index << " of that name ";
0780   } else {
0781     const trigger::Keys& k = triggerObj_->filterKeys(index);
0782     if (!k.empty())
0783       output = true;
0784   }
0785   return output;
0786 }
0787 
0788 #include "FWCore/Framework/interface/MakerMacros.h"
0789 DEFINE_FWK_MODULE(HLTInclusiveVBFSource);