Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-11-19 00:04:53

0001 #include <string>
0002 
0003 #include "DQM/TrackingMonitor/interface/GetLumi.h"
0004 #include "DQMServices/Core/interface/DQMGlobalEDAnalyzer.h"
0005 #include "DataFormats/Luminosity/interface/LumiDetails.h"
0006 #include "DataFormats/Luminosity/interface/LumiSummary.h"
0007 #include "DataFormats/Scalers/interface/LumiScalers.h"
0008 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0009 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0017 #include "FWCore/Utilities/interface/EDGetToken.h"
0018 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0019 
0020 namespace {
0021   typedef dqm::reco::DQMStore DQMStore;
0022   typedef dqm::reco::MonitorElement MonitorElement;
0023   struct MEbinning {
0024     int nbins;
0025     double xmin;
0026     double xmax;
0027   };
0028 
0029   struct Histograms {
0030     dqm::reco::MonitorElement* numberOfPixelClustersVsLS;
0031     dqm::reco::MonitorElement* numberOfPixelClustersVsLumi;
0032     dqm::reco::MonitorElement* lumiVsLS;
0033     dqm::reco::MonitorElement* puVsLS;
0034     dqm::reco::MonitorElement* pixelLumiVsLS;
0035     dqm::reco::MonitorElement* pixelLumiVsLumi;
0036   };
0037 }  // namespace
0038 
0039 //
0040 // class declaration
0041 //
0042 
0043 class LumiMonitor : public DQMGlobalEDAnalyzer<Histograms> {
0044 public:
0045   LumiMonitor(const edm::ParameterSet&);
0046   ~LumiMonitor() override = default;
0047   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0048   static void fillHistoPSetDescription(edm::ParameterSetDescription& pset);
0049   static void fillHistoLSPSetDescription(edm::ParameterSetDescription& pset);
0050 
0051 private:
0052   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&, Histograms&) const override;
0053   void dqmAnalyze(edm::Event const&, edm::EventSetup const&, Histograms const&) const override;
0054 
0055   static MEbinning getHistoPSet(const edm::ParameterSet& pset);
0056   static MEbinning getHistoLSPSet(const edm::ParameterSet& pset);
0057 
0058   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> const trkTopoToken_;
0059 
0060   std::string folderName_;
0061 
0062   edm::EDGetTokenT<LumiScalersCollection> lumiScalersToken_;
0063   MEbinning lumi_binning_;
0064   MEbinning pu_binning_;
0065   MEbinning ls_binning_;
0066 
0067   bool doPixelLumi_;
0068   edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster>> pixelClustersToken_;
0069   bool useBPixLayer1_;
0070   int minNumberOfPixelsPerCluster_;
0071   float minPixelClusterCharge_;
0072   MEbinning pixelCluster_binning_;
0073   MEbinning pixellumi_binning_;
0074 
0075   edm::EDGetTokenT<LumiSummary> lumiSummaryToken_;
0076 
0077   float lumi_factor_per_bx_;
0078 };
0079 
0080 // -----------------------------
0081 //  constructors and destructor
0082 // -----------------------------
0083 
0084 LumiMonitor::LumiMonitor(const edm::ParameterSet& config)
0085     : trkTopoToken_{esConsumes()},
0086       folderName_(config.getParameter<std::string>("FolderName")),
0087       lumiScalersToken_(consumes<LumiScalersCollection>(config.getParameter<edm::InputTag>("scalers"))),
0088       lumi_binning_(getHistoPSet(
0089           config.getParameter<edm::ParameterSet>("histoPSet").getParameter<edm::ParameterSet>("lumiPSet"))),
0090       pu_binning_(
0091           getHistoPSet(config.getParameter<edm::ParameterSet>("histoPSet").getParameter<edm::ParameterSet>("puPSet"))),
0092       ls_binning_(getHistoLSPSet(
0093           config.getParameter<edm::ParameterSet>("histoPSet").getParameter<edm::ParameterSet>("lsPSet"))),
0094       doPixelLumi_(config.getParameter<bool>("doPixelLumi")),
0095       pixelClustersToken_(doPixelLumi_ ? consumes<edmNew::DetSetVector<SiPixelCluster>>(
0096                                              config.getParameter<edm::InputTag>("pixelClusters"))
0097                                        : edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster>>()),
0098       useBPixLayer1_(doPixelLumi_ ? config.getParameter<bool>("useBPixLayer1") : false),
0099       minNumberOfPixelsPerCluster_(doPixelLumi_ ? config.getParameter<int>("minNumberOfPixelsPerCluster") : -1),
0100       minPixelClusterCharge_(doPixelLumi_ ? config.getParameter<double>("minPixelClusterCharge") : -1.),
0101       pixelCluster_binning_(doPixelLumi_ ? getHistoPSet(config.getParameter<edm::ParameterSet>("histoPSet")
0102                                                             .getParameter<edm::ParameterSet>("pixelClusterPSet"))
0103                                          : MEbinning{}),
0104       pixellumi_binning_(doPixelLumi_ ? getHistoPSet(config.getParameter<edm::ParameterSet>("histoPSet")
0105                                                          .getParameter<edm::ParameterSet>("pixellumiPSet"))
0106                                       : MEbinning{}) {
0107   if (useBPixLayer1_) {
0108     lumi_factor_per_bx_ = GetLumi::FREQ_ORBIT * GetLumi::SECONDS_PER_LS / GetLumi::XSEC_PIXEL_CLUSTER;
0109   } else {
0110     lumi_factor_per_bx_ = GetLumi::FREQ_ORBIT * GetLumi::SECONDS_PER_LS / GetLumi::rXSEC_PIXEL_CLUSTER;
0111   }
0112 }
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 LS",
0135                             ls_binning_.nbins,
0136                             ls_binning_.xmin,
0137                             ls_binning_.xmax);
0138     me->setAxisTitle("LS", 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 scal lumi",
0144                             lumi_binning_.nbins,
0145                             lumi_binning_.xmin,
0146                             lumi_binning_.xmax,
0147                             pixelCluster_binning_.xmin,
0148                             pixelCluster_binning_.xmax);
0149     me->setAxisTitle("scal 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 LS",
0155                             ls_binning_.nbins,
0156                             ls_binning_.xmin,
0157                             ls_binning_.xmax,
0158                             pixellumi_binning_.xmin,
0159                             pixellumi_binning_.xmax);
0160     me->setAxisTitle("LS", 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 scal lumi",
0166                             lumi_binning_.nbins,
0167                             lumi_binning_.xmin,
0168                             lumi_binning_.xmax,
0169                             pixellumi_binning_.xmin,
0170                             lumi_binning_.xmax);
0171     me->setAxisTitle("scal 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                                "scal lumi vs LS",
0178                                ls_binning_.nbins,
0179                                ls_binning_.xmin,
0180                                ls_binning_.xmax,
0181                                lumi_binning_.xmin,
0182                                lumi_binning_.xmax);
0183   me->setAxisTitle("LS", 1);
0184   me->setAxisTitle("scal inst lumi E30 [Hz cm^{-2}]", 2);
0185   histograms.lumiVsLS = me;
0186 
0187   me = booker.bookProfile("puVsLS",
0188                           "scal PU vs LS",
0189                           ls_binning_.nbins,
0190                           ls_binning_.xmin,
0191                           ls_binning_.xmax,
0192                           pu_binning_.xmin,
0193                           pu_binning_.xmax);
0194   me->setAxisTitle("LS", 1);
0195   me->setAxisTitle("scal PU", 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 ls = event.id().luminosityBlock();
0203 
0204   float scal_lumi = -1.;
0205   float scal_pu = -1.;
0206   edm::Handle<LumiScalersCollection> lumiScalers;
0207   event.getByToken(lumiScalersToken_, lumiScalers);
0208   if (lumiScalers.isValid() and not lumiScalers->empty()) {
0209     auto scalit = lumiScalers->begin();
0210     scal_lumi = scalit->instantLumi();
0211     scal_pu = scalit->pileup();
0212   } else {
0213     scal_lumi = -1.;
0214     scal_pu = -1.;
0215   }
0216   histograms.lumiVsLS->Fill(ls, scal_lumi);
0217   histograms.puVsLS->Fill(ls, scal_pu);
0218 
0219   if (doPixelLumi_) {
0220     size_t pixel_clusters = 0;
0221     float pixel_lumi = -1.;
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       // Count the number of clusters with at least a minimum
0228       // number of pixels per cluster and at least a minimum charge.
0229       size_t tot = 0;
0230       for (auto pixCluDet = pixelClusters->begin(); pixCluDet != pixelClusters->end(); ++pixCluDet) {
0231         DetId detid = pixCluDet->detId();
0232         size_t subdetid = detid.subdetId();
0233         if (subdetid == (int)PixelSubdetector::PixelBarrel) {
0234           if (tTopo.layer(detid) == 1) {
0235             continue;
0236           }
0237         }
0238 
0239         for (auto pixClu = pixCluDet->begin(); pixClu != pixCluDet->end(); ++pixClu) {
0240           ++tot;
0241           if ((pixClu->size() >= minNumberOfPixelsPerCluster_) and (pixClu->charge() >= minPixelClusterCharge_)) {
0242             ++pixel_clusters;
0243           }
0244         }
0245       }
0246       pixel_lumi = lumi_factor_per_bx_ * pixel_clusters / GetLumi::CM2_TO_NANOBARN;  // ?!?!
0247     } else {
0248       pixel_lumi = -1.;
0249     }
0250 
0251     histograms.numberOfPixelClustersVsLS->Fill(ls, pixel_clusters);
0252     histograms.numberOfPixelClustersVsLumi->Fill(scal_lumi, pixel_clusters);
0253     histograms.pixelLumiVsLS->Fill(ls, pixel_lumi);
0254     histograms.pixelLumiVsLumi->Fill(scal_lumi, pixel_lumi);
0255   }
0256 }
0257 
0258 void LumiMonitor::fillHistoPSetDescription(edm::ParameterSetDescription& pset) {
0259   pset.add<int>("nbins");
0260   pset.add<double>("xmin");
0261   pset.add<double>("xmax");
0262 }
0263 
0264 void LumiMonitor::fillHistoLSPSetDescription(edm::ParameterSetDescription& pset) { pset.add<int>("nbins", 2500); }
0265 
0266 void LumiMonitor::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0267   edm::ParameterSetDescription desc;
0268   desc.add<edm::InputTag>("pixelClusters", edm::InputTag("hltSiPixelClusters"));
0269   desc.add<edm::InputTag>("scalers", edm::InputTag("hltScalersRawToDigi"));
0270   desc.add<std::string>("FolderName", "HLT/LumiMonitoring");
0271   desc.add<bool>("doPixelLumi", false);
0272   desc.add<bool>("useBPixLayer1", false);
0273   desc.add<int>("minNumberOfPixelsPerCluster", 2);  // from DQM/PixelLumi/python/PixelLumiDQM_cfi.py
0274   desc.add<double>("minPixelClusterCharge", 15000.);
0275 
0276   edm::ParameterSetDescription histoPSet;
0277   edm::ParameterSetDescription pixelClusterPSet;
0278   LumiMonitor::fillHistoPSetDescription(pixelClusterPSet);
0279   histoPSet.add("pixelClusterPSet", pixelClusterPSet);
0280 
0281   edm::ParameterSetDescription lumiPSet;
0282   fillHistoPSetDescription(lumiPSet);
0283   histoPSet.add<edm::ParameterSetDescription>("lumiPSet", lumiPSet);
0284 
0285   edm::ParameterSetDescription puPSet;
0286   fillHistoPSetDescription(puPSet);
0287   histoPSet.add<edm::ParameterSetDescription>("puPSet", puPSet);
0288 
0289   edm::ParameterSetDescription pixellumiPSet;
0290   fillHistoPSetDescription(pixellumiPSet);
0291   histoPSet.add<edm::ParameterSetDescription>("pixellumiPSet", pixellumiPSet);
0292 
0293   edm::ParameterSetDescription lsPSet;
0294   fillHistoLSPSetDescription(lsPSet);
0295   histoPSet.add<edm::ParameterSetDescription>("lsPSet", lsPSet);
0296 
0297   desc.add<edm::ParameterSetDescription>("histoPSet", histoPSet);
0298 
0299   descriptions.add("lumiMonitor", desc);
0300 }
0301 
0302 // Define this as a plug-in
0303 #include "FWCore/Framework/interface/MakerMacros.h"
0304 DEFINE_FWK_MODULE(LumiMonitor);