#include <TopQuarkAnalysis/TopTools/plugins/TtJetPartonMatch.h>
Public Member Functions | |
virtual void | produce (edm::Event &, const edm::EventSetup &) |
write jet parton match objects into the event | |
TtJetPartonMatch (const edm::ParameterSet &) | |
default conructor | |
~TtJetPartonMatch () | |
default destructor | |
Private Types | |
typedef std::vector< pat::Jet > | TopJetCollection |
typedef for simplified reading | |
Private Attributes | |
int | algorithm_ |
choice of algorithm (defined in JetPartonMatching) | |
int | event_ |
event counter for internal use with the verbosity level | |
edm::InputTag | jets_ |
jet collection input | |
double | maxDist_ |
threshold for outliers in the case that useMaxDist_ =true | |
int | maxNComb_ |
maximal number of combinationws for which the matching should be stored | |
int | maxNJets_ |
maximal number of jets to be considered for the matching | |
bool | useDeltaR_ |
switch to choose between deltaR/deltaTheta matching | |
bool | useMaxDist_ |
switch to choose whether an outlier rejection should be applied or not | |
int | verbosity_ |
verbolity level |
Definition at line 42 of file TtJetPartonMatch.h.
typedef std::vector<pat::Jet> TtJetPartonMatch< C >::TopJetCollection [private] |
TtJetPartonMatch< C >::TtJetPartonMatch | ( | const edm::ParameterSet & | cfg | ) | [inline, explicit] |
default conructor
Definition at line 86 of file TtJetPartonMatch.h.
00086 : event_(0), 00087 jets_ (cfg.getParameter<edm::InputTag>("jets" )), 00088 maxNJets_ (cfg.getParameter<int> ("maxNJets" )), 00089 maxNComb_ (cfg.getParameter<int> ("maxNComb" )), 00090 algorithm_ (cfg.getParameter<int> ("algorithm" )), 00091 useDeltaR_ (cfg.getParameter<bool> ("useDeltaR" )), 00092 useMaxDist_(cfg.getParameter<bool> ("useMaxDist")), 00093 maxDist_ (cfg.getParameter<double> ("maxDist" )), 00094 verbosity_ (cfg.getParameter<int> ("verbosity" )) 00095 { 00096 // produces a vector of jet/lepton indices in the order of 00097 // * TtSemiLepEventPartons 00098 // * TtFullHadEventPartons 00099 // * TtFullLepEventPartons 00100 // and vectors of the corresponding quality parameters 00101 produces< std::vector<std::vector<int> > >(); 00102 produces< std::vector<double> >("SumPt"); 00103 produces< std::vector<double> >("SumDR"); 00104 }
TtJetPartonMatch< C >::~TtJetPartonMatch | ( | ) | [inline] |
void TtJetPartonMatch< C >::produce | ( | edm::Event & | evt, | |
const edm::EventSetup & | setup | |||
) | [inline, virtual] |
write jet parton match objects into the event
Implements edm::EDProducer.
Definition at line 113 of file TtJetPartonMatch.h.
References TtJetPartonMatch< C >::algorithm_, funct::C, TtJetPartonMatch< C >::event_, TtGenEvtProducer_cfi::genEvt, edm::Event::getByLabel(), JetPartonMatching::getMatchesForPartons(), JetPartonMatching::getNumberOfAvailableCombinations(), JetPartonMatching::getSumDeltaPt(), JetPartonMatching::getSumDeltaR(), jetMatch_cfi::jetPartonMatch, pfTauBenchmarkGeneric_cfi::jets, TtJetPartonMatch< C >::jets_, edm::match(), TtJetPartonMatch< C >::maxDist_, TtJetPartonMatch< C >::maxNComb_, TtJetPartonMatch< C >::maxNJets_, JetPartonMatching::print(), edm::Event::put(), TtJetPartonMatch< C >::useDeltaR_, TtJetPartonMatch< C >::useMaxDist_, and TtJetPartonMatch< C >::verbosity_.
00114 { 00115 edm::Handle<TtGenEvent> genEvt; 00116 evt.getByLabel("genEvt", genEvt); 00117 00118 edm::Handle<TopJetCollection> topJets; 00119 evt.getByLabel(jets_, topJets); 00120 00121 // fill vector of partons in the order of 00122 // * TtFullLepEventPartons 00123 // * TtSemiLepEventPartons 00124 // * TtFullHadEventPartons 00125 C parts; 00126 std::vector<const reco::Candidate*> partons = parts.vec(*genEvt); 00127 00128 // prepare vector of jets 00129 std::vector<pat::Jet> jets; 00130 for(unsigned int ij=0; ij<topJets->size(); ++ij) { 00131 // take all jets if maxNJets_ == -1; otherwise use 00132 // maxNJets_ if maxNJets_ is big enough or use same 00133 // number of jets as partons if maxNJets_ < number 00134 // of partons 00135 if(maxNJets_!=-1) { 00136 if(maxNJets_>=(int)partons.size()) { 00137 if((int)ij==maxNJets_) break; 00138 } 00139 else { 00140 if(ij==partons.size()) break; 00141 } 00142 } 00143 jets.push_back( (*topJets)[ij] ); 00144 } 00145 00146 // do the matching with specified parameters 00147 JetPartonMatching jetPartonMatch(partons, jets, algorithm_, useMaxDist_, useDeltaR_, maxDist_); 00148 00149 // print some info for each event 00150 // if corresponding verbosity level set 00151 if(verbosity_>0 && event_<verbosity_){ 00152 ++event_; 00153 jetPartonMatch.print(); 00154 } 00155 00156 // write 00157 // * parton match 00158 // * sumPt 00159 // * sumDR 00160 // to the event 00161 std::auto_ptr< std::vector<std::vector<int> > > match(new std::vector<std::vector<int> >); 00162 std::auto_ptr< std::vector<double> > sumPt(new std::vector<double>); 00163 std::auto_ptr< std::vector<double> > sumDR(new std::vector<double>); 00164 00165 for(unsigned int ic=0; ic<jetPartonMatch.getNumberOfAvailableCombinations(); ++ic) { 00166 if((int)ic>=maxNComb_ && maxNComb_>=0) break; 00167 match->push_back( jetPartonMatch.getMatchesForPartons(ic) ); 00168 sumPt->push_back( jetPartonMatch.getSumDeltaPt (ic) ); 00169 sumDR->push_back( jetPartonMatch.getSumDeltaR (ic) ); 00170 } 00171 evt.put(match); 00172 evt.put(sumPt, "SumPt"); 00173 evt.put(sumDR, "SumDR"); 00174 }
int TtJetPartonMatch< C >::algorithm_ [private] |
choice of algorithm (defined in JetPartonMatching)
Definition at line 72 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
int TtJetPartonMatch< C >::event_ [private] |
event counter for internal use with the verbosity level
Definition at line 62 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
edm::InputTag TtJetPartonMatch< C >::jets_ [private] |
jet collection input
Definition at line 64 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
double TtJetPartonMatch< C >::maxDist_ [private] |
threshold for outliers in the case that useMaxDist_ =true
Definition at line 80 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
int TtJetPartonMatch< C >::maxNComb_ [private] |
maximal number of combinationws for which the matching should be stored
Definition at line 70 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
int TtJetPartonMatch< C >::maxNJets_ [private] |
maximal number of jets to be considered for the matching
Definition at line 67 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
bool TtJetPartonMatch< C >::useDeltaR_ [private] |
switch to choose between deltaR/deltaTheta matching
Definition at line 74 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
bool TtJetPartonMatch< C >::useMaxDist_ [private] |
switch to choose whether an outlier rejection should be applied or not
Definition at line 77 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().
int TtJetPartonMatch< C >::verbosity_ [private] |
verbolity level
Definition at line 82 of file TtJetPartonMatch.h.
Referenced by TtJetPartonMatch< C >::produce().