Skip to content

[PWGDQ] Add mix event for coherent Jpsi A2 #17033

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
ypwangg wants to merge 9 commits into
base: master
Choose a base branch
from
+99 −2
Open

[PWGDQ] Add mix event for coherent Jpsi A2 #17033

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
0f1bacf
add A2 in MC
ypwangg Jul 3, 2026
ab4fa43
fix bug
ypwangg Jul 3, 2026
6617cfd
fix bug in VM daughter direction
ypwangg Jul 3, 2026
4d53d5f
mix event for coherent jpsi A2
ypwangg Jul 13, 2026
c19f3c3
Please consider the following formatting changes
alibuild Jul 13, 2026
133f9c1
Merge branch 'master' into A2-PR
ypwangg Jul 13, 2026
08c3805
Merge pull request #9 from alibuild/alibot-cleanup-17033
ypwangg Jul 13, 2026
c203e07
Please consider the following formatting changes
alibuild Jul 13, 2026
373e503
Merge pull request #10 from alibuild/alibot-cleanup-17033
ypwangg Jul 13, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
97 changes: 97 additions & 0 deletions PWGDQ/Core/VarManager.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -479,6 +479,14 @@
kQ2X0A2,
kQ2Y0A1,
kQ2Y0A2,
kQ2X0B1,
kQ2X0B2,
kQ2Y0B1,
kQ2Y0B2,
kQ2X0C1,
kQ2X0C2,
kQ2Y0C1,
kQ2Y0C2,
kU2Q2Ev1,
kU2Q2Ev2,
kCos2DeltaPhiEv1,
Expand Down Expand Up @@ -1510,7 +1518,7 @@

static void SetCalibrationType(int type, bool useInterpolation = true)
{
if (type < 0 || type > 2) {

Check failure on line 1521 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
LOG(fatal) << "Invalid calibration type. Must be 0, 1, or 2.";
}
fgCalibrationType = type;
Expand Down Expand Up @@ -1737,7 +1745,7 @@
{
o2::track::TrackParCovFwd fwdtrack = o2::aod::fwdtrackutils::getTrackParCovFwdShift(muon, fgzShiftFwd, muon);
o2::dataformats::GlobalFwdTrack propmuon;
if (static_cast<int>(muon.trackType()) > 2) {

Check failure on line 1748 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
o2::dataformats::GlobalFwdTrack track;
track.setParameters(fwdtrack.getParameters());
track.setZ(fwdtrack.getZ());
Expand All @@ -1762,7 +1770,7 @@
propmuon.setZ(proptrack.getZ());
propmuon.setCovariances(proptrack.getCovariances());

} else if (static_cast<int>(muon.trackType()) < 2) {

Check failure on line 1773 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
std::array<double, 3> dcaInfOrig{999.f, 999.f, 999.f};
fwdtrack.propagateToDCAhelix(fgMagField, {collision.posX(), collision.posY(), collision.posZ()}, dcaInfOrig);
propmuon.setParameters(fwdtrack.getParameters());
Expand Down Expand Up @@ -1824,7 +1832,7 @@

// Redo propagation only for muon tracks
// propagation of MFT tracks alredy done in fwdtrack-extention task
if (static_cast<int>(muon.trackType()) > 2) {

Check failure on line 1835 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
o2::dataformats::GlobalFwdTrack propmuonAtDCA = PropagateMuon(muon, collision, kToDCA);
o2::dataformats::GlobalFwdTrack propmuonAtRabs = PropagateMuon(muon, collision, kToRabs);
float dcaX = (propmuonAtDCA.getX() - collision.posX());
Expand Down Expand Up @@ -2625,7 +2633,7 @@
if (!track.hasTPC()) {
continue; // skip tracks without TPC information
}
if (track.dcaZ() > 998) {

Check failure on line 2636 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
continue; // skip tracks without valid DCAz
}
dcazValues.push_back(track.dcaZ());
Expand Down Expand Up @@ -2713,15 +2721,15 @@
int counter10mm = 0;
for (auto const& d : dcazValues) {
double absD = std::abs(d);
if (absD > 0.01) {

Check failure on line 2724 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
counter100um++;
if (absD > 0.02) {

Check failure on line 2726 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
counter200um++;
if (absD > 0.05) {

Check failure on line 2728 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
counter500um++;
if (absD > 0.1) {

Check failure on line 2730 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
counter1mm++;
if (absD > 0.2) {

Check failure on line 2732 in PWGDQ/Core/VarManager.h

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
counter2mm++;
if (absD > 0.5) {
counter5mm++;
Expand Down Expand Up @@ -2879,13 +2887,33 @@
values[kQ2X0A2] = ev2.q2x0a();
values[kQ2Y0A1] = ev1.q2y0a();
values[kQ2Y0A2] = ev2.q2y0a();

values[kQ2X0B1] = ev1.q2x0b();
values[kQ2X0B2] = ev2.q2x0b();
values[kQ2Y0B1] = ev1.q2y0b();
values[kQ2Y0B2] = ev2.q2y0b();

values[kQ2X0C1] = ev1.q2x0c();
values[kQ2X0C2] = ev2.q2x0c();
values[kQ2Y0C1] = ev1.q2y0c();
values[kQ2Y0C2] = ev2.q2y0c();
}
if constexpr ((fillMap & CollisionQvect) > 0) {
// Tobe used for the calculation of u1q1 and u2q2
values[kQ2X0A1] = (ev1.qvecBPosRe() * ev1.nTrkBPos() + ev1.qvecBNegRe() * ev1.nTrkBNeg()) / (ev1.nTrkBPos() + ev1.nTrkBNeg());
values[kQ2X0A2] = (ev2.qvecBPosRe() * ev2.nTrkBPos() + ev2.qvecBNegRe() * ev2.nTrkBNeg()) / (ev2.nTrkBPos() + ev2.nTrkBNeg());
values[kQ2Y0A1] = (ev1.qvecBPosIm() * ev1.nTrkBPos() + ev1.qvecBNegIm() * ev1.nTrkBNeg()) / (ev1.nTrkBPos() + ev1.nTrkBNeg());
values[kQ2Y0A2] = (ev2.qvecBPosIm() * ev2.nTrkBPos() + ev2.qvecBNegIm() * ev2.nTrkBNeg()) / (ev2.nTrkBPos() + ev2.nTrkBNeg());

values[kQ2X0B1] = ev1.qvecFT0ARe();
values[kQ2X0B2] = ev2.qvecFT0ARe();
values[kQ2Y0B1] = ev1.qvecFT0AIm();
values[kQ2Y0B2] = ev2.qvecFT0AIm();

values[kQ2X0C1] = ev1.qvecFT0CRe();
values[kQ2X0C2] = ev2.qvecFT0CRe();
values[kQ2Y0C1] = ev1.qvecFT0CIm();
values[kQ2Y0C2] = ev2.qvecFT0CIm();
}
}

Expand Down Expand Up @@ -4287,6 +4315,10 @@
// Compute the scalar product UQ for two muon from different event using Q-vector from A, for second and third harmonic
float Psi2A1 = getEventPlane(2, values[kQ2X0A1], values[kQ2Y0A1]);
float Psi2A2 = getEventPlane(2, values[kQ2X0A2], values[kQ2Y0A2]);
float Psi2B1 = getEventPlane(2, values[kQ2X0B1], values[kQ2Y0B1]);
float Psi2B2 = getEventPlane(2, values[kQ2X0B2], values[kQ2Y0B2]);
float Psi2C1 = getEventPlane(2, values[kQ2X0C1], values[kQ2Y0C1]);
float Psi2C2 = getEventPlane(2, values[kQ2X0C2], values[kQ2Y0C2]);
values[kCos2DeltaPhi] = TMath::Cos(2 * (v12.Phi() - Psi2A1)); // WARNING: using the first event EP
values[kCos2DeltaPhiEv1] = TMath::Cos(2 * (v1.Phi() - Psi2A1));
values[kCos2DeltaPhiEv2] = TMath::Cos(2 * (v2.Phi() - Psi2A2));
Expand Down Expand Up @@ -4317,6 +4349,39 @@
values[kV24ME] = (std::isnan(V22ME) || std::isinf(V22ME) || std::isnan(V24ME) || std::isinf(V24ME)) ? 0. : V24ME;
values[kWV24ME] = (std::isnan(V22ME) || std::isinf(V22ME) || std::isnan(V24ME) || std::isinf(V24ME)) ? 0. : 1.0;

// coherent Jpsi A2
bool useCoherentJpsiA2 = fgUsedVars[kA2EP_RP_TPC] || fgUsedVars[kA2EP_RP_FT0A] || fgUsedVars[kA2EP_RP_FT0C];
if (useCoherentJpsiA2) {
ROOT::Math::Boost boostv12{v12.BoostToCM()};
ROOT::Math::PtEtaPhiMVector v_daughter = boostv12(t1.sign() > 0 ? v1 : v2);
float Psi2A = t1.sign() > 0 ? Psi2A1 : Psi2A2;
float Psi2B = t1.sign() > 0 ? Psi2B1 : Psi2B2;
float Psi2C = t1.sign() > 0 ? Psi2C1 : Psi2C2;

// reaction plane
float phi = RecoDecay::constrainAngle(v_daughter.Phi(), -o2::constants::math::PI);
values[kDeltaPhiRP_TPC] = phi > Psi2A ? phi - Psi2A : Psi2A - phi;
values[kDeltaPhiRP_TPC] = values[kDeltaPhiRP_TPC] > TMath::Pi() ? 2. * TMath::Pi() - values[kDeltaPhiRP_TPC] : values[kDeltaPhiRP_TPC];
values[kDeltaPhiRP_FT0A] = phi > Psi2B ? phi - Psi2B : Psi2B - phi;
values[kDeltaPhiRP_FT0A] = values[kDeltaPhiRP_FT0A] > TMath::Pi() ? 2. * TMath::Pi() - values[kDeltaPhiRP_FT0A] : values[kDeltaPhiRP_FT0A];
values[kDeltaPhiRP_FT0C] = phi > Psi2C ? phi - Psi2C : Psi2C - phi;
values[kDeltaPhiRP_FT0C] = values[kDeltaPhiRP_FT0C] > TMath::Pi() ? 2. * TMath::Pi() - values[kDeltaPhiRP_FT0C] : values[kDeltaPhiRP_FT0C];
// fold delta phi to [0, pi/2]
values[kDeltaPhiRP_TPC] = values[kDeltaPhiRP_TPC] > TMath::Pi() / 2. ? TMath::Pi() - values[kDeltaPhiRP_TPC] : values[kDeltaPhiRP_TPC];
values[kDeltaPhiRP_FT0A] = values[kDeltaPhiRP_FT0A] > TMath::Pi() / 2. ? TMath::Pi() - values[kDeltaPhiRP_FT0A] : values[kDeltaPhiRP_FT0A];
values[kDeltaPhiRP_FT0C] = values[kDeltaPhiRP_FT0C] > TMath::Pi() / 2. ? TMath::Pi() - values[kDeltaPhiRP_FT0C] : values[kDeltaPhiRP_FT0C];
values[kCos2DeltaPhiRP_TPC] = TMath::Cos(2. * (phi - Psi2A));
values[kCos2DeltaPhiRP_FT0A] = TMath::Cos(2. * (phi - Psi2B));
values[kCos2DeltaPhiRP_FT0C] = TMath::Cos(2. * (phi - Psi2C));

float A2RP_TPC = values[kCos2DeltaPhiRP_TPC] / values[kR2EP];
float A2RP_FT0A = values[kCos2DeltaPhiRP_FT0A] / values[kR2EP];
float A2RP_FT0C = values[kCos2DeltaPhiRP_FT0C] / values[kR2EP];
values[kA2EP_RP_TPC] = std::isnan(A2RP_TPC) || std::isinf(A2RP_TPC) ? -999. : A2RP_TPC;
values[kA2EP_RP_FT0A] = std::isnan(A2RP_FT0A) || std::isinf(A2RP_FT0A) ? -999. : A2RP_FT0A;
values[kA2EP_RP_FT0C] = std::isnan(A2RP_FT0C) || std::isinf(A2RP_FT0C) ? -999. : A2RP_FT0C;
}

if constexpr ((fillMap & ReducedEventQvectorExtra) > 0) {
std::complex<double> Q21(values[kQ2X0A] * values[kS11A], values[kQ2Y0A] * values[kS11A]);
std::complex<double> Q42(values[kQ42XA], values[kQ42YA]);
Expand Down Expand Up @@ -4388,6 +4453,38 @@
values[kCos2ThetaStarFT0C] = values[kCosThetaStarFT0C] * values[kCosThetaStarFT0C];
}
}

bool useCoherentJpsiA2 = fgUsedVars[kA2EP_RP_TPC] || fgUsedVars[kA2EP_RP_FT0A] || fgUsedVars[kA2EP_RP_FT0C];
if (useCoherentJpsiA2) {
ROOT::Math::Boost boostv12{v12.BoostToCM()};
ROOT::Math::PtEtaPhiMVector v_daughter = boostv12(v1);
float Psi2A = values[kPsi2A];
float Psi2B = values[kPsi2B];
float Psi2C = values[kPsi2C];

// reaction plane
float phi = RecoDecay::constrainAngle(v_daughter.Phi(), -o2::constants::math::PI);
values[kDeltaPhiRP_TPC] = phi > Psi2A ? phi - Psi2A : Psi2A - phi;
values[kDeltaPhiRP_TPC] = values[kDeltaPhiRP_TPC] > TMath::Pi() ? 2. * TMath::Pi() - values[kDeltaPhiRP_TPC] : values[kDeltaPhiRP_TPC];
values[kDeltaPhiRP_FT0A] = phi > Psi2B ? phi - Psi2B : Psi2B - phi;
values[kDeltaPhiRP_FT0A] = values[kDeltaPhiRP_FT0A] > TMath::Pi() ? 2. * TMath::Pi() - values[kDeltaPhiRP_FT0A] : values[kDeltaPhiRP_FT0A];
values[kDeltaPhiRP_FT0C] = phi > Psi2C ? phi - Psi2C : Psi2C - phi;
values[kDeltaPhiRP_FT0C] = values[kDeltaPhiRP_FT0C] > TMath::Pi() ? 2. * TMath::Pi() - values[kDeltaPhiRP_FT0C] : values[kDeltaPhiRP_FT0C];
// fold delta phi to [0, pi/2]
values[kDeltaPhiRP_TPC] = values[kDeltaPhiRP_TPC] > TMath::Pi() / 2. ? TMath::Pi() - values[kDeltaPhiRP_TPC] : values[kDeltaPhiRP_TPC];
values[kDeltaPhiRP_FT0A] = values[kDeltaPhiRP_FT0A] > TMath::Pi() / 2. ? TMath::Pi() - values[kDeltaPhiRP_FT0A] : values[kDeltaPhiRP_FT0A];
values[kDeltaPhiRP_FT0C] = values[kDeltaPhiRP_FT0C] > TMath::Pi() / 2. ? TMath::Pi() - values[kDeltaPhiRP_FT0C] : values[kDeltaPhiRP_FT0C];
values[kCos2DeltaPhiRP_TPC] = TMath::Cos(2. * (phi - Psi2A));
values[kCos2DeltaPhiRP_FT0A] = TMath::Cos(2. * (phi - Psi2B));
values[kCos2DeltaPhiRP_FT0C] = TMath::Cos(2. * (phi - Psi2C));

float A2RP_TPC = values[kCos2DeltaPhiRP_TPC] / values[kR2EP];
float A2RP_FT0A = values[kCos2DeltaPhiRP_FT0A] / values[kR2EP];
float A2RP_FT0C = values[kCos2DeltaPhiRP_FT0C] / values[kR2EP];
values[kA2EP_RP_TPC] = std::isnan(A2RP_TPC) || std::isinf(A2RP_TPC) ? -999. : A2RP_TPC;
values[kA2EP_RP_FT0A] = std::isnan(A2RP_FT0A) || std::isinf(A2RP_FT0A) ? -999. : A2RP_FT0A;
values[kA2EP_RP_FT0C] = std::isnan(A2RP_FT0C) || std::isinf(A2RP_FT0C) ? -999. : A2RP_FT0C;
}
}

template <int pairType, typename T1, typename T2>
Expand Down
4 changes: 2 additions & 2 deletions PWGDQ/Tasks/tableReader_withAssoc_direct.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1694,7 +1694,7 @@ struct AnalysisSameEventPairing {

constexpr bool eventHasQvector = ((TEventFillMap & VarManager::ObjTypes::CollisionQvectCentr) > 0);
constexpr bool trackHasCov = ((TTrackFillMap & VarManager::ObjTypes::TrackCov) > 0);
constexpr bool fillFlowReso = eventHasQvector;
// constexpr bool fillFlowReso = eventHasQvector;

for (auto& event : events) {
if (!event.isEventSelected_bit(0))
Expand All @@ -1711,7 +1711,7 @@ struct AnalysisSameEventPairing {
if (groupedAssocs.size() == 0)
continue;

if (fillFlowReso) {
if (fConfigOptions.useFlowReso) {
if (ResoFlowSP == nullptr || ResoFlowEP == nullptr) {
LOGF(fatal, "Flow resolution histograms are not available, cannot fill flow variables!");
}
Expand Down
Loading