Y. A. Kravtsov, Y. I. Orlov, Geometrical Optics of Inhomogeneous Media (Springer-VerlagBerlin, 1990), Sec. 2.9.

C. H. F. Velzel, J. L. F. de Meijere, “Sensitivity analysis of optical systems using the angle characteristic,” in Current Developments in Optical Engineering and CommercialOptics, R. E. Fischer, H. M. Pollicove, W. J. Smith, eds., Proc. SPIE1168, 164–175 (1989).

[CrossRef]

For a description of characteristic functions see, for example, H. A. Buchdahl, An Introduction to Hamiltonian Optics (Cambridge U. Press, Cambridge, 1970; reprinted by Dover, New York, 1993).

The best choice of reference surface depends on the particular representation of the characteristic function (point or angle, say) and on the application. One possibility is to choose reference surfaces that correspond to the object and the image and use a mixed characteristic. Another, relevant for determining the wave aberration function, would be to use a point representation and take the final reference surface to be a sphere centered on the ideal image location that passes through the center of the exit pupil. In either case, if the object is at infinity, the initial reference surface is chosen arbitrarily.

Note that it turns out that the following derivation is more straightforward when all three coordinates for the intermediate point variable [i.e., (x, y)] are retained. It is possible to repeat the derivations given here for a restricted characteristic function at surface k, but the results are much less compact (even for the first-order correction).

The basic equations of Hamiltonian optics can be found, for example, in Chap. 2 of Ref. 1.

The outer product of two-component vectors, v and w, say, is defined as
v⊗w:=vywyvywzvzwyvzwz.

This can be found by taking the equation of the surface, x=cy2/{1+[1-c2(κ+1)y2]1/2}, making the following replacements: x→x cos θ+y sin θ,y→y cos θ- x sin θ, and solving the resulting equation for x.

For the angle characteristic, the initial coordinate is given by yˆ=∂T/∂pˆ.In generating Fig. 5, I replaced Tby (T0+ θT1) in this equation and numerically solved for the value of p^′ that corresponds to a ray with a given value of ŷ. This procedure was repeated for a variety of values for ŷ. Note that with this procedure the first-order perturbation theory correctly predicts that there are no rays at the extreme bottom of the surface (i.e., values of yˆ approaching -1.3) when the surface is tilted. For example, by plotting ∂(T0+θT1)/∂p^yevaluated at pˆ=0,p^′= (p^y′, 0), as a function of p^y′, it is immediately obvious that this function does not extend below -1.3 for θ greater than or equal to 2°.

Note that if c2(κ+1)y2 is ever negative, the representation of Eq. (40) breaks down. Since the nominal surface used in this example is a prolate spheroid, c2(κ+1)y2 is always positive.

See, for example, A. Walther, The Ray and Wave Theory of Lenses (Cambridge U. Press, Cambridge, 1995), Chap. 17, for the derivation of Eq. (45) and the approximations implicit in that equation.

For each case, T(0, p^′) was sampled on a 128×128 square grid of values of p^′, where the circle defined by |p^′|= n′ cos 30°=0.75 just fit inside the grid [and T(0, p^′) was taken to be zero for |p^′|>0.75]. This array was then embedded in a 256×256 array of zeros, on which the fast Fourier transform was performed. This process was repeated but with T(0, p^′)sampled on a 256×256 square grid that was embedded in a 512×512 array of zeros. Since the results from these two cases were indistinguishable, it is reasonable to assume that the errors associated with the Fast Fourier transform were insignificant. Also, throughout these calculations, A(p^′) for the unperturbed system is used. That is, the perturbations are taken to affect only the phase: exp[ikT(0, p^′)]. Although it is possible to take perturbations into account when one is calculating the amplitude factor A(p^′), this is a slowly varying function that would be affected little by the perturbation.

This follows, for example, by noting that Fermat's principle requires that ∂(W0I+W0II)/∂y=0, where, in W0Iand W0II,x is replaced by f0(y).Combining this equation with Eqs. (25) and solving the resulting expression for y as a function of p̂ and p^′ gives the result of Eq. (26).