Skip to content

Commit 5288bd0

Browse files
njacaziofmazzascFrancesco Mazzaschi
authored
[Common] [PID] add ITS PID (#8419)
Co-authored-by: Francesco Mazzaschi <43742195+fmazzasc@users.noreply.github.com> Co-authored-by: Francesco Mazzaschi <fmazzasc@alipap1.cern.ch>
1 parent 5bc86da commit 5288bd0

File tree

3 files changed

+271
-1
lines changed

3 files changed

+271
-1
lines changed

Common/DataModel/PIDResponseITS.h

Lines changed: 155 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,155 @@
1+
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2+
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3+
// All rights not expressly granted are reserved.
4+
//
5+
// This software is distributed under the terms of the GNU General Public
6+
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7+
//
8+
// In applying this license CERN does not waive the privileges and immunities
9+
// granted to it by virtue of its status as an Intergovernmental Organization
10+
// or submit itself to any jurisdiction.
11+
12+
///
13+
/// \file PIDResponseITS.h
14+
/// \since 2024-11-12
15+
/// \author Nicolò Jacazio nicolo.jacazio@cern.ch
16+
/// \author Francesco Mazzaschi francesco.mazzaschi@cern.ch
17+
/// \brief Set of tables, tasks and utilities to provide the interface between
18+
/// the analysis data model and the PID response of the ITS
19+
///
20+
21+
#ifndef COMMON_DATAMODEL_PIDRESPONSEITS_H_
22+
#define COMMON_DATAMODEL_PIDRESPONSEITS_H_
23+
24+
// O2 includes
25+
#include "Framework/ASoA.h"
26+
#include "Framework/AnalysisDataModel.h"
27+
#include "ReconstructionDataFormats/PID.h"
28+
#include "Framework/Logger.h"
29+
30+
namespace o2::aod
31+
{
32+
33+
struct ITSResponse {
34+
static float averageClusterSize(uint32_t itsClusterSizes)
35+
{
36+
float average = 0;
37+
int nclusters = 0;
38+
39+
for (int layer = 0; layer < 7; layer++) {
40+
if ((itsClusterSizes >> (layer * 4)) & 0xf) {
41+
nclusters++;
42+
average += (itsClusterSizes >> (layer * 4)) & 0xf;
43+
}
44+
}
45+
if (nclusters == 0) {
46+
return 0;
47+
}
48+
return average / nclusters;
49+
};
50+
51+
template <o2::track::PID::ID id>
52+
static float expSignal(const float momentum)
53+
{
54+
static constexpr float inverseMass = 1. / o2::track::pid_constants::sMasses[id];
55+
static constexpr float charge = static_cast<float>(o2::track::pid_constants::sCharges[id]);
56+
const float bg = momentum * inverseMass;
57+
return (mITSRespParams[0] / (std::pow(bg, mITSRespParams[1])) + mITSRespParams[2]) * std::pow(charge, mChargeFactor);
58+
}
59+
60+
template <o2::track::PID::ID id>
61+
static float nSigmaITS(uint32_t itsClusterSizes, float momentum)
62+
{
63+
const float exp = expSignal<id>(momentum);
64+
const float average = averageClusterSize(itsClusterSizes);
65+
const float resolution = mResolution * exp;
66+
return (average - exp) / resolution;
67+
};
68+
69+
static void setParameters(float p0, float p1, float p2, float chargeFactor, float resolution)
70+
{
71+
if (mIsInitialized) {
72+
LOG(fatal) << "ITSResponse parameters already initialized";
73+
}
74+
mIsInitialized = true;
75+
mITSRespParams[0] = p0;
76+
mITSRespParams[1] = p1;
77+
mITSRespParams[2] = p2;
78+
mChargeFactor = chargeFactor;
79+
mResolution = resolution;
80+
}
81+
82+
private:
83+
static std::array<float, 3> mITSRespParams;
84+
static float mChargeFactor;
85+
static float mResolution;
86+
static bool mIsInitialized;
87+
};
88+
89+
std::array<float, 3> ITSResponse::mITSRespParams = {0.903, 2.014, 2.440};
90+
float ITSResponse::mChargeFactor = 2.299999952316284f;
91+
float ITSResponse::mResolution = 0.15f;
92+
bool ITSResponse::mIsInitialized = false;
93+
94+
namespace pidits
95+
{
96+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaElImp, itsNSigmaEl, //! Nsigma separation with the ITS detector for electrons
97+
[](uint32_t itsClusterSizes, float momentum) -> float {
98+
return ITSResponse::nSigmaITS<o2::track::PID::Electron>(itsClusterSizes, momentum);
99+
});
100+
101+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaMuImp, itsNSigmaMu, //! Nsigma separation with the ITS detector for muons
102+
[](uint32_t itsClusterSizes, float momentum) -> float {
103+
return ITSResponse::nSigmaITS<o2::track::PID::Muon>(itsClusterSizes, momentum);
104+
});
105+
106+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaPiImp, itsNSigmaPi, //! Nsigma separation with the ITS detector for pions
107+
[](uint32_t itsClusterSizes, float momentum) -> float {
108+
return ITSResponse::nSigmaITS<o2::track::PID::Pion>(itsClusterSizes, momentum);
109+
});
110+
111+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaKaImp, itsNSigmaKa, //! Nsigma separation with the ITS detector for kaons
112+
[](uint32_t itsClusterSizes, float momentum) -> float {
113+
return ITSResponse::nSigmaITS<o2::track::PID::Kaon>(itsClusterSizes, momentum);
114+
});
115+
116+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaPrImp, itsNSigmaPr, //! Nsigma separation with the ITS detector for protons
117+
[](uint32_t itsClusterSizes, float momentum) -> float {
118+
return ITSResponse::nSigmaITS<o2::track::PID::Proton>(itsClusterSizes, momentum);
119+
});
120+
121+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaDeImp, itsNSigmaDe, //! Nsigma separation with the ITS detector for deuterons
122+
[](uint32_t itsClusterSizes, float momentum) -> float {
123+
return ITSResponse::nSigmaITS<o2::track::PID::Deuteron>(itsClusterSizes, momentum);
124+
});
125+
126+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaTrImp, itsNSigmaTr, //! Nsigma separation with the ITS detector for tritons
127+
[](uint32_t itsClusterSizes, float momentum) -> float {
128+
return ITSResponse::nSigmaITS<o2::track::PID::Triton>(itsClusterSizes, momentum);
129+
});
130+
131+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaHeImp, itsNSigmaHe, //! Nsigma separation with the ITS detector for helium3
132+
[](uint32_t itsClusterSizes, float momentum) -> float {
133+
return ITSResponse::nSigmaITS<o2::track::PID::Helium3>(itsClusterSizes, momentum);
134+
});
135+
136+
DECLARE_SOA_DYNAMIC_COLUMN(ITSNSigmaAlImp, itsNSigmaAl, //! Nsigma separation with the ITS detector for alphas
137+
[](uint32_t itsClusterSizes, float momentum) -> float {
138+
return ITSResponse::nSigmaITS<o2::track::PID::Alpha>(itsClusterSizes, momentum);
139+
});
140+
141+
// Define user friendly names for the columns to join with the tracks
142+
using ITSNSigmaEl = ITSNSigmaElImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
143+
using ITSNSigmaMu = ITSNSigmaMuImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
144+
using ITSNSigmaPi = ITSNSigmaPiImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
145+
using ITSNSigmaKa = ITSNSigmaKaImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
146+
using ITSNSigmaPr = ITSNSigmaPrImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
147+
using ITSNSigmaDe = ITSNSigmaDeImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
148+
using ITSNSigmaTr = ITSNSigmaTrImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
149+
using ITSNSigmaHe = ITSNSigmaHeImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
150+
using ITSNSigmaAl = ITSNSigmaAlImp<o2::aod::track::ITSClusterSizes, o2::aod::track::P>;
151+
152+
} // namespace pidits
153+
} // namespace o2::aod
154+
155+
#endif // COMMON_DATAMODEL_PIDRESPONSEITS_H_

Common/TableProducer/PID/CMakeLists.txt

Lines changed: 6 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -9,8 +9,13 @@
99
# granted to it by virtue of its status as an Intergovernmental Organization
1010
# or submit itself to any jurisdiction.
1111

12-
# TOF
12+
# ITS
13+
o2physics_add_dpl_workflow(pid-its
14+
SOURCES pidITS.cxx
15+
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore
16+
COMPONENT_NAME Analysis)
1317

18+
# TOF
1419
o2physics_add_dpl_workflow(pid-tof-base
1520
SOURCES pidTOFBase.cxx
1621
PUBLIC_LINK_LIBRARIES O2Physics::AnalysisCore O2::TOFBase
Lines changed: 110 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,110 @@
1+
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2+
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3+
// All rights not expressly granted are reserved.
4+
//
5+
// This software is distributed under the terms of the GNU General Public
6+
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7+
//
8+
// In applying this license CERN does not waive the privileges and immunities
9+
// granted to it by virtue of its status as an Intergovernmental Organization
10+
// or submit itself to any jurisdiction.
11+
12+
///
13+
/// \file pidITS.cxx
14+
/// \since 2024-11-12
15+
/// \author Nicolò Jacazio nicolo.jacazio@cern.ch
16+
/// \author Francesco Mazzaschi francesco.mazzaschi@cern.ch
17+
/// \brief Task to produce PID tables for ITS split for each particle.
18+
/// Only the tables for the mass hypotheses requested are filled, the others are sent empty.
19+
///
20+
21+
#include <utility>
22+
#include <vector>
23+
#include <string>
24+
25+
// O2 includes
26+
#include "Framework/runDataProcessing.h"
27+
#include "Framework/AnalysisTask.h"
28+
#include "Framework/HistogramRegistry.h"
29+
#include "ReconstructionDataFormats/Track.h"
30+
#include "CCDB/BasicCCDBManager.h"
31+
#include "TOFBase/EventTimeMaker.h"
32+
33+
// O2Physics includes
34+
#include "Common/DataModel/PIDResponseITS.h"
35+
#include "MetadataHelper.h"
36+
37+
using namespace o2;
38+
using namespace o2::framework;
39+
using namespace o2::framework::expressions;
40+
using namespace o2::track;
41+
42+
MetadataHelper metadataInfo;
43+
44+
static constexpr int nCases = 2;
45+
static constexpr int nParameters = 5;
46+
static const std::vector<std::string> casesNames{"Data", "MC"};
47+
static const std::vector<std::string> parameterNames{"bb1", "bb2", "bb3", "Charge exponent", "Resolution"};
48+
static constexpr float defaultParameters[nCases][nParameters]{{0.903, 2.014, 2.440, 2.299999952316284f, 0.15f},
49+
{0.903, 2.014, 2.440, 2.299999952316284f, 0.15f}};
50+
51+
/// Task to produce the ITS PID information for each particle species
52+
/// The parametrization is: [p0/(bg)**p1 + p2] * pow(q, p3), being bg = p/m and q the charge
53+
struct itsPid {
54+
55+
Configurable<LabeledArray<float>> itsParams{"itsParams",
56+
{defaultParameters[0], nCases, nParameters, casesNames, parameterNames},
57+
"Response parameters"};
58+
Configurable<bool> getFromCCDB{"getFromCCDB", false, "Get the parameters from CCDB"};
59+
60+
Service<o2::ccdb::BasicCCDBManager> ccdb;
61+
Configurable<std::string> paramfile{"param-file", "", "Path to the parametrization object, if empty the parametrization is not taken from file"};
62+
Configurable<std::string> url{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
63+
Configurable<std::string> ccdbPath{"ccdbPath", "Analysis/PID/TPC/Response", "Path of the TPC parametrization on the CCDB"};
64+
Configurable<std::string> recoPass{"recoPass", "", "Reconstruction pass name for CCDB query (automatically takes latest object for timestamp if blank)"};
65+
Configurable<int64_t> ccdbTimestamp{"ccdb-timestamp", 0, "timestamp of the object used to query in CCDB the detector response. Exceptions: -1 gets the latest object, 0 gets the run dependent timestamp"};
66+
67+
void init(o2::framework::InitContext&)
68+
{
69+
if (getFromCCDB) {
70+
ccdb->setURL(url.value);
71+
ccdb->setCaching(true);
72+
ccdb->setLocalObjectValidityChecking();
73+
ccdb->setCreatedNotAfter(std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count());
74+
LOG(fatal) << "Not implemented yet";
75+
} else {
76+
const char* key = metadataInfo.isMC() ? "MC" : "Data";
77+
o2::aod::ITSResponse::setParameters(itsParams->get(key, "bb1"),
78+
itsParams->get(key, "bb2"),
79+
itsParams->get(key, "bb3"),
80+
itsParams->get(key, "Charge exponent"),
81+
itsParams->get(key, "Resolution"));
82+
}
83+
}
84+
85+
/// Dummy process function for BCs, needed in case both Run2 and Run3 process functions are disabled
86+
void process(aod::Timestamps const&) {}
87+
88+
void processTest(o2::soa::Join<aod::TracksIU, aod::TracksExtra> const& tracks)
89+
{
90+
auto tracksWithPid = soa::Attach<o2::soa::Join<aod::TracksIU, aod::TracksExtra>,
91+
aod::pidits::ITSNSigmaEl, aod::pidits::ITSNSigmaMu, aod::pidits::ITSNSigmaPi,
92+
aod::pidits::ITSNSigmaKa, aod::pidits::ITSNSigmaPr, aod::pidits::ITSNSigmaDe,
93+
aod::pidits::ITSNSigmaTr, aod::pidits::ITSNSigmaHe, aod::pidits::ITSNSigmaAl>(tracks);
94+
95+
for (const auto& track : tracksWithPid) {
96+
LOG(info) << track.itsNSigmaEl();
97+
LOG(info) << track.itsNSigmaPi();
98+
LOG(info) << track.itsNSigmaPr();
99+
}
100+
}
101+
PROCESS_SWITCH(itsPid, processTest, "Produce a test", false);
102+
};
103+
104+
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
105+
{
106+
// Parse the metadata
107+
metadataInfo.initMetadata(cfgc);
108+
auto workflow = WorkflowSpec{adaptAnalysisTask<itsPid>(cfgc)};
109+
return workflow;
110+
}

0 commit comments

Comments
 (0)