Skip to content

Commit 4cb452b

Browse files
skundu692skundu692
committed
Add multi-dim ZDC recenter and new hist in spin correlation task
1 parent 2e545a3 commit 4cb452b

2 files changed

Lines changed: 288 additions & 17 deletions

File tree

‎PWGLF/TableProducer/Common/zdc2stagecalibration.cxx‎

Lines changed: 116 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -78,7 +78,9 @@ struct zdc2stagecalibration {
7878
Configurable<bool> deriveSpatialCalib{"deriveSpatialCalib", false, "store spatial regression moments after gain correction"};
7979
Configurable<bool> deriveQRecentering{"deriveQRecentering", false, "store Q-vector recentering regression moments after gain and spatial correction"};
8080
Configurable<bool> useQRecentering{"useQRecentering", false, "apply Q-vector recentering calibration"};
81-
Configurable<std::string> confQRecenteringPath{"confQRecenteringPath", "", "CCDB path to the run-wise Q-vector recentering coefficients"};
81+
Configurable<int> qRecenteringMethod{"qRecenteringMethod", 0, "Q recentering method: 0=current regression, 1=multidimensional cubic regression"};
82+
Configurable<std::string> confQRecenteringPath{"confQRecenteringPath", "", "CCDB path to the current run-wise Q-vector recentering coefficients"};
83+
Configurable<std::string> confQRecenteringMultidimPath{"confQRecenteringMultidimPath", "", "CCDB path to the multidimensional cubic Q-vector recentering coefficients"};
8284

8385
Configurable<int> cfgSpatialCalibBins{
8486
"cfgSpatialCalibBins", 30,
@@ -103,10 +105,16 @@ struct zdc2stagecalibration {
103105
static constexpr int kQRecenteringNMatrixMoments = kQRecenteringNFeatures * (kQRecenteringNFeatures + 1) / 2;
104106
static constexpr int kQRecenteringNComponents = 4;
105107
static constexpr int kQRecenteringNMoments = kQRecenteringNMatrixMoments + kQRecenteringNComponents * kQRecenteringNFeatures;
108+
static constexpr int kQRecenteringMultidimNFeatures = 35;
109+
static constexpr int kQRecenteringMultidimNMatrixMoments = kQRecenteringMultidimNFeatures * (kQRecenteringMultidimNFeatures + 1) / 2;
110+
static constexpr int kQRecenteringMultidimNMoments = kQRecenteringMultidimNMatrixMoments + kQRecenteringNComponents * kQRecenteringMultidimNFeatures;
106111

107112
void init(o2::framework::InitContext&)
108113
{
109114
rctChecker.init(rctCut.cfgEvtRCTFlagCheckerLabel, rctCut.cfgEvtRCTFlagCheckerZDCCheck, rctCut.cfgEvtRCTFlagCheckerLimitAcceptAsBad);
115+
if (qRecenteringMethod.value < 0 || qRecenteringMethod.value > 1) {
116+
LOGF(fatal, "Unsupported qRecenteringMethod=%d. Use 0 for the current regression or 1 for the multidimensional cubic regression.", qRecenteringMethod.value);
117+
}
110118

111119
const int nTimeBins = static_cast<int>(std::ceil(cfgMaxRunHours.value * 60.f / cfgTimeSliceMinutes.value));
112120
AxisSpec timeAxis = {nTimeBins, 0.0, cfgMaxRunHours.value, "time from SOR (h)"};
@@ -119,6 +127,7 @@ struct zdc2stagecalibration {
119127
AxisSpec crossCoordAxis = {2, 0.0, 2.0, "coordinate"};
120128
AxisSpec ratioAxis = {240, 0.0, 2.4, "#Sigma tower/common"};
121129
AxisSpec qRecenteringMomentAxis = {kQRecenteringNMoments, 0.0, static_cast<double>(kQRecenteringNMoments), "regression moment"};
130+
AxisSpec qRecenteringMultidimMomentAxis = {kQRecenteringMultidimNMoments, 0.0, static_cast<double>(kQRecenteringMultidimNMoments), "multidimensional regression moment"};
122131
AxisSpec centralityAxis = {80, 0.0, 80.0, "centrality (%)"};
123132
AxisSpec vertexXYAxis = {100, -0.5, 0.5, "vertex x/y (cm)"};
124133
AxisSpec vertexZAxis = {100, -10.0, 10.0, "vertex z (cm)"};
@@ -187,6 +196,7 @@ struct zdc2stagecalibration {
187196
// For each 1% centrality bin, moments 0..230 are the upper triangle of sum(X_i X_j),
188197
// followed by 4 blocks of 21 sum(X_i Q): QxA, QyA, QxC, QyC.
189198
histos.add("QRecenteringCalibration/hRegressionMoments", "Q recentering regression moments;centrality (%);moment index", kTH2D, {centralityAxis, qRecenteringMomentAxis});
199+
histos.add("QRecenteringCalibration/hRegressionMomentsMultidim", "Multidimensional cubic Q recentering regression moments;centrality (%);moment index", kTH2D, {centralityAxis, qRecenteringMultidimMomentAxis});
190200

191201
histos.add("QRecenteringQA/pQBeforeVsCentrality", "Q before recentering vs centrality;centrality (%);Q component;<Q>", kTProfile2D, {centralityAxis, qComponentAxis});
192202
histos.add("QRecenteringQA/pQBeforeVsTime", "Q before recentering vs time;time from SOR (h);Q component;<Q>", kTProfile2D, {timeAxis, qComponentAxis});
@@ -198,6 +208,22 @@ struct zdc2stagecalibration {
198208
histos.add("QRecenteringQA/pQAfterVsVx", "Q after recentering vs v_{x};v_{x} (cm);Q component;<Q>", kTProfile2D, {vertexXYAxis, qComponentAxis});
199209
histos.add("QRecenteringQA/pQAfterVsVy", "Q after recentering vs v_{y};v_{y} (cm);Q component;<Q>", kTProfile2D, {vertexXYAxis, qComponentAxis});
200210
histos.add("QRecenteringQA/pQAfterVsVz", "Q after recentering vs v_{z};v_{z} (cm);Q component;<Q>", kTProfile2D, {vertexZAxis, qComponentAxis});
211+
histos.add("QRecentering2DQA/pQxAAfterVsCentralityTime", "QxA after recentering;centrality (%);time from SOR (h);<Q_{x}^{A}>", kTProfile2D, {centralityAxis, timeAxis});
212+
histos.add("QRecentering2DQA/pQyAAfterVsCentralityTime", "QyA after recentering;centrality (%);time from SOR (h);<Q_{y}^{A}>", kTProfile2D, {centralityAxis, timeAxis});
213+
histos.add("QRecentering2DQA/pQxCAfterVsCentralityTime", "QxC after recentering;centrality (%);time from SOR (h);<Q_{x}^{C}>", kTProfile2D, {centralityAxis, timeAxis});
214+
histos.add("QRecentering2DQA/pQyCAfterVsCentralityTime", "QyC after recentering;centrality (%);time from SOR (h);<Q_{y}^{C}>", kTProfile2D, {centralityAxis, timeAxis});
215+
histos.add("QRecentering2DQA/pQxAAfterVsCentralityVx", "QxA after recentering;centrality (%);v_{x} (cm);<Q_{x}^{A}>", kTProfile2D, {centralityAxis, vertexXYAxis});
216+
histos.add("QRecentering2DQA/pQyAAfterVsCentralityVx", "QyA after recentering;centrality (%);v_{x} (cm);<Q_{y}^{A}>", kTProfile2D, {centralityAxis, vertexXYAxis});
217+
histos.add("QRecentering2DQA/pQxCAfterVsCentralityVx", "QxC after recentering;centrality (%);v_{x} (cm);<Q_{x}^{C}>", kTProfile2D, {centralityAxis, vertexXYAxis});
218+
histos.add("QRecentering2DQA/pQyCAfterVsCentralityVx", "QyC after recentering;centrality (%);v_{x} (cm);<Q_{y}^{C}>", kTProfile2D, {centralityAxis, vertexXYAxis});
219+
histos.add("QRecentering2DQA/pQxAAfterVsCentralityVy", "QxA after recentering;centrality (%);v_{y} (cm);<Q_{x}^{A}>", kTProfile2D, {centralityAxis, vertexXYAxis});
220+
histos.add("QRecentering2DQA/pQyAAfterVsCentralityVy", "QyA after recentering;centrality (%);v_{y} (cm);<Q_{y}^{A}>", kTProfile2D, {centralityAxis, vertexXYAxis});
221+
histos.add("QRecentering2DQA/pQxCAfterVsCentralityVy", "QxC after recentering;centrality (%);v_{y} (cm);<Q_{x}^{C}>", kTProfile2D, {centralityAxis, vertexXYAxis});
222+
histos.add("QRecentering2DQA/pQyCAfterVsCentralityVy", "QyC after recentering;centrality (%);v_{y} (cm);<Q_{y}^{C}>", kTProfile2D, {centralityAxis, vertexXYAxis});
223+
histos.add("QRecentering2DQA/pQxAAfterVsCentralityVz", "QxA after recentering;centrality (%);v_{z} (cm);<Q_{x}^{A}>", kTProfile2D, {centralityAxis, vertexZAxis});
224+
histos.add("QRecentering2DQA/pQyAAfterVsCentralityVz", "QyA after recentering;centrality (%);v_{z} (cm);<Q_{y}^{A}>", kTProfile2D, {centralityAxis, vertexZAxis});
225+
histos.add("QRecentering2DQA/pQxCAfterVsCentralityVz", "QxC after recentering;centrality (%);v_{z} (cm);<Q_{x}^{C}>", kTProfile2D, {centralityAxis, vertexZAxis});
226+
histos.add("QRecentering2DQA/pQyCAfterVsCentralityVz", "QyC after recentering;centrality (%);v_{z} (cm);<Q_{y}^{C}>", kTProfile2D, {centralityAxis, vertexZAxis});
201227

202228
// A-C Q-vector correlation QA. These histograms are diagnostic only; no Qx/Qy mixing is applied.
203229
histos.add("QCorrelationQA/hQxAQxCBefore", "QxA vs QxC before recentering;Q_{x}^{A};Q_{x}^{C}", kTH2F, {qValueAxis, qValueAxis});
@@ -252,8 +278,11 @@ struct zdc2stagecalibration {
252278
LOGF(warn, "No ZDC spatial calibration found for run %d at timestamp %llu", runNumber, static_cast<unsigned long long>(timestamp));
253279
}
254280
}
255-
if (calibrationStage.value == 2 && useQRecentering.value && !confQRecenteringPath.value.empty()) {
256-
qRecenteringProfile = ccdb->getForTimeStamp<TH3D>(confQRecenteringPath.value, timestamp);
281+
if (calibrationStage.value == 2 && useQRecentering.value) {
282+
const std::string qRecenteringPath = qRecenteringMethod.value == 1 ? confQRecenteringMultidimPath.value : confQRecenteringPath.value;
283+
if (!qRecenteringPath.empty()) {
284+
qRecenteringProfile = ccdb->getForTimeStamp<TH3D>(qRecenteringPath, timestamp);
285+
}
257286
if (!qRecenteringProfile) {
258287
LOGF(warn, "No ZDC Q-recentering calibration found for run %d at timestamp %llu", runNumber, static_cast<unsigned long long>(timestamp));
259288
}
@@ -307,6 +336,26 @@ struct zdc2stagecalibration {
307336
x * y, x * z, y * z};
308337
}
309338

339+
std::array<double, kQRecenteringMultidimNFeatures> makeQRecenteringMultidimFeatures(float timeFromSOR, float vx, float vy, float vz) const
340+
{
341+
const double runHours = (eorTimestamp > sorTimestamp) ? static_cast<double>(eorTimestamp - sorTimestamp) * 1.e-3 / 3600.0 : static_cast<double>(cfgMaxRunHours.value);
342+
const double t = (runHours > 0.0) ? (2.0 * static_cast<double>(timeFromSOR) / runHours - 1.0) : 0.0;
343+
const double x = static_cast<double>(vx) / 0.1;
344+
const double y = static_cast<double>(vy) / 0.1;
345+
const double z = static_cast<double>(vz) / 10.0;
346+
const double t2 = t * t;
347+
const double x2 = x * x;
348+
const double y2 = y * y;
349+
const double z2 = z * z;
350+
351+
return {1.0, t, x, y, z,
352+
t2, x2, y2, z2, t * x, t * y, t * z, x * y, x * z, y * z,
353+
t2 * t, x2 * x, y2 * y, z2 * z,
354+
t2 * x, t2 * y, t2 * z, x2 * t, x2 * y, x2 * z,
355+
y2 * t, y2 * x, y2 * z, z2 * t, z2 * x, z2 * y,
356+
t * x * y, t * x * z, t * y * z, x * y * z};
357+
}
358+
310359
using MyCollisions = o2::soa::Join<o2::aod::Collisions, o2::aod::EvSels, o2::aod::Mults, o2::aod::FT0sCorrected, o2::aod::CentFT0Cs>;
311360

312361
void process(MyCollisions::iterator const& collision, o2::aod::FT0s const& /*ft0s*/, o2::aod::FV0As const& /*fv0s*/, BCsRun3 const& /*bcs*/, o2::aod::Zdcs const& /*zdcs*/)
@@ -632,7 +681,6 @@ struct zdc2stagecalibration {
632681
return;
633682
}
634683

635-
const auto qFeatures = makeQRecenteringFeatures(centrality, timeFromSOR, vx, vy, vz);
636684
std::array<double, 4> qValues = {qxZDCA, qyZDCA, qxZDCC, qyZDCC};
637685

638686
for (int i = 0; i < 4; ++i) {
@@ -654,17 +702,35 @@ struct zdc2stagecalibration {
654702
histos.fill(HIST("QCorrelationQA/pQACorrelationBeforeVsCentrality"), centrality, 3.5, qValues[1] * qValues[2]);
655703

656704
if (calibrationStage.value == 2 && useGainCallib.value && useSpatialCalib.value && deriveQRecentering.value && !useQRecentering.value) {
657-
int moment = 0;
658-
for (int i = 0; i < kQRecenteringNFeatures; ++i) {
659-
for (int j = i; j < kQRecenteringNFeatures; ++j) {
660-
histos.fill(HIST("QRecenteringCalibration/hRegressionMoments"), centrality, moment + 0.5, qFeatures[i] * qFeatures[j]);
661-
++moment;
662-
}
663-
}
664-
for (int component = 0; component < kQRecenteringNComponents; ++component) {
705+
if (qRecenteringMethod.value == 0) {
706+
const auto qFeatures = makeQRecenteringFeatures(centrality, timeFromSOR, vx, vy, vz);
707+
int moment = 0;
665708
for (int i = 0; i < kQRecenteringNFeatures; ++i) {
666-
const int index = kQRecenteringNMatrixMoments + component * kQRecenteringNFeatures + i;
667-
histos.fill(HIST("QRecenteringCalibration/hRegressionMoments"), centrality, index + 0.5, qFeatures[i] * qValues[component]);
709+
for (int j = i; j < kQRecenteringNFeatures; ++j) {
710+
histos.fill(HIST("QRecenteringCalibration/hRegressionMoments"), centrality, moment + 0.5, qFeatures[i] * qFeatures[j]);
711+
++moment;
712+
}
713+
}
714+
for (int component = 0; component < kQRecenteringNComponents; ++component) {
715+
for (int i = 0; i < kQRecenteringNFeatures; ++i) {
716+
const int index = kQRecenteringNMatrixMoments + component * kQRecenteringNFeatures + i;
717+
histos.fill(HIST("QRecenteringCalibration/hRegressionMoments"), centrality, index + 0.5, qFeatures[i] * qValues[component]);
718+
}
719+
}
720+
} else {
721+
const auto qFeatures = makeQRecenteringMultidimFeatures(timeFromSOR, vx, vy, vz);
722+
int moment = 0;
723+
for (int i = 0; i < kQRecenteringMultidimNFeatures; ++i) {
724+
for (int j = i; j < kQRecenteringMultidimNFeatures; ++j) {
725+
histos.fill(HIST("QRecenteringCalibration/hRegressionMomentsMultidim"), centrality, moment + 0.5, qFeatures[i] * qFeatures[j]);
726+
++moment;
727+
}
728+
}
729+
for (int component = 0; component < kQRecenteringNComponents; ++component) {
730+
for (int i = 0; i < kQRecenteringMultidimNFeatures; ++i) {
731+
const int index = kQRecenteringMultidimNMatrixMoments + component * kQRecenteringMultidimNFeatures + i;
732+
histos.fill(HIST("QRecenteringCalibration/hRegressionMomentsMultidim"), centrality, index + 0.5, qFeatures[i] * qValues[component]);
733+
}
668734
}
669735
}
670736
}
@@ -673,7 +739,8 @@ struct zdc2stagecalibration {
673739
if (calibrationStage.value == 2 && useQRecentering.value) {
674740
if (!qRecenteringProfile) {
675741
qRecenteringOK = false;
676-
} else {
742+
} else if (qRecenteringMethod.value == 0) {
743+
const auto qFeatures = makeQRecenteringFeatures(centrality, timeFromSOR, vx, vy, vz);
677744
for (int component = 0; component < kQRecenteringNComponents; ++component) {
678745
double predictedBias = 0.0;
679746
for (int i = 0; i < kQRecenteringNFeatures; ++i) {
@@ -690,6 +757,24 @@ struct zdc2stagecalibration {
690757
}
691758
qValues[component] -= predictedBias;
692759
}
760+
} else {
761+
const auto qFeatures = makeQRecenteringMultidimFeatures(timeFromSOR, vx, vy, vz);
762+
for (int component = 0; component < kQRecenteringNComponents; ++component) {
763+
double predictedBias = 0.0;
764+
for (int i = 0; i < kQRecenteringMultidimNFeatures; ++i) {
765+
const double coefficient = qRecenteringProfile->GetBinContent(qRecenteringProfile->FindBin(centrality, i + 0.5, component + 0.5));
766+
if (!std::isfinite(coefficient)) {
767+
qRecenteringOK = false;
768+
break;
769+
}
770+
predictedBias += coefficient * qFeatures[i];
771+
}
772+
if (!qRecenteringOK || !std::isfinite(predictedBias)) {
773+
qRecenteringOK = false;
774+
break;
775+
}
776+
qValues[component] -= predictedBias;
777+
}
693778
}
694779
}
695780
if (!qRecenteringOK) {
@@ -705,6 +790,22 @@ struct zdc2stagecalibration {
705790
histos.fill(HIST("QRecenteringQA/pQAfterVsVy"), vy, component, qValues[i]);
706791
histos.fill(HIST("QRecenteringQA/pQAfterVsVz"), vz, component, qValues[i]);
707792
}
793+
histos.fill(HIST("QRecentering2DQA/pQxAAfterVsCentralityTime"), centrality, timeFromSOR, qValues[0]);
794+
histos.fill(HIST("QRecentering2DQA/pQyAAfterVsCentralityTime"), centrality, timeFromSOR, qValues[1]);
795+
histos.fill(HIST("QRecentering2DQA/pQxCAfterVsCentralityTime"), centrality, timeFromSOR, qValues[2]);
796+
histos.fill(HIST("QRecentering2DQA/pQyCAfterVsCentralityTime"), centrality, timeFromSOR, qValues[3]);
797+
histos.fill(HIST("QRecentering2DQA/pQxAAfterVsCentralityVx"), centrality, vx, qValues[0]);
798+
histos.fill(HIST("QRecentering2DQA/pQyAAfterVsCentralityVx"), centrality, vx, qValues[1]);
799+
histos.fill(HIST("QRecentering2DQA/pQxCAfterVsCentralityVx"), centrality, vx, qValues[2]);
800+
histos.fill(HIST("QRecentering2DQA/pQyCAfterVsCentralityVx"), centrality, vx, qValues[3]);
801+
histos.fill(HIST("QRecentering2DQA/pQxAAfterVsCentralityVy"), centrality, vy, qValues[0]);
802+
histos.fill(HIST("QRecentering2DQA/pQyAAfterVsCentralityVy"), centrality, vy, qValues[1]);
803+
histos.fill(HIST("QRecentering2DQA/pQxCAfterVsCentralityVy"), centrality, vy, qValues[2]);
804+
histos.fill(HIST("QRecentering2DQA/pQyCAfterVsCentralityVy"), centrality, vy, qValues[3]);
805+
histos.fill(HIST("QRecentering2DQA/pQxAAfterVsCentralityVz"), centrality, vz, qValues[0]);
806+
histos.fill(HIST("QRecentering2DQA/pQyAAfterVsCentralityVz"), centrality, vz, qValues[1]);
807+
histos.fill(HIST("QRecentering2DQA/pQxCAfterVsCentralityVz"), centrality, vz, qValues[2]);
808+
histos.fill(HIST("QRecentering2DQA/pQyCAfterVsCentralityVz"), centrality, vz, qValues[3]);
708809

709810
histos.fill(HIST("QCorrelationQA/hQxAQxCAfter"), qValues[0], qValues[2]);
710811
histos.fill(HIST("QCorrelationQA/hQyAQyCAfter"), qValues[1], qValues[3]);

0 commit comments

Comments
 (0)