Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 #include "FWCore/Framework/interface/ConsumesCollector.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005 #include "FWCore/Framework/interface/stream/EDProducer.h"
0006 #include "RecoParticleFlow/PFClusterProducer/interface/InitialClusteringStepBase.h"
0007 #include "RecoParticleFlow/PFClusterProducer/interface/PFCPositionCalculatorBase.h"
0008 #include "RecoParticleFlow/PFClusterProducer/interface/PFClusterBuilderBase.h"
0009 #include "RecoParticleFlow/PFClusterProducer/interface/PFClusterEnergyCorrectorBase.h"
0010 #include "RecoParticleFlow/PFClusterProducer/interface/RecHitTopologicalCleanerBase.h"
0011 #include "RecoParticleFlow/PFClusterProducer/interface/SeedFinderBase.h"
0012 #include "CondFormats/DataRecord/interface/HcalPFCutsRcd.h"
0013 #include "CondTools/Hcal/interface/HcalPFCutsHandler.h"
0014 
0015 #include <memory>
0016 
0017 class PFMultiDepthClusterProducer : public edm::stream::EDProducer<> {
0018   typedef InitialClusteringStepBase ICSB;
0019   typedef PFClusterBuilderBase PFCBB;
0020   typedef PFCPositionCalculatorBase PosCalc;
0021 
0022 public:
0023   PFMultiDepthClusterProducer(const edm::ParameterSet&);
0024   ~PFMultiDepthClusterProducer() override = default;
0025 
0026   void beginRun(const edm::Run&, const edm::EventSetup&) override;
0027   void produce(edm::Event&, const edm::EventSetup&) override;
0028 
0029 private:
0030   // inputs
0031   edm::EDGetTokenT<reco::PFClusterCollection> _clustersLabel;
0032   edm::ESGetToken<HcalPFCuts, HcalPFCutsRcd> hcalCutsToken_;
0033   // options
0034   // the actual algorithm
0035   std::unique_ptr<PFClusterBuilderBase> _pfClusterBuilder;
0036   std::unique_ptr<PFClusterEnergyCorrectorBase> _energyCorrector;
0037   bool cutsFromDB;
0038   HcalPFCuts const* paramPF = nullptr;
0039 };
0040 
0041 DEFINE_FWK_MODULE(PFMultiDepthClusterProducer);
0042 
0043 #ifdef PFLOW_DEBUG
0044 #define LOGVERB(x) edm::LogVerbatim(x)
0045 #define LOGWARN(x) edm::LogWarning(x)
0046 #define LOGERR(x) edm::LogError(x)
0047 #define LOGDRESSED(x) edm::LogInfo(x)
0048 #else
0049 #define LOGVERB(x) LogTrace(x)
0050 #define LOGWARN(x) edm::LogWarning(x)
0051 #define LOGERR(x) edm::LogError(x)
0052 #define LOGDRESSED(x) LogDebug(x)
0053 #endif
0054 
0055 PFMultiDepthClusterProducer::PFMultiDepthClusterProducer(const edm::ParameterSet& conf) {
0056   _clustersLabel = consumes<reco::PFClusterCollection>(conf.getParameter<edm::InputTag>("clustersSource"));
0057   cutsFromDB = conf.getParameter<bool>("usePFThresholdsFromDB");
0058   const edm::ParameterSet& pfcConf = conf.getParameterSet("pfClusterBuilder");
0059 
0060   edm::ConsumesCollector&& cc = consumesCollector();
0061 
0062   if (cutsFromDB) {
0063     hcalCutsToken_ = esConsumes<HcalPFCuts, HcalPFCutsRcd, edm::Transition::BeginRun>(edm::ESInputTag("", "withTopo"));
0064   }
0065 
0066   if (!pfcConf.empty()) {
0067     const std::string& pfcName = pfcConf.getParameter<std::string>("algoName");
0068     _pfClusterBuilder = PFClusterBuilderFactory::get()->create(pfcName, pfcConf, cc);
0069   }
0070   // see if new need to apply corrections, setup if there.
0071   const edm::ParameterSet& cConf = conf.getParameterSet("energyCorrector");
0072   if (!cConf.empty()) {
0073     const std::string& cName = cConf.getParameter<std::string>("algoName");
0074     _energyCorrector = PFClusterEnergyCorrectorFactory::get()->create(cName, cConf);
0075   }
0076 
0077   produces<reco::PFClusterCollection>();
0078 }
0079 
0080 void PFMultiDepthClusterProducer::beginRun(const edm::Run& run, const edm::EventSetup& es) {
0081   if (cutsFromDB) {
0082     paramPF = &es.getData(hcalCutsToken_);
0083   }
0084   _pfClusterBuilder->update(es);
0085 }
0086 
0087 void PFMultiDepthClusterProducer::produce(edm::Event& e, const edm::EventSetup& es) {
0088   _pfClusterBuilder->reset();
0089 
0090   edm::Handle<reco::PFClusterCollection> inputClusters;
0091   e.getByToken(_clustersLabel, inputClusters);
0092 
0093   std::vector<bool> seedable;
0094 
0095   auto pfClusters = std::make_unique<reco::PFClusterCollection>();
0096   _pfClusterBuilder->buildClusters(*inputClusters, seedable, *pfClusters, paramPF);
0097   LOGVERB("PFMultiDepthClusterProducer::produce()") << *_pfClusterBuilder;
0098 
0099   if (_energyCorrector) {
0100     _energyCorrector->correctEnergies(*pfClusters);
0101   }
0102   e.put(std::move(pfClusters));
0103 }