File indexing completed on 2023-03-17 11:12:09
0001 #include <array>
0002 #include <cassert>
0003 #include <cstring>
0004 #include <memory>
0005 #include <stdexcept>
0006 #include <string>
0007 #include <vector>
0008
0009 #include <fmt/printf.h>
0010
0011 using namespace std::literals;
0012
0013 namespace {
0014
0015 template <class T>
0016 struct Entry {
0017 T value;
0018 const char* tag;
0019 const char* description;
0020 };
0021
0022 template <typename S, typename T, unsigned int N>
0023 std::string build_comment_from_entries(S pre, const Entry<T> (&entries)[N]) {
0024 std::string comment{pre};
0025 size_t length = 0;
0026 for (auto entry : entries)
0027 if (entry.tag)
0028 length = std::max(std::strlen(entry.tag), length);
0029 for (auto entry : entries)
0030 if (entry.tag) {
0031 comment.reserve(comment.size() + length + std::strlen(entry.description) + 8);
0032 comment += "\n \"";
0033 comment += entry.tag;
0034 comment += "\": ";
0035 for (unsigned int i = 0; i < length - std::strlen(entry.tag); ++i)
0036 comment += ' ';
0037 comment += entry.description;
0038 }
0039 return comment;
0040 }
0041
0042 template <typename S1, typename S2, typename T, unsigned int N>
0043 std::string build_comment_from_entries(S1 pre, const Entry<T> (&entries)[N], S2 post) {
0044 std::string comment = build_comment_from_entries(pre, entries);
0045 comment += '\n';
0046 comment += post;
0047 return comment;
0048 }
0049
0050 template <class T>
0051 constexpr T get_enum_value(Entry<T> const* entries, const char* tag) {
0052 for (; entries->tag; ++entries)
0053 if (std::strcmp(entries->tag, tag) == 0)
0054 return entries->value;
0055 throw std::logic_error("invalid tag "s + tag);
0056 }
0057
0058 template <class T>
0059 constexpr T get_enum_value(Entry<T> const* entries, const char* tag, T default_value) {
0060 for (; entries->tag; ++entries)
0061 if (std::strcmp(entries->tag, tag) == 0)
0062 return entries->value;
0063 return default_value;
0064 }
0065
0066 }
0067
0068
0069
0070 #include "FWCore/Framework/interface/Frameworkfwd.h"
0071 #include "FWCore/Framework/interface/Event.h"
0072 #include "FWCore/Framework/interface/one/EDFilter.h"
0073 #include "FWCore/Framework/interface/EventSetup.h"
0074 #include "FWCore/Framework/interface/ESHandle.h"
0075 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0076 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0077 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0078 #include "FWCore/Utilities/interface/EDMException.h"
0079 #include "FWCore/Utilities/interface/ESGetToken.h"
0080 #include "CondFormats/DataRecord/interface/L1TGlobalPrescalesVetosRcd.h"
0081 #include "CondFormats/L1TObjects/interface/L1TGlobalPrescalesVetos.h"
0082 #include "DataFormats/L1TGlobal/interface/GlobalAlgBlk.h"
0083
0084 class L1TGlobalPrescaler : public edm::one::EDFilter<> {
0085 public:
0086 L1TGlobalPrescaler(edm::ParameterSet const& config);
0087
0088 bool filter(edm::Event& event, edm::EventSetup const& setup) override;
0089
0090 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0091
0092 private:
0093 enum class Mode {
0094 ApplyPrescaleValues,
0095 ApplyPrescaleRatios,
0096 ApplyColumnValues,
0097 ApplyColumnRatios,
0098 ForcePrescaleValues,
0099 ForceColumnValues,
0100 Invalid = -1
0101 };
0102
0103 static const constexpr Entry<Mode> s_modes[]{
0104 {Mode::ApplyPrescaleValues, "applyPrescaleValues", "apply the given prescale values"},
0105 {Mode::ApplyPrescaleRatios,
0106 "applyPrescaleRatios",
0107 "apply prescales equal to ratio between the given values and the ones read from the EventSetup"},
0108 {Mode::ApplyColumnValues,
0109 "applyColumnValues",
0110 "apply the prescale values from the EventSetup corresponding to the given column index"},
0111 {Mode::ApplyColumnRatios,
0112 "applyColumnRatios",
0113 "apply prescales equal to ratio between the values corresponsing to the given column index, and the ones read "
0114 "from the EventSetup"},
0115 {Mode::ForcePrescaleValues,
0116 "forcePrescaleValues",
0117 "apply the given prescale values, ignoring the prescales and masks already applied"},
0118 {Mode::ForceColumnValues,
0119 "forceColumnValues",
0120 "apply the prescale values from the EventSetup corresponding to the given column index, ignoring the prescales "
0121 "and masks already applied"},
0122 {Mode::Invalid, nullptr, nullptr}};
0123
0124 const Mode m_mode;
0125 const edm::EDGetTokenT<GlobalAlgBlkBxCollection> m_l1tResultsToken;
0126 const std::array<double, GlobalAlgBlk::maxPhysicsTriggers> m_l1tPrescales;
0127 std::array<double, GlobalAlgBlk::maxPhysicsTriggers> m_prescales;
0128 std::array<unsigned int, GlobalAlgBlk::maxPhysicsTriggers> m_counters;
0129 const int m_l1tPrescaleColumn;
0130 int m_oldIndex;
0131 edm::ESGetToken<L1TGlobalPrescalesVetos, L1TGlobalPrescalesVetosRcd> m_l1tGtPrescalesVetosToken;
0132 };
0133
0134 const constexpr Entry<L1TGlobalPrescaler::Mode> L1TGlobalPrescaler::s_modes[];
0135
0136 L1TGlobalPrescaler::L1TGlobalPrescaler(edm::ParameterSet const& config)
0137 : m_mode(get_enum_value(s_modes, config.getParameter<std::string>("mode").c_str(), Mode::Invalid)),
0138 m_l1tResultsToken(consumes<GlobalAlgBlkBxCollection>(config.getParameter<edm::InputTag>("l1tResults"))),
0139 m_l1tPrescales(m_mode == Mode::ApplyPrescaleValues or m_mode == Mode::ApplyPrescaleRatios or
0140 m_mode == Mode::ForcePrescaleValues
0141 ? config.getParameter<std::array<double, GlobalAlgBlk::maxPhysicsTriggers>>("l1tPrescales")
0142 : std::array<double, GlobalAlgBlk::maxPhysicsTriggers>{}),
0143 m_l1tPrescaleColumn(m_mode == Mode::ApplyColumnValues or m_mode == Mode::ApplyColumnRatios or
0144 m_mode == Mode::ForceColumnValues
0145 ? config.getParameter<uint32_t>("l1tPrescaleColumn")
0146 : 0),
0147 m_oldIndex(-1) {
0148 switch (m_mode) {
0149
0150 case Mode::ApplyPrescaleValues:
0151 case Mode::ForcePrescaleValues:
0152 m_prescales = m_l1tPrescales;
0153 break;
0154
0155
0156 case Mode::ApplyColumnValues:
0157 case Mode::ApplyPrescaleRatios:
0158 case Mode::ApplyColumnRatios:
0159 case Mode::ForceColumnValues:
0160 m_l1tGtPrescalesVetosToken = esConsumes<L1TGlobalPrescalesVetos, L1TGlobalPrescalesVetosRcd>();
0161 break;
0162
0163
0164 case Mode::Invalid:
0165 throw edm::Exception(edm::errors::Configuration)
0166 << "invalid mode \"" << config.getParameter<std::string>("mode") << "\"";
0167 }
0168
0169 m_counters.fill(0);
0170 produces<GlobalAlgBlkBxCollection>();
0171 }
0172
0173 bool L1TGlobalPrescaler::filter(edm::Event& event, edm::EventSetup const& setup) {
0174 edm::Handle<GlobalAlgBlkBxCollection> handle;
0175 event.getByToken(m_l1tResultsToken, handle);
0176
0177
0178
0179 if (handle->isEmpty(0)) {
0180 std::unique_ptr<GlobalAlgBlkBxCollection> result(new GlobalAlgBlkBxCollection());
0181 event.put(std::move(result));
0182 return false;
0183 }
0184
0185
0186 int index = handle->at(0, 0).getPreScColumn();
0187 assert(index >= 0);
0188
0189
0190
0191 if (m_mode == Mode::ApplyPrescaleRatios and m_oldIndex != index) {
0192 edm::ESHandle<L1TGlobalPrescalesVetos> h = setup.getHandle(m_l1tGtPrescalesVetosToken);
0193
0194 auto const& prescaleTable = h->prescale_table_;
0195 if (index >= (int)prescaleTable.size())
0196 throw edm::Exception(edm::errors::LogicError)
0197 << fmt::sprintf("The prescale index %d is invalid, it should be smaller than the prescale table size %d.",
0198 index,
0199 prescaleTable.size());
0200 auto const& prescales = prescaleTable[index];
0201 unsigned long i = 0;
0202 for (; i < std::min(prescales.size(), (unsigned long)GlobalAlgBlk::maxPhysicsTriggers); ++i)
0203 if (m_l1tPrescales[i] == 0) {
0204
0205 m_prescales[i] = 0.;
0206 } else if (prescales[i] == 0) {
0207
0208 m_prescales[i] = 0.;
0209 edm::LogWarning("L1TGlobalPrescaler")
0210 << "Request to enable the trigger " << i << " which was originally disabled\nIt will be kept disabled.";
0211 } else if (m_l1tPrescales[i] < prescales[i]) {
0212
0213 m_prescales[i] = 1.;
0214 edm::LogWarning("L1TGlobalPrescaler")
0215 << "Request to prescale the trigger " << i
0216 << " less than it was originally prescaled\nNo further prescale will be applied.";
0217 } else {
0218
0219 m_prescales[i] = (double)m_l1tPrescales[i] / prescales[i];
0220 }
0221 for (; i < (unsigned long)GlobalAlgBlk::maxPhysicsTriggers; ++i)
0222
0223 m_prescales[i] = 0.;
0224
0225 m_counters.fill(0);
0226 m_oldIndex = index;
0227 }
0228
0229
0230
0231 if ((m_mode == Mode::ApplyColumnValues or m_mode == Mode::ForceColumnValues) and m_oldIndex != m_l1tPrescaleColumn) {
0232 edm::ESHandle<L1TGlobalPrescalesVetos> h = setup.getHandle(m_l1tGtPrescalesVetosToken);
0233 auto const& prescaleTable = h->prescale_table_;
0234 if (m_l1tPrescaleColumn >= (int)prescaleTable.size())
0235 throw edm::Exception(edm::errors::Configuration)
0236 << fmt::sprintf("The prescale index %d is invalid, it should be smaller than the prescale table size %d.",
0237 m_l1tPrescaleColumn,
0238 prescaleTable.size());
0239 auto const& targets = prescaleTable[m_l1tPrescaleColumn];
0240 unsigned long i = 0;
0241 for (; i < std::min(targets.size(), (unsigned long)GlobalAlgBlk::maxPhysicsTriggers); ++i)
0242
0243 m_prescales[i] = targets[i];
0244 for (; i < (unsigned long)GlobalAlgBlk::maxPhysicsTriggers; ++i)
0245
0246 m_prescales[i] = 0.;
0247
0248 m_counters.fill(0);
0249 m_oldIndex = m_l1tPrescaleColumn;
0250 }
0251
0252
0253
0254 if (m_mode == Mode::ApplyColumnRatios and m_oldIndex != index) {
0255 edm::ESHandle<L1TGlobalPrescalesVetos> h = setup.getHandle(m_l1tGtPrescalesVetosToken);
0256 auto const& prescaleTable = h->prescale_table_;
0257 if (index >= (int)prescaleTable.size())
0258 throw edm::Exception(edm::errors::LogicError)
0259 << fmt::sprintf("The prescale index %d is invalid, it should be smaller than the prescale table size %d.",
0260 index,
0261 prescaleTable.size());
0262 if (m_l1tPrescaleColumn >= (int)prescaleTable.size())
0263 throw edm::Exception(edm::errors::Configuration)
0264 << fmt::sprintf("The prescale index %d is invalid, it should be smaller than the prescale table size %d.",
0265 m_l1tPrescaleColumn,
0266 prescaleTable.size());
0267 auto const& prescales = prescaleTable[index];
0268 auto const& targets = prescaleTable[m_l1tPrescaleColumn];
0269 unsigned long i = 0;
0270 for (; i < std::min({prescales.size(), targets.size(), (unsigned long)GlobalAlgBlk::maxPhysicsTriggers}); ++i)
0271 if (prescales[i] == 0)
0272
0273 m_prescales[i] = 0.;
0274 else
0275
0276 m_prescales[i] = targets[i] < prescales[i] ? 1. : (double)targets[i] / prescales[i];
0277 for (; i < (unsigned long)GlobalAlgBlk::maxPhysicsTriggers; ++i)
0278
0279 m_prescales[i] = 0.;
0280
0281 m_counters.fill(0);
0282 m_oldIndex = index;
0283 }
0284
0285
0286 GlobalAlgBlk algoBlock = handle->at(0, 0);
0287
0288 bool finalOr = false;
0289 std::vector<bool> const& decision = (m_mode == Mode::ForceColumnValues or m_mode == Mode::ForcePrescaleValues)
0290 ? algoBlock.getAlgoDecisionInitial()
0291 : algoBlock.getAlgoDecisionFinal();
0292
0293 for (unsigned int i = 0; i < GlobalAlgBlk::maxPhysicsTriggers; ++i) {
0294 if (m_prescales[i] == 0) {
0295
0296 algoBlock.setAlgoDecisionFinal(i, false);
0297 } else if (decision[i]) {
0298
0299 ++m_counters[i];
0300 if (std::fmod(m_counters[i], m_prescales[i]) < 1) {
0301
0302 algoBlock.setAlgoDecisionFinal(i, true);
0303 finalOr = true;
0304 } else {
0305
0306 algoBlock.setAlgoDecisionFinal(i, false);
0307 }
0308 }
0309 }
0310
0311
0312 algoBlock.setFinalORPreVeto(finalOr);
0313 if (algoBlock.getFinalORVeto())
0314 finalOr = false;
0315 algoBlock.setFinalOR(finalOr);
0316
0317
0318 if (m_mode == Mode::ApplyColumnValues or m_mode == Mode::ApplyColumnRatios or m_mode == Mode::ForceColumnValues)
0319 algoBlock.setPreScColumn(m_l1tPrescaleColumn);
0320
0321
0322 std::unique_ptr<GlobalAlgBlkBxCollection> result(new GlobalAlgBlkBxCollection());
0323 result->push_back(0, algoBlock);
0324 event.put(std::move(result));
0325
0326 return finalOr;
0327 }
0328
0329 void L1TGlobalPrescaler::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0330
0331 edm::ParameterDescription<edm::InputTag> l1tResults("l1tResults", edm::InputTag("gtStage2Digis"), true);
0332 l1tResults.setComment("Collection with the original uGT results");
0333
0334
0335 edm::ParameterDescription<std::string> mode("mode", "applyPrescaleValues", true);
0336 mode.setComment(build_comment_from_entries("Define how to apply the prescale values:", s_modes));
0337
0338
0339 edm::ParameterDescription<std::vector<double>> l1tPrescales(
0340 "l1tPrescales", std::vector<double>(GlobalAlgBlk::maxPhysicsTriggers, 1.), true);
0341 l1tPrescales.setComment(
0342 "Target prescale values (for modes \"applyPrescaleValues\", \"applyPrescaleRatios\" or \"forcePrescaleValues\")");
0343
0344
0345 edm::ParameterDescription<uint32_t> l1tPrescaleColumn("l1tPrescaleColumn", 0, true);
0346 l1tPrescaleColumn.setComment(
0347 "Target prescale column (for modes \"applyColumnValues\", \"applyColumnRatios\" or \"forceColumnValues\")");
0348
0349
0350 {
0351 edm::ParameterSetDescription desc;
0352 desc.addNode(l1tResults);
0353 desc.ifValue(
0354 mode,
0355
0356 "applyPrescaleValues" >> l1tPrescales or "applyPrescaleRatios" >> l1tPrescales or
0357 "forcePrescaleValues" >> l1tPrescales or
0358
0359 "applyColumnValues" >> l1tPrescaleColumn or "applyColumnRatios" >> l1tPrescaleColumn or
0360 "forceColumnValues" >> l1tPrescaleColumn);
0361 descriptions.add("l1tGlobalPrescaler", desc);
0362 }
0363
0364
0365 {
0366 edm::ParameterSetDescription desc;
0367 desc.addNode(l1tResults);
0368 desc.add<std::string>("mode", "applyColumnRatios")
0369 ->setComment(
0370 "apply prescales equal to ratio between the values corresponsing to the given column index, and the ones "
0371 "read from the EventSetup");
0372 desc.addNode(l1tPrescaleColumn);
0373 descriptions.add("l1tGlobalPrescalerTargetColumn", desc);
0374 }
0375 }
0376
0377
0378 #include "FWCore/Framework/interface/MakerMacros.h"
0379 DEFINE_FWK_MODULE(L1TGlobalPrescaler);