Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:10:06

0001 // -*- C++ -*-
0002 //
0003 // Package:    DQMDaqInfo
0004 // Class:      DQMDaqInfo
0005 //
0006 /**\class DQMDaqInfo DQMDaqInfo.cc CondCore/DQMDaqInfo/src/DQMDaqInfo.cc
0007    
0008  Description: <one line class summary>
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Ilaria SEGONI
0015 //         Created:  Thu Sep 25 11:17:43 CEST 2008
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 #include <iostream>
0022 #include <fstream>
0023 
0024 // user include files
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026 #include "FWCore/Framework/interface/LuminosityBlock.h"
0027 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0028 #include "FWCore/Framework/interface/Event.h"
0029 #include "FWCore/Framework/interface/MakerMacros.h"
0030 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0031 #include "FWCore/Framework/interface/ESHandle.h"
0032 #include "FWCore/Framework/interface/EventSetup.h"
0033 
0034 //Run Info
0035 #include "CondFormats/DataRecord/interface/RunSummaryRcd.h"
0036 #include "CondFormats/RunInfo/interface/RunSummary.h"
0037 #include "CondFormats/RunInfo/interface/RunInfo.h"
0038 
0039 //DQM
0040 #include "DQMServices/Core/interface/DQMStore.h"
0041 #include "FWCore/ServiceRegistry/interface/Service.h"
0042 
0043 //DataFormats
0044 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0045 
0046 class DQMDaqInfo : public edm::one::EDAnalyzer<> {
0047 public:
0048   typedef dqm::legacy::DQMStore DQMStore;
0049   typedef dqm::legacy::MonitorElement MonitorElement;
0050   explicit DQMDaqInfo(const edm::ParameterSet&);
0051   ~DQMDaqInfo() override = default;
0052 
0053 private:
0054   void beginJob() override;
0055   void beginLuminosityBlock(const edm::LuminosityBlock&, const edm::EventSetup&);
0056   void analyze(const edm::Event&, const edm::EventSetup&) override;
0057 
0058   edm::ESGetToken<RunInfo, RunInfoRcd> runInfoToken_;
0059   DQMStore* dbe_;
0060 
0061   enum subDetList { Pixel, SiStrip, EcalBarrel, EcalEndcap, Hcal, DT, CSC, RPC, L1T };
0062 
0063   MonitorElement* DaqFraction[9];
0064 
0065   std::pair<int, int> PixelRange;
0066   std::pair<int, int> TrackerRange;
0067   std::pair<int, int> CSCRange;
0068   std::pair<int, int> RPCRange;
0069   std::pair<int, int> DTRange;
0070   std::pair<int, int> HcalRange;
0071   std::pair<int, int> ECALBarrRange;
0072   std::pair<int, int> ECALEndcapRangeLow;
0073   std::pair<int, int> ECALEndcapRangeHigh;
0074   std::pair<int, int> L1TRange;
0075 
0076   float NumberOfFeds[9];
0077 };
0078 
0079 DQMDaqInfo::DQMDaqInfo(const edm::ParameterSet& iConfig)
0080     : runInfoToken_{esConsumes<edm::Transition::BeginLuminosityBlock>()} {}
0081 
0082 void DQMDaqInfo::beginLuminosityBlock(const edm::LuminosityBlock& lumiBlock, const edm::EventSetup& iSetup) {
0083   edm::eventsetup::EventSetupRecordKey recordKey(edm::eventsetup::EventSetupRecordKey::TypeTag::findType("RunInfoRcd"));
0084 
0085   if (iSetup.tryToGet<RunInfoRcd>()) {
0086     if (auto sumFED = iSetup.getHandle(runInfoToken_)) {
0087       //const RunInfo* summaryFED=sumFED.product();
0088 
0089       std::vector<int> FedsInIds = sumFED->m_fed_in;
0090 
0091       float FedCount[9] = {0., 0., 0., 0., 0., 0., 0., 0., 0.};
0092 
0093       for (int fedID : FedsInIds) {
0094         if (fedID >= PixelRange.first && fedID <= PixelRange.second)
0095           ++FedCount[Pixel];
0096         if (fedID >= TrackerRange.first && fedID <= TrackerRange.second)
0097           ++FedCount[SiStrip];
0098         if (fedID >= CSCRange.first && fedID <= CSCRange.second)
0099           ++FedCount[CSC];
0100         if (fedID >= RPCRange.first && fedID <= RPCRange.second)
0101           ++FedCount[RPC];
0102         if (fedID >= DTRange.first && fedID <= DTRange.second)
0103           ++FedCount[DT];
0104         if (fedID >= HcalRange.first && fedID <= HcalRange.second)
0105           ++FedCount[Hcal];
0106         if (fedID >= ECALBarrRange.first && fedID <= ECALBarrRange.second)
0107           ++FedCount[EcalBarrel];
0108         if ((fedID >= ECALEndcapRangeLow.first && fedID <= ECALEndcapRangeLow.second) ||
0109             (fedID >= ECALEndcapRangeHigh.first && fedID <= ECALEndcapRangeHigh.second))
0110           ++FedCount[EcalEndcap];
0111         if (fedID >= L1TRange.first && fedID <= L1TRange.second)
0112           ++FedCount[L1T];
0113       }
0114 
0115       for (int detIndex = 0; detIndex < 9; ++detIndex) {
0116         DaqFraction[detIndex]->Fill(FedCount[detIndex] / NumberOfFeds[detIndex]);
0117       }
0118     }
0119   } else {
0120     for (auto& detIndex : DaqFraction)
0121       detIndex->Fill(-1);
0122     return;
0123   }
0124 }
0125 
0126 void DQMDaqInfo::beginJob() {
0127   dbe_ = nullptr;
0128   dbe_ = edm::Service<DQMStore>().operator->();
0129 
0130   std::string commonFolder = "/EventInfo/DAQContents";
0131   std::string subsystFolder;
0132   std::string curentFolder;
0133 
0134   subsystFolder = "Pixel";
0135   curentFolder = subsystFolder + commonFolder;
0136   dbe_->setCurrentFolder(curentFolder);
0137   DaqFraction[Pixel] = dbe_->bookFloat("PixelDaqFraction");
0138 
0139   subsystFolder = "SiStrip";
0140   curentFolder = subsystFolder + commonFolder;
0141   dbe_->setCurrentFolder(curentFolder);
0142   DaqFraction[SiStrip] = dbe_->bookFloat("SiStripDaqFraction");
0143 
0144   subsystFolder = "RPC";
0145   curentFolder = subsystFolder + commonFolder;
0146   dbe_->setCurrentFolder(curentFolder);
0147   DaqFraction[RPC] = dbe_->bookFloat("RPCDaqFraction");
0148 
0149   subsystFolder = "CSC";
0150   curentFolder = subsystFolder + commonFolder;
0151   dbe_->setCurrentFolder(curentFolder);
0152   DaqFraction[CSC] = dbe_->bookFloat("CSCDaqFraction");
0153 
0154   subsystFolder = "DT";
0155   curentFolder = subsystFolder + commonFolder;
0156   dbe_->setCurrentFolder(curentFolder);
0157   DaqFraction[DT] = dbe_->bookFloat("DTDaqFraction");
0158 
0159   subsystFolder = "Hcal";
0160   curentFolder = subsystFolder + commonFolder;
0161   dbe_->setCurrentFolder(curentFolder);
0162   DaqFraction[Hcal] = dbe_->bookFloat("HcalDaqFraction");
0163 
0164   subsystFolder = "EcalBarrel";
0165   curentFolder = subsystFolder + commonFolder;
0166   dbe_->setCurrentFolder(curentFolder);
0167   DaqFraction[EcalBarrel] = dbe_->bookFloat("EcalBarrDaqFraction");
0168 
0169   subsystFolder = "EcalEndcap";
0170   curentFolder = subsystFolder + commonFolder;
0171   dbe_->setCurrentFolder(curentFolder);
0172   DaqFraction[EcalEndcap] = dbe_->bookFloat("EcalEndDaqFraction");
0173 
0174   subsystFolder = "L1T";
0175   curentFolder = subsystFolder + commonFolder;
0176   dbe_->setCurrentFolder(curentFolder);
0177   DaqFraction[L1T] = dbe_->bookFloat("L1TDaqFraction");
0178 
0179   PixelRange.first = FEDNumbering::MINSiPixelFEDID;
0180   PixelRange.second = FEDNumbering::MAXSiPixelFEDID;
0181   TrackerRange.first = FEDNumbering::MINSiStripFEDID;
0182   TrackerRange.second = FEDNumbering::MAXSiStripFEDID;
0183   CSCRange.first = FEDNumbering::MINCSCFEDID;
0184   CSCRange.second = FEDNumbering::MAXCSCFEDID;
0185   RPCRange.first = 790;
0186   RPCRange.second = 792;
0187   DTRange.first = 770;
0188   DTRange.second = 774;
0189   HcalRange.first = FEDNumbering::MINHCALFEDID;
0190   HcalRange.second = FEDNumbering::MAXHCALFEDID;
0191   L1TRange.first = FEDNumbering::MINTriggerGTPFEDID;
0192   L1TRange.second = FEDNumbering::MAXTriggerGTPFEDID;
0193   ECALBarrRange.first = 610;
0194   ECALBarrRange.second = 645;
0195   ECALEndcapRangeLow.first = 601;
0196   ECALEndcapRangeLow.second = 609;
0197   ECALEndcapRangeHigh.first = 646;
0198   ECALEndcapRangeHigh.second = 654;
0199 
0200   NumberOfFeds[Pixel] = PixelRange.second - PixelRange.first + 1;
0201   NumberOfFeds[SiStrip] = TrackerRange.second - TrackerRange.first + 1;
0202   NumberOfFeds[CSC] = CSCRange.second - CSCRange.first + 1;
0203   NumberOfFeds[RPC] = RPCRange.second - RPCRange.first + 1;
0204   NumberOfFeds[DT] = DTRange.second - DTRange.first + 1;
0205   NumberOfFeds[Hcal] = HcalRange.second - HcalRange.first + 1;
0206   NumberOfFeds[EcalBarrel] = ECALBarrRange.second - ECALBarrRange.first + 1;
0207   NumberOfFeds[EcalEndcap] = (ECALEndcapRangeLow.second - ECALEndcapRangeLow.first + 1) +
0208                              (ECALEndcapRangeHigh.second - ECALEndcapRangeHigh.first + 1);
0209   NumberOfFeds[L1T] = L1TRange.second - L1TRange.first + 1;
0210 }
0211 
0212 void DQMDaqInfo::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {}
0213 
0214 #include "FWCore/PluginManager/interface/ModuleDef.h"
0215 #include "FWCore/Framework/interface/MakerMacros.h"
0216 DEFINE_FWK_MODULE(DQMDaqInfo);