Point spread function characteristics analysis of the wavefront coding system

Wenzi Zhang; Zi Ye; Tingyu Zhao; Yanping Chen; Feihong Yu

doi:10.1364/OE.15.001543

1. Introduction

Most of the previous imaging characteristics analysis of the wavefront coding system has been carried out within the frequency domain. The optical transfer function (OTF), mainly the modulation transfer function (MTF), is utilized to analyze the sensitivities to the different aberrations for wavefront coding systems with all kinds of pupil phase masks [1–3]. There are only a few analyses within the spatial domain carried out, and only the on-axis properties are taken into account in most of these analyses, which are due to the implication of the Strehl ratio [4–7]. Besides, phase-space representations, such as the Wigner distribution function, ambiguity function and fractional Fourier transforms, have been employed to analyze the imaging characteristics of the wavefront coding system [8].

In order to get a deeper insight into the imaging characteristics of the wavefront coding system within the spatial domain, the approximate expression of PSF in the presence of defocus aberration is derived for the system with a cubic phase mask. With this expression, the PSF’s boundaries, oscillations, and sensitivities to defocus, astigmatism and coma can be analyzed. Some characteristics unexplained in the previous study are analyzed with approximated expressions in the spatial domain for the first time.

The rest of the paper is organized as follows. In section 2 the approximate expression is derived. Based on this expression, the imaging characteristics are analyzed in section 3. The conclusions are drawn in section 4.

2. Derivation of the PSF

The following analysis is restricted to the case of one-dimensional (1D) optical system. The defocused PSF h(x, W₂₀) of an incoherent imaging system is given by Eq. (1),

h (x, W_{20}) = {∣ \int_{- \infty}^{+ \infty} q (u) \exp (j k W_{20} u^{2} - j 2 πxu) d u ∣}^{2},

where W₂₀ is the traditional defocus aberration constant in unit of wavelength, u is the normalized pupil coordinate, q(u) is the normalized pupil function, and k is the wavenumber 2π/λ. For the wavefront coding system with a cubic phase mask and a bounded aperture, the normalized pupil function can be represented as

q (u) = {\begin{matrix} \frac{1}{\sqrt{2}} \exp (j α u^{3}) ∣ u ∣ \leq 1, \\ 0 ∣ u ∣ > 1 \end{matrix}

where α is the normalized parameter for cubic phase mask in unit of radian. Inserting Eq. (2) into Eq. (1) enables us to rewrite Eq. (1) in the form as shown in Eq. (3),

h (x, W_{20}) = \frac{1}{2} {∣ \int_{- 1}^{+ 1} \exp (j α u^{3} + j k W_{20} u^{2} - j 2 πxu) d u ∣}^{2} .

With the stationary phase method [9], Eq. (3) can be approximated by Eq. (4). Here 3α>2|kW₂₀| is assumed. Please see Appendix A for a detailed derivation of this result.

h (x, W_{20}) = {\begin{matrix} 0 x \leq - \frac{{(k W_{20})}^{2}}{6 πα} \\ \frac{π}{\sqrt{{(k W_{20})}^{2} + 6 π α x}} {1 + \sin [\frac{4 {({(k W_{20})}^{2} + 6 π α x)}^{\frac{3}{2}}}{27 α^{2}}]} - \frac{{(k W_{20})}^{2}}{6 π α} < x \leq \frac{3 α - 2 ∣ k W_{20} ∣}{2 π} \\ \frac{π}{2 \sqrt{{(k W_{20})}^{2} + 6 π α x}} \frac{3 α - 2 ∣ k W_{20} ∣}{2 π} < x \leq \frac{3 α + 2 ∣ k W_{20} ∣}{2 π} \\ 0 x > \frac{3 α + 2 ∣ k W_{20} ∣}{2 π} \end{matrix}

PSFs calculated from Eq. (1) with Fast Fourier Transform (FFT) and PSFs calculated from Eq. (4) for α=30π and W₂₀=0, 2.5λ, 5λ are shown in Figs. 1(a), 1(b) and 1(c), respectively. The boundaries in Eq. (4) and the absolute errors between the FFT and approximate PSFs are also shown in Fig. 1. Though there is a big range in the middle part of approximate PSF with very small absolute errors, large absolute errors occur when the reduced spatial coordinate gets closer to the three boundaries in Eq. (4). This phenomenon can be referred to the stationary phase method and the same problem exists in Dowski and Cathey’s approximate expression of OTF [1], which will not affect the imaging characteristics analysis.

Compared to the FFT approach, the approximate expression needs less sample points calculated to produce an accurate plot of the PSF in the middle part, which can be used to illustrate the PSF characteristics of the wavefront coding system. For a sample point number of N in the sampled pupil range [-L, L], N equally spaced points in the reduced spatial range [-N/(4L), N/(4L)] is calculated within the FFT approach. If d stands for the length of the desired reduced spatial range to plot the PSF, the number of calculated sample points that can be used is m=4dL. In order to acquire an accurate plot of the PSF in the middle part, the sampled pupil range should be enlarged, which leads to an increase of the sample point number N to maintain the FFT accuracy.

Fig. 1. Comparison of the FFT and approximate PSFs in the presence of defocus aberration for α=30π. The defocus values in (a), (b) and (c) are W₂₀=0, 2.5λ and 5λ, respectively.

Download Full Size | PDF

3. PSF characteristics analysis

3.1 Boundaries

Fig. 2. PSF of an unbounded aperture for α=30π and W₂₀=2.5λ.

Download Full Size | PDF

For an unbounded aperture, the defocused PSF of the wavefront coding system with a cubic phase mask can be represented as Eq. (5) and a substitution u= t - kW₂₀/(3α) can be applied,

h (x, W_{20}) = \frac{1}{2} {∣ \int_{- \infty}^{+ \infty} \exp (j α u^{3} + j k W_{20} u^{2} - j 2 π x u) d u ∣}^{2}

where Ai(x) is the Airy function [10]. The PSF of an unbounded aperture system for α=30π and W₂₀=2.5λ. is shown in Fig. 2. A boundary on the left near the origin can be found for both bounded and unbounded aperture, while the boundary on the right only exists for the bounded aperture. The width of the PSF area of the wavefront coding system with a cubic phase mask and a bounded aperture can be approximated by the subtraction between the right and left boundary position in Eq. (4),

D = \frac{3 α + 2 ∣ k W_{20} ∣}{2 π} - [- \frac{{(k W_{20})}^{2}}{6 πα}] = \frac{{(3 α + ∣ k W_{20} ∣)}^{2}}{3 π α f_{0}} \approx \frac{3 α}{π f_{0}} (for small W_{20}) .

where f₀ is the diffraction-limited cutoff frequency of the incoherent imaging system. For the working wavelength λ, f₀ is equal to L/(λd_i), where L is the aperture diameter, and d_i is the distance from the pupil to the diffraction-limited image plane [11]. About 99.79%, 99.77% and 99.68% of the PSF energy are constrained in the areas with the width determined by Eq. (6) for Figs. 1(a), 1(b) and 1(c), respectively. Thus Eq. (6) can be used to evaluate the sharpness of the wavefront-coded intermediate image with adequate accuracies, it is seen that a larger α and a smaller f₀ result in a broader PSF, i.e., a more blurry intermediate image. For an F/5 incoherent imaging system with a working wavelength of 587.6 nm, the diffraction-limited cutoff frequency is f₀ ≈ 340.4 cyc/mm and the PSF area width is D ≈ 0.26 mm when a cubic phase mask of α=30π is applied. Comparing to the width of the PSF area of a traditional imaging system which can be calculated according to 1/f₀, it can be found that the PSF area of the wavefront coding system with a cubic phase mask is much larger and the intermediate image without decoding is deeply blurred. Besides, Eq. (6) can help to prepare for the experiment of PSF measurements [12], which is a key step to construct the kernel filter for the filtering process in the wavefront coding technique.

The boundary limit in the spatial domain always leads to the bandwidth limit in the frequency domain. When checking the stationary points and the boundary limits while applying the stationary phase method to the Fourier transform of Eq. (4), the ratio of the bandwidth between the wavefront coding system with a cubic phase mask and the traditional system can be found to be 1 - |kW₂₀|/(3α). The derivation is not included in the paper because a much easier way of analysis will be presented in the paper. When a small defocus occurs, this ratio is approaching unity. However when a large defocus occurs, the decrease of the bandwidth becomes obvious. The same result can be found in a previous paper [13].

3.2 Oscillations

Fig. 3. PSFs of a 2D system for α=30π and W₂₀=2.5λ. The left image is the FFT PSF, and the right image is the approximate PSF.

Download Full Size | PDF

Figure 1 and Fig. 2 also show that when moving away from the origin in the right side, there is a dramatic decay accompanied by the oscillations with the increase of frequency, no matter the aperture is bounded or unbounded. From Eq. (4), one can find that the PSF in the range of -(kW₂₀)²/(6πα)< x ≤ (3α - 2|kW₂₀|)/(2π) is a sinusoidal decaying function with the increase of the frequency. Hence the bandwidth can be determined by the largest sinusoidal modulation frequency. The sinusoidal modulation frequency within this range can be expressed as

f (x) = \frac{d [\frac{4 {({(k W_{20})}^{2} + 6 π α x)}^{\frac{3}{2}}}{27 α^{2}} \frac{1}{2 π}]}{d x} = \frac{2}{3 α} \sqrt{{(k W_{20})}^{2} + 6 π α x}, - \frac{{(k W_{20})}^{2}}{6 π α} < x \leq \frac{3 α - 2 ∣ k W_{20} ∣}{2 π} .

It’s a monotonously increasing function with the maximal modulation frequency equal to 2[1-|kW₂₀|/(3α)]. The ratio of this bandwidth to the traditional imaging system’s bandwidth, which is 2 when the normalized spatial coordinate for the pupil are used, is equal to the one derived in the subsection 3.1. An FFT PSF and an approximate PSF of a two-dimensional (2D) system with a cubic phase mask in the presence of defocus aberration are shown in Fig. 3 for α=30π and W₂₀=2.5λ. In these 2D images, regions of larger intensity are given by lighter shades. From the previous analysis, it is easy to understand why the grid points with decreasing intensities and spaces exist in the PSF. From Eq. (7) it can be concluded that a larger α leads to a less space between grid points in the PSF, so Eq. (7) can also be used to estimate the strength of the modulation of the cubic phase mask. Compared to the PSF area width approach described in the subsection 3.1, only the middle part of the PSF is needed instead of the whole PSF.

Besides, a tail with few oscillations occurs when the reduced spatial coordinate gets closer to the right boundary, and the tail is longer in the presence of larger defocus, which can be found from Eq. (4) and Fig. 1. This is because there is only one stationary point in this range. However the tail does not exist in the system with an unbounded aperture, as there are two stationary points instead of one.

The oscillations are due to the cubic phase in the pupil function, which leads to two stationary points. The oscillations with increasing frequency and decreasing intensities can also be found in systems with an anti-symmetric phase mask [3], but only the attenuations occur in the system with a symmetric phase mask because there is only one stationary point present there.

3.3 Sensitivities to aberrations

Fig. 4. PSFs calculated from Eq. (4) for α=30π and W₂₀=0, 1λ, 2λ, 3λ.

Download Full Size | PDF

From Eq. (5), one can conclude that for an unbounded aperture there is nothing occurring but a spatial shift of (kW₂₀)²/(6πα) when the wavefront coding system suffers from the defocus aberration, which only leads to a linear phase shift exp[-j(kW₂₀)² u/3α] in OTF. For a bounded aperture, an equation similar to Eq. (5) can be derived from Eq. (3) with a substitution u = t - kW₂₀/(3α),

h (x, W_{20}) = \frac{1}{2} {∣ \int_{- \infty}^{+ \infty} \exp (j α u^{3} + j k W_{20} u^{2} - j 2 π x u) du ∣}^{2}

In addition to the spatial shift, the integral range which corresponds to the aperture shifts by (kW₂₀)/3α from its original position. If the aperture shift is negligible compared to the spatial shift, the same phase shift will take place in the OTF, i.e., the MTF will not be significantly affected. Thus the wavefront coding system with a cubic phase mask is insensitive to the defocus aberration. PSFs calculated from Eq. (4) for α=30π and W₂₀=0, 1λ, 2λ, 3λ are shown in Fig. 4 respectively, where the PSFs with only spatial shifts can be found.

When the astigmatism occurs, an additional wavefront term of exp(jkW₂₂u²) is added into the original wavefront, where W₂₂ is the traditional astigmatism aberration constant in unit of wavelength. It can be noticed that a defocus of kW₂₂ has been introduced, so the wavefront coding system with a cubic phase mask is also insensitive to the astigmatism aberration.

Fig. 5. PSFs calculated from Eq. (4) for α=30π and W₃₁=0, 1λ, 2λ, 3λ.

Download Full Size | PDF

As to the coma aberration, an additional wavefront term of exp[jkW₃₁(u ² + v ²)u] is added into the original wavefront of the 2D system, where W₃₁ is the traditional coma aberration constant in unit of wavelength, u and v are reduced spatial coordinates. For simplicity, v is set to be zero for the case of 1D system, so the PSF in the presence of coma can be given by

h (x, W_{31}) = \frac{1}{2} {∣ \int_{- 1}^{+ 1} \exp (j α u^{3} + j k W_{31} u^{3} - j 2 π x u) d u ∣}^{2}

From Eq. (9), we can see that the normalized parameter α for the cubic phase mask is changed to (α+2πW₃₁/λ). Thus the sensitivity to coma can be derived from Eq. (4). PSFs calculated from Eq. (4) for α=30π and W₃₁=0, 1λ, 2λ, 3λ are shown in Fig. 5. The standard deviation (STD) of PSFs is also shown in Fig. 5 in black line, which indicates large errors between PSFs.

Though the wavefront coding system with a cubic phase mask is sensitive to the coma aberration, another type of phase mask can be employed to eliminate this aberration using the same analysis method applied for the defocus aberration. A spatial shift for a phase mask of order n leads to phase terms of lower order, which can be used not only to produce a phase mask of lower order [14], but also to eliminate the aberrations of lower order. With a small spatial shift of the phase mask, the aberrations of lower order occur without significantly affecting the MTF, thus the system with a phase mask has more or less insensitivities to aberrations of lower order. The PSF of the 1D wavefront coding system with a quartic phase mask in the presence of coma can be written as Eq. (10) with a substitution u = t -β /(4α),

h (x, W_{31}) = \frac{1}{2} {∣ \int_{- 1}^{+ 1} \exp (j α u^{4} + j β u^{3} - j 2 π x u) d u ∣}^{2}

where α is the normalized parameter for quartic phase mask in unit of radian, andβ is equal to kW₃₁. Besides the spatial shift and the aperture shift, coma can be treated as defocus effectively. The system with a quartic phase mask is sensitive to defocus [2], nevertheless the system is insensitive to coma. That is because even if W₃₁ is equal to 6.4λ, only a defocus of -2π is introduced when α is equal to 30π. PSFs of the system with a quartic phase mask in the presence of coma are shown in Fig. 6 for α=30π and W₃₁=0, 2λ, 4λ, 6λ, where small STD errors can be found.

Fig. 6. FFT PSFs of the system with a quartic phase mask for α=30π and W₃₁=0, 2λ, 4λ, 6λ.

Download Full Size | PDF

4. Conclusion

In summary, the approximate expression of PSF in the presence of defocus aberration is derived for the wavefront coding system with a cubic phase mask. Based on this, the width of the PSF area and the bandwidth in the presence of defocus aberration are calculated, and the oscillations that have not been described in the previous study are now analyzed. Both approximate expressions derived for the PSF area width and the oscillations can be applied to evaluate the strength of the cubic phase mask’s modulation. We also analyze the reason within the spatial domain why the system with a cubic phase mask is insensitive to defocus and astigmatism, but sensitive to coma. Besides, phase masks of higher order are proposed to eliminate the lower order aberrations. As an example, a quartic phase mask is used to eliminate the coma.

Appendix A:

Derivation of the PSF for cubic phase mask in the presence of defocus For the integral

{PSF}_{c} (x) = \int_{- 1}^{+ 1} \exp (j α u^{3} + j ψ u^{2} - j 2 πxu) du,

where μ(u) = αu ³ + ψu ² - 2πux, and ψ = kW ₂₀. The first and second derivatives of μ (u) are respectively,

μ' (u) = 3 α u^{2} + 2 ψ u - 2 π x,

μ ′′ (u) = 6 α u + 2 ψ .

The stationary point u₀ exists only when the following conditions are satisfied,

{\begin{matrix} μ' (u_{0}) = 3 α u_{0}^{2} + 2 ψ u_{0} - 2 π x = 0 \\ μ ′′ (u_{0}) = 6 α u_{0} + 2 ψ \neq 0 . \\ ∣ u_{0} ∣ \leq 1 \end{matrix}

So the stationary point u₀ only exists in the range of $- \frac{ψ^{2}}{6 πα} < x \leq \frac{3 α + 2 ∣ ψ ∣}{2 π},$ , the integral outside of which is assumed to be zero.

Two stationary points exists in the range of $- \frac{ψ^{2}}{6 πα} < x \leq \frac{3 α - 2 ∣ ψ ∣}{2 π},$ ,

{\begin{matrix} u_{01} = \frac{- 2 ψ + \sqrt{4 ψ^{2} + 24 π α x}}{6 α} = \frac{- ψ + \sqrt{ψ^{2} + 6 π α x}}{3 α} \\ u_{02} = \frac{- 2 ψ - \sqrt{4 ψ^{2} + 24 π α x}}{6 α} = \frac{- ψ - \sqrt{ψ^{2} + 6 π α x}}{3 α} . \end{matrix}

According to the stationary phase approximation, Eq. (A 1) can be approximated as the sum of the stationary phase approximations evaluated at these two stationary points, so

{PSF}_{c} (x) \approx \sqrt{\frac{2 π}{∣ μ ′′ (u_{01}) ∣}} \exp {jsign [μ ′′ (u_{01})] \frac{π}{4}} \exp [j μ (u_{01})] +

where sign(x) is the signum function defined as

sign (x) = {\begin{matrix} 1 x > 0 \\ 0 x = 0 \cdot \\ - 1 x < 0 \end{matrix}

Only one stationary point exists in the range of $\frac{3 α - 2 ∣ ψ ∣}{2 π} < x \leq \frac{3 α + 2 ∣ ψ ∣}{2 π},$ , so Eq. (A 1) can be approximated as the stationary phase approximation evaluated at the stationary point. When ψ>0, u₀₁ exists, then

{PSF}_{c} (x) \approx \sqrt{\frac{2 π}{∣ μ ′′ (u_{01}) ∣}} \exp {jsign [μ ′′ (u_{01})] \frac{π}{4}} \exp [j μ (u_{01})] .

Otherwise, u₀₂ exists, then

{PSF}_{c} (x) \approx \sqrt{\frac{2 π}{∣ μ ′′ (u_{02}) ∣}} \exp {jsign [μ ′′ (u_{02})] \frac{π}{4}} \exp [j μ (u_{02})] .

Inserting Eq. (A5) into μ(u) and μ″ (u), one can write

{\begin{matrix} μ ′′ (u_{01}) = - μ ′′ (u_{02}) = 2 \sqrt{ψ^{2} + 6 π α x} \\ μ (u_{01}) - μ (u_{02}) = \frac{- 4}{27 α^{2}} {(ψ^{2} + 6 π α x)}^{\frac{3}{2}} \end{matrix} .

With Eq. (A 6), Eq. (A 8), Eq. (A 9) and Eq. (A 10), the modulus of PSFc(x) can be approximated as

∣ {PSE}_{c} (x) ∣ \approx {\begin{matrix} 0 \\ \sqrt{\frac{2 π}{∣ μ^{″} (u_{01}) ∣} {2 + 2 \cos [\frac{π}{2} + μ (u_{01}) - μ (u_{02})]}} \\ \sqrt{\frac{2 π}{∣ μ^{″} (u_{01}) ∣}} \\ 0 \end{matrix} \begin{matrix} x \leq - \frac{ψ^{2}}{6 πα} \\ - \frac{ψ^{2}}{6 πα} < x \leq \frac{3 α - 2 ∣ ψ ∣}{2 π} \\ \frac{3 α - 2 ∣ ψ ∣}{2 π} < x \leq \frac{3 α + 2 ∣ ψ ∣}{2 π} \\ x > \frac{3 α + 2 ∣ ψ ∣}{2 π} \end{matrix}

As to Eq. (3), the PSF h(x, W₂₀) of the incoherent imaging system is equal to the half of the square of modulus of Eq. (A 1), so Eq. (4) can be obtained.

References and links

1. E. R. Dowski and W. T. Cathey, “Extended depth of field through wavefront coding,” Appl. Opt. 34,1859–1866 (1995). [CrossRef] [PubMed]

2. S. Mezouari, G. Muyo, and A. R. Harvey, “Amplitude and phase filters for mitigation of defocus and third-order aberrations,” Proc. SPIE 5249,238–248 (2004). [CrossRef]

3. A. Castro and J. O. Castañeda, “Increased depth of field with phase-only filters: ambiguity function,” Proc. SPIE 5827,1–11 (2005). [CrossRef]

4. S. S. Sherif, E. R. Dowski, and W. T. Cathey, “A logarithmic phase filter to extend the depth of field of incoherent hybrid imaging systems,” Proc. SPIE 4471,272–280 (2001). [CrossRef]

5. S. Mezouari and A. R. Harvey, “Primary aberrations alleviated with phase pupil filters,” Proc. SPIE 4768,21–31 (2002). [CrossRef]

6. N. George and W. Chi, “Extended depth of field using a logarithmic asphere,” J. Opt. A 5,157–163 (2003). [CrossRef]

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

8. G. E. Johnson, P. E. X. Silveira, and E. R. Dowski, “Analysis tools for computational imaging systems,” Proc. SPIE 5817,34–44 (2005). [CrossRef]

9. M. Born and E. Wolf, Principles of Optics (Pergamon, Oxford, 1985).

10. D. L. Marks, R. A. Stack, and D. J. Brady, “Three-dimensional tomography using a cubic-phase plate extended depth-of-field system,” Opt. Lett. 24,253–255 (1999). [CrossRef]

11. J. W. Goodman, Introduction to Fourier optics (McGraw-Hill, 1996), Chap. 6.

12. M. R. Arnison, “Phase control and measurement in digital microscopy” (Sydney Digital Theses, Physics, 2006), http://hdl.handle.net/2123/569.

13. M. Somayji and M. P. Christensen, “Enhancing form factor and light collection of multiplex imaging systems by using cubic phase mask,” Appl. Opt. 45,2911–2923 (2006). [CrossRef]

14. T. Hellmuth, A. Bich, R. Börret, and A. Kelm, “Variable phaseplates for focus invariant optical systems,” in Optical Design and Engineering II, L. Mazuray and R. Wartmann, eds., Proc. SPIE 5962,596215 (2005). [CrossRef]

Point spread function characteristics analysis of the wavefront coding system

Abstract

1. Introduction

2. Derivation of the PSF

3. PSF characteristics analysis

3.1 Boundaries

3.2 Oscillations

3.3 Sensitivities to aberrations

4. Conclusion

Appendix A:

Derivation of the PSF for cubic phase mask in the presence of defocus For the integral

References and links

Cited By

Figures (6)

Equations (32)

Optics Express