diff --git a/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelPhase2DigiToCluster.cc b/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelPhase2DigiToCluster.cc index 9b9a36f207d76..198837c14a9db 100644 --- a/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelPhase2DigiToCluster.cc +++ b/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelPhase2DigiToCluster.cc @@ -13,13 +13,17 @@ #include "DataFormats/SiPixelDigi/interface/PixelDigi.h" #include "DataFormats/SiPixelDigiSoA/interface/alpaka/SiPixelDigiErrorsSoACollection.h" #include "DataFormats/SiPixelDigiSoA/interface/alpaka/SiPixelDigisSoACollection.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/ParameterSet/interface/ParameterSetDescription.h" #include "FWCore/Utilities/interface/ESGetToken.h" #include "FWCore/Utilities/interface/InputTag.h" +#include "Geometry/CommonTopologies/interface/GeomDetEnumerators.h" #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" +#include "Geometry/Records/interface/TrackerTopologyRcd.h" #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" +#include "DataFormats/TrackerCommon/interface/TrackerTopology.h" #include "HeterogeneousCore/AlpakaCore/interface/alpaka/EDPutToken.h" #include "HeterogeneousCore/AlpakaCore/interface/alpaka/Event.h" #include "HeterogeneousCore/AlpakaCore/interface/alpaka/stream/SynchronizingEDProducer.h" @@ -41,8 +45,12 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { private: void acquire(device::Event const& iEvent, device::EventSetup const& iSetup) override; void produce(device::Event& iEvent, device::EventSetup const& iSetup) override; + void beginRun(edm::Run const&, edm::EventSetup const& iSetup) override; const edm::ESGetToken geomToken_; + const edm::ESGetToken geomTokenBeginRun_; // For BeginRun + const edm::ESGetToken topoToken_; + const edm::EDGetTokenT> pixelDigiToken_; const device::EDPutToken digiPutToken_; const device::EDPutToken clusterPutToken_; @@ -51,11 +59,14 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { Algo algo_; uint32_t nDigis_ = 0; std::optional digis_d_; + mutable uint32_t offsetBPIX2_ = pixelTopology::Phase2::layerStart[1]; }; SiPixelPhase2DigiToCluster::SiPixelPhase2DigiToCluster(const edm::ParameterSet& iConfig) : SynchronizingEDProducer(iConfig), geomToken_(esConsumes()), + geomTokenBeginRun_(esConsumes()), + topoToken_(esConsumes()), pixelDigiToken_(consumes>(iConfig.getParameter("InputDigis"))), digiPutToken_(produces()), clusterPutToken_(produces()), @@ -80,6 +91,44 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { desc.add("InputDigis", edm::InputTag("simSiPixelDigis:Pixel")); descriptions.addWithDefaultLabel(desc); } + void SiPixelPhase2DigiToCluster::beginRun(edm::Run const&, edm::EventSetup const& iSetup) { + using namespace pixelTopology; + + auto const& trackerGeometry = iSetup.getData(geomTokenBeginRun_); + auto const& trackerTopology = iSetup.getData(topoToken_); + + auto const& dets = trackerGeometry.detUnits(); + + uint32_t n_modules = 0; + uint32_t oldLayer = std::numeric_limits::max(); + uint32_t layerCount = 0; + uint32_t bpix2Start = 0; + + // Loop over detector modules to find where BPIX2 starts + for (auto& det : dets) { + if (!GeomDetEnumerators::isInnerTracker(det->subDetector())) + continue; + + DetId detId = det->geographicalId(); + auto layer = trackerTopology.layer(detId); + + if (layer != oldLayer) { + if (layerCount == 1) { + // layer 1 is BPIX2 + bpix2Start = n_modules; + } + layerCount++; + oldLayer = layer; + } + n_modules++; + } + + offsetBPIX2_ = bpix2Start; + + LogDebug("SiPixelPhase2DigiToCluster") + << "beginRun: BPIX2 module start = " << offsetBPIX2_ << " (total pixel modules: " << n_modules + << "). Offset from simplePixelTopology = " << pixelTopology::Phase2::layerStart[1] << '\n'; + } void SiPixelPhase2DigiToCluster::acquire(device::Event const& iEvent, device::EventSetup const& iSetup) { auto const& input = iEvent.get(pixelDigiToken_); @@ -117,7 +166,7 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { assert(nDigis == nDigis_); alpaka::memcpy(iEvent.queue(), digis_d_->buffer(), digis_h.buffer()); - algo_.makePhase2ClustersAsync(iEvent.queue(), clusterThresholds_, digis_d_->view(), nDigis_); + algo_.makePhase2ClustersAsync(iEvent.queue(), clusterThresholds_, digis_d_->view(), nDigis_, offsetBPIX2_); } void SiPixelPhase2DigiToCluster::produce(device::Event& iEvent, device::EventSetup const& iSetup) { diff --git a/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.dev.cc b/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.dev.cc index 884cd66643bf7..10cde5a8f3639 100644 --- a/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.dev.cc +++ b/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.dev.cc @@ -661,7 +661,8 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { Queue &queue, const SiPixelClusterThresholds clusterThresholds, SiPixelDigisSoAView &digis_view, - const uint32_t numDigis) { + const uint32_t numDigis, + const uint32_t offsetBPIX2) { using namespace pixelClustering; using pixelTopology::Phase2; nDigis = numDigis; @@ -744,7 +745,7 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { // last element holds the number of all clusters const auto clusModuleStartLastElement = cms::alpakatools::make_device_view( queue, clusters_d->const_view().clusModuleStart().data() + numberOfModules, 1u); - constexpr int startBPIX2 = pixelTopology::Phase2::layerStart[1]; + const int startBPIX2 = offsetBPIX2; // element startBPIX2 hold the number of clusters until BPIX2 const auto bpix2ClusterStart = cms::alpakatools::make_device_view(queue, clusters_d->const_view().clusModuleStart().data() + startBPIX2, 1u); diff --git a/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.h b/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.h index 05cf94a18fac7..e8e56a5c1fa9b 100644 --- a/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.h +++ b/RecoLocalTracker/SiPixelClusterizer/plugins/alpaka/SiPixelRawToClusterKernel.h @@ -168,7 +168,8 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE { void makePhase2ClustersAsync(Queue& queue, const SiPixelClusterThresholds clusterThresholds, SiPixelDigisSoAView& digis_view, - const uint32_t numDigis); + const uint32_t numDigis, + const uint32_t offsetBPIX2); SiPixelDigisSoACollection getDigis() { digis_d->setNModules(nModules_Clusters_h[0]);