Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 // C++ headers
0002 #include <cstring>
0003 #include <iterator>
0004 #include <string>
0005 
0006 // {fmt} headers
0007 #include <fmt/printf.h>
0008 
0009 // CMSSW headers
0010 #include "FWCore/Framework/interface/Frameworkfwd.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/LuminosityBlock.h"
0013 #include "FWCore/Framework/interface/Run.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0017 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0018 #include "FWCore/ParameterSet/interface/Registry.h"
0019 #include "FWCore/ServiceRegistry/interface/Service.h"
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021 #include "DataFormats/Provenance/interface/ProcessHistory.h"
0022 #include "DataFormats/Common/interface/TriggerResults.h"
0023 #include "DataFormats/L1TGlobal/interface/GlobalAlgBlk.h"
0024 #include "CondFormats/DataRecord/interface/L1TUtmTriggerMenuRcd.h"
0025 #include "CondFormats/L1TObjects/interface/L1TUtmTriggerMenu.h"
0026 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0027 #include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h"
0028 
0029 namespace {
0030 
0031   struct RunBasedHistograms {
0032   public:
0033     typedef dqm::reco::MonitorElement MonitorElement;
0034     RunBasedHistograms()
0035         :  // L1T and HLT configuration
0036           hltConfig(),
0037           // L1T and HLT results
0038           tcds_bx_all(nullptr),
0039           l1t_bx_all(nullptr),
0040           hlt_bx_all(nullptr),
0041           tcds_bx(),
0042           l1t_bx(),
0043           hlt_bx(),
0044           tcds_bx_2d(),
0045           l1t_bx_2d(),
0046           hlt_bx_2d() {}
0047 
0048   public:
0049     // HLT configuration
0050     HLTConfigProvider hltConfig;
0051 
0052     // L1T and HLT results
0053     dqm::reco::MonitorElement* tcds_bx_all;
0054     dqm::reco::MonitorElement* l1t_bx_all;
0055     dqm::reco::MonitorElement* hlt_bx_all;
0056     std::vector<dqm::reco::MonitorElement*> tcds_bx;
0057     std::vector<dqm::reco::MonitorElement*> l1t_bx;
0058     std::vector<dqm::reco::MonitorElement*> hlt_bx;
0059     std::vector<dqm::reco::MonitorElement*> tcds_bx_2d;
0060     std::vector<dqm::reco::MonitorElement*> l1t_bx_2d;
0061     std::vector<dqm::reco::MonitorElement*> hlt_bx_2d;
0062   };
0063 
0064 }  // namespace
0065 
0066 class TriggerBxMonitor : public DQMGlobalEDAnalyzer<RunBasedHistograms> {
0067 public:
0068   explicit TriggerBxMonitor(edm::ParameterSet const&);
0069   ~TriggerBxMonitor() override = default;
0070 
0071   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0072 
0073 private:
0074   void dqmBeginRun(edm::Run const&, edm::EventSetup const&, RunBasedHistograms&) const override;
0075   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&, RunBasedHistograms&) const override;
0076   void dqmAnalyze(edm::Event const&, edm::EventSetup const&, RunBasedHistograms const&) const override;
0077 
0078   // number of bunch crossings
0079   static const unsigned int s_bx_range = 3564;
0080 
0081   // TCDS trigger types
0082   // see https://twiki.cern.ch/twiki/bin/viewauth/CMS/TcdsEventRecord
0083   static constexpr const char* s_tcds_trigger_types[] = {
0084       "Empty",          //  0 - No trigger
0085       "Physics",        //  1 - GT trigger
0086       "Calibration",    //  2 - Sequence trigger (calibration)
0087       "Random",         //  3 - Random trigger
0088       "Auxiliary",      //  4 - Auxiliary (CPM front panel NIM input) trigger
0089       nullptr,          //  5 - reserved
0090       nullptr,          //  6 - reserved
0091       nullptr,          //  7 - reserved
0092       "Cyclic",         //  8 - Cyclic trigger
0093       "Bunch-pattern",  //  9 - Bunch-pattern trigger
0094       "Software",       // 10 - Software trigger
0095       "TTS",            // 11 - TTS-sourced trigger
0096       nullptr,          // 12 - reserved
0097       nullptr,          // 13 - reserved
0098       nullptr,          // 14 - reserved
0099       nullptr           // 15 - reserved
0100   };
0101 
0102   // module configuration
0103   const edm::ESGetToken<L1TUtmTriggerMenu, L1TUtmTriggerMenuRcd> m_l1tMenuToken;
0104   const edm::EDGetTokenT<GlobalAlgBlkBxCollection> m_l1t_results;
0105   const edm::EDGetTokenT<edm::TriggerResults> m_hlt_results;
0106   const std::string m_dqm_path;
0107   const bool m_make_1d_plots;
0108   const bool m_make_2d_plots;
0109   const uint32_t m_ls_range;
0110 };
0111 
0112 // definition
0113 constexpr const char* TriggerBxMonitor::s_tcds_trigger_types[];
0114 
0115 void TriggerBxMonitor::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0116   edm::ParameterSetDescription desc;
0117   desc.addUntracked<edm::InputTag>("l1tResults", edm::InputTag("gtStage2Digis"));
0118   desc.addUntracked<edm::InputTag>("hltResults", edm::InputTag("TriggerResults"));
0119   desc.addUntracked<std::string>("dqmPath", "HLT/TriggerBx");
0120   desc.addUntracked<bool>("make1DPlots", true);
0121   desc.addUntracked<bool>("make2DPlots", false);
0122   desc.addUntracked<uint32_t>("lsRange", 4000);
0123   descriptions.add("triggerBxMonitor", desc);
0124 }
0125 
0126 TriggerBxMonitor::TriggerBxMonitor(edm::ParameterSet const& config)
0127     :  // module configuration
0128       m_l1tMenuToken{esConsumes<edm::Transition::BeginRun>()},
0129       m_l1t_results(consumes<GlobalAlgBlkBxCollection>(config.getUntrackedParameter<edm::InputTag>("l1tResults"))),
0130       m_hlt_results(consumes<edm::TriggerResults>(config.getUntrackedParameter<edm::InputTag>("hltResults"))),
0131       m_dqm_path(config.getUntrackedParameter<std::string>("dqmPath")),
0132       m_make_1d_plots(config.getUntrackedParameter<bool>("make1DPlots")),
0133       m_make_2d_plots(config.getUntrackedParameter<bool>("make2DPlots")),
0134       m_ls_range(config.getUntrackedParameter<uint32_t>("lsRange")) {}
0135 
0136 void TriggerBxMonitor::dqmBeginRun(edm::Run const& run,
0137                                    edm::EventSetup const& setup,
0138                                    RunBasedHistograms& histograms) const {
0139   // initialise the TCDS vector
0140   if (m_make_1d_plots) {
0141     histograms.tcds_bx.clear();
0142     histograms.tcds_bx.resize(std::size(s_tcds_trigger_types));
0143   }
0144   if (m_make_2d_plots) {
0145     histograms.tcds_bx_2d.clear();
0146     histograms.tcds_bx_2d.resize(std::size(s_tcds_trigger_types));
0147   }
0148 
0149   // cache the L1 trigger menu
0150   if (m_make_1d_plots) {
0151     histograms.l1t_bx.clear();
0152     histograms.l1t_bx.resize(GlobalAlgBlk::maxPhysicsTriggers);
0153   }
0154   if (m_make_2d_plots) {
0155     histograms.l1t_bx_2d.clear();
0156     histograms.l1t_bx_2d.resize(GlobalAlgBlk::maxPhysicsTriggers);
0157   }
0158 
0159   // initialise the HLTConfigProvider
0160   bool changed = true;
0161   edm::EDConsumerBase::Labels labels;
0162   labelsForToken(m_hlt_results, labels);
0163   if (histograms.hltConfig.init(run, setup, labels.process, changed)) {
0164     if (m_make_1d_plots) {
0165       histograms.hlt_bx.clear();
0166       histograms.hlt_bx.resize(histograms.hltConfig.size());
0167     }
0168     if (m_make_2d_plots) {
0169       histograms.hlt_bx_2d.clear();
0170       histograms.hlt_bx_2d.resize(histograms.hltConfig.size());
0171     }
0172   } else {
0173     // HLTConfigProvider not initialised, skip the the HLT monitoring
0174     edm::LogError("TriggerBxMonitor")
0175         << "failed to initialise HLTConfigProvider, the HLT bx distribution will not be monitored";
0176   }
0177 }
0178 
0179 void TriggerBxMonitor::bookHistograms(DQMStore::IBooker& booker,
0180                                       edm::Run const& run,
0181                                       edm::EventSetup const& setup,
0182                                       RunBasedHistograms& histograms) const {
0183   // TCDS trigger type plots
0184   {
0185     size_t size = std::size(s_tcds_trigger_types);
0186 
0187     // book 2D histogram to monitor all TCDS trigger types in a single plot
0188     booker.setCurrentFolder(m_dqm_path);
0189     histograms.tcds_bx_all = booker.book2D("TCDS Trigger Types",
0190                                            "TCDS Trigger Types vs. bunch crossing",
0191                                            s_bx_range + 1,
0192                                            -0.5,
0193                                            s_bx_range + 0.5,
0194                                            size,
0195                                            -0.5,
0196                                            size - 0.5);
0197 
0198     // book the individual histograms for the known TCDS trigger types
0199     booker.setCurrentFolder(m_dqm_path + "/TCDS");
0200     for (unsigned int i = 0; i < size; ++i) {
0201       if (s_tcds_trigger_types[i]) {
0202         if (m_make_1d_plots) {
0203           histograms.tcds_bx.at(i) =
0204               booker.book1D(s_tcds_trigger_types[i], s_tcds_trigger_types[i], s_bx_range + 1, -0.5, s_bx_range + 0.5);
0205         }
0206         if (m_make_2d_plots) {
0207           std::string const& name_ls = std::string(s_tcds_trigger_types[i]) + " vs LS";
0208           histograms.tcds_bx_2d.at(i) = booker.book2D(
0209               name_ls, name_ls, s_bx_range + 1, -0.5, s_bx_range + 0.5, m_ls_range, 0.5, m_ls_range + 0.5);
0210         }
0211         histograms.tcds_bx_all->setBinLabel(i + 1, s_tcds_trigger_types[i], 2);  // Y axis
0212       }
0213     }
0214   }
0215 
0216   // L1T plots
0217   {
0218     // book 2D histogram to monitor all L1 triggers in a single plot
0219     booker.setCurrentFolder(m_dqm_path);
0220     histograms.l1t_bx_all = booker.book2D("Level 1 Triggers",
0221                                           "Level 1 Triggers vs. bunch crossing",
0222                                           s_bx_range + 1,
0223                                           -0.5,
0224                                           s_bx_range + 0.5,
0225                                           GlobalAlgBlk::maxPhysicsTriggers,
0226                                           -0.5,
0227                                           GlobalAlgBlk::maxPhysicsTriggers - 0.5);
0228 
0229     // book the individual histograms for the L1 triggers that are included in the L1 menu
0230     booker.setCurrentFolder(m_dqm_path + "/L1T");
0231     auto const& l1tMenu = setup.getData(m_l1tMenuToken);
0232     for (auto const& keyval : l1tMenu.getAlgorithmMap()) {
0233       unsigned int bit = keyval.second.getIndex();
0234       std::string const& name = fmt::sprintf("%s (bit %d)", keyval.first, bit);
0235       if (m_make_1d_plots) {
0236         histograms.l1t_bx.at(bit) = booker.book1D(name, name, s_bx_range + 1, -0.5, s_bx_range + 0.5);
0237       }
0238       if (m_make_2d_plots) {
0239         std::string const& name_ls = name + " vs LS";
0240         histograms.l1t_bx_2d.at(bit) =
0241             booker.book2D(name_ls, name_ls, s_bx_range + 1, -0.5, s_bx_range + 0.5, m_ls_range, 0.5, m_ls_range + 0.5);
0242       }
0243       histograms.l1t_bx_all->setBinLabel(bit + 1, keyval.first, 2);  // Y axis
0244     }
0245   }
0246 
0247   // HLT plots
0248   if (histograms.hltConfig.inited()) {
0249     // book 2D histogram to monitor all HLT paths in a single plot
0250     booker.setCurrentFolder(m_dqm_path);
0251     histograms.hlt_bx_all = booker.book2D("High Level Triggers",
0252                                           "High Level Triggers vs. bunch crossing",
0253                                           s_bx_range + 1,
0254                                           -0.5,
0255                                           s_bx_range + 0.5,
0256                                           histograms.hltConfig.size(),
0257                                           -0.5,
0258                                           histograms.hltConfig.size() - 0.5);
0259 
0260     // book the individual HLT triggers histograms
0261     booker.setCurrentFolder(m_dqm_path + "/HLT");
0262     for (unsigned int i = 0; i < histograms.hltConfig.size(); ++i) {
0263       std::string const& name = histograms.hltConfig.triggerName(i);
0264       if (m_make_1d_plots) {
0265         histograms.hlt_bx[i] = booker.book1D(name, name, s_bx_range + 1, -0.5, s_bx_range + 0.5);
0266       }
0267       if (m_make_2d_plots) {
0268         std::string const& name_ls = name + " vs LS";
0269         histograms.hlt_bx_2d[i] =
0270             booker.book2D(name_ls, name_ls, s_bx_range + 1, -0.5, s_bx_range + 0.5, m_ls_range, 0.5, m_ls_range + 0.5);
0271       }
0272       histograms.hlt_bx_all->setBinLabel(i + 1, name, 2);  // Y axis
0273     }
0274   }
0275 }
0276 
0277 void TriggerBxMonitor::dqmAnalyze(edm::Event const& event,
0278                                   edm::EventSetup const& setup,
0279                                   RunBasedHistograms const& histograms) const {
0280   unsigned int bx = event.bunchCrossing();
0281   unsigned int ls = event.luminosityBlock();
0282 
0283   // monitor the bx distribution for the TCDS trigger types
0284   {
0285     size_t size = std::size(s_tcds_trigger_types);
0286     unsigned int type = event.experimentType();
0287     if (type < size) {
0288       if (m_make_1d_plots and histograms.tcds_bx.at(type))
0289         histograms.tcds_bx[type]->Fill(bx);
0290       if (m_make_2d_plots and histograms.tcds_bx_2d.at(type))
0291         histograms.tcds_bx_2d[type]->Fill(bx, ls);
0292     }
0293     histograms.tcds_bx_all->Fill(bx, type);
0294   }
0295 
0296   // monitor the bx distribution for the L1 triggers
0297   {
0298     auto const& bxvector = event.get(m_l1t_results);
0299     if (not bxvector.isEmpty(0)) {
0300       auto const& results = bxvector.at(0, 0);
0301       for (unsigned int i = 0; i < GlobalAlgBlk::maxPhysicsTriggers; ++i)
0302         if (results.getAlgoDecisionFinal(i)) {
0303           if (m_make_1d_plots and histograms.l1t_bx.at(i))
0304             histograms.l1t_bx[i]->Fill(bx);
0305           if (m_make_2d_plots and histograms.l1t_bx_2d.at(i))
0306             histograms.l1t_bx_2d[i]->Fill(bx, ls);
0307           histograms.l1t_bx_all->Fill(bx, i);
0308         }
0309     }
0310   }
0311 
0312   // monitor the bx distribution for the HLT triggers
0313   if (histograms.hltConfig.inited()) {
0314     auto const& hltResults = event.get(m_hlt_results);
0315     for (unsigned int i = 0; i < hltResults.size(); ++i) {
0316       if (hltResults.at(i).accept()) {
0317         if (m_make_1d_plots and histograms.hlt_bx.at(i))
0318           histograms.hlt_bx[i]->Fill(bx);
0319         if (m_make_2d_plots and histograms.hlt_bx_2d.at(i))
0320           histograms.hlt_bx_2d[i]->Fill(bx, ls);
0321         histograms.hlt_bx_all->Fill(bx, i);
0322       }
0323     }
0324   }
0325 }
0326 
0327 //define this as a plug-in
0328 #include "FWCore/Framework/interface/MakerMacros.h"
0329 DEFINE_FWK_MODULE(TriggerBxMonitor);