From f963bf6548bfb35dd186ba91e07e9c644c8986a6 Mon Sep 17 00:00:00 2001 From: Daiki Sekihata Date: Mon, 2 Feb 2026 17:36:48 +0100 Subject: [PATCH 1/3] PWGEM/Dilepton: update muon analyses for dr vs. chi2MatchMCHMFT --- PWGEM/Dilepton/Core/Dilepton.h | 3 + PWGEM/Dilepton/Core/DileptonHadronMPC.h | 3 + PWGEM/Dilepton/Core/DileptonMC.h | 3 + PWGEM/Dilepton/Core/DimuonCut.cxx | 7 + PWGEM/Dilepton/Core/DimuonCut.h | 10 ++ PWGEM/Dilepton/Core/SingleTrackQC.h | 6 + PWGEM/Dilepton/Core/SingleTrackQCMC.h | 6 + .../TableProducer/skimmerPrimaryMuon.cxx | 147 ++++++++++++------ PWGEM/Dilepton/Tasks/matchingMFT.cxx | 95 +++++++---- 9 files changed, 202 insertions(+), 78 deletions(-) diff --git a/PWGEM/Dilepton/Core/Dilepton.h b/PWGEM/Dilepton/Core/Dilepton.h index 7c20f72802f..4f8aac33a70 100644 --- a/PWGEM/Dilepton/Core/Dilepton.h +++ b/PWGEM/Dilepton/Core/Dilepton.h @@ -289,6 +289,8 @@ struct Dilepton { Configurable cfg_max_DPhi_wrt_matchedMCHMID{"cfg_max_DPhi_wrt_matchedMCHMID", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to apply MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; + Configurable cfg_slope_dr_chi2MatchMFTMCH{"cfg_slope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; + Configurable cfg_intercept_dr_chi2MatchMFTMCH{"cfg_intercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; } dimuoncuts; struct : ConfigurableGroup { @@ -809,6 +811,7 @@ struct Dilepton { fDimuonCut.SetMaxPDCARabsDep([&](float rabs) { return (rabs < 26.5 ? 594.f : 324.f); }); fDimuonCut.SetMaxdPtdEtadPhiwrtMCHMID(dimuoncuts.cfg_max_relDPt_wrt_matchedMCHMID, dimuoncuts.cfg_max_DEta_wrt_matchedMCHMID, dimuoncuts.cfg_max_DPhi_wrt_matchedMCHMID); // this is relevant for global muons fDimuonCut.SetMFTHitMap(dimuoncuts.requireMFTHitMap, dimuoncuts.requiredMFTDisks); + fDimuonCut.SetSlopeAndInterceptDRvsChi2MCHMFT(dimuoncuts.cfg_slope_dr_chi2MatchMFTMCH, dimuoncuts.cfg_intercept_dr_chi2MatchMFTMCH); } template diff --git a/PWGEM/Dilepton/Core/DileptonHadronMPC.h b/PWGEM/Dilepton/Core/DileptonHadronMPC.h index 856d5d7f740..ca985aaa1d3 100644 --- a/PWGEM/Dilepton/Core/DileptonHadronMPC.h +++ b/PWGEM/Dilepton/Core/DileptonHadronMPC.h @@ -271,6 +271,8 @@ struct DileptonHadronMPC { Configurable cfg_max_DPhi_wrt_matchedMCHMID{"cfg_max_DPhi_wrt_matchedMCHMID", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to apply MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; + Configurable cfg_slope_dr_chi2MatchMFTMCH{"cfg_slope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; + Configurable cfg_intercept_dr_chi2MatchMFTMCH{"cfg_intercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; } dimuoncuts; EMTrackCut fEMTrackCut; @@ -692,6 +694,7 @@ struct DileptonHadronMPC { fDimuonCut.SetMaxPDCARabsDep([&](float rabs) { return (rabs < 26.5 ? 594.f : 324.f); }); fDimuonCut.SetMaxdPtdEtadPhiwrtMCHMID(dimuoncuts.cfg_max_relDPt_wrt_matchedMCHMID, dimuoncuts.cfg_max_DEta_wrt_matchedMCHMID, dimuoncuts.cfg_max_DPhi_wrt_matchedMCHMID); // this is relevant for global muons fDimuonCut.SetMFTHitMap(dimuoncuts.requireMFTHitMap, dimuoncuts.requiredMFTDisks); + fDimuonCut.SetSlopeAndInterceptDRvsChi2MCHMFT(dimuoncuts.cfg_slope_dr_chi2MatchMFTMCH, dimuoncuts.cfg_intercept_dr_chi2MatchMFTMCH); } void DefineEMTrackCut() diff --git a/PWGEM/Dilepton/Core/DileptonMC.h b/PWGEM/Dilepton/Core/DileptonMC.h index c5b976f0405..3af734dafd3 100644 --- a/PWGEM/Dilepton/Core/DileptonMC.h +++ b/PWGEM/Dilepton/Core/DileptonMC.h @@ -289,6 +289,8 @@ struct DileptonMC { Configurable cfg_max_DPhi_wrt_matchedMCHMID{"cfg_max_DPhi_wrt_matchedMCHMID", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to apply MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; + Configurable cfg_slope_dr_chi2MatchMFTMCH{"cfg_slope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; + Configurable cfg_intercept_dr_chi2MatchMFTMCH{"cfg_intercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; Configurable rejectWrongMatch{"rejectWrongMatch", false, "flag to reject wrong match between MFT and MCH-MID"}; // this is only for MC study, as we don't know correct match in data. } dimuoncuts; @@ -814,6 +816,7 @@ struct DileptonMC { fDimuonCut.SetMaxPDCARabsDep([&](float rabs) { return (rabs < 26.5 ? 594.f : 324.f); }); fDimuonCut.SetMaxdPtdEtadPhiwrtMCHMID(dimuoncuts.cfg_max_relDPt_wrt_matchedMCHMID, dimuoncuts.cfg_max_DEta_wrt_matchedMCHMID, dimuoncuts.cfg_max_DPhi_wrt_matchedMCHMID); // this is relevant for global muons fDimuonCut.SetMFTHitMap(dimuoncuts.requireMFTHitMap, dimuoncuts.requiredMFTDisks); + fDimuonCut.SetSlopeAndInterceptDRvsChi2MCHMFT(dimuoncuts.cfg_slope_dr_chi2MatchMFTMCH, dimuoncuts.cfg_intercept_dr_chi2MatchMFTMCH); } template diff --git a/PWGEM/Dilepton/Core/DimuonCut.cxx b/PWGEM/Dilepton/Core/DimuonCut.cxx index c9cded44408..3371c163129 100644 --- a/PWGEM/Dilepton/Core/DimuonCut.cxx +++ b/PWGEM/Dilepton/Core/DimuonCut.cxx @@ -147,3 +147,10 @@ void DimuonCut::SetMaxdPtdEtadPhiwrtMCHMID(float reldPtMax, float dEtaMax, float LOG(info) << "Dimuon Cut, set max deta between MFT-MCH-MID and associated MCH-MID: " << mMaxdEtawrtMCHMID; LOG(info) << "Dimuon Cut, set max dphi between MFT-MCH-MID and associated MCH-MID: " << mMaxdPhiwrtMCHMID; } +void DimuonCut::SetSlopeAndInterceptDRvsChi2MCHMFT(float slope, float intercept) +{ + mSlope_dr_chi2MatchMFTMCH = slope; + mIntercept_dr_chi2MatchMFTMCH = intercept; + LOG(info) << "Dimuon Cut, set slope between dr and chi2MCHMFT: " << mSlope_dr_chi2MatchMFTMCH; + LOG(info) << "Dimuon Cut, set intercept between dr and chi2MCHMFT: " << mIntercept_dr_chi2MatchMFTMCH; +} diff --git a/PWGEM/Dilepton/Core/DimuonCut.h b/PWGEM/Dilepton/Core/DimuonCut.h index 80968cc9e81..3ffd5198114 100644 --- a/PWGEM/Dilepton/Core/DimuonCut.h +++ b/PWGEM/Dilepton/Core/DimuonCut.h @@ -64,6 +64,7 @@ class DimuonCut : public TNamed kPDCA, kMFTHitMap, kDPtDEtaDPhiwrtMCHMID, + kDr_MatchingChi2MCHMFT_2D, kNCuts }; @@ -172,6 +173,9 @@ class DimuonCut : public TNamed if (track.trackType() == static_cast(o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) && !IsSelectedTrack(track, DimuonCuts::kDPtDEtaDPhiwrtMCHMID)) { return false; } + if (track.trackType() == static_cast(o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) && !IsSelectedTrack(track, DimuonCuts::kDr_MatchingChi2MCHMFT_2D)) { + return false; + } return true; } @@ -232,6 +236,9 @@ class DimuonCut : public TNamed case DimuonCuts::kDPtDEtaDPhiwrtMCHMID: return std::fabs(track.ptMatchedMCHMID() - track.pt()) / track.pt() < mMaxReldPtwrtMCHMID && std::sqrt(std::pow((track.etaMatchedMCHMID() - track.eta()) / mMaxdEtawrtMCHMID, 2) + std::pow((track.phiMatchedMCHMID() - track.phi()) / mMaxdPhiwrtMCHMID, 2)) < 1.f; + case DimuonCuts::kDr_MatchingChi2MCHMFT_2D: + return mSlope_dr_chi2MatchMFTMCH * track.chi2MatchMCHMFT() + mIntercept_dr_chi2MatchMFTMCH < std::sqrt(std::pow(track.etaMatchedMCHMID() - track.eta(), 2) + std::pow(track.phiMatchedMCHMID() - track.phi(), 2)); + default: return false; } @@ -259,6 +266,7 @@ class DimuonCut : public TNamed void SetMaxPDCARabsDep(std::function RabsDepCut); void SetMFTHitMap(bool flag, std::vector hitMap); void SetMaxdPtdEtadPhiwrtMCHMID(float reldPtMax, float dEtaMax, float dPhiMax); // this is relevant for global muons + void SetSlopeAndInterceptDRvsChi2MCHMFT(float slope, float intercept); // this is relevant for global muons private: // pair cuts @@ -290,6 +298,8 @@ class DimuonCut : public TNamed float mMaxReldPtwrtMCHMID{1e10f}, mMaxdEtawrtMCHMID{1e10f}, mMaxdPhiwrtMCHMID{1e10f}; bool mApplyMFTHitMap{false}; std::vector mRequiredMFTDisks{}; + float mSlope_dr_chi2MatchMFTMCH{-0.15 / 30}; + float mIntercept_dr_chi2MatchMFTMCH{1e+10}; ClassDef(DimuonCut, 1); }; diff --git a/PWGEM/Dilepton/Core/SingleTrackQC.h b/PWGEM/Dilepton/Core/SingleTrackQC.h index d9691d54408..7a230d072b9 100644 --- a/PWGEM/Dilepton/Core/SingleTrackQC.h +++ b/PWGEM/Dilepton/Core/SingleTrackQC.h @@ -211,6 +211,8 @@ struct SingleTrackQC { Configurable cfg_max_DPhi_wrt_matchedMCHMID{"cfg_max_DPhi_wrt_matchedMCHMID", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to apply MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; + Configurable cfg_slope_dr_chi2MatchMFTMCH{"cfg_slope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; + Configurable cfg_intercept_dr_chi2MatchMFTMCH{"cfg_intercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; } dimuoncuts; struct : ConfigurableGroup { @@ -320,6 +322,7 @@ struct SingleTrackQC { fRegistry.add("Track/positive/hChi2MatchMCHMID", "chi2 match MCH-MID;chi2", kTH1F, {{100, 0.0f, 100}}, false); fRegistry.add("Track/positive/hChi2MatchMCHMFT", "chi2 match MCH-MFT;chi2", kTH1F, {{100, 0.0f, 100}}, false); fRegistry.add("Track/positive/hMFTClusterMap", "MFT cluster map", kTH1F, {{1024, -0.5, 1023.5}}, false); + fRegistry.add("Track/positive/hdR_Chi2MatchMCHMFT", "dr vs. matching chi2 MCH-MFT;chi2 match MCH-MFT;#DeltaR;", kTH2F, {{200, 0, 50}, {200, 0, 0.5}}, false); fRegistry.addClone("Track/positive/", "Track/negative/"); } } @@ -541,6 +544,7 @@ struct SingleTrackQC { fDimuonCut.SetMaxPDCARabsDep([&](float rabs) { return (rabs < 26.5 ? 594.f : 324.f); }); fDimuonCut.SetMaxdPtdEtadPhiwrtMCHMID(dimuoncuts.cfg_max_relDPt_wrt_matchedMCHMID, dimuoncuts.cfg_max_DEta_wrt_matchedMCHMID, dimuoncuts.cfg_max_DPhi_wrt_matchedMCHMID); // this is relevant for global muons fDimuonCut.SetMFTHitMap(dimuoncuts.requireMFTHitMap, dimuoncuts.requiredMFTDisks); + fDimuonCut.SetSlopeAndInterceptDRvsChi2MCHMFT(dimuoncuts.cfg_slope_dr_chi2MatchMFTMCH, dimuoncuts.cfg_intercept_dr_chi2MatchMFTMCH); } template @@ -675,6 +679,7 @@ struct SingleTrackQC { fRegistry.fill(HIST("Track/positive/hChi2MatchMCHMID"), track.chi2MatchMCHMID()); fRegistry.fill(HIST("Track/positive/hChi2MatchMCHMFT"), track.chi2MatchMCHMFT()); fRegistry.fill(HIST("Track/positive/hMFTClusterMap"), track.mftClusterMap()); + fRegistry.fill(HIST("Track/positive/hdR_Chi2MatchMCHMFT"), track.chi2MatchMCHMFT(), std::sqrt(deta * deta + dphi * dphi)); } else { fRegistry.fill(HIST("Track/negative/hs"), track.pt(), track.eta(), track.phi(), dca_xy, weight); fRegistry.fill(HIST("Track/negative/hEtaPhi_MatchMCHMID"), track.phiMatchedMCHMID(), track.etaMatchedMCHMID(), weight); @@ -698,6 +703,7 @@ struct SingleTrackQC { fRegistry.fill(HIST("Track/negative/hChi2MatchMCHMID"), track.chi2MatchMCHMID()); fRegistry.fill(HIST("Track/negative/hChi2MatchMCHMFT"), track.chi2MatchMCHMFT()); fRegistry.fill(HIST("Track/negative/hMFTClusterMap"), track.mftClusterMap()); + fRegistry.fill(HIST("Track/negative/hdR_Chi2MatchMCHMFT"), track.chi2MatchMCHMFT(), std::sqrt(deta * deta + dphi * dphi)); } } diff --git a/PWGEM/Dilepton/Core/SingleTrackQCMC.h b/PWGEM/Dilepton/Core/SingleTrackQCMC.h index 78fc3481cf2..57c51ccae4e 100644 --- a/PWGEM/Dilepton/Core/SingleTrackQCMC.h +++ b/PWGEM/Dilepton/Core/SingleTrackQCMC.h @@ -218,6 +218,8 @@ struct SingleTrackQCMC { Configurable cfg_max_DPhi_wrt_matchedMCHMID{"cfg_max_DPhi_wrt_matchedMCHMID", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to apply MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; + Configurable cfg_slope_dr_chi2MatchMFTMCH{"cfg_slope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; + Configurable cfg_intercept_dr_chi2MatchMFTMCH{"cfg_intercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; Configurable rejectWrongMatch{"rejectWrongMatch", false, "flag to reject wrong match between MFT and MCH-MID"}; // this is only for MC study, as we don't know correct match in data. } dimuoncuts; @@ -366,6 +368,7 @@ struct SingleTrackQCMC { fRegistry.add("Track/PromptLF/positive/hPtGen_DeltaPtOverPtGen", "muon p_{T} resolution;p_{T}^{gen} (GeV/c);(p_{T}^{rec} - p_{T}^{gen})/p_{T}^{gen}", kTH2F, {{200, 0, 10}, {200, -1.0f, 1.0f}}, true); fRegistry.add("Track/PromptLF/positive/hPtGen_DeltaEta", "muon #eta resolution;p_{T}^{gen} (GeV/c);#eta^{rec} - #eta^{gen}", kTH2F, {{200, 0, 10}, {100, -0.05f, 0.05f}}, true); fRegistry.add("Track/PromptLF/positive/hPtGen_DeltaPhi", "muon #varphi resolution;p_{T}^{gen} (GeV/c);#varphi^{rec} - #varphi^{gen} (rad.)", kTH2F, {{200, 0, 10}, {100, -0.05f, 0.05f}}, true); + fRegistry.add("Track/PromptLF/positive/hdR_Chi2MatchMCHMFT", "dr vs. matching chi2 MCH-MFT;chi2 match MCH-MFT;#DeltaR;", kTH2F, {{200, 0, 50}, {200, 0, 0.5}}, false); } fRegistry.addClone("Track/PromptLF/positive/", "Track/PromptLF/negative/"); fRegistry.addClone("Track/PromptLF/", "Track/NonPromptLF/"); @@ -584,6 +587,7 @@ struct SingleTrackQCMC { fDimuonCut.SetMaxPDCARabsDep([&](float rabs) { return (rabs < 26.5 ? 594.f : 324.f); }); fDimuonCut.SetMaxdPtdEtadPhiwrtMCHMID(dimuoncuts.cfg_max_relDPt_wrt_matchedMCHMID, dimuoncuts.cfg_max_DEta_wrt_matchedMCHMID, dimuoncuts.cfg_max_DPhi_wrt_matchedMCHMID); // this is relevant for global muons fDimuonCut.SetMFTHitMap(dimuoncuts.requireMFTHitMap, dimuoncuts.requiredMFTDisks); + fDimuonCut.SetSlopeAndInterceptDRvsChi2MCHMFT(dimuoncuts.cfg_slope_dr_chi2MatchMFTMCH, dimuoncuts.cfg_intercept_dr_chi2MatchMFTMCH); } template @@ -781,6 +785,7 @@ struct SingleTrackQCMC { fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hChi2MatchMCHMID"), track.chi2MatchMCHMID()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hChi2MatchMCHMFT"), track.chi2MatchMCHMFT()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hMFTClusterMap"), track.mftClusterMap()); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hdR_Chi2MatchMCHMFT"), track.chi2MatchMCHMFT(), std::sqrt(deta * deta + dphi * dphi)); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hPtGen_DeltaPtOverPtGen"), mctrack.pt(), (track.pt() - mctrack.pt()) / mctrack.pt()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hPtGen_DeltaEta"), mctrack.pt(), track.eta() - mctrack.eta()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hPtGen_DeltaPhi"), mctrack.pt(), track.phi() - mctrack.phi()); @@ -813,6 +818,7 @@ struct SingleTrackQCMC { fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hChi2MatchMCHMID"), track.chi2MatchMCHMID()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hChi2MatchMCHMFT"), track.chi2MatchMCHMFT()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hMFTClusterMap"), track.mftClusterMap()); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hdR_Chi2MatchMCHMFT"), track.chi2MatchMCHMFT(), std::sqrt(deta * deta + dphi * dphi)); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hPtGen_DeltaPtOverPtGen"), mctrack.pt(), (track.pt() - mctrack.pt()) / mctrack.pt()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hPtGen_DeltaEta"), mctrack.pt(), track.eta() - mctrack.eta()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hPtGen_DeltaPhi"), mctrack.pt(), track.phi() - mctrack.phi()); diff --git a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx index df73ca2391c..e353da03a88 100644 --- a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx +++ b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx @@ -92,8 +92,11 @@ struct skimmerPrimaryMuon { Configurable minNmuon{"minNmuon", 0, "min number of muon candidates per collision"}; Configurable maxDEta{"maxDEta", 1e+10f, "max. deta between MFT-MCH-MID and MCH-MID"}; Configurable maxDPhi{"maxDPhi", 1e+10f, "max. dphi between MFT-MCH-MID and MCH-MID"}; + Configurable cfgApplyPreselectionInBestMatch{"cfgApplyPreselectionInBestMatch", false, "flag to apply preselection in find best match function"}; Configurable cfgSlope_dr_chi2MatchMFTMCH{"cfgSlope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; Configurable cfgIntercept_dr_chi2MatchMFTMCH{"cfgIntercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; + Configurable cfgPeakPosition_chi2MatchMFTMCH{"cfgPeakPosition_chi2MatchMFTMCH", 0.f, "peak position of chiMatchMCHMFT distribution"}; // 2 + Configurable cfgPeakPosition_dr{"cfgPeakPosition_dr", 0.f, "peak position of dr distribution"}; // 0.01 // for z shift for propagation Configurable cfgApplyZShiftFromCCDB{"cfgApplyZShiftFromCCDB", false, "flag to apply z shift"}; @@ -308,6 +311,9 @@ struct skimmerPrimaryMuon { float etaMatchedMFTatMP = 999.f; float phiMatchedMFTatMP = 999.f; + float deta = 999.f; + float dphi = 999.f; + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { // apply r-absorber cut here to minimize the number of calling propagateMuon. if (fwdtrack.rAtAbsorberEnd() < minRabsGL || maxRabs < fwdtrack.rAtAbsorberEnd()) { @@ -389,6 +395,19 @@ struct skimmerPrimaryMuon { sigma_dcaXY = dcaXY / dcaXYinSigma; } + deta = etaMatchedMCHMID - eta; + dphi = phiMatchedMCHMID - phi; + o2::math_utils::bringToPMPi(dphi); + + if (std::sqrt(std::pow(deta / maxDEta, 2) + std::pow(dphi / maxDPhi, 2)) > 1.f) { + return false; + } + + float dr = std::sqrt(deta * deta + dphi * dphi); + if (cfgSlope_dr_chi2MatchMFTMCH * fwdtrack.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { + return false; + } + if (refitGlobalMuon) { pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); } @@ -424,14 +443,6 @@ struct skimmerPrimaryMuon { return false; } - float deta = etaMatchedMCHMID - eta; - float dphi = phiMatchedMCHMID - phi; - o2::math_utils::bringToPMPi(dphi); - - if (std::sqrt(std::pow(deta / maxDEta, 2) + std::pow(dphi / maxDPhi, 2)) > 1.f) { - return false; - } - if constexpr (fillTable) { float dpt = (ptMatchedMCHMID - pt) / pt; @@ -530,7 +541,8 @@ struct skimmerPrimaryMuon { std::unordered_map map_mfttrackcovs; std::vector> vec_min_chi2MatchMCHMFT; // std::pair -> chi2MatchMCHMFT; - std::vector> vec_min_dr; // std::pair -> dr; + // std::vector> vec_min_dr; // std::pair -> dr; + // std::vector> vec_min_2d; // std::pair -> dr + chi2MatchMCHMFT; template void findBestMatchPerMCHMID(TCollision const& collision, TFwdTrack const& fwdtrack, TFwdTracks const& fwdtracks, TMFTTracks const&) @@ -553,9 +565,10 @@ struct skimmerPrimaryMuon { float phiMatchedMCHMID = propmuonAtPV_Matched.getPhi(); o2::math_utils::bringTo02Pi(phiMatchedMCHMID); - float min_chi2MatchMCHMFT = 1e+10, min_dr = 1e+10; + float min_chi2MatchMCHMFT = 1e+10, min_dr = 1e+10, min_distance_2d = 1e+10; std::tuple tupleIds_at_min_chi2mftmch; std::tuple tupleIds_at_min_dr; + std::tuple tupleIds_at_min_distance_2d; for (const auto& muon_tmp : muons_per_MCHMID) { if (muon_tmp.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { auto tupleId = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId()); @@ -571,30 +584,39 @@ struct skimmerPrimaryMuon { float dphi = phiMatchedMCHMID - phi; o2::math_utils::bringToPMPi(dphi); float dr = std::sqrt(deta * deta + dphi * dphi); - int ndf = 2 * (mchtrack.nClusters() + mfttrack.nClusters()) - 5; + if (cfgApplyPreselectionInBestMatch && cfgSlope_dr_chi2MatchMFTMCH * muon_tmp.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { + continue; + } + fRegistry.fill(HIST("MFTMCHMID/hdR_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), dr); fRegistry.fill(HIST("MFTMCHMID/hdR_Chi2"), muon_tmp.chi2() / ndf, dr); fRegistry.fill(HIST("MFTMCHMID/hChi2_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), muon_tmp.chi2() / ndf); - if (cfgSlope_dr_chi2MatchMFTMCH * muon_tmp.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { - continue; - } // LOGF(info, "muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f", muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId(), muon_tmp.chi2MatchMCHMFT()); - if (0.f < muon_tmp.chi2MatchMCHMFT() && muon_tmp.chi2MatchMCHMFT() < min_chi2MatchMCHMFT) { - min_chi2MatchMCHMFT = muon_tmp.chi2MatchMCHMFT(); - tupleIds_at_min_chi2mftmch = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId()); - } - if (dr < min_dr) { - min_dr = dr; - tupleIds_at_min_dr = tupleId; + if (0.f < muon_tmp.chi2MatchMCHMFT() && std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2)) < min_chi2MatchMCHMFT) { + min_chi2MatchMCHMFT = std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2)); + tupleIds_at_min_chi2mftmch = tupleId; } + + // if (std::sqrt(std::pow(dr - cfgPeakPosition_dr, 2)) < min_dr) { + // min_dr = std::sqrt(std::pow(dr - cfgPeakPosition_dr, 2)); + // tupleIds_at_min_dr = tupleId; + // } + + // float distance_2d = std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2) + std::pow(dr - cfgPeakPosition_dr, 2)); + // if (distance_2d < min_distance_2d) { + // min_distance_2d = distance_2d; + // tupleIds_at_min_distance_2d = tupleId; + // } } } vec_min_chi2MatchMCHMFT.emplace_back(tupleIds_at_min_chi2mftmch); - vec_min_dr.emplace_back(tupleIds_at_min_dr); + // vec_min_dr.emplace_back(tupleIds_at_min_dr); + // vec_min_2d.emplace_back(tupleIds_at_min_distance_2d); + // LOGF(info, "min: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f", std::get<0>(tupleIds_at_min), std::get<1>(tupleIds_at_min), std::get<2>(tupleIds_at_min), min_chi2MatchMCHMFT); } @@ -609,7 +631,8 @@ struct skimmerPrimaryMuon { void processRec_SA(MyCollisions const& collisions, MyFwdTracks const& fwdtracks, aod::MFTTracks const& mfttracks, aod::BCsWithTimestamps const&) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); @@ -680,15 +703,18 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processRec_SA, "process reconstructed info", false); void processRec_TTCA(MyCollisions const& collisions, MyFwdTracks const& fwdtracks, aod::MFTTracks const& mfttracks, aod::BCsWithTimestamps const&, aod::FwdTrackAssoc const& fwdtrackIndices) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -767,8 +793,10 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processRec_TTCA, "process reconstructed info", false); @@ -779,7 +807,8 @@ struct skimmerPrimaryMuon { } vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -858,15 +887,18 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processRec_TTCA_withMFTCov, "process reconstructed info", false); void processRec_SA_SWT(MyCollisionsWithSWT const& collisions, MyFwdTracks const& fwdtracks, aod::MFTTracks const& mfttracks, aod::BCsWithTimestamps const&) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -939,15 +971,18 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processRec_SA_SWT, "process reconstructed info only with standalone", false); void processRec_TTCA_SWT(MyCollisionsWithSWT const& collisions, MyFwdTracks const& fwdtracks, aod::MFTTracks const& mfttracks, aod::BCsWithTimestamps const&, aod::FwdTrackAssoc const& fwdtrackIndices) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1028,8 +1063,10 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processRec_TTCA_SWT, "process reconstructed info", false); @@ -1039,7 +1076,8 @@ struct skimmerPrimaryMuon { map_mfttrackcovs[mfttrackConv.matchMFTTrackId()] = mfttrackConv.globalIndex(); } vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1120,15 +1158,18 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processRec_TTCA_SWT_withMFTCov, "process reconstructed info", false); void processMC_SA(soa::Join const& collisions, MyFwdTracksMC const& fwdtracks, MFTTracksMC const& mfttracks, aod::BCsWithTimestamps const&) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1202,15 +1243,18 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processMC_SA, "process reconstructed and MC info", false); void processMC_TTCA(soa::Join const& collisions, MyFwdTracksMC const& fwdtracks, MFTTracksMC const& mfttracks, aod::BCsWithTimestamps const&, aod::FwdTrackAssoc const& fwdtrackIndices) { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1294,8 +1338,10 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processMC_TTCA, "process reconstructed and MC info", false); @@ -1305,7 +1351,8 @@ struct skimmerPrimaryMuon { map_mfttrackcovs[mfttrackConv.matchMFTTrackId()] = mfttrackConv.globalIndex(); } vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); - vec_min_dr.reserve(fwdtracks.size()); + // vec_min_dr.reserve(fwdtracks.size()); + // vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1389,8 +1436,10 @@ struct skimmerPrimaryMuon { map_mfttrackcovs.clear(); vec_min_chi2MatchMCHMFT.clear(); vec_min_chi2MatchMCHMFT.shrink_to_fit(); - vec_min_dr.clear(); - vec_min_dr.shrink_to_fit(); + // vec_min_dr.clear(); + // vec_min_dr.shrink_to_fit(); + // vec_min_2d.clear(); + // vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(skimmerPrimaryMuon, processMC_TTCA_withMFTCov, "process reconstructed and MC with MFTCov info", false); diff --git a/PWGEM/Dilepton/Tasks/matchingMFT.cxx b/PWGEM/Dilepton/Tasks/matchingMFT.cxx index de8dfc1d722..558d0f418d0 100644 --- a/PWGEM/Dilepton/Tasks/matchingMFT.cxx +++ b/PWGEM/Dilepton/Tasks/matchingMFT.cxx @@ -96,9 +96,12 @@ struct matchingMFT { Configurable requireMFTHitMap{"requireMFTHitMap", false, "flag to require MFT hit map"}; Configurable> requiredMFTDisks{"requiredMFTDisks", std::vector{0}, "hit map on MFT disks [0,1,2,3,4]. logical-OR of each double-sided disk"}; Configurable matchingZ{"matchingZ", -77.5, "z position where matching is performed"}; - Configurable cfgBestMatchFinder{"cfgBestMatchFinder", 0, "matching chi2:0, dr:1"}; + Configurable cfgBestMatchFinder{"cfgBestMatchFinder", 0, "matching chi2:0, dr:1, 2d:2"}; + Configurable cfgApplyPreselectionInBestMatch{"cfgApplyPreselectionInBestMatch", false, "flag to apply preselection in find best match function"}; Configurable cfgSlope_dr_chi2MatchMFTMCH{"cfgSlope_dr_chi2MatchMFTMCH", -0.15 / 30, "slope of chiMatchMCHMFT vs. dR"}; Configurable cfgIntercept_dr_chi2MatchMFTMCH{"cfgIntercept_dr_chi2MatchMFTMCH", 1e+10f, "intercept of chiMatchMCHMFT vs. dR"}; + Configurable cfgPeakPosition_chi2MatchMFTMCH{"cfgPeakPosition_chi2MatchMFTMCH", 2.f, "peak position of chiMatchMCHMFT distribution"}; + Configurable cfgPeakPosition_dr{"cfgPeakPosition_dr", 0.01, "peak position of dr distribution"}; Configurable cfgCentEstimator{"cfgCentEstimator", 2, "FT0M:0, FT0A:1, FT0C:2"}; Configurable cfgCentMin{"cfgCentMin", -1.f, "min. centrality"}; @@ -228,9 +231,7 @@ struct matchingMFT { fRegistry.add("MFTMCHMID/primary/correct/hDCAxyz", "DCA xy vs. z;DCA_{xy} (cm);DCA_{z} (cm);", kTH2F, {{100, 0, 1}, {200, -0.1, 0.1}}, false); fRegistry.add("MFTMCHMID/primary/correct/hMCHBitMap", "MCH bit map;MCH bit map", kTH1F, {{1024, -0.5, 1023.5}}, false); fRegistry.add("MFTMCHMID/primary/correct/hMIDBitMap", "MID bit map;MID bit map", kTH1F, {{256, -0.5, 255.5}}, false); - fRegistry.add("MFTMCHMID/primary/correct/hdR_Chi2MatchMCHMFT", "dr vs. matching chi2 MCH-MFT;chi2 match MCH-MFT;#DeltaR", kTH2F, {{200, 0, 50}, {200, 0, 0.5}}, false); - fRegistry.add("MFTMCHMID/primary/correct/hdR_Chi2", "dr vs. chi2;global chi2/ndf;#DeltaR", kTH2F, {{100, 0, 10}, {200, 0, 0.5}}, false); - fRegistry.add("MFTMCHMID/primary/correct/hChi2_Chi2MatchMCHMFT", "chi2 vs. matching chi2 MCH-MFT;chi2 match MCH-MFT;global chi2/ndf", kTH2F, {{200, 0, 50}, {100, 0, 10}}, false); + fRegistry.add("MFTMCHMID/primary/correct/hdR_Chi2MatchMCHMFT", "dr vs. matching chi2 MCH-MFT;chi2 match MCH-MFT;#DeltaR;", kTH2F, {{200, 0, 50}, {200, 0, 0.5}}, false); fRegistry.add("MFTMCHMID/primary/correct/hCorrectAsocc", "correct fwdtrack-to-collision association", kTH1F, {{2, -0.5, +1.5}}, false); fRegistry.add("MFTMCHMID/primary/correct/hIsCA", "cellular automaton;isCA", kTH1F, {{2, -0.5, 1.5}}, false); fRegistry.add("MFTMCHMID/primary/correct/hProdVtxZ", "prod. vtx Z of muon;V_{z} (cm)", kTH1F, {{200, -100, 100}}, false); @@ -284,6 +285,7 @@ struct matchingMFT { if (chi2_per_ndf < 0.f || maxChi2GL < chi2_per_ndf) { return false; } + } else if (trackType == static_cast(o2::aod::fwdtrack::ForwardTrackTypeEnum::MuonStandaloneTrack)) { if (eta < minEtaSA || maxEtaSA < eta) { return false; @@ -539,6 +541,8 @@ struct matchingMFT { float deta = etaMatchedMCHMID - eta; float dphi = phiMatchedMCHMID - phi; o2::math_utils::bringToPMPi(dphi); + float dr = std::sqrt(deta * deta + dphi * dphi); + if (std::sqrt(std::pow(deta / maxDEta, 2) + std::pow(dphi / maxDPhi, 2)) > 1.f) { return; } @@ -549,6 +553,10 @@ struct matchingMFT { return; } + if (cfgSlope_dr_chi2MatchMFTMCH * fwdtrack.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { + return; + } + if (!isSelected(pt, eta, rAtAbsorberEnd, pDCA, fwdtrack.chi2() / (2.f * (mchtrack.nClusters() + nClustersMFT) - 5.f), fwdtrack.trackType(), dcaXY)) { return; } @@ -746,6 +754,7 @@ struct matchingMFT { std::vector> vec_min_chi2MatchMCHMFT; // std::pair -> chi2MatchMCHMFT; std::vector> vec_min_dr; // std::pair -> deta + dphi; + std::vector> vec_min_2d; // std::pair -> deta + dphi; std::map, bool> mapCorrectMatch; template @@ -760,8 +769,10 @@ struct matchingMFT { std::tuple tupleIds_at_min_chi2mftmch; std::tuple tupleIds_at_min_dr; + std::tuple tupleIds_at_min_distance_2d; float min_chi2MatchMCHMFT = 1e+10; float min_dr = 1e+10; + float min_distance_2d = 1e+10; auto muons_per_MCHMID = fwdtracks.sliceBy(fwdtracksPerMCHTrack, fwdtrack.globalIndex()); // LOGF(info, "muons_per_MCHMID.size() = %d", muons_per_MCHMID.size()); @@ -772,56 +783,53 @@ struct matchingMFT { for (const auto& muon_tmp : muons_per_MCHMID) { if (muon_tmp.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { - auto tupleId = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId()); auto mchtrack = muon_tmp.template matchMCHTrack_as(); // MCH-MID auto mfttrack = muon_tmp.template matchMFTTrack_as(); + // LOGF(info, "muon_tmp.has_mcParticle() = %d, mchtrack.has_mcParticle() = %d, mfttrack.has_mcParticle() = %d", muon_tmp.has_mcParticle(), mchtrack.has_mcParticle(), mfttrack.has_mcParticle()); + if (!muon_tmp.has_mcParticle() || !mchtrack.has_mcParticle() || !mfttrack.has_mcParticle()) { continue; } + // auto mcParticle_MFTMCHMID = muon_tmp.template mcParticle_as(); // this is identical to mcParticle_MCHMID + auto mcParticle_MCHMID = mchtrack.template mcParticle_as(); // this is identical to mcParticle_MFTMCHMID + auto mcParticle_MFT = mfttrack.template mcParticle_as(); + + bool isPrimary = mcParticle_MCHMID.isPhysicalPrimary() || mcParticle_MCHMID.producedByGenerator(); + bool isMatched = (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) && (mcParticle_MFT.mcCollisionId() == mcParticle_MCHMID.mcCollisionId()); + o2::dataformats::GlobalFwdTrack propmuonAtPV = propagateMuon(muon_tmp, muon_tmp, collision, propagationPoint::kToVertex, matchingZ, mBz, mZShift); + // float pt = propmuonAtPV.getPt(); float eta = propmuonAtPV.getEta(); float phi = propmuonAtPV.getPhi(); o2::math_utils::bringTo02Pi(phi); + // if (refitGlobalMuon) { + // pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); + // } float deta = etaMatchedMCHMID - eta; float dphi = phiMatchedMCHMID - phi; o2::math_utils::bringToPMPi(dphi); float dr = std::sqrt(deta * deta + dphi * dphi); + // float chi2ndf = muon_tmp.chi2() / (2.f * (mchtrack.nClusters() + mfttrack.nClusters()) - 5.f); - if (cfgSlope_dr_chi2MatchMFTMCH * muon_tmp.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { + if (cfgApplyPreselectionInBestMatch && cfgSlope_dr_chi2MatchMFTMCH * muon_tmp.chi2MatchMCHMFT() + cfgIntercept_dr_chi2MatchMFTMCH < dr) { continue; } - // auto mcParticle_MFTMCHMID = muon_tmp.template mcParticle_as(); // this is identical to mcParticle_MCHMID - auto mcParticle_MCHMID = mchtrack.template mcParticle_as(); // this is identical to mcParticle_MFTMCHMID - auto mcParticle_MFT = mfttrack.template mcParticle_as(); - float chi2ndf = muon_tmp.chi2() / (2.f * (mchtrack.nClusters() + mfttrack.nClusters()) - 5.f); - - bool isPrimary = mcParticle_MCHMID.isPhysicalPrimary() || mcParticle_MCHMID.producedByGenerator(); - bool isMatched = (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) && (mcParticle_MFT.mcCollisionId() == mcParticle_MCHMID.mcCollisionId()); - if (isPrimary) { if (isMatched) { fRegistry.fill(HIST("MFTMCHMID/primary/correct/hdR_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), dr); - fRegistry.fill(HIST("MFTMCHMID/primary/correct/hdR_Chi2"), chi2ndf, dr); - fRegistry.fill(HIST("MFTMCHMID/primary/correct/hChi2_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), chi2ndf); } else { fRegistry.fill(HIST("MFTMCHMID/primary/wrong/hdR_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), dr); - fRegistry.fill(HIST("MFTMCHMID/primary/wrong/hdR_Chi2"), chi2ndf, dr); - fRegistry.fill(HIST("MFTMCHMID/primary/wrong/hChi2_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), chi2ndf); } } else { if (isMatched) { fRegistry.fill(HIST("MFTMCHMID/secondary/correct/hdR_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), dr); - fRegistry.fill(HIST("MFTMCHMID/secondary/correct/hdR_Chi2"), chi2ndf, dr); - fRegistry.fill(HIST("MFTMCHMID/secondary/correct/hChi2_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), chi2ndf); } else { fRegistry.fill(HIST("MFTMCHMID/secondary/wrong/hdR_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), dr); - fRegistry.fill(HIST("MFTMCHMID/secondary/wrong/hdR_Chi2"), chi2ndf, dr); - fRegistry.fill(HIST("MFTMCHMID/secondary/wrong/hChi2_Chi2MatchMCHMFT"), muon_tmp.chi2MatchMCHMFT(), chi2ndf); } } @@ -833,32 +841,40 @@ struct matchingMFT { // if (std::abs(mcParticle_MCHMID.pdgCode()) == 13 && mcParticle_MCHMID.isPhysicalPrimary()) { // if (mcParticle_MFT.globalIndex() == mcParticle_MCHMID.globalIndex()) { - // LOGF(info, "This is correct match between MFT and MCH-MID: muon_tmp.globalIndex() = %d, chi2/ndf = %f, matching chi2/ndf = %f, mcParticle.pt() = %f, mcParticle.eta() = %f, mcParticle.phi() = %f, reldpt = %f, deta = %f, dphi = %f, dr = %f", muon_tmp.globalIndex(), chi2ndf, muon_tmp.chi2MatchMCHMFT(), mcParticle_MCHMID.pt(), mcParticle_MCHMID.eta(), mcParticle_MCHMID.phi(), reldpt, deta, dphi, dr); + // LOGF(info, "This is correct match between MFT and MCH-MID: muon_tmp.globalIndex() = %d, chi2/ndf = %f, matching chi2/ndf = %f, mcParticle.pt() = %f, mcParticle.eta() = %f, mcParticle.phi() = %f, dr = %f", muon_tmp.globalIndex(), chi2ndf, muon_tmp.chi2MatchMCHMFT(), mcParticle_MCHMID.pt(), mcParticle_MCHMID.eta(), mcParticle_MCHMID.phi(), dr); // } else { - // LOGF(info, "This is wrong match between MFT and MCH-MID: muon_tmp.globalIndex() = %d, chi2/ndf = %f, matching chi2/ndf = %f , mcParticle.pt() = %f, mcParticle.eta() = %f, mcParticle.phi() = %f, reldpt = %f, deta = %f, dphi = %f, dr = %f" , muon_tmp.globalIndex(), chi2ndf, muon_tmp.chi2MatchMCHMFT(), mcParticle_MCHMID.pt(), mcParticle_MCHMID.eta(), mcParticle_MCHMID.phi(), reldpt, deta, dphi, dr); + // LOGF(info, "This is wrong match between MFT and MCH-MID: muon_tmp.globalIndex() = %d, chi2/ndf = %f, matching chi2/ndf = %f , mcParticle.pt() = %f, mcParticle.eta() = %f, mcParticle.phi() = %f, dr = %f" , muon_tmp.globalIndex(), chi2ndf, muon_tmp.chi2MatchMCHMFT(), mcParticle_MCHMID.pt(), mcParticle_MCHMID.eta(), mcParticle_MCHMID.phi(), dr); // } // } - if (0.f < muon_tmp.chi2MatchMCHMFT() && muon_tmp.chi2MatchMCHMFT() < min_chi2MatchMCHMFT) { - min_chi2MatchMCHMFT = muon_tmp.chi2MatchMCHMFT(); + if (0.f < muon_tmp.chi2MatchMCHMFT() && std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2)) < min_chi2MatchMCHMFT) { + min_chi2MatchMCHMFT = std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2)); tupleIds_at_min_chi2mftmch = tupleId; } - if (dr < min_dr) { - min_dr = dr; + if (std::sqrt(std::pow(dr - cfgPeakPosition_dr, 2)) < min_dr) { + min_dr = std::sqrt(std::pow(dr - cfgPeakPosition_dr, 2)); tupleIds_at_min_dr = tupleId; } + float distance_2d = std::sqrt(std::pow(muon_tmp.chi2MatchMCHMFT() - cfgPeakPosition_chi2MatchMFTMCH, 2) + std::pow(dr - cfgPeakPosition_dr, 2)); + if (distance_2d < min_distance_2d) { + min_distance_2d = distance_2d; + tupleIds_at_min_distance_2d = tupleId; + } + } // end of if global muon } // end of candidates loop vec_min_chi2MatchMCHMFT.emplace_back(tupleIds_at_min_chi2mftmch); vec_min_dr.emplace_back(tupleIds_at_min_dr); + vec_min_2d.emplace_back(tupleIds_at_min_distance_2d); // auto mcParticleTMP = fwdtrack.template mcParticle_as(); // this is identical to mcParticle_MFTMCHMID // if (std::abs(mcParticleTMP.pdgCode()) == 13 && mcParticleTMP.isPhysicalPrimary()) { // LOGF(info, "min chi2: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, muon_tmp.chi2MatchMCHMFT() = %f, correct match = %d", std::get<0>(tupleIds_at_min_chi2mftmch), std::get<1>(tupleIds_at_min_chi2mftmch), std::get<2>(tupleIds_at_min_chi2mftmch), min_chi2MatchMCHMFT, mapCorrectMatch[tupleIds_at_min_chi2mftmch]); - // LOGF(info, "min dr: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, dr = %f, correct match = %d", std::get<0>(tupleIds_at_min_dr), std::get<1>(tupleIds_at_min_dr), std::get<2>(tupleIds_at_min_dr), min_dr, mapCorrectMatch[tupleIds_at_min_dr]); + // LOGF(info, "min dr: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, correct match = %d", std::get<0>(tupleIds_at_min_dr), std::get<1>(tupleIds_at_min_dr), std::get<2>(tupleIds_at_min_dr), mapCorrectMatch[tupleIds_at_min_dr]); + // LOGF(info, "min 2d: muon_tmp.globalIndex() = %d, muon_tmp.matchMCHTrackId() = %d, muon_tmp.matchMFTTrackId() = %d, correct match = %d", std::get<0>(tupleIds_at_min_distance_2d), std::get<1>(tupleIds_at_min_distance_2d), std::get<2>(tupleIds_at_min_distance_2d), mapCorrectMatch[tupleIds_at_min_distance_2d]); // } } @@ -1007,6 +1023,7 @@ struct matchingMFT { { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); vec_min_dr.reserve(fwdtracks.size()); + vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1051,6 +1068,10 @@ struct matchingMFT { if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_dr.begin(), vec_min_dr.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_dr.end()) { continue; } + } else if (cfgBestMatchFinder == 2) { // 2d + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_2d.begin(), vec_min_2d.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_2d.end()) { + continue; + } } fillHistograms(collision, fwdtrack, fwdtracks, mfttracks, nullptr); @@ -1061,6 +1082,8 @@ struct matchingMFT { vec_min_chi2MatchMCHMFT.shrink_to_fit(); vec_min_dr.clear(); vec_min_dr.shrink_to_fit(); + vec_min_2d.clear(); + vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(matchingMFT, processWithoutFTTCA, "process without FTTCA", false); @@ -1068,6 +1091,7 @@ struct matchingMFT { { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); vec_min_dr.reserve(fwdtracks.size()); + vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); initCCDB(bc); @@ -1113,6 +1137,10 @@ struct matchingMFT { if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_dr.begin(), vec_min_dr.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_dr.end()) { continue; } + } else if (cfgBestMatchFinder == 2) { // 2d + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_2d.begin(), vec_min_2d.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_2d.end()) { + continue; + } } fillHistograms(collision, fwdtrack, fwdtracks, mfttracks, nullptr); @@ -1123,6 +1151,8 @@ struct matchingMFT { vec_min_chi2MatchMCHMFT.shrink_to_fit(); vec_min_dr.clear(); vec_min_dr.shrink_to_fit(); + vec_min_2d.clear(); + vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(matchingMFT, processWithFTTCA, "process with FTTCA", true); @@ -1136,6 +1166,7 @@ struct matchingMFT { vec_min_chi2MatchMCHMFT.reserve(fwdtracks.size()); vec_min_dr.reserve(fwdtracks.size()); + vec_min_2d.reserve(fwdtracks.size()); for (const auto& collision : collisions) { auto bc = collision.template bc_as(); @@ -1182,6 +1213,10 @@ struct matchingMFT { if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_dr.begin(), vec_min_dr.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_dr.end()) { continue; } + } else if (cfgBestMatchFinder == 2) { // 2d + if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack && std::find(vec_min_2d.begin(), vec_min_2d.end(), std::make_tuple(fwdtrack.globalIndex(), fwdtrack.matchMCHTrackId(), fwdtrack.matchMFTTrackId())) == vec_min_2d.end()) { + continue; + } } fillHistograms(collision, fwdtrack, fwdtracks, mfttracks, mftCovs); } // end of fwdtrack loop @@ -1192,6 +1227,8 @@ struct matchingMFT { vec_min_chi2MatchMCHMFT.shrink_to_fit(); vec_min_dr.clear(); vec_min_dr.shrink_to_fit(); + vec_min_2d.clear(); + vec_min_2d.shrink_to_fit(); } PROCESS_SWITCH(matchingMFT, processWithFTTCA_withMFTCov, "process with FTTCA with MFTCov", false); }; From 8ba1b29e643912a892f8c7de0de4d977797cb35f Mon Sep 17 00:00:00 2001 From: Daiki Sekihata Date: Mon, 2 Feb 2026 19:50:02 +0100 Subject: [PATCH 2/3] Clean up unused variables in skimmerPrimaryMuon.cxx Removed unused variables min_dr and min_distance_2d from skimmerPrimaryMuon.cxx. --- PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx index e353da03a88..64dad7b56bf 100644 --- a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx +++ b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx @@ -565,7 +565,8 @@ struct skimmerPrimaryMuon { float phiMatchedMCHMID = propmuonAtPV_Matched.getPhi(); o2::math_utils::bringTo02Pi(phiMatchedMCHMID); - float min_chi2MatchMCHMFT = 1e+10, min_dr = 1e+10, min_distance_2d = 1e+10; + // float min_chi2MatchMCHMFT = 1e+10, min_dr = 1e+10, min_distance_2d = 1e+10; + float min_chi2MatchMCHMFT = 1e+10; std::tuple tupleIds_at_min_chi2mftmch; std::tuple tupleIds_at_min_dr; std::tuple tupleIds_at_min_distance_2d; From 83f5025d2e661da79930fb8c113350b0c7a66e0d Mon Sep 17 00:00:00 2001 From: Daiki Sekihata Date: Mon, 2 Feb 2026 20:23:25 +0100 Subject: [PATCH 3/3] Remove unused tuple declarations in skimmerPrimaryMuon Comment out unused tuple declarations for minimum distance calculations. --- PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx index 64dad7b56bf..b3e70174527 100644 --- a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx +++ b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx @@ -568,8 +568,8 @@ struct skimmerPrimaryMuon { // float min_chi2MatchMCHMFT = 1e+10, min_dr = 1e+10, min_distance_2d = 1e+10; float min_chi2MatchMCHMFT = 1e+10; std::tuple tupleIds_at_min_chi2mftmch; - std::tuple tupleIds_at_min_dr; - std::tuple tupleIds_at_min_distance_2d; + // std::tuple tupleIds_at_min_dr; + // std::tuple tupleIds_at_min_distance_2d; for (const auto& muon_tmp : muons_per_MCHMID) { if (muon_tmp.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { auto tupleId = std::make_tuple(muon_tmp.globalIndex(), muon_tmp.matchMCHTrackId(), muon_tmp.matchMFTTrackId());