Skip to content

Commit aed038f

Browse files
apalascianoapalasciano
authored
[PWGHF] Fix DCA sel and derived data in hiddenCharm workflow (#16134)
1 parent 47f3014 commit aed038f

2 files changed

Lines changed: 80 additions & 39 deletions

File tree

‎PWGHF/D2H/TableProducer/dataCreatorHiddenCharmReduced.cxx‎

Lines changed: 58 additions & 23 deletions
Original file line numberDiff line numberDiff line change
@@ -76,9 +76,10 @@ struct HfDataCreatorHiddenCharmReduced {
7676
Produces<aod::HcSelTracks> hfTrackLite;
7777

7878
struct : ConfigurableGroup {
79-
// track quality
8079
Configurable<bool> fillHistograms{"fillHistograms", true, "Fill proton QA histograms"};
80+
Configurable<bool> fillSparses{"fillSparses", false, "Flag to enable sparse filling"};
8181
Configurable<bool> selectProtons{"selectProtons", true, "Select protons"};
82+
// track quality
8283
Configurable<int> itsNClsMin{"itsNClsMin", 5, "Minimum number of ITS clusters"};
8384
Configurable<int> tpcNClsFoundMin{"tpcNClsFoundMin", 50, "Minimum number of found TPC clusters"};
8485
Configurable<int> tpcNClsCrossedRowsMin{"tpcNClsCrossedRowsMin", 80, "Minimum number of crossed TPC rows"};
@@ -90,6 +91,7 @@ struct HfDataCreatorHiddenCharmReduced {
9091
Configurable<std::vector<double>> binsPtTrack{"binsPtTrack", std::vector<double>{hf_cuts_single_track::vecBinsPtTrack}, "Track pT bin limits for DCA cuts"};
9192
Configurable<LabeledArray<double>> cutsTrack{"cutsTrack", {hf_cuts_single_track::CutsTrack[0], hf_cuts_single_track::NBinsPtTrack, hf_cuts_single_track::NCutVarsTrack, hf_cuts_single_track::labelsPtTrack, hf_cuts_single_track::labelsCutVarTrack}, "Single-track DCA selections per pT bin"};
9293
// DCA
94+
Configurable<bool> applyPtDependentDcaTrackSelection{"applyPtDependentDcaTrackSelection", false, "Flag to enable pT-dependent DCAxy and DCAz selection on single tracks"};
9395
Configurable<std::vector<float>> paramsDCAxyPtDep{"paramsDCAxyPtDep", std::vector<float>{0.0010, 0.0080, 0.73}, "Parameters for pT-dependent DCAxy cut: p0, p1, p2 for cut = p0 + p1/pt^p2"};
9496
Configurable<std::vector<float>> paramsDCAzPtDep{"paramsDCAzPtDep", std::vector<float>{-0.0044, 0.0152, 0.47}, "Parameters for pT-dependent DCAz cut: p0, p1, p2 for cut = p0 + p1/pt^p2"};
9597
// PID
@@ -147,17 +149,30 @@ struct HfDataCreatorHiddenCharmReduced {
147149
const AxisSpec axisPt{360, 0., 36., "#it{p}_{T}^{proton} (GeV/#it{c})"};
148150
const AxisSpec axisEta{100, -1., 1., "#eta"};
149151
const AxisSpec axisDca{400, -2., 2., "DCA_{xy} to primary vertex (cm)"};
150-
const AxisSpec axisNSigma{100, -5., 5., "n#sigma"};
152+
const AxisSpec axisNSigmaTPC{100, -5., 5., "n#sigma_{TPC}"};
153+
const AxisSpec axisNSigmaTOF{100, -5., 5., "n#sigma_{TOF}"};
154+
const AxisSpec axisMass{100, 2.85, 3.25, "p#overline{p} invariant mass (GeV/#it{c}^{2})"};
151155

152156
registry.add("hPzVtx", "Z position of primary vertex for selected tracks;z_{vtx} (cm);entries", {HistType::kTH1D, {AxisSpec{200, -20., 20., "z_{vtx} (cm)"}}});
153157
registry.add("hPtNoCuts", "All associated tracks;#it{p}_{T}^{track} (GeV/#it{c});entries", {HistType::kTH1D, {axisPt}});
154158
registry.add("hPtCutsProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});entries", {HistType::kTH1D, {axisPt}});
155159
registry.add("hEtaCutsProton", "Selected proton tracks;#eta;entries", {HistType::kTH1D, {axisEta}});
156160
registry.add("hDCAToPrimXYVsPtCutsProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});DCA_{xy} to primary vertex (cm)", {HistType::kTH2D, {axisPt, axisDca}});
157-
registry.add("hNSigmaTPCProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});n#sigma_{TPC}", {HistType::kTH2D, {axisPt, axisNSigma}});
158-
registry.add("hNSigmaTOFProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});n#sigma_{TOF}", {HistType::kTH2D, {axisPt, axisNSigma}});
159-
registry.add("hInvMass", "Invariant mass of selected proton with all other tracks in the event;#it{p}_{T}^{proton} (GeV/#it{c});invariant mass with other tracks (GeV/#it{c}^{2})", {HistType::kTH2D, {axisPt, AxisSpec{100, 2.85, 3.25, "invariant mass with other tracks (GeV/#it{c}^{2})"}}});
161+
registry.add("hNSigmaTPCProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});n#sigma_{TPC}^{p}", {HistType::kTH2D, {axisPt, axisNSigmaTPC}});
162+
registry.add("hNSigmaTOFProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});n#sigma_{TOF}^{p}", {HistType::kTH2D, {axisPt, axisNSigmaTOF}});
163+
registry.add("hNSigmaTPCKaon", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); n#sigma_{TPC}^{K}", {HistType::kTH2D, {axisPt, axisNSigmaTPC}});
164+
registry.add("hNSigmaTOFKaon", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); n#sigma_{TOF}^{K}", {HistType::kTH2D, {axisPt, axisNSigmaTOF}});
165+
registry.add("hNSigmaTPCPion", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); n#sigma_{TPC}^{#pi}", {HistType::kTH2D, {axisPt, axisNSigmaTPC}});
166+
registry.add("hNSigmaTOFPion", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); n#sigma_{TOF}^{#pi}", {HistType::kTH2D, {axisPt, axisNSigmaTOF}});
167+
registry.add("hTpcNClsFound", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); TpcNClsFound", {HistType::kTH2D, {axisPt, {200, 0., 200, "TpcNClsFound"}}});
168+
registry.add("hTpcNClsCrossedRows", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); TpcNClsCrossedRows", {HistType::kTH2D, {axisPt, {200, 0., 200, "TpcNClsCrossedRows"}}});
169+
registry.add("hTpcChi2NCl", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); TpcChi2NCl", {HistType::kTH2D, {axisPt, {10, 0., 10, "TpcChi2NCl"}}});
170+
registry.add("hItsNCls", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});ItsNCls", {HistType::kTH2D, {axisPt, {7, 0, 7, "ItsNCls"}}});
171+
registry.add("hItsChi2NCl", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c}); ItsChi2NCl", {HistType::kTH2D, {axisPt, {10, 0., 10, "ItsChi2NCl"}}});
172+
registry.add("hInvMass", "Invariant mass of selected proton with all other tracks in the event;#it{p}_{T}^{proton} (GeV/#it{c});", {HistType::kTH2D, {axisPt, axisMass}});
173+
registry.add("hInvMassSparse", "Invariant mass of selected proton with all other tracks in the event;#it{p}_{T}^{proton} (GeV/#it{c});", {HistType::kTH2D, {axisPt, axisMass}});
160174
registry.add("hDeDxTPCProton", "Selected proton tracks;#it{p}_{T}^{track} (GeV/#it{c});TPC dE/dx (a.u.)", {HistType::kTH2D, {axisPt, AxisSpec{100, 0., 200., "TPC dE/dx (a.u.)"}}});
175+
registry.add("hMassPtCutVars", "Charmonia candidates;#it{M} (p#overline{p}) (GeV/#it{c}^{2});#it{p}_{T}(cc) (GeV/#it{c}); n#sigma_{TPC}^{p1}; n#sigma_{TOF}^{p1}; dca_{xy}^{p1}; n#sigma_{TPC}^{p2}; n#sigma_{TOF}^{p2}; dca_{xy}^{p2}", {HistType::kTHnSparseF, {axisMass, axisPt, axisNSigmaTPC, axisNSigmaTOF, axisDca, axisNSigmaTPC, axisNSigmaTOF, axisDca}});
161176
}
162177

163178
// init HF event selection helper
@@ -266,8 +281,10 @@ struct HfDataCreatorHiddenCharmReduced {
266281
if (!isSelectedTrackDca(config.binsPtTrack, config.cutsTrack, track.pt(), track.dcaXY(), track.dcaZ())) {
267282
return false;
268283
}
269-
if (dcaSigma(track.pt(), config.paramsDCAxyPtDep.value[0], config.paramsDCAxyPtDep.value[1], config.paramsDCAxyPtDep.value[2]) > std::abs(dcaXY) || dcaSigma(track.pt(), config.paramsDCAzPtDep.value[0], config.paramsDCAzPtDep.value[1], config.paramsDCAzPtDep.value[2]) > std::abs(dcaZ)) {
270-
return false;
284+
if (config.applyPtDependentDcaTrackSelection.value) {
285+
if (std::abs(dcaXY) > dcaSigma(track.pt(), config.paramsDCAxyPtDep.value[0], config.paramsDCAxyPtDep.value[1], config.paramsDCAxyPtDep.value[2]) || std::abs(dcaZ) > dcaSigma(track.pt(), config.paramsDCAzPtDep.value[0], config.paramsDCAzPtDep.value[1], config.paramsDCAzPtDep.value[2])) {
286+
return false;
287+
}
271288
}
272289
return isSelectedPid(track);
273290
}
@@ -304,37 +321,55 @@ struct HfDataCreatorHiddenCharmReduced {
304321
if (!isSelectedTrack(trk)) {
305322
continue;
306323
}
307-
std::array pVecProton{trk.pVector()};
308-
hfTrackLite(trk.globalIndex(), collision.globalIndex(), pVecProton[0], pVecProton[1], pVecProton[2], trk.sign(), static_cast<uint8_t>(TrackType::Proton));
309-
selectedTrackIds.push_back(trk.globalIndex());
324+
selectedTrackIds.push_back(trkId.trackId());
325+
310326
if (config.fillHistograms) {
311327
registry.fill(HIST("hPtCutsProton"), trk.pt());
312328
registry.fill(HIST("hEtaCutsProton"), trk.eta());
313329
registry.fill(HIST("hDCAToPrimXYVsPtCutsProton"), trk.pt(), trk.dcaXY());
314330
registry.fill(HIST("hNSigmaTPCProton"), trk.pt(), trk.tpcNSigmaPr());
315-
registry.fill(HIST("hDeDxTPCProton"), trk.pt(), trk.tpcSignal());
331+
registry.fill(HIST("hNSigmaTPCPion"), trk.pt(), trk.tpcNSigmaPi());
332+
registry.fill(HIST("hNSigmaTPCKaon"), trk.pt(), trk.tpcNSigmaKa());
316333
if (trk.hasTOF()) {
317334
registry.fill(HIST("hNSigmaTOFProton"), trk.pt(), trk.tofNSigmaPr());
335+
registry.fill(HIST("hNSigmaTOFPion"), trk.pt(), trk.tofNSigmaPi());
336+
registry.fill(HIST("hNSigmaTOFKaon"), trk.pt(), trk.tofNSigmaKa());
318337
}
338+
registry.fill(HIST("hTpcNClsFound"), trk.pt(), trk.tpcNClsFound());
339+
registry.fill(HIST("hTpcNClsCrossedRows"), trk.pt(), trk.tpcNClsCrossedRows());
340+
registry.fill(HIST("hTpcChi2NCl"), trk.pt(), trk.tpcChi2NCl());
341+
registry.fill(HIST("hItsNCls"), trk.pt(), trk.itsNCls());
342+
registry.fill(HIST("hItsChi2NCl"), trk.pt(), trk.itsChi2NCl());
319343
}
320344
}
345+
// skip event if not at least a pair of tracks skimming selection
321346
if (selectedTrackIds.size() < NDaughtersCharmMeson) {
322347
continue;
323348
}
324349
hfReducedCollision(collision.posX(), collision.posY(), collision.posZ(), collision.numContrib(), hfRejMap, bz);
325-
for (size_t i = 0; i < selectedTrackIds.size(); ++i) {
326-
auto t1 = tracks.rawIteratorAt(selectedTrackIds[i]);
327-
std::array pVec1{t1.pVector()};
328-
for (size_t j = i + 1; j < selectedTrackIds.size(); ++j) {
329-
auto t2 = tracks.rawIteratorAt(selectedTrackIds[j]);
330-
if (t1.sign() * t2.sign() > 0) {
331-
continue;
332-
}
333-
std::array pVec2{t2.pVector()};
334-
float invMass = RecoDecay::m(std::array{pVec1, pVec2}, std::array{o2::constants::physics::MassProton, o2::constants::physics::MassProton});
335-
float ptEtac = RecoDecay::pt(RecoDecay::sumOfVec(pVec1, pVec2));
336-
if (config.fillHistograms) {
350+
351+
auto reducedEventId = hfReducedCollision.lastIndex();
352+
for (const auto& trkId : selectedTrackIds) {
353+
auto track = tracks.rawIteratorAt(trkId);
354+
std::array pVec{track.pVector()};
355+
hfTrackLite(track.globalIndex(), reducedEventId, pVec[0], pVec[1], pVec[2], track.sign(), static_cast<uint8_t>(TrackType::Proton));
356+
}
357+
358+
if (config.fillHistograms.value) {
359+
for (size_t i = 0; i < selectedTrackIds.size(); ++i) {
360+
auto t1 = tracks.rawIteratorAt(selectedTrackIds[i]);
361+
std::array pVec1{t1.pVector()};
362+
for (size_t j = i + 1; j < selectedTrackIds.size(); ++j) {
363+
auto t2 = tracks.rawIteratorAt(selectedTrackIds[j]);
364+
if (t1.sign() * t2.sign() > 0)
365+
continue;
366+
std::array pVec2{t2.pVector()};
367+
float invMass = RecoDecay::m(std::array{pVec1, pVec2}, std::array{o2::constants::physics::MassProton, o2::constants::physics::MassProton});
368+
float ptEtac = RecoDecay::pt(RecoDecay::sumOfVec(pVec1, pVec2));
337369
registry.fill(HIST("hInvMass"), ptEtac, invMass);
370+
if (config.fillSparses.value) {
371+
registry.fill(HIST("hMassPtCutVars"), invMass, ptEtac, t1.tpcNSigmaPr(), t1.tofNSigmaPr(), t1.dcaXY(), t2.tpcNSigmaPr(), t2.tofNSigmaPr(), t2.dcaXY());
372+
}
338373
}
339374
}
340375
}

‎PWGHF/D2H/Tasks/taskHiddenCharm.cxx‎

Lines changed: 22 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -52,8 +52,8 @@ struct HfTaskHiddenCharm {
5252
Configurable<int> centEstimator{"centEstimator", 0, "Centrality estimation (None: 0, FT0A: 1, FT0C: 2, FT0M: 3, FV0A: 4, NTracksPV: 5, FT0CVariant2: 6)"};
5353
Configurable<float> centralityMin{"centralityMin", 0.f, "Minimum accepted centrality"};
5454
Configurable<float> centralityMax{"centralityMax", 100.f, "Maximum accepted centrality"};
55-
Configurable<bool> fillOnlyUnlikeSign{"fillOnlyUnlikeSign", true, "Fill only unlike-sign proton pairs"};
56-
Configurable<bool> fillOnlyLikeSign{"fillOnlyLikeSign", true, "Fill only like-sign proton pairs"};
55+
Configurable<bool> fillUnlikeSign{"fillUnlikeSign", true, "Fill unlike-sign proton pairs"};
56+
Configurable<bool> fillLikeSign{"fillLikeSign", true, "Fill like-sign proton pairs"};
5757

5858
SliceCache cache;
5959

@@ -94,34 +94,40 @@ struct HfTaskHiddenCharm {
9494
void fillEtac(TCollisions const& collision,
9595
TProtonIds const& protonIds)
9696
{
97-
float cent{-1.f};
97+
float cent{0.f};
9898
if constexpr (CentEstimator != o2::hf_centrality::CentralityEstimator::None) {
9999
cent = o2::hf_centrality::getCentralityColl(collision, centEstimator);
100100
if (cent < centralityMin || cent >= centralityMax) {
101101
return; // skip events outside the centrality range
102102
}
103103
}
104104

105-
for (const auto& proton1 : protonIds) {
106-
for (const auto& proton2 : protonIds) {
107-
if (proton1.trackId() >= proton2.trackId()) {
108-
continue; // avoid double counting and self-pairs
109-
}
105+
for (auto it1 = protonIds.begin(); it1 != protonIds.end(); ++it1) {
106+
for (auto it2 = it1 + 1; it2 != protonIds.end(); ++it2) {
107+
108+
const auto& proton1 = *it1;
109+
const auto& proton2 = *it2;
110110
const int sign = (proton1.sign() * proton2.sign() > 0) ? 1 : -1;
111-
if ((sign == 1 && !fillOnlyLikeSign) || (sign == -1 && !fillOnlyUnlikeSign)) {
111+
112+
if (sign == 1 && !fillLikeSign)
112113
continue;
113-
}
114+
if (sign == -1 && !fillUnlikeSign)
115+
continue;
116+
114117
std::array<float, 3> pVec1{proton1.px(), proton1.py(), proton1.pz()};
115118
std::array<float, 3> pVec2{proton2.px(), proton2.py(), proton2.pz()};
119+
116120
float invMass = RecoDecay::m(std::array{pVec1, pVec2}, std::array{o2::constants::physics::MassProton, o2::constants::physics::MassProton});
117-
float ptEtac = RecoDecay::pt(RecoDecay::sumOfVec(pVec1, pVec2));
118-
registry.fill(HIST("hSparseHiddenCharm"), invMass, ptEtac, sign, cent);
121+
float pt = RecoDecay::pt(RecoDecay::sumOfVec(pVec1, pVec2));
122+
123+
registry.fill(HIST("hSparseHiddenCharm"), invMass, pt, sign, cent);
119124
if (sign == 1) {
120-
registry.fill(HIST("hPtVsInvMassLikeSign"), invMass, ptEtac);
121-
} else if (sign == -1) {
122-
registry.fill(HIST("hPtVsInvMassUnlikeSign"), invMass, ptEtac);
123-
registry.fill(HIST("hPtVsInvMassAllSign"), invMass, ptEtac);
125+
registry.fill(HIST("hPtVsInvMassLikeSign"), invMass, pt);
126+
} else {
127+
registry.fill(HIST("hPtVsInvMassUnlikeSign"), invMass, pt);
124128
}
129+
130+
registry.fill(HIST("hPtVsInvMassAllSign"), invMass, pt);
125131
}
126132
}
127133
}

0 commit comments

Comments
 (0)