Video-rate compressive holographic microscopic tomography

Compressive holography enables 3D reconstruction from a single 2D holographic snapshot for objects that can be sparsely represented in some basis. The snapshot mode enables tomographic imaging of microscopic moving objects. We demonstrate video-rate tomographic image acquisition of two live water cyclopses with 5.2 μm spatial resolution and 60 μm axial resolution. ©2011 Optical Society of America OCIS codes: (090.1995) Digital holography; (070.0070) Fourier optics and signal processing. References and links 1. D. Gabor, “A new microscopic principle,” Nature 161(4098), 777–778 (1948). 2. E. N. Leith, and J. Upatnieks, “Reconstructed wavefronts and communication theory,” J. Opt. Soc. Am. 52(10), 1123–1130 (1962). 3. D. J. Brady, K. Choi, D. L. Marks, R. Horisaki, and S. Lim, “Compressive holography,” Opt. Express 17(15), 13040–13049 (2009). 4. E. J. Candès, J. K. Romberg, and T. Tao, “Stable signal recovery from incomplete and inaccurate measurements,” Commun. Pure Appl. Math. 59(8), 1207–1223 (2006). 5. E. J. Candes, and T. Tao, “Near-optimal signal recovery from random projections: Universal encoding strategies?” IEEE Trans. Inf. Theory 52(12), 5406–5425 (2006). 6. D. L. Donoho, “Compressed sensing,” IEEE Trans. Inf. Theory 52(4), 1289–1306 (2006). 7. L. Denis, D. Lorenz, E. Thiébaut, C. Fournier, and D. Trede, “Inline hologram reconstruction with sparsity constraints,” Opt. Lett. 34(22), 3475–3477 (2009). 8. C. Fournier, L. Denis, and T. Fournel, “On the single point resolution of on-axis digital holography,” J. Opt. Soc. Am. A 27(8), 1856–1862 (2010). 9. A. F. Coskun, I. Sencan, T.-W. Su, and A. Ozcan, “Lensless wide-field fluorescent imaging on a chip using compressive decoding of sparse objects,” Opt. Express 18(10), 10510–10523 (2010). 10. M. M. Marim, M. Atlan, E. Angelini, and J.-C. Olivo-Marin, “Compressed sensing with off-axis frequencyshifting holography,” Opt. Lett. 35(6), 871–873 (2010). 11. Y. Rivenson, A. Stern, and B. Javidi, “Compressive Fresnel holography,” J. Disp. Technol. 6(10), 506–509 (2010). 12. H. Kim, S.-W. Min, B. Lee, and T.-C. Poon, “Optical sectioning for optical scanning holography using phasespace filtering with Wigner distribution functions,” Appl. Opt. 47(19), D164–D175 (2008). 13. E. Y. Lam, X. Zhang, H. Vo, T.-C. Poon, and G. Indebetouw, “Three-dimensional microscopy and sectional image reconstruction using optical scanning holography,” Appl. Opt. 48(34), H113–H119 (2009). 14. J. Garcia-Sucerquia, W. Xu, S. K. Jericho, P. Klages, M. H. Jericho, and H. J. Kreuzer, “Digital in-line holographic microscopy,” Appl. Opt. 45(5), 836–850 (2006). 15. S. Seo, T. W. Su, D. K. Tseng, A. Erlinger, and A. Ozcan, “Lensfree holographic imaging for on-chip cytometry and diagnostics,” Lab Chip 9(6), 777–787 (2009). 16. J. Hahn, S. Lim, K. Choi, R. Horisaki, D. L. Marks, and D. J. Brady, “Compressive Holographic Microscopy,” in Biomedical Optics, OSA Technical Digest (CD) (Optical Society of America, 2010), paper JMA1. 17. R. E. Blahurt, Theory of Remote Image Formation (Cambridge University Press, 2005). 18. L. I. Rudin, S. Osher, and E. Fatemi, “Nonlinear total variation based noise removal algorithms,” Physica D 60(1–4), 259–268 (1992). 19. J. M. Bioucas-Dias, and M. A. T. Figueiredo, “A new twIst: two-step iterative shrinkage/thresholding algorithms for image restoration,” IEEE Trans. Image Process. 16(12), 2992–3004 (2007). 20. D. J. Brady, Optical Imaging and Spectroscopy (Wiley, 2009). 21. D. P. Kelly, B. M. Hennelly, N. Pandey, T. J. Naughton, and W. T. Rhodes, “Resolution limits in practical digital holographic systems,” Opt. Eng. 48(9), 095801 (2009). 22. A. C. Kak, and M. Slaney, Principles of Computerized Tomographic Imaging (SIAM, 2001). #137489 $15.00 USD Received 1 Nov 2010; revised 16 Jan 2011; accepted 3 Feb 2011; published 31 Mar 2011 (C) 2011 OSA 11 April 2011 / Vol. 19, No. 8 / OPTICS EXPRESS 7289 23. D. L. Marks, “A family of approximations spanning the Born and Rytov scattering series,” Opt. Express 14(19), 8837–8848 (2006). 24. J. Chae, and S. Nishida, “Integumental ultrastructure and colour patterns in the iridescent copepods of the family Sapphirinidae (Copepoda: Poecilostomatoida),” Mar. Biol. 119(2), 205–210 (1994).


Introduction
Holography characterizes the optical field, including the amplitude and the phase [1,2], on a 2D surface.Quasi-3D imaging has been achieved using holography via the so-called backpropagation mechanism, which numerically focuses the optical field at several depths, as if lenses focused the field to form an image.However, as discussed in [3], the backpropagation in digital holography recovers merely 2D fields at multiple depths as opposed to the true 3D scattering density that induces the 2D fields.It is possible, however, to form truly tomographic images from holographic data using nonlinear object estimation algorithms.This paper demonstrates the utility of these algorithms for 3D imaging of microscopic moving objects.In particular, we show that tomographic images can be acquired at full frame rate at the longitudinal and transverse resolution limits of the optical system.
Our group has demonstrated that 3D tomographic reconstruction can be achieved from a single holographic projection by using the theory of compressive sensing [4][5][6] provided that the representation of the 3D scattering density in some transformation basis is sufficiently sparse.This method, called compressive holography, has a potential for fast acquisition for 3D imaging.Several papers based on a similar philosophy have been published.In particular, Denis et.al. demonstrated success application of sparse reconstruction techniques to holographic 3D imaging and field-of-view expansion [7,8].Also, Coskun et al. applied a similar sparse estimation method for lenseless wide-field fluorescent imaging on a chip [9].A similar philosophy used in our compressive holography work has been applied to the 2D holographic imaging as well.Marim et al. and Rivenson et al. applied decompressive inference to frequency-shifting hologram [10] and digital Fresnel holography [11].It is also noteworthy that other 3D imaging methods such as optical sectioning [12] and Tikhonovregularized solution [13] are available in the literature even though all these methods are not strictly 3D imaging in that the solutions from the methods do not mathematically correspond to the true 3D scattering densities.
While the resolution of a microscopic system is determined by its numerical aperture (NA), the practical resolving power and reconstruction fidelity of compressive holography largely depend upon the effective NA that is determined by object feature sizes.It is desirable that the diffraction signals fully occupy the detector area such that the NA of an optical system is fully utilized.In practice, the combined effects of these factors may be reasonably estimated from the amounts of diffraction signals produced by the object features.Increasing diffraction signals, and equivalently the effective NA, may be achieved by several ways.For example, a microscopic objective can be used to increase the numerical aperture of a system by magnifying the field entering the objective.Of course, the objective tends to increase system volume and decrease field of view.
Highly compact microscopes may be constructed, as in digital holographic lensless microscopy [14,15], by eliminating optics between the scattering object and the focal plane.This paper focuses on such designs while also obtaining detector limited NA.We illuminate the object with a spherical wave incident field to magnify the diffraction signal on the FPA.Other more sophisticated incident fields may be utilized to further increase the detectable diffraction signals to improve the performance of our methods.Our physical system is similar to the digital in-line holographic microscope [14]; our contribution is to adapt compressive holography to microscopy and to demonstrate its effectiveness in single frame 3D tomography.Building on our first report of this approach in [16], this paper presents a theoretical framework for holographic tomography and demonstrates the utility of the approach for practical imaging.We note that the theoretical framework developed in this paper can be generalized to other microscopic geometries with small modifications which may utilize more complicated incident fields.
This paper is organized as follows.Section 2 describes the mathematical forward model and reconstruction methods for compressive holographic tomography.Section 3 describes an experimental system and presents experimental results.Section 4 concludes with suggestions for future development.

System model
We consider a microscopic system composed of a point source for illuminating the objects and the FPA that measures the Gabor hologram.Figure 1 shows a schematic of such a holographic microscope.The intensity of interference pattern ( , ; ) I x y z measured using the FPA placed at a propagation distance z can be written as where 0 U represents the incident field, and s U represents the complex field at distance FPA z scattered off of the objects.The incident field 0 U is generated by a point source and thus forms a spherical wave.Thus, .
Here,   , xy and   ,, x y z    represent the coordinates in the FPA plane and in the object space, respectively.  ,, x y z     is the scattering density of a 3D object and h is the Huygens- Fresnel point-spread function [17].Since the NA of the system is relatively large, the Fresnel approximation may be inaccurate.Thus, the use of the angular spectrum transfer function for propagation would be desirable.However, the effective NA is smaller than the numerical aperture of the system, which allows us to use the Fresnel approximation without suffering from the loss of numerical accuracy.The concept of the effective NA is discussed in details in Sec. 3.
The incident field   0 ,; U x y z is a spherical field originating from 0 z  and thus can be expressed using   ,; F h x y z .The scattered wave from the objects is represented as , ,

U x y z dx dy dz h x y z x y z h x x y y z z zz x z y z h x y z dx dy dz C z z h x x y y
. Here,

 
,; Considering that the FPA samples the scattered field, the last term in Eq. ( 1) can be accordingly discretized as where  denotes the sampling spacing in the FPA plane and   denotes the sampling spacing in the object space. . Here, ,, . Thus, considering Eq. (1) as in [3], Eq. ( 6) can be algebraically expressed as where f and g denote vectorized versions of Equation ( 7) is solved by finding f that minimizes the total variation (TV) [18].Minimizing the total variation (TV) is equivalent to enforcing the sparsity of f in the variational domain.The minimum TV estimate f can be obtained by solving In this equation, TV f denotes the total variation of f and is defined as , ( ) , where z f denotes a 2D transverse slice of the 3D object datacube, and , () the magnitude of the gradient at location   , xy in the z-th transverse slice.We solve the optimization problem in Eq. ( 9) by adapting the two-step iterative shrinkage/thresholding (TwIST) algorithm [19].
Separation of holographic signal from the background intensity terms ( e in Eq. ( 7)) and from the conjugate image is traditionally challenging for in-line holography.As previously reported in [3], we note that the estimator described in Eq. ( 8) achieves this separation by localizing the conjugate image appropriately outside the reconstruction volume and by localizing the background terms in the 0 z  plane.One expects this separation because the total variation of the propagated background terms is minimized in the 0 z  plane where these terms are real (thus eliminating variation in the complex component).

Experimental results and discussion
We constructed a holographic microscope using a He-Ne laser with 632.8nm wavelength as a light source.A microscope objective with 0.65NA manufactured by LOMO is used to generate a spherical wave for illuminating the sample.While the illumination source is not compact, we note that alternative sources using fiber or holographic components could achieve similar NA.The microscope objective is chosen to have a larger NA (0.65) than the microscope system NA (0.   data are used to avoid wraparound effects.This preprocessed data is defined as g in the algorithm in Eq. ( 8) to reconstruct the object scattering density f (or equivalently,  ).The theoretical limits of the resolution of the proposed system may be expressed as [20] , 2 where the NA is defined as the half width of the FPA over the distance of objects to the FPA plane.The NA is 0.27 in our experiments.From the equations, the theoretical resolution limits are estimated as 1.17 .In practice, however, the resolutions are limited by other factors such as the quantization of FPA, various sources of noise, and object feature sizes.These factors reflect on the hologram as discernible diffraction signals.Practical resolutions may thus be estimated by considering amounts of diffraction signals [21].Note that such resolutions can depend on objects because diffractions depend on object feature sizes.For example, one tail of the water cyclops on the left-hand side in Fig. 3       Note that the Born approximation is valid when the product of the index contrast and object size is less than one-quarter wavelength since the field inside the object is approximated using the incident field [22,23].For example, the platelets of Sapphirinidae that belongs to the same subclass as water cyclopses have a high refractive index as 1.8 [24].So the bodies of cyclopses are not in the range that the Born approximation is valid and its inner structure are not correctly reconstructed by our method based on the Born approximation.The reconstruction along the propagation direction is more inaccurate as in Fig. 6 since the holographic microscope is modeled as a linear system along the propagation direction.
Figure 7 presents an alternative visualization of the same 3D datacube.Figure 7(a) shows an image of the maximum intensity values along the propagation directions in 3D. Figure 7(b) represents a map of the axial positions corresponding to the maximum values in Fig. 7(a).From these two images, a range map is constructed as shown in Fig. 7(c).The map represents the HSV space.The hue (h) represents Fig. 7 Figure 8 shows the videos of the raw measurements and range maps associated with the compressive holography reconstructions.The raw measurements capture the movement of three water cyclopses in a 3mm depth container.In Fig. 8(a), the frame rate for the 2D hologram recording was set to 15 fps whose frames are played at the rate of 3 frames per second for readers' visual convenience. .A water cyclops around the center is swimming down over the range from 1.6mm to 4.5mm .The movement is illustrated by its color changes according to its depths.

Conclusion
We have demonstrated a video-rate compressive holographic microscopic tomography.The method has been applied to microscopy where an illumination wave is generated with a microscope objective to produce a spherical wave incident field.This framework for compressive holographic microscopy can be generalized to any computational holographic microscopy with moderate modifications.We have also analyzed both the theoretical and practical resolutions of the proposed microscope.Our use of the term "video-rate" refers to the acquisition time of the image, processing to reconstruct the tomographic data cube is far from real-time in our current implementation.Future work may focus on the use of application specific software and hardware to speed the reconstruction and visualization process.Alternative sampling strategies using multiple illumination wavelengths or more complex illumination patterns may be effective improving the longitudinal resolution and increasing the range of observable objects.

z
is the sampling spacing in the axial direction in the object space.The sample spacings in the Fourier domain k  satisfy different relations:2 can be rewritten as 27) to ensure uniform illumination intensity on the FPA.A Lumenera CMOS sensor records the hologram.The sensor has 1280 1024  resolution, 5.2 m  pixel pitch, 10bit digitization, and the maximum frame rate 15fps.Figure 2(a) shows a photograph of the microscope instrumentation.The relative positions of the microscope objective and the sample are adjusted by two 3-axis stages.Figure 2(b) shows a photograph of the container in which two live water cyclopses of roughly the same size are floating in water.The depth of the container is 3 mm.The distance from the FPA to the microscope objective for illumination is 10mm .

Fig. 2 .Figure 3
Fig. 2. Photographs of (a) the holographic microscope and (b) a water container in which water cyclopses are swimming.

Figure 4
Figure 4 compares the backpropagation reconstructions with compressive holographic reconstructions obtained with a single 2D measurement of the two live water cyclopses.Figures 4(a) and 4(b) show transverse slices of the backpropagation reconstructions at 3.31 in the reconstruction, the axial sample spacing, namely the distance between two adjacent transverse slices, is set to 60 m  .

Fig. 4 .
Fig. 4. A comparison of reconstructions at chosen axial positions: (a-b) the backpropagation reconstructions and by (c-d) the compressive holographic reconstructions using the data shown in Fig. 3.All the transverse slices of the reconstructions are sequentially shown in (e) obtained by the backpropagation method (Media 1) and in (f) obtained by the compressive holography method (Media 2) for the full range of water container.

Figure 5
Figure 5 compares the magnified tails, marked by rectangles in Figs.4(b) and 4(d), of (a) the backpropagation reconstruction and (b) the compressive holography reconstruction.The detailed comparison shows clear difference between the two reconstructions and their changes in the axial direction.

Fig. 5 .
Fig. 5.A comparison of the magnified tails, marked by rectangles in Figs.4(b) and 4(d), of (a) the backpropagation reconstruction and (b) the compressive holography reconstruction.

Fig. 7 .
Fig. 7. Images of (a) the maximum intensity values of the reconstructed density ( f ) along the propagation directions and (b) a map of the axial positions corresponding to the maximum values in Fig. 7(a).(c) A range colormap represents the HSV space.

Figure 8 (
Figure8shows the videos of the raw measurements and range maps associated with the compressive holography reconstructions.The raw measurements capture the movement of three water cyclopses in a 3mm depth container.In Fig.8(a), the frame rate for the 2D hologram recording was set to 15 fps whose frames are played at the rate of 3 frames per second for readers' visual convenience.Figure 8(b) shows a video of the range maps of the compressive holography reconstructions whose sizes are 1024 1024 35 .This reconstructed video has 109frames and the data processing was performed parallel on a computer cluster of the Scalable Computing Support Center in Duke University.The axial resolution for these reconstructions is set to 90 z