Belle II Software  release-05-01-25
HarmonicMoments.cc
1 /* BASF2 (Belle Analysis Framework 2) *
2  * Copyright(C) 2018 - Belle II Collaboration *
3  * *
4  * Author: The Belle II Collaboration *
5  * Contributors: Umberto Tamponi (tamponi@to.infn.it) *
6  * *
7  * This software is provided "as is" without any warranty. *
8  **************************************************************************/
9 
10 
11 #include <analysis/ContinuumSuppression/HarmonicMoments.h>
12 #include "TVector3.h"
13 
14 using namespace Belle2;
15 
16 
17 
19 {
20  // Loop over the particles' momenta
21  for (auto& p : m_momenta) {
22  // Gets momentum and costheta of the vector
23  double pMag = p.Mag();
24  double cTheta = p.Dot(m_axis) / pMag;
25 
26  // Fills the momenta.
27  // This part is quite ugly, but hard-coding the legendre polinimials makes the code
28  // much faster than using the boost libraries, which are implementing the recursive formulas.
29  // This implementation should also be faster than a switch...case one.
30  double cTheta2 = cTheta * cTheta;
31  double cTheta3 = cTheta2 * cTheta;
32  double cTheta4 = cTheta2 * cTheta2;
33 
34  m_moment[0] += pMag;
35  m_moment[1] += pMag * cTheta;
36  m_moment[2] += pMag * 0.5 * (3.*cTheta2 - 1);
37  m_moment[3] += pMag * 0.5 * (5.*cTheta3 - 3.*cTheta);
38  m_moment[4] += pMag * 0.125 * (35.*cTheta4 - 30.*cTheta2 + 3.);
39  }
40  return;
41 }
42 
43 
45 {
46  // Loop over the particles' momenta
47  for (auto& p : m_momenta) {
48  // gets momentum and costheta of the vector
49  double pMag = p.Mag();
50  double cTheta = p.Dot(m_axis) / pMag;
51 
52  // Fills the momenta.
53  // This part is quite ugly, but hard-coding the legendre polinimials makes the code
54  // much faster than using the boost libraries, which are implementing the recursive formulas.
55  // This implementation should also be faster than a switch...case one.
56  double cTheta2 = cTheta * cTheta;
57  double cTheta3 = cTheta2 * cTheta;
58  double cTheta4 = cTheta2 * cTheta2;
59  double cTheta5 = cTheta4 * cTheta;
60  double cTheta6 = cTheta3 * cTheta3;
61  double cTheta7 = cTheta6 * cTheta;
62  double cTheta8 = cTheta4 * cTheta4;
63 
64  m_moment[0] += pMag;
65  m_moment[1] += pMag * cTheta;
66  m_moment[2] += pMag * 0.5 * (3.*cTheta2 - 1);
67  m_moment[3] += pMag * 0.5 * (5.*cTheta3 - 3.*cTheta);
68  m_moment[4] += pMag * 0.125 * (35.*cTheta4 - 30.*cTheta2 + 3.);
69  m_moment[5] += pMag * 0.125 * (63.*cTheta5 - 70 * cTheta3 + 15.*cTheta);
70  m_moment[6] += pMag * 0.0625 * (231.*cTheta6 - 315 * cTheta4 + 105 * cTheta2 - 5.);
71  m_moment[7] += pMag * 0.0625 * (429.*cTheta7 - 693.*cTheta5 + 315.*cTheta3 - 35.*cTheta);
72  m_moment[8] += pMag * 0.0078125 * (6435.*cTheta8 - 12012.*cTheta6 + 6930.*cTheta4 - 1260.*cTheta2 + 35.);
73  }
74  return;
75 }
Belle2::HarmonicMoments::m_moment
double m_moment[9]
The harmonic moments.
Definition: HarmonicMoments.h:113
Belle2::HarmonicMoments::m_axis
TVector3 m_axis
The reference axis.
Definition: HarmonicMoments.h:115
Belle2::HarmonicMoments::calculateAllMoments
void calculateAllMoments()
Calculates the moments up to order 8.
Definition: HarmonicMoments.cc:44
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::HarmonicMoments::calculateBasicMoments
void calculateBasicMoments()
Calculates the moments up to order 4.
Definition: HarmonicMoments.cc:18
Belle2::HarmonicMoments::m_momenta
std::vector< TVector3 > m_momenta
The list of particles.
Definition: HarmonicMoments.h:114