Phase retrieval for near-field X-ray imaging beyond linearisation or compact support

X-ray phase contrast imaging based on free space propagation relies on phase retrieval to obtain sharp images of micro- and nanoscale objects, with widespread applications in material science and biomedical research. For high resolution synchrotron experiments, phase retrieval is largely based on the single step reconstruction using the contrast transfer function approach (CTF), as introduced almost twenty years ago [Cloetens et al., Appl. Phys. Lett. 75, 2912 (1999)]. Notwithstanding its tremendous merits, this scheme makes stringent assumptions on the optical properties of the object, requiring, in particular, a weakly varying phase. In this work, we show how significant the loss in image quality becomes if these assumption are violated, and how phase retrieval can be easily improved by a simple scheme of alternating projections. Importantly, the approach demonstrated here uses the same input data and constraint sets as the conventional CTF-based phase retrieval, and is particularly well suited for the holographic regime.X-ray phase contrast imaging based on free space propagation relies on phase retrieval to obtain sharp images of micro- and nanoscale objects, with widespread applications in material science and biomedical research. For high resolution synchrotron experiments, phase retrieval is largely based on the single step reconstruction using the contrast transfer function approach (CTF), as introduced almost twenty years ago [Cloetens et al., Appl. Phys. Lett. 75, 2912 (1999)]. Notwithstanding its tremendous merits, this scheme makes stringent assumptions on the optical properties of the object, requiring, in particular, a weakly varying phase. In this work, we show how significant the loss in image quality becomes if these assumption are violated, and how phase retrieval can be easily improved by a simple scheme of alternating projections. Importantly, the approach demonstrated here uses the same input data and constraint sets as the conventional CTF-based phase retrieval, and is particularly well suited for the ...

Phase retrieval for near-field X-ray imaging beyond linearisation or compact support J. Hagemann X-ray phase contrast imaging based on free space propagation relies on phase retrieval to obtain sharp images of micro-and nanoscale objects, with widespread applications in material science and biomedical research. For high resolution synchrotron experiments, phase retrieval is largely based on the single step reconstruction using the contrast transfer function approach (CTF), as introduced almost twenty years ago [Cloetens et al., Appl. Phys. Lett. 75, 2912(1999]. Notwithstanding its tremendous merits, this scheme makes stringent assumptions on the optical properties of the object, requiring, in particular, a weakly varying phase. In this work, we show how significant the loss in image quality becomes if these assumption are violated, and how phase retrieval can be easily improved by a simple scheme of alternating projections. Importantly, the approach demonstrated here uses the same input data and constraint sets as the conventional CTF-based phase retrieval, and is particularly well suited for the holographic regime. V C 2018 Author(s Classic radiography based on absorption contrast has prevailed the first hundred years of radiography and continues to be the most important imaging modality in medical as well as in industrial x-ray imaging. This is despite the fact that x-ray phase contrast has been demonstrated already two decades ago and can be based either on optical elements such as crystals, gratings, and Zernike phase plates, or on free space propagation. 1-3 X-ray phase contrast or near-field holography (NFH) utilizes the self interference of the beam behind an object by free space propagation to turn minute phase shifts of the x-ray wavefront into a measureable intensity pattern. In this way, non-absorbing or weakly absorbing matter such as soft tissues can be visualized. 4 Contrast is particularly strong in the holographic regime, corresponding to small Fresnel numbers Fr ¼ Dx 2 =kz ( 1, where k is the wavelength and z is the propagation distance. Dx is a characteristic length in the object plane; for convenience of numerical propagation, it is defined as the effective pixel size of the reconstruction. A main limitation of NFH is phase retrieval, since simple holographic back-propagation as introduced by Gabor is associated with strong artifacts of the twin image. 5 In most cases, phase-contrast experiments are carried out in the direct-contrast regime Fr ! 1, where linearization along the propagation direction z yields an approximate inversion formula. 6,7 These approaches fail, however, in the holographic regime where the phase-sensitivity of NFH is highest. 8,9 In order to profit from this regime, an inversion based on the known contrast transfer function for an optically weak object has been proposed twenty years ago, 10 which is based on k ¼ 1…K different measurement planes with Fr k . This so-called CTF-phase retrieval has been extremely successful, due to its simple implementation as a deterministic and single step computation, and has found widespread use. However, CTF fails or has limited image quality in many cases, when the basic assumptions of weakly varying phase are not fulfilled. 11 Thus, variations of CTF phase retrieval have been introduced to partly overcome these issues. 12,13 Similar to far-field coherent diffractive imaging (CDI), 14 an alternative phase retrieval approach is iterative projection algorithms, [15][16][17] cycling between the measurement and the object plane, utilizing, for example, the powerful constraint of a compact object support. 18,19 However, many of the more relevant objects are not compactly supported, and hence not amenable to the support constraint. Furthermore, the large (and scalable) field of view (FOV) without scanning of the object as in ptychography 20 is a significant advantage of NFH, for example, for dynamic or tomography studies. Note that by the use of divergent wavefronts behind a focal plane (cone beam geometry), the FOV and magnification M ¼ z 01 =ðz 01 þ z 02 Þ can be easily selected by variation of the source-to-sample distance z 01 and/or equivalently the source-to-detector distance z 02 . We thus face a situation where we lack suitable phase retrieval methods exactly for that situation which would maximize the signal, i.e., in the holographic regime and for wavelength where the object properties cannot be linearized.
In this work, we show that a simple approach of iterative alternative projections (AP) can fill this gap, yielding superior image quality even for samples which do not obey the assumption of a weakly varying phase, which is particularly well suited for the holographic regime. The input data correspond to the same set of measurements which is typically used in the CTF phase retrieval with K ¼ 4. To demonstrate the AP approach for NFH data, we make use of a recently introduced x-ray waveguide (WG) setup which provides a highly coherent, quasi point-like secondary source for holographic illumination, 9,21,22 as sketched in Fig. 1 experiment was carried out at the G€ ottingen instrument for nanoimaging with x-rays (GINIX), installed at the P10 beamline of the PETRA III storage ring (DESY, Hamburg), at a photon energy of 8 keV. A x-ray waveguide was used for spatial and coherence filtering, 21,23 forming a quasi-point source of holographic illumination. In the present case, a 1 mm long channel in silicon fabricated by e-beam lithography 24 was used, with a cross section of 40 nm. The object consisted of a layer of polystyrene spheres with diameter D ¼ 15 lm and was placed at K ¼ 4 a different distance z 01 behind the waveguide, chosen as {1.5, 1.52, 1.6, 1.8} mm to create diversity in the near-field diffraction pattern, as in standard CTF phase retrieval. 1, 25 The holograms were recorded by a fiber coupled scintillator CMOS detector (Photonic Science) with N x Â N y ¼ 2048 Â 2048 pixels and a pixel size of 6.5 lm, placed at z 02 ¼ 5178 mm behind the waveguide exit, for the purpose of geometric magnification. Figure 2 shows an example for the measurements I k , after flat-field correction and alignment (correcting shifts and magnification). The dataset is available online. 26 Note that the computations were carried out in an equivalent parallel beam coordinate frame, using the concept of an effective propagation distance z eff arg1 ¼ z 01 z 12 =z 02 , according to the Fresnel scaling theorem. 4 The reconstruction of the dataset was carried out by two different schemes: (i) CTF and (ii) an iterative reconstruction scheme based on von Neumann's alternating projection algorithm (AP). 27 The schemes are introduced briefly in the following. The CTF phase retrieval is based on a linearization of the image formation also for smaller Fresnel numbers as long as the objects are "weak enough," i.e., as long as the object exhibits a slowly varying phase and weak absorption. 11 In this case, a deterministic phase-retrieval approach based on the known functional form of the CTF can be derived 10,25 for a homogeneous object. A least-square minimization of the CTF leads to 28,29 /ðrÞ ¼ F À1 with an additional frequency-dependent regularization aðk x;y Þ of high spatial frequencies and v ¼ p Fr ðk 2 x þ k 2 y Þ. k x ¼ 2n x =N x and k y ¼ 2n y =N y are spatial frequencies in Fourier space with n x;y 2 ½ÀN x;y =2 … N x;y =2. The ratio of phase shift to absorption is j ¼ d=b. The frequency cut-off of 0.027 for the reconstruction has been determined from the first maximum of the CTF function for the given Fr's.
The AP reconstruction is obtained by simple alternation of the projections P M and P S , which can be written as W nþ1 ¼ P S ðP M ðW n ÞÞ. The projector P M adapts the current iterate to the measured intensities I k D Fr is given by In the current case, the P M was applied in an averaging way. It is also possible to apply P M in a sequential way on the measurements or just select a subset of measurements per iteration. The constraints in real space, i.e., the constraints on the sample are applied by P S . The classical constraint for P S is the support constraint. By defining the area in which the specimen is located, the number of solutions is constrained. This constraint cannot be applied in the case of closely packed specimen or specimen filling the whole FOV. The test sample thus represents a class of extended objects, relevant for most applications. For this study, two variants of P S have been applied. The negative phase constraint P pp S is a range constraint on the phase. It also enforces the assumption of a nonabsorbing object P pp S has a similar effect as a support, i.e., non-object pixels are set to one which leads to a reduction of background variations. The amplitude range constraint P ar S is a relaxation of P pp S in the sense that it permits the complete range of phases. Notably, it allows phase wrapping and it restricts the amplitudes to values smaller than one Figure 3 shows the reconstruction obtained by the two methods. Panels (a) and (b) show the reconstructed phases of the specimen, for CTF and AP, respectively. The insets show a zoom of the region marked in (a) and (b). The artifacts of the CTF reconstruction (a) become clearly visible. The spheres show a cloudy like sub-structure, and the background shows a variation when comparing sphere free regions with regions close to the spheres. The AP reconstruction (b) shows nearly homogeneous spheres, as expected on a flat background. The range of phases exceeds the interval ½Àp; 0, so that phase unwrapping becomes necessary. The reconstruction has been obtained after 2500 iterations of AP with the change of P S , starting with 1000 iterations of P ar S , then 500 iterations with P pp S , and finally 1000 iterations with P ar S . This scheme reduced background artifacts and took the phase wrapping correctly into account. In the supplemental material, we have further analyzed the reconstruction quality and reviewed reasons for the unsatisfying CTF reconstruction: In spite of AP's simplicity, it recovers the best approximation of a sphere lying on a flat surface. Further, we have investigated the effect of reducing the number of measurement planes below four (standard case). We found that, for a single distance as well as for two distances (measurement planes), the reconstruction did not converge to satisfactory image quality.
The reconstructions were carried out on a Nvidia K40m, using MATLAB R2017a, supported by custom CUDA 30 -kernels. The run-time was 295 s, using the same code without further optimizations on a Nvidia P100 resulted in a runtime of 66 s.
Twenty years ago, when CTF phase retrieval was introduced to x-ray phase contrast imaging, numerical computation was in its infancy, and several hundreds to thousand iterations for an equal number of projections in a typical tomography dataset was prohibitive. It is also for this reason that CTF as a single-step phase reconstruction schemes became so successful. 19,31,32 The same is true for the regime of Fr ! 1, the edge enhancement regime, where single step linearized schemes based on the transport-ofintensity (TIE) formalism are equally successful. As we stress here, both cases either violate the respective approximations or use distances and photon energies which are suboptimal for weakly scattering objects. Now, with the increase in computation power, the AP scheme as proposed here is no longer an obstacle. The application of AP to a tomographic dataset consisting of 1000 projections is as well feasible. Figure 4 shows the reconstruction of an Epon embedded Golgi-Cox stained brain slice of the Hippocampus of a wildtype mouse. 33 The measurements have been carried out at the GINIX-experiment using a crossed waveguide with guiding layer thickness of 81 nm at a photon energy of 13.8 keV. The sample was placed at z 01 ¼ 145 mm and the detector (sCMOS, Hamamatsu, Tokyo, Japan, fibre-coupled 20 lm LuAG:Ce single crystal scintillator, pixel-size 6.5 lm) at z 02 ¼ 5040 mm resulting in an effective pixel size of 192 nm for the reconstruction. The acquisition time per projection was 0.4 s. In panel (a) and (b) the reconstructed phases for the same projection of CTF and AP are shown, respectively. The iterative reconstruction for each projection has been obtained after 2000 iterations using P ar S . The CTF reconstruction has been carried out with cut-off frequency of 0.035 and d=b ¼ 6. After the phase reconstruction of all projections the tomographic reconstruction has been carried out with the same corrections and parameters for both sets of reconstructions. All processing steps are listed in the supplemental material. Visual comparison of the x-z tomography orthoslices displayed in (c) and (d) shows that both reconstructions seem to have resolved the same features. Still the CTF exhibits stronger background variations. Overall the iterative reconstruction shows less artifacts on small length scales (reduced fringes). The y-z orthoslices in (e) and (f) show on the other hand the streak artifacts connecting strongly absorbing objects become more visible in the AP reconstruction. The CTF reconstruction has a smaller range of reconstructed phases (from the unmet homogeneity assumption). Contrarily, AP recovers the expected range and boundaries of the metal stained neurons. Given these sharp features the number of tomographic projections is not sufficient and gives rise to the streak artefacts.
In summary, we have demonstrated the failure of CTF phase retrieval even for a simple soft matter material specimen such as polystyrene spheres, where the associated artifacts can be clearly diagnosed, which is not always the case for biological tissues. At the same time, a simple iterative reconstruction based on AP and applied to the same data yields excellent reconstruction quality. With the advent of new computational hardware, especially general purpose graphical processing units, the computation time of iterative algorithms can be drastically reduced. Thus application to large data sets becomes feasible. A common drawback of iterative schemes is the large number of iterations needed to recover the correct phase shift, especially for small spatial frequencies. To address this issue and to improve convergence we have tested other variants of projection schemes, which are more elaborated than AP, such as the relaxed averaged alternating reflections (RAAR). 34 However, the reconstruction obtained by RAAR did not yield a better result at the same number of iterations. We have also tested an iteratively regularized Gauss Newton method 35 using the same set of constraints, which also yielded good reconstructions, but at same or higher cost of computations, as detailed in the supplementary material. The challenge for the future improvement thus lies in finding an optimal preconditioner for iterative schemes, i.e., an appropriate starting guess. A typical choice is the plain holographic reconstruction, i.e., a back propagation of the measured data, which is spoiled by the twin image. Other choices would obviously be the single step phase reconstructions as CTF or Paganin's method. For the data set presented here, the use of the CTF reconstruction did not result in better reconstructions or less iterations needed. The iterative scheme rather needed some iterations "to repair" the starting guess. While we have deliberately chosen a strong phase object with negligible absorption, the work can be generalized to the case of mixed phase and absorption object. CTF cannot treat objects, for which these two properties are uncoupled, unless the object consists of a single material. 1 Hence, CTF requires not only a slowly varying phase, but in addition either negligible absorption or a homogeneous material over the entire FOV. Contrarily, the AP scheme can easily cope with strong phase gradients (see boundaries of the spheres), and in addition it can be implemented without any constraints of the absorption. For the object, even if only the well warranted amplitude and phase range constraints are applied, cf. Eqs. (5) and (4), high image quality can be obtained. 36 See supplementary material for additional comparisons with other iterative algorithms.
We thank Markus Osterhoff and Micheal Sprung for the support during the experiment. We thank Jan Goeman for keeping the number crunching machines up and running. Financial support by CRC 755-Nanoscale Photonic Imaging is greatly appreciated.