Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

SiStripHitEfficiency workflow in the Prompt Calibration Loop #37708

Merged
merged 6 commits into from May 4, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
Expand Up @@ -10,6 +10,10 @@
#include "DQMServices/Core/interface/DQMEDHarvester.h"
#include "FWCore/Framework/interface/ESWatcher.h"
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
#include "FWCore/ParameterSet/interface/ParameterDescription.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
#include "FWCore/ServiceRegistry/interface/Service.h"
#include "Geometry/CommonDetUnit/interface/GeomDet.h"
#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
Expand All @@ -30,11 +34,13 @@ class SiStripHitEfficiencyHarvester : public DQMEDHarvester {
public:
explicit SiStripHitEfficiencyHarvester(const edm::ParameterSet&);
~SiStripHitEfficiencyHarvester() override = default;
static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);

void endRun(edm::Run const&, edm::EventSetup const&) override;
void dqmEndJob(DQMStore::IBooker&, DQMStore::IGetter&) override;

private:
const bool isAtPCL_;
const bool showRings_, autoIneffModTagging_, doStoreOnDB_;
const unsigned int nTEClayers_;
const double threshold_;
Expand All @@ -55,14 +61,16 @@ class SiStripHitEfficiencyHarvester : public DQMEDHarvester {
void writeBadStripPayload(const SiStripQuality& quality) const;
void printTotalStatistics(const std::array<long, 23>& layerFound, const std::array<long, 23>& layerTotal) const;
void printAndWriteBadModules(const SiStripQuality& quality, const SiStripDetInfo& detInfo) const;
bool checkMapsValidity(const std::vector<MonitorElement*>& maps, const std::string& type) const;
void makeSummary(DQMStore::IGetter& getter, TFileService& fs) const;
void makeSummaryVsBX(DQMStore::IGetter& getter, TFileService& fs) const;
void makeSummaryVsLumi(DQMStore::IGetter& getter, TFileService& fs) const;
void makeSummaryVsLumi(DQMStore::IGetter& getter) const;
void makeSummaryVsCM(DQMStore::IGetter& getter, TFileService& fs) const;
};

SiStripHitEfficiencyHarvester::SiStripHitEfficiencyHarvester(const edm::ParameterSet& conf)
: showRings_(conf.getUntrackedParameter<bool>("ShowRings", false)),
: isAtPCL_(conf.getParameter<bool>("isAtPCL")),
showRings_(conf.getUntrackedParameter<bool>("ShowRings", false)),
autoIneffModTagging_(conf.getUntrackedParameter<bool>("AutoIneffModTagging", false)),
doStoreOnDB_(conf.getParameter<bool>("doStoreOnDB")),
nTEClayers_(showRings_ ? 7 : 9), // number of rings or wheels
Expand Down Expand Up @@ -96,28 +104,75 @@ void SiStripHitEfficiencyHarvester::endRun(edm::Run const&, edm::EventSetup cons
}
}

bool SiStripHitEfficiencyHarvester::checkMapsValidity(const std::vector<MonitorElement*>& maps,
const std::string& type) const {
std::vector<bool> isAvailable;
isAvailable.reserve(maps.size());
std::transform(
maps.begin() + 1, maps.end(), std::back_inserter(isAvailable), [](auto& x) { return !(x == nullptr); });

int count{0};
for (const auto& it : isAvailable) {
count++;
LogDebug("SiStripHitEfficiencyHarvester") << " layer: " << count << " " << it << std::endl;
if (it)
LogDebug("SiStripHitEfficiencyHarvester") << "resolving to " << maps[count]->getName() << std::endl;
}

// check on the input TkHistoMap
bool areMapsAvailable{true};
int layerCount{0};
for (const auto& it : isAvailable) {
layerCount++;
if (!it) {
edm::LogError("SiStripHitEfficiencyHarvester")
<< type << " TkHistoMap for layer " << layerCount << " was not found.\n -> Aborting!";
areMapsAvailable = false;
break;
}
}
return areMapsAvailable;
}

void SiStripHitEfficiencyHarvester::dqmEndJob(DQMStore::IBooker& booker, DQMStore::IGetter& getter) {
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.";

edm::Service<TFileService> fs;

auto h_module_total = std::make_unique<TkHistoMap>(tkDetMap_.get());
h_module_total->loadTkHistoMap("SiStrip/HitEfficiency", "perModule_total");
h_module_total->loadTkHistoMap("AlCaReco/SiStripHitEfficiency", "perModule_total");
auto h_module_found = std::make_unique<TkHistoMap>(tkDetMap_.get());
h_module_found->loadTkHistoMap("SiStrip/HitEfficiency", "perModule_found");
h_module_found->loadTkHistoMap("AlCaReco/SiStripHitEfficiency", "perModule_found");

// collect how many layers are missing
const auto& totalMaps = h_module_total->getAllMaps();
const auto& foundMaps = h_module_found->getAllMaps();

LogDebug("SiStripHitEfficiencyHarvester")
<< "totalMaps.size(): " << totalMaps.size() << " foundMaps.size() " << foundMaps.size() << std::endl;

// check on the input TkHistoMaps
bool isTotalMapAvailable = this->checkMapsValidity(totalMaps, std::string("Total"));
bool isFoundMapAvailable = this->checkMapsValidity(foundMaps, std::string("Found"));

LogDebug("SiStripHitEfficiencyHarvester")
<< "isTotalMapAvailable: " << isTotalMapAvailable << " isFoundMapAvailable " << isFoundMapAvailable << std::endl;

// no input TkHistoMaps -> early return
if (!isTotalMapAvailable or !isFoundMapAvailable)
return;

LogDebug("SiStripHitEfficiencyHarvester")
<< "Entries in total TkHistoMap for layer 3: " << h_module_total->getMap(3)->getEntries() << ", found "
<< h_module_found->getMap(3)->getEntries();

std::vector<TH1F*> hEffInLayer(std::size_t(1), nullptr);
std::vector<MonitorElement*> hEffInLayer(std::size_t(1), nullptr);
hEffInLayer.reserve(23);
for (std::size_t i = 1; i != 23; ++i) {
hEffInLayer.push_back(
fs->make<TH1F>(Form("eff_layer%i", int(i)), Form("Module efficiency in layer %i", int(i)), 201, 0, 1.005));
booker.book1D(Form("eff_layer%i", int(i)), Form("Module efficiency in layer %i", int(i)), 201, 0, 1.005));
}
std::array<long, 23> layerTotal{};
std::array<long, 23> layerFound{};
Expand Down Expand Up @@ -176,12 +231,13 @@ void SiStripHitEfficiencyHarvester::dqmEndJob(DQMStore::IBooker& booker, DQMStor
if (autoIneffModTagging_) {
for (Long_t i = 1; i <= 22; i++) {
//Compute threshold to use for each layer
hEffInLayer[i]->GetXaxis()->SetRange(3, hEffInLayer[i]->GetNbinsX() + 1); // Remove from the avg modules below 1%
const double layer_min_eff = hEffInLayer[i]->GetMean() - std::max(2.5 * hEffInLayer[i]->GetRMS(), threshold_);
hEffInLayer[i]->getTH1()->GetXaxis()->SetRange(
3, hEffInLayer[i]->getNbinsX() + 1); // Remove from the avg modules below 1%
const double layer_min_eff = hEffInLayer[i]->getMean() - std::max(2.5 * hEffInLayer[i]->getRMS(), threshold_);
LOGPRINT << "Layer " << i << " threshold for bad modules: <" << layer_min_eff
<< " (layer mean: " << hEffInLayer[i]->GetMean() << " rms: " << hEffInLayer[i]->GetRMS() << ")";
<< " (layer mean: " << hEffInLayer[i]->getMean() << " rms: " << hEffInLayer[i]->getRMS() << ")";

hEffInLayer[i]->GetXaxis()->SetRange(1, hEffInLayer[i]->GetNbinsX() + 1);
hEffInLayer[i]->getTH1()->GetXaxis()->SetRange(1, hEffInLayer[i]->getNbinsX() + 1);

for (auto det : stripDetIds_) {
const auto layer = ::checkLayer(det, tTopo_.get());
Expand Down Expand Up @@ -252,10 +308,14 @@ void SiStripHitEfficiencyHarvester::dqmEndJob(DQMStore::IBooker& booker, DQMStor
// << " lumiBlock " << e.luminosityBlock() << " time " << e.time().value() << "\n-----------------\n";
printAndWriteBadModules(pQuality, detInfo); // TODO

makeSummary(getter, *fs); // TODO
makeSummaryVsBX(getter, *fs); // TODO
makeSummaryVsLumi(getter, *fs); // TODO
makeSummaryVsCM(getter, *fs); // TODO
if (!isAtPCL_) {
edm::Service<TFileService> fs;
makeSummary(getter, *fs); // TODO
makeSummaryVsBX(getter, *fs); // TODO
makeSummaryVsCM(getter, *fs); // TODO
}

makeSummaryVsLumi(getter); // TODO
}

void SiStripHitEfficiencyHarvester::printTotalStatistics(const std::array<long, 23>& layerFound,
Expand Down Expand Up @@ -331,10 +391,24 @@ void SiStripHitEfficiencyHarvester::makeSummaryVsBX(DQMStore::IGetter& getter, T
// use found/totalVsBx_layer%i [0,23)
}

void SiStripHitEfficiencyHarvester::makeSummaryVsLumi(DQMStore::IGetter& getter, TFileService& fs) const {
void SiStripHitEfficiencyHarvester::makeSummaryVsLumi(DQMStore::IGetter& getter) const {
for (unsigned int iLayer = 1; iLayer != (showRings_ ? 20 : 22); ++iLayer) {
auto hfound = getter.get(fmt::format("SiStrip/HitEfficiency/layerfound_vsLumi_layer_{}", iLayer))->getTH1F();
auto htotal = getter.get(fmt::format("SiStrip/HitEfficiency/layertotal_vsLumi_layer_{}", iLayer))->getTH1F();
auto hfound =
getter.get(fmt::format("AlCaReco/SiStripHitEfficiency/layerfound_vsLumi_layer_{}", iLayer))->getTH1F();
auto htotal =
getter.get(fmt::format("AlCaReco/SiStripHitEfficiency/layertotal_vsLumi_layer_{}", iLayer))->getTH1F();

if (hfound == nullptr or htotal == nullptr) {
if (hfound == nullptr)
edm::LogError("SiStripHitEfficiencyHarvester")
<< fmt::format("AlCaReco/SiStripHitEfficiency/layerfound_vsLumi_layer_{}", iLayer) << " was not found!";
if (htotal == nullptr)
edm::LogError("SiStripHitEfficiencyHarvester")
<< fmt::format("AlCaReco/SiStripHitEfficiency/layertotal_vsLumi_layer_{}", iLayer) << " was not found!";
// no input histograms -> continue in the loop
continue;
}

if (!hfound->GetSumw2())
hfound->Sumw2();
if (!htotal->GetSumw2())
Expand All @@ -349,7 +423,6 @@ void SiStripHitEfficiencyHarvester::makeSummaryVsLumi(DQMStore::IGetter& getter,
<< "Total hits for layer " << iLayer << " (vs lumi): " << htotal->GetEntries() << ", found "
<< hfound->GetEntries();
}

// continue
}

Expand Down Expand Up @@ -578,5 +651,19 @@ void SiStripHitEfficiencyHarvester::printAndWriteBadModules(const SiStripQuality
badModules.close();
}

void SiStripHitEfficiencyHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
edm::ParameterSetDescription desc;
desc.add<bool>("isAtPCL", false);
desc.add<bool>("doStoreOnDB", false);
desc.add<std::string>("Record", "SiStripBadStrip");
desc.add<double>("Threshold", 0.1);
desc.add<std::string>("Title", "Hit Efficiency");
desc.add<int>("nModsMin", 5);
desc.addUntracked<bool>("AutoIneffModTagging", false);
desc.addUntracked<double>("TkMapMin", 0.9);
desc.addUntracked<bool>("ShowRings", false);
descriptions.addWithDefaultLabel(desc);
}

#include "FWCore/Framework/interface/MakerMacros.h"
DEFINE_FWK_MODULE(SiStripHitEfficiencyHarvester);
Expand Up @@ -9,8 +9,8 @@

#include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
#include "CalibTracker/Records/interface/SiStripQualityRcd.h"
#include "CalibTracker/SiStripHitEfficiency/interface/TrajectoryAtInvalidHit.h"
#include "CalibTracker/SiStripHitEfficiency/interface/SiStripHitEfficiencyHelpers.h"
#include "CalibTracker/SiStripHitEfficiency/interface/TrajectoryAtInvalidHit.h"
#include "DQM/SiStripCommon/interface/TkHistoMap.h"
#include "DQMServices/Core/interface/DQMEDAnalyzer.h"
#include "DataFormats/Common/interface/DetSetVector.h"
Expand All @@ -35,7 +35,10 @@
#include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
#include "FWCore/ParameterSet/interface/ParameterDescription.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
#include "Geometry/CommonDetUnit/interface/GeomDet.h"
#include "Geometry/CommonDetUnit/interface/GeomDetType.h"
#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
Expand All @@ -58,6 +61,7 @@ class SiStripHitEfficiencyWorker : public DQMEDAnalyzer {
public:
explicit SiStripHitEfficiencyWorker(const edm::ParameterSet& conf);
~SiStripHitEfficiencyWorker() override = default;
static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);

private:
void beginJob(); // TODO remove
Expand Down Expand Up @@ -193,20 +197,19 @@ SiStripHitEfficiencyWorker::SiStripHitEfficiencyWorker(const edm::ParameterSet&
propagatorToken_(esConsumes(edm::ESInputTag{"", "PropagatorWithMaterial"})),
tkDetMapToken_(esConsumes<edm::Transition::BeginRun>()),
layers_(conf.getParameter<int>("Layer")),
DEBUG_(conf.getParameter<bool>("Debug")),
DEBUG_(conf.getUntrackedParameter<bool>("Debug", false)),
addLumi_(conf.getUntrackedParameter<bool>("addLumi", false)),
addCommonMode_(conf.getUntrackedParameter<bool>("addCommonMode", false)),
cutOnTracks_(conf.getUntrackedParameter<bool>("cutOnTracks", false)),
trackMultiplicityCut_(conf.getUntrackedParameter<unsigned int>("trackMultiplicity", 100)),
useFirstMeas_(conf.getUntrackedParameter<bool>("useFirstMeas", false)),
useLastMeas_(conf.getUntrackedParameter<bool>("useLastMeas", false)),
useAllHitsFromTracksWithMissingHits_(
conf.getUntrackedParameter<bool>("useAllHitsFromTracksWithMissingHits", false)),
clusterMatchingMethod_(conf.getUntrackedParameter<int>("ClusterMatchingMethod", 0)),
resXSig_(conf.getUntrackedParameter<double>("ResXSig", -1)),
clusterTracjDist_(conf.getUntrackedParameter<double>("ClusterTrajDist", 64.0)),
stripsApvEdge_(conf.getUntrackedParameter<double>("StripsApvEdge", 10.0)),
useOnlyHighPurityTracks_(conf.getUntrackedParameter<bool>("UseOnlyHighPurityTracks", true)),
cutOnTracks_(conf.getParameter<bool>("cutOnTracks")),
trackMultiplicityCut_(conf.getParameter<unsigned int>("trackMultiplicity")),
useFirstMeas_(conf.getParameter<bool>("useFirstMeas")),
useLastMeas_(conf.getParameter<bool>("useLastMeas")),
useAllHitsFromTracksWithMissingHits_(conf.getParameter<bool>("useAllHitsFromTracksWithMissingHits")),
clusterMatchingMethod_(conf.getParameter<int>("ClusterMatchingMethod")),
resXSig_(conf.getParameter<double>("ResXSig")),
clusterTracjDist_(conf.getParameter<double>("ClusterTrajDist")),
stripsApvEdge_(conf.getParameter<double>("StripsApvEdge")),
useOnlyHighPurityTracks_(conf.getParameter<bool>("UseOnlyHighPurityTracks")),
bunchX_(conf.getUntrackedParameter<int>("BunchCrossing", 0)),
showRings_(conf.getUntrackedParameter<bool>("ShowRings", false)),
showTOB6TEC9_(conf.getUntrackedParameter<bool>("ShowTOB6TEC9", false)) {
Expand Down Expand Up @@ -246,7 +249,7 @@ void SiStripHitEfficiencyWorker::beginJob() {
void SiStripHitEfficiencyWorker::bookHistograms(DQMStore::IBooker& booker,
const edm::Run& run,
const edm::EventSetup& setup) {
const std::string path = "SiStrip/HitEfficiency"; // TODO make this configurable
const std::string path = "AlCaReco/SiStripHitEfficiency"; // TODO make this configurable
booker.setCurrentFolder(path);
h_bx = booker.book1D("bx", "bx", 3600, 0, 3600);
h_instLumi = booker.book1D("instLumi", "inst. lumi.", 250, 0, 25000);
Expand Down Expand Up @@ -904,7 +907,35 @@ void SiStripHitEfficiencyWorker::endJob() {
LogDebug("SiStripHitEfficiencyWorker") << " Number Of Tracked events " << EventTrackCKF;
}

void SiStripHitEfficiencyWorker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
edm::ParameterSetDescription desc;
desc.add<bool>("UseOnlyHighPurityTracks", true);
desc.add<bool>("cutOnTracks", false);
desc.add<bool>("useAllHitsFromTracksWithMissingHits", false);
desc.add<bool>("useFirstMeas", false);
desc.add<bool>("useLastMeas", false);
desc.add<double>("ClusterTrajDist", 64.0);
desc.add<double>("ResXSig", -1);
desc.add<double>("StripsApvEdge", 10.0);
desc.add<edm::InputTag>("combinatorialTracks", edm::InputTag{"generalTracks"});
desc.add<edm::InputTag>("commonMode", edm::InputTag{"siStripDigis", "CommonMode"});
desc.add<edm::InputTag>("lumiScalers", edm::InputTag{"scalersRawToDigi"});
desc.add<edm::InputTag>("siStripClusters", edm::InputTag{"siStripClusters"});
desc.add<edm::InputTag>("siStripDigis", edm::InputTag{"siStripDigis"});
desc.add<edm::InputTag>("trackerEvent", edm::InputTag{"MeasurementTrackerEvent"});
desc.add<edm::InputTag>("trajectories", edm::InputTag{"generalTracks"});
desc.add<int>("ClusterMatchingMethod", 0);
desc.add<int>("Layer", 0);
desc.add<unsigned int>("trackMultiplicity", 100);
desc.addUntracked<bool>("Debug", false);
desc.addUntracked<bool>("ShowRings", false);
desc.addUntracked<bool>("ShowTOB6TEC9", false);
desc.addUntracked<bool>("addCommonMode", false);
desc.addUntracked<bool>("addLumi", true);
desc.addUntracked<int>("BunchCrossing", 0);
desc.addUntracked<std::string>("BadModulesFile", "");
descriptions.addWithDefaultLabel(desc);
}

#include "FWCore/Framework/interface/MakerMacros.h"
DEFINE_FWK_MODULE(SiStripHitEfficiencyWorker);

// TODO next: try to run this
Expand Up @@ -35,9 +35,7 @@
doStoreOnDB = cms.bool(True),
ShowRings = cms.untracked.bool(False), # default False
TkMapMin = cms.untracked.double(0.90), # default 0.90
EffPlotMin = cms.untracked.double(0.90), # default 0.90
Title = cms.string(' Hit Efficiency - run {0:d}'.format(runNumber))
)
Title = cms.string(' Hit Efficiency - run {0:d}'.format(runNumber)))

process.load("DQM.SiStripCommon.TkHistoMap_cff")

Expand Down
22 changes: 11 additions & 11 deletions CalibTracker/SiStripHitEfficiency/test/testHitEffWorker.py
Expand Up @@ -44,28 +44,28 @@
tracks = cms.InputTag("refitTracks")

process.hiteff = cms.EDProducer("SiStripHitEfficiencyWorker",
lumiScalers=cms.InputTag("scalersRawToDigi"),
lumiScalers =cms.InputTag("scalersRawToDigi"),
addLumi = cms.untracked.bool(True),
commonMode=cms.InputTag("siStripDigis", "CommonMode"),
addCommonMode=cms.untracked.bool(False),
commonMode = cms.InputTag("siStripDigis", "CommonMode"),
addCommonMode = cms.untracked.bool(False),
combinatorialTracks = tracks,
trajectories = tracks,
siStripClusters = cms.InputTag("siStripClusters"),
siStripDigis = cms.InputTag("siStripDigis"),
trackerEvent = cms.InputTag("MeasurementTrackerEvent"),
# part 2
Layer = cms.int32(0), # =0 means do all layers
Debug = cms.bool(True),
Debug = cms.untracked.bool(True),
# do not cut on the total number of tracks
cutOnTracks = cms.untracked.bool(True),
trackMultiplicity = cms.untracked.uint32(100),
cutOnTracks = cms.bool(True),
trackMultiplicity = cms.uint32(100),
# use or not first and last measurement of a trajectory (biases), default is false
useFirstMeas = cms.untracked.bool(False),
useLastMeas = cms.untracked.bool(False),
useAllHitsFromTracksWithMissingHits = cms.untracked.bool(False),
useFirstMeas = cms.bool(False),
useLastMeas = cms.bool(False),
useAllHitsFromTracksWithMissingHits = cms.bool(False),
## non-default settings
ClusterMatchingMethod = cms.untracked.int32(4), # default 0 case0,1,2,3,4
ClusterTrajDist = cms.untracked.double(15), # default 64
ClusterMatchingMethod = cms.int32(4), # default 0 case0,1,2,3,4
ClusterTrajDist = cms.double(15), # default 64
)

process.load("DQM.SiStripCommon.TkHistoMap_cff")
Expand Down
@@ -0,0 +1,9 @@
import FWCore.ParameterSet.Config as cms

OutALCARECOPromptCalibProdSiStripHitEfficiency_noDrop = cms.PSet(
SelectEvents = cms.untracked.PSet(
SelectEvents = cms.vstring('pathALCARECOPromptCalibProdSiStripHitEfficiency')),
outputCommands = cms.untracked.vstring('keep *_MEtoEDMConvertSiStripHitEff_*_*'))

OutALCARECOPromptCalibProdSiStripHitEfficiency = OutALCARECOPromptCalibProdSiStripHitEfficiency_noDrop.clone()
OutALCARECOPromptCalibProdSiStripHitEfficiency.outputCommands.insert(0, "drop *")