Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:27:34

0001 #include <array>
0002 #include <memory>
0003 #include <vector>
0004 
0005 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "HeterogeneousCore/AlpakaCore/interface/alpaka/ESProducer.h"
0008 #include "RecoParticleFlow/PFRecHitProducer/interface/PFRecHitParamsHostCollection.h"
0009 #include "RecoParticleFlow/PFRecHitProducer/interface/PFRecHitParamsRecord.h"
0010 #include "CalorimeterDefinitions.h"
0011 
0012 namespace ALPAKA_ACCELERATOR_NAMESPACE {
0013   using namespace particleFlowRecHitProducer;
0014 
0015   class PFRecHitHCALParamsESProducer : public ESProducer {
0016   public:
0017     PFRecHitHCALParamsESProducer(edm::ParameterSet const& iConfig)
0018         : ESProducer(iConfig),
0019           energyThresholdsHB_(iConfig.getParameter<std::array<double, HCAL::kMaxDepthHB>>("energyThresholdsHB")),
0020           energyThresholdsHE_(iConfig.getParameter<std::array<double, HCAL::kMaxDepthHE>>("energyThresholdsHE")) {
0021       setWhatProduced(this);
0022     }
0023 
0024     static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0025       edm::ParameterSetDescription desc;
0026       desc.add<std::vector<double>>("energyThresholdsHB", {0.1, 0.2, 0.3, 0.3});
0027       desc.add<std::vector<double>>("energyThresholdsHE", {0.1, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2});
0028       descriptions.addWithDefaultLabel(desc);
0029     }
0030 
0031     std::unique_ptr<reco::PFRecHitHCALParamsHostCollection> produce(PFRecHitHCALParamsRecord const& iRecord) {
0032       auto product = std::make_unique<reco::PFRecHitHCALParamsHostCollection>(HCAL::kMaxDepthHB + HCAL::kMaxDepthHE,
0033                                                                               cms::alpakatools::host());
0034       for (uint32_t idx = 0; idx < HCAL::kMaxDepthHB; ++idx) {
0035         product->view().energyThresholds()[idx] = energyThresholdsHB_[idx];
0036       }
0037       for (uint32_t idx = 0; idx < HCAL::kMaxDepthHE; ++idx) {
0038         product->view().energyThresholds()[idx + HCAL::kMaxDepthHB] = energyThresholdsHE_[idx];
0039       }
0040       return product;
0041     }
0042 
0043   private:
0044     std::array<double, HCAL::kMaxDepthHB> energyThresholdsHB_;
0045     std::array<double, HCAL::kMaxDepthHE> energyThresholdsHE_;
0046   };
0047 
0048 }  // namespace ALPAKA_ACCELERATOR_NAMESPACE
0049 
0050 #include "HeterogeneousCore/AlpakaCore/interface/alpaka/ModuleFactory.h"
0051 DEFINE_FWK_EVENTSETUP_ALPAKA_MODULE(PFRecHitHCALParamsESProducer);