Skip to content

Commit be7091f

Browse files
omvazqueomvazque
committed
Reduces disk space by using TH2 instead of TH3
1 parent 6d792bb commit be7091f

File tree

1 file changed

+16
-288
lines changed

1 file changed

+16
-288
lines changed

‎PWGLF/Tasks/Nuspex/piKpRAA.cxx‎

Lines changed: 16 additions & 288 deletions
Original file line numberDiff line numberDiff line change
@@ -81,10 +81,10 @@ using TracksMC = soa::Join<aod::Tracks, aod::TracksExtra, aod::TrackSelectionExt
8181

8282
static constexpr int kNEtaHists{8};
8383

84-
std::array<std::shared_ptr<TH3>, kNEtaHists> dEdxPiV0{};
85-
std::array<std::shared_ptr<TH3>, kNEtaHists> dEdxPrV0{};
86-
std::array<std::shared_ptr<TH3>, kNEtaHists> dEdxElV0{};
87-
std::array<std::shared_ptr<TH3>, kNEtaHists> dEdxPiTOF{};
84+
std::array<std::shared_ptr<TH2>, kNEtaHists> dEdxPiV0{};
85+
std::array<std::shared_ptr<TH2>, kNEtaHists> dEdxPrV0{};
86+
std::array<std::shared_ptr<TH2>, kNEtaHists> dEdxElV0{};
87+
std::array<std::shared_ptr<TH2>, kNEtaHists> dEdxPiTOF{};
8888
std::array<std::shared_ptr<TH3>, kNEtaHists> dEdx{};
8989
std::array<std::shared_ptr<TH2>, kNEtaHists> nClVsdEdxPiV0{};
9090
std::array<std::shared_ptr<TH2>, kNEtaHists> nClVsdEdxElV0{};
@@ -218,7 +218,6 @@ struct PiKpRAA {
218218
Configurable<bool> selINELgt0{"selINELgt0", true, "Select INEL > 0?"};
219219
Configurable<bool> isApplyFT0CbasedOccupancy{"isApplyFT0CbasedOccupancy", false, "T0C Occu cut"};
220220
Configurable<bool> applyNchSel{"applyNchSel", false, "Use mid-rapidity-based Nch selection"};
221-
Configurable<bool> skipRecoColGTOne{"skipRecoColGTOne", true, "Remove collisions if reconstructed more than once"};
222221

223222
// Event selection
224223
Configurable<float> posZcut{"posZcut", +10.0, "z-vertex position cut"};
@@ -489,10 +488,10 @@ struct PiKpRAA {
489488
for (int i = 0; i < kNEtaHists; ++i) {
490489
dEdx[i] = registry.add<TH3>(Form("dEdx_%s", endingEta[i]), Form("%s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH3F, {axisPt, axisdEdx, axisCent});
491490
pTVsP[i] = registry.add<TH2>(Form("pTVsP_%s", endingEta[i]), Form("%s;Momentum (GeV/#it{c});#it{p}_{T} (GeV/#it{c});", latexEta[i]), kTH2F, {axisPt, axisPt});
492-
dEdxPiV0[i] = registry.add<TH3>(Form("dEdxPiV0_%s", endingEta[i]), Form("#pi^{+} + #pi^{-}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH3F, {axisPtV0s, axisdEdx, axisCent});
493-
dEdxPrV0[i] = registry.add<TH3>(Form("dEdxPrV0_%s", endingEta[i]), Form("p + #bar{p}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH3F, {axisPtV0s, axisdEdx, axisCent});
494-
dEdxElV0[i] = registry.add<TH3>(Form("dEdxElV0_%s", endingEta[i]), Form("e^{+} + e^{-}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH3F, {axisPtV0s, axisdEdx, axisCent});
495-
dEdxPiTOF[i] = registry.add<TH3>(Form("dEdxPiTOF_%s", endingEta[i]), Form("#pi^{+} + #pi^{-}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH3F, {axisPtV0s, axisdEdx, axisCent});
491+
dEdxPiV0[i] = registry.add<TH2>(Form("dEdxPiV0_%s", endingEta[i]), Form("#pi^{+} + #pi^{-}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH2F, {axisPtV0s, axisdEdx});
492+
dEdxPrV0[i] = registry.add<TH2>(Form("dEdxPrV0_%s", endingEta[i]), Form("p + #bar{p}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH2F, {axisPtV0s, axisdEdx});
493+
dEdxElV0[i] = registry.add<TH2>(Form("dEdxElV0_%s", endingEta[i]), Form("e^{+} + e^{-}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH2F, {axisPtV0s, axisdEdx});
494+
dEdxPiTOF[i] = registry.add<TH2>(Form("dEdxPiTOF_%s", endingEta[i]), Form("#pi^{+} + #pi^{-}, %s;Momentum (GeV/#it{c});dE/dx;", latexEta[i]), kTH2F, {axisPtV0s, axisdEdx});
496495
nClVsdEdxPiV0[i] = registry.add<TH2>(Form("NclVsdEdxPiV0_%s", endingEta[i]), Form("%s;#it{N}_{cl} used for PID;dE/dx;", latexEta[i]), kTH2F, {axisNcl, axisdEdx});
497496
nClVsdEdxElV0[i] = registry.add<TH2>(Form("NclVsdEdxElV0_%s", endingEta[i]), Form("%s;#it{N}_{cl} used for PID;dE/dx;", latexEta[i]), kTH2F, {axisNcl, axisdEdx});
498497
nClVsdEdxPrV0[i] = registry.add<TH2>(Form("NclVsdEdxPrV0_%s", endingEta[i]), Form("%s;#it{N}_{cl} used for PID;dE/dx;", latexEta[i]), kTH2F, {axisNcl, axisdEdx});
@@ -525,7 +524,7 @@ struct PiKpRAA {
525524
xtrkSel->SetBinLabel(12, "Passed all");
526525
}
527526

528-
if (doprocessMC || doprocessSim) {
527+
if (doprocessSim) {
529528
registry.add("zPosMC", "Generated Events With at least One Rec. Collision + Sel. criteria;;Entries;", kTH1F, {axisZpos});
530529
registry.add("dcaVsPtPiDec", "Secondary pions from decays;#it{p}_{T} (GeV/#it{c});DCA_{xy} (cm);Centrality Perc.;", kTH3F, {axisPt, axisDCAxy, axisCent});
531530
registry.add("dcaVsPtPrDec", "Secondary protons from decays;#it{p}_{T} (GeV/#it{c});DCA_{xy} (cm);Centrality Perc.;", kTH3F, {axisPt, axisDCAxy, axisCent});
@@ -534,19 +533,6 @@ struct PiKpRAA {
534533
registry.add("NclVsPhip", Form("#LTNcl#GT used for PID;%s (GeV/#it{c});#varphi", titlePorPt.data()), kTProfile2D, {{{axisXPhiCut}, {350, 0.0, 0.35}}});
535534
}
536535

537-
if (doprocessMC) {
538-
539-
registry.add("EventCounterMC", "Event counter", kTH1F, {axisEvent});
540-
541-
registry.add("PtPiVsCent", "", kTH2F, {axisPt, axisCent});
542-
registry.add("PtKaVsCent", "", kTH2F, {axisPt, axisCent});
543-
registry.add("PtPrVsCent", "", kTH2F, {axisPt, axisCent});
544-
545-
registry.add("PtPiVsCentMC", "", kTH2F, {axisPt, axisCent});
546-
registry.add("PtKaVsCentMC", "", kTH2F, {axisPt, axisCent});
547-
registry.add("PtPrVsCentMC", "", kTH2F, {axisPt, axisCent});
548-
}
549-
550536
if (doprocessSim) {
551537

552538
// MC events passing the TVX requirement
@@ -906,7 +892,7 @@ struct PiKpRAA {
906892
// registry.fill(HIST("TOFExpEl2TOF"), momentum, tExpElTOF / tTOF);
907893

908894
if (std::abs((tExpPiTOF / tTOF) - kOne) < v0Selections.maxExpTOFPi) {
909-
dEdxPiTOF[indexEta]->Fill(momentum, dedx, centrality);
895+
dEdxPiTOF[indexEta]->Fill(momentum, dedx);
910896
}
911897
}
912898
}
@@ -1042,8 +1028,8 @@ struct PiKpRAA {
10421028
nClVsdEdxPiV0[posIndexEta]->Fill(posNcl, posTrkdEdx);
10431029
nClVsdEdxpPiV0[posIndexEta]->Fill(posNcl, posTrkdEdx);
10441030
nClVsPpPiV0[negIndexEta]->Fill(negPorPt, negNcl);
1045-
dEdxPiV0[posIndexEta]->Fill(posTrkP, posTrkdEdx, centrality);
1046-
dEdxPiV0[negIndexEta]->Fill(negTrkP, negTrkdEdx, centrality);
1031+
dEdxPiV0[posIndexEta]->Fill(posTrkP, posTrkdEdx);
1032+
dEdxPiV0[negIndexEta]->Fill(negTrkP, negTrkdEdx);
10471033

10481034
if (posTrkP > kMinPMIP && posTrkP < kMaxPMIP && posTrkdEdx > kMindEdxMIP && posTrkdEdx < kMaxdEdxMIP) {
10491035
registry.fill(HIST("dEdxVsEtaPiMIPV0"), posTrkEta, posTrkdEdx);
@@ -1068,7 +1054,7 @@ struct PiKpRAA {
10681054
nClVsPpPrV0[posIndexEta]->Fill(posPorPt, posNcl);
10691055
nClVsdEdxPrV0[posIndexEta]->Fill(posNcl, posTrkdEdx);
10701056
nClVsdEdxpPrV0[posIndexEta]->Fill(posNcl, posTrkdEdx);
1071-
dEdxPrV0[posIndexEta]->Fill(posTrkP, posTrkdEdx, centrality);
1057+
dEdxPrV0[posIndexEta]->Fill(posTrkP, posTrkdEdx);
10721058
}
10731059
if (negTrackCharge < kZero) {
10741060
registry.fill(HIST("nSigPiFromL"), negTrkPt, negTrack.tpcNSigmaPi());
@@ -1104,7 +1090,7 @@ struct PiKpRAA {
11041090
nClVsPpPrV0[negIndexEta]->Fill(negPorPt, negNcl);
11051091
nClVsdEdxPrV0[negIndexEta]->Fill(negNcl, negTrkdEdx);
11061092
nClVsdEdxpPrV0[negIndexEta]->Fill(negNcl, negTrkdEdx);
1107-
dEdxPrV0[negIndexEta]->Fill(negTrkP, negTrkdEdx, centrality);
1093+
dEdxPrV0[negIndexEta]->Fill(negTrkP, negTrkdEdx);
11081094
}
11091095
}
11101096

@@ -1139,8 +1125,8 @@ struct PiKpRAA {
11391125
registry.fill(HIST("dEdxVsEtaElMIPV0p"), posTrkEta, posTrkdEdx);
11401126
registry.fill(HIST("dEdxVsEtaElMIPV0"), negTrkEta, negTrkdEdx);
11411127
registry.fill(HIST("dEdxVsEtaElMIPV0p"), negTrkEta, negTrkdEdx);
1142-
dEdxElV0[posIndexEta]->Fill(posTrkP, posTrkdEdx, centrality);
1143-
dEdxElV0[negIndexEta]->Fill(negTrkP, negTrkdEdx, centrality);
1128+
dEdxElV0[posIndexEta]->Fill(posTrkP, posTrkdEdx);
1129+
dEdxElV0[negIndexEta]->Fill(negTrkP, negTrkdEdx);
11441130
}
11451131
}
11461132
} // v0s
@@ -1149,264 +1135,6 @@ struct PiKpRAA {
11491135

11501136
Preslice<TracksMC> perCollision = aod::track::collisionId;
11511137
Service<o2::framework::O2DatabasePDG> pdg;
1152-
void processMC(aod::McCollisions::iterator const& mccollision, soa::SmallGroups<ColEvSelsMC> const& collisions, BCsRun3 const& /*bcs*/, aod::FT0s const& /*ft0s*/, aod::McParticles const& mcParticles, TracksMC const& tracksMC)
1153-
{
1154-
for (const auto& collision : collisions) {
1155-
// Event selection
1156-
if (!isEventSelected(collision)) {
1157-
continue;
1158-
}
1159-
// MC collision?
1160-
if (!collision.has_mcCollision()) {
1161-
continue;
1162-
}
1163-
1164-
registry.fill(HIST("EventCounterMC"), EvCutLabel::All);
1165-
1166-
if (std::fabs(mccollision.posZ()) > posZcut)
1167-
continue;
1168-
1169-
registry.fill(HIST("zPos"), collision.posZ());
1170-
registry.fill(HIST("zPosMC"), mccollision.posZ());
1171-
registry.fill(HIST("EventCounterMC"), EvCutLabel::VtxZ);
1172-
1173-
const auto& foundBC = collision.foundBC_as<BCsRun3>();
1174-
uint64_t timeStamp{foundBC.timestamp()};
1175-
const int magField{getMagneticField(timeStamp)};
1176-
1177-
if (v0Selections.applyPhiCut) {
1178-
const int nextRunNumber{foundBC.runNumber()};
1179-
if (currentRunNumberPhiSel != nextRunNumber) {
1180-
loadPhiCutSelections(timeStamp);
1181-
currentRunNumberPhiSel = nextRunNumber;
1182-
LOG(info) << "\tcurrentRunNumberPhiSel= " << currentRunNumberPhiSel << " timeStamp = " << timeStamp;
1183-
}
1184-
1185-
// return if phi cut objects are nullptr
1186-
if (!(phiCut.hPhiCutHigh && phiCut.hPhiCutLow))
1187-
return;
1188-
}
1189-
1190-
// Remove collisions if reconstructed more than once
1191-
if (skipRecoColGTOne && (collisions.size() > kOne))
1192-
continue;
1193-
1194-
const float centrality{isT0Ccent ? collision.centFT0C() : collision.centFT0M()};
1195-
registry.fill(HIST("T0Ccent"), centrality);
1196-
1197-
const auto& groupedTracks{tracksMC.sliceBy(perCollision, collision.globalIndex())};
1198-
1199-
// Track selection with Open DCAxy
1200-
for (const auto& track : groupedTracks) {
1201-
// Track Selection
1202-
if (track.eta() < v0Selections.minEtaDaughter || track.eta() > v0Selections.maxEtaDaughter)
1203-
continue;
1204-
1205-
if (track.pt() < v0Selections.minPt || track.pt() > v0Selections.maxPt)
1206-
continue;
1207-
1208-
if (!trkSelGlobalOpenDCAxy.IsSelected(track))
1209-
continue;
1210-
1211-
if (!track.has_mcParticle())
1212-
continue;
1213-
1214-
// Get the MC particle
1215-
const auto& particle{track.mcParticle()};
1216-
auto charge{0.};
1217-
auto* pdgParticle = pdg->GetParticle(particle.pdgCode());
1218-
if (pdgParticle != nullptr)
1219-
charge = pdgParticle->Charge();
1220-
else
1221-
continue;
1222-
1223-
// Is it a charged particle?
1224-
if (std::abs(charge) < kMinCharge)
1225-
continue;
1226-
1227-
float phiPrime{track.phi()};
1228-
phiPrimeFunc(phiPrime, magField, charge);
1229-
1230-
const float pOrPt{v0Selections.usePinPhiSelection ? track.p() : track.pt()};
1231-
if (v0Selections.applyPhiCut) {
1232-
if (!passesPhiSelection(pOrPt, phiPrime))
1233-
continue;
1234-
}
1235-
1236-
const int16_t nclFound{track.tpcNClsFound()};
1237-
const int16_t nclPID{track.tpcNClsPID()};
1238-
const int16_t ncl = v0Selections.useNclsPID ? nclPID : nclFound;
1239-
if (v0Selections.applyNclSel && ncl < v0Selections.minNcl)
1240-
continue;
1241-
1242-
bool isPrimary{false};
1243-
bool isDecay{false};
1244-
bool isMaterial{false};
1245-
if (particle.isPhysicalPrimary())
1246-
isPrimary = true;
1247-
else if (particle.getProcess() == TMCProcess::kPDecay)
1248-
isDecay = true;
1249-
else
1250-
isMaterial = true;
1251-
1252-
bool isPi{false};
1253-
bool isPr{false};
1254-
if (particle.pdgCode() == PDG_t::kPiPlus || particle.pdgCode() == PDG_t::kPiMinus)
1255-
isPi = true;
1256-
else if (particle.pdgCode() == PDG_t::kProton || particle.pdgCode() == PDG_t::kProtonBar)
1257-
isPr = true;
1258-
else
1259-
continue;
1260-
1261-
if (isPrimary && !isDecay && !isMaterial) {
1262-
if (isPi && !isPr)
1263-
registry.fill(HIST("dcaVsPtPi"), track.pt(), track.dcaXY(), centrality);
1264-
if (isPr && !isPi)
1265-
registry.fill(HIST("dcaVsPtPr"), track.pt(), track.dcaXY(), centrality);
1266-
}
1267-
1268-
if (isDecay && !isPrimary && !isMaterial) {
1269-
if (isPi && !isPr)
1270-
registry.fill(HIST("dcaVsPtPiDec"), track.pt(), track.dcaXY(), centrality);
1271-
if (isPr && !isPi)
1272-
registry.fill(HIST("dcaVsPtPrDec"), track.pt(), track.dcaXY(), centrality);
1273-
}
1274-
1275-
if (isMaterial && !isPrimary && !isDecay) {
1276-
if (isPi && !isPr)
1277-
registry.fill(HIST("dcaVsPtPiMat"), track.pt(), track.dcaXY(), centrality);
1278-
if (isPr && !isPi)
1279-
registry.fill(HIST("dcaVsPtPrMat"), track.pt(), track.dcaXY(), centrality);
1280-
}
1281-
}
1282-
1283-
// Global track + DCAxy selections
1284-
for (const auto& track : groupedTracks) {
1285-
// Track Selection
1286-
if (track.eta() < v0Selections.minEtaDaughter || track.eta() > v0Selections.maxEtaDaughter)
1287-
continue;
1288-
1289-
if (track.pt() < v0Selections.minPt || track.pt() > v0Selections.maxPt)
1290-
continue;
1291-
1292-
if (!trkSelGlobal.IsSelected(track))
1293-
continue;
1294-
1295-
// Has MC particle?
1296-
if (!track.has_mcParticle())
1297-
continue;
1298-
1299-
// Get the MC particle
1300-
const auto& particle{track.mcParticle()};
1301-
auto charge{0.};
1302-
auto* pdgParticle = pdg->GetParticle(particle.pdgCode());
1303-
if (pdgParticle != nullptr)
1304-
charge = pdgParticle->Charge();
1305-
else
1306-
continue;
1307-
1308-
// Is it a charged particle?
1309-
if (std::abs(charge) < kMinCharge)
1310-
continue;
1311-
1312-
float phiPrime{track.phi()};
1313-
phiPrimeFunc(phiPrime, magField, charge);
1314-
1315-
const float pOrPt{v0Selections.usePinPhiSelection ? track.p() : track.pt()};
1316-
if (v0Selections.applyPhiCut) {
1317-
if (!passesPhiSelection(pOrPt, phiPrime))
1318-
continue;
1319-
}
1320-
1321-
const int16_t nclFound{track.tpcNClsFound()};
1322-
const int16_t nclPID{track.tpcNClsPID()};
1323-
const int16_t ncl = v0Selections.useNclsPID ? nclPID : nclFound;
1324-
if (v0Selections.applyNclSel && ncl < v0Selections.minNcl)
1325-
continue;
1326-
1327-
int indexEta{-999};
1328-
const float eta{track.eta()};
1329-
for (int i = 0; i < kNEtaHists; ++i) {
1330-
if (eta >= kLowEta[i] && eta < kHighEta[i]) {
1331-
indexEta = i;
1332-
break;
1333-
}
1334-
}
1335-
1336-
if (indexEta < kZeroInt || indexEta > kSevenInt)
1337-
continue;
1338-
1339-
registry.fill(HIST("NclVsPhip"), pOrPt, phiPrime, ncl);
1340-
registry.fill(HIST("NclVsEtaPID"), eta, ncl);
1341-
registry.fill(HIST("NclVsEtaPIDp"), eta, ncl);
1342-
1343-
bool isPrimary{false};
1344-
if (particle.isPhysicalPrimary())
1345-
isPrimary = true;
1346-
1347-
if (!isPrimary)
1348-
continue;
1349-
1350-
bool isPi{false};
1351-
bool isKa{false};
1352-
bool isPr{false};
1353-
1354-
if (particle.pdgCode() == PDG_t::kPiPlus || particle.pdgCode() == PDG_t::kPiMinus)
1355-
isPi = true;
1356-
else if (particle.pdgCode() == PDG_t::kKPlus || particle.pdgCode() == PDG_t::kKMinus)
1357-
isKa = true;
1358-
else if (particle.pdgCode() == PDG_t::kProton || particle.pdgCode() == PDG_t::kProtonBar)
1359-
isPr = true;
1360-
else
1361-
continue;
1362-
1363-
if (isPi && !isKa && !isPr)
1364-
registry.fill(HIST("PtPiVsCent"), track.pt(), centrality);
1365-
if (isKa && !isPi && !isPr)
1366-
registry.fill(HIST("PtKaVsCent"), track.pt(), centrality);
1367-
if (isPr && !isPi && !isKa)
1368-
registry.fill(HIST("PtPrVsCent"), track.pt(), centrality);
1369-
}
1370-
1371-
// Generated MC
1372-
for (const auto& particle : mcParticles) {
1373-
if (particle.eta() < v0Selections.minEtaDaughter || particle.eta() > v0Selections.maxEtaDaughter)
1374-
continue;
1375-
1376-
if (particle.pt() < v0Selections.minPt || particle.pt() > v0Selections.maxPt)
1377-
continue;
1378-
1379-
auto charge{0.};
1380-
// Get the MC particle
1381-
auto* pdgParticle = pdg->GetParticle(particle.pdgCode());
1382-
if (pdgParticle != nullptr)
1383-
charge = pdgParticle->Charge();
1384-
else
1385-
continue;
1386-
1387-
// Is it a charged particle?
1388-
if (std::abs(charge) < kMinCharge)
1389-
continue;
1390-
1391-
// Is it a primary particle?
1392-
bool isPrimary{true};
1393-
if (!particle.isPhysicalPrimary())
1394-
isPrimary = false;
1395-
1396-
if (isPrimary) {
1397-
if (particle.pdgCode() == PDG_t::kPiPlus || particle.pdgCode() == PDG_t::kPiMinus) // pion
1398-
registry.fill(HIST("PtPiVsCentMC"), particle.pt(), centrality);
1399-
else if (particle.pdgCode() == PDG_t::kKPlus || particle.pdgCode() == PDG_t::kKMinus) // kaon
1400-
registry.fill(HIST("PtKaVsCentMC"), particle.pt(), centrality);
1401-
else if (particle.pdgCode() == PDG_t::kProton || particle.pdgCode() == PDG_t::kProtonBar) // proton
1402-
registry.fill(HIST("PtPrVsCentMC"), particle.pt(), centrality);
1403-
else
1404-
continue;
1405-
}
1406-
}
1407-
} // Collisions
1408-
}
1409-
PROCESS_SWITCH(PiKpRAA, processMC, "Process MC closure", false);
14101138

14111139
void processSim(aod::McCollisions::iterator const& mccollision, soa::SmallGroups<ColEvSelsMC> const& collisions, BCsRun3 const& /*bcs*/, aod::FT0s const& /*ft0s*/, aod::McParticles const& mcParticles, TracksMC const& tracksMC)
14121140
{

0 commit comments

Comments
 (0)