![]() |
![]() |
#include <BTagSkimLeptonJet.h>
Classes | |
class | PtSorter |
Public Member Functions | |
BTagSkimLeptonJet (const edm::ParameterSet &) | |
virtual void | endJob () |
virtual bool | filter (edm::Event &, const edm::EventSetup &) |
~BTagSkimLeptonJet () | |
Private Attributes | |
edm::InputTag | CaloJetInput_ |
edm::InputTag | LeptonInput_ |
std::string | LeptonType_ |
double | MaxCaloJetEta_ |
double | MaxDeltaR_ |
double | MaxLeptonEta_ |
double | MinCaloJetPt_ |
double | MinLeptonPt_ |
int | MinNLeptonJet_ |
double | MinPtRel_ |
unsigned int | nAccepted_ |
unsigned int | nEvents_ |
Definition at line 22 of file BTagSkimLeptonJet.h.
BTagSkimLeptonJet::BTagSkimLeptonJet | ( | const edm::ParameterSet & | iConfig | ) | [explicit] |
Definition at line 37 of file BTagSkimLeptonJet.cc.
References CaloJetInput_, edm::ParameterSet::getParameter(), LeptonInput_, LeptonType_, MaxCaloJetEta_, MaxDeltaR_, MaxLeptonEta_, MinCaloJetPt_, MinLeptonPt_, MinNLeptonJet_, and MinPtRel_.
: nEvents_(0), nAccepted_(0) { CaloJetInput_ = iConfig.getParameter<InputTag>( "CaloJet" ); MinCaloJetPt_ = iConfig.getParameter<double>( "MinimumCaloJetPt" ); MaxCaloJetEta_ = iConfig.getParameter<double>( "MaximumCaloJetEta" ); LeptonType_ = iConfig.getParameter<std::string>("LeptonType"); if ( LeptonType_ != "electron" && LeptonType_ != "muon" ) edm::LogError( "BTagSkimLeptonJet" ) << "unknown lepton type !!"; LeptonInput_ = iConfig.getParameter<InputTag>( "Lepton" ); MinLeptonPt_ = iConfig.getParameter<double>( "MinimumLeptonPt" ); MaxLeptonEta_ = iConfig.getParameter<double>( "MaximumLeptonEta" ); //MinNLepton_ = iConfig.getParameter<int>( "MinimumNLepton" ); MaxDeltaR_ = iConfig.getParameter<double>( "MaximumDeltaR" ); MinPtRel_ = iConfig.getParameter<double>("MinimumPtRel" ); MinNLeptonJet_ = iConfig.getParameter<int>( "MinimumNLeptonJet" ); if ( MinNLeptonJet_ < 1 ) edm::LogError( "BTagSkimLeptonJet" ) << "MinimunNCaloJet < 1 !!"; }
BTagSkimLeptonJet::~BTagSkimLeptonJet | ( | ) |
Definition at line 65 of file BTagSkimLeptonJet.cc.
{}
void BTagSkimLeptonJet::endJob | ( | void | ) | [virtual] |
Reimplemented from edm::EDFilter.
Definition at line 186 of file BTagSkimLeptonJet.cc.
References LeptonType_, nAccepted_, and nEvents_.
{ edm::LogVerbatim( "BTagSkimLeptonJet" ) << "=============================================================================\n" << " Events read: " << nEvents_ << "\n Events accepted by (" << LeptonType_ << ") BTagSkimLeptonJet: " << nAccepted_ << "\n Efficiency: " << (double)(nAccepted_)/(double)(nEvents_) << "\n===========================================================================" << endl; }
bool BTagSkimLeptonJet::filter | ( | edm::Event & | iEvent, |
const edm::EventSetup & | iSetup | ||
) | [virtual] |
Implements edm::EDFilter.
Definition at line 70 of file BTagSkimLeptonJet.cc.
References CaloJetInput_, edm::Event::getByLabel(), LeptonInput_, LeptonType_, MaxCaloJetEta_, MaxDeltaR_, MaxLeptonEta_, MinCaloJetPt_, MinLeptonPt_, MinNLeptonJet_, MinPtRel_, nAccepted_, and nEvents_.
{ nEvents_++; Handle<CaloJetCollection> CaloJetsHandle; //try { iEvent.getByLabel( CaloJetInput_, CaloJetsHandle ); //} //catch ( cms::Exception& ex ) { //edm::LogError( "BTagSkimLeptonJet" ) // << "Unable to get CaloJet collection " // << CaloJetInput_.label(); //return false; //} if ( CaloJetsHandle->empty() ) return false; CaloJetCollection TheCaloJets = *CaloJetsHandle; std::stable_sort( TheCaloJets.begin(), TheCaloJets.end(), PtSorter() ); //int nLepton = 0; Handle<MuonCollection> MuonHandle; MuonCollection TheMuons; if (LeptonType_ == "muon") { //try{ iEvent.getByLabel( LeptonInput_, MuonHandle ); // } // catch ( cms::Exception& ex ) { //edm::LogError( "BTagSkimLeptonJet" ) // << "Unable to get muon collection " // << LeptonInput_.label(); //return false; //} TheMuons = *MuonHandle; std::stable_sort( TheMuons.begin(), TheMuons.end(), PtSorter() ); } Handle<GsfElectronCollection> ElectronHandle; GsfElectronCollection TheElectrons; if (LeptonType_ == "electron") { //try{ iEvent.getByLabel( LeptonInput_, ElectronHandle ); // } // catch ( cms::Exception& ex ) { // edm::LogError( "BTagSkimLeptonJet" ) // << "Unable to get electron collection " // << LeptonInput_.label(); // return false; // } TheElectrons = *ElectronHandle; std::stable_sort( TheElectrons.begin(), TheElectrons.end(), PtSorter() ); } // Jet cuts int nJet = 0; for ( CaloJetCollection::const_iterator ajet = TheCaloJets.begin(); ajet != TheCaloJets.end(); ++ajet ) { if ( (fabs(ajet->eta()) < MaxCaloJetEta_) && (ajet->pt() > MinCaloJetPt_) ) { if (LeptonType_ == "muon" ) { for ( MuonCollection::const_iterator amuon = TheMuons.begin(); amuon != TheMuons.end(); ++amuon ) { // select good muon if ( (amuon->pt() > MinLeptonPt_) && (fabs(amuon->eta()) < MaxLeptonEta_) ) { double deltar = ROOT::Math::VectorUtil::DeltaR(ajet->p4().Vect(), amuon->momentum() ); TVector3 jetvec(ajet->p4().Vect().X(), ajet->p4().Vect().Y(), ajet->p4().Vect().Z() ); TVector3 muvec( amuon->momentum().X(), amuon->momentum().Y(), amuon->momentum().Z() ); jetvec += muvec; double ptrel = muvec.Perp(jetvec); if ( ( deltar < MaxDeltaR_ ) && (ptrel > MinPtRel_) ) nJet++; if ( nJet >= MinNLeptonJet_ ) break; } } } if (LeptonType_ == "electron") { for ( GsfElectronCollection::const_iterator anelectron = TheElectrons.begin(); anelectron != TheElectrons.end(); anelectron++ ) { if ( (anelectron->pt() > MinLeptonPt_) && (fabs(anelectron->eta()) < MaxLeptonEta_ ) ) { double deltar = ROOT::Math::VectorUtil::DeltaR(ajet->p4().Vect(), anelectron->momentum() ); TVector3 jetvec(ajet->p4().Vect().X(), ajet->p4().Vect().Y(), ajet->p4().Vect().Z() ); TVector3 evec( anelectron->momentum().X(), anelectron->momentum().Y(), anelectron->momentum().Z() ); jetvec += evec; double ptrel = evec.Perp(jetvec); if ( ( deltar < MaxDeltaR_ ) && (ptrel > MinPtRel_) ) nJet++; if ( nJet >= MinNLeptonJet_ ) break; } } } }// close jet selection } if ( nJet < MinNLeptonJet_ ) return false; nAccepted_++; return true; }
Definition at line 31 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
edm::InputTag BTagSkimLeptonJet::LeptonInput_ [private] |
Definition at line 36 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
std::string BTagSkimLeptonJet::LeptonType_ [private] |
Definition at line 35 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), endJob(), and filter().
double BTagSkimLeptonJet::MaxCaloJetEta_ [private] |
Definition at line 33 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
double BTagSkimLeptonJet::MaxDeltaR_ [private] |
Definition at line 39 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
double BTagSkimLeptonJet::MaxLeptonEta_ [private] |
Definition at line 38 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
double BTagSkimLeptonJet::MinCaloJetPt_ [private] |
Definition at line 32 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
double BTagSkimLeptonJet::MinLeptonPt_ [private] |
Definition at line 37 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
int BTagSkimLeptonJet::MinNLeptonJet_ [private] |
Definition at line 34 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
double BTagSkimLeptonJet::MinPtRel_ [private] |
Definition at line 40 of file BTagSkimLeptonJet.h.
Referenced by BTagSkimLeptonJet(), and filter().
unsigned int BTagSkimLeptonJet::nAccepted_ [private] |
Definition at line 43 of file BTagSkimLeptonJet.h.
unsigned int BTagSkimLeptonJet::nEvents_ [private] |
Definition at line 42 of file BTagSkimLeptonJet.h.