Skip to content

Commit d769748

Browse files
pbuehlerrolavick
andauthored
[PWGUD,Tutorial] Adding analysis task decayTreeAnalyzer (#9366)
Co-authored-by: rolavick <roman.lavicka@cern.ch>
1 parent baf57de commit d769748

File tree

12 files changed

+2572
-88
lines changed

12 files changed

+2572
-88
lines changed

PWGUD/Core/CMakeLists.txt

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -56,3 +56,11 @@ o2physics_add_library(UPCCutparHolder
5656
o2physics_target_root_dictionary(UPCCutparHolder
5757
HEADERS UPCCutparHolder.h
5858
LINKDEF UPCCutparHolderLinkDef.h)
59+
60+
o2physics_add_library(decayTree
61+
SOURCES decayTree.cxx
62+
PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore ROOT::EG RapidJSON::RapidJSON)
63+
64+
o2physics_target_root_dictionary(decayTree
65+
HEADERS decayTree.h
66+
LINKDEF decayTreeLinkDef.h)

PWGUD/Core/UDHelpers.h

Lines changed: 58 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -18,6 +18,7 @@
1818

1919
#include <vector>
2020
#include <bitset>
21+
2122
#include "TLorentzVector.h"
2223
#include "Framework/Logger.h"
2324
#include "DataFormatsFT0/Digit.h"
@@ -42,7 +43,7 @@ template <bool onlyPV, typename std::enable_if<onlyPV>::type* = nullptr, typenam
4243
int8_t netCharge(TCs tracks)
4344
{
4445
int8_t nch = 0;
45-
for (auto track : tracks) {
46+
for (const auto& track : tracks) {
4647
if (track.isPVContributor()) {
4748
nch += track.sign();
4849
}
@@ -55,7 +56,7 @@ template <bool onlyPV, typename std::enable_if<!onlyPV>::type* = nullptr, typena
5556
int8_t netCharge(TCs tracks)
5657
{
5758
int8_t nch = 0;
58-
for (auto track : tracks) {
59+
for (const auto& track : tracks) {
5960
nch += track.sign();
6061
}
6162
return nch;
@@ -67,7 +68,7 @@ template <bool onlyPV, typename std::enable_if<onlyPV>::type* = nullptr, typenam
6768
float rPVtrwTOF(TCs tracks, int nPVTracks)
6869
{
6970
float rpvrwTOF = 0.;
70-
for (auto& track : tracks) {
71+
for (const auto& track : tracks) {
7172
if (track.isPVContributor() && track.hasTOF()) {
7273
rpvrwTOF += 1.;
7374
}
@@ -83,7 +84,7 @@ template <bool onlyPV, typename std::enable_if<!onlyPV>::type* = nullptr, typena
8384
float rPVtrwTOF(TCs tracks, int nPVTracks)
8485
{
8586
float rpvrwTOF = 0.;
86-
for (auto& track : tracks) {
87+
for (const auto& track : tracks) {
8788
if (track.hasTOF()) {
8889
rpvrwTOF += 1.;
8990
}
@@ -115,14 +116,14 @@ T compatibleBCs(B const& bc, uint64_t const& meanBC, int const& deltaBC, T const
115116
auto bcIter = bcs.iteratorAt(bc.globalIndex());
116117

117118
// range of BCs to consider
118-
uint64_t minBC = (uint64_t)deltaBC < meanBC ? meanBC - (uint64_t)deltaBC : 0;
119-
uint64_t maxBC = meanBC + (uint64_t)deltaBC;
119+
uint64_t minBC = static_cast<uint64_t>(deltaBC) < meanBC ? meanBC - static_cast<uint64_t>(deltaBC) : 0;
120+
uint64_t maxBC = meanBC + static_cast<uint64_t>(deltaBC);
120121
LOGF(debug, " minBC %d maxBC %d bcIterator %d (%d) #BCs %d", minBC, maxBC, bcIter.globalBC(), bcIter.globalIndex(), bcs.size());
121122

122123
// check [min,max]BC to overlap with [bcs.iteratorAt([0,bcs.size() - 1])
123124
if (maxBC < bcs.iteratorAt(0).globalBC() || minBC > bcs.iteratorAt(bcs.size() - 1).globalBC()) {
124125
LOGF(info, "<compatibleBCs> No overlap of [%d, %d] and [%d, %d]", minBC, maxBC, bcs.iteratorAt(0).globalBC(), bcs.iteratorAt(bcs.size() - 1).globalBC());
125-
return T{{bcs.asArrowTable()->Slice(0, 0)}, (uint64_t)0};
126+
return T{{bcs.asArrowTable()->Slice(0, 0)}, static_cast<uint64_t>(0)};
126127
}
127128

128129
// find slice of BCs table with BC in [minBC, maxBC]
@@ -156,7 +157,7 @@ T compatibleBCs(B const& bc, uint64_t const& meanBC, int const& deltaBC, T const
156157
}
157158

158159
// create bc slice
159-
T bcslice{{bcs.asArrowTable()->Slice(minBCId, maxBCId - minBCId + 1)}, (uint64_t)minBCId};
160+
T bcslice{{bcs.asArrowTable()->Slice(minBCId, maxBCId - minBCId + 1)}, static_cast<uint64_t>(minBCId)};
160161
bcs.copyIndexBindings(bcslice);
161162
LOGF(debug, " size of slice %d", bcslice.size());
162163
return bcslice;
@@ -171,7 +172,7 @@ T compatibleBCs(C const& collision, int ndt, T const& bcs, int nMinBCs = 7)
171172

172173
// return if collisions has no associated BC
173174
if (!collision.has_foundBC() || ndt < 0) {
174-
return T{{bcs.asArrowTable()->Slice(0, 0)}, (uint64_t)0};
175+
return T{{bcs.asArrowTable()->Slice(0, 0)}, static_cast<uint64_t>(0)};
175176
}
176177

177178
// get associated BC
@@ -196,7 +197,7 @@ template <typename T>
196197
T compatibleBCs(uint64_t const& meanBC, int const& deltaBC, T const& bcs)
197198
{
198199
// find BC with globalBC ~ meanBC
199-
uint64_t ind = (uint64_t)(bcs.size() / 2);
200+
uint64_t ind = static_cast<uint64_t>(bcs.size() / 2);
200201
auto bcIter = bcs.iteratorAt(ind);
201202

202203
return compatibleBCs(bcIter, meanBC, deltaBC, bcs);
@@ -212,7 +213,7 @@ T MCcompatibleBCs(F const& collision, int const& ndt, T const& bcs, int const& n
212213
// return if collisions has no associated BC
213214
if (!collision.has_foundBC()) {
214215
LOGF(debug, "Collision %i - no BC found!", collision.globalIndex());
215-
return T{{bcs.asArrowTable()->Slice(0, 0)}, (uint64_t)0};
216+
return T{{bcs.asArrowTable()->Slice(0, 0)}, static_cast<uint64_t>(0)};
216217
}
217218

218219
// get associated BC
@@ -245,19 +246,19 @@ bool hasGoodPID(DGCutparHolder diffCuts, TC track)
245246
track.tpcNSigmaPi(),
246247
track.tpcNSigmaKa(),
247248
track.tpcNSigmaPr());
248-
if (TMath::Abs(track.tpcNSigmaEl()) < diffCuts.maxNSigmaTPC()) {
249+
if (std::abs(track.tpcNSigmaEl()) < diffCuts.maxNSigmaTPC()) {
249250
return true;
250251
}
251-
if (TMath::Abs(track.tpcNSigmaMu()) < diffCuts.maxNSigmaTPC()) {
252+
if (std::abs(track.tpcNSigmaMu()) < diffCuts.maxNSigmaTPC()) {
252253
return true;
253254
}
254-
if (TMath::Abs(track.tpcNSigmaPi()) < diffCuts.maxNSigmaTPC()) {
255+
if (std::abs(track.tpcNSigmaPi()) < diffCuts.maxNSigmaTPC()) {
255256
return true;
256257
}
257-
if (TMath::Abs(track.tpcNSigmaKa()) < diffCuts.maxNSigmaTPC()) {
258+
if (std::abs(track.tpcNSigmaKa()) < diffCuts.maxNSigmaTPC()) {
258259
return true;
259260
}
260-
if (TMath::Abs(track.tpcNSigmaPr()) < diffCuts.maxNSigmaTPC()) {
261+
if (std::abs(track.tpcNSigmaPr()) < diffCuts.maxNSigmaTPC()) {
261262
return true;
262263
}
263264

@@ -268,19 +269,19 @@ bool hasGoodPID(DGCutparHolder diffCuts, TC track)
268269
track.tofNSigmaPi(),
269270
track.tofNSigmaKa(),
270271
track.tofNSigmaPr());
271-
if (TMath::Abs(track.tofNSigmaEl()) < diffCuts.maxNSigmaTOF()) {
272+
if (std::abs(track.tofNSigmaEl()) < diffCuts.maxNSigmaTOF()) {
272273
return true;
273274
}
274-
if (TMath::Abs(track.tofNSigmaMu()) < diffCuts.maxNSigmaTOF()) {
275+
if (std::abs(track.tofNSigmaMu()) < diffCuts.maxNSigmaTOF()) {
275276
return true;
276277
}
277-
if (TMath::Abs(track.tofNSigmaPi()) < diffCuts.maxNSigmaTOF()) {
278+
if (std::abs(track.tofNSigmaPi()) < diffCuts.maxNSigmaTOF()) {
278279
return true;
279280
}
280-
if (TMath::Abs(track.tofNSigmaKa()) < diffCuts.maxNSigmaTOF()) {
281+
if (std::abs(track.tofNSigmaKa()) < diffCuts.maxNSigmaTOF()) {
281282
return true;
282283
}
283-
if (TMath::Abs(track.tofNSigmaPr()) < diffCuts.maxNSigmaTOF()) {
284+
if (std::abs(track.tofNSigmaPr()) < diffCuts.maxNSigmaTOF()) {
284285
return true;
285286
}
286287
}
@@ -741,14 +742,14 @@ int64_t sameMCCollision(T tracks, aod::McCollisions, aod::McParticles)
741742
colID = mccol.globalIndex();
742743
} else {
743744
if (colID != mccol.globalIndex()) {
744-
return (int64_t)-1;
745+
return static_cast<int64_t>(-1);
745746
}
746747
}
747748
} else {
748-
return (int64_t)-1;
749+
return static_cast<int64_t>(-1);
749750
}
750751
} else {
751-
return (int64_t)-1;
752+
return static_cast<int64_t>(-1);
752753
}
753754
}
754755

@@ -761,7 +762,7 @@ int64_t sameMCCollision(T tracks, aod::McCollisions, aod::McParticles)
761762
template <typename T>
762763
bool isPythiaCDE(T MCparts)
763764
{
764-
for (auto mcpart : MCparts) {
765+
for (const auto& mcpart : MCparts) {
765766
if (mcpart.pdgCode() == 9900110) {
766767
return true;
767768
}
@@ -780,43 +781,56 @@ bool isSTARLightJPsimumu(T MCparts)
780781
} else {
781782
if (MCparts.iteratorAt(0).pdgCode() != 443013)
782783
return false;
783-
if (abs(MCparts.iteratorAt(1).pdgCode()) != 13)
784+
if (std::abs(MCparts.iteratorAt(1).pdgCode()) != 13)
784785
return false;
785786
if (MCparts.iteratorAt(2).pdgCode() != -MCparts.iteratorAt(1).pdgCode())
786787
return false;
787788
}
788789
return true;
789790
}
790791

792+
// -----------------------------------------------------------------------------
793+
// PbPb di electron production
794+
// [15, 11, 13], [15, 11, 13]
795+
template <typename T>
796+
bool isUpcgen(T MCparts)
797+
{
798+
if (MCparts.size() < 4) {
799+
return false;
800+
} else {
801+
auto pid1 = std::abs(MCparts.iteratorAt(0).pdgCode());
802+
auto pid2 = std::abs(MCparts.iteratorAt(1).pdgCode());
803+
if (pid1 != 11 && pid1 != 13 && pid1 != 15)
804+
return false;
805+
if (pid2 != 11 && pid2 != 13 && pid2 != 15)
806+
return false;
807+
}
808+
return true;
809+
}
810+
791811
// -----------------------------------------------------------------------------
792812
// In pp events produced with GRANIITTI the stack starts with
793-
// 22212/22212/99/22212/2212/99/90
813+
// 22212/22212/22212/2212/[211,321,]/[211,321,]
794814
template <typename T>
795815
bool isGraniittiCDE(T MCparts)
796816
{
797817

798-
for (auto MCpart : MCparts) {
799-
LOGF(debug, " MCpart.pdgCode() %d", MCpart.pdgCode());
800-
}
801-
LOGF(debug, "");
818+
// for (auto MCpart : MCparts) {
819+
// LOGF(info, " MCpart.pdgCode() %d", MCpart.pdgCode());
820+
// }
821+
// LOGF(debug, "");
802822

803-
if (MCparts.size() < 7) {
823+
if (MCparts.size() < 6) {
804824
return false;
805825
} else {
806826
if (MCparts.iteratorAt(0).pdgCode() != 2212)
807827
return false;
808828
if (MCparts.iteratorAt(1).pdgCode() != 2212)
809829
return false;
810-
if (MCparts.iteratorAt(2).pdgCode() != 99)
830+
if (MCparts.iteratorAt(2).pdgCode() != 2212)
811831
return false;
812832
if (MCparts.iteratorAt(3).pdgCode() != 2212)
813833
return false;
814-
if (MCparts.iteratorAt(4).pdgCode() != 2212)
815-
return false;
816-
if (MCparts.iteratorAt(5).pdgCode() != 99)
817-
return false;
818-
if (MCparts.iteratorAt(6).pdgCode() != 90)
819-
return false;
820834
}
821835

822836
return true;
@@ -843,6 +857,11 @@ int isOfInterest(T MCparts)
843857
return 3;
844858
}
845859

860+
// Upcgen
861+
if (isUpcgen(MCparts)) {
862+
return 4;
863+
}
864+
846865
return 0;
847866
}
848867

0 commit comments

Comments
 (0)