SiPixelClusterizer/plugins/SiPixelDigisClustersFromSoA.cc

0001 #include "CUDADataFormats/SiPixelCluster/interface/gpuClusteringConstants.h"
0002 #include "DataFormats/Common/interface/DetSetVector.h"
0003 #include "DataFormats/Common/interface/Handle.h"
0004 #include "DataFormats/DetId/interface/DetId.h"
0005 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0006 #include "DataFormats/SiPixelDigi/interface/PixelDigi.h"
0007 #include "DataFormats/SiPixelDigi/interface/SiPixelDigisSoA.h"
0008 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/EventSetup.h"
0011 #include "FWCore/Framework/interface/MakerMacros.h"
0012 #include "FWCore/Framework/interface/global/EDProducer.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0017 #include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
0018 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0019 #include "RecoLocalTracker/SiPixelClusterizer/interface/SiPixelClusterThresholds.h"
0020
0021 // local include(s)
0022 #include "PixelClusterizerBase.h"
0023
0024 //#define GPU_DEBUG
0025
0026 template <typename TrackerTraits>
0027 class SiPixelDigisClustersFromSoAT : public edm::global::EDProducer<> {
0028 public:
0029   explicit SiPixelDigisClustersFromSoAT(const edm::ParameterSet& iConfig);
0030   ~SiPixelDigisClustersFromSoAT() override = default;
0031
0032   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0033
0034 private:
0035   void produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const override;
0036
0037   const edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topoToken_;
0038
0039   edm::EDGetTokenT<legacy::SiPixelDigisSoA> digiGetToken_;
0040
0041   edm::EDPutTokenT<edm::DetSetVector<PixelDigi>> digiPutToken_;
0042   edm::EDPutTokenT<SiPixelClusterCollectionNew> clusterPutToken_;
0043
0044   const SiPixelClusterThresholds clusterThresholds_;  // Cluster threshold in electrons
0045
0046   const bool produceDigis_;
0047   const bool storeDigis_;
0048 };
0049
0050 template <typename TrackerTraits>
0051 SiPixelDigisClustersFromSoAT<TrackerTraits>::SiPixelDigisClustersFromSoAT(const edm::ParameterSet& iConfig)
0052     : topoToken_(esConsumes()),
0053       digiGetToken_(consumes<legacy::SiPixelDigisSoA>(iConfig.getParameter<edm::InputTag>("src"))),
0054       clusterPutToken_(produces<SiPixelClusterCollectionNew>()),
0055       clusterThresholds_(iConfig.getParameter<int>("clusterThreshold_layer1"),
0056                          iConfig.getParameter<int>("clusterThreshold_otherLayers")),
0057       produceDigis_(iConfig.getParameter<bool>("produceDigis")),
0058       storeDigis_(iConfig.getParameter<bool>("produceDigis") && iConfig.getParameter<bool>("storeDigis")) {
0059   if (produceDigis_)
0060     digiPutToken_ = produces<edm::DetSetVector<PixelDigi>>();
0061 }
0062
0063 template <typename TrackerTraits>
0064 void SiPixelDigisClustersFromSoAT<TrackerTraits>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0065   edm::ParameterSetDescription desc;
0066   desc.add<edm::InputTag>("src", edm::InputTag("siPixelDigisSoA"));
0067   desc.add<int>("clusterThreshold_layer1", gpuClustering::clusterThresholdLayerOne);
0068   desc.add<int>("clusterThreshold_otherLayers", gpuClustering::clusterThresholdOtherLayers);
0069   desc.add<bool>("produceDigis", true);
0070   desc.add<bool>("storeDigis", true);
0071
0072   descriptions.addWithDefaultLabel(desc);
0073 }
0074
0075 template <typename TrackerTraits>
0076 void SiPixelDigisClustersFromSoAT<TrackerTraits>::produce(edm::StreamID,
0077                                                           edm::Event& iEvent,
0078                                                           const edm::EventSetup& iSetup) const {
0079   const auto& digis = iEvent.get(digiGetToken_);
0080   const uint32_t nDigis = digis.size();
0081   const auto& ttopo = iSetup.getData(topoToken_);
0082   constexpr auto maxModules = TrackerTraits::numberOfModules;
0083
0084   std::unique_ptr<edm::DetSetVector<PixelDigi>> collection;
0085   if (produceDigis_)
0086     collection = std::make_unique<edm::DetSetVector<PixelDigi>>();
0087   if (storeDigis_)
0088     collection->reserve(maxModules);
0089   auto outputClusters = std::make_unique<SiPixelClusterCollectionNew>();
0090   outputClusters->reserve(maxModules, nDigis / 2);
0091
0092   edm::DetSet<PixelDigi>* detDigis = nullptr;
0093   uint32_t detId = 0;
0094   for (uint32_t i = 0; i < nDigis; i++) {
0095     // check for uninitialized digis
0096     // this is set in RawToDigi_kernel in SiPixelRawToClusterGPUKernel.cu
0097     if (digis.rawIdArr(i) == 0)
0098       continue;
0099     // check for noisy/dead pixels (electrons set to 0)
0100     if (digis.adc(i) == 0)
0101       continue;
0102
0103     detId = digis.rawIdArr(i);
0104     if (storeDigis_) {
0105       detDigis = &collection->find_or_insert(detId);
0106       if ((*detDigis).empty())
0107         (*detDigis).data.reserve(64);  // avoid the first relocations
0108     }
0109     break;
0110   }
0111
0112   int32_t nclus = -1;
0113   PixelClusterizerBase::AccretionCluster aclusters[TrackerTraits::maxNumClustersPerModules];
0114 #ifdef EDM_ML_DEBUG
0115   auto totClustersFilled = 0;
0116 #endif
0117
0118   auto fillClusters = [&](uint32_t detId) {
0119     if (nclus < 0)
0120       return;  // this in reality should never happen
0121     edmNew::DetSetVector<SiPixelCluster>::FastFiller spc(*outputClusters, detId);
0122     auto layer = (DetId(detId).subdetId() == 1) ? ttopo.pxbLayer(detId) : 0;
0123     auto clusterThreshold = clusterThresholds_.getThresholdForLayerOnCondition(layer == 1);
0124     for (int32_t ic = 0; ic < nclus + 1; ++ic) {
0125       auto const& acluster = aclusters[ic];
0126       // in any case we cannot  go out of sync with gpu...
0127       if (acluster.charge < clusterThreshold)
0128         edm::LogWarning("SiPixelDigisClustersFromSoA") << "cluster below charge Threshold "
0129                                                        << "Layer/DetId/clusId " << layer << '/' << detId << '/' << ic
0130                                                        << " size/charge " << acluster.isize << '/' << acluster.charge;
0131       // sort by row (x)
0132       spc.emplace_back(acluster.isize, acluster.adc, acluster.x, acluster.y, acluster.xmin, acluster.ymin, ic);
0133       aclusters[ic].clear();
0134 #ifdef EDM_ML_DEBUG
0135       ++totClustersFilled;
0136       const auto& cluster{spc.back()};
0137       LogDebug("SiPixelDigisClustersFromSoA")
0138           << "putting in this cluster " << ic << " " << cluster.charge() << " " << cluster.pixelADC().size();
0139 #endif
0140       std::push_heap(spc.begin(), spc.end(), [](SiPixelCluster const& cl1, SiPixelCluster const& cl2) {
0141         return cl1.minPixelRow() < cl2.minPixelRow();
0142       });
0143     }
0144     nclus = -1;
0145     // sort by row (x)
0146     std::sort_heap(spc.begin(), spc.end(), [](SiPixelCluster const& cl1, SiPixelCluster const& cl2) {
0147       return cl1.minPixelRow() < cl2.minPixelRow();
0148     });
0149     if (spc.empty())
0150       spc.abort();
0151   };
0152
0153 #ifdef GPU_DEBUG
0154   std::cout << "Dumping all digis. nDigis = " << nDigis << std::endl;
0155 #endif
0156
0157   for (uint32_t i = 0; i < nDigis; i++) {
0158     // check for uninitialized digis
0159     if (digis.rawIdArr(i) == 0)
0160       continue;
0161     // check for noisy/dead pixels (electrons set to 0)
0162     if (digis.adc(i) == 0)
0163       continue;
0164     if (digis.clus(i) > 9000)
0165       continue;  // not in cluster; TODO add an assert for the size
0166 #ifdef EDM_ML_DEBUG
0167     assert(digis.rawIdArr(i) > 109999);
0168 #endif
0169     if (detId != digis.rawIdArr(i)) {
0170 #ifdef GPU_DEBUG
0171       std::cout << ">> Closed module --" << detId << "; nclus = " << nclus << std::endl;
0172 #endif
0173       // new module
0174       fillClusters(detId);
0175 #ifdef EDM_ML_DEBUG
0176       assert(nclus == -1);
0177 #endif
0178       detId = digis.rawIdArr(i);
0179       if (storeDigis_) {
0180         detDigis = &collection->find_or_insert(detId);
0181         if ((*detDigis).empty())
0182           (*detDigis).data.reserve(64);  // avoid the first relocations
0183         else {
0184           edm::LogWarning("SiPixelDigisClustersFromSoA")
0185               << "Problem det present twice in input! " << (*detDigis).detId();
0186         }
0187       }
0188     }
0189     PixelDigi dig(digis.pdigi(i));
0190
0191 #ifdef GPU_DEBUG
0192     std::cout << i << ";" << digis.rawIdArr(i) << ";" << digis.clus(i) << ";" << digis.pdigi(i) << ";" << digis.adc(i)
0193               << ";" << dig.row() << ";" << dig.column() << std::endl;
0194 #endif
0195
0196     if (storeDigis_)
0197       (*detDigis).data.emplace_back(dig);
0198       // fill clusters
0199 #ifdef EDM_ML_DEBUG
0200     assert(digis.clus(i) >= 0);
0201     assert(digis.clus(i) < static_cast<int32_t>(TrackerTraits::maxNumClustersPerModules));
0202 #endif
0203     nclus = std::max(digis.clus(i), nclus);
0204     auto row = dig.row();
0205     auto col = dig.column();
0206     SiPixelCluster::PixelPos pix(row, col);
0207     aclusters[digis.clus(i)].add(pix, digis.adc(i));
0208   }
0209
0210   // fill final clusters
0211   if (detId > 0)
0212     fillClusters(detId);
0213
0214 #ifdef EDM_ML_DEBUG
0215   LogDebug("SiPixelDigisClustersFromSoA") << "filled " << totClustersFilled << " clusters";
0216 #endif
0217
0218   if (produceDigis_)
0219     iEvent.put(digiPutToken_, std::move(collection));
0220   iEvent.put(clusterPutToken_, std::move(outputClusters));
0221 }
0222
0223 using SiPixelDigisClustersFromSoAPhase1 = SiPixelDigisClustersFromSoAT<pixelTopology::Phase1>;
0224 DEFINE_FWK_MODULE(SiPixelDigisClustersFromSoAPhase1);
0225 using SiPixelDigisClustersFromSoAPhase2 = SiPixelDigisClustersFromSoAT<pixelTopology::Phase2>;
0226 DEFINE_FWK_MODULE(SiPixelDigisClustersFromSoAPhase2);
0227 using SiPixelDigisClustersFromSoAHIonPhase1 = SiPixelDigisClustersFromSoAT<pixelTopology::HIonPhase1>;
0228 DEFINE_FWK_MODULE(SiPixelDigisClustersFromSoAHIonPhase1);