diff --git a/PWGCF/EbyEFluctuations/Tasks/CMakeLists.txt b/PWGCF/EbyEFluctuations/Tasks/CMakeLists.txt index 25d872112b4..f17be4f5601 100644 --- a/PWGCF/EbyEFluctuations/Tasks/CMakeLists.txt +++ b/PWGCF/EbyEFluctuations/Tasks/CMakeLists.txt @@ -88,3 +88,8 @@ o2physics_add_dpl_workflow(strongly-intensive-corr SOURCES stronglyIntensiveCorr.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGCFCore COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(multandpt-fluctuations + SOURCES MultandPtFluctuations.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGCFCore + COMPONENT_NAME Analysis) diff --git a/PWGCF/EbyEFluctuations/Tasks/MultandPtFluctuations.cxx b/PWGCF/EbyEFluctuations/Tasks/MultandPtFluctuations.cxx new file mode 100644 index 00000000000..6736eaa667c --- /dev/null +++ b/PWGCF/EbyEFluctuations/Tasks/MultandPtFluctuations.cxx @@ -0,0 +1,247 @@ +// 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. +/// + +/// \file MultandPtFluctuations.cxx +/// \brief Calculate multiplicity and transverse momentum fluctuations using strongly intensive observables +/// \author Omama Rubza + +#include "Common/Core/TrackSelection.h" +#include "Common/DataModel/Centrality.h" +#include "Common/DataModel/EventSelection.h" +#include "Common/DataModel/FT0Corrected.h" +#include "Common/DataModel/Multiplicity.h" +#include "Common/DataModel/TrackSelectionTables.h" + +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisDataModel.h" +#include "Framework/AnalysisTask.h" +#include "Framework/HistogramRegistry.h" +#include "Framework/HistogramSpec.h" +#include "Framework/runDataProcessing.h" + +#include "TProfile.h" + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::aod; +using namespace o2::framework::expressions; + +struct MultandPtFluctuations { + + // ------ Histogram binning + + Configurable nBinsPt{"nBinsPt", 100, "pT bins"}; + Configurable ndcaXY{"ndcaXY", 100, "DCAxy bins"}; + Configurable ndcaZ{"ndcaZ", 100, "DCAz bins"}; + Configurable nCentBins{"nCentBins", 100, "Number of centrality bins"}; + + // ------ Event Cuts + Configurable vtxZcut{"vtxZcut", 10.0f, "vertex Z"}; + Configurable cfgNoSameBunchPileup{"cfgNoSameBunchPileup", true, "kNoSameBunchPileup"}; + Configurable cfgEvSelUseGoodZvtxFT0vsPV{"cfgEvSelUseGoodZvtxFT0vsPV", true, "Good Zvtx FT0 vs PV"}; + Configurable cfgUseGoodITSLayerAllCut{"cfgUseGoodITSLayerAllCut", true, "Good ITS Layers"}; + + // ----- Centrality estimator + Configurable cFT0M{"cFT0M", false, "FT0M centrality"}; // for pp + Configurable cFT0C{"cFT0C", true, "Use FT0C centrality"}; + + // ------ Track cuts + + Configurable ptMinCut{"ptMinCut", 0.2f, "Minimum pT"}; + Configurable ptMaxCut{"ptMaxCut", 5.0f, "Maximum pT"}; + Configurable etaCut{"etaCut", 0.8f, "Maximum |eta|"}; + Configurable itsNClsCut{"itsNClsCut", 4, "Minimum ITS clusters"}; + Configurable tpcCrossCut{"tpcCrossCut", 70.0f, "TPC crossed rows"}; + Configurable crossedRowsOverFindableCut{"crossedRowsOverFindableCut", 0.8f, "TPC crossed rows over findable"}; + Configurable dcaZCut{"dcaZCut", 2.0f, "Maximum DCAz"}; + Configurable dcaXYCut{"dcaXYCut", 0.2f, "Maximum DCAxy"}; + Configurable tpcChiCut{"tpcChiCut", 4.0f, "TPC chi2/NCls"}; + Configurable itsChiCut{"itsChiCut", 36.0f, "ITS chi2/NCls"}; + Configurable requireITS{"requireITS", true, "Require ITS hit"}; + Configurable requireTPC{"requireTPC", true, "Require TPC hit"}; + Configurable requireInnerITS{"requireInnerITS", true, "At least one hit in ITS layers 0,1,2"}; + + //----- Histogram Registry + HistogramRegistry histos{ + "histos", + {}, + OutputObjHandlingPolicy::AnalysisObject}; + + void init(InitContext const&) + { + const AxisSpec axisEta{30, -1.5, +1.5, "#eta"}; + const AxisSpec axisPt{nBinsPt, 0.15, 10.0, "p_{T}"}; + const AxisSpec axisVtxz{80, -20, 20, "V_{Z} (cm)"}; + const AxisSpec axisPhi{40, -1, 7, "#phi"}; + const AxisSpec axisdcaXY{ndcaXY, -0.3, 0.3, "DCAxy"}; + const AxisSpec axisdcaZ{ndcaZ, -3.0, 3.0, "DCAz"}; + const AxisSpec axisNch{500, 0, 500, "Nch"}; + const AxisSpec axisCent{nCentBins, 0.0, 100.0, "Centrality (%)"}; // ------------- centrality bins + + histos.add("QA/BeforeCut/VtxZ", "Vertex Z", kTH1F, {axisVtxz}); + histos.add("QA/AfterCut/VtxZ", "Vertex Z", kTH1F, {axisVtxz}); + histos.add("QA/BeforeCut/Cent", "Centrality", kTH1F, {axisCent}); + histos.add("QA/AfterCut/Cent", "Centrality", kTH1F, {axisCent}); + histos.add("QA/BeforeCut/Eta", "Eta", kTH1F, {axisEta}); + histos.add("QA/AfterCut/Eta", "Eta", kTH1F, {axisEta}); + histos.add("QA/BeforeCut/Pt", "Pt", kTH1F, {axisPt}); + histos.add("QA/AfterCut/Pt", "Pt", kTH1F, {axisPt}); + histos.add("QA/BeforeCut/Phi", "Phi", kTH1F, {axisPhi}); + histos.add("QA/AfterCut/Phi", "Phi", kTH1F, {axisPhi}); + histos.add("QA/BeforeCut/DcaXY", "DCAxy", kTH1F, {axisdcaXY}); + histos.add("QA/AfterCut/DcaXY", "DCAxy", kTH1F, {axisdcaXY}); + histos.add("QA/BeforeCut/DcaZ", "DCAz", kTH1F, {axisdcaZ}); + histos.add("QA/AfterCut/DcaZ", "DCAz", kTH1F, {axisdcaZ}); + histos.add("hNch", "Nch", kTH1F, {axisNch}); + histos.add("h2_DcaZ", "DCA_{Z}", kTH2D, {{axisPt}, {axisdcaZ}}); + histos.add("h2_DcaXY", "DCA_{XY}", kTH2D, {{axisPt}, {axisdcaXY}}); + histos.add("hEventCounter", "Number of events vs centrality", kTH1F, {axisCent}); + histos.add("p_a", " vs centrality", kTProfile, {axisCent}); + histos.add("p_b", " vs centrality", kTProfile, {axisCent}); + histos.add("p_a2", " vs centrality", kTProfile, {axisCent}); + histos.add("p_b2", " vs centrality", kTProfile, {axisCent}); + histos.add("p_ab", " vs centrality", kTProfile, {axisCent}); + histos.add("p_asumb", " vs centrality", kTProfile, {axisCent}); + } + + //---------- Filters + Filter ptFilter = (aod::track::pt > ptMinCut) && (aod::track::pt < ptMaxCut); + Filter etaFilter = nabs(aod::track::eta) < etaCut; + Filter posZFilter = nabs(aod::collision::posZ) < vtxZcut; + + using myColsData = soa::Join; + using myTracksData = soa::Join; + using myFilteredTracksData = soa::Filtered; + + void process(myColsData::iterator const& col, myFilteredTracksData const& tracks) + { + histos.fill(HIST("QA/BeforeCut/VtxZ"), col.posZ()); + + float cent = -1.0f; + + if (cFT0M) { + cent = col.centFT0M(); + } + + if (cFT0C) { + cent = col.centFT0C(); + // if (!track.isGlobalTrack()) continue; + } + + if (cent < 0) { + return; + } + + histos.fill(HIST("QA/BeforeCut/Cent"), cent); + + if (!col.sel8()) + return; + + if (std::abs(col.posZ()) > vtxZcut) + return; + + if (cfgNoSameBunchPileup && !col.selection_bit(o2::aod::evsel::kNoSameBunchPileup)) + return; + + if (cfgEvSelUseGoodZvtxFT0vsPV && !col.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV)) + return; + + if (cfgUseGoodITSLayerAllCut && !col.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll)) + return; + + histos.fill(HIST("hEventCounter"), cent); + + histos.fill(HIST("QA/AfterCut/VtxZ"), col.posZ()); + histos.fill(HIST("QA/AfterCut/Cent"), cent); + + double a = 0.0; + double b = 0.0; + + // if (!track.isGlobalTrack()) continue; + for (const auto& track : tracks) { + + histos.fill(HIST("QA/BeforeCut/Eta"), track.eta()); + histos.fill(HIST("QA/BeforeCut/Pt"), track.pt()); + histos.fill(HIST("QA/BeforeCut/Phi"), track.phi()); + histos.fill(HIST("QA/BeforeCut/DcaXY"), track.dcaXY()); + histos.fill(HIST("QA/BeforeCut/DcaZ"), track.dcaZ()); + + // Track quality cuts + + if (requireITS && !track.hasITS()) + continue; + if (requireTPC && !track.hasTPC()) + continue; + + if (track.itsNCls() < itsNClsCut) + continue; + if (track.tpcNClsCrossedRows() < tpcCrossCut) + continue; + + if (track.tpcCrossedRowsOverFindableCls() < crossedRowsOverFindableCut) + continue; + + if (std::abs(track.dcaZ()) > dcaZCut) + continue; + if (std::abs(track.dcaXY()) > dcaXYCut) + continue; + + if (track.tpcChi2NCl() > tpcChiCut) + continue; + if (track.itsChi2NCl() > itsChiCut) + continue; + + if (std::abs(track.eta()) > etaCut) + continue; + + if (requireInnerITS) { + + auto itsMap = track.itsClusterMap(); + + if (!(itsMap & (1 << 0)) && + !(itsMap & (1 << 1)) && + !(itsMap & (1 << 2))) { + continue; + } + } + + // After cuts QA + histos.fill(HIST("QA/AfterCut/Eta"), track.eta()); + histos.fill(HIST("QA/AfterCut/Pt"), track.pt()); + histos.fill(HIST("QA/AfterCut/Phi"), track.phi()); + histos.fill(HIST("QA/AfterCut/DcaXY"), track.dcaXY()); + histos.fill(HIST("QA/AfterCut/DcaZ"), track.dcaZ()); + + histos.fill(HIST("h2_DcaZ"), track.pt(), track.dcaZ()); + histos.fill(HIST("h2_DcaXY"), track.pt(), track.dcaXY()); + + a++; // A-nch, B =pt + b += track.pt(); + } + + // ---- Fill TProfiles (once per event) + histos.fill(HIST("p_a"), cent, a); + histos.fill(HIST("p_b"), cent, b); + histos.fill(HIST("p_a2"), cent, a * a); + histos.fill(HIST("p_b2"), cent, b * b); + histos.fill(HIST("p_ab"), cent, a * b); + histos.fill(HIST("p_asumb"), cent, a + b); + } +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{ + adaptAnalysisTask(cfgc)}; +}