1. Introduction and background

Imaging lenses, ranging from microscope objectives to satellite-based cameras, are limited in the total number of features they can resolve. This resolution limit is a function of both the imaging system’s point-spread function (PSF) size and its inherent aberrations across the image plane field of view (FOV). Typically referred to as an optical system’s space-bandwidth product (SBP) [1], this physical resolution limit is uncoupled from the number of pixels at the detector plane. While the SBP scales with the dimensions of the lens [2], it is often on the order of 10–100 megapixels for most applications. For example, a typical microscope objective can offer a sharp intensity PSF (0.5 μm) but only across a narrow FOV (1 mm), while a wide-angle lens can offer a large FOV (10 mm) but at the expense of a blurry intensity PSF (5 μm). Both systems can resolve a maximum of 10⁷ features, which one may ideally sample with a 10 megapixel digital detector.

Often, additional optical elements are used to correct for the aberrations that limit a camera’s SBP. These extra elements simultaneously increase the size, weight and complexity of the physical setup. In this paper, we explore an alternative approach to improving a camera’s resolution performance. Instead of adding corrective optical components, we capture a series of images and apply digital post-processing to measure and remove aberrations. This procedure, termed overlapped Fourier coding (OFC), relaxes the complexity of the optics and digital detector while ensuring a large optical system SBP (defined in [1]). Since maintaining this large SBP now requires multiple image acquisitions, we effectively trade off lens complexity for acquisition speed. Our correction scheme also acquires the phase of the optical field exiting the sample plane, which may in turn be used to digitally refocus the sample post-capture to ensure all image segments are in sharp focus.

OFC first acquires a sequence of images, modulating the camera’s Fourier plane with a shifting aperture function between each snapshot. Here, the shifting aperture function we use is a small square pupil generated on an amplitude-modulating spatial light modulator (SLM), although other aperture function types, including phase-only modulation, may benefit alternative applications. It is important the shifted sub-apertures spatially overlap with one another in the pupil plane by a certain degree (here we use a 75% degree of overlap). Second, OFC uses a unique algorithm to digitally synthesize the sub-aperture images into a high-resolution estimate of the complex optical field at the image plane. This algorithm extends prior work using multiple images captured under varied external illumination provided by an LED array placed beneath a conventional microscope’s sample plane [3]. As such, it is directly related to the well-known methods of phase-retrieval [4, 5] and ptychography [6, 7].

Unlike other holographic [8, 9] or non-holographic [3, 10–12] resolution-enhancing optical setups that measure phase, this work does not use external illumination. So while OFC cannot extend an optical system’s resolution beyond its aperture-defined cutoff, it can increase its SBP via removal of undesired aberrations and misalignments. This sets our goal as distinct from prior methods using aperture-based [13] and SLM amplitude modulation [14] to only acquire phase, without correcting for system aberrations. Although not used by the current OFC setup, phase-based modulation may also lead to accurate sample phase recovery [15–17]. Furthermore, an SLM can help acquire an incoherent light field near a camera’s native spatial resolution [18, 19]. Resolution-improved incoherent intensity images may also be obtained via SLM-generated structured illumination [20–22]. Since incoherent setups do not measure phase, their digital removal of system aberrations must follow an ill-posed intensity PSF deconvolution, previously attempted within both conventional fluorescent [22, 23] and light-field [24] microscopy. By simultaneously measuring a complex sample field and estimating the camera’s complex coherent transfer function, we avoid the common “divide-by-zero” challenges faced when inverting an aberrated optical transfer function [25].

We recover an aberration-free complex sample function using simulated annealing, which iteratively reduces discrepancies between the image’s digital Fourier representation and the camera system’s physical Fourier plane. Annealing is a well-studied optimization procedure [26] that has been previously applied within the area of X-ray ptychography to correct for an optical probe’s unknown shifted location [27, 28]. X-ray ptychography [6, 7] and coherent diffractive imaging [29, 30] setups have also measured and removed the effects of their illumination probe. The proposed algorithm is closely related to this prior work, but instead searches over either a set of unknown Zernike aberration parameters, or a space of Fourier plane misalignments, or both, to increase a camera’s SBP over multiple acquisitions. While many previous systems [31–33] can pre-calibrate for optical aberrations for later removal (assuming they also quantitatively measure sample phase), OFC can continually update an improved estimate of system imperfections from each captured image. This may prove beneficial in microscopy applications where objective lenses are often shifted or replaced, or in conventional cameras where aberrations are a function of object depth and zoom lens position. In a broader sense, this situates the OFC scheme close to the realm of adaptive optics, which utilizes correction schemes targeting and removing optical distortions that can change with each acquired image.

The remainder of this paper is outlined as follows. In Section 2, we present the OFC imaging strategy and discuss its recovery algorithm in the context of an aberration-free setup. In Section 3.1, we introduce the experimentally realistic situation of a camera containing low-order optical aberrations that we model with Zernike circle polynomials [34]. We then demonstrate how OFC uses simulated annealing to estimate and remove these aberrations from its final reconstruction. In Section 3.2, we show how the simulated annealing approach may additionally account for aberrations that do not fall within the standard Zernike polynomial model, such as distortion and lens misalignment. Finally, in Section 4 we experimentally demonstrate OFC with an SLM-based shifting sub-aperture to recover a sample’s image and phase with resolution at the camera’s maximum cutoff spatial frequency. Our annealing algorithm likewise helps remove included aberrations to increase the setup’s SBP.

2. Principle of operation

In this section, we first develop a mathematical model of OFC’s imaging routine. Then, we detail how the OFC algorithm reconstructs a sample’s amplitude and phase at high-resolution from a series of low-resolution measurements, assuming aberration-free optics.

2.1. Image acquisition

As diagrammed in Fig. 1, we implement overlapped Fourier coding in a conventional 4f imaging system with a simple modification: a transmissive amplitude SLM inserted into the camera Fourier plane (i.e., aperture plane). A 4f imaging system is not a requirement for effective operation – placing the SLM at any plane conjugate to the sample plane (i.e., the aperture plane of a compound lens system, or the back-focal plane of a microscope objective) offers similar functionality. We will assume the sample s(x, y) at the object plane is illuminated with spatially coherent and quasi-monochromatic light, as defined in [35]. OFC directly extends to operate with partially coherent illumination. We also assume that our imaging setup contains pixels that are matched to the maximum cutoff frequency of its unmodulated Fourier plane (i.e., each pixel matches the system’s minimum PSF width). This enables us to effectively verify the accuracy of our reconstruction through a direct comparison with an unmodulated image. A more practical OFC setup should match the pixel size to the NA of a single sub-aperture image, which leads to a final reconstructed resolution surpassing the detector pixel count.

 

Fig. 1 Outline of the OFC procedure. (a) We place a transmissive SLM in the Fourier plane of a 4f system to digitally create different sub-aperture functions. (b) We capture a sequence of aberrated images while the SLM displays a sub-aperture shifted to a unique location between each snapshot. (c) We computationally transform the captured image set into a high-resolution amplitude and phase map, as well as an estimation of the camera’s low-order aberrations.

Download Full Size | PDF

The square SLM of width L is configured to display a small square sub-aperture (i.e., optically transparent area) of width and height ℓ. In this initial demonstration, we select a square sub-aperture for its geometric simplicity and x–y separability, which helps minimize SLM pixelization artifacts. Alternative sub-aperture shapes like rectangles or circles may be easily incorporated into the following procedure. Our square sub-aperture is shifted to n² different equally spaced locations arranged on a 2D rectilinear grid. The j^th sub-aperture will be displaced from the origin by a two-dimensional vector c_j = (c_xj, c_yj) for 1 ≤ j < n². By setting n > L/ℓ, we ensure that each sub-aperture overlaps partially with its neighboring sub-apertures. In practice, it is useful to select n such that the sub-aperture images overlap by approximately 70–75% following a similar optimized parameter in ptychography [36]. All simulations and experiments in this work use n = 9, requiring a total of n² = 81 images per capture sequence. We set L/ℓ = 2.5 to ensure each sub-aperture window overlaps with its neighbors by 75%.

An ideal binary amplitude SLM will completely block light when its pixels are switched to opaque and pass 100% of any incident light when switched to transparent. In a practical setup, a transmissive SLM will have a finite optical density b, which we include in our model by defining the j^th SLM sub-aperture’s transmission function W_j as,

Modeling the optical field emerging from the sample surface as s(x, y), we can write the complex field directly before the aperture (i.e., SLM) plane as ℱ[s(x, y)] = Ŝ(k_x, k_y), where ℱ is the Fourier transform operation and we will refer to Ŝ as the sample’s spectrum. After passing through the j^th sub-aperture window, W_j, the exiting field is the product Ŝ(k_x, k_y)W_j(k_x, k_y) assuming the SLM is thin. This modulated field then propagates to the image plane detector, where to a first-order approximation the j^th image measures,

2.2. Aberration-free OFC reconstruction

The basic OFC post-processing algorithm employs the well-known strategy of alternating projections, a common phase retrieval method [4,5]. Our measured images constrain our sample estimate to a known set of amplitudes in the spatial domain, while our shifting SLM sub-aperture uniquely constrains its support in the Fourier domain. Assuming an aberration-free setup, the OFC algorithm reconstructs a complex field exiting the sample plane at the employed imaging system’s maximum resolution. Unlike in the related Fourier ptychographic microscope [3], the sample here may be thick, optically reflective, and under arbitrary quasi-coherent illumination. While the aperture scanning procedure in [13] also follows similar steps, OFC’s SLM modulator requires certain algorithm modifications that we outline below.

The following explains one iteration of the OFC process. To begin, we initialize an estimate of the unknown high-resolution complex sample spectrum, Ŝ₀(k_x, k_y), as shown in Fig. 2(a). A good starting point is an up-sampled version of the low-resolution image recorded with a centered sub-aperture, multiplied with a randomly assigned phase. Then, the iteration first computes a masked spectrum estimate by modulating Ŝ₀ with one of our SLM sub-apertures from Eq. (1). Starting with j = 0 and continuing until j = n² − 1, we compute the spectrum product P̂_j(k_x, k_y) = Ŝ₀(k_x, k_y)W_j(k_x, k_y). For an SLM with perfect contrast, this corresponds to selecting a windowed region of Ŝ₀, as shown in Fig. 2(b). Second, we compute the Fourier transform of this spectrum product to simulate light’s propagation to the image plane: p_j(x, y) = ℱ[P_j((k_x, k_y)]. Third, we replace the amplitude of the resulting Fourier transform p_j with the known amplitude from the j^th detected image $\sqrt{I_j}$ to form p′_j:

 

Fig. 2 One stage of the basic OFC algorithm. For each window position W_j, a segment of the spectrum estimate Ŝ₀(k_x, k_y) is extracted (corresponding to the j^th position of the SLM’s shifting sub-aperture). In the spatial domain, the amplitude associated with this windowed spectrum is constrained with the measured image I_j(x, y) to form p′(x, y). Example detected images are in (d), where we measure the simulation’s spatial frequency “cutoff” at Group 0, Element 4. The result is Fourier transformed back to the aperture plane, where it is used to update the spectrum estimate, Ŝ₀(k_x, k_y).

Download Full Size | PDF

Fifth, we update our sample spectrum estimate with P̂′_j. While we accounted for the SLM’s imperfect modulation when extracting P̂_j from Ŝ₀, we do not use the same support function when re-inserting P̂′_j in this update. Instead, we only update the aperture area that in an ideal case is optically transparent, leaving all other areas unchanged:

2.3. Aberration-free OFC simulation

Figure 3 demonstrates the OFC algorithm’s successful convergence to an accurate amplitude and phase measurement in simulation. With an ideal camera (no aberrations or misalignments, but including Gaussian noise), Section 2.2’s five update steps recover a sample estimate at the system’s full (i.e., unmodulated) resolution. Our ideal simulated 4f camera has a maximum F-number of 7.5 (PSF width= 5μm width at the detector assuming λ = 632 nm illumination), and its detector pixels are 5 μm to match this unmodulated PSF. We shift a square sub-aperture that is ℓ = 4 mm wide across a total distance of L = 10 mm, where each shift is Δc = 1 mm along one dimension. In two dimensions, this leads to a total of 81 images and an aperture overlap percentage of 75%. Simulated intensity images are N = 1000² pixels and have a 12-bit well depth. These setup parameters closely match our experimental setup (see Section 5) and are used for all subsequent simulations.

 

Fig. 3 OFC simulations for an aberration-free noisy imaging system. (a) Air force resolution “target” amplitude and phase (top) and the corresponding OFC reconstruction (bottom), including 2% noise. Note this simulated sample’s slowly-varying phase is poorly reconstructed. We box the reconstruction’s spatial frequency cutoff (Group 1, Element 5). (b) A blood smear “target” sample’s phase is reconstructed with much higher (10X) fidelity than in (a), again including 2% noise. For both simulated samples, reconstruction error grows linearly with system noise while the OFC algorithm continues to converge.

Download Full Size | PDF

In Fig. 3(a), we first reconstruct an Air Force resolution target sample U(x, y) with a 2.5 μm minimum feature size (i.e., half a pixel) with a multiplied cubic phase function: ∠U(x, y) = α(x³ + y³), where ∠ indicates phase and α = 20π. We add 2% random Gaussian noise to each sub-image’s intensities after detection to simulate detector-induced noise. For the sub-image captured with its sub-aperture centered on the optical axis, the minimum resolvable element in the Air Force target is Group 0, Element 4 (using Sparrow’s condition, boxed in Fig. 2(d)). After OFC reconstruction, the intensity’s minimum resolvable element shrinks to Group 1, Element 5 (boxed in Fig. 3(a)), showing an increase in the maximum resolvable spatial frequency by a factor of 2.5. This closely approaches an ideal factor of 3 system resolution gain, as the final synthesized aperture is 3 times wider than a single sub-aperture. The reconstructed phase resembles the target’s original phase, but contains significant artifacts in areas that lack significant phase variation (i.e., towards the image center). As noted in prior work, the ptychography update process has trouble accurately reconstructing such areas of low amplitude and phase variation [37]. This may be especially problematic when rectangular grid scanning is used. However, such a “challenging” sample will offer a useful target to benchmark our simulated annealing technique’s performance in Section 3.

We verify this phenomena with a second OFC simulation of a different target sample containing highly varying phase (a previously measured amplitude and phase map of a monolayer of blood cells, Fig. 3(b)). Here, even in the presence of 2% detector noise, we observe near-perfect phase reconstruction. The normalized mean-squared error (NMSE) E_q between the target sample U(x, y)’s and reconstruction s_q(x, y)’s amplitudes offers a useful metric to examine algorithm performance, as also used in [4]:

3. OFC with simulated annealing

By imaging with an overlapping sub-aperture set, OFC captures slightly redundant data. This redundancy not only allows us to accurately extract sample phase, but also helps us compute a low number of unknown variables that influence the imaging process. A well-known method of searching over a space of unknown variables to maximize or minimize a particular function of merit (i.e., error function) is simulated annealing (SA). Instead of exhaustively searching through all possible unknown variable configurations, SA takes an iterative approach. At one iteration, annealing first randomly searches through a small number of different configurations and selects the configuration that minimizes its error function. Then, it uses this configuration as starting point for the next search. As iteration continues, the algorithm slowly reduces the range over which it randomly searches for error-minimizing states. For many problems, this type of iterative local search is very efficient at seeking out global minima of nonlinear functions [26].

OFC’s space of unknown variables includes any optical system aberrations or misalignments. At each sub-aperture location, our function to minimize is the mean-squared error difference between the recorded intensity image, I_j(x, y), and the corresponding image that would result from detecting our complex windowed spectrum estimate, P̂(k_x, k_y). For each iteration, we make several guesses about how the optical system might be aberrated or misaligned, compute the resulting image captured through each imperfect camera, and then select the imperfections that yield the closest image to our observed data as the starting point for the next iteration. As we will demonstrate, this process is both effective at recovering the correct imperfections and robustly removing these imperfections from a final sample solution. However, we must assume a-priori knowledge about which subspace of imperfections to search through, as search time will scale linearly with the dimensions of this search space.

3.1. Characterization and removal of low-order aberrations

We may account for our 4f setup’s wavefront-based aberrations using a multiplicative phase function, A(k_x, k_y), at its Fourier plane. It is common to decompose A into a sum of weighted Zernike polynomials on the unit disk: A(ρ, θ) = exp(2πi∑_l a_lW_l(ρ, θ)), where W_l are the orthogonal Zernike polynomials and a_l are the associated weights [34]. In this section, we will consider the most significant Zernike polynomial used to describe an optical aberration, defocus. A Zernike polynomial decomposition’s defocus term takes the Cartesian coordinate form,

In addition to an estimated spectrum Ŝ₀, simulated annealing OFC (SA-OFC) also initializes an estimated aberration map A₀ (Fig. 4). We always initialize with A₀(k_x, k_y) = 1. SA-OFC then begins with j = 1 and computes T different candidate aberration functions at the j^th loop:

 

Fig. 4 Schematic of the OFC algorithm with simulated annealing (SA-OFC). We use the same phase retrieval loop outlined in Fig. 2 with three additions. First, an estimate of the system’s aberrations at the Fourier plane A_j is now multiplied with the spectrum estimate Ŝ at each loop. Second, SA is used to compare T different perturbed versions of the j^th windowed spectrum with the j^th intensity measurement. Third, the error-minimizing aberration perturbation $A_j^{t_{\mathit{min}}}$ and corresponding sample estimate ${\psi }_j^{t_{\mathit{min}}}$ are identified with Eq. (9), which then update A_j+1 and Ŝ_j+1 via Eq. (10)–Eq. (11).

Download Full Size | PDF

While Eq. (11) resembles the sample update strategy of earlier ptychography algorithms in [7, 38], SA-OFC must depart from this prior work to overcome a key challenge. All previous ptychographic setups (including Fourier ptychography) acquire a series of images that are similarly altered by the same aberrated probe beam/microscope aperture. Such redundancy enables direct extension of alternating projections methods to jointly solve for the sample and probe/aperture function. In OFC, however, each sub-image is altered by a different windowed area of the aberrated camera aperture. The useful redundancy helping the extended alternating projection algorithms in [7,38] converge to an accurate joint solution is no longer present. OFC must instead iteratively update the global aberration function A(k_x, k_y) in Eq. (10) with partial information from each sub-image through annealing. We now demonstrate this alternative strategy also converges to an accurate aberration measurement.

Figure 5 presents an example simulation of the SA-OFC algorithm using the same 4f setup from Section 2. The target sample matches Fig. 3(a)’s, but now with α = 5 defining its cubic phase envelope. We add defocus aberration by multiplying A(k_x, k_y) in Eq. (6) with each simulated sub-aperture spectrum, setting d = 200. To first test the effect of aberrations without simulated annealing, we run Section 2’s basic OFC algorithm to recover the amplitude and phase maps shown in Fig. 5(a). Because the algorithm incorrectly updates our sample estimate with aberrated low-resolution intensity images, the final solution does not closely resemble the target sample. Switching to Section 3.1’s SA-OFC algorithm enables simultaneous aberration map estimation and removal, as shown in Fig. 5(b). For this example, we use the same initial conditions and noise in Section 2, setting T = 8, q = 20, and r_a = d/2. The recovered aberration map A_q’s defocus coefficient d_q differs from ground truth by Δd = (d_q − d)² = 25.71. Likewise, we now measure the Sparrow limit frequency cutoff as Group 1 Element 3, which is 80% the ideal cutoff due to residual defocus error. Annealing generally performs better with more candidate search functions T, but will require linearly more computation [27]. It is also sensitive to the selected search radius r_a. All simulations here use r_a = d/2 and linearly decrease it to zero by the q^th iteration with α_a = r_a/(q + 1). It may prove useful to optimize over these free parameters, and also randomly offset each sub-aperture center to minimize artifacts in the sample’s reconstructed phase.

 

Fig. 5 Simulation results of the SA-OFC algorithm. (a) Without simulated annealing, OFC cannot combine a set of aberrated sub-images into an accurate full-resolution complex field estimate. The induced defocus aberration is show to right. (b) The annealing OFC algorithm accurately recovers both the sample’s amplitude and phase and the 4f setup’s aberration map A_q(k_x, k_y). Here, q = 20. (c)–(d) Example low-resolution sub-aperture images, for comparison.

Download Full Size | PDF

We plot the SA-OFC algorithm’s NMSE from Eq. (5) versus iteration q in Fig. 6(a). Here, we compare simulated annealing (blue curve) with two alternative post-processing options. We use the same parameters as in Fig. 5 but now average over 5 runs with σ² = 2% noise. First, no annealing (“no SA”, red) leads to a significantly higher sample recovery error (see Fig. 5(a)) but maintains error-decreasing performance, demonstrating algorithm stability. Second, we assume a-priori knowledge of the 4f setup’s aberration map A₀(k_x, k_y) (“known”, green) to recover an almost exact solution. Here, instead of using Eq. (10) to estimate the aberration function each iteration, we simply set A_j+1(k_x, k_y) = A₀(k_x, k_y), the actual aberration map in Eq. (6), for all j. The SA-OFC algorithm closely approaches this exact removal process (within 5% after 20 iterations). Thus, while computational limitations prevent exact recovery, OFC’s captured dataset certainly contains enough redundancy to significantly improve an image’s SBP via defocus aberration removal.

 

Fig. 6 (a) SA-OFC is error reducing (blue) and exhibits much lower error than the regular OFC algorithm (red) in the presence of aberrations (assuming defocus d = 200). The ideal case of performing SA-OFC with an a-priori known aberration initialized and enforced each iteration is plotted in green. (b) As the aberration size increases, SA-OFC, OFC and the ideal case all slowly decrease in MSE performance. SA-OFC exhibits an aberration recovery mean-square error Δd scaling roughly as 10–15% of d, which may be improved with additional fine-tuning of annealing parameters.

Download Full Size | PDF

We once again repeat Fig. 5’s simulation in Fig. 6(b), but now vary the amount of defocus aberration d in the ground-truth aberration A(k_x, k_y). Again, each plot point is an average over 5 independent runs with σ² = 2% noise. As expected, SA-OFC remains bounded below by the case of knowing and directly accounting for a complex aberration map. However, it achieves much lower error than an annealing-free algorithm, even for significantly defocused image sets.

3.2. Characterization and removal of geometric distortion

A number of camera imperfections that negatively impact image fidelity cannot be summarized as a Fourier plane phase-only modification (i.e., with the Zernike polynomial aberration model). Examples of such imperfections include unknown changes in magnification, image distortion, and vignetting caused by system misalignment. In this section, we explore how to correct for one specific form of misalignment that significantly impacts our experiments: the displacement of the Fourier plane from its assumed location, which is directly connected to image distortion. We outline how the SA-OFC algorithm accounts for this type of unknown parameter by updating the estimated size and location of each SLM sub-aperture, leading to an improved-resolution complex sample reconstruction.

It is challenging to construct an imaging system with a perfectly flat, centered Fourier conjugate plane (i.e., aperture plane or back focal plane). Slight curvature across the plane is often encountered, even with the aid of advanced lens design software [39]. Element misalignment during system assembly will additionally shift the Fourier plane away from its ideal location. Misalignments may point perpendicular to the optical axis, which will cause the Fourier plane to laterally shift and add a linear phase ramp across the image plane. Or, they may point axially, which will primarily cause the Fourier plane to scale in magnification. Depending upon the size and linearity of this scaling, the field at the image plane can become magnified, geometrically distorted and also defocused.

The SA-OFC algorithm can help measure and remove this type of misalignment. First, let’s assume that the location of each SLM pixel is known a-priori and we can accurately measure the optical field associated with one sub-aperture image, I_j(x, y). This field’s digital inverse Fourier transform should create a masked spectrum with a clearly visible window function, W_j(k_x, k_y), modulating the sample spectrum via Eq. (2). Any deviation in the position and size of the computed window W_j(k_x, k_y) from the known SLM pixel map will inform us of how the optical setup’s Fourier transformation differs from an exact digital computation. We may computationally account for any such deviation (i.e., system misalignment) by digitally updating the window’s assumed position and size with a more accurate reflection of the imperfect optical setup.

In practice, although we assume accurate knowledge of each sub-aperture’s centered location (c_xj, c_yj) and size ℓ_j when we extract and insert the updated spectra (e.g., Eq. (4)), these variables are not known a-priori. Moreover, each sub-aperture image only measures the complex field’s amplitude, thus preventing us from adopting the above simple misalignment correction scheme. However, just like last section’s unknown aberration function A(k_x, k_y), we can determine the 3 · n representative variables c_jx, c_jy and ℓ_j from measurements of the field’s intensity with simulated annealing. Our “geometric distortion” annealing process proceeds as follows. At the j^th iteration, we construct T candidate window functions $W_j^t$, defined for 1 ≤ t ≤ T as,

Figure 7 presents a simulation of the geometric distortion SA-OFC process. Here, we again use our familiar noisy optical system, now free of Zernike model aberrations but geometrically distorted in a radial manner. Each sub-aperture center is displaced from its original location by δc_j = 0.05c_j, suggesting a Fourier plane that is axially offset from its true location by 250 μm. Figure 7(b) illustrates this distortion, where the blue circles denote the distortion-free rectilinear centers c_j (i.e., the algorithm’s assumed centers) and the x’s mark their actual locations, c_j + δc_j.

 

Fig. 7 SA-OFC removes geometric distortions. (a) SLM sub-aperture centers at the Fourier plane (blue dots) may be distorted by optical misalignments to unknown positions. During algorithm iteration, each estimated center is randomly perturbed by vector (δ_x, δ_y) as we generate and compare trial images to measured data. (b) Ground-truth sub-aperture centers for Fig. 5’s simulation are here radially offset to simulate a Fourier plane displaced by 100 μm axially. The actual center of each displaced sub-aperture is marked with an ‘x’. Initializing each sub-aperture center on an assumed rectilinear grid (blue dots), the simulated annealing process draws these estimates close to the actual centers after 10 iterations. (c) Error decreases with iteration similar to Fig. 6(a).

Download Full Size | PDF

Iterating q = 20 times accurately identifies the unknown misaligned sub-aperture centers. Using parameters T = 8, r_x = 10 pixels and r_w = 0 pixels, we reduce an initial average center offset from 125 μm to 14 μm. Removal of these geometric effects qualitatively improves the final image similar to Fig. 5. As with using SA to correct for misalignments in ptychography [27], we also empirically find that holding ( ${\delta }_x^t$, ${\delta }_Y^t$) to zero for the first few iterations helps improve convergence. Figure 7’s simulation holds the search radius to 0 for the first 5 iterations, as SA’s blue error curve in Fig. 7(c) indicates. As with our aberration correction example in Fig. 6, this curve’s final error is much lower than OFC without simulated annealing (red) and is bounded below by running OFC with the correct geometric distortions known a-priori (green).

3.3. The complete OFC algorithm

We may also simultaneously determine and remove low-order Zernike polynomial-based aberrations and geometric misalignments in one SA loop. By searching over both candidate aberrations with Eq. (7) and sub-aperture positions with Eq. (12) in a parallel manner, we may use the SA process summarized in Eq. (7) – (11) to search over a larger sub-space of unknown camera configurations. However, as we increase the dimension of this parameter search space, computation time will grow exponentially. For example, to search over both T different possible defocus settings and T different sub-aperture positions with the same amount of coverage, each loop must now test T² candidate states. If we would also like to search for the effects of x and y astigmatism, the number of candidate tests jumps to T⁴. For a p-dimensional parameter search space, the number of candidate tests becomes T^p. As our constructed optical setup was primarily influenced by sub-aperture position shift and low-order aberrations, this exponential scaling does not become a major concern, as the following experiment demonstrates. However, for systems suffering from many possible equally-weighted aberrations and geometric imperfections, alternative iterative strategies like conjugate gradient descent [28] will prove more efficient.

4. Experimental results

Our experimental OFC setup closely matches the optical conditions simulated above. We created a 4f camera using two bi-convex lenses (diameter= 25.4 mm, f= 75 mm, Thorlabs AC254-075). Without any additional corrective elements, it is easy to confirm that such a single-element large F-number lens (∼F/3) exhibits significant off-axis aberrations. We use an amplitude-modulating SLM (Epson HDTV LCD, BBS Bildsysteme 1920×1080 pixels placed between two crossed polarizers) to create the Fourier plane’s shifting sub-aperture masks. Used as an amplitude modulator, imperfections in the SLM’s surface flatness or angle-dependent response do not significantly alter detected intensities. One may also add parameters to the SA search process to estimate and correct for such unknown imperfections, albeit at the expense of increased computation time. Furthermore, the SLM’s imperfect pixel fill factor generates multiple diffractive orders. Due to our interest in just the amplitude-modulated zeroth order, we allow all higher orders to exit the optical system.

At the image plane, we place a CMOS detector (Prosilica GX-1920 with 4.54 μm pixels). We vary the central 120² pixels of the SLM, defining the camera’s full aperture width as 2.76 mm with its 23 μm pixels. This finite aperture range turned out to be an optimal compromise between resolution-matching with our detector pixel size and efficiently computing an SAOFC reconstruction. In the presence of large aberrations, such as those created by large-angle rays from our single-element lenses, SA must search over a proportionally increased parameter range, which slowly becomes intractable in a high-dimensional parameter search space.

We illuminate our sample with a collimated 632 nm laser beam (no spatial filtering used). Unlike Fourier ptychographic microscopy [3], the illumination field’s specific shape and phase is not critical and can remain unknown. The spatial coherence length of the illuminating beam must be as wide as the sample at the object plane (approximately 5 mm here). Likewise, the source’s temporal coherence must follow the quasi-monochromatic condition in [35], although a spectral bandwidth full width half maximum of approximately 10 nm has been observed as sufficiently narrow [13,40]. Both conditions may be achieved with an LED placed sufficient far behind the sample. Alternative algorithms exist to incorporate the effects of a partially coherent source and further reduce these illumination requirements [41].

The SLM’s displayed sub-aperture is a square with side length ℓ = 0.92 mm, which we sequentially shift laterally n = 9 times in x and y across an L × L =2.76 mm× 2.76 mm square Fourier plane area. In each step, we shift the aperture by 0.28 mm, raster scanning it until it passes through all n² = 81 unique aperture locations. This leads to a sub-aperture overlap of approximately 70% and an expected resolution gain of 3. To ensure that each of the 81 captured images is properly exposed, we take 3 snapshots at each sub-aperture position with an exposure sequence of .01, .1 and 1 second and combine these 3 images via high-dynamic range (HDR) processing [42]. The HDR process may be omitted by using a high bit-depth detector or if larger sub-apertures are utilized. Finally, we experimentally measure the SLM’s optical density as 22 (i.e., the central sub-aperture region is 22 times as bright as the outer apodizing region), making b = 0.045 in Eq. (1).

Figure 8(a) displays an example sub-aperture image captured with the aperture centered on the optical axis (i.e., c_j = 0). First, we attempt OFC reconstruction without simulated annealing as outlined in Section 2. The reconstructed image intensity after q = 20 iterations is in Fig. 8(b). Note the resolution is improved significantly above a single sub-aperture image, but artifacts remain primarily due to a misaligned Fourier plane. To correct these misalignments, we then implement OFC with simulated annealing. We use Section 3.3’s “complete” SA-OFC algorithm to jointly search over aberrations and geometric misalignments. We adopt an aberration model that contains the first five primary aberrations beyond linear tilt in x and y – defocus, x/y astigmatism and x/y coma: $A(\rho ,\theta )=\text{exp}\left(2\pi i{\sum }_{l=4}^8{a}_l{W}_l(\rho ,\theta )\right)$. Applying Eq. (7)’s search but now over these 5 orthogonal aberration modes, in addition to Eq. (12)’s geometric sub-aperture correction, creates a total search space dimension of p = 6 (T⁶ candidate tests).

 

Fig. 8 Experimental results from an OFC setup. We first image an Air Force resolution target to test the algorithm’s spatial resolution performance. (a) A single sub-aperture image exhibits low spatial resolution. (b) The OFC algorithm without annealing recovers a sharper image, but still contains artifacts. (c) The SA-OFC algorithm further improves the output field s_q(x, y)’s spatial resolution (see text for algorithm parameters), as highlighted by traces in (d).

Download Full Size | PDF

Setting the SA-OFC algorithm parameters to q = 20 iterations, T = 2, r_a = 10, r_x = 8, r_w = 1 and with linearly reducing α’s, we obtain Fig. 8(c)’s intensity image. Figure 8(d) plots a line trace through each image’s group 4, element 6. Both the sub-aperture image (a) and uncorrected output (b) fail to resolve this feature (17.54 μm-width). The aberration-corrected image (expected spatial frequency cutoff at 17.2 μm⁻¹) does resolve this group, confirming SA-OFC’s ability to restore a potentially misaligned camera system close to its native resolution.

The experimental system’s aberration map A_q(k_x, k_y) and Fourier plane misalignments, simultaneously recovered with the image in Fig. 8(c), are shown in Fig. 9(a) and Fig. 9(b), respectively. We note the final aberration function’s modal weights decrease with mode number, as expected for most lens systems. Linear x and y-tilt aberrations are partially corrected for by the geometric Fourier plane realignments and are thus left out of A_q(k_x, k_y) for computational efficiency. Figure 9(b) shows that geometric realignments are primarily a linear function along the vertical direction, as highlighted by its plots of average shift along k_x and k_y (average δc_x and δc_y). This matches observable changes in reconstructed vertical image features without (Fig. 8(b)) and with (Fig. 8(c)) SA. Without SA, the large unknown vertical misalignments cause the OFC algorithm to inaccurately combine each sub-images vertical spatial frequencies into a vertically blurred reconstruction. SA correctly identifies and accounts for this these large misalignments. Finally, we plot the algorithm’s NMSE E_q in Fig. 9(c), calculated via Eq. (5). We also plot E_q after halving the geometric misalignment search radius to r_x = 4, offering an idea of the algorithm’s sensitivity to this type of free parameter. The addition of annealing (SAOFC) improves our reconstruction NMSE by a factor of 3–4, but falls short of perfect recovery. Future efforts can increase reconstruction quality by including higher-order aberration and misalignment terms in the SA-OFC model, which will also increase computational complexity.

 

Fig. 9 (a) Aberration map recovered simultaneously with Fig. 8(c)’s image. (b) Geometric misalignments of our 4f setup’s Fourier plane also simultaneously recovered by the SAOFC algorithm. Average shifts for each row/column are plotted on the side/top. (c) Plot demonstrating the algorithm’s error reduction with iteration. After 10 iterations, SA-OFC’s NMSE is 3 times lower than the direct OFC algorithm without annealing.

Download Full Size | PDF

Next, to demonstrate that the OFC procedure can also acquire quantitative phase, we image a monolayer of polystyrene microspheres (diameter = 117 μm, index n_sphere = 1.594 at λ = 632 nm) coated on a microscope slide and immersed in oil (index n_oil = 1.5915). We use the same SA-OFC algorithm parameters outlined above to converge upon a high-resolution complex sample reconstruction. A phase map containing five microspheres is in Fig. 10(a). We take a line trace through one microsphere’s phase (dashed line) and plot its optical thickness h computed from the measured phase map Δϕ(x) via the following equation: $h=\frac{\lambda }{2\pi }\mathrm{\Delta }\varphi {\left(n_{\mathit{sphere}}-n_{\mathit{oil}}\right)}^{-1}$. This blue curve in Fig. 10(b) closely matches the optical thickness of a perfect sphere of known microsphere diameter, which here in green equals the length of the vertical chord connecting the top and bottom arcs of a 117 μm diameter circle. Minor background phase variations remain (∼0.2 radian peak change), which can be minimized with additional parameter tuning or a modified aperture coding strategy. Figure 10(d) includes the microsphere image’s computed aberration map, indicating successful convergence following two observations. First, the aberration map exhibits a similar structure to the resolution target’s shown in Fig. 9(a), apart from a constant phase offset. Second, microspheres recovered without aberration correction in Fig. 9(e) do not closely match the expected thickness profile, thus pointing towards the necessity of an aberration and misalignment correction strategy to ensure ptychography-based recovery schemes like OFC remain quantitatively accurate.

 

Fig. 10 SA-OFC for quantitative phase. (a) Recovered phase map of 5 microspheres with aberration correction. (b) Line trace through one sphere demonstrates quantitative accuracy. (c) Example captured images. (d) Simultaneously recovered aberration map exhibits a similar structure as Fig. 9(a), as expected for the same optical setup. (e) Phase map recovered by the OFC algorithm without simulated annealing, where low-order aberrations clearly impact the fidelity of the reconstructed phase map.

Download Full Size | PDF

5. Discussion and conclusion

We have demonstrated how a sequence of low-resolution images can be computationally processed into a full-resolution amplitude and phase image while simultaneously extracting camera aberrations. While our experimental imaging setup offers proof-of-concept aberration removal from a series of simple lenses, several steps may help move the OFC concept towards a more practical device. First, the multiple orders generated by the OFC setup’s transmissive SLM are undesirable. A reflective optical modulator such as liquid crystal-on-silicon can help avoid this effect. Second, while we expect minimal variation in the SLM’s angle-dependent amplitude response at our low angles of incidence (< 9°) [43], future increased-numerical aperture OFC setups should consider this effect. Third, the current OFC setup’s detector pixel size oversamples each captured image. A fully optimized OFC setup will match its detector pixel size to the width of each sub-aperture image’s intensity PSF. Since our final OFC reconstruction’s spatial frequency passband is 3 times as wide as each detected sub-image’s passband, an optimized setup will thus reconstruct an image containing 3 times the number of detector pixels along both x and y. This may of greatest use in detector-limited imaging scenarios.

Finally, modifying the selected SLM coding strategy is an important area of future improvement. Due to the limited sub-aperture size, required exposure times are quite long. Moving to an alternative aperture coding strategy to increase light throughput, such as displaying multiple sub-apertures per exposure, can immediately address this shortcoming. Phase-only modulation may also help decrease total capture time, both by increasing per-image light throughput and potentially reducing the total number of required acquisitions following [17]’s demonstration. Likewise, we may minimize artifacts observed in our reconstructed sample phase (e.g., the cubic phase in simulation) with addition of a random offset to each sub-aperture center. This is a direct and promising change to improve quality. Finally, the current coding strategy offers a higher total number of overlapped recording positions at lower spatial frequencies. Alternate sub-aperture arrangements may improve recovery of samples with significant energy in high spatial frequencies.

Two primary shortcomings currently limit the performance of the SA-OFC algorithm. First, as already noted above, computational scaling issues require the annealing search to consider only a small number of Zernike polynomial aberration coefficients and misalignment parameters. Second, instabilities are introduced when solving for multiple parameters that may not be orthogonal. Although each Zernike decomposition polynomial is mutually orthogonal and thus does not confuse the annealing’s movement towards minimal error, the geometric alignments and Zernike polynomials are not necessarily, which might lead to algorithm stagnation. Future work will examine more robust methods to search within a high-dimensional unknown parameter space. Possibilities include gradient descent or a maximum likelihood-based solver, where the log of the aberration function presents itself as a simple linear sum. These alternative strategies may additionally benefit from a new aperture coding strategy.

Besides system and algorithm improvement, several applications may directly benefit from our new imaging strategy. By acquiring an unaberrated quantitative sample phase map, it is possible to digitally refocus through thick or multi-planar objects. With aberration correction, an adaptive optics-inspired approach may additionally remove sample-induced distortions and blurring to help sharpen Fourier ptychography’s performance in this thick and/or turbid object regime. Also, unlike the original Fourier ptychography system, the OFC setup may extend to improve the SBP of fluorescence-based images. These two new capabilities will help extend ptychographic image resolution enhancement into areas that have not yet benefitted much from its recent invention, and will hopefully enable as yet unachieved system space bandwidths.