Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 10:48:43

0001 // -*- C++ -*-
0002 //
0003 // Package:    CondFormats/SiPixelObjects
0004 // Class:      SiPixelFEDChannelContainerFromQualityConverter
0005 //
0006 /**\class SiPixelFEDChannelContainerFromQualityConverter SiPixelFEDChannelContainerFromQualityConverter.cc CondFormats/SiPixelObjects/plugins/SiPixelFEDChannelContainerFromQualityConverter.cc
0007  Description: class to build the SiPixelFEDChannelContainer payloads
0008 */
0009 //
0010 // Original Author:  Marco Musich
0011 //         Created:  Wed, 27 Nov 2018 12:04:36 GMT
0012 //
0013 //
0014 
0015 // system include files
0016 #include <memory>
0017 
0018 // user include files
0019 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
0020 #include "CondFormats/DataRecord/interface/SiPixelFedCablingMapRcd.h"
0021 #include "CondFormats/DataRecord/interface/SiPixelQualityFromDbRcd.h"
0022 #include "CondFormats/SiPixelObjects/interface/CablingPathToDetUnit.h"
0023 #include "CondFormats/SiPixelObjects/interface/PixelROC.h"
0024 #include "CondFormats/SiPixelObjects/interface/SiPixelFEDChannelContainer.h"
0025 #include "CondFormats/SiPixelObjects/interface/SiPixelFedCablingMap.h"
0026 #include "CondFormats/SiPixelObjects/interface/SiPixelFedCablingTree.h"
0027 #include "CondFormats/SiPixelObjects/interface/SiPixelQuality.h"
0028 #include "FWCore/Framework/interface/ESHandle.h"
0029 #include "FWCore/Framework/interface/ESWatcher.h"
0030 #include "FWCore/Framework/interface/Event.h"
0031 #include "FWCore/Framework/interface/Frameworkfwd.h"
0032 #include "FWCore/Framework/interface/MakerMacros.h"
0033 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0034 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0035 #include "FWCore/ServiceRegistry/interface/Service.h"
0036 
0037 //
0038 // class declaration
0039 //
0040 
0041 class SiPixelFEDChannelContainerFromQualityConverter : public edm::one::EDAnalyzer<> {
0042 public:
0043   explicit SiPixelFEDChannelContainerFromQualityConverter(const edm::ParameterSet&);
0044   ~SiPixelFEDChannelContainerFromQualityConverter() override;
0045 
0046   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0047   SiPixelFEDChannelContainer::SiPixelFEDChannelCollection createFromSiPixelQuality(
0048       const SiPixelQuality& theQuality, const SiPixelFedCablingMap& theFedCabling);
0049 
0050 private:
0051   void beginJob() override;
0052   void analyze(const edm::Event&, const edm::EventSetup&) override;
0053   void endJob() override;
0054 
0055   // ----------member data ---------------------------
0056   const edm::ESGetToken<SiPixelQuality, SiPixelQualityFromDbRcd> siPixelQualityToken_;
0057   const edm::ESGetToken<SiPixelFedCablingMap, SiPixelFedCablingMapRcd> siPixelCablingToken_;
0058 
0059   const std::string m_record;
0060   const bool printdebug_;
0061   const bool isMC_;
0062   const bool removeEmptyPayloads_;
0063   std::unique_ptr<SiPixelFEDChannelContainer> myQualities;
0064 
0065   int IOVcount_;
0066   edm::ESWatcher<SiPixelQualityFromDbRcd> SiPixelQualityWatcher_;
0067 };
0068 
0069 //
0070 // constructors and destructor
0071 //
0072 SiPixelFEDChannelContainerFromQualityConverter::SiPixelFEDChannelContainerFromQualityConverter(
0073     const edm::ParameterSet& iConfig)
0074     : siPixelQualityToken_(esConsumes()),
0075       siPixelCablingToken_(esConsumes()),
0076       m_record(iConfig.getParameter<std::string>("record")),
0077       printdebug_(iConfig.getUntrackedParameter<bool>("printDebug", false)),
0078       isMC_(iConfig.getUntrackedParameter<bool>("isMC", true)),
0079       removeEmptyPayloads_(iConfig.getUntrackedParameter<bool>("removeEmptyPayloads", false)) {
0080   //now do what ever initialization is needed
0081   myQualities = std::make_unique<SiPixelFEDChannelContainer>();
0082 }
0083 
0084 SiPixelFEDChannelContainerFromQualityConverter::~SiPixelFEDChannelContainerFromQualityConverter() = default;
0085 
0086 //
0087 // member functions
0088 //
0089 
0090 // ------------ method called for each event  ------------
0091 void SiPixelFEDChannelContainerFromQualityConverter::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0092   using namespace edm;
0093 
0094   unsigned int RunNumber_ = iEvent.eventAuxiliary().run();
0095   unsigned int LuminosityBlockNumber_ = iEvent.eventAuxiliary().luminosityBlock();
0096 
0097   bool hasQualityIOV = SiPixelQualityWatcher_.check(iSetup);
0098 
0099   if (hasQualityIOV) {
0100     //Retrieve the strip quality from conditions
0101     edm::ESHandle<SiPixelQuality> siPixelQuality_ = iSetup.getHandle(siPixelQualityToken_);
0102     edm::ESHandle<SiPixelFedCablingMap> cablingMapHandle = iSetup.getHandle(siPixelCablingToken_);
0103 
0104     std::string scenario = std::to_string(RunNumber_) + "_" + std::to_string(LuminosityBlockNumber_);
0105 
0106     edm::LogInfo("SiPixelFEDChannelContainerFromQualityConverter")
0107         << "Found IOV:" << RunNumber_ << "(" << LuminosityBlockNumber_ << ")" << std::endl;
0108 
0109     auto theSiPixelFEDChannelCollection =
0110         this->createFromSiPixelQuality(*(siPixelQuality_.product()), *(cablingMapHandle.product()));
0111 
0112     if (removeEmptyPayloads_ && theSiPixelFEDChannelCollection.empty())
0113       return;
0114 
0115     myQualities->setScenario(scenario, theSiPixelFEDChannelCollection);
0116 
0117     IOVcount_++;
0118   }
0119 }
0120 
0121 // ------------ method called once each job just before starting event loop  ------------
0122 SiPixelFEDChannelContainer::SiPixelFEDChannelCollection
0123 SiPixelFEDChannelContainerFromQualityConverter::createFromSiPixelQuality(const SiPixelQuality& theQuality,
0124                                                                          const SiPixelFedCablingMap& theFedCabling) {
0125   auto fedid_ = theFedCabling.det2fedMap();
0126 
0127   SiPixelFEDChannelContainer::SiPixelFEDChannelCollection theBadChannelCollection;
0128 
0129   auto theDisabledModules = theQuality.getBadComponentList();
0130   for (const auto& mod : theDisabledModules) {
0131     //mod.DetID, mod.errorType,mod.BadRocs
0132 
0133     int coded_badRocs = mod.BadRocs;
0134     std::vector<PixelFEDChannel> disabledChannelsDetSet;
0135     std::vector<sipixelobjects::CablingPathToDetUnit> path = theFedCabling.pathToDetUnit(mod.DetID);
0136     auto cabling_ = theFedCabling.cablingTree();
0137     unsigned int nrocs_inLink(0);
0138     if (!path.empty()) {
0139       const sipixelobjects::PixelFEDCabling* aFed = cabling_->fed(path.at(0).fed);
0140       const sipixelobjects::PixelFEDLink* link = aFed->link(path.at(0).link);
0141       nrocs_inLink = link->numberOfROCs();
0142     } else {
0143       throw cms::Exception("Inconsistent data") << "could not find CablingPathToDetUnit for detId:" << mod.DetID;
0144     }
0145 
0146     std::bitset<16> bad_rocs(coded_badRocs);
0147     unsigned int n_ch = bad_rocs.size() / nrocs_inLink;
0148 
0149     for (unsigned int i_roc = 0; i_roc < n_ch; ++i_roc) {
0150       unsigned int first_idx = nrocs_inLink * i_roc;
0151       unsigned int sec_idx = nrocs_inLink * (i_roc + 1) - 1;
0152       unsigned int mask = pow(2, nrocs_inLink) - 1;
0153       unsigned int n_setbits = (coded_badRocs >> (i_roc * nrocs_inLink)) & mask;
0154 
0155       if (n_setbits == 0) {
0156         continue;
0157       }
0158 
0159       if (n_setbits != mask) {
0160         edm::LogWarning("SiPixelFEDChannelContainerFromQualityConverter")
0161             << "Mismatch! DetId: " << mod.DetID << " " << n_setbits << " " << mask << std::endl;
0162         continue;
0163       }
0164 
0165       if (printdebug_) {
0166         edm::LogVerbatim("SiPixelFEDChannelContainerFromQualityConverter") << "passed" << std::endl;
0167       }
0168 
0169       unsigned int link_id = 99999;
0170       unsigned int fed_id = 99999;
0171 
0172       for (auto const& p : path) {
0173         const sipixelobjects::PixelFEDCabling* aFed = cabling_->fed(p.fed);
0174         const sipixelobjects::PixelFEDLink* link = aFed->link(p.link);
0175         const sipixelobjects::PixelROC* roc = link->roc(p.roc);
0176         unsigned int first_roc = roc->idInDetUnit();
0177 
0178         if (first_roc == first_idx) {
0179           link_id = p.link;
0180           fed_id = p.fed;
0181           break;
0182         }
0183       }
0184 
0185       if (printdebug_) {
0186         edm::LogVerbatim("SiPixelFEDChannelContainerFromQualityConverter")
0187             << " " << fed_id << " " << link_id << " " << first_idx << "  " << sec_idx << std::endl;
0188       }
0189 
0190       PixelFEDChannel ch = {fed_id, link_id, first_idx, sec_idx};
0191       disabledChannelsDetSet.push_back(ch);
0192 
0193       if (printdebug_) {
0194         edm::LogVerbatim("SiPixelFEDChannelContainerFromQualityConverter")
0195             << i_roc << " " << coded_badRocs << " " << first_idx << " " << sec_idx << std::endl;
0196         edm::LogVerbatim("SiPixelFEDChannelContainerFromQualityConverter")
0197             << "=======================================" << std::endl;
0198       }
0199     }
0200 
0201     if (!disabledChannelsDetSet.empty()) {
0202       theBadChannelCollection[mod.DetID] = disabledChannelsDetSet;
0203     }
0204   }
0205   return theBadChannelCollection;
0206 }
0207 
0208 void SiPixelFEDChannelContainerFromQualityConverter::beginJob() { IOVcount_ = 0; }
0209 
0210 // ------------ method called once each job just after ending the event loop  ------------
0211 void SiPixelFEDChannelContainerFromQualityConverter::endJob() {
0212   edm::LogInfo("SiPixelFEDChannelContainerFromQualityConverter") << "Analyzed " << IOVcount_ << " IOVs" << std::endl;
0213   edm::LogInfo("SiPixelFEDChannelContainerFromQualityConverter")
0214       << "Size of SiPixelFEDChannelContainer object " << myQualities->size() << std::endl
0215       << std::endl;
0216 
0217   if (printdebug_) {
0218     edm::LogInfo("SiPixelFEDChannelContainerFromQualityConverter")
0219         << "Content of SiPixelFEDChannelContainer " << std::endl;
0220 
0221     // use built-in method in the CondFormat
0222     myQualities->printAll();
0223   }
0224 
0225   // Form the data here
0226   edm::Service<cond::service::PoolDBOutputService> poolDbService;
0227   if (poolDbService.isAvailable()) {
0228     cond::Time_t valid_time = poolDbService->currentTime();
0229     // this writes the payload to begin in current run defined in cfg
0230     if (!isMC_) {
0231       poolDbService->writeOneIOV(*myQualities, valid_time, m_record);
0232     } else {
0233       // for MC IOV since=1
0234       poolDbService->writeOneIOV(*myQualities, 1, m_record);
0235     }
0236   }
0237 }
0238 
0239 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0240 void SiPixelFEDChannelContainerFromQualityConverter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0241   edm::ParameterSetDescription desc;
0242   desc.setComment("Writes payloads of type SiPixelFEDChannelContainer");
0243   desc.addUntracked<bool>("printDebug", false);
0244   desc.addUntracked<bool>("removeEmptyPayloads", false);
0245   desc.add<std::string>("record", "SiPixelStatusScenariosRcd");
0246   descriptions.add("SiPixelFEDChannelContainerFromQualityConverter", desc);
0247 }
0248 
0249 //define this as a plug-in
0250 DEFINE_FWK_MODULE(SiPixelFEDChannelContainerFromQualityConverter);