1. Introduction

The amount of residual aberrations is the most important criteria for the quality of an imaging optical design. It is the goal of the optical designers and optical design software to optimize the surface curvatures, element thickness, refractive index and other parameters in such a way as to limit the residual amount of aberrations to an acceptable level. An algorithm for aberration calculation is often the central part during optimization. For a rotational symmetric system, any optical design software such as Zemax or CodeV, allows to apply third order aberration theory and e.g. the Seidel aberrations can easily and fast illustrate the contribution from each surface [1,2]. However, for non-rotational symmetric systems, aberration theory requires much more mathematical effort, such as general vector aberration methods or nodal aberration theory [3,4]. These theories require a deep mathematical treatment of the general aberrations and therefore the mathematical algorithms and the interpretation of the aberrations will be quite complex. Consequently, it is desirable to search for alternate algorithms which can easily compute aberrations for all surfaces just like Seidel aberrations do in rotational system.

The method of phase space is widely employed in both quantum mechanics and classical mechanics [5], where the position and velocity of a particle define the location in phase space. In optics, the position and velocity of a particle correspond to the position and angle of a tracing ray [6,7]. In this view, rays can be represented in the coordinates of phase space. The mathematical similarities between geometrical optics and mechanics consequently result in similar transformation behavior. The employment of phase space has proven to be an effective method for illumination design [8,9] in applications such as LED lighting and solar concentrators [10]. For imaging systems, as described in [11,12], aberrations correspond to non-linear deformations in phase space and phase space can be an alternative access to visualize aberrations in freeform systems [13]. Therefore phase space methods seem to be also promising to quantitatively determine aberration contributions, even for systems without any symmetry. As a first step in this direction, it is the goal of this paper to prove the method for rotational symmetric systems, where a comparison to other methods is possible.

Rotational symmetric elements are the basis of many optical systems, since spherical and rotational symmetric aspherical optical surfaces can be easily fabricated and tested. In terms of aberrations they are relatively simple to study. As mentioned above, third order aberration theory is standard in analysis of such systems. However for rotational symmetric systems, also aberration contributions including all orders can be calculated in an exact way via so-called Aldis calculus. Aldis theorem is published in 1964 by Cox [14] as the exact method to calculate the transverse aberration contribution of each surface for rotational symmetric systems. Brewer and Welford further developed Aldis method in papers [15] and [16], respectively. We will employ these methods as a gold standard to calculate aberrations in rotational symmetric systems.

With the following chapters an alternate approach to surface aberration contributions based on phase space optics is developed. As will be shown, the phase space approach is not so much based on a mathematical treatment but rather developed from a more general analysis of paraxial and non-paraxial ray-propagation behavior within the optical system. We will in depth explain the method and compare to the results of the Aldis calculus in order to verify its correctness.

2. Phase space dynamics for imaging systems

2.1. Basic definitions in phase space

Many well written introductions to the concept of phase space in optics are existing in the literature [6,7], therefore only some important relations and quantities will be explained within this paper. Let us start with a one-dimensional system analysis, i.e. for example only meridional rays are considered. Any meridional ray is completely defined by two quantities, namely the distance to the optical axis x and the angle θ to the optical axis, as illustrated in Fig. 1(a). If we furthermore weight the ray direction with the refractive index n of the propagation medium, we can define$\tilde{u}=n\cdot \sin (\theta )$. In this case x and $\tilde{u}$ represent the phase space trajectory of the ray and the differential area$dH=dx\cdot d\tilde{u}$corresponds to the Lagrange invariant, which is conserved within an optical system. For the general situation of an oblique ray, four parameters$(x,y,\tilde{u},\tilde{v})$ are required, where x and y are the lateral ray height in x and y direction, and $\tilde{u},\tilde{v}$are the index weighted direction cosines. The four-dimensional phase space volume element$dH=dx\cdot dy\cdot d\tilde{u}\cdot d\tilde{v}$in this case corresponds to the differential etendue. In the following however it is required to establish a relation to the ABCD-matrix-formalism, therefore we will use u = n tan(θ) as a slightly modified direction variable for the one-dimensional case, respectively for the general case the direction tangents, which can easily be calculated from the normalized direction vector (L, M, N) of any ray via:

 

Fig. 1 Relation of ray position and direction to phase space coordinates, for (a) 2D phase space and (b) 4D phase space.

Download Full Size | PDF

It is easy to record the phase space behavior from real ray tracing, since all that is required is a monitoring of the ray positions (x, y) and ray direction (L, M, N), respectively (u, v) as the rays propagate through the optical system. Thus the phase space behavior of an optical system, e.g. a simple lens as illustrated in Fig. 2, can be analyzed with the help of any optical design software or ray-tracer. For this paper Synopsis CodeV is employed for the ray-trace analysis and phase space recording.

 

Fig. 2 Ray trace illustration of the ray pair used to calculate aberrations and corresponding phase space coordinates.

Download Full Size | PDF

2.2. Calculating aberrations from phase space

In a short conference contribution [13] we have already sketched the basic idea of extracting aberrations from phase space analysis and also presented results for individual rays. In this paper a more detailed description of the mathematical method is given, and the full ray fan and spot diagrams are analyzed. For simplicity we will first describe the method for meridional rays only, i.e. a two-dimensional phase space (x, u). Later however we will show that the extension to the full four-dimensional case is straightforward and we will use the four-dimensional method to calculate the aberrations.

As in Aldis calculus, which is reviewed later, it is the core of the method to compare the real ray behavior to the linear system behavior. To do so, we slightly modify the optical model by additionally introducing plane dummy surfaces D_i and D_i’ directly before and after each surface of the system, as illustrated in Fig. 2.

By using the ray-tracer we can, for any ray, determine the real ray-data (x_i, u_i) at the dummy surface D_i preceding the optical surface and also the real ray-data (x_i’, u_i’) on the dummy surface D_i’ after the real optical surface.

In the process of calculating aberration contributions the real propagation of the rays in phase space will be compared to the ideal linear aberration-free propagation. The linear propagation is achieved if all refractive powers in the system are replaced by ideal thin elements at the position of the surface vertex, which corresponds to the position of the dummy surfaces. The corresponding linear ray-propagation step for any input ray can be represented by the well-known ABCD matrix formalism:

Here the bar above the variables is used to indicate the linear nature. For spherical surfaces of curvature c_i the relevant matrix S_i can be calculated directly via:

Here n_i-1 and n_i is the refractive index of the medium before and after the surface.

For the free propagation step towards the next system surface i + 1, which is located at some distance d_i, it can be found from simple geometry for any ray under the angle θ to the optical axis:

And since we defined u = n tan(θ) the corresponding ABCD matrix for propagation from surface i to i + 1 is represented by the usual form

Note that the above propagation matrix is exact also for non-paraxial rays and corresponding to real ray-tracing. Therefore for free propagation in-between surfaces there will be no difference between the ideal linear and the real ray. In other words the free propagation does not introduce any aberrations, as expected.

In contrast, the real propagation through the spherical surface i, corresponding to the exact ray-tracing from dummy surface D_i to D_i’, does not coincide with the linear transformation according to Eq. (2) and (3). In consequence every input ray (x_i, u_i) at dummy surface D_i results in two output rays: The linear output ray $\left({\overline{x}}_i\text{'},{\overline{u}}_i\text{'}\right)$after the surface, which can be calculated via the ABCD formalism. And the real ray (x_i’, u_i’), resulting from exact interaction with the spherical surface. Thus for any input ray, two differential rays are resulting at each surface and an aberration vector can be defined as:

This aberration vector describes the angular and lateral difference between an ideal and a real ray introduced by the surface, and thus corresponds to the introduced aberration. Both rays are illustrated in Fig. 2 in real space and in phase space.

2.3. Propagation of the aberration vector to the image plane

However this differential vector describes the aberration within the system. In order to compare to Aldis aberration theory we need to propagate it to the image plane, since the transverse aberrations from Aldis calculus are only calculated at the image. This propagation will transform the aberration vector (δx_i, δu_i), emerging after surface i, to an aberration vector (δx_i^I, δu_i^I) at the image plane. The final lateral position component δx_i^I then corresponds to the transverse aberration introduced by surface i at the image plane.

A simple but approximate method to perform this propagation is to apply a sequence of linear transformations, corresponding to the forward propagation through the rest of the optical system to the image plane, as illustrated in Fig. 3. It only requires the sequence of ABCD matrixes to be applied to the aberration vector, i.e.

 

Fig. 3 Propagation of the differential aberration vector to the image.

Download Full Size | PDF

As mentioned this method is approximate, since only the linear system behavior after the surface i is considered. However it is a very good approximation, since the error made just corresponds to the aberrations on top of the aberration vector; therefore it is of second order. The main effect will be correctly described by the linear propagation.

2.4. Total transverse aberration and oblique rays

The above described process results in the transverse aberration contribution δx_i^I at the image plane resulting from only one surface i. However if we repeatedly apply the method to all of the surfaces in the system, the contributions of all the surfaces can be worked out. The total transverse aberration of the ray on the image δx_i^I is the sum of the transverse aberrations of all the surfaces:

So far the method was explained for the simple meridional case only, since it is much easier for the reader to follow and visualize. However the generalization to the full four-dimensional situation is quite straightforward. In the general case the oblique ray is described by four parameters (x, y, u, v), as defined in chapter 2.1. The real ray-tracing is just the same and is performed by any ray-tracer. Only the paraxial ABCD matrix needs to be increased to 4 × 4 dimension, which for rotational symmetric systems is straightforward and results in:

With these definitions the method can just be performed in the same way as in the above two-dimensional case, but now yielding a four-dimensional aberration vector (δx_i, δy_i, δu_i, δv_i), which needs to be propagated to the image. There the spatial components (δx_i^I, δy_i^I) correspond to saggital and tangential transverse aberration contributions.

3. Aldis calculus

Aldis equations allow an exact computation of the aberration contributions within rotational symmetric systems. The calculus is based on a mathematical treatment of the paraxial ray data, quite similar to the method employed above, in combination with the exact ray-tracing equations. Combining both propagation equations, results in a mathematical exact sum for the transverse ray aberration of an arbitrary ray on the image plane. The individual terms in the sum correspond to the individual surface aberration contributions. A good description of the mathematical derivative of Aldis calculus is given by Brewer [15], and his set of equations is also used now for calculation and comparison. The equations as derived in Brewer’s paper are the following:

Here we have mostly adopted the nomenclature of Brewer and U_i and V_i correspond to the transverse aberration contribution of the i^th surface in x and y direction, (X_i, Y_i, Z_i) represent real ray coordinate at the intersection on the surface, (L_i, M_i, N_i) are the optical direction cosines of the ray after being refracted by surface i, u_i is the angle of axial paraxial ray, h_i means the height of axial paraxial ray, n_i is the refractive index after surface i. Generally the subscript I denotes the image plane and Y_I^’ is the Gaussian image height. These ray tracing parameters are illustrated in Fig. 4, in which the red line represents the real ray, the blue line is the paraxial chief ray related to the Gaussian image height, and the green line corresponds to the paraxial marginal ray.

 

Fig. 4 Ray-data used by Aldis theorem.

Download Full Size | PDF

As mentioned above these equations are exact and include all high order aberration contributions. Same as in the phase space approach the equations require real ray-tracing data, namely the ray intersections with the real surfaces i and the ray-directions. So the required input data are quite similar.

4. Results

4.1. Triplet lens example and calculation methods

Brewer in his paper employed the example of a triplet system [15, 17]. In his paper he accurately describes the system and calculates aberration contributions for a number of rays. Now nearly the same triplet system is used to compare the result of the phase space method to the exact results according to Aldis equations. The layout of the triplet lens is illustrated in Fig. 5 for three field points. The exact lens data of the system are listed in Table 1. The stop of the system is located at the back-side of the second lens, however for ray-tracing we virtually start the rays in the corresponding entrance pupil location, which is listed in the table. The wavelength was set to λ = 587.56nm. For the aberration analysis we use the same setting as Brewer, which are an entrance pupil diameter of 4.348mm and a maximum field of angle of 16.5°. The effective focal length is 10.00mm.

 

Fig. 5 Layout of triplet lens used for the aberration analysis.

Download Full Size | PDF

Table 1. Lens data of triplet lens. The refractive index n is listed behind the glass name.

View Table | View all tables in this article

Aldis theorem and the differential phase space method is now employed on this system to calculate the individual aberration contributions for a number of rays. For both methods paraxial ray-tracing is required, which can be obtained from the ABCD-formalism, as explained in section 2.2. The ABCD matrix values can be directly calculated form the system data above. In addition real ray-tracing is required, where it is necessary to record the ray/surface intersections (X_i, Y_i, Z_i) for the Aldis calculus, and the dummy surface intersections (x_i, y_i) and (x_i’, y_i’) for the phase space calculus. In addition the ray directions (L, M, N), respectively (u, v) need to be recorded at each surface.

From this set of data the resulting surface contributions to transverse ray aberrations is calculated. The result according to Aldis is readily obtained from Eq. (9), whereas the approximate phase space result is obtained from paraxial propagation of the differential aberration vector, as described in section 2.2 to 2.4.

4.2. Transverse aberration contributions of several typical rays

Using both methods the aberration contributions are calculated, surface by surface, for the triplet lens. To compare the results several single rays are chosen, which were also used by Brewer in his paper. These rays are the tangential marginal ray Fig. 6(a), the chief ray Fig. 6(b), a sagittal skew ray with field of angle 16.5° and entrance pupil coordinate (1, 0) for Fig. 6(c) and Fig. 6(d).

 

Fig. 6 Aberration contributions for (a) Marginal ray, (b) Chief ray, (c) Δx of skew ray, (d) Δy of skew ray.

Download Full Size | PDF

In Fig. 6 the resulting transverse aberration contributions at the image plane are plotted for each surface within the triplet lens and also the sum is illustrated. The numerical data of these aberrations are also presented in Tables 2-5, where the label Aldis corresponds to data achieved by Aldis exact equations and PS are data calculated by the approximate differential phase space method. It should be mentioned that the values calculated by Aldis equations are consistent with the results of Brewers paper of 1976 [15].

Table 2. Marginal ray

View Table | View all tables in this article

Table 3. Chief ray

View Table | View all tables in this article

Table 4. Δx of skew ray

View Table | View all tables in this article

Table 5. Δy of skew ray

View Table | View all tables in this article

From Fig. 6 it can be concluded that the individual values obtained from PS and also the sum are quite close to the exact aberration contributions calculated by Aldis. This is also confirmed by the number in Tables 2-5. Consequently, we may conclude that the phase space method can achieve transverse aberration contributions to a very good approximation.

4.3 Transverse aberration curves and spots contributions in the image

The computation and comparison of several rays, as shown above and in [13], has proven the excellent performance of the phase space method. In order to furthermore verify the exactness of the method now a larger number of rays is analyzed. We investigate a meridional ray-fan of an on-axis image point (0°-field angle) and the sagittal and tangential ray-fans of a full-field image point (16.5°-field angle). Figure 7 describes the differences of aberrations achieved by phase space method and Aldis theorem for those rays. Each ray-fan has 201 rays, scanning the full diameter of the entrance pupil from −2.174 to 2.174. In all the four figures the approximation by phase space to Aldis is very good. This proves that phase space method is able to predict the correct final aberrations on the image plane for all rays.

 

Fig. 7 Transverse aberration curves of various ray fans in the image.

Download Full Size | PDF

This analysis is just for verification, since the transverse aberration curves on the image could have been easily obtained with any ray-tracer. However the deeper benefit of Aldis and the phase space method is to be able to calculate the contribution from each surface to the total aberrations, including higher orders. In order to visualize this and to further test the method, we compare spot diagrams in Fig. 8. Within these spot diagrams the transverse aberrations at image plane resulting from only one surface are visualized. The individual dots correspond to the aberrations of a rectangular grid of rays, sampling the full entrance pupil, for this full-field point. The blue circles and red dots indicate the ray aberrations on the image for each surface calculated by phase space and Aldis, respectively. Also the sum of all contributions is computed, which corresponds to the “standard” spot as achieved by the complete triplet system.

 

Fig. 8 Figures of the image space spot diagram, illustrating the contribution by each surface (a-f), and the sum (g).

Download Full Size | PDF

All figures present a good accordance between Aldis and phase space method. Accuracy for rays in the center is better than for those in the corner, which is consistent with the expectations that larger aberrations result in larger deviations.

5. Current limitations and extension to freeform systems

As discussed in section 2.3 the results form phase space are approximate, since the aberration vector δ_i is paraxially propagated to the image plane. A more exact propagation could be achieved if the local transformation in the neighborhood of the differential ray is determined. This requires a more detailed analysis of the phase space transformation between the surface i and the image. Technically this can be derived from a local mapping analysis of a neighboring real ray grid in phase space. Since the real ray data (x_i’, y_i’, u_i’, v_i’) after surface i are recorded, and also the corresponding data at the image plane (x_I, y_I, u_I, v_I) are known, the local linear transformation matrix abcd_I,i in the vicinity of the differential ray of interest can be determined by linear fitting the set of rays in the neighboring grid, e.g. in MatLab, in order to satisfy E.(10).

In consequence this also implies that generally the transformation matrixes used above can be determined from the ray-tracing data itself. The linear (paraxial) surface matrix S_i, can be determined from a ray grid close to the optical axis, whereas the local transformation matrix abcd is determined from rays close to the differential ray of interest. Therefore no assumption about symmetry or spherical shapes, as e.g. in Eq. (3), neither knowledge about the geometrical system data is required. Thus all that is required to transfer the method to freeform systems is a “guiding” ray, defining the optical axis. At each intersection of this guiding ray with a real optical surface, orthogonal dummy planes before and after the surface need to be defined. If ray-tracing is performed on this system, the linear surface behavior S_i is determined from analyzing rays close to the guiding ray. According to Eq. (5) and (2) we can thus calculate the aberration vector. From a ray-analysis of Eq. (10) close to the differential ray of interest, we determine the local propagation matrix abcd and are able to propagate the aberration to the image.

In summary the differential ray approach as explained above at no point requires assumptions about symmetry. Therefore the method should also hold for anamorphic or freeform systems. By analyzing a freeform prism system in [13] we have already evidence that this statement is true, however a full in-depth discussion of freeform systems requires a detailed discussion and verification, which would be beyond this paper.

6. Conclusions

We have developed a method to calculate surface aberration contributions which is based on a differential ray analysis in phase space. The method is based on a comparison of real ray-tracing versus linear ray-transformation according the ABCD formalism. The resulting differential ray is propagated towards the image via a set of linear transformations, resulting in an approximate prediction of individual and total aberrations. By comparing the results to the exact theory of Aldis it was shown that the results of the phase space method are very close to the exact calculus. The comparison is showing very good agreement in terms of ray-aberration curves and also in terms of the individual spot contributions. Further improvement of the method seems possible if the underlying phase space transformations are analyzed in more depth.

So in summary we have found an alternate method to calculate aberration contributions including high orders. Moreover in the novel phase space method the assumption about rotational-symmetry of the system is not required. Therefore we think that method can be easily extended to non-rotational symmetric systems or freeform systems. Future work will investigate these problems for applying phase space in asymmetric systems.