Three Dimensional Polarimetric Neutron Tomography of Magnetic Fields

Through the use of Time-of-Flight Three Dimensional Polarimetric Neutron Tomography (ToF 3DPNT) we have for the first time successfully demonstrated a technique capable of measuring and reconstructing three dimensional magnetic field strengths and directions unobtrusively and non-destructively with the potential to probe the interior of bulk samples which is not amenable otherwise. Using a pioneering polarimetric set-up for ToF neutron instrumentation in combination with a newly developed tailored reconstruction algorithm, the magnetic field generated by a current carrying solenoid has been measured and reconstructed, thereby providing the proof-of-principle of a technique able to reveal hitherto unobtainable information on the magnetic fields in the bulk of materials and devices, due to a high degree of penetration into many materials, including metals, and the sensitivity of neutron polarisation to magnetic fields. The technique puts the potential of the ToF time structure of pulsed neutron sources to full use in order to optimise the recorded information quality and reduce measurement time.

where c = 4.632 × 10 14 T −1 m −2 is the Larmor constant, and B is the magnetic field strength. Using this we can map the strength of a magnetic field along a neutron flight path into a neutron spin precession angle, and repeating this for multiple tomographic projections we can reconstruct the magnetic field probed by the neutrons 2 .
In order to evaluate the potential of the technique we have chosen to measure the magnetic field generated by a current carrying solenoid, the magnetic field of which can be calculated for comparison, thereby providing the possibility for producing a solid proof-of-principle experiment investigating the capabilities of three dimensional magnetic field polarimetric neutron tomography (3DPNT). Previous experiments have successfully used monochromatic polarised neutrons beams for 2D imaging of magnetic fields 3 , as well as 2D time-resolved imaging of periodically changing magnetic fields with a microsecond resolution 4 . 3D reconstructions with an assumption of the sample magnetic field direction to be along a direction perpendicular to the neutron polarisation has been demonstrated as well 5 . In contrast, our technique has been developed in order to measure and reconstruct 3D magnetic fields of arbitrary direction and distribution. This provides a method able to investigate samples without any a priori knowledge of the magnetic field orientation needed. Furthermore, it is able to use the full potential of a polychromatic pulsed neutron beam 6 .

Experimental Set-up
To control the neutron spin direction before and after sample interaction, we use a polarimetric set-up where the initial neutron spin direction can be set to be parallel (or antiparallel) to either the x, y, or z direction (see Fig. 1). The neutron spin component along one of the same three axes can be analysed after the neutron has passed through the magnetic field of the sample 7 .
Our experiments were performed at RADEN, BL22, at J-PARC MLF, Japan 8 , with an instrumental set-up as described by Fig. 1. Four π/2 spin rotators and a π spin flipper are used to select directions of spin polarisation and analysis. The polariser and analyser are polarising supermirrors 9 and a microchannel plate timepix detector 10 with a 512 × 512 array of 55 × 55 μm 2 pixels and a temporal resolution of less than 1 μs was used for neutron detection. We measured with a pulsed polychromatic neutron beam using the time-of-flight (ToF) information to determine the neutron wavelength, with the current in the spin rotators synchronised with the neutron pulse in order to achieve the proper neutron spin rotation for all neutron wavelengths. The time to wavelength conversion was done using a measurement of the Bragg edges of a standard iron sample, which provided the flight path and time delay values required for the conversion of ToF values into neutron wavelength. The images for all the wavelengths were acquired simultaneously and no scanning through energies was required in our set-up utilising the pulsed structure of the neutron beam and high count rate capabilities of our ToF imaging detector. We rebinned our data to have a spatial binning of 10 × 10 pixels -providing a spatial resolution of ~1 mm -and a temporal binning of 0.4992 ms, corresponding to δλ/λ = 3.3% at λ = 3.2 Å.
The sample used was an aluminium solenoid of length L s = 1.55 cm, radius R = 0.55 cm, wire thickness of w t = 0.1 cm, with N = 13.5 windings, and carrying a current of I = 0.6 A.

Reconstruction
In order to reconstruct the measured magnetic field from the recorded polarimetric neutron tomographic data set, we have developed the reconstruction procedure presented in this section. As described in 11,12 tomographic reconstruction of a magnetic field is not as straight forward as standard attenuation tomography since the polarisation of a neutron beam passing through a region of various magnetic field directions and strengths cannot be calculated using a simple line integral because of the non-commuting properties of the neutron spin orientation along the path 3 . Which is the reason we measured projections between 0° and 360°, since neutron paths of opposite direction yield different outcomes. This can be exemplified by imaging a neutron with its spin direction along the y direction passing through a magnetic field region with field direction along x, and of such extent and strength that it will rotate the neutron spin 90° to be along the z direction. A second magnetic field region of same size and strength and with the magnetic field direction along y will further rotate the neutron spin to its final orientation along the x direction. If the order of the two magnetic field regions had been switched, the final neutron spin direction would have been along the z direction (the intermediate neutron spin direction being y). The polarisation was calculated from measured intensities as: The effect of the magnetic field is cumulative along the ray and is governed by an ordinary differential equation 13 . This means that the forward problem mapping the magnetic field, B, to the measurements is not linear and so cannot be inverted with a simple inverse Radon transform. In this paper we apply a transformation of the data and then linearise the problem about B = 0. The measured intensities are of course scalars but the effect of the magnetic field is a rotation matrix P ∈ SO(3) given by 7 : An open beam measurement of I y,y and I −y,y was used to measure the flipping ratio (FR) and correct for the non-perfect polarisation and spin manipulation in the instrumental set-up. The FR was measured to be 22, corresponding to a polarisation of 91%, averaged over the detector at λ = 3.2 Å.
In order to analyse the recorded data P was first transformed to a sample reference system. For each projection, the sample was rotated around the vertical axis, y, and for a given projection angle, θ, R y (−θ) is the rotation matrix: y with which we can calculate the polarisation matrix, P′(θ), in the sample reference system, (x′, y, z′), using: P′ is a rotation matrix that describes the spin rotation of a neutron caused by the magnetic field it travels through.
The matrix exponential map takes skew symmetric matrices (the Lie algebra so (3)) to rotation matrices (the Lie group SO (3)). Geometrically for rotation matrix P′ about an axis k k k k ( , , ) (with Tr denoting the matrix trace) and Working in the logarithmic chart has the advantage that (3) so forms a vector space and the linearized forward problem takes the simple form in each plane of constant y where L is the (rebinned) pixel size and  is the two-dimensional Radon transform in the plane, which can be inverted by the standard filtered back projection methods 15 to obtain a reconstruction of the measured magnetic field described by B = (B x , B y , B z ). A further limitation of the reconstruction algorithm is that it breaks down with phase wrapping, when the neutron spin precession angle gets larger than 180°. Since we normalise by wavelength we can average over all wavelengths where φ ≤ 180°.
The reconstruction procedure has been summed up in Fig. 2, where examples of measured sinograms are shown. Figure 3(b-d) show the x, y, and z components of the reconstructed 3D magnetic field from the measured solenoid. It can be seen that the strongest field region along the solenoid axis is easily reconstructed as well as the weaker magnetic field areas where it "wraps around" the ends of the solenoid. Figure 3(e-h) show selected slices (highlighted in (a)-(d)), where even the field from the current in the wires going to and from the solenoid is reconstructed as seen in (e).

Results
In order to compare our reconstruction to the expected resulting field from a current carrying solenoid, a calculation of the 3D magnetic field was done by dividing a description of the solenoid into 0.1 mm long straight wire segments and calculate the field contribution from each segment in a point cloud surrounding the solenoid using the Bio-Savart law: where ρ is the point where the field is calculated, ρ′ is the vector from the wire segment, dl, to ρ, and μ 0 = 4π × 10 −7 NA −2 is the magnetic constant. For further comparison a simple calculation using Ampère's Law (B = μ 0 IN/L) was done as well. The results are shown in Fig. 4, where the magnetic field strength within the central part of the solenoid is shown. As expected the simple calculation overestimates the magnetic field strength. The same is true to a much smaller extend for the Bio-Savart calculation, though the small discrepancy between this and the reconstructed magnetic field from measurements can be attributed to imperfections in the measured solenoid. Also small imperfections in the instrumental set-up are not taken into account, as they are negligible compared to the polariser, analyser, and π spin flipper efficiencies 9 . In the Supplementary material, Fig. A2, a further comparison between measurements and calculation is shown as well as curves based on ray-tracing simulations.

Discussion
We have with our successful measurements demonstrated the capabilities of a powerful technique for measuring three dimensional magnetic fields using ToF 3DPNT. Using the neutron spin precession in a magnetic field as a probe in combination with complicated reconstruction algorithms extracting information from the recorded data output, our proof-of-principle results are in a good agreement with calculations and serves as an initial demonstration of a novel technique that can extract hitherto unattainable information non-destructively from bulk samples.
In our FR correction we only corrected for the non-perfect polarisation characteristics of the polariser, analyser and π spin flipper. It should be noted that in order to take into account the comparatively much smaller 9 depolarisation in the spin rotators, further open beam measurements could have been performed at the expense of longer measurement time, or instead of using open beam measurements a polar decomposition by scaled Newton iteration could have been used to correct for small depolarisation effects thereby reorthonormalising P′.
The current limitation of our technique is that the reconstruction relies on the assumption that we have linearised around B = 0 and that we have neutron spin precession angles below the phase wrapping limit (φ ≤ 180°). If the phase wrapping limit lies within the measured wavelength band, it can be easily identified by following the progression of φ as the wavelengths increase, however, if the probed magnetic field is of such a strength that there is phase wrapping for even the fastest neutrons, the reconstruction algorithm would have to be expanded to possibly utilise the information contained in the period of φ as a function of wavelength 6,16 . Furthermore, the wavelength band used can be adjusted to stay within the limits of the assumption of linearisation around B = 0. To fully get around the assumption, an iterative reconstruction technique 13,17 with a forward model to approximate the measured field can be considered, as well as using vector field tomographic reconstruction 18 on φk directly (instead of using the assumption of linearisation around B = 0 to break it down to three scalars).
Outlook. The unique information only obtainable with our novel method can be of use in a broad range of fields such as electrical engineering, superconductivity, energy materials, thermoelectrics, etc.
Combining the three dimensional magnetic field information with techniques providing structural information, such as conventional attenuation contrast imaging and more advanced methods 19 like ToF three dimensional neutron diffraction (ToF 3DND) 20 , which can be performed using the same recorded data, 3DPNT provides a straight forward method for investigating the interplay between structural and magnetic sample composition. Data Availability. Data can be obtained from the authors by contacting Morten Sales (msales@fysik.dtu.dk). respectively shows the B x , B y , and B z components of the reconstructed magnetic field for field strengths above the threshold indicated by the colorbars. (e-h) shows different slices through the reconstructed 3D magnetic field volume, with the slice location indicated by the gray planes in the above figures. Note that, as the solenoid was not oriented with its axis exactly along the x-axis, but rather at an ~9° angle in the xz-plane, a central signal can be seen along z for the B z field in (h). (Supplementary animation).

Figure 4.
Comparison between simple calculated magnetic field strength inside solenoid, a more precise calculation using the Biot-Savart law, and the reconstructed magnetic field strength using polarimetric neutron tomography. The field strength shown, as a function of position along the solenoid axis, is for the central cylindric area around the solenoid axis. The mean of the field as well as the standard deviation within this area is shown for the Biot-Savart calculation as well as the reconstructed field averaged over the measured wavelength range. Curves for the reconstructed field using three single wavelengths are included as well.