From 637b7dec02e1a1d9edc2414bce7c2b9d85c4324f Mon Sep 17 00:00:00 2001
From: "s.zharko@gsi.de" <s.zharko@gsi.de>
Date: Fri, 22 Sep 2023 03:00:10 +0200
Subject: [PATCH] CA QA: Automatic plot for tracking setup
---
core/qa/CbmQaTask.cxx | 17 +++++
core/qa/CbmQaTask.h | 13 +++-
macro/mcbm/mcbm_qa.C | 1 +
reco/L1/CbmL1.cxx | 1 +
reco/L1/CbmL1.h | 1 +
reco/L1/CbmL1Performance.cxx | 1 +
reco/L1/qa/CbmCaInputQaSetup.cxx | 119 +++++++++++++++++++++++++------
reco/L1/qa/CbmCaInputQaSetup.h | 18 ++++-
reco/L1/qa/CbmCaOutputQa.cxx | 10 ++-
9 files changed, 156 insertions(+), 25 deletions(-)
diff --git a/core/qa/CbmQaTask.cxx b/core/qa/CbmQaTask.cxx
index ecc17172a5..130638f4f2 100644
--- a/core/qa/CbmQaTask.cxx
+++ b/core/qa/CbmQaTask.cxx
@@ -17,6 +17,8 @@
#include "TAxis.h"
#include "TCanvas.h"
#include "TF1.h"
+#include "TPad.h"
+#include "TPaveText.h"
#include "TString.h"
#include <array>
@@ -137,3 +139,18 @@ bool CbmQaTask::CheckRange(TH1* h, double meanMax, double rmsMin, double rmsMax)
}
return res;
}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void CbmQaTask::PutSetupNameOnPad(double xMin, double yMin, double xMax, double yMax)
+{
+ if (gPad) {
+ gPad->cd();
+ auto* paveText = new TPaveText(xMin, yMin, xMax, yMax, "NDC");
+ paveText->AddText(fsSetupName.c_str());
+ paveText->SetBorderSize(0);
+ paveText->SetTextSize(0.04);
+ paveText->SetFillStyle(0);
+ paveText->Draw("same");
+ }
+}
diff --git a/core/qa/CbmQaTask.h b/core/qa/CbmQaTask.h
index cfd3ea35cc..f09f820141 100644
--- a/core/qa/CbmQaTask.h
+++ b/core/qa/CbmQaTask.h
@@ -58,6 +58,9 @@ public:
CbmQaTask& operator=(const CbmQaTask&) = delete;
CbmQaTask& operator=(CbmQaTask&&) = delete;
+ /// @brief Gets name of the setup
+ const std::string& GetSetupName() const { return fsSetupName; }
+
/// Returns flag, whether MC information is used or not in the task
bool IsMCUsed() const { return fbUseMC; }
@@ -73,6 +76,10 @@ public:
/// FairTask: Defines action of the task in the end of run
void Finish();
+ /// @brief Sets name of the setup
+ void SetSetupName(const char* setup) { fsSetupName = setup; }
+
+protected:
ClassDef(CbmQaTask, 0);
protected:
@@ -104,7 +111,6 @@ protected:
/// Get current event number
int GetEventNumber() const { return fNofEvents.GetVal(); }
-
// ***********************
// ** Utility functions **
// ***********************
@@ -132,12 +138,17 @@ protected:
/// \return False, if variable exits the range
bool CheckRange(TH1* h, double meanMax, double rmsMin, double rmsMax);
+ /// @brief Puts setup title on the canvas
+ void PutSetupNameOnPad(double xMin, double yMin, double xMax, double yMax);
+
private:
/// @brief De-initializes this task
void DeInitBase();
bool fbUseMC = false; ///< Flag, if MC is used
+ std::string fsSetupName = ""; ///< Name of the setup (to draw on the canvases)
+
TParameter<int> fNofEvents {"nEvents", 0}; ///< Number of processed events
};
diff --git a/macro/mcbm/mcbm_qa.C b/macro/mcbm/mcbm_qa.C
index 1a67d7523f..1e1274db87 100644
--- a/macro/mcbm/mcbm_qa.C
+++ b/macro/mcbm/mcbm_qa.C
@@ -251,6 +251,7 @@ void mcbm_qa(Int_t nEvents = 0, TString dataset = "data/mcbm_beam_2020_03_test",
pCaInputQaSetup->SetDetectorFlag(L1DetectorID::kTrd, bUseTrd);
pCaInputQaSetup->SetDetectorFlag(L1DetectorID::kTof, bUseTof);
pCaInputQaSetup->ReadParameters(caParFile.Data());
+ pCaInputQaSetup->SetSetupName(setupName.Data());
run->AddTask(pCaInputQaSetup);
auto* pCaOutputQa = new cbm::ca::OutputQa(verbose, bUseMC);
diff --git a/reco/L1/CbmL1.cxx b/reco/L1/CbmL1.cxx
index 61a2f307a9..5113d5e38d 100644
--- a/reco/L1/CbmL1.cxx
+++ b/reco/L1/CbmL1.cxx
@@ -663,6 +663,7 @@ InitStatus CbmL1::Init()
fMonitor.SetKeyName(EMonitorKey::kEvent, "N events");
fMonitor.SetKeyName(EMonitorKey::kRecoTrack, "N reco tracks");
fMonitor.SetKeyName(EMonitorKey::kRecoHit, "N hits");
+ fMonitor.SetKeyName(EMonitorKey::kGhostTrack, "N ghosts");
fMonitor.SetKeyName(EMonitorKey::kDiscardedHitMvd, "N discarded hits in MVD");
fMonitor.SetKeyName(EMonitorKey::kDiscardedHitSts, "N discarded hits in STS");
fMonitor.SetKeyName(EMonitorKey::kDiscardedHitMuch, "N discarded hits in MUCH");
diff --git a/reco/L1/CbmL1.h b/reco/L1/CbmL1.h
index 8ed28ff156..2da65db31d 100644
--- a/reco/L1/CbmL1.h
+++ b/reco/L1/CbmL1.h
@@ -264,6 +264,7 @@ public:
kEvent, ///< Number of processed events
kRecoTrack, ///< Number of reconstructed tracks
kRecoHit, ///< Number of reconstructed hits
+ kGhostTrack, ///< Number of ghost tracks
kDiscardedHitMvd, ///< Number of discarded hits in MVD
kDiscardedHitSts, ///< Number of discarded hits in STS
kDiscardedHitMuch, ///< Number of discarded hits in MuCh
diff --git a/reco/L1/CbmL1Performance.cxx b/reco/L1/CbmL1Performance.cxx
index 0e73bf2bdc..bd3cf8d243 100644
--- a/reco/L1/CbmL1Performance.cxx
+++ b/reco/L1/CbmL1Performance.cxx
@@ -875,6 +875,7 @@ void CbmL1::HistoPerformance() // TODO: check if works correctly. Change vHitFa
// fill ghost histos
if (prtra->IsGhost()) {
+ fMonitor.Increment(EMonitorKey::kGhostTrack);
h_ghost_purity->Fill(100 * prtra->GetMaxPurity());
if (fabs(prtra->T[4]) > 1.e-10) {
h_ghost_mom->Fill(fabs(1.0 / prtra->T[4]));
diff --git a/reco/L1/qa/CbmCaInputQaSetup.cxx b/reco/L1/qa/CbmCaInputQaSetup.cxx
index f2f8a1844a..2fb7f88ad9 100644
--- a/reco/L1/qa/CbmCaInputQaSetup.cxx
+++ b/reco/L1/qa/CbmCaInputQaSetup.cxx
@@ -13,6 +13,8 @@
#include "FairRootManager.h"
+#include "TAxis.h"
+
#include "L1InitManager.h"
using cbm::ca::InputQaSetup;
@@ -52,11 +54,98 @@ void InputQaSetup::FillHistograms()
if (fvbUseDet[L1DetectorID::kTof]) { this->FillHistogramsDet<L1DetectorID::kTof>(); }
}
+// ---------------------------------------------------------------------------------------------------------------------
+//
+InitStatus InputQaSetup::InitCanvases()
+{
+ // Tracking setup by hits
+ {
+ auto* canv = MakeQaObject<CbmQaCanvas>("c_setup_hits", "Setup by hits", 1000, 1000);
+ canv->Divide(1, 2, 0.000001, 0.000001);
+ canv->cd(1);
+ gPad->SetLogz();
+ gPad->SetBottomMargin(0.000001);
+ gPad->SetRightMargin(0.05);
+ gPad->SetTitle("");
+ auto* ph_hit_xz = dynamic_cast<TH2F*>(fph_hit_xz.back()->Clone("h_hit_xz_clone"));
+ ph_hit_xz->SetTitle(";z [cm];x [cm]");
+ ph_hit_xz->Draw("col");
+ this->PutSetupNameOnPad(0.08, 0.83, 0.50, 0.89);
+ canv->cd(2);
+ gPad->SetLogz();
+ gPad->SetTopMargin(0.000001);
+ gPad->SetRightMargin(0.05);
+ gPad->SetTitle("");
+ auto* ph_hit_yz = dynamic_cast<TH2F*>(fph_hit_yz.back()->Clone("h_hit_yz_clone"));
+ ph_hit_yz->SetTitle(";z [cm];y [cm]");
+ ph_hit_yz->Draw("col");
+ this->PutSetupNameOnPad(0.08, 0.93, 0.50, 0.99);
+
+
+ auto LoBinEdge = [](TAxis* pAxis, double val) { return pAxis->GetBinLowEdge(pAxis->FindBin(val)); };
+ auto UpBinEdge = [](TAxis* pAxis, double val) { return pAxis->GetBinUpEdge(pAxis->FindBin(val)); };
+
+ for (int iDet = 0; iDet < static_cast<int>(fvpDetInterface.size()); ++iDet) {
+ if (!fvbUseDet[iDet]) { continue; }
+ int nSt = fvpDetInterface[iDet]->GetNtrackingStations();
+ for (int iSt = 0; iSt < nSt; ++iSt) {
+ int iStActive = fpParameters->GetStationIndexActive(iSt, static_cast<L1DetectorID>(iDet));
+ Color_t boxColor = (iStActive == -1) ? kGray + 1 : kOrange + 2;
+ const char* detName = Form("%s-%d", fvpDetInterface[iDet]->GetDetectorName().c_str(), iSt);
+ {
+ double zMin = LoBinEdge(fph_hit_xz.back()->GetXaxis(), fvZmin[iDet][iSt]);
+ double zMax = UpBinEdge(fph_hit_xz.back()->GetXaxis(), fvZmax[iDet][iSt]);
+ double xMin = LoBinEdge(fph_hit_xz.back()->GetYaxis(), fvXmin[iDet][iSt]);
+ double xMax = UpBinEdge(fph_hit_xz.back()->GetYaxis(), fvXmax[iDet][iSt]);
+ canv->cd(1);
+ auto* pBox = new TBox(zMin, xMin, zMax, xMax);
+ pBox->SetFillStyle(0);
+ pBox->SetLineWidth(2);
+ pBox->SetLineColor(boxColor);
+ pBox->Draw("same");
+ auto* pText = new TText(zMin, xMax, detName);
+ pText->SetTextColor(boxColor);
+ pText->SetTextSize(0.035);
+ pText->SetTextAngle(45);
+ pText->Draw("same");
+ }
+ {
+ double zMin = LoBinEdge(fph_hit_yz.back()->GetXaxis(), fvZmin[iDet][iSt]);
+ double zMax = UpBinEdge(fph_hit_yz.back()->GetXaxis(), fvZmax[iDet][iSt]);
+ double yMin = LoBinEdge(fph_hit_yz.back()->GetYaxis(), fvYmin[iDet][iSt]);
+ double yMax = UpBinEdge(fph_hit_yz.back()->GetYaxis(), fvYmax[iDet][iSt]);
+ canv->cd(2);
+ auto* pBox = new TBox(zMin, yMin, zMax, yMax);
+ pBox->SetFillStyle(0);
+ pBox->SetLineWidth(2);
+ pBox->SetLineColor(boxColor);
+ pBox->Draw("same");
+ auto* pText = new TText(zMin, yMax, detName);
+ pText->SetTextSize(0.035);
+ pText->SetTextColor(boxColor);
+ pText->SetTextAngle(45);
+ pText->Draw("same");
+ }
+ }
+ }
+ }
+ return kSUCCESS;
+}
+
// ---------------------------------------------------------------------------------------------------------------------
//
InitStatus InputQaSetup::InitDataBranches()
try {
LOG(info) << fName << ": initializing... ";
+ // Tracking parameters
+ {
+ fpParameters = std::make_unique<L1Parameters>();
+
+ L1InitManager manager;
+ manager.ReadParametersObject(fsParametersFilename.c_str());
+ manager.SendParameters(*fpParameters);
+ LOG(info) << "CA tracking parameters object: " << fsParametersFilename << '\n' << fpParameters->ToString(1);
+ }
auto pFairManager = FairRootManager::Instance();
assert(pFairManager);
@@ -68,7 +157,16 @@ try {
fvpBrHits.fill(nullptr);
auto InitDetInterface = [&](CbmTrackingDetectorInterfaceBase* ifs, L1DetectorID detID) {
- if (fvbUseDet[detID]) { fvpDetInterface[detID] = ifs; }
+ if (fvbUseDet[detID]) {
+ fvpDetInterface[detID] = ifs;
+ int nSt = ifs->GetNtrackingStations();
+ fvXmin[detID].resize(nSt, std::numeric_limits<double>::max());
+ fvXmax[detID].resize(nSt, std::numeric_limits<double>::lowest());
+ fvYmin[detID].resize(nSt, std::numeric_limits<double>::max());
+ fvYmax[detID].resize(nSt, std::numeric_limits<double>::lowest());
+ fvZmin[detID].resize(nSt, std::numeric_limits<double>::max());
+ fvZmax[detID].resize(nSt, std::numeric_limits<double>::lowest());
+ }
};
auto InitHitBranch = [&](const char* brName, L1DetectorID detID) {
@@ -153,22 +251,3 @@ InitStatus InputQaSetup::InitHistograms()
return kSUCCESS;
}
-
-
-// ---------------------------------------------------------------------------------------------------------------------
-//
-void InputQaSetup::ReadParameters(const char* filename)
-{
- fpParameters = std::make_unique<L1Parameters>();
-
- L1InitManager manager;
- manager.ReadParametersObject(filename);
- manager.SendParameters(*fpParameters);
- LOG(info) << "CA tracking parameters object: " << filename << '\n' << fpParameters->ToString(1);
-}
-
-
-// ---------------------------------------------------------------------------------------------------------------------
-//
-// ---------------------------------------------------------------------------------------------------------------------
-//
diff --git a/reco/L1/qa/CbmCaInputQaSetup.h b/reco/L1/qa/CbmCaInputQaSetup.h
index ee06e67daf..31dc283356 100644
--- a/reco/L1/qa/CbmCaInputQaSetup.h
+++ b/reco/L1/qa/CbmCaInputQaSetup.h
@@ -50,17 +50,16 @@ namespace cbm::ca
/// @brief Reads defined parameters object from file
/// @param filename Name of parameter file
/// @note TEMPORARY FUNCTION, A SEPARATE PARAMETERS INITIALIZATION CLASS IS TO BE USED
- void ReadParameters(const char* filename);
+ void ReadParameters(const char* filename) { fsParametersFilename = filename; }
/// @brief Sets detector flag
void SetDetectorFlag(L1DetectorID detID, bool flag = true) { fvbUseDet[detID] = flag; }
- protected:
/// @brief Checks results of the QA and returns a success flag
bool Check();
/// @brief Initializes canvases
- InitStatus InitCanvases() { return kSUCCESS; };
+ InitStatus InitCanvases();
/// @brief Initializes data branches
InitStatus InitDataBranches();
@@ -96,7 +95,14 @@ namespace cbm::ca
CbmMCEventList* fpMCEventList = nullptr; ///< Pointer to MC event list
CbmMCDataObject* fpMCEventHeader = nullptr; ///< Pointer to MC event header
std::unique_ptr<L1Parameters> fpParameters = nullptr; ///< Pointer to CA parameters object
+ std::string fsParametersFilename = ""; ///< Filename for the tracking parameters
+ DetIdArr_t<std::vector<double>> fvXmin;
+ DetIdArr_t<std::vector<double>> fvXmax;
+ DetIdArr_t<std::vector<double>> fvYmin;
+ DetIdArr_t<std::vector<double>> fvYmax;
+ DetIdArr_t<std::vector<double>> fvZmin;
+ DetIdArr_t<std::vector<double>> fvZmax;
// ******************
// ** Histograms **
@@ -134,6 +140,12 @@ void cbm::ca::InputQaSetup::FillHistogramsDet()
auto xHit = pHit->GetX();
auto yHit = pHit->GetY();
auto zHit = pHit->GetZ();
+ if (fvXmin[DetID][iStLoc] > xHit) { fvXmin[DetID][iStLoc] = xHit; }
+ if (fvXmax[DetID][iStLoc] < xHit) { fvXmax[DetID][iStLoc] = xHit; }
+ if (fvYmin[DetID][iStLoc] > yHit) { fvYmin[DetID][iStLoc] = yHit; }
+ if (fvYmax[DetID][iStLoc] < yHit) { fvYmax[DetID][iStLoc] = yHit; }
+ if (fvZmin[DetID][iStLoc] > zHit) { fvZmin[DetID][iStLoc] = zHit; }
+ if (fvZmax[DetID][iStLoc] < zHit) { fvZmax[DetID][iStLoc] = zHit; }
fph_hit_xz[iStGeo]->Fill(zHit, xHit);
fph_hit_yz[iStGeo]->Fill(zHit, yHit);
diff --git a/reco/L1/qa/CbmCaOutputQa.cxx b/reco/L1/qa/CbmCaOutputQa.cxx
index 455fc60249..8c114a4fba 100644
--- a/reco/L1/qa/CbmCaOutputQa.cxx
+++ b/reco/L1/qa/CbmCaOutputQa.cxx
@@ -97,6 +97,8 @@ void OutputQa::EnableDebugger(const char* filename)
//
void OutputQa::FillHistograms()
{
+ fMonitor.Increment(EMonitorKey::kEvent);
+
// ************************************************************************
// ** Fill distributions for reconstructed tracks (including ghost ones) **
// ************************************************************************
@@ -611,11 +613,12 @@ bool OutputQa::Check()
aTable->SetRowName(iRow, "N ghosts");
aTable->SetCell(iRow, 0, nGhosts);
aTable->SetRowName(iRow + 1, "Ghost rate");
- aTable->SetCell(iRow + 1, 0, nGhosts / fvpTrackHistograms[ETrackType::kAll]->GetNofMCTracks());
+ aTable->SetCell(iRow + 1, 0, nGhosts / fMonitor.GetValue(EMonitorKey::kTrack));
}
LOG(info) << '\n' << aTable->ToString(3);
}
+ LOG(info) << '\n' << fMonitor.ToString();
return true;
}
@@ -635,6 +638,9 @@ InitStatus OutputQa::InitTimeSlice()
nHits = fvHits.size();
nRecoTracks = fvRecoTracks.size();
+ fMonitor.Increment(EMonitorKey::kTrack, nRecoTracks);
+ fMonitor.Increment(EMonitorKey::kHit, nHits);
+
// Match tracks and points
// Read MC tracks and points
if (IsMCUsed()) {
@@ -642,6 +648,8 @@ InitStatus OutputQa::InitTimeSlice()
nMCPoints = fMCData.GetNofPoints();
nMCTracks = fMCData.GetNofTracks();
fpMCModule->MatchRecoAndMC();
+ fMonitor.Increment(EMonitorKey::kMcPoint, nMCPoints);
+ fMonitor.Increment(EMonitorKey::kMcTrack, nMCTracks);
}
// ****** DEBUG: BEGIN
--
GitLab