Skip to content

[Common] Consider propagation of Tracks in TrackPropagationModule #14114

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
romainschotter wants to merge 10 commits into
base: master
Choose a base branch
from
+101 −19
Draft

[Common] Consider propagation of Tracks in TrackPropagationModule #14114

Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
73c6a33
Consider propagation of Tracks + add method to calculate DCA without …
romainschotter Nov 28, 2025
528c226
Please consider the following formatting changes
alibuild Dec 8, 2025
eb6575a
Merge pull request #62 from alibuild/alibot-cleanup-14114
romainschotter Dec 8, 2025
f0dcd44
Add TrackTuner for prepropagated tracks
romainschotter Dec 9, 2025
59823f3
Please consider the following formatting changes
alibuild Dec 9, 2025
f7320e1
Merge pull request #63 from alibuild/alibot-cleanup-14114
romainschotter Dec 9, 2025
fda90eb
Pass TrackPar by reference in calculateDCA
romainschotter Dec 10, 2025
f2d6b2b
Add propagation to PV for preprogated tracks if track tuner used
romainschotter Dec 12, 2025
23808fc
Please consider the following formatting changes
alibuild Dec 12, 2025
202fb8c
Merge pull request #64 from alibuild/alibot-cleanup-14114
romainschotter Dec 12, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
120 changes: 101 additions & 19 deletions Common/Tools/TrackPropagationModule.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -71,7 +71,7 @@

struct TrackPropagationConfigurables : o2::framework::ConfigurableGroup {
std::string prefix = "trackPropagation";
o2::framework::Configurable<float> minPropagationRadius{"minPropagationDistance", o2::constants::geom::XTPCInnerRef + 0.1, "Only tracks which are at a smaller radius will be propagated, defaults to TPC inner wall"};

Check failure on line 74 in Common/Tools/TrackPropagationModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/configurable] 

Use lowerCamelCase for names of configurables and use the same name for the struct member as for the JSON string. (Declare the type and names on the same line.)
// for TrackTuner only (MC smearing)
o2::framework::Configurable<bool> useTrackTuner{"useTrackTuner", false, "Apply track tuner corrections to MC"};
o2::framework::Configurable<bool> useTrkPid{"useTrkPid", false, "use pid in tracking"};
Expand Down Expand Up @@ -176,7 +176,7 @@
/// to understand whether the TrackTuner::getDcaGraphs function can be called here (input path from string/configurables)
/// or inside the process function, to "auto-detect" the input file based on the run number
const auto& workflows = initContext.services().template get<o2::framework::RunningWorkflowInfo const>();
for (const o2::framework::DeviceSpec& device : workflows.devices) { /// loop over devices

Check failure on line 179 in Common/Tools/TrackPropagationModule.h

View workflow job for this annotation

GitHub Actions / O2 linter

[const-ref-in-for-loop] 

Use constant references for non-modified iterators in range-based for loops.
if (device.name == "propagation-service") {
// loop over the options
// to find the value of TrackTuner::autoDetectDcaCalib
Expand Down Expand Up @@ -271,8 +271,9 @@
// std::array<float, 3> trackPxPyPz;
// std::array<float, 3> trackPxPyPzTuned = {0.0, 0.0, 0.0};
double q2OverPtNew = -9999.;
bool isPropagationOK = true;
// Only propagate tracks which have passed the innermost wall of the TPC (e.g. skipping loopers etc). Others fill unpropagated.
if (track.trackType() == o2::aod::track::TrackIU && track.x() < cGroup.minPropagationRadius.value) {
if (track.x() < cGroup.minPropagationRadius.value) {
if (fillTracksCov) {
if constexpr (isMc) { // checking MC and fillCovMat block begins
// bool hasMcParticle = track.has_mcParticle();
Expand All @@ -287,29 +288,61 @@
}
} // MC and fillCovMat block ends
}
bool isPropagationOK = true;

if (track.has_collision()) {
auto const& collision = collisions.rawIteratorAt(track.collisionId());
if (fillTracksCov) {
mVtx.setPos({collision.posX(), collision.posY(), collision.posZ()});
mVtx.setCov(collision.covXX(), collision.covXY(), collision.covYY(), collision.covXZ(), collision.covYZ(), collision.covZZ());
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz(mVtx, mTrackParCov, 2.f, matCorr, &mDcaInfoCov);

if (track.trackType() == o2::aod::track::TrackIU) {
if (track.has_collision()) {
auto const& collision = collisions.rawIteratorAt(track.collisionId());
if (fillTracksCov) {
mVtx.setPos({collision.posX(), collision.posY(), collision.posZ()});
mVtx.setCov(collision.covXX(), collision.covXY(), collision.covYY(), collision.covXZ(), collision.covYZ(), collision.covZZ());
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz(mVtx, mTrackParCov, 2.f, matCorr, &mDcaInfoCov);
} else {
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz({collision.posX(), collision.posY(), collision.posZ()}, mTrackPar, 2.f, matCorr, &mDcaInfo);
}
} else {
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz({collision.posX(), collision.posY(), collision.posZ()}, mTrackPar, 2.f, matCorr, &mDcaInfo);
if (fillTracksCov) {
mVtx.setPos({ccdbLoader.mMeanVtx->getX(), ccdbLoader.mMeanVtx->getY(), ccdbLoader.mMeanVtx->getZ()});
mVtx.setCov(ccdbLoader.mMeanVtx->getSigmaX() * ccdbLoader.mMeanVtx->getSigmaX(), 0.0f, ccdbLoader.mMeanVtx->getSigmaY() * ccdbLoader.mMeanVtx->getSigmaY(), 0.0f, 0.0f, ccdbLoader.mMeanVtx->getSigmaZ() * ccdbLoader.mMeanVtx->getSigmaZ());
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz(mVtx, mTrackParCov, 2.f, matCorr, &mDcaInfoCov);
} else {
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz({ccdbLoader.mMeanVtx->getX(), ccdbLoader.mMeanVtx->getY(), ccdbLoader.mMeanVtx->getZ()}, mTrackPar, 2.f, matCorr, &mDcaInfo);
}
}
if (isPropagationOK) {
trackType = o2::aod::track::Track;
}
} else {
if (fillTracksCov) {
mVtx.setPos({ccdbLoader.mMeanVtx->getX(), ccdbLoader.mMeanVtx->getY(), ccdbLoader.mMeanVtx->getZ()});
mVtx.setCov(ccdbLoader.mMeanVtx->getSigmaX() * ccdbLoader.mMeanVtx->getSigmaX(), 0.0f, ccdbLoader.mMeanVtx->getSigmaY() * ccdbLoader.mMeanVtx->getSigmaY(), 0.0f, 0.0f, ccdbLoader.mMeanVtx->getSigmaZ() * ccdbLoader.mMeanVtx->getSigmaZ());
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz(mVtx, mTrackParCov, 2.f, matCorr, &mDcaInfoCov);
} else {
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz({ccdbLoader.mMeanVtx->getX(), ccdbLoader.mMeanVtx->getY(), ccdbLoader.mMeanVtx->getZ()}, mTrackPar, 2.f, matCorr, &mDcaInfo);
if (fillTracksCov || fillTracksDCA || fillTracksDCACov) {
if (track.has_collision()) {
auto const& collision = collisions.rawIteratorAt(track.collisionId());
mVtx.setPos({collision.posX(), collision.posY(), collision.posZ()});
if (fillTracksCov || fillTracksDCACov) {
mVtx.setCov(collision.covXX(), collision.covXY(), collision.covYY(), collision.covXZ(), collision.covYZ(), collision.covZZ());
}
} else {
mVtx.setPos({ccdbLoader.mMeanVtx->getX(), ccdbLoader.mMeanVtx->getY(), ccdbLoader.mMeanVtx->getZ()});
if (fillTracksCov || fillTracksDCACov) {
mVtx.setCov(ccdbLoader.mMeanVtx->getSigmaX() * ccdbLoader.mMeanVtx->getSigmaX(), 0.0f, ccdbLoader.mMeanVtx->getSigmaY() * ccdbLoader.mMeanVtx->getSigmaY(), 0.0f, 0.0f, ccdbLoader.mMeanVtx->getSigmaZ() * ccdbLoader.mMeanVtx->getSigmaZ());
}
}
if (fillTracksCov) {
if constexpr (isMc) { // checking MC and fillCovMat block begins
if (cGroup.useTrackTuner.value) {
bool hasMcParticle = track.has_mcParticle();
if (hasMcParticle) {
isPropagationOK = o2::base::Propagator::Instance()->propagateToDCABxByBz(mVtx, mTrackParCov, 2.f, matCorr, &mDcaInfoCov);
}
}
} // MC and fillCovMat block ends
}

if (fillTracksDCACov) {
calculateDCA(mTrackParCov, mVtx, o2::base::Propagator::Instance()->getNominalBz(), &mDcaInfoCov, 999.f);
} else {
calculateDCA(mTrackPar, mVtx, o2::base::Propagator::Instance()->getNominalBz(), &mDcaInfo, 999.f);
}
}
}
if (isPropagationOK) {
trackType = o2::aod::track::Track;
}
// filling some QA histograms for track tuner test purpose
if (fillTracksCov) {
if constexpr (isMc) { // checking MC and fillCovMat block begins
Expand Down Expand Up @@ -356,6 +389,55 @@
}
}
}

//_______________________________________________________________________
// TTrackPar type: either aod::TrackPar or aod::TrackParCov
// TVertex type: either math_utils::Point3D<value_t> or o2::dataformats::VertexBase
// TDCA type: either dim2_t or o2::dataformats::DCA
template <typename TTrackPar, typename TVertex, typename TDCA>
bool calculateDCA(TTrackPar& trackPar, const TVertex& vtx, double b, TDCA* dca, double maxD)
{
// propagate track to DCA to the vertex
float sn, cs, alp = trackPar.getAlpha();
o2::math_utils::detail::sincos(alp, sn, cs);
float x = trackPar.getX(), y = trackPar.getY(), snp = trackPar.getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp));
float xv = vtx.getX() * cs + vtx.getY() * sn, yv = -vtx.getX() * sn + vtx.getY() * cs, zv = vtx.getZ();
x -= xv;
y -= yv;

float crv = trackPar.getCurvature(b);
float tgfv = -(crv * x - snp) / (crv * y + csp);
sn = tgfv / gpu::CAMath::Sqrt(1.f + tgfv * tgfv);
cs = gpu::CAMath::Sqrt((1.f - sn) * (1.f + sn));
cs = (gpu::CAMath::Abs(tgfv) > constants::math::Almost0) ? sn / tgfv : constants::math::Almost1;

x = xv * cs + yv * sn;
yv = -xv * sn + yv * cs;
xv = x;

alp += gpu::CAMath::ASin(sn);
if (!trackPar.rotate(alp)) {
Comment on lines +418 to +419
Copy link
Contributor

@shahor02 shahor02 Dec 12, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Did you simply suppress these 2 lines? If the track would not be rotated to the frame of the PCA after the propagation, the difference would be much larger than what your pdf shows.
Is it possible that you are related the same preprogated track to two different compatible PVs?

Copy link
Collaborator Author

@romainschotter romainschotter Dec 13, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I actually tried different approaches:

  • I considered only the rotation of the PV to the same frame as the track at the PCA (using the exactly the code you suggested 3days ago, i.e.
float sn, cs, alp = trackPar.getAlpha();
o2::math_utils::detail::sincos(alp, sn, cs);
float xv = vtx.getX() * cs + vtx.getY() * sn, yv = -vtx.getX() * sn + vtx.getY() * cs, zv = vtx.getZ();
if constexpr (requires { trackPar.getSigmaY2(); vtx.getSigmaX2(); }) {
  o2::math_utils::detail::sincos(alp, sn, cs);
  auto s2ylocvtx = vtx.getSigmaX2() * sn * sn + vtx.getSigmaY2() * cs * cs - 2. * vtx.getSigmaXY() * cs * sn;
  dca->set(trackPar.getY() - yv, trackPar.getZ() - zv, trackPar.getSigmaY2() + s2ylocvtx, trackPar.getSigmaZY(), trackPar.getSigmaZ2() + vtx.getSigmaZ2());
} else {
  (*dca)[0] = trackPar.getY() - yv;
  (*dca)[1] = trackPar.getZ() - zv;
}

)

  • I also tried removing the following line entirely:
if (!trackPar.rotate(alp)) {

In both cases, I got the same DCA distribution which corresponds to the one shown in open star markers in the PDF attached in my previous message.

I suppose it is possible that I am calculating the DCA of a prepropagated track with respect to another compatible PV. Is there a way to check that?

Copy link
Contributor

@shahor02 shahor02 Dec 13, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see that some analyses are using track.compatibleCollIds().size() > 0 to check if the track is ambiguous, but I don't know if getting this track attribute requires some special. Btw, in my code snippet, which you quoted, I've wrongly left unnecessary sincos call, should be

float sn, cs, alp = trackPar.getAlpha();
o2::math_utils::detail::sincos(alp, sn, cs);
float xv = vtx.getX() * cs + vtx.getY() * sn, yv = -vtx.getX() * sn + vtx.getY() * cs, zv = vtx.getZ();
if constexpr (requires { trackPar.getSigmaY2(); vtx.getSigmaX2(); }) {
  auto s2ylocvtx = vtx.getSigmaX2() * sn * sn + vtx.getSigmaY2() * cs * cs - 2. * vtx.getSigmaXY() * cs * sn;
  dca->set(trackPar.getY() - yv, trackPar.getZ() - zv, trackPar.getSigmaY2() + s2ylocvtx, trackPar.getSigmaZY(), trackPar.getSigmaZ2() + vtx.getSigmaZ2());
} else {
  (*dca)[0] = trackPar.getY() - yv;
  (*dca)[1] = trackPar.getZ() - zv;
}

Copy link
Collaborator Author

@romainschotter romainschotter Dec 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @shahor02 ! Apologies for the late reply!
I tried the correct code snippet this time (without the remaining unnecessary sincos call) but the DCA distribution looked the same as before, meaning that I still obtain the open star markers in the PDF attached in my previous message.

Also, I repeated the comparison focusing only on collisions without any other compatible PV (by requiring track.compatibleCollIds().size() > 0) and now the DCA distribution with your code snipped (without rotation) looks the same as with rotation and the same as without prepropagated tracks.
DCAxy_StdVsPrepropagate_WithoutRotation_withoutCompatibleColl_Strictly0.pdf

It seems that, indeed, I was calculating the DCA of some prepropagated tracks with respect to another compatible PV.
However, I am not entirely sure to understand how that's possible. I am wondering if it is not possible that the tracks are always prepropagated to the mean vertex (https://github.com/AliceO2Group/AliceO2/blob/4090041b401c7aa6c919ca923126fff950cbccd1/Detectors/AOD/src/AODProducerWorkflowSpec.cxx#L2858) instead of the PV of its assigned collision. If I am saying something wrong, please don't hesitate to correct me

Copy link
Contributor

@shahor02 shahor02 Dec 18, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @romainschotter ,
The extra sincos could not change the results, it was simply not needed.

Concerning this:

Also, I repeated the comparison focusing only on collisions without any other compatible PV (by requiring track.compatibleCollIds().size() > 0) and now the DCA distribution with your code snipped (without rotation) looks the same as with rotation and the same as without prepropagated tracks.

Do you mean track.compatibleCollIds().size() == 1 (for the tracks compatible with 1 PV only)?

The track is propagated to the meanvertex only if it is not associated to any real PV. But in case it is associated to >1 PV, the propagation will be done only to 1st one (since the track is saved only once).

But if this is the case, then in your method just rotation would not be enough;
in case of the ambiguous vertex assignment also the propagation would be needed, as it is done in case of trackTuning.
I am wondering if your reference distribution is not suffering from a similar problem.
I think it is worth comparing for ambiguous tracks only (i.e. track.compatibleCollIds().size() > 0) the DCA distribution with old propagate to dca, and your method with rotation only and with rotation + propagation.

// provide default DCA for failed propag
if constexpr (requires { trackPar.getSigmaY2(); vtx.getSigmaX2(); }) {
dca->set(o2::track::DefaultDCA, o2::track::DefaultDCA,
o2::track::DefaultDCACov, o2::track::DefaultDCACov, o2::track::DefaultDCACov);
} else {
(*dca)[0] = o2::track::DefaultDCA;
(*dca)[1] = o2::track::DefaultDCA;
}

return false;
}
if constexpr (requires { trackPar.getSigmaY2(); vtx.getSigmaX2(); }) {
o2::math_utils::detail::sincos(alp, sn, cs);
auto s2ylocvtx = vtx.getSigmaX2() * sn * sn + vtx.getSigmaY2() * cs * cs - 2. * vtx.getSigmaXY() * cs * sn;
dca->set(trackPar.getY() - yv, trackPar.getZ() - zv, trackPar.getSigmaY2() + s2ylocvtx, trackPar.getSigmaZY(), trackPar.getSigmaZ2() + vtx.getSigmaZ2());
} else {
(*dca)[0] = trackPar.getY() - yv;
(*dca)[1] = trackPar.getZ() - zv;
}
return true;
}
};

} // namespace common
Expand Down
Loading