We propose a fast calculation method for diffraction to nonplanar surfaces using the fast-Fourier transform (FFT) algorithm. In this method, the diffracted wavefront on a cylindrical surface is expressed as a convolution between the point response function and the spatial distribution of objects wherein the convolution is calculated using FFT. The principle of the fast calculation and the simulation results are presented.
©2005 Optical Society of America
Computer-generated holography is one of the fields that benefit from the rapid development of computers . Although an enormous amount of time and memory are required for calculating a computer-generated hologram (CGH), present computers enable its synthesis within reasonable time. In addition to computer performance, an algorithm that calculates diffraction should be significant in terms of reducing the calculation time. The fast-Fourier transform (FFT) algorithm enables this reduction. Based on FFT, various types of algorithms for calculating diffraction have been developed . Yoshikawa et al. have proposed a fast calculation method for large size holograms by interpolation .
Most algorithms using FFT are effective only under the condition that both the input and observation surfaces are finite planes that are parallel to each other. Some authors proposed fast calculation methods using FFT that can be applied to the case where the input plane is not parallel to the observation plane [4, 5]. These methods are very useful for the calculation of reconstructed images observed from different points of view [6, 7]. However, since the observation surfaces, in any of the methods, are assumed to be planes, very high-resolution display devices are necessary to enlarge the viewing angles. Moreover, although such devices are available, the reconstructed images cannot be observed from the opposite side of a hologram.
A remarkable technique developed to achieve 360° field of view is 360° holography . If a computer-generated 360° hologram is synthesized in a computer, a numerical simulation of the diffraction on the nonplanar observation surfaces is required. Rosen synthesized a CGH on a spherical observation surface ; however, this method does not yield a 360° field of view because the origin of the object does not correspond to the center of the observation sphere, and considerable computing time is required since the FFT algorithm cannot be applied in this method.
To the best of our knowledge, to data, there exist no reports on the methods that enable the fast calculation of diffraction on spherical or cylindrical observation surfaces. In this paper, we propose such a fast calculation method based on the convolution theorem using the FFT algorithm. In this method, both the object and observation surfaces are cylindrical and concentric about the cylindrical axis. Therefore, this method yields a 360° field of view.
2. Correlation between the object and observation surfaces
To achieve fast calculation, we employ the FFT algorithm to calculate the convolution. Thus, the diffracted wavefront on the observation surface must be expressed in the form of a convolution with the distribution on the object surface. In case of a computer-generated 360° hologram, the shape of the observation surface is usually cylindrical. Therefore, the shape of the object surface must be also cylindrical to form a convolution between them, and both the surfaces must be concentric. The schematic of the geometrical relation between them is shown in Fig. 1, where r and R denote the radii of the object and observation surfaces, respectively. The distributions on the object and observation surfaces are represented by fr(θ,y) and gR(ϕ,y 0), respectively, as shown in Fig. 1. Based on the assumption that the object is the aggregate of many point light sources, the distribution gR(ϕ,y 0) is the summation of the spherical wavefronts that emerge from the point light sources, and it is given by
where k and C denote the wavenumber of the incident light and a constant, respectively, and L represents the distance between two points. Here, the point response function (PRF) is defined as
Then, Eq. (1) is expressed in the form of a two-dimensional convolution as follows:
where * denotes the convolution integral. Therefore, fast acquisition of the diffracted wavefront on the cylindrical observation surface becomes possible by calculating Eq. (4) with an FFT based on the convolution theorem.
However, in case of diffraction for a 3-D object, the 3-D object is defined in cylindrical co-ordinates. Therefore, the radius of the object surface r must be treated not as a constant but a variable, and the PRF expressed in Eq. (3) becomes dependent on r. As a result, the diffracted wavefront on the observation surface is given as the superposition of all the wavefronts diffracted from every cylindrical object surface.
3. Spectral band widths
It is impossible to analytically Fourier transform h(θ,y) in Eq. (3). Therefore, it is necessary to numerically Fourier transform h(θ,y). In order to achieve this, the sampling theorem should be satisfied in both the azimuthal, θ, and the vertical, y, directions during the discretizing process. Therefore, the spectral band width of each direction is discussed here.
3.1. The azimuthal direction
The local spatial frequency νθ (θ,y) in the azimuthal direction at the local point (θ,y) is calculated as follows:
where it is assumed that the variation in (R 2+r 2-2Rr cosθ+y 2)-1/2 is much less than that in exp[ik(R 2+r 2-2Rr cosθ+y 2)]. After simple calculations, the maximum of νθ (θ,y) and the spectral band widthWθ are given by
where ℒ represents the smaller length of r and R.
3.2. The vertical direction
Based on the same assumption as in the case of the azimuthal direction, the local spatial frequency νy(θ,y) in the vertical direction is calculated as follows:
Similarly the maximum of νy(θ,y) and the spectral band width Wy are given by
where Δy represents the height of the object surface, which is the same as that of the observation surface, as shown in Fig. 1.
3.3. Aliasing in the azimuthal direction
fr(θ,y) and h(θ,y) are periodic functions with period 2π in the azimuthal direction. Thus, the analytical convolution between them is also a periodic function with the same period. On the other hand, in case of the numerical convolution between them based on the FFT, aliasing caused by the digitization is significant. However, fortunately, aliasing error can be avoided, but only in the azimuthal direction. Aliasing arises from the superposition of identical signals with different parallel shifts, where the amount of the parallel shift is an integer multiple of the azimuthal size of fr(θ,y) and h(θ,y). Hence, it is possible to match the amount of the parallel shift with integer multiples of the period of fr and h, namely 2π. Therefore, the aliasing error is completely avoided in the azimuthal direction.
4. Numerical experiments
4.1. The spectral band widths of PRF
Figure 2 shows the amplitude distribution of the Fourier spectrum of PRF calculated numerically from Eq. (3). The values of r and R are fixed to 1 cm and 10 cm, respectively and that of Δy is 10 cm. In the range of optical wavelengths, very large pixel sizes are required to satisfy the sampling theorem under these conditions. Hence, the wavelength is set to 300µm, which is in the tera-hertz region. Then, the pixel number in both the azimuthal and vertical directions is 512. Thus, the angular sampling pitch of both θ and ϕ is given by 2π/512. If a CGH were designed with a wavelength of 632.8 nm, the number of pixels required would be N ≅ 2.0×105. On the other hand, according to Eqs. (6) and (8), the spectral band widths in the azimuthal and vertical directions are estimated to be 66.7 rad-1 and 3238 m-1, respectively. These numerical values correspond to the spectral widths shown in Fig. 2.
4.2. Cylindrical observation of Young’s fringe
We have simulated the simplest example of our method, referred to as Young’s experiment, where only two point light sources are distributed on the object surface. The two point sources are located at the same position in the vertical direction, while they are located at θ=-π/32 and π/32, respectively, in the azimuthal direction. The values of r, R, Δy, and the wavelength are the same as those in section 4.1. The interference patterns on the cylindrical observation surface are shown in Fig. 3. Here, the pitches of the fringes vary with the azimuth angles. The calculation time required in our method is 0.9 s. On the other hand, when the same simulation is demonstrated by the direct method used to directly calculate Eq. (1), the calculation time exceeds 9×104 s. Therefore, our method reduces the calculation time by 105 times that achieved by the direct method, thus enabling faster calculation of the diffraction on the cylindrical surface.
4.3. Reconstructed images from the cylindrical CGH
Other examples of our method are presented here. The object used for this simulation is shown in Fig. 4. The parameters used for this simulation are the same as those used for the previous example. According to our method, a CGH consisting of complex values on the cylindrical surface is synthesized from Fig. 4. Then, the reconstructed images on the sectional planes are calculated from the CGH by the direct method. The schematic for the reconstruction is shown in Fig. 5. The position of the light source is virtually set to the origin to produce a spherical wave, which is uniformly illuminated on the cylindrical CGH. Figure 6 shows the reconstructed images in the sectional planes of (a) z=1 cm and (b) -1 cm. Since the radius of the cylindrical object surface is r=1 cm, the character located in the vicinity of θ=0 in Fig. 4, namely “x”, is in focus in the sectional plane of z=1 cm, as shown in Fig. 6(a); further, the character “O”, which is located in the vicinity of θ=π in Fig. 4, is in focus in the sectional plane of z=-1 cm, as shown in Fig. 6(b). Therefore, these results verify our method.
A method for the fast calculation of diffraction on cylindrical surfaces has been proposed in this paper. The algorithm of our method and the analysis of the spectral band widths are discussed. The simulation of Young’s fringe is demonstrated as the simplest example of our method. The reconstructed images in the sectional planes from the cylindrical CGH are also presented. The calculation time required in our method is 105 times less than that of the direct method. Therefore, our method is extremely useful in creating computer-generated 360° holograms.
References and links
2. D. Leseberg, “Sizable Fresnel-type hologram generated by computer,” J. Opt. Soc. Am. A 6, 229–233 (1989). [CrossRef]
3. N. Yoshikawa, M. Itoh, and T. Yatagai, “Interpolation of reconstructed image in Fourier transform computer-generated hologram,” Opt. Commun. 119, 33–40 (1995). [CrossRef]
5. T. Tommasi and B. Bianco, “Computer-generated holograms of tilted planes by a spatial frequency approach,” J. Opt. Soc. Am. A 10, 299–305 (1993). [CrossRef]
6. Y. Takaki and H. Ohzu, “Hybrid Holographic Microscopy: Visualization of Three-Dimensional Object Information by use of Viewing Angles,” Appl. Opt. 39, 5302–5308 (2000). [CrossRef]
7. L. Yu, Y. An, and L. Cai, “Numerical reconstruction of digital holograms with variable viewing angles,” Opt. Express 10, 1250–1257 (2002), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-10-22-1250. [PubMed]
9. J. Rosen, “Computer-generated holograms of images reconstructed on curved surfaces,” Appl. Opt. 38, 6136–6140 (1999). [CrossRef]