Skip to content
Closed
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
18 changes: 8 additions & 10 deletions PWGLF/DataModel/LFAntinCexTables.h
Original file line number Diff line number Diff line change
Expand Up @@ -27,12 +27,12 @@
namespace antin_cex
{
// Metadata
DECLARE_SOA_COLUMN(IsCex, isCex, bool); // 1=CEX (from antin), 0=BG
DECLARE_SOA_COLUMN(MotherPdg, motherPdg, int32_t); // mother PDG
DECLARE_SOA_COLUMN(MotherNHitIB, motherNHitIB, int); // mother IB Hits
DECLARE_SOA_COLUMN(ColId, colId, int32_t); // mcCollisionId
DECLARE_SOA_COLUMN(PId, pId, int32_t); // proton MC id
DECLARE_SOA_COLUMN(AntipId, antipId, int32_t); // antiproton MC id
DECLARE_SOA_COLUMN(IsCex, isCex, bool); // 1=CEX (from antin), 0=BG
DECLARE_SOA_COLUMN(MotherPdg, motherPdg, int32_t); // mother PDG
DECLARE_SOA_COLUMN(MotherP, motherP, float); // mother P
DECLARE_SOA_COLUMN(ColId, colId, int32_t); // mcCollisionId
DECLARE_SOA_COLUMN(PId, pId, int32_t); // proton MC id
DECLARE_SOA_COLUMN(AntipId, antipId, int32_t); // antiproton MC id

// MC (pair)
DECLARE_SOA_COLUMN(McPairP, mcPairP, float);
Expand Down Expand Up @@ -104,13 +104,11 @@
DECLARE_SOA_COLUMN(DPairPz, dPairPz, float);
DECLARE_SOA_COLUMN(DOpenAngle, dOpenAngle, float);

DECLARE_SOA_COLUMN(SVNearestLayerId, svNearestLayerId, int16_t);

Check failure on line 107 in PWGLF/DataModel/LFAntinCexTables.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column]

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.
DECLARE_SOA_COLUMN(SVDeltaRToLayer, svDeltaRToLayer, float);

Check failure on line 108 in PWGLF/DataModel/LFAntinCexTables.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column]

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.

DECLARE_SOA_COLUMN(PTrkItsHitMap, pTrkItsHitMap, uint16_t);
DECLARE_SOA_COLUMN(APTrkItsHitMap, apTrkItsHitMap, uint16_t);

Check failure on line 111 in PWGLF/DataModel/LFAntinCexTables.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column]

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.
DECLARE_SOA_COLUMN(PLayersOk, pLayersOk, int8_t);
DECLARE_SOA_COLUMN(APLayersOk, apLayersOk, int8_t);

DECLARE_SOA_COLUMN(PVtxZ, pVtxZ, float);

Expand All @@ -128,7 +126,7 @@
// Table
DECLARE_SOA_TABLE(AntinCexPairs, "AOD", "ANTINCEX",
antin_cex::IsCex,
antin_cex::MotherPdg, antin_cex::MotherNHitIB, antin_cex::ColId, antin_cex::PId, antin_cex::AntipId,
antin_cex::MotherPdg, antin_cex::MotherP, 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,
Expand All @@ -143,7 +141,7 @@
antin_cex::PairPointingAngleDeg, antin_cex::PvsvThetaDeg, antin_cex::PvsvPhiDeg, 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,
antin_cex::PTrkItsHitMap, antin_cex::APTrkItsHitMap,
antin_cex::PVtxZ,
antin_cex::PTrkItsNSigmaPr, antin_cex::PTrkItsPidValid, antin_cex::PTrkTgl,
antin_cex::AntipTrkItsNSigmaPr, antin_cex::AntipTrkItsPidValid, antin_cex::AntipTrkTgl);
Expand Down
179 changes: 47 additions & 132 deletions PWGLF/TableProducer/Nuspex/nucleiAntineutronCex.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -15,9 +15,8 @@
/// \author Fabiola Lugo
///

#include "PWGLF/DataModel/LFAntinCexTables.h"

#include "Common/DataModel/PIDResponseITS.h"
#include <PWGLF/DataModel/LFAntinCexTables.h>
#include <Common/DataModel/PIDResponseITS.h>

#include <CommonConstants/MathConstants.h>
#include <DCAFitter/DCAFitterN.h>
Expand All @@ -33,7 +32,6 @@
#include <ReconstructionDataFormats/Track.h>
#include <ReconstructionDataFormats/TrackParametrization.h>

#include <TDatabasePDG.h>
#include <TMCProcess.h>
#include <TPDGCode.h>
#include <TVector3.h>
Expand All @@ -51,15 +49,15 @@ using o2::constants::math::Rad2Deg;
struct NucleiAntineutronCex {
// Slicing per colision
Preslice<aod::McParticles> perMcByColl = aod::mcparticle::mcCollisionId;
// Check available tables in the AOD, specifically TracksIU, TracksCovIU

using TracksWCovMc = soa::Join<aod::TracksIU, aod::TracksExtra, aod::McTrackLabels, o2::aod::TracksCovIU>;

// === Cut values ===
static constexpr double kIts2MinR = 4.5; // ITS2 min radius (exluding IB) [cm]
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 = 10.0; // acceptance |vz| [cm]
static constexpr double kAccMaxEta = 1.2; // acceptance |eta| (primary particles)
static constexpr double kAccMaxVz = 10.0; // acceptance |vz| (primary particles) [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 @@ -112,28 +110,13 @@ struct NucleiAntineutronCex {
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.}});
histos.add("cexPairMcPz", "CEX pair p_{z};p_{z} (GeV/c);Entries", kTH1F, {{100, -10., 10.}});
histos.add("cex_pairmcDplane", "CEX pair d_{plane};d_{plane} (cm);Entries", kTH1F, {{100, 0., 10.}});
histos.add("cex_pairmc_angle", "Pair opening angle;Angle (°);Entries", kTH1F, {{180, 0., 180.}});
histos.add("cex_pairmc_vtx", "MC CEX pair vertex;X (cm);Y (cm)", kTH2F, {{100, -50., 50.}, {100, -50., 50.}});
histos.add("cex_pairmc_vtx", "MC CEX pair vertex;X (cm);Y (cm)", kTH2F, {{200, -60., 60.}, {200, -60., 60.}});
histos.add("cex_pairmc_vtxz", "MC secondary vertex Z;Z (cm);Entries", kTH1F, {{200, -60., 60.}});
histos.add("cexPairMcPITScuts", "CEX pair momentum (ITS cuts);|p| (GeV/c);Entries", kTH1F, {{100, 0., 10.}});

Expand All @@ -148,16 +131,12 @@ struct NucleiAntineutronCex {
histos.add("cexbg_pairmc_pz", "Background pair p_{z};p_{z} (GeV/c);Entries", kTH1F, {{100, -10., 10.}});
histos.add("cexbg_pairmcDplane", "Background d_{plane};d_{plane} (cm);Entries", kTH1F, {{100, 0., 10.}});
histos.add("cexbg_pairmc_angle", "Background opening angle;Angle (°);Entries", kTH1F, {{180, 0., 180.}});
histos.add("cexbg_pairmc_vtx", "Background pair vertex;X (cm);Y (cm)", kTH2F, {{100, -50., 50.}, {100, -50., 50.}});
histos.add("cexbg_pairmc_vtx", "Background pair vertex;X (cm);Y (cm)", kTH2F, {{200, -60., 60.}, {200, -60., 60.}});
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 IB hit layers (L0-L2) (CEX);N_{IB layers};Entries", kTH1I, {{5, -1.5, 4.5}});
histos.add("hMotherNHitIB_BG", "Mother IB hit layers (L0-L2) (BG);N_{IB layers};Entries", kTH1I, {{5, -1.5, 4.5}});
// Pi0 events
histos.add("cexPairMcP_pi0", "CEX pair total momentum;|p| (GeV/c);Entries", kTH1F, {{100, 0., 10.}});

// CEX pair from antineutron (TRK)
histos.add("cex_pairtrk_angle", "Pair opening angle (tracks);Angle (°);Entries", kTH1F, {{180, 0., 180.}});
Expand Down Expand Up @@ -329,7 +308,7 @@ struct NucleiAntineutronCex {
double antipE = particle.e();
int antipId = particle.globalIndex();

// Selection conditions: Produced in the ITS
// Selection conditions: Produced in the ITS IB
const double r = std::sqrt(antipVx * antipVx + antipVy * antipVy);
// Config for ITS
// if(3.9<=r && r<=43.0 && std::abs(antipVz)<=48.9){
Expand Down Expand Up @@ -396,6 +375,33 @@ struct NucleiAntineutronCex {
if (pionPlus || pionMinus)
continue;

// Pion 0 count
bool pion0 = false;
for (const auto& particle4 : mcPartsThis) {
if (particle4.mcCollisionId() != colId)
continue;
const auto proc4Enum = particle4.getProcess();
const bool isSecondaryFromMaterial4 = (!particle4.producedByGenerator()) && (proc4Enum == kPHadronic || proc4Enum == kPHInhelastic);
if (particle4.pdgCode() != kPi0 || !isSecondaryFromMaterial4 || particle4.mothersIds().empty())
continue;
bool hasPrimaryMotherPi0 = false;
for (const auto& mother : particle4.mothers_as<aod::McParticles>()) {
if (mother.isPhysicalPrimary()) {
hasPrimaryMotherPi0 = true;
break;
}
}
if (!hasPrimaryMotherPi0)
continue;
double pi0Vx = particle4.vx();
double pi0Vy = particle4.vy();
double pi0Vz = particle4.vz();
if (std::abs(pi0Vx - antipVx) < kVtxTol && std::abs(pi0Vy - antipVy) < kVtxTol && std::abs(pi0Vz - antipVz) < kVtxTol) {
pion0 = true;
break;
}
}

// CEX selection
double dplane = kInitDplane;
double dplaneTmp = 0;
Expand Down Expand Up @@ -453,7 +459,6 @@ struct NucleiAntineutronCex {
continue;

// CEX proton selection
// dplaneTmp = (pPy*antipPz - pPz*antipPy)*(pvtxX-antipVx) + (pPz*antipPx - pPx*antipPz)*(pvtxY-antipVy) + (pPx*antipPy - pPy*antipPx)*(pvtxZ-antipVz);
double nx = (pPy * antipPz - pPz * antipPy);
double ny = (pPz * antipPx - pPx * antipPz);
double nz = (pPx * antipPy - pPy * antipPx);
Expand Down Expand Up @@ -518,6 +523,8 @@ struct NucleiAntineutronCex {
histos.fill(HIST("cexn_pairmc_pt"), cexPairMcPt / motherPt);
if (motherPz != 0)
histos.fill(HIST("cexn_pairmc_pz"), cexPairMcPz / motherPz);
if (motherP != 0 && pion0)
histos.fill(HIST("cexPairMcP_pi0"), cexPairMcP / motherP);
}
// BG mother
if (motherPdg != -kNeutron) {
Expand All @@ -533,8 +540,8 @@ struct NucleiAntineutronCex {
histos.fill(HIST("cexbg_pairmc_pITScuts"), cexPairMcP);
}

// Detector signal
bool antipLayers = false;
// Reconstructed data
bool antipLayers_condition = false;
bool antipHasTrack = false;
double antipTrkPx = 0.;
double antipTrkPy = 0.;
Expand All @@ -549,7 +556,7 @@ struct NucleiAntineutronCex {
int8_t pTrkItsPidValid = 0;
float pTrkTgl = 0.f;

bool pLayers = false;
bool pLayers_condition = false;
bool pHasTrack = false;
double pTrkPx = 0.;
double pTrkPy = 0.;
Expand Down Expand Up @@ -614,9 +621,8 @@ struct NucleiAntineutronCex {
antipTrkItsPidValid = std::isfinite(nsigmaITSantip) ? 1 : 0;
antipHasTrack = true;
apItsMap = static_cast<uint16_t>(track.itsClusterMap());
antipLayers = (apItsMap != 0);
if (layerCondition)
antipLayers = true;
antipLayers_condition = true;
if (motherPdg == -kNeutron) {
histos.fill(HIST("apItsNsigmaPr"), antipTrkItsNSigmaPr);
histos.fill(HIST("apItsPidValid"), antipTrkItsPidValid);
Expand All @@ -643,9 +649,8 @@ struct NucleiAntineutronCex {
pTrkItsPidValid = std::isfinite(nsigmaITSp) ? 1 : 0;
pHasTrack = true;
pItsMap = static_cast<uint16_t>(track.itsClusterMap());
pLayers = (pItsMap != 0);
if (layerCondition)
pLayers = true;
pLayers_condition = true;
if (motherPdg == -kNeutron) {
histos.fill(HIST("pItsNsigmaPr"), pTrkItsNSigmaPr);
histos.fill(HIST("pItsPidValid"), pTrkItsPidValid);
Expand Down Expand Up @@ -723,7 +728,7 @@ struct NucleiAntineutronCex {
if (motherPdg != -kNeutron)
histos.fill(HIST("cexbg_pairtrkVtxfitDcaPair"), dcaPair);

if (!(antipLayers && pLayers))
if (!(antipLayers_condition && pLayers_condition))
continue;
double cexPairTrkP = total_trk_pVec.Mag();
double cexPairTrkPt = total_trk_pVec.Pt();
Expand Down Expand Up @@ -813,105 +818,17 @@ struct NucleiAntineutronCex {
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 {
++vtxNNeut;
}
}

// 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(
isCex,
motherPdg,
motherNHitIB,
motherP,

colId,
pId,
antipId,
Expand Down Expand Up @@ -986,8 +903,6 @@ struct NucleiAntineutronCex {

pItsMap,
apItsMap,
static_cast<int8_t>(pLayers ? 1 : 0),
static_cast<int8_t>(antipLayers ? 1 : 0),

pvtxZ,

Expand Down
Loading