Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:12:25

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