From 915af2f750e35c6a2fbc44b9af7016a4925f3003 Mon Sep 17 00:00:00 2001
From: Dominik Smith <smith@th.physik.uni-frankfurt.de>
Date: Thu, 8 Dec 2022 17:55:35 +0100
Subject: [PATCH] Refactored HitFinderTof. Two micro-optimizations of
 CbmTaskTofHitFinder which save about 10 percent execution time in total.

---
 algo/detectors/tof/HitFinderTof.cxx | 114 +++++++++++++++-------------
 algo/detectors/tof/HitFinderTof.h   |  13 ++++
 reco/tasks/CbmTaskTofHitFinder.cxx  |   8 +-
 3 files changed, 79 insertions(+), 56 deletions(-)

diff --git a/algo/detectors/tof/HitFinderTof.cxx b/algo/detectors/tof/HitFinderTof.cxx
index 630cfae883..8ff2770ffc 100644
--- a/algo/detectors/tof/HitFinderTof.cxx
+++ b/algo/detectors/tof/HitFinderTof.cxx
@@ -16,57 +16,15 @@ namespace cbm::algo
   HitFinderTof::resultType HitFinderTof::operator()(std::vector<CbmTofDigi> digisIn,
                                                     const std::vector<int32_t>& digiIndexIn)
   {
-    // Intermediate storage variables
-    std::vector<std::vector<CbmTofDigi*>> digisExp;  //[nbCh][nDigis]
-    std::vector<std::vector<int32_t>> digisInd;      //[nbCh][nDigis]
-
-    digisExp.resize(fParams.fChanPar.size());
-    digisInd.resize(fParams.fChanPar.size());
-
-    {  // digi calibration
-      const int32_t numClWalkBinX = fParams.numClWalkBinX;
-      const double TOTMax         = fParams.TOTMax;
-      const double TOTMin         = fParams.TOTMin;
-
-      for (size_t iDigi = 0; iDigi < digisIn.size(); iDigi++) {
-        CbmTofDigi* pDigi = &digisIn[iDigi];
-        assert(pDigi);
-
-        // These are doubles in the digi class
-        const double chan   = pDigi->GetChannel();
-        const double side   = pDigi->GetSide();
-        const double charge = pDigi->GetTot() * fParams.fChanPar[chan].fvCPTotGain[side];  // calibrate Digi ToT
-        const double time   = pDigi->GetTime() - fParams.fChanPar[chan].fvCPTOff[side];    // calibrate Digi Time
-
-        digisExp[chan].push_back(pDigi);
-        digisInd[chan].push_back(digiIndexIn[iDigi]);
-
-        pDigi->SetTot(charge);
-
-        {  // walk correction
-          const double totBinSize = (TOTMax - TOTMin) / 2. / numClWalkBinX;
-          int32_t iWx             = (int32_t)((charge - TOTMin / 2.) / totBinSize);
-          if (0 > iWx) iWx = 0;
-          if (iWx > (numClWalkBinX - 1)) iWx = numClWalkBinX - 1;
-
-          std::vector<double>& cpWalk = fParams.fChanPar[chan].fvCPWalk[side];
-          double dWT                  = cpWalk[iWx];
-          const double dDTot          = (charge - TOTMin / 2.) / totBinSize - (double) iWx - 0.5;
+    inputType input = calibrateDigis(digisIn, digiIndexIn);
+    return buildClusters(input);
+  }
 
-          if (dDTot > 0) {
-            if (iWx < numClWalkBinX - 1) {  // linear interpolation to next bin
-              dWT += dDTot * (cpWalk[iWx + 1] - cpWalk[iWx]);
-            }
-          }
-          else {
-            if (0 < iWx) {  // linear interpolation to next bin
-              dWT -= dDTot * (cpWalk[iWx - 1] - cpWalk[iWx]);
-            }
-          }
-          pDigi->SetTime(time - dWT);
-        }
-      }
-    }  // digi calibration
+  HitFinderTof::resultType HitFinderTof::buildClusters(inputType& input)
+  {
+    // Intermediate storage variables
+    std::vector<std::vector<CbmTofDigi*>>& digisExp = input.first;   //[nbCh][nDigis]
+    std::vector<std::vector<int32_t>>& digisInd     = input.second;  //[nbCh][nDigis]
 
     // Hit variables
     TofCluster cluster;
@@ -173,7 +131,61 @@ namespace cbm::algo
       cluster.finalize(*trafoCell, iTrafoCell, fParams);
       clustersOut.push_back(cluster);
     }
-
     return clustersOut;
   }
+
+  HitFinderTof::inputType HitFinderTof::calibrateDigis(std::vector<CbmTofDigi>& digisIn,
+                                                       const std::vector<int32_t>& digiIndexIn)
+  {
+    inputType result;
+    std::vector<std::vector<CbmTofDigi*>>& digisExp = result.first;   //[nbCh][nDigis]
+    std::vector<std::vector<int32_t>>& digisInd     = result.second;  //[nbCh][nDigis]
+
+    digisExp.resize(fParams.fChanPar.size());
+    digisInd.resize(fParams.fChanPar.size());
+
+    const int32_t numClWalkBinX = fParams.numClWalkBinX;
+    const double TOTMax         = fParams.TOTMax;
+    const double TOTMin         = fParams.TOTMin;
+
+    for (size_t iDigi = 0; iDigi < digisIn.size(); iDigi++) {
+      CbmTofDigi* pDigi = &digisIn[iDigi];
+      assert(pDigi);
+
+      // These are doubles in the digi class
+      const double chan   = pDigi->GetChannel();
+      const double side   = pDigi->GetSide();
+      const double charge = pDigi->GetTot() * fParams.fChanPar[chan].fvCPTotGain[side];  // calibrate Digi ToT
+      const double time   = pDigi->GetTime() - fParams.fChanPar[chan].fvCPTOff[side];    // calibrate Digi Time
+
+      digisExp[chan].push_back(pDigi);
+      digisInd[chan].push_back(digiIndexIn[iDigi]);
+
+      pDigi->SetTot(charge);
+
+      {  // walk correction
+        const double totBinSize = (TOTMax - TOTMin) / 2. / numClWalkBinX;
+        int32_t iWx             = (int32_t)((charge - TOTMin / 2.) / totBinSize);
+        if (0 > iWx) iWx = 0;
+        if (iWx > (numClWalkBinX - 1)) iWx = numClWalkBinX - 1;
+
+        std::vector<double>& cpWalk = fParams.fChanPar[chan].fvCPWalk[side];
+        double dWT                  = cpWalk[iWx];
+        const double dDTot          = (charge - TOTMin / 2.) / totBinSize - (double) iWx - 0.5;
+
+        if (dDTot > 0) {
+          if (iWx < numClWalkBinX - 1) {  // linear interpolation to next bin
+            dWT += dDTot * (cpWalk[iWx + 1] - cpWalk[iWx]);
+          }
+        }
+        else {
+          if (0 < iWx) {  // linear interpolation to next bin
+            dWT -= dDTot * (cpWalk[iWx - 1] - cpWalk[iWx]);
+          }
+        }
+        pDigi->SetTime(time - dWT);
+      }
+    }
+    return result;
+  }
 } /* namespace cbm::algo */
diff --git a/algo/detectors/tof/HitFinderTof.h b/algo/detectors/tof/HitFinderTof.h
index cb21482196..5ad658e3db 100644
--- a/algo/detectors/tof/HitFinderTof.h
+++ b/algo/detectors/tof/HitFinderTof.h
@@ -2,6 +2,15 @@
    SPDX-License-Identifier: GPL-3.0-only
    Authors: Dominik Smith [committer], Pierre-Alain Loizeau */
 
+/*
+   This algo was based on CbmTofSimpClusterizer, which can be used only for simulation of the main setup and
+   is as the name implies a simplified solution. 
+   A later step will be the replacement with a version based on CbmTofEventClusterizer, which is the version 
+   currently maintained, based on what we learned from real data at mCBM. 
+   This step will be required to apply the algo to real/online data and to prepare 
+   our simulations for first CBM beam
+*/
+
 #ifndef HITFINDERTOF_H
 #define HITFINDERTOF_H
 
@@ -132,6 +141,7 @@ namespace cbm::algo
   class HitFinderTof {
   public:
     typedef std::vector<TofCluster> resultType;
+    typedef std::pair<std::vector<std::vector<CbmTofDigi*>>, std::vector<std::vector<int32_t>>> inputType;
 
     /**
        ** @brief Constructor.
@@ -156,6 +166,9 @@ namespace cbm::algo
   private:
     HitFinderTofRpcPar fParams = {};  ///< Parameter container
 
+    inputType calibrateDigis(std::vector<CbmTofDigi>& digisIn, const std::vector<int32_t>& digiIndexIn);
+    resultType buildClusters(inputType& input);
+
     int32_t numSameSide;  // Digis quality
   };
 
diff --git a/reco/tasks/CbmTaskTofHitFinder.cxx b/reco/tasks/CbmTaskTofHitFinder.cxx
index a49f1d2329..b6b5e6d5bb 100644
--- a/reco/tasks/CbmTaskTofHitFinder.cxx
+++ b/reco/tasks/CbmTaskTofHitFinder.cxx
@@ -397,8 +397,7 @@ pair<int32_t, int32_t> CbmTaskTofHitFinder::BuildClusters(CbmEvent* event)
   // Loop over the digis array and store the Digis in separate vectors for each RPC modules
   for (int32_t iDigi = 0; iDigi < nDigis; iDigi++) {
     const uint32_t digiIndex = (event ? event->GetIndex(ECbmDataType::kTofDigi, iDigi) : iDigi);
-    assert(fDigiMan->Get<CbmTofDigi>(digiIndex));
-    CbmTofDigi* pDigi = new CbmTofDigi(*(fDigiMan->Get<CbmTofDigi>(digiIndex)));
+    const CbmTofDigi* pDigi  = fDigiMan->Get<CbmTofDigi>(digiIndex);
     assert(pDigi);
 
     // These are doubles in the digi class
@@ -434,12 +433,11 @@ pair<int32_t, int32_t> CbmTaskTofHitFinder::BuildClusters(CbmEvent* event)
             CbmTofHit(cluster.detId, hitpos, hiterr, nHits, cluster.weightedTime, cluster.weightedTimeErr, 0, 0);
           nHits++;
           if (event) event->AddData(ECbmDataType::kTofHit, hitIndex);
-          CbmMatch* digiMatch = new CbmMatch();
+
+          CbmMatch* digiMatch = new ((*fTofDigiMatchColl)[hitIndex]) CbmMatch();
           for (uint32_t i = 0; i < cluster.vDigiIndRef.size(); i++) {
             digiMatch->AddLink(CbmLink(0., cluster.vDigiIndRef.at(i), iEventNr, iInputNr));
           }
-          new ((*fTofDigiMatchColl)[hitIndex]) CbmMatch(*digiMatch);
-          delete digiMatch;
         }
         digiExp.clear();
         digiInd.clear();
-- 
GitLab