## Abstract

The communication modes, which constitute a convenient method for the propagation and information analysis of optical fields, are formulated in the generalized axicon geometry. The transmitting region is the axicon’s annular aperture, and the observation domain is the optical axis containing the focal line segment. We show that in rotational symmetry one may employ the prolate spheroidal wave functions to represent the communication modes. Further, in usual circumstances the modes can be approximated by quadratic waves in the aperture domain and by sinc functions in the image domain. Both the exact communication modes and the approximate technique are confirmed numerically, with linear axicons as examples.

©2004 Optical Society of America

## 1. Introduction

Axicons are optical elements that produce extended line images along the optical axis [1–3]. Due to the narrow focal cross-section and the long depth of field, axicons in various forms have acquired a large number of practical applications [4–8]. The communication modes approach [9–12] is a physically insightful method for (back-)propagating optical waves between two domains and for analyzing the resolution and information content of fields and images. The communication modes have been employed in information transfer between geometries such as two rectangular volumes or planar apertures [11], a pair of parallel (transverse) lines [10,13], and a transverse and an axial line [14, 15]. Characteristic to all these systems are that the geometries are one-dimensional or factorable in Cartesian coordinates and that the communication modes are prolate spheroidal wave functions (PSWFs) [12, 16, 17].

In this paper we consider the communication modes theory in the customary axicon geometry: from a circular or annular transverse aperture to an on-axis image segment. We show that in the case of rotational symmetry, the communication modes can be obtained in terms of the PSWFs. When the number of modes (degrees of freedom) associated with the axicon system is large, the eigenfunctions take on the form of quadratic waves across the aperture plane and sinc functions along the image line. The exact and approximate modes methods are confirmed by numerical simulations.

## 2. Theory and numerical results

The geometry of axicon imaging is shown in Fig. 1. The inner and outer radii of the annular aperture are *R*
_{1} and *R*
_{2}, and the image segment is formed on the *z* axis within the observation region ranging from *d*
_{1} to *d*
_{2}. Consider a field *E*
_{in}(*ρ*), where *ρ* is the radial coordinate in the aperture plane. The field *E*
_{in}(*ρ*) contains both the distribution of the original incident beam and the effect of a phase element inserted into the aperture. The field is assumed monochromatic with time-dependence exp(-*iωt*) (which is suppressed), where *ω*=*ck* is the angular frequency and *k*=2*π*/*λ* is the wave number (*c* is the vacuum speed of light and *λ* is the wavelength).

In cylindrical coordinates and for rotationally symmetric fields, the on-axis optical field produced by the axicon can be found from the Fresnel diffraction integral

On making the change of variables

Eq. (1) becomes

where *Ẽ*
_{in}(*u*)=*E*
_{in}(*ρ*) and *U*=(${R}_{2}^{2}$-${R}_{1}^{2}$)/2. Hence the Green function in Eq. (3) is

i.e., *G*(*z*,*u*) serves as the propagator between the aperture field *Ẽ*
_{in}(*u*) and the image field *E*
_{out}(*z*).

The eigenvalues *g*_{n}
and eigenfunctions *a*_{n}
(*u*) in the transmitting, annular-aperture domain are obtained from the integral equation [11, 13]

where the kernel is

The integral in Eq. (3), and the analysis performed below to find the eigenfunctions, is very similar to the integral and analysis in Appendix B of Ref [14], where a more detailed derivation can be found. Using the definitions ${D}_{\mathrm{sum}}=\frac{1}{2}(1/{d}_{1}+1/{d}_{2})$, ${D}_{\mathrm{diff}}=\frac{1}{2}(1/{d}_{1}-1/{d}_{2})$, and *a*_{n}
(*u*)=exp(-*ikD*
_{sum}
*u*/2)*α*_{n}
(*u*), inserting Eq. (4) into Eqs. (5) and (6), and performing the integration in *z* leads to the expression

The solutions *α*_{n}
(*u*) to this integral equation are the prolate spheroidal wave functions (PSWFs) of scale *U* and bandwidth *kD*
_{diff}/2, while coefficients |*g*_{n}
| are proportional to the square root of the corresponding eigenvalues [12]. The number of degrees of freedom [11] for this geometry is

Hence *N* increases as the axicon’s aperture area and the length of the image section are increased, or the wavelength is decreased.

The eigenfunctions *b*_{n}
(*z*) in the receiving domain, and the phases of the eigenvalues, can be determined from the projection [11]

More specifically, inserting Eq. (4) and the functions *a*_{n}
(*u*) into Eq. (9), and performing the integration using the known, limited Fourier transform of the PSWFs [12], we find the full solutions:

where *α*_{n}
(*u*) are the PSWFs of scale *U* and bandwidth *kD*
_{diff}/2, *β*_{n}
(1/*z* - *D*
_{sum}) are the PSWFs of scale *D*
_{diff} and bandwidth *kU*/2, and *λ*_{n}
are the (real) eigenvalues of the PSWFs [12]. The functions *b*_{n}
(*z*) can also be found from an integral equation with kernel *J*(*z*, *z′*)=${\int}_{\mathit{-}U}^{U}$
*G*(*z*,*u*)*G**(*z′*,*u*)d*u*, yielding identical results.

The communication modes are now the two sets of eigenfunctions *a*_{n}
(*u*) (of *K*) and *b*_{n}
(*z*) (of *J*), and the coefficients *g*_{n}
, where |*g*_{n}
|^{2} are the eigenvalues of both *K*(*u*,*u′*) and *J*(*z*, *z′*). The magnitudes of the coupling coefficients *g*_{n}
are very close to unity for modes up to *N*, and beyond that they decrease rapidly. The Green function can then be expanded as [13]

and on substituting Eq. (13) into Eqs. (3) and (4) the output field is found as

where

i.e., the coefficients *A*_{n}
are the projections of the aperture field onto the eigenfunctions of the transmitting domain. The bi-orthogonal expansion in Eq. (13) is sometimes also referred to as an analytical singular-value decomposition (SVD) of the Green function.

As a physical application we may consider the on-axis resolution. It is found from the eigenfunctions *b*_{n}
(*z*), in analogy with Ref. [14]. The most rapidly oscillating function which gives a significant contribution to *E*
_{out}(*z*) is number N, for which the distance between two consecutive zeros is Δ*s*=2*D*
_{diff}/*N*=2*λ*/(${R}_{2}^{2}$-${R}_{1}^{2}$) where *s*=1/*z*-*D*
_{sum}. The separation on the z-axis is Δ*z*=*z*
^{2}Δ*s*/(1-*z*Δ*s*) [14], and the highest possible spatial frequency therefore is

The axial resolution depends on *z*, and is higher close to the aperture.

Numerical propagation is performed using Eqs. (14) and (15). The PSWFs are easily generated in Matlab [13]. As an example, we consider a linear axicon, which is the best known and most widely used representative of the axicon family. In uniform illumination, it gives rise to the aperture field *E*
_{in}(*ρ*)=exp(-*ikγρ*). The coefficient is *γ*=0.025, the aperture radii are *R*
_{1}=5 mm and *R*
_{2}=10 mm, and the observation domain on the *z* axis extends from *d*
_{1}=150 mm to *d*
_{2}=450 mm. The wavelength is *λ*=600 nm, and thus *N*≈278 from Eq. (8). The result is shown in Fig. 2. As expected, the on-axis intensity *I*
_{cm}(*z*)∝|*E*
_{out}(*z*)|^{2} is confined to the region between *z*=*R*
_{1}/*γ*=200 mm and *z*=*R*
_{2}/*γ*=400 mm. The oscillations are due to diffraction at the inner and outer aperture edges. The intensity distribution is the same as that obtained from a direct evaluation of the Fresnel diffraction integral, within a relative error of 1·10^{-4}, thus confirming the validity of the communication modes theory developed above.

## 3. Approximate solution

It is known that if the kernel of an integral equation of the type of Eq. (5) is stationary, and its effective width is much smaller than the integration domain, the eigenfunctions reduce to harmonic exponential functions and the eigenvalues are samples from the power spectrum [18]. Since the kernel of Eq. (7) depends only on *u′*-*u*, the condition for this approximation is that *U* is much larger than the kernel’s extent, i.e., that *U*≫2*π*/*kD*
_{diff}=*λ*/*D*
_{diff}. Using the number of degrees of freedom *N*=2*UD*
_{diff}/*λ*, the requirement simplifies to

Consequently, if there is a large number of communication modes, the integration limits of Eq. (7) may be extended to infinity and the orthonormal solutions are [18]:

That the harmonic exponential is really a solution is easily verified by inserting Eq. (18) into Eq. (7) with the extended integration limits. The factor *π*/*U* in the exponent follows from the orthogonality condition, and the factor $1/\sqrt{2U}$ from the normalization. The corresponding approximate eigenfunctions *b*_{n}
(*z*) may be found by inserting Eqs. (4) and (18) into the projection formula in Eq. (9), yielding

$$\times {\int}_{-U}^{U}\mathrm{exp}\left\{i\left[\frac{kU}{2}(\frac{1}{z}-{D}_{\mathrm{sum}})-(n-\frac{N}{2})\pi \right]\frac{u}{U}\right\}du.$$

Since the integrand is not limited, the integration limits in Eq. (20) can not be extended. However, the integral is readily evaluated and the result is

This is all the information required to perform both forward and backward propagation. The functions *a*_{n}
(*u*) have quadratic phase factors in *ρ* and, physically, they represent spherical waves converging towards different points on the *z* axis, where each of them yields a sinc field distribution *b*_{n}
(*z*). We note that the functions *b*_{n}
(*z*) no longer are strictly orthogonal, due to the approximation. The longitudinal resolution can be determined by the full width Δ*s* of the sinc functions, found from the first zeros as *kU*Δ*s*/2=*π*. The solution Δ*s*=2*λ*/(${R}_{2}^{2}$-${R}_{1}^{2}$) is the same as for the exact method, and thus the expression in Eq. (16) applies also to the approximate method.

Numerically, the propagation is readily performed using Matlab. The on-axis intensity distribution, with the same parameters as in Fig. 2, is shown in Fig. 3. It contains two graphs, one from the exact mode theory and one from the approximate theory, but the profiles are so alike that the graphs can hardly be told from each other. For *N=*278 the approximation is almost identical to the exact solution. For degrees of freedom *N*≤50 the differences become better visible, but the two methods yield fairly similar results even for as few as 5 modes.

The error of the approximation increases when the number of modes decreases. We define the relative error *ε* of the axial intensity as

where *I*
_{cm}(*z*) is the on-axis intensity according to the full communication mode theory and *I*
_{appr}(*z*) is the approximation. Keeping the axicon geometry and phase function the same, we may vary *N* by changing the wavelength λ. The results are shown in Fig. 4, where the relative error is plotted as a function of the number of degrees of freedom. There are oscillations, but the error steadily decreases on average as *N* increases.

## 4. Conclusions

We have developed the theory of communication modes for axicon geometry, in which light incident on a circular or annular aperture is focused by the (generalized) axicon onto a line segment along the optical axis. A change of variables allows us to find the communication modes that are expressible in terms of the prolate spheroidal wave functions, both across the aperture and on the optical axis. We also introduce an approximation, where the modes in the axicon aperture are complex exponentials in the transformed variable and those on the optical axis are sinc functions. The approximation is valid for axicon systems with a large number of degrees of freedom. Numerical simulations for a linear axicon confirm the accuracy of both the exact and the approximate communication modes methods. While the numerical propagation is not as efficient as e.g. the fast Fourier transform (FFT) method, the communication modes provide physical insight such as resolution analysis, information content of the optical field, and propagation in the presence of noise.

## Acknowledgments

This work was sponsored by the Swedish Research Council. Anna Burvall thanks the Göran Gustafsson Foundation for financial support.

## References and links

**1. **J.H. McLeod, “The axicon: a new type of optical element,” J. Opt. Soc. Am. **44**, 592–597 (1954). [CrossRef]

**2. **L.M. Soroko, *Meso-Optics — Foundations and Applications* (World Scientific, Singapore, 1996), Chap. 2 and references therein.

**3. **Z. Jaroszewicz, *Axicons: Design and Propagation Properties*, Research & Development Treatises, Vol. 5 (SPIE Polish Chapter, Warsaw, 1997).

**4. **J.A. Davis, E. Carcole, and D.M. Cottrell, “Range-finding by triangulation with nondifracting beams,” Appl. Opt. **35**, 2159–2161 (1996). [CrossRef] [PubMed]

**5. **G. Haüsler and W. Heckel, “Light sectioning with large depth and high resolution,” Appl. Opt. **27**, 5165–5169 (1988). [CrossRef] [PubMed]

**6. **V.E. Peet and R.V. Tusbin, “Third-harmonic generation and multiphoton ionization in Bessel beams,” Phys. Rev. A **56**, 1613–1620 (1997). [CrossRef]

**7. **J. Arlt, T. Hitomi, and K. Dholakia, “Atom guiding along Laguerre-Gaussian and Bessel light beams,” Appl. Phys. B **71**, 549–554 (2000). [CrossRef]

**8. **Z. Ding, H. Ren, Y. Zhao, J.S. Nelson, and Z. Chen, “High resolution optical coherence tomography over a large depth range with an axicon lens,” Opt. Lett. **27**, 243–245 (2002). [CrossRef]

**9. **F. Gori and R. Grella, “The converging prolate spheroidal functions and their use in Fresnel optics,” Opt. Commun. **45**, 5–10 (1983). [CrossRef]

**10. **R. Pierri and F. Soldovieri, “On the information content of the radiated fields in the near zone over bounded domains,” Inverse Problems **14**, 321–337 (1998). [CrossRef]

**11. **D.A.B. Miller, “Communicating with waves between volumes: evaluating orthogonal spatial channels and limits on coupling strengths,” Appl. Opt. **39**, 1681–1699 (2000). [CrossRef]

**12. **B.R. Frieden, “Evaluation, design and extrapolation methods for optical signals, based on use of the prolate functions,” in *Progress in Optics*, Vol. IX, E. Wolf, editor (North-Holland, Amsterdam, 1971), pp. 311–407. [CrossRef]

**13. **A. Thaning, P. Martinsson, M. Karelin, and A.T. Friberg, “Limits of diffractive optics by communication modes,” J. Opt. A: Pure Appl. Opt. **5**153–158 (2003). [CrossRef]

**14. **R. Pierri, A. Liseno, F. Soldovieri, and R. Solimene, “In-depth resolution for a strip source in the Fresnel zone,” J. Opt. Soc. Am. A **18**, 352–359 (2001). [CrossRef]

**15. **R. Pierri, A. Liseno, R. Solimene, and F. Tartaglione, “In-depth resolution from multifrequency Born fields scattered by a dielectric strip in the Fresnel zone,” J. Opt. Soc. Am. A **19**, 1234–1238 (2002). [CrossRef]

**16. **D. Slepian and H.O. Pollak, “Prolate spheroidal wave functions, Fourier analysis and uncertainty - I,” Bell Syst. Techn. J. **40**, 43–63 (1961).

**17. **H.J. Landau and H.O. Pollak, “Prolate spheroidal wave functions, Fourier analysis and uncertainty - II,” Bell Syst. Techn. J. **40**, 65–84 (1961).

**18. **B. Saleh, *Photoelectron Statistics* (Springer-Verlag, Heidelberg, 1978), *chap. 2.4.4.*