plugins/cuts/GsfEleEffAreaPFIsoCut.cc

0001 #include "PhysicsTools/SelectorUtils/interface/CutApplicatorWithEventContentBase.h"
0002 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0003 #include "CommonTools/Egamma/interface/EffectiveAreas.h"
0004
0005 class GsfEleEffAreaPFIsoCut : public CutApplicatorWithEventContentBase {
0006 public:
0007   GsfEleEffAreaPFIsoCut(const edm::ParameterSet& c);
0008
0009   result_type operator()(const reco::GsfElectronPtr&) const final;
0010
0011   void setConsumes(edm::ConsumesCollector&) final;
0012   void getEventContent(const edm::EventBase&) final;
0013
0014   double value(const reco::CandidatePtr& cand) const final;
0015
0016   CandidateType candidateType() const final { return ELECTRON; }
0017
0018 private:
0019   // Cut values
0020   const float _isoCutEBLowPt, _isoCutEBHighPt, _isoCutEELowPt, _isoCutEEHighPt;
0021   // Configuration
0022   const float _ptCutOff;
0023   const float _barrelCutOff;
0024   bool _isRelativeIso;
0025   // Effective area constants
0026   EffectiveAreas _effectiveAreas;
0027   // The rho
0028   edm::Handle<double> _rhoHandle;
0029
0030   constexpr static char rhoString_[] = "rho";
0031 };
0032
0033 constexpr char GsfEleEffAreaPFIsoCut::rhoString_[];
0034
0035 DEFINE_EDM_PLUGIN(CutApplicatorFactory, GsfEleEffAreaPFIsoCut, "GsfEleEffAreaPFIsoCut");
0036
0037 GsfEleEffAreaPFIsoCut::GsfEleEffAreaPFIsoCut(const edm::ParameterSet& c)
0038     : CutApplicatorWithEventContentBase(c),
0039       _isoCutEBLowPt(c.getParameter<double>("isoCutEBLowPt")),
0040       _isoCutEBHighPt(c.getParameter<double>("isoCutEBHighPt")),
0041       _isoCutEELowPt(c.getParameter<double>("isoCutEELowPt")),
0042       _isoCutEEHighPt(c.getParameter<double>("isoCutEEHighPt")),
0043       _ptCutOff(c.getParameter<double>("ptCutOff")),
0044       _barrelCutOff(c.getParameter<double>("barrelCutOff")),
0045       _isRelativeIso(c.getParameter<bool>("isRelativeIso")),
0046       _effectiveAreas((c.getParameter<edm::FileInPath>("effAreasConfigFile")).fullPath()) {
0047   edm::InputTag rhoTag = c.getParameter<edm::InputTag>("rho");
0048   contentTags_.emplace(rhoString_, rhoTag);
0049 }
0050
0051 void GsfEleEffAreaPFIsoCut::setConsumes(edm::ConsumesCollector& cc) {
0052   auto rho = cc.consumes<double>(contentTags_[rhoString_]);
0053   contentTokens_.emplace(rhoString_, rho);
0054 }
0055
0056 void GsfEleEffAreaPFIsoCut::getEventContent(const edm::EventBase& ev) {
0057   ev.getByLabel(contentTags_[rhoString_], _rhoHandle);
0058 }
0059
0060 CutApplicatorBase::result_type GsfEleEffAreaPFIsoCut::operator()(const reco::GsfElectronPtr& cand) const {
0061   // Establish the cut value
0062   double absEta = std::abs(cand->superCluster()->eta());
0063   const float isoCut = (cand->pt() < _ptCutOff ? (absEta < _barrelCutOff ? _isoCutEBLowPt : _isoCutEELowPt)
0064                                                : (absEta < _barrelCutOff ? _isoCutEBHighPt : _isoCutEEHighPt));
0065
0066   // Compute the combined isolation with effective area correction
0067   const reco::GsfElectron::PflowIsolationVariables& pfIso = cand->pfIsolationVariables();
0068   const float chad = pfIso.sumChargedHadronPt;
0069   const float nhad = pfIso.sumNeutralHadronEt;
0070   const float pho = pfIso.sumPhotonEt;
0071   const float eA = _effectiveAreas.getEffectiveArea(absEta);
0072   const float rho = _rhoHandle.isValid() ? (float)(*_rhoHandle) : 0;  // std::max likes float arguments
0073   const float iso = chad + std::max(0.0f, nhad + pho - rho * eA);
0074
0075   // Apply the cut and return the result
0076   // Scale by pT if the relative isolation is requested but avoid division by 0
0077   return iso < isoCut * (_isRelativeIso ? cand->pt() : 1.);
0078 }
0079
0080 double GsfEleEffAreaPFIsoCut::value(const reco::CandidatePtr& cand) const {
0081   reco::GsfElectronPtr ele(cand);
0082   // Establish the cut value
0083   double absEta = std::abs(ele->superCluster()->eta());
0084
0085   // Compute the combined isolation with effective area correction
0086   const reco::GsfElectron::PflowIsolationVariables& pfIso = ele->pfIsolationVariables();
0087   const float chad = pfIso.sumChargedHadronPt;
0088   const float nhad = pfIso.sumNeutralHadronEt;
0089   const float pho = pfIso.sumPhotonEt;
0090   float eA = _effectiveAreas.getEffectiveArea(absEta);
0091   float rho = (float)(*_rhoHandle);  // std::max likes float arguments
0092   float iso = chad + std::max(0.0f, nhad + pho - rho * eA);
0093
0094   // Divide by pT if the relative isolation is requested
0095   if (_isRelativeIso)
0096     iso /= ele->pt();
0097
0098   // Apply the cut and return the result
0099   return iso;
0100 }