Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
31 changes: 30 additions & 1 deletion ALICE3/Core/DelphesO2TrackSmearer.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,8 @@

#include "ALICE3/Core/DelphesO2TrackSmearer.h"

#include <Framework/Logger.h>

namespace o2
{
namespace delphes
Expand All @@ -45,11 +47,38 @@

bool TrackSmearer::loadTable(int pdg, const char* filename, bool forceReload)
{
auto ipdg = getIndexPDG(pdg);
if (!filename || filename[0] == '\0') {
LOG(info) << " --- No LUT file provided for PDG " << pdg << ". Skipping load.";
return false;
}
const auto ipdg = getIndexPDG(pdg);
LOGF(info, "Will load %s lut file ..: '%s'", getParticleName(pdg), filename);
if (mLUTHeader[ipdg] && !forceReload) {
std::cout << " --- LUT table for PDG " << pdg << " has been already loaded with index " << ipdg << std::endl;
return false;
}
if (strncmp(filename, "ccdb:", 5) == 0) { // Check if filename starts with "ccdb:"
LOG(info) << " --- LUT file source identified as CCDB.";
std::string path = std::string(filename).substr(5); // Remove "ccdb:" prefix
const std::string outPath = "/tmp/LUTs/";
filename = Form("%s/%s/snapshot.root", outPath.c_str(), path.c_str());
std::ifstream checkFile(filename); // Check if file already exists
if (!checkFile.is_open()) { // File does not exist, retrieve from CCDB
LOG(info) << " --- CCDB source detected for PDG " << pdg << ": " << path;
if (!mCcdbManager) {
LOG(fatal) << " --- CCDB manager not set. Please set it before loading LUT from CCDB.";
}
std::map<std::string, std::string> metadata;
mCcdbManager->getCCDBAccessor().retrieveBlob(path, outPath, metadata, 1);
// Add CCDB handling logic here if needed
LOG(info) << " --- Now retrieving LUT file from CCDB to: " << filename;
} else { // File exists, proceed to load
LOG(info) << " --- LUT file already exists: " << filename << ". Skipping download.";
checkFile.close();
}
return loadTable(pdg, filename, forceReload);
}

mLUTHeader[ipdg] = new lutHeader_t;

std::ifstream lutFile(filename, std::ifstream::binary);
Expand Down Expand Up @@ -139,7 +168,7 @@
}
}
} else {
float comparisonValue = mLUTHeader[ipdg]->nchmap.log ? log10(nch) : nch;

Check failure on line 171 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
if (mWhatEfficiency == 1) {
if (inch > 0 && comparisonValue < mLUTHeader[ipdg]->nchmap.max) {
interpolatedEff = (0.5f + fraction) * mLUTEntry[ipdg][inch][irad][ieta][ipt]->eff + (0.5f - fraction) * mLUTEntry[ipdg][inch - 1][irad][ieta][ipt]->eff;
Expand Down Expand Up @@ -192,7 +221,7 @@
double val = 0.;
for (int j = 0; j < 5; ++j)
val += lutEntry->eigvec[j][i] * o2track.getParam(j);
params_[i] = gRandom->Gaus(val, sqrt(lutEntry->eigval[i]));

Check failure on line 224 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
}
// transform back params vector
for (int i = 0; i < 5; ++i) {
Expand All @@ -202,7 +231,7 @@
o2track.setParam(val, i);
}
// should make a sanity check that par[2] sin(phi) is in [-1, 1]
if (fabs(o2track.getParam(2)) > 1.) {

Check failure on line 234 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
std::cout << " --- smearTrack failed sin(phi) sanity check: " << o2track.getParam(2) << std::endl;
}
// set covariance matrix
Expand All @@ -217,7 +246,7 @@
{

auto pt = o2track.getPt();
if (abs(pdg) == 1000020030) {

Check failure on line 249 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
pt *= 2.f;
}
auto eta = o2track.getEta();
Expand All @@ -234,7 +263,7 @@
{
float dummy = 0.0f;
auto lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
auto val = sqrt(lutEntry->covm[14]) * lutEntry->pt;

Check failure on line 266 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
return val;
}

Expand All @@ -244,8 +273,8 @@
{
float dummy = 0.0f;
auto lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
auto sigmatgl = sqrt(lutEntry->covm[9]); // sigmatgl2

Check failure on line 276 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
auto etaRes = fabs(sin(2.0 * atan(exp(-eta)))) * sigmatgl; // propagate tgl to eta uncertainty

Check failure on line 277 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
etaRes /= lutEntry->eta; // relative uncertainty
return etaRes;
}
Expand All @@ -255,7 +284,7 @@
{
float dummy = 0.0f;
auto lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
auto val = sqrt(lutEntry->covm[14]) * pow(lutEntry->pt, 2);

Check failure on line 287 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
return val;
}

Expand All @@ -265,8 +294,8 @@
{
float dummy = 0.0f;
auto lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
auto sigmatgl = sqrt(lutEntry->covm[9]); // sigmatgl2

Check failure on line 297 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
auto etaRes = fabs(sin(2.0 * atan(exp(-eta)))) * sigmatgl; // propagate tgl to eta uncertainty

Check failure on line 298 in ALICE3/Core/DelphesO2TrackSmearer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[std-prefix]

Use std:: prefix for names from the std namespace.
return etaRes;
}
/*****************************************************************/
Expand Down
47 changes: 38 additions & 9 deletions ALICE3/Core/DelphesO2TrackSmearer.h
Original file line number Diff line number Diff line change
Expand Up @@ -24,12 +24,14 @@
#ifndef ALICE3_CORE_DELPHESO2TRACKSMEARER_H_
#define ALICE3_CORE_DELPHESO2TRACKSMEARER_H_

#include <map>
#include <iostream>
#include <fstream>
#include <CCDB/BasicCCDBManager.h>
#include <ReconstructionDataFormats/Track.h>

#include "TRandom.h"
#include "ReconstructionDataFormats/Track.h"
#include <TRandom.h>

#include <fstream>
#include <iostream>
#include <map>

///////////////////////////////
/// DelphesO2/src/lutCovm.hh //
Expand Down Expand Up @@ -85,7 +87,7 @@ struct map_t {
if (bin > nbins - 1)
return nbins - 1;
return bin;
} //;
} //;
void print() { printf("nbins = %d, min = %f, max = %f, log = %s \n", nbins, min, max, log ? "on" : "off"); } //;
};

Expand Down Expand Up @@ -214,10 +216,34 @@ class TrackSmearer
return 7; // Helium3
default:
return 2; // Default: pion
} //;
} //;
}
}

void setdNdEta(float val) { mdNdEta = val; } //;
const char* getParticleName(int pdg)
{
switch (abs(pdg)) {
case 11:
return "electron";
case 13:
return "muon";
case 211:
return "pion";
case 321:
return "kaon";
case 2212:
return "proton";
case 1000010020:
return "deuteron";
case 1000010030:
return "triton";
case 1000020030:
return "helium3";
default:
return "pion"; // Default: pion
}
}
void setdNdEta(float val) { mdNdEta = val; } //;
void setCcdbManager(o2::ccdb::BasicCCDBManager* mgr) { mCcdbManager = mgr; } //;

protected:
static constexpr unsigned int nLUTs = 8; // Number of LUT available
Expand All @@ -228,6 +254,9 @@ class TrackSmearer
bool mSkipUnreconstructed = true; // don't smear tracks that are not reco'ed
int mWhatEfficiency = 1;
float mdNdEta = 1600.;

private:
o2::ccdb::BasicCCDBManager* mCcdbManager = nullptr;
};

} // namespace delphes
Expand Down
52 changes: 21 additions & 31 deletions ALICE3/TableProducer/OTF/onTheFlyTracker.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -239,42 +239,32 @@ struct OnTheFlyTracker {

// For TGenPhaseSpace seed
TRandom3 rand;
Service<o2::ccdb::BasicCCDBManager> ccdb;

void init(o2::framework::InitContext&)
{

ccdb->setURL("http://alice-ccdb.cern.ch");
ccdb->setTimestamp(-1);

if (enableLUT) {
std::map<int, const char*> mapPdgLut;
const char* lutElChar = lutEl->c_str();
const char* lutMuChar = lutMu->c_str();
const char* lutPiChar = lutPi->c_str();
const char* lutKaChar = lutKa->c_str();
const char* lutPrChar = lutPr->c_str();

LOGF(info, "Will load electron lut file ..: %s", lutElChar);
LOGF(info, "Will load muon lut file ......: %s", lutMuChar);
LOGF(info, "Will load pion lut file ......: %s", lutPiChar);
LOGF(info, "Will load kaon lut file ......: %s", lutKaChar);
LOGF(info, "Will load proton lut file ....: %s", lutPrChar);

mapPdgLut.insert(std::make_pair(11, lutElChar));
mapPdgLut.insert(std::make_pair(13, lutMuChar));
mapPdgLut.insert(std::make_pair(211, lutPiChar));
mapPdgLut.insert(std::make_pair(321, lutKaChar));
mapPdgLut.insert(std::make_pair(2212, lutPrChar));

if (enableNucleiSmearing) {
const char* lutDeChar = lutDe->c_str();
const char* lutTrChar = lutTr->c_str();
const char* lutHe3Char = lutHe3->c_str();
mapPdgLut.insert(std::make_pair(1000010020, lutDeChar));
mapPdgLut.insert(std::make_pair(1000010030, lutTrChar));
mapPdgLut.insert(std::make_pair(1000020030, lutHe3Char));
}
for (const auto& e : mapPdgLut) {
if (!mSmearer.loadTable(e.first, e.second)) {
LOG(fatal) << "Having issue with loading the LUT " << e.first << " " << e.second;
mSmearer.setCcdbManager(ccdb.operator->());

auto loadLUT = [&](int pdg, const std::string& lutFile) {
bool success = mSmearer.loadTable(pdg, lutFile.c_str());
if (!success && !lutFile.empty()) {
LOG(fatal) << "Having issue with loading the LUT " << pdg << " " << lutFile;
}
}
};
loadLUT(11, lutEl.value);
loadLUT(13, lutMu.value);
loadLUT(211, lutPi.value);
loadLUT(321, lutKa.value);
loadLUT(2212, lutPr.value);
loadLUT(1000010020, lutDe.value);
loadLUT(1000010030, lutTr.value);
loadLUT(1000020030, lutHe3.value);

// interpolate efficiencies if requested to do so
mSmearer.interpolateEfficiency(static_cast<bool>(interpolateLutEfficiencyVsNch));

Expand Down
Loading