Shifted angular spectrum method for off-axis numerical propagation

Kyoji Matsushima

doi:10.1364/OE.18.018453

1. Introduction

A variety of numerical methods for free space propagation are required for computational holography, such as digital holography (DH) or computer-generated holograms (CGH). DH can be regarded as analog to digital (A/D) conversion of optical wave fields, whereas CGH is digital to analog (D/A) conversion. Evolution of these conversions may contribute to the feasibility of digital signal processing (DSP) of light, in which numerical simulation of free-space propagation plays a role in the signal processing.

The first category of free-space propagation is propagation between parallel planes [1]. In this category, various methods, such as the single Fourier-transform-based Fresnel method (SFT-FR) [1–5] and the convolution-based Fresnel method (CV-FR) [1,6], are continually being proposed. However, the angular spectrum-based method (AS) [7] is potentially the most powerful, because it is rigorously derived from the Rayleigh-Sommerfeld integral. However, the traditional AS cannot serve as an all-round method in a numerical implementation due to sampling problems. The author recently proposed the band-limited AS (BL-AS) to avoid the sampling problems [8]. This is a simple, yet effective, improvement of the AS that magnifies the range of the effective propagation distance of the AS.

Another category of free space propagation is propagation between non-parallel planes [9–12]. The AS also plays an important role in this category. The rotational transformation of wave fields [9,10], formulated as an expansion of the AS, makes it possible to calculate wave fields in arbitrarily tilted planes from a given source field. This method is used for clear imaging of deeply tilted surfaces [10,13] and forming surface sources of light in the polygon-based method for creating CGHs of surface-modeled objects [14].

Recently, a new category of free space propagation, off-axis numerical propagation, was added to the field. This involves propagation between parallel planes, but with the sampling window of the output destination field shifted from that of the input source field. This is very useful in cases where a field is not paraxial and travels in an off-axis direction. In such cases, the sampling window needs to be expanded in conventional methods, as described in Section 2. This generally leads to a much higher computational effort. The techniques of off-axis numerical propagation can remove this difficulty. Off-axis numerical propagation is also useful for the propagation of large-scale wave fields that cannot be stored simultaneously in main memory due to the size of the data. An extremely high-definition CGH, reconstructing true fine 3D images, can also be calculated using off-axis numerical propagation [15].

The most notable method for off-axis numerical propagation is the shifted Fresnel method (Shift-FR) [5]. This excellent technique is derived from the SFT-FR using a scaled FFT and includes the capabilities of off-axis propagation and a variable sampling interval. However, the Shift-FR has a serious problem in that it causes strong aliasing when applied to short distance propagation for given sampling intervals and field sizes, as discussed in Section 4.

The novel method for off-axis numerical propagation proposed in this paper is a generalization of the BL-AS introduced in [8]. Therefore, the proposed method, called the shifted angular spectrum method (Shift-AS), does not cause any aliasing errors in short distance propagation, unlike the Shift-FR. The formulation thereof and numerical examples are presented in the following sections.

2. Off-axis numerical propagation

In numerical propagation of wave fields, the destination field is sometimes required to be in a region set apart from the optical axis in the destination plane. This situation is depicted in Fig. 1 . The source sampling area is placed around the origin of the coordinate system, whereas the region of interest in the destination plane is located apart from the origin. In this case, the destination field is generally calculated by the following procedure:

Fig. 1 Calculation of the diffracted field using (a) conventional methods and (b) methods for off-axis propagation, where the region of interest is apart from the optical axis in the destination pane.

Download Full Size | PDF

Step 1. Extend the source sampling window so that the sampling window includes the region of interest after numerical propagation.
Step 2. Pad the extended source sampling window with zeros.
Step 3. Numerically propagate the source wave field onto the destination plane.
Step 4. Cut out the region of interest from the destination sampling window.

Here, in Step 3, the source wave field is propagated using conventional numerical methods that do not change the sampling interval, such as the CV-FR, AS, or BL-AS. The diffracted field can be calculated using this procedure. However, the procedure usually requires a huge computational effort, especially in cases where the region of interest is far from the optical axis.

On the other hand, as shown in Fig. 1(b), a numerical method can compute directly the diffracted field in a sampling window set apart from the optical axis without extending the sampling window. This type of propagation is referred to as off-axis numerical propagation in this paper. The most notable method for off-axis numerical propagation is the Shift-FR [5], which makes it possible to change the sampling distance and the center of the sampling area in the destination plane by applying a scaled FFT to the SFT-FR method. Although this is an excellent method, it causes a serious problem in short-distance propagation as discussed in Section 4.

3. Formulation of the shifted angular spectrum method

3.1 The angular spectrum method for shifted coordinates

The coordinate system used for the formulation is shown in Fig. 2 . The source field $g (x, y, 0)$ is given in the source plane $(x, y, 0)$ , while the destination field $g (\hat{x}, \hat{y}, z_{0})$ is given in the destination plane $(\hat{x}, \hat{y}, z_{0})$ , in that the origin of the lateral coordinates $\hat{x}$ and $\hat{y}$ is shifted from the source coordinates as follows:

Fig. 2 The coordinate system used for formulation.

Download Full Size | PDF

\begin{matrix} \hat{x} = x - x_{0} \\ \hat{y} = y - y_{0} \end{matrix}

Without any shifts ( $x_{0} = y_{0} = 0$ ), the spectrum of the destination field is given by a convolution form of the Rayleigh-Sommerfeld integral (see Section 2.1 in [8]).

\begin{matrix} \hat{g} (\hat{x}, \hat{y}, z_{0}) = \hat{g} (x, y, z_{0}) = g (x, y, z_{0}) \\ = \iint g (x^{'}, y^{'}, 0) h (x - x^{'}, y - y^{'}, z_{0}) d x^{'} d y^{'} \end{matrix}

where

h (x, y, z_{0})

is the propagation kernel (the impulse response) of the Rayleigh-Sommerfeld formula. This convolution is rewritten using the convolution theorem as:

G (u, v; z_{0}) = G (u, v; 0) H (u, v; z_{0}),

where u and v are Fourier frequencies with respect to x and y, respectively. The spectrum of the source field is given as:

\begin{matrix} G (u, v; 0) = \iint g (x, y, 0) \exp [- i 2 π (u x + v y)] d x d y \\ = F {g (x, y, 0)} . \end{matrix}

Here, F represents the Fourier transform. The transfer function $H (u, v; z_{0})$ is given by:

\begin{matrix} H (u, v; z_{0}) = \exp [i 2 π w z_{0}], \\ w = w (u, v) = {\begin{array}{l} {(λ^{- 2} - u^{2} - v^{2})}^{1 / 2} & \dots u^{2} + v^{2} \leq λ^{- 2} \\ 0 & \dots otherwise \end{array} \end{matrix}

where λ denotes the wavelength.

In cases with shifted coordinates ( $x_{0} \neq 0$ and/or $y_{0} \neq 0$ ), Eq. (2) is rewritten as:

\begin{matrix} \hat{g} (\hat{x}, \hat{y}, z_{0}) = g (\hat{x} + x_{0}, \hat{y} + y_{0}, z_{0}) \\ = \iint g (x^{'}, y^{'}, 0) h (\hat{x} + x_{0} - x^{'}, \hat{y} + y_{0} - y^{'}, z_{0}) d x^{'} d y^{'} . \end{matrix}

This can be rewritten using the convolution theorem and shift theorem of the Fourier transform as:

\begin{matrix} \hat{G} (u, v; z_{0}) = G (u, v; 0) \hat{H} (u, v; z_{0}), \\ \hat{H} (u, v; z_{0}) = H (u, v; z_{0}) \exp [i 2 π (x_{0} u + y_{0} v)] \\ = \exp [i 2 π (x_{0} u + y_{0} v + z_{0} w)] \end{matrix}

Consequently, we obtain the destination field as follows:

\hat{g} (\hat{x}, \hat{y}, z_{0}) = F^{- 1} {\hat{G} (u, v; z_{0})} .

Nevertheless, simple numerical implementations of the above procedure are most likely to cause considerable sampling problems as described in the next section.

3.2 One-dimensional wave fields

For simplicity, one-dimensional wave fields that are a function of x are discussed in this section. In this case, the transfer function of Eq. (7) is redefined as follows:

\begin{matrix} \hat{H} (u; z_{0}) = \exp [i ϕ (u; z_{0})], \\ ϕ (u; z_{0}) = 2 π [x_{0} u + z_{0} {(λ^{- 2} - u^{2})}^{1 / 2}] . \end{matrix}

The local signal frequency of the transfer function is given by [16]:

\begin{matrix} f_{u} = \frac{1}{2 π} \frac{\partial ϕ}{\partial u} \\ = x_{0} - \frac{u}{{(λ^{- 2} - u^{2})}^{1 / 2}} z_{0}, \end{matrix}

where the local signal frequency denotes neither physical frequencies of time nor space, but the frequency of peaks and valleys in the function

H (u; z_{0})

in a certain period u. Therefore, supposing that the transfer function is sampled at intervals of

Δ u

, the Nyquist theorem requires the following relation to avoid sampling problems:

Δ u^{- 1} \geq 2 | f_{u} | .

By substituting Eq. (10) into Eq. (11), the Nyquist condition is rewritten as:

(x_{0} - \frac{1}{2 Δ u}) \frac{1}{z_{0}} \leq \frac{u}{{(λ^{- 2} - u^{2})}^{1 / 2}} \leq (x_{0} + \frac{1}{2 Δ u}) \frac{1}{z_{0}}

Finally, the conditions are summarized in the following three cases:

{\begin{array}{l} u_{limit}^{(-)} < u < u_{limit}^{(+)} & \dots if & S_{x} < x_{0} \\ - u_{limit}^{(-)} < u < u_{limit}^{(+)} & \dots if & - S_{x} \leq x_{0} < S_{x} \\ - u_{limit}^{(-)} < u < - u_{limit}^{(+)} & \dots if & x_{0} \leq - S_{x} \end{array}

where

S_{x} = {(2 Δ u)}^{- 1}

is the size of the source sampling window. Here, it should be noted that the size of the sampling window is given by

S_{x} = {(2 Δ u)}^{- 1}

, and not by

S_{x} = Δ u^{- 1}

, since the source sampling window is doubled to linearize the discrete convolution, as described in Section 2.2 in [8]. The constant

u_{limit}^{(\pm)}

is given by

u_{limit}^{(\pm)} \equiv {[{(x_{0} \pm \frac{1}{2 Δ u})}^{- 2} z_{0}^{2} + 1]}^{- 1 / 2} λ^{- 1} .

Where

x_{0} = 0

, this agrees with

u_{limit} \equiv {[{(2 Δ u z)}^{2} + 1]}^{- 1 / 2} λ^{- 1}

defined in Eq. (13) in [8].

An example of the sampled function of $\hat{H} (u; z_{0})$ is depicted in Fig. 3 . Here, the center of the destination sampling window is shifted by an amount half the size of the sampling window in this example, i.e., $x_{0} = + S_{x} / 2$ is assumed. Since aliasing errors occur below $- u_{limit}^{(-)}$ or above $u_{limit}^{(+)}$ , in this case the sampled transfer function and the source field should be band-limited within $- u_{limit}^{(-)} < u < u_{limit}^{(+)}$ to avoid aliasing errors.

Fig. 3 Example of a sampled transfer function. Only the real part of $\hat{H} (u; z_{0})$ is depicted in the sampling interval $Δ u = {(2 S_{x})}^{- 1} = {(2 N_{x} Δ x)}^{- 1}$ , where $N_{x} = 1024$ , $Δ x = 2 λ$ , $x_{0} = + S_{x} / 2$ , and $z_{0} = 20 S_{x}$ .

Download Full Size | PDF

Consequently, the transfer function for shifted coordinates is given by:

{\hat{H}}^{'} (u; z_{0}) = \hat{H} (u; z_{0}) rect (\frac{u - u_{0}}{u_{width}}),

where

rect(ξ)

is a rectangular function with unity width. The constants

u_{0}

and

u_{width}

are given in Table 1 .

Table 1. Constants used for the band-limit to avoid aliasing errors

View Table

3.3 Two-dimensional wave fields

In two-dimensional cases, the Nyquist condition is given as two-dimensional regions in the $(u, v)$ plane. The local signal frequencies of the function $\hat{H} (u, v; z_{0}) = \exp [i ϕ (u, v; z_{0})]$ are given by:

f_{u} = \frac{1}{2 π} \frac{\partial ϕ}{\partial u} = x_{0} - \frac{u}{w} z_{0}, f_{v} = \frac{1}{2 π} \frac{\partial ϕ}{\partial v} = y_{0} - \frac{v}{w} z_{0},

where

ϕ (u, v; z_{0}) = 2 π (x_{0} u + y_{0} v + z_{0} w)

. The Nyquist conditions are given by:

Δ u^{- 1} \geq 2 | f_{u} |,

Δ v^{- 1} \geq 2 | f_{v} |,

where

Δ v = {(2 S_{y})}^{- 1}

is the sampling interval and S_y is the size of the source sampling window in the y direction. Note that the sampling window is also doubled in the y direction to linearize the discrete convolution.

Both conditions Eqs. (17) and (18) should be satisfied simultaneously to avoid sampling problems in two-dimensional wave fields. Since these conditions are very similar, let us formulate the region that satisfies only condition Eq. (17) in the $(u, v)$ space as follows:

Case (i): $S_{x} < x_{0}$
$u \geq 0 and \frac{u^{2}}{{[u_{limit}^{(-)}]}^{2}} + \frac{v^{2}}{λ^{- 2}} \geq 1 and \frac{u^{2}}{{[u_{limit}^{(+)}]}^{2}} + \frac{v^{2}}{λ^{- 2}} \leq 1,$
Case (ii): $- S_{x} \leq x_{0} < S_{x}$
$(\frac{u^{2}}{{[u_{limit}^{(-)}]}^{2}} + \frac{v^{2}}{λ^{- 2}} \leq 1 and u \leq 0) or (\frac{u^{2}}{{[u_{limit}^{(+)}]}^{2}} + \frac{v^{2}}{λ^{- 2}} \leq 1 and u > 0),$
Case (iii): $x_{0} \leq - S_{x}$
$u \leq 0 and \frac{u^{2}}{{[u_{limit}^{(+)}]}^{2}} + \frac{v^{2}}{λ^{- 2}} \geq 1 and \frac{u^{2}}{{[u_{limit}^{(-)}]}^{2}} + \frac{v^{2}}{λ^{- 2}} \leq 1.$

Here, the region satisfying the other condition [Eq. (18)] is simply obtained by switching the symbols x and u to y and v in the above relations, respectively. Note that the region under the conditions $x_{0} = 0$ or $y_{0} = 0$ , that is, the special case of (ii) is given by:

\frac{u^{2}}{u_{limit}^{2}} + \frac{v^{2}}{λ^{- 2}} \leq 1, \frac{u^{2}}{λ^{- 2}} + \frac{v^{2}}{v_{limit}^{2}} \leq 1.

These agree with relations Eqs. (18) and (19) in [8].

Relations Eqs. (19)–(21) for shifting $x_{0}$ give a region specified by the combination of vertical ellipsoidal regions with a major diameter 2λ ⁻¹ in the (u, v) plane, while the relations for shifting $y_{0}$ give horizontal ellipsoidal regions. Figure 4 shows these ellipsoidal regions in cases where $x_{0} \neq 0$ and $y_{0} = 0$ . In this figure, the overlap regions that satisfy both Nyquist conditions Eqs. (17) and (18) are schematically depicted as red-hatched areas. The sampled transfer function should be limited in the overlap region to avoid sampling problems. This is equivalent to limiting the bandwidth of the source field within the overlap region.

Fig. 4 Schematic illustrations of frequency regions to avoid aliasing errors of the sampled transfer function. Here, $x_{0} \neq 0$ and $y_{0} = 0$ .

Download Full Size | PDF

Here, note that the inside ellipse in Case (i) is switched to the outside ellipse in Case (iii). This switching occurs when $x_{0}$ changes its sign in Case (ii). Although Eq. (20) seems, in a sense, to be written as a separate form, the mathematical expression of Case (ii) is the same, irrespective of the sign of $x_{0}$ .

3.4 Approximated rectangular region to avoid aliasing errors

Regions that avoid sampling problems are considerably more complicated in two-dimensional wave fields. However, in most cases where $z_{0} > S_{x} and S_{y}$ , the two-dimensional region can be approximated as a simple rectangular region by combining the one-dimensional relations of Eq. (13). An example of the approximated region is illustrated schematically in Fig. 5 .

Fig. 5 Schematic illustrating the approximated rectangular region for the band-limit of a two-dimensional wave field. Here, the region is depicted with $x_{0} > + S_{x}$ and $y_{0} < - S_{y}$ as an example.

Download Full Size | PDF

Finally, the sampled transfer function for two-dimensional wave fields is approximately given by:

{\hat{H}}^{'} (u, v; z_{0}) = \hat{H} (u, v; z_{0}) rect (\frac{u - u_{0}}{u_{width}}) rect (\frac{v - v_{0}}{v_{width}}),

where the additional constants are listed in Table 1 with:

v_{limit}^{(\pm)} \equiv {[{(y_{0} \pm \frac{1}{2 Δ v})}^{- 2} z_{0}^{2} + 1]}^{- 1 / 2} λ^{- 1} .

4. Numerical verification of the shifted angular spectrum method

The optical setup used for numerical verification is shown in Fig. 6(a) . We assume that a plane wave travels at an incident angle of θ. Supposing that the plane wave is diffracted by a circular aperture, we calculate the diffracted field in the plane at position $z_{0}$ . Here, the source field including the plane wave and the aperture is sampled at intervals of $Δ x = Δ y = 8 μ m$ and the number of sampling points is $N_{x} = N_{y} = 1024$ . Thus, the source sampling area is $S_{x} \times S_{y} = 8.192 \times 8.192 {mm}^{2}$ . The incident angle and wavelength of the plane wave are $θ = {1.5}^{\circ}$ and $λ = 532 nm$ , respectively. The amplitude and phase images of the sampled source field are shown in Fig. 6(b).

Fig. 6 Setup for (a) numerical simulation of off-axis propagation, and (b) the sampled source field used in the simulation.

Download Full Size | PDF

The diffracted field is calculated using three different methods: the band-limited angular spectrum (BL-AS), shifted Fresnel (Shift-FR), and shifted angular spectrum method (Shift-AS) proposed in this paper. Since the BL-AS is not a method for off-axis numerical propagation, the number of sampling points in the source field is doubled both on the x- and y-axes so that the diffracted field is included within the destination sampling window. In addition, the approximated transfer function given by Eq. (23) is used for calculating the Shift-AS.

The calculated destination field is shown in Fig. 7(a) –7(c) for three positions of z ₀. The shift of the destination sampling window is x ₀ = 0.4 cm both for z ₀ = 5 cm in (a) and 10 cm in (b), whereas x ₀ = 1 cm for z ₀ = 40 cm in (c). In the case of z ₀ = 5 cm in (a), the position of the destination field is shifted slightly along the x-axis because of the slanted plane wave. In this case, the field calculated by the Shift-AS agrees with that by the BL-AS, whereas the field by the Shift-FR disagrees with these fields because of the aliasing error. The field calculated by the Shift-FR is also not normal in z ₀ = 10 cm in (b), while the BL-AS and Shift-AS again give similar fields.

Fig. 7 Amplitude images of the destination fields calculated for different propagation and shift distances using three different methods.

Download Full Size | PDF

In the case of z ₀ = 40 cm in (c), the destination field lies outside the extended sampling area of the BL-AS. However, the Shift-AS and Shift-FR can calculate the field within the sampling area owing to their off-axis property. The Shift-FR, in this case, gives the same result as that using the Shift-AS, due to the long distance of the propagation.

5. Discussion on the limit frequency

We limit the spatial frequencies of the source fields to a band specified by $u_{limit}^{(\pm)}$ . This limit frequency can be physically interpreted as the local spatial frequency of the field emitted from a point within the source sampling window, as shown in Fig. 8 . The highest spatial frequency is given by the field emitted from the point at the lower end of the source sampling window to the point at the upper end of the destination window. The lowest frequency is given in the same manner. Here, angles $θ_{\max}$ and $θ_{\min}$ in Fig. 8 are given by the geometry of the source and destination window as follows:

Fig. 8 Model for discussing the highest and lowest spatial frequencies required for off-axis numerical propagation.

Download Full Size | PDF

\begin{matrix} \sin θ_{\max} = {(x_{0} + S_{x}) / [z_{0} + {(x_{0} + S_{x})}^{2}]}^{1 / 2}, \\ \sin θ_{\min} = {(x_{0} - S_{x}) / [z_{0} + {(x_{0} - S_{x})}^{2}]}^{1 / 2} . \end{matrix}

By substituting $S_{x} = {(2 Δ u)}^{- 1}$ , the highest and lowest spatial frequencies are given by:

\begin{array}{l} u_{high} = \sin θ_{\max} / λ = u_{limit}^{(+)}, \\ u_{low} = \sin θ_{\min} / λ = u_{limit}^{(-)} . \end{array}

As a result, the limits imposed on the frequencies by the procedure to avoid sampling problems are in agreement with the highest and lowest frequencies. Therefore, we can interpret the condition specified by $u_{limit}^{(\pm)}$ as the frequency band that is required to physically propagate the source wave field onto the shifted destination area.

6. Conclusion

A numerical technique called the shifted angular spectrum method is proposed for off-axis numerical propagation. The Shift-AS is a generalization of the band-limited angular spectrum method, and therefore, the method does not have any defined restrictions on the propagation distance. This new method makes it possible to calculate exact wave fields without any aliasing errors, such as those caused by the shifted Fresnel method.

Acknowledgments

This work was partially supported by the Japan Society for the Promotion of Science (JSPS.KAKENHI) (no. 21500114).

References and links

1. T. M. Kreis, M. Adams, and W. P. O. Jüptner, “Methods of digital holography: A comparison,” Proc. SPIE 3098, 224–233 (1997). [CrossRef]

2. F. Zhang, I. Yamaguchi, and L. P. Yaroslavsky, “Algorithm for reconstruction of digital holograms with adjustable magnification,” Opt. Lett. 29(14), 1668–1670 (2004). [CrossRef] [PubMed]

3. L. Yu and M. K. Kim, “Pixel resolution control in numerical reconstruction of digital holography,” Opt. Lett. 31(7), 897–899 (2006). [CrossRef] [PubMed]

4. D. Wang, J. Zhao, F. Zhang, G. Pedrini, and W. Osten, “High-fidelity numerical realization of multiple-step Fresnel propagation for the reconstruction of digital holograms,” Appl. Opt. 47(19), D12–D20 (2008). [CrossRef] [PubMed]

5. R. P. Muffoletto, J. M. Tyler, and J. E. Tohline, “Shifted Fresnel diffraction for computational holography,” Opt. Express 15(9), 5631–5640 (2007), http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-15-9-5631. [CrossRef] [PubMed]

6. M. Sypek, “Light propagation in the Fresnel region. New numerical approach,” Opt. Commun. 116(1-3), 43–48 (1995). [CrossRef]

7. J. W. Goodman, Introduction to Fourier Optics, 2nd ed. (McGraw-Hill, 1996), chap. 3.10.

8. K. Matsushima and T. Shimobaba, “Band-limited angular spectrum method for numerical simulation of free-space propagation in far and near fields,” Opt. Express 17(22), 19662–19673 (2009), http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-17-22-19662. [CrossRef] [PubMed]

9. K. Matsushima, H. Schimmel, and F. Wyrowski, “Fast calculation method for optical diffraction on tilted planes by use of the angular spectrum of plane waves,” J. Opt. Soc. Am. A 20(9), 1755–1762 (2003), http://www.opticsinfobase.org/josaa/abstract.cfm?URI=josaa-20-9-1755. [CrossRef]

10. K. Matsushima, “Formulation of the rotational transformation of wave fields and their application to digital holography,” Appl. Opt. 47(19), D110–D116 (2008), http://www.opticsinfobase.org/ao/abstract.cfm?URI=ao-47-19-D110. [CrossRef] [PubMed]

11. 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]

12. S. J. Jeong and C. K. Hong, “Pixel-size-maintained image reconstruction of digital holograms on arbitrarily tilted planes by the angular spectrum method,” Appl. Opt. 47(16), 3064–3071 (2008). [CrossRef] [PubMed]

13. S. De Nicola, A. Finizio, G. Pierattini, P. Ferraro, and D. Alfieri, “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]

14. K. Matsushima, “Computer-generated holograms for three-dimensional surface objects with shade and texture,” Appl. Opt. 44(22), 4607–4614 (2005), http://www.opticsinfobase.org/ao/abstract.cfm?URI=ao-44-22-4607. [CrossRef] [PubMed]

15. 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), http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-17-22-19662. [CrossRef] [PubMed]

16. J. W. Goodman, Introduction to Fourier Optics, 2nd ed. (McGraw-Hill, 1996), chap. 2.2.

Case	$u_{0}$	$u_{width}$	Case	$v_{0}$	$v_{width}$
$S_{x} < x_{0}$	$(u_{limit}^{(+)} + u_{limit}^{(-)}) / 2$	$u_{limit}^{(+)} - u_{limit}^{(-)}$	$S_{y} < y_{0}$	$(v_{limit}^{(+)} + v_{limit}^{(-)}) / 2$	$v_{limit}^{(+)} - v_{limit}^{(-)}$
$- S_{x} \leq x_{0} < S_{x}$	$(u_{limit}^{(+)} - u_{limit}^{(-)}) / 2$	$u_{limit}^{(+)} + u_{limit}^{(-)}$	$- S_{y} \leq y_{0} < S_{y}$	$(v_{limit}^{(+)} - v_{limit}^{(-)}) / 2$	$v_{limit}^{(+)} + v_{limit}^{(-)}$
$x_{0} \leq - S_{x}$	$- (u_{limit}^{(+)} + u_{limit}^{(-)}) / 2$	$u_{limit}^{(-)} - u_{limit}^{(+)}$	$y_{0} \leq - S_{y}$	$- (v_{limit}^{(+)} + v_{limit}^{(-)}) / 2$	$v_{limit}^{(-)} - v_{limit}^{(+)}$

Shifted angular spectrum method for off-axis numerical propagation

Abstract

1. Introduction

2. Off-axis numerical propagation

3. Formulation of the shifted angular spectrum method

3.1 The angular spectrum method for shifted coordinates

3.2 One-dimensional wave fields

3.3 Two-dimensional wave fields

3.4 Approximated rectangular region to avoid aliasing errors

4. Numerical verification of the shifted angular spectrum method

5. Discussion on the limit frequency

6. Conclusion

Acknowledgments

References and links

Cited By

Figures (8)

Tables (1)

Equations (26)

Optics Express