Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:28:32

0001 #ifndef RecoTracker_PixelSeeding_plugins_alpaka_CAPixelDoublets_h
0002 #define RecoTracker_PixelSeeding_plugins_alpaka_CAPixelDoublets_h
0003 
0004 #include <type_traits>
0005 
0006 #include <alpaka/alpaka.hpp>
0007 
0008 #include "HeterogeneousCore/AlpakaInterface/interface/config.h"
0009 #include "HeterogeneousCore/AlpakaInterface/interface/traits.h"
0010 #include "HeterogeneousCore/AlpakaInterface/interface/workdivision.h"
0011 
0012 #include "CAPixelDoubletsAlgos.h"
0013 
0014 namespace ALPAKA_ACCELERATOR_NAMESPACE {
0015   using namespace alpaka;
0016   using namespace cms::alpakatools;
0017   namespace caPixelDoublets {
0018 
0019     template <typename TrackerTraits>
0020     class InitDoublets {
0021     public:
0022       template <typename TAcc, typename = std::enable_if_t<alpaka::isAccelerator<TAcc>>>
0023       ALPAKA_FN_ACC void operator()(TAcc const& acc,
0024                                     OuterHitOfCell<TrackerTraits>* isOuterHitOfCell,
0025                                     int nHits,
0026                                     CellNeighborsVector<TrackerTraits>* cellNeighbors,
0027                                     CellNeighbors<TrackerTraits>* cellNeighborsContainer,
0028                                     CellTracksVector<TrackerTraits>* cellTracks,
0029                                     CellTracks<TrackerTraits>* cellTracksContainer) const {
0030         ALPAKA_ASSERT_ACC((*isOuterHitOfCell).container);
0031 
0032         for (auto i : cms::alpakatools::uniform_elements(acc, nHits - isOuterHitOfCell->offset))
0033           (*isOuterHitOfCell).container[i].reset();
0034 
0035         if (cms::alpakatools::once_per_grid(acc)) {
0036           cellNeighbors->construct(TrackerTraits::maxNumOfActiveDoublets, cellNeighborsContainer);
0037           cellTracks->construct(TrackerTraits::maxNumOfActiveDoublets, cellTracksContainer);
0038           [[maybe_unused]] auto i = cellNeighbors->extend(acc);
0039           ALPAKA_ASSERT_ACC(0 == i);
0040           (*cellNeighbors)[0].reset();
0041           i = cellTracks->extend(acc);
0042           ALPAKA_ASSERT_ACC(0 == i);
0043           (*cellTracks)[0].reset();
0044         }
0045       }
0046     };
0047 
0048     // Not used for the moment, see below.
0049     //constexpr auto getDoubletsFromHistoMaxBlockSize = 64;  // for both x and y
0050     //constexpr auto getDoubletsFromHistoMinBlocksPerMP = 16;
0051 
0052     template <typename TrackerTraits>
0053     class GetDoubletsFromHisto {
0054     public:
0055       template <typename TAcc, typename = std::enable_if_t<alpaka::isAccelerator<TAcc>>>
0056       // #ifdef __CUDACC__
0057       //       __launch_bounds__(getDoubletsFromHistoMaxBlockSize, getDoubletsFromHistoMinBlocksPerMP)  // TODO: Alapakify
0058       // #endif
0059       ALPAKA_FN_ACC void operator()(TAcc const& acc,
0060                                     CACellT<TrackerTraits>* cells,
0061                                     uint32_t* nCells,
0062                                     CellNeighborsVector<TrackerTraits>* cellNeighbors,
0063                                     CellTracksVector<TrackerTraits>* cellTracks,
0064                                     HitsConstView<TrackerTraits> hh,
0065                                     OuterHitOfCell<TrackerTraits>* isOuterHitOfCell,
0066                                     uint32_t nActualPairs,
0067                                     const uint32_t maxNumOfDoublets,
0068                                     CellCutsT<TrackerTraits> cuts) const {
0069         doubletsFromHisto<TrackerTraits>(
0070             acc, nActualPairs, maxNumOfDoublets, cells, nCells, cellNeighbors, cellTracks, hh, *isOuterHitOfCell, cuts);
0071       }
0072     };
0073 
0074   }  // namespace caPixelDoublets
0075 
0076 }  // namespace ALPAKA_ACCELERATOR_NAMESPACE
0077 
0078 #endif  // RecoTracker_PixelSeeding_plugins_alpaka_CAPixelDoublets_h