Demonstrate stuck particle in ATLAS

[sethrj]

This illustrates the stuck track found by @drbenmorgan in ATLAS documented in #1902 . Because of how VecGeom 1.x determines boundary crossings, it gets stuck when near a plane at a grazing angle. All of ORANGE, Geant4, and VecGeom surface (current master) handle the crossing without error (though vecgeom surface disagrees with Geant4 and VecGeom on some of the safety distances in the outer polycone).

Call sequence to VecGeom:

• Navigator::LocatePointIn
• Navigator::ComputeStepAndNextVolume
• (update position)
• SetBoundaryState(true)
• Navigator::RelocateToNextVolume

@agheata @SeverinDiederichs This can probably be reproduced in a VecGeom unit test using the included GDML and raytracing from {23.51934584635, 17.141066715148, 344.45000092904} along {0.5784236876658104, 0.8157365000698582, -9.290358099212079e-7}.

I haven't tried vecgeom 2 with surface disabled.

[github-actions]

Test summary

 4 559 files   7 281 suites   13m 18s ⏱️
 1 947 tests  1 934 ✅ 10 💤 3 ❌
25 043 runs  24 942 ✅ 97 💤 4 ❌

For more details on these failures, see this check.

Results for commit baf2576.

♻️ This comment has been updated with latest results.

[mrguilima]

Looks good to me.
Thanks!

[agheata]

@SethJ, using the master of today I've checked:
test/testRaytracing -gdml_name /home/agheata/modeler/surface/atlas-hgtd.gdml -nrays 1 -debug 1 -ongpu 0 -test-bvh 0 -use_surf 0 -point 23.51934584635, 17.141066715148, 344.45000092904 -direction 0.5784236876658104, 0.8157365000698582, -9.290358099212079e-7

which traverses the geometry w/o errors:

info: Tracking single ray with point (23.5193, 17.1411, 344.45) and direction (0.578424, 0.815737, -9.29036e-07)
PropagateRaysSolid debug ray 0 : p{(23.51934584635, 17.141066715148, 344.45000092904)} d{(0.5784236876658105, 0.8157365000698584, -9.290358099212082e-07)}
   start :navInd=s0:11 lastExited=0, id=2, level=1/4,  onBoundary=false, path=</1/11>
     dist = 1120.375825265075
   1 : navInd=s0:1 lastExited=11, id=1, level=0/4,  onBoundary=true, path=</1>
     dist = 13053.03440810621
   2 : navInd=0, id=0, path=outside

so this probably does not hit the volume you expect?

[sethrj]

Thanks for looking @agheata ; the units I posted were cm not mm. Could these be the discrepancy?

[agheata]

Thanks for looking @agheata ; the units I posted were cm not mm. Could these be the discrepancy?

Likely. Running with: test/testRaytracing -gdml_name /home/agheata/modeler/surface/atlas-hgtd.gdml -mmunit 0.1 -nrays 1 -debug 1 -ongpu 0 -test-bvh 0 -use_surf 0 -point 23.51934584635, 17.141066715148, 344.45000092904 -direction 0.5784236876658104, 0.8157365000698582, -9.290358099212079e-7

I get (printing volume names now):

info: Tracking single ray with point (23.5193, 17.1411, 344.45) and direction (0.578424, 0.815737, -9.29036e-07)
PropagateRaysSolid debug ray 0 : p{(23.51934584635, 17.141066715148, 344.45000092904)} d{(0.5784236876658105, 0.8157365000698584, -9.290358099212082e-07)}
   start :navInd=s0:44 lastExited=0, id=4, level=3/4,  onBoundary=false, path=</Atlas0x4174ac60_PV/ITK0x4ffa3e50/HGTD0x65e4d130/SPlate_40x602d7220>
     dist = 1.000004527107722
   1 : navInd=s0:21 lastExited=44, id=3, level=2/4,  onBoundary=true, path=</Atlas0x4174ac60_PV/ITK0x4ffa3e50/HGTD0x65e4d130>
     dist = 81.02188807286879
   2 : navInd=s0:11 lastExited=21, id=2, level=1/4,  onBoundary=true, path=</Atlas0x4174ac60_PV/ITK0x4ffa3e50>
     dist = 4.816418570535086
   3 : navInd=s0:1 lastExited=11, id=1, level=0/4,  onBoundary=true, path=</Atlas0x4174ac60_PV>
     dist = 1305.644686893262
   4 : navInd=0, id=0, path=outside
info: Ray 0 : p{(23.519345846349999, 17.141066715148, 344.45000092904002)} d{(0.57842368766581054, 0.81573650006985843, -9.2903580992120816e-07)} crossed 4 boundaries

Which starts from a different volume, but there are no errors.
The issue you see may be coming from arriving to that point with a different state, so what needs to be debugged is why that happened (i.e. the error happened already before)

[sethrj]

@agheata Does your build have CUDA enabled? Our CI seems to suggest that, somehow, we no longer see the broken behavior on CPU when CUDA is enabled in vecgeom.

[agheata]

@agheata Does your build have CUDA enabled? Our CI seems to suggest that, somehow, we no longer see the broken behavior on CPU when CUDA is enabled in vecgeom.

I tried with CUDA disabled, and I am getting the same. Can you try on your side: test/testRaytracing -gdml_name /home/agheata/modeler/surface/atlas-hgtd.gdml -mmunit 0.1 -nrays 1 -debug 1 -ongpu 0 -test_bvh 1 -use_surf 0 -point 23.51934584635, 17.141066715148, 344.45000092904 -direction 0.5784236876658104, 0.8157365000698582, -9.290358099212079e-7 ? I tried also with -ongpu 1 -use_surf 1 and the surface model gives the same output.

I could not try with CUDA OFF and SURF OFF because the raytracer does not run in this case due to some missing protection.