2 – Fast Methods for Computer-Generated Holography




Abstract




In this chapter, a number of quick methods for generating digital Fresnel holograms are introduced. For an object comprising a small number of depth planes, the layer-based method implemented in the Fourier space is preferred for fast hologram generation. The point-based method is more suitable for generating objects with a large number of object points that are scattered over a wide range of distances from the hologram plane. To enhance the speed of hologram generation with the point-based method, different variants and sizes of the look-up-table (LUT) algorithms to trade-off computation time are described. A number of methods based on the concept of a wavefront recording plane (WPR) are presented. Being different from the LUT approach, the WRP methods speed up the hologram-generation process. Instead of generating the full hologram for each object point, only a small area of fringe patterns is computed on a WRP that is at close proximity to the object space. The computation time is substantially reduced. Further enhancement of the computation speed is attained with the warped wavefront recording plane (WWRP) method. A 3-D object image is decomposed into a 2-D intensity image and a depth map. The intensity image is used to generate an interim hologram on a WRP. Different regions of the hologram fringes on the WRP are resized according to their distances (obtained from the depth map) from the hologram plane to generate a hologram from the WRP.





2 Fast Methods for Computer-Generated Holography




2.1 Introduction


In Chapter 1, computer-generated holography (CGH) was briefly reviewed. Basically, CGH is a simulation of optical holography whereby the interference pattern of light waves diffracted from the surface of a 3-D object is recorded on a light-sensitive medium known as a hologram. The process can be numerically computed by applying the Fresnel diffraction equation on a computer graphic model of the 3-D object. Although CGH is not a complicated process in itself, the amount of computation involved in generating a hologram is overwhelming, especially if the size of the hologram and the object space is large.


To estimate the number of arithmetic operations that are required to generate a hologram with the point-based method, the expression in Eq. (1.14) is rewritten:


H(m,n)=∑pAp×exp[iwn×(m−mp )2×δd2+(n−np)2×δd2+zp2],(2.1)

where wn=2πλ . Suppose a hologram is discretized into N rows and M columns, the total number of hologram pixels is M× N. Assuming the constant δd2 is known in advance, the total number of operations to compute a hologram is in the order of (M×N×P) , where P is the number of object points. Each operation is used to calculate the magnitude and phase of the optical wave of an object point that impinges on a hologram pixel. The arithmetic calculations involved in each operation are listed in Table 2.1.




Table 2.1 List of steps and the associated functions for computing the value of a hologram pixel of a single object point




































Arithmetic functions Quantity
Addition + 2
Subtraction 2
Multiplication × 4
Square a2 4
Square root 1
Exponential function exp() 1

As the pixel size of a hologram is in the order of wavelength of light, the number of pixels in a hologram is generally extremely large. For example, if the pixel size δd is 5 μm, a small hologram of 10 × 10 mm already contains two megapixels (2 × 106), which is equivalent to the number of pixels in a domestic television. If the number of object points is the same as the number of hologram pixels (i.e., PM× N), the number of operations is (2 × 1012). Even with modern computers, generating a hologram at a video rate of 25–30 frames per second is still an extremely challenging task.


Eq. (2.1) shows that the computation of each hologram pixel is independent of the others. A straightforward way of speeding up the hologram-generation process, therefore, is to employ a set of processors to work in parallel. Each processor will be dedicated to the computation of a unique hologram pixel. This kind of parallel operation can be realized with graphical processing units (GPUs). However, the speed enhancement is rather limited, and fast processors also lead to other problems such as size, power consumption, and heat dissipation. Numerous research attempts have been made since the early 1990s to speed up the hologram-generation process with fast algorithms that can be realized with parallel processors.


In this chapter, a number of fast hologram-generation methods that have been developed in the past decades will be reviewed. Some of these algorithms reduce the steps (or arithmetic operations) involved in the numerical realization of the CGH equations, while others are based on simplifying the hologram-generation process. This chapter describes a number of fast algorithms for generating a hologram from a set of object points.



2.2 Realization of CGH with Fourier Transform


As mentioned in Chapter 1, if a 3-D object space is composed of a large number of object points and few layers, the layer-based method is more computationally efficient than the point-based method. In the layer-based method, the object space is decomposed into a sequence of uniformly spaced image planes that are parallel to the hologram plane. Each image plane is located at an axial distance zj from the hologram, and the object points on it are represented with the intensity distribution Ij(m,n) . The distance zj is commonly referred to as the depth, and the image planes as the depth planes. Within the jth depth plane, all the object points have identical depth zj from the hologram, and the hologram Hj(m,n) can be generated as


Hj(m,n)=∑up=0M−1∑vp=0N−1Ij(up,vp)×fzp(m−up,n−vp;zj),(2.2)

where fzp() is the Fresnel zone plate. Equation (2.2) can be rewritten as a convolution operation:


(2.3)Hj(m,n)=Ij(m,n)*fzp(m,n;zj).

Equation (2.3), in turn, can be realized in the frequency space. Suppose Ij˜(ωm,ωn)  and fzp~(ωm,ωn;zj) represent the Fourier transform of Ij(m,n) and fzp(m,n;zj) , respectively, the hologram can be computed as


Hj(m,n)=F−1[Ij˜(ωm,ωn)fzp~(ωm,ωn;zj)]=F−1[F[Ij(m,n)]ℑ[fzp(m,n;zj)]],(2.4)

where ℑ[] is the forward Fourier transform, and ℑ−1[] is the inverse Fourier transform. The hologram of the entire object space is equal to the sum of the holographic signals from all the image planes.


From the software implementation point of view, Eq. (2.4) can be conducted more swiftly than the convolution operation in Eq. (2.3), as it only involves a multiplication operation for each point in the frequency space. Both the forward and the inverse Fourier transform can be realized with fast Fourier transform (FFT), which is more computationally efficient than the convolution operation.


By the same principle, the numerical reconstruction process in each focused plane can also be speeded up with FFT. The reconstructed image IR;j(m,n) at the jth focused plane can be obtained from the hologram as follows:


IR;j(m,n)=H(m,n)*fzp†(m,n;zj),(2.5)

which can be rewritten as


IR;j(m,n)=H(m,n)*fzp†(m,n;zj)=F−1{F[H(m,n)]×F[fzp†(m,n;zj)]}.(2.6)

From the software implementation point of view, realizing Eq. (2.6) with FFT is considerable faster than convolution. However, Eq. (2.4) also shows that the computation can be overwhelming with a large number of depth planes. In the rest of this chapter, some of the contemporary methods for fast generation of digital holograms will be explained.



2.3 Direct Look-Up Table Method


A common practice to decrease the computation time in software programming is to store in the computer memory precomputed results of some tedious functions that have a high chance of being repeatedly used in the future. Such an approach is known as the look-up-table (LUT). The use of LUTs for speeding up the CGH process can be traced back to the works of Lucente in the early 1990s [1]. From Eq. (2.1) and Table 2.1, it can be seen that lots of arithmetic calculations are involved in deriving the expression Ap×fzp(m−mp,n−np;zp) . These calculations can be waived if the values corresponding to all the combinations of the six variables {Ap,m,n,mp,np,zp} have been computed in advance and stored in the computer as a LUT. A hologram pixel at location (m,n) can then be obtained with the flowchart shown in Figure 2.1.





Figure 2.1 Flowchart for generating a hologram pixel with the LUT.


Assuming that the time taken to fetch the data and variables from the computer memory is negligible (which is generally valid for modern computers), generating a hologram pixel only requires one additional calculation (i.e., the shaded box in Figure 2.1) for each object point, which is significantly less than the computation listed in Table 2.1. For a hologram of size M×N , and an object space with P object points, the number of arithmetic operations is reduced to (M×N×P) additions. On the downside, a huge LUT is required to store the precomputed values. An estimation of the size of the LUT is given as follows. Suppose the object space is evenly partitioned into D image planes along the axial direction, and the intensity AP is quantized into L levels. All the image planes are assumed to have identical horizontal and vertical extents as the hologram (i.e., size = M×N ). The total combinations from the six variables {Ap,m,n,mp,np,zj} are L×M2×N2×D . If each entry of the LUT is represented by B bits, the size of the LUT will be (L×M2×N2×D×B) bits. As an example, a small square hologram (i.e., M=N ) of size 1024×1024 representing an object scene of similar dimension, with eight depth planes, L=256 , and B=16 (eight bits for each of the real and imaginary parts), will result in a LUT of 256×10242×10242×8×16 ≈ 36 000 terabits. Even with modern computers it is impractical, if not impossible, to set aside memory with size in the order of astronomical terabits for the LUT.


The memory problem can be alleviated with a hybrid approach, with some of the data being computed instead of storing in the LUT. A smaller LUT is constructed to store the values of the FZP function fzp(m−mp,n−np;zp) for all combinations of values of the five variables {m,n,mp,np,zp} . The hologram pixel of an object point can be computed by multiplying its intensity Ap with the retrieved value from the LUT, as shown in the flowchart in Figure 2.2. The size of the LUT will be decreased dramatically from 35 840 terabits to around 140 terabits, at the expense of increasing the computation loading with (M×N×P) multiplications. Even with the reduction in the size of the LUT, the memory requirement is still astronomical. In the subsequent subsections, different methods for trading off the size of the LUT and the speed enhancement in generating a hologram will be described.





Figure 2.2 Flowchart for generating a hologram pixel with a smaller LUT.



2.4 Novel Look-Up Table Method


As mentioned in Section 2.3, the basic LUT method requires precomputation and storage of the FZP for all the combinations of the six variables {Ap,m,n,mp,np,zp} , resulting in a LUT of formidable size even for a small hologram. Although the size of the LUT can be reduced by discarding the intensity parameter Ap , the amount of computer memory is still enormous. The novel look-up table (N-LUT) method was proposed in [2], with the objective of further reducing the size of the LUT. A closer look at the function fzp(m−mp,n−np;zp) reveals that it is simply a translation of the FZP fzp(0,0;zp) of a point source located at an axial distance zp from the hologram plane. Refer to Figure 1.14 on a layered partitioning of the object space. For each depth plane at zj , a function fzp(0,0;zj) , which is referred to as the principal fringe pattern (PFP) is defined. Based on the PFP, the FZP of object points at other locations, fzp(m−mp,n−np;zp), can be derived by translating the PFP. It is then sufficient to store the PFP for each depth plane in a small look-up table known as the N-LUT.


The construction of the N-LUT and its application in generating a hologram pixel from an object point are shown in Figure 2.3. Suppose the object space is divided into D depth planes within the range [z1,zD] , a set of PFP functions (PFP1 to PFPzD ) is computed and stored in the N-LUT. To compute a hologram for a 3-D object, each depth plane is visited one by one. In each depth plane, the corresponding PFP is retrieved from the N-LUT. For an object point located at (mp,np) on the jth depth plane and having intensity Ap , the PFP of the depth plane at zj is shifted to the position (mp,np) . The translated PFP is then multiplied by Ap and accumulated on the hologram. This process is repeated for all the object points in the set of scan planes, resulting in a hologram of the 3-D object space.





Figure 2.3 The N-LUT method.


Compared with the LUT method, the number of major arithmetic operations is slightly larger as the diffraction pattern of each object point is obtained by translating the PFP along the horizontal and vertical directions. Since the computation loading of these additional processes is negligible in practice, the overall complexity of the N-LUT method is about the same as that of the LUT method, comprising of (M×N×P) multiplications and an equal number of addition operations. However, the size of the N-LUT is substantially reduced to (M×N×D×B) bits, which is M×N times smaller than that in the LUT. Referring to the previous example, the size of the N-LUT for M=N=1024 , D=8 , and B=16 bits (2 bytes) is 1024×1024×8×2 ≈ 16.8 megabytes. Reducing the size of the N-LUT by taking into account the pixel size and the reconstruction distance of the digital hologram was conducted and reported in [3].



2.5 The Line Scanning Method


Although the N-LUT is relatively smaller in size than the LUT, it still requires quite a lot of computer memory to store. To alleviate the memory requirement, Z. Yang has proposed the line scanning (LS) method [4] to further reduce the size of the N-LUT. As mentioned in Section 2.4, an N-LUT has to store a complete image of the PFP for each of the depth planes, resulting in a LUT of (M×N×D×B) bits. A PFP is an FZP in the form of concentric circular rings centered at the origin. Adopting polar coordinates, the location at (x,y) can be expressed as x=rcosΘ , y=rsinΘ , where r and Θ are the radius and azimuth of the hologram plane. The FZP at depth plane zj can be rewritten in the polar form as


fzp(r,Θ;zj)=exp[iwn(rcosΘ)2+(rsinΘ)2+zj2]=exp[iwnr2+zj2]=G(r;zj).(2.7)

From the right-hand side of Eq. (2.7), it can be seen that the value of the function fzp(r,Θ;zj) is only dependent on the radius. If the radius is uniformly discretized into RN intervals, there will only be RN different values for the 1-D function G(r;zj) , with which the PFP on the corresponding depth plane can be generated. A LUT of size (RN×D×B) will be sufficient to encapsulate the PFPs of all the depth planes. Suppose the radius RN of an FZP (which is half the width or height of a square hologram) is 512 , the memory required to store the PFPs of a 3-D object scene with eight depth planes is (RN×D×B)=512×8×2 ≈ 8 k bytes, which is about 2000 times smaller than that required to store the N-LUT.


In computing the hologram, a scan line is retrieved for each object point and used to reconstruct the PFP for Θ=[0,2π) . A straightforward implementation is to compute each point of the PFP as


fzp(m,n;zj)=fzp(r cos Θ,r sin Θ;zj).(2.8)

From Eq. (2.8), it can be seen that each quadrant of the PFP is in fact a mirror image of the others, that is:


fzp(m,n;zj)|m≤0;n≥0=fzp(|m|,n;zj),(2.9)

and


(2.10)fzp(m,n;zj)|n≤0=fzp(m,|n|;zj).

Hence, computation of the PFP from the scan line can be simplified to the computation of the first quadrant of the PFP ( m;n≥0 ), and the remaining quadrants are obtained from the symmetrical relation in Eqs. (2.9) and (2.10). The process can be divided into three steps, as illustrated in Figure 2.4. In step 1, a scan line is generated. In step 2, the first quadrant of the PFP is computed from the scan line. Finally, in step 3, the remaining quadrants of the PFP are duplicated from the mirror images of the first quadrant.





Figure 2.4 The LS method for generating the PFP.


Expanding a radial line into a circular arc can be implemented with Bresenham’s line algorithm [5]. However, some of the pixels in the PFP may be missing. An improved method has been proposed in [6], whereby the “arithmetical circle algorithm” [7] is used to generate the PFP from the scan line.



2.6 The Split-Look-Up-Table (S-LUT) Framework


Another method of reducing the size of the LUT, known as the split-look-up-table method (S-LUT), was proposed by Y. Pan et al. [8]. Previous LUT methods are more or less based on storing the FZP function in its original form. The S-LUT method employs an approximation of the FZP that is formed by the product of a pair of 1-D functions. The hologram is generated by accumulating the approximated FZP function of each object point. Referring to Eq. (2.1), the FZP of the pth point in the object space is a function of five variables {m,n,mp,np,zp} . By grouping the pairs of terms {m,mp} and {n,np} , the FZP can be simplified as a function of three variables


fzp(Δm,Δn;zp)=exp[i2πλΔm2δd2+Δn2δd2+zp2],(2.11)

where Δm=(m−mp) and Δn=(n−np) . The hologram H(m,n) can be rewritten as


(2.12)H(m,n)=∑pAp exp[i2πλ−1Δm2δd2+Δn2δd2+zp2].

As Δm2=|Δm|2 and Δn2=|Δn|2 ,


H(m,n)=∑pAp exp[i2πλ−1|Δm|2δd2+|Δn|2δd2+zp2].(2.13)

Assuming |Δmδd|≪zp and |Δnδd|≪zp , Eq. (2.13) can be approximated as


H(m,n)=∑pAp exp[i2πλ−1|Δm|2δd2+zp2]exp[i2πλ−1|Δn|2δd2+zp2].(2.14)

Eq. (2.14) can be rewritten as


H(m,n)=∑pApOH(|Δm|;zp)OV(|Δn|;zp),(2.15)

where OH(|Δm|;zp)=exp[i2πλ−1|Δm|2δd2+zp2] is the horizontal light modulation factor and OV(|Δn|;zp)=exp[i2πλ−1|Δn|2δd2+zp2] is the vertical light modulation factor. Both factors can be precomputed in advance and stored in the computer memory. In the generation of the hologram, the pair of functions are fetched from the memory based on values of |Δm| , |Δn| , and zp . The remaining calculations will be the product of intensity Ap with the pair of retrieved functions for each object point, and the summing of the holographic signals for all the object points.


Without loss of generality, the pixel size of the object and hologram spaces are assumed to be identical, and both spaces have horizontal and vertical extents of M rows and N columns. The size of the LUTs for a single depth plane are M for OH(|Δm|;zp) and N for OV(|Δn|;zp) (where 0≤|Δm|<M and 0≤|Δn|<N ). However, since the computation of both light modulation factors is the same, only a single LUT is needed. Assuming a square hologram ( M=N ), this results in a LUT with M entries, which is much less than that required in the early LUT method. The size of the S-LUT for multiple plane will be multiplied by the number of depth planes that are used to discretize the 3-D object along the axial direction. The memory size of the S-LUT for each depth plane of a 1024×1024 hologram, with each entry of the S-LUT represented with 16 bits, is 1024×16 ≈ 16 kbits.


Generation of a hologram pixel at location (m,n) for a single object point with the S-LUT method is summarized in Figure 2.5. Given an object point at position (mp,np;zp) , the horizontal and the vertical light modulation factors are retrieved from the LUTs. The pair of light modulation factors is then multiplied with the intensity Ap of the object point to generate the hologram pixel.


Mar 16, 2021 | Posted by in Electrical Engineer | Comments Off on 2 – Fast Methods for Computer-Generated Holography
Premium Wordpress Themes by UFO Themes