1.

Introduction

Fourier ptychographic microscopy (FPM) has been proposed, developed and widely applied since 2013.^1,2 FPM has significant advantages in achieving wide-field and high-resolution imaging.^3–8 In combination with other imaging modalities, it can additionally achieve 3D, high-speed, or subwavelength imaging.^9–13 The key concept in FPM is that synthetic aperture and phase retrieval can be combined to increase the space-bandwidth product of the optical system. In FPM, we apply the idea of synthetic aperture at the Fourier imaging plane to combine the angular spectrum imaging information collected from different oblique angle illuminations. Phase retrieval is used to facilitate the data stitching process by backfilling the missing phase information with the help of data redundancy in the collected image data. Put in a different light, we note that phase information is not directly collected through the FPM imaging system; only intensity data are acquired during the measurement process. The phase information is instead recovered by leveraging redundancies in the intensity measurements with the help of an iterative phase retrieval loop.^1,14–17 With this phase information, we can then combine the raw data and constitute an angular spectrum that is significantly larger than what the collection optics numerical aperture conventionally dictates. Phase retrieval can be applied in FPM imaging systems, despite the fact that the light sources used are incoherent LEDs, because the image patch size used in FPM reconstruction is comparable to or smaller than the spatial coherence length of the illumination. Interested readers are invited to refer to Ref. 18.

A significant feature of the FPM is that the iterative phase retrieval loop can incorporate an aberration correction subroutine to compensate for and eliminate aberrations in the final processed image.^1,3 Among all aberration terms, optical defocus is of most interest and significance. Defocus is also unique among all aberration terms in that it is not attributable to an intrinsic physical deficiency or imperfection in the optical imaging system, unlike coma or astigmatism. An otherwise optically “perfect” system can still present defocus aberration if the sample of interest is displaced from the optical system’s focal plane. FPM can “digitally refocus” images by numerically zeroing out the defocus aberration. Practically, we do this by multiplying a defocus phase term onto the combined pupil function and the sample angular spectrum at each iteration of the FPM reconstruction. This operation allows us to eliminate the defocus aberration in each pupil function at each illumination angle. The phase retrieval process can then proceed and converge on an image solution that is digitally in-focus at the sample plane, even if the sample plane is substantially distant from the system’s actual focal plane. Practically, the distance between the sample plane and the system’s focal plane can be several times the system’s native depth of field (DOF).^1,19 This ability to “digitally refocus” the image is highly attractive in digital pathology and other biomedical microscopy applications.

To date, the “digital refocusing” process has invariably been performed by nesting the aberration correction procedure within the phase-retrieval iteration loop (standard in-iteration digital refocusing). This implies that, to alter the defocus distance (to place the digital focal plane at a different displacement), we need to adjust the aberration correction term and rerun the entire raw dataset through the FPM phase retrieval iteration loop. Because such a rerun is computationally intensive, there is a strong impetus to move the aberration correction out of the phase-retrieval loop. If feasible, this will allow users to significantly reduce the FPM computational load and provide a quick and efficient way to digitally refocus microscopy images on the fly. Recently, Refs. 2021.–22 reported a method of using the linear relationship between defocus distance and image lateral shifts under oblique illuminations to find the amount of defocus. Then, the method mechanically moves the sample to the in-focus position with the calculated defocus value or nests the correct defocus distance in the iteration of FPM optimization. This method solves the focusing issue by finding the correct defocus distance at first and then performing FPM with the correct defocus distance value in the phase-retrieval iteration loop – this approach is different from our current analysis. Our study examines the feasibility of bypassing the iterative optimization process using numerical propagation.

In optical terms, the recovered image solution from a properly executed FPM process (where the aberration correction procedure is embedded within the phase-retrieval iteration loop) can be interpreted as a wavefront emerging from the output surface of the sample. The propagation of such a wavefront along the optical axis can be numerically predicted by simply propagating the wavefront digitally. This mapping of the FPM-processed data into a physically interpretable framework suggests that perhaps postreconstruction digital refocusing can proceed in the following fashion.

This seemingly logical process has surfaced and been discussed on numerous occasions within the computational microscopy community. Yet until recently, no peer-reviewed paper has experimentally explored this avenue in any meaningful way. A recent paper²³ reported experimental findings using a combination of this approach and image processing techniques for a focal range that is less than the DOF of the employed objective. It postulated that an FPM image has a natural DOF that is significantly smaller than the objective’s DOF. As such, applying a reverse defocus would allow an FPM image to be digitally refocused within the objective’s range. More broadly, it is interesting to ask whether postreconstruction digital refocusing can work for any defocus distance – whether within or beyond the objective’s DOF.

This study aims to discuss a subtle but fatal issue in this processing approach. Specifically, without the aberration correction, the FPM raw dataset is actually illconditioned for FPM phase-retrieval iteration in Step 1. Although the algorithm will converge onto a local minimum, the output solution cannot be reasonably assumed to be the wavefront at the system’s focal plane. As such, digitally refocusing the data using Step 2 simply cannot work. This failure can be most clearly understood for the regime in which the defocus distance is large. We performed experiments to verify that the output solution from Step 1 is illmatched to the input raw data, indicating that the solution does not represent the wavefront at the system’s focal plane. We experimentally compared the digital refocusing performance of the method with the standard in-iteration digital refocusing approach and clearly demonstrated the difference.

More interestingly, we examined and demonstrated that postreconstruction digital refocusing is also not feasible in the regime in which the defocus distance is small. In this regime, the FPM phase-retrieval iteration in Step 1 can converge to a solution that approximately resembles the correct image solution – as long as the local defocus phase terms are not severely large. However, this solution lacks the global defocus phase term for which postreconstruction digital refocusing purports to correct. As such, performing postreconstruction digital refocusing is also meaningless in this regime. (The demarcation between the two regimes is arbitrary and largely depends on the user’s tolerance for image quality. We discuss this in detail in Sec. 2.)

Postreconstruction digital refocusing strategies, such as the one proposed in Ref. 23, pivotally assume the existence of a global defocus phase term in the reconstructed FPM image. However, to our knowledge, no research has mathematically or experimentally verified the validity of this assumption. Our study showed that this is a false assumption, and postreconstruction digital refocusing methods that depend on it cannot be expected to work as intended.

In Sec. 2, we first fully describe the data processing steps required for (1) digital refocusing with the focus correction nested in an iterative loop (the standard in-iteration digital refocusing approach) and (2) postreconstruction digital refocusing through numerical wave propagation. We then explain why, for large defocus distances, the postreconstruction digital refocusing method fails to condition the raw data so that it can be processed meaningfully in a phase retrieval algorithm. We also explain why, for small defocus distances, postreconstruction digital refocusing is meaningless because a global defocus phase term is absent from the FPM algorithm output. In totality, this means that postreconstruction digital refocusing does not serve a meaningful purpose for any defocus distance. In Sec. 3, we provide experimental verifications of these considerations. In Sec. 4, we use the mathematical framework from Sec. 2 to suggest a few alternate strategies for digital refocusing of the FPM data. Although each of these methods has its respective disadvantages, they can be useful in specific applications in which the ability to perform some level of dynamic digital refocusing is desired.

2.

Principle

We begin this analysis by first listing the two conditions that underpin our analysis. First, we assume that the sample is sufficiently thin compared with the DOF of the imaging system. Thus, we can treat the sample as a 2D transmissive phase and amplitude object. Upon normal incidence illumination, the transmitted wavefront at the exit surface of the sample is expressed as [Fig. 1(a)]

Fig. 1

A schematic diagram showing the difference in the transmitted wavefront reaching the image plane when the sample is in focus and out of focus. This graphical exercise is aimed at familiarizing the readers with the mathematical functions involved in the imaging process, as well as how defocus is mathematically modeled. The mathematical notations follow the same nomenclature used in the main text. The sample is placed at (a) the focal plane $z=0$ and (b) the defocused plane $z=z_0$. $\mathrm{FT}/\mathrm{IFT}$ indicates (inverse) Fourier transform. $|·{|}^2$ is the amplitude square of the image, which gives the light field intensity distribution.

Second, we assume that the optical imaging systems considered in this analysis are free of all aberrations (optically perfect) except for defocus. This simplification allows us to better see the impact of physical defocus on FPM’s operations. If we can successfully show that postreconstruction digital refocusing fails for an optically perfect system, then by logical extension, we can expect it to fail for more realistic optical systems as well.

Now, consider a simple transmission imaging scenario, as shown in Fig. 1(b), where the sample is placed at a distance $z_0$ from the focal plane of the imaging system. Suppose the sample is illuminated at normal incidence. The wavefront that arrives at the imaging plane (where a camera is placed) is then expressed as the following:

Knowledge of $I_{z_0}(x,y)$ does not provide us with enough information to determine $O_{z_0}(k_x,k_y)$ [Eq. (4)]. However, suppose we gather more information to supplement the knowledge of $I_{z_0}(x,y)$ and reconstruct $O_{z_0}(k_x,k_y)$. Then, digital refocusing can be easily accomplished by simply multiplying $O_{z_0}(k_x,k_y)$ with $H(k_x,k_y;-z_0)$ to recover the original in-focus spectrum $O_0(k_x,k_y)$ (spectrum at $z_0=0$), which is effectively a copy of $O(k_x,k_y)$ truncated by the accepting numerical aperture of the imaging system [Eq. (2)], and this can in turn be used to generate an in-focus image. Put in a slightly different context, we can also say that knowledge of $I_{z_0}(x,y)$ does not allow us to digitally refocus an image, but if we have knowledge of $O_{z_0}(k_x,k_y)$, we can then perform digital refocusing.

From this exercise, we can see that knowledge of ${\mathrm{O}}_{z_0}(k_x,k_y)$ is very useful. With it, a simple multiplication of a suitable radial defocus function will refocus the data to the right focal plane. The primary objective of a postreconstruction digital refocusing method is to derive $O_{z_0}(k_x,k_y)$, so postprocessing digital refocusing can be easily done.

To analyze how this might be accomplished, we first turn our attention to the in-focus case in which the defocus distance of the sample is zero, i.e., $z_0=0$. In this situation, $h$ and $H$ are simplified to a constant ($=1$). ${\mathrm{O}}_{z_0}(k_x,k_y)$ at $z_0=0$ is denoted by $O_0(k_x,k_y)$.

In FPM, we successively vary the incidence angle of illumination and capture the intensity image data associated with each incidence angle (Fig. 2). For the incidence angle ${\beta }_i$, corresponding to the $i$’th LED $(k_{x_i},k_{y_i})$ in Fourier space, the illumination will translate the transmission angular spectrum by the corresponding Fourier space distances $k_{x_i}$ and $k_{y_i}$. The transmitted spectrum is expressed as

Fig. 2

The space and angular spectrum representations of FPM reconstruction for the sample placed at the focal plane. During iteration, the spectrum of each aperture is updated in the iteration with an intensity constraint (blue and yellow double-headed solid arrows). Then, all of the apertures’ spectra are repeatedly updated until they converge to a global minimum solution (green dashed circles). If done correctly, the IFT of each aperture should match the raw FPM images (blue and yellow double headed dashed arrows). The final reconstruction can be obtained by ${|\mathrm{IFT}|}^2$.

The spectrum received at the imaging plane will be the same as the spectrum transmitted, except once again it would be truncated by the numerical aperture of the imaging system

Thus, this allows us to capture a different segment of the sample’s transmitted wavefront $O(k_x,k_y)$. The camera-acquired image is then given as

From Eq. (9), we further note that we are again receiving incomplete information because our camera can only measure intensity information.

Although we receive incomplete information with each raw image frame, the process of capturing intensity image data at varying incidence angles allows us to effectively fill in our knowledge of $O(k_x,k_y)$. In FPM, the information segments are designed to overlap. Thus, we can use the information redundancy of the overlapping data in a phase retrieval algorithm to deduce $O(k_x,k_y)$. In very broad terms, the phase retrieval algorithm makes a guess at $O(k_x,k_y)$ and then verifies to see if each segment of the $O(k_x,k_y)$ guess will give the correct $O_{0i}(k_x,k_y)$ prediction that matches with the $i$’th raw intensity image [Eq. (9)] collected by the FPM system. The process is then iterated to minimize this mismatch across all raw image data. Implicit in FPM’s ability to do such processing is the assumption that there is a contiguous $O(k_x,k_y)$ of which the raw FPM intensity images’ spectra are a part. We will later see that this assumption does not hold for the postreconstruction digital refocusing method.

We now consider imaging of a sample that is placed at a distance $z_0$ from the focal plane in an FPM system (Fig. 3).

Fig. 3

The space and angular spectrum representations of FPM reconstruction for the sample placed at the defocused plane. Due to the presence of the defocus phase term, the optimization result will not converge to a meaningful global minimum. Tellingly, postconvergence, the IFT of each aperture spectrum will not match up with the raw FPM images. Thus, the application of the postreconstruction digital refocusing will not lead to a meaningful in-focus solution (the last row in the figure).

Here, the transmission spectrum will once again be tilted [same expression as Eq. (8)]. Due to the defocus, the spectrum received at the imaging plane will be modified by the radial defocus phase term

Similar to the previous scenario, we can capture a sequence of intensity image data at varying incidence angles. If we follow the same logic, we will then apply the phase retrieval algorithm to recover $O(k_x,k_y;z_0)$, $O(k_x,k_y;z_0)$ and $O_{z_0}(k_x,k_y)$ are different. $O(k_x,k_y;z_0)=a_{0}T(k_x,k_y)H(k_x,k_y;z_0)$ is the spectrum of sample at a defocus distance of $z_0$ without any aperture constraint. Put it in another way, $O_{z_0}(k_x,k_y)=O(k_x,k_y)A(k_x,k_y)$, where $O_{z_0}(k_x,k_y)$ is described in Eq. (4). Unfortunately, this is where the correspondence fails. Although $O_{z_0}(k_x,k_y)$ in Eq. (4) contains a radial defocus phase term centered at the origin of $k_x$, $k_y$ space, the collected raw data does not have such a phase term. Instead, each raw image comes from $O_{z_{0}i}(k_x,k_y)$, where the radial defocus phase term is centered at the middle of each raw image angular spectrum [the middle point being the coordinate point ($k_{x_i}$, $k_{y_i}$)]. In other words, rather than sharing a common radial defocus phase term, each raw image has a different radial defocus phase term. Therefore, the raw FPM intensity images are not part of a contiguous $O(k_x,k_y;z_0)$, and as such, we cannot apply phase retrieval to yield any meaningful estimate of $O(k_x,k_y;z_0)$ (see the green dashed box in Fig. 3 for a graphical representation of this problem). Another way to interpret this in-conciliable breakage is to turn our attention to a point in the spectrum space where data from two raw intensity images overlap. In the original FPM scheme, we can expect the “ground-truth” value of that point to be the same. For this defocus situation, this value would be quite different from the two raw image datasets because they are each modified by a different radial defocus function. Running the phase retrieval algorithm to converge to the common “ground-truth” value is, therefore, futile because there is no common value to converge to.

This difficulty is averted in standard in-iteration digital refocusing. In that processing, the algorithm applies a radial defocus phase term to compensate for the defocus associated with each raw data spectrum segment, as depicted in Fig. 4. The convergence target, in this case, is $O(k_x,k_y)$, rather than $O(k_x,k_y;z_0)$. Once we have compensated for the defocus associated with each raw data spectrum segment, all raw data share the same “ground-truth” – they are all part of $O(k_x,k_y)$. The downside of this method is that, before execution, the defocus must be estimated and input into the algorithm. To retune for another defocus distance, the whole algorithm must be re-executed.

Fig. 4

The standard in-iteration digital refocusing of the FPM diagram. In each iteration, the defocus phase term is corrected for (blue and yellow solid box). That is to say, the sample spectrum is reconstructed at the focal plane. Here, the phase retrieval process will converge to a common “ground-truth” solution. Subsequently, the IFT of each aperture will match the raw FPM images. Finally, the in-focus image can be reconstructed.

The regime in which the defocus distance is small is worth a closer examination. In this regime, the individual radial defocus phase term associated with each raw image is not as drastically large as in the regime in which the defocus distance is substantial. In this situation, though the phase mismatches between the different raw intensity images still exist, the mismatch may be small enough for the FPM algorithm to converge onto something that approximates the “ground-truth.” Exactly how small a defocus distance can be tolerated is very much dependent on the user’s choice of performance benchmarks. Such benchmarks can be based on the structural similarity index measure (SSIM), resolution, or even the expected defocus-associated phase mismatch between each raw data spectrum segment. In this study, we chose to use both SSIM and resolution. To be clear, this is benchmarking the FPM algorithm’s performance for a practically usable but fundamentally wrong local minimum solution because, as long as there is a non-zero defocus distance, we can no longer expect the FPM algorithm to converge onto the correct “ground-truth” solution. In other words, the demarcation between the two regimes is arbitrary – the FPM solution is only an approximation to the “ground-truth” solution, and it is only the degree of error tolerance by which we can accept the solution that sets the demarcation line.

In this regime, when the defocus distance is small, each raw data spectrum segment carries its own radial defocus phase term, which is centered at the middle of each raw image angular spectrum – a radial defocus phase term that is different for each segment. As the totality of these spectra is no longer a contiguous span, the phase retrieval process will only attempt to find an aggregate spectrum phase solution that is, on average, flat across the span. There is no rationale for believing that a global defocus phase term will arise in the eventual titrated FPM image wavefront. Another way of interpreting this is to note that the small NA associated with the objective is simply incapable of sampling the system’s global defocus associated with large $k_x,k_y$ values. If the raw data does not contain this information, it is difficult to imagine how the final processed FPM data will contain the information.

The lack of a global defocus phase term overlaid on top of the converged FPM solution means that any attempt to perform postreconstruction digital refocusing to further sharpen the image would be futile. This is because, in mathematical terms, digital refocusing means applying a reverse radial defocus phase term to compensate for a global defocus phase term. If such a term is absent, applying a reverse radial defocus phase term will either have no significant net impact or a negative impact on image quality.

Perhaps, a more relatable way to put this in context is as follows. The FPM algorithm outputs an image to the best of its ability. In the case in which defocus is present, this image would appear blurry, and the greater the defocus distance is, the blurrier the image is. Although the image contains phase information, the FPM titrated phase information is not a true copy of the sample’s phase wavefront. As such, the quality of the image cannot be improved by applying a reverse defocus phase correction. The only way this FPM image can be improved is by rerunning the FPM algorithm with the correct defocus phase correction function incorporated into the algorithm. There is no possible shortcut of performing postreconstruction digital refocusing by simply applying a reverse defocus phase correction.

In totality, postreconstruction digital refocusing does not perform any meaningful function at any defocus distance.

3.

Experiments

The experimental setup of the FPM system is shown in Fig. 5. The system is a modification of an Olympus IX81 inverted microscope. A $32\times 32$ LED matrix (Adafruit $32\times 32$ RGB LED Matrix panel) with a 4-mm separation distance between two adjacent LEDs illuminates the sample from the top. To help find the region of interest, the sample is placed on the Thorlabs MLS203-1 XY scanning stage. The height from the LED array to the sample is set at 54 mm. A 20× 0.40 NA microscope objective lens (Olympus Plan N 20/0.40) is used in the experiments. A $z$-axis piezoelectric transducer translation stage can move the objective lens with nanometer-level repeatable accuracy. The z-axis translation stage is embedded in the microscope frame and commercially available with Olympus IX81. A CMOS camera (ISVI IC-M50T-CXP-FMB) with a pixel pitch of $4.6\mu \mathrm{m}$ and a resolution of $7920\times 6004$ is used for capturing intensity images, but only the central section of the camera frame is used. In the following experiments, $15\times 15$ LEDs of the LED matrix in the green channel are illuminated sequentially. The green LEDs have a central wavelength of 522.6 nm and a full width at half maximum of 20 nm. For system calibration, one of the LED is aligned with the optical axis of the microscopes with manual adjustments as the normal incidence, and an LED position calibration algorithm²⁴ is used to estimate and calibrate the illumination angle of the remaining used LEDs.

Fig. 5

The experimental setup for (a) FPM system and (b) a schematic illustration.

A Siemens Star is used as the testing object for the performance analysis of the two digital refocusing methods.^12,13,25 The sample is captured at various defocus positions ranging between 0 and $6\mu \mathrm{m}$. The collected data are processed with the gradient-descent-based FPM algorithm³ with the embedded pupil function recovery (EPRY)¹⁵ subalgorithm also in place to correct for aberrations. For the standard in-iteration digital refocusing, the FPM algorithm is given the defocus distance during processing. For postreconstruction digital refocusing, we first execute the FPM algorithm with the defocus distance set to zero, thereby generating a “no defocus correction” FPM image. Then, the postreconstruction digital refocusing FPM image is generated by imposing a reverse global defocus term on the “no defocus correction” FPM image dataset.

3.1.

Experimental Verification that FPM Does Not Converge Correctly for Large Defocus Distance

In our mathematical analysis, we discussed the failure of the FPM algorithm to converge on the correct focal plane wavefront when the defocus distance between the sample and the microscope focal plane is large. To check whether the FPM algorithm has arrived at the correct solution, we can crop out an angular spectrum area that corresponds to the angular spectrum of a raw image. The intensity image rendered from this cropped segment should match the raw intensity image.

Here, we perform such a test. First, we place a Siemens Star target at a defocus distance of $4\mu \mathrm{m}$ from the objective’s focal plane. Then, to remove the defocus aberration introduced in the optical system, we perform two versions of the digital refocusing. The images in the middle column in Fig. 6 are the FPM raw data with normal illumination ${\theta }_{\mathrm{illum}}=0\mathrm{deg}$ and oblique illumination ${\theta }_{\mathrm{illum}}=20\mathrm{deg}$. The images in the left column in Fig. 6 are from the reconstruction result using standard in-iteration digital refocusing. To obtain the intensity image, the corresponding spectrum area is cropped out, and an inverse Fourier transform with a modulus square is performed. The reconstructed images agree well with the FPM raw images. However, for the postreconstruction digital refocusing method (the images in the right column in Fig. 6), the corresponding intensity images significantly deteriorate in quality and do not match the FPM raw images. This clearly showed that, at large defocus distances, the FPM algorithm fails to converge well. As we are not getting a close-to-ground-truth wavefront solution, we cannot expect postreconstruction digital refocusing to serve a meaningful function for large defocus distances.

Fig. 6

FPM raw data comparison with reconstructions using the two digital refocusing methods for Siemens Star. The sample is placed at a defocus distance of $4\mu \mathrm{m}$. The images in the middle column are raw FPM images with illumination angles 0 deg and 20 deg, respectively. The images in the left column are standard in-iteration digital refocusing reconstructions with corresponding spectral areas. Postreconstruction digital refocusing results are in the right column with corresponding spectral areas. FPM arrives at its correct solution when it is able to reproduce the raw FPM input images. The poor matches observed in the postreconstruction digital refocusing results are evidence that the FPM did not converge onto the correct solution.

3.2.

Experimental Verification that, in the Short Defocus Distance Regime, Spectrum of FPM Image Does Not Contain a Global Defocus Phase Term

Now, we examine the performance of FPM in the short defocus distance regime. Here, we only focus on the defocus distances of 1 and $3\mu \mathrm{m}$. We aim to directly verify whether the FPM image contains a global defocus phase term that is vital for postreconstruction digital refocusing.

Figure 7 shows the experimental results for the defocus distance of 1 and $3\mu \mathrm{m}$ of Siemens Star. Figure 7(a) shows the in-focus image of the Siemens Star, which will serve as the “ground-truth” reference. For the defocus distance of $1\mu \mathrm{m}$, we can see that in-iteration digital refocusing FPM processing [Fig. 7(b)] can produce a close match to the “ground-truth” reference. When FPM is used to process the same data without correcting for the defocus, the rendered image (“no defocus correction” FPM image) [Fig. 7(c)] is still fairly close in appearance to the “ground-truth.” When a reverse global defocus phase term is then applied to the data, we see a significantly deteriorated image [Fig. 7(d)] – clearly showing that postreconstruction digital refocusing has failed.

Fig. 7

(a) The FPM reconstructed image for the Siemens Star that is placed at the focal plane (“ground-truth”). For images (b–d), the Siemens Star is placed at $1\text{-}\mu \mathrm{m}$ defocus plane. FPM reconstructed image generated (b) with in-iteration digital refocusing, (c) with no defocus correction, and (d) with postreconstruction digital refocusing. (e) The red curve is the expected global defocus phase term with a $1\text{-}\mu \mathrm{m}$ defocus distance; the blue curve is the spectrum phase of the complex field of (c) averaged over the azimuthal direction. $k_{\text{radial}}$ is the spatial frequency along the radial direction. For better visualization, we plot the entire diameter of the spectrum phase from negative to positive spatial frequencies. (f–j) A repetition of (a–e) for the Siemens Star placed at a $3\text{-}\mu \mathrm{m}$ defocus plane. $\otimes$ denotes the convolution operation. The SSIM numbers are included to show the relative similarity of the images to the “ground-truth.”

This experiment allows us to further determine whether the FPM rendered image contains a global defocus term. To check this, we took the “no defocus correction” FPM image from Fig. 7(c) and computed the average phase as a function of spatial frequency displacement along the radial direction. (Specifically, we averaged the phase values along the azimuthal axis in the spatial frequency domain.) Figure 7(e) shows the result of the computation, along with the predicted global defocus phase term associated with a defocus distance of $1\mu \mathrm{m}$. The global defocus phase term for this defocus distance should look like the red trace in Fig. 7(e), and we can expect the average phase curve to follow this trace if the premise of postreconstruction digital refocusing is correct. On the other hand, if our mathematical analysis is correct, the average phase should be a relatively flat line. We can see from our experimental data that the average phase is consistent with our mathematical analysis – there is no global defocus phase term in the FPM algorithm output. This absence of the global defocus phase term in the FPM algorithm output also explains the deteriorated postreconstruction digital refocusing image – applying a reverse global defocus phase correction to the FPM-processed data, which does not contain any global defocus phase term, will only introduce additional errors into the image.

The same observations can be made for the defocus distance of $3\mu \mathrm{m}$ [Figs. 7(f)–7(j)]. We additionally included the corresponding plots for defocus distances of 1 to $6\mu \mathrm{m}$ in Appendix A. All of these experiments show a consistent result – the FPM-processed image does not contain a global defocus phase term.

3.3.

Experimental Verification that Applying a Reverse Defocus Phase Term in Postreconstruction Digital Refocusing Worsens Image Quality

Here, we report our findings on the image quality performance of the three methods. We performed quantitative analysis on the resolution and the SSIM of the Siemens Star pattern images for different defocus distances. For this experiment, we used a defocus distance that ranged from 1 to $6\mu \mathrm{m}$. The resolution analysis was performed as follows.

Rings of pixels aligned with the center of the Siemens Star are extracted with different radii. Each ring of pixels is a one-dimensional (1D) array containing a 40-line pairs square function pattern. The larger the radius of the ring is, the lower the spatial frequency of the 1D array is. For each 1D array, the modulation transfer function (MTF) as a performance metric is calculated. The computation of MTF follows the methodology described in Ref. ²⁶. The MTFs of the reconstructed intensities for Siemens Star in-focus and defocused at 1 to $6\mu \mathrm{m}$ with a step size of $1\mu \mathrm{m}$ are calculated. The threshold value of MTF is set to be 10%. Figure 8(a) shows the plot of the threshold resolution values for different defocus distances and digital refocusing methods. For defocus distances from 1 to $6\mu \mathrm{m}$, a resolution of $0.73\mu \mathrm{m}$, which is the same as the achieved resolution under in-focus conditions (“ground-truth”), can still be achieved by the standard in-iteration digital refocusing method. This means that the in-iteration digital refocusing method can efficiently correct the defocus aberration. However, the postreconstruction digital refocusing method fails to provide the same resolution for defocus distance beyond $3\mu \mathrm{m}$.

Fig. 8

(a) SSIM and (b) resolution at 10% threshold of MTF of the Siemens Star at defocus distances from 0 to $6\mu \mathrm{m}$. The black plots are for the standard in-iteration digital refocusing method; the red curves depict the metrics for the postreconstruction digital refocusing method; the green lines are under the condition that no defocus corrections are performed; and the magenta vertical lines indicate the half-DOF value ($1.95\mu \mathrm{m}$) of the objective lens.

The comparable resolution for all three approaches for defocus distances up to $2\mu \mathrm{m}$ is consistent with our mathematical analysis. Our analysis in Sec. 2 indicated that, for a sufficiently small defocus distance, the individual radial defocus phase term associated with each raw image would be small enough to ignore and the FPM algorithm can converge onto an approximate “ground-truth” solution. The resolution is therefore not significantly deteriorated for the “no defocus correction” approach in the small defocus distance regime. In the same context, the reverse global defocus phase term applied by the postreconstruction digital refocusing method is small as well, and it does not significantly impact the image resolution performance for small defocus distances.

However, it is interesting to note that the images associated with postreconstruction digital refocusing showed clear signs of image quality degradation as the defocus distance increases [see Fig. 8(a)]. To examine this in greater detail, we plotted the SSIM and compared it with the in-focus image (“ground-truth”) for all three methods [Fig. 8(b)]. The SSIM is more sensitive to the texture of the image and can more sensitively pick up subtle image changes than resolution degradation. From our plots, we can see that the standard in-iteration digital refocusing method shows the slowest SSIM degradation as a function of the defocus distance. The “no defocus correction” approach follows the curve of the standard in-iteration digital refocusing method for small defocus distances before degrading more rapidly for larger defocus distances. This behavior is consistent with our mathematical analysis. For small defocus distances, the individual radial defocus phase term associated with each raw image would be negligible, and the FPM algorithm can converge onto an approximate “ground-truth” solution for the “no defocus correction” approach. In contrast, the postreconstruction digital refocusing shows a steeper SSIM degradation as a function of defocus distance. This is consistent with our mathematical prediction that, because the FPM algorithm outputs a result that does not contain a global defocus phase term, when postreconstruction digital refocusing attempts to apply a reverse global phase term to correct for a nonexistent global defocus phase term, it would only deteriorate the image.

Finally, we note that, for large defocus distances (from 4 to $6\mu \mathrm{m}$), the SSIM values for the postreconstruction digital refocusing method are higher than those for the “no defocus correction” method. We have already established that, at large defocus distances, the FPM solution converges onto a local minimum solution that has no physical meaning or correspondence. These oddities are likewise devoid of physical meanings. Readers can refer to Appendix A for the reconstructed intensity image details at different defocus distances. For large defocus distances, image quality assessment is best done with resolution. In terms of resolution for large defocus distances, both the “no defocus correction” method and the postreconstruction digital refocusing method clearly fare much worse than the standard in-iteration digital refocusing method.

3.4.

Experimental Verification on a Pathology Slide Sample

In our final experiment, we used a blood smear slide as our test target. First, the sample is imaged at the in-focus position [Fig. 9(a)] to serve as a “ground-truth” image. Then, the sample is reimaged at a defocus distance of $2\mu \mathrm{m}$. Figure 9(b) shows the in-iteration digital refocused FPM image of the sample. Unsurprisingly, the image is almost identical to the “ground-truth” image. Figure 9(c) shows the “no defocus correction” FPM image of the sample, for which no defocus correction was done. The image is still a fairly close match with that of the “ground-truth,” but some image quality loss is clearly apparent. By applying a reverse defocus phase correction to this image dataset, we can then derive the postreconstruction digital refocused image [Fig. 9(d)]. This image has poorer quality than the other images. We can also see that the “no defocus correction” FPM image does not contain a global defocus phase term [Fig. 9(e)]. Once again, these findings are consistent with our mathematical analysis that concluded that the FPM algorithm does not output a global defocus phase term.

Fig. 9

The FPM intensity reconstruction for a blood smear pathology slide placed at (a) focal plane and a $2\text{-}\mu \mathrm{m}$ defocus plane with (b) standard in-iteration digital refocusing method; (c) no defocus correction; and (d) postreconstruction digital refocusing method. (e) The red curve is the expected global defocus phase term with a $1\text{-}\mu \mathrm{m}$ defocus distance; the blue curve is the spectrum phase of the complex field of (c) averaged over the azimuthal direction. $k_{\text{radial}}$ is the spatial frequency along the radial direction. For better visualization, we plot the entire diameter of the spectrum phase from negative to positive spatial frequencies. $\otimes$ denotes the convolution operation.

4.

Discussions

Our analysis and experimental results showed that postreconstruction digital refocusing is not feasible with FPM raw data. In the presence of defocus, the angular spectrum data of each raw FPM frame contains a local radial defocus phase term that is centered in the middle of each raw image angular spectrum for large defocus distances (sample is outside the FPM’s DOF). As such, the angular spectra of the frames do not overlap in a continuous fashion and, thus, phase retrieval cannot be applied to generate the correct FPM angular spectrum. Though applying phase retrieval in such a scenario will still lead to a local minimum solution, the solution has no physical meaning or correspondence. Moreover, the local minimum solution cannot reproduce the original raw FPM image data – another indication that the solution is incorrect.

One way to interpret why defocused FPM data cannot undergo easy postreconstruction digital refocusing is to recognize that FPM uses tilted illumination to bring off-axis angular spectrum segments into the accepting numerical aperture of the microscope. In doing so, FPM recenters the selected off-axis angular spectrum segment on the microscope’s optical axis. As the defocus phase transformation is centered along that optical axis, each FPM off-axis spectrum segment will receive a radial defocus phase term centered on their respective origins. In comparison, if we simply illuminate the sample at normal incidence and collect a large NA angular spectrum image in a single shot, the data thus collected will have only a single radial defocus term that impacts all of the spectrum data. The nature of the data for these two situations is clearly quite different.

We project that this finding may seem counter-intuitive to some researchers. There is an often-held notion that the FPM method can be used to recover the focal plane wavefront of the imaging system. As such, digital refocusing should be possible without encountering the issue that we have just raised. However, this notion is incorrect. The FPM method is not an assumption-free method for recovering the focal plane wavefront. Instead, FPM needs the following two conditions to be satisfied so that it can be used to recover the focal plane wavefront:

The consequence of not adhering to these conditions is why postreconstruction digital refocusing fails.

We also examined the situations in which the defocus distance is small, which causes the individual radial defocus term for each raw image to be small. Here, the FPM can converge onto a local minimum solution that can approximate the “ground-truth” solution. However, this solution should be considered to be nominally in-focus already. It does not contain a global defocus phase term that can be further corrected to sharpen the image using the postreconstruction digital refocusing method. We further confirmed that this global defocus phase term is never present in the FPM output in our experiments reported in Secs. 3.2 and 3.4. The absence of this global defocus phase term means that postreconstruction digital refocusing, which consists of applying a reverse global defocus phase term, would serve no meaningful purpose.

In totality, this means that postreconstruction digital refocusing does not serve any meaningful function at all defocus distances.

Although straightforward postreconstruction digital refocusing of FPM data is neither meaningful nor feasible, our analysis clarifies the mathematical framework and suggests a few viable strategies. We discuss these strategies as follows.

The first strategy is already well familiar to the community – specifically, to nest the defocus correction procedure within the phase-retrieval loop (standard in-iteration digital refocusing approach). By flattening the radial defocus term for each raw angular spectrum frame, we ensure that the full angular spectrum undergoing phase retrieval is continuous and amenable to phase retrieval. The appeal of postreconstruction digital refocusing is that it has a much lighter computational load, and the processed dataset is highly compact – effectively 2D in nature. Although the standard in-iteration digital refocusing approach is computationally intensive and likely nothing can be done to lighten the computation load, this standard approach can still be used to render a processed dataset that is highly compact. The strategy would be to compute an exhaustive set of defocus corrected images at different depths. Then, any of the reported all-in-focus algorithms²⁸ can be applied to the stack of images to select patches of in-focus images to generate an effective 2D image in which every point is in focus.

The second approach would be to perform FPM imaging in its aperture scanning modality.^8,27,29–31 In this modality, instead of using a tilted illumination to selectively move different portions of the angular spectrum into the collecting numerical aperture of the microscope, we use a high NA imaging system and illuminate the sample at normal incidence. To select different portions of the angular spectrum for detection, an aperture is then scanned at the Fourier plane of the microscope. In this configuration, the angular spectrum is never shifted off the optical axis, and therefore, the same defocus term is applied uniformly to the optical field (a single global radial defocus phase term). For one portion of the angular spectrum, the defocus aberration for this specific portion is a segment of the global radial defocus phase term. After the phase retrieval loop, the different portions of the defocus aberration will be stitched together, and the same global radial defocus phase term will be recovered in the process. That is to say, the collected raw data can be processed by an FPM algorithm to yield an angular spectrum dataset in which defocus will manifest as a uniform global radial defocus phase term. Because this global defocus term is present in aperture scanning FPM, the digital refocusing can be performed by applying a reverse global phase term after the phase retrieval loop. We note that the aperture scanning FPM is distinctly different from the standard FPM. As such, their respective properties do not generalize between them.

The third approach would be to use the standard FPM system in which the objective is replaced with one with a very low NA. The span of each angular spectrum raw data segment is determined by the system’s NA. Using a low NA, the radial defocus phase term can remain relatively flat across the entire segment. This flatness will ensure that we can still perform phase retrieval across a continuous angular spectrum data space to render a processed FPM image. However, this approach will start to fail for excessively large defocus distances because increased defocus distances give a sharper defocus radial defocus phase term that will eventually break our flat defocus phase assumption. Practically, this approach is worth considering because the processed image would not need digital refocusing, and all planes near the focal plane will be simultaneously in focus. For reference, consider a low NA imaging system with $\mathrm{NA}=0.05$ and an illumination wavelength of 500 nm. If we stipulate that the defocus phase error cannot exceed 1 rad to ensure phase retrieval can proceed effectively, with 225 raw FPM image frames (assuming $\sim 60\%$ overlap), we can generate an FPM image with a synthetic NA of 0.2 and a natural DOF of $44.2\mu \mathrm{m}$. In contrast, a standard $\mathrm{NA}=0.2$ FPM (coherent) microscope would have a DOF of $24.7\mu \mathrm{m}$. (We note that the FPM’s natural DOF is defined differently in different publications.^12,19,23 In this study, we adopted the definition in Refs. 12 and 19 as it is verified with experimental data.)

5.

Appendix A

The additional figures for experimental results of Siemens Star at defocus distances of 1-6 microns are listed below in Figs. 10–12.

Fig. 10

FPM images at defocus distances of 1 and $2\mu \mathrm{m}$ with different digital refocusing methods. For the plots in the right column, the red curves are the expected global defocus phase term with a defocus distance of 1 and $2\mu \mathrm{m}$, respectively; the blue curves are spectrum phases of the complex field for “no defocus correction” averaged over the azimuthal direction. $k_{\text{radial}}$ is the spatial frequency along the radial direction, from negative to positive spatial frequencies for better visualization. $\otimes$ denotes the convolution operation.

Fig. 11

FPM images at defocus distances of 3 and $4\mu \mathrm{m}$ with different digital refocusing methods. For the plots in the right column, the red curves are the expected global defocus phase term with a defocus distance of 3 and $4\mu \mathrm{m}$, respectively; the blue curves are spectrum phases of the complex field of “no defocus correction” averaged over the azimuthal direction. $k_{\text{radial}}$ is the spatial frequency along the radial direction, from negative to positive spatial frequencies for better visualization. $\otimes$ denotes the convolution operation.

Fig. 12

FPM images at defocus distances of 5 and $6\mu \mathrm{m}$ with different digital refocusing methods. For the plots in the right column, the red curves are the expected global defocus phase term with a defocus distance of 5 and $6\mu \mathrm{m}$, respectively; the blue curves are spectrum phases of the complex field of “no correction,” averaged over the azimuthal direction. $k_{\text{radial}}$ is the spatial frequency along the radial direction, from negative to positive spatial frequencies for better visualization. $\otimes$ denotes the convolution operation.