Open Access
Add to BibSonomy
Abstract
This paper presents a method for reconstructing the corneal surface. The proposed method was tested in 56 healthy and 15 post-orthokeratology corneas. The Medmont E300 Corneal Topographer was used to measure the anterior corneal elevation, and custom MATLAB scripts were employed for data analysis, fitting, and other computational processes. The results obtained were compared with the fitting to an ellipsoid and to a biconic, using an alternative method, showing similarities among the different approaches. Additionally, the advantages of this method and the biconic’s generality over the ellipsoid were also demonstrated. In conclusion, the method proposed offers an approach with potential applications in the field of visual and ophthalmic optics related with modeling of the cornea and other optical surfaces.
© 2024 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement
1. INTRODUCTION
The cornea, the clear, outermost layer of the eye, plays a crucial role in the refractive power of the human eye. In fact, the anterior surface of the cornea has the highest refractive power, and due to its importance in ocular refraction, accurate measurement and modeling of its surface are critical in several applications related to visual and ophthalmic optics [1,2].
One of the most common ways to describe the topographic surface of the cornea, $S({x,y})$, is by considering the sum of two contributions:
where $b({x,y})$ is a base function that should appropriately describe the fundamental optical properties of the corneal surface, and $r({x,y})$ represents the residuals after subtracting the base function from the topographic elevation data. These residuals are commonly described using polynomial expansions such as Zernike polynomials, among other methods. A proper description of the base surface in terms of its conic coefficients (apical radius and conic constant) is very useful in visual and ophthalmic optics to infer the cornea’s imaging forming properties, such as its paraxial properties, cylindrical power and axis, spherical aberration, etc. Therefore, it is of great importance to find the base surface that best describes the fundamental optics of the cornea.Different approaches have been used to model the geometrical parameters of the cornea, with different levels of complexity [3]. The most simplistic method is the sphere, described only by its radius of curvature $R$. This simple model is inaccurate to describe the intrinsically aspheric and usually non-rotationally symmetric profile of the cornea, affecting the proper description of corneal aberrations. A conic surface, described by its apical radius and conic constant ($R$ and $Q$), can achieve a more realistic representation of spherical aberration, but still limits the proper description of astigmatic corneas, due to the rotational symmetry of conics. The biconic surface can handle the lack of rotational symmetry by describing the cornea with two principal meridians, with independent radius of curvature and conic constants (${R_x}$, ${Q_x}$, ${R_y}$, ${Q_y}$) [4], assuming that the sagittal depth at different meridians varies with the squared sine of the meridional angle, allowing for the description of corneal astigmatism as well as spherical aberration. Increasing the number of parameters of the model allows for a more precise mathematical representation of the corneal surface and, thus, a more realistic description of its fundamental optical properties. However, the drawback is the increased complexity of its application when considering realistic scenarios, where due to fixation misalignments the corneal and topographical axes do not match. An interesting application can be seen in the conical null-screen proposed recently by Campos-García et al., in which the image is formed by reflection on a biconic surface. [5].
One of the most well-known methods for fitting the corneal surface that tries to merge efficiency and simplicity is the non-revolution second order ellipsoid proposed by Navarro et al. [6]. It consists of a mathematical formalism to obtain the geometrical parameters that allow to describe the optical properties of the corneal surface (apical radius and conic constant) from its principal meridians, as well as a mathematical representation to transform the coordinate system of the topographer to any other system, including the useful canonical representation. However, due to its inherent definition ${R_x} = {a^2}/c$, ${R_y} = {b^2}/c$, ${Q_x} = ({a^2}/{c^2}) - 1$ and ${Q_y} = ({b^2}/{c^2}) - 1$, the ellipsoid is not as general as the biconic due to the correlation among the different parameters, and might limit the representation of less common corneas.
Following a similar approach, this paper presents a new method for fitting the corneal surface that makes use of the orthogonality properties of Zernike polynomials and the linear formulation for a biconic surface. The novel method allows for a precise and numerical efficient reconstruction of the corneal surface to the more general biconic surface and the calculation of its geometrical parameters.
2. METHODS
Anterior corneal elevation maps were acquired with the Placido based system Medmont E300 (Medmont, Australia) from two groups of patients: one consisting of 56 healthy corneas and the other consisting of 15 post-orthokeratology corneas, used to validate the method on topographically altered corneal shapes. The Medmont software allows the output of ASCII files containing 32 rings of corneal height data in 300 equally spaced half-meridians. The “.hgt” files used for the analysis represent the elevation relative to a plane tangent to the vertex normal, while the “.dst” files represent the radial distance from the vertex normal. Custom made software developed in MATLAB (The MathWorks Inc., USA) was employed for data analysis, fitting, and other computational processes. The Levenberg–Marquardt algorithm, implemented in the MATLAB 2023a curve fitting toolbox, was used to perform the nonlinear fittings by means of the “Bisquare” method of weights.
A biconic surface exhibits non-rotational symmetry and is characterized by two principal meridians, Fig. 1, described by two radii of curvature (${R_x}$, ${R_y}$), and two conic constants (${Q_x}$, ${Q_y}$), where $Q = 0$ refers to a sphere. The sagitta of a biconic surface centered and aligned according to its principal meridians in Cartesian coordinates is given by the following expression:
Fig. 1. Biconic surface defined by its radii of curvature $({{R_x},{R_y}})$ and conic constants $({{Q_x},{Q_y}})$ along the two principal meridians.
There are two primary challenges associated with directly applying this fitting approach to corneal elevation data. The first is related to possible misalignments between the optical and the keratometric axes, which introduce tilts and decentrations to the surface and might affect the accurate estimation of its geometric parameters [6]. The second is related to the nonlinear fitting approach that needs to be applied. Nonlinear fitting methods can be highly biased by the initial parameters estimation. A proper selection of the initial values can be crucial to avoid the algorithm to converge to a local minimum, which would result in a different solution from the intended one. Linear algorithms, on the other hand, do not depend on an initial departing value, making it numerically more robust and always produce a unique closed form least squares solution.
The novel method proposed in this work makes use of Zernike polynomials and their orthogonal properties to eliminate the effect of misalignments $({{\theta _x},{\theta _y},{\theta _z}})$ and other surface irregularities. It should be noted that by misalignments we refer to the discrepancies between the orientation of the keratometric axis and the optical axis defined by $({{\theta _x},{\theta _y},{\theta _z}})$. The main reason for such misalignments lies in the fact that topographers typically do not use the visual axis as the main axis of reference; it can also be affected by eye movements or calibration errors in the instrument. This implies that, depending on the degree of misalignment, there is a potential risk of inaccurately estimating the geometric parameters of the cornea. Furthermore, during optical simulations, where the pupil axis is employed, this tilt could impact the assessment of aberrations [7].
First, the elevation map (Fig. 2) is adjusted to a surface defined by Zernike polynomials using Eq. (3):
where $z({x,y})$ represents the sagitta of the surface, and $C_n^m$ is the coefficient associated with the Zernike polynomial $Z_n^m$ of radial frequency $n$ and azimuthal order $m$, defined in Cartesian coordinates $({x,y})$. Due to the orthogonality properties of Zernike polynomials, the biconic properties of the surface will be represented by the spherical and astigmatic Zernike modes, truncated to order 4. The choice of order 4 relates to the fact that when expanding the conic section equation into a Maclaurin series, the fourth order coefficient is the first to introduce the effect of the conic constant: $z\sim{r^2}/({2R}) + ({p/({8{R^3}})}){r^4} \,+\def\LDeqbreak{} ({{p^2}/({16 + {R^5}})}){r^6} + \ldots$, with $p = Q + 1$. Hence, the Zernike surface of order 4 already contains information related to the conical shape, even if some surface deformations are smoothed.The rotation angle around the $z$ axis (${\theta _z}$, astigmatism axis) was obtained from the orientation of the astigmatic modes and subsequently applied to the coordinate system. The expression to obtain the angle at which the surface is rotated is as follows (see Appendix A for the full derivation):
Once the angle has been calculated, it is applied to the inverse rotation to the surface, ${r_z}({- {\theta _z}})$ (Fig. 3):
Fig. 3. Rotation of the surface along the $z$ axis, with the $x$ and $y$ tilts and irregularities removed.
After applying the rotation matrix, the corneal elevation data was fitted to the linear equation of a biconic surface:
The radii of curvature (${R_x}$, ${R_y}$), the conic constants (${Q_x}$, ${Q_y}$), and the rotation angle with respect to the $z$ axis, obtained with the new method were compared to the values obtained with the ellipsoid and biconic methods proposed by Navarro et al. [6,8]. In the method of Navarro, the four parameters derived from the linear ellipsoid for ${R_x}$, ${R_y}$, ${Q_x}$, and ${Q_y}$ are optimized (the biconic assumes the same $z$ axis as the ellipsoid) to account for the increased generality of the biconic shape versus the special case of the ellipsoid. The process is summarized in the diagram of Fig. 4.
Fig. 4. Diagram of the proposed methods for obtaining the radii of curvature and conic constants [6,8].
For the statistical analysis, Shapiro–Wilk and Kruskal–Wallis tests were used to study the normality of the sample and determine whether there were statistically significant differences between the methods.
The goodness of the fit was calculated as the root mean square error (RMSe). All the analyses were conducted for an 8 mm corneal diameter. Finally, illustrative examples are provided to demonstrate situations where the versatility of the biconic surface outweighs that of the ellipsoid.
Table 1. Residual Zernike Coefficients for the Generation of Synthetic Corneas, Showing the Mean of the Coefficients (${\boldsymbol C}$) and Standard Deviation (${{\boldsymbol \sigma}_{\boldsymbol C}}$)
Table 2. Corneal Fitting from 56 Healthy Eyes, Showing the Mean of the Radii of Curvature (${{\boldsymbol R}_{\boldsymbol x}}$, ${{\boldsymbol R}_{\boldsymbol y}}$), Axis Angle, Conic Constants (${{\boldsymbol Q}_{\boldsymbol x}}$, ${{\boldsymbol Q}_{\boldsymbol y}}$), Root Mean Square Error (RMSe), and Standard Deviation (${\boldsymbol \sigma}$)a
Table 3. Corneal Fitting from 15 Post-orthokeratology Treated Eyes, Showing the Mean of the Radii of Curvature (${{\boldsymbol R}_{\boldsymbol x}}$, ${{\boldsymbol R}_{\boldsymbol y}}$), Axis Angle, Conic Constants (${{\boldsymbol Q}_{\boldsymbol x}}$, ${{\boldsymbol Q}_{\boldsymbol y}}$), Root Mean Square Error (RMSe), and Standard Deviation (${\boldsymbol \sigma}$)a
Synthetic corneas with known parameters were generated using Eq. (7), in a ${256} \times {256}$ grid with an 8 mm diameter, to simulate extreme scenarios, and used to demonstrate the greater generality of the biconic surface versus the ellipsoid method:
The biconic surface was used as a basis, introducing the geometric parameters $({{R_x},{R_y},{Q_x},{Q_y}})$ manually and rotating it at an angle ${\theta _z}$. Subsequently, a residual surface of Zernike polynomials of order 4 was added to enhance the realism of the corneas. These residuals, derived from a multivariate Gaussian model based on the mean and covariance matrix of residuals from our 56 healthy corneas, include asymmetric and non-astigmatic terms. Table 1 shows the mean values and the standard deviation of the residuals used.
3. RESULTS
Two groups were evaluated, one consisting of 56 healthy corneas and the other consisting of 15 post-orthokeratology corneas. The resulting average parameters for the 56 normal corneas are shown in Table 2 for the different fitting methods, following the procedures explained above.
It is observed that the geometric parameters obtained are similar in the different methods, as well as the goodness of fit, with the only difference being the lower dispersion in the RMSe of the new method.
In addition to the fits for normal corneas, the same methods were applied to corneas with a less common shapes to demonstrate the robustness of the method [9]. Table 3 shows the fitting results for 15 post-orthokeratology corneas.
For the case of post-orthokeratology corneas, the radii, the axes, and the RMSe obtained are similar, but there were differences in the conic constants obtained with the different methods. Statistical analysis presented in Tables 2 and 3 shows no significant differences between methods.
Finally, multiple extreme scenarios were simulated based on Eq. (3) rotating around the $z$ axis, to investigate the superior generality of the biconic over the ellipsoid. Results are presented in Tables 4–7.
Table 4. Fitting from a Synthetic Regular Corneal Surface, Showing the Radii of Curvature (${{\boldsymbol R}_{\boldsymbol x}}$, ${{\boldsymbol R}_{\boldsymbol y}}$), Axis Angle, Conic Constants (${{\boldsymbol Q}_{\boldsymbol x}}$, ${{\boldsymbol Q}_{\boldsymbol y}}$), and Root Mean Square Error (RMSe)
Table 5. Fitting of a Synthetic Corneal Surface with Highly Different Conic Constants on Each Meridian, Showing the Curvature Radii (${{\boldsymbol R}_{\boldsymbol x}}$, ${{\boldsymbol R}_{\boldsymbol y}}$), Axis Angle, Conic Constants (${{\boldsymbol Q}_{\boldsymbol x}}$, ${{\boldsymbol Q}_{\boldsymbol y}}$), and Root Mean Square Error (RMSe)
Table 6. Fitting of a Synthetic Corneal Surface with Highly Different Radii of Curvature on Each Meridian, Showing the Curvature Radii (${{\boldsymbol R}_{\boldsymbol x}}$, ${{\boldsymbol R}_{\boldsymbol y}}$), Axis Angle, Conic Constants (${{\boldsymbol Q}_{\boldsymbol x}}$, ${{\boldsymbol Q}_{\boldsymbol y}}$), and Root Mean Square Error (RMSe)
Table 7. Fitting from a Synthetic Corneal Surface with Different Conic Constant Signs in Each Meridian, Showing the Curvature Radii (${{\boldsymbol R}_{\boldsymbol x}}$, ${{\boldsymbol R}_{\boldsymbol y}}$), Axis Angle, Conic Constants (${{\boldsymbol Q}_{\boldsymbol x}}$, ${{\boldsymbol Q}_{\boldsymbol y}}$), and Root Mean Square Error (RMSe)
In Table 4, it can be observed that both models reliably recover the geometric parameters of the synthetic cornea.
The next scenario evaluates the fitting of a cornea with significantly different conic constants for each principal meridian. Results are presented in Table 5.
It can be seen that in the case of highly disparate conic constants, the ellipsoid model does not converge to the best-fitting solution. This leads to moderately disparate radii of curvature but completely different conic constants, resulting in an increased RMSe. In contrast, the linear fitting method continues to provide reliable results.
Another possible scenario is the case of a high astigmatic cornea, with highly disparate apical radii of curvature corresponding to about 6D of astigmatism. Results are presented in Table 6.
In this scenario, the ellipsoid reveals to be more robust than in the previous case, but still fails to recover the original conic constants. Once again, the linear fitting of a biconic recovers all the parameters perfectly.
Last, the case of a cornea with similar radii of curvature where one meridian has an oblate shape while the other is prolate is considered. Results are presented in Table 7.
The case of a synthetic cornea with conic constants of different signs does not meet the conditions of the definition of an ellipsoid, which leads to a suboptimal solution that results in a large discrepancy in the paraxial astigmatism predicted by the ellipsoid method.
Table 8. Fit Quality (RMSe) of Healthy (HC) and Post-orthokeratology (OK) Corneas Obtained for the Linear Biconic Using Several Zernike Surfaces of Different Orders
Note that the linear biconic surfaces were fitted with a Zernike polynomial surface of order 4. Table 8 shows the RMSe obtained for the linear biconic fit according to different orders.
As expected, the RMSe obtained for orders higher than four is not significantly different from the RMSe obtained with the fourth order Zernike surface, but it seems to increase for higher orders probably due to the noisy nature of the clinical data.
4. DISCUSSION
To the best of our knowledge, the use of a biconic surface to model the cornea was first proposed by Schwiergerling et al. [4]. These authors introduced an iterative nonlinear least squares fitting method to determine the radii of curvature and conic constants of the principal meridians, from corneal elevation data. While this approach provided a viable solution for describing the corneal surface, it did not consider the misalignments between the keratometric and optical axes of the cornea. Additionally, this proposal suffered from the limitations associated with nonlinear fitting, as mentioned earlier. In 2006, Navarro et al. proposed a potential solution to address these drawbacks [6]. These authors suggested the use of the general expression of a non-revolution ellipsoid to retrieve the intrinsic coordinates of the cornea and recover its geometrical parameters. However, the ellipsoid is not as general as the biconic due to the correlation among the different parameters and might limit the representation of less common corneas. Janunts et al. suggested in 2014 a parametric fitting approach to obtain the radii of curvature and conic constants along the principal meridians [10]. However, this method did not account for surface misalignments. Later in 2019, Navarro et al. proposed the use of the ellipsoid geometrical parameters as initial values for the nonlinear fit of a biconic surface [8]. Although this approach offered improvements over the initial ellipsoid fit, it still requires prior fitting to the ellipsoid and relies on a nonlinear fitting method, which could lead to ambiguous results when dealing with corneas with significantly different radii of curvature and conic constants.
In this work, a new method is presented to address all the aforementioned issues. The method makes use of the orthogonality of Zernike polynomials to correct for misalignments and obtain the axis of the astigmatism. Finally, a linear expression is derived for the biconic surface that allows for a quick and direct determination of the corneal geometric parameters through a closed form linear least squares solution.
From the results described in Tables 2 and 3, both regular and post-orthokeratology corneas are well described with both the biconic linear least squares method and the ellipsoid fitting. The geometric parameters are obtained with similar robustness, with only a small increase in performance for the biconic for the case of regular corneas, expressed by the lower dispersion of the RMSe. The radii of curvature (${R_x}$, ${R_y}$), astigmatism axis, conic constants (${Q_x}$, ${Q_y}$), and RMSe were similar for all methods, demonstrating that the new method proposed is able to account for misalignments of the cornea relative to the keratometric axis, as well as its intrinsic misalignments.
Another key advantage of the new method is the use of a linear least squares fitting method versus a nonlinear method, without the need for initial parameters. When the initial parameters are far from the expected values, the fitting process can converge to a local minimum, resulting in a suboptimal solution. As a result, an alternative option is to pre-fit the elevation data with an ellipsoid and use those parameters as initial estimators for the nonlinear fitting of the biconic surface, in an already aligned coordinate system [8]. However, this approach is less straightforward compared to the proposed one, and in the case of corneas with extreme differences in apical radii of curvature of conic constant between principal meridians, the nonlinear fit may not converge to an optimal solution, as demonstrated by the results summarized in Tables 4–7.
To demonstrate the greater generality of the biconic compared to the ellipsoid, different scenarios were generated using synthetic corneas. In the first scenario, presented in Table 4, both the biconic and the ellipsoid accurately recovered the geometric parameters of the synthetic surface. However, for more extreme scenarios, as the ones presented in Tables 5–7, the ellipsoid fails to recover the original parameters. The reason for these discrepancies lies in the inherent definition of the ellipsoid itself:
where $({a,b,c})$ represent the three semi-axes of the ellipsoid, and the geometric parameters are given by ${R_x} = {a^2}/c$, ${R_y} = {b^2}/c$, ${Q_x} = ({{a^2}/{c^2}}) - 1$, and ${Q_y} = ({{b^2}/{c^2}}) - 1$. It can be seen that these parameters are correlated with each other, unlike in the case of the biconic where all four parameters are independent. Consequently, the ellipsoid is a particular case of a biconic, and cannot accurately represent a biconic surface for the cases when the conic parameters differ significantly between principal meridians. This situation might be present in the case of post-orthokeratology corneas shown in Table 3, where it can be observed that the ellipsoid recovers an average conic constant ${Q_x} = 0.59$, whereas the biconic yields ${Q_x} = 0.45$. It is worth noting that after fitting the elevation data with the biconic surface, starting from the ellipsoid parameters, the conic constant decreases to ${Q_x} = 0.49$, in closer agreement to the linear least squares solution. Another reason may be the changes produced in the anterior surface of the cornea by the orthokeratology treatment means that it can no longer be described correctly with a second-order surface such as an ellipsoid, and requires higher order surfaces, like the biconic, a fourth order surface.To summarize, in both cases of normal and post-orthokeratology corneas, the results for the radii of curvature, conic constants, and astigmatism axes were similar among the different methods, with similar values for the goodness of the fit. An interesting result is the case shown in Table 2, where the mean RMSe and standard deviation for the method presented in this paper are slightly lower than the other two, suggesting that the new method might be numerically more robust. Furthermore, various scenarios have been presented to highlight the advantages of employing a biconic surface in the fitting process versus an ellipsoid. The limitations of the direct fitting of a biconic compared to the ellipsoid, regarding possible misalignments with the keratometric axis, were overcome with the use of Zernike polynomials. In conclusion, the new method proposed offers a precise and robust approach, with potential applications in the field of visual and ophthalmic optics related with modeling of the cornea and other optical surfaces.
APPENDIX A
To enhance the understanding of this study, we present a detailed derivation of the equations used for biconic surfaces. This step-by-step explanation will provide clarity and insight into the mathematical framework employed in our analysis.
Sagitta calculation
We begin by considering the equation of a sphere centered on its sagitta along the $z$ axis at a distance $R$:
If this expression is expanded,
This would be the case for a sphere, where $r^2=x^2+y^2$. However, if we want to generalize it for any conical surface, we need to consider the shape factor, denoted as $p$:
arriving to a similar expression like the one proposed by Baker [11]. Based on this equation, the general sagitta equation for a conic surface will be derived. Then the previous expression is multiplied by ${r^2}$, maintaining the identity, but with those terms dependent on the radius of curvature and the shape factor containing an ${r^2}$:Thus, by means of factoring squares we approach the final expression:
Now we have the general expression to define the sagitta of a conic surface in terms of curvature, $c = 1/R$, and its conic constant, $p = ({1 + Q})$. Therefore, if the coefficients are allowed to vary freely in each of the meridians, the sagitta for a biconic surface is obtained:
There is also another approach in polar coordinates $({r,\theta})$ for the sagitta of a biconic surface:
Then,
$C(\theta)$ and $Q(\theta)$ being the Euler curvature and the asphericity along the meridian at angle $\theta$ [4]:
Next, in order to explore the procedure for obtaining the biconic linear expression, it is simply needed to refer to Eq. (A4):
Considering that both the conic shape factor, $p$, and the curvature radius, $R$, vary differently for each meridian, the linear expression is obtained:
Finally, the following expressions for the sagitta can be encountered depending on the level of generality, in Cartesian coordinates [12]:
Astigmatism axis (Zernike)
In this example we will consider that we fit a Zernike surface up to order 6. Thus, the function containing the astigmatic terms of azimuthal frequency $-2$ and $+2$, $f({r,\theta})$, will be given by the following expression:
To obtain the angle, we must calculate when the derivative of this function, $f({r,\theta})$, equals zero:
Consequently,
Surface tilts (Zernike)
Surface tilts due to poor fixation (${\theta _x},{\theta _y}$) can be defined as rotations (${r_x},{r_y},{r_z}$) in a coordinate system $({x,y,z})$. Therefore, the contributions of tilt in $x$ and $y$, $C_{- 1}^1$ and $C_1^1$, respectively, along with their contribution in the remaining comatic modes, represent surface rotations (${\theta _x},{\theta _y}$). This can be demonstrated and calculated from the Zernike polynomials and their gradients [13], using the following equations, where ${r_{\rm{max}}}$ is the normalization radius:
Funding
H2020 Marie Skłodowska-Curie Actions (956720).
Acknowledgment
This project has received funding from the European Union’s Horizon 2020 Research and Innovation Programme under the Marie Skłodowska-Curie.
Disclosures
The authors have no conflicts of interest to declare.
Data availability
Data and code underlying the results presented in this paper are not available at this time but may be obtained from the authors upon request.
REFERENCES
1. D. Cano, S. Barbero, and S. Marcos, “Comparison of real and computer-simulated outcomes of LASIK refractive surgery,” J. Opt. Soc. Am. A 21, 926–936 (2004). [CrossRef]
2. D. Gatinel, J. Malet, T. Hoang-Xuan, et al., “Corneal elevation topography: best fit sphere, elevation distance, asphericity, toricity, and clinical implications,” Cornea 30, 508–515 (2011). [CrossRef]
3. D. Atchison and G. Smith, Optics of the Human Eye (Butterworth-Heinemann, 2000).
4. J. Schwiegerling and R. W. Snyder, “Custom photorefractive keratectomy ablations for the correction of spherical and cylindrical refractive error and higher-order aberration,” J. Opt. Soc. Am. A 15, 2572–2579 (1998). [CrossRef]
5. M. Campos-García, O. Huerta-Carranza, L. Á. Pantoja-Arredondo, et al., “Conical null-screen design for evaluating a biconical surface using a smartphone-based corneal topographer,” Proc. SPIE 12221, 122210L (2022). [CrossRef]
6. R. Navarro, L. González, and J. L. Hernández, “Optics of the average normal cornea from general and canonical representations of its surface topography,” J. Opt. Soc. Am. A 23, 219–232 (2006). [CrossRef]
7. T. O. Salmon and L. N. Thibos, “Videokeratoscope-line-of-sight misalignment and its effect on measurements of corneal and internal ocular aberrations [published correction appears in J Opt Soc Am A Opt Image Sci Vis. 2003 Jan;20(1):195.],” J. Opt. Soc. Am. A 19, 657–669 (2002). [CrossRef]
8. R. Navarro, J. J. Rozema, M. H. Emamian, et al., “Average biometry of the cornea in a large population of Iranian school children,” J. Opt. Soc. Am. A 36, B85–B92 (2019). [CrossRef]
9. M. Faria-Ribeiro, R. N. Belsue, N. López-Gil, et al., “Morphology, topography, and optics of the orthokeratology cornea,” J. Biomed. Opt. 21, 075011 (2016). [CrossRef]
10. E. Janunts, M. Kannengießer, and A. Langenbucher, “Parametric fitting of corneal height data to a biconic surface,” Z. Med. Phys. 25, 25–35 (2015). [CrossRef]
11. T. Y. Baker, “Ray tracing through non-spherical surfaces,” Proc. Phys. Soc. 55, 361–364 (1943). [CrossRef]
12. T. Houllier, “Optical imaging systems with freeform surfaces: optimization algorithms study and freeform surfaces metrology,” (Université de Lyon, 2021).
13. R. Navarro, “Refractive error sensing from wavefront slopes,” J. Vis. 10(13): 3 (2010). [CrossRef]
