File indexing completed on 2024-04-06 12:31:20
0001 #include "TopQuarkAnalysis/TopJetCombination/interface/TtFullHadHypothesis.h"
0002 #include "AnalysisDataFormats/TopObjects/interface/TtFullHadEvtPartons.h"
0003 #include "DataFormats/Candidate/interface/LeafCandidate.h"
0004
0005 class TtFullHadHypGenMatch : public TtFullHadHypothesis {
0006 public:
0007 explicit TtFullHadHypGenMatch(const edm::ParameterSet& cfg);
0008
0009 private:
0010
0011 void buildKey() override { key_ = TtFullHadronicEvent::kGenMatch; };
0012
0013 void buildHypo(edm::Event& evt,
0014 const edm::Handle<std::vector<pat::Jet> >& jets,
0015 std::vector<int>& match,
0016 const unsigned int iComb) override;
0017
0018 protected:
0019 edm::EDGetTokenT<TtGenEvent> genEvtToken_;
0020 };
0021
0022 TtFullHadHypGenMatch::TtFullHadHypGenMatch(const edm::ParameterSet& cfg)
0023 : TtFullHadHypothesis(cfg), genEvtToken_(consumes<TtGenEvent>(edm::InputTag("genEvt"))) {}
0024
0025 void TtFullHadHypGenMatch::buildHypo(edm::Event& evt,
0026 const edm::Handle<std::vector<pat::Jet> >& jets,
0027 std::vector<int>& match,
0028 const unsigned int iComb) {
0029
0030
0031
0032 edm::Handle<TtGenEvent> genEvt;
0033 evt.getByToken(genEvtToken_, genEvt);
0034
0035
0036
0037
0038 for (unsigned idx = 0; idx < match.size(); ++idx) {
0039 if (isValid(match[idx], jets)) {
0040 switch (idx) {
0041 case TtFullHadEvtPartons::LightQ:
0042 if (std::abs(genEvt->daughterQuarkOfWPlus()->pdgId()) == 4)
0043 lightQ_ = makeCandidate(jets, match[idx], jetCorrectionLevel("cQuark"));
0044 else
0045 lightQ_ = makeCandidate(jets, match[idx], jetCorrectionLevel("udsQuark"));
0046 break;
0047 case TtFullHadEvtPartons::LightQBar:
0048 if (std::abs(genEvt->daughterQuarkBarOfWPlus()->pdgId()) == 4)
0049 lightQBar_ = makeCandidate(jets, match[idx], jetCorrectionLevel("cQuark"));
0050 else
0051 lightQBar_ = makeCandidate(jets, match[idx], jetCorrectionLevel("udsQuark"));
0052 break;
0053 case TtFullHadEvtPartons::B:
0054 b_ = makeCandidate(jets, match[idx], jetCorrectionLevel("bQuark"));
0055 break;
0056 case TtFullHadEvtPartons::LightP:
0057 if (std::abs(genEvt->daughterQuarkOfWMinus()->pdgId()) == 4)
0058 lightP_ = makeCandidate(jets, match[idx], jetCorrectionLevel("cQuark"));
0059 else
0060 lightP_ = makeCandidate(jets, match[idx], jetCorrectionLevel("udsQuark"));
0061 break;
0062 case TtFullHadEvtPartons::LightPBar:
0063 if (std::abs(genEvt->daughterQuarkBarOfWMinus()->pdgId()) == 4)
0064 lightPBar_ = makeCandidate(jets, match[idx], jetCorrectionLevel("cQuark"));
0065 else
0066 lightPBar_ = makeCandidate(jets, match[idx], jetCorrectionLevel("udsQuark"));
0067 break;
0068 case TtFullHadEvtPartons::BBar:
0069 bBar_ = makeCandidate(jets, match[idx], jetCorrectionLevel("bQuark"));
0070 break;
0071 }
0072 }
0073 }
0074 }
0075
0076 #include "FWCore/Framework/interface/MakerMacros.h"
0077 DEFINE_FWK_MODULE(TtFullHadHypGenMatch);