Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:34:07

0001 #ifndef RECOEGAMMA_EGAMMAHLTPRODUCERS_EGAMMAHLTFILTEREDOBJPRODUCER_H
0002 #define RECOEGAMMA_EGAMMAHLTPRODUCERS_EGAMMAHLTFILTEREDOBJPRODUCER_H
0003 
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/stream/EDProducer.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "FWCore/Framework/interface/ConsumesCollector.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 
0014 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0015 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateFwd.h"
0016 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateIsolation.h"
0017 
0018 namespace {
0019   template <typename T>
0020   edm::Handle<T> getHandle(const edm::Event& event, const edm::EDGetTokenT<T>& token) {
0021     edm::Handle<T> handle;
0022     event.getByToken(token, handle);
0023     return handle;
0024   }
0025 
0026 }  // namespace
0027 
0028 template <typename OutCollType>
0029 class EgammaHLTFilteredObjProducer : public edm::stream::EDProducer<> {
0030 public:
0031   class SelectionCut {
0032   public:
0033     SelectionCut(const edm::ParameterSet& pset, edm::ConsumesCollector&& iC)
0034         : ebCut_(pset.getParameter<edm::ParameterSet>("barrelCut")),
0035           eeCut_(pset.getParameter<edm::ParameterSet>("endcapCut")),
0036           varToken_(iC.consumes<reco::RecoEcalCandidateIsolationMap>(pset.getParameter<edm::InputTag>("var"))) {}
0037 
0038     ~SelectionCut() = default;
0039 
0040     bool operator()(const reco::RecoEcalCandidateRef& cand) const {
0041       CutValues cut = std::abs(cand->eta()) < kEBEEEtaBound_ ? ebCut_ : eeCut_;
0042       return cut(*cand, getVar(cand));
0043     }
0044 
0045     float getVar(const reco::RecoEcalCandidateRef& cand) const {
0046       auto res = varHandle_->find(cand);
0047       if (res != varHandle_->end())
0048         return res->val;
0049       else {
0050         //FIX ME: add some provenance info to this
0051         throw cms::Exception("LogicError") << " candidate not found in collection ";
0052       }
0053     }
0054 
0055     void getHandles(const edm::Event& event) { event.getByToken(varToken_, varHandle_); }
0056 
0057   private:
0058     struct CutValues {
0059       float cut;
0060       float cutOverE;
0061       float cutOverE2;
0062       bool useEt;
0063       std::function<bool(float, float)> compFunc;
0064 
0065       CutValues(const edm::ParameterSet& pset)
0066           : cut(pset.getParameter<double>("cut")),
0067             cutOverE(pset.getParameter<double>("cutOverE")),
0068             cutOverE2(pset.getParameter<double>("cutOverE2")),
0069             useEt(pset.getParameter<bool>("useEt")),
0070             compFunc(std::less<float>()) {}
0071 
0072       bool operator()(const reco::RecoEcalCandidate& cand, float value) const {
0073         if (compFunc(value, cut))
0074           return true;
0075         else {
0076           float energyInv = useEt ? 1. / cand.et() : 1. / cand.energy();
0077           return compFunc(value * energyInv, cutOverE) || compFunc(value * energyInv * energyInv, cutOverE2);
0078         }
0079       }
0080     };
0081 
0082     CutValues ebCut_;
0083     CutValues eeCut_;
0084     edm::EDGetTokenT<reco::RecoEcalCandidateIsolationMap> varToken_;
0085     edm::Handle<reco::RecoEcalCandidateIsolationMap> varHandle_;
0086   };
0087 
0088   explicit EgammaHLTFilteredObjProducer(const edm::ParameterSet& pset);
0089   ~EgammaHLTFilteredObjProducer() override = default;
0090   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0091 
0092   void produce(edm::Event&, const edm::EventSetup&) override;
0093 
0094 private:
0095   //different collection types may need to specialise this function
0096   //eg if you wanted superclusters it would output.push_back(cand->superCluster())
0097   static void addObj(const reco::RecoEcalCandidateRef& cand, OutCollType& output) { output.push_back(cand); }
0098 
0099   edm::EDGetTokenT<reco::RecoEcalCandidateCollection> candsToken_;
0100   std::vector<SelectionCut> cuts_;
0101   float minEtCutEB_;
0102   float minEtCutEE_;
0103   static constexpr float kEBEEEtaBound_ = 1.479;
0104 };
0105 
0106 template <typename OutCollType>
0107 EgammaHLTFilteredObjProducer<OutCollType>::EgammaHLTFilteredObjProducer(const edm::ParameterSet& pset)
0108     : candsToken_(consumes<reco::RecoEcalCandidateCollection>(pset.getParameter<edm::InputTag>("cands"))),
0109       minEtCutEB_(pset.getParameter<double>("minEtCutEB")),
0110       minEtCutEE_(pset.getParameter<double>("minEtCutEE")) {
0111   const auto& cutPsets = pset.getParameter<std::vector<edm::ParameterSet> >("cuts");
0112   for (auto& cutPset : cutPsets) {
0113     cuts_.push_back(SelectionCut(cutPset, consumesCollector()));
0114   }
0115 
0116   produces<OutCollType>();
0117 }
0118 
0119 template <typename OutCollType>
0120 void EgammaHLTFilteredObjProducer<OutCollType>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0121   edm::ParameterSetDescription desc;
0122   desc.add<edm::InputTag>("cands", edm::InputTag("hltEgammaCandidates"));
0123   desc.add<double>("minEtCutEB", 0);
0124   desc.add<double>("minEtCutEE", 0);
0125 
0126   edm::ParameterSetDescription cutsDesc;
0127   edm::ParameterSetDescription regionCutsDesc;
0128   regionCutsDesc.add<double>("cut", -1);
0129   regionCutsDesc.add<double>("cutOverE", -1);
0130   regionCutsDesc.add<double>("cutOverE2", -1);
0131   regionCutsDesc.add<bool>("useEt", false);
0132   edm::ParameterSet cutDefaults;
0133   cutDefaults.addParameter<double>("cutOverE", 0.2);
0134   cutDefaults.addParameter<double>("useEt", false);
0135 
0136   cutsDesc.add<edm::ParameterSetDescription>("barrelCut", regionCutsDesc);
0137   cutsDesc.add<edm::ParameterSetDescription>("endcapCut", regionCutsDesc);
0138   cutsDesc.add<edm::InputTag>("var", edm::InputTag("hltEgammaHoverE"));
0139 
0140   edm::ParameterSet defaults;
0141   defaults.addParameter<edm::InputTag>("var", edm::InputTag("hltEgammaHoverE"));
0142   defaults.addParameter<edm::ParameterSet>("barrelCut", cutDefaults);
0143   defaults.addParameter<edm::ParameterSet>("endcapCut", cutDefaults);
0144   desc.addVPSet("cuts", cutsDesc, std::vector<edm::ParameterSet>{defaults});
0145   descriptions.addWithDefaultLabel(desc);
0146 }
0147 
0148 template <typename OutCollType>
0149 void EgammaHLTFilteredObjProducer<OutCollType>::produce(edm::Event& iEvent, const edm::EventSetup&) {
0150   for (auto& cut : cuts_)
0151     cut.getHandles(iEvent);
0152   auto candsHandle = getHandle(iEvent, candsToken_);
0153 
0154   auto output = std::make_unique<OutCollType>();
0155 
0156   for (size_t candNr = 0; candNr < candsHandle->size(); candNr++) {
0157     reco::RecoEcalCandidateRef candRef(candsHandle, candNr);
0158     bool passAllCuts = true;
0159 
0160     float minEtCut = std::abs(candRef->eta()) < kEBEEEtaBound_ ? minEtCutEB_ : minEtCutEE_;
0161     if (candRef->et() < minEtCut)
0162       passAllCuts = false;
0163 
0164     if (passAllCuts) {
0165       for (const auto& cut : cuts_) {
0166         if (!cut(candRef)) {
0167           passAllCuts = false;
0168           break;
0169         }
0170       }
0171     }
0172     if (passAllCuts)
0173       addObj(candRef, *output);
0174   }
0175 
0176   iEvent.put(std::move(output));
0177 }
0178 
0179 #endif