File indexing completed on 2024-04-23 22:56:30
0001 #ifndef RecoTracker_PixelSeeding_plugins_alpaka_CAHitNtupletGeneratorKernels_h
0002 #define RecoTracker_PixelSeeding_plugins_alpaka_CAHitNtupletGeneratorKernels_h
0003
0004
0005
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
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
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
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
0067 ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startAt0(int16_t pid) const {
0068 ALPAKA_ASSERT_ACC(
0069 (pixelTopology::Phase1::layerPairs[pid * 2] == 0) ==
0070 (pid < 3 || pid == 13 || pid == 15 || pid == 16));
0071 return pixelTopology::Phase1::layerPairs[pid * 2] == 0;
0072 }
0073 };
0074
0075 template <typename TrackerTraits>
0076 struct CAParamsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> : public CACommon {
0077 const bool includeFarForwards_;
0078
0079 ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startingLayerPair(int16_t pid) const {
0080 return pid < 33;
0081 }
0082
0083
0084 ALPAKA_FN_ACC ALPAKA_FN_INLINE bool startAt0(int16_t pid) const {
0085 ALPAKA_ASSERT_ACC((pixelTopology::Phase2::layerPairs[pid * 2] == 0) == ((pid < 3) | (pid >= 23 && pid < 28)));
0086 return pixelTopology::Phase2::layerPairs[pid * 2] == 0;
0087 }
0088 };
0089
0090
0091
0092 template <typename TrackerTraits, typename Enable = void>
0093 struct ParamsT : public AlgoParams {
0094
0095
0096 inline uint32_t nPairs() const { return 0; }
0097 };
0098
0099 template <typename TrackerTraits>
0100 struct ParamsT<TrackerTraits, pixelTopology::isPhase1Topology<TrackerTraits>> : public AlgoParams {
0101 using TT = TrackerTraits;
0102 using QualityCuts = ::pixelTrack::QualityCutsT<TT>;
0103 using CellCuts = caPixelDoublets::CellCutsT<TT>;
0104 using CAParams = CAParamsT<TT>;
0105
0106 ParamsT(AlgoParams const& commonCuts,
0107 CellCuts const& cellCuts,
0108 QualityCuts const& cutsCuts,
0109 CAParams const& caParams)
0110 : AlgoParams(commonCuts), cellCuts_(cellCuts), qualityCuts_(cutsCuts), caParams_(caParams) {}
0111
0112 const CellCuts cellCuts_;
0113 const QualityCuts qualityCuts_{
0114 {0.68177776, 0.74609577, -0.08035491, 0.00315399},
0115
0116 10.,
0117
0118 30.,
0119
0120 {
0121 0.3,
0122 0.5,
0123 12.0
0124 },
0125
0126 {
0127 0.5,
0128 0.3,
0129 12.0
0130 }};
0131 const CAParams caParams_;
0132
0133 inline uint32_t nPairs() const {
0134
0135 uint32_t nActualPairs = TT::nPairs;
0136 if (not includeJumpingForwardDoublets_) {
0137
0138 nActualPairs = TT::nPairsForTriplets;
0139 }
0140 if (caParams_.minHitsPerNtuplet_ > 3) {
0141
0142 nActualPairs = TT::nPairsForQuadruplets;
0143 }
0144
0145 return nActualPairs;
0146 }
0147
0148 };
0149
0150 template <typename TrackerTraits>
0151 struct ParamsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> : public AlgoParams {
0152 using TT = TrackerTraits;
0153 using QualityCuts = ::pixelTrack::QualityCutsT<TT>;
0154 using CellCuts = caPixelDoublets::CellCutsT<TT>;
0155 using CAParams = CAParamsT<TT>;
0156
0157 ParamsT(AlgoParams const& commonCuts,
0158 CellCuts const& cellCuts,
0159 QualityCuts const& qualityCuts,
0160 CAParams const& caParams)
0161 : AlgoParams(commonCuts), cellCuts_(cellCuts), qualityCuts_(qualityCuts), caParams_(caParams) {}
0162
0163
0164 const CellCuts cellCuts_;
0165 const QualityCuts qualityCuts_{5.0f, 0.9f, 0.4f, 12.0f };
0166 const CAParams caParams_;
0167
0168 inline uint32_t nPairs() const {
0169
0170 uint32_t nActualPairs = TT::nPairsMinimal;
0171 if (caParams_.includeFarForwards_) {
0172
0173 nActualPairs = TT::nPairsFarForwards;
0174 }
0175 if (includeJumpingForwardDoublets_) {
0176
0177 nActualPairs = TT::nPairs;
0178 }
0179
0180 return nActualPairs;
0181 }
0182
0183 };
0184
0185
0186 struct Counters {
0187 unsigned long long nEvents;
0188 unsigned long long nHits;
0189 unsigned long long nCells;
0190 unsigned long long nTuples;
0191 unsigned long long nFitTracks;
0192 unsigned long long nLooseTracks;
0193 unsigned long long nGoodTracks;
0194 unsigned long long nUsedHits;
0195 unsigned long long nDupHits;
0196 unsigned long long nFishCells;
0197 unsigned long long nKilledCells;
0198 unsigned long long nEmptyCells;
0199 unsigned long long nZeroTrackCells;
0200 };
0201
0202 using Quality = ::pixelTrack::Quality;
0203
0204 }
0205
0206 template <typename TTTraits>
0207 class CAHitNtupletGeneratorKernels {
0208 public:
0209 using TrackerTraits = TTTraits;
0210 using QualityCuts = ::pixelTrack::QualityCutsT<TrackerTraits>;
0211 using CellCuts = caPixelDoublets::CellCutsT<TrackerTraits>;
0212 using Params = caHitNtupletGenerator::ParamsT<TrackerTraits>;
0213 using CAParams = caHitNtupletGenerator::CAParamsT<TrackerTraits>;
0214 using Counters = caHitNtupletGenerator::Counters;
0215
0216 using HitsView = TrackingRecHitSoAView<TrackerTraits>;
0217 using HitsConstView = TrackingRecHitSoAConstView<TrackerTraits>;
0218 using TkSoAView = reco::TrackSoAView<TrackerTraits>;
0219
0220 using HitToTuple = caStructures::template HitToTupleT<TrackerTraits>;
0221 using TupleMultiplicity = caStructures::template TupleMultiplicityT<TrackerTraits>;
0222 struct Testttt {
0223 TupleMultiplicity tm;
0224 };
0225 using CellNeighborsVector = caStructures::CellNeighborsVectorT<TrackerTraits>;
0226 using CellNeighbors = caStructures::CellNeighborsT<TrackerTraits>;
0227 using CellTracksVector = caStructures::CellTracksVectorT<TrackerTraits>;
0228 using CellTracks = caStructures::CellTracksT<TrackerTraits>;
0229 using OuterHitOfCellContainer = caStructures::OuterHitOfCellContainerT<TrackerTraits>;
0230 using OuterHitOfCell = caStructures::OuterHitOfCellT<TrackerTraits>;
0231
0232 using CACell = CACellT<TrackerTraits>;
0233
0234 using Quality = ::pixelTrack::Quality;
0235 using HitContainer = typename reco::TrackSoA<TrackerTraits>::HitContainer;
0236
0237 CAHitNtupletGeneratorKernels(Params const& params, uint32_t nhits, uint32_t offsetBPIX2, Queue& queue);
0238 ~CAHitNtupletGeneratorKernels() = default;
0239
0240 TupleMultiplicity const* tupleMultiplicity() const { return device_tupleMultiplicity_.data(); }
0241
0242 void launchKernels(const HitsConstView& hh, uint32_t offsetBPIX2, TkSoAView& track_view, Queue& queue);
0243
0244 void classifyTuples(const HitsConstView& hh, TkSoAView& track_view, Queue& queue);
0245
0246 void buildDoublets(const HitsConstView& hh, uint32_t offsetBPIX2, Queue& queue);
0247
0248 static void printCounters();
0249
0250 private:
0251
0252 Params const& m_params;
0253 cms::alpakatools::device_buffer<Device, Counters> counters_;
0254
0255
0256 cms::alpakatools::device_buffer<Device, HitToTuple> device_hitToTuple_;
0257 cms::alpakatools::device_buffer<Device, uint32_t[]> device_hitToTupleStorage_;
0258 typename HitToTuple::View device_hitToTupleView_;
0259 cms::alpakatools::device_buffer<Device, TupleMultiplicity> device_tupleMultiplicity_;
0260 cms::alpakatools::device_buffer<Device, CACell[]> device_theCells_;
0261 cms::alpakatools::device_buffer<Device, OuterHitOfCellContainer[]> device_isOuterHitOfCell_;
0262 cms::alpakatools::device_buffer<Device, OuterHitOfCell> isOuterHitOfCell_;
0263 cms::alpakatools::device_buffer<Device, CellNeighborsVector> device_theCellNeighbors_;
0264 cms::alpakatools::device_buffer<Device, CellTracksVector> device_theCellTracks_;
0265 cms::alpakatools::device_buffer<Device, unsigned char[]> cellStorage_;
0266 cms::alpakatools::device_buffer<Device, CellCuts> device_cellCuts_;
0267 CellNeighbors* device_theCellNeighborsContainer_;
0268 CellTracks* device_theCellTracksContainer_;
0269 cms::alpakatools::device_buffer<Device, cms::alpakatools::AtomicPairCounter::DoubleWord[]> device_storage_;
0270 cms::alpakatools::AtomicPairCounter* device_hitTuple_apc_;
0271 cms::alpakatools::AtomicPairCounter* device_hitToTuple_apc_;
0272 cms::alpakatools::device_view<Device, uint32_t> device_nCells_;
0273 };
0274
0275 }
0276
0277 #endif