• Photonics Research
  • Vol. 9, Issue 6, 1003 (2021)
Cheng Shen1、*, Mingshu Liang1, An Pan2, and Changhuei Yang1
Author Affiliations
  • 1Department of Electrical Engineering, California Institute of Technology, Pasadena, California 91125, USA
  • 2Xi’an Institute of Optics and Precision Mechanics (XIOPM), Chinese Academy of Sciences (CAS), Xi’an 710119, China
  • show less
    DOI: 10.1364/PRJ.419886 Cite this Article Set citation alerts
    Cheng Shen, Mingshu Liang, An Pan, Changhuei Yang. Non-iterative complex wave-field reconstruction based on Kramers–Kronig relations[J]. Photonics Research, 2021, 9(6): 1003 Copy Citation Text show less

    Abstract

    A non-iterative and non-interferometric computational imaging method to reconstruct a complex wave field called synthetic aperture imaging based on Kramers–Kronig relations (KKSAI) is reported. By collecting images through a modified microscope system with pupil modulation capability, we show that the phase and amplitude profile of the sample at pupil limited resolution can be extracted from as few as two intensity images by using Kramers–Kronig (KK) relations. It is established that as long as each subaperture’s edge crosses the pupil center, the collected raw images are mathematically analogous to off-axis holograms. This in turn allows us to adapt a recently reported KK-relations-based phase recovery framework in off-axis holography for use in KKSAI. KKSAI is non-iterative, free of parameter tuning, and applicable to a wider range of samples. Simulation and experiment results have proved that it has much lower computational burden and achieves the best reconstruction quality when compared with two existing phase imaging methods.

    1. INTRODUCTION

    Coherent optical field can be described as a 2D complex function under scalar diffraction theory [1]. In most practical situations, we are not able to fully measure the entire complex function because optical detectors are only suited for detecting intensity but not the phase of light. Yet phase image profiles are highly informative and useful in various applications—particularly in life science research, where phase images of living unstained cells can reveal cell structures that are otherwise invisible. To address this need, the optical imaging community has actively worked on and developed methods to infer phase information from pure intensity measurements over the past century or so.

    To date, existing quantitative phase measurement methods can be generally categorized as interferometric and non-interferometric methods. Digital holography [26], phase shifting interferometry [7,8], and optical coherence tomography [9,10] belong in the former group. Non-interferometric methods include iterative phase (diversity) retrieval [1113], (Fourier) ptychography [1417], transport of intensity equation [18,19], and quantitative differential phase contrast [2022]. Non-interferometric methods are inherently attractive, as they are generally simpler to implement and more robust to use. Quantitative phase imaging modalities have found numerous applications in a range of fields, including cellular mechanics, biophysics [23,24], digital pathology [25,26], X-ray crystallography [27], and optical metrology [28]. In this paper, we report a novel computational imaging method that is able to compute the complex wave field from pure intensity measurements. We name the method synthetic aperture imaging based on Kramers–Kronig relations (KKSAI). This method is closely related to pupil modulation Fourier ptychographic microscopy (FPM) [29], pupil modulation quantitative differential phase contrast (DPC) microscopy [22], and digital holographic microscopy [3032]. However, it possesses certain advantages over all these existing methods. Compared to pupil modulation FPM [29], KKSAI is non-iterative and does not require data redundancy for operation. As such, KKSAI requires less raw data as well as processing time and hence has a significantly reduced processing burden and an improved speed advantage. In comparison to pupil modulation DPC [22], KKSAI can work with a wider range of samples, as it does not require samples to adhere to the weak sample assumption (both the absorption and phase variations of the sample must be sufficiently small)—a condition that pupil modulation DPC requires for operation. Compared to in-line holography [30], KKSAI does not have a twin-image issue, as guaranteed by the band-limited signal analyticity. Compared to common-path off-axis holography [31], KKSAI’s space-bandwidth product (SBP) can be notably higher by three- to four-fold.

    This paper is structured as follows. We first present the experimental setup, data acquisition, and reconstruction algorithm of KKSAI. Next, we provide a mathematical proof of KKSAI’s validity and describe two possible scanning schemes. In the third section, we report our findings on simulations and experimental results and show that KKSAI can outperform the pupil modulation FPM and pupil modulation DPC method when given the same measured data. In the discussion section, we clarify two differences between KKSAI and the original off-axis holography paper [32]. Finally, we conclude by addressing the overall efficiency of the KKSAI method compared to other phase imaging methods.

    2. PRINCIPLE

    A. Experimental Setup

    The schematic of our KKSAI system prototype is shown in Fig. 1(a). It is simply a conventional wide-field microscope that has been modified to incorporate pupil modulation capability. As with other phase imaging methods, the objective here is to recover the phase profile of the sample from intensity image measurements.

    Principle of KKSAI. (a) Schematic of experimental setup, where pupil modulation is achieved by an SLM-based module. (b) Simplified 4f system corresponding to (a). (c) Simulated complex-valued sample. (d) Amplitude pupil modulation indicated by the green circle, whose center is (ui,vi). (e) Measured images corresponding to (d). (f) KKSAI reconstruction algorithm flowchart. It finally recovers the pupil-limited sample spectrum. MMF, multimode fiber; CL, collimating lens; RL, relay lens; P, polarizer; BS, beam splitter; SLM, spatial light modulator; TL, tube lens.

    Figure 1.Principle of KKSAI. (a) Schematic of experimental setup, where pupil modulation is achieved by an SLM-based module. (b) Simplified 4f system corresponding to (a). (c) Simulated complex-valued sample. (d) Amplitude pupil modulation indicated by the green circle, whose center is (ui,vi). (e) Measured images corresponding to (d). (f) KKSAI reconstruction algorithm flowchart. It finally recovers the pupil-limited sample spectrum. MMF, multimode fiber; CL, collimating lens; RL, relay lens; P, polarizer; BS, beam splitter; SLM, spatial light modulator; TL, tube lens.

    In the experiment, the pupil plane is relayed outside the objective (10× Mitutoyo Plan Apo infinity corrected objective, 0.28 NA) onto the spatial light modulator (SLM) so that we can perform amplitude modulation of the pupil. The modulation module consists of a reflective mode liquid crystal on silicon (Holoeye LC-R 1080) SLM and a pair of linear polarizers (P1 and P2) with their polarization directions orthogonal to each other in order to maximize the amplitude modulation contrast. The light is reflected off the SLM and finally projected onto the camera. The pixel size of the camera (Allied Vision Prosilica GX 6600) is 5.5 μm. The illumination is provided by a laser diode (Thorlabs DJ532-40) with the central wavelength of 532 nm, coupled into a multimode fiber (Thorlabs FT200emtcustom, 0.39 NA, Ø200  μm). The fiber is vibrated by a motor to wash out the speckle at the camera plane.

    The KKSAI system can be simplified and represented as a 4f imaging system as shown in Fig. 1(b). The coordinates at the sample plane, the pupil plane, and the camera sensor plane are denoted as (x,y), (u,v), and (x,y), respectively. If the complex-valued sample in Fig. 1(c) is described by s(x,y)=a(x,y)ejφ(x,y), its maximum spectrum that can be detected by this system is S(u,v)=F{s}(u,v)·C(u,v),where F{·} is the Fourier transform (FT) operator and C(u,v) is the coherent transfer function (CTF) of imaging system as indicated by the red circle in Fig. 1(d). Its radius is determined by the objective numerical aperture (NA). Our overarching objective is to recover S(u,v), which corresponds to the complex wave field at the imaging system’s pupil-limited resolution.

    B. Data Acquisition

    During KKSAI operation, we scan a binary circular aperture D(u,v) at the pupil plane, and its edge strictly crosses the pupil center, depicted by green circles in Fig. 1(d). Technically, we display the appropriate aperture on the reflective SLM to implement this pupil modulation. We use four scanning steps to fully cover S(u,v), and its offset distance from the pupil center is ρ. At each step, the scanning aperture can be denoted as D(uui,vvi), i=1,2,3,4, and its cropped subregion from sample spectrum is Si(u,v)=S(u,v)·D(uui,vvi).

    Accordingly, the measured intensity image at the camera plane can be expressed as Ii(x,y)=|F1{Si(u,v)}|2.

    These measured intensity images are then fed into the KKSAI reconstruction algorithm (explained in the next subsection). The phase recovery process during reconstruction will be repeated for all four Ii(x,y) to recover the complex field Si(u,v). They are subsequently stitched into a full pupil plane complex-valued spectrum S(u,v). In turn, we can perform an inverse FT to constitute the final amplitude and phase image of sample [see Fig. 1(f) for data flowchart].

    C. Reconstruction Algorithm

    The primary objective of the KKSAI reconstruction algorithm is to recover the complex-valued Si(u,v) from Ii(x,y).

    To best explain the process, we can first express Ii as an intensity result from superposition of an unscattered planar wave and the scattered field from the sample: Ii(x,y)=|F1{Si(u,v)}|2=|F1{Si(u,v)+δ(u,v)}|2=|si(x,y)+ej(0·x+0·y)|2,where si(x,y) is the scattered field exiting from the sample plane and ej(0·x+0·y) represents the ballistic light going through the sample. The ballistic light has a planar wavefront and propagates along the optical axis of our imaging system.

    We can see that the Fourier spectrum of Ii(x,y) is quite similar to that of the off-axis hologram as the FT of I1 in Fig. 2(a) shows. This is because the scattered field si(x,y) is offset to the planar wave ej(0·x+0·y) in the Fourier domain, and they would interfere with each other. [Here in Fig. 2, ρ is decreased so that the scanning aperture will not exceed the support of CTF and S1(u,v) covers a circular region. This ensures that the cross terms in the spectrum correspond well to the cross terms in the traditional off-axis holography case. But all the following deduction holds for any 0<ρρNAape.]

    Analogy between KKSAI measurement and off-axis hologram. (a) Measurement I1 between its corresponding complex-valued spectrum subregion and the amplitude of its FT. (b) Shifted spectrum subregion still brings in the same measurement due to the frequency shifting property of FT and phase loss of the square-law detector. (c) The shifted spectrum subregion can be hypothetically decomposed into a Dirac delta function and the shifted scattered complex-valued function.

    Figure 2.Analogy between KKSAI measurement and off-axis hologram. (a) Measurement I1 between its corresponding complex-valued spectrum subregion and the amplitude of its FT. (b) Shifted spectrum subregion still brings in the same measurement due to the frequency shifting property of FT and phase loss of the square-law detector. (c) The shifted spectrum subregion can be hypothetically decomposed into a Dirac delta function and the shifted scattered complex-valued function.

    To better explain this observation, we can employ the frequency shifting property of FT. If the covered spectrum subregion S1(u,v) is shifted to the pupil plane center, only a phase term will be multiplied on the sensor plane. However, our intensity detector is not able to detect the phase. Therefore, the same intensity image would still be measured as shown in Fig. 2(b). I1(x,y)=|F1{S1(u,v)}|2=|F1{S1(u+u1,v+v1)}|2.

    Then the shifted spectrum subregion S1(u+u1,v+v1) can be hypothetically separated into a Dirac delta function δ(u+u1,v+v1) corresponding to the original planar wave and the shifted scattered field function S˜1(u,v)=S1(u+u1,v+v1)as shown in Fig. 2(c).

    Thus, the measurement I1 can be expressed as I1(x,y)=|F1{S˜1(u,v)+δ(u+u1,v+v1)}|2=|s˜1(x,y)+ej(u1·x+v1·y)|2=|s˜1(x,y)|2+1+s˜1*(x,y)·ej(u1·x+v1·y)+s˜1(x,y)·ej(u1·x+v1·y),where s˜1(x,y) is the inverse FT of S˜1(u,v) and * denotes the complex conjugate operator. Then its FT will be F1{I1}=S˜1(u,v)S˜1(u,v)+δ(u,v)+S˜1*[(u+u1),(v+v1)]+S˜1(uu1,vv1),where represents cross-correlation operator. It can be clearly seen that the first two terms correspond to the self-interference terms in the off-axis hologram, and the other two terms come from cross interference. As such, we can see that our KKSAI measurement is analogous to an off-axis hologram.

    Next we can employ Kramers–Kronig (KK) relations to perform phase recovery. The process is similar to the one reported recently for off-axis holography [32]. As it has been established that our KKSAI measurement shares pertinent similarities to off-axis hologram, we can adapt the mathematical formulation reported in Ref. [32] to recover phase from the KKSAI measurements.

    To clearly explain the process, we use I1 as the starting point. We first generate the hypothetical reference plane wave. It is determined by the offset of S1(u,v) and in turn by the scanning aperture position (u1,v1). For I1, this reference plane wave is expressed as r1(x,y)=F1{δ(u+u1,v+v1)}=ej(u1·x+v1·y).

    Next we specify a Hilbert kernel. This kernel depends on the scanning aperture position. For I1, we can express it as H1(u,v)=jsgn(v1)·sgn(v),where sgn(v) is the sign function. The last step is based on directional Hilbert transform to retrieve the complex field corresponding to the spectrum subregion S1(u,v). We first define an intermediate variable X=ln[F1{S1(u+u1,v+v1)}/r1(x,y)]=ln{[s˜1(x,y)+r1(x,y)]/r1(x,y)}.

    Then, Re{X}=12ln[I1(x,y)/|r1(x,y)|2],Im{X}=F1{F{Re{X}}·H1(u,v)},which will be expanded on in the next subsection. Thus, S1(u+u1,v+v1)=F{eRe{X}+jIm{X}·r1(x,y)}.

    After obtaining the shifted version of four subregions, we can move them back to the correct position and get Si(u,v),i=1,2,3,4.

    Eventually, S(u,v)=i=14Si(u,v)/[i=14D(uui,vvi)+ε],where ε=105 is a small constant for numerical stability in the zero-valued region. As seen in Fig. 1(f), the reconstruction results of KKSAI are in good agreement with the original sample function.

    The validity of this KK-relations-based process has been demonstrated in the original paper [32] for off-axis holography. It leverages the analyticity of band-limited signals. To ensure that it also applies for KKSAI, we need to prove our measured image signal is indeed analytical as well.

    D. Analyticity Proof

    For off-axis holography, the self- and cross-interference terms are required to be separable in Fourier domain such that one cross term can be cropped out and inverse FT can be applied to recover the whole complex field. Such separation is usually achieved by adjusting the reference wave incidence angle. However, the recent study [32] relaxed this restriction and proved that as long as the two cross terms do not overlap, the complex field can be fully determined. We can arrive at this conclusion by using Titchmarsh theorem [33].

    For conciseness, we use l as the vector representative of (x,y) coordinate. Also, without loss of generality, the following variable replacements are used: I1(x,y)I(l)s˜1(x,y)s˜(l)r1(x,y)=ej(u1·x+v1·y)r(l)=ejρr·l,where ρr=(u1,v1). Then Eq. (7) can be expressed as I(l)=|s˜(l)+r(l)|2.

    Then I|r|2=|s˜r+1|2,ln|s˜r+1|=12[ln(I)ln(|r|2)].

    Let X=ln(s˜+rr)=ln(s˜r+1)=Re{X}+jIm{X}. We will have eX=s˜r+1=eRe{X}·ejIm{X},|s˜r+1|=eRe{X},Re{X}=ln|s˜r+1|=12[ln(I)ln(|r|2)].

    Equation (22) shows that Re{X} is a quantity determined by the measurement I(l) and the reference wave amplitude |r|. Since r(l)=ejρr·l is a planar wave, its amplitude |r|=1 is a constant. Thus, Re{X} only depends on I(l) and is known. Now, if X was analytical, KK relations can retrieve Im{X} directly from Re{X}, and the complex field s˜(l) can be reconstructed.

    To find the analyticity condition for X, we define α(l)=s˜r=s˜(l)·ejρr·l.

    Then X=ln(α+1)=n=0(1)nn+1αn+1,where the condition α(l)=s˜/r<1 is satisfied as long as the energy of the DC component is several orders of magnitude higher than the other frequency components in practice. Equation (24) shows that the analyticity of X depends on the analyticity of α.

    Next we introduce the Titchmarsh theorem to verify the analyticity of α. The theorem states that the following conditions are equivalent for a complex-valued function f(l) that is square integrable over the real l axis. The real and imaginary parts of f(l) are Hilbert transforms of each other.Its Fourier transform F{f}(ρ) is 0 or vanishes rapidly for ρ<0.

    An important observation is that s˜(l) is always a band-limited signal under the pupil modulation as shown in Figs. 3(a)–3(c). In our experiments, it is characterized by the scanning aperture radius. For the convenience of discussion, ρ=ρr/|ρr| and its perpendicular direction ρ are defined on the (u,v) plane as a new set of coordinates. Their corresponding axis pair on the (x,y) plane is l and l. If applying 1D FT to both sides of Eq. (23) along l, we have A(ρ,l)=S˜(ρ|ρr|,l),where A(ρ,l) and S˜(ρ,l) are the 1D FTs of α(l) and s˜(l) along l respectively, as shown in Fig. 3(d).

    Titchmarsh theorem applied to a band-limited signal. (a) Amplitude and (b) phase of s˜(l→) with bandwidth of ρNAape. (c) Logarithm of its 2D Fourier amplitude spectrum. (d) Logarithm of its 1D Fourier amplitude spectrum along the l∥ axis and its shifted copies by (e) |ρ→r|<ρNAape and (f) |ρ→r|=ρNAape.

    Figure 3.Titchmarsh theorem applied to a band-limited signal. (a) Amplitude and (b) phase of s˜(l) with bandwidth of ρNAape. (c) Logarithm of its 2D Fourier amplitude spectrum. (d) Logarithm of its 1D Fourier amplitude spectrum along the l axis and its shifted copies by (e) |ρr|<ρNAape and (f) |ρr|=ρNAape.

    Since s˜(l) is band limited, α(l) must be a complex-valued square integrable function along l. Hence the Titchmarsh theorem is valid for α(l). This implies that as long as A(ρ,l)=0 for ρ<0, its real and imaginary parts are Hilbert transforms of each other. Equation (25) shows that A(ρ,l) is a shifted copy of S˜(ρ,l) along ρ. As such, when the shifting distance |ρr| is at least equal to ρNAape as illustrated in Fig. 3(f), the condition is met and α(l) is analytical.

    We further note that, although the analyticity condition |ρr|ρNAape seems 1D along l, it is actually a sufficient condition of 2D analyticity on the plane l in the limit sense, which has been proven in Ref. [34] and discussed in detail in the supplementary material of Ref. [32].

    To sum up, when |ρr|ρNAape, X satisfies the analyticity condition, and as such the KK relations can be applied to recover the complex field s˜(l). However, when |ρr|>ρNAape, the DC component cannot enter the measurement, and the KKSAI’s analogy to the off-axis hologram fails. Therefore, for KKSAI to function properly, we need to operate under the condition where |ρr|=ρNAape. In other words, the scanning aperture edge must cross the pupil origin exactly.

    E. Scanning Scheme

    Figure 4(a) shows the scanning scheme discussed above. It is designed to cover the entire pupil and contains some overlaps. However, since our KKSAI method does not require data redundancy, this scanning is slightly inefficient because of the overlaps. If we stick with scanning a circular aperture to cover the entire pupil, this overlapping issue is inevitable. If we are free to choose any aperture shape, then the scanning scheme with two non-overlapping measurements shown in Fig. 4(b) would be the most efficient. The two scans have aperture edges that exactly cross the pupil origin. In principle, as long as the scanning of any convex aperture covers the whole pupil and its edge crosses the origin, the scanning scheme should work for KKSAI.

    Scanning scheme examples to cover the entire pupil. (a) Four circular apertures; (b) two rectangular apertures. The circled numbers are used to label the measurement sequence.

    Figure 4.Scanning scheme examples to cover the entire pupil. (a) Four circular apertures; (b) two rectangular apertures. The circled numbers are used to label the measurement sequence.

    If not stated, Fig. 4(a) is our default scanning scheme used in the rest of this paper. We chose this scheme because some of the other methods to which we compare KKSAI require overlaps in order to accomplish phase recovery. Finally, we note that the camera sampling rate of our method still needs to satisfy the Nyquist sampling limit 2ρNAape to avoid subsampling aliasing.

    3. SIMULATION

    In this section, we report on a series of simulations that were performed to verify the performance of our proposed method. The results are also compared with two existing imaging modalities, pupil modulation DPC (PM-DPC) [22] and pupil modulation FPM (PM-FPM) [29,35].

    To be more precise, we briefly describe these two existing phase imaging methods here. PM-DPC [22] is a counterpart of original quantitative DPC [20] and is different in that it replaces the asymmetric illumination with asymmetric pupil modulation. Specifically, it captured two or four phase gradient images from complementary half-open pupils. Then, when combined with an additional sample’s intensity image, the method can compute the complete field of the sample. PM-FPM [29,35] is the counterpart of the original FPM [16] and is different in that it replaces sample spectrum translation by oblique illumination with aperture translation on the pupil plane. Its reconstruction usually requires at least 20 raw images to converge. Although we note that there are other variants to the original DPC and FPM [36,37], we chose these two to compare with KKSAI because they all use similar experimental setups.

    All the existing complex field imaging methods have trade-offs between system complexity, sample restriction, measurement data volume, reconstruction time complexity, and reconstruction accuracy. To make our comparison meaningful, we fix the system and measurement data volume so as to examine other aspects. Thus, all the three methods in the following discussions utilize the same dataset Ii for fairness.

    Reconstruction of phase-only samples by two existing methods and our proposed method. (a) Weak phase sample; (b) strong phase sample.

    Figure 5.Reconstruction of phase-only samples by two existing methods and our proposed method. (a) Weak phase sample; (b) strong phase sample.

    Quantitative evaluations show that PM-FPM phase reconstruction is excellent in terms of similarity, but its MSE is higher than that of KKSAI. Also, the PM-FPM iterative computation is quite time consuming. On a computer (i7-7700k) with 64 GB RAM, its run time is 1 order of magnitude higher than the one for the other two non-iterative methods in MATLAB R2018b.

    The PM-DPC method provides reasonable phase rendering for the sample with small phase variations, but it fails when the phase variations are large. This is because its reconstruction algorithm demands that the sample adheres to the weak sample assumption, which is necessary to first-order Taylor expansion approximation.

    Reconstruction of complex-valued sample by two existing methods and our proposed method. (a) Phase; (b) amplitude.

    Figure 6.Reconstruction of complex-valued sample by two existing methods and our proposed method. (a) Phase; (b) amplitude.

    In this simulation study, we conducted an additional evaluation for KKSAI. As the scanning scheme used here collects redundant information in the overlaps of the spectrum [see Fig. 4(a)], we can also compare the reconstruction within the overlap regions using the FSIM metric. The results are summarized in Table 3. We can clearly see that the phase and amplitude recovered from KKSAI are highly consistent—thereby further reinforcing our confidence that KKSAI is providing us with the correct phase estimate.

    Similarity Evaluation of Overlapping Spectrum Regions in Fig. 6

    Overlapping Region
    Phase0.99400.99670.99900.9926
    Amplitude0.99970.99980.99990.9998

    To conclude, we can see that KKSAI outperforms PM-DPC and PM-FPM when we take both reconstruction accuracy and computational load into account.

    Our simulation study allows us to investigate the potential impact of experimental imperfection when KKSAI is implemented practically. Specifically, we are concerned about the way that the aperture edge crosses the pupil center. Ideally, we want the aperture edge to cross the pupil center exactly (zero overlap). If the aperture includes the pupil center (positive overlap), the analyticity condition would be broken. On the other hand, if the aperture excludes the pupil center (negative overlap), Ii(x,y) would not contain the ballistic planar wavefront, and no cross interference would occur—resulting in complete breakdown of the KKSAI processing. We simulated the impact of this overlap issue and showed our results in Fig. 7. We can clearly see that negative overlap results in a catastrophic failure of our KKSAI method, but positive overlap is somewhat tolerable as long as the overlap is small.

    Effect of distance between aperture edge and pupil center on the final reconstruction accuracy.

    Figure 7.Effect of distance between aperture edge and pupil center on the final reconstruction accuracy.

    Our final simulation is focused on the feasibility of the scanning scheme in Fig. 4(b). Here we use an aperture that alternates between two halves of the pupil as shown in Fig. 8(a). We can see from the results reported in Fig. 8(b) that the phase and amplitude recovery result is highly consistent with the ground truth in Fig. 6 as indicated by the quantitative metrics.

    KKSAI based on the scanning scheme with only two measurements. (a) Measurements; (b) reconstruction.

    Figure 8.KKSAI based on the scanning scheme with only two measurements. (a) Measurements; (b) reconstruction.

    4. EXPERIMENT

    With the encouraging simulation results, we next conducted a series of experiments with the setup shown in Fig. 1(a).

    In the first experiment, a plano-convex microlens array with 150 μm pitch (Thorlabs MLA150-7AR-M) was used. It can be regarded as a pure phase test sample. Under the illumination wavelength of 532 nm, four pupil modulated images were acquired and processed by the KKSAI reconstruction algorithm. The reconstructed phase map of a single microlens by KKSAI is shown in Fig. 9(a), together with the reconstructions from the PM-DPC and PM-FPM algorithms using the same measurements. We can clearly see that PM-DPC and KKSAI achieve better recovery of the microlens shape, providing the scanning electron microscope (SEM) image in Fig. 9(b) as reference. To quantify their accuracy, we plotted the radially averaged cross section profile from each result, labeled by green, blue, and red lines in Fig. 9(c). According to the spec sheet of the microlens array, the lens nominal profile is computed and labeled as a black line here. It can be seen that the phase profiles obtained by KKSAI and PM-DPC are in good agreement with the ground truth. The FPM algorithm suffered from the low overlapping rate in the Fourier domain and insufficient redundancy to be noise robust, thus resulting in a poor phase recovery.

    Experimental results for a microlens array. (a) Reconstructions of a single microlens by PM-DPC, PM-FPM, and KKSAI using four measurements. (b) A close-up view of the SEM image of the microlens array (adapted from Thorlabs website). (c) Radial average profile of three phase recoveries compared with the ground truth (GT).

    Figure 9.Experimental results for a microlens array. (a) Reconstructions of a single microlens by PM-DPC, PM-FPM, and KKSAI using four measurements. (b) A close-up view of the SEM image of the microlens array (adapted from Thorlabs website). (c) Radial average profile of three phase recoveries compared with the ground truth (GT).

    Next we imaged a papillary thyroid carcinoma pap smear slide with papanicolaou stain. It is a complex-valued sample. Figure 10(a) displays two raw images acquired in the experiment. Due to the pupil amplitude modulation, the shade effect can be seen in the measurements, and it is clear that their spectrum contains two cross-interference terms similar to those of an off-axis hologram, just as analyzed in the principle section. SLM allows us to carefully align the aperture, such that these two cross terms are tangential to each other. Feeding the data into the KKSAI algorithm, we can then perform the complex field reconstruction, and we present the results in Figs. 10(b) and 10(c).

    Experimental results for a thyroid carcinoma pap smear slide. (a) Two out of four measurements acquired by KKSAI and their Fourier amplitude spectrum. (b) Amplitude reconstruction by PM-FPM and KKSAI compared with ground truth. (c) Phase reconstruction by PM-DPC, PM-FPM, and KKSAI compared with ground truth.

    Figure 10.Experimental results for a thyroid carcinoma pap smear slide. (a) Two out of four measurements acquired by KKSAI and their Fourier amplitude spectrum. (b) Amplitude reconstruction by PM-FPM and KKSAI compared with ground truth. (c) Phase reconstruction by PM-DPC, PM-FPM, and KKSAI compared with ground truth.

    To obtain the ground truth for this complex valued sample, we perform a separate PM-FPM experiment where 47 pupil modulation images with an overlapping rate of about 85% in the Fourier domain were acquired. Its high-resolution reconstruction is taken as the ground truth to evaluate the performance of the three methods. As we can see in Figs. 10(b) and 10(c), PM-DPC cannot recover amplitude, and PM-FPM reconstructions based on four measurements are of poor quality. By directly viewing the images, we can see that both amplitude and phase recovery of KKSAI are the closest to ground truth. To quantify the results, FSIM metrics are calculated and labeled under each reconstruction in red. From the quantifications, we can see that KKSAI clearly outperforms the other two methods.

    Similar to the simulation, we can further evaluate KKSAI’s performance by checking the consistency between reconstructions within the overlapping spectrum region from different scanning apertures. Here the FSIM metric is calculated and summarized in Table 4. As we can see, the recovery from different measurements matches very well in the overlapping regions.

    Similarity Evaluation of Overlapping Spectrum Regions in Fig. 10

    Overlapping Region
    Phase0.98720.99270.98640.9691
    Amplitude0.99910.99980.99860.9991

    Since the proposed KKSAI method is able to reconstruct the whole complex light field, we can conduct many digital postprocessings, for example, digital refocusing. So here we show one possible application based on digital refocusing, axial chromatic aberration correction of an imaging system. To demonstrate this, we image the same pap smear slide under the illumination wavelength of 405 nm (Thorlabs DL5146-101S), 532 nm, and 638 nm (Thorlabs L638P040) sequentially. As shown in Fig. 11(d), the chromatic focus shift causes rainbow fringes along the boundaries in the color composite image. It implies that there are still some residual chromatic aberrations associated with the objective. With the phase retrieved by KKSAI in Fig. 11(b), we can utilize the angular spectrum method [1] to digitally propagate the complex wave field in individual channels and an autofocusing metric based on edge sparsity criterion [39] to find the best focus location. Finally, we can obtain the optimal refocusing distance [labeled under each channel in Fig. 11(c)]. Comparing color composite image before and after digital refocusing in Fig. 11(d), we can clearly see that chromatic aberration artefacts are well corrected.

    Chromatic aberration correction by digital refocusing ability of KKSAI. (a) Reconstructed amplitudes by KKSAI of three channels. (b) Reconstructed phases by KKSAI of three channels. (c) Digitally refocused amplitudes with the corresponding refocusing distance labeled at the bottom. (d) Color composite of three channels before and after digital refocusing with the enlargements showing improved image quality. R, red (638 nm); G, green (532 nm); B, blue (405 nm).

    Figure 11.Chromatic aberration correction by digital refocusing ability of KKSAI. (a) Reconstructed amplitudes by KKSAI of three channels. (b) Reconstructed phases by KKSAI of three channels. (c) Digitally refocused amplitudes with the corresponding refocusing distance labeled at the bottom. (d) Color composite of three channels before and after digital refocusing with the enlargements showing improved image quality. R, red (638 nm); G, green (532 nm); B, blue (405 nm).

    Finally, we conducted an experiment to prove the feasibility of reconstructing a complex wave field from only two measurements by KKSAI. The scanning scheme follows Fig. 4(b). As shown in Fig. 12(a), the spectrum amplitude of two raw images is in agreement with the simulation in Fig. 8(a), although here the two cross terms are much stronger than the autocorrelation term, such that the autocorrelation term is barely discernible in the spectrum.

    Complex wave-field reconstruction by KKSAI based on only two measurements. (a) Scanning scheme, raw measurements, and their spectrum amplitude; (b) and (c) reconstructed amplitude and phase, respectively, by KKSAI from two measurements, four measurements, and PM-FPM with 47 measurements. Here PM-FPM reconstruction is taken as the reference to calculate the FSIM metric for KKSAI reconstruction.

    Figure 12.Complex wave-field reconstruction by KKSAI based on only two measurements. (a) Scanning scheme, raw measurements, and their spectrum amplitude; (b) and (c) reconstructed amplitude and phase, respectively, by KKSAI from two measurements, four measurements, and PM-FPM with 47 measurements. Here PM-FPM reconstruction is taken as the reference to calculate the FSIM metric for KKSAI reconstruction.

    We can perform KKSAI reconstructions based on these two measurements. As shown in Figs. 12(b) and 12(c), its result is compared to the KKSAI result based on four measurements of the scanning scheme in Fig. 4(a) and the aforementioned PM-FPM retrieval with 47 measurements. Taking the PM-FPM result as the reference, the FSIM metric for two KKSAI reconstructions is computed and labeled in red at the bottom of each image. We can see that KKSAI reconstructions with only two measurements show high conformity with the PM-FPM ground truth. To conclude, our proposed KKSAI method can achieve complex wave-field reconstruction with almost the same quality as PM-FPM based on 47 measurements. Furthermore, it is iteration free. As such, it is an attractive potential replacement for the current PM-FPM method. To be clear, PM-FPM is distinct from the standard FPM [16]. PM-FPM does not have the capability to increase resolution beyond the limit of the objective NA, which is a key feature of the standard FPM [16]. By extension, KKSAI with the current setup is not capable of increasing resolution beyond the objective NA.

    5. DISCUSSION

    In terms of reconstruction accuracy, KKSAI should be similar to Ref. [32] because they share nearly identical mathematics. However, our KKSAI method is distinct from off-axis holographic microscopy [32] as it does not actually need a real reference arm. Compared to Ref. [32], our system is simpler to implement. Also, SLM is not necessary for KKSAI. The SLM-based pupil modulation module can be replaced by a simple physical iris mounted and controlled by the motorized stage, once the scanning scheme has been pinned down. By its simplicity, KKSAI is more robust to operate in practical situations. Besides, there actually exists a subtle difference between these two methods.

    Mathematically, the off-axis holography experiment can be interpreted as the addition of a band-limited sample spectrum and a delta function in the frequency domain. Their relative offset will decide the distance between two cross terms. In both KKSAI and off-axis holography [32], the delta function lies exactly on the edge of the sample spectrum. But a key difference between the two methods is that for KKSAI there is a global offset for both the sample spectrum and delta function such that the delta function happens to be located at the origin and is exactly the DC value. This explains the absence of interference fringes in the KKSAI measurements, unlike in off-axis holograms.

    It is also worth mentioning that our pupil modulation mode can be easily adapted into its reciprocal mode—tilted illumination [40]. The multiplication with a phase ramp in the spatial domain by tilted illumination is equivalent to the offset modulation in the frequency domain. However, the advantage of our pupil modulation mode is that the thin sample requirement is circumvented because it only focuses on the exit wave from sample.

    An interesting point is that under tilted illumination mode, our method can be realized by only lighting up LED elements that are located on a ring with the illumination NA matching the objective NA. This is known as an annular illumination scheme. It has been shown that only under this particular condition low-frequency phase components can be transferred into intensity measurements [37]. Although this conclusion is derived from the weak sample assumption, its coincidence with our findings indicates that matching illumination and objective NA or keeping the DC component just at the modulation aperture edge is an excellent operating point for computational imaging. This observation may inspire other new computational imaging methods. For example, the recovery of any coherent signal with a strong DC peak may benefit from the KKSAI concept.

    6. CONCLUSION

    In this paper, a computational imaging method for reconstructing the complex wave field, KKSAI, is reported. Its experimental setup, data acquisition, and reconstruction algorithm are described in detail. By establishing the analogy between KKSAI measurement and off-axis hologram, a recent advance for off-axis holography based on KK relations [32] is adapted here to recover phase from intensity images in a non-interferometric way. The method leverages the analyticity of band-limited signals under pupil modulation to compute phase information from intensity measurements.

    As a computational imaging modality, the KKSAI method codesigns the sensing part and the reconstruction algorithm to optimally operate together. From the perspective of sensing, it requires much fewer measurements than pupil modulation FPM. From the perspective of reconstruction, it is iteration free and does not require any sample priors, and it is thus more generally usable than PM-DPC. Our simulation and experiment results demonstrate that KKSAI has clear advantages over the other two methods.

    Acknowledgment

    Acknowledgment. Cheng Shen thanks Ruizhi Cao and Dr. Baptiste Blochet for helpful discussions on this work.

    References

    [1] J. W. Goodman. Introduction to Fourier Optics(2005).

    [2] P. Marquet, B. Rappaz, P. J. Magistretti, E. Cuche, Y. Emery, T. Colomb, C. Depeursinge. Digital holographic microscopy: a noninvasive contrast imaging technique allowing quantitative visualization of living cells with subwavelength axial accuracy. Opt. Lett., 30, 468-470(2005).

    [3] C. S. Seelamantula, N. Pavillon, C. Depeursinge, M. Unser. Exact complex-wave reconstruction in digital holography. J. Opt. Soc. Am. A, 28, 983-992(2011).

    [4] G. Popescu, T. Ikeda, R. R. Dasari, M. S. Feld. Diffraction phase microscopy for quantifying cell structure and dynamics. Opt. Lett., 31, 775-777(2006).

    [5] Z. Wang, L. Millet, M. Mir, H. Ding, S. Unarunotai, J. Rogers, M. U. Gillette, G. Popescu. Spatial light interference microscopy (SLIM). Opt. Express, 19, 1016-1026(2011).

    [6] A. Anand, V. Chhaniwal, B. Javidi. Tutorial: common path self-referencing digital holographic microscopy. APL Photon., 3, 071101(2018).

    [7] N. T. Shaked, Y. Zhu, M. T. Rinehart, A. Wax. Two-step-only phase-shifting interferometry with optimized detector bandwidth for microscopy of live cells. Opt. Express, 17, 15585-15591(2009).

    [8] P. Gao, B. Yao, J. Min, R. Guo, J. Zheng, T. Ye, I. Harder, V. Nercissian, K. Mantel. Parallel two-step phase-shifting point-diffraction interferometry for microscopy based on a pair of cube beamsplitters. Opt. Express, 19, 1930-1935(2011).

    [9] C. Joo, T. Akkin, B. Cense, B. H. Park, J. F. De Boer. Spectral domain optical coherence phase microscopy for quantitative phase contrast imaging. Opt. Lett., 30, 2131-2133(2005).

    [10] Z. Yaqoob, W. Choi, S. Oh, N. Lue, Y. Park, C. Fang-Yen, R. R. Dasari, K. Badizadegan, M. S. Feld. Improved phase sensitivity in spectral domain phase microscopy using line-field illumination and self phase referencing. Opt. Express, 17, 10681-10687(2009).

    [11] R. W. Gerchberg. A practical algorithm for the determination of phase from image and diffraction plane pictures. Optik, 35, 237-246(1972).

    [12] J. R. Fienup. Phase retrieval algorithms: a comparison. Appl. Opt., 21, 2758-2769(1982).

    [13] C. Shen, X. Bao, J. Tan, S. Liu, Z. Liu. Two noise-robust axial scanning multi-image phase retrieval algorithms based on Pauta criterion and smoothness constraint. Opt. Express, 25, 16235-16249(2017).

    [14] J. M. Rodenburg, H. M. Faulkner. A phase retrieval algorithm for shifting illumination. Appl. Phys. Lett., 85, 4795-4797(2004).

    [15] A. Pan, B. Yao. Three-dimensional space optimization for near-field ptychography. Opt. Express, 27, 5433-5446(2019).

    [16] G. Zheng, R. Horstmeyer, C. Yang. Wide-field, high-resolution Fourier ptychographic microscopy. Nat. Photonics, 7, 739-745(2013).

    [17] A. Pan, Y. Zhang, K. Wen, M. Zhou, J. Min, M. Lei, B. Yao. Subwavelength resolution Fourier ptychography with hemispherical digital condensers. Opt. Express, 26, 23119-23131(2018).

    [18] L. Waller, L. Tian, G. Barbastathis. Transport of intensity phase amplitude imaging with higher order intensity derivatives. Opt. Express, 18, 12552-12561(2010).

    [19] C. Zuo, Q. Chen, Y. Yu, A. Asundi. Transport-of-intensity phase imaging using savitzky-golay differentiation filter-theory and applications. Opt. Express, 21, 5346-5362(2013).

    [20] S. B. Mehta, C. J. Sheppard. Quantitative phase-gradient imaging at high resolution with asymmetric illumination-based differential phase contrast. Opt. Lett., 34, 1924-1926(2009).

    [21] L. Tian, L. Waller. Quantitative differential phase contrast imaging in an LED array microscope. Opt. Express, 23, 11394-11403(2015).

    [22] H. Lu, J. Chung, X. Ou, C. Yang. Quantitative phase imaging and complex field reconstruction by pupil modulation differential phase contrast. Opt. Express, 24, 25345-25361(2016).

    [23] Y. Park, M. Diez-Silva, G. Popescu, G. Lykotrafitis, W. Choi, M. S. Feld, S. Suresh. Refractive index maps and membrane dynamics of human red blood cells parasitized by Plasmodium falciparum. Proc. Natl. Acad. Sci. USA, 105, 13730-13735(2008).

    [24] W. J. Eldridge, S. Ceballos, H. S. Park, A. Wax. Comparing quantitative phase derived cellular mechanical parameters with atomic force microscopy measurements. Proc. SPIE, 10888, 1088803(2019).

    [25] A. Greenbaum, Y. Zhang, A. Feizi, P.-L. Chung, W. Luo, S. R. Kandukuri, A. Ozcan. Wide-field computational imaging of pathology slides using lens-free on-chip microscopy. Sci. Transl. Med., 6, 267ra175(2014).

    [26] R. Horstmeyer, X. Ou, G. Zheng, P. Willems, C. Yang. Digital pathology with Fourier ptychography. Comput. Med. Imaging Graph., 42, 38-43(2015).

    [27] J. Miao, D. Sayre, H. Chapman. Phase retrieval from the magnitude of the Fourier transforms of nonperiodic objects. J. Opt. Soc. Am. A, 15, 1662-1669(1998).

    [28] A. Anand, V. K. Chhaniwal, P. Almoro, G. Pedrini, W. Osten. Shape and deformation measurements of 3D objects using volume speckle field and phase retrieval. Opt. Lett., 34, 1522-1524(2009).

    [29] X. Ou, J. Chung, R. Horstmeyer, C. Yang. Aperture scanning Fourier ptychographic microscopy. Biomed. Opt. Express, 7, 3140-3150(2016).

    [30] W. Zhang, L. Cao, D. J. Brady, H. Zhang, J. Cang, H. Zhang, G. Jin. Twin-image-free holography: a compressive sensing approach. Phys. Rev. Lett., 121, 093902(2018).

    [31] C. Zheng, R. Zhou, C. Kuang, G. Zhao, Z. Yaqoob, P. T. So. Digital micromirror device-based common-path quantitative phase imaging. Opt. Lett., 42, 1448-1451(2017).

    [32] Y. Baek, K. Lee, S. Shin, Y. Park. Kramers–Kronig holographic imaging for high-space-bandwidth product. Optica, 6, 45-51(2019).

    [33] E. C. Titchmarsh. Introduction to the Theory of Fourier Integrals(1948).

    [34] J. P. Havlicek, J. W. Havlicek, A. C. Bovik. The analytic image. Proceedings of International Conference on Image Processing, 2, 446-449(1997).

    [35] S. Dong, R. Horstmeyer, R. Shiradkar, K. Guo, X. Ou, Z. Bian, H. Xin, G. Zheng. Aperture-scanning Fourier ptychography for 3D refocusing and super-resolution macroscopic imaging. Opt. Express, 22, 13586-13599(2014).

    [36] M. Chen, Z. F. Phillips, L. Waller. Quantitative differential phase contrast (DPC) microscopy with computational aberration correction. Opt. Express, 26, 32888-32899(2018).

    [37] J. Sun, C. Zuo, J. Zhang, Y. Fan, Q. Chen. High-speed Fourier ptychographic microscopy based on programmable annular illuminations. Sci. Rep., 8, 7669(2018).

    [38] L. Zhang, L. Zhang, X. Mou, D. Zhang. FSIM: a feature similarity index for image quality assessment. IEEE Trans. Image Process., 20, 2378-2386(2011).

    [39] Y. Zhang, H. Wang, Y. Wu, M. Tamamitsu, A. Ozcan. Edge sparsity criterion for robust holographic autofocusing. Opt. Lett., 42, 3824-3827(2017).

    [40] Y. Baek, Y. Park. Intensity-based holographic imaging via space-domain Kramers–Kronig relations. Nat. Photonics, 15, 354-360(2021).

    Cheng Shen, Mingshu Liang, An Pan, Changhuei Yang. Non-iterative complex wave-field reconstruction based on Kramers–Kronig relations[J]. Photonics Research, 2021, 9(6): 1003
    Download Citation