File indexing completed on 2024-04-06 12:07:37
0001 #include <string>
0002
0003 #include "DQM/TrackingMonitor/interface/GetLumi.h"
0004 #include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h"
0005 #include "DataFormats/OnlineMetaData/interface/OnlineLuminosityRecord.h"
0006 #include "DataFormats/Scalers/interface/LumiScalers.h"
0007 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0008 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "FWCore/Framework/interface/Frameworkfwd.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0016 #include "FWCore/Utilities/interface/EDGetToken.h"
0017 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0018
0019 namespace {
0020 typedef dqm::reco::DQMStore DQMStore;
0021 typedef dqm::reco::MonitorElement MonitorElement;
0022 struct MEbinning {
0023 int nbins;
0024 double xmin;
0025 double xmax;
0026 };
0027
0028 struct Histograms {
0029 dqm::reco::MonitorElement* numberOfPixelClustersVsLS;
0030 dqm::reco::MonitorElement* numberOfPixelClustersVsLumi;
0031 dqm::reco::MonitorElement* lumiVsLS;
0032 dqm::reco::MonitorElement* puVsLS;
0033 dqm::reco::MonitorElement* pixelLumiVsLS;
0034 dqm::reco::MonitorElement* pixelLumiVsLumi;
0035 };
0036 }
0037
0038
0039
0040
0041
0042 class LumiMonitor : public DQMGlobalEDAnalyzer<Histograms> {
0043 public:
0044 LumiMonitor(const edm::ParameterSet&);
0045 ~LumiMonitor() override = default;
0046
0047 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0048
0049 static void fillHistoPSetDescription(edm::ParameterSetDescription& pset,
0050 int const nbins,
0051 double const xmin,
0052 double const xmax);
0053 static void fillHistoLSPSetDescription(edm::ParameterSetDescription& pset, int const nbins);
0054
0055 private:
0056 void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&, Histograms&) const override;
0057 void dqmAnalyze(edm::Event const&, edm::EventSetup const&, Histograms const&) const override;
0058
0059 static MEbinning getHistoPSet(const edm::ParameterSet& pset);
0060 static MEbinning getHistoLSPSet(const edm::ParameterSet& pset);
0061
0062 std::string const folderName_;
0063
0064 edm::EDGetTokenT<LumiScalersCollection> const lumiScalersToken_;
0065 edm::EDGetTokenT<OnlineLuminosityRecord> const onlineMetaDataDigisToken_;
0066 MEbinning const lumi_binning_;
0067 MEbinning const pu_binning_;
0068 MEbinning const ls_binning_;
0069
0070 bool const doPixelLumi_;
0071 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> const trkTopoToken_;
0072 edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster>> const pixelClustersToken_;
0073 bool const useBPixLayer1_;
0074 int const minNumberOfPixelsPerCluster_;
0075 float const minPixelClusterCharge_;
0076 MEbinning const pixelCluster_binning_;
0077 MEbinning const pixellumi_binning_;
0078 float const lumi_factor_per_bx_;
0079 };
0080
0081
0082
0083
0084
0085 LumiMonitor::LumiMonitor(const edm::ParameterSet& config)
0086 : folderName_(config.getParameter<std::string>("folderName")),
0087 lumiScalersToken_(consumes(config.getParameter<edm::InputTag>("scalers"))),
0088 onlineMetaDataDigisToken_(consumes(config.getParameter<edm::InputTag>("onlineMetaDataDigis"))),
0089 lumi_binning_(getHistoPSet(
0090 config.getParameter<edm::ParameterSet>("histoPSet").getParameter<edm::ParameterSet>("lumiPSet"))),
0091 pu_binning_(
0092 getHistoPSet(config.getParameter<edm::ParameterSet>("histoPSet").getParameter<edm::ParameterSet>("puPSet"))),
0093 ls_binning_(getHistoLSPSet(
0094 config.getParameter<edm::ParameterSet>("histoPSet").getParameter<edm::ParameterSet>("lsPSet"))),
0095 doPixelLumi_(config.getParameter<bool>("doPixelLumi")),
0096 trkTopoToken_(doPixelLumi_ ? esConsumes<TrackerTopology, TrackerTopologyRcd>()
0097 : edm::ESGetToken<TrackerTopology, TrackerTopologyRcd>()),
0098 pixelClustersToken_(doPixelLumi_ ? consumes<edmNew::DetSetVector<SiPixelCluster>>(
0099 config.getParameter<edm::InputTag>("pixelClusters"))
0100 : edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster>>()),
0101 useBPixLayer1_(doPixelLumi_ ? config.getParameter<bool>("useBPixLayer1") : false),
0102 minNumberOfPixelsPerCluster_(doPixelLumi_ ? config.getParameter<int>("minNumberOfPixelsPerCluster") : -1),
0103 minPixelClusterCharge_(doPixelLumi_ ? config.getParameter<double>("minPixelClusterCharge") : -1.),
0104 pixelCluster_binning_(doPixelLumi_ ? getHistoPSet(config.getParameter<edm::ParameterSet>("histoPSet")
0105 .getParameter<edm::ParameterSet>("pixelClusterPSet"))
0106 : MEbinning{}),
0107 pixellumi_binning_(doPixelLumi_ ? getHistoPSet(config.getParameter<edm::ParameterSet>("histoPSet")
0108 .getParameter<edm::ParameterSet>("pixellumiPSet"))
0109 : MEbinning{}),
0110 lumi_factor_per_bx_(useBPixLayer1_
0111 ? GetLumi::FREQ_ORBIT * GetLumi::SECONDS_PER_LS / GetLumi::XSEC_PIXEL_CLUSTER
0112 : GetLumi::FREQ_ORBIT * GetLumi::SECONDS_PER_LS / GetLumi::rXSEC_PIXEL_CLUSTER) {}
0113
0114 MEbinning LumiMonitor::getHistoPSet(const edm::ParameterSet& pset) {
0115 return MEbinning{
0116 pset.getParameter<int32_t>("nbins"),
0117 pset.getParameter<double>("xmin"),
0118 pset.getParameter<double>("xmax"),
0119 };
0120 }
0121
0122 MEbinning LumiMonitor::getHistoLSPSet(const edm::ParameterSet& pset) {
0123 return MEbinning{pset.getParameter<int32_t>("nbins"), -0.5, pset.getParameter<int32_t>("nbins") - 0.5};
0124 }
0125
0126 void LumiMonitor::bookHistograms(DQMStore::IBooker& booker,
0127 edm::Run const& run,
0128 edm::EventSetup const& setup,
0129 Histograms& histograms) const {
0130 booker.setCurrentFolder(folderName_);
0131
0132 if (doPixelLumi_) {
0133 auto me = booker.book1D("numberOfPixelClustersVsLS",
0134 "number of pixel clusters vs lumisection",
0135 ls_binning_.nbins,
0136 ls_binning_.xmin,
0137 ls_binning_.xmax);
0138 me->setAxisTitle("lumisection", 1);
0139 me->setAxisTitle("number of pixel clusters", 2);
0140 histograms.numberOfPixelClustersVsLS = me;
0141
0142 me = booker.bookProfile("numberOfPixelClustersVsLumi",
0143 "number of pixel clusters vs online lumi",
0144 lumi_binning_.nbins,
0145 lumi_binning_.xmin,
0146 lumi_binning_.xmax,
0147 pixelCluster_binning_.xmin,
0148 pixelCluster_binning_.xmax);
0149 me->setAxisTitle("online inst lumi E30 [Hz cm^{-2}]", 1);
0150 me->setAxisTitle("number of pixel clusters", 2);
0151 histograms.numberOfPixelClustersVsLumi = me;
0152
0153 me = booker.bookProfile("pixelLumiVsLS",
0154 "pixel-lumi vs lumisection",
0155 ls_binning_.nbins,
0156 ls_binning_.xmin,
0157 ls_binning_.xmax,
0158 pixellumi_binning_.xmin,
0159 pixellumi_binning_.xmax);
0160 me->setAxisTitle("lumisection", 1);
0161 me->setAxisTitle("pixel-based inst lumi E30 [Hz cm^{-2}]", 2);
0162 histograms.pixelLumiVsLS = me;
0163
0164 me = booker.bookProfile("pixelLumiVsLumi",
0165 "pixel-lumi vs online lumi",
0166 lumi_binning_.nbins,
0167 lumi_binning_.xmin,
0168 lumi_binning_.xmax,
0169 pixellumi_binning_.xmin,
0170 lumi_binning_.xmax);
0171 me->setAxisTitle("online inst lumi E30 [Hz cm^{-2}]", 1);
0172 me->setAxisTitle("pixel-based inst lumi E30 [Hz cm^{-2}]", 2);
0173 histograms.pixelLumiVsLumi = me;
0174 }
0175
0176 auto me = booker.bookProfile("lumiVsLS",
0177 "online lumi vs lumisection",
0178 ls_binning_.nbins,
0179 ls_binning_.xmin,
0180 ls_binning_.xmax,
0181 lumi_binning_.xmin,
0182 lumi_binning_.xmax);
0183 me->setAxisTitle("lumisection", 1);
0184 me->setAxisTitle("online inst lumi E30 [Hz cm^{-2}]", 2);
0185 histograms.lumiVsLS = me;
0186
0187 me = booker.bookProfile("puVsLS",
0188 "online pileup vs lumisection",
0189 ls_binning_.nbins,
0190 ls_binning_.xmin,
0191 ls_binning_.xmax,
0192 pu_binning_.xmin,
0193 pu_binning_.xmax);
0194 me->setAxisTitle("lumisection", 1);
0195 me->setAxisTitle("online pileup", 2);
0196 histograms.puVsLS = me;
0197 }
0198
0199 void LumiMonitor::dqmAnalyze(edm::Event const& event,
0200 edm::EventSetup const& setup,
0201 Histograms const& histograms) const {
0202 int const ls = event.id().luminosityBlock();
0203
0204 float online_lumi = -1.f;
0205 float online_pu = -1.f;
0206 auto const lumiScalersHandle = event.getHandle(lumiScalersToken_);
0207 auto const onlineMetaDataDigisHandle = event.getHandle(onlineMetaDataDigisToken_);
0208 if (lumiScalersHandle.isValid() and not lumiScalersHandle->empty()) {
0209 auto const scalit = lumiScalersHandle->begin();
0210 online_lumi = scalit->instantLumi();
0211 online_pu = scalit->pileup();
0212 } else if (onlineMetaDataDigisHandle.isValid()) {
0213 online_lumi = onlineMetaDataDigisHandle->instLumi();
0214 online_pu = onlineMetaDataDigisHandle->avgPileUp();
0215 }
0216 histograms.lumiVsLS->Fill(ls, online_lumi);
0217 histograms.puVsLS->Fill(ls, online_pu);
0218
0219 if (doPixelLumi_) {
0220 size_t pixel_clusters = 0;
0221 float pixel_lumi = -1.f;
0222 edm::Handle<edmNew::DetSetVector<SiPixelCluster>> pixelClusters;
0223 event.getByToken(pixelClustersToken_, pixelClusters);
0224 if (pixelClusters.isValid()) {
0225 auto const& tTopo = setup.getData(trkTopoToken_);
0226
0227
0228
0229 for (auto pixCluDet = pixelClusters->begin(); pixCluDet != pixelClusters->end(); ++pixCluDet) {
0230 DetId detid = pixCluDet->detId();
0231 size_t subdetid = detid.subdetId();
0232 if (subdetid == (int)PixelSubdetector::PixelBarrel) {
0233 if (tTopo.layer(detid) == 1) {
0234 continue;
0235 }
0236 }
0237
0238 for (auto pixClu = pixCluDet->begin(); pixClu != pixCluDet->end(); ++pixClu) {
0239 if ((pixClu->size() >= minNumberOfPixelsPerCluster_) and (pixClu->charge() >= minPixelClusterCharge_)) {
0240 ++pixel_clusters;
0241 }
0242 }
0243 }
0244 pixel_lumi = lumi_factor_per_bx_ * pixel_clusters / GetLumi::CM2_TO_NANOBARN;
0245 } else {
0246 pixel_lumi = -1.;
0247 }
0248
0249 histograms.numberOfPixelClustersVsLS->Fill(ls, pixel_clusters);
0250 histograms.numberOfPixelClustersVsLumi->Fill(online_lumi, pixel_clusters);
0251 histograms.pixelLumiVsLS->Fill(ls, pixel_lumi);
0252 histograms.pixelLumiVsLumi->Fill(online_lumi, pixel_lumi);
0253 }
0254 }
0255
0256 void LumiMonitor::fillHistoPSetDescription(edm::ParameterSetDescription& pset,
0257 int const nbins,
0258 double const xmin,
0259 double const xmax) {
0260 pset.add<int>("nbins", nbins);
0261 pset.add<double>("xmin", xmin);
0262 pset.add<double>("xmax", xmax);
0263 }
0264
0265 void LumiMonitor::fillHistoLSPSetDescription(edm::ParameterSetDescription& pset, int const nbins) {
0266 pset.add<int>("nbins", nbins);
0267 }
0268
0269 void LumiMonitor::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0270 edm::ParameterSetDescription desc;
0271 desc.add<edm::InputTag>("pixelClusters", edm::InputTag("hltSiPixelClusters"));
0272 desc.add<edm::InputTag>("scalers", edm::InputTag("hltScalersRawToDigi"));
0273 desc.add<edm::InputTag>("onlineMetaDataDigis", edm::InputTag("hltOnlineMetaDataDigis"));
0274 desc.add<std::string>("folderName", "HLT/LumiMonitoring");
0275 desc.add<bool>("doPixelLumi", false);
0276 desc.add<bool>("useBPixLayer1", false);
0277 desc.add<int>("minNumberOfPixelsPerCluster", 2);
0278 desc.add<double>("minPixelClusterCharge", 15000.);
0279
0280 edm::ParameterSetDescription histoPSet;
0281
0282 edm::ParameterSetDescription lsPSet;
0283 fillHistoLSPSetDescription(lsPSet, 2500);
0284 histoPSet.add<edm::ParameterSetDescription>("lsPSet", lsPSet);
0285
0286 edm::ParameterSetDescription puPSet;
0287 fillHistoPSetDescription(puPSet, 130, 0, 130);
0288 histoPSet.add<edm::ParameterSetDescription>("puPSet", puPSet);
0289
0290 edm::ParameterSetDescription lumiPSet;
0291 fillHistoPSetDescription(lumiPSet, 5000, 0, 20000);
0292 histoPSet.add<edm::ParameterSetDescription>("lumiPSet", lumiPSet);
0293
0294 edm::ParameterSetDescription pixellumiPSet;
0295 fillHistoPSetDescription(pixellumiPSet, 300, 0, 3);
0296 histoPSet.add<edm::ParameterSetDescription>("pixellumiPSet", pixellumiPSet);
0297
0298 edm::ParameterSetDescription pixelClusterPSet;
0299 fillHistoPSetDescription(pixelClusterPSet, 200, -0.5, 19999.5);
0300 histoPSet.add("pixelClusterPSet", pixelClusterPSet);
0301
0302 desc.add<edm::ParameterSetDescription>("histoPSet", histoPSet);
0303
0304 descriptions.add("lumiMonitor", desc);
0305 }
0306
0307
0308 #include "FWCore/Framework/interface/MakerMacros.h"
0309 DEFINE_FWK_MODULE(LumiMonitor);