File indexing completed on 2024-04-06 12:09:11
0001
0002 #include "DataFormats/BTauReco/interface/JetTag.h"
0003 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0004 #include "DataFormats/Common/interface/Handle.h"
0005 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0006 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0007 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0008 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0009 #include "DataFormats/JetReco/interface/PFJet.h"
0010 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0011 #include "DataFormats/JetReco/interface/PileupJetIdentifier.h"
0012 #include "DataFormats/METReco/interface/PFMET.h"
0013 #include "DataFormats/METReco/interface/PFMETCollection.h"
0014 #include "DataFormats/MuonReco/interface/Muon.h"
0015 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0016 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0017 #include "DataFormats/PatCandidates/interface/Jet.h"
0018 #include "DataFormats/TrackReco/interface/HitPattern.h"
0019 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0020 #include "FWCore/Framework/interface/Event.h"
0021 #include "FWCore/Framework/interface/Frameworkfwd.h"
0022 #include "FWCore/Framework/interface/stream/EDFilter.h"
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0026 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0027
0028
0029 #include "TLorentzVector.h"
0030
0031 class ttbarEventSelector : public edm::stream::EDFilter<> {
0032 public:
0033 explicit ttbarEventSelector(const edm::ParameterSet&);
0034
0035 bool filter(edm::Event&, edm::EventSetup const&) override;
0036 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0037
0038 private:
0039
0040 const edm::InputTag electronTag_;
0041 const edm::InputTag jetsTag_;
0042 const edm::InputTag bjetsTag_;
0043 const edm::InputTag pfmetTag_;
0044 const edm::InputTag muonTag_;
0045 const edm::InputTag bsTag_;
0046 const edm::EDGetTokenT<reco::GsfElectronCollection> electronToken_;
0047 const edm::EDGetTokenT<reco::PFJetCollection> jetsToken_;
0048 const edm::EDGetTokenT<reco::JetTagCollection> bjetsToken_;
0049 const edm::EDGetTokenT<reco::PFMETCollection> pfmetToken_;
0050 const edm::EDGetTokenT<reco::MuonCollection> muonToken_;
0051 const edm::EDGetTokenT<reco::BeamSpot> bsToken_;
0052
0053 const double maxEtaEle_;
0054 const double maxEtaMu_;
0055 const double minPt_;
0056 const double maxDeltaPhiInEB_;
0057 const double maxDeltaEtaInEB_;
0058 const double maxHOEEB_;
0059 const double maxSigmaiEiEEB_;
0060 const double maxDeltaPhiInEE_;
0061 const double maxDeltaEtaInEE_;
0062 const double maxHOEEE_;
0063 const double maxSigmaiEiEEE_;
0064
0065 const double minChambers_;
0066 const double minMatches_;
0067 const double minMatchedStations_;
0068
0069 const double maxEtaHighest_Jets_;
0070 const double maxEta_Jets_;
0071
0072 const double btagFactor_;
0073
0074 const double maxNormChi2_;
0075 const double maxD0_;
0076 const double maxDz_;
0077 const int minPixelHits_;
0078 const int minStripHits_;
0079 const double maxIsoEle_;
0080 const double maxIsoMu_;
0081 const double minPtHighestMu_;
0082 const double minPtHighestEle_;
0083 const double minPtHighest_Jets_;
0084 const double minPt_Jets_;
0085 const double minInvMass_;
0086 const double maxInvMass_;
0087 const double minMet_;
0088 const double maxMet_;
0089 const double minWmass_;
0090 const double maxWmass_;
0091 double getMt(const TLorentzVector& vlep, const reco::PFMET& obj);
0092 int EventCategory(int& nEle, int& nMu, int& nJets, int& nbJets);
0093 };
0094
0095 using namespace std;
0096 using namespace edm;
0097
0098 void ttbarEventSelector::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0099 edm::ParameterSetDescription desc;
0100 desc.addUntracked<edm::InputTag>("electronInputTag", edm::InputTag("gedGsfElectrons"));
0101 desc.addUntracked<edm::InputTag>("jetsInputTag", edm::InputTag("ak4PFJetsCHS"));
0102 desc.addUntracked<edm::InputTag>("bjetsInputTag", edm::InputTag("pfDeepCSVJetTags", "probb"));
0103 desc.addUntracked<edm::InputTag>("pfmetTag", edm::InputTag("pfMet"));
0104 desc.addUntracked<edm::InputTag>("muonInputTag", edm::InputTag("muons"));
0105 desc.addUntracked<edm::InputTag>("offlineBeamSpot", edm::InputTag("offlineBeamSpot"));
0106 desc.addUntracked<double>("maxEtaEle", 2.4);
0107 desc.addUntracked<double>("maxEtaMu", 2.4);
0108 desc.addUntracked<double>("minPt", 5);
0109 desc.addUntracked<double>("maxDeltaPhiInEB", .15);
0110 desc.addUntracked<double>("maxDeltaEtaInEB", .007);
0111 desc.addUntracked<double>("maxHOEEB", .12);
0112 desc.addUntracked<double>("maxSigmaiEiEEB", .01);
0113 desc.addUntracked<double>("maxDeltaPhiInEE", .1);
0114 desc.addUntracked<double>("maxDeltaEtaInEE", .009);
0115 desc.addUntracked<double>("maxHOEEB_", .10);
0116 desc.addUntracked<double>("maxSigmaiEiEEE", .03);
0117 desc.addUntracked<uint32_t>("minChambers", 2);
0118 desc.addUntracked<uint32_t>("minMatches", 2);
0119 desc.addUntracked<double>("minMatchedStations", 2);
0120 desc.addUntracked<double>("maxEtaHighest_Jets", 2.4);
0121 desc.addUntracked<double>("maxEta_Jets", 3.0);
0122 desc.addUntracked<double>("btagFactor", 0.6);
0123 desc.addUntracked<double>("maxNormChi2", 10);
0124 desc.addUntracked<double>("maxD0", 0.02);
0125 desc.addUntracked<double>("maxDz", 20.);
0126 desc.addUntracked<uint32_t>("minPixelHits", 1);
0127 desc.addUntracked<uint32_t>("minStripHits", 8);
0128 desc.addUntracked<double>("maxIsoEle", 0.5);
0129 desc.addUntracked<double>("maxIsoMu", 0.3);
0130 desc.addUntracked<double>("minPtHighestMu", 24);
0131 desc.addUntracked<double>("minPtHighestEle", 32);
0132 desc.addUntracked<double>("minPtHighest_Jets", 30);
0133 desc.addUntracked<double>("minPt_Jets", 20);
0134 desc.addUntracked<double>("minInvMass", 140);
0135 desc.addUntracked<double>("maxInvMass", 200);
0136 desc.addUntracked<double>("minMet", 50);
0137 desc.addUntracked<double>("maxMet", 80);
0138 desc.addUntracked<double>("minWmass", 50);
0139 desc.addUntracked<double>("maxWmass", 130);
0140 descriptions.addWithDefaultLabel(desc);
0141 }
0142
0143 ttbarEventSelector::ttbarEventSelector(const edm::ParameterSet& ps)
0144 : electronTag_(ps.getUntrackedParameter<edm::InputTag>("electronInputTag", edm::InputTag("gedGsfElectrons"))),
0145 jetsTag_(ps.getUntrackedParameter<edm::InputTag>("jetsInputTag", edm::InputTag("ak4PFJetsCHS"))),
0146 bjetsTag_(ps.getUntrackedParameter<edm::InputTag>("bjetsInputTag", edm::InputTag("pfDeepCSVJetTags", "probb"))),
0147 pfmetTag_(ps.getUntrackedParameter<edm::InputTag>("pfmetTag", edm::InputTag("pfMet"))),
0148 muonTag_(ps.getUntrackedParameter<edm::InputTag>("muonInputTag", edm::InputTag("muons"))),
0149 bsTag_(ps.getUntrackedParameter<edm::InputTag>("offlineBeamSpot", edm::InputTag("offlineBeamSpot"))),
0150 electronToken_(consumes<reco::GsfElectronCollection>(electronTag_)),
0151 jetsToken_(consumes<reco::PFJetCollection>(jetsTag_)),
0152 bjetsToken_(consumes<reco::JetTagCollection>(bjetsTag_)),
0153 pfmetToken_(consumes<reco::PFMETCollection>(pfmetTag_)),
0154 muonToken_(consumes<reco::MuonCollection>(muonTag_)),
0155 bsToken_(consumes<reco::BeamSpot>(bsTag_)),
0156 maxEtaEle_(ps.getUntrackedParameter<double>("maxEtaEle", 2.4)),
0157 maxEtaMu_(ps.getUntrackedParameter<double>("maxEtaMu", 2.4)),
0158 minPt_(ps.getUntrackedParameter<double>("minPt", 5)),
0159
0160
0161 maxDeltaPhiInEB_(ps.getUntrackedParameter<double>("maxDeltaPhiInEB", .15)),
0162 maxDeltaEtaInEB_(ps.getUntrackedParameter<double>("maxDeltaEtaInEB", .007)),
0163 maxHOEEB_(ps.getUntrackedParameter<double>("maxHOEEB", .12)),
0164 maxSigmaiEiEEB_(ps.getUntrackedParameter<double>("maxSigmaiEiEEB", .01)),
0165 maxDeltaPhiInEE_(ps.getUntrackedParameter<double>("maxDeltaPhiInEE", .1)),
0166 maxDeltaEtaInEE_(ps.getUntrackedParameter<double>("maxDeltaEtaInEE", .009)),
0167 maxHOEEE_(ps.getUntrackedParameter<double>("maxHOEEB_", .10)),
0168 maxSigmaiEiEEE_(ps.getUntrackedParameter<double>("maxSigmaiEiEEE", .03)),
0169
0170
0171 minChambers_(ps.getUntrackedParameter<uint32_t>("minChambers", 2)),
0172 minMatches_(ps.getUntrackedParameter<uint32_t>("minMatches", 2)),
0173 minMatchedStations_(ps.getUntrackedParameter<double>("minMatchedStations", 2)),
0174
0175
0176 maxEtaHighest_Jets_(ps.getUntrackedParameter<double>("maxEtaHighest_Jets", 2.4)),
0177 maxEta_Jets_(ps.getUntrackedParameter<double>("maxEta_Jets", 3.0)),
0178
0179
0180 btagFactor_(ps.getUntrackedParameter<double>("btagFactor", 0.6)),
0181
0182 maxNormChi2_(ps.getUntrackedParameter<double>("maxNormChi2", 10)),
0183 maxD0_(ps.getUntrackedParameter<double>("maxD0", 0.02)),
0184 maxDz_(ps.getUntrackedParameter<double>("maxDz", 20.)),
0185 minPixelHits_(ps.getUntrackedParameter<uint32_t>("minPixelHits", 1)),
0186 minStripHits_(ps.getUntrackedParameter<uint32_t>("minStripHits", 8)),
0187 maxIsoEle_(ps.getUntrackedParameter<double>("maxIsoEle", 0.5)),
0188 maxIsoMu_(ps.getUntrackedParameter<double>("maxIsoMu", 0.3)),
0189 minPtHighestMu_(ps.getUntrackedParameter<double>("minPtHighestMu", 24)),
0190 minPtHighestEle_(ps.getUntrackedParameter<double>("minPtHighestEle", 32)),
0191 minPtHighest_Jets_(ps.getUntrackedParameter<double>("minPtHighest_Jets", 30)),
0192 minPt_Jets_(ps.getUntrackedParameter<double>("minPt_Jets", 20)),
0193 minInvMass_(ps.getUntrackedParameter<double>("minInvMass", 140)),
0194 maxInvMass_(ps.getUntrackedParameter<double>("maxInvMass", 200)),
0195 minMet_(ps.getUntrackedParameter<double>("minMet", 50)),
0196 maxMet_(ps.getUntrackedParameter<double>("maxMet", 80)),
0197 minWmass_(ps.getUntrackedParameter<double>("minWmass", 50)),
0198 maxWmass_(ps.getUntrackedParameter<double>("maxWmass", 130)) {}
0199
0200 bool ttbarEventSelector::filter(edm::Event& iEvent, edm::EventSetup const& iSetup) {
0201
0202 edm::Handle<reco::BeamSpot> beamSpot;
0203 iEvent.getByToken(bsToken_, beamSpot);
0204
0205 int le = 0, lm = 0, lj = 0, lbj = 0;
0206
0207
0208 edm::Handle<reco::GsfElectronCollection> electronColl;
0209 iEvent.getByToken(electronToken_, electronColl);
0210 std::vector<TLorentzVector> list_ele;
0211 std::vector<int> chrgeList_ele;
0212
0213 if (electronColl.isValid()) {
0214 for (auto const& ele : *electronColl) {
0215 if (!ele.ecalDriven())
0216 continue;
0217 if (ele.pt() < minPt_)
0218 continue;
0219
0220 if (!(ele.isEB() || ele.isEE()))
0221 continue;
0222
0223 double hOverE = ele.hadronicOverEm();
0224 double sigmaee = ele.sigmaIetaIeta();
0225 double deltaPhiIn = ele.deltaPhiSuperClusterTrackAtVtx();
0226 double deltaEtaIn = ele.deltaEtaSuperClusterTrackAtVtx();
0227
0228
0229 if (ele.isEB()) {
0230 if (fabs(deltaPhiIn) >= maxDeltaPhiInEB_ && fabs(deltaEtaIn) >= maxDeltaEtaInEB_ && hOverE >= maxHOEEB_ &&
0231 sigmaee >= maxSigmaiEiEEB_)
0232 continue;
0233 } else if (ele.isEE()) {
0234 if (fabs(deltaPhiIn) >= maxDeltaPhiInEE_ && fabs(deltaEtaIn) >= maxDeltaEtaInEE_ && hOverE >= maxHOEEE_ &&
0235 sigmaee >= maxSigmaiEiEEE_)
0236 continue;
0237 }
0238
0239 reco::GsfTrackRef trk = ele.gsfTrack();
0240 if (!trk.isNonnull())
0241 continue;
0242 double chi2 = trk->chi2();
0243 double ndof = trk->ndof();
0244 double chbyndof = (ndof > 0) ? chi2 / ndof : 0;
0245 if (chbyndof >= maxNormChi2_)
0246 continue;
0247
0248 double trkd0 = trk->d0();
0249 if (beamSpot.isValid()) {
0250 trkd0 = -(trk->dxy(beamSpot->position()));
0251 } else {
0252 edm::LogError("ElectronSelector") << "Error >> Failed to get BeamSpot for label: " << bsTag_;
0253 }
0254 if (std::fabs(trkd0) >= maxD0_)
0255 continue;
0256
0257 const reco::HitPattern& hitp = trk->hitPattern();
0258 int nPixelHits = hitp.numberOfValidPixelHits();
0259 if (nPixelHits < minPixelHits_)
0260 continue;
0261
0262 int nStripHits = hitp.numberOfValidStripHits();
0263 if (nStripHits < minStripHits_)
0264 continue;
0265
0266
0267 reco::GsfElectron::PflowIsolationVariables pfIso = ele.pfIsolationVariables();
0268 const float eiso =
0269 pfIso.sumChargedHadronPt + std::max(0.0, pfIso.sumNeutralHadronEt + pfIso.sumPhotonEt - 0.5 * pfIso.sumPUPt);
0270 if (eiso > maxIsoEle_ * ele.pt())
0271 continue;
0272
0273 TLorentzVector lv_ele;
0274 lv_ele.SetPtEtaPhiE(ele.pt(), ele.eta(), ele.phi(), ele.energy());
0275 list_ele.push_back(lv_ele);
0276 chrgeList_ele.push_back(ele.charge());
0277 }
0278 le = list_ele.size();
0279 } else {
0280 edm::LogError("ElectronSelector") << "Error >> Failed to get ElectronCollection for label: " << electronTag_;
0281 }
0282
0283
0284 edm::Handle<reco::MuonCollection> muonColl;
0285 iEvent.getByToken(muonToken_, muonColl);
0286
0287 std::vector<TLorentzVector> list_mu;
0288 std::vector<int> chrgeList_mu;
0289 if (muonColl.isValid()) {
0290 for (auto const& mu : *muonColl) {
0291 if (!mu.isGlobalMuon())
0292 continue;
0293 if (!mu.isPFMuon())
0294 continue;
0295 if (std::fabs(mu.eta()) >= maxEtaMu_)
0296 continue;
0297 if (mu.pt() < minPt_)
0298 continue;
0299
0300 reco::TrackRef gtk = mu.globalTrack();
0301 double chi2 = gtk->chi2();
0302 double ndof = gtk->ndof();
0303 double chbyndof = (ndof > 0) ? chi2 / ndof : 0;
0304 if (chbyndof >= maxNormChi2_)
0305 continue;
0306
0307 reco::TrackRef tk = mu.innerTrack();
0308 if (beamSpot.isValid()) {
0309 double trkd0 = -(tk->dxy(beamSpot->position()));
0310 if (std::fabs(trkd0) >= maxD0_)
0311 continue;
0312 double trkdz = tk->dz(beamSpot->position());
0313 if (std::fabs(trkdz) >= maxDz_)
0314 continue;
0315 } else {
0316 edm::LogError("MuonSelector") << "Error >> Failed to get BeamSpot for label: " << bsTag_;
0317 }
0318
0319 const reco::HitPattern& hitp = gtk->hitPattern();
0320 if (hitp.numberOfValidPixelHits() < minPixelHits_)
0321 continue;
0322 if (hitp.numberOfValidStripHits() < minStripHits_)
0323 continue;
0324
0325
0326 if (mu.numberOfChambers() < minChambers_)
0327 continue;
0328 if (mu.numberOfMatches() < minMatches_)
0329 continue;
0330 if (mu.numberOfMatchedStations() < minMatchedStations_)
0331 continue;
0332 if (!muon::isGoodMuon(mu, muon::GlobalMuonPromptTight))
0333 continue;
0334
0335
0336 const reco::MuonPFIsolation& pfIso04 = mu.pfIsolationR04();
0337 double absiso = pfIso04.sumChargedHadronPt +
0338 std::max(0.0, pfIso04.sumNeutralHadronEt + pfIso04.sumPhotonEt - 0.5 * pfIso04.sumPUPt);
0339 if (absiso > maxIsoMu_ * mu.pt())
0340 continue;
0341
0342 TLorentzVector lv_mu;
0343 lv_mu.SetPtEtaPhiE(mu.pt(), mu.eta(), mu.phi(), mu.energy());
0344 list_mu.push_back(lv_mu);
0345 chrgeList_mu.push_back(mu.charge());
0346 }
0347 lm = list_mu.size();
0348 } else {
0349 edm::LogError("MuonSelector") << "Error >> Failed to get MuonCollection for label: " << muonTag_;
0350 return false;
0351 }
0352
0353
0354 edm::Handle<reco::PFJetCollection> jetColl;
0355 iEvent.getByToken(jetsToken_, jetColl);
0356 std::vector<TLorentzVector> list_jets;
0357
0358 if (jetColl.isValid()) {
0359 for (const auto& jets : *jetColl) {
0360 if (jets.pt() < minPt_Jets_)
0361 continue;
0362 if (std::fabs(jets.eta()) > maxEta_Jets_)
0363 continue;
0364 TLorentzVector lv_jets;
0365 lv_jets.SetPtEtaPhiE(jets.pt(), jets.eta(), jets.phi(), jets.energy());
0366 list_jets.push_back(lv_jets);
0367 }
0368 lj = list_jets.size();
0369 }
0370
0371 edm::Handle<reco::JetTagCollection> bTagHandle;
0372 iEvent.getByToken(bjetsToken_, bTagHandle);
0373 const reco::JetTagCollection& bTags = *(bTagHandle.product());
0374 std::vector<TLorentzVector> list_bjets;
0375
0376 if (!bTags.empty()) {
0377 for (unsigned bj = 0; bj != bTags.size(); ++bj) {
0378 TLorentzVector lv_bjets;
0379 lv_bjets.SetPtEtaPhiE(
0380 bTags[bj].first->pt(), bTags[bj].first->eta(), bTags[bj].first->phi(), bTags[bj].first->energy());
0381 if ((bTags[bj].second > btagFactor_) && (lv_bjets.Pt() > minPt_Jets_))
0382 list_bjets.push_back(lv_bjets);
0383 }
0384 lbj = list_bjets.size();
0385 }
0386
0387
0388 edm::Handle<reco::PFMETCollection> pfColl;
0389 iEvent.getByToken(pfmetToken_, pfColl);
0390 if (EventCategory(le, lm, lj, lbj) == 11) {
0391 if (list_ele[0].Pt() < minPtHighestEle_)
0392 return false;
0393 if ((list_ele[0].Pt() < list_mu[0].Pt()) || (list_ele[1].Pt() < list_mu[0].Pt()))
0394 return false;
0395 if (chrgeList_ele[0] + chrgeList_ele[1] != 0)
0396 return false;
0397 if (pfColl.isValid()) {
0398 double mt1 = getMt(list_ele[0], pfColl->front());
0399 double mt2 = getMt(list_ele[1], pfColl->front());
0400 double mt = mt1 + mt2;
0401 if (mt < 2 * minMet_ || mt > 2 * maxMet_)
0402 return false;
0403 } else {
0404 edm::LogError("ttbarEventSelector")
0405 << "Error >> Failed to get PFMETCollection in dilepton ele-ele channel for label: " << pfmetTag_;
0406 return false;
0407 }
0408 } else if (EventCategory(le, lm, lj, lbj) == 12) {
0409 if (list_mu[0].Pt() < minPtHighestMu_)
0410 return false;
0411 if ((list_mu[0].Pt() < list_ele[0].Pt()) || (list_mu[1].Pt() < list_ele[0].Pt()))
0412 return false;
0413 if (chrgeList_mu[0] + chrgeList_mu[1] != 0)
0414 return false;
0415 if (pfColl.isValid()) {
0416 double mt1 = getMt(list_mu[0], pfColl->front());
0417 double mt2 = getMt(list_mu[1], pfColl->front());
0418 double mt = mt1 + mt2;
0419 if (mt < 2 * minMet_ || mt > 2 * maxMet_)
0420 return false;
0421 } else {
0422 edm::LogError("ttbarEventSelector")
0423 << "Error >> Failed to get PFMETCollection in dilepton mu-mu channel for label: " << pfmetTag_;
0424 return false;
0425 }
0426 } else if (EventCategory(le, lm, lj, lbj) == 13) {
0427 if ((list_mu[0].Pt() < list_ele[1].Pt()) || (list_ele[0].Pt() < list_mu[1].Pt()))
0428 return false;
0429 if ((list_mu[0].Pt() < minPtHighestMu_) || (list_ele[0].Pt() < minPtHighestEle_))
0430 return false;
0431 if (chrgeList_mu[0] + chrgeList_ele[0] != 0)
0432 return false;
0433 if (pfColl.isValid()) {
0434 double mt1 = getMt(list_mu[0], pfColl->front());
0435 double mt2 = getMt(list_ele[0], pfColl->front());
0436 double mt = mt1 + mt2;
0437 if (mt < 2 * minMet_ || mt > 2 * maxMet_)
0438 return false;
0439 } else {
0440 edm::LogError("ttbarEventSelector")
0441 << "Error >> Failed to get PFMETCollection in dilepton ele-mu channel for label: " << pfmetTag_;
0442 return false;
0443 }
0444 }
0445 if (EventCategory(le, lm, lj, lbj) == 21) {
0446 if (list_jets[0].Pt() < minPtHighest_Jets_)
0447 return false;
0448
0449 if (list_ele[0].Pt() < minPtHighestEle_)
0450 return false;
0451
0452 if ((!list_ele.empty() && list_ele[0].Pt() > minPtHighestEle_) &&
0453 (!list_mu.empty() && list_mu[0].Pt() > minPtHighestMu_))
0454 return false;
0455
0456 return true;
0457 }
0458 if (EventCategory(le, lm, lj, lbj) == 22) {
0459 if (list_jets[0].Pt() < minPtHighest_Jets_)
0460 return false;
0461
0462 if (list_mu[0].Pt() < minPtHighestMu_)
0463 return false;
0464
0465 if ((!list_ele.empty() && list_ele[0].Pt() > minPtHighestEle_) &&
0466 (!list_mu.empty() && list_mu[0].Pt() > minPtHighestMu_))
0467 return false;
0468
0469 return true;
0470 } else if (EventCategory(le, lm, lj, lbj) == 30) {
0471 if (list_jets[0].Pt() < minPtHighest_Jets_)
0472 return false;
0473 for (int i = 0; i < 4; i++) {
0474 TLorentzVector vjet1;
0475 for (int j = i + 1; j < 4; j++) {
0476 TLorentzVector vjet2;
0477 vjet1 = list_jets[i];
0478 vjet2 = list_jets[j];
0479 TLorentzVector vjet = vjet1 + vjet2;
0480 if (vjet.M() < minWmass_ || vjet.M() > maxWmass_)
0481 return false;
0482 }
0483 }
0484 } else if (EventCategory(le, lm, lj, lbj) == 50)
0485 return false;
0486
0487 return false;
0488 }
0489
0490 int ttbarEventSelector::EventCategory(int& nEle, int& nMu, int& nJets, int& nbJets) {
0491 int cat = 0;
0492 if ((nEle >= 2 || nMu >= 2) && nJets > 1 && nbJets > 1) {
0493 if (nEle >= 2 && nJets < 2)
0494 cat = 11;
0495 else if (nMu >= 2 && nJets < 2)
0496 cat = 12;
0497 else if (nEle >= 1 && nMu >= 1 && nJets < 2)
0498 cat = 13;
0499 } else if ((nEle >= 1 || nMu >= 1) && nJets > 3 && nbJets > 2) {
0500 if (nEle >= 1 && nJets > 1)
0501 cat = 21;
0502 else if (nMu >= 1 && nJets > 1)
0503 cat = 22;
0504 } else if ((nEle < 1 && nMu < 1) && nJets > 5 && nbJets > 1)
0505 cat = 30;
0506 else
0507 cat = 50;
0508 return cat;
0509 }
0510
0511 double ttbarEventSelector::getMt(const TLorentzVector& vlep, const reco::PFMET& obj) {
0512 double met = obj.et();
0513 double phi = obj.phi();
0514
0515 TLorentzVector vmet;
0516 double metx = met * std::cos(phi);
0517 double mety = met * std::sin(phi);
0518 vmet.SetPxPyPzE(metx, mety, 0.0, met);
0519
0520
0521 TLorentzVector vw = vlep + vmet;
0522
0523 return std::sqrt(2 * vlep.Et() * met * (1 - std::cos(deltaPhi(vlep.Phi(), phi))));
0524 }
0525
0526
0527 #include "FWCore/Framework/interface/MakerMacros.h"
0528 DEFINE_FWK_MODULE(ttbarEventSelector);