Skip to content

Commit 3aea964

Browse files
deependra170598deependra170598
authored
PWGHF: Add multiplicity dependence v2 (#6957)
* Adding-Cent-Without-New-Table * Implimented suggestion by Anton Alkin on slicing and partition issue. In this commit no new table is created for MC Centrality * formatting * Assigning centrality to McCollision on the basis of Number of Contribution in PV reconstruction * formatting
1 parent e6f05f7 commit 3aea964

1 file changed

Lines changed: 70 additions & 4 deletions

File tree

‎PWGHF/D2H/Tasks/taskDstarToD0Pi.cxx‎

Lines changed: 70 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,7 @@
1717

1818
#include "CommonConstants/PhysicsConstants.h"
1919
#include "Framework/AnalysisTask.h"
20+
#include "Framework/ASoAHelpers.h"
2021
#include "Framework/runDataProcessing.h"
2122

2223
#include "PWGHF/Core/SelectorCuts.h"
@@ -27,6 +28,7 @@ using namespace o2;
2728
using namespace o2::analysis;
2829
using namespace o2::framework;
2930
using namespace o2::framework::expressions;
31+
using namespace o2::soa;
3032

3133
/// Dstar analysis task
3234
struct HfTaskDstarToD0Pi {
@@ -42,12 +44,19 @@ struct HfTaskDstarToD0Pi {
4244
using CandDstarWSelFlagMcRec = soa::Join<aod::HfD0FromDstar, aod::HfCandDstars, aod::HfSelDstarToD0Pi, aod::HfCandDstarMcRec>;
4345
using CandDstarMcGen = soa::Join<aod::McParticles, aod::HfCandDstarMcGen>;
4446

47+
using CollisionsWCent = soa::Join<aod::Collisions, aod::CentFT0Ms>;
48+
using CollisionsWCentMcLabel = soa::Join<aod::Collisions, aod::McCollisionLabels, aod::CentFT0Ms>;
49+
50+
PresliceUnsorted<CollisionsWCentMcLabel> colsPerMcCollision = aod::mccollisionlabel::mcCollisionId;
51+
SliceCache cache;
52+
4553
Partition<CandDstarWSelFlag> rowsSelectedCandDstar = aod::hf_sel_candidate_dstar::isSelDstarToD0Pi == selectionFlagDstarToD0Pi;
4654
Partition<CandDstarWSelFlagMcRec> rowsSelectedCandDstarMcRec = aod::hf_sel_candidate_dstar::isRecoD0Flag == selectionFlagHfD0ToPiK;
4755

4856
ConfigurableAxis binningImpactParam{"binningImpactParam", {1000, 0.1, -0.1}, " Bins of Impact Parameter"};
4957
ConfigurableAxis binningDecayLength{"binningDecayLength", {1000, 0.0, 0.7}, "Bins of Decay Length"};
5058
ConfigurableAxis binningNormDecayLength{"binningNormDecayLength", {1000, 0.0, 40.0}, "Bins of Normalised Decay Length"};
59+
ConfigurableAxis binningCentrality{"binningCentrality", {VARIABLE_WIDTH, 0.0, 10.0, 20.0, 30.0, 60.0, 100.0}, "centrality binning"};
5160

5261
HistogramRegistry registry{
5362
"registry",
@@ -57,7 +66,8 @@ struct HfTaskDstarToD0Pi {
5766
{"QA/hPtProng1D0", "Prong1 of D0 Candidates; Prong1 #it{p}_{T} (GeV/#it{c}); entries", {HistType::kTH1F, {{360, 0., 36.}}}},
5867
{"QA/hPtProng0D0Bar", "Prong0 of D0Bar Candidates; Prong0 #it{p}_{T} (GeV/#it{c}); entries", {HistType::kTH1F, {{360, 0., 36.}}}},
5968
{"QA/hPtProng1D0Bar", "Prong1 of D0Bar Candidates; Prong1 #it{p}_{T} (GeV/#it{c}); entries", {HistType::kTH1F, {{360, 0., 36.}}}},
60-
{"QA/hPtSoftPi", "Soft Pi ; Soft Pi #it{p}_{T} (GeV/#it{c}); entries", {HistType::kTH1F, {{100, 0., 1.}}}}}};
69+
{"QA/hPtSoftPi", "Soft Pi ; Soft Pi #it{p}_{T} (GeV/#it{c}); entries", {HistType::kTH1F, {{100, 0., 1.}}}},
70+
{"QA/hDeltaCentGen", "#{Delta}Cent % distribution of Collisions having same MC Collision;FT0M #{Delta}Cent %; Counts", {HistType::kTH1F, {{100, 0.0, 100.0}}}}}};
6171

6272
void init(InitContext&)
6373
{
@@ -66,7 +76,9 @@ struct HfTaskDstarToD0Pi {
6676
AxisSpec axisImpactParam = {binningImpactParam, "impact parameter (cm)"};
6777
AxisSpec axisDecayLength = {binningDecayLength, " decay length (cm)"};
6878
AxisSpec axisNormDecayLength = {binningNormDecayLength, "normalised decay length (cm)"};
79+
AxisSpec axisCentrality = {binningCentrality, "centrality (%)"};
6980

81+
registry.add("Yield/hDeltaInvMassDstar3D", "#Delta #it{M}_{inv} D* Candidate; inv. mass ((#pi #pi k) - (#pi k)) (GeV/#it{c}^{2});#it{p}_{T} (GeV/#it{c}); FT0M centrality", {HistType::kTH3F, {{100, 0.13, 0.16}, {vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {axisCentrality}}}, true);
7082
registry.add("Yield/hDeltaInvMassDstar2D", "#Delta #it{M}_{inv} D* Candidate; inv. mass ((#pi #pi k) - (#pi k)) (GeV/#it{c}^{2});#it{p}_{T} (GeV/#it{c})", {HistType::kTH2F, {{100, 0.13, 0.16}, {vecPtBins, "#it{p}_{T} (GeV/#it{c})"}}}, true);
7183
registry.add("Yield/hDeltaInvMassDstar1D", "#Delta #it{M}_{inv} D* Candidate; inv. mass ((#pi #pi k) - (#pi k)) (GeV/#it{c}^{2}); entries", {HistType::kTH1F, {{100, 0.13, 0.16}}}, true);
7284
registry.add("Yield/hInvMassDstar", "#Delta #it{M}_{inv} D* Candidate; inv. mass (#pi #pi k) (GeV/#it{c}^{2}); entries", {HistType::kTH1F, {{500, 0., 5.0}}}, true);
@@ -94,7 +106,8 @@ struct HfTaskDstarToD0Pi {
94106
registry.add("QA/hEtaSkimD0RecSig", "MC Matched Skimed D* Candidates at Reconstruction Level; #it{#eta} of D0 Prong", {HistType::kTH1F, {{100, -2., 2.}}});
95107
registry.add("QA/hEtaSkimDstarRecSig", "MC Matched Skimed D* Candidates at Reconstruction Level; #it{#eta} of D* Candidate", {HistType::kTH1F, {{100, -2., 2.}}});
96108
// pt vs y
97-
registry.add("QA/hPtSkimDstarGenSig", "MC Matched Skimed D* Reconstructed Candidates at Generator Level; #it{p}_{T} of D* at Generator Level (GeV/#it{c})", {HistType::kTH1F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}}});
109+
registry.add("QA/hPtSkimDstarGenSig", "MC Matched Skimed D* Reconstructed Candidates at Generator Level; #it{p}_{T} of D* at Generator Level (GeV/#it{c})", {HistType::kTH1F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}}}, true);
110+
registry.add("Efficiency/hPtVsCentSkimDstarGenSig", "MC Matched Skimed D* Reconstructed Candidates at Generator Level; #it{p}_{T} of D* at Generator Level (GeV/#it{c}); Centrality (%)", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {axisCentrality}}}, true);
98111
registry.add("QA/hPtVsYSkimDstarRecSig", "MC Matched Skimed D* Candidates at Reconstruction Level; #it{p}_{T} of D* at Reconstruction Level (GeV/#it{c}); #it{y}", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {100, -5., 5.}}});
99112
registry.add("QA/hPtVsYRecoTopolDstarRecSig", "MC Matched RecoTopol D* Candidates at Reconstruction Level; #it{p}_{T} of D* at Reconstruction Level (GeV/#it{c}); #it{y}", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {100, -5., 5.}}});
100113
registry.add("QA/hPtVsYRecoPidDstarRecSig", "MC Matched RecoPid D* Candidates at Reconstruction Level; #it{p}_{T} of D* at Reconstruction Level (GeV/#it{c}); #it{y}", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {100, -5., 5.}}});
@@ -118,6 +131,7 @@ struct HfTaskDstarToD0Pi {
118131
// MC Matching at Generator level Successful
119132
registry.add("QA/hEtaDstarGen", "MC Matched D* Candidates at Generator Level; #it{#eta}", {HistType::kTH1F, {{100, -2., 2.}}});
120133
registry.add("QA/hPtDstarGen", "MC Matched D* Candidates at Generator Level; #it{p}_{T} of D*", {HistType::kTH1F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}}});
134+
registry.add("Efficiency/hPtVsCentDstarGen", "MC Matched D* Candidates at Generator Level; #it{p}_{T} of D*;Centrality (%) ", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {axisCentrality}}}, true);
121135
registry.add("QA/hPtVsYDstarGen", "MC Matched D* Candidates at Generator Level; #it{p}_{T} of D*; #it{y}", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {100, -5., 5.}}});
122136
// Prompt Gen
123137
registry.add("QA/hPtPromptDstarGen", "MC Matched Prompt D* Candidates at Generator Level; #it{p}_{T} of D*", {HistType::kTH1F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}}});
@@ -127,7 +141,7 @@ struct HfTaskDstarToD0Pi {
127141
registry.add("QA/hPtVsYNonPromptDstarGen", "MC Matched Non-Prompt D* Candidates at Generator Level; #it{p}_{T} of D*; #it{y}", {HistType::kTH2F, {{vecPtBins, "#it{p}_{T} (GeV/#it{c})"}, {100, -5., 5.}}});
128142
}
129143

130-
void process(CandDstarWSelFlag const&)
144+
void process(CollisionsWCent const&, CandDstarWSelFlag const&)
131145
{
132146
for (const auto& candDstar : rowsSelectedCandDstar) {
133147
auto yDstar = candDstar.y(constants::physics::MassDStar);
@@ -162,8 +176,12 @@ struct HfTaskDstarToD0Pi {
162176
auto invD0 = candDstar.invMassD0();
163177
auto invD0Bar = candDstar.invMassD0Bar();
164178

179+
auto collision = candDstar.collision_as<CollisionsWCent>();
180+
auto centrality = collision.centFT0M(); // 0-100%
181+
165182
auto signDstar = candDstar.signSoftPi();
166183
if (signDstar > 0) {
184+
registry.fill(HIST("Yield/hDeltaInvMassDstar3D"), (invDstar - invD0), candDstar.pt(), centrality);
167185
registry.fill(HIST("Yield/hDeltaInvMassDstar2D"), (invDstar - invD0), candDstar.pt());
168186
registry.fill(HIST("Yield/hInvMassD0"), invD0, candDstar.ptD0());
169187
registry.fill(HIST("Yield/hDeltaInvMassDstar1D"), (invDstar - invD0));
@@ -172,6 +190,7 @@ struct HfTaskDstarToD0Pi {
172190
registry.fill(HIST("QA/hPtProng0D0"), candDstar.ptProng0());
173191
registry.fill(HIST("QA/hPtProng1D0"), candDstar.ptProng1());
174192
} else if (signDstar < 0) {
193+
registry.fill(HIST("Yield/hDeltaInvMassDstar3D"), (invAntiDstar - invD0Bar), candDstar.pt(), centrality);
175194
registry.fill(HIST("Yield/hDeltaInvMassDstar2D"), (invAntiDstar - invD0Bar), candDstar.pt());
176195
registry.fill(HIST("Yield/hInvMassD0"), invD0Bar, candDstar.ptD0());
177196
registry.fill(HIST("Yield/hDeltaInvMassDstar1D"), (invAntiDstar - invD0Bar));
@@ -183,10 +202,11 @@ struct HfTaskDstarToD0Pi {
183202
}
184203
}
185204

186-
void processMC(CandDstarWSelFlagMcRec const&,
205+
void processMC(aod::McCollisions const&, CollisionsWCentMcLabel const& collisions, CandDstarWSelFlagMcRec const&,
187206
CandDstarMcGen const& rowsMcPartilces,
188207
aod::TracksWMc const&)
189208
{
209+
rowsSelectedCandDstarMcRec.bindExternalIndices(&collisions);
190210

191211
int8_t signDstar = 0;
192212
// MC at Reconstruction level
@@ -196,11 +216,18 @@ struct HfTaskDstarToD0Pi {
196216
if (yCandDstarRecoMax >= 0. && std::abs(yDstarRecSig) > yCandDstarRecoMax) {
197217
continue;
198218
}
219+
auto collision = candDstarMcRec.collision_as<CollisionsWCentMcLabel>();
220+
auto centrality = collision.centFT0M(); // 0-100%
199221
if (TESTBIT(std::abs(candDstarMcRec.flagMcMatchRec()), aod::hf_cand_dstar::DecayType::DstarToD0Pi)) { // if MC matching is successful at Reconstruction Level
200222
// get MC Mother particle
201223
auto indexMother = RecoDecay::getMother(rowsMcPartilces, candDstarMcRec.prong0_as<aod::TracksWMc>().mcParticle_as<CandDstarMcGen>(), o2::constants::physics::Pdg::kDStar, true, &signDstar, 2);
202224
auto particleMother = rowsMcPartilces.rawIteratorAt(indexMother); // What is difference between rawIterator() or iteratorAt() methods?
203225
registry.fill(HIST("QA/hPtSkimDstarGenSig"), particleMother.pt()); // generator level pt
226+
registry.fill(HIST("Efficiency/hPtVsCentSkimDstarGenSig"), particleMother.pt(), centrality);
227+
228+
// auto recCollision = candDstarMcRec.collision_as<CollisionsWCentMcLabel>();
229+
// float centFT0M = recCollision.centFT0M();
230+
// LOGF(info, "centFT0M: %f", centFT0M);
204231

205232
registry.fill(HIST("QA/hPtVsYSkimDstarRecSig"), ptDstarRecSig, yDstarRecSig); // Skimed at level of trackIndexSkimCreator
206233
if (candDstarMcRec.isRecoTopol()) { // if Topological selection are passed
@@ -253,15 +280,48 @@ struct HfTaskDstarToD0Pi {
253280
// MC at Generator Level
254281
for (const auto& mcParticle : rowsMcPartilces) {
255282
if (TESTBIT(std::abs(mcParticle.flagMcMatchGen()), aod::hf_cand_dstar::DecayType::DstarToD0Pi)) { // MC Matching is successful at Generator Level
283+
256284
auto ptGen = mcParticle.pt();
257285
// auto yGen = mcParticle.y(); // Can we use this definition?
258286
auto yGen = RecoDecay::y(mcParticle.pVector(), o2::constants::physics::MassDStar);
259287
if (yCandDstarGenMax >= 0. && std::abs(yGen) > yCandDstarGenMax) {
260288
continue;
261289
}
290+
291+
auto mcCollision = mcParticle.mcCollision_as<aod::McCollisions>();
292+
auto recCollisions = collisions.sliceBy(colsPerMcCollision, mcCollision.globalIndex());
293+
// auto recCollisions = collisions.sliceByCached(aod::mccollisionlabel::mcCollisionId, mcCollision.globalIndex(), cache);
294+
// auto recCollisions = collisions.sliceByCachedUnsorted(aod::mccollisionlabel::mcCollisionId, mcCollision.globalIndex(), cache);
295+
296+
// looking if a generated collision reconstructed more than a times.
297+
if (recCollisions.size() > 1) {
298+
for (const auto& [c1, c2] : combinations(CombinationsStrictlyUpperIndexPolicy(recCollisions, recCollisions))) {
299+
auto deltaCent = abs(c1.centFT0M() - c2.centFT0M());
300+
registry.fill(HIST("QA/hDeltaCentGen"), deltaCent);
301+
}
302+
}
303+
304+
float centFT0MGen;
305+
// assigning centrality to MC Collision using max FT0M amplitute from Reconstructed collisions
306+
if (recCollisions.size()) {
307+
std::vector<std::pair<soa::Filtered<CollisionsWCentMcLabel>::iterator, int>> tempRecCols;
308+
309+
for (const auto& recCol : recCollisions) {
310+
// if(recCollisions.size()>1) LOGF(info, "cuurent cent: %f",recCol.centFT0M());
311+
tempRecCols.push_back(std::make_pair(recCol, recCol.numContrib()));
312+
}
313+
std::sort(tempRecCols.begin(), tempRecCols.end(), compare);
314+
centFT0MGen = tempRecCols.at(0).first.centFT0M();
315+
// if(recCollisions.size()>1) LOGF(info, "assigned cent: %f",centFT0MGen);
316+
} else {
317+
centFT0MGen = -999.;
318+
}
319+
262320
registry.fill(HIST("QA/hEtaDstarGen"), mcParticle.eta());
263321
registry.fill(HIST("QA/hPtDstarGen"), ptGen);
264322
registry.fill(HIST("QA/hPtVsYDstarGen"), ptGen, yGen);
323+
registry.fill(HIST("Efficiency/hPtVsCentDstarGen"), ptGen, centFT0MGen);
324+
265325
// only promt Dstar candidate at Generator level
266326
if (mcParticle.originMcGen() == RecoDecay::OriginType::Prompt) {
267327
registry.fill(HIST("QA/hPtPromptDstarGen"), ptGen);
@@ -274,6 +334,12 @@ struct HfTaskDstarToD0Pi {
274334
}
275335
}
276336
PROCESS_SWITCH(HfTaskDstarToD0Pi, processMC, "Process MC Data", false);
337+
338+
// Comparator function to sort based on the second argument of a tuple
339+
static bool compare(const std::pair<soa::Filtered<CollisionsWCentMcLabel>::iterator, int>& a, const std::pair<soa::Filtered<CollisionsWCentMcLabel>::iterator, int>& b)
340+
{
341+
return a.second > b.second;
342+
}
277343
};
278344

279345
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)

0 commit comments

Comments
 (0)