Skip to content

Commit

Permalink
[PWGLF] Initial commit for table producer to select events with at le…
Browse files Browse the repository at this point in the history
…ast two phi meson (AliceO2Group#9071)
  • Loading branch information
skundu692 authored Dec 20, 2024
1 parent ccf411b commit d53c7d8
Show file tree
Hide file tree
Showing 3 changed files with 374 additions and 0 deletions.
96 changes: 96 additions & 0 deletions PWGLF/DataModel/ReducedDoublePhiTables.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,96 @@
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.

///
/// \author Sourav Kundu <[email protected]>

#ifndef PWGLF_DATAMODEL_REDUCEDDOUBLEPHITABLES_H_
#define PWGLF_DATAMODEL_REDUCEDDOUBLEPHITABLES_H_

#include <cmath>

#include "Common/DataModel/Centrality.h"
#include "Common/DataModel/Multiplicity.h"
#include "Common/DataModel/PIDResponse.h"
#include "Common/Core/RecoDecay.h"
#include "Common/DataModel/TrackSelectionTables.h"
#include "Framework/AnalysisDataModel.h"
#include "Framework/ASoA.h"

namespace o2::aod
{
namespace redphievent
{
DECLARE_SOA_COLUMN(NumPos, numPos, int); //! Number of positive Kaon
DECLARE_SOA_COLUMN(NumNeg, numNeg, int); //! Number of negative Kaon
} // namespace redphievent
DECLARE_SOA_TABLE(RedPhiEvents, "AOD", "REDPHIEVENT",
o2::soa::Index<>,
bc::GlobalBC,
bc::RunNumber,
timestamp::Timestamp,
collision::PosZ,
collision::NumContrib,
redphievent::NumPos,
redphievent::NumNeg);
using RedPhiEvent = RedPhiEvents::iterator;

namespace phitrack
{
DECLARE_SOA_INDEX_COLUMN(RedPhiEvent, redPhiEvent);
DECLARE_SOA_COLUMN(PhiPx, phiPx, float); //! Phi Px
DECLARE_SOA_COLUMN(PhiPy, phiPy, float); //! Phi Py
DECLARE_SOA_COLUMN(PhiPz, phiPz, float); //! Phi Pz
DECLARE_SOA_COLUMN(Phid1Px, phid1Px, float); //! Phi d1 Px
DECLARE_SOA_COLUMN(Phid1Py, phid1Py, float); //! Phi d1 Py
DECLARE_SOA_COLUMN(Phid1Pz, phid1Pz, float); //! Phi d1 Pz
DECLARE_SOA_COLUMN(Phid2Px, phid2Px, float); //! Phi d2 Px
DECLARE_SOA_COLUMN(Phid2Py, phid2Py, float); //! Phi d2 Py
DECLARE_SOA_COLUMN(Phid2Pz, phid2Pz, float); //! Phi d2 Pz
DECLARE_SOA_COLUMN(PhiMass, phiMass, float); //! Phi Mass
DECLARE_SOA_COLUMN(Phid1Index, phid1Index, int64_t); //! Phi d1 index
DECLARE_SOA_COLUMN(Phid2Index, phid2Index, int64_t); //! Phi d2 index
DECLARE_SOA_COLUMN(Phid1Charge, phid1Charge, float); //! Phi d1 charge
DECLARE_SOA_COLUMN(Phid2Charge, phid2Charge, float); //! Phi d1 charge
DECLARE_SOA_COLUMN(Phid1TPC, phid1TPC, float); //! TPC nsigma d1
DECLARE_SOA_COLUMN(Phid2TPC, phid2TPC, float); //! TPC nsigma d2
DECLARE_SOA_COLUMN(Phid1TOFHit, phid1TOFHit, int); //! TOF hit d1
DECLARE_SOA_COLUMN(Phid2TOFHit, phid2TOFHit, int); //! TOF hit d2
DECLARE_SOA_COLUMN(Phid1TOF, phid1TOF, float); //! TOF nsigma d1
DECLARE_SOA_COLUMN(Phid2TOF, phid2TOF, float); //! TOF nsigma d2
} // namespace phitrack
DECLARE_SOA_TABLE(PhiTracks, "AOD", "PHITRACK",
o2::soa::Index<>,
phitrack::RedPhiEventId,
phitrack::PhiPx,
phitrack::PhiPy,
phitrack::PhiPz,
phitrack::Phid1Px,
phitrack::Phid1Py,
phitrack::Phid1Pz,
phitrack::Phid2Px,
phitrack::Phid2Py,
phitrack::Phid2Pz,
phitrack::PhiMass,
phitrack::Phid1Index,
phitrack::Phid2Index,
phitrack::Phid1Charge,
phitrack::Phid2Charge,
phitrack::Phid1TPC,
phitrack::Phid2TPC,
phitrack::Phid1TOFHit,
phitrack::Phid2TOFHit,
phitrack::Phid1TOF,
phitrack::Phid2TOF);

using PhiTrack = PhiTracks::iterator;
} // namespace o2::aod
#endif // PWGLF_DATAMODEL_REDUCEDDOUBLEPHITABLES_H_
5 changes: 5 additions & 0 deletions PWGLF/TableProducer/Resonances/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,11 @@ o2physics_add_dpl_workflow(f1protonreducedtable
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::DetectorsVertexing O2Physics::EventFilteringUtils
COMPONENT_NAME Analysis)

o2physics_add_dpl_workflow(doublephitable
SOURCES doublephitable.cxx
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::DetectorsVertexing
COMPONENT_NAME Analysis)

o2physics_add_dpl_workflow(filterf1proton
SOURCES filterf1proton.cxx
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::DetectorsVertexing
Expand Down
273 changes: 273 additions & 0 deletions PWGLF/TableProducer/Resonances/doublephitable.cxx
Original file line number Diff line number Diff line change
@@ -0,0 +1,273 @@
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.

/// \file doublephitable.cxx
/// \brief Selection of events with triplets and pairs for femtoscopic studies
///
/// \author Sourav Kundu, [email protected]

#include <Framework/Configurable.h>
#include <Math/GenVector/Boost.h>
#include <Math/Vector4D.h>
#include <TMath.h>
#include <fairlogger/Logger.h>
#include <TDatabasePDG.h> // FIXME
#include <TPDGCode.h> // FIXME

#include <iostream>
#include <iterator>
#include <string>
#include <vector>

#include "EventFiltering/Zorro.h"
#include "EventFiltering/ZorroSummary.h"

#include "PWGLF/DataModel/ReducedDoublePhiTables.h"
#include "Framework/ASoAHelpers.h"
#include "Framework/AnalysisDataModel.h"
#include "Framework/AnalysisTask.h"
#include "Framework/HistogramRegistry.h"
#include "Framework/runDataProcessing.h"
#include "CommonConstants/MathConstants.h"
#include "Common/Core/TrackSelection.h"
#include "Common/DataModel/TrackSelectionTables.h"
#include "Common/DataModel/EventSelection.h"
#include "Common/DataModel/Multiplicity.h"
#include "Common/DataModel/PIDResponse.h"
#include "DataFormatsTPC/BetheBlochAleph.h"
#include "CCDB/BasicCCDBManager.h"
#include "CCDB/CcdbApi.h"

using namespace o2;
using namespace o2::framework;
using namespace o2::framework::expressions;

struct doublephitable {

// Produce derived tables
Produces<aod::RedPhiEvents> redPhiEvents;
Produces<aod::PhiTracks> phiTrack;

// events
Configurable<float> cfgCutVertex{"cfgCutVertex", 10.0f, "Accepted z-vertex range"};
// track
Configurable<bool> useGlobalTrack{"useGlobalTrack", true, "use Global track"};
Configurable<float> cfgCutTOFBeta{"cfgCutTOFBeta", 0.0, "cut TOF beta"};
Configurable<float> cfgCutCharge{"cfgCutCharge", 0.0, "cut on Charge"};
Configurable<float> cfgCutPT{"cfgCutPT", 0.2, "PT cut on daughter track"};
Configurable<float> cfgCutEta{"cfgCutEta", 0.8, "Eta cut on daughter track"};
Configurable<float> cfgCutDCAxy{"cfgCutDCAxy", 2.0f, "DCAxy range for tracks"};
Configurable<float> cfgCutDCAz{"cfgCutDCAz", 2.0f, "DCAz range for tracks"};
Configurable<float> nsigmaCutTPC{"nsigmacutTPC", 3.0, "Value of the TPC Nsigma cut"};
Configurable<float> nsigmaCutTOF{"nsigmaCutTOF", 3.0, "Value of the TOF Nsigma cut"};
Configurable<int> cfgITScluster{"cfgITScluster", 0, "Number of ITS cluster"};
Configurable<int> cfgTPCcluster{"cfgTPCcluster", 70, "Number of TPC cluster"};
Configurable<bool> isDeepAngle{"isDeepAngle", true, "Deep Angle cut"};
Configurable<double> cfgDeepAngle{"cfgDeepAngle", 0.04, "Deep Angle cut value"};
ConfigurableAxis configThnAxisInvMass{"configThnAxisInvMass", {120, 0.98, 1.1}, "#it{M} (GeV/#it{c}^{2})"};
ConfigurableAxis configThnAxisPt{"configThnAxisPt", {100, 0.0, 10.}, "#it{p}_{T} (GeV/#it{c})"};

Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex;
Filter acceptanceFilter = (nabs(aod::track::eta) < cfgCutEta && nabs(aod::track::pt) > cfgCutPT);
Filter DCAcutFilter = (nabs(aod::track::dcaXY) < cfgCutDCAxy) && (nabs(aod::track::dcaZ) < cfgCutDCAz);
Filter PIDcutFilter = nabs(aod::pidtpc::tpcNSigmaKa) < nsigmaCutTPC;

using EventCandidates = soa::Filtered<soa::Join<aod::Collisions, aod::EvSels, aod::Mults>>;
using TrackCandidates = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTOFbeta, aod::pidTPCFullKa, aod::pidTOFFullKa>>;

SliceCache cache;
Partition<TrackCandidates> posTracks = aod::track::signed1Pt > cfgCutCharge;
Partition<TrackCandidates> negTracks = aod::track::signed1Pt < cfgCutCharge;

// Histogram
HistogramRegistry qaRegistry{"QAHistos", {
{"hEventstat", "hEventstat", {HistType::kTH1F, {{2, 0.0f, 2.0f}}}},
{"hInvMassPhi", "hInvMassPhi", {HistType::kTH2F, {{40, 1.0f, 1.04f}, {100, 0.0f, 10.0f}}}},
{"hNsigmaPtkaonTPC", "hNsigmaPtkaonTPC", {HistType::kTH2F, {{200, -10.0f, 10.0f}, {100, 0.0f, 10.0f}}}},
{"hNsigmaPtkaonTOF", "hNsigmaPtkaonTOF", {HistType::kTH2F, {{200, -10.0f, 10.0f}, {100, 0.0f, 10.0f}}}},
},
OutputObjHandlingPolicy::AnalysisObject};

double massKa = o2::constants::physics::MassKPlus;

template <typename T>
bool selectionTrack(const T& candidate)
{
if (useGlobalTrack && !(candidate.isGlobalTrack() && candidate.isPVContributor() && candidate.itsNCls() > cfgITScluster && candidate.tpcNClsFound() > cfgTPCcluster)) {
return false;
}
return true;
}
template <typename T>
bool selectionPID(const T& candidate)
{
if (!candidate.hasTOF() && TMath::Abs(candidate.tpcNSigmaKa()) < nsigmaCutTPC) {
return true;
}
if (candidate.hasTOF() && candidate.beta() > cfgCutTOFBeta && TMath::Abs(candidate.tpcNSigmaKa()) < nsigmaCutTPC && TMath::Abs(candidate.tofNSigmaKa()) < nsigmaCutTOF) {
return true;
}
return false;
}
// deep angle cut on pair to remove photon conversion
template <typename T1, typename T2>
bool selectionPair(const T1& candidate1, const T2& candidate2)
{
double pt1, pt2, pz1, pz2, p1, p2, angle;
pt1 = candidate1.pt();
pt2 = candidate2.pt();
pz1 = candidate1.pz();
pz2 = candidate2.pz();
p1 = candidate1.p();
p2 = candidate2.p();
angle = TMath::ACos((pt1 * pt2 + pz1 * pz2) / (p1 * p2));
if (isDeepAngle && angle < cfgDeepAngle) {
return false;
}
return true;
}

ROOT::Math::PxPyPzMVector KaonPlus, KaonMinus, PhiMesonMother, PhiVectorDummy, Phid1dummy, Phid2dummy;
void processPhiReducedTable(EventCandidates::iterator const& collision, TrackCandidates const&, aod::BCsWithTimestamps const&)
{
bool keepEventDoublePhi = false;
int numberPhi = 0;
auto currentRunNumber = collision.bc_as<aod::BCsWithTimestamps>().runNumber();
auto bc = collision.bc_as<aod::BCsWithTimestamps>();
std::vector<int64_t> Phid1Index = {};
std::vector<int64_t> Phid2Index = {};

std::vector<float> Phid1Charge = {};
std::vector<float> Phid2Charge = {};

std::vector<float> Phid1TPC = {};
std::vector<float> Phid2TPC = {};

std::vector<float> Phid1TOF = {};
std::vector<float> Phid2TOF = {};

std::vector<int> Phid1TOFHit = {};
std::vector<int> Phid2TOFHit = {};

std::vector<ROOT::Math::PtEtaPhiMVector> phiresonance, phiresonanced1, phiresonanced2;

int Npostrack = 0;
int Nnegtrack = 0;

if (collision.sel8() && collision.selection_bit(aod::evsel::kNoTimeFrameBorder) && collision.selection_bit(aod::evsel::kNoITSROFrameBorder) && collision.selection_bit(aod::evsel::kNoSameBunchPileup) && collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) {
auto posThisColl = posTracks->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache);
auto negThisColl = negTracks->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache);
for (auto track1 : posThisColl) {
// track selection
if (!selectionTrack(track1)) {
continue;
}
// PID check
if (!selectionPID(track1)) {
continue;
}
Npostrack = Npostrack + 1;
qaRegistry.fill(HIST("hNsigmaPtkaonTPC"), track1.tpcNSigmaKa(), track1.pt());
if (track1.hasTOF()) {
qaRegistry.fill(HIST("hNsigmaPtkaonTOF"), track1.tofNSigmaKa(), track1.pt());
}
auto track1ID = track1.globalIndex();
for (auto track2 : negThisColl) {
// track selection
if (!selectionTrack(track2)) {
continue;
}
// PID check
if (!selectionPID(track2)) {
continue;
}
if (Npostrack == 1) {
Nnegtrack = Nnegtrack + 1;
}
auto track2ID = track2.globalIndex();
if (track2ID == track1ID) {
continue;
}
if (!selectionPair(track1, track2)) {
continue;
}
KaonPlus = ROOT::Math::PxPyPzMVector(track1.px(), track1.py(), track1.pz(), massKa);
KaonMinus = ROOT::Math::PxPyPzMVector(track2.px(), track2.py(), track2.pz(), massKa);
PhiMesonMother = KaonPlus + KaonMinus;
if (PhiMesonMother.M() > 1.0 && PhiMesonMother.M() < 1.04) {
numberPhi = numberPhi + 1;
ROOT::Math::PtEtaPhiMVector temp1(track1.pt(), track1.eta(), track1.phi(), massKa);
ROOT::Math::PtEtaPhiMVector temp2(track2.pt(), track2.eta(), track2.phi(), massKa);
ROOT::Math::PtEtaPhiMVector temp3(PhiMesonMother.pt(), PhiMesonMother.eta(), PhiMesonMother.phi(), PhiMesonMother.M());

phiresonanced1.push_back(temp1);
phiresonanced2.push_back(temp2);
phiresonance.push_back(temp3);

Phid1Index.push_back(track1.globalIndex());
Phid2Index.push_back(track2.globalIndex());

Phid1Charge.push_back(track1.sign());
Phid2Charge.push_back(track2.sign());

Phid1TPC.push_back(track1.tpcNSigmaKa());
Phid2TPC.push_back(track2.tpcNSigmaKa());
auto d1TOFHit = -1;
auto d2TOFHit = -1;
auto d1TOF = -999.0;
auto d2TOF = -999.0;
if (track1.hasTOF()) {
d1TOFHit = 1;
d1TOF = track1.tofNSigmaKa();
}
if (track2.hasTOF()) {
d2TOFHit = 1;
d2TOF = track2.tofNSigmaKa();
}
Phid1TOF.push_back(d1TOF);
Phid2TOF.push_back(d2TOF);
Phid1TOFHit.push_back(d1TOFHit);
Phid2TOFHit.push_back(d2TOFHit);
qaRegistry.fill(HIST("hInvMassPhi"), PhiMesonMother.M(), PhiMesonMother.Pt());
}
}
}
} // select collision
if (numberPhi >= 2) {
keepEventDoublePhi = true;
}
qaRegistry.fill(HIST("hEventstat"), 0.5);
if (keepEventDoublePhi && numberPhi >= 2 && (phiresonance.size() == phiresonanced1.size()) && (phiresonance.size() == phiresonanced2.size())) {
qaRegistry.fill(HIST("hEventstat"), 1.5);
/////////// Fill collision table///////////////
redPhiEvents(bc.globalBC(), currentRunNumber, bc.timestamp(), collision.posZ(), collision.numContrib(), Npostrack, Nnegtrack);
auto indexEvent = redPhiEvents.lastIndex();
//// Fill track table for Phi//////////////////
for (auto if1 = phiresonance.begin(); if1 != phiresonance.end(); ++if1) {
auto i5 = std::distance(phiresonance.begin(), if1);
PhiVectorDummy = phiresonance.at(i5);
Phid1dummy = phiresonanced1.at(i5);
Phid2dummy = phiresonanced2.at(i5);
phiTrack(indexEvent, PhiVectorDummy.Px(), PhiVectorDummy.Py(), PhiVectorDummy.Pz(), Phid1dummy.Px(), Phid1dummy.Py(), Phid1dummy.Pz(), Phid2dummy.Px(), Phid2dummy.Py(), Phid2dummy.Pz(),
PhiVectorDummy.M(),
Phid1Index.at(i5), Phid2Index.at(i5),
Phid1Charge.at(i5), Phid2Charge.at(i5),
Phid1TPC.at(i5), Phid2TPC.at(i5), Phid1TOFHit.at(i5), Phid2TOFHit.at(i5), Phid1TOF.at(i5), Phid2TOF.at(i5));
}
}
} // process
PROCESS_SWITCH(doublephitable, processPhiReducedTable, "Process table creation for double phi", true);
};

WorkflowSpec defineDataProcessing(ConfigContext const& cfg)
{
return WorkflowSpec{adaptAnalysisTask<doublephitable>(cfg)};
}

0 comments on commit d53c7d8

Please sign in to comment.