Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-11-19 23:59:18

0001 // Note to self: the implementation uses TH1F's to store the L1 and HLT rates.
0002 // Assuming a maximum rate of 100 kHz times a period of 23.31 s, one needs to store counts up to ~2.3e6.
0003 // A "float" has 24 bits of precision, so it can store up to 2**24 ~ 16.7e6 without loss of precision.
0004 
0005 // C++ headers
0006 #include <string>
0007 #include <cstring>
0008 
0009 #include <fmt/printf.h>
0010 
0011 // boost headers
0012 #include <boost/regex.hpp>
0013 
0014 // Root headers
0015 #include <TH1F.h>
0016 
0017 // CMSSW headers
0018 #include "FWCore/Framework/interface/Frameworkfwd.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/LuminosityBlock.h"
0021 #include "FWCore/Framework/interface/Run.h"
0022 #include "FWCore/Framework/interface/EventSetup.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0024 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0025 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0026 #include "FWCore/ParameterSet/interface/Registry.h"
0027 #include "FWCore/ServiceRegistry/interface/Service.h"
0028 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0029 #include "DataFormats/Provenance/interface/ProcessHistory.h"
0030 #include "DataFormats/Common/interface/TriggerResults.h"
0031 #include "DataFormats/L1TGlobal/interface/GlobalAlgBlk.h"
0032 #include "CondFormats/DataRecord/interface/L1TUtmTriggerMenuRcd.h"
0033 #include "CondFormats/L1TObjects/interface/L1TUtmTriggerMenu.h"
0034 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0035 #include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h"
0036 
0037 namespace {
0038 
0039   struct RunBasedHistograms {
0040     typedef dqm::reco::MonitorElement MonitorElement;
0041     // HLT configuration
0042     struct HLTIndices {
0043       unsigned int index_l1_seed;
0044       unsigned int index_prescale;
0045 
0046       HLTIndices() : index_l1_seed((unsigned int)-1), index_prescale((unsigned int)-1) {}
0047     };
0048 
0049     HLTConfigProvider hltConfig;
0050     std::vector<HLTIndices> hltIndices;
0051 
0052     std::vector<std::vector<unsigned int>> datasets;
0053     std::vector<std::vector<unsigned int>> streams;
0054 
0055     // L1T and HLT rate plots
0056 
0057     // per-path HLT plots
0058     struct HLTRatesPlots {
0059       dqm::reco::MonitorElement *pass_l1_seed;
0060       dqm::reco::MonitorElement *pass_prescale;
0061       dqm::reco::MonitorElement *accept;
0062       dqm::reco::MonitorElement *reject;
0063       dqm::reco::MonitorElement *error;
0064     };
0065 
0066     // overall event count and event types
0067     dqm::reco::MonitorElement *events_processed;
0068     std::vector<dqm::reco::MonitorElement *> tcds_counts;
0069 
0070     // L1T triggers
0071     std::vector<dqm::reco::MonitorElement *> l1t_counts;
0072 
0073     // HLT triggers
0074     std::vector<std::vector<HLTRatesPlots>> hlt_by_dataset_counts;
0075 
0076     // datasets
0077     std::vector<dqm::reco::MonitorElement *> dataset_counts;
0078 
0079     // streams
0080     std::vector<dqm::reco::MonitorElement *> stream_counts;
0081 
0082     RunBasedHistograms()
0083         :  // L1T and HLT configuration
0084           hltConfig(),
0085           hltIndices(),
0086           datasets(),
0087           streams(),
0088           // overall event count and event types
0089           events_processed(),
0090           tcds_counts(),
0091           // L1T triggers
0092           l1t_counts(),
0093           // HLT triggers
0094           hlt_by_dataset_counts(),
0095           // datasets
0096           dataset_counts(),
0097           // streams
0098           stream_counts() {}
0099   };
0100 }  // namespace
0101 
0102 class TriggerRatesMonitor : public DQMGlobalEDAnalyzer<RunBasedHistograms> {
0103 public:
0104   explicit TriggerRatesMonitor(edm::ParameterSet const &);
0105   ~TriggerRatesMonitor() override = default;
0106 
0107   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0108 
0109 private:
0110   void dqmBeginRun(edm::Run const &, edm::EventSetup const &, RunBasedHistograms &) const override;
0111   void bookHistograms(DQMStore::IBooker &,
0112                       edm::Run const &,
0113                       edm::EventSetup const &,
0114                       RunBasedHistograms &) const override;
0115   void dqmAnalyze(edm::Event const &, edm::EventSetup const &, RunBasedHistograms const &) const override;
0116 
0117   // TCDS trigger types
0118   // see https://twiki.cern.ch/twiki/bin/viewauth/CMS/TcdsEventRecord
0119   static constexpr const char *const s_tcds_trigger_types[] = {
0120       "Empty",          //  0 - No trigger
0121       "Physics",        //  1 - GT trigger
0122       "Calibration",    //  2 - Sequence trigger (calibration)
0123       "Random",         //  3 - Random trigger
0124       "Auxiliary",      //  4 - Auxiliary (CPM front panel NIM input) trigger
0125       nullptr,          //  5 - reserved
0126       nullptr,          //  6 - reserved
0127       nullptr,          //  7 - reserved
0128       "Cyclic",         //  8 - Cyclic trigger
0129       "Bunch-pattern",  //  9 - Bunch-pattern trigger
0130       "Software",       // 10 - Software trigger
0131       "TTS",            // 11 - TTS-sourced trigger
0132       nullptr,          // 12 - reserved
0133       nullptr,          // 13 - reserved
0134       nullptr,          // 14 - reserved
0135       nullptr           // 15 - reserved
0136   };
0137 
0138   // module configuration
0139   const edm::ESGetToken<L1TUtmTriggerMenu, L1TUtmTriggerMenuRcd> m_l1tMenuToken;
0140   const edm::EDGetTokenT<GlobalAlgBlkBxCollection> m_l1t_results;
0141   const edm::EDGetTokenT<edm::TriggerResults> m_hlt_results;
0142   const std::string m_dqm_path;
0143   const uint32_t m_lumisections_range;
0144 };
0145 
0146 // definition
0147 constexpr const char *const TriggerRatesMonitor::s_tcds_trigger_types[];
0148 
0149 void TriggerRatesMonitor::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0150   edm::ParameterSetDescription desc;
0151   desc.addUntracked<edm::InputTag>("l1tResults", edm::InputTag("gtStage2Digis"));
0152   desc.addUntracked<edm::InputTag>("hltResults", edm::InputTag("TriggerResults"));
0153   desc.addUntracked<std::string>("dqmPath", "HLT/TriggerRates");
0154   desc.addUntracked<uint32_t>("lumisectionRange", 2500);  // ~16 hours
0155   descriptions.add("triggerRatesMonitor", desc);
0156 }
0157 
0158 TriggerRatesMonitor::TriggerRatesMonitor(edm::ParameterSet const &config)
0159     :  // module configuration
0160       m_l1tMenuToken{esConsumes<edm::Transition::BeginRun>()},
0161       m_l1t_results(consumes<GlobalAlgBlkBxCollection>(config.getUntrackedParameter<edm::InputTag>("l1tResults"))),
0162       m_hlt_results(consumes<edm::TriggerResults>(config.getUntrackedParameter<edm::InputTag>("hltResults"))),
0163       m_dqm_path(config.getUntrackedParameter<std::string>("dqmPath")),
0164       m_lumisections_range(config.getUntrackedParameter<uint32_t>("lumisectionRange")) {}
0165 
0166 void TriggerRatesMonitor::dqmBeginRun(edm::Run const &run,
0167                                       edm::EventSetup const &setup,
0168                                       RunBasedHistograms &histograms) const {
0169   histograms.tcds_counts.clear();
0170   histograms.tcds_counts.resize(sizeof(s_tcds_trigger_types) / sizeof(const char *));
0171 
0172   // cache the L1 trigger menu
0173   histograms.l1t_counts.clear();
0174   histograms.l1t_counts.resize(GlobalAlgBlk::maxPhysicsTriggers);
0175 
0176   // initialise the HLTConfigProvider
0177   bool changed = true;
0178   edm::EDConsumerBase::Labels labels;
0179   labelsForToken(m_hlt_results, labels);
0180   if (histograms.hltConfig.init(run, setup, labels.process, changed)) {
0181     histograms.hltIndices.resize(histograms.hltConfig.size());
0182 
0183     unsigned int datasets = histograms.hltConfig.datasetNames().size();
0184     histograms.hlt_by_dataset_counts.clear();
0185     histograms.hlt_by_dataset_counts.resize(datasets);
0186 
0187     histograms.datasets.clear();
0188     histograms.datasets.resize(datasets);
0189     for (unsigned int i = 0; i < datasets; ++i) {
0190       auto const &paths = histograms.hltConfig.datasetContent(i);
0191       histograms.hlt_by_dataset_counts[i].resize(paths.size());
0192       histograms.datasets[i].reserve(paths.size());
0193       for (auto const &path : paths) {
0194         histograms.datasets[i].push_back(histograms.hltConfig.triggerIndex(path));
0195       }
0196     }
0197     histograms.dataset_counts.clear();
0198     histograms.dataset_counts.resize(datasets);
0199 
0200     unsigned int streams = histograms.hltConfig.streamNames().size();
0201     histograms.streams.clear();
0202     histograms.streams.resize(streams);
0203     for (unsigned int i = 0; i < streams; ++i) {
0204       for (auto const &dataset : histograms.hltConfig.streamContent(i)) {
0205         for (auto const &path : histograms.hltConfig.datasetContent(dataset))
0206           histograms.streams[i].push_back(histograms.hltConfig.triggerIndex(path));
0207       }
0208       std::sort(histograms.streams[i].begin(), histograms.streams[i].end());
0209       auto unique_end = std::unique(histograms.streams[i].begin(), histograms.streams[i].end());
0210       histograms.streams[i].resize(unique_end - histograms.streams[i].begin());
0211       histograms.streams[i].shrink_to_fit();
0212     }
0213     histograms.stream_counts.clear();
0214     histograms.stream_counts.resize(streams);
0215   } else {
0216     // HLTConfigProvider not initialised, skip the the HLT monitoring
0217     edm::LogError("TriggerRatesMonitor")
0218         << "failed to initialise HLTConfigProvider, the HLT trigger and datasets rates will not be monitored";
0219   }
0220 }
0221 
0222 void TriggerRatesMonitor::bookHistograms(DQMStore::IBooker &booker,
0223                                          edm::Run const &run,
0224                                          edm::EventSetup const &setup,
0225                                          RunBasedHistograms &histograms) const {
0226   // book the overall event count and event types histograms
0227   booker.setCurrentFolder(m_dqm_path);
0228   histograms.events_processed = booker.book1D(
0229       "events", "Processed events vs. lumisection", m_lumisections_range + 1, -0.5, m_lumisections_range + 0.5);
0230   booker.setCurrentFolder(m_dqm_path + "/TCDS");
0231   for (unsigned int i = 0; i < sizeof(s_tcds_trigger_types) / sizeof(const char *); ++i)
0232     if (s_tcds_trigger_types[i]) {
0233       std::string const &title = fmt::sprintf("%s events vs. lumisection", s_tcds_trigger_types[i]);
0234       histograms.tcds_counts[i] =
0235           booker.book1D(s_tcds_trigger_types[i], title, m_lumisections_range + 1, -0.5, m_lumisections_range + 0.5);
0236     }
0237 
0238   // book the rate histograms for the L1 triggers that are included in the L1 menu
0239   booker.setCurrentFolder(m_dqm_path + "/L1T");
0240   auto const &l1tMenu = setup.getData(m_l1tMenuToken);
0241   for (auto const &keyval : l1tMenu.getAlgorithmMap()) {
0242     unsigned int bit = keyval.second.getIndex();
0243     bool masked = false;  // FIXME read L1 masks once they will be avaiable in the EventSetup
0244     std::string const &name = fmt::sprintf("%s (bit %d)", keyval.first, bit);
0245     std::string const &title =
0246         fmt::sprintf("%s (bit %d)%s vs. lumisection", keyval.first, bit, (masked ? " (masked)" : ""));
0247     histograms.l1t_counts.at(bit) =
0248         booker.book1D(name, title, m_lumisections_range + 1, -0.5, m_lumisections_range + 0.5);
0249   }
0250 
0251   if (histograms.hltConfig.inited()) {
0252     auto const &datasets = histograms.hltConfig.datasetNames();
0253 
0254     // book the rate histograms for the HLT triggers
0255     for (unsigned int d = 0; d < datasets.size(); ++d) {
0256       booker.setCurrentFolder(m_dqm_path + "/HLT/" + datasets[d]);
0257       for (unsigned int i = 0; i < histograms.datasets[d].size(); ++i) {
0258         unsigned int index = histograms.datasets[d][i];
0259         std::string const &name = histograms.hltConfig.triggerName(index);
0260         histograms.hlt_by_dataset_counts[d][i].pass_l1_seed = booker.book1D(name + "_pass_L1_seed",
0261                                                                             name + " pass L1 seed, vs. lumisection",
0262                                                                             m_lumisections_range + 1,
0263                                                                             -0.5,
0264                                                                             m_lumisections_range + 0.5);
0265         histograms.hlt_by_dataset_counts[d][i].pass_prescale = booker.book1D(name + "_pass_prescaler",
0266                                                                              name + " pass prescaler, vs. lumisection",
0267                                                                              m_lumisections_range + 1,
0268                                                                              -0.5,
0269                                                                              m_lumisections_range + 0.5);
0270         histograms.hlt_by_dataset_counts[d][i].accept = booker.book1D(name + "_accept",
0271                                                                       name + " accept, vs. lumisection",
0272                                                                       m_lumisections_range + 1,
0273                                                                       -0.5,
0274                                                                       m_lumisections_range + 0.5);
0275         histograms.hlt_by_dataset_counts[d][i].reject = booker.book1D(name + "_reject",
0276                                                                       name + " reject, vs. lumisection",
0277                                                                       m_lumisections_range + 1,
0278                                                                       -0.5,
0279                                                                       m_lumisections_range + 0.5);
0280         histograms.hlt_by_dataset_counts[d][i].error = booker.book1D(name + "_error",
0281                                                                      name + " error, vs. lumisection",
0282                                                                      m_lumisections_range + 1,
0283                                                                      -0.5,
0284                                                                      m_lumisections_range + 0.5);
0285       }
0286 
0287       //      booker.setCurrentFolder( m_dqm_path + "/HLT/" + datasets[d]);
0288       for (unsigned int i : histograms.datasets[d]) {
0289         // look for the index of the (last) L1 seed and prescale module in each path
0290         histograms.hltIndices[i].index_l1_seed = histograms.hltConfig.size(i);
0291         histograms.hltIndices[i].index_prescale = histograms.hltConfig.size(i);
0292         for (unsigned int j = 0; j < histograms.hltConfig.size(i); ++j) {
0293           std::string const &label = histograms.hltConfig.moduleLabel(i, j);
0294           std::string const &type = histograms.hltConfig.moduleType(label);
0295           if (type == "HLTL1TSeed" or type == "HLTLevel1GTSeed" or type == "HLTLevel1Activity" or
0296               type == "HLTLevel1Pattern") {
0297             // there might be more L1 seed filters in sequence
0298             // keep looking and store the index of the last one
0299             histograms.hltIndices[i].index_l1_seed = j;
0300           } else if (type == "HLTPrescaler") {
0301             // there should be only one prescaler in a path, and it should follow all L1 seed filters
0302             histograms.hltIndices[i].index_prescale = j;
0303             break;
0304           }
0305         }
0306       }
0307     }
0308 
0309     // book the HLT datasets rate histograms
0310     booker.setCurrentFolder(m_dqm_path + "/Datasets");
0311     for (unsigned int i = 0; i < datasets.size(); ++i)
0312       histograms.dataset_counts[i] =
0313           booker.book1D(datasets[i], datasets[i], m_lumisections_range + 1, -0.5, m_lumisections_range + 0.5);
0314 
0315     // book the HLT streams rate histograms
0316     booker.setCurrentFolder(m_dqm_path + "/Streams");
0317     auto const &streams = histograms.hltConfig.streamNames();
0318     for (unsigned int i = 0; i < streams.size(); ++i)
0319       histograms.stream_counts[i] =
0320           booker.book1D(streams[i], streams[i], m_lumisections_range + 1, -0.5, m_lumisections_range + 0.5);
0321   }
0322 }
0323 
0324 void TriggerRatesMonitor::dqmAnalyze(edm::Event const &event,
0325                                      edm::EventSetup const &setup,
0326                                      RunBasedHistograms const &histograms) const {
0327   unsigned int lumisection = event.luminosityBlock();
0328 
0329   // monitor the overall event count and event types rates
0330   histograms.events_processed->Fill(lumisection);
0331   if (histograms.tcds_counts[event.experimentType()])
0332     histograms.tcds_counts[event.experimentType()]->Fill(lumisection);
0333 
0334   // monitor the L1 triggers rates
0335   auto const &bxvector = event.get(m_l1t_results);
0336   if (not bxvector.isEmpty(0)) {
0337     auto const &results = bxvector.at(0, 0);
0338     for (unsigned int i = 0; i < GlobalAlgBlk::maxPhysicsTriggers; ++i)
0339       if (results.getAlgoDecisionFinal(i))
0340         if (histograms.l1t_counts[i])
0341           histograms.l1t_counts[i]->Fill(lumisection);
0342   }
0343 
0344   // monitor the HLT triggers and datsets rates
0345   if (histograms.hltConfig.inited()) {
0346     auto const &hltResults = event.get(m_hlt_results);
0347     if (hltResults.size() != histograms.hltIndices.size()) {
0348       edm::LogWarning("TriggerRatesMonitor")
0349           << "This should never happen: the number of HLT paths has changed since the beginning of the run";
0350     }
0351 
0352     for (unsigned int d = 0; d < histograms.datasets.size(); ++d) {
0353       for (unsigned int i : histograms.datasets[d])
0354         if (hltResults.at(i).accept()) {
0355           histograms.dataset_counts[d]->Fill(lumisection);
0356           // ensure each dataset is incremented only once per event
0357           break;
0358         }
0359       for (unsigned int i = 0; i < histograms.datasets[d].size(); ++i) {
0360         unsigned int index = histograms.datasets[d][i];
0361         edm::HLTPathStatus const &path = hltResults.at(index);
0362 
0363         if (path.index() > histograms.hltIndices[index].index_l1_seed)
0364           histograms.hlt_by_dataset_counts[d][i].pass_l1_seed->Fill(lumisection);
0365         if (path.index() > histograms.hltIndices[index].index_prescale)
0366           histograms.hlt_by_dataset_counts[d][i].pass_prescale->Fill(lumisection);
0367         if (path.accept())
0368           histograms.hlt_by_dataset_counts[d][i].accept->Fill(lumisection);
0369         else if (path.error())
0370           histograms.hlt_by_dataset_counts[d][i].error->Fill(lumisection);
0371         else
0372           histograms.hlt_by_dataset_counts[d][i].reject->Fill(lumisection);
0373       }
0374     }
0375 
0376     for (unsigned int i = 0; i < histograms.streams.size(); ++i)
0377       for (unsigned int j : histograms.streams[i])
0378         if (hltResults.at(j).accept()) {
0379           histograms.stream_counts[i]->Fill(lumisection);
0380           // ensure each stream is incremented only once per event
0381           break;
0382         }
0383   }
0384 }
0385 
0386 //define this as a plug-in
0387 #include "FWCore/Framework/interface/MakerMacros.h"
0388 DEFINE_FWK_MODULE(TriggerRatesMonitor);