Determination of the normalized surface height autocorrelation function of a two-dimensional randomly rough dielectric surface by the inversion of light scattering data

An expression is obtained on the basis of phase perturbation theory for the contribution to the mean differential reflection coefficient from the in-plane co-polarized component of the light scattered diffusely from a two-dimensional randomly rough dielectric surface when the latter is illuminated by s-polarized light. This result forms the basis for an approach to inverting experimental light scattering data to obtain the normalized surface height autocorrelation function of the surface. Several parametrized forms of this correlation function, and the minimization of a cost function with respect to the parameters defining these representations, are used in the inversion scheme. This approach also yields the rms height of the surface roughness, and the dielectric constant of the dielectric substrate if it is not known in advance. The input data used in validating this inversion consists of computer simulation results for surfaces defined by exponential and Gaussian surface height correlation functions, without and with the addition of multiplicative noise, for a single or multiple angles of incidence. The reconstructions obtained by this approach are quite accurate for weakly rough surfaces, and the proposed inversion scheme is computationally efficient.


I. INTRODUCTION
Statistical information about the roughness of a surface is contained in its rms height and in its normalized surface height autocorrelation function. Efforts to obtain these properties of rough surfaces from measurements of light scattered into the far field from them are of interest because of the contactless nature of this approach, and because measurements in the far field are easier to carry out than measurements in the near field. This problem has been studied in the past by several authors. In the case of a one-dimensional randomly rough dielectric surface, Chakrabarti et al. [1] inverted by a Fourier transformation an expression for the contribution to the mean differential reflection coefficient, obtained in the Kirchhoff approximation, from light scattered diffusely from the surface. The incident light was s-polarized and the plane of incidence was perpendicular to the generators of the surface. Good agreement with numerically generated scattering data was obtained for weakly rough surfaces.
The case of a two-dimensional randomly rough surface was studied by Chandley [2], and by Marx and Vorburger [3]. Chandley used scalar diffraction theory and a thin random phase screen approximation [4] to model the interaction of light with the randomly rough surface. A thin phase screen may be regarded as a layer of negligible thickness that alters the phase of the wave scattered from it but does not change its magnitude. It is derived from simple optical path length and geometrical optics arguments. He used the angular dependence of the mean intensity of the scattered light in the far field as the experimental quantity to be inverted. The nature of his scattering model allowed this inversion to be carried out by means of a Fourier transformation. The dielectric constant of the scattering medium does not appear explicitly in Chandley's theory which means that it is impossible to use it to recover the dielectric constant of the scattering medium from the experimental scattering data if it is not known in advance.
In their study of this problem, Marx and Vorburger applied the Kirchhoff approximation for the scattering of a scalar plane wave from a two-dimensional randomly rough perfectly conducting surface to obtain the mean intensity of the scattered field. The determination of the rms height of the surface and the normalized surface height autocorrelation function was achieved by assuming an expression for the latter function of a particular analytic form and by the determination of the parameters defining it by a least squares fit of the theoretical mean intensity to the experimental result.
In contrast to these studies, in this paper we present an approach to the determination of the rms height and the normalized surface height autocorrelation function of a two-dimensional randomly rough penetrable surface, in particular a dielectric surface, from the inversion of optical scattering data. It is based on a vector theory of rough surface scattering rather than on a scalar theory, namely phase perturbation theory [5]. The dielectric constant of the medium is taken into account in this approach. This version of rough surface scattering theory was chosen in this study because in a recent comparison between experimental data and the predictions of three perturbation theories for the scattering of electromagnetic radiation from two-dimensional randomly rough metal surfaces, it produced the best results [5]. We expect it to be equally accurate in describing the scattering of visible light from a two-dimensional randomly rough dielectric surface. Specifically, we use the expression for the contribution to the mean differential reflection coefficient from the in-plane, co-polarized component of the light scattered incoherently when the dielectric surface is illuminated at normal or non-normal incidence by s-polarized light.
This expression is evaluated with the use of an expression for the normalized surface height autocorrelation function that contains adjustable parameters. The values of these parameters are then determined by a least squares fit of the resulting expression to the corresponding experimental scattering data. The reconstruction of the parameters is performed using scattering data for both a single and several angles of incidence, and the sensitivity to noise of the reconstructed parameters is investigated. We note that the contribution to the mean differential reflection coefficient from the in-plane co-polarized component of the light scattered incoherently when the surface is illuminated by ppolarized light can also be used for this purpose. However, the expression one works with to effect this inversion is somewhat simpler in s polarization than in p polarization. In addition, there is no Brewster effect in s polarization, so that a smoother function of the scattering angle is being inverted in s polarization than in p polarization. It is for these reasons that we have chosen to work with s polarization.
x 1 Schematics of the scattering geometry considered in this work.

II. THE PHYSICAL SYSTEM STUDIED
The physical system we study in this paper consists of vacuum in the region x 3 > ζ(x ), and a dielectric medium, characterized by a dielectric constant ε that is real, positive and frequency independent, in the region Fig. 1]. Here x = (x 1 , x 2 , 0) is a position vector in the plane x 3 = 0. The surface profile function ζ(x ) is assumed to be a single-valued function of x that is differentiable with respect to x 1 and x 2 . It is also assumed to constitute a stationary, zero-mean, isotropic, Gaussian random process defined by where the angle brackets denote an average over the ensemble of realizations of ζ(x ), δ is the rms height of the surface, and W (|x |) is the normalized surface height autocorrelation function, with the property that W (0) = 1. The surface profile function has a Fourier integral representation, where Q = (Q 1 , Q 2 , 0) is a two-dimensional wave vector so that We also introduce the notationζ The Fourier coefficientζ(Q ) is also a zero-mean Gaussian random process defined by where g(|Q |), the power spectrum of the surface roughness, is defined by It follows from Eqs. (1) and (5) that g(|Q |) is normalized to unity,

III. SCATTERING THEORY
The surface x 3 = ζ(x ) is illuminated from the vacuum by an electromagnetic field of frequency ω. The electric field in the vacuum above the surface is the sum of an incident and a scattered field, E( The subscripts p and s denote the p-polarized (TM) and s-polarized (TE) components of each of these fields, respectively. The function α 0 (q ) in Eqs. (7) is defined as Maxwell's equations imply linear relations between A α (q ) and B β (q ), which we write in the form (α = p, s, β = p, s) The scattering amplitudes {R αβ (q |k )} play a significant role in the present theory because the mean differential reflection coefficient is given in terms of them. The differential reflection coefficient ∂R αβ (q |k )/∂Ω s is defined such that ∂R αβ (q |k )/∂Ω s dΩ s is the fraction of the total time -averaged flux in an incident field of β polarization the projection of whose wave vector on the mean scattering plane is k , that is scattered into a field of α polarization, the projection of whose wave vector on the mean scattering plane is q , within an element of solid angle dΩ s about the scattering direction defined by the polar and azimuthal angles (θ s , φ s ). It is given by [6,7] ∂R αβ (q |k ) with [see Fig. 1] where (θ 0 , φ 0 ) and (θ s , φ s ) are the polar and azimuthal angles of incidence and scattering, respectively. S is the area of the plane x 3 = 0 covered by the rough surface. As we are dealing with scattering from a randomly rough surface, it is the average of this function over the ensemble of realizations of the surface profile function that we have to calculate. The contribution to this average from the light scattered incoherently is Closely related to the matrix of scattering amplitudes R(q |k ) is the scattering matrix S(q |k ) whose elements {S αβ (q |k )} are given by These elements satisfy the reciprocity relations [8] S pp (q |k which serve as a check on the correctness of their derivation. In terms of the elements of the scattering matrix, Eq. (10) takes the form This is the definition we will work with. In the Appendix it is shown that the ss element of the expression given by Eq. (15) obtained on the basis of second-order phase perturbation theory can be written as In writing this expression we have introduced the functions where The function M (q |k ) is given by [see the Appendix] where with sgn(·) denoting the sign function. We now turn to an evaluation of the ingredients in Eq. (16). We begin with the expression for 2M (q |k ) given by Eqs.
where φ q , φ p , φ k , and φ x are the azimuthal angles of the unit vectorsq ,p ,k , andx , respectively, measured from the positive x 1 axis [see Fig. 1]. On evaluating the angular integrals this result becomes where J n (z) is a Bessel function of the first kind and order n, and we have used the relation x sgn(x) = |x|. Finally, due to the circular symmetry of W (|u |) we obtain the result In the case of normal incidence (k = 0) and in-plane (q k ) scattering, Eq. (16) becomes and Eq. (22) can be written as From Eq. (11) it is noted that for normal incidencek = (cos φ 0 , sin φ 0 , 0), even if k /k is not well-defined in this case. Moreover, for scattering into directions that are normal to the mean surface we haveq =k , so for all directions q corresponding to in-plane scattering |q ·k | = 1. These results were used in arriving at the expression presented in Eqs. (24) and (25).

IV. THE INVERSE PROBLEM
To determine the function W (x ) from scattering data for ∂R ss (θ s )/∂Ω s incoh,input , we assume an analytic form for it that contains adjustable parameters. The values of these parameters, together with the rms height δ, are determined by minimizing a cost function with respect to variations of these parameters. The cost function we use is where P denotes the set of variational parameters used to characterize ∂R ss (θ s )/∂Ω s incoh,calc . The minimization of this function with respect to the elements of P was carried out by the use of the routine "lmdif1" contained in the Fortran package MINPACK which is part of the general purpose mathematical library SLATEC [9]. The routine lmdif1 implements a modified version of the Levenberg-Marquardt algorithm [10,11], and it calculates the Jacobian by a forward-difference approximation. The function ∂R ss (θ s )/∂Ω s incoh,input was obtained from rigorous, nonperturbative, purely numerical solutions [8,12] of the reduced Rayleigh equation for the scattering of polarized light from a two-dimensional randomly rough penetrable dielectric surface [13]. These calculations were carried out for an ensemble of random surfaces generated [12] on the basis of expressions for W (|x |) of either the exponential form or the Gaussian form In Eqs. (27), a denotes the transverse correlation length of the surface roughness. The function ∂R ss (θ s )/∂Ω s incoh,calc was obtained by evaluating the expression for it obtained using phase perturbation theory [Eq. (16)] for the trial function assumed to represent W (|x |). Several forms for this trial function were used in our calculations. In the first set of forms we assumed an exponential or Gaussian trial function, that is, or In this case the variational parameters of the reconstruction are δ , a , and potentially also ε . For the second set of forms for the trial function a stretched exponential was assumed which reduces to the exponential and Gaussian forms when γ = 1 and γ = 2, respectively. In this case the variational parameters of the reconstruction are δ , a , γ , and potentially ε .

V. RESULTS
We will now illustrate the inversion method developed here by applying it to the reconstruction of W (|x |), first by the use of one of the trial functions (28) and then by the use of the more general trial function (29).

A. Exponentially correlated surface roughness
For the first scattering system we consider, it is assumed that the surface height autocorrelation function W (|x |) is exponential, Eq. (27a), and characterized by a transverse correlation length a = 158.20 nm and an rms height of the surface δ = 9.50 nm. The medium above the surface is vacuum and the dielectric constant of the substrate is ε = 2.64 (photoresist). The wavelength (in vacuum) of the s-polarized incident light is λ = 632.8 nm. For this geometry and by the method of Ref. [12], we calculated the mean differential reflection coefficients by averaging the results from 5000 realizations of the surface profile function. For normal incidence, the in-plane, s-to-s co-polarized incoherent component of the mean differential reflection coefficient (DRC) obtained in this way is presented as a function of the scattering angle θ s by open circles in Fig. 2(a) [the same data set also appears in Figs. 3(a)-5(a)]. These data constitute the input function R ss (θ s )/∂Ω s incoh,input for our first set of reconstruction examples.
As our first example of reconstruction based on this data set, we assume that the trial function W (|x |) has the exponential form given by Eq. (28a). The set of variational parameters is therefore P = {δ , a }. The use of a mean differential reflection coefficient generated by the use of a known W (|x |) in our inversion approach enables us to assess the quality of the reconstructions we obtain. By starting the minimization procedure with the values δ = 2.00 nm and a = 75.00 nm, the values of these parameters that minimize the cost function χ 2 (P), Eq. (26), were found to be δ = 9.519 nm and a = 158.565 nm, to be compared with the values δ = 9.50 nm and a = 158.20 nm used to generate the input data. In the minimization procedure we assumed that both δ and a were restricted to positive values. The inversion is quite accurate. The function ∂R ss (θ s )/∂Ω s incoh,calc calculated with the reconstructed values of δ and a by means of second-order phase perturbation theory is plotted as the solid curve in Fig. 2(a) while the reconstructed correlation function W (|x |) is plotted as the solid curve in Fig. 2 Fig. 2(b)]. The shaded region in Fig. 2(b), and in the subsequent plots of W (|x |), represents the magnitude of the difference between the input and reconstructed values of this function. This difference is seen to be very small. In the preceding example it was assumed that the dielectric constant of the scattering medium was known. For our second example we take the input data from our first example, given by the open circles in Fig. 2(a), but now assume that together with the roughness parameters the dielectric constant of the substrate is unknown. Therefore the variational parameter set is now P = {δ , a , ε }. The results of this inversion are shown in Fig. 3, and it is seen that also in this case a rather good reconstruction is obtained. By starting the minimization procedure with the values δ = 2.00 nm, a = 75.00 nm, and ε = 2.00, the values of these parameters that minimize the cost function χ 2 (P) were determined to be δ = 9.272 nm, a = 158.042 nm, and ε = 2.718. These are to be compared with the input values δ = 9.50 nm, a = 158.20 nm, and ε = 2.64. The function ∂R ss (θ s /∂Ω s incoh,calc calculated with the reconstructed values of δ , a , ε by means of second-order phase perturbation theory is plotted as the solid curve in Fig. 3(a), while the reconstructed correlation function W (|x |) is plotted as the solid curve in Fig. 3(b). From a comparison of the results presented in Figs. 2 and 3 it is seen that the addition of a single variational parameter changes the reconstruction of W (|x |) only marginally.
A more stringent test of our inversion scheme is obtained when the trial function for W (|x |) has a functional form that differs from the form assumed in generating the input data that the reconstruction is based on. As our third example, we therefore present results of our calculations when the trial W (|x |) is assumed to have the stretched exponential form given by Eq. (29). The set of variational parameters is now P = {δ , a , γ }. By starting the minimization procedure with the values δ = 2.00 nm, a = 75.00 nm, and γ = 2.00, the values of these parameters that minimize the cost function were found to be δ = 9.425 nm, a = 161.717 nm, and γ = 1.012. These values are fairly close to the input values δ = 9.50 nm, a = 158.20 nm, and γ = 1.0 used in obtaining the simulation data. However, the importance of this example is to show that our minimization procedure is in fact able to distinguish a Gaussian form for the correlation function from an exponential form. In Fig. 4(a) we plot the function ∂R ss (θ s )/∂Ω s incoh,calc calculated by means of the second-order phase perturbation theory for the reconstructed values of δ , a , γ (solid curve), together with a plot of the input function (open circles). The agreement between these two results is quite good. In Fig. 4(b) we present plots of the input (open circles) and reconstructed (solid curve) correlation functions W (|x |). The latter curve very nearly coincides with the former curve. In our final example assuming an exponentially correlated surface, we again use the stretched exponential trial function for W (|x |), but now also assume that the dielectric constant of the scattering medium is unknown. The set of variational parameters is now P = {δ , a , γ , ε }. We start the minimization of the cost function χ 2 (P) with the values δ = 2.00 nm, a = 75.00 nm, and γ = 2.00, and ε = 2.00. The values of these parameters that minimize the cost function are found to be δ = 9.774 nm, a = 166.709 nm, γ = 1.027, and ε = 2.507. The proximity of these values to the input values, that are the same as those used previously, is poorer than for the first three examples, but the reconstructed values are still quite satisfactory. The reconstructed function ∂R ss (θ s )/∂Ω s incoh,calc [ Fig. 5(a)] and the reconstructed correlation function W (|x |) [ Fig. 5(b)], calculated with these reconstructed values are still in good quantitative agreement with the corresponding input functions. I. Summary of the scattering system parameters obtained during the different reconstruction scenarios based on inplane, s-to-s co-polarized scattering data corresponding to an exponentially correlated surface; δ , a , ε , and γ . The scattering system parameters assumed in generating the input data were: δ = 9.50 nm, a = 158.20 nm, ε = 2.64, and θ0 = 0 • . Note that an exponential correlation function corresponds to the exponent γ = 1 for the stretched exponential. The last column indicates the relevant figure where the results of the reconstruction in question is presented. The symbol "-" indicates that the corresponding variable was not reconstructed and instead had the value assumed in the input data (numerical simulations). In the two first reconstructions a trial correlation function of the form (28a) was used, while in the last two the form (29) was assumed. The parameters for the scattering system obtained in the different reconstruction scenarios detailed in this subsection are summarized in Table I.

B. Gaussianly correlated surface roughness
The second scattering system, from which we will use data for the purpose of inversion, is characterized by a Gaussian correlation function instead of the exponential correlation function assumed in the earlier set of examples. Compared to the previous scattering system, in addition to the different form of W (|x |), the only parameters that have changed are the rms-roughness and the dielectric constant of the substrate; they now take the values δ = 15.82 nm (an increase of more than 65% compared to its previous value), and ε = 2.6896, respectively. Except for the angles of incidence, all other parameters characterizing the scattering system remained unchanged, i.e., a = 158.20 nm and λ = 632.8 nm.
For these parameters, a computer simulation approach [12] was used to generate scattering data that were obtained by averaging the results from 24 000 surface realizations. Results obtained this way are presented as open symbols in Fig. 6(a) for the polar angle of incidence θ 0 = 50.2 • . It is these data we will base our inversion on in this subsection, that is, here this data set represents ∂R ss /∂Ω s incoh,input .
Motivated by the reconstruction done in Sec. V A using scattering data obtained for the exponentially correlated surface, we will now perform similar inversions for Gaussianly correlated surfaces using various variational parameter sets, P, that are subsets of {δ , a , ε , γ }. In such cases, the starting values assumed in the minimization will be {2 nm, 75 nm, 1, 2}, respectively, if nothing is said to indicate otherwise.   for W (|x |) is of the Gaussian form, Eq. (28b). In this way, the reconstruction procedure resulted in the numerical values δ = 15.873 nm and a = 158.000 nm. These values agree rather well with the values assumed in generating the scattering data used in the inversion. Moreover, the input and reconstructed correlation functions, as well as the absolute difference between them, are depicted in Fig. 6(b); the solid red line in Fig. 6(a) represents the mean DRC predicted by the inversion. We now continue by reconstructing the same variational parameter sets, P, as were used in Sec. V A for the exponential surface roughness; the only difference now is that we will assume the form of the trial function (28b) where we previously used Eq. (28a). A summary of the reconstructed parameters is presented in Table II. Moreover, in Fig. 7 the input and different reconstructed correlation functions obtained in this way are presented, together with comparisons of them. The mean DRC that result from these reconstructions are visually indistinguishable from those of Fig. 6(a), and such plots have therefore not been presented.
From Figs. 6-7 and Table II we see that the reconstructed results for the Gaussianly correlated surface roughness are in general good, at least for the angle of incidence θ 0 = 50.2 • assumed here. The quality of the results obtained are on a par with the results obtained previously for the exponentially correlated surface roughness, even though in the Gaussian case the rms-roughness is significantly larger and the angle of incidence is non-zero.
It should be mentioned that we have also performed reconstructions based on input data corresponding to other polar angles of incidence than θ 0 = 50.2 • , and it was found that the values of the reconstructed parameters essentially remained unchanged; if there were any changes at all, the quality of the reconstruction seemed to improve for smaller polar angles of incidence. Assuming different initial values of the parameters of the set P seemed not to affect the reconstruction. Hence, our results seem to indicate that reconstruction based on different input data, at least for the scattering data that we considered, influences the numerical values of the reconstructed parameters of the surface height correlation function only to a small extent. Furthermore, the reconstruction seems to be reliable for a wide range of angles of incidence.

C. Sensitivity to noise
Any experimental data set will contain some level of noise. Therefore, it is imperative to have reliable inversion approaches that can be applied successfully to data sets containing noise. Until now, we have used simulated data as the basis for the reconstruction, and the only source of uncertainty (or noise) in such simulations results is the finite number of surface realizations used to obtain them. However, since a sufficiently high number of realizations has been used in generating such data, the uncertainty has been modest. To start investigating the sensitivity of the reconstructed parameters to noise, we will, for reasons of comparison, base it on the data set used in Sec. V B [open symbols in Fig. 6(a)] and add noise to it. Due to the way these simulation results were generated [12], only 36 points exist in this data set. However, experimental angular resolved scattering data sets typically will have significantly better angular resolution (and, thus, more points). Therefore, to better mimic the more relevant experimental situation, we have interpolated by splines these simulation results [solid symbols in Fig. 8(a)] to an angular resolution of ∆θ s = 1 • for θ s in the interval from −90 • to 90 • . Then, to these (locally smooth) interpolated data, we have added multiplicative Gaussian white noise of a standard deviation of 5% and zero mean (gray erratic signal oscillating around zero in Fig. 8(a) resulting in the open blue data in Fig. 8(a). It is this latter data set that will be used as the noisy input signal for the reconstructions to be performed below.
Results for the reconstructed parameters based on this Gaussian white noisy data set, performed in a fashion identical to what was done in the preceding subsection, are listed in Table III. From this table, we observe that the results obtained are rather good, even for the significant noise level assumed. Moreover, one finds that the quality of  Table II. The shaded areas represent the absolute error between the input and reconstructed correlation functions.
the reconstruction is not dramatically degraded compared to what was previously obtained using the non-noisy data set [see Table II]. Not surprisingly, the poorest results of the inversion are obtained for the variational parameter set P of cardinality 4, and this "worst case" is presented as solid red lines in Fig. 8. The results presented in Table III and Fig. 8, which for the non-noisy case should be compared to Table II and Fig. 7(c), support the view that the reconstruction approach presented in this paper is able to produce reliable results also when the input data are noisy. For instance, when reconstructing 4 parameters, the relative error in the reconstructed correlation lengths are about 3.5% and −0.7% for the noisy and non-noisy case, respectively, which is not dramatic given the level of noise that was added to the input data. It should be mentioned that we found that if the parameters were estimated using a noisy version of the original  III. Same as Table II, but now the data set used in the inversion consisted of an interpolated version of the data set used to produce the results of Table II, with added multiplicative Gaussian white noise of a standard deviation of 5% (and mean zero). In Fig. 8, the noisy input data set is depicted as blue open symbols, while the noise appears in gray. scattering data instead of the interpolated data, as done above, the results using the trial function (28b) essentially remained unchanged. However, for the stretched exponential form (29), the results were more sensitive to the initial values used in the minimization, resulting in less robust results than those obtained by the use of the former data set.
We will see in the next subsection that this can also be the case when using multiple angles of incidence. The solid symbols represent the same data set that appears as symbols in Fig. 6(a). This (smooth) data set was first interpolated to an angular resolution of 1 • , and then multiplicative Gaussian white noise of 5% standard deviation was added to it. The open blue symbols represent the resulting noisy signal that the reconstruction is based on; the irregular signal in gray oscillating around zero is the actual noise being added. The (non-interpolated) original data set is indicated by green filled symbols. The solid red line represents the incoherent component of the mean DRC resulting from reconstructing the variational set P = {δ , a , ε , γ } (actual values given by the last line of Table III). (b) The input and reconstructed correlation function W (|x |) for variational parameter set P; the "non-noisy" equivalent of this graph is presented in Fig. 7(c). It should be noted that reconstruction using any subset of P results in more accurate results for the correlation function [see Table III], so what is shown here is indeed the "worst case" situation for the cases that we have considered.

D. Inversion of data obtained from multiple angles of incidence
When the experimental setup is prepared to measure the in-plane scattering of light for one angle of incidence, it is relatively straightforward to perform additional measurements for other angles of incidence. Therefore, it is of interest to study how the reconstructed parameters will depend on using multiple angles of incidence, and thus several data sets, in the inversion. In order to include, say, N data sets into the reconstruction, the cost function used in the minimization is generalized to where χ 2 n (P) is defined by Eq. (26) and corresponds to the contribution to the total cost function χ 2 (P) from a single angle of incidence θ Assuming the Gaussianly correlated surface roughness of the previous subsections, in-plane data for the mean DRCs were obtained from computer simulations for the polar angles of incidence θ 0 = 1.6 • , 25.3 • and 50.2 • [12]. A series of joint inversions were then performed based on the three resulting data sets seen as open symbols in Fig. 9(a). The parameters obtained by such an approach are presented in Table IV and the resulting mean DRCs obtained when reconstructing δ , a and ε are presented as solid lines in Fig. 9(a). A comparison of the results presented in Tables II and IV reveals that including additional data sets into the inversion scheme, at least for the data sets we used in obtaining these results, did not change the estimates of the parameters in any significant way. If there is any change, it may seem as if a multi-angle reconstruction may slightly improve the results when using the trial function (28b), while it becomes slightly worse for the trial function (29). However, before one can draw firm conclusions on this issue, proper estimates of the error bars associated with each of these parameters will have to be obtained, something that is outside the scope of the present work.
If now noise is added to these data sets after first interpolating them to a higher angular resolution, as was done above for the single data set corresponding to θ 0 = 50.2 • , then on performing multi-angle reconstructions based on the resulting data, the values presented in Table V are obtained [see also Fig. 9(b)]. The first thing to observe from Table V is that the multi-angle inversions based on the stretched exponential trial function (29) produce rather inaccurate results compared to inversions based on only one of these data sets [see Table III]. Potentially one could first estimate the exponent of the stretched exponential from a single data set, since we have found that it does not matter which of the data sets we use, and then use this value for γ as a fixed parameter in a multi-angle inversion. However, we will not consider this situation farther here, since we feel that it is probably more fruitful to consider alternative forms of the trial functions.
On the other hand, when using the Gaussian trial function (28b) in the minimization, a comparison of the results presented in Tables III and V shows that the multi-angle reconstructions are producing more accurate results than those obtained when the inversion is based on the single angle scattering data set corresponding to θ 0 = 50.2 • . We have also found this result to be true when any of the other two data sets were used in the single-angle inversion. This is an important result, since it may indicate that including several experimental data sets into the inversion may improve the estimates of the parameters. The results from Table V also hint at the importance of the choice taken for the trial function, since a priori it is not known which form of W (|x |) will result in a cost function, χ(P), that has a deep and well defined minimum, in contrast to many local minima of comparable depths.
It is also interesting to note that the results obtained for the multi-angle noisy case [ Table V] seem to be more accurate than those for the corresponding multi-angle non-noisy case [Table IV]. However, here it is important to recall that significantly more points are used in the inversion in the former than in the latter case, and we speculate that this could be the main reason for the improvement.   Table III, but now the reconstruction is based on several data sets obtained for the polar angles of incidence θ0 = 1.6 • , 25.3 • and 50.2 • . The data sets with the multiplicative 5% standard deviation Gaussian white noise data added to them, on which the reconstructions are based, are depicted as open symbols in Fig. 9 In principle, an inversion scheme, similar to the one proposed in this work, could be based on one of the rigorous simulation approaches that recently have become available to calculate light scattering from two-dimensional randomly rough surfaces [6,7,14]. However, for such an inversion scheme to be practically relevant, it is has to be numerically efficient, since during the inversion process, the forward scattering problem has to be solved for a large set of parameters. From this perspective, the rigorous numerical simulation approaches mentioned previously fall short since they typically require days of computer time, or more, to run just for one set of parameters.
In contrast, the phase perturbation theoretical approach to the forward scattering problem that we base our inversion scheme on, is computationally efficient. For instance, for an angular resolution of 1 • , it takes less than 2 s to obtain the mean DRC in the plane of incidence for one set of surface parameters (and one angle of incidence). Furthermore, the cpu times it took to perform the four inversions whose results are given in Table II were 25 s and 140 s when reconstructing 2 and 4 parameters, respectively [initial values as in Sec. V B]. The remaining two inversions required cpu times somewhere in between the two previously given cpu times. All reported computer cpu times are based on the use of a single computer core on an Intel i7 960 CPU operating at 3.20 GHz.
Even if such cpu times do depend strongly on the angular resolution of the data set used in the inversion, the amount of noise that it contains, and the initial values from which the reconstruction is started, these numbers for the computational cost do illustrate that the proposed inversion method is rather quick to perform; this, together with its robustness, should make it useful in practical situations.

VI. DISCUSSIONS AND CONCLUSIONS
In this paper we have shown that second-order phase perturbation theory can be used to model the incoherent scattering of light by a two-dimensional randomly rough dielectric surface and to calculate the mean differential reflection coefficient for such a surface. As a result it has been expected that it should be an effective tool for the inversion of light scattering data to obtain statistical properties of a random surface on which the mean differential reflection coefficient depends. This expectation has been borne out for weakly rough two-dimensional randomly rough dielectric surfaces. Together with several parametrized forms for the normalized surface height autocorrelation function W (|x |), namely an exponential, a Gaussian, and a stretched exponential, and the minimization of a cost function with respect to the parameters defining these forms, phase perturbation theory has been used in this paper to determine W (|x |), the rms height of the surface, the transverse correlation length of the surface roughness, the dielectric constant of the scattering medium when it was not known in advance, as well as the exponent of the stretched exponential. The function W (|x |) has been reconstructed quite accurately. The agreement between the reconstructed values of δ, a, and ε and the input values of these parameters is gratifyingly satisfactory.
This agreement remains very good when significant multiplicative Gaussian white noise is included in the input data. When simulated data for multiple angles of incidence are used in the inversion scheme, it is found that in the absence of noise in the input data the quality of the reconstructions is slightly poorer than when only a single angle of incidence is used. However, when noise is introduced into the input data the use of results obtained from multiple angles of incidence yields slightly better reconstructions than when data from only a single angle of incidence are used. The reasons for this behavior of the reconstructions is not understood, and deserve further study.
An investigation of the computational cost of the inversion approach developed here shows that it is quite computationally efficient compared to the several orders of magnitude higher computational cost of carrying out the inversions by the use of scattering data obtained by rigorous simulations. The cpu times required for carrying out an inversion calculation using our approach for a given angular resolution, one set of surface parameters, and a single angle of incidence, namely seconds, are short enough that this approach can be useful in practical situations.
The inversion approach developed here needs to be explored to determine ranges of roughness, wavelength, and dielectric parameters for which it gives reliable results. Error estimates for the reconstructed values of the roughness and material parameters sought should be obtained. The use of more flexible trial functions for W (|x |) in reconstructions should be explored, as well as the use of more than one wavelength. Reflectivity data can serve as the basis of an inversion scheme based on phase perturbation theory, and its utility for this purpose should be examined. This issue will be explored in subsequent work.