From c0192bca2e6006e0c419c99202935115b62db93d Mon Sep 17 00:00:00 2001
From: "s.zharko@gsi.de" <s.zharko@gsi.de>
Date: Tue, 13 Dec 2022 19:28:46 +0100
Subject: [PATCH] L1: Updates on external MC module for tracking
---
reco/L1/CbmCaMCModule.cxx | 220 ++++++++++++++++++++++++------
reco/L1/CbmL1ReadEvent.cxx | 11 +-
reco/L1/catools/CaToolsMCData.cxx | 22 +--
reco/L1/catools/CaToolsMCData.h | 86 +++++++-----
reco/L1/catools/CaToolsMCPoint.h | 2 -
5 files changed, 248 insertions(+), 93 deletions(-)
diff --git a/reco/L1/CbmCaMCModule.cxx b/reco/L1/CbmCaMCModule.cxx
index 001fe6a53c..d2e0e03876 100644
--- a/reco/L1/CbmCaMCModule.cxx
+++ b/reco/L1/CbmCaMCModule.cxx
@@ -11,10 +11,12 @@
#include "CbmEvent.h"
#include "CbmL1.h" // for L1DetectorID
+#include "CbmL1Hit.h"
#include "CbmMCDataManager.h"
#include "CbmMCDataObject.h"
#include "CbmMCEventList.h"
#include "CbmMCTrack.h"
+#include "CbmStsHit.h"
#include "CbmTimeSlice.h"
#include "FairEventHeader.h"
@@ -70,6 +72,7 @@ try {
fpTofHitMatches = nullptr;
fpStsClusterMatches = nullptr;
+ fpStsHits = nullptr;
fFileEventIDs.clear();
@@ -84,6 +87,7 @@ try {
fpStsPoints = mcManager->InitBranch("StsPoint");
fpStsHitMatches = dynamic_cast<TClonesArray*>(fairManager->GetObject("StsHitMatch"));
fpStsClusterMatches = dynamic_cast<TClonesArray*>(fairManager->GetObject("StsClusterMatch"));
+ fpStsHits = dynamic_cast<TClonesArray*>(fairManager->GetObject("StsHit"));
}
if (fbUseMuch) {
@@ -126,7 +130,7 @@ catch (const std::logic_error& error) {
// ---------------------------------------------------------------------------------------------------------------------
//
-void CbmCaMCModule::InitEvent(CbmEvent* pEvent)
+void CbmCaMCModule::InitEvent(CbmEvent* /*pEvent*/)
{
std::cout << "\033[1;32mCALL CbmCAMCModule::InitEvent\033[0m\n";
@@ -157,6 +161,172 @@ void CbmCaMCModule::InitEvent(CbmEvent* pEvent)
void CbmCaMCModule::Finish() { std::cout << "\033[1;32mFinishing performance\033[0m\n"; }
+// *****************************************
+// ** Hits and MC-points matching **
+// *****************************************
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+void CbmCaMCModule::MatchPointsWithHits(L1Vector<CbmL1Hit>& vHits)
+{
+ int nHits = vHits.size();
+ for (int iH = 0; iH < nHits; ++iH) {
+ auto& hit = vHits[iH];
+ int iP = fMCData.GetPointIndexOfHit(iH);
+ if (iP > -1) {
+ hit.mcPointIds.push_back_no_warning(iP);
+ fMCData.GetPoint(iP).AddHitID(iH);
+ }
+ } // loop over hit indexes: end
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<>
+int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kMvd>(int iHit) const
+{
+ int iPoint = -1;
+ if (fpMvdHitMatches) {
+ int iHitExt = -(1 + iHit); // TODO: SZh 28.08.2022: this should be replaced with iHitExt = hit.extIdex
+ const auto* hitMatch = dynamic_cast<CbmMatch*>(fpMvdHitMatches->At(iHitExt));
+ assert(hitMatch);
+ if (hitMatch->GetNofLinks() > 0 && hitMatch->GetLink(0).GetIndex() < fvNofPointsUsed[int(L1DetectorID::kMvd)]) {
+ iPoint = hitMatch->GetLink(0).GetIndex();
+ }
+ }
+ return iPoint;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<>
+int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kSts>(int iHit) const
+{
+ int iPoint = -1;
+ const auto* sh = dynamic_cast<CbmStsHit*>(fpStsHits->At(iHit));
+
+ // Match MC point
+ if (fpStsClusterMatches) {
+ const auto* clusterMatchF = static_cast<const CbmMatch*>(fpStsClusterMatches->At(sh->GetFrontClusterId()));
+ const auto* clusterMatchB = static_cast<const CbmMatch*>(fpStsClusterMatches->At(sh->GetBackClusterId()));
+ CbmMatch hitMatch;
+ for (int iLinkF = 0; iLinkF < clusterMatchF->GetNofLinks(); ++iLinkF) {
+ const auto& linkF = clusterMatchF->GetLink(iLinkF);
+ for (int iLinkB = 0; iLinkB < clusterMatchB->GetNofLinks(); ++iLinkB) {
+ const auto& linkB = clusterMatchB->GetLink(iLinkB);
+ if (linkF == linkB) {
+ hitMatch.AddLink(linkF);
+ hitMatch.AddLink(linkB);
+ }
+ }
+ }
+ float bestWeight = 0.f;
+ for (int iLink = 0; iLink < hitMatch.GetNofLinks(); ++iLink) {
+ const CbmLink& link = hitMatch.GetLink(iLink);
+
+ int iIndex = link.GetIndex();
+ if (iIndex < 0) {
+ iPoint = -1;
+ break;
+ }
+ int iFile = link.GetFile();
+ int iEvent = link.GetEntry();
+
+ if (fbLegacyEventMode) {
+ iFile = fFileEventIDs.begin()->first;
+ iEvent = fFileEventIDs.begin()->second;
+
+ // NOTE: Assertions below are temporary for testing consistency between link.GetFile() and
+ // fvFileEvent.begin()->first. It seems like this if-block is not needed, because
+ // assertions do not work out for a single file case (SZh)
+ assert(iFile == link.GetFile());
+ assert(iEvent == link.GetEntry());
+ }
+
+ if (link.GetWeight() > bestWeight) {
+ bestWeight = link.GetWeight();
+ iPoint = fMCData.FindGlobalPointIndex(CalcGlobMCPointIndex(iIndex, L1DetectorID::kSts), iEvent, iFile);
+ assert(iPoint != -1);
+ }
+ }
+ } // Match MC point
+ return iPoint;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<>
+int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kMuch>(int iHit) const
+{
+ int iPoint = -1;
+ const auto* hitMatchMuch = dynamic_cast<CbmMatch*>(fpMuchHitMatches->At(iHit));
+ if (hitMatchMuch) {
+ for (int iLink = 0; iLink < hitMatchMuch->GetNofLinks(); ++iLink) {
+ if (hitMatchMuch->GetLink(iLink).GetIndex() < fvNofPointsUsed[int(L1DetectorID::kMuch)]) {
+ int iMc = hitMatchMuch->GetLink(iLink).GetIndex();
+ if (iMc < 0) {
+ iPoint = -1;
+ break;
+ }
+ int iFile = hitMatchMuch->GetLink(iLink).GetFile();
+ int iEvent = hitMatchMuch->GetLink(iLink).GetEntry();
+ iPoint = fMCData.FindGlobalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kMuch), iEvent, iFile);
+ assert(iPoint != -1);
+ }
+ }
+ }
+ return iPoint;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<>
+int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kTrd>(int iHit) const
+{
+ int iPoint = -1;
+ const auto* hitMatch = dynamic_cast<const CbmMatch*>(fpTrdHitMatches->At(iHit));
+ if (hitMatch) {
+ int iMC = -1;
+ if (hitMatch->GetNofLinks() > 0) {
+ iMC = hitMatch->GetLink(0).GetIndex();
+ assert(iMC >= 0 && iMC < fvNofPointsUsed[int(L1DetectorID::kTrd)]);
+
+ int iMc = hitMatch->GetLink(0).GetIndex();
+ if (iMc < 0) return iPoint;
+ int iFile = hitMatch->GetLink(0).GetFile();
+ int iEvent = hitMatch->GetLink(0).GetEntry();
+
+ iPoint = fMCData.FindGlobalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kTrd), iEvent, iFile);
+ assert(iPoint != -1);
+ }
+ }
+ return iPoint;
+}
+
+// ---------------------------------------------------------------------------------------------------------------------
+//
+template<>
+int CbmCaMCModule::MatchHitWithMc<L1DetectorID::kTof>(int iHit) const
+{
+ int iPoint = -1;
+ const auto* hitMatch = dynamic_cast<const CbmMatch*>(fpTofHitMatches->At(iHit));
+ if (hitMatch) {
+ for (int iLink = 0; iLink < hitMatch->GetNofLinks(); ++iLink) {
+ int iMc = hitMatch->GetLink(iLink).GetIndex();
+ int iFile = hitMatch->GetLink(iLink).GetFile();
+ int iEvent = hitMatch->GetLink(iLink).GetEntry();
+ if (iMc < 0) return iPoint;
+
+ assert(iMc >= 0 && iMc < fpTofPoints->Size(iFile, iEvent));
+ int iPointPrelim = fMCData.FindGlobalPointIndex(CalcGlobMCPointIndex(iMc, L1DetectorID::kTof), iEvent, iFile);
+ if (iPointPrelim == -1) { continue; }
+ iPoint = iPointPrelim;
+ }
+ }
+ return iPoint;
+}
+
+
// ***********************
// ** Accessors **
// ***********************
@@ -174,7 +344,6 @@ void CbmCaMCModule::SetDetector(L1DetectorID detID, bool flag)
}
}
-
// *******************************
// ** Utility functions **
// *******************************
@@ -207,6 +376,7 @@ void CbmCaMCModule::CheckInit() const
if (!fpStsPoints) { throw std::logic_error("MC points unavailable for STS"); }
if (!fpStsHitMatches) { throw std::logic_error("Hit matches unavailable for STS"); }
if (!fpStsClusterMatches) { throw std::logic_error("Cluster matches unavailable for STS"); }
+ if (!fpStsHits) { throw std::logic_error("Hits unavailable for STS"); }
}
if (fbUseMuch) {
@@ -263,28 +433,22 @@ void CbmCaMCModule::ReadMCPointsForDetector<L1DetectorID::kTof>(CbmMCDataArray*
}
int nPointsEvent = pPoints->Size(iFile, iEvent);
- int counter = 0;
- int counter1 = 0;
- // loop over all TOF points in event and file
+ // loop over all TOF points in event and ffvNofPointsUsed[int(L1DetectorID::kMvd)]ile
for (int iP = 0; iP < nPointsEvent; ++iP) {
- if (!vbTofPointMatched[iP]) {
- counter1++;
- continue;
- }
+ if (!vbTofPointMatched[iP]) { continue; } // iP was not matched
auto* pExtPoint = L1_DYNAMIC_CAST<CbmTofPoint*>(pPoints->Get(iFile, iEvent, iP));
LOG_IF(fatal, !pExtPoint) << "NO POINT: TOF: file, event, index = " << iFile << ' ' << iEvent << ' ' << iP
<< '\n';
- if (!pExtPoint) { continue; }
- if (pExtPoint->GetTrackID() < 0) { continue; }
+ if (pExtPoint->GetTrackID() < 0) { continue; } // Point does not have a track
// Cut on time-slice time
if (!fbLegacyEventMode) {
double pointT = pExtPoint->GetTime() + fpMCEventList->GetEventTime(iEvent, iFile);
double startT = fpTimeSlice->GetStartTime();
double endT = fpTimeSlice->GetEndTime();
- if ((startT > 0. && pointT < startT) || (endT > 0. && pointT > endT)) { continue; }
+ if ((startT > 0. && pointT < startT) || (endT > 0. && pointT > endT)) { continue; } // does not fit into TS
}
// Select station index for a point
@@ -293,7 +457,7 @@ void CbmCaMCModule::ReadMCPointsForDetector<L1DetectorID::kTof>(CbmMCDataArray*
int iStSelected = -1; // local geometry index of TOF station
for (int iStLocGeo = 0; iStLocGeo < nTofStationsGeo; ++iStLocGeo) {
int iStActive = fpParameters->GetStationIndexActive(iStLocGeo, L1DetectorID::kTof);
- if (iStActive < 0) { continue; }
+ if (iStActive < 0) { continue; } // station is not used in tracking
const L1Station& station = fpParameters->GetStation(iStActive);
if (fabs(zPos - station.fZ[0]) < bestDist) {
bestDist = fabs(zPos - station.fZ[0]);
@@ -356,36 +520,6 @@ void CbmCaMCModule::ReadMCPoints()
if (fbUseMuch) { this->ReadMCPointsForDetector<L1DetectorID::kMuch>(fpMuchPoints); }
if (fbUseTrd) { this->ReadMCPointsForDetector<L1DetectorID::kTrd>(fpTrdPoints); }
if (fbUseTof) { this->ReadMCPointsForDetector<L1DetectorID::kTof>(fpTofPoints); }
-
-
- // ******* DEBUG: START
- static int nCallsDBG = 0;
- if (nCallsDBG < 1) {
- auto openMode = nCallsDBG ? std::ios_base::app : std::ios_base::out;
- std::ofstream out("mc_points_dbg_new.txt", openMode);
- out << "N of points orig:\n";
- out << "\tMVD: " << fvNofPointsOrig[0] << '\n';
- out << "\tSTS: " << fvNofPointsOrig[1] << '\n';
- out << "\tMuCh: " << fvNofPointsOrig[2] << '\n';
- out << "\tTRD: " << fvNofPointsOrig[3] << '\n';
- out << "\tTOF: " << fvNofPointsOrig[4] << '\n';
- out << "N of points used:\n";
- out << "\tMVD: " << fvNofPointsUsed[0] << '\n';
- out << "\tSTS: " << fvNofPointsUsed[1] << '\n';
- out << "\tMuCh: " << fvNofPointsUsed[2] << '\n';
- out << "\tTRD: " << fvNofPointsUsed[3] << '\n';
- out << "\tTOF: " << fvNofPointsUsed[4] << '\n';
- out << "\ttotal: " << std::accumulate(fvNofPointsUsed.begin(), fvNofPointsUsed.end(), 0) << '\n';
- out << "POINTS:\n";
-
- out << fMCData.GetPoint(0).ToString(2, true) << '\n';
- for (const auto& point : fMCData.GetPointContainer()) {
- out << point.ToString(2) << '\n';
- }
- out.close();
- ++nCallsDBG;
- }
- // ******* DEBUG: END
}
// ---------------------------------------------------------------------------------------------------------------------
diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx
index 817f91c73a..38103534ca 100644
--- a/reco/L1/CbmL1ReadEvent.cxx
+++ b/reco/L1/CbmL1ReadEvent.cxx
@@ -272,12 +272,12 @@ int CbmL1::MatchHitWithMc<L1DetectorID::kTof>(int iHit) const
if (hitMatch) {
for (int iLink = 0; iLink < hitMatch->GetNofLinks(); ++iLink) {
int iMc = hitMatch->GetLink(iLink).GetIndex();
- if (iMc < 0) return iPoint;
int iFile = hitMatch->GetLink(iLink).GetFile();
int iEvent = hitMatch->GetLink(iLink).GetEntry();
-
+ if (iMc < 0) return iPoint;
assert(iMc >= 0 && iMc < fpTofPoints->Size(iFile, iEvent));
int iIndex = iMc + fNpointsMvdAll + fNpointsStsAll + fNpointsMuchAll + fNpointsTrdAll;
+
auto itPoint = fmMCPointsLinksMap.find(CbmL1LinkKey(iIndex, iEvent, iFile));
if (itPoint == fmMCPointsLinksMap.cend()) { continue; }
iPoint = itPoint->second;
@@ -539,6 +539,7 @@ void CbmL1::ReadEvent(CbmEvent* event)
for (int iHit = 0; iHit < fpTofHits->GetEntriesFast(); iHit++) {
CbmMatch* pHitMatch = L1_DYNAMIC_CAST<CbmMatch*>(fpTofHitMatches->At(iHit));
+ assert(pHitMatch);
for (int iLink = 0; iLink < pHitMatch->GetNofLinks(); iLink++) {
Int_t iMC = pHitMatch->GetLink(iLink).GetIndex();
if (iMC < 0) continue;
@@ -1152,7 +1153,7 @@ void CbmL1::ReadEvent(CbmEvent* event)
fIODataManager.ResetInputData();
fIODataManager.ReserveNhits(nHits);
fIODataManager.SetNhitKeys(NStrips);
-
+ if (fPerformance) { fpMCModule->GetMCData()->ReserveNofHits(nHits); }
// ----- Fill
for (int iHit = 0; iHit < nHits; ++iHit) {
@@ -1195,7 +1196,10 @@ void CbmL1::ReadEvent(CbmEvent* event)
fvHitDebugInfo.push_back(s);
fvHitPointIndexes.push_back(th.iMC);
+ fpMCModule->GetMCData()->RegisterPointIndexForHit(iHit, th.iMC);
}
+ if (fPerformance) { HitMatch(); } /// OLD
+ if (fPerformance) { fpMCModule->MatchPointsWithHits(fvExternalHits); } /// NEW
if (fVerbose >= 2) cout << "ReadEvent: mvd and sts are saved." << endl;
@@ -1271,7 +1275,6 @@ void CbmL1::Fill_vMCTracks()
}
// TODO: Add light nuclei (d, t, He3, He4): they are common in tracking but not accounted in TDatabasePDG (S.Zharko)
- //cout << "mc track " << iMCTrack << " pdg " << pdg << " z " << vr.Z() << endl;
Int_t iTrack = fvMCTracks.size(); //or iMCTrack?
CbmL1MCTrack T(mass, q, vr, vp, iTrack, mother_ID, pdg, processID);
diff --git a/reco/L1/catools/CaToolsMCData.cxx b/reco/L1/catools/CaToolsMCData.cxx
index 865d661032..85c1fb0e27 100644
--- a/reco/L1/catools/CaToolsMCData.cxx
+++ b/reco/L1/catools/CaToolsMCData.cxx
@@ -24,8 +24,9 @@ MCData::MCData() {}
MCData::MCData(const MCData& other)
: fvPoints(other.fvPoints)
, fvTracks(other.fvTracks)
- , fmPointsLinksMap(other.fmPointsLinksMap)
- , fmTracksLinksMap(other.fmPointsLinksMap)
+ , fvPointIndexOfHit(other.fvPointIndexOfHit)
+ , fmPointLinkMap(other.fmPointLinkMap)
+ , fmTrackLinkMap(other.fmTrackLinkMap)
{
}
@@ -58,15 +59,16 @@ void MCData::Swap(MCData& other) noexcept
{
std::swap(fvPoints, other.fvPoints);
std::swap(fvTracks, other.fvTracks);
- std::swap(fmPointsLinksMap, other.fmPointsLinksMap);
- std::swap(fmTracksLinksMap, other.fmTracksLinksMap);
+ std::swap(fvPointIndexOfHit, other.fvPointIndexOfHit);
+ std::swap(fmPointLinkMap, other.fmPointLinkMap);
+ std::swap(fmTrackLinkMap, other.fmTrackLinkMap);
}
// ---------------------------------------------------------------------------------------------------------------------
//
void MCData::AddPoint(const MCPoint& point, int index, int event, int file)
{
- fmPointsLinksMap[LinkKey(index, event, file)] = static_cast<int>(fvPoints.size());
+ fmPointLinkMap[LinkKey(index, event, file)] = static_cast<int>(fvPoints.size());
fvPoints.push_back(point);
}
@@ -74,7 +76,7 @@ void MCData::AddPoint(const MCPoint& point, int index, int event, int file)
//
void MCData::AddTrack(const CbmL1MCTrack& track, int index, int event, int file)
{
- fmTracksLinksMap[LinkKey(index, event, file)] = static_cast<int>(fvTracks.size());
+ fmTrackLinkMap[LinkKey(index, event, file)] = static_cast<int>(fvTracks.size());
fvTracks.push_back(track);
}
@@ -84,8 +86,9 @@ void MCData::Clear()
{
fvPoints.clear();
fvTracks.clear();
- fmPointsLinksMap.clear();
- fmTracksLinksMap.clear();
+ fvPointIndexOfHit.clear();
+ fmPointLinkMap.clear();
+ fmTrackLinkMap.clear();
}
// ---------------------------------------------------------------------------------------------------------------------
@@ -94,7 +97,8 @@ std::string MCData::ToString(int verbose) const
{
if (verbose < 1) { return std::string(); }
std::stringstream msg;
- msg << "MCData: " << fvTracks.size() << " tracks, " << fvPoints.size() << " points";
+ msg << "MCData: " << fvTracks.size() << " tracks, " << fvPoints.size() << " points, ";
+ msg << fmTrackLinkMap.size() << " track links, " << fmPointLinkMap.size() << " point links";
if (verbose > 1) {
using std::setfill;
using std::setw;
diff --git a/reco/L1/catools/CaToolsMCData.h b/reco/L1/catools/CaToolsMCData.h
index fba046cbf6..8f4584864b 100644
--- a/reco/L1/catools/CaToolsMCData.h
+++ b/reco/L1/catools/CaToolsMCData.h
@@ -52,6 +52,9 @@ namespace ca::tools
/// Swap method
void Swap(MCData& other) noexcept;
+ /// Adds hit index to MC point
+ /// \param iPoint Index of MC point
+ /// \param
/// Adds an MC point
/// \param point MC point object
@@ -70,60 +73,70 @@ namespace ca::tools
/// Clears contents
void Clear();
- /// Reserves memory for tracks to avoid extra allocations
- void ReserveNofTracks(int nTracks) { fvTracks.reserve(nTracks); }
+ /// Finds an MC point global index by given local index, event and file of a link
+ /// \param index Local index of MC point
+ /// \param event Event of link
+ /// \param file File of link
+ int FindGlobalPointIndex(int index, int event, int file) const
+ {
+ auto it = fmPointLinkMap.find(LinkKey(index, event, file));
+ return (it != fmPointLinkMap.cend()) ? it->second : -1;
+ }
- /// Reserves memory for points to avoid extra allocations
- void ReserveNofPoints(int nPoints) { fvPoints.reserve(nPoints); }
+ /// Finds an MC track global index by given local index, event and file of a link
+ /// \param index Local index of track
+ /// \param event Index of event
+ /// \param file Index of file
+ int FindGlobalTrackIndex(int index, int event, int file) const
+ {
+ auto it = fmTrackLinkMap.find(LinkKey(index, event, file));
+ return (it != fmTrackLinkMap.cend()) ? it->second : -1;
+ }
+ /// Gets number of tracks
+ int GetNofTracks() const { return fvTracks.size(); }
- // **********************
- // ** Getters **
- // **********************
+ /// Gets number of points
+ int GetNofPoints() const { return fvPoints.size(); }
/// Gets a reference to MC point by its index
const auto& GetPoint(int idx) const { return fvPoints[idx]; }
+ /// Gets mutual reference to MC point by its index
+ // TODO: SZh 12.12.2022: Probably, the better solution is to write a specific accessor for
+ // setting indexes to MC points
+ auto& GetPoint(int idx) { return fvPoints[idx]; }
+
/// Gets a reference to the vector of points
const auto& GetPointContainer() const { return fvPoints; }
+ /// Gets point index by global index of hit
+ /// \param iHit Index of hit
+ auto GetPointIndexOfHit(int iHit) const { return fvPointIndexOfHit[iHit]; }
+
/// Gets a reference to MC track by its index
const auto& GetTrack(int idx) const { return fvTracks[idx]; }
/// Gets a reference to the vector of tracks
const auto& GetTrackContainer() const { return fvTracks; }
-
- /// Finds an MC point global index by given local index, event and file of a link
- /// \param index Local index of MC point
- /// \param event Event of link
- /// \param file File of link
- int FindGlobalPointIndex(int index, int event, int file)
- {
- auto it = fmPointsLinksMap.find(LinkKey(index, event, file));
- return (it != fmPointsLinksMap.cend()) ? it->second : -1;
- }
-
- /// Finds an MC track global index by given local index, event and file of a link
- /// \param index Local index of track
- /// \param event Index of event
- /// \param file Index of file
- int FindGlobalTrackIndex(int index, int event, int file)
+ /// Registers index of point for a given index of hit
+ /// \param iHit Index of hit
+ /// \param iPoint Index of point
+ void RegisterPointIndexForHit(int iHit, int iPoint)
{
- auto it = fmTracksLinksMap.find(LinkKey(index, event, file));
- return (it != fmTracksLinksMap.cend()) ? it->second : -1;
+ assert(int(fvPointIndexOfHit.size()) > iHit);
+ fvPointIndexOfHit[iHit] = iPoint;
}
- /// Gets number of tracks
- int GetNofTracks() const { return fvTracks.size(); }
-
- /// Gets number of points
- int GetNofPoints() { return fvPoints.size(); }
+ /// Reserves memory for tracks to avoid extra allocations
+ void ReserveNofTracks(int nTracks) { fvTracks.reserve(nTracks); }
+ /// Reserves memory for points to avoid extra allocations
+ void ReserveNofPoints(int nPoints) { fvPoints.reserve(nPoints); }
- // *********************
- // ** Setters **
- // *********************
+ /// Reserves total number of used hits in the event
+ void ReserveNofHits(int nHits) { fvPointIndexOfHit.reset(nHits); }
/// Prints an example of tracks and points
/// \param verbose Verbose level:
@@ -140,8 +153,11 @@ namespace ca::tools
L1Vector<MCPoint> fvPoints = {"ca::tools::MCData::fvMCPoints"}; ///< Container of points
L1Vector<CbmL1MCTrack> fvTracks = {"ca::tools::MCData::fvMCTracks"}; ///< Container of tracks
- std::unordered_map<LinkKey, int> fmPointsLinksMap = {}; ///< MC point internal index vs. link
- std::unordered_map<LinkKey, int> fmTracksLinksMap = {}; ///< MC track internal index vs. link
+ /// Correspondence of MC point index to the global hit index
+ L1Vector<int> fvPointIndexOfHit = {"ca::tools::MCData::fvPointIndexOfHit"};
+
+ std::unordered_map<LinkKey, int> fmPointLinkMap = {}; ///< MC point internal index vs. link
+ std::unordered_map<LinkKey, int> fmTrackLinkMap = {}; ///< MC track internal index vs. link
};
} // namespace ca::tools
diff --git a/reco/L1/catools/CaToolsMCPoint.h b/reco/L1/catools/CaToolsMCPoint.h
index dedc3ddff2..465b8569c6 100644
--- a/reco/L1/catools/CaToolsMCPoint.h
+++ b/reco/L1/catools/CaToolsMCPoint.h
@@ -20,8 +20,6 @@ enum class L1DetectorID;
/// Class describes a Monte-Carlo point used in CA tracking QA analysis
namespace ca::tools
{
- class MCTrack;
-
class MCPoint {
public:
/// Default constructor
--
GitLab