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