plugins/alpaka/CalibPixel.h

0001 #ifndef RecoLocalTracker_SiPixelClusterizer_plugins_alpaka_CalibPixel_h
0002 #define RecoLocalTracker_SiPixelClusterizer_plugins_alpaka_CalibPixel_h
0003
0004 #include <algorithm>
0005 #include <cstdint>
0006 #include <cstdio>
0007 #include <limits>
0008
0009 #include <alpaka/alpaka.hpp>
0010
0011 #include "CondFormats/SiPixelObjects/interface/SiPixelGainCalibrationForHLTLayout.h"
0012 #include "CondFormats/SiPixelObjects/interface/alpaka/SiPixelGainCalibrationForHLTUtilities.h"
0013 #include "DataFormats/SiPixelClusterSoA/interface/ClusteringConstants.h"
0014 #include "DataFormats/SiPixelClusterSoA/interface/SiPixelClustersSoA.h"
0015 #include "DataFormats/SiPixelDigiSoA/interface/SiPixelDigisSoA.h"
0016 #include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
0017 #include "HeterogeneousCore/AlpakaInterface/interface/workdivision.h"
0018 #include "RecoLocalTracker/SiPixelClusterizer/interface/SiPixelClusterThresholds.h"
0019
0020 //#define GPU_DEBUG
0021
0022 namespace calibPixel {
0023   using namespace cms::alpakatools;
0024
0025   template <bool debug = false>
0026   struct CalibDigis {
0027     template <typename TAcc>
0028     ALPAKA_FN_ACC void operator()(const TAcc& acc,
0029                                   SiPixelClusterThresholds clusterThresholds,
0030                                   SiPixelDigisSoAView view,
0031                                   SiPixelClustersSoAView clus_view,
0032                                   const SiPixelGainCalibrationForHLTSoAConstView gains,
0033                                   int numElements) const {
0034       const float VCaltoElectronGain = clusterThresholds.vCaltoElectronGain;
0035       const float VCaltoElectronGain_L1 = clusterThresholds.vCaltoElectronGain_L1;
0036       const float VCaltoElectronOffset = clusterThresholds.vCaltoElectronOffset;
0037       const float VCaltoElectronOffset_L1 = clusterThresholds.vCaltoElectronOffset_L1;
0038
0039       // zero for next kernels...
0040       if (cms::alpakatools::once_per_grid(acc)) {
0041         clus_view[0].clusModuleStart() = 0;
0042         clus_view[0].moduleStart() = 0;
0043       }
0044       for (auto i : cms::alpakatools::uniform_elements(acc, phase1PixelTopology::numberOfModules)) {
0045         clus_view[i].clusInModule() = 0;
0046       }
0047
0048       for (auto i : cms::alpakatools::uniform_elements(acc, numElements)) {
0049         auto dvgi = view[i];
0050         if (dvgi.moduleId() == ::pixelClustering::invalidModuleId)
0051           continue;
0052
0053         bool isDeadColumn = false, isNoisyColumn = false;
0054         int row = dvgi.xx();
0055         int col = dvgi.yy();
0056         auto ret = SiPixelGainUtilities::getPedAndGain(gains, dvgi.moduleId(), col, row, isDeadColumn, isNoisyColumn);
0057         float pedestal = ret.first;
0058         float gain = ret.second;
0059         if (isDeadColumn | isNoisyColumn) {
0060           if constexpr (debug)
0061             printf("bad pixel at %d in %d\n", i, dvgi.moduleId());
0062           dvgi.moduleId() = ::pixelClustering::invalidModuleId;
0063           dvgi.adc() = 0;
0064         } else {
0065           float vcal = dvgi.adc() * gain - pedestal * gain;
0066
0067           float conversionFactor = dvgi.moduleId() < 96 ? VCaltoElectronGain_L1 : VCaltoElectronGain;
0068           float offset = dvgi.moduleId() < 96 ? VCaltoElectronOffset_L1 : VCaltoElectronOffset;
0069 #ifdef GPU_DEBUG
0070           auto old_adc = dvgi.adc();
0071 #endif
0072           dvgi.adc() = std::max(100, int(vcal * conversionFactor + offset));
0073 #ifdef GPU_DEBUG
0074           if (cms::alpakatools::once_per_grid(acc)) {
0075             printf(
0076                 "module %d pixel %d -> old_adc = %d; vcal = %.2f; conversionFactor = %.2f; offset = %.2f; new_adc = "
0077                 "%d \n",
0078                 dvgi.moduleId(),
0079                 i,
0080                 old_adc,
0081                 vcal,
0082                 conversionFactor,
0083                 offset,
0084                 dvgi.adc());
0085           }
0086 #endif
0087         }
0088       }
0089     }
0090   };
0091
0092   struct CalibDigisPhase2 {
0093     template <typename TAcc>
0094     ALPAKA_FN_ACC void operator()(const TAcc& acc,
0095                                   SiPixelClusterThresholds clusterThresholds,
0096                                   SiPixelDigisSoAView view,
0097                                   SiPixelClustersSoAView clus_view,
0098                                   int numElements) const {
0099       const float ElectronPerADCGain = clusterThresholds.electronPerADCGain;
0100       const int8_t Phase2ReadoutMode = clusterThresholds.phase2ReadoutMode;
0101       const uint16_t Phase2DigiBaseline = clusterThresholds.phase2DigiBaseline;
0102       const uint8_t Phase2KinkADC = clusterThresholds.phase2KinkADC;
0103
0104       // zero for next kernels...
0105       if (cms::alpakatools::once_per_grid(acc)) {
0106         clus_view[0].clusModuleStart() = clus_view[0].moduleStart() = 0;
0107       }
0108
0109       for (uint32_t i : cms::alpakatools::uniform_elements(acc, phase2PixelTopology::numberOfModules)) {
0110         clus_view[i].clusInModule() = 0;
0111       }
0112
0113       for (uint32_t i : cms::alpakatools::uniform_elements(acc, numElements)) {
0114         auto dvgi = view[i];
0115         if (pixelClustering::invalidModuleId != dvgi.moduleId()) {
0116           const int mode = (Phase2ReadoutMode < -1 ? -1 : Phase2ReadoutMode);
0117           int adc_int = dvgi.adc();
0118           if (mode < 0)
0119             adc_int = int(adc_int * ElectronPerADCGain);
0120           else {
0121             if (adc_int < Phase2KinkADC)
0122               adc_int = int((adc_int + 0.5) * ElectronPerADCGain);
0123             else {
0124               const int8_t dspp = (Phase2ReadoutMode < 10 ? Phase2ReadoutMode : 10);
0125               const int8_t ds = int8_t(dspp <= 1 ? 1 : (dspp - 1) * (dspp - 1));
0126               adc_int -= Phase2KinkADC;
0127               adc_int *= ds;
0128               adc_int += Phase2KinkADC;
0129               adc_int = ((adc_int + 0.5 * ds) * ElectronPerADCGain);
0130             }
0131             adc_int += int(Phase2DigiBaseline);
0132           }
0133           dvgi.adc() = std::min(adc_int, int(std::numeric_limits<uint16_t>::max()));
0134         }
0135       }
0136     }
0137   };
0138 }  // namespace calibPixel
0139
0140 #endif  // RecoLocalTracker_SiPixelClusterizer_plugins_alpaka_CalibPixel_h