CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
HGCalBackendLayer2Processor3DClustering.cc
Go to the documentation of this file.
2 
12 
13 #include <utility>
14 
16 public:
18  std::string typeMulticluster(conf.getParameterSet("C3d_parameters").getParameter<std::string>("type_multicluster"));
19  if (typeMulticluster == "dRC3d") {
21  multiclustering_ = std::make_unique<HGCalMulticlusteringImpl>(conf.getParameterSet("C3d_parameters"));
22  } else if (typeMulticluster == "DBSCANC3d") {
24  multiclustering_ = std::make_unique<HGCalMulticlusteringImpl>(conf.getParameterSet("C3d_parameters"));
25  } else if (typeMulticluster == "Histo") {
27  multiclusteringHistoSeeding_ = std::make_unique<HGCalHistoSeedingImpl>(
28  conf.getParameterSet("C3d_parameters").getParameterSet("histoMax_C3d_seeding_parameters"));
29  multiclusteringHistoClustering_ = std::make_unique<HGCalHistoClusteringImpl>(
30  conf.getParameterSet("C3d_parameters").getParameterSet("histoMax_C3d_clustering_parameters"));
31  } else {
32  throw cms::Exception("HGCTriggerParameterError") << "Unknown Multiclustering type '" << typeMulticluster << "'";
33  }
34 
35  for (const auto& interpretationPset : conf.getParameter<std::vector<edm::ParameterSet>>("energy_interpretations")) {
36  std::unique_ptr<HGCalTriggerClusterInterpreterBase> interpreter{
37  HGCalTriggerClusterInterpreterFactory::get()->create(interpretationPset.getParameter<std::string>("type"))};
38  interpreter->initialize(interpretationPset);
39  energy_interpreters_.push_back(std::move(interpreter));
40  }
41  }
42 
44  std::pair<l1t::HGCalMulticlusterBxCollection, l1t::HGCalClusterBxCollection>& be_output) override {
45  if (multiclustering_)
46  multiclustering_->setGeometry(geometry());
48  multiclusteringHistoSeeding_->setGeometry(geometry());
51 
52  auto& collCluster3D = be_output.first;
53  auto& rejectedClusters = be_output.second;
54 
55  /* create a persistent vector of pointers to the trigger-cells */
56  std::vector<edm::Ptr<l1t::HGCalCluster>> clustersPtrs;
57  for (unsigned i = 0; i < collHandle->size(); ++i) {
58  edm::Ptr<l1t::HGCalCluster> ptr(collHandle, i);
59  clustersPtrs.push_back(ptr);
60  }
61 
62  /* create a vector of seed positions and their energy*/
63  std::vector<std::pair<GlobalPoint, double>> seedPositionsEnergy;
64 
65  /* call to multiclustering and compute shower shape*/
66  switch (multiclusteringAlgoType_) {
67  case dRC3d:
68  multiclustering_->clusterizeDR(clustersPtrs, collCluster3D, *geometry());
69  break;
70  case DBSCANC3d:
71  multiclustering_->clusterizeDBSCAN(clustersPtrs, collCluster3D, *geometry());
72  break;
73  case HistoC3d:
74  multiclusteringHistoSeeding_->findHistoSeeds(clustersPtrs, seedPositionsEnergy);
75  multiclusteringHistoClustering_->clusterizeHisto(
76  clustersPtrs, seedPositionsEnergy, *geometry(), collCluster3D, rejectedClusters);
77  break;
78  default:
79  // Should not happen, clustering type checked in constructor
80  break;
81  }
82 
83  // Call all the energy interpretation modules on the cluster collection
84  for (const auto& interpreter : energy_interpreters_) {
85  interpreter->setGeometry(geometry());
86  interpreter->interpret(collCluster3D);
87  }
88  }
89 
90 private:
92 
93  /* algorithms instances */
94  std::unique_ptr<HGCalMulticlusteringImpl> multiclustering_;
95  std::unique_ptr<HGCalHistoSeedingImpl> multiclusteringHistoSeeding_;
96  std::unique_ptr<HGCalHistoClusteringImpl> multiclusteringHistoClustering_;
97 
98  /* algorithm type */
100 
101  std::vector<std::unique_ptr<HGCalTriggerClusterInterpreterBase>> energy_interpreters_;
102 };
103 
106  "HGCalBackendLayer2Processor3DClustering");
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
ParameterSet const & getParameterSet(std::string const &) const
std::unique_ptr< HGCalMulticlusteringImpl > multiclustering_
unsigned size(int bx) const
std::unique_ptr< HGCalHistoSeedingImpl > multiclusteringHistoSeeding_
void run(const edm::Handle< l1t::HGCalClusterBxCollection > &collHandle, std::pair< l1t::HGCalMulticlusterBxCollection, l1t::HGCalClusterBxCollection > &be_output) override
const HGCalTriggerGeometryBase * geometry() const
std::unique_ptr< HGCalHistoClusteringImpl > multiclusteringHistoClustering_
std::vector< std::unique_ptr< HGCalTriggerClusterInterpreterBase > > energy_interpreters_
#define DEFINE_EDM_PLUGIN(factory, type, name)
#define get
def move(src, dest)
Definition: eostools.py:511