We propose a novel design method for a circularly symmetric phase mask to extend the depth of focus. Using the free-form rational function as the solution space, we optimize the profile of the phase mask by analysis of the axial intensity distribution, which can be calculated efficiently by employing the fast Fourier transform algorithm. Numerical comparisons prove the resulting rational phase mask’s superiority to the existing quartic phase mask in intensity distribution and imaging performance.
© 2009 Optical Society of America
Conventional imaging systems often suffer from defocus error, and various optical elements have been proposed to extend the depth of focus (DOF). The apodizer [1, 2], which is an amplitude mask, alleviates the defocus error at the expense of the light transmission and results in a decreased contrast. The odd symmetric phase masks can achieve a large focal depth after digital image restoration without energy loss [3, 4, 5]. However, this approach introduces a phase modulation changing continuously; the resulting artifacts make the image restoration difficult and time-consuming . Also, this kind of odd symmetric phase mask is difficult and expensive to manufacture and test.
Circularly symmetric pure phase elements were also proposed in the literature, such as the logarithmic axicon (LA), the quartic phase mask (QPM), and the refractive and diffractive hybrid lens [7, 8, 9]. The last is fabricated through complicated laser writing on photoresist and subsequent reactive ion etching , while the former two can be more conveniently manufactured by use of traditional techniques such as diamond machining, and the cost can be reduced. Unlike the odd symmetric phase masks whose point-spread functions (PSF) are asymmetric, the PSFs of these circularly symmetric phase masks tend to be circularly symmetric as well. Additionally, the PSFs are also sharp and well localized for a large range of defocus without any shifts in the image sensor. Therefore, the post digital image restoration can be omitted and the circularly symmetric phase masks can satisfy some real-time all-optical applications to extend the DOF.
The LA was derived on the basis of geometrical optics: the geometrical law of energy conservation in  and Fermat’s principle in . The QPM was achieved from the Wigner distribution function in  and the stationary phase method in . They have similar performance in extending the DOF, and comparisons proved that the QPM is even superior to the LA . However, there still exist two main problems with them. First, numerical simulations indicate that the Fresnel diffraction integral shows prominent rapid oscillations of the axial optical intensity, especially near the front and the end of the focal region, which are caused by diffraction and interference effects . Second, near the front of the focal region the intensity distribution across a transverse plane is broad, monotonically decreasing away from the axis, and possesses no zeros. The two techniques also do not preserve uniformity of lateral resolution and energy flow along the focal range . To solve these problems, the apodized annular-aperture approach has been proposed and widely used [14, 15], although it brings the drawbacks of a decrease of signal-to-noise ratio (SNR) and a possible decrease in image resolution. Some other possible approaches, such as employing partially coherent or polychromatic light, have been suggested and analyzed [16, 17], but strict illumination conditions limit their application .
In fact, there is the potential to improve these poor performances by optimizing the profile of the phase mask. In this paper we describe the optimization design method of a circularly symmetric phase mask called the rational phase mask (RPM) and compare it with the well-known QPM derived from analytical methods. We find this numerical approach provides more design freedom and allows obtaining improved performance with wide applications. Section 2 gives the basic theory. Using the free-form rational function as the solution space, we present the optimization design method of the RPM by analysis of the axial intensity distribution in Section 3. Section 4 is the numerical comparisons of intensity distribution and imaging performance between the optimized RPM and the QPM. We conclude the work in Section 5.
We consider a circularly symmetric optical system with the generalized phase pupil function 19]19]4) can be rewritten as5) comes from wave theory and can be calculated efficiently by employing the fast Fourier transform algorithm. Therefore, by analyzing the result of Eq. (5) and with a heuristic search we may obtain a phase profile with the required axial intensity distribution.
To extend the DOF, we usually expect a rectangular axial intensity distribution with a mean value as high as possible in the focal region. However, we observe that from the Parseval theorem , the area under the axial intensity distribution curve is
3. OPTIMIZATION DESIGN
We may take the work of pursuing the best phase mask as the approximation of the unknown optimal phase function . Therefore the search space of the feasible solution should be as large as possible. Usually, a polynomial function is used to approximate an analytical function, but it often suffers from unsatisfactory convergence and a large number of undetermined coefficients, which means an increasing computational cost in the optimization procedure. The most advanced approximant would belong to a class of functions known as the rational function , which can obtain sufficient accuracy with a relatively small number of approximation coefficients and exhibits more flexible behavior than polynomials. Here we present the phase mask design method using the rational function, and the corresponding phase function for this RPM is described as 
Fig. 1 shows a typical axial intensity distribution curve resulting from Eq. (5). We assume the axial intensity oscillation is measured as the peak-to-valley value in the effective focal region and define a constant C as the ratio of the acceptable axial intensity oscillation to the maximum axial intensity value . The effective focal region can be defined as the region where there is no axial intensity value below . Then the merit function is to maximize under required constraints, which can be written as21], and the resulting RPM phase profiles are plotted in Fig. 2 .2. As can be seen, all the slopes of these phase profiles are comparable, which implies these masks may have a similar manufacturability. In Section 4, we will compare them in detail.
4. NUMERICAL RESULTS
4A. Intensity Distribution
From Eqs. (4, 2), the axial intensity distribution and the partial PSF near the optical axis for each of the three phase masks are calculated, respectively, as shown in Figs. 3, 4, 5 . Here, the radial coordinate r at the image plane is measured in units of the radius to the first zero of the CA Airy pattern for the in-focus diffraction-limited case, which is given by 3a and Fig. 4a, the axial intensity oscillations reach up to about 50% of the maximum value, which is the same with most values of α. The axial intensity uniformity for the RPM0.05 and RPM0.1 is improved visibly, which results from the small value of C. If we choose a larger C (e.g., ), the optimized RPMs tend to converge to the corresponding QPMs. Particularly, for the QPM0.25 in Fig. 5a, the oscillations are relatively small, which reach about 13% of the maximum value. However, it is still possible to smooth the axial intensity distribution by optimizing the profile of the RPM. By setting , we achieved the RPM0.25 with an extremely flat top in the focal region.
Simultaneously, the axial intensity uniformity is achieved at a cost of some decrease of the effective focal range compared with the QPM. This phenomenon was also observed when using the apodized annular-aperture approach , where was defined as the transition region between the top and the bottom of the axial intensity distribution curve, and a larger resulted in a better smoothness. In Fig. 5a we also plot the axial intensity distribution of the CA. The extending of the DOF is obvious and Eq. (6) is verified, that is, we cannot achieve a large DOF and a high SNR at the same time.
We also evaluate the off-axis quality of the light. From the partial PSF plotted in Figs. 3, 4, the QPMs’ PSF side lobes or local maxima show some discontinuity along the focal region. Usually a larger central peak corresponds to a smaller second maximum, and a larger third maximum, and so on. The similarity between the PSFs at different image planes is destroyed. For the partial PSF of the RPMs, the local maxima are more continuous along the focal region and the lateral intensity uniformity is preserved better. But for the RPM0.25 in Fig. 5b, this improvement is not so evident. The description above may hold true except for the front of the focal region where the central spot size increases and the local maxima disappear. In the paraxial region of the partial PSF, we calculate the light energy,6 . For the QPM0.1, the nonuniform character of the energy flow along the focal region consists of oscillations and a decline. For the RPM0.1 the main trend is only the decline, and the value between the maximum and the minimum is smaller. The dot curves in Fig. 6 correspond to the results calculated from letting , and the uniformity is improved significantly for both RPM0.1 and QPM0.1. This uniformity of the energy flow may be quite useful in applications such as energy concentration, laser machining, and so on.
4B. Imaging Performance
We have considered the axial and paraxial irradiance, but this common analysis is incomplete for imaging applications that require full information about the PSF or the equivalent optical transfer function (OTF). Next, taking RPM0.1 and QPM0.1, for instance, we concern ourselves with the radius to the first zero of the PSF, which we label as the cutoff lateral position , and we may expect that would be approximately constant along the focal region. However, Fig. 7a illustrates that near the front of the focal region the PSF profiles of QPM0.1 are apt to decrease monotonically and possess no zeros, as pointed out in , while for the RPM0.1, its PSF profiles always possess zeros. Elsewhere, the cutoff lateral positions of the RPM0.1 and QPM0.1 tend to converge to the same value with tiny fluctuations. It is interesting to point out that does not mean the RPM0.1 or QPM0.1 has a higher resolution than the in-focus CA, because the local maxima contain energy no less than the central spot as shown in Fig. 7b. For instance, at the central peak contains only about 1% of the total energy for both the RPM0.1 and QPM0.1; for comparison, at the central peak contains about 16% of the total energy for both the RPM0.1 and QPM0.1. This can also be seen from the modulation transfer functions (MTF) plotted in Fig. 8 , where the cutoff frequencies are very low for . For and , the cutoff frequencies of the RPM0.1 are higher. And for , the cutoff frequencies are almost the same for the two masks. It is notable that all the cutoff frequencies are reduced compared with the those of the in-focus CAs, as shown in Fig. 10a. For applications using pixelated detectors, these masks may be usefully employed as antialiasing filters .
We also compare the simulated imaging performance of the two masks. Figure 9 shows the raw images of a simulated spoke target produced by applying the corresponding OTFs, and there is no other digital image processing. At , , and , the imaging performance of the RPM0.1 is better than that of the QPM0.1, and more details can be seen. At and , the results are comparable, while at the high-frequency components in the center disappear first for the QPM0.1, and at the image degradation becomes serious for both masks, and most details are lost. For comparison, the simulated images of the CA are given in Fig. 10b and the high-frequency cutoff and contrast reversal become evident at and . The RPM0.1 has extended the DOF to a larger range than the QPM0.1, but they both have a contrast decrease resulting from their low MTFs in the focal region. The simulation above is under the noise-free condition, and in future research, the effects of noise should be taken into consideration.
In conclusion, we have proposed a novel design method of the circularly symmetric phase mask to extend the DOF by analysis of the axial intensity distribution, which can be calculated efficiently by employing the fast Fourier transform algorithm. The simple constraints on the axial intensity oscillations result in improved uniformity of the axial intensity, energy flow, lateral resolution, and also improved imaging performance along the focal region. With the free-form rational function as the solution space, the numerical optimization method may present more design flexibility for various applications by more complex statistical analysis of the axial intensity distribution.
1. J. Ojeda-Castenada, P. Andres, and A. Diaz, “Annular apodizers for low sensitivity to defocus and to spherical aberration,” Opt. Lett. 11, 487–489 (1986). [CrossRef]
2. J. Ojeda-Castenada, E. Tepichin, and A. Diaz, “Arbitrary high focal depth with a quasi-optimum real and positive transmittance apodizer,” Appl. Opt. 28, 2666–2670 (1989). [CrossRef]
4. Q. Yang, L. Liu, and J. Sun, “Optimized phase pupil masks for extended depth of field,” Opt. Commun. 272, 56–66 (2007). [CrossRef]
5. S. Bagheri, P. Silveira, and D. Farias, “Analytical optical solution of the extension of the depth of field using cubic-phase wavefront coding. Part II. Design and optimization of the cubic phase,” J. Opt. Soc. Am. A 25, 1064–1074 (2008). [CrossRef]
6. S. Förster, H. Gross, F. Höller, and L. Höring, “Extended depth of focus as a process of pupil manipulation,” Proc. SPIE 5962, 596207 (2005). [CrossRef]
8. D. Zalvidea and E. E. Sicre, “Phase pupil functions for focal depth enhancement derived from a Wigner distribution function,” Appl. Opt. 37, 3623–3627 (1998). [CrossRef]
9. Z. Liu, A. Flores, M. R. Wang, and J. J. Yang, “Diffractive infrared lens with extended depth of focus,” Opt. Eng. 46, 018002 (2007). [CrossRef]
10. W. Chi and N. George, “Electric imaging using a logarithmic asphere,” Opt. Lett. 26, 875–877 (2001). [CrossRef]
12. D. Zalvidea, C. Colautti, and E. E. Sicre, “Quality parameters analysis of optical imaging systems with enhanced focal depth using the Wigner distribution function,” J. Opt. Soc. Am. A 17, 867–873 (2000). [CrossRef]
13. L. R. Staronski, J. Sochacki, Z. Jaroszewicz, and A. Kołodziejczyk, “Lateral distribution and flow of energy in uniform-intensity axicons,” J. Opt. Soc. Am. A 9, 2091–2094 (1992). [CrossRef]
14. J. Sochacki, Z. Jaroszewicz, L. R. Staronski, and A. Kołodziejczyk, “Annular-aperture logarithmic axicons,” J. Opt. Soc. Am. A 10, 1765–1768 (1993). [CrossRef]
15. Z. Jaroszewicz, J. Sochacki, A. Kołodziejczyk, and L. R. Staronski, “Apodized annular-aperture logarithmic axicon: smoothness and uniformity of intensity distributions,” Opt. Lett. 18, 1893–1895 (1993). [CrossRef] [PubMed]
16. S. Y. Popov and A. T. Friberg “Design of diffractive axicons for partially coherent light,” Opt. Lett. 23, 1639–1641 (1998). [CrossRef]
17. Z. Jaroszewicz, J. F. Roman Dopazo, and C. Gomez Reino, “Uniformization of the axial intensity of diffraction axicons by polychromatic illumination,” Appl. Opt. 35, 1025–1031 (1996). [CrossRef] [PubMed]
18. S. Y. Popov and A. T. Friberg, “Apodization of generalized axicons to produce uniform axial line images,” Pure Appl. Opt. 7, 537–548 (1998). [CrossRef]
19. J. W. Goodman, Introduction to Fourier Optics, 3rd ed. (Roberts, 2005).
20. J. H. Mathews and K. D. Fink, Numerical Methods Using MATLAB, 4th ed. (Prentice Hall, 2004).
22. S. Mezouari, G. Muyo, and A. R. Harvey, “Circularly symmetric phase filters for control of primary third-order aberrations: coma and astigmatism,” J. Opt. Soc. Am. A 23, 1058–1062 (2006). [CrossRef]