Optimal spatial bandwidth capacity in multiplexed off-axis holography for rapid quantitative phase reconstruction and visualization

We present new methods for optimizing the spatial bandwidth capacity in off-axis holography using spatial multiplexing. We use optimal spatial multiplexing of off-axis holograms to fill the entire spatial frequency domain, including the space previously occupied by the intensity of the sample. Our approach enables spatial digital compression of eight offaxis holograms into a single real-valued multiplexed hologram, having the same number of pixels as each of the input holograms, but still allowing their full reconstruction without resolution or magnification loss in the reconstructed complex wave fronts. This new method allows 33% improvement in usage of the spatial bandwidth capacity compared to the best available off-axis holography real-value multiplexing method. Since the output multiplexed hologram contains only real values, it can be used for rapid display of eight wave front reconstructions at once, which is useful for real-time visualization, when the hologram display device is slower than the acquiring camera. We further generalize this technique to digital multiplexing of 16 real-valued holograms into a single complex-valued hologram by simple arithmetic operations in the hologram domain. Then, the extraction of the 16 wave fronts includes a single 2-D discrete Fourier transform to access the spatial frequency domain, allowing fast reconstruction, which is useful for real-time processing of off-axis holograms, with improved processing rate compared to current hologram processing algorithms. These new approaches allow the full reconstruction of all compressed data without loss of resolution or magnification, even though the samples are dense such that their frequency content employs the entire range. Both multiplexing architectures are then demonstrated for experimentally-acquired off-axis holograms for quantitative phase imaging of biological cells. © 2017 Optical Society of America under the terms of the OSA Open Access Publishing Agreement OCIS codes: (090.1995) Digital holography; (090.2880) Holographic interferometry; (090.4220) Multiplex holography; (090.5694) Real-time holography. References and links 1. C. M. Vest, Holographic Interferometry (Wiley, 1979). 2. N. T. Shaked, “Quantitative phase microscopy of biological samples using a portable interferometer,” Opt. Lett. 37(11), 2016–2018 (2012). 3. B. A. E. Saleh and M. C. Teich, Fundamentals of Photonics (Wiley, 2007). 4. P. Girshovitz and N. T. Shaked, “Compact and portable low-coherence interferometer with off-axis geometry for quantitative phase microscopy and nanoscopy,” Opt. Express 21(5), 5701–5714 (2013). 5. S. Karepov, N. T. Shaked, and T. Ellenbogen, “Off-axis interferometer with adjustable fringe contrast based on polarization encoding,” Opt. Lett. 40(10), 2273–2276 (2015). 6. A. W. Lohmann, “Reconstruction of vectorial wavefronts,” Appl. Opt. 4(12), 1667–1668 (1965). 7. Y. Ohtsuka and K. Oka, “Contour mapping of the spatiotemporal state of polarization of light,” Appl. Opt. 33(13), 2633–2636 (1994). 8. T. Colomb, F. Dürr, E. Cuche, P. Marquet, H. G. Limberger, R.-P. Salathé, and C. Depeursinge, “Polarization microscopy by use of digital holography: application to optical-fiber birefringence measurements,” Appl. Opt. 44(21), 4461–4469 (2005). Vol. 25, No. 26 | 25 Dec 2017 | OPTICS EXPRESS 33400 #307974 Journal © 2017 https://doi.org/10.1364/OE.25.033400 Received 26 Sep 2017; revised 9 Nov 2017; accepted 21 Nov 2017; published 21 Dec 2017 9. N. A. Turko and N. T. Shaked, “Simultaneous two-wavelength phase unwrapping using an external module for multiplexing off-axis holography,” Opt. Lett. 42(1), 73–76 (2017). 10. X. Wang, H. Zhai, and G. Mu, “Pulsed digital holography system recording ultrafast process of the femtosecond order,” Opt. Lett. 31(11), 1636–1638 (2006). 11. J. Kühn, T. Colomb, F. Montfort, F. Charrière, Y. Emery, E. Cuche, P. Marquet, and C. Depeursinge, “Real-time dual-wavelength digital holographic microscopy with a single hologram acquisition,” Opt. Express 15(12), 7231–7242 (2007). 12. M. Paturzo, P. Memmolo, A. Tulino, A. Finizio, and P. Ferraro, “Investigation of angular multiplexing and demultiplexing of digital holograms recorded in microscope configuration,” Opt. Express 17(11), 8709–8718 (2009). 13. Y. Wu, Y. Yang, H. Zhai, Z. Ma, L. Deng, and Q. Ge, “Single-exposure approach for expanding the sampled area of a dynamic process by digital holography with combined multiplexing,” J. Opt. 15(8), 085402 (2013). 14. P. Girshovitz, I. Frenklach, and N. T. Shaked, “Broadband quantitative phase microscopy with extended field of view using off-axis interferometric multiplexing,” J. Biomed. Opt. 20(11), 111217 (2015). 15. I. Frenklach, P. Girshovitz, and N. T. Shaked, “Off-axis interferometric phase microscopy with tripled imaging area,” Opt. Lett. 39(6), 1525–1528 (2014). 16. M. Rubin, G. Dardikman, S. K. Mirsky, N. A. Turko, and N. T. Shaked, “Six-pack off-axis holography,” Accepted to Opt. Lett. (2017). 17. P. Girshovitz and N. T. Shaked, “Real-time quantitative phase reconstruction in off-axis digital holography using multiplexing,” Opt. Lett. 39(8), 2262–2265 (2014). 18. B. Sha, X. Liu, X. L. Ge, and C. S. Guo, “Fast reconstruction of off-axis digital holograms based on digital spatial multiplexing,” Opt. Express 22(19), 23066–23072 (2014). 19. P. Girshovitz and N. T. Shaked, “Fast phase processing in off-axis holography using multiplexing with complex encoding and live-cell fluctuation map calculation in real-time,” Opt. Express 23(7), 8773–8787 (2015). 20. B. Sha, Y. J. Lu, Y. Y. Xie, Q. Y. Yue, and C. S. Guo, “Fast reconstruction of multiple off-axis holograms based on a combination of complex encoding and digital spatial multiplexing,” Chin. Opt. Lett. 14(6), 060902 (2016). 21. Z. Zhong, H. Bai, M. Shan, Y. Zhang, and L. Guo, “Fast phase retrieval in slightly off-axis holography,” Opt. Lasers Eng. 97, 9–18 (2017). 22. A. V. Zea, J. F. Barrera, and R. Torroba, “Cross-talk free selective reconstruction of individual objects from multiplexed optical field data,” Opt. Lasers Eng. 100, 90–97 (2018). 23. M. Paturzo, P. Memmolo, L. Miccio, A. Finizio, P. Ferraro, A. Tulino, and B. Javidi, “Numerical multiplexing and demultiplexing of digital holographic information for remote reconstruction in amplitude and phase,” Opt. Lett. 33(22), 2629–2631 (2008). 24. C. Rosales-Guzmán, N. Bhebhe, N. Mahonisi, and A. Forbes, “Multiplexing 200 spatial modes with a single hologram,” J. Opt. 19(11), 113501 (2017). 25. T. Tahara, Y. Awatsuji, K. Nishio, S. Ura, O. Matoba, and T. Kubota, “Space-bandwidth capacity-enhanced digital holography,” Appl. Phys. Express 6(2), 022502 (2013). 26. D. Roitshtain, N. A. Turko, B. Javidi, and N. T. Shaked, “Flipping interferometry and its application for quantitative phase microscopy in a micro-channel,” Opt. Lett. 41(10), 2354–2357 (2016). 27. M. A. Herráez, D. R. Burton, M. J. Lalor, and M. A. Gdeisat, “Fast two-dimensional phase-unwrapping algorithm based on sorting by reliability following a noncontinuous path,” Appl. Opt. 41(35), 7437–7444 (2002). 28. T. Tahara, R. Mori, Y. Arai, and Y. Takaki, “Four-step phase-shifting digital holography simultaneously sensing dual-wavelength information using a monochromatic image sensor,” J. Opt. 17, 125707 (2015). 29. N. T. Shaked, Y. Zhu, M. T. Rinehart, and A. Wax, “Two-step-only phase-shifting interferometry with optimized detector bandwidth for microscopy of live cells,” Opt. Express 17(18), 15585–15591 (2009).


Introduction
Digital holography can capture a sample complex wave front (amplitude and phase difference) by acquiring a real-value interference pattern between a sample beam, interacting with the sample, and a reference beam [1,2].In off-axis holography [3][4][5], one of the interfering beams is slightly tilted relative to the other beam, creating a linear shift that allows separation of the field intensity from the two complex-conjugate wave front terms in the spatial-frequency domain (SFD), and thus reconstruction of the complex wave front from a single hologram.This encoding leaves free space in the SFD, into which additional information can be compressed.The compression can be done by optical multiplexing [6][7][8][9][10][11][12][13][14][15][16], by projecting sample and reference beam pairs on the digital camera, each of which creating an off-axis hologram with a different interference fringe direction that positions one wave front in the spatial frequency domain without overlapping other terms.Alternatively, the multiplexing can be performed digitally [16][17][18][19][20][21][22][23], by positioning the wave front matrix in the SFD matrix without overlap to generate the multiplexed hologram matrix.Multiplexing is useful for various applications, with focus on acquisition of fast dynamics using optical multiplexing [6][7][8][9][10][11][12][13][14][15][16], as well as for speeding up hologram reconstruction [17][18][19][20][21] and for compression and storage of holographic data [22,23] using digital multiplexing.Digital multiplexing of hundreds of computer generated holograms can also be performed, e.g., for beam shaping [24]; however, in this case, the input holograms are not acquired optically first, but rather generated digitally.
For optically recorded holograms, the optical setup for the off-axis hologram acquisition dictates the maximal frequency that can be seen on the main axes of the SFD, and thus the number of holograms that can be multiplexed without overlap; a wider bandwidth can be obtained by using a digital camera with higher resolution or by increasing the magnification of the optical system; however, this comes at the cost of wasting camera pixels and decreasing the imaged field of view or, alternatively, the frame rate of the camera.Thus, in off-axis holography in general we require that the bandwidth usage would be minimal.
Maintaining a real-valued multiplexed hologram is crucial for applications including either recording the multiplexed hologram with a camera (in optical hologram multiplexing), or sending the compressed hologram to a display mechanism such as spatial light modulator.Various architectures have been suggested for multiplexing while maintaining both a realvalued multiplexed hologram and the original spatial-frequency content.These include multiplexing two holograms by positioning two complex wave fronts in two orthogonal directions [6][7][8][9]11,14], which can be generalized to multiplexing three [10,15], four [13], or even five holograms [12].We recently presented the six-pack holography (6PH) multiplexing architecture, which positions six complex wave fronts in the SFD without overlap [16]; 6PH represents the optimal bandwidth consumption when the intensity of the sample needs to be acquired together with its wave front in optical hologram multiplexing.
Other than optical recording and displaying applications that mandate real values in the final hologram, spatial digital multiplexing architectures of off-axis holograms were previously suggested for speeding up hologram reconstruction, by using simple arithmetic operations in the hologram domain to compress multiple wave fronts into a single hologram, from which the Fourier transform of all wave fronts can be retrieved with a single 2-D DFT [17][18][19][20].In this case, the multiplexed hologram can contain complex values, where the most recent architecture includes digital multiplexing of eight wave fronts into a single complex hologram [20].
In this paper, we present new methods for optimized spatial bandwidth consumption in off-axis holography, taking advantage of the space that was previously allocated to the sample intensity (zero order or DC terms).We first present the eight-pack holography (8PH) principle, yielding a multiplexed off-axis hologram containing only real values.This new approach optimizes the bandwidth usage in off-axis holography by digitally multiplexing eight complex wave fronts without overlap in the SFD, using the entire bandwidth for wave front encoding.The resulting hologram is still a real-valued matrix as long as spatialfrequency Hermitian symmetry is retained, making it useful for recording and visualization purposes.Specifically, 8PH is useful for displaying fast phenomena in real time, by effectively increasing the maximal display rate of the display device by a factor of eight.We quantitatively compare the SFD usage in the presented approach with all previous ones and present a fast algorithm allowing the digital multiplexing to be applied at a rate of over 50 multiplexed holograms per second.
Next, we generalize the concept of hologram multiplexing to allow multiplexed holograms containing complex values.We present four methods for compressing up to 16 DC-free holograms into a single complex-valued multiplexed hologram.We show that these holograms can be obtained experimentally by subtraction of the two outputs of the interferometer.This multiplexing procedure is useful for rapid phase reconstruction, as the hologram multiplexing can be achieved by simple arithmetic operations in the hologram domain.Then, a single Fourier transform can be used to access the SFD for isolating the 16 compressed sample wave fronts from the multiplexed hologram.

Theory
An off-axis interferogram recorded by a camera can be mathematically expressed as follows: where j denotes the imaginary unit, s E and r E are the sample and reference complex wave fronts, respectively, is the phase difference between them, where λ is the illumination wavelength and OPD is the optical path difference, and o ϕ is the off-axis phase, induced by titling one of the beams in respect to the other beam.The first two terms in Eq. ( 1) represent the sample and reference intensities, also referred to as zero orders, DC terms, or auto-correlation terms, whereas each of the last two terms contain the complex wave front of the sample, also referred to as the complex conjugate cross-correlation (CC) terms.The off-axis phase is given by: [ ] where x θ is the angle between the sample and reference waves in relation to the x axis, and y θ is the angle between the sample and reference waves in relation to the y axis, assuming light propagates in the z direction.This linearly space dependent phase is translated to a linear shift of the CC terms in the SFD to opposite sides of the axis, as shown in Fig. 1(a), respective to the tilt angles between the sample and the reference beams.This holographic encoding leaves empty space in the SFD, into which additional CC terms can be inserted by multiplexing several off-axis holograms, which creates a single multiplexed hologram containing the same number of pixels as each of the original holograms.These CC terms can contain additional information on the sample such as different temporal events, different imaging areas, different color images, etc.The multiplexing has to be designed such that all multiplexed wave fronts can be fully reconstructed from the multiplexed hologram.When the multiplexing is performed optically, we project onto the camera multiple holograms each with different x θ and y θ combinations, where all sample-reference pairs illuminate the camera at once, creating the multiplexed interference pattern on the camera.This, in turn, may create unwanted interference between nonmatching pairs; a phenomenon that can be avoided by either using different wavelengths for each sample and reference beam pair, using coherence gating or using different polarization states [15].This problem is not relevant for digital multiplexing, when the CC terms of the various wave fronts to be multiplexed are positioned in the SFD digitally.Assuming that unwanted cross-terms are not created and that the reference beam has constant intensity, the multiplexed hologram containing N wave fronts in the SFD can be mathematically formulated as follows: where ( , u v ) are the spatial-frequency coordinates and ( , ) u v δ is the Dirac delta function.As can be seen from Eq. (3), the DC terms all remain centered in the frequency domain, while each of the CC pairs are linearly shifted to opposite spatial frequencies respective of their tilt angles, due to the convolution with a shifted delta function.The relation between the tilt angle and the desired shift can be formulated as followed: where 0 u and 0 v express the shift on the u and v axes, respectively.According to the Nyquist-Shannon sampling theorem, the maximal frequency that can be seen in the SFD is half the sampling frequency.Thus, if the detector resolution in the x, y directions is , x y Δ Δ , respectively, the coinciding maximal angular frequencies are given by: max, max, In this analysis, we assume optical acquisition of dense samples (such as biological cells), such that the frequency content of the sample employs the entire range.Assuming that the sample maximal spatial frequency is c ω on both axes, each of the CC terms occupies a bandwidth of 2 c ω (as the reference beam is assumed to contain only the zero frequency), and the auto-correlation terms occupy a double spatial bandwidth of

(d). A comparison between
the various suggested architectures, in terms of sampling bandwidth per unit maximal frequency, number of holograms to bandwidth ratio per unit maximal frequency (representing the efficiency of the bandwidth usage), and relative space occupied by the CC terms in the frequency domain [16], is provided in Table 1.As can be seen, even though the latter method uses more bandwidth, it is still superior to the others by both the number of holograms to bandwidth ratio per unit maximal frequency and the relative space occupied by the CC terms.In this paper, we present the 8PH, a new multiplexing method for obtaining the optimal spatial bandwidth capacity in off-axis holography while using the space previously assigned to the DC terms.We show that the result is still a real-valued multiplexed off-axis hologram containing eight complex wave fronts without an overlap in the SFD, and thus can be reconstructed without magnification or resolution loss.
Optically or digitally removing the DC terms allows using an additional quarter of the bandwidth without loss of information.As long as conjugate (Hermitian) symmetry is retained in the SFD, the multiplexed image space remains real, and two more CC term pairs can be encoded, yielding a total of eight multiplexed holograms, as shown in Fig. 1(e), allowing a 1.33 × improvement factor relative to 6PH, as shown in Table 1.As shown in Fig. 1(e), to avoid negative values in the resulting multiplexed hologram, we still assign a narrow delta zero-order in the origin of the SFD, in the central pixel, which is not used for multiplexing anyway to enable Hermitian symmetry (see DC pixel in center of Fig. 2(d)).
As mentioned above, to obtain a real-value multiplexed hologram, all CC terms have to be located in the SFD in a Hermitian symmetry.Figure 2(a) shows a 1-D even array with Hermitian symmetry, as is received by applying 1-D discrete Fourier transform (DFT) to a real even array.Figure 2(b) shows the same array after applying fftshift, which shifts the zero-frequency component to the center of the spectrum.Figure 2(c) shows a 2-D even array with Hermitian symmetry, as obtained by applying 2-D discrete Fourier transform (DFT) to a real even array.Figure 2(d) shows the same array after applying fftshift.

Fast implementation of digital multiplexing for 8PH
The real-valued multiplexed hologram, containing eight complex wave fronts, can be sent to a display device such as a spatial light modulator (SLM), allowing the simultaneous optical display of eight temporal events by illuminating it with the reference beam.This application can be useful for visualization purposes, if the display device allowing optical reconstruction has lower temporal resolution capabilities than the camera acquiring the holograms but both have the same number of pixels.For example, most amplitude SLMs are limited to an input rate of about 50 fps, but by displaying eight wave fronts at each frame an effective display rate of 400 fps can be achieved.The observer can then decide which of the optical reconstructions should be viewed.
This latter frame rate can be easily obtained by off-the-shelf digital cameras.If both recording and display have to be done in real time, and assuming the recording can also be done at 400 fps, the digital compression has to be done at a rate of at least 400 fps as well, which is 50 compressions of eight holograms per second, in order to not be the time limiting factor.
To demonstrate this, we present an algorithm that can implement the 8PH multiplexing digitally at a rate of over 50 compressions per second for 512 × 512-pixel input holograms using Matlab R2013b on a desktop computer equipped with Intel's core i7-4790 3.60 GHz 16.0 GB RAM 64 bit CPU.We assume square holograms with a bandwidth of 8 c ω , where the CC term shifts are 0 0.75 max u ω = and 0 0 v = (see Fig. 1(a)).In this algorithm, instead of applying the time-costly 2D DFT to each of the eight holograms separately in order to isolate the CC terms from the spatial frequency domain, the 2-D DFT only has to be applied once for each group of four holograms, thus only twice for the entire set of eight holograms compressed into a single multiplexed hologram.This is done by rotating two of the holograms in the group by 90°, and creating a temporary complex synthetic hologram encoding a normal and a rotated hologram in its real part, and additional normal and rotated holograms in its imaginary part, thus taking advantage of the fact the fft2 function used in Matlab assumes the input to be complex.After performing 2-D DFT on the complex multiplexed hologram, 4 different synthetic complex CC terms are established, two vertical and two horizontal.Once all 4 are cropped using a square window, each pair can be used to find the CC terms originating from the initial holograms at the same orientation in the following manner: CC − indicate the vertical/horizontal pair cropped from the Fourier space of the multiplexed complex hologram, the asterisk indicates complex conjugation, j indicates the imaginary unit, and flip() indicated flipping the image on both axes.Note that in order to get all phase reconstructions in their original orientation, the CC terms reconstructed from the vertical pair, originating from the interferograms rotated by 90° at the beginning of the process, need to be rotated back by 90°.For a broader explanation see algorithm D, steps 1-5 in [19].Once all eight CC terms are extracted, they can be placed in the left part of the spectrum (as seen in Fig. 1(e)), and then copied to the right side with conjugation and 180° rotation, to yield Hermitian symmetry, before returning to the image domain by applying an inverse 2D DFT.
Finally, to obtain a digital multiplexed hologram with integer values in the range [0 2 N -1], as needed for an N bit display device, the image plane needs to be subtracted its minimal value, divided by its new maximal value, multiplied by 2 N -1 and rounded.Other than the obvious quantization error due to rounding which is inevitable, this action does not affect the reconstruction process; as mentioned above, the addition of pure DC (minimal value subtraction) simply changes the value of the single DC pixel in the 2-D Fourier space (see the origin of the SFD in Fig. 2(d)), which was left blank anyway to allow Hermitian symmetry, and the multiplication by a constant only influences the absolute value of the wave front, but not its phase, which is usually the relevant quantity for optically transparent biological cells when imaged in vitro.

Experimental results
To demonstrate the proposed 8PH technique, we have acquired off-axis holograms of cancer cells in flow and digitally multiplexed each eight holograms in the sequence into a single dynamic multiplexed hologram containing the entire data in one eighth of the temporal resolution of the original data with no loss of sample dynamics, allowing displaying the reconstruction at eight times the maximal framerate possible by the SLM.
Our experimental setup for capturing the cells during flow is depicted in Fig. 3(a).The setup included a 7.9 mW, 633 nm Helium-Neon laser that illuminates an inverted microscope.The laser beam passes through the sample and is magnified by a 50× , 0.55 NA microscope objective (M Plan Apo, Mitutoyo), followed by a 150 tube lens, causing a total magnification of 37.5×.The beam then passes through a telecentric lens system, causing an additional magnification of 7.5× , meant to prevent undersampling given the large pixel size of the camera (20 μm ).Before the image plane of the system, where the digital high-speed camera is positioned (FASTCAM Mini AX200, Photron), we positioned our external flipping interferometric module, which is based on a modified Michelson interferometer [Fig.3(a)] [26].This interferometer flips half of the optical field of view on top of the other half using a retro-reflector.Assuming that half of the optical field is empty, this half can be considered as the reference beam.The two beams interfere on the digital camera at an angle, creating an off-axis image hologram on the digital camera.To meet the requirement of half empty optical field of view, we positioned the optical field of view on the side of the microfluidic channel (μ-slide VI0.1, IBIDI), so that only half of the optical field of view contains flowing cells and the other one contains the empty channel glass [26].We imaged breast cancer cells from a MDA-MB 468 cell line during flow.800 off-axis 512 × 512 holograms were recorded at a framerate of 400 fps during the cancer cell flow in the microfluidics channel (see Fig. 3(b) and Visualization 1), and then compressed at a rate of 400 fps (50 compressions of 8 holograms per second).Figure 3(c) shows the results of applying digital eight-pack holographic compression to the experimental data using the complex multiplexing algorithm.We then reconstruct the OPD profiles from the multiplexed hologram, where the reconstructions are shown in Fig. 3(d) and Visualization 2, in comparison to the OPD profiles reconstructed without compression, shown in Fig. 3(e) and Visualization 3. The OPD extraction process from the CC terms included applying an inverse 2D-DFT, extracting the phase from the complex field, subtracting a sample-free hologram, applying 2-D phase unwrapping and dividing by the wavenumber [17,27].The MSE of the OPD profiles reconstructed from the compressed, quantized 8PH holograms, when using the fast complex multiplexing algorithm for compression, relative to the OPD profiles reconstructed with no compression/quantization was 31.99 2 nm , resulting from both the quantization and numerical errors in the reconstruction process, mostly in the complex multiplexing algorithm.In comparison, when performing the compression and quantization without the fast complex multiplexing algorithm (but rather by performing eight 2D DFTs), the MSE received is 3.67 2 nm .Note that the problem of quality loss due to quantization stemming from grayscale dynamic range sharing in the multiplexed hologram is expected to be more severe when the amplitude of the sample becomes non-negligible [28], requiring cameras with larger bit depth.This issue was less critical in our case, which deals with almost transparent biological cells.

Generalization to complex-valued multiplexed hologram
For applications not including optical multiplexing or digital multiplexing followed by optical reconstruction for display purposes, the concept of multiplexing can be further generalized to allow multiplexed holograms with complex values.This is beneficial for speeding up phase reconstruction from a sequence of off-axis holograms by compressing multiple holograms from the sequence into a single multiplexed hologram, from which the CC terms of all wave fronts can be retrieved with a single 2-D DFT.This procedure assumes that the hologram multiplexing is done in the image domain by simple arithmetic operations which are less computationally intense than a 2-D DFT.

Algorithms for fast phase reconstruction
The new understanding of the so-far unused SFD space gives rise to four new algorithms, allowing digitally compressing up to 16 wave fronts in a single multiplexed hologram.This digital processing relies on the fact that the input holograms do not have DC terms, and is carried out in the hologram domain by arithmetic operations that are simpler than 2-D DFTs, thus sparing computation time for hologram reconstruction.The former requirement of holograms without DC terms can be achieved even for dynamic samples by using both exits of a Mach-Zehnder interferometer, where the first exit projects onto the first camera .This creates two phase-shifted off-axis holograms encoding the same information, with identical DC terms, but with CC terms of opposite signs.Thus, subtracting one hologram from the other one yields a DC-free hologram, while the CC terms remain [21,29].For all suggested algorithms, the input holograms are assumed to be with a bandwidth of 8 c ω , where the CC term shifts are 0 0.75 max u ω = and 0 0 v = (see Fig. 1(a)).First, we further develop the method suggested by Sha et al. [18] (SM algorithm), which multiplies each hologram with a tilted plane wave, creating a linear shift in the vertical direction in the respective SFD (see the SFD of the resulting hologram in Fig. 4(a)), as following: In the new algorithm suggested here, termed Algorithm G (as a continuation to Algorithms A-F for fast holographic processing previously suggested by us [17,19]), we suggest taking each group of 12 holograms and dividing it into three groups of four holograms each.At step 1, each group undergoes the process suggested by Sha et al. [Eq.( 7)] to create a total of three multiplexed holograms.At step 2, the second and third multiplexed holograms are further multiplied by other tilted plane waves to create the SFDs shown in Figs.4(b) and 4(c), respectively, using the cyclic property of DFT, and then summed up.This can be mathematically formulated as: Note that for efficiency, steps 1 and 2 can be performed as one step, by creating diagonally titled plane waves introducing the necessary shift in a single multiplication.Once the multiplexed hologram is created, we can take the 2-D DFT to obtain the SFD shown in Fig. 4(d).In spite of the overlap of the CC terms created by the summation, we can use the fact the CC terms of each wave front are complex conjugate pairs to retrieve all CC terms with no information loss.The CC terms of holograms 1-4, suffering no overlap, can be simply cropped from the furthest right quarter.The CC terms of holograms 5-8, which are overlapped by the CC terms of holograms 1-4, can be obtained by subtracting the both-axes flipped conjugates of the CC terms of holograms 1-4 from the CC terms in the far left quarter of the same row.The CC terms of holograms 9-12, again with no overlap, can be cropped from the second right quarter.From that point on, the reconstruction can go on the standard manner, by applying an inverse 2-D DFT and extracting the quantitative phase profile.Note that by introducing a third linear phase shift in the image domain, placing the CC terms in the two right columns, and summing with the multiplexed holograms respective of the SFDs shown in Figs.4(a)-4(c), we obtain a multiplexed hologram containing 16 holograms.However, in this case the original CC terms cannot be reconstructed from the overlapping CC terms, as even though the number of unknowns is identical to the number of equations, the resulting equations are non-linear, due to the flipping operation.
We now propose another fast off-axis hologram processing algorithm, termed as Algorithm H.This algorithm takes only three DC-free holograms as an input and skips step 1 of Algorithm G.Then, it uses Eq. ( 8) on the input holograms rather than on multiplexed ones, yielding the SFD illustrated in Fig. 4(e).The result can then be 4 times sampled in the row dimension, prior to applying 1-D DFT on the rows and restoring the CC terms in a similar manner to the one explained for Algorithm G; the CC term of hologram 1 is cropped from the furthest right quarter, the CC term of hologram 2 is acquired by subtracting the x-axis flipped conjugate of the CC term of hologram 1 from the CC term in the far left quarter, and the CC term of hologram 3 is cropped from the second right quarter.Then, we apply a 1-D inverse DFT on the rows and extract the quantitative phase in the standard manner [19].
Next, we further develop a more recent method suggested by Sha et al. [20] (CSM algorithm), which creates a complex hologram from each pair of input real-valued holograms (placing one in the real part and the other in the imaginary part) prior to multiplying each with a tilted plane wave according to Eq. ( 7), allowing multiplexing of 8 holograms .The resultant SFD is illustrated in Fig. 4(f).Although mixing of CC pairs is allowed here, all 8 CC terms can be fully retrieved, because in this complex-hologram encoding the mixed CC terms are not complex conjugates, but rather contain different information.
In the new algorithm suggested here, termed Algorithm I, we suggest taking 16 holograms as an input, and performing the process suggested in Ref [20] on each group of 8.Then, the second multiplexed hologram is multiplied by exp(jπ/(x • Δx)), causing the CC terms to shift to the two inner quarters of the SFD, prior to summing it with the first multiplexed hologram.The coinciding SFD is illustrated in Fig. 4(g).Here too, for efficiency, steps 1 and 2 can be performed together.Once the multiplexed hologram is created by these simple operations in the hologram domain, we can perform a single 2-D DFT to obtain the SFD shown in Fig. 4(g).To retrieve each wave front, we can use the two CC terms containing information about it.For example, to retrieve the CC terms of holograms 1 and 2, we can simply crop the upper right and upper left quarters.Then, the original CC terms can be found using Eq. ( 6).From that point on, the reconstruction can go on the standard manner, by applying an inverse 2-D DFT and extracting the quantitative phase profile.This architecture defines the limit in terms of complex multiplexing capacity.
In the fourth and last algorithm suggested here, termed Algorithm J, only 4 holograms are taken as inputs.As before, a complex hologram is created from each pair, prior to multiplying the second complex hologram by exp(jπ/(x • Δx)), causing the CC terms to be shifted to the two inner quarters of the SFD, and summing it with the first complex hologram to yield the SFD illustrated in Fig. 4(h).This time, however, the result is 4 times sampled in the row dimension.Then, we can apply 1-D DFT on the rows and restore the CC terms in a similar manner to the one explained for Algorithm I, where the only difference is that the flipping operation in Eq. ( 6) is applied only on the x axis.Then, we apply 1-D inverse DFT on the rows and extract the quantitative phase in the standard manner [19].
A comparison of these algorithms relative to previous algorithms suggested in [17][18][19][20], both in terms of reconstruction time and quality, is provided in Table 2 for various input size holograms, using Matlab R2013b on a desktop computer equipped with Intel's core i7-4790 3.60 GHz 16.0 GB RAM 64 bit CPU.
Algorithm B [17] is a basic algorithm, and does not include multiplexing, but rather only a 2-D DFT for each hologram, followed by cropping the CC term and applying inverse 2-D DFT, followed by phase extraction and unwrapping.Algorithm C [17] multiplexes two holograms by summing one hologram with another one which is 90 ֯◌-rotated, before applying a single 2-D DFT to extract both CC terms.Algorithm D [19], briefly described in Section 2.2, multiplexes four holograms by placing two 90 ֯◌-rotated holograms in the real part of a synthetic hologram and additional two in its imaginary part, prior to applying a single 2-D DFT for all four holograms.The SM algorithm [18] also multiplexes four holograms into a single complex hologram, by multiplying each hologram with a tilted plane wave, as explained in Eq. ( 7), and performing 1-D DFT on the rows followed by 1-D DFT on only one quarter of the columns, creating an effective 1.25 DFT.Algorithm E [19] is similar to Algorithm B, only that it resamples the rows four times prior to applying 1-D DFT on the rows only.Algorithm F [19] is a combination of Algorithms D and E, storing one resampled hologram in the real part of a synthetic complex hologram, and a second hologram in its imaginary part.The CSM algorithm [20] is a combination of the SM algorithm and algorithm D, multiplexing eight holograms into a single complex hologram by creating four pairs of complex holograms and multiplying each pair with a tilted plane wave as explained in Eq. (7).All these previous algorithms are compared here to the four new ones, Algorithms G-J, where all algorithms were implemented in the most efficient way possible and on the same computer platform to allow a fair comparison.The timing for Algorithms G-J includes the time needed to subtract the two input holograms.Figure 7 presents the OPD profiles obtained by the various algorithms.As can be seen from Fig. 7, the reconstruction quality obtained using either of the new efficient algorithms is similar to the one obtained using the basic algorithm, further emphasizing that the multiplexing procedure is lossless.Note that in Figs. 6 (a) and 6(b), the raw hologram is shown, supplying amplitude contrast.For this reason, fine details, such as the sperm cell tail, can still be seen.On the other hand, in the phase profiles shown in Fig. 7, the sperm cell tail can hardly be seen due to much lower phase values.In any case, the reconstruction quality is similar for all algorithms.

Discussion and conclusions
We have presented new digital multiplexing architectures for optimally using the SFD in offaxis quantitative phase microscopy, in terms of bandwidth vs. compression potential, by using the space previously occupied by the intensity of the sample (DC), and suggested two main applications benefiting from this observation.
First, we demonstrated spatial digital compression of eight off-axis holograms (8PH) into a single real-valued multiplexed hologram, having the same number of pixels as each of the input holograms, but still allowing their full reconstruction without resolution or magnification loss in the reconstructed complex wave fronts.By comparing this architecture to previous architectures, we showed 33% improvement both in terms of relative space occupied by the CC terms in the frequency domain, and in terms of number of holograms to bandwidth ratio per unit maximal frequency.8PH is suitable for digital multiplexing, not requiring DC-free holograms, and is useful for displaying fast phenomena in real time, by effectively increasing the maximal display rate of the display device by a factor of eight (displaying eight dynamic phenomena in parallel).We provided a fast algorithm for applying the digital multiplexing for 8PH, allowing a rate of 50 compressions of eight off-axis holograms per second for 512 × 512-pixel input holograms, enabling real time display where the recording is done at 400 fps and the SLM is limited to 50 fps.This method was later tested on experimental data, showing no significant quality loss relative to standard reconstruction.In the general case, the new 8PH architecture is limited to digital multiplexing due to the space occupied by the DC terms; nevertheless, it can be applied optically by removing the DC terms, as can be done by subtracting the multiplexed off-axis holograms obtained in the two exits of a Mach-Zehnder interferometer, as explained in Section 3.1 and demonstrated in Section 3.2.In this case, although not demonstrated in this paper, this architecture can be useful for simultaneous optical recording of eight wavelength channels, eight fields of view of the sample, or even eight angular viewpoints of the sample.
Note that in our case we did not experience problems due to the grayscale dynamic range sharing in the multiplexed hologram due to the digital multiplexing procedure.In any case, for phase objects, such as live cells, sharing the dynamic range is not expected to affect the results significantly.In addition, due to the digital multiplexing, we did not experience unwanted cross-terms that might appear in optical multiplexing due to interference of nonmatching beam pairs.In optical multiplexing, this problem can be solved by coherence gating.
Second, we suggested multiplexing either 12 (Algorithm G), 3 resized (Algorithm H), 16 (Algorithm I) or 4 resized (Algorithm J) DC-free holograms into a single complex-valued multiplexed hologram by simple arithmetic operations in the image domain, from which the Fourier transform of the all wave fronts can be retrieved by a single 2-D discrete Fourier transform, allowing fast reconstruction, which is useful for real-time processing of off-axis holograms.Our Algorithms I and J present the maximal multiplexing limit, yielding a 100% compression improvement from the most optimal previous digital multiplexing architecture [20].We showed that these DC-free holograms can be obtained experimentally by subtraction of the two output holograms obtained in the two exits of the interferometer, and provided experimental results demonstrating it.We also compared the reconstruction rates achieved by the new suggested algorithms to previous multiplexing-based fast reconstruction algorithms, showing high variability in efficiency, tightly correlated to the input hologram size.In this analysis, Algorithm G was proved to be most efficient for small input holograms, the SM algorithm was proved to be most efficient for mid-range input holograms, and Algorithms J, H, E and F were proved to be the most efficient for larger input holograms.
We expect the new concept presented here to be useful both for real-time visualization of fast phenomena, and real-time phase reconstruction in off-axis holography.

Funding
Horizon 2020 European Research Council (ERC) Grant No. 678316; Tel Aviv University Center for Light-Matter Interaction.

Fig. 1 .
Fig. 1.Schematic illustrations of the holograms (left) and the coinciding spatial power spectra (right) for various off-axis hologram multiplexing architectures, including bandwidth calculations.(a) Conventional off-axis holography.(b) SPACE[25].(c) Diagonal multiplexing of two channels[25].(d) 6PH[16].(e) 8PH, proposed in the paper.In the spectra images, DC denotes the auto-correlation terms, the numbered circles around it denote the CC terms, where the coinciding complex conjugate CC terms are denoted by a number and an asterisk.Note that the relative size of the CC and DC terms change when the total bandwidth changes, since the number of pixels remains constant.

u
4 c ω .The intensity and power spectrum of an off-axis hologram with a bandwidth of 8 c ω (maximal frequency 4 c v = , are shown in Fig. 1(a).A few architectures have been previously suggested for multiplexing several holograms while maintaining both a real multiplexed image plane and the full frequency content of all samples without overlap.Tahara at el[25].proposed the SPACE architecture for the recording of a single hologram, where each of the CC terms is split on the two sides of the SFD, as shown in Fig.1(b), using the cyclic property of the DFT.This method uses a minimal sampling bandwidth of 5.12 c ω by choosing the CC shifts as 0 ω , the CC terms can be placed on the diagonal, enabling the compression of up to two holograms with the use of a slightly larger bandwidth of 6.24 c ω[25], as shown in Fig.1(c).Most recently, we have suggested multiplexing six holograms (6PH) with the use of a bandwidth of 8 c ω[16], as shown in Fig.1
− indicate the CC terms originating from the initial holograms,

Fig. 3 .
Fig. 3. Quantitative phase imaging of cancer cells during flow.(a) Experimental setup for generating the regular holograms as the input for the digital multiplexing process; MO: microscope objective, RR: retro reflector, TL: tube lens, BS: beam splitter (b) Dynamic offaxis hologram (Visualization 1) (c) Left: multiplexed power spectrum, right: multiplexed 8PH hologram.(d) Eight reconstructed OPD profiles from a single multiplexed hologram (Visualization 2).(e) Eight OPD profiles reconstructed without the 8PH compression (Visualization 3), for comparison.White scale bars represent 5 μm on the sample.Colorbar represents OPD values in nm.

Fig. 4 .
Fig. 4. A schematic illustration of the SFD of various complex multiplexed holograms, matching the different stages of Algorithms G, H, I and J for fast hologram processing presented in this paper.(a) Step 1 in Algorithm G. (b, c) Multiplexed holograms 2 and 3, respectively, in step 2 of Algorithm G.(d) Final multiplexed hologram consisting of 12 wave fronts in Algorithm G, with overlaps of CC terms, but still allowing full reconstruction.(e) Step 1 in Algorithm H, prior to resampling four times in the row dimension.(f) Step 1 in Algorithm I. (g) Step 2 in Algorithm I. (h) Step 1 in Algorithm J, prior to resampling four times in the row dimension.j indicates the imaginary unit.Note that although the coinciding SFDs are shown, all operations are performed in the hologram domain, without the need to transform the holograms to the SFD first.

Fig. 5 .
Fig.5.A schematic drawing of the modified Mach-Zehnder interferometer used for the experimental demonstrations for generating the regular DC-free holograms as the input for the digital multiplexing process.AOTF, acousto optical tunable filter; BS1, BS2, beam splitters; M1, M2, mirrors; RR1, RR2, retroreflectors; S, sample; MO1, MO2, microscope objectives; L1, L2, tube lenses.The two phase-shifted off-axis holograms from the two cameras are subtracted to obtain a DC-free off-axis hologram.

Fig. 6 .
Fig. 6.(a) Hologram from camera 1 (b) Hologram from camera 2 (c) Hologram obtained by subtracting the holograms from (a) and (b).(d-f) Coinciding SFD absolute values for (a-c).(g) Absolute value of multiplexed hologram in Algorithm G. (h) Coinciding SFD absolute value.(i) Absolute value of multiplexed resized hologram in Algorithm H. (j) Coinciding SFD absolute value.(k) Absolute value of multiplexed hologram in Algorithm I. (l) Coinciding SFD absolute value.(m) Absolute value of multiplexed resized hologram in Algorithm J. (n) Coinciding SFD absolute value.