Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:21:29

0001 
0002 #include <vector>
0003 #include <numeric>
0004 #include <string>
0005 
0006 ////////////////////
0007 // FRAMEWORK HEADERS
0008 #include "FWCore/Framework/interface/global/EDProducer.h"
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 
0014 #include "DataFormats/L1TParticleFlow/interface/PFCandidate.h"
0015 #include "DataFormats/L1TParticleFlow/interface/PFJet.h"
0016 #include "L1Trigger/Phase2L1ParticleFlow/interface/corrector.h"
0017 #include "DataFormats/Math/interface/deltaR.h"
0018 
0019 // bitwise emulation headers
0020 #include "L1Trigger/Phase2L1ParticleFlow/interface/jetmet/L1SeedConePFJetEmulator.h"
0021 #include "DataFormats/L1TParticleFlow/interface/gt_datatypes.h"
0022 
0023 class L1SeedConePFJetProducer : public edm::global::EDProducer<> {
0024 public:
0025   explicit L1SeedConePFJetProducer(const edm::ParameterSet&);
0026   ~L1SeedConePFJetProducer() override;
0027   static void fillDescriptions(edm::ConfigurationDescriptions& description);
0028 
0029 private:
0030   /// ///////////////// ///
0031   /// MANDATORY METHODS ///
0032   void produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const override;
0033   /// ///////////////// ///
0034 
0035   const float coneSize;
0036   const unsigned nJets;
0037   const bool HW;
0038   const bool debug;
0039   const bool doCorrections;
0040   L1SCJetEmu emulator;
0041   edm::EDGetTokenT<std::vector<l1t::PFCandidate>> l1PFToken;
0042   l1tpf::corrector corrector;
0043 
0044   std::vector<l1t::PFJet> processEvent_SW(std::vector<edm::Ptr<l1t::PFCandidate>>& parts) const;
0045   std::vector<l1t::PFJet> processEvent_HW(std::vector<edm::Ptr<l1t::PFCandidate>>& parts) const;
0046 
0047   l1t::PFJet makeJet_SW(const std::vector<edm::Ptr<l1t::PFCandidate>>& parts) const;
0048 
0049   std::pair<std::vector<L1SCJetEmu::Particle>, std::unordered_map<const l1t::PFCandidate*, edm::Ptr<l1t::PFCandidate>>>
0050   convertEDMToHW(std::vector<edm::Ptr<l1t::PFCandidate>>& edmParticles) const;
0051 
0052   std::vector<l1t::PFJet> convertHWToEDM(
0053       std::vector<L1SCJetEmu::Jet> hwJets,
0054       std::unordered_map<const l1t::PFCandidate*, edm::Ptr<l1t::PFCandidate>> constituentMap) const;
0055 };
0056 
0057 L1SeedConePFJetProducer::L1SeedConePFJetProducer(const edm::ParameterSet& cfg)
0058     : coneSize(cfg.getParameter<double>("coneSize")),
0059       nJets(cfg.getParameter<unsigned>("nJets")),
0060       HW(cfg.getParameter<bool>("HW")),
0061       debug(cfg.getParameter<bool>("debug")),
0062       doCorrections(cfg.getParameter<bool>("doCorrections")),
0063       emulator(L1SCJetEmu(debug, coneSize, nJets)),
0064       l1PFToken(consumes<std::vector<l1t::PFCandidate>>(cfg.getParameter<edm::InputTag>("L1PFObjects"))) {
0065   produces<l1t::PFJetCollection>();
0066   if (doCorrections) {
0067     corrector = l1tpf::corrector(
0068         cfg.getParameter<std::string>("correctorFile"), cfg.getParameter<std::string>("correctorDir"), -1., debug, HW);
0069   }
0070 }
0071 
0072 void L1SeedConePFJetProducer::produce(edm::StreamID /*unused*/,
0073                                       edm::Event& iEvent,
0074                                       const edm::EventSetup& iSetup) const {
0075   std::unique_ptr<l1t::PFJetCollection> newPFJetCollection(new l1t::PFJetCollection);
0076 
0077   edm::Handle<l1t::PFCandidateCollection> l1PFCandidates;
0078   iEvent.getByToken(l1PFToken, l1PFCandidates);
0079 
0080   std::vector<edm::Ptr<l1t::PFCandidate>> particles;
0081   for (unsigned i = 0; i < (*l1PFCandidates).size(); i++) {
0082     particles.push_back(edm::Ptr<l1t::PFCandidate>(l1PFCandidates, i));
0083   }
0084 
0085   std::vector<l1t::PFJet> jets;
0086   if (HW) {
0087     jets = processEvent_HW(particles);
0088   } else {
0089     jets = processEvent_SW(particles);
0090   }
0091 
0092   std::sort(jets.begin(), jets.end(), [](l1t::PFJet i, l1t::PFJet j) { return (i.pt() > j.pt()); });
0093   newPFJetCollection->swap(jets);
0094   iEvent.put(std::move(newPFJetCollection));
0095 }
0096 
0097 /////////////
0098 // DESTRUCTOR
0099 L1SeedConePFJetProducer::~L1SeedConePFJetProducer() {}
0100 
0101 l1t::PFJet L1SeedConePFJetProducer::makeJet_SW(const std::vector<edm::Ptr<l1t::PFCandidate>>& parts) const {
0102   l1t::PFCandidate seed = *parts.at(0);
0103 
0104   auto sumpt = [](float a, const edm::Ptr<l1t::PFCandidate>& b) { return a + b->pt(); };
0105 
0106   // Sum the pt
0107   float pt = std::accumulate(parts.begin(), parts.end(), 0., sumpt);
0108 
0109   // pt weighted d eta
0110   std::vector<float> pt_deta;
0111   pt_deta.resize(parts.size());
0112   std::transform(parts.begin(), parts.end(), pt_deta.begin(), [&seed, &pt](const edm::Ptr<l1t::PFCandidate>& part) {
0113     return (part->pt() / pt) * (part->eta() - seed.eta());
0114   });
0115   // Accumulate the pt weighted etas. Init to the seed eta, start accumulating at begin()+1 to skip seed
0116   float eta = std::accumulate(pt_deta.begin() + 1, pt_deta.end(), seed.eta());
0117 
0118   // pt weighted d phi
0119   std::vector<float> pt_dphi;
0120   pt_dphi.resize(parts.size());
0121   std::transform(parts.begin(), parts.end(), pt_dphi.begin(), [&seed, &pt](const edm::Ptr<l1t::PFCandidate>& part) {
0122     return (part->pt() / pt) * reco::deltaPhi(part->phi(), seed.phi());
0123   });
0124   // Accumulate the pt weighted phis. Init to the seed phi, start accumulating at begin()+1 to skip seed
0125   float phi = std::accumulate(pt_dphi.begin() + 1, pt_dphi.end(), seed.phi());
0126 
0127   l1t::PFJet jet(pt, eta, phi);
0128   for (auto it = parts.begin(); it != parts.end(); it++) {
0129     jet.addConstituent(*it);
0130   }
0131 
0132   if (doCorrections) {
0133     jet.calibratePt(corrector.correctedPt(jet.pt(), jet.eta()));
0134   }
0135 
0136   return jet;
0137 }
0138 
0139 std::vector<l1t::PFJet> L1SeedConePFJetProducer::processEvent_SW(std::vector<edm::Ptr<l1t::PFCandidate>>& work) const {
0140   // The floating point algorithm simulation
0141   std::stable_sort(work.begin(), work.end(), [](edm::Ptr<l1t::PFCandidate> i, edm::Ptr<l1t::PFCandidate> j) {
0142     return (i->pt() > j->pt());
0143   });
0144   std::vector<l1t::PFJet> jets;
0145   jets.reserve(nJets);
0146   while (!work.empty() && jets.size() < nJets) {
0147     // Take the first (highest pt) candidate as a seed
0148     edm::Ptr<l1t::PFCandidate> seed = work.at(0);
0149     // Get the particles within a coneSize of the seed
0150     std::vector<edm::Ptr<l1t::PFCandidate>> particlesInCone;
0151     std::copy_if(
0152         work.begin(), work.end(), std::back_inserter(particlesInCone), [&](const edm::Ptr<l1t::PFCandidate>& part) {
0153           return reco::deltaR<l1t::PFCandidate, l1t::PFCandidate>(*seed, *part) <= coneSize;
0154         });
0155     jets.push_back(makeJet_SW(particlesInCone));
0156     // remove the clustered particles
0157     work.erase(std::remove_if(work.begin(),
0158                               work.end(),
0159                               [&](const edm::Ptr<l1t::PFCandidate>& part) {
0160                                 return reco::deltaR<l1t::PFCandidate, l1t::PFCandidate>(*seed, *part) <= coneSize;
0161                               }),
0162                work.end());
0163   }
0164   return jets;
0165 }
0166 
0167 std::vector<l1t::PFJet> L1SeedConePFJetProducer::processEvent_HW(std::vector<edm::Ptr<l1t::PFCandidate>>& work) const {
0168   // The fixed point emulator
0169   // Convert the EDM format to the hardware format, and call the standalone emulator
0170   std::pair<std::vector<L1SCJetEmu::Particle>, std::unordered_map<const l1t::PFCandidate*, edm::Ptr<l1t::PFCandidate>>>
0171       particles = convertEDMToHW(work);
0172   std::vector<L1SCJetEmu::Jet> jets = emulator.emulateEvent(particles.first);
0173   return convertHWToEDM(jets, particles.second);
0174 }
0175 
0176 std::pair<std::vector<L1SCJetEmu::Particle>, std::unordered_map<const l1t::PFCandidate*, edm::Ptr<l1t::PFCandidate>>>
0177 L1SeedConePFJetProducer::convertEDMToHW(std::vector<edm::Ptr<l1t::PFCandidate>>& edmParticles) const {
0178   std::vector<l1ct::PuppiObjEmu> hwParticles;
0179   std::unordered_map<const l1t::PFCandidate*, edm::Ptr<l1t::PFCandidate>> candidateMap;
0180   std::for_each(edmParticles.begin(), edmParticles.end(), [&](edm::Ptr<l1t::PFCandidate>& edmParticle) {
0181     l1ct::PuppiObjEmu particle;
0182     particle.initFromBits(edmParticle->encodedPuppi64());
0183     particle.srcCand = edmParticle.get();
0184     candidateMap.insert(std::make_pair(edmParticle.get(), edmParticle));
0185     hwParticles.push_back(particle);
0186   });
0187   return std::make_pair(hwParticles, candidateMap);
0188 }
0189 
0190 std::vector<l1t::PFJet> L1SeedConePFJetProducer::convertHWToEDM(
0191     std::vector<L1SCJetEmu::Jet> hwJets,
0192     std::unordered_map<const l1t::PFCandidate*, edm::Ptr<l1t::PFCandidate>> constituentMap) const {
0193   std::vector<l1t::PFJet> edmJets;
0194   std::for_each(hwJets.begin(), hwJets.end(), [&](L1SCJetEmu::Jet jet) {
0195     if (doCorrections) {
0196       float correctedPt = corrector.correctedPt(jet.floatPt(), jet.floatEta());
0197       jet.hwPt = correctedPt;
0198     }
0199     l1gt::Jet gtJet = jet.toGT();
0200     l1t::PFJet edmJet(l1gt::Scales::floatPt(gtJet.v3.pt),
0201                       l1gt::Scales::floatEta(gtJet.v3.eta),
0202                       l1gt::Scales::floatPhi(gtJet.v3.phi),
0203                       /*mass=*/0.,
0204                       gtJet.v3.pt.V,
0205                       gtJet.v3.eta.V,
0206                       gtJet.v3.phi.V);
0207     edmJet.setEncodedJet(l1t::PFJet::HWEncoding::CT, jet.pack());
0208     edmJet.setEncodedJet(l1t::PFJet::HWEncoding::GT, jet.toGT().pack());
0209     // get back the references to the constituents
0210     std::vector<edm::Ptr<l1t::PFCandidate>> constituents;
0211     std::for_each(jet.constituents.begin(), jet.constituents.end(), [&](auto constituent) {
0212       edmJet.addConstituent(constituentMap[constituent.srcCand]);
0213     });
0214     edmJets.push_back(edmJet);
0215   });
0216   return edmJets;
0217 }
0218 
0219 void L1SeedConePFJetProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0220   edm::ParameterSetDescription desc;
0221   desc.add<edm::InputTag>("L1PFObjects", edm::InputTag("l1tLayer1", "Puppi"));
0222   desc.add<uint32_t>("nJets", 16);
0223   desc.add<double>("coneSize", 0.4);
0224   desc.add<bool>("HW", false);
0225   desc.add<bool>("debug", false);
0226   desc.add<bool>("doCorrections", false);
0227   desc.add<std::string>("correctorFile", "");
0228   desc.add<std::string>("correctorDir", "");
0229   descriptions.addWithDefaultLabel(desc);
0230 }
0231 
0232 DEFINE_FWK_MODULE(L1SeedConePFJetProducer);