diff --git a/PWGCF/GenericFramework/Tasks/flowGenericFramework.cxx b/PWGCF/GenericFramework/Tasks/flowGenericFramework.cxx index 57f9c2c0fea..94314746bdb 100644 --- a/PWGCF/GenericFramework/Tasks/flowGenericFramework.cxx +++ b/PWGCF/GenericFramework/Tasks/flowGenericFramework.cxx @@ -84,37 +84,9 @@ using namespace o2; using namespace o2::framework; using namespace o2::analysis::genericframework; -#define O2_DEFINE_CONFIGURABLE(NAME, TYPE, DEFAULT, HELP) Configurable NAME{#NAME, DEFAULT, HELP}; +#define O2_DEFINE_CONFIGURABLE(NAME, TYPE, DEFAULT, HELP) Configurable NAME{#NAME, (DEFAULT), (HELP)}; // NOLINT(bugprone-macro-parentheses) -namespace o2::analysis::gfw -{ -std::vector ptbinning = {0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2, 2.2, 2.4, 2.6, 2.8, 3, 3.5, 4, 5, 6, 8, 10}; -float ptpoilow = 0.2, ptpoiup = 10.0; -float ptreflow = 0.2, ptrefup = 3.0; -float ptlow = 0.2, ptup = 10.0; -int etabins = 16; -float etalow = -0.8, etaup = 0.8; -int vtxZbins = 40; -float vtxZlow = -10.0, vtxZup = 10.0; -int phibins = 72; -float philow = 0.0; -float phiup = o2::constants::math::TwoPI; -int nchbins = 300; -float nchlow = 0; -float nchup = 3000; -std::vector centbinning(90); -int nBootstrap = 10; -GFWRegions regions; -GFWCorrConfigs configs; -GFWCorrConfigs configsV02; -GFWCorrConfigs configsV0; -std::vector> etagapsPtPt; -std::vector multGlobalCorrCutPars; -std::vector multPVCorrCutPars; -std::vector multGlobalPVCorrCutPars; -std::vector multGlobalV0ACutPars; -std::vector multGlobalT0ACutPars; -} // namespace o2::analysis::gfw +constexpr float DefaultMagneticFieldCut = 99999.f; template auto projectMatrix(Array2D const& mat, std::array& array1, std::array& array2, std::array& array3) @@ -140,9 +112,9 @@ auto readMatrix(Array2D const& mat, P& array) // static constexpr float LongArrayFloat[3][20] = {{1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2}, {2.1, 2.2, 2.3, -2.1, -2.2, -2.3, 1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2}, {3.1, 3.2, 3.3, -3.1, -3.2, -3.3, 1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2, 1.3, -1.1, -1.2, -1.3, 1.1, 1.2}}; // static constexpr int LongArrayInt[3][20] = {{1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1}, {2, 2, 2, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1}, {3, 3, 3, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1}}; -static constexpr float LongArrayFloat[20][3] = {{1.1, 2.1, 3.1}, {1.2, 2.2, 3.2}, {1.3, 2.3, 3.3}, {-1.1, -2.1, -3.1}, {-1.2, -2.2, -3.2}, {-1.3, -2.3, -3.3}, {1.1, 1.1, 1.1}, {1.2, 1.2, 1.2}, {1.3, 1.3, 1.3}, {-1.1, -1.1, -1.1}, {-1.2, -1.2, -1.2}, {-1.3, -1.3, -1.3}, {1.1, 1.1, 1.1}, {1.2, 1.2, 1.2}, {1.3, 1.3, 1.3}, {-1.1, -1.1, -1.1}, {-1.2, -1.2, -1.2}, {-1.3, -1.3, -1.3}, {1.1, 1.1, 1.1}, {1.2, 1.2, 1.2}}; -static constexpr int LongArrayInt[20][3] = {{1, 2, 3}, {1, 2, 3}, {1, 2, 3}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {1, 1, 1}, {1, 1, 1}, {1, 1, 1}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {1, 1, 1}, {1, 1, 1}, {1, 1, 1}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {1, 1, 1}, {1, 1, 1}}; -static constexpr double LongArrayDouble[4][2] = {{-0.8, -0.5}, {0.5, 0.8}, {-2, -2}, {-2, -2}}; +static constexpr std::array, 20> LongArrayFloat = {{{{1.1, 2.1, 3.1}}, {{1.2, 2.2, 3.2}}, {{1.3, 2.3, 3.3}}, {{-1.1, -2.1, -3.1}}, {{-1.2, -2.2, -3.2}}, {{-1.3, -2.3, -3.3}}, {{1.1, 1.1, 1.1}}, {{1.2, 1.2, 1.2}}, {{1.3, 1.3, 1.3}}, {{-1.1, -1.1, -1.1}}, {{-1.2, -1.2, -1.2}}, {{-1.3, -1.3, -1.3}}, {{1.1, 1.1, 1.1}}, {{1.2, 1.2, 1.2}}, {{1.3, 1.3, 1.3}}, {{-1.1, -1.1, -1.1}}, {{-1.2, -1.2, -1.2}}, {{-1.3, -1.3, -1.3}}, {{1.1, 1.1, 1.1}}, {{1.2, 1.2, 1.2}}}}; +static constexpr std::array, 20> LongArrayInt = {{{{1, 2, 3}}, {{1, 2, 3}}, {{1, 2, 3}}, {{0, 0, 0}}, {{0, 0, 0}}, {{0, 0, 0}}, {{1, 1, 1}}, {{1, 1, 1}}, {{1, 1, 1}}, {{0, 0, 0}}, {{0, 0, 0}}, {{0, 0, 0}}, {{1, 1, 1}}, {{1, 1, 1}}, {{1, 1, 1}}, {{0, 0, 0}}, {{0, 0, 0}}, {{0, 0, 0}}, {{1, 1, 1}}, {{1, 1, 1}}}}; +static constexpr std::array, 4> LongArrayDouble = {{{{-0.8, -0.5}}, {{0.5, 0.8}}, {{-2, -2}}, {{-2, -2}}}}; struct FlowGenericFramework { O2_DEFINE_CONFIGURABLE(cfgNbootstrap, int, 10, "Number of subsamples") @@ -156,7 +128,6 @@ struct FlowGenericFramework { O2_DEFINE_CONFIGURABLE(cfgFillQA, bool, false, "Fill QA histograms") O2_DEFINE_CONFIGURABLE(cfgFillV0QA, bool, true, "Fill QA histograms for V0 reconstruction") O2_DEFINE_CONFIGURABLE(cfgFillRunByRunQA, bool, false, "Fill histograms on a run-by-run basis")} cfgFill; - O2_DEFINE_CONFIGURABLE(cfgMultCut, bool, false, "Use additional event cut on mult correlations"); O2_DEFINE_CONFIGURABLE(cfgUseCentralMoments, bool, true, "Use central moments in vn-pt calculations") O2_DEFINE_CONFIGURABLE(cfgUsePID, bool, true, "Enable PID information") O2_DEFINE_CONFIGURABLE(cfgUseGapMethod, bool, false, "Use gap method in vn-pt calculations") @@ -169,16 +140,16 @@ struct FlowGenericFramework { O2_DEFINE_CONFIGURABLE(cfgUsePIDEfficiencies, bool, false, "Use species dependent efficiencies") O2_DEFINE_CONFIGURABLE(cfgUse2DEfficiency, bool, false, "Toggle the use of 2D (pt, centrality) efficiency versus centrality integrated efficiency")} cfgEfficiency; O2_DEFINE_CONFIGURABLE(cfgAcceptance, std::string, "", "CCDB path to acceptance object") - O2_DEFINE_CONFIGURABLE(cfgPtmin, float, 0.2, "minimum pt (GeV/c)"); - O2_DEFINE_CONFIGURABLE(cfgPtmax, float, 10, "maximum pt (GeV/c)"); struct : ConfigurableGroup { Configurable> cfgEtaNptV02{"cfgEtaNptV02", {-0.5, 0.5}, "eta cut for v02 fractional spectra"}; Configurable> cfgEtaNptV0{"cfgEtaNptV0", {-0.5, 0.5}, "eta cut for v0 fractional spectra"}; Configurable> cfgEta{"cfgEta", {-0.8, 0.8}, "eta cut"}; Configurable> cfgEtaNch{"cfgEtaNch", {-0.5, 0.5}, "eta cut for nch selection"}; Configurable> cfgEtaPtPt{"cfgEtaPtPt", {-0.5, 0.5}, "eta for pt-pt correlation"}; + O2_DEFINE_CONFIGURABLE(cfgPtMin, float, 0.2, "minimum pt (GeV/c)"); + O2_DEFINE_CONFIGURABLE(cfgPtmax, float, 10, "maximum pt (GeV/c)"); } cfgKinematics; - Configurable> cfgPtPtGaps{"cfgPtPtGaps", {LongArrayDouble[0], 4, 2, {"subevent 1", "subevent 2", "subevent 3", "subevent 4"}, {"etamin", "etamax"}}, "{etamin,etamax} for all ptpt-subevents"}; + Configurable> cfgPtPtGaps{"cfgPtPtGaps", {LongArrayDouble.front().data(), 4, 2, {"subevent 1", "subevent 2", "subevent 3", "subevent 4"}, {"etamin", "etamax"}}, "{etamin,etamax} for all ptpt-subevents"}; O2_DEFINE_CONFIGURABLE(cfgUsePIDTotal, bool, false, "use fraction of PID total"); O2_DEFINE_CONFIGURABLE(cfgVtxZ, float, 10, "vertex cut (cm)"); struct : ConfigurableGroup { @@ -205,9 +176,12 @@ struct FlowGenericFramework { O2_DEFINE_CONFIGURABLE(cfgTVXinTRD, bool, true, "TVXinTRD - Use TVXinTRD (reject TRD triggered events)"); O2_DEFINE_CONFIGURABLE(cfgIsVertexITSTPC, bool, true, "IsVertexITSTPC - Selects collisions with at least one ITS-TPC track"); } cfgEventCutFlags; - O2_DEFINE_CONFIGURABLE(cfgOccupancySelection, int, 2000, "Max occupancy selection, -999 to disable"); - O2_DEFINE_CONFIGURABLE(cfgDoOccupancySel, bool, true, "Bool for event selection on detector occupancy"); - O2_DEFINE_CONFIGURABLE(cfgMagField, float, 99999, "Configurable magnetic field; default CCDB will be queried"); + struct ConfigurableGroup { + O2_DEFINE_CONFIGURABLE(cfgOccupancySelection, int, 2000, "Max occupancy selection, -999 to disable"); + O2_DEFINE_CONFIGURABLE(cfgDoOccupancySel, bool, true, "Bool for event selection on detector occupancy"); + O2_DEFINE_CONFIGURABLE(cfgMagField, float, 99999, "Configurable magnetic field; default CCDB will be queried"); + O2_DEFINE_CONFIGURABLE(cfgMultCut, bool, false, "Use additional event cut on mult correlations"); + } cfgEventSelection; O2_DEFINE_CONFIGURABLE(cfgUseDensityDependentCorrection, bool, false, "Use density dependent efficiency correction based on Run 2 measurements"); Configurable> cfgTrackDensityP0{"cfgTrackDensityP0", std::vector{0.7217476707, 0.7384792571, 0.7542625668, 0.7640680200, 0.7701951667, 0.7755299053, 0.7805901710, 0.7849446786, 0.7957356586, 0.8113039262, 0.8211968966, 0.8280558878, 0.8329342135}, "parameter 0 for track density efficiency correction"}; Configurable> cfgTrackDensityP1{"cfgTrackDensityP1", std::vector{-2.169488e-05, -2.191913e-05, -2.295484e-05, -2.556538e-05, -2.754463e-05, -2.816832e-05, -2.846502e-05, -2.843857e-05, -2.705974e-05, -2.477018e-05, -2.321730e-05, -2.203315e-05, -2.109474e-05}, "parameter 1 for track density efficiency correction"}; @@ -227,9 +201,9 @@ struct FlowGenericFramework { O2_DEFINE_CONFIGURABLE(cfgGlobalT0AHighSigma, float, 4, "Number of sigma deviations above expected value in global vs T0A correlation"); } cfgMultCorrCuts; struct : ConfigurableGroup { - Configurable> nSigmas{"nSigmas", {LongArrayFloat[0], 6, 3, {"pos_pi", "pos_ka", "pos_pr", "neg_pi", "neg_ka", "neg_pr"}, {"TPC", "TOF", "ITS"}}, "Labeled array for n-sigma values for TPC, TOF, ITS for pions, kaons, protons (positive and negative)"}; - Configurable> resonanceCuts{"resonanceCuts", {LongArrayFloat[0], 14, 3, {"cos_PAs", "massMin", "massMax", "PosTrackPt", "NegTrackPt", "DCAPosToPVMin", "DCANegToPVMin", "DCAxDaughters", "Lifetime", "RadiusMin", "RadiusMax", "Rapidity", "ArmPodMin", "MassRejection"}, {"K0", "Lambda", "Phi"}}, "Labeled array (float) for various cuts on resonances"}; - Configurable> resonanceSwitches{"resonanceSwitches", {LongArrayInt[0], 8, 3, {"UseParticle", "UseCosPA", "NMassBins", "UseDCAxDaughters", "UseProperLifetime", "UseV0Radius", "UseArmPodCut", "UseCompetingMassRejection"}, {"K0", "Lambda", "Phi"}}, "Labeled array (int) for various cuts on resonances"}; + Configurable> nSigmas{"nSigmas", {LongArrayFloat.front().data(), 6, 3, {"pos_pi", "pos_ka", "pos_pr", "neg_pi", "neg_ka", "neg_pr"}, {"TPC", "TOF", "ITS"}}, "Labeled array for n-sigma values for TPC, TOF, ITS for pions, kaons, protons (positive and negative)"}; + Configurable> resonanceCuts{"resonanceCuts", {LongArrayFloat.front().data(), 14, 3, {"cos_PAs", "massMin", "massMax", "PosTrackPt", "NegTrackPt", "DCAPosToPVMin", "DCANegToPVMin", "DCAxDaughters", "Lifetime", "RadiusMin", "RadiusMax", "Rapidity", "ArmPodMin", "MassRejection"}, {"K0", "Lambda", "Phi"}}, "Labeled array (float) for various cuts on resonances"}; + Configurable> resonanceSwitches{"resonanceSwitches", {LongArrayInt.front().data(), 8, 3, {"UseParticle", "UseCosPA", "NMassBins", "UseDCAxDaughters", "UseProperLifetime", "UseV0Radius", "UseArmPodCut", "UseCompetingMassRejection"}, {"K0", "Lambda", "Phi"}}, "Labeled array (int) for various cuts on resonances"}; O2_DEFINE_CONFIGURABLE(cfgUseLsPhi, bool, true, "Use LikeSign for Phi v2") O2_DEFINE_CONFIGURABLE(cfgUseOnlyTPC, bool, true, "Use only TPC PID for daughter selection") O2_DEFINE_CONFIGURABLE(cfgFakeKaonCut, float, 0.1f, "Maximum difference in measured momentum and TPC inner ring momentum of particle") @@ -262,9 +236,42 @@ struct FlowGenericFramework { ConfigurableAxis axisNsigmaTPC{"axisNsigmaTPC", {80, -5, 5}, "nsigmaTPC axis"}; ConfigurableAxis axisNsigmaTOF{"axisNsigmaTOF", {80, -5, 5}, "nsigmaTOF axis"}; + struct GFWMemberCache { + std::vector ptbinning = {0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2, 2.2, 2.4, 2.6, 2.8, 3, 3.5, 4, 5, 6, 8, 10}; + float ptpoilow = 0.2; + float ptpoiup = 10.0; + float ptreflow = 0.2; + float ptrefup = 3.0; + float ptlow = 0.2; + float ptup = 10.0; + int etabins = 16; + float etalow = -0.8; + float etaup = 0.8; + int vtxZbins = 40; + float vtxZlow = -10.0; + float vtxZup = 10.0; + int phibins = 72; + float philow = 0.0; + float phiup = o2::constants::math::TwoPI; + int nchbins = 300; + float nchlow = 0; + float nchup = 3000; + std::vector centbinning = std::vector(90); + GFWRegions regions; + GFWCorrConfigs configs; + GFWCorrConfigs configsV02; + GFWCorrConfigs configsV0; + std::vector> etagapsPtPt; + std::vector multGlobalCorrCutPars; + std::vector multPVCorrCutPars; + std::vector multGlobalPVCorrCutPars; + std::vector multGlobalV0ACutPars; + std::vector multGlobalT0ACutPars; + } gfwMemberCache; + // Connect to ccdb - Service ccdb; - Service pdgDB; + Service ccdb{}; + Service pdgDB{}; o2::aod::ITSResponse itsResponse; @@ -282,11 +289,11 @@ struct FlowGenericFramework { HistogramRegistry registry{"registry"}; HistogramRegistry registryQA{"registryQA"}; - std::array, 14> resoCutVals; - std::array, 8> resoSwitchVals; - std::array tofNsigmaCut; - std::array itsNsigmaCut; - std::array tpcNsigmaCut; + std::array, 14> resoCutVals{}; + std::array, 8> resoSwitchVals{}; + std::array tofNsigmaCut{}; + std::array itsNsigmaCut{}; + std::array tpcNsigmaCut{}; enum FractionSetup { FractionV02 = 0, @@ -370,11 +377,11 @@ struct FlowGenericFramework { SpeciesCount }; enum EfficiencySpecies { - EfficiencyCharged = ChargedID, - EfficiencyPion = PionID, - EfficiencyKaon = KaonID, - EfficiencyProton = ProtonID, - EfficiencyK0 = 4, + EfficiencyCharged = 0, + EfficiencyPion, + EfficiencyKaon, + EfficiencyProton, + EfficiencyK0, EfficiencyLambda, EfficiencySpeciesCount }; @@ -402,10 +409,10 @@ struct FlowGenericFramework { }; enum OutputSpecies { K0 = 0, - Lambda = 1, - PhiMeson = 2, - LambdaBar = 3, - RefParticle = 4, + Lambda, + PhiMeson, + LambdaBar, + RefParticle, OutputSpeciesCount }; enum ParticleCuts { @@ -459,25 +466,24 @@ struct FlowGenericFramework { std::vector corrconfigsV0; TRandom3* fRndm = new TRandom3(0); - TAxis* fPtAxis; + TAxis* fPtAxis = nullptr; int lastRun = -1; std::vector runNumbers; // Density dependent eff correction std::vector funcEff; - TH1D* hFindPtBin; - TF1* funcV2; - TF1* funcV3; - TF1* funcV4; + TH1D* hFindPtBin = nullptr; + TF1* funcV2 = nullptr; + TF1* funcV3 = nullptr; + TF1* funcV4 = nullptr; struct DensityCorr { - double psi2Est; - double psi3Est; - double psi4Est; - double v2; - double v3; - double v4; - int density; - DensityCorr() : psi2Est(0.), psi3Est(0.), psi4Est(0.), v2(0.), v3(0.), v4(0.), density(0) {} + double psi2Est = 0.; + double psi3Est = 0.; + double psi4Est = 0.; + double v2 = 0.; + double v3 = 0.; + double v4 = 0.; + int density = 0; }; // Event selection cuts - multiplicity correlation @@ -502,64 +508,70 @@ struct FlowGenericFramework { void init(InitContext const&) { LOGF(info, "FlowGenericFramework::init()"); - o2::analysis::gfw::regions.SetNames(cfgRegions->GetNames()); - o2::analysis::gfw::regions.SetEtaMin(cfgRegions->GetEtaMin()); - o2::analysis::gfw::regions.SetEtaMax(cfgRegions->GetEtaMax()); - o2::analysis::gfw::regions.SetpTDifs(cfgRegions->GetpTDifs()); - o2::analysis::gfw::regions.SetBitmasks(cfgRegions->GetBitmasks()); - o2::analysis::gfw::configs.SetCorrs(cfgCorrConfig->GetCorrs()); - o2::analysis::gfw::configs.SetHeads(cfgCorrConfig->GetHeads()); - o2::analysis::gfw::configs.SetpTDifs(cfgCorrConfig->GetpTDifs()); - o2::analysis::gfw::configs.SetpTCorrMasks(cfgCorrConfig->GetpTCorrMasks()); - o2::analysis::gfw::regions.Print(); - o2::analysis::gfw::configs.Print(); - o2::analysis::gfw::configsV02.SetCorrs(cfgCorrConfigV02->GetCorrs()); - o2::analysis::gfw::configsV02.SetHeads(cfgCorrConfigV02->GetHeads()); - o2::analysis::gfw::configsV02.SetpTDifs(cfgCorrConfigV02->GetpTDifs()); - o2::analysis::gfw::configsV02.SetpTCorrMasks(cfgCorrConfigV02->GetpTCorrMasks()); - o2::analysis::gfw::configsV02.Print(); - o2::analysis::gfw::configsV0.SetCorrs(cfgCorrConfigV0->GetCorrs()); - o2::analysis::gfw::configsV0.SetHeads(cfgCorrConfigV0->GetHeads()); - o2::analysis::gfw::configsV0.SetpTDifs(cfgCorrConfigV0->GetpTDifs()); - o2::analysis::gfw::configsV0.SetpTCorrMasks(cfgCorrConfigV0->GetpTCorrMasks()); - o2::analysis::gfw::configsV0.Print(); - o2::analysis::gfw::ptbinning = cfgGFWBinning->GetPtBinning(); - o2::analysis::gfw::ptpoilow = cfgGFWBinning->GetPtPOImin(); - o2::analysis::gfw::ptpoiup = cfgGFWBinning->GetPtPOImax(); - o2::analysis::gfw::ptreflow = cfgGFWBinning->GetPtRefMin(); - o2::analysis::gfw::ptrefup = cfgGFWBinning->GetPtRefMax(); - o2::analysis::gfw::ptlow = cfgPtmin; - o2::analysis::gfw::ptup = cfgPtmax; - o2::analysis::gfw::etabins = cfgGFWBinning->GetEtaBins(); - o2::analysis::gfw::vtxZbins = cfgGFWBinning->GetVtxZbins(); - o2::analysis::gfw::phibins = cfgGFWBinning->GetPhiBins(); - o2::analysis::gfw::philow = 0.0f; - o2::analysis::gfw::phiup = o2::constants::math::TwoPI; - o2::analysis::gfw::nchbins = cfgGFWBinning->GetNchBins(); - o2::analysis::gfw::nchlow = cfgGFWBinning->GetNchMin(); - o2::analysis::gfw::nchup = cfgGFWBinning->GetNchMax(); - o2::analysis::gfw::centbinning = cfgGFWBinning->GetCentBinning(); + gfwMemberCache.regions.SetNames(cfgRegions->GetNames()); + gfwMemberCache.regions.SetEtaMin(cfgRegions->GetEtaMin()); + gfwMemberCache.regions.SetEtaMax(cfgRegions->GetEtaMax()); + gfwMemberCache.regions.SetpTDifs(cfgRegions->GetpTDifs()); + gfwMemberCache.regions.SetBitmasks(cfgRegions->GetBitmasks()); + gfwMemberCache.configs.SetCorrs(cfgCorrConfig->GetCorrs()); + gfwMemberCache.configs.SetHeads(cfgCorrConfig->GetHeads()); + gfwMemberCache.configs.SetpTDifs(cfgCorrConfig->GetpTDifs()); + gfwMemberCache.configs.SetpTCorrMasks(cfgCorrConfig->GetpTCorrMasks()); + gfwMemberCache.regions.Print(); + gfwMemberCache.configs.Print(); + gfwMemberCache.configsV02.SetCorrs(cfgCorrConfigV02->GetCorrs()); + gfwMemberCache.configsV02.SetHeads(cfgCorrConfigV02->GetHeads()); + gfwMemberCache.configsV02.SetpTDifs(cfgCorrConfigV02->GetpTDifs()); + gfwMemberCache.configsV02.SetpTCorrMasks(cfgCorrConfigV02->GetpTCorrMasks()); + gfwMemberCache.configsV02.Print(); + gfwMemberCache.configsV0.SetCorrs(cfgCorrConfigV0->GetCorrs()); + gfwMemberCache.configsV0.SetHeads(cfgCorrConfigV0->GetHeads()); + gfwMemberCache.configsV0.SetpTDifs(cfgCorrConfigV0->GetpTDifs()); + gfwMemberCache.configsV0.SetpTCorrMasks(cfgCorrConfigV0->GetpTCorrMasks()); + gfwMemberCache.configsV0.Print(); + gfwMemberCache.ptbinning = cfgGFWBinning->GetPtBinning(); + gfwMemberCache.ptpoilow = cfgGFWBinning->GetPtPOImin(); + gfwMemberCache.ptpoiup = cfgGFWBinning->GetPtPOImax(); + gfwMemberCache.ptreflow = cfgGFWBinning->GetPtRefMin(); + gfwMemberCache.ptrefup = cfgGFWBinning->GetPtRefMax(); + gfwMemberCache.ptlow = cfgKinematics.cfgPtMin; + gfwMemberCache.ptup = cfgKinematics.cfgPtmax; + gfwMemberCache.etabins = cfgGFWBinning->GetEtaBins(); + gfwMemberCache.vtxZbins = cfgGFWBinning->GetVtxZbins(); + gfwMemberCache.phibins = cfgGFWBinning->GetPhiBins(); + gfwMemberCache.philow = 0.0f; + gfwMemberCache.phiup = o2::constants::math::TwoPI; + gfwMemberCache.nchbins = cfgGFWBinning->GetNchBins(); + gfwMemberCache.nchlow = cfgGFWBinning->GetNchMin(); + gfwMemberCache.nchup = cfgGFWBinning->GetNchMax(); + gfwMemberCache.centbinning = cfgGFWBinning->GetCentBinning(); cfgGFWBinning->Print(); LOGF(info, "Eta cuts: Filter [%.1f,%.1f] | Nch [%.1f,%.1f] | Npt v02 [%.1f, %.1f] | Npt v0 [%.1f, %.1f] | Pt-Pt [%.1f, %.1f]", cfgKinematics.cfgEta->first, cfgKinematics.cfgEta->second, cfgKinematics.cfgEtaNch->first, cfgKinematics.cfgEtaNch->second, cfgKinematics.cfgEtaNptV02->first, cfgKinematics.cfgEtaNptV02->second, cfgKinematics.cfgEtaNptV0->first, cfgKinematics.cfgEtaNptV0->second, cfgKinematics.cfgEtaPtPt->first, cfgKinematics.cfgEtaPtPt->second); - o2::analysis::gfw::multGlobalCorrCutPars = cfgMultCorrCuts.cfgMultGlobalCutPars; - o2::analysis::gfw::multPVCorrCutPars = cfgMultCorrCuts.cfgMultPVCutPars; - o2::analysis::gfw::multGlobalPVCorrCutPars = cfgMultCorrCuts.cfgMultGlobalPVCutPars; - o2::analysis::gfw::multGlobalV0ACutPars = cfgMultCorrCuts.cfgMultGlobalV0ACutPars; - o2::analysis::gfw::multGlobalT0ACutPars = cfgMultCorrCuts.cfgMultGlobalT0ACutPars; + gfwMemberCache.multGlobalCorrCutPars = cfgMultCorrCuts.cfgMultGlobalCutPars; + gfwMemberCache.multPVCorrCutPars = cfgMultCorrCuts.cfgMultPVCutPars; + gfwMemberCache.multGlobalPVCorrCutPars = cfgMultCorrCuts.cfgMultGlobalPVCutPars; + gfwMemberCache.multGlobalV0ACutPars = cfgMultCorrCuts.cfgMultGlobalV0ACutPars; + gfwMemberCache.multGlobalT0ACutPars = cfgMultCorrCuts.cfgMultGlobalT0ACutPars; projectMatrix(cfgPIDCuts.nSigmas->getData(), tpcNsigmaCut, tofNsigmaCut, itsNsigmaCut); readMatrix(cfgPIDCuts.resonanceCuts->getData(), resoCutVals); readMatrix(cfgPIDCuts.resonanceSwitches->getData(), resoSwitchVals); printResoCuts(); - int nPtPtSubMax = 4; - for (int i = 0; i < nPtPtSubMax; ++i) { - if (cfgPtPtGaps->getData()[i][0] < -1. || cfgPtPtGaps->getData()[i][1] < -1.) + const auto& ptPtGaps = cfgPtPtGaps->getData(); + + for (uint32_t i = 0; i < ptPtGaps.rows; ++i) { + const auto etaMin = ptPtGaps(i, 0); + const auto etaMax = ptPtGaps(i, 1); + + if (etaMin < -1. || etaMax < -1.) { continue; - o2::analysis::gfw::etagapsPtPt.push_back(std::make_pair(cfgPtPtGaps->getData()[i][0], cfgPtPtGaps->getData()[i][1])); + } + + gfwMemberCache.etagapsPtPt.emplace_back(etaMin, etaMax); } - for (const auto& [etamin, etamax] : o2::analysis::gfw::etagapsPtPt) { + for (const auto& [etamin, etamax] : gfwMemberCache.etagapsPtPt) { LOGF(info, "pt-pt subevent: {%.1f,%.1f}", etamin, etamax); } @@ -577,17 +589,17 @@ struct FlowGenericFramework { LOGF(info, "Flag %d is %senabled", cut.histBin, (cut.enabled) ? "" : "not "); } - AxisSpec phiAxis = {o2::analysis::gfw::phibins, o2::analysis::gfw::philow, o2::analysis::gfw::phiup, "#phi"}; + AxisSpec phiAxis = {gfwMemberCache.phibins, gfwMemberCache.philow, gfwMemberCache.phiup, "#phi"}; AxisSpec phiModAxis = {100, 0, constants::math::PI / 9, "fmod(#varphi,#pi/9)"}; - AxisSpec etaAxis = {o2::analysis::gfw::etabins, cfgKinematics.cfgEta->first, cfgKinematics.cfgEta->second, "#eta"}; - AxisSpec vtxAxis = {o2::analysis::gfw::vtxZbins, -cfgVtxZ, cfgVtxZ, "Vtx_{z} (cm)"}; - AxisSpec ptAxis = {o2::analysis::gfw::ptbinning, "#it{p}_{T} GeV/#it{c}"}; + AxisSpec etaAxis = {gfwMemberCache.etabins, cfgKinematics.cfgEta->first, cfgKinematics.cfgEta->second, "#eta"}; + AxisSpec vtxAxis = {gfwMemberCache.vtxZbins, -cfgVtxZ, cfgVtxZ, "Vtx_{z} (cm)"}; + AxisSpec ptAxis = {gfwMemberCache.ptbinning, "#it{p}_{T} GeV/#it{c}"}; std::string sCentralityEstimator = centNamesMap[cfgCentEstimator] + " centrality (%)"; - AxisSpec centAxis = {o2::analysis::gfw::centbinning, sCentralityEstimator.c_str()}; + AxisSpec centAxis = {gfwMemberCache.centbinning, sCentralityEstimator.c_str()}; std::vector nchbinning; - int nchskip = (o2::analysis::gfw::nchup - o2::analysis::gfw::nchlow) / o2::analysis::gfw::nchbins; - for (int i = 0; i <= o2::analysis::gfw::nchbins; ++i) { - nchbinning.push_back(nchskip * i + o2::analysis::gfw::nchlow + 0.5); + const double nchskip = (gfwMemberCache.nchup - gfwMemberCache.nchlow) / gfwMemberCache.nchbins; + for (int i = 0; i <= gfwMemberCache.nchbins; ++i) { + nchbinning.push_back(nchskip * i + gfwMemberCache.nchlow + 0.5); } AxisSpec nchAxis = {nchbinning, "N_{ch}"}; AxisSpec bAxis = {200, 0, 20, "#it{b}"}; @@ -611,16 +623,17 @@ struct FlowGenericFramework { int64_t now = std::chrono::duration_cast(std::chrono::system_clock::now().time_since_epoch()).count(); ccdb->setCreatedNotAfter(now); - int ptbins = o2::analysis::gfw::ptbinning.size() - 1; - fPtAxis = new TAxis(ptbins, &o2::analysis::gfw::ptbinning[0]); + const int ptbins = static_cast(gfwMemberCache.ptbinning.size() - 1); + fPtAxis = new TAxis(ptbins, gfwMemberCache.ptbinning.data()); if (doprocessMCGen || doprocessOnTheFly) { if (cfgFill.cfgFillQA) { registryQA.add("MCGen/before/pt_gen", "", {HistType::kTH1D, {ptAxis}}); registryQA.add("MCGen/before/phi_eta_vtxZ_gen", "", {HistType::kTH3D, {phiAxis, etaAxis, vtxAxis}}); registryQA.addClone("MCGen/before/", "MCGen/after/"); - if (doprocessOnTheFly) + if (doprocessOnTheFly) { registryQA.add("MCGen/impactParameter", "", {HistType::kTH2D, {{bAxis, nchAxis}}}); + } } } if (doprocessEfficiency) { @@ -639,8 +652,8 @@ struct FlowGenericFramework { registryQA.add("trackQA/before/nTPCCrossedRows", "Number of crossed TPC Rows", {HistType::kTH1D, {{100, 40, 180}}}); registryQA.addClone("trackQA/before/", "trackQA/after/"); - registryQA.add("trackQA/after/pt_ref", "; #it{p}_{T}; Counts", {HistType::kTH1D, {{100, o2::analysis::gfw::ptreflow, o2::analysis::gfw::ptrefup}}}); - registryQA.add("trackQA/after/pt_poi", "; #it{p}_{T}; Counts", {HistType::kTH1D, {{100, o2::analysis::gfw::ptpoilow, o2::analysis::gfw::ptpoiup}}}); + registryQA.add("trackQA/after/pt_ref", "; #it{p}_{T}; Counts", {HistType::kTH1D, {{100, gfwMemberCache.ptreflow, gfwMemberCache.ptrefup}}}); + registryQA.add("trackQA/after/pt_poi", "; #it{p}_{T}; Counts", {HistType::kTH1D, {{100, gfwMemberCache.ptpoilow, gfwMemberCache.ptpoiup}}}); registryQA.add("trackQA/after/Nch_corrected", "; N_{ch}; Counts", {HistType::kTH1D, {nchAxis}}); registryQA.add("trackQA/after/Nch_uncorrected", "; N_{ch}; Counts", {HistType::kTH1D, {nchAxis}}); registryQA.add("trackQA/after/etaNch", "; #eta; Counts", {HistType::kTH1D, {etaAxis}}); @@ -691,7 +704,7 @@ struct FlowGenericFramework { AxisSpec axisLambdaMass = {resoSwitchVals[MassBins][Lambda], resoCutVals[MassMin][Lambda], resoCutVals[MassMax][Lambda]}; // QA histograms for V0s - if (resoSwitchVals[UseParticle][K0]) { + if (resoSwitchVals[UseParticle][K0] != 0) { if (cfgFill.cfgFillV0QA) { registryQA.add("K0/PiPlusTPC_K0", "", {HistType::kTH2D, {{ptAxis, axisNsigmaTPC}}}); registryQA.add("K0/PiMinusTPC_K0", "", {HistType::kTH2D, {{ptAxis, axisNsigmaTPC}}}); @@ -699,7 +712,6 @@ struct FlowGenericFramework { registryQA.add("K0/PiMinusTOF_K0", "", {HistType::kTH2D, {{ptAxis, axisNsigmaTOF}}}); registryQA.add("K0/hK0Phi", "", {HistType::kTH1D, {phiAxis}}); registryQA.add("K0/hK0Eta", "", {HistType::kTH1D, {etaAxis}}); - registryQA.add("K0/hK0Eta", "", {HistType::kTH1D, {etaAxis}}); registryQA.add("K0/hK0AP", "", {HistType::kTH2D, {{100, -1, 1}, {100, 0, 5.}}}); registryQA.add("K0/hK0Mass_sparse", "", {HistType::kTHnSparseF, {{axisK0Mass, ptAxis, nchAxis}}}); registryQA.add("K0/hK0s", "", {HistType::kTH1D, {singleCount}}); @@ -721,7 +733,7 @@ struct FlowGenericFramework { registryQA.get(HIST("K0/hK0Count"))->GetXaxis()->SetBinLabel(FillV0DaughterTrackSelection, "v0 Daughter track selection"); } - if (resoSwitchVals[UseParticle][Lambda]) { + if (resoSwitchVals[UseParticle][Lambda] != 0) { if (cfgFill.cfgFillV0QA) { registryQA.add("Lambda/PrPlusTPC_L", "", {HistType::kTH2D, {{ptAxis, axisNsigmaTPC}}}); registryQA.add("Lambda/PiMinusTPC_L", "", {HistType::kTH2D, {{ptAxis, axisNsigmaTPC}}}); @@ -797,31 +809,35 @@ struct FlowGenericFramework { histosResoNpt[setup][LambdaSideband2] = std::make_unique(Form("npt%sLambdaSB2", setupName), "; #it{p}_{T} (GeV/#it{c}; Count", ptAxis.binEdges.size() - 1, ptAxis.binEdges.data()); } - if (o2::analysis::gfw::regions.GetSize() < 0) + if (gfwMemberCache.regions.GetSize() < 0) { LOGF(error, "Configuration contains vectors of different size - check the GFWRegions configurable"); - for (auto i(0); i < o2::analysis::gfw::regions.GetSize(); ++i) { - fGFW->AddRegion(o2::analysis::gfw::regions.GetNames()[i], o2::analysis::gfw::regions.GetEtaMin()[i], o2::analysis::gfw::regions.GetEtaMax()[i], (o2::analysis::gfw::regions.GetpTDifs()[i]) ? ptbins + 1 : 1, o2::analysis::gfw::regions.GetBitmasks()[i]); } - for (auto i = 0; i < o2::analysis::gfw::configs.GetSize(); ++i) { - corrconfigs.push_back(fGFW->GetCorrelatorConfig(o2::analysis::gfw::configs.GetCorrs()[i], o2::analysis::gfw::configs.GetHeads()[i], o2::analysis::gfw::configs.GetpTDifs()[i])); + for (auto i(0); i < gfwMemberCache.regions.GetSize(); ++i) { + fGFW->AddRegion(gfwMemberCache.regions.GetNames()[i], gfwMemberCache.regions.GetEtaMin()[i], gfwMemberCache.regions.GetEtaMax()[i], (gfwMemberCache.regions.GetpTDifs()[i] != 0) ? ptbins + 1 : 1, gfwMemberCache.regions.GetBitmasks()[i]); } - if (corrconfigs.empty()) + for (auto i = 0; i < gfwMemberCache.configs.GetSize(); ++i) { + corrconfigs.push_back(fGFW->GetCorrelatorConfig(gfwMemberCache.configs.GetCorrs()[i], gfwMemberCache.configs.GetHeads()[i], gfwMemberCache.configs.GetpTDifs()[i] != 0)); + } + if (corrconfigs.empty()) { LOGF(error, "Configuration contains vectors of different size - check the GFWCorrConfig configurable"); + } // Radial flow configs - for (auto i = 0; i < o2::analysis::gfw::configsV02.GetSize(); ++i) { - corrconfigsV02.push_back(fGFW->GetCorrelatorConfig(o2::analysis::gfw::configsV02.GetCorrs()[i], o2::analysis::gfw::configsV02.GetHeads()[i], o2::analysis::gfw::configsV02.GetpTDifs()[i])); + for (auto i = 0; i < gfwMemberCache.configsV02.GetSize(); ++i) { + corrconfigsV02.push_back(fGFW->GetCorrelatorConfig(gfwMemberCache.configsV02.GetCorrs()[i], gfwMemberCache.configsV02.GetHeads()[i], gfwMemberCache.configsV02.GetpTDifs()[i] != 0)); } - if (corrconfigsV02.empty()) + if (corrconfigsV02.empty()) { LOGF(error, "Radial (V02) configuration contains vectors of different size - check the GFWCorrConfig configurable"); - for (auto i = 0; i < o2::analysis::gfw::configsV0.GetSize(); ++i) { - corrconfigsV0.push_back(fGFW->GetCorrelatorConfig(o2::analysis::gfw::configsV0.GetCorrs()[i], o2::analysis::gfw::configsV0.GetHeads()[i], o2::analysis::gfw::configsV0.GetpTDifs()[i])); } - if (corrconfigsV0.empty()) + for (auto i = 0; i < gfwMemberCache.configsV0.GetSize(); ++i) { + corrconfigsV0.push_back(fGFW->GetCorrelatorConfig(gfwMemberCache.configsV0.GetCorrs()[i], gfwMemberCache.configsV0.GetHeads()[i], gfwMemberCache.configsV0.GetpTDifs()[i] != 0)); + } + if (corrconfigsV0.empty()) { LOGF(error, "Radial (V0) configuration contains vectors of different size - check the GFWCorrConfig configurable"); + } fGFW->CreateRegions(); - TObjArray* oba = new TObjArray(); + auto oba = new TObjArray(); addConfigObjectsToObjArray(oba, corrconfigs); addConfigObjectsToObjArray(oba, corrconfigsV02); addConfigObjectsToObjArray(oba, corrconfigsV0); @@ -840,8 +856,8 @@ struct FlowGenericFramework { fFCpt->setUseCentralMoments(cfgUseCentralMoments); fFCpt->setUseGapMethod(cfgUseGapMethod); - fFCpt->initialise(multAxis, cfgMpar, o2::analysis::gfw::configs, cfgNbootstrap); - fFCpt->initialiseSubevent(multAxis, cfgMpar, o2::analysis::gfw::etagapsPtPt.size(), cfgNbootstrap); + fFCpt->initialise(multAxis, cfgMpar, gfwMemberCache.configs, cfgNbootstrap); + fFCpt->initialiseSubevent(multAxis, cfgMpar, gfwMemberCache.etagapsPtPt.size(), cfgNbootstrap); fPtDepDCAxy = new TF1("ptDepDCAxy", Form("[0]*%s", cfgTrackCuts.cfgDCAxyPtDep->c_str()), 0.001, 100); fPtDepDCAxy->SetParameter(0, cfgTrackCuts.cfgDCAxyNSigma / 7.); @@ -852,47 +868,55 @@ struct FlowGenericFramework { } // Multiplicity correlation cuts - if (cfgMultCut) { + if (cfgEventSelection.cfgMultCut) { fMultPVCutLow = new TF1("fMultPVCutLow", cfgMultCorrCuts.cfgMultCorrLowCutFunction->c_str(), 0, 100); - fMultPVCutLow->SetParameters(&(o2::analysis::gfw::multPVCorrCutPars[0])); + fMultPVCutLow->SetParameters(gfwMemberCache.multPVCorrCutPars.data()); fMultPVCutHigh = new TF1("fMultPVCutHigh", cfgMultCorrCuts.cfgMultCorrHighCutFunction->c_str(), 0, 100); - fMultPVCutHigh->SetParameters(&(o2::analysis::gfw::multPVCorrCutPars[0])); + fMultPVCutHigh->SetParameters(gfwMemberCache.multPVCorrCutPars.data()); fMultCutLow = new TF1("fMultCutLow", cfgMultCorrCuts.cfgMultCorrLowCutFunction->c_str(), 0, 100); - fMultCutLow->SetParameters(&(o2::analysis::gfw::multGlobalCorrCutPars[0])); + fMultCutLow->SetParameters(gfwMemberCache.multGlobalCorrCutPars.data()); fMultCutHigh = new TF1("fMultCutHigh", cfgMultCorrCuts.cfgMultCorrHighCutFunction->c_str(), 0, 100); - fMultCutHigh->SetParameters(&(o2::analysis::gfw::multGlobalCorrCutPars[0])); + fMultCutHigh->SetParameters(gfwMemberCache.multGlobalCorrCutPars.data()); fMultPVGlobalCutHigh = new TF1("fMultPVGlobalCutHigh", cfgMultCorrCuts.cfgMultGlobalPVCorrCutFunction->c_str(), 0, nchbinning.back()); - fMultPVGlobalCutHigh->SetParameters(&(o2::analysis::gfw::multGlobalPVCorrCutPars[0])); + fMultPVGlobalCutHigh->SetParameters(gfwMemberCache.multGlobalPVCorrCutPars.data()); LOGF(info, "Global V0A function: %s in range 0-%g", cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), v0aAxis.binEdges.back()); fMultGlobalV0ACutLow = new TF1("fMultGlobalV0ACutLow", cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, v0aAxis.binEdges.back()); - for (std::size_t i = 0; i < o2::analysis::gfw::multGlobalV0ACutPars.size(); ++i) - fMultGlobalV0ACutLow->SetParameter(i, o2::analysis::gfw::multGlobalV0ACutPars[i]); - fMultGlobalV0ACutLow->SetParameter(o2::analysis::gfw::multGlobalV0ACutPars.size(), cfgMultCorrCuts.cfgGlobalV0ALowSigma); - for (int i = 0; i < fMultGlobalV0ACutLow->GetNpar(); ++i) + for (std::size_t i = 0; i < gfwMemberCache.multGlobalV0ACutPars.size(); ++i) { + fMultGlobalV0ACutLow->SetParameter(i, gfwMemberCache.multGlobalV0ACutPars[i]); + } + fMultGlobalV0ACutLow->SetParameter(gfwMemberCache.multGlobalV0ACutPars.size(), cfgMultCorrCuts.cfgGlobalV0ALowSigma); + for (int i = 0; i < fMultGlobalV0ACutLow->GetNpar(); ++i) { LOGF(info, "fMultGlobalV0ACutLow par %d = %g", i, fMultGlobalV0ACutLow->GetParameter(i)); + } fMultGlobalV0ACutHigh = new TF1("fMultGlobalV0ACutHigh", cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, v0aAxis.binEdges.back()); - for (std::size_t i = 0; i < o2::analysis::gfw::multGlobalV0ACutPars.size(); ++i) - fMultGlobalV0ACutHigh->SetParameter(i, o2::analysis::gfw::multGlobalV0ACutPars[i]); - fMultGlobalV0ACutHigh->SetParameter(o2::analysis::gfw::multGlobalV0ACutPars.size(), cfgMultCorrCuts.cfgGlobalV0AHighSigma); - for (int i = 0; i < fMultGlobalV0ACutHigh->GetNpar(); ++i) + for (std::size_t i = 0; i < gfwMemberCache.multGlobalV0ACutPars.size(); ++i) { + fMultGlobalV0ACutHigh->SetParameter(i, gfwMemberCache.multGlobalV0ACutPars[i]); + } + fMultGlobalV0ACutHigh->SetParameter(gfwMemberCache.multGlobalV0ACutPars.size(), cfgMultCorrCuts.cfgGlobalV0AHighSigma); + for (int i = 0; i < fMultGlobalV0ACutHigh->GetNpar(); ++i) { LOGF(info, "fMultGlobalV0ACutHigh par %d = %g", i, fMultGlobalV0ACutHigh->GetParameter(i)); + } LOGF(info, "Global T0A function: %s", cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str()); fMultGlobalT0ACutLow = new TF1("fMultGlobalT0ACutLow", cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, t0aAxis.binEdges.back()); - for (std::size_t i = 0; i < o2::analysis::gfw::multGlobalT0ACutPars.size(); ++i) - fMultGlobalT0ACutLow->SetParameter(i, o2::analysis::gfw::multGlobalT0ACutPars[i]); - fMultGlobalT0ACutLow->SetParameter(o2::analysis::gfw::multGlobalT0ACutPars.size(), cfgMultCorrCuts.cfgGlobalT0ALowSigma); - for (int i = 0; i < fMultGlobalT0ACutLow->GetNpar(); ++i) + for (std::size_t i = 0; i < gfwMemberCache.multGlobalT0ACutPars.size(); ++i) { + fMultGlobalT0ACutLow->SetParameter(i, gfwMemberCache.multGlobalT0ACutPars[i]); + } + fMultGlobalT0ACutLow->SetParameter(gfwMemberCache.multGlobalT0ACutPars.size(), cfgMultCorrCuts.cfgGlobalT0ALowSigma); + for (int i = 0; i < fMultGlobalT0ACutLow->GetNpar(); ++i) { LOGF(info, "fMultGlobalT0ACutLow par %d = %g", i, fMultGlobalT0ACutLow->GetParameter(i)); + } fMultGlobalT0ACutHigh = new TF1("fMultGlobalT0ACutHigh", cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->c_str(), 0, t0aAxis.binEdges.back()); - for (std::size_t i = 0; i < o2::analysis::gfw::multGlobalT0ACutPars.size(); ++i) - fMultGlobalT0ACutHigh->SetParameter(i, o2::analysis::gfw::multGlobalT0ACutPars[i]); - fMultGlobalT0ACutHigh->SetParameter(o2::analysis::gfw::multGlobalT0ACutPars.size(), cfgMultCorrCuts.cfgGlobalT0AHighSigma); - for (int i = 0; i < fMultGlobalT0ACutHigh->GetNpar(); ++i) + for (std::size_t i = 0; i < gfwMemberCache.multGlobalT0ACutPars.size(); ++i) { + fMultGlobalT0ACutHigh->SetParameter(i, gfwMemberCache.multGlobalT0ACutPars[i]); + } + fMultGlobalT0ACutHigh->SetParameter(gfwMemberCache.multGlobalT0ACutPars.size(), cfgMultCorrCuts.cfgGlobalT0AHighSigma); + for (int i = 0; i < fMultGlobalT0ACutHigh->GetNpar(); ++i) { LOGF(info, "fMultGlobalT0ACutHigh par %d = %g", i, fMultGlobalT0ACutHigh->GetParameter(i)); + } } if (cfgTrackCuts.cfgTPCSectorCut) { @@ -902,7 +926,7 @@ struct FlowGenericFramework { // Density dependent corrections if (cfgUseDensityDependentCorrection) { std::vector pTEffBins = {0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.4, 1.8, 2.2, 2.6, 3.0}; - hFindPtBin = new TH1D("hFindPtBin", "hFindPtBin", pTEffBins.size() - 1, &pTEffBins[0]); + hFindPtBin = new TH1D("hFindPtBin", "hFindPtBin", pTEffBins.size() - 1, pTEffBins.data()); funcEff.resize(pTEffBins.size() - 1); // LHC24g3 Eff std::vector f1p0 = cfgTrackDensityP0; @@ -920,7 +944,7 @@ struct FlowGenericFramework { } } - static constexpr std::string_view FillTimeName[] = {"before/", "after/"}; + static constexpr std::array FillTimeName = {"before/", "after/"}; void printResoCuts() { @@ -944,15 +968,17 @@ struct FlowGenericFramework { // Determine row label width size_t rowLabelWidth = 0; - for (const auto& rowName : lbl.labels_rows) + for (const auto& rowName : lbl.labels_rows) { rowLabelWidth = std::max(rowLabelWidth, rowName.size()); + } rowLabelWidth += 2; // for ": " // Print header std::ostringstream header; header << std::setw(rowLabelWidth) << " "; - for (size_t c = 0; c < nCols; ++c) + for (size_t c = 0; c < nCols; ++c) { header << std::setw(colWidths[c]) << lbl.labels_cols[c] << " "; + } LOGF(info, "%s", header.str().c_str()); // Print rows @@ -978,14 +1004,16 @@ struct FlowGenericFramework { // ----- Resonance Cuts ----- std::vector> resoCutsVals(resoCutVals.size()); - for (size_t r = 0; r < resoCutVals.size(); ++r) + for (size_t r = 0; r < resoCutVals.size(); ++r) { resoCutsVals[r] = std::vector(resoCutVals[r].begin(), resoCutVals[r].end()); + } printTable(cfgPIDCuts.resonanceCuts.value, resoCutsVals, "Resonance Cuts"); // ----- Resonance Switches ----- std::vector> resoSwitchValsF(resoSwitchVals.size()); - for (size_t r = 0; r < resoSwitchVals.size(); ++r) + for (size_t r = 0; r < resoSwitchVals.size(); ++r) { resoSwitchValsF[r] = std::vector(resoSwitchVals[r].begin(), resoSwitchVals[r].end()); + } printTable(cfgPIDCuts.resonanceSwitches.value, resoSwitchValsF, "Resonance Switches"); } enum QAFillTime { @@ -1000,7 +1028,7 @@ struct FlowGenericFramework { std::string suffix = "_ptDiff"; for (auto i = 0; i < fPtAxis->GetNbins(); ++i) { std::string index = Form("_pt_%i", i + 1); - oba->Add(new TNamed(it->Head.c_str() + index, it->Head.c_str() + suffix)); + oba->Add(new TNamed(it->Head + index, it->Head + suffix)); } } else { oba->Add(new TNamed(it->Head.c_str(), it->Head.c_str())); @@ -1028,8 +1056,9 @@ struct FlowGenericFramework { void loadCorrections(aod::BCsWithTimestamps::iterator const& bc) { uint64_t timestamp = bc.timestamp(); - if (!cfgRunByRun && cfg.correctionsLoaded) + if (!cfgRunByRun && cfg.correctionsLoaded) { return; + } if (!cfgAcceptance.value.empty()) { std::string runstr = (cfgRunByRun) ? "RunByRun/" : ""; cfg.mAcceptance.clear(); @@ -1044,14 +1073,16 @@ struct FlowGenericFramework { } } // Run-by-run efficiencies are not supported at the moment - if (cfg.correctionsLoaded) + if (cfg.correctionsLoaded) { return; + } if (!cfgEfficiency.cfgEfficiencyPath.value.empty()) { if (!cfgEfficiency.cfgUsePIDEfficiencies) { - if (cfgEfficiency.cfgUse2DEfficiency) + if (cfgEfficiency.cfgUse2DEfficiency) { cfg.mEfficiency = ccdb->getForTimeStamp(cfgEfficiency.cfgEfficiencyPath, timestamp); - else + } else { cfg.mEfficiency = ccdb->getForTimeStamp(cfgEfficiency.cfgEfficiencyPath, timestamp); + } if (cfg.mEfficiency == nullptr) { LOGF(fatal, "Could not load efficiency histogram from %s", cfgEfficiency.cfgEfficiencyPath.value.c_str()); } @@ -1059,12 +1090,14 @@ struct FlowGenericFramework { } else { std::vector species = {"ch", "pi", "ka", "pr", "k0", "lambda"}; for (const auto& sp : species) { - if (cfgEfficiency.cfgUse2DEfficiency) + if (cfgEfficiency.cfgUse2DEfficiency) { cfg.mPIDEfficiencies.push_back(ccdb->getForTimeStamp(cfgEfficiency.cfgEfficiencyPath.value + "/" + sp, timestamp)); - else + } else { cfg.mPIDEfficiencies.push_back(ccdb->getForTimeStamp(cfgEfficiency.cfgEfficiencyPath.value + "/" + sp, timestamp)); - if (cfg.mPIDEfficiencies.back() == nullptr) + } + if (cfg.mPIDEfficiencies.back() == nullptr) { LOGF(fatal, "Could not load PID efficiency histograms from %s", cfgEfficiency.cfgEfficiencyPath.value + "/" + sp); + } LOGF(info, "Loaded PID efficiency histogram from %s", cfgEfficiency.cfgEfficiencyPath.value + "/" + sp); } } @@ -1073,41 +1106,46 @@ struct FlowGenericFramework { } template - double getAcceptance(TTrack track, const double& vtxz, int index) + double getAcceptance(const TTrack& track, const double& vtxz, int index) { // 0 ref, 1 ch, 2 pi, 3 ka, 4 pr double wacc = 1; - if (!cfg.mAcceptance.empty()) + if (!cfg.mAcceptance.empty()) { wacc = cfg.mAcceptance[index]->getNUA(track.phi(), track.eta(), vtxz); + } return wacc; } template - double getEfficiency(TTrack track, const float& centrality, int pidIndex = 0) + double getEfficiency(const TTrack& track, const float& centrality, int pidIndex = 0) { //-1 ref, 0 ch, 1 pi, 2 ka, 3 pr, 4 k0, 5 lambda double eff = 1.; if (!cfgEfficiency.cfgUsePIDEfficiencies) { - if (!cfg.mEfficiency) + if (!cfg.mEfficiency) { return eff; - if (cfgEfficiency.cfgUse2DEfficiency) + } + if (cfgEfficiency.cfgUse2DEfficiency) { eff = dynamic_cast(cfg.mEfficiency)->GetBinContent(dynamic_cast(cfg.mEfficiency)->FindBin(track.pt(), centrality)); - else + } else { eff = dynamic_cast(cfg.mEfficiency)->GetBinContent(dynamic_cast(cfg.mEfficiency)->FindBin(track.pt())); + } } else { - if (cfg.mPIDEfficiencies.empty()) + if (cfg.mPIDEfficiencies.empty()) { return eff; - if (cfgEfficiency.cfgUse2DEfficiency) + } + if (cfgEfficiency.cfgUse2DEfficiency) { eff = dynamic_cast(cfg.mPIDEfficiencies[pidIndex])->GetBinContent(dynamic_cast(cfg.mPIDEfficiencies[pidIndex])->FindBin(track.pt(), centrality)); - else + } else { eff = dynamic_cast(cfg.mPIDEfficiencies[pidIndex])->GetBinContent(dynamic_cast(cfg.mPIDEfficiencies[pidIndex])->FindBin(track.pt())); + } } - if (eff == 0) + if (eff == 0) { return -1.; - else - return 1. / eff; + } + return 1. / eff; } template - int getNsigmaPID(TTrack track) + int getNsigmaPID(const TTrack& track) { // Computing Nsigma arrays for pion, kaon, and protons std::array nSigmaTPC = {track.tpcNSigmaPi(), track.tpcNSigmaKa(), track.tpcNSigmaPr()}; @@ -1117,21 +1155,24 @@ struct FlowGenericFramework { // Choose which nSigma to use std::array nSigmaToUse = (track.pt() > cfgPIDCuts.cfgTofPtCut && track.hasTOF()) ? nSigmaCombined : nSigmaTPC; - if (track.pt() > cfgPIDCuts.cfgTofPtCut && !track.hasTOF()) + if (track.pt() > cfgPIDCuts.cfgTofPtCut && !track.hasTOF()) { return 0; + } const int numSpecies = 3; int pidCount = 0; // Select particle with the lowest nsigma for (int i = 0; i < numSpecies; ++i) { if (std::abs(nSigmaToUse[i]) < nsigma) { - if (pidCount > 0 && cfgPIDCuts.cfgUseStrictPID) + if (pidCount > 0 && cfgPIDCuts.cfgUseStrictPID) { return 0; // more than one particle with low nsigma + } pidCount++; pid = i; - if (!cfgPIDCuts.cfgUseStrictPID) + if (!cfgPIDCuts.cfgUseStrictPID) { nsigma = std::abs(nSigmaToUse[i]); + } } } @@ -1139,7 +1180,7 @@ struct FlowGenericFramework { } template - int getNsigmaPIDAssymmetric(TTrack track) + int getNsigmaPIDAssymmetric(const TTrack& track) { // Computing Nsigma arrays for pion, kaon, and protons std::array nSigmaTPC = {track.tpcNSigmaPi(), track.tpcNSigmaKa(), track.tpcNSigmaPr()}; @@ -1150,7 +1191,7 @@ struct FlowGenericFramework { std::array nSigmaToUse = cfgPIDCuts.cfgUseItsPID ? nSigmaITS : nSigmaTPC; // Choose which nSigma to use: TPC or ITS std::array detectorNsigmaCut = cfgPIDCuts.cfgUseItsPID ? itsNsigmaCut : tpcNsigmaCut; // Choose which nSigma to use: TPC or ITS - bool isPion, isKaon, isProton; + bool isPion = false, isKaon = false, isProton = false; bool isDetectedPion = nSigmaToUse[0] < detectorNsigmaCut[0] && nSigmaToUse[0] > detectorNsigmaCut[0 + 3]; bool isDetectedKaon = nSigmaToUse[1] < detectorNsigmaCut[1] && nSigmaToUse[1] > detectorNsigmaCut[1 + 3]; bool isDetectedProton = nSigmaToUse[2] < detectorNsigmaCut[2] && nSigmaToUse[2] > detectorNsigmaCut[2 + 3]; @@ -1161,7 +1202,9 @@ struct FlowGenericFramework { if (track.pt() > cfgPIDCuts.cfgTofPtCut && !track.hasTOF()) { return 0; - } else if (track.pt() > cfgPIDCuts.cfgTofPtCut && track.hasTOF()) { + } + + if (track.pt() > cfgPIDCuts.cfgTofPtCut && track.hasTOF()) { isPion = isTofPion && isDetectedPion; isKaon = isTofKaon && isDetectedKaon; isProton = isTofProton && isDetectedProton; @@ -1189,7 +1232,7 @@ struct FlowGenericFramework { } template - bool eventSelected(TCollision collision, const int multTrk, const float& centrality, const int run) + bool eventSelected(const TCollision& collision, const int multTrk, const float& centrality, const int run) { // Cut on trigger alias if (cfgEventCutFlags.cfgTVXinTRD) { @@ -1199,32 +1242,38 @@ struct FlowGenericFramework { return false; } registryQA.fill(HIST("eventQA/eventSel"), TVXinTRD); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(TVXinTRD); + } } // Cut on event selection flags for (const auto& cut : eventcutflags) { - if (!cut.enabled) + if (!cut.enabled) { continue; - if (!collision.selection_bit(cut.flag)) + } + if (!collision.selection_bit(cut.flag)) { return false; + } registryQA.fill(HIST("eventQA/eventSel"), cut.histBin); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(cut.histBin); + } } // Cut on vertex - if (!selectVertex(collision)) + if (!selectVertex(collision)) { return false; + } // Cut on multiplicity correlations - data driven - if (cfgMultCut) { - if (!selectMultiplicityCorrelation(collision, multTrk, centrality, run)) + if (cfgEventSelection.cfgMultCut) { + if (!selectMultiplicityCorrelation(collision, multTrk, centrality, run)) { return false; + } } return true; } template - bool selectVertex(TCollision collision) + bool selectVertex(const TCollision& collision) { float vtxz = -999; if (collision.numContrib() > 1) { @@ -1232,81 +1281,103 @@ struct FlowGenericFramework { float zRes = std::sqrt(collision.covZZ()); float minZRes = 0.25; int minNContrib = 20; - if (zRes > minZRes && collision.numContrib() < minNContrib) + if (zRes > minZRes && collision.numContrib() < minNContrib) { vtxz = -999; + } } - if (vtxz > o2::analysis::gfw::vtxZup || vtxz < o2::analysis::gfw::vtxZlow) + if (vtxz > gfwMemberCache.vtxZup || vtxz < gfwMemberCache.vtxZlow) { return false; - else - return true; + } + + return true; } template - bool selectMultiplicityCorrelation(TCollision collision, const int multTrk, const float& centrality, const int run) + bool selectMultiplicityCorrelation(const TCollision& collision, const int multTrk, const float& centrality, const int run) { auto multNTracksPV = collision.multNTracksPV(); - if (multNTracksPV < fMultPVCutLow->Eval(centrality)) + if (multNTracksPV < fMultPVCutLow->Eval(centrality)) { return false; - if (multNTracksPV > fMultPVCutHigh->Eval(centrality)) + } + if (multNTracksPV > fMultPVCutHigh->Eval(centrality)) { return false; - if (multTrk < fMultCutLow->Eval(centrality)) + } + if (multTrk < fMultCutLow->Eval(centrality)) { return false; - if (multTrk > fMultCutHigh->Eval(centrality)) + } + if (multTrk > fMultCutHigh->Eval(centrality)) { return false; - if (multTrk > fMultPVGlobalCutHigh->Eval(collision.multNTracksPV())) + } + if (multTrk > fMultPVGlobalCutHigh->Eval(collision.multNTracksPV())) { return false; + } - if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFV0A()) < fMultGlobalV0ACutLow->Eval(multTrk)) + if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFV0A()) < fMultGlobalV0ACutLow->Eval(multTrk)) { return false; - if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFV0A()) > fMultGlobalV0ACutHigh->Eval(multTrk)) + } + if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFV0A()) > fMultGlobalV0ACutHigh->Eval(multTrk)) { return false; - if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFT0A()) < fMultGlobalT0ACutLow->Eval(multTrk)) + } + if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFT0A()) < fMultGlobalT0ACutLow->Eval(multTrk)) { return false; - if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFT0A()) > fMultGlobalT0ACutHigh->Eval(multTrk)) + } + if (!(cfgMultCorrCuts.cfgMultGlobalASideCorrCutFunction->empty()) && static_cast(collision.multFT0A()) > fMultGlobalT0ACutHigh->Eval(multTrk)) { return false; + } registryQA.fill(HIST("eventQA/eventSel"), MultCuts); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(MultCuts); + } return true; } template - bool trackSelected(TTrack track, const int field) + bool trackSelected(const TTrack& track, const int field) { if (cfgTrackCuts.cfgTPCSectorCut) { double phimodn = track.phi(); - if (field < 0) // for negative polarity field + if (field < 0) { // for negative polarity field phimodn = o2::constants::math::TwoPI - phimodn; - if (track.sign() < 0) // for negative charge + } + if (track.sign() < 0) { // for negative charge phimodn = o2::constants::math::TwoPI - phimodn; - if (phimodn < 0) + } + if (phimodn < 0) { LOGF(warning, "phi < 0: %g", phimodn); + } phimodn += o2::constants::math::PI / 18.0; // to center gap in the middle phimodn = fmod(phimodn, o2::constants::math::PI / 9.0); - if (cfgFill.cfgFillQA) + if (cfgFill.cfgFillQA) { registryQA.fill(HIST("trackQA/before/pt_phi"), track.pt(), phimodn); - if (phimodn < fPhiCutHigh->Eval(track.pt()) && phimodn > fPhiCutLow->Eval(track.pt())) + } + if (phimodn < fPhiCutHigh->Eval(track.pt()) && phimodn > fPhiCutLow->Eval(track.pt())) { return false; // reject track - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { registryQA.fill(HIST("trackQA/after/pt_phi"), track.pt(), phimodn); + } } - if (cfgTrackCuts.cfgDCAxyNSigma && (std::fabs(track.dcaXY()) > fPtDepDCAxy->Eval(track.pt()))) + if (cfgTrackCuts.cfgDCAxyNSigma && (std::fabs(track.dcaXY()) > fPtDepDCAxy->Eval(track.pt()))) { return false; - if (!cfgTrackCuts.cfgDCAzPtDep.value.empty() && std::fabs(track.dcaZ() > fPtDepDCAz->Eval(track.pt()))) + } + if (!cfgTrackCuts.cfgDCAzPtDep.value.empty() && std::fabs(track.dcaZ() > fPtDepDCAz->Eval(track.pt()))) { return false; + } return ((track.tpcNClsCrossedRows() >= cfgTrackCuts.cfgNTPCXrows) && (track.tpcNClsFound() >= cfgTrackCuts.cfgNTPCCls) && (track.itsNCls() >= cfgTrackCuts.cfgMinNITSCls)); } template - bool nchSelected(TTrack track) + bool nchSelected(const TTrack& track) { // Renormalise to default cut const float defaultNsigma = 7; - if (cfgTrackCuts.cfgDCAxyNSigma && (std::fabs(track.dcaXY()) > defaultNsigma / cfgTrackCuts.cfgDCAxyNSigma * fPtDepDCAxy->Eval(track.pt()))) + if (cfgTrackCuts.cfgDCAxyNSigma && (std::fabs(track.dcaXY()) > defaultNsigma / cfgTrackCuts.cfgDCAxyNSigma * fPtDepDCAxy->Eval(track.pt()))) { return false; - if (!cfgTrackCuts.cfgDCAzPtDep.value.empty() && std::fabs(track.dcaZ() > fPtDepDCAz->Eval(track.pt()))) + } + if (!cfgTrackCuts.cfgDCAzPtDep.value.empty() && std::fabs(track.dcaZ() > fPtDepDCAz->Eval(track.pt()))) { return false; + } int tpcNClsCrossedRowsDefault = 70; int tpcNClsFoundDefault = 50; int itsNclsDefault = 5; @@ -1319,21 +1390,24 @@ struct FlowGenericFramework { }; template - void fillWeights(const TTrack track, const double vtxz, const int pid_index, const int run) + void fillWeights(const TTrack& track, const double vtxz, const int pid_index, const int run) { if (cfgUsePID) { - double ptpidmins[] = {o2::analysis::gfw::ptpoilow, o2::analysis::gfw::ptpoilow, 0.3, 0.5}; // min pt for ch, pi, ka, pr - double ptpidmaxs[] = {o2::analysis::gfw::ptpoiup, o2::analysis::gfw::ptpoiup, 6.0, 6.0}; // max pt for ch, pi, ka, pr - bool withinPtPOI = (ptpidmins[pid_index] < track.pt()) && (track.pt() < ptpidmaxs[pid_index]); // within POI pT range - bool withinPtRef = (o2::analysis::gfw::ptreflow < track.pt()) && (track.pt() < o2::analysis::gfw::ptrefup); // within RF pT range + std::array ptpidmins = {gfwMemberCache.ptpoilow, gfwMemberCache.ptpoilow, 0.3, 0.5}; // min pt for ch, pi, ka, pr + std::array ptpidmaxs = {gfwMemberCache.ptpoiup, gfwMemberCache.ptpoiup, 6.0, 6.0}; // max pt for ch, pi, ka, pr + bool withinPtPOI = (ptpidmins[pid_index] < track.pt()) && (track.pt() < ptpidmaxs[pid_index]); // within POI pT range + bool withinPtRef = (gfwMemberCache.ptreflow < track.pt()) && (track.pt() < gfwMemberCache.ptrefup); // within RF pT range if (cfgFill.cfgFillRunByRunQA) { - if (withinPtRef && !pid_index) + if (withinPtRef && !pid_index) { th3sList[run][NUAref]->Fill(track.phi(), track.eta(), vtxz); // pt-subset of charged particles for ref flow - if (withinPtPOI) + } + if (withinPtPOI) { th3sList[run][NUAch + pid_index]->Fill(track.phi(), track.eta(), vtxz); // charged and id'ed particle weights + } } else { - if (withinPtRef && !pid_index) + if (withinPtRef && !pid_index) { registryQA.fill(HIST("phi_eta_vtxz_ref"), track.phi(), track.eta(), vtxz); // pt-subset of charged particles for ref flow + } if (withinPtPOI) { switch (pid_index) { case 0: @@ -1348,28 +1422,30 @@ struct FlowGenericFramework { case 3: registryQA.fill(HIST("phi_eta_vtxz_pr"), track.phi(), track.eta(), vtxz); // proton weights break; + default: + break; } } } } else { - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th3sList[run][NUAref]->Fill(track.phi(), track.eta(), vtxz); - else + } else { registryQA.fill(HIST("phi_eta_vtxz_ref"), track.phi(), track.eta(), vtxz); + } } - return; } void createRunByRunHistograms(const int run) { LOGF(info, "Creating histograms for run %d", run); - AxisSpec phiAxis = {o2::analysis::gfw::phibins, o2::analysis::gfw::philow, o2::analysis::gfw::phiup, "#phi"}; + AxisSpec phiAxis = {gfwMemberCache.phibins, gfwMemberCache.philow, gfwMemberCache.phiup, "#phi"}; AxisSpec phiModAxis = {100, 0, constants::math::PI / 9, "fmod(#varphi,#pi/9)"}; - AxisSpec etaAxis = {o2::analysis::gfw::etabins, cfgKinematics.cfgEta->first, cfgKinematics.cfgEta->second, "#eta"}; - AxisSpec vtxAxis = {o2::analysis::gfw::vtxZbins, -cfgVtxZ, cfgVtxZ, "Vtx_{z} (cm)"}; - AxisSpec nchAxis = {o2::analysis::gfw::nchbins, o2::analysis::gfw::nchlow, o2::analysis::gfw::nchup, "N_{ch}"}; - AxisSpec centAxis = {o2::analysis::gfw::centbinning, "Centrality (%)"}; - AxisSpec ptAxis = {o2::analysis::gfw::ptbinning, "#it{p}_{T} GeV/#it{c}"}; + AxisSpec etaAxis = {gfwMemberCache.etabins, cfgKinematics.cfgEta->first, cfgKinematics.cfgEta->second, "#eta"}; + AxisSpec vtxAxis = {gfwMemberCache.vtxZbins, -cfgVtxZ, cfgVtxZ, "Vtx_{z} (cm)"}; + AxisSpec nchAxis = {gfwMemberCache.nchbins, gfwMemberCache.nchlow, gfwMemberCache.nchup, "N_{ch}"}; + AxisSpec centAxis = {gfwMemberCache.centbinning, "Centrality (%)"}; + AxisSpec ptAxis = {gfwMemberCache.ptbinning, "#it{p}_{T} GeV/#it{c}"}; std::vector> histos(TH1NameCount); histos[Phi] = registryQA.add(Form("%d/phi", run), "", {HistType::kTH1D, {phiAxis}}); histos[Eta] = registryQA.add(Form("%d/eta", run), "", {HistType::kTH1D, {etaAxis}}); @@ -1405,7 +1481,6 @@ struct FlowGenericFramework { } histos3d[PtPhiMult] = registryQA.add(Form("%d/pt_phi_mult", run), "", {HistType::kTH3D, {ptAxis, phiModAxis, (cfgUseNch) ? nchAxis : centAxis}}); th3sList.insert(std::make_pair(run, histos3d)); - return; } struct AcceptedTracks { @@ -1429,48 +1504,61 @@ struct FlowGenericFramework { void fillNptRegistry(FractionSetup setup, const float& centmult, const NptHistos& nptHistos, const NptDenominators& dns) { - if (dns[ChargedID] <= 0) + if (dns[ChargedID] <= 0) { return; + } for (int i = 1; i <= fPtAxis->GetNbins(); ++i) { if (setup == FractionV02) { registry.fill(HIST("npt_v02_ch"), fPtAxis->GetBinCenter(i), centmult, nptHistos[ChargedID]->GetBinContent(i) / dns[ChargedID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[ChargedID] : 1.0); - if (dns[PionID] > 0) + if (dns[PionID] > 0) { registry.fill(HIST("npt_v02_pi"), fPtAxis->GetBinCenter(i), centmult, nptHistos[PionID]->GetBinContent(i) / dns[PionID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[PionID] : 1.0); - if (dns[KaonID] > 0) + } + if (dns[KaonID] > 0) { registry.fill(HIST("npt_v02_ka"), fPtAxis->GetBinCenter(i), centmult, nptHistos[KaonID]->GetBinContent(i) / dns[KaonID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[KaonID] : 1.0); - if (dns[ProtonID] > 0) + } + if (dns[ProtonID] > 0) { registry.fill(HIST("npt_v02_pr"), fPtAxis->GetBinCenter(i), centmult, nptHistos[ProtonID]->GetBinContent(i) / dns[ProtonID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[ProtonID] : 1.0); + } } else { registry.fill(HIST("npt_v0_ch"), fPtAxis->GetBinCenter(i), centmult, nptHistos[ChargedID]->GetBinContent(i) / dns[ChargedID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[ChargedID] : 1.0); - if (dns[PionID] > 0) + if (dns[PionID] > 0) { registry.fill(HIST("npt_v0_pi"), fPtAxis->GetBinCenter(i), centmult, nptHistos[PionID]->GetBinContent(i) / dns[PionID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[PionID] : 1.0); - if (dns[KaonID] > 0) + } + if (dns[KaonID] > 0) { registry.fill(HIST("npt_v0_ka"), fPtAxis->GetBinCenter(i), centmult, nptHistos[KaonID]->GetBinContent(i) / dns[KaonID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[KaonID] : 1.0); - if (dns[ProtonID] > 0) + } + if (dns[ProtonID] > 0) { registry.fill(HIST("npt_v0_pr"), fPtAxis->GetBinCenter(i), centmult, nptHistos[ProtonID]->GetBinContent(i) / dns[ProtonID], cfgEventWeight.cfgUseMultiplicityFractionWeights ? dns[ProtonID] : 1.0); + } } } } void fillNptHistos(NptHistos& nptHistos, const float pt, const int pidIndex, const double weightCh, const double weightPid) { - if (weightCh > 0) + if (weightCh > 0) { nptHistos[ChargedID]->Fill(pt, weightCh); - if (pidIndex == PionID && weightPid > 0) + } + if (pidIndex == PionID && weightPid > 0) { nptHistos[PionID]->Fill(pt, weightPid); - if (pidIndex == KaonID && weightPid > 0) + } + if (pidIndex == KaonID && weightPid > 0) { nptHistos[KaonID]->Fill(pt, weightPid); - if (pidIndex == ProtonID && weightPid > 0) + } + if (pidIndex == ProtonID && weightPid > 0) { nptHistos[ProtonID]->Fill(pt, weightPid); + } } void fillNptHistosForEta(const float eta, const float pt, const int pidIndex, const double weightCh, const double weightPid) { - if (eta > cfgKinematics.cfgEtaNptV02->first && eta < cfgKinematics.cfgEtaNptV02->second) + if (eta > cfgKinematics.cfgEtaNptV02->first && eta < cfgKinematics.cfgEtaNptV02->second) { fillNptHistos(histosNpt[FractionV02], pt, pidIndex, weightCh, weightPid); - if (eta > cfgKinematics.cfgEtaNptV0->first && eta < cfgKinematics.cfgEtaNptV0->second) + } + if (eta > cfgKinematics.cfgEtaNptV0->first && eta < cfgKinematics.cfgEtaNptV0->second) { fillNptHistos(histosNpt[FractionV0], pt, pidIndex, weightCh, weightPid); + } } template @@ -1482,32 +1570,37 @@ struct FlowGenericFramework { fFCpt->fillSubeventPtProfiles(centmult, rndm); fFCpt->fillCMProfiles(centmult, rndm); fFCpt->fillCMSubeventProfiles(centmult, rndm); - if (!cfgUseGapMethod) + if (!cfgUseGapMethod) { fFCpt->fillVnPtStdProfiles(centmult, rndm); + } for (uint l_ind = 0; l_ind < corrconfigs.size(); ++l_ind) { if (!corrconfigs.at(l_ind).pTDif) { auto dnx = fGFW->Calculate(corrconfigs.at(l_ind), 0, kTRUE).real(); - if (dnx == 0) + if (dnx == 0) { continue; + } auto val = fGFW->Calculate(corrconfigs.at(l_ind), 0, kFALSE).real() / dnx; if (std::abs(val) < 1) { - if (corrconfigs.at(l_ind).Head.find("3pcW") != std::string::npos && cfgEventWeight.cfgUseMultiplicityFractionWeights) + if (corrconfigs.at(l_ind).Head.find("3pcW") != std::string::npos && cfgEventWeight.cfgUseMultiplicityFractionWeights) { dnx *= histosNpt[FractionV02][ChargedID]->Integral(); + } (dt == Gen) ? fFCgen->FillProfile(corrconfigs.at(l_ind).Head.c_str(), centmult, val, cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm) : fFC->FillProfile(corrconfigs.at(l_ind).Head.c_str(), centmult, val, cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm); if (cfgUseGapMethod) { - fFCpt->fillVnPtProfiles(centmult, val, dnx, rndm, o2::analysis::gfw::configs.GetpTCorrMasks()[l_ind]); + fFCpt->fillVnPtProfiles(centmult, val, dnx, rndm, gfwMemberCache.configs.GetpTCorrMasks()[l_ind]); } } continue; } for (int i = 1; i <= fPtAxis->GetNbins(); i++) { auto dnx = fGFW->Calculate(corrconfigs.at(l_ind), i - 1, kTRUE).real(); - if (dnx == 0) + if (dnx == 0) { continue; + } auto val = fGFW->Calculate(corrconfigs.at(l_ind), i - 1, kFALSE).real() / dnx; - if (std::abs(val) < 1) + if (std::abs(val) < 1) { (dt == Gen) ? fFCgen->FillProfile(Form("%s_pt_%i", corrconfigs.at(l_ind).Head.c_str(), i), centmult, val, cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm) : fFC->FillProfile(Form("%s_pt_%i", corrconfigs.at(l_ind).Head.c_str(), i), centmult, val, cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm); + } } } @@ -1519,65 +1612,73 @@ struct FlowGenericFramework { fillNptRegistry(FractionV02, centmult, nptV02, dnsV02); fillNptRegistry(FractionV0, centmult, nptV0, dnsV0); - if (corrconfigsV02.size() < SpeciesCount) + if (corrconfigsV02.size() < SpeciesCount) { return; + } if (dnsV02[ChargedID] > 0) { for (uint l_ind = 0; l_ind < SpeciesCount; ++l_ind) { for (int i = 1; i <= fPtAxis->GetNbins(); i++) { auto dnx = fGFW->Calculate(corrconfigsV02.at(l_ind), i - 1, kTRUE).real(); - if (dnx == 0) + if (dnx == 0) { continue; + } auto val = fGFW->Calculate(corrconfigsV02.at(l_ind), i - 1, kFALSE).real() / dnx; if (std::abs(val) < 1 && dnsV02[l_ind] > 0) { - if (cfgEventWeight.cfgUseMultiplicityFractionWeights) + if (cfgEventWeight.cfgUseMultiplicityFractionWeights) { dnx *= dnsV02[l_ind]; + } (dt == Gen) ? fFCgen->FillProfile(Form("%s_pt_%i", corrconfigsV02.at(l_ind).Head.c_str(), i), centmult, val * nptV02[l_ind]->GetBinContent(i) / dnsV02[l_ind], cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm) : fFC->FillProfile(Form("%s_pt_%i", corrconfigsV02.at(l_ind).Head.c_str(), i), centmult, val * nptV02[l_ind]->GetBinContent(i) / dnsV02[l_ind], cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm); } } } } - if (corrconfigsV0.size() < SpeciesCount) + if (corrconfigsV0.size() < SpeciesCount) { return; + } double mpt = 0; double dnx = 0; if (cfgKinematics.cfgEtaPtPt->first * cfgKinematics.cfgEtaPtPt->second >= 0) { - if (fFCpt->corrDen[1] == 0.) + if (fFCpt->corrDen[1] == 0.) { return; + } dnx = fFCpt->corrDen[1]; mpt = fFCpt->corrNum[1] / dnx; } else { - if (fFCpt->corrDenSub[0][1] == 0. || fFCpt->corrDenSub[1][1] == 0.) + if (fFCpt->corrDenSub[0][1] == 0. || fFCpt->corrDenSub[1][1] == 0.) { return; + } double mptSub1 = fFCpt->corrNumSub[0][1] / fFCpt->corrDenSub[0][1]; double mptSub2 = fFCpt->corrNumSub[1][1] / fFCpt->corrDenSub[1][1]; dnx = 0.5 * (fFCpt->corrDenSub[0][1] + fFCpt->corrDenSub[1][1]); mpt = 0.5 * (mptSub1 + mptSub2); } - if (std::isnan(mpt)) + if (std::isnan(mpt)) { return; + } for (uint l_ind = 0; l_ind < SpeciesCount; ++l_ind) { for (int i = 1; i <= fPtAxis->GetNbins(); i++) { if (dnsV0[l_ind] > 0) { double profileWeight = dnx; - if (cfgEventWeight.cfgUseMultiplicityFractionWeights) + if (cfgEventWeight.cfgUseMultiplicityFractionWeights) { profileWeight *= dnsV0[l_ind]; + } (dt == Gen) ? fFCgen->FillProfile(Form("%s_pt_%i", corrconfigsV0.at(l_ind).Head.c_str(), i), centmult, mpt * nptV0[l_ind]->GetBinContent(i) / dnsV0[l_ind], cfgEventWeight.cfgUseMultiplicityFlowWeights ? profileWeight : 1.0, rndm) : fFC->FillProfile(Form("%s_pt_%i", corrconfigsV0.at(l_ind).Head.c_str(), i), centmult, mpt * nptV0[l_ind]->GetBinContent(i) / dnsV0[l_ind], cfgEventWeight.cfgUseMultiplicityFlowWeights ? profileWeight : 1.0, rndm); } } } - return; } template void fillResonanceOutput(FractionSetup setup, const float& centmult, const double& rndm) { if (setup == FractionV02) { - if (histosNpt[FractionV02][ChargedID]->Integral() <= 0) + if (histosNpt[FractionV02][ChargedID]->Integral() <= 0) { return; + } double dnK0SB1 = histosNpt[FractionV02][ChargedID]->Integral(); double dnK0Sig = histosNpt[FractionV02][ChargedID]->Integral(); @@ -1596,18 +1697,24 @@ struct FlowGenericFramework { } for (int i = 1; i <= fPtAxis->GetNbins(); ++i) { - if (dnK0SB1 > 0) + if (dnK0SB1 > 0) { registry.fill(HIST("npt_v02_K0_sb1"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV02][K0Sideband1]->GetBinContent(i) / dnK0SB1, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnK0SB1 : 1.); - if (dnK0Sig > 0) + } + if (dnK0Sig > 0) { registry.fill(HIST("npt_v02_K0_sig"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV02][K0Signal]->GetBinContent(i) / dnK0Sig, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnK0Sig : 1.); - if (dnK0SB2 > 0) + } + if (dnK0SB2 > 0) { registry.fill(HIST("npt_v02_K0_sb2"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV02][K0Sideband2]->GetBinContent(i) / dnK0SB2, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnK0SB2 : 1.); - if (dnLambdaSB1 > 0) + } + if (dnLambdaSB1 > 0) { registry.fill(HIST("npt_v02_Lambda_sb1"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV02][LambdaSideband1]->GetBinContent(i) / dnLambdaSB1, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnLambdaSB1 : 1.); - if (dnLambdaSig > 0) + } + if (dnLambdaSig > 0) { registry.fill(HIST("npt_v02_Lambda_sig"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV02][LambdaSignal]->GetBinContent(i) / dnLambdaSig, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnLambdaSig : 1.); - if (dnLambdaSB2 > 0) + } + if (dnLambdaSB2 > 0) { registry.fill(HIST("npt_v02_Lambda_sb2"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV02][LambdaSideband2]->GetBinContent(i) / dnLambdaSB2, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnLambdaSB2 : 1.); + } } std::vector dns = {dnK0SB1, dnK0Sig, dnK0SB2, dnLambdaSB1, dnLambdaSig, dnLambdaSB2}; @@ -1615,19 +1722,22 @@ struct FlowGenericFramework { for (uint l_ind = 4; l_ind < corrconfigsV02.size(); ++l_ind) { for (int i = 1; i <= fPtAxis->GetNbins(); i++) { auto dnx = fGFW->Calculate(corrconfigsV02.at(l_ind), i - 1, kTRUE).real(); - if (dnx == 0) + if (dnx == 0) { continue; + } auto val = fGFW->Calculate(corrconfigsV02.at(l_ind), i - 1, kFALSE).real() / dnx; if (std::abs(val) < 1 && dns[l_ind - 4] > 0) { - if (cfgEventWeight.cfgUseMultiplicityFractionWeights) + if (cfgEventWeight.cfgUseMultiplicityFractionWeights) { dnx *= dns[l_ind - 4]; + } (dt == Gen) ? fFCgen->FillProfile(Form("%s_pt_%i", corrconfigsV02.at(l_ind).Head.c_str(), i), centmult, val * histosResoNpt[FractionV02][l_ind - 4]->GetBinContent(i) / dns[l_ind - 4], cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm) : fFC->FillProfile(Form("%s_pt_%i", corrconfigsV02.at(l_ind).Head.c_str(), i), centmult, val * histosResoNpt[FractionV02][l_ind - 4]->GetBinContent(i) / dns[l_ind - 4], cfgEventWeight.cfgUseMultiplicityFlowWeights ? dnx : 1.0, rndm); } } } } else { - if (histosNpt[FractionV0][ChargedID]->Integral() <= 0) + if (histosNpt[FractionV0][ChargedID]->Integral() <= 0) { return; + } double dnK0SB1 = histosNpt[FractionV0][ChargedID]->Integral(); double dnK0Sig = histosNpt[FractionV0][ChargedID]->Integral(); @@ -1646,49 +1756,61 @@ struct FlowGenericFramework { } for (int i = 1; i <= fPtAxis->GetNbins(); ++i) { - if (dnK0SB1 > 0) + if (dnK0SB1 > 0) { registry.fill(HIST("npt_v0_K0_sb1"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV0][K0Sideband1]->GetBinContent(i) / dnK0SB1, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnK0SB1 : 1.); - if (dnK0Sig > 0) + } + if (dnK0Sig > 0) { registry.fill(HIST("npt_v0_K0_sig"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV0][K0Signal]->GetBinContent(i) / dnK0Sig, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnK0Sig : 1.); - if (dnK0SB2 > 0) + } + if (dnK0SB2 > 0) { registry.fill(HIST("npt_v0_K0_sb2"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV0][K0Sideband2]->GetBinContent(i) / dnK0SB2, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnK0SB2 : 1.); - if (dnLambdaSB1 > 0) + } + if (dnLambdaSB1 > 0) { registry.fill(HIST("npt_v0_Lambda_sb1"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV0][LambdaSideband1]->GetBinContent(i) / dnLambdaSB1, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnLambdaSB1 : 1.); - if (dnLambdaSig > 0) + } + if (dnLambdaSig > 0) { registry.fill(HIST("npt_v0_Lambda_sig"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV0][LambdaSignal]->GetBinContent(i) / dnLambdaSig, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnLambdaSig : 1.); - if (dnLambdaSB2 > 0) + } + if (dnLambdaSB2 > 0) { registry.fill(HIST("npt_v0_Lambda_sb2"), fPtAxis->GetBinCenter(i), centmult, histosResoNpt[FractionV0][LambdaSideband2]->GetBinContent(i) / dnLambdaSB2, cfgEventWeight.cfgUseMultiplicityFractionWeights ? dnLambdaSB2 : 1.); + } } std::vector dns = {dnK0SB1, dnK0Sig, dnK0SB2, dnLambdaSB1, dnLambdaSig, dnLambdaSB2}; - if (fFCpt->corrDenSub[0][1] == 0. || fFCpt->corrDenSub[1][1] == 0.) + if (fFCpt->corrDenSub[0][1] == 0. || fFCpt->corrDenSub[1][1] == 0.) { return; + } double mpt = 0; double dnx = 0; if (cfgKinematics.cfgEtaPtPt->first * cfgKinematics.cfgEtaPtPt->second >= 0) { - if (fFCpt->corrDen[1] == 0.) + if (fFCpt->corrDen[1] == 0.) { return; + } dnx = fFCpt->corrDen[1]; mpt = fFCpt->corrNum[1] / dnx; } else { - if (fFCpt->corrDenSub[0][1] == 0. || fFCpt->corrDenSub[1][1] == 0.) + if (fFCpt->corrDenSub[0][1] == 0. || fFCpt->corrDenSub[1][1] == 0.) { return; + } double mptSub1 = fFCpt->corrNumSub[0][1] / fFCpt->corrDenSub[0][1]; double mptSub2 = fFCpt->corrNumSub[1][1] / fFCpt->corrDenSub[1][1]; dnx = 0.5 * (fFCpt->corrDenSub[0][1] + fFCpt->corrDenSub[1][1]); mpt = 0.5 * (mptSub1 + mptSub2); } - if (std::isnan(mpt)) + if (std::isnan(mpt)) { return; + } for (uint l_ind = 4; l_ind < corrconfigsV0.size(); ++l_ind) { for (int i = 1; i <= fPtAxis->GetNbins(); i++) { - if (dns[l_ind - 4] > 0) - if (cfgEventWeight.cfgUseMultiplicityFractionWeights) + if (dns[l_ind - 4] > 0) { + if (cfgEventWeight.cfgUseMultiplicityFractionWeights) { dnx *= dns[l_ind - 4]; + } + } (dt == Gen) ? fFCgen->FillProfile(Form("%s_pt_%i", corrconfigsV0.at(l_ind).Head.c_str(), i), centmult, mpt * histosResoNpt[FractionV0][l_ind - 4]->GetBinContent(i) / dns[l_ind - 4], 1.0, rndm) : fFC->FillProfile(Form("%s_pt_%i", corrconfigsV0.at(l_ind).Head.c_str(), i), centmult, mpt * histosResoNpt[FractionV0][l_ind - 4]->GetBinContent(i) / dns[l_ind - 4], 1.0, rndm); } } @@ -1696,16 +1818,19 @@ struct FlowGenericFramework { } template - void processCollision(TCollision collision, TTracks tracks, TV0s v0s, const float& centrality, const int field, const int run) + void processCollision(const TCollision& collision, const TTracks& tracks, const TV0s& v0s, const float& centrality, const int field, const int run) { - if (tracks.size() < 1) + if (tracks.size() < 1) { return; - if (dt != Gen && (centrality < o2::analysis::gfw::centbinning.front() || centrality > o2::analysis::gfw::centbinning.back())) + } + if (dt != Gen && (centrality < gfwMemberCache.centbinning.front() || centrality > gfwMemberCache.centbinning.back())) { return; + } if (dt != Gen) { registryQA.fill(HIST("eventQA/eventSel"), TrackCent); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(TrackCent); + } } float vtxz = collision.posZ(); if (dt != Gen && cfgFill.cfgFillRunByRunQA) { @@ -1727,7 +1852,7 @@ struct FlowGenericFramework { double q3x = 0, q3y = 0; double q4x = 0, q4y = 0; for (const auto& track : tracks) { - bool withinPtRef = (o2::analysis::gfw::ptreflow < track.pt()) && (track.pt() < o2::analysis::gfw::ptrefup); // within RF pT rang + bool withinPtRef = (gfwMemberCache.ptreflow < track.pt()) && (track.pt() < gfwMemberCache.ptrefup); // within RF pT rang if (withinPtRef) { q2x += std::cos(2 * track.phi()); q2y += std::sin(2 * track.phi()); @@ -1755,8 +1880,9 @@ struct FlowGenericFramework { AcceptedTracks acceptedTracks; // Reset fraction histograms per event for (const auto& vec : histosNpt) { - for (const auto& h : vec) + for (const auto& h : vec) { h->Reset("ICESM"); + } } for (const auto& track : tracks) { processTrack(track, vtxz, field, run, densitycorrections, centrality, acceptedTracks); @@ -1766,7 +1892,7 @@ struct FlowGenericFramework { registryQA.fill(HIST("trackQA/after/Nch_uncorrected"), acceptedTracks.totaluncorr); } - int multiplicity = 0; + float multiplicity = 0.f; switch (cfgUseNchCorrection) { case 0: multiplicity = tracks.size(); @@ -1782,23 +1908,28 @@ struct FlowGenericFramework { break; } - if (cfgFill.cfgFillWeights) + if (cfgFill.cfgFillWeights) { return; + } fillOutputContainers
((cfgUseNch) ? multiplicity : centrality, lRandom); // Reset fraction histograms per event - for (const auto& vec : histosResoNpt) - for (const auto& h : vec) + for (const auto& vec : histosResoNpt) { + for (const auto& h : vec) { h->Reset("ICESM"); + } + } // Process V0s only for reconstructed-track workflows. if constexpr (dt != Gen) { - for (const auto& v0 : v0s) + for (const auto& v0 : v0s) { processV0(v0, tracks, collision, centrality); + } } else { - for (const auto& particle : tracks) + for (const auto& particle : tracks) { processV0(particle, tracks, collision, centrality); + } } for (auto setup = 0; setup < FractionSetupCount; ++setup) { @@ -1808,38 +1939,48 @@ struct FlowGenericFramework { } template - inline void processTrack(TTrack const& track, const float& vtxz, const int field, const int run, DensityCorr densitycorrections, const float& centrality, AcceptedTracks& acceptedTracks) + inline void processTrack(const TTrack& track, const float& vtxz, const int field, const int run, const DensityCorr& densitycorrections, const float& centrality, AcceptedTracks& acceptedTracks) { if constexpr (framework::has_type_v) { - if (track.mcParticleId() < 0 || !(track.has_mcParticle())) + if (track.mcParticleId() < 0 || !(track.has_mcParticle())) { return; + } auto mcParticle = track.mcParticle(); - if (!mcParticle.isPhysicalPrimary()) + if (!mcParticle.isPhysicalPrimary()) { return; - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { fillTrackQA(track, vtxz); - if (mcParticle.eta() < o2::analysis::gfw::etalow || mcParticle.eta() > o2::analysis::gfw::etaup || mcParticle.pt() < o2::analysis::gfw::ptlow || mcParticle.pt() > o2::analysis::gfw::ptup) + } + if (mcParticle.eta() < gfwMemberCache.etalow || mcParticle.eta() > gfwMemberCache.etaup || mcParticle.pt() < gfwMemberCache.ptlow || mcParticle.pt() > gfwMemberCache.ptup) { return; + } // Select tracks with nominal cuts always - if (!nchSelected(track)) + if (!nchSelected(track)) { return; + } double weffCh = getEfficiency(track, centrality, 0); if (track.eta() > cfgKinematics.cfgEtaNch->first && track.eta() < cfgKinematics.cfgEtaNch->second) { - if (weffCh > 0) + if (weffCh > 0) { acceptedTracks.total += (cfgUseNchCorrection) ? weffCh : 1.0; + } ++acceptedTracks.totaluncorr; } - if (!trackSelected(track, field)) + if (!trackSelected(track, field)) { return; + } int pidIndex = 0; - if (std::abs(mcParticle.pdgCode()) == kPiPlus) + if (std::abs(mcParticle.pdgCode()) == kPiPlus) { pidIndex = PionID; - if (std::abs(mcParticle.pdgCode()) == kKPlus) + } + if (std::abs(mcParticle.pdgCode()) == kKPlus) { pidIndex = KaonID; - if (std::abs(mcParticle.pdgCode()) == kProton) + } + if (std::abs(mcParticle.pdgCode()) == kProton) { pidIndex = ProtonID; + } double weff = getEfficiency(track, centrality, pidIndex); double nptWeightCh = (weffCh > 0) ? ((cfgUseNchCorrection) ? weffCh : 1.0) : -1.0; double nptWeightPid = (weff > 0) ? ((cfgUseNchCorrection) ? weff : 1.0) : -1.0; @@ -1859,20 +2000,26 @@ struct FlowGenericFramework { } } else if constexpr (framework::has_type_v) { - if (!track.isPhysicalPrimary()) + if (!track.isPhysicalPrimary()) { return; - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { fillTrackQA(track, vtxz); - if (track.eta() < o2::analysis::gfw::etalow || track.eta() > o2::analysis::gfw::etaup || track.pt() < o2::analysis::gfw::ptlow || track.pt() > o2::analysis::gfw::ptup) + } + if (track.eta() < gfwMemberCache.etalow || track.eta() > gfwMemberCache.etaup || track.pt() < gfwMemberCache.ptlow || track.pt() > gfwMemberCache.ptup) { return; + } int pidIndex = 0; - if (std::abs(track.pdgCode()) == kPiPlus) + if (std::abs(track.pdgCode()) == kPiPlus) { pidIndex = PionID; - if (std::abs(track.pdgCode()) == kKPlus) + } + if (std::abs(track.pdgCode()) == kKPlus) { pidIndex = KaonID; - if (std::abs(track.pdgCode()) == kProton) + } + if (std::abs(track.pdgCode()) == kProton) { pidIndex = ProtonID; + } if (track.eta() > cfgKinematics.cfgEtaNch->first && track.eta() < cfgKinematics.cfgEtaNch->second) { ++acceptedTracks.total; @@ -1882,23 +2029,28 @@ struct FlowGenericFramework { fillPtSums(track, centrality, vtxz); fillGFW(track, centrality, vtxz, pidIndex, densitycorrections); - if (cfgFill.cfgFillQA) + if (cfgFill.cfgFillQA) { fillTrackQA(track, vtxz); + } } else { - if (cfgFill.cfgFillQA) + if (cfgFill.cfgFillQA) { fillTrackQA(track, vtxz); + } // Select tracks with nominal cuts always - if (!nchSelected(track)) + if (!nchSelected(track)) { return; + } double weffCh = getEfficiency(track, centrality, 0); if (track.eta() > cfgKinematics.cfgEtaNch->first && track.eta() < cfgKinematics.cfgEtaNch->second) { - if (weffCh > 0) + if (weffCh > 0) { acceptedTracks.total += (cfgUseNchCorrection) ? weffCh : 1.0; + } ++acceptedTracks.totaluncorr; } - if (!trackSelected(track, field)) + if (!trackSelected(track, field)) { return; + } // int pidIndex = 0; // if (cfgUsePID) Need PID for v02 int pidIndex = getNsigmaPID(track); @@ -1924,108 +2076,128 @@ struct FlowGenericFramework { } template - inline void processV0(TV0 v0, TTracks tracks, TCollision collision, const float& centrality) + inline void processV0(const TV0& v0, const TTracks& tracks, const TCollision& collision, const float& centrality) { if constexpr (dt == Reco) { using V0TrackTable = std::decay_t; auto postrack = v0.template posTrack_as(); auto negtrack = v0.template negTrack_as(); - if (resoSwitchVals[UseParticle][K0]) { + if (resoSwitchVals[UseParticle][K0] != 0) { double weff = 1; bool fillK0 = true; if (cfgEfficiency.cfgUsePIDEfficiencies) { weff = getEfficiency(v0, centrality, 4); - if (weff < 0) + if (weff < 0) { fillK0 = false; + } } if (fillK0) { auto k0Selection = selectK0(collision, v0, tracks); if (k0Selection.selected) { for (auto setup = 0; setup < FractionSetupCount; ++setup) { auto fractionSetup = static_cast(setup); - if (!selectionV0DaughterEta(postrack, fractionSetup) || !selectionV0DaughterEta(negtrack, fractionSetup)) + if (!selectionV0DaughterEta(postrack, fractionSetup) || !selectionV0DaughterEta(negtrack, fractionSetup)) { continue; + } - if (cfgFill.cfgFillV0QA && fractionSetup == FractionV02) + if (cfgFill.cfgFillV0QA && fractionSetup == FractionV02) { fillV0QA(k0Selection, v0, postrack, negtrack, centrality, weff); + } registryQA.fill(HIST("K0/hK0Count"), (fractionSetup == FractionV02) ? FillV02DaughterTrackSelection : FillV0DaughterTrackSelection); - if (v0.mK0Short() > cfgPIDCuts.cfgK0SideBand1Min && v0.mK0Short() < cfgPIDCuts.cfgK0SideBand1Max) + if (v0.mK0Short() > cfgPIDCuts.cfgK0SideBand1Min && v0.mK0Short() < cfgPIDCuts.cfgK0SideBand1Max) { histosResoNpt[fractionSetup][K0Sideband1]->Fill(v0.pt(), (cfgUseNchCorrection) ? weff : 1.0); - if (v0.mK0Short() > cfgPIDCuts.cfgK0SignalMin && v0.mK0Short() < cfgPIDCuts.cfgK0SignalMax) + } + if (v0.mK0Short() > cfgPIDCuts.cfgK0SignalMin && v0.mK0Short() < cfgPIDCuts.cfgK0SignalMax) { histosResoNpt[fractionSetup][K0Signal]->Fill(v0.pt(), (cfgUseNchCorrection) ? weff : 1.0); - if (v0.mK0Short() > cfgPIDCuts.cfgK0SideBand2Min && v0.mK0Short() < cfgPIDCuts.cfgK0SideBand2Max) + } + if (v0.mK0Short() > cfgPIDCuts.cfgK0SideBand2Min && v0.mK0Short() < cfgPIDCuts.cfgK0SideBand2Max) { histosResoNpt[fractionSetup][K0Sideband2]->Fill(v0.pt(), (cfgUseNchCorrection) ? weff : 1.0); + } } } } } - if (resoSwitchVals[UseParticle][Lambda]) { + if (resoSwitchVals[UseParticle][Lambda] != 0) { double weff = 1.; bool fillLambda = true; if (cfgEfficiency.cfgUsePIDEfficiencies) { weff = getEfficiency(v0, centrality, 5); - if (weff < 0) + if (weff < 0) { fillLambda = false; + } } if (fillLambda) { auto lambdaSelection = selectLambda(collision, v0, tracks); if (lambdaSelection.selected) { for (auto setup = 0; setup < FractionSetupCount; ++setup) { auto fractionSetup = static_cast(setup); - if (!selectionLambdaDaughterEta(postrack, negtrack, lambdaSelection, fractionSetup)) + if (!selectionLambdaDaughterEta(postrack, negtrack, lambdaSelection, fractionSetup)) { continue; + } registryQA.fill(HIST("Lambda/hLambdaCount"), (fractionSetup == FractionV02) ? FillV02DaughterTrackSelection : FillV0DaughterTrackSelection); - if (cfgFill.cfgFillV0QA && fractionSetup == FractionV02) + if (cfgFill.cfgFillV0QA && fractionSetup == FractionV02) { fillV0QA(lambdaSelection, v0, postrack, negtrack, centrality, weff); + } - if (v0.mLambda() > cfgPIDCuts.cfgLambdaSideBand1Min && v0.mLambda() < cfgPIDCuts.cfgLambdaSideBand1Max) + if (v0.mLambda() > cfgPIDCuts.cfgLambdaSideBand1Min && v0.mLambda() < cfgPIDCuts.cfgLambdaSideBand1Max) { histosResoNpt[fractionSetup][LambdaSideband1]->Fill(v0.pt(), (cfgUseNchCorrection) ? weff : 1.0); - if (v0.mLambda() > cfgPIDCuts.cfgLambdaSignalMin && v0.mLambda() < cfgPIDCuts.cfgLambdaSignalMax) + } + if (v0.mLambda() > cfgPIDCuts.cfgLambdaSignalMin && v0.mLambda() < cfgPIDCuts.cfgLambdaSignalMax) { histosResoNpt[fractionSetup][LambdaSignal]->Fill(v0.pt(), (cfgUseNchCorrection) ? weff : 1.0); - if (v0.mLambda() > cfgPIDCuts.cfgLambdaSideBand2Min && v0.mLambda() < cfgPIDCuts.cfgLambdaSideBand2Max) + } + if (v0.mLambda() > cfgPIDCuts.cfgLambdaSideBand2Min && v0.mLambda() < cfgPIDCuts.cfgLambdaSideBand2Max) { histosResoNpt[fractionSetup][LambdaSideband2]->Fill(v0.pt(), (cfgUseNchCorrection) ? weff : 1.0); + } } } } } } else { - if (!v0.isPhysicalPrimary()) + if (!v0.isPhysicalPrimary()) { return; - if (!v0.has_daughters()) + } + if (!v0.has_daughters()) { return; + } bool isK0 = (v0.pdgCode() == PDG_t::kK0Short); bool isLambda = (v0.pdgCode() == PDG_t::kLambda0); - if (!isK0 && !isLambda) + if (!isK0 && !isLambda) { return; + } // To match reconstructed - check that daughters are within v02/v0 eta acceptance bool isWithinV0Acceptance = true; bool isWithinV02Acceptance = true; for (const auto& d : v0.template daughters_as()) { - if (d.eta() < cfgKinematics.cfgEtaNptV02->first || d.eta() > cfgKinematics.cfgEtaNptV02->second) + if (d.eta() < cfgKinematics.cfgEtaNptV02->first || d.eta() > cfgKinematics.cfgEtaNptV02->second) { isWithinV02Acceptance = false; + } if (d.eta() < cfgKinematics.cfgEtaNptV0->first || d.eta() > cfgKinematics.cfgEtaNptV0->second) { isWithinV0Acceptance = false; } } if (isK0) { - if (isWithinV02Acceptance) + if (isWithinV02Acceptance) { histosResoNpt[FractionV02][K0Signal]->Fill(v0.pt(), 1.0); - if (isWithinV0Acceptance) + } + if (isWithinV0Acceptance) { histosResoNpt[FractionV0][K0Signal]->Fill(v0.pt(), 1.0); + } } if (isLambda) { - if (isWithinV02Acceptance) + if (isWithinV02Acceptance) { histosResoNpt[FractionV02][LambdaSignal]->Fill(v0.pt(), 1.0); - if (isWithinV0Acceptance) + } + if (isWithinV0Acceptance) { histosResoNpt[FractionV0][LambdaSignal]->Fill(v0.pt(), 1.0); + } } } } @@ -2043,61 +2215,74 @@ struct FlowGenericFramework { }; template - bool selectionV0Daughter(TTrack const& track, int pid) + bool selectionV0Daughter(const TTrack& track, int pid) { - if (!(track.itsNCls() > cfgTrackCuts.cfgMinNITSCls)) + if (!(track.itsNCls() > cfgTrackCuts.cfgMinNITSCls)) { return 0; - if (!track.hasTPC()) + } + if (!track.hasTPC()) { return false; - if (track.tpcNClsFound() < cfgTrackCuts.cfgNTPCCls) + } + if (track.tpcNClsFound() < cfgTrackCuts.cfgNTPCCls) { return false; - if (!(track.tpcNClsCrossedRows() > cfgTrackCuts.cfgNTPCXrows)) + } + if (!(track.tpcNClsCrossedRows() > cfgTrackCuts.cfgNTPCXrows)) { return 0; + } if (cfgPIDCuts.cfgUseOnlyTPC) { - if (pid == Pions && std::abs(track.tpcNSigmaPi()) > cfgPIDCuts.cfgTPCNsigmaCut) + if (pid == Pions && std::abs(track.tpcNSigmaPi()) > cfgPIDCuts.cfgTPCNsigmaCut) { return false; - if (pid == Kaons && std::abs(track.tpcNSigmaKa()) > cfgPIDCuts.cfgTPCNsigmaCut) + } + if (pid == Kaons && std::abs(track.tpcNSigmaKa()) > cfgPIDCuts.cfgTPCNsigmaCut) { return false; - if (pid == Protons && std::abs(track.tpcNSigmaPr()) > cfgPIDCuts.cfgTPCNsigmaCut) + } + if (pid == Protons && std::abs(track.tpcNSigmaPr()) > cfgPIDCuts.cfgTPCNsigmaCut) { return false; + } } else { int partIndex = cfgPIDCuts.cfgUseAsymmetricPID ? getNsigmaPIDAssymmetric(track) : getNsigmaPID(track); int pidIndex = partIndex - 1; // 0 = pion, 1 = kaon, 2 = proton - if (pidIndex != pid) + if (pidIndex != pid) { return false; + } } return true; } template - bool selectionV0DaughterEta(TTrack const& track, FractionSetup setup) + bool selectionV0DaughterEta(const TTrack& track, FractionSetup setup) { const auto [etaMin, etaMax] = getFractionEtaRange(setup); - if (track.eta() < etaMin || track.eta() > etaMax) + if (track.eta() < etaMin || track.eta() > etaMax) { return false; + } if (cfgFill.cfgFillV0QA) { - if (setup == FractionV02) + if (setup == FractionV02) { registryQA.fill(HIST("trackQA/after/etaV02"), track.eta()); - if (setup == FractionV0) + } + if (setup == FractionV0) { registryQA.fill(HIST("trackQA/after/etaV0"), track.eta()); + } } return true; } template - bool selectionLambdaDaughterEta(TTrack const& postrack, TTrack const& negtrack, const V0Selection& lambdaSelection, FractionSetup setup) + bool selectionLambdaDaughterEta(const TTrack& postrack, const TTrack& negtrack, const V0Selection& lambdaSelection, FractionSetup setup) { - if (lambdaSelection.isL && (!selectionV0DaughterEta(postrack, setup) || !selectionV0DaughterEta(negtrack, setup))) + if (lambdaSelection.isL && (!selectionV0DaughterEta(postrack, setup) || !selectionV0DaughterEta(negtrack, setup))) { return false; - if (lambdaSelection.isAL && (!selectionV0DaughterEta(postrack, setup) || !selectionV0DaughterEta(negtrack, setup))) + } + if (lambdaSelection.isAL && (!selectionV0DaughterEta(postrack, setup) || !selectionV0DaughterEta(negtrack, setup))) { return false; + } return true; } template - V0Selection selectK0(TCollision const& collision, TV0 const& v0, TTracks const&) + V0Selection selectK0(const TCollision& collision, const TV0& v0, const TTracks&) { using V0TrackTable = std::decay_t; @@ -2109,48 +2294,60 @@ struct FlowGenericFramework { auto negtrack = v0.template negTrack_as(); registryQA.fill(HIST("K0/hK0Count"), FillCandidate); - if (postrack.pt() < resoCutVals[PosTrackPt][K0] || negtrack.pt() < resoCutVals[NegTrackPt][K0]) + if (postrack.pt() < resoCutVals[PosTrackPt][K0] || negtrack.pt() < resoCutVals[NegTrackPt][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillDaughterPt); - if (massK0s < resoCutVals[MassMin][K0] && massK0s > resoCutVals[MassMax][K0]) + if (massK0s < resoCutVals[MassMin][K0] && massK0s > resoCutVals[MassMax][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillMassCut); // Rapidity correction - if (v0.yK0Short() > resoCutVals[Rapidity][K0]) + if (v0.yK0Short() > resoCutVals[Rapidity][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillRapidityCut); // DCA cuts for K0short - if (std::abs(v0.dcapostopv()) < resoCutVals[DCAPosToPVMin][K0] || std::abs(v0.dcanegtopv()) < resoCutVals[DCANegToPVMin][K0]) + if (std::abs(v0.dcapostopv()) < resoCutVals[DCAPosToPVMin][K0] || std::abs(v0.dcanegtopv()) < resoCutVals[DCANegToPVMin][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillDCAtoPV); - if (resoSwitchVals[UseDCAxDaughters][K0] && std::abs(v0.dcaV0daughters()) > resoCutVals[DCAxDaughters][K0]) + if (resoSwitchVals[UseDCAxDaughters][K0] != 0 && std::abs(v0.dcaV0daughters()) > resoCutVals[DCAxDaughters][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillDCAxDaughters); // v0 radius cuts - if (resoSwitchVals[UseV0Radius][K0] && (v0.v0radius() < resoCutVals[RadiusMin][K0] || v0.v0radius() > resoCutVals[RadiusMax][K0])) + if (resoSwitchVals[UseV0Radius][K0] != 0 && (v0.v0radius() < resoCutVals[RadiusMin][K0] || v0.v0radius() > resoCutVals[RadiusMax][K0])) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillV0Radius); // cosine pointing angle cuts - if (v0.v0cosPA() < resoCutVals[CosPA][K0]) + if (v0.v0cosPA() < resoCutVals[CosPA][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillCosPA); // Proper lifetime - if (resoSwitchVals[UseProperLifetime][K0] && v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * massK0Short > resoCutVals[LifeTime][K0]) + if (resoSwitchVals[UseProperLifetime][K0] != 0 && v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * massK0Short > resoCutVals[LifeTime][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillProperLifetime); // ArmenterosPodolanskiCut - if (resoSwitchVals[UseArmPodCut][K0] && (v0.qtarm() / std::abs(v0.alpha())) < resoCutVals[ArmPodMin][K0]) + if (resoSwitchVals[UseArmPodCut][K0] != 0 && (v0.qtarm() / std::abs(v0.alpha())) < resoCutVals[ArmPodMin][K0]) { return selection; + } registryQA.fill(HIST("K0/hK0Count"), FillArmPodCut); - if (resoSwitchVals[UseCompetingMassRejection][K0]) { - if (std::abs(v0.mLambda() - o2::constants::physics::MassLambda0) < resoCutVals[MassRejection][K0]) + if (resoSwitchVals[UseCompetingMassRejection][K0] != 0) { + if (std::abs(v0.mLambda() - o2::constants::physics::MassLambda0) < resoCutVals[MassRejection][K0]) { return selection; - if (std::abs(v0.mAntiLambda() - o2::constants::physics::MassLambda0) < resoCutVals[MassRejection][K0]) + } + if (std::abs(v0.mAntiLambda() - o2::constants::physics::MassLambda0) < resoCutVals[MassRejection][K0]) { return selection; + } } registryQA.fill(HIST("K0/hK0Count"), FillCompetingMass); - if (!selectionV0Daughter(postrack, Pions) || !selectionV0Daughter(negtrack, Pions)) + if (!selectionV0Daughter(postrack, Pions) || !selectionV0Daughter(negtrack, Pions)) { return selection; + } selection.selected = true; selection.isK0 = true; @@ -2159,7 +2356,7 @@ struct FlowGenericFramework { } template - V0Selection selectLambda(TCollision const& collision, TV0 const& v0, TTracks const&) + V0Selection selectLambda(const TCollision& collision, const TV0& v0, const TTracks&) { using V0TrackTable = std::decay_t; V0Selection selection; @@ -2172,14 +2369,17 @@ struct FlowGenericFramework { auto negtrack = v0.template negTrack_as(); registryQA.fill(HIST("Lambda/hLambdaCount"), FillCandidate); - if (postrack.pt() < resoCutVals[PosTrackPt][Lambda] || negtrack.pt() < resoCutVals[NegTrackPt][Lambda]) + if (postrack.pt() < resoCutVals[PosTrackPt][Lambda] || negtrack.pt() < resoCutVals[NegTrackPt][Lambda]) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillDaughterPt); - if (mlambda > resoCutVals[MassMin][Lambda] && mlambda < resoCutVals[MassMax][Lambda]) + if (mlambda > resoCutVals[MassMin][Lambda] && mlambda < resoCutVals[MassMax][Lambda]) { selection.isL = true; - if (mantilambda > resoCutVals[MassMin][Lambda] && mantilambda < resoCutVals[MassMax][Lambda]) + } + if (mantilambda > resoCutVals[MassMin][Lambda] && mantilambda < resoCutVals[MassMax][Lambda]) { selection.isAL = true; + } if (!selection.isL && !selection.isAL) { return selection; @@ -2187,41 +2387,50 @@ struct FlowGenericFramework { registryQA.fill(HIST("Lambda/hLambdaCount"), FillMassCut); // Rapidity correction - if (v0.yLambda() > resoCutVals[Rapidity][Lambda]) + if (v0.yLambda() > resoCutVals[Rapidity][Lambda]) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillRapidityCut); // DCA cuts for lambda and antilambda if (selection.isL) { - if (std::abs(v0.dcapostopv()) < resoCutVals[DCAPosToPVMin][Lambda] || std::abs(v0.dcanegtopv()) < resoCutVals[DCANegToPVMin][Lambda]) + if (std::abs(v0.dcapostopv()) < resoCutVals[DCAPosToPVMin][Lambda] || std::abs(v0.dcanegtopv()) < resoCutVals[DCANegToPVMin][Lambda]) { return selection; + } } if (selection.isAL) { - if (std::abs(v0.dcapostopv()) < resoCutVals[DCANegToPVMin][Lambda] || std::abs(v0.dcanegtopv()) < resoCutVals[DCAPosToPVMin][Lambda]) + if (std::abs(v0.dcapostopv()) < resoCutVals[DCANegToPVMin][Lambda] || std::abs(v0.dcanegtopv()) < resoCutVals[DCAPosToPVMin][Lambda]) { return selection; + } } registryQA.fill(HIST("Lambda/hLambdaCount"), FillDCAtoPV); - if (resoSwitchVals[UseDCAxDaughters][Lambda] && std::abs(v0.dcaV0daughters()) > resoCutVals[DCAxDaughters][Lambda]) + if (resoSwitchVals[UseDCAxDaughters][Lambda] != 0 && std::abs(v0.dcaV0daughters()) > resoCutVals[DCAxDaughters][Lambda]) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillDCAxDaughters); // v0 radius cuts - if (resoSwitchVals[UseV0Radius][Lambda] && (v0.v0radius() < resoCutVals[RadiusMin][Lambda] || v0.v0radius() > resoCutVals[RadiusMax][Lambda])) + if (resoSwitchVals[UseV0Radius][Lambda] != 0 && (v0.v0radius() < resoCutVals[RadiusMin][Lambda] || v0.v0radius() > resoCutVals[RadiusMax][Lambda])) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillV0Radius); // cosine pointing angle cuts - if (v0.v0cosPA() < resoCutVals[CosPA][Lambda]) + if (v0.v0cosPA() < resoCutVals[CosPA][Lambda]) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillCosPA); // Proper lifetime - if (resoSwitchVals[UseProperLifetime][Lambda] && v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * massLambda > resoCutVals[LifeTime][Lambda]) + if (resoSwitchVals[UseProperLifetime][Lambda] != 0 && v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * massLambda > resoCutVals[LifeTime][Lambda]) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillProperLifetime); // ArmenterosPodolanskiCut - if (resoSwitchVals[UseArmPodCut][Lambda] && (v0.qtarm() / std::abs(v0.alpha())) < resoCutVals[ArmPodMin][Lambda]) + if (resoSwitchVals[UseArmPodCut][Lambda] != 0 && (v0.qtarm() / std::abs(v0.alpha())) < resoCutVals[ArmPodMin][Lambda]) { return selection; + } registryQA.fill(HIST("Lambda/hLambdaCount"), FillArmPodCut); - if (resoSwitchVals[UseCompetingMassRejection][Lambda]) { - if (std::abs(v0.mK0Short() - o2::constants::physics::MassK0Short) < resoCutVals[MassRejection][Lambda]) + if (resoSwitchVals[UseCompetingMassRejection][Lambda] != 0) { + if (std::abs(v0.mK0Short() - o2::constants::physics::MassK0Short) < resoCutVals[MassRejection][Lambda]) { return selection; + } } registryQA.fill(HIST("Lambda/hLambdaCount"), FillCompetingMass); if (selection.isL) { @@ -2242,20 +2451,22 @@ struct FlowGenericFramework { } template - inline void fillPtSums(TTrack track, const float& centrality, const double& vtxz) + inline void fillPtSums(const TTrack& track, const float& centrality, const double& vtxz) { double wacc = (dt == Gen) ? 1. : getAcceptance(track, vtxz, 0); double weff = (dt == Gen) ? 1. : getEfficiency(track, centrality); - if (weff < 0) + if (weff < 0) { return; + } // Fill the nominal sums - if (track.eta() > cfgKinematics.cfgEtaPtPt->first && track.eta() < cfgKinematics.cfgEtaPtPt->second) + if (track.eta() > cfgKinematics.cfgEtaPtPt->first && track.eta() < cfgKinematics.cfgEtaPtPt->second) { fFCpt->fill(weff, track.pt()); + } // Fill the subevent sums std::size_t index = 0; - for (const auto& [etamin, etamax] : o2::analysis::gfw::etagapsPtPt) { + for (const auto& [etamin, etamax] : gfwMemberCache.etagapsPtPt) { if (etamin < track.eta() && track.eta() < etamax) { fFCpt->fillSub(weff, track.pt(), index); } @@ -2270,40 +2481,49 @@ struct FlowGenericFramework { } template - inline void fillGFW(TTrack track, const float& centrality, const double& vtxz, int pid_index, DensityCorr densitycorrections) + inline void fillGFW(const TTrack& track, const float& centrality, const double& vtxz, int pid_index, const DensityCorr& densitycorrections) { if (cfgUsePID) { // Analysing POI flow with id'ed particles - double ptmins[] = {o2::analysis::gfw::ptpoilow, o2::analysis::gfw::ptpoilow, 0.3, 0.5}; - double ptmaxs[] = {o2::analysis::gfw::ptpoiup, o2::analysis::gfw::ptpoiup, 6.0, 6.0}; - bool withinPtRef = (track.pt() > o2::analysis::gfw::ptreflow && track.pt() < o2::analysis::gfw::ptrefup); + std::array ptmins = {gfwMemberCache.ptpoilow, gfwMemberCache.ptpoilow, 0.3, 0.5}; + std::array ptmaxs = {gfwMemberCache.ptpoiup, gfwMemberCache.ptpoiup, 6.0, 6.0}; + bool withinPtRef = (track.pt() > gfwMemberCache.ptreflow && track.pt() < gfwMemberCache.ptrefup); bool withinPtPOI = (track.pt() > ptmins[pid_index] && track.pt() < ptmaxs[pid_index]); bool withinPtNch = (track.pt() > ptmins[0] && track.pt() < ptmaxs[0]); - if (!withinPtPOI && !withinPtRef) + if (!withinPtPOI && !withinPtRef) { return; + } double waccRef = (dt == Gen) ? 1. : getAcceptance(track, vtxz, 0); double waccPOI = (dt == Gen) ? 1. : withinPtPOI ? getAcceptance(track, vtxz, pid_index + 1) : getAcceptance(track, vtxz, 0); // - if (withinPtRef && withinPtPOI && pid_index) + if (withinPtRef && withinPtPOI && pid_index) { waccRef = waccPOI; // if particle is both (then it's overlap), override ref with POI - if (withinPtRef) + } + if (withinPtRef) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), waccRef, 1); - if (withinPtPOI && pid_index) + } + if (withinPtPOI && pid_index) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), waccPOI, (1 << (pid_index + 1))); - if (withinPtNch) + } + if (withinPtNch) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), waccPOI, 2); - if (withinPtPOI && withinPtRef && pid_index) + } + if (withinPtPOI && withinPtRef && pid_index) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), waccPOI, (1 << (pid_index + 5))); - if (withinPtNch && withinPtRef) + } + if (withinPtNch && withinPtRef) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), waccPOI, 32); + } } else { // Analysing only integrated flow - bool withinPtRef = (track.pt() > o2::analysis::gfw::ptreflow && track.pt() < o2::analysis::gfw::ptrefup); - bool withinPtPOI = (track.pt() > o2::analysis::gfw::ptpoilow && track.pt() < o2::analysis::gfw::ptpoiup); - if (!withinPtPOI && !withinPtRef) + bool withinPtRef = (track.pt() > gfwMemberCache.ptreflow && track.pt() < gfwMemberCache.ptrefup); + bool withinPtPOI = (track.pt() > gfwMemberCache.ptpoilow && track.pt() < gfwMemberCache.ptpoiup); + if (!withinPtPOI && !withinPtRef) { return; + } double weff = (dt == Gen) ? 1. : getEfficiency(track, centrality, 0); - if (weff < 0) + if (weff < 0) { return; + } if (cfgUseDensityDependentCorrection && withinPtRef && dt != Gen) { double fphi = densitycorrections.v2 * std::cos(2 * (track.phi() - densitycorrections.psi2Est)) + densitycorrections.v3 * std::cos(3 * (track.phi() - densitycorrections.psi3Est)) + densitycorrections.v4 * std::cos(4 * (track.phi() - densitycorrections.psi4Est)); fphi = (1 + 2 * fphi); @@ -2317,18 +2537,20 @@ struct FlowGenericFramework { } } double wacc = (dt == Gen) ? 1. : getAcceptance(track, vtxz, 0); - if (withinPtRef) + if (withinPtRef) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), weff * wacc, 1); - if (withinPtPOI) + } + if (withinPtPOI) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), weff * wacc, 2); - if (withinPtRef && withinPtPOI) + } + if (withinPtRef && withinPtPOI) { fGFW->Fill(track.eta(), fPtAxis->FindBin(track.pt()) - 1, track.phi(), weff * wacc, 4); + } } - return; } template - inline void fillTrackQA(TTrack track, const float vtxz) + inline void fillTrackQA(const TTrack& track, const float vtxz) { if constexpr (dt == Gen) { registryQA.fill(HIST("MCGen/") + HIST(FillTimeName[ft]) + HIST("phi_eta_vtxZ_gen"), track.phi(), track.eta(), vtxz); @@ -2348,16 +2570,18 @@ struct FlowGenericFramework { registryQA.fill(HIST("trackQA/") + HIST(FillTimeName[ft]) + HIST("pt_ref"), track.pt()); registryQA.fill(HIST("trackQA/") + HIST(FillTimeName[ft]) + HIST("pt_poi"), track.pt()); - if (track.eta() > cfgKinematics.cfgEtaNch->first && track.eta() < cfgKinematics.cfgEtaNch->second) + if (track.eta() > cfgKinematics.cfgEtaNch->first && track.eta() < cfgKinematics.cfgEtaNch->second) { registryQA.fill(HIST("trackQA/after/etaNch"), track.eta()); - if (track.eta() > cfgKinematics.cfgEtaPtPt->first && track.eta() < cfgKinematics.cfgEtaPtPt->second && cfgFill.cfgFillQA) + } + if (track.eta() > cfgKinematics.cfgEtaPtPt->first && track.eta() < cfgKinematics.cfgEtaPtPt->second && cfgFill.cfgFillQA) { registryQA.fill(HIST("trackQA/after/etaPtPt"), track.eta()); + } } } } template - inline void fillV0QA(const V0Selection& selection, TV0 const& v0, TTrack postrack, TTrack negtrack, const float& centrality, const float& weff) + inline void fillV0QA(const V0Selection& selection, TV0 const& v0, const TTrack& postrack, const TTrack& negtrack, const float& centrality, const float& weff) { if (selection.isK0) { registryQA.fill(HIST("K0/hK0Mass_sparse"), v0.mK0Short(), v0.pt(), centrality); @@ -2369,8 +2593,9 @@ struct FlowGenericFramework { registryQA.fill(HIST("K0/PiMinusTOF_K0"), negtrack.pt(), negtrack.tofNSigmaKa()); registryQA.fill(HIST("K0/hK0s"), 0.5, 1); - if (cfgEfficiency.cfgUsePIDEfficiencies) + if (cfgEfficiency.cfgUsePIDEfficiencies) { registryQA.fill(HIST("K0/hK0s_corrected"), 0.5, weff); + } } if (selection.isL) { registryQA.fill(HIST("Lambda/hLambdaMass_sparse"), v0.mLambda(), v0.pt(), centrality); @@ -2392,13 +2617,14 @@ struct FlowGenericFramework { } if (selection.isL || selection.isAL) { registryQA.fill(HIST("Lambda/hLambdas"), 0.5, 1); - if (cfgEfficiency.cfgUsePIDEfficiencies) + if (cfgEfficiency.cfgUsePIDEfficiencies) { registryQA.fill(HIST("Lambda/hLambdas_corrected"), 0.5, weff); + } } } template - float getCentrality(TCollision collision) + float getCentrality(const TCollision& collision) { switch (cfgCentEstimator) { case CentFT0C: @@ -2420,8 +2646,8 @@ struct FlowGenericFramework { } } - template - inline void fillEventQA(CollisionObject collision, TracksObject tracks) + template + inline void fillEventQA(const TCollision& collision, const TTracks& tracks) { registryQA.fill(HIST("eventQA/") + HIST(FillTimeName[ft]) + HIST("globalTracks_centT0C"), collision.centFT0C(), tracks.size()); registryQA.fill(HIST("eventQA/") + HIST(FillTimeName[ft]) + HIST("PVTracks_centT0C"), collision.centFT0C(), collision.multNTracksPV()); @@ -2430,11 +2656,10 @@ struct FlowGenericFramework { registryQA.fill(HIST("eventQA/") + HIST(FillTimeName[ft]) + HIST("globalTracks_multV0A"), collision.multFV0A(), tracks.size()); registryQA.fill(HIST("eventQA/") + HIST(FillTimeName[ft]) + HIST("multV0A_multT0A"), collision.multFT0A(), collision.multFV0A()); registryQA.fill(HIST("eventQA/") + HIST(FillTimeName[ft]) + HIST("multT0C_centT0C"), collision.centFT0C(), collision.multFT0C()); - return; } o2::framework::expressions::Filter collisionFilter = nabs(aod::collision::posZ) < cfgVtxZ; - o2::framework::expressions::Filter trackFilter = (aod::track::eta > cfgKinematics.cfgEta->first) && (aod::track::eta < cfgKinematics.cfgEta->second) && (aod::track::pt > cfgPtmin) && (aod::track::pt < cfgPtmax) && ((requireGlobalTrackInFilter()) || (aod::track::isGlobalTrackSDD == (uint8_t)true)) && (aod::track::itsChi2NCl < cfgTrackCuts.cfgChi2PrITSCls) && (aod::track::tpcChi2NCl < cfgTrackCuts.cfgChi2PrTPCCls) && nabs(aod::track::dcaZ) < cfgTrackCuts.cfgDCAz; + o2::framework::expressions::Filter trackFilter = (aod::track::eta > cfgKinematics.cfgEta->first) && (aod::track::eta < cfgKinematics.cfgEta->second) && (aod::track::pt > cfgKinematics.cfgPtMin) && (aod::track::pt < cfgKinematics.cfgPtmax) && ((requireGlobalTrackInFilter()) || (aod::track::isGlobalTrackSDD == (uint8_t)true)) && (aod::track::itsChi2NCl < cfgTrackCuts.cfgChi2PrITSCls) && (aod::track::tpcChi2NCl < cfgTrackCuts.cfgChi2PrTPCCls) && nabs(aod::track::dcaZ) < cfgTrackCuts.cfgDCAz; using GFWCollisions = soa::Filtered>; using GFWMCCollisions = soa::Join; @@ -2459,52 +2684,64 @@ struct FlowGenericFramework { LOGF(info, "run = %d", run); if (cfgRunByRun) { if (std::find(runNumbers.begin(), runNumbers.end(), run) == runNumbers.end()) { - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { createRunByRunHistograms(run); + } runNumbers.push_back(run); } else { LOGF(info, "run %d already in runNumbers", run); } - if (!cfgFill.cfgFillWeights) + if (!cfgFill.cfgFillWeights) { loadCorrections(bc); + } } } - if (!cfgFill.cfgFillWeights && !cfgRunByRun) + if (!cfgFill.cfgFillWeights && !cfgRunByRun) { loadCorrections(bc); + } registryQA.fill(HIST("eventQA/eventSel"), 0.5); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(0.5); - if (!collision.sel8()) + } + if (!collision.sel8()) { return; + } registryQA.fill(HIST("eventQA/eventSel"), 1.5); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(1.5); + } float centrality = getCentrality(collision); - if (cfgDoOccupancySel) { + if (cfgEventSelection.cfgDoOccupancySel) { int occupancy = collision.trackOccupancyInTimeRange(); - if (cfgFill.cfgFillQA) + if (cfgFill.cfgFillQA) { registryQA.fill(HIST("eventQA/before/occ_mult_cent"), occupancy, tracks.size(), centrality); - if (occupancy < 0 || occupancy > cfgOccupancySelection) + } + if (occupancy < 0 || occupancy > cfgEventSelection.cfgOccupancySelection) { return; - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { registryQA.fill(HIST("eventQA/after/occ_mult_cent"), occupancy, tracks.size(), centrality); + } } registryQA.fill(HIST("eventQA/eventSel"), 2.5); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(2.5); - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { fillEventQA(collision, tracks); + } registryQA.fill(HIST("eventQA/before/centrality"), centrality); registryQA.fill(HIST("eventQA/before/multiplicity"), tracks.size()); - if (!eventSelected(collision, tracks.size(), centrality, run)) + if (!eventSelected(collision, tracks.size(), centrality, run)) { return; - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { fillEventQA(collision, tracks); + } registryQA.fill(HIST("eventQA/after/centrality"), centrality); registryQA.fill(HIST("eventQA/after/multiplicity"), tracks.size()); // Get magnetic field polarity - int defaultMagCut = 99999; - auto field = (cfgMagField == defaultMagCut) ? getMagneticField(bc.timestamp()) : cfgMagField; + auto field = (cfgEventSelection.cfgMagField == DefaultMagneticFieldCut) ? getMagneticField(bc.timestamp()) : static_cast(cfgEventSelection.cfgMagField); processCollision(collision, tracks, v0s, centrality, field, run); } @@ -2516,45 +2753,54 @@ struct FlowGenericFramework { int run = bc.runNumber(); if (run != lastRun) { lastRun = run; - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { createRunByRunHistograms(run); + } } registryQA.fill(HIST("eventQA/eventSel"), 0.5); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(0.5); + } - if (!collision.sel8()) + if (!collision.sel8()) { return; + } registryQA.fill(HIST("eventQA/eventSel"), 1.5); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(1.5); + } const auto centrality = getCentrality(collision); - if (cfgDoOccupancySel) { + if (cfgEventSelection.cfgDoOccupancySel) { int occupancy = collision.trackOccupancyInTimeRange(); registryQA.fill(HIST("eventQA/before/occ_mult_cent"), occupancy, tracks.size(), centrality); - if (occupancy < 0 || occupancy > cfgOccupancySelection) + if (occupancy < 0 || occupancy > cfgEventSelection.cfgOccupancySelection) { return; + } registryQA.fill(HIST("eventQA/after/occ_mult_cent"), occupancy, tracks.size(), centrality); } registryQA.fill(HIST("eventQA/eventSel"), 2.5); - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { th1sList[run][EventSel]->Fill(2.5); + } - if (cfgFill.cfgFillQA) + if (cfgFill.cfgFillQA) { fillEventQA(collision, tracks); - if (!eventSelected(collision, tracks.size(), centrality, run)) + } + if (!eventSelected(collision, tracks.size(), centrality, run)) { return; - if (cfgFill.cfgFillQA) + } + if (cfgFill.cfgFillQA) { fillEventQA(collision, tracks); + } - if (!cfgFill.cfgFillWeights) + if (!cfgFill.cfgFillWeights) { loadCorrections(bc); - int defaultMagCut = 99999; - auto field = (cfgMagField == defaultMagCut) ? getMagneticField(bc.timestamp()) : cfgMagField; + } + auto field = (cfgEventSelection.cfgMagField == DefaultMagneticFieldCut) ? getMagneticField(bc.timestamp()) : static_cast(cfgEventSelection.cfgMagField); processCollision(collision, tracks, v0s, centrality, field, run); } PROCESS_SWITCH(FlowGenericFramework, processMCReco, "Process analysis for MC reconstructed events", false); @@ -2562,8 +2808,9 @@ struct FlowGenericFramework { o2::framework::expressions::Filter mcCollFilter = nabs(aod::mccollision::posZ) < cfgVtxZ; void processMCGen(soa::Filtered::iterator const& mcCollision, soa::SmallGroups> const& collisions, aod::McParticles const& particles, aod::V0Datas const& v0s) { - if (collisions.size() != 1) + if (collisions.size() != 1) { return; + } float centrality = -1; for (const auto& collision : collisions) { centrality = getCentrality(collision); @@ -2582,87 +2829,100 @@ struct FlowGenericFramework { PROCESS_SWITCH(FlowGenericFramework, processOnTheFly, "Process analysis for MC on-the-fly generated events", false); template - bool isGeneratedEfficiencyV0(TParticle const& particle, int pdgCode, int resonanceIndex) + bool isGeneratedEfficiencyV0(const TParticle& particle, int pdgCode, int resonanceIndex) { return particle.isPhysicalPrimary() && std::abs(particle.pdgCode()) == pdgCode && std::abs(particle.y()) < resoCutVals[Rapidity][resonanceIndex]; } template - bool isGeneratedEfficiencyTrack(TParticle const& particle) + bool isGeneratedEfficiencyTrack(const TParticle& particle) { auto pdgParticle = pdgDB->GetParticle(particle.pdgCode()); - return pdgParticle && pdgParticle->Charge() != 0 && particle.isPhysicalPrimary() && particle.eta() > o2::analysis::gfw::etalow && particle.eta() < o2::analysis::gfw::etaup && particle.pt() > o2::analysis::gfw::ptlow && particle.pt() < o2::analysis::gfw::ptup; + return pdgParticle && pdgParticle->Charge() != 0 && particle.isPhysicalPrimary() && particle.eta() > gfwMemberCache.etalow && particle.eta() < gfwMemberCache.etaup && particle.pt() > gfwMemberCache.ptlow && particle.pt() < gfwMemberCache.ptup; } template - int getEfficiencyTruthPID(TParticle const& particle) + int getEfficiencyTruthPID(const TParticle& particle) { - if (std::abs(particle.pdgCode()) == kPiPlus) + if (std::abs(particle.pdgCode()) == kPiPlus) { return EfficiencyPion; - if (std::abs(particle.pdgCode()) == kKPlus) + } + if (std::abs(particle.pdgCode()) == kKPlus) { return EfficiencyKaon; - if (std::abs(particle.pdgCode()) == kProton) + } + if (std::abs(particle.pdgCode()) == kProton) { return EfficiencyProton; + } return EfficiencyCharged; } template - void fillGeneratedEfficiencyTrack(TParticle const& particle, float centrality) + void fillGeneratedEfficiencyTrack(const TParticle& particle, float centrality) { - if (!isGeneratedEfficiencyTrack(particle)) + if (!isGeneratedEfficiencyTrack(particle)) { return; - if (particle.eta() <= cfgKinematics.cfgEta->first || particle.eta() >= cfgKinematics.cfgEta->second) + } + if (particle.eta() <= cfgKinematics.cfgEta->first || particle.eta() >= cfgKinematics.cfgEta->second) { return; + } const int pidIndex = getEfficiencyTruthPID(particle); registry.fill(HIST("Efficiency/efficiencyHist"), particle.pt(), centrality, EfficiencyCharged, EfficiencyGenerated); - if (pidIndex != EfficiencyCharged) + if (pidIndex != EfficiencyCharged) { registry.fill(HIST("Efficiency/efficiencyHist"), particle.pt(), centrality, pidIndex, EfficiencyGenerated); + } } template - void fillGeneratedEfficiencyV0(TParticle const& particle, int particleIndex, float centrality) + void fillGeneratedEfficiencyV0(const TParticle& particle, int particleIndex, float centrality) { - if (!particle.has_daughters()) + if (!particle.has_daughters()) { return; + } for (const auto& daughter : particle.template daughters_as()) { - if (daughter.eta() <= cfgKinematics.cfgEta->first || daughter.eta() >= cfgKinematics.cfgEta->second) + if (daughter.eta() <= cfgKinematics.cfgEta->first || daughter.eta() >= cfgKinematics.cfgEta->second) { return; + } } registry.fill(HIST("Efficiency/efficiencyHist"), particle.pt(), centrality, particleIndex, EfficiencyGenerated); } template - bool isEfficiencyEventSelected(TCollision const& collision, const int multTrk, float& centrality) + bool isEfficiencyEventSelected(const TCollision& collision, const int multTrk, float& centrality) { - if (!collision.sel8()) + if (!collision.sel8()) { return false; + } centrality = getCentrality(collision); - if (centrality < o2::analysis::gfw::centbinning.front() || centrality > o2::analysis::gfw::centbinning.back()) + if (centrality < gfwMemberCache.centbinning.front() || centrality > gfwMemberCache.centbinning.back()) { return false; - if (cfgDoOccupancySel) { + } + if (cfgEventSelection.cfgDoOccupancySel) { int occupancy = collision.trackOccupancyInTimeRange(); - if (occupancy < 0 || occupancy > cfgOccupancySelection) + if (occupancy < 0 || occupancy > cfgEventSelection.cfgOccupancySelection) { return false; + } } const int run = 0; return eventSelected(collision, multTrk, centrality, run); } template - bool isWithinEfficiencyEtaAcceptance(TTrack const& track) + bool isWithinEfficiencyEtaAcceptance(const TTrack& track) { return track.eta() > cfgKinematics.cfgEta->first && track.eta() < cfgKinematics.cfgEta->second; } template - bool areGeneratedDaughtersWithinEfficiencyEtaAcceptance(TParticle const& particle) + bool areGeneratedDaughtersWithinEfficiencyEtaAcceptance(const TParticle& particle) { - if (!particle.has_daughters()) + if (!particle.has_daughters()) { return false; + } for (const auto& daughter : particle.template daughters_as()) { - if (daughter.eta() <= cfgKinematics.cfgEta->first || daughter.eta() >= cfgKinematics.cfgEta->second) + if (daughter.eta() <= cfgKinematics.cfgEta->first || daughter.eta() >= cfgKinematics.cfgEta->second) { return false; + } } return true; } @@ -2673,35 +2933,42 @@ struct FlowGenericFramework { } template - void fillEnabledV0EfficiencyRecoDebug(TV0 const& v0, float centrality, int stage) + void fillEnabledV0EfficiencyRecoDebug(const TV0& v0, float centrality, int stage) { - if (resoSwitchVals[UseParticle][K0]) + if (resoSwitchVals[UseParticle][K0] != 0) { fillV0EfficiencyRecoDebug(v0.pt(), centrality, EfficiencyK0, stage); - if (resoSwitchVals[UseParticle][Lambda]) + } + if (resoSwitchVals[UseParticle][Lambda] != 0) { fillV0EfficiencyRecoDebug(v0.pt(), centrality, EfficiencyLambda, stage); + } } template - void fillEfficiencyRecoTrack(TTrack const& track, int field, float centrality) + void fillEfficiencyRecoTrack(const TTrack& track, int field, float centrality) { - if (track.mcParticleId() < 0 || !track.has_mcParticle()) + if (track.mcParticleId() < 0 || !track.has_mcParticle()) { return; + } auto mcParticle = track.mcParticle(); - if (!isGeneratedEfficiencyTrack(mcParticle)) + if (!isGeneratedEfficiencyTrack(mcParticle)) { return; - if (!trackSelected(track, field)) + } + if (!trackSelected(track, field)) { return; - if (!isWithinEfficiencyEtaAcceptance(mcParticle) || !isWithinEfficiencyEtaAcceptance(track)) + } + if (!isWithinEfficiencyEtaAcceptance(mcParticle) || !isWithinEfficiencyEtaAcceptance(track)) { return; + } const int pidIndex = getEfficiencyTruthPID(mcParticle); registry.fill(HIST("Efficiency/efficiencyHist"), mcParticle.pt(), centrality, EfficiencyCharged, EfficiencyReconstructed); - if (pidIndex != EfficiencyCharged) + if (pidIndex != EfficiencyCharged) { registry.fill(HIST("Efficiency/efficiencyHist"), mcParticle.pt(), centrality, pidIndex, EfficiencyReconstructed); + } } template - void fillEfficiencyRecoV0(TV0 const& v0, TCollision const& collision, TTracks const& tracks, float centrality) + void fillEfficiencyRecoV0(const TV0& v0, const TCollision& collision, const TTracks& tracks, float centrality) { using V0TrackTable = std::decay_t; fillEnabledV0EfficiencyRecoDebug(v0, centrality, V0EfficiencyCandidate); @@ -2709,33 +2976,39 @@ struct FlowGenericFramework { auto postrack = v0.template posTrack_as(); auto negtrack = v0.template negTrack_as(); - if (!postrack.has_mcParticle() || !negtrack.has_mcParticle()) + if (!postrack.has_mcParticle() || !negtrack.has_mcParticle()) { return; + } fillEnabledV0EfficiencyRecoDebug(v0, centrality, V0EfficiencyDaughterMCLabels); auto posMcParticle = postrack.template mcParticle_as(); auto negMcParticle = negtrack.template mcParticle_as(); - if (!posMcParticle.has_mothers() || !negMcParticle.has_mothers()) + if (!posMcParticle.has_mothers() || !negMcParticle.has_mothers()) { return; + } fillEnabledV0EfficiencyRecoDebug(v0, centrality, V0EfficiencyDaughterMothers); for (const auto& posMother : posMcParticle.template mothers_as()) { for (const auto& negMother : negMcParticle.template mothers_as()) { - if (posMother.globalIndex() != negMother.globalIndex() || !posMother.isPhysicalPrimary()) + if (posMother.globalIndex() != negMother.globalIndex() || !posMother.isPhysicalPrimary()) { continue; + } int pdgCode = std::abs(posMother.pdgCode()); - if (pdgCode == PDG_t::kK0Short && resoSwitchVals[UseParticle][K0]) { + if (pdgCode == PDG_t::kK0Short && resoSwitchVals[UseParticle][K0] != 0) { fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyK0, V0EfficiencyCommonPrimaryMother); fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyK0, V0EfficiencyMotherSpecies); - if (std::abs(posMother.y()) >= resoCutVals[Rapidity][K0]) + if (std::abs(posMother.y()) >= resoCutVals[Rapidity][K0]) { return; + } fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyK0, V0EfficiencyMotherRapidity); - if (!areGeneratedDaughtersWithinEfficiencyEtaAcceptance(posMother)) + if (!areGeneratedDaughtersWithinEfficiencyEtaAcceptance(posMother)) { return; + } fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyK0, V0EfficiencyGeneratedDaughterEta); - if (!isWithinEfficiencyEtaAcceptance(postrack) || !isWithinEfficiencyEtaAcceptance(negtrack)) + if (!isWithinEfficiencyEtaAcceptance(postrack) || !isWithinEfficiencyEtaAcceptance(negtrack)) { return; + } fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyK0, V0EfficiencyRecoDaughterEta); auto selection = selectK0(collision, v0, tracks); if (selection.selected) { @@ -2744,17 +3017,22 @@ struct FlowGenericFramework { fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyK0, V0EfficiencyFilled); } return; - } else if (pdgCode == PDG_t::kLambda0 && resoSwitchVals[UseParticle][Lambda]) { + } + + if (pdgCode == PDG_t::kLambda0 && resoSwitchVals[UseParticle][Lambda] != 0) { fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyLambda, V0EfficiencyCommonPrimaryMother); fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyLambda, V0EfficiencyMotherSpecies); - if (std::abs(posMother.y()) >= resoCutVals[Rapidity][Lambda]) + if (std::abs(posMother.y()) >= resoCutVals[Rapidity][Lambda]) { return; + } fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyLambda, V0EfficiencyMotherRapidity); - if (!areGeneratedDaughtersWithinEfficiencyEtaAcceptance(posMother)) + if (!areGeneratedDaughtersWithinEfficiencyEtaAcceptance(posMother)) { return; + } fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyLambda, V0EfficiencyGeneratedDaughterEta); - if (!isWithinEfficiencyEtaAcceptance(postrack) || !isWithinEfficiencyEtaAcceptance(negtrack)) + if (!isWithinEfficiencyEtaAcceptance(postrack) || !isWithinEfficiencyEtaAcceptance(negtrack)) { return; + } fillV0EfficiencyRecoDebug(posMother.pt(), centrality, EfficiencyLambda, V0EfficiencyRecoDaughterEta); auto selection = selectLambda(collision, v0, tracks); if (selection.selected) { @@ -2778,12 +3056,14 @@ struct FlowGenericFramework { for (const auto& collision : collisions) { int tracksThisCollision = 0; for (const auto& track : tracks) { - if (track.collisionId() == collision.globalIndex()) + if (track.collisionId() == collision.globalIndex()) { ++tracksThisCollision; + } } float centrality = -1.f; - if (!isEfficiencyEventSelected(collision, tracksThisCollision, centrality)) + if (!isEfficiencyEventSelected(collision, tracksThisCollision, centrality)) { continue; + } if (!hasSelectedCollision || collision.numContrib() > bestNumContrib) { hasSelectedCollision = true; bestNumContrib = collision.numContrib(); @@ -2791,33 +3071,39 @@ struct FlowGenericFramework { selectedCentrality = centrality; } } - if (!hasSelectedCollision) + if (!hasSelectedCollision) { return; + } for (const auto& particle : mcParticles) { - if (particle.mcCollisionId() != mcCollision.globalIndex()) + if (particle.mcCollisionId() != mcCollision.globalIndex()) { continue; + } fillGeneratedEfficiencyTrack(particle, selectedCentrality); - if (isGeneratedEfficiencyV0(particle, PDG_t::kK0Short, K0) && resoSwitchVals[UseParticle][K0]) + if (isGeneratedEfficiencyV0(particle, PDG_t::kK0Short, K0) && resoSwitchVals[UseParticle][K0] != 0) { fillGeneratedEfficiencyV0(particle, EfficiencyK0, selectedCentrality); - if (isGeneratedEfficiencyV0(particle, PDG_t::kLambda0, Lambda) && resoSwitchVals[UseParticle][Lambda]) + } + if (isGeneratedEfficiencyV0(particle, PDG_t::kLambda0, Lambda) && resoSwitchVals[UseParticle][Lambda] != 0) { fillGeneratedEfficiencyV0(particle, EfficiencyLambda, selectedCentrality); + } } for (const auto& collision : collisions) { - if (collision.globalIndex() != bestCollisionIndex) + if (collision.globalIndex() != bestCollisionIndex) { continue; + } auto bc = collision.template bc_as(); - int defaultMagCut = 99999; - auto field = (cfgMagField == defaultMagCut) ? getMagneticField(bc.timestamp()) : cfgMagField; + auto field = (cfgEventSelection.cfgMagField == DefaultMagneticFieldCut) ? getMagneticField(bc.timestamp()) : static_cast(cfgEventSelection.cfgMagField); for (const auto& track : tracks) { - if (track.collisionId() != bestCollisionIndex) + if (track.collisionId() != bestCollisionIndex) { continue; + } fillEfficiencyRecoTrack(track, field, selectedCentrality); } for (const auto& v0 : v0s) { - if (v0.collisionId() != bestCollisionIndex) + if (v0.collisionId() != bestCollisionIndex) { continue; + } fillEfficiencyRecoV0(v0, collision, tracks, selectedCentrality); } break; @@ -2831,16 +3117,18 @@ struct FlowGenericFramework { int run = bc.runNumber(); if (run != lastRun) { lastRun = run; - if (cfgFill.cfgFillRunByRunQA) + if (cfgFill.cfgFillRunByRunQA) { createRunByRunHistograms(run); + } } - if (!collision.sel7()) + if (!collision.sel7()) { return; + } const auto centrality = collision.centRun2V0M(); - if (!cfgFill.cfgFillWeights) + if (!cfgFill.cfgFillWeights) { loadCorrections(bc); - int defaultMagCut = 99999; - auto field = (cfgMagField == defaultMagCut) ? getMagneticField(bc.timestamp()) : cfgMagField; + } + auto field = (cfgEventSelection.cfgMagField == DefaultMagneticFieldCut) ? getMagneticField(bc.timestamp()) : static_cast(cfgEventSelection.cfgMagField); processCollision(collision, tracks, v0s, centrality, field, run); } PROCESS_SWITCH(FlowGenericFramework, processRun2, "Process analysis for Run 2 converted data", false);