File indexing completed on 2021-02-14 14:20:45
0001
0002 #include <algorithm>
0003 #include <chrono>
0004 #include <ctime>
0005
0006
0007 #include <fmt/printf.h>
0008
0009
0010 #include "DQMServices/Core/interface/DQMStore.h"
0011 #include "FWCore/ParameterSet/interface/EmptyGroupDescription.h"
0012 #include "FWCore/Utilities/interface/TimeOfDay.h"
0013 #include "HLTrigger/Timer/interface/processor_model.h"
0014 #include "ThroughputService.h"
0015
0016
0017 void ThroughputService::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0018 edm::ParameterSetDescription desc;
0019 desc.addUntracked<uint32_t>("eventRange", 10000)->setComment("Preallocate a buffer for N events");
0020 desc.addUntracked<uint32_t>("eventResolution", 1)->setComment("Sample the processing time every N events");
0021 desc.addUntracked<bool>("printEventSummary", false);
0022 desc.ifValue(edm::ParameterDescription<bool>("enableDQM", true, false),
0023
0024 true >> (edm::ParameterDescription<bool>("dqmPathByProcesses", false, false) and
0025 edm::ParameterDescription<std::string>("dqmPath", "HLT/Throughput", false) and
0026 edm::ParameterDescription<double>("timeRange", 60000.0, false) and
0027 edm::ParameterDescription<double>("timeResolution", 10.0, false)) or
0028
0029 false >> edm::EmptyGroupDescription());
0030 descriptions.add("ThroughputService", desc);
0031 }
0032
0033 ThroughputService::ThroughputService(const edm::ParameterSet& config, edm::ActivityRegistry& registry)
0034 :
0035 m_startup(std::chrono::system_clock::now()),
0036
0037 m_resolution(config.getUntrackedParameter<uint32_t>("eventResolution")),
0038 m_counter(0),
0039 m_events(config.getUntrackedParameter<uint32_t>("eventRange") / m_resolution),
0040 m_print_event_summary(config.getUntrackedParameter<bool>("printEventSummary")),
0041 m_enable_dqm(config.getUntrackedParameter<bool>("enableDQM")),
0042 m_dqm_bynproc(m_enable_dqm ? config.getUntrackedParameter<bool>("dqmPathByProcesses") : false),
0043 m_dqm_path(m_enable_dqm ? config.getUntrackedParameter<std::string>("dqmPath") : ""),
0044 m_time_range(m_enable_dqm ? config.getUntrackedParameter<double>("timeRange") : 0.),
0045 m_time_resolution(m_enable_dqm ? config.getUntrackedParameter<double>("timeResolution") : 0.) {
0046 m_events.clear();
0047 registry.watchPreGlobalBeginRun(this, &ThroughputService::preGlobalBeginRun);
0048 registry.watchPreSourceEvent(this, &ThroughputService::preSourceEvent);
0049 registry.watchPostEvent(this, &ThroughputService::postEvent);
0050 registry.watchPostEndJob(this, &ThroughputService::postEndJob);
0051 }
0052
0053 void ThroughputService::preallocate(edm::service::SystemBounds const& bounds) {
0054 auto concurrent_streams = bounds.maxNumberOfStreams();
0055 auto concurrent_threads = bounds.maxNumberOfThreads();
0056
0057 if (m_enable_dqm and m_dqm_bynproc)
0058 m_dqm_path += fmt::sprintf(
0059 "/Running on %s with %d streams on %d threads", processor_model, concurrent_streams, concurrent_threads);
0060 }
0061
0062 void ThroughputService::preGlobalBeginRun(edm::GlobalContext const& gc) {
0063
0064
0065 if (m_enable_dqm and not edm::Service<dqm::legacy::DQMStore>().isAvailable()) {
0066
0067 m_enable_dqm = false;
0068 edm::LogWarning("ThroughputService") << "The DQMStore is not avalable, the DQM plots will not be generated";
0069 }
0070
0071 if (m_enable_dqm) {
0072 std::string y_axis_title = fmt::sprintf("events / %g s", m_time_resolution);
0073 unsigned int bins = std::round(m_time_range / m_time_resolution);
0074 double range = bins * m_time_resolution;
0075
0076
0077 auto bookTransactionCallback = [&, this](DQMStore::IBooker& booker, DQMStore::IGetter&) {
0078 booker.setCurrentFolder(m_dqm_path);
0079 m_sourced_events = booker.book1D("throughput_sourced", "Throughput (sourced events)", bins, 0., range);
0080 m_sourced_events->setXTitle("time [s]");
0081 m_sourced_events->setYTitle(y_axis_title);
0082 m_retired_events = booker.book1D("throughput_retired", "Throughput (retired events)", bins, 0., range);
0083 m_retired_events->setXTitle("time [s]");
0084 m_retired_events->setYTitle(y_axis_title);
0085 };
0086
0087
0088 edm::Service<DQMStore>()->meBookerGetter(bookTransactionCallback);
0089 } else {
0090 m_sourced_events = nullptr;
0091 m_retired_events = nullptr;
0092 }
0093 }
0094
0095 void ThroughputService::preSourceEvent(edm::StreamID sid) {
0096 auto timestamp = std::chrono::system_clock::now();
0097 auto interval = std::chrono::duration_cast<std::chrono::duration<double>>(timestamp - m_startup).count();
0098 if (m_enable_dqm) {
0099 m_sourced_events->Fill(interval);
0100 }
0101 }
0102
0103 void ThroughputService::postEvent(edm::StreamContext const& sc) {
0104 auto timestamp = std::chrono::system_clock::now();
0105 auto interval = std::chrono::duration_cast<std::chrono::duration<double>>(timestamp - m_startup).count();
0106 if (m_enable_dqm) {
0107 m_retired_events->Fill(interval);
0108 }
0109 ++m_counter;
0110 if (m_counter % m_resolution == 0) {
0111 m_events.push_back(timestamp);
0112 }
0113 }
0114
0115 void ThroughputService::postEndJob() {
0116 if (m_counter < 2 * m_resolution) {
0117
0118 edm::LogWarning("ThroughputService") << "Not enough events to measure the throughput with a resolution of "
0119 << m_resolution << " events";
0120 return;
0121 }
0122
0123 edm::LogInfo info("ThroughputService");
0124
0125 if (m_print_event_summary) {
0126 for (uint32_t i = 0; i < m_events.size(); ++i) {
0127 info << std::setw(8) << (i + 1) * m_resolution << ", " << std::setprecision(6) << edm::TimeOfDay(m_events[i])
0128 << "\n";
0129 }
0130 info << '\n';
0131 }
0132
0133
0134 uint32_t blocks = m_counter / m_resolution - 1;
0135 std::vector<double> delta(blocks);
0136 for (uint32_t i = 0; i < blocks; ++i) {
0137 delta[i] = std::chrono::duration_cast<std::chrono::duration<double>>(m_events[i + 1] - m_events[i]).count();
0138 }
0139
0140 double time_avg = TMath::Mean(delta.begin(), delta.begin() + blocks);
0141 double time_dev = TMath::StdDev(delta.begin(), delta.begin() + blocks);
0142
0143 double throughput_avg = double(m_resolution) / time_avg;
0144 double throughput_dev = double(m_resolution) * time_dev / time_avg / time_avg;
0145
0146 info << "Average throughput: " << throughput_avg << " ± " << throughput_dev << " ev/s";
0147 }
0148
0149
0150 #include "FWCore/ServiceRegistry/interface/ServiceMaker.h"
0151 DEFINE_FWK_SERVICE(ThroughputService);