Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:10:53

0001 #include "DQM/Physics/src/B2GDQM.h"
0002 
0003 #include <memory>
0004 
0005 // DQM
0006 #include "DQMServices/Core/interface/DQMStore.h"
0007 
0008 // Framework
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "FWCore/Framework/interface/LuminosityBlock.h"
0012 #include "FWCore/ServiceRegistry/interface/Service.h"
0013 #include "FWCore/ParameterSet/interface/FileInPath.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Framework/interface/Run.h"
0016 #include "DataFormats/Provenance/interface/EventID.h"
0017 
0018 // Candidate handling
0019 #include "DataFormats/Candidate/interface/Candidate.h"
0020 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0021 #include "DataFormats/Candidate/interface/OverlapChecker.h"
0022 #include "DataFormats/Candidate/interface/CompositeCandidate.h"
0023 #include "DataFormats/Candidate/interface/CompositeCandidateFwd.h"
0024 #include "DataFormats/Candidate/interface/CandMatchMap.h"
0025 
0026 // Vertex utilities
0027 #include "DataFormats/VertexReco/interface/Vertex.h"
0028 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0029 
0030 // Other
0031 #include "DataFormats/TrackReco/interface/Track.h"
0032 #include "DataFormats/DetId/interface/DetId.h"
0033 #include "DataFormats/Common/interface/RefToBase.h"
0034 
0035 // Math
0036 #include "DataFormats/Math/interface/Vector3D.h"
0037 #include "DataFormats/Math/interface/LorentzVector.h"
0038 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0039 #include "DataFormats/Common/interface/AssociationVector.h"
0040 #include "DataFormats/Math/interface/Point3D.h"
0041 #include "DataFormats/Math/interface/deltaR.h"
0042 #include "DataFormats/Math/interface/deltaPhi.h"
0043 
0044 // vertexing
0045 
0046 // Transient tracks
0047 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0048 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0049 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0050 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
0051 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0052 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0053 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0054 
0055 // JetCorrection
0056 #include "JetMETCorrections/Objects/interface/JetCorrector.h"
0057 
0058 // ROOT
0059 #include "TLorentzVector.h"
0060 
0061 // STDLIB
0062 #include <iostream>
0063 #include <iomanip>
0064 #include <cstdio>
0065 #include <string>
0066 #include <sstream>
0067 #include <cmath>
0068 
0069 using namespace edm;
0070 using namespace std;
0071 using namespace reco;
0072 using namespace trigger;
0073 
0074 //
0075 // -- Constructor
0076 //
0077 B2GDQM::B2GDQM(const edm::ParameterSet& ps) {
0078   edm::LogInfo("B2GDQM") << " Starting B2GDQM "
0079                          << "\n";
0080 
0081   typedef std::vector<edm::InputTag> vtag;
0082 
0083   // Get parameters from configuration file
0084   // Trigger
0085   theTriggerResultsCollection = ps.getParameter<InputTag>("triggerResultsCollection");
0086   triggerToken_ = consumes<edm::TriggerResults>(theTriggerResultsCollection);
0087 
0088   // Jets
0089   jetLabels_ = ps.getParameter<std::vector<edm::InputTag> >("jetLabels");
0090   for (std::vector<edm::InputTag>::const_iterator jetlabel = jetLabels_.begin(), jetlabelEnd = jetLabels_.end();
0091        jetlabel != jetlabelEnd;
0092        ++jetlabel) {
0093     jetTokens_.push_back(consumes<edm::View<reco::Jet> >(*jetlabel));
0094   }
0095   sdjetLabel_ = ps.getParameter<edm::InputTag>("sdjetLabel");
0096   sdjetToken_ = consumes<edm::View<reco::BasicJet> >(sdjetLabel_);
0097 
0098   muonToken_ = consumes<edm::View<reco::Muon> >(ps.getParameter<edm::InputTag>("muonSrc"));
0099   electronToken_ = consumes<edm::View<reco::GsfElectron> >(ps.getParameter<edm::InputTag>("electronSrc"));
0100 
0101   jetPtMins_ = ps.getParameter<std::vector<double> >("jetPtMins");
0102   allHadPtCut_ = ps.getParameter<double>("allHadPtCut");
0103   allHadRapidityCut_ = ps.getParameter<double>("allHadRapidityCut");
0104   allHadDeltaPhiCut_ = ps.getParameter<double>("allHadDeltaPhiCut");
0105 
0106   semiMu_HadJetPtCut_ = ps.getParameter<double>("semiMu_HadJetPtCut");
0107   semiMu_LepJetPtCut_ = ps.getParameter<double>("semiMu_LepJetPtCut");
0108   semiMu_dphiHadCut_ = ps.getParameter<double>("semiMu_dphiHadCut");
0109   semiMu_dRMin_ = ps.getParameter<double>("semiMu_dRMin");
0110   semiMu_ptRel_ = ps.getParameter<double>("semiMu_ptRel");
0111   muonSelect_ = std::make_shared<StringCutObjectSelector<reco::Muon> >(
0112 
0113       ps.getParameter<std::string>("muonSelect"));
0114 
0115   semiE_HadJetPtCut_ = ps.getParameter<double>("semiE_HadJetPtCut");
0116   semiE_LepJetPtCut_ = ps.getParameter<double>("semiE_LepJetPtCut");
0117   semiE_dphiHadCut_ = ps.getParameter<double>("semiE_dphiHadCut");
0118   semiE_dRMin_ = ps.getParameter<double>("semiE_dRMin");
0119   semiE_ptRel_ = ps.getParameter<double>("semiE_ptRel");
0120   elecSelect_ = std::make_shared<StringCutObjectSelector<reco::GsfElectron> >(
0121 
0122       ps.getParameter<std::string>("elecSelect"));
0123 
0124   PFJetCorService_ = ps.getParameter<std::string>("PFJetCorService");
0125 
0126   // MET
0127   PFMETLabel_ = ps.getParameter<InputTag>("pfMETCollection");
0128   PFMETToken_ = consumes<std::vector<reco::PFMET> >(PFMETLabel_);
0129 }
0130 
0131 //
0132 // -- Destructor
0133 //
0134 B2GDQM::~B2GDQM() {
0135   edm::LogInfo("B2GDQM") << " Deleting B2GDQM "
0136                          << "\n";
0137 }
0138 
0139 //
0140 //  -- Book histograms
0141 //
0142 void B2GDQM::bookHistograms(DQMStore::IBooker& bei, edm::Run const&, edm::EventSetup const&) {
0143   bei.setCurrentFolder("Physics/B2G");
0144 
0145   //--- Jets
0146 
0147   for (unsigned int icoll = 0; icoll < jetLabels_.size(); ++icoll) {
0148     std::stringstream ss;
0149     ss << "Physics/B2G/" << jetLabels_[icoll].label();
0150     bei.setCurrentFolder(ss.str());
0151     pfJet_pt.push_back(bei.book1D("pfJet_pt", "Pt of PFJet (GeV)", 50, 0.0, 1000));
0152     pfJet_y.push_back(bei.book1D("pfJet_y", "Rapidity of PFJet", 60, -6.0, 6.0));
0153     pfJet_phi.push_back(bei.book1D("pfJet_phi", "#phi of PFJet (radians)", 60, -3.14159, 3.14159));
0154     pfJet_m.push_back(bei.book1D("pfJet_m", "Mass of PFJet (GeV)", 50, 0.0, 500));
0155     pfJet_chef.push_back(bei.book1D("pfJet_pfchef", "PFJetID CHEF", 50, 0.0, 1.0));
0156     pfJet_nhef.push_back(bei.book1D("pfJet_pfnhef", "PFJetID NHEF", 50, 0.0, 1.0));
0157     pfJet_cemf.push_back(bei.book1D("pfJet_pfcemf", "PFJetID CEMF", 50, 0.0, 1.0));
0158     pfJet_nemf.push_back(bei.book1D("pfJet_pfnemf", "PFJetID NEMF", 50, 0.0, 1.0));
0159 
0160     boostedJet_subjetPt.push_back(bei.book1D("boostedJet_subjetPt", "Pt of subjets (GeV)", 50, 0.0, 500));
0161     boostedJet_subjetY.push_back(bei.book1D("boostedJet_subjetY", "Rapidity of subjets", 60, -6.0, 6.0));
0162     boostedJet_subjetPhi.push_back(
0163         bei.book1D("boostedJet_subjetPhi", "#phi of subjets (radians)", 60, -3.14159, 3.14159));
0164     boostedJet_subjetM.push_back(bei.book1D("boostedJet_subjetM", "Mass of subjets (GeV)", 50, 0.0, 250.));
0165     boostedJet_subjetN.push_back(bei.book1D("boostedJet_subjetN", "Number of subjets", 10, 0, 10));
0166     boostedJet_massDrop.push_back(bei.book1D("boostedJet_massDrop", "Mass drop for W-like jets", 50, 0.0, 1.0));
0167     boostedJet_wMass.push_back(bei.book1D("boostedJet_wMass", "W Mass for top-like jets", 50, 0.0, 250.0));
0168   }
0169 
0170   bei.setCurrentFolder("Physics/B2G/MET");
0171   pfMet_pt = bei.book1D("pfMet_pt", "Pf Missing p_{T}; GeV", 50, 0.0, 500);
0172   pfMet_phi = bei.book1D("pfMet_phi", "Pf Missing p_{T} #phi;#phi (radians)", 35, -3.5, 3.5);
0173 
0174   //--- Mu+Jets
0175   bei.setCurrentFolder("Physics/B2G/SemiMu");
0176   semiMu_muPt = bei.book1D("semiMu_muPt", "Pt of Muon in #mu+Jets Channel (GeV)", 50, 0.0, 1000);
0177   semiMu_muEta = bei.book1D("semiMu_muEta", "#eta of Muon in #mu+Jets Channel", 60, -6.0, 6.0);
0178   semiMu_muPhi = bei.book1D("semiMu_muPhi", "#phi of Muon in #mu+Jets Channel (radians)", 60, -3.14159, 3.14159);
0179   semiMu_muDRMin = bei.book1D("semiMu_muDRMin", "#Delta R(E,nearest jet) in #mu+Jets Channel", 50, 0, 10.0);
0180   semiMu_muPtRel = bei.book1D("semiMu_muPtRel", "p_{T}^{REL} in #mu+Jets Channel", 60, 0, 300.);
0181   semiMu_hadJetDR = bei.book1D("semiMu_hadJetDR", "#Delta R(E,had jet) in #mu+Jets Channel", 50, 0, 10.0);
0182   semiMu_hadJetPt =
0183       bei.book1D("semiMu_hadJetPt", "Pt of Leading Hadronic Jet in #mu+Jets Channel (GeV)", 50, 0.0, 1000);
0184   semiMu_hadJetY = bei.book1D("semiMu_hadJetY", "Rapidity of Leading Hadronic Jet in #mu+Jets Channel", 60, -6.0, 6.0);
0185   semiMu_hadJetPhi = bei.book1D(
0186       "semiMu_hadJetPhi", "#phi of Leading Hadronic Jet in #mu+Jets Channel (radians)", 60, -3.14159, 3.14159);
0187   semiMu_hadJetMass =
0188       bei.book1D("semiMu_hadJetMass", "Mass of Leading Hadronic Jet in #mu+Jets Channel (GeV)", 50, 0.0, 500);
0189   semiMu_hadJetWMass =
0190       bei.book1D("semiMu_hadJetwMass", "W Mass for Leading Hadronic Jet in #mu+Jets Channel (GeV)", 50, 0.0, 250.0);
0191   semiMu_mttbar = bei.book1D("semiMu_mttbar", "Mass of #mu+Jets ttbar Candidate", 100, 0., 5000.);
0192 
0193   //--- E+Jets
0194   bei.setCurrentFolder("Physics/B2G/SemiE");
0195   semiE_ePt = bei.book1D("semiE_ePt", "Pt of Electron in e+Jets Channel (GeV)", 50, 0.0, 1000);
0196   semiE_eEta = bei.book1D("semiE_eEta", "#eta of Electron in e+Jets Channel", 60, -6.0, 6.0);
0197   semiE_ePhi = bei.book1D("semiE_ePhi", "#phi of Electron in e+Jets Channel (radians)", 60, -3.14159, 3.14159);
0198   semiE_eDRMin = bei.book1D("semiE_eDRMin", "#Delta R(E,nearest jet) in e+Jets Channel", 50, 0, 10.0);
0199   semiE_ePtRel = bei.book1D("semiE_ePtRel", "p_{T}^{REL} in e+Jets Channel", 60, 0, 300.);
0200   semiE_hadJetDR = bei.book1D("semiE_hadJetDR", "#Delta R(E,had jet) in e+Jets Channel", 50, 0, 10.0);
0201   semiE_hadJetPt = bei.book1D("semiE_hadJetPt", "Pt of Leading Hadronic Jet in e+Jets Channel (GeV)", 50, 0.0, 1000);
0202   semiE_hadJetY = bei.book1D("semiE_hadJetY", "Rapidity of Leading Hadronic Jet in e+Jets Channel", 60, -6.0, 6.0);
0203   semiE_hadJetPhi =
0204       bei.book1D("semiE_hadJetPhi", "#phi of Leading Hadronic Jet in e+Jets Channel (radians)", 60, -3.14159, 3.14159);
0205   semiE_hadJetMass =
0206       bei.book1D("semiE_hadJetMass", "Mass of Leading Hadronic Jet in e+Jets Channel (GeV)", 50, 0.0, 500);
0207   semiE_hadJetWMass =
0208       bei.book1D("semiE_hadJetwMass", "W Mass for Leading Hadronic Jet in e+Jets Channel (GeV)", 50, 0.0, 250.0);
0209   semiE_mttbar = bei.book1D("semiE_mttbar", "Mass of e+Jets ttbar Candidate", 100, 0., 5000.);
0210 
0211   //--- All-hadronic
0212   bei.setCurrentFolder("Physics/B2G/AllHad");
0213   allHad_pt0 = bei.book1D("allHad_pt0", "Pt of Leading All-Hadronic PFJet (GeV)", 50, 0.0, 1000);
0214   allHad_y0 = bei.book1D("allHad_y0", "Rapidity of Leading All-Hadronic PFJet", 60, -6.0, 6.0);
0215   allHad_phi0 = bei.book1D("allHad_phi0", "#phi of Leading All-Hadronic PFJet (radians)", 60, -3.14159, 3.14159);
0216   allHad_mass0 = bei.book1D("allHad_mass0", "Mass of Leading All-Hadronic PFJet (GeV)", 50, 0.0, 500);
0217   allHad_wMass0 = bei.book1D("allHad_wMass0", "W Mass for Leading All-Hadronic PFJet (GeV)", 50, 0.0, 250.0);
0218   allHad_pt1 = bei.book1D("allHad_pt1", "Pt of Subleading All-Hadronic PFJet (GeV)", 50, 0.0, 1000);
0219   allHad_y1 = bei.book1D("allHad_y1", "Rapidity of Subleading All-Hadronic PFJet", 60, -6.0, 6.0);
0220   allHad_phi1 = bei.book1D("allHad_phi1", "#phi of Subleading All-Hadronic PFJet (radians)", 60, -3.14159, 3.14159);
0221   allHad_mass1 = bei.book1D("allHad_mass1", "Mass of Subleading All-Hadronic PFJet (GeV)", 50, 0.0, 500);
0222   allHad_wMass1 = bei.book1D("allHad_wMass1", "W Mass for Subleading All-Hadronic PFJet (GeV)", 50, 0.0, 250.0);
0223   allHad_mttbar = bei.book1D("allHad_mttbar", "Mass of All-Hadronic ttbar Candidate", 100, 0., 5000.);
0224 }
0225 
0226 //
0227 //  -- Analyze
0228 //
0229 void B2GDQM::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0230   analyzeJets(iEvent, iSetup);
0231   analyzeSemiMu(iEvent, iSetup);
0232   analyzeSemiE(iEvent, iSetup);
0233   analyzeAllHad(iEvent, iSetup);
0234 }
0235 
0236 void B2GDQM::analyzeJets(const Event& iEvent, const edm::EventSetup& iSetup) {
0237   // Loop over the different types of jets,
0238   //   Loop over the jets in that collection,
0239   //     fill PF jet information as well as substructure
0240   //     information for boosted jets.
0241   // Utilizes the CMS top-tagging algorithm and the "mass drop" W-tagger.
0242   for (unsigned int icoll = 0; icoll < jetLabels_.size(); ++icoll) {
0243     edm::Handle<edm::View<reco::Jet> > pfJetCollection;
0244     bool ValidPFJets = iEvent.getByToken(jetTokens_[icoll], pfJetCollection);
0245     if (!ValidPFJets)
0246       continue;
0247     edm::View<reco::Jet> const& pfjets = *pfJetCollection;
0248 
0249     // Jet Correction
0250     int countJet = 0;
0251     // const JetCorrector* pfcorrector =
0252     // JetCorrector::getJetCorrector(PFJetCorService_,iSetup);
0253 
0254     for (edm::View<reco::Jet>::const_iterator jet = pfjets.begin(), jetEnd = pfjets.end(); jet != jetEnd; ++jet) {
0255       if (jet->pt() < jetPtMins_[icoll])
0256         continue;
0257       pfJet_pt[icoll]->Fill(jet->pt());
0258       pfJet_y[icoll]->Fill(jet->rapidity());
0259       pfJet_phi[icoll]->Fill(jet->phi());
0260       pfJet_m[icoll]->Fill(jet->mass());
0261 
0262       // Dynamic cast the base class (reco::Jet) to the derived class (PFJet)
0263       // to access the PFJet information
0264       reco::PFJet const* pfjet = dynamic_cast<reco::PFJet const*>(&*jet);
0265 
0266       if (pfjet != nullptr) {
0267         pfJet_chef[icoll]->Fill(pfjet->chargedHadronEnergyFraction());
0268         pfJet_nhef[icoll]->Fill(pfjet->neutralHadronEnergyFraction());
0269         pfJet_cemf[icoll]->Fill(pfjet->chargedEmEnergyFraction());
0270         pfJet_nemf[icoll]->Fill(pfjet->neutralEmEnergyFraction());
0271       }
0272 
0273       // Dynamic cast the base class (reco::Jet) to the derived class (BasicJet)
0274       // to access the substructure information
0275       reco::BasicJet const* basicjet = dynamic_cast<reco::BasicJet const*>(&*jet);
0276 
0277       if (basicjet != nullptr) {
0278         boostedJet_subjetN[icoll]->Fill(jet->numberOfDaughters());
0279 
0280         for (unsigned int ida = 0; ida < jet->numberOfDaughters(); ++ida) {
0281           reco::Candidate const* subjet = jet->daughter(ida);
0282           boostedJet_subjetPt[icoll]->Fill(subjet->pt());
0283           boostedJet_subjetY[icoll]->Fill(subjet->rapidity());
0284           boostedJet_subjetPhi[icoll]->Fill(subjet->phi());
0285           boostedJet_subjetM[icoll]->Fill(subjet->mass());
0286         }
0287         // Check the various tagging algorithms
0288         if ((jetLabels_[icoll].label() == "ak8PFJetsPuppiSoftdrop")) {
0289           if (jet->numberOfDaughters() > 1) {
0290             reco::Candidate const* da0 = jet->daughter(0);
0291             reco::Candidate const* da1 = jet->daughter(1);
0292             if (da0->mass() > da1->mass()) {
0293               boostedJet_wMass[icoll]->Fill(da0->mass());
0294               boostedJet_massDrop[icoll]->Fill(da0->mass() / jet->mass());
0295             } else {
0296               boostedJet_wMass[icoll]->Fill(da1->mass());
0297               boostedJet_massDrop[icoll]->Fill(da1->mass() / jet->mass());
0298             }
0299 
0300           } else {
0301             boostedJet_massDrop[icoll]->Fill(-1.0);
0302           }
0303 
0304         }  // end if collection is AK8 PFJets Puppi soft-drop
0305 
0306       }  // end if basic jet != 0
0307       countJet++;
0308     }
0309   }
0310 
0311   // PFMETs
0312   edm::Handle<std::vector<reco::PFMET> > pfMETCollection;
0313   bool ValidPFMET = iEvent.getByToken(PFMETToken_, pfMETCollection);
0314   if (!ValidPFMET)
0315     return;
0316 
0317   pfMet_pt->Fill((*pfMETCollection)[0].pt());
0318   pfMet_phi->Fill((*pfMETCollection)[0].phi());
0319 }
0320 
0321 void B2GDQM::analyzeAllHad(const Event& iEvent, const edm::EventSetup& iSetup) {
0322   edm::Handle<edm::View<reco::BasicJet> > jetCollection;
0323   bool validJets = iEvent.getByToken(sdjetToken_, jetCollection);
0324   if (!validJets)
0325     return;
0326 
0327   // Require two back-to-back jets at high pt with |delta y| < 1.0
0328   if (jetCollection->size() < 2)
0329     return;
0330   edm::Ptr<reco::BasicJet> jet0 = jetCollection->ptrAt(0);
0331   edm::Ptr<reco::BasicJet> jet1 = jetCollection->ptrAt(1);
0332   if (jet0.isAvailable() == false || jet1.isAvailable() == false)
0333     return;
0334   if (jet0->pt() < allHadPtCut_ || jet1->pt() < allHadPtCut_)
0335     return;
0336   if (std::abs(jet0->rapidity() - jet1->rapidity()) > allHadRapidityCut_)
0337     return;
0338   if (std::abs(reco::deltaPhi(jet0->phi(), jet1->phi())) < M_PI * 0.5)
0339     return;
0340 
0341   allHad_pt0->Fill(jet0->pt());
0342   allHad_y0->Fill(jet0->rapidity());
0343   allHad_phi0->Fill(jet0->phi());
0344   allHad_mass0->Fill(jet0->mass());
0345   if (jet0->numberOfDaughters() > 2) {
0346     double wMass =
0347         jet0->daughter(0)->mass() >= jet0->daughter(1)->mass() ? jet0->daughter(0)->mass() : jet0->daughter(1)->mass();
0348     allHad_wMass0->Fill(wMass);
0349   } else {
0350     allHad_wMass0->Fill(-1.0);
0351   }
0352 
0353   allHad_pt1->Fill(jet1->pt());
0354   allHad_y1->Fill(jet1->rapidity());
0355   allHad_phi1->Fill(jet1->phi());
0356   allHad_mass1->Fill(jet1->mass());
0357   if (jet1->numberOfDaughters() > 2) {
0358     double wMass =
0359         jet1->daughter(0)->mass() >= jet1->daughter(1)->mass() ? jet1->daughter(0)->mass() : jet1->daughter(1)->mass();
0360     allHad_wMass1->Fill(wMass);
0361   } else {
0362     allHad_wMass1->Fill(-1.0);
0363   }
0364 
0365   auto p4cand = (jet0->p4() + jet1->p4());
0366   allHad_mttbar->Fill(p4cand.mass());
0367 }
0368 
0369 void B2GDQM::analyzeSemiMu(const Event& iEvent, const edm::EventSetup& iSetup) {
0370   edm::Handle<edm::View<reco::Muon> > muonCollection;
0371   bool validMuons = iEvent.getByToken(muonToken_, muonCollection);
0372 
0373   if (!validMuons)
0374     return;
0375   if (muonCollection->empty())
0376     return;
0377   reco::Muon const& muon = muonCollection->at(0);
0378   if (!(*muonSelect_)(muon))
0379     return;
0380 
0381   edm::Handle<edm::View<reco::BasicJet> > jetCollection;
0382   bool validJets = iEvent.getByToken(sdjetToken_, jetCollection);
0383   if (!validJets)
0384     return;
0385   if (jetCollection->size() < 2)
0386     return;
0387 
0388   double pt0 = -1.0;
0389   double dRMin = 999.0;
0390   edm::Ptr<reco::BasicJet> hadJet;  // highest pt jet with dphi(lep,jet) > pi/2
0391   edm::Ptr<reco::BasicJet> lepJet;  // closest jet to lepton with pt > ptMin
0392 
0393   for (auto ijet = jetCollection->begin(), ijetBegin = ijet, ijetEnd = jetCollection->end(); ijet != ijetEnd; ++ijet) {
0394     // Hadronic jets
0395     if (std::abs(reco::deltaPhi(muon, *ijet)) > M_PI * 0.5) {
0396       if (ijet->pt() > pt0 && ijet->p() > semiMu_HadJetPtCut_) {
0397         hadJet = jetCollection->ptrAt(ijet - ijetBegin);
0398         pt0 = hadJet->pt();
0399       }
0400     }
0401     // Leptonic jets
0402     else if (ijet->pt() > semiMu_LepJetPtCut_) {
0403       auto idRMin = reco::deltaR(muon, *ijet);
0404       if (idRMin < dRMin) {
0405         lepJet = jetCollection->ptrAt(ijet - ijetBegin);
0406         dRMin = idRMin;
0407       }
0408     }
0409   }
0410   if (hadJet.isAvailable() == false || lepJet.isAvailable() == false)
0411     return;
0412 
0413   auto lepJetP4 = lepJet->p4();
0414   const auto& muonP4 = muon.p4();
0415 
0416   double tot = lepJetP4.mag2();
0417   double ss = muonP4.Dot(lepJet->p4());
0418   double per = muonP4.mag2();
0419   if (tot > 0.0)
0420     per -= ss * ss / tot;
0421   if (per < 0)
0422     per = 0;
0423   double ptRel = per;
0424   bool pass2D = dRMin > semiMu_dRMin_ || ptRel > semiMu_ptRel_;
0425 
0426   if (!pass2D)
0427     return;
0428 
0429   semiMu_muPt->Fill(muon.pt());
0430   semiMu_muEta->Fill(muon.eta());
0431   semiMu_muPhi->Fill(muon.phi());
0432   semiMu_muDRMin->Fill(dRMin);
0433   semiMu_muPtRel->Fill(ptRel);
0434 
0435   semiMu_hadJetDR->Fill(reco::deltaR(muon, *hadJet));
0436   semiMu_mttbar->Fill(0.0);
0437 
0438   semiMu_hadJetPt->Fill(hadJet->pt());
0439   semiMu_hadJetY->Fill(hadJet->rapidity());
0440   semiMu_hadJetPhi->Fill(hadJet->phi());
0441   semiMu_hadJetMass->Fill(hadJet->mass());
0442   if (hadJet->numberOfDaughters() > 2) {
0443     double wMass = hadJet->daughter(0)->mass() >= hadJet->daughter(1)->mass() ? hadJet->daughter(0)->mass()
0444                                                                               : hadJet->daughter(1)->mass();
0445     semiMu_hadJetWMass->Fill(wMass);
0446   } else {
0447     semiMu_hadJetWMass->Fill(-1.0);
0448   }
0449 }
0450 
0451 void B2GDQM::analyzeSemiE(const Event& iEvent, const edm::EventSetup& iSetup) {
0452   edm::Handle<edm::View<reco::GsfElectron> > electronCollection;
0453   bool validElectrons = iEvent.getByToken(electronToken_, electronCollection);
0454 
0455   if (!validElectrons)
0456     return;
0457   if (electronCollection->empty())
0458     return;
0459   reco::GsfElectron const& electron = electronCollection->at(0);
0460   if (!(*elecSelect_)(electron))
0461     return;
0462 
0463   edm::Handle<edm::View<reco::BasicJet> > jetCollection;
0464   bool validJets = iEvent.getByToken(sdjetToken_, jetCollection);
0465   if (!validJets)
0466     return;
0467   if (jetCollection->size() < 2)
0468     return;
0469 
0470   double pt0 = -1.0;
0471   double dRMin = 999.0;
0472   edm::Ptr<reco::BasicJet> hadJet;  // highest pt jet with dphi(lep,jet) > pi/2
0473   edm::Ptr<reco::BasicJet> lepJet;  // closest jet to lepton with pt > ptMin
0474 
0475   for (auto ijet = jetCollection->begin(), ijetBegin = ijet, ijetEnd = jetCollection->end(); ijet != ijetEnd; ++ijet) {
0476     // Hadronic jets
0477     if (std::abs(reco::deltaPhi(electron, *ijet)) > M_PI * 0.5) {
0478       if (ijet->pt() > pt0 && ijet->p() > semiE_HadJetPtCut_) {
0479         hadJet = jetCollection->ptrAt(ijet - ijetBegin);
0480         pt0 = hadJet->pt();
0481       }
0482     }
0483     // Leptonic jets
0484     else if (ijet->pt() > semiE_LepJetPtCut_) {
0485       auto idRMin = reco::deltaR(electron, *ijet);
0486       if (idRMin < dRMin) {
0487         lepJet = jetCollection->ptrAt(ijet - ijetBegin);
0488         dRMin = idRMin;
0489       }
0490     }
0491   }
0492   if (hadJet.isAvailable() == false || lepJet.isAvailable() == false)
0493     return;
0494 
0495   auto lepJetP4 = lepJet->p4();
0496   const auto& electronP4 = electron.p4();
0497 
0498   double tot = lepJetP4.mag2();
0499   double ss = electronP4.Dot(lepJet->p4());
0500   double per = electronP4.mag2();
0501   if (tot > 0.0)
0502     per -= ss * ss / tot;
0503   if (per < 0)
0504     per = 0;
0505   double ptRel = per;
0506   bool pass2D = dRMin > semiE_dRMin_ || ptRel > semiE_ptRel_;
0507 
0508   if (!pass2D)
0509     return;
0510 
0511   semiE_ePt->Fill(electron.pt());
0512   semiE_eEta->Fill(electron.eta());
0513   semiE_ePhi->Fill(electron.phi());
0514   semiE_eDRMin->Fill(dRMin);
0515   semiE_ePtRel->Fill(ptRel);
0516 
0517   semiE_hadJetDR->Fill(reco::deltaR(electron, *hadJet));
0518   semiE_mttbar->Fill(0.0);
0519 
0520   semiE_hadJetPt->Fill(hadJet->pt());
0521   semiE_hadJetY->Fill(hadJet->rapidity());
0522   semiE_hadJetPhi->Fill(hadJet->phi());
0523   semiE_hadJetMass->Fill(hadJet->mass());
0524   if (hadJet->numberOfDaughters() > 2) {
0525     double wMass = hadJet->daughter(0)->mass() >= hadJet->daughter(1)->mass() ? hadJet->daughter(0)->mass()
0526                                                                               : hadJet->daughter(1)->mass();
0527     semiE_hadJetWMass->Fill(wMass);
0528   } else {
0529     semiE_hadJetWMass->Fill(-1.0);
0530   }
0531 }