Bug in trapezium rule integrations in calculate_bfield_phase_error function in field_generator.py

### Problem description

- Trapezium rule to compute the second integral of motion from the bfield tensor.
- Bfield tensor is shape `(eval_X, eval_Z, eval_S, 3)` where `3` refers to slices for the measured field strength in `X`, `Z`, and `S` axes at each coordinate of the evaluation mesh sampled at points `x,z,s in (eval_X, eval_Z, eval_S)`.
- **L169** - Rolling by `1` index in the `X` axis constructs trapezoidal regions in the `XZ` plane staggered along the `S` axis.
- This should be roll by `1` index in the `S` axis (axis `2`).
https://github.com/DiamondLightSource/Opt-ID/blob/3c7e44951ddb73360cadb58f883bc5366ed2c9c7/IDSort/src/field_generator.py#L169
- **L174-L175** - The two cumulative sums are performed correctly along the S axis but the `X` and `Z` integrals are swapped in the order they are stored, and the `Z` integral is negated.
https://github.com/DiamondLightSource/Opt-ID/blob/3c7e44951ddb73360cadb58f883bc5366ed2c9c7/IDSort/src/field_generator.py#L174-L175
- Neither the flip or negation are needed in any downstream use-case, but the fact they have always been in the code means that expected output data files used by the smoke tests expect to see tensors packed in that order and negation.

- Trapezium rule to compute the first integral of motion from the second integral tensor.
- Second integral tensor already has `X` and `Z` integrals swapped and `Z` integral negated, so these properties persist transparently through this integration for the first integral of motion.
- **L177** - Rolling by `4` indices in the `X` axis constructs extra wide trapezoidal regions in the `XZ` plane staggered along the `S` axis.
- This should be roll by `1` index in the `S` axis (axis `2`) so that trapezoidal regions are constructed in the `ZS`.
https://github.com/DiamondLightSource/Opt-ID/blob/3c7e44951ddb73360cadb58f883bc5366ed2c9c7/IDSort/src/field_generator.py#L177

---

### Explanation for why results are still usable and are in the ballpark of the correct answer

- In the insertion devices the fields oscillate the electron beam strongly in the `Z` axis along the length of the `S` axis.
- Divergence in the `X` axis is small, especially for linearly polarized devices such as `Hybrid Symmetric CPMU` and `PPM Anti Symmetric`.
- This means that the results will be incorrect by a more significant factor for a bfield tensor representing an `APPLE-II Symmetric` device **where the field was measured / simulated in helical mode for circular polarization.**
- In practice we simulate `APPLE-II Symmetric` devices in standard mode for linear polarization so the results would still be close in this case.
- Because divergence in field strength is small across the `X` axis along the `S` axis, the computation of the integrals of motion along `S` will be dominated by the field strength `Fx` or `Fz`.
- **This means the result will be somewhat close to the correct integral.**

- Current code with bug that integrates the area of slices taken in the `XF` plane integrated along `S`. Where `F` is the field strength `Fx` or `Fz` at that location in `XZS` (There is also the measured field strength `Fs` at each point but we only consider integrals of motion in `X` and `Z`).

- Proposed fix that integrates the area of slices taken in the `SF` plane integrated along `S`. Where `F` is the field strength `Fx` or `Fz` at that location in `XZS` (There is also the measured field strength `Fs` at each point but we only consider integrals of motion in `X` and `Z`).

---

### Proposed fix

PR #60 merges in the following code that refactors and documents the `calculate_bfield_phase_error` function. This refactor retains the current bugs intentionally and lays out a strategy to remove them cleanly.

https://github.com/JossWhittle/Opt-ID/blob/0cacbc275876197f86c28b5a390d4c2ed64ff15f/IDSort/src/field_generator.py#L160-L234

- [L181](https://github.com/JossWhittle/Opt-ID/blob/0cacbc275876197f86c28b5a390d4c2ed64ff15f/IDSort/src/field_generator.py#L181)
- Refactor to become `trap_bfield = np.roll(bfield[...,:2], shift=1, axis=2)`
- [L191](https://github.com/JossWhittle/Opt-ID/blob/0cacbc275876197f86c28b5a390d4c2ed64ff15f/IDSort/src/field_generator.py#L191)
- Refactor to become `trap_traj_2nd_integral = np.roll(traj_2nd_integral, shift=1, axis=2)`
- [L203](https://github.com/JossWhittle/Opt-ID/blob/0cacbc275876197f86c28b5a390d4c2ed64ff15f/IDSort/src/field_generator.py#L203)
- Refactor to become `trajectories = np.concatenate([traj_1st_integral, traj_2nd_integral], axis=-1)`

Applying these fixes has the expected effect of lowering the observed phase error from passing in the field simulated from a `Hybrid Symmetric CPMU` device built with perfect magnets. There is always some phase error due to the inclusion of differently sized end `HE` and kicker `HT` magnets.
- Phase error with bugs `0.35...`
- Phase error with fixes `0.30...`

Merging these fixes will require regeneration and validation of expected output data for the smoke tests of multiple test suites and will likely alter the deterministic behavior of the `mpi_runner_test` and `mpi_runner_shim_opt_test`.

However, as we are on the cusp of re-modularization that will certainly break these smoke tests from being perfect consistent anyway, this might be the best time to apply this fix and correct for it before that process begins.