Refractive error sensing from wavefront slopes

Abstract

The problem of measuring the objective refractive error with an aberrometer has shown to be more elusive than expected. Here, the formalism of differential geometry is applied to develop a theoretical framework of refractive error sensing. At each point of the pupil, the local refractive error is given by the wavefront curvature, which is a 2 × 2 symmetric matrix, whose elements are directly related to sphere, cylinder, and axis. Aberrometers usually measure the local gradient of the wavefront. Then refractive error sensing consists of differentiating the gradient, instead of integrating as in wavefront sensing. A statistical approach is proposed to pass from the local to the global (clinically meaningful) refractive error, in which the best correction is assumed to be the maximum likelihood estimation. In the practical implementation, this corresponds to the mode of the joint histogram of the 3 different elements of the curvature matrix. Results obtained both in computer simulations and with real data provide a close agreement and consistency with the main optical image quality metrics such as the Strehl ratio.

Introduction

The problem of obtaining the objective refractive error from the wave aberration of an eye is still subject of study (Thibos, Hong, Bradley, & Applegate, 2004). It is possible to predict refractive errors from wave aberration data (Guirao & Williams, 2003) using image quality metrics (Cheng, Bradley, & Thibos, 2004), and several image quality metrics show a good correlation with subjective visual acuity (Applegate, Marsack, & Thibos, 2006). Accurate visual acuity predictions can also be computed directly through schematic models of early visual processing (Dalimier, Pailos, Rivera, & Navarro, 2009; Nestares, Navarro, & Antona, 2003; Watson & Ahumada, 2008). Even though these different metrics and models show a high predictive capability, they present some problems still unsolved. 

A major difficulty is the different nature of the aberrometric measurements (i.e., wavefront or wavefront derivatives) and image quality metrics. Most image quality metrics are derived either from the point spread function (PSF) or from its Fourier transform (the optical transfer function, OTF). Passing from the wavefront W to the PSF (or OTF) involves two nonlinear stages: (i) the complex exponential ( e ^ikW ) to form the complex pupil function and (ii) the squared modulus (∣ ∣ ²) of the far field diffracted amplitude. As a result, in order to find the prescription (best correction) from image quality metrics, one has to solve a nonlinear optimization problem. This nonlinear search is not trivial since it involves a 3-dimensional space of three unknown variables: sphere S, cylinder C, and axis θ ₀ (or S, C ₀ and C ₄₅.) It is well known that nonlinear optimization methods in a multidimensional space can stagnate in local minima and usually require departing from an initial guess not too far from the solution (global minimum.) It is always possible to adopt a sequential strategy (in a similar way as in clinical refraction), searching first for the best sphere, and then for axis and magnitude of cylinder; later adding the spherical equivalent, etc. 

Nevertheless, an even more challenging issue is that different image quality metrics may eventually give different results. Thus there is not a single unified criterion even for objective refraction metrics. Furthermore, subjective refraction might depend on the visual task (Jansonius, 2010). In particular, contrast is a key factor for optimal visual performance, and hence the visual Strehl ratio, VSR (Guirao & Williams, 2003), may be a convenient image quality metric. However, if one wants to optimize visual acuity, then the cut-off frequency of the visual MTF seems to be an even more appropriate metric (Dalimier et al., 2009). In fact, contrast and/or resolution are the two main criteria for image quality and are on the basis of most useful metrics. 

For objective refraction (pure image quality criterion for no particular visual task), the Strehl ratio (SR) seems an especially interesting metric. It has a twofold meaning as the peak intensity of the PSF and as the volume under the MTF. Roughly speaking, volume is proportional both to the covered area of spatial frequency plane (resolution) and to the mean height (contrast), so that SR seems a good compromise of the two main image quality criteria of contrast and resolution. The problem is that the computation of SR from the wave aberration ( W) requires nonlinear operations from the aberrometric data. These nonlinear computations are hardly invertible, so that it is difficult to retrieve the W from the PSF, which precludes a direct feedback. The goal of the present work is to provide a theoretical formulation of the problem of measuring the objective refractive error with an aberrometer. A truly satisfactory method to obtain the refractive error from aberrometric data should meet two main criteria: On the one hand, it should be based on optical (aberrometric) rather than image magnitudes and use linear calculus (as in wavefront sensing where W is obtained by integration of the gradient) and linear algebra (as in conventional refraction where diopters follow additive rules). On the other hand, it is expected that the refractive error measured with an aberrometer should provide a reasonable approximation to that predicted by the best image quality metrics. This would guarantee that both types of metrics are consistent. 

Iskander, Davis, Collins, and Franklin ( 2007) proposed a method for objective refraction from monochromatic wavefront aberrations via Zernike power polynomials, using wavefront vergence as the metric for refractive error. This method was later improved using orthogonal polynomials (Nam, Thibos, & Iskander, 2009a, 2009b). In a previous work (Navarro, 2009), a different definition of vergence was proposed. Vergence error, in diopters, of a single ray normal to the wavefront passing through the exit pupil was formulated as a symmetric 2 × 2 matrix V with three independent elements. These elements are given by the ratios of image/pupil coordinates of the rays, where the image coordinates are transverse aberrations proportional to wavefront slopes. The V matrix was derived using finite geometry of rays, and it represents a conical (ellipse, hyperbola, etc.) elementary beam associated to each ray. The conical beam is characterized by the three components of refractive errors S, C ₀, and C ₄₅, i.e., the beam associated to each ray is equivalent to that produced by the combination of spherical and cylindrical lenses. That formulation of refractive and vergence error has several advantages; an important one is that the vergence error matrix for each pupil sample (ray) can be directly obtained from the raw aberrometric data. The matrix formulation has some similarities with the well-known power vectors notation (Thibos, Wheeler, & Horner, 1997), but as it is shown below, the matrix formulation can be computationally powerful and mathematically rigorous as it permits to establish a direct link with the fundamental quadratic forms of differential geometry (wavefront curvature). 

Once we have a formulation for the refractive error of a single ray (sample), a key issue is how to determine the global refractive error from that set of individual samples. In the previous work (Navarro, 2009), the proposed metric was the simple pupil average, but the mean has the potential problem of being strongly affected by outliers. In other words, a few highly aberrated rays (high refractive errors) can potentially bias the pupil average, thus yielding a suboptimal estimation. 

In the present work, differential geometry is applied to determine the local curvature of the wavefront. This formalism adds further support to the idea that locally, at a given point, the wavefront can be described as an infinitesimal conic wavefront, represented by a 2 × 2 symmetric matrix (second fundamental quadratic form of differential geometry) representing the wavefront curvature W″. The elements of V and W″ are proportional to quotients between image and pupil coordinates, but the difference is that these coordinates are finite in V and differential in W″. These two alternative metrics are compared both with numerical examples and with real data. 

Differential geometry of the wavefront and refractive error

Let us start reviewing some basic definitions and properties of the differential geometry of wavefronts. The wave aberration can be defined as a surface formed by points with coordinates x = ( x, y, z = W( x, y)) (see Figure 1). The geometrical aberration of a ray passing through point x, y is given by the normal to that surface. In 3 dimensions (3D), the unitary normal vector to the surface W( x, y) − z = 0 is the normalized gradient:  

$\mathbf{N}=\frac{\nabla W}{\parallel \nabla W\parallel }=\frac{(W_X^{\prime },W_Y^{\prime },-1)}{\sqrt{1+W_X^{\prime 2}+W_Y^{\prime 2}}}\approx (W_X^{\prime },W_Y^{\prime },-1),$

(1)

where

$W_X^{\prime }$

= ∂ W/∂ x,

$W_Y^{\prime }$

= ∂ W/∂ y, and ∇ W is the gradient. For moderate aberrations, the wavefront slopes (Born & Wolf, 1993) are small (of the order of 10 ⁻³ or 10 ⁻²) so that

$W_X^{\prime 2}$

+

$W_Y^{\prime 2}$

≪ 1 (usually below 10 ⁻⁴) and hence

$\sqrt{1+W_X^{\prime 2}+W_Y^{\prime 2}}$

≈ 1. For this reason, it is more common to work with the 2D gradient (

$W_X^{\prime }$

,

$W_Y^{\prime }$

). The wavefront curvatures are determined from the two fundamental quadratic forms of differential geometry (Kreyszig, 1991):  

$\mathbf{I}=\left(\begin{array}{cc}1+W_X^{\prime 2}& W_X^{\prime }{W}_Y^{\prime }\\ W_X^{\prime }{W}_Y^{\prime }& 1+W_Y^{\prime 2}\end{array}\right),$

(2a)

and  

$\mathbf{I}\mathbf{I}=\frac{1}{\sqrt{1+W_X^{\prime 2}+W_Y^{\prime 2}}}\left(\begin{array}{cc}{W}_{XX}^{″}& W_{XY}^{″}\\ W_{XY}^{″}& W_{YY}^{″}\end{array}\right)\approx \left(\begin{array}{cc}{W}_{XX}^{″}& W_{XY}^{″}\\ W_{XY}^{″}& W_{YY}^{″}\end{array}\right)={\mathbf{W}}^{″}.$

(2b)

The same approximation

$\sqrt{1+W_X^{\prime 2}+W_Y^{\prime 2}}$

≈ 1 is valid here. W″ is symmetric since

$W_{YX}^{\prime \prime }$

=

$W_{XY}^{\prime \prime }$

. The Gaussian curvature, the product of the two main curvatures, is given by the ratio of the determinants of these two matrices:  

$K=k_1{k}_2=\frac{Det(\mathbf{I}\mathbf{I})}{Det\left(\mathbf{I}\right)}=\frac{W_{XX}^{″}{W}_{YY}^{″}-W_{XY}^{″2}}{{(1+W_X^{\prime 2}+W_Y^{\prime 2})}^2}\approx W_{XX}^{″}{W}_{YY}^{″}-W_{XY}^{″2}=Det\left({\mathbf{W}}^{″}\right).$

(3)

This is an intrinsic and invariant property of the surface. Then, at any point of the wavefront we have a normal vector ( Equation 1), which defines the direction of the trajectory of the ray associated to the wavefront at that point. At the same point, we also have an associated infinitesimal conic wavefront. The sign of the Gaussian curvature determines the type of conic: elliptical for K > 0; hyperbolic for K < 0; and for K = 0, it can be either flat or parabolic. This conic wavefront is approximately defined by W″. Its orientation and main curvatures can be found by straightforward diagonalization of the W″ matrix. It is worth remarking that W″ is an approximated expression, valid as far as the norm of the gradient is close to 1. This is true for the usual magnitude of aberrations in normal human eyes, but for extremely aberrated eyes the exact expressions should be used instead. The Gaussian curvature provides a nice intrinsic metric of the local refractive error, but it has the drawback of being quadratic (nonlinear).

Figure 1

 

Diagram of an aberrated wavefront (red) and the ideal reference sphere (dashed black). At any point x, there are both a normal ( N) ray and two principal curvatures defining an infinitesimal conic wavefront. The vergence error of the ray can be estimated from the proportionality of the two right triangles.

 

Refractive error of an infinitesimal sample of a wavefront

The correction of refractive error consists of finding the combination of prisms, spherical and cylindrical lenses, which modify the wavefront to cancel that error. If we add these correcting elements, then the new wavefront will be W _corrected = W + P + S + C. For now on, we shall assume that the wavefront can be expressed as a Zernike polynomial (

$Z_n^m$

) expansion  

$W(\rho ,\theta )\approx \sum _{n,m}{c}_n^m{Z}_n^m(\rho ,\theta ),$

(4)

with coefficients

$c_n^m$

, and that the correcting elements (lenses and prisms) are ideal, free from higher order aberrations. It is well known that a perfect prism P only contributes with first-order terms (tilts) P =

$p_1^{-1}{Z}_1^{-1}$

+

$p_1^1{Z}_1^1$

, whereas perfect lenses (spheres S or cylinders C) only add second-order terms L = S + C =

$l_2^0{Z}_2^0$

+

$l_2^2{Z}_2^2$

+

$l_2^{-2}{Z}_2^{-2}$

to the wavefront. This means that the refractive correction can be described in terms of constant deviations (slopes or gradient) for prisms, and constant curvatures (second derivatives) for lenses. Thus, the gradient produced by a prism P′ can be obtained by computing the first derivatives of

$Z_1^1$

and

$Z_1^{-1}$

:  

${\mathbf{P}}^{\prime }=\left(\begin{array}{c}{P}_X^{\prime }\\ P_Y^{\prime }\end{array}\right)=\frac{1}{R}\left(\begin{array}{c}2p_1^1\\ 2p_1^{-1}\end{array}\right),$

(5)

where R is the pupil radius. Gradient is given in tangent (dimensionless) units. (In clinics, the convention is to pass to prismatic diopters by expressing tangents as percent units.) In a similar way, the second fundamental form II ( Equation 2b) provides an unambiguous way to represent the refractive correction with lenses L:  

${\mathbf{L}}^{″}=\left(\begin{array}{cc}{L}_{XX}^{″}& L_{XY}^{″}\\ L_{XY}^{″}& L_{YY}^{″}\end{array}\right)=\frac{1}{R^2}\left(\begin{array}{cc}4\sqrt{3}{l}_2^0+2\sqrt{6}{l}_2^2& 2\sqrt{6}{l}_2^{-2}\\ 2\sqrt{6}{l}_2^{-2}& 4\sqrt{3}{l}_2^0-2\sqrt{6}{l}_2^2\end{array}\right)\approx \left(\begin{array}{cc}S+C_0& C_{45}\\ C_{45}& S-C_0\end{array}\right)=\mathbf{D},$

(6)

where S is sphere, C and θ ₀ are magnitude and axis of astigmatism, C ₀ =

$\frac{1}{2}$

Ccos2 θ ₀, and C ₄₅ =

$\frac{1}{2}$

Csin2 θ ₀ (Navarro, 2009). Equation 6 is a fundamental and direct result of the differential geometry of the wavefront, which links the elements of the dioptric power matrix D to the curvature of the wavefront for ideal ophthalmic lenses. When the coefficients

$l_2^m$

are in micrometers and the pupil radius is in millimeters, then D and L″ are in diopters. The main curvatures of the resulting conic wavefront are k ₁ = S + C/2 and k ₂ = S − C/2, the mean curvature is S, and the Gaussian curvature is K = S ² −

$C^2/4$

. 

In the presence of higher order aberrations (HOAs), both gradient and curvature of the wavefront vary across the pupil and Equations 5 and 6 do not apply anymore. Nevertheless, the Gauss theorem implies that for an infinitesimal area around a given point ( x, y) the wavefront is free from HOA, and then a local correction is possible making P′ = −∇ W( x, y) and L″ ≈ − W″( x, y). This means that the local refractive error is not constant but may change across the pupil (Charman & Walsh, 1989). The prism cancels the gradient and the lens cancels the curvature of the wavefront, but only locally for a given wavefront sample around an infinitesimal area around point ( x, y). If we substitute L″ by − W″( x, y) in Equation 6, then we obtain the power of the correcting spherical and cross cylinder lenses:  

$S=-\frac{1}{2}(W_{XX}^{″}(x,y)+W_{YY}^{″}(x,y));C_0=-\frac{1}{2}(W_{XX}^{″}(x,y)-W_{YY}^{″}(x,y))\mathrm{a}\mathrm{n}\mathrm{d}{C}_{45}=-W_{XY}^{″}(x,y).$

(7)

These are the basic expressions of refractive error sensing. Note that this notation differs from that used in clinics. Here S is the pure spherical component (average power of the lens across meridians), which is the spherical equivalent in optometry. The clinical sphere is Sc = S − C/2. In addition, here C is positive, whereas negative cylinders are more frequent in clinics. In most aberrometric methods, the raw data are the components of the 2D gradient (

$W_X^{\prime }$

,

$W_Y^{\prime }$

) (i.e., transverse aberrations). Then, wavefront sensing consists of integrating the gradient to obtain the wavefront W. Under this theoretical framework, refractive error sensing consists of differentiating (partial derivatives). Therefore, an aberrometer can be a double sensor of wavefront (integral) and refractive error (partial derivatives). Note that Equation 7 describes the local refractive error, so that the elements of D are not constant over the pupil but are functions of the coordinates S( x, y), C ₀( x, y), and C ₄₅( x, y). When the wavefront is expressed in terms of Zernike coefficients

$c_n^m$

, then the particular case of refractive correction canceling second-order aberrations (

$l_2^0$

= −

$c_2^0$

;

$l_2^2$

= −

$c_2^2$

:

$l_2^{-2}$

= −

$c_2^{-2}$

) may be far from being optimal. This correction is optimal only when HOAs are small; under the Marèchal criterion, this occurs when the root mean square (RMS) wavefront error is below λ/14. However, this is not the case for human eyes, which usually show higher values of HOA (Porter, Guirao, Cox, & Williams, 2001). 

Equation 6 permits to obtain an expression for the RMS metric of refractive error. For matrices, the RMS of the matrix elements is called the Frobenius norm:  

$\parallel {\mathbf{L}}^{″}\parallel =\frac{4\sqrt{3}}{R^2}\sqrt{2l_2^{0^2}+l_2^{-2^2}+l_2^{2^2}}\approx \parallel \mathbf{D}\parallel =\sqrt{2}\sqrt{S^2+C_0^2+C_{45}^2}.$

(8)

This expression can be useful to obtain the RMS refractive error, which is not proportional to the second-order RMS wavefront error, due to the

$\sqrt{2}$

weight of the defocus coefficient. 

Wavefront curvature versus vergence error

In a previous work, the vergence error of a ray passing through the point x, y at the exit pupil plane and intercepting the image (retina) at X′, Y′ was formulated as a 2 × 2 matrix (Navarro, 2009). The geometries of the wavefront and a ray normal to it are illustrated in Figure 1. The formal derivation of V was described in detail in (Navarro, 2009; Equations 7, 8, and 9; now x, y are used instead of ξ, η and slopes

$W_x^{\prime }$

,

$W_y^{\prime }$

instead of X, Y) which is given by the following equation:  

$\mathbf{V}=\left(\begin{array}{cc}\frac{W_X^{\prime }}{x}& \frac{1}{2}(\frac{W_X^{\prime }}{y}+\frac{W_Y^{\prime }}{x})\\ \frac{1}{2}(\frac{W_X^{\prime }}{y}+\frac{W_Y^{\prime }}{x})& \frac{W_Y^{\prime }}{y}\end{array}\right).$

(9)

The formulation applies directly to ray tracing type of aberrometers, but the formulation for the great majority of aberrometers, such as the popular Hartmann–Shack (H–S) wavefront sensors, is similar (Moreno-Barriuso & Navarro, 2000). 

Equations 6 and 9 suggest that the refractive error (except prism) can be expressed as a 2 × 2 symmetrical matrix, which represents a conic wavefront produced by a combination of ideal spherical and cylindrical lenses. When pupil coordinates ( x, y) are normalized to the pupil radius, then the refractive error is given by V/ R ² or W″/ R ². It is worth remarking the parallelism existing between vergence error of a finite ray V and the local curvature W″ of the wavefront. The elements of V are quotients between finite magnitudes (

$\frac{W_X^{\prime }}{x}$

, etc.) and those of W″ are quotients between differential versions of the same magnitudes

$\frac{\partial W_X^{\prime }}{\partial x}$

=

$W_{XX}^{\prime \prime }$

. For the particular case of pure second-order aberrations W =

$Z_2^m$

, then curvatures and vergences are equal because

$W_{XX}^{\prime \prime }$

=

$\frac{W_X^{\prime }}{x}$

,

$W_{YY}^{\prime \prime }$

=

$\frac{W_Y^{\prime }}{y}$

, etc. One potential problem of using V as a metric for vergence error is that first-order terms (tilts) have a nonzero contribution, and hence the presence of prism could bias the estimation of refractive error. Therefore, it is necessary to cancel first-order contributions before computing V. 

In the presence of HOA, when the wavefront is expressed as a Taylor series of monomials, then for each monomial, vergence error and wavefront curvature are proportional (but not equal). The coefficient of proportionality is n − 1 (order of the first derivative of the monomial). In general, when W is given in terms of a Zernike polynomial expansion, it is straightforward to show that vergence error and curvature elements are related through a sparse triangular matrix. Therefore, under the presence of HOA, matrices W″ and V provide different but related metrics of refractive error, one based on differential geometry and one based on finite geometry, respectively. Vergence error has two important drawbacks: it has singularities at the coordinate axes ( x = 0 or y = 0) and it is not shift-invariant. It is possible to remove singularities (Nam et al., 2009a, 2009b) and subtract prism, but these bad features can compromise robustness in practice. One advantage is that it can be computed directly from the (centroided) raw data (Navarro, 2009). On the other hand, wavefront curvatures are well defined for all points in the wavefront and are directly derived from the rigorous formalism of differential geometry. Nevertheless, both estimations of refractive error, W″ and V, involve different approximations that explain why they tend to differ as HOAs increase (see Discussion section). Both metrics could be well suited to analyze old experimental measures of longitudinal aberrations, which basically consisted of measuring longitudinal refractive errors (in diopters) for different pupil radius (to obtain the longitudinal spherical aberration, LSA; Ivanoff, 1953; Koomen, Tousey, & Scolnik, 1949) or for different wavelengths (longitudinal chromatic aberration, LCA; Bedford & Wyszecki, 1957; Charman & Jennings, 1976). The main difference between early and modern aberrometry (Howland & Howland, 1977; Liang, Grimm, Goelz, & Bille, 1994; Navarro & Losada, 1997) is that today's aberrometers measure transverse displacements of spots, instead of longitudinal shifts (LAs), and then compute the wavefront error (by some numerical integration). Both V and W″ permit us to link transverse to longitudinal aberrations, through the gradient or curvature of the wavefront, respectively. 

Local and global refractive errors

Vergence error and wavefront curvature provide local measures of refractive error under differential and finite geometry approaches, respectively. In general, the refractive error will not be homogeneous across the pupil due to the effect of HOA. Thus, some criterion is needed to pass from that inhomogeneous distribution of refractive error to a global single representative value. Another potential problem is that, in the presence of HOA, V and W″ are nearly proportional but not identical, thus they could eventually give different measures of refractive error. Both issues are analyzed next. 

Refractive correction of the chief ray

An especially interesting case is the refractive correction necessary to have the paraxial image focused onto the retina and free from astigmatism. Speaking more rigorously, this is the particular case of canceling the refractive error of the chief ray, i.e., correcting for the very central area of the pupil. In this context, paraxial means that the chief ray is considered to be the axis connecting pupil center and image, but such axis can differ from the optical or visual axes. Obviously, vergence error V is not properly defined for the chief ray characterized by x = y = 0, but we can use the curvature. Thus, correcting the refractive error of the chief ray means that both gradient (prism) and curvature (lens) must be zero: ∇ W _corrected(0, 0) = 0 and

${\mathbf{W}}_{\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{r}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}}^{\prime \prime }$

(0, 0) = 0. The later condition also implies canceling the intrinsic (Gauss) curvature K(0, 0) = Det(

${\mathbf{W}}_{\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{r}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}}^{\prime \prime }$

(0, 0)) = 0. This “paraxial” correction is an especially interesting case as it may be a good first approximation to the refractive error for normal eyes under full photopic conditions (pupil diameters ≤3 mm) when HOAs are expected to be low or moderate. In addition, as we will see next, both gradient and curvature of the chief ray are unaffected by the presence of aberrations with angular frequencies m > 2 (trefoil, tetrafoil, etc.). Furthermore, we will see that each element of the refractive correction is only affected by its specific angular frequency m. 

It is straightforward to see that there are two ways to compute the prism correction for the chief (or central) ray, either by directly using P′ ⁰ = −∇ W(0, 0), or as the pupil average P′ ⁰ = −〈∇ W( x, y)〉; both give the same value:  

${\mathbf{P}}^{\prime 0}=\frac{1}{R}\left(\begin{array}{c}-2c_1^1+4\sqrt{2}{c}_3^1-6\sqrt{3}{c}_5^1+\dots \\ -2c_1^{-1}+4\sqrt{2}{c}_3^{-1}-6\sqrt{3}{c}_5^{-1}+\dots \end{array}\right).$

(10)

The lens correction for the chief ray will be L″ ⁰ = − W″(0, 0), and after applying Equation 7, we obtain the powers of the correcting lenses:  

$S^0=-\frac{1}{R^2}(4\sqrt{3}{c}_2^0-12\sqrt{5}{c}_4^0+24\sqrt{7}{c}_6^0+\dots ),$

(11a)

 

$C_0^0=-\frac{1}{R^2}(2\sqrt{6}{c}_2^2-6\sqrt{10}{c}_4^2+12\sqrt{14}{c}_6^2+\dots ),$

(11b)

 

$C_{45}^0=-\frac{1}{R^2}(2\sqrt{6}{c}_2^{-2}-6\sqrt{10}{c}_4^{-2}+12\sqrt{14}{c}_6^{-2}+\dots ).$

(11c)

Therefore, only terms with m = 0 contribute to sphere, only m = ±1 terms contribute to horizontal/vertical prism, and only m = ±2 terms contribute to cross cylinders, 0° and 45°, respectively, which is a nice invariant property. This correction is equivalent to the “paraxial” metric used by Cheng et al. ( 2004). In other words, the chief ray correction has several equivalent interpretations: Cancelation of the wavefront curvature at the pupil center, or cancelation of the Seidel (monomial expansion) defocus, or cancelation of the second-order Taylor series (around the origin) approximation to the wavefront. In general, central (either chief ray or paraxial) correction might be far from being optimal, except for prism. Prism is a particular case as it produces a global displacement of the image and affects to all rays (samples) equally, so that central value is equal to the pupil average. 

Histogram of refractive errors

In presence of HOA, refractive errors are different for each ray or pupil sample. Therefore, a complete correction is not possible with standard lenses, unless one applies a specific local value for each pupil position. Certainly, this would require deformable mirrors (Liang, Williams, & Miller, 1997), phase plates (Navarro, Moreno-Barriuso, Bará, & Mancebo, 2000), or some kind of inhomogeneous microlens array. This is not the case of clinical refraction where one can only apply a single constant correction for the whole pupil. There are different possible strategies to obtain the best correction. For instance, in optical design it is common to correct spherical aberration either for zonal ( R/

$\sqrt{2}$

) or marginal ( R) rays. Analogously, one might choose a central, zonal, marginal, or any pupil area and correct refractive error for that particular sample, but all these possible choices seem arbitrary. A more efficient strategy is to choose a value for the correction that minimizes the resulting set of refractive errors. When the probability distribution of errors is symmetric (zero skewness), the average is a good strategy as it maximizes the number of samples corrected. If one subtracts that average, the resulting distribution of refractive error over the pupil will have zero mean. That criterion (pupil average) was proposed before for the vergence error 〈 V〉 (Navarro, 2009). However, for skewed distributions, the mean is not necessarily the value most likely. In that case, the mean is highly sensitive to outliers (i.e., a few samples with high values may strongly bias the mean) and needs a high number of samples to have an accurate estimation. 

Let us start analyzing the one-dimensional histogram of refractive error S in the presence of spherical aberration, SA. Different possible combinations of fourth- and sixth-order SAs are compared in Figure 2, for the two proposed metrics, wavefront curvature (upper panel) and vergence error (lower panel). The first case corresponds to pure positive

$c_4^0$

= +0.5 μm fourth-order SA (plain green line). In the other two cases, sixth order is negative,

$c_6^0$

= −0.125 μm, but the fourth-order SA can be either negative,

$c_4^0$

= −0.5 μm (blue open circles), or positive,

$c_4^0$

= +0.5 μm (red filled circles). The last case corresponds to the typical signs (fourth, positive and sixth, negative) that one finds in human eyes. The other possibilities (

$c_6^0$

= +0.125 μm) are mirror symmetric and are not included. The chosen total amount (RMS) of spherical aberration is about 1 wavelength to have a strong enough effect. The histograms were computed from 3209 samples on a square grid over the pupil and correspond to the upper left element of the matrix (

$W_{XX}^{\prime \prime }$

and V _XX , respectively) scaled to micrometers of defocus coefficient

$c_2^0$

. In this example with rotational symmetry, both diagonal elements are equal, whereas cross terms are expected to have no contribution to the global error. In the case of wavefront curvature (upper panel), we can see three different types of distribution. In the presence of

$c_6^0$

< 0, when the fourth order is positive (red), then the histogram is fairly symmetric, basically free from skewness. However, when fourth and sixth orders have the same sign, then the histogram is strongly asymmetric, then showing an exponential-like shape (maximum skewness). In the case of vergence error, both histograms show an almost identical asymmetric distribution. This suggests that W″ is discriminating between both cases, whereas V seems insensitive. Another difference is that the tails of W can be much longer (mainly for the blue asymmetric histogram). This reflects the fact that curvatures for highly aberrated rays can be up to 5 times higher than vergence errors. When the sign of fourth and sixth orders is the same, then the higher order dominates for marginal rays so that W″ ≈ ( n − 1) V (with n = 6), as expected. Despite those differences, in the presence of sixth-order SA, the peaks (statistical mode) of the four histograms are practically coincident. The exact peak positions as well as other statistical moments are listed in Table 1 and also scaled to micrometers of defocus

$c_2^0$

. The case of pure fourth-order SA (green) is totally different. For V (lower panel), the histogram shows a rippled but essentially flat distribution between −2 μm and +2 μm, with no clear peak. The mean is near zero except for a small positive bias of 0.138 μm. The curvature histogram is also flat within the same interval but shows a long tail toward positive values. This affects the mean, which is 1.93 μm.

Figure 2

 

Histogram of wavefront curvature (upper panel) and vergence error (lower panel) for the two cases when fourth- and sixth-order spherical aberrations have the same (blue) or opposite (red) sign.

Table 1

 

 

Statistical moments of refractive errors for the three examples of Figure 2 (in micrometers of wavefront defocus). For a pupil diameter of 5.264 mm, these values also correspond to diopters.

Table 1

 

Statistical moments of refractive errors for the three examples of Figure 2 (in micrometers of wavefront defocus). For a pupil diameter of 5.264 mm, these values also correspond to diopters.

	Chief	Mean	Median	Mode	Standard deviation	Skewness	Kurtosis
Wavefront curvature
$c_4^0$ = 0.5; $c_6^0$ = 0	−1.94	1.93	1.38	–	2.73	0.70	−0.33
$c_4^0$ = 0.5; $c_6^0$ = −0.125	−3.08	0.80	0.89	0.97	1.72	−0.17	−0.54
$c_4^0$ = −0.5; $c_6^0$ = −0.125	0.79	−3.06	−1.01	0.95	4.97	−1.52	1.66
Vergence error
$c_4^0$ = 0.5; $c_6^0$ = 0	–	0.02	0.03	–	1.11	−0.015	−1.19
$c_4^0$ = 0.5; $c_6^0$ = −0.125	–	−0.16	0.32	0.92	1.17	−0.89	−0.48
$c_4^0$ = −0.5; $c_6^0$ = −0.125	–	−0.2	0.26	0.92	1.19	−0.84	−0.58

 

From these examples, it follows that we can have a wide variety of histograms of refractive error; they may be flat, or present strong skewness, or even show a totally asymmetric exponential-like distribution. In the last case, the mean is a poor estimate of the best refractive correction, as it is strongly influenced by the long tail of the skewed distribution, i.e., a few marginal rays with high refractive error may bias the estimation. In general, there are four possible estimators of the most likely value, including the value at the origin (chief ray, only available for the wavefront curvature), all of them listed in Table 1. For low skewness (second and fourth rows in Table 1), mean, median, and mode (when available) are close as expected. In that case, the median (close to the average of these three estimators) is probably the more likely estimate. For the other four cases, skewness is high, and then the mean is quite different from mode (the median is always in between). In that case, the most probable value, the statistical mode, seems a more robust and efficient estimate than the mean. In the presence of sixth-order SA, the peak is always clearly defined and contrasted, high and narrow. For pure fourth-order SA, however, rather than a peak there is a wide (from about −2 to +2 μm) rippled plateau, and hence it is difficult to localize the mode. The width of the (thin) peak or (wide) plateau is probably related to the depth of focus, as it is discussed below. On the other hand, the central (chief ray) value does not seem to be a good estimator in general, but in some cases (see the third row in Table 1) it may be not too far from the peak. In conclusion, the mode is a good candidate to estimate the best global refractive correction. When one applies that correction (by subtracting the mode peak), then the effect is to shift the histogram placing the peak at zero. This is equivalent to maximizing the number of points (or pupil area) corrected. When instead of a narrow peak there is a plateau, one could try to find the center (or mean) of the plateau, but obviously the correction will be less precise (higher uncertainty) and less effective (lower number of rays corrected). Probably, this uncertainty reflects a real physical limit. When the histogram shows a wide plateau, there is a wide range of equally likely best corrections, which suggests that in those cases aberrations are increasing the depth of focus, which is analyzed next. 

Through focus image quality and refractive error metrics

In the above examples ( Table 1), different metrics may give estimations of refractive error with differences as high as 6–7 wavelengths (same number of diopters for a 5.3-mm pupil diameter) when SA is about 1 λ. As an initial test of the above metrics for refractive error, let us compare them against the Strehl ratio (SR) as a standard reference metric for objective image quality. The through focus SR is a function of sphere S, but it would be equally possible to compute the two through cylinder SR as well. As clinical refraction involves three variables, in general, we will have a function SR( S, C₀, C₄₅). Nevertheless, in the above examples of pure SA we only have the variable S. Figure 3 shows the through focus SR( S) for the same cases, but now the fourth-order SA is always positive and the sixth-order coefficient can be negative, zero, or positive. Interestingly, one can find some parallelism between these curves and the histograms of Figure 2. The case of pure fourth-order SA (continuous green line) shows ripples on a mirror symmetric sort of plateau. There are two main symmetric peaks at about −1.2 μm and +1.2 μm. Thus, there is no unique solution, but two opposite best SR corrections. Furthermore, such a waving plateau suggests a sort of focal segment with a broad depth of focus of more than ±1.2 μm. That focal segment, and its associated ambiguity, is consistent with the plateau of the green histograms in Figure 2. When sixth-order SA is added, then the symmetry is broken and there is a marked (high and narrow) single SR peak at minus or plus 0.88 μm, which suggests a much narrower depth of focus. For this amount (about 1 λ) of total SA, the focal zone spans more than 4 μm (about 7 λ), but the full width at half height of the SR peak (depth of focus) is 10 times lower, about 0.4 μm. The two curves for positive (dotted red line) and negative (dashed blue line) sixth-order SAs are mirror symmetric. On the contrary, reversing the sign of the fourth-order SA has no effect and the resulting curve is almost identical. In other words, the sign of the fourth order has little effect on SR( S), i.e., the through focus SR is basically insensitive to the sign of

$c_4^0$

. However, as shown in Figure 4, the corresponding PSFs are different in the two cases of positive and negative fourth-order SAs (and negative sixth order).

Figure 3

 

Through focus Strehl ratio for three cases of pure spherical aberration.

Figure 4

 

Through focus PSF for the cases when fourth- and sixth-order spherical aberrations have the same sign (upper panel) and opposite sign (lower panel).

 

We can compare the PSFs in the upper and lower panels of Figure 4, which correspond to the cases of negative and positive fourth-order SAs, when sixth order is negative. These two cases correspond to the same dashed blue line in Figure 3, but the PSFs look quite different. Each PSF corresponds to the optimal defocus correction for different criteria as indicated by the labels. Label Wrms, minimum RMS wavefront error, corresponds to

$c_2^0$

= 0. Note that the histograms in Figure 2 and the origin of the through focus curves in Figure 3 correspond to that case. Both SR peaks are obtained for the same defocus coefficient, −0.88 μm, despite the fact that the corresponding PSFs (Strehl label) in the upper and lower panels of Figure 4 look quite different. The common feature of the Strehl PSFs is that both display the highest intensity of the central peak but also the dimmer surrounding disk (maximum contrast). The PSFs corresponding to the mode of the histogram (Mode label) look quite similar to the Strehl ones, since the defocus (see Table 1) is close to the SR peak. The minimum RMS wavefront error (Wrms label),

$c_2^0$

= 0, provides a somewhat smaller disk but a lower intensity of the central peak (lower SR). The minimum disk corresponds to the circle of least confusion, given as the minimum RMS geometrical radius of the spot diagram ( R( W′)rms label) commonly used in optical design. The defocus values are 2.05 μm and −0.80 μm for the upper and lower panel, respectively. This metric seems to have a poor consistency with the SR. The chief ray (paraxial) curvature provides, in one case (upper panel), a good PSF close to the Strehl and Mode ones, but in the other case (lower panel), the PSF shows a poor optical quality and the defocus coefficient differs by almost 4 μm from the optimal value.

Table 2

 

 

Subjective (upper row) and objective refractions for different metrics, subject RN. The column RMS contains the distance (rms) between aberrometric and SR peak. The two right columns are the resulting image (SR) and optical (statistical mode K and percent frequency) quality metrics.

Table 2

 

Subjective (upper row) and objective refractions for different metrics, subject RN. The column RMS contains the distance (rms) between aberrometric and SR peak. The two right columns are the resulting image (SR) and optical (statistical mode K and percent frequency) quality metrics.

Metrics	S	C	Axis	RMS	SR	( K) %
Subjective	0.25*	0.5	0°	0.19	4.4%	(0.21) 53%
SR peak	0.12	0.5	169°	0	6.7%	(−0.14) 45%
Mode( W)	0.22	0.53	1°	0.17	4.0%	(0.36) 49%
Mode( V)	0.29	0.62	179°	0.24	3.5%	(0.44) 43%
W″(0, 0) (parax.)	0.13	0.92	167°	0.39	2.7%	(0.25) 32%
Wrms	−0.11	0.39	8°	0.46	2.7%	(0.29) 34%

 

Note: *Spherical equivalent.

 

The PSFs in the upper panel ( Figure 4), when both fourth and sixth orders are negative, are characterized by having a narrow central core (peak) plus a dim and extended surrounding disk. There is a sort of analogous pattern in the long tails of the blue histogram in the upper panel of Figure 2, i.e., points with high curvatures at the pupil correspond to points placed far from the center in the PSF. As the frequency of these highly curved points is low in the histogram, then the corresponding contribution to the PSF is small and hence the disk they form is dim. In the lower panel (

$c_4^0$

> 0), the PSFs tend to be more compact with a narrower and brighter surrounding disk. This corresponds to the red histogram in the upper panel of Figure 2, which is nearly symmetric and with much shorter tails. 

In summary, there are three metrics in Table 1, SR peak, Mode( W), and Mode( V), which provide similar values of best focus. They seem nearly invariant under sign reversal of the fourth-order SA. The standard W _rms is basically insensitive to the presence of SA due to the orthogonality of Zernike polynomials. The other metrics, chief ray curvature ( W″(0, 0)), average of curvature and vergence (〈 W″〉 and 〈 V〉), or the disk of least confusion ( R( W′)rms), show a low robustness. It is interesting to note that the average and central curvatures show a sort of mirror symmetry. The absolute value of the difference of these two metrics seems to provide a rough estimate of the focal segment (total width of the dashed blue line in Figure 3). 

Implementation and results

The proposed formulation was applied to the experimental data of 10 subjects taken from a recent study on the experimental validation of a custom Bayesian (ideal observer) model of visual acuity (Dalimier et al., 2009). For this group, the range of pupil diameters goes from 4 mm to 6.5 mm, with an average of 5.15 ± 0.87 mm; the RMS HOA varies from 0.052 μm to 0.60 μm, with an average of 0.25 ± 0.19 μm. The implementation details and results are analyzed first for the right eye of subject RN (the author), for a 6.5-mm pupil diameter, RMS HOA of 0.50 μm, and clinical subjective refraction of −1 + 0.5 × 0°. Note that wavefront sensing and refractive error sensing are similar in practice. The Zernike coefficients (up to seventh order here) are obtained by least square fitting of the raw data. Then, one can reconstruct either the wavefront, or its gradient vector or its curvature matrix. The distribution of the three elements of the curvature matrix W″ is plotted as a 3D cloud of points in Figure 5. This distribution shows a region with a higher density of dots, even though that region has an odd (curved star-like) shape. There are several possible approaches to obtain the most likely combination of the three variables. A direct way is to compute the point with the highest density (3D mode) volume. A somewhat simpler strategy is to compute the three statistical modes of the marginal 1D distribution of each variable, as shown in Figure 6. The two approaches (single 3D or three 1D marginal distributions) give the same result only when the variables are statistically independent, but this is not the case for the three second derivatives. The current implementation consists of two stages: A first initial estimate is obtained by computing the peaks of the three 1D histograms of 3209 samples with 200 1D bins ( Figure 6). Then, to save computing time, the 3D histogram is obtained but only for a 9 × 9 × 9 square grid around the point defined by the three 1D modes. For the 3D histogram, the sampling interval was coarser, typically δ = 0.15 D along the three dimensions, and even higher, 0.2 D, when the RMS HOA was greater than 0.4 μm. The resulting density inside that 3D grid was typically low. For this reason, overlapping bins were used, computing the number of points within boxes of double length, ± δD, so that the effective volume (and the average frequency) is 8 times higher. In this way, the mode bin always contained 100 points at least (for low RMS HOA the mode bin could contain more than 1000 samples, that is more than 30%). In addition, to keep enough samples within the central boxes, bin overlapping helps to smooth ripples, which is especially important in 3D histograms.

Figure 5

 

Frequency distribution of the 3 different elements of the wavefront curvature (in diopters) for subject RN, 6.5-mm pupil.

Figure 6

 

Marginal 1D histograms of the elements of the wavefront curvature (upper panel) and vergence error (lower panel) matrices. Subject RN, 6.5-mm pupil.

 

Figure 7 shows an alternative way to view these refractive errors, now as maps (Charman & Walsh, 1989) representing the spatial distribution of S, C ₀, and C ₄₅ across the pupil. The little arrows signal the value with the highest frequency in each map. The resulting refractive errors for this subject, obtained with different metrics, are listed in Table 2. The first row corresponds to the subjective (clinical) refraction. The value of S (*) is set to the spherical equivalent of astigmatism, because sphere was compensated during the aberrometric measurements by means of a Badal lens. The second row corresponds to the Strehl ratio peak that is our standard reference for objective image quality. The other rows correspond to four different aberrometric metrics. The second to fourth columns include S, C and axis (diopters and degrees), respectively. The fifth column contains the RMS difference of refractive error with respect to the reference (SR peak), computed as ∥ D∥/

$\sqrt{2}$

( Equation 8). This RMS refraction difference was used to sort the different metrics (except for the subjective one). The SR value and the Gaussian curvature, K, are listed in the sixth and seventh columns, as metrics of image and optical (wavefront) qualities, respectively. The K column contains two numbers, value (in square diopters) and percentage of pupil area contained within the mode peak (height) of K. It is interesting that after correcting S and C, the percentage of homogeneous area of the pupil having the same K is relatively high (>40% in most cases) despite the poor image quality (SR < 7%). Interestingly, the maximum homogeneity (percentage) of K is obtained for the subjective refraction. However, the lowest (absolute) value of K is obtained for the SR peak correction. Attending to these two objective quality metrics SR and K, the SR peak seems the best correction. Subjective refraction and mode ( W″) give similar quality values, probably because the RMS difference is also similar. (Nevertheless, a direct comparison of subjective and objective metrics is risky in this context due to the different nature of objective and subjective measurements.) The next metric, the mode of the vergence error, has a somewhat higher RMS difference, 0.24 D, and hence it provides a lower performance but also yields a high homogeneity of Gaussian curvature (43%). The other two metrics, chief ray curvature W″(0, 0) and minimum RMS wavefront error, show a higher difference (0.39 and 0.46 D, respectively), thus providing a worse performance in terms of both SR and K.

Figure 7

 

Distribution of refractive error across the pupil computed from wavefront curvature W″ for subject RN, 6.5-mm pupil. The arrows signal the statistical mode.

 

The results for the group of 10 eyes are summarized in Figure 8. The bars represent the average difference between the refraction corresponding to the four aberrometric metrics and the reference, the SR peak; for the three variables S, C ₀, and C ₄₅ and for the global RMS distance. The error bars correspond to standard deviations for this group of subjects. All metrics give low average differences with the reference, always below 0.07 D for S, C ₀, and C ₄₅. However, the standard deviations (error bars) are always several times greater than the means, and hence, none of these metrics seems to have a statistically significant systematic bias. In this context, error bars are especially relevant, since a small error bar indicates a high consistency between a given metric and the reference, and conversely. In this sense, all these metrics except the mode of wavefront curvature ( W″) show a low consistency with image quality (SR). The average RMS difference (Euclidean distance) is similar for the RMS wavefront error and for the curvature of the chief ray W″(0, 0), about 0.25 D, but the error bar indicates that in some subjects, the RMS difference can be around 0.5 D. The mode of the vergence error performs somewhat worse, especially for C ₄₅. Perhaps the computation of the cross term of vergence as an average, ( V _xy + V _yx )/2, might be contributing to amplify errors. On the contrary, the mode of W″ provides the closest agreement with the reference for all variables and for the global RMS difference, average below 0.15 D. The global error bar, about 0.1 D, suggests a high degree of consistency between wavefront curvature and Strehl ratio. It is worth remarking that most of this difference can be explained by the quantization of the variables used in the implementation. In fact, the sampling intervals were 0.05 D for the SR and 0.15 D for the 3D mode of W″, respectively, which are similar to the average difference between these two metrics.

Figure 8

 

Average difference between the aberrometric metrics and the SR peak for the three variables and for the global RMS difference.

 

Discussion

In the preceding sections, refractive error was formulated in terms of the local curvature of an infinitesimal area of the wavefront. (The particular case of prism was also formulated in terms of wavefront gradient.) W″ was compared to vergence (or vergence error) V proposed previously (Iskander et al., 2007; Navarro, 2009). These two alternative metrics V and W″ are equivalent for low-order aberrations, but HOAs make that they tend to differ at points with higher curvatures. One relevant result of this study is that the differential geometry approach, W″, is superior, both in theoretical and practical aspects. Its formulation is more rigorous (straightforward derivation from the fundamental forms of Equations 2a and 2b), more exact, and robust. Another important difference becomes patent if we observe the PSFs of Figure 4 and then compare the two panels of Figure 2. (Note that these histograms correspond to the PSFs with label Wrms.) The histograms of vergence (lower panel) are basically equal for positive and negative fourth-order SAs, whereas the PSFs look totally different. On the contrary, the histograms of W″ adequately describe such important difference: The long tail of the blue histogram in the upper panel ( Figure 2), with high negative curvatures, corresponds to the more peripheral area of the dim wide disk of the PSF in the upper panel of Figure 4. Conversely, the nearly symmetric red (close circles) histogram shows much shorter tails in accordance with the narrower (and brighter) disk of the PSF in the bottom panel of Figure 4. Therefore, there seems to be a clear correspondence between the histogram of W″ and the PSF. This is an essential property, which seems lost in the case of vergence histograms. To finish this comparison, the computation of the average of cross vergence, ( V _xy + V _yx )/2, also involves a practical implementation issue. The green histogram in the lower panel of Figure 6 has a peak broader than that in the upper panel. This broadening tends to increase the uncertainty and potential bias in the estimation of the more likely value. In fact, this problem seems to contribute to the bias and the long error bar for the mode( V) metric for C ₄₅ in Figure 8. 

In conclusion, the wavefront curvature W″ shows to be superior to the vergence error, as it provides a more complete, robust, and exact metric of the refractive error. In the context of wavefront sensing, W″ provides a natural formulation for refractive error sensing, simply computing partial derivatives instead of (or in addition to) the integral. Another advantage of W″ is that it contains other metrics used before in the literature (Cheng et al., 2004). In particular, the central or chief ray curvature W″(0, 0) is equivalent to the Seidel (monomial) refractive errors (defocus and astigmatism) or to the second-order Taylor series approximation to the wavefront. Globally, its histogram provides a representation of optical quality, which keeps a high fidelity to the image quality (PSF). In fact, W″ is a Hessian matrix and its determinant (Gaussian curvature, K) can be used to estimate the geometrical PSF, as it provides a measure of the ratio between the densities of rays (per unit of area) at the pupil and at the image (Smith, 1966). Gaussian curvature is possibly a good metric of the optical quality of a wavefront (see the right column in Table 2). Anyway, the correspondence between wavefront curvature and image quality is an important property, because W″ is a pure aberrometric (optical quality) quantity, which can be obtained through linear operations (derivatives) from the wavefront, or directly from the aberrometric measurements (wavefront gradient). The formalism of differential geometry is widely applied in corneal topography, but the gradient and curvature of the corneal surface are much higher, and then, the approximations in Equations 2b and 3 do not hold. Thus, the computation of the two main curvatures from I and II becomes harder and time consuming. For this reason, maps of curvature along some predetermined directions (axial, tangential) are usually given instead of the main curvatures. Some topographers can also provide the Gaussian curvature (González, Hernández-Matamoros, & Navarro, 2008), which is especially interesting because it is an intrinsic invariant property of the surface. 

The proposed refractive error sensing method is local for each pupil sample. Thus, some additional analysis is needed to pass from the set of local values to a single global estimate of the best correction. A statistical approach was proposed here in which the most likely best correction is estimated as the most frequent among all the wavefront samples across the pupil, i.e., the statistical mode. As discussed in the Local and global refractive errors section, the estimation of the most likely value in a skewed irregular 3D distribution is not trivial. The current implementation involves two stages. An initial estimate is obtained by computing the modes of the three marginal histograms, and then the mode of the joint 3D histogram is computed using a limited number of overlapping bins around that initial estimate. The distribution and size of the cubic bins were adjusted so that the mode bin contained a significant percentage of samples and moderate to low ripples. The results improved by adapting the sampling interval, δ, to the amount of RMS wavefront error. Both the differential geometry formalism and the statistical analysis of local refractive error are general and powerful, and hence further improvements and generalizations (postsurgical eyes, bifocals, etc.) are possible, but this will be subject of future work. There are also other potentially interesting alternative implementations that may deserve further study. It is possible to obtain the refractive error by direct differentiation of the raw data, thus avoiding data fitting. Differentiation (or integration) can strongly amplify the measurement noise, but the use of a statistical maximum likelihood estimator (mode) of the best correction may alleviate the effect of noise amplification. 

An interesting result is the totally different effect of pure fourth-order spherical aberration or when it is combined with sixth-order SA (see Local and global refractive errors section). Pure

$c_4^0$

causes an ambiguous solution for the best objective correction, both in terms of the histograms of refractive errors and in terms of through focus Strehl ratio. As shown in Figure 3, pure

$c_4^0$

seems to increase the depth of focus. The through focus SR shows lower values and presents two mirror symmetric peaks (bimodality), but only when higher order coefficients of SA are zero. Any mixture of fourth- and sixth-order SAs seems to break that symmetry providing a single, high, and narrow peak (best objective focus). Recent experiments in which HOAs are modified with adaptive optics show that HOAs, and

$c_4^0$

in particular, expand the depth of focus (Rocha, Vabre, Chateau, & Krueger, 2009). The examples and real cases analyzed in the present study suggest a general trend: when HOAs are high and most orders n and angular frequencies m are present, histograms tend to be broader and flatter, and the mode peak becomes less pronounced, or even one can observe secondary peaks, i.e., bimodal or multimodal histograms. In these cases, the through focus SR also tends to show two or more peaks and an increased depth of focus. As a result, there is a higher uncertainty in the position of the best image plane. Indeed, the attainable precision in finding the most likely refractive correction seems limited by the depth of focus, which can be expanded by the presence of HOA (Rocha et al., 2009). Recent studies with (Chen, Kruger, Hofer, Singer, & Williams, 2006; Gambra, Sawides, Dorronsoro, & Marcos, 2009) and without (López-Gil, Fernández-Sánchez, Thibos, & Montés-Micó, 2009) adaptive optics suggest that HOA, and SA in particular, could play a crucial role in accommodation. Accommodation response seems to improve by canceling HOA (Gambra et al., 2009), which suggests that finding the best image plane becomes easier for the visual system. These studies when combined with appropriate theoretical framework and models may help find the type of information and metrics used or preferred by the visual system. Experiments by Rocha et al. ( 2009) suggest that SA (fourth spherical aberration) induces a defocus of about 2.62 D/ μm of induced SA, which means that the visual system is apparently choosing one (the positive) of the two peaks (see green curve in Figure 3). There are two mechanisms (at least) that may help the visual system to choose between these two options. One is the Stiles–Crawford effect (SCE), optical apodization. It is straightforward to implement the SCE in the objective metrics and verify that it can break the symmetry and disambiguate, but the SCE has a significant influence only for big pupils (>5–6 mm). The SCE alone cannot explain experimental results for small pupils. A possibly more important factor is the neural transfer function (NTF). The band-pass neural response helps attenuate the contrast loss due to wide and dim surrounding disks. (Thus, between the two best PSFs (Strehl) in Figure 4, both having the same Strehl ratio, the visual system is expected to choose the one in the upper panel, which has a more spread, but dimmer, disk and hence higher visual Strehl ratio (VSR).) Neural response is what distinguishes performance in a visual task from objective optical quality. However, in the particular problem of measuring objective refraction with an aberrometer, the subject's NTF is usually unknown, and some nominal NTF is used to compute visual performance metrics such as the VSR. As a consequence, these visual metrics have a predictive value but cannot be considered as true measurements. On the contrary, from aberrometric raw measurements W′ (gradient), one can compute the second derivatives W″ for local refractive error sensing, in addition to the integral W, for standard wavefront sensing. In this context, objective refraction appears as a new application for curvature sensors (Roddier, 1988). 

Acknowledgments

This research has been supported by the Spanish CICyT, Grant FIS2008-00697, and Red Española de Optometría (Ministerio de Ciencia e Innovación, SAF2008-01114-E). 

Commercial relationships: none. 

Corresponding author: Rafael Navarro. 

Email: [e-mail]. 

Address: ICMA, Consejo Superior de Investigaciones Científicas and Facultad de Ciencias, Universidad de Zaragoza, Plaza San Francisco s/n, 50 009 Zaragoza, Spain.