Skip to content

Commit d0dc689

Browse files
Luzhiyonggalibuild
Luzhiyongg
and
alibuild
authored
[PWGCF] add manual ck calculations (#18342)
Co-authored-by: ALICE Action Bot <alibuild@cern.ch>
1 parent a4ecfca commit d0dc689

1 file changed

Lines changed: 58 additions & 19 deletions

File tree

‎PWGCF/Flow/Tasks/flowTask.cxx‎

Lines changed: 58 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -243,6 +243,7 @@ struct FlowTask {
243243
std::vector<float> nPt; // fraction of particles in each pT bin in one event
244244
std::vector<std::shared_ptr<TProfile2D>> listNptV2;
245245
std::vector<std::shared_ptr<TH2>> listPtX;
246+
O2_DEFINE_CONFIGURABLE(cfgManualCk, bool, false, "Produce pt-pt variance")
246247
} cfgAdditionObs;
247248

248249
ConfigurableAxis axisPt{"axisPt", {VARIABLE_WIDTH, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2, 2.2, 2.4, 2.6, 2.8, 3, 3.5, 4, 5, 6, 8, 10}, "pt axis for histograms"};
@@ -287,7 +288,7 @@ struct FlowTask {
287288
std::vector<GFW::CorrConfig> corrconfigs;
288289
GFWCorrConfigs gfwConfigs;
289290
std::vector<GFW::CorrConfig> corrconfigsPtVn;
290-
TAxis* fPtAxis;
291+
TAxis* fPtAxis = nullptr;
291292
TRandom3* fRndm = new TRandom3(0);
292293
std::vector<std::pair<double, double>> etagapsPtPt;
293294
std::vector<std::vector<std::shared_ptr<TProfile>>> bootstrapArray;
@@ -324,15 +325,15 @@ struct FlowTask {
324325
double mMinSeconds{-1.};
325326
std::unordered_map<int, TH2*> gHadronicRate;
326327
ctpRateFetcher mRateFetcher;
327-
TH2* gCurrentHadronicRate;
328+
TH2* gCurrentHadronicRate = nullptr;
328329
RCTFlagsChecker rctChecker{"CBT"};
329330

330331
// phi-EP correction
331332
std::vector<TF1*> funcEff;
332-
TH1D* hFindPtBin;
333-
TF1* funcV2;
334-
TF1* funcV3;
335-
TF1* funcV4;
333+
TH1D* hFindPtBin = nullptr;
334+
TF1* funcV2 = nullptr;
335+
TF1* funcV3 = nullptr;
336+
TF1* funcV4 = nullptr;
336337

337338
void init(InitContext const&)
338339
{
@@ -439,6 +440,13 @@ struct FlowTask {
439440
registry.add("hPtA_afterDptCut", "p_{T} distribution", {HistType::kTH1D, {axisPt}});
440441
registry.add("hPtB_afterDptCut", "p_{T} distribution", {HistType::kTH1D, {axisPt}});
441442
}
443+
if (cfgAdditionObs.cfgManualCk) {
444+
registry.add("ptpt_1sub1_1sub3_2pc", "", {HistType::kTProfile, {axisIndependent}});
445+
registry.add("pt_1sub1_1pc", "", {HistType::kTProfile, {axisIndependent}});
446+
registry.add("pt_1sub1_2pc", "", {HistType::kTProfile, {axisIndependent}});
447+
registry.add("pt_1sub3_1pc", "", {HistType::kTProfile, {axisIndependent}});
448+
registry.add("pt_1sub3_2pc", "", {HistType::kTProfile, {axisIndependent}});
449+
}
442450
if (doprocessMCGen) {
443451
registry.add("MCGen/MChPhi", "#phi distribution", {HistType::kTH1D, {axisPhi}});
444452
registry.add("MCGen/MChEta", "#eta distribution", {HistType::kTH1D, {axisEta}});
@@ -910,45 +918,45 @@ struct FlowTask {
910918
if (cfgUserIO.cfgAcceptance.value.empty() == false) {
911919
mAcceptance = ccdb->getForTimeStamp<GFWWeights>(cfgUserIO.cfgAcceptance, timestamp);
912920
if (mAcceptance)
913-
LOGF(info, "Loaded acceptance weights from %s (%p)", cfgUserIO.cfgAcceptance.value.c_str(), (void*)mAcceptance);
921+
LOGF(info, "Loaded acceptance weights from %s (%p)", cfgUserIO.cfgAcceptance.value.c_str(), static_cast<void*>(mAcceptance));
914922
else
915-
LOGF(warning, "Could not load acceptance weights from %s (%p)", cfgUserIO.cfgAcceptance.value.c_str(), (void*)mAcceptance);
923+
LOGF(warning, "Could not load acceptance weights from %s (%p)", cfgUserIO.cfgAcceptance.value.c_str(), static_cast<void*>(mAcceptance));
916924
}
917925
if (cfgUserIO.cfgEfficiency.value.empty() == false) {
918926
mEfficiency = ccdb->getForTimeStamp<TH1D>(cfgUserIO.cfgEfficiency, timestamp);
919927
if (mEfficiency == nullptr) {
920928
LOGF(fatal, "Could not load efficiency histogram for trigger particles from %s", cfgUserIO.cfgEfficiency.value.c_str());
921929
}
922-
LOGF(info, "Loaded efficiency histogram from %s (%p)", cfgUserIO.cfgEfficiency.value.c_str(), (void*)mEfficiency);
930+
LOGF(info, "Loaded efficiency histogram from %s (%p)", cfgUserIO.cfgEfficiency.value.c_str(), static_cast<void*>(mEfficiency));
923931
}
924932
if (cfgUserIO.cfgEfficiencyForNch.value.empty() == false) {
925933
mEfficiencyForNch = ccdb->getForTimeStamp<TH1D>(cfgUserIO.cfgEfficiencyForNch, timestamp);
926934
if (mEfficiencyForNch == nullptr) {
927935
LOGF(fatal, "Could not load efficiency histogram for trigger particles from %s", cfgUserIO.cfgEfficiencyForNch.value.c_str());
928936
}
929-
LOGF(info, "Loaded efficiency histogram from %s (%p)", cfgUserIO.cfgEfficiencyForNch.value.c_str(), (void*)mEfficiencyForNch);
937+
LOGF(info, "Loaded efficiency histogram from %s (%p)", cfgUserIO.cfgEfficiencyForNch.value.c_str(), static_cast<void*>(mEfficiencyForNch));
930938
}
931939
if (cfgAdditionObs.cfgDptDisEnable && cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.empty() == false) {
932940
cfgAdditionObs.hEvAvgMeanPt = ccdb->getForTimeStamp<TH1D>(cfgAdditionObs.cfgDptDishEvAvgMeanPt, timestamp);
933941
if (cfgAdditionObs.hEvAvgMeanPt == nullptr) {
934942
LOGF(fatal, "Could not load mean pT histogram from %s", cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.c_str());
935943
}
936-
LOGF(info, "Loaded mean pT histogram from %s (%p)", cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.c_str(), (void*)cfgAdditionObs.hEvAvgMeanPt);
944+
LOGF(info, "Loaded mean pT histogram from %s (%p)", cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.c_str(), static_cast<void*>(cfgAdditionObs.hEvAvgMeanPt));
937945
}
938946
if (cfgAdditionObs.cfgDptDisEnable && cfgAdditionObs.cfgDptDisSelectionSwitch > kNoDptCut) {
939947
if (cfgAdditionObs.cfgDptDisCutLow.value.empty() == false) {
940948
cfgAdditionObs.fDptDisCutLow = ccdb->getForTimeStamp<TH1D>(cfgAdditionObs.cfgDptDisCutLow, timestamp);
941949
if (cfgAdditionObs.fDptDisCutLow == nullptr) {
942950
LOGF(fatal, "Could not load dptDis low cut histogram from %s", cfgAdditionObs.cfgDptDisCutLow.value.c_str());
943951
}
944-
LOGF(info, "Loaded dptDis low cut histogram from %s (%p)", cfgAdditionObs.cfgDptDisCutLow.value.c_str(), (void*)cfgAdditionObs.fDptDisCutLow);
952+
LOGF(info, "Loaded dptDis low cut histogram from %s (%p)", cfgAdditionObs.cfgDptDisCutLow.value.c_str(), static_cast<void*>(cfgAdditionObs.fDptDisCutLow));
945953
}
946954
if (cfgAdditionObs.cfgDptDisCutHigh.value.empty() == false) {
947955
cfgAdditionObs.fDptDisCutHigh = ccdb->getForTimeStamp<TH1D>(cfgAdditionObs.cfgDptDisCutHigh, timestamp);
948956
if (cfgAdditionObs.fDptDisCutHigh == nullptr) {
949957
LOGF(fatal, "Could not load dptDis high cut histogram from %s", cfgAdditionObs.cfgDptDisCutHigh.value.c_str());
950958
}
951-
LOGF(info, "Loaded dptDis high cut histogram from %s (%p)", cfgAdditionObs.cfgDptDisCutHigh.value.c_str(), (void*)cfgAdditionObs.fDptDisCutHigh);
959+
LOGF(info, "Loaded dptDis high cut histogram from %s (%p)", cfgAdditionObs.cfgDptDisCutHigh.value.c_str(), static_cast<void*>(cfgAdditionObs.fDptDisCutHigh));
952960
}
953961
}
954962
correctionsLoaded = true;
@@ -1289,6 +1297,10 @@ struct FlowTask {
12891297
double weffEventWithinGap08 = 0., weffEventSquareWithinGap08 = 0.;
12901298
double sumPtsquareWsquareWithinGap08 = 0., sumPtWsquareWithinGap08 = 0.;
12911299
double nTracksCorrected = 0;
1300+
double mpt_sub1 = 0;
1301+
double mpt_sub1_weffEv = 0;
1302+
double mpt_sub3 = 0;
1303+
double mpt_sub3_weffEv = 0;
12921304
int magnetfield = 0;
12931305
float independent = cent;
12941306
float nTracksUncorrected = 0;
@@ -1405,6 +1417,16 @@ struct FlowTask {
14051417
weffEventWithinGap08 += weff;
14061418
weffEventSquareWithinGap08 += weff * weff;
14071419
}
1420+
if (cfgAdditionObs.cfgManualCk) {
1421+
if (track.eta() >= etagapsPtPt[0].first && track.eta() <= etagapsPtPt[0].second) {
1422+
mpt_sub1 += weff * track.pt();
1423+
mpt_sub1_weffEv += weff;
1424+
}
1425+
if (track.eta() >= etagapsPtPt[2].first && track.eta() <= etagapsPtPt[2].second) {
1426+
mpt_sub3 += weff * track.pt();
1427+
mpt_sub3_weffEv += weff;
1428+
}
1429+
}
14081430
}
14091431
if (cfgUserIO.cfgConsistentEventFlag && consistentEventVector.size() == 6) { // o2-linter: disable=magic-number (size match)
14101432
acceptedTracks.nFull += 1;
@@ -1477,13 +1499,12 @@ struct FlowTask {
14771499
if (weffEvent) {
14781500
registry.fill(HIST("hMeanPt"), independent, ptSum / weffEvent, weffEvent);
14791501
}
1480-
if (weffEventWithinGap08)
1502+
if (weffEventWithinGap08) {
14811503
registry.fill(HIST("hMeanPtWithinGap08"), independent, ptSum_Gap08 / weffEventWithinGap08, cfgUserIO.cfgUserPtVnEvWeightEnabled ? weffEventWithinGap08 : 1.);
1482-
if (weffEventWithinGap08)
14831504
bootstrapArray[sampleIndex][kMeanPtWithinGap08]->Fill(independent, ptSum_Gap08 / weffEventWithinGap08, 1.0);
1484-
// c22_gap8 * pt_withGap8
1485-
if (weffEventWithinGap08)
1505+
// c22_gap8 * pt_withGap8
14861506
fillpTvnProfile(corrconfigs.at(7), ptSum_Gap08, weffEventWithinGap08, HIST("c22_gap08_Weff"), HIST("c22_gap08_trackMeanPt"), independent);
1507+
}
14871508
// PtVariance
14881509
if (weffEventDiffWithGap08) {
14891510
registry.fill(HIST("PtVariance_partA_WithinGap08"), independent,
@@ -1494,6 +1515,25 @@ struct FlowTask {
14941515
weffEventDiffWithGap08);
14951516
}
14961517

1518+
if (cfgAdditionObs.cfgManualCk) {
1519+
// 2pc event weight:
1520+
double ev2pc = mpt_sub1_weffEv * mpt_sub3_weffEv;
1521+
// manual pt-pt variance
1522+
if (mpt_sub1_weffEv) {
1523+
mpt_sub1 = mpt_sub1 / mpt_sub1_weffEv;
1524+
registry.fill(HIST("pt_1sub1_1pc"), independent, mpt_sub1, mpt_sub1_weffEv);
1525+
}
1526+
if (mpt_sub3_weffEv) {
1527+
mpt_sub3 = mpt_sub3 / mpt_sub3_weffEv;
1528+
registry.fill(HIST("pt_1sub3_1pc"), independent, mpt_sub3, mpt_sub3_weffEv);
1529+
}
1530+
if (mpt_sub1_weffEv && mpt_sub3_weffEv) {
1531+
registry.fill(HIST("ptpt_1sub1_1sub3_2pc"), independent, mpt_sub1 * mpt_sub3, ev2pc);
1532+
registry.fill(HIST("pt_1sub1_2pc"), independent, mpt_sub1, ev2pc);
1533+
registry.fill(HIST("pt_1sub3_2pc"), independent, mpt_sub3, ev2pc);
1534+
}
1535+
}
1536+
14971537
// v02
14981538
if (cfgAdditionObs.cfgV02Enabled && weffEventFrac > 0) {
14991539
for (auto ipt = 0; ipt < fPtAxis->GetNbins(); ipt++) {
@@ -1569,9 +1609,8 @@ struct FlowTask {
15691609
ptSum_Gap08 += mcParticle.pt();
15701610
count_Gap08 += 1.;
15711611
}
1572-
}
1573-
if (withinPtRef)
15741612
fGFW->Fill(mcParticle.eta(), fPtAxis->FindBin(mcParticle.pt()) - 1, mcParticle.phi(), 1., 1);
1613+
}
15751614
if (withinPtPOI)
15761615
fGFW->Fill(mcParticle.eta(), fPtAxis->FindBin(mcParticle.pt()) - 1, mcParticle.phi(), 1., 2);
15771616
if (withinPtPOI && withinPtRef)

0 commit comments

Comments
 (0)