-
Notifications
You must be signed in to change notification settings - Fork 615
[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.
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
base: master
Are you sure you want to change the base?
Changes from 6 commits
73c6a33
528c226
eb6575a
f0dcd44
59823f3
f7320e1
fda90eb
f2d6b2b
23808fc
202fb8c
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -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
|
||
| // 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"}; | ||
|
|
@@ -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 | ||
| if (device.name == "propagation-service") { | ||
| // loop over the options | ||
| // to find the value of TrackTuner::autoDetectDcaCalib | ||
|
|
@@ -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(); | ||
|
|
@@ -287,29 +288,51 @@ | |
| } | ||
| } // 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 (fillTracksDCA || fillTracksDCACov) { | ||
| if (track.has_collision()) { | ||
| auto const& collision = collisions.rawIteratorAt(track.collisionId()); | ||
| mVtx.setPos({collision.posX(), collision.posY(), collision.posZ()}); | ||
| if (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 (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 (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 | ||
|
|
@@ -356,6 +379,67 @@ | |
| } | ||
| } | ||
| } | ||
|
|
||
| //_______________________________________________________________________ | ||
| // 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 | ||
| double sn, cs, alp = trackPar.getAlpha(); | ||
| o2::math_utils::detail::sincos(alp, sn, cs); | ||
| double x = trackPar.getX(), y = trackPar.getY(), snp = trackPar.getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); | ||
| double xv = vtx.getX() * cs + vtx.getY() * sn, yv = -vtx.getX() * sn + vtx.getY() * cs, zv = vtx.getZ(); | ||
| x -= xv; | ||
| y -= yv; | ||
| // Estimate the impact parameter neglecting the track curvature | ||
| double d = gpu::CAMath::Abs(x * snp - y * csp); | ||
| if (d > maxD) { | ||
| // 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; | ||
| } | ||
|
||
| double crv = trackPar.getCurvature(b); | ||
| double 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
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe 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.
Collaborator
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I actually tried different approaches:
)
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?
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I see that some analyses are using
Collaborator
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Hi @shahor02 ! Apologies for the late reply! Also, I repeated the comparison focusing only on collisions without any other compatible PV (by requiring It seems that, indeed, I was calculating the DCA of some prepropagated tracks with respect to another compatible PV.
Contributor
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Hi @romainschotter , Concerning this:
Do you mean 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; |
||
| // 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 | ||
|
|
||
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
What is the reason to pass the trackPar by copy rather than by reference?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Hi @shahor02 ! That's a mistake from my side, thank you for pointing that out! I will change it tomorrow morning