CalibTracker/plugins/SiStripQualityHotStripIdentifierRoot.cc

0001 #include "DQMOffline/CalibTracker/plugins/SiStripQualityHotStripIdentifierRoot.h"
0002
0003 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0004
0005 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0006 #include "DataFormats/TrackerRecHit2D/interface/ProjectedSiStripRecHit2D.h"
0007 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0008 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0009 #include "DataFormats/TrackReco/interface/Track.h"
0010 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0011
0012 #include <fstream>
0013 #include <iostream>
0014 #include <memory>
0015 #include <sstream>
0016
0017 #include "TH1F.h"
0018
0019 //Insert here the include to the algos
0020 #include "CalibTracker/SiStripQuality/interface/SiStripHotStripAlgorithmFromClusterOccupancy.h"
0021 #include "CalibTracker/SiStripQuality/interface/SiStripBadAPVAlgorithmFromClusterOccupancy.h"
0022 #include "CalibTracker/SiStripQuality/interface/SiStripBadAPVandHotStripAlgorithmFromClusterOccupancy.h"
0023 #include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h"
0024
0025 SiStripQualityHotStripIdentifierRoot::SiStripQualityHotStripIdentifierRoot(const edm::ParameterSet& iConfig)
0026     : ConditionDBWriter<SiStripBadStrip>(iConfig),
0027       UseInputDB_(iConfig.getUntrackedParameter<bool>("UseInputDB", false)),
0028       conf_(iConfig),
0029       tracker_(nullptr),
0030       tTopo(nullptr),
0031       filename(iConfig.getUntrackedParameter<std::string>("rootFilename", "CondDB_TKCC_20X_v3_hlt_50822.root")),
0032       dirpath(iConfig.getUntrackedParameter<std::string>("rootDirPath", "")),
0033       TotNumberOfEvents(0),
0034       MeanNumberOfCluster(0),
0035       calibrationthreshold(iConfig.getUntrackedParameter<uint32_t>("CalibrationThreshold", 10000)),
0036       tTopoToken_(esConsumes<edm::Transition::BeginRun>()),
0037       tkGeomToken_(esConsumes<edm::Transition::BeginRun>()),
0038       stripQualityToken_(UseInputDB_
0039                              ? decltype(stripQualityToken_)(esConsumes<edm::Transition::BeginRun>(
0040                                    edm::ESInputTag{"", iConfig.getUntrackedParameter<std::string>("dataLabel", "")}))
0041                              : decltype(stripQualityToken_){}) {
0042   dqmStore_ = edm::Service<DQMStore>().operator->();
0043
0044   if (!filename.empty()) {
0045     edm::LogInfo("SiStripQualityHotStripIdentifierRoot") << " before opening file " << filename.c_str();
0046     dqmStore_->open(filename, false);
0047     edm::LogInfo("SiStripQualityHotStripIdentifierRoot") << " after opening file ";
0048   }
0049 }
0050
0051 SiStripQualityHotStripIdentifierRoot::~SiStripQualityHotStripIdentifierRoot() {}
0052
0053 std::unique_ptr<SiStripBadStrip> SiStripQualityHotStripIdentifierRoot::getNewObject() {
0054   edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0055       << "SiStripQualityHotStripIdentifierRoot::getNewObject called" << std::endl;
0056
0057   auto obj = std::make_unique<SiStripBadStrip>();
0058
0059   edm::ParameterSet parameters = conf_.getParameter<edm::ParameterSet>("AlgoParameters");
0060   std::string AlgoName = parameters.getParameter<std::string>("AlgoName");
0061
0062   bookHistos();
0063
0064   if (TotNumberOfEvents >= calibrationthreshold) {
0065     edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0066         << " [SiStripQualityHotStripIdentifierRoot::getNewObject] Total number of events is " << TotNumberOfEvents
0067         << "  which is bigger than  " << calibrationthreshold << ". Calibration is preparing to be launched."
0068         << std::endl;
0069     if (MeanNumberOfCluster > 0.01) {
0070       edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0071           << "Mean Number of Cluster is " << MeanNumberOfCluster << ". Calibration is launched." << std::endl;
0072
0073       if (AlgoName == "SiStripHotStripAlgorithmFromClusterOccupancy") {
0074         edm::LogInfo("SiStripQualityHotStripIdentifierRoot") << " [SiStripQualityHotStripIdentifierRoot::getNewObject] "
0075                                                                 "call to SiStripHotStripAlgorithmFromClusterOccupancy"
0076                                                              << std::endl;
0077
0078         theIdentifier = new SiStripHotStripAlgorithmFromClusterOccupancy(conf_, tTopo);
0079         theIdentifier->setProbabilityThreshold(parameters.getUntrackedParameter<double>("ProbabilityThreshold", 1.E-7));
0080         theIdentifier->setMinNumEntries(parameters.getUntrackedParameter<uint32_t>("MinNumEntries", 100));
0081         theIdentifier->setMinNumEntriesPerStrip(parameters.getUntrackedParameter<uint32_t>("MinNumEntriesPerStrip", 5));
0082         theIdentifier->setOccupancyThreshold(parameters.getUntrackedParameter<double>("OccupancyThreshold", 1.E-5));
0083         theIdentifier->setNumberOfEvents(TotNumberOfEvents);
0084         theIdentifier->setOutputFileName(
0085             conf_.getUntrackedParameter<std::string>("OccupancyRootFile", "Occupancy.root"),
0086             conf_.getUntrackedParameter<bool>("WriteOccupancyRootFile", false));
0087         theIdentifier->setTrackerGeometry(tracker_);
0088
0089         const auto detInfo =
0090             SiStripDetInfoFileReader::read(edm::FileInPath{SiStripDetInfoFileReader::kDefaultFile}.fullPath());
0091         SiStripQuality* qobj = new SiStripQuality(detInfo);
0092         theIdentifier->extractBadStrips(
0093             qobj,
0094             ClusterPositionHistoMap,
0095             SiStripQuality_);  //here I insert SiStripQuality as input and get qobj as output
0096
0097         //----------
0098
0099         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0100             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] copy SiStripObject in SiStripBadStrip"
0101             << std::endl;
0102
0103         std::stringstream ss;
0104
0105         SiStripBadStrip::RegistryIterator rIter = qobj->getRegistryVectorBegin();
0106         SiStripBadStrip::RegistryIterator rIterEnd = qobj->getRegistryVectorEnd();
0107         for (; rIter != rIterEnd; ++rIter) {
0108           SiStripBadStrip::Range range(qobj->getDataVectorBegin() + rIter->ibegin,
0109                                        qobj->getDataVectorBegin() + rIter->iend);
0110           if (!obj->put(rIter->detid, range))
0111             edm::LogError("SiStripQualityHotStripIdentifierRoot")
0112                 << "[SiStripQualityHotStripIdentifierRoot::getNewObject] detid already exists" << std::endl;
0113         }
0114         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0115             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] " << ss.str() << std::endl;
0116       } else if (AlgoName == "SiStripBadAPVAlgorithmFromClusterOccupancy") {
0117         edm::LogInfo("SiStripQualityHotStripIdentifierRoot") << " [SiStripQualityHotStripIdentifierRoot::getNewObject] "
0118                                                                 "call to SiStripBadAPVAlgorithmFromClusterOccupancy"
0119                                                              << std::endl;
0120
0121         theIdentifier2 = new SiStripBadAPVAlgorithmFromClusterOccupancy(conf_, tTopo);
0122         theIdentifier2->setLowOccupancyThreshold(parameters.getUntrackedParameter<double>("LowOccupancyThreshold", 5));
0123         theIdentifier2->setHighOccupancyThreshold(
0124             parameters.getUntrackedParameter<double>("HighOccupancyThreshold", 10));
0125         theIdentifier2->setAbsoluteLowThreshold(parameters.getUntrackedParameter<double>("AbsoluteLowThreshold", 0));
0126         theIdentifier2->setNumberIterations(parameters.getUntrackedParameter<uint32_t>("NumberIterations", 2));
0127         theIdentifier2->setAbsoluteOccupancyThreshold(
0128             parameters.getUntrackedParameter<double>("OccupancyThreshold", 1.E-5));
0129         theIdentifier2->setNumberOfEvents(TotNumberOfEvents);
0130         theIdentifier2->setMinNumOfEvents();
0131         theIdentifier2->setOutputFileName(
0132             conf_.getUntrackedParameter<std::string>("OccupancyRootFile", "Occupancy.root"),
0133             conf_.getUntrackedParameter<bool>("WriteOccupancyRootFile", false));
0134         theIdentifier2->setTrackerGeometry(tracker_);
0135
0136         const auto detInfo =
0137             SiStripDetInfoFileReader::read(edm::FileInPath{SiStripDetInfoFileReader::kDefaultFile}.fullPath());
0138         SiStripQuality* qobj = new SiStripQuality(detInfo);
0139         theIdentifier2->extractBadAPVs(qobj, ClusterPositionHistoMap, SiStripQuality_);
0140
0141         //----------
0142
0143         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0144             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] copy SiStripObject in SiStripBadStrip"
0145             << std::endl;
0146
0147         std::stringstream ss;
0148
0149         SiStripBadStrip::RegistryIterator rIter = qobj->getRegistryVectorBegin();
0150         SiStripBadStrip::RegistryIterator rIterEnd = qobj->getRegistryVectorEnd();
0151         for (; rIter != rIterEnd; ++rIter) {
0152           SiStripBadStrip::Range range(qobj->getDataVectorBegin() + rIter->ibegin,
0153                                        qobj->getDataVectorBegin() + rIter->iend);
0154           if (!obj->put(rIter->detid, range))
0155             edm::LogError("SiStripQualityHotStripIdentifierRoot")
0156                 << "[SiStripQualityHotStripIdentifierRoot::getNewObject] detid already exists" << std::endl;
0157         }
0158         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0159             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] " << ss.str() << std::endl;
0160
0161       } else if (AlgoName == "SiStripBadAPVandHotStripAlgorithmFromClusterOccupancy") {
0162         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0163             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] call to "
0164                "SiStripBadAPVandHotStripAlgorithmFromClusterOccupancy"
0165             << std::endl;
0166
0167         theIdentifier3 = new SiStripBadAPVandHotStripAlgorithmFromClusterOccupancy(conf_, tTopo);
0168         theIdentifier3->setProbabilityThreshold(
0169             parameters.getUntrackedParameter<double>("ProbabilityThreshold", 1.E-7));
0170         theIdentifier3->setMinNumEntries(parameters.getUntrackedParameter<uint32_t>("MinNumEntries", 100));
0171         theIdentifier3->setMinNumEntriesPerStrip(
0172             parameters.getUntrackedParameter<uint32_t>("MinNumEntriesPerStrip", 5));
0173         theIdentifier3->setNumberOfEvents(TotNumberOfEvents);
0174         theIdentifier3->setOutputFileName(
0175             conf_.getUntrackedParameter<std::string>("OccupancyRootFile", "Occupancy.root"),
0176             conf_.getUntrackedParameter<bool>("WriteOccupancyRootFile", false),
0177             conf_.getUntrackedParameter<std::string>("DQMHistoOutputFile", "DQMHistos.root"),
0178             conf_.getUntrackedParameter<bool>("WriteDQMHistoOutputFile", false));
0179         theIdentifier3->setTrackerGeometry(tracker_);
0180         theIdentifier3->setLowOccupancyThreshold(parameters.getUntrackedParameter<double>("LowOccupancyThreshold", 5));
0181         theIdentifier3->setHighOccupancyThreshold(
0182             parameters.getUntrackedParameter<double>("HighOccupancyThreshold", 10));
0183         theIdentifier3->setAbsoluteLowThreshold(parameters.getUntrackedParameter<double>("AbsoluteLowThreshold", 0));
0184         theIdentifier3->setNumberIterations(parameters.getUntrackedParameter<uint32_t>("NumberIterations", 2));
0185         theIdentifier3->setAbsoluteOccupancyThreshold(
0186             parameters.getUntrackedParameter<double>("OccupancyThreshold", 1.E-5));
0187         theIdentifier3->setMinNumOfEvents();
0188
0189         const auto detInfo =
0190             SiStripDetInfoFileReader::read(edm::FileInPath{SiStripDetInfoFileReader::kDefaultFile}.fullPath());
0191         SiStripQuality* qobj = new SiStripQuality(detInfo);
0192         theIdentifier3->extractBadAPVSandStrips(
0193             qobj,
0194             ClusterPositionHistoMap,
0195             SiStripQuality_);  //here I insert SiStripQuality as input and get qobj as output
0196
0197         //----------
0198
0199         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0200             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] copy SiStripObject in SiStripBadStrip"
0201             << std::endl;
0202
0203         std::stringstream ss;
0204
0205         SiStripBadStrip::RegistryIterator rIter = qobj->getRegistryVectorBegin();
0206         SiStripBadStrip::RegistryIterator rIterEnd = qobj->getRegistryVectorEnd();
0207         for (; rIter != rIterEnd; ++rIter) {
0208           SiStripBadStrip::Range range(qobj->getDataVectorBegin() + rIter->ibegin,
0209                                        qobj->getDataVectorBegin() + rIter->iend);
0210           if (!obj->put(rIter->detid, range))
0211             edm::LogError("SiStripQualityHotStripIdentifierRoot")
0212                 << "[SiStripQualityHotStripIdentifierRoot::getNewObject] detid already exists" << std::endl;
0213         }
0214         edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0215             << " [SiStripQualityHotStripIdentifierRoot::getNewObject] " << ss.str() << std::endl;
0216       } else {
0217         edm::LogError("SiStripQualityHotStripIdentifierRoot") << " [SiStripQualityHotStripIdentifierRoot::getNewObject]"
0218                                                                  " call for a unknow HotStrip identification algoritm"
0219                                                               << std::endl;
0220
0221         std::vector<uint32_t> a;
0222         SiStripBadStrip::Range range(a.begin(), a.end());
0223         if (!obj->put(0xFFFFFFFF, range))
0224           edm::LogError("SiStripQualityHotStripIdentifierRoot")
0225               << "[SiStripQualityHotStripIdentifierRoot::getNewObject] detid already exists" << std::endl;
0226       }
0227     }
0228
0229     else {
0230       edm::LogWarning("SiStripQualityHotStripIdentifierRoot")
0231           << " [SiStripQualityHotStripIdentifierRoot::getNewObject] :: Empty Events :: "
0232           << "Mean Number of Cluster is  " << MeanNumberOfCluster << " Calibration NOT launched" << std::endl;
0233       setDoStore(false);  // Don't put anything in the sqlite-file!
0234     }
0235   } else {
0236     edm::LogWarning("SiStripQualityHotStripIdentifierRoot")
0237         << " [SiStripQualityHotStripIdentifierRoot::getNewObject] :: Not Enough Events :: "
0238         << "Total number of events is  " << TotNumberOfEvents << " Calibration NOT launched" << std::endl;
0239     setDoStore(false);  // Don't put anything in the sqlite-file!
0240   }
0241
0242   return obj;
0243 }
0244
0245 void SiStripQualityHotStripIdentifierRoot::algoBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0246   tTopo = &iSetup.getData(tTopoToken_);
0247   tracker_ = &iSetup.getData(tkGeomToken_);
0248   if (UseInputDB_ && stripQualityWatcher_.check(iSetup)) {
0249     SiStripQuality_ = &iSetup.getData(stripQualityToken_);
0250   }
0251 }
0252
0253 void SiStripQualityHotStripIdentifierRoot::algoEndJob() {
0254   //Clear map
0255   ClusterPositionHistoMap.clear();
0256 }
0257
0258 void SiStripQualityHotStripIdentifierRoot::bookHistos() {
0259   edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0260       << " [SiStripQualityHotStripIdentifierRoot::bookHistos] " << dirpath << std::endl;
0261   std::vector<MonitorElement*> MEs = dqmStore_->getAllContents(dirpath);
0262   //"DQMData/Run 50908/SiStrip/MechanicalView");
0263   //                                   "/DQMData/Run 50908/SiStrip/Run summary/MechanicalView/TID/side_2/wheel_3/ring_2/mono_modules/module_402676874");
0264   edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0265       << " [SiStripQualityHotStripIdentifierRoot::bookHistos] vector size " << MEs.size() << std::endl;
0266
0267   std::vector<MonitorElement*>::const_iterator iter = MEs.begin();
0268   std::vector<MonitorElement*>::const_iterator iterEnd = MEs.end();
0269
0270   edm::ParameterSet parameters = conf_.getParameter<edm::ParameterSet>("AlgoParameters");
0271   bool gotNentries = true;
0272   if (parameters.getUntrackedParameter<uint32_t>("NumberOfEvents", 0) == 0 &&
0273       parameters.getUntrackedParameter<double>("OccupancyThreshold", 0) != 0)
0274     gotNentries = false;
0275   edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0276       << "[SiStripQualityHotStripIdentifierRoot::bookHistos]  gotNentries flag " << gotNentries << " number of "
0277       << parameters.getUntrackedParameter<uint32_t>("NumberOfEvents", 0) << " occup "
0278       << parameters.getUntrackedParameter<double>("OccupancyThreshold", 0) << " OccupancyHisto"
0279       << parameters.getUntrackedParameter<std::string>("OccupancyHisto") << std::endl;
0280
0281   // Check Number of Events
0282   for (; iter != iterEnd; ++iter) {
0283     std::string me_name = (*iter)->getName();
0284
0285     if (!gotNentries && strstr(me_name.c_str(), "TotalNumberOfCluster__T") != nullptr &&
0286         strstr(me_name.c_str(), "Profile") == nullptr) {
0287       TotNumberOfEvents = ((TH1F*)(*iter)->getTH1F())->GetEntries();
0288       MeanNumberOfCluster = ((TH1F*)(*iter)->getTH1F())->GetMean();
0289       edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0290           << "Total Number of Events: " << TotNumberOfEvents << std::endl;
0291
0292       gotNentries = true;
0293       edm::LogInfo("SiStripQualityHotStripIdentifierRoot")
0294           << "[SiStripQualityHotStripIdentifierRoot::bookHistos]  gotNentries flag " << gotNentries << std::endl;
0295       break;
0296     }
0297   }
0298   if (!gotNentries) {
0299     edm::LogWarning("SiStripQualityHotStripIdentifierRoot")
0300         << " [SiStripQualityHotStripIdentifierRoot::bookHistos] :: Histogram with to check # of events missing"
0301         << std::endl;
0302   }
0303   for (; iter != iterEnd; ++iter) {
0304     std::string me_name = (*iter)->getName();
0305
0306     if (strstr(me_name.c_str(), (parameters.getUntrackedParameter<std::string>("OccupancyHisto")).c_str()) == nullptr)
0307       continue;
0308
0309     unsigned int detid = 0;
0310     char title[128];
0311     sprintf(title, "%s", me_name.c_str());
0312     char* saveptr;
0313     char* ptr = strtok_r(title, "__", &saveptr);
0314     int c = 0;
0315     while (ptr != nullptr) {
0316       if (c == 2) {
0317         detid = atol(ptr);
0318         break;
0319       }
0320       ptr = strtok_r(nullptr, "_", &saveptr);
0321       c++;
0322     }
0323     LogDebug("SiStripQualityHotStripIdentifierRoot")
0324         << " [SiStripQualityHotStripIdentifierRoot::bookHistos] detid " << detid << std::endl;
0325
0326     ClusterPositionHistoMap[detid] = std::make_shared<TH1F>(*(*iter)->getTH1F());
0327   }
0328 }