## Abstract

In several optical systems, a specific Point Spread Function (PSF) needs to be generated. This can be achieved by shaping the complex field at the pupil. The Extended Nijboer-Zernike (ENZ) theory relates complex Zernike modes on the pupil directly to functions in the focal region. In this paper, we introduce a method to engineer a PSF using the ENZ theory. In particular, we present an optimization algorithm to design an extended depth of focus with high lateral resolution, while keeping the transmission of light high (over 60%). We also have demonstrated three outcomes of the algorithm using a Spatial Light Modulator (SLM).

© 2014 Optical Society of America

## 1. Introduction

Engineering an elongated focal spot by manipulating the incident wavefront has been done for a long time. It has been known that by putting an annulus in the lens pupil the depth of focus is increased [1]. In fact, it turned out that if the annulus is made extremely narrow, a diffraction-free Bessel beam is approximated [2]. Another method which has been proposed to extend the depth of focus is by dividing the lens pupil in rings and modulate the phase and amplitude of each ring [3]. These and other methods are discussed in more detail in [4–10].

In this article, a new method of designing pupil function is proposed. It relies on a result of the Extended Nijboer-Zernike theory, as has been suggested by [11]. By exploiting a one-toone correspondence between the Zernike polynomials which compose the pupil function and functions which compose the focal field [12], the number of degrees of freedom reduces to the number of polynomials with which the pupil is constructed, while the computation of the focal field is done by a quick matrix multiplication, rather than a time consuming diffraction integral. In particular, we present in this work an algorithm to find a pupil function which gives rise to an extended depth of focus for the case of a 0.4 numerical aperture, with low loss (less than 40%) of light intensity.

The relevance of such a focal field with an extended depth of focus is apparent from the many applications of Bessel beams, which range from optical tweezers to certain forms of microscopy to barcode scanners. The solution which is presented here however, has the advantage that light passes through the pupil in more areas than just the outer ring, which is a standard method to create a Bessel beam, and thereby achieving higher transmittance.

Finally, three outcomes of the algrotihm have been obtained experimentally using a phase-only Spatial Light Modulator (SLM) to shape the pupil field [13].

## 2. Theory

In this section we briefly discuss theoretical results required for the optimization algorithm that is later presented in this paper. First we give the definition of the complex Zernike polynomials, in which any pupil function can be decomposed. Then, we explain a result from ENZ-theory [12,14], which relates the pupil Zernike modes ${Z}_{n}^{m}$ to basic functions ${V}_{n}^{m}$ in the field in the focal region.

#### 2.1. Complex Zernike polynomials

Any complex function defined on the unit disk can be expanded in terms of complex Zernike polynomials

where the coefficients ${\beta}_{n}^{m}$ can be complex and where the complex Zernike polynomials are defined as where#### 2.2. Result from ENZ-theory

In [12] it is shown that there exists a one-to-one correspondence between the complex Zernike polynomials
${Z}_{n}^{m}(\rho ,\theta )$ which compose the pupil field *E*_{0}(*ρ*, *θ*) and functions
${V}_{n}^{\left|m\right|}(r,f)$ which compose the electric field in the focal region *E*(*r*, *ϕ*, *f*):

*ρ*,

*θ*) are polar coordinates, (

*r*,

*ϕ*,

*f*) cylindrical coordinates, and ${\beta}_{n}^{m}$ complex coefficients. An expression for the ${V}_{n}^{m}$-functions is given in [15] and is included in the Appendix.

The importance of this result lies in the fact that now the intensity distribution can be computed directly using the Zernike coefficients
${\beta}_{n}^{m}$ and the precalculated
${V}_{n}^{m}$ rather than using the input field *E*_{0}(*ρ*, *θ*) and the diffraction integral [16].

#### 2.3. The optimization algorithm

In order to find good solutions efficiently, we rely on Eq. (4). This equation states that once one has found *β* such that

*m*= 0 implies that we only consider circularly symmetric pupil functions. Since the ${V}_{n}^{m}$’s only need to be precalculated once, the intensity

*I*can computed this way significantly faster than using the diffraction integral.

In the algorithm, only the intensity distributions along the lateral axis (*x*-axis) and along the optical axis (*z*-axis) are considered. More specifically, let (0, 0) be the Gaussian focal point. Then *I _{x}* =

*I*(

*x*,

*z*= 0) is the intensity along the lateral axis in the focal plane and

*I*=

_{z}*I*(

*x*= 0,

*z*) is the intensity along the optical axis. Using a least squares method, the ${\beta}_{n}^{m}$ coefficients are found such that

*I*and

_{x}*I*approximate two target functions

_{z}*I*

_{x}_{,target}and

*I*

_{z}_{,target}as shown in Fig. 1. After that, the transmittance of the pupil function is optimized using the Nelder-Mead algorithm. The algorithm is explained in more detail in the Appendix.

In Fig. 2 the cross sections of the phase and amplitude of the obtained pupil function are shown, accompanied by the intensity profiles in the focal region along the *z*-axis and *x*-axis. The focal spot for this pupil is longer in the axial direction than in the aberration-free case, meaning that the depth of focus is increased. Since this pupil involves a quite structured amplitude illumination, we have made some approximations in such a way that the pupil could more easily be implemented in an experimental setup. This is treated in detail in the following section.

#### 2.4. Creating binary pupil functions

In order to simplify the experimental implementation, the pupil functions found by the algorithm (for example the one shown in Fig. 2) are made binary. Consider the ${\beta}_{n}^{m}$ shown in Table 1. According to theory these coefficients should give pupil functions which give rise to an extended depth of focus and a lateral resolution below the diffraction limit. To construct the binary pupil functions we first define

where ${Z}_{n}^{m}$ are the complex Zernike polynomials as defined in Eq. (2). The pupil function for ‘*Phase*’ (named so because only the phase is modulated) is obtained using Note that this is well-defined, because using only Zernike modes with

*m*= 0 means that

*Ẽ*

_{0}(

*ρ*) is real. The pupil functions for ‘

*Complex 1*’ and ‘

*Complex 2*’ (named so because both amplitude and phase are modulated) are obtained using

*t*= 0.57 is a good choice for ‘

*Complex 1*’ and

*t*= 0.55 for ‘

*Complex 2*’. The resulting pupil functions which will be assigned to the SLM are shown in Figs. 3(a)–3(c). Note that in the regions where the amplitude should be 0, a phase ramp is added so that light hitting that region of the SLM is tilted away. The corresponding theoretical intensity distributions in the focal region produced by these pupil functions are shown in Figs. 3(d)–3(f). For ‘

*Complex 2*’, the pupil function that gives the longest focal depth, 61% of the light is transmitted.

## 3. Experiment

In the experiment we used a Holoeye PLUTO liquid crystal Spatial Light Modulator (SLM) to assign the pattern to the laser beam. The specifications of the SLM are given in [17]. A schematic of the experimental setup is shown in Fig. 4. A He-Ne laser with a wavelength of 633nm is used as the source. Its light is coupled to a single-mode optical fiber after which it emerges as a point source (by approximation). The light is then collimated using a lens with a focal length of 80cm and cut off by an aperture. The collimated light passes through the beamsplitter, hits the SLM perpendicularly, and is reflected by the same beamsplitter. The plane of the SLM is conjugated with the plane of the entrance pupil of a 0.4 NA microscope objective using two lenses with a focal length of 30cm. The light in the focal region of the 0.4 NA microscope objective is then imaged onto a CCD camera using a 0.9 NA microscope objective and a tube lens. The 0.4 NA microscope objective can be translated back and forth using a piezo-stage, which allows us to scan different planes of the focal field. An interferometer is used to determine the exact position of the piezo-stage [18]. Several important details of the experimental setup are summarized in Table 2.

#### 3.1. Measurement method

First, we need to realize that not all the light that hits the SLM gets modulated by it due to interpixel space and a ‘wrong’ polarization. Thus in order to separate the modulated from the unmodulated beam an additional phase ramp is added on top of the images shown in Fig. 3 [19].

To view different focal planes, a lens phase is added to the image assigned to the SLM. The strength of this lens phase can be expressed as the Zernike defocus coefficient ${\alpha}_{2}^{0}$ (note that the corresponding Zernike polynomial is real and describes the phase only, as opposed to the coefficients and Zernike polynomials used in the optimization which are complex).

However, if we want to compare measurement results to the theoretical predictions, we first need to find out how the Zernike defocus coefficient
${\alpha}_{2}^{0}$ relates to the Rayleigh unit *R _{E}*, i.e. which value for
${\alpha}_{2}^{0}$ corresponds to a shift of 1

*R*. This has been done by scanning through focus without a phase mask, and observing for what value for ${\alpha}_{2}^{0}$ the intensity has dropped to 1/

_{E}*e*times the maximum intensity. This value of ${\alpha}_{2}^{0}$ corresponds to 1

*R*. In Fig. 5 the result of this measurement is shown.

_{E}#### 3.2. Results and discussion

Once the calibrations are done, we obtained the focal field intensity distributions for the three pupil functions (as in Table 1) and compared them with the theoretical predictions as shown in Fig. 6. In Fig. 7 the lateral width (along the *x*-axis) of the focal spot for the aberration-free case and for the pupil function ‘*Complex 2*’ are compared. Although the measurement results clearly indicate an extended depth of focus of a length similar to the prediction, there are still some discrepancies. To find out the possible causes, we compared the cases where the through focus scans were made by adding a lens phase, or by moving the 0.4 NA microscope objective. The comparison is shown in Fig. 8. It turns out that these measurements give somewhat different results, but it still does not account for the mismatch with the theory. A possible explanation may be that because the 0.9 NA objective used for imaging is infinity corrected, so whenever a plane other than the focal plane is imaged, aberrations may come into play. This is particularly relevant in our case since the focal depth is around 18 *R _{E}*. Nevertheless, the experimental results have demonstrated that pupil masks indeed create elongated focal spots and that the lateral resolution at the focal plane is below the diffraction limit.

## 4. Conclusion

In conclusion, we used the Extended Nijboer-Zernike theory to obtain a focused field with elongated focal length (up to 18 Rayleigh distances) with diffraction-limited spot size, while keeping the transmittance over 60%. Since this theory is based on functions that can be pre-calculated, the optimisation parameters are reduced to a limited number of Zernike coefficients that compose the pupil, allowing fast calculation of the focal field. We have found not only a way to create elongated focal spots, but also we have shown more generally that ENZ theory has the potential to be used for pupil engineering, as had been suggested by [11]. We also demonstrated experimentally three pupil functions generated by an SLM that produce elongated focal spots.

## Appendix A: An expression for the ${V}_{n}^{m}$ functions

An expression of ${V}_{n}^{m}$ is given in [15]:

*p*= (

*n*− |

*m*|)/2,

*q*= (

*n*+ |

*m*|)/2, and the minimum bound of

*l*in the summation is given as:

*l*∈ [max(0,

*k*−

*q*,

*p*−

*k*),

*k*+

*p*]. When

*m*= 0, we have:

*m*= 0:

*δ*is Kronecker’s delta.

## Appendix B: The optimization algorithm

## B.1. Formulating the problem

The goal is to construct a radially symmetric pupil function *E*_{0}(*ρ*) defined in the exit pupil (which is the unit disk) such that the intensity distribution in the focal field *I*(*x*, *y*, *z*) = |*E*(*x*, *y*, *z*)|^{2} has a sharp peak near (*x*, *y*) = (0, 0) for a long range along the optical *z*-axis. This means that the lateral resolution is high and the depth of focus long. In order to give a more precise formulation, we introduce the following simplifications and conventions:

- ★ We only consider pupil functions which can be expressed as a linear combination of
*N*Zernike polynomials of the form ${Z}_{n}^{0}$,*N*being finite: ${E}_{0}(\rho )={\sum}_{k=0}^{N}{\beta}_{2k}^{0}{Z}_{2k}^{0}$. The pupil function then depends on*ρ*only and therefore has circular symmetry, and the electric field*E*(*r*,*ϕ*,*z*) in the focal region does not depend on*ϕ*. - ★ We only consider real Zernike coefficients
*β*. This implies that*I*(*x*,*y*,*z*) is symmetric in*z*. - ★ The cross-sections
*I*(*x*,*z*) =*I*(*x*, 0) and*I*(*x*,*z*) =*I*(0,*z*) are a sufficiently good representation of*I*(*x*,*z*) in the entire positive*xz*-plane, hence we only consider the cross-sections. - ★ The transverse coordinates in the focal region are normalized by the Airy radius ${R}_{A}=\frac{0.61\lambda}{\text{NA}}$, and
*z*is normalized by the Rayleigh unit ${R}_{E}=\frac{n\lambda}{{\text{NA}}^{2}}$.*R*denotes the lateral width of the focal spot of a focused plane wave, while_{A}*R*denotes the length of the spot along the optical axis for this field._{E}

We are now able to define our decision variables:

- Define $\overrightarrow{\beta}={\left[{\beta}_{0}^{0},{\beta}_{2}^{0}\cdots {\beta}_{N}^{0}\right]}^{T}$.
- Define
*x*_{0}to be the location of the local minimum of*I*(*x*, 0) closest to the optical axis. - Define
*z*_{0}to be the location of the point with smallest value where*I*(0,*z*) < 0.9. - Define
*M*to be $M=\underset{x>{x}_{0}}{\text{max}}I(x,0)$.

Find *β⃗* such that:

## B.2. Finding solutions

In order to find good solutions efficiently, we rely on Eq. (4). This equation states that once one has found *β* such that

#### B.2.1. Precalculations

We start by defining grid points on *x* and *z* where we want to calculate the focal field intensity, and the points *ρ* where we want to find the pupil function

*V*

_{2}

*and ${Z}_{2k}^{0}$ for those points*

_{k}**X**,

**Z**and

**R**as

*I*(

*x*, 0),

*I*(0,

*z*) and |

*E*

_{0}|

^{2}are calculated and normalized according to

#### B.2.2. The algorithm

We start the algorithm by defining the target functions *I*_{target}(*x*, 0) and *I*_{target}(0, *z*) which roughly resemble our desired intensity distribution as in Fig. 1. This is shown in the flowchart in Fig. 9. Then we find *β⃗* which closely match these targets using the built-in Matlab algorithm
`lsqnonlin`. The merit function *f* is defined as:

`lsqnonlin`a random starting

*β⃗*

_{0}, it will apply the

*trust-region-reflective algorithm*using

*f*and returns a

*β⃗*

_{1}. We keep creating random starting

*β⃗*

_{0}until the output is satisfactory.

Having found an initial solution we proceed by improving *x*_{0}, *z*_{0}, *T* and *D* using another Matlab function,
`fminsearch`. Given a starting point, it applies the *Nelder-Mead algorithm*, a heuristic search method, to find minima. Of course, while improving *T*, we do not want *x*_{0} or *z*_{0} to worsen. Therefore we use the following scheme:

- Given a starting solution
*β⃗*_{1}, we calculate ${x}_{0}^{*}$ and ${z}_{0}^{*}$. - We define a function
*MT*(*β⃗*) which, given a certain*β⃗*, calculates*x*_{0},*z*_{0}, and*T*, and returns$$MT(\overrightarrow{\beta})=\{\begin{array}{ll}\frac{1}{T(\overrightarrow{\beta})}\hfill & \text{if}\hspace{0.17em}\hspace{0.17em}{x}_{0}(\overrightarrow{\beta})\le {x}_{0}^{*},{z}_{0}(\overrightarrow{\beta})\ge {z}_{0}^{*}\hspace{0.17em}\hspace{0.17em}\text{and}\hspace{0.17em}\hspace{0.17em}M<0.2\hfill \\ \infty \hfill & \text{otherwise}\hfill \end{array}$$So basically ${x}_{0}^{*}$ and ${z}_{0}^{*}$ are constraints for*x*_{0}and*z*_{0}. - We apply the Nelder-Mead algorithm to the function
*MT*(*β⃗*), using the starting solution*β*_{1}. It then returns a new*β⃗*, say*β⃗*_{2}, with improved*T*.

The algorithm is summarized in Fig. 9.

## Acknowledgments

A.P. Konijnenberg would like to thank J.J.B Teuwen for insightful discussions about the theoretical aspect of this research. L. Wei would like to thank Bernd Kleemann especially for discussions and supports. We would like to thank Axel Wiegmann for his contributions to the initial design of the experimental setup.

## References and links

**1. **G. B. Airy, “On the diffraction of an annular aperture,” Philos. Mag. **18**, 1–10 (1841).

**2. **J. Durnin and J. J. Miceli Jr., “Diffraction-free beams,” Phys. Rev. Lett. **58**, 1499–1501 (1987). [CrossRef] [PubMed]

**3. **Y. Xu, J. Singh, C. Sheppard, and N. Chen, “Ultra long high resolution beam by multi-zone rotationally symmetrical complex pupil filter,” Opt. Express **15**, 6409–6413 (2007). [CrossRef] [PubMed]

**4. **H. Wang, C. J. R. Sheppard, K. Ravi, S. T. Ho, and G. Vienne, “Fighting against diffraction: apodization and near field diffraction structures,” Laser Photonics Rev. **6**, 354–392 (2012). [CrossRef]

**5. **B. J. Thompson, “Diffraction by semitransparent and phase annuli,” J. Opt. Soc. Am. **55**, 145–148 (1965). [CrossRef]

**6. **G. Toraldo di Francia, “Super-gain antennas and optical resolving power,” Nuovo Cimento **9**, 426–438 (1952). [CrossRef]

**7. **H. Wang and F. Gan, “High focal depth with a pure-phase apodizer,” Appl. Opt. **40**, 5658–5662 (2001). [CrossRef]

**8. **C. J. R. Sheppard, “Synthesis of filters for specified axial properties,” J. Mod. Opt. **43**, 525–536 (1996). [CrossRef]

**9. **J. Ojeda-Castañeda and L. R. Berriel-Valdos, “Arbitrarily high focal depth with finite apertures,” Opt. Lett. **13**, 183–185 (1988). [CrossRef] [PubMed]

**10. **H. Wang, G. Yuan, W. Tan, L. Shi, and T. Chong, “Spot size and depth of focus in optical data storage system,” Opt. Eng. **46**, 065201 (2007). [CrossRef]

**11. **L. Wei, “Evaluation of Extended Nijboer-Zernike theory as a tool for complex point spread function calculation,” Master thesis (Delft University of Technology, 2012).

**12. **J. J. M. Braat, A. J. E. M. Janssen, P. Dirksen, and S. van Haver, “Assessment of optical systems by means of point spread functions,” Prog. Opt. **51**, 349–468 (2008).

**13. **C. Maurer, A. Jesacher, S. Bernet, and M. Ritsch-Marte, “What spatial light modulators can do for optical microscopy,” Laser Photonics Rev. **5**, 81–101 (2011). [CrossRef]

**14. **A. J. E. M. Janssen, “ENZ approach for the computation of optical point spread function,” J. Opt. Soc. Am. A **19**, 849–857 (2002). [CrossRef]

**15. **A. J. E. M. Janssen, J. J. M. Braat, and P. Dirksen, “On the computation of the Nijboer-Zernike aberration integrals at arbitrary defocus,” J. Mod. Opt. **51**, 687–703 (2004). [CrossRef]

**16. **S. van Haver and A. J. E. M. Janssen, “Advanced analytic treatment and efficient computation of the diffraction integrals in the Extended Nijboer-Zernike theory,” J. Eur. Opt. Soc. Rap. Publ. **8**, 13044 (2013). [CrossRef]

**17. ** Holoeye, “PLUTO Phase Only Spatial Light Modulator (Reflective),” http://holoeye.com/spatial-light-modulators/slm-pluto-phase-only/, June 24, 2013.

**18. **O. El Gawhary, A. Wiegmann, N. Kumar, S. F. Pereira, and H. P. Urbach, “Through-focus phase retrieval and its connection to the spatial correlation for propagating fields,” Opt. Express **21**, 5550–5560 (2013). [CrossRef] [PubMed]

**19. **H. Zhang, J. Xie, J. Liu, and Y. Wang, “Elimination of a zero-order beam induced by a pixelated spatial light modulator for holographic projection,” Appl. Opt. **48**, 5834–5841 (2009). [CrossRef] [PubMed]