Skip to content

[Common] setting up irSource automatically using the CollisionHelper #13282

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
ddobrigk merged 7 commits into from
Oct 28, 2025
Merged

[Common] setting up irSource automatically using the CollisionHelper #13282

Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
c7c8dc3
fixing errors
Oct 22, 2025
ebac66f
follow rules
amaringarcia Oct 22, 2025
0a147e5
follow rules
amaringarcia Oct 22, 2025
ceb9654
follow rules
amaringarcia Oct 22, 2025
41ea9ad
follow rules
amaringarcia Oct 22, 2025
1c0d40e
follow rules
amaringarcia Oct 22, 2025
f3a5a2c
follow rules
Oct 23, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
166 changes: 116 additions & 50 deletions Common/Tools/PID/pidTPCModule.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,12 +22,14 @@
#define COMMON_TOOLS_PID_PIDTPCMODULE_H_

#include "Common/CCDB/ctpRateFetcher.h"
#include "Common/Core/CollisionTypeHelper.h"
#include "Common/Core/PID/TPCPIDResponse.h"
#include "Common/Core/TableHelper.h"
#include "Common/DataModel/PIDResponseTPC.h"
#include "Common/TableProducer/PID/pidTPCBase.h"
#include "Tools/ML/model.h"

#include <DataFormatsParameters/GRPLHCIFData.h>
#include <Framework/AnalysisDataModel.h>
#include <Framework/AnalysisHelpers.h>
#include <Framework/Configurable.h>
Expand Down Expand Up @@ -82,10 +84,10 @@

struct pidTPCConfigurables : o2::framework::ConfigurableGroup {
std::string prefix = "pidTPC";
o2::framework::Configurable<std::string> paramfile{"param-file", "", "Path to the parametrization object, if empty the parametrization is not taken from file"};

Check failure on line 87 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<std::string> ccdbPath{"ccdbPath", "Analysis/PID/TPC/Response", "Path of the TPC parametrization on the CCDB"};
o2::framework::Configurable<std::string> recoPass{"recoPass", "", "Reconstruction pass name for CCDB query (automatically takes latest object for timestamp if blank)"};
o2::framework::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"};

Check failure on line 90 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
// Parameters for loading network from a file / downloading the file
o2::framework::Configurable<bool> useNetworkCorrection{"useNetworkCorrection", 0, "(bool) Wether or not to use the network correction for the TPC dE/dx signal"};
o2::framework::Configurable<bool> autofetchNetworks{"autofetchNetworks", 1, "(bool) Automatically fetches networks from CCDB for the correct run number"};
Expand All @@ -98,14 +100,14 @@
o2::framework::Configurable<int> savedEdxsCorrected{"savedEdxsCorrected", -1, {"Save table with corrected dE/dx calculated on the spot. 0: off, 1: on, -1: auto"}};
o2::framework::Configurable<bool> useCorrecteddEdx{"useCorrecteddEdx", false, "(bool) If true, use corrected dEdx value in Nsigma calculation instead of the one in the AO2D"};

o2::framework::Configurable<int> pidFullEl{"pid-full-el", -1, {"Produce PID information for the Electron mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 103 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullMu{"pid-full-mu", -1, {"Produce PID information for the Muon mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 104 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullPi{"pid-full-pi", -1, {"Produce PID information for the Pion mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 105 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullKa{"pid-full-ka", -1, {"Produce PID information for the Kaon mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 106 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullPr{"pid-full-pr", -1, {"Produce PID information for the Proton mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 107 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullDe{"pid-full-de", -1, {"Produce PID information for the Deuterons mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 108 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullTr{"pid-full-tr", -1, {"Produce PID information for the Triton mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 109 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullHe{"pid-full-he", -1, {"Produce PID information for the Helium3 mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};

Check failure on line 110 in Common/Tools/PID/pidTPCModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
o2::framework::Configurable<int> pidFullAl{"pid-full-al", -1, {"Produce PID information for the Alpha mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};
o2::framework::Configurable<int> pidTinyEl{"pid-tiny-el", -1, {"Produce PID information for the Electron mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};
o2::framework::Configurable<int> pidTinyMu{"pid-tiny-mu", -1, {"Produce PID information for the Muon mass hypothesis, overrides the automatic setup: the corresponding table can be set off (0) or on (1)"}};
Expand All @@ -127,13 +129,13 @@
o2::framework::Configurable<int> useNetworkHe{"useNetworkHe", 1, {"Switch for applying neural network on the helium3 mass hypothesis (if network enabled) (set to 0 to disable)"}};
o2::framework::Configurable<int> useNetworkAl{"useNetworkAl", 1, {"Switch for applying neural network on the alpha mass hypothesis (if network enabled) (set to 0 to disable)"}};
o2::framework::Configurable<float> networkBetaGammaCutoff{"networkBetaGammaCutoff", 0.45, {"Lower value of beta-gamma to override the NN application"}};
o2::framework::Configurable<std::string> irSource{"irSource", "ZNC hadronic", "Estimator of the interaction rate (Recommended: pp --> T0VTX, Pb-Pb --> ZNC hadronic)"};
o2::framework::Configurable<std::string> cfgPathGrpLhcIf{"ccdb-path-grplhcif", "GLO/Config/GRPLHCIF", "Path on the CCDB for the GRPLHCIF object"};
};

// helper getter - FIXME should be separate
int getPIDIndex(const int pdgCode) // Get O2 PID index corresponding to MC PDG code
{
switch (abs(pdgCode)) {
switch (std::abs(pdgCode)) {
case 11:
return o2::track::PID::Electron;
case 13:
Expand Down Expand Up @@ -214,6 +216,10 @@
std::vector<int> speciesNetworkFlags = std::vector<int>(9);
std::string networkVersion;

// To get automatically the proper Hadronic Rate
std::string irSource = "";
o2::common::core::CollisionSystemType::collType collsys = o2::common::core::CollisionSystemType::kCollSysUndef;

// Parametrization configuration
bool useCCDBParam = false;

Expand Down Expand Up @@ -323,6 +329,14 @@
}
LOG(info) << "Successfully retrieved TPC PID object from CCDB for timestamp " << time << ", period " << headers["LPMProductionTag"] << ", recoPass " << headers["RecoPassName"];
metadata["RecoPassName"] = headers["RecoPassName"]; // Force pass number for NN request to match retrieved BB
o2::parameters::GRPLHCIFData* grpo = ccdb->template getForTimeStamp<o2::parameters::GRPLHCIFData>(pidTPCopts.cfgPathGrpLhcIf.value, time);
LOG(info) << " collision type::" << CollisionSystemType::getCollisionTypeFromGrp(grpo);
collsys = CollisionSystemType::getCollisionTypeFromGrp(grpo);
if (collsys == CollisionSystemType::kCollSyspp) {
irSource = std::string("T0VTX");
} else {
irSource = std::string("ZNC hadronic");
}
response->PrintAll();
}
}
Expand Down Expand Up @@ -368,8 +382,8 @@
} // end init

//__________________________________________________
template <typename TCCDB, typename TCCDBApi, typename C, typename M, typename T, typename B>
std::vector<float> createNetworkPrediction(TCCDB& ccdb, TCCDBApi& ccdbApi, C const& collisions, M const& mults, T const& tracks, B const& bcs, const size_t size)
template <typename TCCDB, typename TCCDBApi, typename M, typename T, typename B>
std::vector<float> createNetworkPrediction(TCCDB& ccdb, TCCDBApi& ccdbApi, soa::Join<aod::Collisions, aod::EvSels> const& collisions, M const& mults, T const& tracks, B const& bcs, const size_t size)
{

std::vector<float> network_prediction;
Expand Down Expand Up @@ -397,6 +411,14 @@
}
LOG(info) << "Successfully retrieved TPC PID object from CCDB for timestamp " << bc.timestamp() << ", period " << headers["LPMProductionTag"] << ", recoPass " << headers["RecoPassName"];
metadata["RecoPassName"] = headers["RecoPassName"]; // Force pass number for NN request to match retrieved BB
o2::parameters::GRPLHCIFData* grpo = ccdb->template getForTimeStamp<o2::parameters::GRPLHCIFData>(pidTPCopts.cfgPathGrpLhcIf.value, bc.timestamp());
LOG(info) << "Collision type::" << CollisionSystemType::getCollisionTypeFromGrp(grpo);
collsys = CollisionSystemType::getCollisionTypeFromGrp(grpo);
if (collsys == CollisionSystemType::kCollSyspp) {
irSource = std::string("T0VTX");
} else {
irSource = std::string("ZNC hadronic");
}
response->PrintAll();
}

Expand Down Expand Up @@ -430,11 +452,26 @@
uint64_t counter_track_props = 0;
int loop_counter = 0;

// To load the Hadronic rate once for each collision
float hadronicRateBegin = 0.;
std::vector<float> hadronicRateForCollision(collisions.size(), 0.0f);
size_t i = 0;
for (const auto& collision : collisions) {
const auto& bc = collision.template bc_as<B>();
hadronicRateForCollision[i] = mRateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), irSource) * 1.e-3;
i++;
}
auto bc = bcs.begin();
hadronicRateBegin = mRateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), irSource) * 1.e-3; // kHz

// Filling a std::vector<float> to be evaluated by the network
// Evaluation on single tracks brings huge overhead: Thus evaluation is done on one large vector
for (int i = 0; i < 9; i++) { // Loop over particle number for which network correction is used
float hadronicRate = 0.;
uint64_t timeStamp_bcOld = 0;
static constexpr int NParticleTypes = 9;
constexpr int ExpectedInputDimensionsNNV2 = 7;
constexpr int ExpectedInputDimensionsNNV3 = 8;
constexpr auto NetworkVersionV2 = "2";
constexpr auto NetworkVersionV3 = "3";
for (int i = 0; i < NParticleTypes; i++) { // Loop over particle number for which network correction is used
for (auto const& trk : tracks) {
if (!trk.hasTPC()) {
continue;
Expand All @@ -450,20 +487,24 @@
track_properties[counter_track_props + 3] = o2::track::pid_constants::sMasses[i];
track_properties[counter_track_props + 4] = trk.has_collision() ? mults[trk.collisionId()] / 11000. : 1.;
track_properties[counter_track_props + 5] = std::sqrt(nNclNormalization / trk.tpcNClsFound());
if (input_dimensions == 7 && networkVersion == "2") {
if (input_dimensions == ExpectedInputDimensionsNNV2 && networkVersion == NetworkVersionV2) {
track_properties[counter_track_props + 6] = trk.has_collision() ? collisions.iteratorAt(trk.collisionId()).ft0cOccupancyInTimeRange() / 60000. : 1.;
}
if (input_dimensions == 8 && networkVersion == "3") {
if (input_dimensions == ExpectedInputDimensionsNNV3 && networkVersion == NetworkVersionV3) {
track_properties[counter_track_props + 6] = trk.has_collision() ? collisions.iteratorAt(trk.collisionId()).ft0cOccupancyInTimeRange() / 60000. : 1.;
if (trk.has_collision()) {
auto trk_bc = (collisions.iteratorAt(trk.collisionId())).template bc_as<B>();
if (trk_bc.timestamp() != timeStamp_bcOld) {
hadronicRate = mRateFetcher.fetch(ccdb.service, trk_bc.timestamp(), trk_bc.runNumber(), pidTPCopts.irSource.value) * 1.e-3;
if (collsys == CollisionSystemType::kCollSyspp) {
track_properties[counter_track_props + 7] = hadronicRateForCollision[trk.collisionId()] / 1500.;
} else {
track_properties[counter_track_props + 7] = hadronicRateForCollision[trk.collisionId()] / 50.;
}
timeStamp_bcOld = trk_bc.timestamp();
track_properties[counter_track_props + 7] = hadronicRate / 50.;
} else {
track_properties[counter_track_props + 7] = 1;
// asign Hadronic Rate at beginning of run if track does not belong to a collision
if (collsys == CollisionSystemType::kCollSyspp) {
track_properties[counter_track_props + 7] = hadronicRateBegin / 1500.;
} else {
track_properties[counter_track_props + 7] = hadronicRateBegin / 50.;
}
}
}
counter_track_props += input_dimensions;
Expand Down Expand Up @@ -493,7 +534,7 @@

//__________________________________________________
template <typename T, typename NSF, typename NST>
void makePidTables(const int flagFull, NSF& tableFull, const int flagTiny, NST& tableTiny, const o2::track::PID::ID pid, const float tpcSignal, const T& trk, const long multTPC, const std::vector<float>& network_prediction, const int& count_tracks, const int& tracksForNet_size)
void makePidTables(const int flagFull, NSF& tableFull, const int flagTiny, NST& tableTiny, const o2::track::PID::ID pid, const float tpcSignal, const T& trk, const int64_t multTPC, const std::vector<float>& network_prediction, const int& count_tracks, const int& tracksForNet_size)
{
if (flagFull != 1 && flagTiny != 1) {
return;
Expand Down Expand Up @@ -526,22 +567,24 @@

float nSigma = -999.f;
float bg = trk.tpcInnerParam() / o2::track::pid_constants::sMasses[pid]; // estimated beta-gamma for network cutoff
constexpr int NumOutputNodesSymmetricSigma = 2;
constexpr int NumOutputNodesAsymmetricSigma = 3;
if (pidTPCopts.useNetworkCorrection && speciesNetworkFlags[pid] && trk.has_collision() && bg > pidTPCopts.networkBetaGammaCutoff) {

// Here comes the application of the network. The output--dimensions of the network determine the application: 1: mean, 2: sigma, 3: sigma asymmetric
// For now only the option 2: sigma will be used. The other options are kept if there would be demand later on
if (network.getNumOutputNodes() == 1) { // Expected mean correction; no sigma correction
nSigma = (tpcSignal - network_prediction[count_tracks + tracksForNet_size * pid] * expSignal) / expSigma;
} else if (network.getNumOutputNodes() == 2) { // Symmetric sigma correction
expSigma = (network_prediction[2 * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[2 * (count_tracks + tracksForNet_size * pid)]) * expSignal;
nSigma = (tpcSignal / expSignal - network_prediction[2 * (count_tracks + tracksForNet_size * pid)]) / (network_prediction[2 * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[2 * (count_tracks + tracksForNet_size * pid)]);
} else if (network.getNumOutputNodes() == 3) { // Asymmetric sigma corection
if (tpcSignal / expSignal >= network_prediction[3 * (count_tracks + tracksForNet_size * pid)]) {
expSigma = (network_prediction[3 * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[3 * (count_tracks + tracksForNet_size * pid)]) * expSignal;
nSigma = (tpcSignal / expSignal - network_prediction[3 * (count_tracks + tracksForNet_size * pid)]) / (network_prediction[3 * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[3 * (count_tracks + tracksForNet_size * pid)]);
} else if (network.getNumOutputNodes() == NumOutputNodesSymmetricSigma) { // Symmetric sigma correction
expSigma = (network_prediction[NumOutputNodesSymmetricSigma * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[NumOutputNodesSymmetricSigma * (count_tracks + tracksForNet_size * pid)]) * expSignal;
nSigma = (tpcSignal / expSignal - network_prediction[NumOutputNodesSymmetricSigma * (count_tracks + tracksForNet_size * pid)]) / (network_prediction[NumOutputNodesSymmetricSigma * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[NumOutputNodesSymmetricSigma * (count_tracks + tracksForNet_size * pid)]);
} else if (network.getNumOutputNodes() == NumOutputNodesAsymmetricSigma) { // Asymmetric sigma corection
if (tpcSignal / expSignal >= network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)]) {
expSigma = (network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)]) * expSignal;
nSigma = (tpcSignal / expSignal - network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)]) / (network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid) + 1] - network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)]);
} else {
expSigma = (network_prediction[3 * (count_tracks + tracksForNet_size * pid)] - network_prediction[3 * (count_tracks + tracksForNet_size * pid) + 2]) * expSignal;
nSigma = (tpcSignal / expSignal - network_prediction[3 * (count_tracks + tracksForNet_size * pid)]) / (network_prediction[3 * (count_tracks + tracksForNet_size * pid)] - network_prediction[3 * (count_tracks + tracksForNet_size * pid) + 2]);
expSigma = (network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)] - network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid) + 2]) * expSignal;
nSigma = (tpcSignal / expSignal - network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)]) / (network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid)] - network_prediction[NumOutputNodesAsymmetricSigma * (count_tracks + tracksForNet_size * pid) + 2]);
}
} else {
LOGF(fatal, "Network output-dimensions incompatible!");
Expand All @@ -556,8 +599,8 @@
};

//__________________________________________________
template <typename TCCDB, typename TCCDBApi, typename TBCs, typename TCollisions, typename TTracks, typename TTracksQA, typename TProducts>
void process(TCCDB& ccdb, TCCDBApi& ccdbApi, TBCs const& bcs, TCollisions const& cols, TTracks const& tracks, TTracksQA const& tracksQA, TProducts& products)
template <typename TCCDB, typename TCCDBApi, typename TBCs, typename TTracks, typename TTracksQA, typename TProducts>
void process(TCCDB& ccdb, TCCDBApi& ccdbApi, TBCs const& bcs, soa::Join<aod::Collisions, aod::EvSels> const& cols, TTracks const& tracks, TTracksQA const& tracksQA, TProducts& products)
{
if (tracks.size() == 0) {
return; // empty protection
Expand All @@ -568,9 +611,9 @@
}

// preparatory step: we need the multiplicities for each collision
std::vector<int> pidmults;
long totalTPCtracks = 0;
long totalTPCnotStandalone = 0;
std::vector<int64_t> pidmults;
int64_t totalTPCtracks = 0;
int64_t totalTPCnotStandalone = 0;
pidmults.resize(cols.size(), 0);

// faster counting
Expand Down Expand Up @@ -627,19 +670,33 @@

//_______________________________________
// process tracksQA in case present
std::vector<int64_t> indexTrack2TrackQA;
std::vector<int64_t> indexTrack2TrackQA(outTable_size, -1);
if constexpr (soa::is_table<TTracksQA>) {
for (const auto& trackQA : tracksQA) {
indexTrack2TrackQA[trackQA.trackId()] = trackQA.globalIndex();
}
}
//_______________________________________

// Fill Hadronic rate per collision in case CorrectedDEdx is requested
std::vector<float> hadronicRateForCollision(cols.size(), 0.0f);
float hadronicRateBegin = 0.0f;
if (pidTPCopts.useCorrecteddEdx) {
size_t i = 0;
for (const auto& collision : cols) {
const auto& bc = collision.template bc_as<aod::BCsWithTimestamps>();
hadronicRateForCollision[i] = mRateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), irSource) * 1.e-3;
i++;
}
auto bc = bcs.begin();
hadronicRateBegin = mRateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), irSource) * 1.e-3; // kHz
}

for (auto const& trk : tracks) {
// get the TPC signal to be used in the PID
float tpcSignalToEvaluatePID = trk.tpcSignal();

int multTPC = 0;
int64_t multTPC = 0;
if (trk.has_collision()) {
multTPC = pidmults[trk.collisionId()];
}
Expand All @@ -662,24 +719,25 @@
int occupancy;
if (trk.has_collision()) {
auto collision = cols.iteratorAt(trk.collisionId());
auto bc = collision.template bc_as<aod::BCsWithTimestamps>();
const int runnumber = bc.runNumber();
hadronicRate = mRateFetcher.fetch(ccdb.service, bc.timestamp(), runnumber, "ZNC hadronic") * 1.e-3; // kHz
hadronicRate = hadronicRateForCollision[trk.collisionId()];
occupancy = collision.trackOccupancyInTimeRange();
} else {
auto bc = bcs.begin();
const int runnumber = bc.runNumber();
hadronicRate = mRateFetcher.fetch(ccdb.service, bc.timestamp(), runnumber, "ZNC hadronic") * 1.e-3; // kHz
hadronicRate = hadronicRateBegin;
occupancy = 0;
}

constexpr float kExpectedTPCSignalMIP = 50.0f;
constexpr float kMaxAllowedRatio = 1.05f;
constexpr float kMinAllowedRatio = 0.05f;
constexpr float kMaxAllowedOcc = 12.0f;

float fTPCSignal = tpcSignalToEvaluatePID;
float fNormMultTPC = multTPC / 11000.;

float fTrackOccN = occupancy / 1000.;
float fOccTPCN = fNormMultTPC * 10; //(fNormMultTPC*10).clip(0,12)
if (fOccTPCN > 12)
fOccTPCN = 12;
if (fOccTPCN > kMaxAllowedOcc)
fOccTPCN = kMaxAllowedOcc;
else if (fOccTPCN < 0)
fOccTPCN = 0;

Expand All @@ -688,11 +746,11 @@
float a1pt = std::abs(trk.signed1Pt());
float a1pt2 = a1pt * a1pt;
float atgl = std::abs(trk.tgl());
float mbb0R = 50 / fTPCSignal;
if (mbb0R > 1.05)
mbb0R = 1.05;
else if (mbb0R < 0.05)
mbb0R = 0.05;
float mbb0R = kExpectedTPCSignalMIP / fTPCSignal;
if (mbb0R > kMaxAllowedRatio)
mbb0R = kMaxAllowedRatio;
else if (mbb0R < kMinAllowedRatio)
mbb0R = kMinAllowedRatio;
// float mbb0R = max(0.05, min(50 / fTPCSignal, 1.05));
float a1ptmbb0R = a1pt * mbb0R;
float atglmbb0R = atgl * mbb0R;
Expand All @@ -702,11 +760,11 @@

float fTPCSignalN_CR0 = str_dedx_correction.fReal_fTPCSignalN(vec_occu, vec_track);

float mbb0R1 = 50 / (fTPCSignal / fTPCSignalN_CR0);
if (mbb0R1 > 1.05)
mbb0R1 = 1.05;
else if (mbb0R1 < 0.05)
mbb0R1 = 0.05;
float mbb0R1 = kExpectedTPCSignalMIP / (fTPCSignal / fTPCSignalN_CR0);
if (mbb0R1 > kMaxAllowedRatio)
mbb0R1 = kMaxAllowedRatio;
else if (mbb0R1 < kMinAllowedRatio)
mbb0R1 = kMinAllowedRatio;

std::vector<float> vec_track1 = {mbb0R1, a1pt, atgl, atgl * mbb0R1, a1pt * mbb0R1, side, a1pt2};
float fTPCSignalN_CR1 = str_dedx_correction.fReal_fTPCSignalN(vec_occu, vec_track1);
Expand Down Expand Up @@ -738,6 +796,14 @@
}
}
LOG(info) << "Successfully retrieved TPC PID object from CCDB for timestamp " << bc.timestamp() << ", period " << headers["LPMProductionTag"] << ", recoPass " << headers["RecoPassName"];
o2::parameters::GRPLHCIFData* grpo = ccdb->template getForTimeStamp<o2::parameters::GRPLHCIFData>(pidTPCopts.cfgPathGrpLhcIf.value, bc.timestamp());
LOG(info) << "Collisions type::" << CollisionSystemType::getCollisionTypeFromGrp(grpo);
collsys = CollisionSystemType::getCollisionTypeFromGrp(grpo);
if (collsys == CollisionSystemType::kCollSyspp) {
irSource = std::string("T0VTX");
} else {
irSource = std::string("ZNC hadronic");
}
response->PrintAll();
}

Expand Down
Loading