Noniterative spatially partially coherent diffractive imaging using pinhole array mask

Abstract. We propose and experimentally demonstrate a noniterative diffractive imaging method for reconstructing the complex-valued transmission function of an object illuminated by spatially partially coherent light from the far-field diffraction pattern. Our method is based on a pinhole array mask, which is specially designed such that the correlation function in the mask plane can be obtained directly by inverse Fourier transforming the diffraction pattern. Compared to the traditional iterative diffractive imaging methods using spatially partially coherent illumination, our method is noniterative and robust to the degradation of the spatial coherence of the illumination. In addition to diffractive imaging, the proposed method can also be applied to spatial coherence property characterization, e.g., free-space optical communication and optical coherence singularity measurement.


Introduction
Coherent diffractive imaging (CDI) is an important tool for reconstructing the complex-valued transmission function of an object from the far-field diffraction pattern and has been widely applied in material and biological sciences. 1,2 Miao et al. 3 first experimentally realized imaging of submicrometer sized noncrystalline specimen using CDI. Many CDI approaches have been developed in the past decades; they can be divided into two types: the iterative methods [4][5][6][7][8][9] and the noniterative methods. [10][11][12][13][14] However, most of the traditional CDI methods cannot be directly applied to using spatially partially coherent illumination without a proper modification, which hence have limited applications at short wavelengths, e.g., in the x-ray and electron regime, or in the unstable experimental environment. For example, the degradation of the spatial coherence may also be caused by the mechanical movement of the sample and the experimental setup, or by the fluctuation of the ambient medium, e.g., atmospheric turbulence. [15][16][17] Iterative algorithms retrieve the complex-valued transmission function of an object by propagating the field back-and-forth between the object plane and the far-field diffraction plane, and imposing constraints on the field in both planes. Gerchberg and Saxton pioneered the iterative algorithms in 1972 by proposing a straightforward method using two intensities measured in the object plane and in the far-field, respectively. 4 Iterative algorithms using only one intensity measurement of the far-field diffraction pattern, as proposed by Fienup,5,6 require prior knowledge of the support of the object for imposing the support constraint in the object plane. Recently, ptychographic algorithms have become an essential technique for imaging nanoscale objects using short-wavelength sources. 8 The key feature is that illumination areas at the neighboring shift positions must overlap, and this overlap improves the convergence of the ptychographic algorithms. 7,8 Recent developments in ptychography allow for simultaneous reconstruction of the probe. This significantly reduces the complexity of experimental setup compared to other CDI methods. 9 For spatially partially coherent illumination, the propagation of light is described using the mutual coherence function (MCF) instead of the field. Many efforts have been spent to adapt CDI methods to spatially partially coherent illumination. [15][16][17][18][19][20][21][22][23][24] The modification of the iterative algorithms was first reported by Whitehead et al. 15 using mode decomposition of the MCF. 16 Later, ptychographic algorithms have also been modified to work for spatially partially coherent illumination by decomposing the MCF into orthogonal modes. [17][18][19][20][21] Thibault and Menzel 17 proposed a mixed state model from the quantum perspective, which effectively deals with a series of multistate mixing problems including partially coherent illumination and enables more applications of ptychography, such as continuous-scan ptychography 18 and dynamic imaging of a vibrating sample. 19 Furthermore, ptychographic imaging with the simultaneous presence of both multiple probe and multiple object states was also demonstrated. 21 However, the accuracy of mode decomposition relies on the number of modes for accurately representing the MCF, and it increases as the spatial coherence of the illumination decreases.
Compared to iterative methods, noniterative methods 10,11,14 do not suffer from issues such as stagnation or nonuniqueness of the solution to the diffractive imaging problem, especially when the illumination becomes spatially partially coherent. 21,24 In holographic methods, the field transmitted by the object is perturbed such that the object's transmission function can be directly extracted from the inverse Fourier transform of the diffraction pattern. 14 This perturbation can be achieved by introducing a pinhole, e.g., Fourier transform holography (FTH), [25][26][27] or by changing the transmission function at a point of the object, e.g., Zernike quantitative phase imaging. 28 The performance of applying holographic methods to spatially partially coherent illumination has been discussed in Ref. 14. Alternative methods extract the object information from the three-dimensional autocorrelation functions obtained by inverse Fourier transforming the three-dimensional data set (e.g., the data set measured by varying focus 12 or another optical parameter 13 ).
It has been demonstrated that using noniterative methods can avoid errors due to truncating the number of the modes for representing the MCF. [12][13][14] However, the degree of the spatial coherence of the illumination limits the field of view (FOV) of the reconstructed object. To be precise, what is reconstructed is a product of the object's transmission function and a correlation function of the illumination. This correlation function has a maximum at the perturbation point, and its value decreases at a rate that depends on the degree of spatial coherence, as the distance between the perturbation point and the observation point increases. Therefore, the lower the illumination's degree of spatial coherence is, the smaller the region of the object that can be reconstructed reliably.
In this paper, we propose a noniterative method based on a pinhole array mask (PAM). We place the PAM between the object and the detector, and we measure the far-field diffraction pattern of the spatially partially coherent field transmitted by the PAM. The PAM consists of a periodic array of measurement pinholes and an extra reference pinhole, which is analogous to the perturbation point in FTH. 14 It forms an interference between fields transmitted by the reference pinhole and by the measurement pinholes, and we can directly retrieve the correlation function of the incident light with respect to the reference pinhole at the locations of the measurement pinholes from the interference pattern.
In FTH, since the reference pinhole is placed far from the measurement window, the FOV is rather small due to the low correlation. Our method is advantageous compared to FTH, because splitting the measurement window into a periodic array of measurement pinholes keeps the reference pinhole close to all measurement pinholes and thus maintains a high correlation that results in a large FOV.
Our method places the object at a certain distance before the PAM, instead of superposing the object with the PAM. In practice, this not only offers flexibility for arranging the experimental setup but also allows us to adjust the sampling of the reconstructed object. When the propagation distance satisfies the condition for Fresnel or Fraunhofer diffraction, an object with finite support can be reconstructed from the retrieved correlation function, and its sampling is related to the sampling of the PAM by the Shannon-Nyquist sampling theorem.
Because our method reconstructs the product of the object's transmission function and the illumination's correlation function and hence can be used not only for object reconstruction but also for characterization of the spatial coherence structure, it is useful for a broad range of applications in coherent optics, such as the measurement of optical coherence singularity 29,30 and free-space optical communication through turbulent media. 31,32

Methods
The schematic plot of the experimental setup of our method is shown in Fig. 1, which shows that a transmissive object is illuminated by spatially partially coherent light. In our method, unlike traditional CDI algorithms, we place a PAM between the object and the detector. The location of the PAM is chosen such that light propagation from the object to the PAM and from the Fig. 1 Schematic plot of the experimental setup and the concept of the noniterative diffractive imaging method. (a) A PAM is placed between the object and the camera. The PAM is specially designed such that we can retrieve the correlation function of the incident light by inverse Fourier transform of the measured diffraction pattern. (b) The PAM consists of a reference pinhole at the origin (gray square) and a periodic array of the measurement pinholes with certain offset (white squares). (c) In the inverse Fourier transform of the diffraction pattern, the autocorrelation terms (gray squares) and the two interference terms (red and blue squares) are separated by the offset. Each interference term directly contains the correlation between the fields at the reference pinhole and at the measurement pinholes.
PAM to detector obeys either Fresnel or Fraunhofer propagation and hence can be described using Fourier transforms. By doing so, we can divide our method into two steps: (1) retrieving the correlation function of incident light in the PAM plane using a noniterative approach and (2) reconstructing the product of the object's transmission function and the illumination's correlation function using a differential method, which requires two diffraction patterns corresponding to the object with and without transmission perturbation, respectively. It is worth noting that the use of this differential method is not necessary for completely spatially coherent illumination.

Retrieving the Correlation Function in the PAM Plane
Let the coordinate of the PAM plane be denoted by r ¼ ðx; yÞ.
The PAM consists of a reference pinhole at the origin shown by the gray square in Fig. 1(b), and a periodic array of measurement pinholes around the reference pinhole is shown by the white squares in Fig. 1(b). The center of the periodic array is shifted relative to the reference pinhole by certain offset and is depicted by the white spot at the corner of the reference pinhole. We assume that the reference pinhole and the measurement pinholes are identical and are so small that each pinhole can be approximated by a delta function. This assumption allows us to write the transmission function of the PAM by where 0 ¼ ð0; 0Þ denotes the location of the origin, r m;n ¼ ðΔx þ mp x ; Δy þ np y Þ denotes the location of the measurement pinhole in the periodic array, where ðm; nÞ is the index of the measurement pinhole, ðp x ; p y Þ is the pitch of the periodic array, and ðΔx; ΔyÞ is the offset of the periodic array relative to the reference pinhole. The incident light transmitted by the PAM generates a diffraction pattern in the detector plane. We denote the MCF of the incident light in the PAM plane by Wðr 1 ; r 2 Þ, which describes the correlation between the fields at r 1 and r 2 . Because the light propagation from the PAM to the detector satisfies the condition for either Fresnel or Fraunhofer propagation, we can express the diffraction pattern measured by the detector using the Fourier transform as follows: Wð0;r m 2 ;n 2 Þexpðþi2πk · r m 2 ;n 2 Þ; where k denotes the detector plane coordinate. Denoting the original sampling grid of the detector plane by k 0 , we have k ¼ k 0 ∕ðλzÞ, where λ is the wavelength of illumination and z is the propagation distance, and k is conjugated with the coordinate of the PAM plane r. By taking the inverse Fourier transform of the diffraction pattern [Eq.
(2)], we obtain Wð0;r m 2 ;n 2 Þδðrþr m 2 ;n 2 Þ; where F −1 denotes the operation of the inverse Fourier transform. We can observe that Eq. (3) consists of four terms: 1. The first term is Wð0; 0Þ, which is a constant multiplied by a delta function that appears only at the origin of the coordinate system rð0; 0Þ. 2. The second term is Wðr m 1 ;n 1 ; r m 2 ;n 2 Þ located on the periodic array defined by r ¼ r m 1 ;n 1 − r m 2 ;n 2 ¼ ½ðm 1 − m 2 Þp x ; ðn 1 − n 2 Þp y , which is depicted by the gray squares in Fig. 1(c). This periodic array has the same pitch as the periodic array of the measurement pinholes ðp x ; p y Þ but zero offset relative to the origin. 3. The third term is Wðr m;n ; 0Þ located on the periodic array defined by r ¼ r m;n , which is depicted by the blue squares in Fig. 1(c). 4. The fourth term is Wð0; r m;n Þ ¼ Wðr m;n ; 0Þ Ã located on the periodic array defined by r ¼ −r m;n , which is depicted by the red squares in Fig. 1(c).
The role of the reference pinhole of the PAM is analogous to the perturbation point in FTH, 14 namely to create interference between the incident light transmitted by the reference pinhole and by the measurement pinholes. The interference induces the third term and the fourth term in Eq. (3), which are referred to as the "interference terms." The first term and the second term in Eq. (3) are the autocorrelation of the reference pinhole and of the measurement pinholes, respectively, and hence are referred to as the "autocorrelation terms." The layout of the four terms of Eq. (3) is illustrated by the schematic plot in Fig. 1(c). Figure 1(c) shows that the periodic arrays of the autocorrelation term (the gray squares) and the two interference terms (the blue squares and the red squares) have the same pitch but different offset. This allows us to separate the two interference terms and the autocorrelation term by multiplying a spatial filter to the inverse Fourier transform of the diffraction pattern. The expression for this spatial filter is given by We can multiply F −1 ½IðkÞðrÞ either by M 0 ðrÞ to obtain Wðr m;n ; 0Þ or by M 0 ð−rÞ to obtain Wðr m;n ; 0Þ Ã . As a consequence, we can retrieve the correlation function Wðr m;n ; 0Þ of the incident light in the PAM plane with respect to the location of the reference pinhole r ¼ 0 at the locations of the measurement pinholes r ¼ r m;n . Note that Wðr m;n ; 0Þ is a function of the locations of the measurement pinholes on the periodic array defined by r ¼ r m;n . The sampling interval of r m;n is given by the pitch of the PAM ðp x ; p y Þ. For rectangular pinholes with size ðw x ; w y Þ, the sampling interval should be ðp x ≥ 3w x ; p y ≥ 3w y Þ such that the autocorrelation term and the two interference terms do not overlap as shown in Fig. 1(c). However, this sampling interval ðp x ; p y Þ is usually larger than the diffraction-limited sampling interval according to the Shannon-Nyquist criterion: (p x;0 ¼ λz∕l x , p y;0 ¼ λz∕l y ), where ðl x ; l y Þ is the size of the detector.

Reconstructing the Transmission Function of the Object and the Correlation Function of the Illumination
The sampling interval of the reconstructed object will be higher in the plane at distance z 0 before the PAM plane than exactly in the PAM plane. We denote the coordinate of the object plane by ρ, which is original sampling grid ρ 0 normalized by λz 0 , and the complex-valued transmission function of the object by OðρÞ.
In our method, the MCF in the object plane and the MCF in the PAM plane are related to the Fourier transform in the case of either Fresnel or Fraunhofer propagation. Here we give the example of Fraunhofer propagation as follows: where W 0 ðρ 1 ; ρ 2 Þ is the MCF of the incident beam that describes the correlation between fields at ρ 1 and ρ 2 . According to the Shannon-Nyquist criterion, the pitch of the PAM ðp x ; p y Þ determines the size of the object: (l x ≤ λz 0 ∕p x , l y ≤ λz 0 ∕p y ). Therefore, we need to find the propagation distance z 0 such that the size of the reconstructed object matches the pitch of the PAM. In Eq. (5), by setting r 1 ¼ r m;n and r 2 ¼ 0, we can obtain an expression for computing the retrieved correlation function in the PAM plane as follows: Wðr m;n ; 0Þ ¼
We shall note that the modulation TðρÞ depends on not only the MCF of the illumination W 0 ðρ 1 ; ρ 2 Þ but also the transmission function of the object OðρÞ. To eliminate this modulation, we use the differential approach. This approach requires two measurements: one with a point perturbation to the transmission function of the object at ρ ¼ ρ 0 and the other without perturbing the object. The perturbation is achieved by the change of either the amplitude or the phase of the transmission for ρ 0 inside the object or by introducing an extra spot that lets light pass through for ρ 0 outside the object. It is worth mentioning that for completely spatially coherent illumination, the MCF of the illumination, W 0 ðρ 1 ; ρ 2 Þ, becomes a constant, and there is no need for using the differential approach to eliminate the modulation TðρÞ.
Substituting the transmission function of the perturbed object OðρÞ þ Cδðρ − ρ 0 Þ, where C is a complex-valued constant representing the perturbation, into Eq. (6), we obtain W δ ðr m;n ; 0Þ ¼ Expanding the brackets in this expression leads to W δ ðr m;n ; 0Þ ¼ Wðr m;n ; 0Þ þ jCj Using the property of the delta function, we can derive the expression of the retrieved correlation function in the PAM plane for the perturbed object as follows: W δ ðr m;n ;0Þ¼ Wðr m;n ;0ÞþjCj 2 W 0 ðρ 0 ;ρ 0 Þ þ Z C·Oðρ 2 Þ Ã W 0 ðρ 0 ;ρ 2 Þd 2 ρ 2 expð−i2πρ 0 ·r m;n Þ þ Z Oðρ 1 ÞC Ã W 0 ðρ 1 ;ρ 0 Þexpð−i2πρ 1 ·r m;n Þd 2 ρ 1 : Finally, we take the inverse Fourier transform of the difference between the retrieved correlations in the PAM plane for the unperturbed object and for the perturbed object, and the result yields F −1 ½W δ ðr m;n ; 0Þ − Wðr m;n ; 0Þ Neglecting the delta functions in Eq. (12), we can obtain the product of the transmission function of the object OðρÞ and the correlation function W 0 ðρ; ρ 0 Þ, which describes the correlation between the fields at the point of perturbation ρ 0 and at all other points W 0 ðρ; ρ 0 Þ. The correlation function is determined by the MCF of the illumination but also depends on the location of the perturbation ρ ¼ ρ 0 . Usually, W 0 ðρ; ρ 0 Þ is simply Gaussian distributed, e.g., when using the Gaussian-Schell-model (GSM) beam for illumination. However, when the MCF of the illumination is complicated and not known, we need to calibrate W 0 ðρ; ρ 0 Þ by performing a reconstruction using an empty window as object and then divide the reconstructed product OðρÞW 0 ðρ; ρ 0 Þ by W 0 ðρ; ρ 0 Þ.

Results and Discussions
In the experiment, we use the GSM beam and the Laguerre-Gaussian-Schell-model (LGSM) beam as illumination to validate our method. In the object plane, the MCF of the GSM beam can be expressed as where w 0 and σ are the width of the Gaussian distribution of the intensity distribution and the correlation function, respectively. The experimental setup for generating the GSM beam is shown in Fig. 2. We expand a coherent laser beam at wavelength λ ¼ 532 nm using a beam expander (BE) and then focus it on a rotating ground-glass disk (RGGD) using a lens L1. Because the focal spot follows a Gaussian distribution, the spatially partially coherent light generated due to the scattering by the RGGD satisfies Gaussian statistics, namely, the correlation between the fields in any pair of points follows a Gaussian distribution. We then collimate the spatially partially coherent light by lens L2. By passing the collimated beam through a Gaussian amplitude filter (GAF), we can obtain the GSM beam whose intensity distribution also follows the Gaussian distribution. The MCF of the LGSM beam in the object plane is described by where L m n ½ is the Laguerre polynomial of order n and m ¼ 0. The experimental generation of LGSM beam has been reported in Refs. 33 and 34. In the experimental setup shown in Fig. 2, we need to insert a spiral phase plate between the BE and the focusing lens L1, which produces a dark hollow focal spot on the RGGD. The order n of the LGSM beam is determined by the topological charge of the spiral phase plate. When n ¼ 0, the spiral phase plate has a constant phase and the LGSM beam becomes the GSM beam. However, when n ≠ 0, the MCFs of the LGSM beam and the GSM beam have the same amplitude but different phases.
In the experiment, the beam width w 0 of the Gaussian distribution of the intensity distribution is determined by the GAF and is set to be 0.85 mm, whereas the coherence width σ of the Gaussian distribution of the correlation function, also known as the degree of coherence, is determined by the size of the focal spot on the RGGD. We can control σ by translating back-andforth the focusing lens L1, which determines the size of the focal spot. The degree of spatial coherence σ is calibrated using the method proposed in Refs. 35-37. For the diffractive imaging experiment, we use a phase object, whose amplitude is flat and phase has a binary distribution (0.1π and 0.9π) in the shape of a panda, displayed on the phase spatial light modulator (SLM) (Pluto, Holoeye Inc., with resolution 1920 × 1080 pixel size 8 μm, and frame rate 60 Hz). The pixels inside the support of the object reflect the incident beam back to the beam splitter, whereas the pixels outside the support direct the incident beam to other directions. The beam reflected by the SLM propagates to the PAM. Finally, we measure the farfield diffraction pattern of the light transmitted by the PAM using a charge-coupled device (CCD) (Eco445MVGe, SVS-Vistek Inc. with resolution 1296 × 964 pixel size 3.75 μm, and frame rate 30 fps) camera, which is placed at the focal plane of the Fourier transform lens L3. In the experiment, we set the pitch of the PAM and the size of the pinhole to be p x ¼ p y ¼ 270 μm and w x ¼ w y ¼ 54 μm, respectively. The object size is l u ¼ l v ¼ 1.92 mm, which requires the propagation distance between the object and the PAM to be z 0 ¼ 1170 mm.

Experimental Results Using GSM Beam Illumination
Equation (12) indicates that for spatially partially coherent illumination, to reconstruct the product of the object's transmission function OðρÞ and the illumination's correlation function W 0 ðρ; ρ 0 Þ, our method needs two measurements of the diffraction pattern, one without the perturbation and the other with the perturbation of the object's transmission at ρ 0 . W 0 ðρ; ρ 0 Þ describes the correlation between the fields at the perturbation point ρ 0 and other points ρ, which decreases as the distance between ρ 0 and ρ increases. As a consequence, the reconstructed object's transmission function OðρÞ has a limited FOV since the OðρÞ cannot be reconstructed at locations, where W 0 ðρ; ρ 0 Þ is corrupted by the noise.
In the experiment, we place the perturbation point at the head of the panda by −0.3π. We show the object's transmission function with and without the perturbation point in Figs. 3(a) and 3(b), and the amplitude and the phase of the reconstructed product for various degrees of spatial coherence in Figs. 3(c1)-3(c3) and 3(d1)-3(d3). Because the MCF of the GSM beam has a uniform phase, the phase of the reconstructed product is only given by the phase of the object. The amplitude of the reconstructed product follows the Gaussian distribution of the MCF of the GSM beam. The panda shape in the amplitude is due to the discontinuity of the phase and low-pass filtering. As mentioned in Ref. 38, it is the phase jump between the inside and outside area that enables destructive interference along the outline, therefore, leading to the observation of the dark panda shape in a very bright background. In addition, the finite boundaries of PAM, Fourier transform lens (L3), and CCD constitute low-pass filters, which result in the disappearance of the panda contour that corresponds to high spatial frequencies. 21 We can observe in Fig. 3 that for a lower degree of spatial coherence σ, the amplitude of the correlation function W 0 ðρ; ρ 0 Þ decreases more Fig. 2 The experimental setup for generating the GSM beam and for diffractive imaging. BE, beam expander; RGGD, rotating ground-glass disk; GAF, Gaussian amplitude filter; BS, beam splitter; SLM, spatial light modulator; PAM, pinhole array mask; L1, L2, and L3, thin lenses; and CCD, charge-coupled device.
rapidly as the distance ρ 0 − ρ increases, and hence the FOV of object's transmission function OðρÞ is smaller. Figure 3 shows that, to increase the FOV, we can either increase the degree of spatial coherence σ or decrease the noise level. In Fig. 4, we demonstrate that using more than only one perturbation point, placed at different locations of the object, the object's transmission function can still be reconstructed in the whole FOV in the case of the lowest degree of spatial coherence (σ ¼ 0.21 mm). This requires us to repeat the measurement and the reconstruction procedure for each perturbation point to reconstruct different parts of the object. By combining the different parts reconstructed using low σ illumination together, we can obtain the object's transmission function as if using high σ illumination.

Experimental Results Using LGSM Beam Illumination
In Figs. 3 and 4, the phase of the reconstructed product is given by only the object since the MCF of the GSM beam has a uniform phase. However, for illumination using LGSM beam, the phase of its MCF is not uniform. Therefore, we need to calibrate W 0 ðρ; ρ 0 Þ so that we can divide the reconstructed product by W 0 ðρ; ρ 0 Þ to obtain the object's transmission function OðρÞ alone. We show the amplitude and the phase of the reconstructed product using LGSM illumination beam in Fig. 5(a). Compared to the reconstructed results using GSM illumination beam in Figs. 3 and 4, we can see that now the phase of the object's transmission function OðρÞ is modulated by the phase of the  correlation function W 0 ðρ; ρ 0 Þ of the LGSM beam, and we cannot see the panda in the phase of the reconstructed product. We show the amplitude and the phase of the correlation function W 0 ðρ; ρ 0 Þ calibrated using an empty window as the object in Fig. 5(b). In Fig. 5(c), we demonstrate the object's transmission function OðρÞ obtained by dividing the reconstructed product OðρÞW 0 ðρ; ρ 0 Þ by the calibrated W 0 ðρ; ρ 0 Þ. The panda in the phase of the reconstructed object can clearly be seen. This example verifies that our method can be applied to object reconstruction in cases using an illumination beam whose MCF is not known prior.

Conclusion
In summary, we develop and validate a noniterative method to reconstruct the complex-valued transmission function of an object illuminated by spatially partially coherent beam using a PAM placed between the object and the detector. Our method overcomes several challenges of conventional iterative CDI algorithms and holographic methods. In particular, our method does not depend on the mode decomposition of the MCF of the spatially partially coherent light and has the freedom to choose the location of the point where the transmission function of the object is perturbed, which is particularly beneficial for achieving large FOV when using a low degree of spatial coherence illumination. Moreover, we also demonstrate that our method can be used to calibrate the MCF of an arbitrary spatially partially coherent beam. This calibration allows us to reconstruct the object's transmission function almost as accurately as if using complete illumination. The calibration itself can also be used for spatial coherence property characterization, which is needed as an approach for applications like the measurement of optical coherence singularity. 29,30 Therefore, in addition to diffractive imaging, our method also provides an approach. Finally, our method is wavelength independent and can be applied to a wide range of wavelengths, from x-rays to infrared light.