Utils.hpp
Go to the documentation of this file.
188 void filtereff(std::vector<const HEPUtils::Particle*>& particles, double eff, bool do_delete=false); 191 void filtereff(std::vector<const HEPUtils::Particle*>& particles, std::function<double(const HEPUtils::Particle*)> eff_fn, bool do_delete=false); 194 void filtereff_pt(std::vector<const HEPUtils::Particle*>& particles, const HEPUtils::BinnedFn1D<double>& eff_pt, bool do_delete=false); 197 void filtereff_etapt(std::vector<const HEPUtils::Particle*>& particles, const HEPUtils::BinnedFn2D<double>& eff_etapt, bool do_delete=false); 216 inline std::map<const HEPUtils::Jet*,bool> generateBTagsMap(const std::vector<const HEPUtils::Jet*>& jets, 246 inline size_t binIndex(NUM1 val, const std::vector<NUM2>& binedges, bool allow_overflow=false) { 249 return std::distance(binedges.begin(), --std::upper_bound(binedges.begin(), binedges.end(), val)); 269 inline std::vector<std::shared_ptr<HEPUtils::Jet>> get_jets(const std::vector<MOM*>& moms, double R, 278 for (const FJNS::PseudoJet& j : jets) rtn.push_back(std::make_shared<HEPUtils::Jet>(HEPUtils::mk_p4(j))); 284 inline bool object_in_cone(const HEPUtils::Event& e, const HEPUtils::P4& p4, double ptmin, double rmax, double rmin=0.05) { 286 if (p->pT() > ptmin && HEPUtils::in_range(HEPUtils::deltaR_eta(p4, *p), rmin, rmax)) return true; 288 if (j->pT() > ptmin && HEPUtils::in_range(HEPUtils::deltaR_eta(p4, *j), rmin, rmax)) return true; 299 void removeOverlap(MOMPTRS1& momstofilter, const MOMPTRS2& momstocmp, double deltaRMax, bool use_rapidity=false, double pTmax = DBL_MAX, double btageff = 0) 304 const double dR = (use_rapidity) ? deltaR_rap(mom1->mom(), mom2->mom()) : deltaR_eta(mom1->mom(), mom2->mom()); 305 if (dR < deltaRMax && mom1->pT() < pTmax && ( !amIaBJet(mom1) || !random_bool(btageff) ) ) return true; 318 void removeOverlap(MOMPTRS1& momstofilter, const MOMPTRS2& momstocmp, double (*deltaRMax)(const double), bool use_rapidity=false, double pTmax = DBL_MAX, double btageff = 0) 323 const double dR = (use_rapidity) ? deltaR_rap(mom1->mom(), mom2->mom()) : deltaR_eta(mom1->mom(), mom2->mom()); 324 if (dR < deltaRMax(mom1->pT()) && mom1->pT() < pTmax && ( !amIaBJet(mom1) || !random_bool(btageff) ) ) return true; 335 void removeOverlapIfBjet(MOMPTRS1& momstofilter, std::vector<const HEPUtils::Jet*>& jets, double deltaRMax, bool use_rapidity=false, double pTmax = DBL_MAX) 340 const double dR = (use_rapidity) ? deltaR_rap(mom1->mom(), jet->mom()) : deltaR_eta(mom1->mom(), jet->mom()); 368 std::vector<std::vector<const HEPUtils::Particle*>> getSFOSpairs(std::vector<const HEPUtils::Particle*> particles); 371 std::vector<std::vector<const HEPUtils::Particle*>> getOSpairs(std::vector<const HEPUtils::Particle*> particles); 374 std::vector<std::vector<const HEPUtils::Particle*>> getSSpairs(std::vector<const HEPUtils::Particle*> particles); 381 inline void sortBy(ParticlePtrs& particles, std::function<bool(const Particle*, const Particle*)> cmpfn) { 386 inline bool cmpParticlesByPt(const HEPUtils::Particle* lep1, const HEPUtils::Particle* lep2) { return lep1->pT() > lep2->pT(); } 398 inline bool cmpJetsByPt(const HEPUtils::Jet* jet1, const HEPUtils::Jet* jet2) { return jet1->pT() > jet2->pT(); } void filtereff_etapt(std::vector< const 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 void removeOverlapIfBjet(MOMPTRS1 &momstofilter, std::vector< const HEPUtils::Jet *> &jets, double deltaRMax, bool use_rapidity=false, double pTmax=DBL_MAX) Overlap removal for checking against b-jets – discard from first list if within deltaRMax of a b-jet... Definition: Utils.hpp:335 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:207 bool object_in_cone(const HEPUtils::Event &e, const HEPUtils::P4 &p4, double ptmin, double rmax, double rmin=0.05) Check if there's a physics object above ptmin in an annulus rmin..rmax around the given four-momentum... Definition: Utils.hpp:284 void removeOverlap(MOMPTRS1 &momstofilter, const MOMPTRS2 &momstocmp, double deltaRMax, bool use_rapidity=false, double pTmax=DBL_MAX, double btageff=0) Overlap removal – discard from first list if within deltaRMax of any from the second list Optional a... Definition: Utils.hpp:299 double scalarSumPt(const std::vector< const Particle *> &particles, double pTlim=0.) Scalar sum pT of particles with pT > pTlim (default pTlim = 0) Definition: Utils.hpp:437 std::vector< const HEPUtils::Jet * > JetPtrs Typedef for a vector of Jet pointers. Definition: Utils.hpp:57 ParticlePtrs filter_reject(const ParticlePtrs &particles, std::function< bool(const Particle *)> rejfn, bool do_delete=true) Filter a supplied particle vector by rejecting those which fail a supplied cut. Definition: Utils.hpp:105 ParticlePtrs filter_select(const ParticlePtrs &particles, std::function< bool(const Particle *)> selfn, bool do_delete=true) Filter a supplied particle vector by keeping those which pass a supplied cut. Definition: Utils.hpp:113 FJNS::PseudoJet mk_pseudojet(const Vec4T &p) Definition: Py8Utils.hpp:39 std::vector< double > makeBinValues(const std::vector< double > &binEdgeValues) Alias. Definition: Utils.hpp:262 bool none_of(const CONTAINER &c, const FN &f) Non-iterator version of std::none_of. Definition: Utils.hpp:362 std::vector< std::vector< const HEPUtils::Particle * > > getOSpairs(std::vector< const HEPUtils::Particle *> particles) Utility function for returning a collection of oppsosite-sign particle pairs. Definition: Utils.cpp:86 std::vector< const HEPUtils::Particle * > ParticlePtrs Typedef for a vector of Particle pointers. Definition: Utils.hpp:54 std::map< const HEPUtils::Jet *, bool > generateBTagsMap(const std::vector< const HEPUtils::Jet *> &jets, double bTagEff, double cMissTagEff, double otherMissTagEff, double pTmin=0., double absEtaMax=DBL_MAX) Return a map<Jet*,bool> containing a generated b-tag for every jet in the input vector. Definition: Utils.hpp:216 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:254 bool cmpJetsByPt(const HEPUtils::Jet *jet1, const HEPUtils::Jet *jet2) Comparison function to give a jets sorting order decreasing by pT. Definition: Utils.hpp:398 int countPt(const std::vector< const Particle *> &particles, double pTlim) Count number of particles that have pT > pTlim. Definition: Utils.hpp:409 size_t binIndex(NUM1 val, const std::vector< NUM2 > &binedges, bool allow_overflow=false) Definition: Utils.hpp:246 void ifilter_reject(ParticlePtrs &particles, std::function< bool(const Particle *)> rejfn, bool do_delete=true) In-place filter a supplied particle vector by rejecting those which fail a supplied cut... Definition: Utils.hpp:87 std::vector< std::vector< const HEPUtils::Particle * > > getSFOSpairs(std::vector< const HEPUtils::Particle *> particles) Utility function for returning a collection of same-flavour, oppsosite-sign particle pairs... Definition: Utils.cpp:69 void sortBy(ParticlePtrs &particles, std::function< bool(const Particle *, const Particle *)> cmpfn) Particle-sorting function. Definition: Utils.hpp:381 bool cmpParticlesByPt(const HEPUtils::Particle *lep1, const HEPUtils::Particle *lep2) Comparison function to give a particles sorting order decreasing by pT. Definition: Utils.hpp:386 void filtereff_pt(std::vector< const 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 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:269 bool any_of(const CONTAINER &c, const FN &f) Non-iterator version of std::any_of. Definition: Utils.hpp:356 bool all_of(const CONTAINER &c, const FN &f) Non-iterator version of std::all_of. Definition: Utils.hpp:350 void iremoveerase(CONTAINER &c, const RMFN &fn) Convenience combination of remove_if and erase. Definition: Utils.hpp:81 std::vector< std::vector< const HEPUtils::Particle * > > getSSpairs(std::vector< const HEPUtils::Particle *> particles) Utility function for returning a collection of same-sign particle pairs. Definition: Utils.cpp:103 void ifilter_select(ParticlePtrs &particles, std::function< bool(const Particle *)> selfn, bool do_delete=true) In-place filter a supplied particle vector by keeping those which pass a supplied cut... Definition: Utils.hpp:97 void filtereff(std::vector< const 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 bool random_bool(double eff) Return a random true/false at a success rate given by a number. Definition: Utils.cpp:10 |