diff --git a/PWGCF/FemtoDream/TableProducer/femtoDreamProducerReducedTask.cxx b/PWGCF/FemtoDream/TableProducer/femtoDreamProducerReducedTask.cxx index dc5fcc8a862..819ca1b4ca6 100644 --- a/PWGCF/FemtoDream/TableProducer/femtoDreamProducerReducedTask.cxx +++ b/PWGCF/FemtoDream/TableProducer/femtoDreamProducerReducedTask.cxx @@ -55,6 +55,7 @@ namespace o2::aod using FemtoFullCollision = soa::Join::iterator; using FemtoFullCollisionMC = soa::Join::iterator; using FemtoFullCollision_noCent_MC = soa::Join::iterator; +using FemtoFullCollision_CentPbPb = soa::Join::iterator; using FemtoFullTracks = soa::Join outputCollision; + Produces outputExtQnCollision; Produces outputParts; Produces outputPartsMC; Produces outputDebugParts; @@ -104,14 +106,39 @@ struct femtoDreamProducerReducedTask { Configurable> ConfTrkITSnclsIbMin{FemtoDreamTrackSelection::getSelectionName(femtoDreamTrackSelection::kITSnClsIbMin, "ConfTrk"), std::vector{-1.f, 1.f}, FemtoDreamTrackSelection::getSelectionHelper(femtoDreamTrackSelection::kITSnClsIbMin, "Track selection: ")}; Configurable> ConfTrkDCAxyMax{FemtoDreamTrackSelection::getSelectionName(femtoDreamTrackSelection::kDCAxyMax, "ConfTrk"), std::vector{0.1f, 0.5f}, FemtoDreamTrackSelection::getSelectionHelper(femtoDreamTrackSelection::kDCAxyMax, "Track selection: ")}; /// here we need an open cut to do the DCA fits later on! Configurable> ConfTrkDCAzMax{FemtoDreamTrackSelection::getSelectionName(femtoDreamTrackSelection::kDCAzMax, "ConfTrk"), std::vector{0.2f, 0.5f}, FemtoDreamTrackSelection::getSelectionHelper(femtoDreamTrackSelection::kDCAzMax, "Track selection: ")}; - Configurable> ConfTrkPIDnSigmaMax{FemtoDreamTrackSelection::getSelectionName(femtoDreamTrackSelection::kPIDnSigmaMax, "Conf"), std::vector{3.5f, 3.f, 2.5f}, FemtoDreamTrackSelection::getSelectionHelper(femtoDreamTrackSelection::kPIDnSigmaMax, "Track selection: ")}; + Configurable> ConfTrkPIDnSigmaMax{FemtoDreamTrackSelection::getSelectionName(femtoDreamTrackSelection::kPIDnSigmaMax, "ConfTrk"), std::vector{3.5f, 3.f, 2.5f}, FemtoDreamTrackSelection::getSelectionHelper(femtoDreamTrackSelection::kPIDnSigmaMax, "Track selection: ")}; // off set the center of the nsigma distribution to deal with bad TPC/TOF calibration Configurable ConfTrkPIDnSigmaOffsetTPC{"ConfTrkPIDnSigmaOffsetTPC", 0., "Offset for TPC nSigma because of bad calibration"}; Configurable ConfTrkPIDnSigmaOffsetTOF{"ConfTrkPIDnSigmaOffsetTOF", 0., "Offset for TOF nSigma because of bad calibration"}; Configurable> ConfTrkPIDspecies{"ConfTrkPIDspecies", std::vector{o2::track::PID::Pion, o2::track::PID::Kaon, o2::track::PID::Proton, o2::track::PID::Deuteron}, "Trk sel: Particles species for PID"}; + struct : o2::framework::ConfigurableGroup { + Configurable ConfgFlowCalculate{"ConfgFlowCalculate", false, "Evt sel: Cumulant of flow"}; // To do + Configurable ConfgQnSeparation{"ConfgQnSeparation", false, "Evt sel: Qn of event"}; + Configurable> ConfQnBinSeparator{"ConfQnBinSeparator", std::vector{-999.f, -999.f, -999.f}, "Qn bin separator"}; + Configurable ConfCentralityMax{"ConfCentralityMax", 80.f, "Evt sel: Maximum Centrality cut"}; + Configurable ConfCentBinWidth{"ConfCentBinWidth", 1.f, "Centrality bin length for qn separator"}; + Configurable ConfQnBinMin{"ConfQnBinMin", 0, "Minimum qn bin"}; + Configurable ConfQnBinMax{"ConfQnBinMax", 10, "Maximum qn bin"}; + } qnCal; + + struct : o2::framework::ConfigurableGroup { + Configurable ConfIsPbPb{"ConfIsPbPb", false, "Running on Run3 or Run2"}; // Choose if running on PbPb data + Configurable ConfIsUsePileUp{"ConfIsUsePileUp", false, "Required for choosing whether to run the pile-up cuts"}; + Configurable ConfEvNoSameBunchPileup{"ConfEvNoSameBunchPileup", false, "Require kNoSameBunchPileup selection on Events."}; + Configurable ConfEvIsGoodZvtxFT0vsPV{"ConfEvIsGoodZvtxFT0vsPV", false, "Require kIsGoodZvtxFT0vsPV selection on Events."}; + Configurable ConfEvIsGoodITSLayersAll{"ConfEvIsGoodITSLayersAll", false, "Require kIsGoodITSLayersAll selection on Events."}; + Configurable ConfEvNoCollInRofStandard{"ConfEvNoCollInRofStandard", false, "Require kNoCollInRofStandard selection on Events."}; + Configurable ConfEvNoHighMultCollInPrevRof{"ConfEvNoHighMultCollInPrevRof", false, "Require kNoHighMultCollInPrevRof selection on Events."}; + Configurable ConfEvNoCollInTimeRangeStandard{"ConfEvNoCollInTimeRangeStandard", false, "Require kNoCollInTimeRangeStandard selection on Events."}; + Configurable ConfEvIsVertexITSTPC{"ConfEvIsVertexITSTPC", false, "Require kIsVertexITSTPC selection on Events"}; + Configurable ConfTPCOccupancyMin{"ConfTPCOccupancyMin", 0, "Minimum value for TPC Occupancy selection"}; + Configurable ConfTPCOccupancyMax{"ConfTPCOccupancyMax", 1000, "Maximum value for TPC Occupancy selection"}; + } evtSel_PbPb; + HistogramRegistry qaRegistry{"QAHistos", {}, OutputObjHandlingPolicy::AnalysisObject}; HistogramRegistry Registry{"Tracks", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry FlowRegistry{"QnandFlowInfo", {}, OutputObjHandlingPolicy::AnalysisObject}; int mRunNumber; float mMagField; @@ -144,6 +171,11 @@ struct femtoDreamProducerReducedTask { trackCuts.init(&qaRegistry, &Registry); + + if (qnCal.ConfgFlowCalculate) { + colCuts.initFlow(&FlowRegistry, qnCal.ConfgQnSeparation); + } + mRunNumber = 0; mMagField = 0.0; /// Initializing CCDB @@ -330,8 +362,166 @@ struct femtoDreamProducerReducedTask { } } - void - processData(aod::FemtoFullCollision const& col, aod::BCsWithTimestamps const&, aod::FemtoFullTracks const& tracks) + // Centrality (Multiplicity percentile) obtained from FT0C + // Pile-up rejection involved + template + void fillCollisionsAndTracks_PbPb(CollisionType const& col, TrackType const& tracks) + { + const auto vtxZ = col.posZ(); + const auto spher = colCuts.computeSphericity(col, tracks); + int mult = 0; + int multNtr = 0; + if (ConfIsRun3) { + if constexpr (useCentrality) { + mult = col.centFT0C(); + } else { + mult = 0.; + } + multNtr = col.multNTracksPV(); + } else { + mult = 1; // multiplicity percentile is known in Run 2 + multNtr = col.multTracklets(); + } + if (ConfEvtUseTPCmult) { + multNtr = col.multTPC(); + } + colCuts.fillQA(col, mult); + + /// First thing to do is to check whether the basic event selection criteria are fulfilled + /// That includes checking if there are any usable tracks in a collision + if (!colCuts.isSelectedCollision(col)) { + return; + } + if (colCuts.isEmptyCollision(col, tracks, trackCuts)) { + return; + } + + // Pileup rejection in PbPb data + if (evtSel_PbPb.ConfIsPbPb && evtSel_PbPb.ConfIsUsePileUp && + !colCuts.isPileUpCollisionPbPb(col, evtSel_PbPb.ConfEvNoSameBunchPileup, evtSel_PbPb.ConfEvIsGoodZvtxFT0vsPV, + evtSel_PbPb.ConfEvIsGoodITSLayersAll, evtSel_PbPb.ConfEvNoCollInRofStandard, + evtSel_PbPb.ConfEvNoHighMultCollInPrevRof, evtSel_PbPb.ConfEvNoCollInTimeRangeStandard, + evtSel_PbPb.ConfEvIsVertexITSTPC, + evtSel_PbPb.ConfTPCOccupancyMin, evtSel_PbPb.ConfTPCOccupancyMax)) { + return; + } + // now the table is filled + outputCollision(vtxZ, mult, multNtr, spher, mMagField); + + // these IDs are necessary to keep track of the children + // since this producer only produces the tables for tracks, there are no children + std::vector childIDs = {0, 0}; + for (auto& track : tracks) { + /// if the most open selection criteria are not fulfilled there is no point looking further at the track + if (!trackCuts.isSelectedMinimal(track)) { + continue; + } + trackCuts.fillQA(track); + // an array of two bit-wise containers of the systematic variations is obtained + // one container for the track quality cuts and one for the PID cuts + auto cutContainer = trackCuts.getCutContainer(track, track.pt(), track.eta(), sqrtf(powf(track.dcaXY(), 2.f) + powf(track.dcaZ(), 2.f))); + + // now the table is filled + outputParts(outputCollision.lastIndex(), + track.pt(), + track.eta(), + track.phi(), + aod::femtodreamparticle::ParticleType::kTrack, + cutContainer.at(femtoDreamTrackSelection::TrackContainerPosition::kCuts), + cutContainer.at(femtoDreamTrackSelection::TrackContainerPosition::kPID), + track.dcaXY(), childIDs, 0, 0); + if constexpr (isMC) { + fillMCParticle(col, track, o2::aod::femtodreamparticle::ParticleType::kTrack); + } + + if (ConfIsDebug) { + outputDebugParts(track.sign(), + (uint8_t)track.tpcNClsFound(), + track.tpcNClsFindable(), + (uint8_t)track.tpcNClsCrossedRows(), + track.tpcNClsShared(), + track.tpcInnerParam(), + track.itsNCls(), + track.itsNClsInnerBarrel(), + track.dcaXY(), + track.dcaZ(), + track.tpcSignal(), + track.tpcNSigmaEl(), + track.tpcNSigmaPi(), + track.tpcNSigmaKa(), + track.tpcNSigmaPr(), + track.tpcNSigmaDe(), + track.tpcNSigmaTr(), + track.tpcNSigmaHe(), + track.tofNSigmaEl(), + track.tofNSigmaPi(), + track.tofNSigmaKa(), + track.tofNSigmaPr(), + track.tofNSigmaDe(), + track.tofNSigmaTr(), + track.tofNSigmaHe(), + -1, + track.itsNSigmaEl(), + track.itsNSigmaPi(), + track.itsNSigmaKa(), + track.itsNSigmaPr(), + track.itsNSigmaDe(), + track.itsNSigmaTr(), + track.itsNSigmaHe(), + -999., -999., -999., -999., -999., -999., + -999., -999., -999., -999., -999., -999., -999.); + } + } + } + + // Calculate and separate qn bins + // Do and fill cumulant in qn bins + template + void fillCollisionsFlow(CollisionType const& col, TrackType const& tracks) + { + // get magnetic field for run + + const auto spher = colCuts.computeSphericity(col, tracks); + int multNtr = 0; + if (ConfIsRun3) { + multNtr = col.multNTracksPV(); + } else { + multNtr = col.multTracklets(); + } + if (ConfEvtUseTPCmult) { + multNtr = col.multTPC(); + } + + /// First thing to do is to check whether the basic event selection criteria are fulfilled + /// That includes checking if there are any usable tracks in a collision + if (!colCuts.isSelectedCollision(col)) { + return; + } + if (colCuts.isEmptyCollision(col, tracks, trackCuts)) { + return; + } + + // Pileup rejection in PbPb data + if (evtSel_PbPb.ConfIsPbPb && evtSel_PbPb.ConfIsUsePileUp && + !colCuts.isPileUpCollisionPbPb(col, evtSel_PbPb.ConfEvNoSameBunchPileup, evtSel_PbPb.ConfEvIsGoodZvtxFT0vsPV, + evtSel_PbPb.ConfEvIsGoodITSLayersAll, evtSel_PbPb.ConfEvNoCollInRofStandard, + evtSel_PbPb.ConfEvNoHighMultCollInPrevRof, evtSel_PbPb.ConfEvNoCollInTimeRangeStandard, + evtSel_PbPb.ConfEvIsVertexITSTPC, + evtSel_PbPb.ConfTPCOccupancyMin, evtSel_PbPb.ConfTPCOccupancyMax)) { + return; + } + + // Calculate and fill qnBins + auto qnBin = colCuts.myqnBin(col, qnCal.ConfCentralityMax, qnCal.ConfQnBinSeparator, spher, multNtr, 1.f); + if (qnBin < qnCal.ConfQnBinMin || qnBin > qnCal.ConfQnBinMax) { + qnBin = -999; + } + colCuts.fillCumulants(col, tracks); + colCuts.doCumulants(col, qnCal.ConfgQnSeparation, qnBin); + outputExtQnCollision(qnBin); + } + + void processData(aod::FemtoFullCollision const& col, aod::BCsWithTimestamps const&, aod::FemtoFullTracks const& tracks) { // get magnetic field for run getMagneticFieldTesla(col.bc_as()); @@ -370,8 +560,24 @@ struct femtoDreamProducerReducedTask { fillCollisionsAndTracks(col, tracksWithItsPid); } PROCESS_SWITCH(femtoDreamProducerReducedTask, processMC_noCentrality, "Provide MC data", false); -}; + void processData_FlowCalc(aod::FemtoFullCollision_CentPbPb const& col, + aod::BCsWithTimestamps const&, + aod::FemtoFullTracks const& tracks) + { + // get magnetic field for run + getMagneticFieldTesla(col.bc_as()); + auto tracksWithItsPid = soa::Attach(tracks); + // fill the tables + fillCollisionsAndTracks_PbPb(col, tracksWithItsPid); + if (qnCal.ConfgQnSeparation) { + fillCollisionsFlow(col, tracksWithItsPid); + } + } + PROCESS_SWITCH(femtoDreamProducerReducedTask, processData_FlowCalc, + "Provide experimental data with cumulant flow calculation", false); +}; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { WorkflowSpec workflow{adaptAnalysisTask(cfgc)};