Skip to content

Commit 8abbdaf

Browse files
amatyjaamatyja
committed
Avoid TLorentzVector and pdg magic numbers
1 parent 2b699a6 commit 8abbdaf

File tree

1 file changed

+83
-48
lines changed

1 file changed

+83
-48
lines changed

‎PWGUD/Tasks/upcTauTau13topo.cxx‎

Lines changed: 83 additions & 48 deletions
Original file line numberDiff line numberDiff line change
@@ -24,7 +24,8 @@
2424
// #include "TDatabasePDG.h" // not recommended in o2
2525
#include "Framework/O2DatabasePDGPlugin.h"
2626

27-
#include "TLorentzVector.h"
27+
// #include "TLorentzVector.h"
28+
#include "Math/Vector4D.h"
2829
// #include "Common/DataModel/EventSelection.h"
2930
// #include "Common/DataModel/TrackSelectionTables.h"
3031

@@ -1635,7 +1636,14 @@ struct TauTau13topo {
16351636
return std::sqrt(E * E - px * px - py * py - pz * pz);
16361637
}
16371638

1638-
float calculateDeltaPhi(TLorentzVector p, TLorentzVector p1)
1639+
float energy(float px, float py, float pz, float mass)
1640+
// Just a simple function to return energy
1641+
{
1642+
return std::sqrt(px * px + py * py + pz * pz + mass * mass);
1643+
}
1644+
1645+
// float calculateDeltaPhi(TLorentzVector p, TLorentzVector p1)
1646+
float calculateDeltaPhi(ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> p, ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> p1)
16391647
{
16401648
// float delta = p.Phi();
16411649
float delta = RecoDecay::constrainAngle(p.Phi());
@@ -2124,7 +2132,8 @@ struct TauTau13topo {
21242132
int nEtaIn15 = 0;
21252133
int nITSbits = -1;
21262134
int npT100 = 0;
2127-
TLorentzVector p;
2135+
// TLorentzVector p;
2136+
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> p;
21282137
// auto const pionMass = MassPiPlus;
21292138
// auto const electronMass = MassElectron;
21302139
bool flagGlobalCheck = true;
@@ -2133,7 +2142,8 @@ struct TauTau13topo {
21332142
// loop over PV contributors
21342143
for (const auto& trk : PVContributors) {
21352144
qtot += trk.sign();
2136-
p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2145+
// p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2146+
p.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassPiPlus));
21372147
registry.get<TH1>(HIST("global/hTrackPtPV"))->Fill(p.Pt());
21382148
if (std::abs(p.Eta()) < trkEtacut)
21392149
nEtaIn15++; // 1.5 is a default
@@ -2207,7 +2217,8 @@ struct TauTau13topo {
22072217
registry.get<TH1>(HIST("global1/hNTracks"))->Fill(dgtracks.size());
22082218
registry.get<TH1>(HIST("global1/hNTracksPV"))->Fill(PVContributors.size());
22092219
for (const auto& trk : PVContributors) {
2210-
p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2220+
// p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2221+
p.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassPiPlus));
22112222
registry.get<TH1>(HIST("global1/hTrackPtPV"))->Fill(p.Pt());
22122223
registry.get<TH2>(HIST("global1/hTrackEtaPhiPV"))->Fill(p.Eta(), p.Phi());
22132224
}
@@ -2367,23 +2378,27 @@ struct TauTau13topo {
23672378
// 12 34 | 01 23 |//1 //6 | 0 5 |counter<3?counter:5-counter counter<3?0:1
23682379
// 13 24 | 02 13 |//2 //5 | 1 4 |
23692380
// 14 23 | 03 12 |//3 //4 | 2 3 |
2370-
TLorentzVector p1, p2;
2371-
TLorentzVector gammaPair[3][2];
2372-
float invMass2El[3][2];
2381+
// TLorentzVector p1, p2;
2382+
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> p1, p2;
2383+
// TLorentzVector gammaPair[3][2];
2384+
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> gammaPair[3][2];
2385+
float invMass2El[3][2]; // inv. mass squared
23732386
int counterTmp = 0;
23742387
bool flagIMGam2ePV[4] = {true, true, true, true};
23752388

23762389
for (const auto& trk : PVContributors) {
2377-
p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassElectron);
2390+
// p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassElectron);
2391+
p.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassElectron));
23782392
for (const auto& trk1 : PVContributors) {
23792393
if (trk.index() >= trk1.index())
23802394
continue;
23812395
if (trk1.hasTPC())
23822396
nPiHasTPC[trk.index()]++;
2383-
p1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), MassElectron);
2384-
invMass2El[(counterTmp < 3 ? counterTmp : 5 - counterTmp)][(counterTmp < 3 ? 0 : 1)] = (p + p1).Mag2();
2397+
// p1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), MassElectron);
2398+
p1.SetXYZT(trk1.px(), trk1.py(), trk1.pz(), energy(trk1.px(), trk1.py(), trk1.pz(), MassElectron));
2399+
invMass2El[(counterTmp < 3 ? counterTmp : 5 - counterTmp)][(counterTmp < 3 ? 0 : 1)] = (p + p1).mag2();
23852400
gammaPair[(counterTmp < 3 ? counterTmp : 5 - counterTmp)][(counterTmp < 3 ? 0 : 1)] = (p + p1);
2386-
registry.get<TH1>(HIST("control/cut0/hInvMass2ElAll"))->Fill((p + p1).Mag2());
2401+
registry.get<TH1>(HIST("control/cut0/hInvMass2ElAll"))->Fill((p + p1).mag2());
23872402
counterTmp++;
23882403
if ((p + p1).M() < 0.015) {
23892404
flagIMGam2ePV[trk.index()] = false;
@@ -2393,15 +2408,17 @@ struct TauTau13topo {
23932408
} // end of loop over PVContributors
23942409

23952410
// first loop to add all the tracks together
2396-
p = TLorentzVector(0., 0., 0., 0.);
2411+
// p = TLorentzVector(0., 0., 0., 0.);
2412+
p.SetXYZT(0., 0., 0., 0.);
23972413
for (const auto& trk : PVContributors) {
2398-
p1.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2414+
// p1.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2415+
p1.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassPiPlus));
23992416
p += p1;
24002417
scalarPtsum += trk.pt();
24012418
} // end of loop over PVContributors
24022419

24032420
float pttot = p.Pt();
2404-
float mass4pi = p.Mag();
2421+
float mass4pi = p.mag();
24052422

24062423
TVector3 v1(0, 0, 0);
24072424
TVector3 vtmp(0, 0, 0);
@@ -2423,9 +2440,11 @@ struct TauTau13topo {
24232440
registry.get<TH1>(HIST("global/hTrkCheck"))->Fill(tmpTrkCheck);
24242441

24252442
// inv mass of 3pi + 1e
2426-
p1.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2427-
p2.SetXYZM(trk.px(), trk.py(), trk.pz(), MassElectron);
2428-
mass3pi1e[counterTmp] = (p - p1 + p2).Mag();
2443+
// p1.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2444+
p1.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassPiPlus));
2445+
// p2.SetXYZM(trk.px(), trk.py(), trk.pz(), MassElectron);
2446+
p2.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassElectron));
2447+
mass3pi1e[counterTmp] = (p - p1 + p2).mag();
24292448

24302449
v1.SetXYZ(trk.px(), trk.py(), trk.pz());
24312450
for (const auto& trk1 : PVContributors) {
@@ -2472,14 +2491,15 @@ struct TauTau13topo {
24722491
trkTime[counterTmp] = trk.trackTime();
24732492
trkTimeRes[counterTmp] = trk.trackTimeRes();
24742493

2475-
p1.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2494+
// p1.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
2495+
p1.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(),MassPiPlus));
24762496
tmpMomentum[counterTmp] = p1.P();
24772497
tmpPt[counterTmp] = p1.Pt();
24782498
tmpDedx[counterTmp] = trk.tpcSignal();
24792499
tmpTofNsigmaEl[counterTmp] = trk.tofNSigmaEl();
24802500

24812501
deltaPhiTmp = calculateDeltaPhi(p - p1, p1);
2482-
pi3invMass[counterTmp] = (p - p1).Mag();
2502+
pi3invMass[counterTmp] = (p - p1).mag();
24832503
pi3pt[counterTmp] = (p - p1).Pt();
24842504
pi3deltaPhi[counterTmp] = deltaPhiTmp;
24852505
pi3assymav[counterTmp] = (p1.Pt() - (scalarPtsum - p1.Pt()) / 3.) / (p1.Pt() + (scalarPtsum - p1.Pt()) / 3.);
@@ -3342,7 +3362,8 @@ struct TauTau13topo {
33423362
countGen++;
33433363
if (mcParticle.isPhysicalPrimary()) {
33443364
countBoth++;
3345-
if (mcParticle.pdgCode() != 22 && std::abs(mcParticle.pdgCode()) != 12 && std::abs(mcParticle.pdgCode()) != 14 && std::abs(mcParticle.pdgCode()) != 16 && mcParticle.pdgCode() != 130 && mcParticle.pdgCode() != 111) {
3365+
// if (mcParticle.pdgCode() != 22 && std::abs(mcParticle.pdgCode()) != 12 && std::abs(mcParticle.pdgCode()) != 14 && std::abs(mcParticle.pdgCode()) != 16 && mcParticle.pdgCode() != 130 && mcParticle.pdgCode() != 111) {
3366+
if (mcParticle.pdgCode() != kGamma && std::abs(mcParticle.pdgCode()) != kNuE && std::abs(mcParticle.pdgCode()) != kNuMu && std::abs(mcParticle.pdgCode()) != kNuTau && mcParticle.pdgCode() != kK0Long && mcParticle.pdgCode() != kPi0) {
33463367
countCharged++;
33473368

33483369
registryMC.get<TH1>(HIST("globalMC/hMCetaGen"))->Fill(mcParticle.eta());
@@ -3352,7 +3373,8 @@ struct TauTau13topo {
33523373

33533374
if (mcParticle.has_mothers()) {
33543375
auto const& mother = mcParticle.mothers_first_as<aod::McParticles>();
3355-
if (std::abs(mother.pdgCode()) == 15) {
3376+
// if (std::abs(mother.pdgCode()) == 15) {
3377+
if (std::abs(mother.pdgCode()) == kTauMinus) {
33563378
countChargedFromTau++;
33573379
} // mother is tau
33583380
} // mc particle has mother
@@ -3363,7 +3385,8 @@ struct TauTau13topo {
33633385
//
33643386
// tau+/-
33653387
//
3366-
if (std::abs(mcParticle.pdgCode()) == 15) { // tau+/-
3388+
// if (std::abs(mcParticle.pdgCode()) == 15) { // tau+/-
3389+
if (std::abs(mcParticle.pdgCode()) == kTauMinus) { // 15 = tau+/-
33673390
countTau++;
33683391
if (countTau <= 2) {
33693392
etaTau[countTau - 1] = mcParticle.eta();
@@ -3380,15 +3403,15 @@ struct TauTau13topo {
33803403
if (mcParticle.has_daughters()) {
33813404
for (const auto& daughter : mcParticle.daughters_as<aod::McParticles>()) {
33823405
// pions from tau
3383-
if (std::abs(daughter.pdgCode()) == 211) { // 211 = pi+
3406+
if (std::abs(daughter.pdgCode()) == kPiPlus) { // 211 = pi+
33843407
pionCounter++;
33853408
tmpPionIndex = daughter.index(); // returns index of daughter of tau, not in the event, not in the MC particles
33863409
if (std::abs(daughter.eta()) > 0.9)
33873410
partFromTauInEta = false;
33883411
} // end of pion check
33893412
// electron from tau
3390-
if (std::abs(daughter.pdgCode()) == 11) { // 11 = electron
3391-
if (daughter.pdgCode() == 11)
3413+
if (std::abs(daughter.pdgCode()) == kElectron) { // 11 = electron
3414+
if (daughter.pdgCode() == kElectron)
33923415
flagElPlusElMinus = true;
33933416
registryMC.get<TH1>(HIST("electronMC/hMCeta"))->Fill(daughter.eta());
33943417
registryMC.get<TH1>(HIST("electronMC/hMCphi"))->Fill(daughter.phi());
@@ -3403,8 +3426,8 @@ struct TauTau13topo {
34033426
partFromTauInEta = false;
34043427
} // end of electron check
34053428
// muon from tau
3406-
if (std::abs(daughter.pdgCode()) == 13) {
3407-
if (daughter.pdgCode() == 13)
3429+
if (std::abs(daughter.pdgCode()) == kMuonMinus) { // 13
3430+
if (daughter.pdgCode() == kMuonMinus) // 13
34083431
flagMuPlusMuMinus = true;
34093432
muonFound = !muonFound;
34103433
partPt = static_cast<float>(daughter.pt());
@@ -3420,7 +3443,7 @@ struct TauTau13topo {
34203443
singlePionFound = true;
34213444
singlePionIndex = tmpPionIndex;
34223445
auto mcPartTmp = mcParticle.daughters_as<aod::McParticles>().begin() + singlePionIndex;
3423-
if (mcPartTmp.pdgCode() == -211)
3446+
if (mcPartTmp.pdgCode() == kPiMinus) // -211
34243447
flagPiPlusPiMinus = true;
34253448
partPt = static_cast<float>(mcPartTmp.pt());
34263449
// motherOfSinglePionIndex = mcParticle.index();
@@ -3561,7 +3584,8 @@ struct TauTau13topo {
35613584
bool tauInRapidity = true;
35623585
bool partFromTauInEta = true;
35633586

3564-
TLorentzVector tmp(0., 0., 0., 0.);
3587+
// TLorentzVector tmp(0., 0., 0., 0.);
3588+
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> tmp(0., 0., 0., 0.);
35653589

35663590
// get reconstructed collisions associated to generated collision
35673591
// auto const& collisions = collisionsFull.sliceBy(colPerMcCollision, mcCollision.globalIndex());
@@ -3574,23 +3598,24 @@ struct TauTau13topo {
35743598
for (const auto& mcParticle : mcParticles) {
35753599
// LOGF(info, "<tautau13topo_MC> mcParticle pdg %d", mcParticle.pdgCode());
35763600
if (mcParticle.isPhysicalPrimary()) {
3577-
if (mcParticle.pdgCode() != 22 && std::abs(mcParticle.pdgCode()) != 12 && std::abs(mcParticle.pdgCode()) != 14 && std::abs(mcParticle.pdgCode()) != 16 && mcParticle.pdgCode() != 130 && mcParticle.pdgCode() != 111) {
3601+
// if (mcParticle.pdgCode() != 22 && std::abs(mcParticle.pdgCode()) != 12 && std::abs(mcParticle.pdgCode()) != 14 && std::abs(mcParticle.pdgCode()) != 16 && mcParticle.pdgCode() != 130 && mcParticle.pdgCode() != 111) {
3602+
if (mcParticle.pdgCode() != kGamma && std::abs(mcParticle.pdgCode()) != kNuE && std::abs(mcParticle.pdgCode()) != kNuMu && std::abs(mcParticle.pdgCode()) != kNuTau && mcParticle.pdgCode() != kK0Long && mcParticle.pdgCode() != kPi0) {
35783603
if (mcParticle.has_mothers()) {
35793604
auto const& mother = mcParticle.mothers_first_as<aod::UDMcParticles>();
35803605
// LOGF(info, "<tautau13topo_MC> mcParticle has mother %d",mother.pdgCode());
35813606
tmp.SetPxPyPzE(mother.px(), mother.py(), mother.pz(), mother.e());
3582-
if (std::abs(mother.pdgCode()) == 15) {
3607+
if (std::abs(mother.pdgCode()) == kTauMinus) { // 15
35833608
if (std::abs(rapidity(mother.e(), mother.pz())) > 0.9)
35843609
// if (std::abs(tmp.Rapidity()) > 0.9)
35853610
tauInRapidity = false;
35863611
if (std::abs(eta(mcParticle.px(), mcParticle.py(), mcParticle.pz())) > 0.9)
35873612
// if (std::abs(tmp.Eta()) > 0.9)
35883613
partFromTauInEta = false;
35893614

3590-
if (std::abs(mcParticle.pdgCode()) == 11) {
3615+
if (std::abs(mcParticle.pdgCode()) == kElectron) { // 11
35913616
index1ProngMC = mcParticle.index();
35923617
is1ProngElectronMC = true;
3593-
} else if (std::abs(mcParticle.pdgCode()) == 13) {
3618+
} else if (std::abs(mcParticle.pdgCode()) == kMuonMinus) { // 13
35943619
index1ProngMC = mcParticle.index();
35953620
// is1ProngMuonMC = true;
35963621
}
@@ -3638,10 +3663,10 @@ struct TauTau13topo {
36383663
// // motherIndex[i] = (tmpMC.mothers_first_as<aod::UDMcParticles>()).globalIndex();
36393664
// }
36403665
int motherIndex1Pi = motherIndex[0];
3641-
int motherIndexNew = -1;
3666+
int motherIndexNew = 3; // was -1
36423667
int nDifferences = 0;
36433668
for (int i = 1; i < 4; i++) {
3644-
if (motherIndex1Pi != motherIndex[i]) { // the same mother index
3669+
if (motherIndex1Pi != motherIndex[i]) { // the same mother (tau) index
36453670
nDifferences++;
36463671
motherIndexNew = i;
36473672
}
@@ -3654,7 +3679,7 @@ struct TauTau13topo {
36543679
// if (!onlyPi) LOGF(info, "ERROR: should be 4 pions, but they are not!");
36553680
} // end of special check for pi + 3pi
36563681

3657-
int index3ProngMC[3];
3682+
int index3ProngMC[3] = {0, 0, 0}; // initialised of request
36583683
if (index1ProngMC > 0) { // electron or muon case + 3pi
36593684
int index3pi = 0;
36603685
for (int i = 0; i < 4; i++) {
@@ -3673,7 +3698,8 @@ struct TauTau13topo {
36733698
auto const& tmpPion2MC = mcParticles.begin() + index3ProngMC[1];
36743699
auto const& tmpPion3MC = mcParticles.begin() + index3ProngMC[2];
36753700

3676-
if (std::abs(tmpPion1MC.pdgCode()) == 211 && std::abs(tmpPion2MC.pdgCode()) == 211 && std::abs(tmpPion3MC.pdgCode()) == 211)
3701+
//if (std::abs(tmpPion1MC.pdgCode()) == 211 && std::abs(tmpPion2MC.pdgCode()) == 211 && std::abs(tmpPion3MC.pdgCode()) == 211)
3702+
if (std::abs(tmpPion1MC.pdgCode()) == kPiPlus && std::abs(tmpPion2MC.pdgCode()) == kPiPlus && std::abs(tmpPion3MC.pdgCode()) == kPiPlus) // 211 211 211
36773703
is3prong3PiMC = true;
36783704

36793705
//
@@ -4003,8 +4029,10 @@ struct TauTau13topo {
40034029
registryMC.get<TH1>(HIST("global1MCrec/hNTracks"))->Fill(groupedTracks.size());
40044030
registryMC.get<TH1>(HIST("global1MCrec/hNTracksPV"))->Fill(nPVTracks);
40054031

4006-
TLorentzVector p, p1;
4007-
p.SetXYZM(0., 0., 0., 0.);
4032+
// TLorentzVector p, p1;
4033+
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> p, p1;
4034+
// p.SetXYZM(0., 0., 0., 0.);
4035+
p.SetXYZT(0., 0., 0., 0.);
40084036
TVector3 v1(0, 0, 0);
40094037
TVector3 v2(0, 0, 0);
40104038
float scalarPtsum = 0;
@@ -4038,10 +4066,13 @@ struct TauTau13topo {
40384066
float tmpPhiData = phi(tmptrack.px(), tmptrack.py());
40394067
registryMC.get<TH2>(HIST("global1MCrec/hTrackEtaPhiPV"))->Fill(tmpEtaData, tmpPhiData);
40404068
registryMC.get<TH1>(HIST("global1MCrec/hTrackPtPV"))->Fill(tmptrack.pt());
4041-
p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), MassPiPlus); // in case of ghost
4069+
// p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), MassPiPlus); // in case of ghost
4070+
p1.SetXYZT(v1.X(), v1.Y(), v1.Z(), energy(v1.X(), v1.Y(), v1.Z(), MassPiPlus)); // in case of ghost
40424071

40434072
if (trackMCId[i] >= 0) {
4044-
p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), (std::abs(tmptrack.udMcParticle().pdgCode()) == 211 ? MassPiPlus : MassElectron));
4073+
// p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), (std::abs(tmptrack.udMcParticle().pdgCode()) == 211 ? MassPiPlus : MassElectron));
4074+
// p1.SetXYZT(v1.X(), v1.Y(), v1.Z(), energy(v1.X(), v1.Y(), v1.Z(), (std::abs(tmptrack.udMcParticle().pdgCode()) == 211 ? MassPiPlus : MassElectron)));
4075+
p1.SetXYZT(v1.X(), v1.Y(), v1.Z(), energy(v1.X(), v1.Y(), v1.Z(), (std::abs(tmptrack.udMcParticle().pdgCode()) == kPiPlus ? MassPiPlus : MassElectron))); // 211
40454076
float tmpPt = pt(tmptrack.udMcParticle().px(), tmptrack.udMcParticle().py());
40464077
float tmpEta = eta(tmptrack.udMcParticle().px(), tmptrack.udMcParticle().py(), tmptrack.udMcParticle().pz());
40474078
float tmpPhi = phi(tmptrack.udMcParticle().px(), tmptrack.udMcParticle().py());
@@ -4156,7 +4187,7 @@ struct TauTau13topo {
41564187
//
41574188
// temporary control variables per event with combinatorics
41584189
float pttot = p.Pt();
4159-
float mass4pi = p.Mag();
4190+
float mass4pi = p.mag();
41604191
int counterTmp = 0;
41614192

41624193
float nSigmaEl[4];
@@ -4199,9 +4230,11 @@ struct TauTau13topo {
41994230
auto const tmptrack = groupedTracks.begin() + trackId[i];
42004231
// if (tmptrack.hasTOF()) trkHasTof[i] = true;
42014232
v1.SetXYZ(tmptrack.px(), tmptrack.py(), tmptrack.pz());
4202-
p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), MassPiPlus); // in case of ghost
4233+
// p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), MassPiPlus); // in case of ghost
4234+
p1.SetXYZT(v1.X(), v1.Y(), v1.Z(), energy(v1.X(), v1.Y(), v1.Z(), MassPiPlus)); // in case of ghost
42034235
if (trackMCId[i] >= 0) {
4204-
p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), (i == matchedElIndexToData ? MassElectron : MassPiPlus));
4236+
// p1.SetXYZM(v1.X(), v1.Y(), v1.Z(), (i == matchedElIndexToData ? MassElectron : MassPiPlus));
4237+
p1.SetXYZT(v1.X(), v1.Y(), v1.Z(), energy(v1.X(), v1.Y(), v1.Z(), (i == matchedElIndexToData ? MassElectron : MassPiPlus)));
42054238
}
42064239

42074240
nSigmaEl[counterTmp] = tmptrack.tpcNSigmaEl();
@@ -4243,7 +4276,7 @@ struct TauTau13topo {
42434276
tmpTofNsigmaEl[counterTmp] = tmptrack.tofNSigmaEl();
42444277

42454278
deltaPhiTmp = calculateDeltaPhi(p - p1, p1);
4246-
pi3invMass[counterTmp] = (p - p1).Mag();
4279+
pi3invMass[counterTmp] = (p - p1).mag();
42474280
pi3pt[counterTmp] = (p - p1).Pt();
42484281
pi3deltaPhi[counterTmp] = deltaPhiTmp;
42494282
pi3assymav[counterTmp] = (p1.Pt() - (scalarPtsum - p1.Pt()) / 3.) / (p1.Pt() + (scalarPtsum - p1.Pt()) / 3.);
@@ -4994,10 +5027,12 @@ struct TauTau13topo {
49945027
int npT100 = 0;
49955028
// int qtot = 0;
49965029
int8_t qtot = 0;
4997-
TLorentzVector p;
5030+
// TLorentzVector p;
5031+
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double>> p;
49985032
for (const auto& trk : PVContributors) {
49995033
qtot += trk.sign();
5000-
p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
5034+
// p.SetXYZM(trk.px(), trk.py(), trk.pz(), MassPiPlus);
5035+
p.SetXYZT(trk.px(), trk.py(), trk.pz(), energy(trk.px(), trk.py(), trk.pz(), MassPiPlus));
50015036
if (std::abs(p.Eta()) < trkEtacut)
50025037
nEtaIn15++; // 0.9 is a default
50035038
if (trk.pt() > 0.1)

0 commit comments

Comments
 (0)