AliceO2Group · mfasDa · Jul 19, 2024 · Jul 12, 2024
@@ -14,6 +14,7 @@
 #include <iosfwd>
 #include "SimulationDataFormat/BaseHits.h"
 #include "CommonUtils/ShmAllocator.h"
+#include <boost/container_hash/hash.hpp>

 namespace o2::focal
 {
@@ -35,6 +36,39 @@ class Hit : public o2::BasicXYZEHit<float>
    UNKNOWN  ///< Undefined
  };

+  /// \struct HitID
+  /// \brief Mapped information of a channel
+  struct HitID {
+    int mParentID;   ///< parentID of the track creating this hit
+    uint8_t mRow;    ///< Row of the hit in the calorimeter
+    uint8_t mColumn; ///< Column of the hit in the calorimeter
+    uint8_t mLayer;  ///< Layer the was hit
+
+    bool operator==(const HitID& other) const
+    {
+      return mParentID == other.mParentID && mRow == other.mRow && mColumn == other.mColumn && mLayer == other.mLayer;
+    }
+    friend std::ostream& operator<<(std::ostream& stream, const Hit::HitID& channel);
+  };
+
+  /// \struct HitIDHasher
+  /// \brief Hash functor for hit ID
+  struct HitIDHasher {
+
+    /// \brief Functor implementation
+    /// \param s Hit for which to determine a hash value
+    /// \return hash value for channel ID
+    size_t operator()(const HitID& s) const
+    {
+      std::size_t seed = 0;
+      boost::hash_combine(seed, s.mParentID);
+      boost::hash_combine(seed, s.mRow);
+      boost::hash_combine(seed, s.mColumn);
+      boost::hash_combine(seed, s.mLayer);
+      return seed;
+    }
+  };
+
  /// \brief Dummy constructor
  Hit() = default;


@@ -184,13 +184,14 @@ class Detector : public o2::base::DetImpl<Detector>

  Geometry* mGeometry; //!<! Geometry pointer

-  int mMedSensHCal = 0; //!<! Sensitive Medium for HCal
+  int mMedSensHCal = 0;    //!<! Sensitive Medium for HCal
  int mMedSensECalPad = 0; //!<! Sensitive Medium for ECal Pads
  int mMedSensECalPix = 0; //!<! Sensitive Medium for ECal Pixels

  std::vector<const Composition*> mGeoCompositions; //!<! list of FOCAL compositions

-  std::vector<o2::focal::Hit>* mHits; ///< Container with hits
+  std::vector<o2::focal::Hit>* mHits;                                              ///< Container with hits
+  std::unordered_map<Hit::HitID, unsigned int, Hit::HitIDHasher> mHitIndexMapping; ///< Mapping the hits to a cell in the detector

  std::unordered_map<int, int> mSuperParentsIndices; //!<! Super parent indices (track index - superparent index)
  std::unordered_map<int, Parent> mSuperParents;     //!<! Super parent kine info (superparent index - superparent object)

@@ -28,6 +28,7 @@ using namespace o2::focal;
 Detector::Detector(bool active, std::string geofilename)
  : o2::base::DetImpl<Detector>("FOC", active),
    mHits(o2::utils::createSimVector<Hit>()),
+    mHitIndexMapping(),
    mGeometry(nullptr),
    mMedSensHCal(-1),
    mMedSensECalPad(-1),
@@ -166,23 +167,23 @@ Hit* Detector::AddHit(int trackID, int primary, double initialEnergy, int detID,
  LOG(debug3) << "Adding hit for track " << trackID << " with position (" << pos.X() << ", "
              << pos.Y() << ", " << pos.Z() << ")  with energy " << initialEnergy << " loosing " << eLoss;
  mHits->emplace_back(primary, trackID, detID, subsystem, initialEnergy, pos, time, eLoss);
+
+  auto [isin, col, row, layer, segment] = mGeometry->getVirtualInfo(pos.X(), pos.Y(), pos.Z());
+  mHitIndexMapping.insert(std::pair<Hit::HitID, unsigned int>({trackID, uint8_t(row), uint8_t(col), uint8_t(layer)}, static_cast<unsigned int>(mHits->size() - 1)));
+
  return &(mHits->back());
 }

 Hit* Detector::FindHit(int parentID, int col, int row, int layer)
 {
+  Hit::HitID hitToFind{parentID, uint8_t(row), uint8_t(col), uint8_t(layer)};
+  auto found = mHitIndexMapping.find(hitToFind);

-  auto HitComparison = [&](const Hit& hit) {
-    auto information = mGeometry->getVirtualInfo(hit.GetX(), hit.GetY(), hit.GetZ());
-    // FIXME Should we compare segments instead of layers ???
-    return hit.GetTrackID() == parentID && col == std::get<1>(information) && row == std::get<2>(information) && layer == std::get<3>(information);
-  };
-
-  auto result = std::find_if(mHits->begin(), mHits->end(), HitComparison);
-  if (result == mHits->end()) {
+  if (found == mHitIndexMapping.end()) {
    return nullptr;
  }
-  return &(*result);
+
+  return &((*mHits)[found->second]);
 }

 Parent* Detector::AddSuperparent(int trackID, int pdg, double energy)
@@ -205,6 +206,7 @@ void Detector::Reset()
  if (!o2::utils::ShmManager::Instance().isOperational()) {
    mHits->clear();
  }
+  mHitIndexMapping.clear();

  mSuperParentsIndices.clear();
  mSuperParents.clear();
@@ -928,12 +930,12 @@ void Detector::CreateECALGeometry()

  if (geom->getTowerGapMaterial() == "Cu") { // Copper
    // if (contains(geom->getTowerGapMaterial(), "Cu")) { // Copper
-    volumeColdPlate = new TGeoVolume("volColdPlate", coldPlateBox, gGeoManager->GetMedium("FOCAL_Cu$"));
+    volumeColdPlate = new TGeoVolume("volColdPlate", coldPlateBox, gGeoManager->GetMedium(getMediumID(ID_COPPER)));
  } else if (geom->getTowerGapMaterial() == "Al") { // Aluminium
    // else if (contains(geom->getTowerGapMaterial(), "Al")) {   // Aluminium
-    volumeColdPlate = new TGeoVolume("volColdPlate", coldPlateBox, gGeoManager->GetMedium("FOCAL_AlPlate"));
+    volumeColdPlate = new TGeoVolume("volColdPlate", coldPlateBox, gGeoManager->GetMedium(getMediumID(ID_ALUMINIUM)));
  } else {
-    volumeColdPlate = new TGeoVolume("volColdPlate", coldPlateBox, gGeoManager->GetMedium("FOCAL_AirGaps$"));
+    volumeColdPlate = new TGeoVolume("volColdPlate", coldPlateBox, gGeoManager->GetMedium(getMediumID(ID_AIR)));
  }
  // mSensitiveECALPad.push_back(volumeColdPlate->GetName());
  mSensitive.push_back(volumeColdPlate->GetName());
@@ -957,7 +959,7 @@ void Detector::CreateECALGeometry()
  // Create SiPad box for the two sensitive layers to be placed in front of ECAL
  TGeoBBox* siPadBox = new TGeoBBox("SiPadBox", geom->getTowerSizeX() / 2. + geom->getTowerGapSizeX() / 2.,
                                    geom->getTowerSizeY() / 2. + geom->getTowerGapSizeY() / 2., 0.03 / 2.0);
-  TGeoVolume* volumeSiPad = new TGeoVolume("volSiPad", siPadBox, gGeoManager->GetMedium("FOCAL_SiSens$"));
+  TGeoVolume* volumeSiPad = new TGeoVolume("volSiPad", siPadBox, gGeoManager->GetMedium(getMediumID(ID_SILICON)));
  volumeSiPad->SetLineColor(kOrange + 7);
  // mSensitiveECALPad.push_back(volumeSiPad->GetName());
  mSensitive.push_back(volumeSiPad->GetName());