File indexing completed on 2024-04-06 12:18:25
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/Framework/interface/LuminosityBlock.h"
0011 #include "FWCore/ServiceRegistry/interface/PathContext.h"
0012 #include "FWCore/ServiceRegistry/interface/PlaceInPathContext.h"
0013 #include "FWCore/ServiceRegistry/interface/ModuleCallingContext.h"
0014 #include "FWCore/ServiceRegistry/interface/Service.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0017 #include "HLTrigger/HLTcore/interface/HLTPrescaler.h"
0018
0019
0020
0021
0022
0023 const unsigned int HLTPrescaler::prescaleSeed_ = 65537;
0024
0025
0026
0027
0028
0029
0030 HLTPrescaler::HLTPrescaler(edm::ParameterSet const& iConfig, const trigger::Efficiency* efficiency)
0031 : prescaleSet_(0),
0032 prescaleFactor_(1),
0033 eventCount_(0),
0034 acceptCount_(0),
0035 offsetCount_(0),
0036 offsetPhase_(iConfig.getParameter<unsigned int>("offset")),
0037 prescaleService_(nullptr),
0038 newLumi_(true),
0039 gtDigiTag_(iConfig.getParameter<edm::InputTag>("L1GtReadoutRecordTag")),
0040 gtDigi1Token_(consumes<L1GlobalTriggerReadoutRecord>(gtDigiTag_)),
0041 gtDigi2Token_(consumes<GlobalAlgBlkBxCollection>(gtDigiTag_)) {
0042 if (edm::Service<edm::service::PrescaleService>().isAvailable())
0043 prescaleService_ = edm::Service<edm::service::PrescaleService>().operator->();
0044 else
0045 LogDebug("NoPrescaleService") << "PrescaleService unavailable, prescaleFactor=1!";
0046 }
0047
0048
0049 HLTPrescaler::~HLTPrescaler() = default;
0050
0051
0052
0053
0054
0055 void HLTPrescaler::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0056 edm::ParameterSetDescription desc;
0057 desc.add<unsigned int>("offset", 0);
0058 desc.add<edm::InputTag>("L1GtReadoutRecordTag", edm::InputTag("hltGtStage2Digis"));
0059 descriptions.add("hltPrescaler", desc);
0060 }
0061
0062
0063 void HLTPrescaler::beginLuminosityBlock(edm::LuminosityBlock const& lb, edm::EventSetup const& iSetup) {
0064 newLumi_ = true;
0065 }
0066
0067
0068 bool HLTPrescaler::filter(edm::Event& iEvent, const edm::EventSetup&) {
0069
0070
0071 if (newLumi_) {
0072 newLumi_ = false;
0073
0074 bool needsInit(eventCount_ == 0);
0075
0076 if (prescaleService_) {
0077 std::string const& pathName = iEvent.moduleCallingContext()->placeInPathContext()->pathContext()->pathName();
0078 const unsigned int oldSet(prescaleSet_);
0079 const unsigned int oldPrescale(prescaleFactor_);
0080
0081 edm::Handle<GlobalAlgBlkBxCollection> handle2;
0082 iEvent.getByToken(gtDigi2Token_, handle2);
0083 if (handle2.isValid()) {
0084 if (not handle2->isEmpty(0)) {
0085 prescaleSet_ = static_cast<unsigned int>(handle2->begin(0)->getPreScColumn());
0086 prescaleFactor_ = prescaleService_->getPrescale(prescaleSet_, pathName);
0087 } else {
0088 edm::LogWarning("HLT")
0089 << "Cannot read prescale column index from GT2 data: using default as defined by configuration or DAQ";
0090 prescaleFactor_ = prescaleService_->getPrescale(pathName);
0091 }
0092 } else {
0093 edm::Handle<L1GlobalTriggerReadoutRecord> handle1;
0094 iEvent.getByToken(gtDigi1Token_, handle1);
0095 if (handle1.isValid()) {
0096 prescaleSet_ = handle1->gtFdlWord().gtPrescaleFactorIndexAlgo();
0097
0098
0099 prescaleFactor_ = prescaleService_->getPrescale(prescaleSet_, pathName);
0100 } else {
0101 edm::LogWarning("HLT")
0102 << "Cannot read prescale column index from GT1 data: using default as defined by configuration or DAQ";
0103 prescaleFactor_ = prescaleService_->getPrescale(pathName);
0104 }
0105 }
0106
0107 if (prescaleSet_ != oldSet) {
0108 edm::LogInfo("ChangedPrescale") << "lumiBlockNb = " << iEvent.getLuminosityBlock().id().luminosityBlock()
0109 << ", set = " << prescaleSet_ << " [" << oldSet << "]"
0110 << ", path = " << pathName << ": " << prescaleFactor_ << " [" << oldPrescale
0111 << "]";
0112
0113 needsInit = true;
0114 }
0115 }
0116
0117 if (needsInit && (prescaleFactor_ != 0)) {
0118
0119 offsetCount_ = ((uint64_t)(iEvent.id().event() + offsetPhase_) * prescaleSeed_) % prescaleFactor_;
0120 }
0121 }
0122
0123 const bool result((prescaleFactor_ == 0) ? false : ((eventCount_ + offsetCount_) % prescaleFactor_ == 0));
0124
0125 ++eventCount_;
0126 if (result)
0127 ++acceptCount_;
0128 return result;
0129 }
0130
0131
0132 void HLTPrescaler::endStream() {
0133
0134
0135 globalCache()->eventCount_ += eventCount_;
0136 globalCache()->acceptCount_ += acceptCount_;
0137 return;
0138 }
0139
0140
0141 void HLTPrescaler::globalEndJob(const trigger::Efficiency* efficiency) {
0142 unsigned int accept(efficiency->acceptCount_);
0143 unsigned int event(efficiency->eventCount_);
0144 edm::LogInfo("PrescaleSummary") << accept << "/" << event << " ("
0145 << 100. * accept / static_cast<double>(std::max(1u, event))
0146 << "% of events accepted).";
0147 return;
0148 }
0149
0150 #include "FWCore/Framework/interface/MakerMacros.h"
0151 DEFINE_FWK_MODULE(HLTPrescaler);