## Abstract

The phase retrieval algorithm has been used in this paper for whole reconstruction of the optical wave fields. The quantitative information of the phase distribution as well as the intensity distribution of the reconstruction field at different locations along the propagation direction has been achieved from double or multi in-line holograms. Numerical reconstructions of the wave fields from experimentally recorded in-line holograms are presented. This technique can be potentially applied for aberrated wave front analyzing and 3D imaging.

©2003 Optical Society of America

## 1. Introduction

With the development of the CCD and the computer techniques, the digital holography (DH) technique has been widely used in different fields of optics such as microscopy [1], deformation analysis [2], object contouring [3,4], and particles sizing [5]. In the DH, the hologram is recorded by a CCD and the reconstruction is carried out in a computer. Compared to the conventional holographic technique, the numerical reconstruction allows accessing not only the intensity distribution but also the phase distribution of the reconstructed field, which is impossible by using traditional optical reconstruction. Therefore, it can be used to compensate for aberrations [6] and for correcting image reconstruction in the presence of anamorphism [7], furthermore, the DH can also be used in whole field reconstruction [8], i.e. for determination of the intensity and of the phase distribution of the wave field at any arbitrary plane located between the object and the recording planes. Quantitative determination of the complex amplitude of the field propagating away from the object allows three-dimensional imaging, investigation of propagation of light in the phase-distorting media, and other applications.

Grilli [8] *et al.* studied the principle of the method for numerical reconstruction of the complex amplitude and have shown that this technique can be used for simultaneously determining the intensity and phase distributions at any location along the propagation. However, they just gave the intensity reconstructions from digitized experimental holograms, the zero order diffractions and the ghost images have not been removed. The emergence of the ghost image will blur the real object and affects the following analysis. In this paper, we will combine the in-line holography with the phase retrieval algorithm to reconstruct the whole optical field. Although the phase retrieval algorithm has already been used to reconstruct the in-line hologram by Liu *et al.* [9], the object is assumed to be a pure real or pure imaginary one, therefore, this approach is limited for wider applications. Here we reconstruct both the intensity and the phase distributions of an object from double or multi holograms, present the numerical reconstruction of the whole wave field from digitally recorded in-line holograms.

## 2. Description of the in-line holography

Holography is a method for reconstruction of optical wave fields. In the DH, the hologram is recorded onto a high resolution CCD, the numerical method is used to reconstruct the original object digitally. Although the reconstruction of off-axis holograms can separate the real image and the ghost image effectively, the size of the object is limited in order to ensure the separation between the real and twin images. Phase-shifting method [10] can derive the phase and the amplitude information of the object wave from four holograms, however, the precisely controlled subwavelength translations of the object or of the reference mirror are required. The in-line holography is a simple way for high-resolution imaging. Due to the fact that no optical elements such as a lens or a beam splitter are needed, it is suitable for short wavelengths imaging.

The typical setup for recording in-line holograms is shown in Fig. 1. The input plane wave serves not only as illumining light but also the reference light. The complex wave field in the hologram plane can be described by the Rayleigh-Sommerfield diffraction formula as

where *h*(*x*
_{1}, *y*
_{1}, *d*) is the impulse response of the free space propagation, d is the distance between the object and the hologram planes, “⊗” expresses the convolution operation. *U*
_{1} and *U*
_{2} are complex wave field distribution in the object and hologram planes, (*x*
_{1}, *y*
_{1}) and (*x*
_{2}, *y*
_{2}) are coordinates in the object and hologram planes, respectively. *U*
_{1}(*x*
_{1}, *y*
_{1})=1-*a*(*x*
_{1}, *y*
_{1}) and *a*(*x*
_{1}, *y*
_{1}) is the opacity function of the object. In this paper, we use Fresnel approximation of the Rayleigh-Sommerfield formula, i. e. the impulse response function is

Due to the CCD is only sensitive to the intensity of the field in the hologram plane, the distribution of the recorded hologram will be

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}={\mid 1-a({x}_{1},{y}_{1})\otimes h({x}_{1},{y}_{1},d)\mid}^{2}$$

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}=1-{a}^{*}({x}_{1},{y}_{1})\otimes {h}^{*}({x}_{1},{y}_{1},d)-a({x}_{1},{y}_{1})\otimes h({x}_{1},{y}_{1},d)$$

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}+{\mid a({x}_{1},{y}_{1})\otimes h({x}_{1},{y}_{1},d)\mid}^{2}.$$

Then the reconstructed field *U _{r}*(

*x*,

_{r}*y*,

_{r}*d*′) in any plane between the object and hologram planes can be obtained by using inverse Fresnel transformation as

where *d*′ is the distance backward measured from the hologram plane to the image plane. When *d*′=*d*, the reconstruction obtained will be

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}+{\mid a({x}_{1},{y}_{1})\otimes h({x}_{1},{y}_{1},d)\mid}^{2}\otimes h({x}_{1},{y}_{1},-d),$$

the first term is a uniform component, the second and third terms lead to the real and ghost images of the object, respectively, and the last term lead to a nonlinear component. In order to keep the size of the reconstructed image constant [11], the wave field *U _{r}*(

*x*,

_{r}*y*

_{r},

*d*′) can be alternatively computed by

where *F*{ } expresses the Fourier transform, *H*(ξ, η, *d*′) is the Fourier transform of the impulse response function *h*(*x*
_{2}, *y*
_{2}, -*d*′)

With this method, the size of the reconstructed image does not change when modifying the reconstruction distance. It will be easier to compare the reconstructed images in different planes with different distances *d*′. The fast Fourier transform can also be used to reduce the computing time. Furthermore, the exact solution of the diffraction integral can be archived as far as the sampling (Nyquist) theorem is not violated.

Experiments are carried out by using the setup shown in Fig. 1. The wavelength of the collimating light is 532nm. The used CCD has 1300×1030 pixels with pixel size of 6.7µ*m*×6.7µ*m*. Only 1024×1024 pixels are used to reconstruct images. The pure phase sample used in experiment is the logo of the Institut für Technische Optik shown in Fig. 2. The fabricated logo has a step height of ~330nm at the location of underlying glass substrate, the corresponding phase difference is ~1.8rad. The widths of lines are 40~60µ*m* depending on the different parts in the logo. Parts of recorded holograms with distances *d*=7.40*cm*, 8.15*cm*, 8.90*cm* are shown in Fig. 3(a), (b), and (c), respectively. Using the hologram shown in Fig. 3(c) and Eq. (5), the wave fields from *d*′=8.90*cm* to *d*′=2.00*cm* with step of Δ*d*′=1*mm* are reconstructed. Fig. 4 presents the clip videos of the intensity (a) and phase (b) distributions of these reconstructed filed, the phase distribution is wrapped in the range of [-π, π]. Due to the in-line geometry of the setup and lost the phase information in the hologram plane, the real and the ghost images are superposed. When *d*′=8.90*cm*, the original pure phase object should be reconstructed. However, the intensity distribution of the reconstructed field shows clear diffractive patterns of the object instead of the uniform background. The phase distribution is also distorted by a ghost image. In order to analyze the wave fields quantitatively, not only the intensity but also the phase distribution of the wave fields in the hologram plane should be provided. Many methods have been proposed to eliminate the ghost image. Superposition can be avoided by using off-axis reconstruction of the in-line holograms [12], however, this method is limited to small discrete objects due to losing low frequency components. The linear filtering method [13] requires the object should be a purely absorptive one. Nugent [14] proposed a method which can reconstruct complex object, however, the object should obey a certain size constraint. In the next Section, we will use phase retrieval algorithm to recover the necessary phase distribution from double or multi in-line holograms.

## 3. Whole field reconstruction from double holograms

Due to the fact of the lost phase information in the hologram plane and the superposition of the reference and illuminating lights, the twin image appears in the reconstructed field which will blur the real image. If we could recover the phase distribution in the hologram plane, we could reconstruct the whole wave field exactly. Many approaches have been proposed to recover phase distribution from intensity measurements [15–19]. Tiller *et al.* [19] developed a deterministic algorithm based on the transport of intensity equation, however, this method requires the distance between two closed spaced planes must be sufficiently small that the intensity varies approximately linearly with the distance between the planes. This distance must be large enough to compensate for spurious signals in the presence of noise; therefore, it will be difficult to select this distance without a priori information. In this paper, we will employ Gerchberg-Saxton-Fienup algorithm [16] to recover the phase information from two recorded holograms and then reconstruct the whole optical field. The distance between two recorded planes can be selected freely, which does not depend on the sample structure, noise, and other experimental parameters.

We use the holograms *H*
_{0} and *H*
_{1} to recover the phase distribution in the *H*
_{1} plane. After 500 iterations, the corresponding sum squared error [18] (SSE) of the predicted intensity relative to the known image intensity in the *H*
_{1} plane is SSE=0.2956. Using the recovered phase and known intensity distribution in the *H*
_{1} plane, the wave field distributions from *d*′=8.90*cm* to *d*′=2.00*cm* with the step of Δ*d*′=1*mm* have also been reconstructed. Fig. 5 (a) shows the intensity distribution and Fig. 5(b) describes the phase change along the wave propagation. Comparing with Fig. 4(b), Fig. 5(b) presents a clearer phase distribution of the original object when *d*′=8.90*cm*. The variations of the intensity and phase of the wave field with the propagation distance can be clearly seen. This character can be used to analyze the propagation of light in the phase-distorting media. However, the intensity distribution at *d*′=8.90*cm* is not an ideal uniform plane wave, some characters of original phase distribution can be found. This is because of the low resolution of the used CCD, where the high frequency components have been filtered out. Furthermore, the noise in recorded holograms and errors in distances measurement will also bring mistakes in the reconstruction. It should be pointed out that although this algorithm works in practice, it is sensitive to the initial value of phase, sometimes it may fail to find the true phase distribution in the hologram plane [17]. In order to overcome this drawback, we will try to reconstruct the whole wave field using multi in-line holograms.

## 4. Whole field reconstruction from multi holograms

Unlike the GS algorithm used between an assumed pure phase or pure amplitude object plane and hologram plane, the GS algorithm used between the two hologram planes plane seems to be under constrained and susceptible to uniqueness concerns. This problem may be solved by using multi holograms. Due to more information (multi holograms) are used, the recovered phase distribution should satisfy more constraints, and therefore, this extended algorithm is expected to provide a deterministic phase solution. Here, we use three holograms as a demonstration for whole optical field reconstruction, and the problem of uniqueness will be the subject of further work. As shown in Fig. 6, the iteration process begins from the central (*H*
_{0}) plane, uses the amplitude measured from *H*
_{0} and arbitrarily initial phase to construct a wave function. This wave function propagates to the plane *H*
_{1}, then replaces the amplitude with the amplitude measured from *H*
_{1} and keeps the phase unchanged, propagates back to the plane *H*
_{0}, replaces the amplitude with the amplitude with the amplitude measured from *H*
_{0} and leaving the phase unchanged, and continually propagates back to the plane *H*
_{-1}, the amplitude is again corrected to the value measured from *H*
_{-1}, and the corrected wave function was propagated to the central (*H*
_{0}) plane. After a complete cycle, the SSE will be checked. The iteration continues until the SSE is sufficient small or the SSE no longer decreases.

The reconstruction has been performed by using the three holograms shown in Fig. 3. After 75 iterations, the SSE arrives at its minimal value 15.96. Using the recovered phase distribution in the *H*
_{1} plane, the whole wave field from *d*′=8.90*cm* to *d*′=2.00*cm* with the step of Δ*d*′=1*mm* has been reconstructed. Fig. 7(a) and (b) show the intensity and phase distributions, respectively. Although we used three holograms, more information has been used, the SSE has not been reduced and the reconstruction quality has not been improved compares with double holograms. The noise in holograms and errors in distance measurements could be the main reasons. If the distances between the holograms are not correct, the corresponding Fresnel transforms cannot be calculated exactly, thus, these errors will be transferred into the recovered phase distribution after using iterative phase retrieval algorithm. The quality of the whole optical field reconstructed based on this phase distribution will be effected.

## 5. Summary

In this paper, we have investigated a new approach for a complete reconstruction of wave fields from double or multi holograms, the assumption of a purely absorptive or pure phase object can be avoided. The phase retrieval algorithms have been used to recover the phase distribution in the hologram plane and then the whole wave field in any plane between the object and hologram planes can be achieved. The intensity and phase distributions at different locations along the propagation direction can be determined simultaneously. The numerical reconstructed whole wave fields have been presented from experimentally recorded in-line holograms. This novel approach can be used in 3D imaging and for the analysis of aberrated wave front.

Acknowledgments

The authors thank Thomas Schoder for sample fabrication. One of the authors Yan Zhang acknowledges the financial support of the Alexander von Humboldt foundation.

## References and Links

**1. **U. Takaki and H. Ohzu, “Fast numerical reconstruction technique for high-resolutions hybrid holographic microscopy,” Appl. Opt. **38**2204–2211 (1999). [CrossRef]

**2. **S. Schedin, G. Pedrini, H. Tiziani, A. K. Aggarwal, and M. E. Gusev, “Highly sensitive pulsed digital holography for built-in defect analysis with laser excitation,” Appl. Opt. **40**100–117 (2001). [CrossRef]

**3. **G. Pedrini, P. Froning, H. Tiziani, and F. Santoyo, “Shape measurement of microscopic structures using digital holograms,” Opt. Commun. **164**257–268 (1999). [CrossRef]

**4. **C. Wagner, W. Osten, and S. Seebacher, “Direct shape measurement by digital wavefront reconstruction and multiwavelength contouring,” Opt. Eng. **39**79–85 (2000). [CrossRef]

**5. **S. Murata and N. Yasuda, “Potential of digital holography in particle measurement,” Optics & Laser Technology **32**567–574 (2000). [CrossRef]

**6. **A. Stadelmaier and H. H. Massing, “Compensation of lens aberrations in digital holography,” Opt. Lett. **25**1630–1632 (2000). [CrossRef]

**7. **S. De Nicola, P. Ferraro, A. Finizio, and G. Pierattini, “Correct-image reconstruction in the presence of severe anamorphism by means of digital holography,” Opt. Lett. **26**974–976 (2001). [CrossRef]

**8. **S. Grilli, P. Ferraro, S. D. Nicola, A. Finizio, and G. Pierattini, “Whole optical wavefields reconstruction by digital holography,” Opt. Express **9**294–302 (2001), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-9-6-294 [CrossRef] [PubMed]

**9. **G. Liu and P. D. Scott, “Phase retrieval and twin-image elimination for in-line Fresnel holograms,” J. Opt. Soc. Am. A **4**159–165 (1987). [CrossRef]

**10. **T. Zhang and I. Yamaguchi, “Three-dimensional microscopywith phase-shifting digital holography,” Opt. Lett. **23**1221–1223 (1998). [CrossRef]

**11. **J. W. Goodman, *Introduction to Fourier Optics*, (McGraw-Hill, San Francisco, Calif., 1968) Cap. 5.

**12. **S. Lai, B. Kemper, and G. v. Bally, “Off-axis reconstructions of in-line holograms for twin-image elimination,” Opt. Commun. **169**37–43 (1999). [CrossRef]

**13. **L. Onural and P. D. Scott, “Digital decoding of in-line holograms,” Opt. Eng. 26, 1124–1132 (1987).

**14. **K. A. Nugent, “Twin-image elimination in Gabor holography,” Opt. Commun. **78**293–299 (1990). [CrossRef]

**15. **R. W. Gerchberg and W. O. Saxton, “A practical algorithm for the determination of phase from image and diffraction plane pictures,” Optik (Stuttgart) **35**237–246 (1978).

**16. **J. R. Fienup, “Phase retrieval algorithms: a comparision,” Appl. Opt. **21**2758- 27 (1982). [CrossRef] [PubMed]

**17. **J. R. Fienup, “Phase-retrieval stagnation problems and solutions,” J. Opt. Soc. Am. A **3 1**897–1907 (1986).

**18. **G. Z. Yang, B. Z. Dong, B. Y. Gu, J. Zhuang, and O. K. Ersoy, “Gerchberg-Saxton and Yang-Gu algorithms for phase retrieval in a nonunitary transform system: a comparision,” Appl. Opt. **33**209–218. (1994). [CrossRef] [PubMed]

**19. **J. B. Tiller, A. Barty, D. Paganin, and K. A. Nugent, “The holographic twin image problem: a deterministic phase solution,” Opt. Commun. **183**7–14 (2000). [CrossRef]