Skip to content

Commit

Permalink
Merge pull request #43893 from mmusich/mm_dev_misc_improv_forPhase2
Browse files Browse the repository at this point in the history
[14.0.X]  miscellaneous updates to improve validation of phase-2 Tracker Alignment samples
  • Loading branch information
cmsbuild authored Feb 9, 2024
2 parents df18e5b + 0cea335 commit 0e0543d
Show file tree
Hide file tree
Showing 7 changed files with 272 additions and 96 deletions.
21 changes: 13 additions & 8 deletions Alignment/OfflineValidation/macros/loopAndPlot.C
Original file line number Diff line number Diff line change
Expand Up @@ -32,12 +32,12 @@ void MakeNiceProfile(TProfile *prof);
//void MakeNicePlotStyle(TH1 *hist);
void plot2Histograms(TH1 *h1, TH1 *h2, const TString &label1, const TString &label2);
void plot2Profiles(TProfile *h1, TProfile *h2, const TString &label1, const TString &label2);
void recurseOverKeys(TDirectory *target1, const std::vector<TString> &labels);
void recurseOverKeys(TDirectory *target1, const std::vector<TString> &labels, bool isNorm);
void recurseOverKeys(TDirectory *target1, const TString &label1, const TString &label2);
void plotHistograms(std::vector<TH1 *> histos, const std::vector<TString> &labels);
void plotHistograms(std::vector<TH1 *> histos, const std::vector<TString> &labels, bool isNormalized = false);

/************************************************/
void loopAndPlot(TString namesandlabels)
void loopAndPlot(TString namesandlabels, bool doNormalize = false)
/************************************************/
{
std::vector<TString> labels;
Expand All @@ -56,15 +56,15 @@ void loopAndPlot(TString namesandlabels)
}
}

recurseOverKeys(sourceFiles[0], labels);
recurseOverKeys(sourceFiles[0], labels, doNormalize);

for (const auto &file : sourceFiles) {
file->Close();
}
}

/************************************************/
void recurseOverKeys(TDirectory *target1, const std::vector<TString> &labels)
void recurseOverKeys(TDirectory *target1, const std::vector<TString> &labels, bool isNorm)
/************************************************/
{
// Figure out where we are
Expand Down Expand Up @@ -107,7 +107,7 @@ void recurseOverKeys(TDirectory *target1, const std::vector<TString> &labels)

//outputFilename=histName;
//plot2Histograms(htemp1, htemp2, outputFolder+path+"/"+outputFilename+"."+imageType);
plotHistograms(histos, labels);
plotHistograms(histos, labels, isNorm);

} else if (obj->IsA()->InheritsFrom("TDirectory")) {
// it's a subdirectory
Expand All @@ -122,14 +122,14 @@ void recurseOverKeys(TDirectory *target1, const std::vector<TString> &labels)
if ((TString(obj->GetName())).Contains("Residuals"))
continue;

recurseOverKeys((TDirectory *)obj, labels);
recurseOverKeys((TDirectory *)obj, labels, isNorm);

} // end of IF a TDriectory
}
}

/************************************************/
void plotHistograms(std::vector<TH1 *> histos, const std::vector<TString> &labels) {
void plotHistograms(std::vector<TH1 *> histos, const std::vector<TString> &labels, bool isNormalized) {
/************************************************/

TGaxis::SetMaxDigits(3);
Expand All @@ -143,6 +143,11 @@ void plotHistograms(std::vector<TH1 *> histos, const std::vector<TString> &label
for (const auto &histo : histos) {
MakeNicePlotStyle<TH1>(histo);

if (isNormalized) {
Double_t scale = 1. / histo->Integral();
histo->Scale(scale);
}

histo->SetLineColor(def_colors[index]);
histo->SetMarkerColor(def_colors[index]);
histo->SetMarkerStyle(20);
Expand Down
73 changes: 42 additions & 31 deletions Alignment/OfflineValidation/plugins/DMRChecker.cc
Original file line number Diff line number Diff line change
Expand Up @@ -154,10 +154,16 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
runInfoToken_(esConsumes<edm::Transition::BeginRun>()),
magFieldToken_(esConsumes<edm::Transition::BeginRun>()),
topoToken_(esConsumes<edm::Transition::BeginRun>()),
latencyToken_(esConsumes<edm::Transition::BeginRun>()),
isCosmics_(pset.getParameter<bool>("isCosmics")) {
isCosmics_(pset.getParameter<bool>("isCosmics")),
doLatencyAnalysis_(pset.getParameter<bool>("doLatencyAnalysis")) {
usesResource(TFileService::kSharedResource);

if (doLatencyAnalysis_) {
latencyToken_ = esConsumes<edm::Transition::BeginRun>();
} else {
latencyToken_ = edm::ESGetToken<SiStripLatency, SiStripLatencyRcd>();
}

TkTag_ = pset.getParameter<edm::InputTag>("TkTag");
theTrackCollectionToken_ = consumes<reco::TrackCollection>(TkTag_);

Expand Down Expand Up @@ -231,7 +237,8 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
const edm::ESGetToken<RunInfo, RunInfoRcd> runInfoToken_;
const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magFieldToken_;
const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoToken_;
const edm::ESGetToken<SiStripLatency, SiStripLatencyRcd> latencyToken_;
// not const
edm::ESGetToken<SiStripLatency, SiStripLatencyRcd> latencyToken_;

const MagneticField *magneticField_;
const TrackerGeometry *trackerGeometry_;
Expand Down Expand Up @@ -463,12 +470,12 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
int ievt;
int itrks;
int mode;
bool firstEvent_;

SiPixelPI::phase phase_;
float etaMax_;

const bool isCosmics_;
const bool doLatencyAnalysis_;

edm::InputTag TkTag_;
edm::InputTag TriggerResultsTag_;
Expand Down Expand Up @@ -506,19 +513,6 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh

edm::Handle<reco::TrackCollection> trackCollection = event.getHandle(theTrackCollectionToken_);

if (firstEvent_) {
if (trackerGeometry_->isThere(GeomDetEnumerators::P2PXB) ||
trackerGeometry_->isThere(GeomDetEnumerators::P2PXEC)) {
phase_ = SiPixelPI::phase::two;
} else if (trackerGeometry_->isThere(GeomDetEnumerators::P1PXB) ||
trackerGeometry_->isThere(GeomDetEnumerators::P1PXEC)) {
phase_ = SiPixelPI::phase::one;
} else {
phase_ = SiPixelPI::phase::zero;
}
firstEvent_ = false;
}

GlobalPoint zeroPoint(0, 0, 0);
if (DEBUG)
edm::LogVerbatim("DMRChecker") << "event #" << ievt << " Event ID = " << event.id()
Expand Down Expand Up @@ -589,7 +583,8 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
continue; // is it a single physical module?

if ((*iHit)->isValid() && (subid > PixelSubdetector::PixelEndcap)) {
tmap->fill(detid_db, 1);
if (phase_ != SiPixelPI::phase::two)
tmap->fill(detid_db, 1);

//LocalPoint lp = (*iHit)->localPosition();
//LocalError le = (*iHit)->localPositionError();
Expand Down Expand Up @@ -1105,20 +1100,38 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
edm::LogVerbatim("DMRChecker") << "time difference: " << seconds << " s" << std::endl;
timeMap_[run.run()] = seconds;

//SiStrip Latency
const SiStripLatency *apvlat = &setup.getData(latencyToken_);
if (apvlat->singleReadOutMode() == 1) {
mode = 1; // peak mode
} else if (apvlat->singleReadOutMode() == 0) {
mode = -1; // deco mode
// set geometry and topology
trackerGeometry_ = &setup.getData(geomToken_);
if (trackerGeometry_->isThere(GeomDetEnumerators::P2PXB) || trackerGeometry_->isThere(GeomDetEnumerators::P2PXEC)) {
phase_ = SiPixelPI::phase::two;
} else if (trackerGeometry_->isThere(GeomDetEnumerators::P1PXB) ||
trackerGeometry_->isThere(GeomDetEnumerators::P1PXEC)) {
phase_ = SiPixelPI::phase::one;
} else {
phase_ = SiPixelPI::phase::zero;
}

trackerTopology_ = &setup.getData(topoToken_);

// if it's a phase-2 geometry there are no phase-1 conditions
if (phase_ == SiPixelPI::phase::two) {
mode = 0;
} else {
if (doLatencyAnalysis_) {
//SiStrip Latency
const SiStripLatency *apvlat = &setup.getData(latencyToken_);
if (apvlat->singleReadOutMode() == 1) {
mode = 1; // peak mode
} else if (apvlat->singleReadOutMode() == 0) {
mode = -1; // deco mode
}
} else {
mode = 0.;
}
}

conditionsMap_[run.run()].first = mode;
conditionsMap_[run.run()].second = B_;

// set geometry and topology
trackerGeometry_ = &setup.getData(geomToken_);
trackerTopology_ = &setup.getData(topoToken_);
}

//*************************************************************
Expand Down Expand Up @@ -1378,11 +1391,8 @@ class DMRChecker : public edm::one::EDAnalyzer<edm::one::WatchRuns, edm::one::Sh
vTrack2DHistos_.push_back(book<TH2F>("h2_kappa_vs_phi", "#kappa vs. #phi;#phi_{Track};#kappa", 100, -M_PI, M_PI, 100, .0, .05));
vTrack2DHistos_.push_back(book<TH2F>("h2_kappa_vs_eta", "#kappa vs. #eta;#eta_{Track};#kappa", 100, -etaMax_, etaMax_, 100, .0, .05));
vTrack2DHistos_.push_back(book<TH2F>("h2_normchi2_vs_kappa", "#kappa vs. #chi^{2}/ndof;#chi^{2}/ndof;#kappa", 100, 0., 10, 100, -.03, .03));

// clang-format on

firstEvent_ = true;

// create the full maps
fullPixelmapXDMR->createTrackerBaseMap();
fullPixelmapYDMR->createTrackerBaseMap();
Expand Down Expand Up @@ -2102,6 +2112,7 @@ void DMRChecker::fillDescriptions(edm::ConfigurationDescriptions &descriptions)
desc.add<edm::InputTag>("BeamSpotTag", edm::InputTag("offlineBeamSpot"));
desc.add<edm::InputTag>("VerticesTag", edm::InputTag("offlinePrimaryVertices"));
desc.add<bool>("isCosmics", false);
desc.add<bool>("doLatencyAnalysis", true);
descriptions.addWithDefaultLabel(desc);
}

Expand Down
86 changes: 52 additions & 34 deletions Alignment/OfflineValidation/plugins/DiMuonVertexValidation.cc
Original file line number Diff line number Diff line change
Expand Up @@ -79,18 +79,17 @@ class DiMuonVertexValidation : public edm::one::EDAnalyzer<edm::one::SharedResou
const bool useReco_;
const bool useClosestVertex_;
std::vector<double> pTthresholds_;
float maxSVdist_;
const float maxSVdist_;

// plot configurations

edm::ParameterSet CosPhiConfiguration_;
edm::ParameterSet CosPhi3DConfiguration_;
edm::ParameterSet VtxProbConfiguration_;
edm::ParameterSet VtxDistConfiguration_;
edm::ParameterSet VtxDist3DConfiguration_;
edm::ParameterSet VtxDistSigConfiguration_;
edm::ParameterSet VtxDist3DSigConfiguration_;
edm::ParameterSet DiMuMassConfiguration_;
const edm::ParameterSet CosPhiConfiguration_;
const edm::ParameterSet CosPhi3DConfiguration_;
const edm::ParameterSet VtxProbConfiguration_;
const edm::ParameterSet VtxDistConfiguration_;
const edm::ParameterSet VtxDist3DConfiguration_;
const edm::ParameterSet VtxDistSigConfiguration_;
const edm::ParameterSet VtxDist3DSigConfiguration_;
const edm::ParameterSet DiMuMassConfiguration_;

// control plots

Expand Down Expand Up @@ -127,13 +126,14 @@ class DiMuonVertexValidation : public edm::one::EDAnalyzer<edm::one::SharedResou

const edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> ttbESToken_;

edm::EDGetTokenT<reco::TrackCollection> tracksToken_; //used to select what tracks to read from configuration file
edm::EDGetTokenT<reco::VertexCollection> vertexToken_; //used to select what vertices to read from configuration file
//used to select what tracks to read from configuration file
edm::EDGetTokenT<reco::TrackCollection> tracksToken_;
//used to select what vertices to read from configuration file
const edm::EDGetTokenT<reco::VertexCollection> vertexToken_;

// either on or the other!
edm::EDGetTokenT<reco::MuonCollection> muonsToken_; //used to select what tracks to read from configuration file
edm::EDGetTokenT<reco::TrackCollection>
alcaRecoToken_; //used to select what muon tracks to read from configuration file
edm::EDGetTokenT<reco::MuonCollection> muonsToken_; // used to select tracks to read from configuration file
edm::EDGetTokenT<reco::TrackCollection> alcaRecoToken_; //used to select muon tracks to read from configuration file
};

//
Expand Down Expand Up @@ -350,7 +350,7 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
// fill the VtxProb plots
VtxProbPlots.fillPlots(SVProb, tktk_p4);

math::XYZPoint MainVertex(0, 0, 0);
math::XYZPoint mainVtxPos(0, 0, 0);
const reco::Vertex* theClosestVertex = nullptr;
// get collection of reconstructed vertices from event
edm::Handle<reco::VertexCollection> vertexHandle = iEvent.getHandle(vertexToken_);
Expand All @@ -362,35 +362,35 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
return;
}

reco::Vertex TheMainVertex;
reco::Vertex theMainVertex;
if (!useClosestVertex_ || theClosestVertex == nullptr) {
// if the closest vertex is not available, or explicitly not chosen
TheMainVertex = vertexHandle.product()->front();
theMainVertex = vertexHandle.product()->front();
} else {
TheMainVertex = *theClosestVertex;
theMainVertex = *theClosestVertex;
}

MainVertex.SetXYZ(TheMainVertex.position().x(), TheMainVertex.position().y(), TheMainVertex.position().z());
mainVtxPos.SetXYZ(theMainVertex.position().x(), theMainVertex.position().y(), theMainVertex.position().z());
const math::XYZPoint myVertex(
aTransientVertex.position().x(), aTransientVertex.position().y(), aTransientVertex.position().z());
const math::XYZPoint deltaVtx(
MainVertex.x() - myVertex.x(), MainVertex.y() - myVertex.y(), MainVertex.z() - myVertex.z());
mainVtxPos.x() - myVertex.x(), mainVtxPos.y() - myVertex.y(), mainVtxPos.z() - myVertex.z());

#ifdef EDM_ML_DEBUG
edm::LogVerbatim("DiMuonVertexValidation")
<< "mm vertex position:" << aTransientVertex.position().x() << "," << aTransientVertex.position().y() << ","
<< aTransientVertex.position().z();

edm::LogVerbatim("DiMuonVertexValidation") << "main vertex position:" << TheMainVertex.position().x() << ","
<< TheMainVertex.position().y() << "," << TheMainVertex.position().z();
edm::LogVerbatim("DiMuonVertexValidation") << "main vertex position:" << theMainVertex.position().x() << ","
<< theMainVertex.position().y() << "," << theMainVertex.position().z();
#endif

if (TheMainVertex.isValid()) {
if (theMainVertex.isValid()) {
// Z Vertex distance in the xy plane

VertexDistanceXY vertTool;
double distance = vertTool.distance(aTransientVertex, TheMainVertex).value();
double dist_err = vertTool.distance(aTransientVertex, TheMainVertex).error();
double distance = vertTool.distance(aTransientVertex, theMainVertex).value();
double dist_err = vertTool.distance(aTransientVertex, theMainVertex).error();

hSVDist_->Fill(distance * cmToum);
hSVDistSig_->Fill(distance / dist_err);
Expand All @@ -404,8 +404,8 @@ void DiMuonVertexValidation::analyze(const edm::Event& iEvent, const edm::EventS
// Z Vertex distance in 3D

VertexDistance3D vertTool3D;
double distance3D = vertTool3D.distance(aTransientVertex, TheMainVertex).value();
double dist3D_err = vertTool3D.distance(aTransientVertex, TheMainVertex).error();
double distance3D = vertTool3D.distance(aTransientVertex, theMainVertex).value();
double dist3D_err = vertTool3D.distance(aTransientVertex, theMainVertex).error();

hSVDist3D_->Fill(distance3D * cmToum);
hSVDist3DSig_->Fill(distance3D / dist3D_err);
Expand Down Expand Up @@ -459,14 +459,32 @@ void DiMuonVertexValidation::beginJob() {
TH1F::SetDefaultSumw2(kTRUE);
hSVProb_ = fs->make<TH1F>("VtxProb", ";#mu^{+}#mu^{-} vertex probability;N(#mu#mu pairs)", 100, 0., 1.);

hSVDist_ = fs->make<TH1F>("VtxDist", ";PV-#mu^{+}#mu^{-} vertex xy distance [#mum];N(#mu#mu pairs)", 100, 0., 300.);
hSVDistSig_ = fs->make<TH1F>("VtxDistSig", ";PV-#mu^{+}#mu^{-} vertex xy distance signficance;N(#mu#mu pairs)", 100, 0., 5.);
auto extractRangeValues = [](const edm::ParameterSet& PSetConfiguration_) -> std::pair<double, double> {
double min = PSetConfiguration_.getParameter<double>("ymin");
double max = PSetConfiguration_.getParameter<double>("ymax");
return {min, max};
};

// take the range from the 2D histograms
const auto& svDistRng = extractRangeValues(VtxDistConfiguration_);
hSVDist_ = fs->make<TH1F>("VtxDist", ";PV-#mu^{+}#mu^{-} vertex xy distance [#mum];N(#mu#mu pairs)", 100, svDistRng.first, svDistRng.second);

// take the range from the 2D histograms
const auto& svDistSigRng = extractRangeValues(VtxDistSigConfiguration_);
hSVDistSig_ = fs->make<TH1F>("VtxDistSig", ";PV-#mu^{+}#mu^{-} vertex xy distance signficance;N(#mu#mu pairs)", 100, svDistSigRng.first, svDistRng.second);

// take the range from the 2D histograms
const auto& svDist3DRng = extractRangeValues(VtxDist3DConfiguration_);
hSVDist3D_ = fs->make<TH1F>("VtxDist3D", ";PV-#mu^{+}#mu^{-} vertex 3D distance [#mum];N(#mu#mu pairs)", 100, svDist3DRng.first, svDist3DRng.second);

hSVDist3D_ = fs->make<TH1F>("VtxDist3D", ";PV-#mu^{+}#mu^{-} vertex 3D distance [#mum];N(#mu#mu pairs)", 100, 0., 300.);
hSVDist3DSig_ = fs->make<TH1F>("VtxDist3DSig", ";PV-#mu^{+}#mu^{-} vertex 3D distance signficance;N(#mu#mu pairs)", 100, 0., 5.);
// take the range from the 2D histograms
const auto& svDist3DSigRng = extractRangeValues(VtxDist3DSigConfiguration_);
hSVDist3DSig_ = fs->make<TH1F>("VtxDist3DSig", ";PV-#mu^{+}#mu^{-} vertex 3D distance signficance;N(#mu#mu pairs)", 100, svDist3DSigRng.first, svDist3DSigRng.second);

hInvMass_ = fs->make<TH1F>("InvMass", ";M(#mu#mu) [GeV];N(#mu#mu pairs)", 70., 50., 120.);
hTrackInvMass_ = fs->make<TH1F>("TkTkInvMass", ";M(tk,tk) [GeV];N(tk tk pairs)", 70., 50., 120.);
// take the range from the 2D histograms
const auto& massRng = extractRangeValues(DiMuMassConfiguration_);
hInvMass_ = fs->make<TH1F>("InvMass", ";M(#mu#mu) [GeV];N(#mu#mu pairs)", 70., massRng.first, massRng.second);
hTrackInvMass_ = fs->make<TH1F>("TkTkInvMass", ";M(tk,tk) [GeV];N(tk tk pairs)", 70., massRng.first, massRng.second);

hCosPhi_ = fs->make<TH1F>("CosPhi", ";cos(#phi_{xy});N(#mu#mu pairs)", 50, -1., 1.);
hCosPhi3D_ = fs->make<TH1F>("CosPhi3D", ";cos(#phi_{3D});N(#mu#mu pairs)", 50, -1., 1.);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -869,7 +869,7 @@ class GeneralPurposeTrackAnalyzer : public edm::one::EDAnalyzer<edm::one::WatchR
hEta = book<TH1D>("h_Eta", "Track pseudorapidity; track #eta;tracks", 100, -etaMax_, etaMax_);
hPhi = book<TH1D>("h_Phi", "Track azimuth; track #phi;tracks", 100, -M_PI, M_PI);

hPhiBarrel = book<TH1D>("h_PhiBarrel", "hPhiBarrel (0<|#eta|<0.8);track #Phi;tracks", 100, -M_PI, M_PI);
hPhiBarrel = book<TH1D>("h_PhiBarrel", "hPhiBarrel (0<|#eta|<0.8);track #phi;tracks", 100, -M_PI, M_PI);
hPhiOverlapPlus =
book<TH1D>("h_PhiOverlapPlus", "hPhiOverlapPlus (0.8<#eta<1.4);track #phi;tracks", 100, -M_PI, M_PI);
hPhiOverlapMinus =
Expand All @@ -880,7 +880,7 @@ class GeneralPurposeTrackAnalyzer : public edm::one::EDAnalyzer<edm::one::WatchR
h_BSx0 = book<TH1F>("h_BSx0", "x-coordinate of reco beamspot;x^{BS}_{0};n_{events}", 100, -0.1, 0.1);
h_BSy0 = book<TH1F>("h_BSy0", "y-coordinate of reco beamspot;y^{BS}_{0};n_{events}", 100, -0.1, 0.1);
h_BSz0 = book<TH1F>("h_BSz0", "z-coordinate of reco beamspot;z^{BS}_{0};n_{events}", 100, -1., 1.);
h_Beamsigmaz = book<TH1F>("h_Beamsigmaz", "z-coordinate beam width;#sigma_{Z}^{beam};n_{events}", 100, 0., 1.);
h_Beamsigmaz = book<TH1F>("h_Beamsigmaz", "z-coordinate beam width;#sigma_{Z}^{beam};n_{events}", 100, 0., 7.);
h_BeamWidthX = book<TH1F>("h_BeamWidthX", "x-coordinate beam width;#sigma_{X}^{beam};n_{events}", 100, 0., 0.01);
h_BeamWidthY = book<TH1F>("h_BeamWidthY", "y-coordinate beam width;#sigma_{Y}^{beam};n_{events}", 100, 0., 0.01);
h_BSdxdz = book<TH1F>("h_BSdxdz", "BeamSpot dxdz;beamspot dx/dz;n_{events}", 100, -0.0003, 0.0003);
Expand Down Expand Up @@ -1126,6 +1126,12 @@ class GeneralPurposeTrackAnalyzer : public edm::one::EDAnalyzer<edm::one::WatchR
fieldByRun_->GetXaxis()->SetBinLabel((the_r - theRuns_.front()) + 1, std::to_string(the_r).c_str());
}

static const int kappadiffindex = this->index(vTrackHistos_, "h_diff_curvature");
vTrackHistos_[kappadiffindex]->Add(vTrackHistos_[this->index(vTrackHistos_, "h_curvature_neg")],
vTrackHistos_[this->index(vTrackHistos_, "h_curvature_pos")],
-1,
1);

if (phase_ < SiPixelPI::phase::two) {
if (phase_ == SiPixelPI::phase::zero) {
pmap->save(true, 0, 0, "PixelHitMap.pdf", 600, 800);
Expand Down
Loading

0 comments on commit 0e0543d

Please sign in to comment.