Skip to content

[PWGDQ] Code cleanup and logic improvements for muon processing in tableMaker(MC)_withAssoc #17893

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

Merged
mcoquet642 merged 9 commits into from
Sep 13, 2026
+119 −121
Merged

[PWGDQ] Code cleanup and logic improvements for muon processing in tableMaker(MC)_withAssoc #17893

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
110f7bc
processPPMuonRefit in tableMaker_withAssoc
mcoquet642 Sep 12, 2026
d0091ae
Moving computation of muon dca to FillTrackCollision
mcoquet642 Sep 12, 2026
bd1e8e9
Removing obsolete FillMuonPDca
mcoquet642 Sep 12, 2026
67773a6
Reorganizing computation of dca and pdca
mcoquet642 Sep 12, 2026
5dcad1e
logic fix
mcoquet642 Sep 12, 2026
3841990
removing old tests
mcoquet642 Sep 13, 2026
0ec5d32
removing deprecated PV shift options
mcoquet642 Sep 13, 2026
fbd1553
fixing few errors
mcoquet642 Sep 13, 2026
3a3348d
dummy commit
mcoquet642 Sep 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
92 changes: 39 additions & 53 deletions PWGDQ/Core/VarManager.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1449,8 +1449,6 @@
template <typename T, typename C>
static o2::track::TrackParCovFwd PropagateFwd(const T& track, const C& cov, float z);
template <uint32_t fillMap, typename T, typename C>
static void FillMuonPDca(const T& muon, const C& collision, float* values = nullptr);
template <uint32_t fillMap, typename T, typename C>
static void FillPropagateMuon(const T& muon, const C& collision, float* values = nullptr);
template <typename T>
static void FillBC(T const& bc, float* values = nullptr);
Expand Down Expand Up @@ -1593,7 +1591,7 @@

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

Check failure on line 1594 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 @@ -1836,12 +1834,8 @@
GetFwdShiftForY(muon.y(), xShift, yShift, zShift);
o2::track::TrackParCovFwd fwdtrack = o2::aod::fwdtrackutils::getTrackParCovFwd3DShift(muon, xShift, yShift, zShift, muon);
o2::dataformats::GlobalFwdTrack propmuon;
if (static_cast<int>(muon.trackType()) > 2) {

Check failure on line 1837 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());
track.setCovariances(fwdtrack.getCovariances());
auto mchTrack = mMatching.FwdtoMCH(track);
auto mchTrack = mMatching.FwdtoMCH(fwdtrack);

if (endPoint == kToVertex) {
o2::mch::TrackExtrap::extrapToVertex(mchTrack, collision.posX(), collision.posY(), collision.posZ(), collision.covXX(), collision.covYY());
Expand All @@ -1856,17 +1850,12 @@
o2::mch::TrackExtrap::extrapToVertexWithoutBranson(mchTrack, fgzMatching);
}

auto proptrack = mMatching.MCHtoFwd(mchTrack);
propmuon.setParameters(proptrack.getParameters());
propmuon.setZ(proptrack.getZ());
propmuon.setCovariances(proptrack.getCovariances());
propmuon = mMatching.MCHtoFwd(mchTrack);

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

Check failure on line 1855 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());
propmuon.setZ(fwdtrack.getZ());
propmuon.setCovariances(fwdtrack.getCovariances());
propmuon = fwdtrack;
}
return propmuon;
}
Expand All @@ -1879,25 +1868,6 @@
return fwdtrack;
}

template <uint32_t fillMap, typename T, typename C>
void VarManager::FillMuonPDca(const T& muon, const C& collision, float* values)
{
if (!values) {
values = fgValues;
}

if constexpr ((fillMap & MuonCov) > 0 || (fillMap & ReducedMuonCov) > 0) {

o2::dataformats::GlobalFwdTrack propmuon = PropagateMuon(muon, collision);
o2::dataformats::GlobalFwdTrack propmuonAtDCA = PropagateMuon(muon, collision, kToDCA);

float dcaX = (propmuonAtDCA.getX() - collision.posX());
float dcaY = (propmuonAtDCA.getY() - collision.posY());
float dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
values[kMuonPDca] = muon.p() * dcaXY;
}
}

template <uint32_t fillMap, typename T, typename C>
void VarManager::FillPropagateMuon(const T& muon, const C& collision, float* values)
{
Expand All @@ -1921,20 +1891,6 @@
values[kTgl] = propmuon.getTgl();
values[kPhi] = propmuon.getPhi();

// Redo propagation only for muon tracks
// propagation of MFT tracks alredy done in fwdtrack-extention task
if (static_cast<int>(muon.trackType()) > 2) {
o2::dataformats::GlobalFwdTrack propmuonAtDCA = PropagateMuon(muon, collision, kToDCA);
o2::dataformats::GlobalFwdTrack propmuonAtRabs = PropagateMuon(muon, collision, kToRabs);
float dcaX = (propmuonAtDCA.getX() - collision.posX());
float dcaY = (propmuonAtDCA.getY() - collision.posY());
values[kMuonDCAx] = dcaX;
values[kMuonDCAy] = dcaY;
double xAbs = propmuonAtRabs.getX();
double yAbs = propmuonAtRabs.getY();
values[kMuonRAtAbsorberEnd] = std::sqrt(xAbs * xAbs + yAbs * yAbs);
}

const SMatrix55& cov = propmuon.getCovariances();
values[kMuonCXX] = cov(0, 0);
values[kMuonCXY] = cov(1, 0);
Expand Down Expand Up @@ -1971,6 +1927,7 @@
double pz = propmuon.getP() * std::cos(o2::constants::math::PIHalf - std::atan(mfttrack.tgl()));
double pt = std::sqrt(std::pow(px, 2) + std::pow(py, 2));
auto mftprop = o2::aod::fwdtrackutils::getTrackParCovFwd3DShift(mfttrack, xShift, yShift, zShift);
mftprop.setInvQPt(static_cast<double>(muontrack.sign()) / pt);
values[kX] = mftprop.getX();
values[kY] = mftprop.getY();
values[kZ] = mftprop.getZ();
Expand All @@ -1979,6 +1936,12 @@
values[kPz] = pz;
values[kEta] = mftprop.getEta();
values[kPhi] = mftprop.getPhi();

// Helix DCA of the refitted global track w.r.t. the associated collision
std::array<double, 3> dca{999., 999., 999.};
mftprop.propagateToDCAhelix(fgMagField, {collision.posX(), collision.posY(), collision.posZ()}, dca);
values[kMuonDCAx] = static_cast<float>(dca[0]);
values[kMuonDCAy] = static_cast<float>(dca[1]);
}
}

Expand Down Expand Up @@ -2006,6 +1969,12 @@
values[kPz] = globalRefit.getPz();
values[kEta] = globalRefit.getEta();
values[kPhi] = globalRefit.getPhi();

// Helix DCA of the covariance-refitted global track w.r.t. the associated collision
std::array<double, 3> dca{999., 999., 999.};
globalRefit.propagateToDCAhelix(fgMagField, {collision.posX(), collision.posY(), collision.posZ()}, dca);
values[kMuonDCAx] = static_cast<float>(dca[0]);
values[kMuonDCAy] = static_cast<float>(dca[1]);
}
}
}
Expand Down Expand Up @@ -2732,7 +2701,7 @@
if (!track.hasTPC()) {
continue; // skip tracks without TPC information
}
if (track.dcaZ() > 998) {

Check failure on line 2704 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 @@ -2820,17 +2789,17 @@
int counter10mm = 0;
for (auto const& d : dcazValues) {
double absD = std::abs(d);
if (absD > 0.01) {

Check failure on line 2792 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 2794 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 2796 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 2798 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 2800 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) {

Check failure on line 2802 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.
counter5mm++;
if (absD > 1.0) {
counter10mm++;
Expand Down Expand Up @@ -3477,13 +3446,30 @@
}
}
if constexpr ((fillMap & MuonCov) > 0 || (fillMap & MuonCovRealign) > 0 || (fillMap & ReducedMuonCov) > 0) {
float dcaX = 999.f;
float dcaY = 999.f;
if (static_cast<int>(track.trackType()) <= 2) {
// Global / MCH-MID: helix DCA only (for globals, kMuonPDca is filled from the matched MCH in skimMuons)
float xShift = 0.f;
float yShift = 0.f;
float zShift = 0.f;
GetFwdShiftForY(track.y(), xShift, yShift, zShift);
o2::track::TrackParCovFwd fwdtrack = o2::aod::fwdtrackutils::getTrackParCovFwd3DShift(track, xShift, yShift, zShift, track);
std::array<double, 3> dca{999., 999., 999.};
fwdtrack.propagateToDCAhelix(fgMagField, {collision.posX(), collision.posY(), collision.posZ()}, dca);
dcaX = static_cast<float>(dca[0]);
dcaY = static_cast<float>(dca[1]);
} else {
// MCH standalone: DCA and pDCA from MCH extrapolation
o2::dataformats::GlobalFwdTrack propmuonAtDCA = PropagateMuon(track, collision, kToDCA);
dcaX = propmuonAtDCA.getX() - collision.posX();
dcaY = propmuonAtDCA.getY() - collision.posY();
float dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
values[kMuonPDca] = track.p() * dcaXY;
}

o2::dataformats::GlobalFwdTrack propmuonAtDCA = PropagateMuon(track, collision, kToDCA);

float dcaX = (propmuonAtDCA.getX() - collision.posX());
float dcaY = (propmuonAtDCA.getY() - collision.posY());
float dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
values[kMuonPDca] = track.p() * dcaXY;
values[kMuonDCAx] = dcaX;
values[kMuonDCAy] = dcaY;
}
}

Expand Down
82 changes: 41 additions & 41 deletions PWGDQ/TableProducer/tableMakerMC_withAssoc.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -144,7 +144,7 @@ template <typename TMap>
void PrintBitMap(TMap map, int nbits)
{
for (int i = 0; i < nbits; i++) {
cout << ((map & (TMap(1) << i)) > 0 ? "1" : "0");
LOG(info) << ((map & (TMap(1) << i)) > 0 ? "1" : "0");
}
}
*/
Expand Down Expand Up @@ -234,9 +234,6 @@ struct TableMakerMC {
Configurable<std::string> fConfigCcdbUrl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
Configurable<std::string> fGeoPath{"geoPath", "GLO/Config/GeometryAligned", "Path of the geometry file"};
Configurable<std::string> fGrpMagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"};
Configurable<std::string> fZShiftPath{"zShiftPath", "Users/m/mcoquet/ZShift", "CCDB path for z shift to apply to forward tracks"};
Configurable<bool> fUseRemoteZShift{"cfgUseRemoteZShift", false, "Enable getting Zshift from ccdb"};
Configurable<float> fManualZShift{"cfgManualZShift", 0.f, "Manual value for the Zshift for muons."};
Configurable<std::string> fGrpMagPathRun2{"grpmagPathRun2", "GLO/GRP/GRP", "CCDB path of the GRPObject (Usage for Run 2)"};
Configurable<int64_t> timestampCCDB{"timestampCCDB", -1, "timestamp of the ONNX file for ML model used to query in CCDB"};
} fConfigCCDB;
Expand Down Expand Up @@ -580,17 +577,16 @@ struct TableMakerMC {
/*if ((std::abs(mctrack.pdgCode())>400 && std::abs(mctrack.pdgCode())<599) ||
(std::abs(mctrack.pdgCode())>4000 && std::abs(mctrack.pdgCode())<5999) ||
(mcflags > 0)) {
cout << ">>>>>>>>>>>>>>>>>>>>>>> track idx / pdg / process / status code / HEPMC status / primary : "
<< mctrack.globalIndex() << " / " << mctrack.pdgCode() << " / "
<< mctrack.getProcess() << " / " << mctrack.getGenStatusCode() << " / " << mctrack.getHepMCStatusCode() << " / " << mctrack.isPhysicalPrimary() << endl;
cout << ">>>>>>>>>>>>>>>>>>>>>>> track bitmap: ";
LOG(info) << ">>>>>>>>>>>>>>>>>>>>>> track idx / pdg / process / status code / HEPMC status / primary : "
<< mctrack.globalIndex() << " / " << mctrack.pdgCode() << " / "
<< mctrack.getProcess() << " / " << mctrack.getGenStatusCode() << " / " << mctrack.getHepMCStatusCode() << " / " << mctrack.isPhysicalPrimary();
LOG(info) << ">>>>>>>>>>>>>>>>>>>>>> track bitmap: ";
PrintBitMap(mcflags, 16);
cout << endl;
if (mctrack.has_mothers()) {
for (const auto& m : mctrack.mothersIds()) {
if (m < mcTracks.size()) { // protect against bad mother indices
auto aMother = mcTracks.rawIteratorAt(m);
cout << "<<<<<< mother idx / pdg: " << m << " / " << aMother.pdgCode() << endl;
LOG(info) << "<<<<<< mother idx / pdg: " << m << " / " << aMother.pdgCode();
}
}
}
Expand All @@ -600,7 +596,7 @@ struct TableMakerMC {

if (d < mcTracks.size()) { // protect against bad daughter indices
auto aDaughter = mcTracks.rawIteratorAt(d);
cout << "<<<<<< daughter idx / pdg: " << d << " / " << aDaughter.pdgCode() << endl;
LOG(info) << "<<<<<< daughter idx / pdg: " << d << " / " << aDaughter.pdgCode();
}
}
}
Expand Down Expand Up @@ -1133,22 +1129,29 @@ struct TableMakerMC {
VarManager::FillTrack<TMuonFillMap>(muon);
// NOTE: If a muon is associated to multiple collisions, depending on the selections,
// it may be accepted for some associations and rejected for other
if (fConfigVariousOptions.fPropMuon) {
if (static_cast<int>(muon.trackType()) > 2 && fConfigVariousOptions.fPropMuon) {
VarManager::FillPropagateMuon<TMuonFillMap>(muon, collision);
}
// recalculte pDca and global muon kinematics
if (static_cast<int>(muon.trackType()) < 2 && fConfigVariousOptions.fRefitGlobalMuon) {
// recalculate pDca / DCA and global muon kinematics
// kMuonPDca is always taken from MCH (standalone or the MCH matched to a global)
if (static_cast<int>(muon.trackType()) <= 2) {
auto muontrack = muon.template matchMCHTrack_as<TMuons>();
if (muontrack.eta() < fConfigVariousOptions.fMuonMatchEtaMin || muontrack.eta() > fConfigVariousOptions.fMuonMatchEtaMax) {
continue;
}
auto mfttrack = muon.template matchMFTTrack_as<TMFTTracks>();
VarManager::FillTrackCollision<TMuonFillMap>(muontrack, collision);
if constexpr (static_cast<bool>(TMFTFillMap & VarManager::ObjTypes::MFTCov)) {
auto const& mfttrackcov = mfCovs.rawIteratorAt(map_mfttrackcovs[mfttrack.globalIndex()]);
VarManager::FillGlobalMuonRefitCov<TMuonFillMap, TMFTFillMap>(muontrack, mfttrack, collision, mfttrackcov);
if (fConfigVariousOptions.fRefitGlobalMuon) {
if (muontrack.eta() < fConfigVariousOptions.fMuonMatchEtaMin || muontrack.eta() > fConfigVariousOptions.fMuonMatchEtaMax) {
continue;
}
auto mfttrack = muon.template matchMFTTrack_as<TMFTTracks>();
// Helix DCA (kMuonDCAx/y) is filled from the refitted parameters inside FillGlobalMuonRefit(Cov)
if constexpr (static_cast<bool>(TMFTFillMap & VarManager::ObjTypes::MFTCov)) {
auto const& mfttrackcov = mfCovs.rawIteratorAt(map_mfttrackcovs[mfttrack.globalIndex()]);
VarManager::FillGlobalMuonRefitCov<TMuonFillMap, TMFTFillMap>(muontrack, mfttrack, collision, mfttrackcov);
} else {
VarManager::FillGlobalMuonRefit<TMuonFillMap>(muontrack, mfttrack, collision);
}
} else {
VarManager::FillGlobalMuonRefit<TMuonFillMap>(muontrack, mfttrack, collision);
// Helix DCA of the global track; leaves kMuonPDca from the matched MCH above
VarManager::FillTrackCollision<TMuonFillMap>(muon, collision);
}
} else {
VarManager::FillTrackCollision<TMuonFillMap>(muon, collision);
Expand Down Expand Up @@ -1262,21 +1265,28 @@ struct TableMakerMC {
}

VarManager::FillTrack<TMuonFillMap>(muon);
if (fConfigVariousOptions.fPropMuon) {
if (static_cast<int>(muon.trackType()) > 2 && fConfigVariousOptions.fPropMuon) {
VarManager::FillPropagateMuon<TMuonFillMap>(muon, collision);
}
// recalculte pDca and global muon kinematics
// recalculate pDca / DCA and global muon kinematics
// kMuonPDca is always taken from MCH (standalone or the MCH matched to a global)
int globalClusters = muon.nClusters();
if (static_cast<int>(muon.trackType()) < 2 && fConfigVariousOptions.fRefitGlobalMuon) {
if (static_cast<int>(muon.trackType()) <= 2) {
auto muontrack = muon.template matchMCHTrack_as<TMuons>();
auto mfttrack = muon.template matchMFTTrack_as<TMFTTracks>();
globalClusters += mfttrack.nClusters();
VarManager::FillTrackCollision<TMuonFillMap>(muontrack, collision);
if constexpr (static_cast<bool>(TMFTFillMap & VarManager::ObjTypes::MFTCov)) {
auto const& mfttrackcov = mfCovs.rawIteratorAt(map_mfttrackcovs[mfttrack.globalIndex()]);
VarManager::FillGlobalMuonRefitCov<TMuonFillMap, TMFTFillMap>(muontrack, mfttrack, collision, mfttrackcov);
if (fConfigVariousOptions.fRefitGlobalMuon) {
auto mfttrack = muon.template matchMFTTrack_as<TMFTTracks>();
globalClusters += mfttrack.nClusters();
// Helix DCA (kMuonDCAx/y) is filled from the refitted parameters inside FillGlobalMuonRefit(Cov)
if constexpr (static_cast<bool>(TMFTFillMap & VarManager::ObjTypes::MFTCov)) {
auto const& mfttrackcov = mfCovs.rawIteratorAt(map_mfttrackcovs[mfttrack.globalIndex()]);
VarManager::FillGlobalMuonRefitCov<TMuonFillMap, TMFTFillMap>(muontrack, mfttrack, collision, mfttrackcov);
} else {
VarManager::FillGlobalMuonRefit<TMuonFillMap>(muontrack, mfttrack, collision);
}
} else {
VarManager::FillGlobalMuonRefit<TMuonFillMap>(muontrack, mfttrack, collision);
// Helix DCA of the global track; leaves kMuonPDca from the matched MCH above
VarManager::FillTrackCollision<TMuonFillMap>(muon, collision);
}
} else {
VarManager::FillTrackCollision<TMuonFillMap>(muon, collision);
Expand Down Expand Up @@ -1323,16 +1333,6 @@ struct TableMakerMC {
o2::base::Propagator::initFieldFromGRP(fGrpMag);
VarManager::SetMagneticField(fGrpMag->getNominalL3Field());
}
if (fConfigCCDB.fUseRemoteZShift) {
auto* fZShift = fCCDB->getForTimeStamp<std::vector<float>>(fConfigCCDB.fZShiftPath, bcs.begin().timestamp());
if (fZShift != nullptr && !fZShift->empty()) {
VarManager::SetZShift((*fZShift)[0]);
} else {
LOG(fatal) << "Could not retrieve Z-shift value from CCDB";
}
} else {
VarManager::SetZShift(fConfigCCDB.fManualZShift.value);
}
if (fConfigVariousOptions.fPropMuon) {
VarManager::SetupMuonMagField();
}
Expand Down
Loading
Loading