File indexing completed on 2024-04-06 12:23:24
0001 #include "FWCore/Framework/interface/global/EDProducer.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0006 #include "DataFormats/Common/interface/View.h"
0007 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0008 #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
0009 #include "DataFormats/JetReco/interface/GenJet.h"
0010 #include "DataFormats/JetReco/interface/GenJetCollection.h"
0011 #include "PhysicsTools/JetMCUtils/interface/JetMCTag.h"
0012 #include "DataFormats/TauReco/interface/PFTau.h"
0013 #include "DataFormats/Math/interface/deltaR.h"
0014 #include "CommonTools/Utils/interface/PtComparator.h"
0015
0016 #include <vector>
0017 #include <iostream>
0018
0019 class GenVisTauProducer : public edm::global::EDProducer<> {
0020 public:
0021 GenVisTauProducer(const edm::ParameterSet& params)
0022 : src_(consumes<reco::GenJetCollection>(params.getParameter<edm::InputTag>("src"))),
0023 srcGenParticles_(consumes<reco::GenParticleCollection>(params.getParameter<edm::InputTag>("srcGenParticles"))),
0024 pTComparator_() {
0025 produces<reco::GenParticleCollection>();
0026 }
0027
0028 ~GenVisTauProducer() override {}
0029
0030 void produce(edm::StreamID id, edm::Event& evt, const edm::EventSetup& es) const override {
0031 edm::Handle<reco::GenJetCollection> genTauJets;
0032 evt.getByToken(src_, genTauJets);
0033
0034 edm::Handle<reco::GenParticleCollection> genParticles;
0035 evt.getByToken(srcGenParticles_, genParticles);
0036 size_t numGenParticles = genParticles->size();
0037
0038 auto genVisTaus = std::make_unique<reco::GenParticleCollection>();
0039
0040 for (const auto& genTauJet : *genTauJets) {
0041 std::string decayMode_string = JetMCTagUtils::genTauDecayMode(genTauJet);
0042
0043 if (decayMode_string == "electron" || decayMode_string == "muon")
0044 continue;
0045 int decayMode = reco::PFTau::kNull;
0046 if (decayMode_string == "oneProng0Pi0")
0047 decayMode = reco::PFTau::kOneProng0PiZero;
0048 else if (decayMode_string == "oneProng1Pi0")
0049 decayMode = reco::PFTau::kOneProng1PiZero;
0050 else if (decayMode_string == "oneProng2Pi0")
0051 decayMode = reco::PFTau::kOneProng2PiZero;
0052 else if (decayMode_string == "threeProng0Pi0")
0053 decayMode = reco::PFTau::kThreeProng0PiZero;
0054 else if (decayMode_string == "threeProng1Pi0")
0055 decayMode = reco::PFTau::kThreeProng1PiZero;
0056 else
0057 decayMode = reco::PFTau::kRareDecayMode;
0058
0059 int pdgId = (genTauJet.charge() > 0) ? -15 : +15;
0060
0061
0062 reco::GenParticle genVisTau(genTauJet.charge(), genTauJet.p4(), genTauJet.vertex(), pdgId, decayMode, true);
0063
0064
0065 for (size_t idxGenParticle = 0; idxGenParticle < numGenParticles; ++idxGenParticle) {
0066 const reco::GenParticle& genTau = (*genParticles)[idxGenParticle];
0067 if (abs(genTau.pdgId()) == 15 && genTau.status() == 2) {
0068 reco::Candidate::LorentzVector daughterVisP4;
0069 for (const reco::GenParticleRef& daughter : genTau.daughterRefVector()) {
0070 int abs_pdgId = abs(daughter->pdgId());
0071
0072 if (abs_pdgId == 12 || abs_pdgId == 14 || abs_pdgId == 16)
0073 continue;
0074 daughterVisP4 += daughter->p4();
0075 }
0076 double dR2 = deltaR2(daughterVisP4, genVisTau);
0077 if (dR2 < 1.e-4) {
0078 genVisTau.addMother(reco::GenParticleRef(genParticles, idxGenParticle));
0079 break;
0080 }
0081 }
0082 }
0083
0084 genVisTaus->push_back(genVisTau);
0085 }
0086
0087 std::sort(genVisTaus->begin(), genVisTaus->end(), pTComparator_);
0088 evt.put(std::move(genVisTaus));
0089 }
0090
0091 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0092 edm::ParameterSetDescription desc;
0093 desc.add<edm::InputTag>("src")->setComment("collection of visible gen taus (as reco::GenJetCollection)");
0094 desc.add<edm::InputTag>("srcGenParticles")->setComment("collections of gen particles");
0095 descriptions.add("genVisTaus", desc);
0096 }
0097
0098 private:
0099 const edm::EDGetTokenT<reco::GenJetCollection> src_;
0100 const edm::EDGetTokenT<reco::GenParticleCollection> srcGenParticles_;
0101 const GreaterByPt<reco::GenParticle> pTComparator_;
0102 };
0103
0104 #include "FWCore/Framework/interface/MakerMacros.h"
0105 DEFINE_FWK_MODULE(GenVisTauProducer);