ZToMuMuIsoDepositSelector.cc

Go to the documentation of this file.

/* \class ZToMuMuIsolationSelector
 *
 * \author Luca Lista, INFN
 *
 */

struct IsolatedSelector {
  IsolatedSelector(double cut) : cut_(cut) {}
  bool operator()(double i1, double i2) const { return i1 < cut_ && i2 < cut_; }
  double cut() const { return cut_; }

private:
  double cut_;
};

struct NonIsolatedSelector {
  NonIsolatedSelector(double cut) : isolated_(cut) {}
  bool operator()(double i1, double i2) const { return !isolated_(i1, i2); }
  double cut() const { return isolated_.cut(); }

private:
  IsolatedSelector isolated_;
};

struct OneNonIsolatedSelector {
  OneNonIsolatedSelector(double cut) : cut_(cut) {}
  bool operator()(double i1, double i2) const { return (i1 < cut_ && i2 >= cut_) || (i1 >= cut_ && i2 < cut_); }
  double cut() const { return cut_; }

private:
  double cut_;
};

struct TwoNonIsolatedSelector {
  TwoNonIsolatedSelector(double cut) : cut_(cut) {}
  bool operator()(double i1, double i2) const { return i1 >= cut_ && i2 >= cut_; }
  double cut() const { return cut_; }

private:
  double cut_;
};

#include "FWCore/Framework/interface/ConsumesCollector.h"
#include "DataFormats/PatCandidates/interface/Muon.h"
#include "DataFormats/PatCandidates/interface/GenericParticle.h"
#include "FWCore/Utilities/interface/EDMException.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "DataFormats/RecoCandidate/interface/IsoDepositVetos.h"
#include "DataFormats/RecoCandidate/interface/IsoDepositDirection.h"
#include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
#include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
#include "DataFormats/PatCandidates/interface/Isolation.h"

using namespace reco;
using namespace isodeposit;

template <typename Isolator>
class ZToMuMuIsoDepositSelector {
public:
  ZToMuMuIsoDepositSelector(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC)
      : isolator_(cfg.template getParameter<double>("isoCut")),
        ptThreshold(cfg.getUntrackedParameter<double>("ptThreshold")),
        etEcalThreshold(cfg.getUntrackedParameter<double>("etEcalThreshold")),
        etHcalThreshold(cfg.getUntrackedParameter<double>("etHcalThreshold")),
        dRVetoTrk(cfg.getUntrackedParameter<double>("deltaRVetoTrk")),
        dRTrk(cfg.getUntrackedParameter<double>("deltaRTrk")),
        dREcal(cfg.getUntrackedParameter<double>("deltaREcal")),
        dRHcal(cfg.getUntrackedParameter<double>("deltaRHcal")),
        alpha(cfg.getUntrackedParameter<double>("alpha")),
        beta(cfg.getUntrackedParameter<double>("beta")),
        relativeIsolation(cfg.template getParameter<bool>("relativeIsolation")) {}

  template <typename T>
  double isolation(const T* t) const {
    const pat::IsoDeposit* trkIso = t->isoDeposit(pat::TrackIso);
    const pat::IsoDeposit* ecalIso = t->isoDeposit(pat::EcalIso);
    const pat::IsoDeposit* hcalIso = t->isoDeposit(pat::HcalIso);

    Direction dir = Direction(t->eta(), t->phi());

    IsoDeposit::AbsVetos vetosTrk;
    vetosTrk.push_back(new ConeVeto(dir, dRVetoTrk));
    vetosTrk.push_back(new ThresholdVeto(ptThreshold));

    IsoDeposit::AbsVetos vetosEcal;
    vetosEcal.push_back(new ConeVeto(dir, 0.));
    vetosEcal.push_back(new ThresholdVeto(etEcalThreshold));

    IsoDeposit::AbsVetos vetosHcal;
    vetosHcal.push_back(new ConeVeto(dir, 0.));
    vetosHcal.push_back(new ThresholdVeto(etHcalThreshold));

    double isovalueTrk = (trkIso->sumWithin(dRTrk, vetosTrk));
    double isovalueEcal = (ecalIso->sumWithin(dREcal, vetosEcal));
    double isovalueHcal = (hcalIso->sumWithin(dRHcal, vetosHcal));

    double iso =
        alpha * (((1 + beta) / 2 * isovalueEcal) + ((1 - beta) / 2 * isovalueHcal)) + ((1 - alpha) * isovalueTrk);
    if (relativeIsolation)
      iso /= t->pt();
    return iso;
  }

  double candIsolation(const reco::Candidate* c) const {
    const pat::Muon* mu = dynamic_cast<const pat::Muon*>(c);
    if (mu != nullptr)
      return isolation(mu);
    const pat::GenericParticle* trk = dynamic_cast<const pat::GenericParticle*>(c);
    if (trk != nullptr)
      return isolation(trk);
    throw edm::Exception(edm::errors::InvalidReference)
        << "Candidate daughter #0 is neither pat::Muons nor pat::GenericParticle\n";
    return -1;
  }
  bool operator()(const reco::Candidate& z) const {
    if (z.numberOfDaughters() != 2)
      throw edm::Exception(edm::errors::InvalidReference)
          << "Candidate has " << z.numberOfDaughters() << " daughters, 2 expected\n";
    const reco::Candidate* dau0 = z.daughter(0);
    const reco::Candidate* dau1 = z.daughter(1);
    if (!(dau0->hasMasterClone() && dau1->hasMasterClone()))
      throw edm::Exception(edm::errors::InvalidReference) << "Candidate daughters have no master clone\n";
    const reco::Candidate *m0 = &*dau0->masterClone(), *m1 = &*dau1->masterClone();
    return isolator_(candIsolation(m0), candIsolation(m1));
  }

private:
  Isolator isolator_;
  double ptThreshold, etEcalThreshold, etHcalThreshold, dRVetoTrk, dRTrk, dREcal, dRHcal, alpha, beta;
  bool relativeIsolation;
};

#include "CommonTools/UtilAlgos/interface/SingleObjectSelector.h"
#include "CommonTools/UtilAlgos/interface/AndSelector.h"
#include "CommonTools/UtilAlgos/interface/StringCutObjectSelector.h"

typedef SingleObjectSelector<
    reco::CandidateView,
    AndSelector<ZToMuMuIsoDepositSelector<IsolatedSelector>, StringCutObjectSelector<reco::Candidate> > >
    ZToMuMuIsolatedIDSelector;

typedef SingleObjectSelector<
    reco::CandidateView,
    AndSelector<ZToMuMuIsoDepositSelector<NonIsolatedSelector>, StringCutObjectSelector<reco::Candidate> > >
    ZToMuMuNonIsolatedIDSelector;

typedef SingleObjectSelector<
    reco::CandidateView,
    AndSelector<ZToMuMuIsoDepositSelector<OneNonIsolatedSelector>, StringCutObjectSelector<reco::Candidate> > >
    ZToMuMuOneNonIsolatedIDSelector;

typedef SingleObjectSelector<
    reco::CandidateView,
    AndSelector<ZToMuMuIsoDepositSelector<TwoNonIsolatedSelector>, StringCutObjectSelector<reco::Candidate> > >
    ZToMuMuTwoNonIsolatedIDSelector;

#include "FWCore/Framework/interface/MakerMacros.h"

DEFINE_FWK_MODULE(ZToMuMuIsolatedIDSelector);
DEFINE_FWK_MODULE(ZToMuMuNonIsolatedIDSelector);
DEFINE_FWK_MODULE(ZToMuMuOneNonIsolatedIDSelector);
DEFINE_FWK_MODULE(ZToMuMuTwoNonIsolatedIDSelector);