Skip to content

Commit d2d4b91

Browse files
committed
mv Cascade AP-variables to Datamodel; bugfix qtarm eval
1 parent 395f55a commit d2d4b91

File tree

2 files changed

+63
-36
lines changed

2 files changed

+63
-36
lines changed

PWGDQ/Tasks/v0selector.cxx

Lines changed: 7 additions & 32 deletions
Original file line numberDiff line numberDiff line change
@@ -174,39 +174,13 @@ struct v0selector {
174174
return kUndef;
175175
}
176176

177-
std::pair<float, float> evalArmenterosPodolanskiCascade(aod::CascData const& casc)
178-
{
179-
const float pxv0 = casc.pxpos() + casc.pxneg();
180-
const float pyv0 = casc.pypos() + casc.pyneg();
181-
const float pzv0 = casc.pzpos() + casc.pzneg();
182-
183-
const float pxbach = casc.pxbach();
184-
const float pybach = casc.pybach();
185-
const float pzbach = casc.pzbach();
186-
187-
const double momTot = RecoDecay::p2(pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach);
188-
189-
const float lQlv0 = RecoDecay::dotProd(std::array{pxv0, pyv0, pzv0}, std::array{pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach}) / momTot;
190-
const float lQlbach = RecoDecay::dotProd(std::array{pxbach, pybach, pzbach}, std::array{pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach}) / momTot;
191-
const float dp = RecoDecay::dotProd(std::array{pxbach, pybach, pzbach}, std::array{pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach});
192-
193-
float alpha = (lQlv0 - lQlbach) / (lQlv0 + lQlbach);
194-
const float qtarm = std::sqrt(RecoDecay::p2(pxbach, pybach, pzbach) - dp * dp / momTot / momTot);
195-
196-
if (casc.sign() > 0) {
197-
alpha = -alpha;
198-
}
199-
200-
return std::make_pair(alpha, qtarm);
201-
}
202-
203177
int checkCascade(float alpha, float qt)
204178
{
205179
const bool isAlphaPos = alpha > 0;
206180
alpha = std::fabs(alpha);
207181

208-
const float qUp = cutAPOmegaUp1 * std::sqrt(std::abs(1.0f - ((alpha - cutAPOmegaUp2) * (alpha - cutAPOmegaUp2)) / (cutAPOmegaUp3 * cutAPOmegaUp3)));
209-
const float qDown = cutAPOmegaDown1 * std::sqrt(std::abs(1.0f - ((alpha - cutAPOmegaDown2) * (alpha - cutAPOmegaDown2)) / (cutAPOmegaDown3 * cutAPOmegaDown3)));
182+
const float qUp = std::abs(alpha - cutAPOmegaUp2) > std::abs(cutAPOmegaUp3) ? 0. : cutAPOmegaUp1 * std::sqrt(1.0f - ((alpha - cutAPOmegaUp2) * (alpha - cutAPOmegaUp2)) / (cutAPOmegaUp3 * cutAPOmegaUp3));
183+
const float qDown = std::abs(alpha - cutAPOmegaDown2) > std::abs(cutAPOmegaDown3) ? 0. : cutAPOmegaDown1 * std::sqrt(1.0f - ((alpha - cutAPOmegaDown2) * (alpha - cutAPOmegaDown2)) / (cutAPOmegaDown3 * cutAPOmegaDown3));
210184

211185
if (alpha < cutAlphaOmegaLow || alpha > cutAlphaOmegaHigh || qt < qDown || qt > qUp) {
212186
return kUndef;
@@ -262,8 +236,8 @@ struct v0selector {
262236
registry.add("hV0APplot", "hV0APplot", HistType::kTH2F, {{200, -1.0f, +1.0f}, {250, 0.0f, 0.25f}});
263237
registry.add("hV0APplotSelected", "hV0APplotSelected", HistType::kTH2F, {{200, -1.0f, +1.0f}, {250, 0.0f, 0.25f}});
264238
registry.add("hV0Psi", "hV0Psi", HistType::kTH2F, {{100, 0, TMath::PiOver2()}, {100, 0, 0.1}});
265-
registry.add("hCascAPplot", "hCascAPplot", HistType::kTH2F, {{200, -1.0f, +1.0f}, {500, 0.0f, 0.5f}});
266-
registry.add("hCascAPplotSelected", "hCascAPplotSelected", HistType::kTH2F, {{200, -1.0f, +1.0f}, {500, 0.0f, 0.5f}});
239+
registry.add("hCascAPplot", "hCascAPplot", HistType::kTH2F, {{200, -1.0f, +1.0f}, {300, 0.0f, 0.3f}});
240+
registry.add("hCascAPplotSelected", "hCascAPplotSelected", HistType::kTH2F, {{200, -1.0f, +1.0f}, {300, 0.0f, 0.3f}});
267241
registry.add("hMassOmega", "hMassOmega", HistType::kTH2F, {{900, 0.0f, 90.0f}, {100, 1.62f, 1.72f}});
268242
registry.add("hMassAntiOmega", "hMassAntiOmega", HistType::kTH2F, {{900, 0.0f, 90.0f}, {100, 1.62f, 1.72f}});
269243
}
@@ -283,7 +257,7 @@ struct v0selector {
283257
std::vector<int8_t> cascpidmap;
284258
cascpidmap.clear();
285259
if (produceCascID.value) {
286-
cascpidmap.resize(V0s.size(), kUndef);
260+
cascpidmap.resize(Cascs.size(), kUndef);
287261
}
288262
for (auto& V0 : V0s) {
289263
// if (!(V0.posTrack_as<FullTracksExt>().trackType() & o2::aod::track::TPCrefit)) {
@@ -467,7 +441,8 @@ struct v0selector {
467441

468442
const float mOmega = Casc.mOmega();
469443

470-
auto [alpha, qt] = evalArmenterosPodolanskiCascade(Casc);
444+
const float alpha = Casc.alpha();
445+
const float qt = Casc.qtarm();
471446

472447
if (fillhisto) {
473448
registry.fill(HIST("hCascAPplot"), alpha, qt);

PWGLF/DataModel/LFStrangenessTables.h

Lines changed: 56 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -713,21 +713,21 @@ DECLARE_SOA_DYNAMIC_COLUMN(MLambda, mLambda, //! mass under lambda hypothesis
713713
[](int v0type, float pxpos, float pypos, float pzpos, float pxneg, float pyneg, float pzneg) -> float {
714714
if ((v0type & (0x1 << o2::dataformats::V0Index::kPhotonOnly)) != 0) {
715715
return 0.0f; // provide mass only if NOT a photon with TPC-only tracks (special handling)
716-
};
716+
}
717717
return RecoDecay::m(std::array{std::array{pxpos, pypos, pzpos}, std::array{pxneg, pyneg, pzneg}}, std::array{o2::constants::physics::MassProton, o2::constants::physics::MassPionCharged});
718718
});
719719
DECLARE_SOA_DYNAMIC_COLUMN(MAntiLambda, mAntiLambda, //! mass under antilambda hypothesis
720720
[](int v0type, float pxpos, float pypos, float pzpos, float pxneg, float pyneg, float pzneg) -> float {
721721
if ((v0type & (0x1 << o2::dataformats::V0Index::kPhotonOnly)) != 0) {
722722
return 0.0f; // provide mass only if NOT a photon with TPC-only tracks (special handling)
723-
};
723+
}
724724
return RecoDecay::m(std::array{std::array{pxpos, pypos, pzpos}, std::array{pxneg, pyneg, pzneg}}, std::array{o2::constants::physics::MassPionCharged, o2::constants::physics::MassProton});
725725
});
726726
DECLARE_SOA_DYNAMIC_COLUMN(MK0Short, mK0Short, //! mass under K0short hypothesis
727727
[](int v0type, float pxpos, float pypos, float pzpos, float pxneg, float pyneg, float pzneg) -> float {
728728
if ((v0type & (0x1 << o2::dataformats::V0Index::kPhotonOnly)) != 0) {
729729
return 0.0f; // provide mass only if NOT a photon with TPC-only tracks (special handling)
730-
};
730+
}
731731
return RecoDecay::m(std::array{std::array{pxpos, pypos, pzpos}, std::array{pxneg, pyneg, pzneg}}, std::array{o2::constants::physics::MassPionCharged, o2::constants::physics::MassPionCharged});
732732
});
733733
DECLARE_SOA_DYNAMIC_COLUMN(MLambda_unchecked, mLambda_unchecked, //! mass under lambda hypothesis without v0 type check (will include TPC only and potentially duplicates! use with care)
@@ -1417,6 +1417,53 @@ DECLARE_SOA_DYNAMIC_COLUMN(Phi, phi, //! Cascade phi in the range [0, 2pi)
14171417
[](float px, float py) -> float { return RecoDecay::phi(px, py); });
14181418
DECLARE_SOA_DYNAMIC_COLUMN(Eta, eta, //! Cascade pseudorapidity
14191419
[](float px, float py, float pz) -> float { return RecoDecay::eta(std::array{px, py, pz}); });
1420+
1421+
// Armenteros-Podolanski variables
1422+
DECLARE_SOA_DYNAMIC_COLUMN(Alpha, alpha, //! Armenteros Alpha
1423+
[](float pxpos, float pypos, float pzpos, float pxneg, float pyneg, float pzneg, float pxbach, float pybach, float pzbach, int sign) {
1424+
const float pxv0 = pxpos + pxneg;
1425+
const float pyv0 = pypos + pyneg;
1426+
const float pzv0 = pzpos + pzneg;
1427+
1428+
// No need to divide by momentum of the cascade (as in the v0data namespace) since the ratio of lQl is evaluated
1429+
const float lQlBach = RecoDecay::dotProd(std::array{pxbach, pybach, pzbach}, std::array{pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach});
1430+
const float lQlV0 = RecoDecay::dotProd(std::array{pxv0, pyv0, pzv0}, std::array{pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach});
1431+
float alpha = (lQlBach - lQlV0) / (lQlV0 + lQlBach); // alphacascade
1432+
if (sign < 0) {
1433+
alpha = -alpha;
1434+
}
1435+
return alpha;
1436+
});
1437+
1438+
DECLARE_SOA_DYNAMIC_COLUMN(QtArm, qtarm, //! Armenteros Qt
1439+
[](float pxpos, float pypos, float pzpos, float pxneg, float pyneg, float pzneg, float pxbach, float pybach, float pzbach) {
1440+
const float pxv0 = pxpos + pxneg;
1441+
const float pyv0 = pypos + pyneg;
1442+
const float pzv0 = pzpos + pzneg;
1443+
1444+
const float momTot2 = RecoDecay::p2(pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach);
1445+
const float dp = RecoDecay::dotProd(std::array{pxbach, pybach, pzbach}, std::array{pxv0 + pxbach, pyv0 + pybach, pzv0 + pzbach});
1446+
return std::sqrt(RecoDecay::p2(pxbach, pybach, pzbach) - dp * dp / momTot2); // qtarm
1447+
});
1448+
1449+
// Psi pair angle: angle between the plane defined by the v0 and bachelor momenta and the xy plane
1450+
DECLARE_SOA_DYNAMIC_COLUMN(PsiPair, psipair, //! psi pair angle
1451+
[](float pxpos, float pypos, float pzpos, float pxneg, float pyneg, float pzneg, float pxbach, float pybach, float pzbach, int sign) {
1452+
auto clipToPM1 = [](float x) { return x < -1.f ? -1.f : (x > 1.f ? 1.f : x); };
1453+
const float pxv0 = pxpos + pxneg;
1454+
const float pyv0 = pypos + pyneg;
1455+
const float pzv0 = pzpos + pzneg;
1456+
1457+
const float ptot2 = RecoDecay::p2(pxbach, pybach, pzbach) * RecoDecay::p2(pxv0, pyv0, pzv0);
1458+
const float argcos = RecoDecay::dotProd(std::array{pxbach, pybach, pzbach}, std::array{pxv0, pyv0, pzv0}) / std::sqrt(ptot2);
1459+
const float thetaV0 = std::atan2(RecoDecay::sqrtSumOfSquares(pxv0, pyv0), pzv0);
1460+
const float thetaBach = std::atan2(RecoDecay::sqrtSumOfSquares(pxbach, pybach), pzbach);
1461+
float argsin = (thetaV0 - thetaBach) / std::acos(clipToPM1(argcos));
1462+
if (sign < 0) {
1463+
argsin = -argsin;
1464+
}
1465+
return std::asin(clipToPM1(argsin));
1466+
});
14201467
} // namespace cascdata
14211468

14221469
//______________________________________________________
@@ -1493,7 +1540,12 @@ DECLARE_SOA_TABLE(StoredCascCores, "AOD", "CASCCORE", //! core information about
14931540
cascdata::PositiveEta<cascdata::PxPos, cascdata::PyPos, cascdata::PzPos>,
14941541
cascdata::PositivePhi<cascdata::PxPos, cascdata::PyPos>,
14951542
cascdata::BachelorEta<cascdata::PxBach, cascdata::PyBach, cascdata::PzBach>,
1496-
cascdata::BachelorPhi<cascdata::PxBach, cascdata::PyBach>);
1543+
cascdata::BachelorPhi<cascdata::PxBach, cascdata::PyBach>,
1544+
1545+
// Armenteros-Podolanski
1546+
cascdata::Alpha<cascdata::PxPos, cascdata::PyPos, cascdata::PzPos, cascdata::PxNeg, cascdata::PyNeg, cascdata::PzNeg, cascdata::PxBach, cascdata::PyBach, cascdata::PzBach, cascdata::Sign>,
1547+
cascdata::QtArm<cascdata::PxPos, cascdata::PyPos, cascdata::PzPos, cascdata::PxNeg, cascdata::PyNeg, cascdata::PzNeg, cascdata::PxBach, cascdata::PyBach, cascdata::PzBach>,
1548+
cascdata::PsiPair<cascdata::PxPos, cascdata::PyPos, cascdata::PzPos, cascdata::PxNeg, cascdata::PyNeg, cascdata::PzNeg, cascdata::PxBach, cascdata::PyBach, cascdata::PzBach, cascdata::Sign>);
14971549

14981550
DECLARE_SOA_TABLE(StoredKFCascCores, "AOD", "KFCASCCORE", //!
14991551
cascdata::Sign, cascdata::MXi, cascdata::MOmega,

0 commit comments

Comments
 (0)