2020#include < chrono>
2121#include < cmath>
2222#include < cstdlib>
23+ #include < gsl/span>
2324#include < limits>
2425#include < numeric>
2526#include < queue>
2829#include < string>
2930#include < utility>
3031
32+ #include " CommonConstants/MathConstants.h"
3133#include " Framework/Logger.h"
3234#include " GPUCommonMath.h"
3335#include " ITSMFTTracking/BoundedAllocator.h"
3436#include " ITSMFTTracking/IndexTableConfiguration.h"
3537#include " ITSMFTTracking/MaterialPhysics.h"
36- #include " ITSMFTTracking/detail/TrackerTraversalPreparation .h"
38+ #include " ITSMFTTracking/MathUtils .h"
3739
3840namespace o2 ::itsmft::tracking
3941{
@@ -181,6 +183,50 @@ void prepareIterationConfiguration(const DetectorConfiguration& detector,
181183 validateSparsePlan (configuration, iteration, topology);
182184}
183185
186+ float diskLayerMultipleScatteringAngle (float layerxX0, float layerRadius, float referenceCoordinate, float trackletMinPt)
187+ {
188+ const float invP = 1 .f / trackletMinPt;
189+ const float tanlRef = (std::abs (layerRadius) > 1e-6f )
190+ ? referenceCoordinate / layerRadius
191+ : 0 .f ;
192+ const float absTanl = std::abs (tanlRef);
193+ const float cscLambda = (absTanl > 1e-6f )
194+ ? std::sqrt (1 .f + tanlRef * tanlRef) / absTanl
195+ : 1e6f;
196+ return 0 .0136f * invP * std::sqrt (layerxX0 * cscLambda);
197+ }
198+
199+ float clampEdgeCurvature (float oneOverR, float outerRadius) noexcept
200+ {
201+ return (outerRadius > 0 .f && 0 .5f * oneOverR >= 1 .f / outerRadius)
202+ ? (2 .f / outerRadius) - o2::constants::math::Almost0
203+ : oneOverR;
204+ }
205+
206+ struct EdgeScatteringBendingPrep {
207+ float msAngle;
208+ float phiCut;
209+ };
210+
211+ EdgeScatteringBendingPrep prepareEdgeScatteringAndBending (
212+ gsl::span<const float > perLayerMSAngle, int fromLayer, int toLayer,
213+ float r1, float r2, float clampedOneOverR, float res1, float res2) noexcept
214+ {
215+ float ms2 = 0 .f ;
216+ for (int layer = fromLayer; layer < toLayer; ++layer) {
217+ ms2 += o2::its::math_utils::Sq (perLayerMSAngle[layer]);
218+ }
219+ const float msAngle = o2::gpu::CAMath::Sqrt (ms2);
220+ const float cosTheta1half = o2::gpu::CAMath::Sqrt (1 .f - o2::its::math_utils::Sq (0 .5f * r1 * clampedOneOverR));
221+ const float cosTheta2half = o2::gpu::CAMath::Sqrt (1 .f - o2::its::math_utils::Sq (0 .5f * r2 * clampedOneOverR));
222+ const float x = (r2 * cosTheta1half) - (r1 * cosTheta2half);
223+ const float delta = o2::gpu::CAMath::Sqrt (1 .f / (1 .f - 0 .25f * o2::its::math_utils::Sq (x * clampedOneOverR)) *
224+ (o2::its::math_utils::Sq ((0 .25f * r1 * r2 * o2::its::math_utils::Sq (clampedOneOverR) / cosTheta2half) + cosTheta1half) * o2::its::math_utils::Sq (res1) +
225+ o2::its::math_utils::Sq ((0 .25f * r1 * r2 * o2::its::math_utils::Sq (clampedOneOverR) / cosTheta1half) + cosTheta2half) * o2::its::math_utils::Sq (res2)));
226+ const float phiCut = o2::gpu::CAMath::Min (o2::gpu::CAMath::ASin (0 .5f * x * clampedOneOverR) + 2 .f * msAngle + delta, o2::constants::math::PI * 0 .5f );
227+ return {msAngle, phiCut};
228+ }
229+
184230void prepareTraversalEdgeTolerances (
185231 IterationContext& context,
186232 int iteration)
@@ -194,15 +240,13 @@ void prepareTraversalEdgeTolerances(
194240 std::array<float , MaxLayoutSurfaces> msAngles{};
195241 for (int iLayer{0 }; iLayer < layerCount; ++iLayer) {
196242 const auto surface = LayerId{static_cast <uint16_t >(iLayer)};
197- if ( topology.getSurface (surface). kind == SurfaceKind::Cylinder) {
198- msAngles[iLayer] = cylinderLayerMultipleScatteringAngle (
199- CylinderLayerScatteringInputs{topology. getSurface (surface). material . xOverX0 } , trkParam.TrackletMinPt );
243+ const auto & descriptor = topology.getSurface (surface);
244+ if (descriptor. kind == SurfaceKind::Cylinder) {
245+ msAngles[iLayer] = o2::its::math_utils::MSangle ( 0 . 14f , trkParam.TrackletMinPt , descriptor. material . xOverX0 );
200246 } else {
201247 msAngles[iLayer] = diskLayerMultipleScatteringAngle (
202- DiskLayerScatteringInputs{topology.getSurface (surface).material .xOverX0 ,
203- context.detectorConfiguration .getRepresentativeRadius (surface),
204- topology.getSurface (surface).referenceCoordinate },
205- trkParam.TrackletMinPt );
248+ descriptor.material .xOverX0 , context.detectorConfiguration .getRepresentativeRadius (surface),
249+ descriptor.referenceCoordinate , trkParam.TrackletMinPt );
206250 }
207251 }
208252
@@ -225,7 +269,7 @@ void prepareTraversalEdgeTolerances(
225269 const float edgeOneOverR = clampEdgeCurvature (oneOverR, r2);
226270 const float res1 = o2::gpu::CAMath::Hypot (trkParam.PVres , context.detectorConfiguration .positionResolutions [fromLayer]);
227271 const float res2 = o2::gpu::CAMath::Hypot (trkParam.PVres , context.detectorConfiguration .positionResolutions [toLayer]);
228- const auto prep = :: o2::itsmft::tracking:: prepareEdgeScatteringAndBending (
272+ const auto prep = prepareEdgeScatteringAndBending (
229273 gsl::span<const float >(msAngles.data (), static_cast <std::size_t >(layerCount)), fromLayer, toLayer, r1, r2, edgeOneOverR, res1, res2);
230274 edgeMSAngles[*edgeSlot] = prep.msAngle ;
231275 edgePhiCuts[*edgeSlot] = prep.phiCut ;
0 commit comments