|
14 | 14 | // Class to handle analysis variables |
15 | 15 | // |
16 | 16 |
|
17 | | -#ifndef VarManager_H |
18 | | -#define VarManager_H |
19 | | - |
20 | | -#include <TObject.h> |
21 | | -#include <TString.h> |
22 | | -#include <Math/Vector4D.h> |
23 | | -#include "Math/Vector3D.h" |
24 | | -#include "Math/GenVector/Boost.h" |
25 | | -#include <TRandom.h> |
| 17 | +#ifndef PWGDQ_CORE_VARMANAGER_H_ |
| 18 | +#define PWGDQ_CORE_VARMANAGER_H_ |
26 | 19 |
|
27 | 20 | #include <vector> |
28 | 21 | #include <map> |
29 | 22 | #include <cmath> |
30 | 23 | #include <iostream> |
| 24 | +#include <utility> |
| 25 | + |
| 26 | +#include <TObject.h> |
| 27 | +#include <TString.h> |
| 28 | +#include "TRandom.h" |
| 29 | +#include "Math/Vector4D.h" |
| 30 | +#include "Math/Vector3D.h" |
| 31 | +#include "Math/GenVector/Boost.h" |
31 | 32 |
|
32 | 33 | #include "Framework/DataTypes.h" |
33 | 34 | #include "ReconstructionDataFormats/Track.h" |
@@ -315,6 +316,7 @@ class VarManager : public TObject |
315 | 316 | kPairPtDau, |
316 | 317 | kPairEta, |
317 | 318 | kPairPhi, |
| 319 | + kPairPhiv, |
318 | 320 | kDeltaEta, |
319 | 321 | kDeltaPhi, |
320 | 322 | kDeltaPhiSym, |
@@ -671,13 +673,13 @@ void VarManager::FillTrack(T const& track, float* values) |
671 | 673 | values[kIsGlobalTrack] = track.filteringFlags() & (uint64_t(1) << 0); |
672 | 674 | values[kIsGlobalTrackSDD] = track.filteringFlags() & (uint64_t(1) << 1); |
673 | 675 |
|
674 | | - values[kIsLegFromGamma] = bool(track.filteringFlags() & (uint64_t(1) << 2)); |
675 | | - values[kIsLegFromK0S] = bool(track.filteringFlags() & (uint64_t(1) << 3)); |
676 | | - values[kIsLegFromLambda] = bool(track.filteringFlags() & (uint64_t(1) << 4)); |
677 | | - values[kIsLegFromAntiLambda] = bool(track.filteringFlags() & (uint64_t(1) << 5)); |
678 | | - values[kIsLegFromOmega] = bool(track.filteringFlags() & (uint64_t(1) << 6)); |
| 676 | + values[kIsLegFromGamma] = static_cast<bool>(track.filteringFlags() & (uint64_t(1) << 2)); |
| 677 | + values[kIsLegFromK0S] = static_cast<bool>(track.filteringFlags() & (uint64_t(1) << 3)); |
| 678 | + values[kIsLegFromLambda] = static_cast<bool>(track.filteringFlags() & (uint64_t(1) << 4)); |
| 679 | + values[kIsLegFromAntiLambda] = static_cast<bool>(track.filteringFlags() & (uint64_t(1) << 5)); |
| 680 | + values[kIsLegFromOmega] = static_cast<bool>(track.filteringFlags() & (uint64_t(1) << 6)); |
679 | 681 |
|
680 | | - values[kIsProtonFromLambdaAndAntiLambda] = bool((values[kIsLegFromLambda] * track.sign() > 0) || (values[kIsLegFromAntiLambda] * (-track.sign()) > 0)); |
| 682 | + values[kIsProtonFromLambdaAndAntiLambda] = static_cast<bool>((values[kIsLegFromLambda] * track.sign() > 0) || (values[kIsLegFromAntiLambda] * (-track.sign()) > 0)); |
681 | 683 | } |
682 | 684 | } |
683 | 685 |
|
@@ -922,6 +924,70 @@ void VarManager::FillPair(T1 const& t1, T2 const& t2, float* values) |
922 | 924 | values[kQuadDCAsigXY] = std::sqrt((dca1sig * dca1sig + dca2sig * dca2sig) / 2); |
923 | 925 | } |
924 | 926 | } |
| 927 | + if (fgUsedVars[kPairPhiv]) { |
| 928 | + // cos(phiv) = w*a /|w||a| |
| 929 | + // with w = u x v |
| 930 | + // and a = u x z / |u x z| , unit vector perpendicular to v12 and z-direction (magnetic field) |
| 931 | + // u = v12 / |v12| , the unit vector of v12 |
| 932 | + // v = v1 x v2 / |v1 x v2| , unit vector perpendicular to v1 and v2 |
| 933 | + |
| 934 | + float bz = fgFitterTwoProngBarrel.getBz(); |
| 935 | + |
| 936 | + // ordering of tracks, so v1 has larger momentum |
| 937 | + if (v1.P() < v2.P()) { |
| 938 | + ROOT::Math::PtEtaPhiMVector v3 = v1; |
| 939 | + v1 = v2; |
| 940 | + v2 = v3; |
| 941 | + } |
| 942 | + |
| 943 | + // momentum of e+ and e- in (ax,ay,az) axis. Note that az=0 by definition. |
| 944 | + // vector product of pep X pem |
| 945 | + float vpx = 0, vpy = 0, vpz = 0; |
| 946 | + if (t1.sign() * t2.sign() > 0) { // Like Sign |
| 947 | + if (bz * t1.sign() < 0) { |
| 948 | + vpx = v1.Py() * v2.Pz() - v1.Pz() * v2.Py(); |
| 949 | + vpy = v1.Pz() * v2.Px() - v1.Px() * v2.Pz(); |
| 950 | + vpz = v1.Px() * v2.Py() - v1.Py() * v2.Px(); |
| 951 | + } else { |
| 952 | + vpx = v2.Py() * v1.Pz() - v2.Pz() * v1.Py(); |
| 953 | + vpy = v2.Pz() * v1.Px() - v2.Px() * v1.Pz(); |
| 954 | + vpz = v2.Px() * v1.Py() - v2.Py() * v1.Px(); |
| 955 | + } |
| 956 | + } else { // Unlike Sign |
| 957 | + if (bz * t1.sign() > 0) { |
| 958 | + vpx = v1.Py() * v2.Pz() - v1.Pz() * v2.Py(); |
| 959 | + vpy = v1.Pz() * v2.Px() - v1.Px() * v2.Pz(); |
| 960 | + vpz = v1.Px() * v2.Py() - v1.Py() * v2.Px(); |
| 961 | + } else { |
| 962 | + vpx = v2.Py() * v1.Pz() - v2.Pz() * v1.Py(); |
| 963 | + vpy = v2.Pz() * v1.Px() - v2.Px() * v1.Pz(); |
| 964 | + vpz = v2.Px() * v1.Py() - v2.Py() * v1.Px(); |
| 965 | + } |
| 966 | + } |
| 967 | + |
| 968 | + // unit vector of pep X pem |
| 969 | + float vx = vpx / TMath::Sqrt(vpx * vpx + vpy * vpy + vpz * vpz); |
| 970 | + float vy = vpy / TMath::Sqrt(vpx * vpx + vpy * vpy + vpz * vpz); |
| 971 | + float vz = vpz / TMath::Sqrt(vpx * vpx + vpy * vpy + vpz * vpz); |
| 972 | + |
| 973 | + float px = v12.Px(); |
| 974 | + float py = v12.Py(); |
| 975 | + float pz = v12.Pz(); |
| 976 | + |
| 977 | + // unit vector of (pep+pem) |
| 978 | + float ux = px / TMath::Sqrt(px * px + py * py + pz * pz); |
| 979 | + float uy = py / TMath::Sqrt(px * px + py * py + pz * pz); |
| 980 | + float uz = pz / TMath::Sqrt(px * px + py * py + pz * pz); |
| 981 | + float ax = uy / TMath::Sqrt(ux * ux + uy * uy); |
| 982 | + float ay = -ux / TMath::Sqrt(ux * ux + uy * uy); |
| 983 | + |
| 984 | + // The third axis defined by vector product (ux,uy,uz)X(vx,vy,vz) |
| 985 | + float wx = uy * vz - uz * vy; |
| 986 | + float wy = uz * vx - ux * vz; |
| 987 | + // by construction, (wx,wy,wz) must be a unit vector. Measure angle between (wx,wy,wz) and (ax,ay,0). |
| 988 | + // The angle between them should be small if the pair is conversion. This function then returns values close to pi! |
| 989 | + values[kPairPhiv] = TMath::ACos(wx * ax + wy * ay); // phiv in [0,pi] //cosPhiV = wx * ax + wy * ay; |
| 990 | + } |
925 | 991 | } |
926 | 992 |
|
927 | 993 | template <int pairType, typename T1, typename T2> |
@@ -1399,4 +1465,4 @@ void VarManager::FillDileptonHadron(T1 const& dilepton, T2 const& hadron, float* |
1399 | 1465 | values[kDeltaEta] = dilepton.eta() - hadron.eta(); |
1400 | 1466 | } |
1401 | 1467 | } |
1402 | | -#endif |
| 1468 | +#endif // PWGDQ_CORE_VARMANAGER_H_ |
0 commit comments