From f089221b2b8c6a98ce199c8ad53f6c0ceb3eb5b1 Mon Sep 17 00:00:00 2001
From: "se.gorbunov" <se.gorbunov@gsi.de>
Date: Wed, 1 Mar 2023 00:31:23 +0000
Subject: [PATCH] L1: new subslice splitting algorithm that searches in
partially disordered data
---
reco/L1/CbmL1ReadEvent.cxx | 9 ----
reco/L1/L1Algo/L1Algo.h | 4 +-
reco/L1/L1Algo/L1CaTrackFinder.cxx | 79 ++++++++++++++++++++----------
3 files changed, 57 insertions(+), 35 deletions(-)
diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx
index 7c8762b088..a4173fb6c7 100644
--- a/reco/L1/CbmL1ReadEvent.cxx
+++ b/reco/L1/CbmL1ReadEvent.cxx
@@ -1178,15 +1178,6 @@ void CbmL1::ReadEvent(CbmEvent* event)
* Measured hits gathering: END
*/
-
- //
- // Hits sorting
- //
- // Two hits are compared as follows. If the hits are measured with two different stations, the smallest hit has the smallest
- // station ID. If the hits are measured within one station, the smallest hit has the smallest y position coordinate.
- //
- sort(tmpHits.begin(), tmpHits.end(), TmpHit::Compare);
-
if (fVerbose >= 10) { cout << "ReadEvent: strips are read." << endl; }
diff --git a/reco/L1/L1Algo/L1Algo.h b/reco/L1/L1Algo/L1Algo.h
index db73196b6c..71f0568377 100644
--- a/reco/L1/L1Algo/L1Algo.h
+++ b/reco/L1/L1Algo/L1Algo.h
@@ -75,6 +75,8 @@ class L1AlgoEfficiencyPerformance;
struct L1HitTimeInfo {
fscal fEventTimeMin {-std::numeric_limits<fscal>::max() / 2.};
fscal fEventTimeMax {std::numeric_limits<fscal>::max() / 2.};
+ fscal fMaxTimeBeforeHit {0.}; //< max event time for hits [0 .. hit] in the station hit array
+ fscal fMinTimeAfterHit {0.}; //< min event time for hits [hit ... ] in the station hit array
};
/// Main class of L1 CA track finder algorithm
@@ -463,7 +465,7 @@ public:
L1Vector<omp_lock_t> fStripToTrackLock {"L1Algo::fStripToTrackLock"};
#endif
- L1Vector<int> fStripToTrack {"L1Algo::fStripToTrack"}; // front strip to track pointers
+ L1Vector<int> fStripToTrack {"L1Algo::fStripToTrack"}; // front strip to track pointers
// L1Vector<int> fStripToTrack1B {"L1Algo::fStripToTrackB"}; // back strip to track pointers
int fNThreads {0};
diff --git a/reco/L1/L1Algo/L1CaTrackFinder.cxx b/reco/L1/L1Algo/L1CaTrackFinder.cxx
index e5a33fb6e8..ba9d808422 100644
--- a/reco/L1/L1Algo/L1CaTrackFinder.cxx
+++ b/reco/L1/L1Algo/L1CaTrackFinder.cxx
@@ -63,34 +63,48 @@ void L1Algo::CaTrackFinder()
for (int iS = 0; iS < fParameters.GetNstationsActive(); ++iS) {
const L1Station& st = fParameters.GetStation(iS);
- L1HitIndex_t nStationHits = fInputData.GetNhits(iS);
- for (L1HitIndex_t ih = 0; ih < nStationHits; ++ih) {
- int caHitId = fInputData.GetStartHitIndex(iS) + ih;
- const L1Hit& h = fInputData.GetHit(caHitId);
+ int nStationHits = fInputData.GetNhits(iS);
+
+ fscal maxTimeBeforeHit = std::numeric_limits<fscal>::min();
+
+ // loop in the reverse order to fill L1HitTimeInfo::fMinTimeAfterHit fields
+
+ for (int ih = 0; ih < nStationHits; ++ih) {
+ L1HitIndex_t caHitId = fInputData.GetStartHitIndex(iS) + ih;
+ const L1Hit& h = fInputData.GetHit(caHitId);
+
+ auto [x, y] = st.ConvUVtoXY(h.u, h.v);
+ fscal dx = x - targX;
+ fscal dy = y - targY;
+ fscal dz = h.z - targZ;
+ fscal l = sqrt(dx * dx + dy * dy + dz * dz);
+
+ fscal timeOfFlightMin = l * L1TrackPar::kClightNsInv;
+ fscal timeOfFlightMax =
+ 1.5 * l * sqrt(1. + L1TrackPar::kProtonMass * L1TrackPar::kProtonMass / minProtonMomentum / minProtonMomentum)
+ * L1TrackPar::kClightNsInv;
+ fscal dt = 3.5 * sqrt(h.dt2);
+
L1HitTimeInfo& info = fHitTimeInfo[caHitId];
- if (0 && !st.timeInfo) {
- info.fEventTimeMin = -std::numeric_limits<fscal>::max() / 2.;
- info.fEventTimeMax = std::numeric_limits<fscal>::max() / 2.;
- }
- else {
- auto [x, y] = st.ConvUVtoXY(h.u, h.v);
- fscal dx = x - targX;
- fscal dy = y - targY;
- fscal dz = h.z - targZ;
- fscal l = sqrt(dx * dx + dy * dy + dz * dz);
-
- fscal timeOfFlightMin = l * L1TrackPar::kClightNsInv;
- fscal timeOfFlightMax =
- 1.5 * l * sqrt(1. + L1TrackPar::kProtonMass * L1TrackPar::kProtonMass / minProtonMomentum / minProtonMomentum)
- * L1TrackPar::kClightNsInv;
- fscal dt = 3.5 * sqrt(h.dt2);
- info.fEventTimeMin = h.t - dt - timeOfFlightMax;
- info.fEventTimeMax = h.t + dt - timeOfFlightMin;
- }
+ info.fEventTimeMin = h.t - dt - timeOfFlightMax;
+ info.fEventTimeMax = h.t + dt - timeOfFlightMin;
+
+ if (maxTimeBeforeHit < info.fEventTimeMax) { maxTimeBeforeHit = info.fEventTimeMax; }
+ info.fMaxTimeBeforeHit = maxTimeBeforeHit;
+
if (tsStart > info.fEventTimeMax - 5) {
tsStart = info.fEventTimeMax - 5; // 5 ns margin
}
}
+
+ fscal minTimeAfterHit = std::numeric_limits<fscal>::max();
+
+ for (int ih = nStationHits - 1; ih >= 0; --ih) {
+ L1HitIndex_t caHitId = fInputData.GetStartHitIndex(iS) + ih;
+ L1HitTimeInfo& info = fHitTimeInfo[caHitId];
+ if (minTimeAfterHit > info.fEventTimeMin) { minTimeAfterHit = info.fEventTimeMin; }
+ info.fMinTimeAfterHit = minTimeAfterHit;
+ }
}
// cut data into sub-timeslices and process them one by one
@@ -98,6 +112,12 @@ void L1Algo::CaTrackFinder()
bool areDataLeft = true; // is the whole TS processed
int nSubSlices = 0;
+ L1HitIndex_t sliceFirstHit[L1Constants::size::kMaxNstations] {0};
+
+ for (int iS = 0; iS < fParameters.GetNstationsActive(); ++iS) {
+ sliceFirstHit[iS] = fInputData.GetStartHitIndex(iS);
+ }
+
while (areDataLeft) {
nSubSlices++;
@@ -111,8 +131,8 @@ void L1Algo::CaTrackFinder()
for (int iS = 0; iS < fParameters.GetNstationsActive(); ++iS) {
fSliceHitIdsStartIndex[iS] = fSliceHitIds.size();
- for (L1HitIndex_t caHitId = fInputData.GetStartHitIndex(iS); caHitId < fInputData.GetStopHitIndex(iS);
- ++caHitId) {
+
+ for (L1HitIndex_t caHitId = sliceFirstHit[iS]; caHitId < fInputData.GetStopHitIndex(iS); ++caHitId) {
L1HitTimeInfo& info = fHitTimeInfo[caHitId];
const L1Hit& h = fInputData.GetHit(caHitId);
if (fvHitKeyFlags[h.f] || fvHitKeyFlags[h.b]) { // the hit is already reconstructed
@@ -120,13 +140,22 @@ void L1Algo::CaTrackFinder()
}
if (info.fEventTimeMin > tsStart + tsLength) { // the hit is too late for the sub slice
areDataLeft = true;
+ if (info.fMinTimeAfterHit > tsStart + tsLength) {
+ // this hit and all later hits are out of the sub-slice
+ break;
+ }
}
else {
if (tsStart <= info.fEventTimeMax) { // the hit belongs to the sub-slice
fSliceHitIds.push_back(caHitId);
+ if (info.fMaxTimeBeforeHit < tsStart + tsLength - tsOverlap) {
+ // this hit and all hits before are before the overlap
+ sliceFirstHit[iS] = caHitId + 1;
+ }
}
} // else the hit has been alread processed in previous sub-slices
}
+
fSliceHitIdsStopIndex[iS] = fSliceHitIds.size();
}
--
GitLab