plugins/cuts/GsfEleCalPFClusterIsoCut.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 GsfEleCalPFClusterIsoCut : public CutApplicatorWithEventContentBase {
0006 public:
0007   enum IsoType { UNDEF = -1, ISO_ECAL = 0, ISO_HCAL = 1 };
0008
0009   GsfEleCalPFClusterIsoCut(const edm::ParameterSet& c);
0010
0011   result_type operator()(const reco::GsfElectronPtr&) const final;
0012
0013   void setConsumes(edm::ConsumesCollector&) final;
0014   void getEventContent(const edm::EventBase&) final;
0015
0016   double value(const reco::CandidatePtr& cand) const final;
0017
0018   CandidateType candidateType() const final { return ELECTRON; }
0019
0020 private:
0021   // Cut values
0022   const float isoCutEBLowPt_, isoCutEBHighPt_, isoCutEELowPt_, isoCutEEHighPt_;
0023   // Configuration
0024   const int isoType_;
0025   const float ptCutOff_;
0026   const float barrelCutOff_;
0027   bool isRelativeIso_;
0028   // Effective area constants
0029   EffectiveAreas effectiveAreas_;
0030   // The rho
0031   edm::Handle<double> rhoHandle_;
0032
0033   constexpr static char rhoString_[] = "rho";
0034 };
0035
0036 constexpr char GsfEleCalPFClusterIsoCut::rhoString_[];
0037
0038 DEFINE_EDM_PLUGIN(CutApplicatorFactory, GsfEleCalPFClusterIsoCut, "GsfEleCalPFClusterIsoCut");
0039
0040 GsfEleCalPFClusterIsoCut::GsfEleCalPFClusterIsoCut(const edm::ParameterSet& c)
0041     : CutApplicatorWithEventContentBase(c),
0042       isoCutEBLowPt_(c.getParameter<double>("isoCutEBLowPt")),
0043       isoCutEBHighPt_(c.getParameter<double>("isoCutEBHighPt")),
0044       isoCutEELowPt_(c.getParameter<double>("isoCutEELowPt")),
0045       isoCutEEHighPt_(c.getParameter<double>("isoCutEEHighPt")),
0046       isoType_(c.getParameter<int>("isoType")),
0047       ptCutOff_(c.getParameter<double>("ptCutOff")),
0048       barrelCutOff_(c.getParameter<double>("barrelCutOff")),
0049       isRelativeIso_(c.getParameter<bool>("isRelativeIso")),
0050       effectiveAreas_((c.getParameter<edm::FileInPath>("effAreasConfigFile")).fullPath()) {
0051   edm::InputTag rhoTag = c.getParameter<edm::InputTag>("rho");
0052   contentTags_.emplace(rhoString_, rhoTag);
0053 }
0054
0055 void GsfEleCalPFClusterIsoCut::setConsumes(edm::ConsumesCollector& cc) {
0056   auto rho = cc.consumes<double>(contentTags_[rhoString_]);
0057   contentTokens_.emplace(rhoString_, rho);
0058 }
0059
0060 void GsfEleCalPFClusterIsoCut::getEventContent(const edm::EventBase& ev) {
0061   ev.getByLabel(contentTags_[rhoString_], rhoHandle_);
0062 }
0063
0064 CutApplicatorBase::result_type GsfEleCalPFClusterIsoCut::operator()(const reco::GsfElectronPtr& cand) const {
0065   // Establish the cut value
0066   double absEta = std::abs(cand->superCluster()->eta());
0067
0068   const pat::Electron* elPat = dynamic_cast<const pat::Electron*>(cand.get());
0069   if (!elPat) {
0070     throw cms::Exception("ERROR: this VID selection is meant to be run on miniAOD/PAT only")
0071         << std::endl
0072         << "Change input format to PAT/miniAOD or contact Egamma experts" << std::endl
0073         << std::endl;
0074   }
0075
0076   const float isoCut = (cand->pt() < ptCutOff_ ? (absEta < barrelCutOff_ ? isoCutEBLowPt_ : isoCutEELowPt_)
0077                                                : (absEta < barrelCutOff_ ? isoCutEBHighPt_ : isoCutEEHighPt_));
0078
0079   const float eA = effectiveAreas_.getEffectiveArea(absEta);
0080   const float rho = rhoHandle_.isValid() ? (float)(*rhoHandle_) : 0;  // std::max likes float arguments
0081
0082   float isoValue = -999;
0083   if (isoType_ == ISO_ECAL) {
0084     isoValue = elPat->ecalPFClusterIso();
0085   } else if (isoType_ == ISO_HCAL) {
0086     isoValue = elPat->hcalPFClusterIso();
0087   } else {
0088     throw cms::Exception("ERROR: unknown type requested for PF cluster isolation.")
0089         << std::endl
0090         << "Check VID configuration." << std::endl;
0091   }
0092   float isoValueCorr = std::max(0.0f, isoValue - rho * eA);
0093
0094   // Apply the cut and return the result
0095   // Scale by pT if the relative isolation is requested but avoid division by 0
0096   return isoValueCorr < isoCut * (isRelativeIso_ ? cand->pt() : 1.);
0097 }
0098
0099 double GsfEleCalPFClusterIsoCut::value(const reco::CandidatePtr& cand) const {
0100   reco::GsfElectronPtr ele(cand);
0101   // Establish the cut value
0102   double absEta = std::abs(ele->superCluster()->eta());
0103
0104   const pat::Electron* elPat = dynamic_cast<const pat::Electron*>(ele.get());
0105   if (!elPat) {
0106     throw cms::Exception("ERROR: this VID selection is meant to be run on miniAOD/PAT only")
0107         << std::endl
0108         << "Change input format to PAT/miniAOD or contact Egamma experts" << std::endl
0109         << std::endl;
0110   }
0111
0112   const float eA = effectiveAreas_.getEffectiveArea(absEta);
0113   const float rho = rhoHandle_.isValid() ? (float)(*rhoHandle_) : 0;  // std::max likes float arguments
0114
0115   float isoValue = -999;
0116   if (isoType_ == ISO_ECAL) {
0117     isoValue = elPat->ecalPFClusterIso();
0118   } else if (isoType_ == ISO_HCAL) {
0119     isoValue = elPat->hcalPFClusterIso();
0120   } else {
0121     throw cms::Exception("ERROR: unknown type requested for PF cluster isolation.")
0122         << std::endl
0123         << "Check VID configuration." << std::endl;
0124   }
0125   float isoValueCorr = std::max(0.0f, isoValue - rho * eA);
0126
0127   // Divide by pT if the relative isolation is requested
0128   if (isRelativeIso_)
0129     isoValueCorr /= ele->pt();
0130
0131   // Apply the cut and return the result
0132   return isoValueCorr;
0133 }