From c8a96c09e7d526d7e4ec196889a979d48713773d Mon Sep 17 00:00:00 2001
From: "s.zharko@gsi.de" <s.zharko@gsi.de>
Date: Thu, 2 Mar 2023 03:01:31 +0100
Subject: [PATCH] CA: output qa, updated MC initialization
---
reco/L1/CbmCaMCModule.cxx | 109 +++++++++++++----------------
reco/L1/CbmCaMCModule.h | 115 ++++++++++++++-----------------
reco/L1/CbmCaTimeSliceReader.cxx | 11 ---
reco/L1/CbmCaTimeSliceReader.h | 8 +--
reco/L1/qa/CbmCaOutputQa.cxx | 49 +++++++++++--
reco/L1/qa/CbmCaOutputQa.h | 10 ++-
6 files changed, 151 insertions(+), 151 deletions(-)
diff --git a/reco/L1/CbmCaMCModule.cxx b/reco/L1/CbmCaMCModule.cxx
index c81d4e7b11..126cc20bda 100644
--- a/reco/L1/CbmCaMCModule.cxx
+++ b/reco/L1/CbmCaMCModule.cxx
@@ -122,23 +122,6 @@ try {
LOG(info) << "\tTRD: " << (fbUseTrd ? "ON" : "OFF");
LOG(info) << "\tTOF: " << (fbUseTof ? "ON" : "OFF");
- // Init QA module
- fpQaModule = std::make_unique<ca::tools::Qa>();
- fpQaModule->SetParameters(fpParameters);
- fpQaModule->SetMCData(&fMCData);
- fpQaModule->SetRecoTrackContainer(fpvRecoTracks);
- fpQaModule->SetHitContainer(fpvHits);
-
- // Setup output names
- // Names are defined using the name for main reconstruction output file, defined from the FairRunAna class
- boost::filesystem::path sOutPath = FairRunAna::Instance()->GetUserOutputFileName().Data();
- std::string sOutDir = sOutPath.parent_path().string();
- if (sOutDir.empty()) { sOutDir = "."; }
- std::string sOutHistoQa = sOutDir + "/CaHistoQa." + sOutPath.filename().string();
- fpQaModule->SetOutputFilename(sOutHistoQa.c_str());
-
- fpQaModule->InitHistograms();
-
return true;
}
catch (const std::logic_error& error) {
@@ -170,20 +153,22 @@ void CbmCaMCModule::InitEvent(CbmEvent* /*pEvent*/)
}
// ------ Read data
- fMCData.Clear();
+ fpMCData->Clear();
this->ReadMCTracks();
this->ReadMCPoints();
+ L1_SHOW(fpMCData->GetNofTracks());
+ L1_SHOW(fpMCData->GetNofPoints());
// ------ Prepare tracks: set point indexes and redefine indexes from external to internal containers
- for (auto& aTrk : fMCData.GetTrackContainer()) {
+ for (auto& aTrk : fpMCData->GetTrackContainer()) {
aTrk.SortPointIndexes(
- [&](const int& iPl, const int& iPr) { return fMCData.GetPoint(iPl).GetZ() < fMCData.GetPoint(iPr).GetZ(); });
+ [&](const int& iPl, const int& iPr) { return fpMCData->GetPoint(iPl).GetZ() < fpMCData->GetPoint(iPr).GetZ(); });
}
}
// ---------------------------------------------------------------------------------------------------------------------
//
-void CbmCaMCModule::ProcessEvent(CbmEvent*) { fpQaModule->FillHistograms(); }
+void CbmCaMCModule::ProcessEvent(CbmEvent*) {}
// ---------------------------------------------------------------------------------------------------------------------
//
@@ -191,11 +176,11 @@ void CbmCaMCModule::InitTrackInfo(const L1Vector<CbmL1HitDebugInfo>& vHits)
{
LOG(info) << "\033[1;32m!!!! FLAG 1\033[0m";
// ----- Initialize stations arrangement and hit indexes
- fMCData.InitTrackInfo(vHits);
+ fpMCData->InitTrackInfo(vHits);
LOG(info) << "\033[1;32m!!!! FLAG 2\033[0m";
// ----- Define reconstructable and additional flags
- for (auto& aTrk : fMCData.GetTrackContainer()) {
+ for (auto& aTrk : fpMCData->GetTrackContainer()) {
bool isRec = true; // is track reconstructable
// Cut on momentum
@@ -219,7 +204,7 @@ void CbmCaMCModule::InitTrackInfo(const L1Vector<CbmL1HitDebugInfo>& vHits)
isAdd &= aTrk.GetNofConsStationsWithPoint() == aTrk.GetTotNofStationsWithHit();
isAdd &= aTrk.GetTotNofStationsWithHit() == aTrk.GetTotNofStationsWithPoint();
isAdd &= aTrk.GetNofConsStationsWithHit() >= 3;
- isAdd &= fMCData.GetPoint(aTrk.GetPointIndexes()[0]).GetStationId() == 0;
+ isAdd &= fpMCData->GetPoint(aTrk.GetPointIndexes()[0]).GetStationId() == 0;
isAdd &= !isRec;
}
else {
@@ -233,11 +218,7 @@ void CbmCaMCModule::InitTrackInfo(const L1Vector<CbmL1HitDebugInfo>& vHits)
// ---------------------------------------------------------------------------------------------------------------------
//
-void CbmCaMCModule::Finish()
-{
- std::cout << "\033[1;32mFinishing performance\033[0m\n";
- fpQaModule->SaveHistograms();
-}
+void CbmCaMCModule::Finish() {}
// **********************************
@@ -246,20 +227,20 @@ void CbmCaMCModule::Finish()
// ---------------------------------------------------------------------------------------------------------------------
//
-void CbmCaMCModule::MatchPointsWithHits(L1Vector<CbmL1HitDebugInfo>& vInputExtHits)
+void CbmCaMCModule::MatchPointsWithHits()
{
// FIXME: Cleaning up makes it safer, but probably is not needed after old code removal
- for (auto& hit : vInputExtHits) {
+ for (auto& hit : (*fpvHitIds)) {
hit.mcPointIds.clear();
}
- int nHits = vInputExtHits.size();
+ int nHits = fpvHitIds->size();
for (int iH = 0; iH < nHits; ++iH) {
- auto& hit = vInputExtHits[iH];
- int iP = fMCData.GetPointIndexOfHit(iH);
+ auto& hit = (*fpvHitIds)[iH];
+ int iP = fpMCData->GetPointIndexOfHit(iH);
if (iP > -1) {
hit.mcPointIds.push_back_no_warning(iP);
- fMCData.GetPoint(iP).AddHitID(iH);
+ fpMCData->GetPoint(iP).AddHitID(iH);
}
} // loop over hit indexes: end
}
@@ -329,7 +310,7 @@ int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kSts>(int iHit) const
if (link.GetWeight() > bestWeight) {
bestWeight = link.GetWeight();
- iPoint = fMCData.FindInternalPointIndex(CalcGlobMCPointIndex(iIndex, L1DetectorID::kSts), iEvent, iFile);
+ iPoint = fpMCData->FindInternalPointIndex(CalcGlobMCPointIndex(iIndex, L1DetectorID::kSts), iEvent, iFile);
assert(iPoint != -1);
}
}
@@ -354,7 +335,7 @@ int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kMuch>(int iHit) const
}
int iFile = hitMatchMuch->GetLink(iLink).GetFile();
int iEvent = hitMatchMuch->GetLink(iLink).GetEntry();
- iPoint = fMCData.FindInternalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kMuch), iEvent, iFile);
+ iPoint = fpMCData->FindInternalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kMuch), iEvent, iFile);
assert(iPoint != -1);
}
}
@@ -380,7 +361,7 @@ int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kTrd>(int iHit) const
int iFile = hitMatch->GetLink(0).GetFile();
int iEvent = hitMatch->GetLink(0).GetEntry();
- iPoint = fMCData.FindInternalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kTrd), iEvent, iFile);
+ iPoint = fpMCData->FindInternalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kTrd), iEvent, iFile);
assert(iPoint != -1);
}
}
@@ -402,7 +383,7 @@ int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kTof>(int iHit) const
if (iMc < 0) return iPoint;
assert(iMc >= 0 && iMc < fpTofPoints->Size(iFile, iEvent));
- int iPointPrelim = fMCData.FindInternalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kTof), iEvent, iFile);
+ int iPointPrelim = fpMCData->FindInternalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kTof), iEvent, iFile);
if (iPointPrelim == -1) { continue; }
iPoint = iPointPrelim;
}
@@ -412,20 +393,18 @@ int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kTof>(int iHit) const
// ---------------------------------------------------------------------------------------------------------------------
//
-void CbmCaMCModule::MatchRecoAndMCTracks(L1Vector<CbmL1Track>& vRecoTracks, L1Vector<CbmL1HitDebugInfo>& vInputExtHits)
+void CbmCaMCModule::MatchRecoAndMCTracks()
{
- // NOTE: Potentially there is a case, when
-
- for (int iTre = 0; iTre < int(vRecoTracks.size()); ++iTre) {
- auto& trkRe = vRecoTracks[iTre];
+ for (int iTre = 0; iTre < int(fpvRecoTracks->size()); ++iTre) {
+ auto& trkRe = (*fpvRecoTracks)[iTre];
// ----- Count number of hits from each MC track belonging to this reconstructed track
auto& mNofHitsVsMCTrkID = trkRe.hitMap;
mNofHitsVsMCTrkID.clear();
for (int iH : trkRe.Hits) {
- for (int iP : vInputExtHits[iH].mcPointIds) {
+ for (int iP : (*fpvHitIds)[iH].mcPointIds) {
assert(iP > -1); // Should never be triggered
- int iTmc = fMCData.GetPoint(iP).GetTrackId();
+ int iTmc = fpMCData->GetPoint(iP).GetTrackId();
if (mNofHitsVsMCTrkID.find(iTmc) == mNofHitsVsMCTrkID.end()) { mNofHitsVsMCTrkID[iTmc] = 1; }
else {
mNofHitsVsMCTrkID[iTmc] += 1;
@@ -449,7 +428,7 @@ void CbmCaMCModule::MatchRecoAndMCTracks(L1Vector<CbmL1Track>& vRecoTracks, L1Ve
if (double(nHitsTrkMc) > double(nHitsTrkRe) * maxPurity) { maxPurity = double(nHitsTrkMc) / double(nHitsTrkRe); }
- ca::tools::MCTrack& trkMc = fMCData.GetTrack(iTmc);
+ ca::tools::MCTrack& trkMc = fpMCData->GetTrack(iTmc);
// Match reconstructed and MC tracks, if purity with this MC track passes the threshold
if (double(nHitsTrkMc) >= CbmL1Constants::MinPurity * double(nHitsTrkRe)) {
@@ -494,7 +473,15 @@ void CbmCaMCModule::SetDetector(L1DetectorID detID, bool flag)
void CbmCaMCModule::CheckInit() const
{
// ----- Check parameters
- if (!fpParameters) { throw std::logic_error("Tracking parameters object was not defined"); }
+ if (!fpParameters.get()) { throw std::logic_error("Tracking parameters object was not defined"); }
+
+ if (!fpMCData) { throw std::logic_error("MC data object was not registered"); }
+
+ if (!fpvRecoTracks) { throw std::logic_error("Reconstructed track container was not registered"); }
+
+ if (!fpvHitIds) { throw std::logic_error("Hit index container was not registered"); }
+
+ if (!fpInputData) { throw std::logic_error("Input data object was not registered"); }
// Check event list
if (!fbLegacyEventMode && !fpMCEventList) { throw std::logic_error("MC event list was not found"); }
@@ -551,8 +538,8 @@ void CbmCaMCModule::ReadMCPointsForDetector<L1DetectorID::kTof>(CbmMCDataArray*
std::vector<std::vector<int>> vTofPointToTrack(nTofStationsGeo); // TODO
std::vector<std::vector<float>> vTofPointToTrackDZ(nTofStationsGeo); // TODO
for (int iStLocGeo = 0; iStLocGeo < nTofStationsGeo; ++iStLocGeo) {
- vTofPointToTrack[iStLocGeo].resize(fMCData.GetNofTracks(), -1);
- vTofPointToTrackDZ[iStLocGeo].resize(fMCData.GetNofTracks(), 1.e+5f);
+ vTofPointToTrack[iStLocGeo].resize(fpMCData->GetNofTracks(), -1);
+ vTofPointToTrackDZ[iStLocGeo].resize(fpMCData->GetNofTracks(), 1.e+5f);
}
// Array of flags, if a given TOF point is matched
@@ -605,7 +592,7 @@ void CbmCaMCModule::ReadMCPointsForDetector<L1DetectorID::kTof>(CbmMCDataArray*
}
}
- int iTrack = fMCData.FindInternalTrackIndex(pExtPoint->GetTrackID(), iEvent, iFile);
+ int iTrack = fpMCData->FindInternalTrackIndex(pExtPoint->GetTrackID(), iEvent, iFile);
assert(iTrack > -1);
if (iStSelected != -1) {
int iStActive = fpParameters->GetStationIndexActive(iStSelected, L1DetectorID::kTof);
@@ -626,9 +613,9 @@ void CbmCaMCModule::ReadMCPointsForDetector<L1DetectorID::kTof>(CbmMCDataArray*
if (FillMCPoint<L1DetectorID::kTof>(vTofPointToTrack[iStLocGeo][iTrk], iEvent, iFile, aPoint)) {
aPoint.SetStationId(iStActive);
aPoint.SetExternalId(CalcGlobMCPointIndex(vTofPointToTrack[iStLocGeo][iTrk], L1DetectorID::kTof));
- int iPInt = fMCData.GetNofPoints(); // assign an index of point in internal container
- if (aPoint.GetTrackId() > -1) { fMCData.GetTrack(aPoint.GetTrackId()).AddPointIndex(iPInt); }
- fMCData.AddPoint(aPoint);
+ int iPInt = fpMCData->GetNofPoints(); // assign an index of point in internal container
+ if (aPoint.GetTrackId() > -1) { fpMCData->GetTrack(aPoint.GetTrackId()).AddPointIndex(iPInt); }
+ fpMCData->AddPoint(aPoint);
++fvNofPointsUsed[int(L1DetectorID::kTof)];
}
} // loop over tracks: end
@@ -640,8 +627,8 @@ void CbmCaMCModule::ReadMCPointsForDetector<L1DetectorID::kTof>(CbmMCDataArray*
//
void CbmCaMCModule::ReadMCPoints()
{
- int nPointsEstimated = 5 * fMCData.GetNofTracks() * fpParameters->GetNstationsActive();
- fMCData.ReserveNofPoints(nPointsEstimated);
+ int nPointsEstimated = 5 * fpMCData->GetNofTracks() * fpParameters->GetNstationsActive();
+ fpMCData->ReserveNofPoints(nPointsEstimated);
// ----- Initialise the number of points
std::fill(fvNofPointsOrig.begin(), fvNofPointsOrig.end(), 0);
@@ -673,7 +660,7 @@ void CbmCaMCModule::ReadMCTracks()
for (const auto& key : fFileEventIDs) {
nTracksTot += fpMCTracks->Size(key.first, key.second); /// iFile, iEvent
}
- fMCData.ReserveNofTracks(nTracksTot);
+ fpMCData->ReserveNofTracks(nTracksTot);
// ----- Loop over MC events
for (const auto& key : fFileEventIDs) {
@@ -695,8 +682,8 @@ void CbmCaMCModule::ReadMCTracks()
// Create a CbmL1MCTrack
ca::tools::MCTrack aTrk;
- aTrk.SetId(fMCData.GetNofTracks()); // assign current number of tracks read so far as an ID of a new track
- aTrk.SetExternalId(iTrkExt); // external index of track is its index from CbmMCTrack objects container
+ aTrk.SetId(fpMCData->GetNofTracks()); // assign current number of tracks read so far as an ID of a new track
+ aTrk.SetExternalId(iTrkExt); // external index of track is its index from CbmMCTrack objects container
aTrk.SetEventId(iEvent);
aTrk.SetFileId(iFile);
@@ -728,7 +715,7 @@ void CbmCaMCModule::ReadMCTracks()
}
else {
// This is a secondary track, mother ID should be recalculated for the internal array
- int motherId = fMCData.FindInternalTrackIndex(extMotherId, iEvent, iFile);
+ int motherId = fpMCData->FindInternalTrackIndex(extMotherId, iEvent, iFile);
assert(motherId > -1);
aTrk.SetMotherId(motherId);
}
@@ -736,7 +723,7 @@ void CbmCaMCModule::ReadMCTracks()
aTrk.SetProcessId(pExtMCTrk->GetGeantProcessId());
aTrk.SetPdgCode(pExtMCTrk->GetPdgCode());
- fMCData.AddTrack(aTrk);
+ fpMCData->AddTrack(aTrk);
} // Loop over MC tracks: end
} // Loop over MC events: end
}
diff --git a/reco/L1/CbmCaMCModule.h b/reco/L1/CbmCaMCModule.h
index da9658ecfe..224fa60031 100644
--- a/reco/L1/CbmCaMCModule.h
+++ b/reco/L1/CbmCaMCModule.h
@@ -34,8 +34,8 @@
#include "CaToolsMCData.h"
#include "CaToolsMCPoint.h"
-#include "CaToolsQa.h"
#include "L1Constants.h"
+#include "L1InputData.h"
#include "L1Parameters.h"
#include "L1Undef.h"
@@ -52,6 +52,7 @@ class CbmCaMCModule {
public:
/// @brief Constructor
/// @param verbosity Verbosity level
+ /// @param perfMode Performance mode (defines cut on number of consecutive stations with hit or point)
CbmCaMCModule(int verb = 0, int perfMode = 3) : fVerbose(verb), fPerformanceMode(perfMode) {}
/// @brief Destructor
@@ -84,14 +85,6 @@ public:
return std::accumulate(fvNofPointsUsed.cbegin(), fvNofPointsUsed.cbegin() + int(detID) + 1, 0) - 1;
}
- /// @brief Gets MC data object
- /// @return Constant reference to MC data object
- const ca::tools::MCData& GetMCData() const { return fMCData; }
-
- /// @brief Gets mutable pointer to MC data object
- /// @return Mutable reference to MC data object
- ca::tools::MCData* GetMCData() { return &fMCData; }
-
/// @brief Defines performance action in the beginning of each event or time slice
/// @note This function should be called before hits initialization
/// @param pEvent Pointer to a current CbmEvent
@@ -126,12 +119,9 @@ public:
int MatchHitWithMc(int iHit) const;
/// @brief Assigns MC point indexes to each hit and hit indexes to each MC point
- /// @param vInputExtHits Vector of external input hits
- void MatchPointsWithHits(L1Vector<CbmL1HitDebugInfo>& vInputExtHits);
+ void MatchPointsWithHits();
/// @brief Matches reconstructed tracks with MC tracks
- /// @param vRecoTracks Vector of reconstructed tracks
- /// @param vInputExtHits Vector of external input hits
///
/// In the procedure, the maps of associated MC track indexes to corresponding number of hits are filled out and the
/// max purity is calculated for each reconstructed track in the TS/event. If the value of purity for a given MC track
@@ -139,7 +129,7 @@ public:
/// time the index of the reconstructed track object is saved to this MC track object. If purity for the MC track does
/// not pass the threshold, its index will not be accounted as a matched track, and the index of reconstructed track
/// will be added as an index of "touched" track.
- void MatchRecoAndMCTracks(L1Vector<CbmL1Track>& vRecoTracks, L1Vector<CbmL1HitDebugInfo>& vInputExtHits);
+ void MatchRecoAndMCTracks();
/// @brief Sets used detector subsystems
/// @param detID Id of detector
@@ -147,20 +137,36 @@ public:
/// @note Should be called before this->Init()
void SetDetector(L1DetectorID detID, bool flag);
- /// @brief Sets pointer to the vector of hits from reference
- void SetHitContainer(const L1Vector<CbmL1HitDebugInfo>& vHits) { fpvHits = &vHits; }
-
/// @brief Sets legacy event mode:
/// @param flag Flag:
/// - true: runs on events base,
/// - false: runs on time-slice base
void SetLegacyEventMode(bool flag) { fbLegacyEventMode = flag; }
- /// @brief Registers pointer to the tracking parameters object
- void SetParameters(const L1Parameters* pParameters) { fpParameters = pParameters; }
+ /// @brief Registers MC data object
+ /// @param mcData Instance of MC data
+ void RegisterMCData(ca::tools::MCData& mcData) { fpMCData = &mcData; }
+
+ /// @brief Registers reconstructed track container
+ /// @param vRecoTrack Reference to reconstructed track container
+ void RegisterRecoTrackContainer(L1Vector<CbmL1Track>& vRecoTracks) { fpvRecoTracks = &vRecoTracks; }
+
+ /// @brief Registers hit index container
+ /// @param vHitIds Reference to hit index container
+ void RegisterHitIndexContainer(L1Vector<CbmL1HitDebugInfo>& vHitIds) { fpvHitIds = &vHitIds; }
+
+ /// @brief Registers input data object
+ /// @param inputData Instance of the input data object
+ void RegisterInputData(L1InputData& inputData) { fpInputData = &inputData; }
+
+ // /// @brief Registers debug hit container
+ // /// @param vDbgHits Reference to debug hit container
+ // void RegisterHitDebugInfoContainer(L1Vector<CbmL1HitStore>& vDbgHits) { fpvDbgHits = &vDbgHits; }
+
+ /// @brief Registers CA parameters object
+ /// @param pParameters A shared pointer to the parameters object
+ void RegisterParameters(std::shared_ptr<L1Parameters>& pParameters) { fpParameters = pParameters; }
- /// @brief Sets pointer to the vector of reconstructed tracks from reference
- void SetRecoTrackContainer(const L1Vector<CbmL1Track>& vRecoTracks) { fpvRecoTracks = &vRecoTracks; }
private:
/// @brief Calculates global index of MC point
@@ -200,14 +206,9 @@ private:
bool FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MCPoint& point);
- // *******************************
- // ** Utility variables **
- // *******************************
-
- const L1Parameters* fpParameters = nullptr; ///< Pointer to tracking parameters object
-
- // ** Flags **
+ std::shared_ptr<L1Parameters> fpParameters = nullptr; ///< Pointer to tracking parameters object
+ // ------ Flags
bool fbUseMvd = false; ///< is MVD used
bool fbUseSts = false; ///< is STS used
bool fbUseMuch = false; ///< is MuCh used
@@ -218,18 +219,15 @@ private:
int fVerbose = 0; ///< Verbosity level
int fPerformanceMode = -1; ///< Mode of performance
- // *********************************
- // ** Input data branches **
- // *********************************
-
+ // ------ Input data branches
const CbmTimeSlice* fpTimeSlice = nullptr; ///< Current time slice
- // ------ Mc-event
+ // Mc-event
CbmMCEventList* fpMCEventList = nullptr; ///< MC event list
CbmMCDataObject* fpMCEventHeader = nullptr; ///< MC event header
CbmMCDataArray* fpMCTracks = nullptr; ///< MC tracks input
- // ------ Mc-points
+ // Mc-points
CbmMCDataArray* fpMvdPoints = nullptr; ///< MVD MC-points input container
CbmMCDataArray* fpStsPoints = nullptr; ///< STS MC-points input container
CbmMCDataArray* fpMuchPoints = nullptr; ///< MuCh MC-points input container
@@ -239,7 +237,7 @@ private:
std::array<int, L1Constants::size::kMaxNdetectors> fvNofPointsOrig = {}; ///< Number of points by detector provided
std::array<int, L1Constants::size::kMaxNdetectors> fvNofPointsUsed = {}; ///< Number of points used in performance
- // ------ Matches
+ // Matches
TClonesArray* fpMvdHitMatches = nullptr; ///< Array of MVD hit matches
TClonesArray* fpStsHitMatches = nullptr; ///< Array of STS hit matches
TClonesArray* fpMuchHitMatches = nullptr; ///< Array of MuCh hit matches
@@ -249,23 +247,16 @@ private:
TClonesArray* fpStsHits = nullptr; ///< Array of STS hits (currently needed for matching)
TClonesArray* fpStsClusterMatches = nullptr; ///< Array of STS cluster matches
- // ------ Matching information
+ // Matching information
std::set<std::pair<int, int>> fFileEventIDs; ///< Set of file and event indexes: first - iFile, second - iEvent
- // ************************************
- // ** Current MC data object **
- // ************************************
-
- ca::tools::MCData fMCData {}; ///< MC-information in CA tracking internal format (tracks and points)
-
- std::unique_ptr<ca::tools::Qa> fpQaModule = nullptr; ///< Pointer to CA QA module
+ ca::tools::MCData* fpMCData = nullptr; ///< MC-information in CA tracking internal format (tracks and points)
- // ********************************************
- // ** Pointers to reconstructed data **
- // ********************************************
-
- const L1Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to reconstructed track container
- const L1Vector<CbmL1HitDebugInfo>* fpvHits = nullptr; ///< Pointer to hits container
+ // Reconstructed data
+ L1Vector<CbmL1Track>* fpvRecoTracks = nullptr; ///< Pointer to reconstructed track container
+ L1Vector<CbmL1HitDebugInfo>* fpvHitIds = nullptr; ///< Pointer to hit index container
+ //L1Vector<CbmL1HitStore>* fpvDbgHits = nullptr; ///< Pointer to debug hit container
+ L1InputData* fpInputData = nullptr; ///< Pointer to input data object
};
// **********************************************
@@ -286,7 +277,7 @@ bool CbmCaMCModule::FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MC
int iTmcExt = undef::kI32; // Track ID in the external container
// MVD
if constexpr (L1DetectorID::kMvd == DetID) {
- auto* pExtPoint = L1_DYNAMIC_CAST<CbmMvdPoint*>(fpMvdPoints->Get(iFile, iEvent, iExtId));
+ auto* pExtPoint = dynamic_cast<CbmMvdPoint*>(fpMvdPoints->Get(iFile, iEvent, iExtId));
if (!pExtPoint) {
LOG(warn) << "CbmCaMCModule: MVD MC point with iExtId = " << iExtId << ", iEvent = " << iEvent
<< ", iFile = " << iFile << " does not exist";
@@ -301,7 +292,7 @@ bool CbmCaMCModule::FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MC
}
// STS
else if constexpr (L1DetectorID::kSts == DetID) {
- auto* pExtPoint = L1_DYNAMIC_CAST<CbmStsPoint*>(fpStsPoints->Get(iFile, iEvent, iExtId));
+ auto* pExtPoint = dynamic_cast<CbmStsPoint*>(fpStsPoints->Get(iFile, iEvent, iExtId));
if (!pExtPoint) {
LOG(warn) << "CbmCaMCModule: STS MC point with iExtId = " << iExtId << ", iEvent = " << iEvent
<< ", iFile = " << iFile << " does not exist";
@@ -316,7 +307,7 @@ bool CbmCaMCModule::FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MC
}
// MuCh
else if constexpr (L1DetectorID::kMuch == DetID) {
- auto* pExtPoint = L1_DYNAMIC_CAST<CbmMuchPoint*>(fpMuchPoints->Get(iFile, iEvent, iExtId));
+ auto* pExtPoint = dynamic_cast<CbmMuchPoint*>(fpMuchPoints->Get(iFile, iEvent, iExtId));
if (!pExtPoint) {
LOG(warn) << "CbmCaMCModule: MuCh MC point with iExtId = " << iExtId << ", iEvent = " << iEvent
<< ", iFile = " << iFile << " does not exist";
@@ -331,7 +322,7 @@ bool CbmCaMCModule::FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MC
}
// TRD
else if constexpr (L1DetectorID::kTrd == DetID) {
- auto* pExtPoint = L1_DYNAMIC_CAST<CbmTrdPoint*>(fpTrdPoints->Get(iFile, iEvent, iExtId));
+ auto* pExtPoint = dynamic_cast<CbmTrdPoint*>(fpTrdPoints->Get(iFile, iEvent, iExtId));
if (!pExtPoint) {
LOG(warn) << "CbmCaMCModule: TRD MC point with iExtId = " << iExtId << ", iEvent = " << iEvent
<< ", iFile = " << iFile << " does not exist";
@@ -346,7 +337,7 @@ bool CbmCaMCModule::FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MC
}
// TOF
else if constexpr (L1DetectorID::kTof == DetID) {
- auto* pExtPoint = L1_DYNAMIC_CAST<CbmTofPoint*>(fpTofPoints->Get(iFile, iEvent, iExtId));
+ auto* pExtPoint = dynamic_cast<CbmTofPoint*>(fpTofPoints->Get(iFile, iEvent, iExtId));
if (!pExtPoint) {
LOG(warn) << "CbmCaMCModule: TOF MC point with iExtId = " << iExtId << ", iEvent = " << iEvent
<< ", iFile = " << iFile << " does not exist";
@@ -418,34 +409,28 @@ bool CbmCaMCModule::FillMCPoint(int iExtId, int iEvent, int iFile, ca::tools::MC
point.SetEventId(iEvent);
point.SetFileId(iFile);
point.SetTime(time);
-
point.SetX(posMid.X());
point.SetY(posMid.Y());
point.SetZ(posMid.Z());
-
point.SetXIn(posIn.X());
point.SetYIn(posIn.Y());
point.SetZIn(posIn.Z());
-
point.SetXOut(posOut.X());
point.SetYOut(posOut.Y());
point.SetZOut(posOut.Z());
-
point.SetPx(momMid.X());
point.SetPy(momMid.Y());
point.SetPz(momMid.Z());
-
point.SetPxIn(momIn.X());
point.SetPyIn(momIn.Y());
point.SetPzIn(momIn.Z());
-
point.SetPxOut(momOut.X());
point.SetPyOut(momOut.Y());
point.SetPzOut(momOut.Z());
- int iTmcInt = fMCData.FindInternalTrackIndex(iTmcExt, iEvent, iFile);
+ int iTmcInt = fpMCData->FindInternalTrackIndex(iTmcExt, iEvent, iFile);
- point.SetId(fMCData.GetNofPoints()); // select current number of points as a local id of points
+ point.SetId(fpMCData->GetNofPoints()); // select current number of points as a local id of points
point.SetTrackId(iTmcInt);
point.SetStationId(stationID);
point.SetDetectorId(DetID);
@@ -481,9 +466,9 @@ void CbmCaMCModule::ReadMCPointsForDetector(CbmMCDataArray* pPoints)
ca::tools::MCPoint aPoint;
if (FillMCPoint<DetID>(iP, iEvent, iFile, aPoint)) {
aPoint.SetExternalId(CalcGlobMCPointIndex(iP, DetID));
- int iPInt = fMCData.GetNofPoints(); // assign an index of point in internal container
- if (aPoint.GetTrackId() > -1) { fMCData.GetTrack(aPoint.GetTrackId()).AddPointIndex(iPInt); }
- fMCData.AddPoint(aPoint);
+ int iPInt = fpMCData->GetNofPoints(); // assign an index of point in internal container
+ if (aPoint.GetTrackId() > -1) { fpMCData->GetTrack(aPoint.GetTrackId()).AddPointIndex(iPInt); }
+ fpMCData->AddPoint(aPoint);
++fvNofPointsUsed[int(DetID)];
}
} // iP: end
diff --git a/reco/L1/CbmCaTimeSliceReader.cxx b/reco/L1/CbmCaTimeSliceReader.cxx
index e4f7fa6b59..8c85c555b1 100644
--- a/reco/L1/CbmCaTimeSliceReader.cxx
+++ b/reco/L1/CbmCaTimeSliceReader.cxx
@@ -59,9 +59,6 @@ void TimeSliceReader::Clear()
fpvDbgHits = nullptr;
fpIODataManager = nullptr;
fpvTracks = nullptr;
-
- // Other
- fpParameters = nullptr;
}
// ---------------------------------------------------------------------------------------------------------------------
@@ -261,14 +258,6 @@ void TimeSliceReader::SetDetector(L1DetectorID detID, bool flag)
}
}
-// ---------------------------------------------------------------------------------------------------------------------
-//
-void TimeSliceReader::SetParameters(const L1Parameters* pParameters)
-{
- LOG_IF(fatal, !pParameters) << "TimeSliceReader: passed null pointer as L1IODataManager instance";
- fpParameters = pParameters;
-}
-
// ---------------------------------------------------------------------------------------------------------------------
//
void TimeSliceReader::StoreHitRecord(const HitRecord& hitRecord)
diff --git a/reco/L1/CbmCaTimeSliceReader.h b/reco/L1/CbmCaTimeSliceReader.h
index 808bd6ada2..83a4e548e6 100644
--- a/reco/L1/CbmCaTimeSliceReader.h
+++ b/reco/L1/CbmCaTimeSliceReader.h
@@ -126,10 +126,6 @@ namespace cbm::ca
/// @note Should be called before this->Init()
void SetDetector(L1DetectorID detID, bool flag = true);
- /// @brief Sets the instance of parameters object
- /// @param pParameters Pointer to the parameters instance
- void SetParameters(const L1Parameters* pParameters);
-
private:
/// @brief Reads hits for a given detector subsystem
/// @tparam Detector ID
@@ -184,10 +180,8 @@ namespace cbm::ca
L1Vector<CbmL1Track>* fpvTracks = nullptr; ///< Pointer to array of reconstructed tracks
std::shared_ptr<L1IODataManager> fpIODataManager = nullptr; ///< Pointer to input data manager
-
// Additional
- const L1Parameters* fpParameters = nullptr; ///< Pointer to the defined parameters instance
- ECbmTrackingMode fTrackingMode; ///< Tracking mode
+ ECbmTrackingMode fTrackingMode; ///< Tracking mode
// Indexes
int fNofHits = 0; ///< Stored number of hits
diff --git a/reco/L1/qa/CbmCaOutputQa.cxx b/reco/L1/qa/CbmCaOutputQa.cxx
index 72df21eb8a..b8a7fd83cd 100644
--- a/reco/L1/qa/CbmCaOutputQa.cxx
+++ b/reco/L1/qa/CbmCaOutputQa.cxx
@@ -32,13 +32,16 @@ void OutputQa::EnableDebugger(const char* filename)
InitStatus OutputQa::InitDataBranches()
{
LOG(info) << fName << ": Initializing data branches";
-
- if (!fpTSReader.get()) { fpTSReader = std::make_unique<TimeSliceReader>(fTrackingMode); }
+ LOG_IF(fatal, !fpParameters.get())
+ << fName << ": CA parameters object is not defined. Please, provide initializer or read parameters from binary "
+ << "via OutputQa::ReadParameters(filename) from the qa macro";
+ // Initialize IO data manager
if (!fpDataManager.get()) { fpDataManager = std::make_shared<L1IODataManager>(); }
// Initialize time slice reader instance
+ if (!fpTSReader.get()) { fpTSReader = std::make_unique<TimeSliceReader>(fTrackingMode); }
fpTSReader->SetDetector(L1DetectorID::kMvd, fbUseMvd);
fpTSReader->SetDetector(L1DetectorID::kSts, fbUseSts);
fpTSReader->SetDetector(L1DetectorID::kMuch, fbUseMuch);
@@ -48,9 +51,29 @@ InitStatus OutputQa::InitDataBranches()
fpTSReader->RegisterIODataManager(fpDataManager);
fpTSReader->RegisterTracksContainer(fvRecoTracks);
fpTSReader->RegisterHitDebugInfoContainer(fvDbgHits);
- fpTSReader->RegisterHitIndexContainer(fvHitIds);
+ //fpTSReader->RegisterHitIndexContainer(fvHitIds);
fpTSReader->InitRun();
+
+ // Initialize MC module
+ if (IsMCUsed()) {
+ if (!fpMCModule.get()) { fpMCModule = std::make_unique<CbmCaMCModule>(fVerbose, fPerformanceMode); }
+ fpMCModule->SetDetector(L1DetectorID::kMvd, fbUseMvd);
+ fpMCModule->SetDetector(L1DetectorID::kSts, fbUseSts);
+ fpMCModule->SetDetector(L1DetectorID::kMuch, fbUseMuch);
+ fpMCModule->SetDetector(L1DetectorID::kTrd, fbUseTrd);
+ fpMCModule->SetDetector(L1DetectorID::kTof, fbUseTof);
+
+ fpMCModule->RegisterInputData(fInputData);
+ fpMCModule->RegisterMCData(fMCData);
+ fpMCModule->RegisterRecoTrackContainer(fvRecoTracks);
+ fpMCModule->RegisterHitIndexContainer(fvHitIds);
+ fpMCModule->RegisterParameters(fpParameters);
+
+ fpMCModule->InitRun();
+ }
+
+
return kSUCCESS;
}
@@ -58,14 +81,30 @@ InitStatus OutputQa::InitDataBranches()
//
InitStatus OutputQa::InitTimeSlice()
{
+ int nMCTracks = -1;
+ int nMCPoints = -1;
+ int nHits = -1;
+ int nRecoTracks = -1;
+
// Read hits: fill fvHitIds, fvDbgHits and fpDataManager
fpTSReader->ReadHits();
fpDataManager->SendInputData(fInputData);
// Read reconstructed tracks
fpTSReader->ReadRecoTracks();
- LOG(info) << fName << ": Time slice was read: " << fInputData.GetNhits() << " hits, " << fvRecoTracks.size()
- << " tracks";
+
+ nHits = fInputData.GetNhits();
+ nRecoTracks = fvRecoTracks.size();
+
+ if (IsMCUsed()) {
+ // Read MC information
+ fpMCModule->InitEvent(nullptr);
+ nMCTracks = fMCData.GetNofTracks();
+ nMCPoints = fMCData.GetNofPoints();
+ }
+
+ LOG(info) << fName << ": Event or time slice consists from " << nHits << " hits, " << nRecoTracks << " reco tracks, "
+ << nMCTracks << " MC tracks, " << nMCPoints << " MC points";
return kSUCCESS;
}
diff --git a/reco/L1/qa/CbmCaOutputQa.h b/reco/L1/qa/CbmCaOutputQa.h
index 09c2efbe96..eafee38a12 100644
--- a/reco/L1/qa/CbmCaOutputQa.h
+++ b/reco/L1/qa/CbmCaOutputQa.h
@@ -66,6 +66,9 @@ namespace cbm::ca
/// @brief Sets mCBM global tracking mode
void SetMcbmTrackingMode() { fTrackingMode = ECbmTrackingMode::kMCBM; }
+ /// @brief Sets performance mode
+ /// @param mode Performance mode (3 is default: TODO - test)
+ void SetPerformanceMode(int mode) { fPerformanceMode = mode; }
ClassDef(OutputQa, 0);
@@ -102,6 +105,8 @@ namespace cbm::ca
bool fbUseTrd = false; ///< is TRD used
bool fbUseTof = false; ///< is TOF used
+ int fPerformanceMode = 3; ///< Performance mode
+
ECbmTrackingMode fTrackingMode = ECbmTrackingMode::kSTS; ///< Tracking mode
std::unique_ptr<TimeSliceReader> fpTSReader = nullptr; ///< Reader of the time slice
@@ -110,10 +115,11 @@ namespace cbm::ca
std::shared_ptr<::ca::tools::Debugger> fpDebugger = nullptr; ///< Debugger
std::shared_ptr<L1Parameters> fpParameters = nullptr; ///< Tracking parameters object
- L1Vector<CbmL1HitId> fvHitIds {"CbmCaOutputQa::fvHitIds"};
+ L1Vector<CbmL1HitDebugInfo> fvHitIds {"CbmCaOutputQa::fvHitIds"};
L1Vector<CbmL1HitDebugInfo> fvDbgHits {"CbmCaOutputQa::fvDbgHits"};
L1Vector<CbmL1Track> fvRecoTracks {"CbmCaOutputQa::fvRecoTracks"};
- L1InputData fInputData; ///< Input hits to L1 algo
+ L1InputData fInputData; ///< Input hits to CA algo
+ ::ca::tools::MCData fMCData; ///< Input MC data (points and tracks)
};
} // namespace cbm::ca
--
GitLab