File indexing completed on 2023-03-17 10:45:24
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021 #include <vector>
0022
0023
0024 #include "FWCore/Framework/interface/Frameworkfwd.h"
0025 #include "FWCore/Framework/interface/stream/EDProducer.h"
0026
0027 #include "FWCore/Framework/interface/Event.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029
0030 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0031 #include "DataFormats/Common/interface/View.h"
0032 #include "DataFormats/Candidate/interface/Candidate.h"
0033 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0034 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0035 #include "DataFormats/Common/interface/ValueMap.h"
0036 #include "fastjet/contrib/SoftKiller.hh"
0037
0038
0039
0040
0041
0042 class SoftKillerProducer : public edm::stream::EDProducer<> {
0043 public:
0044 typedef math::XYZTLorentzVector LorentzVector;
0045 typedef std::vector<LorentzVector> LorentzVectorCollection;
0046 typedef std::vector<reco::PFCandidate> PFOutputCollection;
0047
0048 explicit SoftKillerProducer(const edm::ParameterSet&);
0049 ~SoftKillerProducer() override;
0050
0051 private:
0052 void produce(edm::Event&, const edm::EventSetup&) override;
0053
0054 edm::EDGetTokenT<reco::CandidateView> tokenPFCandidates_;
0055
0056 double Rho_EtaMax_;
0057 double rParam_;
0058 };
0059
0060
0061
0062
0063
0064
0065
0066
0067
0068
0069
0070
0071 SoftKillerProducer::SoftKillerProducer(const edm::ParameterSet& iConfig)
0072 : Rho_EtaMax_(iConfig.getParameter<double>("Rho_EtaMax")), rParam_(iConfig.getParameter<double>("rParam")) {
0073 tokenPFCandidates_ = consumes<reco::CandidateView>(iConfig.getParameter<edm::InputTag>("PFCandidates"));
0074
0075 produces<edm::ValueMap<LorentzVector> >("SoftKillerP4s");
0076 produces<PFOutputCollection>();
0077 }
0078
0079 SoftKillerProducer::~SoftKillerProducer() {}
0080
0081
0082
0083
0084
0085
0086 void SoftKillerProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0087 std::unique_ptr<PFOutputCollection> pOutput(new PFOutputCollection);
0088
0089
0090 edm::Handle<reco::CandidateView> pfCandidates;
0091 iEvent.getByToken(tokenPFCandidates_, pfCandidates);
0092
0093 std::vector<fastjet::PseudoJet> fjInputs;
0094 for (auto i = pfCandidates->begin(), ibegin = pfCandidates->begin(), iend = pfCandidates->end(); i != iend; ++i) {
0095 fjInputs.push_back(fastjet::PseudoJet(i->px(), i->py(), i->pz(), i->energy()));
0096 fjInputs.back().set_user_index(i - ibegin);
0097 }
0098
0099
0100 fastjet::contrib::SoftKiller soft_killer(Rho_EtaMax_, rParam_);
0101
0102 double pt_threshold = 0.;
0103 std::vector<fastjet::PseudoJet> soft_killed_event;
0104 soft_killer.apply(fjInputs, soft_killed_event, pt_threshold);
0105
0106 std::unique_ptr<edm::ValueMap<LorentzVector> > p4SKOut(new edm::ValueMap<LorentzVector>());
0107 LorentzVectorCollection skP4s;
0108
0109 static const reco::PFCandidate dummySinceTranslateIsNotStatic;
0110
0111
0112
0113 for (auto j = fjInputs.begin(), jend = fjInputs.end(); j != jend; ++j) {
0114 const reco::Candidate& cand = pfCandidates->at(j->user_index());
0115 auto id = dummySinceTranslateIsNotStatic.translatePdgIdToType(cand.pdgId());
0116 const reco::PFCandidate* pPF = dynamic_cast<const reco::PFCandidate*>(&cand);
0117 reco::PFCandidate pCand(pPF ? *pPF : reco::PFCandidate(cand.charge(), cand.p4(), id));
0118 auto val = j->user_index();
0119 auto skmatch = find_if(soft_killed_event.begin(), soft_killed_event.end(), [&val](fastjet::PseudoJet const& i) {
0120 return i.user_index() == val;
0121 });
0122 LorentzVector pVec;
0123 if (skmatch != soft_killed_event.end()) {
0124 pVec.SetPxPyPzE(skmatch->px(), skmatch->py(), skmatch->pz(), skmatch->E());
0125 } else {
0126 pVec.SetPxPyPzE(0., 0., 0., 0.);
0127 }
0128 pCand.setP4(pVec);
0129 skP4s.push_back(pVec);
0130 pOutput->push_back(pCand);
0131 }
0132
0133
0134 edm::ValueMap<LorentzVector>::Filler p4SKFiller(*p4SKOut);
0135 p4SKFiller.insert(pfCandidates, skP4s.begin(), skP4s.end());
0136 p4SKFiller.fill();
0137
0138 iEvent.put(std::move(p4SKOut), "SoftKillerP4s");
0139 iEvent.put(std::move(pOutput));
0140 }
0141
0142
0143 DEFINE_FWK_MODULE(SoftKillerProducer);