Extended depth-of-field imaging and ranging in a snapshot

Traditional approaches to imaging require that an increase in depth of field is associated with a reduction in numerical aperture, and hence with a reduction in resolution and optical throughput. In their seminal work, Dowski and Cathey reported how the asymmetric point-spread function generated by a cubic-phase aberration encodes the detected image such that digital recovery can yield images with an extended depth of field without sacrificing resolution [Appl. Opt. 34, 1859 (1995)]. Unfortunately recovered images are generally visibly degraded by artifacts arising from subtle variations in point-spread functions with defocus. We report a technique that involves determination of the spatially variant translation of image components that accompanies defocus to enable determination of spatially variant defocus. This in turn enables recovery of artifact-free, extended depth-of-field images together with a two-dimensional defocus and range map of the imaged scene. We demonstrate the technique for high-quality macroscopic and microscopic imaging of scenes presenting an extended defocus of up to two waves, and for generation of defocus maps with an uncertainty of 0.036 waves. © 2014 Optical Society of America


INTRODUCTION
We report here a new modality of hybrid imaging that provides three-dimensional image ranging with an extended depth-offield (DOF) and with enhanced image quality.Perhaps the earliest demonstration of the benefit of hybrid imaging was the report by Häusler in 1972 [1] of combining time-sequential, swept-focus imaging of a deep object with post-detection image recovery to yield a sharp image with extended DOF.The average modulation-transfer function (MTF) of sweptfocus imaging exhibits no nulls and the absence of phase effects enabled simple recovery using coherent optical processing based on a photographic transparency, which provided more convenient recovery than digital computation at that time.The seminal 1995 paper by Dowski and Cathey [2] showed that a cubic optical phase function introduced into the exit pupil of an imaging system yields a point-spread function (PSF) that is approximately invariant to defocus, exhibits no nulls in the MTF, and therefore enables digital recovery of a high-quality image for an extended DOF [3].Control of focus-related aberrations is a major challenge in lens design, and so hybrid imaging using cubic or other antisymmetric phase functions such as trefoil [4] has been exploited for simplification of lens design and manufacture in miniaturization of zoom lenses [5] and thermal imaging [6], particularly for infinite-conjugate imaging.
Although the MTF for a cubic phase function is approximately invariant with defocus, there is translation of the PSF and strong phase modulation on the optical-transfer function (OTF) [7,8].Image recovery with a single kernel therefore introduces phase mismatches between the coding optical phasetransfer function (PTF) and the PTF of the digital filter used for image recovery, leading to range-dependent translation and image-replication artifacts [9].It is perhaps for these reasons that practical exploitation of this so-called wavefront coding (WC) technique appears to have slowed in recent years.
Attempts have been made to reduce the impact of these phase-induced artifacts: algorithms based on wavelet transforms have been used to estimate the magnitude of the imagereplication artifacts for a single image and to parametrically estimate the optimal kernel for artifact-free image recovery [10].For practical three-dimensional scenes, however, this requires accurate segmentation of scene components corresponding to objects at different ranges.Mo et al. proposed to use a deconvolution kernel based on the average optical PTF over the selected DOF together with a PTF-based optimization technique [11]; however, artifacts remain in simulated recovered images.Similarly, nonlinear filtering instead of the conventional Wiener filter has also succeeded in suppressing artifacts, but not eradicating them [12].
Radially symmetric phase functions have also been proposed for increase of DOF, including quartic and logarithmic phase functions [13] and the logarithmic asphere [14], and these enable recovery of images without phase-induced artifacts but with higher levels of noise amplification than for the antisymmetric functions [15] for reasons that appear to be fundamentally associated with the antisymmetry [7].
We report here theory and experimental results for a generalized technique we call Complementary Kernel Matching (CKM), which exploits the range-dependent translation of two images recorded with dissimilar imaging PSFs obtained from cubic phase functions, followed by inference of defocus and recovery of an image with the correct kernel and hence without translation or image-replication artifacts.By measuring translation, rather than artifact magnitude as used in [10], more robust determination of range is accomplished and furthermore this is achieved on a per-pixel basis, without the need for segmentation of image components.It is, however, necessary to record two images with dissimilar optical modulations and we present methods to achieve this in both time-sequential and snapshot modes.
The CKM technique involves both recovery of an extended-DOF image and simultaneously determination of a depth map.Estimation of a depth map from variations in image quality associated with a range-variant PSF has previously been demonstrated using shape from defocus of a conventional PSF using a sequence of defocused images [16]; Quirin and Piestun proposed a technique to estimate depth from two snapshots taken with two different engineered PSFs [17], where one PSF yields depth information and the other yields an extended-DOF image; Blanchard and Greenaway extracted depth information using a diffractive optical element to generate multiple defocused images on a single camera array and solved the intensity-transport equation to determine depth [18].When the scene consists of point sources, the form of the PSF is directly measurable and this enables range measurement in particle image velocimetry [19,20] and three-dimensional stochastic optical reconstruction microscopy (STORM) [21,22].The principle of the CKM technique we describe differs in that it exploits range-dependent translation of image pixels rather than blurring, offers the unique advantages of combining extended DOF with snapshot operation, is inherently achromatic, and is relatively simple to implement.
In the following sections, the theory of CKM is developed, simulation results are presented, and finally, an experimental evaluation of the technique and the obtained results are presented and analyzed.Although we consider here only the cubic function for generation of the encoding function, the principle is expected to be applicable to various other masks.CKM is based on the use of phase masks that produce images exhibiting detectable changes with defocus.The robustness of the technique benefits from the image-domain translation associated with some asymmetric phase functions, such as cubic functions, which is not present for radially symmetric or trefoil phase functions.

THEORY AND SIMULATION
Complementary Kernel matching imaging is based on recording two independent images with two phase masks providing complementary PSFs, that is, imaging kernels.We restrict our analysis initially to complex-conjugate pairs of the cubic phase function, ψ and ψ , where ψ exp2πiϕ and ϕ αx 3 y 3 which, when combined with defocus W 20 , provides the two pupil functions P exp2πiW 20 x 2 y 2 ϕ; (1) where α determines the strengths of the phase mask and x; y are the normalized pupil coordinates.This choice is motivated solely by the fact that the complex conjugate is easily generated by simply rotating the CPM half a revolution in the pupil plane.However, CKM requires only two PSFs which respond differently to defocus, and we present initially only a specific and idealized example of a more general principle.As will be seen, the technique is robust even when P and P − are not simply related by complex conjugation.Simple image recovery, using a Wiener filter for example, of the intensity distribution recorded using the phase functions P and P − , yields images exhibiting both artifacts and translation when the optical defocus W 20 is dissimilar to the defocus, W 20 , used for the recovery kernel.Translation, ρ, is parallel to the ξ η direction in the image plane and is proportional to W 2 20 − W 2 20 ∕α for the phase function ψ and proportionate to − W 2 20 − W 2 20 ∕α for ψ [2,7,9,23].The value of optical defocus, W 20 , used to record an image can thus be determined by identifying the matching W 20 that produces no displacement between the recovered images corresponding to ψ and ψ .In these circumstances, recovered images will also be free of phase-error artifacts.
More generally we write W 20 ξ; η to indicate that image defocus varies according to the range of scene components, and this leads to a spatially varying shift ρξ; η.Estimation of ρξ; η, based on a map of the disparity (as used to characterize stereopsis for example) between the images recorded with the phase masks ψ and ψ enables determination of a spatially varying W 20 ξ; η that enables recovery of an artifact-free image.Furthermore, W 20 ξ; η yields an estimate of the depth map of the scene.
If r W 20 ;W 20 and r W 20 ;W 20 − are the images restored using the deconvolving PSF corresponding to P and P − , with uniform defocus W 20 and image defocus W 20 ξ; η, the spatial-shifting properties lead to a disparity between them such that where ρξ; η is the scalar value describing a spatially variant displacement of scene components in the ξ η direction according to the mismatch between W 20 and W 20 ξ; η.
The effect of a spatially variant displacement of the recovered scene is illustrated for a simulated image of a spoke target in Fig. 1, for which W 20 ξ; η spatially varies with segment from zero to three waves as the orientations of the segment spokes varies from zero through 2π clockwise.Red arrows indicating disparities in r W 20 0;W 20 and r W 20 0;W 20 − are superimposed on the recovered images.
Additional to the shifts, phase-mismatch artifacts are apparent from mismatches between W 20 and W 20 .In consequence, Eq. ( 2) does not strictly hold.Nevertheless, calculation of the disparity map 2ρξ; η is possible since the shift is by far the dominant effect.Conversely, and more interesting, for spatially variant deconvolution with W 20 ξ; η W 20 ξ; η, ρξ; η vanishes and therefore phase artifacts are not present in Eq. ( 2).The recovered image is then, in the absence of noise, free of phase-induced artifacts.
Rather than calculating the disparity map which results from the recovery using uniform W 20 , it is more convenient to seek the spatially variant W 20 ξ; η that produces a null disparity map, since the artifacts will be suppressed as the solution is approached, that is, we calculate W 20 ξ; η that minimizes the difference between r W 20 ;W 20 ξ;η ξ; η and r W 20 ;W 20 ξ;η − ξ; η.Calculation of W 20 may be performed on a per-pixel basis, but to improve robustness, evaluation in a small neighborhood of each pixel is desirable to reduce sensitivity to noise.We perform this neighborhood evaluation using a smoothed representation of the difference image where G σ is a Gaussian low-pass filter.The standard deviation σ is a compromise between noise robustness and the size of the response function that tends to segment the image according to W 20 ξ; η.When W σ 20 ξ; η W 20 ξ; η ∀ ξ; η, r W 20 ;W 20 r W 20 ;W 20 − and each recovered image will be artifactfree.If σ is large, evaluation of Eq. ( 3) becomes less sensitive to noise and to scene features, but at the expense of lower lateral resolution for detection of spatial changes in defocus; for example, for discriminating between one object in front of another.
We reconstruct the deconvolved image by Wiener filtering based on the spatially varying PSF corresponding to the values of W 20 ξ; η that minimize the metric in Eq. ( 3).We denote this image obtained by r to indicate that a neighborhood of size σ is used and to highlight the field dependence on both the imaging and recovery kernels, W 20 ξ; η and W σ 20 ξ; η, respectively.In practice, this reconstruction strategy can result in local changes in contrast which manifest as slightly brighter or darker patches in the recovered image.This effect is removed by averaging the images obtained for a range of values of σ.These changes in contrast are quite abrupt and hence introduce nonphysical frequency components beyond the optical cutoff frequency, and these we attenuate using a low-pass filter with negligible effect on the image quality.The final reconstructed image can thus be written as where for clarity the ξ; η dependence of r and r − has been suppressed but is implicit, F •; ν c is an ideal low-pass filter, and ν c is the optical cutoff frequency.Note that reconstruction in Eq. ( 4) follows from the set of images r Simulation results obtained from the application of this algorithm for the recovery of simulated images are shown in Fig. 2. A reference diffraction-limited image is shown in Fig. 2(a) and blurring by uniform defocus W 20 4 is shown in Fig. 2(b).The image in Fig. 2(c) is obtained by conventional image recovery with W 20 0 as typically employed for WC, while the image in Fig. 2(d) was recovered using the CKM method.Detected images included zero-mean, white Gaussian noise with 46 dB signal-to-noise ratio.The variation with W 20 of the argument of the argmin function in Eq. ( 3) is depicted in Fig. 2(e), where the minimum indicates the correct detected defocus at W 20 W 20 4. By comparing Figs.2(c) and 2(d), it is apparent that through CKM, the phase-modulation artifacts that are generally observed in WC image recovery [9] have been eliminated, thus enabling a higher quality, artifactfree recovery to be obtained.Similar results were robustly obtained for a range of scene types and noise levels.
We illustrate the ability to recover images with varying defocus across the field of view in Fig. 3. Conventional in-focus images of a spoke target and fishing boat are shown in column (a); simulations of the effects of spatially varying defocus W 20 ξ; η for conventional imaging are shown in column (b), where W 20 ξ; η varies with orientation of the spokes for the spoke target (as for Fig. 1) and varies linearly with elevation from positive to negative defocus for the fishing boat; conventional image recovery of images recorded with a cubic phase mask and W 20 0 is shown in column (c) and with CKM recovery in column (d).The presence of artifacts for conventional WC image recovery and their absence in CKM is clearly evident.Furthermore, the CKM algorithm also yields a defocus map and hence also a range map.Although the raw defocus map W σ 20 ξ; η can be noisy where regions without texture give rise to ambiguity, restored image are nevertheless free of phase-induced artifacts.
The two images corresponding to ψ and ψ can be recorded time-sequentially, for example, by rotating a phase mask through an angle of π, or in a snapshot by using one of the configurations shown in Fig. 4. In Fig. 4(a), the pupil and its conjugate are implemented in two distinct optical paths and using one detector for each image, while in (b) a minor modification enables a single detector array to be used.It is important to note that CKM image recovery involves individual calibration of the PSFs for P and P − , and hence is un-affected if practical implementations, such as those depicted in Fig. 4, yield two PSFs that do not correspond to the ideal phase conjugation between the two phase functions.This will arise, for example, for small dissimilarities in the aberrations between  the imaging systems.Figure 4(c) depicts a scheme similar to that reported by Blanchard and Greenaway for producing conjugate wavefronts [18]: a dislocated grating located in the pupil of an imaging system yields three diffraction orders that are captured by a single imaging detector (the light from higher orders is lost).Encoding of ψ into the grating spatial modulation yields replicated images at the detector encoded with PSFs due to ψ and ψ in the 1 and −1 diffraction orders.The zero-order image corresponds to a nonmodulated pupil and is not used.The 1 diffraction orders thus yield equivalent images to those provided by Figs.4(a) and 4(b) but with the advantage of a common lens for both images recorded for ψ and ψ , thus simplifying the calibration process.This method offers a convenient simultaneous recording of images for ψ and ψ phase masks with a single detector array, but is inefficient in its use of light and of detector pixels (the zero order is not used), and in its simplest form requires imaging with narrow-band light, although some achromatization is possible [24].
For our proof-of-concept demonstration, we report here a time-sequential implementation of ψ and ψ using a Spatial Light Modulator (SLM) [25], which proves convenient for calibration.

EXPERIMENTAL DEMONSTRATION A. Experimental Setup
We demonstrate the CKM technique using a conventional finite-conjugate imaging configuration employing an SLM to time-sequentially implement the phase functions ψ and ψ and indeed the more general phase functions P and P − that include defocus, W 20 .This is intended as a proof of principle, and we will report a snapshot implementation in the near future.A mechanical time-sequential implementation is also possible: for example, by rotation of a refractive cubic-phase mask through π radians about the optic axis to switch between ψ and ψ ; however, practical implementation introduces significant error into the measurement of ρ associated with uncertainties in the deviation of the images by the phase mask.Furthermore, accurate calibration of the variation of ρ and PSF with W 20 based on defocus of the system PSF can be more readily achieved with an SLM than by mechanical defocus of, for example, a pinhole.
The experimental setup is shown in Fig. 5. Using object and image distances of 350 and 1650 mm, an f/15 singlet lens of focal length 300 mm forms an in-focus image with a nominal magnification of approximately 5 on a Hamamatsu Orca CCD detector.A polarizer, quarter-wave plate, and analyzer are used with the SLM [25,26], yielding a maximum phase modulation of 3π∕2 with total amplitude modulation of <4% using illumination at a wavelength of 543 nm.An iris placed adjacent to the SLM ensures that the SLM is in the aperture stop of the system, so that phase coding is independent of field angle.
As described below, the CKM method is demonstrated for imaging of three-dimensional test targets involving significant defocus.Calibration of the variation of PSF intensity distribution and displacements ρ with W 20 , as required for accurate image recovery and estimation of W 20 ξ; η, was achieved by recording the image of a 1 μm pinhole located in the object plane of the imaging system.Two hundred PSFs were recorded for pupil functions P and P − implemented with the SLM, that is, for encoding functions for ψ and ψ , and for defocus varying equidistantly in the range −3 ≤ W 20 ≤ 3. The variations of ρ with W 20 , as determined by correlation of each PSF with the reference (in-focus) PSF for P and P − together with least-square fits of quadratic functions, are shown in Fig. 6.The quadratic fits were used to improve the estimate of the noise-free PSF position that was subsequently used in image recovery.
The best-fit quadratic functions correspond to an α that is approximately 21% greater for P than for P − .Measurement with a Shack-Hartman sensor of the phase fronts produced by the SLM yields best-fit cubic wavefronts for P and P − that correspond to values of α comparable to this measured asymmetry.Although this asymmetry could be removed by a pixel-wise calibration of the SLM phase function, it is more interesting to demonstrate that the technique is robust to such aberrations.

B. Results-Artifact Removal
As a pertinent demonstration of this technique, we present the application of CKM to imaging of microscope slides.Shown in Fig. 7 is an image of a section of the petiole of a leaf and in Fig. 8, a sample of seeds and a pine-leaf section.To provide appreciable range of defocus, the object in Fig. 7(a) was tilted to provide a linear variation in defocus in the vertical direction, and the two samples in Fig. 8(a) were separated by a glass slide of constant thickness.In each case, the defocus is W 20 ≈ 1.6.
Figures 7(b) and 8(b) show the corresponding images captured by means of a WC system with W 20 0 in the recovery, where the out-of-focus regions of the scene exhibit clear phase-error artifacts.For the images captured and reconstructed with a CKM system shown in Figs.7(c) and 8(c), however, the artifacts are effectively eliminated, which in turn enables image details to be more readily discerned.In addition, Fig. 8(d) shows the image recovered using a commercial focusstack algorithm (Helicon focus V5.3) generated using 201 images taken over a defocus range of −3 ≤ W 20 ≤ 3 at steps of 0.03 waves.Comparing Figs.8(c) and 8(d) one can observe that, using just two recorded images, the CKM technique demonstrates comparable image quality to that obtained with the 201-image focus-stack algorithm.In the near future, we will demonstrate this using a single snapshot.
The line plots shown in Fig. 8

C. Results-Depth Estimation
As stated previously, the CKM technique involves evaluation of defocus in a small region of an image so that it is possible to reconstruct a scene in 3D provided sufficient texture is present.We assessed the use of CKM image reconstruction for threedimensional ranging by imaging a calibration target consisting of a regularly spaced array of disks tilted at an angle of approximately 65°with respect to the nominal image plane, introducing defocus of 0.3 ≤ W 20 ≤ 2.0.Images of the calibration target captured with a conventional imaging system, WC system, and CKM system are shown in Figs.9(a), 9(b), and 9(c), respectively.
To reconstruct the scene in three dimensions, the defocus map was averaged over the region of each disk, and textureless  Red is for the WC system, blue is for the CKM system, and green is for the focus stack.areas which do not provide defocus information were masked.The calibration target was aligned such that there is effectively uniform defocus along each row of disks and a linear variation in defocus along each column.A ground truth slope was calculated by centroiding several disks and taking the ratio between their average horizontal to vertical separation.This was found to be 1.976 0.002, that is, approximately tan65° as expected.
Measurements of the range of each disk using the snapshot CKM measurement of defocus and based on a time-sequentially recorded focus stack were performed.For the focus-stack measurement, an image-sharpness criterion was used to select the plane of best focus, as the SLM was used to vary the focus of the image in steps of 0.03 waves.The range of each disk as computed by both the CKM and the focus stack are shown plotted in Fig. 10(b).Additionally, for illustration purposes, a three-dimensional reconstruction of the calibration target as computed by the CKM method is shown in Fig. 10(a).
The gradients of the least-square linear fits are 2.00 0.01 for the CKM technique and 1.95 0.02 for the focus-stack measurement.That is, the CKM measurement is 1.2% greater than the ground truth while the focus-stack measurement is 1.3% smaller.These results suggest that the ranging accuracy of CKM is comparable to the accuracy for focus stacking.For this specific target, the uncertainty in range for the CKM technique is 40 μm, which corresponds to an uncertainty in W 20 of 0.036 waves of defocus (slightly larger than the quantization step) and to 1.84% of the depth range.The corresponding uncertainty for the focus stack was found to be approximately twice as much.If implemented on a typical microscope with a numerical aperture of 0.65 and a typical DOF of 1 μm, this would correspond to a depth resolution of approximately 5% of the conventional DOF for the CKM and 10% for the focus stack.This implies that the performance of CKM is comparable to that of a focus stack.We therefore conclude that the CKM technique can be employed to capture threedimensional range-resolved images with extended DOF using just two data acquisitions and, in principle, using one of the techniques shown in Fig. 4 in a single snapshot.This is pertinent to a wide range of time-resolved imaging applications, including extended DOF microscopy or particle image velocimetry [19].

CONCLUSIONS
We report to our knowledge the first demonstration of artifactfree, extended-DOF imaging with simultaneous ranging using hybrid-imaging techniques.In previously reported approaches to hybrid imaging, a range of phase functions have been reported that tend to fall into two classes: either the antisymmetric cubic and trefoil masks (or qualitatively similar shapes) or symmetric masks [14,13], and these provide complementary advantages.While the former offers a superior trade-off between enhanced DOF and noise amplification [15], the spatial-phase effects introduced by the asymmetry introduce highly problematic artifacts and also range-dependent image shifts [9] that are absent for symmetric phase functions.Here we demonstrate that, by recording images with complementary OTF characteristics, it is possible to benefit from the enhanced DOF of an antisymmetric mask but without introducing image artifacts, and combined with the enhancement of three-dimensional ranging.This so-called complementary kernel imaging constitutes a new paradigm for hybrid imaging that builds on the early use of multiple kernels by Hausler [1].As reported here, CKM exploits translation of the PSF with defocus, and this translation provides robust image recovery.Translation does not occur with the radially symmetric phase masks [14,13] or the trefoil mask [4]; however, the presence of defocus-related phase-induced image artifacts for the trefoil and other masks does suggest CKM is applicable provided the image recovery is sensitive to defocus.In this case, CKM would not rely on translation but on a distortion and blurring combination.Indeed, two different phase mask types might be implemented, benefiting from complementary advantages.
The CKM method has been verified both by simulation and by experiment.Simulation results show suppression of restoration artifacts even at relatively low signal-to-noise ratios.Experimental results show that high-quality, artifact-free images can be obtained even for large defocus.In addition, CKM provides a means of 3D image reconstruction and, following proper calibration, range measurement.Future work will transfer this technique to high-numerical-aperture imaging, use of broadband refractive phase masks, and snapshot operation.

Fig. 3 .
Fig. 3. Images with varying amounts of defocus up to three waves; spoke target shows angular step change in defocus and boat image shows continuous vertical defocus change.(a) Diffraction-limited reference.(b) Conventional optics imaging.(c) WC with α 4. (d) CKM with α 4.
(e) correspond to a mid-height intensity horizontal profile taken along the dashed lines shown in Figs.8(b), 8(c), and 8(d).While strong artifactual oscillations are evident on the traditional WC profile, these are absent in the CKM images.The commercial focus-stack reconstruction algorithm employs some form of smoothing which we do not use in the CKM recovery, as can be observed by comparing Figs.8(d) to 8(c).This explains why some peaks are shallower than the corresponding peaks in the CKM recovery.

Fig. 8 .
Fig. 8. Step-defocus pine-leaf section and seeds.(a) Conventional imaging system, (b) WC system, (c) CKM system, and (d) focus stack.(e) Line profile taken along the lines with the corresponding color in (b)-(d).Red is for the WC system, blue is for the CKM system, and green is for the focus stack.