We present an imaging procedure that simultaneously optimizes a camera’s resolution and retrieves a sample’s phase over a sequence of snapshots. The technique, termed overlapped Fourier coding (OFC), first digitally pans a small aperture across a camera’s pupil plane with a spatial light modulator. At each aperture location, a unique image is acquired. The OFC algorithm then fuses these low-resolution images into a full-resolution estimate of the complex optical field incident upon the detector. Simultaneously, the algorithm utilizes redundancies within the acquired dataset to computationally estimate and remove unknown optical aberrations and system misalignments via simulated annealing. The result is an imaging system that can computationally overcome its optical imperfections to offer enhanced resolution, at the expense of taking multiple snapshots over time.
© 2014 Optical Society of America
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) , 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 , 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 107 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 ). 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 . 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  and SLM amplitude modulation  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  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 .
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  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 . 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 . 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.
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 n2 different equally spaced locations arranged on a 2D rectilinear grid. The jth sub-aperture will be displaced from the origin by a two-dimensional vector cj = (cxj, cyj) for 1 ≤ j < n2. 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 . All simulations and experiments in this work use n = 9, requiring a total of n2 = 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 jth SLM sub-aperture’s transmission function Wj as,Equation (1) states the SLM’s modulation is a biased rect function of width ℓ and center cj in the Fourier plane.
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)] = Ŝ(kx, ky), where ℱ is the Fourier transform operation and we will refer to Ŝ as the sample’s spectrum. After passing through the jth sub-aperture window, Wj, the exiting field is the product Ŝ(kx, ky)Wj(kx, ky) assuming the SLM is thin. This modulated field then propagates to the image plane detector, where to a first-order approximation the jth 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 , the sample here may be thick, optically reflective, and under arbitrary quasi-coherent illumination. While the aperture scanning procedure in  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, Ŝ0(kx, ky), 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 Ŝ0 with one of our SLM sub-apertures from Eq. (1). Starting with j = 0 and continuing until j = n2 − 1, we compute the spectrum product P̂j(kx, ky) = Ŝ0(kx, ky)Wj(kx, ky). For an SLM with perfect contrast, this corresponds to selecting a windowed region of Ŝ0, 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: pj(x, y) = ℱ[Pj((kx, ky)]. Third, we replace the amplitude of the resulting Fourier transform pj with the known amplitude from the jth detected image to form p′j:
Fifth, we update our sample spectrum estimate with P̂′j. While we accounted for the SLM’s imperfect modulation when extracting P̂j from Ŝ0, 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:4] and iterated projections  that do not strictly enforce an optical system’s support constraint every iteration. While alternative update strategies are certainly possible, we empirically determined that Eq. (4) leads to quick and accurate algorithm convergence. After updating all n2 overlapping sub-aperture areas of Ŝ0, we proceed to a second iteration by advancing our spectrum estimate: Ŝ1 ← Ŝ0 and resetting j = 0. After q iterations, we inverse Fourier transform the final recovered spectrum Ŝq(kx, ky) to recover our high-resolution sample estimate, sq(x, y).
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 = 10002 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.
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) = α(x3 + y3), 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 . 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) Eq between the target sample U(x, y)’s and reconstruction sq(x, y)’s amplitudes offers a useful metric to examine algorithm performance, as also used in :Figure 3 includes plots of the NMSE as a function of 15 iterations, where one iteration updates all 9 × 9 overlapping sub-aperture regions. Each curve represents a different amount of zero-mean Gaussian noise (variance σ2 = 0, 0.5%, and 2% of the maximum signal value). As expected, the algorithm is error-reducing and convergence accuracy decreases with increased noise for both samples. In addition, the blood cell sample’s highly varying phase leads to approximately 10 times lower NMSE than the slowly-varying cubic phase sample. Alternative aperture coding strategies may help reduce this performance gap , including addition of a random offset to the center of each sub-aperture to break scan pattern symmetry. But as with most phase-retrieval based solvers, we expect OFC’s convergence to remain sample dependent.
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 .
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, Ij(x, y), and the corresponding image that would result from detecting our complex windowed spectrum estimate, P̂(kx, ky). 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(kx, ky), 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 alWl(ρ, θ)), where Wl are the orthogonal Zernike polynomials and al are the associated weights . 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,Eq. (6) as our aberration model calls upon prior knowledge that this simple example system is primarily susceptible to defocus. In the following section, we discuss how to correct for an aberration function composed of multiple higher-order Zernike polynomials.
In addition to an estimated spectrum Ŝ0, simulated annealing OFC (SA-OFC) also initializes an estimated aberration map A0 (Fig. 4). We always initialize with A0(kx, ky) = 1. SA-OFC then begins with j = 1 and computes T different candidate aberration functions at the jth loop:Eq. (9)’s merit function only considers each ’s amplitude. We use the error-minimizing aberration perturbation as our annealing search’s starting point for the next (j + 1)th sub-aperture image: Eq. (8) to find , our optimal aberrated image field estimate. We then constrain with our measured intensities following Eq. (3): . After Fourier transforming ψ′j into Ψ′j, we are then ready to update our unaberrated sample spectrum estimate, Ŝj. To remove the effects of aberrations, we adopt a strategy common to prior algorithms like ePIE  and EPRY  and effectively divide out the aberration function estimate Aj from Ψ′j: Eq. (11)’s modified notation includes the iteration number in parenthesis and indicates update by the jth sub-image with a subscript. This random search and update process is repeated for all 0 ≤ j ≤ n2 − 1 sub-apertures to complete one iteration. To encourage the SA-OFC algorithm’s convergence, we linearly reduce the search radius ra after each iteration: ra ← ra − αa, where typically αa ≈ ra/q. Iteration continues for q loops to jointly solve for our final camera aberration map A(q) and sample spectrum solution Ŝ(q).
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(kx, ky) 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(kx, ky) 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 ra = d/2. The recovered aberration map Aq’s defocus coefficient dq differs from ground truth by Δd = (dq − d)2 = 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 . It is also sensitive to the selected search radius ra. All simulations here use ra = d/2 and linearly decrease it to zero by the qth iteration with αa = ra/(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.
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 = 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 A0(kx, ky) (“known”, green) to recover an almost exact solution. Here, instead of using Eq. (10) to estimate the aberration function each iteration, we simply set Aj+1(kx, ky) = A0(kx, ky), 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.
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(kx, ky). Again, each plot point is an average over 5 independent runs with σ2 = 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 . 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, Ij(x, y). This field’s digital inverse Fourier transform should create a masked spectrum with a clearly visible window function, Wj(kx, ky), modulating the sample spectrum via Eq. (2). Any deviation in the position and size of the computed window Wj(kx, ky) 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 (cxj, cyj) 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(kx, ky), we can determine the 3 · n representative variables cjx, cjy and ℓj from measurements of the field’s intensity with simulated annealing. Our “geometric distortion” annealing process proceeds as follows. At the jth iteration, we construct T candidate window functions , defined for 1 ≤ t ≤ T as,Eq. (9). We save a map of all updated centers, which serve as the starting window positions and sizes during the next cycle through all n windowed images (i.e., as the algorithm proceeds from the 1st to qth iteration).
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 δcj = 0.05cj, 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 cj (i.e., the algorithm’s assumed centers) and the x’s mark their actual locations, cj + δcj.
Iterating q = 20 times accurately identifies the unknown misaligned sub-aperture centers. Using parameters T = 8, rx = 10 pixels and rw = 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 , we also empirically find that holding ( , ) 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 T2 candidate states. If we would also like to search for the effects of x and y astigmatism, the number of candidate tests jumps to T4. For a p-dimensional parameter search space, the number of candidate tests becomes Tp. 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  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 1202 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 , 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 , 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 .
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 n2 = 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 . 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., cj = 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: . 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 (T6 candidate tests).
Setting the SA-OFC algorithm parameters to q = 20 iterations, T = 2, ra = 10, rx = 8, rw = 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−1) 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 Aq(kx, ky) 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 Aq(kx, ky) 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 kx and ky (average δcx and δcy). 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 Eq in Fig. 9(c), calculated via Eq. (5). We also plot Eq after halving the geometric misalignment search radius to rx = 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.
Next, to demonstrate that the OFC procedure can also acquire quantitative phase, we image a monolayer of polystyrene microspheres (diameter = 117 μm, index nsphere = 1.594 at λ = 632 nm) coated on a microscope slide and immersed in oil (index noil = 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: . 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.
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°) , 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 ’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.
The authors acknowledge funding support from the National Institutes of Health (grant no. 1DP2OD007307-01) and Clearbridge Biophotonics Pte Ltd., Singapore (Agency Award no. Clearbridge 1).
References and links
1. A. W. Lohmann, R. G. Dorsch, D. Mendlovic, Z. Zalevsky, and C. Ferreira, “Space-bandwidth product of optical signals and systems,” J. Opt. Soc. Am. A 13, 470–473 (1996). [CrossRef]
2. A. W. Lohmann, “Scaling laws for lens systems,” Appl. Optics 28, 4996–4998 (1989). [CrossRef]
3. G. Zheng, R. Horstmeyer, and C. Yang, “Wide-field, high-resolution Fourier ptychographic microscopy,” Nature Photon. 7, 739–745 (2013). [CrossRef]
4. J. R. Fienup, “Phase retrieval algorithms: a comparison,” Appl. Optics 21, 2758–2769 (1982). [CrossRef]
5. V. Elser, “Phase retrieval by iterated projections,” J. Opt. Soc. Am. A 20, 40–56 (2003). [CrossRef]
9. T. R. Hillman, T. Gutzler, S. A. Alexandrov, and D. D. Sampson, “High-resolution, wide-field object reconstruction with synthetic aperture Fourier holographic optical microscopy,” Opt. Express 17, 7873–7892 (2009). [CrossRef] [PubMed]
10. S. Chowdhury and J. Izatt, “Structured illumination diffraction phase microscopy for broadband, subdiffraction resolution, quantitative phase imaging,” Opt. Lett. 39, 1015–1018 (2014). [CrossRef] [PubMed]
13. S. Dong, R. Horstmeyer, R. Shiradkar, K. Guo, X. Ou, Z. Bian, H. Xin, and G. Zheng, “Aperture-scanning Fourier ptychography for 3D refocusing and super-resolution macroscopic imaging,” Opt. Express 22, 13586–13599 (2014). [CrossRef] [PubMed]
16. M. P. Lee, G. M. Gibson, R. Bowman, S. Bernet, M. Ritsch-Marte, D. B. Phillips, and M. J. Padgett, “A multi-modal stereo microscope based on a spatial light modulator,” Opt. Express 21, 16541–16551 (2013). [CrossRef] [PubMed]
17. F. Zhang and J. M. Rodenburg, “Phase retrieval based on wave-front relay and modulation,” Phys. Rev. B 82, 121104(R) (2010). [CrossRef]
18. C. Liang, T. Lin, B. Wong, C. Liu, and H. Chen, “Programmable aperture photography: multiplexed light field acquisition,” ACM Trans. Graph. 27(3), 55 (2008). [CrossRef]
19. K. Marwah, G. Wetzstein, Y. Bando, and R. Raskar, “Compressive light field photography using overcomplete dictionaries and optimized projections,” ACM Trans. Graph. 32(4), 46 (2013). [CrossRef]
20. P. Kner, B. B. Chhun, E. R. Griffis, L. Winoto, and M. G. L. Gustafsson, “Super-resolution video microscopy of live cells by structured illumination,” Nat. Methods 6(5), 339–342 (2009). [CrossRef] [PubMed]
21. R. Fiolka, M. Beck, and A. Stemmer, “Structured illumination in total internal reflection fluorescence microscopy using a spatial light modulator,” Opt. Lett. 33(14), 1629–1631 (2008). [CrossRef] [PubMed]
22. S. Quirin, D. S. Peterka, and R. Yuste, “Instantaneous three-dimensional sensing using spatial light modulator illumination with extended depth of field imaging,” Opt. Express 21(13), 16007–16021 (2013). [CrossRef] [PubMed]
23. S. Sarder and A. Nehorai, “Deconvolution methods for 3-D fluorescence microscopy images,” IEEE Sig. Proc. Mag. 23(3), 32–45 (2006). [CrossRef]
24. M. Broxton, L. Grosenick, S. Yang, N. Cohen, A. Andalman, K. Deisseroth, and M. Levoy, “Wave optics theory and 3-D deconvolution for the light field microscope,” Opt. Express 21(21), 25418–25439 (2013). [CrossRef] [PubMed]
25. J. Goodman, Introduction to Fourier Optics (McGraw-Hill, 1996).
27. A. M. Maiden, M. J. Humphry, M. C. Sarahan, B. Kraus, and J. M. Rodenburg, “An annealing algorithm to correct positioning errors in ptychography,” Ultramicroscopy 120, 64–72 (2012). [CrossRef] [PubMed]
28. A. Tripathi, I. McNulty, and O. G. Shpyrko, “Ptychographic overlap constraint errors and the limits of their numerical recovery using conjugate gradient descent methods,” Opt. Express 22, 1452–1466 (2013). [CrossRef]
29. B. Abbey, K. A. Nugent, G. J. Williams, J. N. Clark, A. G. Peele, M. A. Pfeifer, M. De Jonge, and I. McNulty, “Keyhole coherent diffractive imaging,” Nature Phys. 4, 394–398 (2008). [CrossRef]
30. I. Johnson, K. Jefimovs, O. Bunk, C. David, M. Dierolf, J. Gray, D. Renker, and F. Pfeiffer, “Coherent diffractive imaging using phase front modifications,” Phys. Rev. Lett. 100, 155503 (2008). [CrossRef] [PubMed]
31. R. G. Paxman, T. J. Schulz, and J. R. Fienup, “Joint estimation of object and aberrations by using phase diversity,” J. Opt. Soc. Am. A 9(7), 10721085 (1992). [CrossRef]
35. A. W. Lohmann, “Matched filtering with self-luminous objects,” Appl. Opt. 7(3), 561–563 (1968). [CrossRef]
36. O. Bunk, M. Dierolf, S. Kynde, I. Johnson, O. Marti, and F. Pfeiffer, “Influence of the overlap parameter on the convergence of the ptychographical iterative engine,” Ultramicroscopy 108, 481–487 (2008). [CrossRef]
37. M. Dierolf, P. Thibault, A. Menzel, C. M. Kewish, K. Jefimovs, I. Schlichting, K. von Koning, O. Bunk, and F. Pfeiffer, “Ptychographic coherent diffractive imaging of weakly scattering specimens,” New J. Phys. 12, 035017 (2010). [CrossRef]
39. R. Kingslake and R. B. Johnson, Lens Design Fundamentals: Second Edition (Elsevier and SPIE, 2010).
42. P. E. Debevec and J. Malik, “Recovering high dynamic range radiance maps from photographs,” In Proc. SIGGRAPH 97, ACM SIGGRAPH/Addison Wesley Computer Graphics Proceedings, Annual Conference Series, 369–378 (1997).
43. A. Lizana, N. Martin, M Estape, E. Fernandez, I. Moreno, A. Marquez, C. Iemmi, J. Campos, and M. J. Yzuel, “Influence of the incident angle in the performance of Liquid Crystal on Silicon displays”, Opt. Express 17(10), 8491–8505 (2009). [CrossRef] [PubMed]