Skip to content
Merged
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
40 changes: 20 additions & 20 deletions PWGLF/DataModel/LFPhotonDeuteronTables.h
Original file line number Diff line number Diff line change
Expand Up @@ -23,30 +23,30 @@ namespace o2::aod
{
namespace photondeuteron
{
DECLARE_SOA_COLUMN(PhotonPt, photonPt, float); //! Photon transverse momentum
DECLARE_SOA_COLUMN(PhotonEta, photonEta, float); //! Photon pseudorapidity
DECLARE_SOA_COLUMN(PhotonPhi, photonPhi, float); //! Photon azimuthal angle
DECLARE_SOA_COLUMN(PhotonMass, photonMass, float); //! Photon invariant mass
DECLARE_SOA_COLUMN(PhotonPosPt, photonPosPt, float); //! Positive daughter (e+) transverse momentum
DECLARE_SOA_COLUMN(PhotonPosEta, photonPosEta, float); //! Positive daughter (e+) pseudorapidity
DECLARE_SOA_COLUMN(PhotonPosPhi, photonPosPhi, float); //! Positive daughter (e+) azimuthal angle
DECLARE_SOA_COLUMN(PhotonPt, photonPt, float); //! Photon transverse momentum
DECLARE_SOA_COLUMN(PhotonEta, photonEta, float); //! Photon pseudorapidity
DECLARE_SOA_COLUMN(PhotonPhi, photonPhi, float); //! Photon azimuthal angle
DECLARE_SOA_COLUMN(PhotonMass, photonMass, float); //! Photon invariant mass
DECLARE_SOA_COLUMN(PhotonPosPt, photonPosPt, float); //! Positive daughter (e+) transverse momentum
DECLARE_SOA_COLUMN(PhotonPosEta, photonPosEta, float); //! Positive daughter (e+) pseudorapidity
DECLARE_SOA_COLUMN(PhotonPosPhi, photonPosPhi, float); //! Positive daughter (e+) azimuthal angle
DECLARE_SOA_COLUMN(PhotonPosNSigmaElTPC, photonPosNSigmaElTPC, float); //! Positive daughter (e+) NSigma electron TPC
DECLARE_SOA_COLUMN(PhotonNegPt, photonNegPt, float); //! Negative daughter (e-) transverse momentum
DECLARE_SOA_COLUMN(PhotonNegEta, photonNegEta, float); //! Negative daughter (e-) pseudorapidity
DECLARE_SOA_COLUMN(PhotonNegPhi, photonNegPhi, float); //! Negative daughter (e-) azimuthal angle
DECLARE_SOA_COLUMN(PhotonNegPt, photonNegPt, float); //! Negative daughter (e-) transverse momentum
DECLARE_SOA_COLUMN(PhotonNegEta, photonNegEta, float); //! Negative daughter (e-) pseudorapidity
DECLARE_SOA_COLUMN(PhotonNegPhi, photonNegPhi, float); //! Negative daughter (e-) azimuthal angle
DECLARE_SOA_COLUMN(PhotonNegNSigmaElTPC, photonNegNSigmaElTPC, float); //! Negative daughter (e-) NSigma electron TPC
DECLARE_SOA_COLUMN(DeuteronPt, deuteronPt, float); //! Deuteron transverse momentum
DECLARE_SOA_COLUMN(DeuteronEta, deuteronEta, float); //! Deuteron pseudorapidity
DECLARE_SOA_COLUMN(DeuteronPhi, deuteronPhi, float); //! Deuteron azimuthal angle
DECLARE_SOA_COLUMN(DeuteronNSigmaTPC, deuteronNSigmaTPC, float); //! Deuteron NSigma TPC
DECLARE_SOA_COLUMN(DeuteronNSigmaTOF, deuteronNSigmaTOF, float); //! Deuteron NSigma TOF
DECLARE_SOA_COLUMN(DeltaPhi, deltaPhi, float); //! Delta phi between photon and deuteron
DECLARE_SOA_COLUMN(DeltaEta, deltaEta, float); //! Delta eta between photon and deuteron
DECLARE_SOA_COLUMN(InvMass, invMass, float); //! Photon-deuteron invariant mass
DECLARE_SOA_COLUMN(RelativeMomentum, relativeMomentum, float); //! Relative momentum k*_pn
DECLARE_SOA_COLUMN(DeuteronPt, deuteronPt, float); //! Deuteron transverse momentum
DECLARE_SOA_COLUMN(DeuteronEta, deuteronEta, float); //! Deuteron pseudorapidity
DECLARE_SOA_COLUMN(DeuteronPhi, deuteronPhi, float); //! Deuteron azimuthal angle
DECLARE_SOA_COLUMN(DeuteronNSigmaTPC, deuteronNSigmaTPC, float); //! Deuteron NSigma TPC
DECLARE_SOA_COLUMN(DeuteronNSigmaTOF, deuteronNSigmaTOF, float); //! Deuteron NSigma TOF
DECLARE_SOA_COLUMN(DeltaPhi, deltaPhi, float); //! Delta phi between photon and deuteron
DECLARE_SOA_COLUMN(DeltaEta, deltaEta, float); //! Delta eta between photon and deuteron
DECLARE_SOA_COLUMN(InvMass, invMass, float); //! Photon-deuteron invariant mass
DECLARE_SOA_COLUMN(RelativeMomentum, relativeMomentum, float); //! Relative momentum k*_pn
} // namespace photondeuteron

DECLARE_SOA_TABLE(PhotonDeuteronPairs, "AOD", "PHOTONDPAIRS", //! Table for photon-deuteron pairs
DECLARE_SOA_TABLE(PhotonDeuteronPairs, "AOD", "PHOTONDPAIRS", //! Table for photon-deuteron pairs
photondeuteron::PhotonPt,
photondeuteron::PhotonEta,
photondeuteron::PhotonPhi,
Expand Down
90 changes: 46 additions & 44 deletions PWGLF/TableProducer/Nuspex/photonDeuteron.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -13,18 +13,20 @@
/// \brief Task to analyze photon-deuteron correlations
/// \author Arvind Khuntia <arvind.khuntia@cern.ch> and Francesco Noferini <francesco.noferini@cern.ch>

#include <TLorentzVector.h>
#include <TVector3.h>
#include "PWGLF/DataModel/LFPhotonDeuteronTables.h"
#include "PWGLF/DataModel/LFStrangenessTables.h"

#include "Common/DataModel/EventSelection.h"
#include "Common/DataModel/PIDResponse.h"
#include "Common/DataModel/TrackSelectionTables.h"
#include "Framework/AnalysisTask.h"

#include "Framework/ASoAHelpers.h"
#include "Framework/AnalysisTask.h"
#include "Framework/HistogramRegistry.h"
#include "Framework/runDataProcessing.h"
#include "PWGLF/DataModel/LFStrangenessTables.h"
#include "PWGLF/DataModel/LFPhotonDeuteronTables.h"

#include <TLorentzVector.h>
#include <TVector3.h>

using namespace o2;
using namespace o2::framework;
Expand All @@ -43,13 +45,13 @@ struct PhotonDeuteronCorrelation {
Configurable<float> cfgPhotonMaxPt{"cfgPhotonMaxPt", 10.0, "Maximum photon pT (GeV/c)"};
Configurable<float> cfgPhotonMaxEta{"cfgPhotonMaxEta", 0.8, "Maximum photon eta"};
Configurable<float> cfgPhotonDaughterNSigmaElTPCMax{"cfgPhotonDaughterNSigmaElTPCMax", 5.0, "Maximum NSigma TPC electron for photon daughters (999=off)"};

Configurable<float> cfgDeuteronMinPt{"cfgDeuteronMinPt", 0.2, "Minimum deuteron pT (GeV/c)"};
Configurable<float> cfgDeuteronMaxPt{"cfgDeuteronMaxPt", 5.0, "Maximum deuteron pT (GeV/c)"};
Configurable<float> cfgDeuteronMaxEta{"cfgDeuteronMaxEta", 0.8, "Maximum deuteron eta"};
Configurable<float> cfgDeuteronNSigmaTPCMax{"cfgDeuteronNSigmaTPCMax", 3.0, "Maximum NSigma TPC for deuteron"};
Configurable<float> cfgDeuteronNSigmaTOFMax{"cfgDeuteronNSigmaTOFMax", 3.0, "Maximum NSigma TOF for deuteron"};

Configurable<float> cfgV0CosPA{"cfgV0CosPA", 0.995, "Minimum V0 cosine of pointing angle"};
Configurable<float> cfgV0Radius{"cfgV0Radius", 5.0, "Minimum V0 radius (cm)"};
Configurable<bool> cfgUsePhotonDaughterPIDTPCOnly{"cfgUsePhotonDaughterPIDTPCOnly", true, "Use TPC-only PID for photon daughters"};
Expand All @@ -75,7 +77,7 @@ struct PhotonDeuteronCorrelation {

// Event histograms
histos.add("hEventCounter", "Event counter", kTH1F, {{5, 0.5, 5.5}});

// V0 (photon) histograms
histos.add("hV0Mass", "V0 invariant mass", kTH1F, {axisPhotonMass});
histos.add("hV0Pt", "V0 p_{T}", kTH1F, {axisPt});
Expand All @@ -84,7 +86,7 @@ struct PhotonDeuteronCorrelation {
histos.add("hV0CosPA", "V0 cos(PA)", kTH1F, {axisV0CosPA});
histos.add("hV0Radius", "V0 radius", kTH1F, {axisV0Radius});
histos.add("hPhotonPtEta", "Photon p_{T} vs #eta", kTH2F, {axisPt, axisEta});

// Deuteron histograms
histos.add("hDeuteronPt", "Deuteron p_{T}", kTH1F, {axisPt});
histos.add("hDeuteronEta", "Deuteron #eta", kTH1F, {axisEta});
Expand All @@ -93,7 +95,7 @@ struct PhotonDeuteronCorrelation {
histos.add("hDeuteronNSigmaTOF", "Deuteron n#sigma TOF", kTH1F, {axisNSigma});
histos.add("hDeuteronNSigmaTPCvsPt", "Deuteron n#sigma TPC vs p_{T}", kTH2F, {axisPt, axisNSigma});
histos.add("hDeuteronNSigmaTOFvsPt", "Deuteron n#sigma TOF vs p_{T}", kTH2F, {axisPt, axisNSigma});

// Correlation histograms
histos.add("hPhotonDeuteronDeltaPhi", "Photon-Deuteron #Delta#phi", kTH1F, {axisDeltaPhi});
histos.add("hPhotonDeuteronDeltaEta", "Photon-Deuteron #Delta#eta", kTH1F, {axisDeltaEta});
Expand All @@ -115,33 +117,33 @@ struct PhotonDeuteronCorrelation {
if (v0.v0radius() < cfgV0Radius) {
return false;
}

// Photon mass window
if (std::abs(v0.mGamma()) > cfgPhotonMassWindow) {
return false;
}

// Kinematic cuts
if (v0.pt() < cfgPhotonMinPt || v0.pt() > cfgPhotonMaxPt) {
return false;
}
if (std::abs(v0.eta()) > cfgPhotonMaxEta) {
return false;
}

// Optional electron PID cuts for daughter tracks (TPC-only)
if (cfgUsePhotonDaughterPIDTPCOnly) {
auto posTrack = v0.template posTrack_as<TTracks>();
auto negTrack = v0.template negTrack_as<TTracks>();

if (std::abs(posTrack.tpcNSigmaEl()) > cfgPhotonDaughterNSigmaElTPCMax) {
return false;
}
if (std::abs(negTrack.tpcNSigmaEl()) > cfgPhotonDaughterNSigmaElTPCMax) {
return false;
}
}

return true;
}

Expand All @@ -156,21 +158,21 @@ struct PhotonDeuteronCorrelation {
if (std::abs(track.eta()) > cfgDeuteronMaxEta) {
return false;
}

// TPC PID
if (std::abs(track.tpcNSigmaDe()) > cfgDeuteronNSigmaTPCMax) {
return false;
}

// TOF PID (if available)
if (track.hasTOF() && std::abs(track.tofNSigmaDe()) > cfgDeuteronNSigmaTOFMax) {
return false;
}

return true;
}

//range [-pi/2, 3pi/2]
// range [-pi/2, 3pi/2]
float getDeltaPhi(float phi1, float phi2)
{
float dphi = phi1 - phi2;
Expand All @@ -189,21 +191,21 @@ struct PhotonDeuteronCorrelation {
if (invMass <= 0.0) {
return -1.0;
}

float M = invMass;
float M2 = M * M;
float M4 = M2 * M2;
float mn2 = massNeutron * massNeutron;
float mp2 = massProton * massProton;
float deltaMass2 = mn2 - mp2;
float sumMass2 = mn2 + mp2;

float term = M4 + deltaMass2 * deltaMass2 - 2.0 * sumMass2 * M2;

if (term < 0.0) {
return -1.0;
}

float kpn = 0.5 / M * std::sqrt(term);
return kpn;
}
Expand All @@ -217,28 +219,29 @@ struct PhotonDeuteronCorrelation {
{
histos.fill(HIST("hEventCounter"), 1);


// Loop over V0s to find photons
std::vector<int> photonIndices;
for (auto const& v0 : V0s) {
// Fill V0 QA histograms
histos.fill(HIST("hV0Mass"), v0.mGamma());
histos.fill(HIST("hV0CosPA"), v0.v0cosPA());
histos.fill(HIST("hV0Radius"), v0.v0radius());

// Select photons
if (!selectPhoton(v0, tracks)) {
continue;
}

// Fill photon histograms
histos.fill(HIST("hV0Pt"), v0.pt());
histos.fill(HIST("hV0Eta"), v0.eta());
histos.fill(HIST("hV0Phi"), v0.phi());
histos.fill(HIST("hPhotonPtEta"), v0.pt(), v0.eta());

if (v0.isPhotonTPConly()) photonIndices.push_back(v0.index());
if(v0.isPhotonTPConly()) std::cout<<" [main] global index photon: "<<v0.globalIndex()<<" v0 id: "<<v0.index()<<" pt "<<v0.pt()<<std::endl;

if (v0.isPhotonTPConly())
photonIndices.push_back(v0.index());
if (v0.isPhotonTPConly())
std::cout << " [main] global index photon: " << v0.globalIndex() << " v0 id: " << v0.index() << " pt " << v0.pt() << std::endl;
}

// Loop over tracks to find deuterons
Expand All @@ -248,19 +251,19 @@ struct PhotonDeuteronCorrelation {
if (!track.isGlobalTrack()) {
continue;
}

// Select deuterons
if (!selectDeuteron(track)) {
continue;
}

// Fill deuteron histograms
histos.fill(HIST("hDeuteronPt"), track.pt());
histos.fill(HIST("hDeuteronEta"), track.eta());
histos.fill(HIST("hDeuteronPhi"), track.phi());
histos.fill(HIST("hDeuteronNSigmaTPC"), track.tpcNSigmaDe());
histos.fill(HIST("hDeuteronNSigmaTPCvsPt"), track.pt(), track.tpcNSigmaDe());

if (track.hasTOF()) {
histos.fill(HIST("hDeuteronNSigmaTOF"), track.tofNSigmaDe());
histos.fill(HIST("hDeuteronNSigmaTOFvsPt"), track.pt(), track.tofNSigmaDe());
Expand All @@ -269,11 +272,11 @@ struct PhotonDeuteronCorrelation {
}
// Calculate correlations between photons and deuterons
for (auto const& photonIdx : photonIndices) {
const auto& photon = V0s.iteratorAt(photonIdx);
const auto& photon = V0s.iteratorAt(photonIdx);

for (auto const& deuteronIdx : deuteronIndices) {
const auto& deuteron = tracks.iteratorAt(deuteronIdx);
const auto& deuteron = tracks.iteratorAt(deuteronIdx);

// Calculate angular correlations
float deltaPhi = getDeltaPhi(photon.phi(), deuteron.phi());
float deltaEta = photon.eta() - deuteron.eta();
Expand All @@ -283,27 +286,27 @@ struct PhotonDeuteronCorrelation {
histos.fill(HIST("hPhotonDeuteronDeltaEta"), deltaEta);
histos.fill(HIST("hPhotonDeuteronCorrelation"), deltaPhi, deltaEta);
histos.fill(HIST("hPhotonDeuteronPtCorr"), photon.pt(), deuteron.pt());

// Calculate invariant mass
TLorentzVector photonVec, deuteronVec;
photonVec.SetPtEtaPhiM(photon.pt(), photon.eta(), photon.phi(), 0.0); // Photon-mass = 0
photonVec.SetPtEtaPhiM(photon.pt(), photon.eta(), photon.phi(), 0.0); // Photon-mass = 0
deuteronVec.SetPtEtaPhiM(deuteron.pt(), deuteron.eta(), deuteron.phi(), massDeuteron); // Deuteron-mass

TLorentzVector combinedVec = photonVec + deuteronVec;
float invMass = combinedVec.M();
histos.fill(HIST("hPhotonDeuteronInvMass"), invMass);

// Calculate relative momentum using Equation 4.6
float kpn = calculateRelativeMomentum(invMass);
if (kpn >= 0.0) {
histos.fill(HIST("hRelativeMomentum"), kpn);
histos.fill(HIST("hRelativeMomentumVsInvMass"), invMass, kpn);
}

// Fill the output table
auto posTrack = photon.posTrack_as<soa::Join<aod::Tracks, aod::TracksExtra, aod::TrackSelection, aod::pidTPCFullDe, aod::pidTOFFullDe, aod::pidTPCFullEl, aod::pidTOFFullEl>>();
auto negTrack = photon.negTrack_as<soa::Join<aod::Tracks, aod::TracksExtra, aod::TrackSelection, aod::pidTPCFullDe, aod::pidTOFFullDe, aod::pidTPCFullEl, aod::pidTOFFullEl>>();

photonDeuteronTable(
photon.pt(),
photon.eta(),
Expand All @@ -325,8 +328,7 @@ struct PhotonDeuteronCorrelation {
deltaPhi,
deltaEta,
invMass,
kpn
);
kpn);
}
}
}
Expand Down