Skip to content

Commit 9fc286b

Browse files
njacazionjacazio
authored
Update DelphesO2TrackSmearer.cxx
1 parent b76bdaf commit 9fc286b

File tree

1 file changed

+77
-39
lines changed

1 file changed

+77
-39
lines changed

‎ALICE3/Core/DelphesO2TrackSmearer.cxx‎

Lines changed: 77 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -13,27 +13,13 @@
1313
/// \file DelphesO2TrackSmearer.cxx
1414
/// \author Roberto Preghenella
1515
/// \brief Porting to O2Physics of DelphesO2 code.
16-
/// Minimal changes have been made to the original code for adaptation purposes, formatting and commented parts have been considered.
1716
/// Relevant sources:
1817
/// DelphesO2/src/lutCovm.hh https://github.com/AliceO2Group/DelphesO2/blob/master/src/lutCovm.hh
1918
/// DelphesO2/src/TrackSmearer.cc https://github.com/AliceO2Group/DelphesO2/blob/master/src/TrackSmearer.cc
2019
/// DelphesO2/src/TrackSmearer.hh https://github.com/AliceO2Group/DelphesO2/blob/master/src/TrackSmearer.hh
2120
/// @email: preghenella@bo.infn.it
2221
///
2322

24-
//////////////////////////////
25-
// DelphesO2/src/lutCovm.cc //
26-
//////////////////////////////
27-
28-
/// @author: Roberto Preghenella
29-
/// @email: preghenella@bo.infn.it
30-
31-
// #include "TrackSmearer.hh"
32-
// #include "TrackUtils.hh"
33-
// #include "TRandom.h"
34-
// #include <iostream>
35-
// #include <fstream>
36-
3723
#include "ALICE3/Core/DelphesO2TrackSmearer.h"
3824

3925
#include "ALICE3/Core/GeometryContainer.h"
@@ -51,7 +37,59 @@ namespace delphes
5137

5238
/*****************************************************************/
5339

54-
bool TrackSmearer::loadTable(int pdg, const char* filename, bool forceReload)
40+
int DelphesO2TrackSmearer::getIndexPDG(const int pdg)
41+
{
42+
switch (abs(pdg)) {
43+
case 11:
44+
return 0; // Electron
45+
case 13:
46+
return 1; // Muon
47+
case 211:
48+
return 2; // Pion
49+
case 321:
50+
return 3; // Kaon
51+
case 2212:
52+
return 4; // Proton
53+
case 1000010020:
54+
return 5; // Deuteron
55+
case 1000010030:
56+
return 6; // Triton
57+
case 1000020030:
58+
return 7; // Helium3
59+
case 1000020040:
60+
return 8; // Alphas
61+
default:
62+
return 2; // Default: pion
63+
}
64+
}
65+
66+
const char* DelphesO2TrackSmearer::getParticleName(int pdg)
67+
{
68+
switch (abs(pdg)) {
69+
case 11:
70+
return "electron";
71+
case 13:
72+
return "muon";
73+
case 211:
74+
return "pion";
75+
case 321:
76+
return "kaon";
77+
case 2212:
78+
return "proton";
79+
case 1000010020:
80+
return "deuteron";
81+
case 1000010030:
82+
return "triton";
83+
case 1000020030:
84+
return "helium3";
85+
case 1000020040:
86+
return "alpha";
87+
default:
88+
return "pion"; // Default: pion
89+
}
90+
}
91+
92+
bool DelphesO2TrackSmearer::loadTable(int pdg, const char* filename, bool forceReload)
5593
{
5694
if (!filename || filename[0] == '\0') {
5795
LOG(info) << " --- No LUT file provided for PDG " << pdg << ". Skipping load.";
@@ -109,17 +147,17 @@ bool TrackSmearer::loadTable(int pdg, const char* filename, bool forceReload)
109147
const int nrad = mLUTHeader[ipdg]->radmap.nbins;
110148
const int neta = mLUTHeader[ipdg]->etamap.nbins;
111149
const int npt = mLUTHeader[ipdg]->ptmap.nbins;
112-
mLUTEntry[ipdg] = new lutEntry_t****[nnch];
150+
mLUTEntry[ipdg] = new DelphesO2TrackSmearer::lutEntry_t****[nnch];
113151
for (int inch = 0; inch < nnch; ++inch) {
114-
mLUTEntry[ipdg][inch] = new lutEntry_t***[nrad];
152+
mLUTEntry[ipdg][inch] = new DelphesO2TrackSmearer::lutEntry_t***[nrad];
115153
for (int irad = 0; irad < nrad; ++irad) {
116-
mLUTEntry[ipdg][inch][irad] = new lutEntry_t**[neta];
154+
mLUTEntry[ipdg][inch][irad] = new DelphesO2TrackSmearer::lutEntry_t**[neta];
117155
for (int ieta = 0; ieta < neta; ++ieta) {
118-
mLUTEntry[ipdg][inch][irad][ieta] = new lutEntry_t*[npt];
156+
mLUTEntry[ipdg][inch][irad][ieta] = new DelphesO2TrackSmearer::lutEntry_t*[npt];
119157
for (int ipt = 0; ipt < npt; ++ipt) {
120-
mLUTEntry[ipdg][inch][irad][ieta][ipt] = new lutEntry_t;
121-
lutFile.read(reinterpret_cast<char*>(mLUTEntry[ipdg][inch][irad][ieta][ipt]), sizeof(lutEntry_t));
122-
if (lutFile.gcount() != sizeof(lutEntry_t)) {
158+
mLUTEntry[ipdg][inch][irad][ieta][ipt] = new DelphesO2TrackSmearer::lutEntry_t;
159+
lutFile.read(reinterpret_cast<char*>(mLUTEntry[ipdg][inch][irad][ieta][ipt]), sizeof(DelphesO2TrackSmearer::lutEntry_t));
160+
if (lutFile.gcount() != sizeof(DelphesO2TrackSmearer::lutEntry_t)) {
123161
LOG(info) << " --- troubles reading covariance matrix entry for PDG " << pdg << ": " << localFilename << std::endl;
124162
LOG(info) << " --- expected/detected " << sizeof(lutHeader_t) << "/" << lutFile.gcount() << std::endl;
125163
return false;
@@ -137,7 +175,7 @@ bool TrackSmearer::loadTable(int pdg, const char* filename, bool forceReload)
137175

138176
/*****************************************************************/
139177

140-
lutEntry_t* TrackSmearer::getLUTEntry(const int pdg, const float nch, const float radius, const float eta, const float pt, float& interpolatedEff)
178+
DelphesO2TrackSmearer::lutEntry_t* DelphesO2TrackSmearer::getLUTEntry(const int pdg, const float nch, const float radius, const float eta, const float pt, float& interpolatedEff)
141179
{
142180
const int ipdg = getIndexPDG(pdg);
143181
if (!mLUTHeader[ipdg]) {
@@ -210,7 +248,7 @@ lutEntry_t* TrackSmearer::getLUTEntry(const int pdg, const float nch, const floa
210248

211249
/*****************************************************************/
212250

213-
bool TrackSmearer::smearTrack(O2Track& o2track, lutEntry_t* lutEntry, float interpolatedEff)
251+
bool DelphesO2TrackSmearer::smearTrack(o2::track::TrackParCov& o2track, DelphesO2TrackSmearer::lutEntry_t* lutEntry, float interpolatedEff)
214252
{
215253
bool isReconstructed = true;
216254
// generate efficiency
@@ -263,7 +301,7 @@ bool TrackSmearer::smearTrack(O2Track& o2track, lutEntry_t* lutEntry, float inte
263301

264302
/*****************************************************************/
265303

266-
bool TrackSmearer::smearTrack(O2Track& o2track, int pdg, float nch)
304+
bool DelphesO2TrackSmearer::smearTrack(o2::track::TrackParCov& o2track, int pdg, float nch)
267305
{
268306
auto pt = o2track.getPt();
269307
switch (pdg) {
@@ -274,67 +312,67 @@ bool TrackSmearer::smearTrack(O2Track& o2track, int pdg, float nch)
274312
}
275313
auto eta = o2track.getEta();
276314
float interpolatedEff = 0.0f;
277-
lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, interpolatedEff);
315+
DelphesO2TrackSmearer::lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, interpolatedEff);
278316
if (!lutEntry || !lutEntry->valid)
279317
return false;
280318
return smearTrack(o2track, lutEntry, interpolatedEff);
281319
}
282320

283321
/*****************************************************************/
284322
// relative uncertainty on pt
285-
double TrackSmearer::getPtRes(const int pdg, const float nch, const float eta, const float pt)
323+
double DelphesO2TrackSmearer::getPtRes(const int pdg, const float nch, const float eta, const float pt)
286324
{
287325
float dummy = 0.0f;
288-
lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
326+
DelphesO2TrackSmearer::lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
289327
auto val = std::sqrt(lutEntry->covm[14]) * lutEntry->pt;
290328
return val;
291329
}
292330

293331
/*****************************************************************/
294332
// relative uncertainty on eta
295-
double TrackSmearer::getEtaRes(const int pdg, const float nch, const float eta, const float pt)
333+
double DelphesO2TrackSmearer::getEtaRes(const int pdg, const float nch, const float eta, const float pt)
296334
{
297335
float dummy = 0.0f;
298-
lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
336+
DelphesO2TrackSmearer::lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
299337
auto sigmatgl = std::sqrt(lutEntry->covm[9]); // sigmatgl2
300338
auto etaRes = std::fabs(std::sin(2.0 * std::atan(std::exp(-eta)))) * sigmatgl; // propagate tgl to eta uncertainty
301339
etaRes /= lutEntry->eta; // relative uncertainty
302340
return etaRes;
303341
}
304342
/*****************************************************************/
305343
// absolute uncertainty on pt
306-
double TrackSmearer::getAbsPtRes(const int pdg, const float nch, const float eta, const float pt)
344+
double DelphesO2TrackSmearer::getAbsPtRes(const int pdg, const float nch, const float eta, const float pt)
307345
{
308346
float dummy = 0.0f;
309-
lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
347+
DelphesO2TrackSmearer::lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
310348
auto val = std::sqrt(lutEntry->covm[14]) * lutEntry->pt * lutEntry->pt;
311349
return val;
312350
}
313351

314352
/*****************************************************************/
315353
// absolute uncertainty on eta
316-
double TrackSmearer::getAbsEtaRes(const int pdg, const float nch, const float eta, const float pt)
354+
double DelphesO2TrackSmearer::getAbsEtaRes(const int pdg, const float nch, const float eta, const float pt)
317355
{
318356
float dummy = 0.0f;
319-
lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
357+
DelphesO2TrackSmearer::lutEntry_t* lutEntry = getLUTEntry(pdg, nch, 0., eta, pt, dummy);
320358
auto sigmatgl = std::sqrt(lutEntry->covm[9]); // sigmatgl2
321359
auto etaRes = std::fabs(std::sin(2.0 * std::atan(std::exp(-eta)))) * sigmatgl; // propagate tgl to eta uncertainty
322360
return etaRes;
323361
}
324362
/*****************************************************************/
325363
// efficiency
326-
double TrackSmearer::getEfficiency(const int pdg, const float nch, const float eta, const float pt)
364+
double DelphesO2TrackSmearer::getEfficiency(const int pdg, const float nch, const float eta, const float pt)
327365
{
328366
float efficiency = 0.0f;
329367
getLUTEntry(pdg, nch, 0., eta, pt, efficiency);
330368
return efficiency;
331369
}
332370
/*****************************************************************/
333371
// Only in DelphesO2
334-
// bool TrackSmearer::smearTrack(Track& track, bool atDCA)
372+
// bool DelphesO2TrackSmearer::smearTrack(Track& track, bool atDCA)
335373
// {
336374

337-
// O2Track o2track;
375+
// o2::track::TrackParCov o2track;
338376
// TrackUtils::convertTrackToO2Track(track, o2track, atDCA);
339377
// int pdg = track.PID;
340378
// float nch = mdNdEta; // use locally stored dNch/deta for the time being
@@ -344,11 +382,11 @@ double TrackSmearer::getEfficiency(const int pdg, const float nch, const float e
344382
// return true;
345383

346384
// #if 0
347-
// lutEntry_t* lutEntry = getLUTEntry(track.PID, 0., 0., track.Eta, track.PT);
385+
// DelphesO2TrackSmearer::lutEntry_t* lutEntry = getLUTEntry(track.PID, 0., 0., track.Eta, track.PT);
348386
// if (!lutEntry)
349387
// return;
350388

351-
// O2Track o2track;
389+
// o2::track::TrackParCov o2track;
352390
// TrackUtils::convertTrackToO2Track(track, o2track, atDCA);
353391
// smearTrack(o2track, lutEntry);
354392
// TrackUtils::convertO2TrackToTrack(o2track, track, atDCA);

0 commit comments

Comments
 (0)