From a7a427072f3dc8e495ff7a328f2b1570e026dae2 Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Thu, 11 Jun 2026 12:45:31 +0200 Subject: [PATCH 1/9] add 2 and 4 particle correlation and fix errors --- .../Tasks/multiparticleCumulants.cxx | 1366 ++++++++++++----- 1 file changed, 963 insertions(+), 403 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index 0b5ceafdd47..1be1f6e0f67 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -1,6 +1,6 @@ // Copyright 2019-2020 CERN and copyright holders of ALICE O2. -// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. -// All rights not expressly granted are reserved. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright +// holders. All rights not expressly granted are reserved. // // This software is distributed under the terms of the GNU General Public // License v3 (GPL Version 3), copied verbatim in the file "COPYING". @@ -30,6 +30,7 @@ #include #include +#include #include #include #include @@ -51,26 +52,34 @@ using namespace o2; using namespace o2::framework; // Definitions of join tables for Run 3 analysis: -using EventSelection = soa::Join; -using CollisionRec = soa::Join::iterator; // use in json "isMC": "true" for "event-selection-task" -using CollisionRecSim = soa::Join::iterator; +using EventSelection = + soa::Join; // no FV0Ms for cent +using CollisionRec = soa::Join:: + iterator; // use in json "isMC": "true" for "event-selection-task" +using CollisionRecSim = soa::Join::iterator; using CollisionSim = aod::McCollision; -using TracksRec = soa::Join; -using TrackRec = soa::Join::iterator; -using TracksRecSim = soa::Join; // + use in json "isMC" : "true" -using TrackRecSim = soa::Join::iterator; +using TracksRec = soa::Join; +using TrackRec = soa::Join::iterator; +using TracksRecSim = + soa::Join; // + use in json "isMC" : "true" +using TrackRecSim = + soa::Join::iterator; using TracksSim = aod::McParticles; using TrackSim = aod::McParticles::iterator; using namespace std; // *) Define enums: -enum EnRlMc { eRl = 0, - eMc }; +enum EnRlMc { eRl = 0, eMc }; -enum EnRecSim { eRec = 0, - eSim, - eRecAndSim }; +enum EnRecSim { eRec = 0, eSim, eRecAndSim }; enum EnProcess { eProcessRec = 0, // Run 3, only reconstructed @@ -89,100 +98,181 @@ enum EnEventHistograms { eEventHistograms_N }; -const char* eventHistNames[eEventHistograms_N] = { - "Centrality", - "Multiplicity", - "VertexX", - "VertexY", - "VertexZ", - "NumContrib"}; - -enum EnParticleHistograms { - ePt, - ePhi, - eParticleHistograms_N +const char *eventHistNames[eEventHistograms_N] = {"Centrality", "Multiplicity", + "VertexX", "VertexY", + "VertexZ", "NumContrib"}; + +enum EnParticleHistograms { ePt, ePhi, eParticleHistograms_N }; + +const char *particleHistNames[eParticleHistograms_N] = {"Pt", "Phi"}; + +enum EnQAHistograms { eQACent, eQAMultNumContrib, eQAHistograms_N }; + +enum EnCorrHistograms { eCorrCent, eCorrMult, eCorrHistograms_N }; + +const char *corrHistNames[eCorrHistograms_N] = { + "Centrality", + "Multiplicity", +}; + +enum EnCentEstm { eCentFT0C, eCentFT0M, eCentFV0A, eCentEstm_N }; + +const char *centEstmNames[eCentEstm_N] = { + "FT0C", + "FT0M", + "FV0A", }; -const char* particleHistNames[eParticleHistograms_N] = { - "Pt", - "Phi"}; +enum EnMultEstm { eMultFT0C, eMultFT0M, eMultFV0A, eMultEstm_N }; + +const char *multEstmNames[eMultEstm_N] = {"FT0C", "FT0M", "FV0A"}; -enum EnQAHistograms { - eQACent, - eQAMultNumContrib, - eQAHistograms_N +enum EnCutBeforeAfter { eBefore, eAfter, eCutBeforeAfter_N }; + +const char *cutBeforeAfterNames[eCutBeforeAfter_N] = { + "before", + "after", }; // *) Main task: -struct MultiparticleCumulants { // this name is used in lower-case format to name the TDirectoryFile in AnalysisResults.root +struct MultiparticleCumulants { // this name is used in lower-case format to + // name the TDirectoryFile in + // AnalysisResults.root // *) Base TList to hold all output objects: - TString sBaseListName = "Default list name"; // yes, I declare it separately, because I need it also later in BailOut() function + TString sBaseListName = + "Default list name"; // yes, I declare it separately, because I need it + // also later in BailOut() function OutputObj fBaseList{sBaseListName.Data(), OutputObjHandlingPolicy::AnalysisObject, OutputObjSourceType::OutputObjSource}; // *) Service: - Service ccdb; // support for offline callibration data base + Service + ccdb; // support for offline callibration data base Service pdg; // *) Define configurables: - Configurable cfDryRun{"cfDryRun", false, "book all histos and run without filling and calculating anything"}; // example for built-in type (float, string, etc.) - Configurable cfCentEstm{"cfCentEstm", "FT0M", "centrality estimator: FT0M, FV0A, NTPV, FT0C"}; - Configurable cfMultEstm{"cfMultEstm", "FT0A", "multiplicity estimator: FT0A, FT0C, FT0M, FV0A, FV0C, NTracksPV"}; + Configurable cfDryRun{ + "cfDryRun", false, + "book all histos and run without filling and calculating anything"}; + Configurable cfCentEstm{"cfCentEstm", "FT0M", + "centrality estimator: TBD"}; + Configurable cfMultEstm{"cfMultEstm", "FV0A", + "multiplicity estimator: TBD"}; Configurable cfQASwitch{"cfQASwitch", true, "quality assurance switch"}; Configurable cfWeightSwitch{"cfWeightSwitch", true, "weight switch"}; - Configurable cfPrintSwitch{"cfPrintSwitch", true, "printing result switch"}; + Configurable cfPrintSwitch{"cfPrintSwitch", true, + "printing result switch"}; // *) Event cut switches - Configurable cfVertexZCutSwitch{"cfVertexZCutSwitch", true, "vertex z cut switch"}; - Configurable cfSel8CutSwitch{"cfSel8CutSwitch", true, "Sel8 cut switch"}; - Configurable cfCentCutSwitch{"cfCentCutSwitch", true, "centrality cut switch"}; - Configurable cfNumContribCutSwitch{"cfNumContribCutSwitch", true, "NContribution cut switch"}; + Configurable cfVertexZCutSwitch{"cfVertexZCutSwitch", true, + "vertex z cut switch"}; + Configurable cfSel8CutSwitch{"cfSel8CutSwitch", true, + "Sel8 cut switch"}; + Configurable cfCentCutSwitch{"cfCentCutSwitch", true, + "centrality cut switch"}; + Configurable cfNumContribCutSwitch{"cfNumContribCutSwitch", true, + "NContribution cut switch"}; + Configurable cfCentCorrCutSwitch{"cfCentCorrCutSwitch", true, + "centrality correlation cut switch"}; + Configurable cfMultCorrCutSwitch{"cfMultCorrCutSwitch", true, + "multiplicity correlation cut switch"}; // *) Particle cut switches Configurable cfPtCutSwitch{"cfPtCutSwitch", true, "Pt cut switch"}; Configurable cfEtaCutSwitch{"cfEtaCutSwitch", true, "Eta cut switch"}; - Configurable cfSignCutSwitch{"cfSignCutSwitch", true, "Charge cut switch"}; - Configurable cfTpcNClsFoundCutSwitch{"cfTpcNClsFoundCutSwitch", true, ""}; + Configurable cfSignCutSwitch{"cfSignCutSwitch", true, + "Charge cut switch"}; + Configurable cfTpcNClsFoundCutSwitch{"cfTpcNClsFoundCutSwitch", true, + ""}; Configurable cfDCAXYCutSwitch{"cfDCAXYCutSwitch", true, ""}; Configurable cfDCAZCutSwitch{"cfDCAZCutSwitch", true, ""}; // *) Event cut - Configurable> cfVertexZCut{"cfVertexZCut", {-10., 10.}, "vertex z position range: {min, max}[cm]"}; - Configurable> cfCentCut{"cfCentCut", {10., 20.}, "centrality range: {min, max}[%]"}; - Configurable> cfNumContribCut{"cfNumContribCut", {0, 3000.}, "NContribution range: {min, max}"}; + Configurable> cfVertexZCut{ + "cfVertexZCut", {-10., 10.}, "vertex z position range: {min, max}[cm]"}; + Configurable> cfCentCut{ + "cfCentCut", {10., 20.}, "centrality range: {min, max}[%]"}; + Configurable> cfNumContribCut{ + "cfNumContribCut", {0, 3000.}, "NContribution range: {min, max}"}; + Configurable cfCentCorrCut{"cfCentCorrCut", 5, + "centrality difference maximum"}; + Configurable cfMultCorrCut{"cfMultCorrCut", 0.2, + "multiplicity difference maximum"}; // *) Particle cut - Configurable> cfPtCut{"cfPtCut", {0.2, 5.0}, "Pt range: {min, max}[GeV], with convention: min <= Pt < max"}; - Configurable> cfEtaCut{"cfEtaCut", {-0.8, 0.8}, "Eta range: {min, max}, with convention: min <= Eta < max"}; - Configurable> cfSignCut{"cfSignCut", {1, 0, 1}, "sign of charge, 1 to keep and 0 to discard, {negative, neutral, positive}"}; - Configurable> cfTpcNClsFoundCut{"cfTpcNClsFoundCut", {70., 160.}, "range of found TPC clusters for this track geometry: {min, max}"}; - Configurable> cfDCAXYCut{"cfDCAXYCut", {-3.2, 3.2}, "range of distance-of-closest-approach (DCA) of the extrapolated track to the primary position in XY-direction: {min, max}[cm]"}; - Configurable> cfDCAZCut{"cfDCAZCut", {-2.4, 2.4}, "range of distance-of-closest-approach (DCA) of the extrapolated track to the primary position in Z-direction: {min, max}[cm]"}; + Configurable> cfPtCut{ + "cfPtCut", + {0.2, 5.0}, + "Pt range: {min, max}[GeV], with convention: min <= Pt < max"}; + Configurable> cfEtaCut{ + "cfEtaCut", + {-0.8, 0.8}, + "Eta range: {min, max}, with convention: min <= Eta < max"}; + Configurable> cfSignCut{ + "cfSignCut", + {1, 0, 1}, + "sign of charge, 1 to keep and 0 to discard, {negative, neutral, " + "positive}"}; + Configurable> cfTpcNClsFoundCut{ + "cfTpcNClsFoundCut", + {70., 160.}, + "range of found TPC clusters for this track geometry: {min, max}"}; + Configurable> cfDCAXYCut{ + "cfDCAXYCut", + {-3.2, 3.2}, + "range of distance-of-closest-approach (DCA) of the extrapolated track " + "to the primary position in XY-direction: {min, max}[cm]"}; + Configurable> cfDCAZCut{ + "cfDCAZCut", + {-2.4, 2.4}, + "range of distance-of-closest-approach (DCA) of the extrapolated track " + "to the primary position in Z-direction: {min, max}[cm]"}; // *) - Configurable cfFileWithWeights{"cfFileWithWeights", "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root", "path to external ROOT file which holds all particle weights in O2 format"}; - Configurable cfRunNumber{"cfRunNumber", "000123456", "run number"}; + Configurable cfFileWithWeights{ + "cfFileWithWeights", + "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root", + "path to external ROOT file which holds all particle weights in O2 " + "format"}; + Configurable cfRunNumber{"cfRunNumber", "000123456", + "run number"}; // *) Bins - Configurable> cfPtBins{"cfPtBins", {1000, 0., 100.}, "nPtBins, ptMin, ptMax"}; - Configurable> cfPhiBins{"cfPhiBins", {1000, 0., o2::constants::math::TwoPI}, "nPhiBins, phiMin, phiMax"}; - Configurable> cfCentBins{"cfCentBins", {100, 0., 100.}, "nCenBins, cenMin, cenMax"}; - Configurable> cfMultBins{"cfMultBins", {100, 0., 5000.}, "nMultBins, MultMin, MultMax"}; - Configurable> cfVerXBins{"cfVerXBins", {100, -0.05, 0.05}, "nVerXBins, VerXMin, VerXMax"}; - Configurable> cfVerYBins{"cfVerYBins", {100, -0.05, 0.05}, "nVerYBins, VerYMin, VerYMax"}; - Configurable> cfVerZBins{"cfVerZBins", {100, -50., 50.}, "nVerZBins, VerZMin, VerZMax"}; - Configurable> cfNumContribBins{"cfNumContribBins", {100, 0., 5000.}, "nNumContribBins, NumContribMin, NumContribMax"}; - - // *) Define and initialize all data members to be called in the main process* functions: + Configurable> cfPtBins{ + "cfPtBins", {1000, 0., 100.}, "nPtBins, ptMin, ptMax"}; + Configurable> cfPhiBins{ + "cfPhiBins", + {1000, 0., o2::constants::math::TwoPI}, + "nPhiBins, phiMin, phiMax"}; + Configurable> cfCentBins{ + "cfCentBins", {100, 0., 100.}, "nCenBins, cenMin, cenMax"}; + Configurable> cfMultBins{ + "cfMultBins", {100, 0., 5000.}, "nMultBins, MultMin, MultMax"}; + Configurable> cfVerXBins{ + "cfVerXBins", {100, -0.05, 0.05}, "nVerXBins, VerXMin, VerXMax"}; + Configurable> cfVerYBins{ + "cfVerYBins", {100, -0.05, 0.05}, "nVerYBins, VerYMin, VerYMax"}; + Configurable> cfVerZBins{ + "cfVerZBins", {100, -50., 50.}, "nVerZBins, VerZMin, VerZMax"}; + Configurable> cfNumContribBins{ + "cfNumContribBins", + {100, 0., 5000.}, + "nNumContribBins, NumContribMin, NumContribMax"}; + + // *) Define and initialize all data members to be called in the main process* + // functions: // **) Task configuration: struct TaskConfiguration { - bool fProcess[eProcess_N] = {kFALSE}; // Set what to process. See enum EnProcess for full description. Set via implicit variables within a PROCESS_SWITCH clause. - bool fDryRun = kFALSE; // book all histos and run without filling and calculating anything + bool fProcess[eProcess_N] = { + false}; // Set what to process. See enum EnProcess for full description. + // Set via implicit variables within a PROCESS_SWITCH clause. + bool fDryRun = false; // book all histos and run without filling and + // calculating anything std::string fCentEstm = "FT0M"; - std::string fMultEstm = "FT0A"; + std::string fMultEstm = "FV0A"; bool fPrintSwitch = true; @@ -190,6 +280,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam bool fSel8CutSwitch = true; bool fCentCutSwitch = true; bool fNumContribCutSwitch = true; + bool fCentCorrCutSwitch = true; + bool fMultCorrCutSwitch = true; bool fPtCutSwitch = true; bool fEtaCutSwitch = true; @@ -201,6 +293,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam std::vector fVertexZCut = {-10., 10.}; std::vector fCentCut = {10., 20.}; std::vector fNumContribCut = {0, 3000.}; + float fCentCorrCut = 5; + float fMultCorrCut = 0.2; std::vector fPtCut = {0.2, 5.0}; std::vector fEtaCut = {-0.8, 0.8}; @@ -219,32 +313,47 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam std::vector fVerZBins = {-50., 50.}; std::vector fNumContribBins = {0, 5000}; - std::string fFileWithWeights = "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root"; + std::string fFileWithWeights = + "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root"; std::string fRunNumber = "000123456"; - } tc; // you have to prepend "tc." for all objects name in this group later in the code + + } tc; // you have to prepend "tc." for all objects name in this group later in + // the code // **) Particle histograms: struct ParticleHistograms { - TList* fParticleHistogramsList = NULL; //! fWeightHistograms; + TList *fWeightHistogramsList = NULL; + std::vector fWeightHistograms; } wt; struct EventByEventQuantities { @@ -253,33 +362,88 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam float fCentralitySim = 0.; float fImpactParameter = 0.; float fNumContrib = 0.; + + float fTwoParticleCorrelationEbye[eCutBeforeAfter_N][3] = { + {0., 0., 0.}, {0., 0., 0.}}; //<2> [before, after][v2^2, v3^2, v4^2] + float fFourParticleCorrelationEbye[eCutBeforeAfter_N][2] = { + {0., 0.}, {0., 0.}}; //<4> [before, after][v3^2 v2^2, v4^2 v2^2] } ebye; - template - bool ctEventCuts(T1 const& collision) - { + struct MultiparticleCorrelationProfile { + TList *fMultiparticleCorrelationProfilesList = NULL; + TProfile *fTwoParticleCorrelationProfiles[eCutBeforeAfter_N] = {NULL}; + TProfile *fFourParticleCorrelationProfiles[eCutBeforeAfter_N] = {NULL}; + } mc; + + struct MultiparticleCorrelationCalculation { + int h1, h2, h3, h4, h5, h6, h7, h8; + // Book Q-vector components: + static constexpr int MaxCorrelator = 2; //<> + static constexpr int MaxHarmonic = 5; // n+1 in vn + static constexpr int MaxPower = MaxCorrelator + 1; + TComplex fQvectorBefore[MaxHarmonic][MaxPower]; + TComplex fQvectorAfter[MaxHarmonic] + [MaxPower]; // All needed Q-vector components + // Store the final results here: + // Remark: [2][MaxCorrelator] => [Cos,Sin][<2>,<3>,<4>,<5>,<6>,<7>,<8>] + } mcc; + + template + bool ctEventCuts(T1 const &collision, T2 rlCollisionCentAll, + T3 rlCollisionMultAll) { bool pass = true; bool bVertexZCut = true; bool bSel8Cut = true; bool bCentCut = true; bool bNumContribCut = true; + bool bCentCorrCut = true; + bool bMultCorrCut = true; // *) For real event and MC event - bVertexZCut = collision.posZ() < tc.fVertexZCut[1] && collision.posZ() > tc.fVertexZCut[0]; + bVertexZCut = collision.posZ() < tc.fVertexZCut[1] && + collision.posZ() > tc.fVertexZCut[0]; // *) For real event only if constexpr (rm == eRl) { + // *) Sel8Cut bSel8Cut = collision.sel8(); - bCentCut = ebye.fCentrality < tc.fCentCut[1] && ebye.fCentrality > tc.fCentCut[0]; - bNumContribCut = ebye.fNumContrib < tc.fNumContribCut[1] && ebye.fNumContrib > tc.fNumContribCut[0]; + // *) CentCut + bCentCut = ebye.fCentrality < tc.fCentCut[1] && + ebye.fCentrality > tc.fCentCut[0]; + // *) NumContribCut + bNumContribCut = ebye.fNumContrib < tc.fNumContribCut[1] && + ebye.fNumContrib > tc.fNumContribCut[0]; + // *) CentCorrCut + float iCent = 0.; + float jCent = 0.; + for (int i = 0; i < eCentEstm_N; i++) { + iCent = rlCollisionCentAll[i]; + for (int j = i + 1; j < eCentEstm_N; j++) { + jCent = rlCollisionCentAll[j]; + bCentCorrCut = std::abs((iCent - jCent)) < tc.fCentCorrCut; + } + } + // *) MultCorrCut + float iMult = 0.; + float jMult = 0.; + for (int i = 0; i < eMultEstm_N; i++) { + iMult = rlCollisionMultAll[i]; + for (int j = i + 1; j < eMultEstm_N; j++) { + jMult = rlCollisionMultAll[j]; + bMultCorrCut = std::abs((iMult - jMult) / iMult) < tc.fMultCorrCut; + } + } } // *) For MC event only if constexpr (rm == eMc) { - bCentCut = ebye.fCentralitySim < tc.fCentCut[1] && ebye.fCentralitySim > tc.fCentCut[0]; + // *) CentCut + bCentCut = ebye.fCentralitySim < tc.fCentCut[1] && + ebye.fCentralitySim > tc.fCentCut[0]; } + // *) Combine all switches if (tc.fVertexZCutSwitch) { pass = pass && bVertexZCut; } @@ -292,13 +456,17 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (tc.fNumContribCutSwitch) { pass = pass && bNumContribCut; } + if (tc.fCentCorrCutSwitch) { + pass = pass && bCentCorrCut; + } + if (tc.fMultCorrCutSwitch) { + pass = pass && bMultCorrCut; + } return pass; } - template - bool ctParticleCuts(T1 const& track) - { + template bool ctParticleCuts(T1 const &track) { bool pass = true; bool bPtCut = true; @@ -314,24 +482,28 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // *) For real event only if constexpr (rm == eRl) { - bSignCut = (track.sign() == -1 && tc.fSignCut[0]) || (track.sign() == 0 && tc.fSignCut[1]) || (track.sign() == 1 && tc.fSignCut[2]); - bTpcNClsFoundCut = track.tpcNClsFound() < tc.fTpcNClsFoundCut[1] && track.tpcNClsFound() > tc.fTpcNClsFoundCut[0]; - bDCAXYCut = track.dcaXY() < tc.fDCAXYCut[1] && track.dcaXY() > tc.fDCAXYCut[0]; + bSignCut = (track.sign() == -1 && tc.fSignCut[0]) || + (track.sign() == 0 && tc.fSignCut[1]) || + (track.sign() == 1 && tc.fSignCut[2]); + bTpcNClsFoundCut = track.tpcNClsFound() < tc.fTpcNClsFoundCut[1] && + track.tpcNClsFound() > tc.fTpcNClsFoundCut[0]; + bDCAXYCut = + track.dcaXY() < tc.fDCAXYCut[1] && track.dcaXY() > tc.fDCAXYCut[0]; bDCAZCut = track.dcaZ() < tc.fDCAZCut[1] && track.dcaZ() > tc.fDCAZCut[0]; } // *) For mc event only if constexpr (rm == eMc) { - // TDatabasePDG *db= TDatabasePDG::Instance(); - TParticlePDG* particle = pdg->GetParticle(track.pdgCode()); - + TParticlePDG *particle = pdg->GetParticle(track.pdgCode()); if (!particle) { // LOGF(warning, "PDG code %d not found", track.pdgCode()); bSignCut = false; } else { // LOGF(info, "PDG code %d found", track.pdgCode()); float charge = particle->Charge(); - bSignCut = (charge < 0 && tc.fSignCut[0]) || (charge == 0 && tc.fSignCut[1]) || (charge > 0 && tc.fSignCut[2]); + bSignCut = (charge < 0 && tc.fSignCut[0]) || + (charge == 0 && tc.fSignCut[1]) || + (charge > 0 && tc.fSignCut[2]); } } @@ -357,15 +529,74 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam return pass; } - TObject* getObjectFromList(TList* list, const char* objectName) - { - // Get TObject pointer from TList, even if it's in some nested TList. Foreseen - // to be used to fetch histograms or profiles from files directly. + TComplex mccQ(int n, int p, EnCutBeforeAfter eba) { + // Using the fact that Q{-n,p} = Q{n,p}^*. + + if (eba == eBefore) { + if (n >= 0) { + return mcc.fQvectorBefore[n][p]; + } + return TComplex::Conjugate(mcc.fQvectorBefore[-n][p]); + } else if (eba == eAfter) { + if (n >= 0) { + return mcc.fQvectorAfter[n][p]; + } + return TComplex::Conjugate(mcc.fQvectorAfter[-n][p]); + } + + } // TComplex Q(int n, int p) + + TComplex mccRecursion(int n, int *harmonic, EnCutBeforeAfter eba, + int mult = 1, int skip = 0) { + // Calculate multi-particle correlators by using recursion (an improved + // faster version) originally developed by Kristjan Gulbrandsen + // (gulbrand@nbi.dk). + + int nm1 = n - 1; + TComplex c(Q(harmonic[nm1], mult, eba)); + if (nm1 == 0) + return c; + c *= mccRecursion(nm1, harmonic, eba); + if (nm1 == skip) + return c; + + int multp1 = mult + 1; + int nm2 = n - 2; + int counter1 = 0; + int hhold = harmonic[counter1]; + harmonic[counter1] = harmonic[nm2]; + harmonic[nm2] = hhold + harmonic[nm1]; + TComplex c2(mccRecursion(nm1, harmonic, eba, multp1, nm2)); + int counter2 = n - 3; + while (counter2 >= skip) { + harmonic[nm2] = harmonic[counter1]; + harmonic[counter1] = hhold; + ++counter1; + hhold = harmonic[counter1]; + harmonic[counter1] = harmonic[nm2]; + harmonic[nm2] = hhold + harmonic[nm1]; + c2 += mccRecursion(nm1, harmonic, eba, multp1, counter2); + --counter2; + } + harmonic[nm2] = harmonic[counter1]; + harmonic[counter1] = hhold; + + if (mult == 1) + return c - c2; + return c - double(mult) * c2; + + } // TComplex AliFlowAnalysisWithMultiparticleCorrelations::mccRecursion(int + // n, int* harmonic, int mult = 1, int skip = 0) + + TObject *getObjectFromList(TList *list, const char *objectName) { + // Get TObject pointer from TList, even if it's in some nested TList. + // Foreseen to be used to fetch histograms or profiles from files directly. // Some ideas taken from TCollection::ls() - // If you have added histograms directly to files (without TList's), then you can fetch them directly with - // file->Get("hist-name"). + // If you have added histograms directly to files (without TList's), then + // you can fetch them directly with file->Get("hist-name"). - // Usage: TH1D *hist = (TH1D*) getObjectFromList("some-valid-TList-pointer","some-object-name"); + // Usage: TH1D *hist = (TH1D*) + // getObjectFromList("some-valid-TList-pointer","some-object-name"); // Insanity checks: if (!list) { @@ -379,18 +610,21 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // The object is in the current base list: - TObject* objectFinal = list->FindObject(objectName); // the final object I am after + TObject *objectFinal = + list->FindObject(objectName); // the final object I am after if (objectFinal) { return objectFinal; } // Otherwise, search for the object recursively in the nested lists: - TObject* objectIter; // iterator object in the loop below + TObject *objectIter; // iterator object in the loop below TIter next(list); - while ((objectIter = next())) // double round braces are to silence the warnings + while ((objectIter = + next())) // double round braces are to silence the warnings { if (TString(objectIter->ClassName()).EqualTo("TList")) { - objectFinal = getObjectFromList(reinterpret_cast(objectIter), objectName); + objectFinal = getObjectFromList(reinterpret_cast(objectIter), + objectName); if (objectFinal) return objectFinal; } @@ -400,17 +634,23 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // TObject* getObjectFromList(TList *list, char *objectName) - std::vector getHistogramsWithWeights(const char* filePath, const char* runNumber) - { + std::vector getHistogramsWithWeights(const char *filePath, + const char *runNumber) { // a) Return value: - std::vector histograms; - TList* baseList = NULL; // base top-level list in the TFile, e.g. named "ccdb_object" - TList* listWithRuns = NULL; // nested list with run-wise TList's holding run-specific weights - - // c) Determine from filePath if the file in on a local machine, or in home dir AliEn, or in CCDB: - // Algorithm: If filePath begins with "/alice/data/CCDB/" then it's in home dir AliEn. - // If filePath begins with "/alice-ccdb.cern.ch/" then it's in CCDB. Therefore, files in AliEn and CCDB must be specified with abs path, - // for local files both abs and relative paths are just fine. + std::vector histograms; + TList *baseList = + NULL; // base top-level list in the TFile, e.g. named "ccdb_object" + TList *listWithRuns = + NULL; // nested list with run-wise TList's holding run-specific weights + + // c) Determine from filePath if the file in on a local machine, or in home + // dir AliEn, or in CCDB: + // Algorithm: If filePath begins with "/alice/data/CCDB/" then it's in + // home dir AliEn. + // If filePath begins with "/alice-ccdb.cern.ch/" then + // it's in CCDB. Therefore, files in AliEn and CCDB must + // be specified with abs path, for local files both abs + // and relative paths are just fine. bool bFileIsInAliEn = false; bool bFileIsInCCDB = false; @@ -424,11 +664,17 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (bFileIsInAliEn) { // File you want to access is in your home dir in AliEn: - TGrid* alien = TGrid::Connect("alien", gSystem->Getenv("USER"), "", ""); // do not forget to add #include to the preamble of your analysis task + TGrid *alien = + TGrid::Connect("alien", gSystem->Getenv("USER"), "", + ""); // do not forget to add #include to the + // preamble of your analysis task if (!alien) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - TFile* weightsFile = TFile::Open(Form("alien://%s", filePath), "READ"); // yes, ROOT can open a file transparently, even if it's sitting in AliEn, with this specific syntax + TFile *weightsFile = TFile::Open( + Form("alien://%s", filePath), + "READ"); // yes, ROOT can open a file transparently, even if it's + // sitting in AliEn, with this specific syntax if (!weightsFile) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -438,49 +684,67 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - // Finally, from the top-level TList, get the desired nested TList => the technical problem here is that it can be nested at any level, - // for thare there is a helper utility function getObjectFromList(...) , see its implementation further below - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); + // Finally, from the top-level TList, get the desired nested TList => the + // technical problem here is that it can be nested at any level, for thare + // there is a helper utility function getObjectFromList(...) , see its + // implementation further below + listWithRuns = + reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; - runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + runNumberWithLeadingZeroes += + runNumber; // another try, with "000" prepended to run number + listWithRuns = reinterpret_cast( + getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } } } else if (bFileIsInCCDB) { // File you want to access is in your home dir in CCDB: - // Remember that here I do not access the file; instead, I directly access the object in that file. - // My home dir in CCDB: https://alice-ccdb.cern.ch/browse/Users/a/abilandz/ => adapt for your case - ccdb->setURL("https://alice-ccdb.cern.ch"); // to be able to use "ccdb" this object in your analysis task, see 4b/ below - baseList = reinterpret_cast(ccdb->get(TString(filePath).ReplaceAll("/alice-ccdb.cern.ch/", "").Data())); + // Remember that here I do not access the file; instead, I directly access + // the object in that file. My home dir in CCDB: + // https://alice-ccdb.cern.ch/browse/Users/a/abilandz/ => adapt for + // your case + ccdb->setURL( + "https://alice-ccdb.cern.ch"); // to be able to use "ccdb" this object + // in your analysis task, see 4b/ below + baseList = reinterpret_cast(ccdb->get( + TString(filePath).ReplaceAll("/alice-ccdb.cern.ch/", "").Data())); baseList->ls(); if (!baseList) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); + listWithRuns = + reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; - runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + runNumberWithLeadingZeroes += + runNumber; // another try, with "000" prepended to run number + listWithRuns = reinterpret_cast( + getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } } - // OK, we got the desired TList with efficiency corrections, after that we can use the common code for all 3 cases (local, AliEn, CCDB, that common code is below) + // OK, we got the desired TList with efficiency corrections, after that we + // can use the common code for all 3 cases (local, AliEn, CCDB, that + // common code is below) } else { // this is the local case, please handle this one now: // Check if the external ROOT file exists at specified path: if (gSystem->AccessPathName(filePath, kFileExists)) { - LOGF(info, "\033[1;33m if(gSystem->AccessPathName(filePath,kFileExists)), filePath = %s \033[0m", filePath); + LOGF(info, + "\033[1;33m if(gSystem->AccessPathName(filePath,kFileExists)), " + "filePath = %s \033[0m", + filePath); LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - TFile* weightsFile = TFile::Open(filePath, "READ"); + TFile *weightsFile = TFile::Open(filePath, "READ"); if (!weightsFile) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -492,20 +756,27 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); // baseList->FindObject(runNumber) + listWithRuns = reinterpret_cast(getObjectFromList( + baseList, runNumber)); // baseList->FindObject(runNumber) if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; - runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + runNumberWithLeadingZeroes += + runNumber; // another try, with "000" prepended to run number + listWithRuns = reinterpret_cast( + getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { baseList->ls(); - LOGF(fatal, "\033[1;31m%s at line %d : this crash can happen if in the output file there is no list with weights for the current run number = %s\033[0m", __FUNCTION__, __LINE__, runNumber); + LOGF(fatal, + "\033[1;31m%s at line %d : this crash can happen if in the " + "output file there is no list with weights for the current run " + "number = %s\033[0m", + __FUNCTION__, __LINE__, runNumber); } } } TIter next(listWithRuns); - TObject* object = nullptr; + TObject *object = nullptr; while (true) { @@ -514,13 +785,13 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam break; } - auto* hist = dynamic_cast(object); + auto *hist = dynamic_cast(object); if (!hist) { continue; } hist->SetDirectory(0); - auto* histClone = dynamic_cast(hist->Clone()); + auto *histClone = dynamic_cast(hist->Clone()); if (!histClone) { LOGF(fatal, "Failed to clone histogram %s", hist->GetName()); } @@ -533,42 +804,58 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // *) Define all member functions to be called in the main process* functions: template - void Steer(T1 const& collision, T2 const& tracks) - { + void Steer(T1 const &collision, T2 const &tracks) { // Dry run: if (tc.fDryRun) { return; } + TH1F *histAcceptanceWeight = wt.fWeightHistograms[1]; + + for (int h = 0; h < mcc.MaxHarmonic; h++) { + for (int p = 0; p < mcc.MaxPower; p++) { + mcc.fQvectorBefore[h][p] = TComplex(0., 0.); + mcc.fQvectorAfter[h][p] = TComplex(0., 0.); + } // for(int p=0;p(collision.multFT0C()), + static_cast(collision.multFT0M()), + static_cast(collision.multFV0A())}; float rlCollisionMult = 0.; - if (tc.fMultEstm == "FT0A") - rlCollisionMult = collision.multFT0A(); - else if (tc.fMultEstm == "FT0C") - rlCollisionMult = collision.multFT0C(); - else if (tc.fMultEstm == "FT0M") - rlCollisionMult = collision.multFT0M(); - else if (tc.fMultEstm == "FV0A") - rlCollisionMult = collision.multFV0A(); - else if (tc.fMultEstm == "FV0C") - rlCollisionMult = collision.multFV0C(); - else if (tc.fMultEstm == "NTracksPV") - rlCollisionMult = collision.multNTracksPV(); + + for (int i = 0; i < eMultEstm_N; i++) { + if (tc.fPrintSwitch) { + LOGF(info, "%s Multiplicity: %f", multEstmNames[i], + rlCollisionMultAll[i]); + } + if (tc.fMultEstm == multEstmNames[i]) { + rlCollisionMult = rlCollisionMultAll[i]; + } + } float rlCollisionNumContrib = 0.; rlCollisionNumContrib = static_cast(collision.numContrib()); @@ -586,31 +873,53 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam LOGF(info, "Vertex Z position: %f", collision.posZ()); // Print NContributors - LOGF(info, "NContributors: %f", static_cast(rlCollisionNumContrib)); + LOGF(info, "NContributors: %f", + static_cast(rlCollisionNumContrib)); } ebye.fCentrality = rlCollisionCent; ebye.fReferenceMultiplicity = rlCollisionMult; ebye.fNumContrib = rlCollisionNumContrib; if constexpr (rs == eRec || rs == eRecAndSim) { - ev.fEventHistograms[eCent][eRec][0]->Fill(rlCollisionCent); - ev.fEventHistograms[eMult][eRec][0]->Fill(rlCollisionMult); - ev.fEventHistograms[eVertexX][eRec][0]->Fill(collision.posX()); - ev.fEventHistograms[eVertexY][eRec][0]->Fill(collision.posY()); - ev.fEventHistograms[eVertexZ][eRec][0]->Fill(collision.posZ()); - ev.fEventHistograms[eNumContrib][eRec][0]->Fill(rlCollisionNumContrib); - - if (ctEventCuts(collision)) { - ev.fEventHistograms[eCent][eRec][1]->Fill(rlCollisionCent); - ev.fEventHistograms[eMult][eRec][1]->Fill(rlCollisionMult); - ev.fEventHistograms[eVertexX][eRec][1]->Fill(collision.posX()); - ev.fEventHistograms[eVertexY][eRec][1]->Fill(collision.posY()); - ev.fEventHistograms[eVertexZ][eRec][1]->Fill(collision.posZ()); - ev.fEventHistograms[eNumContrib][eRec][1]->Fill(rlCollisionNumContrib); + ev.fEventHistograms[eCent][eRec][eBefore]->Fill(rlCollisionCent); + ev.fEventHistograms[eMult][eRec][eBefore]->Fill(rlCollisionMult); + ev.fEventHistograms[eVertexX][eRec][eBefore]->Fill(collision.posX()); + ev.fEventHistograms[eVertexY][eRec][eBefore]->Fill(collision.posY()); + ev.fEventHistograms[eVertexZ][eRec][eBefore]->Fill(collision.posZ()); + ev.fEventHistograms[eNumContrib][eRec][eBefore]->Fill( + rlCollisionNumContrib); + + // if (tc.fCentCorrCutSwitch) { + for (int i = 0; i < eCentEstm_N; i++) { + for (int j = i + 1; j < eCentEstm_N; j++) { + auto *h = cr.fCorrHistograms[eCorrCent][i][j][eBefore]; + if (!h) { + LOGF(fatal, + "Missing histogram " + "cr.fCorrHistograms[eCorrCent][%d][%d][eBefore]", + i, j); + } + h->Fill(rlCollisionCentAll[i], rlCollisionCentAll[j]); + } + } + // } + + // if (tc.fMultCorrCutSwitch) { + for (int i = 0; i < eMultEstm_N; i++) { + for (int j = i + 1; j < eMultEstm_N; j++) { + auto *h = cr.fCorrHistograms[eCorrMult][i][j][eBefore]; + if (!h) { + LOGF(fatal, + "Missing histogram " + "cr.fCorrHistograms[eCorrMult][%d][%d][eBefore]", + i, j); + } + h->Fill(rlCollisionMultAll[i], rlCollisionMultAll[j]); + } } + // } if constexpr (rs == eRecAndSim) { - if (!collision.has_mcCollision()) { if (tc.fPrintSwitch) { LOGF(warning, " No MC collision for this collision, skip..."); @@ -622,8 +931,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam float mcCollisionCent = 0.; - float b = mccollision.impactParameter() * std::pow(10, -15); // convert fm to m - float xs = 7.71 * std::pow(10, -28); // convert barn to m^2 + float b = mccollision.impactParameter() * + std::pow(10, -15); // convert fm to m + float xs = 7.71 * std::pow(10, -28); // convert barn to m^2 mcCollisionCent = o2::constants::math::PI * b * b / xs * 100; ebye.fCentralitySim = mcCollisionCent; @@ -634,31 +944,90 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam LOGF(info, "mc centrality: %f", mcCollisionCent); } - ev.fEventHistograms[eCent][eSim][0]->Fill(mcCollisionCent); - ev.fEventHistograms[eVertexX][eSim][0]->Fill(mccollision.posX()); - ev.fEventHistograms[eVertexY][eSim][0]->Fill(mccollision.posY()); - ev.fEventHistograms[eVertexZ][eSim][0]->Fill(mccollision.posZ()); - - if (ctEventCuts(mccollision)) { - ev.fEventHistograms[eCent][eSim][1]->Fill(mcCollisionCent); - ev.fEventHistograms[eVertexX][eSim][1]->Fill(mccollision.posX()); - ev.fEventHistograms[eVertexY][eSim][1]->Fill(mccollision.posY()); - ev.fEventHistograms[eVertexZ][eSim][1]->Fill(mccollision.posZ()); + ev.fEventHistograms[eCent][eSim][eBefore]->Fill(mcCollisionCent); + ev.fEventHistograms[eVertexX][eSim][eBefore]->Fill(mccollision.posX()); + ev.fEventHistograms[eVertexY][eSim][eBefore]->Fill(mccollision.posY()); + ev.fEventHistograms[eVertexZ][eSim][eBefore]->Fill(mccollision.posZ()); + + if (ctEventCuts(mccollision, rlCollisionCentAll, + rlCollisionMultAll)) { + ev.fEventHistograms[eCent][eSim][eAfter]->Fill(mcCollisionCent); + ev.fEventHistograms[eVertexX][eSim][eAfter]->Fill(mccollision.posX()); + ev.fEventHistograms[eVertexY][eSim][eAfter]->Fill(mccollision.posY()); + ev.fEventHistograms[eVertexZ][eSim][eAfter]->Fill(mccollision.posZ()); } if (qa.fQASwitch) { - qa.fQAHistograms[eQACent][0]->Fill(rlCollisionCent, mcCollisionCent); - qa.fQAHistograms[eQAMultNumContrib][0]->Fill(rlCollisionMult, rlCollisionNumContrib); - if (ctEventCuts(collision) && ctEventCuts(mccollision)) { - qa.fQAHistograms[eQACent][1]->Fill(rlCollisionCent, mcCollisionCent); - qa.fQAHistograms[eQAMultNumContrib][1]->Fill(rlCollisionMult, rlCollisionNumContrib); + qa.fQAHistograms[eQACent][eBefore]->Fill(rlCollisionCent, + mcCollisionCent); + qa.fQAHistograms[eQAMultNumContrib][eBefore]->Fill( + rlCollisionMult, rlCollisionNumContrib); + if (ctEventCuts(collision, rlCollisionCentAll, + rlCollisionMultAll) && + ctEventCuts(mccollision, rlCollisionCentAll, + rlCollisionMultAll)) { + qa.fQAHistograms[eQACent][eAfter]->Fill(rlCollisionCent, + mcCollisionCent); + qa.fQAHistograms[eQAMultNumContrib][eAfter]->Fill( + rlCollisionMult, rlCollisionNumContrib); } } } + + if (ctEventCuts(collision, rlCollisionCentAll, rlCollisionMultAll)) { + ev.fEventHistograms[eCent][eRec][eAfter]->Fill(rlCollisionCent); + ev.fEventHistograms[eMult][eRec][eAfter]->Fill(rlCollisionMult); + ev.fEventHistograms[eVertexX][eRec][eAfter]->Fill(collision.posX()); + ev.fEventHistograms[eVertexY][eRec][eAfter]->Fill(collision.posY()); + ev.fEventHistograms[eVertexZ][eRec][eAfter]->Fill(collision.posZ()); + ev.fEventHistograms[eNumContrib][eRec][eAfter]->Fill( + rlCollisionNumContrib); + + if (tc.fCentCorrCutSwitch) { + for (int i = 0; i < eCentEstm_N; i++) { + for (int j = i + 1; j < eCentEstm_N; j++) { + auto *h = cr.fCorrHistograms[eCorrCent][i][j][eAfter]; + if (!h) { + LOGF(fatal, + "Missing histogram " + "cr.fCorrHistograms[eCorrCent][%d][%d][eAfter]", + i, j); + } + h->Fill(rlCollisionCentAll[i], rlCollisionCentAll[j]); + } + } + } + + if (tc.fMultCorrCutSwitch) { + + for (int i = 0; i < eMultEstm_N; i++) { + for (int j = i + 1; j < eMultEstm_N; j++) { + auto *h = cr.fCorrHistograms[eCorrMult][i][j][eAfter]; + if (!h) { + LOGF(fatal, + "Missing histogram " + "cr.fCorrHistograms[eCorrMult][%d][%d][eAfter]", + i, j); + } + h->Fill(rlCollisionMultAll[i], rlCollisionMultAll[j]); + } + } + } + + } else + return; } + int nTracksBefore = tracks.size(); + int nTracksAfter = 0; + + // Calculate Q-vectors for available angles and weights: + double dPhi = 0.; // particle angle + double wPhi = 1.; // particle weight + double wPhiToPowerP = 1.; // particle weight raised to power p + // Main loop over particles: - for (auto const& track : tracks) { + for (auto const &track : tracks) { // LOGF(info, "Track azimuthal angle: %f", track.phi()); // LOGF(info, "Transverse momentum: %f", track.pt()); @@ -666,18 +1035,51 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if constexpr (rs == eRec || rs == eRecAndSim) { // Fill track pt distribution: - pc.fParticleHistograms[ePt][eRec][0]->Fill(track.pt()); - pc.fParticleHistograms[ePhi][eRec][0]->Fill(track.phi()); + pc.fParticleHistograms[ePt][eRec][eBefore]->Fill(track.pt()); + pc.fParticleHistograms[ePhi][eRec][eBefore]->Fill(track.phi()); + + dPhi = track.phi(); + if (wt.fWeightSwitch) { + auto *hist = wt.fWeightHistograms[1]; + wPhi = hist->GetBinContent( + wt.fWeightHistograms[1]->GetXaxis()->FindBin(dPhi)); + } + + for (int h = 0; h < mcc.MaxHarmonic; h++) { + for (int p = 0; p < mcc.MaxPower; p++) { + if (wt.fWeightSwitch) { + wPhiToPowerP = std::pow(wPhi, p); + } + mcc.fQvectorBefore[h][p] += + TComplex(wPhiToPowerP * o2::constants::math::Cos(h * dPhi), + wPhiToPowerP * o2::constants::math::Sin(h * dPhi)); + } // for(int p=0;p(track)) { - pc.fParticleHistograms[ePt][eRec][1]->Fill(track.pt()); - pc.fParticleHistograms[ePhi][eRec][1]->Fill(track.phi()); + pc.fParticleHistograms[ePt][eRec][eAfter]->Fill(track.pt()); + pc.fParticleHistograms[ePhi][eRec][eAfter]->Fill(track.phi()); + + for (int h = 0; h < mcc.MaxHarmonic; h++) { + for (int p = 0; p < mcc.MaxPower; p++) { + if (wt.fWeightSwitch) { + wPhiToPowerP = std::pow(wPhi, p); + } + mcc.fQvectorAfter[h][p] += + TComplex(wPhiToPowerP * std::cos(h * dPhi), + wPhiToPowerP * std::sin(h * dPhi)); + } // for(int p=0;pFill(mcparticle.pt()); - pc.fParticleHistograms[ePhi][eSim][0]->Fill(mcparticle.phi()); + auto mcparticle = + track.mcParticle(); // corresponding MC truth simulated particle + pc.fParticleHistograms[ePt][eSim][eBefore]->Fill(mcparticle.pt()); + pc.fParticleHistograms[ePhi][eSim][eBefore]->Fill(mcparticle.phi()); if (ctParticleCuts(mcparticle)) { - pc.fParticleHistograms[ePt][eSim][1]->Fill(mcparticle.pt()); - pc.fParticleHistograms[ePhi][eSim][1]->Fill(mcparticle.phi()); + pc.fParticleHistograms[ePt][eSim][eAfter]->Fill(mcparticle.pt()); + pc.fParticleHistograms[ePhi][eSim][eAfter]->Fill(mcparticle.phi()); } } // end of if constexpr (rs == eRecAndSim) { @@ -698,216 +1101,318 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // end of for (int64_t i = 0; i < tracks.size(); i++) { - } // end of template void Steer(T1 const& collision, T2 const& tracks) { + for (int i = 0; i < 3; i++) { + mcc.h1 = -(i + 2); + mcc.h2 = i + 2; + int harmonicsTwoNum[2] = {mcc.h1, mcc.h2}; + int harmonicsTwoDen[2] = {0, 0}; + + TComplex twoRecursionBefore = + mccRecursion(2, harmonicsTwoNum, eBefore) / + mccRecursion(2, harmonicsTwoDen, eBefore).Re(); + double wTwoRecursionBefore = + mccRecursion(2, harmonicsTwoDen, eBefore).Re(); + ebye.fTwoParticleCorrelationEbye[eBefore][i] = twoRecursionBefore.Re(); + + TComplex twoRecursionAfter = + mccRecursion(2, harmonicsTwoNum, eAfter) / + mccRecursion(2, harmonicsTwoDen, eAfter).Re(); + double wTwoRecursionAfter = mccRecursion(2, harmonicsTwoDen, eAfter).Re(); + ebye.fTwoParticleCorrelationEbye[eAfter][i] = twoRecursionAfter.Re(); + + if (nTracksBefore > 1) { + mc.fTwoParticleCorrelationProfiles[eBefore]->Fill( + i + 2.5, ebye.fTwoParticleCorrelationEbye[eBefore][i], + wTwoRecursionBefore); + } + if (nTracksAfter > 1) { + mc.fTwoParticleCorrelationProfiles[eAfter]->Fill( + i + 2.5, ebye.fTwoParticleCorrelationEbye[eAfter][i], + wTwoRecursionAfter); + } + } + + for (int i = 0; i < 2; i++) { + // 4-p correlations: + mcc.h1 = -(i + 3); + mcc.h2 = -2; + mcc.h3 = 2; + mcc.h4 = i + 3; + int harmonicsFourNum[4] = {mcc.h1, mcc.h2, mcc.h3, mcc.h4}; + int harmonicsFourDen[4] = {0, 0, 0, 0}; + + TComplex fourRecursionBefore = + mccRecursion(4, harmonicsFourNum, eBefore) / + mccRecursion(4, harmonicsFourDen, eBefore).Re(); + double wFourRecursionBefore = + mccRecursion(4, harmonicsFourDen, eBefore).Re(); + ebye.fFourParticleCorrelationEbye[eBefore][i] = fourRecursionBefore.Re(); + + TComplex fourRecursionAfter = + mccRecursion(4, harmonicsFourNum, eAfter) / + mccRecursion(4, harmonicsFourDen, eAfter).Re(); + double wFourRecursionAfter = + mccRecursion(4, harmonicsFourDen, eAfter).Re(); + ebye.fFourParticleCorrelationEbye[eAfter][i] = fourRecursionAfter.Re(); + + if (nTracksBefore > 3) { + mc.fFourParticleCorrelationProfiles[eBefore]->Fill( + i + 3.5, ebye.fFourParticleCorrelationEbye[eBefore][i], + wFourRecursionBefore); + } + if (nTracksAfter > 3) { + mc.fFourParticleCorrelationProfiles[eAfter]->Fill( + i + 3.5, ebye.fFourParticleCorrelationEbye[eAfter][i], + wFourRecursionAfter); + } + } + + } // end of template void Steer(T1 + // const& collision, T2 const& tracks) { template - void BookParticleHistograms(T1 const& lPcBins, ParticleHistograms& pc) - { - // *) structure: - // hists without cut - // - name - // - new hists - // └ rec - // └ sim - // hists with cut - // - name - // - new hists - // └ rec - // └ sim - - std::vector lPtBins = lPcBins[histType]; // define local array and initialize it from an array set in the configurables + void BookParticleHistograms(T1 const &lPcBins, ParticleHistograms &pc) { + std::vector lPtBins = + lPcBins[histType]; // define local array and initialize it from an array + // set in the configurables int nBinsPt = static_cast(lPtBins[0]); float minPt = lPtBins[1]; float maxPt = lPtBins[2]; - std::string nameRecNocut = std::string("fHist") + particleHistNames[histType] + std::string("[eRec][before cut]"); - std::string nameSimNocut = std::string("fHist") + particleHistNames[histType] + std::string("[eSim][before cut]"); - std::string nameRecNocutfull = particleHistNames[histType] + std::string(" distribution for reconstructed particles"); - std::string nameSimNocutfull = particleHistNames[histType] + std::string(" distribution for simulated particles"); - - if (doprocessRec || doprocessRecSim) { - pc.fParticleHistograms[histType][eRec][0] = new TH1F(nameRecNocut.c_str(), nameRecNocutfull.c_str(), nBinsPt, minPt, maxPt); - pc.fParticleHistograms[histType][eRec][0]->GetXaxis()->SetTitle(particleHistNames[histType]); - pc.fParticleHistogramsList->Add(pc.fParticleHistograms[histType][eRec][0]); - } - - if (doprocessSim || doprocessRecSim) { - pc.fParticleHistograms[histType][eSim][0] = new TH1F(nameSimNocut.c_str(), nameSimNocutfull.c_str(), nBinsPt, minPt, maxPt); - pc.fParticleHistograms[histType][eSim][0]->GetXaxis()->SetTitle(particleHistNames[histType]); - pc.fParticleHistogramsList->Add(pc.fParticleHistograms[histType][eSim][0]); - } - - std::string nameRecCut = std::string("fHist") + particleHistNames[histType] + std::string("[eRec][after cut]"); - std::string nameSimCut = std::string("fHist") + particleHistNames[histType] + std::string("[eSim][after cut]"); - std::string nameRecCutfull = particleHistNames[histType] + std::string(" distribution for reconstructed particles"); - std::string nameSimCutfull = particleHistNames[histType] + std::string(" distribution for simulated particles"); - - if (doprocessRec || doprocessRecSim) { - pc.fParticleHistograms[histType][eRec][1] = new TH1F(nameRecCut.c_str(), nameRecCutfull.c_str(), nBinsPt, minPt, maxPt); - pc.fParticleHistograms[histType][eRec][1]->GetXaxis()->SetTitle(particleHistNames[histType]); - pc.fParticleHistogramsList->Add(pc.fParticleHistograms[histType][eRec][1]); - } + for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { + + std::string nameRec = + Form("fHist%s[eRec][%s cut]", particleHistNames[histType], + cutBeforeAfterNames[ba]); + std::string nameSim = + Form("fHist%s[eSim][%s cut]", particleHistNames[histType], + cutBeforeAfterNames[ba]); + std::string nameRecfull = + Form("%s distribution for reconstructed particles", + particleHistNames[histType]); + std::string nameSimfull = Form("%s distribution for simulated particles", + particleHistNames[histType]); + + if (doprocessRec || doprocessRecSim) { + pc.fParticleHistograms[histType][eRec][ba] = new TH1F( + nameRec.c_str(), nameRecfull.c_str(), nBinsPt, minPt, maxPt); + pc.fParticleHistograms[histType][eRec][ba]->GetXaxis()->SetTitle( + particleHistNames[histType]); + pc.fParticleHistogramsList->Add( + pc.fParticleHistograms[histType][eRec][ba]); + } - if (doprocessSim || doprocessRecSim) { - pc.fParticleHistograms[histType][eSim][1] = new TH1F(nameSimCut.c_str(), nameSimCutfull.c_str(), nBinsPt, minPt, maxPt); - pc.fParticleHistograms[histType][eSim][1]->GetXaxis()->SetTitle(particleHistNames[histType]); - pc.fParticleHistogramsList->Add(pc.fParticleHistograms[histType][eSim][1]); + if (doprocessSim || doprocessRecSim) { + pc.fParticleHistograms[histType][eSim][ba] = new TH1F( + nameSim.c_str(), nameSimfull.c_str(), nBinsPt, minPt, maxPt); + pc.fParticleHistograms[histType][eSim][ba]->GetXaxis()->SetTitle( + particleHistNames[histType]); + pc.fParticleHistogramsList->Add( + pc.fParticleHistograms[histType][eSim][ba]); + } } } template - void BookEventHistograms(T1 const& lEvBins, EventHistograms& ev) - { - // *) structure: - // hists without cut - // - name - // └ centrality + estimator - // └ multiplicity + estimator - // └ others - // - new hists - // └ rec - // └ sim (without multiplicity and nContrib) - // hists with cut - // - name - // └ centrality + estimator - // └ multiplicity + estimator - // └ others - // - new hists - // └ rec - // └ sim (without multiplicity and nContrib) - - std::vector lCentBins = lEvBins[histType]; // define local array and initialize it from an array set in the configurables + void BookEventHistograms(T1 const &lEvBins, EventHistograms &ev) { + std::vector lCentBins = + lEvBins[histType]; // define local array and initialize it from an array + // set in the configurables int nBinsCent = static_cast(lCentBins[0]); float minCent = lCentBins[1]; float maxCent = lCentBins[2]; - std::string nameRecNocut = std::string("fHist") + eventHistNames[histType] + std::string("[eRec][before cut]"); - std::string nameSimNocut = std::string("fHist") + eventHistNames[histType] + std::string("[eSim][before cut]"); - std::string nameRecNocutfull; - std::string nameSimNocutfull; - - if constexpr (histType == eCent) { - std::string nameRecNocutfull = tc.fCentEstm + eventHistNames[histType] + std::string(" distribution for reconstructed events"); - std::string nameSimNocutfull = tc.fCentEstm + eventHistNames[histType] + std::string(" distribution for simulated events"); - } else if constexpr (histType == eMult) { - std::string nameRecNocutfull = tc.fMultEstm + eventHistNames[histType] + std::string(" distribution for reconstructed events"); - std::string nameSimNocutfull = tc.fMultEstm + eventHistNames[histType] + std::string(" distribution for simulated events"); - } else { - std::string nameRecNocutfull = eventHistNames[histType] + std::string(" distribution for reconstructed events"); - std::string nameSimNocutfull = eventHistNames[histType] + std::string(" distribution for simulated events"); - } + for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { + + std::string nameRec = + Form("fHist%s[eRec][%s cut]", eventHistNames[histType], + cutBeforeAfterNames[ba]); + std::string nameSim = + Form("fHist%s[eSim][%s cut]", eventHistNames[histType], + cutBeforeAfterNames[ba]); + std::string nameRecfull; + std::string nameSimfull; + + if constexpr (histType == eCent) { + std::string nameRecfull = + Form("%s %s distribution for reconstructed events", + tc.fCentEstm.c_str(), eventHistNames[histType]); + std::string nameSimfull = + Form("%s %s distribution for simulated events", + tc.fCentEstm.c_str(), eventHistNames[histType]); + } else if constexpr (histType == eMult) { + std::string nameRecfull = + Form("%s %s distribution for reconstructed events", + tc.fMultEstm.c_str(), eventHistNames[histType]); + std::string nameSimfull = + Form("%s %s distribution for simulated events", + tc.fMultEstm.c_str(), eventHistNames[histType]); + } else { + std::string nameRecfull = + Form("%s distribution for reconstructed events", + eventHistNames[histType]); + std::string nameSimfull = Form("%s distribution for simulated events", + eventHistNames[histType]); + } - if (doprocessRec || doprocessRecSim) { - ev.fEventHistograms[histType][eRec][0] = new TH1F(nameRecNocut.c_str(), nameRecNocutfull.c_str(), nBinsCent, minCent, maxCent); - ev.fEventHistograms[histType][eRec][0]->GetXaxis()->SetTitle(eventHistNames[histType]); - ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eRec][0]); - } + if (doprocessRec || doprocessRecSim) { + ev.fEventHistograms[histType][eRec][ba] = new TH1F( + nameRec.c_str(), nameRecfull.c_str(), nBinsCent, minCent, maxCent); + ev.fEventHistograms[histType][eRec][ba]->GetXaxis()->SetTitle( + eventHistNames[histType]); + ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eRec][ba]); + } - if (doprocessSim || doprocessRecSim) { - if constexpr (histType != eNumContrib && histType != eMult) { - ev.fEventHistograms[histType][eSim][0] = new TH1F(nameSimNocut.c_str(), nameSimNocutfull.c_str(), nBinsCent, minCent, maxCent); - ev.fEventHistograms[histType][eSim][0]->GetXaxis()->SetTitle(eventHistNames[histType]); - ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eSim][0]); - } // No nContrib and multiplicity for processSim + if (doprocessSim || doprocessRecSim) { + if constexpr (histType != eNumContrib && histType != eMult) { + ev.fEventHistograms[histType][eSim][ba] = + new TH1F(nameSim.c_str(), nameSimfull.c_str(), nBinsCent, minCent, + maxCent); + ev.fEventHistograms[histType][eSim][ba]->GetXaxis()->SetTitle( + eventHistNames[histType]); + ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eSim][ba]); + } // No nContrib and multiplicity for processSim + } } + } - std::string nameRecCut = std::string("fHist") + eventHistNames[histType] + std::string("[eRec][after cut]"); - std::string nameSimCut = std::string("fHist") + eventHistNames[histType] + std::string("[eSim][after cut]"); - std::string nameRecCutfull; - std::string nameSimCutfull; - - if constexpr (histType == eCent) { - std::string nameRecCutfull = tc.fCentEstm + eventHistNames[histType] + std::string(" distribution for reconstructed events"); - std::string nameSimCutfull = tc.fCentEstm + eventHistNames[histType] + std::string(" distribution for simulated events"); - } else if constexpr (histType == eMult) { - std::string nameRecCutfull = tc.fMultEstm + eventHistNames[histType] + std::string(" distribution for reconstructed events"); - std::string nameSimCutfull = tc.fMultEstm + eventHistNames[histType] + std::string(" distribution for simulated events"); + template + void BookQAHistograms(T1 const &lQABins, QAHistograms &qa) { + int nBinsCentX = 0; + float minCentX = 0.; + float maxCentX = 0.; + int nBinsCentY = 0; + float minCentY = 0.; + float maxCentY = 0.; + int nBinsCent = 0; + float minCent = 0.; + float maxCent = 0.; + + if (histType == eMult) { + std::vector lCentBinsX = lQABins[1]; // MultBins + nBinsCentX = static_cast(lCentBinsX[0]); + minCentX = lCentBinsX[1]; + maxCentX = lCentBinsX[2]; + std::vector lCentBinsY = lQABins[2]; // nContribBins + nBinsCentY = static_cast(lCentBinsY[0]); + minCentY = lCentBinsY[1]; + maxCentY = lCentBinsY[2]; } else { - std::string nameRecCutfull = eventHistNames[histType] + std::string(" distribution for reconstructed events"); - std::string nameSimCutfull = eventHistNames[histType] + std::string(" distribution for simulated events"); + std::vector lCentBins = lQABins[histType]; + nBinsCent = static_cast(lCentBins[0]); + minCent = lCentBins[1]; + maxCent = lCentBins[2]; } - if (doprocessRec || doprocessRecSim) { - ev.fEventHistograms[histType][eRec][1] = new TH1F(nameRecCut.c_str(), nameRecCutfull.c_str(), nBinsCent, minCent, maxCent); - ev.fEventHistograms[histType][eRec][1]->GetXaxis()->SetTitle(eventHistNames[histType]); - ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eRec][1]); - } + for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { + std::string name = Form("fHist%s[%s cut]", eventHistNames[histType], + cutBeforeAfterNames[ba]); + std::string namefull; + + if constexpr (histType == eCent) { + std::string namefull = + Form("Quality assurance of %s %s", tc.fCentEstm.c_str(), + eventHistNames[histType]); + } else if constexpr (histType == eMult) { + std::string namefull = + Form("Quality assurance of %s %s vs. NContributors", + tc.fMultEstm.c_str(), eventHistNames[histType]); + } else { + std::string namefull = + Form("Quality assurance of %s", eventHistNames[histType]); + } - if (doprocessSim || doprocessRecSim) { - if constexpr (histType != eNumContrib && histType != eMult) { - ev.fEventHistograms[histType][eSim][1] = new TH1F(nameSimCut.c_str(), nameSimCutfull.c_str(), nBinsCent, minCent, maxCent); - ev.fEventHistograms[histType][eSim][1]->GetXaxis()->SetTitle(eventHistNames[histType]); - ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eSim][1]); + if constexpr (histType == eMult) { + qa.fQAHistograms[histType][ba] = + new TH2F(name.c_str(), namefull.c_str(), nBinsCentX, minCentX, + maxCentX, nBinsCentY, minCentY, maxCentY); + qa.fQAHistograms[histType][ba]->GetYaxis()->SetTitle("NContributors"); + qa.fQAHistograms[histType][ba]->GetXaxis()->SetTitle( + "Reference multiplicity"); + } else { + qa.fQAHistograms[histType][ba] = + new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, + maxCent, nBinsCent, minCent, maxCent); + qa.fQAHistograms[histType][ba]->GetYaxis()->SetTitle( + Form("Simulated %s", eventHistNames[histType])); + qa.fQAHistograms[histType][ba]->GetXaxis()->SetTitle( + Form("Reconstructed %s", eventHistNames[histType])); } + qa.fQAHistogramsList->Add(qa.fQAHistograms[histType][ba]); } } - template - void BookQAHistograms(T1 const& lQABins, QAHistograms& qa) - { - // *) structure: - // hists without cut - // └ simulated centrality vs. reconstructed centrality - // └ nContrib vs. multiplicity - // └ others - // hists with cut - // └ simulated centrality vs. reconstructed centrality - // └ nContrib vs. multiplicity - // └ others - - std::vector lCentBins = lQABins[histType]; + template + void BookCorrHistograms(T1 const &lCrBins, CorrHistograms &cr) { + + std::vector lCentBins = + lCrBins[histType]; // define local array and initialize it from an array + // set in the configurables int nBinsCent = static_cast(lCentBins[0]); float minCent = lCentBins[1]; float maxCent = lCentBins[2]; - std::string nameNocut = std::string("fHist") + eventHistNames[histType] + std::string("[before cut]"); - std::string nameNocutfull; - if constexpr (histType == eCent) { - std::string nameNocutfull = std::string("Quality assurance of ") + tc.fCentEstm + eventHistNames[histType]; - } else if constexpr (histType == eNumContrib) { - std::string nameNocutfull = std::string("Quality assurance of ") + tc.fMultEstm + eventHistNames[histType] + std::string(" vs. NContributors"); - } else { - std::string nameNocutfull = std::string("Quality assurance of ") + eventHistNames[histType]; - } + for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { + std::string name = ""; + std::string namefull = + Form("%s correlation %s cut", corrHistNames[histType], + cutBeforeAfterNames[ba]); + + int N = 0; + if constexpr (histType == eCorrCent) { + N = eCentEstm_N; + } else if constexpr (histType == eCorrMult) { + N = eMultEstm_N; + } + for (int i = 0; i < N; i++) { + + std::string titleX = ""; + if (histType == eCorrCent) { + titleX = Form("%s %s", centEstmNames[i], corrHistNames[histType]); + } else if (histType == eCorrMult) { + titleX = Form("%s %s", multEstmNames[i], corrHistNames[histType]); + // x axis bin... + } - qa.fQAHistograms[histType][0] = new TH2F(nameNocut.c_str(), nameNocutfull.c_str(), nBinsCent, minCent, maxCent, nBinsCent, minCent, maxCent); - if constexpr (histType == eNumContrib) { - qa.fQAHistograms[histType][0]->GetYaxis()->SetTitle("NContributors"); - qa.fQAHistograms[histType][0]->GetXaxis()->SetTitle("Reference multiplicity"); - } else { - qa.fQAHistograms[histType][0]->GetYaxis()->SetTitle(Form("Simulated %s", eventHistNames[histType])); - qa.fQAHistograms[histType][0]->GetXaxis()->SetTitle(Form("Reconstructed %s", eventHistNames[histType])); - } - qa.fQAHistogramsList->Add(qa.fQAHistograms[histType][0]); - - std::string nameCut = std::string("fHist") + eventHistNames[histType] + std::string("[after cut]"); - std::string nameCutfull; - if constexpr (histType == eCent) { - std::string nameCutfull = std::string("Quality assurance of ") + tc.fCentEstm + eventHistNames[histType]; - } else if constexpr (histType == eNumContrib) { - std::string nameCutfull = std::string("Quality assurance of ") + tc.fMultEstm + eventHistNames[histType]; - } else { - std::string nameCutfull = std::string("Quality assurance of ") + eventHistNames[histType]; - } + for (int j = i + 1; j < N; j++) { + + std::string titleY = ""; + if (histType == eCorrCent) { + name = Form("fHist%s[%s vs. %s][%s cut]", corrHistNames[histType], + centEstmNames[i], centEstmNames[j], + cutBeforeAfterNames[ba]); + titleY = Form("%s %s", centEstmNames[j], corrHistNames[histType]); + cr.fCorrHistograms[histType][i][j][ba] = + new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, + maxCent, nBinsCent, minCent, maxCent); + } else if (histType == eCorrMult) { + name = Form("fHist%s[%s vs. %s][%s cut]", corrHistNames[histType], + multEstmNames[i], multEstmNames[j], + cutBeforeAfterNames[ba]); + titleY = Form("%s %s", multEstmNames[j], corrHistNames[histType]); + // y axis bins... + cr.fCorrHistograms[histType][i][j][ba] = + new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, + maxCent, nBinsCent, minCent, maxCent); + } - qa.fQAHistograms[histType][1] = new TH2F(nameCut.c_str(), nameCutfull.c_str(), nBinsCent, minCent, maxCent, nBinsCent, minCent, maxCent); - if constexpr (histType == eNumContrib) { - qa.fQAHistograms[histType][1]->GetYaxis()->SetTitle("NContributors"); - qa.fQAHistograms[histType][1]->GetXaxis()->SetTitle("Reference multiplicity"); - } else { - qa.fQAHistograms[histType][1]->GetYaxis()->SetTitle(Form("Simulated %s", eventHistNames[histType])); - qa.fQAHistograms[histType][1]->GetXaxis()->SetTitle(Form("Reconstructed %s", eventHistNames[histType])); + cr.fCorrHistograms[histType][i][j][ba]->GetYaxis()->SetTitle( + titleY.c_str()); + cr.fCorrHistograms[histType][i][j][ba]->GetXaxis()->SetTitle( + titleX.c_str()); + cr.fCorrHistogramsList->Add(cr.fCorrHistograms[histType][i][j][ba]); + } + } } - qa.fQAHistogramsList->Add(qa.fQAHistograms[histType][1]); } // *) Initialize and book all objects: - void init(InitContext&) - { + void init(InitContext &) { // ... code to book and initialize all analysis objects ... - // *) Set automatically what to process, from an implicit variable "doprocessSomeProcessName" within a PROCESS_SWITCH clause: + // *) Set automatically what to process, from an implicit variable + // "doprocessSomeProcessName" within a PROCESS_SWITCH clause: tc.fProcess[eProcessRec] = doprocessRec; tc.fProcess[eProcessRecSim] = doprocessRecSim; tc.fProcess[eProcessSim] = doprocessSim; @@ -921,6 +1426,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam tc.fSel8CutSwitch = cfSel8CutSwitch; tc.fCentCutSwitch = cfCentCutSwitch; tc.fNumContribCutSwitch = cfNumContribCutSwitch; + tc.fCentCorrCutSwitch = cfCentCorrCutSwitch; + tc.fMultCorrCutSwitch = cfMultCorrCutSwitch; tc.fPtCutSwitch = cfPtCutSwitch; tc.fEtaCutSwitch = cfEtaCutSwitch; @@ -932,6 +1439,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam tc.fVertexZCut = cfVertexZCut; tc.fCentCut = cfCentCut; tc.fNumContribCut = cfNumContribCut; + tc.fCentCorrCut = cfCentCorrCut; + tc.fMultCorrCut = cfMultCorrCut; tc.fPtCut = cfPtCut; tc.fEtaCut = cfEtaCut; @@ -958,34 +1467,51 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam wt.fWeightSwitch = cfWeightSwitch; // *) Book base list: - TList* temp = new TList(); + TList *temp = new TList(); temp->SetOwner(kTRUE); fBaseList.setObject(temp); // *) Book and nest all other TLists: pc.fParticleHistogramsList = new TList(); pc.fParticleHistogramsList->SetName("ParticleHistograms"); - pc.fParticleHistogramsList->SetOwner(kTRUE); - fBaseList->Add(pc.fParticleHistogramsList); // any nested TList in the base TList appears as a subdir in the output ROOT file + pc.fParticleHistogramsList->SetOwner(kFALSE); + fBaseList->Add(pc.fParticleHistogramsList); // any nested TList in the base + // TList appears as a subdir in + // the output ROOT file ev.fEventHistogramsList = new TList(); ev.fEventHistogramsList->SetName("EventHistograms"); - ev.fEventHistogramsList->SetOwner(kTRUE); + ev.fEventHistogramsList->SetOwner(kFALSE); fBaseList->Add(ev.fEventHistogramsList); qa.fQAHistogramsList = new TList(); qa.fQAHistogramsList->SetName("QualityAssuranceHistograms"); - qa.fQAHistogramsList->SetOwner(kTRUE); + qa.fQAHistogramsList->SetOwner(kFALSE); fBaseList->Add(qa.fQAHistogramsList); wt.fWeightHistogramsList = new TList(); wt.fWeightHistogramsList->SetName("WeightHistograms"); - wt.fWeightHistogramsList->SetOwner(kTRUE); + wt.fWeightHistogramsList->SetOwner(kFALSE); fBaseList->Add(wt.fWeightHistogramsList); + cr.fCorrHistogramsList = new TList(); + cr.fCorrHistogramsList->SetName("CorrelationHistograms"); + cr.fCorrHistogramsList->SetOwner(kFALSE); + fBaseList->Add(cr.fCorrHistogramsList); + + mc.fMultiparticleCorrelationProfilesList = new TList(); + mc.fMultiparticleCorrelationProfilesList->SetName( + "MultiparticleCorrelationProfiles"); + mc.fMultiparticleCorrelationProfilesList->SetOwner(kFALSE); + fBaseList->Add(mc.fMultiparticleCorrelationProfilesList); + std::vector> lPcBins = {tc.fPtBins, tc.fPhiBins}; - std::vector> lEvBins = {tc.fCentBins, tc.fMultBins, tc.fVerXBins, tc.fVerYBins, tc.fVerZBins, tc.fNumContribBins}; - std::vector> lQABins = {tc.fCentBins, tc.fMultBins, tc.fNumContribBins}; + std::vector> lEvBins = { + tc.fCentBins, tc.fMultBins, tc.fVerXBins, + tc.fVerYBins, tc.fVerZBins, tc.fNumContribBins}; + std::vector> lQABins = {tc.fCentBins, tc.fMultBins, + tc.fNumContribBins}; + std::vector> lCrBins = {tc.fCentBins, tc.fMultBins}; BookParticleHistograms(lPcBins, pc); BookParticleHistograms(lPcBins, pc); @@ -997,50 +1523,84 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam BookEventHistograms(lEvBins, ev); BookQAHistograms(lQABins, qa); BookQAHistograms(lQABins, qa); + BookCorrHistograms(lCrBins, cr); // if switch on ... + BookCorrHistograms(lCrBins, cr); if (wt.fWeightSwitch) { - wt.fWeightHistograms = getHistogramsWithWeights(tc.fFileWithWeights.c_str(), tc.fRunNumber.c_str()); - for (auto* hist : wt.fWeightHistograms) { + wt.fWeightHistograms = getHistogramsWithWeights( + tc.fFileWithWeights.c_str(), tc.fRunNumber.c_str()); + for (auto *hist : wt.fWeightHistograms) { wt.fWeightHistogramsList->Add(hist); } } + mc.fTwoParticleCorrelationProfiles[eBefore] = + new TProfile("prof2Before", "2-p correlation before cut", 3, 2., 5.); + mc.fTwoParticleCorrelationProfiles[eAfter] = + new TProfile("prof2After", "2-p correlation after cut", 3, 2., 5.); + mc.fTwoParticleCorrelationProfiles[eBefore]->Sumw2(); + mc.fTwoParticleCorrelationProfiles[eAfter]->Sumw2(); + mc.fMultiparticleCorrelationProfilesList->Add( + mc.fTwoParticleCorrelationProfiles[eBefore]); + mc.fMultiparticleCorrelationProfilesList->Add( + mc.fTwoParticleCorrelationProfiles[eAfter]); + + mc.fFourParticleCorrelationProfiles[eBefore] = + new TProfile("prof4Before", "4-p correlation before cut", 2, 3., 5.); + mc.fFourParticleCorrelationProfiles[eAfter] = + new TProfile("prof4After", "4-p correlation after cut", 2, 3., 5.); + mc.fFourParticleCorrelationProfiles[eBefore]->Sumw2(); + mc.fFourParticleCorrelationProfiles[eAfter]->Sumw2(); + mc.fMultiparticleCorrelationProfilesList->Add( + mc.fFourParticleCorrelationProfiles[eBefore]); + mc.fMultiparticleCorrelationProfilesList->Add( + mc.fFourParticleCorrelationProfiles[eAfter]); + } // end of void init(InitContext&) { // A) Process only reconstructed data: - void processRec(CollisionRec const& collision, aod::BCs const&, TracksRec const& tracks) - { + void processRec(CollisionRec const &collision, aod::BCs const &, + TracksRec const &tracks) { // ... // *) Steer all analysis steps: Steer(collision, tracks); } - PROCESS_SWITCH(MultiparticleCumulants, processRec, "process only reconstructed data", true); // yes, keep always one process switch "true", so that there is default running version + PROCESS_SWITCH(MultiparticleCumulants, processRec, + "process only reconstructed data", + true); // yes, keep always one process switch "true", so that + // there is default running version // ------------------------------------------- // B) Process both reconstructed and corresponding MC truth simulated data: - void processRecSim(CollisionRecSim const& collision, aod::BCs const&, TracksRecSim const& tracks, aod::McParticles const&, aod::McCollisions const&) - { + void processRecSim(CollisionRecSim const &collision, aod::BCs const &, + TracksRecSim const &tracks, aod::McParticles const &, + aod::McCollisions const &) { Steer(collision, tracks); } - PROCESS_SWITCH(MultiparticleCumulants, processRecSim, "process both reconstructed and corresponding MC truth simulated data", false); + PROCESS_SWITCH( + MultiparticleCumulants, processRecSim, + "process both reconstructed and corresponding MC truth simulated data", + false); // ------------------------------------------- // C) Process only simulated data: - void processSim(CollisionSim const& /*collision*/, aod::BCs const&, TracksSim const& /*tracks*/) - { - // Steer(collision, tracks); // TBI 20241105 not ready yet, but I do not really need this one urgently, since RecSim is working, and I need that one for efficiencies... + void processSim(CollisionSim const & /*collision*/, aod::BCs const &, + TracksSim const & /*tracks*/) { + // Steer(collision, tracks); // TBI 20241105 not ready yet, but I do + // not really need this one urgently, since RecSim is working, and I need + // that one for efficiencies... } - PROCESS_SWITCH(MultiparticleCumulants, processSim, "process only simulated data", false); + PROCESS_SWITCH(MultiparticleCumulants, processSim, + "process only simulated data", false); }; // struct MultiparticleCumulants { // *) The final touch: -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) -{ +WorkflowSpec defineDataProcessing(ConfigContext const &cfgc) { return WorkflowSpec{ - adaptAnalysisTask(cfgc), + adaptAnalysisTask(cfgc), }; } // WorkflowSpec... From 2ca97c8734dcea1110ac689f6ba143d9f75fdc2a Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Thu, 11 Jun 2026 16:57:23 +0200 Subject: [PATCH 2/9] fix O2linter and errors --- .../Tasks/multiparticleCumulants.cxx | 840 +++++++----------- 1 file changed, 307 insertions(+), 533 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index 1be1f6e0f67..6fdfea1851e 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -52,34 +52,30 @@ using namespace o2; using namespace o2::framework; // Definitions of join tables for Run 3 analysis: -using EventSelection = - soa::Join; // no FV0Ms for cent -using CollisionRec = soa::Join:: - iterator; // use in json "isMC": "true" for "event-selection-task" -using CollisionRecSim = soa::Join::iterator; +using EventSelection = soa::Join; // no FV0Ms for cent +using CollisionRec = soa::Join::iterator; // use in json "isMC": "true" for "event-selection-task" +using CollisionRecSim = soa::Join::iterator; using CollisionSim = aod::McCollision; -using TracksRec = soa::Join; -using TrackRec = soa::Join::iterator; -using TracksRecSim = - soa::Join; // + use in json "isMC" : "true" -using TrackRecSim = - soa::Join::iterator; +using TracksRec = soa::Join; +using TrackRec = soa::Join::iterator; +using TracksRecSim = soa::Join; // + use in json "isMC" : "true" +using TrackRecSim = soa::Join::iterator; using TracksSim = aod::McParticles; using TrackSim = aod::McParticles::iterator; using namespace std; // *) Define enums: -enum EnRlMc { eRl = 0, eMc }; +enum EnRlMc { + eRl = 0, + eMc +}; -enum EnRecSim { eRec = 0, eSim, eRecAndSim }; +enum EnRecSim { + eRec = 0, + eSim, + eRecAndSim +}; enum EnProcess { eProcessRec = 0, // Run 3, only reconstructed @@ -98,179 +94,151 @@ enum EnEventHistograms { eEventHistograms_N }; -const char *eventHistNames[eEventHistograms_N] = {"Centrality", "Multiplicity", - "VertexX", "VertexY", - "VertexZ", "NumContrib"}; +const char *eventHistNames[eEventHistograms_N] = { + "Centrality", + "Multiplicity", + "VertexX", + "VertexY", + "VertexZ", + "NumContrib" +}; -enum EnParticleHistograms { ePt, ePhi, eParticleHistograms_N }; +enum EnParticleHistograms { + ePt, + ePhi, + eParticleHistograms_N +}; -const char *particleHistNames[eParticleHistograms_N] = {"Pt", "Phi"}; +const char *particleHistNames[eParticleHistograms_N] = { + "Pt", + "Phi" +}; -enum EnQAHistograms { eQACent, eQAMultNumContrib, eQAHistograms_N }; +enum EnQAHistograms { + eQACent, + eQAMultNumContrib, + eQAHistograms_N +}; -enum EnCorrHistograms { eCorrCent, eCorrMult, eCorrHistograms_N }; +enum EnCorrHistograms { + eCorrCent, + eCorrMult, + eCorrHistograms_N +}; const char *corrHistNames[eCorrHistograms_N] = { - "Centrality", - "Multiplicity", + "Centrality", + "Multiplicity", }; -enum EnCentEstm { eCentFT0C, eCentFT0M, eCentFV0A, eCentEstm_N }; +enum EnCentEstm { + eCentFT0C, + eCentFT0M, + eCentFV0A, + eCentEstm_N +}; const char *centEstmNames[eCentEstm_N] = { - "FT0C", - "FT0M", - "FV0A", + "FT0C", + "FT0M", + "FV0A", }; -enum EnMultEstm { eMultFT0C, eMultFT0M, eMultFV0A, eMultEstm_N }; +enum EnMultEstm { + eMultFT0C, + eMultFT0M, + eMultFV0A, + eMultEstm_N +}; -const char *multEstmNames[eMultEstm_N] = {"FT0C", "FT0M", "FV0A"}; +const char *multEstmNames[eMultEstm_N] = { + "FT0C", + "FT0M", + "FV0A" +}; -enum EnCutBeforeAfter { eBefore, eAfter, eCutBeforeAfter_N }; +enum EnCutBeforeAfter { + eBefore, + eAfter, + eCutBeforeAfter_N +}; const char *cutBeforeAfterNames[eCutBeforeAfter_N] = { - "before", - "after", + "before", + "after", }; // *) Main task: -struct MultiparticleCumulants { // this name is used in lower-case format to - // name the TDirectoryFile in - // AnalysisResults.root +struct MultiparticleCumulants { // this name is used in lower-case format to name the TDirectoryFile in AnalysisResults.root // *) Base TList to hold all output objects: - TString sBaseListName = - "Default list name"; // yes, I declare it separately, because I need it - // also later in BailOut() function - OutputObj fBaseList{sBaseListName.Data(), - OutputObjHandlingPolicy::AnalysisObject, - OutputObjSourceType::OutputObjSource}; + TString sBaseListName = "Default list name"; // yes, I declare it separately, because I need it also later in BailOut() function + OutputObj fBaseList{sBaseListName.Data(), OutputObjHandlingPolicy::AnalysisObject, OutputObjSourceType::OutputObjSource}; // *) Service: - Service - ccdb; // support for offline callibration data base + Service ccdb; // support for offline callibration data base Service pdg; // *) Define configurables: - Configurable cfDryRun{ - "cfDryRun", false, - "book all histos and run without filling and calculating anything"}; - Configurable cfCentEstm{"cfCentEstm", "FT0M", - "centrality estimator: TBD"}; - Configurable cfMultEstm{"cfMultEstm", "FV0A", - "multiplicity estimator: TBD"}; + Configurable cfDryRun{"cfDryRun", false, "book all histos and run without filling and calculating anything"}; + Configurable cfCentEstm{"cfCentEstm", "FT0M", "centrality estimator: TBD"}; + Configurable cfMultEstm{"cfMultEstm", "FV0A", "multiplicity estimator: TBD"}; Configurable cfQASwitch{"cfQASwitch", true, "quality assurance switch"}; Configurable cfWeightSwitch{"cfWeightSwitch", true, "weight switch"}; - Configurable cfPrintSwitch{"cfPrintSwitch", true, - "printing result switch"}; + Configurable cfPrintSwitch{"cfPrintSwitch", true, "printing result switch"}; // *) Event cut switches - Configurable cfVertexZCutSwitch{"cfVertexZCutSwitch", true, - "vertex z cut switch"}; - Configurable cfSel8CutSwitch{"cfSel8CutSwitch", true, - "Sel8 cut switch"}; - Configurable cfCentCutSwitch{"cfCentCutSwitch", true, - "centrality cut switch"}; - Configurable cfNumContribCutSwitch{"cfNumContribCutSwitch", true, - "NContribution cut switch"}; - Configurable cfCentCorrCutSwitch{"cfCentCorrCutSwitch", true, - "centrality correlation cut switch"}; - Configurable cfMultCorrCutSwitch{"cfMultCorrCutSwitch", true, - "multiplicity correlation cut switch"}; + Configurable cfVertexZCutSwitch{"cfVertexZCutSwitch", true, "vertex z cut switch"}; + Configurable cfSel8CutSwitch{"cfSel8CutSwitch", true, "Sel8 cut switch"}; + Configurable cfCentCutSwitch{"cfCentCutSwitch", true, "centrality cut switch"}; + Configurable cfNumContribCutSwitch{"cfNumContribCutSwitch", true, "NContribution cut switch"}; + Configurable cfCentCorrCutSwitch{"cfCentCorrCutSwitch", true, "centrality correlation cut switch"}; + Configurable cfMultCorrCutSwitch{"cfMultCorrCutSwitch", true, "multiplicity correlation cut switch"}; // *) Particle cut switches Configurable cfPtCutSwitch{"cfPtCutSwitch", true, "Pt cut switch"}; Configurable cfEtaCutSwitch{"cfEtaCutSwitch", true, "Eta cut switch"}; - Configurable cfSignCutSwitch{"cfSignCutSwitch", true, - "Charge cut switch"}; - Configurable cfTpcNClsFoundCutSwitch{"cfTpcNClsFoundCutSwitch", true, - ""}; + Configurable cfSignCutSwitch{"cfSignCutSwitch", true, "Charge cut switch"}; + Configurable cfTpcNClsFoundCutSwitch{"cfTpcNClsFoundCutSwitch", true, ""}; Configurable cfDCAXYCutSwitch{"cfDCAXYCutSwitch", true, ""}; Configurable cfDCAZCutSwitch{"cfDCAZCutSwitch", true, ""}; // *) Event cut - Configurable> cfVertexZCut{ - "cfVertexZCut", {-10., 10.}, "vertex z position range: {min, max}[cm]"}; - Configurable> cfCentCut{ - "cfCentCut", {10., 20.}, "centrality range: {min, max}[%]"}; - Configurable> cfNumContribCut{ - "cfNumContribCut", {0, 3000.}, "NContribution range: {min, max}"}; - Configurable cfCentCorrCut{"cfCentCorrCut", 5, - "centrality difference maximum"}; - Configurable cfMultCorrCut{"cfMultCorrCut", 0.2, - "multiplicity difference maximum"}; + Configurable> cfVertexZCut{"cfVertexZCut", {-10., 10.}, "vertex z position range: {min, max}[cm]"}; + Configurable> cfCentCut{"cfCentCut", {10., 20.}, "centrality range: {min, max}[%]"}; + Configurable> cfNumContribCut{"cfNumContribCut", {0, 3000.}, "NContribution range: {min, max}"}; + Configurable cfCentCorrCut{"cfCentCorrCut", 5, "centrality difference maximum"}; + Configurable cfMultCorrCut{"cfMultCorrCut", 0.2, "multiplicity difference maximum"}; // *) Particle cut - Configurable> cfPtCut{ - "cfPtCut", - {0.2, 5.0}, - "Pt range: {min, max}[GeV], with convention: min <= Pt < max"}; - Configurable> cfEtaCut{ - "cfEtaCut", - {-0.8, 0.8}, - "Eta range: {min, max}, with convention: min <= Eta < max"}; - Configurable> cfSignCut{ - "cfSignCut", - {1, 0, 1}, - "sign of charge, 1 to keep and 0 to discard, {negative, neutral, " - "positive}"}; - Configurable> cfTpcNClsFoundCut{ - "cfTpcNClsFoundCut", - {70., 160.}, - "range of found TPC clusters for this track geometry: {min, max}"}; - Configurable> cfDCAXYCut{ - "cfDCAXYCut", - {-3.2, 3.2}, - "range of distance-of-closest-approach (DCA) of the extrapolated track " - "to the primary position in XY-direction: {min, max}[cm]"}; - Configurable> cfDCAZCut{ - "cfDCAZCut", - {-2.4, 2.4}, - "range of distance-of-closest-approach (DCA) of the extrapolated track " - "to the primary position in Z-direction: {min, max}[cm]"}; - - // *) - Configurable cfFileWithWeights{ - "cfFileWithWeights", - "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root", - "path to external ROOT file which holds all particle weights in O2 " - "format"}; - Configurable cfRunNumber{"cfRunNumber", "000123456", - "run number"}; + Configurable> cfPtCut{"cfPtCut", {0.2, 5.0}, "Pt range: {min, max}[GeV], with convention: min <= Pt < max"}; + Configurable> cfEtaCut{"cfEtaCut", {-0.8, 0.8}, "Eta range: {min, max}, with convention: min <= Eta < max"}; + Configurable> cfSignCut{"cfSignCut", {1, 0, 1}, "sign of charge, 1 to keep and 0 to discard, {negative, neutral, positive}"}; + Configurable> cfTpcNClsFoundCut{"cfTpcNClsFoundCut", {70., 160.}, "range of found TPC clusters for this track geometry: {min, max}"}; + Configurable> cfDCAXYCut{"cfDCAXYCut", {-3.2, 3.2}, "range of distance-of-closest-approach (DCA) of the extrapolated track to the primary position in XY-direction: {min, max}[cm]"}; + Configurable> cfDCAZCut{"cfDCAZCut", {-2.4, 2.4}, "range of distance-of-closest-approach (DCA) of the extrapolated track to the primary position in Z-direction: {min, max}[cm]"}; + + // *) Others + Configurable cfFileWithWeights{"cfFileWithWeights", "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root", "path to external ROOT file which holds all particle weights in O2 format"}; + Configurable cfRunNumber{"cfRunNumber", "000123456", "run number"}; // *) Bins - Configurable> cfPtBins{ - "cfPtBins", {1000, 0., 100.}, "nPtBins, ptMin, ptMax"}; - Configurable> cfPhiBins{ - "cfPhiBins", - {1000, 0., o2::constants::math::TwoPI}, - "nPhiBins, phiMin, phiMax"}; - Configurable> cfCentBins{ - "cfCentBins", {100, 0., 100.}, "nCenBins, cenMin, cenMax"}; - Configurable> cfMultBins{ - "cfMultBins", {100, 0., 5000.}, "nMultBins, MultMin, MultMax"}; - Configurable> cfVerXBins{ - "cfVerXBins", {100, -0.05, 0.05}, "nVerXBins, VerXMin, VerXMax"}; - Configurable> cfVerYBins{ - "cfVerYBins", {100, -0.05, 0.05}, "nVerYBins, VerYMin, VerYMax"}; - Configurable> cfVerZBins{ - "cfVerZBins", {100, -50., 50.}, "nVerZBins, VerZMin, VerZMax"}; - Configurable> cfNumContribBins{ - "cfNumContribBins", - {100, 0., 5000.}, - "nNumContribBins, NumContribMin, NumContribMax"}; - - // *) Define and initialize all data members to be called in the main process* - // functions: + Configurable> cfPtBins{"cfPtBins", {1000, 0., 100.}, "nPtBins, ptMin, ptMax"}; + Configurable> cfPhiBins{"cfPhiBins", {1000, 0., o2::constants::math::TwoPI}, "nPhiBins, phiMin, phiMax"}; + Configurable> cfCentBins{"cfCentBins", {100, 0., 100.}, "nCenBins, cenMin, cenMax"}; + Configurable> cfMultBins{"cfMultBins", {100, 0., 5000.}, "nMultBins, MultMin, MultMax"}; + Configurable> cfVerXBins{"cfVerXBins", {100, -0.05, 0.05}, "nVerXBins, VerXMin, VerXMax"}; + Configurable> cfVerYBins{"cfVerYBins", {100, -0.05, 0.05}, "nVerYBins, VerYMin, VerYMax"}; + Configurable> cfVerZBins{"cfVerZBins", {100, -50., 50.}, "nVerZBins, VerZMin, VerZMax"}; + Configurable> cfNumContribBins{"cfNumContribBins", {100, 0., 5000.}, "nNumContribBins, NumContribMin, NumContribMax"}; + + // *) Define and initialize all data members to be called in the main process* functions: // **) Task configuration: struct TaskConfiguration { - bool fProcess[eProcess_N] = { - false}; // Set what to process. See enum EnProcess for full description. - // Set via implicit variables within a PROCESS_SWITCH clause. - bool fDryRun = false; // book all histos and run without filling and - // calculating anything + bool fProcess[eProcess_N] = {false}; // Set what to process. See enum EnProcess for full description. Set via implicit variables within a PROCESS_SWITCH clause. + bool fDryRun = false; // book all histos and run without filling and calculating anything + std::string fCentEstm = "FT0M"; std::string fMultEstm = "FV0A"; @@ -313,41 +281,30 @@ struct MultiparticleCumulants { // this name is used in lower-case format to std::vector fVerZBins = {-50., 50.}; std::vector fNumContribBins = {0, 5000}; - std::string fFileWithWeights = - "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root"; + std::string fFileWithWeights = "/scratch3/go52dab/O2tutorial/tutorial3-6/weights.root"; std::string fRunNumber = "000123456"; - } tc; // you have to prepend "tc." for all objects name in this group later in - // the code + } tc; - // **) Particle histograms: struct ParticleHistograms { - TList *fParticleHistogramsList = - NULL; //! [before, after][v2^2, v3^2, v4^2] - float fFourParticleCorrelationEbye[eCutBeforeAfter_N][2] = { - {0., 0.}, {0., 0.}}; //<4> [before, after][v3^2 v2^2, v4^2 v2^2] + float fTwoParticleCorrelationEbye[eCutBeforeAfter_N][3] = {{0., 0., 0.}, {0., 0., 0.}}; //<2> [before, after][v2^2, v3^2, v4^2] + float fFourParticleCorrelationEbye[eCutBeforeAfter_N][2] = {{0., 0.}, {0., 0.}}; //<4> [before, after][v3^2 v2^2, v4^2 v2^2] } ebye; struct MultiparticleCorrelationProfile { @@ -378,19 +333,15 @@ struct MultiparticleCumulants { // this name is used in lower-case format to struct MultiparticleCorrelationCalculation { int h1, h2, h3, h4, h5, h6, h7, h8; // Book Q-vector components: - static constexpr int MaxCorrelator = 2; //<> - static constexpr int MaxHarmonic = 5; // n+1 in vn + static constexpr int MaxCorrelator = 4; //<> + static constexpr int MaxHarmonic = 5; // n+1 in vn, in this case n=4 as we need v2, v3, v4 static constexpr int MaxPower = MaxCorrelator + 1; + static constexpr int NumSC = 2; //need SC(3,2) and SC(4,2) TComplex fQvectorBefore[MaxHarmonic][MaxPower]; - TComplex fQvectorAfter[MaxHarmonic] - [MaxPower]; // All needed Q-vector components - // Store the final results here: - // Remark: [2][MaxCorrelator] => [Cos,Sin][<2>,<3>,<4>,<5>,<6>,<7>,<8>] + TComplex fQvectorAfter[MaxHarmonic][MaxPower]; // All needed Q-vector components } mcc; - template - bool ctEventCuts(T1 const &collision, T2 rlCollisionCentAll, - T3 rlCollisionMultAll) { + template bool ctEventCuts(T1 const& collision, T2 rlCollisionCentAll, T3 rlCollisionMultAll) { bool pass = true; bool bVertexZCut = true; @@ -466,7 +417,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to return pass; } - template bool ctParticleCuts(T1 const &track) { + template bool ctParticleCuts(T1 const& track) { bool pass = true; bool bPtCut = true; @@ -487,9 +438,10 @@ struct MultiparticleCumulants { // this name is used in lower-case format to (track.sign() == 1 && tc.fSignCut[2]); bTpcNClsFoundCut = track.tpcNClsFound() < tc.fTpcNClsFoundCut[1] && track.tpcNClsFound() > tc.fTpcNClsFoundCut[0]; - bDCAXYCut = - track.dcaXY() < tc.fDCAXYCut[1] && track.dcaXY() > tc.fDCAXYCut[0]; - bDCAZCut = track.dcaZ() < tc.fDCAZCut[1] && track.dcaZ() > tc.fDCAZCut[0]; + bDCAXYCut = track.dcaXY() < tc.fDCAXYCut[1] && + track.dcaXY() > tc.fDCAXYCut[0]; + bDCAZCut = track.dcaZ() < tc.fDCAZCut[1] && + track.dcaZ() > tc.fDCAZCut[0]; } // *) For mc event only @@ -544,16 +496,13 @@ struct MultiparticleCumulants { // this name is used in lower-case format to return TComplex::Conjugate(mcc.fQvectorAfter[-n][p]); } - } // TComplex Q(int n, int p) + } - TComplex mccRecursion(int n, int *harmonic, EnCutBeforeAfter eba, - int mult = 1, int skip = 0) { - // Calculate multi-particle correlators by using recursion (an improved - // faster version) originally developed by Kristjan Gulbrandsen - // (gulbrand@nbi.dk). + TComplex mccRecursion(int n, int *harmonic, EnCutBeforeAfter eba, int mult = 1, int skip = 0) { + // Calculate multi-particle correlators by using recursion (an improved faster version) originally developed by Kristjan Gulbrandsen (gulbrand@nbi.dk). int nm1 = n - 1; - TComplex c(Q(harmonic[nm1], mult, eba)); + TComplex c(mccQ(harmonic[nm1], mult, eba)); if (nm1 == 0) return c; c *= mccRecursion(nm1, harmonic, eba); @@ -585,8 +534,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to return c - c2; return c - double(mult) * c2; - } // TComplex AliFlowAnalysisWithMultiparticleCorrelations::mccRecursion(int - // n, int* harmonic, int mult = 1, int skip = 0) + } TObject *getObjectFromList(TList *list, const char *objectName) { // Get TObject pointer from TList, even if it's in some nested TList. @@ -619,12 +567,10 @@ struct MultiparticleCumulants { // this name is used in lower-case format to // Otherwise, search for the object recursively in the nested lists: TObject *objectIter; // iterator object in the loop below TIter next(list); - while ((objectIter = - next())) // double round braces are to silence the warnings + while ((objectIter = next())) // double round braces are to silence the warnings { if (TString(objectIter->ClassName()).EqualTo("TList")) { - objectFinal = getObjectFromList(reinterpret_cast(objectIter), - objectName); + objectFinal = getObjectFromList(reinterpret_cast(objectIter), objectName); if (objectFinal) return objectFinal; } @@ -634,23 +580,17 @@ struct MultiparticleCumulants { // this name is used in lower-case format to } // TObject* getObjectFromList(TList *list, char *objectName) - std::vector getHistogramsWithWeights(const char *filePath, - const char *runNumber) { + std::vector getHistogramsWithWeights(const char *filePath, const char *runNumber) { // a) Return value: std::vector histograms; - TList *baseList = - NULL; // base top-level list in the TFile, e.g. named "ccdb_object" - TList *listWithRuns = - NULL; // nested list with run-wise TList's holding run-specific weights - - // c) Determine from filePath if the file in on a local machine, or in home - // dir AliEn, or in CCDB: - // Algorithm: If filePath begins with "/alice/data/CCDB/" then it's in - // home dir AliEn. - // If filePath begins with "/alice-ccdb.cern.ch/" then - // it's in CCDB. Therefore, files in AliEn and CCDB must - // be specified with abs path, for local files both abs - // and relative paths are just fine. + TList *baseList = NULL; // base top-level list in the TFile, e.g. named "ccdb_object" + TList *listWithRuns = NULL; // nested list with run-wise TList's holding run-specific weights + + // c) Determine from filePath if the file in on a local machine, or in home dir AliEn, or in CCDB: + // Algorithm: + // If filePath begins with "/alice/data/CCDB/" then it's in home dir AliEn. + // If filePath begins with "/alice-ccdb.cern.ch/" then it's in CCDB. + // Therefore, files in AliEn and CCDB must be specified with abs path, for local files both abs and relative paths are just fine. bool bFileIsInAliEn = false; bool bFileIsInCCDB = false; @@ -664,17 +604,11 @@ struct MultiparticleCumulants { // this name is used in lower-case format to if (bFileIsInAliEn) { // File you want to access is in your home dir in AliEn: - TGrid *alien = - TGrid::Connect("alien", gSystem->Getenv("USER"), "", - ""); // do not forget to add #include to the - // preamble of your analysis task + TGrid *alien = TGrid::Connect("alien", gSystem->Getenv("USER"), "", ""); // do not forget to add #include to the preamble of your analysis task if (!alien) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - TFile *weightsFile = TFile::Open( - Form("alien://%s", filePath), - "READ"); // yes, ROOT can open a file transparently, even if it's - // sitting in AliEn, with this specific syntax + TFile *weightsFile = TFile::Open(Form("alien://%s", filePath), "READ"); // yes, ROOT can open a file transparently, even if it's sitting in AliEn, with this specific syntax if (!weightsFile) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -684,31 +618,22 @@ struct MultiparticleCumulants { // this name is used in lower-case format to LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - // Finally, from the top-level TList, get the desired nested TList => the - // technical problem here is that it can be nested at any level, for thare - // there is a helper utility function getObjectFromList(...) , see its - // implementation further below - listWithRuns = - reinterpret_cast(getObjectFromList(baseList, runNumber)); + // Finally, from the top-level TList, get the desired nested TList => the technical problem here is that it can be nested at any level, + // for thare there is a helper utility function GetObjectFromList(...) , see its implementation further below + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; - runNumberWithLeadingZeroes += - runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast( - getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } } } else if (bFileIsInCCDB) { // File you want to access is in your home dir in CCDB: - // Remember that here I do not access the file; instead, I directly access - // the object in that file. My home dir in CCDB: - // https://alice-ccdb.cern.ch/browse/Users/a/abilandz/ => adapt for - // your case - ccdb->setURL( - "https://alice-ccdb.cern.ch"); // to be able to use "ccdb" this object - // in your analysis task, see 4b/ below + // Remember that here I do not access the file; instead, I directly access the object in that file. + // My home dir in CCDB: https://alice-ccdb.cern.ch/browse/Users/a/abilandz/ => adapt for your case + ccdb->setURL("https://alice-ccdb.cern.ch"); // to be able to use "ccdb" this object in your analysis task, see 4b/ below baseList = reinterpret_cast(ccdb->get( TString(filePath).ReplaceAll("/alice-ccdb.cern.ch/", "").Data())); baseList->ls(); @@ -716,14 +641,11 @@ struct MultiparticleCumulants { // this name is used in lower-case format to LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - listWithRuns = - reinterpret_cast(getObjectFromList(baseList, runNumber)); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; - runNumberWithLeadingZeroes += - runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast( - getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -738,9 +660,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to if (gSystem->AccessPathName(filePath, kFileExists)) { LOGF(info, - "\033[1;33m if(gSystem->AccessPathName(filePath,kFileExists)), " - "filePath = %s \033[0m", - filePath); + "\033[1;33m if(gSystem->AccessPathName(filePath,kFileExists)), filePath = %s \033[0m", filePath); LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -756,21 +676,15 @@ struct MultiparticleCumulants { // this name is used in lower-case format to LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - listWithRuns = reinterpret_cast(getObjectFromList( - baseList, runNumber)); // baseList->FindObject(runNumber) + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; - runNumberWithLeadingZeroes += - runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast( - getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { baseList->ls(); LOGF(fatal, - "\033[1;31m%s at line %d : this crash can happen if in the " - "output file there is no list with weights for the current run " - "number = %s\033[0m", - __FUNCTION__, __LINE__, runNumber); + "\033[1;31m%s at line %d : this crash can happen if in the output file there is no list with weights for the current runnumber = %s\033[0m", __FUNCTION__, __LINE__, runNumber); } } } @@ -804,7 +718,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to // *) Define all member functions to be called in the main process* functions: template - void Steer(T1 const &collision, T2 const &tracks) { + void runLoop(T1 const& collision, T2 const& tracks) { // Dry run: if (tc.fDryRun) { @@ -817,8 +731,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to for (int p = 0; p < mcc.MaxPower; p++) { mcc.fQvectorBefore[h][p] = TComplex(0., 0.); mcc.fQvectorAfter[h][p] = TComplex(0., 0.); - } // for(int p=0;p(collision.multFT0C()), - static_cast(collision.multFT0M()), - static_cast(collision.multFV0A())}; + static_cast(collision.multFT0C()), + static_cast(collision.multFT0M()), + static_cast(collision.multFV0A()) + }; float rlCollisionMult = 0.; for (int i = 0; i < eMultEstm_N; i++) { if (tc.fPrintSwitch) { - LOGF(info, "%s Multiplicity: %f", multEstmNames[i], - rlCollisionMultAll[i]); + LOGF(info, "%s Multiplicity: %f", multEstmNames[i], rlCollisionMultAll[i]); } if (tc.fMultEstm == multEstmNames[i]) { rlCollisionMult = rlCollisionMultAll[i]; @@ -886,38 +803,27 @@ struct MultiparticleCumulants { // this name is used in lower-case format to ev.fEventHistograms[eVertexX][eRec][eBefore]->Fill(collision.posX()); ev.fEventHistograms[eVertexY][eRec][eBefore]->Fill(collision.posY()); ev.fEventHistograms[eVertexZ][eRec][eBefore]->Fill(collision.posZ()); - ev.fEventHistograms[eNumContrib][eRec][eBefore]->Fill( - rlCollisionNumContrib); + ev.fEventHistograms[eNumContrib][eRec][eBefore]->Fill(rlCollisionNumContrib); - // if (tc.fCentCorrCutSwitch) { for (int i = 0; i < eCentEstm_N; i++) { for (int j = i + 1; j < eCentEstm_N; j++) { auto *h = cr.fCorrHistograms[eCorrCent][i][j][eBefore]; if (!h) { - LOGF(fatal, - "Missing histogram " - "cr.fCorrHistograms[eCorrCent][%d][%d][eBefore]", - i, j); + LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrCent][%d][%d][eBefore]", i, j); } h->Fill(rlCollisionCentAll[i], rlCollisionCentAll[j]); } } - // } - // if (tc.fMultCorrCutSwitch) { for (int i = 0; i < eMultEstm_N; i++) { for (int j = i + 1; j < eMultEstm_N; j++) { auto *h = cr.fCorrHistograms[eCorrMult][i][j][eBefore]; if (!h) { - LOGF(fatal, - "Missing histogram " - "cr.fCorrHistograms[eCorrMult][%d][%d][eBefore]", - i, j); + LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrMult][%d][%d][eBefore]", i, j); } h->Fill(rlCollisionMultAll[i], rlCollisionMultAll[j]); } } - // } if constexpr (rs == eRecAndSim) { if (!collision.has_mcCollision()) { @@ -927,12 +833,11 @@ struct MultiparticleCumulants { // this name is used in lower-case format to return; } - auto mccollision = collision.mcCollision(); // McCollisionLabels + auto mccollision = collision.mcCollision(); float mcCollisionCent = 0.; - float b = mccollision.impactParameter() * - std::pow(10, -15); // convert fm to m + float b = mccollision.impactParameter() * std::pow(10, -15); // convert fm to m float xs = 7.71 * std::pow(10, -28); // convert barn to m^2 mcCollisionCent = o2::constants::math::PI * b * b / xs * 100; @@ -949,8 +854,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to ev.fEventHistograms[eVertexY][eSim][eBefore]->Fill(mccollision.posY()); ev.fEventHistograms[eVertexZ][eSim][eBefore]->Fill(mccollision.posZ()); - if (ctEventCuts(mccollision, rlCollisionCentAll, - rlCollisionMultAll)) { + if (ctEventCuts(mccollision, rlCollisionCentAll, rlCollisionMultAll)) { ev.fEventHistograms[eCent][eSim][eAfter]->Fill(mcCollisionCent); ev.fEventHistograms[eVertexX][eSim][eAfter]->Fill(mccollision.posX()); ev.fEventHistograms[eVertexY][eSim][eAfter]->Fill(mccollision.posY()); @@ -958,18 +862,12 @@ struct MultiparticleCumulants { // this name is used in lower-case format to } if (qa.fQASwitch) { - qa.fQAHistograms[eQACent][eBefore]->Fill(rlCollisionCent, - mcCollisionCent); - qa.fQAHistograms[eQAMultNumContrib][eBefore]->Fill( - rlCollisionMult, rlCollisionNumContrib); - if (ctEventCuts(collision, rlCollisionCentAll, - rlCollisionMultAll) && - ctEventCuts(mccollision, rlCollisionCentAll, - rlCollisionMultAll)) { - qa.fQAHistograms[eQACent][eAfter]->Fill(rlCollisionCent, - mcCollisionCent); - qa.fQAHistograms[eQAMultNumContrib][eAfter]->Fill( - rlCollisionMult, rlCollisionNumContrib); + qa.fQAHistograms[eQACent][eBefore]->Fill(rlCollisionCent, mcCollisionCent); + qa.fQAHistograms[eQAMultNumContrib][eBefore]->Fill(rlCollisionMult, rlCollisionNumContrib); + if (ctEventCuts(collision, rlCollisionCentAll, rlCollisionMultAll) && + ctEventCuts(mccollision, rlCollisionCentAll, rlCollisionMultAll)) { + qa.fQAHistograms[eQACent][eAfter]->Fill(rlCollisionCent, mcCollisionCent); + qa.fQAHistograms[eQAMultNumContrib][eAfter]->Fill(rlCollisionMult, rlCollisionNumContrib); } } } @@ -980,18 +878,14 @@ struct MultiparticleCumulants { // this name is used in lower-case format to ev.fEventHistograms[eVertexX][eRec][eAfter]->Fill(collision.posX()); ev.fEventHistograms[eVertexY][eRec][eAfter]->Fill(collision.posY()); ev.fEventHistograms[eVertexZ][eRec][eAfter]->Fill(collision.posZ()); - ev.fEventHistograms[eNumContrib][eRec][eAfter]->Fill( - rlCollisionNumContrib); + ev.fEventHistograms[eNumContrib][eRec][eAfter]->Fill(rlCollisionNumContrib); if (tc.fCentCorrCutSwitch) { for (int i = 0; i < eCentEstm_N; i++) { for (int j = i + 1; j < eCentEstm_N; j++) { auto *h = cr.fCorrHistograms[eCorrCent][i][j][eAfter]; if (!h) { - LOGF(fatal, - "Missing histogram " - "cr.fCorrHistograms[eCorrCent][%d][%d][eAfter]", - i, j); + LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrCent][%d][%d][eAfter]", i, j); } h->Fill(rlCollisionCentAll[i], rlCollisionCentAll[j]); } @@ -1004,10 +898,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to for (int j = i + 1; j < eMultEstm_N; j++) { auto *h = cr.fCorrHistograms[eCorrMult][i][j][eAfter]; if (!h) { - LOGF(fatal, - "Missing histogram " - "cr.fCorrHistograms[eCorrMult][%d][%d][eAfter]", - i, j); + LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrMult][%d][%d][eAfter]", i, j); } h->Fill(rlCollisionMultAll[i], rlCollisionMultAll[j]); } @@ -1027,7 +918,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to double wPhiToPowerP = 1.; // particle weight raised to power p // Main loop over particles: - for (auto const &track : tracks) { + for (auto const& track : tracks) { // LOGF(info, "Track azimuthal angle: %f", track.phi()); // LOGF(info, "Transverse momentum: %f", track.pt()); @@ -1041,8 +932,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to dPhi = track.phi(); if (wt.fWeightSwitch) { auto *hist = wt.fWeightHistograms[1]; - wPhi = hist->GetBinContent( - wt.fWeightHistograms[1]->GetXaxis()->FindBin(dPhi)); + wPhi = hist->GetBinContent(wt.fWeightHistograms[1]->GetXaxis()->FindBin(dPhi)); } for (int h = 0; h < mcc.MaxHarmonic; h++) { @@ -1050,11 +940,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to if (wt.fWeightSwitch) { wPhiToPowerP = std::pow(wPhi, p); } - mcc.fQvectorBefore[h][p] += - TComplex(wPhiToPowerP * o2::constants::math::Cos(h * dPhi), - wPhiToPowerP * o2::constants::math::Sin(h * dPhi)); - } // for(int p=0;p(track)) { pc.fParticleHistograms[ePt][eRec][eAfter]->Fill(track.pt()); @@ -1065,11 +953,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to if (wt.fWeightSwitch) { wPhiToPowerP = std::pow(wPhi, p); } - mcc.fQvectorAfter[h][p] += - TComplex(wPhiToPowerP * std::cos(h * dPhi), - wPhiToPowerP * std::sin(h * dPhi)); - } // for(int p=0;pFill(mcparticle.pt()); pc.fParticleHistograms[ePhi][eSim][eBefore]->Fill(mcparticle.phi()); if (ctParticleCuts(mcparticle)) { @@ -1101,38 +986,29 @@ struct MultiparticleCumulants { // this name is used in lower-case format to } // end of for (int64_t i = 0; i < tracks.size(); i++) { - for (int i = 0; i < 3; i++) { + for (int i = 0; i < mcc.MaxHarmonic-2 ; i++) { mcc.h1 = -(i + 2); mcc.h2 = i + 2; int harmonicsTwoNum[2] = {mcc.h1, mcc.h2}; int harmonicsTwoDen[2] = {0, 0}; - TComplex twoRecursionBefore = - mccRecursion(2, harmonicsTwoNum, eBefore) / - mccRecursion(2, harmonicsTwoDen, eBefore).Re(); - double wTwoRecursionBefore = - mccRecursion(2, harmonicsTwoDen, eBefore).Re(); + TComplex twoRecursionBefore = mccRecursion(2, harmonicsTwoNum, eBefore) / mccRecursion(2, harmonicsTwoDen, eBefore).Re(); + double wTwoRecursionBefore = mccRecursion(2, harmonicsTwoDen, eBefore).Re(); ebye.fTwoParticleCorrelationEbye[eBefore][i] = twoRecursionBefore.Re(); - TComplex twoRecursionAfter = - mccRecursion(2, harmonicsTwoNum, eAfter) / - mccRecursion(2, harmonicsTwoDen, eAfter).Re(); + TComplex twoRecursionAfter = mccRecursion(2, harmonicsTwoNum, eAfter) / mccRecursion(2, harmonicsTwoDen, eAfter).Re(); double wTwoRecursionAfter = mccRecursion(2, harmonicsTwoDen, eAfter).Re(); ebye.fTwoParticleCorrelationEbye[eAfter][i] = twoRecursionAfter.Re(); if (nTracksBefore > 1) { - mc.fTwoParticleCorrelationProfiles[eBefore]->Fill( - i + 2.5, ebye.fTwoParticleCorrelationEbye[eBefore][i], - wTwoRecursionBefore); + mc.fTwoParticleCorrelationProfiles[eBefore]->Fill(i + 2.5, ebye.fTwoParticleCorrelationEbye[eBefore][i], wTwoRecursionBefore); } if (nTracksAfter > 1) { - mc.fTwoParticleCorrelationProfiles[eAfter]->Fill( - i + 2.5, ebye.fTwoParticleCorrelationEbye[eAfter][i], - wTwoRecursionAfter); + mc.fTwoParticleCorrelationProfiles[eAfter]->Fill(i + 2.5, ebye.fTwoParticleCorrelationEbye[eAfter][i], wTwoRecursionAfter); } } - for (int i = 0; i < 2; i++) { + for (int i = 0; i < mcc.NumSC; i++) { // 4-p correlations: mcc.h1 = -(i + 3); mcc.h2 = -2; @@ -1141,143 +1017,94 @@ struct MultiparticleCumulants { // this name is used in lower-case format to int harmonicsFourNum[4] = {mcc.h1, mcc.h2, mcc.h3, mcc.h4}; int harmonicsFourDen[4] = {0, 0, 0, 0}; - TComplex fourRecursionBefore = - mccRecursion(4, harmonicsFourNum, eBefore) / - mccRecursion(4, harmonicsFourDen, eBefore).Re(); - double wFourRecursionBefore = - mccRecursion(4, harmonicsFourDen, eBefore).Re(); + TComplex fourRecursionBefore = mccRecursion(4, harmonicsFourNum, eBefore) / mccRecursion(4, harmonicsFourDen, eBefore).Re(); + double wFourRecursionBefore = mccRecursion(4, harmonicsFourDen, eBefore).Re(); ebye.fFourParticleCorrelationEbye[eBefore][i] = fourRecursionBefore.Re(); - TComplex fourRecursionAfter = - mccRecursion(4, harmonicsFourNum, eAfter) / - mccRecursion(4, harmonicsFourDen, eAfter).Re(); - double wFourRecursionAfter = - mccRecursion(4, harmonicsFourDen, eAfter).Re(); + TComplex fourRecursionAfter = mccRecursion(4, harmonicsFourNum, eAfter) / mccRecursion(4, harmonicsFourDen, eAfter).Re(); + double wFourRecursionAfter = mccRecursion(4, harmonicsFourDen, eAfter).Re(); ebye.fFourParticleCorrelationEbye[eAfter][i] = fourRecursionAfter.Re(); - if (nTracksBefore > 3) { - mc.fFourParticleCorrelationProfiles[eBefore]->Fill( - i + 3.5, ebye.fFourParticleCorrelationEbye[eBefore][i], - wFourRecursionBefore); + if (nTracksBefore > mcc.MaxCorrelator) { + mc.fFourParticleCorrelationProfiles[eBefore]->Fill(i + 3.5, ebye.fFourParticleCorrelationEbye[eBefore][i], wFourRecursionBefore); } - if (nTracksAfter > 3) { - mc.fFourParticleCorrelationProfiles[eAfter]->Fill( - i + 3.5, ebye.fFourParticleCorrelationEbye[eAfter][i], - wFourRecursionAfter); + if (nTracksAfter > mcc.MaxCorrelator) { + mc.fFourParticleCorrelationProfiles[eAfter]->Fill(i + 3.5, ebye.fFourParticleCorrelationEbye[eAfter][i], wFourRecursionAfter); } } - } // end of template void Steer(T1 - // const& collision, T2 const& tracks) { + } template - void BookParticleHistograms(T1 const &lPcBins, ParticleHistograms &pc) { - std::vector lPtBins = - lPcBins[histType]; // define local array and initialize it from an array - // set in the configurables + void bookParticleHistograms(T1 const& lPcBins, ParticleHistograms &pc) { + std::vector lPtBins = lPcBins[histType]; // define local array and initialize it from an array set in the configurables int nBinsPt = static_cast(lPtBins[0]); float minPt = lPtBins[1]; float maxPt = lPtBins[2]; for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { - std::string nameRec = - Form("fHist%s[eRec][%s cut]", particleHistNames[histType], - cutBeforeAfterNames[ba]); - std::string nameSim = - Form("fHist%s[eSim][%s cut]", particleHistNames[histType], - cutBeforeAfterNames[ba]); - std::string nameRecfull = - Form("%s distribution for reconstructed particles", - particleHistNames[histType]); - std::string nameSimfull = Form("%s distribution for simulated particles", - particleHistNames[histType]); + std::string nameRec = Form("fHist%s[eRec][%s cut]", particleHistNames[histType], cutBeforeAfterNames[ba]); + std::string nameSim = Form("fHist%s[eSim][%s cut]", particleHistNames[histType], cutBeforeAfterNames[ba]); + std::string nameRecfull = Form("%s distribution for reconstructed particles", particleHistNames[histType]); + std::string nameSimfull = Form("%s distribution for simulated particles", particleHistNames[histType]); if (doprocessRec || doprocessRecSim) { - pc.fParticleHistograms[histType][eRec][ba] = new TH1F( - nameRec.c_str(), nameRecfull.c_str(), nBinsPt, minPt, maxPt); - pc.fParticleHistograms[histType][eRec][ba]->GetXaxis()->SetTitle( - particleHistNames[histType]); - pc.fParticleHistogramsList->Add( - pc.fParticleHistograms[histType][eRec][ba]); + pc.fParticleHistograms[histType][eRec][ba] = new TH1F(nameRec.c_str(), nameRecfull.c_str(), nBinsPt, minPt, maxPt); + pc.fParticleHistograms[histType][eRec][ba]->GetXaxis()->SetTitle(particleHistNames[histType]); + pc.fParticleHistogramsList->Add(pc.fParticleHistograms[histType][eRec][ba]); } if (doprocessSim || doprocessRecSim) { - pc.fParticleHistograms[histType][eSim][ba] = new TH1F( - nameSim.c_str(), nameSimfull.c_str(), nBinsPt, minPt, maxPt); - pc.fParticleHistograms[histType][eSim][ba]->GetXaxis()->SetTitle( - particleHistNames[histType]); - pc.fParticleHistogramsList->Add( - pc.fParticleHistograms[histType][eSim][ba]); + pc.fParticleHistograms[histType][eSim][ba] = new TH1F(nameSim.c_str(), nameSimfull.c_str(), nBinsPt, minPt, maxPt); + pc.fParticleHistograms[histType][eSim][ba]->GetXaxis()->SetTitle(particleHistNames[histType]); + pc.fParticleHistogramsList->Add(pc.fParticleHistograms[histType][eSim][ba]); } } } template - void BookEventHistograms(T1 const &lEvBins, EventHistograms &ev) { - std::vector lCentBins = - lEvBins[histType]; // define local array and initialize it from an array - // set in the configurables + void bookEventHistograms(T1 const& lEvBins, EventHistograms &ev) { + std::vector lCentBins = lEvBins[histType]; // define local array and initialize it from an array set in the configurables int nBinsCent = static_cast(lCentBins[0]); float minCent = lCentBins[1]; float maxCent = lCentBins[2]; for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { - std::string nameRec = - Form("fHist%s[eRec][%s cut]", eventHistNames[histType], - cutBeforeAfterNames[ba]); - std::string nameSim = - Form("fHist%s[eSim][%s cut]", eventHistNames[histType], - cutBeforeAfterNames[ba]); + std::string nameRec = Form("fHist%s[eRec][%s cut]", eventHistNames[histType], cutBeforeAfterNames[ba]); + std::string nameSim = Form("fHist%s[eSim][%s cut]", eventHistNames[histType], cutBeforeAfterNames[ba]); std::string nameRecfull; std::string nameSimfull; if constexpr (histType == eCent) { - std::string nameRecfull = - Form("%s %s distribution for reconstructed events", - tc.fCentEstm.c_str(), eventHistNames[histType]); - std::string nameSimfull = - Form("%s %s distribution for simulated events", - tc.fCentEstm.c_str(), eventHistNames[histType]); + std::string nameRecfull = Form("%s %s distribution for reconstructed events", tc.fCentEstm.c_str(), eventHistNames[histType]); + std::string nameSimfull = Form("%s %s distribution for simulated events", tc.fCentEstm.c_str(), eventHistNames[histType]); } else if constexpr (histType == eMult) { - std::string nameRecfull = - Form("%s %s distribution for reconstructed events", - tc.fMultEstm.c_str(), eventHistNames[histType]); - std::string nameSimfull = - Form("%s %s distribution for simulated events", - tc.fMultEstm.c_str(), eventHistNames[histType]); + std::string nameRecfull = Form("%s %s distribution for reconstructed events", tc.fMultEstm.c_str(), eventHistNames[histType]); + std::string nameSimfull = Form("%s %s distribution for simulated events", tc.fMultEstm.c_str(), eventHistNames[histType]); } else { - std::string nameRecfull = - Form("%s distribution for reconstructed events", - eventHistNames[histType]); - std::string nameSimfull = Form("%s distribution for simulated events", - eventHistNames[histType]); + std::string nameRecfull = Form("%s distribution for reconstructed events", eventHistNames[histType]); + std::string nameSimfull = Form("%s distribution for simulated events", eventHistNames[histType]); } if (doprocessRec || doprocessRecSim) { - ev.fEventHistograms[histType][eRec][ba] = new TH1F( - nameRec.c_str(), nameRecfull.c_str(), nBinsCent, minCent, maxCent); - ev.fEventHistograms[histType][eRec][ba]->GetXaxis()->SetTitle( - eventHistNames[histType]); + ev.fEventHistograms[histType][eRec][ba] = new TH1F(nameRec.c_str(), nameRecfull.c_str(), nBinsCent, minCent, maxCent); + ev.fEventHistograms[histType][eRec][ba]->GetXaxis()->SetTitle(eventHistNames[histType]); ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eRec][ba]); } if (doprocessSim || doprocessRecSim) { if constexpr (histType != eNumContrib && histType != eMult) { - ev.fEventHistograms[histType][eSim][ba] = - new TH1F(nameSim.c_str(), nameSimfull.c_str(), nBinsCent, minCent, - maxCent); - ev.fEventHistograms[histType][eSim][ba]->GetXaxis()->SetTitle( - eventHistNames[histType]); + ev.fEventHistograms[histType][eSim][ba] = new TH1F(nameSim.c_str(), nameSimfull.c_str(), nBinsCent, minCent, maxCent); + ev.fEventHistograms[histType][eSim][ba]->GetXaxis()->SetTitle(eventHistNames[histType]); ev.fEventHistogramsList->Add(ev.fEventHistograms[histType][eSim][ba]); } // No nContrib and multiplicity for processSim } } } - template - void BookQAHistograms(T1 const &lQABins, QAHistograms &qa) { + template void bookQAHistograms(T1 const& lQABins, QAHistograms &qa) { int nBinsCentX = 0; float minCentX = 0.; float maxCentX = 0.; @@ -1305,66 +1132,49 @@ struct MultiparticleCumulants { // this name is used in lower-case format to } for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { - std::string name = Form("fHist%s[%s cut]", eventHistNames[histType], - cutBeforeAfterNames[ba]); + std::string name = Form("fHist%s[%s cut]", eventHistNames[histType], cutBeforeAfterNames[ba]); std::string namefull; if constexpr (histType == eCent) { - std::string namefull = - Form("Quality assurance of %s %s", tc.fCentEstm.c_str(), - eventHistNames[histType]); + std::string namefull = Form("Quality assurance of %s %s", tc.fCentEstm.c_str(), eventHistNames[histType]); } else if constexpr (histType == eMult) { - std::string namefull = - Form("Quality assurance of %s %s vs. NContributors", - tc.fMultEstm.c_str(), eventHistNames[histType]); + std::string namefull = Form("Quality assurance of %s %s vs. NContributors", tc.fMultEstm.c_str(), eventHistNames[histType]); } else { - std::string namefull = - Form("Quality assurance of %s", eventHistNames[histType]); + std::string namefull = Form("Quality assurance of %s", eventHistNames[histType]); } if constexpr (histType == eMult) { - qa.fQAHistograms[histType][ba] = - new TH2F(name.c_str(), namefull.c_str(), nBinsCentX, minCentX, - maxCentX, nBinsCentY, minCentY, maxCentY); + qa.fQAHistograms[histType][ba] = new TH2F(name.c_str(), namefull.c_str(), nBinsCentX, minCentX, maxCentX, nBinsCentY, minCentY, maxCentY); qa.fQAHistograms[histType][ba]->GetYaxis()->SetTitle("NContributors"); - qa.fQAHistograms[histType][ba]->GetXaxis()->SetTitle( - "Reference multiplicity"); + qa.fQAHistograms[histType][ba]->GetXaxis()->SetTitle("Reference multiplicity"); } else { - qa.fQAHistograms[histType][ba] = - new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, - maxCent, nBinsCent, minCent, maxCent); - qa.fQAHistograms[histType][ba]->GetYaxis()->SetTitle( - Form("Simulated %s", eventHistNames[histType])); - qa.fQAHistograms[histType][ba]->GetXaxis()->SetTitle( - Form("Reconstructed %s", eventHistNames[histType])); + qa.fQAHistograms[histType][ba] = new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, maxCent, nBinsCent, minCent, maxCent); + qa.fQAHistograms[histType][ba]->GetYaxis()->SetTitle(Form("Simulated %s", eventHistNames[histType])); + qa.fQAHistograms[histType][ba]->GetXaxis()->SetTitle(Form("Reconstructed %s", eventHistNames[histType])); } qa.fQAHistogramsList->Add(qa.fQAHistograms[histType][ba]); } } template - void BookCorrHistograms(T1 const &lCrBins, CorrHistograms &cr) { + void bookCorrHistograms(T1 const& lCrBins, CorrHistograms &cr) { - std::vector lCentBins = - lCrBins[histType]; // define local array and initialize it from an array - // set in the configurables + std::vector lCentBins = lCrBins[histType]; // define local array and initialize it from an array set in the configurables int nBinsCent = static_cast(lCentBins[0]); float minCent = lCentBins[1]; float maxCent = lCentBins[2]; for (int ba = 0; ba < eCutBeforeAfter_N; ba++) { std::string name = ""; - std::string namefull = - Form("%s correlation %s cut", corrHistNames[histType], - cutBeforeAfterNames[ba]); + std::string namefull = Form("%s correlation %s cut", corrHistNames[histType], cutBeforeAfterNames[ba]); - int N = 0; + int nEstm = 0; if constexpr (histType == eCorrCent) { - N = eCentEstm_N; + nEstm = eCentEstm_N; } else if constexpr (histType == eCorrMult) { - N = eMultEstm_N; + nEstm = eMultEstm_N; } - for (int i = 0; i < N; i++) { + for (int i = 0; i < nEstm; i++) { std::string titleX = ""; if (histType == eCorrCent) { @@ -1374,32 +1184,22 @@ struct MultiparticleCumulants { // this name is used in lower-case format to // x axis bin... } - for (int j = i + 1; j < N; j++) { + for (int j = i + 1; j < nEstm; j++) { std::string titleY = ""; if (histType == eCorrCent) { - name = Form("fHist%s[%s vs. %s][%s cut]", corrHistNames[histType], - centEstmNames[i], centEstmNames[j], - cutBeforeAfterNames[ba]); + name = Form("fHist%s[%s vs. %s][%s cut]", corrHistNames[histType], centEstmNames[i], centEstmNames[j], cutBeforeAfterNames[ba]); titleY = Form("%s %s", centEstmNames[j], corrHistNames[histType]); - cr.fCorrHistograms[histType][i][j][ba] = - new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, - maxCent, nBinsCent, minCent, maxCent); + cr.fCorrHistograms[histType][i][j][ba] = new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, maxCent, nBinsCent, minCent, maxCent); } else if (histType == eCorrMult) { - name = Form("fHist%s[%s vs. %s][%s cut]", corrHistNames[histType], - multEstmNames[i], multEstmNames[j], - cutBeforeAfterNames[ba]); + name = Form("fHist%s[%s vs. %s][%s cut]", corrHistNames[histType], multEstmNames[i], multEstmNames[j], cutBeforeAfterNames[ba]); titleY = Form("%s %s", multEstmNames[j], corrHistNames[histType]); // y axis bins... - cr.fCorrHistograms[histType][i][j][ba] = - new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, - maxCent, nBinsCent, minCent, maxCent); + cr.fCorrHistograms[histType][i][j][ba] = new TH2F(name.c_str(), namefull.c_str(), nBinsCent, minCent, maxCent, nBinsCent, minCent, maxCent); } - cr.fCorrHistograms[histType][i][j][ba]->GetYaxis()->SetTitle( - titleY.c_str()); - cr.fCorrHistograms[histType][i][j][ba]->GetXaxis()->SetTitle( - titleX.c_str()); + cr.fCorrHistograms[histType][i][j][ba]->GetYaxis()->SetTitle(titleY.c_str()); + cr.fCorrHistograms[histType][i][j][ba]->GetXaxis()->SetTitle(titleX.c_str()); cr.fCorrHistogramsList->Add(cr.fCorrHistograms[histType][i][j][ba]); } } @@ -1475,9 +1275,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to pc.fParticleHistogramsList = new TList(); pc.fParticleHistogramsList->SetName("ParticleHistograms"); pc.fParticleHistogramsList->SetOwner(kFALSE); - fBaseList->Add(pc.fParticleHistogramsList); // any nested TList in the base - // TList appears as a subdir in - // the output ROOT file + fBaseList->Add(pc.fParticleHistogramsList); // any nested TList in the base TList appears as a subdir in the output ROOT file ev.fEventHistogramsList = new TList(); ev.fEventHistogramsList->SetName("EventHistograms"); @@ -1500,36 +1298,32 @@ struct MultiparticleCumulants { // this name is used in lower-case format to fBaseList->Add(cr.fCorrHistogramsList); mc.fMultiparticleCorrelationProfilesList = new TList(); - mc.fMultiparticleCorrelationProfilesList->SetName( - "MultiparticleCorrelationProfiles"); + mc.fMultiparticleCorrelationProfilesList->SetName("MultiparticleCorrelationProfiles"); mc.fMultiparticleCorrelationProfilesList->SetOwner(kFALSE); fBaseList->Add(mc.fMultiparticleCorrelationProfilesList); std::vector> lPcBins = {tc.fPtBins, tc.fPhiBins}; - std::vector> lEvBins = { - tc.fCentBins, tc.fMultBins, tc.fVerXBins, - tc.fVerYBins, tc.fVerZBins, tc.fNumContribBins}; - std::vector> lQABins = {tc.fCentBins, tc.fMultBins, - tc.fNumContribBins}; + std::vector> lEvBins = {tc.fCentBins, tc.fMultBins, tc.fVerXBins, tc.fVerYBins, tc.fVerZBins, tc.fNumContribBins}; + std::vector> lQABins = {tc.fCentBins, tc.fMultBins, tc.fNumContribBins}; std::vector> lCrBins = {tc.fCentBins, tc.fMultBins}; - BookParticleHistograms(lPcBins, pc); - BookParticleHistograms(lPcBins, pc); - BookEventHistograms(lEvBins, ev); - BookEventHistograms(lEvBins, ev); - BookEventHistograms(lEvBins, ev); - BookEventHistograms(lEvBins, ev); - BookEventHistograms(lEvBins, ev); - BookEventHistograms(lEvBins, ev); - BookQAHistograms(lQABins, qa); - BookQAHistograms(lQABins, qa); - BookCorrHistograms(lCrBins, cr); // if switch on ... - BookCorrHistograms(lCrBins, cr); + bookParticleHistograms(lPcBins, pc); + bookParticleHistograms(lPcBins, pc); + bookEventHistograms(lEvBins, ev); + bookEventHistograms(lEvBins, ev); + bookEventHistograms(lEvBins, ev); + bookEventHistograms(lEvBins, ev); + bookEventHistograms(lEvBins, ev); + bookEventHistograms(lEvBins, ev); + bookQAHistograms(lQABins, qa); + bookQAHistograms(lQABins, qa); + bookCorrHistograms(lCrBins, cr); // if switch on ... + bookCorrHistograms(lCrBins, cr); if (wt.fWeightSwitch) { wt.fWeightHistograms = getHistogramsWithWeights( tc.fFileWithWeights.c_str(), tc.fRunNumber.c_str()); - for (auto *hist : wt.fWeightHistograms) { + for (auto* const& hist : wt.fWeightHistograms) { //o2-linter: disable=const-ref-in-for-loop wt.fWeightHistogramsList->Add(hist); } } @@ -1540,66 +1334,46 @@ struct MultiparticleCumulants { // this name is used in lower-case format to new TProfile("prof2After", "2-p correlation after cut", 3, 2., 5.); mc.fTwoParticleCorrelationProfiles[eBefore]->Sumw2(); mc.fTwoParticleCorrelationProfiles[eAfter]->Sumw2(); - mc.fMultiparticleCorrelationProfilesList->Add( - mc.fTwoParticleCorrelationProfiles[eBefore]); - mc.fMultiparticleCorrelationProfilesList->Add( - mc.fTwoParticleCorrelationProfiles[eAfter]); - - mc.fFourParticleCorrelationProfiles[eBefore] = - new TProfile("prof4Before", "4-p correlation before cut", 2, 3., 5.); - mc.fFourParticleCorrelationProfiles[eAfter] = - new TProfile("prof4After", "4-p correlation after cut", 2, 3., 5.); + mc.fMultiparticleCorrelationProfilesList->Add(mc.fTwoParticleCorrelationProfiles[eBefore]); + mc.fMultiparticleCorrelationProfilesList->Add(mc.fTwoParticleCorrelationProfiles[eAfter]); + + mc.fFourParticleCorrelationProfiles[eBefore] = new TProfile("prof4Before", "4-p correlation before cut", 2, 3., 5.); + mc.fFourParticleCorrelationProfiles[eAfter] = new TProfile("prof4After", "4-p correlation after cut", 2, 3., 5.); mc.fFourParticleCorrelationProfiles[eBefore]->Sumw2(); mc.fFourParticleCorrelationProfiles[eAfter]->Sumw2(); - mc.fMultiparticleCorrelationProfilesList->Add( - mc.fFourParticleCorrelationProfiles[eBefore]); - mc.fMultiparticleCorrelationProfilesList->Add( - mc.fFourParticleCorrelationProfiles[eAfter]); + mc.fMultiparticleCorrelationProfilesList->Add(mc.fFourParticleCorrelationProfiles[eBefore]); + mc.fMultiparticleCorrelationProfilesList->Add(mc.fFourParticleCorrelationProfiles[eAfter]); } // end of void init(InitContext&) { // A) Process only reconstructed data: - void processRec(CollisionRec const &collision, aod::BCs const &, - TracksRec const &tracks) { + void processRec(CollisionRec const& collision, aod::BCs const&, TracksRec const& tracks) { // ... - // *) Steer all analysis steps: - Steer(collision, tracks); + runLoop(collision, tracks); } - PROCESS_SWITCH(MultiparticleCumulants, processRec, - "process only reconstructed data", - true); // yes, keep always one process switch "true", so that - // there is default running version + PROCESS_SWITCH(MultiparticleCumulants, processRec, "process only reconstructed data", true); // yes, keep always one process switch "true", so that there is default running version // ------------------------------------------- // B) Process both reconstructed and corresponding MC truth simulated data: - void processRecSim(CollisionRecSim const &collision, aod::BCs const &, - TracksRecSim const &tracks, aod::McParticles const &, - aod::McCollisions const &) { - Steer(collision, tracks); + void processRecSim(CollisionRecSim const& collision, aod::BCs const&, TracksRecSim const& tracks, aod::McParticles const&, aod::McCollisions const&) { + runLoop(collision, tracks); } - PROCESS_SWITCH( - MultiparticleCumulants, processRecSim, - "process both reconstructed and corresponding MC truth simulated data", - false); + PROCESS_SWITCH(MultiparticleCumulants, processRecSim, "process both reconstructed and corresponding MC truth simulated data", false); // ------------------------------------------- // C) Process only simulated data: - void processSim(CollisionSim const & /*collision*/, aod::BCs const &, - TracksSim const & /*tracks*/) { - // Steer(collision, tracks); // TBI 20241105 not ready yet, but I do - // not really need this one urgently, since RecSim is working, and I need - // that one for efficiencies... + void processSim(CollisionSim const& /*collision*/, aod::BCs const&, TracksSim const& /*tracks*/) { + // runLoop(collision, tracks); // TBI 20241105 not ready yet, but I do not really need this one urgently, since RecSim is working, and I need that one for efficiencies... } - PROCESS_SWITCH(MultiparticleCumulants, processSim, - "process only simulated data", false); + PROCESS_SWITCH(MultiparticleCumulants, processSim, "process only simulated data", false); }; // struct MultiparticleCumulants { // *) The final touch: -WorkflowSpec defineDataProcessing(ConfigContext const &cfgc) { +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{ adaptAnalysisTask(cfgc), }; From f63f16ee9f886d9465b5d57479f320b3ce6c8ea4 Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Thu, 11 Jun 2026 18:25:19 +0200 Subject: [PATCH 3/9] fix clang format --- .../Tasks/multiparticleCumulants.cxx | 215 +++++++++--------- 1 file changed, 110 insertions(+), 105 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index 6fdfea1851e..a9b0bf57637 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -1,6 +1,6 @@ // Copyright 2019-2020 CERN and copyright holders of ALICE O2. -// See https://alice-o2.web.cern.ch/copyright for details of the copyright -// holders. All rights not expressly granted are reserved. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. // // This software is distributed under the terms of the GNU General Public // License v3 (GPL Version 3), copied verbatim in the file "COPYING". @@ -16,7 +16,7 @@ #include "Common/DataModel/Centrality.h" #include "Common/DataModel/EventSelection.h" #include "Common/DataModel/Multiplicity.h" -#include "Common/DataModel/TrackSelectionTables.h" // needed for aod::TracksDCA table +#include "Common/DataModel/TrackSelectionTables.h" #include #include @@ -52,8 +52,8 @@ using namespace o2; using namespace o2::framework; // Definitions of join tables for Run 3 analysis: -using EventSelection = soa::Join; // no FV0Ms for cent -using CollisionRec = soa::Join::iterator; // use in json "isMC": "true" for "event-selection-task" +using EventSelection = soa::Join; +using CollisionRec = soa::Join::iterator; using CollisionRecSim = soa::Join::iterator; using CollisionSim = aod::McCollision; using TracksRec = soa::Join; @@ -68,12 +68,12 @@ using namespace std; // *) Define enums: enum EnRlMc { eRl = 0, - eMc + eMc }; enum EnRecSim { eRec = 0, - eSim, + eSim, eRecAndSim }; @@ -94,14 +94,13 @@ enum EnEventHistograms { eEventHistograms_N }; -const char *eventHistNames[eEventHistograms_N] = { +const char* eventHistNames[eEventHistograms_N] = { "Centrality", "Multiplicity", "VertexX", "VertexY", "VertexZ", - "NumContrib" -}; + "NumContrib"}; enum EnParticleHistograms { ePt, @@ -109,10 +108,9 @@ enum EnParticleHistograms { eParticleHistograms_N }; -const char *particleHistNames[eParticleHistograms_N] = { +const char* particleHistNames[eParticleHistograms_N] = { "Pt", - "Phi" -}; + "Phi"}; enum EnQAHistograms { eQACent, @@ -126,7 +124,7 @@ enum EnCorrHistograms { eCorrHistograms_N }; -const char *corrHistNames[eCorrHistograms_N] = { +const char* corrHistNames[eCorrHistograms_N] = { "Centrality", "Multiplicity", }; @@ -138,11 +136,10 @@ enum EnCentEstm { eCentEstm_N }; -const char *centEstmNames[eCentEstm_N] = { +const char* centEstmNames[eCentEstm_N] = { "FT0C", "FT0M", - "FV0A", -}; + "FV0A",}; enum EnMultEstm { eMultFT0C, @@ -151,11 +148,10 @@ enum EnMultEstm { eMultEstm_N }; -const char *multEstmNames[eMultEstm_N] = { +const char* multEstmNames[eMultEstm_N] = { "FT0C", "FT0M", - "FV0A" -}; + "FV0A"}; enum EnCutBeforeAfter { eBefore, @@ -163,7 +159,7 @@ enum EnCutBeforeAfter { eCutBeforeAfter_N }; -const char *cutBeforeAfterNames[eCutBeforeAfter_N] = { +const char* cutBeforeAfterNames[eCutBeforeAfter_N] = { "before", "after", }; @@ -237,8 +233,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // **) Task configuration: struct TaskConfiguration { bool fProcess[eProcess_N] = {false}; // Set what to process. See enum EnProcess for full description. Set via implicit variables within a PROCESS_SWITCH clause. - bool fDryRun = false; // book all histos and run without filling and calculating anything - + bool fDryRun = false; // book all histos and run without filling and calculating anything + std::string fCentEstm = "FT0M"; std::string fMultEstm = "FV0A"; @@ -287,30 +283,30 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } tc; struct ParticleHistograms { - TList *fParticleHistogramsList = NULL; //! fWeightHistograms; + TList* fWeightHistogramsList = NULL; + std::vector fWeightHistograms; } wt; struct EventByEventQuantities { @@ -321,13 +317,13 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam float fNumContrib = 0.; float fTwoParticleCorrelationEbye[eCutBeforeAfter_N][3] = {{0., 0., 0.}, {0., 0., 0.}}; //<2> [before, after][v2^2, v3^2, v4^2] - float fFourParticleCorrelationEbye[eCutBeforeAfter_N][2] = {{0., 0.}, {0., 0.}}; //<4> [before, after][v3^2 v2^2, v4^2 v2^2] + float fFourParticleCorrelationEbye[eCutBeforeAfter_N][2] = {{0., 0.}, {0., 0.}}; //<4> [before, after][v3^2 v2^2, v4^2 v2^2] } ebye; struct MultiparticleCorrelationProfile { - TList *fMultiparticleCorrelationProfilesList = NULL; - TProfile *fTwoParticleCorrelationProfiles[eCutBeforeAfter_N] = {NULL}; - TProfile *fFourParticleCorrelationProfiles[eCutBeforeAfter_N] = {NULL}; + TList* fMultiparticleCorrelationProfilesList = NULL; + TProfile* fTwoParticleCorrelationProfiles[eCutBeforeAfter_N] = {NULL}; + TProfile* fFourParticleCorrelationProfiles[eCutBeforeAfter_N] = {NULL}; } mc; struct MultiparticleCorrelationCalculation { @@ -336,12 +332,14 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam static constexpr int MaxCorrelator = 4; //<> static constexpr int MaxHarmonic = 5; // n+1 in vn, in this case n=4 as we need v2, v3, v4 static constexpr int MaxPower = MaxCorrelator + 1; - static constexpr int NumSC = 2; //need SC(3,2) and SC(4,2) + static constexpr int NumSC = 2; // need SC(3,2) and SC(4,2) TComplex fQvectorBefore[MaxHarmonic][MaxPower]; TComplex fQvectorAfter[MaxHarmonic][MaxPower]; // All needed Q-vector components } mcc; - template bool ctEventCuts(T1 const& collision, T2 rlCollisionCentAll, T3 rlCollisionMultAll) { + template + bool ctEventCuts(T1 const& collision, T2 rlCollisionCentAll, T3 rlCollisionMultAll) + { bool pass = true; bool bVertexZCut = true; @@ -417,7 +415,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam return pass; } - template bool ctParticleCuts(T1 const& track) { + template + bool ctParticleCuts(T1 const& track) + { bool pass = true; bool bPtCut = true; @@ -446,7 +446,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // *) For mc event only if constexpr (rm == eMc) { - TParticlePDG *particle = pdg->GetParticle(track.pdgCode()); + TParticlePDG* particle = pdg->GetParticle(track.pdgCode()); if (!particle) { // LOGF(warning, "PDG code %d not found", track.pdgCode()); bSignCut = false; @@ -481,7 +481,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam return pass; } - TComplex mccQ(int n, int p, EnCutBeforeAfter eba) { + TComplex mccQ(int n, int p, EnCutBeforeAfter eba) + { // Using the fact that Q{-n,p} = Q{n,p}^*. if (eba == eBefore) { @@ -495,10 +496,10 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } return TComplex::Conjugate(mcc.fQvectorAfter[-n][p]); } - } - TComplex mccRecursion(int n, int *harmonic, EnCutBeforeAfter eba, int mult = 1, int skip = 0) { + TComplex mccRecursion(int n, int* harmonic, EnCutBeforeAfter eba, int mult = 1, int skip = 0) + { // Calculate multi-particle correlators by using recursion (an improved faster version) originally developed by Kristjan Gulbrandsen (gulbrand@nbi.dk). int nm1 = n - 1; @@ -532,18 +533,19 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (mult == 1) return c - c2; - return c - double(mult) * c2; + return c - double(mult)* c2; } - TObject *getObjectFromList(TList *list, const char *objectName) { + TObject* getObjectFromList(TList* list, const char* objectName) + { // Get TObject pointer from TList, even if it's in some nested TList. // Foreseen to be used to fetch histograms or profiles from files directly. // Some ideas taken from TCollection::ls() // If you have added histograms directly to files (without TList's), then // you can fetch them directly with file->Get("hist-name"). - // Usage: TH1D *hist = (TH1D*) + // Usage: TH1D* hist = (TH1D*) // getObjectFromList("some-valid-TList-pointer","some-object-name"); // Insanity checks: @@ -558,19 +560,18 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // The object is in the current base list: - TObject *objectFinal = - list->FindObject(objectName); // the final object I am after + TObject* objectFinal = list->FindObject(objectName); // the final object I am after if (objectFinal) { return objectFinal; } // Otherwise, search for the object recursively in the nested lists: - TObject *objectIter; // iterator object in the loop below + TObject* objectIter; // iterator object in the loop below TIter next(list); while ((objectIter = next())) // double round braces are to silence the warnings { if (TString(objectIter->ClassName()).EqualTo("TList")) { - objectFinal = getObjectFromList(reinterpret_cast(objectIter), objectName); + objectFinal = getObjectFromList(reinterpret_cast(objectIter), objectName); if (objectFinal) return objectFinal; } @@ -578,18 +579,19 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam return NULL; - } // TObject* getObjectFromList(TList *list, char *objectName) + } // TObject* getObjectFromList(TList* list, char* objectName) - std::vector getHistogramsWithWeights(const char *filePath, const char *runNumber) { + std::vector getHistogramsWithWeights(const char* filePath, const char* runNumber) + { // a) Return value: - std::vector histograms; - TList *baseList = NULL; // base top-level list in the TFile, e.g. named "ccdb_object" - TList *listWithRuns = NULL; // nested list with run-wise TList's holding run-specific weights + std::vector histograms; + TList* baseList = NULL; // base top-level list in the TFile, e.g. named "ccdb_object" + TList* listWithRuns = NULL; // nested list with run-wise TList's holding run-specific weights // c) Determine from filePath if the file in on a local machine, or in home dir AliEn, or in CCDB: - // Algorithm: + // Algorithm: // If filePath begins with "/alice/data/CCDB/" then it's in home dir AliEn. - // If filePath begins with "/alice-ccdb.cern.ch/" then it's in CCDB. + // If filePath begins with "/alice-ccdb.cern.ch/" then it's in CCDB. // Therefore, files in AliEn and CCDB must be specified with abs path, for local files both abs and relative paths are just fine. bool bFileIsInAliEn = false; bool bFileIsInCCDB = false; @@ -604,11 +606,11 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (bFileIsInAliEn) { // File you want to access is in your home dir in AliEn: - TGrid *alien = TGrid::Connect("alien", gSystem->Getenv("USER"), "", ""); // do not forget to add #include to the preamble of your analysis task + TGrid* alien = TGrid::Connect("alien", gSystem->Getenv("USER"), "", ""); // do not forget to add #include to the preamble of your analysis task if (!alien) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - TFile *weightsFile = TFile::Open(Form("alien://%s", filePath), "READ"); // yes, ROOT can open a file transparently, even if it's sitting in AliEn, with this specific syntax + TFile* weightsFile = TFile::Open(Form("alien://%s", filePath), "READ"); // yes, ROOT can open a file transparently, even if it's sitting in AliEn, with this specific syntax if (!weightsFile) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -620,11 +622,11 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // Finally, from the top-level TList, get the desired nested TList => the technical problem here is that it can be nested at any level, // for thare there is a helper utility function GetObjectFromList(...) , see its implementation further below - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -634,18 +636,17 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // Remember that here I do not access the file; instead, I directly access the object in that file. // My home dir in CCDB: https://alice-ccdb.cern.ch/browse/Users/a/abilandz/ => adapt for your case ccdb->setURL("https://alice-ccdb.cern.ch"); // to be able to use "ccdb" this object in your analysis task, see 4b/ below - baseList = reinterpret_cast(ccdb->get( - TString(filePath).ReplaceAll("/alice-ccdb.cern.ch/", "").Data())); + baseList = reinterpret_cast(ccdb->get(TString(filePath).ReplaceAll("/alice-ccdb.cern.ch/", "").Data())); baseList->ls(); if (!baseList) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -664,7 +665,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - TFile *weightsFile = TFile::Open(filePath, "READ"); + TFile* weightsFile = TFile::Open(filePath, "READ"); if (!weightsFile) { LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } @@ -676,11 +677,11 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam LOGF(fatal, "\033[1;31m%s at line %d\033[0m", __FUNCTION__, __LINE__); } - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumber)); if (!listWithRuns) { TString runNumberWithLeadingZeroes = "000"; runNumberWithLeadingZeroes += runNumber; // another try, with "000" prepended to run number - listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); + listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { baseList->ls(); LOGF(fatal, @@ -690,7 +691,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } TIter next(listWithRuns); - TObject *object = nullptr; + TObject* object = nullptr; while (true) { @@ -699,13 +700,13 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam break; } - auto *hist = dynamic_cast(object); + auto* hist = dynamic_cast(object); if (!hist) { continue; } hist->SetDirectory(0); - auto *histClone = dynamic_cast(hist->Clone()); + auto* histClone = dynamic_cast(hist->Clone()); if (!histClone) { LOGF(fatal, "Failed to clone histogram %s", hist->GetName()); } @@ -716,16 +717,17 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam return histograms; } - // *) Define all member functions to be called in the main process* functions: + //* ) Define all member functions to be called in the main process* functions: template - void runLoop(T1 const& collision, T2 const& tracks) { + void runLoop(T1 const& collision, T2 const& tracks) + { // Dry run: if (tc.fDryRun) { return; } - TH1F *histAcceptanceWeight = wt.fWeightHistograms[1]; + TH1F* histAcceptanceWeight = wt.fWeightHistograms[1]; for (int h = 0; h < mcc.MaxHarmonic; h++) { for (int p = 0; p < mcc.MaxPower; p++) { @@ -740,11 +742,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } float rlCollisionCentAll[eCentEstm_N] = { - // collision.centFT0A(), collision.centFT0C(), collision.centFT0M(), - collision.centFV0A() - }; + collision.centFV0A()}; float rlCollisionCent = 0.; @@ -761,8 +761,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam float rlCollisionMultAll[eMultEstm_N] = { static_cast(collision.multFT0C()), static_cast(collision.multFT0M()), - static_cast(collision.multFV0A()) - }; + static_cast(collision.multFV0A())}; float rlCollisionMult = 0.; for (int i = 0; i < eMultEstm_N; i++) { @@ -807,7 +806,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam for (int i = 0; i < eCentEstm_N; i++) { for (int j = i + 1; j < eCentEstm_N; j++) { - auto *h = cr.fCorrHistograms[eCorrCent][i][j][eBefore]; + auto* h = cr.fCorrHistograms[eCorrCent][i][j][eBefore]; if (!h) { LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrCent][%d][%d][eBefore]", i, j); } @@ -817,7 +816,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam for (int i = 0; i < eMultEstm_N; i++) { for (int j = i + 1; j < eMultEstm_N; j++) { - auto *h = cr.fCorrHistograms[eCorrMult][i][j][eBefore]; + auto* h = cr.fCorrHistograms[eCorrMult][i][j][eBefore]; if (!h) { LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrMult][%d][%d][eBefore]", i, j); } @@ -838,7 +837,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam float mcCollisionCent = 0.; float b = mccollision.impactParameter() * std::pow(10, -15); // convert fm to m - float xs = 7.71 * std::pow(10, -28); // convert barn to m^2 + float xs = 7.71 * std::pow(10, -28); // convert barn to m^2 mcCollisionCent = o2::constants::math::PI * b * b / xs * 100; ebye.fCentralitySim = mcCollisionCent; @@ -883,7 +882,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (tc.fCentCorrCutSwitch) { for (int i = 0; i < eCentEstm_N; i++) { for (int j = i + 1; j < eCentEstm_N; j++) { - auto *h = cr.fCorrHistograms[eCorrCent][i][j][eAfter]; + auto* h = cr.fCorrHistograms[eCorrCent][i][j][eAfter]; if (!h) { LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrCent][%d][%d][eAfter]", i, j); } @@ -896,7 +895,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam for (int i = 0; i < eMultEstm_N; i++) { for (int j = i + 1; j < eMultEstm_N; j++) { - auto *h = cr.fCorrHistograms[eCorrMult][i][j][eAfter]; + auto* h = cr.fCorrHistograms[eCorrMult][i][j][eAfter]; if (!h) { LOGF(fatal, "Missing histogram cr.fCorrHistograms[eCorrMult][%d][%d][eAfter]", i, j); } @@ -931,7 +930,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam dPhi = track.phi(); if (wt.fWeightSwitch) { - auto *hist = wt.fWeightHistograms[1]; + auto* hist = wt.fWeightHistograms[1]; wPhi = hist->GetBinContent(wt.fWeightHistograms[1]->GetXaxis()->FindBin(dPhi)); } @@ -986,7 +985,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // end of for (int64_t i = 0; i < tracks.size(); i++) { - for (int i = 0; i < mcc.MaxHarmonic-2 ; i++) { + for (int i = 0; i < mcc.MaxHarmonic - 2 ; i++) { mcc.h1 = -(i + 2); mcc.h2 = i + 2; int harmonicsTwoNum[2] = {mcc.h1, mcc.h2}; @@ -1036,7 +1035,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } template - void bookParticleHistograms(T1 const& lPcBins, ParticleHistograms &pc) { + void bookParticleHistograms(T1 const& lPcBins, ParticleHistograms& pc) + { std::vector lPtBins = lPcBins[histType]; // define local array and initialize it from an array set in the configurables int nBinsPt = static_cast(lPtBins[0]); float minPt = lPtBins[1]; @@ -1064,7 +1064,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } template - void bookEventHistograms(T1 const& lEvBins, EventHistograms &ev) { + void bookEventHistograms(T1 const& lEvBins, EventHistograms& ev) + { std::vector lCentBins = lEvBins[histType]; // define local array and initialize it from an array set in the configurables int nBinsCent = static_cast(lCentBins[0]); float minCent = lCentBins[1]; @@ -1104,7 +1105,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } } - template void bookQAHistograms(T1 const& lQABins, QAHistograms &qa) { + template + void bookQAHistograms(T1 const& lQABins, QAHistograms& qa) + { int nBinsCentX = 0; float minCentX = 0.; float maxCentX = 0.; @@ -1157,7 +1160,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } template - void bookCorrHistograms(T1 const& lCrBins, CorrHistograms &cr) { + void bookCorrHistograms(T1 const& lCrBins, CorrHistograms& cr) + { std::vector lCentBins = lCrBins[histType]; // define local array and initialize it from an array set in the configurables int nBinsCent = static_cast(lCentBins[0]); @@ -1207,7 +1211,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // *) Initialize and book all objects: - void init(InitContext &) { + void init(InitContext&) + { // ... code to book and initialize all analysis objects ... @@ -1267,7 +1272,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam wt.fWeightSwitch = cfWeightSwitch; // *) Book base list: - TList *temp = new TList(); + TList* temp = new TList(); temp->SetOwner(kTRUE); fBaseList.setObject(temp); @@ -1323,15 +1328,13 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (wt.fWeightSwitch) { wt.fWeightHistograms = getHistogramsWithWeights( tc.fFileWithWeights.c_str(), tc.fRunNumber.c_str()); - for (auto* const& hist : wt.fWeightHistograms) { //o2-linter: disable=const-ref-in-for-loop + for (auto* const& hist : wt.fWeightHistograms) { // o2-linter: disable=const-ref-in-for-loop wt.fWeightHistogramsList->Add(hist); } } - mc.fTwoParticleCorrelationProfiles[eBefore] = - new TProfile("prof2Before", "2-p correlation before cut", 3, 2., 5.); - mc.fTwoParticleCorrelationProfiles[eAfter] = - new TProfile("prof2After", "2-p correlation after cut", 3, 2., 5.); + mc.fTwoParticleCorrelationProfiles[eBefore] = new TProfile("prof2Before", "2-p correlation before cut", 3, 2., 5.); + mc.fTwoParticleCorrelationProfiles[eAfter] = new TProfile("prof2After", "2-p correlation after cut", 3, 2., 5.); mc.fTwoParticleCorrelationProfiles[eBefore]->Sumw2(); mc.fTwoParticleCorrelationProfiles[eAfter]->Sumw2(); mc.fMultiparticleCorrelationProfilesList->Add(mc.fTwoParticleCorrelationProfiles[eBefore]); @@ -1347,7 +1350,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // end of void init(InitContext&) { // A) Process only reconstructed data: - void processRec(CollisionRec const& collision, aod::BCs const&, TracksRec const& tracks) { + void processRec(CollisionRec const& collision, aod::BCs const&, TracksRec const& tracks) + { // ... // *) Steer all analysis steps: runLoop(collision, tracks); @@ -1357,7 +1361,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // ------------------------------------------- // B) Process both reconstructed and corresponding MC truth simulated data: - void processRecSim(CollisionRecSim const& collision, aod::BCs const&, TracksRecSim const& tracks, aod::McParticles const&, aod::McCollisions const&) { + void processRecSim(CollisionRecSim const& collision, aod::BCs const&, TracksRecSim const& tracks, aod::McParticles const&, aod::McCollisions const&) + { runLoop(collision, tracks); } PROCESS_SWITCH(MultiparticleCumulants, processRecSim, "process both reconstructed and corresponding MC truth simulated data", false); @@ -1365,7 +1370,8 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // ------------------------------------------- // C) Process only simulated data: - void processSim(CollisionSim const& /*collision*/, aod::BCs const&, TracksSim const& /*tracks*/) { + void processSim(CollisionSim const& /*collision*/, aod::BCs const&, TracksSim const& /*tracks*/) + { // runLoop(collision, tracks); // TBI 20241105 not ready yet, but I do not really need this one urgently, since RecSim is working, and I need that one for efficiencies... } PROCESS_SWITCH(MultiparticleCumulants, processSim, "process only simulated data", false); @@ -1373,8 +1379,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam }; // struct MultiparticleCumulants { // *) The final touch: -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { - return WorkflowSpec{ - adaptAnalysisTask(cfgc), - }; +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{adaptAnalysisTask(cfgc),}; } // WorkflowSpec... From 4cf96bf1c74ed639bc166ad105f94b3e1a6222f6 Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Thu, 11 Jun 2026 18:52:33 +0200 Subject: [PATCH 4/9] try to fix clang format --- .../Tasks/multiparticleCumulants.cxx | 24 +++++++++---------- 1 file changed, 11 insertions(+), 13 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index a9b0bf57637..f6275ddd5f1 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -139,7 +139,7 @@ enum EnCentEstm { const char* centEstmNames[eCentEstm_N] = { "FT0C", "FT0M", - "FV0A",}; + "FV0A"}; enum EnMultEstm { eMultFT0C, @@ -283,19 +283,19 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } tc; struct ParticleHistograms { - TList* fParticleHistogramsList = NULL; //!> + static constexpr int MaxCorrelator = 4; // <> static constexpr int MaxHarmonic = 5; // n+1 in vn, in this case n=4 as we need v2, v3, v4 static constexpr int MaxPower = MaxCorrelator + 1; static constexpr int NumSC = 2; // need SC(3,2) and SC(4,2) @@ -533,7 +533,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (mult == 1) return c - c2; - return c - double(mult)* c2; + return c - double(mult) * c2; } @@ -585,7 +585,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam { // a) Return value: std::vector histograms; - TList* baseList = NULL; // base top-level list in the TFile, e.g. named "ccdb_object" + TList* baseList = NULL; // base top-level list in the TFile, e.g. named "ccdb_object" TList* listWithRuns = NULL; // nested list with run-wise TList's holding run-specific weights // c) Determine from filePath if the file in on a local machine, or in home dir AliEn, or in CCDB: @@ -985,7 +985,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } // end of for (int64_t i = 0; i < tracks.size(); i++) { - for (int i = 0; i < mcc.MaxHarmonic - 2 ; i++) { + for (int i = 0; i < mcc.MaxHarmonic - 2; i++) { mcc.h1 = -(i + 2); mcc.h2 = i + 2; int harmonicsTwoNum[2] = {mcc.h1, mcc.h2}; @@ -1031,7 +1031,6 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam mc.fFourParticleCorrelationProfiles[eAfter]->Fill(i + 3.5, ebye.fFourParticleCorrelationEbye[eAfter][i], wFourRecursionAfter); } } - } template @@ -1326,8 +1325,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam bookCorrHistograms(lCrBins, cr); if (wt.fWeightSwitch) { - wt.fWeightHistograms = getHistogramsWithWeights( - tc.fFileWithWeights.c_str(), tc.fRunNumber.c_str()); + wt.fWeightHistograms = getHistogramsWithWeights(tc.fFileWithWeights.c_str(), tc.fRunNumber.c_str()); for (auto* const& hist : wt.fWeightHistograms) { // o2-linter: disable=const-ref-in-for-loop wt.fWeightHistogramsList->Add(hist); } @@ -1381,5 +1379,5 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam // *) The final touch: WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { - return WorkflowSpec{adaptAnalysisTask(cfgc),}; + return WorkflowSpec{adaptAnalysisTask(cfgc)}; } // WorkflowSpec... From 266f5056e361f9be336274a5aff9328abc44e824 Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Thu, 11 Jun 2026 21:15:28 +0200 Subject: [PATCH 5/9] try to fix clang format --- .../Tasks/multiparticleCumulants.cxx | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index f6275ddd5f1..baf061d5e83 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -533,8 +533,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (mult == 1) return c - c2; - return c - double(mult) * c2; - + return c - static_cast(mult) * c2; } TObject* getObjectFromList(TList* list, const char* objectName) @@ -684,8 +683,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam listWithRuns = reinterpret_cast(getObjectFromList(baseList, runNumberWithLeadingZeroes.Data())); if (!listWithRuns) { baseList->ls(); - LOGF(fatal, - "\033[1;31m%s at line %d : this crash can happen if in the output file there is no list with weights for the current runnumber = %s\033[0m", __FUNCTION__, __LINE__, runNumber); + LOGF(fatal, "\033[1;31m%s at line %d : this crash can happen if in the output file there is no list with weights for the current runnumber = %s\033[0m", __FUNCTION__, __LINE__, runNumber); } } } @@ -904,8 +902,9 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam } } - } else + } else { return; + } } int nTracksBefore = tracks.size(); From f72d4d78900f34ae8cb22f4bf014c4cfbaaf0323 Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Fri, 12 Jun 2026 12:42:31 +0200 Subject: [PATCH 6/9] fix error --- .../Tasks/multiparticleCumulants.cxx | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index baf061d5e83..30f47f1976e 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -489,12 +489,16 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (n >= 0) { return mcc.fQvectorBefore[n][p]; } - return TComplex::Conjugate(mcc.fQvectorBefore[-n][p]); - } else if (eba == eAfter) { + else { + return TComplex::Conjugate(mcc.fQvectorBefore[-n][p]); + } + } else { if (n >= 0) { return mcc.fQvectorAfter[n][p]; } - return TComplex::Conjugate(mcc.fQvectorAfter[-n][p]); + else { + return TComplex::Conjugate(mcc.fQvectorAfter[-n][p]); + } } } @@ -725,7 +729,7 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam return; } - TH1F* histAcceptanceWeight = wt.fWeightHistograms[1]; + // TH1F* histAcceptanceWeight = wt.fWeightHistograms[1]; for (int h = 0; h < mcc.MaxHarmonic; h++) { for (int p = 0; p < mcc.MaxPower; p++) { From 8b513bd06732f933879c022c3f46835562cabfad Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Fri, 12 Jun 2026 12:58:18 +0200 Subject: [PATCH 7/9] fix format --- .../MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index 30f47f1976e..aa66af40ccb 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -1383,4 +1383,4 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{adaptAnalysisTask(cfgc)}; -} // WorkflowSpec... +} // WorkflowSpec... \ No newline at end of file From 73b9d468b7ff4cd28e9a96eb236354ec5612f655 Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Fri, 12 Jun 2026 13:03:21 +0200 Subject: [PATCH 8/9] fix format --- .../Tasks/multiparticleCumulants.cxx | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index aa66af40ccb..db4b3e0d187 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -488,15 +488,13 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam if (eba == eBefore) { if (n >= 0) { return mcc.fQvectorBefore[n][p]; - } - else { + } else { return TComplex::Conjugate(mcc.fQvectorBefore[-n][p]); } } else { if (n >= 0) { return mcc.fQvectorAfter[n][p]; - } - else { + } else { return TComplex::Conjugate(mcc.fQvectorAfter[-n][p]); } } From 32965e59317451b8bc9dcff8f93dda455d27e86a Mon Sep 17 00:00:00 2001 From: Pei-Ying Kuan Date: Fri, 12 Jun 2026 13:09:56 +0200 Subject: [PATCH 9/9] add new line --- .../MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx index db4b3e0d187..4abfd65310a 100644 --- a/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx +++ b/PWGCF/MultiparticleCorrelations/Tasks/multiparticleCumulants.cxx @@ -1381,4 +1381,4 @@ struct MultiparticleCumulants { // this name is used in lower-case format to nam WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{adaptAnalysisTask(cfgc)}; -} // WorkflowSpec... \ No newline at end of file +} // WorkflowSpec...