Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-04 00:15:52

0001 #ifndef RecoPixelVertexing_PixelTriplets_plugins_gpuFishbone_h
0002 #define RecoPixelVertexing_PixelTriplets_plugins_gpuFishbone_h
0003 
0004 #include <algorithm>
0005 #include <cmath>
0006 #include <cstdint>
0007 #include <cstdio>
0008 #include <limits>
0009 
0010 #include "DataFormats/Math/interface/approx_atan2.h"
0011 #include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
0012 #include "HeterogeneousCore/CUDAUtilities/interface/VecArray.h"
0013 #include "HeterogeneousCore/CUDAUtilities/interface/cuda_assert.h"
0014 
0015 #include "GPUCACell.h"
0016 
0017 namespace gpuPixelDoublets {
0018 
0019   __global__ void fishbone(GPUCACell::Hits const* __restrict__ hhp,
0020                            GPUCACell* cells,
0021                            uint32_t const* __restrict__ nCells,
0022                            GPUCACell::OuterHitOfCell const isOuterHitOfCellWrap,
0023                            int32_t nHits,
0024                            bool checkTrack) {
0025     constexpr auto maxCellsPerHit = GPUCACell::maxCellsPerHit;
0026 
0027     auto const& hh = *hhp;
0028 
0029     auto const isOuterHitOfCell = isOuterHitOfCellWrap.container;
0030     int32_t offset = isOuterHitOfCellWrap.offset;
0031 
0032     // x runs faster...
0033     auto firstY = threadIdx.y + blockIdx.y * blockDim.y;
0034     auto firstX = threadIdx.x;
0035 
0036     float x[maxCellsPerHit], y[maxCellsPerHit], z[maxCellsPerHit], n[maxCellsPerHit];
0037     uint32_t cc[maxCellsPerHit];
0038     uint16_t d[maxCellsPerHit];
0039     uint8_t l[maxCellsPerHit];
0040 
0041     for (int idy = firstY, nt = nHits - offset; idy < nt; idy += gridDim.y * blockDim.y) {
0042       auto const& vc = isOuterHitOfCell[idy];
0043       auto size = vc.size();
0044       if (size < 2)
0045         continue;
0046       // if alligned kill one of the two.
0047       // in principle one could try to relax the cut (only in r-z?) for jumping-doublets
0048       auto const& c0 = cells[vc[0]];
0049       auto xo = c0.outer_x(hh);
0050       auto yo = c0.outer_y(hh);
0051       auto zo = c0.outer_z(hh);
0052       auto sg = 0;
0053       for (int32_t ic = 0; ic < size; ++ic) {
0054         auto& ci = cells[vc[ic]];
0055         if (ci.unused())
0056           continue;  // for triplets equivalent to next
0057         if (checkTrack && ci.tracks().empty())
0058           continue;
0059         cc[sg] = vc[ic];
0060         l[sg] = ci.layerPairId();
0061         d[sg] = ci.inner_detIndex(hh);
0062         x[sg] = ci.inner_x(hh) - xo;
0063         y[sg] = ci.inner_y(hh) - yo;
0064         z[sg] = ci.inner_z(hh) - zo;
0065         n[sg] = x[sg] * x[sg] + y[sg] * y[sg] + z[sg] * z[sg];
0066         ++sg;
0067       }
0068       if (sg < 2)
0069         continue;
0070       // here we parallelize
0071       for (int32_t ic = firstX; ic < sg - 1; ic += blockDim.x) {
0072         auto& ci = cells[cc[ic]];
0073         for (auto jc = ic + 1; jc < sg; ++jc) {
0074           auto& cj = cells[cc[jc]];
0075           // must be different detectors
0076           //        if (d[ic]==d[jc]) continue;
0077           auto cos12 = x[ic] * x[jc] + y[ic] * y[jc] + z[ic] * z[jc];
0078           if (d[ic] != d[jc] && cos12 * cos12 >= 0.99999f * (n[ic] * n[jc])) {
0079             // alligned:  kill farthest (prefer consecutive layers)
0080             // if same layer prefer farthest (longer level arm) and make space for intermediate hit
0081             bool sameLayer = l[ic] == l[jc];
0082             if (n[ic] > n[jc]) {
0083               if (sameLayer) {
0084                 cj.kill();  // closest
0085                 ci.setFishbone(cj.inner_hit_id(), cj.inner_z(hh), hh);
0086               } else {
0087                 ci.kill();  // farthest
0088                 // break;  // removed to improve reproducibility. keep it for reference and tests
0089               }
0090             } else {
0091               if (!sameLayer) {
0092                 cj.kill();  // farthest
0093               } else {
0094                 ci.kill();  // closest
0095                 cj.setFishbone(ci.inner_hit_id(), ci.inner_z(hh), hh);
0096                 // break;  // removed to improve reproducibility. keep it for reference    and tests
0097               }
0098             }
0099           }
0100         }  // cj
0101       }    // ci
0102     }      // hits
0103   }
0104 
0105 }  // namespace gpuPixelDoublets
0106 
0107 #endif  // RecoPixelVertexing_PixelTriplets_plugins_gpuFishbone_h