From 9f5bf631e0aa330c8cde79f99b1cd809f386f2e0 Mon Sep 17 00:00:00 2001 From: Daiki Sekihata Date: Tue, 21 Jul 2026 16:31:14 +0200 Subject: [PATCH 1/2] [PWGEM/Dilepton] update muon cov --- Common/Core/fwdtrackUtilities.h | 62 +++++++++ PWGEM/Dilepton/Core/Dilepton.h | 8 +- PWGEM/Dilepton/Core/DileptonHadronMPC.h | 8 +- PWGEM/Dilepton/Core/DileptonSV.h | 16 +-- PWGEM/Dilepton/Core/SingleTrackQC.h | 12 +- PWGEM/Dilepton/Core/SingleTrackQCMC.h | 12 +- PWGEM/Dilepton/DataModel/dileptonTables.h | 38 ++++- .../TableProducer/Converters/CMakeLists.txt | 5 + .../Converters/muonConverter4.cxx | 67 +++++++++ .../TableProducer/skimmerPrimaryMuon.cxx | 130 ++++++++++++------ PWGEM/Dilepton/Utils/EMFwdTrack.h | 6 +- PWGEM/Dilepton/Utils/EMTrackUtilities.h | 6 +- 12 files changed, 292 insertions(+), 78 deletions(-) create mode 100644 PWGEM/Dilepton/TableProducer/Converters/muonConverter4.cxx diff --git a/Common/Core/fwdtrackUtilities.h b/Common/Core/fwdtrackUtilities.h index c5576f3d638..7265f1fe581 100644 --- a/Common/Core/fwdtrackUtilities.h +++ b/Common/Core/fwdtrackUtilities.h @@ -270,6 +270,68 @@ o2::dataformats::GlobalFwdTrack refitGlobalMuonCov(TFwdTrack const& muon, TMFTTr return globalTrack; } +template +float getFwdChi2IP(TTrackParCovFwd const& inputTrk, TCollision const& collision, const float bz) +{ + // this function returns imcompatibility of fwdtrack respect to a given PV. + // fwdtracks are never PV contributors in ALICE. + // chi2IP is defined as chi2^{PV}_{with fwdtrack} - chi2^{PV}_{without fwdtrack}. https://arxiv.org/abs/2604.11574 + // chi2IP cannot be used to decide the best fwdtrack-to-collision match or MFT-MCH match, because it gives biases toward small muon impact parameter. + // chi2IP should be used only after the best fwdtrack-to-collision association and the best MFT-MCH match are defined. + + auto trk = inputTrk; // mutable copy + + if (std::abs(bz) < 1e-12) { + trk.propagateToZlinear(collision.posZ()); + } else { + trk.propagateToZhelix(collision.posZ(), bz); + } + + float x = trk.getX(); + float y = trk.getY(); + float phi = trk.getPhi(); + float tgl = trk.getTanl(); + + if (!std::isfinite(x) || !std::isfinite(y) || !std::isfinite(phi) || !std::isfinite(tgl) || std::abs(tgl) < 1e-8) { + return -1.f; + } + + float dx = x - collision.posX(); + float dy = y - collision.posY(); + + // dx/dz and dy/dz for TrackParFwd. + float tx = std::cos(phi) / tgl; + float ty = std::sin(phi) / tgl; + + const auto& ct = trk.getCovariances(); + + float trkXX = ct(0, 0); + float trkXY = ct(0, 1); + float trkYY = ct(1, 1); + + float pvXX = collision.covXX() - 2.0 * tx * collision.covXZ() + tx * tx * collision.covZZ(); + float pvXY = collision.covXY() - ty * collision.covXZ() - tx * collision.covYZ() + tx * ty * collision.covZZ(); + float pvYY = collision.covYY() - 2.0 * ty * collision.covYZ() + ty * ty * collision.covZZ(); + + float sXX = trkXX + pvXX; + float sXY = trkXY + pvXY; + float sYY = trkYY + pvYY; + + float det = sXX * sYY - sXY * sXY; + + if (!std::isfinite(det) || det <= 0.0) { + return -1.f; + } + + float chi2 = (sYY * dx * dx - 2.0 * sXY * dx * dy + sXX * dy * dy) / det; + + if (!std::isfinite(chi2) || chi2 < -1e-8) { + return -1.f; + } + + return std::max(0.f, chi2); +} + template float getFwdChi2IP(const TFullFwdTrack& fwdtrack, const TCollision& collision, const float bz, const float zShift) { diff --git a/PWGEM/Dilepton/Core/Dilepton.h b/PWGEM/Dilepton/Core/Dilepton.h index 1d3e2a75864..148ed2d99c6 100644 --- a/PWGEM/Dilepton/Core/Dilepton.h +++ b/PWGEM/Dilepton/Core/Dilepton.h @@ -1221,9 +1221,9 @@ struct Dilepton { used_trackIds_per_col.emplace_back(t1.globalIndex()); if (cfgDoMix) { if (t1.sign() > 0) { - emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY())); + emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA())); } else { - emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY())); + emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA())); } } } @@ -1231,9 +1231,9 @@ struct Dilepton { used_trackIds_per_col.emplace_back(t2.globalIndex()); if (cfgDoMix) { if (t2.sign() > 0) { - emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY())); + emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA())); } else { - emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY())); + emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA())); } } } diff --git a/PWGEM/Dilepton/Core/DileptonHadronMPC.h b/PWGEM/Dilepton/Core/DileptonHadronMPC.h index 5ce486d89c0..ea1740e3ff6 100644 --- a/PWGEM/Dilepton/Core/DileptonHadronMPC.h +++ b/PWGEM/Dilepton/Core/DileptonHadronMPC.h @@ -827,9 +827,9 @@ struct DileptonHadronMPC { used_trackIds_per_col.emplace_back(t1.globalIndex()); if (cfgDoMix) { if (t1.sign() > 0) { - emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY())); + emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA())); } else { - emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY())); + emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA())); } } } @@ -837,9 +837,9 @@ struct DileptonHadronMPC { used_trackIds_per_col.emplace_back(t2.globalIndex()); if (cfgDoMix) { if (t2.sign() > 0) { - emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY())); + emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA())); } else { - emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY())); + emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA())); } } } diff --git a/PWGEM/Dilepton/Core/DileptonSV.h b/PWGEM/Dilepton/Core/DileptonSV.h index 1df85e1d81f..4acaf0f6028 100644 --- a/PWGEM/Dilepton/Core/DileptonSV.h +++ b/PWGEM/Dilepton/Core/DileptonSV.h @@ -1433,16 +1433,16 @@ struct DileptonSV { used_trackIds_per_col.emplace_back(t1.globalIndex()); if (mixingGroup.cfgDoMix) { if (t1.sign() > 0) { - emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt1, candidate.eta1, candidate.phi1, leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY())); - // emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY(), + emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt1, candidate.eta1, candidate.phi1, leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA())); + // emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA(), // t1.x(), t1.y(), t1.z(), t1.tgl(), // t1.cPhiX(), t1.cPhiY(), t1.cPhiPhi(), // t1.cTglX(), t1.cTglY(), t1.cTglPhi(), t1.cTglTgl(), // t1.c1PtX(), t1.c1PtY(), t1.c1PtPhi(), t1.c1PtTgl(), t1.c1Pt21Pt2(), t1.chi2())); } else { - emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt1, candidate.eta1, candidate.phi1, leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY())); - // emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXX(), t1.cXY(), t1.cYY(), + emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt1, candidate.eta1, candidate.phi1, leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA())); + // emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t1.pt(), t1.eta(), t1.phi(), leptonM1, t1.sign(), t1.fwdDcaX(), t1.fwdDcaY(), t1.cXXatDCA(), t1.cXYatDCA(), t1.cYYatDCA(), // t1.x(), t1.y(), t1.z(), t1.tgl(), // t1.cPhiX(), t1.cPhiY(), t1.cPhiPhi(), // t1.cTglX(), t1.cTglY(), t1.cTglPhi(), t1.cTglTgl(), @@ -1454,15 +1454,15 @@ struct DileptonSV { used_trackIds_per_col.emplace_back(t2.globalIndex()); if (mixingGroup.cfgDoMix) { if (t2.sign() > 0) { - emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt2, candidate.eta2, candidate.phi2, leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY())); - // emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY(), + emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt2, candidate.eta2, candidate.phi2, leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA())); + // emh_pos->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA(), // t2.x(), t2.y(), t2.z(), t2.tgl(), // t2.cPhiX(), t2.cPhiY(), t2.cPhiPhi(), // t2.cTglX(), t2.cTglY(), t2.cTglPhi(), t2.cTglTgl(), // t2.c1PtX(), t2.c1PtY(), t2.c1PtPhi(), t2.c1PtTgl(), t2.c1Pt21Pt2(), t2.chi2())); } else { - emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt2, candidate.eta2, candidate.phi2, leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY())); - // emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXX(), t2.cXY(), t2.cYY(), + emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrack(candidate.pt2, candidate.eta2, candidate.phi2, leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA())); + // emh_neg->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMFwdTrackWithCov(t2.pt(), t2.eta(), t2.phi(), leptonM2, t2.sign(), t2.fwdDcaX(), t2.fwdDcaY(), t2.cXXatDCA(), t2.cXYatDCA(), t2.cYYatDCA(), // t2.x(), t2.y(), t2.z(), t2.tgl(), // t2.cPhiX(), t2.cPhiY(), t2.cPhiPhi(), // t2.cTglX(), t2.cTglY(), t2.cTglPhi(), t2.cTglTgl(), diff --git a/PWGEM/Dilepton/Core/SingleTrackQC.h b/PWGEM/Dilepton/Core/SingleTrackQC.h index 792c376a5e4..0f0fe143a4b 100644 --- a/PWGEM/Dilepton/Core/SingleTrackQC.h +++ b/PWGEM/Dilepton/Core/SingleTrackQC.h @@ -664,9 +664,9 @@ struct SingleTrackQC { fRegistry.fill(HIST("Track/positive/hQoverPt"), track.sign() / track.pt()); fRegistry.fill(HIST("Track/positive/hTrackType"), track.trackType()); fRegistry.fill(HIST("Track/positive/hDCAxy"), track.fwdDcaX(), track.fwdDcaY()); - fRegistry.fill(HIST("Track/positive/hDCAxySigma"), track.fwdDcaX() / std::sqrt(track.cXX()), track.fwdDcaY() / std::sqrt(track.cYY())); - fRegistry.fill(HIST("Track/positive/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXX()) * 1e+4); - fRegistry.fill(HIST("Track/positive/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYY()) * 1e+4); + fRegistry.fill(HIST("Track/positive/hDCAxySigma"), track.fwdDcaX() / std::sqrt(track.cXXatDCA()), track.fwdDcaY() / std::sqrt(track.cYYatDCA())); + fRegistry.fill(HIST("Track/positive/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXXatDCA()) * 1e+4); + fRegistry.fill(HIST("Track/positive/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYYatDCA()) * 1e+4); fRegistry.fill(HIST("Track/positive/hDCAxyRes_Pt"), track.pt(), o2::aod::pwgem::dilepton::utils::emtrackutil::sigmaFwdDcaXY(track) * 1e+4); fRegistry.fill(HIST("Track/positive/hDCAx_PosZ"), collision.posZ(), track.fwdDcaX()); fRegistry.fill(HIST("Track/positive/hDCAy_PosZ"), collision.posZ(), track.fwdDcaY()); @@ -688,9 +688,9 @@ struct SingleTrackQC { fRegistry.fill(HIST("Track/negative/hQoverPt"), track.sign() / track.pt()); fRegistry.fill(HIST("Track/negative/hTrackType"), track.trackType()); fRegistry.fill(HIST("Track/negative/hDCAxy"), track.fwdDcaX(), track.fwdDcaY()); - fRegistry.fill(HIST("Track/negative/hDCAxySigma"), track.fwdDcaX() / std::sqrt(track.cXX()), track.fwdDcaY() / std::sqrt(track.cYY())); - fRegistry.fill(HIST("Track/negative/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXX()) * 1e+4); - fRegistry.fill(HIST("Track/negative/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYY()) * 1e+4); + fRegistry.fill(HIST("Track/negative/hDCAxySigma"), track.fwdDcaX() / std::sqrt(track.cXXatDCA()), track.fwdDcaY() / std::sqrt(track.cYYatDCA())); + fRegistry.fill(HIST("Track/negative/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXXatDCA()) * 1e+4); + fRegistry.fill(HIST("Track/negative/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYYatDCA()) * 1e+4); fRegistry.fill(HIST("Track/negative/hDCAxyRes_Pt"), track.pt(), o2::aod::pwgem::dilepton::utils::emtrackutil::sigmaFwdDcaXY(track) * 1e+4); fRegistry.fill(HIST("Track/negative/hDCAx_PosZ"), collision.posZ(), track.fwdDcaX()); fRegistry.fill(HIST("Track/negative/hDCAy_PosZ"), collision.posZ(), track.fwdDcaY()); diff --git a/PWGEM/Dilepton/Core/SingleTrackQCMC.h b/PWGEM/Dilepton/Core/SingleTrackQCMC.h index c0b119333b1..b2d1243df6d 100644 --- a/PWGEM/Dilepton/Core/SingleTrackQCMC.h +++ b/PWGEM/Dilepton/Core/SingleTrackQCMC.h @@ -776,9 +776,9 @@ struct SingleTrackQCMC { fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hTrackType"), track.trackType()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxy"), std::sqrt(std::pow(track.fwdDcaX(), 2) + std::pow(track.fwdDcaY(), 2))); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxy2D"), track.fwdDcaX(), track.fwdDcaY()); - fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxy2DinSigma"), track.fwdDcaX() / std::sqrt(track.cXX()), track.fwdDcaY() / std::sqrt(track.cYY())); - fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXX()) * 1e+4); - fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYY()) * 1e+4); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxy2DinSigma"), track.fwdDcaX() / std::sqrt(track.cXXatDCA()), track.fwdDcaY() / std::sqrt(track.cYYatDCA())); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXXatDCA()) * 1e+4); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYYatDCA()) * 1e+4); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAxyRes_Pt"), track.pt(), o2::aod::pwgem::dilepton::utils::emtrackutil::sigmaFwdDcaXY(track) * 1e+4); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAx_PosZ"), collision.posZ(), track.fwdDcaX()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("positive/hDCAy_PosZ"), collision.posZ(), track.fwdDcaY()); @@ -809,9 +809,9 @@ struct SingleTrackQCMC { fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hTrackType"), track.trackType()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxy"), std::sqrt(std::pow(track.fwdDcaX(), 2) + std::pow(track.fwdDcaY(), 2))); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxy2D"), track.fwdDcaX(), track.fwdDcaY()); - fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxy2DinSigma"), track.fwdDcaX() / std::sqrt(track.cXX()), track.fwdDcaY() / std::sqrt(track.cYY())); - fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXX()) * 1e+4); - fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYY()) * 1e+4); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxy2DinSigma"), track.fwdDcaX() / std::sqrt(track.cXXatDCA()), track.fwdDcaY() / std::sqrt(track.cYYatDCA())); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxRes_Pt"), track.pt(), std::sqrt(track.cXXatDCA()) * 1e+4); + fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAyRes_Pt"), track.pt(), std::sqrt(track.cYYatDCA()) * 1e+4); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAxyRes_Pt"), track.pt(), o2::aod::pwgem::dilepton::utils::emtrackutil::sigmaFwdDcaXY(track) * 1e+4); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAx_PosZ"), collision.posZ(), track.fwdDcaX()); fRegistry.fill(HIST("Track/") + HIST(lepton_source_types[lepton_source_id]) + HIST("negative/hDCAy_PosZ"), collision.posZ(), track.fwdDcaY()); diff --git a/PWGEM/Dilepton/DataModel/dileptonTables.h b/PWGEM/Dilepton/DataModel/dileptonTables.h index fe484802594..6074a1e66f8 100644 --- a/PWGEM/Dilepton/DataModel/dileptonTables.h +++ b/PWGEM/Dilepton/DataModel/dileptonTables.h @@ -1131,7 +1131,31 @@ DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuons_003, "AOD", "EMPRIMARYMU", 3, //! emprimarymuon::Py, emprimarymuon::Pz); -using EMPrimaryMuons = EMPrimaryMuons_003; +DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuons_004, "AOD", "EMPRIMARYMU", 4, //! + o2::soa::Index<>, emprimarymuon::CollisionId, + emprimarymuon::FwdTrackId, emprimarymuon::MFTTrackId, emprimarymuon::MCHTrackId, fwdtrack::TrackType, + fwdtrack::Pt, fwdtrack::Eta, fwdtrack::Phi, emprimarymuon::Sign, + fwdtrack::FwdDcaX, fwdtrack::FwdDcaY, emprimarymuon::CXXatDCA, emprimarymuon::CYYatDCA, emprimarymuon::CXYatDCA, + emprimarymuon::PtMatchedMCHMID, emprimarymuon::EtaMatchedMCHMID, emprimarymuon::PhiMatchedMCHMID, + + fwdtrack::NClusters, fwdtrack::PDca, fwdtrack::RAtAbsorberEnd, + fwdtrack::Chi2, fwdtrack::Chi2MatchMCHMID, fwdtrack::Chi2MatchMCHMFT, emprimarymuon::DiffChi2MatchMCHMFT, + fwdtrack::MCHBitMap, fwdtrack::MIDBitMap, fwdtrack::MIDBoards, + fwdtrack::MFTClusterSizesAndTrackFlags, emprimarymuon::Chi2MFT, emprimarymuon::IsAssociatedToMPC, emprimarymuon::IsAmbiguous, + + // dynamic column + emprimarymuon::Signed1Pt, + emprimarymuon::Tgl, + emprimarymuon::NClustersMFT, + fwdtrack::IsCA, + emprimarymuon::MFTClusterMap, + emprimarymuon::P, + emprimarymuon::Px, + emprimarymuon::Py, + emprimarymuon::Pz); + + +using EMPrimaryMuons = EMPrimaryMuons_004; // iterators using EMPrimaryMuon = EMPrimaryMuons::iterator; @@ -1152,7 +1176,17 @@ DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuonsCov_002, "AOD", "EMPRIMARYMUCOV", 2, / o2::aod::fwdtrack::C1PtPhi, o2::aod::fwdtrack::C1PtTgl, o2::aod::fwdtrack::C1Pt21Pt2); -using EMPrimaryMuonsCov = EMPrimaryMuonsCov_002; + + +DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuonsCov_003, "AOD", "EMPRIMARYMUCOV", 3, //! at PV. Signed1Pt, Tgl and Phi are in EMPrimaryMuons table. + fwdtrack::X, fwdtrack::Y, fwdtrack::Z, + o2::aod::fwdtrack::CXX, + o2::aod::fwdtrack::CXY, o2::aod::fwdtrack::CYY, + o2::aod::fwdtrack::CPhiX, o2::aod::fwdtrack::CPhiY, o2::aod::fwdtrack::CPhiPhi, + o2::aod::fwdtrack::CTglX, o2::aod::fwdtrack::CTglY, o2::aod::fwdtrack::CTglPhi, o2::aod::fwdtrack::CTglTgl, + o2::aod::fwdtrack::C1PtX, o2::aod::fwdtrack::C1PtY, o2::aod::fwdtrack::C1PtPhi, o2::aod::fwdtrack::C1PtTgl, o2::aod::fwdtrack::C1Pt21Pt2); + +using EMPrimaryMuonsCov = EMPrimaryMuonsCov_003; // iterators using EMPrimaryMuonCov = EMPrimaryMuonsCov::iterator; diff --git a/PWGEM/Dilepton/TableProducer/Converters/CMakeLists.txt b/PWGEM/Dilepton/TableProducer/Converters/CMakeLists.txt index 5ac1908c3ff..3185ecb31af 100644 --- a/PWGEM/Dilepton/TableProducer/Converters/CMakeLists.txt +++ b/PWGEM/Dilepton/TableProducer/Converters/CMakeLists.txt @@ -75,6 +75,11 @@ o2physics_add_dpl_workflow(muon-converter3 PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore COMPONENT_NAME Analysis) +o2physics_add_dpl_workflow(muon-converter4 + SOURCES muonConverter4.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore + COMPONENT_NAME Analysis) + o2physics_add_dpl_workflow(muon-selfid-converter1 SOURCES muonSelfIdConverter1.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore diff --git a/PWGEM/Dilepton/TableProducer/Converters/muonConverter4.cxx b/PWGEM/Dilepton/TableProducer/Converters/muonConverter4.cxx new file mode 100644 index 00000000000..98035efe8e9 --- /dev/null +++ b/PWGEM/Dilepton/TableProducer/Converters/muonConverter4.cxx @@ -0,0 +1,67 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +// +// ======================== +// +// This code produces muon table 004 from 003. +// Please write to: daiki.sekihata@cern.ch + +#include "PWGEM/Dilepton/DataModel/dileptonTables.h" + +#include +#include +#include + +using namespace o2; +using namespace o2::aod; +using namespace o2::framework; +using namespace o2::framework::expressions; +using namespace o2::soa; + +struct muonConverter4 { + Produces muon_004; + Produces muoncov_003; + + void processMuon(o2::aod::EMPrimaryMuons_003 const& muons) + { + for (const auto& muon : muons) { + muon_004( + muon.collisionId(), + muon.fwdtrackId(), muon.mfttrackId(), muon.mchtrackId(), muon.trackType(), + muon.pt(), muon.eta(), muon.phi(), muon.sign(), + muon.fwdDcaX(), muon.fwdDcaY(), muon.cXX(), muon.cYY(), muon.cXY(), + muon.ptMatchedMCHMID(), muon.etaMatchedMCHMID(), muon.phiMatchedMCHMID(), + muon.nClusters(), muon.pDca(), muon.rAtAbsorberEnd(), + muon.chi2(), muon.chi2MatchMCHMID(), muon.chi2MatchMCHMFT(), muon.diffChi2MatchingMCHMFT(), + muon.mchBitMap(), muon.midBitMap(), muon.midBoards(), + muon.mftClusterSizesAndTrackFlags(), muon.chi2MFT(), muon.isAssociatedToMPC(), muon.isAmbiguous()); + } // end of muon loop + } + PROCESS_SWITCH(muonConverter4, processMuon, "convert primary muon 003 into 004", true); + + void processMuonCov(soa::Join const& muons) + { + for (const auto& muon : muons) { + muoncov_003(muon.x(), muon.y(), muon.y(), + muon.cXX(), + muon.cYY(), muon.cXY(), + muon.cPhiX(), muon.cPhiY(), muon.cPhiPhi(), + muon.cTglX(), muon.cTglY(), muon.cTglPhi(), muon.cTglTgl(), + muon.c1PtX(), muon.c1PtY(), muon.c1PtPhi(), muon.c1PtTgl(), muon.c1Pt21Pt2()); + } // end of muon loop + } + PROCESS_SWITCH(muonConverter4, processMuonCov, "convert primary muon 002 into 003", false); +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{adaptAnalysisTask(cfgc, TaskName{"muon-converter4"})}; +} diff --git a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx index b77586d63d6..210e41414cd 100644 --- a/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx +++ b/PWGEM/Dilepton/TableProducer/skimmerPrimaryMuon.cxx @@ -54,8 +54,6 @@ #include #include -#include - using namespace o2; using namespace o2::soa; using namespace o2::framework; @@ -205,9 +203,10 @@ struct skimmerPrimaryMuon { fRegistry.add("MFTMCHMID/hChi2MatchMCHMFT_Pt", "chi2 match MCH-MFT;p_{T} (GeV/c);matching chi2 MCH-MFT", kTH2F, {{100, 0, 10}, {200, 0.0f, 100}}, false); fRegistry.add("MFTMCHMID/hDCAxy2D", "DCA x vs. y;DCA_{x} (cm);DCA_{y} (cm)", kTH2F, {{200, -1, 1}, {200, -1, +1}}, false); fRegistry.add("MFTMCHMID/hDCAxy2DinSigma", "DCA x vs. y in sigma;DCA_{x} (#sigma);DCA_{y} (#sigma)", kTH2F, {{200, -10, 10}, {200, -10, +10}}, false); - fRegistry.add("MFTMCHMID/hDCAxy", "DCAxy;DCA_{xy} (cm);", kTH1F, {{100, 0, 1}}, false); - fRegistry.add("MFTMCHMID/hDCAxyz", "DCA xy vs. z;DCA_{xy} (cm);DCA_{z} (cm)", kTH2F, {{100, 0, 1}, {200, -0.1, 0.1}}, false); + fRegistry.add("MFTMCHMID/hDCAxy", "DCAxy;DCA_{xy} (cm);", kTH1F, {{1000, 0, 1}}, false); fRegistry.add("MFTMCHMID/hDCAxyinSigma", "DCAxy in sigma;DCA_{xy} (#sigma);", kTH1F, {{100, 0, 10}}, false); + fRegistry.add("MFTMCHMID/hLog10Chi2IP", "chi2IP;log_{10}(#chi^{2}_{IP})", kTH1F, {{100, -5, 5}}, false); + fRegistry.add("MFTMCHMID/hSqrtChi2IP", "chi2IP;#sqrt{#chi^{2}_{IP}}", kTH1F, {{100, 0, 10}}, false); fRegistry.add("MFTMCHMID/hDCAx_PosZ", "DCAx vs. posZ;Z_{vtx} (cm);DCA_{x} (cm)", kTH2F, {{200, -10, +10}, {400, -0.2, +0.2}}, false); fRegistry.add("MFTMCHMID/hDCAy_PosZ", "DCAy vs. posZ;Z_{vtx} (cm);DCA_{y} (cm)", kTH2F, {{200, -10, +10}, {400, -0.2, +0.2}}, false); fRegistry.add("MFTMCHMID/hDCAx_Phi", "DCAx vs. #varphi;#varphi (rad.);DCA_{x} (cm)", kTH2F, {{180, -M_PI, M_PI}, {400, -0.2, +0.2}}, false); @@ -284,21 +283,22 @@ struct skimmerPrimaryMuon { float phi = propmuonAtPV.getPhi(); o2::math_utils::bringTo02Pi(phi); - float dcaX = propmuonAtPV.getX() - collision.posX(); - float dcaY = propmuonAtPV.getY() - collision.posY(); - float dcaZ = propmuonAtPV.getZ() - collision.posZ(); + o2::dataformats::GlobalFwdTrack propmuonAtDCA = propagateMuon(fwdtrack, fwdtrack, collision, propagationPoint::kToDCA, matchingZ, mBz, mZShift); + + float dcaX = propmuonAtDCA.getX() - collision.posX(); + float dcaY = propmuonAtDCA.getY() - collision.posY(); float dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY); float rAtAbsorberEnd = fwdtrack.rAtAbsorberEnd(); // this works only for GlobalMuonTrack - float cXX = propmuonAtPV.getSigma2X(); - float cYY = propmuonAtPV.getSigma2Y(); - float cXY = propmuonAtPV.getSigmaXY(); + float cXXatDCA = propmuonAtDCA.getSigma2X(); + float cYYatDCA = propmuonAtDCA.getSigma2Y(); + float cXYatDCA = propmuonAtDCA.getSigmaXY(); - float det = cXX * cYY - cXY * cXY; // determinanat + float det = cXXatDCA * cYYatDCA - cXYatDCA * cXYatDCA; // determinanat float dcaXYinSigma = 999.f; if (det < 0) { dcaXYinSigma = 999.f; } else { - dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYY + dcaY * dcaY * cXX - 2.f * dcaX * dcaY * cXY) / det / 2.f)); // dca xy in sigma + dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYYatDCA + dcaY * dcaY * cXXatDCA - 2.f * dcaX * dcaY * cXYatDCA) / det / 2.f)); // dca xy in sigma } float sigma_dcaXY = dcaXY / dcaXYinSigma; @@ -326,6 +326,8 @@ struct skimmerPrimaryMuon { float dphi = 999.f; float xMFT = 0.f, yMFT = 0.f; + float chi2IP = -1; + float sfPt = 1.f; // scaling factor for pT of global muons if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::GlobalMuonTrack) { if (fwdtrack.chi2MatchMCHMFT() < 0.f || maxMatchingChi2MCHMFT < fwdtrack.chi2MatchMCHMFT()) { @@ -404,33 +406,79 @@ struct skimmerPrimaryMuon { phi = globalMuon.getPhi(); o2::math_utils::bringTo02Pi(phi); - cXX = globalMuon.getSigma2X(); - cYY = globalMuon.getSigma2Y(); - cXY = globalMuon.getSigmaXY(); + cXXatDCA = globalMuon.getSigma2X(); + cYYatDCA = globalMuon.getSigma2Y(); + cXYatDCA = globalMuon.getSigmaXY(); dcaX = globalMuon.getX() - collision.posX(); dcaY = globalMuon.getY() - collision.posY(); - dcaZ = globalMuon.getZ() - collision.posZ(); dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY); - det = cXX * cYY - cXY * cXY; // determinanat + det = cXXatDCA * cYYatDCA - cXYatDCA * cXYatDCA; // determinanat dcaXYinSigma = 999.f; if (det < 0) { dcaXYinSigma = 999.f; } else { - dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYY + dcaY * dcaY * cXX - 2.f * dcaX * dcaY * cXY) / det / 2.f)); // dca xy in sigma + dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYYatDCA + dcaY * dcaY * cXXatDCA - 2.f * dcaX * dcaY * cXYatDCA) / det / 2.f)); // dca xy in sigma } sigma_dcaXY = dcaXY / dcaXYinSigma; - } + } // end of withMFTCov deta = etaMatchedMCHMID - eta; dphi = phiMatchedMCHMID - phi; o2::math_utils::bringToPMPi(dphi); + chi2IP = getFwdChi2IP(fwdtrack, collision, mBz, mZShift); + if (std::sqrt(std::pow(deta / maxDEta, 2) + std::pow(dphi / maxDPhi, 2)) > 1.f) { return false; } - if (refitGlobalMuon) { - pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); + if (refitGlobalMuon) { // if pT of global muon is modified, one has to redo propagations with updated covariance matrix. + sfPt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))) / propmuonAtPV.getPt(); + pt = propmuonAtPV_Matched.getP() * std::sin(2.f * std::atan(std::exp(-eta))); // sfPt * propmuonAtPV.getPt(); + + const auto& fwdcov = propmuonAtPV.getCovariances(); // covatiance matrix at PV + + auto globalMuonManual = getTrackParCovFwdShiftManual( + propmuonAtPV.getX(), propmuonAtPV.getY(), propmuonAtPV.getPhi(), propmuonAtPV.getTgl(), propmuonAtPV.getInvQPt() / sfPt, + fwdcov(0, 0), + fwdcov(1, 0), fwdcov(1, 1), + fwdcov(2, 0), fwdcov(2, 1), fwdcov(2, 2), + fwdcov(3, 0), fwdcov(3, 1), fwdcov(3, 2), fwdcov(3, 3), + fwdcov(4, 0) / sfPt, fwdcov(4, 1) / sfPt, fwdcov(4, 2) / sfPt, fwdcov(4, 3) / sfPt, fwdcov(4, 4) / sfPt / sfPt, + propmuonAtPV.getZ(), 0.0, fwdtrack.chi2()); + + auto globalMuonManualAtZPV = o2::aod::fwdtrackutils::propagateTrackParCovFwd(globalMuonManual, fwdtrack.trackType(), collision, propagationPoint::kToDCA, matchingZ, mBz); + + dcaX = globalMuonManualAtZPV.getX() - collision.posX(); + dcaY = globalMuonManualAtZPV.getY() - collision.posY(); + cXXatDCA = globalMuonManualAtZPV.getSigma2X(); + cYYatDCA = globalMuonManualAtZPV.getSigma2Y(); + cXYatDCA = globalMuonManualAtZPV.getSigmaXY(); + // add convariance matrix of vertex + + float tx = std::cos(globalMuonManualAtZPV.getPhi()) / globalMuonManualAtZPV.getTgl(); // dx/dz + float ty = std::sin(globalMuonManualAtZPV.getPhi()) / globalMuonManualAtZPV.getTgl(); // dy/dz + + // PV covariance projected into the (dx, dy) residual space. + float resPVXX = collision.covXX() - 2.0 * tx * collision.covXZ() + tx * tx * collision.covZZ(); + float resPVXY = collision.covXY() - ty * collision.covXZ() - tx * collision.covYZ() + tx * ty * collision.covZZ(); + float resPVYY = collision.covYY() - 2.0 * ty * collision.covYZ() + ty * ty * collision.covZZ(); + + // Full residual covariance. + cXXatDCA = cXXatDCA + resPVXX; + cXYatDCA = cXYatDCA + resPVXY; + cYYatDCA = cYYatDCA + resPVYY; + + det = cXXatDCA * cYYatDCA - cXYatDCA * cXYatDCA; // determinanat + dcaXYinSigma = 999.f; + if (det < 0) { + dcaXYinSigma = 999.f; + } else { + dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYYatDCA + dcaY * dcaY * cXXatDCA - 2.f * dcaX * dcaY * cXYatDCA) / det / 2.f)); // dca xy in sigma + } + sigma_dcaXY = dcaXY / dcaXYinSigma; + + chi2IP = getFwdChi2IP(globalMuonManual, collision, mBz); } } else if (fwdtrack.trackType() == o2::aod::fwdtrack::ForwardTrackTypeEnum::MuonStandaloneTrack) { o2::dataformats::GlobalFwdTrack propmuonAtRabs = propagateMuon(fwdtrack, fwdtrack, collision, propagationPoint::kToRabs, matchingZ, mBz, mZShift); // this is necessary only for MuonStandaloneTrack @@ -439,21 +487,20 @@ struct skimmerPrimaryMuon { rAtAbsorberEnd = std::sqrt(xAbs * xAbs + yAbs * yAbs); // Redo propagation only for muon tracks // propagation of MFT tracks alredy done in reconstruction o2::dataformats::GlobalFwdTrack propmuonAtDCA = propagateMuon(fwdtrack, fwdtrack, collision, propagationPoint::kToDCA, matchingZ, mBz, mZShift); - cXX = propmuonAtDCA.getSigma2X(); - cYY = propmuonAtDCA.getSigma2Y(); - cXY = propmuonAtDCA.getSigmaXY(); + cXXatDCA = propmuonAtDCA.getSigma2X(); + cYYatDCA = propmuonAtDCA.getSigma2Y(); + cXYatDCA = propmuonAtDCA.getSigmaXY(); dcaX = propmuonAtDCA.getX() - collision.posX(); dcaY = propmuonAtDCA.getY() - collision.posY(); - dcaZ = propmuonAtDCA.getZ() - collision.posZ(); dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY); pDCA = fwdtrack.p() * dcaXY; - det = cXX * cYY - cXY * cXY; // determinanat + det = cXXatDCA * cYYatDCA - cXYatDCA * cXYatDCA; // determinanat dcaXYinSigma = 999.f; if (det < 0) { dcaXYinSigma = 999.f; } else { - dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYY + dcaY * dcaY * cXX - 2.f * dcaX * dcaY * cXY) / det / 2.f)); // dca xy in sigma + dcaXYinSigma = std::sqrt(std::fabs((dcaX * dcaX * cYYatDCA + dcaY * dcaY * cXXatDCA - 2.f * dcaX * dcaY * cXYatDCA) / det / 2.f)); // dca xy in sigma } sigma_dcaXY = dcaXY / dcaXYinSigma; } else { @@ -475,19 +522,18 @@ struct skimmerPrimaryMuon { // LOGF(info, "isAmbiguous = %d, isAssociatedToMPC = %d, fwdtrack.globalIndex() = %d, fwdtrack.collisionId() = %d, collision.globalIndex() = %d", isAmbiguous, isAssociatedToMPC, fwdtrack.globalIndex(), fwdtrack.collisionId(), collision.globalIndex()); emprimarymuons(collision.globalIndex(), fwdtrack.globalIndex(), fwdtrack.matchMFTTrackId(), fwdtrack.matchMCHTrackId(), fwdtrack.trackType(), - pt, eta, phi, fwdtrack.sign(), dcaX, dcaY, cXX, cYY, cXY, ptMatchedMCHMID, etaMatchedMCHMID, phiMatchedMCHMID, - // etaMatchedMCHMIDatMP, phiMatchedMCHMIDatMP, etaMatchedMFTatMP, phiMatchedMFTatMP, + pt, eta, phi, fwdtrack.sign(), dcaX, dcaY, cXXatDCA, cYYatDCA, cXYatDCA, ptMatchedMCHMID, etaMatchedMCHMID, phiMatchedMCHMID, fwdtrack.nClusters(), pDCA, rAtAbsorberEnd, fwdtrack.chi2(), fwdtrack.chi2MatchMCHMID(), fwdtrack.chi2MatchMCHMFT(), map_diff_chi2MatchMCHMFT[fwdtrack.globalIndex()], fwdtrack.mchBitMap(), fwdtrack.midBitMap(), fwdtrack.midBoards(), mftClusterSizesAndTrackFlags, chi2mft, isAssociatedToMPC, isAmbiguous); const auto& fwdcov = propmuonAtPV.getCovariances(); // covatiance matrix at PV emprimarymuonscov( propmuonAtPV.getX(), propmuonAtPV.getY(), propmuonAtPV.getZ(), - // fwdcov(0, 0), - // fwdcov(0, 1), fwdcov(1, 1), + fwdcov(0, 0), + fwdcov(1, 0), fwdcov(1, 1), fwdcov(2, 0), fwdcov(2, 1), fwdcov(2, 2), fwdcov(3, 0), fwdcov(3, 1), fwdcov(3, 2), fwdcov(3, 3), - fwdcov(4, 0), fwdcov(4, 1), fwdcov(4, 2), fwdcov(4, 3), fwdcov(4, 4)); + fwdcov(4, 0) / sfPt, fwdcov(4, 1) / sfPt, fwdcov(4, 2) / sfPt, fwdcov(4, 3) / sfPt, fwdcov(4, 4) / sfPt / sfPt); // See definition DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackFwd.h // Covariance matrix of track parameters, ordered as follows: @@ -516,13 +562,14 @@ struct skimmerPrimaryMuon { fRegistry.fill(HIST("MFTMCHMID/hChi2MatchMCHMID_Pt"), pt, fwdtrack.chi2MatchMCHMID()); fRegistry.fill(HIST("MFTMCHMID/hChi2MatchMCHMFT_Pt"), pt, fwdtrack.chi2MatchMCHMFT()); fRegistry.fill(HIST("MFTMCHMID/hDCAxy2D"), dcaX, dcaY); - fRegistry.fill(HIST("MFTMCHMID/hDCAxy2DinSigma"), dcaX / std::sqrt(cXX), dcaY / std::sqrt(cYY)); + fRegistry.fill(HIST("MFTMCHMID/hDCAxy2DinSigma"), dcaX / std::sqrt(cXXatDCA), dcaY / std::sqrt(cYYatDCA)); fRegistry.fill(HIST("MFTMCHMID/hDCAxy"), dcaXY); - fRegistry.fill(HIST("MFTMCHMID/hDCAxyz"), dcaXY, dcaZ); fRegistry.fill(HIST("MFTMCHMID/hDCAxyinSigma"), dcaXYinSigma); - fRegistry.fill(HIST("MFTMCHMID/hDCAxResolutionvsPt"), pt, std::sqrt(cXX) * 1e+4); // convert cm to um - fRegistry.fill(HIST("MFTMCHMID/hDCAyResolutionvsPt"), pt, std::sqrt(cYY) * 1e+4); // convert cm to um - fRegistry.fill(HIST("MFTMCHMID/hDCAxyResolutionvsPt"), pt, sigma_dcaXY * 1e+4); // convert cm to um + fRegistry.fill(HIST("MFTMCHMID/hDCAxResolutionvsPt"), pt, std::sqrt(cXXatDCA) * 1e+4); // convert cm to um + fRegistry.fill(HIST("MFTMCHMID/hDCAyResolutionvsPt"), pt, std::sqrt(cYYatDCA) * 1e+4); // convert cm to um + fRegistry.fill(HIST("MFTMCHMID/hDCAxyResolutionvsPt"), pt, sigma_dcaXY * 1e+4); // convert cm to um + fRegistry.fill(HIST("MFTMCHMID/hLog10Chi2IP"), std::log10(chi2IP)); + fRegistry.fill(HIST("MFTMCHMID/hSqrtChi2IP"), std::sqrt(chi2IP)); fRegistry.fill(HIST("MFTMCHMID/hDCAx_PosZ"), collision.posZ(), dcaX); fRegistry.fill(HIST("MFTMCHMID/hDCAy_PosZ"), collision.posZ(), dcaY); fRegistry.fill(HIST("MFTMCHMID/hDCAx_Phi"), std::atan2(yMFT, xMFT), dcaX); @@ -546,13 +593,12 @@ struct skimmerPrimaryMuon { fRegistry.fill(HIST("MCHMID/hChi2MatchMCHMID_Pt"), pt, fwdtrack.chi2MatchMCHMID()); fRegistry.fill(HIST("MCHMID/hChi2MatchMCHMFT_Pt"), pt, fwdtrack.chi2MatchMCHMFT()); fRegistry.fill(HIST("MCHMID/hDCAxy2D"), dcaX, dcaY); - fRegistry.fill(HIST("MCHMID/hDCAxy2DinSigma"), dcaX / std::sqrt(cXX), dcaY / std::sqrt(cYY)); + fRegistry.fill(HIST("MCHMID/hDCAxy2DinSigma"), dcaX / std::sqrt(cXXatDCA), dcaY / std::sqrt(cYYatDCA)); fRegistry.fill(HIST("MCHMID/hDCAxy"), dcaXY); - fRegistry.fill(HIST("MCHMID/hDCAxyz"), dcaXY, dcaZ); fRegistry.fill(HIST("MCHMID/hDCAxyinSigma"), dcaXYinSigma); - fRegistry.fill(HIST("MCHMID/hDCAxResolutionvsPt"), pt, std::sqrt(cXX) * 1e+4); // convert cm to um - fRegistry.fill(HIST("MCHMID/hDCAyResolutionvsPt"), pt, std::sqrt(cYY) * 1e+4); // convert cm to um - fRegistry.fill(HIST("MCHMID/hDCAxyResolutionvsPt"), pt, sigma_dcaXY * 1e+4); // convert cm to um + fRegistry.fill(HIST("MCHMID/hDCAxResolutionvsPt"), pt, std::sqrt(cXXatDCA) * 1e+4); // convert cm to um + fRegistry.fill(HIST("MCHMID/hDCAyResolutionvsPt"), pt, std::sqrt(cYYatDCA) * 1e+4); // convert cm to um + fRegistry.fill(HIST("MCHMID/hDCAxyResolutionvsPt"), pt, sigma_dcaXY * 1e+4); // convert cm to um } } } diff --git a/PWGEM/Dilepton/Utils/EMFwdTrack.h b/PWGEM/Dilepton/Utils/EMFwdTrack.h index 72b384560ef..a4ad5c43f14 100644 --- a/PWGEM/Dilepton/Utils/EMFwdTrack.h +++ b/PWGEM/Dilepton/Utils/EMFwdTrack.h @@ -51,9 +51,9 @@ class EMFwdTrack float pz() const { return fPt * std::sinh(fEta); } float signed1Pt() const { return fSign / fPt; } - float cXX() const { return fCXX; } - float cXY() const { return fCXY; } - float cYY() const { return fCYY; } + float cXXatDCA() const { return fCXX; } + float cXYatDCA() const { return fCXY; } + float cYYatDCA() const { return fCYY; } protected: float fPt; diff --git a/PWGEM/Dilepton/Utils/EMTrackUtilities.h b/PWGEM/Dilepton/Utils/EMTrackUtilities.h index c06903b1eee..09cdd18f820 100644 --- a/PWGEM/Dilepton/Utils/EMTrackUtilities.h +++ b/PWGEM/Dilepton/Utils/EMTrackUtilities.h @@ -108,9 +108,9 @@ float dcaZinSigma(T const& track) template float fwdDcaXYinSigma(T const& track) { - float cXX = track.cXX(); // in cm^2 - float cYY = track.cYY(); // in cm^2 - float cXY = track.cXY(); // in cm^2 + float cXX = track.cXXatDCA(); // in cm^2 + float cYY = track.cYYatDCA(); // in cm^2 + float cXY = track.cXYatDCA(); // in cm^2 float dcaX = track.fwdDcaX(); // in cm float dcaY = track.fwdDcaY(); // in cm float det = cXX * cYY - cXY * cXY; // determinant From 30f17dbc5861d633a3c07ca3e6243c945feab10f Mon Sep 17 00:00:00 2001 From: ALICE Action Bot Date: Tue, 21 Jul 2026 15:21:37 +0000 Subject: [PATCH 2/2] Please consider the following formatting changes --- PWGEM/Dilepton/DataModel/dileptonTables.h | 2 -- 1 file changed, 2 deletions(-) diff --git a/PWGEM/Dilepton/DataModel/dileptonTables.h b/PWGEM/Dilepton/DataModel/dileptonTables.h index 6074a1e66f8..78ebebbed7f 100644 --- a/PWGEM/Dilepton/DataModel/dileptonTables.h +++ b/PWGEM/Dilepton/DataModel/dileptonTables.h @@ -1154,7 +1154,6 @@ DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuons_004, "AOD", "EMPRIMARYMU", 4, //! emprimarymuon::Py, emprimarymuon::Pz); - using EMPrimaryMuons = EMPrimaryMuons_004; // iterators using EMPrimaryMuon = EMPrimaryMuons::iterator; @@ -1177,7 +1176,6 @@ DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuonsCov_002, "AOD", "EMPRIMARYMUCOV", 2, / o2::aod::fwdtrack::C1PtTgl, o2::aod::fwdtrack::C1Pt21Pt2); - DECLARE_SOA_TABLE_VERSIONED(EMPrimaryMuonsCov_003, "AOD", "EMPRIMARYMUCOV", 3, //! at PV. Signed1Pt, Tgl and Phi are in EMPrimaryMuons table. fwdtrack::X, fwdtrack::Y, fwdtrack::Z, o2::aod::fwdtrack::CXX,