IsolationEstimators.hh
// -*- C++ -*-
#ifndef RIVET_IsolationEstimators_HH
#define RIVET_IsolationEstimators_HH
#include "Rivet/Math/MathUtils.hh"
#include "Rivet/Particle.hh"
#include "Rivet/Jet.hh"
#include <vector>
#include <typeinfo>
namespace Rivet {
template < typename T, typename C >
class IsolationEstimator {
public:
virtual ~IsolationEstimator(){};
virtual double estimate(const T & t, const C & c) const = 0;
virtual int compare(const IsolationEstimator < T, C > *other) const = 0;
virtual bool before(const IsolationEstimator * other) const {
const std::type_info & thisid = typeid(*this);
const std::type_info & otherid = typeid(*other);
if (thisid == otherid) {
return compare(other) < 0;
} else {
return thisid.before(otherid);
}
}
bool operator < (const IsolationEstimator& other) const{
return this->before(&other);
}
};
// An estimator for the sum of the pt of the particles in collection C
// being within radius from t
template < class T, class C >
class PtInConeEstimator : public IsolationEstimator < T, C > {
public:
PtInConeEstimator(double radius, double ptmin = 0.0)
: _radius(radius), _ptmin(ptmin) { }
virtual double estimate(const T & t, const C & c) const {
double ptsum = 0;
for (typename C::const_iterator ic = c.begin(); ic != c.end(); ++ic) {
if (ic->momentum().pT() < _ptmin)
continue;
if (deltaR(t.momentum(), ic->momentum()) < _radius) {
ptsum += ic->momentum().pT();
}
}
return ptsum / t.momentum().pT();
}
virtual int compare(const IsolationEstimator < T, C > *other) const {
const PtInConeEstimator *concreteother = dynamic_cast<const PtInConeEstimator*>(other);
int ptmincmp = cmp(_ptmin, concreteother->getPtMin());
if (ptmincmp != 0)
return ptmincmp;
int radcmp = cmp(_radius, concreteother->getRadius());
if (radcmp != 0)
return radcmp;
return 0;
}
double radius() const {
return _radius;
}
double ptMin() const {
return _ptmin;
}
private:
double _radius;
double _ptmin;
};
// An estimator for the number of particles in collection C
// being within radius from t
template < class T, class C >
class MultiplicityInConeEstimator : public IsolationEstimator < T, C > {
public:
MultiplicityInConeEstimator(double radius, double ptmin = 0.0)
: _radius(radius), _ptmin(ptmin) { }
virtual double estimate(const T & t, const C & c) const {
double npart = 0;
for (typename C::const_iterator ic = c.begin(); ic != c.end(); ++ic) {
if (ic->momentum().pT() < _ptmin)
continue;
if (deltaR(t.momentum(), ic->momentum()) < _radius) {
npart++;
}
}
return npart;
}
virtual int compare(const IsolationEstimator < T, C > *other) const {
const MultiplicityInConeEstimator *concreteother = dynamic_cast < const MultiplicityInConeEstimator * >(other);
int ptmincmp = cmp(_ptmin, concreteother->getPtMin());
if (ptmincmp != 0)
return ptmincmp;
int radcmp = cmp(_radius, concreteother->getRadius());
if (radcmp != 0)
return radcmp;
return 0;
}
double radius() const {
return _radius;
}
double ptMin() const {
return _ptmin;
}
private:
double _radius;
double _ptmin;
};
typedef MultiplicityInConeEstimator < Jet, std::vector < Jet > >JetIsoEstimatorByMultiplicity;
typedef MultiplicityInConeEstimator < Particle, std::vector < Jet > >ParticleFromJetIsoEstimatorByMultiplicity;
typedef MultiplicityInConeEstimator < Particle, std::vector < Particle > >ParticleIsoEstimatorByMultiplicity;
typedef PtInConeEstimator < Jet, std::vector < Jet > >JetIsoEstimatorByPt;
typedef PtInConeEstimator < Particle, std::vector < Jet > >ParticleFromJetIsoEstimatorByPt;
typedef PtInConeEstimator < Particle, std::vector < Particle > >ParticleIsoEstimatorByPt;
template <typename TYPE1, typename TYPE2> struct isohelper {
typedef IsolationEstimator<TYPE1, TYPE2> estimatorhelper;
};
}
#endif
Rivet::cmp
Cmp< T > cmp(const T &t1, const T &t2)
Global helper function for easy creation of Cmp objects.
Definition: Cmp.hh:285