fix: Remove incomplete initialisation from InterpolatedBFieldMap

[timadye]

• Remove unneeded/misleading partial initialisation.
  • Note std::array<Vector3, 8> neighbors{Vector3{0.0, 0.0, 0.0}}; only initialises the first element of the std::array.
• assert if neighbors array not fully filled.

Originally "fix: B-field interpolation errors", but this PR no longer provides a fix to the following problem.

• Remove InterpolatedBFieldMap boundary check. Instead, use the under/overflow bins to return the B-field outside the grid.
  • Was generating ~259 errors per event (full_chain_itk.py ttbar+PU200), with each error generating 2 error messages (one quite verbose) and 4 warnings.
  • We don't need to propagate these track candidates further, but we should leave that decision up to the caller and not flood the logfile.

[timadye]

Label WiP because it may take some time for deal with the fallout from a unit test failure. The ATLAS B-field mostly works fine because the under/overflow bins all have B=(0,0,0).

The unit test expects the error, so can be fixed by expecting ok. But adding some tests for B-field outside the range hit issues with the Grid methods that mostly do not support the user/overflow bins. In particular, the neighboring bin finding doesn't work for underflow bins (can't return bin=-1 in an unsigned int, so doesn't fully initialise the array of neighbors). Also the interpolation doesn't do extrapolation.

I think I can work round these issues and improve the unit test for these cases. But maybe this isn't the right way to go about it.

Stepping back for a moment, it seems to me that InterpolatedBFieldMap should be able to cope with requests outside its range, either by returning B=(0,0,0) or the edge values. The former would be simpler, but the latter fits the binned map model.

[paulgessinger]

Thanks first of all.

I don't think the magnetic field interpolation should fail silently in this way, to be honest. If the particle propagation goes out of bounds of the interpolation domain, that's a serious issue. It either means that your field map is not sufficiently large or there's a navigation error for some reason.

I would personally prefer encouraging users to fix either the geometry (if that's the source of the navigation issues), fix a bug in the navigation, or explicitly extent the interpolation domain.

Thoughts?

[timadye]

I don't think the magnetic field interpolation should fail silently in this way, to be honest.

InterpolatedBFieldMap defines underflow and overflow bins, so there is a clear definition for values outside the range. So I think it is a matter of interpretation whether this is a failure or just a limited domain of implementation. Other field map types (constant, solenoid) don't impose a limited domain except where they make no sense (eg. R<0).

If the particle propagation goes out of bounds of the interpolation domain, that's a serious issue.

For me, "if it goes out of the bounds of the detector, that's a serious issue". The magnetic field goes on forever! (even if it is 0). In this case, it may have identified an error, but that's using the field map to cross-check the propagation. Perhaps useful, but a bit ugly.

Sorry for the philosophy. For my bug (#1484), I'd be happy with anything that gets rid of the GBytes of logfile messages.

[andiwand]

if the overflow was designed for this case I would agree to drop the error in this case. but we should make sure that we still error if there is no overflow value.

the underlying problem #1484 looks like a navigation failure to me so if that would be fixed the magnetic field should be fine right?

[paulgessinger]

While the InterpolatedBField indeed has under and overflow values, we had issues in the past where it went out of bounds and we didn't notice and only after non-negligible debugging did we figure out that was the case.
I'll try to dig up the issues.

[paulgessinger]

Got it, it's here #784 and then I introduced this error mode in #825.

I also had some discussion with @HadrienG2 about this in #821.

One technical issue is that for the clamped B field to be correct, you'd have to synthesize grid points out of bounds of the grid itself, since the interpolation uses the lower left corner of the grid bin, and the normal neighbor mechanism doesn't return neighbors beyond the overflow bin I don't think.

Overall I don't think this is the way to go: in my opinion this should be an explicit and loud error, and either the navigation or the field map should be fixed. If this requires changing the CKF to behave correctly in this case, I think that's what we'll have to do.

If we move ahead with this, I think at the very least this PR should roll back all the changes I made back in #825.

[andiwand]

makes sense @paulgessinger but if we do not intend to extrapolate the B field at all the overflow functionality sounds kind of useless to me?

if we want to have a defined B fields in case of overflow the more appropriate point to check for navigation failures would be in the navigator IMO. we could require an end of world bounds in the config and check at each iteration if we stepped into the void and error there which would be more explicit

[paulgessinger]

We use the overflow functionality of Grid in the SurfaceArray. If we can make it optional, that would be a possibility.

To your other point: fully agree that ideally the navigation should catch this. Isn't that what the EndOfWorld aborter is supposed to do?

[andiwand]

you are right EndOfWorld does the job usually only CKF sneaks in because you cannot use standard aborters there. but I guess that is a detail of CKF and has to be dealt with there specifically

[timadye]

Got it, it's here #784 and then I introduced this error mode in #825.
I also had some discussion with @HadrienG2 about this in #821.

Sorry, I didn't know that this had already been widely discussed and the current behavior was (fairly) recently introduced. I might have argued for your options 1 or 2, but not sufficiently strongly to warrant reopening the issue.

I notice that #825 introduced the initialisation in InterpolatedBFieldMap.hpp line 197:

    std::array<Vector3, nCorners> neighbors{Vector3{0.0, 0.0, 0.0}};

I think that's a mistake, which I tried to correct in this PR (not yet quite right there). It only initialises the first element of neighbors. The array will be filled later, as long as the closestPointsIndices iterator yields nCorners values. I think it should always do that, except in the MagneticFieldError::OutOfBounds case (actually, IIRC, I found only in the lower or left bounds, while the upper and right bounds still give the correct number of indices).

Anyway, I think it would be faster to remove the initialisation, or perhaps safer to use std::array::fill to fill all elements. Either would be better than a deeply misleading first-value initialisation. I checked on Godbolt, and fill (with Vector3{0,0,0}) after declaration is just as fast as initialising all elements, but much easier to do for constexpr nCorners.

While we are looking here, I see that there may be an optimisation opportunity to replace the isInsideLocal, which compares the coordinate position with the double bounds, with comparisons on the indices. I don't know if that's also safer in edge cases.

[paulgessinger]

Hey @timadye, that sounds fair. I think I didn't realize back then that this only initializes the first value.

I guess we can just drop the initialization if it's basically guaranteed to be initialized later on.

[stale]

This issue/PR has been automatically marked as stale because it has not had recent activity. The stale label will be removed if any interaction occurs.

[andiwand]

how do we move forward with this @paulgessinger @timadye ?

[timadye]

how do we move forward with this @paulgessinger @timadye ?

Good question. On this PR, we can take just the first changed line, removing the incomplete but unneeded initialisation. Should I also add an assert to check all elements are actually filled in as they are supposed to be by the static iterator?

Probably that is enough for this PR. There was also a potential optimisation opportunity with isInsideLocal I mentioned above. Is that worth following up separately? Probably would need someone more expert than me.

The point of this PR was to remove all the error messages. I haven't checked recently if they are still there. I guess we can follow that up in #1484, but that issue is also stale.

[codecov]

Codecov Report

Merging #1467 (f596fd2) into main (e7d0ce9) will decrease coverage by 0.01%.
The diff coverage is 0.00%.

❗ Current head f596fd2 differs from pull request most recent head a9a4042. Consider uploading reports for the commit a9a4042 to get more accurate results

@@            Coverage Diff             @@
##             main    #1467      +/-   ##
==========================================
- Coverage   49.32%   49.31%   -0.01%     
==========================================
  Files         441      441              
  Lines       25208    25209       +1     
  Branches    11626    11627       +1     
==========================================
  Hits        12433    12433              
  Misses       4513     4513              
- Partials     8262     8263       +1     

Impacted Files	Coverage Δ
...clude/Acts/MagneticField/InterpolatedBFieldMap.hpp	55.07% <0.00%> (-0.81%)	⬇️

📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more

[github-actions]

📊 Physics performance monitoring for a9a4042

Summary
Full report
Seeding: seeded, truth estimated, orthogonal
CKF: seeded, truth smeared, truth estimated, orthogonal
IVF: seeded, truth smeared, truth estimated, orthogonal
AMVF: seeded, truth smeared, truth estimated, orthogonal
Ambiguity resolution: seeded, orthogonal
Truth tracking
Truth tracking (GSF)

Vertexing

Vertexing vs. mu

IVF seeded

IVF truth_smeared

IVF truth_estimated

IVF orthogonal

AMVF seeded

AMVF truth_smeared

AMVF truth_estimated

AMVF orthogonal

Seeding

Seeding seeded

Seeding truth_estimated

Seeding orthogonal

CKF

CKF seeded

CKF truth_smeared

CKF truth_estimated

CKF orthogonal

Ambiguity resolution

seeded

Truth tracking (Kalman Filter)

Truth tracking

Truth tracking (GSF)

Truth tracking