File indexing completed on 2023-03-17 10:40:54
0001
0002
0003
0004 #include "FWCore/Utilities/interface/EDMException.h"
0005 #include "CommonTools/CandUtils/interface/pdgIdUtils.h"
0006 #include "AnalysisDataFormats/TopObjects/interface/TtGenEvent.h"
0007 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0008
0009
0010 TtGenEvent::TtGenEvent(reco::GenParticleRefProd& decaySubset, reco::GenParticleRefProd& initSubset) {
0011 parts_ = decaySubset;
0012 initPartons_ = initSubset;
0013 if (top() && topBar())
0014 topPair_ = math::XYZTLorentzVector(top()->p4() + topBar()->p4());
0015 }
0016
0017 bool TtGenEvent::fromGluonFusion() const {
0018 const reco::GenParticleCollection& initPartsColl = *initPartons_;
0019 if (initPartsColl.size() == 2)
0020 if (initPartsColl[0].pdgId() == 21 && initPartsColl[1].pdgId() == 21)
0021 return true;
0022 return false;
0023 }
0024
0025 bool TtGenEvent::fromQuarkAnnihilation() const {
0026 const reco::GenParticleCollection& initPartsColl = *initPartons_;
0027 if (initPartsColl.size() == 2)
0028 if (std::abs(initPartsColl[0].pdgId()) < TopDecayID::tID && initPartsColl[0].pdgId() == -initPartsColl[1].pdgId())
0029 return true;
0030 return false;
0031 }
0032
0033 WDecay::LeptonType TtGenEvent::semiLeptonicChannel() const {
0034 WDecay::LeptonType type = WDecay::kNone;
0035 if (isSemiLeptonic() && singleLepton()) {
0036 if (std::abs(singleLepton()->pdgId()) == TopDecayID::elecID)
0037 type = WDecay::kElec;
0038 if (std::abs(singleLepton()->pdgId()) == TopDecayID::muonID)
0039 type = WDecay::kMuon;
0040 if (std::abs(singleLepton()->pdgId()) == TopDecayID::tauID)
0041 type = WDecay::kTau;
0042 }
0043 return type;
0044 }
0045
0046 std::pair<WDecay::LeptonType, WDecay::LeptonType> TtGenEvent::fullLeptonicChannel() const {
0047 WDecay::LeptonType typeA = WDecay::kNone, typeB = WDecay::kNone;
0048 if (isFullLeptonic()) {
0049 if (lepton()) {
0050 if (std::abs(lepton()->pdgId()) == TopDecayID::elecID)
0051 typeA = WDecay::kElec;
0052 if (std::abs(lepton()->pdgId()) == TopDecayID::muonID)
0053 typeA = WDecay::kMuon;
0054 if (std::abs(lepton()->pdgId()) == TopDecayID::tauID)
0055 typeA = WDecay::kTau;
0056 }
0057 if (leptonBar()) {
0058 if (std::abs(leptonBar()->pdgId()) == TopDecayID::elecID)
0059 typeB = WDecay::kElec;
0060 if (std::abs(leptonBar()->pdgId()) == TopDecayID::muonID)
0061 typeB = WDecay::kMuon;
0062 if (std::abs(leptonBar()->pdgId()) == TopDecayID::tauID)
0063 typeB = WDecay::kTau;
0064 }
0065 }
0066 return (std::pair<WDecay::LeptonType, WDecay::LeptonType>(typeA, typeB));
0067 }
0068
0069 const reco::GenParticle* TtGenEvent::lepton(bool excludeTauLeptons) const {
0070 const reco::GenParticle* cand = nullptr;
0071 const reco::GenParticleCollection& partsColl = *parts_;
0072 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0073 if (reco::isLepton(partsColl[i]) && partsColl[i].mother() &&
0074 std::abs(partsColl[i].mother()->pdgId()) == TopDecayID::WID) {
0075 if (reco::flavour(partsColl[i]) > 0) {
0076 cand = &partsColl[i];
0077 }
0078 }
0079 }
0080 return cand;
0081 }
0082
0083 const reco::GenParticle* TtGenEvent::leptonBar(bool excludeTauLeptons) const {
0084 const reco::GenParticle* cand = nullptr;
0085 const reco::GenParticleCollection& partsColl = *parts_;
0086 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0087 if (reco::isLepton(partsColl[i]) && partsColl[i].mother() &&
0088 std::abs(partsColl[i].mother()->pdgId()) == TopDecayID::WID) {
0089 if (reco::flavour(partsColl[i]) < 0) {
0090 cand = &partsColl[i];
0091 }
0092 }
0093 }
0094 return cand;
0095 }
0096
0097 const reco::GenParticle* TtGenEvent::singleLepton(bool excludeTauLeptons) const {
0098 const reco::GenParticle* cand = nullptr;
0099 if (isSemiLeptonic(excludeTauLeptons)) {
0100 const reco::GenParticleCollection& partsColl = *parts_;
0101 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0102 if (reco::isLepton(partsColl[i]) && partsColl[i].mother() &&
0103 std::abs(partsColl[i].mother()->pdgId()) == TopDecayID::WID) {
0104 cand = &partsColl[i];
0105 }
0106 }
0107 }
0108 return cand;
0109 }
0110
0111 const reco::GenParticle* TtGenEvent::neutrino(bool excludeTauLeptons) const {
0112 const reco::GenParticle* cand = nullptr;
0113 const reco::GenParticleCollection& partsColl = *parts_;
0114 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0115 if (reco::isNeutrino(partsColl[i]) && partsColl[i].mother() &&
0116 std::abs(partsColl[i].mother()->pdgId()) == TopDecayID::WID) {
0117 if (reco::flavour(partsColl[i]) > 0) {
0118 cand = &partsColl[i];
0119 }
0120 }
0121 }
0122 return cand;
0123 }
0124
0125 const reco::GenParticle* TtGenEvent::neutrinoBar(bool excludeTauLeptons) const {
0126 const reco::GenParticle* cand = nullptr;
0127 const reco::GenParticleCollection& partsColl = *parts_;
0128 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0129 if (reco::isNeutrino(partsColl[i]) && partsColl[i].mother() &&
0130 std::abs(partsColl[i].mother()->pdgId()) == TopDecayID::WID) {
0131 if (reco::flavour(partsColl[i]) < 0) {
0132 cand = &partsColl[i];
0133 }
0134 }
0135 }
0136 return cand;
0137 }
0138
0139 const reco::GenParticle* TtGenEvent::singleNeutrino(bool excludeTauLeptons) const {
0140 const reco::GenParticle* cand = nullptr;
0141 if (isSemiLeptonic(excludeTauLeptons)) {
0142 const reco::GenParticleCollection& partsColl = *parts_;
0143 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0144 if (reco::isNeutrino(partsColl[i]) && partsColl[i].mother() &&
0145 std::abs(partsColl[i].mother()->pdgId()) == TopDecayID::WID) {
0146 cand = &partsColl[i];
0147 }
0148 }
0149 }
0150 return cand;
0151 }
0152
0153 const reco::GenParticle* TtGenEvent::hadronicDecayQuark(bool invertFlavor) const {
0154 const reco::GenParticle* cand = nullptr;
0155
0156
0157
0158
0159 if (singleLepton(false)) {
0160 for (reco::GenParticleCollection::const_iterator w = parts_->begin(); w != parts_->end(); ++w) {
0161 if (std::abs(w->pdgId()) == TopDecayID::WID) {
0162
0163 for (reco::GenParticle::const_iterator wd = w->begin(); wd != w->end(); ++wd) {
0164
0165 if (std::abs(wd->pdgId()) < TopDecayID::tID) {
0166 if (invertFlavor ? reco::flavour(*wd) < 0 : reco::flavour(*wd) > 0) {
0167 cand = dynamic_cast<const reco::GenParticle*>(&(*wd));
0168 if (cand == nullptr) {
0169 throw edm::Exception(edm::errors::InvalidReference, "Not a GenParticle");
0170 }
0171 break;
0172 }
0173 }
0174 }
0175 }
0176 }
0177 }
0178 return cand;
0179 }
0180
0181 const reco::GenParticle* TtGenEvent::hadronicDecayB(bool excludeTauLeptons) const {
0182 const reco::GenParticle* cand = nullptr;
0183 if (singleLepton(excludeTauLeptons)) {
0184 const reco::GenParticleCollection& partsColl = *parts_;
0185 const reco::GenParticle& singleLep = *singleLepton(excludeTauLeptons);
0186 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0187 if (std::abs(partsColl[i].pdgId()) == TopDecayID::bID &&
0188 reco::flavour(singleLep) == reco::flavour(partsColl[i])) {
0189 cand = &partsColl[i];
0190 }
0191 }
0192 }
0193 return cand;
0194 }
0195
0196 const reco::GenParticle* TtGenEvent::hadronicDecayW(bool excludeTauLeptons) const {
0197 const reco::GenParticle* cand = nullptr;
0198 if (singleLepton(excludeTauLeptons)) {
0199 const reco::GenParticleCollection& partsColl = *parts_;
0200 const reco::GenParticle& singleLep = *singleLepton(excludeTauLeptons);
0201 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0202 if (std::abs(partsColl[i].pdgId()) == TopDecayID::WID &&
0203 reco::flavour(singleLep) != -reco::flavour(partsColl[i])) {
0204
0205 cand = &partsColl[i];
0206 }
0207 }
0208 }
0209 return cand;
0210 }
0211
0212 const reco::GenParticle* TtGenEvent::hadronicDecayTop(bool excludeTauLeptons) const {
0213 const reco::GenParticle* cand = nullptr;
0214 if (singleLepton(excludeTauLeptons)) {
0215 const reco::GenParticleCollection& partsColl = *parts_;
0216 const reco::GenParticle& singleLep = *singleLepton(excludeTauLeptons);
0217 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0218 if (std::abs(partsColl[i].pdgId()) == TopDecayID::tID &&
0219 reco::flavour(singleLep) == reco::flavour(partsColl[i])) {
0220 cand = &partsColl[i];
0221 }
0222 }
0223 }
0224 return cand;
0225 }
0226
0227 const reco::GenParticle* TtGenEvent::leptonicDecayB(bool excludeTauLeptons) const {
0228 const reco::GenParticle* cand = nullptr;
0229 if (singleLepton(excludeTauLeptons)) {
0230 const reco::GenParticleCollection& partsColl = *parts_;
0231 const reco::GenParticle& singleLep = *singleLepton(excludeTauLeptons);
0232 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0233 if (std::abs(partsColl[i].pdgId()) == TopDecayID::bID &&
0234 reco::flavour(singleLep) != reco::flavour(partsColl[i])) {
0235 cand = &partsColl[i];
0236 }
0237 }
0238 }
0239 return cand;
0240 }
0241
0242 const reco::GenParticle* TtGenEvent::leptonicDecayW(bool excludeTauLeptons) const {
0243 const reco::GenParticle* cand = nullptr;
0244 if (singleLepton(excludeTauLeptons)) {
0245 const reco::GenParticleCollection& partsColl = *parts_;
0246 const reco::GenParticle& singleLep = *singleLepton(excludeTauLeptons);
0247 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0248 if (std::abs(partsColl[i].pdgId()) == TopDecayID::WID &&
0249 reco::flavour(singleLep) == -reco::flavour(partsColl[i])) {
0250
0251 cand = &partsColl[i];
0252 }
0253 }
0254 }
0255 return cand;
0256 }
0257
0258 const reco::GenParticle* TtGenEvent::leptonicDecayTop(bool excludeTauLeptons) const {
0259 const reco::GenParticle* cand = nullptr;
0260 if (singleLepton(excludeTauLeptons)) {
0261 const reco::GenParticleCollection& partsColl = *parts_;
0262 const reco::GenParticle& singleLep = *singleLepton(excludeTauLeptons);
0263 for (unsigned int i = 0; i < partsColl.size(); ++i) {
0264 if (std::abs(partsColl[i].pdgId()) == TopDecayID::tID &&
0265 reco::flavour(singleLep) != reco::flavour(partsColl[i])) {
0266 cand = &partsColl[i];
0267 }
0268 }
0269 }
0270 return cand;
0271 }
0272
0273 std::vector<const reco::GenParticle*> TtGenEvent::leptonicDecayTopRadiation(bool excludeTauLeptons) const {
0274 if (leptonicDecayTop(excludeTauLeptons)) {
0275 return (leptonicDecayTop(excludeTauLeptons)->pdgId() > 0 ? radiatedGluons(TopDecayID::tID)
0276 : radiatedGluons(-TopDecayID::tID));
0277 }
0278 std::vector<const reco::GenParticle*> rad;
0279 return (rad);
0280 }
0281
0282 std::vector<const reco::GenParticle*> TtGenEvent::hadronicDecayTopRadiation(bool excludeTauLeptons) const {
0283 if (hadronicDecayTop(excludeTauLeptons)) {
0284 return (hadronicDecayTop(excludeTauLeptons)->pdgId() > 0 ? radiatedGluons(TopDecayID::tID)
0285 : radiatedGluons(-TopDecayID::tID));
0286 }
0287 std::vector<const reco::GenParticle*> rad;
0288 return (rad);
0289 }