From 98d431e8bd3cb70df6f3b4014cb0dafcfeb0a058 Mon Sep 17 00:00:00 2001 From: Gleb Romanenko Date: Fri, 19 Apr 2024 18:34:59 +0200 Subject: [PATCH 1/4] adding table with pileup flags and new analysis distributions --- PWGCF/Femto3D/DataModel/singletrackselector.h | 9 + PWGCF/Femto3D/TableProducer/CMakeLists.txt | 5 + .../TableProducer/singleTrackSelector.cxx | 92 +++--- .../singleTrackSelectorExtra.cxx | 49 ++++ PWGCF/Femto3D/Tasks/femto3dPairTask.cxx | 100 +++++-- PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx | 262 ++++++++++-------- PWGCF/Femto3D/Tasks/femto3dQA.cxx | 31 ++- 7 files changed, 362 insertions(+), 186 deletions(-) create mode 100755 PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx diff --git a/PWGCF/Femto3D/DataModel/singletrackselector.h b/PWGCF/Femto3D/DataModel/singletrackselector.h index c1217abb080..e91dafc7bbd 100755 --- a/PWGCF/Femto3D/DataModel/singletrackselector.h +++ b/PWGCF/Femto3D/DataModel/singletrackselector.h @@ -104,6 +104,10 @@ DECLARE_SOA_COLUMN(MultPercentile, multPerc, float); // Percentiles of multiplic DECLARE_SOA_COLUMN(PosZ, posZ, float); // Vertex of the collision DECLARE_SOA_COLUMN(MagField, magField, float); // Magnetic field corresponding to a collision (in T) +DECLARE_SOA_COLUMN(IsNoSameBunchPileup, isNoSameBunchPileup, bool); +DECLARE_SOA_COLUMN(IsGoodZvtxFT0vsPV, isGoodZvtxFT0vsPV, bool); +DECLARE_SOA_COLUMN(IsVertexITSTPC, isVertexITSTPC, bool); + } // namespace singletrackselector DECLARE_SOA_TABLE(SingleCollSels, "AOD", "SINGLECOLLSEL", // Table of the variables for single track selection. @@ -113,6 +117,11 @@ DECLARE_SOA_TABLE(SingleCollSels, "AOD", "SINGLECOLLSEL", // Table of the variab singletrackselector::PosZ, singletrackselector::MagField); +DECLARE_SOA_TABLE(SingleCollExtras, "AOD", "SINGLECOLLEXTRA", // Joinable collision table with Pile-Up flags + singletrackselector::IsNoSameBunchPileup, + singletrackselector::IsGoodZvtxFT0vsPV, + singletrackselector::IsVertexITSTPC); + namespace singletrackselector { DECLARE_SOA_INDEX_COLUMN(SingleCollSel, singleCollSel); // Index to the collision diff --git a/PWGCF/Femto3D/TableProducer/CMakeLists.txt b/PWGCF/Femto3D/TableProducer/CMakeLists.txt index 331407e6537..48750decf2e 100755 --- a/PWGCF/Femto3D/TableProducer/CMakeLists.txt +++ b/PWGCF/Femto3D/TableProducer/CMakeLists.txt @@ -12,4 +12,9 @@ o2physics_add_dpl_workflow(single-track-selector SOURCES singleTrackSelector.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGCFCore + COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(single-track-selector-extra + SOURCES singleTrackSelectorExtra.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGCFCore COMPONENT_NAME Analysis) \ No newline at end of file diff --git a/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx b/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx index f61dd67ec1f..3f52e009cf4 100755 --- a/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx +++ b/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx @@ -52,13 +52,14 @@ struct singleTrackSelector { Configurable centTableToUse{"centTableToUse", 1, "Flag to choose cent./mult.perc. estimator (Run3 only [FTOC for PbPb; FTOM for pp], for Run2 the V0M is used): 0 -> CentFV0As, 1 -> CentFT0Ms, 2 -> CentFT0As, 3 -> CentFT0Cs, 4 -> CentFDDMs, 5 -> CentNTPVs"}; Configurable multTableToUse{"multTableToUse", 1, "Flag to choose mult. estimator (Run3 only): 0 -> TPCMults, 1 -> MultNTracksPV, 2 -> MultNTracksPVeta1"}; Configurable rejectNotPropagatedTrks{"rejectNotPropagatedTrks", true, "rejects tracks that are not propagated to the primary vertex"}; - Configurable removeTFBorder{"removeTFBorder", false, "Remove TF border"}; Configurable> _particlesToKeep{"particlesToKeepPDGs", std::vector{2212, 1000010020}, "PDG codes of perticles for which the 'singletrackselector' tables will be created (only proton and deurton are supported now)"}; Configurable> keepWithinNsigmaTPC{"keepWithinNsigmaTPC", std::vector{-4.0f, 4.0f}, "TPC range for preselection of particles specified with PDG"}; Configurable> _particlesToReject{"particlesToRejectPDGs", std::vector{211, 321}, "PDG codes of particles that will be rejected with TOF (only pion, kaon, proton and deurton are supported now)"}; Configurable> rejectWithinNsigmaTOF{"rejectWithinNsigmaTOF", std::vector{-5.0f, 5.0f}, "TOF rejection Nsigma range for particles specified with PDG to be rejected"}; + Configurable _pRemoveTofOutOfRange{"pRemoveTofOutOfRange", 100.f, "momentum starting from which request TOF nSigma to be within the stored range (-10 < Nsigma < 10)"}; + Configurable _min_P{"min_P", 0.f, "lower mometum limit"}; Configurable _max_P{"max_P", 100.f, "upper mometum limit"}; Configurable _eta{"eta", 100.f, "abs eta value limit"}; @@ -78,8 +79,9 @@ struct singleTrackSelector { using CollRun2 = soa::Join; using CollRun3 = soa::Join; - Produces tableRow; Produces tableRowColl; + Produces tableRowCollExtra; + Produces tableRow; Produces tableRowExtra; Produces tableRowMC; @@ -100,7 +102,7 @@ struct singleTrackSelector { std::vector particlesToKeep; std::vector particlesToReject; - void init(InitContext&) + void init(InitContext& context) { particlesToKeep = _particlesToKeep; @@ -169,6 +171,7 @@ struct singleTrackSelector { for (auto ii : particlesToKeep) if (o2::aod::singletrackselector::TPCselection(track, std::make_pair(ii, keepWithinNsigmaTPC))) { + if(track.p() > _pRemoveTofOutOfRange && !o2::aod::singletrackselector::TOFselection(track, std::make_pair(ii, std::vector{-10.0, 10.0}), 10.0)) continue; tableRow(tableRowColl.lastIndex(), track.p(), @@ -261,58 +264,59 @@ struct singleTrackSelector { auto bc = collision.bc_as(); initCCDB(bc); - if (removeTFBorder && collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { + float centValue = -100.0f; - float centValue = -100.0f; + switch (centTableToUse) { + case 0: + centValue = collision.centFV0A(); + break; + case 1: + centValue = collision.centFT0M(); + break; + case 2: + centValue = collision.centFT0A(); + break; + case 3: + centValue = collision.centFT0C(); + break; + case 4: + centValue = collision.centFDDM(); + break; + case 5: + centValue = collision.centNTPV(); + break; + default: + LOGF(fatal, "Invalid flag for cent./mult.perc. estimator has been choosen. Please check."); + break; + } + if (centValue >= _centCut.value.first && centValue <= _centCut.value.second) { + int multValue = -1; - switch (centTableToUse) { + switch (multTableToUse) { case 0: - centValue = collision.centFV0A(); + multValue = collision.multTPC(); break; case 1: - centValue = collision.centFT0M(); + multValue = collision.multNTracksPV(); break; case 2: - centValue = collision.centFT0A(); - break; - case 3: - centValue = collision.centFT0C(); - break; - case 4: - centValue = collision.centFDDM(); - break; - case 5: - centValue = collision.centNTPV(); + multValue = collision.multNTracksPVeta1(); break; default: - LOGF(fatal, "Invalid flag for cent./mult.perc. estimator has been choosen. Please check."); + LOGF(fatal, "Invalid flag for mult. estimator has been choosen. Please check."); break; } - if (centValue >= _centCut.value.first && centValue <= _centCut.value.second) { - int multValue = -1; - - switch (multTableToUse) { - case 0: - multValue = collision.multTPC(); - break; - case 1: - multValue = collision.multNTracksPV(); - break; - case 2: - multValue = collision.multNTracksPVeta1(); - break; - default: - LOGF(fatal, "Invalid flag for mult. estimator has been choosen. Please check."); - break; - } - tableRowColl(multValue, - centValue, - collision.posZ(), - d_bz); + tableRowColl(multValue, + centValue, + collision.posZ(), + d_bz); - fillTrackTables(tracks); - } + tableRowCollExtra(collision.selection_bit(aod::evsel::kNoSameBunchPileup), + collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV), + collision.selection_bit(aod::evsel::kIsVertexITSTPC)); + + fillTrackTables(tracks); } } PROCESS_SWITCH(singleTrackSelector, processDataRun3, "process data Run3", true); @@ -405,6 +409,10 @@ struct singleTrackSelector { collision.posZ(), d_bz); + tableRowCollExtra(collision.selection_bit(aod::evsel::kNoSameBunchPileup), + collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV), + collision.selection_bit(aod::evsel::kIsVertexITSTPC)); + fillTrackTables(tracks); } } diff --git a/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx b/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx new file mode 100755 index 00000000000..f04f9eaced0 --- /dev/null +++ b/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx @@ -0,0 +1,49 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +// +/// \brief create a table applying some basic cuts on the ITS and DCA. +/// \author Sofia Tomassini, Gleb Romanenko, Nicolò Jacazio +/// \since 31 May 2023 + + +// this task produces a dummy "SingleCollExtras" table that is now required in the analysis tasks. Needed to have a compatibility with old der. data + +#include +#include + +#include "PWGCF/Femto3D/DataModel/singletrackselector.h" + +#include "Framework/AnalysisTask.h" +#include "Framework/runDataProcessing.h" + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; +using namespace o2::track; +using namespace o2::aod; +//::singletrackselector; // the namespace defined in .h + +struct singleTrackSelectorDummy { + + Produces tableRowCollExtra; + + void process(aod::SingleCollSels::iterator const&) + { + tableRowCollExtra(true, + true, + true); + } +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{adaptAnalysisTask(cfgc)}; +} diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx index f9b85d33978..396eec4f7ff 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx @@ -47,6 +47,10 @@ struct FemtoCorrelations { /// Construct a registry object with direct declaration HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject}; + Configurable _removeSameBunchPileup{"removeSameBunchPileup", false, ""}; + Configurable _requestGoodZvtxFT0vsPV{"requestGoodZvtxFT0vsPV", false, ""}; + Configurable _requestVertexITSTPC{"requestVertexITSTPC", false, ""}; + Configurable _min_P{"min_P", 0.0, "lower mometum limit"}; Configurable _max_P{"max_P", 100.0, "upper mometum limit"}; Configurable _eta{"eta", 100.0, "abs eta value limit"}; @@ -88,6 +92,7 @@ struct FemtoCorrelations { ConfigurableAxis CFkStarBinning{"CFkStarBinning", {500, 0.005, 5.005}, "k* binning of the CF (Nbins, lowlimit, uplimit)"}; Configurable _fill3dCF{"fill3dCF", false, "flag for filling 3D LCMS histos: true -- fill; false -- not"}; + Configurable _fillDetaDphi{"fillDetaDphi", false, "flag for filling dEta(dPhi*) histos: true -- fill; false -- not; note: they're filled before the double track cut"}; ConfigurableAxis CF3DqLCMSBinning{"CF3DqLCMSBinning", {60, -0.3, 0.3}, "q_out/side/long binning of the CF 3D in LCMS (Nbins, lowlimit, uplimit)"}; // the next configarable is responsible for skipping (pseudo)randomly chosen ($value -1) pairs of events in the mixing process // migth be useful (for the sake of execution time and the output file size ...) in case of too many events per DF since in the SE thacks are mixed in N_ev and in the ME 0.5*N_ev*(N_ev - 1) @@ -106,7 +111,7 @@ struct FemtoCorrelations { std::pair> TPCcuts_2; std::pair> TOFcuts_2; - using FilteredCollisions = aod::SingleCollSels; + using FilteredCollisions = soa::Join; using FilteredTracks = aod::SingleTrackSels; typedef std::shared_ptr::iterator> trkType; @@ -142,6 +147,9 @@ struct FemtoCorrelations { std::vector>> MEhistos_3D; std::vector>> qLCMSvskStar; + std::vector>> DoubleTrack_SE_histos; + std::vector>> DoubleTrack_ME_histos; + void init(o2::framework::InitContext&) { @@ -205,6 +213,21 @@ struct FemtoCorrelations { MEhistos_3D.push_back(std::move(MEperMult_3D)); qLCMSvskStar.push_back(std::move(qLCMSvskStarperMult)); } + + if (_fillDetaDphi) { + std::vector> DoubleTrack_SE_histos_perMult; + std::vector> DoubleTrack_ME_histos_perMult; + + for (int j = 0; j < _kTbins.value.size() - 1; j++) { + auto hDblTrk_SE = registry.add(Form("Cent%i/DoubleTrackEffects_SE_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_SE_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); + auto hDblTrk_ME = registry.add(Form("Cent%i/DoubleTrackEffects_ME_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_ME_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); + + DoubleTrack_SE_histos_perMult.push_back(std::move(hDblTrk_SE)); + DoubleTrack_ME_histos_perMult.push_back(std::move(hDblTrk_ME)); + } + DoubleTrack_SE_histos.push_back(std::move(DoubleTrack_SE_histos_perMult)); + DoubleTrack_ME_histos.push_back(std::move(DoubleTrack_ME_histos_perMult)); + } } registry.add("p_first", Form("p_%i", static_cast(_particlePDG_1)), kTH1F, {{100, 0., 5., "p"}}); @@ -248,16 +271,17 @@ struct FemtoCorrelations { LOGF(fatal, "kTbin value obtained for a pair is less than 0"); } - if (!Pair->IsClosePair(_deta, _dphi, _radiusTPC)) { - kThistos[multBin][kTbin]->Fill(pair_kT); - mThistos[multBin][kTbin]->Fill(Pair->GetMt()); // test - SEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); // close pair rejection and fillig the SE histo + if (_fillDetaDphi) DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + if (Pair->IsClosePair(_deta, _dphi, _radiusTPC)) continue; - if (_fill3dCF) { - std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); - TVector3 qLCMS = std::pow(-1, (mt() % 2)) * Pair->GetQLCMS(); // introducing randomness to the pair order ([first, second]); important only for 3D because if there are any sudden order/correlation in the tables, it could couse unwanted asymmetries in the final 3d rel. momentum distributions; irrelevant in 1D case because the absolute value of the rel.momentum is taken - SEhistos_3D[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z()); - } + kThistos[multBin][kTbin]->Fill(pair_kT); + mThistos[multBin][kTbin]->Fill(Pair->GetMt()); // test + SEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); // close pair rejection and fillig the SE histo + + if (_fill3dCF) { + std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); + TVector3 qLCMS = std::pow(-1, (mt() % 2)) * Pair->GetQLCMS(); // introducing randomness to the pair order ([first, second]); important only for 3D because if there are any sudden order/correlation in the tables, it could couse unwanted asymmetries in the final 3d rel. momentum distributions; irrelevant in 1D case because the absolute value of the rel.momentum is taken + SEhistos_3D[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z()); } Pair->ResetPair(); } @@ -295,32 +319,37 @@ struct FemtoCorrelations { LOGF(fatal, "kTbin value obtained for a pair is less than 0"); } - if (!Pair->IsClosePair(_deta, _dphi, _radiusTPC)) { - if (!SE_or_ME) { - SEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); - kThistos[multBin][kTbin]->Fill(pair_kT); - mThistos[multBin][kTbin]->Fill(Pair->GetMt()); // test + if (_fillDetaDphi) { + if (!SE_or_ME) DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + else DoubleTrack_ME_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + } - if (_fill3dCF) { - std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); - TVector3 qLCMS = std::pow(-1, (mt() % 2)) * Pair->GetQLCMS(); // introducing randomness to the pair order ([first, second]); important only for 3D because if there are any sudden order/correlation in the tables, it could couse unwanted asymmetries in the final 3d rel. momentum distributions; irrelevant in 1D case because the absolute value of the rel.momentum is taken - SEhistos_3D[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z()); - } - } else { - MEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); + if (Pair->IsClosePair(_deta, _dphi, _radiusTPC)) continue; - if (_fill3dCF) { - std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); - TVector3 qLCMS = std::pow(-1, (mt() % 2)) * Pair->GetQLCMS(); // introducing randomness to the pair order ([first, second]); important only for 3D because if there are any sudden order/correlation in the tables, it could couse unwanted asymmetries in the final 3d rel. momentum distributions; irrelevant in 1D case because the absolute value of the rel.momentum is taken - MEhistos_3D[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z()); - qLCMSvskStar[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z(), Pair->GetKstar()); - } + if (!SE_or_ME) { + SEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); + kThistos[multBin][kTbin]->Fill(pair_kT); + mThistos[multBin][kTbin]->Fill(Pair->GetMt()); // test + + if (_fill3dCF) { + std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); + TVector3 qLCMS = std::pow(-1, (mt() % 2)) * Pair->GetQLCMS(); // introducing randomness to the pair order ([first, second]); important only for 3D because if there are any sudden order/correlation in the tables, it could couse unwanted asymmetries in the final 3d rel. momentum distributions; irrelevant in 1D case because the absolute value of the rel.momentum is taken + SEhistos_3D[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z()); + } + } else { + MEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); + + if (_fill3dCF) { + std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); + TVector3 qLCMS = std::pow(-1, (mt() % 2)) * Pair->GetQLCMS(); // introducing randomness to the pair order ([first, second]); important only for 3D because if there are any sudden order/correlation in the tables, it could couse unwanted asymmetries in the final 3d rel. momentum distributions; irrelevant in 1D case because the absolute value of the rel.momentum is taken + MEhistos_3D[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z()); + qLCMSvskStar[multBin][kTbin]->Fill(qLCMS.X(), qLCMS.Y(), qLCMS.Z(), Pair->GetKstar()); } } Pair->ResetPair(); } } - } + } void process(soa::Filtered const& collisions, soa::Filtered const& tracks) { @@ -330,6 +359,12 @@ struct FemtoCorrelations { for (auto track : tracks) { if (abs(track.template singleCollSel_as>().posZ()) > _vertexZ) continue; + if (_removeSameBunchPileup && !track.template singleCollSel_as>().isNoSameBunchPileup()) + continue; + if (_requestGoodZvtxFT0vsPV && !track.template singleCollSel_as>().isGoodZvtxFT0vsPV()) + continue; + if (_requestVertexITSTPC && !track.template singleCollSel_as>().isVertexITSTPC()) + continue; if (track.tpcFractionSharedCls() > _tpcFractionSharedCls || track.itsNCls() < _itsNCls) continue; if (track.template singleCollSel_as>().multPerc() < *_centBins.value.begin() || track.template singleCollSel_as>().multPerc() >= *(_centBins.value.end() - 1)) @@ -386,6 +421,13 @@ struct FemtoCorrelations { if (collision.multPerc() < *_centBins.value.begin() || collision.multPerc() >= *(_centBins.value.end() - 1)) continue; + if (_removeSameBunchPileup && !collision.isNoSameBunchPileup()) + continue; + if (_requestGoodZvtxFT0vsPV && !collision.isGoodZvtxFT0vsPV()) + continue; + if (_requestVertexITSTPC && !collision.isVertexITSTPC()) + continue; + if (selectedtracks_1.find(collision.globalIndex()) == selectedtracks_1.end()) { if (IsIdentical) continue; diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx index c14a97fcc06..7c26c74612f 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx @@ -75,10 +75,11 @@ struct FemtoCorrelationsMC { Configurable _radiusTPC{"radiusTPC", 1.2, "TPC radius to calculate phi_star for"}; - ConfigurableAxis CFkStarBinning{"CFkStarBinning", {500, 0.005, 5.005}, "k* binning of the res. matrix (Nbins, lowlimit, uplimit)"}; - Configurable _vertexNbinsToMix{"vertexNbinsToMix", 10, "Number of vertexZ bins for the mixing"}; - Configurable _multNsubBins{"multSubBins", 10, "number of sub-bins to perform the mixing within"}; + Configurable> _centBins{"multBins", std::vector{0.0f, 100.0f}, "multiplicity percentile/centrality binning (min:0, max:100)"}; + Configurable _multNsubBins{"multSubBins", 1, "number of sub-bins to perform the mixing within"}; + Configurable> _kTbins{"kTbins", std::vector{0.0f, 100.0f}, "pair transverse momentum kT binning"}; + ConfigurableAxis CFkStarBinning{"CFkStarBinning", {500, 0.005, 5.005}, "k* binning of the res. matrix (Nbins, lowlimit, uplimit)"}; bool IsIdentical; @@ -111,6 +112,17 @@ struct FemtoCorrelationsMC { Filter vertexFilter = nabs(o2::aod::singletrackselector::posZ) < _vertexZ; + std::vector>> DCA_histos_1; // key -- origin; origin = 0 - primiry, 1 - weak, 2 - material; + std::vector>> DCA_histos_2; // key -- origin; origin = 0 - primiry, 1 - weak, 2 - material; + + std::vector>> Purity_histos_1; // key -- PDG; PDG == 0 -> all selected tracks + std::vector>> Purity_histos_2; // key -- PDG; PDG == 0 -> all selected tracks + + std::vector>> kThistos; + std::vector>> Resolution_histos; + std::vector>> DoubleTrack_SE_histos; + std::vector>> DoubleTrack_ME_histos; + void init(o2::framework::InitContext&) { @@ -125,77 +137,147 @@ struct FemtoCorrelationsMC { TPCcuts_2 = std::make_pair(_particlePDG_2, _tpcNSigma_2); TOFcuts_2 = std::make_pair(_particlePDG_2, _tofNSigma_2); - registry.add("FirstParticle/dcaxyz_vs_pt_primary", "dcaxyz_vs_pt_primary", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) primary"}, {250, -1., 1., "DCA_Z(pt) primary"}}); - registry.add("FirstParticle/dcaxyz_vs_pt_weakdecay", "dcaxyz_vs_pt_weakdecay", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) weakdecay"}, {250, -1., 1., "DCA_Z(pt) weakdecay"}}); - registry.add("FirstParticle/dcaxyz_vs_pt_material", "dcaxyz_vs_pt_material", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) material"}, {250, -1., 1., "DCA_Z(pt) material"}}); + for (int i = 0; i < _centBins.value.size() - 1; i++) { + + std::map> DCA_histos_1_perMult; + DCA_histos_1_perMult[0] = registry.add(Form("Cent%i/FirstParticle/dcaxyz_vs_pt_primary", i), "dcaxyz_vs_pt_primary", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) primary"}, {250, -1., 1., "DCA_Z(pt) primary"}}); + DCA_histos_1_perMult[1] = registry.add(Form("Cent%i/FirstParticle/dcaxyz_vs_pt_weakdecay", i), "dcaxyz_vs_pt_weakdecay", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) weakdecay"}, {250, -1., 1., "DCA_Z(pt) weakdecay"}}); + DCA_histos_1_perMult[2] = registry.add(Form("Cent%i/FirstParticle/dcaxyz_vs_pt_material", i), "dcaxyz_vs_pt_material", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) material"}, {250, -1., 1., "DCA_Z(pt) material"}}); + + std::map> Purity_histos_1_perMult; + Purity_histos_1_perMult[11] = registry.add(Form("Cent%i/FirstParticle/pSpectraEl", i), "pSpectraEl", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_1_perMult[13] = registry.add(Form("Cent%i/FirstParticle/pSpectraMu", i), "pSpectraMu", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_1_perMult[211] = registry.add(Form("Cent%i/FirstParticle/pSpectraPi", i), "pSpectraPi", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_1_perMult[321] = registry.add(Form("Cent%i/FirstParticle/pSpectraKa", i), "pSpectraKa", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_1_perMult[2212] = registry.add(Form("Cent%i/FirstParticle/pSpectraPr", i), "pSpectraPr", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_1_perMult[1000010020] = registry.add(Form("Cent%i/FirstParticle/pSpectraDe", i), "pSpectraDe", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_1_perMult[0] = registry.add(Form("Cent%i/FirstParticle/pSpectraAll", i), "pSpectrAll", kTH1F, {{100, 0., 5., "p"}}); + + DCA_histos_1.push_back(std::move(DCA_histos_1_perMult)); + Purity_histos_1.push_back(std::move(Purity_histos_1_perMult)); + + if (!IsIdentical) { + std::map> DCA_histos_2_perMult; + DCA_histos_2_perMult[0] = registry.add(Form("Cent%i/SecondParticle/dcaxyz_vs_pt_primary", i), "dcaxyz_vs_pt_primary", kTH3F, {{100, 0., 5., "pt"}, {200, -1., 1., "DCA_XY(pt) primary"}, {200, -1., 1., "DCA_Z(pt) primary"}}); + DCA_histos_2_perMult[1] = registry.add(Form("Cent%i/SecondParticle/dcaxyz_vs_pt_weakdecay", i), "dcaxyz_vs_pt_weakdecay", kTH3F, {{100, 0., 5., "pt"}, {200, -1., 1., "DCA_XY(pt) weakdecay"}, {200, -1., 1., "DCA_Z(pt) weakdecay"}}); + DCA_histos_2_perMult[2] = registry.add(Form("Cent%i/SecondParticle/dcaxyz_vs_pt_material", i), "dcaxyz_vs_pt_material", kTH3F, {{100, 0., 5., "pt"}, {200, -1., 1., "DCA_XY(pt) material"}, {200, -1., 1., "DCA_Z(pt) material"}}); + + std::map> Purity_histos_2_perMult; + Purity_histos_2_perMult[11] = registry.add(Form("Cent%i/SecondParticle/pSpectraEl", i), "pSpectraEl", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_2_perMult[13] = registry.add(Form("Cent%i/SecondParticle/pSpectraMu", i), "pSpectraMu", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_2_perMult[211] = registry.add(Form("Cent%i/SecondParticle/pSpectraPi", i), "pSpectraPi", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_2_perMult[321] = registry.add(Form("Cent%i/SecondParticle/pSpectraKa", i), "pSpectraKa", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_2_perMult[2212] = registry.add(Form("Cent%i/SecondParticle/pSpectraPr", i), "pSpectraPr", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_2_perMult[1000010020] = registry.add(Form("Cent%i/SecondParticle/pSpectraDe", i), "pSpectraDe", kTH1F, {{100, 0., 5., "p"}}); + Purity_histos_2_perMult[0] = registry.add(Form("Cent%i/SecondParticle/pSpectraAll", i), "pSpectrAll", kTH1F, {{100, 0., 5., "p"}}); + + DCA_histos_2.push_back(std::move(DCA_histos_2_perMult)); + Purity_histos_2.push_back(std::move(Purity_histos_2_perMult)); + } - registry.add("FirstParticle/pSpectraEl", "pSpectraEl", kTH1F, {{100, 0., 5., "p"}}); - registry.add("FirstParticle/pSpectraMu", "pSpectraMu", kTH1F, {{100, 0., 5., "p"}}); - registry.add("FirstParticle/pSpectraPi", "pSpectraPi", kTH1F, {{100, 0., 5., "p"}}); - registry.add("FirstParticle/pSpectraKa", "pSpectraKa", kTH1F, {{100, 0., 5., "p"}}); - registry.add("FirstParticle/pSpectraPr", "pSpectraPr", kTH1F, {{100, 0., 5., "p"}}); - registry.add("FirstParticle/pSpectraDe", "pSpectraDe", kTH1F, {{100, 0., 5., "p"}}); - registry.add("FirstParticle/pSpectraAll", "pSpectrAll", kTH1F, {{100, 0., 5., "p"}}); + std::vector> kThistos_perMult; + std::vector> Resolution_histos_perMult; + std::vector> DoubleTrack_SE_histos_perMult; + std::vector> DoubleTrack_ME_histos_perMult; + + for (int j = 0; j < _kTbins.value.size() - 1; j++) { + auto kT_tmp = registry.add(Form("Cent%i/kT_cent%i_kT%i", i, i, j), Form("kT_cent%i_kT%i", i, j), kTH1F, {{500, 0., 5., "kT"}}); + auto Res_tmp = registry.add(Form("Cent%i/ResolutionMatrix_cent%i_kT%i", i, i, j), Form("ResolutionMatrix_rec(gen)_cent%i_kT%i", i, j), kTH2F, {{CFkStarBinning, "k*_gen (GeV/c)"}, {CFkStarBinning, "k*_rec (GeV/c)"}}); + auto DblTrk_SE_tmp = registry.add(Form("Cent%i/DoubleTrackEffects_SE_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_SE_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); + auto DblTrk_ME_tmp = registry.add(Form("Cent%i/DoubleTrackEffects_ME_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_ME_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); + kThistos_perMult.push_back(std::move(kT_tmp)); + Resolution_histos_perMult.push_back(std::move(Res_tmp)); + DoubleTrack_SE_histos_perMult.push_back(std::move(DblTrk_SE_tmp)); + DoubleTrack_ME_histos_perMult.push_back(std::move(DblTrk_ME_tmp)); + } - if (!IsIdentical) { - registry.add("SecondParticle/dcaxyz_vs_pt_primary", "dcaxyz_vs_pt_primary", kTH3F, {{100, 0., 5., "pt"}, {200, -1., 1., "DCA_XY(pt) primary"}, {200, -1., 1., "DCA_Z(pt) primary"}}); - registry.add("SecondParticle/dcaxyz_vs_pt_weakdecay", "dcaxyz_vs_pt_weakdecay", kTH3F, {{100, 0., 5., "pt"}, {200, -1., 1., "DCA_XY(pt) weakdecay"}, {200, -1., 1., "DCA_Z(pt) weakdecay"}}); - registry.add("SecondParticle/dcaxyz_vs_pt_material", "dcaxyz_vs_pt_material", kTH3F, {{100, 0., 5., "pt"}, {200, -1., 1., "DCA_XY(pt) material"}, {200, -1., 1., "DCA_Z(pt) material"}}); - - registry.add("SecondParticle/pSpectraEl", "pSpectraEl", kTH1F, {{100, 0., 5., "p"}}); - registry.add("SecondParticle/pSpectraMu", "pSpectraMu", kTH1F, {{100, 0., 5., "p"}}); - registry.add("SecondParticle/pSpectraPi", "pSpectraPi", kTH1F, {{100, 0., 5., "p"}}); - registry.add("SecondParticle/pSpectraKa", "pSpectraKa", kTH1F, {{100, 0., 5., "p"}}); - registry.add("SecondParticle/pSpectraPr", "pSpectraPr", kTH1F, {{100, 0., 5., "p"}}); - registry.add("SecondParticle/pSpectraDe", "pSpectraDe", kTH1F, {{100, 0., 5., "p"}}); - registry.add("SecondParticle/pSpectraAll", "pSpectrAll", kTH1F, {{100, 0., 5., "p"}}); + kThistos.push_back(std::move(kThistos_perMult)); + Resolution_histos.push_back(std::move(Resolution_histos_perMult)); + DoubleTrack_SE_histos.push_back(std::move(DoubleTrack_SE_histos_perMult)); + DoubleTrack_ME_histos.push_back(std::move(DoubleTrack_ME_histos_perMult)); } - - registry.add("ResolutionMatrix", "ResolutionMatrix_rec(gen)", kTH2F, {{CFkStarBinning, "k*_gen (GeV/c)"}, {CFkStarBinning, "k*_rec (GeV/c)"}}); - registry.add("DoubleTrackEffects", "DoubleTrackEffects_deta(dphi*)", kTH2F, {{200, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); } template - void fillEtaPhi(Type const& tracks) + void fillEtaPhi(Type const& tracks, int centBin) { // template for particles from the same collision identical for (int ii = 0; ii < tracks.size(); ii++) { // nested loop for all the combinations for (int iii = ii + 1; iii < tracks.size(); iii++) { Pair->SetPair(tracks[ii], tracks[iii]); + float pair_kT = Pair->GetKt(); + + if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) + continue; + + int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin >= 0) { + if (kTbin > DoubleTrack_SE_histos[centBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); + } else { + LOGF(fatal, "kTbin value obtained for a pair is less than 0"); + } - registry.fill(HIST("DoubleTrackEffects"), Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + kThistos[centBin][kTbin]->Fill(pair_kT); + DoubleTrack_SE_histos[centBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); Pair->ResetPair(); } } } template - void fillEtaPhi(Type const& tracks1, Type const& tracks2) + void fillEtaPhi(Type const& tracks1, Type const& tracks2, int centBin) { // template for particles from the same collision non-identical for (auto ii : tracks1) { for (auto iii : tracks2) { Pair->SetPair(ii, iii); + float pair_kT = Pair->GetKt(); - registry.fill(HIST("DoubleTrackEffects"), Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) + continue; + + int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin >= 0) { + if (kTbin > DoubleTrack_SE_histos[centBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); + } else { + LOGF(fatal, "kTbin value obtained for a pair is less than 0"); + } + + kThistos[centBin][kTbin]->Fill(pair_kT); + DoubleTrack_SE_histos[centBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); Pair->ResetPair(); } } } template - void fillResMatrix(Type const& tracks1, Type const& tracks2) + void fillResMatrix(Type const& tracks1, Type const& tracks2, int centBin) { // template for ME for (auto ii : tracks1) { for (auto iii : tracks2) { Pair->SetPair(ii, iii); + float pair_kT = Pair->GetKt(); + + if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) + continue; + + int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin >= 0) { + if (kTbin > Resolution_histos[centBin].size() || kTbin > DoubleTrack_ME_histos[centBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); + } else { + LOGF(fatal, "kTbin value obtained for a pair is less than 0"); + } TLorentzVector first4momentumGen; first4momentumGen.SetPtEtaPhiM(ii->pt_MC(), ii->eta_MC(), ii->phi_MC(), particle_mass(_particlePDG_1)); TLorentzVector second4momentumGen; second4momentumGen.SetPtEtaPhiM(iii->pt_MC(), iii->eta_MC(), iii->phi_MC(), particle_mass(_particlePDG_2)); - registry.fill(HIST("ResolutionMatrix"), o2::aod::singletrackselector::GetKstarFrom4vectors(first4momentumGen, second4momentumGen, IsIdentical), Pair->GetKstar()); + Resolution_histos[centBin][kTbin]->Fill(o2::aod::singletrackselector::GetKstarFrom4vectors(first4momentumGen, second4momentumGen, IsIdentical), Pair->GetKstar()); + DoubleTrack_ME_histos[centBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); Pair->ResetPair(); } } @@ -213,109 +295,61 @@ struct FemtoCorrelationsMC { continue; if (track.tpcFractionSharedCls() > _tpcFractionSharedCls || track.itsNCls() < _itsNCls) continue; + if (track.template singleCollSel_as>().multPerc() < *_centBins.value.begin() || track.template singleCollSel_as>().multPerc() >= *(_centBins.value.end() - 1)) + continue; + + int centBin = o2::aod::singletrackselector::getBinIndex(track.template singleCollSel_as>().multPerc(), _centBins); if (track.sign() == _sign_1 && (track.p() < _PIDtrshld_1 ? o2::aod::singletrackselector::TPCselection(track, TPCcuts_1) : o2::aod::singletrackselector::TOFselection(track, TOFcuts_1, _tpcNSigmaResidual_1))) { + trackOrigin = track.origin(); - switch (trackOrigin) { - case 0: - registry.fill(HIST("FirstParticle/dcaxyz_vs_pt_primary"), track.pt(), track.dcaXY(), track.dcaZ()); - break; - case 1: - registry.fill(HIST("FirstParticle/dcaxyz_vs_pt_weakdecay"), track.pt(), track.dcaXY(), track.dcaZ()); - break; - case 2: - registry.fill(HIST("FirstParticle/dcaxyz_vs_pt_material"), track.pt(), track.dcaXY(), track.dcaZ()); - break; - } + + if(trackOrigin > -1 && trackOrigin < 3) DCA_histos_1[centBin][track.origin()]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); + if (abs(track.dcaXY()) > _dcaXY || abs(track.dcaZ()) > _dcaZ) continue; - selectedtracks_1[track.singleCollSelId()].push_back(std::make_shared(track)); // filling the map: eventID <-> selected particles1 + trackPDG = track.pdgCode(); - trackPDG = abs(track.pdgCode()); - - registry.fill(HIST("FirstParticle/pSpectraAll"), track.p_MC()); - - switch (trackPDG) { - case 11: - registry.fill(HIST("FirstParticle/pSpectraEl"), track.p_MC()); - break; - case 13: - registry.fill(HIST("FirstParticle/pSpectraMu"), track.p_MC()); - break; - case 211: - registry.fill(HIST("FirstParticle/pSpectraPi"), track.p_MC()); - break; - case 321: - registry.fill(HIST("FirstParticle/pSpectraKa"), track.p_MC()); - break; - case 2212: - registry.fill(HIST("FirstParticle/pSpectraPr"), track.p_MC()); - break; - case 1000010020: - registry.fill(HIST("FirstParticle/pSpectraDe"), track.p_MC()); - break; - } + Purity_histos_1[centBin][0]->Fill(track.p()); + if(trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) Purity_histos_1[centBin][trackPDG]->Fill(track.p()); + + selectedtracks_1[track.singleCollSelId()].push_back(std::make_shared(track)); // filling the map: eventID <-> selected particles1 } if (IsIdentical) { continue; } else if (track.sign() != _sign_2 && !TOFselection(track, std::make_pair(_particlePDGtoReject, _rejectWithinNsigmaTOF)) && (track.p() < _PIDtrshld_2 ? o2::aod::singletrackselector::TPCselection(track, TPCcuts_2) : o2::aod::singletrackselector::TOFselection(track, TOFcuts_2, _tpcNSigmaResidual_2))) { // filling the map: eventID <-> selected particles2 if (see condition above ^) + trackOrigin = track.origin(); - switch (trackOrigin) { - case 0: - registry.fill(HIST("SecondParticle/dcaxyz_vs_pt_primary"), track.pt(), track.dcaXY(), track.dcaZ()); - break; - case 1: - registry.fill(HIST("SecondParticle/dcaxyz_vs_pt_weakdecay"), track.pt(), track.dcaXY(), track.dcaZ()); - break; - case 2: - registry.fill(HIST("SecondParticle/dcaxyz_vs_pt_material"), track.pt(), track.dcaXY(), track.dcaZ()); - break; - } + + if(trackOrigin > -1 && trackOrigin < 3) DCA_histos_2[centBin][track.origin()]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); + if (abs(track.dcaXY()) > _dcaXY || abs(track.dcaZ()) > _dcaZ) continue; - selectedtracks_2[track.singleCollSelId()].push_back(std::make_shared(track)); // filling the map: eventID <-> selected particles2 + Purity_histos_2[centBin][0]->Fill(track.p()); + if(trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) Purity_histos_2[centBin][trackPDG]->Fill(track.p()); - trackPDG = abs(track.pdgCode()); - - registry.fill(HIST("SecondParticle/pSpectraAll"), track.p_MC()); - - switch (trackPDG) { - case 11: - registry.fill(HIST("SecondParticle/pSpectraEl"), track.p_MC()); - break; - case 13: - registry.fill(HIST("SecondParticle/pSpectraMu"), track.p_MC()); - break; - case 211: - registry.fill(HIST("SecondParticle/pSpectraPi"), track.p_MC()); - break; - case 321: - registry.fill(HIST("SecondParticle/pSpectraKa"), track.p_MC()); - break; - case 2212: - registry.fill(HIST("SecondParticle/pSpectraPr"), track.p_MC()); - break; - case 1000010020: - registry.fill(HIST("SecondParticle/pSpectraDe"), track.p_MC()); - break; - } + selectedtracks_2[track.singleCollSelId()].push_back(std::make_shared(track)); // filling the map: eventID <-> selected particles2 } } for (auto collision : collisions) { + if (collision.multPerc() < *_centBins.value.begin() || collision.multPerc() >= *(_centBins.value.end() - 1)) + continue; + if (selectedtracks_1.find(collision.globalIndex()) == selectedtracks_1.end()) { if (IsIdentical) continue; else if (selectedtracks_2.find(collision.globalIndex()) == selectedtracks_2.end()) continue; } + int vertexBinToMix = std::floor((collision.posZ() + _vertexZ) / (2 * _vertexZ / _vertexNbinsToMix)); - int centBinToMix = std::floor(collision.multPerc() / (100.0 / _multNsubBins)); + float centBinToMix = o2::aod::singletrackselector::getBinIndex(collision.multPerc(), _centBins, _multNsubBins); - mixbins[std::pair{vertexBinToMix, centBinToMix}].push_back(std::make_shared(collision)); + mixbins[std::pair{vertexBinToMix, centBinToMix}].push_back(std::make_shared(collision)); } //====================================== filling deta(dphi*) & res. matrix starts here ====================================== @@ -331,14 +365,16 @@ struct FemtoCorrelationsMC { Pair->SetMagField1(col1->magField()); Pair->SetMagField2(col1->magField()); - fillEtaPhi(selectedtracks_1[col1->index()]); // filling deta(dphi*) -- SE identical + int centBin = std::floor((i->first).second); + + fillEtaPhi(selectedtracks_1[col1->index()], centBin); // filling deta(dphi*) -- SE identical for (int indx2 = indx1 + 1; indx2 < (i->second).size(); indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin auto col2 = (i->second)[indx2]; Pair->SetMagField2(col2->magField()); - fillResMatrix(selectedtracks_1[col1->index()], selectedtracks_1[col2->index()]); // filling res. matrix -- ME identical + fillResMatrix(selectedtracks_1[col1->index()], selectedtracks_1[col2->index()], centBin); // filling res. matrix -- ME identical } } } @@ -354,14 +390,16 @@ struct FemtoCorrelationsMC { Pair->SetMagField1(col1->magField()); Pair->SetMagField2(col1->magField()); - fillEtaPhi(selectedtracks_1[col1->index()], selectedtracks_2[col1->index()]); // filling deta(dphi*) -- SE non-identical + int centBin = std::floor((i->first).second); + + fillEtaPhi(selectedtracks_1[col1->index()], selectedtracks_2[col1->index()], centBin); // filling deta(dphi*) -- SE non-identical for (int indx2 = indx1 + 1; indx2 < (i->second).size(); indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin auto col2 = (i->second)[indx2]; Pair->SetMagField2(col2->magField()); - fillResMatrix(selectedtracks_1[col1->index()], selectedtracks_2[col2->index()]); // filling res. matrix -- ME non-identical + fillResMatrix(selectedtracks_1[col1->index()], selectedtracks_2[col2->index()], centBin); // filling res. matrix -- ME non-identical } } } diff --git a/PWGCF/Femto3D/Tasks/femto3dQA.cxx b/PWGCF/Femto3D/Tasks/femto3dQA.cxx index cb4a6f5680a..e695c8d5a5a 100755 --- a/PWGCF/Femto3D/Tasks/femto3dQA.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dQA.cxx @@ -40,6 +40,10 @@ struct QAHistograms { /// Construct a registry object with direct declaration HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject}; + Configurable _removeSameBunchPileup{"removeSameBunchPileup", false, ""}; + Configurable _requestGoodZvtxFT0vsPV{"requestGoodZvtxFT0vsPV", false, ""}; + Configurable _requestVertexITSTPC{"requestVertexITSTPC", false, ""}; + Configurable _sign{"sign", 1, "sign of a track"}; Configurable _vertexZ{"VertexZ", 10.0, "abs vertexZ value limit"}; Configurable _min_P{"min_P", 0.0, "lower mometum limit"}; @@ -62,6 +66,8 @@ struct QAHistograms { Configurable _particlePDGtoReject{"particlePDGtoReject", 211, "PDG codes of perticles that will be rejected with TOF (only pion, kaon, proton and deurton are supported now)"}; Configurable> _rejectWithinNsigmaTOF{"rejectWithinNsigmaTOF", std::vector{-0.0f, 0.0f}, "TOF rejection Nsigma range for particles specified with PDG to be rejected"}; + Configurable> _centCut{"centCut", std::pair{0.f, 100.f}, "[min., max.] centrality range to keep tracks within"}; + std::pair> TPCcuts; std::pair> TOFcuts; @@ -132,20 +138,39 @@ struct QAHistograms { } registry.add("posZ", "posZ", kTH1F, {{300, -16., 16., "posZ"}}); - registry.add("mult", "mult", kTH1F, {{500, 0., 500., "mult"}}); + registry.add("mult", "mult", kTH1F, {{5001, -0.5, 5000.5, "mult."}}); } template void fillHistograms(ColsType const& collisions, TracksType const& tracks) { for (auto& collision : collisions) { + if (_removeSameBunchPileup && !collision.isNoSameBunchPileup()) + continue; + if (_requestGoodZvtxFT0vsPV && !collision.isGoodZvtxFT0vsPV()) + continue; + if (_requestVertexITSTPC && !collision.isVertexITSTPC()) + continue; + if (collision.multPerc() < _centCut.value.first || collision.multPerc() >= _centCut.value.second) + continue; + registry.fill(HIST("posZ"), collision.posZ()); registry.fill(HIST("mult"), collision.mult()); } for (auto& track : tracks) { + + if (_removeSameBunchPileup && !track.template singleCollSel_as().isNoSameBunchPileup()) + continue; + if (_requestGoodZvtxFT0vsPV && !track.template singleCollSel_as().isGoodZvtxFT0vsPV()) + continue; + if (_requestVertexITSTPC && !track.template singleCollSel_as().isVertexITSTPC()) + continue; + if (abs(track.template singleCollSel_as().posZ()) > _vertexZ) continue; + if (track.template singleCollSel_as().multPerc() < _centCut.value.first || track.template singleCollSel_as().multPerc() >= _centCut.value.second) + continue; if ((track.tpcFractionSharedCls()) > _tpcFractionSharedCls || (track.itsNCls()) < _itsNCls) continue; @@ -203,13 +228,13 @@ struct QAHistograms { } } - void processDefault(soa::Filtered const& collisions, soa::Filtered const& tracks) + void processDefault(soa::Filtered> const& collisions, soa::Filtered const& tracks) { fillHistograms(collisions, tracks); } PROCESS_SWITCH(QAHistograms, processDefault, "process default", true); - void processExtra(soa::Filtered const& collisions, soa::Filtered> const& tracks) + void processExtra(soa::Filtered> const& collisions, soa::Filtered> const& tracks) { fillHistograms(collisions, tracks); } From 06d357146b1550f4ef018439312d7fddb1d030fc Mon Sep 17 00:00:00 2001 From: Gleb Romanenko Date: Fri, 19 Apr 2024 18:50:48 +0200 Subject: [PATCH 2/4] fixing formatting --- .../TableProducer/singleTrackSelector.cxx | 3 ++- .../singleTrackSelectorExtra.cxx | 7 +++---- PWGCF/Femto3D/Tasks/femto3dPairTask.cxx | 19 ++++++++++++------- PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx | 12 ++++++++---- 4 files changed, 25 insertions(+), 16 deletions(-) diff --git a/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx b/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx index 3f52e009cf4..800cf120a60 100755 --- a/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx +++ b/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx @@ -171,7 +171,8 @@ struct singleTrackSelector { for (auto ii : particlesToKeep) if (o2::aod::singletrackselector::TPCselection(track, std::make_pair(ii, keepWithinNsigmaTPC))) { - if(track.p() > _pRemoveTofOutOfRange && !o2::aod::singletrackselector::TOFselection(track, std::make_pair(ii, std::vector{-10.0, 10.0}), 10.0)) continue; + if (track.p() > _pRemoveTofOutOfRange && !o2::aod::singletrackselector::TOFselection(track, std::make_pair(ii, std::vector{-10.0, 10.0}), 10.0)) + continue; tableRow(tableRowColl.lastIndex(), track.p(), diff --git a/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx b/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx index f04f9eaced0..eec6197c98f 100755 --- a/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx +++ b/PWGCF/Femto3D/TableProducer/singleTrackSelectorExtra.cxx @@ -13,7 +13,6 @@ /// \author Sofia Tomassini, Gleb Romanenko, Nicolò Jacazio /// \since 31 May 2023 - // this task produces a dummy "SingleCollExtras" table that is now required in the analysis tasks. Needed to have a compatibility with old der. data #include @@ -37,9 +36,9 @@ struct singleTrackSelectorDummy { void process(aod::SingleCollSels::iterator const&) { - tableRowCollExtra(true, - true, - true); + tableRowCollExtra(true, + true, + true); } }; diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx index 396eec4f7ff..aa0873f3503 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx @@ -221,7 +221,7 @@ struct FemtoCorrelations { for (int j = 0; j < _kTbins.value.size() - 1; j++) { auto hDblTrk_SE = registry.add(Form("Cent%i/DoubleTrackEffects_SE_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_SE_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); auto hDblTrk_ME = registry.add(Form("Cent%i/DoubleTrackEffects_ME_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_ME_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); - + DoubleTrack_SE_histos_perMult.push_back(std::move(hDblTrk_SE)); DoubleTrack_ME_histos_perMult.push_back(std::move(hDblTrk_ME)); } @@ -271,8 +271,10 @@ struct FemtoCorrelations { LOGF(fatal, "kTbin value obtained for a pair is less than 0"); } - if (_fillDetaDphi) DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); - if (Pair->IsClosePair(_deta, _dphi, _radiusTPC)) continue; + if (_fillDetaDphi) + DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + if (Pair->IsClosePair(_deta, _dphi, _radiusTPC)) + continue; kThistos[multBin][kTbin]->Fill(pair_kT); mThistos[multBin][kTbin]->Fill(Pair->GetMt()); // test @@ -320,11 +322,14 @@ struct FemtoCorrelations { } if (_fillDetaDphi) { - if (!SE_or_ME) DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); - else DoubleTrack_ME_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + if (!SE_or_ME) + DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); + else + DoubleTrack_ME_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); } - if (Pair->IsClosePair(_deta, _dphi, _radiusTPC)) continue; + if (Pair->IsClosePair(_deta, _dphi, _radiusTPC)) + continue; if (!SE_or_ME) { SEhistos_1D[multBin][kTbin]->Fill(Pair->GetKstar()); @@ -349,7 +354,7 @@ struct FemtoCorrelations { Pair->ResetPair(); } } - } + } void process(soa::Filtered const& collisions, soa::Filtered const& tracks) { diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx index 7c26c74612f..3768972019d 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx @@ -304,7 +304,8 @@ struct FemtoCorrelationsMC { trackOrigin = track.origin(); - if(trackOrigin > -1 && trackOrigin < 3) DCA_histos_1[centBin][track.origin()]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); + if (trackOrigin > -1 && trackOrigin < 3) + DCA_histos_1[centBin][track.origin()]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); if (abs(track.dcaXY()) > _dcaXY || abs(track.dcaZ()) > _dcaZ) continue; @@ -312,7 +313,8 @@ struct FemtoCorrelationsMC { trackPDG = track.pdgCode(); Purity_histos_1[centBin][0]->Fill(track.p()); - if(trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) Purity_histos_1[centBin][trackPDG]->Fill(track.p()); + if (trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) + Purity_histos_1[centBin][trackPDG]->Fill(track.p()); selectedtracks_1[track.singleCollSelId()].push_back(std::make_shared(track)); // filling the map: eventID <-> selected particles1 } @@ -323,13 +325,15 @@ struct FemtoCorrelationsMC { trackOrigin = track.origin(); - if(trackOrigin > -1 && trackOrigin < 3) DCA_histos_2[centBin][track.origin()]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); + if (trackOrigin > -1 && trackOrigin < 3) + DCA_histos_2[centBin][track.origin()]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); if (abs(track.dcaXY()) > _dcaXY || abs(track.dcaZ()) > _dcaZ) continue; Purity_histos_2[centBin][0]->Fill(track.p()); - if(trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) Purity_histos_2[centBin][trackPDG]->Fill(track.p()); + if (trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) + Purity_histos_2[centBin][trackPDG]->Fill(track.p()); selectedtracks_2[track.singleCollSelId()].push_back(std::make_shared(track)); // filling the map: eventID <-> selected particles2 } From 23ca772d3fbe00463b889a90ddab0855b476cbb7 Mon Sep 17 00:00:00 2001 From: Gleb Romanenko Date: Sat, 20 Apr 2024 20:49:10 +0200 Subject: [PATCH 3/4] fixing warnings --- .../TableProducer/singleTrackSelector.cxx | 2 +- PWGCF/Femto3D/Tasks/femto3dPairTask.cxx | 82 ++++++++----------- PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx | 64 ++++++--------- 3 files changed, 61 insertions(+), 87 deletions(-) diff --git a/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx b/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx index 800cf120a60..f65cb137188 100755 --- a/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx +++ b/PWGCF/Femto3D/TableProducer/singleTrackSelector.cxx @@ -102,7 +102,7 @@ struct singleTrackSelector { std::vector particlesToKeep; std::vector particlesToReject; - void init(InitContext& context) + void init(InitContext&) { particlesToKeep = _particlesToKeep; diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx index aa0873f3503..efbb623befa 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx @@ -171,7 +171,7 @@ struct FemtoCorrelations { TPCcuts_2 = std::make_pair(_particlePDG_2, _tpcNSigma_2); TOFcuts_2 = std::make_pair(_particlePDG_2, _tofNSigma_2); - for (int i = 0; i < _centBins.value.size() - 1; i++) { + for (unsigned int i = 0; i < _centBins.value.size() - 1; i++) { std::vector> SEperMult_1D; std::vector> MEperMult_1D; std::vector> kTperMult; @@ -180,7 +180,7 @@ struct FemtoCorrelations { auto hMult = registry.add(Form("Cent%i/TPCMult_cent%i", i, i), Form("TPCMult_cent%i", i), kTH1F, {{5001, -0.5, 5000.5, "Mult."}}); MultHistos.push_back(std::move(hMult)); - for (int j = 0; j < _kTbins.value.size() - 1; j++) { + for (unsigned int j = 0; j < _kTbins.value.size() - 1; j++) { auto hSE_1D = registry.add(Form("Cent%i/SE_1D_cent%i_kT%i", i, i, j), Form("SE_1D_cent%i_kT%i", i, j), kTH1F, {{CFkStarBinning, "k* (GeV/c)"}}); auto hME_1D = registry.add(Form("Cent%i/ME_1D_cent%i_kT%i", i, i, j), Form("ME_1D_cent%i_kT%i", i, j), kTH1F, {{CFkStarBinning, "k* (GeV/c)"}}); auto hkT = registry.add(Form("Cent%i/kT_cent%i_kT%i", i, i, j), Form("kT_cent%i_kT%i", i, j), kTH1F, {{500, 0., 5., "kT"}}); @@ -201,7 +201,7 @@ struct FemtoCorrelations { std::vector> MEperMult_3D; std::vector> qLCMSvskStarperMult; - for (int j = 0; j < _kTbins.value.size() - 1; j++) { + for (unsigned int j = 0; j < _kTbins.value.size() - 1; j++) { auto hSE_3D = registry.add(Form("Cent%i/SE_3D_cent%i_kT%i", i, i, j), Form("SE_3D_cent%i_kT%i", i, j), kTH3F, {{CF3DqLCMSBinning, "q_out (GeV/c)"}, {CF3DqLCMSBinning, "q_side (GeV/c)"}, {CF3DqLCMSBinning, "q_long (GeV/c)"}}); auto hME_3D = registry.add(Form("Cent%i/ME_3D_cent%i_kT%i", i, i, j), Form("ME_3D_cent%i_kT%i", i, j), kTH3F, {{CF3DqLCMSBinning, "q_out (GeV/c)"}, {CF3DqLCMSBinning, "q_side (GeV/c)"}, {CF3DqLCMSBinning, "q_long (GeV/c)"}}); auto hqLCMSvskStar = registry.add(Form("Cent%i/qLCMSvskStar_cent%i_kT%i", i, i, j), Form("qLCMSvskStar_cent%i_kT%i", i, j), kTH3F, {{CF3DqLCMSBinning, "q_out (GeV/c)"}, {CF3DqLCMSBinning, "q_side (GeV/c)"}, {CF3DqLCMSBinning, "q_long (GeV/c)"}}); @@ -218,7 +218,7 @@ struct FemtoCorrelations { std::vector> DoubleTrack_SE_histos_perMult; std::vector> DoubleTrack_ME_histos_perMult; - for (int j = 0; j < _kTbins.value.size() - 1; j++) { + for (unsigned int j = 0; j < _kTbins.value.size() - 1; j++) { auto hDblTrk_SE = registry.add(Form("Cent%i/DoubleTrackEffects_SE_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_SE_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); auto hDblTrk_ME = registry.add(Form("Cent%i/DoubleTrackEffects_ME_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_ME_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); @@ -241,19 +241,15 @@ struct FemtoCorrelations { } template - void mixTracks(Type const& tracks, int multBin) + void mixTracks(Type const& tracks, unsigned int multBin) { // template for identical particles from the same collision - if (multBin >= 0) { - if (multBin > SEhistos_1D.size()) - LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (1D)"); - if (_fill3dCF && multBin > SEhistos_3D.size()) - LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (3D)"); - } else { - LOGF(fatal, "multBin value passed to the mixTracks function is less than 0"); - } + if (multBin > SEhistos_1D.size()) + LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (1D)"); + if (_fill3dCF && multBin > SEhistos_3D.size()) + LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (3D)"); - for (int ii = 0; ii < tracks.size(); ii++) { // nested loop for all the combinations - for (int iii = ii + 1; iii < tracks.size(); iii++) { + for (unsigned int ii = 0; ii < tracks.size(); ii++) { // nested loop for all the combinations + for (unsigned int iii = ii + 1; iii < tracks.size(); iii++) { Pair->SetPair(tracks[ii], tracks[iii]); float pair_kT = Pair->GetKt(); @@ -261,15 +257,11 @@ struct FemtoCorrelations { if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) continue; - int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); - if (kTbin >= 0) { - if (kTbin > SEhistos_1D[multBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (1D)"); - if (_fill3dCF && kTbin > SEhistos_3D[multBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (3D)"); - } else { - LOGF(fatal, "kTbin value obtained for a pair is less than 0"); - } + unsigned int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin > SEhistos_1D[multBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (1D)"); + if (_fill3dCF && kTbin > SEhistos_3D[multBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (3D)"); if (_fillDetaDphi) DoubleTrack_SE_histos[multBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); @@ -291,16 +283,12 @@ struct FemtoCorrelations { } template - void mixTracks(Type const& tracks1, Type const& tracks2, int multBin) + void mixTracks(Type const& tracks1, Type const& tracks2, unsigned int multBin) { // last value: 0 -- SE; 1 -- ME - if (multBin >= 0) { - if (multBin > SEhistos_1D.size()) - LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (1D)"); - if (_fill3dCF && multBin > SEhistos_3D.size()) - LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (3D)"); - } else { - LOGF(fatal, "multBin value passed to the mixTracks function is less than 0"); - } + if (multBin > SEhistos_1D.size()) + LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (1D)"); + if (_fill3dCF && multBin > SEhistos_3D.size()) + LOGF(fatal, "multBin value passed to the mixTracks function exceeds the configured number of Cent. bins (3D)"); for (auto ii : tracks1) { for (auto iii : tracks2) { @@ -311,15 +299,11 @@ struct FemtoCorrelations { if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) continue; - int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); - if (kTbin >= 0) { - if (kTbin > SEhistos_1D[multBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (1D)"); - if (_fill3dCF && kTbin > SEhistos_3D[multBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (3D)"); - } else { - LOGF(fatal, "kTbin value obtained for a pair is less than 0"); - } + unsigned int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin > SEhistos_1D[multBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (1D)"); + if (_fill3dCF && kTbin > SEhistos_3D[multBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins (3D)"); if (_fillDetaDphi) { if (!SE_or_ME) @@ -452,19 +436,19 @@ struct FemtoCorrelations { for (auto i = mixbins.begin(); i != mixbins.end(); i++) { // iterating over all vertex&mult bins int EvPerBin = (i->second).size(); - for (int indx1 = 0; indx1 < EvPerBin; indx1++) { // loop over all the events in each vertex&mult bin + for (unsigned int indx1 = 0; indx1 < EvPerBin; indx1++) { // loop over all the events in each vertex&mult bin auto col1 = (i->second)[indx1]; Pair->SetMagField1(col1->magField()); Pair->SetMagField2(col1->magField()); - int centBin = std::floor((i->first).second); + unsigned int centBin = std::floor((i->first).second); MultHistos[centBin]->Fill(col1->mult()); mixTracks(selectedtracks_1[col1->index()], centBin); // mixing SE identical - for (int indx2 = indx1 + 1; indx2 < EvPerBin; indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin + for (unsigned int indx2 = indx1 + 1; indx2 < EvPerBin; indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin if (_MEreductionFactor.value > 1) { std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); if ((mt() % (_MEreductionFactor.value + 1)) < _MEreductionFactor.value) @@ -482,21 +466,21 @@ struct FemtoCorrelations { } else { //====================================== mixing non-identical ====================================== for (auto i = mixbins.begin(); i != mixbins.end(); i++) { // iterating over all vertex&mult bins - int EvPerBin = (i->second).size(); + unsigned int EvPerBin = (i->second).size(); - for (int indx1 = 0; indx1 < EvPerBin; indx1++) { // loop over all the events in each vertex&mult bin + for (unsigned int indx1 = 0; indx1 < EvPerBin; indx1++) { // loop over all the events in each vertex&mult bin auto col1 = (i->second)[indx1]; Pair->SetMagField1(col1->magField()); Pair->SetMagField2(col1->magField()); - int centBin = std::floor((i->first).second); + unsigned int centBin = std::floor((i->first).second); MultHistos[centBin]->Fill(col1->mult()); mixTracks<0>(selectedtracks_1[col1->index()], selectedtracks_2[col1->index()], centBin); // mixing SE non-identical, in <> brackets: 0 -- SE; 1 -- ME - for (int indx2 = indx1 + 1; indx2 < EvPerBin; indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin + for (unsigned int indx2 = indx1 + 1; indx2 < EvPerBin; indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin if (_MEreductionFactor.value > 1) { std::mt19937 mt(std::chrono::steady_clock::now().time_since_epoch().count()); if (mt() % (_MEreductionFactor.value + 1) < _MEreductionFactor.value) diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx index 3768972019d..596c6039933 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTaskMC.cxx @@ -137,7 +137,7 @@ struct FemtoCorrelationsMC { TPCcuts_2 = std::make_pair(_particlePDG_2, _tpcNSigma_2); TOFcuts_2 = std::make_pair(_particlePDG_2, _tofNSigma_2); - for (int i = 0; i < _centBins.value.size() - 1; i++) { + for (unsigned int i = 0; i < _centBins.value.size() - 1; i++) { std::map> DCA_histos_1_perMult; DCA_histos_1_perMult[0] = registry.add(Form("Cent%i/FirstParticle/dcaxyz_vs_pt_primary", i), "dcaxyz_vs_pt_primary", kTH3F, {{100, 0., 5., "pt"}, {250, -1., 1., "DCA_XY(pt) primary"}, {250, -1., 1., "DCA_Z(pt) primary"}}); @@ -180,7 +180,7 @@ struct FemtoCorrelationsMC { std::vector> DoubleTrack_SE_histos_perMult; std::vector> DoubleTrack_ME_histos_perMult; - for (int j = 0; j < _kTbins.value.size() - 1; j++) { + for (unsigned int j = 0; j < _kTbins.value.size() - 1; j++) { auto kT_tmp = registry.add(Form("Cent%i/kT_cent%i_kT%i", i, i, j), Form("kT_cent%i_kT%i", i, j), kTH1F, {{500, 0., 5., "kT"}}); auto Res_tmp = registry.add(Form("Cent%i/ResolutionMatrix_cent%i_kT%i", i, i, j), Form("ResolutionMatrix_rec(gen)_cent%i_kT%i", i, j), kTH2F, {{CFkStarBinning, "k*_gen (GeV/c)"}, {CFkStarBinning, "k*_rec (GeV/c)"}}); auto DblTrk_SE_tmp = registry.add(Form("Cent%i/DoubleTrackEffects_SE_cent%i_kT%i", i, i, j), Form("DoubleTrackEffects_deta(dphi*)_SE_cent%i_kT%i", i, j), kTH2F, {{600, -M_PI, M_PI, "dphi*"}, {200, -0.5, 0.5, "deta"}}); @@ -199,10 +199,10 @@ struct FemtoCorrelationsMC { } template - void fillEtaPhi(Type const& tracks, int centBin) - { // template for particles from the same collision identical - for (int ii = 0; ii < tracks.size(); ii++) { // nested loop for all the combinations - for (int iii = ii + 1; iii < tracks.size(); iii++) { + void fillEtaPhi(Type const& tracks, unsigned int centBin) + { // template for particles from the same collision identical + for (unsigned int ii = 0; ii < tracks.size(); ii++) { // nested loop for all the combinations + for (unsigned int iii = ii + 1; iii < tracks.size(); iii++) { Pair->SetPair(tracks[ii], tracks[iii]); float pair_kT = Pair->GetKt(); @@ -210,13 +210,9 @@ struct FemtoCorrelationsMC { if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) continue; - int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); - if (kTbin >= 0) { - if (kTbin > DoubleTrack_SE_histos[centBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); - } else { - LOGF(fatal, "kTbin value obtained for a pair is less than 0"); - } + unsigned int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin > DoubleTrack_SE_histos[centBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); kThistos[centBin][kTbin]->Fill(pair_kT); DoubleTrack_SE_histos[centBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); @@ -226,7 +222,7 @@ struct FemtoCorrelationsMC { } template - void fillEtaPhi(Type const& tracks1, Type const& tracks2, int centBin) + void fillEtaPhi(Type const& tracks1, Type const& tracks2, unsigned int centBin) { // template for particles from the same collision non-identical for (auto ii : tracks1) { for (auto iii : tracks2) { @@ -237,13 +233,9 @@ struct FemtoCorrelationsMC { if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) continue; - int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); - if (kTbin >= 0) { - if (kTbin > DoubleTrack_SE_histos[centBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); - } else { - LOGF(fatal, "kTbin value obtained for a pair is less than 0"); - } + unsigned int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin > DoubleTrack_SE_histos[centBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); kThistos[centBin][kTbin]->Fill(pair_kT); DoubleTrack_SE_histos[centBin][kTbin]->Fill(Pair->GetPhiStarDiff(_radiusTPC), Pair->GetEtaDiff()); @@ -253,7 +245,7 @@ struct FemtoCorrelationsMC { } template - void fillResMatrix(Type const& tracks1, Type const& tracks2, int centBin) + void fillResMatrix(Type const& tracks1, Type const& tracks2, unsigned int centBin) { // template for ME for (auto ii : tracks1) { for (auto iii : tracks2) { @@ -263,13 +255,9 @@ struct FemtoCorrelationsMC { if (pair_kT < *_kTbins.value.begin() || pair_kT >= *(_kTbins.value.end() - 1)) continue; - int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); - if (kTbin >= 0) { - if (kTbin > Resolution_histos[centBin].size() || kTbin > DoubleTrack_ME_histos[centBin].size()) - LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); - } else { - LOGF(fatal, "kTbin value obtained for a pair is less than 0"); - } + unsigned int kTbin = o2::aod::singletrackselector::getBinIndex(pair_kT, _kTbins); + if (kTbin > Resolution_histos[centBin].size() || kTbin > DoubleTrack_ME_histos[centBin].size()) + LOGF(fatal, "kTbin value obtained for a pair exceeds the configured number of kT bins"); TLorentzVector first4momentumGen; first4momentumGen.SetPtEtaPhiM(ii->pt_MC(), ii->eta_MC(), ii->phi_MC(), particle_mass(_particlePDG_1)); @@ -298,7 +286,7 @@ struct FemtoCorrelationsMC { if (track.template singleCollSel_as>().multPerc() < *_centBins.value.begin() || track.template singleCollSel_as>().multPerc() >= *(_centBins.value.end() - 1)) continue; - int centBin = o2::aod::singletrackselector::getBinIndex(track.template singleCollSel_as>().multPerc(), _centBins); + unsigned int centBin = o2::aod::singletrackselector::getBinIndex(track.template singleCollSel_as>().multPerc(), _centBins); if (track.sign() == _sign_1 && (track.p() < _PIDtrshld_1 ? o2::aod::singletrackselector::TPCselection(track, TPCcuts_1) : o2::aod::singletrackselector::TOFselection(track, TOFcuts_1, _tpcNSigmaResidual_1))) { @@ -310,7 +298,7 @@ struct FemtoCorrelationsMC { if (abs(track.dcaXY()) > _dcaXY || abs(track.dcaZ()) > _dcaZ) continue; - trackPDG = track.pdgCode(); + trackPDG = abs(track.pdgCode()); Purity_histos_1[centBin][0]->Fill(track.p()); if (trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) @@ -331,6 +319,8 @@ struct FemtoCorrelationsMC { if (abs(track.dcaXY()) > _dcaXY || abs(track.dcaZ()) > _dcaZ) continue; + trackPDG = abs(track.pdgCode()); + Purity_histos_2[centBin][0]->Fill(track.p()); if (trackPDG == 11 || trackPDG == 13 || trackPDG == 211 || trackPDG == 321 || trackPDG == 2212 || trackPDG == 1000010020) Purity_histos_2[centBin][trackPDG]->Fill(track.p()); @@ -362,18 +352,18 @@ struct FemtoCorrelationsMC { for (auto i = mixbins.begin(); i != mixbins.end(); i++) { // iterating over all vertex&mult bins - for (int indx1 = 0; indx1 < (i->second).size(); indx1++) { // iterating over all selected collisions with selected tracks + for (unsigned int indx1 = 0; indx1 < (i->second).size(); indx1++) { // iterating over all selected collisions with selected tracks auto col1 = (i->second)[indx1]; Pair->SetMagField1(col1->magField()); Pair->SetMagField2(col1->magField()); - int centBin = std::floor((i->first).second); + unsigned int centBin = std::floor((i->first).second); fillEtaPhi(selectedtracks_1[col1->index()], centBin); // filling deta(dphi*) -- SE identical - for (int indx2 = indx1 + 1; indx2 < (i->second).size(); indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin + for (unsigned int indx2 = indx1 + 1; indx2 < (i->second).size(); indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin auto col2 = (i->second)[indx2]; @@ -387,18 +377,18 @@ struct FemtoCorrelationsMC { for (auto i = mixbins.begin(); i != mixbins.end(); i++) { // iterating over all vertex&mult bins - for (int indx1 = 0; indx1 < (i->second).size(); indx1++) { // iterating over all selected collisions with selected tracks1 + for (unsigned int indx1 = 0; indx1 < (i->second).size(); indx1++) { // iterating over all selected collisions with selected tracks1 auto col1 = (i->second)[indx1]; Pair->SetMagField1(col1->magField()); Pair->SetMagField2(col1->magField()); - int centBin = std::floor((i->first).second); + unsigned int centBin = std::floor((i->first).second); fillEtaPhi(selectedtracks_1[col1->index()], selectedtracks_2[col1->index()], centBin); // filling deta(dphi*) -- SE non-identical - for (int indx2 = indx1 + 1; indx2 < (i->second).size(); indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin + for (unsigned int indx2 = indx1 + 1; indx2 < (i->second).size(); indx2++) { // nested loop for all the combinations of collisions in a chosen mult/vertex bin auto col2 = (i->second)[indx2]; From bb4a8d5c45b2b3415797f9174022b3b42a76af6d Mon Sep 17 00:00:00 2001 From: Gleb Romanenko Date: Mon, 22 Apr 2024 01:40:04 +0200 Subject: [PATCH 4/4] fixing warnings x2 --- PWGCF/Femto3D/Core/femto3dPairTask.h | 4 ++-- PWGCF/Femto3D/Tasks/femto3dPairTask.cxx | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/PWGCF/Femto3D/Core/femto3dPairTask.h b/PWGCF/Femto3D/Core/femto3dPairTask.h index 43baeb6a9bc..2bda63575a8 100755 --- a/PWGCF/Femto3D/Core/femto3dPairTask.h +++ b/PWGCF/Femto3D/Core/femto3dPairTask.h @@ -42,8 +42,8 @@ namespace o2::aod::singletrackselector template Type getBinIndex(float const& value, std::vector const& binning, int const& NsubBins = 1) { - Type res = -100; - for (int i = 0; i < binning.size() - 1; i++) { + Type res = 10e6; + for (unsigned int i = 0; i < binning.size() - 1; i++) { if (value >= binning[i] && binning[i + 1] > value) { if (NsubBins < 2) { res = (Type)i; diff --git a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx index efbb623befa..6d9fd003c1b 100755 --- a/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx +++ b/PWGCF/Femto3D/Tasks/femto3dPairTask.cxx @@ -101,7 +101,7 @@ struct FemtoCorrelations { // P.P.S. the chosen way of optimizing the mixing midgt not be the correct one -- feel free to propose the right one! // P.P.P.S. choose wisely.... // P.P.P.P.S this way is still being testing i might be reconsidered; might change in the future, keep looking at the source code - Configurable _MEreductionFactor{"MEreductionFactor", 1, "only one (pseudo)randomly choosen event out per pair $value events will be processed and contribute to the final mixing (if < 1 -> all the possible event pairs (per vertex¢ bin) will be processed); implemented for the sake of efficiency; look at the source code;"}; + Configurable _MEreductionFactor{"MEreductionFactor", 1, "only one (pseudo)randomly choosen event out per pair $value events will be processed and contribute to the final mixing (if < 2 -> all the possible event pairs (per vertex¢ bin) will be processed); implemented for the sake of efficiency; look at the source code;"}; bool IsIdentical; @@ -434,7 +434,7 @@ struct FemtoCorrelations { if (IsIdentical) { //====================================== mixing identical ====================================== for (auto i = mixbins.begin(); i != mixbins.end(); i++) { // iterating over all vertex&mult bins - int EvPerBin = (i->second).size(); + unsigned int EvPerBin = (i->second).size(); for (unsigned int indx1 = 0; indx1 < EvPerBin; indx1++) { // loop over all the events in each vertex&mult bin