Calculating the Fresnel diffraction of light from a shifted and tilted plane

Kenji Yamamoto; Yasuyuki Ichihashi; Takanori Senoh; Ryutaro Oi; Taiichiro Kurita

doi:10.1364/OE.20.012949

1. Introduction

Since holography is a technique that enables light to be recorded and played back with the precision of lightwave, it is expected to be used as a high-precision 3D imaging technique. As a result, various kinds of research has been conducted including research on sensitive materials for still-image holograms, research concerning dynamic image display, research for obtaining holograms using dynamic image sensors, and research in converting computer graphics data to holograms (computer generated holograms). Research concerning computer generated holograms conventionally had represented the subject by a set of point light sources, but recently, there has been an increase in research that represents the subject by a set of planes in a similar manner as is used in computer graphics. In this type of research, the wavefield on the plane eventually ends up as a collection of point light sources. However, stipulating that there are point light sources on a grid on the plane is different than research that deals with the subject using conventional point light sources. Stipulating this makes methods available that enable the diffraction to be calculated quickly and allow radiation characteristics to be assigned to the light that is emitted from the plane. In addition, the researches that do not split planes into point light sources have emerged [1].

Since calculating the diffraction of light from one plane to another is a basic technique when the subject is represented by a set of planes, this has been investigated many times in the past. For example, Tommasi, Matsushima and Nicola proposed methods of calculating the diffraction from tilted planes by taking the angular spectrum into consideration in the References [2–5]. Delen proposed two methods, one for shifted planes and the other for tilted planes, and these can be combined sequentially for a shifted and tilted plane [6]. In addition, his approach can be applied to wide angle diffraction since it is based on Rayleigh-Sommerfeld diffraction. However, these sorts of methods require the angular spectrum to be calculated in addition to the diffraction of light. Muffoletto proposed a method of calculating the diffraction from a shifted plane, but this method cannot handle tilted planes [7, 8]. Yu proposed a method of calculating the diffraction from a tilted plane, but the rotation was only around one axis and the diffraction could not be calculated from a shifted plane [9]. Miura proposed a method of calculating the diffraction from a tilted and shifted plane, but the rotation was only around one axis and the error increased for planes tilted at large angles [10]. Using the concept of impulse functions over a surface, Onural systematically analyzes what we are now researching in these topics [11].

In this paper, we describe a single technique for calculating the diffraction from a tilted and shifted plane. This technique deals with rotation around three axes. In addition, angular spectrum calculations are unnecessary. This paper is organized as follows. Chapter 2 describes the theoretical basis of the proposed technique, Chapter 3 describes the numerical experiment and shows the calculation times, and Chapter 4 presents conclusions.

2. Theory

Figure 1 shows the setup of the two planes that are dealt with in this paper. The x-, y-, and zaxes are in the world coordinate system, and coordinates are represented by [x,y,z]^T. Note that [ · ]^T indicates the transpose vector. All coordinates are assumed to be represented in the world coordinate system below. The center P₀ of the source plane is [x₀, y₀, z₀]^T. The vectors $\vec{Δ s}$ and $\vec{Δ t}$ are defined to represent the position P_st on this plane. $\vec{Δ s}$ and $\vec{Δ t}$ , which are represented as $\vec{Δ s} = {[Δ s_{x}, Δ s_{y}, Δ s_{z}]}^{T}$ and $\vec{Δ t} = {[Δ t_{x}, Δ t_{y}, Δ t_{z}]}^{T}$ , are assumed to be orthogonal. In other words, Δs_xΔt_x + Δs_yΔt_y + Δs_zΔt_z = 0. The center P₁ of the plane where the diffracted light is to be obtained (destination plane) is [x₁, y₁, z₁]^T. The vectors $\vec{Δ u}$ and $\vec{Δ v}$ , which are parallel to the x- and y-axes respectively, are defined to represent the position P_uv on this plane. $\vec{Δ u} = {[Δ u, 0, 0]}^{T}$ and $\vec{Δ v} = {[0, Δ v, 0]}^{T}$ . According to these definitions, the rotation of the source plane is represented by $\vec{Δ s}$ and $\vec{Δ t}$ , and the translation is represented by the difference between P₁ and P₀.

Fig. 1 Setup of source plane and destination plane

Download Full Size | PDF

The coordinates [x_st, y_st, z_st]^T of P_st and the coordinates [x_uv, y_uv, z_uv]^T of P_uv are defined as follows.

\begin{array}{l} [\begin{matrix} x_{s t} \\ y_{s t} \\ z_{s t} \end{matrix}] & = & P_{0} + s \vec{Δ s} + t \vec{Δ t} \\ = & [\begin{matrix} x_{0} \\ y_{0} \\ z_{0} \end{matrix}] + s [\begin{matrix} Δ s_{x} \\ Δ s_{y} \\ Δ s_{z} \end{matrix}] + t [\begin{matrix} Δ t_{x} \\ Δ t_{y} \\ Δ t_{z} \end{matrix}] \end{array}

\begin{array}{l} [\begin{matrix} x_{u v} \\ y_{u v} \\ z_{u v} \end{matrix}] & = & P_{1} + u \vec{Δ u} + v \vec{Δ v} \\ = & [\begin{matrix} x_{1} \\ y_{1} \\ z_{1} \end{matrix}] + u [\begin{matrix} Δ u \\ 0 \\ 0 \end{matrix}] + v [\begin{matrix} 0 \\ Δ v \\ 0 \end{matrix}] \end{array}

where s and t are real numbers for representing P_st and u and v are real numbers for representing P_uv. However, in the actual calculations, these are all handled as integers to make the calculations discrete.

The source plane wavefield U₀ diffracts to create the wavefield U₁ on the destination plane. The diffracted result U₁ can be represented as follows.

U_{1} (P_{u v}) = \iint \frac{U_{0} (P_{s t})}{j λ | z_{u v} - z_{s t} |} exp {j k \sqrt{{(x_{u v} - x_{s t})}^{2} + {(y_{u v} - y_{s t})}^{2} + {(z_{u v} - z_{s t})}^{2}}} d s d t

where j represents the imaginary unit, λ represents the wavelength of the light, and k represents the wave number (k = 2π/λ). Also, the two planes are assumed to be separated by a distance for which the Fresnel diffraction holds. In other words, Eq. (4) is assumed to hold at all P_st and P_uv [12].

{(z_{u v} - z_{s t})}^{3} ≫ \frac{π}{4 λ} {[{(x_{u v} - x_{s t})}^{2} + {(y_{u v} - y_{s t})}^{2}]}^{2}

The diffraction into the destination plane is also assumed to be limited by the angle for which the Fresnel diffraction holds.

We can approximate Eq. (3) by rewriting it as Eq. (15) if we define the new variables shown in Eqs. (5) to (12), approximate the square root and drop some terms under the case of Eqs. (13) and (14).

z_{01} = z_{1} - z_{0}

Δ z = - s Δ s_{z} - t Δ t_{z}

x_{01} = x_{1} - x_{0}

x_{s t}^{'} = {(s Δ s_{x})}^{2} + {(t Δ t_{x})}^{2} - 2 x_{01} s Δ s_{x} - 2 x_{01} t Δ t_{x} + 2 Δ s_{x} Δ t_{x} s t

x_{u v}^{'} = x_{01}^{2} + {(u Δ u)}^{2} + 2 x_{01} u Δ u

y_{01} = y_{1} - y_{0}

y_{s t}^{'} = {(s Δ s_{y})}^{2} + {(t Δ t_{y})}^{2} - 2 y_{01} s Δ s_{y} - 2 y_{01} t Δ t_{y} + 2 Δ s_{y} Δ t_{y} s t

y_{u v}^{'} = y_{01}^{2} + {(v Δ v)}^{2} + 2 y_{01} v Δ v

{(Δ z)}^{2} ≪ z_{01}^{2}

{(Δ z)}^{2} ≪ z_{01}

\begin{array}{l} U_{1} (P_{u v}) & = & \frac{1}{j λ} exp {j k \frac{x_{u v}^{'} + y_{u v}^{'}}{2 z_{01}}} \\ \iint U_{0} (P_{s t}) \frac{| z_{01} - Δ z |}{z_{01}^{2}} exp {j k (z_{01} + Δ z)} exp {j k \frac{x_{s t}^{'} + y_{s t}^{'}}{2 z_{01}}} \\ exp {- j k \frac{(Δ u Δ s_{x} u + Δ v Δ s_{y} v) s + (Δ u Δ t_{x} u + Δ v Δ t_{y} v) t}{z_{01}}} d s d t \end{array}

It is apparent that the integral part of Eq. (15) contains the scaled Fourier transform shown in Eq. (16) [7, 13]. As a result, the calculation speed can be increased since a scaled Fourier transform can be used in part of the calculation.

F (u) = \int f (s) exp {- j 2 π a u s} d s

where a is the scale parameter.

The s,t,u,v are all handled as integers to make the calculations discrete in actuality. When you discretize the Eq. (15) and expand it, you can find the terms of discretized scaled Fourier transform.

3. Numerical experiment and calculation times

First, we performed two types of numerical experiments using the optical system shown in Fig. 2 and Table 1. In these numerical experiments, the light emitted from the surface plane passes through the lens and forms an image near the image plane. The intensity of the light at the image plane is calculated as the experimental result. If this technique is right, the intensity forms an image which corresponds to the location of surface plane, lens and image plane.

Fig. 2 Optical system used in numerical experiments

Download Full Size | PDF

Table 1. Numerical experiment setup

View Table | View all tables in this article

In experiment 1 and experiment 2 described below, the surface plane, which consists of a set of points on a 256 × 256 grid, had the intensity shown in Fig. 3. The phase at this time was assumed to be random. Also, the lens and image plane consisted of sets of points on a 768×768 grid. Since the numbers of points on the surface plane and lens differed, the calculation of the diffraction from the surface plane to the lens was performed a total of 9 times with a size of 256 × 256 while moving P₁. In other words, the calculation was performed 3 times in the x-axis direction and 3 times in the y-axis direction for a total of 9 times. The diffraction from the lens to the image plane was performed 1 time with a size of 768 × 768. To reduce aliasing, the scaled Fourier transform was calculated with a doubled size in both the horizontal and vertical directions. In other words, the calculation was performed, for example, using a scaled Fourier transform with a size of 512 × 512 for the set of 256 × 256 points. Since the surface plane had random phases, speckle noise was conspicuous at the image plane. Therefore, we performed the experiment 30 times while changing the random phase for each condition and let the accumulated value be the experimental result.

Fig. 3 Intensity distribution of surface plane

Download Full Size | PDF

For experiment 1, we performed a numerical experiment to confirm that the calculations could be performed appropriately with respect to rotation. In this experiment, we fixed the locations of the lens and image plane and changed the inclination of the surface plane by facing it forward or tilting it toward the top or bottom or to the left or right. Figure 4(a) shows the setup of the entire optical system, and Fig. 4(b) shows the inclinations of the surface plane within the optical system and the experimental results. The experimental results are inverted left to right and top to bottom relative to Fig. 3 because of the lens. It is apparent from Fig. 4(b) that appropriate experimental results correspond to the inclinations of the surface plane.

Fig. 4 Experiment 1: experimental results for rotation. (a) Movement of surface plane. (b) Experimental results.

Download Full Size | PDF

For experiment 2, we performed a numerical experiment to confirm that the calculations could be performed appropriately with respect to translation. In this experiment, we fixed the location of the surface plane and shifted the lens and image plane. We performed the experiment when the surface plane was facing forward, when it was rotated around the y-axis, and when it was tilted. Figures 5, 6, and 7 show the various experimental results. In Fig. 5 to 7, (a) shows the setup of the entire optical system, and (b) shows the locations of the lens and image plane within the optical system and the experimental results. The experimental results are inverted left to right and top to bottom relative to Fig. 3 because of the lens. It is apparent from Figs. 5(b), 6(b), and 7(b) that appropriate experimental results correspond to the translations of the locations of the lens and image plane. It is also apparent that the various inclinations of the surface plane are reflected. We prepared a video of this experiment for your reference ( Media 1).

Fig. 5 Experiment 2: experimental results for translation with no rotation ( Media 1). (a) Position of surface plane and movement of image plane. (b) Experimental results.

Download Full Size | PDF

Fig. 6 Experiment 2: experimental results for translation with y-axis positive rotation ( Media 1). (a) Position of surface plane and movement of image plane. (b) Experimental results.

Download Full Size | PDF

Fig. 7 Experiment 2: experimental results of translation with x-axis negative rotation and y-axis positive rotation ( Media 1). (a) Position of surface plane and movement of image plane. (b) Experimental results.

Download Full Size | PDF

It is apparent from the results of experiments 1 and 2 that this technique can be used to calculate the diffraction of light from a tilted and shifted plane.

Second, we compared the time required for this calculation with the time required using a Fresnel transform. Specifically, we measured the calculation times when the setup shown in Fig. 1 was used to sequentially calculate the diffraction from each point of the source plane to each point of the destination plane using a Fresnel transform and when the proposed technique was used to calculate the diffraction from the source plane to the destination plane. The measurements were conducted using six sizes of grids from 32 × 32 (N = 32) to 1024 × 1024 (N = 1024). Here, N(N = 32, 64, 128, 256, 512, 1024) is the number of points along an axis of the grid. Note that we doubled the vertical and horizontal sizes for the scaled Fourier transform as we did in experiments 1 and 2. The PC used for the measurements had an Intel Core i7 CPU with Microsoft Windows 7 64-bit edition operating system running single threaded. Table 2 shows the measurement results. From these measurement results, it is apparent that the proposed technique was able shorten the calculation time t_N to approximately 1/9 (≒ 193/1, 734) that of the Fresnel transform for a 64 × 64 grid and to approximately 1/70 (≒ 99, 918/7, 065, 681) for a 512 × 512 grid.

Table 2. Calculation time t_N

View Table | View all tables in this article

In order to confirm the computational complexity, we defined a_N and b_N as follows.

a_{N} = t_{N} / N^{4}

b_{N} = t_{N} / N^{2} \log_{2} N

If the index log₁₀(a_N/a₃₂) is around 0 for any N, it suggests the computational complexity is O(N⁴). In the same way, if the index log₁₀(b_N/b₃₂) is around 0, it suggests O(N² log₂ N). Figures 8(a) and 8(b) show the indexes at Fresnel transform and at proposed technique respectively. As you can see, the indexes suggest O(N⁴) in (a); however, they do not suggest O(N⁴) in (b).

Fig. 8 Computational complexity. (a) Index at Fresnel transform. (b) Index at proposed technique.

Download Full Size | PDF

The indexes do not suggest O(N² log₂ N) so much in (b). It would be the reason that the Eq. (15) is not a two-dimensional scaled Fourier transform, but we are not be confident about it. We have to find out definite reason in the future.

4. Conclusions

In this paper, we described a technique for calculating the diffraction of light from a tilted and shifted plane. This technique can be used to calculate the diffraction of light from a source plane that is rotated and translated relative to three axes. In addition, it can perform the calculations quickly since a scaled Fourier transform is used for part of the calculation. We confirmed through numerical experiments that the diffraction could be calculated appropriately. We also measured actual calculation times to confirm that the diffraction could be calculated faster than using a Fresnel transform at each point.

In this paper, to provide degrees of freedom for the sampling points on the source plane, we assumed that $\vec{Δ s} = {[Δ s_{x}, Δ s_{y}, Δ s_{z}]}^{T}$ and $\vec{Δ t} = {[Δ t_{x}, Δ t_{y}, Δ t_{z}]}^{T}$ . If this assumption were restricted to $\vec{Δ s} = {[Δ s_{x}, 0, Δ s_{z}]}^{T}$ and $\vec{Δ t} = {[0, Δ t_{y}, Δ t_{z}]}^{T}$ (that is, Δs_y = 0 and Δt_x = 0), then the calculations could be performed even faster since Eq. (15) would become a two-dimensional scaled Fourier transform. We plan to investigate the advantages and disadvantages of this restriction as well as the computational complexity in the future. In addition, we also plan to confirm the calculation when the rotation is large, for example in the case that the angle between surface plane and destination plane is nearly 90 degree.

References and links

1. L. Ahrenberg, P. Benzie, M. Magnor, and J. Watson, “Computer generated holograms from three dimensional meshes using an analytic light transport model,” Appl. Opt. 47(10), 1567–1574 (2008). [CrossRef] [PubMed]

2. T. Tommasi and B. Bianco, “Computer-generated holograms of tilted planes by a spatial frequency approach,” J. Opt. Soc. Am. A 10(2), 299–305 (1993). [CrossRef]

3. K. Matsushima, H. Schimmel, and F. Wyrowski, “Fast calculation method for optical diffraction on titled planes by use of the angular spectrum of plane waves,” J. Opt. Soc. Am. A 20(9), 1755–1762 (2003). [CrossRef]

4. S. D. Nicola, A. Finizio, and G. Pierattini, “Angular spectrum method with correction of anamorphism for numerical reconstruction of digital holograms on tilted planes,” Opt. Express 13(24), 9935–9940 (2005). [CrossRef] [PubMed]

5. K. Matsushima, “Shifted angular spectrum method for off-axis numerical propagation,” Opt. Express 18(17), 18453–18463 (2010). [CrossRef] [PubMed]

6. N. Delen and B. Hooker, “Free-space beam propagation between arbitrarily oriented planes based on full diffraction theory: a fast Fourier transform approach,” J. Opt. Soc. Am. A 15(4), 857–867 (1998). [CrossRef]

7. R. P. Muffoletto, “Numerical techniques for fresnel diffraction in computational holography,” PhD thesis (Louisiana State University, 2006).

8. R. P. Muffoletto, J. M. Tyler, and J. E. Tohline, “Shifted Fresnel diffraction for computational holography,” Opt. Express 15(9), 5631–5640 (2007). [CrossRef] [PubMed]

9. L. Yu, U. An, and L. Cai, “Numerical reconstruction of digital holograms with variable viewing angles,” Opt. Express 10(22), 1250–1257 (2002). [PubMed]

10. J. Miura and T. Shimobaba, “Shifted-fresnel diffraction between two tilted planes (in Japanese),” Watake Seminar in Tohoku YS–6–52 (2008).

11. L. Onural, “Exact solution for scalar diffraction between tilted and translated planes using impulse functions over a surface,” J. Opt. Soc. Am. A 28(3), 290–295 (2011). [CrossRef]

12. J. W. Goodman, Introduction to Fourier Optics, 3rd ed. (Roberts and Company Publishers, Englewood, CO, 2005), Chap. 4.

13. D. H. Bailey and P. N. Swarztrauber, “The fractional fourier transform and applications,” SIAM Rev . 33(3), 389–404 (1991). [CrossRef]

λ	633 [nm]
$\| \vec{Δ s} \|$ of surface plane	0.707 [μm]
$\| \vec{Δ t} \|$ of surface plane	0.707 [μm]
Δu of incident side of lens	1 [μm]
Δv of incident side of lens	1 [μm]
$\| \vec{Δ s} \|$ of exiting side of lens	1 [μm]
$\| \vec{Δ t} \|$ of exiting side of lens	1 [μm]
Δu of image plane	1 [μm]
Δv of image plane	1 [μm]
Distance between center of surface plane and lens	4 [mm]
Distance between lens and center of image plane	8 [mm]
Focal length of lens	2.667 [mm]

N	Grid size (N × N)	t_N at Fresnel transform	t_N at proposed technique
32	32 × 32	112 [ms]	52 [ms]
64	64 × 64	1,734 [ms]	193 [ms]
128	128 × 128	27,647 [ms]	1,524 [ms]
256	256 × 256	441,904 [ms]	12,305 [ms]
512	512 × 512	7,065,681 [ms]	99,918 [ms]
1024	1024 × 1024	112,998,993 [ms]	805,339 [ms]

λ	633 [nm]
$\| \vec{Δ s} \|$ of surface plane	0.707 [μm]
$\| \vec{Δ t} \|$ of surface plane	0.707 [μm]
Δu of incident side of lens	1 [μm]
Δv of incident side of lens	1 [μm]
$\| \vec{Δ s} \|$ of exiting side of lens	1 [μm]
$\| \vec{Δ t} \|$ of exiting side of lens	1 [μm]
Δu of image plane	1 [μm]
Δv of image plane	1 [μm]
Distance between center of surface plane and lens	4 [mm]
Distance between lens and center of image plane	8 [mm]
Focal length of lens	2.667 [mm]

N	Grid size (N × N)	t_N at Fresnel transform	t_N at proposed technique
32	32 × 32	112 [ms]	52 [ms]
64	64 × 64	1,734 [ms]	193 [ms]
128	128 × 128	27,647 [ms]	1,524 [ms]
256	256 × 256	441,904 [ms]	12,305 [ms]
512	512 × 512	7,065,681 [ms]	99,918 [ms]
1024	1024 × 1024	112,998,993 [ms]	805,339 [ms]

Calculating the Fresnel diffraction of light from a shifted and tilted plane

Abstract

1. Introduction

2. Theory

3. Numerical experiment and calculation times

4. Conclusions

References and links

Supplementary Material (1)

Cited By

Figures (8)

Tables (2)

Equations (18)

Optics Express