Rivet  3.1.3
ParticleBaseUtils.hh
1 #ifndef RIVET_PARTICLEBASEUTILS_HH
2 #define RIVET_PARTICLEBASEUTILS_HH
3 
4 #include "Rivet/ParticleBase.hh"
5 
6 namespace Rivet {
7 
8 
11 
17 
19  using ParticleBaseSelector = function<bool(const ParticleBase&)>;
21  using ParticleBaseSorter = function<bool(const ParticleBase&, const ParticleBase&)>;
22 
23 
26  virtual bool operator()(const ParticleBase& p) const = 0;
27  virtual ~BoolParticleBaseFunctor() {}
28  };
29 
30 
32  struct PtGtr : public BoolParticleBaseFunctor {
33  PtGtr(double pt) : ptcut(pt) { }
34  PtGtr(const FourMomentum& p) : ptcut(p.pT()) { }
35  bool operator()(const ParticleBase& p) const { return p.pT() > ptcut; }
36  double ptcut;
37  };
38  using pTGtr = PtGtr;
39  using ptGtr = PtGtr;
40 
42  struct PtLess : public BoolParticleBaseFunctor {
43  PtLess(const FourMomentum& p) : ptcut(p.pT()) { }
44  PtLess(double pt) : ptcut(pt) { }
45  bool operator()(const ParticleBase& p) const { return p.pT() < ptcut; }
46  double ptcut;
47  };
48  using pTLess = PtLess;
49  using ptLess = PtLess;
50 
53  PtInRange(pair<double, double> ptcuts) : ptcut(ptcuts) { }
54  PtInRange(double ptlow, double pthigh) : PtInRange(make_pair(ptlow, pthigh)) { }
55  PtInRange(const FourMomentum& p1, const FourMomentum& p2) : PtInRange(p1.pT(), p2.pT()) { }
56  bool operator()(const ParticleBase& p) const { return p.pT() >= ptcut.first && p.pT() < ptcut.second; }
57  pair<double,double> ptcut;
58  };
59  using pTInRange = PtInRange;
60  using ptInRange = PtInRange;
61 
62 
64  struct EtaGtr : public BoolParticleBaseFunctor {
65  EtaGtr(double eta) : etacut(eta) { }
66  EtaGtr(const FourMomentum& p) : etacut(p.eta()) { }
67  bool operator()(const ParticleBase& p) const { return p.eta() > etacut; }
68  double etacut;
69  };
70  using etaGtr = EtaGtr;
71 
73  struct EtaLess : public BoolParticleBaseFunctor {
74  EtaLess(double eta) : etacut(eta) { }
75  EtaLess(const FourMomentum& p) : etacut(p.eta()) { }
76  bool operator()(const ParticleBase& p) const { return p.eta() < etacut; }
77  double etacut;
78  };
79  using etaLess = EtaLess;
80 
83  EtaInRange(pair<double, double> etacuts) : etacut(etacuts) { }
84  EtaInRange(double etalow, double etahigh) : EtaInRange(make_pair(etalow, etahigh)) { }
85  EtaInRange(const FourMomentum& p1, const FourMomentum& p2) : EtaInRange(p1.eta(), p2.eta()) { }
86  bool operator()(const ParticleBase& p) const { return p.eta() >= etacut.first && p.eta() < etacut.second; }
87  pair<double,double> etacut;
88  };
89  using etaInRange = EtaInRange;
90 
91 
94  AbsEtaGtr(double abseta) : absetacut(abseta) { }
95  AbsEtaGtr(const FourMomentum& p) : absetacut(p.abseta()) { }
96  bool operator()(const ParticleBase& p) const { return p.abseta() > absetacut; }
97  double absetacut;
98  };
99  using absEtaGtr = AbsEtaGtr;
100  using absetaGtr = AbsEtaGtr;
101 
104  AbsEtaLess(double abseta) : absetacut(abseta) { }
105  AbsEtaLess(const FourMomentum& p) : absetacut(p.abseta()) { }
106  bool operator()(const ParticleBase& p) const { return p.abseta() < absetacut; }
107  double absetacut;
108  };
109  using absEtaLess = AbsEtaLess;
110  using absetaLess = AbsEtaLess;
111 
114  AbsEtaInRange(const pair<double, double>& absetacuts) : absetacut(absetacuts) { }
115  AbsEtaInRange(double absetalow, double absetahigh) : AbsEtaInRange(make_pair(absetalow, absetahigh)) { }
116  AbsEtaInRange(const FourMomentum& p1, const FourMomentum& p2) : AbsEtaInRange(p1.abseta(), p2.abseta()) { }
117  bool operator()(const ParticleBase& p) const { return p.abseta() >= absetacut.first && p.abseta() < absetacut.second; }
118  pair<double,double> absetacut;
119  };
122 
123 
125  struct RapGtr : public BoolParticleBaseFunctor {
126  RapGtr(double rap) : rapcut(rap) { }
127  RapGtr(const FourMomentum& p) : rapcut(p.rap()) { }
128  bool operator()(const ParticleBase& p) const { return p.rap() > rapcut; }
129  double rapcut;
130  };
131  using rapGtr = RapGtr;
132 
135  RapLess(double rap) : rapcut(rap) { }
136  RapLess(const FourMomentum& p) : rapcut(p.rap()) { }
137  bool operator()(const ParticleBase& p) const { return p.rap() < rapcut; }
138  double rapcut;
139  };
140  using rapLess = RapLess;
141 
144  RapInRange(const pair<double, double>& rapcuts) : rapcut(rapcuts) { }
145  RapInRange(double raplow, double raphigh) : RapInRange(make_pair(raplow, raphigh)) { }
146  RapInRange(const FourMomentum& p1, const FourMomentum& p2) : RapInRange(p1.rap(), p2.rap()) { }
147  bool operator()(const ParticleBase& p) const { return p.rap() >= rapcut.first && p.rap() < rapcut.second; }
148  pair<double,double> rapcut;
149  };
150  using rapInRange = RapInRange;
151 
152 
155  AbsRapGtr(double absrap) : absrapcut(absrap) { }
156  AbsRapGtr(const FourMomentum& p) : absrapcut(p.absrap()) { }
157  bool operator()(const ParticleBase& p) const { return p.absrap() > absrapcut; }
158  double absrapcut;
159  };
160  using absRapGtr = AbsRapGtr;
161  using absrapGtr = AbsRapGtr;
162 
165  AbsRapLess(double absrap) : absrapcut(absrap) { }
166  AbsRapLess(const FourMomentum& p) : absrapcut(p.absrap()) { }
167  bool operator()(const ParticleBase& p) const { return p.absrap() < absrapcut; }
168  double absrapcut;
169  };
170  using absRapLess = AbsRapLess;
171  using absrapLess = AbsRapLess;
172 
175  AbsRapInRange(const pair<double, double>& absrapcuts) : absrapcut(absrapcuts) { }
176  AbsRapInRange(double absraplow, double absraphigh) : AbsRapInRange(make_pair(absraplow, absraphigh)) { }
177  AbsRapInRange(const FourMomentum& p1, const FourMomentum& p2) : AbsRapInRange(p1.absrap(), p2.absrap()) { }
178  bool operator()(const ParticleBase& p) const { return p.absrap() >= absrapcut.first && p.absrap() < absrapcut.second; }
179  pair<double,double> absrapcut;
180  };
183 
184 
186 
187 
190  DeltaRGtr(const ParticleBase& vec, double dr, RapScheme scheme=PSEUDORAPIDITY)
191  : refvec(vec.mom()), drcut(dr), rapscheme(scheme) { }
192  DeltaRGtr(const FourMomentum& vec, double dr, RapScheme scheme=PSEUDORAPIDITY)
193  : refvec(vec), drcut(dr), rapscheme(scheme) { }
194  DeltaRGtr(const Vector3& vec, double dr)
195  : drcut(dr), rapscheme(PSEUDORAPIDITY) { refvec.setPx(vec.x()); refvec.setPy(vec.y()); refvec.setPz(vec.z()); }
196  bool operator()(const ParticleBase& p) const { return deltaR(p, refvec, rapscheme) > drcut; }
197  FourMomentum refvec;
198  double drcut;
199  RapScheme rapscheme;
200  };
201  using deltaRGtr = DeltaRGtr;
202 
205  DeltaRLess(const ParticleBase& vec, double dr, RapScheme scheme=PSEUDORAPIDITY)
206  : refvec(vec.mom()), drcut(dr), rapscheme(scheme) { }
207  DeltaRLess(const FourMomentum& vec, double dr, RapScheme scheme=PSEUDORAPIDITY)
208  : refvec(vec), drcut(dr), rapscheme(scheme) { }
209  DeltaRLess(const Vector3& vec, double dr)
210  : drcut(dr), rapscheme(PSEUDORAPIDITY) { refvec.setPx(vec.x()); refvec.setPy(vec.y()); refvec.setPz(vec.z()); }
211  bool operator()(const ParticleBase& p) const { return deltaR(p, refvec, rapscheme) < drcut; }
212  FourMomentum refvec;
213  double drcut;
214  RapScheme rapscheme;
215  };
216  using deltaRLess = DeltaRLess;
217 
220  DeltaRInRange(const ParticleBase& vec, const pair<double,double>& dr, RapScheme scheme=PSEUDORAPIDITY)
221  : refvec(vec.mom()), drcut(dr), rapscheme(scheme) { }
222  DeltaRInRange(const ParticleBase& vec, double drmin, double drmax, RapScheme scheme=PSEUDORAPIDITY)
223  : DeltaRInRange(vec, make_pair(drmin, drmax), scheme) { }
224  DeltaRInRange(const FourMomentum& vec, const pair<double,double>& dr, RapScheme scheme=PSEUDORAPIDITY)
225  : refvec(vec), drcut(dr), rapscheme(scheme) { }
226  DeltaRInRange(const FourMomentum& vec, double drmin, double drmax, RapScheme scheme=PSEUDORAPIDITY)
227  : DeltaRInRange(vec, make_pair(drmin, drmax), scheme) { }
228  DeltaRInRange(const Vector3& vec, const pair<double,double>& dr)
229  : drcut(dr), rapscheme(PSEUDORAPIDITY) { refvec.setPx(vec.x()); refvec.setPy(vec.y()); refvec.setPz(vec.z()); }
230  DeltaRInRange(const Vector3& vec, double drmin, double drmax)
231  : DeltaRInRange(vec, make_pair(drmin, drmax)) { }
232  bool operator()(const ParticleBase& p) const {
233  const double dR = deltaR(p, refvec, rapscheme);
234  return dR >= drcut.first && dR < drcut.second;
235  }
236  FourMomentum refvec;
237  pair<double,double> drcut;
238  RapScheme rapscheme;
239  };
241 
242 
245  DeltaPhiGtr(const ParticleBase& vec, double dphi)
246  : refvec(vec.p3()), dphicut(dphi) { }
247  DeltaPhiGtr(const FourMomentum& vec, double dphi)
248  : refvec(vec.p3()), dphicut(dphi) { }
249  DeltaPhiGtr(const Vector3& vec, double dphi)
250  : refvec(vec), dphicut(dphi) { }
251  bool operator()(const ParticleBase& p) const { return deltaPhi(p, refvec) > dphicut; }
252  Vector3 refvec;
253  double dphicut;
254  };
255  using deltaPhiGtr = DeltaPhiGtr;
256 
259  DeltaPhiLess(const ParticleBase& vec, double dphi)
260  : refvec(vec.p3()), dphicut(dphi) { }
261  DeltaPhiLess(const FourMomentum& vec, double dphi)
262  : refvec(vec.p3()), dphicut(dphi) { }
263  DeltaPhiLess(const Vector3& vec, double dphi)
264  : refvec(vec), dphicut(dphi) { }
265  bool operator()(const ParticleBase& p) const { return deltaPhi(p, refvec) < dphicut; }
266  Vector3 refvec;
267  double dphicut;
268  };
269  using deltaPhiLess = DeltaPhiLess;
270 
273  DeltaPhiInRange(const ParticleBase& vec, const pair<double,double>& dphi)
274  : refvec(vec.mom()), dphicut(dphi) { }
275  DeltaPhiInRange(const ParticleBase& vec, double dphimin, double dphimax)
276  : DeltaPhiInRange(vec, make_pair(dphimin, dphimax)) { }
277  DeltaPhiInRange(const FourMomentum& vec, const pair<double,double>& dphi)
278  : refvec(vec), dphicut(dphi) { }
279  DeltaPhiInRange(const FourMomentum& vec, double dphimin, double dphimax)
280  : DeltaPhiInRange(vec, make_pair(dphimin, dphimax)) { }
281  DeltaPhiInRange(const Vector3& vec, const pair<double,double>& dphi)
282  : refvec(vec), dphicut(dphi) { }
283  DeltaPhiInRange(const Vector3& vec, double dphimin, double dphimax)
284  : DeltaPhiInRange(vec, make_pair(dphimin, dphimax)) { }
285  bool operator()(const ParticleBase& p) const {
286  const double dphi = deltaPhi(p, refvec);
287  return dphi >= dphicut.first && dphi < dphicut.second;
288  }
289  Vector3 refvec;
290  pair<double,double> dphicut;
291  };
293 
294 
297  DeltaEtaGtr(const ParticleBase& vec, double deta)
298  : refvec(vec.p3()), detacut(deta) { }
299  DeltaEtaGtr(const FourMomentum& vec, double deta)
300  : refvec(vec.p3()), detacut(deta) { }
301  DeltaEtaGtr(const Vector3& vec, double deta)
302  : refvec(vec), detacut(deta) { }
303  bool operator()(const ParticleBase& p) const { return std::abs(deltaEta(p, refvec)) > detacut; }
304  Vector3 refvec;
305  double detacut;
306  };
307  using deltaEtaGtr = DeltaEtaGtr;
308 
311  DeltaEtaLess(const ParticleBase& vec, double deta)
312  : refvec(vec.p3()), detacut(deta) { }
313  DeltaEtaLess(const FourMomentum& vec, double deta)
314  : refvec(vec.p3()), detacut(deta) { }
315  DeltaEtaLess(const Vector3& vec, double deta)
316  : refvec(vec), detacut(deta) { }
317  bool operator()(const ParticleBase& p) const { return std::abs(deltaEta(p, refvec)) < detacut; }
318  Vector3 refvec;
319  double detacut;
320  };
321  using deltaEtaLess = DeltaEtaLess;
322 
325  DeltaEtaInRange(const ParticleBase& vec, const pair<double,double>& deta)
326  : refvec(vec.mom()), detacut(deta) { }
327  DeltaEtaInRange(const ParticleBase& vec, double detamin, double detamax)
328  : DeltaEtaInRange(vec, make_pair(detamin, detamax)) { }
329  DeltaEtaInRange(const FourMomentum& vec, const pair<double,double>& deta)
330  : refvec(vec), detacut(deta) { }
331  DeltaEtaInRange(const FourMomentum& vec, double detamin, double detamax)
332  : DeltaEtaInRange(vec, make_pair(detamin, detamax)) { }
333  DeltaEtaInRange(const Vector3& vec, const pair<double,double>& deta)
334  : refvec(vec), detacut(deta) { }
335  DeltaEtaInRange(const Vector3& vec, double detamin, double detamax)
336  : DeltaEtaInRange(vec, make_pair(detamin, detamax)) { }
337  bool operator()(const ParticleBase& p) const {
338  const double deta = deltaEta(p, refvec);
339  return deta >= detacut.first && deta < detacut.second;
340  }
341  Vector3 refvec;
342  pair<double,double> detacut;
343  };
345 
346 
349  DeltaRapGtr(const ParticleBase& vec, double drap)
350  : refvec(vec.mom()), drapcut(drap) { }
351  DeltaRapGtr(const FourMomentum& vec, double drap)
352  : refvec(vec), drapcut(drap) { }
353  bool operator()(const ParticleBase& p) const { return std::abs(deltaRap(p, refvec)) > drapcut; }
354  FourMomentum refvec;
355  double drapcut;
356  };
357  using deltaRapGtr = DeltaRapGtr;
358 
361  DeltaRapLess(const ParticleBase& vec, double drap)
362  : refvec(vec.mom()), drapcut(drap) { }
363  DeltaRapLess(const FourMomentum& vec, double drap)
364  : refvec(vec), drapcut(drap) { }
365  bool operator()(const ParticleBase& p) const { return std::abs(deltaRap(p, refvec)) < drapcut; }
366  FourMomentum refvec;
367  double drapcut;
368  };
369  using deltaRapLess = DeltaRapLess;
370 
373  DeltaRapInRange(const ParticleBase& vec, const pair<double,double>& drap)
374  : refvec(vec.mom()), drapcut(drap) { }
375  DeltaRapInRange(const ParticleBase& vec, double drapmin, double drapmax)
376  : DeltaRapInRange(vec, make_pair(drapmin, drapmax)) { }
377  DeltaRapInRange(const FourMomentum& vec, const pair<double,double>& drap)
378  : refvec(vec), drapcut(drap) { }
379  DeltaRapInRange(const FourMomentum& vec, double drapmin, double drapmax)
380  : DeltaRapInRange(vec, make_pair(drapmin, drapmax)) { }
381  bool operator()(const ParticleBase& p) const {
382  const double drap = deltaRap(p, refvec);
383  return drap >= drapcut.first && drap < drapcut.second;
384  }
385  FourMomentum refvec;
386  pair<double,double> drapcut;
387  };
389 
391 
392 
398 
401  virtual double operator()(const ParticleBase& p) const = 0;
402  virtual ~DoubleParticleBaseFunctor() {}
403  };
404 
407  DeltaRWRT(const ParticleBase& pb, RapScheme scheme=PSEUDORAPIDITY) : p(pb.mom()), rapscheme(scheme) {}
408  DeltaRWRT(const FourMomentum& p4, RapScheme scheme=PSEUDORAPIDITY) : p(p4), rapscheme(scheme) {}
409  DeltaRWRT(const Vector3& p3) : p(p3.mod(), p3.x(), p3.y(), p3.z()), rapscheme(PSEUDORAPIDITY) {}
410  double operator()(const ParticleBase& pb) const { return deltaR(p, pb, rapscheme); }
411  double operator()(const FourMomentum& p4) const { return deltaR(p, p4, rapscheme); }
412  double operator()(const Vector3& p3) const { return deltaR(p, p3); }
413  const FourMomentum p;
414  RapScheme rapscheme;
415  };
416  using deltaRWRT = DeltaRWRT;
417 
420  DeltaPhiWRT(const ParticleBase& pb) : p(pb.mom().vector3()) {}
421  DeltaPhiWRT(const FourMomentum& p4) : p(p4.vector3()) {}
422  DeltaPhiWRT(const Vector3& p3) : p(p3) {}
423  double operator()(const ParticleBase& pb) const { return deltaPhi(p, pb); }
424  double operator()(const FourMomentum& p4) const { return deltaPhi(p, p4); }
425  double operator()(const Vector3& p3) const { return deltaPhi(p, p3); }
426  const Vector3 p;
427  };
428  using deltaPhiWRT = DeltaPhiWRT;
429 
432  DeltaEtaWRT(const ParticleBase& pb) : p(pb.mom().vector3()) {}
433  DeltaEtaWRT(const FourMomentum& p4) : p(p4.vector3()) {}
434  DeltaEtaWRT(const Vector3& p3) : p(p3) {}
435  double operator()(const ParticleBase& pb) const { return deltaEta(p, pb); }
436  double operator()(const FourMomentum& p4) const { return deltaEta(p, p4); }
437  double operator()(const Vector3& p3) const { return deltaEta(p, p3); }
438  const Vector3 p;
439  };
440  using deltaEtaWRT = DeltaEtaWRT;
441 
444  AbsDeltaEtaWRT(const ParticleBase& pb) : p(pb.mom().vector3()) {}
445  AbsDeltaEtaWRT(const FourMomentum& p4) : p(p4.vector3()) {}
446  AbsDeltaEtaWRT(const Vector3& p3) : p(p3) {}
447  double operator()(const ParticleBase& pb) const { return fabs(deltaEta(p, pb)); }
448  double operator()(const FourMomentum& p4) const { return fabs(deltaEta(p, p4)); }
449  double operator()(const Vector3& p3) const { return fabs(deltaEta(p, p3)); }
450  const Vector3 p;
451  };
453 
456  DeltaRapWRT(const ParticleBase& pb) : p(pb.mom()) {}
457  DeltaRapWRT(const FourMomentum& p4) : p(p4) {}
458  double operator()(const ParticleBase& pb) const { return deltaRap(p, pb); }
459  double operator()(const FourMomentum& p4) const { return deltaRap(p, p4); }
460  const FourMomentum p;
461  };
462  using deltaRapWRT = DeltaRapWRT;
463 
466  AbsDeltaRapWRT(const ParticleBase& pb) : p(pb.mom()) {}
467  AbsDeltaRapWRT(const FourMomentum& p4) : p(p4) {}
468  double operator()(const ParticleBase& pb) const { return fabs(deltaRap(p, pb)); }
469  double operator()(const FourMomentum& p4) const { return fabs(deltaRap(p, p4)); }
470  const FourMomentum p;
471  };
473 
475 
476 
479 
480  template<typename PBCONTAINER1, typename PBCONTAINER2>
481  inline void idiscardIfAny(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
482  typename std::function<bool(const typename PBCONTAINER1::value_type&,
483  const typename PBCONTAINER2::value_type&)> fn) {
484  for (const auto& pbcmp : tocompare) {
485  ifilter_discard(tofilter, [&](const typename PBCONTAINER1::value_type& pbfilt){ return fn(pbfilt, pbcmp); });
486  }
487  }
488 
489  template<typename PBCONTAINER1, typename PBCONTAINER2>
490  inline PBCONTAINER1 discardIfAny(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
491  typename std::function<bool(const typename PBCONTAINER1::value_type&,
492  const typename PBCONTAINER2::value_type&)> fn) {
493  PBCONTAINER1 tmp{tofilter};
494  idiscardIfAny(tmp, tocompare, fn);
495  return tmp;
496  }
497 
498 
499  template<typename PBCONTAINER1, typename PBCONTAINER2>
500  inline PBCONTAINER1 selectIfAny(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
501  typename std::function<bool(const typename PBCONTAINER1::value_type&,
502  const typename PBCONTAINER2::value_type&)> fn) {
503  PBCONTAINER1 selected;
504  for (const auto& pbfilt : tofilter) {
505  if (any(tocompare, [&](const typename PBCONTAINER2::value_type& pbcmp){ return fn(pbfilt, pbcmp); })) {
506  selected += pbfilt;
507  }
508  }
509  return selected;
510  }
511 
512  template<typename PBCONTAINER1, typename PBCONTAINER2>
513  inline void iselectIfAny(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
514  typename std::function<bool(const typename PBCONTAINER1::value_type&,
515  const typename PBCONTAINER2::value_type&)> fn) {
516  tofilter = selectIfAny(tofilter, tocompare, fn);
517  }
518 
519 
520 
521  template<typename PBCONTAINER1, typename PBCONTAINER2>
522  inline PBCONTAINER1 discardIfAll(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
523  typename std::function<bool(const typename PBCONTAINER1::value_type&,
524  const typename PBCONTAINER2::value_type&)> fn) {
525  PBCONTAINER1 selected;
526  for (const auto& pbfilt : tofilter) {
527  if (!all(tocompare, [&](const typename PBCONTAINER2::value_type& pbcmp){ return fn(pbfilt, pbcmp); })) {
528  selected += pbfilt;
529  }
530  }
531  return selected;
532  }
533 
534  template<typename PBCONTAINER1, typename PBCONTAINER2>
535  inline void idiscardIfAll(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
536  typename std::function<bool(const typename PBCONTAINER1::value_type&,
537  const typename PBCONTAINER2::value_type&)> fn) {
538  tofilter = discardIfAll(tofilter, tocompare, fn);
539  }
540 
541 
542  template<typename PBCONTAINER1, typename PBCONTAINER2>
543  inline PBCONTAINER1 selectIfAll(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
544  typename std::function<bool(const typename PBCONTAINER1::value_type&,
545  const typename PBCONTAINER2::value_type&)> fn) {
546  PBCONTAINER1 selected;
547  for (const auto& pbfilt : tofilter) {
548  if (all(tocompare, [&](const typename PBCONTAINER2::value_type& pbcmp){ return fn(pbfilt, pbcmp); })) {
549  selected += pbfilt;
550  }
551  }
552  return selected;
553  }
554 
555  template<typename PBCONTAINER1, typename PBCONTAINER2>
556  inline void iselectIfAll(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare,
557  typename std::function<bool(const typename PBCONTAINER1::value_type&,
558  const typename PBCONTAINER2::value_type&)> fn) {
559  tofilter = selectIfAll(tofilter, tocompare, fn);
560  }
561 
563 
564 
567 
568  template<typename PBCONTAINER1, typename PBCONTAINER2>
569  inline void idiscardIfAnyDeltaRLess(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dR) {
570  for (const typename PBCONTAINER2::value_type& pb : tocompare) {
571  ifilter_discard(tofilter, deltaRLess(pb, dR));
572  }
573  }
574 
575  template<typename PBCONTAINER1, typename PBCONTAINER2>
576  inline PBCONTAINER1 discardIfAnyDeltaRLess(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dR) {
577  PBCONTAINER1 tmp{tofilter};
578  idiscardIfAnyDeltaRLess(tmp, tocompare, dR);
579  return tmp;
580  }
581 
582  template<typename PBCONTAINER1, typename PBCONTAINER2>
583  inline void idiscardIfAnyDeltaPhiLess(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dphi) {
584  for (const typename PBCONTAINER2::value_type& pb : tocompare) {
585  ifilter_discard(tofilter, deltaPhiLess(pb, dphi));
586  }
587  }
588 
589  template<typename PBCONTAINER1, typename PBCONTAINER2>
590  inline PBCONTAINER1 discardIfAnyDeltaPhiLess(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dphi) {
591  PBCONTAINER1 tmp{tofilter};
592  idiscardIfAnyDeltaPhiLess(tmp, tocompare, dphi);
593  return tmp;
594  }
595 
596 
597 
598  template<typename PBCONTAINER1, typename PBCONTAINER2>
599  inline PBCONTAINER1 selectIfAnyDeltaRLess(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dR) {
600  PBCONTAINER1 selected;
601  for (const typename PBCONTAINER1::value_type& f : tofilter) {
602  if (any(tocompare, deltaRLess(f, dR))) selected.push_back(f);
603  }
604  return selected;
605  }
606 
607  template<typename PBCONTAINER1, typename PBCONTAINER2>
608  inline void iselectIfAnyDeltaRLess(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dR) {
609  tofilter = selectIfAnyDeltaRLess(tofilter, tocompare, dR);
610  }
611 
612 
613  template<typename PBCONTAINER1, typename PBCONTAINER2>
614  inline PBCONTAINER1 selectIfAnyDeltaPhiLess(const PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dphi) {
615  PBCONTAINER1 selected;
616  for (const typename PBCONTAINER1::value_type& f : tofilter) {
617  if (any(tocompare, deltaPhiLess(f, dphi))) selected.push_back(f);
618  }
619  return selected;
620  }
621 
622  template<typename PBCONTAINER1, typename PBCONTAINER2>
623  inline void iselectIfAnyDeltaPhiLess(PBCONTAINER1& tofilter, const PBCONTAINER2& tocompare, double dphi) {
624  tofilter = selectIfAnyDeltaPhiLess(tofilter, tocompare, dphi);
625  }
626 
627 
629 
631 
632 
633 
639  namespace Kin {
640 
642  inline FourMomentum mom(const ParticleBase& p) { return p.mom(); }
644  inline FourMomentum p4(const ParticleBase& p) { return p.mom(); }
645 
647  inline Vector3 p3(const ParticleBase& p) { return p.p3(); }
648 
650  inline Vector3 pTvec(const ParticleBase& p) { return p.pTvec(); }
651 
653  inline double p(const ParticleBase& p) { return p.p(); }
654 
656  inline double pT(const ParticleBase& p) { return p.pT(); }
657 
659  inline double Et(const ParticleBase& p) { return p.Et(); }
660 
662  inline double eta(const ParticleBase& p) { return p.eta(); }
663 
665  inline double abseta(const ParticleBase& p) { return p.abseta(); }
666 
668  inline double rap(const ParticleBase& p) { return p.rap(); }
669 
671  inline double absrap(const ParticleBase& p) { return p.absrap(); }
672 
674  inline double mass(const ParticleBase& p) { return p.mass(); }
675 
676 
678  inline double pairPt(const ParticleBase& p1, const ParticleBase& p2) { return (p1.mom() + p2.mom()).pT(); }
679 
681  inline double pairMass(const ParticleBase& p1, const ParticleBase& p2) { return (p1.mom() + p2.mom()).mass(); }
682 
683  }
684 
685  // Import Kin namespace into Rivet
686  using namespace Kin;
687 
689 
690 
692 
693 }
694 
695 #endif
Definition: MC_Cent_pPb.hh:10
double abseta() const
Get the directly (alias).
Definition: Vector4.hh:159
double pT() const
Calculate the transverse momentum .
Definition: Vector4.hh:628
bool any(const CONTAINER &c)
Return true if x is true for any x in container c, otherwise false.
Definition: Utils.hh:323
Abs pseudorapidity greater-than functor.
Definition: ParticleBaseUtils.hh:93
double eta() const
Synonym for pseudorapidity.
Definition: Vector4.hh:152
Pseudorapidity less-than functor.
Definition: ParticleBaseUtils.hh:73
Base type for Particle -> double functors.
Definition: ParticleBaseUtils.hh:400
double absrap() const
Absolute rapidity.
Definition: Vector4.hh:605
(with respect to another momentum, vec) greater-than functor
Definition: ParticleBaseUtils.hh:348
(with respect to another momentum, vec) greater-than functor
Definition: ParticleBaseUtils.hh:296
function< bool(const ParticleBase &, const ParticleBase &)> ParticleBaseSorter
std::function instantiation for functors taking two ParticleBase and returning a bool ...
Definition: ParticleBaseUtils.hh:21
Abs pseudorapidity momentum less-than functor.
Definition: ParticleBaseUtils.hh:103
const FourMomentum & mom() const
Get equivalent single momentum four-vector (const) (alias).
Definition: ParticleBase.hh:39
double abseta() const
Get the directly (alias).
Definition: ParticleBase.hh:91
Vector3 p3() const
Get 3-momentum part, .
Definition: Vector4.hh:578
Pseudorapidity in-range functor.
Definition: ParticleBaseUtils.hh:82
Abs rapidity in-range functor.
Definition: ParticleBaseUtils.hh:174
function< bool(const ParticleBase &)> ParticleBaseSelector
std::function instantiation for functors taking a ParticleBase and returning a bool ...
Definition: ParticleBaseUtils.hh:19
Base type for Particle -> bool functors.
Definition: ParticleBaseUtils.hh:25
double mass() const
Get the mass directly.
Definition: ParticleBase.hh:80
(with respect to another 4-momentum, vec) greater-than functor
Definition: ParticleBaseUtils.hh:189
double deltaEta(double eta1, double eta2, bool sign=false)
Definition: MathUtils.hh:637
(with respect to another 4-momentum, vec) in-range functor
Definition: ParticleBaseUtils.hh:272
Calculator of with respect to a given momentum.
Definition: ParticleBaseUtils.hh:419
Transverse momentum in-range functor.
Definition: ParticleBaseUtils.hh:52
double Et() const
Get the directly.
Definition: ParticleBase.hh:75
(with respect to another momentum, vec) less-than functor
Definition: ParticleBaseUtils.hh:360
double deltaPhi(double phi1, double phi2, bool sign=false)
Calculate the difference between two angles in radians.
Definition: MathUtils.hh:629
double rap() const
Get the directly (alias).
Definition: ParticleBase.hh:96
Transverse momentum greater-than functor.
Definition: ParticleBaseUtils.hh:32
Calculator of with respect to a given momentum.
Definition: ParticleBaseUtils.hh:431
Rapidity in-range functor.
Definition: ParticleBaseUtils.hh:143
Vector3 pTvec() const
Get the transverse 3-momentum directly.
Definition: ParticleBase.hh:117
Pseudorapidity greater-than functor.
Definition: ParticleBaseUtils.hh:64
Abs rapidity greater-than functor.
Definition: ParticleBaseUtils.hh:154
Vector3 p3() const
Get the 3-momentum directly.
Definition: ParticleBase.hh:108
double pT() const
Get the directly (alias).
Definition: ParticleBase.hh:63
Calculator of with respect to a given momentum.
Definition: ParticleBaseUtils.hh:465
Calculator of with respect to a given momentum.
Definition: ParticleBaseUtils.hh:406
(with respect to another 4-momentum, vec) less-than functor
Definition: ParticleBaseUtils.hh:204
Base class for particle-like things like Particle and Jet.
Definition: ParticleBase.hh:13
double deltaRap(double y1, double y2, bool sign=false)
Definition: MathUtils.hh:645
(with respect to another momentum, vec) less-than functor
Definition: ParticleBaseUtils.hh:310
double deltaR(double rap1, double phi1, double rap2, double phi2)
Definition: MathUtils.hh:659
Abs rapidity momentum less-than functor.
Definition: ParticleBaseUtils.hh:164
double absrap() const
Get the directly (alias).
Definition: ParticleBase.hh:100
RapScheme
Enum for rapidity variable to be used in calculating , applying rapidity cuts, etc.
Definition: MathConstants.hh:46
Vector3 vector3() const
Get the spatial part of the 4-vector as a 3-vector.
Definition: Vector4.hh:162
Abs pseudorapidity in-range functor.
Definition: ParticleBaseUtils.hh:113
Calculator of with respect to a given momentum.
Definition: ParticleBaseUtils.hh:455
Transverse momentum less-than functor.
Definition: ParticleBaseUtils.hh:42
double eta() const
Get the directly (alias).
Definition: ParticleBase.hh:87
Rapidity greater-than functor.
Definition: ParticleBaseUtils.hh:125
(with respect to another 4-momentum, vec) in-range functor
Definition: ParticleBaseUtils.hh:219
Three-dimensional specialisation of Vector.
Definition: Vector3.hh:26
(with respect to another momentum, vec) less-than functor
Definition: ParticleBaseUtils.hh:258
bool all(const CONTAINER &c)
Return true if x is true for all x in container c, otherwise false.
Definition: Utils.hh:345
(with respect to another momentum, vec) greater-than functor
Definition: ParticleBaseUtils.hh:244
Calculator of with respect to a given momentum.
Definition: ParticleBaseUtils.hh:443
Specialized version of the FourVector with momentum/energy functionality.
Definition: Vector4.hh:301
Jets & ifilter_discard(Jets &jets, const Cut &c)
Filter a jet collection in-place to the subset that fails the supplied Cut.
double p() const
Get the 3-momentum magnitude directly.
Definition: ParticleBase.hh:110
double rap() const
Alias for rapidity.
Definition: Vector4.hh:596
double mod() const
Calculate the modulus of a vector. .
Definition: VectorN.hh:95
(with respect to another 4-momentum, vec) in-range functor
Definition: ParticleBaseUtils.hh:324
(with respect to another 4-momentum, vec) in-range functor
Definition: ParticleBaseUtils.hh:372
Rapidity momentum less-than functor.
Definition: ParticleBaseUtils.hh:134