An optical analysis is developed to separate forward light scatter of the human eye from the conventional wavefront aberrations in a double pass optical system. To quantify the separate contributions made by these micro- and macro-aberrations, respectively, to the spot image blur in the Shark-Hartmann aberrometer, we develop a metric called radial variance for spot blur. We prove an additivity property for radial variance that allows us to distinguish between spot blurs from macro-aberrations and micro-aberrations. When the method is applied to tear break-up in the human eye, we find that micro-aberrations in the second pass accounts for about 87% of the double pass image blur in the Shack-Hartmann wavefront aberrometer under our experimental conditions.
© 2011 OSA
The anterior corneal surface provides the main optical power for the eye, but its optical quality relies on the integrity of the pre-corneal tear film. Normally the tear film varies a few microns in thickness across the cornea  and provides a smooth optical surface covering the microscopically rough corneal surface . The surface topography of the normal tear film, therefore, mirrors the macroscopic structure of the corneal surface but hides its microscopic structure. Evaporation between blinks causes tear film thinning [3,4] that has the potential to change the wavefront aberrations of the eye. Several recent studies using wavefront sensing have shown that higher order aberrations (HOA) increase shortly after a blink and, if blinking is suppressed, they can become significantly elevated relative to those observed during the normal post-blink period [5–10]. Elevated post-blink increases in HOAs have also been reported in cases of dry eye [11–14] and in eyes wearing soft contact lenses . These post-blink changes in tear film quality are associated with a reduction in optical quality [16–19] and loss of visual performance [18,20–22]. Eventually, thinning may cause a localized break in the tear film that exposes the irregular, mucin-covered corneal epithelial surface , introducing small scale optical disturbances that scatter light and further degrade retinal image quality .
Experimental measurements of tear film aberrations have relied primarily upon Shack-Hartmann (SH) wavefront aberrometers. However, there is reason to be concerned that SH aberrometry lacks the necessary resolution to fully capture the optical disturbances associated with the rough corneal surface exposed by tear break-up (TBU) [10,23]. Conventional analysis of Shack-Hartmann (SH) images quantifies aberrations based on displacement of a multitude of spots formed by an array of lenslets, each of which is the image of a retinal beacon of reflected light . Typically the measured spot displacements are subjected to modal analysis in which Zernike polynomials are fit to the data by the method of least-squares. This modal method of wavefront reconstruction tends to smooth irregularities in the data associated with local disturbances in the tear film and may fail to capture most of the loss in image quality associated with TBU . Moreover, it is clear from earlier reports of TBU [10,23,25], that some SH spots are blurred as well as displaced. These degraded SH spots are indicative of poor optical quality within the local region of the eye’s optics being sampled by an individual lenslet in the aberrometer. Conventional SH wavefront analysis ignores these local optical imperfections associated with TBU.
Our report aims to remedy this shortcoming by introducing a new method of data analysis that extracts additional information from the SH spots to quantify wavefront quality locally over the domain of each lenslet sampling the eye’s pupil. We argue that blurring of SH spots is evidence of unresolved micro-aberrations that exist on a fine spatial scale smaller than the diameter of the lenslet. This situation is analogous to the inverse scattering problem of determining surface characteristics from scattered light measurements . We recover information about the magnitude of micro-aberrations from measurements of the size of individual SH spots. In this way we extend the utility of the SH aberrometer as a spatially-resolved method that also quantifies the micro-aberrations and light scatter at each location in the eye’s pupil. We note that although TBU provides a motivation for the current study, the scope of its applications is not limited to TBU, but is applicable whenever fine micro-aberrations are present, for example in cataract .
Measuring the blur of individual SH spot images is straightforward using a variety of image quality metrics . However, in addition to the micro-aberrations present in the region of optics represented by each lenslet, these spot images are also blurred by image degradation on the first pass through the eye’s optics and the resolved macro-aberrations on the second pass. Our goal is to apportion these different contributions to spot blur and thereby to quantify the magnitude of micro-aberrations present in the eye. To accomplish this goal we develop a mathematical model of optical blur produced by a thin random screen using the concept of the scattered light field (section 2). To quantify spot blur, we develop a metric called radial variance which exhibits a convenient additivity property that enables a simple algorithm to estimate the blur produced by each pass in a double-pass SH aberrometer (section 3). We validated the algorithm theoretically using local perturbations to the wavefront aberrations and demonstrate its feasibility for human eyes (section 4) and investigations of tear film breakup (section 5).
2. Terminology and theoretical background
In order to describe ocular aberrations with a wavefront error map spanning the full extent of the eye's pupil, we sample the emerging wavefront at numerous locations using an array of non-overlapping sub-apertures. For an aberrometer employing a Shack-Hartmann wavefront sensor (SHWFS), the subapertures are defined by the faces of individual lenslets in the array. If the wavefront over individual subapertures is not flat, local small scale aberrations are present which will cause blur in the spot image produced by the lenslet. The classical design of the SHWFS assumes the lenslets have sufficiently small diameter that the wavefront is approximately planar over each subaperture. In this case, all of the spot images are diffraction-limited Airy-disks free of aberration-induced blur and spot displacement is a measurement of wavefront tilt over each lenslet. In effect, the wavefront for the whole pupil is being approximated by a surface tessellated by flat tiles. The SHWFS determines the slope of each tile by measuring the displacement of each spot from the optical axis of the corresponding lenslet. The wavefront is then reconstructed from slope measurements mathematically by algorithms that either retain all features of the slope data (i.e. zonal reconstruction) or that smooth the wavefront to deemphasize local irregularities (i.e. modal reconstruction) (Fig. 1 ).
Modal reconstructions using Zernike basis functions describe wavefront aberrations by a vector of Zernike coefficients associated with a series of Zernike polynomials of increasing order. The low-order aberrations (LOA) refer to Zernike polynomials of orders 0-2. The high order aberrations (HOA) could, in principle, extend from 3 to N where N is the largest order needed to fit the data exactly (i.e. a fully determined least-squares regression). However, in practice, HOA are usually specified for a smaller range 3 to R where R < N and frequently R ≤ 10. Thus a modal reconstruction of wavefront aberrations based on Zernike coefficients of LOA and HOA typically omits the very high order aberrations (VHOA) of orders R to N. The total number P of Zernike polynomials contained in orders 0 to R can be determined from the recursive sequence P(R) = 1 + R + P(R-1), where P(0) = 1. For example, P(2) = 6. By contrast, a zonal reconstruction such as integration by many orbits  includes the VHOA. Thus the difference between wavefronts reconstructed using modal and zonal methods represents the contribution of VHOA, which are known to be indicative of tear film anomalies [23,30].
In this paper we use the term “macro-aberrations” to describe wavefronts reconstructed from slope measurements. Such aberrations exist on a spatial scale that is large compared to lenslet diameter. By comparison, aberrations with a spatial periodicity of oscillation that is small compared to lenslet diameter are called “micro-aberrations”. Micro-aberrations are assumed to be square integrable and may contain random perturbations. They are too fine to be resolved by conventional reconstruction methods but nevertheless increase the size of SH spots. This division between macro- and micro-aberrations is based on the sampling density of the wavefront sensor, analogous to the Nyquist frequency of communication theory that separates the frequency spectrum into a band of well-sampled signals and a band of under-sampled signals. In principle, the range of macro-aberrations can be expanded further into the domain of micro-aberrations by reducing lenslet diameter. In practice, however, safety considerations and increased diffraction by the lenslet aperture limit the density of lenslet arrays that can be used with the human eye. For this reason we quantify micro-aberrations by taking advantage of the fact that each lenslet in a SHWFS is a miniature optical system that produces a blurred image of the retinal beacon. Thus the size of SH spots is a measure of wavefront quality over an individual lenslet that is related to the magnitude of micro-aberrations.
To measure micro-aberrations with a SHWFS, we envision fine ripples riding on the waves of macro-aberrations (Fig. 2 ). These ripples (micro-aberrations) add texture to the wavefront without affecting curvature or mean slope over the lenslet domain. From a geometrical optics perspective, the optical effect of the ripples over any given lenslet is to scatter the rays of light in a disorganized plethora of directions, thereby blurring the spot image of the retinal beacon . To model such phenomenon mathematically, we include multiplicative perturbations in the conventional complex-valued pupil function . This extra light field in the traditional pupil function due to micro-aberrations is called the Rytov approximation for light scatter in standard references . Since macro-aberrations can also blur the SH spot (if they produce sufficient wavefront curvature over individual lenslets), our goal was to develop a method for isolating the micro-aberration contribution to blurring of individual SH spots.
Our method of analysis is motivated by the statistical optical description by Goodman  of imaging through random media. We consider each lenslet of the SHWFS as a small imaging system. Phase disturbances over the face of each lenslet are caused by two independent sources: gross, structural errors of the eye that cause macro-aberrations and fine, micro-aberrations modeled below as a random phase screen. We first briefly review Goodman’s formulation in the context of the imaging system for human eyes. We then develop a metric of radial variance to quantify the blurring of SH spots due to the forward propagation of the probe beam that produces the retinal beacon, followed by additional blurring due to the reverse propagation of reflected light captured by the SHWFS.
2.1 Modeling imaging of the retinal beacon by an individual SHWFS lenslet
In the imaging optics for the human eye, random phase disturbance can occur from various sources. The optical effects of this random disturbance can be statistically modeled, and therefore changes of the optical transfer function (OTF) can be modeled by the presence of a thin random screen in the optical system . If this thin random screen is located in the pupil of the human eye, light fields propagate in a disorganized way and the scattered light fields will modify the pupil function by some complex-valued random disturbance We assume the screen is thin enough that the light arrives and leaves the pupil at the same location. The new pupil function can be written as the product The optical transfer function for is computed as the autocorrelation of the pupil function,
Unfortunately, due to lack of knowledge of the transmittance it is in general impossible to analyze Eq. (1) deterministically and so an alternative approach is needed. One way to proceed is to assume a certain physical model for and then derive properties of Alternatively, by assuming that the phase of follows Gaussian statistics for 200 eyes, we previously calculated the phase structure function of and used it to estimate the resolution limit (the Fried number) of the eye as an imaging system in terms of the volume under the OTF . However, that work used slope measurements from SHWFS and took into account only the classical macro-aberrations measured by the SH aberrometer, thereby ignoring the micro-aberration contribution over each lenslet.
Goodman’s approach  is to look for the ensemble average OTF, among all possible realizations of the random process. The expected value of Eq. (1) is a straightforward definition of but since the numerator and the denominator are correlated random variables, the computation may not be easy to carry out. Instead, we consider an alternative definition for as a ratio of expected values,Eq. (1) is a real number and the expected values from the two definitions in Eq. (1) and (2) are equal. As a second example, if the spatial periodicity of the amplitude of the pupil function in the denominator is much larger than the correlation width of then the denominator of Eq. (1) is almost constant and independent of each specific realization of the screen.Eq. (4). We use in the lower case for the Fourier transform pairs of respectively,Eq. (5): the combined PSF for any given lenslet is the convolution of the PSF due to the micro-aberrations with the PSF due to the macro-aberrations. For simplicity, we drop the overline above the average OTF and PSF and use the plain letters hereafter.
Our goal is to obtain information about micro-aberrations from the blurred SH spots in a double pass aberrometer in which scatter from micro aberrations occurs when light from a narrow laser beam enters the eye, and again within each sub-aperture when light is reflected out of the eye. Any additional scattering of light within the retinal tissue that enlarges the retinal beacon is bundled with forward scatter in our analysis. Thus the single-pass Eq. (5) is employed twice, with possibly imperfect knowledge of the pupil entry point of the laser beam on the forward pass. Recovering the full function might be possible either using Eq. (4) or Eq. (5). However, a method based on Eq. (4) involves dividing H by could be problematic when has zero values. Also, a method based on Eq. (5) requires a stable deconvolution algorithm that is currently lacking. Therefore, rather than attempting the difficult problem of recovering the full function we chose instead to recover a scalar measure of the width of this spread function. Our measure of choice is the radial variance, which has the convenient property of converting convolution into addition, making the computation more feasible while seizing a characteristic of the scatter PSF useful for characterizing micro-aberrations.
The model presented above assumes that is a stationary random process. This means that the exact pattern of fine ripples in the wavefront (Fig. 2) are unpredictable, yet their statistical properties are stable from one blink to another. Although the magnitude of micro-aberrations increases significantly as TBU evolves, and can vary dramatically from one part of the cornea to another (see Section 5), it is not unreasonable to assume statistical stability at any particular corneal location at some fixed time after blinks occur .
2.2 Additivity of radial variance
Although the detailed shape of blurred images can be computed by the convolution relationship developed in section 2.1, for the more limited goal of quantifying how much larger an image is compared to the object, it is sufficient to use a scalar measure of size that is additive under convolution. This condition is satisfied for univariate density functions by variance, defined as the mean-square deviation about the centroid . For bivariate distributions of light in polar coordinates, the corresponding measure of size is radial variance defined as the second moment of light about the centroid in the radial direction In Cartesian coordinates, radial variance is computed as
In our optics problem, convolution and variance additivity arise in two contexts. One is the standard imaging problem in which the image is the convolution of an object with the PSF of the system. Another is the scatter problem in which the average PSF of the whole system is computed by convolving the PSF due to macro-aberrations with the average PSF due to a thin random screen. Thus we have two variance relations for our imaging systems. First, the radial variance of the image is equal to the radial variance of an object plus the radial variance of the PSF of the imaging system
3. Algorithm for scatter analysis
In our double-pass optical model, the SH spots are blurred for five reasons: macro and micro aberrations encountered on the inward path through a single unknown sub-aperture, retinal scatter, and macro and micro aberrations encountered on the second pass through an array of lenslets. The macro-aberrations of the eye over each lenslet are known from conventional analysis of wavefront slopes but the other factors are unknown. Our goal is to estimate the blur due to micro-aberrations on the second pass over each lenslet, despite our incomplete knowledge of pupil entry point of the probe beam or the magnitude of retinal scatter. To do this we first estimate upper and lower bounds for the radial variance of the retinal beacon. We then use those estimates to derive lower and upper bounds for the blur contributed on the second pass by micro-aberrations.
3.1 Upper and lower bounds for the retinal image variance
In the first pass, the narrow beam of light comes into the eye and forms an image R on the retina (Fig. 3 ). Although the incoming beam does not pass through the lenslet array, we make the simplifying assumption that the beam diameter is the same as the diameter of one lenslet. In case the beam diameter is different from the lenslet diameter, a telescope system with an appropriate magnification can be used to meet the assumption. A unit-magnification afocal telescope makes the lenslet array optically conjugate to the eye’s pupil. Thus the pupil entry point of the laser beam can be specified in the coordinate system of the lenslet array. We let the PSF over the specific subaperture that the probe beam entered be Since the object is a point source, By using Eq. (7), we conclude that the variance of the retinal beacon R is equal to the variance of the PSF for the eye over the sub-aperture defined by the entering laser beam plus some additional, unknown amount of variance S contributed by retinal scatter,36]. However, in what follows we assume the exact location of the beam entering the eye is unknown. We note that even if the pupil entry point of the incoming beam is known, is unknown due to the lack of information about
In the second pass of light through the eye, the retinal beacon becomes a common object to be imaged by all the lenslets in the aberrometer. Each spot image in the SH image is the convolution of the object Rwith the PSF for the eye + -th lenslet that quantifies the combined effects of macro- and micro-aberrations. Thus, the variance of each spot is the sum of the variances of the retinal beacon and the corresponding eye + lenslet PSF (Eqs. (7) and (8)). Specifically, the radial variance of is the sum of the radial variance of and the radial variance ofEq. (10) for requires an estimate of We do this by averaging the lower and upper bounds for estimated using the following rationale.
Although we cannot measure the retinal beacon directly, the size of the smallest SH spot constrains the size of the retinal beacon. For example, since objects are always smaller than images, the retinal beacon must be smaller than the smallest SH spot, which means V(R) < min[V(grc)], where (r,c) is the location of the smallest SH spot. In fact, the retinal beacon must be even smaller because some of the variance in SH spots is caused by aberrations and diffraction in the exiting path. This blurring in the second pass is minimized when aberrations are absent, which means the maximum possible size of the retinal beacon is found by assuming the smallest SH spot was produced by a diffraction-limited second pass. This argument establishes the upper bound U for : V(R) < min[V(grc)] - V(diffraction-limited PSF) = :
A lower bound L for V(R) is established by a similar argument. The smallest retinal beacon is produced when the probe beam passes through the least aberrated part of the pupil. The image of the beacon formed by the lenslet located at the same optimal pupil location will produce the smallest SH spot. Thus the smallest SH spot will have been blurred twice by the best part of the eye’s optical system. Both passes contribute the same amount of radial variance to the SH spot, which means the retinal beacon cannot be smaller than half the size of the best SH spot. This argument establishes the lower bound L for: V(R) > 0.5*min[V(grc)] = :
3.2 Lower and Upper bounds for the PSF blur due to micro-aberrations
Given these lower bound L and upper bound U for we can estimate lower and upper bounds for of the PSF created by the second pass through each lenslet. Since the following analysis applies to each lenslet, we omit the subscripts for simplicity. Using Eq. (10), our estimates for are as follows:
We recall Eq. (8) that the radial variance of the PSF h of an optical system is given by a sum of the radial variance of the PSF for macro-aberrations and the radial variance of the PSF for micro-aberrations:
Finally, in light of Eq. (11), we conclude that the radial variance V(hs) over each lenslet is bounded below by and above by
3.3 SH simulation program
Our SH simulation program was designed to match the specifications of a laboratory SH aberrometer used to collect data from human eyes. We use this simulation program to calculate the PSF for each portion of the eye’s pupil defined by a lenslet, Only macro-aberrations are used to calculate the radial variances, and The input for this SH simulation program is a high-resolution wavefront map interpolated from measurements of wavefront slope. Over each lenslet, we extract the local wavefront phase and use it to compute a simulated spot image that has no measurement noise. As an example, in the diffraction limited case, our program faithfully represents the Airy pattern in the simulated SH image for a point-source retinal beacon.
Macro-aberrations were represented by wavefronts constructed by either modal or zonal methods (Fig. 4 ). For modal wavefronts, the macro-aberrations are represented as a series of Zernike coefficients. The wavefront aberrations can be sampled at any arbitrary rate. The local wavefront over each lenslet is unlikely to be planar, and instead contains varying wavefront slopes. In the case of the zonal wavefronts, each lenslet produces a single value of wavefront slopes. We interpolate these wavefront slopes to generate a planar wavefront surface over each lenslet (bi-linear interpolation). A collection of those local planar wavefront surfaces provides an approximation for the wavefront map. Therefore, if two neighboring local wavefront surfaces happen to have opposite wavefront slopes forming a cusp between them, a finite sampling of such wavefronts may lead double spots in one lenslet (Fig. 4).
4. Validation using simulated aberrometer images
In this section, we test the assertion that radial variances of object and PSF can be added together to produce the radial variance of the image by using the concept of a localized Gaussian disturbance of the wavefront. This validation is placed in the context of tear film breakup by modeling the localized thinning of the tear film as an application of a small drop of dryness called a “xerop” (from the Greek word xeros (ξερός) for “dry”) to the tear layer. The result is a localized shortening of the optical path length from retina to SH wavefront sensor that perturbs the wavefront aberration function. For demonstration purposes, we assume this perturbation is small enough to fit inside a lenslet face. For simplicity, we choose a Gaussian xerop,Fig. 5 ).
From these simulated images we computed in Table 1 the radial variances of object (Fig. 5(a)), PSF (Fig. 5(c)) and images (Fig. 5(d)) associated with Gaussian xerops (Fig. 5(b)). At each location, the parameters C and σ for Gaussian xerops are also shown in Table 1. For all test cases, the error is less than 4% of The presence of small errors is due to finite discretization of the related functions.
4.1 Resolution of the SHWFS detector
The additive property of the radial variances is proved analytically in Appendix and is valid over the infinite domain. In practice, the additive rule holds only approximately for several reasons. The apertures of the optical systems are of finite size whereas the proof in Appendix assumes the domain to be of infinite size. The resolution of the SHWFS detector imposes a limit since the individual spots are sampled in a discrete grid. Furthermore, the aberrometer samples wavefront slope at a finite number of pupil locations. Although each of the factors contributes to the total error of the additivity, in this section we focus on detector resolution as it applies to the additive property of the radial variances.
For the simulations shown in Table 1, we assumed that the angular resolution of the detector of spot images cast by individual lenslets is 0.11 arcmin per pixel. For example, if the focal length of the lenslet is 24 mm, this angular resolution is equivalent to a linear resolution of 0.76 micron per pixel. Figure 6 shows that the error of the additive rule is about 4% of the radial variance V(image) of the image. Since V(image) = V(object) + V(PSF), the error relative to V(object) or V(PSF) will be in fact doubled. If the detector resolution decreases, the spot image will be more coarsely sampled and the numerical value of V(image) will be less accurate. Thus the difference between the numerical and true V(image) will increase, causing the error of the additivity to increase. Computer simulations show that as detector resolution decreases to 2.7 arcmin per pixel, the error increases to 16.5% of V(image). The error relative to V(object) or V(PSF) will be double this value. Furthermore, we observed that the computation becomes less reliable.
5. Feasibility analysis for human eyes
5.1 Radial variances of the SHWFS spot images
The raw SH spot images g were acquired from human eyes by a 633 nm laboratory aberrometer . The lenslet diameter was 400 microns referenced in the eye’s pupil with the telescope magnification of 1. The camera resolution was 6.8 micron per pixel producing a 59x59 pixel window for each lenslet image. When the tear film forms a smooth surface, individual raw spot images are well focused and appear in a regular array as shown in Fig. 7(a) . When the tear film disrupts locally, the spot images become not only displaced from the regular pattern, but also blurred due to the light scatter as shown in Fig. 7(b). Tear break-up (TBU) was achieved by asking subjects to refrain from blinking for prolonged periods . The pupil radius was set at 3 mm for the analysis. Since the lenslet diameter in the SH detector is 400 micro meters, the lenslet array contains 15 rows and 15 columns of the spot images. The standard aberration analysis from the SH aberrometry with these two arrays of the spot images produces two sets of the Zernike aberration coefficients, which are shown in Fig. 8 . The TBU data exhibits larger higher order aberrations than the base line data, which is consistent with the previous literatures .
The radial variance maps of these double pass images are shown in Fig. 9(a) and 9(b). The original radial variance maps, which hold one value for each lenslet, are bilinearly interpolated so that the patterns are more easily readable. Redder colors indicate larger radial variances (arcmin2). In general, the spot images possess larger radial variances during the TBU (Fig. 9(b)) than for the baseline data (Fig. 9(a)). We argue that since any blur from the first pass combined with the retinal scatter will form a common object for the second pass, the spatial differences shown in Fig. 9(a) and 9(b) must originate from the second pass. Radial variance maps of the simulated spot images from second pass macro-aberrations are shown in Fig. 9(c) and 9(d) for the baseline data and the TBU data, respectively. We use the zonal wavefront to quantify V from macro-aberrations. In contrast to the radial variances of the double pass spot images, the radial variances of the simulated spot images for the second pass macro-aberrations are in general much smaller in magnitude (median = 13% of the radial variance for the double pass blur for both the baseline and TBU data), which implies that macro-aberrations contribute only a small amount to the blur of the raw SH spot images. Their overall values and spatial variations across the pupil are small.
5.2 Light scatter for the second pass in a double pass optical system
The spot blur analysis algorithm in section 3 enables us to isolate the light scattering effect of micro-aberrations for the second pass of a double-pass optical system. To demonstrate the idea, we apply the algorithm to the data set from section 5.1. A lower bound L and an upper bound U for the retinal image quality for the first pass will be computed from image quality of the best spot, and therefore include the retinal scatter. The bounds L and U are then averaged to calculate the light scatter contributed by the second pass, following Eq. (15). Although any of the three criteria - lower bound maps, upper bound maps, or the average maps - for the second pass can be equally suitable to quantify the light scatter, we plot the upper and lower bound average scatter only for simplicity. The radial variance of the SH spots due to the micro-aberrations is isolated by subtracting the blur introduced by the first pass and by the macro-aberrations on the second pass from the total blur observed experimentally:
The values of the radial variance maps for the micro-aberrations are, on average, always larger in the TBU data than in the baseline data (Fig. 10 ). Over the lenslets with large image blur, the average radial variances contributed by the micro-aberrations are considerably greater after TBU than before TBU. The radial variances for the micro-aberrations are smaller in scale when compared with the radial variances for the double pass total blur as shown in Fig. 9(a) and 9(b), but their spatial patterns are very similar because contributions from the radial variances for the macro-aberrations are almost uniform across the pupil. Although the raw spot images for the baseline data (Fig. 7(a)) appear to be well focused, there still remains light scatter as indicated in Fig. 10(a). The confidence interval for the scatter can be quantified by the difference of the lower and upper bounds for For each data set, although the lower and upper bounds vary across the pupil, as our algorithm in Section 3 indicates, the confidence interval is uniform. The confidence interval is only slightly larger with the TBU data (75.06 arcmin2) than with the base line data (63.3 arcmin2), which indicates that the best spots for the two data sets have similar blur in both pre and post TBU.
The post-TBU radial variance maps in Fig. 10(b) can include the scatter of the eye that may not be related to TBU . To see the change in light scatter attributable to TBU, we compute the Post-Pre TBU difference maps (Fig. 11 ). The pattern of blurred spots visible in the raw data (Fig, 9b) is similar to the pattern seen in the difference scatter map (Fig. 11). This implies that the relative scales of the radial variance maps are independent of the topographic characteristics of the patterns. Interestingly, over some regions, the radial variances become smaller after TBU (Fig. 11).
6.1 Correlation between the topographic distribution of scatter and the TBU visible in the raw SH image
The spatial distribution of forward light scatter produced by optical perturbations finer than can be resolved by our SHWFS (micro-aberrations), was characterized by quantifying spot blur in the SHWS images. Since the spot images were acquired in a double-pass optical system, the blur due to the first pass and the blur due to the macro-aberrations on the second pass have both been subtracted from the measured blur of the raw spot images in order to quantify the 2nd pass blur caused by unresolved micro-aberrations (Fig. 10). The result is not zero. Although blur from micro-aberrations increases dramatically after TBU, the blur from the macro-aberrations on the second pass contributes about 29% of the total blur on the second pass for both the baseline and TBU data.
Similar spatial patterns of the SH spot image blur were reported by Himebaugh et. al  in cases of TBU. They report that fluorescence and retro-illumination images are indicative of the increase in microscopic surface irregularities in regions of TBU. They also show SHWFS spots become very blurred in the corresponding TBU regions. By simulation they show that the SH spot blur cannot be explained by the second pass macro-aberrations (similar to our comparison of Fig. 9(a) and 9(b) to Fig. 9(c) and 9(d)). Therefore, they interpret this as evidence that TBU introduces micro-aberrations to the eye’s optics. Our study develops a metric to quantify the image blur caused by these micro-aberrations. We find that the majority of the double pass optical aberrations during TBU are indeed described by the unresolved micro-aberrations. From these two studies, we can deduce that the spatial patterns of micro-aberrations observed during TBU are likely to be correlated with the spatial distribution of exposed rough surface in the tear film.
To see the changes of the light scatter during TBU, we subtracted the scatter map for the baseline data from that for the TBU data (Fig. 11). This difference of the two average scatter maps is our preferred measure of the optical change for the following reason: Since the lower bound L and the upper bound U for the first pass were estimated based on the quality of the best spot in the raw image, both of the bounds include the retinal scatter. The average (L + U)/2 is subtracted from the radial variance for the double pass PSF blur. The radial variance for the second pass macro-aberrations is further subtracted leaving the radial variance for the second pass micro-aberrations. However, it is possible that this scatter map for the second pass micro-aberrations might contain the blur from other sources such as the crystalline lens or aqueous or vitreous humour, some of which are varying slowly and others of which are varying dynamically during TBU. The difference of the two average scatter maps records the optical changes that have occurred dynamically during the intervening time. For a similar reason, the difference in the lower bound maps or in the upper bound maps can be used to bound such changes. Furthermore, the average maps differ by only a small amount, (U−L)/2, from either of the boundary maps, which imply that indeed our bounds are sharp. Comparing the difference map in Fig. 11 with the raw image in Fig. 9(b) shows that the features in the original spot images are also present in the radial variance map. This is a promising optical method for understanding time dynamics of micro-aberrations during TBU. The current method for spot blur analysis can also be applied to other problems where fine scale aberrations are important, e.g. cataract.
6.2 Resolution of SH aberrometer
The ability of the analytical method described in this paper to quantify scatter produced by micro-aberrations is limited by two instrument-specific parameters. First, the cut-off order for the macro-aberrations depends on the resolution of the aberrometer (Fig. 1). Theoretically, as we reduce the size of the lenslets, the zonal wavefront may be able to quantify finer scale aberrations present in the eye by increasing the cut-off order and therefore expanding the range of square integrable functions.
One currently available high-resolution aberrometer (COAS HD, AMO) has square lenslets of size 108 microns and of a focal length 2.421 mm. Its camera captures the spot images formed by each lenslet with a 16x16 pixel array. The angular resolution of this system for each pixel is therefore 0.108mm/2.421mm x 1/16 x 180°/x 60 = 9.58 arcmin/pixel. Our analysis from Section 4 suggests that in order to quantify blurring of SH spots, an angular resolution of ≤ 1 armin/pixel is required. This is about 1/10 of the angular resolution of the high resolution aberrometer that is currently available. The aberrometer used in the current study was built in our lab and has a resolution of 0.9740 arcmin/pixel. That is, we have used a low resolution aberrometer but with a high sampling density of the SH image, whereas newer aberrometers might have higher pupil sampling density, but lower sampling density of the image.
The second important instrument limitation is introduced by the commonly employed “range limiting aperture” that prevents highly aberrated rays from being imaged by the lenslets. These limiting apertures, located at the common focal point of the relay lenses shown in Fig. 3, constrain the angular range of scattered light that will be included within the blurred images of the SHWFS spots. Thus largely aberrated rays and scatter with large angles are prevented from reaching the lenslets. The current method of spot blur analysis is likely to underestimate the scatter of the eye in that case.
6.3 Micro-aberrations, Strehl ratio and radial variance
Although the current algorithm quantifies the blur created by the micro-aberrations unresolved by the SH aberrometer, it does not quantify the characteristics of these micro-aberrations. Some approximate formulas have been suggested by Marechal  or Mahajan  by linking the Strehl ratio to the variance of the total phase aberrations. If the aberrations are weak (e.g. RMS of wavefront error < wavelength of light), the overall magnitude of wavefront RMS predicts (and can be derived from) the blurred Strehl ratio.
Alternatively, it is possible to link the Strehl ratio and the micro-aberrations. When we decompose the total phase aberrations into the deterministic macro-aberrations and the stochastic micro-aberrations, the RMS of the phase aberrations will also be statistically varying depending on each realization of the micro-aberrations. The average RMS of the phase aberrations can be expressed in terms of the RMS of the macro-aberrations and the average RMS of the micro-aberrations. Since the RMS of the macro-aberrations can be estimated, we can derive the RMS of the micro-aberrations from the measured Strehl ratio. We are currently investigating this approach.
Some linkage between the micro-aberrations and the radial variance can also be possible: If the phase aberrations are large, the phase difference between two points in the pupil can be approximated by the gradient of the phase. Since OTF is the autocorrelation of the pupil function, the average OTF can be approximated by the joint characteristic function of the two wavefront slopes. The PSF can be represented as the joint probability density function of the x- and y- slopes. Since the radial variance of the PSF is simply a sum of two second moments, by using the identity in Fourier analysis, it will be equal to the Laplacian of the OTF evaluated at the origin. In particular, if we assume that the x- and y- slopes of the wavefronts follow the Gaussian distributions with zero means and variances, and respectively, then the radial variance of the PSF will be proportional to the sum . This analysis implies that when the phase variance is large, the wavefront slopes undergo large fluctuation. It is the slopes of the wavefront that determine the distribution of light energy in the PSF. Applying this calculation to the radial variance of the PSF contributed by micro-aberrations, we can extract estimates for the variance of the wavefront slopes for the random phase screen, details of which will be reported in near future.
6.4 The whole eye PSF
Conventional techniques to measure the optical performance of the human eye have been dominated by the resolved macro-aberrations. The modal wavefront aberrations are often favored. The zonal wavefront aberrations have recently gained popularity among the clinical scientists studying tear film optics . However, our results show that across individual lenslet apertures the resolved macro-aberrations comprise only 13% of the total aberrations on the double pass (or 29% of the total blur created by the second pass) and that much of the blur of SH spots, therefore, was caused by the micro-aberrations. This implies that the conventional SH analysis is not sufficient to explain all of the optical complications of the eye and can be improved by including the micro-aberrations and their blurring of the retinal image in order to explain visual performance. Deducing the PSFs of the whole eye from blurring of SH spots produced by wave aberrations over individual lenslet apertures is challenging. It is unclear whether the PSF for the whole eye can be built from these local PSFs. In order to reveal the optical and visual impact of these micro-aberrations, a method for combining these spatially uncorrelated (statistical) micro-aberrations from each sub-aperture needs to be developed.
Because of the inherent limitations of SH aberrometer spatial resolution, it is not possible to quantify aberrations deterministically beyond the instrument resolution limit. However, we have developed a statistical description of the optical contributions of these micro-aberrations. Applying this approach to the human eye, we showed that high levels of forward light scatter are correlated with the regions of TBU, which we attribute to micro-aberrations.
Our attention has so far been focused on demonstrating the feasibility of scatter analysis algorithm for human eyes. An extended analysis for a large population of eyes or for dry eyes could be a useful future study. This method could also be useful to quantify the optical qualities of eyes with cataract or refractive surgery or to test the quality of intraocular lenses, where fine micro-aberrations are major concerns. Our data confirm that classic SHWFS analysis may fail to reveal visually significant optical perturbations, but can be extended successfully into the domain of micro-aberrations by the analysis of spot quality.
Appendix: Radial variance of bivariate random variables and point spread functions
The mapping between optics and probability theory enables us to explore the radial variances both for random variables and for point spread functions (PSF) at the same time. We extend below the optical analysis for the 1 dimensional case by Bracewell  to the 2 dimensional case without imposing additional assumptions on random variables.
Let X, Y, and Z be an object, an imaging system, and an image, respectively, in two dimensional space.
We denote the point spread function (PSF) for X, Y, Z by respectively. Using standard imaging optics, the PSF for the image is a convolution of the PSF for the object and the PSF for the system,
A natural choice for generalizing the second moment would be
Since measures the spread of in the radial direction, we call the radial variance. The definitions for the radial variances for and are the same by replacing with and respectively. We claim that those radial variances add under convolution,
We note that
and are the second order joint moments of the PSF The similar relation holds for pX and
Our proof utilizes optical transfer function (OTF). We let be the OTFs of and respectively. Since the OTFs and the PSFs are Fourier transform pairs, from Eq. (16), we can deduce
The OTFs can be expanded in terms of and moments of the corresponding PSFs. The second derivatives of the OTF can be found using the second moments of
On the other hand, by applying the product rule to Eq. (18), the second derivatives of are given in terms of the first and second derivatives of and
Evaluating them at
Similarly, the first moments can be expressed in terms of the first moments for and
Finally, the radial variance for
We thank Dr. Ray Applegate, Dr. Tom Van den Berg, Dr. Thomas Applegate, and Dr. Erez Ribak for helpful comments to improve our manuscript. Our research is supported by NEI/NIH grant R01EY05109 to LNT.
References and links
1. P. E. King-Smith, B. A. Fink, J. J. Nichols, K. K. Nichols, and R. M. Hill, “Interferometric imaging of the full thickness of the precorneal tear film,” J. Opt. Soc. Am. A 23(9), 2097–2104 (2006). [CrossRef]
2. S. D. Klyce and R. W. Beuerman, “Structure and Function of the Cornea,” in The Cornea, H. E. Kaufman, B. A. Barron, and M. B. McDonald, eds. (Butterworth-Heinemann, 1998), pp. 3–50.
4. P. E. King-Smith, J. J. Nichols, K. K. Nichols, B. A. Fink, and R. J. Braun, “Contributions of evaporation and other mechanisms to tear film thinning and break-up,” Optom. Vis. Sci. 85(8), 623–630 (2008). [CrossRef] [PubMed]
6. S. Koh, N. Maeda, Y. Hirohara, T. Mihashi, S. Ninomiya, K. Bessho, H. Watanabe, T. Fujikado, and Y. Tano, “Serial measurements of higher-order aberrations after blinking in normal subjects,” Invest. Ophthalmol. Vis. Sci. 47(8), 3318–3324 (2006). [CrossRef] [PubMed]
7. S. Koh, N. Maeda, T. Kuroda, Y. Hori, H. Watanabe, T. Fujikado, Y. Tano, Y. Hirohara, and T. Mihashi, “Effect of tear film break-up on higher-order aberrations measured with wavefront sensor,” Am. J. Ophthalmol. 134(1), 115–117 (2002). [CrossRef] [PubMed]
8. T. Mihashi, Y. Hirohara, S. Koh, S. Ninomiya, N. Maeda, and T. Fujikado, “Tear film break-up time evaluated by real-time Hartmann-Shack wavefront sensing,” Jpn. J. Ophthalmol. 50(2), 85–89 (2006). [CrossRef] [PubMed]
9. R. Montés-Micó, J. L. Alió, G. Muñoz, J. J. Pérez-Santonja, and W. N. Charman, “Postblink changes in total and corneal ocular aberrations,” Ophthalmology 111(4), 758–767 (2004). [CrossRef] [PubMed]
11. S. Koh, N. Maeda, Y. Hirohara, T. Mihashi, K. Bessho, Y. Hori, T. Inoue, H. Watanabe, T. Fujikado, and Y. Tano, “Serial measurements of higher-order aberrations after blinking in patients with dry eye,” Invest. Ophthalmol. Vis. Sci. 49(1), 133–138 (2008). [CrossRef] [PubMed]
12. S. Koh, N. Maeda, Y. Hori, T. Inoue, H. Watanabe, Y. Hirohara, T. Mihashi, T. Fujikado, and Y. Tano, “Effects of suppression of blinking on quality of vision in borderline cases of evaporative dry eye,” Cornea 27(3), 275–278 (2008). [CrossRef] [PubMed]
13. R. Montés-Micó, A. Cáliz, and J. L. Alió, “Wavefront analysis of higher order aberrations in dry eye patients,” J. Refract. Surg. 20(3), 243–247 (2004). [PubMed]
14. R. Montés-Micó, A. Cáliz, and J. L. Alió, “Changes in ocular aberrations after instillation of artificial tears in dry-eye patients,” J. Cataract Refract. Surg. 30(8), 1649–1652 (2004). [CrossRef] [PubMed]
17. T. Ferrer-Blasco, S. García-Lázaro, R. Montés-Micó, A. Cerviño, and J. M. González-Méijome, “Dynamic changes in the air-tear film interface modulation transfer function,” Graefes Arch. Clin. Exp. Ophthalmol. 248(1), 127–132 (2010). [CrossRef]
18. R. Tutt, A. Bradley, C. Begley, and L. N. Thibos, “Optical and visual impact of tear break-up in human eyes,” Invest. Ophthalmol. Vis. Sci. 41(13), 4117–4123 (2000). [PubMed]
19. N. L. Himebaugh, L. N. Thibos, A. Bradley, G. Wilson, and C. G. Begley, “Predicting optical effects of tear film break up on retinal image quality using the Shack-Hartmann aberrometer and computational optical modeling,” Adv. Exp. Med. Biol. 506(Pt B), 1141–1147 (2002). [CrossRef]
21. G. T. Timberlake, M. G. Doane, and J. H. Bertera, “Short-term, low-contrast visual acuity reduction associated with in vivo contact lens drying,” Optom. Vis. Sci. 69(10), 755–760 (1992). [CrossRef] [PubMed]
23. H. Liu, L. Thibos, C. G. Begley, and A. Bradley, “Measurement of the time course of optical quality and visual deterioration during tear break-up,” Invest. Ophthalmol. Vis. Sci. 51(6), 3318–3326 (2010). [CrossRef] [PubMed]
24. J. Liang, B. Grimm, S. Goelz, and J. F. Bille, “Objective measurement of wave aberrations of the human eye with the use of a Hartmann-Shack wave-front sensor,” J. Opt. Soc. Am. A 11(7), 1949–1957 (1994). [CrossRef]
25. N. H. Himebaugh, K. E. Thorn, L. N. Thibos, A. Bradley, and C. G. Begley, “Temporal changes in visual acuity, optical aberrations and light scatter secondary to tear film break-up in soft contact lens wearers,” Optom. Vis. Sci. 80(12s), 5 (2003).
26. J. E. Harvey and R. V. Shack, Light-Scattering Characteristics of Optical Surfaces, (University of Arizona, 1975).
27. W. J. Donnelly 3rd, K. Pesudovs, J. D. Marsack, E. J. Sarver, and R. A. Applegate, “Quantifying scatter in Shack-Hartmann images to evaluate nuclear cataract,” J. Refract. Surg. 20(5), S515–S522 (2004). [PubMed]
29. J. Nam and J. Rubinstein, “Weighted Zernike expansion with applications to the optical aberration of the human eye,” J. Opt. Soc. Am. A 22(9), 1709–1716 (2005). [CrossRef]
30. D. R. Neal, C. D. Baer, and D. M. Topa, “Errors in Zernike transformations and non-modal reconstruction methods,” J. Refract. Surg. 21(5), S558–S562 (2005). [PubMed]
31. J. W. Goodman, Statistical Optics, Wiley classics library ed., Wiley classics library (Wiley, 2000).
32. M. Born and E. Wolf, Principles of Optics: Electromagnetic Theory of Propagation, Interference and Diffraction of Light, 7th (expanded) ed. (Cambridge University Press, 1999).
33. H. Zhao, D. T. Miller, L. N. Thibos, X. Hong, A. Bradley, X. Cheng, and N. Himebaugh, “A Fried's parameter for the human eye?,” Abstract for the annual Optical Society of America meeting (2000).
35. R. N. Bracewell, The Fourier Transform and its Applications, 3rd ed., MCGraw-Hill series in electrical and computer engineering (McGraw Hill, 2000).
36. J. Nam, J. Rubinstein, and L. N. Thibos, “Wavelength adjustment using an eye model from aberrometry data,” J. Opt. Soc. Am. A 27(7), 1561–1574 (2010). [CrossRef]
37. N. L. Himebaugh, C. G. Begley, L. N. Thibos, and A. Bradley, “The spatial orrelation between tear break-up and optical aberrations,” in preparation (2011).
38. W. Gao, R. S. Jonnal, B. Cense, O. P. Kocaoglu, Q. Wang, and D. T. Miller, “Measuring directionality of the retinal reflection with a Shack-Hartmann wavefront sensor,” Opt. Express 17(25), 23085–23097 (2009). [CrossRef]
39. A. Marechal, “Etude des effets combines de la diffraction et des aberrations geometriques sur l'image d'un point lumineux,” Rev. d'Opt 26, 257–277 (1947).
40. V. N. Mahajan, “Strehl ratio for primary aberrations in terms of their aberration variance,” J. Opt. Soc. Am. 73(6), 860–861 (1983). [CrossRef]