Skip to content

Commit 270f2bd

Browse files
cnkostercnkoster
committed
Add possibility to recenter for timestamp intervals and add cutanalysis for different centrality estimators
1 parent c985c8c commit 270f2bd

File tree

1 file changed

+119
-26
lines changed

1 file changed

+119
-26
lines changed

‎PWGCF/Flow/TableProducer/zdcQVectors.cxx‎

Lines changed: 119 additions & 26 deletions
Original file line numberDiff line numberDiff line change
@@ -74,6 +74,7 @@ int counter = 0;
7474
// Energy calibration:
7575
std::vector<TString> namesEcal(10, "");
7676
std::vector<std::vector<TString>> names(5, std::vector<TString>()); //(1x 4d 4x 1d)
77+
std::vector<TString> namesTS; // for timestamo recentering
7778
std::vector<TString> vnames = {"hvertex_vx", "hvertex_vy"};
7879

7980
// https://alice-notes.web.cern.ch/system/files/notes/analysis/620/017-May-31-analysis_note-ALICE_analysis_note_v2.pdf
@@ -97,6 +98,7 @@ int lastRunNumber = 0;
9798
std::vector<float> v(3, 0); // vx, vy, vz
9899
bool isSelected = true;
99100
std::vector<float> cents; // centrality estimaters
101+
uint64_t timestamp = 0;
100102

101103
} // namespace o2::analysis::qvectortask
102104

@@ -122,6 +124,11 @@ struct ZdcQVectors {
122124
O2_DEFINE_CONFIGURABLE(cfgCentMax, float, 90, "Maximum cenrality for selected events");
123125
} EvSel;
124126

127+
struct : ConfigurableGroup {
128+
O2_DEFINE_CONFIGURABLE(cfgRecenterForTimestamp, bool, false, "Add 1D recentering for timestamp");
129+
O2_DEFINE_CONFIGURABLE(cfgCCDBdir_Timestamp, std::string, "Users/c/ckoster/ZDC/LHC23_PbPb_pass5/Timestamp", "CCDB directory for Timestamp recentering");
130+
} extraTS;
131+
125132
ConfigurableAxis axisCent{"axisCent", {90, 0, 90}, "Centrality axis in 1% bins"};
126133
ConfigurableAxis axisCent10{"axisCent10", {9, 0, 90}, "Centrality axis in 10% bins"};
127134
ConfigurableAxis axisQ{"axisQ", {100, -2, 2}, "Q vector (xy) in ZDC"};
@@ -131,6 +138,7 @@ struct ZdcQVectors {
131138
ConfigurableAxis axisVx{"axisVx", {100, -0.01, 0.01}, "for Pos X of collision"};
132139
ConfigurableAxis axisVy{"axisVy", {100, -0.01, 0.01}, "for Pos Y of collision"};
133140
ConfigurableAxis axisVz{"axisVz", {100, -10, 10}, "for vz of collision"};
141+
ConfigurableAxis axisTimestamp{"axisTimestamp", {100, 0, 100}, "for timestamp of collision in (TSi - TS_start) / (TS_eind - TS_start) x 100% "};
134142

135143
// Centrality Estimators -> standard is FT0C
136144
O2_DEFINE_CONFIGURABLE(cfgFT0Cvariant1, bool, false, "Set centrality estimator to cfgFT0Cvariant1");
@@ -178,7 +186,8 @@ struct ZdcQVectors {
178186
enum CalibModes {
179187
kEnergyCal,
180188
kMeanv,
181-
kRec
189+
kRec,
190+
kTimestamp
182191
};
183192

184193
// Define output
@@ -188,8 +197,8 @@ struct ZdcQVectors {
188197

189198
// keep track of calibration histos for each given step and iteration
190199
struct Calib {
191-
std::vector<TList*> calibList = std::vector<TList*>(3, nullptr); // [0] Enerfy cal, [1] vmean, [2] recentering
192-
std::vector<bool> calibfilesLoaded = std::vector<bool>(3, false);
200+
std::vector<TList*> calibList = std::vector<TList*>(4, nullptr); // [0] Enerfy cal, [1] vmean, [2] recentering, [3] timestamp
201+
std::vector<bool> calibfilesLoaded = std::vector<bool>(4, false);
193202
int atStep = 0;
194203
int atIteration = 0;
195204

@@ -252,6 +261,7 @@ struct ZdcQVectors {
252261
names[2].push_back(TString::Format("hQ%s%s_mean_vx_run", coord, side));
253262
names[3].push_back(TString::Format("hQ%s%s_mean_vy_run", coord, side));
254263
names[4].push_back(TString::Format("hQ%s%s_mean_vz_run", coord, side));
264+
namesTS.push_back(TString::Format("hQ%s%s_mean_timestamp_run", coord, side));
255265
} // end of capCOORDS
256266
}
257267

@@ -271,6 +281,7 @@ struct ZdcQVectors {
271281
registry.add<TProfile>(Form("QA/before/hQ%sA_Q%sC_vs_vx", COORD1, COORD2), Form("hQ%sA_Q%sC_vs_vx", COORD1, COORD2), kTProfile, {axisVx});
272282
registry.add<TProfile>(Form("QA/before/hQ%sA_Q%sC_vs_vy", COORD1, COORD2), Form("hQ%sA_Q%sC_vs_vy", COORD1, COORD2), kTProfile, {axisVy});
273283
registry.add<TProfile>(Form("QA/before/hQ%sA_Q%sC_vs_vz", COORD1, COORD2), Form("hQ%sA_Q%sC_vs_vz", COORD1, COORD2), kTProfile, {axisVz});
284+
registry.add<TProfile>(Form("QA/before/hQ%sA_Q%sC_vs_timestamp", COORD1, COORD2), Form("hQ%sA_Q%sC_vs_timestamp", COORD1, COORD2), kTProfile, {axisTimestamp});
274285
}
275286
}
276287

@@ -282,6 +293,7 @@ struct ZdcQVectors {
282293
registry.add(Form("QA/before/hQ%s%s_vs_vx", coord, side), Form("hQ%s%s_vs_vx", coord, side), {HistType::kTProfile, {axisVx}});
283294
registry.add(Form("QA/before/hQ%s%s_vs_vy", coord, side), Form("hQ%s%s_vs_vy", coord, side), {HistType::kTProfile, {axisVy}});
284295
registry.add(Form("QA/before/hQ%s%s_vs_vz", coord, side), Form("hQ%s%s_vs_vz", coord, side), {HistType::kTProfile, {axisVz}});
296+
registry.add(Form("QA/before/hQ%s%s_vs_timestamp", coord, side), Form("hQ%s%s_vs_timestamp", coord, side), {HistType::kTProfile, {axisTimestamp}});
285297
registry.add(Form("QA/Q%s%s_vs_iteration", coord, side), Form("hQ%s%s_vs_iteration", coord, side), {HistType::kTH2D, {{25, 0, 25}, axisQ}});
286298
} // end of capCOORDS
287299
} // end of sides
@@ -361,16 +373,24 @@ struct ZdcQVectors {
361373
// Tower mean energies vs. centrality used for tower gain equalisation
362374
int totalTowers = 10;
363375
for (int tower = 0; tower < totalTowers; tower++) {
364-
registry.add<TProfile2D>(Form("CutAnalysis/%s", namesEcal[tower].Data()), Form("%s", namesEcal[tower].Data()), kTProfile2D, {axisCent, {nEventSelections, 0, nEventSelections}});
376+
registry.add<TProfile2D>(Form("CutAnalysis/%s", namesEcal[tower].Data()), Form("%s", namesEcal[tower].Data()), kTProfile2D, {axisCent, {nEventSelections + 5, 0, nEventSelections + 5}});
365377
}
366378
// recentered q-vectors (to check what steps are finished in the end)
367-
registry.add<TProfile2D>("CutAnalysis/hvertex_vx", "hvertex_vx", kTProfile2D, {{1, 0., 1.}, {nEventSelections, 0, nEventSelections}});
368-
registry.add<TProfile2D>("CutAnalysis/hvertex_vy", "hvertex_vy", kTProfile2D, {{1, 0., 1.}, {nEventSelections, 0, nEventSelections}});
369-
registry.add<TProfile2D>("CutAnalysis/hvertex_vz", "hvertex_vz", kTProfile2D, {{1, 0., 1.}, {nEventSelections, 0, nEventSelections}});
379+
registry.add<TProfile2D>("CutAnalysis/hvertex_vx", "hvertex_vx", kTProfile2D, {{1, 0., 1.}, {nEventSelections + 5, 0, nEventSelections + 5}});
380+
registry.add<TProfile2D>("CutAnalysis/hvertex_vy", "hvertex_vy", kTProfile2D, {{1, 0., 1.}, {nEventSelections + 5, 0, nEventSelections + 5}});
381+
registry.add<TProfile2D>("CutAnalysis/hvertex_vz", "hvertex_vz", kTProfile2D, {{1, 0., 1.}, {nEventSelections + 5, 0, nEventSelections + 5}});
370382
}
371383
}
372384
}
373385

386+
double rescaleTimestamp(uint64_t timestamp, int runnumber) {
387+
auto& ccdb = o2::ccdb::BasicCCDBManager::instance();
388+
auto duration = ccdb.getRunDuration(runnumber);
389+
double ts = (static_cast<double>(timestamp - duration.first) / static_cast<double> (duration.second - duration.first)) * 100.0;
390+
391+
return ts;
392+
}
393+
374394
template <typename TCollision, typename TZdc>
375395
inline void fillCutAnalysis(TCollision collision, TZdc zdcBC, int evSel)
376396
{
@@ -379,6 +399,7 @@ struct ZdcQVectors {
379399

380400
if (!cfgFillCutAnalysis || cfgFillNothing)
381401
return;
402+
// Add default with different centrality estimators as well
382403
// Here we fill the Energy and mean vx, vy vz histograms with an extra dimension for all the event selections used.
383404
registry.get<TProfile2D>(HIST("CutAnalysis/hvertex_vx"))->Fill(Form("%d", runnumber), evSel, collision.posX());
384405
registry.get<TProfile2D>(HIST("CutAnalysis/hvertex_vy"))->Fill(Form("%d", runnumber), evSel, collision.posY());
@@ -394,6 +415,34 @@ struct ZdcQVectors {
394415
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t2_cent"))->Fill(centrality, evSel, zdcBC.energySectorZNC()[1], 1);
395416
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t3_cent"))->Fill(centrality, evSel, zdcBC.energySectorZNC()[2], 1);
396417
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t4_cent"))->Fill(centrality, evSel, zdcBC.energySectorZNC()[3], 1);
418+
419+
420+
if(evSel == nEventSelections) {
421+
int centCounter = 0;
422+
423+
std::vector<float> cents = {
424+
collision.centFT0C(),
425+
collision.centFT0CVariant1(),
426+
collision.centFT0M(),
427+
collision.centFV0A(),
428+
collision.centNGlobal()
429+
};
430+
431+
for(auto& cent : cents) {
432+
LOGF(info, "Filling centrality: %f for event selection: %d", cent, evSel + centCounter);
433+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNA_mean_t0_cent"))->Fill(cent, evSel + centCounter, zdcBC.energyCommonZNA(), 1);
434+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNA_mean_t1_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNA()[0], 1);
435+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNA_mean_t2_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNA()[1], 1);
436+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNA_mean_t3_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNA()[2], 1);
437+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNA_mean_t4_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNA()[3], 1);
438+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t0_cent"))->Fill(cent, evSel + centCounter, zdcBC.energyCommonZNC(), 1);
439+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t1_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNC()[0], 1);
440+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t2_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNC()[1], 1);
441+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t3_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNC()[2], 1);
442+
registry.get<TProfile2D>(HIST("CutAnalysis/hZNC_mean_t4_cent"))->Fill(cent, evSel + centCounter, zdcBC.energySectorZNC()[3], 1);
443+
centCounter++;
444+
}
445+
}
397446
}
398447

399448
template <typename TCollision, typename TBunchCrossing>
@@ -474,6 +523,10 @@ struct ZdcQVectors {
474523
fillCutAnalysis(collision, bunchCrossing, evSel_RCTFlagsZDC);
475524
}
476525

526+
// Fill for centrality estimators!
527+
fillCutAnalysis(collision, bunchCrossing, nEventSelections);
528+
529+
477530
return selectionBits;
478531
}
479532

@@ -529,6 +582,16 @@ struct ZdcQVectors {
529582
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQYA_QXC_vs_vz"), v[2], qya * qxc);
530583
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQXA_QYC_vs_vz"), v[2], qxa * qyc);
531584

585+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQXA_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qxa);
586+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQYA_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qya);
587+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQXC_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qxc);
588+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQYC_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qyc);
589+
590+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQXA_QXC_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qxa * qxc);
591+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQYA_QYC_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qya * qyc);
592+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQYA_QXC_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qya * qxc);
593+
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/hQXA_QYC_vs_timestamp"), rescaleTimestamp(timestamp, runnumber), qxa * qyc);
594+
532595
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/ZNA_Qx_vs_Centrality"), centrality, qxa);
533596
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/ZNA_Qy_vs_Centrality"), centrality, qya);
534597
registry.fill(HIST("QA/") + HIST(Time[ft]) + HIST("/ZNC_Qx_vs_Centrality"), centrality, qxc);
@@ -544,7 +607,7 @@ struct ZdcQVectors {
544607
}
545608

546609
template <CalibModes cm>
547-
void loadCalibrations(uint64_t timestamp, std::string ccdb_dir)
610+
void loadCalibrations(std::string ccdb_dir)
548611
{
549612
// iteration = 0 (Energy calibration) -> step 0 only
550613
// iteration 1,2,3,4,5 = recentering -> 5 steps per iteration (1x 4D + 4x 1D)
@@ -570,11 +633,11 @@ struct ZdcQVectors {
570633
T* hist = nullptr;
571634
double calibConstant{0};
572635

573-
if (cm == kEnergyCal) {
636+
if (cm == kEnergyCal || cm == kMeanv) {
574637
TList* list = cal.calibList[cm];
575638
hist = reinterpret_cast<T*>(list->FindObject(Form("%s", objName)));
576-
} else if (cm == kMeanv) {
577-
TList* list = cal.calibList[cm];
639+
} else if (cm == kTimestamp) {
640+
TList* list = reinterpret_cast<TList*>(cal.calibList[cm]->FindObject(Form("it%i_step%i", iteration, step)));
578641
hist = reinterpret_cast<T*>(list->FindObject(Form("%s", objName)));
579642
} else if (cm == kRec) {
580643
TList* list = reinterpret_cast<TList*>(cal.calibList[cm]->FindObject(Form("it%i_step%i", iteration, step)));
@@ -588,10 +651,10 @@ struct ZdcQVectors {
588651
cal.atStep = step;
589652
cal.atIteration = iteration;
590653
}
591-
592654
if (!hist) {
593655
LOGF(fatal, "%s not available.. Abort..", objName);
594656
}
657+
595658

596659
if (hist->InheritsFrom("TProfile2D")) {
597660
// needed for energy calibration!
@@ -604,16 +667,24 @@ struct ZdcQVectors {
604667
TProfile* h = reinterpret_cast<TProfile*>(hist);
605668
TString name = h->GetName();
606669
int bin{};
607-
if (name.Contains("mean_vx"))
670+
if (name.Contains("mean_vx")) {
608671
bin = h->GetXaxis()->FindBin(v[0]);
609-
if (name.Contains("mean_vy"))
672+
}
673+
if (name.Contains("mean_vy")) {
610674
bin = h->GetXaxis()->FindBin(v[1]);
611-
if (name.Contains("mean_vz"))
675+
}
676+
if (name.Contains("mean_vz")) {
612677
bin = h->GetXaxis()->FindBin(v[2]);
613-
if (name.Contains("mean_cent"))
678+
}
679+
if (name.Contains("mean_cent")) {
614680
bin = h->GetXaxis()->FindBin(centrality);
615-
if (name.Contains("vertex"))
681+
}
682+
if (name.Contains("vertex")) {
616683
bin = h->GetXaxis()->FindBin(TString::Format("%i", runnumber));
684+
}
685+
if (name.Contains("timestamp")) {
686+
bin = h->GetXaxis()->FindBin(rescaleTimestamp(timestamp, runnumber));
687+
}
617688
calibConstant = h->GetBinContent(bin);
618689
} else if (hist->InheritsFrom("THnSparse")) {
619690
std::vector<int> sparsePars;
@@ -693,6 +764,8 @@ struct ZdcQVectors {
693764

694765
registry.fill(HIST("hEventCount"), evSel_FilteredEvent);
695766

767+
timestamp = foundBC.timestamp();
768+
696769
if (!foundBC.has_zdc()) {
697770
isSelected = false;
698771
spTableZDC(runnumber, cents, v, foundBC.timestamp(), 0, 0, 0, 0, isSelected, 0);
@@ -761,21 +834,24 @@ struct ZdcQVectors {
761834
cal.calibfilesLoaded[2] = false;
762835
cal.calibList[2] = nullptr;
763836

837+
cal.calibfilesLoaded[3] = false;
838+
cal.calibList[3] = nullptr;
839+
764840
cal.isShiftProfileFound = false;
765841
cal.shiftprofileC = nullptr;
766842
cal.shiftprofileA = nullptr;
767843
}
768844

769845
// load the calibration histos for iteration 0 step 0 (Energy Calibration)
770-
loadCalibrations<kEnergyCal>(foundBC.timestamp(), cfgEnergyCal.value);
846+
loadCalibrations<kEnergyCal>(cfgEnergyCal.value);
771847

772848
if (!cal.calibfilesLoaded[0]) {
773849
if (counter < 1) {
774850
LOGF(info, " --> No Energy calibration files found.. -> Only Energy calibration will be done. ");
775851
}
776852
}
777853
// load the calibrations for the mean v
778-
loadCalibrations<kMeanv>(foundBC.timestamp(), cfgMeanv.value);
854+
loadCalibrations<kMeanv>(cfgMeanv.value);
779855

780856
if (!cfgFillNothing) {
781857
registry.get<TProfile>(HIST("vmean/hvertex_vx"))->Fill(Form("%d", runnumber), v[0]);
@@ -928,7 +1004,11 @@ struct ZdcQVectors {
9281004
return;
9291005
}
9301006

931-
loadCalibrations<kRec>(foundBC.timestamp(), cfgRec.value);
1007+
loadCalibrations<kRec>(cfgRec.value);
1008+
1009+
if(extraTS.cfgRecenterForTimestamp) {
1010+
loadCalibrations<kTimestamp>(extraTS.cfgCCDBdir_Timestamp.value);
1011+
}
9321012

9331013
std::vector<double> qRec(q);
9341014

@@ -981,15 +1061,28 @@ struct ZdcQVectors {
9811061
corrQyC.push_back(getCorrection<TProfile, kRec>(names[step - 1][3].Data(), it, step));
9821062
pb++;
9831063
}
1064+
1065+
if(extraTS.cfgRecenterForTimestamp) {
1066+
corrQxA.push_back(getCorrection<TProfile, kTimestamp>(namesTS[0].Data(), it, 6));
1067+
corrQyA.push_back(getCorrection<TProfile, kTimestamp>(namesTS[1].Data(), it, 6));
1068+
corrQxC.push_back(getCorrection<TProfile, kTimestamp>(namesTS[2].Data(), it, 6));
1069+
corrQyC.push_back(getCorrection<TProfile, kTimestamp>(namesTS[3].Data(), it, 6));
1070+
pb++;
1071+
}
9841072
}
9851073

986-
for (int cor = 0; cor < pb; cor++) {
987-
qRec[0] -= corrQxA[cor];
988-
qRec[1] -= corrQyA[cor];
989-
qRec[2] -= corrQxC[cor];
990-
qRec[3] -= corrQyC[cor];
1074+
double totalCorrectionQxA = std::accumulate(corrQxA.begin(), corrQxA.end(), 0.0);
1075+
double totalCorrectionQyA = std::accumulate(corrQyA.begin(), corrQyA.end(), 0.0);
1076+
double totalCorrectionQxC = std::accumulate(corrQxC.begin(), corrQxC.end(), 0.0);
1077+
double totalCorrectionQyC = std::accumulate(corrQyC.begin(), corrQyC.end(), 0.0);
1078+
1079+
qRec[0] -= totalCorrectionQxA;
1080+
qRec[1] -= totalCorrectionQyA;
1081+
qRec[2] -= totalCorrectionQxC;
1082+
qRec[3] -= totalCorrectionQyC;
9911083

992-
if (cfgFillHistRegistry && !cfgFillNothing) {
1084+
if (cfgFillHistRegistry && !cfgFillNothing) {
1085+
for (int cor = 0; cor < pb; cor++) {
9931086
registry.get<TH2>(HIST("QA/QXA_vs_iteration"))->Fill(cor, qRec[0]);
9941087
registry.get<TH2>(HIST("QA/QYA_vs_iteration"))->Fill(cor, qRec[1]);
9951088
registry.get<TH2>(HIST("QA/QXC_vs_iteration"))->Fill(cor, qRec[2]);

0 commit comments

Comments
 (0)