*Article* **Scaling of Three-Dimensional Computer-Generated Holograms with Layer-Based Shifted Fresnel Di**ff**raction**

### **Hao Zhang \*, Liangcai Cao and Guofan Jin**

State Key Laboratory of Precision Measurement Technology and Instruments, Department of Precision Instrument, Tsinghua University, Beijing 100084, China; clc@tsinghua.edu.cn (L.C.); jgf-dpi@tsinghua.edu.cn (G.J.)

**\*** Correspondence: haozhang274@tsinghua.edu.cn

Received: 27 March 2019; Accepted: 22 May 2019; Published: 24 May 2019

**Abstract:** Holographic three-dimensional (3D) displays can reconstruct a whole wavefront of a 3D scene and provide rich depth information for the human eyes. Computer-generated holographic techniques offer an efficient way for reconstructing holograms without complicated interference recording systems. In this work, we present a technique for generating 3D computer-generated holograms (CGHs) with scalable samplings, by using layer-based diffraction calculations. The 3D scene is partitioned into multiple layers according to its depth image. Shifted Fresnel diffraction is used for calculating the wave diffractions from the partitioned layers to the CGH plane with adjustable sampling rates, while maintaining the depth information. The algorithm provides an effective method for scaling 3D CGHs without an optical zoom module in the holographic display system. Experiments have been performed, demonstrating that the proposed method can reconstruct quality 3D images at different scale factors.

**Keywords:** holography; computer-generated hologram; holographic display; 3D display

#### **1. Introduction**

The holographic three-dimensional (3D) display is a powerful way of providing depth information for the human eyes, since it can reconstruct the entire whole optical wavefront of the 3D scene [1]. With the development of computing technology and spatial light modulators (SLMs), 3D computer-generated holograms (CGHs) can be reconstructed dynamically without the complicated interference recording systems [2]. The CGH algorithms are used to encode the mathematical representations of the 3D scenes into the holograms, which are closely related to image reconstruction quality and computing efficiency.

Point-based and polygon-based algorithms are commonly used when synthesizing CGHs for various types of 3D scenes [3–6]. These algorithms simulate the wave propagation process from the 3D scene to the CGH plane. The 3D objects are often segmented into many point sources or polygons, and they provide accurate geometric relationships of the 3D scenes. Since the computational load increases with the number of primitives. The computing efficiency is aggravated drastically when calculating the CGH of a complicated 3D scene.

Stereogram-based algorithms can use computer graphics-rendering techniques to process different formats of 3D objects efficiently [7–9]. During the calculations, each holographic element (hogel) of the holographic stereogram is calculated, based on a two-dimensional (2D) parallax image; hence, depth performances are often limited in these stereogram-based CGHs. Recently, several techniques have been developed to improve the depth performances of the stereogram-based algorithms [10–12]. These algorithms have integrated physically-based algorithms and stereogram-based algorithms, in order to reconstruct the accurate depth information of the 3D scenes. However, in stereogram-based algorithms, each hogel corresponds to a viewpoint during the rendering procedure, and this means that we need to render multiple times before implementing the CGH algorithm, which would increase the rendering time. Also, the continuity of the motion parallax and occlusion would be affected by the segmentation setup of the CGH.

Layer-based techniques are efficient in calculating the CGH, by partitioning the 3D scene in multiple layers, according to the depth information [13–16]. Fast Fourier transform (FFT)-based diffraction can be used in simulating the wave propagation processes, from the layers to the CGH plane. However, there are sampling restrictions in the current layer-based 3D CGH algorithms. In FFT-based Fresnel diffraction calculations, the sampling interval on the observation plane is proportional to the propagation distance, while the sampling interval would remain unchanged in convolution-based Fresnel algorithm and angular spectrum method [17]. Hence, the flexibility of the holographic 3D display technique, based on the layer-based algorithm, is limited. Although zero-padding and resampling methods can be used to change the sampling rates, this would take a larger amount of memory and a longer computation time, due to more sampling points that need to be processed.

A double Fresnel transform algorithm was developed, to calculate the wave propagation, with variable magnification, by splitting the propagation distance into two partial steps, which has been used in generating layer-based 3D CGHs with amplitude modulation [18,19]. Since there are minimum propagation distances for the two propagation steps, the propagation distance and the adjustable range of the sampling rates are limited. The wavefront recording plane (WRP) method with non-uniform fast Fourier transform (NUFFT) was used to simplify the 3D CGH calculation [20]. Since the WRP needs to be located close to the object points, the overall depth range is limited during reconstruction.

Chirp Z-transform with Bluestein's algorithm provides a powerful way to perform discrete Fourier transforms with different scaling factors [21,22], which has been used in the Fourier transform-based Fresnel diffraction calculations. Shifted Fresnel diffraction was introduced to simulate the scalable Fresnel diffraction between shifted parallel planes in computational holography [23]. Later, aliasing-reduced shifted and scaled (ARSS) Fresnel diffraction was introduced in the algorithm, to reduce aliasing by a band-limiting function [24], which was also used in lenseless holographic projections [25]. Similar algorithms was also introduced in digital holography to perform magnified reconstructions of the digital holograms [26]. However, the above applications, based on the scaling diffractions, were mainly focused on 2D reconstructions. Also, the sampling technique in 3D CGH still needs further investigation.

In this work, layer-based shifted Fresnel diffraction is implemented to calculate the wave propagation from the partitioned layers to the CGH plane, with adjustable sampling rates. Shifted Fresnel diffraction can simulate the wave propagation between parallel planes with scalable sampling, by using a Bluestein transform. This layer-based processing with scalable diffraction calculation provides an effective way of generating a CGH, with flexible sampling of the 3D scene. By implementing a scalable diffraction calculation into layer-based processing for 3D CGH generation, the reconstructed 3D images are zoomable without an additional optical zoom module. The sampling parameters of the reconstructed 3D images can be adjusted accordingly. Numerical and optical experiments are performed, to demonstrate the effectiveness of our proposed method.

#### **2. Layer-Based Shifted Fresnel Di**ff**raction**

The Fresnel diffraction of the 2D optical wave field can be expressed by:

$$\mathcal{U}(u,v) = \frac{\exp(jkz)}{j\lambda z} \exp\left[j\frac{k}{2z}(u^2+v^2)\right] \left[\int\_{\mathcal{O}} \mathcal{U}(\mathbf{x},y) \exp\left[j\frac{k}{2z}(\mathbf{x}^2+y^2)\right] \exp\left[-j2\pi(f\_x\mathbf{x}+f\_y y)\right] d\mathbf{x}dyd,\tag{1}$$

where (*x*, *y*) and (*u*, *v*) denote the coordinates in the source plane and observation plane, respectively. *z* is the wave propagation distance, *k* = 2/λ is the wave number, and *fx* and *fy* are calculated as:

$$f\_x = \frac{u}{\lambda z'} \ f\_y = \frac{v}{\lambda z}. \tag{2}$$

By imposing a 2D Fourier transform, Equation (1) can be represented as:

$$\mathcal{U}(\boldsymbol{u},\boldsymbol{v}) = \frac{\exp(j k \boldsymbol{z})}{j \lambda \boldsymbol{z}} \exp\left[j \frac{k}{2\boldsymbol{z}} (\boldsymbol{u}^2 + \boldsymbol{v}^2)\right] \mathcal{F}\left\{\mathcal{U}(\boldsymbol{x},\boldsymbol{y}) \exp\left[j \frac{k}{2\boldsymbol{z}} (\boldsymbol{x}^2 + \boldsymbol{y}^2)\right]\right\},\tag{3}$$

where *F* denotes the Fourier transform. According to the sampling rule of the discrete Fourier transform, the following relation exists:

$$
\Delta f\_{\mathbf{x}} = \frac{1}{N\_{\mathbf{x}} \Delta \mathbf{x}'} \,\,\Delta f\_{\mathbf{y}} = \frac{1}{N\_{\mathbf{y}} \Delta \mathbf{y}'} \,\tag{4}
$$

where *Nx* and *Ny* are the sampling numbers in *x* and *y* coordinates, Δ*x* and Δ*y* are the sampling intervals in the *x* and *y* directions, and Δ*fx* and Δ*fy* are the frequency sampling intervals, respectively. Hence, the sampling distances in the observation plane can be deduced as:

$$
\Delta \mu = \frac{\lambda z}{N\_{\text{X}} \Delta \mathbf{x}}, \ \Delta v = \frac{\lambda z}{N\_{\text{y}} \Delta y}. \tag{5}
$$

where *Nx*Δ*x* and *Ny*Δ*y* denote the size of the sampling window in the source plane. From Equation (5), we can see that the sampling interval in the observation plane is proportional to the wave propagation distance, and inversely proportional to the size of the sampling window. The sampling ranges in the observation plane can be given by:

$$L\_{\rm u} = N\_{\rm x} \Delta \mu = \frac{\lambda z}{\Delta x} L\_{\rm v} = N\_{\rm y} \Delta v = \frac{\lambda z}{\Delta y'} \tag{6}$$

where *Lu* and *Lv* are the sampling ranges along the *u* and *v* axes in the observation plane. We can see that the sampling parameters in the observation plane cannot be adjusted freely in the Fresnel diffraction calculation, which are determined by the sampling parameters of the source plane and the wave propagation distance.

Shifted Fresnel diffraction can overcome the limitations of traditional Fourier-based Fresnel diffraction, by allowing for arbitrary sampling intervals in the source plane and the observation plane. Figure 1 is the geometric setup of the shifted Fresnel diffraction. (*x*0,*y*0) is the center of the sampling grid in the source plane (*x*,*y*), with a resolution of *P* × *Q*, and (*u*0,*v*0) is the center of the sampling grid in the observation plane (*u*,*v*), with a resolution of *M* × *N*. Δ*x* and Δ*y* are the sampling intervals in the source plane, and Δ*u* and Δ*v* are the sampling intervals in the observation plane. The coordinates of these two grids are determined by:

$$\begin{aligned} \mathbf{x}\_{\mathcal{P}} &= \mathbf{x}\_0 + p\Delta \mathbf{x}\_{\prime} \ y\_{\mathcal{I}} = y\_0 + q\Delta y\_{\prime} \\ \mathbf{u}\_m &= \mathbf{u}\_0 + m\Delta \mathbf{u}\_{\prime} \ v\_n = v\_0 + n\Delta v\_{\prime} \end{aligned} \tag{7}$$

where *p*, *q*, *m*, *n* are the integer indices:

$$\begin{aligned} -\frac{P}{2} &\le p \le \frac{P}{2} - 1, -\frac{Q}{2} \le q \le \frac{Q}{2} - 1, \\\ -\frac{M}{2} &\le m \le \frac{M}{2} - 1, -\frac{N}{2} \le n \le \frac{N}{2} - 1. \end{aligned} \tag{8}$$

**Figure 1.** Geometry of shifted Fresnel diffraction.

The sampled 2D Fresnel diffraction can be expressed as:

$$\begin{split} \mathcal{U}(m,n) &= \frac{\exp(jkz)}{j\lambda z} \exp\left[j\frac{k}{2\pi}(\mu\_{m}^{2}+\upsilon\_{n}^{2})\right] \sum\_{p} \sum\_{q} \mathcal{U}(p,q) \exp\left[j\frac{k}{2\pi}(\upsilon\_{p}^{2}+y\_{q}^{2})\right] \exp\left[-j\frac{2\pi}{\lambda z}(\mathbf{x}\_{p}\mu\_{m}+y\_{q}\upsilon\_{n})\right] \\ &= \frac{\exp(jkz)}{j\lambda z} \exp\left[j\frac{k}{2\pi}(\mu\_{m}^{2}+\upsilon\_{n}^{2})\right] \exp\left[-j\frac{2\pi}{\lambda z}(\mathbf{x}\_{0}m\Delta u + y\_{0}m\Delta v)\right] \\ &\sum\_{p} \sum\_{q} \mathcal{U}(p,q) \exp\left[j\frac{k}{2\pi}(\upsilon\_{p}^{2}+y\_{q}^{2})\right] \exp\left[-j\frac{2\pi}{\lambda z}(\mathbf{x}\_{p}\mu\_{0}+y\_{q}\upsilon\_{0})\right] \exp\left[-j\frac{2\pi}{\lambda z}(pm\Delta x\Delta u + qn\Delta y \Delta v)\right]. \end{split} \tag{9}$$

Let:

$$\begin{aligned} 2pm &= p^2 + m^2 - \left( m - p \right)^2, \\ 2qn &= q^2 + n^2 - \left( n - q \right)^2. \end{aligned} \tag{10}$$

According to the convolution theorem, the optical wave field *U*(*m*,*n*) can be calculated as:

$$\mathcal{L}I(m,n) = \mathbb{C} \cdot F^{-1}\{F[a(p,q)]F[b(p,q)]\},\tag{11}$$

where:

$$\mathbf{C} = \frac{\exp\left(jkz\right)}{j\lambda z} \exp\left[j\frac{k}{2z}\left(\mu\_{tt}\,^2 + \upsilon\_{\mathbf{z}}\,^2\right)\right] \exp\left[-j\frac{2\pi}{\lambda z}(\mathbf{x}\_0 m\Delta u + y\_0 n\Delta v)\right] \exp\left[-j\pi\left(\frac{\Delta x\Delta u}{\lambda z}m^2 + \frac{\Delta y\Delta v}{\lambda z}n^2\right)\right] \tag{12}$$

$$a(p,q) = \mathcal{U}(p,q) \exp\left[j\frac{k}{2\pi}(\mathbf{x}\_{\parallel}\mathbf{}^2 + y\_q\mathbf{1})\right] \exp\left[-j\frac{2\pi}{\lambda z}(\mathbf{x}\_{\parallel}\mathbf{u}\_0 + y\_q\mathbf{u}\_0)\right] \exp\left[-j\pi\left(\frac{\Delta x\Delta u}{\lambda z}p^2 + \frac{\Delta y\Delta v}{\lambda z}q^2\right)\right] \tag{13}$$

$$b(p,q) = \exp\left[j\pi \left(\frac{\Delta x \Delta u}{\lambda z} p^2 + \frac{\Delta y \Delta v}{\lambda z} q^2\right)\right].\tag{14}$$

By imposing the required sampling rules in the source plane and the observation plane, diffraction patterns with different scale factors can be reconstructed with the help of shifted Fresnel diffraction. Figure 2 numerically demonstrates the shifted Fresnel diffraction of a square aperture, with three different scale factors. The size of the aperture is 2.4 mm × 2.4 mm. The sampling number in the source plane is 1024 × 1024, and the sampling interval is 8 μm. The wave propagation distance is 300 mm, and the wavelength used in the simulation is 633 nm. The diffraction patterns in the observation plane are computed with different scale factors. Figure 2a–c shows the diffraction patterns, with Δ*u* = 6 μm, Δ*u* = 8 μm, and Δ*u* = 10 μm, respectively.

**Figure 2.** Shifted Fresnel diffraction with different scale factors: (**a**) Δ*u* = 6 μm, (**b**) Δ*u* = 8 μm, (**c**) Δ*u* = 10 μm.

In layer-based 3D CGH synthesis, the 3D scene is firstly partitioned into multiple layers, which are parallel to the CGH plane. Then, the wave fields diffracted from the partitioned layers to the CGH plane are calculated. By using shifted Fresnel diffraction in the diffraction simulations from each layer to the CGH plane, the sampling parameters addressed on the 3D scene can be adjusted accordingly. The diagram of the layer-based shifted Fresnel diffraction is shown in Figure 3. The 3D model is partitioned into different parallel layers, based on its depth information. For each layer, amplitude distribution is extracted according to the rendered image. In order to simulate the diffusive effect of the object surface, random phase distribution is added to each layer. The complex distribution of each layer can be calculated by using shifted Fresnel diffractions with preset sampling rules. The complex distribution of the CGH plane can be acquired after adding the contributions from all the layers [27]. After obtaining the field distribution of the hologram plane, the phase distribution is then extracted for uploading to the phase SLM.

**Figure 3.** Diagram of the layer-based shifted Fresnel diffraction for 3D computer-generated holograms (CGHs) calculations.

#### **3. Reconstruction Results**

To demonstrate the effectiveness of the layer-based shifted Fresnel diffraction for 3D reconstruction, optical experiments are implemented. Figure 4 illustrates the optical diagram of the holographic reconstruction experiment. In the experimental work, the PLUTO phase-only SLM was used for optical reconstructions. The pixel number of the SLM was 1920 × 1080. The pixel pitch was 8 μm, and the SLM was addressed with 8-bit grayscale levels. The wavelength of the laser used in our experiment was 633 nm. A 4-f system was used in the reconstruction path, for filtering out the zero-order interruption and unwanted orders [28–30].

**Figure 4.** Optical diagram of the holographic reconstruction experiment.

When generating CGHs by using a layer-based algorithm with a shifted Fresnel diffraction, the source planes are the partitioned layers of the 3D scene, and the observation plane is the hologram plane. The original 3D scene used for generating the CGH was a dragon model located 200 mm behind the hologram plane. The dragon model was partitioned into 50 layers, according to its depth information. Since shifted Fresnel diffraction is used in the CGH calculation, the sampling intervals of the object layers can be adjusted. In order to scale the reconstructed 3D scene, the lateral sampling distances of the object layers were set to 6 μm, 8 μm, and 10 μm, respectively. Figure 5a,c,e presents the numerical reconstruction results when focusing on the head of the model. Figure 5b,d,f presents the numerical reconstruction results when focusing on the tail of the model.

**Figure 5.** Numerical reconstructions of CGHs with different sampling rates. (**a**,**c**,**e**) focus on the head; (**b**,**d**,**f**) focus on the tail.

Figure 6 demonstrates the optical reconstruction results of CGHs with different sampling rates, which are consistent with the numerical results shown in Figure 5. According to the numerical and optical demonstrations, we can see that the depth information of the CGHs with different scale factors could be reconstructed with high quality, which demonstrates that the reconstructed images can be scaled with precise control of the sampling rates. Since the scaled sampling was addressed on the object layers of the original 3D object, and because the reconstruction processes of different CGHs are based on the unified sampling rules, smaller sampling interval would lead to a smaller size for the reconstruction image.

**Figure 6.** Optical reconstructions of CGHs at different sampling rates. (**a**,**c**,**e**) focus on the head; (**b**,**d**,**f**) focus on the tail.

#### **4. Accuracy Analysis with Di**ff**erent Sampling Parameters**

To evaluate the accuracy of the proposed algorithm, a single slit with a size of 2.4 mm was placed at the center of the source plane, to calculate the diffraction fields with different parameters. The sampling interval on the source plane was 8 μm, and the sampling number was 1024. The wavelength was 633 nm. The accuracy was evaluated with a signal-to-noise ratio (SNR) of the wavefield that was compared to the Rayleigh–Sommerfeld diffraction integral.

Figure 7 illustrates the accuracies with different sampling rates in the observation plane as a function of the propagation distance. Due to the paraxial approximation, the accuracies of the shifted-Fresnel diffraction were low for the short propagation distances. The SNR obviously improved with an increase in the propagation distance. On the other hand, the sampling rate in the observation plane also affected the calculation accuracy. Larger sampling intervals would lead to a larger sampling range, according to Equation (6), which induces a calculation that satisfies the paraxial approximation at a longer propagation distance. Additionally, all of the SNRs are larger than 30 dB when the propagation distance is larger than 150 mm.

**Figure 7.** Accuracies of the diffraction calculation with different sampling rates in the observation plane.

Figure 8 illustrates the amplitude profiles of the diffraction fields with different scaling factors, calculated by the proposed method and the Rayleigh–Sommerfeld diffraction integral. In the simulation, the propagation distance was set to 300 mm. From the figure, we can see that the calculation results were in excellent agreement.

**Figure 8.** Amplitude profiles of the diffraction fields at different propagation distances: (**a**) Shifted Fresnel diffraction, (**b**) Rayleigh–Sommerfeld diffraction integral.

#### **5. Conclusions**

In conclusion, we propose an effective method for scaling 3D CGHs, using layer-based shifted Fresnel diffraction. During the calculation, the wave propagation from the partitioned layers to the CGH plane can be calculated with a flexible sampling rate; hence, the scaling display of the 3D image can be realized by changing the sampling rates of the object layers. A phase-only SLM is used for optical reconstruction, and it demonstrates that the proposed system can scale a 3D scene using different scale factors, without the need for an optical zoom module.

**Author Contributions:** Conceptualization, H.Z. and L.C.; methodology, H.Z.; software, H.Z.; validation, L.C. and G.J.; writing—original draft preparation, H.Z; writing—review and editing, H.Z., L.C., and G.J.; funding acquisition, H.Z. and L.C.

**Funding:** This research was funded by the National Natural Science Foundation of China (NSFC), grant number 61875105.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
