From eed1fae40594d9d6dedea0cdf7554b48804a94bb Mon Sep 17 00:00:00 2001
From: "se.gorbunov" <se.gorbunov@gsi.de>
Date: Mon, 11 Sep 2023 17:28:39 +0000
Subject: [PATCH] CA: fix the numeration of Much tracking stations

---
 core/base/CbmTrackingDetectorInterfaceBase.h  |  8 +--
 .../much/CbmMuchTrackingInterface.cxx         | 15 ++++-
 .../detectors/much/CbmMuchTrackingInterface.h | 47 +++++----------
 core/detectors/mvd/CbmMvdTrackingInterface.h  |  8 +--
 core/detectors/sts/CbmStsTrackingInterface.h  |  8 +--
 core/detectors/tof/CbmTofTrackingInterface.h  |  4 +-
 core/detectors/trd/CbmTrdTrackingInterface.h  |  8 +--
 reco/L1/CbmL1ReadEvent.cxx                    | 57 ++++++++-----------
 reco/L1/qa/CbmCaInputQaTof.cxx                |  4 +-
 9 files changed, 73 insertions(+), 86 deletions(-)

diff --git a/core/base/CbmTrackingDetectorInterfaceBase.h b/core/base/CbmTrackingDetectorInterfaceBase.h
index ebbc8ce740..48807029ec 100644
--- a/core/base/CbmTrackingDetectorInterfaceBase.h
+++ b/core/base/CbmTrackingDetectorInterfaceBase.h
@@ -20,7 +20,7 @@
 
 #include <cmath>
 
-class CbmPixelHit;
+class CbmHit;
 
 /// @class CbmTrackingDetectorInterfaceBase
 /// @brief Abstract class, which should be inherited by every detecting subsystem tracking interface class
@@ -45,10 +45,10 @@ public:
   /// @return [phiU, phiV] - Stereo angles [rad]
   virtual std::tuple<double, double> GetStereoAnglesSensor(int address) const = 0;
 
-  /// @brief  Gets a tracking station of a CbmPixelHit
-  /// @param  hit  A pointer to CbmPixelHit
+  /// @brief  Gets a tracking station of a CbmHit
+  /// @param  hit  A pointer to CbmHit
   /// @return Local index of the tracking station
-  virtual int GetTrackingStationIndex(const CbmPixelHit* hit) const = 0;
+  virtual int GetTrackingStationIndex(const CbmHit* hit) const = 0;
 
   /// @brief  Gets a tracking station by the address of element
   /// @param  address  Unique element address
diff --git a/core/detectors/much/CbmMuchTrackingInterface.cxx b/core/detectors/much/CbmMuchTrackingInterface.cxx
index 465a4b898a..0c29a6cd81 100644
--- a/core/detectors/much/CbmMuchTrackingInterface.cxx
+++ b/core/detectors/much/CbmMuchTrackingInterface.cxx
@@ -78,8 +78,9 @@ InitStatus CbmMuchTrackingInterface::Init()
   }
 
   // Keep Track of how many layers are present in each station (may be heterogenous)
-  for (int iMuchSt = 0; iMuchSt < fGeoScheme->GetNStations(); ++iMuchSt) {
-    fNbLayersPerStation.push_back(fGeoScheme->GetStation(iMuchSt)->GetNLayers());
+  for (int iMuchSt = 0, iFirstLayer = 0; iMuchSt < fGeoScheme->GetNStations(); ++iMuchSt) {
+    fFirstTrackingStation.push_back(iFirstLayer);
+    iFirstLayer += fGeoScheme->GetStation(iMuchSt)->GetNLayers();
   }
 
   // Check the validity of the parameters
@@ -106,3 +107,13 @@ void CbmMuchTrackingInterface::SetParContainers() {}
 
 //-------------------------------------------------------------------------------------------------------------------------------------
 //
+
+std::pair<int, int> CbmMuchTrackingInterface::ConvTrackingStationId2MuchId(int traStationId) const
+{
+  int muchStation = fFirstTrackingStation.size() - 1;
+  while (muchStation > 0 && (fFirstTrackingStation[muchStation] > traStationId)) {
+    muchStation--;
+  }
+  int muchLayer = traStationId - fFirstTrackingStation[muchStation];
+  return std::pair(muchStation, muchLayer);
+}
diff --git a/core/detectors/much/CbmMuchTrackingInterface.h b/core/detectors/much/CbmMuchTrackingInterface.h
index ee8ec96393..8989e4719c 100644
--- a/core/detectors/much/CbmMuchTrackingInterface.h
+++ b/core/detectors/much/CbmMuchTrackingInterface.h
@@ -12,11 +12,11 @@
 #ifndef CbmMuchTrackingInterface_h
 #define CbmMuchTrackingInterface_h 1
 
+#include "CbmHit.h"
 #include "CbmMuchGeoScheme.h"
 #include "CbmMuchModuleGem.h"
 #include "CbmMuchPad.h"
 #include "CbmMuchStation.h"
-#include "CbmPixelHit.h"
 #include "CbmTrackingDetectorInterfaceBase.h"
 
 #include "FairTask.h"
@@ -85,15 +85,19 @@ public:
   /// @return max Z of the station [cm]
   double GetZmax(int stationId) const { return GetZref(stationId) + 0.5 * GetMuchLayer(stationId)->GetDz(); }
 
-  /// @brief  Gets a tracking station of a CbmPixelHit
-  /// @param  hit  A pointer to CbmPixelHit
+  /// @brief  Gets a tracking station of a CbmHit
+  /// @param  hit  A pointer to CbmHit
   /// @return Local index of the tracking station
-  int GetTrackingStationIndex(const CbmPixelHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
+  int GetTrackingStationIndex(const CbmHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
 
   /// @brief  Gets a tracking station by the address of element
   /// @param  address  Unique element address
   /// @return Local index of the tracking station
-  int GetTrackingStationIndex(int address) const { return CbmMuchGeoScheme::GetStationIndex(address); }
+  int GetTrackingStationIndex(int address) const
+  {
+    return fFirstTrackingStation.at(CbmMuchGeoScheme::GetStationIndex(address))
+           + CbmMuchGeoScheme::GetLayerIndex(address);
+  }
 
   // TODO: SZh 30.08.2023: Provide automatic definintion of bounds
   /// @brief  Gets upper bound of a station along the X-axis
@@ -137,45 +141,24 @@ private:
   /// @brief  Gets pointer to the TRD module
   /// @param  stationId  Tracking staton ID
   /// @return Pointer to the particular CbmTrdParModDigi object
-  __attribute__((always_inline)) CbmMuchLayer* GetMuchLayer(int stationId) const
+  __attribute__((always_inline)) CbmMuchLayer* GetMuchLayer(int traStationId) const
   {
-    return fGeoScheme->GetLayer(GetMuchStationId(stationId), GetMuchLayerId(stationId));
-  }
-
-  /// @brief Calculates MuCh layer ID from tracker station ID
-  /// @param stationId  Index of the tracking station
-  __attribute__((always_inline)) int GetMuchLayerId(int stationId) const
-  {
-    int layersCount  = 0;
-    auto itStatNbLay = fNbLayersPerStation.cbegin();
-    while ((layersCount + *itStatNbLay) <= stationId) {
-      layersCount += *itStatNbLay;
-      itStatNbLay++;
-    }
-    return stationId - layersCount;
+    auto [muchSta, muchLayer] = ConvTrackingStationId2MuchId(traStationId);
+    return fGeoScheme->GetLayer(muchSta, muchLayer);
   }
 
   /// @brief Calculates MuCh station ID from tracker station ID
   /// @param stationId  Index of the tracking station
-  __attribute__((always_inline)) int GetMuchStationId(int stationId) const
-  {
-    int layersCount   = 0;
-    int stationsCount = 0;
-    // Assume that never station ID as input bigger than total number of stations
-    while ((layersCount + fNbLayersPerStation[stationsCount]) <= stationId) {
-      layersCount += fNbLayersPerStation[stationsCount];
-      stationsCount++;
-    }
-    return stationsCount;
-  }
+  std::pair<int, int> ConvTrackingStationId2MuchId(int traStationId) const;
 
   static CbmMuchTrackingInterface* fpInstance;  ///< Instance of the class
 
   CbmMuchGeoScheme* fGeoScheme {nullptr};  ///< MuCh geometry scheme instance
 
-  std::vector<uint16_t> fNbLayersPerStation = {};
+  std::vector<int> fFirstTrackingStation {};
 
   ClassDef(CbmMuchTrackingInterface, 0);
 };
 
+
 #endif  // CbmMuchTrackingInterface
diff --git a/core/detectors/mvd/CbmMvdTrackingInterface.h b/core/detectors/mvd/CbmMvdTrackingInterface.h
index c2854be056..5c75038dcf 100644
--- a/core/detectors/mvd/CbmMvdTrackingInterface.h
+++ b/core/detectors/mvd/CbmMvdTrackingInterface.h
@@ -12,10 +12,10 @@
 #ifndef CbmMvdTrackingInterface_h
 #define CbmMvdTrackingInterface_h 1
 
+#include "CbmHit.h"                            // for CbmHit
 #include "CbmMvdDetectorId.h"                  // for CbmMvdDetectorId
 #include "CbmMvdHit.h"                         // for CbmMvdHit
 #include "CbmMvdStationPar.h"                  // for CbmMvdStationPar
-#include "CbmPixelHit.h"                       // for CbmPixelHit
 #include "CbmTrackingDetectorInterfaceBase.h"  // for CbmTrackingDetectorInt...
 
 #include <FairTask.h>  // for InitStatus, FairTask
@@ -79,10 +79,10 @@ public:
   /// @return Station thickness [cm]
   [[deprecated]] double GetSensorThickness(int stationId) const { return fMvdStationPar->GetZThickness(stationId); }
 
-  /// @brief  Gets a tracking station of a CbmPixelHit
-  /// @param  hit  A pointer to CbmPixelHit
+  /// @brief  Gets a tracking station of a CbmHit
+  /// @param  hit  A pointer to CbmHit
   /// @return Local index of the tracking station
-  int GetTrackingStationIndex(const CbmPixelHit* hit) const
+  int GetTrackingStationIndex(const CbmHit* hit) const
   {
     auto hitMvd = [&] {
       if constexpr (kUseDynamicCast) { return dynamic_cast<const CbmMvdHit*>(hit); }
diff --git a/core/detectors/sts/CbmStsTrackingInterface.h b/core/detectors/sts/CbmStsTrackingInterface.h
index 5cc2203fd0..941a597ece 100644
--- a/core/detectors/sts/CbmStsTrackingInterface.h
+++ b/core/detectors/sts/CbmStsTrackingInterface.h
@@ -12,7 +12,7 @@
 #ifndef CbmStsTrackingInterface_h
 #define CbmStsTrackingInterface_h 1
 
-#include "CbmPixelHit.h"
+#include "CbmHit.h"
 #include "CbmStsParSetModule.h"
 #include "CbmStsParSetSensor.h"
 #include "CbmStsParSetSensorCond.h"
@@ -73,10 +73,10 @@ public:
   /// @return [phiU, phiV] - Stereo angles [rad]
   std::tuple<double, double> GetStereoAnglesSensor(int address) const;
 
-  /// @brief  Gets a tracking station of a CbmPixelHit
-  /// @param  hit  A pointer to CbmPixelHit
+  /// @brief  Gets a tracking station of a CbmHit
+  /// @param  hit  A pointer to CbmHit
   /// @return Local index of the tracking station
-  int GetTrackingStationIndex(const CbmPixelHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
+  int GetTrackingStationIndex(const CbmHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
 
   /// @brief  Gets a tracking station by the address of element
   /// @param  address  Unique element address
diff --git a/core/detectors/tof/CbmTofTrackingInterface.h b/core/detectors/tof/CbmTofTrackingInterface.h
index 060fee3695..1ae4ff16a3 100644
--- a/core/detectors/tof/CbmTofTrackingInterface.h
+++ b/core/detectors/tof/CbmTofTrackingInterface.h
@@ -12,7 +12,7 @@
 #ifndef CbmTofTrackingInterface_h
 #define CbmTofTrackingInterface_h 1
 
-#include "CbmPixelHit.h"
+#include "CbmHit.h"
 #include "CbmTofAddress.h"
 #include "CbmTofCell.h"
 #include "CbmTofDigiBdfPar.h"
@@ -67,7 +67,7 @@ public:
   /// @brief  Gets a tracking station of a ToF hit
   /// @param  hit  A pointer to ToF hit
   /// @return Local index of the tracking station
-  int GetTrackingStationIndex(const CbmPixelHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
+  int GetTrackingStationIndex(const CbmHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
 
   /// @brief  Gets a tracking station by the address of element
   /// @param  address  Unique element address
diff --git a/core/detectors/trd/CbmTrdTrackingInterface.h b/core/detectors/trd/CbmTrdTrackingInterface.h
index be853c115f..9f6670d826 100644
--- a/core/detectors/trd/CbmTrdTrackingInterface.h
+++ b/core/detectors/trd/CbmTrdTrackingInterface.h
@@ -12,7 +12,7 @@
 #ifndef CbmTrdTrackingInterface_h
 #define CbmTrdTrackingInterface_h 1
 
-#include "CbmPixelHit.h"
+#include "CbmHit.h"
 #include "CbmTrackingDetectorInterfaceBase.h"
 #include "CbmTrdAddress.h"
 #include "CbmTrdParModDigi.h"
@@ -63,10 +63,10 @@ public:
   /// @return [phiU, phiV] - Stereo angles [rad]
   std::tuple<double, double> GetStereoAnglesSensor(int /*address*/) const { return std::tuple(0., TMath::Pi() / 2.); }
 
-  /// @brief  Gets a tracking station of a CbmPixelHit
-  /// @param  hit  A pointer to CbmPixelHit
+  /// @brief  Gets a tracking station of a CbmHit
+  /// @param  hit  A pointer to CbmHit
   /// @return Local index of the tracking station
-  int GetTrackingStationIndex(const CbmPixelHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
+  int GetTrackingStationIndex(const CbmHit* hit) const { return GetTrackingStationIndex(hit->GetAddress()); }
 
   /// @brief  Gets a tracking station by the address
   /// @param  address  Unique element address
diff --git a/reco/L1/CbmL1ReadEvent.cxx b/reco/L1/CbmL1ReadEvent.cxx
index 5d90973d45..46ade00051 100644
--- a/reco/L1/CbmL1ReadEvent.cxx
+++ b/reco/L1/CbmL1ReadEvent.cxx
@@ -23,21 +23,23 @@
 #include "CbmL1.h"
 #include "CbmMCDataObject.h"
 #include "CbmMatch.h"
-#include "CbmMuchGeoScheme.h"
 #include "CbmMuchPixelHit.h"
 #include "CbmMuchPoint.h"
 #include "CbmMuchTrackingInterface.h"
+#include "CbmMvdTrackingInterface.h"
 #include "CbmStsAddress.h"
 #include "CbmStsCluster.h"
 #include "CbmStsDigi.h"
 #include "CbmStsHit.h"
 #include "CbmStsPoint.h"
 #include "CbmStsSetup.h"
+#include "CbmStsTrackingInterface.h"
 #include "CbmTofHit.h"
 #include "CbmTofPoint.h"
 #include "CbmTofTrackingInterface.h"
 #include "CbmTrdHit.h"
 #include "CbmTrdPoint.h"
+#include "CbmTrdTrackingInterface.h"
 
 #include "FairMCEventHeader.h"
 
@@ -173,8 +175,8 @@ void CbmL1::ReadEvent(CbmEvent* event)
     }
 
     for (DFSET::iterator set_it = fvSelectedMcEvents.begin(); set_it != fvSelectedMcEvents.end(); ++set_it) {
-      Int_t iFile   = set_it->first;
-      Int_t iEvent  = set_it->second;
+      Int_t iFile  = set_it->first;
+      Int_t iEvent = set_it->second;
       if (fUseMVD && fpMvdPoints) {
         Int_t fNpointsMvdInEvent = fpMvdPoints->Size(iFile, iEvent);
         double maxDeviation      = 0;
@@ -440,13 +442,13 @@ void CbmL1::ReadEvent(CbmEvent* event)
       CbmMvdHit* h = L1_DYNAMIC_CAST<CbmMvdHit*>(fpMvdHits->At(hitIndex));
       {
         th.ExtIndex = hitIndex;
-        th.iStation = h->GetStationNr();
+        th.iStation = fpAlgo->GetParameters()->GetStationIndexActive(
+          CbmMvdTrackingInterface::Instance()->GetTrackingStationIndex(h), L1DetectorID::kMvd);
+
+        if (th.iStation < 0) continue;
 
-        int stIdx = fpAlgo->GetParameters()->GetStationIndexActive(h->GetStationNr(), L1DetectorID::kMvd);
-        if (stIdx == -1) continue;
-        th.iStation = stIdx;
-        th.iStripF  = firstDetStrip + j;
-        th.iStripB  = th.iStripF;
+        th.iStripF = firstDetStrip + j;
+        th.iStripB = th.iStripF;
         if (NStrips <= th.iStripF) { NStrips = th.iStripF + 1; }
 
         TVector3 pos, err;
@@ -501,15 +503,14 @@ void CbmL1::ReadEvent(CbmEvent* event)
       {
         th.ExtIndex = hitIndex;
         th.Det      = 1;
-        int stIdx   = fpAlgo->GetParameters()->GetStationIndexActive(
-          CbmStsSetup::Instance()->GetStationNumber(h->GetAddress()), L1DetectorID::kSts);
 
-        if (stIdx == -1) continue;
+        th.iStation = fpAlgo->GetParameters()->GetStationIndexActive(
+          CbmStsTrackingInterface::Instance()->GetTrackingStationIndex(h), L1DetectorID::kSts);
 
+        if (th.iStation < 0) continue;
 
-        th.iStation = stIdx;
-        th.iStripF  = firstDetStrip + h->GetFrontClusterId();
-        th.iStripB  = firstDetStrip + h->GetBackClusterId();
+        th.iStripF = firstDetStrip + h->GetFrontClusterId();
+        th.iStripB = firstDetStrip + h->GetBackClusterId();
 
         if (NStrips <= th.iStripF) { NStrips = th.iStripF + 1; }
         if (NStrips <= th.iStripB) { NStrips = th.iStripB + 1; }
@@ -557,14 +558,10 @@ void CbmL1::ReadEvent(CbmEvent* event)
         th.ExtIndex = hitIndex;
         th.Det      = 2;
 
-        Int_t stationNumber = CbmMuchGeoScheme::GetStationIndex(h->GetAddress());
-        Int_t layerNumber   = CbmMuchGeoScheme::GetLayerIndex(h->GetAddress());
-
-        int DetId = stationNumber * 3 + layerNumber;
+        th.iStation = fpAlgo->GetParameters()->GetStationIndexActive(
+          CbmMuchTrackingInterface::Instance()->GetTrackingStationIndex(h), L1DetectorID::kMuch);
 
-        int stIdx = fpAlgo->GetParameters()->GetStationIndexActive(DetId, L1DetectorID::kMuch);
-        if (stIdx == -1) continue;
-        th.iStation = stIdx;  //h->GetStationNr() - 1;
+        if (th.iStation < 0) continue;
 
         //Get time
         th.time = h->GetTime();
@@ -619,11 +616,10 @@ void CbmL1::ReadEvent(CbmEvent* event)
       th.ExtIndex = hitIndex;
       th.Det      = 3;
 
-      int stIdx = fpAlgo->GetParameters()->GetStationIndexActive(h->GetPlaneId(), L1DetectorID::kTrd);
-      //assert(stIdx != -1);
-      if (stIdx == -1) continue;  // Station was disabled and thus should not be proceede
+      th.iStation = fpAlgo->GetParameters()->GetStationIndexActive(
+        CbmTrdTrackingInterface::Instance()->GetTrackingStationIndex(h), L1DetectorID::kTrd);
 
-      th.iStation = stIdx;
+      if (th.iStation < 0) continue;
 
       //  if (h->GetPlaneId()==0) continue;
       //  if (h->GetPlaneId()==1) continue;
@@ -683,9 +679,10 @@ void CbmL1::ReadEvent(CbmEvent* event)
 
       if (5 == CbmTofAddress::GetSmType(h->GetAddress())) continue;  // Skip T0 hits from TOF
 
-      int sttof = CbmTofTrackingInterface::Instance()->GetTrackingStationIndex(h);
+      th.iStation = fpAlgo->GetParameters()->GetStationIndexActive(
+        CbmTofTrackingInterface::Instance()->GetTrackingStationIndex(h), L1DetectorID::kTof);
 
-      if (sttof < 0) continue;
+      if (th.iStation < 0) continue;
 
       th.time = h->GetTime();
       th.dt   = h->GetTimeError();
@@ -711,10 +708,6 @@ void CbmL1::ReadEvent(CbmEvent* event)
       //TODO:  is it still needed here? (S.Zharko)
       if (L1Algo::TrackingMode::kMcbm == fTrackingMode && th.z > 400) continue;
 
-      int stIdx = fpAlgo->GetParameters()->GetStationIndexActive(sttof, L1DetectorID::kTof);
-      if (stIdx == -1) continue;
-      th.iStation = stIdx;
-
       th.iMC = fPerformance ? MatchHitWithMc<L1DetectorID::kTof>(hitIndex) : -1;
 
       th.fDataStream = (static_cast<int64_t>(th.Det) << 60) | h->GetAddress();
diff --git a/reco/L1/qa/CbmCaInputQaTof.cxx b/reco/L1/qa/CbmCaInputQaTof.cxx
index 7554768a3c..cfe02e1e97 100644
--- a/reco/L1/qa/CbmCaInputQaTof.cxx
+++ b/reco/L1/qa/CbmCaInputQaTof.cxx
@@ -236,12 +236,12 @@ InitStatus CbmCaInputQaTof::InitHistograms()
         title += ";x_{hit} [cm];y_{hit} [cm]";
         fvph_hit_xy_vs_cell[iSmType][iSm][iRpc] =
           MakeQaObject<TH2F>(name, title, fNbXo, fLoXo, fUpXo, fNbYo, fLoYo, fUpYo);
-        name  = Form("%s/occup_zx_smt%d_sm%d_rpc%d_ch%d", dir, iSmType, iSm, iRpc);
+        name  = Form("%s/occup_zx_smt%d_sm%d_rpc%d", dir, iSmType, iSm, iRpc);
         title = Form("Hit Occupancy in zx-Plane for iSmType = %d, iSm = %d, iRpc = %d", iSmType, iSm, iRpc);
         title += ";z_{hit} [cm];x_{hit} [cm]";
         fvph_hit_zx_vs_cell[iSmType][iSm][iRpc] =
           MakeQaObject<TH2F>(name, title, fNbinsZ, frZmin.back(), frZmax.back(), fNbXo, fLoXo, fUpXo);
-        name  = Form("%s/occup_zy_smt%d_sm%d_rpc%d_ch%d", dir, iSmType, iSm, iRpc);
+        name  = Form("%s/occup_zy_smt%d_sm%d_rpc%d", dir, iSmType, iSm, iRpc);
         title = Form("Hit Occupancy in zy-Plane for iSmType = %d, iSm = %d, iRpc = %d", iSmType, iSm, iRpc);
         title += ";z_{hit} [cm];y_{hit} [cm]";
         fvph_hit_zy_vs_cell[iSmType][iSm][iRpc] =
-- 
GitLab