Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-03-14 23:57:40

0001 #include <cuda_runtime.h>
0002 
0003 #include "CUDADataFormats/BeamSpot/interface/BeamSpotCUDA.h"
0004 #include "CUDADataFormats/SiPixelCluster/interface/SiPixelClustersCUDA.h"
0005 #include "CUDADataFormats/SiPixelDigi/interface/SiPixelDigisCUDA.h"
0006 #include "CUDADataFormats/TrackingRecHit/interface/TrackingRecHit2DHeterogeneous.h"
0007 #include "CUDADataFormats/Common/interface/HostProduct.h"
0008 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0009 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0012 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 #include "FWCore/Framework/interface/MakerMacros.h"
0016 #include "FWCore/Framework/interface/global/EDProducer.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/Records/interface/TrackerDigiGeometryRecord.h"
0022 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0023 #include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
0024 #include "RecoLocalTracker/Records/interface/TkPixelCPERecord.h"
0025 #include "RecoLocalTracker/SiPixelRecHits/interface/PixelCPEBase.h"
0026 #include "RecoLocalTracker/SiPixelRecHits/interface/PixelCPEFast.h"
0027 
0028 #include "gpuPixelRecHits.h"
0029 
0030 class SiPixelRecHitSoAFromLegacy : public edm::global::EDProducer<> {
0031 public:
0032   explicit SiPixelRecHitSoAFromLegacy(const edm::ParameterSet& iConfig);
0033   ~SiPixelRecHitSoAFromLegacy() override = default;
0034 
0035   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0036 
0037   using HitModuleStart = std::array<uint32_t, gpuClustering::maxNumModules + 1>;
0038   using HMSstorage = HostProduct<uint32_t[]>;
0039 
0040 private:
0041   void produce(edm::StreamID streamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const override;
0042 
0043   const edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geomToken_;
0044   const edm::ESGetToken<PixelClusterParameterEstimator, TkPixelCPERecord> cpeToken_;
0045   const edm::EDGetTokenT<reco::BeamSpot> bsGetToken_;
0046   const edm::EDGetTokenT<SiPixelClusterCollectionNew> clusterToken_;  // Legacy Clusters
0047   const edm::EDPutTokenT<TrackingRecHit2DCPU> tokenHit_;
0048   const edm::EDPutTokenT<HMSstorage> tokenModuleStart_;
0049   const bool convert2Legacy_;
0050   const bool isPhase2_;
0051 };
0052 
0053 SiPixelRecHitSoAFromLegacy::SiPixelRecHitSoAFromLegacy(const edm::ParameterSet& iConfig)
0054     : geomToken_(esConsumes()),
0055       cpeToken_(esConsumes(edm::ESInputTag("", iConfig.getParameter<std::string>("CPE")))),
0056       bsGetToken_{consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamSpot"))},
0057       clusterToken_{consumes<SiPixelClusterCollectionNew>(iConfig.getParameter<edm::InputTag>("src"))},
0058       tokenHit_{produces<TrackingRecHit2DCPU>()},
0059       tokenModuleStart_{produces<HMSstorage>()},
0060       convert2Legacy_(iConfig.getParameter<bool>("convertToLegacy")),
0061       isPhase2_(iConfig.getParameter<bool>("isPhase2")) {
0062   if (convert2Legacy_)
0063     produces<SiPixelRecHitCollectionNew>();
0064 }
0065 
0066 void SiPixelRecHitSoAFromLegacy::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0067   edm::ParameterSetDescription desc;
0068 
0069   desc.add<edm::InputTag>("beamSpot", edm::InputTag("offlineBeamSpot"));
0070   desc.add<edm::InputTag>("src", edm::InputTag("siPixelClustersPreSplitting"));
0071   desc.add<std::string>("CPE", "PixelCPEFast");
0072   desc.add<bool>("convertToLegacy", false);
0073   desc.add<bool>("isPhase2", false);
0074   descriptions.addWithDefaultLabel(desc);
0075 }
0076 
0077 void SiPixelRecHitSoAFromLegacy::produce(edm::StreamID streamID, edm::Event& iEvent, const edm::EventSetup& es) const {
0078   const TrackerGeometry* geom_ = &es.getData(geomToken_);
0079   PixelCPEFast const* fcpe = dynamic_cast<const PixelCPEFast*>(&es.getData(cpeToken_));
0080   if (not fcpe) {
0081     throw cms::Exception("Configuration") << "SiPixelRecHitSoAFromLegacy can only use a CPE of type PixelCPEFast";
0082   }
0083   auto const& cpeView = fcpe->getCPUProduct();
0084 
0085   const reco::BeamSpot& bs = iEvent.get(bsGetToken_);
0086 
0087   BeamSpotPOD bsHost;
0088   bsHost.x = bs.x0();
0089   bsHost.y = bs.y0();
0090   bsHost.z = bs.z0();
0091 
0092   edm::Handle<SiPixelClusterCollectionNew> hclusters;
0093   iEvent.getByToken(clusterToken_, hclusters);
0094   auto const& input = *hclusters;
0095 
0096   const int nMaxModules = isPhase2_ ? phase2PixelTopology::numberOfModules : phase1PixelTopology::numberOfModules;
0097   const int startBPIX2 = isPhase2_ ? phase2PixelTopology::layerStart[1] : phase1PixelTopology::layerStart[1];
0098 
0099   assert(nMaxModules < gpuClustering::maxNumModules);
0100   assert(startBPIX2 < nMaxModules);
0101 
0102   // allocate a buffer for the indices of the clusters
0103   auto hmsp = std::make_unique<uint32_t[]>(nMaxModules + 1);
0104   // hitsModuleStart is a non-owning pointer to the buffer
0105   auto hitsModuleStart = hmsp.get();
0106   // wrap the buffer in a HostProduct
0107   auto hms = std::make_unique<HMSstorage>(std::move(hmsp));
0108   // move the HostProduct to the Event, without reallocating the buffer or affecting hitsModuleStart
0109   iEvent.put(tokenModuleStart_, std::move(hms));
0110 
0111   // legacy output
0112   auto legacyOutput = std::make_unique<SiPixelRecHitCollectionNew>();
0113 
0114   // storage
0115   std::vector<uint16_t> xx;
0116   std::vector<uint16_t> yy;
0117   std::vector<uint16_t> adc;
0118   std::vector<uint16_t> moduleInd;
0119   std::vector<int32_t> clus;
0120 
0121   std::vector<edm::Ref<edmNew::DetSetVector<SiPixelCluster>, SiPixelCluster>> clusterRef;
0122 
0123   constexpr uint32_t maxHitsInModule = gpuClustering::maxHitsInModule();
0124 
0125   HitModuleStart moduleStart_;  // index of the first pixel of each module
0126   HitModuleStart clusInModule_;
0127   memset(&clusInModule_, 0, sizeof(HitModuleStart));  // needed??
0128   memset(&moduleStart_, 0, sizeof(HitModuleStart));
0129   assert(gpuClustering::maxNumModules + 1 == clusInModule_.size());
0130   assert(0 == clusInModule_[gpuClustering::maxNumModules]);
0131   uint32_t moduleId_;
0132   moduleStart_[1] = 0;  // we run sequentially....
0133 
0134   SiPixelClustersCUDA::SiPixelClustersCUDASOAView clusterView{
0135       moduleStart_.data(), clusInModule_.data(), &moduleId_, hitsModuleStart};
0136 
0137   // fill cluster arrays
0138   int numberOfClusters = 0;
0139   for (auto const& dsv : input) {
0140     unsigned int detid = dsv.detId();
0141     DetId detIdObject(detid);
0142     const GeomDetUnit* genericDet = geom_->idToDetUnit(detIdObject);
0143     auto gind = genericDet->index();
0144     assert(gind < nMaxModules);
0145     auto const nclus = dsv.size();
0146     clusInModule_[gind] = nclus;
0147     numberOfClusters += nclus;
0148   }
0149   hitsModuleStart[0] = 0;
0150 
0151   for (int i = 1, n = nMaxModules + 1; i < n; ++i)
0152     hitsModuleStart[i] = hitsModuleStart[i - 1] + clusInModule_[i - 1];
0153 
0154   assert(numberOfClusters == int(hitsModuleStart[nMaxModules]));
0155 
0156   // output SoA
0157   // element 96 is the start of BPIX2 (i.e. the number of clusters in BPIX1)
0158 
0159   auto output = std::make_unique<TrackingRecHit2DCPU>(
0160       numberOfClusters, isPhase2_, hitsModuleStart[startBPIX2], &cpeView, hitsModuleStart, nullptr);
0161   assert(output->nMaxModules() == uint32_t(nMaxModules));
0162 
0163   if (0 == numberOfClusters) {
0164     iEvent.put(std::move(output));
0165     if (convert2Legacy_)
0166       iEvent.put(std::move(legacyOutput));
0167     return;
0168   }
0169 
0170   if (convert2Legacy_)
0171     legacyOutput->reserve(nMaxModules, numberOfClusters);
0172 
0173   int numberOfDetUnits = 0;
0174   int numberOfHits = 0;
0175   for (auto const& dsv : input) {
0176     numberOfDetUnits++;
0177     unsigned int detid = dsv.detId();
0178     DetId detIdObject(detid);
0179     const GeomDetUnit* genericDet = geom_->idToDetUnit(detIdObject);
0180     auto const gind = genericDet->index();
0181     assert(gind < nMaxModules);
0182     const PixelGeomDetUnit* pixDet = dynamic_cast<const PixelGeomDetUnit*>(genericDet);
0183     assert(pixDet);
0184     auto const nclus = dsv.size();
0185     assert(clusInModule_[gind] == nclus);
0186     if (0 == nclus)
0187       continue;  // is this really possible?
0188 
0189     auto const fc = hitsModuleStart[gind];
0190     auto const lc = hitsModuleStart[gind + 1];
0191     assert(lc > fc);
0192     LogDebug("SiPixelRecHitSoAFromLegacy") << "in det " << gind << ": conv " << nclus << " hits from " << dsv.size()
0193                                            << " legacy clusters" << ' ' << fc << ',' << lc;
0194     assert((lc - fc) == nclus);
0195     if (nclus > maxHitsInModule)
0196       printf(
0197           "WARNING: too many clusters %d in Module %d. Only first %d Hits converted\n", nclus, gind, maxHitsInModule);
0198 
0199     // fill digis
0200     xx.clear();
0201     yy.clear();
0202     adc.clear();
0203     moduleInd.clear();
0204     clus.clear();
0205     clusterRef.clear();
0206     moduleId_ = gind;
0207     uint32_t ic = 0;
0208     uint32_t ndigi = 0;
0209     for (auto const& clust : dsv) {
0210       assert(clust.size() > 0);
0211       for (int i = 0, nd = clust.size(); i < nd; ++i) {
0212         auto px = clust.pixel(i);
0213         xx.push_back(px.x);
0214         yy.push_back(px.y);
0215         adc.push_back(px.adc);
0216         moduleInd.push_back(gind);
0217         clus.push_back(ic);
0218         ++ndigi;
0219       }
0220 
0221       if (convert2Legacy_)
0222         clusterRef.emplace_back(edmNew::makeRefTo(hclusters, &clust));
0223       ic++;
0224     }
0225     assert(nclus == ic);
0226     assert(clus.size() == ndigi);
0227     numberOfHits += nclus;
0228     // filled creates view
0229     SiPixelDigisCUDASOAView digiView;
0230     digiView.xx_ = xx.data();
0231     digiView.yy_ = yy.data();
0232     digiView.adc_ = adc.data();
0233     digiView.moduleInd_ = moduleInd.data();
0234     digiView.clus_ = clus.data();
0235     digiView.pdigi_ = nullptr;
0236     digiView.rawIdArr_ = nullptr;
0237     assert(digiView.adc(0) != 0);
0238     // we run on blockId.x==0
0239     gpuPixelRecHits::getHits(&cpeView, &bsHost, digiView, ndigi, &clusterView, output->view());
0240     for (auto h = fc; h < lc; ++h)
0241       if (h - fc < maxHitsInModule)
0242         assert(gind == output->view()->detectorIndex(h));
0243       else
0244         assert(gpuClustering::invalidModuleId == output->view()->detectorIndex(h));
0245     if (convert2Legacy_) {
0246       SiPixelRecHitCollectionNew::FastFiller recHitsOnDetUnit(*legacyOutput, detid);
0247       for (auto h = fc; h < lc; ++h) {
0248         auto ih = h - fc;
0249 
0250         if (ih >= maxHitsInModule)
0251           break;
0252         assert(ih < clusterRef.size());
0253         LocalPoint lp(output->view()->xLocal(h), output->view()->yLocal(h));
0254         LocalError le(output->view()->xerrLocal(h), 0, output->view()->yerrLocal(h));
0255         SiPixelRecHitQuality::QualWordType rqw = 0;
0256         SiPixelRecHit hit(lp, le, rqw, *genericDet, clusterRef[ih]);
0257         recHitsOnDetUnit.push_back(hit);
0258       }
0259     }
0260   }
0261 
0262   assert(numberOfHits == numberOfClusters);
0263 
0264   // fill data structure to support CA
0265   const auto nLayers = isPhase2_ ? phase2PixelTopology::numberOfLayers : phase1PixelTopology::numberOfLayers;
0266   for (auto i = 0U; i < nLayers + 1; ++i) {
0267     output->hitsLayerStart()[i] = hitsModuleStart[cpeView.layerGeometry().layerStart[i]];
0268     LogDebug("SiPixelRecHitSoAFromLegacy")
0269         << "Layer n." << i << " - starting at module: " << cpeView.layerGeometry().layerStart[i]
0270         << " - starts ad cluster: " << output->hitsLayerStart()[i] << "\n";
0271   }
0272 
0273   cms::cuda::fillManyFromVector(output->phiBinner(),
0274                                 nLayers,
0275                                 output->iphi(),
0276                                 output->hitsLayerStart(),
0277                                 numberOfHits,
0278                                 256,
0279                                 output->phiBinnerStorage());
0280 
0281   LogDebug("SiPixelRecHitSoAFromLegacy") << "created HitSoa for " << numberOfClusters << " clusters in "
0282                                          << numberOfDetUnits << " Dets";
0283   iEvent.put(std::move(output));
0284   if (convert2Legacy_)
0285     iEvent.put(std::move(legacyOutput));
0286 }
0287 
0288 DEFINE_FWK_MODULE(SiPixelRecHitSoAFromLegacy);