Update indexing scheme - fix bug

Author: scannito

Detectors/Upgrades/ALICE3/TRK/base/include/TRKBase/GeometryTGeo.h

@@ -83,14 +83,14 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache
 
   int getSubDetID(int index) const;
   int getPetalCase(int index) const;
+  int getDisk(int index) const;
   int getLayer(int index) const;
   int getStave(int index) const;
   int getHalfStave(int index) const;
-  int getDisk(int index) const;
   int getModule(int index) const;
   int getChip(int index) const;
 
-  /// This routine computes the chip index number from the subDetID, petal, disk, layer, stave /// TODO: retrieve also from chip when chips will be available
+  /// This routine computes the chip index number from the subDetID, petal, disk, layer, stave, half stave, module, chip
   /// \param int subDetID The subdetector ID, 0 for VD, 1 for MLOT
   /// \param int petalcase The petal case number for VD, from 0 to 3
   /// \param int disk The disk number for VD, from 0 to 5
@@ -101,7 +101,7 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache
   /// \param int chip The chip number for MLOT, from 0 to 8
   int getChipIndex(int subDetID, int petalcase, int disk, int lay, int stave, int halfstave, int mod, int chip) const;
 
-  /// This routine computes the chip index number from the subDetID, volume, layer, stave /// TODO: retrieve also from chip when chips will be available
+  /// This routine computes the chip index number from the subDetID, volume, layer, stave, half stave, module, chip
   /// \param int subDetID The subdetector ID, 0 for VD, 1 for MLOT
   /// \param int volume is needed only with the current configuration for VD where each single element is a volume. // TODO: when the geometry naming scheme will be changed, change this method
   /// \param int lay The layer number for the MLOT. In the current configuration for VD this is not needed. // TODO: when the geometry naming scheme will be changed, change this method
@@ -111,7 +111,7 @@ class GeometryTGeo : public o2::detectors::DetMatrixCache
   /// \param int chip The chip number for MLOT, from 0 to 8
   int getChipIndex(int subDetID, int volume, int lay, int stave, int halfstave, int mod, int chip) const;
 
-  /// This routine computes subDetID, petal, disk, layer, stave given the chip index number /// TODO: copute also from chip when chips will be available
+  /// This routine computes subDetID, petal, disk, layer, stave, half stave, module, chip, given the chip index number
   /// \param int index The chip index number, starting from 0
   /// \param int subDetID The subdetector ID, 0 for VD, 1 for MLOT
   /// \param int petalcase The petal case number for VD, from 0 to 3

Detectors/Upgrades/ALICE3/TRK/base/src/GeometryTGeo.cxx

@@ -148,6 +148,22 @@ int GeometryTGeo::getPetalCase(int index) const
   return petalcase;
 }
 
+//__________________________________________________________________________
+int GeometryTGeo::getDisk(int index) const
+{
+  int subDetID = getSubDetID(index);
+  int petalcase = getPetalCase(index);
+
+  if (subDetID == 0) { /// VD
+    if (index % mNumberOfChipsPerPetalVD[petalcase] < mNumberOfLayersVD) {
+      return -1; /// layers
+    }
+    return (index % mNumberOfChipsPerPetalVD[petalcase]) - mNumberOfLayersVD;
+  }
+
+  return -1; /// not found or ML/OT
+}
+
 //__________________________________________________________________________
 int GeometryTGeo::getLayer(int index) const
 {
@@ -180,9 +196,23 @@ int GeometryTGeo::getStave(int index) const
   } else if (subDetID == 1) { /// MLOT
     int lay = getLayer(index);
     index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer
-    return index / mNumberOfHalfStaves[lay];
+
+    const int Nhs = mNumberOfHalfStaves[lay];
+    const int Nmod = mNumberOfModules[lay];
+    const int Nchip = mNumberOfChips[lay];
+
+    if (Nhs == 2) {
+      int chipsPerModule = Nchip;
+      int chipsPerHalfStave = Nmod * chipsPerModule;
+      int chipsPerStave = Nhs * chipsPerHalfStave;
+      return index / chipsPerStave;
+    } else if (Nhs == 1) {
+      int chipsPerModule = Nchip;
+      int chipsPerStave = Nmod * chipsPerModule;
+      return index / chipsPerStave;
+    }
   }
-  return -1; /// not found
+  return -1;
 }
 
 //__________________________________________________________________________
@@ -191,48 +221,85 @@ int GeometryTGeo::getHalfStave(int index) const
   int subDetID = getSubDetID(index);
   int lay = getLayer(index);
   int petalcase = getPetalCase(index);
-  int stave = getStave(index);
 
   if (subDetID == 0) { /// VD
     return -1;
   } else if (subDetID == 1) { /// MLOT
     int lay = getLayer(index);
     index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer
-    return index % 2;                                     /// 0 = half stave left, 1 = half stave right, as geometry is filled /// TODO: generalize once chips will be in place. Can it be working also with chips?
+
+    const int Nhs = mNumberOfHalfStaves[lay];
+    const int Nmod = mNumberOfModules[lay];
+    const int Nchip = mNumberOfChips[lay];
+
+    int chipsPerModule = Nchip;
+    int chipsPerHalfStave = Nmod * chipsPerModule;
+    int chipsPerStave = Nhs * chipsPerHalfStave;
+
+    int rem = index % chipsPerStave;
+    return rem / chipsPerHalfStave; // 0 = left, 1 = right
   }
-  return -1; /// not found
+  return -1;
 }
 
 //__________________________________________________________________________
-int GeometryTGeo::getDisk(int index) const
+int GeometryTGeo::getModule(int index) const
 {
   int subDetID = getSubDetID(index);
+  int lay = getLayer(index);
   int petalcase = getPetalCase(index);
 
   if (subDetID == 0) { /// VD
-    if (index % mNumberOfChipsPerPetalVD[petalcase] < mNumberOfLayersVD) {
-      return -1; /// layers
+    return -1;
+  } else if (subDetID == 1) { /// MLOT
+    int lay = getLayer(index);
+    index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer
+
+    const int Nhs = mNumberOfHalfStaves[lay];
+    const int Nmod = mNumberOfModules[lay];
+    const int Nchip = mNumberOfChips[lay];
+
+    if (Nhs == 2) {
+      int chipsPerModule = Nchip;
+      int chipsPerHalfStave = Nmod * chipsPerModule;
+      int rem = index % (Nhs * chipsPerHalfStave);
+      rem = rem % chipsPerHalfStave;
+      return rem / chipsPerModule;
+    } else if (Nhs == 1) {
+      int chipsPerModule = Nchip;
+      int rem = index % (Nmod * chipsPerModule);
+      return rem / chipsPerModule;
     }
-    return (index % mNumberOfChipsPerPetalVD[petalcase]) - mNumberOfLayersVD;
   }
-
-  return -1; /// not found or ML/OT
+  return -1;
 }
 
 //__________________________________________________________________________
-int GeometryTGeo::getModule(int index) const
+int GeometryTGeo::getChip(int index) const
 {
   int subDetID = getSubDetID(index);
+  int lay = getLayer(index);
+  int petalcase = getPetalCase(index);
 
-  return -1; /// not implemented yet
-}
+  if (subDetID == 0) { /// VD
+    return -1;
+  } else if (subDetID == 1) { /// MLOT
+    int lay = getLayer(index);
+    index -= getFirstChipIndex(lay, petalcase, subDetID); // get the index of the sensing element in the layer
 
-//__________________________________________________________________________
-int GeometryTGeo::getChip(int index) const
-{
-  int subDetID = getSubDetID(index);
+    const int Nhs = mNumberOfHalfStaves[lay];
+    const int Nmod = mNumberOfModules[lay];
+    const int Nchip = mNumberOfChips[lay];
 
-  return -1; /// not implemented yet
+    if (Nhs == 2) {
+      int chipsPerModule = Nchip;
+      return index % chipsPerModule;
+    } else if (Nhs == 1) {
+      int chipsPerModule = Nchip;
+      return index % chipsPerModule;
+    }
+  }
+  return -1;
 }
 
 //__________________________________________________________________________
@@ -244,11 +311,20 @@ int GeometryTGeo::getChipIndex(int subDetID, int petalcase, int disk, int lay, i
     } else { // layer
       return getFirstChipIndex(lay, petalcase, subDetID) + lay;
     }
-  } else if (subDetID == 1) {            // MLOT
-    if (mNumberOfHalfStaves[lay] == 2) { // staggered geometry
-      return getFirstChipIndex(lay, petalcase, subDetID) + stave * mNumberOfHalfStaves[lay] + halfstave * mNumberOfModules[lay] + mod * mNumberOfChips[lay] + chip;
-    } else if (mNumberOfHalfStaves[lay] == 1) { // turbo geometry
-      return getFirstChipIndex(lay, petalcase, subDetID) + stave * mNumberOfModules[lay] + mod * mNumberOfChips[lay] + chip;
+  } else if (subDetID == 1) {                 // MLOT
+    const int Nhs = mNumberOfHalfStaves[lay]; // 1 or 2
+    const int Nmod = mNumberOfModules[lay];   // module per half-stave (per stave if Nhs==1)
+    const int Nchip = mNumberOfChips[lay];    // chips per module
+
+    if (Nhs == 2) { // staggered geometry: layer -> stave -> halfstave -> mod -> chip
+      int chipsPerModule = Nchip;
+      int chipsPerHalfStave = Nmod * chipsPerModule;
+      int chipsPerStave = Nhs * chipsPerHalfStave;
+      return getFirstChipIndex(lay, petalcase, subDetID) + stave * chipsPerStave + halfstave * chipsPerHalfStave + mod * chipsPerModule + chip;
+    } else if (Nhs == 1) { // turbo geometry: layer -> stave -> mod -> chip (no halfstave)
+      int chipsPerModule = Nchip;
+      int chipsPerStave = Nmod * chipsPerModule;
+      return getFirstChipIndex(lay, petalcase, subDetID) + stave * chipsPerStave + mod * chipsPerModule + chip;
     }
   }
   return -1; // not found
@@ -260,11 +336,20 @@ int GeometryTGeo::getChipIndex(int subDetID, int volume, int lay, int stave, int
   if (subDetID == 0) { // VD
     return volume;     /// In the current configuration for VD, each volume is the sensor element = chip. // TODO: when the geometry naming scheme will be changed, change this method
 
-  } else if (subDetID == 1) {            // MLOT
-    if (mNumberOfHalfStaves[lay] == 2) { // staggered geometry
-      return getFirstChipIndex(lay, -1, subDetID) + stave * mNumberOfHalfStaves[lay] + halfstave * mNumberOfModules[lay] + mod * mNumberOfChips[lay] + chip;
-    } else if (mNumberOfHalfStaves[lay] == 1) { // turbo geometry
-      return getFirstChipIndex(lay, -1, subDetID) + stave * mNumberOfModules[lay] + mod * mNumberOfChips[lay] + chip;
+  } else if (subDetID == 1) {                 // MLOT
+    const int Nhs = mNumberOfHalfStaves[lay]; // 1 or 2
+    const int Nmod = mNumberOfModules[lay];   // module per half-stave (per stave if Nhs==1)
+    const int Nchip = mNumberOfChips[lay];    // chips per module
+
+    if (Nhs == 2) { // staggered geometry: layer -> stave -> halfstave -> mod -> chip
+      int chipsPerModule = Nchip;
+      int chipsPerHalfStave = Nmod * chipsPerModule;
+      int chipsPerStave = Nhs * chipsPerHalfStave;
+      return getFirstChipIndex(lay, -1, subDetID) + stave * chipsPerStave + halfstave * chipsPerHalfStave + mod * chipsPerModule + chip;
+    } else if (Nhs == 1) { // turbo geometry: layer -> stave -> mod -> chip (no halfstave)
+      int chipsPerModule = Nchip;
+      int chipsPerStave = Nmod * chipsPerModule;
+      return getFirstChipIndex(lay, -1, subDetID) + stave * chipsPerStave + mod * chipsPerModule + chip;
     }
   }
   return -1; // not found
@@ -278,6 +363,11 @@ bool GeometryTGeo::getChipID(int index, int& subDetID, int& petalcase, int& disk
   disk = getDisk(index);
   lay = getLayer(index);
   stave = getStave(index);
+  if (mNumberOfHalfStaves[lay] == 2) {
+    halfstave = getHalfStave(index);
+  } else {
+    halfstave = 0; // if not staggered geometry, return 0
+  }
   halfstave = getHalfStave(index);
   mod = getModule(index);
   chip = getChip(index);

Detectors/Upgrades/ALICE3/TRK/simulation/src/Detector.cxx

@@ -104,9 +104,6 @@ void Detector::configDefault()
   mLayers.clear();
 
   LOGP(warning, "Loading Scoping Document configuration for ALICE3 TRK");
-  // mLayers.emplace_back(0, GeometryTGeo::getTRKLayerPattern() + std::to_string(0), 0.5f, 4, 100.e-4);
-  // mLayers.emplace_back(1, GeometryTGeo::getTRKLayerPattern() + std::to_string(1), 1.2f, 4, 100.e-4);
-  // mLayers.emplace_back(2, GeometryTGeo::getTRKLayerPattern() + std::to_string(2), 2.5f, 4, 100.e-4);
   mLayers.emplace_back(0, GeometryTGeo::getTRKLayerPattern() + std::to_string(0), 3.78f, 10, 100.e-3);
   mLayers.emplace_back(1, GeometryTGeo::getTRKLayerPattern() + std::to_string(1), 7.f, 10, 100.e-3);
   mLayers.emplace_back(2, GeometryTGeo::getTRKLayerPattern() + std::to_string(2), 12.f, 10, 100.e-3);