gambit is hosted by Hepforge, IPPP Durham
GAMBIT  v1.5.0-252-gf9a3f78
a Global And Modular Bsm Inference Tool
Utils.hpp
Go to the documentation of this file.
1 #pragma once
2 #include "HEPUtils/MathUtils.h"
3 #include "HEPUtils/BinnedFn.h"
4 #include "HEPUtils/Event.h"
5 #include "HEPUtils/FastJet.h"
6 #include <functional>
7 #include <memory>
8 
9 namespace Gambit
10 {
11 
12  namespace ColliderBit
13  {
14 
15 
17  static const double GeV = 1, MeV = 1e-3, TeV = 1e3;
18 
19 
21  using HEPUtils::Event;
23  using HEPUtils::Particle;
25  using HEPUtils::Jet;
27  using HEPUtils::P4;
28 
29 
31  typedef std::vector<HEPUtils::Particle*> ParticlePtrs;
33  typedef std::vector<const HEPUtils::Particle*> ConstParticlePtrs;
34 
35 
37  typedef std::vector<HEPUtils::Jet*> JetPtrs;
39  typedef std::vector<const HEPUtils::Jet*> ConstJetPtrs;
40 
41 
43 
44 
46  template <typename CONTAINER, typename RMFN>
47  inline void iremoveerase(CONTAINER& c, const RMFN& fn) {
48  auto newend = std::remove_if(c.begin(), c.end(), fn);
49  c.erase(newend, c.end());
50  }
51 
53  inline void ifilter_reject(ParticlePtrs& particles,
54  std::function<bool(Particle*)> rejfn, bool do_delete=true) {
55  iremoveerase(particles, [&](Particle* p) {
56  const bool rm = rejfn(p);
57  if (rm && do_delete) delete p;
58  return rm;
59  });
60  }
61 
63  inline void ifilter_select(ParticlePtrs& particles,
64  std::function<bool(Particle*)> selfn, bool do_delete=true) {
65  ifilter_reject(particles, [&](Particle* p) { return !selfn(p); }, do_delete);
66  }
67 
68 
71  inline ParticlePtrs filter_reject(const ParticlePtrs& particles,
72  std::function<bool(Particle*)> rejfn, bool do_delete=true) {
73  ParticlePtrs rtn = particles;
74  ifilter_reject(rtn, rejfn, do_delete);
75  return rtn;
76  }
77 
79  inline ParticlePtrs filter_select(const ParticlePtrs& particles,
80  std::function<bool(Particle*)> selfn, bool do_delete=true) {
81  return filter_reject(particles, [&](Particle* p) { return !selfn(p); }, do_delete);
82  }
83 
85 
86 
87 
89 
90 
92  inline void ifilter_reject(JetPtrs& jets,
93  std::function<bool(Jet*)> rejfn, bool do_delete=true) {
94  iremoveerase(jets, [&](Jet* j) {
95  const bool rm = rejfn(j);
96  if (rm && do_delete) delete j;
97  return rm;
98  });
99  }
100 
102  inline void ifilter_select(JetPtrs& jets,
103  std::function<bool(Jet*)> selfn, bool do_delete=true) {
104  ifilter_reject(jets, [&](Jet* j) { return !selfn(j); }, do_delete);
105  }
106 
107 
110  inline JetPtrs filter_reject(const JetPtrs& jets,
111  std::function<bool(Jet*)> rejfn, bool do_delete=true) {
112  JetPtrs rtn = jets;
113  ifilter_reject(rtn, rejfn, do_delete);
114  return rtn;
115  }
116 
118  inline JetPtrs filter_select(const JetPtrs& jets,
119  std::function<bool(Jet*)> selfn, bool do_delete=true) {
120  return filter_reject(jets, [&](Jet* j) { return !selfn(j); }, do_delete);
121  }
122 
124 
125 
126 
128 
129 
130 
132 
133 
135  bool random_bool(double eff);
136 
138  inline bool random_bool(const HEPUtils::BinnedFn1D<double>& effmap, double x) {
139  return random_bool( effmap.get_at(x) );
140  }
141 
143  inline bool random_bool(const HEPUtils::BinnedFn2D<double>& effmap, double x, double y) {
144  return random_bool( effmap.get_at(x, y) );
145  }
146 
148 
149 
151 
152 
154  void filtereff(std::vector<HEPUtils::Particle*>& particles, double eff, bool do_delete=false);
155 
157  void filtereff(std::vector<HEPUtils::Particle*>& particles, std::function<double(HEPUtils::Particle*)> eff_fn, bool do_delete=false);
158 
160  void filtereff_pt(std::vector<HEPUtils::Particle*>& particles, const HEPUtils::BinnedFn1D<double>& eff_pt, bool do_delete=false);
161 
163  void filtereff_etapt(std::vector<HEPUtils::Particle*>& particles, const HEPUtils::BinnedFn2D<double>& eff_etapt, bool do_delete=false);
164 
166 
167 
169 
170 
173  inline bool has_tag(const HEPUtils::BinnedFn2D<double>& effmap, double eta, double pt) {
174  try {
175  return random_bool(effmap, fabs(eta), pt);
176  } catch (...) {
177  return false; // No tag if outside lookup range... be careful!
178  }
179  }
180 
181  template <typename NUM1, typename NUM2>
182  inline size_t binIndex(NUM1 val, const std::vector<NUM2>& binedges, bool allow_overflow=false) {
183  if (val < binedges.front()) return -1;
184  if (val >= binedges.back()) return allow_overflow ? int(binedges.size())-1 : -1;
185  return std::distance(binedges.begin(), --std::upper_bound(binedges.begin(), binedges.end(), val));
186  }
187 
188 
190  inline std::vector<double> mk_bin_values(const std::vector<double>& binEdgeValues) {
191  std::vector<double> results;
192  results.reserve(binEdgeValues.size()-1);
193  for (size_t i = 0; i < binEdgeValues.size()-1; ++i)
194  results.push_back( (binEdgeValues[i] + binEdgeValues[i+1])/2.0 );
195  return results;
196  }
198  inline std::vector<double> makeBinValues(const std::vector<double>& binEdgeValues) {
199  return mk_bin_values(binEdgeValues);
200  }
201 
202 
204  template <typename MOM>
205  inline std::vector<std::shared_ptr<HEPUtils::Jet>> get_jets(const std::vector<MOM*>& moms, double R,
206  double ptmin=0*GeV, FJNS::JetAlgorithm alg=FJNS::antikt_algorithm) {
207  // Make PseudoJets
208  std::vector<FJNS::PseudoJet> constituents;
209  for (const MOM* p : moms) constituents.push_back(HEPUtils::mk_pseudojet(*p));
210  // Run clustering
211  std::vector<FJNS::PseudoJet> jets = HEPUtils::get_jets(constituents, R, ptmin, alg);
212  // Make newly-allocated Jets
213  std::vector<std::shared_ptr<HEPUtils::Jet>> rtn;
214  for (const FJNS::PseudoJet& j : jets) rtn.push_back(std::make_shared<HEPUtils::Jet>(HEPUtils::mk_p4(j)));
215  return rtn;
216  }
217 
218 
220  inline bool object_in_cone(const HEPUtils::Event& e, const HEPUtils::P4& p4, double ptmin, double rmax, double rmin=0.05) {
221  for (const HEPUtils::Particle* p : e.visible_particles())
222  if (p->pT() > ptmin && HEPUtils::in_range(HEPUtils::deltaR_eta(p4, *p), rmin, rmax)) return true;
223  for (const HEPUtils::Jet* j : e.jets())
224  if (j->pT() > ptmin && HEPUtils::in_range(HEPUtils::deltaR_eta(p4, *j), rmin, rmax)) return true;
225  return false;
226  }
227 
228 
230  template<typename MOMPTRS1, typename MOMPTRS2>
231  void removeOverlap(MOMPTRS1& momstofilter, const MOMPTRS2& momstocmp, double deltaRMax, bool use_rapidity=false) {
232  ifilter_reject(momstofilter, [&](const typename MOMPTRS1::value_type& mom1) {
233  for (const typename MOMPTRS2::value_type& mom2 : momstocmp) {
234  const double dR = (use_rapidity) ? deltaR_rap(mom1->mom(), mom2->mom()) : deltaR_eta(mom1->mom(), mom2->mom());
235  if (dR < deltaRMax) return true;
236  }
237  return false;
238  }, false);
239  }
240 
241 
243  template <typename CONTAINER, typename FN>
244  inline bool all_of(const CONTAINER& c, const FN& f) {
245  return std::all_of(std::begin(c), std::end(c), f);
246  }
247 
249  template <typename CONTAINER, typename FN>
250  inline bool any_of(const CONTAINER& c, const FN& f) {
251  return std::any_of(std::begin(c), std::end(c), f);
252  }
253 
255  template <typename CONTAINER, typename FN>
256  inline bool none_of(const CONTAINER& c, const FN& f) {
257  return std::none_of(std::begin(c), std::end(c), f);
258  }
259 
260 
262  std::vector<std::vector<HEPUtils::Particle*>> getSFOSpairs(std::vector<HEPUtils::Particle*> particles);
263 
265  std::vector<std::vector<HEPUtils::Particle*>> getOSpairs(std::vector<HEPUtils::Particle*> particles);
266 
268  std::vector<std::vector<HEPUtils::Particle*>> getSSpairs(std::vector<HEPUtils::Particle*> particles);
269 
270 
272 
273 
275  inline void sortBy(ParticlePtrs& particles, std::function<bool(const Particle*, const Particle*)> cmpfn) {
276  std::sort(particles.begin(), particles.end(), cmpfn);
277  }
278 
280  inline bool cmpParticlesByPt(const HEPUtils::Particle* lep1, const HEPUtils::Particle* lep2) { return lep1->pT() > lep2->pT(); }
281 
282  // Sort a particles list by decreasing pT
283  inline void sortByPt(ParticlePtrs& particles) { sortBy(particles, cmpParticlesByPt); }
284 
285 
287  inline void sortBy(JetPtrs& jets, std::function<bool(const Jet*, const Jet*)> cmpfn) {
288  std::sort(jets.begin(), jets.end(), cmpfn);
289  }
290 
292  inline bool cmpJetsByPt(const HEPUtils::Jet* jet1, const HEPUtils::Jet* jet2) { return jet1->pT() > jet2->pT(); }
293 
294  // Sort a jets list by decreasing pT
295  inline void sortByPt(JetPtrs& jets) { sortBy(jets, cmpJetsByPt); }
296 
298 
299 
300  }
301 }
ParticlePtrs filter_select(const ParticlePtrs &particles, std::function< bool(Particle *)> selfn, bool do_delete=true)
Filter a supplied particle vector by keeping those which pass a supplied cut.
Definition: Utils.hpp:79
bool has_tag(const HEPUtils::BinnedFn2D< double > &effmap, double eta, double pt)
Randomly get a tag result (can be anything) from a 2D |eta|-pT efficiency map.
Definition: Utils.hpp:173
std::vector< const HEPUtils::Jet * > ConstJetPtrs
Typedef for a vector of const Jet pointers.
Definition: Utils.hpp:39
bool object_in_cone(const HEPUtils::Event &e, const HEPUtils::P4 &p4, double ptmin, double rmax, double rmin=0.05)
Check if there&#39;s a physics object above ptmin in an annulus rmin..rmax around the given four-momentum...
Definition: Utils.hpp:220
void filtereff(std::vector< HEPUtils::Particle *> &particles, double eff, bool do_delete=false)
Utility function for filtering a supplied particle vector by sampling wrt an efficiency scalar...
Definition: Utils.cpp:16
std::vector< HEPUtils::Particle * > ParticlePtrs
Typedef for a vector of Particle pointers.
Definition: Utils.hpp:31
FJNS::PseudoJet mk_pseudojet(const Vec4T &p)
Definition: Py8Utils.hpp:35
std::vector< double > makeBinValues(const std::vector< double > &binEdgeValues)
Alias.
Definition: Utils.hpp:198
HEPUtils::P4 mk_p4(const Vec4T &p)
Definition: Py8Utils.hpp:41
bool none_of(const CONTAINER &c, const FN &f)
Non-iterator version of std::none_of.
Definition: Utils.hpp:256
std::vector< std::vector< HEPUtils::Particle * > > getSSpairs(std::vector< HEPUtils::Particle *> particles)
Utility function for returning a collection of same-sign particle pairs.
Definition: Utils.cpp:103
std::vector< double > mk_bin_values(const std::vector< double > &binEdgeValues)
Make a vector of central bin values from a vector of bin edge values using linear interpolation...
Definition: Utils.hpp:190
bool cmpJetsByPt(const HEPUtils::Jet *jet1, const HEPUtils::Jet *jet2)
Comparison function to give a jets sorting order decreasing by pT.
Definition: Utils.hpp:292
void filtereff_pt(std::vector< HEPUtils::Particle *> &particles, const HEPUtils::BinnedFn1D< double > &eff_pt, bool do_delete=false)
Utility function for filtering a supplied particle vector by sampling wrt a binned 1D efficiency map ...
Definition: Utils.cpp:43
void ifilter_reject(ParticlePtrs &particles, std::function< bool(Particle *)> rejfn, bool do_delete=true)
In-place filter a supplied particle vector by rejecting those which fail a supplied cut...
Definition: Utils.hpp:53
size_t binIndex(NUM1 val, const std::vector< NUM2 > &binedges, bool allow_overflow=false)
Definition: Utils.hpp:182
void sortByPt(ParticlePtrs &particles)
Definition: Utils.hpp:283
std::vector< const HEPUtils::Particle * > ConstParticlePtrs
Typedef for a vector of const Particle pointers.
Definition: Utils.hpp:33
void ifilter_select(ParticlePtrs &particles, std::function< bool(Particle *)> selfn, bool do_delete=true)
In-place filter a supplied particle vector by keeping those which pass a supplied cut...
Definition: Utils.hpp:63
std::vector< std::vector< HEPUtils::Particle * > > getOSpairs(std::vector< HEPUtils::Particle *> particles)
Utility function for returning a collection of oppsosite-sign particle pairs.
Definition: Utils.cpp:86
void sortBy(ParticlePtrs &particles, std::function< bool(const Particle *, const Particle *)> cmpfn)
Particle-sorting function.
Definition: Utils.hpp:275
bool cmpParticlesByPt(const HEPUtils::Particle *lep1, const HEPUtils::Particle *lep2)
Comparison function to give a particles sorting order decreasing by pT.
Definition: Utils.hpp:280
void filtereff_etapt(std::vector< HEPUtils::Particle *> &particles, const HEPUtils::BinnedFn2D< double > &eff_etapt, bool do_delete=false)
Utility function for filtering a supplied particle vector by sampling wrt a binned 2D efficiency map ...
Definition: Utils.cpp:56
std::vector< std::shared_ptr< HEPUtils::Jet > > get_jets(const std::vector< MOM *> &moms, double R, double ptmin=0 *GeV, FJNS::JetAlgorithm alg=FJNS::antikt_algorithm)
Run jet clustering from any P4-compatible momentum type.
Definition: Utils.hpp:205
bool any_of(const CONTAINER &c, const FN &f)
Non-iterator version of std::any_of.
Definition: Utils.hpp:250
std::vector< std::vector< HEPUtils::Particle * > > getSFOSpairs(std::vector< HEPUtils::Particle *> particles)
Utility function for returning a collection of same-flavour, oppsosite-sign particle pairs...
Definition: Utils.cpp:69
void removeOverlap(MOMPTRS1 &momstofilter, const MOMPTRS2 &momstocmp, double deltaRMax, bool use_rapidity=false)
Overlap removal – discard from first list if within deltaRMax of any from the second list...
Definition: Utils.hpp:231
ParticlePtrs filter_reject(const ParticlePtrs &particles, std::function< bool(Particle *)> rejfn, bool do_delete=true)
Filter a supplied particle vector by rejecting those which fail a supplied cut.
Definition: Utils.hpp:71
bool all_of(const CONTAINER &c, const FN &f)
Non-iterator version of std::all_of.
Definition: Utils.hpp:244
void iremoveerase(CONTAINER &c, const RMFN &fn)
Convenience combination of remove_if and erase.
Definition: Utils.hpp:47
TODO: see if we can use this one:
Definition: Analysis.hpp:33
std::vector< HEPUtils::Jet * > JetPtrs
Typedef for a vector of Jet pointers.
Definition: Utils.hpp:37
bool random_bool(double eff)
Return a random true/false at a success rate given by a number.
Definition: Utils.cpp:10