Tracker Traits and Enabling Phase2 for Inner Tracker Reconstruction on GPU

CUDADataFormats/Common/interface/HeterogeneousSoA.h

@@ -92,6 +92,11 @@ namespace cms {
         return cms::cuda::make_host_unique<T>(stream);
       }
 
+      template <typename T>
+      static auto make_unique(size_t size, cudaStream_t stream) {
+        return cms::cuda::make_host_unique<T>(size, stream);
+      }
+
       template <typename T>
       static auto make_host_unique(cudaStream_t stream) {
         return cms::cuda::make_host_unique<T>(stream);

CUDADataFormats/SiPixelDigi/interface/SiPixelDigisCUDASOAView.h

@@ -12,7 +12,10 @@
 class SiPixelDigisCUDASOAView {
 public:
   friend class SiPixelDigisCUDA;
-  friend class SiPixelRecHitSoAFromLegacy;
+
+  template <typename TrackerTraits>
+  friend class SiPixelRecHitSoAFromLegacyT;
+
   enum class StorageLocation {
     kCLUS = 0,
     kPDIGI = 2,

CUDADataFormats/Track/interface/PixelTrackHeterogeneous.h

@@ -3,7 +3,9 @@
 
 #include "CUDADataFormats/Common/interface/HeterogeneousSoA.h"
 #include "CUDADataFormats/Track/interface/TrackSoAHeterogeneousT.h"
+#include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
 
-using PixelTrackHeterogeneous = HeterogeneousSoA<pixelTrack::TrackSoA>;
+template <typename TrackerTraits>
+using PixelTrackHeterogeneousT = HeterogeneousSoA<pixelTrack::TrackSoAT<TrackerTraits>>;
 
-#endif  // #ifndef CUDADataFormats_Track_PixelTrackHeterogeneous_h
+#endif  // #ifndef CUDADataFormats_Track_PixelTrackHeterogeneous_h

CUDADataFormats/Track/interface/TrackSoAHeterogeneousT.h

@@ -5,33 +5,38 @@
 #include <algorithm>
 
 #include "CUDADataFormats/Track/interface/TrajectoryStateSoAT.h"
-#include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
 #include "HeterogeneousCore/CUDAUtilities/interface/HistoContainer.h"
-
+#include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
 #include "CUDADataFormats/Common/interface/HeterogeneousSoA.h"
+#include "DataFormats/Common/interface/CMS_CLASS_VERSION.h"
 
 namespace pixelTrack {
+
   enum class Quality : uint8_t { bad = 0, edup, dup, loose, strict, tight, highPurity, notQuality };
   constexpr uint32_t qualitySize{uint8_t(Quality::notQuality)};
   const std::string qualityName[qualitySize]{"bad", "edup", "dup", "loose", "strict", "tight", "highPurity"};
   inline Quality qualityByName(std::string const &name) {
     auto qp = std::find(qualityName, qualityName + qualitySize, name) - qualityName;
     return static_cast<Quality>(qp);
   }
+
 }  // namespace pixelTrack
 
-template <int32_t S>
+template <typename TrackerTraits>
 class TrackSoAHeterogeneousT {
 public:
+  static constexpr int32_t S = TrackerTraits::maxNumberOfTuples;
+  static constexpr int32_t H = TrackerTraits::maxHitsOnTrack;  // Average hits rather than max?
   static constexpr int32_t stride() { return S; }
 
+  using hindex_type = uint32_t;  //TrackerTraits::hindex_type ?
+
   using Quality = pixelTrack::Quality;
-  using hindex_type = uint32_t;
-  using HitContainer = cms::cuda::OneToManyAssoc<hindex_type, S + 1, 5 * S>;
+  using HitContainer = cms::cuda::OneToManyAssoc<hindex_type, S + 1, H * S>;
 
   // Always check quality is at least loose!
   // CUDA does not support enums  in __lgc ...
-private:
+protected:
   eigenSoA::ScalarSoA<uint8_t, S> quality_;
 
 public:
@@ -56,9 +61,9 @@ class TrackSoAHeterogeneousT {
     // layers are in order and we assume tracks are either forward or backward
     auto pdet = detIndices.begin(i);
     int nl = 1;
-    auto ol = phase1PixelTopology::getLayer(*pdet);
+    auto ol = pixelTopology::getLayer<TrackerTraits>(*pdet);
     for (; pdet < detIndices.end(i); ++pdet) {
-      auto il = phase1PixelTopology::getLayer(*pdet);
+      auto il = pixelTopology::getLayer<TrackerTraits>(*pdet);
       if (il != ol)
         ++nl;
       ol = il;
@@ -90,17 +95,100 @@ class TrackSoAHeterogeneousT {
 
 namespace pixelTrack {
 
-#ifdef GPU_SMALL_EVENTS
-  // kept for testing and debugging
-  constexpr uint32_t maxNumber() { return 2 * 1024; }
-#else
-  // tested on MC events with 55-75 pileup events
-  constexpr uint32_t maxNumber() { return 32 * 1024; }
-#endif
+  template <typename TrackerTraits>
+  using TrackSoAT = TrackSoAHeterogeneousT<TrackerTraits>;
+
+  template <typename TrackerTraits>
+  using HitContainerT = typename TrackSoAHeterogeneousT<TrackerTraits>::HitContainer;
+
+  //Used only to ease classes definitions
+  using TrackSoAPhase1 = TrackSoAHeterogeneousT<pixelTopology::Phase1>;
+  using TrackSoAPhase2 = TrackSoAHeterogeneousT<pixelTopology::Phase2>;
+
+  template <typename TrackerTraits, typename Enable = void>
+  struct QualityCutsT {};
+
+  template <typename TrackerTraits>
+  struct QualityCutsT<TrackerTraits, pixelTopology::isPhase1Topology<TrackerTraits>> {
+    // chi2 cut = chi2Scale * (chi2Coeff[0] + pT/GeV * (chi2Coeff[1] + pT/GeV * (chi2Coeff[2] + pT/GeV * chi2Coeff[3])))
+    float chi2Coeff[4];
+    float chi2MaxPt;  // GeV
+    float chi2Scale;
+
+    struct Region {
+      float maxTip;  // cm
+      float minPt;   // GeV
+      float maxZip;  // cm
+    };
+
+    Region triplet;
+    Region quadruplet;
+
+    __device__ __forceinline__ bool isHP(TrackSoAHeterogeneousT<TrackerTraits> const *__restrict__ tracks,
+                                         int nHits,
+                                         int it) const {
+      // impose "region cuts" based on the fit results (phi, Tip, pt, cotan(theta)), Zip)
+      // default cuts:
+      //   - for triplets:    |Tip| < 0.3 cm, pT > 0.5 GeV, |Zip| < 12.0 cm
+      //   - for quadruplets: |Tip| < 0.5 cm, pT > 0.3 GeV, |Zip| < 12.0 cm
+      // (see CAHitNtupletGeneratorGPU.cc)
+      auto const &region = (nHits > 3) ? quadruplet : triplet;
+      return (std::abs(tracks->tip(it)) < region.maxTip) and (tracks->pt(it) > region.minPt) and
+             (std::abs(tracks->zip(it)) < region.maxZip);
+    }
 
-  using TrackSoA = TrackSoAHeterogeneousT<maxNumber()>;
-  using TrajectoryState = TrajectoryStateSoAT<maxNumber()>;
-  using HitContainer = TrackSoA::HitContainer;
+    __device__ __forceinline__ bool strictCut(TrackSoAHeterogeneousT<TrackerTraits> const *__restrict__ tracks,
+                                              int it) const {
+      auto roughLog = [](float x) {
+        // max diff [0.5,12] at 1.25 0.16143
+        // average diff  0.0662998
+        union IF {
+          uint32_t i;
+          float f;
+        };
+        IF z;
+        z.f = x;
+        uint32_t lsb = 1 < 21;
+        z.i += lsb;
+        z.i >>= 21;
+        auto f = z.i & 3;
+        int ex = int(z.i >> 2) - 127;
+
+        // log2(1+0.25*f)
+        // averaged over bins
+        const float frac[4] = {0.160497f, 0.452172f, 0.694562f, 0.901964f};
+        return float(ex) + frac[f];
+      };
+
+      float pt = std::min<float>(tracks->pt(it), chi2MaxPt);
+      float chi2Cut = chi2Scale * (chi2Coeff[0] + roughLog(pt) * chi2Coeff[1]);
+      if (tracks->chi2(it) >= chi2Cut) {
+#ifdef NTUPLE_FIT_DEBUG
+        printf("Bad chi2 %d pt %f eta %f chi2 %f\n", it, tracks->pt(it), tracks->eta(it), tracks->chi2(it));
+#endif
+        return true;
+      }
+      return false;
+    }
+  };
+
+  template <typename TrackerTraits>
+  struct QualityCutsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> {
+    float maxChi2;
+    float minPt;
+    float maxTip;
+    float maxZip;
+
+    __device__ __forceinline__ bool isHP(TrackSoAHeterogeneousT<TrackerTraits> const *__restrict__ tracks,
+                                         int nHits,
+                                         int it) const {
+      return (std::abs(tracks->tip(it)) < maxTip) and (tracks->pt(it) > minPt) and (std::abs(tracks->zip(it)) < maxZip);
+    }
+    __device__ __forceinline__ bool strictCut(TrackSoAHeterogeneousT<TrackerTraits> const *__restrict__ tracks,
+                                              int it) const {
+      return tracks->chi2(it) >= maxChi2;
+    }
+  };
 
 }  // namespace pixelTrack
 

CUDADataFormats/Track/src/classes_def.xml

@@ -1,6 +1,15 @@
 <lcgdict>
-  <class name="cms::cuda::Product<HeterogeneousSoA<pixelTrack::TrackSoA>>" persistent="false"/>
-  <class name="edm::Wrapper<cms::cuda::Product<HeterogeneousSoA<pixelTrack::TrackSoA>>>" persistent="false"/>
-  <class name="HeterogeneousSoA<pixelTrack::TrackSoA>" persistent="false"/>
-  <class name="edm::Wrapper<HeterogeneousSoA<pixelTrack::TrackSoA>>" persistent="false"/>
+
+  <class name="pixelTrack::TrackSoAPhase1" persistent="false"/>
+  <class name="edm::Wrapper<HeterogeneousSoA<pixelTrack::TrackSoAPhase1>>" persistent="false" />
+  <class name="cms::cuda::Product<HeterogeneousSoA<pixelTrack::TrackSoAPhase1>>" persistent="false"/>
+  <class name="edm::Wrapper<cms::cuda::Product<HeterogeneousSoA<pixelTrack::TrackSoAPhase1>>>" persistent="false"/>
+  <class name="HeterogeneousSoA<pixelTrack::TrackSoAPhase1>" persistent="false" />
+
+  <class name="pixelTrack::TrackSoAPhase2" persistent="false"/>
+  <class name="edm::Wrapper<HeterogeneousSoA<pixelTrack::TrackSoAPhase2>>" persistent="false" />
+  <class name="cms::cuda::Product<HeterogeneousSoA<pixelTrack::TrackSoAPhase2>>" persistent="false"/>
+  <class name="edm::Wrapper<cms::cuda::Product<HeterogeneousSoA<pixelTrack::TrackSoAPhase2>>>" persistent="false"/>
+  <class name="HeterogeneousSoA<pixelTrack::TrackSoAPhase2>" persistent="false" />
+
 </lcgdict>