Accelerated synthesis algorithm of polygon computer-generated holograms

For the real-time computation of computer-generated holograms (CGHs), various accelerated algorithms have been actively investigated. This paper proposes a novel concept of sparse computation of polygon CGH, which is inspired by an observation of the sparsity in the angular spectrum of a unit triangular polygon and present the accelerated algorithm using the intrinsic sparsity in the polygon CGH pattern for the enhancement of computational efficiency effectively. It is shown with numerical results that computation efficiency can be greatly improved without degrading the quality of holographic image. ©2015 Optical Society of America OCIS codes: (090.2870) Holographic display; (090.1995) Digital holography References and links 1. J. Hong, Y. Kim, H.-J. Choi, J. Hahn, J.-H. Park, H. Kim, S.-W. Min, N. Chen, and B. Lee, “Three-dimensional display technologies of recent interest: principles, status, and issues [Invited],” Appl. Opt. 50(34), H87–H115 (2011). 2. R. Haussler, A. Schwerdtner, and N. Leister, “Large holographic displays as an alternative to stereoscopic displays,” Proc. SPIE 6803, 68030M (2008). 3. R. Stahl, V. Rochus, X. Rottenberg, S. Cosemans, L. Haspeslagh, S. Severi, G. V. Plas, G. Lafruit, and S. Donnay, “Modular sub-wavelength diffractive light modulator for high-definition holographic displays,” J. Phys. Conf. Ser. 415, 012057 (2013). 4. K. Aoshima, N. Funabashi, K. Machida, Y. Miyamoto, K. Kuga, T. Ishibashi, N. Shimidzu, and F. Sato, “Submicron magneto-optical spatial light modulation device for holographic displays driven by spin-polarized electrons,” J. Disp. Tech. 6(9), 374–380 (2010). 5. J. Sun, E. Timurdogan, A. Yaacobi, E. S. Hosseini, and M. R. Watts, “Large-scale nanophotonic phased array,” Nature 493(7431), 195–199 (2013). 6. S. Reichelt, R. Häussler, G. Fütterer, N. Leister, H. Kato, N. Usukura, and Y. Kanbayashi, “Full-range, complex spatial light modulator for real-time holography,” Opt. Lett. 37(11), 1955–1957 (2012). 7. H. Kim, C.-Y. Hwang, K.-S. Kim, J. Roh, W. Moon, S. Kim, B.-R. Lee, S. Oh, and J. Hahn, “Anamorphic optical transformation of an amplitude spatial light modulator to a complex spatial light modulator with square pixels [invited],” Appl. Opt. 53(27), G139–G146 (2014). 8. N. Takada, T. Shimobaba, H. Nakayama, A. Shiraki, N. Okada, M. Oikawa, N. Masuda, and T. Ito, “Fast highresolution computer-generated hologram computation using multiple graphics processing unit cluster system,” Appl. Opt. 51(30), 7303–7307 (2012). 9. H. Nakayama, N. Takada, Y. Ichihashi, S. Awazu, T. Shimobaba, N. Masuda, and T. Ito, “Real-time color electroholography using multiple graphics processing units and multiple high-definition liquid-crystal display panels,” Appl. Opt. 49(31), 5993–5996 (2010). 10. T. Shimobaba, T. Ito, N. Masuda, Y. Ichihashi, and N. Takada, “Fast calculation of computer-generatedhologram on AMD HD5000 series GPU and OpenCL,” Opt. Express 18(10), 9955–9960 (2010). 11. Y.-Z. Liu, J.-W. Dong, Y.-Y. Pu, B.-C. Chen, H.-X. He, and H.-Z. Wang, “High-speed full analytical holographic computations for true-life scenes,” Opt. Express 18(4), 3345–3351 (2010). 12. J. Cho, J. Hahn, and H. Kim, “Fast reconfiguration algorithm of computer generated holograms for adaptive view direction change in holographic three-dimensional display,” Opt. Express 20(27), 28282–28291 (2012). 13. Y. Pan, Y. Wang, J. Liu, X. Li, and J. Jia, “Fast polygon-based method for calculating computer-generated holograms in three-dimensional display,” Appl. Opt. 52(1), A290–A299 (2013). 14. W. Lee, D. Im, J. Paek, J. Hahn, and H. Kim, “Semi-analytic texturing algorithm for polygon computergenerated holograms,” Opt. Express 22(25), 31180–31191 (2014). 15. H. Kim, J. Hahn, and B. Lee, “Mathematical modeling of triangle-mesh-modeled three-dimensional surface objects for digital holography,” Appl. Opt. 47(19), D117–D127 (2008). #228323 $15.00 USD Received 20 Nov 2014; revised 15 Jan 2015; accepted 23 Jan 2015; published 30 Jan 2015 © 2015 OSA 9 Feb 2015 | Vol. 23, No. 3 | DOI:10.1364/OE.23.002863 | OPTICS EXPRESS 2863 16. D. Im, E. Moon, Y. Park, D. Lee, J. Hahn, and H. Kim, “Phase-regularized polygon computer-generated holograms,” Opt. Lett. 39(12), 3642–3645 (2014). 17. D. Leseberg and C. Frère, “Computer-generated holograms of 3-D objects composed of tilted planar segments,” Appl. Opt. 27(14), 3020–3024 (1988). 18. K. Matsushima, “Computer-generated holograms for three-dimensional surface objects with shade and texture,” Appl. Opt. 44(22), 4607–4614 (2005). 19. K. Matsushima and S. Nakahara, “Extremely high-definition full-parallax computer-generated hologram created by the polygon-based method,” Appl. Opt. 48(34), H54–H63 (2009). 20. K. Matsushima, M. Nakamura, and S. Nakahara, “Silhouette method for hidden surface removal in computer holography and its acceleration using the switch-back technique,” Opt. Express 22(20), 24450–24465 (2014).


Introduction
The recent increase in public interest in holographic three-dimensional (3D) displays is attributable to their promising possibility of unprecedented 3D image quality and the delivery of natural 3D cues without the problems commonly associated with commercial stereoscopic 3D displays, such as accommodation-vergence conflict, eye strain, and fatigue [1,2].In principle, the holographic 3D display is therefore the ultimate form of 3D display, but several technological challenges remain to be solved in order to develop commercializable products and achieve commercial success.Significant research efforts have been made to overcome these technological challenges in recent years.One hardware challenge is the development of spatial light modulators (SLMs) with wavelength-scale pixels.Several studies have actively pursued, wavelength-scale SLMs based on novel technologies such as nanoelectromechanical systems (NEMS) [3], spin transfer switching magneto-optical (STS-MO) devices [4] and the active thermo-optic phase tuning [5] in addition to advanced liquid crystal (LC) SLMs [6,7].The development of efficient computational algorithms for the seamless real-time generation of large-scale photorealistic CGHs is a software challenge.In particular, recent research in this area has tended to focus on computational enhancement using parallel computing technology based on specific parallel computing hardware modules such as graphic processing units (GPUs) and parallel computing architectures [8][9][10][11].More fundamental innovations in mathematical modeling for computational efficiency enhancement are, however, relatively rare, although some theoretical studies to improve the computational efficiency of the polygon CGH were reported in some previous papers [12][13][14].
In this paper, we propose a simple and effective method for the enhancement of the computational efficiency of the polygon CGH synthesis algorithm [15][16][17][18][19][20].The idea for this was obtained from the observation that the two-dimensional (2D) angular spectrum profile of a triangular unit, i.e. the elementary unit of the polygon CGH, is inherently sparse.In other words, very large percentage of the angular spectrum domain area is filled with negligibly small values that can be approximated to zero, while the angular spectrum values in only small rest of the area determine the whole angular spectrum of the triangular unit.The efficient use of this sparsity in polygon CGH computation is expected to allow for increased computational efficiency.The proposed algorithm is not a parallel processing algorithm based on conventional mathematical formulas, but rather a novel acceleration algorithm with mathematically improved features that enhance its effectiveness even for single-computingcore platforms.Of course, extending the proposed acceleration algorithm to accelerated parallel processing frameworks is straightforward, so it is expected to lead to unprecedented computational acceleration in the polygon CGH generation.
In this study, the validity and effectiveness of the proposed algorithm will be tested with numerical simulations.In Section 2, the proposed algorithm is described and the analysis of the computation time enhancement is presented.In Section 3, the numerical reconstruction of the polygon CGH of a complex 3D object is presented to validate our approach and the possible influence on image reconstruction quality is analyzed, followed by concluding remarks.

Accelerated synthesis algorithm for polygon computer-generated holograms
In order to describe the proposed algorithm, a brief review of the numerical framework of polygon CGH synthesis algorithms is needed.In the polygon CGH method, 3D objects are depicted by triangular meshes.The light field distribution of a single triangle facet is an elementary unit of full 3D light field distribution, and the holographic 3D light field is synthesized by the total superposition of the elementary light fields as shown in Fig. 1.Here a triangle-meshed 3D lobster is chosen as an example target 3D object.In Fig. 2(a), a unit triangular facet in a global coordinate system is shown.On the k th triangular facet of the polygon object, its local coordinate system is depicted, where the center of mass and the normal vector of the triangle facet are assumed to be ( ) ( ) cos sin ,sin sin , cos , respectively, and the center of mass is set to the origin of the local coordinate system.The spatial coordinates ( ) , , x y z ′ ′ ′ in the local coordinate system and ( ) , , x y z in the global coordinate system are related by the rotational transform, cos cos cos sin sin sin cos 0 .sin cos sin sin cos In the local coordinate system, the normal vector k n is equal to the z′ axis, so the triangular facet is placed on the x y ′ ′ plane.
where ( ) , α β are the spatial frequency and γ is given by . λ is the wavelength of light in free space.The angular spectrum in the global coordinate system ( ) and k η are the spatial frequency in the local coordinate system and the complex amplitude. ( ) ( ) ( )  The main idea for achieving the computational acceleration is that the angular spectrum profile can be approximated to a sparse pattern by excluding the small angular spectrum values smaller than a threshold set near to zero in the calculation.For this, is set the region of interest (ROI) that indicates the area with angular spectrum values used in the CGH calculation.In the conventional computation schemes of polygon CGH, this sparse property within the angular spectrum has not been perceived.However, because significant proportion of computation time is spent on the calculation of near zero values in the angular spectrum, it is expected that the selective calculation of the angular spectrum in this ROI, where nonnegligible values are compactly distributed, can greatly save computation time.The ROI is essential for accelerating the CGH synthesis algorithm.Our finding is that the ROI can be very efficiently constructed via a simple mathematical formulation.
Let us consider the triangular facet in the local coordinate system in Fig. 3(a) with the three apex points of the triangle 1 2 3 P P P Δ , where ( ) , P x y ′ ′ = and ( ) The perpendicular bisectors of the triangle, 12 t′ , 23 t′ , and 31 t′ are obtained, respectively, by ) where FT means the Fourier transform.In Fig. 5, the reconstructed images of the sparse angular spectrum and the conventional dense angular spectrum are compared.The sparsity, i.e., the ratio of the ROI area to the whole area of the angular spectrum domain, is approximately 5.69% for the sparse angular spectrum.Figures 5(c) and 5(d) show the reconstruction results of the tilted triangular facet image from which no perceptible discrimination is observed.Consequently, the sparse angular spectrum obtained by a simple cut-off masking process where only a small percentage in the entire angular spectrum area is filled with non-zero data and the rest of the area that is approximated to be zero does not greatly influence the reconstructed image but can significantly reduce computation.

Numerical reconstruction of polygon CGHs
To verify the proposed method, the numerical simulation of CGH synthesis and reconstruction is performed.In order to clear and accurate holographic 3D images, the polygon CGHs have been calculated using the phase regularization technique [16].The simulation setup is illustrated with parameter specifications in Fig. 6.The optical system schematic of the eye observation of a holographic 3D image is shown.The CGH plane ( ) , the eye lens plane ( ) , u v , and the retina plane ( ) Figs. 7(a), 7(b), and 7(c), the ROI is set to 1.2%, 4.0% and 100% of total area, respectively.For the ROI of 1.2%, only a minor difference in the object surface is discernible from the fullcalculation result.After expanding the ROI to that of 4.0%, it is seen that there is very little difference between the images displayed in Figs.7(b) and 7(c).On average, it took 0.0569 seconds and 0.2052 seconds per triangle facet to calculate the ROIs of 1.2% and 4.0%, respectively, while, the conventional full-calculation required 5.1562 seconds.Thus, the computation times for the first and second cases are approximately 91 times and 25 times faster than conventional calculation, without a significant image quality degradation.This comparison of computation times supports our proposal.The lobster object used in the simulation example is a very general shape object.From this, we expect that even 2% ROI will produce quite good holographic images for other 3D objects.In Fig. 8   The ROI build-up process is analyzed in four steps.As explained in Section 2, the first step calculates the bisectors of the triangular facet in the local coordinate system and maps the bisectors to the ROI indicators in the local angular spectrum domain.In the second step, the ROI indicators are translated from the local to the global coordinate system.The third step is to construct the ROI so that it maintain a constant thickness using the convolution method.In the fourth step, the angular spectrum is calculated in the constructed ROI region.As shown in Fig. 8(b), the first, second and fourth steps are not time-consuming.The time-consuming step is the third one; it is approximately ten times slower than the other steps because of the convolution operation.Since the time required in the third step is determined by the resolution, the computational efficiency cannot improve on this.

Conclusion
In conclusion, we have proposed an accelerated polygon CGH algorithm using the inherent sparisity in the polygon CGH synthesis that allows the effective size reduction of the computation domain.The proposed method was verified through the numerical simulation of the target 3D object showing approximately 90-times higher computation efficiency enhancement without a loss of quality in the reconstructed image.We believe that the performance of the polygon CGH synthesis algorithm will be further improved if the proposed algorithm is combined within parallel processing.

Fig. 1 .
Fig. 1. Construction of a holographic 3D image and its CGH pattern by the superposition of (a) the light field elements of unit triangular facets and (b) the elementary sparse angular spectrums of these facets

Fig. 2 .
Fig. 2. Elementary units of the polygon CGH synthesis algorithm: (a) unit triangular facet and (b) its angular spectrum profile As detailed in ref [15], the light field distribution in the global coordinate system, ( ) , , k W x y z , is described by the angular spectrum representation, ⋅ are the local coordinate angular spectrum of the k th triangular facet and the Heaviside unit step function [15].The spatial frequency variables ( ) , α β ′ ′ in the local coordinate system and ( ) , α β in the global coordinate system are related by the same rotation transform as in Eq. (4).The spatial frequency ( ) , α β in the global coordinate system is given by the matrix equation #228323 -$15.00USD

Figure 2 (
Figure 2(b) shows the angular spectrum profile in the global spatial frequency domain ( ) ,

Fig. 3 .
Fig. 3. (a) A unit triangle facet in a local coordinate system and its angular spectrum profile in (b) the local coordinate system and (c) in the global coordinate system Here, a finding is that the bisectors of the triangle are exactly matched to the angular spectrum profile with the center ( ) 0 0 , α β ′ ′ as shown in Fig. 3(b).The line indicators of the ROI in the spatial frequency domain ( - α β

,
, displayed in Fig. 4(a) with the circular function ( ) Circ α β , seen in Fig. 4(b), in the spatial frequency domain.This is represented as Figure 4(c) shows the convolution result of the indicator and the circular function with coarse sampling.The area of the ROI can be adjusted with the radius of the circular function.

. 6 .Fig. 5 .Fig. 6 .
Fig. 5. (a) Sparse and (b) dense angular spectrums and their reconstruction images of (c) the sparse and (d) dense angular spectrums Figs. 7(a), 7(b), and 7(c), the ROI is set to 1.2%, 4.0% and 100% of total area, respectively.For the ROI of 1.2%, only a minor difference in the object surface is discernible from the fullcalculation result.After expanding the ROI to that of 4.0%, it is seen that there is very little difference between the images displayed in Figs.7(b) and 7(c).On average, it took 0.0569 seconds and 0.2052 seconds per triangle facet to calculate the ROIs of 1.2% and 4.0%, respectively, while, the conventional full-calculation required 5.1562 seconds.Thus, the computation times for the first and second cases are approximately 91 times and 25 times faster than conventional calculation, without a significant image quality degradation.This comparison of computation times supports our proposal.The lobster object used in the simulation example is a very general shape object.From this, we expect that even 2% ROI will produce quite good holographic images for other 3D objects.In Fig.8(a), the total computation times according to the sparsity of the ROI region is presented.The red line represents the total computation time for the full area computation, set as the upper limit for the computation time.The blue line in this figure indicates the computation time of the sparse polygon CGH, which shows that the computation time increases proportionally with the ROI area.This linear relationship between the computation region and the computations speed is as expected.The computation time of the blue line includes both the ROI identification process and the angular spectrum computation, which was performed by the convolution process of the sampled ROI indicator and the circular function with a sizeable radius.The black line is the estimated computing time after the ROI identification time has been removed.In Fig.8(b), we analyze in more detail the cause of the time consumption in the ROI build-up stage.

Fig. 7 .
Fig. 7.The reconstructed images and CGH patterns for the ROIs with various sparsity of (a) 1.2014%, (b) 4.0222% and (c) 100%.The average computation times of a unit triangle in the CGH patterns are estimated to (a) 0.0569, (b) 0.2052 and (c) 5.1562 (seconds), respectively, for those cases.

Fig. 8 .
Fig. 8. (a) Computation time versus the ROI computation region (%) and (b) the computing time used up in the ROI build-up stage.
are the traces of the central point of the ROI, we need to further extend the effective ROI area to cover a range of area around the curved indicator.The necessary ROI construction process, in order to construct an ROI with a finite area, is performed by the convolution of the sampled ROI indicator l