Skip to content

Commit a727717

Browse files
rolavickalibuild
rolavick
and
alibuild
authored
[PWGUD] Modification to global muon producer (#15697)
Co-authored-by: ALICE Action Bot <alibuild@cern.ch>
1 parent e908014 commit a727717

File tree

1 file changed

+69
-52
lines changed

1 file changed

+69
-52
lines changed

‎PWGUD/TableProducer/upcCandProducerGlobalMuon.cxx‎

Lines changed: 69 additions & 52 deletions
Original file line numberDiff line numberDiff line change
@@ -113,6 +113,12 @@ struct UpcCandProducerGlobalMuon {
113113
static constexpr double fCcenterMFT[3] = {0, 0, -61.4}; // Field evaluation point at center of MFT
114114
float fZShift{0}; // z-vertex shift for forward track propagation
115115

116+
// Named constants (avoid magic numbers in expressions)
117+
static constexpr double kInvalidDCA = 999.; // Sentinel for "no valid DCA computed yet"
118+
static constexpr double kBcTimeRoundingOffset = 1.; // Offset used when rounding trackTime to BC units
119+
static constexpr uint16_t kMinTracksForPair = 2; // Minimum tracks required to compute a pair invariant mass
120+
static constexpr uint16_t kMinTracksForCandidate = 1; // Minimum contributors required to save a candidate
121+
116122
void init(InitContext&)
117123
{
118124
fUpcCuts = (UPCCutparHolder)fUpcCutsConf;
@@ -132,30 +138,28 @@ struct UpcCandProducerGlobalMuon {
132138
histRegistry.get<TH1>(HIST("MuonsSelCounter"))->GetXaxis()->SetBinLabel(upchelpers::kFwdSelChi2 + 1, "Chi2");
133139

134140
// NEW: Add histograms for global track monitoring
135-
if (fEnableMFT) {
136-
const AxisSpec axisTrackType{5, -0.5, 4.5, "Track Type"};
137-
histRegistry.add("hTrackTypes", "Track type distribution", kTH1F, {axisTrackType});
138-
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(1, "MuonStandalone");
139-
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(2, "MCHStandalone");
140-
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(3, "GlobalMuon");
141-
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(4, "GlobalFwd");
142-
143-
const AxisSpec axisEta{100, -4.0, -2.0, "#eta"};
144-
histRegistry.add("hEtaGlobal", "Global track eta", kTH1F, {axisEta});
145-
146-
const AxisSpec axisDCAxy{200, 0., 10., "DCA_{xy} (cm)"};
147-
const AxisSpec axisDCAz{200, -10., 10., "DCA_{z} (cm)"};
148-
histRegistry.add("hDCAxyGlobal", "DCAxy of global tracks to best collision", kTH1F, {axisDCAxy});
149-
histRegistry.add("hDCAzGlobal", "DCAz of global tracks to best collision", kTH1F, {axisDCAz});
150-
histRegistry.add("hNCompatColls", "Number of compatible collisions per global track", kTH1F, {{21, -0.5, 20.5}});
151-
152-
const AxisSpec axisChi2Match{200, 0., 100., "#chi^{2}_{MCH-MFT}"};
153-
histRegistry.add("hChi2MatchMCHMFT", "Chi2 of MCH-MFT matching (before cut)", kTH1F, {axisChi2Match});
154-
155-
const AxisSpec axisMass{500, 0., 10., "m_{inv} (GeV/c^{2})"};
156-
histRegistry.add("hMassGlobalMuon", "Invariant mass from MCH-MID-MFT tracks only", kTH1F, {axisMass});
157-
histRegistry.add("hMassGlobalMuonWithMCHMFT", "Invariant mass from MCH-MID-MFT + MCH-MFT tracks", kTH1F, {axisMass});
158-
}
141+
const AxisSpec axisTrackType{5, -0.5, 4.5, "Track Type"};
142+
histRegistry.add("hTrackTypes", "Track type distribution", kTH1F, {axisTrackType});
143+
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(1, "MuonStandalone");
144+
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(2, "MCHStandalone");
145+
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(3, "GlobalMuon");
146+
histRegistry.get<TH1>(HIST("hTrackTypes"))->GetXaxis()->SetBinLabel(4, "GlobalFwd");
147+
148+
const AxisSpec axisEta{100, -4.0, -2.0, "#eta"};
149+
histRegistry.add("hEtaGlobal", "Global track eta", kTH1F, {axisEta});
150+
151+
const AxisSpec axisDCAxy{200, 0., 10., "DCA_{xy} (cm)"};
152+
const AxisSpec axisDCAz{200, -10., 10., "DCA_{z} (cm)"};
153+
histRegistry.add("hDCAxyGlobal", "DCAxy of global tracks to best collision", kTH1F, {axisDCAxy});
154+
histRegistry.add("hDCAzGlobal", "DCAz of global tracks to best collision", kTH1F, {axisDCAz});
155+
histRegistry.add("hNCompatColls", "Number of compatible collisions per global track", kTH1F, {{21, -0.5, 20.5}});
156+
157+
const AxisSpec axisChi2Match{200, 0., 100., "#chi^{2}_{MCH-MFT}"};
158+
histRegistry.add("hChi2MatchMCHMFT", "Chi2 of MCH-MFT matching (before cut)", kTH1F, {axisChi2Match});
159+
160+
const AxisSpec axisMass{500, 0., 10., "m_{inv} (GeV/c^{2})"};
161+
histRegistry.add("hMassGlobalMuon", "Invariant mass from MCH-MID-MFT tracks only", kTH1F, {axisMass});
162+
histRegistry.add("hMassGlobalMuonWithMCHMFT", "Invariant mass from MCH-MID-MFT + MCH-MFT tracks", kTH1F, {axisMass});
159163
}
160164

161165
bool cut(const o2::dataformats::GlobalFwdTrack& pft, const ForwardTracks::iterator& fwdTrack)
@@ -393,13 +397,11 @@ struct UpcCandProducerGlobalMuon {
393397
float px, py, pz;
394398
int sign;
395399

396-
// NEW: Fill track type histogram if MFT enabled
397-
if (fEnableMFT) {
398-
histRegistry.fill(HIST("hTrackTypes"), track.trackType());
399-
if (track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack ||
400-
track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalForwardTrack) {
401-
histRegistry.fill(HIST("hEtaGlobal"), track.eta());
402-
}
400+
// NEW: Fill track type histogram
401+
histRegistry.fill(HIST("hTrackTypes"), track.trackType());
402+
if (track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack ||
403+
track.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalForwardTrack) {
404+
histRegistry.fill(HIST("hEtaGlobal"), track.eta());
403405
}
404406

405407
if (fUpcCuts.getUseFwdCuts()) {
@@ -476,11 +478,11 @@ struct UpcCandProducerGlobalMuon {
476478
}
477479

478480
int newId = 0;
479-
for (auto trackId : trackIds) {
481+
for (const auto& trackId : trackIds) {
480482
auto it = clustersPerTrack.find(trackId);
481483
if (it != clustersPerTrack.end()) {
482484
const auto& clusters = it->second;
483-
for (auto clsId : clusters) {
485+
for (const auto& clsId : clusters) {
484486
const auto& clsInfo = fwdTrkCls.iteratorAt(clsId);
485487
udFwdTrkClusters(newId, clsInfo.x(), clsInfo.y(), clsInfo.z(), clsInfo.clInfo());
486488
}
@@ -502,6 +504,23 @@ struct UpcCandProducerGlobalMuon {
502504
return {dcaXY, dcaOrig[2], dcaOrig[0], dcaOrig[1]};
503505
}
504506

507+
// Dispatch DCA computation: use MFT helix propagation when fEnableMFT is true,
508+
// otherwise fall back to MCH extrapolation to (0,0,0) via propagateToZero.
509+
// Returns {DCAxy, DCAz, DCAx, DCAy}
510+
std::array<double, 4> propagateFwdToDCA(ForwardTracks::iterator const& track,
511+
double colX, double colY, double colZ)
512+
{
513+
if (fEnableMFT) {
514+
return propagateGlobalToDCA(track, colX, colY, colZ);
515+
}
516+
auto pft = propagateToZero(track);
517+
double dcaX = pft.getX() - colX;
518+
double dcaY = pft.getY() - colY;
519+
double dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
520+
double dcaZ = pft.getZ() - colZ;
521+
return {dcaXY, dcaZ, dcaX, dcaY};
522+
}
523+
505524
void createCandidates(ForwardTracks const& fwdTracks,
506525
o2::aod::FwdTrkCls const& fwdTrkCls,
507526
o2::aod::AmbiguousFwdTracks const& ambFwdTracks,
@@ -610,18 +629,18 @@ struct UpcCandProducerGlobalMuon {
610629

611630
auto trackId = fwdTrack.globalIndex();
612631
int64_t indexBC = vAmbFwdTrackIndex[trackId] < 0 ? vColIndexBCs[fwdTrack.collisionId()] : vAmbFwdTrackIndexBCs[vAmbFwdTrackIndex[trackId]];
613-
auto globalBC = vGlobalBCs[indexBC] + TMath::FloorNint(fwdTrack.trackTime() / o2::constants::lhc::LHCBunchSpacingNS + 1.);
632+
auto globalBC = vGlobalBCs[indexBC] + TMath::FloorNint(fwdTrack.trackTime() / o2::constants::lhc::LHCBunchSpacingNS + kBcTimeRoundingOffset);
614633

615634
if (trackType == MuonStandaloneTrack) { // MCH-MID
616635
mapGlobalBcsWithMCHMIDTrackIds[globalBC].push_back(trackId);
617636
} else if (trackType == MCHStandaloneTrack) { // MCH-only
618637
mapGlobalBcsWithMCHTrackIds[globalBC].push_back(trackId);
619-
} else if (fEnableMFT && trackType == GlobalMuonTrack) { // MCH-MID-MFT: good timing, used as anchor
638+
} else if (trackType == GlobalMuonTrack) { // MCH-MID-MFT: good timing, used as anchor
620639
histRegistry.fill(HIST("hChi2MatchMCHMFT"), fwdTrack.chi2MatchMCHMFT());
621640
if (fwdTrack.chi2MatchMCHMFT() > 0 && fwdTrack.chi2MatchMCHMFT() < fMaxChi2MatchMCHMFT) {
622641
mapGlobalBcsWithGlobalMuonTrackIds[globalBC].push_back(trackId);
623642
}
624-
} else if (fEnableMFT && trackType == GlobalForwardTrack) { // MCH-MFT: poor timing, matched to anchors
643+
} else if (trackType == GlobalForwardTrack) { // MCH-MFT: poor timing, matched to anchors
625644
histRegistry.fill(HIST("hChi2MatchMCHMFT"), fwdTrack.chi2MatchMCHMFT());
626645
if (fwdTrack.chi2MatchMCHMFT() > 0 && fwdTrack.chi2MatchMCHMFT() < fMaxChi2MatchMCHMFT) {
627646
mapGlobalBcsWithMCHMFTTrackIds[globalBC].push_back(trackId);
@@ -631,11 +650,9 @@ struct UpcCandProducerGlobalMuon {
631650

632651
// Map global BC to collisions for DCA-based vertex assignment
633652
std::map<uint64_t, std::vector<int64_t>> mapGlobalBCtoCollisions;
634-
if (fEnableMFT) {
635-
for (const auto& col : collisions) {
636-
uint64_t gbc = vGlobalBCs[col.bcId()];
637-
mapGlobalBCtoCollisions[gbc].push_back(col.globalIndex());
638-
}
653+
for (const auto& col : collisions) {
654+
uint64_t gbc = vGlobalBCs[col.bcId()];
655+
mapGlobalBCtoCollisions[gbc].push_back(col.globalIndex());
639656
}
640657

641658
std::vector<int64_t> selTrackIds{}; // NEW: For cluster saving
@@ -645,7 +662,7 @@ struct UpcCandProducerGlobalMuon {
645662
// Process global tracks: MCH-MID-MFT anchors + MCH-MFT in BC window
646663
// MCH-MID-MFT tracks have good timing (from MID) and serve as anchors.
647664
// MCH-MFT tracks have poor timing and are searched in a BC window around anchors.
648-
if (fEnableMFT && !mapGlobalBcsWithGlobalMuonTrackIds.empty()) {
665+
if (!mapGlobalBcsWithGlobalMuonTrackIds.empty()) {
649666
for (const auto& gbc_anchorids : mapGlobalBcsWithGlobalMuonTrackIds) {
650667
uint64_t globalBcAnchor = gbc_anchorids.first;
651668
auto itFv0Id = mapGlobalBcWithV0A.find(globalBcAnchor);
@@ -675,7 +692,7 @@ struct UpcCandProducerGlobalMuon {
675692

676693
// Step 1: Find best collision vertex using DCA-based propagation
677694
float bestVtxX = 0., bestVtxY = 0., bestVtxZ = 0.;
678-
double bestAvgDCA = 999.;
695+
double bestAvgDCA = kInvalidDCA;
679696
bool hasVertex = false;
680697
int nCompatColls = 0;
681698

@@ -684,18 +701,18 @@ struct UpcCandProducerGlobalMuon {
684701
auto itCol = mapGlobalBCtoCollisions.find(searchBC);
685702
if (itCol == mapGlobalBCtoCollisions.end())
686703
continue;
687-
for (auto colIdx : itCol->second) {
704+
for (const auto& colIdx : itCol->second) {
688705
nCompatColls++;
689706
const auto& col = collisions.iteratorAt(colIdx);
690707
double sumDCAxy = 0.;
691708
int nTracks = 0;
692709
for (const auto& iglobal : allTrackIds) {
693710
const auto& trk = fwdTracks.iteratorAt(iglobal);
694-
auto dca = propagateGlobalToDCA(trk, col.posX(), col.posY(), col.posZ());
711+
auto dca = propagateFwdToDCA(trk, col.posX(), col.posY(), col.posZ());
695712
sumDCAxy += dca[0];
696713
nTracks++;
697714
}
698-
double avgDCA = nTracks > 0 ? sumDCAxy / nTracks : 999.;
715+
double avgDCA = nTracks > 0 ? sumDCAxy / nTracks : kInvalidDCA;
699716
if (!hasVertex || avgDCA < bestAvgDCA) {
700717
bestAvgDCA = avgDCA;
701718
bestVtxX = col.posX();
@@ -715,7 +732,7 @@ struct UpcCandProducerGlobalMuon {
715732
for (const auto& ianchor : vAnchorIds) {
716733
if (hasVertex) {
717734
const auto& trk = fwdTracks.iteratorAt(ianchor);
718-
auto dca = propagateGlobalToDCA(trk, bestVtxX, bestVtxY, bestVtxZ);
735+
auto dca = propagateFwdToDCA(trk, bestVtxX, bestVtxY, bestVtxZ);
719736
histRegistry.fill(HIST("hDCAxyGlobal"), dca[0]);
720737
histRegistry.fill(HIST("hDCAzGlobal"), dca[1]);
721738
if (dca[0] > static_cast<double>(fMaxDCAxy))
@@ -735,7 +752,7 @@ struct UpcCandProducerGlobalMuon {
735752

736753
// Fill invariant mass from MCH-MID-MFT anchors only
737754
uint16_t numContribAnch = numContrib;
738-
if (numContribAnch >= 2) {
755+
if (numContribAnch >= kMinTracksForPair) {
739756
double mass2 = sumE * sumE - sumPx * sumPx - sumPy * sumPy - sumPz * sumPz;
740757
histRegistry.fill(HIST("hMassGlobalMuon"), mass2 > 0. ? std::sqrt(mass2) : 0.);
741758
}
@@ -744,7 +761,7 @@ struct UpcCandProducerGlobalMuon {
744761
for (const auto& [imchMft, gbc] : mapMchMftIdBc) {
745762
if (hasVertex) {
746763
const auto& trk = fwdTracks.iteratorAt(imchMft);
747-
auto dca = propagateGlobalToDCA(trk, bestVtxX, bestVtxY, bestVtxZ);
764+
auto dca = propagateFwdToDCA(trk, bestVtxX, bestVtxY, bestVtxZ);
748765
histRegistry.fill(HIST("hDCAxyGlobal"), dca[0]);
749766
histRegistry.fill(HIST("hDCAzGlobal"), dca[1]);
750767
if (dca[0] > static_cast<double>(fMaxDCAxy))
@@ -763,12 +780,12 @@ struct UpcCandProducerGlobalMuon {
763780
}
764781

765782
// Fill invariant mass including MCH-MFT tracks (only if MCH-MFT tracks were added)
766-
if (numContrib > numContribAnch && numContrib >= 2) {
783+
if (numContrib > numContribAnch && numContrib >= kMinTracksForPair) {
767784
double mass2 = sumE * sumE - sumPx * sumPx - sumPy * sumPy - sumPz * sumPz;
768785
histRegistry.fill(HIST("hMassGlobalMuonWithMCHMFT"), mass2 > 0. ? std::sqrt(mass2) : 0.);
769786
}
770787

771-
if (numContrib < 1)
788+
if (numContrib < kMinTracksForCandidate)
772789
continue;
773790

774791
eventCandidates(globalBcAnchor, runNumber, bestVtxX, bestVtxY, bestVtxZ, 0, numContrib, 0, 0);
@@ -817,7 +834,7 @@ struct UpcCandProducerGlobalMuon {
817834
numContrib++;
818835
selTrackIds.push_back(imuon);
819836
}
820-
if (numContrib < 1) // didn't save any MCH-MID tracks
837+
if (numContrib < kMinTracksForCandidate) // didn't save any MCH-MID tracks
821838
continue;
822839
std::map<int64_t, uint64_t> mapMchIdBc{};
823840
getMchTrackIds(globalBcMid, mapGlobalBcsWithMCHTrackIds, fBcWindowMCH, mapMchIdBc);

0 commit comments

Comments
 (0)