diff --git a/DPG/Tasks/AOTTrack/PID/HMPID/CMakeLists.txt b/DPG/Tasks/AOTTrack/PID/HMPID/CMakeLists.txt index d33a73eb12f..748576f57eb 100644 --- a/DPG/Tasks/AOTTrack/PID/HMPID/CMakeLists.txt +++ b/DPG/Tasks/AOTTrack/PID/HMPID/CMakeLists.txt @@ -20,7 +20,13 @@ o2physics_add_dpl_workflow(hmpid-table-producer PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::HMPIDBase COMPONENT_NAME Analysis) -o2physics_add_dpl_workflow(hmpid-deuteron - SOURCES hmpidDeuteron.cxx +o2physics_add_dpl_workflow(hmpid-abs-cross-section + SOURCES hmpidAbsCrossSection.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2::ReconstructionDataFormats O2Physics::AnalysisCore COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(hmpid-cross-section-corrections + SOURCES hmpidCrossSectionCorrections.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2::ReconstructionDataFormats O2Physics::AnalysisCore + COMPONENT_NAME Analysis) + diff --git a/DPG/Tasks/AOTTrack/PID/HMPID/hmpidAbsCrossSection.cxx b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidAbsCrossSection.cxx new file mode 100644 index 00000000000..aa8a79e6bc2 --- /dev/null +++ b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidAbsCrossSection.cxx @@ -0,0 +1,389 @@ +// 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 Analysis task for the absorption cross section using HMPID (de,pi,ka,pr) +// \authors rocco.liotino@cern.ch + +#include "tableHMPID.h" + +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +using namespace o2; +using namespace o2::framework; + +struct HmpidAbsCrossSection { + + HistogramRegistry registryDA{"registryDA", {}, OutputObjHandlingPolicy::AnalysisObject}; + + Configurable minReqClusterITS{"minReqClusterITS", 4.0, "min number of clusters required in ITS"}; + Configurable minTPCnClsFound{"minTPCnClsFound", 50.0f, "minTPCnClsFound"}; + Configurable minNCrossedRowsTPC{"minNCrossedRowsTPC", 70.0f, "min number of crossed rows TPC"}; + Configurable maxChi2ITS{"maxChi2ITS", 36.0f, "max chi2 per cluster ITS"}; + Configurable maxChi2TPC{"maxChi2TPC", 4.0f, "max chi2 per cluster TPC"}; + Configurable maxDCAxy{"maxDCAxy", 0.5f, "maxDCAxy"}; + Configurable maxDCAz{"maxDCAz", 0.5f, "maxDCAz"}; + + // Configurables for PbPb collisions + Configurable collisionsPbPb{"collisionsPbPb", false, "Require PbPb collisions"}; + Configurable centralityCut{"centralityCut", 30.0f, "centrality cut for PbPb collisions"}; + + // Configurables for PID cuts + Configurable nsigmaTPCMin{"nsigmaTPCMin", -3.0, "nsigmaTPCMin"}; + Configurable nsigmaTPCMax{"nsigmaTPCMax", +3.0, "nsigmaTPCMax"}; + Configurable nsigmaTOFMin{"nsigmaTOFMin", -3.0, "nsigmaTOFMin"}; + Configurable nsigmaTOFMax{"nsigmaTOFMax", +3.0, "nsigmaTOFMax"}; + + // Particle species selection + Configurable particleSpecies{"particleSpecies", "de", "Particle species: de (deuteron), pr (proton), pi (pion), ka (kaon)"}; + + // Configurable offset for residuals distributions + Configurable offsetX{"offsetX", 0.696712f, "Residuals offset - X direction"}; // mean value (old) + Configurable offsetY{"offsetY", 2.32849f, "Residuals offset - Y direction"}; // mean value (old) + + static constexpr float UnsetOffset = -999.f; + static constexpr int MaxOffsetBins = 35; + static constexpr int OffsetColumns = 8; + + static std::vector makeOffsetRowLabels() + { + std::vector labels; + labels.reserve(MaxOffsetBins); + for (int bin = 0; bin < MaxOffsetBins; ++bin) { + labels.emplace_back("pBin" + std::to_string(bin)); + } + return labels; + } + + Configurable offsetMomMin{"offsetMomMin", 0.5f, "Minimum momentum for residual offset correction"}; + Configurable offsetMomMax{"offsetMomMax", 4.0f, "Maximum momentum for residual offset correction"}; + Configurable offsetMomStep{"offsetMomStep", 0.5f, "Momentum bin width for residual offset correction"}; + + Configurable> residualOffsets{"residualOffsets", LabeledArray{Array2D{std::vector(MaxOffsetBins * OffsetColumns, UnsetOffset), MaxOffsetBins, OffsetColumns}, makeOffsetRowLabels(), {"xPosR2", "yPosR2", "xNegR2", "yNegR2", "xPosR4", "yPosR4", "xNegR4", "yNegR4"}}, "Residual offsets by momentum bin, chamber and charge"}; + + // Configurables for plotting + Configurable momNBins{"momNBins", 450, "Momentum bins"}; + Configurable momMin{"momMin", 0.5f, "Momentum min"}; + Configurable momMax{"momMax", 5.0f, "Momentum max"}; + + Configurable deltaRNBins{"deltaRNBins", 100, "DeltaR bins"}; + Configurable deltaRMin{"deltaRMin", 0.0f, "DeltaR min"}; + Configurable deltaRMax{"deltaRMax", 50.0f, "DeltaR max"}; + + Configurable nSigmaTPCNBins{"nSigmaTPCNBins", 100, "nSigmaTPC bins"}; + Configurable nSigmaTPCMinPlot{"nSigmaTPCMinPlot", -5.0f, "nSigmaTPC min"}; + Configurable nSigmaTPCMaxPlot{"nSigmaTPCMaxPlot", 5.0f, "nSigmaTPC max"}; + + Configurable nSigmaTOFNBins{"nSigmaTOFNBins", 100, "nSigmaTOF bins"}; + Configurable nSigmaTOFMinPlot{"nSigmaTOFMinPlot", -5.0f, "nSigmaTOF min"}; + Configurable nSigmaTOFMaxPlot{"nSigmaTOFMaxPlot", 5.0f, "nSigmaTOF max"}; + + static constexpr int Rich2 = 2, Rich4 = 4; + // static constexpr Rich0 = 0, Rich1 = 1, Rich3 = 3, Rich5 = 5, Rich6 = 6; + // static constexpr int Nch = 7; + + void init(o2::framework::InitContext const&) + { + // Axes for histograms + AxisSpec axisMom{momNBins, momMin, momMax, "#it{p} (GeV/#it{c})", "axisMom"}; + AxisSpec axisDeltaR{deltaRNBins, deltaRMin, deltaRMax, "#Delta R (cm)", "axisDeltaR"}; + AxisSpec axisNsigmaTOF{nSigmaTOFNBins, nSigmaTOFMinPlot, nSigmaTOFMaxPlot, "n#sigma_{TOF}", "axisNsigmaTOF"}; + AxisSpec axisNsigmaTPC{nSigmaTPCNBins, nSigmaTPCMinPlot, nSigmaTPCMaxPlot, "n#sigma_{TPC}", "axisNsigmaTPC"}; + AxisSpec axisChamber{10, -1.5, 8.5, "Chamber number", "axisChamber"}; + + // 1D integrated histograms (quick view) + registryDA.add("incoming_Pos_8cm", "incoming_Pos_8cm", HistType::kTH1F, {axisMom}); + registryDA.add("incoming_Pos_4cm", "incoming_Pos_4cm", HistType::kTH1F, {axisMom}); + + registryDA.add("DeltaR_Pos_8cm", "DeltaR_Pos_8cm", HistType::kTH1F, {axisDeltaR}); + registryDA.add("DeltaR_Pos_4cm", "DeltaR_Pos_4cm", HistType::kTH1F, {axisDeltaR}); + + registryDA.add("incoming_Neg_8cm", "incoming_Neg_8cm", HistType::kTH1F, {axisMom}); + registryDA.add("incoming_Neg_4cm", "incoming_Neg_4cm", HistType::kTH1F, {axisMom}); + + registryDA.add("DeltaR_Neg_8cm", "DeltaR_Neg_8cm", HistType::kTH1F, {axisDeltaR}); + registryDA.add("DeltaR_Neg_4cm", "DeltaR_Neg_4cm", HistType::kTH1F, {axisDeltaR}); + + // Positive candidates histograms + registryDA.add("nSigmaTPC_vs_momHMPID_PreCut_Pos_8cm", "nSigmaTPC_vs_momHMPID_PreCut_Pos_8cm", HistType::kTH2F, {axisMom, axisNsigmaTPC}); + registryDA.add("nSigmaTOF_vs_momHMPID_PreCut_Pos_8cm", "nSigmaTOF_vs_momHMPID_PreCut_Pos_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + registryDA.add("nSigmaTPC_vs_momHMPID_PreCut_Pos_4cm", "nSigmaTPC_vs_momHMPID_PreCut_Pos_4cm", HistType::kTH2F, {axisMom, axisNsigmaTPC}); + registryDA.add("nSigmaTOF_vs_momHMPID_PreCut_Pos_4cm", "nSigmaTOF_vs_momHMPID_PreCut_Pos_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + + // TPC only candidates histograms + registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_Pos_8cm", "nSigmaTOF_vs_momHMPID_TPCOnly_Pos_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_Pos_4cm", "nSigmaTOF_vs_momHMPID_TPCOnly_Pos_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + + // after TPC+TOF cuts candidates histograms + registryDA.add("nSigmaTOF_vs_momHMPID_Pos_8cm", "nSigmaTOF_vs_momHMPID_Pos_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + registryDA.add("nSigmaTOF_vs_momHMPID_Pos_4cm", "nSigmaTOF_vs_momHMPID_Pos_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + + registryDA.add("DeltaR_vs_momHMPID_Pos_8cm", "DeltaR_vs_momHMPID_Pos_8cm", HistType::kTH2F, {axisMom, axisDeltaR}); + registryDA.add("DeltaR_vs_momHMPID_Pos_4cm", "DeltaR_vs_momHMPID_Pos_4cm", HistType::kTH2F, {axisMom, axisDeltaR}); + + // corrected residuals distributions + registryDA.add("residualsX_Vs_momentum_Pos_8cm", "residualsX_Vs_momentum_Pos_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); + registryDA.add("residualsY_Vs_momentum_Pos_8cm", "residualsY_Vs_momentum_Pos_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); + registryDA.add("residualsX_Vs_momentum_Pos_4cm", "residualsX_Vs_momentum_Pos_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); + registryDA.add("residualsY_Vs_momentum_Pos_4cm", "residualsY_Vs_momentum_Pos_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); + + // Negative candidates histograms + registryDA.add("nSigmaTPC_vs_momHMPID_PreCut_Neg_8cm", "nSigmaTPC_vs_momHMPID_PreCut_Neg_8cm", HistType::kTH2F, {axisMom, axisNsigmaTPC}); + registryDA.add("nSigmaTOF_vs_momHMPID_PreCut_Neg_8cm", "nSigmaTOF_vs_momHMPID_PreCut_Neg_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + registryDA.add("nSigmaTPC_vs_momHMPID_PreCut_Neg_4cm", "nSigmaTPC_vs_momHMPID_PreCut_Neg_4cm", HistType::kTH2F, {axisMom, axisNsigmaTPC}); + registryDA.add("nSigmaTOF_vs_momHMPID_PreCut_Neg_4cm", "nSigmaTOF_vs_momHMPID_PreCut_Neg_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + + // TPC only candidates histograms + registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_Neg_8cm", "nSigmaTOF_vs_momHMPID_TPCOnly_Neg_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_Neg_4cm", "nSigmaTOF_vs_momHMPID_TPCOnly_Neg_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + + // after TPC+TOF cuts candidates histograms + registryDA.add("nSigmaTOF_vs_momHMPID_Neg_8cm", "nSigmaTOF_vs_momHMPID_Neg_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + registryDA.add("nSigmaTOF_vs_momHMPID_Neg_4cm", "nSigmaTOF_vs_momHMPID_Neg_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); + + registryDA.add("DeltaR_vs_momHMPID_Neg_8cm", "DeltaR_vs_momHMPID_Neg_8cm", HistType::kTH2F, {axisMom, axisDeltaR}); + registryDA.add("DeltaR_vs_momHMPID_Neg_4cm", "DeltaR_vs_momHMPID_Neg_4cm", HistType::kTH2F, {axisMom, axisDeltaR}); + + // corrected residuals distributions + registryDA.add("residualsX_Vs_momentum_Neg_8cm", "residualsX_Vs_momentum_Neg_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); + registryDA.add("residualsY_Vs_momentum_Neg_8cm", "residualsY_Vs_momentum_Neg_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); + registryDA.add("residualsX_Vs_momentum_Neg_4cm", "residualsX_Vs_momentum_Neg_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); + registryDA.add("residualsY_Vs_momentum_Neg_4cm", "residualsY_Vs_momentum_Neg_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); + + // nsigma TPC vs nsigma TOF for candidates + registryDA.add("pidTPC_TOF_p_ch_Pos", "pidTPC_TOF_p_ch_Pos", HistType::kTHnSparseF, {axisMom, axisNsigmaTPC, axisNsigmaTOF, axisChamber}); + registryDA.add("pidTPC_TOF_p_ch_Neg", "pidTPC_TOF_p_ch_Neg", HistType::kTHnSparseF, {axisMom, axisNsigmaTPC, axisNsigmaTOF, axisChamber}); + } // end init + + void process(aod::HmpidAnalysis const& hmpidTable) + { + for (const auto& hmpid : hmpidTable) { + + // ------------------------- + // track filters + // ------------------------- + if (hmpid.itsNCluster() < minReqClusterITS) { + continue; + } + if (hmpid.tpcNCluster() < minTPCnClsFound) { + continue; + } + if (hmpid.tpcNClsCrossedRows() < minNCrossedRowsTPC) { + continue; + } + if (hmpid.tpcChi2() > maxChi2TPC) { + continue; + } + if (hmpid.itsChi2() > maxChi2ITS) { + continue; + } + if (std::abs(hmpid.dcaXY()) > maxDCAxy) { + continue; + } + if (std::abs(hmpid.dcaZ()) > maxDCAz) { + continue; + } + + // ------------------------- + // derived quantities + // ------------------------- + const float momHmpid = hmpid.momentumHmpid(); + const float momAbs = std::abs(momHmpid); + + const int chamber = static_cast(hmpid.chamber()); + const bool isPos = momHmpid > 0; + const bool isNeg = momHmpid < 0; + + // Select chambers with absorber + if (chamber != Rich2 && chamber != Rich4) { + continue; + } + + // Default residual offsets + float usedOffsetX = offsetX.value; + float usedOffsetY = offsetY.value; + + // Momentum-dependent offset correction + if (momAbs >= offsetMomMin.value && momAbs < offsetMomMax.value && (isPos || isNeg)) { + const int pBin = static_cast((momAbs - offsetMomMin.value) / offsetMomStep.value); + if (pBin >= 0 && pBin < static_cast(residualOffsets->rows())) { + const int col = (chamber == Rich2 ? 0 : 4) + (isPos ? 0 : 2); + const float binnedOffsetX = residualOffsets->get(pBin, col); + const float binnedOffsetY = residualOffsets->get(pBin, col + 1); + + if (binnedOffsetX != UnsetOffset) { + usedOffsetX = binnedOffsetX; + } + if (binnedOffsetY != UnsetOffset) { + usedOffsetY = binnedOffsetY; + } + } + } + + // Corrected residuals + const float dxCorr = (hmpid.xMip() - hmpid.xTrack()) - usedOffsetX; + const float dyCorr = (hmpid.yMip() - hmpid.yTrack()) - usedOffsetY; + const float drCorr = std::hypot(dxCorr, dyCorr); + + // check only HMPID chamber with an absorber (only RICH 2 (4cm abs) and RICH4 (8cm abs)) + if (chamber != Rich2 && chamber != Rich4) { + continue; + } + + // select particle species + + float nSigmaTPC = 0.f; + float nSigmaTOF = 0.f; + + // Select particle species + if (particleSpecies.value == "de") { + nSigmaTPC = hmpid.tpcNSigmaDe(); + nSigmaTOF = hmpid.tofNSigmaDe(); + } else if (particleSpecies.value == "pr") { + nSigmaTPC = hmpid.tpcNSigmaPr(); + nSigmaTOF = hmpid.tofNSigmaPr(); + } else if (particleSpecies.value == "pi") { + nSigmaTPC = hmpid.tpcNSigmaPi(); + nSigmaTOF = hmpid.tofNSigmaPi(); + } else if (particleSpecies.value == "ka") { + nSigmaTPC = hmpid.tpcNSigmaKa(); + nSigmaTOF = hmpid.tofNSigmaKa(); + } else { + LOGF(fatal, "Invalid particle species: %s", + particleSpecies.value.c_str()); + } + + // ------------------------- + // nsigma distributions pre-cuts + // ------------------------- + if (isPos) { + registryDA.fill(HIST("pidTPC_TOF_p_ch_Pos"), momAbs, nSigmaTPC, nSigmaTOF, chamber); + if (chamber == Rich2) { + registryDA.fill(HIST("nSigmaTPC_vs_momHMPID_PreCut_Pos_4cm"), momAbs, nSigmaTPC); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_PreCut_Pos_4cm"), momAbs, nSigmaTOF); + } + if (chamber == Rich4) { + registryDA.fill(HIST("nSigmaTPC_vs_momHMPID_PreCut_Pos_8cm"), momAbs, nSigmaTPC); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_PreCut_Pos_8cm"), momAbs, nSigmaTOF); + } + } + if (isNeg) { + registryDA.fill(HIST("pidTPC_TOF_p_ch_Neg"), momAbs, nSigmaTPC, nSigmaTOF, chamber); + if (chamber == Rich2) { + registryDA.fill(HIST("nSigmaTPC_vs_momHMPID_PreCut_Neg_4cm"), momAbs, nSigmaTPC); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_PreCut_Neg_4cm"), momAbs, nSigmaTOF); + } + if (chamber == Rich4) { + registryDA.fill(HIST("nSigmaTPC_vs_momHMPID_PreCut_Neg_8cm"), momAbs, nSigmaTPC); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_PreCut_Neg_8cm"), momAbs, nSigmaTOF); + } + } + + // ------------------------- + // deuteron candidate cuts - TPC + // ------------------------- + if (nSigmaTPC < nsigmaTPCMin || nSigmaTPC > nsigmaTPCMax) { + continue; + } + + // ------------------------- + // post TPC cut + // ------------------------- + if (isPos) { + if (chamber == Rich2) { + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_Pos_8cm"), momAbs, nSigmaTOF); + } + if (chamber == Rich4) { + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_Pos_4cm"), momAbs, nSigmaTOF); + } + } + if (isNeg) { + if (chamber == Rich2) { + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_Neg_8cm"), momAbs, nSigmaTOF); + } + if (chamber == Rich4) { + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_Neg_4cm"), momAbs, nSigmaTOF); + } + } + + // ------------------------- + // deuteron candidate cuts - TOF + // ------------------------- + if (nSigmaTOF < nsigmaTOFMin || nSigmaTOF > nsigmaTOFMax) { + continue; + } + + // Pb-Pb collisions centrality requirement + if (collisionsPbPb) { + if (hmpid.centrality() < centralityCut) { + continue; + } + } + + // Fill distributions for absorption cross section measurement + if (isPos && chamber == Rich4) { + registryDA.fill(HIST("incoming_Pos_8cm"), momAbs); + registryDA.fill(HIST("DeltaR_Pos_8cm"), drCorr); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_Pos_8cm"), momAbs, nSigmaTOF); + registryDA.fill(HIST("DeltaR_vs_momHMPID_Pos_8cm"), momAbs, drCorr); + registryDA.fill(HIST("residualsX_Vs_momentum_Pos_8cm"), momAbs, dxCorr); + registryDA.fill(HIST("residualsY_Vs_momentum_Pos_8cm"), momAbs, dyCorr); + } + + if (isNeg && chamber == Rich4) { + registryDA.fill(HIST("incoming_Neg_8cm"), momAbs); + registryDA.fill(HIST("DeltaR_Neg_8cm"), drCorr); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_Neg_8cm"), momAbs, nSigmaTOF); + registryDA.fill(HIST("DeltaR_vs_momHMPID_Neg_8cm"), momAbs, drCorr); + + registryDA.fill(HIST("residualsX_Vs_momentum_Neg_8cm"), momAbs, dxCorr); + registryDA.fill(HIST("residualsY_Vs_momentum_Neg_8cm"), momAbs, dyCorr); + } + + if (isPos && chamber == Rich2) { + registryDA.fill(HIST("incoming_Pos_4cm"), momAbs); + registryDA.fill(HIST("DeltaR_Pos_4cm"), drCorr); + registryDA.fill(HIST("DeltaR_vs_momHMPID_Pos_4cm"), momAbs, drCorr); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_Pos_4cm"), momAbs, nSigmaTOF); + + registryDA.fill(HIST("residualsX_Vs_momentum_Pos_4cm"), momAbs, dxCorr); + registryDA.fill(HIST("residualsY_Vs_momentum_Pos_4cm"), momAbs, dyCorr); + } + + if (isNeg && chamber == Rich2) { + registryDA.fill(HIST("incoming_Neg_4cm"), momAbs); + registryDA.fill(HIST("DeltaR_Neg_4cm"), drCorr); + registryDA.fill(HIST("DeltaR_vs_momHMPID_Neg_4cm"), momAbs, drCorr); + registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_Neg_4cm"), momAbs, nSigmaTOF); + + registryDA.fill(HIST("residualsX_Vs_momentum_Neg_4cm"), momAbs, dxCorr); + registryDA.fill(HIST("residualsY_Vs_momentum_Neg_4cm"), momAbs, dyCorr); + } + } // end of loop over hmpidTable + } // end of process +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{ + adaptAnalysisTask(cfgc), + }; +} diff --git a/DPG/Tasks/AOTTrack/PID/HMPID/hmpidCrossSectionCorrections.cxx b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidCrossSectionCorrections.cxx new file mode 100644 index 00000000000..8ca94444d60 --- /dev/null +++ b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidCrossSectionCorrections.cxx @@ -0,0 +1,374 @@ +// 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 Analysis task for the (anti-)deuteron absorption cross section using HMPID +// \authors rocco.liotino@cern.ch + +#include "tableHMPID.h" + +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + +using namespace o2; +using namespace o2::framework; + +struct HmpidCrossSectionCorrections { + + HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject}; + + // Configurables for main track quality criteria + Configurable minReqClusterITS{"minReqClusterITS", 4.0, "min number of clusters required in ITS"}; + Configurable minTPCnClsFound{"minTPCnClsFound", 50.0f, "minTPCnClsFound"}; + Configurable minNCrossedRowsTPC{"minNCrossedRowsTPC", 70.0f, "min number of crossed rows TPC"}; + Configurable maxChi2ITS{"maxChi2ITS", 36.0f, "max chi2 per cluster ITS"}; + Configurable maxChi2TPC{"maxChi2TPC", 4.0f, "max chi2 per cluster TPC"}; + Configurable maxDCAxy{"maxDCAxy", 0.5f, "maxDCAxy"}; + Configurable maxDCAz{"maxDCAz", 0.5f, "maxDCAz"}; + + // Configurables for QA selection + Configurable qMipMin{"qMipMin", 200.0f, "minimum MIP charge for QA selection"}; + Configurable qMipMax{"qMipMax", 2000.0f, "maximum MIP charge for QA selection"}; + Configurable cutDistanceMipTrack{"cutDistanceMipTrack", 10.0f, "maximum raw MIP-track distance for QA selection"}; + + // Configurables for plotting + Configurable momNBins{"momNBins", 490, "Momentum bins"}; + Configurable momMin{"momMin", 0.1f, "Momentum min"}; + Configurable momMax{"momMax", 5.0f, "Momentum max"}; + + Configurable etaNBins{"etaNBins", 100, "Eta bins"}; + Configurable etaMin{"etaMin", -5.0f, "Eta min"}; + Configurable etaMax{"etaMax", 5.0f, "Eta max"}; + + Configurable phiNBins{"phiNBins", 100, "Phi bins"}; + Configurable phiMin{"phiMin", -10.0f, "Phi min"}; + Configurable phiMax{"phiMax", 10.0f, "Phi max"}; + + Configurable mapNBins{"mapNBins", 135, "Map bins (x or y)"}; + Configurable mapMin{"mapMin", 0.0f, "Map min (x or y) in cm"}; + Configurable mapMax{"mapMax", 135.0f, "Map max (x or y) in cm"}; + + Configurable residualsNBins{"residualsNBins", 400, "Residuals bins"}; + Configurable residualsMin{"residualsMin", -20.0f, "Residuals min (x or y) in cm"}; + Configurable residualsMax{"residualsMax", 20.0f, "Residuals max (x or y) in cm"}; + + static constexpr int Rich0 = 0, Rich1 = 1, Rich2 = 2, Rich3 = 3, Rich4 = 4, Rich5 = 5, Rich6 = 6; + static constexpr int Nch = 7; + + void init(o2::framework::InitContext const&) + { + // Axes for histograms + AxisSpec axisHmpidMom{momNBins, momMin, momMax, "#it{p} (GeV/#it{c})", "axisHmpidMom"}; + AxisSpec axisGlobalMom{momNBins, momMin, momMax, "#it{p} (GeV/#it{c})", "axisGlobalMom"}; + AxisSpec axisMipX{mapNBins, mapMin, mapMax, "X_{MIP} (cm)", "axisMipX"}; + AxisSpec axisMipY{mapNBins, mapMin, mapMax, "Y_{MIP} (cm)", "axisMipY"}; + AxisSpec axisGlobalTrackX{mapNBins, mapMin, mapMax, "X_{track} (cm)", "axisGlobalTrackX"}; + AxisSpec axisGlobalTrackY{mapNBins, mapMin, mapMax, "Y_{track} (cm)", "axisGlobalTrackY"}; + AxisSpec axisEta{etaNBins, etaMin, etaMax, "#eta", "axisEta"}; + AxisSpec axisPhi{phiNBins, phiMin, phiMax, "#phi", "axisPhi"}; + AxisSpec axisResiduals{residualsNBins, residualsMin, residualsMax, "Residuals (cm)", "axisResiduals"}; + AxisSpec axisChamber{10, -1.5, 8.5, "Chamber number", "axisChamber"}; + + for (int iCh = 0; iCh < Nch; iCh++) { + registry.add(Form("hEtaRich%d", iCh), Form("hEtaRich%d", iCh), HistType::kTH1F, {axisEta}); + registry.add(Form("hPhiRich%d", iCh), Form("hPhiRich%d", iCh), HistType::kTH1F, {axisPhi}); + + registry.add(Form("hmpidXYHmpidMipRich%dPos", iCh), "MIP XY - Q_{MIP} cut;X_{MIP} (cm);Y_{MIP} (cm)", HistType::kTHnSparseF, {axisHmpidMom, axisMipX, axisMipY}); + registry.add(Form("hmpidXYGlobalTrackRich%dPos", iCh), "Global track XY - all;X_{track} (cm);Y_{track} (cm)", HistType::kTHnSparseF, {axisGlobalMom, axisGlobalTrackX, axisGlobalTrackY}); + registry.add(Form("hmpidXYHmpidMipRich%dNeg", iCh), "MIP XY - Q_{MIP} cut;X_{MIP} (cm);Y_{MIP} (cm)", HistType::kTHnSparseF, {axisHmpidMom, axisMipX, axisMipY}); + registry.add(Form("hmpidXYGlobalTrackRich%dNeg", iCh), "Global track XY - all;X_{track} (cm);Y_{track} (cm)", HistType::kTHnSparseF, {axisGlobalMom, axisGlobalTrackX, axisGlobalTrackY}); + + registry.add(Form("residualsXVsMomentumRich%dPos", iCh), Form("residualsXVsMomentumRich%dPos", iCh), HistType::kTH2F, {axisHmpidMom, axisResiduals}); + registry.add(Form("residualsYVsMomentumRich%dPos", iCh), Form("residualsYVsMomentumRich%dPos", iCh), HistType::kTH2F, {axisHmpidMom, axisResiduals}); + registry.add(Form("residualsXVsMomentumRich%dNeg", iCh), Form("residualsXVsMomentumRich%dNeg", iCh), HistType::kTH2F, {axisHmpidMom, axisResiduals}); + registry.add(Form("residualsYVsMomentumRich%dNeg", iCh), Form("residualsYVsMomentumRich%dNeg", iCh), HistType::kTH2F, {axisHmpidMom, axisResiduals}); + + // Momentum and Position correlations + registry.add(Form("hMomTrackVsHmpidRich%dPos", iCh), Form("hMomTrackVsHmpidRich%dPos", iCh), HistType::kTH2F, {axisGlobalMom, axisHmpidMom}); + registry.add(Form("hXTrackVsXMipRich%dPos", iCh), Form("hXTrackVsXMipRich%dPos", iCh), HistType::kTH2F, {axisGlobalTrackX, axisMipX}); + registry.add(Form("hYTrackVsYMipRich%dPos", iCh), Form("hYTrackVsYMipRich%dPos", iCh), HistType::kTH2F, {axisGlobalTrackY, axisMipY}); + registry.add(Form("hMomTrackVsHmpidRich%dNeg", iCh), Form("hMomTrackVsHmpidRich%dNeg", iCh), HistType::kTH2F, {axisGlobalMom, axisHmpidMom}); + registry.add(Form("hXTrackVsXMipRich%dNeg", iCh), Form("hXTrackVsXMipRich%dNeg", iCh), HistType::kTH2F, {axisGlobalTrackX, axisMipX}); + registry.add(Form("hYTrackVsYMipRich%dNeg", iCh), Form("hYTrackVsYMipRich%dNeg", iCh), HistType::kTH2F, {axisGlobalTrackY, axisMipY}); + } // end loop over chambers + + } // end init + + void process(aod::HmpidAnalysis const& hmpidTable) + { + for (const auto& hmpid : hmpidTable) { + if (hmpid.itsNCluster() < minReqClusterITS) { + continue; + } + if (hmpid.tpcNCluster() < minTPCnClsFound) { + continue; + } + if (hmpid.tpcNClsCrossedRows() < minNCrossedRowsTPC) { + continue; + } + if (hmpid.tpcChi2() > maxChi2TPC) { + continue; + } + if (hmpid.itsChi2() > maxChi2ITS) { + continue; + } + if (std::abs(hmpid.dcaXY()) > maxDCAxy) { + continue; + } + if (std::abs(hmpid.dcaZ()) > maxDCAz) { + continue; + } + + // ------------------------- + // derived quantities + // ------------------------- + const int chamber = static_cast(hmpid.chamber()); + + const float momHmpid = hmpid.momentumHmpid(); + const float momTrackAbs = hmpid.momentumTrack(); + const float momHmpidAbs = std::abs(momHmpid); + + const bool isPos = momHmpid > 0; + const bool isNeg = momHmpid < 0; + const bool passQ = hmpid.chargeMip() > qMipMin && hmpid.chargeMip() < qMipMax; + + switch (chamber) { + case Rich0: + registry.fill(HIST("hEtaRich0"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich0"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich0Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich0Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich0Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich0Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich0Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich0Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich0Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich0Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich0Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich0Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich0Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich0Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich0Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich0Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + case Rich1: + registry.fill(HIST("hEtaRich1"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich1"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich1Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich1Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich1Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich1Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich1Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich1Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich1Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich1Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich1Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich1Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich1Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich1Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich1Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich1Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + case Rich2: + registry.fill(HIST("hEtaRich2"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich2"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich2Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich2Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich2Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich2Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich2Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich2Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich2Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich2Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich2Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich2Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich2Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich2Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich2Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich2Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + case Rich3: + registry.fill(HIST("hEtaRich3"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich3"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich3Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich3Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich3Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich3Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich3Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich3Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich3Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich3Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich3Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich3Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich3Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich3Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich3Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich3Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + case Rich4: + registry.fill(HIST("hEtaRich4"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich4"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich4Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich4Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich4Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich4Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich4Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich4Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich4Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich4Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich4Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich4Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich4Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich4Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich4Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich4Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + case Rich5: + registry.fill(HIST("hEtaRich5"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich5"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich5Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich5Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich5Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich5Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich5Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich5Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich5Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich5Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich5Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich5Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich5Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich5Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich5Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich5Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + case Rich6: + registry.fill(HIST("hEtaRich6"), hmpid.etaTrack()); + registry.fill(HIST("hPhiRich6"), hmpid.phiTrack()); + if (isPos) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich6Pos"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich6Pos"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich6Pos"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich6Pos"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich6Pos"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich6Pos"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich6Pos"), hmpid.yTrack(), hmpid.yMip()); + } + if (isNeg) { + if (passQ) { + registry.fill(HIST("hmpidXYHmpidMipRich6Neg"), momHmpidAbs, hmpid.xMip(), hmpid.yMip()); + } + registry.fill(HIST("hmpidXYGlobalTrackRich6Neg"), momHmpidAbs, hmpid.xTrack(), hmpid.yTrack()); + + registry.fill(HIST("residualsXVsMomentumRich6Neg"), momHmpidAbs, hmpid.xTrack() - hmpid.xMip()); + registry.fill(HIST("residualsYVsMomentumRich6Neg"), momHmpidAbs, hmpid.yTrack() - hmpid.yMip()); + + registry.fill(HIST("hMomTrackVsHmpidRich6Neg"), momTrackAbs, momHmpidAbs); + registry.fill(HIST("hXTrackVsXMipRich6Neg"), hmpid.xTrack(), hmpid.xMip()); + registry.fill(HIST("hYTrackVsYMipRich6Neg"), hmpid.yTrack(), hmpid.yMip()); + } + break; + default: + break; + } // end switch on chamber + } // end of loop over hmpidTable + } // end process +}; // end struct + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{ + adaptAnalysisTask(cfgc), + }; +} diff --git a/DPG/Tasks/AOTTrack/PID/HMPID/hmpidDeuteron.cxx b/DPG/Tasks/AOTTrack/PID/HMPID/hmpidDeuteron.cxx deleted file mode 100644 index af8f96dac9a..00000000000 --- a/DPG/Tasks/AOTTrack/PID/HMPID/hmpidDeuteron.cxx +++ /dev/null @@ -1,620 +0,0 @@ -// 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 Analysis task for the (anti-)deuteron absorption cross section using HMPID -// \authors rocco.liotino@cern.ch - -#include "tableHMPID.h" - -#include -#include -#include -#include -#include -#include -#include - -#include - -#include - -using namespace o2; -using namespace o2::framework; - -struct HmpidDeuteron { - HistogramRegistry registryDA{"registryDA", {}, OutputObjHandlingPolicy::AnalysisObject}; - - Configurable nsigmaTPCMin{"nsigmaTPCMin", -3.0, "nsigmaTPCMin"}; - Configurable nsigmaTPCMax{"nsigmaTPCMax", +3.0, "nsigmaTPCMax"}; - Configurable nsigmaTOFMin{"nsigmaTOFMin", -3.0, "nsigmaTOFMin"}; - Configurable nsigmaTOFMax{"nsigmaTOFMax", +3.0, "nsigmaTOFMax"}; - Configurable minReqClusterITS{"minReqClusterITS", 4.0, "min number of clusters required in ITS"}; - Configurable minTPCnClsFound{"minTPCnClsFound", 50.0f, "minTPCnClsFound"}; - Configurable minNCrossedRowsTPC{"minNCrossedRowsTPC", 70.0f, "min number of crossed rows TPC"}; - Configurable maxChi2ITS{"maxChi2ITS", 36.0f, "max chi2 per cluster ITS"}; - Configurable maxChi2TPC{"maxChi2TPC", 4.0f, "max chi2 per cluster TPC"}; - Configurable maxDCAxy{"maxDCAxy", 0.5f, "maxDCAxy"}; - Configurable maxDCAz{"maxDCAz", 0.5f, "maxDCAz"}; - - Configurable qMipThreshold{"qMipThreshold", 120.0f, "minimum MIP charge required for a real signal"}; - Configurable cutDistanceMipTrack{"cutDistanceMipTrack", 6.0f, "cut distance between MIP and track"}; - - Configurable offsetX{"offsetX", 0.696712f, "Residuals offset - X direction"}; - Configurable offsetY{"offsetY", 2.32849f, "Residuals offset - Y direction"}; - - Configurable collisionsPbPb{"collisionsPbPb", false, "Require PbPb collisions"}; - - Configurable centralityCut{"centralityCut", 20.0f, "centrality cut for PbPb collisions"}; - - // Configurables for plotting - // Momentum - Configurable momNBins{"momNBins", 490, "Momentum bins"}; - Configurable momMin{"momMin", 0.1f, "Momentum min"}; - Configurable momMax{"momMax", 5.0f, "Momentum max"}; - - // DeltaR - Configurable deltaRNBins{"deltaRNBins", 50, "DeltaR bins"}; - Configurable deltaRMin{"deltaRMin", 0.0f, "DeltaR min"}; - Configurable deltaRMax{"deltaRMax", 50.0f, "DeltaR max"}; - - // nSigma TPC - Configurable nSigmaTPCNBins{"nSigmaTPCNBins", 100, "nSigmaTPC bins"}; - Configurable nSigmaTPCMinPlot{"nSigmaTPCMinPlot", -5.0f, "nSigmaTPC min"}; - Configurable nSigmaTPCMaxPlot{"nSigmaTPCMaxPlot", 5.0f, "nSigmaTPC max"}; - - // nSigma TOF - Configurable nSigmaTOFNBins{"nSigmaTOFNBins", 100, "nSigmaTOF bins"}; - Configurable nSigmaTOFMinPlot{"nSigmaTOFMinPlot", -5.0f, "nSigmaTOF min"}; - Configurable nSigmaTOFMaxPlot{"nSigmaTOFMaxPlot", 5.0f, "nSigmaTOF max"}; - - // Charge - Configurable chargeNBins{"chargeNBins", 1000, "Charge bins"}; - Configurable chargeMin{"chargeMin", 200.0f, "Charge min"}; - Configurable chargeMax{"chargeMax", 2200.0f, "Charge max"}; - - // Eta - Configurable etaNBins{"etaNBins", 100, "Eta bins"}; - Configurable etaMin{"etaMin", -5.0f, "Eta min"}; - Configurable etaMax{"etaMax", 5.0f, "Eta max"}; - - // Phi - Configurable phiNBins{"phiNBins", 100, "Phi bins"}; - Configurable phiMin{"phiMin", -10.0f, "Phi min"}; - Configurable phiMax{"phiMax", 10.0f, "Phi max"}; - - // Mass - Configurable massNBins{"massNBins", 1000, "Mass bins"}; - Configurable massMin{"massMin", -10.0f, "Mass min"}; - Configurable massMax{"massMax", 10.0f, "Mass max"}; - - // variables for chamber_number and HVs/PCs - static constexpr int Rich0 = 0, Rich1 = 1, Rich2 = 2, Rich3 = 3, Rich4 = 4, Rich5 = 5, Rich6 = 6; - static constexpr int Nch = 7; - static constexpr float Nmean = 1.28947; // meanIdxRad(); method from param not working, using value from HMPIDBase/Param.cxx - - void init(o2::framework::InitContext const&) - { - - // Axes for histograms - AxisSpec axisMom{momNBins, momMin, momMax, "#it{p} (GeV/#it{c})", "axisMom"}; - AxisSpec axisDeltaR{deltaRNBins, deltaRMin, deltaRMax, "#Delta R (cm)", "axisDeltaR"}; - AxisSpec axisNsigmaTOF{nSigmaTOFNBins, nSigmaTOFMinPlot, nSigmaTOFMaxPlot, "n#sigma_{TOF}", "axisNsigmaTOF"}; - AxisSpec axisNsigmaTPC{nSigmaTPCNBins, nSigmaTPCMinPlot, nSigmaTPCMaxPlot, "n#sigma_{TPC}", "axisNsigmaTPC"}; - AxisSpec axisCharge{chargeNBins, chargeMin, chargeMax, "Charge (ADC)", "axisCharge"}; - AxisSpec axisEta{etaNBins, etaMin, etaMax, "#eta", "axisEta"}; - AxisSpec axisPhi{phiNBins, phiMin, phiMax, "#phi", "axisPhi"}; - AxisSpec axisMass{massNBins, massMin, massMax, "mass^{2} (GeV/#it{c}^{2})", "axisMass"}; - AxisSpec axisChamber{10, -1.5, 8.5, "Chamber number", "axisChamber"}; - - // Deuteron Pos - registryDA.add("incomingDe_Pos_8cm", "incomingDe_Pos_8cm", HistType::kTH1F, {axisMom}); - registryDA.add("incomingDe_Pos_4cm", "incomingDe_Pos_4cm", HistType::kTH1F, {axisMom}); - - registryDA.add("De_Pos_deltaR_8cm", "De_Pos_deltaR_8cm", HistType::kTH1F, {axisDeltaR}); - registryDA.add("De_Pos_deltaR_4cm", "De_Pos_deltaR_4cm", HistType::kTH1F, {axisDeltaR}); - - registryDA.add("nSigmaTOF_vs_momHMPID_DePos_8cm", "nSigmaTOF_vs_momHMPID_DePos_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - registryDA.add("nSigmaTOF_vs_momHMPID_DePos_4cm", "nSigmaTOF_vs_momHMPID_DePos_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - - registryDA.add("survivingDe_Pos_8cm", "survivingDe_Pos_8cm", HistType::kTH2F, {axisMom, axisDeltaR}); - registryDA.add("survivingDe_Pos_4cm", "survivingDe_Pos_4cm", HistType::kTH2F, {axisMom, axisDeltaR}); - - registryDA.add("De_Pos_Q_8cm", "De_Pos_Q_8cm", HistType::kTH2F, {axisMom, axisCharge}); - registryDA.add("De_Pos_Q_4cm", "De_Pos_Q_4cm", HistType::kTH2F, {axisMom, axisCharge}); - - // nsigma plots - registryDA.add("nSigmaTPC_vs_momHMPID_Cut_DePos", "nSigmaTPC_vs_momHMPID_Cut_DePos", HistType::kTH2F, {axisMom, axisNsigmaTPC}); - registryDA.add("nSigmaTOF_vs_momHMPID_Cut_DePos", "nSigmaTOF_vs_momHMPID_Cut_DePos", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - - // check residuals distributions - registryDA.add("residualsX_Vs_momentum_Neg_8cm", "residualsX_Vs_momentum_Neg_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); - registryDA.add("residualsY_Vs_momentum_Neg_8cm", "residualsY_Vs_momentum_Neg_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); - registryDA.add("residualsX_Vs_momentum_Pos_8cm", "residualsX_Vs_momentum_Pos_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); - registryDA.add("residualsY_Vs_momentum_Pos_8cm", "residualsY_Vs_momentum_Pos_8cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); - - // Deuteron Neg - registryDA.add("incomingDe_Neg_8cm", "incomingDe_Neg_8cm", HistType::kTH1F, {axisMom}); - registryDA.add("incomingDe_Neg_4cm", "incomingDe_Neg_4cm", HistType::kTH1F, {axisMom}); - - registryDA.add("De_Neg_deltaR_8cm", "De_Neg_deltaR_8cm", HistType::kTH1F, {axisDeltaR}); - registryDA.add("De_Neg_deltaR_4cm", "De_Neg_deltaR_4cm", HistType::kTH1F, {axisDeltaR}); - - registryDA.add("nSigmaTOF_vs_momHMPID_DeNeg_8cm", "nSigmaTOF_vs_momHMPID_DeNeg_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - registryDA.add("nSigmaTOF_vs_momHMPID_DeNeg_4cm", "nSigmaTOF_vs_momHMPID_DeNeg_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - - registryDA.add("survivingDe_Neg_8cm", "survivingDe_Neg_8cm", HistType::kTH2F, {axisMom, axisDeltaR}); - registryDA.add("survivingDe_Neg_4cm", "survivingDe_Neg_4cm", HistType::kTH2F, {axisMom, axisDeltaR}); - registryDA.add("De_Neg_Q_8cm", "De_Neg_Q_8cm", HistType::kTH2F, {axisMom, axisCharge}); - registryDA.add("De_Neg_Q_4cm", "De_Neg_Q_4cm", HistType::kTH2F, {axisMom, axisCharge}); - registryDA.add("De_Neg_momentum", "De_Neg_momentum", HistType::kTH2F, {{100, 0.0, 5.0, "#it{p}_{vtx} (GeV/#it{c})"}, {100, 0.0, 5.0, "#it{p}_{hmpid} (GeV/#it{c})"}}); - - // check residuals distributions - registryDA.add("residualsX_Vs_momentum_Neg_4cm", "residualsX_Vs_momentum_Neg_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); - registryDA.add("residualsY_Vs_momentum_Neg_4cm", "residualsY_Vs_momentum_Neg_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); - registryDA.add("residualsX_Vs_momentum_Pos_4cm", "residualsX_Vs_momentum_Pos_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "x_{mip} - x_{track} (cm)"}}); - registryDA.add("residualsY_Vs_momentum_Pos_4cm", "residualsY_Vs_momentum_Pos_4cm", HistType::kTH2F, {axisMom, {400, -20.0f, 20.0f, "y_{mip} - y_{track} (cm)"}}); - - // nsigma plots - registryDA.add("nSigmaTPC_vs_momHMPID_Cut_DeNeg", "nSigmaTPC_vs_momHMPID_Cut_DeNeg", HistType::kTH2F, {axisMom, axisNsigmaTPC}); - registryDA.add("nSigmaTOF_vs_momHMPID_Cut_DeNeg", "nSigmaTOF_vs_momHMPID_Cut_DeNeg", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - - registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_DePos", "nSigmaTOF_vs_momHMPID_TPCOnly_DePos", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg", "nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - - registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_DePos_4cm", "nSigmaTOF_vs_momHMPID_TPCOnly_DePos_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_DePos_8cm", "nSigmaTOF_vs_momHMPID_TPCOnly_DePos_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg_4cm", "nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg_4cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - registryDA.add("nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg_8cm", "nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg_8cm", HistType::kTH2F, {axisMom, axisNsigmaTOF}); - - registryDA.add("pidTPC_TOF_p_ch_DePos", "pidTPC_TOF_p_ch_DePos", HistType::kTHnSparseF, {axisMom, axisNsigmaTPC, axisNsigmaTOF, axisChamber}); - registryDA.add("pidTPC_TOF_p_ch_DeNeg", "pidTPC_TOF_p_ch_DeNeg", HistType::kTHnSparseF, {axisMom, axisNsigmaTPC, axisNsigmaTOF, axisChamber}); - - // general plots - registryDA.add("hEta", "hEta", kTH1F, {axisEta}); - registryDA.add("hPhi", "hPhi", kTH1F, {axisPhi}); - registryDA.add("hMomentumHmpid", "hMomentumHmpid", kTH1F, {{1000, -10.f, 10.f, "#it{p}_{hmpid} (GeV/#it{c})"}}); - - registryDA.add("hMass", "hMass", kTH1F, {axisMass}); - registryDA.add("hMass_postDeuteron", "hMass_postDeuteron", kTH1F, {axisMass}); - registryDA.add("hMass_Vs_centrality", "hMass_Vs_centrality", kTH2F, {axisMass, {100, 0.0f, 100.0f, "centrality (%)"}}); - - for (int iCh = 0; iCh < Nch; iCh++) { - registryDA.add(Form("De_Pos_deltaR_%d", iCh), Form("De_Pos_deltaR_%d", iCh), HistType::kTH1F, {axisDeltaR}); - registryDA.add(Form("De_Neg_deltaR_%d", iCh), Form("De_Neg_deltaR_%d", iCh), HistType::kTH1F, {axisDeltaR}); - registryDA.add(Form("hEta_%d", iCh), Form("hEta_%d", iCh), HistType::kTH1F, {axisEta}); - registryDA.add(Form("hPhi_%d", iCh), Form("hPhi_%d", iCh), HistType::kTH1F, {axisPhi}); - - // tracking efficiency hmpid - registryDA.add(Form("hAllTracks_%d", iCh), Form("hAllTracks_%d", iCh), HistType::kTH1F, {axisMom}); - registryDA.add(Form("hCutQandD_%d", iCh), Form("hCutQandD_%d", iCh), HistType::kTH1F, {axisMom}); - registryDA.add(Form("hAllTracks_Pos_%d", iCh), Form("hAllTracks_Pos_%d", iCh), HistType::kTH1F, {axisMom}); - registryDA.add(Form("hAllTracks_Neg_%d", iCh), Form("hAllTracks_Neg_%d", iCh), HistType::kTH1F, {axisMom}); - registryDA.add(Form("hCutQandD_Pos_%d", iCh), Form("hCutQandD_Pos_%d", iCh), HistType::kTH1F, {axisMom}); - registryDA.add(Form("hCutQandD_Neg_%d", iCh), Form("hCutQandD_Neg_%d", iCh), HistType::kTH1F, {axisMom}); - } - } - - void process(aod::HmpidAnalysis const& hmpidTable) - { - for (const auto& hmpid : hmpidTable) { - - // ------------------------- - // track filters - // ------------------------- - if (hmpid.itsNCluster() < minReqClusterITS) { - continue; - } - if (hmpid.tpcNCluster() < minTPCnClsFound) { - continue; - } - if (hmpid.tpcNClsCrossedRows() < minNCrossedRowsTPC) { - continue; - } - if (hmpid.tpcChi2() > maxChi2TPC) { - continue; - } - if (hmpid.itsChi2() > maxChi2ITS) { - continue; - } - if (std::abs(hmpid.dcaXY()) > maxDCAxy) { - continue; - } - if (std::abs(hmpid.dcaZ()) > maxDCAz) { - continue; - } - - // ------------------------- - // derived quantities - // ------------------------- - const float momHmpid = hmpid.momentumHmpid(); - const float momAbs = std::abs(momHmpid); - - const float dx = hmpid.xMip() - hmpid.xTrack(); - const float dy = hmpid.yMip() - hmpid.yTrack(); - const float dr = std::hypot(dx - offsetX, dy - offsetY); - const float dMipTrack = std::hypot(dx, dy); - - const float mass = - std::pow(Nmean * momAbs * std::cos(hmpid.chAngle()), 2) - - std::pow(momAbs, 2); - - const int chamber = static_cast(hmpid.chamber()); - const bool isPos = momHmpid > 0; - const bool isNeg = momHmpid < 0; - - const bool passQandD = hmpid.chargeMip() > qMipThreshold && dMipTrack < cutDistanceMipTrack; - - // ------------------------- - // global distributions - // ------------------------- - registryDA.fill(HIST("hEta"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi"), hmpid.phiTrack()); - registryDA.fill(HIST("hMomentumHmpid"), momHmpid); - - if (hmpid.chAngle() > 0) { - registryDA.fill(HIST("hMass"), mass); - registryDA.fill(HIST("hMass_Vs_centrality"), mass, hmpid.centrality()); - } - - // ------------------------- - // Fill plots per chamber - // ------------------------- - - switch (chamber) { - - case Rich0: - registryDA.fill(HIST("hEta_0"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_0"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_0"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_0"), momAbs); - } - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_0"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_0"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_0"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_0"), momAbs); - } - } - break; - - case Rich1: - registryDA.fill(HIST("hEta_1"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_1"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_1"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_1"), momAbs); - } - - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_1"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_1"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_1"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_1"), momAbs); - } - } - break; - - case Rich2: - registryDA.fill(HIST("hEta_2"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_2"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_2"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_2"), momAbs); - } - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_2"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_2"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_2"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_2"), momAbs); - } - } - break; - - case Rich3: - registryDA.fill(HIST("hEta_3"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_3"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_3"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_3"), momAbs); - } - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_3"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_3"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_3"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_3"), momAbs); - } - } - break; - - case Rich4: - registryDA.fill(HIST("hEta_4"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_4"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_4"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_4"), momAbs); - } - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_4"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_4"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_4"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_4"), momAbs); - } - } - break; - - case Rich5: - registryDA.fill(HIST("hEta_5"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_5"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_5"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_5"), momAbs); - } - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_5"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_5"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_5"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_5"), momAbs); - } - } - break; - - case Rich6: - registryDA.fill(HIST("hEta_6"), hmpid.etaTrack()); - registryDA.fill(HIST("hPhi_6"), hmpid.phiTrack()); - registryDA.fill(HIST("hAllTracks_6"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_6"), momAbs); - } - if (isPos) { - registryDA.fill(HIST("hAllTracks_Pos_6"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Pos_6"), momAbs); - } - } - if (isNeg) { - registryDA.fill(HIST("hAllTracks_Neg_6"), momAbs); - if (passQandD) { - registryDA.fill(HIST("hCutQandD_Neg_6"), momAbs); - } - } - break; - - default: - break; - } - - // ------------------------- - // TPC + TOF PID distributions - // ------------------------- - if (isPos) { - registryDA.fill(HIST("pidTPC_TOF_p_ch_DePos"), momAbs, hmpid.tpcNSigmaDe(), hmpid.tofNSigmaDe(), chamber); - } - if (isNeg) { - registryDA.fill(HIST("pidTPC_TOF_p_ch_DeNeg"), momAbs, hmpid.tpcNSigmaDe(), hmpid.tofNSigmaDe(), chamber); - } - - // ------------------------- - // deuteron candidate cuts - TPC - // ------------------------- - if (hmpid.tpcNSigmaDe() < nsigmaTPCMin || - hmpid.tpcNSigmaDe() > nsigmaTPCMax) { - continue; - } - - // ------------------------- - // post TPC cut - // ------------------------- - if (isPos) { - registryDA.fill(HIST("nSigmaTPC_vs_momHMPID_Cut_DePos"), momAbs, hmpid.tpcNSigmaDe()); - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_DePos"), momAbs, hmpid.tofNSigmaDe()); - } - if (isNeg) { - registryDA.fill(HIST("nSigmaTPC_vs_momHMPID_Cut_DeNeg"), momAbs, hmpid.tpcNSigmaDe()); - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg"), momAbs, hmpid.tofNSigmaDe()); - } - - if (isPos && chamber == Rich2) { - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_DePos_4cm"), momAbs, hmpid.tofNSigmaDe()); - } - if (isPos && chamber == Rich4) { - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_DePos_8cm"), momAbs, hmpid.tofNSigmaDe()); - } - if (isNeg && chamber == Rich2) { - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg_4cm"), momAbs, hmpid.tofNSigmaDe()); - } - if (isNeg && chamber == Rich4) { - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_TPCOnly_DeNeg_8cm"), momAbs, hmpid.tofNSigmaDe()); - } - - registryDA.fill(HIST("hMass_postDeuteron"), mass); - - // ------------------------- - // deuteron candidate cuts - TOF - // ------------------------- - if (hmpid.tofNSigmaDe() < nsigmaTOFMin || - hmpid.tofNSigmaDe() > nsigmaTOFMax) { - continue; - } - - // ------------------------- - // TOF nsigma - // ------------------------- - if (isPos) { - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_Cut_DePos"), - momAbs, hmpid.tofNSigmaDe()); - } - - if (isNeg) { - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_Cut_DeNeg"), - momAbs, hmpid.tofNSigmaDe()); - } - - // ------------------------- - // DEUTERON CANDIDATES AFTER ALL CUTS - // ------------------------- - - // Pb-Pb collisions centrality requirement - if (collisionsPbPb) { - if (hmpid.centrality() < centralityCut) { - continue; - } - } - - // ------------------------- - // deltaR per chamber (Pos/Neg) - // ------------------------- - if (chamber == Rich0) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_0"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_0"), dr); - } - } - if (chamber == Rich1) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_1"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_1"), dr); - } - } - if (chamber == Rich2) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_2"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_2"), dr); - } - } - if (chamber == Rich3) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_3"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_3"), dr); - } - } - if (chamber == Rich4) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_4"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_4"), dr); - } - } - if (chamber == Rich5) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_5"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_5"), dr); - } - } - if (chamber == Rich6) { - if (isPos) { - registryDA.fill(HIST("De_Pos_deltaR_6"), dr); - } - if (isNeg) { - registryDA.fill(HIST("De_Neg_deltaR_6"), dr); - } - } - - // ------------------------- - // absorbers - // ------------------------- - const int abs4cm = 2; - const int abs8cm = 4; - - if (isPos && chamber == abs8cm) { - registryDA.fill(HIST("incomingDe_Pos_8cm"), momAbs); - registryDA.fill(HIST("survivingDe_Pos_8cm"), momAbs, dr); - registryDA.fill(HIST("De_Pos_Q_8cm"), momAbs, hmpid.chargeMip()); - registryDA.fill(HIST("De_Pos_deltaR_8cm"), dr); - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_DePos_8cm"), momAbs, hmpid.tofNSigmaDe()); - - registryDA.fill(HIST("residualsX_Vs_momentum_Pos_8cm"), momAbs, dx); - registryDA.fill(HIST("residualsY_Vs_momentum_Pos_8cm"), momAbs, dy); - } - - if (isNeg && chamber == abs8cm) { - registryDA.fill(HIST("incomingDe_Neg_8cm"), momAbs); - registryDA.fill(HIST("survivingDe_Neg_8cm"), momAbs, dr); - registryDA.fill(HIST("De_Neg_Q_8cm"), momAbs, hmpid.chargeMip()); - registryDA.fill(HIST("De_Neg_deltaR_8cm"), dr); - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_DeNeg_8cm"), momAbs, hmpid.tofNSigmaDe()); - - registryDA.fill(HIST("residualsX_Vs_momentum_Neg_8cm"), momAbs, dx); - registryDA.fill(HIST("residualsY_Vs_momentum_Neg_8cm"), momAbs, dy); - } - - if (isPos && chamber == abs4cm) { - registryDA.fill(HIST("incomingDe_Pos_4cm"), momAbs); - registryDA.fill(HIST("survivingDe_Pos_4cm"), momAbs, dr); - registryDA.fill(HIST("De_Pos_Q_4cm"), momAbs, hmpid.chargeMip()); - registryDA.fill(HIST("De_Pos_deltaR_4cm"), dr); - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_DePos_4cm"), momAbs, hmpid.tofNSigmaDe()); - - registryDA.fill(HIST("residualsX_Vs_momentum_Pos_4cm"), momAbs, dx); - registryDA.fill(HIST("residualsY_Vs_momentum_Pos_4cm"), momAbs, dy); - } - - if (isNeg && chamber == abs4cm) { - registryDA.fill(HIST("incomingDe_Neg_4cm"), momAbs); - registryDA.fill(HIST("survivingDe_Neg_4cm"), momAbs, dr); - registryDA.fill(HIST("De_Neg_Q_4cm"), momAbs, hmpid.chargeMip()); - registryDA.fill(HIST("De_Neg_deltaR_4cm"), dr); - registryDA.fill(HIST("nSigmaTOF_vs_momHMPID_DeNeg_4cm"), momAbs, hmpid.tofNSigmaDe()); - - registryDA.fill(HIST("residualsX_Vs_momentum_Neg_4cm"), momAbs, dx); - registryDA.fill(HIST("residualsY_Vs_momentum_Neg_4cm"), momAbs, dy); - } - } // end of loop over hmpidTable - } // end of process -}; - -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) -{ - return WorkflowSpec{ - adaptAnalysisTask(cfgc), - }; -}