File indexing completed on 2025-03-06 03:07:34
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 "RecoEgamma/EgammaIsolationAlgos/interface/HcalPFClusterIsolation.h"
0015 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0016
0017 #include "FWCore/Framework/interface/Frameworkfwd.h"
0018 #include "FWCore/Framework/interface/global/EDProducer.h"
0019
0020 #include "FWCore/Framework/interface/Event.h"
0021 #include "FWCore/Framework/interface/MakerMacros.h"
0022
0023 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0024
0025 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0026 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateIsolation.h"
0027
0028 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0029 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateIsolation.h"
0030
0031 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0032 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
0033
0034 template <typename T1>
0035 class HLTHcalPFClusterIsolationProducer : public edm::global::EDProducer<> {
0036 typedef std::vector<T1> T1Collection;
0037 typedef edm::Ref<T1Collection> T1Ref;
0038 typedef edm::AssociationMap<edm::OneToValue<std::vector<T1>, float>> T1IsolationMap;
0039
0040 public:
0041 explicit HLTHcalPFClusterIsolationProducer(const edm::ParameterSet&);
0042 ~HLTHcalPFClusterIsolationProducer() override;
0043
0044 void produce(edm::StreamID sid, edm::Event&, const edm::EventSetup&) const override;
0045 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0046
0047 private:
0048 edm::EDGetTokenT<T1Collection> recoCandidateProducer_;
0049 const edm::EDGetTokenT<reco::PFClusterCollection> pfClusterProducerHCAL_;
0050 const edm::EDGetTokenT<double> rhoProducer_;
0051 const edm::EDGetTokenT<reco::PFClusterCollection> pfClusterProducerHFEM_;
0052 const edm::EDGetTokenT<reco::PFClusterCollection> pfClusterProducerHFHAD_;
0053
0054 const bool useHF_;
0055
0056 const double drMax_;
0057 const double drVetoBarrel_;
0058 const double drVetoEndcap_;
0059 const double etaStripBarrel_;
0060 const double etaStripEndcap_;
0061 const double energyBarrel_;
0062 const double energyEndcap_;
0063 const bool useEt_;
0064
0065 const bool doRhoCorrection_;
0066 const double rhoMax_;
0067 const double rhoScale_;
0068 const std::vector<double> effectiveAreas_;
0069 const std::vector<double> absEtaLowEdges_;
0070
0071 const bool doEffAreaCorrection_;
0072 const std::vector<double> effectiveAreasCorr_;
0073 const std::vector<double> effectiveAreasThres_;
0074 };
0075
0076 template <typename T1>
0077 HLTHcalPFClusterIsolationProducer<T1>::HLTHcalPFClusterIsolationProducer(const edm::ParameterSet& config)
0078 : pfClusterProducerHCAL_(
0079 consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducerHCAL"))),
0080 rhoProducer_(consumes<double>(config.getParameter<edm::InputTag>("rhoProducer"))),
0081 pfClusterProducerHFEM_(
0082 consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducerHFEM"))),
0083 pfClusterProducerHFHAD_(
0084 consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducerHFHAD"))),
0085 useHF_(config.getParameter<bool>("useHF")),
0086 drMax_(config.getParameter<double>("drMax")),
0087 drVetoBarrel_(config.getParameter<double>("drVetoBarrel")),
0088 drVetoEndcap_(config.getParameter<double>("drVetoEndcap")),
0089 etaStripBarrel_(config.getParameter<double>("etaStripBarrel")),
0090 etaStripEndcap_(config.getParameter<double>("etaStripEndcap")),
0091 energyBarrel_(config.getParameter<double>("energyBarrel")),
0092 energyEndcap_(config.getParameter<double>("energyEndcap")),
0093 useEt_(config.getParameter<bool>("useEt")),
0094 doRhoCorrection_(config.getParameter<bool>("doRhoCorrection")),
0095 rhoMax_(config.getParameter<double>("rhoMax")),
0096 rhoScale_(config.getParameter<double>("rhoScale")),
0097 effectiveAreas_(config.getParameter<std::vector<double>>("effectiveAreas")),
0098 absEtaLowEdges_(config.getParameter<std::vector<double>>("absEtaLowEdges")),
0099 doEffAreaCorrection_(config.getParameter<bool>("doEffAreaCorrection")),
0100 effectiveAreasCorr_(config.getParameter<std::vector<double>>("effectiveAreasCorr")),
0101 effectiveAreasThres_(config.getParameter<std::vector<double>>("effectiveAreasThres")) {
0102 if (doRhoCorrection_) {
0103 if (absEtaLowEdges_.size() != effectiveAreas_.size())
0104 throw cms::Exception("IncompatibleVects") << "absEtaLowEdges and effectiveAreas should be of the same size. \n";
0105
0106 if (absEtaLowEdges_.at(0) != 0.0)
0107 throw cms::Exception("IncompleteCoverage") << "absEtaLowEdges should start from 0. \n";
0108
0109 for (unsigned int aIt = 0; aIt < absEtaLowEdges_.size() - 1; aIt++) {
0110 if (!(absEtaLowEdges_.at(aIt) < absEtaLowEdges_.at(aIt + 1)))
0111 throw cms::Exception("ImproperBinning") << "absEtaLowEdges entries should be in increasing order. \n";
0112 }
0113 }
0114
0115 std::string recoCandidateProducerName = "recoCandidateProducer";
0116 if ((typeid(HLTHcalPFClusterIsolationProducer<T1>) ==
0117 typeid(HLTHcalPFClusterIsolationProducer<reco::RecoEcalCandidate>)))
0118 recoCandidateProducerName = "recoEcalCandidateProducer";
0119 recoCandidateProducer_ = consumes<T1Collection>(config.getParameter<edm::InputTag>(recoCandidateProducerName));
0120
0121 produces<T1IsolationMap>();
0122 }
0123
0124 template <typename T1>
0125 HLTHcalPFClusterIsolationProducer<T1>::~HLTHcalPFClusterIsolationProducer() {}
0126
0127 template <typename T1>
0128 void HLTHcalPFClusterIsolationProducer<T1>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0129 std::string recoCandidateProducerName = "recoCandidateProducer";
0130 if ((typeid(HLTHcalPFClusterIsolationProducer<T1>) ==
0131 typeid(HLTHcalPFClusterIsolationProducer<reco::RecoEcalCandidate>)))
0132 recoCandidateProducerName = "recoEcalCandidateProducer";
0133
0134 edm::ParameterSetDescription desc;
0135 desc.add<edm::InputTag>(recoCandidateProducerName, edm::InputTag("hltL1SeededRecoEcalCandidatePF"));
0136 desc.add<edm::InputTag>("pfClusterProducerHCAL", edm::InputTag("hltParticleFlowClusterHCAL"));
0137 desc.ifValue(
0138 edm::ParameterDescription<bool>("useHF", false, true),
0139 true >> (edm::ParameterDescription<edm::InputTag>(
0140 "pfClusterProducerHFEM", edm::InputTag("hltParticleFlowClusterHFEM"), true) and
0141 edm::ParameterDescription<edm::InputTag>(
0142 "pfClusterProducerHFHAD", edm::InputTag("hltParticleFlowClusterHFHAD"), true)) or
0143 false >> (edm::ParameterDescription<edm::InputTag>("pfClusterProducerHFEM", edm::InputTag(""), true) and
0144 edm::ParameterDescription<edm::InputTag>("pfClusterProducerHFHAD", edm::InputTag(""), true)));
0145 desc.add<edm::InputTag>("rhoProducer", edm::InputTag("fixedGridRhoFastjetAllCalo"));
0146 desc.add<bool>("doRhoCorrection", false);
0147 desc.add<double>("rhoMax", 9.9999999E7);
0148 desc.add<double>("rhoScale", 1.0);
0149 desc.add<double>("drMax", 0.3);
0150 desc.add<double>("drVetoBarrel", 0.0);
0151 desc.add<double>("drVetoEndcap", 0.0);
0152 desc.add<double>("etaStripBarrel", 0.0);
0153 desc.add<double>("etaStripEndcap", 0.0);
0154 desc.add<double>("energyBarrel", 0.0);
0155 desc.add<double>("energyEndcap", 0.0);
0156 desc.add<bool>("useEt", true);
0157 desc.add<std::vector<double>>("effectiveAreas", {0.2, 0.25});
0158 desc.add<std::vector<double>>("absEtaLowEdges", {0.0, 1.479});
0159 desc.add<bool>("doEffAreaCorrection", false);
0160 desc.add<std::vector<double>>("effectiveAreasCorr", {0.0, 0.0});
0161 desc.add<std::vector<double>>("effectiveAreasThres", {0.0, 0.0});
0162 descriptions.add(defaultModuleLabel<HLTHcalPFClusterIsolationProducer<T1>>(), desc);
0163 }
0164
0165 template <typename T1>
0166 void HLTHcalPFClusterIsolationProducer<T1>::produce(edm::StreamID sid,
0167 edm::Event& iEvent,
0168 const edm::EventSetup&) const {
0169 edm::Handle<double> rhoHandle;
0170 double rho = 0.0;
0171 if (doRhoCorrection_) {
0172 iEvent.getByToken(rhoProducer_, rhoHandle);
0173 rho = *(rhoHandle.product());
0174 }
0175
0176 if (rho > rhoMax_)
0177 rho = rhoMax_;
0178
0179 rho = rho * rhoScale_;
0180
0181 edm::Handle<T1Collection> recoCandHandle;
0182
0183 std::vector<edm::Handle<reco::PFClusterCollection>> clusterHandles;
0184 edm::Handle<reco::PFClusterCollection> clusterHcalHandle;
0185 edm::Handle<reco::PFClusterCollection> clusterHfemHandle;
0186 edm::Handle<reco::PFClusterCollection> clusterHfhadHandle;
0187
0188 iEvent.getByToken(recoCandidateProducer_, recoCandHandle);
0189 iEvent.getByToken(pfClusterProducerHCAL_, clusterHcalHandle);
0190 clusterHandles.push_back(clusterHcalHandle);
0191
0192 if (useHF_) {
0193 iEvent.getByToken(pfClusterProducerHFEM_, clusterHfemHandle);
0194 clusterHandles.push_back(clusterHfemHandle);
0195 iEvent.getByToken(pfClusterProducerHFHAD_, clusterHfhadHandle);
0196 clusterHandles.push_back(clusterHfhadHandle);
0197 }
0198
0199 T1IsolationMap recoCandMap(recoCandHandle);
0200 HcalPFClusterIsolation<T1> isoAlgo(
0201 drMax_, drVetoBarrel_, drVetoEndcap_, etaStripBarrel_, etaStripEndcap_, energyBarrel_, energyEndcap_, useEt_);
0202
0203 for (unsigned int iReco = 0; iReco < recoCandHandle->size(); iReco++) {
0204 T1Ref candRef(recoCandHandle, iReco);
0205 float sum = isoAlgo.getSum(candRef, clusterHandles);
0206
0207 if (doRhoCorrection_) {
0208 int iEA = -1;
0209 auto cEta = std::abs(candRef->eta());
0210 for (int bIt = absEtaLowEdges_.size() - 1; bIt > -1; bIt--) {
0211 if (cEta >= absEtaLowEdges_[bIt]) {
0212 iEA = bIt;
0213 break;
0214 }
0215 }
0216 if (doEffAreaCorrection_) {
0217 float correction = (rho > effectiveAreasThres_[iEA]) ? (1 + effectiveAreasCorr_[iEA]) : 1;
0218 sum -= rho * correction * effectiveAreas_[iEA];
0219 } else {
0220 sum -= rho * effectiveAreas_[iEA];
0221 }
0222 }
0223
0224 recoCandMap.insert(candRef, sum);
0225 }
0226
0227 iEvent.put(std::make_unique<T1IsolationMap>(recoCandMap));
0228 }
0229
0230 typedef HLTHcalPFClusterIsolationProducer<reco::RecoEcalCandidate> EgammaHLTHcalPFClusterIsolationProducer;
0231 typedef HLTHcalPFClusterIsolationProducer<reco::RecoChargedCandidate> MuonHLTHcalPFClusterIsolationProducer;
0232
0233 DEFINE_FWK_MODULE(EgammaHLTHcalPFClusterIsolationProducer);
0234 DEFINE_FWK_MODULE(MuonHLTHcalPFClusterIsolationProducer);
