Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-19 02:16:54

0001 #ifndef RecoTracker_PixelSeeding_plugins_alpaka_CAHitNtupletGeneratorKernels_h
0002 #define RecoTracker_PixelSeeding_plugins_alpaka_CAHitNtupletGeneratorKernels_h
0003 
0004 //#define GPU_DEBUG
0005 //#define DUMP_GPU_TK_TUPLES
0006 
0007 #include <cstdint>
0008 
0009 #include <alpaka/alpaka.hpp>
0010 
0011 #include "DataFormats/TrackSoA/interface/TrackDefinitions.h"
0012 #include "DataFormats/TrackSoA/interface/TracksHost.h"
0013 #include "DataFormats/TrackSoA/interface/alpaka/TrackUtilities.h"
0014 #include "DataFormats/TrackingRecHitSoA/interface/TrackingRecHitsSoA.h"
0015 #include "HeterogeneousCore/AlpakaInterface/interface/AtomicPairCounter.h"
0016 #include "HeterogeneousCore/AlpakaInterface/interface/HistoContainer.h"
0017 #include "HeterogeneousCore/AlpakaInterface/interface/config.h"
0018 #include "HeterogeneousCore/AlpakaInterface/interface/memory.h"
0019 
0020 #include "CACell.h"
0021 #include "CAPixelDoublets.h"
0022 #include "CAStructures.h"
0023 
0024 namespace ALPAKA_ACCELERATOR_NAMESPACE {
0025   namespace caHitNtupletGenerator {
0026 
0027     //Configuration params common to all topologies, for the algorithms
0028     struct AlgoParams {
0029       const uint32_t minHitsForSharingCut_;
0030       const bool useRiemannFit_;
0031       const bool fitNas4_;
0032       const bool includeJumpingForwardDoublets_;
0033       const bool earlyFishbone_;
0034       const bool lateFishbone_;
0035       const bool doStats_;
0036       const bool doSharedHitCut_;
0037       const bool dupPassThrough_;
0038       const bool useSimpleTripletCleaner_;
0039     };
0040 
0041     //CAParams
0042     struct CACommon {
0043       const uint32_t maxNumberOfDoublets_;
0044       const uint32_t minHitsPerNtuplet_;
0045       const float ptmin_;
0046       const float CAThetaCutBarrel_;
0047       const float CAThetaCutForward_;
0048       const float hardCurvCut_;
0049       const float dcaCutInnerTriplet_;
0050       const float dcaCutOuterTriplet_;
0051     };
0052 
0053     template <typename TrackerTraits, typename Enable = void>
0054     struct CAParamsT : public CACommon {
0055       ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startingLayerPair(int16_t pid) const { return false; };
0056       ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startAt0(int16_t pid) const { return false; };
0057     };
0058 
0059     template <typename TrackerTraits>
0060     struct CAParamsT<TrackerTraits, pixelTopology::isPhase1Topology<TrackerTraits>> : public CACommon {
0061       /// Is is a starting layer pair?
0062       ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startingLayerPair(int16_t pid) const {
0063         return minHitsPerNtuplet_ > 3 ? pid < 3 : pid < 8 || pid > 12;
0064       }
0065 
0066       /// Is this a pair with inner == 0?
0067       ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startAt0(int16_t pid) const {
0068         assert((pixelTopology::Phase1::layerPairs[pid * 2] == 0) ==
0069                (pid < 3 || pid == 13 || pid == 15 || pid == 16));  // to be 100% sure it's working, may be removed
0070         return pixelTopology::Phase1::layerPairs[pid * 2] == 0;
0071       }
0072     };
0073 
0074     template <typename TrackerTraits>
0075     struct CAParamsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> : public CACommon {
0076       const bool includeFarForwards_;
0077       /// Is is a starting layer pair?
0078       ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startingLayerPair(int16_t pid) const {
0079         return pid < 33;  // in principle one could remove 5,6,7 23, 28 and 29
0080       }
0081 
0082       /// Is this a pair with inner == 0
0083       ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startAt0(int16_t pid) const {
0084         assert((pixelTopology::Phase2::layerPairs[pid * 2] == 0) == ((pid < 3) | (pid >= 23 && pid < 28)));
0085         return pixelTopology::Phase2::layerPairs[pid * 2] == 0;
0086       }
0087     };
0088 
0089     //Full list of params = algo params + ca params + cell params + quality cuts
0090     //Generic template
0091     template <typename TrackerTraits, typename Enable = void>
0092     struct ParamsT : public AlgoParams {
0093       // one should define the params for its own pixelTopology
0094       // not defining anything here
0095       inline uint32_t nPairs() const { return 0; }
0096     };
0097 
0098     template <typename TrackerTraits>
0099     struct ParamsT<TrackerTraits, pixelTopology::isPhase1Topology<TrackerTraits>> : public AlgoParams {
0100       using TT = TrackerTraits;
0101       using QualityCuts = ::pixelTrack::QualityCutsT<TT>;  //track quality cuts
0102       using CellCuts = caPixelDoublets::CellCutsT<TT>;     //cell building cuts
0103       using CAParams = CAParamsT<TT>;                      //params to be used on device
0104 
0105       ParamsT(AlgoParams const& commonCuts,
0106               CellCuts const& cellCuts,
0107               QualityCuts const& cutsCuts,
0108               CAParams const& caParams)
0109           : AlgoParams(commonCuts), cellCuts_(cellCuts), qualityCuts_(cutsCuts), caParams_(caParams) {}
0110 
0111       const CellCuts cellCuts_;
0112       const QualityCuts qualityCuts_{// polynomial coefficients for the pT-dependent chi2 cut
0113                                      {0.68177776, 0.74609577, -0.08035491, 0.00315399},
0114                                      // max pT used to determine the chi2 cut
0115                                      10.,
0116                                      // chi2 scale factor: 30 for broken line fit, 45 for Riemann fit
0117                                      30.,
0118                                      // regional cuts for triplets
0119                                      {
0120                                          0.3,  // |Tip| < 0.3 cm
0121                                          0.5,  // pT > 0.5 GeV
0122                                          12.0  // |Zip| < 12.0 cm
0123                                      },
0124                                      // regional cuts for quadruplets
0125                                      {
0126                                          0.5,  // |Tip| < 0.5 cm
0127                                          0.3,  // pT > 0.3 GeV
0128                                          12.0  // |Zip| < 12.0 cm
0129                                      }};
0130       const CAParams caParams_;
0131       /// Compute the number of pairs
0132       inline uint32_t nPairs() const {
0133         // take all layer pairs into account
0134         uint32_t nActualPairs = TT::nPairs;
0135         if (not includeJumpingForwardDoublets_) {
0136           // exclude forward "jumping" layer pairs
0137           nActualPairs = TT::nPairsForTriplets;
0138         }
0139         if (caParams_.minHitsPerNtuplet_ > 3) {
0140           // for quadruplets, exclude all "jumping" layer pairs
0141           nActualPairs = TT::nPairsForQuadruplets;
0142         }
0143 
0144         return nActualPairs;
0145       }
0146 
0147     };  // Params Phase1
0148 
0149     template <typename TrackerTraits>
0150     struct ParamsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> : public AlgoParams {
0151       using TT = TrackerTraits;
0152       using QualityCuts = ::pixelTrack::QualityCutsT<TT>;
0153       using CellCuts = caPixelDoublets::CellCutsT<TT>;
0154       using CAParams = CAParamsT<TT>;
0155 
0156       ParamsT(AlgoParams const& commonCuts,
0157               CellCuts const& cellCuts,
0158               QualityCuts const& qualityCuts,
0159               CAParams const& caParams)
0160           : AlgoParams(commonCuts), cellCuts_(cellCuts), qualityCuts_(qualityCuts), caParams_(caParams) {}
0161 
0162       // quality cuts
0163       const CellCuts cellCuts_;
0164       const QualityCuts qualityCuts_{5.0f, /*chi2*/ 0.9f, /* pT in Gev*/ 0.4f, /*zip in cm*/ 12.0f /*tip in cm*/};
0165       const CAParams caParams_;
0166 
0167       inline uint32_t nPairs() const {
0168         // take all layer pairs into account
0169         uint32_t nActualPairs = TT::nPairsMinimal;
0170         if (caParams_.includeFarForwards_) {
0171           // considera far forwards (> 11 & > 23)
0172           nActualPairs = TT::nPairsFarForwards;
0173         }
0174         if (includeJumpingForwardDoublets_) {
0175           // include jumping forwards
0176           nActualPairs = TT::nPairs;
0177         }
0178 
0179         return nActualPairs;
0180       }
0181 
0182     };  // Params Phase1
0183 
0184     // counters
0185     struct Counters {
0186       unsigned long long nEvents;
0187       unsigned long long nHits;
0188       unsigned long long nCells;
0189       unsigned long long nTuples;
0190       unsigned long long nFitTracks;
0191       unsigned long long nLooseTracks;
0192       unsigned long long nGoodTracks;
0193       unsigned long long nUsedHits;
0194       unsigned long long nDupHits;
0195       unsigned long long nFishCells;
0196       unsigned long long nKilledCells;
0197       unsigned long long nEmptyCells;
0198       unsigned long long nZeroTrackCells;
0199     };
0200 
0201     using Quality = ::pixelTrack::Quality;
0202 
0203   }  // namespace caHitNtupletGenerator
0204 
0205   template <typename TTTraits>
0206   class CAHitNtupletGeneratorKernels {
0207   public:
0208     using TrackerTraits = TTTraits;
0209     using QualityCuts = ::pixelTrack::QualityCutsT<TrackerTraits>;
0210     using CellCuts = caPixelDoublets::CellCutsT<TrackerTraits>;
0211     using Params = caHitNtupletGenerator::ParamsT<TrackerTraits>;
0212     using CAParams = caHitNtupletGenerator::CAParamsT<TrackerTraits>;
0213     using Counters = caHitNtupletGenerator::Counters;
0214 
0215     using HitsView = TrackingRecHitSoAView<TrackerTraits>;
0216     using HitsConstView = TrackingRecHitSoAConstView<TrackerTraits>;
0217     using TkSoAView = reco::TrackSoAView<TrackerTraits>;
0218 
0219     using HitToTuple = caStructures::template HitToTupleT<TrackerTraits>;
0220     using TupleMultiplicity = caStructures::template TupleMultiplicityT<TrackerTraits>;
0221     struct Testttt {
0222       TupleMultiplicity tm;
0223     };
0224     using CellNeighborsVector = caStructures::CellNeighborsVectorT<TrackerTraits>;
0225     using CellNeighbors = caStructures::CellNeighborsT<TrackerTraits>;
0226     using CellTracksVector = caStructures::CellTracksVectorT<TrackerTraits>;
0227     using CellTracks = caStructures::CellTracksT<TrackerTraits>;
0228     using OuterHitOfCellContainer = caStructures::OuterHitOfCellContainerT<TrackerTraits>;
0229     using OuterHitOfCell = caStructures::OuterHitOfCellT<TrackerTraits>;
0230 
0231     using CACell = CACellT<TrackerTraits>;
0232 
0233     using Quality = ::pixelTrack::Quality;
0234     using HitContainer = typename reco::TrackSoA<TrackerTraits>::HitContainer;
0235 
0236     CAHitNtupletGeneratorKernels(Params const& params, uint32_t nhits, uint32_t offsetBPIX2, Queue& queue);
0237     ~CAHitNtupletGeneratorKernels() = default;
0238 
0239     TupleMultiplicity const* tupleMultiplicity() const { return device_tupleMultiplicity_.data(); }
0240 
0241     void launchKernels(const HitsConstView& hh, uint32_t offsetBPIX2, TkSoAView& track_view, Queue& queue);
0242 
0243     void classifyTuples(const HitsConstView& hh, TkSoAView& track_view, Queue& queue);
0244 
0245     void buildDoublets(const HitsConstView& hh, uint32_t offsetBPIX2, Queue& queue);
0246 
0247     static void printCounters();
0248 
0249   private:
0250     // params
0251     Params const& m_params;
0252     cms::alpakatools::device_buffer<Device, Counters> counters_;
0253 
0254     // workspace
0255     cms::alpakatools::device_buffer<Device, HitToTuple> device_hitToTuple_;
0256     cms::alpakatools::device_buffer<Device, uint32_t[]> device_hitToTupleStorage_;
0257     typename HitToTuple::View device_hitToTupleView_;
0258     cms::alpakatools::device_buffer<Device, TupleMultiplicity> device_tupleMultiplicity_;
0259     cms::alpakatools::device_buffer<Device, CACell[]> device_theCells_;
0260     cms::alpakatools::device_buffer<Device, OuterHitOfCellContainer[]> device_isOuterHitOfCell_;
0261     cms::alpakatools::device_buffer<Device, OuterHitOfCell> isOuterHitOfCell_;
0262     cms::alpakatools::device_buffer<Device, CellNeighborsVector> device_theCellNeighbors_;
0263     cms::alpakatools::device_buffer<Device, CellTracksVector> device_theCellTracks_;
0264     cms::alpakatools::device_buffer<Device, unsigned char[]> cellStorage_;
0265     cms::alpakatools::device_buffer<Device, CellCuts> device_cellCuts_;
0266     CellNeighbors* device_theCellNeighborsContainer_;
0267     CellTracks* device_theCellTracksContainer_;
0268     cms::alpakatools::device_buffer<Device, cms::alpakatools::AtomicPairCounter::DoubleWord[]> device_storage_;
0269     cms::alpakatools::AtomicPairCounter* device_hitTuple_apc_;
0270     cms::alpakatools::AtomicPairCounter* device_hitToTuple_apc_;
0271     cms::alpakatools::device_view<Device, uint32_t> device_nCells_;
0272   };
0273 
0274 }  // namespace ALPAKA_ACCELERATOR_NAMESPACE
0275 
0276 #endif  // RecoTracker_PixelSeeding_plugins_alpaka_CAHitNtupletGeneratorKernels_h