Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-12-24 02:22:17

0001 #include <cuda_runtime.h>
0002 
0003 #include <fmt/printf.h>
0004 
0005 #include "CUDADataFormats/Common/interface/HostProduct.h"
0006 #include "CUDADataFormats/Common/interface/Product.h"
0007 #include "CUDADataFormats/TrackingRecHit/interface/TrackingRecHit2DHeterogeneous.h"
0008 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0009 #include "DataFormats/Common/interface/Handle.h"
0010 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0011 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/Framework/interface/MakerMacros.h"
0015 #include "FWCore/Framework/interface/stream/EDProducer.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0020 #include "FWCore/Utilities/interface/InputTag.h"
0021 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0022 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0023 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0024 #include "HeterogeneousCore/CUDACore/interface/ScopedContext.h"
0025 #include "RecoLocalTracker/SiPixelRecHits/interface/pixelCPEforGPU.h"
0026 
0027 class SiPixelRecHitFromCUDA : public edm::stream::EDProducer<edm::ExternalWork> {
0028 public:
0029   explicit SiPixelRecHitFromCUDA(const edm::ParameterSet& iConfig);
0030   ~SiPixelRecHitFromCUDA() override = default;
0031 
0032   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0033 
0034   using HMSstorage = HostProduct<uint32_t[]>;
0035 
0036 private:
0037   void acquire(edm::Event const& iEvent,
0038                edm::EventSetup const& iSetup,
0039                edm::WaitingTaskWithArenaHolder waitingTaskHolder) override;
0040   void produce(edm::Event& iEvent, edm::EventSetup const& iSetup) override;
0041 
0042   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0043   const edm::EDGetTokenT<cms::cuda::Product<TrackingRecHit2DGPU>> hitsToken_;  // CUDA hits
0044   const edm::EDGetTokenT<SiPixelClusterCollectionNew> clusterToken_;           // legacy clusters
0045   const edm::EDPutTokenT<SiPixelRecHitCollection> rechitsPutToken_;            // legacy rechits
0046   const edm::EDPutTokenT<HMSstorage> hostPutToken_;
0047 
0048   uint32_t nHits_;
0049   uint32_t nMaxModules_;
0050   cms::cuda::host::unique_ptr<float[]> store32_;
0051   cms::cuda::host::unique_ptr<uint32_t[]> hitsModuleStart_;
0052 };
0053 
0054 SiPixelRecHitFromCUDA::SiPixelRecHitFromCUDA(const edm::ParameterSet& iConfig)
0055     : geomToken_(esConsumes()),
0056       hitsToken_(
0057           consumes<cms::cuda::Product<TrackingRecHit2DGPU>>(iConfig.getParameter<edm::InputTag>("pixelRecHitSrc"))),
0058       clusterToken_(consumes<SiPixelClusterCollectionNew>(iConfig.getParameter<edm::InputTag>("src"))),
0059       rechitsPutToken_(produces<SiPixelRecHitCollection>()),
0060       hostPutToken_(produces<HMSstorage>()) {}
0061 
0062 void SiPixelRecHitFromCUDA::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0063   edm::ParameterSetDescription desc;
0064   desc.add<edm::InputTag>("pixelRecHitSrc", edm::InputTag("siPixelRecHitsPreSplittingCUDA"));
0065   desc.add<edm::InputTag>("src", edm::InputTag("siPixelClustersPreSplitting"));
0066   descriptions.addWithDefaultLabel(desc);
0067 }
0068 
0069 void SiPixelRecHitFromCUDA::acquire(edm::Event const& iEvent,
0070                                     edm::EventSetup const& iSetup,
0071                                     edm::WaitingTaskWithArenaHolder waitingTaskHolder) {
0072   cms::cuda::Product<TrackingRecHit2DGPU> const& inputDataWrapped = iEvent.get(hitsToken_);
0073   cms::cuda::ScopedContextAcquire ctx{inputDataWrapped, std::move(waitingTaskHolder)};
0074   auto const& inputData = ctx.get(inputDataWrapped);
0075 
0076   nHits_ = inputData.nHits();
0077   nMaxModules_ = inputData.nMaxModules();
0078   LogDebug("SiPixelRecHitFromCUDA") << "converting " << nHits_ << " Hits";
0079 
0080   if (0 == nHits_)
0081     return;
0082   store32_ = inputData.localCoordToHostAsync(ctx.stream());
0083   hitsModuleStart_ = inputData.hitsModuleStartToHostAsync(ctx.stream());
0084 }
0085 
0086 void SiPixelRecHitFromCUDA::produce(edm::Event& iEvent, edm::EventSetup const& es) {
0087   // allocate a buffer for the indices of the clusters
0088   auto hmsp = std::make_unique<uint32_t[]>(nMaxModules_ + 1);
0089 
0090   SiPixelRecHitCollection output;
0091   output.reserve(nMaxModules_, nHits_);
0092 
0093   if (0 == nHits_) {
0094     iEvent.emplace(rechitsPutToken_, std::move(output));
0095     iEvent.emplace(hostPutToken_, std::move(hmsp));
0096     return;
0097   }
0098   output.reserve(nMaxModules_, nHits_);
0099 
0100   std::copy(hitsModuleStart_.get(), hitsModuleStart_.get() + nMaxModules_ + 1, hmsp.get());
0101   // wrap the buffer in a HostProduct, and move it to the Event, without reallocating the buffer or affecting hitsModuleStart
0102   iEvent.emplace(hostPutToken_, std::move(hmsp));
0103 
0104   auto xl = store32_.get();
0105   auto yl = xl + nHits_;
0106   auto xe = yl + nHits_;
0107   auto ye = xe + nHits_;
0108 
0109   const TrackerGeometry* geom = &es.getData(geomToken_);
0110 
0111   edm::Handle<SiPixelClusterCollectionNew> hclusters = iEvent.getHandle(clusterToken_);
0112   auto const& input = *hclusters;
0113 
0114   constexpr uint32_t maxHitsInModule = gpuClustering::maxHitsInModule();
0115 
0116   int numberOfDetUnits = 0;
0117   int numberOfClusters = 0;
0118   for (auto const& dsv : input) {
0119     numberOfDetUnits++;
0120     unsigned int detid = dsv.detId();
0121     DetId detIdObject(detid);
0122     const GeomDetUnit* genericDet = geom->idToDetUnit(detIdObject);
0123     auto gind = genericDet->index();
0124     const PixelGeomDetUnit* pixDet = dynamic_cast<const PixelGeomDetUnit*>(genericDet);
0125     assert(pixDet);
0126     SiPixelRecHitCollection::FastFiller recHitsOnDetUnit(output, detid);
0127     auto fc = hitsModuleStart_[gind];
0128     auto lc = hitsModuleStart_[gind + 1];
0129     auto nhits = lc - fc;
0130 
0131     assert(lc > fc);
0132     LogDebug("SiPixelRecHitFromCUDA") << "in det " << gind << ": conv " << nhits << " hits from " << dsv.size()
0133                                       << " legacy clusters" << ' ' << fc << ',' << lc;
0134     if (nhits > maxHitsInModule)
0135       edm::LogWarning("SiPixelRecHitFromCUDA") << fmt::sprintf(
0136           "Too many clusters %d in module %d. Only the first %d hits will be converted", nhits, gind, maxHitsInModule);
0137     nhits = std::min(nhits, maxHitsInModule);
0138 
0139     LogDebug("SiPixelRecHitFromCUDA") << "in det " << gind << "conv " << nhits << " hits from " << dsv.size()
0140                                       << " legacy clusters" << ' ' << lc << ',' << fc;
0141 
0142     if (0 == nhits)
0143       continue;
0144     auto jnd = [&](int k) { return fc + k; };
0145     assert(nhits <= dsv.size());
0146     if (nhits != dsv.size()) {
0147       edm::LogWarning("GPUHits2CPU") << "nhits!= nclus " << nhits << ' ' << dsv.size();
0148     }
0149     for (auto const& clust : dsv) {
0150       assert(clust.originalId() >= 0);
0151       assert(clust.originalId() < dsv.size());
0152       if (clust.originalId() >= nhits)
0153         continue;
0154       auto ij = jnd(clust.originalId());
0155       LocalPoint lp(xl[ij], yl[ij]);
0156       LocalError le(xe[ij], 0, ye[ij]);
0157       SiPixelRecHitQuality::QualWordType rqw = 0;
0158 
0159       numberOfClusters++;
0160 
0161       /* cpu version....  (for reference)
0162       std::tuple<LocalPoint, LocalError, SiPixelRecHitQuality::QualWordType> tuple = cpe_->getParameters( clust, *genericDet );
0163       LocalPoint lp( std::get<0>(tuple) );
0164       LocalError le( std::get<1>(tuple) );
0165       SiPixelRecHitQuality::QualWordType rqw( std::get<2>(tuple) );
0166       */
0167 
0168       // Create a persistent edm::Ref to the cluster
0169       edm::Ref<edmNew::DetSetVector<SiPixelCluster>, SiPixelCluster> cluster = edmNew::makeRefTo(hclusters, &clust);
0170       // Make a RecHit and add it to the DetSet
0171       recHitsOnDetUnit.emplace_back(lp, le, rqw, *genericDet, cluster);
0172       // =============================
0173 
0174       LogDebug("SiPixelRecHitFromCUDA") << "cluster " << numberOfClusters << " at " << lp << ' ' << le;
0175 
0176     }  //  <-- End loop on Clusters
0177 
0178     //  LogDebug("SiPixelRecHitGPU")
0179     LogDebug("SiPixelRecHitFromCUDA") << "found " << recHitsOnDetUnit.size() << " RecHits on " << detid;
0180 
0181   }  //    <-- End loop on DetUnits
0182 
0183   LogDebug("SiPixelRecHitFromCUDA") << "found " << numberOfDetUnits << " dets, " << numberOfClusters << " clusters";
0184 
0185   iEvent.emplace(rechitsPutToken_, std::move(output));
0186 }
0187 
0188 DEFINE_FWK_MODULE(SiPixelRecHitFromCUDA);