File indexing completed on 2024-04-06 12:07:38
0001
0002 #include <cstring>
0003 #include <iterator>
0004 #include <string>
0005
0006
0007 #include <fmt/printf.h>
0008
0009
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 :
0036 hltConfig(),
0037
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
0050 HLTConfigProvider hltConfig;
0051
0052
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 }
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
0079 static const unsigned int s_bx_range = 3564;
0080
0081
0082
0083 static constexpr const char* s_tcds_trigger_types[] = {
0084 "Empty",
0085 "Physics",
0086 "Calibration",
0087 "Random",
0088 "Auxiliary",
0089 nullptr,
0090 nullptr,
0091 nullptr,
0092 "Cyclic",
0093 "Bunch-pattern",
0094 "Software",
0095 "TTS",
0096 nullptr,
0097 nullptr,
0098 nullptr,
0099 nullptr
0100 };
0101
0102
0103 const edm::ESGetToken<L1TUtmTriggerMenu, L1TUtmTriggerMenuRcd> m_l1tMenu_token;
0104 const edm::InputTag m_l1t_results_inputTag;
0105 const edm::EDGetTokenT<GlobalAlgBlkBxCollection> m_l1t_results_token;
0106 const edm::EDGetTokenT<edm::TriggerResults> m_hlt_results_token;
0107 const std::string m_dqm_path;
0108 const bool m_make_1d_plots;
0109 const bool m_make_2d_plots;
0110 const uint32_t m_ls_range;
0111 };
0112
0113
0114 constexpr const char* TriggerBxMonitor::s_tcds_trigger_types[];
0115
0116 void TriggerBxMonitor::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0117 edm::ParameterSetDescription desc;
0118 desc.addUntracked<edm::InputTag>("l1tResults", edm::InputTag("gtStage2Digis"));
0119 desc.addUntracked<edm::InputTag>("hltResults", edm::InputTag("TriggerResults"));
0120 desc.addUntracked<std::string>("dqmPath", "HLT/TriggerBx");
0121 desc.addUntracked<bool>("make1DPlots", true);
0122 desc.addUntracked<bool>("make2DPlots", false);
0123 desc.addUntracked<uint32_t>("lsRange", 4000);
0124 descriptions.add("triggerBxMonitor", desc);
0125 }
0126
0127 TriggerBxMonitor::TriggerBxMonitor(edm::ParameterSet const& config)
0128 :
0129 m_l1tMenu_token{esConsumes<edm::Transition::BeginRun>()},
0130 m_l1t_results_inputTag{config.getUntrackedParameter<edm::InputTag>("l1tResults")},
0131 m_l1t_results_token{consumes(m_l1t_results_inputTag)},
0132 m_hlt_results_token{consumes(config.getUntrackedParameter<edm::InputTag>("hltResults"))},
0133 m_dqm_path{config.getUntrackedParameter<std::string>("dqmPath")},
0134 m_make_1d_plots{config.getUntrackedParameter<bool>("make1DPlots")},
0135 m_make_2d_plots{config.getUntrackedParameter<bool>("make2DPlots")},
0136 m_ls_range{config.getUntrackedParameter<uint32_t>("lsRange")} {}
0137
0138 void TriggerBxMonitor::dqmBeginRun(edm::Run const& run,
0139 edm::EventSetup const& setup,
0140 RunBasedHistograms& histograms) const {
0141
0142 if (m_make_1d_plots) {
0143 histograms.tcds_bx.clear();
0144 histograms.tcds_bx.resize(std::size(s_tcds_trigger_types));
0145 }
0146 if (m_make_2d_plots) {
0147 histograms.tcds_bx_2d.clear();
0148 histograms.tcds_bx_2d.resize(std::size(s_tcds_trigger_types));
0149 }
0150
0151
0152 if (m_make_1d_plots) {
0153 histograms.l1t_bx.clear();
0154 histograms.l1t_bx.resize(GlobalAlgBlk::maxPhysicsTriggers);
0155 }
0156 if (m_make_2d_plots) {
0157 histograms.l1t_bx_2d.clear();
0158 histograms.l1t_bx_2d.resize(GlobalAlgBlk::maxPhysicsTriggers);
0159 }
0160
0161
0162 bool changed = true;
0163 edm::EDConsumerBase::Labels labels;
0164 labelsForToken(m_hlt_results_token, labels);
0165 if (histograms.hltConfig.init(run, setup, labels.process, changed)) {
0166 if (m_make_1d_plots) {
0167 histograms.hlt_bx.clear();
0168 histograms.hlt_bx.resize(histograms.hltConfig.size());
0169 }
0170 if (m_make_2d_plots) {
0171 histograms.hlt_bx_2d.clear();
0172 histograms.hlt_bx_2d.resize(histograms.hltConfig.size());
0173 }
0174 } else {
0175
0176 edm::LogError("TriggerBxMonitor")
0177 << "failed to initialise HLTConfigProvider, the HLT bx distribution will not be monitored";
0178 }
0179 }
0180
0181 void TriggerBxMonitor::bookHistograms(DQMStore::IBooker& booker,
0182 edm::Run const& run,
0183 edm::EventSetup const& setup,
0184 RunBasedHistograms& histograms) const {
0185
0186 {
0187 size_t size = std::size(s_tcds_trigger_types);
0188
0189
0190 booker.setCurrentFolder(m_dqm_path);
0191 histograms.tcds_bx_all = booker.book2D("TCDS Trigger Types",
0192 "TCDS Trigger Types vs. bunch crossing",
0193 s_bx_range + 1,
0194 -0.5,
0195 s_bx_range + 0.5,
0196 size,
0197 -0.5,
0198 size - 0.5);
0199
0200
0201 booker.setCurrentFolder(m_dqm_path + "/TCDS");
0202 for (unsigned int i = 0; i < size; ++i) {
0203 if (s_tcds_trigger_types[i]) {
0204 if (m_make_1d_plots) {
0205 histograms.tcds_bx.at(i) =
0206 booker.book1D(s_tcds_trigger_types[i], s_tcds_trigger_types[i], s_bx_range + 1, -0.5, s_bx_range + 0.5);
0207 }
0208 if (m_make_2d_plots) {
0209 std::string const& name_ls = std::string(s_tcds_trigger_types[i]) + " vs LS";
0210 histograms.tcds_bx_2d.at(i) = booker.book2D(
0211 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);
0212 }
0213 histograms.tcds_bx_all->setBinLabel(i + 1, s_tcds_trigger_types[i], 2);
0214 }
0215 }
0216 }
0217
0218
0219 {
0220
0221 booker.setCurrentFolder(m_dqm_path);
0222 histograms.l1t_bx_all = booker.book2D("Level 1 Triggers",
0223 "Level 1 Triggers vs. bunch crossing",
0224 s_bx_range + 1,
0225 -0.5,
0226 s_bx_range + 0.5,
0227 GlobalAlgBlk::maxPhysicsTriggers,
0228 -0.5,
0229 GlobalAlgBlk::maxPhysicsTriggers - 0.5);
0230
0231
0232 booker.setCurrentFolder(m_dqm_path + "/L1T");
0233 auto const& l1tMenu = setup.getData(m_l1tMenu_token);
0234 for (auto const& keyval : l1tMenu.getAlgorithmMap()) {
0235 unsigned int bit = keyval.second.getIndex();
0236 std::string const& name = fmt::sprintf("%s (bit %d)", keyval.first, bit);
0237 if (m_make_1d_plots) {
0238 histograms.l1t_bx.at(bit) = booker.book1D(name, name, s_bx_range + 1, -0.5, s_bx_range + 0.5);
0239 }
0240 if (m_make_2d_plots) {
0241 std::string const& name_ls = name + " vs LS";
0242 histograms.l1t_bx_2d.at(bit) =
0243 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);
0244 }
0245 histograms.l1t_bx_all->setBinLabel(bit + 1, keyval.first, 2);
0246 }
0247 }
0248
0249
0250 if (histograms.hltConfig.inited()) {
0251
0252 booker.setCurrentFolder(m_dqm_path);
0253 histograms.hlt_bx_all = booker.book2D("High Level Triggers",
0254 "High Level Triggers vs. bunch crossing",
0255 s_bx_range + 1,
0256 -0.5,
0257 s_bx_range + 0.5,
0258 histograms.hltConfig.size(),
0259 -0.5,
0260 histograms.hltConfig.size() - 0.5);
0261
0262
0263 booker.setCurrentFolder(m_dqm_path + "/HLT");
0264 for (unsigned int i = 0; i < histograms.hltConfig.size(); ++i) {
0265 std::string const& name = histograms.hltConfig.triggerName(i);
0266 if (m_make_1d_plots) {
0267 histograms.hlt_bx[i] = booker.book1D(name, name, s_bx_range + 1, -0.5, s_bx_range + 0.5);
0268 }
0269 if (m_make_2d_plots) {
0270 std::string const& name_ls = name + " vs LS";
0271 histograms.hlt_bx_2d[i] =
0272 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);
0273 }
0274 histograms.hlt_bx_all->setBinLabel(i + 1, name, 2);
0275 }
0276 }
0277 }
0278
0279 void TriggerBxMonitor::dqmAnalyze(edm::Event const& event,
0280 edm::EventSetup const& setup,
0281 RunBasedHistograms const& histograms) const {
0282 unsigned int bx = event.bunchCrossing();
0283 unsigned int ls = event.luminosityBlock();
0284
0285
0286 {
0287 size_t size = std::size(s_tcds_trigger_types);
0288 unsigned int type = event.experimentType();
0289 if (type < size) {
0290 if (m_make_1d_plots and histograms.tcds_bx.at(type))
0291 histograms.tcds_bx[type]->Fill(bx);
0292 if (m_make_2d_plots and histograms.tcds_bx_2d.at(type))
0293 histograms.tcds_bx_2d[type]->Fill(bx, ls);
0294 }
0295 histograms.tcds_bx_all->Fill(bx, type);
0296 }
0297
0298
0299 {
0300 auto const& algBlkBxVecHandle = event.getHandle(m_l1t_results_token);
0301 if (not algBlkBxVecHandle.isValid()) {
0302 edm::LogError("TriggerBxMonitor")
0303 << "L1 trigger results with label [" << m_l1t_results_inputTag.encode()
0304 << "] not present or invalid. MonitorElements of L1T results not filled for this event.";
0305 } else if (algBlkBxVecHandle->isEmpty(0)) {
0306 edm::LogError("TriggerBxMonitor")
0307 << "L1 trigger results with label [" << m_l1t_results_inputTag.encode()
0308 << "] empty for BX=0. MonitorElements of L1T results not filled for this event.";
0309 } else {
0310 auto const& results = algBlkBxVecHandle->at(0, 0);
0311 for (unsigned int i = 0; i < GlobalAlgBlk::maxPhysicsTriggers; ++i) {
0312 if (results.getAlgoDecisionFinal(i)) {
0313 if (m_make_1d_plots and histograms.l1t_bx.at(i))
0314 histograms.l1t_bx[i]->Fill(bx);
0315 if (m_make_2d_plots and histograms.l1t_bx_2d.at(i))
0316 histograms.l1t_bx_2d[i]->Fill(bx, ls);
0317 histograms.l1t_bx_all->Fill(bx, i);
0318 }
0319 }
0320 }
0321 }
0322
0323
0324 if (histograms.hltConfig.inited()) {
0325 auto const& hltResults = event.get(m_hlt_results_token);
0326 for (unsigned int i = 0; i < hltResults.size(); ++i) {
0327 if (hltResults.at(i).accept()) {
0328 if (m_make_1d_plots and histograms.hlt_bx.at(i))
0329 histograms.hlt_bx[i]->Fill(bx);
0330 if (m_make_2d_plots and histograms.hlt_bx_2d.at(i))
0331 histograms.hlt_bx_2d[i]->Fill(bx, ls);
0332 histograms.hlt_bx_all->Fill(bx, i);
0333 }
0334 }
0335 }
0336 }
0337
0338
0339 #include "FWCore/Framework/interface/MakerMacros.h"
0340 DEFINE_FWK_MODULE(TriggerBxMonitor);