Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 9 additions & 1 deletion PWGLF/DataModel/LFAntinCexTables.h
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,12 @@ DECLARE_SOA_COLUMN(McAngleDeg, mcAngleDeg, float);
DECLARE_SOA_COLUMN(McVtxX, mcVtxX, float);
DECLARE_SOA_COLUMN(McVtxY, mcVtxY, float);
DECLARE_SOA_COLUMN(McVtxZ, mcVtxZ, float);
DECLARE_SOA_COLUMN(VtxNAll, vtxNAll, int16_t);
DECLARE_SOA_COLUMN(VtxNCh, vtxNCh, int16_t);
DECLARE_SOA_COLUMN(VtxNNeut, vtxNNeut, int16_t);
DECLARE_SOA_COLUMN(VtxNPi0, vtxNPi0, int16_t);
DECLARE_SOA_COLUMN(VtxNGamma, vtxNGamma, int16_t);
DECLARE_SOA_COLUMN(VtxNN, vtxNN, int16_t);

// Tracks (pair, fitter)
DECLARE_SOA_COLUMN(TrkPairP, trkPairP, float);
Expand Down Expand Up @@ -85,6 +91,7 @@ DECLARE_SOA_COLUMN(AntipTrkNClsIts, antipTrkNClsIts, int16_t);
DECLARE_SOA_COLUMN(SelMask, selMask, uint32_t);

DECLARE_SOA_COLUMN(PairPointingAngleDeg, pairPointingAngleDeg, float);
DECLARE_SOA_COLUMN(PvsvAngleZDeg, pvsvAngleZDeg, float);
DECLARE_SOA_COLUMN(PairPBalance, pairPBalance, float);
DECLARE_SOA_COLUMN(PairPtBalance, pairPtBalance, float);
DECLARE_SOA_COLUMN(PairQ, pairQ, float);
Expand Down Expand Up @@ -121,6 +128,7 @@ DECLARE_SOA_TABLE(AntinCexPairs, "AOD", "ANTINCEX",
antin_cex::MotherPdg, antin_cex::ColId, antin_cex::PId, antin_cex::AntipId,
antin_cex::McPairP, antin_cex::McPairPt, antin_cex::McPairPz,
antin_cex::McDplane, antin_cex::McAngleDeg, antin_cex::McVtxX, antin_cex::McVtxY, antin_cex::McVtxZ,
antin_cex::VtxNAll, antin_cex::VtxNCh, antin_cex::VtxNNeut, antin_cex::VtxNPi0, antin_cex::VtxNGamma, antin_cex::VtxNN,
antin_cex::TrkPairP, antin_cex::TrkPairPt, antin_cex::TrkPairPz, antin_cex::TrkAngleDeg,
antin_cex::TrkVtxfitDcaPair, antin_cex::TrkVtxfitR, antin_cex::TrkVtxfitDistToPv,
antin_cex::TrkVtxfitSecVtxX, antin_cex::TrkVtxfitSecVtxY, antin_cex::TrkVtxfitSecVtxZ,
Expand All @@ -129,7 +137,7 @@ DECLARE_SOA_TABLE(AntinCexPairs, "AOD", "ANTINCEX",
antin_cex::PTrkP, antin_cex::PTrkPx, antin_cex::PTrkPy, antin_cex::PTrkPz, antin_cex::PTrkEta, antin_cex::PTrkTpcSignal, antin_cex::PTrkNClsIts,
antin_cex::AntipTrkP, antin_cex::AntipTrkPx, antin_cex::AntipTrkPy, antin_cex::AntipTrkPz, antin_cex::AntipTrkEta, antin_cex::AntipTrkTpcSignal, antin_cex::AntipTrkNClsIts,
antin_cex::SelMask,
antin_cex::PairPointingAngleDeg, antin_cex::PairPBalance, antin_cex::PairPtBalance, antin_cex::PairQ,
antin_cex::PairPointingAngleDeg, antin_cex::PvsvAngleZDeg, antin_cex::PairPBalance, antin_cex::PairPtBalance, antin_cex::PairQ,
antin_cex::DPairP, antin_cex::DPairPt, antin_cex::DPairPz, antin_cex::DOpenAngle,
antin_cex::SVNearestLayerId, antin_cex::SVDeltaRToLayer,
antin_cex::PTrkItsHitMap, antin_cex::APTrkItsHitMap, antin_cex::PLayersOk, antin_cex::APLayersOk,
Expand Down
151 changes: 150 additions & 1 deletion PWGLF/TableProducer/Nuspex/nucleiAntineutronCex.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,7 @@
#include <Framework/runDataProcessing.h>
#include <ReconstructionDataFormats/TrackParametrization.h>

#include <TDatabasePDG.h>
#include <TMCProcess.h>
#include <TMath.h>
#include <TPDGCode.h>
Expand All @@ -53,7 +54,7 @@ struct NucleiAntineutronCex {
static constexpr double kIts2MaxR = 48.0; // ITS2 max radius [cm]
static constexpr double kIts2MaxVz = 39.0; // ITS2 max |vz| [cm]
static constexpr double kAccMaxEta = 1.2; // acceptance |eta|
static constexpr double kAccMaxVz = 5.3; // acceptance |vz| [cm]
static constexpr double kAccMaxVz = 10.0; // acceptance |vz| [cm]
static constexpr double kStrictEta = 0.9; // tighter eta cut
static constexpr double kInitDplane = 10.0; // init dplane
static constexpr double kHuge = 1e9; // fallback for bad denom
Expand Down Expand Up @@ -105,6 +106,25 @@ struct NucleiAntineutronCex {
// test (MC)
histos.add("antip_test", "Secondary antiprotons;|p| (GeV/c);Entries", kTH1F, {{100, 0., 10.}});

// Process enum breakdown (secondary antiproton that anchors the SV)
histos.add("hProcEnumAP_CEX", "procEnum of secondary #bar{p} (CEX);procEnum;Entries", kTH1I, {{100, -0.5, 99.5}});
histos.add("hProcEnumAP_BG", "procEnum of secondary #bar{p} (BG);procEnum;Entries", kTH1I, {{100, -0.5, 99.5}});

// Multiplicity/composition at the SV (MC truth, for FINAL selected candidates)
histos.add("hVtxNAll_CEX", "N(all) secondaries at SV (CEX);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNAll_BG", "N(all) secondaries at SV (BG);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNCh_CEX", "N(charged) secondaries at SV (CEX);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNCh_BG", "N(charged) secondaries at SV (BG);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNNeut_CEX", "N(neutral) secondaries at SV (CEX);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNNeut_BG", "N(neutral) secondaries at SV (BG);N;Entries", kTH1I, {{60, -0.5, 59.5}});

histos.add("hVtxNPi0_CEX", "N(#pi^{0}) at SV (CEX);N;Entries", kTH1I, {{40, -0.5, 39.5}});
histos.add("hVtxNPi0_BG", "N(#pi^{0}) at SV (BG);N;Entries", kTH1I, {{40, -0.5, 39.5}});
histos.add("hVtxNGamma_CEX", "N(#gamma) at SV (CEX);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNGamma_BG", "N(#gamma) at SV (BG);N;Entries", kTH1I, {{60, -0.5, 59.5}});
histos.add("hVtxNN_CEX", "N(n) at SV (CEX);N;Entries", kTH1I, {{40, -0.5, 39.5}});
histos.add("hVtxNN_BG", "N(n) at SV (BG);N;Entries", kTH1I, {{40, -0.5, 39.5}});

// CEX pair from antineutron (MC)
histos.add("cexPairMcP", "CEX pair total momentum;|p| (GeV/c);Entries", kTH1F, {{100, 0., 10.}});
histos.add("cexPairMcPt", "CEX pair p_{T};p_{T} (GeV/c);Entries", kTH1F, {{100, 0., 10.}});
Expand All @@ -130,6 +150,13 @@ struct NucleiAntineutronCex {
histos.add("cexbg_pairmc_vtxz", "Background secondary vertex Z;Z (cm);Entries", kTH1F, {{200, -60., 60.}});
histos.add("cexbg_pairmc_pITScuts", "Background momentum (ITS cuts);|p| (GeV/c);Entries", kTH1F, {{100, 0., 10.}});

histos.add("hDeltaP_CEX", "|p_{mother}-Σp_{SV}| (CEX);Δp (GeV/c);Entries", kTH1F, {{200, 0., 10.}});
histos.add("hDeltaP_BG", "|p_{mother}-Σp_{SV}| (BG);Δp (GeV/c);Entries", kTH1F, {{200, 0., 10.}});

// Mother IB hits
histos.add("hMotherNHitIB_CEX", "Mother #bar{p} IB hit layers (L0-L2) (CEX);N_{IB layers};Entries", kTH1I, {{4, -0.5, 3.5}});
histos.add("hMotherNHitIB_BG", "Mother #bar{p} IB hit layers (L0-L2) (BG);N_{IB layers};Entries", kTH1I, {{4, -0.5, 3.5}});

// CEX pair from antineutron (TRK)
histos.add("cex_pairtrk_angle", "Pair opening angle (tracks);Angle (°);Entries", kTH1F, {{180, 0., 180.}});
histos.add("cexPairTrkP", "Pair momentum (tracks);|p| (GeV/c);Entries", kTH1F, {{120, 0., 12.}});
Expand Down Expand Up @@ -265,21 +292,27 @@ struct NucleiAntineutronCex {
// Primary mother
bool hasPrimaryMotherAntip = false;
double motherPt = 0.0;
double motherPx = 0.0;
double motherPy = 0.0;
double motherPz = 0.0;
double motherVz = 0.0;
double motherP = 0.0;
double motherEta = 0.0;
int motherPdg = 0;
int motherId = -1;

for (const auto& mother : particle.mothers_as<aod::McParticles>()) {
if (mother.isPhysicalPrimary()) {
hasPrimaryMotherAntip = true;
motherPt = mother.pt();
motherPx = mother.px();
motherPy = mother.py();
motherPz = mother.pz();
motherVz = mother.vz();
motherP = mother.p();
motherEta = mother.eta();
motherPdg = mother.pdgCode();
motherId = mother.globalIndex();
break;
}
}
Expand Down Expand Up @@ -530,6 +563,10 @@ struct NucleiAntineutronCex {
int8_t antipTrkItsPidValid = 0;
float antipTrkTgl = 0.f;

bool motherHasTrack = false;
uint16_t motherItsMap = 0;
int motherNHitIB = 0; // number of hits in IB (L0-L2)

o2::aod::ITSResponse itsResponse;

for (const auto& track : tracks) {
Expand Down Expand Up @@ -557,6 +594,12 @@ struct NucleiAntineutronCex {
int nITS = track.itsNCls();
bool layerCondition = (!hitIB) && hitOuter && (nITS >= kMinItsHits);

if (mc.globalIndex() == motherId) {
motherHasTrack = true;
motherItsMap = static_cast<uint16_t>(track.itsClusterMap());
motherNHitIB = static_cast<int>(hitL0) + static_cast<int>(hitL1) + static_cast<int>(hitL2);
}

if (mc.globalIndex() == antipId) {
antipTrkP = track.p();
antipTrkPx = track.px();
Expand Down Expand Up @@ -695,6 +738,9 @@ struct NucleiAntineutronCex {
const TVector3 pv2sv(secX - pvtxX, secY - pvtxY, secZ - pvtxZ);
const double pairPointingAngleDeg = pv2sv.Angle(total_trk_pVec) * Rad2Deg;

const TVector3 zAxis(0., 0., 1.);
const double pvsvAngleZDeg = pv2sv.Angle(zAxis) * Rad2Deg;

const double pP = pVecProton_trk.Mag();
const double pAP = AntipVecProton_trk.Mag();
const double ptP = pVecProton_trk.Pt();
Expand Down Expand Up @@ -757,11 +803,106 @@ struct NucleiAntineutronCex {

const bool isCex = (motherPdg == -kNeutron);

// Nature of the process
if (isCex) {
histos.fill(HIST("hProcEnumAP_CEX"), static_cast<int>(procEnum));
} else {
histos.fill(HIST("hProcEnumAP_BG"), static_cast<int>(procEnum));
}

// Count material secondaries produced at the same SV as the selected secondary antiproton.
int vtxNAll = 0;
int vtxNCh = 0;
int vtxNNeut = 0;
int vtxNPi0 = 0;
int vtxNGamma = 0;
int vtxNN = 0;
double sumPx_vtx = 0.0;
double sumPy_vtx = 0.0;
double sumPz_vtx = 0.0;
auto* pdgDB = TDatabasePDG::Instance();

for (const auto& particle5 : mcPartsThis) {
if (particle5.mcCollisionId() != colId) {
continue;
}
// Same SV (use the SV of the selected secondary antiproton)
if (std::abs(particle5.vx() - antipVx) >= kVtxTol || std::abs(particle5.vy() - antipVy) >= kVtxTol || std::abs(particle5.vz() - antipVz) >= kVtxTol) {
continue;
}
const auto proc5Enum = particle5.getProcess();
const bool isSecondaryFromMaterial5 =
(!particle5.producedByGenerator()) && (proc5Enum == kPHadronic || proc5Enum == kPHInhelastic);
if (!isSecondaryFromMaterial5) {
continue;
}
++vtxNAll;
sumPx_vtx += particle5.px();
sumPy_vtx += particle5.py();
sumPz_vtx += particle5.pz();
const int pdg = particle5.pdgCode();
if (pdg == kPi0) {
++vtxNPi0;
}
if (pdg == kGamma) {
++vtxNGamma;
}
if (pdg == kNeutron) {
++vtxNN;
}
// Charged vs neutral via PDG database (Charge() is in units of e/3)
double q = 0.0;
if (auto* part = pdgDB->GetParticle(pdg)) {
q = part->Charge() / 3.0;
}
if (std::abs(q) > 0.0) {
++vtxNCh;
} else {
++vtxNN;
}
}

// Fill histos (final selected candidates only)
if (isCex) {
histos.fill(HIST("hVtxNAll_CEX"), vtxNAll);
histos.fill(HIST("hVtxNCh_CEX"), vtxNCh);
histos.fill(HIST("hVtxNNeut_CEX"), vtxNNeut);
histos.fill(HIST("hVtxNPi0_CEX"), vtxNPi0);
histos.fill(HIST("hVtxNGamma_CEX"), vtxNGamma);
histos.fill(HIST("hVtxNN_CEX"), vtxNN);
} else {
histos.fill(HIST("hVtxNAll_BG"), vtxNAll);
histos.fill(HIST("hVtxNCh_BG"), vtxNCh);
histos.fill(HIST("hVtxNNeut_BG"), vtxNNeut);
histos.fill(HIST("hVtxNPi0_BG"), vtxNPi0);
histos.fill(HIST("hVtxNGamma_BG"), vtxNGamma);
histos.fill(HIST("hVtxNN_BG"), vtxNN);
}

const float vtxfitDX = secX - antipVx;
const float vtxfitDY = secY - antipVy;
const float vtxfitDZ = secZ - antipVz;
const float vtxfitD3D = std::sqrt(vtxfitDX * vtxfitDX + vtxfitDY * vtxfitDY + vtxfitDZ * vtxfitDZ);

const double dPx = motherPx - sumPx_vtx;
const double dPy = motherPy - sumPy_vtx;
const double dPz = motherPz - sumPz_vtx;
const double deltaP = std::sqrt(dPx * dPx + dPy * dPy + dPz * dPz);

if (isCex) {
histos.fill(HIST("hDeltaP_CEX"), deltaP);
} else {
histos.fill(HIST("hDeltaP_BG"), deltaP);
}

if (motherHasTrack) {
if (isCex) {
histos.fill(HIST("hMotherNHitIB_CEX"), motherNHitIB);
} else {
histos.fill(HIST("hMotherNHitIB_BG"), motherNHitIB);
}
}

const uint32_t selMask = 0u;

outPairs(
Expand All @@ -780,6 +921,13 @@ struct NucleiAntineutronCex {
antipVy,
antipVz,

static_cast<int16_t>(vtxNAll),
static_cast<int16_t>(vtxNCh),
static_cast<int16_t>(vtxNNeut),
static_cast<int16_t>(vtxNPi0),
static_cast<int16_t>(vtxNGamma),
static_cast<int16_t>(vtxNN),

cexPairTrkP,
cexPairTrkPt,
cexPairTrkPz,
Expand Down Expand Up @@ -818,6 +966,7 @@ struct NucleiAntineutronCex {
selMask,

pairPointingAngleDeg,
pvsvAngleZDeg,
pairPBalance,
pairPtBalance,
pairQ,
Expand Down
Loading