Optimisation of phase imaging geometry

Certain phase retrieval methods use knowledge about the free space propagation of a wave to phase a paraxial beam passing through one or more measurement planes. This approach has been widely applied and has been shown to quantitatively retrieve the refractive index profile of a sample. The quality of the phase retrieval will depend on a range of factors including sample feature size, propagation distance, measurement plane separation, wavelength and noise. Here we describe an optimisation study for two-plane phase retrieval using a laboratory-based X-ray source that considers all of these factors. We discuss our results in the context of a three-dimensional reconstruction of a sample refractive index profile. ©2010 Optical Society of America OCIS codes: (340.0340) X-ray optics; (340.7440) X-ray imaging; (100.5070) Phase retrieval. References and links 1. D. Paganin, and K. A. Nugent, “Noninterferometric phase imaging with partially coherent light,” Phys. Rev. Lett. 80(12), 2586–2589 (1998). 2. P. Cloetens, W. Ludwig, J. Baruchel, D. Van Dyck, J. Van Landuyt, J. P. Guigay, and M. Schlenker, “Holotomography: Quantitative phase tomography with micrometer resolution using hard synchrotron radiation x-rays,” Appl. Phys. Lett. 75(19), 2912–2914 (1999). 3. J. N. Clark, G. J. Williams, H. M. Quiney, L. Whitehead, M. D. de Jonge, E. Hanssen, M. Altissimo, K. A. Nugent, and A. G. Peele, “Quantitative phase measurement in coherent diffraction imaging,” Opt. Express 16(5), 3342–3348 (2008). 4. W. Coene, G. Janssen, M. Op de Beeck, and D. Van Dyck, “Phase retrieval through focus variation for ultraresolution in field-emission transmission electron microscopy,” Phys. Rev. Lett. 69(26), 3743–3746 (1992). 5. K. A. Nugent, T. E. Gureyev, D. J. Cookson, D. Paganin, and Z. Barnea, “Quantitative phase imaging using hard X rays,” Phys. Rev. Lett. 77(14), 2961–2964 (1996). 6. D. Paganin, S. C. Mayo, T. E. Gureyev, P. R. Miller, and S. W. Wilkins, “Simultaneous phase and amplitude extraction from a single defocused image of a homogeneous object,” J. Microsc. 206(1), 33–40 (2002). 7. M. E. Friese, T. A. Nieminen, N. R. Heckenberg, and H. Rubinsztein-Dunlop, “Optical torque controlled by elliptical polarization,” Opt. Lett. 23(1), 1–3 (1998). 8. T. E. Gureyev, A. W. Stevenson, D. Paganin, S. C. Mayo, A. Pogany, D. Gao, and S. W. Wilkins, “Quantitative methods in phase-contrast x-ray imaging,” J. Digit. Imaging 13(S1), 121–126 (2000). 9. B. E. Allman, P. J. McMahon, K. A. Nugent, D. Paganin, D. L. Jacobson, M. Arif, and S. A. Werner, “Phase radiography with neutrons,” Nature 408(6809), 158–159 (2000). 10. S. Bajt, A. Barty, K. A. Nugent, M. McCartney, M. Wall, and D. Paganin, “Quantitative phase-sensitive imaging in a transmission electron microscope,” Ultramicroscopy 83(1-2), 67–73 (2000). 11. D. Paganin, A. Barty, P. J. McMahon, and K. A. Nugent, “Quantitative phase-amplitude microscopy. III. The effects of noise,” J. Microsc. 214(1), 51–61 (2004). 12. M. R. Teague, “Deterministic phase retrieval: a Green's function solution,” J. Opt. Soc. Am. 73(11), 1434–1441 (1983). 13. B. D. Arhatari, W. P. Gates, N. Eshtiaghi, and A. G. Peele, “Phase retrieval tomography in the presence of noise,” J. Appl. Phys. 107(034904), 1–7 (2010). 14. B. D. Arhatari, K. Hannah, E. Balaur, and A. G. Peele, “Phase imaging using a polychromatic x-ray laboratory source,” Opt. Express 16(24), 19950–19956 (2008). 15. T. E. Gureyev, A. Pogany, D. M. Paganin, and S. W. Wilkins, “Linear algorithms for phase retrieval in the Fresnel region,” Opt. Commun. 231(1-6), 53–70 (2004). 16. L. D. Turner, B. Dhal, J. Hayes, A. P. Mancuso, K. Nugent, D. Paterson, R. Scholten, C. Tran, and A. Peele, “Xray phase imaging: Demonstration of extended conditions for homogeneous objects,” Opt. Express 12(13), 2960– 2965 (2004). 17. T. E. Gureyev, C. Raven, A. Snigirev, I. Snigireva, and S. W. Wilkins, “Hard X-ray quantitative noninterferometric phase-contrast microscopy,” J. Phys. D Appl. Phys. 32(5), 563–567 (1999). #133365 $15.00 USD Received 13 Aug 2010; revised 5 Oct 2010; accepted 24 Oct 2010; published 27 Oct 2010 (C) 2010 OSA 8 November 2010 / Vol. 18, No. 23 / OPTICS EXPRESS 23727 18. T. E. Gureyev, and S. W. Wilkins, “On x-ray phase imaging with a point source,” J. Opt. Soc. Am. A 15(3), 579– 585 (1998). 19. D. M. Paganin, Coherent X-Ray Optics (Oxford University Press, 2006). 20. D. Paganin, “Studies in phase retrieval,” PhD thesis (University of Melbourne, 1999). 21. A. N. Tikhonov, “Solution of incorrectly formulated problems and the regularization method,” Sov. Math. Dokl. 4, 1035–1038 (1963). 22. B. D. Arhatari, A. P. Mancuso, A. G. Peele, and K. A. Nugent, “Phase contrast radiography: Image modelling and optimization,” Rev. Sci. Instrum. 75(12), 5271–5276 (2004). 23. K. A. Nugent, “Partially coherent diffraction patterns and coherence measurement,” J. Opt. Soc. Am. A 8(10), 1574–1579 (1991). 24. T. E. Gureyev, S. C. Mayo, D. E. Myers, Y. Nesterets, D. M. Paganin, A. Pogany, A. W. Stevenson, and S. W. Wilkins, “Refracting Röntgen‟s rays: Propagation-based x-ray phase contrast for biomedical imaging,” J. Appl. Phys. 105(10), 102005 (2009).


Introduction
Two-plane based phase retrieval relying on free space propagation [1] is one of several existing methods [2][3][4] for quantitative phase imaging.Quantitative methods are important as they can provide useful information as to the materials property of a sample.In the propagation based methods diffraction patterns are recorded in one or more planes.The choice of retrieval algorithm then depends upon the range of the propagation distances, the x-ray energy in the experiment and the properties of the sample.The two-plane based method investigated here uses a solution to a paraxial form of the free-space wave equation known as the Transport of Intensity Equation (TIE).The non-iterative phase retrieval algorithms [5,6], of which the two-plane method is one, find a solution, which is a mapping between the intensity data and the retrieved phase.This non-iterative solution can provide significant computational advantages over iterative methods.TIE-based phase retrieval imaging has been demonstrated for a range of samples using visible light [7], x-rays [8], neutrons [9] and electrons [10].On the whole, these demonstrations have been typified by a variety of choices as to propagation distance and, where two-plane measurements have been used, inter-plane separation.There have been studies that examine the qualitative effects of noise [11,12] and which have shown that extending the separation distances can reduce the noise effect [13].Other studies have examined the effect of x-ray illumination that has a broad spectrum, such as that found in a laboratory source [14].In this paper we investigate the optimisation of TIE-based two-plane phase retrieval when all of these conditions are incorporated simultaneously.We examine conditions that appertain to a particular laboratory set-up but our approach is general enough to be applied broadly.In the following sections we describe the basis for our simulation study and the results of our experiments confirming the simulation results.

Two-plane TIE-based phase retrieval
The TIE has the form [5]: where r is the position vector in a plane perpendicular to the propagation direction, , I/z is the derivative of the intensity with respect to distance along the propagation direction, k = 2π/λ, where λ is the wavelength, and I is the intensity, and φ the phase of the measured wavefield.Within the conditions for applicability of the TIE [15,16] the recovered phase can be taken as a good approximation to the object phase [17,18].The solution for the phase from Eq. ( 1) can be obtained as [19]: In the Fourier domain, the Laplacian has a simple expression given by is the Fourier transform operator and u is the Fourier variable conjugate to the position coordinates.
According to Eq. ( 2) to obtain the phase, we need to have two intensity images taken at different propagation distances with a separation distance δz instead in order to obtain the intensity derivative, I/z.The distance δz needs to be optimised.It has to be small to keep valid and it has to be large enough to reduce the effect of noise [11,12].Hence phase retrieval based on Eq. ( 2) is referred to as two-plane TIE-based phase retrieval.

Simulation studies
The set up of the system can be seen in Fig. 1.
In the simulation we used the parameters of our experimental work as the starting basis and expanded the range in order to explore optimisation.We used 45 μm polystyrene spheres as a sample.The experiment was performed with a laboratory-based x-ray source (Xradia Inc, microXCT) located in the Physics Department, La Trobe University.The source has a closed x-ray tube containing a Tungsten target.A tube voltage of 40 kV was used.In the simulation we use an effective monochromatic x-ray energy of 11.5 keV according to the prescription described in Section 4 [14].The sourcesample distance (z 12 ) used was between 0.01 -1.0 m for the simulations while for the experimental work the physical range was limited to 0.02 -0.125 m.The sampledetector distance (z 23 , z 24 ) used in both simulations and the experimental work was between 0.02 -0.54 m.The geometric magnification used in the experimental work was about 1.2 to 14.In this study we used one detector by taking images at each distance z 23 and z 24 respectively.The imaging detector is a CCD camera coupled with a scintillator system and objective lenses.The objective lens used will determine the field of view and the resolution of the image.The scintillator set in front of the lens is made from CsI(Tl) crystal.The CCD camera (Andor Technology) has 2048x2048 pixels with a physical pixel size of 13.5 μm.A 20x objective lens is used to image the scintillator.

Noise regularization
In Eq. ( 2), low spatial frequencies in the intensity derivative and in the filtered intensity derivative divided by the intensity will be enhanced by the 1/u 2 filtering term in the solution to the phase.This will enhance the low spatial frequency features and therefore create low frequency artefacts [13] for data containing significant amounts of noise.Such low spatial frequency artefacts can be reduced using regularization.This is here done by introducing a Tikhonof regularization parameter, α [6,20,21]: The regularization parameter in Eq. ( 5) acts as a noise suppression filter as it suppresses the effect of pixel to pixel variations in the measured data.In order to choose a value for the regularization parameter for the remainder of this study we simulated our experiment and then performed two-plane TIE phase retrieval on the simulated intensity at two measurement planes for three different noise levels.
In order to model the noise, we use pseudo-random numbers with a normal (Gaussian) distribution with a mean of zero and a standard deviation expressed as a percentage of the mean normalized intensity.In the lens coupled imaging detector, the intensity and noise distribution are not uniform in the entire field of view.The intensity is high in the centre and drops off to the side of the field of view.The noise at the side of the field of view is therefore larger (about 1.5 to 2 times) than the noise in the centre.The reported noise percentage in the simulations in this work is the noise at the centre of the field of view and so underestimates the overall noise slightly.For each pixel location, the noise is added to each raw intensity value and the retrieved phase then compared with the known phase to calculate χ 2 , as shown in Eq. (6).
where the sums are performed over all pixels, i.The process was repeated 10 times with a different set of random numbers each time.The error bars for the χ 2 noise curves are the standard deviations of the resulting distribution of χ 2 values.Figure 2 shows the simulation result for varied values of the regularization parameter for three different noise levels.The higher noise level of 1.4% produces larger estimations of the χ 2 value than the lower noise level of 0.7%.The uncertainty is produced randomly in the retrieved result by the low frequency artefacts due to noise.Therefore, there is no uncertainty for zero noise level.From Fig. 2, we see that the effects of noise can be minimised for certain values of α. Figure 2 also shows that a choice of α can be made for a given noise level that will be close to the optimal result.Here for a noise level at 1.4%, α = 2 is viable, for 0.7%, α = 1 and for 0%, α = 0. Too high a choice of the regularisation parameter will therefore, introduce artefacts in the retrieved shape of an object due to the suppression of low frequency terms [13].

Propagation distance optimization
In practice we cannot avoid the presence of noise in the detected intensities.Noise will change the intensity from pixel to pixel in the detector.This will be similar to the effect of the varying intensity seen from just-resolved Fresnel diffraction fringes for a certain feature size and will lead to errors in the retrieved phase.In order to accurately estimate the intensity derivative in Eq. ( 2) we therefore want the change in the fringes to be big enough so that the signal obtained by subtracting the intensity measured at one detector from the other is significant compared to the noise [11].This suggests a large interplanar separation is desirable.We simulated intensities at two measurement planes for three different noise levels respectively and then performed two-plane TIE phase retrieval.We use a regularisation parameter of σ = 2 for 1.4% noise, σ = 1 for 0.7% noise and σ = 0 for 0% noise.We evaluated the retrieved phase by calculating the χ 2 values for a range of distances, z 12 , z 23 and z 24 under different noise conditions.The simulation result when varying z 23 and z 24 respectively is shown in Fig. 3. Figure 3(a) shows that the shortest values of z 23 produce the best retrieved phase for a given noise level and fixed z 12 and z 24 .The simulation result when varying z 24 for fixed z 12 and z 23 is shown in Fig. 3(b).As expected maximising the inter-plane separation, z 34 , improves the phase retrieval.Figure 3(c) shows the effect of maintaining a fixed inter-plane separation and co-varying z 23 and z 24 .In this case the trends are for shorter propagation distances for both noise and no-noise conditions.Accordingly, in this work we use the smallest value of z 23 (0.02 m) that is compatible with allowing the sample stage to rotate with enough space for a tomography scan.

System resolution on phase retrieval
In spherical beam geometry, such as for a laboratory x-ray source, the system resolution can be limited by two options.Firstly, under small geometric magnification, it is limited by the detector system [22].We use a high resolution detector in our experimental system.Secondly, under large geometric magnification, the system resolution is limited by the source size.The x-ray source size in the system is about 6 μm full-width at half maximum.The system resolution for both conditions is shown in Fig. 6.In a noise free environment with no limitation on detector resolution, the best condition for phase recovery is when both images are taken in the small propagation distances.In this case source blurring is negligible and the retrieved phase will be perfect for all feature sizes in the object.Accordingly, the shortest propagation distance for which the Fresnel diffraction fringes for a given feature can be resolved is one criterion.In a system with noise, we need to have a sufficient separation distance for both images to reduce the effect of noise.Accordingly, one image needs to be taken in large magnification where the source blurring will have an effect.
To evaluate this effect on the phase retrieval, we simulated a kapton grating sample with variations in the object sizes.Maximum thickness of the sample is 20 μm.We consider an object with a one dimensional sinusoidal variation.The results of evaluating the χ 2 values for the actual and retrieved profiles as a function of object periods for a single wavelength (1 Å) are plotted in Fig. 7 with a comparison between 0% noise and 0.4% noise.We use a regularisation parameter of σ = 0 for 0% noise and σ = 2 for 0.4% noise.Simulations have been made for various set-ups.The corresponding magnified or demagnified sources due to the configuration setup are presented in the right column.We assume the source has a Gaussian distribution characterised by a full-width at half maximum of 6 μm.The recorded intensity in the planes z 23 and z 24 respectively will be blurred by the convolution of the appropriately scaled source distribution [23] which will then affect the phase retrieval results.
It can be seen from the top row of Fig. 7 that the inherent filtering action of the phase retrieval step degrades also the resolution of the result.If we define an acceptable retrieval as when χ 2 < 0.2, then for the zero noise case in the top row of Fig. 7 we see that the smallest object that can be retrieved acceptably is approximately 7 μm.Similarly, for zero noise case, larger effective source sizes result in further resolution loss as shown in the middle and bottom rows (9 μm and 12 μm smallest acceptably retrieved objects respectively).When resolution levels of noise are included it can also be seen that the system resolution is not significantly degraded (middle and bottom rows).However, noise can degrade the resolution level overall in the set-up of the top row where the inter-plane separation distance is not large enough to act as noise suppression.

Polychromatic effects
If quantitative phase retrieval is desired then, when a broad-band source is used, such as a laboratory x-ray source, a way must be found to deal with the fact that the measured intensity is a weighted sum of individual wavelength contributions.The weighting function will include contributions from the source spectrum, the detector response function and phase and absorption effects in the geometric path, including from the sample.In some cases the sample effects can be isolated so that quantitative information about the sample can be retrieved [14].Such a technique is also applicable for two-plane TIE-based phase retrieval because the wavelength dependent term in Eq. ( 2) can be removed by division with k on both sides and by using     kt    rr , where δ is the decrement of the real part refractive index and t is the projected thickness of a homogenous sample.For the monochromatic case Eq. ( 2) can be written as: For the polychromatic case, we need to incorporate the spectral response function, S(λ) where we assume that for a low absorption object     0 I r I r  and we can replace the right hand side of Eq. ( 8) with that of Eq. ( 7) with the substitution poly II   where is I poly defined in the same way as δ poly .This expression can easily be calculated as in an experiment we measure I poly and dI poly /dz.Retrieved δ poly t vs its real value for kapton sinusoidal sample.
The sample for the polychromatic case simulation is a sinusoidal profile sample with a period of 45 μm made from a polyimide film (Kapton, C 22 H 10 N 2 O 4 ) with a density of 1.45 g/cm 3 .The chosen period is big enough so that no blurring affects the retrieved phaseas discussed in Section 3.3.The amplitude of the sample is 20 μm.The simulation geometry is z 12 = 0.1 m, z 23 = 0.02 m, z 24 = 0.125 m.We measured the spectral response function, S(λ), which takes into account factors such as the source spectrum, the detection efficiency of the energy sensitive detector used, absorption in the beryllium exit window, and the scintillator efficiency of the image detector.The spectral response function at 40 kV tube voltage is shown in Fig. 8(a) [14].We then calculate δ poly using the known spectral response function S(λ) according to Eq. ( 8) and obtain a value of 2.36x10 6 .This corresponds to the value for δ for kapton at an energy of 11.5 keV.We then calculate the phase retrieval value for δ poly t for zero noise using the simulated intensities.The result is shown in Fig. 8(b) and is in excellent agreement with the values calculated.

Experimental results
We used a sinusoidal kapton sample with 42.5 μm period and a maximum amplitude of 20 ± 2 μm as measured using an optical microscope.From the results from the previous sections we select the following condition to obtain a known "good" phase retrieval.Two images were taken at z 23 = 0.02 m and z 24 = 0.125 m, as shown in Fig. 9(a) and 9(b) at fixed z 12 = 0.1 m.The retrieved δ poly t obtained using the method discussed in the previous section, is shown in Fig. 9(c).We used the regularization parameter of α = 2.The experimental result was averaged in the row direction to increase the signal to noise ratio and the resulting average δ poly t is shown in Fig. 10. Figure 10 shows excellent agreement between the experimentally phase retrieved results and the calculated value for δ poly t.We also evaluated the phase retrieval by comparing simulations with experimental data for 45 μm diameter polystyrene spheres, as shown in Fig. 11.The noise level in each simulation result has been set according to the measured noise levels from the experimental data for each configuration.Figure 11(a) used the same experimental configuration as before and in Fig. 11(b) the results for z 12 = 0.12 m, z 23 = 0.02 m, z 24 = 0.54 m, are shown.From Fig. 4 we expected that both cases will provide reasonable phase retrieval with the configuration for Fig. 11(b) being better, in the sense that the χ 2 for the simulation is less.Figure 11(c) shows the result from the first case shown in Fig. 7 (z 12 = 0.1 m, z 23 = 0.02m, z 24 = 0.05 m) where the shorter z 24 is expected to produce poorer noise suppression.Again the expectation is borne out.Finally Fig. 11(d) shows the result when the set-up used in Fig. 11(a) is changed so that z 23 = 0.09 m.As expected the decreased inter-plane separation distance degrades the result.The cloudy background in all the images is due to the noise in the measurement which produces low spatial frequency artefacts in the retrieved phase.As the noise suppression is worse in the set-up for Fig. 11(d  The two-plane TIE-based phase retrieval images can easily be extended to three dimensional imaging once the corresponding phase retrieved images have been obtained for each projection angle.Here, we used polystyrene spheres mixed with plastic fibre as a sample.The same 45μm polystyrene diameter has been used.A tube voltage of 40 kV was used.The optimal parameters for the setup, as discussed above, were used z 12 = 0.1 m, z 23 = 0.02 m and z 24 = 0.125 m.A total of 361 projections were acquired over an angular range of 180°.A reconstructed slice is shown in Fig. 12(a) and a three-dimensional rendering is shown in Fig. 12(b).As with the two-dimensional case, the three-dimensional result also shows the effect of noise as low spatial frequency artefacts, or "clouds" in the image.

Conclusion
For the fairly generic objects described here, a large distance between detector 1 and detector 2, when detector1 is placed as close as possible to the sample, is required to have good quality retrieved phase.When detector 1 is placed at a larger distance then the separation distance must be much larger again.Simulations showed that a source sample distance larger than 0.1m will produce a good retrieved phase.Increasing the source sample distance will help to reduce the effect of noise in the phase retrieved image in the case when the two plane separation distances are big enough.These results were demonstrated in simulation and validated for a sinusoidal sample in an experimental result.Further experimental results were also consistent with the expectations from simulation and were able to show that the same benefits can be achieved in three-dimensional imaging.Increasing the detector separation distances will help to further improve the retrieved phase due to noise but the effect of source blurring will start to degrade the resolution.If distances are increased too much exposure times will have to be increased significantly to maintain comparable noise characteristics compared to shorter detector distances.
We have repeated this investigation for different image maps, compositions and energies in the simulation work and have reproduced the same broad findings described above.It is however, difficult to be quantitative for arbitrary samples.As the precise trade-off between sample feature size, noise and source and detector distances will vary we have therefore presented here a procedure that can be used to explore the relevant parameter space.Such a procedure can be used to refine the qualitative rules for "good" imaging as demonstrated in section 3 and as discussed elsewhere [24].It is expected that a good set of operating conditions can be found that will apply to a broad range of samples.

Fig. 1 .
Fig. 1.The set up of the system used in this study: z 12 is the source to sample distance, z 23 is the sample to detector1 distance, z 24 is the sample to detector2 distance and z 34 is the inter plane separation between the two detector positions.Nominal images recorded at the detector positions are shown in the insets.In the projection magnification case used here what has become referred to as the Fresnel scaling theorem is used [19].In the TIE solution the measured intensity is scaled by the geometric magnification and the propagation distance (z 23 , z 24 ) is replaced by an effective propagation distance, as shown in Eq. (3).

Fig. 3 .
Fig. 3. Simulation result for (a) fixed z 12 = 0.12 m and z 24 = 0.54 m with varied distance z 23 .(b) fixed z 12 = 0.12 m and z 23 = 0.02 m with varied distance z 24 .(c) fixed z 12 = 0.12 m and z 34 = 0.25 m with co-varied distance z 23 and z 24 .Note the result for 0% noise is multiplied by 10x to make it visible.

Figure 4
Figure 4 shows the effect of varying z 12 for different fixed z 23 and z 24 .The results are consistent with there being a reduction in χ 2 with increased z 34eff .Similarly, for a given z 12 Fig.4(c) shows a lower χ 2 than Fig. 4(a) as z 34eff is larger.The same result holds between Fig. 4(d) and 4(b).The relevant values of z 34eff can be found in Eq. (4) and are plotted in Fig. 5.

Fig. 5 .
Fig. 5. Plot of z 34eff as a function of source sample distance variation, z 12 , for the set-ups in Fig. 4.

Fig. 6 .
Fig. 6.System resolution for (a) small geometric magnification (b).large geometric magnification.Nominal images recorded by the detector are shown in the insets.

Fig. 7 .
Fig. 7. χ 2 values are plotted as a function of object periods for various set-ups at 0% noise (left column) and 0.4% noise (middle column).The corresponding scaled Gaussian sources are plotted in the right column.

Fig. 8 .
Fig. 8. (a).The spectral response function of the system taken at 40kV tube voltage.(b).Retrieved δ poly t vs its real value for kapton sinusoidal sample.

Fig. 10 .
Fig. 10.Comparison between the retrieved δ poly t from experimental data and the real value using δ poly of 2.36x10 6 .

Fig. 11 .
Fig. 11.Simulation vs experimental result for varied distances as explained in the text.