Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-06-12 23:30:06

0001 #ifndef RecoTracker_LSTCore_src_alpaka_Kernels_h
0002 #define RecoTracker_LSTCore_src_alpaka_Kernels_h
0003 
0004 #include "HeterogeneousCore/AlpakaInterface/interface/workdivision.h"
0005 #include "FWCore/Utilities/interface/CMSUnrollLoop.h"
0006 
0007 #include "RecoTracker/LSTCore/interface/alpaka/Common.h"
0008 #include "RecoTracker/LSTCore/interface/ModulesSoA.h"
0009 #include "RecoTracker/LSTCore/interface/ObjectRangesSoA.h"
0010 #include "RecoTracker/LSTCore/interface/MiniDoubletsSoA.h"
0011 #include "RecoTracker/LSTCore/interface/PixelQuintupletsSoA.h"
0012 #include "RecoTracker/LSTCore/interface/PixelTripletsSoA.h"
0013 #include "RecoTracker/LSTCore/interface/PixelSegmentsSoA.h"
0014 #include "RecoTracker/LSTCore/interface/QuintupletsSoA.h"
0015 #include "RecoTracker/LSTCore/interface/SegmentsSoA.h"
0016 #include "RecoTracker/LSTCore/interface/TripletsSoA.h"
0017 
0018 namespace ALPAKA_ACCELERATOR_NAMESPACE::lst {
0019   ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmQuintupletFromMemory(Quintuplets quintuplets,
0020                                                              unsigned int quintupletIndex,
0021                                                              bool secondpass = false) {
0022     quintuplets.isDup()[quintupletIndex] |= 1 + secondpass;
0023   }
0024 
0025   ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmPixelTripletFromMemory(PixelTriplets pixelTriplets,
0026                                                                unsigned int pixelTripletIndex) {
0027     pixelTriplets.isDup()[pixelTripletIndex] = true;
0028   }
0029 
0030   ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmPixelQuintupletFromMemory(PixelQuintuplets pixelQuintuplets,
0031                                                                   unsigned int pixelQuintupletIndex) {
0032     pixelQuintuplets.isDup()[pixelQuintupletIndex] = true;
0033   }
0034 
0035   ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmPixelSegmentFromMemory(PixelSegments pixelSegments,
0036                                                                unsigned int pixelSegmentArrayIndex,
0037                                                                bool secondpass = false) {
0038     pixelSegments.isDup()[pixelSegmentArrayIndex] |= 1 + secondpass;
0039   }
0040 
0041   ALPAKA_FN_ACC ALPAKA_FN_INLINE int checkHitsT5(unsigned int ix, unsigned int jx, QuintupletsConst quintuplets) {
0042     unsigned int hits1[Params_T5::kHits];
0043     unsigned int hits2[Params_T5::kHits];
0044 
0045     for (int i = 0; i < Params_T5::kHits; i++) {
0046       hits1[i] = quintuplets.hitIndices()[ix][i];
0047       hits2[i] = quintuplets.hitIndices()[jx][i];
0048     }
0049 
0050     int nMatched = 0;
0051     for (int i = 0; i < Params_T5::kHits; i++) {
0052       bool matched = false;
0053       for (int j = 0; j < Params_T5::kHits; j++) {
0054         if (hits1[i] == hits2[j]) {
0055           matched = true;
0056           break;
0057         }
0058       }
0059       if (matched) {
0060         nMatched++;
0061       }
0062     }
0063     return nMatched;
0064   }
0065 
0066   ALPAKA_FN_ACC ALPAKA_FN_INLINE int checkHitspT5(unsigned int ix,
0067                                                   unsigned int jx,
0068                                                   PixelQuintupletsConst pixelQuintuplets) {
0069     unsigned int hits1[Params_pT5::kHits];
0070     unsigned int hits2[Params_pT5::kHits];
0071 
0072     for (int i = 0; i < Params_pT5::kHits; i++) {
0073       hits1[i] = pixelQuintuplets.hitIndices()[ix][i];
0074       hits2[i] = pixelQuintuplets.hitIndices()[jx][i];
0075     }
0076 
0077     int nMatched = 0;
0078     for (int i = 0; i < Params_pT5::kHits; i++) {
0079       bool matched = false;
0080       for (int j = 0; j < Params_pT5::kHits; j++) {
0081         if (hits1[i] == hits2[j]) {
0082           matched = true;
0083           break;
0084         }
0085       }
0086       if (matched) {
0087         nMatched++;
0088       }
0089     }
0090     return nMatched;
0091   }
0092 
0093   ALPAKA_FN_ACC ALPAKA_FN_INLINE void checkHitspT3(unsigned int ix,
0094                                                    unsigned int jx,
0095                                                    PixelTripletsConst pixelTriplets,
0096                                                    int* matched) {
0097     int phits1[Params_pLS::kHits];
0098     int phits2[Params_pLS::kHits];
0099 
0100     for (int i = 0; i < Params_pLS::kHits; i++) {
0101       phits1[i] = pixelTriplets.hitIndices()[ix][i];
0102       phits2[i] = pixelTriplets.hitIndices()[jx][i];
0103     }
0104 
0105     int npMatched = 0;
0106     for (int i = 0; i < Params_pLS::kHits; i++) {
0107       bool pmatched = false;
0108       for (int j = 0; j < Params_pLS::kHits; j++) {
0109         if (phits1[i] == phits2[j]) {
0110           pmatched = true;
0111           break;
0112         }
0113       }
0114       if (pmatched) {
0115         npMatched++;
0116       }
0117     }
0118 
0119     int hits1[Params_T3::kHits];
0120     int hits2[Params_T3::kHits];
0121 
0122     for (int i = 0; i < Params_T3::kHits; i++) {
0123       hits1[i] = pixelTriplets.hitIndices()[ix][i + 4];  // Omitting the pLS hits
0124       hits2[i] = pixelTriplets.hitIndices()[jx][i + 4];  // Omitting the pLS hits
0125     }
0126 
0127     int nMatched = 0;
0128     for (int i = 0; i < Params_T3::kHits; i++) {
0129       bool tmatched = false;
0130       for (int j = 0; j < Params_T3::kHits; j++) {
0131         if (hits1[i] == hits2[j]) {
0132           tmatched = true;
0133           break;
0134         }
0135       }
0136       if (tmatched) {
0137         nMatched++;
0138       }
0139     }
0140 
0141     matched[0] = npMatched;
0142     matched[1] = nMatched;
0143   }
0144 
0145   struct RemoveDupQuintupletsAfterBuild {
0146     ALPAKA_FN_ACC void operator()(Acc3D const& acc,
0147                                   ModulesConst modules,
0148                                   Quintuplets quintuplets,
0149                                   QuintupletsOccupancyConst quintupletsOccupancy,
0150                                   ObjectRangesConst ranges) const {
0151       for (unsigned int lowmod : cms::alpakatools::uniform_elements_z(acc, modules.nLowerModules())) {
0152         unsigned int nQuintuplets_lowmod = quintupletsOccupancy.nQuintuplets()[lowmod];
0153         int quintupletModuleIndices_lowmod = ranges.quintupletModuleIndices()[lowmod];
0154 
0155         for (unsigned int ix1 : cms::alpakatools::uniform_elements_y(acc, nQuintuplets_lowmod)) {
0156           unsigned int ix = quintupletModuleIndices_lowmod + ix1;
0157           float eta1 = __H2F(quintuplets.eta()[ix]);
0158           float phi1 = __H2F(quintuplets.phi()[ix]);
0159           float score_rphisum1 = __H2F(quintuplets.score_rphisum()[ix]);
0160 
0161           for (unsigned int jx1 : cms::alpakatools::uniform_elements_x(acc, ix1 + 1, nQuintuplets_lowmod)) {
0162             unsigned int jx = quintupletModuleIndices_lowmod + jx1;
0163 
0164             float eta2 = __H2F(quintuplets.eta()[jx]);
0165             float phi2 = __H2F(quintuplets.phi()[jx]);
0166             float dEta = alpaka::math::abs(acc, eta1 - eta2);
0167             float dPhi = cms::alpakatools::deltaPhi(acc, phi1, phi2);
0168             float score_rphisum2 = __H2F(quintuplets.score_rphisum()[jx]);
0169 
0170             if (dEta > 0.1f)
0171               continue;
0172 
0173             if (alpaka::math::abs(acc, dPhi) > 0.1f)
0174               continue;
0175 
0176             int nMatched = checkHitsT5(ix, jx, quintuplets);
0177             const int minNHitsForDup_T5 = 7;
0178             if (nMatched >= minNHitsForDup_T5) {
0179               if (score_rphisum1 >= score_rphisum2) {
0180                 rmQuintupletFromMemory(quintuplets, ix);
0181               } else {
0182                 rmQuintupletFromMemory(quintuplets, jx);
0183               }
0184             }
0185           }
0186         }
0187       }
0188     }
0189   };
0190 
0191   struct RemoveDupQuintupletsBeforeTC {
0192     ALPAKA_FN_ACC void operator()(Acc2D const& acc,
0193                                   Quintuplets quintuplets,
0194                                   QuintupletsOccupancyConst quintupletsOccupancy,
0195                                   ObjectRangesConst ranges) const {
0196       for (unsigned int lowmodIdx1 : cms::alpakatools::uniform_elements_y(acc, ranges.nEligibleT5Modules())) {
0197         uint16_t lowmod1 = ranges.indicesOfEligibleT5Modules()[lowmodIdx1];
0198         unsigned int nQuintuplets_lowmod1 = quintupletsOccupancy.nQuintuplets()[lowmod1];
0199         if (nQuintuplets_lowmod1 == 0)
0200           continue;
0201 
0202         unsigned int quintupletModuleIndices_lowmod1 = ranges.quintupletModuleIndices()[lowmod1];
0203 
0204         for (unsigned int lowmodIdx2 :
0205              cms::alpakatools::uniform_elements_x(acc, lowmodIdx1, ranges.nEligibleT5Modules())) {
0206           uint16_t lowmod2 = ranges.indicesOfEligibleT5Modules()[lowmodIdx2];
0207           unsigned int nQuintuplets_lowmod2 = quintupletsOccupancy.nQuintuplets()[lowmod2];
0208           if (nQuintuplets_lowmod2 == 0)
0209             continue;
0210 
0211           unsigned int quintupletModuleIndices_lowmod2 = ranges.quintupletModuleIndices()[lowmod2];
0212 
0213           for (unsigned int ix1 = 0; ix1 < nQuintuplets_lowmod1; ix1 += 1) {
0214             unsigned int ix = quintupletModuleIndices_lowmod1 + ix1;
0215             if (quintuplets.isDup()[ix] & 1)
0216               continue;
0217 
0218             bool isPT5_ix = quintuplets.partOfPT5()[ix];
0219 
0220             for (unsigned int jx1 = 0; jx1 < nQuintuplets_lowmod2; jx1++) {
0221               unsigned int jx = quintupletModuleIndices_lowmod2 + jx1;
0222               if (ix == jx)
0223                 continue;
0224 
0225               if (quintuplets.isDup()[jx] & 1)
0226                 continue;
0227 
0228               bool isPT5_jx = quintuplets.partOfPT5()[jx];
0229 
0230               if (isPT5_ix && isPT5_jx)
0231                 continue;
0232 
0233               float eta1 = __H2F(quintuplets.eta()[ix]);
0234               float phi1 = __H2F(quintuplets.phi()[ix]);
0235               float score_rphisum1 = __H2F(quintuplets.score_rphisum()[ix]);
0236 
0237               float eta2 = __H2F(quintuplets.eta()[jx]);
0238               float phi2 = __H2F(quintuplets.phi()[jx]);
0239               float score_rphisum2 = __H2F(quintuplets.score_rphisum()[jx]);
0240 
0241               float dEta = alpaka::math::abs(acc, eta1 - eta2);
0242               float dPhi = cms::alpakatools::deltaPhi(acc, phi1, phi2);
0243 
0244               if (dEta > 0.1f)
0245                 continue;
0246 
0247               if (alpaka::math::abs(acc, dPhi) > 0.1f)
0248                 continue;
0249 
0250               float dR2 = dEta * dEta + dPhi * dPhi;
0251               int nMatched = checkHitsT5(ix, jx, quintuplets);
0252               const int minNHitsForDup_T5 = 5;
0253 
0254               float d2 = 0.f;
0255               CMS_UNROLL_LOOP
0256               for (unsigned int k = 0; k < Params_T5::kEmbed; ++k) {
0257                 float diff = quintuplets.t5Embed()[ix][k] - quintuplets.t5Embed()[jx][k];
0258                 d2 += diff * diff;
0259               }
0260 
0261               if (((dR2 < 0.001f || nMatched >= minNHitsForDup_T5) && d2 < 1.0f) || (dR2 < 0.02f && d2 < 0.1f)) {
0262                 if (isPT5_jx || score_rphisum1 > score_rphisum2) {
0263                   rmQuintupletFromMemory(quintuplets, ix, true);
0264                 } else if (isPT5_ix || score_rphisum1 < score_rphisum2) {
0265                   rmQuintupletFromMemory(quintuplets, jx, true);
0266                 } else {
0267                   rmQuintupletFromMemory(quintuplets, (ix < jx ? ix : jx), true);
0268                 }
0269               }
0270             }
0271           }
0272         }
0273       }
0274     }
0275   };
0276 
0277   struct RemoveDupPixelTripletsFromMap {
0278     ALPAKA_FN_ACC void operator()(Acc2D const& acc, PixelTriplets pixelTriplets) const {
0279       for (unsigned int ix : cms::alpakatools::uniform_elements_y(acc, pixelTriplets.nPixelTriplets())) {
0280         for (unsigned int jx : cms::alpakatools::uniform_elements_x(acc, pixelTriplets.nPixelTriplets())) {
0281           if (ix == jx)
0282             continue;
0283 
0284           int nMatched[2];
0285           checkHitspT3(ix, jx, pixelTriplets, nMatched);
0286           const int minNHitsForDup_pT3 = 5;
0287           if ((nMatched[0] + nMatched[1]) >= minNHitsForDup_pT3) {
0288             // Check the layers
0289             if (pixelTriplets.logicalLayers()[jx][2] < pixelTriplets.logicalLayers()[ix][2]) {
0290               rmPixelTripletFromMemory(pixelTriplets, ix);
0291               break;
0292             } else if (pixelTriplets.logicalLayers()[ix][2] == pixelTriplets.logicalLayers()[jx][2] &&
0293                        __H2F(pixelTriplets.score()[ix]) > __H2F(pixelTriplets.score()[jx])) {
0294               rmPixelTripletFromMemory(pixelTriplets, ix);
0295               break;
0296             } else if (pixelTriplets.logicalLayers()[ix][2] == pixelTriplets.logicalLayers()[jx][2] &&
0297                        (__H2F(pixelTriplets.score()[ix]) == __H2F(pixelTriplets.score()[jx])) && (ix < jx)) {
0298               rmPixelTripletFromMemory(pixelTriplets, ix);
0299               break;
0300             }
0301           }
0302         }
0303       }
0304     }
0305   };
0306 
0307   struct RemoveDupPixelQuintupletsFromMap {
0308     ALPAKA_FN_ACC void operator()(Acc2D const& acc, PixelQuintuplets pixelQuintuplets) const {
0309       unsigned int nPixelQuintuplets = pixelQuintuplets.nPixelQuintuplets();
0310       for (unsigned int ix : cms::alpakatools::uniform_elements_y(acc, nPixelQuintuplets)) {
0311         float score1 = __H2F(pixelQuintuplets.score()[ix]);
0312         for (unsigned int jx : cms::alpakatools::uniform_elements_x(acc, nPixelQuintuplets)) {
0313           if (ix == jx)
0314             continue;
0315 
0316           int nMatched = checkHitspT5(ix, jx, pixelQuintuplets);
0317           float score2 = __H2F(pixelQuintuplets.score()[jx]);
0318           const int minNHitsForDup_pT5 = 7;
0319           if (nMatched >= minNHitsForDup_pT5) {
0320             if (score1 > score2 or ((score1 == score2) and (ix > jx))) {
0321               rmPixelQuintupletFromMemory(pixelQuintuplets, ix);
0322               break;
0323             }
0324           }
0325         }
0326       }
0327     }
0328   };
0329 
0330   struct CheckHitspLS {
0331     ALPAKA_FN_ACC void operator()(Acc2D const& acc,
0332                                   ModulesConst modules,
0333                                   SegmentsOccupancyConst segmentsOccupancy,
0334                                   PixelSeedsConst pixelSeeds,
0335                                   PixelSegments pixelSegments,
0336                                   bool secondpass) const {
0337       int pixelModuleIndex = modules.nLowerModules();
0338       unsigned int nPixelSegments = segmentsOccupancy.nSegments()[pixelModuleIndex];
0339 
0340       if (nPixelSegments > n_max_pixel_segments_per_module)
0341         nPixelSegments = n_max_pixel_segments_per_module;
0342 
0343       for (unsigned int ix : cms::alpakatools::uniform_elements_y(acc, nPixelSegments)) {
0344         if (secondpass && (!pixelSeeds.isQuad()[ix] || (pixelSegments.isDup()[ix] & 1)))
0345           continue;
0346 
0347         auto const& phits1 = pixelSegments.pLSHitsIdxs()[ix];
0348         float eta_pix1 = pixelSeeds.eta()[ix];
0349         float phi_pix1 = pixelSeeds.phi()[ix];
0350 
0351         for (unsigned int jx : cms::alpakatools::uniform_elements_x(acc, ix + 1, nPixelSegments)) {
0352           float eta_pix2 = pixelSeeds.eta()[jx];
0353           float phi_pix2 = pixelSeeds.phi()[jx];
0354 
0355           if (alpaka::math::abs(acc, eta_pix2 - eta_pix1) > 0.1f)
0356             continue;
0357 
0358           if (secondpass && (!pixelSeeds.isQuad()[jx] || (pixelSegments.isDup()[jx] & 1)))
0359             continue;
0360 
0361           int8_t quad_diff = pixelSeeds.isQuad()[ix] - pixelSeeds.isQuad()[jx];
0362           float score_diff = pixelSegments.score()[ix] - pixelSegments.score()[jx];
0363           // Always keep quads over trips. If they are the same, we want the object with better score
0364           int idxToRemove;
0365           if (quad_diff > 0)
0366             idxToRemove = jx;
0367           else if (quad_diff < 0)
0368             idxToRemove = ix;
0369           else if (score_diff < 0)
0370             idxToRemove = jx;
0371           else if (score_diff > 0)
0372             idxToRemove = ix;
0373           else
0374             idxToRemove = ix;
0375 
0376           auto const& phits2 = pixelSegments.pLSHitsIdxs()[jx];
0377 
0378           int npMatched = 0;
0379           for (int i = 0; i < Params_pLS::kHits; i++) {
0380             bool pmatched = false;
0381             for (int j = 0; j < Params_pLS::kHits; j++) {
0382               if (phits1[i] == phits2[j]) {
0383                 pmatched = true;
0384                 break;
0385               }
0386             }
0387             if (pmatched) {
0388               npMatched++;
0389               // Only one hit is enough
0390               if (secondpass)
0391                 break;
0392             }
0393           }
0394           const int minNHitsForDup_pLS = 3;
0395           if (npMatched >= minNHitsForDup_pLS) {
0396             rmPixelSegmentFromMemory(pixelSegments, idxToRemove, secondpass);
0397           }
0398           if (secondpass) {
0399             float dEta = alpaka::math::abs(acc, eta_pix1 - eta_pix2);
0400             float dPhi = cms::alpakatools::deltaPhi(acc, phi_pix1, phi_pix2);
0401 
0402             float dR2 = dEta * dEta + dPhi * dPhi;
0403             if ((npMatched >= 1) || (dR2 < 1e-5f)) {
0404               rmPixelSegmentFromMemory(pixelSegments, idxToRemove, secondpass);
0405             }
0406           }
0407         }
0408       }
0409     }
0410   };
0411 }  // namespace ALPAKA_ACCELERATOR_NAMESPACE::lst
0412 #endif