Performance Estimation of Intensity Accumulation Display by Computer-Generated Holograms

: The image generated by binary computer-generated holograms (CGHs) always suffers from serious speckle noise. Thanks to the fast frame rate of the binary spatial light modulator, the speckle can be signiﬁcantly suppressed by intensity accumulation, i.e., the sequential display of multiple CGHs of the same scene. If enough randomness is added to the CGHs, the speckle noise can be mostly averaged out. Intuitively, the quality of the reconstructed image should be proportional to the number of intensity accumulation. However, there is no simple method to predict the dependence of the average noise and accumulation number, and we can only know the results after ﬁnishing the full computation. In this paper, we propose an empirical formula of the average noise based on the speckle phenomenon in a laser projector. Using this model, we have conﬁrmed that the randomness induced by random phase is equivalent to that induced by random down-sampling for the generation of binary CGHs. In addition, if the computational efﬁciency is a concern, the CGH calculated with iterations is not recommended for intensity accumulation display. Finally, there is an upper-quality limit of the reconstructed image by intensity accumulation. Thus, a strategy for efﬁcient intensity accumulation is suggested.


Introduction
The use of a computer-generated hologram (CGH) is a promising technique for naked-eye three-dimensional (3D) display without the issue of vergence-accommodation conflict [1]. To realize the display, the CGHs are applied to address a spatial light modulator (SLM), generating the light field of a 3D scene. Usually, the SLM is a phase-modulation (PM) SLM for the high diffraction efficiency. Alternatively, the 3D display can also be realized by using a binary SLM; i.e., the modulation is 0, π in phase or 0, 1 in amplitude [2][3][4][5][6]. The display quality of a binary SLM is usually low because of the limited degrees of freedom. Although optimization algorithms, such as binary search [7][8][9] or simulated annealing [10], can be applied to improve the quality of binary CGH display, the quality of a single binary CGH is still limited. Error diffusion [11][12][13] can generate a binary CGH with good quality only when the object distance, the normalization factor, and the weighting factors are optimized. Super-pixel display schemes [14][15][16] can achieve high-quality display by binary CGHs. However, its optical filtering setup significantly increases the form factor, thereby limiting its application potential. On the other hand, it is easier to drive a binary SLM, and thus the display frame rate of a binary SLM can be ten to a thousand times faster than an eight-bit PM SLM. Therefore, it is capable of quickly displaying different CGHs of the same scene, which is called intensity accumulation (IA) [17][18][19][20]. In this way, the speckle noise can be averaged out and thus the quality of a reconstructed image can be significantly improved. However, there is no simple method to predict the dependence of average noise and accumulation number, and we can only know the results after finishing the full computation. In this paper, our goal is to find a suitable model to describe the dependence between the average noise and the accumulation number. This model is confirmed to fit the properties of binary CGHs generated by the random phase (RP) method, the modified iterative Fresnel algorithm (MIFA), and the localized random down-sampling (LRDS) method. With the help of the noise model, we can optimize the design strategy for various application circumstances. Figure 1 illustrates the setup for calculating binary CGHs. First, the complex amplitude of the object light at the hologram plane is calculated by

Binary Computer-Generated Holograms
where O(x, y) is the object function (the amplitude transmittance of object), which is assumed to be two-dimensional (2D) for simplicity; T( f x , f y ; z) is the free-space transfer function, ( f x , f y ) are spatial frequencies corresponding to spatial coordinates x and y; z 0 is the distance from the object to the hologram plane; and finally F{·} and F −1 {·} denote the 2D Fourier transform and 2D inverse Fourier transform, respectively. In the digital calculation of free-space propagation, the workspace of both the object plane and the hologram plane must be much larger than the object; otherwise, aliasing occurs [21]. Subsequently, the light field inside the window region of the hologram plane is cropped and binarized by where w(x, y) is the window function, and B{·} is the binarization operator, defined by [17][18][19][20]. In this way, the speckle noise can be averaged out and thus the quality of a reconstructed image can be significantly improved. However, there is no simple method to predict the dependence of average noise and accumulation number, and we can only know the results after finishing the full computation. In this paper, our goal is to find a suitable model to describe the dependence between the average noise and the accumulation number. This model is confirmed to fit the properties of binary CGHs generated by the random phase (RP) method, the modified iterative Fresnel algorithm (MIFA), and the localized random down-sampling (LRDS) method. With the help of the noise model, we can optimize the design strategy for various application circumstances. Figure 1 illustrates the setup for calculating binary CGHs. First, the complex amplitude of the object light at the hologram plane is calculated by

Binary Computer-Generated Holograms
where ( , ) O x y is the object function (the amplitude transmittance of object), which is assumed to be two-dimensional (2D) for simplicity; denote the 2D Fourier transform and 2D inverse Fourier transform, respectively. In the digital calculation of free-space propagation, the workspace of both the object plane and the hologram plane must be much larger than the object; otherwise, aliasing occurs [21]. Subsequently, the light field inside the window region of the hologram plane is cropped and binarized by   In reconstruction, the binary hologram is backpropagated from the hologram plane to the object plane. There are two problems in the calculation of binary CGHs. First, if the object function is spatially smooth, the binary CGH will reconstruct a high-pass filtered version of the image and the low spatial-frequency information is lost [22]. Therefore, the object function should be pre-processed before the forward propagation. The processing methods will be introduced in the next section. The second problem is that both the field Appl. Sci. 2021, 11, 7729 3 of 10 truncation due to the window function and binarization will result in speckle noise. In this paper, our focus is on speckle noise suppression by intensity accumulation. The image reconstructed by intensity accumulation is expressed as where E n r (x, y) is the reconstructed field of the n-th binary hologram H n b (x, y), and N is the total number of holograms for the same object. Finally, we use peak-signal-to-noise ratio (PSNR) as the metric to estimate the quality of the reconstructed image. PSNR is defined as where I MAX and I RMSE are, respectively, the maximum signal and the root mean square error of I r (x, y) in comparison with |O(x, y)| 2 .

Processing of the Object Function
The typical method to process the object function is the random phase (RP) method, in which a random phase function is multiplied to an object function, giving a new function, where φ n RP (x, y) is a random value in the range (−π, π). Because there will be numerous distributions of the random phase, it is easy to generate N different random phase functions for N binary CGHs. A variation of the RP method is the modified iterative Fresnel algorithm (MIFA), in which the object field and the binary hologram are iteratively obtained with each other [23,24]. The object-amplitude constraint and the binary-amplitude constraint are applied on the retrieved object and the hologram, respectively. Therefore, each iteration in MIFA needs one forward propagation and one backward propagation. After multiple iterations, the object can be regarded as being modulated by an optimized phase mask, and the phase distribution of the object still looks random.
The second method of object processing is localized random down-sampling (LRDS). The down-sampled object function is expressed as where M m LRDS (x, y) is a binary mask with the same size as the object function. The mask M m LRDS (x, y) is segmented into many tiles with size k × k pixels. For each tile, one pixel is randomly selected as unity (on), and the others are zero (off). Therefore, the down-sampling rate of the mask is 1/k 2 , and many LRDS masks with different sampling positions can be generated. Based on LRDS, there is an advanced method called adaptive intensity accumulation (AIA) [18]. In the AIA method, the object will not only be down-sampled but will also be modified by an adaptive binary mask. The n-th adaptive mask is determined from the reconstruction of previous n CGHs to enhance the image quality. Therefore, every CGH for the AIA method needs one forward propagation and one backward propagation.

Noise Model of Intensity Accumulation
Before comparing different methods of binary CGHs, we need to find a model to quantitatively describe the dependence of the speckle noise on the number of intensity accumulation. Because this is a statistical problem and needs complicated calculation, we intend not to derive the exact solution. By contrast, we only build an empirical formula according to the simulation data. The candidate formula comes from an analogous speckle phenomenon, i.e., the summation of multiple speckled images in laser projection [25,26]. This postulation is reasonable because the holographic recording is similar to the step of light from a micro displayer to the entrance pupil of the projector; the holographic reconstruction is comparable to the procedure of light from the exit pupil of the projector to the screen. In short, the speckle contrast of the displayed image is expressed as where N is the number of summation, and K is the ratio of imaging resolution (the minimum imaging spot size) over the speckle size. In our problem, i.e., the intensity accumulation of CGHs, K depends on object distance, the pixel pitch, and the hologram size. According to our simulation condition, K ranges from 15 to 25 normally. However, direct calculation of K is difficult because it is also affected by the distribution of object function and the processing method. Finally, by assuming that all the noise in reconstruction comes from speckle, the PSNR in this model is expressed as where S is the peak signal, which depends on the distribution of the objective function and the object processing method. It should be noted that Equation (8) is originally derived under the assumption that the signal is uniform. In practice, the signal, i.e., the original object, contains different brightness in various locations. Therefore, we can only find the mean values of S and K from simulation or experimental data.

Simulation
We performed simulations to confirm the noise model. The simulation specifications are as follows. The wavelength and the pixel pitch are 0.65 µm and 7.56 µm, respectively. The computing workspace is 2000 × 2000 pixels, while both the object size and the hologram size are 512 × 512 pixels. The object distance is 190 mm. First, we used the RP method to generate 100 binary CGHs and calculated the PSNR of the reconstructed image with different accumulation numbers. The ground truth image (object) applied in the simulation is shown in Figure 2a. Figure 2b shows a typical reconstructed image, which is from the intensity accumulation of images reconstructed by 100 RP binary CGHs. The resolution of the reconstructed image is apparently less than the ground truth image because of the natural lowpass filtering of the limited size of the hologram (Figure 1). The PSNR of the reconstructed image is calculated in the region of interest (the same size as the object) at every number of intensity accumulation. As we obtain the PSNR curve from simulations, we can fit the simulation curve using Equations (8) and (9) to find the values of S and K. Both the PSNR from simulations (in circle) and the fitting PSNR C (in solid line) are shown in Figure 3. Apparently, the proposed noise formula matches the simulation very well. In this case, the values of the parameters (S, K) are found to be 2.04 and 17.36, respectively. The values can be also regarded as a reference of comparison. The same procedures for MIFA CGHs (five iterations and ten iterations), the LRDS CGHS, and the LRDS+AIA CGHs were repeated. The PSNR data, together with the fitting curves, are shown in Figures 4-7, respectively. The model fit well most types of CGHs, except for the AIA method. These simulation results are discussed in the next section.

Discussions
The important parameters obtained from the empirical formula are listed in Table 1. In particular, we consider the PSNR at different propagations (P) or different numbers of intensity accumulation (N). The concern of different propagations is for the evaluation of the computing-capacity-limited case, in which the binary CGHs must be generated using a limited computing power (e.g., the real-time display). The concern of different numbers of intensity accumulation is for the evaluation of the frame-rate-limited case, in which the frame rate is limited but the computing power is not a constraint (i.e., non-real-time applications). From Equations (8) and (9) and for N → ∞ , we can also define a maximum achievable PSNR as 20 log 10 S √ K . In other words, the image quality cannot be infinitely improved via intensity accumulation. This is very different from previous studies [19,20,27], in which the noise was assumed to be proportional to 1/ √ N and thus the quality could be continuously improved. For the MIFA method (Figure 4), the S value raises to above 3, but K decreases, in comparison with the RP method. For MIFA with 10 iterations (Figure 5), the S value increases further than that of MIFA, with five iterations. It is postulated that a greater iteration number can increase the maximum achievable PSNR. On the contrary, the maximum achievable PSNR decreases slightly in the case of more iterations. This counterintuitive result can be explained by the fact that the iteration in producing one hologram is a procedure of signal enhancement. Meanwhile, the same procedure also increases the correlation between different holograms, which effectively reduces the value of K in the average process. The significantly reduced K also results in another problem. The efficiency (PSNR at 100 propagations) of intensity accumulation by MIFA CGHs is the worst among all discussed methods. Accordingly, MIFA is not recommended for the intensity accumulation, provided the computing cost is the main concern. The characteristics of LRDS CGH ( Figure 6) are almost the same as those of RP CGH. This implies that the randomness generated by random phase and by random down-sampling can be regarded to be identical in practice. Therefore, for a pointwise CGH algorithm (e.g., [28]), the use of random down-sampling can significantly reduce the computing load without degrading the display quality. The fitting model is not applicable for the LRDS+AIA method (Figure 7). The reason is that AIA is an adaptive method to improve the quality of reconstruction globally. Therefore, the statistical randomness of each binary hologram has been broken in this method. According to the above discussion and Table 1, we can conclude that the RP method and the LRDS method are better for the computing-capacity-limited case, while MIFA method can provide a better PSNR for the frame-rate-limited case. The LRDS + AIA method always exhibits the best PSNR, except for the case of a small number of propagations (P = 10).
Finally, we will consider the design strategy of intensity accumulation. Since there is a maximum achievable PSNR, the display by using many numbers of intensity accumulation is a waste of computing power. Our goal is to find a reasonable number of intensity accumulation to achieve a PSNR only slightly lower than the maximum achievable PSNR. The difference between the maximum achievable PSNR and PSNR at a specific number N ∆ is calculated as ∆ = 20 log 10 S √ K − 20 log 10 It is noted that ∆ is independent of S. From Equation (10) and for a given ∆ and K, we find N ∆ N ∆ = K + 1 10 0.1∆ − 1 .
The needed accumulation numbers for different K's and ∆'s were calculated and are listed in Table 2. Normally, we can set ∆ = 1 dB to achieve a balance between quality and computing cost.

Conclusions
In this paper, we have proposed a model to estimate the dependence of PSNR on the number of intensity accumulation by binary CGHs. The model is confirmed to fit the properties of binary CGHs generated by the random phase (RP) method, the modified iterative Fresnel algorithm (MIFA), and the localized random down-sampling (LRDS) method. In summary, MIFA can achieve good quality but is worse in its computing efficiency. The randomness characteristic of RP CGHs and LRDS CGHs is almost the same. Therefore, for a pointwise CGH algorithm, the LRDS method is a better object processing method because the computing load can be significantly reduced without degrading the display quality. For any method, there is a maximum achievable PSNR for the reconstructed image, and thus the number of accumulations must be accurately determined in order to avoid the waste of computing power. In this paper, we only considered the speckle due to wave truncation and quantization. In practice, the speckle due to scattered light from optical elements cannot be suppressed by intensity accumulation, and thus the image quality is usually worse than the simulation. On the other hand, the coherence of light source is another important factor of image quality. If a low-coherence source is applied, the speckle can be better suppressed. This aspect has not been considered in the present paper, and it is worthy of further study.