pixelClustering::ClusterChargeCut< TrackerTraits > Struct Template Reference

#include <ClusterChargeCut.h>

Public Member Functions
template<typename TAcc >
ALPAKA_FN_ACC void	operator() (TAcc const &acc, SiPixelDigisSoAView digi_view, SiPixelClustersSoAView clus_view, SiPixelClusterThresholds clusterThresholds, const uint32_t numElements) const

Detailed Description

template<typename TrackerTraits>
struct pixelClustering::ClusterChargeCut< TrackerTraits >

Definition at line 20 of file ClusterChargeCut.h.

Member Function Documentation

operator()()

template<typename TrackerTraits >

template<typename TAcc >

ALPAKA_FN_ACC void pixelClustering::ClusterChargeCut< TrackerTraits >::operator() ( TAcc const &  acc, SiPixelDigisSoAView  digi_view, SiPixelClustersSoAView  clus_view, SiPixelClusterThresholds  clusterThresholds, const uint32_t  numElements  ) const

inline

Definition at line 22 of file ClusterChargeCut.h.

References gpuClustering::adc, ALPAKA_ACCELERATOR_NAMESPACE::caPixelDoublets::ALPAKA_ASSERT_ACC(), cms::cudacompat::atomicAdd(), cms::alpakatools::blockPrefixScan(), ALCARECOTkAlJpsiMuMu_cff::charge, DMR_cfg::chargeCut, ALPAKA_ACCELERATOR_NAMESPACE::brokenline::constexpr(), gpuClustering::endModule, SiPixelClusterThresholds::getThresholdForLayerOnCondition(), caHitNtupletGeneratorKernels::good, mps_fire::i, cms::alpakatools::independent_group_elements(), cms::alpakatools::independent_groups(), pixelClustering::invalidModuleId, phase1PixelTopology::layerStart, pixelClustering::maxNumClustersPerModules, pixelClustering::maxNumModules, SiStripPI::min, callgraph::module, gpuClustering::moduleId, gpuClustering::newclusId, phase1PixelTopology::numberOfModules, gpuClustering::numElements, hltrates_dqm_sourceclient-live_cfg::offset, convertSQLiteXML::ok, cms::alpakatools::once_per_block(), gpuClustering::startBPIX2, and ALPAKA_ACCELERATOR_NAMESPACE::vertexFinder::ws.

                                                                    {
      constexpr int32_t maxNumClustersPerModules = TrackerTraits::maxNumClustersPerModules;

      auto& charge = alpaka::declareSharedVar<int32_t[maxNumClustersPerModules], __COUNTER__>(acc);
      auto& ok = alpaka::declareSharedVar<uint8_t[maxNumClustersPerModules], __COUNTER__>(acc);
      auto& newclusId = alpaka::declareSharedVar<uint16_t[maxNumClustersPerModules], __COUNTER__>(acc);

      constexpr int startBPIX2 = TrackerTraits::layerStart[1];

      ALPAKA_ASSERT_ACC(TrackerTraits::numberOfModules < maxNumModules);
      ALPAKA_ASSERT_ACC(startBPIX2 < TrackerTraits::numberOfModules);

      auto endModule = clus_view[0].moduleStart();
      for (auto module : cms::alpakatools::independent_groups(acc, endModule)) {
        auto firstPixel = clus_view[1 + module].moduleStart();
        auto thisModuleId = digi_view[firstPixel].moduleId();
        while (thisModuleId == invalidModuleId and firstPixel < numElements) {
          // skip invalid or duplicate pixels
          ++firstPixel;
          thisModuleId = digi_view[firstPixel].moduleId();
        }
        if (firstPixel >= numElements) {
          // reached the end of the input while skipping the invalid pixels, nothing left to do
          break;
        }
        if (thisModuleId != clus_view[module].moduleId()) {
          // reached the end of the module while skipping the invalid pixels, skip this module
          continue;
        }
        ALPAKA_ASSERT_ACC(thisModuleId < TrackerTraits::numberOfModules);

        uint32_t nclus = clus_view[thisModuleId].clusInModule();
        if (nclus == 0)
          return;

        if (cms::alpakatools::once_per_block(acc) && nclus > maxNumClustersPerModules)
          printf("Warning: too many clusters in module %u in block %u: %u > %d\n",
                 thisModuleId,
                 module,
                 nclus,
                 maxNumClustersPerModules);

        if (nclus > maxNumClustersPerModules) {
          // remove excess  FIXME find a way to cut charge first....
          for (auto i : cms::alpakatools::independent_group_elements(acc, firstPixel, numElements)) {
            if (digi_view[i].moduleId() == invalidModuleId)
              continue;  // not valid
            if (digi_view[i].moduleId() != thisModuleId)
              break;  // end of module
            if (digi_view[i].clus() >= maxNumClustersPerModules) {
              digi_view[i].moduleId() = invalidModuleId;
              digi_view[i].clus() = invalidModuleId;
            }
          }
          nclus = maxNumClustersPerModules;
        }

#ifdef GPU_DEBUG
        if (thisModuleId % 100 == 1)
          if (cms::alpakatools::once_per_block(acc))
            printf("start cluster charge cut for module %d in block %d\n", thisModuleId, module);
#endif

        ALPAKA_ASSERT_ACC(nclus <= maxNumClustersPerModules);
        for (auto i : cms::alpakatools::independent_group_elements(acc, nclus)) {
          charge[i] = 0;
        }
        alpaka::syncBlockThreads(acc);

        for (auto i : cms::alpakatools::independent_group_elements(acc, firstPixel, numElements)) {
          if (digi_view[i].moduleId() == invalidModuleId)
            continue;  // not valid
          if (digi_view[i].moduleId() != thisModuleId)
            break;  // end of module
          alpaka::atomicAdd(acc,
                            &charge[digi_view[i].clus()],
                            static_cast<int32_t>(digi_view[i].adc()),
                            alpaka::hierarchy::Threads{});
        }
        alpaka::syncBlockThreads(acc);

        auto chargeCut = clusterThresholds.getThresholdForLayerOnCondition(thisModuleId < startBPIX2);

        bool good = true;
        for (auto i : cms::alpakatools::independent_group_elements(acc, nclus)) {
          newclusId[i] = ok[i] = (charge[i] >= chargeCut) ? 1 : 0;
          if (0 == ok[i])
            good = false;
        }

        // if all clusters are above threshold, do nothing
        if (alpaka::syncBlockThreadsPredicate<alpaka::BlockAnd>(acc, good))
          continue;

        // renumber
        auto& ws = alpaka::declareSharedVar<uint16_t[32], __COUNTER__>(acc);
        // FIXME this value should come from cms::alpakatools::blockPrefixScan itself
        constexpr uint32_t maxThreads = 1024;
        auto minClust = std::min(nclus, maxThreads);

        cms::alpakatools::blockPrefixScan(acc, newclusId, minClust, ws);
        if constexpr (maxNumClustersPerModules > maxThreads)  //only if needed
        {
          for (uint32_t offset = maxThreads; offset < nclus; offset += maxThreads) {
            cms::alpakatools::blockPrefixScan(acc, newclusId + offset, nclus - offset, ws);
            for (uint32_t i : cms::alpakatools::independent_group_elements(acc, offset, nclus)) {
              uint32_t prevBlockEnd = (i / maxThreads) * maxThreads - 1;
              newclusId[i] += newclusId[prevBlockEnd];
            }
            alpaka::syncBlockThreads(acc);
          }
        }
        ALPAKA_ASSERT_ACC(nclus >= newclusId[nclus - 1]);

        clus_view[thisModuleId].clusInModule() = newclusId[nclus - 1];

        // reassign id
        for (auto i : cms::alpakatools::independent_group_elements(acc, firstPixel, numElements)) {
          if (digi_view[i].moduleId() == invalidModuleId)
            continue;  // not valid
          if (digi_view[i].moduleId() != thisModuleId)
            break;  // end of module
          if (0 == ok[digi_view[i].clus()])
            digi_view[i].moduleId() = digi_view[i].clus() = invalidModuleId;
          else
            digi_view[i].clus() = newclusId[digi_view[i].clus()] - 1;
        }

        // done
        alpaka::syncBlockThreads(acc);
      }
    }