Tomographic reconstruction of a three-dimensional magnetization vector field

Using x-ray magnetic nanotomography the internal magnetization structure within extended samples can be determined with high spatial resolution and element specificity, without the need for assumptions or prior knowledge of the magnetic properties of a sample. Here we present the details of a new algorithm for the reconstruction of a three-dimensional magnetization vector field, discussing both the mathematical description of the problem, and details of the gradient-based iterative reconstruction routine. To test the accuracy of the algorithm the method is demonstrated for a complex simulated magnetization configuration obtained from micromagnetic simulations. The reconstruction of the complex three-dimensional magnetic nanostructure, including the surroundings of magnetic singularities (or Bloch points), exhibits an excellent qualitative and quantitative agreement with the simulated magnetic structure. This method provides a robust route for the reconstruction of internal three-dimensional magnetization structures obtained from x-ray magnetic tomographic datasets, which can be acquired with either hard or soft x-rays, and can be applied to a wide variety of three-dimensional magnetic systems.


Introduction
Magnetic materials play an important role in modern technological and engineering applications, from data storage and sensors to electric motors and energy harvesting. In each case, their function is closely related to the details of the magnetic structure. Soft magnets for sensors can be more efficient, for example, if they are singledomain, while electric cores rely on the presence of domain walls and their displacement. In permanent magnets, the influence of grain boundaries on the local magnetization directly affects the coercivity of the magnet. In extended, bulk-like systems the magnetostatic interaction can often lead to the formation of complex three-dimensional magnetic structures.
In order to characterise the behaviour of these magnetic systems, and increase their efficiency, the internal magnetic structure needs to be determined. However, these bulk magnets have been historically difficult to investigate, relying mainly on indirect probing techniques, and it is only recently that tomographic techniques for imaging of three-dimensional magnetization configurations have been developed. In first demonstrations of magnetic tomography, neutron imaging was used to visualise the magnetic fields [1] and magnetic domain walls [2] within bulk magnetic systems with spatial resolutions of tens to hundreds of micrometres. Higher spatial resolutions have recently been achieved with electron [3,4] and soft x-ray [5,6] tomography. Imaging of magnetization structures with these methods is, however, limited to samples with thicknesses below 200 nm. This is particularly useful for the investigation of relatively thin three-dimensional magnetic structures [7,8], which may support complex magnetic configurations [9], exhibit magnetochiral effects [10,11], enable high domain wall velocities [12,13] and lead to asymmetric spin wave dispersion [14]. Recently we have developed an x-ray based tomography method with high spatial resolution and large penetration depth to probe bulk samples that can also be adapted to investigate thin films, making possible the investigation of a wide variety of extended magnetic systems [9].
One of the main challenges in the development of magnetic tomography has been the development of an appropriate reconstruction algorithm to recover all three components of the magnetization. Previously, the reconstruction of the three-dimensional magnetic structure has been achieved through additional constraints, such as the incorporation of prior information [5,15], or by analysing the angular dependence of the magnetic signal [6].
Here we present an iterative technique for the tomographic reconstruction of three-dimensional magnetization structures measured using x-ray magnetic circular dichroism (XMCD) that does not require prior information about the magnetic properties of the sample and can be applied to the reconstruction of complex magnetic structures found in extended magnetic systems. This technique is a significant improvement on the reconstruction method that was used in our experimental demonstration in [9], where the complex internal magnetic configuration of a m 5 m diameter pillar containing a number of fundamental magnetic structures was recovered with a spatial resolution of 100 nm. In that work, the three-dimensional magnetization was reconstructed in two steps. First the magnetization in two planes perpendicular to the axis of rotation was recovered, one with the sample untilted and the other with the sample tilted, and then the three-dimensional magnetization was determined from the two components of the magnetization in each plane by solving a set of simultaneous equations. The reconstruction algorithm used in [9] is described in detail in [16] and a schematic overview of the reconstruction algorithm is given in figure 1(a). We now present a more versatile reconstruction algorithm, schematically shown in figure 1(b), where all tomographic projections are combined to recover the three components of the magnetization in a single iterative reconstruction routine. Not only is this method more streamlined and generalised to an arbitrary measuring geometry but, by simultaneously combining data from multiple sample orientations, we circumvent known limitations in the reconstruction of two-dimensional divergent magnetization vector fields when using single rotation axis measurements. In this way we approach the correct solution with a more direct route, which results in significantly lower errors in the reconstructed magnetization.
In this article, we describe the single-step magnetic reconstruction algorithm that we have used to reconstruct the magnetization within a simulated mesoscopic pillar that contains a number of fundamental magnetic structures, which are presented in section 5. In addition, a comparison with the previous 2-step reconstruction algorithm used in [9] is given in section 5.3.

XMCD projections
In magnetic tomography, it is necessary to take into account the angular dependence of the magnetic signal. To do this, we first consider the nature of the XMCD, which is used to probe the magnetization of a specific element  [9], first the magnetization in the plane perpendicular to the rotation axis is reconstructed using an iterative algorithm for two measurement geometries. In a second step the two datasets are combined as described in [16] to obtain the three-dimensional magnetic structure. (b) For the single-step reconstruction, which is described in detail in this article, the magnetic projections measured for different sample orientations are combined to reconstruct the three-dimensional magnetization in a single step.
by tuning the energy of the x-rays to the resonant absorption edge of the element. For an x-ray beam propagating along theẑ direction, the change in the scattering factor due to the XMCD signal is proportional to the component of the magnetization parallel to the direction of propagation of the x-ray beam,ẑ: where f c and ( ) f m 1 are the electron density and XMCD scattering factors, respectively, m is the magnetization vector, = ( ) /M m M s , which is a function of the Cartesian coordinate vector ¢ = ¢ ¢ ¢ ( ) x y z r , , , where ¢ ¢ ¢ x y z , , are the object coordinates, and the two scattering channels corresponding to two circular polarisations C R -to-C R and C L -to-C L result in either a positive or negative contribution of the magnetic signal to the overall scattering factor [17][18][19]. This implies that, when rotating the sample about a tomographic axis by an angle θ, with the axis of rotation perpendicular to the direction of x-ray propagation, one probes only the magnetization in the plane perpendicular to the axis of rotation. In our dual-axis experimental setup, such tomographic measurements are performed for the sample at different orientations, tilted by an angle f with respect to the rotation axis, in order to probe the components of the magnetization in multiple planes, and thus all components of the magnetization vector, as depicted schematically in figure 2. The experimental procedure, however, can be extended to measure projections at arbitrary sample orientations, for example using an Euler cradle. To generalise our discussion to arbitrary orientations of the sample, we therefore consider the magnetic signal as the sample is rotated using a rotation matrix ( ) R n for the nth tomographic projection. The nth XMCD projection, P n (x, y), which for the case of absorption measurements is the projection of the imaginary part of the refractive index, β, can be expressed in terms of the magnetization scattering factor f as: where (x, y, z) are the laboratory coordinates, and ¢ ¢ ¢ ( ) x y z , , the object coordinates, respectively, and the integral with respect to dz represents the projection of the structure along the direction of the x-rays. Here r e is the classical electron radius, λ is the x-ray wavelength, n k at is the atomic density of the kth element, and we sum over all chemical elements present. Once we take into account that the x-rays are tuned to the absorption edge of a specific magnetic element, the equation becomes: where the magnetic contribution ( ) f m 1 originates only from the resonant element, whilst the electronic contribution of the magnetic element is absorbed into the sum over all elements in the first term, n at mag is the atomic density of the magnetic element, and ( ) R n and ( ) † R n are the rotation matrices defining the sample orientation, and its adjoint, respectively. Here, the magnetic scattering is dependent on the rotated magnetization vector represents the sample coordinates ¢ r , that are rotated from the laboratory coordinates by ( ) R n . The measured amplitude A n (x, y) is related to the XMCD projections P n (x, y) by: For coherent diffractive imaging techniques including ptychography, A n (x, y) is the measured quantity, as opposed to other methods such as full field transmission x-ray microscopy, for example, where | ( )| A x y , n 2 is measured. During the reconstruction, an estimated dataset of the projections, represented byP, is calculated from the reconstructed components of the magnetization, which is given bym. In this case, for discrete computation, we approximate the integrals by sums, and the estimated projections are then given by: whereô is the estimated non-magnetic structure, and c is a constant that relates the XMCD signal to the magnetization, which is given by: Here, for the sake of simplicity, we set the constant c=1 and take this factor as an overall scaling value in the reconstruction.

Requirements for the reconstruction of a vector field
Usually in tomography one reconstructs a scalar value for each voxel within the object, meaning that the problem is well posed due to the matching number of unknowns and equations, and that one scalar value for each pixel, for example the density or refractive index, is reconstructed from one tomographic dataset. For a three-dimensional vector field, however, rather than a single scalar value, all three components of the magnetization need to be recovered, giving rise to particular requirements in the measured data or on the constraints that must be applied during the reconstruction in order for the problem to be well posed. Following the discussion by Norton in the 1980s of the reconstruction of 2D vector fields [20], in which the reconstruction of a divergence-less field was shown to be unique, a more generalised discussion of the tomographic reconstruction of the components of arbitrary vector fields was given by Prince in the 1990s [21,22]. He demonstrated that, for the reconstruction of an arbitrary n-dimensional vector field, n tomographic projection datasets in which the probe is sensitive to n different directions of the vector field are required. This implies that three tomographic datasets are sufficient for the reconstruction of a three-dimensional vector field, such as the magnetization.
Measuring tomographic datasets that are sensitive to different components of the magnetization are, however, difficult to realise experimentally. In a tomographic measurement, the probe is normally sensitive to a particular component, or components, of the vector field with respect to the direction of propagation of the probe. Indeed, Braun and Heuck considered the scenario of vector tomography with different types of probe: a 'longitudinal' probe and a 'transverse' probe, which are sensitive to the components of the vector field parallel and perpendicular to the direction of propagation of the radiation, respectively [23]. They made use of the fact that a vector field such as the magnetization can be expressed as a combination of a source-free (solenoidal) component q S and a curl-free (irrotational) component q I : . For a vector field that is either irrotational or solenoidal, the requirements for a complete reconstruction are reduced, and fewer datasets are required [21,22]. Braun and Hauck [23] found that longitudinal and transverse probe measurements taken around a single axis of rotation can be used to reconstruct the solenoidal and irrotational components of a vector field, respectively, and they demonstrated the successful reconstruction of the solenoidal case of vortices within a liquid flow [23]. While knowledge of the properties of the vector field provides an opportunity to reduce the complexity of the problem, this work made clear the limitation of vector field tomography with a single type of probe.
For magnetic samples, XMCD measurements are sensitive to the component of magnetization parallel to the x-ray beam and therefore correspond to the longitudinal case, while measurements of x-ray magnetic linear dichroism (XMLD) and electron microscopy measurements are sensitive to the components of the magnetization and the magnetic field, respectively, perpendicular to the probe beam [24], and are therefore related to the transverse case. An experimental setup in which the XMCD signal is probed [9,19] in principle allows us to reconstruct a solenoidal magnetization vector field in the rotation plane by using the closed-form reconstruction approach in [23]. To reconstruct the full vector field, including the irrotational component, a different set or type of projections is required, or additional a priori information is needed. This would require an independent dataset, which could be obtained with tomographic measurements with two tilt orientations of the sample with respect to the axis of rotation. Indeed, measuring the sample response at diverse orientations has already been shown to provide enough information for higher dimensional reconstructions [25,26].

Reconstruction of the three-dimensional magnetization
The reconstruction of the magnetization is based on gradient-based optimisation. First, an estimate of the object is made, and projections of the current estimate object are calculated according to equation (5), thus creating an estimated projection dataset,ˆ( ) P x y , n . The calculated estimate projections are compared with the measured data, P n (x, y), through an error metric, which quantifies the difference between the measured and estimated data. Next, the current estimate of the object is updated in such a way that the error metric is reduced, directed by a gradient that computes how the object in each pixel should be changed to reduce the error metric. An initial estimate of the object can consist of an empty matrix, or an estimate stemming from prior information such as a first reconstruction obtained with a filtered back projection method, or a simulation. The process is then iterated, updating the current estimate of the object at each iteration, until a pre-specified convergence criterion, or a maximum number of iterations, is reached.
The error metric ò, is defined as: The gradient of the error metric is calculated from an analytical expression, rather than finite differences, for the sake of computational efficiency. For an arbitrary variable, α, the gradient is therefore given by: Using equation (5), the gradients of the error metric with respect to each variable are given by:

Numerical simulations
To validate the effectiveness of the reconstruction algorithm presented in this article, we performed numerical simulations of x-ray magnetic tomography using a model of a complex magnetic structure within a mesoscopic GdCo 2 pillar that was calculated using micromagnetic simulations [27].
The micromagnetic simulations were performed for a GdCo 2 cylindrical pillar of diameter m 1 m and a length of m 2 m, with a spatial resolution corresponding to a minimum feature size of approximately 6 nm. The magnetic configuration was obtained after relaxation from a state in which the magnetization was saturated perpendicular to the long axis of the pillar and consists of a complex magnetic structure similar to the one experimentally observed in our first demonstration of this technique [9] consisting of two magnetic domains along the height of the pillar, separated by a domain wall, as shown in figure 3. Vortices are equally present in the top and bottom regions of the pillar and, owing to inversion symmetry, one possesses a clockwise circulation, while the other has a counterclockwise circulation [28]. The core of each one of these vortices intersects the domain wall, giving rise to a Bloch point/anti-Bloch point pair. Such a complex magnetic structure provides a realistic challenge for the reconstruction algorithm, testing both how well a complex internal structure within an extended sample can be recovered, as well as whether the magnetization structure surrounding inherently divergent structures such as Bloch points can be correctly reconstructed.
Dual-axis magnetic tomography was simulated with similar parameters to our experimental demonstration [9]. In particular, single circular left polarisation projections with both electronic and magnetic contrast were calculated, as detailed in section 2, for tomographic measurements with two tilt orientations of the sample with respect to the axis of rotation, as shown schematically in figure 2. For each tilt axis, single circular polarisation projections were calculated over  360 with an angular spacing of  2 . For this geometry, the rotation matrix ( ) R n is defined as: where ( ) R z n tilts the sample away from the rotation axis, setting the sample orientation, while ( ) R y n rotates the sample about the tomographic rotation axis.
The three components of the magnetization were then reconstructed using 250 iterations of the gradientbased optimisation routine detailed in section 4. During the reconstruction, the magnetization was constrained to the location of the magnetic material using a mask that corresponds to the boundaries of the simulated sample. In experiments, the scalar electron density tomogram can be used to generate a mask that defines the location of the magnetic material. No constraint on the absolute value of the magnetization was applied. We note that from this point on, when referring to the magnetization, for simplicity we consider the object in an untilted orientation, and let the object coordinates be (x, y, z).
To evaluate the magnetic reconstruction, we compare the original micromagnetic simulation to the reconstructed magnetization structure. In particular, we identify specific complex magnetic textures that are present in the simulated structure, and determine the ability of the magnetic reconstruction to correctly recover both the nanoscale features in the magnetic structure, as well as their position within the pillar. An overview of the magnetic structure within the pillar is given in figure 3(a).
The magnetization within a horizontal slice of the pillar is shown in figure 3(b), where the direction of the magnetization is represented by arrows, and their colour indicates theŷ component of the magnetization. In the slice shown, the magnetization forms a vortex and, in addition, from the colour of the arrows we can identify two domains where theŷ component of the magnetization points in opposite directions, indicated by the red and blue regions in figure 3(b). The reconstruction of the magnetization distribution ( figure 3(c)) can be directly compared with the magnetic state within the slice in figure 3(b). One can see that not only are both the vortex structure and the m y domains correctly recovered, but an accurate reconstruction of the position of the center of the vortex as well as of the location of the domain wall separating the two domains alongŷ is achieved.
We now consider the validity of the reconstruction of smaller features exhibiting a strong divergence of the magnetization, such as Bloch points. These pose a challenge since the magnetization is maximally divergent around the singularities, within a radius of the order of the exchange length of the magnetic material, which is approximately 5 nm for GdCo 2 . While micromagnetic simulations cannot accurately describe the abrupt change in the magnetization at such a singularity, the magnetization can nevertheless be adequately described on a sphere surrounding the Bloch point, especially when its radius is greater than the exchange length, resulting in a region with a non-zero divergence of the magnetization around each Bloch point. Given that the spatial resolution currently available with magnetic tomography at light sources such as the Swiss Light Source is currently above the ferromagnetic exchange length of the material, only the somewhat smoother magnetization distribution surrounding the Bloch point at a larger radius can currently be imaged experimentally.
The magnetic structure surrounding a Bloch point at a radius of 50 nm is shown in figure 3(d). We have identified this structure as a circulating Bloch point [29]. The reconstructed magnetic structure surrounding the Bloch point shown in figure 3(e) exhibits an excellent agreement with the original structure ( figure 3(d)).
In addition, the reconstruction of the structure surrounding an anti-Bloch point [29] given in figure 3(g) is compared with its magnetic structure in figure 3(f). The twisted antivortex-like structure of the anti-Bloch point is recognisable in the reconstruction in 3(g), and also exhibits a very good qualitative and quantitative agreement with the micromagnetic simulations.
As well as obtaining a correct reconstruction of the magnetic structures present within the sample, one can also obtain an estimate of the precision with which one can locate them. We first map the centre of the domain wall separating domains of ±m y by plotting the isosurface of m y =0 for the simulations in figure 3(h) and for the reconstructed magnetization in 3(i). One can see that, apart from some slight artefacts at the top of the structure, the magnetic domain wall is correctly reconstructed with high precision throughout the pillar.
We next map the central region of the vortices. In thin films, within the core of the vortex, the magnetization rotates out of the plane of the sample, meaning that the core can be mapped by plotting an isosurface representing a small in-plane magnetization, such as, for example, . This isosurface is shown in figure 3(j) where one can see that this description also holds for the two vortices present in the system, despite their more complex, three-dimensional character. At the middle height of the pillar the isosurface expands to form two horizontal surfaces. These correspond to a transition region in which the magnetization is almost completely oriented along theŷ axis, as shown in figure 4, and thus are not trivially related to the location of the vortex cores. When the isosurface is plotted for the reconstructed magnetic structure in figure 3(k), one can see a very good agreement between the simulation and the reconstruction of the position of the defined vortex core regions (compare figure 3(j) and (k)). A quantitative estimate of the error in the core position is obtained by using the minimum of the value + m m x z 2 2 within an axial plane. The error in the vortex core position of the reconstruction, δ core , is therefore calculated as: where x sim (x rec ) and z sim (z rec ) are the x and z coordinates of the vortex cores for the original simulation (reconstruction). For the reconstruction of the magnetization, the error was found to range between 0 and 1 voxel, with an RMS error of 0.15 voxels. With the presented reconstruction algorithm and sufficient signal-tonoise ratio, one can therefore expect a precision of less than 1 voxel, which in this case corresponds to approximately 6 nm, in locating specific structures within the sample.

3D reconstruction of the solenoidal and irrotational part of the magnetization
The dual-axis measurement approach presented in [9] and used here provides additional information that helps counter the expected limitations of single-axis tomographic measurements, in particular in the reconstruction of the irrotational component of the magnetization. In [16] (chapter 5, section 5.2.2), we showed that, for the successful reconstruction of the irrotational part of 2D magnetization structures with a single-axis measurement, additional constraints, such as limiting the absolute value of the magnetization, were required.
Here we evaluate how effective the recovery of the irrotational component of the micromagnetic structure is via the dual-axis tomography measurement. In figure 5(a), the divergence of the magnetization within an example horizontal slice of the micromagnetic simulations is given, which can be compared with the divergence of the xz components of the magnetization of the single rotation axis reconstruction in figure 5(b) and the dual rotation axis reconstruction presented in this article in figure 5(c).
In the single axis reconstruction, shown in figure 5(b), the reconstruction has no internal divergence, which is consistent with the findings of Braun and Hauck [23]. In contrast, in the dual-axis reconstruction in figure 5(c) a non-zero divergence in the x-z plane is recovered, which exhibits a good agreement with the in-plane divergence of the micromagnetic simulation, recovering the divergence with a standard error of D  = -( · ) m 20% x z . We can therefore conclude that the measurement about the second rotation axis provides additional information that enables the reconstruction of the irrotational component of the magnetization in the plane. We perform the same analysis for the internal divergence of the vector field in two other planes, the x-y and y-z planes, in figures 5(d)-(g), where the irrotational component of the in-plane magnetization is accurately reconstructed, with a calculated error of D  = -( · ) m 10% x y and D  = -( · ) m 15% y z , respectively. We note that, with magnetic tomography using XMCD contrast, it is possible to recover the two-dimensional irrotational component of the magnetization with a dual-axis measurement. However, a reconstruction of the three-dimensional irrotational component of the vector field,  · m, has only been demonstrated using additional types of probes, such as a transverse probe, [23] or by using additional constraints such as on the absolute value of the magnetization [16].
For completeness, we also consider the reconstruction of the solenoidal part of the magnetization vector field. In particular, the component of the curl of the magnetization in different perpendicular planes is plotted in figure 6. We see that, as predicted for tomographic measurements with a longitudinal probe [23], all three components of the curl of the magnetization are recovered, with an error of D ´=  ( ) m 9% x , Δ(∇×m) y =±5% and Δ(∇×m) z =±15%.

Quantitative analysis of the error of the dual-axis magnetic reconstruction
To confirm the validity of the reconstruction of the magnetization, a quantative analysis of the error in the reconstruction is performed. In particular, the errors in the angle of the reconstructed magnetization and the magnitude of the reconstructed magnetization with respect to the original micromagnetic simulations are calculated. The errors are shown in a bivariate histogram in figure 7, where one can see that the error is in general limited to small angles with 95% of the voxels having an error of less than 1% of the magnitude of the magnetization, and less than 1°error in the reconstructed angle of the magnetization in three dimensions.
To visualise the spatial distribution of the errors, we plot the error in the angle and the magnitude of the single-step reconstructed magnetization for two vertical slices through the pillar that contain the anti-Bloch point, the Bloch point as well as the core of the vortex in figure 8. Within the volume of the pillar we observe an Figure 5. Comparison of the 2D divergence of the simulated and the reconstructed magnetization in different perpendicular planes, schematically illustrated in the insets. The x-z divergence of the magnetization is given for a horizontal slice of the pillar for (a) the micromagnetic simulations, (b) a single axis measurement and (c) the dual-axis measurement. As found in the 2D single-axis tomographic simulations presented in [16], for a single rotation axis, no internal divergence is recovered, and only the magnetic charges at the surface are reconstructed. With the dual-axis measurement, however, the divergence in the x-zplane is recovered. A good agreement is observed between the simulated and the reconstructed divergence for the dual-axis reconstruction for both the x-y plane (d), (e) and the y-z plane (f), (g) for both the structure and the magnitude of the divergence.   increase of the error in the reconstruction of the magnetization in the vicinity of the Bloch points. In particular, the error in the direction of the magnetization increases to approximately   5 and   10 in the close vicinity of the Bloch point ( figure 8(d)) and the anti-Bloch point ( figure 8(b)), respectively. Indeed, the error in both the magnitude and direction of the magnetization is greater for the anti-Bloch point, as can be seen in figures 8(a), (b), than for the Bloch point or the vortex core, as seen in figures 8(c), (d). This increase in the error of the reconstruction is due to the fact that the magnetization structure surrounding the anti-Bloch point is more inhomogeneous over a larger distance than for a Bloch point, and is therefore more divergent. This relationship between the error and divergence of the magnetization reflects the challenges of reconstructing divergent structures, such as the Bloch points, as discussed in section 3. Nevertheless, the fact that the magnetization is reconstructed well to within 10°, even in the close vicinity of such strongly inhomogeneous structures, gives us confidence in the robustness of our magnetic reconstruction algorithm.
So far in this study we have simulated dual-axis tomography, as we demonstrated experimentally in [9]. However, in principle there is no reason to limit ourselves to a dual-axis experimental setup, and it is has been suggested that more than two independent axes of rotation may provide an advantage for the reconstruction of a three-dimensional vector field [22]. In fact, since we have demonstrated that the second orientation provides additional non-redundant information, it is interesting to determine what further improvement can be expected from more than two tilt orientations. To investigate the optimal imaging conditions, we have therefore simulated magnetic tomography as above for up to 5 axes of rotation, keeping the total number of projections constant at 360. This is achieved by performing tomography with a single rotation axis for multiple tilt orientations of the sample f (see figure 2). Reconstructions are performed with 250 iterations of the gradientbased iterative reconstruction algorithm presented in section 4 in this paper, and the different combinations of parameters are summarised in table 1. We note that for these simulations, purely magnetic projections "XMCD" were measured over  180 , which is equivalent to measuring single circular polarisation projections over  360 .
To determine the influence of the number of tilt axes on the quality and validity of the reconstruction, we first consider the evolution of the error metric during the 250 iterations of the reconstruction in figure 9(a), which is plotted for 2-5 tilt axes. We see immediately that the error metric continuously decreases with the number of iterations for all scenarios. However, the rate at which it decreases depends greatly on the number of tilt axes over which the projections are distributed. In particular, as the number of tilt axes is increased from 2 to Table 1. Parameters for magnetic tomography with multiple tilt axes. The tilt angle f refers to tilt orientations of the sample for which tomography is performed, as defined in equation (11) and shown schematically in figure 2. Note that the total number of projections is kept constant.

No.
No. projections Angular Tilt tilt axes per axis spacing angle (f)  3 axes, a significant increase in the rate of convergence is obtained. On increasing the number of tilt axes to 4, the convergence is further improved. However a further increase to 5 tilt axes does not lead to a significant further improvement of the convergence.
To determine the accuracy of the reconstructions, we calculate the error in the direction and magnitude of the reconstructed magnetization for each number of tilt axes after 250 iterations of the single-step reconstruction algorithm. The errors are plotted as a function of the number of tilt axes in figure 9(b) and one can see that, with an increasing number of tilt axes, the error in both the magnitude and direction of the magnetization decreases. Again, an increase from 2 to 3 and 4 tilt axes provides significant improvements in the reconstruction, while a further increase to 5 tilt axes does not result in significant improvements. We therefore conclude that having the projections equally distributed over 3 or 4 axes of rotation is optimal for x-ray magnetic tomography.

Comparison between the single-step and 2-step algorithms
The single-step reconstruction algorithm presented in this article represents a more direct route to the reconstruction of the three-dimensional magnetization when compared to the algorithm that we have used in the first experimental demonstration of the technique in [9]. To determine the effectiveness of the previous 2-step reconstruction algorithm, in addition to reconstructing the simulated tomographic data with the singlestep algorithm presented in this article, we used the 2-step algorithm to reconstruct the magnetization configuration, and therefore to obtain a comparison between the two algorithms.
In the 2-step reconstruction method, schematically shown in figure 1(a), first the magnetization in the plane perpendicular to the rotation axis is reconstructed for each of the sample orientations. In a second step, the three components of the magnetization are recovered from the various planes using an iterative reconstruction routine that is described in detail in [16]. Here we reconstruct the magnetization using 50 iterations in the first step, and 50 iterations in the second step. More iterations were not necessary as the reconstruction normally converged within this number of iterations.
We first evaluate the 2-step magnetic reconstruction with a qualitative comparison between the simulated and reconstructed magnetic structures of the features shown in figure 3. The magnetic structure within the slice in figure 10(a) is directly compared with the corresponding magnetic structure reconstructed using the 2-step reconstruction algorithm implemented in [9] in figure 10(b), and with the single-step reconstruction algorithm presented in this article, in figure 10(c). While it is clear that the in-plane vortex structure and the out-of-plane two-domain state (red and blue coloured regions) are adequately reconstructed in both cases, the single-step algorithm leads to a more accurate reconstruction of the position of the center of the vortex as well as of the location of the domain wall separating the two domains alongŷ .
In addition to magnetic domains, we can consider the validity of the reconstruction of structures such as Bloch points. The magnetic structure surrounding a circulating Bloch point in the micromagnetic simulations is shown in figure 10(d), which can be directly compared with the magnetic structure of the Bloch point reconstructed with the 2-step algorithm shown in figure 10(e). The reconstruction exhibits a good agreement with the original structure. However, the magnetic structure obtained with the single-step reconstruction algorithm, shown in figure 10(f), is once again noticeably more accurate. In addition, the magnetic structure of the anti-Bloch point in figure 10(g) is compared with the 2-step and single-step reconstructions in figure 10(h) and 10(i), respectively. One can see that the twisted structure of the anti-Bloch point, recovered with the 2-step reconstruction algorithm, is recognisable in the reconstruction in 10(h) but is noticeably distorted. The deviation in the reconstructed structure from the true magnetic structure is mainly attributed to the fact that the anti-Bloch point is located close to the edge of the pillar, which decreases the quality and the resolution of the reconstruction in the case of the 2-step algorithm. Indeed, during the rotation of the structure to tilt the sample in the second stage of the 2-step reconstruction algorithm, some artefacts are introduced at the edge of the pillar and, as a result, noise is locally introduced into the reconstruction. In contrast, with the single-step reconstruction algorithm we observe a very good agreement between the simulated and reconstructed magnetic structure of the anti-Bloch point (compare figure 10(g) and 10(i)), showing not only that the new reconstruction algorithm offers an improvement in the quantitative reconstruction of the magnetic components, but also that it appears to be more robust to edge artefacts compared to the 2-step algorithm.
The errors for the 2-step reconstruction are shown in a bivariate histogram in figure 11(a), where one can see that the error is in general limited to small angles with 95% of the voxels (of a total 1612221 in the pillar) having an error of less than 2.4% in the magnitude of the magnetization, and less than  15.3 error in the reconstructed angle of the magnetization in three dimensions. When we compare the error of the reconstruction with the 2-step algorithm in figure 11(a) to the error of the single-step reconstruction in figure 11(b), we see that, while the error in the magnitude of the magnetization is comparable, the error in the direction of the reconstructed magnetization is significantly reduced. In fact, 95% of the voxels have an error of less than 1°in the direction of  Figure 11. A bivariate histogram of the error in the angle θ m and the magnitude | | m of the reconstructed magnetization is given for (a) the 2-step reconstruction algorithm used in [9] and for (b) the new single-step reconstruction algorithm presented in this paper. the reconstructed magnetization, indicating a large improvement in the quantitative reconstruction of the orientation of the three-dimensional magnetization.
It is likely that the larger error in the 2-step reconstruction originates from the limitation in the single rotation axis reconstruction of the divergence of an in-plane vector field as shown in figure 5(b). In the 2-step routine, two planes are reconstructed independently using a single rotation axis reconstruction, and are subsequently combined to retrieve the three-dimensional magnetization. With the single-step reconstruction algorithm, the single axis reconstruction and its associated limitations are avoided, which leads to a more accurate solution. We note that, as the spatial resolution of the reconstruction is limited by the signal-to-noise ratio as described in [9], the higher accuracy of the new algorithm does not necessarily provide an improvement in the reconstruction for data acquired at third generation synchrotron facilities. Nevertheless, with the advent of diffraction-limited synchrotron sources and the subsequent increase in spatial resolution, this new algorithm will be very useful for magnetic tomography as well as alternative new experimental geometries.
Computationally, the single-step reconstruction is more time consuming than the 2-step algorithm due to the requirement to calculate projections at arbitrary directions, which requires a more computationally demanding interpolation. In particular, with a CPU-based implementation, the single-step reconstruction and the 2-step reconstruction required approximately 8 h and 5 h, respectively, for 50 iterations and a sample size of 1612221 voxels on an Intel Ceon E5-2690v3 with 24 cores, 2.6 GHz, and 256 GB of RAM. With a recent implementation of our algorithm using the parallel computing toolbox for graphics processing unit (GPU) calculations and a modified ASTRA toolkit [30][31][32], the time required for 50 iterations of the single-step reconstruction is reduced to approximately 8 min using a GPU card with NVIDIA Quadro K4200 and 4 GB of memory.

Conclusion
We have presented an iterative reconstruction algorithm for the tomographic reconstruction of a threedimensional magnetization vector field from tomographic data measured with XMCD. Unlike the 2-step approach presented in [9], this single-step algorithm simultaneously combines all available XMCD projections to obtain a three-dimensional map of the magnetization vector field. Numerical simulations of a dual-axis measurement show a significant improvement in the accuracy of the reconstruction compared to the earlier approach for a complex simulated three-dimensional magnetic configuration that contains a variety of fundamental structures including vortices, domain walls and micromagnetic singularities. Comparisons between the reconstructions and the simulated structure demonstrate that the magnetic configurations are in general successfully recovered with both reconstruction algorithms, even in the vicinity of strongly diverging magnetic textures, leading to the successful identification and reconstruction of the magnetization distribution around Bloch points. With the presented single-step reconstruction algorithm, we observe a significant improvement in the error of the reconstructed magnetization, and a quantitative analysis of the reconstruction reveals that for 95% of the voxels, the magnitude of the magnetization is successfully reconstructed to within 1%, and the direction of the three-dimensional magnetic vector to within 1°. This new technique offers a route to determining arbitrary three-dimensional magnetization structures with resonant x-ray measurements, and thus opens the door to improved three-dimensional magnetic investigations.
Additionally, with our numerical study we show that, although with our technique it is not possible to recover the in-plane divergence of a magnetic structure with a single axis of rotation, in agreement with the work of Braun and Hauck [23], probing additional orientations of the sample provides the non-redundant information needed for a successful reconstruction of both the solenoidal and irrotational in-plane components of the vector field. In particular, we validate the dual-axis experimental approach presented in [9].
Although x-ray ptychography was used in the experimental demonstration of magnetic nanotomography in [9], x-ray ptychography is not required for the reconstruction of the magnetization, and magnetic tomography can be used with a number of different x-ray imaging techniques. Furthermore, in combination with lower energy soft x-rays [33], where the magnetic signal is significantly higher, high spatial resolution reconstructions of the three-dimensional magnetization within thin films and magnetic nanostructures with a spatial resolution below 10 nm are now feasible.
Here we have presented a generic framework that allows reconstructions of the three-dimensional magnetization vector field based on XMCD projections from arbitrary orientations. Based on this study, a number of possible improvements and extensions for future magnetic tomography experiments and reconstruction techniques have become clear. For example, we have observed further improvements in the accuracy of the reconstruction when tomographic projections are measured for three or four different sample orientations with respect to the axis of rotation, rather than two. In addition, the reconstruction algorithm can be directly implemented with alternative measurement geometries such as laminography [34], a three-dimensional imaging technique in which the rotation axis is not perpendicular to the direction of propagation of the x-rays, and is therefore particularly suited for the study of flat sample geometries such as magnetic thin films. In addition, the gradients in equation (10a) can be used directly in Simultaneous Iterative Reconstruction Technique [35] or Simultaneous Algebraic Reconstruction Technique [36] iterative approaches as alternatives to the gradient descent optimisation. We also note that the presented reconstruction method is valid for the geometries proposed for soft x-ray tomographic measurements in [37].
Finally, the dual-axis tomographic technique could be extended to measurements of XMLD, which may allow for the determination of the three-dimensional magnetic structure of antiferromagnetic materials. With recent developments in the manipulation of antiferromagnets with electrical switching protocols [38], there is a growing interest in their technological applications since antiferromagnets are robust against external magnetic fields, and can be used for devices requiring switching at ultrafast timescales.