From 3eb0b2f79b1c0461d3f33ea9246a2a17f4fd4475 Mon Sep 17 00:00:00 2001 From: wuctlby Date: Mon, 17 Feb 2025 12:25:24 +0100 Subject: [PATCH 1/3] resolution with occupancy --- PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx | 71 +++++++++++++++++------- 1 file changed, 51 insertions(+), 20 deletions(-) diff --git a/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx b/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx index a154b416b2b..857d32156b5 100644 --- a/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx +++ b/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx @@ -65,6 +65,7 @@ struct HfTaskFlowCharmHadrons { Configurable centralityMax{"centralityMax", 100., "Maximum centrality accepted in SP/EP computation (not applied in resolution process)"}; Configurable storeEP{"storeEP", false, "Flag to store EP-related axis"}; Configurable storeMl{"storeMl", false, "Flag to store ML scores"}; + Configurable storeResoOccu{"storeResoOccu", false, "Flag to store Occupancy in resolution ThnSparse"}; Configurable occEstimator{"occEstimator", 0, "Occupancy estimation (0: None, 1: ITS, 2: FT0C)"}; Configurable saveEpResoHisto{"saveEpResoHisto", false, "Flag to save event plane resolution histogram"}; Configurable ccdbUrl{"ccdbUrl", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; @@ -85,6 +86,9 @@ struct HfTaskFlowCharmHadrons { ConfigurableAxis thnConfigAxisNoCollInTimeRangeNarrow{"thnConfigAxisNoCollInTimeRangeNarrow", {2, 0, 2}, ""}; ConfigurableAxis thnConfigAxisNoCollInTimeRangeStandard{"thnConfigAxisNoCollInTimeRangeStandard", {2, 0, 2}, ""}; ConfigurableAxis thnConfigAxisNoCollInRofStandard{"thnConfigAxisNoCollInRofStandard", {2, 0, 2}, ""}; + ConfigurableAxis thnConfigAxisResoFT0cFV0a{"thnConfigAxisResoFT0cFV0a", {160, -8, -8}, ""}; + ConfigurableAxis thnConfigAxisResoFT0cTPCtot{"thnConfigAxisResoFT0cTPCtot", {160, -8, -8}, ""}; + ConfigurableAxis thnConfigAxisResoFV0aTPCtot{"thnConfigAxisResoFV0aTPCtot", {160, -8, -8}, ""}; using CandDsDataWMl = soa::Filtered>; using CandDsData = soa::Filtered>; @@ -124,6 +128,9 @@ struct HfTaskFlowCharmHadrons { void init(InitContext&) { + if (storeResoOccu && occEstimator == 0) { + LOGP(fatal, "Occupancy estimation must be enabled to store resolution THnSparse! Please check your configuration!"); + } const AxisSpec thnAxisInvMass{thnConfigAxisInvMass, "Inv. mass (GeV/#it{c}^{2})"}; const AxisSpec thnAxisPt{thnConfigAxisPt, "#it{p}_{T} (GeV/#it{c})"}; const AxisSpec thnAxisCent{thnConfigAxisCent, "Centrality"}; @@ -139,6 +146,10 @@ struct HfTaskFlowCharmHadrons { const AxisSpec thnAxisNoCollInTimeRangeNarrow{thnConfigAxisNoCollInTimeRangeNarrow, "NoCollInTimeRangeNarrow"}; const AxisSpec thnAxisNoCollInTimeRangeStandard{thnConfigAxisNoCollInTimeRangeStandard, "NoCollInTimeRangeStandard"}; const AxisSpec thnAxisNoCollInRofStandard{thnConfigAxisNoCollInRofStandard, "NoCollInRofStandard"}; + // TODO: currently only the Q vector of FT0c FV0a and TPCtot are considered + const AxisSpec thnAxisResoFT0cFV0a{thnConfigAxisResoFT0cFV0a, "Q_{FT0c} #bullet Q_{FV0a}"}; + const AxisSpec thnAxisResoFT0cTPCtot{thnConfigAxisResoFT0cTPCtot, "Q_{FT0c} #bullet Q_{TPCtot}"}; + const AxisSpec thnAxisResoFV0aTPCtot{thnConfigAxisResoFV0aTPCtot, "Q_{FV0a} #bullet Q_{TPCtot}"}; std::vector axes = {thnAxisInvMass, thnAxisPt, thnAxisCent, thnAxisScalarProd}; if (storeEP) { @@ -199,6 +210,18 @@ struct HfTaskFlowCharmHadrons { registry.add("epReso/hEpResoFV0aTPCneg", "hEpResoFV0aTPCneg; centrality; #Delta#Psi_{sub}", {HistType::kTH2F, {thnAxisCent, thnAxisCosNPhi}}); registry.add("epReso/hEpResoFV0aTPCtot", "hEpResoFV0aTPCtot; centrality; #Delta#Psi_{sub}", {HistType::kTH2F, {thnAxisCent, thnAxisCosNPhi}}); registry.add("epReso/hEpResoTPCposTPCneg", "hEpResoTPCposTPCneg; centrality; #Delta#Psi_{sub}", {HistType::kTH2F, {thnAxisCent, thnAxisCosNPhi}}); + + if (storeResoOccu) { + std::vector axes_reso = {thnAxisCent, thnAxisResoFT0cFV0a, thnAxisResoFT0cTPCtot, thnAxisResoFV0aTPCtot}; + if (occEstimator == 1) { + axes_reso.push_back(thnAxisOccupancyITS, thnAxisNoSameBunchPileup, thnAxisOccupancy, + thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard); + } else { + axes_reso.push_back(thnAxisOccupancyFT0C, thnAxisNoSameBunchPileup, thnAxisOccupancy, + thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard); + } + registry.add("hSparseReso", "THn for resolution with occupancy", HistType::kTHnSparseF, axes_reso); + } } hfEvSel.addHistograms(registry); // collision monitoring @@ -260,6 +283,19 @@ struct HfTaskFlowCharmHadrons { return deltaPsi; } + /// Get the event selection flags + /// \param hfevselflag is the event selection flag + std::vector getEventSelectionFlags(uint16_t hfevselflag) + { + return { + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoSameBunchPileup), + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)}; + }; + } + /// Fill THnSparse /// \param mass is the invariant mass of the candidate /// \param pt is the transverse momentum of the candidate @@ -281,37 +317,22 @@ struct HfTaskFlowCharmHadrons { uint16_t& hfevselflag) { if (occEstimator != 0) { + std::vector evtSelFlags = getEventSelectionFlags(hfevselflag); if (storeMl) { if (storeEP) { registry.fill(HIST("hSparseFlowCharm"), mass, pt, cent, sp, cosNPhi, cosDeltaPhi, outputMl[0], outputMl[1], occupancy, - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoSameBunchPileup), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)); + evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); } else { registry.fill(HIST("hSparseFlowCharm"), mass, pt, cent, sp, outputMl[0], outputMl[1], occupancy, - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoSameBunchPileup), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)); + evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); } } else { if (storeEP) { registry.fill(HIST("hSparseFlowCharm"), mass, pt, cent, sp, cosNPhi, cosDeltaPhi, occupancy, - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoSameBunchPileup), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)); + evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); } else { registry.fill(HIST("hSparseFlowCharm"), mass, pt, cent, sp, occupancy, - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoSameBunchPileup), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)); + evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); } } } else { @@ -671,6 +692,16 @@ struct HfTaskFlowCharmHadrons { registry.fill(HIST("epReso/hEpResoFV0aTPCneg"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFV0a, epBNegs))); registry.fill(HIST("epReso/hEpResoFV0aTPCtot"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFV0a, epBTots))); registry.fill(HIST("epReso/hEpResoTPCposTPCneg"), centrality, std::cos(harmonic * getDeltaPsiInRange(epBPoss, epBNegs))); + + if (storeResoOccu) { + float occupancy = 0.; + uint16_t hfevflag{}; + occupancy = getOccupancyColl(collision, occEstimator); + registry.fill(HIST("trackOccVsFT0COcc"), collision.trackOccupancyInTimeRange(), collision.ft0cOccupancyInTimeRange()); + hfevflag = hfEvSel.getHfCollisionRejectionMask(collision, cent, ccdb, registry); + std::vector evtSelFlags = getEventSelectionFlags(hfevselflag); + registry.fill(HIST("hSparseReso"), centrality, xQVecFT0c * xQVecFV0a + yQVecFT0c * yQVecFV0a, xQVexFT0c * xQVecTPCtot + yQVecFT0c * yQVecTPCtot, xQVecFV0a * xQVecTPCtot + yQVecFV0a * yQVecTPCtot, evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); + } } } PROCESS_SWITCH(HfTaskFlowCharmHadrons, processResolution, "Process resolution", false); From 45cde62c4e897718a01900fc64024db2047edb50 Mon Sep 17 00:00:00 2001 From: wuctlby Date: Mon, 17 Feb 2025 14:13:06 +0100 Subject: [PATCH 2/3] Debug --- PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx | 49 +++++++++++++----------- 1 file changed, 26 insertions(+), 23 deletions(-) diff --git a/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx b/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx index 857d32156b5..f5399d04d88 100644 --- a/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx +++ b/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx @@ -86,9 +86,9 @@ struct HfTaskFlowCharmHadrons { ConfigurableAxis thnConfigAxisNoCollInTimeRangeNarrow{"thnConfigAxisNoCollInTimeRangeNarrow", {2, 0, 2}, ""}; ConfigurableAxis thnConfigAxisNoCollInTimeRangeStandard{"thnConfigAxisNoCollInTimeRangeStandard", {2, 0, 2}, ""}; ConfigurableAxis thnConfigAxisNoCollInRofStandard{"thnConfigAxisNoCollInRofStandard", {2, 0, 2}, ""}; - ConfigurableAxis thnConfigAxisResoFT0cFV0a{"thnConfigAxisResoFT0cFV0a", {160, -8, -8}, ""}; - ConfigurableAxis thnConfigAxisResoFT0cTPCtot{"thnConfigAxisResoFT0cTPCtot", {160, -8, -8}, ""}; - ConfigurableAxis thnConfigAxisResoFV0aTPCtot{"thnConfigAxisResoFV0aTPCtot", {160, -8, -8}, ""}; + ConfigurableAxis thnConfigAxisResoFT0cFV0a{"thnConfigAxisResoFT0cFV0a", {160, -8, 8}, ""}; + ConfigurableAxis thnConfigAxisResoFT0cTPCtot{"thnConfigAxisResoFT0cTPCtot", {160, -8, 8}, ""}; + ConfigurableAxis thnConfigAxisResoFV0aTPCtot{"thnConfigAxisResoFV0aTPCtot", {160, -8, 8}, ""}; using CandDsDataWMl = soa::Filtered>; using CandDsData = soa::Filtered>; @@ -210,18 +210,18 @@ struct HfTaskFlowCharmHadrons { registry.add("epReso/hEpResoFV0aTPCneg", "hEpResoFV0aTPCneg; centrality; #Delta#Psi_{sub}", {HistType::kTH2F, {thnAxisCent, thnAxisCosNPhi}}); registry.add("epReso/hEpResoFV0aTPCtot", "hEpResoFV0aTPCtot; centrality; #Delta#Psi_{sub}", {HistType::kTH2F, {thnAxisCent, thnAxisCosNPhi}}); registry.add("epReso/hEpResoTPCposTPCneg", "hEpResoTPCposTPCneg; centrality; #Delta#Psi_{sub}", {HistType::kTH2F, {thnAxisCent, thnAxisCosNPhi}}); + } - if (storeResoOccu) { - std::vector axes_reso = {thnAxisCent, thnAxisResoFT0cFV0a, thnAxisResoFT0cTPCtot, thnAxisResoFV0aTPCtot}; - if (occEstimator == 1) { - axes_reso.push_back(thnAxisOccupancyITS, thnAxisNoSameBunchPileup, thnAxisOccupancy, - thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard); - } else { - axes_reso.push_back(thnAxisOccupancyFT0C, thnAxisNoSameBunchPileup, thnAxisOccupancy, - thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard); - } - registry.add("hSparseReso", "THn for resolution with occupancy", HistType::kTHnSparseF, axes_reso); + if (storeResoOccu) { + std::vector axes_reso = {thnAxisCent, thnAxisResoFT0cFV0a, thnAxisResoFT0cTPCtot, thnAxisResoFV0aTPCtot}; + if (occEstimator == 1) { + axes_reso.insert(axes_reso.end(), {thnAxisOccupancyITS, thnAxisNoSameBunchPileup, thnAxisOccupancy, + thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard}); + } else { + axes_reso.insert(axes_reso.end(), {thnAxisOccupancyFT0C, thnAxisNoSameBunchPileup, thnAxisOccupancy, + thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard}); } + registry.add("spReso/hSparseReso", "THn for resolution with occupancy", HistType::kTHnSparseF, axes_reso); } hfEvSel.addHistograms(registry); // collision monitoring @@ -292,7 +292,7 @@ struct HfTaskFlowCharmHadrons { TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)}; + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard) }; } @@ -692,16 +692,19 @@ struct HfTaskFlowCharmHadrons { registry.fill(HIST("epReso/hEpResoFV0aTPCneg"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFV0a, epBNegs))); registry.fill(HIST("epReso/hEpResoFV0aTPCtot"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFV0a, epBTots))); registry.fill(HIST("epReso/hEpResoTPCposTPCneg"), centrality, std::cos(harmonic * getDeltaPsiInRange(epBPoss, epBNegs))); + } - if (storeResoOccu) { - float occupancy = 0.; - uint16_t hfevflag{}; - occupancy = getOccupancyColl(collision, occEstimator); - registry.fill(HIST("trackOccVsFT0COcc"), collision.trackOccupancyInTimeRange(), collision.ft0cOccupancyInTimeRange()); - hfevflag = hfEvSel.getHfCollisionRejectionMask(collision, cent, ccdb, registry); - std::vector evtSelFlags = getEventSelectionFlags(hfevselflag); - registry.fill(HIST("hSparseReso"), centrality, xQVecFT0c * xQVecFV0a + yQVecFT0c * yQVecFV0a, xQVexFT0c * xQVecTPCtot + yQVecFT0c * yQVecTPCtot, xQVecFV0a * xQVecTPCtot + yQVecFV0a * yQVecTPCtot, evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); - } + if (storeResoOccu) { + float occupancy = 0.; + uint16_t hfevflag{}; + occupancy = getOccupancyColl(collision, occEstimator); + registry.fill(HIST("trackOccVsFT0COcc"), collision.trackOccupancyInTimeRange(), collision.ft0cOccupancyInTimeRange()); + hfevflag = hfEvSel.getHfCollisionRejectionMask(collision, centrality, ccdb, registry); + std::vector evtSelFlags = getEventSelectionFlags(hfevflag); + registry.fill(HIST("spReso/hSparseReso"), centrality, xQVecFT0c * xQVecFV0a + yQVecFT0c * yQVecFV0a, + xQVecFT0c * xQVecBTot + yQVecFT0c * yQVecBTot, + xQVecFV0a * xQVecBTot + yQVecFV0a * yQVecBTot, + occupancy, evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); } } PROCESS_SWITCH(HfTaskFlowCharmHadrons, processResolution, "Process resolution", false); From 7d105d66bb36872fc24484b9ecb89095d387ff32 Mon Sep 17 00:00:00 2001 From: ALICE Action Bot Date: Mon, 17 Feb 2025 16:26:22 +0000 Subject: [PATCH 3/3] Please consider the following formatting changes --- PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx | 15 +++++++-------- 1 file changed, 7 insertions(+), 8 deletions(-) diff --git a/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx b/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx index f5399d04d88..ae06cf65024 100644 --- a/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx +++ b/PWGHF/D2H/Tasks/taskFlowCharmHadrons.cxx @@ -216,10 +216,10 @@ struct HfTaskFlowCharmHadrons { std::vector axes_reso = {thnAxisCent, thnAxisResoFT0cFV0a, thnAxisResoFT0cTPCtot, thnAxisResoFV0aTPCtot}; if (occEstimator == 1) { axes_reso.insert(axes_reso.end(), {thnAxisOccupancyITS, thnAxisNoSameBunchPileup, thnAxisOccupancy, - thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard}); - } else { + thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard}); + } else { axes_reso.insert(axes_reso.end(), {thnAxisOccupancyFT0C, thnAxisNoSameBunchPileup, thnAxisOccupancy, - thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard}); + thnAxisNoCollInTimeRangeNarrow, thnAxisNoCollInTimeRangeStandard, thnAxisNoCollInRofStandard}); } registry.add("spReso/hSparseReso", "THn for resolution with occupancy", HistType::kTHnSparseF, axes_reso); } @@ -292,8 +292,7 @@ struct HfTaskFlowCharmHadrons { TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::Occupancy), TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeNarrow), TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInTimeRangeStandard), - TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard) - }; + TESTBIT(hfevselflag, o2::hf_evsel::EventRejection::NoCollInRofStandard)}; } /// Fill THnSparse @@ -702,9 +701,9 @@ struct HfTaskFlowCharmHadrons { hfevflag = hfEvSel.getHfCollisionRejectionMask(collision, centrality, ccdb, registry); std::vector evtSelFlags = getEventSelectionFlags(hfevflag); registry.fill(HIST("spReso/hSparseReso"), centrality, xQVecFT0c * xQVecFV0a + yQVecFT0c * yQVecFV0a, - xQVecFT0c * xQVecBTot + yQVecFT0c * yQVecBTot, - xQVecFV0a * xQVecBTot + yQVecFV0a * yQVecBTot, - occupancy, evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); + xQVecFT0c * xQVecBTot + yQVecFT0c * yQVecBTot, + xQVecFV0a * xQVecBTot + yQVecFV0a * yQVecBTot, + occupancy, evtSelFlags[0], evtSelFlags[1], evtSelFlags[2], evtSelFlags[3], evtSelFlags[4]); } } PROCESS_SWITCH(HfTaskFlowCharmHadrons, processResolution, "Process resolution", false);