Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-30 02:33:25

0001 #include "DataFormats/ParticleFlowReco/interface/PFRecHitFraction.h"
0002 #include "RecoParticleFlow/PFClusterProducer/interface/InitialClusteringStepBase.h"
0003 
0004 class Basic2DClusterForEachSeed : public InitialClusteringStepBase {
0005 public:
0006   Basic2DClusterForEachSeed(const edm::ParameterSet& conf, edm::ConsumesCollector& cc)
0007       : InitialClusteringStepBase(conf, cc) {}
0008   ~Basic2DClusterForEachSeed() override = default;
0009 
0010   void buildClusters(const edm::Handle<reco::PFRecHitCollection>&,
0011                      const std::vector<bool>&,
0012                      const std::vector<bool>&,
0013                      reco::PFClusterCollection&) override;
0014 };
0015 
0016 DEFINE_EDM_PLUGIN(InitialClusteringStepFactory, Basic2DClusterForEachSeed, "Basic2DClusterForEachSeed");
0017 
0018 void Basic2DClusterForEachSeed::buildClusters(const edm::Handle<reco::PFRecHitCollection>& input,
0019                                               const std::vector<bool>& rechitMask,
0020                                               const std::vector<bool>& seedable,
0021                                               reco::PFClusterCollection& output) {
0022   auto const& hits = *input;
0023 
0024   // loop over seeds and make clusters
0025   reco::PFCluster cluster;
0026   for (unsigned int hit = 0; hit < hits.size(); ++hit) {
0027     if (!rechitMask[hit] || !seedable[hit])
0028       continue;  // if not seed, ignore.
0029     cluster.reset();
0030 
0031     // seed
0032     auto refhit = makeRefhit(input, hit);
0033     auto rhf = reco::PFRecHitFraction(refhit, 1.0);  // entire rechit energy should go to a cluster
0034 
0035     // add the hit to the cluster
0036     cluster.addRecHitFraction(rhf);
0037 
0038     // extract
0039     const auto rh_energy = refhit->energy();
0040 
0041     // fill cluster information
0042     cluster.setSeed(refhit->detId());
0043     cluster.setEnergy(rh_energy);
0044     cluster.setTime(refhit->time());
0045     cluster.setLayer(refhit->layer());
0046     cluster.setPosition(math::XYZPoint(refhit->position().x(), refhit->position().y(), refhit->position().z()));
0047     cluster.calculatePositionREP();
0048     cluster.setDepth(refhit->depth());
0049 
0050     output.push_back(cluster);
0051 
0052   }  // looping over seeds ends
0053 }