File indexing completed on 2025-02-05 03:15:14
0001
0002 #include <memory>
0003 #include <optional>
0004 #include <stdexcept>
0005 #include <string>
0006 #include <utility>
0007 #include <vector>
0008
0009 #include <alpaka/alpaka.hpp>
0010
0011 #include "DataFormats/Common/interface/DetSetVector.h"
0012 #include "DataFormats/SiPixelClusterSoA/interface/alpaka/SiPixelClustersSoACollection.h"
0013 #include "DataFormats/SiPixelDigi/interface/PixelDigi.h"
0014 #include "DataFormats/SiPixelDigiSoA/interface/alpaka/SiPixelDigiErrorsSoACollection.h"
0015 #include "DataFormats/SiPixelDigiSoA/interface/alpaka/SiPixelDigisSoACollection.h"
0016 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0019 #include "FWCore/Utilities/interface/ESGetToken.h"
0020 #include "FWCore/Utilities/interface/InputTag.h"
0021 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0022 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0023 #include "HeterogeneousCore/AlpakaCore/interface/alpaka/EDPutToken.h"
0024 #include "HeterogeneousCore/AlpakaCore/interface/alpaka/Event.h"
0025 #include "HeterogeneousCore/AlpakaCore/interface/alpaka/stream/SynchronizingEDProducer.h"
0026 #include "HeterogeneousCore/AlpakaInterface/interface/config.h"
0027 #include "RecoLocalTracker/SiPixelClusterizer/interface/SiPixelClusterThresholds.h"
0028
0029 #include "SiPixelRawToClusterKernel.h"
0030
0031 namespace ALPAKA_ACCELERATOR_NAMESPACE {
0032
0033 class SiPixelPhase2DigiToCluster : public stream::SynchronizingEDProducer<> {
0034 public:
0035 explicit SiPixelPhase2DigiToCluster(const edm::ParameterSet& iConfig);
0036 ~SiPixelPhase2DigiToCluster() override = default;
0037
0038 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0039 using Algo = pixelDetails::SiPixelRawToClusterKernel<pixelTopology::Phase2>;
0040
0041 private:
0042 void acquire(device::Event const& iEvent, device::EventSetup const& iSetup) override;
0043 void produce(device::Event& iEvent, device::EventSetup const& iSetup) override;
0044
0045 const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0046 const edm::EDGetTokenT<edm::DetSetVector<PixelDigi>> pixelDigiToken_;
0047 const device::EDPutToken<SiPixelDigisSoACollection> digiPutToken_;
0048 const device::EDPutToken<SiPixelClustersSoACollection> clusterPutToken_;
0049 const SiPixelClusterThresholds clusterThresholds_;
0050
0051 Algo algo_;
0052 uint32_t nDigis_ = 0;
0053 std::optional<SiPixelDigisSoACollection> digis_d_;
0054 };
0055
0056 SiPixelPhase2DigiToCluster::SiPixelPhase2DigiToCluster(const edm::ParameterSet& iConfig)
0057 : SynchronizingEDProducer(iConfig),
0058 geomToken_(esConsumes()),
0059 pixelDigiToken_(consumes<edm::DetSetVector<PixelDigi>>(iConfig.getParameter<edm::InputTag>("InputDigis"))),
0060 digiPutToken_(produces()),
0061 clusterPutToken_(produces()),
0062 clusterThresholds_{iConfig.getParameter<int32_t>("clusterThreshold_layer1"),
0063 iConfig.getParameter<int32_t>("clusterThreshold_otherLayers"),
0064 static_cast<float>(iConfig.getParameter<double>("ElectronPerADCGain")),
0065 static_cast<int8_t>(iConfig.getParameter<int>("Phase2ReadoutMode")),
0066 static_cast<uint16_t>(iConfig.getParameter<uint32_t>("Phase2DigiBaseline")),
0067 static_cast<uint8_t>(iConfig.getParameter<uint32_t>("Phase2KinkADC"))} {}
0068
0069 void SiPixelPhase2DigiToCluster::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0070 edm::ParameterSetDescription desc;
0071
0072 desc.add<bool>("IncludeErrors", true);
0073 desc.add<int32_t>("clusterThreshold_layer1",
0074 pixelClustering::clusterThresholdPhase2LayerOne);
0075 desc.add<int32_t>("clusterThreshold_otherLayers", pixelClustering::clusterThresholdPhase2OtherLayers);
0076 desc.add<double>("ElectronPerADCGain", 1500.);
0077 desc.add<int32_t>("Phase2ReadoutMode", 3);
0078 desc.add<uint32_t>("Phase2DigiBaseline", 1000);
0079 desc.add<uint32_t>("Phase2KinkADC", 8);
0080 desc.add<edm::InputTag>("InputDigis", edm::InputTag("simSiPixelDigis:Pixel"));
0081 descriptions.addWithDefaultLabel(desc);
0082 }
0083
0084 void SiPixelPhase2DigiToCluster::acquire(device::Event const& iEvent, device::EventSetup const& iSetup) {
0085 auto const& input = iEvent.get(pixelDigiToken_);
0086
0087 const TrackerGeometry* geom_ = &iSetup.getData(geomToken_);
0088
0089 nDigis_ = 0;
0090 for (const auto& det : input) {
0091 nDigis_ += det.size();
0092 }
0093 digis_d_ = SiPixelDigisSoACollection(nDigis_, iEvent.queue());
0094
0095 if (nDigis_ == 0)
0096 return;
0097
0098 SiPixelDigisHost digis_h(nDigis_, iEvent.queue());
0099
0100 uint32_t nDigis = 0;
0101 for (const auto& det : input) {
0102 unsigned int detid = det.detId();
0103 DetId detIdObject(detid);
0104 const GeomDetUnit* genericDet = geom_->idToDetUnit(detIdObject);
0105 auto const gind = genericDet->index();
0106 for (auto const& px : det) {
0107 digis_h.view()[nDigis].moduleId() = uint16_t(gind);
0108 digis_h.view()[nDigis].xx() = uint16_t(px.row());
0109 digis_h.view()[nDigis].yy() = uint16_t(px.column());
0110 digis_h.view()[nDigis].adc() = uint16_t(px.adc());
0111 digis_h.view()[nDigis].clus() = 0;
0112 digis_h.view()[nDigis].pdigi() = uint32_t(px.packedData());
0113 digis_h.view()[nDigis].rawIdArr() = uint32_t(detid);
0114 ++nDigis;
0115 }
0116 }
0117 assert(nDigis == nDigis_);
0118
0119 alpaka::memcpy(iEvent.queue(), digis_d_->buffer(), digis_h.buffer());
0120 algo_.makePhase2ClustersAsync(iEvent.queue(), clusterThresholds_, digis_d_->view(), nDigis_);
0121 }
0122
0123 void SiPixelPhase2DigiToCluster::produce(device::Event& iEvent, device::EventSetup const& iSetup) {
0124 if (nDigis_ == 0) {
0125 iEvent.emplace(digiPutToken_, std::move(*digis_d_));
0126 iEvent.emplace(clusterPutToken_, pixelTopology::Phase2::numberOfModules, iEvent.queue());
0127 } else {
0128 digis_d_->setNModules(algo_.nModules());
0129 iEvent.emplace(digiPutToken_, std::move(*digis_d_));
0130 iEvent.emplace(clusterPutToken_, algo_.getClusters());
0131 }
0132 digis_d_.reset();
0133 }
0134
0135 }
0136
0137
0138 #include "HeterogeneousCore/AlpakaCore/interface/alpaka/MakerMacros.h"
0139 DEFINE_FWK_ALPAKA_MODULE(SiPixelPhase2DigiToCluster);