From f4fd6718453077e8f25c5279bacb2f4eeef81ba7 Mon Sep 17 00:00:00 2001 From: ptak Date: Wed, 5 Aug 2026 17:10:51 +0200 Subject: [PATCH] [PWGCF] Femto: Store SH Ylm histograms with centrality and kT axes --- PWGCF/Femto/Core/pairHistManager.h | 133 +++++++++++++++-------------- 1 file changed, 69 insertions(+), 64 deletions(-) diff --git a/PWGCF/Femto/Core/pairHistManager.h b/PWGCF/Femto/Core/pairHistManager.h index a917fdb3579..99d66198ad7 100644 --- a/PWGCF/Femto/Core/pairHistManager.h +++ b/PWGCF/Femto/Core/pairHistManager.h @@ -234,7 +234,7 @@ struct ConfPairBinning : o2::framework::ConfigurableGroup { o2::framework::Configurable shUseCent{"shUseCent", false, "SH: bin by centrality instead of multiplicity"}; o2::framework::ConfigurableAxis shCentBins{"shCentBins", {o2::framework::VARIABLE_WIDTH, 0.0f, 200.0f}, "SH: multiplicity/centrality bin edges (like FemtoUniverse confMultKstarBins)"}; o2::framework::ConfigurableAxis shKtBins{"shKtBins", {o2::framework::VARIABLE_WIDTH, 0.1f, 0.2f, 0.3f, 0.4f}, "SH: kT bin edges (like FemtoUniverse confKtKstarBins)"}; - o2::framework::Configurable shPlot1D{"shPlot1D", false, "(SH) Also fill 1D qinv/k* numerator/denominator (h1D) per (mult,kT) bin"}; + o2::framework::Configurable shPlot1D{"shPlot1D", false, "(SH) Also fill the 1D qinv/k* distribution (h1D) and the bin occupancy (BinCount) per (mult,kT) bin"}; }; struct ConfPairCuts : o2::framework::ConfigurableGroup { @@ -586,6 +586,8 @@ class PairHistManager // copy bin edges, stripping the leading VARIABLE_WIDTH (0) marker mShCentEdges.assign(ConfPairBinning.shCentBins.value.begin() + 1, ConfPairBinning.shCentBins.value.end()); mShKtEdges.assign(ConfPairBinning.shKtBins.value.begin() + 1, ConfPairBinning.shKtBins.value.end()); + mShCentSpec = {ConfPairBinning.shCentBins, mShUseCent ? "centrality (%)" : "multiplicity"}; + mShKtSpec = {ConfPairBinning.shKtBins, "k_{T} (GeV/#it{c})"}; } // transverse mass type @@ -1029,80 +1031,80 @@ class PairHistManager const int nKt = static_cast(mShKtEdges.size()) - 1; mShYlmBuffer.assign(nJM, {}); - mShReal.resize(nCent); - mShImag.resize(nCent); + mShReal.resize(nJM); + mShImag.resize(nJM); mShCov.resize(nCent); mSh1D.resize(nCent); mShBinCount.resize(nCent); + + const std::string dir = std::string(prefix) + std::string(AnalysisDir) + "SH/"; + int ihist = 0; + for (int l = 0; l <= mShLMax; ++l) { + for (int m = -l; m <= l; ++m) { + std::string lm = std::to_string(l); + lm += (m < 0) ? std::to_string(l - m) : std::to_string(m); + std::string nameRe = dir; + nameRe += "ReYlm"; + nameRe += lm; + std::string nameIm = dir; + nameIm += "ImYlm"; + nameIm += lm; + // shared "Y_{l}^{m}" suffix for both titles + std::string ylmLabel = "Y_{"; + ylmLabel += std::to_string(l); + ylmLabel += "}^{"; + ylmLabel += std::to_string(m); + ylmLabel += "}"; + std::string titleRe = "Re "; + titleRe += ylmLabel; + titleRe += "; k* (GeV/#it{c}); mult/cent; k_{T} (GeV/#it{c})"; + std::string titleIm = "Im "; + titleIm += ylmLabel; + titleIm += "; k* (GeV/#it{c}); mult/cent; k_{T} (GeV/#it{c})"; + mShReal[ihist] = mHistogramRegistry->add(nameRe.c_str(), titleRe.c_str(), o2::framework::kTH3D, {mShKstarSpec, mShCentSpec, mShKtSpec}); + mShImag[ihist] = mHistogramRegistry->add(nameIm.c_str(), titleIm.c_str(), o2::framework::kTH3D, {mShKstarSpec, mShCentSpec, mShKtSpec}); + mShReal[ihist]->Sumw2(); + mShImag[ihist]->Sumw2(); + ++ihist; + } + } + + const int nAxisLM = 2 * nJM; + const o2::framework::AxisSpec covLmAxis{nAxisLM, -0.5, static_cast(nAxisLM) - 0.5, "l,m #times (re,im)"}; + for (int iCent = 0; iCent < nCent; ++iCent) { - mShReal[iCent].resize(nKt); - mShImag[iCent].resize(nKt); mShCov[iCent].resize(nKt); mSh1D[iCent].resize(nKt); mShBinCount[iCent].resize(nKt); - // folder name: mult_{low}_{high} - const std::string centFolder = "mult_" + std::to_string(static_cast(mShCentEdges[iCent])) + - "_" + std::to_string(static_cast(mShCentEdges[iCent + 1])); + // name suffix: mult_{low}_{high} + std::string centSuffix = "_mult_"; + centSuffix += std::to_string(static_cast(mShCentEdges[iCent])); + centSuffix += "_"; + centSuffix += std::to_string(static_cast(mShCentEdges[iCent + 1])); for (int iKt = 0; iKt < nKt; ++iKt) { - mShReal[iCent][iKt].resize(nJM); - mShImag[iCent][iKt].resize(nJM); - // folder name: kT_{low*100}_{high*100} - std::string ktFolder = "kT_"; - ktFolder += std::to_string(static_cast(mShKtEdges[iKt] * 100.0)); - ktFolder += "_"; - ktFolder += std::to_string(static_cast(mShKtEdges[iKt + 1] * 100.0)); - std::string dir = std::string(prefix) + std::string(AnalysisDir) + "SH/"; - dir += centFolder; - dir += "/"; - dir += ktFolder; - dir += "/"; - - int ihist = 0; - for (int l = 0; l <= mShLMax; ++l) { - for (int m = -l; m <= l; ++m) { - std::string lm = std::to_string(l); - lm += (m < 0) ? std::to_string(l - m) : std::to_string(m); - std::string nameRe = dir; - nameRe += "ReYlm"; - nameRe += lm; - std::string nameIm = dir; - nameIm += "ImYlm"; - nameIm += lm; - // shared "Y_{l}^{m}" suffix for both titles - std::string ylmLabel = "Y_{"; - ylmLabel += std::to_string(l); - ylmLabel += "}^{"; - ylmLabel += std::to_string(m); - ylmLabel += "}"; - std::string titleRe = "Re "; - titleRe += ylmLabel; - titleRe += "; k* (GeV/#it{c}); Re[A_{l}^{m}]"; - std::string titleIm = "Im "; - titleIm += ylmLabel; - titleIm += "; k* (GeV/#it{c}); Im[A_{l}^{m}]"; - mShReal[iCent][iKt][ihist] = mHistogramRegistry->add(nameRe.c_str(), titleRe.c_str(), o2::framework::kTH1D, {mShKstarSpec}); - mShImag[iCent][iKt][ihist] = mHistogramRegistry->add(nameIm.c_str(), titleIm.c_str(), o2::framework::kTH1D, {mShKstarSpec}); - mShReal[iCent][iKt][ihist]->Sumw2(); - mShImag[iCent][iKt][ihist]->Sumw2(); - ++ihist; - } - } + // name suffix: _mult_{low}_{high}_kT_{low*100}_{high*100} + std::string cellSuffix = centSuffix; + cellSuffix += "_kT_"; + cellSuffix += std::to_string(static_cast(mShKtEdges[iKt] * 100.0)); + cellSuffix += "_"; + cellSuffix += std::to_string(static_cast(mShKtEdges[iKt + 1] * 100.0)); // SH covariance TH3D - const int nAxisLM = 2 * nJM; - const o2::framework::AxisSpec covLmAxis{nAxisLM, -0.5, static_cast(nAxisLM) - 0.5, "l,m #times (re,im)"}; std::string nameCov = dir; nameCov += "Cov"; + nameCov += cellSuffix; mShCov[iCent][iKt] = mHistogramRegistry->add(nameCov.c_str(), "SH covariance; k* (GeV/#it{c}); l,m; l,m", o2::framework::kTH3D, {mShKstarSpec, covLmAxis, covLmAxis}); mShCov[iCent][iKt]->Sumw2(); - std::string nameBinCount = dir; - nameBinCount += "BinCount"; - mShBinCount[iCent][iKt] = mHistogramRegistry->add(nameBinCount.c_str(), "SH bin occupancy; k* (GeV/#it{c}); Entries", o2::framework::kTH1D, {mShKstarSpec}); - if (mShPlot1D) { + std::string nameBinCount = dir; + nameBinCount += "BinCount"; + nameBinCount += cellSuffix; + mShBinCount[iCent][iKt] = mHistogramRegistry->add(nameBinCount.c_str(), "SH bin occupancy; k* (GeV/#it{c}); Entries", o2::framework::kTH1D, {mShKstarSpec}); + std::string name1D = dir; name1D += "h1D"; + name1D += cellSuffix; mSh1D[iCent][iKt] = mHistogramRegistry->add(name1D.c_str(), "1D distribution; k* (GeV/#it{c}); Entries", o2::framework::kTH1D, {mShKstarSpec}); mSh1D[iCent][iKt]->Sumw2(); } @@ -1304,13 +1306,14 @@ class PairHistManager mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kQoutQsideQlong, HistTable)), mQout, mQside, mQlong); } if (mPlotSH) { - const int iCent = findShBin(mShUseCent ? mCent : mMult, mShCentEdges); + const float shCentValue = mShUseCent ? mCent : mMult; + const int iCent = findShBin(shCentValue, mShCentEdges); const int iKt = findShBin(mKt, mShKtEdges); if (iCent >= 0 && iKt >= 0) { mYlm.doYlmUpToL(mShLMax, mShOut, mShSide, mShLong, mShYlmBuffer.data()); for (std::size_t i = 0; i < mShYlmBuffer.size(); ++i) { - mShReal[iCent][iKt][i]->Fill(mShKv, std::real(mShYlmBuffer[i])); - mShImag[iCent][iKt][i]->Fill(mShKv, -std::imag(mShYlmBuffer[i])); + mShReal[i]->Fill(mShKv, shCentValue, mKt, std::real(mShYlmBuffer[i])); + mShImag[i]->Fill(mShKv, shCentValue, mKt, -std::imag(mShYlmBuffer[i])); } // covariance: outer product of the (re, -im) Ylm vector packed on 2*nJM axes // (each Ylm contributes two consecutive axis bins: even = real, odd = -imag) @@ -1324,8 +1327,8 @@ class PairHistManager } } - mShBinCount[iCent][iKt]->Fill(mShKv, 1.0); if (mShPlot1D) { + mShBinCount[iCent][iKt]->Fill(mShKv, 1.0); // FemtoUniverse h1D = f3d[0]: qinv (=2k*) for identical-LCMS, else k*. const float sh1DValue = (mShFrame == ShFrameLcmsIdentical) ? (2.0f * mKstar) : mKstar; mSh1D[iCent][iKt]->Fill(sh1DValue); @@ -1693,6 +1696,8 @@ class PairHistManager static constexpr int ShFramePrf = 2; o2::framework::AxisSpec mShKstarSpec{{60, 0.0f, 0.3f}, "k* (GeV/#it{c})"}; // set in init() + o2::framework::AxisSpec mShCentSpec{{1, 0.0f, 200.0f}, "mult/cent"}; + o2::framework::AxisSpec mShKtSpec{{3, 0.1f, 0.4f}, "k_{T} (GeV/#it{c})"}; // kinematics computed in setPair(): axis value + 3 components feeding Ylm float mShKv = 0.f; // kstar (non-identical) or qinv (identical) @@ -1700,9 +1705,9 @@ class PairHistManager float mShSide = 0.f; float mShLong = 0.f; - // SH histograms binned in [iCent][iKt][ihist]; ihist = l*(l+1)+m - std::vector>>> mShReal; - std::vector>>> mShImag; + // SH histograms per [ihist] (ihist = l*(l+1)+m); TH3: k* on X, mult/cent on Y, kT on Z + std::vector> mShReal; + std::vector> mShImag; // SH covariance matrix per [iCent][iKt]; TH3d: k* on X, 2*nJM (l,m x re/im) std::vector>> mShCov; bool mShPlot1D = false;