Abstract

We describe a numerical method that can be used to calculate the propagation of light in a medium of constant (possibly complex) index of refraction n. The method integrates the Rayleigh–Sommerfeld diffraction integral numerically. After an appropriate change of integration variables, the integrand of the diffraction integral is split into a slowly varying and an (often fast) oscillating quadratic factor. The slowly varying factor is approximated by a spline fit, and the resulting Fresnel integrals are subsequently integrated exactly. Although the method is not as fast as methods involving a fast Fourier transform, such as plane-wave propagation or Fresnel approximation, it is accurate over a greater range than these methods.

© 2005 Optical Society of America

Full Article  |  PDF Article

References

  • View by:
  • |
  • |
  • |

  1. J. J. Stamnes, Waves in Focal Regions (Institute of Physics, Bristol, UK, 1986).
  2. M. Mansuripur, “Distribution of light at and near the focus of high-numerical-aperture objectives,” J. Opt. Soc. Am. A 3, 2086–2093 (1986).
    [CrossRef]
  3. P. Török, S. J. Hewlett, P. Varga, “On the series expansion of high-aperture, vectorial diffraction integrals,” J. Mod. Opt. 44, 493–503 (1997).
    [CrossRef]
  4. C. J. R. Sheppard, P. Török, “Efficient calculation of electromagnetic diffraction in optical systems using a multipole expansion,” J. Mod. Opt. 44, 803–818 (1997).
    [CrossRef]
  5. S. Stallinga, “Axial birefringence in high-numerical-aperture optical systems and the light distribution close to focus,” J. Opt. Soc. Am. A 18, 2846–2859 (2001).
    [CrossRef]
  6. P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
    [CrossRef]
  7. J. M. Brok, H. P. Urbach, “Rigorous model of the scat-tering of a focused spot by a grating and its application in optical recording,” J. Opt. Soc. Am. A 20, 256–272 (2003).
    [CrossRef]
  8. T. Gravelsaeter, J. J. Stamnes, “Diffraction by circular apertures. 1: Method of linear phase and amplitude approximation,” Appl. Opt. 21, 3644–3651 (1982).
    [CrossRef] [PubMed]
  9. J. J. Stamnes, “Hybrid integration technique for efficient and accurate computation of diffraction integrals,” J. Opt. Soc. Am. A 6, 1330–1342 (1989).
    [CrossRef]
  10. L. d’Arcio, J. J. M. Braat, H. J. Frankema, “Numerical evaluation of diffraction integrals for apertures of complicated shape,” J. Opt. Soc. Am. A 11, 2664–2674 (1994).
    [CrossRef]
  11. See, e.g., J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, New York, 1996).
  12. J. A. C. Veerman, “Evaluation of spot propagation calculations in 1D,” Internal Note, 2001.
  13. J. Stoer, Einführung in die Numerische Mathematik (Springer Verlag, Berlin, 1979).
  14. See http://www.tgs.com .

2003 (1)

2001 (2)

1997 (3)

P. Török, S. J. Hewlett, P. Varga, “On the series expansion of high-aperture, vectorial diffraction integrals,” J. Mod. Opt. 44, 493–503 (1997).
[CrossRef]

C. J. R. Sheppard, P. Török, “Efficient calculation of electromagnetic diffraction in optical systems using a multipole expansion,” J. Mod. Opt. 44, 803–818 (1997).
[CrossRef]

P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
[CrossRef]

1994 (1)

1989 (1)

1986 (1)

1982 (1)

Braat, J. J. M.

Brok, J. M.

d’Arcio, L.

Frankema, H. J.

Goodman, J. W.

See, e.g., J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, New York, 1996).

Gravelsaeter, T.

Hewlett, S. J.

P. Török, S. J. Hewlett, P. Varga, “On the series expansion of high-aperture, vectorial diffraction integrals,” J. Mod. Opt. 44, 493–503 (1997).
[CrossRef]

Mansuripur, M.

Morton, R. D.

P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
[CrossRef]

Put, P. L. M.

P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
[CrossRef]

Rusch, J. J.

P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
[CrossRef]

Sheppard, C. J. R.

C. J. R. Sheppard, P. Török, “Efficient calculation of electromagnetic diffraction in optical systems using a multipole expansion,” J. Mod. Opt. 44, 803–818 (1997).
[CrossRef]

Stallinga, S.

Stamnes, J. J.

Stoer, J.

J. Stoer, Einführung in die Numerische Mathematik (Springer Verlag, Berlin, 1979).

Török, P.

P. Török, S. J. Hewlett, P. Varga, “On the series expansion of high-aperture, vectorial diffraction integrals,” J. Mod. Opt. 44, 493–503 (1997).
[CrossRef]

C. J. R. Sheppard, P. Török, “Efficient calculation of electromagnetic diffraction in optical systems using a multipole expansion,” J. Mod. Opt. 44, 803–818 (1997).
[CrossRef]

Urbach, H. P.

J. M. Brok, H. P. Urbach, “Rigorous model of the scat-tering of a focused spot by a grating and its application in optical recording,” J. Opt. Soc. Am. A 20, 256–272 (2003).
[CrossRef]

P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
[CrossRef]

Varga, P.

P. Török, S. J. Hewlett, P. Varga, “On the series expansion of high-aperture, vectorial diffraction integrals,” J. Mod. Opt. 44, 493–503 (1997).
[CrossRef]

Veerman, J. A. C.

J. A. C. Veerman, “Evaluation of spot propagation calculations in 1D,” Internal Note, 2001.

Appl. Opt. (1)

Internal Note (1)

J. A. C. Veerman, “Evaluation of spot propagation calculations in 1D,” Internal Note, 2001.

J. Mod. Opt. (2)

P. Török, S. J. Hewlett, P. Varga, “On the series expansion of high-aperture, vectorial diffraction integrals,” J. Mod. Opt. 44, 493–503 (1997).
[CrossRef]

C. J. R. Sheppard, P. Török, “Efficient calculation of electromagnetic diffraction in optical systems using a multipole expansion,” J. Mod. Opt. 44, 803–818 (1997).
[CrossRef]

J. Opt. Soc. Am. A (5)

Jpn. J. Appl. Phys. (1)

P. L. M. Put, H. P. Urbach, R. D. Morton, J. J. Rusch, “Resolution limit of optical disc mastering,” Jpn. J. Appl. Phys. 36, 539–548 (1997).
[CrossRef]

Other (4)

See, e.g., J. W. Goodman, Introduction to Fourier Optics (McGraw-Hill, New York, 1996).

J. Stoer, Einführung in die Numerische Mathematik (Springer Verlag, Berlin, 1979).

See http://www.tgs.com .

J. J. Stamnes, Waves in Focal Regions (Institute of Physics, Bristol, UK, 1986).

Cited By

OSA participates in CrossRef's Cited-By Linking service. Citing articles from OSA journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (10)

Fig. 1
Fig. 1

Important circles for the point ( x 0 ,   y 0 ) in the case of a rectangular diaphragm Ω. The integral over ϕ runs over the overlap region of Ω and the area enclosed in two adjacent circles. Note that when the point ( x 0 ,   y 0 ) lies within Ω the first circle has radius r = 0 .

Fig. 2
Fig. 2

Absolute value of the field for the Gauss spot for i z = 0   ( z 0 = 0 ) (the initial field).

Fig. 3
Fig. 3

Absolute value of the field for the bloc spot for i z = 0   ( z 0 = 0 ) (the initial field). Owing to the coarseness of the grid, the spot is not entirely circular.

Fig. 4
Fig. 4

Absolute value of the field for the Airy spot for i z = 0   ( z 0 = 0 ) (the initial field). The central spot is well contained within the box. The sidelobes fall partially outside the box.

Fig. 5
Fig. 5

Absolute value of the field for the Gauss spot for i z = 4   ( z 0 = 64 λ ) as calculated by the PWP method. In this case the spot shows aliasing effects when the edge of the box is approached.

Fig. 6
Fig. 6

Absolute value of the field for the Gauss spot for i z = 4   ( z 0 = 64 λ ) as calculated by the RS method. Although the spot nears the edge of the calculational box, the spot remains smooth.

Fig. 7
Fig. 7

Absolute value of the field for the block spot for i z = 2   ( z 0 = 32 λ ) calculated by the PWP method. The squarelike symmetry remains present.

Fig. 8
Fig. 8

Absolute value of the field for the bloc spot for i z = 2   ( z 0 = 32 λ ) calculated by the RS method. The spot has propagated and becomes smoother than the initial spot.

Fig. 9
Fig. 9

Absolute value of the field for the Airy spot for i z = 2   ( z 0 = 32 λ ) , calculated by the PWP method. The central spot agrees well with that of the RS method. For the sidelobes, a clear difference is visible.

Fig. 10
Fig. 10

Absolute value of the field for the Airy spot for i z = 2   ( z 0 = 32 λ ) , calculated by the RS method. Owing to the cutoff in the calculational box, the sidelobes are not exactly circular.

Equations (85)

Equations on this page are rendered with MathJax. Learn more.

I = 0 s 1 G ( s ) exp ( iks 2 ) d s ,
F ( U ) ( k x ,   k y ) = 1 4 π 2   - - U ( x ,   y ,   z = 0 ) × exp ( - ik x x - ik y y ) d x d y .
k z = ( 2 π n ) 2 λ 2 - k x 2 - k y 2 1 / 2 .
U ( x ,   y ,   z = z 0 ) = - - F ( U ) ( k x ,   k y ) × exp ( ik z z 0 ) exp ( ik x x + ik y y ) d k x d k y .
U ( x ,   y ,   z 0 ) - i   exp 2 π i   z 0 n λ n z 0 λ × exp π in z 0 λ   ( x 2 + y 2 ) × F U ( x ,   y ,   z = 0 ) exp π in z 0 λ   ( x 2 + y 2 ) exp 2 π n xz 0   ( x 0 x + y 0 y ) xn z 0 λ ,   yn z 0 λ d x d y ,
n = n + in
z 0 3 π n 4 λ   [ ( x - x ) 2 + ( y - y ) 2 ] max 2 ,
U ( x 0 ,   y 0 ,   z 0 ) = - Ω U ( x ,   y ,   0 )   2 z 0 | r 0 - r | × ik - 1 | r 0 - r |   exp ( ik | r 0 - r | ) 4 π | r 0 - r |   d x d y ,
k = ω 0 μ 0 n .
Ω ˜ = Ω - ( x 0 ,   y 0 ,   0 ) ,
U ( x 0 ,   y 0 ,   z 0 ) = - Ω ˜ U ( x 0 + x ,   y 0 + y ,   0 )   2 z 0 ( x 2 + y 2 + z 0 2 ) 1 / 2 × ik - 1 ( x 2 + y 2 + z 0 2 ) 1 / 2 × exp [ ik ( x 2 + y 2 + z 0 2 ) 1 / 2 ] 4 π ( x 2 + y 2 + z 0 2 ) 1 / 2   d x d y .
r = ( x 2 + y 2 ) 1 / 2 ,
x = r   cos   ϕ , y = r   sin   ϕ
ρ = ( x 2 + y 2 + z 0 2 ) 1 / 2 - z 0 = ( r 2 + z 0 2 ) 1 / 2 - z 0 .
r = ρ 1 + 2 z 0 ρ 1 / 2 .
U ( x 0 ,   y 0 ,   z 0 ) = Ω ˜ ˜ f ( ρ ,   ϕ )   exp ( ik ρ ) ρ + z 0   d ϕ d ρ ,
f ( ρ ,   ϕ ) = z 0 n   exp ( ikz 0 ) i λ   U { x 0 + [ ( ρ + z 0 ) 2 - z 0 2 ] 1 / 2 × cos   ϕ ,   y 0 + [ ( ρ + z 0 ) 2 - z 0 2 ] 1 / 2 × sin   ϕ ,   0 } 1 - 1 ik ( ρ + z 0 )
U ( x 0 ,   y 0 ,   z 0 ) = Ω ˜ ˜ f ( ρ ,   ϕ )   exp ( ik ρ ) ρ + z 0   d ϕ d ρ = l = 1 N ρ 1 ( l ) ρ 2 ( l )   exp ( ik ρ ) F l ( ρ ) d ρ ,
F l ( ρ ) = 1 ρ + z 0   j = 1 M ( l ) ϕ 2 j - 1 ( ρ ) ϕ 2 j ( ρ ) f ( ρ ,   ϕ ) d ϕ .
U ( x 0 ,   y 0 ,   z 0 ) = l = 1 N ρ 1 ( l ) ρ 2 ( l ) F l ( ρ ) exp ( ik ρ ) d ρ .
ρ ( ϕ ) ρ 1 ( l ) × C 2   [ ϕ - ϕ 2 j - 1 ( ρ ) ] 2 ,
F l ( ρ ) = a 0 + a 1 ( ρ - ρ 1 ) 1 / 2 + a 2 ( ρ - ρ 1 ) + a 3 ( ρ - ρ 1 ) 3 / 2 + ,
F l ( ρ ) = b 0 + b 1 ( ρ 2 - ρ ) 1 / 2 + b 2 ( ρ 2 - ρ ) + b 3 ( ρ 2 - ρ ) 3 / 2 +   .
ρ 1 ρ 2   exp ( ik ρ ) F l ( ρ ) d ρ = ρ 1 ( ρ 1 + ρ 2 ) / 2   exp ( ik ρ ) F l ( ρ ) d ρ + ( ρ 1 + ρ 2 ) / 2 ρ 2   exp ( ik ρ ) F l ( ρ ) d ρ = exp ( ik ρ 1 ) 0 ( ρ 2 - ρ 1 ) / 2   exp ( iks 2 ) F left ( s ) d s + exp ( ik ρ 2 ) 0 ( ρ 2 - ρ 1 ) / 2 exp ( iks 2 ) F right ( s ) d s * ,
s = ρ - ρ 1 ,
s = ρ 2 - ρ ,
F left ( s ) = 2 s F l ( ρ 1 + s 2 ) ,
F right ( s ) = 2 s F l ( ρ 2 - s 2 ) * .
I = 0 ( ρ 2 - ρ 1 ) / 2   exp ( iks 2 ) G ( s ) d s .
0 ( ρ 2 - ρ 1 ) / 2 exp ( iks 2 ) G ( s ) d s - 0 ( ρ 2 - ρ 1 ) / 2 exp ( iks 2 ) σ ( s ) d s 0 ( ρ 2 - ρ 1 ) / 2 exp ( iks 2 ) G ( s ) d s C   h 2 ( h + 1 ) k ,
n = 1 .
λ = 0.4579 × 10 - 6     m .
L x = L y = L = 0.011424 × 10 - 3     m .
z 0 ( i z ) = 0.0075 × 10 - 3 i z     m ( i z = 1 , , 4 ) ,
z 0 ( 0 ) = 1 × 10 - 8     m ;
N x = N y = 40 ,
x ( l ,   j ) = l - N x 2 L , l = 0 , ,   N x - 1 ,
j = 0 , N y - 1 ,
y ( l ,   j ) = j - N y 2 L , l = 0 , ,   N x - 1 ,
j = 0 , N y - 1 ,
x FA ( l ,   j ,   i z ) = k x ( l ,   j ) z 0 ( i z ) λ ,
y FA ( l ,   j ,   i z ) = k y ( l ,   j ) z 0 ( i z ) λ .
rel ( x ,   y ,   z ) = | f 2 ( x ,   y ,   z ) - f 1 ( x ,   y ,   z ) | f 1 ,
f Gauss = π   a ( 140 ) 2   exp - a 2 x 2 L 2 + y 2 L 2 ,
a = 256 .
f bloc = 1 if x 2 L 2 + x 2 L 2 < 100 N x 2 ,
f bloc = 0 elsewhere .
f Airy = 2 J 1 ( arg ) arg ,
arg = 10 π N x x 2 L 2 + y 2 L 2 1 / 2 ,
I = 0 ( ρ 1 - ρ 2 ) / 2 exp ( iks 2 ) G ( s ) d s
0 < s 0 < s 1 < < s ν = ( ρ 1 - ρ 2 ) / 2 .
σ ( s ) = D i ( s - s i ) 3 + C i ( s - s i ) 2 + B i ( s - s i ) + A i ,
lim s s i   σ ( s ) = lim s s i   σ ( s ) = m i ; i = 1 , ,   ν - 1 ,
lim s 0   σ ( s ) = m 0 ; lim s s ν   σ ( s ) = m ν .
lim s 0   σ ( s ) = G 0 ; lim s s ν   σ ( s ) = G ν ,
lim s s i   σ ( s ) = lim s s i   σ ( s ) = G i , i = 1 , ,   ν - 1 .
h i = s i - s i - 1 ,
σ ( s ) = m i + 1 - m i 6 h i + 1   s 3 + m i s i + 1 - m i + 1 s i 2 h i + 1   s 2 + - m i s i + 1 2 - m i + 1 s i 2 2 h i + 1 + G i + 1 - G i h i + 1 - h i + 1 ( m i + 1 - m i ) 6 s + m i s i + 1 3 - m i + 1 s i 3 6 h i + 1 - G i + 1 - G i h i + 1   s i + h i + 1 ( m i + 1 - m i ) s i 6 + G i - h i + 1 2 6   m i .
σ ( s 0 ) = G ( s 0 ) ; σ ( s ν ) = G ( s ν ) .
I i μ = s i s i + 1 exp ( iks 2 ) s μ d s , μ = 0 , 1 ,   .
I i 0 = λ 2 C 2 s i + 1 λ - C 2 s i λ + i S 2 s i + 1 λ - S 2 s i λ ,
I i 1 = exp ( iks i + 1 2 ) - exp ( iks i 2 ) 2 ik ,
C ( x ) = 0 x cos π t 2 2 d t , S ( x ) = 0 x sin π t 2 2 d t .
I i μ = s i + 1 μ - 1 exp ( iks i + 1 2 ) - s i μ - 1 exp ( iks i 2 ) 2 ik - μ - 1 2 ik   I i μ - 2 .
I = i = 0 ν - 1 m i + 1 - m i 6 h i + 1   I i 3 + m i s i + 1 - m i + 1 s i 2 h i + 1   I i 2 + - m i s i + 1 2 - m i + 1 s i 2 2 h i + 1 + G i + 1 - G i h i + 1 - h i + 1 ( m i + 1 - m i ) 6 I i 1 + m i s i + 1 3 - m i + 1 s i 3 6 h i + 1 - G i + 1 - G i h i + 1   s i + h i + 1 ( m i + 1 - m i ) s i 6 + G i - h i + 1 2 6   m i I i 0 .
I = 0 ( ρ 1 - ρ 2 ) / 2 exp ( iks 2 ) G ( s ) d s .
I G 0 1 exp ( iks 2 ) G ( s ) d s
I σ 0 1 exp ( iks 2 ) σ ( s ) d s ,
I f 0 1 exp ( iks 2 ) f ( s ) d s = 1 2 ik 0 1 f ( s ) s d d s exp ( iks 2 ) d s = 1 2 ik f ( s ) s exp ( iks 2 ) | 0 1 - 1 2 ik 0 1 f ( s ) s - f ( s ) s 2 exp ( iks 2 ) d s = 1 2 ik   [ f ( 1 ) exp ( ik ) - f ( 0 ) ] + 1 2 ik 0 1 f ( s ) - sf ( s ) s 2 exp ( iks 2 ) d s .
I g 0 1 exp ( iks 2 ) g ( s ) d s = g ( 0 ) 0 1 exp ( iks 2 ) d s + 0 1 exp ( iks 2 ) [ g ( s ) - g ( 0 ) ] d s = g ( 0 ) π 2 k 1 / 2 C 2 k π 1 / 2 + iS 2 k π 1 / 2 + g ( 1 ) - g ( 0 ) 2 ik exp ( ik ) - g ( 0 ) 2 ik + 1 2 ik 0 1 g ( s ) - g ( 0 ) - sg ( s ) s 2 exp ( iks 2 ) d s ,
C 2 k π 1 / 2 + iS 2 k π 1 / 2
exp [ i ( π / 4 ) ] 2 - ( 1 + i ) 2 π k exp ( ik ) + O 1 k 3 / 2 .
g ( s ) = G ( s ) ,
f ( s ) = G ( s ) - σ ( s ) ,
max 0 s 1 | f ( j ) ( s ) |   C G 4 , h 4 - j ; j = 0 ,   1 ,   2 ,   3 ,
G 4 , = max 0 s 1 [ | G ( s ) | + | G ( 1 ) ( s ) | + | G ( 2 ) ( s ) | + | G ( 3 ) ( s ) | + | G ( 4 ) ( s ) | ] .
| I G - I σ | = | I G - σ | = 0 1   exp ( iks 2 ) [ G ( s ) - σ ( s ) ] d s G - σ 4 , h 4 ,
I G G ( 0 ) 2 π k 1 / 2 exp i   π 4 + O 1 k ,
I G - σ σ ( 0 ) - G ( 0 ) 2 ik + 1 2 ik 0 1 G ( s ) - σ ( s ) - s [ G ( s ) - σ ( s ) ] s 2   exp ( iks 2 ) d s .
p ( s ) = G ( s ) - σ ( s ) - s [ G ( s ) - σ ( s ) ] .
G ( s ) - σ ( s ) - s [ G ( s ) - σ ( s ) ] = [ σ ( t ( s ) ) - G ( t ( s ) ) ] t ( s ) s .
| G ( s ) - σ ( s ) - s [ G ( s ) - σ ( s ) ] | s 2 C G 4 , h 2 .
| G ( 0 ) - σ ( 0 ) |   C G 4 , h 3 ,
| I G - σ |   C G 4 , h 2 ( h + 1 ) 2 k .
| I G - σ | | I G | C G 4 , π G ( 0 ) h 2 ( h + 1 ) k ; k .

Metrics