TrackingMonitorSource/plugins/WtoLNuSelector.cc

0001 // W->lnu Event Selector
0002
0003 // user includes
0004 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0005 #include "DataFormats/Common/interface/Handle.h"
0006 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0007 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0008 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0009 #include "DataFormats/METReco/interface/PFMET.h"
0010 #include "DataFormats/METReco/interface/PFMETCollection.h"
0011 #include "DataFormats/METReco/interface/PFMETFwd.h"
0012 #include "DataFormats/Math/interface/deltaPhi.h"
0013 #include "DataFormats/Math/interface/deltaR.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/TrackReco/interface/HitPattern.h"
0018 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/Frameworkfwd.h"
0021 #include "FWCore/Framework/interface/stream/EDFilter.h"
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0024 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0026
0027 // ROOT includes
0028 #include "TLorentzVector.h"
0029 #include "TMath.h"
0030
0031 class WtoLNuSelector : public edm::stream::EDFilter<> {
0032 public:
0033   explicit WtoLNuSelector(const edm::ParameterSet&);
0034   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0035
0036 private:
0037   bool filter(edm::Event&, edm::EventSetup const&) override;
0038   double getMt(const TLorentzVector& vlep, const reco::PFMET& obj);
0039
0040 private:
0041   // module config parameters
0042   const edm::InputTag electronTag_;
0043   const edm::InputTag bsTag_;
0044   const edm::InputTag muonTag_;
0045   const edm::InputTag pfmetTag_;
0046   const edm::EDGetTokenT<reco::GsfElectronCollection> electronToken_;
0047   const edm::EDGetTokenT<reco::BeamSpot> bsToken_;
0048   const edm::EDGetTokenT<reco::MuonCollection> muonToken_;
0049   const edm::EDGetTokenT<reco::PFMETCollection> pfmetToken_;
0050 };
0051
0052 void WtoLNuSelector::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0053   edm::ParameterSetDescription desc;
0054   desc.addUntracked<edm::InputTag>("electronInputTag", edm::InputTag("gedGsfElectrons"));
0055   desc.addUntracked<edm::InputTag>("offlineBeamSpot", edm::InputTag("offlineBeamSpot"));
0056   desc.addUntracked<edm::InputTag>("muonInputTag", edm::InputTag("muons"));
0057   desc.addUntracked<edm::InputTag>("pfmetTag", edm::InputTag("pfMetT1T2Txy"));
0058   descriptions.addWithDefaultLabel(desc);
0059 }
0060
0061 WtoLNuSelector::WtoLNuSelector(const edm::ParameterSet& ps)
0062     : electronTag_(ps.getUntrackedParameter<edm::InputTag>("electronInputTag", edm::InputTag("gedGsfElectrons"))),
0063       bsTag_(ps.getUntrackedParameter<edm::InputTag>("offlineBeamSpot", edm::InputTag("offlineBeamSpot"))),
0064       muonTag_(ps.getUntrackedParameter<edm::InputTag>("muonInputTag", edm::InputTag("muons"))),
0065       pfmetTag_(ps.getUntrackedParameter<edm::InputTag>("pfmetTag", edm::InputTag("pfMetT1T2Txy"))),
0066       electronToken_(consumes<reco::GsfElectronCollection>(electronTag_)),
0067       bsToken_(consumes<reco::BeamSpot>(bsTag_)),
0068       muonToken_(consumes<reco::MuonCollection>(muonTag_)),
0069       pfmetToken_(consumes<reco::PFMETCollection>(pfmetTag_)) {}
0070
0071 bool WtoLNuSelector::filter(edm::Event& iEvent, edm::EventSetup const& iSetup) {
0072   // Read Electron Collection
0073   edm::Handle<reco::GsfElectronCollection> electronColl;
0074   iEvent.getByToken(electronToken_, electronColl);
0075
0076   edm::Handle<reco::BeamSpot> beamSpot;
0077   iEvent.getByToken(bsToken_, beamSpot);
0078
0079   std::vector<TLorentzVector> eleList;
0080   if (electronColl.isValid()) {
0081     for (auto const& ele : *electronColl) {
0082       if (!ele.ecalDriven())
0083         continue;
0084
0085       double hOverE = ele.hadronicOverEm();
0086       double sigmaee = ele.sigmaIetaIeta();
0087       double deltaPhiIn = ele.deltaPhiSuperClusterTrackAtVtx();
0088       double deltaEtaIn = ele.deltaEtaSuperClusterTrackAtVtx();
0089
0090       // separate cut for barrel and endcap
0091       if (ele.isEB()) {
0092         if (std::fabs(deltaPhiIn) >= .15 && std::fabs(deltaEtaIn) >= .007 && hOverE >= .12 && sigmaee >= .01)
0093           continue;
0094       } else if (ele.isEE()) {
0095         if (std::fabs(deltaPhiIn) >= .10 && std::fabs(deltaEtaIn) >= .009 && hOverE >= .10 && sigmaee >= .03)
0096           continue;
0097       }
0098
0099       reco::GsfTrackRef trk = ele.gsfTrack();
0100       if (!trk.isNonnull())
0101         continue;  // only electrons wd tracks
0102       double chi2 = trk->chi2();
0103       double ndof = trk->ndof();
0104       double chbyndof = (ndof > 0) ? chi2 / ndof : 0;
0105       double trkd0 = trk->d0();
0106       double trkdz = trk->dz();
0107       if (beamSpot.isValid()) {
0108         trkd0 = -(trk->dxy(beamSpot->position()));
0109         trkdz = trk->dz(beamSpot->position());
0110       } else {
0111         edm::LogError("WtoLNuSelector") << "Error >> Failed to get BeamSpot for label: " << bsTag_;
0112       }
0113       if (chbyndof >= 10 || std::fabs(trkd0) >= 0.02 || std::fabs(trkdz) >= 20)
0114         continue;
0115       const reco::HitPattern& hitp = trk->hitPattern();
0116       int nPixelHits = hitp.numberOfValidPixelHits();
0117       int nStripHits = hitp.numberOfValidStripHits();
0118       if (nPixelHits < 1 || nStripHits < 8)
0119         continue;
0120
0121       // PF Isolation
0122       reco::GsfElectron::PflowIsolationVariables pfIso = ele.pfIsolationVariables();
0123       float absiso =
0124           pfIso.sumChargedHadronPt + std::max(0.0, pfIso.sumNeutralHadronEt + pfIso.sumPhotonEt - 0.5 * pfIso.sumPUPt);
0125       float eiso = absiso / ele.pt();
0126       if (eiso > 0.3)
0127         continue;
0128
0129       TLorentzVector le;
0130       le.SetPtEtaPhiE(ele.pt(), ele.eta(), ele.phi(), ele.energy());
0131       eleList.push_back(le);
0132     }
0133   } else {
0134     edm::LogError("WtoLNuSelector") << "Error >> Failed to get ElectronCollection for label: " << electronTag_;
0135   }
0136
0137   // Read Muon Collection
0138   edm::Handle<reco::MuonCollection> muonColl;
0139   iEvent.getByToken(muonToken_, muonColl);
0140
0141   std::vector<TLorentzVector> muList;
0142   if (muonColl.isValid()) {
0143     for (auto const& mu : *muonColl) {
0144       if (!mu.isGlobalMuon() || !mu.isPFMuon() || std::fabs(mu.eta()) > 2.1 || mu.pt() <= 5)
0145         continue;
0146
0147       reco::TrackRef gtrkref = mu.globalTrack();
0148       if (!gtrkref.isNonnull())
0149         continue;
0150       const reco::Track* gtk = &(*gtrkref);
0151       double chi2 = gtk->chi2();
0152       double ndof = gtk->ndof();
0153       double chbyndof = (ndof > 0) ? chi2 / ndof : 0;
0154
0155       const reco::HitPattern& hitp = gtk->hitPattern();
0156       int nPixelHits = hitp.numberOfValidPixelHits();
0157       int nStripHits = hitp.numberOfValidStripHits();
0158
0159       reco::TrackRef itrkref = mu.innerTrack();  // tracker segment only
0160       if (!itrkref.isNonnull())
0161         continue;
0162       const reco::Track* tk = &(*itrkref);
0163       double trkd0 = tk->d0();
0164       double trkdz = tk->dz();
0165       if (beamSpot.isValid()) {
0166         trkd0 = -(tk->dxy(beamSpot->position()));
0167         trkdz = tk->dz(beamSpot->position());
0168       }
0169       // Hits/section in the muon chamber
0170       int nChambers = mu.numberOfChambers();
0171       int nMatches = mu.numberOfMatches();
0172       int nMatchedStations = mu.numberOfMatchedStations();
0173
0174       // PF Isolation
0175       const reco::MuonPFIsolation& pfIso04 = mu.pfIsolationR04();
0176       double absiso = pfIso04.sumChargedParticlePt +
0177                       std::max(0.0, pfIso04.sumNeutralHadronEt + pfIso04.sumPhotonEt - 0.5 * pfIso04.sumPUPt);
0178
0179       // Final selection
0180       if (chbyndof < 10 && std::fabs(trkd0) < 0.02 && std::fabs(trkdz) < 20.0 && nPixelHits > 1 && nStripHits > 8 &&
0181           nChambers > 2 && nMatches > 2 && nMatchedStations > 2 && absiso / mu.pt() < 0.3) {
0182         TLorentzVector lm;
0183         lm.SetPtEtaPhiE(mu.pt(), mu.eta(), mu.phi(), mu.energy());
0184         muList.push_back(lm);
0185       }
0186     }
0187   } else {
0188     edm::LogError("WtoLNuSelector") << "Error >> Failed to get MuonCollection for label: " << muonTag_;
0189   }
0190
0191   // Require either a high pt electron or muon
0192   if (eleList.empty() && muList.empty())
0193     return false;
0194
0195   // Both should not be present at the same time
0196   if ((!eleList.empty() && eleList[0].Pt() > 20) && (!muList.empty() && muList[0].Pt() > 20))
0197     return false;
0198
0199   // find the high pt lepton
0200   TLorentzVector vlep;
0201   if (!eleList.empty() && !muList.empty()) {
0202     vlep = (eleList[0].Pt() > muList[0].Pt()) ? eleList[0] : muList[0];
0203   } else if (!eleList.empty()) {
0204     vlep = eleList[0];
0205   } else {
0206     vlep = muList[0];
0207   }
0208   if (vlep.Pt() < 20)
0209     return false;
0210
0211   edm::Handle<reco::PFMETCollection> pfColl;
0212   iEvent.getByToken(pfmetToken_, pfColl);
0213
0214   if (pfColl.isValid()) {
0215     double mt = getMt(vlep, pfColl->front());
0216     if (mt < 60 || mt > 80)
0217       return false;
0218   } else {
0219     edm::LogError("WtoLNuSelector") << "Error >> Failed to get PFMETCollection for label: " << pfmetTag_;
0220     return false;
0221   }
0222
0223   return true;
0224 }
0225 double WtoLNuSelector::getMt(const TLorentzVector& vlep, const reco::PFMET& obj) {
0226   double met = obj.et();
0227   double phi = obj.phi();
0228
0229   TLorentzVector vmet;
0230   double metx = met * std::cos(phi);
0231   double mety = met * std::sin(phi);
0232   vmet.SetPxPyPzE(metx, mety, 0.0, met);
0233
0234   // transverse mass
0235   TLorentzVector vw = vlep + vmet;
0236
0237   return std::sqrt(2 * vlep.Et() * met * (1 - std::cos(deltaPhi(vlep.Phi(), phi))));
0238 }
0239 // Define this as a plug-in
0240 #include "FWCore/Framework/interface/MakerMacros.h"
0241 DEFINE_FWK_MODULE(WtoLNuSelector);