Skip to content

Commit 8e908a9

Browse files
LuzhiyonggLuzhiyongg
committed
modify corrfig log info and add V02
1 parent e7f4405 commit 8e908a9

File tree

1 file changed

+92
-34
lines changed

1 file changed

+92
-34
lines changed

‎PWGCF/Flow/Tasks/flowTask.cxx‎

Lines changed: 92 additions & 34 deletions
Original file line numberDiff line numberDiff line change
@@ -164,6 +164,13 @@ struct FlowTask {
164164
O2_DEFINE_CONFIGURABLE(cfgTPCPhiCutPtMin, float, 2.0f, "start point of phi-pt cut")
165165
TF1* fPhiCutLow = nullptr;
166166
TF1* fPhiCutHigh = nullptr;
167+
// Functional form of pt-dependent DCAxy cut
168+
TF1* fPtDepDCAxy = nullptr;
169+
O2_DEFINE_CONFIGURABLE(cfgDCAxyNSigma, float, 0, "0: disable; Cut on number of sigma deviations from expected DCA in the transverse direction, nsigma=7 is the same with global track");
170+
O2_DEFINE_CONFIGURABLE(cfgDCAxy, std::string, "(0.0026+0.005/(x^1.01))", "Functional form of pt-dependent DCAxy cut");
171+
} cfgFuncParas;
172+
173+
struct : ConfigurableGroup {
167174
// for deltaPt/<pT> vs centrality
168175
O2_DEFINE_CONFIGURABLE(cfgDptDisEnable, bool, false, "Produce deltaPt/meanPt vs centrality")
169176
O2_DEFINE_CONFIGURABLE(cfgDptDisSelectionSwitch, int, 0, "0: disable, 1: use low cut, 2:use high cut")
@@ -174,13 +181,14 @@ struct FlowTask {
174181
O2_DEFINE_CONFIGURABLE(cfgDptDisCutLow, std::string, "", "CCDB path to dpt lower boundary")
175182
O2_DEFINE_CONFIGURABLE(cfgDptDisCutHigh, std::string, "", "CCDB path to dpt higher boundary")
176183
ConfigurableAxis cfgDptDisAxisNormal{"cfgDptDisAxisNormal", {200, -1., 1.}, "normalized axis"};
177-
// Functional form of pt-dependent DCAxy cut
178-
TF1* fPtDepDCAxy = nullptr;
179-
O2_DEFINE_CONFIGURABLE(cfgDCAxyNSigma, float, 0, "0: disable; Cut on number of sigma deviations from expected DCA in the transverse direction, nsigma=7 is the same with global track");
180-
O2_DEFINE_CONFIGURABLE(cfgDCAxy, std::string, "(0.0026+0.005/(x^1.01))", "Functional form of pt-dependent DCAxy cut");
181-
} cfgFuncParas;
184+
// for v02(pT)
185+
O2_DEFINE_CONFIGURABLE(cfgV02Enabled, bool, false, "Produce v02(pT)")
186+
O2_DEFINE_CONFIGURABLE(cfgV02CorrIndex, int, 8, "correlation index for n(pT)*corr")
187+
std::vector<float> nPt; // fraction of particles in each pT bin in one event
188+
std::vector<std::shared_ptr<TProfile2D>> listNptV2;
189+
std::vector<std::shared_ptr<TH2>> listPtX;
190+
} cfgAdditionObs;
182191

183-
ConfigurableAxis axisPtHist{"axisPtHist", {100, 0., 10.}, "pt axis for histograms"};
184192
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"};
185193
ConfigurableAxis axisIndependent{"axisIndependent", {VARIABLE_WIDTH, 0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90}, "X axis for histograms"};
186194
ConfigurableAxis axisNch{"axisNch", {4000, 0, 4000}, "N_{ch}"};
@@ -338,8 +346,8 @@ struct FlowTask {
338346
registry.add("hPhi", "#phi distribution", {HistType::kTH1D, {axisPhi}});
339347
registry.add("hPhiWeighted", "corrected #phi distribution", {HistType::kTH1D, {axisPhi}});
340348
registry.add("hEta", "#eta distribution", {HistType::kTH1D, {axisEta}});
341-
registry.add("hPt", "p_{T} distribution before cut", {HistType::kTH1D, {axisPtHist}});
342-
registry.add("hPtRef", "p_{T} distribution after cut", {HistType::kTH1D, {axisPtHist}});
349+
registry.add("hPt", "p_{T} distribution", {HistType::kTH1D, {axisPt}});
350+
registry.add("hPtRef", "p_{T} distribution for refernce particle", {HistType::kTH1D, {axisPt}});
343351
registry.add("pt_phi_bef", "before cut;p_{T};#phi_{modn}", {HistType::kTH2D, {axisPt, cfgFuncParas.axisPhiMod}});
344352
registry.add("pt_phi_aft", "after cut;p_{T};#phi_{modn}", {HistType::kTH2D, {axisPt, cfgFuncParas.axisPhiMod}});
345353
registry.add("hChi2prTPCcls", "#chi^{2}/cluster for the TPC track segment", {HistType::kTH1D, {{100, 0., 5.}}});
@@ -364,13 +372,13 @@ struct FlowTask {
364372
registry.add("c22_gap08_trackMeanPt", "", {HistType::kTProfile, {axisIndependent}});
365373
registry.add("PtVariance_partA_WithinGap08", "", {HistType::kTProfile, {axisIndependent}});
366374
registry.add("PtVariance_partB_WithinGap08", "", {HistType::kTProfile, {axisIndependent}});
367-
if (cfgFuncParas.cfgDptDisEnable) {
368-
registry.add("hNormDeltaPt_X", "; #delta p_{T}/[p_{T}]; X", {HistType::kTH2D, {cfgFuncParas.cfgDptDisAxisNormal, axisIndependent}});
375+
if (cfgAdditionObs.cfgDptDisEnable) {
376+
registry.add("hNormDeltaPt_X", "; #delta p_{T}/[p_{T}]; X", {HistType::kTH2D, {cfgAdditionObs.cfgDptDisAxisNormal, axisIndependent}});
369377
}
370378
if (doprocessMCGen) {
371379
registry.add("MCGen/MChPhi", "#phi distribution", {HistType::kTH1D, {axisPhi}});
372380
registry.add("MCGen/MChEta", "#eta distribution", {HistType::kTH1D, {axisEta}});
373-
registry.add("MCGen/MChPt", "p_{T} distribution after cut", {HistType::kTH1D, {axisPtHist}});
381+
registry.add("MCGen/MChPt", "p_{T} distribution after cut", {HistType::kTH1D, {axisPt}});
374382
registry.add("hMeanPtWithinGap08_MC", "mean p_{T}", {HistType::kTProfile, {axisIndependent}});
375383
for (auto i = 0; i < cfgNbootstrap; i++) {
376384
bootstrapArray[i][kMeanPtWithinGap08_MC] = registry.add<TProfile>(Form("BootstrapContainer_%d/hMeanPtWithinGap08_MC", i), "", {HistType::kTProfile, {axisIndependent}});
@@ -386,6 +394,17 @@ struct FlowTask {
386394
fWeights->setPtBins(nPtBins, ptBins);
387395
fWeights->init(true, false);
388396
}
397+
if (cfgAdditionObs.cfgV02Enabled) {
398+
cfgAdditionObs.nPt.resize(fPtAxis->GetNbins());
399+
cfgAdditionObs.listNptV2.resize(cfgNbootstrap + 1);
400+
cfgAdditionObs.listPtX.resize(cfgNbootstrap + 1);
401+
cfgAdditionObs.listNptV2[0] = registry.add<TProfile2D>(Form("V02/hNptV2"), "", {HistType::kTProfile2D, {axisIndependent, axisPt}});
402+
cfgAdditionObs.listPtX[0] = registry.add<TH2>(Form("V02/hPt_X"), "", {HistType::kTH2D, {axisIndependent, axisPt}});
403+
for (auto ibo = 0; ibo < cfgNbootstrap; ibo++) {
404+
cfgAdditionObs.listNptV2[ibo + 1] = registry.add<TProfile2D>(Form("V02/hNptV2_%d", ibo), "", {HistType::kTProfile2D, {axisIndependent, axisPt}});
405+
cfgAdditionObs.listPtX[ibo + 1] = registry.add<TH2>(Form("V02/hPt_X_%d", ibo), ";p_{T};X", {HistType::kTH2D, {axisIndependent, axisPt}});
406+
}
407+
}
389408

390409
// add in FlowContainer to Get boostrap sample automatically
391410
TObjArray* oba = new TObjArray();
@@ -548,6 +567,10 @@ struct FlowTask {
548567
corrconfigs.push_back(fGFW->GetCorrelatorConfig("refN10 {4 2} refP10 {-4 -2}", "Ch10Gap4242", kFALSE));
549568
corrconfigs.push_back(fGFW->GetCorrelatorConfig("refN10 {2 2} refP10 {-2 -2}", "Ch10Gap24", kFALSE));
550569
corrconfigs.push_back(fGFW->GetCorrelatorConfig("poiN10 refN10 | olN10 {2 2} refP10 {-2 -2}", "Ch10Gap24", kTRUE));
570+
LOGF(info, "Embedded GFW CorrelatorConfig:");
571+
for (auto icorr = 0; icorr < corrconfigs.size(); icorr++) {
572+
LOGF(info, "corrconfigs.at(%d): %s", icorr, corrconfigs.at(icorr).Head.c_str());
573+
}
551574
if (!userDefineGFWCorr.empty() && !userDefineGFWName.empty()) {
552575
LOGF(info, "User adding GFW CorrelatorConfig:");
553576
// attentaion: here we follow the index of cfgUserDefineGFWCorr
@@ -557,11 +580,11 @@ struct FlowTask {
557580
break;
558581
}
559582
if (userDefineGFWCorr.at(i).find("poi") != std::string::npos) {
560-
LOGF(info, "%d: enable pt-Diff for %s %s", i, userDefineGFWCorr.at(i).c_str(), userDefineGFWName.at(i).c_str());
561583
corrconfigs.push_back(fGFW->GetCorrelatorConfig(userDefineGFWCorr.at(i).c_str(), userDefineGFWName.at(i).c_str(), kTRUE));
584+
LOGF(info, "corrconfigs.at(%d): enable pt-Diff for %s %s", corrconfigs.size()-1, userDefineGFWCorr.at(i).c_str(), userDefineGFWName.at(i).c_str());
562585
} else {
563-
LOGF(info, "%d: %s %s", i, userDefineGFWCorr.at(i).c_str(), userDefineGFWName.at(i).c_str());
564586
corrconfigs.push_back(fGFW->GetCorrelatorConfig(userDefineGFWCorr.at(i).c_str(), userDefineGFWName.at(i).c_str(), kFALSE));
587+
LOGF(info, "corrconfigs.at(%d): %s %s", corrconfigs.size()-1, userDefineGFWCorr.at(i).c_str(), userDefineGFWName.at(i).c_str());
565588
}
566589
}
567590
}
@@ -702,6 +725,25 @@ struct FlowTask {
702725
return;
703726
}
704727

728+
void fillNptV2Profile(const GFW::CorrConfig& corrconf, const std::vector<float>& npt, const double& cent, const int& isample)
729+
{
730+
double dnx, val;
731+
dnx = fGFW->Calculate(corrconf, 0, kTRUE).real();
732+
if (dnx == 0)
733+
return;
734+
if (!corrconf.pTDif) {
735+
val = fGFW->Calculate(corrconf, 0, kFALSE).real() / dnx;
736+
if (std::fabs(val) < 1) {
737+
for (auto iPt = 0; iPt < fPtAxis->GetNbins(); iPt++) {
738+
auto pt = fPtAxis->GetBinCenter(iPt + 1);
739+
cfgAdditionObs.listNptV2[isample]->Fill(cent, pt, val * npt[iPt], dnx);
740+
}
741+
}
742+
return;
743+
}
744+
return;
745+
}
746+
705747
template <DataType dt>
706748
void fillFC(const GFW::CorrConfig& corrconf, const double& cent, const double& rndm)
707749
{
@@ -796,27 +838,27 @@ struct FlowTask {
796838
}
797839
LOGF(info, "Loaded efficiency histogram from %s (%p)", cfgUserIO.cfgEfficiencyForNch.value.c_str(), (void*)mEfficiencyForNch);
798840
}
799-
if (cfgFuncParas.cfgDptDisEnable && cfgFuncParas.cfgDptDishEvAvgMeanPt.value.empty() == false) {
800-
cfgFuncParas.hEvAvgMeanPt = ccdb->getForTimeStamp<TH1D>(cfgFuncParas.cfgDptDishEvAvgMeanPt, timestamp);
801-
if (cfgFuncParas.hEvAvgMeanPt == nullptr) {
802-
LOGF(fatal, "Could not load mean pT histogram from %s", cfgFuncParas.cfgDptDishEvAvgMeanPt.value.c_str());
841+
if (cfgAdditionObs.cfgDptDisEnable && cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.empty() == false) {
842+
cfgAdditionObs.hEvAvgMeanPt = ccdb->getForTimeStamp<TH1D>(cfgAdditionObs.cfgDptDishEvAvgMeanPt, timestamp);
843+
if (cfgAdditionObs.hEvAvgMeanPt == nullptr) {
844+
LOGF(fatal, "Could not load mean pT histogram from %s", cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.c_str());
803845
}
804-
LOGF(info, "Loaded mean pT histogram from %s (%p)", cfgFuncParas.cfgDptDishEvAvgMeanPt.value.c_str(), (void*)cfgFuncParas.hEvAvgMeanPt);
846+
LOGF(info, "Loaded mean pT histogram from %s (%p)", cfgAdditionObs.cfgDptDishEvAvgMeanPt.value.c_str(), (void*)cfgAdditionObs.hEvAvgMeanPt);
805847
}
806-
if (cfgFuncParas.cfgDptDisEnable && cfgFuncParas.cfgDptDisSelectionSwitch > kNoDptCut) {
807-
if (cfgFuncParas.cfgDptDisCutLow.value.empty() == false) {
808-
cfgFuncParas.fDptDisCutLow = ccdb->getForTimeStamp<TH1D>(cfgFuncParas.cfgDptDisCutLow, timestamp);
809-
if (cfgFuncParas.fDptDisCutLow == nullptr) {
810-
LOGF(fatal, "Could not load dptDis low cut histogram from %s", cfgFuncParas.cfgDptDisCutLow.value.c_str());
848+
if (cfgAdditionObs.cfgDptDisEnable && cfgAdditionObs.cfgDptDisSelectionSwitch > kNoDptCut) {
849+
if (cfgAdditionObs.cfgDptDisCutLow.value.empty() == false) {
850+
cfgAdditionObs.fDptDisCutLow = ccdb->getForTimeStamp<TH1D>(cfgAdditionObs.cfgDptDisCutLow, timestamp);
851+
if (cfgAdditionObs.fDptDisCutLow == nullptr) {
852+
LOGF(fatal, "Could not load dptDis low cut histogram from %s", cfgAdditionObs.cfgDptDisCutLow.value.c_str());
811853
}
812-
LOGF(info, "Loaded dptDis low cut histogram from %s (%p)", cfgFuncParas.cfgDptDisCutLow.value.c_str(), (void*)cfgFuncParas.fDptDisCutLow);
854+
LOGF(info, "Loaded dptDis low cut histogram from %s (%p)", cfgAdditionObs.cfgDptDisCutLow.value.c_str(), (void*)cfgAdditionObs.fDptDisCutLow);
813855
}
814-
if (cfgFuncParas.cfgDptDisCutHigh.value.empty() == false) {
815-
cfgFuncParas.fDptDisCutHigh = ccdb->getForTimeStamp<TH1D>(cfgFuncParas.cfgDptDisCutHigh, timestamp);
816-
if (cfgFuncParas.fDptDisCutHigh == nullptr) {
817-
LOGF(fatal, "Could not load dptDis high cut histogram from %s", cfgFuncParas.cfgDptDisCutHigh.value.c_str());
856+
if (cfgAdditionObs.cfgDptDisCutHigh.value.empty() == false) {
857+
cfgAdditionObs.fDptDisCutHigh = ccdb->getForTimeStamp<TH1D>(cfgAdditionObs.cfgDptDisCutHigh, timestamp);
858+
if (cfgAdditionObs.fDptDisCutHigh == nullptr) {
859+
LOGF(fatal, "Could not load dptDis high cut histogram from %s", cfgAdditionObs.cfgDptDisCutHigh.value.c_str());
818860
}
819-
LOGF(info, "Loaded dptDis high cut histogram from %s (%p)", cfgFuncParas.cfgDptDisCutHigh.value.c_str(), (void*)cfgFuncParas.fDptDisCutHigh);
861+
LOGF(info, "Loaded dptDis high cut histogram from %s (%p)", cfgAdditionObs.cfgDptDisCutHigh.value.c_str(), (void*)cfgAdditionObs.fDptDisCutHigh);
820862
}
821863
}
822864
correctionsLoaded = true;
@@ -1090,12 +1132,14 @@ struct FlowTask {
10901132
return;
10911133
registry.fill(HIST("hEventCount"), 3.5);
10921134
float lRandom = fRndm->Rndm();
1135+
int sampleIndex = static_cast<int>(cfgNbootstrap * lRandom);
10931136
float vtxz = collision.posZ();
10941137
registry.fill(HIST("hVtxZ"), vtxz);
10951138
registry.fill(HIST("hMult"), tracks.size());
10961139
registry.fill(HIST("hCent"), cent);
10971140
fGFW->Clear();
10981141
fFCpt->clearVector();
1142+
cfgAdditionObs.nPt.clear();
10991143
if (cfgGetInteractionRate) {
11001144
initHadronicRate(bc);
11011145
double hadronicRate = mRateFetcher.fetch(ccdb.service, bc.timestamp(), mRunNumber, "ZNC hadronic") * 1.e-3; //
@@ -1201,6 +1245,12 @@ struct FlowTask {
12011245
}
12021246
}
12031247
registry.fill(HIST("hPt"), track.pt());
1248+
if (cfgAdditionObs.cfgV02Enabled) {
1249+
cfgAdditionObs.listPtX[0]->Fill(independent, track.pt(), weff);
1250+
cfgAdditionObs.listPtX[sampleIndex + 1]->Fill(independent, track.pt(), weff);
1251+
if ((fPtAxis->FindBin(track.pt()) - 1) < fPtAxis->GetNbins()) // sanity check to avoid out of range error
1252+
cfgAdditionObs.nPt[fPtAxis->FindBin(track.pt()) - 1] += weff;
1253+
}
12041254
if (withinPtRef) {
12051255
registry.fill(HIST("hPhi"), track.phi());
12061256
registry.fill(HIST("hPhiWeighted"), track.phi(), wacc);
@@ -1240,15 +1290,15 @@ struct FlowTask {
12401290
independent = nTracksCorrected;
12411291
}
12421292

1243-
if (cfgFuncParas.cfgDptDisEnable) {
1293+
if (cfgAdditionObs.cfgDptDisEnable) {
12441294
double meanPt = ptSum / weffEvent;
1245-
double deltaPt = meanPt - cfgFuncParas.hEvAvgMeanPt->GetBinContent(cfgFuncParas.hEvAvgMeanPt->FindBin(independent));
1295+
double deltaPt = meanPt - cfgAdditionObs.hEvAvgMeanPt->GetBinContent(cfgAdditionObs.hEvAvgMeanPt->FindBin(independent));
12461296
double normDeltaPt = deltaPt / meanPt;
12471297
registry.fill(HIST("hNormDeltaPt_X"), normDeltaPt, independent, weffEvent);
1248-
if (cfgFuncParas.cfgDptDisSelectionSwitch == kLowDptCut && normDeltaPt > cfgFuncParas.fDptDisCutLow->GetBinContent(cfgFuncParas.fDptDisCutLow->FindBin(independent))) {
1298+
if (cfgAdditionObs.cfgDptDisSelectionSwitch == kLowDptCut && normDeltaPt > cfgAdditionObs.fDptDisCutLow->GetBinContent(cfgAdditionObs.fDptDisCutLow->FindBin(independent))) {
12491299
// only keep low 10% dpt event
12501300
return;
1251-
} else if (cfgFuncParas.cfgDptDisSelectionSwitch == kHighDptCut && normDeltaPt < cfgFuncParas.fDptDisCutHigh->GetBinContent(cfgFuncParas.fDptDisCutHigh->FindBin(independent))) {
1301+
} else if (cfgAdditionObs.cfgDptDisSelectionSwitch == kHighDptCut && normDeltaPt < cfgAdditionObs.fDptDisCutHigh->GetBinContent(cfgAdditionObs.fDptDisCutHigh->FindBin(independent))) {
12521302
// only keep high 10% dpt event
12531303
return;
12541304
}
@@ -1261,7 +1311,6 @@ struct FlowTask {
12611311
}
12621312
if (weffEventWithinGap08)
12631313
registry.fill(HIST("hMeanPtWithinGap08"), independent, ptSum_Gap08 / weffEventWithinGap08, 1.0);
1264-
int sampleIndex = static_cast<int>(cfgNbootstrap * lRandom);
12651314
if (weffEventWithinGap08)
12661315
bootstrapArray[sampleIndex][kMeanPtWithinGap08]->Fill(independent, ptSum_Gap08 / weffEventWithinGap08, 1.0);
12671316
// c22_gap8 * pt_withGap8
@@ -1277,6 +1326,15 @@ struct FlowTask {
12771326
weffEventDiffWithGap08);
12781327
}
12791328

1329+
// v02
1330+
if (cfgAdditionObs.cfgV02Enabled && weffEvent > 0) {
1331+
for (auto ipt = 0; ipt < fPtAxis->GetNbins(); ipt++) {
1332+
cfgAdditionObs.nPt[ipt] /= weffEvent;
1333+
}
1334+
fillNptV2Profile(corrconfigs.at(cfgAdditionObs.cfgV02CorrIndex), cfgAdditionObs.nPt, independent, 0);
1335+
fillNptV2Profile(corrconfigs.at(cfgAdditionObs.cfgV02CorrIndex), cfgAdditionObs.nPt, independent, sampleIndex + 1);
1336+
}
1337+
12801338
// Filling Flow Container
12811339
for (uint l_ind = 0; l_ind < corrconfigs.size(); l_ind++) {
12821340
fillFC<kReco>(corrconfigs.at(l_ind), independent, lRandom);

0 commit comments

Comments
 (0)