Skip to content

Commit 5d092e1

Browse files
rolavickalibuild
andauthored
[PWGUD] GlobalMuon matching based on chi2matching (#15938)
Co-authored-by: ALICE Action Bot <alibuild@cern.ch>
1 parent c90afcd commit 5d092e1

1 file changed

Lines changed: 70 additions & 17 deletions

File tree

PWGUD/TableProducer/upcCandProducerGlobalMuon.cxx

Lines changed: 70 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -18,6 +18,8 @@
1818
#include "PWGUD/Core/UPCHelpers.h"
1919
#include "PWGUD/DataModel/UDTables.h"
2020

21+
#include "Common/Core/fwdtrackUtilities.h"
22+
2123
#include <CCDB/BasicCCDBManager.h>
2224
#include <CCDB/CcdbApi.h>
2325
#include <CommonConstants/LHCConstants.h>
@@ -55,6 +57,8 @@
5557
#include <map>
5658
#include <numeric>
5759
#include <string>
60+
#include <unordered_map>
61+
#include <utility>
5862
#include <vector>
5963

6064
using namespace o2::framework;
@@ -64,6 +68,7 @@ struct UpcCandProducerGlobalMuon {
6468
bool fDoMC{false};
6569

6670
std::map<int32_t, int32_t> fNewPartIDs;
71+
std::map<uint32_t, bool> fBestMuonMatch;
6772

6873
Produces<o2::aod::UDMcCollisions> udMCCollisions;
6974
Produces<o2::aod::UDMcParticles> udMCParticles;
@@ -97,6 +102,7 @@ struct UpcCandProducerGlobalMuon {
97102
Configurable<bool> fUseAbsDCA{"fUseAbsDCA", true, "Use absolute DCA in FwdDCAFitter"};
98103
Configurable<float> fMaxChi2MatchMCHMFT{"fMaxChi2MatchMCHMFT", 4.f, "Maximum chi2 for MCH-MFT matching quality filter"};
99104
Configurable<int> fBcWindowMCHMFT{"fBcWindowMCHMFT", 20, "BC window for searching MCH-MFT tracks around MCH-MID-MFT anchors"};
105+
Configurable<bool> fKeepBestMuonMatch{"fKeepBestMuonMatch", true, "Keep only the best MCH-MFT match per MCH track (lowest chi2)"};
100106

101107
using ForwardTracks = o2::soa::Join<o2::aod::FwdTracks, o2::aod::FwdTracksCov>;
102108

@@ -111,9 +117,10 @@ struct UpcCandProducerGlobalMuon {
111117
o2::vertexing::FwdDCAFitterN<2> fFwdFitter;
112118

113119
// Named constants (avoid magic numbers in expressions)
114-
static constexpr double kBcTimeRoundingOffset = 1.; // Offset used when rounding trackTime to BC units
115-
static constexpr uint16_t kMinTracksForPair = 2; // Minimum tracks required to compute a pair invariant mass
116-
static constexpr uint16_t kMinTracksForCandidate = 1; // Minimum contributors required to save a candidate
120+
static constexpr double kBcTimeRoundingOffset = 1.; // Offset used when rounding trackTime to BC units
121+
static constexpr uint16_t kMinTracksForPair = 2; // Minimum tracks required to compute a pair invariant mass
122+
static constexpr uint16_t kMinTracksForCandidate = 1; // Minimum contributors required to save a candidate
123+
static constexpr int kUpperBoundaryToTrackTypeEnum = 2; // Make sure you use MFT tracks
117124

118125
void init(InitContext&)
119126
{
@@ -384,13 +391,15 @@ struct UpcCandProducerGlobalMuon {
384391
return propmuon;
385392
}
386393

387-
bool addToFwdTable(int64_t candId, int64_t trackId, uint64_t gbc, float trackTime, ForwardTracks const& fwdTracks, const o2::aod::McFwdTrackLabels* mcFwdTrackLabels)
394+
bool addToFwdTable(int64_t candId, int64_t trackId, uint64_t gbc, float trackTime,
395+
ForwardTracks const& fwdTracks, o2::aod::MFTTracks const& /*mftTracks*/,
396+
const o2::aod::McFwdTrackLabels* mcFwdTrackLabels)
388397
{
389398
const auto& track = fwdTracks.iteratorAt(trackId);
390399
float px, py, pz;
391400
int sign;
392401

393-
// NEW: Fill track type histogram
402+
// Fill track type histogram
394403
histRegistry.fill(HIST("hTrackTypes"), track.trackType());
395404
if (track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack ||
396405
track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalForwardTrack) {
@@ -402,12 +411,16 @@ struct UpcCandProducerGlobalMuon {
402411
track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalForwardTrack;
403412
o2::dataformats::GlobalFwdTrack pft;
404413
if (isGlobal && fEnableMFT) {
405-
// For global tracks, use raw kinematics for cuts;
406-
// the FwdDCAFitter handles proper vertex propagation
407-
auto trackPar = fwdToTrackPar(track);
408-
pft.setParameters(trackPar.getParameters());
409-
pft.setZ(trackPar.getZ());
410-
pft.setCovariances(trackPar.getCovariances());
414+
// Refit global muon: propagate MCH component to vertex, combine with MFT spatial info
415+
auto mchTrack = track.matchMCHTrack_as<ForwardTracks>();
416+
auto propMuon = propagateToZero(mchTrack);
417+
auto mfttrack = track.matchMFTTrack_as<o2::aod::MFTTracks>();
418+
using SMatrix5 = ROOT::Math::SVector<double, 5>;
419+
using SMatrix55 = ROOT::Math::SMatrix<double, 5, 5, ROOT::Math::MatRepSym<double, 5>>;
420+
SMatrix5 tpars(mfttrack.x(), mfttrack.y(), mfttrack.phi(), mfttrack.tgl(), mfttrack.signed1Pt());
421+
SMatrix55 tcovs{};
422+
o2::track::TrackParCovFwd mft{mfttrack.z(), tpars, tcovs, mfttrack.chi2()};
423+
pft = o2::aod::fwdtrackutils::refitGlobalMuonCov(propMuon, mft);
411424
} else {
412425
pft = propagateToZero(track);
413426
}
@@ -495,7 +508,33 @@ struct UpcCandProducerGlobalMuon {
495508
return o2::track::TrackParCovFwd{track.z(), tpars, tcovs, track.chi2()};
496509
}
497510

511+
// Select the best MCH-MFT match per MCH track based on lowest chi2MatchMCHMFT.
512+
// Multiple global tracks can share the same MCH track with different MFT matches;
513+
// this function keeps only the best one to reduce combinatorial background.
514+
void selectBestMuonMatches(ForwardTracks const& fwdTracks)
515+
{
516+
fBestMuonMatch.clear();
517+
std::unordered_map<int, std::pair<float, int>> mCandidates;
518+
for (const auto& muon : fwdTracks) {
519+
if (static_cast<int>(muon.trackType()) < kUpperBoundaryToTrackTypeEnum) {
520+
auto muonID = muon.matchMCHTrackId();
521+
auto chi2 = muon.chi2MatchMCHMFT();
522+
if (mCandidates.find(muonID) == mCandidates.end()) {
523+
mCandidates[muonID] = {chi2, muon.globalIndex()};
524+
} else {
525+
if (chi2 < mCandidates[muonID].first) {
526+
mCandidates[muonID] = {chi2, muon.globalIndex()};
527+
}
528+
}
529+
}
530+
}
531+
for (const auto& pairCand : mCandidates) {
532+
fBestMuonMatch[pairCand.second.second] = true;
533+
}
534+
}
535+
498536
void createCandidates(ForwardTracks const& fwdTracks,
537+
o2::aod::MFTTracks const& mftTracks,
499538
o2::aod::FwdTrkCls const& fwdTrkCls,
500539
o2::aod::AmbiguousFwdTracks const& ambFwdTracks,
501540
o2::aod::BCs const& bcs,
@@ -580,6 +619,11 @@ struct UpcCandProducerGlobalMuon {
580619
vAmbFwdTrackIndexBCs[ambTr.globalIndex()] = ambTr.bcIds()[0];
581620
}
582621

622+
// Select best MCH-MFT match per MCH track before sorting into BC maps
623+
if (fKeepBestMuonMatch) {
624+
selectBestMuonMatches(fwdTracks);
625+
}
626+
583627
std::map<uint64_t, std::vector<int64_t>> mapGlobalBcsWithMCHMIDTrackIds;
584628
std::map<uint64_t, std::vector<int64_t>> mapGlobalBcsWithMCHTrackIds;
585629
std::map<uint64_t, std::vector<int64_t>> mapGlobalBcsWithGlobalMuonTrackIds; // MCH-MID-MFT (good timing from MID)
@@ -595,6 +639,13 @@ struct UpcCandProducerGlobalMuon {
595639
trackType != GlobalForwardTrack)
596640
continue;
597641

642+
// For global tracks, skip if not the best match for this MCH track
643+
if (fKeepBestMuonMatch && static_cast<int>(trackType) < kUpperBoundaryToTrackTypeEnum) {
644+
if (fBestMuonMatch.find(fwdTrack.globalIndex()) == fBestMuonMatch.end()) {
645+
continue;
646+
}
647+
}
648+
598649
auto trackId = fwdTrack.globalIndex();
599650
int64_t indexBC = vAmbFwdTrackIndex[trackId] < 0 ? vColIndexBCs[fwdTrack.collisionId()] : vAmbFwdTrackIndexBCs[vAmbFwdTrackIndex[trackId]];
600651
auto globalBC = vGlobalBCs[indexBC] + TMath::FloorNint(fwdTrack.trackTime() / o2::constants::lhc::LHCBunchSpacingNS + kBcTimeRoundingOffset);
@@ -675,7 +726,7 @@ struct UpcCandProducerGlobalMuon {
675726
uint16_t numContrib = 0;
676727
double sumPx = 0., sumPy = 0., sumPz = 0., sumE = 0.;
677728
for (const auto& ianchor : vAnchorIds) {
678-
if (!addToFwdTable(candId, ianchor, globalBcAnchor, 0., fwdTracks, mcFwdTrackLabels))
729+
if (!addToFwdTable(candId, ianchor, globalBcAnchor, 0., fwdTracks, mftTracks, mcFwdTrackLabels))
679730
continue;
680731
const auto& trk = fwdTracks.iteratorAt(ianchor);
681732
double p2 = trk.px() * trk.px() + trk.py() * trk.py() + trk.pz() * trk.pz();
@@ -696,7 +747,7 @@ struct UpcCandProducerGlobalMuon {
696747

697748
// Step 3: Write matched MCH-MFT tracks with adjusted track time
698749
for (const auto& [imchMft, gbc] : mapMchMftIdBc) {
699-
if (!addToFwdTable(candId, imchMft, gbc, (gbc - globalBcAnchor) * o2::constants::lhc::LHCBunchSpacingNS, fwdTracks, mcFwdTrackLabels))
750+
if (!addToFwdTable(candId, imchMft, gbc, (gbc - globalBcAnchor) * o2::constants::lhc::LHCBunchSpacingNS, fwdTracks, mftTracks, mcFwdTrackLabels))
700751
continue;
701752
const auto& trk = fwdTracks.iteratorAt(imchMft);
702753
double p2 = trk.px() * trk.px() + trk.py() * trk.py() + trk.pz() * trk.pz();
@@ -758,7 +809,7 @@ struct UpcCandProducerGlobalMuon {
758809
auto& vMuonIds = gbc_muids.second;
759810
// writing MCH-MID tracks
760811
for (const auto& imuon : vMuonIds) {
761-
if (!addToFwdTable(candId, imuon, globalBcMid, 0., fwdTracks, mcFwdTrackLabels))
812+
if (!addToFwdTable(candId, imuon, globalBcMid, 0., fwdTracks, mftTracks, mcFwdTrackLabels))
762813
continue;
763814
numContrib++;
764815
selTrackIds.push_back(imuon);
@@ -769,7 +820,7 @@ struct UpcCandProducerGlobalMuon {
769820
getMchTrackIds(globalBcMid, mapGlobalBcsWithMCHTrackIds, fBcWindowMCH, mapMchIdBc);
770821
// writing MCH-only tracks
771822
for (const auto& [imch, gbc] : mapMchIdBc) {
772-
if (!addToFwdTable(candId, imch, gbc, (gbc - globalBcMid) * o2::constants::lhc::LHCBunchSpacingNS, fwdTracks, mcFwdTrackLabels))
823+
if (!addToFwdTable(candId, imch, gbc, (gbc - globalBcMid) * o2::constants::lhc::LHCBunchSpacingNS, fwdTracks, mftTracks, mcFwdTrackLabels))
773824
continue;
774825
numContrib++;
775826
selTrackIds.push_back(imch);
@@ -816,6 +867,7 @@ struct UpcCandProducerGlobalMuon {
816867
}
817868

818869
void processFwd(ForwardTracks const& fwdTracks,
870+
o2::aod::MFTTracks const& mftTracks,
819871
o2::aod::FwdTrkCls const& fwdTrkCls,
820872
o2::aod::AmbiguousFwdTracks const& ambFwdTracks,
821873
o2::aod::BCs const& bcs,
@@ -824,10 +876,11 @@ struct UpcCandProducerGlobalMuon {
824876
o2::aod::Zdcs const& zdcs)
825877
{
826878
fDoMC = false;
827-
createCandidates(fwdTracks, fwdTrkCls, ambFwdTracks, bcs, collisions, fv0s, zdcs, (o2::aod::McFwdTrackLabels*)nullptr);
879+
createCandidates(fwdTracks, mftTracks, fwdTrkCls, ambFwdTracks, bcs, collisions, fv0s, zdcs, (o2::aod::McFwdTrackLabels*)nullptr);
828880
}
829881

830882
void processFwdMC(ForwardTracks const& fwdTracks,
883+
o2::aod::MFTTracks const& mftTracks,
831884
o2::aod::FwdTrkCls const& fwdTrkCls,
832885
o2::aod::AmbiguousFwdTracks const& ambFwdTracks,
833886
o2::aod::BCs const& bcs,
@@ -840,7 +893,7 @@ struct UpcCandProducerGlobalMuon {
840893
{
841894
fDoMC = true;
842895
skimMCInfo(mcCollisions, mcParticles);
843-
createCandidates(fwdTracks, fwdTrkCls, ambFwdTracks, bcs, collisions, fv0s, zdcs, &mcFwdTrackLabels);
896+
createCandidates(fwdTracks, mftTracks, fwdTrkCls, ambFwdTracks, bcs, collisions, fv0s, zdcs, &mcFwdTrackLabels);
844897
fNewPartIDs.clear();
845898
}
846899

0 commit comments

Comments
 (0)