Rivet  1.8.3
IsolationEstimators.hh
1 // -*- C++ -*-
2 #ifndef RIVET_IsolationEstimators_HH
3 #define RIVET_IsolationEstimators_HH
4 
5 #include "Rivet/Math/MathUtils.hh"
6 #include "Rivet/Particle.hh"
7 #include "Rivet/Jet.hh"
8 
9 #include <vector>
10 #include <typeinfo>
11 
12 namespace Rivet {
13 
15 
16  template < typename T, typename C >
17  class IsolationEstimator {
18 
19  public:
20 
21  virtual ~IsolationEstimator(){};
22 
23  virtual double estimate(const T & t, const C & c) const = 0;
24 
25  virtual int compare(const IsolationEstimator < T, C > *other) const = 0;
26 
27  virtual bool before(const IsolationEstimator * other) const {
28  const std::type_info & thisid = typeid(*this);
29  const std::type_info & otherid = typeid(*other);
30  if (thisid == otherid) {
31  return compare(other) < 0;
32  } else {
33  return thisid.before(otherid);
34  }
35  }
36 
37  bool operator < (const IsolationEstimator& other) const{
38  return this->before(&other);
39  }
40  };
41 
42 
43  // An estimator for the sum of the pt of the particles in collection C
44  // being within radius from t
45  template < class T, class C >
46  class PtInConeEstimator : public IsolationEstimator < T, C > {
47  public:
48  PtInConeEstimator(double radius, double ptmin = 0.0)
49  : _radius(radius), _ptmin(ptmin) { }
50 
51  virtual double estimate(const T & t, const C & c) const {
52  double ptsum = 0;
53  for (typename C::const_iterator ic = c.begin(); ic != c.end(); ++ic) {
54  if (ic->momentum().pT() < _ptmin)
55  continue;
56  if (deltaR(t.momentum(), ic->momentum()) < _radius) {
57  ptsum += ic->momentum().pT();
58  }
59  }
60  return ptsum / t.momentum().pT();
61  }
62 
63  virtual int compare(const IsolationEstimator < T, C > *other) const {
64  const PtInConeEstimator *concreteother = dynamic_cast<const PtInConeEstimator*>(other);
65  int ptmincmp = cmp(_ptmin, concreteother->getPtMin());
66  if (ptmincmp != 0)
67  return ptmincmp;
68  int radcmp = cmp(_radius, concreteother->getRadius());
69  if (radcmp != 0)
70  return radcmp;
71  return 0;
72  }
73 
74  double radius() const {
75  return _radius;
76  }
77 
78  double ptMin() const {
79  return _ptmin;
80  }
81  private:
82  double _radius;
83  double _ptmin;
84  };
85 
86 
87  // An estimator for the number of particles in collection C
88  // being within radius from t
89  template < class T, class C >
90  class MultiplicityInConeEstimator : public IsolationEstimator < T, C > {
91  public:
92  MultiplicityInConeEstimator(double radius, double ptmin = 0.0)
93  : _radius(radius), _ptmin(ptmin) { }
94 
95  virtual double estimate(const T & t, const C & c) const {
96  double npart = 0;
97  for (typename C::const_iterator ic = c.begin(); ic != c.end(); ++ic) {
98  if (ic->momentum().pT() < _ptmin)
99  continue;
100  if (deltaR(t.momentum(), ic->momentum()) < _radius) {
101  npart++;
102  }
103  }
104  return npart;
105  }
106 
107  virtual int compare(const IsolationEstimator < T, C > *other) const {
108  const MultiplicityInConeEstimator *concreteother = dynamic_cast < const MultiplicityInConeEstimator * >(other);
109  int ptmincmp = cmp(_ptmin, concreteother->getPtMin());
110  if (ptmincmp != 0)
111  return ptmincmp;
112  int radcmp = cmp(_radius, concreteother->getRadius());
113  if (radcmp != 0)
114  return radcmp;
115  return 0;
116  }
117 
118  double radius() const {
119  return _radius;
120  }
121 
122  double ptMin() const {
123  return _ptmin;
124  }
125 
126  private:
127  double _radius;
128  double _ptmin;
129  };
130 
131 
133  typedef MultiplicityInConeEstimator < Jet, std::vector < Jet > >JetIsoEstimatorByMultiplicity;
134  typedef MultiplicityInConeEstimator < Particle, std::vector < Jet > >ParticleFromJetIsoEstimatorByMultiplicity;
135  typedef MultiplicityInConeEstimator < Particle, std::vector < Particle > >ParticleIsoEstimatorByMultiplicity;
136 
137  typedef PtInConeEstimator < Jet, std::vector < Jet > >JetIsoEstimatorByPt;
138  typedef PtInConeEstimator < Particle, std::vector < Jet > >ParticleFromJetIsoEstimatorByPt;
139  typedef PtInConeEstimator < Particle, std::vector < Particle > >ParticleIsoEstimatorByPt;
140  template <typename TYPE1, typename TYPE2> struct isohelper {
141  typedef IsolationEstimator<TYPE1, TYPE2> estimatorhelper;
142  };
143 
144 
146 
147 }
148 
149 #endif
Cmp< T > cmp(const T &t1, const T &t2)
Global helper function for easy creation of Cmp objects.
Definition: Cmp.hh:285