CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
HGCalSoALayerClustersProducer.cc
Go to the documentation of this file.
16 
23 
25 
27  public:
29  : getTokenDeviceRecHits_{consumes(config.getParameter<edm::InputTag>("hgcalRecHitsSoA"))},
30  getTokenDeviceClusters_{consumes(config.getParameter<edm::InputTag>("hgcalRecHitsLayerClustersSoA"))},
32  thresholdW0_(config.getParameter<double>("thresholdW0")),
33  positionDeltaRho2_(config.getParameter<double>("positionDeltaRho2")) {}
34 
35  ~HGCalSoALayerClustersProducer() override = default;
36 
37  void acquire(device::Event const& iEvent, device::EventSetup const& iSetup) override {
38  // Get LayerClusters almost-SoA on device: this has still the same
39  // cardinality as the RecHitsSoA, but has all the required information
40  // to assemble the clusters, i.e., it has the cluster index assigned to
41  // each rechit.
42  auto const& deviceInputClusters = iEvent.get(getTokenDeviceClusters_);
43  auto const inputClusters_v = deviceInputClusters.view();
44  //
45  // Allocate output SoA for the clusters, one entry for each cluster
46  auto device_numclusters = cms::alpakatools::make_device_view<const unsigned int>(
47  alpaka::getDev(iEvent.queue()), inputClusters_v.numberOfClustersScalar());
48  auto host_numclusters = cms::alpakatools::make_host_view<unsigned int>(num_clusters_);
49  alpaka::memcpy(iEvent.queue(), host_numclusters, device_numclusters);
50  }
51 
52  void produce(device::Event& iEvent, device::EventSetup const& iSetup) override {
53  // Get RecHitsSoA on the device
54  auto const& deviceInputRecHits = iEvent.get(getTokenDeviceRecHits_);
55  auto const inputRechits_v = deviceInputRecHits.view();
56 
57  // Get LayerClusters almost-SoA on device: this has still the same
58  // cardinality as the RecHitsSoA, but has all the required information
59  // to assemble the clusters, i.e., it has the cluster index assigned to
60  // each rechit.
61  auto const& deviceInputClusters = iEvent.get(getTokenDeviceClusters_);
62  auto const inputClusters_v = deviceInputClusters.view();
63 
65  auto output_v = output.view();
66  // Allocate workspace SoA cluster
68  auto output_workspace_v = outputWorkspace.view();
69 
70  algo_.run(iEvent.queue(),
74  inputRechits_v,
75  inputClusters_v,
76  output_v,
77  output_workspace_v);
79  }
80 
81  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
83  desc.add<edm::InputTag>("hgcalRecHitsLayerClustersSoA", edm::InputTag("TO BE DEFINED"));
84  desc.add<edm::InputTag>("hgcalRecHitsSoA", edm::InputTag("TO BE DEFINED"));
85  desc.add<double>("thresholdW0", 2.9);
86  desc.add<double>("positionDeltaRho2", 1.69);
87  descriptions.addWithDefaultLabel(desc);
88  }
89 
90  private:
95  unsigned int num_clusters_;
96  float thresholdW0_;
98  };
99 
100 } // namespace ALPAKA_ACCELERATOR_NAMESPACE
101 
103 DEFINE_FWK_ALPAKA_MODULE(HGCalSoALayerClustersProducer);
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
device::EDGetToken< HGCalSoARecHitsDeviceCollection > const getTokenDeviceRecHits_
device::EDPutToken< HGCalSoAClustersDeviceCollection > const deviceTokenSoAClusters_
Definition: config.py:1
int iEvent
Definition: GenABIO.cc:224
void produce(device::Event &iEvent, device::EventSetup const &iSetup) override
void run(Queue &queue, const unsigned int numer_of_clusters, float thresholdW0, float positionDeltaRho2, const HGCalSoARecHitsDeviceCollection::ConstView input_rechits_soa, const HGCalSoARecHitsExtraDeviceCollection::ConstView input_clusters_soa, HGCalSoAClustersDeviceCollection::View outputs, HGCalSoAClustersExtraDeviceCollection::View outputs_service) const
PortableCollection< HGCalSoAClustersExtra > HGCalSoAClustersExtraDeviceCollection
PortableCollection< HGCalSoAClusters > HGCalSoAClustersDeviceCollection
device::EDGetToken< HGCalSoARecHitsExtraDeviceCollection > const getTokenDeviceClusters_
void acquire(device::Event const &iEvent, device::EventSetup const &iSetup) override
auto produces(std::string instanceName) noexcept
declare what type of product will make and with which optional label
#define DEFINE_FWK_ALPAKA_MODULE(name)
Definition: MakerMacros.h:16
Definition: output.py:1
def move(src, dest)
Definition: eostools.py:511