Skip to content

Commit c39720a

Browse files
[PWGLF] Modifications in sigma0 tasks (AliceO2Group#9014)
Co-authored-by: ALICE Action Bot <[email protected]>
1 parent b8b7493 commit c39720a

File tree

3 files changed

+261
-173
lines changed

3 files changed

+261
-173
lines changed

PWGLF/DataModel/LFSigmaTables.h

Lines changed: 3 additions & 20 deletions
Original file line numberDiff line numberDiff line change
@@ -25,19 +25,6 @@
2525
// Creating output TTree for sigma analysis
2626
namespace o2::aod
2727
{
28-
DECLARE_SOA_TABLE(Sigma0Collisions, "AOD", "SIGMA0COLLISION", //! basic collision properties: position
29-
o2::soa::Index<>, collision::PosX, collision::PosY, collision::PosZ,
30-
cent::CentFT0M, cent::CentFT0A, cent::CentFT0C, cent::CentFV0A);
31-
32-
using Sigma0Collision = Sigma0Collisions::iterator;
33-
34-
namespace sigma0Core
35-
{
36-
//______________________________________________________
37-
// REGULAR COLUMNS FOR INDEXING
38-
// FOR DERIVED
39-
DECLARE_SOA_INDEX_COLUMN(Sigma0Collision, sigma0Collision); //!
40-
} // namespace sigma0Core
4128

4229
// for real data
4330
namespace sigma0Core
@@ -46,20 +33,16 @@ DECLARE_SOA_COLUMN(SigmapT, sigmapT, float);
4633
DECLARE_SOA_COLUMN(SigmaMass, sigmaMass, float);
4734
DECLARE_SOA_COLUMN(SigmaRapidity, sigmaRapidity, float);
4835
DECLARE_SOA_COLUMN(SigmaOPAngle, sigmaOPAngle, float);
49-
DECLARE_SOA_COLUMN(SigmaDeltaEta, sigmaDeltaEta, float);
50-
DECLARE_SOA_COLUMN(SigmaDeltaPhi, sigmaDeltaPhi, float);
36+
DECLARE_SOA_COLUMN(SigmaCentrality, sigmaCentrality, float);
37+
5138
} // namespace sigma0Core
5239

5340
DECLARE_SOA_TABLE(Sigma0Cores, "AOD", "SIGMA0CORES",
5441
sigma0Core::SigmapT,
5542
sigma0Core::SigmaMass,
5643
sigma0Core::SigmaRapidity,
5744
sigma0Core::SigmaOPAngle,
58-
sigma0Core::SigmaDeltaEta,
59-
sigma0Core::SigmaDeltaPhi);
60-
61-
DECLARE_SOA_TABLE(Sigma0CollRefs, "AOD", "SIGMA0COLLREF", //! optional table to refer back to a collision
62-
o2::soa::Index<>, sigma0Core::Sigma0CollisionId);
45+
sigma0Core::SigmaCentrality);
6346

6447
// For Photon extra info
6548
namespace sigmaPhotonExtra

PWGLF/TableProducer/Strangeness/sigma0builder.cxx

Lines changed: 135 additions & 34 deletions
Original file line numberDiff line numberDiff line change
@@ -62,8 +62,6 @@ using V0StandardDerivedDatas = soa::Join<aod::V0Cores, aod::V0CollRefs, aod::V0E
6262
struct sigma0builder {
6363
SliceCache cache;
6464

65-
Produces<aod::Sigma0Collision> sigma0Coll; // characterises collisions
66-
Produces<aod::Sigma0CollRefs> sigma0CollRefs; // characterises collisions
6765
Produces<aod::Sigma0Cores> sigma0cores; // save sigma0 candidates for analysis
6866
Produces<aod::SigmaPhotonExtras> sigmaPhotonExtras; // save sigma0 candidates for analysis
6967
Produces<aod::SigmaLambdaExtras> sigmaLambdaExtras; // save sigma0 candidates for analysis
@@ -104,6 +102,20 @@ struct sigma0builder {
104102
Configurable<float> Sigma0Window{"Sigma0Window", 0.1, "Mass window around expected (in GeV/c2)"};
105103
Configurable<float> SigmaMaxRap{"SigmaMaxRap", 0.8, "Max sigma0 rapidity"};
106104

105+
//// Extras:
106+
Configurable<bool> doPi0QA{"doPi0QA", true, "Flag to fill QA histos for pi0 rejection study."};
107+
Configurable<float> Pi0PhotonMinDCADauToPv{"Pi0PhotonMinDCADauToPv", 0.0, "Min DCA daughter To PV (cm)"};
108+
Configurable<float> Pi0PhotonMaxDCAV0Dau{"Pi0PhotonMaxDCAV0Dau", 3.5, "Max DCA V0 Daughters (cm)"};
109+
Configurable<int> Pi0PhotonMinTPCCrossedRows{"Pi0PhotonMinTPCCrossedRows", 0, "Min daughter TPC Crossed Rows"};
110+
Configurable<int> Pi0PhotonMaxTPCNSigmas{"Pi0PhotonMaxTPCNSigmas", 7, "Max TPC NSigmas for daughters"};
111+
Configurable<float> Pi0PhotonMaxEta{"Pi0PhotonMaxEta", 0.8, "Max photon rapidity"};
112+
Configurable<float> Pi0PhotonMinRadius{"Pi0PhotonMinRadius", 3.0, "Min photon conversion radius (cm)"};
113+
Configurable<float> Pi0PhotonMaxRadius{"Pi0PhotonMaxRadius", 115, "Max photon conversion radius (cm)"};
114+
Configurable<float> Pi0PhotonMaxQt{"Pi0PhotonMaxQt", 0.05, "Max photon qt value (AP plot) (GeV/c)"};
115+
Configurable<float> Pi0PhotonMaxAlpha{"Pi0PhotonMaxAlpha", 0.95, "Max photon alpha absolute value (AP plot)"};
116+
Configurable<float> Pi0PhotonMinV0cospa{"Pi0PhotonMinV0cospa", 0.80, "Min V0 CosPA"};
117+
Configurable<float> Pi0PhotonMaxMass{"Pi0PhotonMaxMass", 0.10, "Max photon mass (GeV/c^{2})"};
118+
107119
// Axis
108120
// base properties
109121
ConfigurableAxis axisPt{"axisPt", {VARIABLE_WIDTH, 0.0f, 0.1f, 0.2f, 0.3f, 0.4f, 0.5f, 0.6f, 0.7f, 0.8f, 0.9f, 1.0f, 1.1f, 1.2f, 1.3f, 1.4f, 1.5f, 1.6f, 1.7f, 1.8f, 1.9f, 2.0f, 2.2f, 2.4f, 2.6f, 2.8f, 3.0f, 3.2f, 3.4f, 3.6f, 3.8f, 4.0f, 4.4f, 4.8f, 5.2f, 5.6f, 6.0f, 6.5f, 7.0f, 7.5f, 8.0f, 9.0f, 10.0f, 11.0f, 12.0f, 13.0f, 14.0f, 15.0f, 17.0f, 19.0f, 21.0f, 23.0f, 25.0f, 30.0f, 35.0f, 40.0f, 50.0f}, "pt axis for analysis"};
@@ -176,6 +188,8 @@ struct sigma0builder {
176188
histos.add("GeneralQA/h2dMassGammaVsK0SAfterMassSel", "h2dMassGammaVsK0SAfterMassSel", kTH2D, {axisPhotonMass, {200, 0.4f, 0.6f}});
177189
histos.add("GeneralQA/h2dMassLambdaVsK0SAfterMassSel", "h2dMassLambdaVsK0SAfterMassSel", kTH2D, {axisLambdaMass, {200, 0.4f, 0.6f}});
178190
histos.add("GeneralQA/h2dMassGammaVsLambdaAfterMassSel", "h2dMassGammaVsLambdaAfterMassSel", kTH2D, {axisPhotonMass, axisLambdaMass});
191+
histos.add("GeneralQA/h2dPtVsMassPi0BeforeSel_Candidates", "h2dPtVsMassPi0BeforeSel_Candidates", kTH2D, {axisPt, {500, 0.08f, 0.18f}});
192+
histos.add("GeneralQA/h2dPtVsMassPi0AfterSel_Candidates", "h2dPtVsMassPi0AfterSel_Candidates", kTH2D, {axisPt, {500, 0.08f, 0.18f}});
179193

180194
// MC
181195
histos.add("MC/h2dPtVsCentrality_GammaBeforeSel", "h2dPtVsCentrality_GammaBeforeSel", kTH2D, {axisCentrality, axisPt});
@@ -210,9 +224,110 @@ struct sigma0builder {
210224
histos.add("MC/h2dFakeGammaTrueLambdaMatrix", "h2dFakeGammaTrueLambdaMatrix", kTHnSparseD, {{10001, -5000.5f, +5000.5f}, {10001, -5000.5f, +5000.5f}});
211225
histos.add("MC/h2dFakeDaughtersMatrix", "h2dFakeDaughtersMatrix", kTHnSparseD, {{10001, -5000.5f, +5000.5f}, {10001, -5000.5f, +5000.5f}});
212226

227+
// For Pi0 QA
228+
histos.add("MC/h2dPtVsMassPi0BeforeSel_SignalOnly", "h2dPtVsMassPi0BeforeSel_SignalOnly", kTH2D, {axisPt, {500, 0.08f, 0.18f}});
229+
histos.add("MC/h2dPtVsMassPi0AfterSel_SignalOnly", "h2dPtVsMassPi0AfterSel_SignalOnly", kTH2D, {axisPt, {500, 0.08f, 0.18f}});
230+
213231
histos.add("h3dMassSigmasBeforeSel", "h3dMassSigmasBeforeSel", kTH3F, {axisCentrality, axisPt, axisSigmaMass});
214232
histos.add("h3dMassSigmasAfterSel", "h3dMassSigmasAfterSel", kTH3F, {axisCentrality, axisPt, axisSigmaMass});
215233
}
234+
template <typename TV0Object>
235+
void runPi0QA(TV0Object const& gamma1, TV0Object const& gamma2)
236+
{
237+
238+
// Check if both V0s are made of the same tracks
239+
if (gamma1.posTrackExtraId() == gamma2.posTrackExtraId() ||
240+
gamma1.negTrackExtraId() == gamma2.negTrackExtraId() ||
241+
gamma1.posTrackExtraId() == gamma2.negTrackExtraId() ||
242+
gamma1.negTrackExtraId() == gamma2.posTrackExtraId()) {
243+
return;
244+
}
245+
246+
// Calculate pi0 properties
247+
std::array<float, 3> pVecGamma1{gamma1.px(), gamma1.py(), gamma1.pz()};
248+
std::array<float, 3> pVecGamma2{gamma2.px(), gamma2.py(), gamma2.pz()};
249+
std::array arrpi0{pVecGamma1, pVecGamma2};
250+
float pi0Mass = RecoDecay::m(arrpi0, std::array{o2::constants::physics::MassPhoton, o2::constants::physics::MassPhoton});
251+
float pi0Pt = RecoDecay::pt(std::array{gamma1.px() + gamma2.px(), gamma1.py() + gamma2.py()});
252+
float pi0Y = RecoDecay::y(std::array{gamma1.px() + gamma2.px(), gamma1.py() + gamma2.py(), gamma1.pz() + gamma2.pz()}, o2::constants::physics::MassPi0);
253+
254+
// MC-specific variables
255+
bool fIsPi0 = false, fIsMC = false;
256+
257+
// Check if MC data and populate fIsMC, fIsPi0
258+
if constexpr (requires { gamma1.pdgCode(); gamma2.pdgCode(); }) {
259+
fIsMC = true;
260+
if (gamma1.pdgCode() == 22 && gamma2.pdgCode() == 22 &&
261+
gamma1.pdgCodeMother() == 111 && gamma2.pdgCodeMother() == 111 &&
262+
gamma1.motherMCPartId() == gamma2.motherMCPartId()) {
263+
fIsPi0 = true;
264+
histos.fill(HIST("MC/h2dPtVsMassPi0BeforeSel_SignalOnly"), pi0Pt, pi0Mass);
265+
}
266+
} else {
267+
histos.fill(HIST("GeneralQA/h2dPtVsMassPi0BeforeSel_Candidates"), pi0Pt, pi0Mass);
268+
}
269+
270+
// Photon-specific selections
271+
auto posTrackGamma1 = gamma1.template posTrackExtra_as<dauTracks>();
272+
auto negTrackGamma1 = gamma1.template negTrackExtra_as<dauTracks>();
273+
auto posTrackGamma2 = gamma2.template posTrackExtra_as<dauTracks>();
274+
auto negTrackGamma2 = gamma2.template negTrackExtra_as<dauTracks>();
275+
276+
// Gamma1 Selection
277+
bool passedTPCGamma1 = (posTrackGamma1.tpcNSigmaEl() == -999.f || TMath::Abs(posTrackGamma1.tpcNSigmaEl()) < Pi0PhotonMaxTPCNSigmas) &&
278+
(negTrackGamma1.tpcNSigmaEl() == -999.f || TMath::Abs(negTrackGamma1.tpcNSigmaEl()) < Pi0PhotonMaxTPCNSigmas);
279+
280+
if (TMath::Abs(gamma1.mGamma()) > Pi0PhotonMaxMass ||
281+
gamma1.qtarm() >= Pi0PhotonMaxQt ||
282+
TMath::Abs(gamma1.alpha()) >= Pi0PhotonMaxAlpha ||
283+
TMath::Abs(gamma1.dcapostopv()) < Pi0PhotonMinDCADauToPv ||
284+
TMath::Abs(gamma1.dcanegtopv()) < Pi0PhotonMinDCADauToPv ||
285+
TMath::Abs(gamma1.dcaV0daughters()) > Pi0PhotonMaxDCAV0Dau ||
286+
TMath::Abs(gamma1.negativeeta()) >= Pi0PhotonMaxEta ||
287+
TMath::Abs(gamma1.positiveeta()) >= Pi0PhotonMaxEta ||
288+
gamma1.v0cosPA() <= Pi0PhotonMinV0cospa ||
289+
gamma1.v0radius() <= Pi0PhotonMinRadius ||
290+
gamma1.v0radius() >= Pi0PhotonMaxRadius ||
291+
posTrackGamma1.tpcCrossedRows() < Pi0PhotonMinTPCCrossedRows ||
292+
negTrackGamma1.tpcCrossedRows() < Pi0PhotonMinTPCCrossedRows ||
293+
!passedTPCGamma1) {
294+
return;
295+
}
296+
297+
// Gamma2 Selection
298+
bool passedTPCGamma2 = (posTrackGamma2.tpcNSigmaEl() == -999.f || TMath::Abs(posTrackGamma2.tpcNSigmaEl()) < Pi0PhotonMaxTPCNSigmas) &&
299+
(negTrackGamma2.tpcNSigmaEl() == -999.f || TMath::Abs(negTrackGamma2.tpcNSigmaEl()) < Pi0PhotonMaxTPCNSigmas);
300+
301+
if (TMath::Abs(gamma2.mGamma()) > Pi0PhotonMaxMass ||
302+
gamma2.qtarm() >= Pi0PhotonMaxQt ||
303+
TMath::Abs(gamma2.alpha()) >= Pi0PhotonMaxAlpha ||
304+
TMath::Abs(gamma2.dcapostopv()) < Pi0PhotonMinDCADauToPv ||
305+
TMath::Abs(gamma2.dcanegtopv()) < Pi0PhotonMinDCADauToPv ||
306+
TMath::Abs(gamma2.dcaV0daughters()) > Pi0PhotonMaxDCAV0Dau ||
307+
TMath::Abs(gamma2.negativeeta()) >= Pi0PhotonMaxEta ||
308+
TMath::Abs(gamma2.positiveeta()) >= Pi0PhotonMaxEta ||
309+
gamma2.v0cosPA() <= Pi0PhotonMinV0cospa ||
310+
gamma2.v0radius() <= Pi0PhotonMinRadius ||
311+
gamma2.v0radius() >= Pi0PhotonMaxRadius ||
312+
posTrackGamma2.tpcCrossedRows() < Pi0PhotonMinTPCCrossedRows ||
313+
negTrackGamma2.tpcCrossedRows() < Pi0PhotonMinTPCCrossedRows ||
314+
!passedTPCGamma2) {
315+
return;
316+
}
317+
318+
// Pi0-specific selections:
319+
if (TMath::Abs(pi0Y) > 0.5) {
320+
return;
321+
}
322+
323+
// Fill histograms
324+
if (fIsMC) {
325+
if (fIsPi0)
326+
histos.fill(HIST("MC/h2dPtVsMassPi0AfterSel_SignalOnly"), pi0Pt, pi0Mass);
327+
} else {
328+
histos.fill(HIST("GeneralQA/h2dPtVsMassPi0AfterSel_Candidates"), pi0Pt, pi0Mass);
329+
}
330+
}
216331

217332
// Process sigma candidate and store properties in object
218333
template <typename TV0Object>
@@ -350,19 +465,10 @@ struct sigma0builder {
350465
histos.fill(HIST("hCandidateBuilderSelection"), 13.);
351466
return true;
352467
}
353-
// Helper struct to pass v0 information
354-
struct {
355-
float mass;
356-
float pT;
357-
float Rapidity;
358-
float OPAngle;
359-
float DeltaEta;
360-
float DeltaPhi;
361-
} sigmaCandidate;
362468

363469
// Fill tables with reconstructed sigma0 candidate
364-
template <typename TV0Object>
365-
void fillTables(TV0Object const& lambda, TV0Object const& gamma)
470+
template <typename TV0Object, typename TCollision>
471+
void fillTables(TV0Object const& lambda, TV0Object const& gamma, TCollision const& coll)
366472
{
367473

368474
float GammaBDTScore = -1;
@@ -464,23 +570,15 @@ struct sigma0builder {
464570
TVector3 v1(gamma.px(), gamma.py(), gamma.pz());
465571
TVector3 v2(lambda.px(), lambda.py(), lambda.pz());
466572

467-
sigmaCandidate.mass = RecoDecay::m(arrMom, std::array{o2::constants::physics::MassPhoton, o2::constants::physics::MassLambda0});
468-
sigmaCandidate.pT = RecoDecay::pt(array{gamma.px() + lambda.px(), gamma.py() + lambda.py()});
469-
sigmaCandidate.Rapidity = RecoDecay::y(std::array{gamma.px() + lambda.px(), gamma.py() + lambda.py(), gamma.pz() + lambda.pz()}, o2::constants::physics::MassSigma0);
470-
sigmaCandidate.OPAngle = v1.Angle(v2);
471-
sigmaCandidate.DeltaEta = fLambdaEta - fPhotonEta;
472-
sigmaCandidate.DeltaPhi = fLambdaPhi - fPhotonPhi;
473-
474573
// Sigma related
475-
float fSigmapT = sigmaCandidate.pT;
476-
float fSigmaMass = sigmaCandidate.mass;
477-
float fSigmaRap = sigmaCandidate.Rapidity;
478-
float fSigmaOPAngle = sigmaCandidate.OPAngle;
479-
float fSigmaDeltaEta = sigmaCandidate.DeltaEta;
480-
float fSigmaDeltaPhi = sigmaCandidate.DeltaPhi;
574+
float fSigmapT = RecoDecay::pt(array{gamma.px() + lambda.px(), gamma.py() + lambda.py()});
575+
float fSigmaMass = RecoDecay::m(arrMom, std::array{o2::constants::physics::MassPhoton, o2::constants::physics::MassLambda0});
576+
float fSigmaRap = RecoDecay::y(std::array{gamma.px() + lambda.px(), gamma.py() + lambda.py(), gamma.pz() + lambda.pz()}, o2::constants::physics::MassSigma0);
577+
float fSigmaOPAngle = v1.Angle(v2);
578+
float fSigmaCentrality = coll.centFT0C();
481579

482580
// Filling TTree for ML analysis
483-
sigma0cores(fSigmapT, fSigmaMass, fSigmaRap, fSigmaOPAngle, fSigmaDeltaEta, fSigmaDeltaPhi);
581+
sigma0cores(fSigmapT, fSigmaMass, fSigmaRap, fSigmaOPAngle, fSigmaCentrality);
484582

485583
sigmaPhotonExtras(fPhotonPt, fPhotonMass, fPhotonQt, fPhotonAlpha, fPhotonRadius,
486584
fPhotonCosPA, fPhotonDCADau, fPhotonDCANegPV, fPhotonDCAPosPV, fPhotonZconv,
@@ -501,7 +599,7 @@ struct sigma0builder {
501599
fLambdaV0Type, LambdaBDTScore, AntiLambdaBDTScore);
502600
}
503601

504-
void processMonteCarlo(soa::Join<aod::StraCollisions, aod::StraCents> const& collisions, V0DerivedMCDatas const& V0s)
602+
void processMonteCarlo(soa::Join<aod::StraCollisions, aod::StraCents> const& collisions, V0DerivedMCDatas const& V0s, dauTracks const&)
505603
{
506604
for (const auto& coll : collisions) {
507605
// Do analysis with collision-grouped V0s, retain full collision information
@@ -549,6 +647,9 @@ struct sigma0builder {
549647
}
550648

551649
for (auto& lambda : V0Table_thisCollision) { // selecting lambdas from Sigma0
650+
if (doPi0QA) // Pi0 QA study
651+
runPi0QA(gamma, lambda);
652+
552653
// Sigma0 candidate properties
553654
std::array<float, 3> pVecPhotons{gamma.px(), gamma.py(), gamma.pz()};
554655
std::array<float, 3> pVecLambda{lambda.px(), lambda.py(), lambda.pz()};
@@ -607,11 +708,14 @@ struct sigma0builder {
607708
auto V0Table_thisCollision = V0s.sliceBy(perCollisionSTDDerived, collIdx);
608709

609710
histos.fill(HIST("hEventCentrality"), coll.centFT0C());
610-
sigma0Coll(coll.posX(), coll.posY(), coll.posZ(), coll.centFT0M(), coll.centFT0A(), coll.centFT0C(), coll.centFV0A());
611711

612712
// V0 table sliced
613713
for (auto& gamma : V0Table_thisCollision) { // selecting photons from Sigma0
614714
for (auto& lambda : V0Table_thisCollision) { // selecting lambdas from Sigma0
715+
if (doPi0QA) // Pi0 QA study
716+
runPi0QA(gamma, lambda);
717+
718+
// Sigma0 candidate properties
615719
std::array<float, 3> pVecPhotons{gamma.px(), gamma.py(), gamma.pz()};
616720
std::array<float, 3> pVecLambda{lambda.px(), lambda.py(), lambda.pz()};
617721
auto arrMom = std::array{pVecPhotons, pVecLambda};
@@ -625,8 +729,7 @@ struct sigma0builder {
625729

626730
histos.fill(HIST("h3dMassSigmasAfterSel"), coll.centFT0C(), SigmapT, SigmaMass);
627731

628-
sigma0CollRefs(collIdx);
629-
fillTables(lambda, gamma); // filling tables with accepted candidates
732+
fillTables(lambda, gamma, coll); // filling tables with accepted candidates
630733

631734
nSigmaCandidates++;
632735
if (nSigmaCandidates % 5000 == 0) {
@@ -645,7 +748,6 @@ struct sigma0builder {
645748
auto V0Table_thisCollision = V0s.sliceBy(perCollisionMLDerived, collIdx);
646749

647750
histos.fill(HIST("hEventCentrality"), coll.centFT0C());
648-
sigma0Coll(coll.posX(), coll.posY(), coll.posZ(), coll.centFT0M(), coll.centFT0A(), coll.centFT0C(), coll.centFV0A());
649751

650752
// V0 table sliced
651753
for (auto& gamma : V0Table_thisCollision) { // selecting photons from Sigma0
@@ -657,8 +759,7 @@ struct sigma0builder {
657759
if (nSigmaCandidates % 5000 == 0) {
658760
LOG(info) << "Sigma0 Candidates built: " << nSigmaCandidates;
659761
}
660-
sigma0CollRefs(collIdx);
661-
fillTables(lambda, gamma); // filling tables with accepted candidates
762+
fillTables(lambda, gamma, coll); // filling tables with accepted candidates
662763
}
663764
}
664765
}

0 commit comments

Comments
 (0)