diff --git a/PWGLF/DataModel/LFKinkDecayTables.h b/PWGLF/DataModel/LFKinkDecayTables.h index 5f8798856d5..9533bd66116 100644 --- a/PWGLF/DataModel/LFKinkDecayTables.h +++ b/PWGLF/DataModel/LFKinkDecayTables.h @@ -15,10 +15,11 @@ /// \author Francesco Mazzaschi /// -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" #include "Common/Core/RecoDecay.h" +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisDataModel.h" + #ifndef PWGLF_DATAMODEL_LFKINKDECAYTABLES_H_ #define PWGLF_DATAMODEL_LFKINKDECAYTABLES_H_ @@ -77,6 +78,13 @@ DECLARE_SOA_DYNAMIC_COLUMN(MSigmaPlus, mSigmaPlus, //! mass under sigma plus hyp float pzneut = pzmoth - pzch; return RecoDecay::m(std::array{std::array{pxch, pych, pzch}, std::array{pxneut, pyneut, pzneut}}, std::array{o2::constants::physics::MassProton, o2::constants::physics::MassPionNeutral}); }); +DECLARE_SOA_DYNAMIC_COLUMN(MXiMinus, mXiMinus, //! mass under Xi minus hypothesis + [](float pxmoth, float pymoth, float pzmoth, float pxch, float pych, float pzch) -> float { + float pxneut = pxmoth - pxch; + float pyneut = pymoth - pych; + float pzneut = pzmoth - pzch; + return RecoDecay::m(std::array{std::array{pxch, pych, pzch}, std::array{pxneut, pyneut, pzneut}}, std::array{o2::constants::physics::MassPionCharged, o2::constants::physics::MassLambda}); }); + } // namespace kinkcand DECLARE_SOA_TABLE(KinkCands, "AOD", "KINKCANDS", @@ -93,7 +101,24 @@ DECLARE_SOA_TABLE(KinkCands, "AOD", "KINKCANDS", kinkcand::PtMoth, kinkcand::PtDaug, kinkcand::MSigmaMinus, - kinkcand::MSigmaPlus); + kinkcand::MSigmaPlus, + kinkcand::MXiMinus); + +DECLARE_SOA_TABLE(KinkCandsUnbound, "AOD", "UBKINKCANDS", + o2::soa::Index<>, kinkcand::XDecVtx, kinkcand::YDecVtx, kinkcand::ZDecVtx, + kinkcand::MothSign, kinkcand::PxMoth, kinkcand::PyMoth, kinkcand::PzMoth, + kinkcand::PxDaug, kinkcand::PyDaug, kinkcand::PzDaug, + kinkcand::DcaMothPv, kinkcand::DcaDaugPv, kinkcand::DcaKinkTopo, + + // dynamic columns + kinkcand::PxDaugNeut, + kinkcand::PyDaugNeut, + kinkcand::PzDaugNeut, + kinkcand::PtMoth, + kinkcand::PtDaug, + kinkcand::MSigmaMinus, + kinkcand::MSigmaPlus, + kinkcand::MXiMinus); } // namespace o2::aod diff --git a/PWGLF/TableProducer/Common/kinkBuilder.cxx b/PWGLF/TableProducer/Common/kinkBuilder.cxx index 51be861741a..64605ec9897 100644 --- a/PWGLF/TableProducer/Common/kinkBuilder.cxx +++ b/PWGLF/TableProducer/Common/kinkBuilder.cxx @@ -13,30 +13,31 @@ /// \brief Builder task for kink decay topologies using ITS standalone tracks for the mother /// \author Francesco Mazzaschi -#include -#include -#include -#include -#include +#include "PWGLF/DataModel/LFKinkDecayTables.h" +#include "PWGLF/DataModel/LFParticleIdentification.h" +#include "PWGLF/Utils/svPoolCreator.h" -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" -#include "ReconstructionDataFormats/Track.h" +#include "Common/DataModel/Centrality.h" #include "Common/DataModel/EventSelection.h" #include "Common/DataModel/Multiplicity.h" -#include "Common/DataModel/Centrality.h" -#include "DetectorsBase/Propagator.h" -#include "DetectorsBase/GeometryManager.h" -#include "DataFormatsParameters/GRPObject.h" -#include "DataFormatsParameters/GRPMagField.h" + #include "CCDB/BasicCCDBManager.h" #include "DCAFitter/DCAFitterN.h" +#include "DataFormatsParameters/GRPMagField.h" +#include "DataFormatsParameters/GRPObject.h" +#include "DetectorsBase/GeometryManager.h" +#include "DetectorsBase/Propagator.h" +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisDataModel.h" +#include "Framework/AnalysisTask.h" +#include "Framework/runDataProcessing.h" +#include "ReconstructionDataFormats/Track.h" -#include "PWGLF/DataModel/LFParticleIdentification.h" -#include "PWGLF/Utils/svPoolCreator.h" -#include "PWGLF/DataModel/LFKinkDecayTables.h" +#include +#include +#include +#include +#include using namespace o2; using namespace o2::framework; @@ -93,9 +94,12 @@ struct kinkBuilder { kHyperhelium4sigma }; Produces outputDataTable; + Produces outputDataTableUB; + Service ccdb; Configurable hypoMoth{"hypoMoth", kSigmaMinus, "Mother particle hypothesis"}; + Configurable fillDebugTable{"fillDebugTable", false, "If true, fill the debug table with all candidates unbound"}; // Selection criteria Configurable maxDCAMothToPV{"maxDCAMothToPV", 0.1, "Max DCA of the mother to the PV"}; Configurable minDCADaugToPV{"minDCADaugToPV", 0., "Min DCA of the daughter to the PV"}; @@ -430,11 +434,18 @@ struct kinkBuilder { std::sort(kinkCandidates.begin(), kinkCandidates.end(), [](const kinkCandidate& a, const kinkCandidate& b) { return a.collisionID < b.collisionID; }); for (const auto& kinkCand : kinkCandidates) { - outputDataTable(kinkCand.collisionID, kinkCand.mothTrackID, kinkCand.daugTrackID, - kinkCand.decVtx[0], kinkCand.decVtx[1], kinkCand.decVtx[2], - kinkCand.mothSign, kinkCand.momMoth[0], kinkCand.momMoth[1], kinkCand.momMoth[2], - kinkCand.momDaug[0], kinkCand.momDaug[1], kinkCand.momDaug[2], - kinkCand.dcaXYmoth, kinkCand.dcaXYdaug, kinkCand.dcaKinkTopo); + if (fillDebugTable) { + outputDataTableUB(kinkCand.decVtx[0], kinkCand.decVtx[1], kinkCand.decVtx[2], + kinkCand.mothSign, kinkCand.momMoth[0], kinkCand.momMoth[1], kinkCand.momMoth[2], + kinkCand.momDaug[0], kinkCand.momDaug[1], kinkCand.momDaug[2], + kinkCand.dcaXYmoth, kinkCand.dcaXYdaug, kinkCand.dcaKinkTopo); + } else { + outputDataTable(kinkCand.collisionID, kinkCand.mothTrackID, kinkCand.daugTrackID, + kinkCand.decVtx[0], kinkCand.decVtx[1], kinkCand.decVtx[2], + kinkCand.mothSign, kinkCand.momMoth[0], kinkCand.momMoth[1], kinkCand.momMoth[2], + kinkCand.momDaug[0], kinkCand.momDaug[1], kinkCand.momDaug[2], + kinkCand.dcaXYmoth, kinkCand.dcaXYdaug, kinkCand.dcaKinkTopo); + } } } }; diff --git a/PWGLF/Tasks/Strangeness/sigmaminustask.cxx b/PWGLF/Tasks/Strangeness/sigmaminustask.cxx index 8b309476913..f7332656cb2 100644 --- a/PWGLF/Tasks/Strangeness/sigmaminustask.cxx +++ b/PWGLF/Tasks/Strangeness/sigmaminustask.cxx @@ -13,11 +13,13 @@ /// \brief Example of a simple task for the analysis of the Sigma-minus /// \author Francesco Mazzaschi -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Common/DataModel/EventSelection.h" #include "PWGLF/DataModel/LFKinkDecayTables.h" + +#include "Common/DataModel/EventSelection.h" #include "Common/DataModel/PIDResponse.h" + +#include "Framework/AnalysisTask.h" +#include "Framework/runDataProcessing.h" #include "ReconstructionDataFormats/PID.h" using namespace o2; @@ -26,6 +28,7 @@ using namespace o2::framework::expressions; using TracksFull = soa::Join; using CollisionsFull = soa::Join; +using CollisionsFullMC = soa::Join; struct sigmaminustask { // Histograms are defined with HistogramRegistry @@ -42,18 +45,24 @@ struct sigmaminustask { const AxisSpec ptAxis{50, -10, 10, "#it{p}_{T} (GeV/#it{c})"}; const AxisSpec nSigmaPiAxis{100, -5, 5, "n#sigma_{#pi}"}; const AxisSpec sigmaMassAxis{100, 1.1, 1.4, "m (GeV/#it{c}^{2})"}; + const AxisSpec xiMassAxis{100, 1.2, 1.6, "m_{#Xi} (GeV/#it{c}^{2})"}; + const AxisSpec pdgAxis{10001, -5000, 5000, "PDG code"}; const AxisSpec vertexZAxis{100, -15., 15., "vrtx_{Z} [cm]"}; // Event selection rEventSelection.add("hVertexZRec", "hVertexZRec", {HistType::kTH1F, {vertexZAxis}}); - // Sigma-minus reconstruction rSigmaMinus.add("h2MassSigmaMinusPt", "h2MassSigmaMinusPt", {HistType::kTH2F, {ptAxis, sigmaMassAxis}}); + rSigmaMinus.add("h2SigmaMassVsXiMass", "h2SigmaMassVsXiMass", {HistType::kTH2F, {xiMassAxis, sigmaMassAxis}}); rSigmaMinus.add("h2NSigmaPiPt", "h2NSigmaPiPt", {HistType::kTH2F, {ptAxis, nSigmaPiAxis}}); + + if (doprocessMC) { + // Add MC histograms if needed + rSigmaMinus.add("h2MassSigmaMinusPtMC", "h2MassSigmaMinusPtMC", {HistType::kTH2F, {ptAxis, sigmaMassAxis}}); + } } - void process(soa::Join::iterator const& collision, - aod::KinkCands const& KinkCands, TracksFull const&) + void processData(CollisionsFull::iterator const& collision, aod::KinkCands const& KinkCands, TracksFull const&) { if (std::abs(collision.posZ()) > cutzvertex || !collision.sel8()) { return; @@ -65,9 +74,54 @@ struct sigmaminustask { continue; } rSigmaMinus.fill(HIST("h2MassSigmaMinusPt"), kinkCand.mothSign() * kinkCand.ptMoth(), kinkCand.mSigmaMinus()); + rSigmaMinus.fill(HIST("h2SigmaMassVsXiMass"), kinkCand.mXiMinus(), kinkCand.mSigmaMinus()); rSigmaMinus.fill(HIST("h2NSigmaPiPt"), kinkCand.mothSign() * kinkCand.ptMoth(), dauTrack.tpcNSigmaPi()); } } + PROCESS_SWITCH(sigmaminustask, processData, "Data processing", true); + + void processMC(CollisionsFullMC::iterator const& collision, aod::KinkCands const& KinkCands, aod::McTrackLabels const& trackLabelsMC, aod::McParticles const&, TracksFull const&) + { + if (std::abs(collision.posZ()) > cutzvertex || !collision.sel8()) { + return; + } + rEventSelection.fill(HIST("hVertexZRec"), collision.posZ()); + for (const auto& kinkCand : KinkCands) { + auto dauTrack = kinkCand.trackDaug_as(); + auto mothTrack = kinkCand.trackMoth_as(); + if (dauTrack.sign() != mothTrack.sign()) { + LOG(info) << "Skipping kink candidate with opposite sign daughter and mother: " << kinkCand.globalIndex(); + continue; // Skip if the daughter has the opposite sign as the mother + } + if (abs(dauTrack.tpcNSigmaPi()) > cutNSigmaPi) { + continue; + } + + rSigmaMinus.fill(HIST("h2MassSigmaMinusPt"), kinkCand.mothSign() * kinkCand.ptMoth(), kinkCand.mSigmaMinus()); + rSigmaMinus.fill(HIST("h2SigmaMassVsXiMass"), kinkCand.mXiMinus(), kinkCand.mSigmaMinus()); + rSigmaMinus.fill(HIST("h2NSigmaPiPt"), kinkCand.mothSign() * kinkCand.ptMoth(), dauTrack.tpcNSigmaPi()); + // do MC association + auto mcLabSigma = trackLabelsMC.rawIteratorAt(mothTrack.globalIndex()); + auto mcLabPiDau = trackLabelsMC.rawIteratorAt(dauTrack.globalIndex()); + if (mcLabSigma.has_mcParticle() && mcLabPiDau.has_mcParticle()) { + auto mcTrackSigma = mcLabSigma.mcParticle_as(); + auto mcTrackPiDau = mcLabPiDau.mcParticle_as(); + if (!mcTrackPiDau.has_mothers()) { + continue; + } + for (auto& piMother : mcTrackPiDau.mothers_as()) { + if (piMother.globalIndex() != mcTrackSigma.globalIndex()) { + continue; + } + if (std::abs(mcTrackSigma.pdgCode()) != 3112 || std::abs(mcTrackPiDau.pdgCode()) != 211) { + continue; + } + rSigmaMinus.fill(HIST("h2MassSigmaMinusPtMC"), kinkCand.mothSign() * kinkCand.ptMoth(), kinkCand.mSigmaMinus()); + } + } + } + } + PROCESS_SWITCH(sigmaminustask, processMC, "MC processing", false); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)