Skip to content

Commit 96ef627

Browse files
author
Maxim Virta
committed
[PWGCF] Added multiplicity study & pt correlations
1 parent c8c6eb7 commit 96ef627

1 file changed

Lines changed: 105 additions & 9 deletions

File tree

‎PWGCF/GenericFramework/Tasks/flowGfwV02.cxx‎

Lines changed: 105 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -92,7 +92,7 @@ float philow = 0.0;
9292
float phiup = o2::constants::math::TwoPI;
9393
int nchbins = 300;
9494
float nchlow = 0;
95-
float nchup = 3000;
95+
float nchup = 300;
9696
std::vector<double> centbinning(90);
9797
int nBootstrap = 10;
9898
std::vector<std::pair<double, double>> etagapsPtPt;
@@ -104,6 +104,7 @@ std::vector<double> multGlobalPVCorrCutPars;
104104
} // namespace o2::analysis::gfw
105105

106106
struct FlowGfwV02 {
107+
107108
O2_DEFINE_CONFIGURABLE(cfgNbootstrap, int, 10, "Number of subsamples")
108109
O2_DEFINE_CONFIGURABLE(cfgMpar, int, 4, "Highest order of pt-pt correlations")
109110
O2_DEFINE_CONFIGURABLE(cfgCentEstimator, int, 0, "0:FT0C; 1:FT0CVariant1; 2:FT0M; 3:FT0A")
@@ -124,6 +125,7 @@ struct FlowGfwV02 {
124125
O2_DEFINE_CONFIGURABLE(cfgNormalizeByCharged, bool, true, "Enable or disable the normalization by charged particles");
125126
O2_DEFINE_CONFIGURABLE(cfgConsistentEventFlag, int, 15, "Flag for consistent event selection");
126127
O2_DEFINE_CONFIGURABLE(cfgMultCut, bool, true, "Use additional event cut on mult correlations");
128+
O2_DEFINE_CONFIGURABLE(cfgUseV0, bool, false, "Use V0 analysis");
127129

128130
// Event selection cuts
129131
struct : ConfigurableGroup {
@@ -139,6 +141,16 @@ struct FlowGfwV02 {
139141
O2_DEFINE_CONFIGURABLE(cfgIsVertexITSTPC, bool, true, "kIsVertexITSTPC - Selects collisions with at least one ITS-TPC track");
140142
} cfgEventCutFlags;
141143

144+
// Event selection cuts
145+
struct : ConfigurableGroup {
146+
O2_DEFINE_CONFIGURABLE(cfgEtaSubAMin, float, -0.8, "Minimum eta for subevent A");
147+
O2_DEFINE_CONFIGURABLE(cfgEtaSubAMax, float, -0.5, "Maximum eta for subevent A");
148+
O2_DEFINE_CONFIGURABLE(cfgEtaSubBMin, float, 0.5, "Minimum eta for subevent B");
149+
O2_DEFINE_CONFIGURABLE(cfgEtaSubBMax, float, 0.8, "Maximum eta for subevent B");
150+
O2_DEFINE_CONFIGURABLE(cfgEtaSubCMin, float, -0.4, "Minimum eta for subevent C");
151+
O2_DEFINE_CONFIGURABLE(cfgEtaSubCMax, float, 0.4, "Maximum eta for subevent C");
152+
} cfgSubeventCuts;
153+
142154
struct : ConfigurableGroup {
143155
Configurable<std::vector<double>> cfgMultGlobalCutPars{"cfgMultGlobalCutPars", std::vector<double>{2272.16, -76.6932, 1.01204, -0.00631545, 1.59868e-05, 136.336, -4.97006, 0.121199, -0.0015921, 7.66197e-06}, "Global vs FT0C multiplicity cut parameter values"};
144156
Configurable<std::vector<double>> cfgMultPVCutPars{"cfgMultPVCutPars", std::vector<double>{3074.43, -106.192, 1.46176, -0.00968364, 2.61923e-05, 182.128, -7.43492, 0.193901, -0.00256715, 1.22594e-05}, "PV vs FT0C multiplicity cut parameter values"};
@@ -253,6 +265,8 @@ struct FlowGfwV02 {
253265
std::array<float, 6> itsNsigmaCut;
254266
std::array<float, 6> tpcNsigmaCut;
255267
std::array<std::unique_ptr<TH1D>, 4> hPtMid{};
268+
std::array<std::unique_ptr<TH1D>, 4> hPtForward{};
269+
std::array<std::unique_ptr<TH1D>, 4> hPtBackward{};
256270
};
257271
PIDState pidStates;
258272

@@ -275,7 +289,8 @@ struct FlowGfwV02 {
275289
PidCharged = 0,
276290
PidPions,
277291
PidKaons,
278-
PidProtons
292+
PidProtons,
293+
PidTotal
279294
};
280295
enum PiKpArrayIndex {
281296
IndPionUp = 0,
@@ -378,6 +393,25 @@ struct FlowGfwV02 {
378393
pidStates.hPtMid[PidKaons]->SetDirectory(nullptr);
379394
pidStates.hPtMid[PidProtons]->SetDirectory(nullptr);
380395

396+
pidStates.hPtForward[PidCharged] = std::make_unique<TH1D>("hPtForward_charged", "hPtForward_charged", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
397+
pidStates.hPtForward[PidPions] = std::make_unique<TH1D>("hPtForward_pions", "hPtForward_pions", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
398+
pidStates.hPtForward[PidKaons] = std::make_unique<TH1D>("hPtForward_kaons", "hPtForward_kaons", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
399+
pidStates.hPtForward[PidProtons] = std::make_unique<TH1D>("hPtForward_protons", "hPtForward_protons", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
400+
pidStates.hPtForward[PidCharged]->SetDirectory(nullptr);
401+
pidStates.hPtForward[PidPions]->SetDirectory(nullptr);
402+
pidStates.hPtForward[PidKaons]->SetDirectory(nullptr);
403+
pidStates.hPtForward[PidProtons]->SetDirectory(nullptr);
404+
405+
pidStates.hPtBackward[PidCharged] = std::make_unique<TH1D>("hPtBackward_charged", "hPtBackward_charged", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
406+
pidStates.hPtBackward[PidPions] = std::make_unique<TH1D>("hPtBackward_pions", "hPtBackward_pions", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
407+
pidStates.hPtBackward[PidKaons] = std::make_unique<TH1D>("hPtBackward_kaons", "hPtBackward_kaons", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
408+
pidStates.hPtBackward[PidProtons] = std::make_unique<TH1D>("hPtBackward_protons", "hPtBackward_protons", o2::analysis::gfw::ptbinning.size() - 1, &o2::analysis::gfw::ptbinning[0]);
409+
pidStates.hPtBackward[PidCharged]->SetDirectory(nullptr);
410+
pidStates.hPtBackward[PidPions]->SetDirectory(nullptr);
411+
pidStates.hPtBackward[PidKaons]->SetDirectory(nullptr);
412+
pidStates.hPtBackward[PidProtons]->SetDirectory(nullptr);
413+
414+
381415
AxisSpec phiAxis = {o2::analysis::gfw::phibins, o2::analysis::gfw::philow, o2::analysis::gfw::phiup, "#phi"};
382416
AxisSpec etaAxis = {o2::analysis::gfw::etabins, -cfgTrackCuts.cfgEtaMax, cfgTrackCuts.cfgEtaMax, "#eta"};
383417
AxisSpec vtxAxis = {o2::analysis::gfw::vtxZbins, -cfgEventCuts.cfgZvtxMax, cfgEventCuts.cfgZvtxMax, "Vtx_{z} (cm)"};
@@ -398,9 +432,20 @@ struct FlowGfwV02 {
398432
AxisSpec multpvAxis = {600, 0, 600, "N_{ch} (PV)"};
399433
AxisSpec dcaZAxis = {200, -2, 2, "DCA_{z} (cm)"};
400434
AxisSpec dcaXYAxis = {200, -0.5, 0.5, "DCA_{xy} (cm)"};
435+
AxisSpec pidAxis = {4, -0.5, 3.5, "PID"}; // 0 = not identified, 1 = pion, 2 = kaon, 3 = proton
401436

402437
registry.add("v02pt", "", {HistType::kTProfile2D, {ptAxis, centAxis}});
403-
registry.add("nchMid", "", {HistType::kTProfile2D, {ptAxis, centAxis}});
438+
registry.add("nchMid", "", {HistType::kTProfile3D, {ptAxis, centAxis, nchAxis}});
439+
registry.add("v02centmult", "", {HistType::kTProfile2D, {centAxis, nchAxis}});
440+
441+
442+
registry.add("analysis/v0AB", "", {HistType::kTProfile3D, {pidAxis,ptAxis, centAxis}});
443+
registry.add("analysis/v0BA", "", {HistType::kTProfile3D, {pidAxis,ptAxis, centAxis}});
444+
registry.add("analysis/nchA", "", {HistType::kTProfile3D, {pidAxis,ptAxis, centAxis}});
445+
registry.add("analysis/nchB", "", {HistType::kTProfile3D, {pidAxis,ptAxis, centAxis}});
446+
registry.add("analysis/ptA", "", {HistType::kTProfile3D, {pidAxis, centAxis, nchAxis}});
447+
registry.add("analysis/ptB", "", {HistType::kTProfile3D, {pidAxis, centAxis, nchAxis}});
448+
registry.add("analysis/ptAB", "", {HistType::kTProfile3D, {pidAxis, centAxis, nchAxis}});
404449

405450
ccdb->setURL("http://alice-ccdb.cern.ch");
406451
ccdb->setCaching(true);
@@ -614,7 +659,7 @@ struct FlowGfwV02 {
614659
}
615660
if (cfgPIDEfficiency) {
616661
const std::array<std::string, 4> pidStrings = {"ch", "pi", "ka", "pr"};
617-
for (int i = 1; i < 4; i++) {
662+
for (int i = 1; i < PidTotal; i++) {
618663

619664
cfg.mEfficiency[i] = ccdb->getForTimeStamp<TH1D>(cfgEfficiency.value + pidStrings[i], timestamp);
620665
if (cfg.mEfficiency[i] == nullptr) {
@@ -806,7 +851,7 @@ struct FlowGfwV02 {
806851
}
807852

808853
template <DataType dt>
809-
void fillOutputContainers(const float& centmult, const double& rndm, const int& /*run*/ = 0)
854+
void fillOutputContainers(const float& centmult, const int& multiplicity, const double& rndm, const int& /*run*/ = 0)
810855
{
811856
for (uint l_ind = 0; l_ind < corrconfigs.size(); ++l_ind) {
812857
if (!corrconfigs.at(l_ind).pTDif) {
@@ -846,6 +891,32 @@ struct FlowGfwV02 {
846891
}
847892
}
848893
}
894+
895+
if (cfgUseV0) {
896+
double v0corrAB = 0;
897+
double v0corrBA = 0;
898+
double ptMeanForward = pidStates.hPtForward[PidCharged]->GetMean();
899+
double ptMeanBackward = pidStates.hPtBackward[PidCharged]->GetMean();
900+
double ptFractionForward = 0.;
901+
double ptFractionBackward = 0.;
902+
for (int pid = 0; pid < PidTotal; pid++) {
903+
int normIndex = (cfgNormalizeByCharged) ? PidCharged : pid;
904+
for (int i = 1; i <= fSecondAxis->GetNbins(); i++) {
905+
ptFractionForward = pidStates.hPtForward[pid]->GetBinContent(i) / pidStates.hPtForward[normIndex]->Integral();
906+
ptFractionBackward = pidStates.hPtBackward[pid]->GetBinContent(i) / pidStates.hPtBackward[normIndex]->Integral();
907+
v0corrAB = ptFractionForward * ptMeanBackward;
908+
v0corrBA = ptFractionBackward * ptMeanForward;
909+
registry.fill(HIST("analysis/v0AB"), pid, fSecondAxis->GetBinCenter(i), centmult, v0corrAB);
910+
registry.fill(HIST("analysis/v0BA"), pid, fSecondAxis->GetBinCenter(i), centmult, v0corrBA);
911+
registry.fill(HIST("analysis/nchA"), pid, fSecondAxis->GetBinCenter(i), centmult, ptFractionForward);
912+
registry.fill(HIST("analysis/nchB"), pid, fSecondAxis->GetBinCenter(i), centmult, ptFractionBackward);
913+
}
914+
registry.fill(HIST("analysis/ptA"), pid, centmult, multiplicity, ptMeanForward);
915+
registry.fill(HIST("analysis/ptB"), pid, centmult, multiplicity, ptMeanBackward);
916+
registry.fill(HIST("analysis/ptAB"), pid, centmult, multiplicity, ptMeanForward * ptMeanBackward);
917+
}
918+
}
919+
849920
// Fill the profiles for each pT bin
850921
auto dnx = fGFW->Calculate(corrconfigs.at(0), 0, kTRUE).real();
851922
if (dnx == 0)
@@ -857,9 +928,10 @@ struct FlowGfwV02 {
857928
ptFraction = pidStates.hPtMid[PidCharged]->GetBinContent(i) / pidStates.hPtMid[PidCharged]->Integral();
858929
if (std::abs(val) < 1)
859930
registry.fill(HIST("v02pt"), fSecondAxis->GetBinCenter(i), centmult, val * ptFraction, (cfgUseMultiplicityFlowWeights) ? dnx : 1.0);
860-
registry.fill(HIST("nchMid"), fSecondAxis->GetBinCenter(i), centmult, ptFraction);
931+
registry.fill(HIST("nchMid"), fSecondAxis->GetBinCenter(i), centmult, multiplicity, ptFraction);
861932
}
862933
}
934+
registry.fill(HIST("v02centmult"), centmult, multiplicity, val);
863935
return;
864936
}
865937

@@ -891,23 +963,47 @@ struct FlowGfwV02 {
891963
pidStates.hPtMid[PidPions]->Reset();
892964
pidStates.hPtMid[PidKaons]->Reset();
893965
pidStates.hPtMid[PidProtons]->Reset();
966+
pidStates.hPtBackward[PidCharged]->Reset();
967+
pidStates.hPtBackward[PidPions]->Reset();
968+
pidStates.hPtBackward[PidKaons]->Reset();
969+
pidStates.hPtBackward[PidProtons]->Reset();
970+
pidStates.hPtForward[PidCharged]->Reset();
971+
pidStates.hPtForward[PidPions]->Reset();
972+
pidStates.hPtForward[PidKaons]->Reset();
973+
pidStates.hPtForward[PidProtons]->Reset();
894974

895975
float lRandom = fRndm->Rndm();
896976

897977
// Loop over tracks and check if they are accepted
898978
AcceptedTracks acceptedTracks{0, 0, 0, 0};
899979
for (const auto& track : tracks) {
900980
processTrack(track, vtxz, xaxis.multiplicity, run, acceptedTracks);
901-
if (track.eta() > -0.4 && track.eta() < 0.4)
981+
if (track.eta() > cfgSubeventCuts.cfgEtaSubCMin && track.eta() < cfgSubeventCuts.cfgEtaSubCMax)
902982
pidStates.hPtMid[PidCharged]->Fill(track.pt(), getEfficiency(track, PidCharged));
983+
if (track.eta() > cfgSubeventCuts.cfgEtaSubAMin && track.eta() < cfgSubeventCuts.cfgEtaSubAMax) // add mean pT
984+
pidStates.hPtBackward[PidCharged]->Fill(track.pt(), getEfficiency(track, PidCharged));
985+
if (track.eta() > cfgSubeventCuts.cfgEtaSubBMin && track.eta() < cfgSubeventCuts.cfgEtaSubBMax) // add mean pT
986+
pidStates.hPtForward[PidCharged]->Fill(track.pt(), getEfficiency(track, PidCharged));
903987
// If PID is identified, fill pt spectrum for the corresponding particle
904988
int pidInd = getNsigmaPID(track);
905-
if (pidInd != -1 && track.eta() > -0.4 && track.eta() < 0.4) {
989+
if (pidInd != -1 && track.eta() > cfgSubeventCuts.cfgEtaSubCMin && track.eta() < cfgSubeventCuts.cfgEtaSubCMax) {
906990
if (cfgPIDEfficiency)
907991
pidStates.hPtMid[pidInd]->Fill(track.pt(), getEfficiency(track, pidInd));
908992
else
909993
pidStates.hPtMid[pidInd]->Fill(track.pt(), getEfficiency(track, PidCharged)); // Default to charged particles if PID efficiency is not used
910994
}
995+
if (pidInd != -1 && track.eta() > cfgSubeventCuts.cfgEtaSubAMin && track.eta() < cfgSubeventCuts.cfgEtaSubAMax) {
996+
if (cfgPIDEfficiency)
997+
pidStates.hPtBackward[pidInd]->Fill(track.pt(), getEfficiency(track, pidInd));
998+
else
999+
pidStates.hPtBackward[pidInd]->Fill(track.pt(), getEfficiency(track, PidCharged)); // Default to charged particles if PID efficiency is not used
1000+
}
1001+
if (pidInd != -1 && track.eta() > cfgSubeventCuts.cfgEtaSubBMin && track.eta() < cfgSubeventCuts.cfgEtaSubBMax) {
1002+
if (cfgPIDEfficiency)
1003+
pidStates.hPtForward[pidInd]->Fill(track.pt(), getEfficiency(track, pidInd));
1004+
else
1005+
pidStates.hPtForward[pidInd]->Fill(track.pt(), getEfficiency(track, PidCharged)); // Default to charged particles if PID efficiency is not used
1006+
}
9111007
}
9121008
if (cfgConsistentEventFlag & 1)
9131009
if (!acceptedTracks.nPos || !acceptedTracks.nNeg)
@@ -922,7 +1018,7 @@ struct FlowGfwV02 {
9221018
if (acceptedTracks.nPos < 2 || acceptedTracks.nMid < 2 || acceptedTracks.nNeg < 2) // o2-linter: disable=magic-number (at least two tracks in all three subevents)
9231019
return;
9241020
// Fill output containers
925-
fillOutputContainers<dt>(xaxis.centrality, lRandom, run);
1021+
fillOutputContainers<dt>(xaxis.centrality, xaxis.multiplicity, lRandom, run);
9261022
}
9271023

9281024
template <typename TTrack>

0 commit comments

Comments
 (0)