File indexing completed on 2024-04-06 12:28:33
0001 #ifndef RecoTracker_PixelSeeding_plugins_CAHitNtupletGeneratorKernels_h
0002 #define RecoTracker_PixelSeeding_plugins_CAHitNtupletGeneratorKernels_h
0003
0004
0005
0006
0007 #include "CUDADataFormats/Common/interface/HeterogeneousSoA.h"
0008 #include "CUDADataFormats/Track/interface/PixelTrackUtilities.h"
0009 #include "CUDADataFormats/Track/interface/TrackSoAHeterogeneousHost.h"
0010 #include "CUDADataFormats/TrackingRecHit/interface/TrackingRecHitSoADevice.h"
0011 #include "CUDADataFormats/TrackingRecHit/interface/TrackingRecHitsUtilities.h"
0012
0013 #include "GPUCACell.h"
0014 #include "gpuPixelDoublets.h"
0015
0016 namespace caHitNtupletGenerator {
0017
0018
0019 struct AlgoParams {
0020 const bool onGPU_;
0021 const uint32_t minHitsForSharingCut_;
0022 const bool useRiemannFit_;
0023 const bool fitNas4_;
0024 const bool includeJumpingForwardDoublets_;
0025 const bool earlyFishbone_;
0026 const bool lateFishbone_;
0027 const bool doStats_;
0028 const bool doSharedHitCut_;
0029 const bool dupPassThrough_;
0030 const bool useSimpleTripletCleaner_;
0031 };
0032
0033
0034 struct CACommon {
0035 const uint32_t maxNumberOfDoublets_;
0036 const uint32_t minHitsPerNtuplet_;
0037 const float ptmin_;
0038 const float CAThetaCutBarrel_;
0039 const float CAThetaCutForward_;
0040 const float hardCurvCut_;
0041 const float dcaCutInnerTriplet_;
0042 const float dcaCutOuterTriplet_;
0043 };
0044
0045 template <typename TrackerTraits, typename Enable = void>
0046 struct CAParamsT : public CACommon {
0047 __device__ __forceinline__ bool startingLayerPair(int16_t pid) const { return false; };
0048 __device__ __forceinline__ bool startAt0(int16_t pid) const { return false; };
0049 };
0050
0051 template <typename TrackerTraits>
0052 struct CAParamsT<TrackerTraits, pixelTopology::isPhase1Topology<TrackerTraits>> : public CACommon {
0053
0054 __device__ __forceinline__ bool startingLayerPair(int16_t pid) const {
0055 return minHitsPerNtuplet_ > 3 ? pid < 3 : pid < 8 || pid > 12;
0056 }
0057
0058
0059 __device__ __forceinline__ bool startAt0(int16_t pid) const {
0060 assert((pixelTopology::Phase1::layerPairs[pid * 2] == 0) ==
0061 (pid < 3 || pid == 13 || pid == 15 || pid == 16));
0062 return pixelTopology::Phase1::layerPairs[pid * 2] == 0;
0063 }
0064 };
0065
0066 template <typename TrackerTraits>
0067 struct CAParamsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> : public CACommon {
0068 const bool includeFarForwards_;
0069
0070 __device__ __forceinline__ bool startingLayerPair(int16_t pid) const {
0071 return pid < 33;
0072 }
0073
0074
0075 __device__ __forceinline__ bool startAt0(int16_t pid) const {
0076 assert((pixelTopology::Phase2::layerPairs[pid * 2] == 0) == ((pid < 3) | (pid >= 23 && pid < 28)));
0077 return pixelTopology::Phase2::layerPairs[pid * 2] == 0;
0078 }
0079 };
0080
0081
0082
0083 template <typename TrackerTraits, typename Enable = void>
0084 struct ParamsT : public AlgoParams {
0085
0086
0087 inline uint32_t nPairs() const { return 0; }
0088 };
0089
0090 template <typename TrackerTraits>
0091 struct ParamsT<TrackerTraits, pixelTopology::isPhase1Topology<TrackerTraits>> : public AlgoParams {
0092 using TT = TrackerTraits;
0093 using QualityCuts = pixelTrack::QualityCutsT<TT>;
0094 using CellCuts = gpuPixelDoublets::CellCutsT<TT>;
0095 using CAParams = CAParamsT<TT>;
0096
0097 ParamsT(AlgoParams const& commonCuts,
0098 CellCuts const& cellCuts,
0099 QualityCuts const& cutsCuts,
0100 CAParams const& caParams)
0101 : AlgoParams(commonCuts), cellCuts_(cellCuts), qualityCuts_(cutsCuts), caParams_(caParams) {}
0102
0103 const CellCuts cellCuts_;
0104 const QualityCuts qualityCuts_{
0105 {0.68177776, 0.74609577, -0.08035491, 0.00315399},
0106
0107 10.,
0108
0109 30.,
0110
0111 {
0112 0.3,
0113 0.5,
0114 12.0
0115 },
0116
0117 {
0118 0.5,
0119 0.3,
0120 12.0
0121 }};
0122 const CAParams caParams_;
0123
0124 inline uint32_t nPairs() const {
0125
0126 uint32_t nActualPairs = TT::nPairs;
0127 if (not includeJumpingForwardDoublets_) {
0128
0129 nActualPairs = TT::nPairsForTriplets;
0130 }
0131 if (caParams_.minHitsPerNtuplet_ > 3) {
0132
0133 nActualPairs = TT::nPairsForQuadruplets;
0134 }
0135
0136 return nActualPairs;
0137 }
0138
0139 };
0140
0141 template <typename TrackerTraits>
0142 struct ParamsT<TrackerTraits, pixelTopology::isPhase2Topology<TrackerTraits>> : public AlgoParams {
0143 using TT = TrackerTraits;
0144 using QualityCuts = pixelTrack::QualityCutsT<TT>;
0145 using CellCuts = gpuPixelDoublets::CellCutsT<TT>;
0146 using CAParams = CAParamsT<TT>;
0147
0148 ParamsT(AlgoParams const& commonCuts,
0149 CellCuts const& cellCuts,
0150 QualityCuts const& qualityCuts,
0151 CAParams const& caParams)
0152 : AlgoParams(commonCuts), cellCuts_(cellCuts), qualityCuts_(qualityCuts), caParams_(caParams) {}
0153
0154
0155 const CellCuts cellCuts_;
0156 const QualityCuts qualityCuts_{5.0f, 0.9f, 0.4f, 12.0f };
0157 const CAParams caParams_;
0158
0159 inline uint32_t nPairs() const {
0160
0161 uint32_t nActualPairs = TT::nPairsMinimal;
0162 if (caParams_.includeFarForwards_) {
0163
0164 nActualPairs = TT::nPairsFarForwards;
0165 }
0166 if (includeJumpingForwardDoublets_) {
0167
0168 nActualPairs = TT::nPairs;
0169 }
0170
0171 return nActualPairs;
0172 }
0173
0174 };
0175
0176
0177 struct Counters {
0178 unsigned long long nEvents;
0179 unsigned long long nHits;
0180 unsigned long long nCells;
0181 unsigned long long nTuples;
0182 unsigned long long nFitTracks;
0183 unsigned long long nLooseTracks;
0184 unsigned long long nGoodTracks;
0185 unsigned long long nUsedHits;
0186 unsigned long long nDupHits;
0187 unsigned long long nFishCells;
0188 unsigned long long nKilledCells;
0189 unsigned long long nEmptyCells;
0190 unsigned long long nZeroTrackCells;
0191 };
0192
0193 using Quality = pixelTrack::Quality;
0194
0195 }
0196
0197 template <typename TTraits, typename TTTraits>
0198 class CAHitNtupletGeneratorKernels {
0199 public:
0200 using Traits = TTraits;
0201 using TrackerTraits = TTTraits;
0202 using QualityCuts = pixelTrack::QualityCutsT<TrackerTraits>;
0203 using Params = caHitNtupletGenerator::ParamsT<TrackerTraits>;
0204 using CAParams = caHitNtupletGenerator::CAParamsT<TrackerTraits>;
0205 using CellCuts = gpuPixelDoublets::CellCutsT<TrackerTraits>;
0206 using Counters = caHitNtupletGenerator::Counters;
0207
0208 template <typename T>
0209 using unique_ptr = typename Traits::template unique_ptr<T>;
0210
0211 using HitsView = TrackingRecHitSoAView<TrackerTraits>;
0212 using HitsConstView = TrackingRecHitSoAConstView<TrackerTraits>;
0213 using TkSoAView = TrackSoAView<TrackerTraits>;
0214
0215 using HitToTuple = caStructures::HitToTupleT<TrackerTraits>;
0216 using TupleMultiplicity = caStructures::TupleMultiplicityT<TrackerTraits>;
0217 using CellNeighborsVector = caStructures::CellNeighborsVectorT<TrackerTraits>;
0218 using CellNeighbors = caStructures::CellNeighborsT<TrackerTraits>;
0219 using CellTracksVector = caStructures::CellTracksVectorT<TrackerTraits>;
0220 using CellTracks = caStructures::CellTracksT<TrackerTraits>;
0221 using OuterHitOfCellContainer = caStructures::OuterHitOfCellContainerT<TrackerTraits>;
0222 using OuterHitOfCell = caStructures::OuterHitOfCellT<TrackerTraits>;
0223
0224 using CACell = GPUCACellT<TrackerTraits>;
0225
0226 using Quality = pixelTrack::Quality;
0227 using HitContainer = typename TrackSoA<TrackerTraits>::HitContainer;
0228
0229 CAHitNtupletGeneratorKernels(Params const& params)
0230 : params_(params), paramsMaxDoubletes3Quarters_(3 * params.caParams_.maxNumberOfDoublets_ / 4) {}
0231
0232 ~CAHitNtupletGeneratorKernels() = default;
0233
0234 TupleMultiplicity const* tupleMultiplicity() const { return device_tupleMultiplicity_.get(); }
0235
0236 void launchKernels(const HitsConstView& hh, TkSoAView& track_view, cudaStream_t cudaStream);
0237
0238 void classifyTuples(const HitsConstView& hh, TkSoAView& track_view, cudaStream_t cudaStream);
0239
0240 void buildDoublets(const HitsConstView& hh, cudaStream_t stream);
0241 void allocateOnGPU(int32_t nHits, cudaStream_t stream);
0242 void cleanup(cudaStream_t cudaStream);
0243
0244 static void printCounters(Counters const* counters);
0245 void setCounters(Counters* counters) { counters_ = counters; }
0246
0247 protected:
0248 Counters* counters_ = nullptr;
0249
0250
0251 unique_ptr<unsigned char[]> cellStorage_;
0252 unique_ptr<CellNeighborsVector> device_theCellNeighbors_;
0253 CellNeighbors* device_theCellNeighborsContainer_;
0254 unique_ptr<CellTracksVector> device_theCellTracks_;
0255 CellTracks* device_theCellTracksContainer_;
0256
0257 unique_ptr<CACell[]> device_theCells_;
0258 unique_ptr<OuterHitOfCellContainer[]> device_isOuterHitOfCell_;
0259 OuterHitOfCell isOuterHitOfCell_;
0260 uint32_t* device_nCells_ = nullptr;
0261
0262 unique_ptr<HitToTuple> device_hitToTuple_;
0263 unique_ptr<uint32_t[]> device_hitToTupleStorage_;
0264 typename HitToTuple::View hitToTupleView_;
0265
0266 unique_ptr<CellCuts> device_cellCuts_;
0267
0268 cms::cuda::AtomicPairCounter* device_hitToTuple_apc_ = nullptr;
0269
0270 cms::cuda::AtomicPairCounter* device_hitTuple_apc_ = nullptr;
0271
0272 unique_ptr<TupleMultiplicity> device_tupleMultiplicity_;
0273
0274 unique_ptr<cms::cuda::AtomicPairCounter::c_type[]> device_storage_;
0275
0276
0277 Params params_;
0278
0279 const uint32_t paramsMaxDoubletes3Quarters_;
0280
0281 inline uint32_t nDoubletBlocks(uint32_t blockSize) {
0282
0283 return (paramsMaxDoubletes3Quarters_ + blockSize - 1) / blockSize;
0284 }
0285
0286
0287 inline uint32_t nQuadrupletBlocks(uint32_t blockSize) {
0288
0289 return (3 * TrackerTraits::maxNumberOfQuadruplets / 4 + blockSize - 1) / blockSize;
0290 }
0291 };
0292
0293 template <typename TrackerTraits>
0294 class CAHitNtupletGeneratorKernelsGPU : public CAHitNtupletGeneratorKernels<cms::cudacompat::GPUTraits, TrackerTraits> {
0295 using CAHitNtupletGeneratorKernels<cms::cudacompat::GPUTraits, TrackerTraits>::CAHitNtupletGeneratorKernels;
0296
0297 using Counters = caHitNtupletGenerator::Counters;
0298 using CAParams = caHitNtupletGenerator::CAParamsT<TrackerTraits>;
0299
0300 using HitContainer = typename TrackSoA<TrackerTraits>::HitContainer;
0301
0302 using CellNeighborsVector = caStructures::CellNeighborsVectorT<TrackerTraits>;
0303 using HitToTuple = caStructures::HitToTupleT<TrackerTraits>;
0304 using CellTracksVector = caStructures::CellTracksVectorT<TrackerTraits>;
0305 using TupleMultiplicity = caStructures::TupleMultiplicityT<TrackerTraits>;
0306
0307 using HitsConstView = TrackingRecHitSoAConstView<TrackerTraits>;
0308 using TkSoAView = TrackSoAView<TrackerTraits>;
0309
0310 using Params = caHitNtupletGenerator::ParamsT<TrackerTraits>;
0311
0312 public:
0313 void launchKernels(const HitsConstView& hh, TkSoAView& track_view, cudaStream_t cudaStream);
0314 void classifyTuples(const HitsConstView& hh, TkSoAView& track_view, cudaStream_t cudaStream);
0315 void buildDoublets(const HitsConstView& hh, int32_t offsetBPIX2, cudaStream_t stream);
0316 void allocateOnGPU(int32_t nHits, cudaStream_t stream);
0317 static void printCounters(Counters const* counters);
0318 };
0319
0320 template <typename TrackerTraits>
0321 class CAHitNtupletGeneratorKernelsCPU : public CAHitNtupletGeneratorKernels<cms::cudacompat::CPUTraits, TrackerTraits> {
0322 using CAHitNtupletGeneratorKernels<cms::cudacompat::CPUTraits, TrackerTraits>::CAHitNtupletGeneratorKernels;
0323
0324 using Counters = caHitNtupletGenerator::Counters;
0325 using CAParams = caHitNtupletGenerator::CAParamsT<TrackerTraits>;
0326
0327 using HitContainer = typename TrackSoA<TrackerTraits>::HitContainer;
0328
0329 using CellNeighborsVector = caStructures::CellNeighborsVectorT<TrackerTraits>;
0330 using HitToTuple = caStructures::HitToTupleT<TrackerTraits>;
0331 using CellTracksVector = caStructures::CellTracksVectorT<TrackerTraits>;
0332 using TupleMultiplicity = caStructures::TupleMultiplicityT<TrackerTraits>;
0333
0334 using HitsConstView = TrackingRecHitSoAConstView<TrackerTraits>;
0335 using TkSoAView = TrackSoAView<TrackerTraits>;
0336
0337 using Params = caHitNtupletGenerator::ParamsT<TrackerTraits>;
0338
0339 public:
0340 void launchKernels(const HitsConstView& hh, TkSoAView& track_view, cudaStream_t cudaStream);
0341 void classifyTuples(const HitsConstView& hh, TkSoAView& track_view, cudaStream_t cudaStream);
0342 void buildDoublets(const HitsConstView& hh, int32_t offsetBPIX2, cudaStream_t stream);
0343 void allocateOnGPU(int32_t nHits, cudaStream_t stream);
0344 static void printCounters(Counters const* counters);
0345 };
0346
0347 #endif