Make exp filter alpha a parameter

Author: fweig

GPU/GPUTracking/Definitions/GPUSettingsList.h

@@ -114,6 +114,7 @@ AddOptionRTC(maxConsecTimeBinAboveThreshold, uint16_t, 200, "", 0, "Except pad f
 AddOptionRTC(noisyPadSaturationThreshold, uint16_t, 700, "", 0, "Threshold where a timebin is considered saturated, disabling the noisy pad check for that pad")
 AddOptionRTC(hipTailFilter, uint8_t, 0, "", 0, "Enable Highly Ionising Particle tail filter in CheckPadBaseline (0 = disable, 1 = filter tails)")
 AddOptionRTC(hipTailFilterThreshold, uint16_t, 150, "", 0, "Threshold that must be exceeded for a timebin to be counted towards Highly Ionising Particle tail")
+AddOptionRTC(hipTailFilterAlpha, float, 0.25f, "", 0, "Smoothing factor for the exponential Highly Ionising Particle tail filter")
 AddOptionRTC(occupancyMapTimeBins, uint16_t, 16, "", 0, "Number of timebins per histogram bin of occupancy map (0 = disable occupancy map)")
 AddOptionRTC(occupancyMapTimeBinsAverage, uint16_t, 0, "", 0, "Number of timebins +/- to use for the averaging")
 AddOptionRTC(trackFitCovLimit, uint16_t, 1000, "", 0, "Abort fit when y/z cov exceed the limit")

GPU/GPUTracking/TPCClusterFinder/GPUTPCCFCheckPadBaseline.cxx

@@ -41,6 +41,11 @@ using namespace o2::gpu::tpccf;
 
 using Kernel = GPUTPCCFCheckPadBaseline;
 
+static GPUdi() Charge UpdateHIPTailFilter(Charge filteredCharge, Charge charge, Charge alpha)
+{
+  return filteredCharge + alpha * (charge - filteredCharge);
+}
+
 // Collect tails marked for closing across the workgroup using a prefix scan,
 // then cooperatively zero the charge map entries for each closed tail.
 // Caller must set acc.activeHIPTail.end before calling if the tail is open.
@@ -74,6 +79,7 @@ static GPUdi() uint16_t CloseHIPTails(
       }
 
       acc.tailQTot = 0;
+      acc.tailFilterCharge = 0;
       acc.activeHIPTail.Reset();
     }
 
@@ -94,7 +100,7 @@ static GPUdi() uint16_t CloseHIPTails(
 }
 
 template <bool CheckHIPTrigger, bool CheckHIPTailEnd>
-static GPUdi() void ScanCachedCharges(Kernel::GPUSharedMemory& smem, uint16_t timeOffset, uint16_t pad, Charge hipTailThreshold, Kernel::PadChargeAccu& acc)
+static GPUdi() void ScanCachedCharges(Kernel::GPUSharedMemory& smem, uint16_t timeOffset, uint16_t pad, Charge hipTailThreshold, Charge hipTailFilterAlpha, Kernel::PadChargeAccu& acc)
 {
   for (int32_t i = 0; i < Kernel::NumOfCachedTBs; i++) {
     const Charge qs = smem.charges[i][pad];
@@ -123,7 +129,8 @@ static GPUdi() void ScanCachedCharges(Kernel::GPUSharedMemory& smem, uint16_t ti
     if constexpr (CheckHIPTailEnd) {
       if (acc.activeHIPTail.IsOpen()) {
         acc.tailQTot += qs;
-        if (qs < hipTailThreshold) {
+        acc.tailFilterCharge = UpdateHIPTailFilter(acc.tailFilterCharge, qs, hipTailFilterAlpha);
+        if (acc.tailFilterCharge < hipTailThreshold) {
           acc.activeHIPTail.end = curTB;
         }
       }
@@ -157,6 +164,7 @@ GPUd() void GPUTPCCFCheckPadBaseline::CheckBaselineGPU(int32_t nBlocks, int32_t
   const CfFragment& fragment = clusterer.mPmemory->fragment;
   const bool hipFilterOn = clusterer.Param().rec.tpc.hipTailFilter;
   const Charge hipTailThreshold = clusterer.Param().rec.tpc.hipTailFilterThreshold;
+  const Charge hipTailFilterAlpha = clusterer.Param().rec.tpc.hipTailFilterAlpha;
   CfArray2D<PackedCharge> chargeMap(reinterpret_cast<PackedCharge*>(clusterer.mPchargeMap));
 
   const auto iRow = iBlock;
@@ -205,15 +213,15 @@ GPUd() void GPUTPCCFCheckPadBaseline::CheckBaselineGPU(int32_t nBlocks, int32_t
       // Why is the old version so much slower, when we just add short branches to the loop???
       if (!hasHIPTrigger) [[likely]] {
         if (!acc.activeHIPTail.IsOpen()) {
-          ScanCachedCharges<false, false>(smem, t, iPadHandle, hipTailThreshold, acc);
+          ScanCachedCharges<false, false>(smem, t, iPadHandle, hipTailThreshold, hipTailFilterAlpha, acc);
         } else {
-          ScanCachedCharges<false, true>(smem, t, iPadHandle, hipTailThreshold, acc);
+          ScanCachedCharges<false, true>(smem, t, iPadHandle, hipTailThreshold, hipTailFilterAlpha, acc);
         }
       } else {
         if (!acc.activeHIPTail.IsOpen()) {
-          ScanCachedCharges<true, false>(smem, t, iPadHandle, hipTailThreshold, acc);
+          ScanCachedCharges<true, false>(smem, t, iPadHandle, hipTailThreshold, hipTailFilterAlpha, acc);
         } else {
-          ScanCachedCharges<true, true>(smem, t, iPadHandle, hipTailThreshold, acc);
+          ScanCachedCharges<true, true>(smem, t, iPadHandle, hipTailThreshold, hipTailFilterAlpha, acc);
         }
       }
     }
@@ -255,6 +263,7 @@ GPUd() void GPUTPCCFCheckPadBaseline::CheckBaselineGPU(int32_t nBlocks, int32_t
       if (acc.HIPtb > -1) {
         DPRINT("%d: start = %d\n", iThread, acc.HIPtb);
         acc.activeHIPTail.SetOpen(acc.HIPtb);
+        acc.tailFilterCharge = Charge(MaxADC);
       }
 
       // Clear smem between iterations to prevent stale entries

GPU/GPUTracking/TPCClusterFinder/GPUTPCCFCheckPadBaseline.h

@@ -18,8 +18,8 @@
 ///
 /// Optionally detects Highly Ionising Particle (HIP) tails: when a saturated
 /// ADC value (1023) is found, the tail region on the triggering pad and its
-/// neighbors is zeroed in the charge map until charges drop below a
-/// configurable threshold.
+/// neighbors is zeroed in the charge map until an exponential charge filter
+/// drops below a configurable threshold.
 
 #ifndef O2_GPU_GPU_TPC_CF_CHECK_PAD_BASELINE_H
 #define O2_GPU_GPU_TPC_CF_CHECK_PAD_BASELINE_H
@@ -107,6 +107,7 @@ class GPUTPCCFCheckPadBaseline : public GPUKernelTemplate
     int16_t HIPtb = -1;
     int16_t aboveThresholdStart = -1; // first TB of current above-hipTailThreshold streak; used to extend the tail back over the rising edge before saturation
     HipTailRange activeHIPTail{-1, -1};
+    tpccf::Charge tailFilterCharge = 0;
     float tailQTot = 0;
   };