Merge branch 'AliceO2Group:master' into master

Author: mherzer28

.git-blame-ignore-revs

@@ -38,7 +38,7 @@ jobs:
         id: ml
         # You can override MegaLinter flavor used to have faster performances
         # More info at https://megalinter.io/flavors/
-        uses: oxsecurity/megalinter@v8.5.0
+        uses: oxsecurity/megalinter@v8.7.0
         env:
           # All available variables are described in documentation:
           # https://megalinter.io/configuration/

.github/workflows/mega-linter.yml

@@ -39,4 +39,5 @@ PYTHON_PYRIGHT_CONFIG_FILE: pyproject.toml
 PYTHON_RUFF_CONFIG_FILE: pyproject.toml
 CPP_CPPLINT_FILE_EXTENSIONS: [".C", ".c", ".c++", ".cc", ".cl", ".cpp", ".cu", ".cuh", ".cxx", ".cxx.in", ".h", ".h++", ".hh", ".h.in", ".hpp", ".hxx", ".inc", ".inl", ".macro"]
 CPP_CLANG_FORMAT_FILE_EXTENSIONS: [".C", ".c", ".c++", ".cc", ".cl", ".cpp", ".cu", ".cuh", ".cxx", ".cxx.in", ".h", ".h++", ".hh", ".h.in", ".hpp", ".hxx", ".inc", ".inl", ".macro"]
+CPP_CPPCHECK_ARGUMENTS: --language=c++ --std=c++20 --check-level=exhaustive --suppressions-list=cppcheck_config
 REPOSITORY_GITLEAKS_PR_COMMITS_SCAN: true

.mega-linter.yml

@@ -12,17 +12,22 @@
 o2physics_add_library(ALICE3Core
                       SOURCES TOFResoALICE3.cxx
                               DelphesO2TrackSmearer.cxx
-                      PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore)
+                      PUBLIC_LINK_LIBRARIES O2::Framework
+                                            O2Physics::AnalysisCore)
 
 o2physics_target_root_dictionary(ALICE3Core
                       HEADERS TOFResoALICE3.h
+                              TrackUtilities.h
                               DelphesO2TrackSmearer.h
                       LINKDEF ALICE3CoreLinkDef.h)
 
 o2physics_add_library(FastTracker
                       SOURCES FastTracker.cxx
-                      PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore)
+                              DetLayer.cxx
+                      PUBLIC_LINK_LIBRARIES O2::Framework
+                                            O2Physics::AnalysisCore)
 
 o2physics_target_root_dictionary(FastTracker
                       HEADERS FastTracker.h
+                              DetLayer.h
                       LINKDEF FastTrackerLinkDef.h)

ALICE3/Core/CMakeLists.txt

@@ -0,0 +1,19 @@
+// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
+// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
+// All rights not expressly granted are reserved.
+//
+// This software is distributed under the terms of the GNU General Public
+// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
+//
+// In applying this license CERN does not waive the privileges and immunities
+// granted to it by virtue of its status as an Intergovernmental Organization
+// or submit itself to any jurisdiction.
+
+///
+/// \file   DetLayer.cxx
+/// \author David Dobrigkeit Chinellato
+/// \since  11/03/2021
+/// \brief  Basic struct to hold information regarding a detector layer to be used in fast simulation
+///
+
+#include "DetLayer.h"

ALICE3/Core/DetLayer.cxx

@@ -10,6 +10,7 @@
 // or submit itself to any jurisdiction.
 
 #include <vector>
+#include <string>
 #include "TMath.h"
 #include "TMatrixD.h"
 #include "TRandom.h"
@@ -60,7 +61,7 @@ void FastTracker::AddLayer(TString name, float r, float z, float x0, float xrho,
   layers.push_back(newLayer);
 }
 
-DetLayer FastTracker::GetLayer(int layer, bool ignoreBarrelLayers)
+DetLayer FastTracker::GetLayer(int layer, bool ignoreBarrelLayers) const
 {
   int layerIdx = layer;
   if (ignoreBarrelLayers) {
@@ -74,6 +75,18 @@ DetLayer FastTracker::GetLayer(int layer, bool ignoreBarrelLayers)
   return layers[layerIdx];
 }
 
+int FastTracker::GetLayerIndex(std::string name) const
+{
+  int i = 0;
+  for (const auto& layer : layers) {
+    if (layer.name == name) {
+      return i;
+    }
+    i++;
+  }
+  return -1;
+}
+
 void FastTracker::Print()
 {
   // print out layer setup
@@ -100,18 +113,18 @@ void FastTracker::AddSiliconALICE3v4(std::vector<float> pixelResolution)
   float resRPhiOT = pixelResolution[2];
   float resZOT = pixelResolution[3];
 
-  layers.push_back(DetLayer{"bpipe0", 0.48, 250, 0.00042, 2.772e-02, 0.0f, 0.0f, 0.0f, 0}); // 150 mum Be
-  layers.push_back(DetLayer{"ddd0", 0.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1});
-  layers.push_back(DetLayer{"ddd1", 1.2, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1});
-  layers.push_back(DetLayer{"ddd2", 2.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1});
-  layers.push_back(DetLayer{"bpipe1", 5.7, 250, 0.0014, 9.24e-02, 0.0f, 0.0f, 0.0f, 0}); // 500 mum Be
-  layers.push_back(DetLayer{"ddd3", 7., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"ddd4", 10., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"ddd5", 13., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"ddd6", 16., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"ddd7", 25., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"ddd8", 40., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"ddd9", 45., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
+  AddLayer("bpipe0", 0.48, 250, 0.00042, 2.772e-02, 0.0f, 0.0f, 0.0f, 0); // 150 mum Be
+  AddLayer("ddd0", 0.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1);
+  AddLayer("ddd1", 1.2, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1);
+  AddLayer("ddd2", 2.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1);
+  AddLayer("bpipe1", 5.7, 250, 0.0014, 9.24e-02, 0.0f, 0.0f, 0.0f, 0); // 500 mum Be
+  AddLayer("ddd3", 7., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("ddd4", 10., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("ddd5", 13., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("ddd6", 16., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("ddd7", 25., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("ddd8", 40., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("ddd9", 45., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
 }
 
 void FastTracker::AddSiliconALICE3v2(std::vector<float> pixelResolution)
@@ -128,19 +141,19 @@ void FastTracker::AddSiliconALICE3v2(std::vector<float> pixelResolution)
   float resRPhiOT = pixelResolution[2];
   float resZOT = pixelResolution[3];
 
-  layers.push_back(DetLayer{"bpipe0", 0.48, 250, 0.00042, 2.772e-02, 0.0f, 0.0f, 0.0f, 0}); // 150 mum Be
-  layers.push_back(DetLayer{"B00", 0.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1});
-  layers.push_back(DetLayer{"B01", 1.2, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1});
-  layers.push_back(DetLayer{"B02", 2.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1});
-  layers.push_back(DetLayer{"bpipe1", 3.7, 250, 0.0014, 9.24e-02, 0.0f, 0.0f, 0.0f, 0}); // 500 mum Be
-  layers.push_back(DetLayer{"B03", 3.75, 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B04", 7., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B05", 12., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B06", 20., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B07", 30., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B08", 45., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B09", 60., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
-  layers.push_back(DetLayer{"B10", 80., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1});
+  AddLayer("bpipe0", 0.48, 250, 0.00042, 2.772e-02, 0.0f, 0.0f, 0.0f, 0); // 150 mum Be
+  AddLayer("B00", 0.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1);
+  AddLayer("B01", 1.2, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1);
+  AddLayer("B02", 2.5, 250, x0IT, xrhoIB, resRPhiIT, resZIT, eff, 1);
+  AddLayer("bpipe1", 3.7, 250, 0.0014, 9.24e-02, 0.0f, 0.0f, 0.0f, 0); // 500 mum Be
+  AddLayer("B03", 3.75, 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B04", 7., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B05", 12., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B06", 20., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B07", 30., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B08", 45., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B09", 60., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
+  AddLayer("B10", 80., 250, x0OT, xrhoOT, resRPhiOT, resZOT, eff, 1);
 }
 
 void FastTracker::AddTPC(float phiResMean, float zResMean)
@@ -202,7 +215,7 @@ float FastTracker::Dist(float z, float r)
     if (i == nSteps - 1)
       index = 1;
     z0 = -4 * sigmaD + i * dz0;
-    dist += index * (dz0 / 3.) * (1 / o2::math_utils::sqrt(o2::constants::math::TwoPI) / sigmaD) * exp(-z0 * z0 / 2. / sigmaD / sigmaD) * (1 / o2::math_utils::sqrt((z - z0) * (z - z0) + r * r));
+    dist += index * (dz0 / 3.) * (1 / o2::math_utils::sqrt(o2::constants::math::TwoPI) / sigmaD) * std::exp(-z0 * z0 / 2. / sigmaD / sigmaD) * (1 / o2::math_utils::sqrt((z - z0) * (z - z0) + r * r));
     if (index != 4)
       index = 4;
     else
@@ -283,7 +296,7 @@ float FastTracker::ProbGoodChiSqHit(float radius, float searchRadiusRPhi, float
 
 // function to provide a reconstructed track from a perfect input track
 // returns number of intercepts (generic for now)
-int FastTracker::FastTrack(o2::track::TrackParCov inputTrack, o2::track::TrackParCov& outputTrack, float nch)
+int FastTracker::FastTrack(o2::track::TrackParCov inputTrack, o2::track::TrackParCov& outputTrack, const float nch)
 {
   hits.clear();
   nIntercepts = 0;
@@ -333,7 +346,6 @@ int FastTracker::FastTrack(o2::track::TrackParCov inputTrack, o2::track::TrackPa
     // was there a problem on this layer?
     if (!ok && il > 0) { // may fail to reach target layer due to the eloss
       float rad2 = inputTrack.getX() * inputTrack.getX() + inputTrack.getY() * inputTrack.getY();
-      float fMinRadTrack = 132.;
       float maxR = layers[il - 1].r + kTrackingMargin * 2;
       float minRad = (fMinRadTrack > 0 && fMinRadTrack < maxR) ? fMinRadTrack : maxR;
       if (rad2 - minRad * minRad < kTrackingMargin * kTrackingMargin) { // check previously reached layer

ALICE3/Core/FastTracker.cxx

@@ -14,6 +14,7 @@
 
 #include <fairlogger/Logger.h> // not a system header but megalinter thinks so
 #include <vector>
+#include <string>
 #include "DetLayer.h"
 #include "ReconstructionDataFormats/Track.h"
 
@@ -34,15 +35,38 @@ class FastTracker
   FastTracker();
   virtual ~FastTracker() {}
 
+  // Layer and layer configuration
   void AddLayer(TString name, float r, float z, float x0, float xrho, float resRPhi = 0.0f, float resZ = 0.0f, float eff = 0.0f, int type = 0);
-  DetLayer GetLayer(const int layer, bool ignoreBarrelLayers = true);
+  DetLayer GetLayer(const int layer, bool ignoreBarrelLayers = true) const;
+  int GetLayerIndex(const std::string name) const;
+  void SetRadiationLength(const std::string layerName, float x0) { layers[GetLayerIndex(layerName)].x0 = x0; }
+  void SetRadius(const std::string layerName, float r) { layers[GetLayerIndex(layerName)].r = r; }
+  void SetResolutionRPhi(const std::string layerName, float resRPhi) { layers[GetLayerIndex(layerName)].resRPhi = resRPhi; }
+  void SetResolutionZ(const std::string layerName, float resZ) { layers[GetLayerIndex(layerName)].resZ = resZ; }
+  void SetResolution(const std::string layerName, float resRPhi, float resZ)
+  {
+    SetResolutionRPhi(layerName, resRPhi);
+    SetResolutionZ(layerName, resZ);
+  }
 
   void AddSiliconALICE3v4(std::vector<float> pixelResolution);
   void AddSiliconALICE3v2(std::vector<float> pixelResolution);
   void AddTPC(float phiResMean, float zResMean);
 
   void Print();
-  int FastTrack(o2::track::TrackParCov inputTrack, o2::track::TrackParCov& outputTrack, float nch);
+
+  /**
+   * @brief Performs fast tracking on the input track parameters.
+   *
+   * Propagates the given input track through the detector layers, applying
+   * relevant corrections and updates, and stores the result in outputTrack.
+   *
+   * @param inputTrack The input track parameters and covariance (const, by value).
+   * @param outputTrack Reference to the output track parameters and covariance, to be filled.
+   * @param nch Charged particle multiplicity (used for hit density calculations).
+   * @return int i.e. number of intercepts (implementation-defined).
+   */
+  int FastTrack(o2::track::TrackParCov inputTrack, o2::track::TrackParCov& outputTrack, const float nch);
 
   // For efficiency calculation
   float Dist(float z, float radius);
@@ -52,6 +76,32 @@ class FastTracker
   float HitDensity(float radius);
   float ProbGoodChiSqHit(float radius, float searchRadiusRPhi, float searchRadiusZ);
 
+  // Setters and getters for configuration
+  void SetIntegrationTime(float t) { integrationTime = t; }
+  void SetMaxRadiusOfSlowDetectors(float r) { maxRadiusSlowDet = r; }
+  void SetAvgRapidity(float y) { avgRapidity = y; }
+  void SetdNdEtaCent(float d) { dNdEtaCent = d; }
+  void SetLhcUPCscale(float s) { lhcUPCScale = s; }
+  void SetBField(float b) { magneticField = b; }
+  void SetMinRadTrack(float r) { fMinRadTrack = r; }
+  void SetMagneticField(float b) { magneticField = b; }
+  void SetApplyZacceptance(bool b) { applyZacceptance = b; }
+  void SetApplyMSCorrection(bool b) { applyMSCorrection = b; }
+  void SetApplyElossCorrection(bool b) { applyElossCorrection = b; }
+
+  // Getters for the last track
+  int GetNIntercepts() const { return nIntercepts; }
+  int GetNSiliconPoints() const { return nSiliconPoints; }
+  int GetNGasPoints() const { return nGasPoints; }
+  float GetGoodHitProb(int layer) const { return goodHitProbability[layer]; }
+  std::size_t GetNHits() const { return hits.size(); }
+  float GetHitX(const int i) const { return hits[i][0]; }
+  float GetHitY(const int i) const { return hits[i][1]; }
+  float GetHitZ(const int i) const { return hits[i][2]; }
+  uint64_t GetCovMatOK() const { return covMatOK; }
+  uint64_t GetCovMatNotOK() const { return covMatNotOK; }
+
+ private:
   // Definition of detector layers
   std::vector<DetLayer> layers;
   std::vector<std::vector<float>> hits; // bookkeep last added hits
@@ -75,6 +125,7 @@ class FastTracker
   float avgRapidity;
   float lhcUPCScale;
   float upcBackgroundMultiplier;
+  float fMinRadTrack = 132.;
 
   uint64_t covMatOK;    // cov mat has negative eigenvals
   uint64_t covMatNotOK; // cov mat has negative eigenvals

ALICE3/Core/FastTracker.h

@@ -0,0 +1,107 @@
+// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
+// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
+// All rights not expressly granted are reserved.
+//
+// This software is distributed under the terms of the GNU General Public
+// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
+//
+// In applying this license CERN does not waive the privileges and immunities
+// granted to it by virtue of its status as an Intergovernmental Organization
+// or submit itself to any jurisdiction.
+///
+/// \file TrackUtilities.h
+///
+/// \brief Set of utilities for the ALICE3 track handling
+///
+/// \since  May 21, 2025
+///
+
+#ifndef ALICE3_CORE_TRACKUTILITIES_H_
+#define ALICE3_CORE_TRACKUTILITIES_H_
+
+#include <vector>
+
+#include "ReconstructionDataFormats/Track.h"
+#include "Framework/O2DatabasePDGPlugin.h"
+#include "Framework/AnalysisHelpers.h"
+#include "TLorentzVector.h"
+
+namespace o2::upgrade
+{
+
+/// Function to convert a TLorentzVector into a perfect Track
+/// \param charge particle charge (integer)
+/// \param particle the particle to convert (TLorentzVector)
+/// \param productionVertex where the particle was produced
+/// \param o2track the address of the resulting TrackParCov
+void convertTLorentzVectorToO2Track(const int charge,
+                                    const TLorentzVector particle,
+                                    const std::vector<double> productionVertex,
+                                    o2::track::TrackParCov& o2track)
+{
+  std::array<float, 5> params;
+  std::array<float, 15> covm = {0.};
+  float s, c, x;
+  o2::math_utils::sincos(static_cast<float>(particle.Phi()), s, c);
+  o2::math_utils::rotateZInv(static_cast<float>(productionVertex[0]), static_cast<float>(productionVertex[1]), x, params[0], s, c);
+  params[1] = static_cast<float>(productionVertex[2]);
+  params[2] = 0.; // since alpha = phi
+  const auto theta = 2. * std::atan(std::exp(-particle.PseudoRapidity()));
+  params[3] = 1. / std::tan(theta);
+  params[4] = charge / particle.Pt();
+
+  // Initialize TrackParCov in-place
+  new (&o2track)(o2::track::TrackParCov)(x, particle.Phi(), params, covm);
+}
+
+/// Function to convert a TLorentzVector into a perfect Track
+/// \param pdgCode particle pdg
+/// \param particle the particle to convert (TLorentzVector)
+/// \param productionVertex where the particle was produced
+/// \param o2track the address of the resulting TrackParCov
+/// \param pdg the pdg service
+void convertTLorentzVectorToO2Track(int pdgCode,
+                                    TLorentzVector particle,
+                                    std::vector<double> productionVertex,
+                                    o2::track::TrackParCov& o2track,
+                                    const o2::framework::Service<o2::framework::O2DatabasePDG>& pdg)
+{
+  const auto pdgInfo = pdg->GetParticle(pdgCode);
+  int charge = 0;
+  if (pdgInfo != nullptr) {
+    charge = pdgInfo->Charge() / 3;
+  }
+  convertTLorentzVectorToO2Track(charge, particle, productionVertex, o2track);
+}
+
+/// Function to convert a McParticle into a perfect Track
+/// \param particle the particle to convert (mcParticle)
+/// \param o2track the address of the resulting TrackParCov
+/// \param pdg the pdg service
+template <typename McParticleType>
+void convertMCParticleToO2Track(McParticleType& particle,
+                                o2::track::TrackParCov& o2track,
+                                const o2::framework::Service<o2::framework::O2DatabasePDG>& pdg)
+{
+  static TLorentzVector tlv;
+  tlv.SetPxPyPzE(particle.px(), particle.py(), particle.pz(), particle.e());
+  tlv.SetXYZT(particle.vx(), particle.vy(), particle.vz(), particle.vt());
+  convertTLorentzVectorToO2Track(particle.pdgCode(), tlv, {particle.vx(), particle.vy(), particle.vz()}, o2track, pdg);
+}
+
+/// Function to convert a McParticle into a perfect Track
+/// \param particle the particle to convert (mcParticle)
+/// \param o2track the address of the resulting TrackParCov
+/// \param pdg the pdg service
+template <typename McParticleType>
+o2::track::TrackParCov convertMCParticleToO2Track(McParticleType& particle,
+                                                  const o2::framework::Service<o2::framework::O2DatabasePDG>& pdg)
+{
+  o2::track::TrackParCov o2track;
+  convertMCParticleToO2Track(particle, o2track, pdg);
+  return o2track;
+}
+
+} // namespace o2::upgrade
+
+#endif // ALICE3_CORE_TRACKUTILITIES_H_