Rivet  3.1.3
MathUtils.hh
1 // -*- C++ -*-
2 #ifndef RIVET_MathUtils_HH
3 #define RIVET_MathUtils_HH
4 
5 #include "Rivet/Math/MathConstants.hh"
6 #include <type_traits>
7 #include <cassert>
8 
9 namespace Rivet {
10 
11 
13 
14 
16 
17 
22  template <typename NUM>
23  inline typename std::enable_if<std::is_floating_point<NUM>::value, bool>::type
24  isZero(NUM val, double tolerance=1e-8) {
25  return fabs(val) < tolerance;
26  }
27 
32  template <typename NUM>
33  inline typename std::enable_if<std::is_integral<NUM>::value, bool>::type
34  isZero(NUM val, double=1e-5) { //< NB. unused tolerance parameter for ints, still needs a default value!
35  return val == 0;
36  }
37 
39  template <typename NUM>
40  inline typename std::enable_if<std::is_floating_point<NUM>::value, bool>::type
41  isNaN(NUM val) { return std::isnan(val); }
42 
44  template <typename NUM>
45  inline typename std::enable_if<std::is_floating_point<NUM>::value, bool>::type
46  notNaN(NUM val) { return !std::isnan(val); }
47 
53  template <typename N1, typename N2>
54  inline typename std::enable_if<
55  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value &&
56  (std::is_floating_point<N1>::value || std::is_floating_point<N2>::value), bool>::type
57  fuzzyEquals(N1 a, N2 b, double tolerance=1e-5) {
58  const double absavg = (std::abs(a) + std::abs(b))/2.0;
59  const double absdiff = std::abs(a - b);
60  const bool rtn = (isZero(a) && isZero(b)) || absdiff < tolerance*absavg;
61  return rtn;
62  }
63 
68  template <typename N1, typename N2>
69  inline typename std::enable_if<
70  std::is_integral<N1>::value && std::is_integral<N2>::value, bool>::type
71  fuzzyEquals(N1 a, N2 b, double) { //< NB. unused tolerance parameter for ints, still needs a default value!
72  return a == b;
73  }
74 
75 
79  template <typename N1, typename N2>
80  inline typename std::enable_if<
81  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value, bool>::type
82  fuzzyGtrEquals(N1 a, N2 b, double tolerance=1e-5) {
83  return a > b || fuzzyEquals(a, b, tolerance);
84  }
85 
86 
90  template <typename N1, typename N2>
91  inline typename std::enable_if<
92  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value, bool>::type
93  fuzzyLessEquals(N1 a, N2 b, double tolerance=1e-5) {
94  return a < b || fuzzyEquals(a, b, tolerance);
95  }
96 
98  template <typename N1, typename N2>
99  inline typename std::enable_if<
100  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value,
101  typename std::common_type<N1,N2>::type >::type
102  min(N1 a, N2 b) {
103  return a > b ? b : a;
104  }
105 
107  template <typename N1, typename N2>
108  inline typename std::enable_if<
109  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value,
110  typename std::common_type<N1,N2>::type >::type
111  max(N1 a, N2 b) {
112  return a > b ? a : b;
113  }
114 
116 
117 
119 
120 
125  enum RangeBoundary { OPEN=0, SOFT=0, CLOSED=1, HARD=1 };
126 
130  template <typename N1, typename N2, typename N3>
131  inline typename std::enable_if<
132  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value && std::is_arithmetic<N3>::value, bool>::type
133  inRange(N1 value, N2 low, N3 high,
134  RangeBoundary lowbound=CLOSED, RangeBoundary highbound=OPEN) {
135  if (lowbound == OPEN && highbound == OPEN) {
136  return (value > low && value < high);
137  } else if (lowbound == OPEN && highbound == CLOSED) {
138  return (value > low && value <= high);
139  } else if (lowbound == CLOSED && highbound == OPEN) {
140  return (value >= low && value < high);
141  } else { // if (lowbound == CLOSED && highbound == CLOSED) {
142  return (value >= low && value <= high);
143  }
144  }
145 
150  template <typename N1, typename N2, typename N3>
151  inline typename std::enable_if<
152  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value && std::is_arithmetic<N3>::value, bool>::type
153  fuzzyInRange(N1 value, N2 low, N3 high,
154  RangeBoundary lowbound=CLOSED, RangeBoundary highbound=OPEN) {
155  if (lowbound == OPEN && highbound == OPEN) {
156  return (value > low && value < high);
157  } else if (lowbound == OPEN && highbound == CLOSED) {
158  return (value > low && fuzzyLessEquals(value, high));
159  } else if (lowbound == CLOSED && highbound == OPEN) {
160  return (fuzzyGtrEquals(value, low) && value < high);
161  } else { // if (lowbound == CLOSED && highbound == CLOSED) {
162  return (fuzzyGtrEquals(value, low) && fuzzyLessEquals(value, high));
163  }
164  }
165 
167  template <typename N1, typename N2, typename N3>
168  inline typename std::enable_if<
169  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value && std::is_arithmetic<N3>::value, bool>::type
170  inRange(N1 value, pair<N2, N3> lowhigh,
171  RangeBoundary lowbound=CLOSED, RangeBoundary highbound=OPEN) {
172  return inRange(value, lowhigh.first, lowhigh.second, lowbound, highbound);
173  }
174 
175 
176  // Alternative forms, with snake_case names and boundary types in names rather than as args -- from MCUtils
177 
181  template <typename N1, typename N2, typename N3>
182  inline typename std::enable_if<
183  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value && std::is_arithmetic<N3>::value, bool>::type
184  in_range(N1 val, N2 low, N3 high) {
185  return inRange(val, low, high, CLOSED, OPEN);
186  }
187 
191  template <typename N1, typename N2, typename N3>
192  inline typename std::enable_if<
193  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value && std::is_arithmetic<N3>::value, bool>::type
194  in_closed_range(N1 val, N2 low, N3 high) {
195  return inRange(val, low, high, CLOSED, CLOSED);
196  }
197 
201  template <typename N1, typename N2, typename N3>
202  inline typename std::enable_if<
203  std::is_arithmetic<N1>::value && std::is_arithmetic<N2>::value && std::is_arithmetic<N3>::value, bool>::type
204  in_open_range(N1 val, N2 low, N3 high) {
205  return inRange(val, low, high, OPEN, OPEN);
206  }
207 
209 
211 
212 
214 
215 
217  template <typename NUM>
218  inline typename std::enable_if<std::is_arithmetic<NUM>::value, NUM>::type
219  sqr(NUM a) {
220  return a*a;
221  }
222 
227  // template <typename N1, typename N2>
228  template <typename NUM>
229  inline typename std::enable_if<std::is_arithmetic<NUM>::value, NUM>::type
230  //std::common_type<N1, N2>::type
231  add_quad(NUM a, NUM b) {
232  return sqrt(a*a + b*b);
233  }
234 
239  // template <typename N1, typename N2>
240  template <typename NUM>
241  inline typename std::enable_if<std::is_arithmetic<NUM>::value, NUM>::type
242  //std::common_type<N1, N2, N3>::type
243  add_quad(NUM a, NUM b, NUM c) {
244  return sqrt(a*a + b*b + c*c);
245  }
246 
249  inline double safediv(double num, double den, double fail=0.0) {
250  return (!isZero(den)) ? num/den : fail;
251  }
252 
254  template <typename NUM>
255  inline typename std::enable_if<std::is_arithmetic<NUM>::value, NUM>::type
256  intpow(NUM val, unsigned int exp) {
257  assert(exp >= 0);
258  if (exp == 0) return (NUM) 1;
259  else if (exp == 1) return val;
260  return val * intpow(val, exp-1);
261  }
262 
264  template <typename NUM>
265  inline typename std::enable_if<std::is_arithmetic<NUM>::value, int>::type
266  sign(NUM val) {
267  if (isZero(val)) return ZERO;
268  const int valsign = (val > 0) ? PLUS : MINUS;
269  return valsign;
270  }
271 
273 
274 
276 
277 
279  inline double cdfBW(double x, double mu, double gamma) {
280  // normalize to (0;1) distribution
281  const double xn = (x - mu)/gamma;
282  return std::atan(xn)/M_PI + 0.5;
283  }
284 
286  inline double invcdfBW(double p, double mu, double gamma) {
287  const double xn = std::tan(M_PI*(p-0.5));
288  return gamma*xn + mu;
289  }
290 
292 
293 
295 
296 
303  inline vector<double> linspace(size_t nbins, double start, double end, bool include_end=true) {
304  assert(end >= start);
305  assert(nbins > 0);
306  vector<double> rtn;
307  const double interval = (end-start)/static_cast<double>(nbins);
308  for (size_t i = 0; i < nbins; ++i) {
309  rtn.push_back(start + i*interval);
310  }
311  assert(rtn.size() == nbins);
312  if (include_end) rtn.push_back(end); //< exact end, not result of n * interval
313  return rtn;
314  }
315 
316 
326  inline vector<double> aspace(double step, double start, double end, bool include_end=true, double tol=1e-2) {
327  assert(end >= start);
328  assert(step > 0);
329  vector<double> rtn;
330  double next = start;
331  while (true) {
332  if (next > end) break;
333  rtn.push_back(next);
334  next += step;
335  }
336  if (include_end) {
337  if (end - rtn[rtn.size()-1] > tol*step) rtn.push_back(end);
338  }
339  return rtn;
340  }
341 
342 
350  inline vector<double> logspace(size_t nbins, double start, double end, bool include_end=true) {
351  assert(end >= start);
352  assert(start > 0);
353  assert(nbins > 0);
354  const double logstart = std::log(start);
355  const double logend = std::log(end);
356  const vector<double> logvals = linspace(nbins, logstart, logend, false);
357  assert(logvals.size() == nbins);
358  vector<double> rtn; rtn.reserve(nbins+1);
359  rtn.push_back(start); //< exact start, not exp(log(start))
360  for (size_t i = 1; i < logvals.size(); ++i) {
361  rtn.push_back(std::exp(logvals[i]));
362  }
363  assert(rtn.size() == nbins);
364  if (include_end) rtn.push_back(end); //< exact end, not exp(n * loginterval)
365  return rtn;
366  }
367 
368 
370 
371 
379  inline vector<double> bwspace(size_t nbins, double start, double end, double mu, double gamma) {
380  assert(end >= start);
381  assert(nbins > 0);
382  const double pmin = cdfBW(start, mu, gamma);
383  const double pmax = cdfBW(end, mu, gamma);
384  const vector<double> edges = linspace(nbins, pmin, pmax);
385  assert(edges.size() == nbins+1);
386  vector<double> rtn;
387  for (double edge : edges) {
388  rtn.push_back(invcdfBW(edge, mu, gamma));
389  }
390  assert(rtn.size() == nbins+1);
391  return rtn;
392  }
393 
394 
396  template <typename NUM, typename CONTAINER>
397  inline typename std::enable_if<std::is_arithmetic<NUM>::value && std::is_arithmetic<typename CONTAINER::value_type>::value, int>::type
398  _binIndex(NUM val, const CONTAINER& binedges, bool allow_overflow=false) {
399  if (val < *begin(binedges)) return -1;
400  // CONTAINER::iterator_type itend =
401  if (val >= *(end(binedges)-1)) return allow_overflow ? int(binedges.size())-1 : -1;
402  auto it = std::upper_bound(begin(binedges), end(binedges), val);
403  return std::distance(begin(binedges), --it);
404  }
405 
414  template <typename NUM1, typename NUM2>
415  inline typename std::enable_if<std::is_arithmetic<NUM1>::value && std::is_arithmetic<NUM2>::value, int>::type
416  binIndex(NUM1 val, std::initializer_list<NUM2> binedges, bool allow_overflow=false) {
417  return _binIndex(val, binedges, allow_overflow);
418  }
419 
428  template <typename NUM, typename CONTAINER>
429  inline typename std::enable_if<std::is_arithmetic<NUM>::value && std::is_arithmetic<typename CONTAINER::value_type>::value, int>::type
430  binIndex(NUM val, const CONTAINER& binedges, bool allow_overflow=false) {
431  return _binIndex(val, binedges, allow_overflow);
432  }
433 
435 
436 
438 
439 
442  template <typename NUM>
443  inline typename std::enable_if<std::is_arithmetic<NUM>::value, NUM>::type
444  median(const vector<NUM>& sample) {
445  if (sample.empty()) throw RangeError("Can't compute median of an empty set");
446  vector<NUM> tmp = sample;
447  std::sort(tmp.begin(), tmp.end());
448  const size_t imid = tmp.size()/2; // len1->idx0, len2->idx1, len3->idx1, len4->idx2, ...
449  if (sample.size() % 2 == 0) return (tmp.at(imid-1) + tmp.at(imid)) / 2.0;
450  else return tmp.at(imid);
451  }
452 
453 
456  template <typename NUM>
457  inline typename std::enable_if<std::is_arithmetic<NUM>::value, double>::type
458  mean(const vector<NUM>& sample) {
459  if (sample.empty()) throw RangeError("Can't compute mean of an empty set");
460  double mean = 0.0;
461  for (size_t i = 0; i < sample.size(); ++i) {
462  mean += sample[i];
463  }
464  return mean/sample.size();
465  }
466 
467  // Calculate the error on the mean, assuming Poissonian errors
469  template <typename NUM>
470  inline typename std::enable_if<std::is_arithmetic<NUM>::value, double>::type
471  mean_err(const vector<NUM>& sample) {
472  if (sample.empty()) throw RangeError("Can't compute mean_err of an empty set");
473  double mean_e = 0.0;
474  for (size_t i = 0; i < sample.size(); ++i) {
475  mean_e += sqrt(sample[i]);
476  }
477  return mean_e/sample.size();
478  }
479 
480 
483  template <typename NUM>
484  inline typename std::enable_if<std::is_arithmetic<NUM>::value, double>::type
485  covariance(const vector<NUM>& sample1, const vector<NUM>& sample2) {
486  if (sample1.empty() || sample2.empty()) throw RangeError("Can't compute covariance of an empty set");
487  if (sample1.size() != sample2.size()) throw RangeError("Sizes of samples must be equal for covariance calculation");
488  const double mean1 = mean(sample1);
489  const double mean2 = mean(sample2);
490  const size_t N = sample1.size();
491  double cov = 0.0;
492  for (size_t i = 0; i < N; i++) {
493  const double cov_i = (sample1[i] - mean1)*(sample2[i] - mean2);
494  cov += cov_i;
495  }
496  if (N > 1) return cov/(N-1);
497  else return 0.0;
498  }
499 
502  template <typename NUM>
503  inline typename std::enable_if<std::is_arithmetic<NUM>::value, double>::type
504  covariance_err(const vector<NUM>& sample1, const vector<NUM>& sample2) {
505  if (sample1.empty() || sample2.empty()) throw RangeError("Can't compute covariance_err of an empty set");
506  if (sample1.size() != sample2.size()) throw RangeError("Sizes of samples must be equal for covariance_err calculation");
507  const double mean1 = mean(sample1);
508  const double mean2 = mean(sample2);
509  const double mean1_e = mean_err(sample1);
510  const double mean2_e = mean_err(sample2);
511  const size_t N = sample1.size();
512  double cov_e = 0.0;
513  for (size_t i = 0; i < N; i++) {
514  const double cov_i = (sqrt(sample1[i]) - mean1_e)*(sample2[i] - mean2) +
515  (sample1[i] - mean1)*(sqrt(sample2[i]) - mean2_e);
516  cov_e += cov_i;
517  }
518  if (N > 1) return cov_e/(N-1);
519  else return 0.0;
520  }
521 
522 
525  template <typename NUM>
526  inline typename std::enable_if<std::is_arithmetic<NUM>::value, double>::type
527  correlation(const vector<NUM>& sample1, const vector<NUM>& sample2) {
528  const double cov = covariance(sample1, sample2);
529  const double var1 = covariance(sample1, sample1);
530  const double var2 = covariance(sample2, sample2);
531  const double correlation = cov/sqrt(var1*var2);
532  const double corr_strength = correlation*sqrt(var2/var1);
533  return corr_strength;
534  }
535 
538  template <typename NUM>
539  inline typename std::enable_if<std::is_arithmetic<NUM>::value, double>::type
540  correlation_err(const vector<NUM>& sample1, const vector<NUM>& sample2) {
541  const double cov = covariance(sample1, sample2);
542  const double var1 = covariance(sample1, sample1);
543  const double var2 = covariance(sample2, sample2);
544  const double cov_e = covariance_err(sample1, sample2);
545  const double var1_e = covariance_err(sample1, sample1);
546  const double var2_e = covariance_err(sample2, sample2);
547 
548  // Calculate the correlation
549  const double correlation = cov/sqrt(var1*var2);
550  // Calculate the error on the correlation
551  const double correlation_err = cov_e/sqrt(var1*var2) -
552  cov/(2*pow(3./2., var1*var2)) * (var1_e * var2 + var1 * var2_e);
553 
554  // Calculate the error on the correlation strength
555  const double corr_strength_err = correlation_err*sqrt(var2/var1) +
556  correlation/(2*sqrt(var2/var1)) * (var2_e/var1 - var2*var1_e/pow(2, var2));
557 
558  return corr_strength_err;
559  }
560 
562 
563 
565 
566 
571  inline double _mapAngleM2PITo2Pi(double angle) {
572  double rtn = fmod(angle, TWOPI);
573  if (isZero(rtn)) return 0;
574  assert(rtn >= -TWOPI && rtn <= TWOPI);
575  return rtn;
576  }
577 
579  inline double mapAngleMPiToPi(double angle) {
580  double rtn = _mapAngleM2PITo2Pi(angle);
581  if (isZero(rtn)) return 0;
582  if (rtn > PI) rtn -= TWOPI;
583  if (rtn <= -PI) rtn += TWOPI;
584  assert(rtn > -PI && rtn <= PI);
585  return rtn;
586  }
587 
589  inline double mapAngle0To2Pi(double angle) {
590  double rtn = _mapAngleM2PITo2Pi(angle);
591  if (isZero(rtn)) return 0;
592  if (rtn < 0) rtn += TWOPI;
593  if (rtn == TWOPI) rtn = 0;
594  assert(rtn >= 0 && rtn < TWOPI);
595  return rtn;
596  }
597 
599  inline double mapAngle0ToPi(double angle) {
600  double rtn = fabs(mapAngleMPiToPi(angle));
601  if (isZero(rtn)) return 0;
602  assert(rtn > 0 && rtn <= PI);
603  return rtn;
604  }
605 
607  inline double mapAngle(double angle, PhiMapping mapping) {
608  switch (mapping) {
609  case MINUSPI_PLUSPI:
610  return mapAngleMPiToPi(angle);
611  case ZERO_2PI:
612  return mapAngle0To2Pi(angle);
613  case ZERO_PI:
614  return mapAngle0To2Pi(angle);
615  default:
616  throw Rivet::UserError("The specified phi mapping scheme is not implemented");
617  }
618  }
619 
621 
622 
624 
625 
629  inline double deltaPhi(double phi1, double phi2, bool sign=false) {
630  const double x = mapAngleMPiToPi(phi1 - phi2);
631  return sign ? x : fabs(x);
632  }
633 
637  inline double deltaEta(double eta1, double eta2, bool sign=false) {
638  const double x = eta1 - eta2;
639  return sign ? x : fabs(x);
640  }
641 
645  inline double deltaRap(double y1, double y2, bool sign=false) {
646  const double x = y1 - y2;
647  return sign? x : fabs(x);
648  }
649 
652  inline double deltaR2(double rap1, double phi1, double rap2, double phi2) {
653  const double dphi = deltaPhi(phi1, phi2);
654  return sqr(rap1-rap2) + sqr(dphi);
655  }
656 
659  inline double deltaR(double rap1, double phi1, double rap2, double phi2) {
660  return sqrt(deltaR2(rap1, phi1, rap2, phi2));
661  }
662 
664  inline double rapidity(double E, double pz) {
665  if (isZero(E - pz)) {
666  throw std::runtime_error("Divergent positive rapidity");
667  return DBL_MAX;
668  }
669  if (isZero(E + pz)) {
670  throw std::runtime_error("Divergent negative rapidity");
671  return -DBL_MAX;
672  }
673  return 0.5*log((E+pz)/(E-pz));
674  }
675 
677 
678 
681  inline double mT(double pT1, double pT2, double dphi) {
682  return sqrt(2*pT1*pT2 * (1 - cos(dphi)) );
683  }
684 
685 
686 }
687 
688 
689 #endif
Definition: MC_Cent_pPb.hh:10
std::enable_if< std::is_floating_point< NUM >::value, bool >::type notNaN(NUM val)
Check if a number is non-NaN.
Definition: MathUtils.hh:46
std::enable_if< std::is_arithmetic< NUM1 >::value &&std::is_arithmetic< NUM2 >::value, int >::type binIndex(NUM1 val, std::initializer_list< NUM2 > binedges, bool allow_overflow=false)
Return the bin index of the given value, val, given a vector of bin edges.
Definition: MathUtils.hh:416
double rapidity(double E, double pz)
Calculate a rapidity value from the supplied energy E and longitudinal momentum pz.
Definition: MathUtils.hh:664
double mapAngleMPiToPi(double angle)
Map an angle into the range (-PI, PI].
Definition: MathUtils.hh:579
double safediv(double num, double den, double fail=0.0)
Definition: MathUtils.hh:249
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value &&std::is_arithmetic< N3 >::value, bool >::type in_range(N1 val, N2 low, N3 high)
Boolean function to determine if value is within the given range.
Definition: MathUtils.hh:184
static const double PI
Definition: MathConstants.hh:13
Error specialisation for where the problem is between the chair and the computer. ...
Definition: Exceptions.hh:55
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value, typename std::common_type< N1, N2 >::type >::type min(N1 a, N2 b)
Get the minimum of two numbers.
Definition: MathUtils.hh:102
std::enable_if< std::is_floating_point< NUM >::value, bool >::type isNaN(NUM val)
Check if a number is NaN.
Definition: MathUtils.hh:41
double mapAngle(double angle, PhiMapping mapping)
Map an angle into the enum-specified range.
Definition: MathUtils.hh:607
std::enable_if< std::is_arithmetic< NUM >::value, NUM >::type add_quad(NUM a, NUM b)
Named number-type addition in quadrature operation.
Definition: MathUtils.hh:231
vector< double > linspace(size_t nbins, double start, double end, bool include_end=true)
Make a list of nbins + 1 values equally spaced between start and end inclusive.
Definition: MathUtils.hh:303
double cdfBW(double x, double mu, double gamma)
CDF for the Breit-Wigner distribution.
Definition: MathUtils.hh:279
vector< double > logspace(size_t nbins, double start, double end, bool include_end=true)
Make a list of nbins + 1 values exponentially spaced between start and end inclusive.
Definition: MathUtils.hh:350
double deltaEta(double eta1, double eta2, bool sign=false)
Definition: MathUtils.hh:637
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value, typename std::common_type< N1, N2 >::type >::type max(N1 a, N2 b)
Get the maximum of two numbers.
Definition: MathUtils.hh:111
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value &&std::is_arithmetic< N3 >::value, bool >::type in_open_range(N1 val, N2 low, N3 high)
Boolean function to determine if value is within the given range.
Definition: MathUtils.hh:204
PhiMapping
Enum for range of to be mapped into.
Definition: MathConstants.hh:49
std::enable_if< std::is_arithmetic< NUM >::value, double >::type covariance(const vector< NUM > &sample1, const vector< NUM > &sample2)
Definition: MathUtils.hh:485
double deltaPhi(double phi1, double phi2, bool sign=false)
Calculate the difference between two angles in radians.
Definition: MathUtils.hh:629
vector< double > aspace(double step, double start, double end, bool include_end=true, double tol=1e-2)
Make a list of values equally spaced by step between start and end inclusive.
Definition: MathUtils.hh:326
std::enable_if< std::is_arithmetic< NUM >::value, double >::type mean_err(const vector< NUM > &sample)
Definition: MathUtils.hh:471
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value &&std::is_arithmetic< N3 >::value, bool >::type inRange(N1 value, N2 low, N3 high, RangeBoundary lowbound=CLOSED, RangeBoundary highbound=OPEN)
Determine if value is in the range low to high, for floating point numbers.
Definition: MathUtils.hh:133
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value &&std::is_arithmetic< N3 >::value, bool >::type in_closed_range(N1 val, N2 low, N3 high)
Boolean function to determine if value is within the given range.
Definition: MathUtils.hh:194
double mT(double pT1, double pT2, double dphi)
Definition: MathUtils.hh:681
double invcdfBW(double p, double mu, double gamma)
Inverse CDF for the Breit-Wigner distribution.
Definition: MathUtils.hh:286
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value &&(std::is_floating_point< N1 >::value||std::is_floating_point< N2 >::value), bool >::type fuzzyEquals(N1 a, N2 b, double tolerance=1e-5)
Compare two numbers for equality with a degree of fuzziness.
Definition: MathUtils.hh:57
static const double TWOPI
A pre-defined value of .
Definition: MathConstants.hh:16
double deltaRap(double y1, double y2, bool sign=false)
Definition: MathUtils.hh:645
double deltaR2(double rap1, double phi1, double rap2, double phi2)
Definition: MathUtils.hh:652
double deltaR(double rap1, double phi1, double rap2, double phi2)
Definition: MathUtils.hh:659
RangeBoundary
Definition: MathUtils.hh:125
double mapAngle0To2Pi(double angle)
Map an angle into the range [0, 2PI).
Definition: MathUtils.hh:589
double mapAngle0ToPi(double angle)
Map an angle into the range [0, PI].
Definition: MathUtils.hh:599
std::enable_if< std::is_arithmetic< NUM >::value, double >::type covariance_err(const vector< NUM > &sample1, const vector< NUM > &sample2)
Definition: MathUtils.hh:504
std::enable_if< std::is_floating_point< NUM >::value, bool >::type isZero(NUM val, double tolerance=1e-8)
Compare a number to zero.
Definition: MathUtils.hh:24
double angle(const Vector2 &a, const Vector2 &b)
Angle (in radians) between two 2-vectors.
Definition: Vector2.hh:175
std::enable_if< std::is_arithmetic< NUM >::value, double >::type correlation_err(const vector< NUM > &sample1, const vector< NUM > &sample2)
Definition: MathUtils.hh:540
std::enable_if< std::is_arithmetic< NUM >::value, NUM >::type median(const vector< NUM > &sample)
Definition: MathUtils.hh:444
std::enable_if< std::is_arithmetic< NUM >::value, double >::type mean(const vector< NUM > &sample)
Definition: MathUtils.hh:458
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value, bool >::type fuzzyGtrEquals(N1 a, N2 b, double tolerance=1e-5)
Compare two numbers for >= with a degree of fuzziness.
Definition: MathUtils.hh:82
std::enable_if< std::is_arithmetic< NUM >::value, NUM >::type sqr(NUM a)
Named number-type squaring operation.
Definition: MathUtils.hh:219
std::enable_if< std::is_arithmetic< NUM >::value, NUM >::type intpow(NUM val, unsigned int exp)
A more efficient version of pow for raising numbers to integer powers.
Definition: MathUtils.hh:256
Error for e.g. use of invalid bin ranges.
Definition: Exceptions.hh:22
vector< double > bwspace(size_t nbins, double start, double end, double mu, double gamma)
Make a list of nbins + 1 values spaced for equal area Breit-Wigner binning between start and end incl...
Definition: MathUtils.hh:379
std::enable_if< std::is_arithmetic< NUM >::value, int >::type sign(NUM val)
Find the sign of a number.
Definition: MathUtils.hh:266
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value &&std::is_arithmetic< N3 >::value, bool >::type fuzzyInRange(N1 value, N2 low, N3 high, RangeBoundary lowbound=CLOSED, RangeBoundary highbound=OPEN)
Determine if value is in the range low to high, for floating point numbers.
Definition: MathUtils.hh:153
std::enable_if< std::is_arithmetic< N1 >::value &&std::is_arithmetic< N2 >::value, bool >::type fuzzyLessEquals(N1 a, N2 b, double tolerance=1e-5)
Compare two floating point numbers for <= with a degree of fuzziness.
Definition: MathUtils.hh:93
std::enable_if< std::is_arithmetic< NUM >::value, double >::type correlation(const vector< NUM > &sample1, const vector< NUM > &sample2)
Definition: MathUtils.hh:527