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