diff --git a/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx b/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx index 4de102546a0..f66844cb1e2 100644 --- a/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx +++ b/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx @@ -1263,7 +1263,7 @@ struct FemtoUniverseProducerTask { } } - template + template void fillTracksFullNsigma(TrackType const& tracks) { std::vector childIDs = {0, 0}; @@ -1318,6 +1318,12 @@ struct FemtoUniverseProducerTask { track.tofNSigmaDe(), track.tpcNSigmaEl(), track.tpcNSigmaPi(), track.tpcNSigmaKa(), track.tpcNSigmaPr(), track.tpcNSigmaDe()); + + if constexpr (fillITS) { + outputDebugITSParts(track.itsNSigmaEl(), track.itsNSigmaPi(), + track.itsNSigmaKa(), track.itsNSigmaPr(), + track.itsNSigmaDe()); + } } } @@ -2995,10 +3001,19 @@ struct FemtoUniverseProducerTask { soa::Filtered> const& tracksNsigma) { + + auto tracksWithItsPid = soa::Attach>, + aod::pidits::ITSNSigmaEl, aod::pidits::ITSNSigmaPi, aod::pidits::ITSNSigmaKa, aod::pidits::ITSNSigmaPr, aod::pidits::ITSNSigmaDe, aod::pidits::ITSNSigmaTr, aod::pidits::ITSNSigmaHe>(tracksNsigma); + // fill the tables const auto colcheck = fillCollisionsCentRun3(col); if (colcheck) { - fillTracksFullNsigma(tracksNsigma); + if (!confFillITSPid) { + fillTracksFullNsigma(tracksNsigma); + } else { + fillTracksFullNsigma(tracksWithItsPid); + } } } PROCESS_SWITCH(FemtoUniverseProducerTask, processTrackCentRun3DataFullNSigma, "Provide experimental data for Run 3 with centrality for track track with full NSigma information", false); diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackTrackMultKtExtended.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackTrackMultKtExtended.cxx index 24cd5ab2610..19b0e9255da 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackTrackMultKtExtended.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackTrackMultKtExtended.cxx @@ -74,39 +74,50 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { struct : o2::framework::ConfigurableGroup { Configurable isKaonNsigma{"isKaonNsigma", false, "Enable a strict cut selection for K+ and K-"}; Configurable isNsigmaRectangular{"isNsigmaRectangular", false, "Apply rectangular TPC and TOF nSigma cut, instead of a combined one (TPC confIsAddITSNsigma{"confIsAddITSNsigma", false, "Apply cuts on ITS nSigma"}; + + Configurable confNsigmaTPCITSCombined{"confNsigmaTPCITSCombined", 3.0f, "TPC and ITS nSigma (combined) for momentum < confTOFpMin"}; Configurable confNsigmaCombined{"confNsigmaCombined", 3.0f, "TPC and TOF Pion Sigma (combined) for momentum > confTOFpMin"}; Configurable confNsigmaTPC{"confNsigmaTPC", 3.0f, "TPC Pion Sigma for momentum < confTOFpMin"}; Configurable confNsigmaKaonRejection{"confNsigmaKaonRejection", 0, "In proton selection reject tracks with combined Nsigma confNsigmaPionRejection{"confNsigmaPionRejection", 0, "In proton selection reject tracks with combined Nsigma confTPCNsigmaKaonRejection{"confTPCNsigmaKaonRejection", 0, "In proton selection reject tracks with TPC Nsigma confTPCNsigmaPionRejection{"confTPCNsigmaPionRejection", 0, "In proton selection reject tracks with TPC Nsigma confTOFNsigmaKaonRejection{"confTOFNsigmaKaonRejection", 0, "In proton selection reject tracks with TOF Nsigma confTOFNsigmaPionRejection{"confTOFNsigmaPionRejection", 0, "In proton selection reject tracks with TOF Nsigma confTOFpMinRej{"confTOFpMinRej", 0.5f, "Min. momentum for which TOF is required for PID for rejection."}; + Configurable isNsigmaRejectionCombined{"isNsigmaRejectionCombined", true, "True: reject based on combined nSigma, false: reject based on TPC or TOF nSigma (confTOFpMinRej)."}; Configurable confTOFpMin{"confTOFpMin", 0.5f, "Min. momentum for which TOF is required for PID."}; Configurable confEtaMax{"confEtaMax", 0.8f, "Higher limit for |Eta| (the same for both particles)"}; Configurable confIsMC{"confIsMC", false, "Enable additional Histogramms in the case of a MonteCarlo Run"}; - Configurable> confTrkPIDnSigmaMax{"confTrkPIDnSigmaMax", std::vector{4.f, 3.f, 2.f}, "This configurable needs to be the same as the one used in the producer task"}; Configurable confUse3D{"confUse3D", false, "Enable three dimensional histogramms (to be used only for analysis with high statistics): k* vs mT vs multiplicity"}; Configurable confOnlyPrimaryMCPair{"confOnlyPrimaryMCPair", false, "Fill MC pair histograms only with primary particles"}; + Configurable confTrkTPCnclsMin{"confTrkTPCnclsMin", 0, "Min. TPC clusters (if stricter than in the Producer)"}; + Configurable confTrkTPCcRowsMin{"confTrkTPCcRowsMin", 0, "Min. TPC crossed rows (if stricter than in the Producer)"}; + Configurable confDcaXYCustom1FilterCut{"confDcaXYCustom1FilterCut", 10, "Value for [1] custom DCAxy cut -> |DCAxy| < [1] + [2]/pT (if stricter than in the Producer)"}; + Configurable confDcaXYCustom2FilterCut{"confDcaXYCustom2FilterCut", 0, "Value for [2] custom DCAxy cut -> |DCAxy| < [1] + [2]/pT (if stricter than in the Producer)"}; } twotracksconfigs; - struct : o2::framework::ConfigurableGroup { - Configurable confTrkTPCnclsMin{"confTrkTPCnclsMin", 0, "Min. TPC clusters"}; - Configurable confTrkTPCcRowsMin{"confTrkTPCcRowsMin", 0, "Min. TPC crossed rows"}; - Configurable confDcaXYCustom1FilterCut{"confDcaXYCustom1FilterCut", 10, "Value for [1] custom DCAxy cut -> |DCAxy| < [1] + [2]/pT"}; - Configurable confDcaXYCustom2FilterCut{"confDcaXYCustom2FilterCut", 0, "Value for [2] custom DCAxy cut -> |DCAxy| < [1] + [2]/pT"}; - } additionalcuts; - using FemtoFullParticles = soa::Join; + using FemtoFullParticlesITS = soa::Join; using FemtoFullSigmaParticles = soa::Join; + using FemtoFullSigmaParticlesITS = soa::Join; // Filters for selecting particles (both p1 and p2) Filter trackAdditionalfilter = (nabs(aod::femtouniverseparticle::eta) < twotracksconfigs.confEtaMax); // example filtering on configurable using FilteredFemtoFullParticles = soa::Filtered; using FilteredFemtoFullSigmaParticles = soa::Filtered; + using FilteredFemtoFullParticlesITS = soa::Filtered; + using FilteredFemtoFullSigmaParticlesITS = soa::Filtered; // using FilteredFemtoFullParticles = FemtoFullParticles; //if no filtering is applied uncomment this option SliceCache cache; Preslice perCol = aod::femtouniverseparticle::fdCollisionId; Preslice perColSigma = aod::femtouniverseparticle::fdCollisionId; + Preslice perColITS = aod::femtouniverseparticle::fdCollisionId; + Preslice perColSigmaITS = aod::femtouniverseparticle::fdCollisionId; using FemtoTruthParticles = soa::Filtered>; Preslice perColMCTruth = aod::femtouniverseparticle::fdCollisionId; @@ -114,6 +125,9 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { using FemtoRecoParticles = soa::Filtered>; Preslice perColMC = aod::femtouniverseparticle::fdCollisionId; + using FemtoRecoParticlesITS = soa::Filtered>; + Preslice perColMCITS = aod::femtouniverseparticle::fdCollisionId; + /// Particle 1 struct : o2::framework::ConfigurableGroup { Configurable confPDGCodePartOne{"confPDGCodePartOne", 211, "Particle 1 -- PDG code"}; @@ -123,11 +137,14 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { } trackonefilter; /// Partition for particle 1 - Partition partsOne = (nabs(aod::track::dcaXY) < (additionalcuts.confDcaXYCustom1FilterCut + additionalcuts.confDcaXYCustom2FilterCut / aod::femtouniverseparticle::pt)) && (aod::femtouniverseparticle::tpcNClsFound) > as(additionalcuts.confTrkTPCnclsMin) && (aod::femtouniverseparticle::tpcNClsCrossedRows) > as(additionalcuts.confTrkTPCcRowsMin) && (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt trackonefilter.confpLowPart1; + Partition partsOne = (nabs(aod::track::dcaXY) < (twotracksconfigs.confDcaXYCustom1FilterCut + twotracksconfigs.confDcaXYCustom2FilterCut / aod::femtouniverseparticle::pt)) && (aod::femtouniverseparticle::tpcNClsFound) > as(twotracksconfigs.confTrkTPCnclsMin) && (aod::femtouniverseparticle::tpcNClsCrossedRows) > as(twotracksconfigs.confTrkTPCcRowsMin) && (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt trackonefilter.confpLowPart1; + Partition partsOneITS = (nabs(aod::track::dcaXY) < (twotracksconfigs.confDcaXYCustom1FilterCut + twotracksconfigs.confDcaXYCustom2FilterCut / aod::femtouniverseparticle::pt)) && (aod::femtouniverseparticle::tpcNClsFound) > as(twotracksconfigs.confTrkTPCnclsMin) && (aod::femtouniverseparticle::tpcNClsCrossedRows) > as(twotracksconfigs.confTrkTPCcRowsMin) && (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt trackonefilter.confpLowPart1; Partition partsOneSigma = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt < trackonefilter.confPtHighPart1 && aod::femtouniverseparticle::pt > trackonefilter.confpLowPart1; + Partition partsOneSigmaITS = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt < trackonefilter.confPtHighPart1 && aod::femtouniverseparticle::pt > trackonefilter.confpLowPart1; Partition partsOneMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt < trackonefilter.confPtHighPart1 && aod::femtouniverseparticle::pt > trackonefilter.confpLowPart1; + Partition partsOneMCITS = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(trackonefilter.confChargePart1) && aod::femtouniverseparticle::pt < trackonefilter.confPtHighPart1 && aod::femtouniverseparticle::pt > trackonefilter.confpLowPart1; Partition partsOneMCTruth = aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kMCTruthTrack) && aod::femtouniverseparticle::pt < trackonefilter.confPtHighPart1 && aod::femtouniverseparticle::pt > trackonefilter.confpLowPart1; @@ -143,11 +160,14 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { } tracktwofilter; /// Partition for particle 2 - Partition partsTwo = (nabs(aod::track::dcaXY) < (additionalcuts.confDcaXYCustom1FilterCut + additionalcuts.confDcaXYCustom2FilterCut / aod::femtouniverseparticle::pt)) && (aod::femtouniverseparticle::tpcNClsFound > as(additionalcuts.confTrkTPCnclsMin)) && aod::femtouniverseparticle::tpcNClsCrossedRows > as(additionalcuts.confTrkTPCcRowsMin) && (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2) && aod::femtouniverseparticle::pt tracktwofilter.confpLowPart2; + Partition partsTwo = (nabs(aod::track::dcaXY) < (twotracksconfigs.confDcaXYCustom1FilterCut + twotracksconfigs.confDcaXYCustom2FilterCut / aod::femtouniverseparticle::pt)) && (aod::femtouniverseparticle::tpcNClsFound > as(twotracksconfigs.confTrkTPCnclsMin)) && aod::femtouniverseparticle::tpcNClsCrossedRows > as(twotracksconfigs.confTrkTPCcRowsMin) && (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2) && aod::femtouniverseparticle::pt tracktwofilter.confpLowPart2; + Partition partsTwoITS = (nabs(aod::track::dcaXY) < (twotracksconfigs.confDcaXYCustom1FilterCut + twotracksconfigs.confDcaXYCustom2FilterCut / aod::femtouniverseparticle::pt)) && (aod::femtouniverseparticle::tpcNClsFound > as(twotracksconfigs.confTrkTPCnclsMin)) && aod::femtouniverseparticle::tpcNClsCrossedRows > as(twotracksconfigs.confTrkTPCcRowsMin) && (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2) && aod::femtouniverseparticle::pt tracktwofilter.confpLowPart2; Partition partsTwoSigma = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && (aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2)) && aod::femtouniverseparticle::pt < tracktwofilter.confPtHighPart2 && aod::femtouniverseparticle::pt > tracktwofilter.confpLowPart2; + Partition partsTwoSigmaITS = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && (aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2)) && aod::femtouniverseparticle::pt < tracktwofilter.confPtHighPart2 && aod::femtouniverseparticle::pt > tracktwofilter.confpLowPart2; Partition partsTwoMC = aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack) && (aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2)) && aod::femtouniverseparticle::pt < tracktwofilter.confPtHighPart2 && aod::femtouniverseparticle::pt > tracktwofilter.confpLowPart2; + Partition partsTwoMCITS = aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack) && (aod::femtouniverseparticle::sign == as(tracktwofilter.confChargePart2)) && aod::femtouniverseparticle::pt < tracktwofilter.confPtHighPart2 && aod::femtouniverseparticle::pt > tracktwofilter.confpLowPart2; Partition partsTwoMCTruth = aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kMCTruthTrack) && aod::femtouniverseparticle::pt < tracktwofilter.confPtHighPart2 && aod::femtouniverseparticle::pt > tracktwofilter.confpLowPart2; @@ -170,33 +190,38 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { // Filter trackAdditionalfilter = (nabs(aod::femtouniverseparticle::eta) < twotracksconfigs.confEtaMax); // example filtering on configurable /// Particle part - ConfigurableAxis confTempFitVarBins{"confTempFitVarBins", {300, -0.15, 0.15}, "binning of the TempFitVar in the pT vs. TempFitVar plot"}; - ConfigurableAxis confTempFitVarpTBins{"confTempFitVarpTBins", {20, 0.5, 4.05}, "pT binning of the pT vs. TempFitVar plot"}; - ConfigurableAxis confTOFnSigmaBins{"confTOFnSigmaBins", {200, -4.975, 5.025}, "Binning of the TOF Nsigma vs. pT vs. centrality plot"}; + struct : o2::framework::ConfigurableGroup { + ConfigurableAxis confTempFitVarBins{"confTempFitVarBins", {300, -0.15, 0.15}, "binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confTempFitVarpTBins{"confTempFitVarpTBins", {20, 0.5, 4.05}, "pT binning of the pT vs. TempFitVar plot"}; + ConfigurableAxis confTOFnSigmaBins{"confTOFnSigmaBins", {200, -4.975, 5.025}, "Binning of the TOF Nsigma vs. pT vs. centrality plot"}; + } partBinning; /// Correlation part ConfigurableAxis confMultBins{"confMultBins", {VARIABLE_WIDTH, 0.0f, 4.0f, 8.0f, 12.0f, 16.0f, 20.0f, 24.0f, 28.0f, 32.0f, 36.0f, 40.0f, 44.0f, 48.0f, 52.0f, 56.0f, 60.0f, 64.0f, 68.0f, 72.0f, 76.0f, 80.0f, 84.0f, 88.0f, 92.0f, 96.0f, 100.0f, 200.0f, 99999.f}, "Mixing bins - multiplicity or centrality"}; // \todo to be obtained from the hash task ConfigurableAxis confMultKstarBins{"confMultKstarBins", {VARIABLE_WIDTH, 0.0f, 13.0f, 20.0f, 30.0f, 40.0f, 50.0f, 100.0f, 99999.f}, "Bins for kstar analysis in multiplicity or centrality bins (10 is maximum)"}; ConfigurableAxis confKtKstarBins{"confKtKstarBins", {VARIABLE_WIDTH, 0.0f, 0.2f, 0.4f, 0.6f, 0.8f, 1.0f, 2.0f, 99999.f}, "Bins for kstar analysis in kT bins (10 is maximum)"}; ConfigurableAxis confVtxBins{"confVtxBins", {VARIABLE_WIDTH, -10.0f, -8.f, -6.f, -4.f, -2.f, 0.f, 2.f, 4.f, 6.f, 8.f, 10.f}, "Mixing bins - z-vertex"}; - + ConfigurableAxis confkstarBins{"confkstarBins", {1500, 0., 6.}, "binning kstar"}; + ConfigurableAxis confkTBins{"confkTBins", {150, 0., 9.}, "binning kT"}; + ConfigurableAxis confmTBins{"confmTBins", {225, 0., 7.5}, "binning mT"}; + Configurable confNEventsMix{"confNEventsMix", 5, "Number of events for mixing"}; ConfigurableAxis confmTBins3D{"confmTBins3D", {VARIABLE_WIDTH, 1.02f, 1.14f, 1.20f, 1.26f, 1.38f, 1.56f, 1.86f, 4.50f}, "mT Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; ConfigurableAxis confmultBins3D{"confmultBins3D", {VARIABLE_WIDTH, 0.0f, 20.0f, 30.0f, 40.0f, 99999.0f}, "multiplicity Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; ColumnBinningPolicy colBinning{{confVtxBins, confMultBins}, true}; - ConfigurableAxis confkstarBins{"confkstarBins", {1500, 0., 6.}, "binning kstar"}; - ConfigurableAxis confkTBins{"confkTBins", {150, 0., 9.}, "binning kT"}; - ConfigurableAxis confmTBins{"confmTBins", {225, 0., 7.5}, "binning mT"}; - Configurable confNEventsMix{"confNEventsMix", 5, "Number of events for mixing"}; - Configurable confIsCPR{"confIsCPR", true, "Close Pair Rejection"}; - Configurable confCPRPlotPerRadii{"confCPRPlotPerRadii", false, "Plot CPR per radii"}; - Configurable confCPRdeltaPhiCutMax{"confCPRdeltaPhiCutMax", 0.0, "Delta Phi max cut for Close Pair Rejection"}; - Configurable confCPRdeltaPhiCutMin{"confCPRdeltaPhiCutMin", 0.0, "Delta Phi min cut for Close Pair Rejection"}; - Configurable confCPRdeltaEtaCutMax{"confCPRdeltaEtaCutMax", 0.0, "Delta Eta max cut for Close Pair Rejection"}; - Configurable confCPRdeltaEtaCutMin{"confCPRdeltaEtaCutMin", 0.0, "Delta Eta min cut for Close Pair Rejection"}; - Configurable confCPRChosenRadii{"confCPRChosenRadii", 0.80, "Chosen CPR radii"}; - Configurable confRandomizeCPR{"confRandomizeCPR", false, "Fill Dete-Dphi histograms in random order: p1, p2 or p2, p1"}; + struct : o2::framework::ConfigurableGroup { + Configurable confIsCPR{"confIsCPR", true, "Close Pair Rejection"}; + Configurable confCPRPlotPerRadii{"confCPRPlotPerRadii", false, "Plot CPR per radii"}; + Configurable confCPRdeltaPhiCutMax{"confCPRdeltaPhiCutMax", 0.0, "Delta Phi max cut for Close Pair Rejection"}; + Configurable confCPRdeltaPhiCutMin{"confCPRdeltaPhiCutMin", 0.0, "Delta Phi min cut for Close Pair Rejection"}; + Configurable confCPRdeltaEtaCutMax{"confCPRdeltaEtaCutMax", 0.0, "Delta Eta max cut for Close Pair Rejection"}; + Configurable confCPRdeltaEtaCutMin{"confCPRdeltaEtaCutMin", 0.0, "Delta Eta min cut for Close Pair Rejection"}; + Configurable confCPRChosenRadii{"confCPRChosenRadii", 0.80, "Chosen CPR radii"}; + Configurable confRandomizeCPR{"confRandomizeCPR", false, "Fill Dete-Dphi histograms in random order: p1, p2 or p2, p1"}; + ConfigurableAxis confDeltaEtaAxis{"confDeltaEtaAxis", {100, -0.15, 0.15}, "DeltaEta"}; + ConfigurableAxis confDeltaPhiStarAxis{"confDeltaPhiStarAxis", {100, -0.15, 0.15}, "DeltaPhiStar"}; + } cprConfigs; Configurable isPairIdentical{"isPairIdentical", true, "'true' for identical particles, 'false' for non-identical particles"}; struct : o2::framework::ConfigurableGroup { @@ -213,9 +238,6 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { Configurable confMaxEEMinv{"confMaxEEMinv", 0.002, "Max. minv of e-e+ pair for gamma pair rejection"}; Configurable confMaxDTheta{"confMaxDTheta", 0.008, "Max. DeltaTheta of pair for gamma pair rejection"}; - ConfigurableAxis confDeltaEtaAxis{"confDeltaEtaAxis", {100, -0.15, 0.15}, "DeltaEta"}; - ConfigurableAxis confDeltaPhiStarAxis{"confDeltaPhiStarAxis", {100, -0.15, 0.15}, "DeltaPhiStar"}; - FemtoUniverseFemtoContainer sameEventCont; FemtoUniverseFemtoContainer mixedEventCont; @@ -261,47 +283,64 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { HistogramRegistry sphericityRegistry{"SphericityHisto", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; /// TPC Pion/Kaon/Proton Sigma selection (general) - bool isNSigma(float mom, float nsigmaTPC, float nsigmaTOF) + bool isNSigma(float mom, float nsigmaTPC, float nsigmaTOF, float nsigmaITS = -999) { // |nsigma_TPC| < 3 for p < 0.5 GeV/c // |nsigma_combined| < 3 for p > 0.5 // using configurables: // confTOFpMin -- momentum value when we start using TOF; set to 1000 if TOF not needed + // confNsigmaTPCITSCombined -> TPC and ITS nSigma combined for momentum TPC Sigma for momentum < confTOFpMin // confNsigmaCombined -> TPC and TOF Sigma (combined) for momentum > confTOFpMin if (twotracksconfigs.isNsigmaRectangular) { if (mom < twotracksconfigs.confTOFpMin) { + if (twotracksconfigs.confIsAddITSNsigma) { + return (std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC && std::abs(nsigmaITS) < twotracksconfigs.confNsigmaTPCITSCombined); + } return std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC; } return (std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC && std::abs(nsigmaTOF) < twotracksconfigs.confNsigmaCombined); } if (mom < twotracksconfigs.confTOFpMin) { + if (twotracksconfigs.confIsAddITSNsigma) { + return std::hypot(nsigmaITS, nsigmaTPC) < twotracksconfigs.confNsigmaTPCITSCombined; + } return std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC; } return std::hypot(nsigmaTOF, nsigmaTPC) < twotracksconfigs.confNsigmaCombined; } - bool isNSigmaProton(float mom, float nsigmaTPC, float nsigmaTOF, float nsigmaTPCkaon, float nsigmaTOFkaon, float nsigmaTPCpion, float nsigmaTOFpion) + bool isNSigmaProton(float mom, float nsigmaTPC, float nsigmaTOF, float nsigmaTPCkaon, float nsigmaTOFkaon, float nsigmaTPCpion, float nsigmaTOFpion, float nsigmaITS = -999) { // with additional nSigma kaon and pion rejection - if (mom < twotracksconfigs.confTOFpMin) { - if (std::abs(nsigmaTPCkaon) < twotracksconfigs.confNsigmaKaonRejection || std::abs(nsigmaTPCpion) < twotracksconfigs.confNsigmaPionRejection) { + if (mom < twotracksconfigs.confTOFpMinRej) { // rejection based on TPC + if (std::abs(nsigmaTPCkaon) < twotracksconfigs.confTPCNsigmaKaonRejection || std::abs(nsigmaTPCpion) < twotracksconfigs.confTPCNsigmaPionRejection) { return false; } - } else if (std::hypot(nsigmaTOFkaon, nsigmaTPCkaon) < twotracksconfigs.confNsigmaKaonRejection || std::hypot(nsigmaTOFpion, nsigmaTPCpion) < twotracksconfigs.confNsigmaPionRejection) { - return false; + } else if (twotracksconfigs.isNsigmaRejectionCombined) { // rejection based on TPC and TOF combined + if (std::hypot(nsigmaTOFkaon, nsigmaTPCkaon) < twotracksconfigs.confNsigmaKaonRejection || std::hypot(nsigmaTOFpion, nsigmaTPCpion) < twotracksconfigs.confNsigmaPionRejection) { + return false; + } + } else if (std::abs(nsigmaTOFkaon) < twotracksconfigs.confTOFNsigmaKaonRejection || std::abs(nsigmaTOFpion) < twotracksconfigs.confTOFNsigmaPionRejection) { + return false; // rejection based on TOF } if (twotracksconfigs.isNsigmaRectangular) { if (mom < twotracksconfigs.confTOFpMin) { + if (twotracksconfigs.confIsAddITSNsigma) { + return (std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC && std::abs(nsigmaITS) < twotracksconfigs.confNsigmaTPCITSCombined); + } return std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC; } return (std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC && std::abs(nsigmaTOF) < twotracksconfigs.confNsigmaCombined); } if (mom < twotracksconfigs.confTOFpMin) { + if (twotracksconfigs.confIsAddITSNsigma) { + return std::hypot(nsigmaITS, nsigmaTPC) < twotracksconfigs.confNsigmaTPCITSCombined; + } return std::abs(nsigmaTPC) < twotracksconfigs.confNsigmaTPC; } return std::hypot(nsigmaTOF, nsigmaTPC) < twotracksconfigs.confNsigmaCombined; @@ -332,25 +371,53 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { return isNSigma(mom, nsigmaTPCK, nsigmaTOFK); } - bool isParticleNSigma(int8_t particle_number, float mom, float nsigmaTPCPr, float nsigmaTOFPr, float nsigmaTPCPi, float nsigmaTOFPi, float nsigmaTPCK, float nsigmaTOFK) + template + bool isParticleNSigma(int8_t particle_number, PartType const& part) { + float mom = part.p(); + float nsigmaTPCPr = trackCuts.getNsigmaTPC(part, o2::track::PID::Proton); + float nsigmaTOFPr = trackCuts.getNsigmaTOF(part, o2::track::PID::Proton); + float nsigmaITSPr = -999; + + float nsigmaTPCK = trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon); + float nsigmaTOFK = trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon); + float nsigmaITSK = -999; + + float nsigmaTPCPi = trackCuts.getNsigmaTPC(part, o2::track::PID::Pion); + float nsigmaTOFPi = trackCuts.getNsigmaTOF(part, o2::track::PID::Pion); + float nsigmaITSPi = -999; + + if constexpr (requires(PartType& part) { part.tpcFullNSigmaPr(); }) { + nsigmaTPCPr = part.tpcFullNSigmaPr(); + nsigmaTOFPr = part.tofFullNSigmaPr(); + nsigmaTPCK = part.tpcFullNSigmaKa(); + nsigmaTOFK = part.tofFullNSigmaKa(); + nsigmaTPCPi = part.tpcFullNSigmaPi(); + nsigmaTOFPi = part.tofFullNSigmaPi(); + } + + if constexpr (requires(PartType& part) { part.itsNSigmaPr(); }) { + nsigmaITSPr = trackCuts.getNsigmaITS(part, o2::track::PID::Proton); + nsigmaITSK = trackCuts.getNsigmaITS(part, o2::track::PID::Kaon); + nsigmaITSPi = trackCuts.getNsigmaITS(part, o2::track::PID::Pion); + } + enum ParticleNumber { partOne = 1, partTwo = 2 }; if (particle_number == partOne) { switch (trackonefilter.confPDGCodePartOne) { case kProton: case kProtonBar: - return isNSigmaProton(mom, nsigmaTPCPr, nsigmaTOFPr, nsigmaTPCK, nsigmaTOFK, nsigmaTPCPi, nsigmaTOFPi); + return isNSigmaProton(mom, nsigmaTPCPr, nsigmaTOFPr, nsigmaTPCK, nsigmaTOFK, nsigmaTPCPi, nsigmaTOFPi, nsigmaITSPr); case kPiPlus: case kPiMinus: case kPi0: - return isNSigma(mom, nsigmaTPCPi, nsigmaTOFPi); + return isNSigma(mom, nsigmaTPCPi, nsigmaTOFPi, nsigmaITSPi); case kKPlus: case kKMinus: - return isKaonNsigma(mom, nsigmaTPCK, nsigmaTOFK); case kK0Long: case kK0Short: - return isNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isNSigma(mom, nsigmaTPCK, nsigmaTOFK, nsigmaITSK); default: return false; } @@ -358,17 +425,16 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { switch (tracktwofilter.confPDGCodePartTwo) { case kProton: case kProtonBar: - return isNSigmaProton(mom, nsigmaTPCPr, nsigmaTOFPr, nsigmaTPCK, nsigmaTOFK, nsigmaTPCPi, nsigmaTOFPi); + return isNSigmaProton(mom, nsigmaTPCPr, nsigmaTOFPr, nsigmaTPCK, nsigmaTOFK, nsigmaTPCPi, nsigmaTOFPi, nsigmaITSPr); case kPiPlus: case kPiMinus: case kPi0: - return isNSigma(mom, nsigmaTPCPi, nsigmaTOFPi); + return isNSigma(mom, nsigmaTPCPi, nsigmaTOFPi, nsigmaITSPi); case kKPlus: case kKMinus: - return isKaonNsigma(mom, nsigmaTPCK, nsigmaTOFK); case kK0Long: case kK0Short: - return isNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isNSigma(mom, nsigmaTPCK, nsigmaTOFK, nsigmaITSK); default: return false; } @@ -395,7 +461,27 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { case kK0Short: return part.tofFullNSigmaKa(); default: - return 0; + return -999; + } + } + template + float whichTPCNSigma(int PDGcode, PartType const& part) + { + switch (PDGcode) { + case kProton: + case kProtonBar: + return part.tpcFullNSigmaPr(); + case kPiPlus: + case kPiMinus: + case kPi0: + return part.tpcFullNSigmaPi(); + case kKPlus: + case kKMinus: + case kK0Long: + case kK0Short: + return part.tpcFullNSigmaKa(); + default: + return -999; } } @@ -405,12 +491,17 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { eventHisto.init(&qaRegistry); - trackHistoPartOne.init(&qaRegistry, confTempFitVarpTBins, confTempFitVarBins, twotracksconfigs.confIsMC, trackonefilter.confPDGCodePartOne, true); + trackHistoPartOne.init(&qaRegistry, partBinning.confTempFitVarpTBins, partBinning.confTempFitVarBins, twotracksconfigs.confIsMC, trackonefilter.confPDGCodePartOne, true); - trackHistoPartTwo.init(&qaRegistry, confTempFitVarpTBins, confTempFitVarBins, twotracksconfigs.confIsMC, tracktwofilter.confPDGCodePartTwo, true); - qaRegistry.add("Tracks_one/nSigmaTOF", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}; Centrality", kTH3F, {confTempFitVarpTBins, confTOFnSigmaBins, confMultBins}); - qaRegistry.add("Tracks_two/nSigmaTOF", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}; Centrality", kTH3F, {confTempFitVarpTBins, confTOFnSigmaBins, confMultBins}); + trackHistoPartTwo.init(&qaRegistry, partBinning.confTempFitVarpTBins, partBinning.confTempFitVarBins, twotracksconfigs.confIsMC, tracktwofilter.confPDGCodePartTwo, true); + if (doprocessPurityQA || doprocessPurityQAWithITS) { + qaRegistry.add("Tracks_one/nSigmaTOF", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}; Centrality", kTH3F, {partBinning.confTempFitVarpTBins, partBinning.confTOFnSigmaBins, confMultBins}); + qaRegistry.add("Tracks_two/nSigmaTOF", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}; Centrality", kTH3F, {partBinning.confTempFitVarpTBins, partBinning.confTOFnSigmaBins, confMultBins}); + + qaRegistry.add("Tracks_one/nSigmaTPC", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}; Centrality", kTH3F, {partBinning.confTempFitVarpTBins, partBinning.confTOFnSigmaBins, confMultBins}); + qaRegistry.add("Tracks_two/nSigmaTPC", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}; Centrality", kTH3F, {partBinning.confTempFitVarpTBins, partBinning.confTOFnSigmaBins, confMultBins}); + } if (confFillDebug) { sphericityRegistry.add("sphericity", ";Sphericity;Entries", kTH1F, {{150, 0.0, 3, "Sphericity"}}); mixQaRegistry.add("MixingQA/hSECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); @@ -482,8 +573,8 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { } pairCleaner.init(&qaRegistry); - if (confIsCPR.value) { - pairCloseRejection.init(&resultRegistry, &qaRegistry, confDeltaEtaAxis, confDeltaPhiStarAxis, confCPRdeltaPhiCutMin.value, confCPRdeltaPhiCutMax.value, confCPRdeltaEtaCutMin.value, confCPRdeltaEtaCutMax.value, confCPRChosenRadii.value, confCPRPlotPerRadii.value); + if (cprConfigs.confIsCPR.value) { + pairCloseRejection.init(&resultRegistry, &qaRegistry, cprConfigs.confDeltaEtaAxis, cprConfigs.confDeltaPhiStarAxis, cprConfigs.confCPRdeltaPhiCutMin.value, cprConfigs.confCPRdeltaPhiCutMax.value, cprConfigs.confCPRdeltaEtaCutMin.value, cprConfigs.confCPRdeltaEtaCutMax.value, cprConfigs.confCPRChosenRadii.value, cprConfigs.confCPRPlotPerRadii.value); } } @@ -517,7 +608,7 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { /// Histogramming same event if ((pairType == type1 || pairType == type2) && fillQA) { for (const auto& part : groupPartsOne) { - if (!isParticleNSigma((int8_t)1, part.p(), trackCuts.getNsigmaTPC(part, o2::track::PID::Proton), trackCuts.getNsigmaTOF(part, o2::track::PID::Proton), trackCuts.getNsigmaTPC(part, o2::track::PID::Pion), trackCuts.getNsigmaTOF(part, o2::track::PID::Pion), trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)1, part)) { continue; } trackHistoPartOne.fillQA(part); @@ -526,7 +617,7 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { if ((pairType == type1 || pairType == type3) && fillQA) { for (const auto& part : groupPartsTwo) { - if (!isParticleNSigma((int8_t)2, part.p(), trackCuts.getNsigmaTPC(part, o2::track::PID::Proton), trackCuts.getNsigmaTOF(part, o2::track::PID::Proton), trackCuts.getNsigmaTPC(part, o2::track::PID::Pion), trackCuts.getNsigmaTOF(part, o2::track::PID::Pion), trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)2, part)) { continue; } trackHistoPartTwo.fillQA(part); @@ -538,11 +629,11 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { /// Now build the combinations for non-identical particle pairs for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { - if (!isParticleNSigma((int8_t)1, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)1, p1)) { continue; } - if (!isParticleNSigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)2, p2)) { continue; } @@ -555,12 +646,12 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { continue; } - if (confIsCPR.value) { + if (cprConfigs.confIsCPR.value) { double rand = 0; auto part1 = p1; auto part2 = p2; - if (confRandomizeCPR) { + if (cprConfigs.confRandomizeCPR) { randgen = new TRandom2(0); rand = randgen->Rndm(); @@ -587,20 +678,20 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { /// Now build the combinations for identical particles pairs for (const auto& [p1, p2] : combinations(CombinationsStrictlyUpperIndexPolicy(groupPartsOne, groupPartsOne))) { - if (!isParticleNSigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)2, p1)) { continue; } - if (!isParticleNSigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)2, p2)) { continue; } - if (confIsCPR.value) { + if (cprConfigs.confIsCPR.value) { double rand = 0; auto part1 = p1; auto part2 = p2; - if (confRandomizeCPR) { + if (cprConfigs.confRandomizeCPR) { randgen = new TRandom2(0); rand = randgen->Rndm(); @@ -687,27 +778,68 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { float const multCol = col.multV0M(); for (const auto& part : thegroupPartsOne) { - if (!isParticleNSigma((int8_t)1, part.p(), part.tpcFullNSigmaPr(), part.tofFullNSigmaPr(), part.tpcFullNSigmaPi(), part.tofFullNSigmaPi(), part.tpcFullNSigmaKa(), part.tofFullNSigmaKa())) { + if (!isParticleNSigma((int8_t)1, part)) { continue; } trackHistoPartOne.fillQA(part); qaRegistry.fill(HIST("Tracks_one/nSigmaTOF"), part.pt(), whichTOFNSigma(trackonefilter.confPDGCodePartOne, part), multCol); + qaRegistry.fill(HIST("Tracks_one/nSigmaTPC"), part.pt(), whichTPCNSigma(trackonefilter.confPDGCodePartOne, part), multCol); } if (!isPairIdentical) { // fill second for (const auto& part : thegroupPartsTwo) { - if (!isParticleNSigma((int8_t)2, part.p(), part.tpcFullNSigmaPr(), part.tofFullNSigmaPr(), part.tpcFullNSigmaPi(), part.tofFullNSigmaPi(), part.tpcFullNSigmaKa(), part.tofFullNSigmaKa())) { + if (!isParticleNSigma((int8_t)2, part)) { continue; } trackHistoPartTwo.fillQA(part); qaRegistry.fill(HIST("Tracks_two/nSigmaTOF"), part.pt(), whichTOFNSigma(tracktwofilter.confPDGCodePartTwo, part), multCol); + qaRegistry.fill(HIST("Tracks_two/nSigmaTPC"), part.pt(), whichTPCNSigma(trackonefilter.confPDGCodePartOne, part), multCol); } } } PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processPurityQA, "Enable processing QA for purity estimation", false); + /// \param col subscribe to the collision table (Data) + /// \param parts subscribe to the femtoUniverseParticleTable + void processPurityQAWithITS(soa::Filtered::iterator const& col, + FilteredFemtoFullSigmaParticlesITS const&) + { + fillCollision(col); + if (confFillDebug) { + sphericityRegistry.fill(HIST("sphericity"), col.sphericity()); + } + + auto thegroupPartsOne = partsOneSigmaITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + auto thegroupPartsTwo = partsTwoSigmaITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + + float const multCol = col.multV0M(); + + for (const auto& part : thegroupPartsOne) { + if (!isParticleNSigma((int8_t)1, part)) { + continue; + } + + trackHistoPartOne.fillQA(part); + qaRegistry.fill(HIST("Tracks_one/nSigmaTOF"), part.pt(), whichTOFNSigma(trackonefilter.confPDGCodePartOne, part), multCol); + qaRegistry.fill(HIST("Tracks_one/nSigmaTPC"), part.pt(), whichTPCNSigma(trackonefilter.confPDGCodePartOne, part), multCol); + } + + if (!isPairIdentical) { + // fill second + for (const auto& part : thegroupPartsTwo) { + if (!isParticleNSigma((int8_t)2, part)) { + continue; + } + trackHistoPartTwo.fillQA(part); + qaRegistry.fill(HIST("Tracks_two/nSigmaTOF"), part.pt(), whichTOFNSigma(tracktwofilter.confPDGCodePartTwo, part), multCol); + qaRegistry.fill(HIST("Tracks_two/nSigmaTPC"), part.pt(), whichTPCNSigma(trackonefilter.confPDGCodePartOne, part), multCol); + } + } + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processPurityQAWithITS, "Enable processing QA for purity estimation", false); + /// process function for to call doSameEvent with Data /// \param col subscribe to the collision table (Data) /// \param parts subscribe to the femtoUniverseParticleTable @@ -737,6 +869,35 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { } PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processSameEvent, "Enable processing same event", true); + /// process function for to call doSameEvent with Data + /// \param col subscribe to the collision table (Data) + /// \param parts subscribe to the femtoUniverseParticleTable + void processSameEventWithITS(soa::Filtered::iterator const& col, + FilteredFemtoFullParticlesITS const& parts) + { + fillCollision(col); + if (confFillDebug) { + sphericityRegistry.fill(HIST("sphericity"), col.sphericity()); + } + + auto thegroupPartsOne = partsOneITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + auto thegroupPartsTwo = partsTwoITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + + bool fillQA = true; + + if (processPair.cfgProcessPM) { + doSameEvent(thegroupPartsOne, thegroupPartsTwo, parts, col.magField(), col.multV0M(), 1, fillQA); + fillQA = false; + } + if (processPair.cfgProcessPP) { + doSameEvent(thegroupPartsOne, thegroupPartsOne, parts, col.magField(), col.multV0M(), 2, fillQA); + } + if (processPair.cfgProcessMM) { + doSameEvent(thegroupPartsTwo, thegroupPartsTwo, parts, col.magField(), col.multV0M(), 3, fillQA); + } + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processSameEventWithITS, "Enable processing same event", false); + /// process function for to call doSameEvent with Monte Carlo /// \param col subscribe to the collision table (Monte Carlo Reconstructed reconstructed) /// \param parts subscribe to joined table FemtoUniverseParticles and FemtoUniverseMCLabels to access Monte Carlo truth @@ -780,11 +941,11 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { { for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { - if (!isParticleNSigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)2, p1)) { continue; } - if (!isParticleNSigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon))) { + if (!isParticleNSigma((int8_t)2, p2)) { continue; } @@ -792,12 +953,12 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { continue; } - if (confIsCPR.value) { + if (cprConfigs.confIsCPR.value) { double rand = 0; auto part1 = p1; auto part2 = p2; - if (confRandomizeCPR) { + if (cprConfigs.confRandomizeCPR) { randgen = new TRandom2(0); rand = randgen->Rndm(); @@ -912,6 +1073,45 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { } PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processMixedEvent, "Enable processing mixed events", true); + /// process function for to call doMixedEvent with Data + /// \param cols subscribe to the collisions table (Data) + /// \param parts subscribe to the femtoUniverseParticleTable + void processMixedEventWithITS(soa::Filtered const& cols, + FilteredFemtoFullParticlesITS const& parts) + { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinning, confNEventsMix, -1, cols, cols)) { + + const float multiplicityCol = collision1.multV0M(); + if (confFillDebug) { + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinning.getBin({collision1.posZ(), multiplicityCol})); + } + + const auto& magFieldTesla1 = collision1.magField(); + const auto& magFieldTesla2 = collision2.magField(); + + if (magFieldTesla1 != magFieldTesla2) { + continue; + } + + if (processPair.cfgProcessPM) { + auto groupPartsOne = partsOneITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwoITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + doMixedEvent(groupPartsOne, groupPartsTwo, parts, magFieldTesla1, multiplicityCol, 1); + } + if (processPair.cfgProcessPP) { + auto groupPartsOne = partsOneITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsOneITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + doMixedEvent(groupPartsOne, groupPartsTwo, parts, magFieldTesla1, multiplicityCol, 2); + } + if (processPair.cfgProcessMM) { + auto groupPartsOne = partsTwoITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwoITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + doMixedEvent(groupPartsOne, groupPartsTwo, parts, magFieldTesla1, multiplicityCol, 3); + } + } + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processMixedEventWithITS, "Enable processing mixed events", false); + /// brief process function for to call doMixedEvent with Monte Carlo /// \param cols subscribe to the collisions table (Monte Carlo Reconstructed reconstructed) /// \param parts subscribe to joined table FemtoUniverseParticles and FemtoUniverseMCLables to access Monte Carlo truth @@ -954,6 +1154,76 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { } PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processMixedEventMC, "Enable processing mixed events MC reco", false); + template + void doFractions(PartType const& p1, PartType const& p2, PartsType const& parts, PartsMcType const& mparts, float magFieldTesla, int partType) + { + if (!isParticleNSigma((int8_t)2, p1) || !isParticleNSigma((int8_t)2, p2)) { + return; + } + + if (cprConfigs.confIsCPR.value) { + double rand = 0; + auto part1 = p1; + auto part2 = p2; + + if (cprConfigs.confRandomizeCPR) { + randgen = new TRandom2(0); + rand = randgen->Rndm(); + + if (rand > randCompareValue) { + part1 = p2; + part2 = p1; + } + delete randgen; + } + if (pairCloseRejection.isClosePair(part1, part2, parts, magFieldTesla, femto_universe_container::EventType::mixed)) { + return; + } + } + switch (partType) { + case 1: + mixedMultRegistryPM.fill(HIST("MCreco/motherParticle"), p1.motherPDG(), p2.motherPDG()); + break; + case 2: + mixedMultRegistryPP.fill(HIST("MCreco/motherParticle"), p1.motherPDG(), p2.motherPDG()); + break; + case 3: + mixedMultRegistryMM.fill(HIST("MCreco/motherParticle"), p1.motherPDG(), p2.motherPDG()); + break; + default: + break; + } + auto mcPartId1 = p1.fdMCParticleId(); + auto mcPartId2 = p2.fdMCParticleId(); + if (mcPartId1 == -1 || mcPartId2 == -1) { + return; + } + const auto& mcParticle1 = mparts.iteratorAt(mcPartId1); + const auto& mcParticle2 = mparts.iteratorAt(mcPartId2); + switch (partType) { + case 1: + if ((trackonefilter.confPDGCodePartOne != mcParticle1.pdgMCTruth()) || (tracktwofilter.confPDGCodePartTwo != mcParticle2.pdgMCTruth())) { + return; + } + mixedMultRegistryPM.fill(HIST("MCreco/motherParticlePDGCheck"), p1.motherPDG(), p2.motherPDG()); + break; + case 2: + if ((trackonefilter.confPDGCodePartOne != mcParticle1.pdgMCTruth()) || (trackonefilter.confPDGCodePartOne != mcParticle2.pdgMCTruth())) { + return; + } + mixedMultRegistryPP.fill(HIST("MCreco/motherParticlePDGCheck"), p1.motherPDG(), p2.motherPDG()); + break; + case 3: + if ((tracktwofilter.confPDGCodePartTwo != mcParticle1.pdgMCTruth()) || (tracktwofilter.confPDGCodePartTwo != mcParticle2.pdgMCTruth())) { + return; + } + mixedMultRegistryMM.fill(HIST("MCreco/motherParticlePDGCheck"), p1.motherPDG(), p2.motherPDG()); + break; + default: + break; + } + } + /// brief process function for to call process fractions with Monte Carlo reco /// \param cols subscribe to the collisions table (Monte Carlo Reconstructed reconstructed) /// \param parts subscribe to joined table FemtoUniverseParticles and FemtoUniverseMCLables to access Monte Carlo truth @@ -962,81 +1232,48 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { FemtoRecoParticles const& parts, o2::aod::FdMCParticles const& mparts) { - auto doFractions = [&](auto& p1, auto& p2, auto& magFieldTesla, int partType) -> void { - if (!isParticleNSigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon))) { - return; + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinning, confNEventsMix, -1, cols, cols)) { + const float multiplicityCol = collision1.multV0M(); + if (confFillDebug) { + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinning.getBin({collision1.posZ(), multiplicityCol})); } + const auto& magFieldTesla1 = collision1.magField(); + const auto& magFieldTesla2 = collision2.magField(); - if (!isParticleNSigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon))) { - return; + if (magFieldTesla1 != magFieldTesla2) { + continue; } - if (confIsCPR.value) { - double rand = 0; - auto part1 = p1; - auto part2 = p2; - - if (confRandomizeCPR) { - randgen = new TRandom2(0); - rand = randgen->Rndm(); + auto groupPartsOne = partsOneMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); - if (rand > randCompareValue) { - part1 = p2; - part2 = p1; - } - delete randgen; - } - if (pairCloseRejection.isClosePair(part1, part2, parts, magFieldTesla, femto_universe_container::EventType::mixed)) { - return; + if (processPair.cfgProcessPM) { + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + doFractions(p1, p2, parts, mparts, magFieldTesla1, 1); } } - switch (partType) { - case 1: - mixedMultRegistryPM.fill(HIST("MCreco/motherParticle"), p1.motherPDG(), p2.motherPDG()); - break; - case 2: - mixedMultRegistryPP.fill(HIST("MCreco/motherParticle"), p1.motherPDG(), p2.motherPDG()); - break; - case 3: - mixedMultRegistryMM.fill(HIST("MCreco/motherParticle"), p1.motherPDG(), p2.motherPDG()); - break; - default: - break; - } - auto mcPartId1 = p1.fdMCParticleId(); - if (mcPartId1 == -1) { - return; - } - auto mcPartId2 = p2.fdMCParticleId(); - if (mcPartId2 == -1) { - return; + if (processPair.cfgProcessPP) { + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsOne))) { + doFractions(p1, p2, parts, mparts, magFieldTesla1, 2); + } } - const auto& mcParticle1 = mparts.iteratorAt(mcPartId1); - const auto& mcParticle2 = mparts.iteratorAt(mcPartId2); - switch (partType) { - case 1: - if ((trackonefilter.confPDGCodePartOne != mcParticle1.pdgMCTruth()) || (tracktwofilter.confPDGCodePartTwo != mcParticle2.pdgMCTruth())) { - return; - } - mixedMultRegistryPM.fill(HIST("MCreco/motherParticlePDGCheck"), p1.motherPDG(), p2.motherPDG()); - break; - case 2: - if ((trackonefilter.confPDGCodePartOne != mcParticle1.pdgMCTruth()) || (trackonefilter.confPDGCodePartOne != mcParticle2.pdgMCTruth())) { - return; - } - mixedMultRegistryPP.fill(HIST("MCreco/motherParticlePDGCheck"), p1.motherPDG(), p2.motherPDG()); - break; - case 3: - if ((tracktwofilter.confPDGCodePartTwo != mcParticle1.pdgMCTruth()) || (tracktwofilter.confPDGCodePartTwo != mcParticle2.pdgMCTruth())) { - return; - } - mixedMultRegistryMM.fill(HIST("MCreco/motherParticlePDGCheck"), p1.motherPDG(), p2.motherPDG()); - break; - default: - break; + if (processPair.cfgProcessMM) { + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsTwo, groupPartsTwo))) { + doFractions(p1, p2, parts, mparts, magFieldTesla1, 3); + } } - }; + } + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processFractions, "Enable processing fractions for MC reco", false); + /// brief process function for to call process fractions with Monte Carlo reco + /// \param cols subscribe to the collisions table (Monte Carlo Reconstructed reconstructed) + /// \param parts subscribe to joined table FemtoUniverseParticles and FemtoUniverseMCLables to access Monte Carlo truth + /// \param FemtoUniverseMCParticles subscribe to the Monte Carlo truth table + void processFractionsWithITS(o2::aod::FdCollisions const& cols, + FemtoRecoParticlesITS const& parts, + o2::aod::FdMCParticles const& mparts) + { for (const auto& [collision1, collision2] : soa::selfCombinations(colBinning, confNEventsMix, -1, cols, cols)) { const float multiplicityCol = collision1.multV0M(); if (confFillDebug) { @@ -1049,27 +1286,27 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended { continue; } - auto groupPartsOne = partsOneMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); - auto groupPartsTwo = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + auto groupPartsOne = partsOneMCITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwoMCITS->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); if (processPair.cfgProcessPM) { for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { - doFractions(p1, p2, magFieldTesla1, 1); + doFractions(p1, p2, parts, mparts, magFieldTesla1, 1); } } if (processPair.cfgProcessPP) { for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsOne))) { - doFractions(p1, p2, magFieldTesla1, 2); + doFractions(p1, p2, parts, mparts, magFieldTesla1, 2); } } if (processPair.cfgProcessMM) { for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsTwo, groupPartsTwo))) { - doFractions(p1, p2, magFieldTesla1, 3); + doFractions(p1, p2, parts, mparts, magFieldTesla1, 3); } } } } - PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processFractions, "Enable processing fractions for MC reco", false); + PROCESS_SWITCH(FemtoUniversePairTaskTrackTrackMultKtExtended, processFractionsWithITS, "Enable processing fractions for MC reco", false); /// This function processes fills histograms for fractions analysis /// \todo the trivial loops over the collisions and tracks should be factored out since they will be common to all combinations of T-T, T-V0, V0-V0, ...