File indexing completed on 2021-07-30 02:33:25
0001 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0002 #include "DataFormats/ParticleFlowReco/interface/PFRecHitFraction.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "RecoParticleFlow/PFClusterProducer/interface/InitialClusteringStepBase.h"
0006 #include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
0007 #include "SimDataFormats/CaloAnalysis/interface/SimClusterFwd.h"
0008
0009 class GenericSimClusterMapper : public InitialClusteringStepBase {
0010 typedef GenericSimClusterMapper B2DGT;
0011
0012 public:
0013 GenericSimClusterMapper(const edm::ParameterSet& conf, edm::ConsumesCollector& cc)
0014 : InitialClusteringStepBase(conf, cc) {
0015 _simClusterToken = cc.consumes<SimClusterCollection>(conf.getParameter<edm::InputTag>("simClusterSrc"));
0016 }
0017 ~GenericSimClusterMapper() override = default;
0018 GenericSimClusterMapper(const B2DGT&) = delete;
0019 B2DGT& operator=(const B2DGT&) = delete;
0020
0021 void updateEvent(const edm::Event&) final;
0022
0023 void buildClusters(const edm::Handle<reco::PFRecHitCollection>&,
0024 const std::vector<bool>&,
0025 const std::vector<bool>&,
0026 reco::PFClusterCollection&) override;
0027
0028 private:
0029 edm::EDGetTokenT<SimClusterCollection> _simClusterToken;
0030 edm::Handle<SimClusterCollection> _simClusterH;
0031 };
0032
0033 DEFINE_EDM_PLUGIN(InitialClusteringStepFactory, GenericSimClusterMapper, "GenericSimClusterMapper");
0034
0035 #ifdef PFLOW_DEBUG
0036 #define LOGVERB(x) edm::LogVerbatim(x)
0037 #define LOGWARN(x) edm::LogWarning(x)
0038 #define LOGERR(x) edm::LogError(x)
0039 #define LOGDRESSED(x) edm::LogInfo(x)
0040 #else
0041 #define LOGVERB(x) LogTrace(x)
0042 #define LOGWARN(x) edm::LogWarning(x)
0043 #define LOGERR(x) edm::LogError(x)
0044 #define LOGDRESSED(x) LogDebug(x)
0045 #endif
0046
0047 void GenericSimClusterMapper::updateEvent(const edm::Event& ev) { ev.getByToken(_simClusterToken, _simClusterH); }
0048
0049 void GenericSimClusterMapper::buildClusters(const edm::Handle<reco::PFRecHitCollection>& input,
0050 const std::vector<bool>& rechitMask,
0051 const std::vector<bool>& seedable,
0052 reco::PFClusterCollection& output) {
0053 const SimClusterCollection& simClusters = *_simClusterH;
0054 auto const& hits = *input;
0055
0056
0057 std::unordered_map<uint32_t, size_t> detIdToIndex(hits.size());
0058 for (uint32_t i = 0; i < hits.size(); ++i) {
0059 detIdToIndex[hits[i].detId()] = i;
0060 auto ref = makeRefhit(input, i);
0061 }
0062
0063 for (const auto& sc : simClusters) {
0064 output.emplace_back();
0065 reco::PFCluster& back = output.back();
0066 edm::Ref<std::vector<reco::PFRecHit> > seed;
0067 double energy = 0.0, highest_energy = 0.0;
0068 auto hitsAndFractions = sc.hits_and_fractions();
0069 for (const auto& hAndF : hitsAndFractions) {
0070 auto itr = detIdToIndex.find(hAndF.first);
0071 if (itr == detIdToIndex.end())
0072 continue;
0073 auto ref = makeRefhit(input, itr->second);
0074 const double hit_energy = hAndF.second * ref->energy();
0075 energy += hit_energy;
0076 back.addRecHitFraction(reco::PFRecHitFraction(ref, hAndF.second));
0077 if (hit_energy > highest_energy || highest_energy == 0.0) {
0078 highest_energy = hit_energy;
0079 seed = ref;
0080 }
0081 }
0082 if (!back.hitsAndFractions().empty()) {
0083 back.setSeed(seed->detId());
0084 back.setEnergy(energy);
0085 back.setCorrectedEnergy(energy);
0086 } else {
0087 back.setSeed(-1);
0088 back.setEnergy(0.f);
0089 }
0090 }
0091 }