From 76c3ff3ded5f8e6220ed8180f2610bd57fd9e9f1 Mon Sep 17 00:00:00 2001 From: mmusich Date: Mon, 6 Mar 2023 17:29:48 +0100 Subject: [PATCH 1/9] fix segmentation violation in SiStripHitEffFromCalibTree in ASAN_X --- .../plugins/SiStripHitEffFromCalibTree.cc | 30 +++++++++++-------- 1 file changed, 18 insertions(+), 12 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc index d51f2f36adb7b..7619a84bd51d0 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc @@ -119,6 +119,7 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { // to be used everywhere static constexpr int SiStripLayers = 22; + static constexpr double nBxInAnOrbit = 3565; edm::Service fs; SiStripDetInfo _detInfo; @@ -179,6 +180,8 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { vector layertotal_vsPU; vector layerfound_vsCM; vector layertotal_vsCM; + vector layerfound_vsBX; + vector layertotal_vsBX; int goodlayertotal[35]; int goodlayerfound[35]; int alllayertotal[35]; @@ -284,6 +287,11 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve layertotal_vsPU.push_back( fs->make(Form("layertotal_vsPU_layer_%i", (int)(ilayer)), lyrName.c_str(), 45, 0, 90)); + layerfound_vsBX.push_back(fs->make( + Form("foundVsBx_layer%i", (int)ilayer), Form("layer %i", (int)ilayer), nBxInAnOrbit, 0, nBxInAnOrbit)); + layertotal_vsBX.push_back(fs->make( + Form("totalVsBx_layer%i", (int)ilayer), Form("layer %i", (int)ilayer), nBxInAnOrbit, 0, nBxInAnOrbit)); + if (_useCM) { layerfound_vsCM.push_back( fs->make(Form("layerfound_vsCM_layer_%i", (int)(ilayer)), lyrName.c_str(), 20, 0, 400)); @@ -1376,27 +1384,25 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { nLayers = 20; for (unsigned int ilayer = 1; ilayer < nLayers; ilayer++) { - TH1F* hfound = fs->make(Form("foundVsBx_layer%i", ilayer), Form("layer %i", ilayer), 3565, 0, 3565); - TH1F* htotal = fs->make(Form("totalVsBx_layer%i", ilayer), Form("layer %i", ilayer), 3565, 0, 3565); - for (unsigned int ibx = 0; ibx < 3566; ibx++) { - hfound->SetBinContent(ibx, 1e-6); - htotal->SetBinContent(ibx, 1); + layerfound_vsBX[ilayer]->SetBinContent(ibx, 1e-6); + layertotal_vsBX[ilayer]->SetBinContent(ibx, 1); } + map >::iterator iterMapvsBx; for (iterMapvsBx = layerfound_perBx.begin(); iterMapvsBx != layerfound_perBx.end(); ++iterMapvsBx) - hfound->SetBinContent(iterMapvsBx->first, iterMapvsBx->second[ilayer]); + layerfound_vsBX[ilayer]->SetBinContent(iterMapvsBx->first, iterMapvsBx->second[ilayer]); for (iterMapvsBx = layertotal_perBx.begin(); iterMapvsBx != layertotal_perBx.end(); ++iterMapvsBx) if (iterMapvsBx->second[ilayer] > 0) - htotal->SetBinContent(iterMapvsBx->first, iterMapvsBx->second[ilayer]); + layertotal_vsBX[ilayer]->SetBinContent(iterMapvsBx->first, iterMapvsBx->second[ilayer]); - hfound->Sumw2(); - htotal->Sumw2(); + layerfound_vsBX[ilayer]->Sumw2(); + layertotal_vsBX[ilayer]->Sumw2(); TGraphAsymmErrors* geff = fs->make(3564); geff->SetName(Form("effVsBx_layer%i", ilayer)); geff->SetTitle(fmt::format("Hit Efficiency vs bx - {}", ::layerName(ilayer, _showRings, nTEClayers)).c_str()); - geff->BayesDivide(hfound, htotal); + geff->BayesDivide(layerfound_vsBX[ilayer], layertotal_vsBX[ilayer]); //Average over trains TGraphAsymmErrors* geff_avg = fs->make(); @@ -1432,8 +1438,8 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { firstbx = ibx; } sum_bx += ibx; - found += hfound->GetBinContent(ibx); - total += htotal->GetBinContent(ibx); + found += layerfound_vsBX[ilayer]->GetBinContent(ibx); + total += layertotal_vsBX[ilayer]->GetBinContent(ibx); nbx++; previous_bx = ibx; From cf0a20d1b3bc89a5c802b02eef1a18a74825a0e7 Mon Sep 17 00:00:00 2001 From: mmusich Date: Sat, 25 Feb 2023 13:15:58 +0100 Subject: [PATCH 2/9] minor fixes to SiStrip HitEfficiency codes --- CalibTracker/SiStripHitEfficiency/plugins/HitEff.cc | 7 ++++--- CalibTracker/SiStripHitEfficiency/plugins/HitEff.h | 2 +- CalibTracker/SiStripHitEfficiency/test/testHitEffWorker.py | 2 +- 3 files changed, 6 insertions(+), 5 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/HitEff.cc b/CalibTracker/SiStripHitEfficiency/plugins/HitEff.cc index 82eb3ff03737a..c366f9966726b 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/HitEff.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/HitEff.cc @@ -92,6 +92,7 @@ HitEff::HitEff(const edm::ParameterSet& conf) chi2MeasurementEstimatorToken_(esConsumes(edm::ESInputTag("", "Chi2"))), propagatorToken_(esConsumes(edm::ESInputTag("", "PropagatorWithMaterial"))), conf_(conf) { + usesResource(TFileService::kSharedResource); compSettings = conf_.getUntrackedParameter("CompressionSettings", -1); layers = conf_.getParameter("Layer"); DEBUG = conf_.getParameter("Debug"); @@ -243,7 +244,7 @@ void HitEff::analyze(const edm::Event& e, const edm::EventSetup& es) { //e.getByLabel("siStripDigis", fedErrorIds ); e.getByToken(digis_token_, fedErrorIds); - ESHandle measurementTrackerHandle = es.getHandle(measurementTkToken_); + edm::ESHandle measurementTrackerHandle = es.getHandle(measurementTkToken_); edm::Handle measurementTrackerEvent; //e.getByLabel("MeasurementTrackerEvent", measurementTrackerEvent); @@ -298,7 +299,7 @@ void HitEff::analyze(const edm::Event& e, const edm::EventSetup& es) { #ifdef ExtendedCALIBTree //get dEdx info if available - Handle > dEdxUncalibHandle; + edm::Handle > dEdxUncalibHandle; if (e.getByLabel("dedxMedianCTF", dEdxUncalibHandle)) { const ValueMap dEdxTrackUncalib = *dEdxUncalibHandle.product(); @@ -311,7 +312,7 @@ void HitEff::analyze(const edm::Event& e, const edm::EventSetup& es) { } //get muon and ecal timing info if available - Handle muH; + edm::Handle muH; if (e.getByLabel("muonsWitht0Correction", muH)) { const MuonCollection& muonsT0 = *muH.product(); if (!muonsT0.empty()) { diff --git a/CalibTracker/SiStripHitEfficiency/plugins/HitEff.h b/CalibTracker/SiStripHitEfficiency/plugins/HitEff.h index bdf4e176ed44f..0d0188232c456 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/HitEff.h +++ b/CalibTracker/SiStripHitEfficiency/plugins/HitEff.h @@ -49,7 +49,7 @@ class TrackerTopology; -class HitEff : public edm::one::EDAnalyzer<> { +class HitEff : public edm::one::EDAnalyzer { public: explicit HitEff(const edm::ParameterSet& conf); ~HitEff() override = default; diff --git a/CalibTracker/SiStripHitEfficiency/test/testHitEffWorker.py b/CalibTracker/SiStripHitEfficiency/test/testHitEffWorker.py index 47ba605514de7..745431d77939c 100644 --- a/CalibTracker/SiStripHitEfficiency/test/testHitEffWorker.py +++ b/CalibTracker/SiStripHitEfficiency/test/testHitEffWorker.py @@ -10,7 +10,7 @@ options.parseArguments() process = cms.Process("HitEff") -process.load("Configuration/StandardSequences/MagneticField_cff") +process.load("Configuration.StandardSequences.MagneticField_cff") process.load("Configuration.StandardSequences.GeometryRecoDB_cff") process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff') process.load('FWCore.MessageService.MessageLogger_cfi') From 13a7f1b5fd20d9407be710db57851326b8b4e11f Mon Sep 17 00:00:00 2001 From: mmusich Date: Sat, 25 Feb 2023 15:29:17 +0100 Subject: [PATCH 3/9] apply CMS code rules in SiStripHitEffFromCalibTree --- .../plugins/SiStripHitEffFromCalibTree.cc | 347 +++++++++--------- 1 file changed, 175 insertions(+), 172 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc index 7619a84bd51d0..f1e400ad5d975 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc @@ -98,9 +98,13 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { ~SiStripHitEffFromCalibTree() override = default; private: + // overridden from ConditionDBWriter void algoBeginJob(const edm::EventSetup&) override; void algoEndJob() override; void algoAnalyze(const edm::Event& e, const edm::EventSetup& c) override; + std::unique_ptr getNewObject() override; + + // native methods void setBadComponents(int i, int component, SiStripQuality::BadComponent& BC, @@ -118,40 +122,38 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { TString getLayerSideName(Long_t k); // to be used everywhere - static constexpr int SiStripLayers = 22; - static constexpr double nBxInAnOrbit = 3565; + static constexpr int SiStripLayers_ = 22; + static constexpr double nBxInAnOrbit_ = 3565; edm::Service fs; - SiStripDetInfo _detInfo; - edm::FileInPath FileInPath_; + SiStripDetInfo detInfo_; + edm::FileInPath fileInPath_; SiStripQuality* quality_; - std::unique_ptr getNewObject() override; - TTree* CalibTree; - vector CalibTreeFilenames; - float threshold; - unsigned int nModsMin; - unsigned int doSummary; - string _badModulesFile; - bool _autoIneffModTagging; - unsigned int _clusterMatchingMethod; - float _ResXSig; - float _clusterTrajDist; - float _stripsApvEdge; - bool _useOnlyHighPurityTracks; - unsigned int _bunchx; - unsigned int _spaceBetweenTrains; - bool _useCM; - bool _showEndcapSides; - bool _showRings; - bool _showTOB6TEC9; - bool _showOnlyGoodModules; - float _tkMapMin; - float _effPlotMin; - TString _title; - - edm::ESGetToken _tkGeomToken; - edm::ESGetToken _tTopoToken; + const edm::ESGetToken tkGeomToken_; + const edm::ESGetToken tTopoToken_; + + TTree* calibTree_; + vector calibTreeFileNames_; + float threshold_; + unsigned int nModsMin_; + string badModulesFile_; + bool autoIneffModTagging_; + unsigned int clusterMatchingMethod_; + float resXSig_; + float clusterTrajDist_; + float stripsApvEdge_; + bool useOnlyHighPurityTracks_; + unsigned int bunchX_; + unsigned int spaceBetweenTrains_; + bool useCM_; + bool showEndcapSides_; + bool showRings_; + bool showTOB6TEC9_; + bool showOnlyGoodModules_; + float tkMapMin_; + float effPlotMin_; + TString title_; unsigned int nTEClayers; @@ -190,38 +192,38 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { }; SiStripHitEffFromCalibTree::SiStripHitEffFromCalibTree(const edm::ParameterSet& conf) - : ConditionDBWriter(conf), FileInPath_(SiStripDetInfoFileReader::kDefaultFile) { + : ConditionDBWriter(conf), + fileInPath_(SiStripDetInfoFileReader::kDefaultFile), + tkGeomToken_(esConsumes()), + tTopoToken_(esConsumes()) { usesResource(TFileService::kSharedResource); - CalibTreeFilenames = conf.getUntrackedParameter >("CalibTreeFilenames"); - threshold = conf.getParameter("Threshold"); - nModsMin = conf.getParameter("nModsMin"); - doSummary = conf.getParameter("doSummary"); - _badModulesFile = conf.getUntrackedParameter("BadModulesFile", ""); - _autoIneffModTagging = conf.getUntrackedParameter("AutoIneffModTagging", false); - _clusterMatchingMethod = conf.getUntrackedParameter("ClusterMatchingMethod", 0); - _ResXSig = conf.getUntrackedParameter("ResXSig", -1); - _clusterTrajDist = conf.getUntrackedParameter("ClusterTrajDist", 64.0); - _stripsApvEdge = conf.getUntrackedParameter("StripsApvEdge", 10.0); - _useOnlyHighPurityTracks = conf.getUntrackedParameter("UseOnlyHighPurityTracks", true); - _bunchx = conf.getUntrackedParameter("BunchCrossing", 0); - _spaceBetweenTrains = conf.getUntrackedParameter("SpaceBetweenTrains", 25); - _useCM = conf.getUntrackedParameter("UseCommonMode", false); - _showEndcapSides = conf.getUntrackedParameter("ShowEndcapSides", true); - _showRings = conf.getUntrackedParameter("ShowRings", false); - _showTOB6TEC9 = conf.getUntrackedParameter("ShowTOB6TEC9", false); - _showOnlyGoodModules = conf.getUntrackedParameter("ShowOnlyGoodModules", false); - _tkMapMin = conf.getUntrackedParameter("TkMapMin", 0.9); - _effPlotMin = conf.getUntrackedParameter("EffPlotMin", 0.9); - _title = conf.getParameter("Title"); - _tkGeomToken = esConsumes(); - _tTopoToken = esConsumes(); - _detInfo = SiStripDetInfoFileReader::read(FileInPath_.fullPath()); + calibTreeFileNames_ = conf.getUntrackedParameter >("CalibTreeFilenames"); + threshold_ = conf.getParameter("Threshold"); + nModsMin_ = conf.getParameter("nModsMin"); + badModulesFile_ = conf.getUntrackedParameter("BadModulesFile", ""); + autoIneffModTagging_ = conf.getUntrackedParameter("AutoIneffModTagging", false); + clusterMatchingMethod_ = conf.getUntrackedParameter("ClusterMatchingMethod", 0); + resXSig_ = conf.getUntrackedParameter("ResXSig", -1); + clusterTrajDist_ = conf.getUntrackedParameter("ClusterTrajDist", 64.0); + stripsApvEdge_ = conf.getUntrackedParameter("StripsApvEdge", 10.0); + useOnlyHighPurityTracks_ = conf.getUntrackedParameter("UseOnlyHighPurityTracks", true); + bunchX_ = conf.getUntrackedParameter("BunchCrossing", 0); + spaceBetweenTrains_ = conf.getUntrackedParameter("SpaceBetweenTrains", 25); + useCM_ = conf.getUntrackedParameter("UseCommonMode", false); + showEndcapSides_ = conf.getUntrackedParameter("ShowEndcapSides", true); + showRings_ = conf.getUntrackedParameter("ShowRings", false); + showTOB6TEC9_ = conf.getUntrackedParameter("ShowTOB6TEC9", false); + showOnlyGoodModules_ = conf.getUntrackedParameter("ShowOnlyGoodModules", false); + tkMapMin_ = conf.getUntrackedParameter("TkMapMin", 0.9); + effPlotMin_ = conf.getUntrackedParameter("EffPlotMin", 0.9); + title_ = conf.getParameter("Title"); + detInfo_ = SiStripDetInfoFileReader::read(fileInPath_.fullPath()); nTEClayers = 9; // number of wheels - if (_showRings) + if (showRings_) nTEClayers = 7; // number of rings - quality_ = new SiStripQuality(_detInfo); + quality_ = new SiStripQuality(detInfo_); } void SiStripHitEffFromCalibTree::algoBeginJob(const edm::EventSetup&) {} @@ -229,14 +231,14 @@ void SiStripHitEffFromCalibTree::algoBeginJob(const edm::EventSetup&) {} void SiStripHitEffFromCalibTree::algoEndJob() {} void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::EventSetup& c) { - const auto& tkgeom = c.getData(_tkGeomToken); - const auto& tTopo = c.getData(_tTopoToken); + const auto& tkgeom = c.getData(tkGeomToken_); + const auto& tTopo = c.getData(tTopoToken_); // read bad modules to mask ifstream badModules_file; set badModules_list; - if (!_badModulesFile.empty()) { - badModules_file.open(_badModulesFile.c_str()); + if (!badModulesFile_.empty()) { + badModules_file.open(badModulesFile_.c_str()); uint32_t badmodule_detid; int mods, fiber1, fiber2, fiber3; if (badModules_file.is_open()) { @@ -273,7 +275,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve TH1F* resolutionPlots[23]; for (Long_t ilayer = 0; ilayer < 23; ilayer++) { - std::string lyrName = ::layerName(ilayer, _showRings, nTEClayers); + std::string lyrName = ::layerName(ilayer, showRings_, nTEClayers); resolutionPlots[ilayer] = fs->make(Form("resol_layer_%i", (int)(ilayer)), lyrName.c_str(), 125, -125, 125); resolutionPlots[ilayer]->GetXaxis()->SetTitle("trajX-clusX [strip unit]"); @@ -288,11 +290,11 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve fs->make(Form("layertotal_vsPU_layer_%i", (int)(ilayer)), lyrName.c_str(), 45, 0, 90)); layerfound_vsBX.push_back(fs->make( - Form("foundVsBx_layer%i", (int)ilayer), Form("layer %i", (int)ilayer), nBxInAnOrbit, 0, nBxInAnOrbit)); + Form("foundVsBx_layer%i", (int)ilayer), Form("layer %i", (int)ilayer), nBxInAnOrbit_, 0, nBxInAnOrbit_)); layertotal_vsBX.push_back(fs->make( - Form("totalVsBx_layer%i", (int)ilayer), Form("layer %i", (int)ilayer), nBxInAnOrbit, 0, nBxInAnOrbit)); + Form("totalVsBx_layer%i", (int)ilayer), Form("layer %i", (int)ilayer), nBxInAnOrbit_, 0, nBxInAnOrbit_)); - if (_useCM) { + if (useCM_) { layerfound_vsCM.push_back( fs->make(Form("layerfound_vsCM_layer_%i", (int)(ilayer)), lyrName.c_str(), 20, 0, 400)); layertotal_vsCM.push_back( @@ -302,19 +304,19 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve layerfound[ilayer] = 0; } - if (!_autoIneffModTagging) - LOGPRINT << "A module is bad if efficiency < " << threshold << " and has at least " << nModsMin << " nModsMin."; + if (!autoIneffModTagging_) + LOGPRINT << "A module is bad if efficiency < " << threshold_ << " and has at least " << nModsMin_ << " nModsMin."; else - LOGPRINT << "A module is bad if the upper limit on the efficiency is < to the avg in the layer - " << threshold - << " and has at least " << nModsMin << " nModsMin."; + LOGPRINT << "A module is bad if the upper limit on the efficiency is < to the avg in the layer - " << threshold_ + << " and has at least " << nModsMin_ << " nModsMin."; unsigned int run, evt, bx{0}; double instLumi, PU; //Open the ROOT Calib Tree - for (unsigned int ifile = 0; ifile < CalibTreeFilenames.size(); ifile++) { - LOGPRINT << "Loading file: " << CalibTreeFilenames[ifile]; - TFile* CalibTreeFile = TFile::Open(CalibTreeFilenames[ifile].c_str(), "READ"); + for (unsigned int ifile = 0; ifile < calibTreeFileNames_.size(); ifile++) { + LOGPRINT << "Loading file: " << calibTreeFileNames_[ifile]; + TFile* CalibTreeFile = TFile::Open(calibTreeFileNames_[ifile].c_str(), "READ"); // Get event infos bool foundEventInfos = false; @@ -362,39 +364,39 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve // Get hit infos CalibTreeFile->cd("anEff"); - CalibTree = (TTree*)(gDirectory->Get("traj")); - - runLf = CalibTree->GetLeaf("run"); - evtLf = CalibTree->GetLeaf("event"); - TLeaf* BadLf = CalibTree->GetLeaf("ModIsBad"); - TLeaf* sistripLf = CalibTree->GetLeaf("SiStripQualBad"); - TLeaf* idLf = CalibTree->GetLeaf("Id"); - TLeaf* acceptLf = CalibTree->GetLeaf("withinAcceptance"); - TLeaf* layerLf = CalibTree->GetLeaf("layer"); - //TLeaf* nHitsLf = CalibTree->GetLeaf("nHits"); - TLeaf* highPurityLf = CalibTree->GetLeaf("highPurity"); - TLeaf* xLf = CalibTree->GetLeaf("TrajGlbX"); - TLeaf* yLf = CalibTree->GetLeaf("TrajGlbY"); - TLeaf* zLf = CalibTree->GetLeaf("TrajGlbZ"); - TLeaf* ResXSigLf = CalibTree->GetLeaf("ResXSig"); - TLeaf* TrajLocXLf = CalibTree->GetLeaf("TrajLocX"); - TLeaf* TrajLocYLf = CalibTree->GetLeaf("TrajLocY"); - TLeaf* ClusterLocXLf = CalibTree->GetLeaf("ClusterLocX"); - BunchLf = CalibTree->GetLeaf("bunchx"); - InstLumiLf = CalibTree->GetLeaf("instLumi"); - PULf = CalibTree->GetLeaf("PU"); + calibTree_ = (TTree*)(gDirectory->Get("traj")); + + runLf = calibTree_->GetLeaf("run"); + evtLf = calibTree_->GetLeaf("event"); + TLeaf* BadLf = calibTree_->GetLeaf("ModIsBad"); + TLeaf* sistripLf = calibTree_->GetLeaf("SiStripQualBad"); + TLeaf* idLf = calibTree_->GetLeaf("Id"); + TLeaf* acceptLf = calibTree_->GetLeaf("withinAcceptance"); + TLeaf* layerLf = calibTree_->GetLeaf("layer"); + //TLeaf* nHitsLf = calibTree_->GetLeaf("nHits"); + TLeaf* highPurityLf = calibTree_->GetLeaf("highPurity"); + TLeaf* xLf = calibTree_->GetLeaf("TrajGlbX"); + TLeaf* yLf = calibTree_->GetLeaf("TrajGlbY"); + TLeaf* zLf = calibTree_->GetLeaf("TrajGlbZ"); + TLeaf* ResXSigLf = calibTree_->GetLeaf("ResXSig"); + TLeaf* TrajLocXLf = calibTree_->GetLeaf("TrajLocX"); + TLeaf* TrajLocYLf = calibTree_->GetLeaf("TrajLocY"); + TLeaf* ClusterLocXLf = calibTree_->GetLeaf("ClusterLocX"); + BunchLf = calibTree_->GetLeaf("bunchx"); + InstLumiLf = calibTree_->GetLeaf("instLumi"); + PULf = calibTree_->GetLeaf("PU"); TLeaf* CMLf = nullptr; - if (_useCM) - CMLf = CalibTree->GetLeaf("commonMode"); + if (useCM_) + CMLf = calibTree_->GetLeaf("commonMode"); - int nevents = CalibTree->GetEntries(); + int nevents = calibTree_->GetEntries(); LOGPRINT << "Successfully loaded analyze function with " << nevents << " events!\n"; map, array >::iterator itEventInfos; //Loop through all of the events for (int j = 0; j < nevents; j++) { - CalibTree->GetEntry(j); + calibTree_->GetEntry(j); run = (unsigned int)runLf->GetValue(); evt = (unsigned int)evtLf->GetValue(); unsigned int isBad = (unsigned int)BadLf->GetValue(); @@ -403,7 +405,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve unsigned int accept = (unsigned int)acceptLf->GetValue(); unsigned int layer_wheel = (unsigned int)layerLf->GetValue(); unsigned int layer = layer_wheel; - if (_showRings && layer > 10) { // use rings instead of wheels + if (showRings_ && layer > 10) { // use rings instead of wheels if (layer < 14) layer = 10 + ((id >> 9) & 0x3); //TID 3 disks and also 3 rings -> use the same container else @@ -435,7 +437,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve PU = PULf->GetValue(); // branch not filled by default } int CM = -100; - if (_useCM) + if (useCM_) CM = CMLf->GetValue(); // Get infos from eventInfos if they exist @@ -451,19 +453,19 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve //We have two things we want to do, both an XY color plot, and the efficiency measurement //First, ignore anything that isn't in acceptance and isn't good quality - if (_bunchx > 0 && _bunchx != bx) + if (bunchX_ > 0 && bunchX_ != bx) continue; //if(quality == 1 || accept != 1 || nHits < 8) continue; if (accept != 1) continue; - if (_useOnlyHighPurityTracks && !highPurity) + if (useOnlyHighPurityTracks_ && !highPurity) continue; if (quality == 1) badquality = true; // don't compute efficiencies in modules from TOB6 and TEC9 - if (!_showTOB6TEC9 && (layer_wheel == 10 || layer_wheel == SiStripLayers)) + if (!showTOB6TEC9_ && (layer_wheel == 10 || layer_wheel == SiStripLayers_)) continue; // don't use bad modules given in the bad module list @@ -477,11 +479,11 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve bool badflag = false; // By default uses the old matching method - if (_ResXSig < 0) { + if (resXSig_ < 0) { if (isBad == 1) badflag = true; // isBad set to false in the tree when resxsig<999.0 } else { - if (isBad == 1 || resxsig > _ResXSig) + if (isBad == 1 || resxsig > resXSig_) badflag = true; } @@ -534,24 +536,24 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve int capv = -9; float stripInAPV = 64.; - if (_clusterMatchingMethod >= 1) { + if (clusterMatchingMethod_ >= 1) { badflag = false; // reset if (resxsig == 1000.0) { // default value when no cluster found in the module badflag = true; // consider the module inefficient in this case } else { - if (_clusterMatchingMethod == 2 || - _clusterMatchingMethod == 4) { // check the distance between cluster and trajectory position - if (abs(stripCluster - stripTrajMid) > _clusterTrajDist) + if (clusterMatchingMethod_ == 2 || + clusterMatchingMethod_ == 4) { // check the distance between cluster and trajectory position + if (abs(stripCluster - stripTrajMid) > clusterTrajDist_) badflag = true; } - if (_clusterMatchingMethod == 3 || - _clusterMatchingMethod == + if (clusterMatchingMethod_ == 3 || + clusterMatchingMethod_ == 4) { // cluster and traj have to be in the same APV (don't take edges into accounts) tapv = (int)stripTrajMid / 128; capv = (int)stripCluster / 128; stripInAPV = stripTrajMid - tapv * 128; - if (stripInAPV < _stripsApvEdge || stripInAPV > 128 - _stripsApvEdge) + if (stripInAPV < stripsApvEdge_ || stripInAPV > 128 - stripsApvEdge_) continue; if (tapv != capv) badflag = true; @@ -598,7 +600,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve layerfound_vsPU[layer]->Fill(PU); layertotal_vsPU[layer]->Fill(PU); - if (_useCM) { + if (useCM_) { if (!badflag) layerfound_vsCM[layer]->Fill(CM); layertotal_vsCM[layer]->Fill(CM); @@ -619,7 +621,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve goodlayerfound[layer + 3]++; goodlayertotal[layer + 3]++; } - } else if (layer > 13 && layer <= SiStripLayers) { + } else if (layer > 13 && layer <= SiStripLayers_) { if (((id >> 18) & 0x3) == 1) { if (!badflag) goodlayerfound[layer + 3]++; @@ -646,7 +648,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve alllayerfound[layer + 3]++; alllayertotal[layer + 3]++; } - } else if (layer > 13 && layer <= SiStripLayers) { + } else if (layer > 13 && layer <= SiStripLayers_) { if (((id >> 18) & 0x3) == 1) { if (!badflag) alllayerfound[layer + 3]++; @@ -662,13 +664,13 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve } // go to next CalibTreeFile makeHotColdMaps(); - makeTKMap(_autoIneffModTagging); + makeTKMap(autoIneffModTagging_); makeSQLite(); totalStatistics(); makeSummary(); makeSummaryVsBx(); makeSummaryVsLumi(); - if (_useCM) + if (useCM_) makeSummaryVsCM(); //////////////////////////////////////////////////////////////////////// @@ -785,7 +787,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve percentage += range; } if (percentage != 0) - percentage /= 128. * _detInfo.getNumberOfApvsAndStripLength(detid).first; + percentage /= 128. * detInfo_.getNumberOfApvsAndStripLength(detid).first; if (percentage > 1) edm::LogError("SiStripQualityStatistics") << "PROBLEM detid " << detid << " value " << percentage << std::endl; } @@ -888,7 +890,7 @@ void SiStripHitEffFromCalibTree::makeHotColdMaps() { //Already have access to the data as a private variable //Create all of the histograms in the TFileService TH2F* temph2; - for (Long_t maplayer = 1; maplayer <= SiStripLayers; maplayer++) { + for (Long_t maplayer = 1; maplayer <= SiStripLayers_; maplayer++) { //Initialize all of the histograms if (maplayer > 0 && maplayer <= 4) { //We are in the TIB @@ -962,7 +964,7 @@ void SiStripHitEffFromCalibTree::makeHotColdMaps() { HotColdMaps.push_back(temph2); } } - for (Long_t mylayer = 1; mylayer <= SiStripLayers; mylayer++) { + for (Long_t mylayer = 1; mylayer <= SiStripLayers_; mylayer++) { //Determine what kind of plot we want to write out //Loop through the entirety of each layer //Create an array of the histograms @@ -1009,14 +1011,14 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { LOGPRINT << "Entering TKMap generation!\n"; tkmap = new TrackerMap(" Detector Inefficiency "); tkmapbad = new TrackerMap(" Inefficient Modules "); - tkmapeff = new TrackerMap(_title.Data()); + tkmapeff = new TrackerMap(title_.Data()); tkmapnum = new TrackerMap(" Detector numerator "); tkmapden = new TrackerMap(" Detector denominator "); double myeff, mynum, myden, myeff_up; double layer_min_eff = 0; - for (Long_t i = 1; i <= SiStripLayers; i++) { + for (Long_t i = 1; i <= SiStripLayers_; i++) { //Loop over every layer, extracting the information from //the map of the efficiencies layertotal[i] = 0; @@ -1037,21 +1039,21 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { hEffInLayer->Fill(myeff); if (!autoTagging) { - if ((myden >= nModsMin) && (myeff < threshold)) { + if ((myden >= nModsMin_) && (myeff < threshold_)) { //We have a bad module, put it in the list! BadModules[(*ih).first] = myeff; tkmapbad->fillc((*ih).first, 255, 0, 0); - LOGPRINT << "Layer " << i << " (" << ::layerName(i, _showRings, nTEClayers) << ") module " << (*ih).first + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first << " efficiency: " << myeff << " , " << mynum << "/" << myden; } else { //Fill the bad list with empty results for every module tkmapbad->fillc((*ih).first, 255, 255, 255); } - if (myeff < threshold) - LOGPRINT << "Layer " << i << " (" << ::layerName(i, _showRings, nTEClayers) << ") module " << (*ih).first + if (myeff < threshold_) + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first << " efficiency: " << myeff << " , " << mynum << "/" << myden; - if (myden < nModsMin) { - LOGPRINT << "Layer " << i << " (" << ::layerName(i, _showRings, nTEClayers) << ") module " << (*ih).first + if (myden < nModsMin_) { + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first << " is under occupancy at " << myden; } } @@ -1072,8 +1074,8 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { hEffInLayer->GetXaxis()->SetRange(3, hEffInLayer->GetNbinsX() + 1); // Remove from the avg modules below 1% layer_min_eff = hEffInLayer->GetMean() - 2.5 * hEffInLayer->GetRMS(); // uses RMS in case the distribution is wide - if (threshold > 2.5 * hEffInLayer->GetRMS()) - layer_min_eff = hEffInLayer->GetMean() - threshold; // otherwise uses the parameter 'threshold' + if (threshold_ > 2.5 * hEffInLayer->GetRMS()) + layer_min_eff = hEffInLayer->GetMean() - threshold_; // otherwise uses the parameter 'threshold' LOGPRINT << "Layer " << i << " threshold for bad modules: <" << layer_min_eff << " (layer mean: " << hEffInLayer->GetMean() << " rms: " << hEffInLayer->GetRMS() << ")"; @@ -1089,7 +1091,7 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { myeff = 0; // upper limit on the efficiency myeff_up = TEfficiency::Bayesian(myden, mynum, .99, 1, 1, true); - if ((myden >= nModsMin) && (myeff_up < layer_min_eff)) { + if ((myden >= nModsMin_) && (myeff_up < layer_min_eff)) { //We have a bad module, put it in the list! BadModules[(*ih).first] = myeff; tkmapbad->fillc((*ih).first, 255, 0, 0); @@ -1098,10 +1100,10 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { tkmapbad->fillc((*ih).first, 255, 255, 255); } if (myeff_up < layer_min_eff + 0.08) // printing message also for modules slighly above (8%) the limit - LOGPRINT << "Layer " << i << " (" << ::layerName(i, _showRings, nTEClayers) << ") module " << (*ih).first + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first << " efficiency: " << myeff << " , " << mynum << "/" << myden << " , upper limit: " << myeff_up; - if (myden < nModsMin) { - LOGPRINT << "Layer " << i << " (" << ::layerName(i, _showRings, nTEClayers) << ") module " << (*ih).first + if (myden < nModsMin_) { + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first << " layer " << i << " is under occupancy at " << myden; } } @@ -1109,7 +1111,7 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { } tkmap->save(true, 0, 0, "SiStripHitEffTKMap.png"); tkmapbad->save(true, 0, 0, "SiStripHitEffTKMapBad.png"); - tkmapeff->save(true, _tkMapMin, 1., "SiStripHitEffTKMapEff.png"); + tkmapeff->save(true, tkMapMin_, 1., "SiStripHitEffTKMapEff.png"); tkmapnum->save(true, 0, 0, "SiStripHitEffTKMapNum.png"); tkmapden->save(true, 0, 0, "SiStripHitEffTKMapDen.png"); LOGPRINT << "Finished TKMap Generation\n"; @@ -1121,7 +1123,7 @@ void SiStripHitEffFromCalibTree::makeSQLite() { std::vector BadStripList; unsigned short NStrips; unsigned int id1; - std::unique_ptr pQuality = std::make_unique(_detInfo); + std::unique_ptr pQuality = std::make_unique(detInfo_); //This is the list of the bad strips, use to mask out entire APVs //Now simply go through the bad hit list and mask out things that //are bad! @@ -1129,7 +1131,7 @@ void SiStripHitEffFromCalibTree::makeSQLite() { for (it = BadModules.begin(); it != BadModules.end(); it++) { //We need to figure out how many strips are in this particular module //To Mask correctly! - NStrips = _detInfo.getNumberOfApvsAndStripLength((*it).first).first * 128; + NStrips = detInfo_.getNumberOfApvsAndStripLength((*it).first).first * 128; LOGPRINT << "Number of strips module " << (*it).first << " is " << NStrips; BadStripList.push_back(pQuality->encode(0, NStrips, 0)); //Now compact into a single bad module @@ -1157,9 +1159,9 @@ void SiStripHitEffFromCalibTree::totalStatistics() { subdettotal[i] = 0; } - for (Long_t i = 1; i <= SiStripLayers; i++) { + for (Long_t i = 1; i <= SiStripLayers_; i++) { layereff = double(layerfound[i]) / double(layertotal[i]); - LOGPRINT << "Layer " << i << " (" << ::layerName(i, _showRings, nTEClayers) << ") has total efficiency " << layereff + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") has total efficiency " << layereff << " " << layerfound[i] << "/" << layertotal[i]; totalfound += layerfound[i]; totaltotal += layertotal[i]; @@ -1196,11 +1198,11 @@ void SiStripHitEffFromCalibTree::makeSummary() { //setTDRStyle(); int nLayers = 34; - if (_showRings) + if (showRings_) nLayers = 30; - if (!_showEndcapSides) { - if (!_showRings) - nLayers = SiStripLayers; + if (!showEndcapSides_) { + if (!showRings_) + nLayers = SiStripLayers_; else nLayers = 20; } @@ -1226,7 +1228,7 @@ void SiStripHitEffFromCalibTree::makeSummary() { c7->SetGrid(); int nLayers_max = nLayers + 1; // barrel+endcap - if (!_showEndcapSides) + if (!showEndcapSides_) nLayers_max = 11; // barrel for (Long_t i = 1; i < nLayers_max; ++i) { LOGPRINT << "Fill only good modules layer " << i << ": S = " << goodlayerfound[i] @@ -1244,7 +1246,7 @@ void SiStripHitEffFromCalibTree::makeSummary() { } // endcap - merging sides - if (!_showEndcapSides) { + if (!showEndcapSides_) { for (Long_t i = 11; i < 14; ++i) { // TID disks LOGPRINT << "Fill only good modules layer " << i << ": S = " << goodlayerfound[i] + goodlayerfound[i + 3] << " B = " << goodlayertotal[i] + goodlayertotal[i + 3]; @@ -1301,12 +1303,12 @@ void SiStripHitEffFromCalibTree::makeSummary() { gr->SetLineColor(2); gr->SetLineWidth(4); gr->SetMarkerStyle(20); - gr->SetMinimum(_effPlotMin); + gr->SetMinimum(effPlotMin_); gr->SetMaximum(1.001); gr->GetYaxis()->SetTitle("Efficiency"); gStyle->SetTitleFillColor(0); gStyle->SetTitleBorderSize(0); - gr->SetTitle(_title); + gr->SetTitle(title_); gr2->GetXaxis()->SetLimits(0, nLayers); gr2->SetMarkerColor(1); @@ -1314,27 +1316,27 @@ void SiStripHitEffFromCalibTree::makeSummary() { gr2->SetLineColor(1); gr2->SetLineWidth(4); gr2->SetMarkerStyle(21); - gr2->SetMinimum(_effPlotMin); + gr2->SetMinimum(effPlotMin_); gr2->SetMaximum(1.001); gr2->GetYaxis()->SetTitle("Efficiency"); - gr2->SetTitle(_title); + gr2->SetTitle(title_); for (Long_t k = 1; k < nLayers + 1; k++) { TString label; - if (_showEndcapSides) + if (showEndcapSides_) label = getLayerSideName(k); else - label = ::layerName(k, _showRings, nTEClayers); - if (!_showTOB6TEC9) { + label = ::layerName(k, showRings_, nTEClayers); + if (!showTOB6TEC9_) { if (k == 10) label = ""; - if (!_showRings && k == nLayers) + if (!showRings_ && k == nLayers) label = ""; - if (!_showRings && _showEndcapSides && k == 25) + if (!showRings_ && showEndcapSides_ && k == 25) label = ""; } - if (!_showRings) { - if (_showEndcapSides) { + if (!showRings_) { + if (showEndcapSides_) { gr->GetXaxis()->SetBinLabel(((k + 1) * 100 + 2) / (nLayers)-4, label); gr2->GetXaxis()->SetBinLabel(((k + 1) * 100 + 2) / (nLayers)-4, label); } else { @@ -1342,7 +1344,7 @@ void SiStripHitEffFromCalibTree::makeSummary() { gr2->GetXaxis()->SetBinLabel((k + 1) * 100 / (nLayers)-6, label); } } else { - if (_showEndcapSides) { + if (showEndcapSides_) { gr->GetXaxis()->SetBinLabel((k + 1) * 100 / (nLayers)-4, label); gr2->GetXaxis()->SetBinLabel((k + 1) * 100 / (nLayers)-4, label); } else { @@ -1362,12 +1364,12 @@ void SiStripHitEffFromCalibTree::makeSummary() { overlay->SetFrameFillStyle(4000); overlay->Draw("same"); overlay->cd(); - if (!_showOnlyGoodModules) + if (!showOnlyGoodModules_) gr2->Draw("AP"); TLegend* leg = new TLegend(0.70, 0.27, 0.88, 0.40); leg->AddEntry(gr, "Good Modules", "p"); - if (!_showOnlyGoodModules) + if (!showOnlyGoodModules_) leg->AddEntry(gr2, "All Modules", "p"); leg->SetTextSize(0.020); leg->SetFillColor(0); @@ -1379,8 +1381,8 @@ void SiStripHitEffFromCalibTree::makeSummary() { void SiStripHitEffFromCalibTree::makeSummaryVsBx() { LOGPRINT << "Computing efficiency vs bx"; - unsigned int nLayers = SiStripLayers; - if (_showRings) + unsigned int nLayers = SiStripLayers_; + if (showRings_) nLayers = 20; for (unsigned int ilayer = 1; ilayer < nLayers; ilayer++) { @@ -1401,13 +1403,14 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { TGraphAsymmErrors* geff = fs->make(3564); geff->SetName(Form("effVsBx_layer%i", ilayer)); - geff->SetTitle(fmt::format("Hit Efficiency vs bx - {}", ::layerName(ilayer, _showRings, nTEClayers)).c_str()); + + geff->SetTitle(fmt::format("Hit Efficiency vs bx - {}", ::layerName(ilayer, showRings_, nTEClayers)).c_str()); geff->BayesDivide(layerfound_vsBX[ilayer], layertotal_vsBX[ilayer]); //Average over trains TGraphAsymmErrors* geff_avg = fs->make(); geff_avg->SetName(Form("effVsBxAvg_layer%i", ilayer)); - geff_avg->SetTitle(fmt::format("Hit Efficiency vs bx - {}", ::layerName(ilayer, _showRings, nTEClayers)).c_str()); + geff_avg->SetTitle(fmt::format("Hit Efficiency vs bx - {}", ::layerName(ilayer, showRings_, nTEClayers)).c_str()); geff_avg->SetMarkerStyle(20); int ibx = 0; int previous_bx = -80; @@ -1423,7 +1426,7 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { ibx = iterMapvsBx->first; delta_bx = ibx - previous_bx; // consider a new train - if (delta_bx > (int)_spaceBetweenTrains && nbx > 0 && total > 0) { + if (delta_bx > (int)spaceBetweenTrains_ && nbx > 0 && total > 0) { eff = found / (float)total; //LOGPRINT<<"new train "<SetPoint(ipt, sum_bx / nbx, eff); @@ -1455,8 +1458,8 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { } void SiStripHitEffFromCalibTree::computeEff(vector& vhfound, vector& vhtotal, string name) { - unsigned int nLayers = SiStripLayers; - if (_showRings) + unsigned int nLayers = SiStripLayers_; + if (showRings_) nLayers = 20; TH1F* hfound; @@ -1482,12 +1485,12 @@ void SiStripHitEffFromCalibTree::computeEff(vector& vhfound, vectorBayesDivide(hfound, htotal); if (name == "effVsLumi") geff->SetTitle( - fmt::format("Hit Efficiency vs inst. lumi. - {}", ::layerName(ilayer, _showRings, nTEClayers)).c_str()); + fmt::format("Hit Efficiency vs inst. lumi. - {}", ::layerName(ilayer, showRings_, nTEClayers)).c_str()); if (name == "effVsPU") - geff->SetTitle(fmt::format("Hit Efficiency vs pileup - {}", ::layerName(ilayer, _showRings, nTEClayers)).c_str()); + geff->SetTitle(fmt::format("Hit Efficiency vs pileup - {}", ::layerName(ilayer, showRings_, nTEClayers)).c_str()); if (name == "effVsCM") geff->SetTitle( - fmt::format("Hit Efficiency vs common Mode - {}", ::layerName(ilayer, _showRings, nTEClayers)).c_str()); + fmt::format("Hit Efficiency vs common Mode - {}", ::layerName(ilayer, showRings_, nTEClayers)).c_str()); geff->SetMarkerStyle(20); } } @@ -1501,8 +1504,8 @@ void SiStripHitEffFromCalibTree::makeSummaryVsLumi() { else { // from infos per hit - unsigned int nLayers = SiStripLayers; - if (_showRings) + unsigned int nLayers = SiStripLayers_; + if (showRings_) nLayers = 20; unsigned int nLayersForAvg = 0; float layerLumi = 0; @@ -1538,7 +1541,7 @@ void SiStripHitEffFromCalibTree::makeSummaryVsCM() { TString SiStripHitEffFromCalibTree::getLayerSideName(Long_t k) { TString layername = ""; TString ringlabel = "D"; - if (_showRings) + if (showRings_) ringlabel = "R"; if (k > 0 && k < 5) { layername = TString("TIB L") + k; @@ -1578,7 +1581,7 @@ std::unique_ptr SiStripHitEffFromCalibTree::getNewObject() { void SiStripHitEffFromCalibTree::setBadComponents( int i, int component, SiStripQuality::BadComponent& BC, std::stringstream ssV[4][19], int NBadComponent[4][19][4]) { - int napv = _detInfo.getNumberOfApvsAndStripLength(BC.detid).first; + int napv = detInfo_.getNumberOfApvsAndStripLength(BC.detid).first; ssV[i][component] << "\n\t\t " << BC.detid << " \t " << BC.BadModule << " \t " << ((BC.BadFibers) & 0x1) << " "; if (napv == 4) From 4ad02441cff5752b684d9a0b28e873fd9cb59e0f Mon Sep 17 00:00:00 2001 From: mmusich Date: Sat, 25 Feb 2023 15:38:28 +0100 Subject: [PATCH 4/9] remove unused methods --- .../plugins/SiStripHitEffFromCalibTree.cc | 6 ------ 1 file changed, 6 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc index f1e400ad5d975..ce7854d5566e9 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc @@ -99,8 +99,6 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { private: // overridden from ConditionDBWriter - void algoBeginJob(const edm::EventSetup&) override; - void algoEndJob() override; void algoAnalyze(const edm::Event& e, const edm::EventSetup& c) override; std::unique_ptr getNewObject() override; @@ -226,10 +224,6 @@ SiStripHitEffFromCalibTree::SiStripHitEffFromCalibTree(const edm::ParameterSet& quality_ = new SiStripQuality(detInfo_); } -void SiStripHitEffFromCalibTree::algoBeginJob(const edm::EventSetup&) {} - -void SiStripHitEffFromCalibTree::algoEndJob() {} - void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::EventSetup& c) { const auto& tkgeom = c.getData(tkGeomToken_); const auto& tTopo = c.getData(tTopoToken_); From b40891c216bc65b8f0e64bdd87218b8b28cdfba9 Mon Sep 17 00:00:00 2001 From: mmusich Date: Sat, 25 Feb 2023 16:12:07 +0100 Subject: [PATCH 5/9] range based loops in StripHitEffFromCalibTree --- .../plugins/SiStripHitEffFromCalibTree.cc | 85 ++++++++++--------- 1 file changed, 43 insertions(+), 42 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc index ce7854d5566e9..9eb2817e3aefb 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc @@ -251,8 +251,9 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve if (!badModules_list.empty()) LOGPRINT << "Remove additionnal bad modules from the analysis: "; set::iterator itBadMod; - for (itBadMod = badModules_list.begin(); itBadMod != badModules_list.end(); ++itBadMod) - LOGPRINT << " " << *itBadMod; + for (const auto& badMod : badModules_list) { + LOGPRINT << " " << badMod; + } // initialze counters and histos @@ -308,9 +309,9 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve double instLumi, PU; //Open the ROOT Calib Tree - for (unsigned int ifile = 0; ifile < calibTreeFileNames_.size(); ifile++) { - LOGPRINT << "Loading file: " << calibTreeFileNames_[ifile]; - TFile* CalibTreeFile = TFile::Open(calibTreeFileNames_[ifile].c_str(), "READ"); + for (const auto& calibTreeFileName : calibTreeFileNames_) { + LOGPRINT << "Loading file: " << calibTreeFileName; + TFile* CalibTreeFile = TFile::Open(calibTreeFileName.c_str(), "READ"); // Get event infos bool foundEventInfos = false; @@ -1020,12 +1021,11 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { TH1F* hEffInLayer = fs->make(Form("eff_layer%i", int(i)), Form("Module efficiency in layer %i", int(i)), 201, 0, 1.005); - map >::const_iterator ih; - for (ih = modCounter[i].begin(); ih != modCounter[i].end(); ih++) { + for (const auto& ih : modCounter[i]) { //We should be in the layer in question, and looping over all of the modules in said layer //Generate the list for the TKmap, and the bad module list - mynum = (double)(((*ih).second).second); - myden = (double)(((*ih).second).first); + mynum = (double)((ih.second).second); + myden = (double)((ih.second).first); if (myden > 0) myeff = mynum / myden; else @@ -1035,28 +1035,28 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { if (!autoTagging) { if ((myden >= nModsMin_) && (myeff < threshold_)) { //We have a bad module, put it in the list! - BadModules[(*ih).first] = myeff; - tkmapbad->fillc((*ih).first, 255, 0, 0); - LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first + BadModules[ih.first] = myeff; + tkmapbad->fillc(ih.first, 255, 0, 0); + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << ih.first << " efficiency: " << myeff << " , " << mynum << "/" << myden; } else { //Fill the bad list with empty results for every module - tkmapbad->fillc((*ih).first, 255, 255, 255); + tkmapbad->fillc(ih.first, 255, 255, 255); } if (myeff < threshold_) - LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << ih.first << " efficiency: " << myeff << " , " << mynum << "/" << myden; if (myden < nModsMin_) { - LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << ih.first << " is under occupancy at " << myden; } } //Put any module into the TKMap - tkmap->fill((*ih).first, 1. - myeff); - tkmapeff->fill((*ih).first, myeff); - tkmapnum->fill((*ih).first, mynum); - tkmapden->fill((*ih).first, myden); + tkmap->fill(ih.first, 1. - myeff); + tkmapeff->fill(ih.first, myeff); + tkmapnum->fill(ih.first, mynum); + tkmapden->fill(ih.first, myden); //Add the number of hits in the layer layertotal[i] += long(myden); @@ -1075,10 +1075,10 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { hEffInLayer->GetXaxis()->SetRange(1, hEffInLayer->GetNbinsX() + 1); - for (ih = modCounter[i].begin(); ih != modCounter[i].end(); ih++) { + for (const auto& ih : modCounter[i]) { // Second loop over modules to tag inefficient ones - mynum = (double)(((*ih).second).second); - myden = (double)(((*ih).second).first); + mynum = (double)((ih.second).second); + myden = (double)((ih.second).first); if (myden > 0) myeff = mynum / myden; else @@ -1087,17 +1087,17 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { myeff_up = TEfficiency::Bayesian(myden, mynum, .99, 1, 1, true); if ((myden >= nModsMin_) && (myeff_up < layer_min_eff)) { //We have a bad module, put it in the list! - BadModules[(*ih).first] = myeff; - tkmapbad->fillc((*ih).first, 255, 0, 0); + BadModules[ih.first] = myeff; + tkmapbad->fillc(ih.first, 255, 0, 0); } else { //Fill the bad list with empty results for every module - tkmapbad->fillc((*ih).first, 255, 255, 255); + tkmapbad->fillc(ih.first, 255, 255, 255); } if (myeff_up < layer_min_eff + 0.08) // printing message also for modules slighly above (8%) the limit - LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << ih.first << " efficiency: " << myeff << " , " << mynum << "/" << myden << " , upper limit: " << myeff_up; if (myden < nModsMin_) { - LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << (*ih).first + LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") module " << ih.first << " layer " << i << " is under occupancy at " << myden; } } @@ -1121,15 +1121,14 @@ void SiStripHitEffFromCalibTree::makeSQLite() { //This is the list of the bad strips, use to mask out entire APVs //Now simply go through the bad hit list and mask out things that //are bad! - map::const_iterator it; - for (it = BadModules.begin(); it != BadModules.end(); it++) { + for (const auto& it : BadModules) { //We need to figure out how many strips are in this particular module //To Mask correctly! - NStrips = detInfo_.getNumberOfApvsAndStripLength((*it).first).first * 128; - LOGPRINT << "Number of strips module " << (*it).first << " is " << NStrips; + NStrips = detInfo_.getNumberOfApvsAndStripLength(it.first).first * 128; + LOGPRINT << "Number of strips module " << it.first << " is " << NStrips; BadStripList.push_back(pQuality->encode(0, NStrips, 0)); //Now compact into a single bad module - id1 = (unsigned int)(*it).first; + id1 = (unsigned int)it.first; LOGPRINT << "ID1 shoudl match list of modules above " << id1; quality_->compact(id1, BadStripList); SiStripQuality::Range range(BadStripList.begin(), BadStripList.end()); @@ -1380,22 +1379,24 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { nLayers = 20; for (unsigned int ilayer = 1; ilayer < nLayers; ilayer++) { - for (unsigned int ibx = 0; ibx < 3566; ibx++) { + for (unsigned int ibx = 0; ibx <= nBxInAnOrbit_; ibx++) { layerfound_vsBX[ilayer]->SetBinContent(ibx, 1e-6); layertotal_vsBX[ilayer]->SetBinContent(ibx, 1); } - map >::iterator iterMapvsBx; - for (iterMapvsBx = layerfound_perBx.begin(); iterMapvsBx != layerfound_perBx.end(); ++iterMapvsBx) - layerfound_vsBX[ilayer]->SetBinContent(iterMapvsBx->first, iterMapvsBx->second[ilayer]); - for (iterMapvsBx = layertotal_perBx.begin(); iterMapvsBx != layertotal_perBx.end(); ++iterMapvsBx) - if (iterMapvsBx->second[ilayer] > 0) - layertotal_vsBX[ilayer]->SetBinContent(iterMapvsBx->first, iterMapvsBx->second[ilayer]); + for (const auto& iterMapvsBx : layerfound_perBx) { + layerfound_vsBX[ilayer]->SetBinContent(iterMapvsBx.first, iterMapvsBx.second[ilayer]); + } + for (const auto& iterMapvsBx : layertotal_perBx) { + if (iterMapvsBx.second[ilayer] > 0) { + layertotal_vsBX[ilayer]->SetBinContent(iterMapvsBx.first, iterMapvsBx.second[ilayer]); + } + } layerfound_vsBX[ilayer]->Sumw2(); layertotal_vsBX[ilayer]->Sumw2(); - TGraphAsymmErrors* geff = fs->make(3564); + TGraphAsymmErrors* geff = fs->make(nBxInAnOrbit_ - 1); geff->SetName(Form("effVsBx_layer%i", ilayer)); geff->SetTitle(fmt::format("Hit Efficiency vs bx - {}", ::layerName(ilayer, showRings_, nTEClayers)).c_str()); @@ -1416,8 +1417,8 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { int ipt = 0; float low, up, eff; int firstbx = 0; - for (iterMapvsBx = layertotal_perBx.begin(); iterMapvsBx != layertotal_perBx.end(); ++iterMapvsBx) { - ibx = iterMapvsBx->first; + for (const auto& iterMapvsBx : layertotal_perBx) { + ibx = iterMapvsBx.first; delta_bx = ibx - previous_bx; // consider a new train if (delta_bx > (int)spaceBetweenTrains_ && nbx > 0 && total > 0) { From ec62881a1d1d18cff5b956e495cab8ec66e94839 Mon Sep 17 00:00:00 2001 From: mmusich Date: Sun, 26 Feb 2023 14:05:11 +0100 Subject: [PATCH 6/9] SiStripHitEfficiencyHarvester: add generic function for summary plots of efficiency trends per layer vs variable (PU, inst. lumi, bx number) --- .../interface/SiStripHitEfficiencyHelpers.h | 14 ++ .../plugins/SiStripHitEfficiencyHarvester.cc | 120 +++++++++++------- 2 files changed, 91 insertions(+), 43 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h b/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h index 5f401cad531f1..76bdb9bd1de95 100644 --- a/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h +++ b/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h @@ -23,6 +23,20 @@ namespace { k_END_OF_LAYS_AND_RINGS = 35 }; + /* + * for the trend plots of efficiency vs some variable + */ + enum projections { k_vs_LUMI = 0, k_vs_PU = 1, k_vs_BX = 2, k_SIZE = 3 }; + + const std::array projFolder = {{"VsLumi", "VsPu", "VsBx"}}; + const std::array projFoundHisto = { + {"layerfound_vsLumi_layer_", "layerfound_vsPU_layer_", "foundVsBx_layer"}}; + const std::array projTotalHisto = { + {"layertotal_vsLumi_layer_", "layertotal_vsPU_layer_", "totalVsBx_layer"}}; + const std::array projTitle = {{"Inst Lumi", "Pile-Up", "Bunch Crossing"}}; + const std::array projXtitle = { + {"instantaneous luminosity [Hz/cm^{2}]", "Pile-Up events", "Bunch Crossing number"}}; + inline void replaceInString(std::string& str, const std::string& from, const std::string& to) { if (from.empty()) return; diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc index bca6d0de07e18..6f063234271c2 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc @@ -23,16 +23,17 @@ #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" //system includes -#include #include +#include #include // for std::accumulate +#include // ROOT includes -#include "TGraphAsymmErrors.h" #include "TCanvas.h" +#include "TEfficiency.h" +#include "TGraphAsymmErrors.h" #include "TLegend.h" #include "TStyle.h" -#include "TEfficiency.h" #include "TTree.h" // custom made printout @@ -92,9 +93,7 @@ class SiStripHitEfficiencyHarvester : public DQMEDHarvester { void makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker) const; template void setEffBinLabels(const T gr, const T gr2, const unsigned int nLayers) const; - void makeSummaryVsBX(DQMStore::IGetter& getter, TFileService& fs) const; - void makeSummaryVsLumi(DQMStore::IGetter& getter) const; - void makeSummaryVsCM(DQMStore::IGetter& getter, TFileService& fs) const; + void makeSummaryVsVariable(DQMStore::IGetter& getter, DQMStore::IBooker& booker, ::projections theProj) const; }; SiStripHitEfficiencyHarvester::SiStripHitEfficiencyHarvester(const edm::ParameterSet& conf) @@ -384,7 +383,7 @@ void SiStripHitEfficiencyHarvester::dqmEndJob(DQMStore::IBooker& booker, DQMStor } // loop on DetIds if (autoIneffModTagging_) { - for (Long_t i = 1; i <= k_LayersAtTECEnd; i++) { + for (unsigned int i = 1; i <= k_LayersAtTECEnd; i++) { //Compute threshold to use for each layer hEffInLayer[i]->getTH1()->GetXaxis()->SetRange( 3, hEffInLayer[i]->getNbinsX() + 1); // Remove from the avg modules below 1% @@ -489,14 +488,9 @@ void SiStripHitEfficiencyHarvester::dqmEndJob(DQMStore::IBooker& booker, DQMStor // make summary plots makeSummary(getter, booker); - makeSummaryVsLumi(getter); // TODO - - /* - if (!isAtPCL_) { - makeSummaryVsBX(getter, *fs); // TODO - makeSummaryVsCM(getter, *fs); // TODO - } - */ + makeSummaryVsVariable(getter, booker, projections::k_vs_LUMI); + makeSummaryVsVariable(getter, booker, projections::k_vs_PU); + makeSummaryVsVariable(getter, booker, projections::k_vs_BX); } void SiStripHitEfficiencyHarvester::printTotalStatistics( @@ -509,12 +503,12 @@ void SiStripHitEfficiencyHarvester::printTotalStatistics( int subdetfound[5]; int subdettotal[5]; - for (Long_t i = 1; i < 5; i++) { + for (unsigned int i = 1; i < 5; i++) { subdetfound[i] = 0; subdettotal[i] = 0; } - for (Long_t i = 1; i <= bounds::k_LayersAtTECEnd; i++) { + for (unsigned int i = 1; i <= bounds::k_LayersAtTECEnd; i++) { layereff = double(layerFound[i]) / double(layerTotal[i]); LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers_) << ") has total efficiency " << layereff << " " << layerFound[i] << "/" << layerTotal[i]; @@ -567,7 +561,6 @@ void SiStripHitEfficiencyHarvester::writeBadStripPayload(const SiStripQuality& q void SiStripHitEfficiencyHarvester::makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker) const { // use goodlayer_total/found and alllayer_total/found, collapse side and/or ring if needed - unsigned int nLayers{34}; // default if (showRings_) nLayers = 30; @@ -774,7 +767,7 @@ void SiStripHitEfficiencyHarvester::setEffBinLabels(const T gr, const T gr2, con << " showEndcapSides: " << showEndcapSides_ << " type of object is " << boost::typeindex::type_id().pretty_name(); - for (Long_t k = 1; k < nLayers + 1; k++) { + for (unsigned int k = 1; k < nLayers + 1; k++) { std::string label{}; if (showEndcapSides_) label = ::layerSideName(k, showRings_, nTEClayers_); @@ -816,45 +809,86 @@ void SiStripHitEfficiencyHarvester::setEffBinLabels(const T gr, const T gr2, con } } -// not yet implemented -void SiStripHitEfficiencyHarvester::makeSummaryVsBX(DQMStore::IGetter& getter, TFileService& fs) const { - // use found/totalVsBx_layer%i [0,23) -} +void SiStripHitEfficiencyHarvester::makeSummaryVsVariable(DQMStore::IGetter& getter, + DQMStore::IBooker& booker, + ::projections theProj) const { + std::vector effVsVariable; + effVsVariable.reserve(showRings_ ? 20 : 22); -// not yet impelement -void SiStripHitEfficiencyHarvester::makeSummaryVsCM(DQMStore::IGetter& getter, TFileService& fs) const {} + const auto& folderString = ::projFolder[theProj]; + const auto& foundHistoString = ::projFoundHisto[theProj]; + const auto& totalHistoString = ::projTotalHisto[theProj]; + const auto& titleString = ::projTitle[theProj]; + const auto& titleXString = ::projXtitle[theProj]; + + LogDebug("SiStripHitEfficiencyHarvester") + << " inside" << __PRETTY_FUNCTION__ << " from " << ::projFolder[theProj] << " " << __LINE__ << std::endl; -void SiStripHitEfficiencyHarvester::makeSummaryVsLumi(DQMStore::IGetter& getter) const { for (unsigned int iLayer = 1; iLayer != (showRings_ ? 20 : 22); ++iLayer) { - auto hfound = getter.get(fmt::format("{}/VsLumi/layerfound_vsLumi_layer_{}", inputFolder_, iLayer))->getTH1F(); - auto htotal = getter.get(fmt::format("{}/VsLumi/layertotal_vsLumi_layer_{}", inputFolder_, iLayer))->getTH1F(); + LogDebug("SiStripHitEfficiencyHarvester") + << "iLayer " << iLayer << " " << fmt::format("{}/{}/{}{}", inputFolder_, folderString, foundHistoString, iLayer) + << std::endl; + + const auto lyrName = ::layerName(iLayer, showRings_, nTEClayers_); + auto hfound = getter.get(fmt::format("{}/{}/{}{}", inputFolder_, folderString, foundHistoString, iLayer)); + auto htotal = getter.get(fmt::format("{}/{}/{}{}", inputFolder_, folderString, totalHistoString, iLayer)); if (hfound == nullptr or htotal == nullptr) { if (hfound == nullptr) edm::LogError("SiStripHitEfficiencyHarvester") - << fmt::format("{}/VsLumi/layerfound_vsLumi_layer_{}", inputFolder_, iLayer) << " was not found!"; + << fmt::format("{}/{}/{}{}", inputFolder_, folderString, foundHistoString, iLayer) << " was not found!"; if (htotal == nullptr) edm::LogError("SiStripHitEfficiencyHarvester") - << fmt::format("{}/VsLumi/layertotal_vsLumi_layer_{}", inputFolder_, iLayer) << " was not found!"; + << fmt::format("{}/{}/{}{}", inputFolder_, folderString, totalHistoString, iLayer) << " was not found!"; // no input histograms -> continue in the loop continue; } - if (!hfound->GetSumw2()) - hfound->Sumw2(); - if (!htotal->GetSumw2()) - htotal->Sumw2(); - for (Long_t i = 0; i != hfound->GetNbinsX() + 1; ++i) { - if (hfound->GetBinContent(i) == 0) - hfound->SetBinContent(i, 1e-6); - if (htotal->GetBinContent(i) == 0) - htotal->SetBinContent(i, 1); + // in order to display correct errors when taking the ratio + if (!hfound->getTH1F()->GetSumw2()) + hfound->getTH1F()->Sumw2(); + if (!htotal->getTH1F()->GetSumw2()) + htotal->getTH1F()->Sumw2(); + + // prevent dividing by 0 + for (int i = 0; i != hfound->getNbinsX() + 1; ++i) { + if (hfound->getBinContent(i) == 0) + hfound->setBinContent(i, 1e-6); + if (htotal->getBinContent(i) == 0) + htotal->setBinContent(i, 1); } + LogDebug("SiStripHitEfficiencyHarvester") << "Total hits for layer " << iLayer << " (" << folderString + << "): " << htotal->getEntries() << ", found " << hfound->getEntries(); + + booker.setCurrentFolder(fmt::format("{}/EfficiencySummary{}", inputFolder_, folderString)); + effVsVariable[iLayer] = booker.book1D( + fmt::sprintf("eff%sLayer%s", folderString, lyrName), + fmt::sprintf("Efficiency vs %s for layer %s;%s;SiStrip Hit efficiency", titleString, lyrName, titleXString), + hfound->getNbinsX(), + hfound->getAxisMin(), + hfound->getAxisMax()); + LogDebug("SiStripHitEfficiencyHarvester") - << "Total hits for layer " << iLayer << " (vs lumi): " << htotal->GetEntries() << ", found " - << hfound->GetEntries(); - } - // continue + << " bin 0 " << hfound->getAxisMin() << " bin last: " << hfound->getAxisMax() << std::endl; + + for (int i = 0; i != hfound->getNbinsX() + 1; ++i) { + const auto& den = htotal->getBinContent(i); + const auto& num = hfound->getBinContent(i); + + // fill all modules efficiency + // error on efficiency computed with binomial approximation + if (den > 0.) { + float eff = num / den; + float err_eff = (eff * (1 - eff)) / den; + effVsVariable[iLayer]->setBinContent(i, eff); + effVsVariable[iLayer]->setBinError(i, err_eff); + } + } + + // graphics adjustment + effVsVariable[iLayer]->getTH1F()->SetMinimum(tkMapMin_); + + } // loop on layers } namespace { From b227118b85ec0a2b2192e785dcb81362cc36cb21 Mon Sep 17 00:00:00 2001 From: mmusich Date: Wed, 15 Mar 2023 14:39:58 +0100 Subject: [PATCH 7/9] use uniform convention for data-members in SiStripHitEffFromCalibTree --- .../plugins/SiStripHitEffFromCalibTree.cc | 24 +++++++++---------- 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc index 9eb2817e3aefb..cf0bb46a04362 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEffFromCalibTree.cc @@ -120,7 +120,7 @@ class SiStripHitEffFromCalibTree : public ConditionDBWriter { TString getLayerSideName(Long_t k); // to be used everywhere - static constexpr int SiStripLayers_ = 22; + static constexpr int siStripLayers_ = 22; static constexpr double nBxInAnOrbit_ = 3565; edm::Service fs; @@ -460,7 +460,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve badquality = true; // don't compute efficiencies in modules from TOB6 and TEC9 - if (!showTOB6TEC9_ && (layer_wheel == 10 || layer_wheel == SiStripLayers_)) + if (!showTOB6TEC9_ && (layer_wheel == 10 || layer_wheel == siStripLayers_)) continue; // don't use bad modules given in the bad module list @@ -616,7 +616,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve goodlayerfound[layer + 3]++; goodlayertotal[layer + 3]++; } - } else if (layer > 13 && layer <= SiStripLayers_) { + } else if (layer > 13 && layer <= siStripLayers_) { if (((id >> 18) & 0x3) == 1) { if (!badflag) goodlayerfound[layer + 3]++; @@ -643,7 +643,7 @@ void SiStripHitEffFromCalibTree::algoAnalyze(const edm::Event& e, const edm::Eve alllayerfound[layer + 3]++; alllayertotal[layer + 3]++; } - } else if (layer > 13 && layer <= SiStripLayers_) { + } else if (layer > 13 && layer <= siStripLayers_) { if (((id >> 18) & 0x3) == 1) { if (!badflag) alllayerfound[layer + 3]++; @@ -885,7 +885,7 @@ void SiStripHitEffFromCalibTree::makeHotColdMaps() { //Already have access to the data as a private variable //Create all of the histograms in the TFileService TH2F* temph2; - for (Long_t maplayer = 1; maplayer <= SiStripLayers_; maplayer++) { + for (Long_t maplayer = 1; maplayer <= siStripLayers_; maplayer++) { //Initialize all of the histograms if (maplayer > 0 && maplayer <= 4) { //We are in the TIB @@ -959,7 +959,7 @@ void SiStripHitEffFromCalibTree::makeHotColdMaps() { HotColdMaps.push_back(temph2); } } - for (Long_t mylayer = 1; mylayer <= SiStripLayers_; mylayer++) { + for (Long_t mylayer = 1; mylayer <= siStripLayers_; mylayer++) { //Determine what kind of plot we want to write out //Loop through the entirety of each layer //Create an array of the histograms @@ -1013,7 +1013,7 @@ void SiStripHitEffFromCalibTree::makeTKMap(bool autoTagging = false) { double myeff, mynum, myden, myeff_up; double layer_min_eff = 0; - for (Long_t i = 1; i <= SiStripLayers_; i++) { + for (Long_t i = 1; i <= siStripLayers_; i++) { //Loop over every layer, extracting the information from //the map of the efficiencies layertotal[i] = 0; @@ -1152,7 +1152,7 @@ void SiStripHitEffFromCalibTree::totalStatistics() { subdettotal[i] = 0; } - for (Long_t i = 1; i <= SiStripLayers_; i++) { + for (Long_t i = 1; i <= siStripLayers_; i++) { layereff = double(layerfound[i]) / double(layertotal[i]); LOGPRINT << "Layer " << i << " (" << ::layerName(i, showRings_, nTEClayers) << ") has total efficiency " << layereff << " " << layerfound[i] << "/" << layertotal[i]; @@ -1195,7 +1195,7 @@ void SiStripHitEffFromCalibTree::makeSummary() { nLayers = 30; if (!showEndcapSides_) { if (!showRings_) - nLayers = SiStripLayers_; + nLayers = siStripLayers_; else nLayers = 20; } @@ -1374,7 +1374,7 @@ void SiStripHitEffFromCalibTree::makeSummary() { void SiStripHitEffFromCalibTree::makeSummaryVsBx() { LOGPRINT << "Computing efficiency vs bx"; - unsigned int nLayers = SiStripLayers_; + unsigned int nLayers = siStripLayers_; if (showRings_) nLayers = 20; @@ -1453,7 +1453,7 @@ void SiStripHitEffFromCalibTree::makeSummaryVsBx() { } void SiStripHitEffFromCalibTree::computeEff(vector& vhfound, vector& vhtotal, string name) { - unsigned int nLayers = SiStripLayers_; + unsigned int nLayers = siStripLayers_; if (showRings_) nLayers = 20; @@ -1499,7 +1499,7 @@ void SiStripHitEffFromCalibTree::makeSummaryVsLumi() { else { // from infos per hit - unsigned int nLayers = SiStripLayers_; + unsigned int nLayers = siStripLayers_; if (showRings_) nLayers = 20; unsigned int nLayersForAvg = 0; From 70466342ae1b803c852c21f9d752d0590cfd0980 Mon Sep 17 00:00:00 2001 From: mmusich Date: Wed, 15 Mar 2023 14:41:44 +0100 Subject: [PATCH 8/9] add provision to create TProfiles of efficiency and use Clopper-Pearson errors in trend plots --- .../interface/SiStripHitEfficiencyHelpers.h | 58 +++++++++++++++- .../plugins/SiStripHitEfficiencyHarvester.cc | 68 +++++++++++++++---- 2 files changed, 109 insertions(+), 17 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h b/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h index 76bdb9bd1de95..b25afcd32803c 100644 --- a/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h +++ b/CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h @@ -4,13 +4,22 @@ // A bunch of helper functions to deal with menial tasks in the // hit efficiency computation for the PCL workflow -#include "TString.h" -#include +// system includes #include +#include + +// user includes #include "DataFormats/TrackerCommon/interface/TrackerTopology.h" +#include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" #include "RecoLocalTracker/ClusterParameterEstimator/interface/StripClusterParameterEstimator.h" #include "TrackingTools/PatternTools/interface/TrajectoryMeasurement.h" +// ROOT includes +#include "TEfficiency.h" +#include "TProfile.h" +#include "TString.h" + namespace { enum bounds { @@ -215,5 +224,50 @@ namespace { return phi; } + inline TProfile* computeEff(const TH1F* num, const TH1F* denum, const std::string nameHist) { + std::string name = "eff_" + nameHist; + std::string title = "SiStrip Hit Efficiency" + std::string(num->GetTitle()); + TProfile* efficHist = new TProfile(name.c_str(), + title.c_str(), + denum->GetXaxis()->GetNbins(), + denum->GetXaxis()->GetXmin(), + denum->GetXaxis()->GetXmax()); + + for (int i = 1; i <= denum->GetNbinsX(); i++) { + double nNum = num->GetBinContent(i); + double nDenum = denum->GetBinContent(i); + if (nDenum == 0 || nNum == 0) { + continue; + } + if (nNum > nDenum) { + edm::LogWarning("SiStripHitEfficiencyHelpers") + << "Alert! specific bin's num is bigger than denum " << i << " " << nNum << " " << nDenum; + nNum = nDenum; // set the efficiency to 1 + } + const double effVal = nNum / nDenum; + efficHist->SetBinContent(i, effVal); + efficHist->SetBinEntries(i, 1); + const double errLo = TEfficiency::ClopperPearson((int)nDenum, (int)nNum, 0.683, false); + const double errUp = TEfficiency::ClopperPearson((int)nDenum, (int)nNum, 0.683, true); + const double errVal = (effVal - errLo > errUp - effVal) ? effVal - errLo : errUp - effVal; + efficHist->SetBinError(i, sqrt(effVal * effVal + errVal * errVal)); + + LogDebug("SiStripHitEfficiencyHelpers") << __PRETTY_FUNCTION__ << " " << nameHist << " bin:" << i + << " err:" << sqrt(effVal * effVal + errVal * errVal); + } + return efficHist; + } + + inline Measurement1D computeCPEfficiency(const double num, const double den) { + if (den > 0) { + const double effVal = num / den; + const double errLo = TEfficiency::ClopperPearson((int)den, (int)num, 0.683, false); + const double errUp = TEfficiency::ClopperPearson((int)den, (int)num, 0.683, true); + const double errVal = (effVal - errLo > errUp - effVal) ? effVal - errLo : errUp - effVal; + return Measurement1D(effVal, errVal); + } else { + return Measurement1D(0., 0.); + } + } } // namespace #endif diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc index 6f063234271c2..fd9ffc9eab6f3 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc @@ -90,7 +90,7 @@ class SiStripHitEfficiencyHarvester : public DQMEDHarvester { void printAndWriteBadModules(const SiStripQuality& quality, const SiStripDetInfo& detInfo) const; bool checkMapsValidity(const std::vector& maps, const std::string& type) const; unsigned int countTotalHits(const std::vector& maps); /* to check if TK was ON */ - void makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker) const; + void makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker, bool doProfiles = false) const; template void setEffBinLabels(const T gr, const T gr2, const unsigned int nLayers) const; void makeSummaryVsVariable(DQMStore::IGetter& getter, DQMStore::IBooker& booker, ::projections theProj) const; @@ -559,7 +559,9 @@ void SiStripHitEfficiencyHarvester::writeBadStripPayload(const SiStripQuality& q } } -void SiStripHitEfficiencyHarvester::makeSummary(DQMStore::IGetter& getter, DQMStore::IBooker& booker) const { +void SiStripHitEfficiencyHarvester::makeSummary(DQMStore::IGetter& getter, + DQMStore::IBooker& booker, + bool doProfiles) const { // use goodlayer_total/found and alllayer_total/found, collapse side and/or ring if needed unsigned int nLayers{34}; // default if (showRings_) @@ -657,6 +659,23 @@ void SiStripHitEfficiencyHarvester::makeSummary(DQMStore::IGetter& getter, DQMSt MonitorElement* h_eff_good = booker.book1D("eff_good", "Strip hit efficiency for good modules", nLayers + 1, 0, nLayers + 1); + if (doProfiles) { + // now do the profile + TProfile* profile_all = ::computeEff(found2->getTH1F(), all2->getTH1F(), "all"); + profile_all->SetMinimum(tkMapMin_); + profile_all->SetTitle("Strip hit efficiency for all modules"); + booker.bookProfile(profile_all->GetName(), profile_all); + + TProfile* profile_good = ::computeEff(found->getTH1F(), all->getTH1F(), "good"); + profile_good->SetMinimum(tkMapMin_); + profile_good->SetTitle("Strip hit efficiency for good modules"); + booker.bookProfile(profile_good->GetName(), profile_good); + + // clean the house + delete profile_all; + delete profile_good; + } + for (int i = 1; i < found->getNbinsX(); i++) { const auto& den_all = all2->getBinContent(i); const auto& num_all = found2->getBinContent(i); @@ -665,18 +684,26 @@ void SiStripHitEfficiencyHarvester::makeSummary(DQMStore::IGetter& getter, DQMSt // fill all modules efficiency if (den_all > 0.) { - float eff_all = num_all / den_all; - float err_eff_all = (eff_all * (1 - eff_all)) / den_all; - h_eff_all->setBinContent(i, eff_all); - h_eff_all->setBinError(i, err_eff_all); + // naive binomial errors + //float eff_all = num_all / den_all; + //float err_eff_all = (eff_all * (1 - eff_all)) / den_all; + + // use Clopper-Pearson errors + const auto& effPair_all = ::computeCPEfficiency(num_all, den_all); + h_eff_all->setBinContent(i, effPair_all.value()); + h_eff_all->setBinError(i, effPair_all.error()); } // fill good modules efficiency if (den_good > 0.) { - float eff_good = num_good / den_good; - float err_eff_good = (eff_good * (1 - eff_good)) / den_good; - h_eff_good->setBinContent(i, eff_good); - h_eff_good->setBinError(i, err_eff_good); + // naive binomial errors + //float eff_good = num_good / den_good; + //float err_eff_good = (eff_good * (1 - eff_good)) / den_good; + + // use Clopper-Pearson errors + const auto& effPair_good = ::computeCPEfficiency(num_good, den_good); + h_eff_good->setBinContent(i, effPair_good.value()); + h_eff_good->setBinError(i, effPair_good.error()); } } @@ -876,18 +903,29 @@ void SiStripHitEfficiencyHarvester::makeSummaryVsVariable(DQMStore::IGetter& get const auto& num = hfound->getBinContent(i); // fill all modules efficiency - // error on efficiency computed with binomial approximation if (den > 0.) { - float eff = num / den; - float err_eff = (eff * (1 - eff)) / den; - effVsVariable[iLayer]->setBinContent(i, eff); - effVsVariable[iLayer]->setBinError(i, err_eff); + const auto& effPair = ::computeCPEfficiency(num, den); + effVsVariable[iLayer]->setBinContent(i, effPair.value()); + effVsVariable[iLayer]->setBinError(i, effPair.error()); + + LogDebug("SiStripHitEfficiencyHarvester") + << __PRETTY_FUNCTION__ << " " << lyrName << " bin:" << i << " err:" << effPair.error() << std::endl; } } // graphics adjustment effVsVariable[iLayer]->getTH1F()->SetMinimum(tkMapMin_); + // now do the profile + TProfile* profile = ::computeEff(hfound->getTH1F(), htotal->getTH1F(), lyrName); + TString title = + fmt::sprintf("Efficiency vs %s for layer %s;%s;SiStrip Hit efficiency", titleString, lyrName, titleXString); + profile->SetMinimum(tkMapMin_); + + profile->SetTitle(title.Data()); + booker.bookProfile(profile->GetName(), profile); + + delete profile; } // loop on layers } From 8307e38118ad0e2f37e352fb3bb32ea5e40324ee Mon Sep 17 00:00:00 2001 From: mmusich Date: Wed, 15 Mar 2023 21:23:16 +0100 Subject: [PATCH 9/9] direct initialization of subdetfound and subdettotal --- .../plugins/SiStripHitEfficiencyHarvester.cc | 9 ++------- 1 file changed, 2 insertions(+), 7 deletions(-) diff --git a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc index fd9ffc9eab6f3..6758675c647ea 100644 --- a/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc +++ b/CalibTracker/SiStripHitEfficiency/plugins/SiStripHitEfficiencyHarvester.cc @@ -500,13 +500,8 @@ void SiStripHitEfficiencyHarvester::printTotalStatistics( int totalfound = 0; int totaltotal = 0; double layereff; - int subdetfound[5]; - int subdettotal[5]; - - for (unsigned int i = 1; i < 5; i++) { - subdetfound[i] = 0; - subdettotal[i] = 0; - } + int subdetfound[5] = {0, 0, 0, 0, 0}; + int subdettotal[5] = {0, 0, 0, 0, 0}; for (unsigned int i = 1; i <= bounds::k_LayersAtTECEnd; i++) { layereff = double(layerFound[i]) / double(layerTotal[i]);