ITS: simplify constants + mathutils (#14383)

Author: f3sch

Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TrackingKernels.h

@@ -22,7 +22,11 @@ class CellSeed;
 class ExternalAllocator;
 namespace gpu
 {
+
 #ifdef GPUCA_GPUCODE // GPUg() global kernels must only when compiled by GPU compiler
+
+GPUdi() int4 getEmptyBinsRect() { return int4{0, 0, 0, 0}; }
+
 GPUd() bool fitTrack(TrackITSExt& track,
                      int start,
                      int end,

Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/VertexerTraitsGPU.h

@@ -18,7 +18,6 @@
 #define ITSTRACKINGGPU_VERTEXERTRAITSGPU_H_
 
 #include <vector>
-#include <array>
 
 #include "ITStracking/VertexerTraits.h"
 #include "ITStracking/Configuration.h"
@@ -29,13 +28,8 @@
 
 #include "ITStrackingGPU/TimeFrameGPU.h"
 
-namespace o2
+namespace o2::its
 {
-namespace its
-{
-class ROframe;
-
-using constants::its2::InversePhiBinSize;
 
 class VertexerTraitsGPU final : public VertexerTraits
 {
@@ -63,6 +57,6 @@ inline void VertexerTraitsGPU::adoptTimeFrame(TimeFrame<7>* tf) noexcept
   mTimeFrame = static_cast<TimeFrame<7>*>(tf);
 }
 
-} // namespace its
-} // namespace o2
+} // namespace o2::its
+
 #endif

Detectors/ITSMFT/ITS/tracking/GPU/cuda/TrackingKernels.cu

@@ -52,14 +52,8 @@ using namespace o2::track;
 
 namespace o2::its
 {
-using namespace constants::its2;
 using Vertex = o2::dataformats::Vertex<o2::dataformats::TimeStamp<int>>;
 
-GPUdii() float Sq(float v)
-{
-  return v * v;
-}
-
 namespace gpu
 {
 
@@ -99,9 +93,9 @@ GPUd() const int4 getBinsRect(const Cluster& currentCluster, const int layerInde
                               const float z1, const float z2, float maxdeltaz, float maxdeltaphi)
 {
   const float zRangeMin = o2::gpu::CAMath::Min(z1, z2) - maxdeltaz;
-  const float phiRangeMin = (maxdeltaphi > constants::math::Pi) ? 0.f : currentCluster.phi - maxdeltaphi;
+  const float phiRangeMin = (maxdeltaphi > o2::constants::math::PI) ? 0.f : currentCluster.phi - maxdeltaphi;
   const float zRangeMax = o2::gpu::CAMath::Max(z1, z2) + maxdeltaz;
-  const float phiRangeMax = (maxdeltaphi > constants::math::Pi) ? constants::math::TwoPi : currentCluster.phi + maxdeltaphi;
+  const float phiRangeMax = (maxdeltaphi > o2::constants::math::PI) ? o2::constants::math::TwoPI : currentCluster.phi + maxdeltaphi;
 
   if (zRangeMax < -utils.getLayerZ(layerIndex) ||
       zRangeMin > utils.getLayerZ(layerIndex) || zRangeMin > zRangeMax) {
@@ -129,7 +123,7 @@ GPUd() bool fitTrack(TrackITSExt& track,
                      o2::base::PropagatorF::MatCorrType matCorrType)
 {
   for (int iLayer{start}; iLayer != end; iLayer += step) {
-    if (track.getClusterIndex(iLayer) == constants::its::UnusedIndex) {
+    if (track.getClusterIndex(iLayer) == constants::UnusedIndex) {
       continue;
     }
     const TrackingFrameInfo& trackingHit = tfInfos[iLayer][track.getClusterIndex(iLayer)];
@@ -316,7 +310,7 @@ GPUg() void fitTrackSeedsKernel(
     temporaryTrack.setChi2(0);
     int* clusters = seed.getClusters();
     for (int iL{0}; iL < 7; ++iL) {
-      temporaryTrack.setExternalClusterIndex(iL, clusters[iL], clusters[iL] != constants::its::UnusedIndex);
+      temporaryTrack.setExternalClusterIndex(iL, clusters[iL], clusters[iL] != constants::UnusedIndex);
     }
     bool fitSuccess = fitTrack(temporaryTrack,               // TrackITSExt& track,
                                0,                            // int lastLayer,
@@ -422,8 +416,6 @@ GPUg() void computeLayerCellsKernel(
   const float cellDeltaTanLambdaSigma,
   const float nSigmaCut)
 {
-  constexpr float radl = 9.36f;                                                           // Radiation length of Si [cm].
-  constexpr float rho = 2.33f;                                                            // Density of Si [g/cm^3].
   constexpr float layerxX0[7] = {5.e-3f, 5.e-3f, 5.e-3f, 1.e-2f, 1.e-2f, 1.e-2f, 1.e-2f}; // Hardcoded here for the moment.
   for (int iCurrentTrackletIndex = blockIdx.x * blockDim.x + threadIdx.x; iCurrentTrackletIndex < nTrackletsCurrent; iCurrentTrackletIndex += blockDim.x * gridDim.x) {
     const Tracklet& currentTracklet = tracklets[layer][iCurrentTrackletIndex];
@@ -462,7 +454,7 @@ GPUg() void computeLayerCellsKernel(
             break;
           }
 
-          if (!track.correctForMaterial(layerxX0[layer + iC], layerxX0[layer] * radl * rho, true)) {
+          if (!track.correctForMaterial(layerxX0[layer + iC], layerxX0[layer] * constants::Radl * constants::Rho, true)) {
             break;
           }
 
@@ -548,12 +540,12 @@ GPUg() void computeLayerTrackletsMultiROFKernel(
         if (primaryVertex.isFlagSet(2) && iteration != 3) {
           continue;
         }
-        const float resolution = o2::gpu::CAMath::Sqrt(Sq(resolutionPV) / primaryVertex.getNContributors() + Sq(positionResolution));
+        const float resolution = o2::gpu::CAMath::Sqrt(math_utils::Sq(resolutionPV) / primaryVertex.getNContributors() + math_utils::Sq(positionResolution));
         const float tanLambda{(currentCluster.zCoordinate - primaryVertex.getZ()) * inverseR0};
         const float zAtRmin{tanLambda * (minR - currentCluster.radius) + currentCluster.zCoordinate};
         const float zAtRmax{tanLambda * (maxR - currentCluster.radius) + currentCluster.zCoordinate};
-        const float sqInverseDeltaZ0{1.f / (Sq(currentCluster.zCoordinate - primaryVertex.getZ()) + 2.e-8f)}; /// protecting from overflows adding the detector resolution
-        const float sigmaZ{o2::gpu::CAMath::Sqrt(Sq(resolution) * Sq(tanLambda) * ((Sq(inverseR0) + sqInverseDeltaZ0) * Sq(meanDeltaR) + 1.f) + Sq(meanDeltaR * MSAngle))};
+        const float sqInverseDeltaZ0{1.f / (math_utils::Sq(currentCluster.zCoordinate - primaryVertex.getZ()) + 2.e-8f)}; /// protecting from overflows adding the detector resolution
+        const float sigmaZ{o2::gpu::CAMath::Sqrt(math_utils::Sq(resolution) * math_utils::Sq(tanLambda) * ((math_utils::Sq(inverseR0) + sqInverseDeltaZ0) * math_utils::Sq(meanDeltaR) + 1.f) + math_utils::Sq(meanDeltaR * MSAngle))};
         const int4 selectedBinsRect{getBinsRect(currentCluster, layerIndex + 1, *utils, zAtRmin, zAtRmax, sigmaZ * NSigmaCut, phiCut)};
         if (selectedBinsRect.x == 0 && selectedBinsRect.y == 0 && selectedBinsRect.z == 0 && selectedBinsRect.w == 0) {
           continue;
@@ -587,7 +579,7 @@ GPUg() void computeLayerTrackletsMultiROFKernel(
               const float deltaPhi{o2::gpu::CAMath::Abs(currentCluster.phi - nextCluster.phi)};
               const float deltaZ{o2::gpu::CAMath::Abs(tanLambda * (nextCluster.radius - currentCluster.radius) + currentCluster.zCoordinate - nextCluster.zCoordinate)};
               const int nextSortedIndex{ROFClusters[layerIndex + 1][rof1] + nextClusterIndex};
-              if (deltaZ / sigmaZ < NSigmaCut && (deltaPhi < phiCut || o2::gpu::CAMath::Abs(deltaPhi - constants::math::TwoPi) < phiCut)) {
+              if (deltaZ / sigmaZ < NSigmaCut && (deltaPhi < phiCut || o2::gpu::CAMath::Abs(deltaPhi - o2::constants::math::TwoPI) < phiCut)) {
                 if constexpr (initRun) {
                   trackletsLUT[layerIndex][currentSortedIndex]++; // we need l0 as well for usual exclusive sums.
                 } else {
@@ -634,8 +626,6 @@ GPUg() void processNeighboursKernel(const int layer,
                                     const o2::base::Propagator* propagator,
                                     const o2::base::PropagatorF::MatCorrType matCorrType)
 {
-  constexpr float radl = 9.36f;                                                           // Radiation length of Si [cm].
-  constexpr float rho = 2.33f;                                                            // Density of Si [g/cm^3].
   constexpr float layerxX0[7] = {5.e-3f, 5.e-3f, 5.e-3f, 1.e-2f, 1.e-2f, 1.e-2f, 1.e-2f}; // Hardcoded here for the moment.
   for (unsigned int iCurrentCell = blockIdx.x * blockDim.x + threadIdx.x; iCurrentCell < nCurrentCells; iCurrentCell += blockDim.x * gridDim.x) {
     int foundSeeds{0};
@@ -678,7 +668,7 @@ GPUg() void processNeighboursKernel(const int layer,
       }
 
       if (matCorrType == o2::base::PropagatorF::MatCorrType::USEMatCorrNONE) {
-        if (!seed.correctForMaterial(layerxX0[layer - 1], layerxX0[layer - 1] * radl * rho, true)) {
+        if (!seed.correctForMaterial(layerxX0[layer - 1], layerxX0[layer - 1] * constants::Radl * constants::Rho, true)) {
           continue;
         }
       }

Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cxx

@@ -40,7 +40,7 @@ void VertexerTraitsGPU::updateVertexingParameters(const std::vector<VertexingPar
   mTfGPUParams = tfPar;
   mIndexTableUtils.setTrackingParameters(vrtPar[0]);
   for (auto& par : mVrtParams) {
-    par.phiSpan = static_cast<int>(std::ceil(mIndexTableUtils.getNphiBins() * par.phiCut / constants::math::TwoPi));
+    par.phiSpan = static_cast<int>(std::ceil(mIndexTableUtils.getNphiBins() * par.phiCut / o2::constants::math::TwoPI));
     par.zSpan = static_cast<int>(std::ceil(par.zCut * mIndexTableUtils.getInverseZCoordinate(0)));
   }
 }

Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexingKernels.cu

@@ -20,10 +20,7 @@ namespace o2
 {
 namespace its
 {
-using constants::its::VertexerHistogramVolume;
-using constants::math::TwoPi;
 using math_utils::getNormalizedPhi;
-using namespace constants::its2;
 
 namespace gpu
 {
@@ -58,11 +55,6 @@ void trackletFinderHandler(const Cluster* clustersNextLayer,    // 0 2
     maxTrackletsPerCluster); // const unsigned int maxTrackletsPerCluster = 1e2
 }
 /*
-GPUd() float smallestAngleDifference(float a, float b)
-{
-  float diff = fmod(b - a + constants::math::Pi, constants::math::TwoPi) - constants::math::Pi;
-  return (diff < -constants::math::Pi) ? diff + constants::math::TwoPi : ((diff > constants::math::Pi) ? diff - constants::math::TwoPi : diff);
-}
 
 GPUd() const int4 getBinsRect(const Cluster& currentCluster, const int layerIndex,
                               const float z1, float maxdeltaz, float maxdeltaphi)

Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h

@@ -58,9 +58,9 @@ class Configuration : public Param
 };
 
 struct TrackingParameters {
-  int CellMinimumLevel() { return MinTrackLength - constants::its::ClustersPerCell + 1; }
-  int CellsPerRoad() const { return NLayers - 2; }
-  int TrackletsPerRoad() const { return NLayers - 1; }
+  int CellMinimumLevel() const noexcept { return MinTrackLength - constants::ClustersPerCell + 1; }
+  int CellsPerRoad() const noexcept { return NLayers - 2; }
+  int TrackletsPerRoad() const noexcept { return NLayers - 1; }
   std::string asString() const;
 
   int NLayers = 7;

Detectors/ITSMFT/ITS/tracking/include/ITStracking/Constants.h

@@ -17,112 +17,25 @@
 #define TRACKINGITSU_INCLUDE_CONSTANTS_H_
 
 #include "ITStracking/Definitions.h"
-#include "CommonConstants/MathConstants.h"
 
-#include "GPUCommonMath.h"
-#include "GPUCommonDef.h"
-
-#ifndef GPUCA_GPUCODE_DEVICE
-#include <climits>
-#include <vector>
-#include <array>
-#endif
-
-namespace o2
-{
-namespace its
-{
-
-namespace constants
+namespace o2::its::constants
 {
 constexpr float MB = 1024.f * 1024.f;
 constexpr float GB = 1024.f * 1024.f * 1024.f;
 constexpr bool DoTimeBenchmarks = true;
 constexpr bool SaveTimeBenchmarks = false;
 
-namespace math
-{
-constexpr float Pi{3.14159265359f};
-constexpr float TwoPi{2.0f * Pi};
-constexpr float FloatMinThreshold{1e-20f};
-} // namespace math
-
-namespace its
-{
-constexpr int LayersNumberVertexer{3};
+constexpr float Tolerance{1e-12}; // numerical tolerance
 constexpr int ClustersPerCell{3};
 constexpr int UnusedIndex{-1};
 constexpr float Resolution{0.0005f};
-
-GPUhdi() constexpr std::array<float, 3> VertexerHistogramVolume()
-{
-  return std::array<float, 3>{{1.98, 1.98, 40.f}};
-}
-} // namespace its
-
-namespace its2
-{
-constexpr int LayersNumber{7};
-constexpr int TrackletsPerRoad{LayersNumber - 1};
-constexpr int CellsPerRoad{LayersNumber - 2};
-
-GPUhdi() constexpr std::array<float, LayersNumber> LayersZCoordinate()
-{
-  constexpr double s = 1.; // safety margin
-  return std::array<float, LayersNumber>{16.333f + s, 16.333f + s, 16.333f + s, 42.140f + s, 42.140f + s, 73.745f + s, 73.745f + s};
-}
-
-GPUhdi() constexpr std::array<float, LayersNumber> LayersRCoordinate()
+constexpr float Radl = 9.36f; // Radiation length of Si [cm]
+constexpr float Rho = 2.33f;  // Density of Si [g/cm^3]
+namespace its                 // to be removed
 {
-  return std::array<float, LayersNumber>{{2.33959f, 3.14076f, 3.91924f, 19.6213f, 24.5597f, 34.388f, 39.3329f}};
-}
-
-constexpr int ZBins{256};
-constexpr int PhiBins{128};
-constexpr float InversePhiBinSize{PhiBins / constants::math::TwoPi};
-GPUhdi() constexpr std::array<float, LayersNumber> InverseZBinSize()
-{
-  constexpr auto zSize = LayersZCoordinate();
-  return std::array<float, LayersNumber>{0.5f * ZBins / (zSize[0]), 0.5f * ZBins / (zSize[1]), 0.5f * ZBins / (zSize[2]),
-                                         0.5f * ZBins / (zSize[3]), 0.5f * ZBins / (zSize[4]), 0.5f * ZBins / (zSize[5]),
-                                         0.5f * ZBins / (zSize[6])};
-}
-
-GPUhdi() constexpr float getInverseZCoordinate(const int layerIndex)
-{
-  return 0.5f * ZBins / LayersZCoordinate()[layerIndex];
-}
-
-GPUhdi() int getZBinIndex(const int layerIndex, const float zCoordinate)
-{
-  return (zCoordinate + LayersZCoordinate()[layerIndex]) *
-         InverseZBinSize()[layerIndex];
-}
-
-GPUhdi() int getPhiBinIndex(const float currentPhi)
-{
-  return (currentPhi * InversePhiBinSize);
-}
-
-GPUhdi() int getBinIndex(const int zIndex, const int phiIndex)
-{
-  return o2::gpu::GPUCommonMath::Min(phiIndex * ZBins + zIndex,
-                                     ZBins * PhiBins - 1);
-}
-
-GPUhdi() constexpr int4 getEmptyBinsRect() { return int4{0, 0, 0, 0}; }
-
-} // namespace its2
-
-namespace pdgcodes
-{
-constexpr int PionCode{211};
-}
-} // namespace constants
-#ifndef GPUCA_GPUCODE_DEVICE
-typedef std::vector<std::vector<int>> index_table_t;
-#endif
+constexpr int UnusedIndex{-1};
+constexpr float Resolution{0.0005f};
 } // namespace its
-} // namespace o2
+} // namespace o2::its::constants
 
 #endif /* TRACKINGITSU_INCLUDE_CONSTANTS_H_ */

Detectors/ITSMFT/ITS/tracking/include/ITStracking/IndexTableUtils.h

@@ -16,9 +16,9 @@
 #ifndef TRACKINGITSU_INCLUDE_INDEXTABLEUTILS_H_
 #define TRACKINGITSU_INCLUDE_INDEXTABLEUTILS_H_
 
-#include "ITStracking/Constants.h"
 #include "ITStracking/Configuration.h"
 #include "ITStracking/Definitions.h"
+#include "CommonConstants/MathConstants.h"
 #include "GPUCommonMath.h"
 #include "GPUCommonDef.h"
 
@@ -55,7 +55,7 @@ class IndexTableUtils
 template <class T>
 inline void IndexTableUtils::setTrackingParameters(const T& params)
 {
-  mInversePhiBinSize = params.PhiBins / constants::math::TwoPi;
+  mInversePhiBinSize = params.PhiBins / o2::constants::math::TwoPI;
   mNzBins = params.ZBins;
   mNphiBins = params.PhiBins;
   for (int iLayer{0}; iLayer < params.LayerZ.size(); ++iLayer) {