Skip to content

Commit 83c6e4d

Browse files
lubynetslubynets
authored
[PWGHF] Lc decay time: fix definition and provide additional info (#17518)
1 parent 7dbf78e commit 83c6e4d

2 files changed

Lines changed: 83 additions & 60 deletions

File tree

‎PWGHF/D2H/Tasks/taskLc.cxx‎

Lines changed: 67 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -63,6 +63,7 @@ using namespace o2::framework::expressions;
6363
using namespace o2::hf_centrality;
6464
using namespace o2::hf_occupancy;
6565
using namespace o2::hf_evsel;
66+
using namespace o2::constants::physics;
6667

6768
/// Λc± → p± K∓ π± analysis task
6869
struct HfTaskLc {
@@ -74,7 +75,7 @@ struct HfTaskLc {
7475
Configurable<bool> fillTHn{"fillTHn", false, "fill THn"};
7576
Configurable<bool> storeOccupancy{"storeOccupancy", true, "Flag to store occupancy information"};
7677
Configurable<int> occEstimator{"occEstimator", 2, "Occupancy estimation (None: 0, ITS: 1, FT0C: 2)"};
77-
Configurable<bool> storeProperLifetime{"storeProperLifetime", false, "Flag to store proper lifetime"};
78+
Configurable<bool> storeProperDecayTime{"storeProperDecayTime", false, "Flag to store proper decay time"};
7879
// CCDB configuration
7980
Configurable<std::string> ccdbUrl{"ccdbUrl", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
8081
Configurable<std::string> ccdbPathGrp{"ccdbPathGrp", "GLO/GRP/GRP", "Path of the grp file (Run 2)"};
@@ -117,16 +118,9 @@ struct HfTaskLc {
117118
ConfigurableAxis thnConfigAxisGenPtB{"thnConfigAxisGenPtB", {1000, 0, 100}, "Gen Pt B"};
118119
ConfigurableAxis thnConfigAxisNumPvContr{"thnConfigAxisNumPvContr", {200, -0.5, 199.5}, "Number of PV contributors"};
119120
ConfigurableAxis thnConfigAxisOccupancy{"thnConfigAxisOccupancy", {14, 0, 14000}, "axis for centrality"};
120-
ConfigurableAxis thnConfigAxisProperLifetime{"thnConfigAxisProperLifetime", {200, 0, 2}, "Proper lifetime, ps"};
121+
ConfigurableAxis thnConfigAxisProperDecayTime{"thnConfigAxisProperDecayTime", {200, 0, 2}, "Proper decay time, ps"};
121122
HistogramRegistry registry{"registry", {}};
122123

123-
// Factors for conversion between units
124-
constexpr static float CtToProperLifetimePs = 1.f / o2::constants::physics::LightSpeedCm2PS;
125-
constexpr static float NanoToPico = 1000.f;
126-
// Names of folders and suffixes for MC signal histograms
127-
constexpr static std::string_view SignalFolders[] = {"signal", "prompt", "nonprompt"};
128-
constexpr static std::string_view SignalSuffixes[] = {"", "Prompt", "NonPrompt"};
129-
130124
enum MlClasses : int {
131125
MlClassBackground = 0,
132126
MlClassPrompt,
@@ -137,9 +131,16 @@ struct HfTaskLc {
137131
enum SignalClasses : int {
138132
Signal = 0,
139133
Prompt,
140-
NonPrompt
134+
NonPrompt,
135+
NumberOfSignalClasses
141136
};
142137

138+
// Factor for conversion between units
139+
constexpr static float CtToProperDecayTimePs = 1.f / o2::constants::physics::LightSpeedCm2PS;
140+
// Names of folders and suffixes for MC signal histograms
141+
constexpr static std::array<std::string_view, NumberOfSignalClasses> SignalFolders = {"signal", "prompt", "nonprompt"};
142+
constexpr static std::array<std::string_view, NumberOfSignalClasses> SignalSuffixes = {"", "Prompt", "NonPrompt"};
143+
143144
void init(InitContext&)
144145
{
145146
const std::array<bool, 12> doprocess{doprocessDataStd, doprocessDataStdWithFT0C, doprocessDataStdWithFT0M, doprocessDataWithMl, doprocessDataWithMlWithFT0C, doprocessDataWithMlWithFT0M, doprocessMcStd, doprocessMcStdWithFT0C, doprocessMcStdWithFT0M, doprocessMcWithMl, doprocessMcWithMlWithFT0C, doprocessMcWithMlWithFT0M};
@@ -186,8 +187,8 @@ struct HfTaskLc {
186187
addHistogramsRec("hDecLength", "decay length (cm)", "entries", {HistType::kTH1F, {{400, 0., 1.}}});
187188
/// decay length xy candidate
188189
addHistogramsRec("hDecLengthxy", "decay length xy (cm)", "entries", {HistType::kTH1F, {{400, 0., 1.}}});
189-
/// proper lifetime
190-
addHistogramsRec("hCt", "proper lifetime (#Lambda_{c}) * #it{c} (cm)", "entries", {HistType::kTH1F, {{100, 0., 0.2}}});
190+
/// proper decay time
191+
addHistogramsRec("hCt", "proper decay time (#Lambda_{c}) * #it{c} (cm)", "entries", {HistType::kTH1F, {{100, 0., 0.2}}});
191192
/// cosine of pointing angle
192193
addHistogramsRec("hCPA", "cosine of pointing angle", "entries", {HistType::kTH1F, {{110, -1.1, 1.1}}});
193194
/// cosine of pointing angle xy
@@ -219,8 +220,8 @@ struct HfTaskLc {
219220
/// decay length xy candidate
220221
addHistogramsRec("hDecLengthxyVsPt", "decay length xy (cm)", "#it{p}_{T} (GeV/#it{c})", {HistType::kTH2F, {{400, 0., 1.}, {vbins}}});
221222

222-
/// proper lifetime
223-
addHistogramsRec("hCtVsPt", "proper lifetime (#Lambda_{c}) * #it{c} (cm)", "#it{p}_{T} (GeV/#it{c})", {HistType::kTH2F, {{100, 0., 0.2}, {vbins}}});
223+
/// proper decay time
224+
addHistogramsRec("hCtVsPt", "proper decay time (#Lambda_{c}) * #it{c} (cm)", "#it{p}_{T} (GeV/#it{c})", {HistType::kTH2F, {{100, 0., 0.2}, {vbins}}});
224225

225226
/// cosine of pointing angle
226227
addHistogramsRec("hCPAVsPt", "cosine of pointing angle", "#it{p}_{T} (GeV/#it{c})", {HistType::kTH2F, {{110, -1.1, 1.1}, {vbins}}});
@@ -268,7 +269,8 @@ struct HfTaskLc {
268269
const AxisSpec thnAxisPtB{thnConfigAxisGenPtB, "#it{p}_{T}^{B} (GeV/#it{c})"};
269270
const AxisSpec thnAxisTracklets{thnConfigAxisNumPvContr, "Number of PV contributors"};
270271
const AxisSpec thnAxisOccupancy{thnConfigAxisOccupancy, "Occupancy"};
271-
const AxisSpec thnAxisProperLifetime{thnConfigAxisProperLifetime, "T_{proper} (ps)"};
272+
const AxisSpec thnAxisProperDecayTime{thnConfigAxisProperDecayTime, "#it{t}_{proper} (ps)"};
273+
const AxisSpec thnAxisProperDecayTimeGen{thnConfigAxisProperDecayTime, "#it{t}_{proper, gen} (ps)"};
272274

273275
bool const isDataWithMl = doprocessDataWithMl || doprocessDataWithMlWithFT0C || doprocessDataWithMlWithFT0M;
274276
bool const isMcWithMl = doprocessMcWithMl || doprocessMcWithMlWithFT0C || doprocessMcWithMlWithFT0M;
@@ -300,10 +302,13 @@ struct HfTaskLc {
300302
}
301303
}
302304
}
303-
if (storeProperLifetime) {
305+
if (storeProperDecayTime) {
304306
for (const auto& axes : std::array<std::vector<AxisSpec>*, 3>{&axesWithBdt, &axesStd, &axesGen}) {
305307
if (!axes->empty()) {
306-
axes->push_back(thnAxisProperLifetime);
308+
axes->push_back(thnAxisProperDecayTime);
309+
if (!isData && axes != &axesGen) {
310+
axes->push_back(thnAxisProperDecayTimeGen);
311+
}
307312
}
308313
}
309314
}
@@ -334,13 +339,34 @@ struct HfTaskLc {
334339
return o2::hf_centrality::getCentralityColl<Coll>(collision);
335340
}
336341

342+
/// Evaluate decay time of generated particle
343+
/// \param mcParticleProng0 one of generated particle's daughters
344+
/// \param motherParticle generated particle
345+
/// \return decay time in picoseconds. For nonprompt particles it is evaluated as if it was prompt (as it is calculated for data)
346+
float evaluateMcGenDecayTime(const McParticles3ProngMatched::iterator& mcParticleProng0, const McParticles3ProngMatched::iterator& motherParticle)
347+
{
348+
const auto mcCollision = motherParticle.template mcCollision_as<aod::McCollisions>();
349+
const float pMother = motherParticle.p();
350+
const float pvX = mcCollision.posX();
351+
const float pvY = mcCollision.posY();
352+
const float pvZ = mcCollision.posZ();
353+
const float svX = mcParticleProng0.vx();
354+
const float svY = mcParticleProng0.vy();
355+
const float svZ = mcParticleProng0.vz();
356+
357+
const float decayLength = static_cast<float>(RecoDecay::distance(std::array<float, 3>{svX, svY, svZ}, std::array<float, 3>{pvX, pvY, pvZ}));
358+
const float properDecayTime = decayLength * static_cast<float>(MassLambdaCPlus) / LightSpeedCm2PS / pMother;
359+
360+
return properDecayTime;
361+
}
362+
337363
/// Helper function for filling MC reconstructed histograms for prompt, nonpromt and common (signal)
338364
/// \param candidate is a reconstructed candidate
339365
/// \tparam SignalType is an enum defining which histogram in which folder (signal, prompt or nonpromt) to fill
340366
template <int SignalType, typename CandidateType>
341367
void fillHistogramsRecSig(CandidateType const& candidate)
342368
{
343-
const auto& mcParticleProng0 = candidate.template prong0_as<aod::TracksWMc>().template mcParticle_as<soa::Join<aod::McParticles, aod::HfCand3ProngMcGen>>();
369+
const auto& mcParticleProng0 = candidate.template prong0_as<aod::TracksWMc>().template mcParticle_as<McParticles3ProngMatched>();
344370
const auto pdgCodeProng0 = std::abs(mcParticleProng0.pdgCode());
345371
if ((candidate.isSelLcToPKPi() >= selectionFlagLc) && pdgCodeProng0 == kProton) {
346372
registry.fill(HIST("MC/reconstructed/") + HIST(SignalFolders[SignalType]) + HIST("/hMassRecSig") + HIST(SignalSuffixes[SignalType]), HfHelper::invMassLcToPKPi(candidate));
@@ -385,8 +411,8 @@ struct HfTaskLc {
385411

386412
/// Fill MC histograms at reconstruction level
387413
/// \tparam FillMl switch to fill ML histograms
388-
template <bool FillMl, typename CollType, typename CandLcMcRec, typename CandLcMcGen>
389-
void fillHistosMcRec(CollType const& collision, CandLcMcRec const& candidates, CandLcMcGen const& mcParticles)
414+
template <bool FillMl, typename CollType, typename CandLcMcRec>
415+
void fillHistosMcRec(CollType const& collision, CandLcMcRec const& candidates, McParticles3ProngMatched const& mcParticles)
390416
{
391417
const auto thisCollId = collision.globalIndex();
392418
const auto& groupedLcCandidates = candidates.sliceBy(candLcPerCollision, thisCollId);
@@ -403,7 +429,7 @@ struct HfTaskLc {
403429

404430
if (std::abs(candidate.flagMcMatchRec()) == hf_decay::hf_cand_3prong::DecayChannelMain::LcToPKPi) {
405431
// Get the corresponding MC particle.
406-
const auto& mcParticleProng0 = candidate.template prong0_as<aod::TracksWMc>().template mcParticle_as<soa::Join<aod::McParticles, aod::HfCand3ProngMcGen>>();
432+
const auto& mcParticleProng0 = candidate.template prong0_as<aod::TracksWMc>().template mcParticle_as<McParticles3ProngMatched>();
407433
const auto pdgCodeProng0 = std::abs(mcParticleProng0.pdgCode());
408434
const auto indexMother = RecoDecay::getMother(mcParticles, mcParticleProng0, o2::constants::physics::Pdg::kLambdaCPlus, true);
409435
const auto particleMother = mcParticles.rawIteratorAt(indexMother);
@@ -420,6 +446,8 @@ struct HfTaskLc {
420446
const auto numPvContributors = collision.numContrib();
421447
const auto ptRecB = candidate.ptBhadMotherPart();
422448

449+
const float properDecayTimeGen = evaluateMcGenDecayTime(mcParticleProng0, particleMother);
450+
423451
/// MC reconstructed signal
424452
fillHistogramsRecSig<Signal>(candidate);
425453

@@ -438,7 +466,7 @@ struct HfTaskLc {
438466
occ = o2::hf_occupancy::getOccupancyColl(collision, occEstimator);
439467
}
440468
double outputBkg(-1), outputPrompt(-1), outputFD(-1);
441-
const float properLifetime = HfHelper::ctLc(candidate) * CtToProperLifetimePs;
469+
const float properDecayTime = HfHelper::ctLc(candidate) * CtToProperDecayTimePs;
442470

443471
auto fillTHnRecSig = [&](bool isPKPi) {
444472
const auto massLc = isPKPi ? HfHelper::invMassLcToPKPi(candidate) : HfHelper::invMassLcToPiKP(candidate);
@@ -461,8 +489,9 @@ struct HfTaskLc {
461489
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
462490
valuesToFill.push_back(occ);
463491
}
464-
if (storeProperLifetime) {
465-
valuesToFill.push_back(properLifetime);
492+
if (storeProperDecayTime) {
493+
valuesToFill.push_back(properDecayTime);
494+
valuesToFill.push_back(properDecayTimeGen);
466495
}
467496
if constexpr (FillMl) {
468497
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(valuesToFill.data());
@@ -485,8 +514,8 @@ struct HfTaskLc {
485514
/// Helper function for filling MC generated histograms for prompt, nonpromt and common (signal)
486515
/// \param particle is a generated particle
487516
/// \tparam SignalType is an enum defining which histogram in which folder (signal, prompt or nonpromt) to fill
488-
template <int SignalType, typename ParticleType>
489-
void fillHistogramsGen(ParticleType const& particle)
517+
template <int SignalType>
518+
void fillHistogramsGen(McParticles3ProngMatched::iterator const& particle)
490519
{
491520
registry.fill(HIST("MC/generated/") + HIST(SignalFolders[SignalType]) + HIST("/hPtGen") + HIST(SignalSuffixes[SignalType]), particle.pt());
492521
registry.fill(HIST("MC/generated/") + HIST(SignalFolders[SignalType]) + HIST("/hEtaGen") + HIST(SignalSuffixes[SignalType]), particle.eta());
@@ -498,8 +527,8 @@ struct HfTaskLc {
498527
}
499528

500529
/// Fill MC histograms at generated level
501-
template <typename CandLcMcGen, typename Coll>
502-
void fillHistosMcGen(CandLcMcGen const& mcParticles, Coll const& recoCollisions)
530+
template <typename Coll>
531+
void fillHistosMcGen(McParticles3ProngMatched const& mcParticles, Coll const& recoCollisions)
503532
{
504533
// MC gen.
505534
for (const auto& particle : mcParticles) {
@@ -522,10 +551,9 @@ struct HfTaskLc {
522551
occ = o2::hf_occupancy::getOccupancyGenColl(recoCollsPerMcColl, occEstimator);
523552
}
524553

525-
const auto& mcDaughter0 = particle.template daughters_as<soa::Join<aod::McParticles, aod::HfCand3ProngMcGen>>().begin();
526-
const float p2m = particle.p() / o2::constants::physics::MassLambdaCPlus;
527-
const float gamma = std::sqrt(1 + p2m * p2m); // mother's particle Lorentz factor
528-
const float properLifetime = mcDaughter0.vt() * NanoToPico / gamma; // from ns to ps * from lab time to proper time
554+
const auto mcDaughter0 = particle.template daughters_as<McParticles3ProngMatched>().begin();
555+
556+
const float properDecayTime = evaluateMcGenDecayTime(mcDaughter0, particle);
529557

530558
fillHistogramsGen<Signal>(particle);
531559

@@ -537,8 +565,8 @@ struct HfTaskLc {
537565
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
538566
valuesToFill.push_back(occ);
539567
}
540-
if (storeProperLifetime) {
541-
valuesToFill.push_back(properLifetime);
568+
if (storeProperDecayTime) {
569+
valuesToFill.push_back(properDecayTime);
542570
}
543571
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(valuesToFill.data());
544572
}
@@ -631,7 +659,7 @@ struct HfTaskLc {
631659
occ = o2::hf_occupancy::getOccupancyColl(collision, occEstimator);
632660
}
633661
double outputBkg(-1), outputPrompt(-1), outputFD(-1);
634-
const float properLifetime = HfHelper::ctLc(candidate) * CtToProperLifetimePs;
662+
const float properDecayTime = HfHelper::ctLc(candidate) * CtToProperDecayTimePs;
635663

636664
auto fillTHnData = [&](bool isPKPi) {
637665
const auto massLc = isPKPi ? HfHelper::invMassLcToPKPi(candidate) : HfHelper::invMassLcToPiKP(candidate);
@@ -654,8 +682,8 @@ struct HfTaskLc {
654682
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
655683
valuesToFill.push_back(occ);
656684
}
657-
if (storeProperLifetime) {
658-
valuesToFill.push_back(properLifetime);
685+
if (storeProperDecayTime) {
686+
valuesToFill.push_back(properDecayTime);
659687
}
660688
if constexpr (FillMl) {
661689
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(valuesToFill.data());
@@ -687,10 +715,10 @@ struct HfTaskLc {
687715

688716
/// Run the analysis on MC data
689717
/// \tparam FillMl switch to fill ML histograms
690-
template <bool FillMl, typename CollType, typename CandType, typename CandLcMcGen>
718+
template <bool FillMl, typename CollType, typename CandType>
691719
void runAnalysisPerCollisionMc(CollType const& collisions,
692720
CandType const& candidates,
693-
CandLcMcGen const& mcParticles)
721+
McParticles3ProngMatched const& mcParticles)
694722
{
695723
for (const auto& collision : collisions) {
696724
// MC Rec.

0 commit comments

Comments
 (0)