File indexing completed on 2024-04-06 12:31:19
0001 #ifndef TtSemiLepHypothesis_h
0002 #define TtSemiLepHypothesis_h
0003
0004 #include <memory>
0005 #include <vector>
0006
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "FWCore/Framework/interface/stream/EDProducer.h"
0009 #include "FWCore/Framework/interface/Frameworkfwd.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "DataFormats/PatCandidates/interface/Jet.h"
0013 #include "DataFormats/PatCandidates/interface/MET.h"
0014 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0015 #include "DataFormats/Candidate/interface/ShallowClonePtrCandidate.h"
0016
0017 #include "AnalysisDataFormats/TopObjects/interface/TtSemiLeptonicEvent.h"
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031 class TtSemiLepHypothesis : public edm::stream::EDProducer<> {
0032 public:
0033
0034 explicit TtSemiLepHypothesis(const edm::ParameterSet&);
0035
0036 protected:
0037
0038 void produce(edm::Event&, const edm::EventSetup&) override;
0039
0040 void resetCandidates();
0041
0042
0043 std::string jetCorrectionLevel(const std::string& quarkType);
0044
0045 template <typename C>
0046 std::unique_ptr<reco::ShallowClonePtrCandidate> makeCandidate(const edm::Handle<C>& handle, const int& idx);
0047
0048 std::unique_ptr<reco::ShallowClonePtrCandidate> makeCandidate(const edm::Handle<std::vector<pat::Jet> >& handle,
0049 const int& idx,
0050 const std::string& correctionLevel);
0051
0052 void setNeutrino(const edm::Handle<std::vector<pat::MET> >& met,
0053 const edm::Handle<edm::View<reco::RecoCandidate> >& leps,
0054 const int& idx,
0055 const int& type);
0056
0057
0058 void buildHypo(const edm::Handle<edm::View<reco::RecoCandidate> >& leps,
0059 const edm::Handle<std::vector<pat::MET> >& mets,
0060 const edm::Handle<std::vector<pat::Jet> >& jets,
0061 std::vector<int>& jetPartonAssociation);
0062 int key() const { return key_; };
0063
0064 reco::CompositeCandidate hypo();
0065
0066 bool isValid(const int& idx, const edm::Handle<std::vector<pat::Jet> >& jets) {
0067 return (0 <= idx && idx < (int)jets->size());
0068 };
0069
0070 WDecay::LeptonType leptonType(const reco::RecoCandidate* cand);
0071
0072
0073
0074
0075
0076
0077
0078 virtual void buildKey() = 0;
0079
0080 virtual void buildHypo(edm::Event& event,
0081 const edm::Handle<edm::View<reco::RecoCandidate> >& lepton,
0082 const edm::Handle<std::vector<pat::MET> >& neutrino,
0083 const edm::Handle<std::vector<pat::Jet> >& jets,
0084 std::vector<int>& jetPartonAssociation,
0085 const unsigned int iComb) = 0;
0086
0087 protected:
0088
0089
0090 bool getMatch_;
0091
0092 edm::EDGetTokenT<std::vector<pat::Jet> > jetsToken_;
0093 edm::EDGetTokenT<edm::View<reco::RecoCandidate> > lepsToken_;
0094 edm::EDGetTokenT<std::vector<pat::MET> > metsToken_;
0095 edm::EDGetTokenT<std::vector<std::vector<int> > > matchToken_;
0096 edm::EDGetTokenT<int> nJetsConsideredToken_;
0097
0098
0099 std::string jetCorrectionLevel_;
0100
0101 int key_;
0102
0103 int neutrinoSolutionType_;
0104
0105
0106 int numberOfRealNeutrinoSolutions_;
0107
0108
0109 std::unique_ptr<reco::ShallowClonePtrCandidate> lightQ_;
0110 std::unique_ptr<reco::ShallowClonePtrCandidate> lightQBar_;
0111 std::unique_ptr<reco::ShallowClonePtrCandidate> hadronicB_;
0112 std::unique_ptr<reco::ShallowClonePtrCandidate> leptonicB_;
0113 std::unique_ptr<reco::ShallowClonePtrCandidate> neutrino_;
0114 std::unique_ptr<reco::ShallowClonePtrCandidate> lepton_;
0115 };
0116
0117
0118
0119 template <typename C>
0120 std::unique_ptr<reco::ShallowClonePtrCandidate> TtSemiLepHypothesis::makeCandidate(const edm::Handle<C>& handle,
0121 const int& idx) {
0122 typedef typename C::value_type O;
0123 edm::Ptr<O> ptr = edm::Ptr<O>(handle, idx);
0124 return std::make_unique<reco::ShallowClonePtrCandidate>(ptr, ptr->charge(), ptr->p4(), ptr->vertex());
0125 }
0126
0127 #endif