1919// / @email: preghenella@bo.infn.it
2020// /
2121
22- #include < cstdio >
22+ #include " ALICE3/Core/DelphesO2LutWriter.h "
2323
2424#include " ALICE3/Core/DelphesO2TrackSmearer.h"
25- #include " ALICE3/Core/DelphesO2LutWriter.h"
26- #include < iostream>
27- #include " TMatrixD.h"
28- #include " TVectorD.h"
25+ #include " ALICE3/Core/FastTracker.h"
26+ #include " ALICE3/Core/TrackUtilities.h"
27+
2928#include " TAxis.h"
30- #include " TMatrixDSymEigen.h"
3129#include " TDatabasePDG.h"
3230#include " TLorentzVector.h"
33- #include " ALICE3/Core/FastTracker.h"
34- #include " ALICE3/Core/TrackUtilities.h"
31+ #include " TMatrixD.h"
32+ #include " TMatrixDSymEigen.h"
33+ #include " TVectorD.h"
34+
35+ #include < cstdio>
36+ #include < iostream>
37+ #include < string>
3538
3639// #define USE_FWD_PARAM
3740#ifdef USE_FWD_PARAM
4144namespace o2 ::fastsim
4245{
4346
44- void DelphesO2LutWriter::printLutWriterConfiguration ()
47+ void DelphesO2LutWriter::Print () const
4548{
46- std::cout << " --- Printing configuration of LUT writer --- " << std::endl;
47- std::cout << " -> etaMaxBarrel = " << etaMaxBarrel << std::endl;
48- std::cout << " -> usePara = " << usePara << std::endl;
49- std::cout << " -> useDipole = " << useDipole << std::endl;
50- std::cout << " -> useFlatDipole = " << useFlatDipole << std::endl;
49+ LOG (info) << " --- Printing configuration of LUT writer --- " ;
50+ LOG (info) << " -> etaMaxBarrel = " << etaMaxBarrel;
51+ LOG (info) << " -> usePara = " << usePara;
52+ LOG (info) << " -> useDipole = " << useDipole;
53+ LOG (info) << " -> useFlatDipole = " << useFlatDipole;
54+ LOG (info) << " -> mAtLeastHits = " << mAtLeastHits ;
55+ LOG (info) << " -> mAtLeastCorr = " << mAtLeastCorr ;
56+ LOG (info) << " -> mAtLeastFake = " << mAtLeastFake ;
57+ LOG (info) << " -> Nch Binning: = " << mNchBinning .toString ();
58+ LOG (info) << " -> Radius Binning: = " << mRadiusBinning .toString ();
59+ LOG (info) << " -> Eta Binning: = " << mEtaBinning .toString ();
60+ LOG (info) << " -> Pt Binning: = " << mPtBinning .toString ();
61+ LOG (info) << " --- End of configuration --- " ;
62+ }
63+
64+ std::string DelphesO2LutWriter::LutBinning::toString () const
65+ {
66+ std::string str = " " ;
67+ str.append (log ? " log" : " lin" );
68+ str.append (" nbins: " );
69+ str.append (std::to_string (nbins));
70+ str.append (" min: " );
71+ str.append (std::to_string (min));
72+ str.append (" max: " );
73+ str.append (std::to_string (max));
74+ return str;
5175}
5276
5377bool DelphesO2LutWriter::fatSolve (lutEntry_t& lutEntry,
5478 float pt,
5579 float eta,
5680 const float mass,
57- int itof,
58- int otof,
81+ size_t itof,
82+ size_t otof,
5983 int q,
6084 const float nch)
6185{
@@ -65,13 +89,19 @@ bool DelphesO2LutWriter::fatSolve(lutEntry_t& lutEntry,
6589 tlv.SetPtEtaPhiM (pt, eta, 0 ., mass);
6690 o2::track::TrackParCov trkIn;
6791 o2::upgrade::convertTLorentzVectorToO2Track (q, tlv, {0 ., 0 ., 0 .}, trkIn);
92+ // tlv.Print();
93+ // return fmt::format("X:{:+.4e} Alp:{:+.3e} Par: {:+.4e} {:+.4e} {:+.4e} {:+.4e} {:+.4e} |Q|:{:d} {:s}\n",
94+ // getX(), getAlpha(), getY(), getZ(), getSnp(), getTgl(), getQ2Pt(), getAbsCharge(), getPID().getName());
95+ // trkIn.print();
6896 o2::track::TrackParCov trkOut;
6997 const int status = fat.FastTrack (trkIn, trkOut, nch);
70- if (status <= 0 ) {
98+ if (status <= mAtLeastHits ) {
7199 Printf (" --- fatSolve: FastTrack failed --- \n " );
72- tlv.Print ();
100+ // tlv.Print();
73101 return false ;
74102 }
103+ Printf (" --- fatSolve: FastTrack succeeded %d --- \n " , status);
104+ // trkOut.print();
75105 lutEntry.valid = true ;
76106 lutEntry.itof = fat.GetGoodHitProb (itof);
77107 lutEntry.otof = fat.GetGoodHitProb (otof);
@@ -81,13 +111,15 @@ bool DelphesO2LutWriter::fatSolve(lutEntry_t& lutEntry,
81111 // define the efficiency
82112 auto totfake = 0 .;
83113 lutEntry.eff = 1 .;
84- for (int i = 1 ; i < 20 ; ++i) {
114+ for (size_t i = 1 ; i < fat.GetNLayers (); ++i) {
115+ if (fat.IsLayerInert (i))
116+ continue ; // skip inert layers
85117 auto igoodhit = fat.GetGoodHitProb (i);
86118 if (igoodhit <= 0 . || i == itof || i == otof)
87119 continue ;
88120 lutEntry.eff *= igoodhit;
89121 auto pairfake = 0 .;
90- for (int j = i + 1 ; j < 20 ; ++j) {
122+ for (size_t j = i + 1 ; j < fat. GetNLayers () ; ++j) {
91123 auto jgoodhit = fat.GetGoodHitProb (j);
92124 if (jgoodhit <= 0 . || j == itof || j == otof)
93125 continue ;
@@ -180,7 +212,7 @@ bool DelphesO2LutWriter::fwdPara(lutEntry_t& lutEntry, float pt, float eta, floa
180212 return true ;
181213}
182214
183- void DelphesO2LutWriter::lutWrite (const char * filename, int pdg, float field, int itof, int otof)
215+ void DelphesO2LutWriter::lutWrite (const char * filename, int pdg, float field, size_t itof, size_t otof)
184216{
185217
186218 if (useFlatDipole && useDipole) {
@@ -206,26 +238,21 @@ void DelphesO2LutWriter::lutWrite(const char* filename, int pdg, float field, in
206238 return ;
207239 }
208240 lutHeader.field = field;
241+ auto setMap = [](map_t & map, LutBinning b) {
242+ map.log = b.log ;
243+ map.nbins = b.nbins ;
244+ map.min = b.min ;
245+ map.max = b.max ;
246+ };
209247 // nch
210- lutHeader.nchmap .log = true ;
211- lutHeader.nchmap .nbins = 20 ;
212- lutHeader.nchmap .min = 0.5 ;
213- lutHeader.nchmap .max = 3.5 ;
248+ setMap (lutHeader.nchmap , mNchBinning );
214249 // radius
215- lutHeader.radmap .log = false ;
216- lutHeader.radmap .nbins = 1 ;
217- lutHeader.radmap .min = 0 .;
218- lutHeader.radmap .max = 100 .;
250+ setMap (lutHeader.radmap , mRadiusBinning );
219251 // eta
220- lutHeader.etamap .log = false ;
221- lutHeader.etamap .nbins = 80 ;
222- lutHeader.etamap .min = -4 .;
223- lutHeader.etamap .max = 4 .;
252+ setMap (lutHeader.etamap , mEtaBinning );
224253 // pt
225- lutHeader.ptmap .log = true ;
226- lutHeader.ptmap .nbins = 200 ;
227- lutHeader.ptmap .min = -2 ;
228- lutHeader.ptmap .max = 2 .;
254+ setMap (lutHeader.ptmap , mPtBinning );
255+
229256 lutFile.write (reinterpret_cast <char *>(&lutHeader), sizeof (lutHeader));
230257
231258 // entries
@@ -247,11 +274,11 @@ void DelphesO2LutWriter::lutWrite(const char* filename, int pdg, float field, in
247274 for (int irad = 0 ; irad < nrad; ++irad) {
248275 Printf (" --- writing irad = %d/%d" , irad, nrad);
249276 for (int ieta = 0 ; ieta < neta; ++ieta) {
250- nCalls++;
251277 Printf (" --- writing ieta = %d/%d" , ieta, neta);
252278 auto eta = lutHeader.etamap .eval (ieta);
253279 lutEntry.eta = lutHeader.etamap .eval (ieta);
254280 for (int ipt = 0 ; ipt < npt; ++ipt) {
281+ nCalls++;
255282 Printf (" --- writing ipt = %d/%d" , ipt, npt);
256283 lutEntry.pt = lutHeader.ptmap .eval (ipt);
257284 lutEntry.valid = true ;
@@ -325,7 +352,7 @@ void DelphesO2LutWriter::diagonalise(lutEntry_t& lutEntry)
325352 }
326353 }
327354
328- m.Print ();
355+ // m.Print();
329356 TMatrixDSymEigen eigen (m);
330357 // eigenvalues vector
331358 TVectorD eigenVal = eigen.GetEigenValues ();
@@ -345,7 +372,7 @@ void DelphesO2LutWriter::diagonalise(lutEntry_t& lutEntry)
345372
346373TGraph* DelphesO2LutWriter::lutRead (const char * filename, int pdg, int what, int vs, float nch, float radius, float eta, float pt)
347374{
348- Printf ( " --- reading LUT file %s" , filename);
375+ LOGF (info, " --- reading LUT file %s" , filename);
349376 // vs
350377 static const int kNch = 0 ;
351378 static const int kEta = 1 ;
@@ -363,6 +390,7 @@ TGraph* DelphesO2LutWriter::lutRead(const char* filename, int pdg, int what, int
363390 o2::delphes::DelphesO2TrackSmearer smearer;
364391 smearer.loadTable (pdg, filename);
365392 auto lutHeader = smearer.getLUTHeader (pdg);
393+ lutHeader->print ();
366394 map_t lutMap;
367395 switch (vs) {
368396 case kNch :
0 commit comments