Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2022-06-10 01:53:53 

0001 // system include files
0002 #include <memory>
0003 #include <algorithm>
0004 #include <fstream>
0005 #include <cstdio>
0006 
0007 // user include files
0008 #include "FWCore/Framework/interface/Frameworkfwd.h"
0009 #include "FWCore/Framework/interface/stream/EDProducer.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/ParameterSet/interface/FileInPath.h"
0013 
0014 #include "DataFormats/Common/interface/View.h"
0015 #include "DataFormats/Common/interface/RefToPtr.h"
0016 #include "DataFormats/L1TParticleFlow/interface/PFCandidate.h"
0017 #include "DataFormats/L1Trigger/interface/Vertex.h"
0018 #include "DataFormats/L1Trigger/interface/VertexWord.h"
0019 
0020 #include "DataFormats/Math/interface/deltaR.h"
0021 
0022 #include "L1Trigger/Phase2L1ParticleFlow/interface/l1-converters/tkinput_ref.h"
0023 #include "L1Trigger/Phase2L1ParticleFlow/interface/l1-converters/muonGmtToL1ct_ref.h"
0024 #include "L1Trigger/Phase2L1ParticleFlow/interface/regionizer/regionizer_base_ref.h"
0025 #include "L1Trigger/Phase2L1ParticleFlow/interface/regionizer/multififo_regionizer_ref.h"
0026 #include "L1Trigger/Phase2L1ParticleFlow/interface/regionizer/tdr_regionizer_ref.h"
0027 #include "L1Trigger/Phase2L1ParticleFlow/interface/pf/pfalgo2hgc_ref.h"
0028 #include "L1Trigger/Phase2L1ParticleFlow/interface/pf/pfalgo3_ref.h"
0029 #include "L1Trigger/Phase2L1ParticleFlow/interface/pf/pfalgo_dummy_ref.h"
0030 #include "L1Trigger/Phase2L1ParticleFlow/interface/puppi/linpuppi_ref.h"
0031 #include "L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegalgo_ref.h"
0032 #include "L1Trigger/Phase2L1ParticleFlow/interface/pf/pfalgo_common_ref.h"
0033 #include "L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_ref.h"
0034 
0035 #include "DataFormats/L1TCorrelator/interface/TkElectron.h"
0036 #include "DataFormats/L1TCorrelator/interface/TkElectronFwd.h"
0037 #include "DataFormats/L1Trigger/interface/EGamma.h"
0038 #include "DataFormats/L1TCorrelator/interface/TkEm.h"
0039 #include "DataFormats/L1TCorrelator/interface/TkEmFwd.h"
0040 
0041 //--------------------------------------------------------------------------------------------------
0042 class L1TCorrelatorLayer1Producer : public edm::stream::EDProducer<> {
0043 public:
0044   explicit L1TCorrelatorLayer1Producer(const edm::ParameterSet &);
0045   ~L1TCorrelatorLayer1Producer() override;
0046 
0047 private:
0048   edm::ParameterSet config_;
0049 
0050   bool hasTracks_;
0051   edm::EDGetTokenT<l1t::PFTrackCollection> tkCands_;
0052   float trkPt_;
0053   bool emuTkVtx_;
0054   edm::EDGetTokenT<std::vector<l1t::Vertex>> extTkVtx_;
0055   edm::EDGetTokenT<std::vector<l1t::VertexWord>> tkVtxEmu_;
0056   int nVtx_;
0057 
0058   edm::EDGetTokenT<l1t::SAMuonCollection> muCands_;  // standalone muons
0059 
0060   std::vector<edm::EDGetTokenT<l1t::PFClusterCollection>> emCands_;
0061   std::vector<edm::EDGetTokenT<l1t::PFClusterCollection>> hadCands_;
0062 
0063   float emPtCut_, hadPtCut_;
0064 
0065   l1ct::Event event_;
0066   std::unique_ptr<l1ct::TrackInputEmulator> trackInput_;
0067   std::unique_ptr<l1ct::GMTMuonDecoderEmulator> muonInput_;
0068   std::unique_ptr<l1ct::RegionizerEmulator> regionizer_;
0069   std::unique_ptr<l1ct::PFAlgoEmulatorBase> l1pfalgo_;
0070   std::unique_ptr<l1ct::LinPuppiEmulator> l1pualgo_;
0071   std::unique_ptr<l1ct::PFTkEGAlgoEmulator> l1tkegalgo_;
0072   std::unique_ptr<l1ct::PFTkEGSorterEmulator> l1tkegsorter_;
0073 
0074   bool writeEgSta_;
0075   // Region dump
0076   const std::string regionDumpName_;
0077   bool writeRawHgcalCluster_;
0078   std::fstream fRegionDump_;
0079 
0080   // region of interest debugging
0081   float debugEta_, debugPhi_, debugR_;
0082 
0083   // these are used to link items back
0084   std::unordered_map<const l1t::PFCluster *, l1t::PFClusterRef> clusterRefMap_;
0085   std::unordered_map<const l1t::PFTrack *, l1t::PFTrackRef> trackRefMap_;
0086   std::unordered_map<const l1t::SAMuon *, l1t::PFCandidate::MuonRef> muonRefMap_;
0087 
0088   // main methods
0089   void beginStream(edm::StreamID) override;
0090   void produce(edm::Event &, const edm::EventSetup &) override;
0091   void addUInt(unsigned int value, std::string iLabel, edm::Event &iEvent);
0092 
0093   void initSectorsAndRegions(const edm::ParameterSet &iConfig);
0094   void initEvent(const edm::Event &e);
0095   // add object, tracking references
0096   void addTrack(const l1t::PFTrack &t, l1t::PFTrackRef ref);
0097   void addMuon(const l1t::SAMuon &t, l1t::PFCandidate::MuonRef ref);
0098   void addHadCalo(const l1t::PFCluster &t, l1t::PFClusterRef ref);
0099   void addEmCalo(const l1t::PFCluster &t, l1t::PFClusterRef ref);
0100   // add objects in already-decoded format
0101   void addDecodedTrack(l1ct::DetectorSector<l1ct::TkObjEmu> &sec, const l1t::PFTrack &t);
0102   void addDecodedMuon(l1ct::DetectorSector<l1ct::MuObjEmu> &sec, const l1t::SAMuon &t);
0103   void addDecodedHadCalo(l1ct::DetectorSector<l1ct::HadCaloObjEmu> &sec, const l1t::PFCluster &t);
0104   void addDecodedEmCalo(l1ct::DetectorSector<l1ct::EmCaloObjEmu> &sec, const l1t::PFCluster &t);
0105 
0106   void addRawHgcalCluster(l1ct::DetectorSector<ap_uint<256>> &sec, const l1t::PFCluster &c);
0107 
0108   // fetching outputs
0109   std::unique_ptr<l1t::PFCandidateCollection> fetchHadCalo() const;
0110   std::unique_ptr<l1t::PFCandidateCollection> fetchEmCalo() const;
0111   std::unique_ptr<l1t::PFCandidateCollection> fetchTracks() const;
0112   std::unique_ptr<l1t::PFCandidateCollection> fetchPF() const;
0113   void putPuppi(edm::Event &iEvent) const;
0114 
0115   void putEgStaObjects(edm::Event &iEvent,
0116                        const std::string &egLablel,
0117                        std::vector<edm::Ref<BXVector<l1t::EGamma>>> &egsta_refs);
0118   void putEgObjects(edm::Event &iEvent,
0119                     const bool writeEgSta,
0120                     const std::vector<edm::Ref<BXVector<l1t::EGamma>>> &egsta_refs,
0121                     const std::string &tkEmLabel,
0122                     const std::string &tkEmPerBoardLabel,
0123                     const std::string &tkEleLabel,
0124                     const std::string &tkElePerBoardLabel) const;
0125 
0126   template <typename T>
0127   void setRefs_(l1t::PFCandidate &pf, const T &p) const;
0128 
0129   void doVertexings(std::vector<float> &pvdz) const;
0130   // for multiplicities
0131   enum InputType { caloType = 0, emcaloType = 1, trackType = 2, l1muType = 3 };
0132   static constexpr const char *inputTypeName[l1muType + 1] = {"Calo", "EmCalo", "TK", "Mu"};
0133   std::unique_ptr<std::vector<unsigned>> vecSecInput(InputType i) const;
0134   std::unique_ptr<std::vector<unsigned>> vecRegInput(InputType i) const;
0135   typedef l1ct::OutputRegion::ObjType OutputType;
0136   std::unique_ptr<std::vector<unsigned>> vecOutput(OutputType i, bool usePuppi) const;
0137   std::pair<unsigned int, unsigned int> totAndMax(const std::vector<unsigned> &perRegion) const;
0138 };
0139 
0140 //
0141 // constructors and destructor
0142 //
0143 L1TCorrelatorLayer1Producer::L1TCorrelatorLayer1Producer(const edm::ParameterSet &iConfig)
0144     : config_(iConfig),
0145       hasTracks_(!iConfig.getParameter<edm::InputTag>("tracks").label().empty()),
0146       tkCands_(hasTracks_ ? consumes<l1t::PFTrackCollection>(iConfig.getParameter<edm::InputTag>("tracks"))
0147                           : edm::EDGetTokenT<l1t::PFTrackCollection>()),
0148       trkPt_(iConfig.getParameter<double>("trkPtCut")),
0149       muCands_(consumes<l1t::SAMuonCollection>(iConfig.getParameter<edm::InputTag>("muons"))),
0150       emPtCut_(iConfig.getParameter<double>("emPtCut")),
0151       hadPtCut_(iConfig.getParameter<double>("hadPtCut")),
0152       regionizer_(nullptr),
0153       l1pfalgo_(nullptr),
0154       l1pualgo_(nullptr),
0155       l1tkegalgo_(nullptr),
0156       l1tkegsorter_(nullptr),
0157       regionDumpName_(iConfig.getUntrackedParameter<std::string>("dumpFileName", "")),
0158       writeRawHgcalCluster_(iConfig.getUntrackedParameter<bool>("writeRawHgcalCluster", false)),
0159       debugEta_(iConfig.getUntrackedParameter<double>("debugEta", 0)),
0160       debugPhi_(iConfig.getUntrackedParameter<double>("debugPhi", 0)),
0161       debugR_(iConfig.getUntrackedParameter<double>("debugR", -1)) {
0162   produces<l1t::PFCandidateCollection>("PF");
0163   produces<l1t::PFCandidateCollection>("Puppi");
0164   produces<l1t::PFCandidateRegionalOutput>("PuppiRegional");
0165 
0166   produces<l1t::PFCandidateCollection>("EmCalo");
0167   produces<l1t::PFCandidateCollection>("Calo");
0168   produces<l1t::PFCandidateCollection>("TK");
0169 #if 0  // LATER
0170   produces<l1t::PFCandidateCollection>("TKVtx");
0171 #endif
0172 
0173   for (const auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("emClusters")) {
0174     emCands_.push_back(consumes<l1t::PFClusterCollection>(tag));
0175   }
0176   for (const auto &tag : iConfig.getParameter<std::vector<edm::InputTag>>("hadClusters")) {
0177     hadCands_.push_back(consumes<l1t::PFClusterCollection>(tag));
0178   }
0179 
0180   if (hasTracks_) {
0181     const std::string &tkInAlgo = iConfig.getParameter<std::string>("trackInputConversionAlgo");
0182     if (tkInAlgo == "Emulator") {
0183       trackInput_ = std::make_unique<l1ct::TrackInputEmulator>(
0184           iConfig.getParameter<edm::ParameterSet>("trackInputConversionParameters"));
0185     } else if (tkInAlgo != "Ideal")
0186       throw cms::Exception("Configuration", "Unsupported trackInputConversionAlgo");
0187   }
0188 
0189   const std::string &muInAlgo = iConfig.getParameter<std::string>("muonInputConversionAlgo");
0190   if (muInAlgo == "Emulator") {
0191     muonInput_ = std::make_unique<l1ct::GMTMuonDecoderEmulator>(
0192         iConfig.getParameter<edm::ParameterSet>("muonInputConversionParameters"));
0193   } else if (muInAlgo != "Ideal")
0194     throw cms::Exception("Configuration", "Unsupported muonInputConversionAlgo");
0195 
0196   const std::string &regalgo = iConfig.getParameter<std::string>("regionizerAlgo");
0197   if (regalgo == "Ideal") {
0198     regionizer_ =
0199         std::make_unique<l1ct::RegionizerEmulator>(iConfig.getParameter<edm::ParameterSet>("regionizerAlgoParameters"));
0200   } else if (regalgo == "Multififo") {
0201     regionizer_ = std::make_unique<l1ct::MultififoRegionizerEmulator>(
0202         iConfig.getParameter<edm::ParameterSet>("regionizerAlgoParameters"));
0203   } else if (regalgo == "TDR") {
0204     regionizer_ = std::make_unique<l1ct::TDRRegionizerEmulator>(
0205         iConfig.getParameter<edm::ParameterSet>("regionizerAlgoParameters"));
0206   } else
0207     throw cms::Exception("Configuration", "Unsupported regionizerAlgo");
0208 
0209   const std::string &algo = iConfig.getParameter<std::string>("pfAlgo");
0210   if (algo == "PFAlgo3") {
0211     l1pfalgo_ = std::make_unique<l1ct::PFAlgo3Emulator>(iConfig.getParameter<edm::ParameterSet>("pfAlgoParameters"));
0212   } else if (algo == "PFAlgo2HGC") {
0213     l1pfalgo_ = std::make_unique<l1ct::PFAlgo2HGCEmulator>(iConfig.getParameter<edm::ParameterSet>("pfAlgoParameters"));
0214   } else if (algo == "PFAlgoDummy") {
0215     l1pfalgo_ =
0216         std::make_unique<l1ct::PFAlgoDummyEmulator>(iConfig.getParameter<edm::ParameterSet>("pfAlgoParameters"));
0217   } else
0218     throw cms::Exception("Configuration", "Unsupported pfAlgo");
0219 
0220   const std::string &pualgo = iConfig.getParameter<std::string>("puAlgo");
0221   if (pualgo == "LinearizedPuppi") {
0222     l1pualgo_ = std::make_unique<l1ct::LinPuppiEmulator>(iConfig.getParameter<edm::ParameterSet>("puAlgoParameters"));
0223   } else
0224     throw cms::Exception("Configuration", "Unsupported puAlgo");
0225 
0226   l1tkegalgo_ = std::make_unique<l1ct::PFTkEGAlgoEmulator>(
0227       l1ct::PFTkEGAlgoEmuConfig(iConfig.getParameter<edm::ParameterSet>("tkEgAlgoParameters")));
0228 
0229   l1tkegsorter_ =
0230       std::make_unique<l1ct::PFTkEGSorterEmulator>(iConfig.getParameter<edm::ParameterSet>("tkEgSorterParameters"));
0231 
0232   if (l1tkegalgo_->writeEgSta())
0233     produces<BXVector<l1t::EGamma>>("L1Eg");
0234   produces<l1t::TkElectronCollection>("L1TkEle");
0235   produces<l1t::TkElectronRegionalOutput>("L1TkElePerBoard");
0236   produces<l1t::TkEmCollection>("L1TkEm");
0237   produces<l1t::TkEmRegionalOutput>("L1TkEmPerBoard");
0238 
0239   emuTkVtx_ = iConfig.getParameter<bool>("vtxCollectionEmulation");
0240   if (emuTkVtx_) {
0241     tkVtxEmu_ = consumes<std::vector<l1t::VertexWord>>(iConfig.getParameter<edm::InputTag>("vtxCollection"));
0242   } else {
0243     extTkVtx_ = consumes<std::vector<l1t::Vertex>>(iConfig.getParameter<edm::InputTag>("vtxCollection"));
0244   }
0245   nVtx_ = iConfig.getParameter<int>("nVtx");
0246 
0247   const char *iprefix[4] = {"totNReg", "maxNReg", "totNSec", "maxNSec"};
0248   for (int i = 0; i <= l1muType; ++i) {
0249     for (int ip = 0; ip < 4; ++ip) {
0250       produces<unsigned int>(std::string(iprefix[ip]) + inputTypeName[i]);
0251     }
0252     produces<std::vector<unsigned>>(std::string("vecNReg") + inputTypeName[i]);
0253     produces<std::vector<unsigned>>(std::string("vecNSec") + inputTypeName[i]);
0254   }
0255   const char *oprefix[4] = {"totNPF", "maxNPF", "totNPuppi", "maxNPuppi"};
0256   for (int i = 0; i < l1ct::OutputRegion::nPFTypes; ++i) {
0257     for (int ip = 0; ip < 4; ++ip) {
0258       produces<unsigned int>(std::string(oprefix[ip]) + l1ct::OutputRegion::objTypeName[i]);
0259     }
0260     produces<std::vector<unsigned>>(std::string("vecNPF") + l1ct::OutputRegion::objTypeName[i]);
0261     produces<std::vector<unsigned>>(std::string("vecNPuppi") + l1ct::OutputRegion::objTypeName[i]);
0262   }
0263 
0264   initSectorsAndRegions(iConfig);
0265 }
0266 
0267 L1TCorrelatorLayer1Producer::~L1TCorrelatorLayer1Producer() {}
0268 
0269 void L1TCorrelatorLayer1Producer::beginStream(edm::StreamID id) {
0270   if (!regionDumpName_.empty()) {
0271     if (id == 0) {
0272       fRegionDump_.open(regionDumpName_.c_str(), std::ios::out | std::ios::binary);
0273     } else {
0274       edm::LogWarning("L1TCorrelatorLayer1Producer")
0275           << "Job running with multiple streams, but dump file will have only events on stream zero.";
0276     }
0277   }
0278 }
0279 
0280 // ------------ method called to produce the data  ------------
0281 void L1TCorrelatorLayer1Producer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0282   // clear the regions also at the beginning, in case one event didn't complete but the job continues on
0283   initEvent(iEvent);
0284 
0285   /// ------ READ TRACKS ----
0286   if (hasTracks_) {
0287     edm::Handle<l1t::PFTrackCollection> htracks;
0288     iEvent.getByToken(tkCands_, htracks);
0289     const auto &tracks = *htracks;
0290     for (unsigned int itk = 0, ntk = tracks.size(); itk < ntk; ++itk) {
0291       const auto &tk = tracks[itk];
0292       // adding objects to PF
0293       if (debugR_ > 0 && deltaR(tk.eta(), tk.phi(), debugEta_, debugPhi_) > debugR_)
0294         continue;
0295       if (tk.pt() > trkPt_) {
0296         addTrack(tk, l1t::PFTrackRef(htracks, itk));
0297       }
0298     }
0299   }
0300 
0301   /// ------ READ MUONS ----
0302   edm::Handle<l1t::SAMuonCollection> muons;
0303   iEvent.getByToken(muCands_, muons);
0304   for (unsigned int i = 0, n = muons->size(); i < n; ++i) {
0305     const l1t::SAMuon &mu = (*muons)[i];
0306     if (debugR_ > 0 && deltaR(mu.eta(), mu.phi(), debugEta_, debugPhi_) > debugR_)
0307       continue;
0308     addMuon(mu, l1t::PFCandidate::MuonRef(muons, i));
0309   }
0310   // ------ READ CALOS -----
0311   edm::Handle<l1t::PFClusterCollection> caloHandle;
0312   for (const auto &tag : emCands_) {
0313     iEvent.getByToken(tag, caloHandle);
0314     const auto &calos = *caloHandle;
0315     for (unsigned int ic = 0, nc = calos.size(); ic < nc; ++ic) {
0316       const auto &calo = calos[ic];
0317       if (debugR_ > 0 && deltaR(calo.eta(), calo.phi(), debugEta_, debugPhi_) > debugR_)
0318         continue;
0319       if (calo.pt() > emPtCut_)
0320         addEmCalo(calo, l1t::PFClusterRef(caloHandle, ic));
0321     }
0322   }
0323   for (const auto &tag : hadCands_) {
0324     iEvent.getByToken(tag, caloHandle);
0325     const auto &calos = *caloHandle;
0326     for (unsigned int ic = 0, nc = calos.size(); ic < nc; ++ic) {
0327       const auto &calo = calos[ic];
0328       if (debugR_ > 0 && deltaR(calo.eta(), calo.phi(), debugEta_, debugPhi_) > debugR_)
0329         continue;
0330       if (calo.pt() > hadPtCut_)
0331         addHadCalo(calo, l1t::PFClusterRef(caloHandle, ic));
0332     }
0333   }
0334 
0335   regionizer_->run(event_.decoded, event_.pfinputs);
0336 
0337   // First, get a copy of the discretized and corrected inputs, and write them out
0338   iEvent.put(fetchEmCalo(), "EmCalo");
0339   iEvent.put(fetchHadCalo(), "Calo");
0340   iEvent.put(fetchTracks(), "TK");
0341 
0342   // Then do the vertexing, and save it out
0343   std::vector<float> z0s;
0344   std::vector<std::pair<float, float>> ptsums;
0345   float z0 = 0;
0346   double ptsum = 0;
0347   l1t::VertexWord pvwd;
0348   // FIXME: collections seem to be already sorted
0349   if (emuTkVtx_) {
0350     edm::Handle<std::vector<l1t::VertexWord>> vtxEmuHandle;
0351     iEvent.getByToken(tkVtxEmu_, vtxEmuHandle);
0352     for (const auto &vtx : *vtxEmuHandle) {
0353       ptsums.push_back(std::pair<float, float>(vtx.pt(), vtx.z0()));
0354       if (ptsum == 0 || vtx.pt() > ptsum) {
0355         ptsum = vtx.pt();
0356         pvwd = vtx;
0357       }
0358     }
0359   } else {
0360     edm::Handle<std::vector<l1t::Vertex>> vtxHandle;
0361     iEvent.getByToken(extTkVtx_, vtxHandle);
0362     for (const auto &vtx : *vtxHandle) {
0363       ptsums.push_back(std::pair<float, float>(vtx.pt(), vtx.z0()));
0364       if (ptsum == 0 || vtx.pt() > ptsum) {
0365         ptsum = vtx.pt();
0366         z0 = vtx.z0();
0367       }
0368     }
0369     pvwd = l1t::VertexWord(1, z0, 1, ptsum, 1, 1, 1);
0370   }
0371   l1ct::PVObjEmu hwpv;
0372   hwpv.hwZ0 = l1ct::Scales::makeZ0(pvwd.z0());
0373   event_.pvs.push_back(hwpv);
0374   event_.pvs_emu.push_back(pvwd.vertexWord());
0375   //Do a quick histogram vertexing to get multiple vertices (Hack for the taus)
0376   if (nVtx_ > 1) {
0377     std::stable_sort(ptsums.begin(), ptsums.end(), [](const auto &a, const auto &b) { return a.first > b.first; });
0378     for (int i0 = 0; i0 < std::min(int(ptsums.size()), int(nVtx_)); i0++) {
0379       z0s.push_back(ptsums[i0].second);
0380     }
0381     for (unsigned int i = 1; i < z0s.size(); ++i) {
0382       l1ct::PVObjEmu hwpv;
0383       hwpv.hwZ0 = l1ct::Scales::makeZ0(z0s[i]);
0384       event_.pvs.push_back(hwpv);  //Skip emu
0385     }
0386   }
0387   // Then also save the tracks with a vertex cut
0388 #if 0
0389   iEvent.put(l1regions_.fetchTracks(/*ptmin=*/0.0, /*fromPV=*/true), "TKVtx");
0390 #endif
0391 
0392   // Then run PF in each region
0393   event_.out.resize(event_.pfinputs.size());
0394   for (unsigned int ir = 0, nr = event_.pfinputs.size(); ir < nr; ++ir) {
0395     l1pfalgo_->run(event_.pfinputs[ir], event_.out[ir]);
0396     l1pfalgo_->mergeNeutrals(event_.out[ir]);
0397     l1tkegalgo_->run(event_.pfinputs[ir], event_.out[ir]);
0398     l1tkegalgo_->runIso(event_.pfinputs[ir], event_.pvs, event_.out[ir]);
0399   }
0400 
0401   // Then run puppi (regionally)
0402   for (unsigned int ir = 0, nr = event_.pfinputs.size(); ir < nr; ++ir) {
0403     l1pualgo_->run(event_.pfinputs[ir], event_.pvs, event_.out[ir]);
0404     //l1pualgo_->runNeutralsPU(l1region, z0, -1., puGlobals);
0405   }
0406 
0407   // NOTE: This needs to happen before the EG sorting per board so that the EG objects
0408   // get a global reference to the EGSta before being mixed among differente regions
0409   std::vector<edm::Ref<BXVector<l1t::EGamma>>> egsta_refs;
0410   if (l1tkegalgo_->writeEgSta()) {
0411     putEgStaObjects(iEvent, "L1Eg", egsta_refs);
0412   }
0413 
0414   // l1tkegsorter_->setDebug(true);
0415   for (auto &board : event_.board_out) {
0416     l1tkegsorter_->run(event_.pfinputs, event_.out, board.region_index, board.egphoton);
0417     l1tkegsorter_->run(event_.pfinputs, event_.out, board.region_index, board.egelectron);
0418   }
0419 
0420   // save PF into the event
0421   iEvent.put(fetchPF(), "PF");
0422 
0423   // and save puppi
0424   putPuppi(iEvent);
0425 
0426   // save the EG objects
0427   putEgObjects(iEvent, l1tkegalgo_->writeEgSta(), egsta_refs, "L1TkEm", "L1TkEmPerBoard", "L1TkEle", "L1TkElePerBoard");
0428 
0429   // Then go do the multiplicities
0430   for (int i = 0; i <= l1muType; ++i) {
0431     auto vecInputs = vecSecInput(InputType(i));
0432     auto tm = totAndMax(*vecInputs);
0433     addUInt(tm.first, std::string("totNSec") + inputTypeName[i], iEvent);
0434     addUInt(tm.second, std::string("maxNSec") + inputTypeName[i], iEvent);
0435     iEvent.put(std::move(vecInputs), std::string("vecNSec") + inputTypeName[i]);
0436   }
0437   for (int i = 0; i <= l1muType; ++i) {
0438     auto vecInputs = vecRegInput(InputType(i));
0439     auto tm = totAndMax(*vecInputs);
0440     addUInt(tm.first, std::string("totNReg") + inputTypeName[i], iEvent);
0441     addUInt(tm.second, std::string("maxNReg") + inputTypeName[i], iEvent);
0442     iEvent.put(std::move(vecInputs), std::string("vecNReg") + inputTypeName[i]);
0443   }
0444   for (int i = 0; i < l1ct::OutputRegion::nPFTypes; ++i) {
0445     auto vecPF = vecOutput(OutputType(i), false);
0446     auto tmPF = totAndMax(*vecPF);
0447     addUInt(tmPF.first, std::string("totNPF") + l1ct::OutputRegion::objTypeName[i], iEvent);
0448     addUInt(tmPF.second, std::string("maxNPF") + l1ct::OutputRegion::objTypeName[i], iEvent);
0449     iEvent.put(std::move(vecPF), std::string("vecNPF") + l1ct::OutputRegion::objTypeName[i]);
0450     auto vecPuppi = vecOutput(OutputType(i), true);
0451     auto tmPuppi = totAndMax(*vecPuppi);
0452     addUInt(tmPuppi.first, std::string("totNPuppi") + l1ct::OutputRegion::objTypeName[i], iEvent);
0453     addUInt(tmPuppi.second, std::string("maxNPuppi") + l1ct::OutputRegion::objTypeName[i], iEvent);
0454     iEvent.put(std::move(vecPuppi), std::string("vecNPuppi") + l1ct::OutputRegion::objTypeName[i]);
0455   }
0456 
0457   if (fRegionDump_.is_open()) {
0458     event_.write(fRegionDump_);
0459   }
0460 
0461   // finally clear the regions
0462   event_.clear();
0463 }
0464 
0465 void L1TCorrelatorLayer1Producer::addUInt(unsigned int value, std::string iLabel, edm::Event &iEvent) {
0466   iEvent.put(std::make_unique<unsigned>(value), iLabel);
0467 }
0468 
0469 void L1TCorrelatorLayer1Producer::initSectorsAndRegions(const edm::ParameterSet &iConfig) {
0470   // the track finder geometry is fixed
0471   unsigned int TF_phiSlices = 9;
0472   float TF_phiWidth = 2 * M_PI / TF_phiSlices;
0473   event_.decoded.track.clear();
0474   for (unsigned int ieta = 0, neta = 2; ieta < neta; ++ieta) {
0475     for (unsigned int iphi = 0; iphi < TF_phiSlices; ++iphi) {
0476       float phiCenter = reco::reduceRange(iphi * TF_phiWidth);
0477       event_.decoded.track.emplace_back((ieta ? 0. : -2.5), (ieta ? 2.5 : 0.0), phiCenter, TF_phiWidth);
0478       event_.raw.track.emplace_back((ieta ? 0. : -2.5), (ieta ? 2.5 : 0.0), phiCenter, TF_phiWidth);
0479     }
0480   }
0481 
0482   event_.decoded.emcalo.clear();
0483   event_.decoded.hadcalo.clear();
0484   event_.raw.hgcalcluster.clear();
0485 
0486   for (const edm::ParameterSet &preg : iConfig.getParameter<std::vector<edm::ParameterSet>>("caloSectors")) {
0487     std::vector<double> etaBoundaries = preg.getParameter<std::vector<double>>("etaBoundaries");
0488     if (!std::is_sorted(etaBoundaries.begin(), etaBoundaries.end()))
0489       throw cms::Exception("Configuration", "caloSectors.etaBoundaries not sorted\n");
0490     unsigned int phiSlices = preg.getParameter<uint32_t>("phiSlices");
0491     float phiWidth = 2 * M_PI / phiSlices;
0492     if (phiWidth > 2 * l1ct::Scales::maxAbsPhi())
0493       throw cms::Exception("Configuration", "caloSectors phi range too large for phi_t data type");
0494     double phiZero = preg.getParameter<double>("phiZero");
0495     for (unsigned int ieta = 0, neta = etaBoundaries.size() - 1; ieta < neta; ++ieta) {
0496       float etaWidth = etaBoundaries[ieta + 1] - etaBoundaries[ieta];
0497       if (etaWidth > 2 * l1ct::Scales::maxAbsEta())
0498         throw cms::Exception("Configuration", "caloSectors eta range too large for eta_t data type");
0499       for (unsigned int iphi = 0; iphi < phiSlices; ++iphi) {
0500         float phiCenter = reco::reduceRange(iphi * phiWidth + phiZero);
0501         event_.decoded.hadcalo.emplace_back(etaBoundaries[ieta], etaBoundaries[ieta + 1], phiCenter, phiWidth);
0502         event_.decoded.emcalo.emplace_back(etaBoundaries[ieta], etaBoundaries[ieta + 1], phiCenter, phiWidth);
0503         event_.raw.hgcalcluster.emplace_back(etaBoundaries[ieta], etaBoundaries[ieta + 1], phiCenter, phiWidth);
0504       }
0505     }
0506   }
0507 
0508   event_.decoded.muon.region = l1ct::PFRegionEmu(0., 0.);  // centered at (0,0)
0509   event_.raw.muon.region = l1ct::PFRegionEmu(0., 0.);      // centered at (0,0)
0510 
0511   event_.pfinputs.clear();
0512   for (const edm::ParameterSet &preg : iConfig.getParameter<std::vector<edm::ParameterSet>>("regions")) {
0513     std::vector<double> etaBoundaries = preg.getParameter<std::vector<double>>("etaBoundaries");
0514     if (!std::is_sorted(etaBoundaries.begin(), etaBoundaries.end()))
0515       throw cms::Exception("Configuration", "regions.etaBoundaries not sorted\n");
0516     unsigned int phiSlices = preg.getParameter<uint32_t>("phiSlices");
0517     float etaExtra = preg.getParameter<double>("etaExtra");
0518     float phiExtra = preg.getParameter<double>("phiExtra");
0519     float phiWidth = 2 * M_PI / phiSlices;
0520     for (unsigned int ieta = 0, neta = etaBoundaries.size() - 1; ieta < neta; ++ieta) {
0521       for (unsigned int iphi = 0; iphi < phiSlices; ++iphi) {
0522         float phiCenter = reco::reduceRange(iphi * phiWidth);  //align with L1 TrackFinder phi sector indexing
0523         event_.pfinputs.emplace_back(
0524             etaBoundaries[ieta], etaBoundaries[ieta + 1], phiCenter, phiWidth, etaExtra, phiExtra);
0525       }
0526     }
0527   }
0528 
0529   event_.board_out.clear();
0530   const std::vector<edm::ParameterSet> &board_params = iConfig.getParameter<std::vector<edm::ParameterSet>>("boards");
0531   event_.board_out.resize(board_params.size());
0532   for (unsigned int bidx = 0; bidx < board_params.size(); bidx++) {
0533     event_.board_out[bidx].region_index = board_params[bidx].getParameter<std::vector<unsigned int>>("regions");
0534     float etaBoard = 0.;
0535     float phiBoard = 0.;
0536     for (auto ridx : event_.board_out[bidx].region_index) {
0537       etaBoard += event_.pfinputs[ridx].region.floatEtaCenter();
0538       phiBoard += event_.pfinputs[ridx].region.floatPhiCenter();
0539     }
0540     event_.board_out[bidx].eta = etaBoard / event_.board_out[bidx].region_index.size();
0541     event_.board_out[bidx].phi = phiBoard / event_.board_out[bidx].region_index.size();
0542   }
0543 }
0544 
0545 void L1TCorrelatorLayer1Producer::initEvent(const edm::Event &iEvent) {
0546   event_.clear();
0547   event_.run = iEvent.id().run();
0548   event_.lumi = iEvent.id().luminosityBlock();
0549   event_.event = iEvent.id().event();
0550   clusterRefMap_.clear();
0551   trackRefMap_.clear();
0552   muonRefMap_.clear();
0553 }
0554 
0555 void L1TCorrelatorLayer1Producer::addTrack(const l1t::PFTrack &t, l1t::PFTrackRef ref) {
0556   auto &rawsectors = event_.raw.track;
0557   auto &sectors = event_.decoded.track;
0558   assert(sectors.size() == 18 && rawsectors.size() == 18);
0559   int isec = t.track()->phiSector() + (t.eta() >= 0 ? 9 : 0);
0560   rawsectors[isec].obj.push_back(t.trackWord().getTrackWord());
0561   addDecodedTrack(sectors[isec], t);
0562   trackRefMap_[&t] = ref;
0563 }
0564 void L1TCorrelatorLayer1Producer::addMuon(const l1t::SAMuon &mu, l1t::PFCandidate::MuonRef ref) {
0565   event_.raw.muon.obj.emplace_back(mu.word());
0566   addDecodedMuon(event_.decoded.muon, mu);
0567   muonRefMap_[&mu] = ref;
0568 }
0569 void L1TCorrelatorLayer1Producer::addHadCalo(const l1t::PFCluster &c, l1t::PFClusterRef ref) {
0570   int sidx = 0;
0571   for (auto &sec : event_.decoded.hadcalo) {
0572     if (sec.region.contains(c.eta(), c.phi())) {
0573       addDecodedHadCalo(sec, c);
0574       if (writeRawHgcalCluster_)
0575         addRawHgcalCluster(event_.raw.hgcalcluster[sidx], c);
0576     }
0577     sidx++;
0578   }
0579   clusterRefMap_[&c] = ref;
0580 }
0581 void L1TCorrelatorLayer1Producer::addEmCalo(const l1t::PFCluster &c, l1t::PFClusterRef ref) {
0582   for (auto &sec : event_.decoded.emcalo) {
0583     if (sec.region.contains(c.eta(), c.phi())) {
0584       addDecodedEmCalo(sec, c);
0585     }
0586   }
0587   clusterRefMap_[&c] = ref;
0588 }
0589 
0590 void L1TCorrelatorLayer1Producer::addDecodedTrack(l1ct::DetectorSector<l1ct::TkObjEmu> &sec, const l1t::PFTrack &t) {
0591   std::pair<l1ct::TkObjEmu, bool> tkAndSel;
0592   if (trackInput_) {
0593     tkAndSel = trackInput_->decodeTrack(t.trackWord().getTrackWord(), sec.region);
0594   } else {
0595     tkAndSel.first.hwPt = l1ct::Scales::makePtFromFloat(t.pt());
0596     tkAndSel.first.hwEta =
0597         l1ct::Scales::makeGlbEta(t.caloEta()) -
0598         sec.region.hwEtaCenter;  // important to enforce that the region boundary is on a discrete value
0599     tkAndSel.first.hwPhi = l1ct::Scales::makePhi(sec.region.localPhi(t.caloPhi()));
0600     tkAndSel.first.hwCharge = t.charge() > 0;
0601     tkAndSel.first.hwQuality = t.quality();
0602     tkAndSel.first.hwDEta = l1ct::Scales::makeEta(t.eta() - t.caloEta());
0603     tkAndSel.first.hwDPhi = l1ct::Scales::makePhi(std::abs(reco::deltaPhi(t.phi(), t.caloPhi())));
0604     tkAndSel.first.hwZ0 = l1ct::Scales::makeZ0(t.vertex().Z());
0605     tkAndSel.first.hwDxy = 0;
0606     tkAndSel.second = t.quality() > 0;
0607   }
0608   // CMSSW-only extra info
0609   tkAndSel.first.hwChi2 = round(t.chi2() * 10);
0610   tkAndSel.first.hwStubs = t.nStubs();
0611   tkAndSel.first.simPt = t.pt();
0612   tkAndSel.first.simCaloEta = t.caloEta();
0613   tkAndSel.first.simCaloPhi = t.caloPhi();
0614   tkAndSel.first.simVtxEta = t.eta();
0615   tkAndSel.first.simVtxPhi = t.phi();
0616   tkAndSel.first.simZ0 = t.vertex().Z();
0617   tkAndSel.first.simD0 = t.vertex().Rho();
0618   tkAndSel.first.src = &t;
0619   // If the track fails, we set its pT to zero, so that the decoded tracks are still aligned with the raw tracks
0620   // Downstream, the regionizer will just ignore zero-momentum tracks
0621   if (!tkAndSel.second)
0622     tkAndSel.first.hwPt = 0;
0623   sec.obj.push_back(tkAndSel.first);
0624 }
0625 
0626 void L1TCorrelatorLayer1Producer::addDecodedMuon(l1ct::DetectorSector<l1ct::MuObjEmu> &sec, const l1t::SAMuon &t) {
0627   l1ct::MuObjEmu mu;
0628   if (muonInput_) {
0629     mu = muonInput_->decode(t.word());
0630   } else {
0631     mu.hwPt = l1ct::Scales::makePtFromFloat(t.pt());
0632     mu.hwEta = l1ct::Scales::makeGlbEta(t.eta());  // IMPORTANT: input is in global coordinates!
0633     mu.hwPhi = l1ct::Scales::makeGlbPhi(t.phi());
0634     mu.hwCharge = !t.hwCharge();
0635     mu.hwQuality = t.hwQual() / 2;
0636     mu.hwDEta = 0;
0637     mu.hwDPhi = 0;
0638     mu.hwZ0 = l1ct::Scales::makeZ0(t.vertex().Z());
0639     mu.hwDxy = 0;  // Dxy not defined yet
0640   }
0641   mu.src = &t;
0642   sec.obj.push_back(mu);
0643 }
0644 
0645 void L1TCorrelatorLayer1Producer::addDecodedHadCalo(l1ct::DetectorSector<l1ct::HadCaloObjEmu> &sec,
0646                                                     const l1t::PFCluster &c) {
0647   l1ct::HadCaloObjEmu calo;
0648   calo.hwPt = l1ct::Scales::makePtFromFloat(c.pt());
0649   calo.hwEta = l1ct::Scales::makeGlbEta(c.eta()) -
0650                sec.region.hwEtaCenter;  // important to enforce that the region boundary is on a discrete value
0651   calo.hwPhi = l1ct::Scales::makePhi(sec.region.localPhi(c.phi()));
0652   calo.hwEmPt = l1ct::Scales::makePtFromFloat(c.emEt());
0653   calo.hwEmID = c.hwEmID();
0654   calo.src = &c;
0655   sec.obj.push_back(calo);
0656 }
0657 
0658 void L1TCorrelatorLayer1Producer::addRawHgcalCluster(l1ct::DetectorSector<ap_uint<256>> &sec, const l1t::PFCluster &c) {
0659   ap_uint<256> cwrd = 0;
0660   ap_uint<14> w_pt = round(c.pt() / 0.25);
0661   ap_uint<14> w_empt = round(c.emEt() / 0.25);
0662   constexpr float ETAPHI_LSB = M_PI / 720;
0663   ap_int<9> w_eta = round(sec.region.localEta(c.eta()) / ETAPHI_LSB);
0664   ap_int<9> w_phi = round(sec.region.localPhi(c.phi()) / ETAPHI_LSB);
0665   ap_uint<10> w_qual = c.hwQual();
0666 
0667   cwrd(13, 0) = w_pt;
0668   cwrd(27, 14) = w_empt;
0669   cwrd(72, 64) = w_eta;
0670   cwrd(81, 73) = w_phi;
0671   cwrd(115, 106) = w_qual;
0672 
0673   sec.obj.push_back(cwrd);
0674 }
0675 
0676 void L1TCorrelatorLayer1Producer::addDecodedEmCalo(l1ct::DetectorSector<l1ct::EmCaloObjEmu> &sec,
0677                                                    const l1t::PFCluster &c) {
0678   l1ct::EmCaloObjEmu calo;
0679   calo.hwPt = l1ct::Scales::makePtFromFloat(c.pt());
0680   calo.hwEta = l1ct::Scales::makeGlbEta(c.eta()) -
0681                sec.region.hwEtaCenter;  // important to enforce that the region boundary is on a discrete value
0682   calo.hwPhi = l1ct::Scales::makePhi(sec.region.localPhi(c.phi()));
0683   calo.hwPtErr = l1ct::Scales::makePtFromFloat(c.ptError());
0684   calo.hwEmID = c.hwEmID();
0685   calo.src = &c;
0686   sec.obj.push_back(calo);
0687 }
0688 
0689 template <typename T>
0690 void L1TCorrelatorLayer1Producer::setRefs_(l1t::PFCandidate &pf, const T &p) const {
0691   if (p.srcCluster) {
0692     auto match = clusterRefMap_.find(p.srcCluster);
0693     if (match == clusterRefMap_.end()) {
0694       throw cms::Exception("CorruptData") << "Invalid cluster pointer in PF candidate id " << p.intId() << " pt "
0695                                           << p.floatPt() << " eta " << p.floatEta() << " phi " << p.floatPhi();
0696     }
0697     pf.setPFCluster(match->second);
0698   }
0699   if (p.srcTrack) {
0700     auto match = trackRefMap_.find(p.srcTrack);
0701     if (match == trackRefMap_.end()) {
0702       throw cms::Exception("CorruptData") << "Invalid track pointer in PF candidate id " << p.intId() << " pt "
0703                                           << p.floatPt() << " eta " << p.floatEta() << " phi " << p.floatPhi();
0704     }
0705     pf.setPFTrack(match->second);
0706   }
0707   if (p.srcMu) {
0708     auto match = muonRefMap_.find(p.srcMu);
0709     if (match == muonRefMap_.end()) {
0710       throw cms::Exception("CorruptData") << "Invalid muon pointer in PF candidate id " << p.intId() << " pt "
0711                                           << p.floatPt() << " eta " << p.floatEta() << " phi " << p.floatPhi();
0712     }
0713     pf.setMuon(match->second);
0714   }
0715 }
0716 
0717 template <>
0718 void L1TCorrelatorLayer1Producer::setRefs_<l1ct::PFNeutralObjEmu>(l1t::PFCandidate &pf,
0719                                                                   const l1ct::PFNeutralObjEmu &p) const {
0720   if (p.srcCluster) {
0721     auto match = clusterRefMap_.find(p.srcCluster);
0722     if (match == clusterRefMap_.end()) {
0723       throw cms::Exception("CorruptData") << "Invalid cluster pointer in PF candidate id " << p.intId() << " pt "
0724                                           << p.floatPt() << " eta " << p.floatEta() << " phi " << p.floatPhi();
0725     }
0726     pf.setPFCluster(match->second);
0727   }
0728 }
0729 
0730 template <>
0731 void L1TCorrelatorLayer1Producer::setRefs_<l1ct::HadCaloObjEmu>(l1t::PFCandidate &pf,
0732                                                                 const l1ct::HadCaloObjEmu &p) const {
0733   if (p.src) {
0734     auto match = clusterRefMap_.find(p.src);
0735     if (match == clusterRefMap_.end()) {
0736       throw cms::Exception("CorruptData") << "Invalid cluster pointer in hadcalo candidate  pt " << p.floatPt()
0737                                           << " eta " << p.floatEta() << " phi " << p.floatPhi();
0738     }
0739     pf.setPFCluster(match->second);
0740   }
0741 }
0742 
0743 template <>
0744 void L1TCorrelatorLayer1Producer::setRefs_<l1ct::EmCaloObjEmu>(l1t::PFCandidate &pf,
0745                                                                const l1ct::EmCaloObjEmu &p) const {
0746   if (p.src) {
0747     auto match = clusterRefMap_.find(p.src);
0748     if (match == clusterRefMap_.end()) {
0749       throw cms::Exception("CorruptData") << "Invalid cluster pointer in emcalo candidate  pt " << p.floatPt()
0750                                           << " eta " << p.floatEta() << " phi " << p.floatPhi();
0751     }
0752     pf.setPFCluster(match->second);
0753   }
0754 }
0755 
0756 template <>
0757 void L1TCorrelatorLayer1Producer::setRefs_<l1ct::TkObjEmu>(l1t::PFCandidate &pf, const l1ct::TkObjEmu &p) const {
0758   if (p.src) {
0759     auto match = trackRefMap_.find(p.src);
0760     if (match == trackRefMap_.end()) {
0761       throw cms::Exception("CorruptData") << "Invalid track pointer in track candidate  pt " << p.floatPt() << " eta "
0762                                           << p.floatEta() << " phi " << p.floatPhi();
0763     }
0764     pf.setPFTrack(match->second);
0765   }
0766 }
0767 
0768 std::unique_ptr<l1t::PFCandidateCollection> L1TCorrelatorLayer1Producer::fetchHadCalo() const {
0769   auto ret = std::make_unique<l1t::PFCandidateCollection>();
0770   for (const auto &r : event_.pfinputs) {
0771     const auto &reg = r.region;
0772     for (const auto &p : r.hadcalo) {
0773       if (p.hwPt == 0 || !reg.isFiducial(p))
0774         continue;
0775       reco::Particle::PolarLorentzVector p4(p.floatPt(), reg.floatGlbEtaOf(p), reg.floatGlbPhiOf(p), 0.13f);
0776       l1t::PFCandidate::ParticleType type = p.hwIsEM() ? l1t::PFCandidate::Photon : l1t::PFCandidate::NeutralHadron;
0777       ret->emplace_back(type, 0, p4, 1, p.intPt(), p.intEta(), p.intPhi());
0778       ret->back().setHwEmID(p.hwEmID);
0779       setRefs_(ret->back(), p);
0780     }
0781   }
0782   return ret;
0783 }
0784 std::unique_ptr<l1t::PFCandidateCollection> L1TCorrelatorLayer1Producer::fetchEmCalo() const {
0785   auto ret = std::make_unique<l1t::PFCandidateCollection>();
0786   for (const auto &r : event_.pfinputs) {
0787     const auto &reg = r.region;
0788     for (const auto &p : r.emcalo) {
0789       if (p.hwPt == 0 || !reg.isFiducial(p))
0790         continue;
0791       reco::Particle::PolarLorentzVector p4(p.floatPt(), reg.floatGlbEtaOf(p), reg.floatGlbPhiOf(p), 0.13f);
0792       ret->emplace_back(l1t::PFCandidate::Photon, 0, p4, 1, p.intPt(), p.intEta(), p.intPhi());
0793       ret->back().setHwEmID(p.hwEmID);
0794       setRefs_(ret->back(), p);
0795     }
0796   }
0797   return ret;
0798 }
0799 std::unique_ptr<l1t::PFCandidateCollection> L1TCorrelatorLayer1Producer::fetchTracks() const {
0800   auto ret = std::make_unique<l1t::PFCandidateCollection>();
0801   for (const auto &r : event_.pfinputs) {
0802     const auto &reg = r.region;
0803     for (const auto &p : r.track) {
0804       if (p.hwPt == 0 || !reg.isFiducial(p))
0805         continue;
0806       reco::Particle::PolarLorentzVector p4(
0807           p.floatPt(), reg.floatGlbEta(p.hwVtxEta()), reg.floatGlbPhi(p.hwVtxPhi()), 0.13f);
0808       ret->emplace_back(l1t::PFCandidate::ChargedHadron, p.intCharge(), p4, 1, p.intPt(), p.intEta(), p.intPhi());
0809       setRefs_(ret->back(), p);
0810     }
0811   }
0812   return ret;
0813 }
0814 
0815 std::unique_ptr<l1t::PFCandidateCollection> L1TCorrelatorLayer1Producer::fetchPF() const {
0816   auto ret = std::make_unique<l1t::PFCandidateCollection>();
0817   for (unsigned int ir = 0, nr = event_.pfinputs.size(); ir < nr; ++ir) {
0818     const auto &reg = event_.pfinputs[ir].region;
0819     for (const auto &p : event_.out[ir].pfcharged) {
0820       if (p.hwPt == 0 || !reg.isFiducial(p))
0821         continue;
0822       reco::Particle::PolarLorentzVector p4(
0823           p.floatPt(), reg.floatGlbEta(p.hwVtxEta()), reg.floatGlbPhi(p.hwVtxPhi()), 0.13f);
0824       l1t::PFCandidate::ParticleType type = l1t::PFCandidate::ChargedHadron;
0825       if (p.hwId.isMuon())
0826         type = l1t::PFCandidate::Muon;
0827       else if (p.hwId.isElectron())
0828         type = l1t::PFCandidate::Electron;
0829       ret->emplace_back(type, p.intCharge(), p4, 1, p.intPt(), p.intEta(), p.intPhi());
0830       ret->back().setZ0(p.floatZ0());
0831       ret->back().setDxy(p.floatDxy());
0832       ret->back().setHwZ0(p.hwZ0);
0833       ret->back().setHwDxy(p.hwDxy);
0834       ret->back().setHwTkQuality(p.hwTkQuality);
0835       setRefs_(ret->back(), p);
0836     }
0837     for (const auto &p : event_.out[ir].pfneutral) {
0838       if (p.hwPt == 0 || !reg.isFiducial(p))
0839         continue;
0840       reco::Particle::PolarLorentzVector p4(p.floatPt(), reg.floatGlbEtaOf(p), reg.floatGlbPhiOf(p), 0.13f);
0841       l1t::PFCandidate::ParticleType type =
0842           p.hwId.isPhoton() ? l1t::PFCandidate::Photon : l1t::PFCandidate::NeutralHadron;
0843       ret->emplace_back(type, 0, p4, 1, p.intPt(), p.intEta(), p.intPhi());
0844       ret->back().setHwEmID(p.hwEmID);
0845       setRefs_(ret->back(), p);
0846     }
0847   }
0848   return ret;
0849 }
0850 
0851 void L1TCorrelatorLayer1Producer::putPuppi(edm::Event &iEvent) const {
0852   auto refprod = iEvent.getRefBeforePut<l1t::PFCandidateCollection>("Puppi");
0853   auto coll = std::make_unique<l1t::PFCandidateCollection>();
0854   auto reg = std::make_unique<l1t::PFCandidateRegionalOutput>(refprod);
0855   std::vector<int> nobj;
0856   for (unsigned int ir = 0, nr = event_.pfinputs.size(); ir < nr; ++ir) {
0857     nobj.clear();
0858     for (const auto &p : event_.out[ir].puppi) {
0859       if (p.hwPt == 0)
0860         continue;
0861       // note: Puppi candidates are already in global coordinates & fiducial-only!
0862       l1t::PFCandidate::ParticleType type;
0863       float mass = 0.13f;
0864       if (p.hwId.charged()) {
0865         if (p.hwId.isMuon()) {
0866           type = l1t::PFCandidate::Muon;
0867           mass = 0.105;
0868         } else if (p.hwId.isElectron()) {
0869           type = l1t::PFCandidate::Electron;
0870           mass = 0.005;
0871         } else
0872           type = l1t::PFCandidate::ChargedHadron;
0873       } else {
0874         type = p.hwId.isPhoton() ? l1t::PFCandidate::Photon : l1t::PFCandidate::NeutralHadron;
0875         mass = p.hwId.isPhoton() ? 0.0 : 0.5;
0876       }
0877       reco::Particle::PolarLorentzVector p4(p.floatPt(), p.floatEta(), p.floatPhi(), mass);
0878       coll->emplace_back(type, p.intCharge(), p4, p.floatPuppiW(), p.intPt(), p.intEta(), p.intPhi());
0879       if (p.hwId.charged()) {
0880         coll->back().setZ0(p.floatZ0());
0881         coll->back().setDxy(p.floatDxy());
0882         coll->back().setHwZ0(p.hwZ0());
0883         coll->back().setHwDxy(p.hwDxy());
0884         coll->back().setHwTkQuality(p.hwTkQuality());
0885       } else {
0886         coll->back().setHwPuppiWeight(p.hwPuppiW());
0887         coll->back().setHwEmID(p.hwEmID());
0888       }
0889       coll->back().setEncodedPuppi64(p.pack().to_uint64());
0890       setRefs_(coll->back(), p);
0891       nobj.push_back(coll->size() - 1);
0892     }
0893     reg->addRegion(nobj, event_.pfinputs[ir].region.floatEtaCenter(), event_.pfinputs[ir].region.floatPhiCenter());
0894   }
0895   iEvent.put(std::move(coll), "Puppi");
0896   iEvent.put(std::move(reg), "PuppiRegional");
0897 }
0898 
0899 // NOTE: as a side effect we change the "sta_idx" of TkEle and TkEm objects to an index of the
0900 // vector of refs, for this reason this is not const. We could make this more explicit via arguments
0901 void L1TCorrelatorLayer1Producer::putEgStaObjects(edm::Event &iEvent,
0902                                                   const std::string &egLablel,
0903                                                   std::vector<edm::Ref<BXVector<l1t::EGamma>>> &egsta_refs) {
0904   auto egs = std::make_unique<BXVector<l1t::EGamma>>();
0905   edm::RefProd<BXVector<l1t::EGamma>> ref_egs = iEvent.getRefBeforePut<BXVector<l1t::EGamma>>(egLablel);
0906 
0907   edm::Ref<BXVector<l1t::EGamma>>::key_type idx = 0;
0908   // FIXME: in case more BXes are introduced shuld probably use egs->key(egs->end(bx));
0909 
0910   for (unsigned int ir = 0, nr = event_.pfinputs.size(); ir < nr; ++ir) {
0911     const auto &reg = event_.pfinputs[ir].region;
0912 
0913     std::vector<unsigned int> ref_pos(event_.out[ir].egsta.size());
0914 
0915     // EG standalone objects
0916     for (unsigned int ieg = 0, neg = event_.out[ir].egsta.size(); ieg < neg; ++ieg) {
0917       const auto &p = event_.out[ir].egsta[ieg];
0918       if (p.hwPt == 0 || !reg.isFiducial(p))
0919         continue;
0920       l1t::EGamma eg(
0921           reco::Candidate::PolarLorentzVector(p.floatPt(), reg.floatGlbEta(p.hwEta), reg.floatGlbPhi(p.hwPhi), 0.));
0922       eg.setHwQual(p.hwQual);
0923       egs->push_back(0, eg);
0924       egsta_refs.push_back(edm::Ref<BXVector<l1t::EGamma>>(ref_egs, idx++));
0925       ref_pos[ieg] = egsta_refs.size() - 1;
0926     }
0927 
0928     for (auto &egiso : event_.out[ir].egphoton) {
0929       if (egiso.hwPt == 0)
0930         continue;
0931       egiso.sta_idx = ref_pos[egiso.sta_idx];
0932     }
0933 
0934     for (auto &egele : event_.out[ir].egelectron) {
0935       if (egele.hwPt == 0)
0936         continue;
0937       egele.sta_idx = ref_pos[egele.sta_idx];
0938     }
0939   }
0940 
0941   iEvent.put(std::move(egs), egLablel);
0942 }
0943 
0944 void L1TCorrelatorLayer1Producer::putEgObjects(edm::Event &iEvent,
0945                                                const bool writeEgSta,
0946                                                const std::vector<edm::Ref<BXVector<l1t::EGamma>>> &egsta_refs,
0947                                                const std::string &tkEmLabel,
0948                                                const std::string &tkEmPerBoardLabel,
0949                                                const std::string &tkEleLabel,
0950                                                const std::string &tkElePerBoardLabel) const {
0951   auto tkems = std::make_unique<l1t::TkEmCollection>();
0952   auto tkemRefProd = iEvent.getRefBeforePut<l1t::TkEmCollection>(tkEmLabel);
0953   auto tkemPerBoard = std::make_unique<l1t::TkEmRegionalOutput>(tkemRefProd);
0954   auto tkeles = std::make_unique<l1t::TkElectronCollection>();
0955   auto tkeleRefProd = iEvent.getRefBeforePut<l1t::TkElectronCollection>(tkEleLabel);
0956   auto tkelePerBoard = std::make_unique<l1t::TkElectronRegionalOutput>(tkeleRefProd);
0957 
0958   // TkEG objects are written out after the per-board sorting.
0959   // The mapping to each board is saved into the regionalmap for further (stage-2 consumption)
0960   std::vector<int> nele_obj;
0961   std::vector<int> npho_obj;
0962 
0963   for (const auto &board : event_.board_out) {
0964     npho_obj.clear();
0965     for (const auto &egiso : board.egphoton) {
0966       if (egiso.hwPt == 0)
0967         continue;
0968 
0969       edm::Ref<BXVector<l1t::EGamma>> ref_egsta;
0970       if (writeEgSta) {
0971         ref_egsta = egsta_refs[egiso.sta_idx];
0972       } else {
0973         auto egptr = egiso.srcCluster->constituentsAndFractions()[0].first;
0974         ref_egsta =
0975             edm::Ref<BXVector<l1t::EGamma>>(egptr.id(), dynamic_cast<const l1t::EGamma *>(egptr.get()), egptr.key());
0976       }
0977 
0978       reco::Candidate::PolarLorentzVector mom(egiso.floatPt(), egiso.floatEta(), egiso.floatPhi(), 0.);
0979 
0980       l1t::TkEm tkem(reco::Candidate::LorentzVector(mom),
0981                      ref_egsta,
0982                      egiso.floatRelIso(l1ct::EGIsoObjEmu::IsoType::TkIso),
0983                      egiso.floatRelIso(l1ct::EGIsoObjEmu::IsoType::TkIsoPV));
0984       tkem.setHwQual(egiso.hwQual);
0985       tkem.setPFIsol(egiso.floatRelIso(l1ct::EGIsoObjEmu::IsoType::PfIso));
0986       tkem.setPFIsolPV(egiso.floatRelIso(l1ct::EGIsoObjEmu::IsoType::PfIsoPV));
0987       tkem.setEgBinaryWord(egiso.pack());
0988       tkems->push_back(tkem);
0989       npho_obj.push_back(tkems->size() - 1);
0990     }
0991     tkemPerBoard->addRegion(npho_obj, board.eta, board.phi);
0992 
0993     nele_obj.clear();
0994     for (const auto &egele : board.egelectron) {
0995       if (egele.hwPt == 0)
0996         continue;
0997 
0998       edm::Ref<BXVector<l1t::EGamma>> ref_egsta;
0999       if (writeEgSta) {
1000         ref_egsta = egsta_refs[egele.sta_idx];
1001       } else {
1002         auto egptr = egele.srcCluster->constituentsAndFractions()[0].first;
1003         ref_egsta =
1004             edm::Ref<BXVector<l1t::EGamma>>(egptr.id(), dynamic_cast<const l1t::EGamma *>(egptr.get()), egptr.key());
1005       }
1006 
1007       reco::Candidate::PolarLorentzVector mom(egele.floatPt(), egele.floatEta(), egele.floatPhi(), 0.);
1008 
1009       l1t::TkElectron tkele(reco::Candidate::LorentzVector(mom),
1010                             ref_egsta,
1011                             edm::refToPtr(egele.srcTrack->track()),
1012                             egele.floatRelIso(l1ct::EGIsoEleObjEmu::IsoType::TkIso));
1013       tkele.setHwQual(egele.hwQual);
1014       tkele.setPFIsol(egele.floatRelIso(l1ct::EGIsoEleObjEmu::IsoType::PfIso));
1015       tkele.setEgBinaryWord(egele.pack());
1016       tkeles->push_back(tkele);
1017       nele_obj.push_back(tkeles->size() - 1);
1018     }
1019     tkelePerBoard->addRegion(nele_obj, board.eta, board.phi);
1020   }
1021 
1022   iEvent.put(std::move(tkems), tkEmLabel);
1023   iEvent.put(std::move(tkemPerBoard), tkEmPerBoardLabel);
1024   iEvent.put(std::move(tkeles), tkEleLabel);
1025   iEvent.put(std::move(tkelePerBoard), tkElePerBoardLabel);
1026 }
1027 
1028 std::unique_ptr<std::vector<unsigned>> L1TCorrelatorLayer1Producer::vecSecInput(InputType t) const {
1029   auto v = std::make_unique<std::vector<unsigned>>();
1030   {
1031     switch (t) {
1032       case caloType:
1033         for (const auto &s : event_.decoded.hadcalo)
1034           v->push_back(s.size());
1035         break;
1036       case emcaloType:
1037         for (const auto &s : event_.decoded.emcalo)
1038           v->push_back(s.size());
1039         break;
1040       case trackType:
1041         for (const auto &s : event_.decoded.track)
1042           v->push_back(s.size());
1043         break;
1044       case l1muType:
1045         v->push_back(event_.decoded.muon.size());
1046         break;
1047     }
1048   }
1049   return v;
1050 }
1051 
1052 std::unique_ptr<std::vector<unsigned>> L1TCorrelatorLayer1Producer::vecRegInput(InputType t) const {
1053   auto v = std::make_unique<std::vector<unsigned>>();
1054   for (const auto &reg : event_.pfinputs) {
1055     switch (t) {
1056       case caloType:
1057         v->push_back(reg.hadcalo.size());
1058         break;
1059       case emcaloType:
1060         v->push_back(reg.emcalo.size());
1061         break;
1062       case trackType:
1063         v->push_back(reg.track.size());
1064         break;
1065       case l1muType:
1066         v->push_back(reg.muon.size());
1067         break;
1068     }
1069   }
1070   return v;
1071 }
1072 
1073 std::unique_ptr<std::vector<unsigned>> L1TCorrelatorLayer1Producer::vecOutput(OutputType i, bool usePuppi) const {
1074   auto v = std::make_unique<std::vector<unsigned>>();
1075   for (const auto &reg : event_.out) {
1076     v->push_back(reg.nObj(i, usePuppi));
1077   }
1078   return v;
1079 }
1080 std::pair<unsigned int, unsigned int> L1TCorrelatorLayer1Producer::totAndMax(
1081     const std::vector<unsigned> &perRegion) const {
1082   unsigned int ntot = 0, nmax = 0;
1083   for (unsigned ni : perRegion) {
1084     ntot += ni;
1085     nmax = std::max(nmax, ni);
1086   }
1087   return std::make_pair(ntot, nmax);
1088 }
1089 //define this as a plug-in
1090 #include "FWCore/Framework/interface/MakerMacros.h"
1091 DEFINE_FWK_MODULE(L1TCorrelatorLayer1Producer);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!