File indexing completed on 2024-04-19 02:16:50
0001 #include <iostream>
0002 #include <vector>
0003 #include <memory>
0004
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "DataFormats/Common/interface/Handle.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 #include "FWCore/Utilities/interface/Exception.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0013
0014 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0015 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateIsolation.h"
0016
0017 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0018 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
0019
0020 #include "DataFormats/Math/interface/deltaR.h"
0021
0022 #include "RecoEgamma/EgammaIsolationAlgos/interface/EcalPFClusterIsolation.h"
0023 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0024
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026 #include "FWCore/Framework/interface/stream/EDProducer.h"
0027
0028 #include "FWCore/Framework/interface/Event.h"
0029 #include "FWCore/Framework/interface/MakerMacros.h"
0030
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032
0033 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0034 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateIsolation.h"
0035
0036 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0037 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateIsolation.h"
0038
0039 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0040 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
0041
0042 template <typename T1>
0043 class HLTEcalPFClusterIsolationProducer : public edm::stream::EDProducer<> {
0044 typedef std::vector<T1> T1Collection;
0045 typedef edm::Ref<T1Collection> T1Ref;
0046 typedef edm::AssociationMap<edm::OneToValue<std::vector<T1>, float>> T1IsolationMap;
0047
0048 public:
0049 explicit HLTEcalPFClusterIsolationProducer(const edm::ParameterSet&);
0050 ~HLTEcalPFClusterIsolationProducer() override;
0051
0052 void produce(edm::Event&, const edm::EventSetup&) override;
0053 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0054
0055 private:
0056 bool computedRVeto(T1Ref candRef, reco::PFClusterRef pfclu);
0057
0058 edm::EDGetTokenT<T1Collection> recoCandidateProducer_;
0059 const edm::EDGetTokenT<reco::PFClusterCollection> pfClusterProducer_;
0060 const edm::EDGetTokenT<double> rhoProducer_;
0061
0062 const double drMax_;
0063 const double drVetoBarrel_;
0064 const double drVetoEndcap_;
0065 const double etaStripBarrel_;
0066 const double etaStripEndcap_;
0067 const double energyBarrel_;
0068 const double energyEndcap_;
0069
0070 const bool doRhoCorrection_;
0071 const double rhoMax_;
0072 const double rhoScale_;
0073 const std::vector<double> effectiveAreas_;
0074 const std::vector<double> absEtaLowEdges_;
0075 };
0076
0077 template <typename T1>
0078 HLTEcalPFClusterIsolationProducer<T1>::HLTEcalPFClusterIsolationProducer(const edm::ParameterSet& config)
0079 : pfClusterProducer_(consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducer"))),
0080 rhoProducer_(consumes<double>(config.getParameter<edm::InputTag>("rhoProducer"))),
0081 drMax_(config.getParameter<double>("drMax")),
0082 drVetoBarrel_(config.getParameter<double>("drVetoBarrel")),
0083 drVetoEndcap_(config.getParameter<double>("drVetoEndcap")),
0084 etaStripBarrel_(config.getParameter<double>("etaStripBarrel")),
0085 etaStripEndcap_(config.getParameter<double>("etaStripEndcap")),
0086 energyBarrel_(config.getParameter<double>("energyBarrel")),
0087 energyEndcap_(config.getParameter<double>("energyEndcap")),
0088 doRhoCorrection_(config.getParameter<bool>("doRhoCorrection")),
0089 rhoMax_(config.getParameter<double>("rhoMax")),
0090 rhoScale_(config.getParameter<double>("rhoScale")),
0091 effectiveAreas_(config.getParameter<std::vector<double>>("effectiveAreas")),
0092 absEtaLowEdges_(config.getParameter<std::vector<double>>("absEtaLowEdges")) {
0093 if (doRhoCorrection_) {
0094 if (absEtaLowEdges_.size() != effectiveAreas_.size())
0095 throw cms::Exception("IncompatibleVects") << "absEtaLowEdges and effectiveAreas should be of the same size. \n";
0096
0097 if (absEtaLowEdges_.at(0) != 0.0)
0098 throw cms::Exception("IncompleteCoverage") << "absEtaLowEdges should start from 0. \n";
0099
0100 for (unsigned int aIt = 0; aIt < absEtaLowEdges_.size() - 1; aIt++) {
0101 if (!(absEtaLowEdges_.at(aIt) < absEtaLowEdges_.at(aIt + 1)))
0102 throw cms::Exception("ImproperBinning") << "absEtaLowEdges entries should be in increasing order. \n";
0103 }
0104 }
0105
0106 std::string recoCandidateProducerName = "recoCandidateProducer";
0107 if ((typeid(HLTEcalPFClusterIsolationProducer<T1>) ==
0108 typeid(HLTEcalPFClusterIsolationProducer<reco::RecoEcalCandidate>)))
0109 recoCandidateProducerName = "recoEcalCandidateProducer";
0110
0111 recoCandidateProducer_ = consumes<T1Collection>(config.getParameter<edm::InputTag>(recoCandidateProducerName));
0112 produces<T1IsolationMap>();
0113 }
0114
0115 template <typename T1>
0116 HLTEcalPFClusterIsolationProducer<T1>::~HLTEcalPFClusterIsolationProducer() {}
0117
0118 template <typename T1>
0119 void HLTEcalPFClusterIsolationProducer<T1>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0120 std::string recoCandidateProducerName = "recoCandidateProducer";
0121 if ((typeid(HLTEcalPFClusterIsolationProducer<T1>) ==
0122 typeid(HLTEcalPFClusterIsolationProducer<reco::RecoEcalCandidate>)))
0123 recoCandidateProducerName = "recoEcalCandidateProducer";
0124
0125 edm::ParameterSetDescription desc;
0126 desc.add<edm::InputTag>(recoCandidateProducerName, edm::InputTag("hltL1SeededRecoEcalCandidatePF"));
0127 desc.add<edm::InputTag>("pfClusterProducer", edm::InputTag("hltParticleFlowClusterECAL"));
0128 desc.add<edm::InputTag>("rhoProducer", edm::InputTag("fixedGridRhoFastjetAllCalo"));
0129 desc.add<bool>("doRhoCorrection", false);
0130 desc.add<double>("rhoMax", 9.9999999E7);
0131 desc.add<double>("rhoScale", 1.0);
0132 desc.add<double>("drMax", 0.3);
0133 desc.add<double>("drVetoBarrel", 0.0);
0134 desc.add<double>("drVetoEndcap", 0.0);
0135 desc.add<double>("etaStripBarrel", 0.0);
0136 desc.add<double>("etaStripEndcap", 0.0);
0137 desc.add<double>("energyBarrel", 0.0);
0138 desc.add<double>("energyEndcap", 0.0);
0139 desc.add<std::vector<double>>("effectiveAreas", {0.29, 0.21});
0140 desc.add<std::vector<double>>("absEtaLowEdges", {0.0, 1.479});
0141 descriptions.add(defaultModuleLabel<HLTEcalPFClusterIsolationProducer<T1>>(), desc);
0142 }
0143
0144 template <typename T1>
0145 void HLTEcalPFClusterIsolationProducer<T1>::produce(edm::Event& iEvent, const edm::EventSetup&) {
0146 edm::Handle<double> rhoHandle;
0147 double rho = 0.0;
0148 if (doRhoCorrection_) {
0149 iEvent.getByToken(rhoProducer_, rhoHandle);
0150 rho = *(rhoHandle.product());
0151 }
0152
0153 if (rho > rhoMax_)
0154 rho = rhoMax_;
0155
0156 rho = rho * rhoScale_;
0157
0158 edm::Handle<T1Collection> recoCandHandle;
0159 edm::Handle<reco::PFClusterCollection> clusterHandle;
0160
0161 iEvent.getByToken(recoCandidateProducer_, recoCandHandle);
0162 iEvent.getByToken(pfClusterProducer_, clusterHandle);
0163
0164 EcalPFClusterIsolation<T1> isoAlgo(
0165 drMax_, drVetoBarrel_, drVetoEndcap_, etaStripBarrel_, etaStripEndcap_, energyBarrel_, energyEndcap_);
0166 T1IsolationMap recoCandMap(recoCandHandle);
0167
0168 for (unsigned int iReco = 0; iReco < recoCandHandle->size(); iReco++) {
0169 T1Ref candRef(recoCandHandle, iReco);
0170
0171 float sum = isoAlgo.getSum(candRef, clusterHandle);
0172
0173 if (doRhoCorrection_) {
0174 int iEA = -1;
0175 auto cEta = std::abs(candRef->eta());
0176 for (int bIt = absEtaLowEdges_.size() - 1; bIt > -1; bIt--) {
0177 if (cEta >= absEtaLowEdges_[bIt]) {
0178 iEA = bIt;
0179 break;
0180 }
0181 }
0182 sum = sum - rho * effectiveAreas_[iEA];
0183 }
0184
0185 recoCandMap.insert(candRef, sum);
0186 }
0187
0188 iEvent.put(std::make_unique<T1IsolationMap>(recoCandMap));
0189 }
0190
0191 typedef HLTEcalPFClusterIsolationProducer<reco::RecoEcalCandidate> EgammaHLTEcalPFClusterIsolationProducer;
0192 typedef HLTEcalPFClusterIsolationProducer<reco::RecoChargedCandidate> MuonHLTEcalPFClusterIsolationProducer;
0193
0194 DEFINE_FWK_MODULE(EgammaHLTEcalPFClusterIsolationProducer);
0195 DEFINE_FWK_MODULE(MuonHLTEcalPFClusterIsolationProducer);