Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-05-18 03:27:07

0001 // Package:    RPCDqmClient
0002 // Original Author:  Anna Cimmino
0003 #include "DQM/RPCMonitorClient/interface/RPCDqmClient.h"
0004 #include "DQM/RPCMonitorClient/interface/RPCNameHelper.h"
0005 #include "DQM/RPCMonitorClient/interface/RPCBookFolderStructure.h"
0006 //include client headers
0007 #include "DQM/RPCMonitorClient/interface/RPCDeadChannelTest.h"
0008 #include "DQM/RPCMonitorClient/interface/RPCMultiplicityTest.h"
0009 #include "DQM/RPCMonitorClient/interface/RPCClusterSizeTest.h"
0010 #include "DQM/RPCMonitorClient/interface/RPCOccupancyTest.h"
0011 #include "DQM/RPCMonitorClient/interface/RPCNoisyStripTest.h"
0012 //Geometry
0013 #include "Geometry/RPCGeometry/interface/RPCGeomServ.h"
0014 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0015 //Framework
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/Framework/interface/ESHandle.h"
0018 
0019 #include <fmt/format.h>
0020 
0021 RPCDqmClient::RPCDqmClient(const edm::ParameterSet& pset) {
0022   edm::LogVerbatim("rpcdqmclient") << "[RPCDqmClient]: Constructor";
0023 
0024   offlineDQM_ = pset.getUntrackedParameter<bool>("OfflineDQM", true);
0025   useRollInfo_ = pset.getUntrackedParameter<bool>("UseRollInfo", false);
0026   //check enabling
0027   enableDQMClients_ = pset.getUntrackedParameter<bool>("EnableRPCDqmClient", true);
0028   minimumEvents_ = pset.getUntrackedParameter<int>("MinimumRPCEvents", 10000);
0029   numberOfDisks_ = pset.getUntrackedParameter<int>("NumberOfEndcapDisks", 4);
0030   numberOfRings_ = pset.getUntrackedParameter<int>("NumberOfEndcapRings", 2);
0031 
0032   std::string subsystemFolder = pset.getUntrackedParameter<std::string>("RPCFolder", "RPC");
0033   std::string recHitTypeFolder = pset.getUntrackedParameter<std::string>("RecHitTypeFolder", "AllHits");
0034   std::string summaryFolder = pset.getUntrackedParameter<std::string>("SummaryFolder", "SummaryHistograms");
0035 
0036   prefixDir_ = subsystemFolder + "/" + recHitTypeFolder;
0037   globalFolder_ = subsystemFolder + "/" + recHitTypeFolder + "/" + summaryFolder;
0038 
0039   //get prescale factor
0040   prescaleGlobalFactor_ = pset.getUntrackedParameter<int>("DiagnosticGlobalPrescale", 5);
0041 
0042   //make default client list
0043   clientList_ = {{"RPCMultiplicityTest", "RPCDeadChannelTest", "RPCClusterSizeTest"}};
0044   clientList_ = pset.getUntrackedParameter<std::vector<std::string> >("RPCDqmClientList", clientList_);
0045 
0046   //get all the possible RPC DQM clients
0047   this->makeClientMap(pset);
0048 
0049   //clear counters
0050   lumiCounter_ = 0;
0051 
0052   rpcGeomToken_ = esConsumes<edm::Transition::EndLuminosityBlock>();
0053 }
0054 
0055 void RPCDqmClient::beginJob() {
0056   if (!enableDQMClients_) {
0057     return;
0058   };
0059   edm::LogVerbatim("rpcdqmclient") << "[RPCDqmClient]: Begin Job";
0060 
0061   //Do whatever the begin jobs of all client modules do
0062   for (auto& module : clientModules_) {
0063     module->beginJob(globalFolder_);
0064   }
0065 }
0066 
0067 void RPCDqmClient::dqmEndLuminosityBlock(DQMStore::IBooker& ibooker,
0068                                          DQMStore::IGetter& igetter,
0069                                          edm::LuminosityBlock const& lumiSeg,
0070                                          edm::EventSetup const& c) {
0071   if (!enableDQMClients_) {
0072     return;
0073   }
0074   edm::LogVerbatim("rpcdqmclient") << "[RPCDqmClient]: End DQM LB";
0075 
0076   if (myDetIds_.empty()) {
0077     //Get RPCdetId...
0078 
0079     this->getRPCdetId(c);
0080 
0081     //...book summary histograms
0082     for (auto& module : clientModules_) {
0083       module->myBooker(ibooker);
0084     }
0085   }
0086 
0087   if (!offlineDQM_) {  //Do this only for the online
0088 
0089     if (lumiCounter_ == 0) {  //only for the first lumi section do this...
0090       // ...get chamber based histograms and pass them to the client modules
0091       this->getMonitorElements(igetter);
0092     }
0093 
0094     //Do not perform client oparations every lumi block
0095     ++lumiCounter_;
0096     if (lumiCounter_ % prescaleGlobalFactor_ != 0) {
0097       return;
0098     }
0099 
0100     //Check if there's enough statistics
0101     float rpcevents = minimumEvents_;
0102     if (RPCEvents_) {
0103       rpcevents = RPCEvents_->getBinContent(1);
0104     }
0105     if (rpcevents < minimumEvents_) {
0106       return;
0107     }
0108 
0109     edm::LogVerbatim("rpcdqmclient") << "[RPCDqmClient]: Client operations";
0110     for (auto& module : clientModules_) {
0111       module->clientOperation();
0112     }
0113   }  //end of online operations
0114 }
0115 
0116 void RPCDqmClient::dqmEndJob(DQMStore::IBooker& ibooker, DQMStore::IGetter& igetter) {
0117   if (!enableDQMClients_) {
0118     return;
0119   }
0120 
0121   edm::LogVerbatim("rpcdqmclient") << "[RPCDqmClient]: End DQM Job";
0122 
0123   if (offlineDQM_) {  // ...get chamber based histograms and pass them to the client modules
0124     this->getMonitorElements(igetter);
0125   }
0126 
0127   float rpcevents = minimumEvents_;
0128   if (RPCEvents_) {
0129     rpcevents = RPCEvents_->getBinContent(1);
0130   }
0131   if (rpcevents < minimumEvents_) {
0132     return;
0133   }
0134 
0135   edm::LogVerbatim("rpcdqmclient") << "[RPCDqmClient]: Client operations";
0136   for (auto& module : clientModules_) {
0137     module->clientOperation();
0138   }
0139 }
0140 
0141 void RPCDqmClient::getMonitorElements(DQMStore::IGetter& igetter) {
0142   std::vector<MonitorElement*> myMeVect;
0143   std::vector<RPCDetId> myDetIds;
0144 
0145   //loop on all geometry and get all histos
0146   for (auto& detId : myDetIds_) {
0147     //Get name
0148     const std::string rollName = RPCNameHelper::name(detId, useRollInfo_);
0149     const std::string folder = RPCBookFolderStructure::folderStructure(detId);
0150 
0151     //loop on clients
0152     for (unsigned int cl = 0, nCL = clientModules_.size(); cl < nCL; ++cl) {
0153       if (clientHisto_[cl] == "ClusterSize")
0154         continue;
0155 
0156       MonitorElement* myMe = igetter.get(fmt::format("{}/{}/{}_{}", prefixDir_, folder, clientHisto_[cl], rollName));
0157       if (!myMe)
0158         continue;
0159 
0160       myMeVect.push_back(myMe);
0161       myDetIds.push_back(detId);
0162 
0163     }  //end loop on clients
0164   }    //end loop on all geometry and get all histos
0165 
0166   //Clustersize
0167   std::vector<MonitorElement*> myMeVectCl;
0168   const std::array<std::string, 4> chNames = {{"CH01-CH09", "CH10-CH18", "CH19-CH27", "CH28-CH36"}};
0169 
0170   //Retrieve barrel clustersize
0171   for (int wheel = -2; wheel <= 2; wheel++) {
0172     for (int sector = 1; sector <= 12; sector++) {
0173       MonitorElement* myMeCl = igetter.get(fmt::format(
0174           "{}/Barrel/Wheel_{}/SummaryBySectors/ClusterSize_Wheel_{}_Sector_{}", prefixDir_, wheel, wheel, sector));
0175       myMeVectCl.push_back(myMeCl);
0176     }
0177   }
0178   //Retrieve endcaps clustersize
0179   for (int region = -1; region <= 1; region++) {
0180     if (region == 0)
0181       continue;
0182 
0183     std::string regionName = "Endcap-";
0184     if (region == 1)
0185       regionName = "Endcap+";
0186 
0187     for (int disk = 1; disk <= numberOfDisks_; disk++) {
0188       for (int ring = numberOfRings_; ring <= 3; ring++) {
0189         for (unsigned int ich = 0; ich < chNames.size(); ich++) {
0190           MonitorElement* myMeCl =
0191               igetter.get(fmt::format("{}/{}/Disk_{}/SummaryByRings/ClusterSize_Disk_{}_Ring_{}_{}",
0192                                       prefixDir_,
0193                                       regionName,
0194                                       (region * disk),
0195                                       (region * disk),
0196                                       ring,
0197                                       chNames[ich]));
0198           myMeVectCl.push_back(myMeCl);
0199         }
0200       }
0201     }
0202   }
0203 
0204   RPCEvents_ = igetter.get(prefixDir_ + "/RPCEvents");
0205   for (unsigned int cl = 0; cl < clientModules_.size(); ++cl) {
0206     if (clientHisto_[cl] == "ClusterSize")
0207       clientModules_[cl]->getMonitorElements(myMeVectCl, myDetIds, clientHisto_[cl]);
0208     else
0209       clientModules_[cl]->getMonitorElements(myMeVect, myDetIds, clientHisto_[cl]);
0210   }
0211 }
0212 
0213 void RPCDqmClient::getRPCdetId(const edm::EventSetup& eventSetup) {
0214   myDetIds_.clear();
0215 
0216   auto rpcGeo = eventSetup.getHandle(rpcGeomToken_);
0217 
0218   for (auto& det : rpcGeo->dets()) {
0219     const RPCChamber* ch = dynamic_cast<const RPCChamber*>(det);
0220     if (!ch)
0221       continue;
0222 
0223     //Loop on rolls in given chamber
0224     for (auto& r : ch->rolls()) {
0225       myDetIds_.push_back(r->id());
0226     }
0227   }
0228 }
0229 
0230 void RPCDqmClient::makeClientMap(const edm::ParameterSet& pset) {
0231   for (unsigned int i = 0; i < clientList_.size(); i++) {
0232     if (clientList_[i] == "RPCMultiplicityTest") {
0233       clientHisto_.push_back("Multiplicity");
0234       // clientTag_.push_back(rpcdqm::MULTIPLICITY);
0235       clientModules_.emplace_back(new RPCMultiplicityTest(pset));
0236     } else if (clientList_[i] == "RPCDeadChannelTest") {
0237       clientHisto_.push_back("Occupancy");
0238       clientModules_.emplace_back(new RPCDeadChannelTest(pset));
0239       // clientTag_.push_back(rpcdqm::OCCUPANCY);
0240     } else if (clientList_[i] == "RPCClusterSizeTest") {
0241       clientHisto_.push_back("ClusterSize");
0242       clientModules_.emplace_back(new RPCClusterSizeTest(pset));
0243       // clientTag_.push_back(rpcdqm::CLUSTERSIZE);
0244     } else if (clientList_[i] == "RPCOccupancyTest") {
0245       clientHisto_.push_back("Occupancy");
0246       clientModules_.emplace_back(new RPCOccupancyTest(pset));
0247       // clientTag_.push_back(rpcdqm::OCCUPANCY);
0248     } else if (clientList_[i] == "RPCNoisyStripTest") {
0249       clientHisto_.push_back("Occupancy");
0250       clientModules_.emplace_back(new RPCNoisyStripTest(pset));
0251       //clientTag_.push_back(rpcdqm::OCCUPANCY);
0252     }
0253   }
0254 
0255   return;
0256 }