Calculates eigenvalues and eigenvectors.
18{
19 Eigen::Matrix3f sphericityTensor;
20
21
22
23
24
25
26
27
28
29
30 double elements[9] = {0.};
31
32
33 double norm = 0;
34
35
37 B2WARNING("The particle list has less than 2 elements. The sphericity matrix will not be calculated");
38 return;
39 }
40
42 elements[0] += p.X() * p.X();
43 elements[1] += p.X() * p.Y();
44 elements[2] += p.X() * p.Z();
45
46 elements[3] += p.Y() * p.X();
47 elements[4] += p.Y() * p.Y();
48 elements[5] += p.Y() * p.Z();
49
50 elements[6] += p.Z() * p.X();
51 elements[7] += p.Z() * p.Y();
52 elements[8] += p.Z() * p.Z();
53 norm += p.P2();
54 }
55
56 for (short i = 0; i < 3; i++) {
57 for (short j = 0; j < 3; j++) {
58 sphericityTensor(i, j) = elements[i * 3 + j] / norm;
59 }
60 }
61
62 auto eigenVals = sphericityTensor.eigenvalues();
63 Eigen::ComplexEigenSolver<Eigen::MatrixXcf> ces(sphericityTensor);
64
65
66
67
68
69 short order[3] = {0, 1, 2};
70
71 std::vector<float> tmpLambda;
72 for (short i = 0; i < 3; i++) {
73 tmpLambda.push_back(eigenVals[i].real());
74 }
75
76
77 order[0] = std::distance(tmpLambda.begin(), std::max_element(tmpLambda.begin(), tmpLambda.end()));
78
79 order[2] = std::distance(tmpLambda.begin(), std::min_element(tmpLambda.begin(), tmpLambda.end()));
80
81 order[1] = (short)(3.1 - (order[0] + order[2]));
82
83 for (short i = 0; i < 3; i++) {
84 short n = order[i];
86 auto eigenVector = ces.eigenvectors().col(n);
90 }
91
92 return;
93}
std::vector< ROOT::Math::PxPyPzEVector > m_momenta
The particles' momenta.
ROOT::Math::XYZVector m_eVector[3]
The eigenvectors.
double m_lambda[3]
The eigenvalues.
1.13.2