Abstract
A common approach to non-uniformity is to assume that the local thicknesses inside the light spot are distributed according to a certain distribution, such as the uniform distribution or the Wigner semicircle distribution. A model considered in this work uses a different approach in which the local thicknesses are given by a polynomial in the coordinates x and y along the surface of the film. An approach using the Gaussian quadrature is very efficient for including the influence of the non-uniformity on the measured ellipsometric quantities. However, the nodes and weights for the Gaussian quadrature must be calculated numerically if the non-uniformity is parameterized by the second or higher degree polynomial. A method for calculating these nodes and weights which is both efficient and numerically stable is presented. The presented method with a model using a second-degree polynomial is demonstrated on the sample of highly non-uniform polymer-like thin film characterized using variable-angle spectroscopic ellipsometry. The results are compared with those obtained using a model assuming the Wigner semicircle distribution.
© 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
Corrections
25 February 2020: A typographical correction was made to the funding section.
1. Introduction
Thin films prepared by various deposition techniques are almost never ideal films without any defects. One of the defects that is encountered very often is area non-uniformity of the films. For example, in plasma chemical techniques it is difficult to ensure that the deposition conditions are the same along the entire area of the sample. Even if the film is relatively uniform in the center of the substrate, it is almost always non-uniform near the edges and corners of the substrates because the electric fields are disrupted in these regions.
Thin films can exhibit non-uniformity in thickness as well as non-uniformity in the optical constants. However, the influence of non-uniformity is mostly through the optical thickness, thus, these two effects are highly correlated. This means that the model assuming only thickness non-uniformity is often sufficient to describe the effects of area non-uniformity on the measured optical quantities.
Variable-angle spectroscopic ellipsometry is a tool, which is often used for the optical characterization of various thin films. This technique is well suited for the study of area non-uniformity. Firstly, the area non-uniformity is manifested as depolarization and can be distinguished from many other frequently encountered defects, such as overlayers, transition layers, inhomogeneity, etc., which do not cause depolarization. Secondly, the changes in the incidence angle affect the size and shape of the light spot, i.e. the area probed by the ellipsometer is different for each incidence angle. This allows us to obtain additional information about the area non-uniformity.
The formulae needed to describe the area non-uniformity can be easily formulated using the Stokes–Mueller formalism. The Mueller matrix of the sample non-uniform along its surface can be calculated as [1]
where the integration is performed over the area of the light spot, $S$ is the size of this area and $M(x,y)$ is the local Mueller matrix corresponding to the given point on the surface of the sample. The local Mueller matrices can be calculated from the local reflection or transmission coefficients using formulae presented, for example, in [2,3]. We should note that the above formula was derived under the assumption that all points inside the light spot contribute with the same weight, or in other words, the illumination by the incident beam is assumed to be constant over the entire area of the light spot. In the case of thin films non-uniform in thickness, the local Mueller matrix depends only on the local thickness and it is possible to write where $h(x,y)$ is the local thickness at the given point on the surface of the sample. Instead of using this formula, which requires evaluation of a double integral, it is possible to use a formula which uses the probability density function for the distribution of local thicknesses inside the light spot $\rho (h)$ and requires only evaluation of a single integral Throughout this paper it will be assumed that the light beam of the ellipsometer has circular cross section, thus, the light spots on the sample are ellipses with the lengths of the semi-major and semi-minor axes given as where $\theta$ is the incidence angle in vacuum and $R$ is the radius of the light beam. The coordinates $x$ and $y$ along the surface of the sample are chosen such that the $x$-axis is parallel to the plane of incidence, the $y$-axis is perpendicular to this plane and the origin is located at the center of the elliptical light spot.The probability density function $\rho (h)$ can be calculated using the formula [4]
where the integration is performed over the curve $\mathcal {C}_h$ with constant thickness $h$ lying inside the light spot and $|\textrm {grad}\, h|$ denotes the size of the gradient at the given point on the curve. While this formula is very useful for deriving the distribution of thicknesses for some simple models, it becomes impractical if finding the contour lines is a difficult problem.In most cases, the thickness non-uniformity is very small and a model of wedge-shaped non-uniformity can be used. The probability density function can be easily derived using the formula (5) and it is given as [4]
A modification of the model suitable for describing thickness non-uniformity that slightly differs from the ideal wedge was presented in paper [11]. In this model it is assumed that for a given angle of incidence the distribution of thicknesses has the same form as for the wedge-shaped thickness non-uniformity, i.e. it is given by the Wigner semicircle distribution
The influence of the area non-uniformity on the ellipsometric measurements can be minimized by reducing the spot size of the ellipsometer or by orienting the sample such that the non-uniformity is largest along the direction where the light spot is shortest (i.e. along the semi-minor axis). However, it is not always possible to change the orientation of the light spot on the sample and the reduction of the spot size reduces the intensity of the signal registered by the detector and lowers the signal to noise ratio. It is often easier to incorporate the area non-uniformity into the processing of the experimental data than to change the arrangement of the experiment in the hope of eliminating its influence.
2. Polynomial thickness non-uniformity
The local thicknesses in the model of non-uniformity considered in this work are given by a second-degree polynomial in the coordinates $x$ and $y$
The inclusion of thickness non-uniformity into formulae for the optical quantities requires evaluation of either the integral (2) or the integral (3). The first form (2) requires the evaluation of a double integral, which could be numerically demanding. The second form (3) requires evaluation of a single integral, however, it requires the knowledge of the probability density function for the distribution of thicknesses, which could be difficult to determine for a general shape of non-uniformity.
The integral (3) can be evaluated very efficiently using the Gaussian quadrature. The $n$ point Gaussian quadrature approximates the integral by a sum
where the weights $w_k$ and the nodes $h_k$ are chosen such that the quadrature gives exact results if the elements of the Mueller matrix $M(h)$ are polynomials of degree up to $2n-1$. Of course, the values of the weights and the positions of the nodes depend on the course of the function $\rho (h)$.For a probability density function corresponding to a wedge-shaped non-uniformity, the Gaussian quadrature is known as the Chebyshev–Gauss quadrature and the nodes and weights can be in this case calculated by means of the closed-form formulae [7,8]. For a general shape of non-uniformity the nodes and weights must be calculated numerically. Although this is more computationally demanding than the use of the closed-form formulae, it could still be used as a basis for a very efficient algorithm. Because the nodes and weights depend only on the shape of the non-uniformity and incidence angle, they could be calculated once and then reused for all points in the studied spectra.
A general algorithm for determining the nodes and weights of the Gaussian quadrature defines a sequence of orthogonal polynomials with the orthogonality relation given by the inner product
where the symbols $f(h)$ and $g(h)$ are used to represent functions of variable $h$. The function $\rho (h)$, which is in this context called the weight function, must be non-negative. This is not a problem since $\rho (h)$ is a probability density function, thus, it cannot be negative. The nodes $h_k$ for the $n$-point Gaussian quadrature rule are then chosen as the roots of the $n$-th degree polynomial in this sequence of orthogonal polynomials. Once we know the positions of the nodes it is not difficult to find the values of the weights using simple linear algebra.The inner product can also be written as the surface integral over the area of the light spot
This form is advantageous since it does not require the knowledge of the distribution function $\rho (h)$. The surface integral in the above formula with the function $h(x,y)$ given by the polynomial (9) can be evaluated using the method described in appendix A. The numerical procedure for finding the nodes and weights for the Gaussian quadrature using this method to evaluate the surface integrals and the Golub–Welsch algorithm is presented in appendix B.Although it was assumed that the shape of non-uniformity is described by a polynomial with at most quadratic terms, it is not difficult to generalize the results to models describing the non-uniformity by polynomials of higher degrees. This generalization is described in appendix C.
3. Experiment
The method will be demonstrated on a sample of a polymer-like thin film exhibiting relatively large thickness non-uniformity deposited onto a silicon substrate.
3.1 Sample preparation
The film was deposited by the plasma enhanced chemical vapor deposition (PECVD) onto a double-side polished silicon single crystal wafer. The deposition was performed inside the parallel-plate PECVD reactor using capacitively coupled radiofrequency glow discharge in a mixture of oxygen and trimethylsilyl acetate (TMSA – C$_5$H$_{12}$O$_2$Si). Prior to the deposition, the wafer was cleaned in the argon plasma to increase the adhesion of the deposited film. In order to create a film with large thickness non-uniformity, a small part of the silicon wafer was covered by an obstacle (a small steel plate) during the deposition. The obstacle disrupted the homogeneity of the discharge and the parts of the film located near the obstacle then exhibited large non-uniformity.
3.2 Experimental data
The sample was measured using the Horiba Jobin Yvon UVISEL phase modulated ellipsometer for five incidence angles 55$^{\circ }$, 60$^{\circ }$, 65$^{\circ }$, 70$^{\circ }$, 75$^{\circ }$ in the spectral range 0.6–6.3 eV (197–2066 nm). The ellipsometric experimental data are represented by the quantities $I_\textrm {s}$, $I_\textrm {c}$ and $I_\textrm {n}$, which correspond to independent elements of the normalized Mueller matrix of isotropic systems. In particular, for reflected light, the sample is characterized by the Mueller matrix
3.3 Structural model
The polymer-like film was modeled as a homogeneous film non-uniform in thickness. The thickness non-uniformity was described using the model introduced in section 2. It was also necessary to assume that there is a thin transition layer between the polymer-like film and the silicon substrate. The transition layer was also modeled as a homogeneous film. It was assumed that the thickness of the transition layer is the same in the entire area measured by the ellipsometer, i.e. it is unaffected by the thickness non-uniformity. Thus, the local reflection coefficients needed to calculate the Mueller matrices in (10) correspond to a two-layer system with the thickness of the layer representing the polymer-like film given by the function $h(x,y)$ and the thickness of the transition layer fixed in value $h_\textrm {tl}$.
Since silicon becomes transparent in the infrared region, it was necessary to take into account reflections on the back side of the double-side polished silicon wafer. This effect, which is an additional source of depolarization in the region below $\approx$ 1.2 eV, was incorporated into the formulae for the ellipsometric quantities using the method described in [3].
The finite spectral resolution and the beam divergence of the used ellipsometer are another effects contributing to depolarization. Their influence was taken into account by averaging the Mueller matrices over the distribution of frequencies around the nominal value due to the limited resolving power of the monochromator and the distribution of incidence angles in the divergent beam [3]. However, the depolarization due to these effects is much smaller than the depolarization caused by the thickness non-uniformity of the film.
3.4 Dispersion models
The optical constants of the polymer-like film and the transition layer were modeled using the universal dispersion model (UDM) [12] implemented in the newAD2 software [13]. The dielectric response in the UDM is described by Kramers–Kronig consistent models derived on the basis of the parameterization of the joint density of states. The UDM allows us to describe the dielectric response of a wide range of materials by combining contributions representing various elementary excitations in condensed matter, e.g. the interband transitions, absorption on localized states etc.
The dielectric function describing the polymer-like film is given by two contributions. The first contribution describes the interband transitions and the second contribution models the weak absorption below the fundamental absorption edge (band gap energy). The determined spectral dependencies of the optical constants of the non-uniform polymer-like film are plotted in Fig. 2(a).
The Campi–Coriasso model [14] is used to model the interband transitions. Since the fundamental absorption edge is located almost at the end of the measured spectral range, the choice of the model of interband transitions is not very important. Other models combining the Lorentz model and the Tauc’s law [15] or the model in which the interband transitions are modeled by a broad absorption band [12,16] could also be used instead of the Campi–Coriasso model. The fact that these models contain enough parameters to accurately model the dispersion below the fundamental absorption edge is more important than the concrete shape of the absorption structure representing interband transitions.
The weak absorption below the fundamental absorption edge is described using a contribution that models the imaginary part of the dielectric function as decaying exponential in the region below the fundamental absorption edge. Although the absorption below the fundamental absorption edge is weak (the extinction coefficients is only several thousandths), it affects the value of the refractive index through the Kramers–Kronig relations. The small changes in the value of the refractive index are more important than almost negligible attenuation of the electromagnetic waves. If this part of the model was omitted, it would not be possible to correctly match the positions of interference extrema, which are sensitive to the value of the refractive index. The used model is very similar to the model of the Urbach tail presented in [12,16]. However, in contrast to the model of the Urbach tail which assumes zero absorption below half the band gap energy, it models nonzero absorption down to zero frequency. The details concerning this model of the exponential tail have not been published so far, they will be presented in our upcoming paper.
The dielectric response of the transition layer was described using the model combining two Campi–Coriasso terms and the contribution representing the exponential tail. The determined spectral dependencies of the optical constants of the transition layer are plotted in Fig. 2(b). The thickness of the transition layer was determined in the value $12.0 \pm 0.9$ nm. The transition layer probably represents the upper layer of the silicon wafer perturbed by the cleaning in the argon plasma. Of course, there could also be an inhomogeneous part of the polymer-like film created in the initial phase of the deposition which would also manifest in this transition layer. The spectral dependencies of the optical constants of the transition layer are similar to those of the amorphous silicon (see e.g. [17]), however, the absolute value of the refractive index is lower and the band gap energy is slightly shifted to higher energies. The two peaks visible in Fig. 2(b) evidently correspond to the peaks at $E_1$ and $E_2$ critical points energies in crystalline silicon (see e.g. [18]) but they are much less pronounced than in crystalline silicon.
The optical constants of the silicon single crystal substrate were not sought within the data processing, it was assumed that they are the same as those determined in [19].
3.5 Results
The experimental ellipsometric data were processed using the least-squares method. The agreement between the theoretically calculated curves and the experimental data is shown in Fig. 1. The values of the parameters describing the thickness non-uniformity are presented in Table 1(a). The uncertainties of the parameters introduced in this table are the uncertainties determined by the least-squares method. In order to compare the quality of the fits corresponding to different models of thickness non-uniformity, it is convenient to define a quantity $\chi$, which assesses goodness of the fits
where $\Sigma$ is the weighted sum of squared residuals and $N$ is the number of experimental values. Lower values of this quantity correspond to better fit of the experimental data. The value of this quantity for the results presented in Figs. 1, 2 and Table 1(a) is $\chi =1.76$.The value of the coefficient $h_\textrm {xx}$ is negative while the value of the coefficient $h_\textrm {yy}$ is positive. This means that the surface of the film is concave along the $x$-axis and convex along the $y$-axis, which suggests that it exhibits a saddle-like shape.
The parameters $h_\textrm {y}$, $h_\textrm {xy}$ and $h_\textrm {yy}$ were determined with much larger uncertainties than the parameters $h_\textrm {x}$ and $h_\textrm {xx}$. This is because the light spot of the ellipsometer is much longer along the $x$-axis than along the $y$-axis and, moreover, the dependence of the thickness distribution on the angle of incidence is determined mostly by the shape of the non-uniformity along the $x$-axis. In fact, if a model assuming only non-uniformity along the $x$-axis described by only two parameters $h_\textrm {x}$ and $h_\textrm {xx}$ is used (i.e. the other parameters $h_\textrm {y}$, $h_\textrm {xy}$ and $h_\textrm {yy}$ are set to zero), then the quality of the fit of the experimental data corresponds to $\chi =1.97$. This is not much worse fit than in the model using all five parameters $h_\textrm {x}$, $h_\textrm {y}$, $h_\textrm {xx}$, $h_\textrm {xy}$ and $h_\textrm {yy}$ which results in $\chi =1.76$.
The probability density functions for the distribution of thicknesses inside the light spot of the ellipsometer at three selected incidence angles are plotted in Fig. 3(a). The calculated nodes and weights for the Gaussian quadrature corresponding to these distributions are displayed in Fig. 3(b).
In order to estimate the accuracy of the Gaussian quadrature achieved for a given number of points, the values of the ellipsometric quantities calculated with integrals (3) evaluated using the Gaussian quadrature were compared with values calculated using method which directly evaluated the integrals (1) using the Romberg method. The maximal $\Delta _\textrm {max}$ and RMS $\Delta _\textrm {RMS}$ errors of the ellipsometric quantities $I_\textrm {s}$, $I_\textrm {c}$, $I_\textrm {n}$ achieved for a given number of points used in the Gaussian quadrature are introduced in Table 1(b). These errors were calculated using five incidence angles 55$^{\circ }$, 60$^{\circ }$, 65$^{\circ }$, 70$^{\circ }$, 75$^{\circ }$ and 229 spectral points in the range 0.6–6.3 eV, i.e. the same collection of points that was used when the ellipsometric quantities were experimentally measured. The calculation of the errors was performed with the values of the structural and dispersion parameters fixed at values corresponding to the best fit of the experimental data. The results presented in Figs. 1, 2 and Table 1(a) were obtained using the Gaussian quadrature rule with 15 points.
4. Comparison with other models of thickness non-uniformity
In this section, the results achieved using the model describing the non-uniformity using the polynomial (9) will be compared with the results achieved by another three models. These models are:
The probability density functions for the distribution of heights corresponding to the empirical model using the series in $1/\cos ^{2}\theta$ and the model of a wedge-shaped non-uniformity are plotted in Fig. 4(a). The values of the mean thickness $\bar {h}(\theta )$ and the RMS value of thickness deviations from the mean value $\sigma (\theta )$ corresponding to various models of non-uniformity are presented in Table 2. The last row in this table shows the quantity $\chi$ defined in (15), which assesses the goodness of the fits. From the values of the quantity $\chi$ it is evident that the model of a wedge-shaped non-uniformity gives only slightly worse result than the model using the series in $1/\cos ^{2}\theta$ and the model in which $\bar {h}(\theta )$ and $\sigma (\theta )$ are sought independently for each incidence angle. The fits obtained by these three models are considerably worse than the fit obtained using the model describing the shape of non-uniformity by the polynomial (9). An unsatisfactory agreement between the experimentally determined values of the degree of polarization and the values calculated using the empirical model using the series in $1/\cos ^{2}\theta$ is also visible in Fig. 4(b).From Table 2 it is evident that the values corresponding to the model which determines $\bar {h}(\theta )$ and $\sigma (\theta )$ using the series (8) in $1/\cos ^{2}\theta$ are very close to those determined if these values are sought independently for each incidence angle. This means that the model using the series in $1/\cos ^{2}\theta$ contains enough parameters to model the shifts in the mean thickness and changes in the RMS value of deviations. Thus, the worse fit of the experimental data in the models using the Wigner semicircle distribution is not caused by the inability to take into account the changes in these quantities. In order to obtain better fits of the experimental data, it is necessary to correctly model the course of the probability density function for the distribution of local thicknesses.
5. Conclusion
The model of thickness non-uniformity in which the local thicknesses are given by a polynomial in the coordinates $x$ and $y$ along the surface of the film was considered. The formulae used to include the effect of thickness non-uniformity on the ellipsometric quantities uses the Gaussian quadrature with the weight function given as the probability density function for the thickness distribution inside the light spot of the instrument. An elliptical light spot, which corresponds to ellipsometer with a circular cross section of the beam, was considered in this work. If the local thicknesses are given using the second or higher degree polynomial then it is difficult to find the corresponding probability density function. Fortunately, it is possible to find the nodes and weights for the Gaussian quadrature without explicitly knowing this probability density function. An efficient method for calculating these nodes and weights was presented. This method is also numerically stable, meaning that it behaves well with respect to the propagation of the rounding errors. Although most of the results are presented for the model using the second-degree polynomial, it was shown how these results could be generalized to polynomials of higher degrees.
The usability of the method was demonstrated on the optical characterization of a highly non-uniform polymer-like film. The studied sample was measured using the variable-angle spectroscopic ellipsometry. The variable-angle spectroscopic ellipsometry is suitable for studying the non-uniformity since the size of the light spot, and as a consequence the distribution of thicknesses in this light spot, changes with the angle of incidence.
Most of the papers dealing with thickness non-uniformity assume that the local thicknesses inside the light spot are distributed according to a certain distribution (see e.g. [5,6,9–11]). A typical examples are the uniform distribution or the Wigner semicircle distribution, which corresponds to a wedge-shaped non-uniformity. The advantage of the model using the Wigner semicircle distribution is that the nodes and weights for the Gaussian quadrature are given by the closed-form formulae. A model of thickness non-uniformity using this distribution presented in [11] was also used to process the experimental data. While the model using the distribution corresponding to the non-uniformity specified by the second-degree polynomial gave excellent results, the model using the Wigner semicircle distribution was unable to describe the non-uniformity even if several parameters modeling the dependence of the mean value and variance of this distribution on the angle of incidence were used.
In most cases, the thickness non-uniformity represents a small defect and the choice of the probability density function representing the distribution of thicknesses inside the light spot of the instrument is not crucial. However, as was demonstrated in this work, if the shape of the non-uniformity is more complicated, it is better to use the distribution of local thicknesses derived on the basis of the real shape of non-uniformity than to assume some chosen distribution function. Although the nodes and weights for the Gaussian quadrature must be calculated numerically in this case, the resulting method is not that much complicated. Moreover, since the calculated nodes and weights are reused for all spectral points, the method is also very efficient.
Appendix A. Evaluation of the inner product
An efficient method for evaluatiing the surface integral defining the inner product (12) will be presented in this appendix. As will be shown in appendix B., in order to find the nodes and weights for the Gaussian quadrature using $n$ points, it is necessary to evaluate scalar products between polynomials $f(h)$ and $g(h)$ with the sum of their degrees at most $2n - 1$ (i.e. the same degree that can be integrated exactly using the Gaussian quadrature). Thus, it will be assumed that the integrand $F(h(x,y))=f(h(x,y))g(h(x,y))$ in (12) is a polynomial of degree at most $2n - 1$ in $h(x,y)$.
The ellipse representing the domain of the surface integral can be parameterized by scaling the $x$ and $y$ axes such that we obtain the unit circle and then using the polar coordinates. In this parameterization it holds
where $r\in [0,1]$ and $\phi \in [0,2\pi ]$ are the coordinates parameterizing the unit circle. Using this substitution, the integral (12) can be written asThe integral over the angular coordinate in (17) can be replaced by a sum over the discrete set of equally spaced values. If $N_\phi$ is a positive integer and $\vert q\vert\;<\;N_\phi$, then it holds
The results for the integrals over the radial and angular coordinates can be combined into the following formula for the inner product
Appendix B. Calculation of the nodes and weights for the Gaussian quadrature
It can be shown that the orthogonal polynomials satisfy the recurrence relation [7]
where the coefficients $\alpha _j$, $\beta _{j-1}^{2}$ can be calculated using the inner product asThe inner products in (25) can be calculated using the formula (23) provided we know the values of $p_j(h(r,\phi ))$ on the net consisting of points $r_k$ and $\phi _m$. In order to make our formulae shorter, these values will be denoted as
If we combine the formula (23) for the inner product with (24) and (25), we can formulate the recurrence procedure for the calculation of the coefficients $\alpha _j$ and $\beta _{j-1}$.At the beginning we start with
Note that the Golub–Welsch algorithm can also be used to calculate the weights $\eta _k$ and the nodes $r_k$ for the Gaussian quadrature used in (22) and (23) to calculate the integral over the radial coordinate. However, in this case the tridiagonal matrix (26) has dimension $N_\textrm {r} \times N_\textrm {r}$, the elements of this matrix are given as
and the zeroth moment is equal to $\mu _0=1/2$.Appendix C. Generalization to polynomials of higher degrees
It is also straightforward to generalize the presented formulae to polynomials $h(x,y)$ of higher than the second degree. All that is needed is to adjust the numbers $N_\textrm {r}$ and $N_\phi$ such that the formula (23) gives exact results and change the formula (29) for $\alpha _1$ so that it takes into account higher order terms. By repeating the considerations concerning the number of points required for the angular coordinate in (21) and for the Gauss-Jacobi quadrature in (22) with general degree of polynomial $h(x,y)$, we conclude that it is necessary to choose
Funding
Ministerstvo Školství, Mládeže a Tělovýchovy (LM2018097); Grantová Agentura České Republiky (GACR 19-15240S).
Disclosures
The authors declare no conflicts of interest.
References
1. R. Ossikovski, M. Kildemo, M. Stchakovsky, and M. Mooney, “Anisotropic Incoherent Reflection Model for Spectroscopic Ellipsometry of a Thick Semitransparent Anisotropic Substrate,” Appl. Opt. 39(13), 2071–2077 (2000). [CrossRef]
2. E. Garcia-Caurel, R. Ossikovski, M. Foldyna, A. Pierangelo, B. Drévillon, and A. De Martino, “Advanced mueller ellipsometry instrumentation and data analysis,” in Ellipsometry at the Nanoscale, M. Losurdo and K. Hingerl, eds. (Springer, 2013).
3. I. Ohlídal, J. Vohánka, M. Čermák, and D. Franta, “Ellipsometry of layered systems,” in Optical Characterization of Thin Solid Films, O. Stenzel and M. Ohlídal, eds. (Springer International Publishing, 2018).
4. D. Nečas, I. Ohlídal, and D. Franta, “The reflectance of non-uniform thin films,” J. Opt. A: Pure Appl. Opt. 11(4), 045202 (2009). [CrossRef]
5. I. Ohlídal, D. Nečas, D. Franta, and V. Buršíková, “Characterization of non-uniform diamond-like carbon films by spectroscopic ellipsometry,” Diamond Relat. Mater. 18(2-3), 364–367 (2009). [CrossRef]
6. I. Ohlídal, D. Franta, and D. Nečas, “Ellipsometric and reflectometric characterization of thin films exhibiting thickness non-uniformity and boundary roughness,” Appl. Surf. Sci. 421, 687–696 (2017). [CrossRef]
7. M. Abramowitz and I. A. Stegun, eds., Handbook of mathematical functions with Formulas, Graphs, and Mathematical Tables (National Bureau of Standards, 1964).
8. A. Ralston and P. Rabinowitz, A First Course in Numerical Analysis (2nd ed.) (Dover Publications, 2001).
9. T. Pisarkiewicz, “Reflection spectrum for a thin film with non-uniform thickness,” J. Phys. D: Appl. Phys. 27(1), 160–164 (1994). [CrossRef]
10. U. Richter, “Application of the degree of polarization to film thickness gradients,” Thin Solid Films 313-314, 102–107 (1998). [CrossRef]
11. D. Nečas, I. Ohlídal, and D. Franta, “Variable-angle spectroscopic ellipsometry of considerably non-uniform thin films,” J. Opt. 13(8), 085705 (2011). [CrossRef]
12. D. Franta, J. Vohánka, and M. Čermák, “Universal dispersion model for characterisation of thin films over wide spectral range,” in Optical Characterization of Thin Solid Films, O. Stenzel and M. Ohlídal, eds. (Springer International Publishing, 2018).
13. D. Franta, D. Nečas, and J. Vohánka, “Software for optical characterization newAD2,” http://newad.physics.muni.cz/ (2014).
14. D. Campi and C. Coriasso, “Prediction of Optical Properties of Amorphous Tetrahedrally Bounded Materials,” J. Appl. Phys. 64(8), 4128–4134 (1988). [CrossRef]
15. D. Franta, M. Čermák, J. Vohánka, and I. Ohlídal, “Dispersion models describing interband electronic transitions combining Tauc’s law and Lorentz model,” Thin Solid Films 631, 12–22 (2017). [CrossRef]
16. D. Franta, D. Nečas, and I. Ohlídal, “Universal dispersion model for characterization of optical thin films over wide spectral range: Application to hafnia,” Appl. Opt. 54(31), 9108–9119 (2015). [CrossRef]
17. D. T. Pierce and W. E. Spicer, “Electronic structure of amorphous Si from photoemission and optical studies,” Phys. Rev. B 5(8), 3017–3029 (1972). [CrossRef]
18. H. R. Philipp and E. A. Taft, “Optical constants of silicon in the region 1 to 10 eV,” Phys. Rev. 120(1), 37–38 (1960). [CrossRef]
19. D. Franta, A. Dubroka, C. Wang, A. Giglia, J. Vohánka, P. Franta, and I. Ohlídal, “Temperature-dependent dispersion model of float zone crystalline silicon,” Appl. Surf. Sci. 421, 405–419 (2017). [CrossRef]
20. G. H. Golub and J. H. Welsch, “Calculation of Gauss quadrature rules,” Math. Comp. 23(106), 221 (1969). [CrossRef]