## Abstract

A novel fast double-phase retrieval algorithm for lensless optical security systems based on the Fresnel domain is presented in this paper. Two phase-only masks are efficiently determined by using a modified Gerchberg-Saxton algorithm, in which two cascaded Fresnel transforms are replaced by one Fourier transform with compensations to reduce the consumed computations. Simulation results show that the proposed algorithm substantially speeds up the iterative process, while keeping the reconstructed image highly correlated with the original one.

©2009 Optical Society of America

## 1. Introduction

Optical security systems with digital transmission have gained increasing attention [1–24] due to the convenience in transmitting enormous data through the Internet or World Wide Web. There are two kinds of optical security techniques. First is the optical encryption, which is a process of scrambling perceptible information into a ciphered one by using a secret key. The ciphered information can only be perceived by those who have the correct key for decryption. Second is the optical authentication, a process of verifying the integrity of information. Most of the previous optical information security techniques (optical encryption [1–12, 19] and optical authentication [14–18, 20]) utilize the conventional 4-f correlator architecture based on the optical Fourier Transform (FT) [1–4, 14–17, 18, 19, 23]. Some modified schemes [15–17] were proposed to improve the 4-f correlator architecture, but they still need the use of lenses.

Situ and Zhang [20] proposed an algorithm which uses a novel iterative scheme to retrieve two phase-only masks (POMs) in the Fresnel domain. Differing from conventional lost-phase reconstruction, phase retrieval here is meant to encode a target image into two POMs. Figure 1 shows its optical setup, which is lensless and thus easy to be implemented. There are three planes defined: (1) the input plane where POM_{2} is located, (2) the transform plane where POM_{1} is located at a distance of *z*
_{2} from the input plane, and (3) the output plane which is at a distance of *z*
_{1} from the transform plane. The distance parameters, *z*
_{1} and *z*
_{2}, are determined according to the size of the POMs to satisfy the Fresnel approximation [21]. The target image can be obtained at the output plane when the system is directly illuminated with a plane wave of a specific wavelength. The phase retrieval in their system is performed digitally, while the target image reconstruction can be implemented optically or digitally. For optical implementation, the phase distributions can be produced in a form of POMs with micro-optics fabrication techniques [22].

Situ and Zhang’s iterative algorithm is based on a modified POCSA (Projection Onto Constraint Sets Algorithm) which adopts the Fresnel transform (FrT) [20, 28] instead of the conventional FT [29, 30]. The algorithm is summarized as follows.

Step 1: POM_{1} and POM_{2} are first generated with random phase functions, which will be updated in each of following iterations (Steps 2~6).

Step 2: The above-generated POM_{2} is Fresnel-transformed to the transform plane.

Step 3: The result of Step 2 is multiplied by the POM_{1} function.

Step 4: An approximated image is obtained in the output plane by performing FrT on the result of Step 3. If the amplitude of the approximated image satisfies a convergence criterion, then stop the process; otherwise, modify the approximated image according to a certain constraint.

Step 5: The modified amplitudes and phases at the output plane are inversely Fresnel-transformed back to the transform plane for updating the phase function of POM_{1}.

Step 6: The result of Step 5 is inversely Fresnel-transformed back to the input plane for updating the phase function of POM_{2}.

Referring to the above algorithm, it is obvious that each pass of iteration includes two FrTs and two Inverse FrTs (IFrTs). As shown in Fig. 1, the forward sub-pass performs two stages of FrTs: the first is from the input plane to the transform plane and the second is from the transform plane to the output plane. Similarly, the backward sub-pass also includes two stages of IFrTs. It is clear that Situ and Zhang’s iterative algorithm based on FrT requires the same number of transforms in an iteration as that in the conventional POCSA based on FT [29, 30]. If the FrT and IFrT operations are simplified or replaced by other more efficient transform, a faster phase retrieval can be achieved.

A fast phase-retrieval algorithm is proposed in this work, which is also advantageous of lensless optical setup, high security and digital implementation. Our proposed method modifies the traditional Gerchberg-Saxton algorithm (GSA) based on the FTs [26, 27] to solve the above-mentioned problem. It is more efficient since only a half number of operations is required with respect to POCSA. Based on FrT, POM_{1} and POM_{2} can be generated straightforwardly and the reconstructed image is highly correlated to the target one. However, as shown later in Section 3, the requirements on the precision of the wavelength and the POMs location are of an order of 10^{0} nm and 10^{2} µm, respectively, which is difficult to practice for today’s optical engineering. Thus the proposed algorithm is more suitable for digital implementation, which can be performed simply with a digital computer. This, however, gives another advantage: the retrieved phase data can be directly transmitted over a digital communication channel and the target (original) image can be faithfully recovered only in case that correct keys (the positions of POMs and the wavelength) are available for authorized receivers.

The rest of this paper is organized as follows. The proposed phase-retrieval algorithm is presented in Section 2. In Section 3, the algorithm is simulated and verified digitally, with impact analyses on the wavelength and POM position imprecision. Finally, the conclusions are drawn in Section 4.

## 2. Proposed phase retrieval algorithm

In our lensless security system, two phase masks are also required to reconstruct a target image. Let the target image be denoted as *g*(*x*
_{0}, *y*
_{0}) and the two phase distributions be denoted as *ψ*
_{2}(*x*
_{2}, *y*
_{2}) (POM_{2}) and *ψ*
_{1}(*x*
_{1}, *y*
_{1}) (POM_{1}). Figure 2 shows the systematic diagram of the proposed modified Gerchberg-Saxton algorithm (MGSA). First, the given target image is inversely Fourier-transformed to obtain its amplitude function *G*(*x*
_{1},*y*
_{1}) and phase function *ψ*
_{1}(*x*
_{1},*y*
_{1}). Then, the amplitude function *G*(*x*
_{1},*y*
_{1}) is inversely Fourier-transformed with an amplitude constraint to get another phase function *ψ _{G}*(

*x*

_{1},

*y*

_{1}) (see the details about MGSA below), which, together with

*ψ*(

_{t}*x*

_{1},

*y*

_{1}), are encoded into POM

_{1}and POM

_{2}, respectively. Since the FT is used during the phase retrieval process, but the FrT is adopted for target image reconstruction, some quadratic phase terms need to be manipulated to account for this mismatch, which will be explained in details below.

Generally, phase-retrieval adopting Fresnel diffraction can be expressed as

$$\genfrac{}{}{0.1ex}{}{\mathrm{exp}\left(\genfrac{}{}{0.1ex}{}{j2\pi z}{\lambda}\right)}{j\lambda z}\int \int f({\mu}_{0},{v}_{0})\mathrm{exp}\left\{\genfrac{}{}{0.1ex}{}{j\pi}{\lambda z}\left[{\left({\mu}_{1}-{\mu}_{0}\right)}^{2}+{\left({v}_{1}-{v}_{0}\right)}^{2}\right]\right\}d{\mu}_{0}d{v}_{0},$$

where *λ* denotes the wavelength of the plane wave and *z* denotes the distance between the input (*µ*
_{0}, *ν*
_{0}) and the output (*µ*
_{1}, *ν*
_{1}) plane. When a unit plane wave is incident on a pattern function *f*(*µ*
_{0}, *ν*
_{0}) in the (*µ*
_{0}, *ν*
_{0})domain, a Fresnel diffraction pattern *F*(*µ*
_{1}, *ν*
_{1}) can be obtained in the (*µ*
_{1}, *ν*
_{1})domain. Inversely Fourier- transforming (IFTing) the target image *g*(*x*
_{0}, *y*
_{0}) at a distance *z*
_{1} from the output plane obtains

$$=G({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{t}({x}_{1},{y}_{1})\right],$$

where exp[*jψ _{t}*(

*x*

_{1},

*y*

_{1})] is one of the components of POM

_{1}whose derivation will be shown below. In order to obtain a pure phase distribution to approximate

*G*(

*x*

_{1},

*y*

_{1}) in Eq. (2), an MGSA based on the Fourier domain is proposed, whose block diagram is shown in Fig. 2. The algorithm starts to perform IFT on

*G*(

*x*

_{1},

*y*

_{1}) and retains the resulting phase function

*ψ*′

_{2}(

*x*

_{2},

*y*

_{2}). Next, the phase function

*ψ*′

_{2}(

*x*

_{2},

*y*

_{2}) is constrained with unity amplitudes and then Fourier-transformed. From the FT result, the approximation

*G*̂(

*x*

_{1},

*y*

_{1}) and the phase function

*ψ*(

_{G}*x*

_{1},

*y*

_{1}) can be obtained. So far, only iterations between the input and transform planes is required. Again, perform IFT on

*G*(

*x*

_{1},

*y*

_{1}) with an updated phase function

*ψ*(

_{G}*x*

_{1},

*y*

_{1}) and continue the above iterative process. When the required correlation (similarity) between the exact

*G*(

*x*

_{1},

*y*

_{1}) and the approximation

*G*̂(

*x*

_{1},

*y*

_{1}) is achieved (e. g., a correlation coefficient

*ρ*=0.999), the resultant pure phase distribution

*ψ*′

_{2}(

*x*

_{2},

*y*

_{2}) can be determined. This iterative process can converge efficiently. As to the processing between the transform and the output planes, no iterations are needed since

*G*(

*x*

_{1},

*y*

_{1})can be obtained directly by one IFT as in Eq. (2). After the phase function

*ψ*′

_{2}(

*x*

_{2},

*y*

_{2}) is determined,

*G*̂(

*x*

_{1},

*y*

_{1}) can be obtained when

*ψ*′

_{2}(

*x*

_{2},

*y*

_{2}) is Fourier-transformed as in Eq. (3). Differing from the traditional GSA [26,27] which is generally applied to solve the problem of retrieving the lost phase in the Fourier domain, our proposed MGSA is suitable for generating a pure phase distribution

*ψ*′

_{2}(

*x*

_{2},

*y*

_{2}) from the amplitude function

*G*(

*x*

_{1},

*y*

_{1}) with only a few iterations, as shown in Fig. 2. The relationship is

$$=\hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{G}({x}_{1},{y}_{1})\right],$$

where exp[*jψ _{G}*(

*x*

_{1},

*y*

_{1})] is an accompanied phase term when we perform forward FT on exp[

*jψ*′

_{2}(

*x*

_{2},

*y*

_{2})]. The amplitude

*G*̂(

*x*

_{1},

*y*

_{1}) is then used to produce the approximated target image

*g*̂(

*x*

_{0},

*y*

_{0}), as in Eq. (2):

In the sequel, we first encode *ψ*′_{2}(*x*
_{2}, *y*
_{2}) into POM_{2}(*ψ*
_{2}(*x*
_{2},*y*
_{2})). Since the amplitude *G*̂(*x*
_{1},*y*
_{1}) is related to the phase *ψ*′_{2}(*x*
_{2},*y*
_{2}) via FT (see Eq. (3)) and needs to be reproduced from the FrT of POM_{2}(see Eq. (5) below), the quadratic and the propagation phase terms arising from the mismatch between these two transforms need to be compensated. Let

${\psi}_{2}({x}_{2},{y}_{2})=-\genfrac{}{}{0.1ex}{}{2\pi {z}_{2}}{\lambda}+{\psi}_{2}^{\prime}({x}_{2},{y}_{2})-\genfrac{}{}{0.1ex}{}{\pi \left({x}_{2}^{2}+{y}_{2}^{2}\right)}{\lambda {z}_{2}}.$.

Performing FrT of the phase function *ψ*
_{2}(*x*
_{2},*y*
_{2}) will get *G*̂(*x*
_{1},*y*
_{1}) with some accompanying terms:

$$=\genfrac{}{}{0.1ex}{}{\mathrm{exp}\left(\genfrac{}{}{0.1ex}{}{j2\pi {z}_{2}}{\lambda}\right)}{j\lambda {z}_{2}}\int \int \mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{-j2\pi {z}_{2}}{\lambda}+j{\psi}_{2}^{\prime}({x}_{2},{y}_{2})-\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{2}^{2}+{y}_{2}^{2}\right)}{\lambda {z}_{2}}\right]$$

$$\times \mathrm{exp}[\genfrac{}{}{0.1ex}{}{j\pi}{\lambda {z}_{2}}{({x}_{2}-{x}_{1})}^{2}+\genfrac{}{}{0.1ex}{}{j\pi}{\lambda {z}_{2}}{\left({y}_{2}-{y}_{1}\right)}^{2}]{\mathrm{dx}}_{2}{\mathrm{dy}}_{2}$$

$$=\genfrac{}{}{0.1ex}{}{1}{j\lambda {z}_{2}}\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{1}^{2}+{y}_{1}^{2}\right)}{\lambda {z}_{2}}\right]\int \int \mathrm{exp}\left[j{\psi}_{2}^{\prime}({x}_{2},{y}_{2})\right]\mathrm{exp}\left\{\genfrac{}{}{0.1ex}{}{-j2\pi \left({x}_{2}{x}_{1}+{y}_{2}{y}_{1}\right)}{\lambda {z}_{2}}\right\}{\mathrm{dx}}_{2}{\mathrm{dy}}_{2}\text{}$$

$$=\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{1}^{2}+{y}_{1}^{2}\right)}{\lambda {z}_{2}}\right]\mathrm{FT}\left\{\mathrm{exp}\left[j{\psi}_{2}^{\prime}({x}_{2},{y}_{2})\right];\lambda ;{z}_{2}\right\}$$

$$=\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{1}^{2}+{y}_{1}^{2}\right)}{\lambda {z}_{2}}\right]\hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{G}({x}_{1},{y}_{1})\right].$$

After rearrangement of the above equation, we get

Referring to Eq. (4), POM_{1} must contain the information of exp[*jψ _{t}*(

*x*

_{1},

*y*

_{1})], which, together with

*G*̂(

*x*

_{1},

*y*

_{1}), is Fourier-transformed to produce the approximation

*g*̂(

*x*

_{0},

*y*

_{0}) in the Fresnel domain. Let

It can be shown below from Eq. (8) that the approximation image *g*̂(*x*
_{0},*y*
_{0}) can be produced when the phase functions POM_{2} and POM_{1} are Fresnel-transformed with a placement of *z*
_{1} and *z*
_{2} distances as indicated in Fig. 1.

$$=\mid \mathrm{FrT}\left\{\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{1}^{2}+{y}_{1}^{2}\right)}{\lambda {z}_{2}}\right]\hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{G}({x}_{1},{y}_{1})\right]\times \mathrm{exp}\left[j{\psi}_{1}({x}_{1},{y}_{1})\right];\lambda ;{z}_{1}\right\}\mid $$

$$=\mid \mathrm{FrT}\left\{\mathrm{exp}\left(\genfrac{}{}{0.1ex}{}{-j2\pi {z}_{1}}{\lambda}\right)\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{-j\pi \left({x}_{1}^{2}+{y}_{1}^{2}\right)}{\lambda {z}_{1}}\right]\hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{t}({x}_{1},{y}_{1})\right];\lambda ;{z}_{1}\right\}\mid $$

$$=\mid \genfrac{}{}{0.1ex}{}{\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j2\pi {z}_{1}}{\lambda}\right]}{j\lambda {z}_{1}}\int \int \mathrm{exp}\left(\genfrac{}{}{0.1ex}{}{-j2\pi {z}_{1}}{\lambda}\right)\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{-j\pi \left({x}_{1}^{2}+{y}_{1}^{2}\right)}{\lambda {z}_{1}}\right]\hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{t}({x}_{1},{y}_{1})\right]\times \mathrm{exp}[\genfrac{}{}{0.1ex}{}{j\pi}{\lambda {z}_{1}}{\left({x}_{1}-{x}_{0}\right)}^{2}+\genfrac{}{}{0.1ex}{}{j\pi}{\lambda {z}_{1}}{\left({y}_{1}-{y}_{0}\right)}^{2}]{\mathrm{dx}}_{1}{\mathrm{dy}}_{1}\mid $$

$$=\mid \genfrac{}{}{0.1ex}{}{1}{j\lambda {z}_{1}}\mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{0}^{2}+{y}_{0}^{2}\right)}{\lambda {z}_{1}}\right]\int \int \hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{t}({x}_{1},{y}_{1})\right]\times \mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{-j2\pi \left({x}_{1}{x}_{0}+{y}_{1}{y}_{0}\right)}{\lambda {z}_{1}}\right]\text{d}{x}_{1}\text{d}{y}_{1}\mid $$

$$=\mid \mathrm{exp}\left[\genfrac{}{}{0.1ex}{}{j\pi \left({x}_{0}^{2}+{y}_{0}^{2}\right)}{\lambda {z}_{1}}\right]\mathrm{FT}\left\{\hat{G}({x}_{1},{y}_{1})\mathrm{exp}\left[j{\psi}_{t}({x}_{1},{y}_{1})\right]\right\}\mid $$

$$=\hat{g}({x}_{0},{y}_{0}).$$

Correlation coefficient *ρ*[20] is used to evaluate the quality of the approximation image *g*̂(*x*
_{0}, *y*
_{0}) and the sensitivity of the secret keys (i.e., *λ*, *z*
_{1}, and *z*
_{2}):

where *E*[·] denotes the expectation value operator. The value of *ρ* achieves a maximum of unity if *g*̂(*x*
_{0},*y*
_{0}) is perfectly matched to *g*(*x*
_{0}, *y*
_{0}). Besides, *ρ* can be also used to measure the correlation between *G*(*x*
_{1},*y*
_{1}) and *G*̂(*x*
_{1},*y*
_{1}).

For convenience, the formulas in the proposed phase retrieval algorithm are summarized as follows:

Also the target image reconstruction can be expressed as

Compared to the previous literatures [18, 20, 28–30], the proposed MGSA performs only a single stage of FTs, rather than two stages of FrTs in each iteration. Clearly, only the iterative process between the input and the transform planes is required. On the contrary, the previous literatures perform POCSA not only between the input and the transform planes, but also between the transform and the output planes. Substantially, the proposed phase retrieval scheme is more efficient with a shorter iterative process.

## 3. Simulation and analysis

The computer simulation of the proposed phase retrieval algorithm is performed to validate its correctness. Also, sensitivity on the wavelength and the axial shift is studied. Figures 3(a) and 3(b) show two target images, “Baboon” and “Lena”, respectively, which are of 256×256 pixels in dimensions and 256 grayscales in intensity. The sizes of both POMs are the same as that of the target image. Assume that the wavelength *λ* of the incident unit plane wave is 632.8 nm, the position parameters *z*
_{1} and *z*
_{2} are 50 and 30 mm, respectively, and the aperture of the output plane is 3.2×3.2 mm^{2}. In order to implement the proposed scheme digitally, we need to perform discrete FrT. The discrete Fresnel transform (DFrT) and its inverse transform IDFrT defined in Ref. [31] are used in the following process. Based on Eqs. (10a)–(10c), POM_{2} can be calculated and the phase distributions *ψ*
_{2}(*x*
_{2}, *y*
_{2}) for both target images are shown in Figs. 4(a) and 4(d). Based on Eqs. (10d) and (10e), POM_{1} can be obtained and the noise-like phase distributions *ψ*
_{1}(*x*
_{1}, *y*
_{1}) are shown in Fig. 4(b) and 4(e). The target image is reconstructed based on these two POMs via Eq. (11). Figures 4(c) and 4(f) show the reconstructed images which own a correlation coefficient of at least 0.999 with respect to the true target image.

Furthermore, the proposed MGSA is faster and more efficient than the traditional POCSA [18, 20, 28–30]. Tables 1 and 2 show the results of correlation coefficient *ρ* under different numbers of fast Fourier transform (FFT) and inverse fast Fourier transform (IFFT) operations used in the MGSA and the POCSA processes for “Baboon” and “Lena”, respectively. Each iteration in POCSA contains two stages of FFT and two stages of IFFT [18, 20, 28–30] (four stages totally). However, only two stages of FFT/IFFT are needed in the MGSA (referring to Fig. 2). At the same number of FFT and IFFT operations, MGSA can obtain higher *ρ*‘s than POCSA. According to Tables 1&2, requesting *ρ*>0.999 will take about 1024 FFT/IFFT operations for POCSA, but only 512 for the proposed MGSA. In this view, the proposed phase-retrieval algorithm is more efficient than that of the previous literatures [18, 20, 28–30].

Consider the phase distributions shown in Figs. 4(a) and (d), which are not random enough and could be vulnerable to known-plaintext attack. In order to further enhance the security of the proposed MGSA, the term *g*(*x*
_{0}, *y*
_{0}) in Eq. (2) can be replaced with *g*(*x*
_{0},*y*
_{0})exp[*jψ _{r}*(

*x*

_{0},

*y*

_{0})], where

*ψ*(

_{r}*x*

_{0},

*y*

_{0}) is a computer-generated random number within the interval (-

*π*,+

*π*]. We call this

*phase pre-scrambling*. Based on Eqs. (2)–(8) and the fact that only intensities are detected by optical detectors, it is clear that the approximated image thus derived (denoted as

*g*̃(

*x*

_{0},

*y*

_{0})) can still be highly correlated with

*g*(

*x*

_{0},

*y*

_{0}) when correct keys are available, where |

*g*̃(

*x*

_{0},

*y*

_{0})|≈|

*g*̂(

*x*

_{0},

*y*

_{0})exp[

*jψ*̂

*(*

_{r}*x*

_{0},

*y*

_{0})]|≈|

*g*̂(

*x*

_{0},

*y*

_{0})| and

*ψ*̂

*(*

_{r}*x*

_{0},

*y*

_{0}) is the approximation of the phase term

*ψ*(

_{r}*x*

_{0},

*y*

_{0}). Based on such a phase pre-scrambling by

*ψ*(

_{r}*x*

_{0},

*y*

_{0}), both POM

_{2}and POM

_{1}are produced in random noise-like phase distributions. Figures 5(a), 5(b), 5(d), and 5(e) show the derived phase distributions based on this modification, in comparison to Figs. 4(a), 4(b), 4(d), and 4(e), respectively. This apparently enhances the security of the proposed method. Furthermore, the correlation coefficients between the original and the decrypted images, shown in Figs. 5(c) and 5(f), are also higher than 0.999.

The generated phase functions, *ψ*
_{2}(*x*
_{2}, *y*
_{2}) and *ψ*
_{1}(*x*
_{1}, *y*
_{1}), can be transmitted over digital communication channels. Authorized receivers with correct keys, *z*
_{1}, *z*
_{2}, and *λ*, can transform *ψ*
_{2}(*x*
_{2}, *y*
_{2}) and *ψ*
_{1}(*x*
_{1}, *y*
_{1}) based on the proposed MGSA to obtain the target image. It is hard for intruders in the Internet to decipher the intercepted phase functions by randomly selecting *z*
_{1}, *z*
_{2}, and *λ*.

## 4. Conclusion

An improved phase retrieval algorithm, MGSA, for lensless optical security systems based on computer generated POMs is proposed. It is featured of less iteration without degrading the reconstruction quality. This algorithm is both straightforward and efficient compared with previous works [18, 20, 28–30]. Mathematical derivation of the related formulas is presented and the simulation results show that it works very well. Also the sensitivity of the reconstructed image on the parameters such as the wavelength and the axial positions of the POMs are investigated. Numerical results indicate that precisions of the wavelength and the POMs positions (i.e., the keys) are crucial to the reconstruction quality, which then ensures the security of the system. Another advantage is that the POM data generated can be transmitted via digital communication channels and then recovered optically or digitally.

## Acknowledgements

This research was supported in part of National Science Council under contract number NSC 97-2221-E-235 -003 -.

## References and links

**1. **H. T. Chang, “Image encryption using separable amplitude-based virtual image and iteratively retrieved phase information,” Opt. Eng. **40**, 2165–2171 (2001) [CrossRef]

**2. **G.-S. Lin, H. T. Chang, W.-N. Lie, and C.-H. Chuang, “A public-key-based optical image cryptosystem based on data embedding techniques,” Opt. Eng. **42**, 2331–2339 (2003) [CrossRef]

**3. **P. Réfrégier and B. Javidi, “Optical image encryption based on input plane and Fourier plane random encoding,” Opt. Lett. **20**, 767–769 (1995) [CrossRef] [PubMed]

**4. **O. Matoba and B. Javidi, “Encrypted Optical Memory Using Multi-Dimensional Keys,” Opt. Lett. **24**, 762–765 (1999) [CrossRef]

**5. **O. Matoba and B. Javidi, “Encrypted optical storage with wavelength-key and random phase codes,” Appl. Opt. **38**, 6785–6790 (1999) [CrossRef]

**6. **E. Tajahuerce, O. Matoba, S. C. Verrall, and B. Javidi, “Optoelectronic information encryption with phase shifting interferometry,” Appl. Opt. **39**, 2313–2320 (2000) [CrossRef]

**7. **G. Unnikrishnan, J. Joseph, and K. Singh, “Optical encryption by double-random phase encoding in the fractional Fourier domain,” Opt. Lett. **25**, 887–889 (2000) [CrossRef]

**8. **B. Zhu, S. Liu, and Q. Ran, “Optical image encryption based on multifractional Fourier transforms,” Opt. Lett. **25**, 1159–1161 (2000) [CrossRef]

**9. **Y. Zhang, C.-H. Zheng, and N. Tanno, “Optical encryption based on iterative fractional Fourier transform,” Opt. Commun. **202**, 277–285 (2002) [CrossRef]

**10. **B. Hennelly and J. T. Sheridan, “Optical image encryption by random shifting in fractional Fourier domains,” Opt. Lett. **28**, 269–271 (2003) [CrossRef] [PubMed]

**11. **P. C. Mogensen and J. Glückstad, “Phase-only optical encryption,” Opt. Lett. **25**, 566–568 (2000) [CrossRef]

**12. **X. Peng, L. Yu, and L. Cai, “Double-lock for image encryption with virtual optical wavelength,” Opt. Express **10**, 41–45 (2002) [PubMed]

**13. **S. Sinzinger, “Microoptically integrated correlators for security applications,” Opt. Commun. **209**, 69–74 (2002) [CrossRef]

**14. **R. K. Wang, I. A. Watson, and C. Chatwin, “Random phase encoding for optical secutity,” Opt. Eng. **35**, 2464–2469 (1996) [CrossRef]

**15. **Y. Li, K. Kreske, and J. Rosen, “Security and encryption optical systems based on a correlator with significant output image,” Appl. Opt. **39**, 5295–5301 (2000) [CrossRef]

**16. **H. T. Chang, W. C. Lu, and C. J. Kuo, “Multiple-phase retrieval for optical security systems by use of random-phase encoding,” Appl. Opt. **41**, 4825–4834 (2002) [CrossRef] [PubMed]

**17. **Y. C. Chang, Hsuan T. Chang, and C. J. Kuo, “Hybrid image cryptosystem based on dyadic phase displacement in the Fourier domain,” Opt. Commun. **236**, 245–257 (2004) [CrossRef]

**18. **G. Situ and J. Zhang, “A cascaded iterative Fourier transform algorithm for optical security applications,” Optik **114**, 473–477 (2003) [CrossRef]

**19. **L. G. Neto and Y. Sheng, “Optical implementation of image encryption using random phase encoding,” Opt. Eng. **35**, 2459–2463 (1996) [CrossRef]

**20. **G. Situ and J. Zhang, “A lensless optical security system based on computer-generated phase only masks,” Opt. Commun. **232**, 115–122 (2004) [CrossRef]

**21. **J. W. Goodman, *Introduction to Fourier Optics* (McGraw-Hill, New York, 1968)

**22. **H.P. Herzig (Ed.), *Micro-Optics* (Taylor & Francis, London, UK, 1996)

**23. **B. Javidi, G. Zhang, and J. Li, “Encrypted optical memory using double-random phase encoding,” Appl. Opt. **36**, 1054–1058 (1997) [CrossRef] [PubMed]

**24. **X. Tan, O. Matoba, T. Shimura, K. Kuroda, and B. Javidi, “Secure optical storage that uses fully phase encryption,” Appl. Opt. **39**, 6689–6694 (2000) [CrossRef]

**25. **X. Tan, O. Matoba, T. Shimura, and K. Kuroda, “Improvement in holographic storage capacity by use of double-random phase encryption,” Appl. Opt. **40**, 4721–4727 (2001) [CrossRef]

**26. **R. W. Gerchberg and W. O. Saxton, “Phase determination for image and diffraction plane pictures in the electron microscope,” Optik **34**, 275–284 (1971)

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

**28. **S. Deng, L. Liu, H. Lang, W. Pan, and D. Zhao, “Hiding an image in cascaded Fresnel digital holograms,” Chin. Opt. Lett. **4**, 268–271 (2006)

**29. **D. Abookasis, O. Montal, O. Abramson, and J. Rosen, “Watermarks encrypted in a concealogram and deciphered by a modified joint-transform correlator,” Appl. Opt. **44**, 3019–3023 (2005) [CrossRef] [PubMed]

**30. **J. Rosen, “Learning in correlators based on projections onto constraint sets,” Opt. Lett. **18**, 1183–1815 (1993) [CrossRef] [PubMed]

**31. **W. X. Cong, N. X. Chen, and B.Y. Gu, “Phase retrieval in the Fresnel transform system: a recursive algorithm,” J. Opt. Soc. Am. A **16**, 1827–1830 (1999) [CrossRef]