GPU TPC: Make workarounds for cyclic merge graphs optional (to be checked, but should be removed eventually)

Author: davidrohr

GPU/GPUTracking/Definitions/GPUSettingsList.h

@@ -199,6 +199,7 @@ AddOptionRTC(fitInProjections, int8_t, -1, "", 0, "Fit in projection, -1 to enab
 AddOptionRTC(fitPropagateBzOnly, int8_t, -1, "", 0, "Propagate using Bz only for n passes")
 AddOptionRTC(useMatLUT, int8_t, 0, "", 0, "Use material lookup table for TPC refit")
 AddOptionRTC(trackingRefitGPUModel, int8_t, 1, "", 0, "Use GPU track model for the Global Track Refit")
+AddOptionRTC(enableCyclicGraphWorkarounds, int8_t, 1, "", 0, "Apply workarounds to avoid cyclic merge graphs, should not be needed")
 AddCustomCPP(void SetMinTrackPtB5(float v) { maxTrackQPtB5 = v > 0.001f ? (1.f / v) : (1.f / 0.001f); })
 AddSubConfig(GPUSettingsRecTPC, tpc)
 AddSubConfig(GPUSettingsRecTRD, trd)

GPU/GPUTracking/Merger/GPUTPCGMMergedTrack.h

@@ -50,23 +50,23 @@ class GPUTPCGMMergedTrack
   GPUd() bool MergedLooper() const { return mFlags & 0x30; }
   GPUd() int32_t PrevSegment() const { return mPrevSegment; }
   template <class T>
-  GPUd() static T* GetFirstSegment_helper(T* me, T* base)
+  GPUd() static T* GetFirstSegment_helper(T* me, T* base, bool workaround)
   {
     if (me->mPrevSegment < 0) {
       return me;
     }
     T* cur = &base[me->mPrevSegment];
     while (cur->mPrevSegment >= 0) {
       T* next = &base[cur->mPrevSegment];
-      if (next == me) {
+      if (workaround && next == me) {
         return cur;
       }
       cur = next;
     }
     return cur;
   }
-  GPUd() GPUTPCGMMergedTrack* GetFirstSegment(GPUTPCGMMergedTrack* base) { return GetFirstSegment_helper<GPUTPCGMMergedTrack>(this, base); }
-  GPUd() const GPUTPCGMMergedTrack* GetFirstSegment(const GPUTPCGMMergedTrack* base) const { return GetFirstSegment_helper<const GPUTPCGMMergedTrack>(this, base); }
+  GPUd() GPUTPCGMMergedTrack* GetFirstSegment(GPUTPCGMMergedTrack* base, bool workaround) { return GetFirstSegment_helper<GPUTPCGMMergedTrack>(this, base, workaround); }
+  GPUd() const GPUTPCGMMergedTrack* GetFirstSegment(const GPUTPCGMMergedTrack* base, bool workaround) const { return GetFirstSegment_helper<const GPUTPCGMMergedTrack>(this, base, workaround); }
   GPUd() uint8_t Leg() const { return mLeg; }
   GPUd() uint8_t Flags() const { return mFlags; }
 

GPU/GPUTracking/Merger/GPUTPCGMMerger.cxx

@@ -1489,7 +1489,7 @@ GPUd() void GPUTPCGMMerger::CollectMergedTracks(int32_t nBlocks, int32_t nThread
       if (jtr >= 0) {
         int32_t lasttr = itr;
         while (jtr >= 0) { // --------------- count segments ---------------
-          if (&mSectorTrackInfos[jtr] == trbase) {
+          if (Param().rec.enableCyclicGraphWorkarounds && &mSectorTrackInfos[jtr] == trbase) {
             break; // Break cyclic graph
           }
           lasttr = jtr;
@@ -1512,7 +1512,7 @@ GPUd() void GPUTPCGMMerger::CollectMergedTracks(int32_t nBlocks, int32_t nThread
               mainT = t;
             }
             int32_t next = trchk->NextSegmentNeighbour();
-            if (next < 0 || next == ichk) {
+            if (next < 0 || (Param().rec.enableCyclicGraphWorkarounds && next == ichk)) {
               break; // Breaks also cycles
             }
             trchk = &mSectorTrackInfos[next];
@@ -1533,7 +1533,7 @@ GPUd() void GPUTPCGMMerger::CollectMergedTracks(int32_t nBlocks, int32_t nThread
               length = trchk->OrigTrack()->NHits();
             }
             int32_t next = trchk->NextSegmentNeighbour();
-            if (next < 0 || next == ichk) {
+            if (next < 0 || (Param().rec.enableCyclicGraphWorkarounds && next == ichk)) {
               break; // Breaks also cycles
             }
             trchk = &mSectorTrackInfos[next];
@@ -1797,7 +1797,7 @@ GPUd() void GPUTPCGMMerger::PrepareForFit1(int32_t nBlocks, int32_t nThreads, in
         CAMath::AtomicAdd(&mSharedCount[mClusters[trk.FirstClusterRef() + j].num], 1u);
       }
       if (!trk.CCE() && !trk.MergedLooper()) {
-        GPUTPCGMMergedTrack* updTrk = trk.GetFirstSegment(mMergedTracks);
+        GPUTPCGMMergedTrack* updTrk = trk.GetFirstSegment(mMergedTracks, Param().rec.enableCyclicGraphWorkarounds);
         const auto &cl0 = mClusters[trk.FirstClusterRef()], &cln = mClusters[updTrk->FirstClusterRef() + updTrk->NClusters() - 1];
         const auto& GPUrestrict() cls = GetConstantMem()->ioPtrs.clustersNative->clustersLinear;
         float z0 = cls[cl0.num].getTime(), zn = cls[cln.num].getTime();
@@ -1806,7 +1806,7 @@ GPUd() void GPUTPCGMMerger::PrepareForFit1(int32_t nBlocks, int32_t nThreads, in
         updTrk = &trk;
         while (updTrk->PrevSegment() >= 0) {
           auto next = &mMergedTracks[updTrk->PrevSegment()];
-          if (next == &trk) {
+          if (Param().rec.enableCyclicGraphWorkarounds && next == &trk) {
             break;
           }
           updTrk = next;
@@ -1966,7 +1966,7 @@ GPUd() void GPUTPCGMMerger::MergeLoopersMain(int32_t nBlocks, int32_t nThreads,
       const GPUTPCGMMergedTrack* trkI = &mMergedTracks[candidates[i].id];
       float refZI = candidates[i].refz;
       {
-        const auto* tmp = trkI->GetFirstSegment(mMergedTracks);
+        const auto* tmp = trkI->GetFirstSegment(mMergedTracks, Param().rec.enableCyclicGraphWorkarounds);
         if (tmp != trkI && tmp->CSide() == trkI->CSide() && CAMath::Abs(tmp->GetParam().GetZ()) > CAMath::Abs(trkI->GetParam().GetZ())) {
           float tmpRefZ = refZI + tmp->GetParam().GetZ() - trkI->GetParam().GetZ();
           if (CAMath::Abs(tmpRefZ) < CAMath::Abs(candidates[j].refz) && CAMath::Abs(tmpRefZ) > CAMath::Abs(refZI)) {

GPU/GPUTracking/Merger/GPUTPCGMO2Output.cxx

@@ -212,7 +212,7 @@ GPUdii() void GPUTPCGMO2Output::Thread<GPUTPCGMO2Output::output>(int32_t nBlocks
     }
 
     if (track.PrevSegment() >= 0) {
-      const GPUTPCGMMergedTrack* chkTrk = track.GetFirstSegment(tracks);
+      const GPUTPCGMMergedTrack* chkTrk = track.GetFirstSegment(tracks, merger.Param().rec.enableCyclicGraphWorkarounds);
       const auto& firstPrevCluster = trackClusters[chkTrk->FirstClusterRef()];
       t1 = clusters->clustersLinear[firstPrevCluster.num].getTime();
       sector1 = firstPrevCluster.sector;