Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-10-01 22:40:58

0001 #ifndef RecoPixelVertexing_PixelTriplets_plugins_gpuPixelDoublets_h
0002 #define RecoPixelVertexing_PixelTriplets_plugins_gpuPixelDoublets_h
0003 
0004 #include "RecoPixelVertexing/PixelTriplets/plugins/gpuPixelDoubletsAlgos.h"
0005 
0006 #define CONSTANT_VAR __constant__
0007 
0008 namespace gpuPixelDoublets {
0009 
0010   constexpr int nPairsForQuadruplets = 13;                     // quadruplets require hits in all layers
0011   constexpr int nPairsForTriplets = nPairsForQuadruplets + 2;  // include barrel "jumping" layer pairs
0012   constexpr int nPairs = nPairsForTriplets + 4;                // include forward "jumping" layer pairs
0013   static_assert(nPairs <= caConstants::maxNumberOfLayerPairs);
0014 
0015   // start constants
0016   // clang-format off
0017 
0018   CONSTANT_VAR const uint8_t layerPairs[2 * nPairs] = {
0019       0, 1, 0, 4, 0, 7,              // BPIX1 (3)
0020       1, 2, 1, 4, 1, 7,              // BPIX2 (6)
0021       4, 5, 7, 8,                    // FPIX1 (8)
0022       2, 3, 2, 4, 2, 7, 5, 6, 8, 9,  // BPIX3 & FPIX2 (13)
0023       0, 2, 1, 3,                    // Jumping Barrel (15)
0024       0, 5, 0, 8,                    // Jumping Forward (BPIX1,FPIX2)
0025       4, 6, 7, 9                     // Jumping Forward (19)
0026   };
0027 
0028   constexpr int16_t phi0p05 = 522;  // round(521.52189...) = phi2short(0.05);
0029   constexpr int16_t phi0p06 = 626;  // round(625.82270...) = phi2short(0.06);
0030   constexpr int16_t phi0p07 = 730;  // round(730.12648...) = phi2short(0.07);
0031 
0032   CONSTANT_VAR const int16_t phicuts[nPairs]{phi0p05,
0033                                              phi0p07,
0034                                              phi0p07,
0035                                              phi0p05,
0036                                              phi0p06,
0037                                              phi0p06,
0038                                              phi0p05,
0039                                              phi0p05,
0040                                              phi0p06,
0041                                              phi0p06,
0042                                              phi0p06,
0043                                              phi0p05,
0044                                              phi0p05,
0045                                              phi0p05,
0046                                              phi0p05,
0047                                              phi0p05,
0048                                              phi0p05,
0049                                              phi0p05,
0050                                              phi0p05};
0051   //   phi0p07, phi0p07, phi0p06,phi0p06, phi0p06,phi0p06};  // relaxed cuts
0052 
0053   CONSTANT_VAR float const minz[nPairs] = {
0054       -20., 0., -30., -22., 10., -30., -70., -70., -22., 15., -30, -70., -70., -20., -22., 0, -30., -70., -70.};
0055   CONSTANT_VAR float const maxz[nPairs] = {
0056       20., 30., 0., 22., 30., -10., 70., 70., 22., 30., -15., 70., 70., 20., 22., 30., 0., 70., 70.};
0057   CONSTANT_VAR float const maxr[nPairs] = {
0058       20., 9., 9., 20., 7., 7., 5., 5., 20., 6., 6., 5., 5., 20., 20., 9., 9., 9., 9.};
0059 
0060   // end constants
0061   // clang-format on
0062 
0063   using CellNeighbors = caConstants::CellNeighbors;
0064   using CellTracks = caConstants::CellTracks;
0065   using CellNeighborsVector = caConstants::CellNeighborsVector;
0066   using CellTracksVector = caConstants::CellTracksVector;
0067 
0068   __global__ void initDoublets(GPUCACell::OuterHitOfCell isOuterHitOfCell,
0069                                int nHits,
0070                                CellNeighborsVector* cellNeighbors,
0071                                CellNeighbors* cellNeighborsContainer,
0072                                CellTracksVector* cellTracks,
0073                                CellTracks* cellTracksContainer) {
0074     assert(isOuterHitOfCell.container);
0075     int first = blockIdx.x * blockDim.x + threadIdx.x;
0076     for (int i = first; i < nHits - isOuterHitOfCell.offset; i += gridDim.x * blockDim.x)
0077       isOuterHitOfCell.container[i].reset();
0078 
0079     if (0 == first) {
0080       cellNeighbors->construct(caConstants::maxNumOfActiveDoublets, cellNeighborsContainer);
0081       cellTracks->construct(caConstants::maxNumOfActiveDoublets, cellTracksContainer);
0082       auto i = cellNeighbors->extend();
0083       assert(0 == i);
0084       (*cellNeighbors)[0].reset();
0085       i = cellTracks->extend();
0086       assert(0 == i);
0087       (*cellTracks)[0].reset();
0088     }
0089   }
0090 
0091   constexpr auto getDoubletsFromHistoMaxBlockSize = 64;  // for both x and y
0092   constexpr auto getDoubletsFromHistoMinBlocksPerMP = 16;
0093 
0094   __global__
0095 #ifdef __CUDACC__
0096   __launch_bounds__(getDoubletsFromHistoMaxBlockSize, getDoubletsFromHistoMinBlocksPerMP)
0097 #endif
0098       void getDoubletsFromHisto(GPUCACell* cells,
0099                                 uint32_t* nCells,
0100                                 CellNeighborsVector* cellNeighbors,
0101                                 CellTracksVector* cellTracks,
0102                                 TrackingRecHit2DSOAView const* __restrict__ hhp,
0103                                 GPUCACell::OuterHitOfCell isOuterHitOfCell,
0104                                 int nActualPairs,
0105                                 bool ideal_cond,
0106                                 bool doClusterCut,
0107                                 bool doZ0Cut,
0108                                 bool doPtCut,
0109                                 uint32_t maxNumOfDoublets) {
0110     auto const& __restrict__ hh = *hhp;
0111     doubletsFromHisto(layerPairs,
0112                       nActualPairs,
0113                       cells,
0114                       nCells,
0115                       cellNeighbors,
0116                       cellTracks,
0117                       hh,
0118                       isOuterHitOfCell,
0119                       phicuts,
0120                       minz,
0121                       maxz,
0122                       maxr,
0123                       ideal_cond,
0124                       doClusterCut,
0125                       doZ0Cut,
0126                       doPtCut,
0127                       maxNumOfDoublets);
0128   }
0129 
0130 }  // namespace gpuPixelDoublets
0131 
0132 #endif  // RecoPixelVertexing_PixelTriplets_plugins_gpuPixelDoublets_h