Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-12-05 03:12:13

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