ATLAS Offline Software
AnalysisMisc.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #ifndef ANALYSISUTILS_ANALYSISMISC_H
6 #define ANALYSISUTILS_ANALYSISMISC_H
7 
14 #include <math.h>
15 #include <vector>
16 #include <algorithm>
17 
18 #include "CxxUtils/fpcompare.h"
19 #include "AthContainers/DataVector.h"
20 #include "TruthUtils/MagicNumbers.h"
21 #include <type_traits>
22 
23 
24 #ifndef XAOD_STANDALONE
25 #include "CLHEP/Vector/LorentzVector.h"
26 #include "EventKernel/INavigable4Momentum.h"
27 #include "McParticleEvent/TruthParticleContainer.h"
28 #endif
29 
30 namespace AnalysisUtils {
31 
32 #ifndef XAOD_STANDALONE
33 
36  namespace Delta {
37 
40  inline double phi (const INavigable4Momentum *p1, const INavigable4Momentum *p2) {
41  double phi1 = (p1->phi()>M_PI) ? p1->phi()-2*M_PI : p1->phi();
42  double phi2 = (p2->phi()>M_PI) ? p2->phi()-2*M_PI : p2->phi();
43  double dphi = fabs(phi1-phi2);
44  if(dphi>M_PI) dphi = 2*M_PI-dphi;
45  return dphi;
46  }
47 
50  inline double R (const INavigable4Momentum *p1, const double v_eta, const double v_phi) {
51  double phi1 = (p1->phi()>M_PI) ? p1->phi()-2*M_PI : p1->phi();
52  double phi2 = (v_phi>M_PI) ? v_phi-2*M_PI : v_phi;
53  double dphi = fabs(phi1-phi2);
54  if(dphi>M_PI) dphi = 2*M_PI-dphi;
55  double deta = p1->eta() - v_eta;
56  return sqrt(dphi*dphi+deta*deta);
57  }
58 
61  inline double R (const INavigable4Momentum *p1, const INavigable4Momentum *p2) {
62  return R (p1, p2->eta(), p2->phi());
63  }
64  } // end of Delta namespace
65 
67 
70  namespace Imass {
71 
72  inline double two (const INavigable4Momentum *p1, const INavigable4Momentum *p2) {
73  return (p1->hlv()+p2->hlv()).m();
74  }
75 
76  inline double four (const INavigable4Momentum *p1, const INavigable4Momentum *p2,
77  const INavigable4Momentum *p3, const INavigable4Momentum *p4) {
78  return (p1->hlv()+p2->hlv()+p3->hlv()+p4->hlv()).m();
79  }
80  } // end of Imass namespace
81 
83 
86  namespace Match {
87 
93  template <class COLL>
94  inline bool R (const double eta, const double phi, COLL *coll, int &index, double &deltaR, const int pdg)
95  {
96  deltaR = 10000.; // big value
97  bool l_return = false;
98  int l_idx = 0;
99  typename COLL::const_iterator it = coll->begin();
100  typename COLL::const_iterator itE = coll->end();
101  for (; it != itE; ++it)
102  {
103  if (((*it)->pdgId()==pdg) || pdg==0 )
104  {
105  double rtu = Delta::R(*it,eta,phi);
106  if ( rtu < deltaR )
107  {
108  index = l_idx;
109  deltaR = rtu;
110  l_return = true;
111  }
112  }
113  ++l_idx;
114  }
115  return l_return;
116  }
117 
118  inline bool R (const double eta, const double phi, const TruthParticleContainer *coll, int &index,
119  double &deltaR, const int pdg, const bool genOnly=true)
120  {
121  deltaR = 10000.; // big value
122  bool l_return = false;
123  int l_idx = 0;
124  TruthParticleContainer::const_iterator it = coll->begin();
125  TruthParticleContainer::const_iterator itE = coll->end();
126  if ( genOnly )
127  {
128  for (; it != itE; ++it)
129  {
130  if ( !HepMC::is_simulation_particle(*it) ) // only generator particles
131  {
132  if (((*it)->pdgId()==pdg) || pdg==0 )
133  {
134  double rtu = Delta::R(*it,eta,phi);
135  if ( rtu < deltaR )
136  {
137  index = l_idx;
138  deltaR = rtu;
139  l_return = true;
140  }
141  }
142  }
143  ++l_idx;
144  }
145  return l_return;
146  }
147  else return R( eta, phi, coll, index, deltaR, pdg);
148  }
149 
155  template <class COLL>
156  inline bool R (const double eta, const double phi, const double e, COLL *coll,
157  int &index, double &deltaR, const int pdg, double &deltaE)
158  {
159  deltaR = 1.0e+20; // big value
160  deltaE = 1.0e+20;
161  bool l_return = false;
162  int l_idx = 0;
163  typename COLL::const_iterator it = coll->begin();
164  typename COLL::const_iterator itE = coll->end();
165  for (; it != itE; ++it)
166  {
167  if (((*it)->pdgId()==pdg) || pdg==0 )
168  {
169  double rtu = Delta::R(*it,eta,phi);
170  double dE = fabs( e-(*it)->e() );
171  if ( rtu < deltaR && dE < deltaE)
172  {
173  index = l_idx;
174  deltaR = rtu;
175  deltaE = dE;
176  l_return = true;
177  }
178  }
179  ++l_idx;
180  }
181  return l_return;
182  }
183 
184  inline bool R (const double eta, const double phi, const double e, const TruthParticleContainer *coll,
185  int &index, double &deltaR, const int pdg, double &deltaE, const bool genOnly=true)
186  {
187  deltaR = 1.0e+20; // big value
188  deltaE = 1.0e+20;
189  bool l_return = false;
190  int l_idx = 0;
191  TruthParticleContainer::const_iterator it = coll->begin();
192  TruthParticleContainer::const_iterator itE = coll->end();
193  if ( genOnly )
194  {
195  for (; it != itE; ++it)
196  {
197  if ( !HepMC::is_simulation_particle(*it) ) // only generator particles
198  {
199  if (((*it)->pdgId()==pdg) || pdg==0 )
200  {
201  double rtu = Delta::R(*it,eta,phi);
202  double dE = fabs( e-(*it)->e() );
203  if ( rtu < deltaR && dE < deltaE)
204  {
205  index = l_idx;
206  deltaR = rtu;
207  deltaE = dE;
208  l_return = true;
209  }
210  }
211  }
212  ++l_idx;
213  }
214  return l_return;
215  }
216  else return R(eta, phi, e, coll, index, deltaR, pdg, deltaE);
217  }
218 
224  template <class COLL>
225  inline bool R (const INavigable4Momentum *t, COLL *coll, int &index, double &deltaR, const int pdg)
226  {
227  return R (t->eta(), t->phi(), coll, index, deltaR, pdg);
228  }
229 
230  inline bool R (const INavigable4Momentum *t, const TruthParticleContainer *coll, int &index, double &deltaR,
231  const int pdg, const bool genOnly=true)
232  {
233  return R (t->eta(), t->phi(), coll, index, deltaR, pdg, genOnly);
234  }
235 
241  template <class COLL>
242  inline bool R (const INavigable4Momentum *t, COLL *coll, int &index, double &deltaR,
243  const int pdg, double &deltaE)
244  {
245  return R (t->eta(), t->phi(), t->e(), coll, index, deltaR, pdg, deltaE);
246  }
247 
248  inline bool R (const INavigable4Momentum *t, const TruthParticleContainer *coll, int &index,
249  double &deltaR, const int pdg, double &deltaE, const bool genOnly=true)
250  {
251  return R (t->eta(), t->phi(), t->e(), coll, index, deltaR, pdg, deltaE, genOnly);
252  }
253 
257  template <class COLL>
258  inline bool R (const double eta, const double phi, COLL *coll, int &index, double &deltaR)
259  {
260  deltaR = 10000.; // big value
261  bool l_return = false;
262  int l_idx = 0;
263  typename COLL::const_iterator it = coll->begin();
264  typename COLL::const_iterator itE = coll->end();
265  for (; it != itE; ++it)
266  {
267  double rtu = Delta::R(*it,eta,phi);
268  if ( rtu < deltaR )
269  {
270  index = l_idx;
271  deltaR = rtu;
272  l_return = true;
273  }
274  ++l_idx;
275  }
276  return l_return;
277  }
278 
279  inline bool R (const double eta, const double phi, const TruthParticleContainer *coll,
280  int &index, double &deltaR, const bool genOnly=true)
281  {
282  deltaR = 10000.; // big value
283  bool l_return = false;
284  int l_idx = 0;
285  TruthParticleContainer::const_iterator it = coll->begin();
286  TruthParticleContainer::const_iterator itE = coll->end();
287  if ( genOnly )
288  {
289  for (; it != itE; ++it)
290  {
291  if ( !HepMC::is_simulation_particle(*it) ) // only generator particles
292  {
293  double rtu = Delta::R(*it,eta,phi);
294  if ( rtu < deltaR )
295  {
296  index = l_idx;
297  deltaR = rtu;
298  l_return = true;
299  }
300  }
301  ++l_idx;
302  }
303  return l_return;
304  }
305  else return R(eta,phi,coll, index, deltaR);
306  }
307 
311  template <class COLL>
312  inline bool R (const double eta, const double phi, const double e, COLL *coll, int &index,
313  double &deltaR, double &deltaE)
314  {
315  deltaR = 1.0e+20; // big value
316  deltaE = 1.0e+20;
317  bool l_return = false;
318  int l_idx = 0;
319  typename COLL::const_iterator it = coll->begin();
320  typename COLL::const_iterator itE = coll->end();
321  for (; it != itE; ++it)
322  {
323  double rtu = Delta::R(*it,eta,phi);
324  double dE = fabs( e-(*it)->e() );
325  if ( rtu < deltaR && dE < deltaE)
326  {
327  index = l_idx;
328  deltaR = rtu;
329  deltaE = dE;
330  l_return = true;
331  }
332  ++l_idx;
333  }
334  return l_return;
335  }
336 
337  inline bool R (const double eta, const double phi, const double e, const TruthParticleContainer *coll, int &index,
338  double &deltaR, double &deltaE, const bool genOnly=true)
339  {
340  deltaR = 1.0e+20; // big value
341  deltaE = 1.0e+20;
342  bool l_return = false;
343  int l_idx = 0;
344  TruthParticleContainer::const_iterator it = coll->begin();
345  TruthParticleContainer::const_iterator itE = coll->end();
346  if ( genOnly )
347  {
348  for (; it != itE; ++it)
349  {
350  if ( !HepMC::is_simulation_particle(*it) ) // only generator particles
351  {
352  double rtu = Delta::R(*it,eta,phi);
353  double dE = fabs( e-(*it)->e() );
354  if ( rtu < deltaR && dE < deltaE)
355  {
356  index = l_idx;
357  deltaR = rtu;
358  deltaE = dE;
359  l_return = true;
360  }
361  }
362  ++l_idx;
363  }
364  return l_return;
365  }
366  else return R(eta, phi, e, coll, index, deltaR, deltaE);
367  }
368 
372  template <class COLL>
373  inline bool R (const INavigable4Momentum *t, COLL *coll, int &index, double &deltaR)
374  {
375  return R (t->eta(), t->phi(), coll, index, deltaR);
376  }
377 
378  inline bool R (const INavigable4Momentum *t, const TruthParticleContainer *coll, int &index, double &deltaR,
379  const bool genOnly=true)
380  {
381  return R (t->eta(), t->phi(), coll, index, deltaR, genOnly);
382  }
383 
387  template <class COLL>
388  inline bool R (const INavigable4Momentum *t, COLL *coll, int &index, double &deltaR, double &deltaE)
389  {
390  return R (t->eta(), t->phi(), t->e(), coll, index, deltaR, deltaE);
391  }
392 
393  inline bool R (const INavigable4Momentum *t, const TruthParticleContainer *coll, int &index,
394  double &deltaR, double &deltaE, const bool genOnly=true)
395  {
396  return R (t->eta(), t->phi(), t->e(), coll, index, deltaR, deltaE, genOnly);
397  }
398 
399  } // end of Math namespace
400 
402 
403 #endif
404 
407  namespace Sort {
408 
409  namespace Private {
410 
411  using CxxUtils::fpcompare::greater;
412 
413  // function object for pT sorting
414  template <class T> inline bool compPt (T a, T b)
415  {
416  return greater (a->pt(), b->pt());
417  }
418 
419  // function object for e sorting
420  template <class T> inline bool compE (T a, T b)
421  {
422  return greater (a->e(), b->e());
423  }
424 
425  // function object for eta sorting
426  template <class T> inline bool compEta (T a, T b)
427  {
428  return greater (a->eta(), b->eta());
429  }
430 
431  // function object for phi sorting
432  template <class T> inline bool compPhi (T a, T b)
433  {
434  return greater (a->phi(), b->phi());
435  }
436 
437  template <bool b>
438  struct dosort_imp
439  {
440  template <class COLL, class COMP>
441  static void dosort (COLL& coll, COMP comp)
442  { std::sort (coll.begin(), coll.end(), comp); }
443  };
444 
445  template <>
446  struct dosort_imp<true>
447  {
448  template <class COLL, class COMP>
449  static void dosort (COLL& coll, COMP comp)
450  { coll.sort(comp); }
451  };
452 
453  template <class COLL, class COMP>
454  inline void dosort (COLL& coll, COMP comp)
455  {
456  typedef typename std::remove_pointer<typename COLL::value_type>::type valtype;
457  typedef DataVector<valtype> dvtype;
458  const bool is_dv = std::is_convertible<COLL&, dvtype&>::value;
459  dosort_imp<is_dv>::dosort (coll, comp);
460  }
461 
462  } // end of Private namespace
463 
468  template <class COLL> inline void pT (COLL *coll)
469  {
470  // sort
471  Private::dosort (*coll, Private::compPt<typename COLL::value_type>);
472  }
473 
478  template <class COLL> inline void e (COLL *coll)
479  {
480  // sort
481  Private::dosort (*coll, Private::compE<typename COLL::value_type>);
482  }
483 
488  template <class COLL> inline void eta (COLL *coll)
489  {
490  // sort
491  Private::dosort (*coll, Private::compEta<typename COLL::value_type>);
492  }
493 
498  template <class COLL> inline void phi (COLL *coll)
499  {
500  // sort
501  Private::dosort (*coll, Private::compPhi<typename COLL::value_type>);
502  }
503  } // end of Sort namespace
504 
505 #ifndef XAOD_STANDALONE
506 
510  namespace Classify {
511 
518  template <class COLL> inline void charge (COLL *coll,
519  std::vector<typename COLL::value_type> &pos,
520  std::vector<typename COLL::value_type> &neg)
521  {
522  pos.clear();
523  neg.clear();
524 
525  // classify
526  typename COLL::const_iterator it = coll->begin();
527  typename COLL::const_iterator itE = coll->end();
528  for (; it != itE; ++it)
529  {
530  if ((*it)->pdgId() > 0)
531  pos.push_back(*it);
532  else if ((*it)->pdgId() < 0)
533  neg.push_back(*it);
534  }
535  }
536  } // end of Classify namespace
537 
538 #endif
539 
540 } // end of AnalysisUtils namespace
541 
542 #endif
plotting.yearwise_luminosity_vs_mu.comp
comp
Definition: yearwise_luminosity_vs_mu.py:24
DataModel_detail::const_iterator
Const iterator class for DataVector/DataList.
Definition: DVLIterator.h:82
AnalysisUtils::Classify::charge
void charge(COLL *coll, std::vector< typename COLL::value_type > &pos, std::vector< typename COLL::value_type > &neg)
classify by charge
Definition: AnalysisMisc.h:518
phi
Scalar phi() const
phi method
Definition: AmgMatrixBasePlugin.h:64
INavigable4Momentum.h
TruthParticleContainer.h
eta
Scalar eta() const
pseudorapidity method
Definition: AmgMatrixBasePlugin.h:79
index
Definition: index.py:1
AnalysisUtils::Sort::Private::dosort_imp
Definition: AnalysisMisc.h:439
skel.it
it
Definition: skel.GENtoEVGEN.py:423
M_PI
#define M_PI
Definition: ActiveFraction.h:11
AnalysisUtils::Sort::Private::dosort_imp< true >::dosort
static void dosort(COLL &coll, COMP comp)
Definition: AnalysisMisc.h:449
AnalysisUtils::Delta::phi
double phi(const INavigable4Momentum *p1, const INavigable4Momentum *p2)
Definition: AnalysisMisc.h:40
athena.value
value
Definition: athena.py:122
AnalysisUtils
utility class to select combination of elements in a collection
Definition: AnalysisCombination.h:16
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
TruthParticleContainer
Definition: PhysicsAnalysis/TruthParticleID/McParticleEvent/McParticleEvent/TruthParticleContainer.h:42
AnalysisUtils::Sort::Private::dosort_imp::dosort
static void dosort(COLL &coll, COMP comp)
Definition: AnalysisMisc.h:441
I4Momentum::hlv
virtual CLHEP::HepLorentzVector hlv() const =0
CLHEP HepLorentzVector.
CxxUtils::fpcompare::greater
bool greater(double a, double b)
Compare two FP numbers, working around x87 precision issues.
Definition: fpcompare.h:140
TruthTest.itE
itE
Definition: TruthTest.py:25
ParticleJetTools::p3
Amg::Vector3D p3(const xAOD::TruthVertex *p)
Definition: ParticleJetLabelCommon.cxx:55
AnalysisUtils::Imass::two
double two(const INavigable4Momentum *p1, const INavigable4Momentum *p2)
Definition: AnalysisMisc.h:72
AnalysisUtils::Sort::eta
void eta(COLL *coll)
sort by eta
Definition: AnalysisMisc.h:488
HepMC::is_simulation_particle
bool is_simulation_particle(const T &p)
Method to establish if a particle (or barcode) was created during the simulation (TODO update to be s...
Definition: MagicNumbers.h:299
python.StandardJetMods.Sort
Sort
Define the simple modifier setups here – those defined in JetRec.
Definition: StandardJetMods.py:35
AnalysisUtils::Sort::e
void e(COLL *coll)
sort by e
Definition: AnalysisMisc.h:478
AnalysisUtils::Match::R
bool R(const double eta, const double phi, COLL *coll, int &index, double &deltaR, const int pdg)
find the closest element in R
Definition: AnalysisMisc.h:94
fpcompare.h
Workaround x86 precision issues for FP inequality comparisons.
AnalysisUtils::Sort::pT
void pT(COLL *coll)
sort by pT
Definition: AnalysisMisc.h:468
I4Momentum::eta
virtual double eta() const =0
pseudo rapidity
AnalysisUtils::Sort::Private::compPt
bool compPt(T a, T b)
Definition: AnalysisMisc.h:414
AnalysisUtils::Delta::R
double R(const INavigable4Momentum *p1, const double v_eta, const double v_phi)
Definition: AnalysisMisc.h:50
DataVector
Derived DataVector<T>.
Definition: DataVector.h:581
I4Momentum::phi
virtual double phi() const =0
phi in [-pi,pi[
AnalysisUtils::Imass::four
double four(const INavigable4Momentum *p1, const INavigable4Momentum *p2, const INavigable4Momentum *p3, const INavigable4Momentum *p4)
Definition: AnalysisMisc.h:76
AnalysisUtils::Sort::phi
void phi(COLL *coll)
sort by phi
Definition: AnalysisMisc.h:498
AnalysisUtils::Sort::Private::compEta
bool compEta(T a, T b)
Definition: AnalysisMisc.h:426
MagicNumbers.h
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:77
AnalysisUtils::Sort::Private::compE
bool compE(T a, T b)
Definition: AnalysisMisc.h:420
python.LumiBlobConversion.pos
pos
Definition: LumiBlobConversion.py:18
DataVector::end
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
DataVector.h
An STL vector of pointers that by default owns its pointed-to elements.
a
TList * a
Definition: liststreamerinfos.cxx:10
INavigable4Momentum
Definition: INavigable4Momentum.h:21
python.CaloScaleNoiseConfig.type
type
Definition: CaloScaleNoiseConfig.py:78
AnalysisUtils::Sort::Private::compPhi
bool compPhi(T a, T b)
Definition: AnalysisMisc.h:432
AnalysisUtils::Sort::Private::dosort
void dosort(COLL &coll, COMP comp)
Definition: AnalysisMisc.h:454
makeComparison.deltaR
float deltaR
Definition: makeComparison.py:36
DataVector::begin
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
Generated on Tue Jun 25 2024 21:10:05 for ATLAS Offline Software by   doxygen 1.8.18