Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-08-30 02:10:49

0001 // -*- C++ -*-
0002 //
0003 // Package:    L1Trigger/L1TCaloSummary
0004 // Class:      L1TCaloSummary
0005 //
0006 /**\class L1TCaloSummary L1TCaloSummary.cc L1Trigger/L1TCaloSummary/plugins/L1TCaloSummary.cc
0007 
0008    Description: The package L1Trigger/L1TCaloSummary is prepared for monitoring the CMS Layer-1 Calorimeter Trigger.
0009 
0010    Implementation:
0011    It prepares region objects and puts them in the event
0012 */
0013 //
0014 // Original Author:  Sridhara Dasu
0015 //         Created:  Sat, 14 Nov 2015 14:18:27 GMT
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 
0022 // user include files
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/stream/EDProducer.h"
0025 
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028 
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 
0031 #include "DataFormats/EcalDigi/interface/EcalDigiCollections.h"
0032 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0033 
0034 #include "L1Trigger/L1TCaloLayer1/src/UCTLayer1.hh"
0035 #include "L1Trigger/L1TCaloLayer1/src/UCTCrate.hh"
0036 #include "L1Trigger/L1TCaloLayer1/src/UCTCard.hh"
0037 #include "L1Trigger/L1TCaloLayer1/src/UCTRegion.hh"
0038 #include "L1Trigger/L1TCaloLayer1/src/UCTTower.hh"
0039 #include "L1Trigger/L1TCaloLayer1/src/UCTGeometry.hh"
0040 
0041 #include "L1Trigger/L1TCaloLayer1/src/UCTObject.hh"
0042 #include "L1Trigger/L1TCaloLayer1/src/UCTSummaryCard.hh"
0043 #include "L1Trigger/L1TCaloLayer1/src/UCTGeometryExtended.hh"
0044 
0045 #include "DataFormats/L1Trigger/interface/L1EmParticle.h"
0046 #include "DataFormats/L1Trigger/interface/L1EmParticleFwd.h"
0047 #include "DataFormats/L1Trigger/interface/L1JetParticle.h"
0048 #include "DataFormats/L1Trigger/interface/L1JetParticleFwd.h"
0049 #include "DataFormats/L1Trigger/interface/L1EtMissParticle.h"
0050 #include "DataFormats/L1Trigger/interface/L1EtMissParticleFwd.h"
0051 
0052 #include "DataFormats/L1CaloTrigger/interface/L1CaloCollections.h"
0053 #include "DataFormats/L1CaloTrigger/interface/L1CaloRegion.h"
0054 #include "DataFormats/L1CaloTrigger/interface/CICADA.h"
0055 
0056 #include "DataFormats/Math/interface/LorentzVector.h"
0057 
0058 #include "L1Trigger/L1TCaloLayer1/src/UCTLogging.hh"
0059 #include <bitset>
0060 
0061 #include <string>
0062 #include <sstream>
0063 
0064 //Anomaly detection includes
0065 #include "ap_fixed.h"
0066 #include "hls4ml/emulator.h"
0067 
0068 using namespace l1tcalo;
0069 using namespace l1extra;
0070 using namespace std;
0071 
0072 //
0073 // class declaration
0074 //
0075 
0076 template <class INPUT, class OUTPUT>
0077 class L1TCaloSummary : public edm::stream::EDProducer<> {
0078 public:
0079   explicit L1TCaloSummary(const edm::ParameterSet&);
0080   ~L1TCaloSummary() override = default;
0081 
0082   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0083 
0084 private:
0085   //void beginJob() override;
0086   void produce(edm::Event&, const edm::EventSetup&) override;
0087   //void endJob() override;
0088 
0089   void print();
0090 
0091   // ----------member data ---------------------------
0092 
0093   uint32_t nPumBins;
0094 
0095   std::vector<std::vector<std::vector<uint32_t>>> pumLUT;
0096 
0097   double caloScaleFactor;
0098 
0099   uint32_t jetSeed;
0100   uint32_t tauSeed;
0101   float tauIsolationFactor;
0102   uint32_t eGammaSeed;
0103   double eGammaIsolationFactor;
0104   double boostedJetPtFactor;
0105 
0106   bool verbose;
0107   int fwVersion;
0108 
0109   edm::EDGetTokenT<L1CaloRegionCollection> regionToken;
0110   edm::EDGetTokenT<L1CaloRegionCollection> backupRegionToken;
0111 
0112   UCTLayer1* layer1;
0113 
0114   hls4mlEmulator::ModelLoader loader;
0115   std::shared_ptr<hls4mlEmulator::Model> model;
0116 
0117   bool overwriteWithTestPatterns;
0118   std::vector<edm::ParameterSet> testPatterns;
0119 };
0120 
0121 //
0122 // constants, enums and typedefs
0123 //
0124 
0125 //
0126 // static data member definitions
0127 //
0128 
0129 //
0130 // constructors and destructor
0131 //
0132 template <class INPUT, class OUTPUT>
0133 L1TCaloSummary<INPUT, OUTPUT>::L1TCaloSummary(const edm::ParameterSet& iConfig)
0134     : nPumBins(iConfig.getParameter<unsigned int>("nPumBins")),
0135       pumLUT(nPumBins, std::vector<std::vector<uint32_t>>(2, std::vector<uint32_t>(13))),
0136       caloScaleFactor(iConfig.getParameter<double>("caloScaleFactor")),
0137       jetSeed(iConfig.getParameter<unsigned int>("jetSeed")),
0138       tauSeed(iConfig.getParameter<unsigned int>("tauSeed")),
0139       tauIsolationFactor(iConfig.getParameter<double>("tauIsolationFactor")),
0140       eGammaSeed(iConfig.getParameter<unsigned int>("eGammaSeed")),
0141       eGammaIsolationFactor(iConfig.getParameter<double>("eGammaIsolationFactor")),
0142       boostedJetPtFactor(iConfig.getParameter<double>("boostedJetPtFactor")),
0143       verbose(iConfig.getParameter<bool>("verbose")),
0144       fwVersion(iConfig.getParameter<int>("firmwareVersion")),
0145       regionToken(consumes<L1CaloRegionCollection>(iConfig.getParameter<edm::InputTag>("caloLayer1Regions"))),
0146       backupRegionToken(consumes<L1CaloRegionCollection>(edm::InputTag("simCaloStage2Layer1Digis"))),
0147       loader(hls4mlEmulator::ModelLoader(iConfig.getParameter<string>("CICADAModelVersion"))),
0148       overwriteWithTestPatterns(iConfig.getParameter<bool>("useTestPatterns")),
0149       testPatterns(iConfig.getParameter<std::vector<edm::ParameterSet>>("testPatterns")) {
0150   std::vector<double> pumLUTData;
0151   char pumLUTString[10];
0152   for (uint32_t pumBin = 0; pumBin < nPumBins; pumBin++) {
0153     for (uint32_t side = 0; side < 2; side++) {
0154       if (side == 0)
0155         sprintf(pumLUTString, "pumLUT%2.2dp", pumBin);
0156       else
0157         sprintf(pumLUTString, "pumLUT%2.2dn", pumBin);
0158       pumLUTData = iConfig.getParameter<std::vector<double>>(pumLUTString);
0159       for (uint32_t iEta = 0; iEta < std::max((uint32_t)pumLUTData.size(), MaxUCTRegionsEta); iEta++) {
0160         pumLUT[pumBin][side][iEta] = (uint32_t)round(pumLUTData[iEta] / caloScaleFactor);
0161       }
0162       if (pumLUTData.size() != (MaxUCTRegionsEta))
0163         edm::LogError("L1TCaloSummary") << "PUM LUT Data size integrity check failed; Expected size = "
0164                                         << MaxUCTRegionsEta << "; Provided size = " << pumLUTData.size()
0165                                         << "; Will use what is provided :(" << std::endl;
0166     }
0167   }
0168   produces<L1JetParticleCollection>("Boosted");
0169 
0170   //anomaly trigger loading
0171   model = loader.load_model();
0172   produces<l1t::CICADABxCollection>("CICADAScore");
0173 }
0174 
0175 //
0176 // member functions
0177 //
0178 
0179 // ------------ method called to produce the data  ------------
0180 template <class INPUT, class OUTPUT>
0181 void L1TCaloSummary<INPUT, OUTPUT>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0182   using namespace edm;
0183 
0184   std::unique_ptr<L1JetParticleCollection> bJetCands(new L1JetParticleCollection);
0185 
0186   std::unique_ptr<l1t::CICADABxCollection> CICADAScore = std::make_unique<l1t::CICADABxCollection>();
0187   CICADAScore->setBXRange(-2, 2);
0188 
0189   UCTGeometry g;
0190 
0191   // Here we read region data from the region collection created by L1TCaloLayer1 instead of
0192   // independently creating regions from TPGs for processing by the summary card. This results
0193   // in a single region vector of size 252 whereas from independent creation we had 3*6 vectors
0194   // of size 7*2. Indices are mapped in UCTSummaryCard accordingly.
0195   UCTSummaryCard summaryCard =
0196       UCTSummaryCard(&pumLUT, jetSeed, tauSeed, tauIsolationFactor, eGammaSeed, eGammaIsolationFactor);
0197   std::vector<UCTRegion*> inputRegions;
0198   inputRegions.clear();
0199   edm::Handle<std::vector<L1CaloRegion>> regionCollection;
0200   if (!iEvent.getByToken(regionToken, regionCollection))
0201     edm::LogError("L1TCaloSummary") << "UCT: Failed to get regions from region collection!";
0202   iEvent.getByToken(regionToken, regionCollection);
0203 
0204   if (regionCollection->size() == 0) {
0205     iEvent.getByToken(backupRegionToken, regionCollection);
0206     edm::LogWarning("L1TCaloSummary") << "Switched to emulated regions since data regions was empty.\n";
0207   }
0208 
0209   //Model input
0210   //This is done as a flat vector input, but future versions may involve 2D input
0211   //This will have to be handled later
0212   INPUT modelInput[252];
0213   for (const L1CaloRegion& i : *regionCollection) {
0214     UCTRegionIndex r = g.getUCTRegionIndexFromL1CaloRegion(i.gctEta(), i.gctPhi());
0215     UCTTowerIndex t = g.getUCTTowerIndexFromL1CaloRegion(r, i.raw());
0216     uint32_t absCaloEta = std::abs(t.first);
0217     uint32_t absCaloPhi = std::abs(t.second);
0218     bool negativeEta = false;
0219     if (t.first < 0)
0220       negativeEta = true;
0221     uint32_t crate = g.getCrate(t.first, t.second);
0222     uint32_t card = g.getCard(t.first, t.second);
0223     uint32_t region = g.getRegion(absCaloEta, absCaloPhi);
0224     UCTRegion* test = new UCTRegion(crate, card, negativeEta, region, fwVersion);
0225     test->setRegionSummary(i.raw());
0226     inputRegions.push_back(test);
0227     //This *should* fill the tensor in the proper order to be fed to the anomaly model
0228     //We take 4 off of the GCT eta/iEta.
0229     //iEta taken from this ranges from 4-17, (I assume reserving lower and higher for forward regions)
0230     //So our first index, index 0, is technically iEta=4, and so-on.
0231     //CICADA reads this as a flat vector
0232     modelInput[14 * i.gctPhi() + (i.gctEta() - 4)] = i.et();
0233   }
0234   // Check if we're using test patterns. If so, we overwrite the inputs with a test pattern
0235   if (overwriteWithTestPatterns) {
0236     unsigned int evt = iEvent.id().event();
0237     unsigned int totalTestPatterns = testPatterns.size();
0238     unsigned int patternElement = evt % totalTestPatterns;
0239     const edm::ParameterSet& element = testPatterns.at(patternElement);
0240     std::stringstream inputStream;
0241     std::string PhiRowString;
0242 
0243     edm::LogWarning("L1TCaloSummary") << "Overwriting existing CICADA input with test pattern!\n";
0244 
0245     for (unsigned short int iPhi = 1; iPhi <= 18; ++iPhi) {
0246       PhiRowString = "";
0247       std::stringstream PhiRowStringStream;
0248       PhiRowStringStream << "iPhi_" << iPhi;
0249       PhiRowString = PhiRowStringStream.str();
0250       std::vector<unsigned int> phiRow = element.getParameter<std::vector<unsigned int>>(PhiRowString);
0251       for (unsigned short int iEta = 1; iEta <= 14; ++iEta) {
0252         modelInput[14 * (iPhi - 1) + (iEta - 1)] = phiRow.at(iEta - 1);
0253         inputStream << phiRow.at(iEta - 1) << " ";
0254       }
0255       inputStream << "\n";
0256     }
0257     edm::LogInfo("L1TCaloSummary") << "Input Stream:\n" << inputStream.str();
0258   }
0259 
0260   //Extract model output
0261   OUTPUT modelResult[1] = {
0262       OUTPUT("0.0", 10)};  //the 10 here refers to the fact that we read in "0.0" as a decimal number
0263   model->prepare_input(modelInput);
0264   model->predict();
0265   model->read_result(modelResult);
0266 
0267   CICADAScore->push_back(0, modelResult[0].to_float());
0268 
0269   if (overwriteWithTestPatterns)
0270     edm::LogInfo("L1TCaloSummary") << "Test Pattern Output: " << CICADAScore->at(0, 0);
0271 
0272   summaryCard.setRegionData(inputRegions);
0273 
0274   if (!summaryCard.process()) {
0275     edm::LogError("L1TCaloSummary") << "UCT: Failed to process summary card" << std::endl;
0276     exit(1);
0277   }
0278 
0279   double pt = 0;
0280   double eta = -999.;
0281   double phi = -999.;
0282   double mass = 0;
0283 
0284   std::list<UCTObject*> boostedJetObjs = summaryCard.getBoostedJetObjs();
0285   for (std::list<UCTObject*>::const_iterator i = boostedJetObjs.begin(); i != boostedJetObjs.end(); i++) {
0286     const UCTObject* object = *i;
0287     pt = ((double)object->et()) * caloScaleFactor * boostedJetPtFactor;
0288     eta = g.getUCTTowerEta(object->iEta());
0289     phi = g.getUCTTowerPhi(object->iPhi());
0290     bitset<3> activeRegionEtaPattern = 0;
0291     for (uint32_t iEta = 0; iEta < 3; iEta++) {
0292       bool activeStrip = false;
0293       for (uint32_t iPhi = 0; iPhi < 3; iPhi++) {
0294         if (object->boostedJetRegionET()[3 * iEta + iPhi] > 30 &&
0295             object->boostedJetRegionET()[3 * iEta + iPhi] > object->et() * 0.0625)
0296           activeStrip = true;
0297       }
0298       if (activeStrip)
0299         activeRegionEtaPattern |= (0x1 << iEta);
0300     }
0301     bitset<3> activeRegionPhiPattern = 0;
0302     for (uint32_t iPhi = 0; iPhi < 3; iPhi++) {
0303       bool activeStrip = false;
0304       for (uint32_t iEta = 0; iEta < 3; iEta++) {
0305         if (object->boostedJetRegionET()[3 * iEta + iPhi] > 30 &&
0306             object->boostedJetRegionET()[3 * iEta + iPhi] > object->et() * 0.0625)
0307           activeStrip = true;
0308       }
0309       if (activeStrip)
0310         activeRegionPhiPattern |= (0x1 << iPhi);
0311     }
0312     string regionEta = activeRegionEtaPattern.to_string<char, std::string::traits_type, std::string::allocator_type>();
0313     string regionPhi = activeRegionPhiPattern.to_string<char, std::string::traits_type, std::string::allocator_type>();
0314 
0315     bool centralHighest = object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[0] &&
0316                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[1] &&
0317                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[2] &&
0318                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[3] &&
0319                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[5] &&
0320                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[6] &&
0321                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[7] &&
0322                           object->boostedJetRegionET()[4] >= object->boostedJetRegionET()[8];
0323 
0324     if (abs(eta) < 2.5 && ((regionEta == "101" && (regionPhi == "110" || regionPhi == "101" || regionPhi == "010")) ||
0325                            ((regionEta == "110" || regionEta == "101" || regionEta == "010") && regionPhi == "101") ||
0326                            (regionEta == "111" && (regionPhi == "110" || regionPhi == "010")) ||
0327                            ((regionEta == "110" || regionEta == "010") && regionPhi == "111") ||
0328                            ((regionEta == "010" || regionPhi == "010" || regionEta == "110" || regionPhi == "110" ||
0329                              regionEta == "011" || regionPhi == "011") &&
0330                             centralHighest)))
0331       bJetCands->push_back(L1JetParticle(math::PtEtaPhiMLorentzVector(pt, eta, phi, mass), L1JetParticle::kCentral));
0332   }
0333 
0334   iEvent.put(std::move(bJetCands), "Boosted");
0335   //Write out anomaly score
0336   iEvent.put(std::move(CICADAScore), "CICADAScore");
0337 }
0338 
0339 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0340 template <class INPUT, class OUTPUT>
0341 void L1TCaloSummary<INPUT, OUTPUT>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0342   //The following says we do not know what parameters are allowed so do no validation
0343   // Please change this to state exactly what you do use, even if it is no parameters
0344   edm::ParameterSetDescription desc;
0345   desc.setUnknown();
0346   descriptions.addDefault(desc);
0347 }
0348 
0349 // Initial version, X.0.0, input/output typing
0350 typedef L1TCaloSummary<ap_ufixed<10, 10>, ap_fixed<11, 5>> L1TCaloSummary_CICADA_v1p0p0;
0351 typedef L1TCaloSummary<ap_uint<10>, ap_ufixed<16, 8>> L1TCaloSummary_CICADA_v2p0p0;
0352 DEFINE_FWK_MODULE(L1TCaloSummary_CICADA_v1p0p0);
0353 DEFINE_FWK_MODULE(L1TCaloSummary_CICADA_v2p0p0);
0354 // X.1.0 version input.output typing
0355 typedef L1TCaloSummary<ap_ufixed<10, 10>, ap_fixed<11, 5>> L1TCaloSummary_CICADA_v1p1p0;
0356 typedef L1TCaloSummary<ap_uint<10>, ap_ufixed<16, 8>> L1TCaloSummary_CICADA_v2p1p0;
0357 DEFINE_FWK_MODULE(L1TCaloSummary_CICADA_v1p1p0);
0358 DEFINE_FWK_MODULE(L1TCaloSummary_CICADA_v2p1p0);
0359 // X.1.1 version input/output typing
0360 typedef L1TCaloSummary<ap_uint<10>, ap_ufixed<16, 8, AP_RND, AP_SAT, AP_SAT>> L1TCaloSummary_CICADA_vXp1p1;
0361 DEFINE_FWK_MODULE(L1TCaloSummary_CICADA_vXp1p1);
0362 // X.1.2 version input/output typing
0363 typedef L1TCaloSummary<ap_uint<10>, ap_ufixed<16, 8, AP_RND_CONV, AP_SAT>> L1TCaloSummary_CICADA_vXp1p2;
0364 DEFINE_FWK_MODULE(L1TCaloSummary_CICADA_vXp1p2);