File indexing completed on 2021-02-14 13:34:09
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
0072 template <typename T1>
0073 HLTHcalPFClusterIsolationProducer<T1>::HLTHcalPFClusterIsolationProducer(const edm::ParameterSet& config)
0074 : pfClusterProducerHCAL_(
0075 consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducerHCAL"))),
0076 rhoProducer_(consumes<double>(config.getParameter<edm::InputTag>("rhoProducer"))),
0077 pfClusterProducerHFEM_(
0078 consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducerHFEM"))),
0079 pfClusterProducerHFHAD_(
0080 consumes<reco::PFClusterCollection>(config.getParameter<edm::InputTag>("pfClusterProducerHFHAD"))),
0081 useHF_(config.getParameter<bool>("useHF")),
0082 drMax_(config.getParameter<double>("drMax")),
0083 drVetoBarrel_(config.getParameter<double>("drVetoBarrel")),
0084 drVetoEndcap_(config.getParameter<double>("drVetoEndcap")),
0085 etaStripBarrel_(config.getParameter<double>("etaStripBarrel")),
0086 etaStripEndcap_(config.getParameter<double>("etaStripEndcap")),
0087 energyBarrel_(config.getParameter<double>("energyBarrel")),
0088 energyEndcap_(config.getParameter<double>("energyEndcap")),
0089 useEt_(config.getParameter<bool>("useEt")),
0090 doRhoCorrection_(config.getParameter<bool>("doRhoCorrection")),
0091 rhoMax_(config.getParameter<double>("rhoMax")),
0092 rhoScale_(config.getParameter<double>("rhoScale")),
0093 effectiveAreas_(config.getParameter<std::vector<double>>("effectiveAreas")),
0094 absEtaLowEdges_(config.getParameter<std::vector<double>>("absEtaLowEdges")) {
0095 if (doRhoCorrection_) {
0096 if (absEtaLowEdges_.size() != effectiveAreas_.size())
0097 throw cms::Exception("IncompatibleVects") << "absEtaLowEdges and effectiveAreas should be of the same size. \n";
0098
0099 if (absEtaLowEdges_.at(0) != 0.0)
0100 throw cms::Exception("IncompleteCoverage") << "absEtaLowEdges should start from 0. \n";
0101
0102 for (unsigned int aIt = 0; aIt < absEtaLowEdges_.size() - 1; aIt++) {
0103 if (!(absEtaLowEdges_.at(aIt) < absEtaLowEdges_.at(aIt + 1)))
0104 throw cms::Exception("ImproperBinning") << "absEtaLowEdges entries should be in increasing order. \n";
0105 }
0106 }
0107
0108 std::string recoCandidateProducerName = "recoCandidateProducer";
0109 if ((typeid(HLTHcalPFClusterIsolationProducer<T1>) ==
0110 typeid(HLTHcalPFClusterIsolationProducer<reco::RecoEcalCandidate>)))
0111 recoCandidateProducerName = "recoEcalCandidateProducer";
0112 recoCandidateProducer_ = consumes<T1Collection>(config.getParameter<edm::InputTag>(recoCandidateProducerName));
0113
0114 produces<T1IsolationMap>();
0115 }
0116
0117 template <typename T1>
0118 HLTHcalPFClusterIsolationProducer<T1>::~HLTHcalPFClusterIsolationProducer() {}
0119
0120 template <typename T1>
0121 void HLTHcalPFClusterIsolationProducer<T1>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0122 std::string recoCandidateProducerName = "recoCandidateProducer";
0123 if ((typeid(HLTHcalPFClusterIsolationProducer<T1>) ==
0124 typeid(HLTHcalPFClusterIsolationProducer<reco::RecoEcalCandidate>)))
0125 recoCandidateProducerName = "recoEcalCandidateProducer";
0126
0127 edm::ParameterSetDescription desc;
0128 desc.add<edm::InputTag>(recoCandidateProducerName, edm::InputTag("hltL1SeededRecoEcalCandidatePF"));
0129 desc.add<edm::InputTag>("pfClusterProducerHCAL", edm::InputTag("hltParticleFlowClusterHCAL"));
0130 desc.ifValue(
0131 edm::ParameterDescription<bool>("useHF", false, true),
0132 true >> (edm::ParameterDescription<edm::InputTag>(
0133 "pfClusterProducerHFEM", edm::InputTag("hltParticleFlowClusterHFEM"), true) and
0134 edm::ParameterDescription<edm::InputTag>(
0135 "pfClusterProducerHFHAD", edm::InputTag("hltParticleFlowClusterHFHAD"), true)) or
0136 false >> (edm::ParameterDescription<edm::InputTag>("pfClusterProducerHFEM", edm::InputTag(""), true) and
0137 edm::ParameterDescription<edm::InputTag>("pfClusterProducerHFHAD", edm::InputTag(""), true)));
0138 desc.add<edm::InputTag>("rhoProducer", edm::InputTag("fixedGridRhoFastjetAllCalo"));
0139 desc.add<bool>("doRhoCorrection", false);
0140 desc.add<double>("rhoMax", 9.9999999E7);
0141 desc.add<double>("rhoScale", 1.0);
0142 desc.add<double>("drMax", 0.3);
0143 desc.add<double>("drVetoBarrel", 0.0);
0144 desc.add<double>("drVetoEndcap", 0.0);
0145 desc.add<double>("etaStripBarrel", 0.0);
0146 desc.add<double>("etaStripEndcap", 0.0);
0147 desc.add<double>("energyBarrel", 0.0);
0148 desc.add<double>("energyEndcap", 0.0);
0149 desc.add<bool>("useEt", true);
0150 desc.add<std::vector<double>>("effectiveAreas", {0.2, 0.25});
0151 desc.add<std::vector<double>>("absEtaLowEdges", {0.0, 1.479});
0152 descriptions.add(defaultModuleLabel<HLTHcalPFClusterIsolationProducer<T1>>(), desc);
0153 }
0154
0155 template <typename T1>
0156 void HLTHcalPFClusterIsolationProducer<T1>::produce(edm::StreamID sid,
0157 edm::Event& iEvent,
0158 const edm::EventSetup&) const {
0159 edm::Handle<double> rhoHandle;
0160 double rho = 0.0;
0161 if (doRhoCorrection_) {
0162 iEvent.getByToken(rhoProducer_, rhoHandle);
0163 rho = *(rhoHandle.product());
0164 }
0165
0166 if (rho > rhoMax_)
0167 rho = rhoMax_;
0168
0169 rho = rho * rhoScale_;
0170
0171 edm::Handle<T1Collection> recoCandHandle;
0172
0173 std::vector<edm::Handle<reco::PFClusterCollection>> clusterHandles;
0174 edm::Handle<reco::PFClusterCollection> clusterHcalHandle;
0175 edm::Handle<reco::PFClusterCollection> clusterHfemHandle;
0176 edm::Handle<reco::PFClusterCollection> clusterHfhadHandle;
0177
0178 iEvent.getByToken(recoCandidateProducer_, recoCandHandle);
0179 iEvent.getByToken(pfClusterProducerHCAL_, clusterHcalHandle);
0180
0181 clusterHandles.push_back(clusterHcalHandle);
0182
0183 if (useHF_) {
0184 iEvent.getByToken(pfClusterProducerHFEM_, clusterHfemHandle);
0185 clusterHandles.push_back(clusterHfemHandle);
0186 iEvent.getByToken(pfClusterProducerHFHAD_, clusterHfhadHandle);
0187 clusterHandles.push_back(clusterHfhadHandle);
0188 }
0189
0190 T1IsolationMap recoCandMap(recoCandHandle);
0191 HcalPFClusterIsolation<T1> isoAlgo(
0192 drMax_, drVetoBarrel_, drVetoEndcap_, etaStripBarrel_, etaStripEndcap_, energyBarrel_, energyEndcap_, useEt_);
0193
0194 for (unsigned int iReco = 0; iReco < recoCandHandle->size(); iReco++) {
0195 T1Ref candRef(recoCandHandle, iReco);
0196
0197 float sum = isoAlgo.getSum(candRef, clusterHandles);
0198
0199 if (doRhoCorrection_) {
0200 int iEA = -1;
0201 auto cEta = std::abs(candRef->eta());
0202 for (int bIt = absEtaLowEdges_.size() - 1; bIt > -1; bIt--) {
0203 if (cEta > absEtaLowEdges_.at(bIt)) {
0204 iEA = bIt;
0205 break;
0206 }
0207 }
0208
0209 sum = sum - rho * effectiveAreas_.at(iEA);
0210 }
0211
0212 recoCandMap.insert(candRef, sum);
0213 }
0214
0215 iEvent.put(std::make_unique<T1IsolationMap>(recoCandMap));
0216 }
0217
0218 typedef HLTHcalPFClusterIsolationProducer<reco::RecoEcalCandidate> EgammaHLTHcalPFClusterIsolationProducer;
0219 typedef HLTHcalPFClusterIsolationProducer<reco::RecoChargedCandidate> MuonHLTHcalPFClusterIsolationProducer;
0220
0221 DEFINE_FWK_MODULE(EgammaHLTHcalPFClusterIsolationProducer);
0222 DEFINE_FWK_MODULE(MuonHLTHcalPFClusterIsolationProducer);
