Understanding the effects of groove structures on the MTF

Kevin Liang; Miguel A. Alonso

doi:10.1364/OE.25.018827

1. Introduction

The ability to more finely control subaperture tools used for the creation and finishing/polishing of optical surfaces has been linked to the possibility of generating surface shapes of greater complexity. Such freeform surfaces provide the opportunity of improvements in system volume and performance. However, subaperture tools are inevitably accompanied by mid-spatial frequency (MSF) error structures. With characteristic spatial frequencies greater than those of the typical low-order Zernike aberrations and lower than those of surface roughness, MSF errors pose a range of problems that include complicating the ability to set effective tolerances for surface quality [1–4].

In this work, we provide a new approach for understanding the effects of MSF structures and other errors on the performance of an optical system, characterized by its optical transfer function (OTF) and its magnitude, the modulation transfer function (MTF). We start by following similar ideas to those presented previously for another performance metric, the Strehl ratio [5], under the standard perturbation model [6]. The statistics of MSF structures on a given optical surface are gathered in what we refer to as the pupil-difference probability density (PDPD) function, which is the probability density for the aberration difference between any two pupil points with prescribed separation. The link between the PDPD and the OTF, for all wavelengths, is then given by a Fourier relation. This framework is used for analyzing several simple one-dimensional periodic groove patterns. These ideas are then extended to two-dimensional structures, and approximate results are found for the cases of parallel grooves (like those resulting from diamond milling) and concentric circular grooves (such as those resulting from diamond turning). This framework is novel in that it provides mathematical intuition on how surface errors resulting from manufacturing parameters such as tool profile and path affect the degradation of an optical system’s performance. We also comment on the effect of concatenating several types of error.

2. Pupil Difference Probability Density

For a given wavenumber, k = 2π/λ, the optical transfer function (OTF) of a system is defined as the autocorrelation of the system’s field at the exit pupil [7] due to an object point source:

OTF (k, ρ) = \frac{\iint A (q - ρ / 2) e^{i k W (q - ρ / 2)} A (q + ρ / 2) e^{- i k W (q + ρ / 2)} d^{2} q}{\iint A^{2} (q) d^{2} q},

where W (q) is the aberration function at the pupil point q (in dimensionless pupil units), A(q) ∈ [0, 1] is the pupil amplitude function (which is normally binary, although the treatment we present remains valid for an apodized pupil), and ρ is a pupil coordinate displacement. We now rewrite Eq. (1) by adopting the shorthand notation:

OTF (k, ρ) = \frac{1}{a} \iint_{O (ρ)} e^{i k Δ W (q, ρ)} d^{2} q,

where ΔW (q, ρ) = W (q − ρ/2) − W (q + ρ/2 and a = ∫∫ A²(q) d²q. The region of integration, O(ρ), refers to the (potentially weighted) area of the overlap region of the autocorrelating pupils. Note that we write for convenience the OTF as a function of separation ρ in the plane of the pupil, although the results will also be expressed in terms of spatial frequency f, proportional to ρ.

We now define the PDPD, denoted by P(η, ρ), as a probability density for the error ΔW with given pupil displacement ρ. That is, for given ρ, P(η, ρ)dη is the fraction of values of q over the overlap region for which the aberration difference ΔW (q, ρ) falls within the infinitesimal range [η, η + dη]. This function can be expressed mathematically as

P (η, ρ) = \frac{\iint_{O (ρ)} δ [η - Δ W (q, ρ)] d^{2} q}{a {OTF}_{perf} (ρ)},

where δ is the Dirac delta distribution and OTF_perf(ρ) = ∫∫_O₍_ρ₎ d²q/a is the system’s optical transfer function in the absence of errors (W = 0). We can see from this definition that ∫ P(η, ρ) dη = 1, for any ρ, in accordance with its interpretations as a probability density function. Note that, by construction, P(−η, ρ) = P(η, ρ), and so P(η, ρ) has a zero mean in η.

By inserting unity in the form of ∫ δ[η − ΔW (q, ρ)] dη = 1 into the integrand of Eq. (2) and reversing the order of integration, the OTF can be rewritten as

OTF (k, ρ) = \int \frac{1}{a} \iint_{O (ρ)} δ [η - Δ W (q, ρ)] d^{2} q e^{i k η} d η = {OTF}_{perf} (ρ) \tilde{P} (k, ρ),

where

\tilde{P}

is the Fourier transform of the PDPD:

\tilde{P} (k, ρ) = \int P (η, ρ) e^{i k η} d η .

The corresponding result for the MTF is

MTF (k, ρ) = {MTF}_{perf} (ρ) | \tilde{P} (k, ρ) |, {MTF}_{perf} (ρ) = | {OTF}_{perf} (ρ) | .

Equations (4) and (6) state the first main result of this work: for systems where W can be regarded as independent of wavelength (e.g. reflective systems), the OTF as a function of wavenumber k is simply the Fourier transform of the PDPD times the OTF of the unaberrated system. Note that for systems in which W does depend on wavelength, this formulation is still valid, i.e., Eq. (5) still applies. However, because P(η, ρ) then also depends on k, the inverse process of recovering the PDPD from the OTF by inverse Fourier transformation is not possible.

3. One-dimensional periodic groove structures

To illustrate how the PDPD gives a simple picture of MTF performance, we now consider four different periodic error patterns in one dimension that allow a simple analytic treatment: a binary rectangular groove with fill factor 1/2, a symmetric triangular groove, a sinusoidal groove, and a piecewise parabolic groove. These groove profiles are chosen because their simplicity helps illustrate the concepts introduced here and allow analytic solutions. In particular, the piecewise parabolic groove structure mimics the profile left behind by a diamond tool. For all patterns, the spatial period is denoted by T and the peak-to-valley height by h. Throughout, the dependence on the one-dimensional separation ρ is expressed in terms of the shorthands ${\hat{ρ}}_{T} = (ρ \mod T) / T$ and $z = 1 - | 1 - 2 {\hat{ρ}}_{T} |$ . The groove patterns we examine, and the corresponding PDPD and its Fourier transform are presented in Table 1. Illustrations of these results are shown in Fig. 1 and Visualization 1. The plots for the PDPD within the second column are shown on their side so that they can be related to the error differences shown in purple in the first column.

Table 1. Four groove shapes with their PDPD (second column) and the PDPD Fourier transform (third column).

View Table | View all tables in this article

Fig. 1 The PDPD functions for one-dimensional (a) rectangular, (b) triangular, (c) sinusoidal, and (d) piecewise parabolic periodic groove patterns. The left column shows two copies of the groove pattern (blue and orange) separated by ρ, and their difference (purple). The center and right columns are the corresponding PDPD and the modulus its Fourier transform (solid green). The dashed and dotted green lines show the approximations 1 − k²D(ρ)²/2 and exp[−k²D(ρ)²/2], described in Section 6. Here, we used ρ = T/6. Visualization 1 shows the same information for ρ varing over a complete period.

Download Full Size | PDF

4. Parallel groove patterns

The one-dimensional cases in Section 3 provide a way of understanding the two-dimensional case corresponding to milled surfaces, namely of surfaces that exhibit periodic groove patterns in one direction and are constant in the other. If the surface in question is located near a plane conjugate to the aperture stop, the error W at the pupil is approximately proportional to the surface error. Consider, for example, a circular pupil of radius R across which the error has the shape of grooves aligned with the y-axis. Under the assumption that R is much larger than the groove width T we can use the approximation P(η, ρ) ≈ P₁(η, ρ_x), where ρ_x is the x-component of ρ. This means that the MTF is given by

MTF (k, ρ) \approx {MTF}_{perf} (ρ) | {\tilde{P}}_{1} (k, ρ_{x}) | .

Figure 2 illustrates this result for the case of 5 parabolic groove profiles across the pupil.

Fig. 2 (a) Circular pupil whose error presents five parabolic grooves of height h across its diameter, aligned in the vertical direction. (b-d) MTF as a function of both kh and spatial frequency for different radial slices, corresponding to the directions of the orange, green, and blue arrows shown in (a). The axis variable ρ/R equals the normalized spatial frequency f/f_c, where f is the spatial frequency in the direction of the radial slice and f_c is the coherent cutoff frequency for the circular pupil.

Download Full Size | PDF

5. Circular groove patterns

It is also insightful to see how this framework can be used for estimating the MTF of a system for which the errors are composed of concentric circular grooves with identical cross sections, such as those generated by diamond turning. That is, we assume that the error is only a function of the radial pupil coordinate q according to W (q), with radial period T. The rotational symmetry of the pupil error can be exploited to derive approximate analytic results that build upon the one-dimensional analysis presented in Section 3.

Consider again a circular pupil of radius R. Once again, ρ denotes the relative displacement between the two shifted copies of the pupil, but it is treated as a scalar given the problem’s rotational symmetry. Let a point q within the overlap region be specified by r₁ and r₂, which are the distances from the centers of the two shifted pupil copies, as shown in Fig. 3. Note that, due to the problem’s symmetry, it is sufficient to consider points above the line joining these centers. In fact, we can just consider points for which r₁ > r₂. That is, only a quarter of the overlap region is sufficient. Let us define the centroid and difference coordinates:

r = \frac{r_{1} + r_{2}}{2}, \bar{ρ} = r_{2} - r_{1} .

Fig. 3 Diagrams representing the various coordinates, for ρ < R (a,b) and ρ ≥ R (c,d). In (a) and (c), two copies of the pupil with a rotationally symmetric periodic groove pattern are shown, mutually displaced by ρ. The distances to a given point q from the centers of the pupils are r₁ and r₂. The overlap region is mapped onto the $(\bar{ρ}, r)$ plane in (b) and (d), taking the shape of a polygon whose boundaries are given in Eqs. (10) and (11). The elliptic and hyperbolic (red and blue) regions in the (q_x, q_y) space over which both sets of grooves are nearly parallel map onto straight bands (red and blue) in the $(\bar{ρ}, r)$ space. Visualization 2 shows the same information for varying ρ.

Download Full Size | PDF

The change of variables from q to r and $\bar{ρ}$ has a Jacobian given by

J (r, \bar{ρ}) = \sqrt{\frac{r^{2} - {(ρ / 2)}^{2}}{ρ^{2} - {\bar{ρ}}^{2}}} + \frac{1}{4} \sqrt{\frac{ρ^{2} - {\bar{ρ}}^{2}}{r^{2} - {(ρ / 2)}^{2}}} .

The region of integration O(ρ) is defined by the bounds of r and $\bar{ρ}$ , given either by

\frac{ρ}{2} \leq r \leq R - \frac{| \bar{ρ} |}{2} and 0 \leq \bar{ρ} \leq R - | R - ρ |,

or

0 \leq \bar{ρ} \leq Min [ρ, 2 (R - r)] and \frac{ρ}{2} \leq r \leq R,

depending on whether we want the range in r to depend on

\bar{ρ}

or vice-versa. Of course, integration of the Jacobian over these bounds gives back F(ρ) = OTF_perf (ρ)/4, which is the area of a quadrant of the overlap region. The regions defined by these bounds in the spaces of q and

(\bar{ρ}, r)

are shown in Fig. 3 for the two distinct regimes ρ < R and ρ > R. In terms of these new variables, Eq. (3) can be written as

P (η, ρ) = \frac{1}{F (ρ)} \iint_{O (ρ)} J (r, \bar{ρ}) \bar{δ} (η, r, \bar{ρ}) d r d \bar{ρ},

where

\bar{δ} (η, r, \bar{ρ}) = δ [η - W (r_{1}) + W (r_{2})] = δ [η - W (r - \bar{ρ} / 2) + W (r + \bar{ρ} / 2)]

, which has period 2T in

\bar{ρ}

and period T in r.

It is easy to see that there are two regions where the Jacobian in Eq. (9) becomes large, i.e. where the mapping from the area in the pupil overlap region bunches up in the space of $\bar{ρ}$ and r. The first of these regions is near the boundary r = ρ/2 and the second is near the boundary $\bar{ρ} = ρ$ . The first corresponds to the line segment between the two centers, and the second to the line segments between the centers and the edges of the overlap region (only when ρ < R). Note from Fig. 3 that these are regions where there is significant alignment of the grooves from both pupil copies (the relative displacement of the grooves of the two copies is at most one period more than at the axis). Over the rest of the area, the grooves of both copies intersect at a considerable angle.

The approximate results that we set out to obtain are expressed in terms of three quantities. The first is precisely the PDPD for the one-dimensional groove structure in question, like those summarized in Table 1. This one-dimensional PDPD can be written as an integral over a single period of $\bar{δ} (η, r, \bar{ρ})$ in r:

P_{1} (η, \bar{ρ}) = \frac{1}{T} \int_{0}^{T} \bar{δ} (η, r, \bar{ρ}) d r .

Although its integrand has period 2T in $\bar{ρ}$ , it can be shown that $P_{1} (η, \bar{ρ})$ has period T (as do the PDPDs in Table 1). The second is defined similarly, but the integration is over the difference:

{\bar{P}}_{1} (η, r) = \frac{1}{2 T} \int_{0}^{2 T} \bar{δ} (η, r, \bar{ρ}) d \bar{ρ} .

It can be shown that ${\bar{P}}_{1} (η, r)$ has period T/2 even though the integrand has period T in r. Notice that if the error is symmetric around some pivot point r_c, i.e. W (r) = W (r_c − r) (as are all groove structures in Table 1), then ${\bar{P}}_{1} (η, r) = P_{1} (η, 2 r - r_{c})$ . The third quantity is defined as

P_{0} (η) = \frac{1}{T} \int_{0}^{T} P_{1} (η, \bar{ρ}) d \bar{ρ} = \frac{2}{T} \int_{0}^{T / 2} \bar{P} (η, r) d r,

which is the average value of

P_{1} (η, \bar{ρ})

or

{\bar{P}}_{1} (η, r)

over a period. Note that this is the probability density that was used in [5] to describe the Strehl ratio.

We show in Appendix A that an approximate expression for the PDPD is given by

P (η, ρ) \approx \frac{1}{F (ρ)} [\frac{A (ρ)}{2} {\hat{Q}}_{t} P_{1} (η, ρ) + \frac{B (ρ)}{\sqrt{2}} {\hat{Q}}_{τ} {\bar{P}}_{1} (η, ρ / 2) + E (ρ) P_{0} (η)],

and the corresponding OTF estimate is found as 4F times the Fourier transform in η of this result:

OTF (k, ρ) \approx 4 [\frac{A (ρ)}{2} {\hat{Q}}_{t} {\tilde{P}}_{1} (k, ρ) + \frac{B (ρ)}{\sqrt{2}} {\hat{Q}}_{τ} {\tilde{\bar{P}}}_{1} (k, ρ / 2) + E (ρ) {\tilde{P}}_{0} (k)] .

In these results, A(ρ), B(ρ), and E(ρ) are the areas shown in blue, red, and green in Fig. 3, respectively, such that A(ρ) + B(ρ) + E(ρ) = F(ρ). Their explicit functional forms are given in Appendix A for the case considered here of a circular aperture. Note, however, that were the pupil to have a different shape or not be centered at the center of the grooves (due, for example, to an off-axis beam footprint), the approximate result above would still be valid, as long as the areas A, B, E, and F are calculated with the appropriate boundaries. Also, t = Min(ρ, T) and τ = −Min(2R − ρ, T), and the integral operator ${\hat{Q}}_{x}$ applied to a function f(ρ) is defined as

{\hat{Q}}_{x} f (ρ) = \int_{0}^{1} \frac{f (ρ - x v)}{\sqrt{v}} d v .

Since this operator has the form of a convolution and is applied to periodic functions with period T, the easiest way to calculate its effect is through a Fourier series expansion of the function it acts on. That is, if the periodic function (assumed here to be even) is written as

f (ρ) = \frac{1}{2} a_{0} + \sum_{m = 1}^{\infty} a_{m} \cos (\frac{2 π m ρ}{T}),

then it can be shown that this operator simply modifies the Fourier series as

{\hat{Q}}_{x} f (ρ) = a_{0} + \sqrt{\frac{T}{| x |}} \sum_{m = 1}^{\infty} \frac{a_{m}}{\sqrt{m}} [C (2 \sqrt{\frac{m | x |}{T}}) \cos (\frac{2 π m ρ}{T}) + sgn (x) S (2 \sqrt{\frac{m | x |}{T}}) \sin (\frac{2 π m ρ}{T})],

where C and S are the Fresnel cosine and sine integrals, respectively. Hence, in order to solve the integrals in Eq. (17), we have to find the Fourier coefficients a_m for

{\tilde{P}}_{1}

. These coefficients are summarized in Table 2 for the 1D PDPD Fourier transforms in Table 1. Details of their derivation are given in Appendix B. (It turns out that the coefficients in the third row were found previously in [3].) Note that, for computational purposes, only a finite number of terms in the Fourier series are used, and the number required to achieve a given level of error depends on kh. The last column on Table 2, labelled as m_10%, gives a simple estimate of the terms m ≤ Max(0, m_10%) to be used to achieve an error of under 10%. In all cases, a single Fourier coefficient (m = 1) is sufficient to achieve this level of error for peak-to-valley aberration heights up to about half a wavelength. Finally, notice that

{\tilde{P}}_{0} (k)

in Eq. (17) corresponds simply to a₀/2, given explicitly in Table 2.

Table 2. Fourier expansion coefficients (second and third columns) for the functions $\tilde{P} (k, ρ)$ in Table 1, corresponding to the groove structures shown in the first column. The fourth column shows an estimate of how many Fourier terms are needed (m ≤ Max(0, m_10%)) for the truncated Fourier series to achieve an error below 10%.

View Table | View all tables in this article

Figure 4 shows the resulting MTF plots of the four groove shapes as calculated by Eq. (17), in comparison to those generated numerically. One can see that, even with the approximations used in the analytic result and despite the fact that there are only ten grooves across the pupil radius, the analytic MTF plots match the numerical ones very well, giving results that are consistent to within about 5%. Figure 5 shows, for the same structures, plots of these MTF estimates, as well as of each of the three individual contributions to the estimate in Eq. (17), for kh = 1.6, and Visualization 3 shows the same data for varying kh. Notice that for kh ± 0 the MTF presents an initial peak for ρ < T, followed by oscillations around a baseline, shown in gray in Fig. 5; this general behavior for the MTF has been described before by employing a different approach [1]. It turns out that the equation for this baseline is very simple, and gives a basic estimate of the MTF for ρ > T. Note from Eq. (20) and the fact that $a_{0} = 2 {\tilde{P}}_{0} (k)$ that the first and second terms in Eq. (17) include parts proportional to ${\tilde{P}}_{0} (k)$ , and that these, added to the third term of Eq. (17), give simply ${MTF}_{perf} (ρ) {\tilde{P}}_{0} (k)$ . Therefore, this baseline is simply the MTF of the unaberrated system scaled by ${\tilde{P}}_{0} (k)$ . The number of oscillations around the baseline equals the number of grooves across the diameter, and their amplitude scales not only with kh but also with $\sqrt{T / R}$ , since the areas A(ρ) and B(ρ) are roughly proportional to the square root of the period.

Fig. 4 MTF for four concentric circular groove shapes with R = 5T: (a) rectangular with a fill factor of 1/2, (b) triangular, (c) sinusoidal, and (d) piecewise parabolic. For each, the upper plot was generated from Eq. (17), and the lower plot resulted from a numerical calculation.

Download Full Size | PDF

Fig. 5 MTF for the same circular groove shapes used in Fig. 4, for kh = 1.6. In each figure, the orange dots correspond to a numerical computation, of the MTF, the black curve corresponds to the estimate in Eq. (17), and the blue, red and green curves correspond to the contributions from each of the terms in Eq. (17), in that order. The gray curves show the baseline ${MTF}_{perf} (ρ) {\tilde{P}}_{0} (k)$ . Visualization 3 shows the same information for varying kh.

Download Full Size | PDF

6. Concatenation of errors, randomness, and relation to other measures

The use of the PDPD as a probability density provides a simple way of understanding the concatenation of independent errors, and allows the study of not only deterministic error shapes like those described so far, but also of stochastic errors that follow certain statistics. There are often multiple error contributions within an optical system, coming from separate surfaces, from different manufacturing processes for a single surface, or even from fluctuating media such as the atmosphere. For the cases where these contributions are approximately uncorrelated, the total PDPD is constructed following the standard concatenation rule for independent probability densities, as a convolution of the PDPDs P_i from each contribution:

P (η, ρ) = (P_{1} * P_{2} * \dots * P_{n}) (η, ρ)

Then, from the Fourier convolution theorem, the MTF of the total system is simply the MTF of the unaberrated system times the product of the Fourier transforms of the individual PDPDs [8]:

MTF (k, ρ) = {MTF}_{perf} (ρ) \prod_{i = 1}^{n} | {\tilde{P}}_{i} (k, ρ) | .

Note that, for this relation to hold, the two different error contributions do not need to be uncorrelated due to independent random fluctuations, but simply from the point of view of the PDPD; for example, two deterministic groove structures in different directions or with different periods have PDPDs that are effectively uncorrelated. This is because, for any given ρ, if we were to pick a large series of points q spanning the overlap region, the values of ΔW for the two contributions would be largely uncorrelated.

Let us now discuss the connection between the PDPD and another error measure that has also been used to study effects on optical performance, the structure function [9–11]. The structure function D(ρ) is defined as the integral of ΔW²(q, ρ) over the overlap region [12]. It is easy to see that the structure function is simply the variance of the PDPD, that is,

D (ρ) = \iint_{O (ρ)} Δ W^{2} (q, ρ) d^{2} q = \int η^{2} P (η, ρ) d η .

While the PDPD contains more information than the structure function, it is worth mentioning two cases where D(ρ) gives an approximate full description of both the PDPD and MTF. The first is that where P(η, ρ) is localized within an interval in η smaller than the wavelength. An approximate expression for the MTF can then be found by expanding the exponential in Eq. (5) in a Taylor series up to second order, so that Eq. (6) becomes

MTF (k, ρ) \approx {MTF}_{perf} (ρ) | \int P (η, ρ) (1 + i k η - k^{2} \frac{η^{2}}{2}) d η | = {MTF}_{perf} (ρ) [1 - \frac{k^{2} D (ρ)}{2}],

where we used the fact that the first moment of P(η, ρ) in η vanishes. For the second case, suppose that the PDPD is approximately Gaussian, e.g., due to the concatenation of several statistically independent errors discussed earlier, following the Central Limit Theorem. Then, since the Fourier transform of a Gaussian is simply another Gaussian with the reciprocal standard deviation (the structure function), we arrive at the standard approximation for the MTF used, for example, in atmospheric studies [13]:

MTF (k, ρ) \approx {MTF}_{perf} (ρ) \exp [- \frac{k^{2} D (ρ)}{2}] .

Note that both limiting results above have counterparts in the case of the Strehl ratio, as described in [5], corresponding to the expressions by Maréchal [14] and Mahajan [15]. The factors multiplying MTF_perf in these two simple estimates are shown for the one-dimensional groove structures in the fourth column of Fig. 1 as dashed and dotted green lines, where we see that they only resemble the true results for small kh.

7. Conclusion

We showed how the MTF of an optical system for all wavelengths can be understood in terms of the PDPD, a function that characterizes the MSF phase errors (whether these are periodic or not). Once the PDPD is known for a given wavefront error in an optical system, it can be Fourier transformed to obtain the MTF. The implementation of this framework is illustrated for several simple 1D periodic groove structures. Further, it is shown that for 2D errors consisting of parallel or concentric periodic grooves over the pupil, the MTF can be easily calculated from the PDPD of the corresponding one-dimensional pattern. These surfaces are important because they are characteristic of the MSF errors left behind on optical surfaces under the diamond milling and turning processes of freeform manufacturing, and when these surfaces are nearly conjugate to the aperture stop, the wavefront error they produce on the pupil mimics their surface error [6].

Several simple generalizations can be made to this framework. First, the groove patterns could differ from the linear and concentric ones examined explicitly in this paper. Changing the groove pattern would change the effective areas A(ρ) and B(ρ) over which the grooves from the pupil replicas align to within a period. A similar modification can be made to account for fluctuations in the periodicity, which blur the oscillations in the MTF. In fact, the case of non-periodic surface errors can be thought of as the limit when A(ρ) and B(ρ) tend to zero. Therefore, aside from the initial peak of the MTF for |ρ| smaller than the typical feature size, we can see from Eq. (17) that the MTF is in this case dominated by ${\tilde{P}}_{0} (k)$ . That is, except for the central peak, the MTF for non-periodic errors is well approximated by the baseline ${MTF}_{perf} (ρ) {\tilde{P}}_{0} (k)$ .

Finally, while we focused on the MTF, these results can be extended to another well-known metric of optical performance, the encircled energy, which can be relatede to the OTF as [16]:

EE (k, ρ) = \frac{1}{OTF (k, 0) S_{O} (0)} \int \int_{- \infty}^{\infty} OTF (k, ρ) S_{O} (ρ) S_{D} (ρ) d^{2} ρ,

where S_O and S_D are the Fourier transforms of the object and detector functions, respectively. By using Eq. (17), one can find an expression for the encircled energy explicitly in terms of the PDPD, as we have done with the MTF. Note that for a point source and a point detector, Eq. (26) reduces to the expression for the Strehl ratio [5].

A. Derivation of PDPD for circular grooves

We start by inserting Eq. (9) into Eq. (12) and separating the result into two parts as P(η, ρ) ≈ (I₁ + I₂)/F, where

I_{1} = \int_{0}^{R - | R - ρ |} \frac{1}{\sqrt{ρ^{2} - {\bar{ρ}}^{2}}} \int_{ρ / 2}^{R - \bar{ρ} / 2} \bar{δ} (η, r, \bar{ρ}) \sqrt{r^{2} - {(ρ / 2)}^{2}} d r d \bar{ρ},

I_{2} = \frac{1}{4} \int_{ρ / 2}^{R} \frac{1}{\sqrt{r^{2} - {(ρ / 2)}^{2}}} \int_{0}^{Min [ρ, 2 (R - r)]} \bar{δ} (η, r, \bar{ρ}) \sqrt{ρ^{2} - {\bar{ρ}}^{2}} d \bar{ρ} d r .

Note that for I₁ and I₂ the limits of integration were specified according to Eqs. (10) and (11), respectively. In particular, let us look at the integral in r within I₁, and break it as a sum of many integrals over a period each. Due to the periodicity of $\bar{δ}$ we can write

\int_{ρ / 2}^{R - \bar{ρ} / 2} \bar{δ} (η, r, \bar{ρ}) \sqrt{r^{2} - {(ρ / 2)}^{2}} d r \approx \int_{0}^{T} \bar{δ} (η, ρ / 2 + v, \bar{ρ}) \sum_{n = 0}^{N - 1} \sqrt{r_{n}^{2} - {(ρ / 2)}^{2}} d v,

where

N = ⌊ (R - | \bar{ρ} / 2 | - ρ / 2) / T ⌋

, r_n = ρ/2+nT+v, and we used

\bar{δ} (η, r_{n}, \bar{ρ}) = \bar{δ} (η, ρ / 2 + v, \bar{ρ})

. Equation (29) is an approximation because the integration range is generally not an integer multiple of T. If T is small compared to ρ, we can approximate the summation as another integral:

\sum_{n = 0}^{N - 1} \sqrt{r_{n}^{2} - {(ρ / 2)}^{2}} \approx \frac{1}{T} \int_{ρ / 2}^{R - \bar{ρ} / 2} \sqrt{{r^{'}}^{2} - {(ρ / 2)}^{2}} d r^{'} .

Substituting Eqs. (29) and (30) into Eq. (27) and then using Eq. (13) gives

I_{1} \approx \int_{0}^{R - | R - ρ |} P_{1} (η, \bar{ρ}) \frac{α (ρ, \bar{ρ})}{\sqrt{ρ - \bar{ρ}}} d \bar{ρ}, α (ρ, \bar{ρ}) = \int_{ρ / 2}^{R - \bar{ρ} / 2} \sqrt{\frac{{r^{'}}^{2} - {(ρ / 2)}^{2}}{ρ + \bar{ρ}}} d r^{'} .

We can use an analogous procedure for I₂ in Eq. (28) to find

I_{2} \approx \int_{ρ / 2}^{R} {\bar{P}}_{1} (η, r) \frac{β (ρ, r)}{\sqrt{r - ρ / 2}} d r, β (ρ, r) = \frac{1}{4} \int_{0}^{Min [ρ, 2 (R - r)]} \sqrt{\frac{ρ^{2} - {\bar{ρ}}^{'}^{2}}{r + ρ / 2}} d {\bar{ρ}}^{'} .

To simplify these expressions, we discretize again the outermost integrals. Let us this time start with I₂ (and then apply a similar procedure to I₁). The integral in r in Eq. (32) is broken into a sum of integrals, each over an interval of T/2 (the period of ${\bar{P}}_{1}$ ). By writing r_m = ρ/2 + mT/2 + u with M = ⌊R − ρ/2⌋, and using ${\bar{P}}_{1} (η, r_{m}) = {\bar{P}}_{1} (η, ρ / 2 + u)$ , we find

I_{2} \approx \int_{0}^{T / 2} {\bar{P}}_{1} (η, ρ / 2 + u) \sum_{m = 0}^{M - 1} \frac{β (ρ, r_{m})}{\sqrt{r_{m} - ρ / 2}} d u .

Notice that for ρ > 0, β is a smooth, nonsingular function of its second argument. However, we cannot approximate the sum in Eq. (33) as an integral, given the strong dependence of its first term (m = 0) on u. It can be shown, however, that if this term is separated the remaining sum can be approximated as an integral:

\begin{array}{l} I_{2} \approx \int_{0}^{T / 2} \frac{{\bar{P}}_{1} (η, ρ / 2 + u)}{\sqrt{u}} β (ρ, ρ / 2 + u) d u + \int_{0}^{T / 2} {\bar{P}}_{1} (η, ρ / 2 + u) \sum_{m = 1}^{M - 1} \frac{β (ρ, r_{m})}{\sqrt{r_{m} - ρ / 2}} d u \\ \approx β (ρ, ρ / 2) \int_{0}^{T / 2} \frac{{\bar{P}}_{1} (η, ρ / 2 + u)}{\sqrt{u}} d u + \int_{0}^{T / 2} {\bar{P}}_{1} (η, ρ / 2 + u) \frac{2}{T} \int_{ρ / 2 + T / 2}^{R} \frac{β (ρ, r^{'})}{\sqrt{r^{'} - ρ / 2}} d r^{'} d u \\ = β (ρ, ρ / 2) \int_{0}^{T / 2} \frac{{\bar{P}}_{1} (η, ρ / 2 + u)}{\sqrt{u}} d u + P_{0} (η) \int_{ρ / 2 + T / 2}^{R} \frac{β (ρ, r^{'})}{\sqrt{r^{'} - ρ / 2}} d r^{'}, \end{array}

where we used β(ρ, ρ/2+u) ≈ β(ρ, ρ/2) in the second step. The error of the approximation in Eq. (34) is surprisingly small, even for small M.

An analogous result can be found for I₁ by breaking the integral in $\bar{ρ}$ into a sum of integrals over a period each. However, it is important to note that the singularity of integrand is contained in the range of integration only when ρ ≤ R. Therefore, only in this case must we separate one term of the sum (the last one) before approximating the remaining sum as an integral. The result is then

I_{1} \approx {\begin{array}{l} α (ρ, ρ) \int_{0}^{T} \frac{P_{1} (η, ρ - u)}{\sqrt{u}} d u + P_{0} (η) \int_{0}^{ρ - T} \frac{α (ρ, \bar{ρ})}{\sqrt{ρ, \bar{ρ}}} d \bar{ρ}, & ρ \leq R, \\ P_{0} (η) \int_{0}^{2 R - ρ} \frac{α (ρ, \bar{ρ})}{\sqrt{ρ, \bar{ρ}}} d \bar{ρ}, & ρ > R, \end{array}

To finish the derivation, note that

\begin{array}{l} β (ρ, ρ / 2) = \frac{1}{4} \int_{0}^{Min [ρ, 2 (R - r)]} \sqrt{\frac{ρ^{2} - {\bar{ρ}}^{2}}{ρ}} d \bar{ρ} \\ \approx \frac{1}{\sqrt{2 T}} \int_{ρ / 2}^{ρ / 2 + T / 2} \int_{0}^{Min [ρ, 2 (R - r)]} J (r, \bar{ρ}) d \bar{ρ} d r = \frac{B (ρ)}{\sqrt{2 T}}, \end{array}

where B(ρ) is the area shown in red in Fig. 3. Similarly, for ρ ≤ R we get

α (ρ, ρ) \approx \frac{A (ρ)}{2 \sqrt{T}},

where A(ρ) is the area shown in blue in Fig. 3. Finally, the sum of the coefficients of P₀ in Eqs. (34) and (35) is approximately equal to E(ρ) = F(ρ) − A(ρ) − B(ρ), the remaining area shown in green in Fig. 3.

The approximations made in this appendix depend on the fact that T is small compared to ρ. This allowed the approximations in Eqs. (30) and (34). However, when considering ρ < T, adjustments are necessary in the expression for α(ρ, ρ). In particular, the “average value” over the integrals would not be over an entire period. That is, the first integral in Eq. (35) would go from 0 to ρ, rather than T. As a result,

α (ρ, ρ) \approx \frac{A (ρ)}{2 \sqrt{ρ}} .

Similar adjustments are necessary for β(ρ, ρ/2) when ρ > 2R − T. That is, in the region where B(ρ) is the dominant share of the overlap area F(ρ). In this case, the first integral in the last line of Eq. (34) would go from 0 to R − ρ/2, rather than T/2. As a result,

β (ρ, ρ / 2) \approx \frac{B (ρ)}{\sqrt{2 (R - ρ / 2)}} .

These changes are taken into account in Eq. (16).

To finish this appendix, we give the explicit expressions for the areas A, B, E and F for a circular pupil of radius R. For ρ < R we have

A (ρ) = \frac{1}{2} [\frac{V_{-}}{U_{-}} X_{-} \sqrt{X_{-}^{2} - U_{-}^{2}} - U_{-} V_{-} arccosh (\frac{X_{-}}{U_{-}}) + R^{2} arccos (\frac{X_{-}^{'}}{R}) - X_{-}^{'} \sqrt{R^{2} - X_{-}^{' 2}}],

B (ρ) = \frac{π}{4} U_{+} V_{+},

while for ρ ≥ R we have

A (ρ) = 0,

B (ρ) = \frac{1}{2} [\frac{V_{+}}{U_{+}} X_{+} \sqrt{U_{+}^{2} - X_{+}^{2}} + U_{+} V_{+} arcsin (\frac{X_{+}}{U_{+}}) + R^{2} arccos (\frac{X_{+}^{'}}{R}) - X_{-}^{'} \sqrt{R^{2} - X_{+}^{' 2}}],

where

U_{\pm} = \frac{ρ \pm T}{2}, V_{\pm} = \frac{\sqrt{T (2 ρ \pm T)}}{2}, X_{\pm} = \frac{[2 R - (ρ \pm T)] (ρ \pm T)}{2 ρ}, X_{\pm}^{'} = X_{\pm} + \frac{ρ}{2} .

Finally, for all ρ ∈ [0, 2R], we have

F (ρ) = \frac{1}{2} [R^{2} arccos (\frac{ρ}{2 R}) - \frac{ρ}{2} \sqrt{R^{2} - {(\frac{ρ}{2})}^{2}}],

E (ρ) = F (ρ) - A (ρ) - B (ρ) .

B. Fourier coefficients for ${\tilde{P}}_{1}$

Here we calculate the Fourier series coefficients of the functions ${\tilde{P}}_{1}$ given in the third column of Table 1, according to

a_{m} = 2 \int_{0}^{1} {\tilde{P}}_{1} (k, u T) \cos (2 π m u) d u,

for integer m ≥ 0. The cases corresponding to the first two rows are trivial to calculate after using the substitution z = 1 − |1 − 2u|. The result for the third row also follows straightforwardly from using a standard integral property of Bessel functions [17]. We therefore focus on the case within the fourth row, corresponding to piecewise parabolic groove structures. The substitution of

\tilde{P}

in Eq. (47) gives

a_{m} = \int_{0}^{1} \frac{\sin [4 k h u (1 - u)]}{2 k h u (1 - u)} \cos (2 π m u) d u,

where

u = {\hat{\bar{ρ}}}_{T}

. To begin, we separate the integrand using partial fraction decomposition, noting that

\frac{1}{u (1 - u)} = \frac{1}{u} + \frac{1}{1 - u} .

Then we have

a_{m} = \int_{0}^{1} \frac{\sin [4 k h u (1 - u)]}{2 k h u} \cos (2 π m u) d u + \int_{0}^{1} \frac{\sin [4 k h u (1 - u)]}{2 k h (1 - u)} \cos (2 π m u) d u .

Notice that we can make the change of variables 1 − u → u in the second integral without changing the limits of integration, which gives

a_{m} = \int_{0}^{1} \frac{\sin [4 k h u (1 - u)]}{2 k h u} \cos (2 π m u) d u + \int_{0}^{1} \frac{\sin [4 k h u (1 - u)]}{2 k h u} \cos [2 π m (1 - u) d u .

However, cos(2πmu) = cos[2πm(1 − u)], so the two integrals are identical. We can then write

a_{m} = \int_{0}^{1} \frac{\sin [4 k h u (1 - u)]}{k h u} \cos (2 π m u) d u = \int_{- 1}^{1} \frac{\sin [4 k h u (1 - | u |)]}{2 k h u} \cos (2 π m u) d u .

Note that the ratio within the integrand can now be written as another integral in the form

a_{m} = \int_{- 1}^{1} \int_{- 1 + | u |}^{1 - | u |} \exp (i 4 k h u u^{'}) d u^{'} \cos (2 π m u) d u .

We can now change variables of integration according to u = (u₁ + u₂)/2, u′ = (u₁ − u₂)/2, so that du du′ = du₁ du₂/2, giving

a_{m} = \frac{1}{2} \int_{- 1}^{1} \int_{- 1}^{1} \exp [i k h (u_{1}^{2} - u_{2}^{2})] \cos [π m (u_{1} + u_{2})] d u_{1} d u_{2} .

Using the cosine angle addition formula, we have

\begin{array}{l} a_{m} = \frac{1}{2} {\int \int_{- 1}^{1} \exp [i k h (u_{1}^{2} - u_{2}^{2})] \cos (π m u_{1}) \cos (π m u_{2}) d u_{1} d u_{2} \\ - \int \int_{- 1}^{1} \exp [i k h (u_{1}^{2} - u_{2}^{2})] \sin (π m u_{1}) \sin (π m u_{2}) d u_{1} d u_{2}} \\ = \frac{1}{2} [{| \int_{- 1}^{1} \exp (i k h u^{2}) \cos (π m u) d u |}^{2} - {| \int_{- 1}^{1} \exp (i k h u^{2}) \sin (π m u) d u |}^{2}] \\ = \frac{1}{2} {| \int_{- 1}^{1} \exp (i k h u^{2}) \cos (π m u) d u |}^{2}, \end{array}

where we used the fact that the second integral in the third line vanishes due to parity. The resulting integral can be solved to yield

a_{m} = \frac{π}{8 | k h |} {| Erfi [\frac{(1 + i) (2 k h - m π)}{4 \sqrt{k h / 2}}] + Erfi [\frac{(1 + i) (2 k h + m π)}{4 \sqrt{k h / 2}}] |}^{2},

where Erfi is the imaginary error function. Note that a₀ can be written in the more compact form

a_{0} = \frac{π}{2 k h} [C^{2} (\sqrt{\frac{2 k h}{π}}) + S^{2} (\sqrt{\frac{2 k h}{π}})] .

Funding

National Science Foundation (NSF) I/UCRC Center for Freeform Optics (IIP-1338877)

Acknowledgments

We thank the members of CeFO for their support and insights, and acknowledge in particular useful discussions with Thomas Suleski, Reza Aryan, Christoph Menke, Thomas Köhler, Flemming Tinker, James Fienup, Angela Davies, Chris Evans, Jannick Rolland, and Greg Forbes.

References and links

1. R. J. Noll, “Effect of Mid- and High-Spatial Frequencies on Optical Performance,” Opt. Eng. 18(2), 182137 (1979). [CrossRef]

2. D. Aikens, J. E. DeGroote, and R. N. Youngworth, “Specification and Control of Mid-Spatial Frequency Wavefront Errors in Optical Systems,” in Frontiers in Optics 2008/Laser Science XXIV/Plasmonics and Metamaterials/Optical Fabrication and Testing, OSA Technical Digest (CD) (Optical Society of America, 2008), paper OTuA1. [CrossRef]

3. J. M. Tamkin, T. D. Milster, and W. Dallas, “Theory of modulation transfer function artifacts due to mid-spatial-frequency errors and its application to optical tolerancing,” Appl. Opt. 49, 4895 (2010).

4. G. W. Forbes, “Never-ending struggles with mid-spatial frequencies,” Optical Measurement Systems for Industrial Inspection IX, Proc. of SPIE9525, 95251B (2015).

5. M. A. Alonso and G. W. Forbes, “The Strehl ratio as the Fourier transform of a probability density,” Opt. Lett. 41, 3735–3738 (2016). [CrossRef] [PubMed]

6. R. N. Youngsworth and B. D. Stone, “Simple estimates for the effects of mid-spatial-frequency surface errors on image quality,” Appl. Opt. 39, 2198–2209 (2000). [CrossRef]

7. J. W. Goodman, Introduction to Fourier Optics (Roberts & Company Publishers, 2005), Chap. 6.

8. G. D. Boremann, Modulation Transfer Function in Optical and Electro-Optical Systems (SPIE, 2001), pp. 85–88. [CrossRef]

9. L. He, C. J. Evans, and A. Davies, “Two-quadrant area structure function analysis for optical surface characterization,” Opt. Express 20, 23275–23280 (2012). [CrossRef] [PubMed]

10. L. He, C. J. Evans, and A. Davies, “Optical surface characterization with the area structure function,” CIRP Annals -Manufacturing Technology , 62(1), 539–542 (2013). [CrossRef]

11. R. E. Parks and M. T. Tuell, “The structure function as a metric for roughness and figure,” Proc. of SPIE 995199510 (2016).

12. J. W. Goodman, Statistical Optics (John Wiley & Sons, 2015), Chap. 8.

13. R. E. Hufnagel and N. R. Stanley, “Modulation Transfer Function associated with Image Transmission through Turbulent Media,” J. Opt. Soc. Am. 54, 52–61 (1964). [CrossRef]

14. A. Maréchal, “Etude des effets combinés de la diffraction et des aberrations géométriques sur l’image d’un point lumineux,” Rev. Opt. 2, 257–277 (1947).

15. V. N. Mahajan, “Strehl ratio for primary aberrations in terms of their aberration variance,” J. Opt. Soc. Am. 73, 860–861 (1983). [CrossRef]

16. J. V. Baliga and B. D. Cohn, “Simplified method for calculating encircled energy,” Proc. SPIE 892, 152 (1988). [CrossRef]

17. I. S. Gradshteyn and I. M. Ryxhik, Table of Integrals, Series, and Products (Academic, 1980), p. 739, ET II 360(14).

Name	Description
Visualization 1	This shows the autocorrelating groove structures (and their difference), the PDPD, and the magnitude of the Fourier transform of the PDPD in the right, center, and left columns, respectively.
Visualization 2	This shows the three different regions in the overlap region of two autocorrelating pupils, and the corresponding areas in the space with the change of coordinates.
Visualization 3	This shows the comparison between the MTFs generated analytically and numerically. Also shown are the components to the analytic MTF, which relate to the three regions on the overlap region.

Groove shape	P(η, ρ)	$\tilde{P} (k, ρ)$
	$\frac{z}{2} [δ (η - h) + δ (η + h)] + (1 - z) δ (η)$	$z \cos (k h) + (1 - z)$
	$\frac{1 - z}{2} [δ (η - z h) + δ (η + z h)] + \frac{1}{2 h} rect (\frac{η}{2 z h})$	$(1 - z) \cos (z k h) + \frac{\sin (z k h)}{k h}$
	$Re (\frac{1}{π \sqrt{{[h \sin (π {\hat{ρ}}_{T})]}^{2} - η^{2}}})$	$J_{0} [k h \sin (π {\hat{ρ}}_{T})]$
	$\frac{1}{8 h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})} rect [\frac{η}{8 h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})}]$	$\frac{\sin [4 k h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})]}{4 k h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})}$

${\tilde{P}}_{0} (k) = a_{0} / 2$	a_m	m_10%
$\frac{1 + \cos (k h)}{2}$	$[1 - \cos (k h)] \frac{2 [1 - {(- 1)}^{m}]}{m^{2} π^{2}}$	1
$\frac{2 - 2 \cos (k h)}{k^{2} h^{2}}$	$\frac{4 k^{2} h^{2} [1 - {(- 1)}^{m} \cos (k h)]}{{(k^{2} h^{2} - m^{2} π^{2})}^{2}}$	−0.375 +0.330kh
$J_{0}^{2} (k h / 2)$	$2 J_{m}^{2} (k h / 2)$	−0:851 +0.48kh
$\frac{π}{2 k h} [C^{2} (\sqrt{\frac{2 k h}{π}}) + S^{2} (\sqrt{\frac{2 k h}{π}})]$	$\frac{π}{8 \| k h \|} {\| Erfi [\frac{(1 + i) (2 k h - m π)}{4 \sqrt{k h / 2}}] + Erfi [\frac{(1 + i) (2 k h + m π)}{4 \sqrt{k h / 2}}] \|}^{2}$	−1.184 +0.521kh

Groove shape	P(η, ρ)	$\tilde{P} (k, ρ)$
	$\frac{z}{2} [δ (η - h) + δ (η + h)] + (1 - z) δ (η)$	$z \cos (k h) + (1 - z)$
	$\frac{1 - z}{2} [δ (η - z h) + δ (η + z h)] + \frac{1}{2 h} rect (\frac{η}{2 z h})$	$(1 - z) \cos (z k h) + \frac{\sin (z k h)}{k h}$
	$Re (\frac{1}{π \sqrt{{[h \sin (π {\hat{ρ}}_{T})]}^{2} - η^{2}}})$	$J_{0} [k h \sin (π {\hat{ρ}}_{T})]$
	$\frac{1}{8 h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})} rect [\frac{η}{8 h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})}]$	$\frac{\sin [4 k h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})]}{4 k h {\hat{ρ}}_{T} (1 - {\hat{ρ}}_{T})}$

${\tilde{P}}_{0} (k) = a_{0} / 2$	a_m	m_10%
$\frac{1 + \cos (k h)}{2}$	$[1 - \cos (k h)] \frac{2 [1 - {(- 1)}^{m}]}{m^{2} π^{2}}$	1
$\frac{2 - 2 \cos (k h)}{k^{2} h^{2}}$	$\frac{4 k^{2} h^{2} [1 - {(- 1)}^{m} \cos (k h)]}{{(k^{2} h^{2} - m^{2} π^{2})}^{2}}$	−0.375 +0.330kh
$J_{0}^{2} (k h / 2)$	$2 J_{m}^{2} (k h / 2)$	−0:851 +0.48kh
$\frac{π}{2 k h} [C^{2} (\sqrt{\frac{2 k h}{π}}) + S^{2} (\sqrt{\frac{2 k h}{π}})]$	$\frac{π}{8 \| k h \|} {\| Erfi [\frac{(1 + i) (2 k h - m π)}{4 \sqrt{k h / 2}}] + Erfi [\frac{(1 + i) (2 k h + m π)}{4 \sqrt{k h / 2}}] \|}^{2}$	−1.184 +0.521kh

Understanding the effects of groove structures on the MTF

Abstract

1. Introduction

2. Pupil Difference Probability Density

3. One-dimensional periodic groove structures

4. Parallel groove patterns

5. Circular groove patterns

6. Concatenation of errors, randomness, and relation to other measures

7. Conclusion

A. Derivation of PDPD for circular grooves

B. Fourier coefficients for ${\tilde{P}}_{1}$

Funding

Acknowledgments

References and links

Supplementary Material (3)

Cited By

Figures (5)

Tables (2)

Equations (57)

Optics Express

Abstract

1. Introduction

2. Pupil Difference Probability Density

3. One-dimensional periodic groove structures

4. Parallel groove patterns

5. Circular groove patterns

6. Concatenation of errors, randomness, and relation to other measures

7. Conclusion

A. Derivation of PDPD for circular grooves

B. Fourier coefficients for P˜1

Funding

Acknowledgments

References and links

Supplementary Material (3)

Cited By

Figures (5)

Tables (2)

Equations (57)

Optics Express

B. Fourier coefficients for ${\tilde{P}}_{1}$