File indexing completed on 2022-12-09 23:46:17
0001
0002 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0003 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
0004 #include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h"
0005 #include "CalibTracker/SiStripHitEfficiency/interface/SiStripHitEffData.h"
0006 #include "CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h"
0007 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0008 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
0009 #include "DQM/SiStripCommon/interface/TkHistoMap.h"
0010 #include "DQMServices/Core/interface/DQMEDHarvester.h"
0011 #include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h" /* for STRIPS_PER_APV*/
0012 #include "FWCore/Framework/interface/ESWatcher.h"
0013 #include "FWCore/Framework/interface/Frameworkfwd.h"
0014 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0015 #include "FWCore/ParameterSet/interface/ParameterDescription.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0018 #include "FWCore/ServiceRegistry/interface/Service.h"
0019 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0020 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0021 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0022 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0023 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0024
0025
0026 #include <sstream>
0027 #include <boost/type_index.hpp>
0028 #include <numeric> // for std::accumulate
0029
0030
0031 #include "TGraphAsymmErrors.h"
0032 #include "TCanvas.h"
0033 #include "TLegend.h"
0034 #include "TStyle.h"
0035 #include "TEfficiency.h"
0036 #include "TTree.h"
0037
0038
0039 #define LOGPRINT edm::LogPrint("SiStripHitEfficiencyHarvester")
0040
0041 class SiStripHitEfficiencyHarvester : public DQMEDHarvester {
0042 public:
0043 explicit SiStripHitEfficiencyHarvester(const edm::ParameterSet&);
0044 ~SiStripHitEfficiencyHarvester() override = default;
0045 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0046
0047 void endRun(edm::Run const&, edm::EventSetup const&) override;
0048 void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override;
0049
0050 private:
0051 SiStripHitEffData calibData_;
0052 const std::string inputFolder_;
0053 const bool isAtPCL_;
0054 const bool autoIneffModTagging_, doStoreOnDB_;
0055 const bool doStoreOnTree_;
0056 const bool showRings_, showEndcapSides_, showTOB6TEC9_, showOnlyGoodModules_;
0057 const unsigned int nTEClayers_;
0058 const double threshold_;
0059 const int nModsMin_;
0060 const float effPlotMin_;
0061 const double tkMapMin_;
0062 const std::string title_, record_;
0063
0064 const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0065 const edm::ESGetToken<TkDetMap, TrackerTopologyRcd> tkDetMapToken_;
0066 const edm::ESGetToken<SiStripQuality, SiStripQualityRcd> stripQualityToken_;
0067 const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkGeomToken_;
0068
0069 std::unique_ptr<TrackerTopology> tTopo_;
0070 std::unique_ptr<TkDetMap> tkDetMap_;
0071 std::unique_ptr<SiStripQuality> stripQuality_;
0072 std::vector<DetId> stripDetIds_;
0073
0074 int goodlayertotal[bounds::k_END_OF_LAYS_AND_RINGS];
0075 int goodlayerfound[bounds::k_END_OF_LAYS_AND_RINGS];
0076 int alllayertotal[bounds::k_END_OF_LAYS_AND_RINGS];
0077 int alllayerfound[bounds::k_END_OF_LAYS_AND_RINGS];
0078
0079
0080 TTree* tree;
0081 unsigned int t_DetId, t_found, t_total;
0082 unsigned char t_layer;
0083 bool t_isTaggedIneff;
0084 float t_threshold;
0085
0086 void writeBadStripPayload(const SiStripQuality& quality) const;
0087 void printTotalStatistics(const std::array<long, bounds::k_END_OF_LAYERS>& layerFound,
0088 const std::array<long, bounds::k_END_OF_LAYERS>& layerTotal) const;
0089 void printAndWriteBadModules(const SiStripQuality& quality, const SiStripDetInfo& detInfo) const;
0090 bool checkMapsValidity(const std::vector<MonitorElement*>& maps, const std::string& type) const;
0091 unsigned int countTotalHits(const std::vector<MonitorElement*>& maps);
0092 void makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker) const;
0093 template <typename T>
0094 void setEffBinLabels(const T gr, const T gr2, const unsigned int nLayers) const;
0095 void makeSummaryVsBX(DQMStore::IGetter& getter, TFileService& fs) const;
0096 void makeSummaryVsLumi(DQMStore::IGetter& getter) const;
0097 void makeSummaryVsCM(DQMStore::IGetter& getter, TFileService& fs) const;
0098 };
0099
0100 SiStripHitEfficiencyHarvester::SiStripHitEfficiencyHarvester(const edm::ParameterSet& conf)
0101 : inputFolder_(conf.getParameter<std::string>("inputFolder")),
0102 isAtPCL_(conf.getParameter<bool>("isAtPCL")),
0103 autoIneffModTagging_(conf.getUntrackedParameter<bool>("AutoIneffModTagging", false)),
0104 doStoreOnDB_(conf.getParameter<bool>("doStoreOnDB")),
0105 doStoreOnTree_(conf.getUntrackedParameter<bool>("doStoreOnTree")),
0106 showRings_(conf.getUntrackedParameter<bool>("ShowRings", false)),
0107 showEndcapSides_(conf.getUntrackedParameter<bool>("ShowEndcapSides", true)),
0108 showTOB6TEC9_(conf.getUntrackedParameter<bool>("ShowTOB6TEC9", false)),
0109 showOnlyGoodModules_(conf.getUntrackedParameter<bool>("ShowOnlyGoodModules", false)),
0110 nTEClayers_(showRings_ ? 7 : 9),
0111 threshold_(conf.getParameter<double>("Threshold")),
0112 nModsMin_(conf.getParameter<int>("nModsMin")),
0113 effPlotMin_(conf.getUntrackedParameter<double>("EffPlotMin", 0.9)),
0114 tkMapMin_(conf.getUntrackedParameter<double>("TkMapMin", 0.9)),
0115 title_(conf.getParameter<std::string>("Title")),
0116 record_(conf.getParameter<std::string>("Record")),
0117 tTopoToken_(esConsumes<edm::Transition::EndRun>()),
0118 tkDetMapToken_(esConsumes<edm::Transition::EndRun>()),
0119 stripQualityToken_(esConsumes<edm::Transition::EndRun>()),
0120 tkGeomToken_(esConsumes<edm::Transition::EndRun>()) {
0121
0122 for (int l = 0; l < bounds::k_END_OF_LAYS_AND_RINGS; l++) {
0123 goodlayertotal[l] = 0;
0124 goodlayerfound[l] = 0;
0125 alllayertotal[l] = 0;
0126 alllayerfound[l] = 0;
0127 }
0128 }
0129
0130 void SiStripHitEfficiencyHarvester::endRun(edm::Run const&, edm::EventSetup const& iSetup) {
0131 if (!tTopo_) {
0132 tTopo_ = std::make_unique<TrackerTopology>(iSetup.getData(tTopoToken_));
0133 }
0134 if (!tkDetMap_) {
0135 tkDetMap_ = std::make_unique<TkDetMap>(iSetup.getData(tkDetMapToken_));
0136 }
0137 if (!stripQuality_) {
0138 stripQuality_ = std::make_unique<SiStripQuality>(iSetup.getData(stripQualityToken_));
0139 }
0140 if (stripDetIds_.empty()) {
0141 const auto& tkGeom = iSetup.getData(tkGeomToken_);
0142 for (const auto& det : tkGeom.detUnits()) {
0143 if (dynamic_cast<const StripGeomDetUnit*>(det)) {
0144 stripDetIds_.push_back(det->geographicalId());
0145 }
0146 }
0147 }
0148 }
0149
0150 bool SiStripHitEfficiencyHarvester::checkMapsValidity(const std::vector<MonitorElement*>& maps,
0151 const std::string& type) const {
0152 std::vector<bool> isThere;
0153 isThere.reserve(maps.size());
0154 std::transform(maps.begin() + 1, maps.end(), std::back_inserter(isThere), [](auto& x) { return !(x == nullptr); });
0155
0156 int count{0};
0157 for (const auto& it : isThere) {
0158 count++;
0159 LogDebug("SiStripHitEfficiencyHarvester") << " layer: " << count << " " << it << std::endl;
0160 if (it)
0161 LogDebug("SiStripHitEfficiencyHarvester") << "resolving to " << maps[count]->getName() << std::endl;
0162 }
0163
0164
0165 bool areMapsAvailable{true};
0166 int layerCount{0};
0167 for (const auto& it : isThere) {
0168 layerCount++;
0169 if (!it) {
0170 edm::LogError("SiStripHitEfficiencyHarvester")
0171 << type << " TkHistoMap for layer " << layerCount << " was not found.\n -> Aborting!";
0172 areMapsAvailable = false;
0173 break;
0174 }
0175 }
0176 return areMapsAvailable;
0177 }
0178
0179 unsigned int SiStripHitEfficiencyHarvester::countTotalHits(const std::vector<MonitorElement*>& maps) {
0180 return std::accumulate(maps.begin() + 1, maps.end(), 0, [](unsigned int total, MonitorElement* item) {
0181 return total + item->getEntries();
0182 });
0183 }
0184
0185 void SiStripHitEfficiencyHarvester::dqmEndJob(DQMStore::IBooker& booker, DQMStore::IGetter& getter) {
0186 if (!isAtPCL_) {
0187 edm::Service<TFileService> fs;
0188 if (!fs.isAvailable()) {
0189 throw cms::Exception("BadConfig") << "TFileService unavailable: "
0190 << "please add it to config file";
0191 }
0192
0193 if (doStoreOnTree_) {
0194
0195 tree = fs->make<TTree>("ModEff", "ModEff");
0196 tree->Branch("DetId", &t_DetId, "DetId/i");
0197 tree->Branch("Layer", &t_layer, "Layer/b");
0198 tree->Branch("FoundHits", &t_found, "FoundHits/i");
0199 tree->Branch("AllHits", &t_total, "AllHits/i");
0200 tree->Branch("IsTaggedIneff", &t_isTaggedIneff, "IsTaggedIneff/O");
0201 tree->Branch("TagThreshold", &t_threshold, "TagThreshold/F");
0202 }
0203 }
0204
0205 if (!autoIneffModTagging_)
0206 LOGPRINT << "A module is bad if efficiency < " << threshold_ << " and has at least " << nModsMin_ << " nModsMin.";
0207 else
0208 LOGPRINT << "A module is bad if the upper limit on the efficiency is < to the avg in the layer - " << threshold_
0209 << " and has at least " << nModsMin_ << " nModsMin.";
0210
0211 auto h_module_total = std::make_unique<TkHistoMap>(tkDetMap_.get());
0212 h_module_total->loadTkHistoMap(fmt::format("{}/TkDetMaps", inputFolder_), "perModule_total");
0213 auto h_module_found = std::make_unique<TkHistoMap>(tkDetMap_.get());
0214 h_module_found->loadTkHistoMap(fmt::format("{}/TkDetMaps", inputFolder_), "perModule_found");
0215
0216
0217 const auto& totalMaps = h_module_total->getAllMaps();
0218 const auto& foundMaps = h_module_found->getAllMaps();
0219
0220 LogDebug("SiStripHitEfficiencyHarvester")
0221 << "totalMaps.size(): " << totalMaps.size() << " foundMaps.size() " << foundMaps.size() << std::endl;
0222
0223
0224 bool isTotalMapAvailable = this->checkMapsValidity(totalMaps, std::string("Total"));
0225 bool isFoundMapAvailable = this->checkMapsValidity(foundMaps, std::string("Found"));
0226
0227 LogDebug("SiStripHitEfficiencyHarvester")
0228 << "isTotalMapAvailable: " << isTotalMapAvailable << " isFoundMapAvailable " << isFoundMapAvailable << std::endl;
0229
0230
0231 if (!isTotalMapAvailable or !isFoundMapAvailable)
0232 return;
0233
0234 LogDebug("SiStripHitEfficiencyHarvester")
0235 << "Entries in total TkHistoMap for layer 3: " << h_module_total->getMap(3)->getEntries() << ", found "
0236 << h_module_found->getMap(3)->getEntries();
0237
0238
0239 const unsigned int totalHits = this->countTotalHits(totalMaps);
0240
0241
0242 for (size_t i = 1; i < totalMaps.size(); i++) {
0243 h_module_total->getMap(i)->setOption("colz");
0244 h_module_found->getMap(i)->setOption("colz");
0245 }
0246
0247
0248 booker.setCurrentFolder(inputFolder_);
0249
0250 std::vector<MonitorElement*> hEffInLayer(std::size_t(1), nullptr);
0251 hEffInLayer.reserve(bounds::k_END_OF_LAYERS);
0252 for (std::size_t i = 1; i != bounds::k_END_OF_LAYERS; ++i) {
0253 const auto lyrName = ::layerName(i, showRings_, nTEClayers_);
0254 hEffInLayer.push_back(booker.book1D(
0255 Form("eff_layer%i", int(i)), Form("Module efficiency in layer %s", lyrName.c_str()), 201, 0, 1.005));
0256 }
0257 std::array<long, bounds::k_END_OF_LAYERS> layerTotal{};
0258 std::array<long, bounds::k_END_OF_LAYERS> layerFound{};
0259 layerTotal.fill(0);
0260 layerFound.fill(0);
0261
0262
0263
0264
0265
0266 TrackerMap tkMap{" Detector Inefficiency "};
0267 TrackerMap tkMapBad{" Inefficient Modules "};
0268 TrackerMap tkMapEff{title_};
0269 TrackerMap tkMapNum{" Detector numerator "};
0270 TrackerMap tkMapDen{" Detector denominator "};
0271 std::map<unsigned int, double> badModules;
0272
0273
0274 const auto& EventStats = getter.get(fmt::format("{}/EventInfo/EventStats", inputFolder_));
0275 const int totalEvents = EventStats->getBinContent(1., 1.);
0276 calibData_.FEDErrorOccupancy = std::make_unique<TkHistoMap>(tkDetMap_.get());
0277 calibData_.FEDErrorOccupancy->loadTkHistoMap(fmt::format("{}/FEDErrorTkDetMaps", inputFolder_),
0278 "perModule_FEDErrors");
0279
0280
0281 calibData_.fillMapFromTkMap(totalEvents, 0.75, stripDetIds_);
0282
0283 for (const auto& [badId, fraction] : calibData_.fedErrorCounts) {
0284 LogDebug("SiStripHitEfficiencyHarvester")
0285 << __PRETTY_FUNCTION__ << " bad module from FEDError " << badId << "," << fraction << std::endl;
0286 }
0287
0288 for (auto det : stripDetIds_) {
0289 auto layer = ::checkLayer(det, tTopo_.get());
0290 const auto num = h_module_found->getValue(det);
0291 const auto denom = h_module_total->getValue(det);
0292 if (denom) {
0293
0294 if (stripQuality_->getBadApvs(det) == 0 && calibData_.checkFedError(det)) {
0295 const auto eff = num / denom;
0296 hEffInLayer[layer]->Fill(eff);
0297 if (!autoIneffModTagging_) {
0298 if ((denom >= nModsMin_) && (eff < threshold_)) {
0299
0300 badModules[det] = eff;
0301 tkMapBad.fillc(det, 255, 0, 0);
0302 LOGPRINT << "Layer " << layer << " (" << ::layerName(layer, showRings_, nTEClayers_) << ") module "
0303 << det.rawId() << " efficiency: " << eff << " , " << num << "/" << denom;
0304 } else {
0305
0306 tkMapBad.fillc(det, 255, 255, 255);
0307 }
0308 if (eff < threshold_)
0309 LOGPRINT << "Layer " << layer << " (" << ::layerName(layer, showRings_, nTEClayers_) << ") module "
0310 << det.rawId() << " efficiency: " << eff << " , " << num << "/" << denom;
0311
0312 if (denom < nModsMin_) {
0313 LOGPRINT << "Layer " << layer << " (" << ::layerName(layer, showRings_, nTEClayers_) << ") module "
0314 << det.rawId() << " is under occupancy at " << denom;
0315 }
0316
0317 if (doStoreOnTree_ && !isAtPCL_) {
0318 t_DetId = det.rawId();
0319 t_layer = layer;
0320 t_found = num;
0321 t_total = denom;
0322 t_isTaggedIneff = false;
0323 t_threshold = 0;
0324 tree->Fill();
0325 }
0326 }
0327
0328
0329 tkMap.fill(det, 1. - eff);
0330 tkMapEff.fill(det, eff);
0331 tkMapNum.fill(det, num);
0332 tkMapDen.fill(det, denom);
0333
0334 layerTotal[layer] += denom;
0335 layerFound[layer] += num;
0336
0337
0338
0339 if (layer <= bounds::k_LayersAtTOBEnd) {
0340 goodlayerfound[layer] += num;
0341 goodlayertotal[layer] += denom;
0342 } else if (layer > bounds::k_LayersAtTOBEnd && layer <= bounds::k_LayersAtTIDEnd) {
0343 if (tTopo_->tidSide(det) == 1) {
0344 goodlayerfound[layer] += num;
0345 goodlayertotal[layer] += denom;
0346 } else if (tTopo_->tidSide(det) == 2) {
0347 goodlayerfound[layer + 3] += num;
0348 goodlayertotal[layer + 3] += denom;
0349 }
0350 } else if (layer > bounds::k_LayersAtTIDEnd && layer <= bounds::k_LayersAtTECEnd) {
0351 if (tTopo_->tecSide(det) == 1) {
0352 goodlayerfound[layer + 3] += num;
0353 goodlayertotal[layer + 3] += denom;
0354 } else if (tTopo_->tecSide(det) == 2) {
0355 goodlayerfound[layer + 3 + nTEClayers_] += num;
0356 goodlayertotal[layer + 3 + nTEClayers_] += denom;
0357 }
0358 }
0359 }
0360
0361
0362 if (layer <= bounds::k_LayersAtTOBEnd) {
0363 alllayerfound[layer] += num;
0364 alllayertotal[layer] += denom;
0365 } else if (layer > bounds::k_LayersAtTOBEnd && layer <= bounds::k_LayersAtTIDEnd) {
0366 if (tTopo_->tidSide(det) == 1) {
0367 alllayerfound[layer] += num;
0368 alllayertotal[layer] += denom;
0369 } else if (tTopo_->tidSide(det) == 2) {
0370 alllayerfound[layer + 3] += num;
0371 alllayertotal[layer + 3] += denom;
0372 }
0373 } else if (layer > bounds::k_LayersAtTIDEnd && layer <= bounds::k_LayersAtTECEnd) {
0374 if (tTopo_->tecSide(det) == 1) {
0375 alllayerfound[layer + 3] += num;
0376 alllayertotal[layer + 3] += denom;
0377 } else if (tTopo_->tecSide(det) == 2) {
0378 alllayerfound[layer + 3 + nTEClayers_] += num;
0379 alllayertotal[layer + 3 + nTEClayers_] += denom;
0380 }
0381 }
0382
0383 }
0384 }
0385
0386 if (autoIneffModTagging_) {
0387 for (Long_t i = 1; i <= k_LayersAtTECEnd; i++) {
0388
0389 hEffInLayer[i]->getTH1()->GetXaxis()->SetRange(
0390 3, hEffInLayer[i]->getNbinsX() + 1);
0391 const double layer_min_eff = hEffInLayer[i]->getMean() - std::max(2.5 * hEffInLayer[i]->getRMS(), threshold_);
0392 LOGPRINT << "Layer " << i << " threshold for bad modules: <" << layer_min_eff
0393 << " (layer mean: " << hEffInLayer[i]->getMean() << " rms: " << hEffInLayer[i]->getRMS() << ")";
0394
0395 hEffInLayer[i]->getTH1()->GetXaxis()->SetRange(1, hEffInLayer[i]->getNbinsX() + 1);
0396
0397 for (auto det : stripDetIds_) {
0398
0399 if (stripQuality_->getBadApvs(det) == 0 && calibData_.checkFedError(det)) {
0400 const auto layer = ::checkLayer(det, tTopo_.get());
0401 if (layer == i) {
0402 const auto num = h_module_found->getValue(det);
0403 const auto denom = h_module_total->getValue(det);
0404 if (denom) {
0405 const auto eff = num / denom;
0406 const auto eff_up = TEfficiency::Bayesian(denom, num, .99, 1, 1, true);
0407
0408 if ((denom >= nModsMin_) && (eff_up < layer_min_eff)) {
0409
0410 badModules[det] = eff;
0411 tkMapBad.fillc(det, 255, 0, 0);
0412 if (!isAtPCL_ && doStoreOnTree_) {
0413 t_isTaggedIneff = true;
0414 }
0415 } else {
0416
0417 tkMapBad.fillc(det, 255, 255, 255);
0418 if (!isAtPCL_ && doStoreOnTree_) {
0419 t_isTaggedIneff = false;
0420 }
0421 }
0422 if (eff_up < layer_min_eff + 0.08) {
0423
0424 LOGPRINT << "Layer " << layer << " (" << ::layerName(layer, showRings_, nTEClayers_) << ") module "
0425 << det.rawId() << " efficiency: " << eff << " , " << num << "/" << denom
0426 << " , upper limit: " << eff_up;
0427 }
0428 if (denom < nModsMin_) {
0429 LOGPRINT << "Layer " << layer << " (" << ::layerName(layer, showRings_, nTEClayers_) << ") module "
0430 << det.rawId() << " layer " << layer << " is under occupancy at " << denom;
0431 }
0432
0433 if (!isAtPCL_ && doStoreOnTree_) {
0434 t_DetId = det.rawId();
0435 t_layer = layer;
0436 t_found = num;
0437 t_total = denom;
0438 t_threshold = layer_min_eff;
0439 tree->Fill();
0440 }
0441 }
0442 }
0443 }
0444 }
0445 }
0446 }
0447
0448 tkMap.save(true, 0, 0, "SiStripHitEffTKMap_NEW.png");
0449 tkMapBad.save(true, 0, 0, "SiStripHitEffTKMapBad_NEW.png");
0450 tkMapEff.save(true, tkMapMin_, 1., "SiStripHitEffTKMapEff_NEW.png");
0451 tkMapNum.save(true, 0, 0, "SiStripHitEffTKMapNum_NEW.png");
0452 tkMapDen.save(true, 0, 0, "SiStripHitEffTKMapDen_NEW.png");
0453
0454 const auto detInfo =
0455 SiStripDetInfoFileReader::read(edm::FileInPath{SiStripDetInfoFileReader::kDefaultFile}.fullPath());
0456 SiStripQuality pQuality{detInfo};
0457
0458
0459
0460 for (const auto it : badModules) {
0461 const auto det = it.first;
0462 std::vector<unsigned int> badStripList;
0463
0464
0465 const auto nStrips = detInfo.getNumberOfApvsAndStripLength(det).first * sistrip::STRIPS_PER_APV;
0466 LOGPRINT << "Number of strips module " << det << " is " << nStrips;
0467 badStripList.push_back(pQuality.encode(0, nStrips, 0));
0468
0469 LOGPRINT << "ID1 shoudl match list of modules above " << det;
0470 pQuality.compact(det, badStripList);
0471 pQuality.put(det, SiStripQuality::Range(badStripList.begin(), badStripList.end()));
0472 }
0473 pQuality.fillBadComponents();
0474 if (doStoreOnDB_) {
0475 if (totalHits > 0u) {
0476 writeBadStripPayload(pQuality);
0477 } else {
0478 edm::LogPrint("SiStripHitEfficiencyHarvester")
0479 << __PRETTY_FUNCTION__ << " There are no SiStrip hits for a valid measurement, skipping!";
0480 }
0481 } else {
0482 edm::LogInfo("SiStripHitEfficiencyHarvester") << "Will not produce payload!";
0483 }
0484
0485 printTotalStatistics(layerFound, layerTotal);
0486
0487
0488 printAndWriteBadModules(pQuality, detInfo);
0489
0490
0491 makeSummary(getter, booker);
0492 makeSummaryVsLumi(getter);
0493
0494
0495
0496
0497
0498
0499
0500 }
0501
0502 void SiStripHitEfficiencyHarvester::printTotalStatistics(
0503 const std::array<long, bounds::k_END_OF_LAYERS>& layerFound,
0504 const std::array<long, bounds::k_END_OF_LAYERS>& layerTotal) const {
0505
0506 int totalfound = 0;
0507 int totaltotal = 0;
0508 double layereff;
0509 int subdetfound[5];
0510 int subdettotal[5];
0511
0512 for (Long_t i = 1; i < 5; i++) {
0513 subdetfound[i] = 0;
0514 subdettotal[i] = 0;
0515 }
0516
0517 for (Long_t i = 1; i <= bounds::k_LayersAtTECEnd; i++) {
0518 layereff = double(layerFound[i]) / double(layerTotal[i]);
0519 LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers_) << ") has total efficiency "
0520 << layereff << " " << layerFound[i] << "/" << layerTotal[i];
0521 totalfound += layerFound[i];
0522 totaltotal += layerTotal[i];
0523 if (i <= bounds::k_LayersAtTIBEnd) {
0524 subdetfound[1] += layerFound[i];
0525 subdettotal[1] += layerTotal[i];
0526 }
0527 if (i > bounds::k_LayersAtTIBEnd && i <= bounds::k_LayersAtTOBEnd) {
0528 subdetfound[2] += layerFound[i];
0529 subdettotal[2] += layerTotal[i];
0530 }
0531 if (i > bounds::k_LayersAtTOBEnd && i <= bounds::k_LayersAtTIDEnd) {
0532 subdetfound[3] += layerFound[i];
0533 subdettotal[3] += layerTotal[i];
0534 }
0535 if (i > bounds::k_LayersAtTIDEnd) {
0536 subdetfound[4] += layerFound[i];
0537 subdettotal[4] += layerTotal[i];
0538 }
0539 }
0540
0541 LOGPRINT << "The total efficiency is " << double(totalfound) / double(totaltotal);
0542 LOGPRINT << " TIB: " << double(subdetfound[1]) / subdettotal[1] << " " << subdetfound[1] << "/"
0543 << subdettotal[1];
0544 LOGPRINT << " TOB: " << double(subdetfound[2]) / subdettotal[2] << " " << subdetfound[2] << "/"
0545 << subdettotal[2];
0546 LOGPRINT << " TID: " << double(subdetfound[3]) / subdettotal[3] << " " << subdetfound[3] << "/"
0547 << subdettotal[3];
0548 LOGPRINT << " TEC: " << double(subdetfound[4]) / subdettotal[4] << " " << subdetfound[4] << "/"
0549 << subdettotal[4];
0550 }
0551
0552 void SiStripHitEfficiencyHarvester::writeBadStripPayload(const SiStripQuality& quality) const {
0553 SiStripBadStrip pBadStrip{};
0554 const auto pQdvBegin = quality.getDataVectorBegin();
0555 for (auto rIt = quality.getRegistryVectorBegin(); rIt != quality.getRegistryVectorEnd(); ++rIt) {
0556 const auto range = SiStripBadStrip::Range(pQdvBegin + rIt->ibegin, pQdvBegin + rIt->iend);
0557 if (!pBadStrip.put(rIt->detid, range))
0558 edm::LogError("SiStripHitEfficiencyHarvester") << "detid already exists in SiStripBadStrip";
0559 }
0560 edm::Service<cond::service::PoolDBOutputService> poolDbService;
0561 if (poolDbService.isAvailable()) {
0562 poolDbService->writeOneIOV(pBadStrip, poolDbService->currentTime(), record_);
0563 } else {
0564 throw cms::Exception("PoolDBService required");
0565 }
0566 }
0567
0568 void SiStripHitEfficiencyHarvester::makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker) const {
0569
0570
0571 unsigned int nLayers{34};
0572 if (showRings_)
0573 nLayers = 30;
0574 if (!showEndcapSides_) {
0575 if (!showRings_)
0576 nLayers = 22;
0577 else
0578 nLayers = 20;
0579 }
0580
0581
0582 booker.setCurrentFolder(fmt::format("{}/EfficiencySummary", inputFolder_));
0583 MonitorElement* found = booker.book1D("found", "found", nLayers + 1, 0, nLayers + 1);
0584 MonitorElement* all = booker.book1D("all", "all", nLayers + 1, 0, nLayers + 1);
0585 MonitorElement* found2 = booker.book1D("found2", "found", nLayers + 1, 0, nLayers + 1);
0586 MonitorElement* all2 = booker.book1D("all2", "all2", nLayers + 1, 0, nLayers + 1);
0587
0588
0589 found->setBinContent(0, -1);
0590 all->setBinContent(0, 1);
0591
0592
0593 for (unsigned int i = 1; i < nLayers + 2; ++i) {
0594 found->setBinContent(i, 1e-6);
0595 all->setBinContent(i, 1);
0596 found2->setBinContent(i, 1e-6);
0597 all2->setBinContent(i, 1);
0598 }
0599
0600 TCanvas* c7 = new TCanvas("c7", " test ", 10, 10, 800, 600);
0601 c7->SetFillColor(0);
0602 c7->SetGrid();
0603
0604 unsigned int nLayers_max = nLayers + 1;
0605 if (!showEndcapSides_)
0606 nLayers_max = 11;
0607 for (unsigned int i = 1; i < nLayers_max; ++i) {
0608 LOGPRINT << "Fill only good modules layer " << i << ": S = " << goodlayerfound[i]
0609 << " B = " << goodlayertotal[i];
0610 if (goodlayertotal[i] > 5) {
0611 found->setBinContent(i, goodlayerfound[i]);
0612 all->setBinContent(i, goodlayertotal[i]);
0613 }
0614
0615 LOGPRINT << "Filling all modules layer " << i << ": S = " << alllayerfound[i] << " B = " << alllayertotal[i];
0616 if (alllayertotal[i] > 5) {
0617 found2->setBinContent(i, alllayerfound[i]);
0618 all2->setBinContent(i, alllayertotal[i]);
0619 }
0620 }
0621
0622
0623 if (!showEndcapSides_) {
0624 for (unsigned int i = 11; i < 14; ++i) {
0625 LOGPRINT << "Fill only good modules layer " << i << ": S = " << goodlayerfound[i] + goodlayerfound[i + 3]
0626 << " B = " << goodlayertotal[i] + goodlayertotal[i + 3];
0627 if (goodlayertotal[i] + goodlayertotal[i + 3] > 5) {
0628 found->setBinContent(i, goodlayerfound[i] + goodlayerfound[i + 3]);
0629 all->setBinContent(i, goodlayertotal[i] + goodlayertotal[i + 3]);
0630 }
0631 LOGPRINT << "Filling all modules layer " << i << ": S = " << alllayerfound[i] + alllayerfound[i + 3]
0632 << " B = " << alllayertotal[i] + alllayertotal[i + 3];
0633 if (alllayertotal[i] + alllayertotal[i + 3] > 5) {
0634 found2->setBinContent(i, alllayerfound[i] + alllayerfound[i + 3]);
0635 all2->setBinContent(i, alllayertotal[i] + alllayertotal[i + 3]);
0636 }
0637 }
0638 for (unsigned int i = 17; i < 17 + nTEClayers_; ++i) {
0639 LOGPRINT << "Fill only good modules layer " << i - 3
0640 << ": S = " << goodlayerfound[i] + goodlayerfound[i + nTEClayers_]
0641 << " B = " << goodlayertotal[i] + goodlayertotal[i + nTEClayers_];
0642 if (goodlayertotal[i] + goodlayertotal[i + nTEClayers_] > 5) {
0643 found->setBinContent(i - 3, goodlayerfound[i] + goodlayerfound[i + nTEClayers_]);
0644 all->setBinContent(i - 3, goodlayertotal[i] + goodlayertotal[i + nTEClayers_]);
0645 }
0646 LOGPRINT << "Filling all modules layer " << i - 3
0647 << ": S = " << alllayerfound[i] + alllayerfound[i + nTEClayers_]
0648 << " B = " << alllayertotal[i] + alllayertotal[i + nTEClayers_];
0649 if (alllayertotal[i] + alllayertotal[i + nTEClayers_] > 5) {
0650 found2->setBinContent(i - 3, alllayerfound[i] + alllayerfound[i + nTEClayers_]);
0651 all2->setBinContent(i - 3, alllayertotal[i] + alllayertotal[i + nTEClayers_]);
0652 }
0653 }
0654 }
0655
0656 found->getTH1F()->Sumw2();
0657 all->getTH1F()->Sumw2();
0658
0659 found2->getTH1F()->Sumw2();
0660 all2->getTH1F()->Sumw2();
0661
0662 MonitorElement* h_eff_all =
0663 booker.book1D("eff_all", "Strip hit efficiency for all modules", nLayers + 1, 0, nLayers + 1);
0664 MonitorElement* h_eff_good =
0665 booker.book1D("eff_good", "Strip hit efficiency for good modules", nLayers + 1, 0, nLayers + 1);
0666
0667 for (int i = 1; i < found->getNbinsX(); i++) {
0668 const auto& den_all = all2->getBinContent(i);
0669 const auto& num_all = found2->getBinContent(i);
0670 const auto& den_good = all->getBinContent(i);
0671 const auto& num_good = found->getBinContent(i);
0672
0673
0674 if (den_all > 0.) {
0675 float eff_all = num_all / den_all;
0676 float err_eff_all = (eff_all * (1 - eff_all)) / den_all;
0677 h_eff_all->setBinContent(i, eff_all);
0678 h_eff_all->setBinError(i, err_eff_all);
0679 }
0680
0681
0682 if (den_good > 0.) {
0683 float eff_good = num_good / den_good;
0684 float err_eff_good = (eff_good * (1 - eff_good)) / den_good;
0685 h_eff_good->setBinContent(i, eff_good);
0686 h_eff_good->setBinError(i, err_eff_good);
0687 }
0688 }
0689
0690 h_eff_all->getTH1F()->SetMinimum(effPlotMin_);
0691 h_eff_good->getTH1F()->SetMinimum(effPlotMin_);
0692
0693
0694 this->setEffBinLabels(h_eff_all->getTH1F(), h_eff_good->getTH1F(), nLayers);
0695
0696 if (!isAtPCL_) {
0697
0698 edm::Service<TFileService> fs;
0699 if (!fs.isAvailable()) {
0700 throw cms::Exception("BadConfig") << "TFileService unavailable: "
0701 << "please add it to config file";
0702 }
0703
0704 TGraphAsymmErrors* gr = (*fs).make<TGraphAsymmErrors>(nLayers + 1);
0705 gr->SetName("eff_good");
0706 gr->BayesDivide(found->getTH1F(), all->getTH1F());
0707
0708 TGraphAsymmErrors* gr2 = (*fs).make<TGraphAsymmErrors>(nLayers + 1);
0709 gr2->SetName("eff_all");
0710 gr2->BayesDivide(found2->getTH1F(), all2->getTH1F());
0711
0712 for (unsigned int j = 0; j < nLayers + 1; j++) {
0713 gr->SetPointError(j, 0., 0., gr->GetErrorYlow(j), gr->GetErrorYhigh(j));
0714 gr2->SetPointError(j, 0., 0., gr2->GetErrorYlow(j), gr2->GetErrorYhigh(j));
0715 }
0716
0717 this->setEffBinLabels(gr, gr2, nLayers);
0718
0719 gr->GetXaxis()->SetLimits(0, nLayers);
0720 gr->SetMarkerColor(2);
0721 gr->SetMarkerSize(1.2);
0722 gr->SetLineColor(2);
0723 gr->SetLineWidth(4);
0724 gr->SetMarkerStyle(20);
0725 gr->SetMinimum(effPlotMin_);
0726 gr->SetMaximum(1.001);
0727 gr->GetYaxis()->SetTitle("Efficiency");
0728 gStyle->SetTitleFillColor(0);
0729 gStyle->SetTitleBorderSize(0);
0730 gr->SetTitle("SiStripHitEfficiency by Layer");
0731
0732 gr2->GetXaxis()->SetLimits(0, nLayers);
0733 gr2->SetMarkerColor(1);
0734 gr2->SetMarkerSize(1.2);
0735 gr2->SetLineColor(1);
0736 gr2->SetLineWidth(4);
0737 gr2->SetMarkerStyle(21);
0738 gr2->SetMinimum(effPlotMin_);
0739 gr2->SetMaximum(1.001);
0740 gr2->GetYaxis()->SetTitle("Efficiency");
0741 gr2->SetTitle("SiStripHitEfficiency by Layer");
0742
0743 gr->Draw("AP");
0744 gr->GetXaxis()->SetNdivisions(36);
0745
0746 c7->cd();
0747 TPad* overlay = new TPad("overlay", "", 0, 0, 1, 1);
0748 overlay->SetFillStyle(4000);
0749 overlay->SetFillColor(0);
0750 overlay->SetFrameFillStyle(4000);
0751 overlay->Draw("same");
0752 overlay->cd();
0753 if (!showOnlyGoodModules_)
0754 gr2->Draw("AP");
0755
0756 TLegend* leg = new TLegend(0.70, 0.27, 0.88, 0.40);
0757 leg->AddEntry(gr, "Good Modules", "p");
0758 if (!showOnlyGoodModules_)
0759 leg->AddEntry(gr2, "All Modules", "p");
0760 leg->SetTextSize(0.020);
0761 leg->SetFillColor(0);
0762 leg->Draw("same");
0763
0764 c7->SaveAs("Summary.png");
0765 c7->SaveAs("Summary.root");
0766 }
0767 }
0768
0769 template <typename T>
0770 void SiStripHitEfficiencyHarvester::setEffBinLabels(const T gr, const T gr2, const unsigned int nLayers) const {
0771 LogDebug("SiStripHitEfficiencyHarvester")
0772 << "nLayers = " << nLayers << " number of bins, gr1: " << gr->GetXaxis()->GetNbins()
0773 << " number of bins, gr2: " << gr2->GetXaxis()->GetNbins() << " showRings: " << showRings_
0774 << " showEndcapSides: " << showEndcapSides_ << " type of object is "
0775 << boost::typeindex::type_id<T>().pretty_name();
0776
0777 for (Long_t k = 1; k < nLayers + 1; k++) {
0778 std::string label{};
0779 if (showEndcapSides_)
0780 label = ::layerSideName(k, showRings_, nTEClayers_);
0781 else
0782 label = ::layerName(k, showRings_, nTEClayers_);
0783 if (!showTOB6TEC9_) {
0784 if (k == 10)
0785 label = "";
0786 if (!showRings_ && k == nLayers)
0787 label = "";
0788 if (!showRings_ && showEndcapSides_ && k == 25)
0789 label = "";
0790 }
0791
0792 int bin{-1};
0793 if constexpr (std::is_same_v<T, TGraphAsymmErrors*>) {
0794 edm::LogInfo("SiStripHitEfficiencyHarvester")
0795 << "class name: " << gr->ClassName() << " expected TGraphAsymErrors" << std::endl;
0796 if (!showRings_) {
0797 if (showEndcapSides_) {
0798 bin = (((k + 1) * 100 + 2) / (nLayers)-4);
0799 } else {
0800 bin = ((k + 1) * 100 / (nLayers)-6);
0801 }
0802 } else {
0803 if (showEndcapSides_) {
0804 bin = ((k + 1) * 100 / (nLayers)-4);
0805 } else {
0806 bin = ((k + 1) * 100 / (nLayers)-7);
0807 }
0808 }
0809 } else {
0810 edm::LogInfo("SiStripHitEfficiencyHarvester")
0811 << "class name: " << gr->ClassName() << " expected TH1F" << std::endl;
0812 bin = k;
0813 }
0814 gr->GetXaxis()->SetBinLabel(bin, label.data());
0815 gr2->GetXaxis()->SetBinLabel(bin, label.data());
0816 }
0817 }
0818
0819
0820 void SiStripHitEfficiencyHarvester::makeSummaryVsBX(DQMStore::IGetter& getter, TFileService& fs) const {
0821
0822 }
0823
0824
0825 void SiStripHitEfficiencyHarvester::makeSummaryVsCM(DQMStore::IGetter& getter, TFileService& fs) const {}
0826
0827 void SiStripHitEfficiencyHarvester::makeSummaryVsLumi(DQMStore::IGetter& getter) const {
0828 for (unsigned int iLayer = 1; iLayer != (showRings_ ? 20 : 22); ++iLayer) {
0829 auto hfound = getter.get(fmt::format("{}/VsLumi/layerfound_vsLumi_layer_{}", inputFolder_, iLayer))->getTH1F();
0830 auto htotal = getter.get(fmt::format("{}/VsLumi/layertotal_vsLumi_layer_{}", inputFolder_, iLayer))->getTH1F();
0831
0832 if (hfound == nullptr or htotal == nullptr) {
0833 if (hfound == nullptr)
0834 edm::LogError("SiStripHitEfficiencyHarvester")
0835 << fmt::format("{}/VsLumi/layerfound_vsLumi_layer_{}", inputFolder_, iLayer) << " was not found!";
0836 if (htotal == nullptr)
0837 edm::LogError("SiStripHitEfficiencyHarvester")
0838 << fmt::format("{}/VsLumi/layertotal_vsLumi_layer_{}", inputFolder_, iLayer) << " was not found!";
0839
0840 continue;
0841 }
0842
0843 if (!hfound->GetSumw2())
0844 hfound->Sumw2();
0845 if (!htotal->GetSumw2())
0846 htotal->Sumw2();
0847 for (Long_t i = 0; i != hfound->GetNbinsX() + 1; ++i) {
0848 if (hfound->GetBinContent(i) == 0)
0849 hfound->SetBinContent(i, 1e-6);
0850 if (htotal->GetBinContent(i) == 0)
0851 htotal->SetBinContent(i, 1);
0852 }
0853 LogDebug("SiStripHitEfficiencyHarvester")
0854 << "Total hits for layer " << iLayer << " (vs lumi): " << htotal->GetEntries() << ", found "
0855 << hfound->GetEntries();
0856 }
0857
0858 }
0859
0860 namespace {
0861 void setBadComponents(int i,
0862 int comp,
0863 const SiStripQuality::BadComponent& bc,
0864 std::stringstream ssV[4][19],
0865 int nBad[4][19][4],
0866 int nAPV) {
0867 ssV[i][comp] << "\n\t\t " << bc.detid << " \t " << bc.BadModule << " \t " << ((bc.BadFibers) & 0x1) << " ";
0868 if (nAPV == 4)
0869 ssV[i][comp] << "x " << ((bc.BadFibers >> 1) & 0x1);
0870
0871 if (nAPV == 6)
0872 ssV[i][comp] << ((bc.BadFibers >> 1) & 0x1) << " " << ((bc.BadFibers >> 2) & 0x1);
0873 ssV[i][comp] << " \t " << ((bc.BadApvs) & 0x1) << " " << ((bc.BadApvs >> 1) & 0x1) << " ";
0874 if (nAPV == 4)
0875 ssV[i][comp] << "x x " << ((bc.BadApvs >> 2) & 0x1) << " " << ((bc.BadApvs >> 3) & 0x1);
0876 if (nAPV == 6)
0877 ssV[i][comp] << ((bc.BadApvs >> 2) & 0x1) << " " << ((bc.BadApvs >> 3) & 0x1) << " " << ((bc.BadApvs >> 4) & 0x1)
0878 << " " << ((bc.BadApvs >> 5) & 0x1) << " ";
0879
0880 if (bc.BadApvs) {
0881 nBad[i][0][2] += ((bc.BadApvs >> 5) & 0x1) + ((bc.BadApvs >> 4) & 0x1) + ((bc.BadApvs >> 3) & 0x1) +
0882 ((bc.BadApvs >> 2) & 0x1) + ((bc.BadApvs >> 1) & 0x1) + ((bc.BadApvs) & 0x1);
0883 nBad[i][comp][2] += ((bc.BadApvs >> 5) & 0x1) + ((bc.BadApvs >> 4) & 0x1) + ((bc.BadApvs >> 3) & 0x1) +
0884 ((bc.BadApvs >> 2) & 0x1) + ((bc.BadApvs >> 1) & 0x1) + ((bc.BadApvs) & 0x1);
0885 }
0886 if (bc.BadFibers) {
0887 nBad[i][0][1] += ((bc.BadFibers >> 2) & 0x1) + ((bc.BadFibers >> 1) & 0x1) + ((bc.BadFibers) & 0x1);
0888 nBad[i][comp][1] += ((bc.BadFibers >> 2) & 0x1) + ((bc.BadFibers >> 1) & 0x1) + ((bc.BadFibers) & 0x1);
0889 }
0890 if (bc.BadModule) {
0891 nBad[i][0][0]++;
0892 nBad[i][comp][0]++;
0893 }
0894 }
0895 }
0896
0897 void SiStripHitEfficiencyHarvester::printAndWriteBadModules(const SiStripQuality& quality,
0898 const SiStripDetInfo& detInfo) const {
0899
0900
0901
0902 int nTkBadComp[4];
0903 int nBadComp[4][19][4];
0904
0905
0906
0907 std::stringstream ssV[4][19];
0908
0909 for (int i = 0; i < 4; ++i) {
0910 nTkBadComp[i] = 0;
0911 for (int j = 0; j < 19; ++j) {
0912 ssV[i][j].str("");
0913 for (int k = 0; k < 4; ++k)
0914 nBadComp[i][j][k] = 0;
0915 }
0916 }
0917
0918 for (const auto& bc : quality.getBadComponentList()) {
0919
0920 if (bc.BadModule)
0921 nTkBadComp[0]++;
0922 if (bc.BadFibers)
0923 nTkBadComp[1] += ((bc.BadFibers >> 2) & 0x1) + ((bc.BadFibers >> 1) & 0x1) + ((bc.BadFibers) & 0x1);
0924 if (bc.BadApvs)
0925 nTkBadComp[2] += ((bc.BadApvs >> 5) & 0x1) + ((bc.BadApvs >> 4) & 0x1) + ((bc.BadApvs >> 3) & 0x1) +
0926 ((bc.BadApvs >> 2) & 0x1) + ((bc.BadApvs >> 1) & 0x1) + ((bc.BadApvs) & 0x1);
0927
0928 DetId det(bc.detid);
0929 if ((det.subdetId() >= SiStripSubdetector::TIB) && (det.subdetId() <= SiStripSubdetector::TEC)) {
0930 const auto nAPV = detInfo.getNumberOfApvsAndStripLength(det).first;
0931 switch (det.subdetId()) {
0932 case SiStripSubdetector::TIB:
0933 setBadComponents(0, tTopo_->tibLayer(det), bc, ssV, nBadComp, nAPV);
0934 break;
0935 case SiStripSubdetector::TID:
0936 setBadComponents(1,
0937 (tTopo_->tidSide(det) == 2 ? tTopo_->tidWheel(det) : tTopo_->tidWheel(det) + 3),
0938 bc,
0939 ssV,
0940 nBadComp,
0941 nAPV);
0942 break;
0943 case SiStripSubdetector::TOB:
0944 setBadComponents(2, tTopo_->tobLayer(det), bc, ssV, nBadComp, nAPV);
0945 break;
0946 case SiStripSubdetector::TEC:
0947 setBadComponents(3,
0948 (tTopo_->tecSide(det) == 2 ? tTopo_->tecWheel(det) : tTopo_->tecWheel(det) + 9),
0949 bc,
0950 ssV,
0951 nBadComp,
0952 nAPV);
0953 break;
0954 default:
0955 break;
0956 }
0957 }
0958 }
0959
0960 for (auto rp = quality.getRegistryVectorBegin(); rp != quality.getRegistryVectorEnd(); ++rp) {
0961 DetId det{rp->detid};
0962 int subdet = -999;
0963 int component = -999;
0964 switch (det.subdetId()) {
0965 case SiStripSubdetector::TIB:
0966 subdet = 0;
0967 component = tTopo_->tibLayer(det);
0968 break;
0969 case SiStripSubdetector::TID:
0970 subdet = 1;
0971 component = tTopo_->tidSide(det) == 2 ? tTopo_->tidWheel(det) : tTopo_->tidWheel(det) + 3;
0972 break;
0973 case SiStripSubdetector::TOB:
0974 subdet = 2;
0975 component = tTopo_->tobLayer(det);
0976 break;
0977 case SiStripSubdetector::TEC:
0978 subdet = 3;
0979 component = tTopo_->tecSide(det) == 2 ? tTopo_->tecWheel(det) : tTopo_->tecWheel(det) + 9;
0980 break;
0981 default:
0982 break;
0983 }
0984
0985 const auto pQdvBegin = quality.getDataVectorBegin();
0986 const auto sqrange = SiStripQuality::Range(pQdvBegin + rp->ibegin, pQdvBegin + rp->iend);
0987 float percentage = 0;
0988 for (int it = 0; it < sqrange.second - sqrange.first; it++) {
0989 unsigned int range = quality.decode(*(sqrange.first + it)).range;
0990 nTkBadComp[3] += range;
0991 nBadComp[subdet][0][3] += range;
0992 nBadComp[subdet][component][3] += range;
0993 percentage += range;
0994 }
0995 if (percentage != 0)
0996 percentage /= (sistrip::STRIPS_PER_APV * detInfo.getNumberOfApvsAndStripLength(det).first);
0997 if (percentage > 1)
0998 edm::LogError("SiStripHitEfficiencyHarvester") << "PROBLEM detid " << det.rawId() << " value " << percentage;
0999 }
1000
1001
1002 std::ostringstream ss;
1003 ss << "\n-----------------\nGlobal Info\n-----------------";
1004 ss << "\nBadComp \t Modules \tFibers "
1005 "\tApvs\tStrips\n----------------------------------------------------------------";
1006 ss << "\nTracker:\t\t" << nTkBadComp[0] << "\t" << nTkBadComp[1] << "\t" << nTkBadComp[2] << "\t" << nTkBadComp[3];
1007 ss << "\nTIB:\t\t\t" << nBadComp[0][0][0] << "\t" << nBadComp[0][0][1] << "\t" << nBadComp[0][0][2] << "\t"
1008 << nBadComp[0][0][3];
1009 ss << "\nTID:\t\t\t" << nBadComp[1][0][0] << "\t" << nBadComp[1][0][1] << "\t" << nBadComp[1][0][2] << "\t"
1010 << nBadComp[1][0][3];
1011 ss << "\nTOB:\t\t\t" << nBadComp[2][0][0] << "\t" << nBadComp[2][0][1] << "\t" << nBadComp[2][0][2] << "\t"
1012 << nBadComp[2][0][3];
1013 ss << "\nTEC:\t\t\t" << nBadComp[3][0][0] << "\t" << nBadComp[3][0][1] << "\t" << nBadComp[3][0][2] << "\t"
1014 << nBadComp[3][0][3];
1015 ss << "\n";
1016
1017 for (int i = 1; i < 5; ++i)
1018 ss << "\nTIB Layer " << i << " :\t\t" << nBadComp[0][i][0] << "\t" << nBadComp[0][i][1] << "\t" << nBadComp[0][i][2]
1019 << "\t" << nBadComp[0][i][3];
1020 ss << "\n";
1021 for (int i = 1; i < 4; ++i)
1022 ss << "\nTID+ Disk " << i << " :\t\t" << nBadComp[1][i][0] << "\t" << nBadComp[1][i][1] << "\t" << nBadComp[1][i][2]
1023 << "\t" << nBadComp[1][i][3];
1024 for (int i = 4; i < 7; ++i)
1025 ss << "\nTID- Disk " << i - 3 << " :\t\t" << nBadComp[1][i][0] << "\t" << nBadComp[1][i][1] << "\t"
1026 << nBadComp[1][i][2] << "\t" << nBadComp[1][i][3];
1027 ss << "\n";
1028 for (int i = 1; i < 7; ++i)
1029 ss << "\nTOB Layer " << i << " :\t\t" << nBadComp[2][i][0] << "\t" << nBadComp[2][i][1] << "\t" << nBadComp[2][i][2]
1030 << "\t" << nBadComp[2][i][3];
1031 ss << "\n";
1032 for (int i = 1; i < 10; ++i)
1033 ss << "\nTEC+ Disk " << i << " :\t\t" << nBadComp[3][i][0] << "\t" << nBadComp[3][i][1] << "\t" << nBadComp[3][i][2]
1034 << "\t" << nBadComp[3][i][3];
1035 for (int i = 10; i < 19; ++i)
1036 ss << "\nTEC- Disk " << i - 9 << " :\t\t" << nBadComp[3][i][0] << "\t" << nBadComp[3][i][1] << "\t"
1037 << nBadComp[3][i][2] << "\t" << nBadComp[3][i][3];
1038 ss << "\n";
1039
1040 ss << "\n----------------------------------------------------------------\n\t\t Detid \tModules Fibers "
1041 "Apvs\n----------------------------------------------------------------";
1042 for (int i = 1; i < 5; ++i)
1043 ss << "\nTIB Layer " << i << " :" << ssV[0][i].str();
1044 ss << "\n";
1045 for (int i = 1; i < 4; ++i)
1046 ss << "\nTID+ Disk " << i << " :" << ssV[1][i].str();
1047 for (int i = 4; i < 7; ++i)
1048 ss << "\nTID- Disk " << i - 3 << " :" << ssV[1][i].str();
1049 ss << "\n";
1050 for (int i = 1; i < 7; ++i)
1051 ss << "\nTOB Layer " << i << " :" << ssV[2][i].str();
1052 ss << "\n";
1053 for (int i = 1; i < 10; ++i)
1054 ss << "\nTEC+ Disk " << i << " :" << ssV[3][i].str();
1055 for (int i = 10; i < 19; ++i)
1056 ss << "\nTEC- Disk " << i - 9 << " :" << ssV[3][i].str();
1057
1058 LOGPRINT << ss.str();
1059
1060
1061 std::ofstream badModules;
1062 badModules.open("BadModules_NEW.log");
1063 badModules << "\n----------------------------------------------------------------\n\t\t Detid \tModules Fibers "
1064 "Apvs\n----------------------------------------------------------------";
1065 for (int i = 1; i < 5; ++i)
1066 badModules << "\nTIB Layer " << i << " :" << ssV[0][i].str();
1067 badModules << "\n";
1068 for (int i = 1; i < 4; ++i)
1069 badModules << "\nTID+ Disk " << i << " :" << ssV[1][i].str();
1070 for (int i = 4; i < 7; ++i)
1071 badModules << "\nTID- Disk " << i - 3 << " :" << ssV[1][i].str();
1072 badModules << "\n";
1073 for (int i = 1; i < 7; ++i)
1074 badModules << "\nTOB Layer " << i << " :" << ssV[2][i].str();
1075 badModules << "\n";
1076 for (int i = 1; i < 10; ++i)
1077 badModules << "\nTEC+ Disk " << i << " :" << ssV[3][i].str();
1078 for (int i = 10; i < 19; ++i)
1079 badModules << "\nTEC- Disk " << i - 9 << " :" << ssV[3][i].str();
1080 badModules.close();
1081 }
1082
1083 void SiStripHitEfficiencyHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
1084 edm::ParameterSetDescription desc;
1085 desc.add<std::string>("inputFolder", "AlCaReco/SiStripHitEfficiency");
1086 desc.add<bool>("isAtPCL", false);
1087 desc.add<bool>("doStoreOnDB", false);
1088 desc.add<std::string>("Record", "SiStripBadStrip");
1089 desc.add<double>("Threshold", 0.1);
1090 desc.add<std::string>("Title", "Hit Efficiency");
1091 desc.add<int>("nModsMin", 5);
1092 desc.addUntracked<bool>("doStoreOnTree", false);
1093 desc.addUntracked<bool>("AutoIneffModTagging", false);
1094 desc.addUntracked<double>("TkMapMin", 0.9);
1095 desc.addUntracked<double>("EffPlotMin", 0.9);
1096 desc.addUntracked<bool>("ShowRings", false);
1097 desc.addUntracked<bool>("ShowEndcapSides", true);
1098 desc.addUntracked<bool>("ShowTOB6TEC9", false);
1099 desc.addUntracked<bool>("ShowOnlyGoodModules", false);
1100 descriptions.addWithDefaultLabel(desc);
1101 }
1102
1103 #include "FWCore/Framework/interface/MakerMacros.h"
1104 DEFINE_FWK_MODULE(SiStripHitEfficiencyHarvester);