Accepted Manuscript Spectroscopic ellipsometry data analysis using penalized splines representation for the dielectric function
D.V. Likhachev PII: DOI: Reference:
S0040-6090(18)30741-7 https://doi.org/10.1016/j.tsf.2018.10.057 TSF 36974
To appear in:
Thin Solid Films
Received date: Revised date: Accepted date:
3 May 2018 21 September 2018 31 October 2018
Please cite this article as: D.V. Likhachev , Spectroscopic ellipsometry data analysis using penalized splines representation for the dielectric function. Tsf (2018), https://doi.org/ 10.1016/j.tsf.2018.10.057
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT
Spectroscopic ellipsometry data analysis using penalized splines representation for the dielectric
PT
function
RI
D.V. Likhachev*
SC
GLOBALFOUNDRIES Dresden Module One LLC & Co. KG, Wilschdorfer Landstr. 101,
NU
D-01109 Dresden, Germany
MA
Keywords Spectroscopic ellipsometry, data analysis, optical modeling, dielectric function, parameterization, penalized splines, optical metrology
ED
Accuracy of the spectroscopic ellipsometry data analysis strongly depends on appropriate
EP T
modeling of the complex dielectric function ε (or the complex index of refraction N) over required spectral range. In this paper, we outline penalized B-spline (P-spline) formulation for an arbitrary dielectric function ε modeling in spectroscopic ellipsometry data analysis. The
AC C
main idea is to first use a generous number of equally-spaced knots (up to the number of spectral points) to describe even highly complicated structure in the imaginary part ε2 of the dielectric function and then provide a certain penalty on the coefficients of adjacent B-splines to tune the smoothness of the ε curve. Two real-data applications have been provided to evaluate the practical performance and effectiveness of the proposed penalized formulation for dielectric function representation. A comparison of obtained results with the findings of our previous studies without penalization demonstrates a good agreement between the proposed
* Electronic mail:
[email protected]. 1
ACCEPTED MANUSCRIPT method and ordinary B-spline representation in an equally-spaced fashion with optimal number of knots.
1. Introduction Accuracy of the spectroscopic ellipsometry data analysis strongly depends on appropriate
PT
modeling of the complex dielectric function ε = ε1 − iε2 (or the complex index of refraction N = n − ik) over required spectral range. There are several common ways of representing the
RI
dielectric function ε (or the index N) in the data analysis such as tabulated data for the n&k’s
SC
of each layer, effective medium approximation or physics-based analytical models [1–5]. A
NU
less common, but an elegant and advantageous approach to describe an arbitrary dielectric function is to use the polynomial-based models [6–12]. One of the polynomial-based models
MA
is the dielectric function representation by means of B-splines [13], a purely mathematical way of dealing with the optical properties modeling which uses a class of fairly easy to
ED
construct spline functions called B-splines (or Basis-splines) [14]. B-splines of various orders are piecewise polynomial functions with local support designed for approximation of
EP T
arbitrary functions and the order of the B-spline resembles the degree of the polynomials plus one. These polynomials are smoothly and continuously connected together at the joining
AC C
points the abscissas of which are called “knots”. By summing different numbers of shifted Bspline basis functions with unequal weights we can obtain the curves of various shapes and complexity (see Fig. 1 as an example).
To establish a model for the dielectric function one can parametrize the imaginary part ε2 (or k) by B-spline curve and then the real part ε1 (or n) can be derived using the Kramers-Kronig transformation. In our case, a B-spline curve S(x), which represents ε2 , is uniquely
2
ACCEPTED MANUSCRIPT constructed as a linear combination (weighted sum) of localized B-spline basis functions Bip of degree p:
S ( x ) bi Bip ( x ),
(1)
PT
i
where bi denotes the weight (spline coefficient) for the ith basis function of pth degree which
RI
can be recursively calculated from the lower-order functions using the Cox-de Boor recursion
NU
SC
formula [14]:
ti p1 x p 1 x ti p 1 Bi ( x ) B ( x ), p 1. ti p ti ti p1 ti 1 i 1
(2)
ED
Bip ( x )
MA
1, if ti x ti 1 Bi0 ( x ) , 0, otherwise
EP T
Here t = (t i) (i = 0, …, l) is the knot vector specified in wavelength or photon energy, i.e., a non-decreasing sequence of real numbers of length at least p + 2 (t 0 ≤ t 1 ≤ … ≤ t i ≤ t i+1 ≤ …
AC C
t m ) (possible divisions by zero are resolved by convention that “anything divided by zero is zero”, i.e., 0/0 = 0). If the number and the placement of the knots (the knot vector) were chosen somehow, then the weights bi are the only parameters which need to be estimated from a least-squares fitting of appropriate optical model to the spectroscopic ellipsometry data. Also, the weights bi are the ordinates of so-called “control points” or “de Boor points”.
3
SC
RI
PT
ACCEPTED MANUSCRIPT
NU
Fig. 1. Example of shifted B-spline basis functions of degree p = 3 over the interval [0, 8] with equally-spaced internal knots (uniform B-splines) and different weights (black dashed
MA
and blue solid lines) and their linear combination (red solid line). The knots are indicated by black dots. Such equidistant knot arrangement can be justified if the sample size, i.e., the
ED
number of data points, is not small comparing to the size of the knot vector and the data
EP T
points are distributed uniformly in the considered spectral range.
Some useful properties of B-splines are the following. First of all, the basis functions of
AC C
degree p have support only over p + 1 intervals, which contains p + 2 consecutive knots, and they are exactly zero outside the finite range [t i, t i+p+1 ]. In other words, each of the basis functions has compact support, as in the example illustrated in Fig. 1: the B-spline (marked by blue color) is non-zero over the interval ]t 3 = 1, t 7 = 3[ and is equal to zero outside it. As a direct consequence of this compact support property, changes in the spline coefficients (or/and knot positions, in more general case of non-uniform B-splines) modify the splinecurve shape only locally and this fact directly relates to the numerical stability and efficiency of the B-splines evaluation. Another attractive property, the “convex hull” property [15],
4
ACCEPTED MANUSCRIPT offers a simple and well-defined control over the sign of the imaginary part of the dielectric function ε2 : in order to ensure that ε2 ≥ 0, all spline coefficients bi must be non-negative. This property,
together with existing Kramers-Kronig consistent B-spline formulation [13],
ensures physical validity of the (ε1 ,ε2 )-curve shape in the considered spectral range. A comprehensive description of B-splines and their properties can be found in Ref. [14]. Ref.
PT
[13] also provides a glimpse into many of the practical issues that accompany use of B-spline
RI
parameterization of the dielectric functions including, in particular, knot termination and
EP T
ED
MA
NU
SC
possible issues with presence of absorption features outside of the spectral range.
AC C
Fig. 2. Influence of the number of knots on the B-spline representations of the TiN complex index of refraction (the knots are equally spaced in the wavelength range of 190–1000 nm). One can clearly see that the B-spline model with 35 knots is too flexible and overparameterized the experimental data resulting in a few numerical artifacts at ~200, 245, 320, 440, 730 nm due to small knot spacing. Reduction of the number of knots suppresses this “wiggly” overfitting behavior.
5
ACCEPTED MANUSCRIPT B-splines have already shown their
effectiveness for material modeling in multiple
spectroscopic ellipsometry applications [16–27]. But the B-spline representation of the dielectric functions also has two main disadvantages. First of all, since the dielectric function representation by B-splines is not a physics-based model, an extraction of some physical information about the material under study is not straightforward. As a consequence,
PT
additional effort is required to determine, for instance, bandgap values, just as it has been
RI
carried out for the point-by-point extracted imaginary part of the dielectric function (ε2 ),
SC
extinction coefficient (k) and the absorption coefficient (α) using various linear extrapolation methods [28]. Another important problem arises from the fact that the actual accuracy of B-
NU
spline parameterization in ellipsometric data analysis strongly depends on the number of knots (and, in general, their locations) used to describe the dielectric function. Using fewer
MA
knots in the considered wavelength range can lead to data underfitting and some relevant spectral features will be neglected. On the other hand, too many knots within the wavelength
ED
range beyond the optimal value put us in danger to overfit the data and, hence, to catch random noise in the data or unrealistic curve bending. We illustrate this issue by B-spline
EP T
representation of the complex index of refraction of thin titanium nitride film (Fig. 2). Therefore, to restrict flexibility of the B-spline representation and prevent overfitting, a
AC C
choice of the optimal number of knots is one of the main challenges in practical implementation of the B-spline approach and a very subjective and empirically-driven task. In previous papers [29,30], we used well-established and useful Akaike and Bayesian information criteria (AIC, BIC) for choosing a suitable B-spline model with optimal number of knots to balance fidelity to the measured data (goodness of fit) and parametric complexity of the fit.
6
ACCEPTED MANUSCRIPT An alternative method to avoid overfitting by B-spline parameterization is to use one of the variants of smoothing splines called “penalized splines”, or simply “P-splines”, which include a penalty on B-spline coefficients bi to ensure that adjoining coefficients do not differ too much from each other and produce less pliant fitted curve [31–34]. This paper focuses on application of the P-spline approach to dielectric function representation in spectroscopic
PT
ellipsometry data analysis. First, since many ellipsometry users are not very familiar with the
RI
P-spline method, an outline of penalized smoothing with B-splines will be provided.
SC
Secondly, some complications in practical implementation of the penalization will be briefly discussed. Finally, two real-data applications will be provided to evaluate the practical
NU
performance and effectiveness of the proposed penalized fitting of the spectroscopic
MA
ellipsometry data.
2. General overview of ellipsometric data fitting with P-splines
ED
In this section, the dielectric function representation by P-splines will be briefly described, covering all key aspects of our approach. A comprehensive description of penalized spline
EP T
smoothing methods is given in Refs. [35–38].
AC C
2.1. P-spline parameterization in a nutshell Although a penalized optimization, that balances a quality-of-fit measure (for example, the mean squared error, MSE) with a roughness penalty, along with selection of optimal degree of penalization have been developed for a long time (see, e.g., Refs. [39–41]), the P-spline approach itself has become a popular technique for smoothing in the last two decades [34]. The amount of penalty, which based on the second- (or higher-) order finite differences of the coefficients of adjacent B-splines, is easily controlled by non-negative smoothing (or penalty)
7
ACCEPTED MANUSCRIPT parameter. For example, the smoothing B-spline function S(x) for the dielectric function ε is the function that minimizes the penalized least squares error
2
(3)
i
PT
j
S ( x j ) k bi , 2
j
where the εj’s are the values of the non-smoothed dielectric function at the spectral points x j,
RI
is the non-negative smoothing parameter which controls the influence of the penalty and
SC
may be chosen by appropriate selection strategy, k is a positive integer, k is the k-order
NU
difference operator defined by means of 1bi bi bi 1 , 2bi 1 1bi bi 2bi 1 bi 2 ,
MA
k bi 1 k 1bi and the first sum is over all measured spectral points. A second-order
difference penalty has become a common choice. Other orders of the differences in the penalty can be also used to enforce shape constraints. The choice of the order might be very
ED
important in some applications for optimal extrapolation and the details can be found in Ref.
EP T
[34]. However, in our case, the order of the penalty should have no marked effect on the fit, especially if the number of knots is sufficiently large. Commonly, in statistical literature the Greek letter λ is used as a generic symbol for the smoothing parameter. However, since in
AC C
spectroscopic ellipsometry it is usually used to denote wavelength of light, we use the letter
to represent the smoothing parameter. P-spline smoothing typically uses a relatively large number of equally-spaced knots and this fact greatly simplifies the construction of penalties. P-splines have a very interesting property that makes them very attractive for the dielectric function representation, namely, the “power of the penalty” property: “The number of Bsplines can be (much) larger than the number of observations. The penalty makes the fitting procedure well-conditioned. This should be taken literally: even a thousand splines will fit ten observations without problems” [34]. Therefore, the elaborate problem of knot optimization 8
ACCEPTED MANUSCRIPT can be circumvented and replaced by a task of choosing the smoothing parameter for P-spline dielectric function representation with rather large number of knots (even equal to the number of spectral points) [42]. Generally speaking, the number of knots should be sufficiently large
MA
NU
SC
RI
PT
to ensure that ε can be represented with enough flexibility to get proper fit to the data.
Fig. 3. Penalized B-spline parameterizations of the imaginary part of the complex dielectric
ED
function for bulk aluminum with 188 equidistant knots: using no penalty ( = 0) on the spline coefficients (red solid line) and with a strong penalty ( = 0.1) to constrain the
EP T
coefficients of B-splines (green solid line). Restricted wavelength range was used to exclude
AC C
the Drude term contribution since it is not efficiently represented by polynomial functions.
It is perhaps easiest to illustrate P-splines by means of an example. Fig. 3 shows two B-spline parameterizations of the aluminum optical properties obtained from the same set of measured data. When setting a large number of equidistant knots, the B-spline parameterization without penalty tends to overfit the data and, although we get a perfect fit to our particular set of data, this fit is perfectly useless since it also reproduces an existing measurement noise. Not only a noise but, for instance, various random defects (point or line defects, vacancies, etc.) may be present in the samples and give their contribution to the dielectric function curve. The P-
9
ACCEPTED MANUSCRIPT spline representation with proper smoothing value will remove such unnecessary spectral features and enhance one of the main characteristics of any (dispersion) model, namely, prediction performance. In case of a strong penalty, the resulting penalized spline estimate of the imaginary part of the complex dielectric function has a nice smooth shape and approaches polynomial fit of some degree. Here an analogy can be drawn between P-spline
PT
representation and model selection approach [29,30]: the first term in Eq. (3) is a measure of
RI
the B-spline model’s quality of fit to the non-smoothed dielectric function while the second
SC
term penalizes a roughness of the resulted curve and, therefore, the smoothing parameter sets a compromise between two terms. Minimization of the penalized least squares error via Eq.
NU
(3) is used in conjunction with fitting of complete optical model, which consists of the material optical properties and thickness (or other topographic parameters) for each layer, to
MA
experimental data.
ED
The main idea of the P-spline parameterization can be summarized as follows: Represent the unknown dielectric function ε by means of ordinary B-splines [13] using a
EP T
large number of equally-spaced knots (close or even equal to the number of spectral points x j) to ensure that ε can be described with excessive flexibility (overfitted). The
AC C
result of this step, to some extent, reminds the direct data inversion when all existing measurement noise will be directly translated into the dielectric function curve although the B-spline formulation, unlike to direct inversion, is Kramers-Kronig consistent.
Introduce an additional penalty term (Eq. (3)) with certain smoothing parameter to tune the smoothness of the ε curve.
Perform usual data fitting with an appropriate optical model which include the smoothed dielectric function ε to minimize the ordinary least squares criterion
10
ACCEPTED MANUSCRIPT
min
1 W
m
y j 1 M
2
j ,M
yˆ j ,M ( S , | x j ) .
(4)
Here y j ,M and yˆ j ,M ( S , | x j ) denote the measured and calculated ellipsometric quantities, i.e., the ellipsometric angles Ψ and ∆, m is the number of the spectral points, the
PT
summation over subscript M corresponds to a particular y-quantity, W is the weighting
SC
RI
factor.
We complete this section with the following important remark. The main advantage of the
NU
penalization in the dielectric function representation by B-splines is that the smoothness of a spline curve is now controlled by single smoothing parameter and does not depend on the
MA
number and positions of knots. Since we now use a generous number of equidistant knots, an
ED
“optimized” knot positioning [13] is not a key determining factor.
2.2. Choosing the optimal degree of smoothing
EP T
Selection of 0 corresponds to unpenalized B-spline representation for the dielectric function and in case of a large number of knots will potentially lead to data overfitting (cf.
AC C
Fig. 2). The larger the value of the smoothing parameter and, therefore, the weight of the second term in Eq. (3), the smoother the result one can get, - a straight line in the limit
(or 1 if the smoothing parameter is defined on the interval [0, 1]). Therefore, we face the challenging problem of finding an optimal or, at least, acceptable value for the smoothing parameter to make a trade-off between data underfitting (or oversmoothing) and overfitting, i.e., we want to control the smoothness of the fitted curve to be neither “too smooth” nor “too rough”.
11
ACCEPTED MANUSCRIPT In statistics and data analysis when considering usual problem of scattered data smoothing the choice of α is typically made using some selection criteria such as (generalized) crossvalidation (e.g., Refs. [35–38]) or AIC (BIC), which require knowledge of the effective degree of freedom (DF). Eilers and Marx in Ref. [31], for instance, proposed to use the
PT
following form of AIC:
(5)
SC
RI
AIC( ) dev( y; b, a) 2 dim(b, ),
where dev(y; b, α) is the deviance, i.e., a goodness-of-fit measure for a model, and dim(b,α) is
NU
the effective dimension of the vector of parameters, also called the effective number of parameters or the effective DF. In other words, dim(b,α) can be interpreted by analogy with
MA
the number of fitted parameters in common AIC model selection as a peculiar measure of the model’s parametric complexity. However, we should bear in mind that, in general, DF(α) is
EP T
ED
not an integer in contrast to the true number of model parameters.
Our suggested application of P-splines in ellipsometric data fitting is different in the sense that we do not simply smooth the ellipsometric data but assign the smoothed optical function
AC C
ε, associated with particular smoothing parameter, to one of the layers in optical model and then perform data fitting process. Therefore, the smoothing parameter appears as an integral part of the model. A series of optical models with the dielectric function of the layer under test, represented by P-splines with different smoothing parameters, and other model variables (layer thicknesses, the optical functions of other layers, etc.) forms a set of candidate models which can be fitted to measured data one at a time to evaluate the quality of fit (e.g., MSE). Then one of the mentioned above “formal” model selection techniques can be used to select the “best” P-spline model with optimal amount of smoothing. However, rigorous estimation
12
ACCEPTED MANUSCRIPT of DF(α) in this case may be somewhat complex and some sort of empirical alternative for the smoothing parameter selection can be used, - for instance, by careful examination of the MSE behavior as a function of the smoothing parameter (or its derivative which, under smoothness assumption, can be approximated by finite differences). Such delicate ad hoc approach, of course, is not free from subjectivity and arbitrariness. Therefore, it might be
PT
possible that there is a certain range of reasonable smoothing parameter values and in practice
SC
RI
it would be acceptable to put the upper and lower bounds on α to alleviate the problem.
3. Real-data applications
NU
This section reports the results of the proposed P-spline approach to dielectric function parameterization. In the first application, we approximate the TiN dielectric function by P-
MA
splines with a set of the smoothing parameter values. In the second application, we applied the proposed methodology to represent the NiSi dielectric function. A commercially available
ED
tool for metal physical vapor deposition (PVD) was used to deposit approximately 200 Å of TiNx on standard 300 mm silicon (100) substrate with ~ 4000 Å of thermal silicon dioxide.
EP T
The PVD process uses titanium as a target material and working (Ar) and reactive (N 2 ) gases at certain partial pressures which exact values cannot be disclosed due to confidentiality
AC C
reasons. To prepare a ~ 250 Å-thick nickel-monosilicide (NiSi) film, an approximately 110 Å-thick Ni blanket film was deposited on 300 mm Si substrate by PVD process and then the sample was annealed by using a conventional two-step rapid thermal annealing with required silicidation processing temperature. Since we focus only on the P-spline representation of the TiN and NiSi optical properties, the details of all deposition processes (power settings, working pressure, deposition time, temperature and other process parameters) used to fabricate these metal films are beyond the scope of this paper. Smoothing of unpenalized Bspline dielectric functions for TiN and NiSi was done with relevant functions provided by the
13
ACCEPTED MANUSCRIPT Curve Fitting Toolbox 3.5 available in MATLAB® high-level language and interactive environment (version 8.4 (R2014b), The MathWorks, Inc., Natick, MA, U.S.A.). The measured ellipsometric data were fitted to the optical models containing the penalized (smoothed) TiN and NiSi optical properties by using CompleteEASE® software (version
PT
4.86) from J.A. Woollam Co., Inc. (Lincoln, NE, U.S.A.).
RI
3.1. Application 1: TiN
SC
In this first example, the ellipsometric measurements were performed in the spectral range of 1.25–6.45 eV (192–992 nm; 105 data points) at three angles of incidence (65°, 70°, and 75°)
NU
on an unpatterned test wafer with ~ 200 Å-thick TiN layer, deposited on thick SiO 2 film on Si substrate. The presence of the transparent oxide layer creates different optical paths for
MA
different angles of incidence hereby reducing the correlation between the TiN optical properties and thickness (see, for instance, Ref. [5], pp. 142–143). The optical properties of
ED
the TiN layer were represented by the degree 3 B-splines with 105 equidistant internal knots and tunable smoothing parameter. A surface roughness layer with variable thickness on the
EP T
TiN was also included in the model and its optical properties were obtained as a mixture of voids (50%) and underlying TiN film (50%) via the Bruggeman effective medium The silicon substrate and
AC C
approximation.
silicon dioxide optical constants, originally
published in Ref. [43], were taken from the J.A. Woollam Co., Inc. database of materials.
Resulting MSE values and log-values of the MSE average rate of change are plotted in Fig. 4 as functions of the smoothing parameter α. In our case, a choice of 0 corresponds to unpenalized dielectric function and, therefore, the lowest MSE value. The MSE curve itself represents monotonically increasing function. Note that the MSE curve is fairly smooth with a high-curvature point where the rapidly increasing MSE function of α transforms into slowly
14
ACCEPTED MANUSCRIPT increasing one. Although such point is clearly seen, it is still not obvious how to choose an appropriate value of the smoothing parameter α beyond which the quality of fit ceases to deteriorate significantly. However, on the MSE average rate of change curve we see more clearly an “elbow/knee point” (maximum curvature) at α = 0.02 and this value of the smoothing parameter can be taken as the best. Yet, we wish to use a more formal method to
PT
locate the elbow point of the curve and determine the optimal α. One of the simplest and
RI
intuitive, though somewhat naive, methods to identify the elbow point is known as “L-
SC
method” [44,45]. The L-method fits one straight line to the points on the left side of the graph and another straight line to the points on the right side and the intersection of these lines will
NU
approximately indicate the required smoothing parameter value. The line fitting should be performed for all possible pairs of lines with different sequences of points by minimizing the
MA
total root-mean-squared error, i.e., the combined MSE for the lines at the left and right sides of the graph. A more detailed description can be found in Refs. [44,45]. Using the L-method
ED
the point of maximum curvature in the MSE average rate of change curve was estimated to be α = 0.016. A similar value can be obtained by applying the Kneedle algorithm proposed in
EP T
Ref. [46]. Thus, in current example all these methods indicate that α ≈ 0.02 is the best value
AC C
of the smoothing parameter.
15
ACCEPTED MANUSCRIPT
Fig. 4. Thin-TiN-film example with 105 equally spaced knots: plot of MSE and log-values of the MSE average rate of change (difference quotient) against the smoothing parameter α. Vertical bar indicates points of maximum curvature and the abscissas of which identify the
ED
MA
NU
SC
RI
PT
optimal α.
EP T
Fig. 5. Refractive index (n) and extinction coefficient (k) plotted vs. wavelength for thin TiN film for different values of the smoothing parameter. To enhance visibility, the curves for
0 and
optimally penalized optical functions are shifted up vertically by 0.5
AC C
unpenalized
and 0.25 units, respectively.
As an example, in Fig. 5 three curves for the real, n, and imaginary, k, parts of the complex index of refraction obtained with different values of the smoothing parameter are shown. It can be clearly seen that the curve for α = 0 (unpenalized), which yields the lowest MSE value of 2.355, reveals noticeable “wiggly” behavior caused by very small knot spacing. At the same time, the curve with the highest value of parameter α (≃ 1.00), used in the analysis, produces higher misfit (MSE = 6.130) and, therefore, misses some essential spectral features. 16
ACCEPTED MANUSCRIPT The curve generated with selected by L-method optimal value of α ≈ 0.02 displays the absence of pronounced “wiggliness” and yields acceptable MSE value of 3.849. It is perhaps worth to note that this MSE value is comparatively close to the MSE value of 3.817 for the same application example obtained in our previous study [29] without penalization but by
PT
using optimization of the number of knots instead.
RI
In Fig. 5 the n and k curves, even after optimal smoothing with α = 0.02, still reveal a
SC
discernible local artifact at ~ 450 nm, possibly a systematic error, although all random errors (a.k.a. measurement noise) were wiped out completely. The resulting degree of smoothness
NU
simply depends on a value of the smoothing parameter and “amplitude” of error. Thereby, the fitting with increased parameter value of 0.1 practically removes this slight non-random flaw
MA
although the MSE value increases a little from 3.849 to 4.277 due to the fact that one does not fit the local artefact in the measurement spectra anymore. This example illustrates an obvious
ED
condition to have correct and optimal P-spline representation for the dielectric function, namely, the systematic errors (biases) should be identified and eliminated after careful
EP T
inspection of the experimental data. However, unlike random errors, the systematic errors are quite difficult to spot: “…there is a well-developed mathematical theory to deal with random
AC C
errors while there is seemingly nothing to say about systematic errors other than “Be clever and catch them.” [47]. We believe that the best way of avoiding systematic biases by practitioners of ellipsometry is to use various verification techniques, for instance, via multisample, multi-ellipsometer examination.
3.2. Application 2: NiSi In the second example, we investigate a thin metal layer of NiSi on a Si substrate. A NiSi single film was formed by silicidation at elevated temperatures of a thin Ni layer deposited
17
ACCEPTED MANUSCRIPT onto a Si substrate using commercially available tool. The ellipsometric data were acquired in the spectral range of 1.24–6.17 eV (200–1000 nm; 505 data points) at the angle of incidence of 65°. The optical properties of the NiSi layer were parametrized by means of the degree 3 B-splines with 504 equidistant internal knots and adjustable smoothing parameter. Good starting optical properties for the layer to set up the B-spline model have been obtained using
PT
various oscillator models using the global fit option to avoid false solutions. The silicon
RI
substrate optical constants were also taken from the J.A. Woollam Co., Inc. database of
EP T
ED
MA
NU
SC
materials.
AC C
Fig. 6. Thin-NiSi-film example with 504 equally spaced knots: plot of MSE and log-values of the MSE average rate of change (difference quotient) against the smoothing parameter α. Vertical bar indicates points of maximum curvature and the abscissas of which identify the optimal α.
In this example, the resulting MSE values and log-values of the MSE average rate of change are plotted in Fig. 6 as functions of the smoothing parameter α. Examining obtained results, one can find the point of maximum curvature in the MSE average rate of change curve at α =
18
ACCEPTED MANUSCRIPT 0.0019 (from the L-method). Then we can accept this value as an appropriate choice of the smoothing parameter. As above, the corresponding real and imaginary parts of the complex
MA
NU
SC
RI
PT
index of refraction for different values of the smoothing parameter are illustrated in Fig. 7.
ED
Fig. 7. Refractive index (n) and extinction coefficient (k) plotted vs. wavelength for thin NiSi film for different values of the smoothing parameter. To enhance visibility, the curves for
0 and
optimally penalized optical functions are shifted up vertically by 0.5
EP T
unpenalized
and 0.25 units, respectively. The insets show the enlarged views of the n and k spectra in the
AC C
spectral range of 300–450 nm.
The curve corresponding to the unpenalized dielectric function of NiSi generates, as expected, the lowest MSE value of 0.866 and visible “wiggliness” can be clearly seen between 300 and 450 nm indicating that the measurement noise is not filtered out (see insets in Fig. 7). The curve obtained with the highest used value of parameter α = 1.00 leads to pretty high MSE value of 14.52, thereby denoting a significant misfit. Finally, use of optimally smoothed dielectric function of NiSi with α ≈ 0.002 results in reasonable MSE value of 1.355. It is worth mentioning that the NiSi dielectric function obtained from 19
ACCEPTED MANUSCRIPT choosing the optimal number of knots by the information criteria (AIC, BIC), as it has earlier been suggested in Refs. [29,30], yields quite close MSE value of 1.519.
Summarizing, we can say that penalized B-splines provide a simple and yet effective method for representing the dielectric functions of materials. In fact, P-splines eliminate an ill-posed
PT
problem of choosing the optimal number and placement of B-spline knots, thereby
SC
RI
simplifying the ellipsometric data analysis process.
4. Conclusions
NU
In this paper, we outline penalized B-spline (P-spline) formulation for an arbitrary dielectric function ε modeling in spectroscopic ellipsometry data analysis. The main idea is to first use
MA
a generous number of equally-spaced knots (up to the number of spectral points) to describe even highly complicated structure in the the imaginary part ε2 of the dielectric function and
ED
then provide a certain penalty on the coefficients of adjacent B-splines to tune the smoothness of the ε curve. Current implementation of this approach requires a small amount of user
EP T
intervention to select the right value of penalization for each application or, at least, indicate a certain range of reasonable smoothing parameter values. Two real-data applications have
AC C
been provided to evaluate the practical performance and effectiveness of the proposed penalized formulation for dielectric function representation. A comparison of obtained results with the findings of our previous studies without penalization [29,30] demonstrates a good agreement between the proposed method and ordinary B-spline representation in an equallyspaced fashion with optimal number of knots.
Finally, we would like to remark that only the simplest case of P-spline formulation with a global smoothing parameter, when the degree of smoothness of the dielectric function stays
20
ACCEPTED MANUSCRIPT about the same across the entire spectral region, has been considered in this work. If the dielectric function curve has a complicated shape and contains abrupt as well as smooth curvature changes, then a single value of α might not be efficient enough for an entire spectral region. Then we should apply, for example, more complicated approach called highly adaptive smoothing [35] wherein α varies as a function of some independent parameter (like
PT
a gradually changing smoothness or derivative information). Another approach, which still
RI
requires more investigation, is to introduce certain data-dependent weights on the differences
SC
in the penalty term in Eq. (3) (similarly to introduction of the adaptive penalty function for
NU
the smoothing splines [48]).
Acknowledgements
MA
The author wish to acknowledge Dr. Martin Weisheit (GLOBALFOUNDRIES, Dresden,
ED
Germany) for providing TiN ellipsometric data used in present work.
References
EP T
[1] H.G. Tompkins, W.A. McGahan, Spectroscopic Ellipsometry and Reflectometry: A User's Guide, John Wiley & Sons Inc., New York NY, 1999.
AC C
[2] H. Fujiwara, Spectroscopic Ellipsometry: Principles and Applications, John Wiley & Sons Ltd., Chichester, West Sussex, U.K., 2007. [3] R.W. Collins, A.S. Ferlauto, Optical physics of materials, in: H.G. Tompkins, E.A. Irene (Eds.), Handbook of Ellipsometry, William Andrew Publishing/Noyes, Norwich NY 2005, pp. 93–235. [4] J. Humlicek, Data analysis for nanomaterials: Effective medium approximation, its limits and implementations, in: M. Losurdo, K. Hingerl (Eds.), Ellipsometry at the Nanoscale, Springer-Verlag, Berlin, Heidelberg 2013, pp. 145–178.
21
ACCEPTED MANUSCRIPT [5] H.G. Tompkins, J.N. Hilfiker, Spectroscopic Ellipsometry: Practical Application to Thin Film Characterization, Momentum Press, New York NY, 2016. [6] C.C. Kim, J.W. Garland, H. Abad, P.M. Raccah, Modeling the optical dielectric function
of
semiconductors:
extension
of
the
critical-point
parabolic-band
approximation, Phys. Rev. B 45 (1992) 11749–11767.
PT
[7] C.C. Kim, J.W. Garland, P.M. Raccah, Modeling the optical dielectric function of the
RI
alloy system AlxGa1−xAs, Phys. Rev. B 47 (1993) 1876–1888.
SC
[8] J.-Th. Zettler, Th. Trepk, L. Spanos, Y.-Z. Hu, W. Richter, High precision UVvisible-near-IR Stokes vector spectroscopy, Thin Solid Films 234 (1993) 402–407.
NU
[9] M. Zorn, T. Trepk, J.-T. Zettler, B. Junno, C. Meyne, K. Knorr, T. Wethkamp, M. Klein, M. Miller, W. Richter, L. Samuelson, Temperature dependence of the InP(001)
MA
bulk and surface dielectric function, Appl. Phys. A 65 (1997) 333–339. [10] B. Johs, C.M. Herzinger, J.H. Dinan, A. Cornfeld, J.D. Benson, Development of a
ED
parametric optical constant model for Hg1−xCdxTe for control of composition by
137–142.
EP T
spectroscopic ellipsometry during MBE growth, Thin Solid Films 313–314 (1998)
[11] M. Gilliot, A. Hadjadj, M. Stchakovsky, Spectroscopic ellipsometry data inversion
AC C
using constrained splines and application to characterization of ZnO with various morphologies”, Appl. Surf. Sci., 421 (2016) 453–459. [12] M. Gilliot, Inversion of ellipsometry data using constrained spline analysis, Appl. Opt. 56 (2017) 1173–1182. [13] B. Johs, J.S. Hale, Dielectric function representation by B-splines, Phys. Status Solidi A 205 (2008) 715–719. [14] C. de Boor, A Practical Guide to Splines, Revised Edition, Springer-Verlag, New York NY, 2001.
22
ACCEPTED MANUSCRIPT [15] P. Dierckx, Curve and Surface Fitting with Splines, Oxford University Press, Inc., New York NY, 1993. [16] J.W. Weber, T.A.R. Hansen, M.C.M. van de Sanden, R. Engeln, B-spline parametrization of the dielectric function applied to spectroscopic ellipsometry on amorphous carbon, J. Appl. Phys. 106 (2009) 123503.
PT
[17] J.W. Weber, V.E. Calado, M.C.M. van de Sanden, Optical constants of graphene
RI
measured by spectroscopic ellipsometry, Appl. Phys. Lett. 97 (2010) 091904.
SC
[18] S.G. Choi, J. Zúñiga-Pérez, V. Muñoz-Sanjosé, A.G. Norman, C.L. Perkins, D.H. Levi, Complex dielectric function and refractive index spectra of epitaxial CdO thin
NU
film grown on r-plane sapphire from 0.74 to 6.45 eV, J. Vac. Sci. Technol. B 28 (2010) 1120–1124.
MA
[19] H.T. Beyene, J.W. Weber, M.A. Verheijen, M.C.M. van de Sanden, M. Creatore, Real time in situ spectroscopic ellipsometry of the growth and plasmonic properties of Au
ED
nanoparticles on SiO 2 , Nano Res. 5 (2012) 513–520. [20] S.G. Choi, H.Y. Zhao, C. Persson, C.L. Perkins, A.L. Donohue, B. To, A.G. Norman, Li,
I.L.
Repins,
Dielectric function spectra and critical-point energies of
EP T
J.
Cu2 ZnSnSe4 from 0.5 to 9.0 eV, J. Appl. Phys. 111 (2012) 033506.
AC C
[21] A. Jimenez, D. Lepage, J. Beauvais, J.J. Dubowski, Study of surface morphology and refractive index of dielectric and metallic films used for the fabrication of monolithically
integrated
surface
plasmon
resonance
biosensing
devices,
Microelectron. Eng. 93 (2012) 91–94. [22] E. Agocs, B. Fodor, B. Pollakowski, B. Beckhoff, A. Nutsch, M. Jank, P. Petrik, Approaches to calculate the dielectric function of ZnO around the band gap, Thin Solid Films 571 (2014) 684–688.
23
ACCEPTED MANUSCRIPT [23] L. Kőrösi, A. Scarpellini, P. Petrik, S. Papp, I. Dékány, Sol–gel synthesis of nanostructured indium tin oxide with controlled morphology and porosity, Appl. Surf. Sci. 320 (2014) 725–731. [24] L.S. Abdallah, S. Zollner, C. Lavoie, A. Ozcan, M. Raymond, Compositional dependence of the optical conductivity of Ni1
− x Ptx
alloys (0 < x < 0.25) determined
PT
by spectroscopic ellipsometry, Thin Solid Films 571 (2014) 484–489.
RI
[25] V. Zviagin, P. Richter, T. Böntgen, M. Lorenz, M. Ziese, D.R.T. Zahn, G. Salvan, M.
SC
Grundmann, R. Schmidt-Grund, Comparative study of optical and magneto‐optical properties of normal, disordered, and inverse spinel‐type oxides, Phys. Status Solidi B
NU
253 (2016) 429–436.
[26] J.N. Hilfiker, J.S. Hale, C.M. Herzinger, T. Tiwald, N. Hong, S. Schöche, H. Arwin,
MA
Estimating depolarization with the Jones matrix quality factor, Appl. Surf. Sci. (2017) 494–499.
421
ED
[27] N. Hong, R.A. Synowicki, J.N. Hilfiker, Mueller matrix characterization of flexible plastic substrates, Appl. Surf. Sci. 421 (2017) 518–528.
EP T
[28] M. Di, E. Bersch, A.C. Diebold, S. Consiglio, R.D. Clark, G.J. Leusink, T. Kaack, Comparison of methods to determine bandgaps of ultrathin HfO 2 films using
AC C
spectroscopic ellipsometry, J. Vac. Sci. Technol. A 29 (2011) 041001. [29] D.V. Likhachev, Selecting the right number of knots for B-spline parameterization of the dielectric functions in spectroscopic ellipsometry data analysis, Thin Solid Films 636 (2017) 519–526. [30] D.V. Likhachev, B-spline parameterization of the dielectric function and information criteria: The craft of non-overfitting, in: B. Bodermann, K. Frenner, R.M. Silver (Eds.), Modeling Aspects in Optical Metrology VI, Munich, Germany, June 25-29, 2017, SPIE Proc. 10330 2017, 103300B.
24
ACCEPTED MANUSCRIPT [31] P.H.C. Eilers, B.D. Marx, Flexible smoothing with B-splines and penalties, Statist. Sci. 11 (1996) 89–102. [32] B.D. Marx, P.H.C. Eilers, Generalized linear regression on sampled signals and curves: a P-spline approach, Technometrics 41 (1999) 1–13. [33] P.H.C. Eilers, B.D. Marx, Splines, knots, and penalties, WIREs Comp. Stat. 2 (2010)
PT
637–653.
RI
[34] P.H.C. Eilers, B.D. Marx, M. Durbán, Twenty years of P-splines, SORT-Stat. Oper.
SC
Res. Trans. 39 (2015) 149–186.
[35] D. Ruppert, M.P. Wand, R.J. Carroll, Semiparametric Regression, Cambridge
NU
University Press, Cambridge, U.K., 2003.
[36] D. Ruppert, Statistics and Finance: An Introduction, Springer-Verlag, New York NY,
MA
2004.
[37] L. Fahrmeir, T. Kneib, S. Lang, B. Marx, Regression: Methods, Models and
ED
Applications, Springer-Verlag, Berlin, Heidelberg, 2013. [38] S.N. Wood, Generalized Additive Models: An Introduction with R, 2nd Ed., CRC
EP T
Press, Boca Raton FL, 2017.
[39] A.N. Tikhonov, V.Y. Arsenin, Solutions of Ill-Posed Problems, V.H. Winston &
AC C
Sons, Washington D.C., 1977. [40] A.N. Tikhonov, A.S. Leonov, A.G. Yagola, Nonlinear Ill-Posed Problems, Vols. 1 & 2, Chapman & Hall, London, U.K., 1998. [41] L. Tenorio, Statistical regularization of inverse problems, SIAM Rev. 43 (2001) 347– 366. [42] D.V. Likhachev, Dielectric function parameterization by penalized splines, in: B. Bodermann, K. Frenner, R.M. Silver (Eds.), Modeling Aspects in Optical Metrology VI, Munich, Germany, June 25-29, 2017, SPIE Proc. 10330 2017, 103301J.
25
ACCEPTED MANUSCRIPT [43] C.M. Herzinger, B. Johs, W.A. McGahan, J.A. Woollam, W. Paulson, Ellipsometric determination of optical constants for silicon and thermally grown silicon dioxide via a multi-sample, multi-wavelength, multi-angle investigation, J. Appl. Phys. 83 (1998) 3323–3336. [44] S. Salvador, P. Chan, Determining the number of clusters/segments in hierarchical
PT
clustering/segmentation algorithms, in: Proceedings of the 16th IEEE International
SC
U.S.A., November 15-17, 2004, pp. 576–584.
RI
Conference on Tools with Artificial Intelligence (ICTAI 2004), Boca Raton FL,
[45] S. Salvador, P. Chan, Learning states and rules for detecting anomalies in time series,
NU
Appl. Intell. 23 (2005) 241–255.
[46] V. Satopää, J. Albrecht, D. Irwin, B. Raghavan, Finding a "kneedle" in a haystack:
Conference on Distributed
MA
Detecting knee points in system behavior, in: Proceedings of the 31st International Computing Systems Workshops (ICDCSW 2011),
ED
Minneapolis MN, U.S.A. June 20-24, 2011, pp. 166–171. [47] J. Bechhoefer, Curve fits in the presence of random and systematic error, Am. J. Phys.
EP T
68 (2000) 424–429.
[48] X. Wang, P. Du, J Shen., Smoothing splines with varying smoothing parameter,
AC C
Biometrika 100 (2013) 955–970.
26
ACCEPTED MANUSCRIPT Highlights D.V. Likhachev, Spectroscopic ellipsometry data analysis using penalized splines representation for the dielectric function (Ref. No.: TSF-D-18-00749)
P-spline approach for representing the dielectric functions of materials is proposed
It eliminates an ill-posed problem of knot selection in B-spline parameterization
It uses a generous number of equally-spaced knots to ensure enough flexibility
A roughness penalty on the coefficients of adjacent B-splines is applied
Effectiveness of the proposed approach is illustrated by real-data examples
AC C
EP T
ED
MA
NU
SC
RI
PT
27
ACCEPTED MANUSCRIPT Figure Captions D.V. Likhachev, Spectroscopic ellipsometry data analysis using penalized splines representation for the dielectric function (Ref. No.: TSF-D-18-00749) Fig. 1. Example of shifted B-spline basis functions of degree p = 3 over the interval [0, 8] with equally-spaced internal knots (uniform B-splines) and different weights (black dashed
PT
and blue solid lines) and their linear combination (red solid line). The knots are indicated by black dots. Such equidistant knot arrangement can be justified if the sample size, i.e., the
RI
number of data points, is not small comparing to the size of the knot vector and the data
SC
points are distributed uniformly in the considered spectral range.
NU
Fig. 2. Influence of the number of knots on the B-spline representations of the TiN complex index of refraction (the knots are equally spaced in the wavelength range of 190–1000 nm).
MA
One can clearly see that the B-spline model with 35 knots is too flexible and overparameterized the experimental data resulting in a few numerical artifacts at ~200, 245, 320,
“wiggly” overfitting behavior.
ED
440, 730 nm due to small knot spacing. Reduction of the number of knots suppresses this
EP T
Fig. 3. Penalized B-spline parameterizations of the imaginary part of the complex dielectric function for bulk aluminum with 188 equidistant knots: using no penalty ( = 0) on the
AC C
spline coefficients (red solid line) and with a strong penalty ( = 0.1) to constrain the coefficients of B-splines (green solid line). Restricted wavelength range was used to exclude the Drude term contribution since it is not efficiently represented by polynomial functions. Fig. 4. Thin-TiN-film example with 105 equally spaced knots: plot of MSE and log-values of the MSE average rate of change (difference quotient) against the smoothing parameter α. Vertical bar indicates points of maximum curvature and the abscissas of which identify the optimal α. Fig. 5. Refractive index (n) and extinction coefficient (k) plotted vs. wavelength for thin TiN film for different values of the smoothing parameter. To enhance visibility, the curves for 28
ACCEPTED MANUSCRIPT unpenalized (α = 0) and optimally penalized optical functions are shifted up vertically by 0.5 and 0.25 units, respectively. Fig. 6. Thin-NiSi-film example with 504 equally spaced knots: plot of MSE and log-values of the MSE average rate of change (difference quotient) against the smoothing parameter α. Vertical bar indicates points of maximum curvature and the abscissas of which identify the
PT
optimal α.
RI
Fig. 7. Refractive index (n) and extinction coefficient (k) plotted vs. wavelength for thin NiSi
SC
film for different values of the smoothing parameter. To enhance visibility, the curves for unpenalized (α = 0) and optimally penalized optical functions are shifted up vertically by 0.5
NU
and 0.25 units, respectively. The insets show the enlarged views of the n and k spectra in the
AC C
EP T
ED
MA
spectral range of 300–450 nm.
29