## Abstract

The optical constants of ruthenium in the spectral range 8 nm – 23.75 nm are determined with their corresponding uncertainties from the reflectance of a sputtered ruthenium thin film, measured using monochromatized synchrotron radiation. This work emphasizes the correlation between structure modelling and the determined optical parameters in a robust inverse-problem solving strategy. Complementary X-ray Reflectivity (XRR) measurements are coupled with Markov chain Monte Carlo (MCMC) based Bayesian inferences and a quasi-model-independent method to create a model factoring the sample’s oxidation, contamination, and interfacial imperfections. The robustness of the modelling scheme against contamination and oxidation is tested and verified by measurements after hydrogen-radical cleaning of the sample’s surface.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. Introduction

Extreme Ultraviolet Lithography (EUVL) systems are now used in high-volume manufacturing at the 7-nm-node. Its capabilities need to be further improved and current challenges mitigated to support future technology nodes. The characteristics of Extreme Ultraviolet (EUV) radiation, with materials being highly absorbing requires the utilization of reflective optics and creates specific challenges for EUVL. To further enhance the capabilities of EUVL and to investigate the possible improvement through material optimization, simulations are carried out [1]. The optical constants of the materials involved, n(λ) = 1-δ(λ) + iβ(λ) are a key parameter for the simulations. δ and β relate to the phase velocity and to the absorption of electromagnetic radiation with a wavelength λ, respectively. A hindrance for EUVL development is that the optical constants for many materials in the EUV spectral range are not adequately determined. In addition to the lack of calculated uncertainties for the optical constants of most materials, inconsistencies are reported between different references for some materials [2,3].

A technologically relevant material is ruthenium. It is used in many critical optical components of EUVL systems [4]. Due to its high chemical stability and oxidation resistance, high reflectivity in the EUV spectral range and high etch selectivity to the materials of mask absorber layers (e.g.: TaN), Ru thin films are used as capping layers of the Mo/Si multilayer mirrors (MLMs) of the masks, to extend the lifetime of the MLM by protecting it during the pattern repair and cleaning processes [4,5]. In addition to the use of Ru capping layers in EUVL masks, Ru is also used as capping or reflective layer (grazing incidence) throughout the optical path of EUVL machines and due to its high thermal emissivity Ru coated pellicles are reported to have better performance characteristics [4]. Ru-based alloys are promising materials to mitigate the induced 3D mask imaging effects that are crucial for the development of next-generation EUVL technology. This could be achieved by using thinner EUV mask absorbers than the currently used tantalum-based absorbers [6]. However, the study presented here could also be useful for the development of EUV astronomy optics, since the spectral range studied covers the emission lines of some ionized states of common elements, that are considered in solar weather research and are of interest for the observation of planetary magnetospheres [7,8].

The two most widely used methods to determine optical constants in the EUV spectral range are angle-dependent reflectometry and transmission-mode measurements. The latter method usually requires free-standing thin films, of well-defined thickness, from which the absorption part β can be determined and the dispersive part δ is then calculated using Kramers–Kronig analysis [9–11]. This requires a comprehensive compilation of optical data, spanning the X-ray range to the infrared which is simply missing for many materials particularly alloys. However, transmission-mode measurements are relatively insensitive to the surface state of the sample, and effects of contamination can be reliably accounted for [11].

The reflectometry method allows non-destructive characterization of the thin films and substrates with minimal preparation requirements, enabling the simultaneous determination of δ and β [12]. To circumvent the limitations associated with determining optical constants using transmission-mode measurements, and since many optical components of EUVL systems are rather stratified systems, angle-dependent reflectometry was selected here. In this work, the optical constants of ruthenium in the spectral range 8 nm – 23.75 nm (52.2 eV – 154.9 eV) are determined from angle-dependent EUV reflectometry (EUVR) measurements conducted in the radiometry laboratory of the Physikalisch-Technische Bundesanstalt (PTB) at the electron storage ring of Berliner Elektronenspeicherring-Gesellschaft für Synchrotronstrahlung (BESSY II).

Another advantage for reflectometry is that it allows the characterization of thin films with similar settings of their intended application. Concurrently, this advantage accompanies a challenge, since the optical constants of a thin film in certain spectral ranges are correlated with its morphology and thickness [13]. Even the deposition settings could alter the optical characteristics of a thin film, although this is not expected to be the case of ruthenium and other platinum-group metals, due to their chemical inertness [14].

Reflectometry is highly sensitive to interfacial imperfections, surface contamination and oxidation. Additionally, the determination of optical constants is achieved indirectly via an inverse-problem solving strategy [15]. The selection of a proper model of a sample in a reflectometry inverse-problem is often not straightforward. Even single thin films on substrates often have to be modeled as multilayer systems due to oxidation, contamination or interdiffusion processes. This requires a proper initialization of the reflectometry inverse-problem, especially regarding the sample’s geometrical structure.

To aide in the determination of the optical constants, complementary X-ray Reflectivity (XRR) data was used here to construct a geometrical model of the sample’s structure. Given the high sensitivity of XRR for the structural characteristics and interfacial imperfections of thin films, a reliable model of the sample can be derived for the optical constants’ determination in the EUV spectral range.

In this work, the selection of the assumed model of the sample is supported by time-frequency analysis, a quasi-model independent approach for XRR data analysis. That is crucial since a deficient or flawed geometrical model ultimately influences the determined optical constants. We put special emphasis on the determination of the parameter uncertainties using Markov chain Monte Carlo (MCMC)-based Bayesian inferences, with focus on the correlations between the determined optical constants and the presumed sample’s structure.

Our strategy to determine optical constants using reflectometry is capable of factoring the presence of contamination and oxidation. One direct way to avoid contamination and oxidation of thin films when determining optical constants would be to conduct the measurements in-situ (i.e., to deposit the films and conduct measurements without exposing the sample to air) [16]. However, EUVR and deposition of thin films both require sophisticated and complex set-ups which cannot easily be combined. Also, molecular contamination which is critical here may also happen in the experimental vacuum environment [17]. Therefore, the aim of this work is to present a comprehensive approach that is robust against the effects of contamination and oxidation.

The robustness of the presented approach is substantiated by a comparison between the results of analyzing datasets collected from a sample considering significantly different surface conditions, where the data was collected before and after hydrogen radical cleaning. These results are also highly relevant for the EUVL technology where optical elements with ruthenium capping are routinely currently used in EUVL systems, and are exposed to similar cleaning processes to mitigate the effects of tin debris due to the laser-produced plasma [18].

## 2. Experimental setup

#### 2.1 Sample preparation

The ruthenium thin film was deposited using DC magnetron sputtering onto a 300 mm Si wafer with a surface roughness of less than 0.1 nm root-mean-square (RMS). After evacuating the sputtering chamber (p_{0} < 10^{−7} mbar), argon with high purity (6N) was used as a sputtering gas. The native oxide of the substrate was not removed before deposition of the ruthenium film. Afterwards, the coated wafer was cleaved into smaller square coupons with a side length of 25 mm. Using a laboratory-based X-ray diffractometer (D8 Discover by Bruker) in the laboratories of OptiX fab GmbH (Jena), an XRR profile was measured at Cu-Kα from a witness sample. By analyzing the peaks of the XRR profile, the critical angle was calculated and the approximated mass density from the electron density of the ruthenium film was found to be close to that tabulated of the bulk of 12.1 g/ cm^{3} [19–21].

#### 2.2 EUV reflectivity measurements

The EUVR data was collected at the soft X-ray radiometry beamline (SX700) in the EUV radiometry laboratory of PTB at the electron storage ring facility BESSY II [22–24]. The main characteristics of the SX700 beamline are shown in Table 1. A vacuum tank housing a reflectometer setup serves as an experimental endstation for the beamline to allow the necessary handling of the samples.

The reflectometer’s sample holder has 6-degrees of freedom placing the sample in any orientation to the incident beam, particularly incidence angles between 1.5° and 90° relative to the surface’s normal in S- and P- polarization directions with respect to the linear polarized incoming radiation can be realized. The radiation is detected with GaAsP photodiodes. The machineries are lubricant-free, and the vacuum tank is housed in a clean room environment to reduce the potential contamination of the probed optical surfaces.

Additionally, the SX700 reflectometer is equipped with a load-lock chamber to allow in vacuo transfer of samples to an auxiliary process chamber where hydrogen-radical cleaning is feasible [25]. In this chamber’s setup, hydrogen radicals are produced by thermal dissociation of hydrogen molecules using a heated tungsten filament [18].

#### 2.3 XRR Measurement

An additional XRR measurement at a photon energy of 1750eV was performed at another beamline, the plane grating monochromator (PGM) beamline of PTB, which provides soft X-ray radiation of high spectral purity [26]. Few of the main characteristics of the PGM beamline are shown in Table 2.

The experimental end station used is an ultra-high vacuum chamber equipped with a multi-axis manipulator that also allows XRR experiments in θ-2θ geometry over a wide angular range [27].

## 3. Bayesian framework for optical constants determination

The determination of optical constants from reflectance data is possible by the optimization of an inverse problem. This can be achieved by iteratively minimizing an objective function that compares the measured dataset against a simulated one, based on a presupposed model. In the optimization process, the parameters of the model are refined iteratively and are assumed to be descriptive of the probed sample once the stopping criterion of the optimization is met.

Preliminary results for the optical constants of ruthenium, in the spectral range 10 nm – 20 nm, based on the aforementioned optimization scheme but without uncertainties have been presented before [28].

Optimizing reflectivity data using the previously addressed formalism does not provide the uncertainties and the correlations of the determined parameters. The anticipated vast solution space is not guaranteed to be fully sampled. Moreover, even if multiple optimization trials were performed, potential multimodalities are hardly detected.

A framework mitigating the previously addressed complications and enabling an assertive approach for the determination of optical constants, can be realized combining MCMC algorithms and Bayesian inferences.

Applying Bayesian inferences to extract information from reflectivity data is known, and has been particularly useful for selecting a model to solve the relevant reflectometry inverse-problem [29]. Bayesian inferences have also been used for the design of MLM [30]. In this work, the model is selected using time-frequency analysis of XRR data. The gain from applying a Bayesian based framework here is the extended investigation of the posterior distribution of the model parameters. This investigation can easily indicate if there are problems with the presupposed selected model. It could also qualitatively indicate convergence and it is very informative about the occurrence of multimodalities and possible correlations of the parameters.

More details on the applicability of Bayesian inferences to reflectometry data, in addition to MCMC algorithms are discussed elsewhere [29–34].

Earlier, the determination of optical constants with their corresponding uncertainties from magnetron sputtered boron carbide thin films using MCMC, in the spectral range 40 nm – 80 nm has been conducted by treating each wavelength separately [12]. The attempt here regards a global simulation approach, presumably, this approach would yield higher accuracy. A global simulation approach is also expected to be highly instructional considering the presupposed model, because then the tolerance for the model’s characteristics is firmer. When each wavelength scan is taken separately, the model’s deficiencies can be compensated on the expense of biasing different degrees of freedom. For example, a lower density of the film could be mixed with a slightly higher thickness value, but such an issue is less likely to occur when treating numerous wavelength scans simultaneously.

MCMC algorithms allow exploring the entire parameter space and Bayes’ inferences aide in correlating the probabilities updated from the MCMC sampled space to the determined parameters [30,35]. For this work, the MCMC sampler proposals are generated using the “Stretch Move” ensemble method proposed by Goodman and Weare [35,36].

From Bayes’ theorem, the probability distribution of the parameters array **p** to represent the characteristics of the sample given the measured data (the posterior) is given by [30]:

where the set of target parameters to be determined from the measured reflectivity data are presented in the following array:

where d_{m} and s_{m} denote the mth layer thickness and the mth interfacial roughness, respectively, and $\Pr ({{{\mathbf R}_{measured}}\textrm{|}{\mathbf p}} )$ is the probability distribution of obtaining the measured reflectivity data given the parameters array **p**. $\Pr ({\mathbf p} )$ represents our (prior) knowledge about the sample. The evidence $\Pr ({{{\mathbf R}_{measured}}} )$ is considered here as a normalization constant and hence [Eq. (1)] can be reformulated as [30]:

Assuming that the measurement errors have a Gaussian distribution and are not correlated, the likelihood function comparing between the measured reflectivity and the simulated, supposing that the sample characteristics are represented by the parameters array **p** is given by [37]:

where ${\sigma _n}$ denotes the error in the nth measured data point and it is given here as a summation [37]:

*a*and

*b*are supposed to partially account for the model’s deficiencies in describing the sample’s real structure.

*a*and

*b*are included in the parameters array

**p**and to be sampled in the MCMC-based evaluation. In [Eq. (5)] we assume that the modelling error is normally distributed and has the same properties as the known measurement error.

To restrict the posterior from updating with probabilities from an overwhelmingly large range we define the prior as a function for the set of parameters in **p** [35], so the likelihood search space is restricted to maintain higher efficiency for the MCMC algorithm:

The desired result in the MCMC-based Bayesian inferences approach is the stationary distribution of the approximated posterior regarding all the parameters in the array **p** [35]. A major drawback of this method then is the required computational effort since it increases with the number of parameters. The required number of MCMC proposals to reach the anticipated target distribution is usually considerable, especially with the large number of sampled parameters.

## 4. Hybrid model parametrization

The proper determination of optical constants regarding thin films from reflectivity data is only attainable with the simultaneous determination of the geometrical characteristics of the probed sample precisely. This poses a conspicuous obstacle since EUV radiation is highly attenuated by almost all materials. This limitation is exacerbated by the presence of inevitable contamination and possible surface oxidation concurrent with interfacial roughness and interdiffusion. To aide in the parameterization of the EUVR inverse problem, XRR measurements were opted for their non-destructiveness and high penetration depth enabling high sensitivity even for buried interfaces [38,39].

Deducing information regarding a multilayer from XRR data is not straightforward. Usually, extracting information is demonstrated with a time-consuming model-dependent approach utilizing trial and error optimization scheme. This becomes more laborious if none of the sample’s characteristics is properly known a priori. Although model selection in the realm of reflectometry can be demonstrated via Bayesian inferences [31], it is rarely done this way. On many occasions, researchers employ Occam's razor in selecting a model to simulate reflectivity data [40]. The interpretation for experimental data simulation is to model the data with the smallest possible number of parameters [41]. Because a “smallest possible number” is hard to define, model- independent methods for directly inferring information from reflectivity data are required. Such an approach is also interesting for avoiding overfitting and underfitting by defining the number of layers needed. Methods like Fourier analysis of XRR data can retrieve the density profile of a multilayer system, aiding in setting a preliminary optimization model or to verify the results from one [20].

The Fourier transform (FT) of XRR profiles is a reliable complementary approach to aid in the model-independent analysis. The peaks of a FT spectrum of an XRR profile indicate the depths of the interfaces present in the sample, from which the layers thicknesses can be estimated [40], As adapted from Smigiel et al. [42], a limitation of FT is the loss of information regarding the temporal evolution in the case of a signal progressing in a time domain, which translates to the loss of information regarding the stratification order for XRR profiles with the intensity modulations progressing with increasing grazing incidence angle. Also, the FT resolution is vulnerable to the windowing procedure of the experimental data [20]. Time-frequency representations (TFRs) allow investigating signals in both time and frequency domains simultaneously, and can be used to analyze XRR profiles [42–44] . The angle (or the wave vector) and the interfacial depth (which translates to the layers thickness) replace the pair of conjugate variables of time and frequency, respectively [42–44]. Smigiel et al. [42] demonstrated the use of a TFR for a study of XRR profiles, via the short-time Fourier transform (STFT), enabling a quasi-model-independent approach to infer valuable information from phaseless XRR data. Similarly, Smigiel et al. [43] further elaborated the application of continuous wavelet transform (CWT) on XRR data. Where the properties of CWT produce a scalogram with higher quality and less noise enhancing the interpretability over the spectrogram obtained from STFT. That enabled a better TFR, additionally allowing a model-independent estimation of the RMS roughness.

Further details on the basics and applicability of time-frequency analysis on XRR data can be found in the Refs. [42–44].

The approach utilized here to extract information from XRR data is hierarchical beginning with analyzing a TFR based on the Wigner-Ville distribution (WVD) to infer the layers thicknesses and the stratification order (assuming a multilayer system) [45,46]. This is followed by a fit using the global optimization metaheuristic algorithm Differential Evolution (DE) [47], to arrive at refined values for the target parameters in particular the interfacial roughness. The WVD for a signal can be written as [48,49]:

where τ is the lag variable, s(t) is the analytic signal derived from the real signal x(t) and s*(t) is the complex analytic, i.e.:

with H{x(t)} being the Hilbert transform of the real signal.The WVD based TFR is used here due to its desirable properties as a bilinear TFR that decomposes a wave based on the density of its spectral content (energy), theoretically enabling high resolution readability when extracting features from the TFRs regarding non-stationary signals [50–52], that in some cases surpass the resolution obtained from CWT and STFT [52]. The bilinearity of the WVD introduces cross-terms in the TFRs that usually appear between the signals’ resolved features worsening the resolution, yet there are methods to eliminate those artefacts [50].

To illustrate the potential of a WVD-based TFR in analyzing signals, an example clarifies for the case of a time-dependent superposition in the appendices.

## 5. Data acquisition and analysis

#### 5.1 EUVR Data

From two different positions (Position 1 and Position 2) from the sample’s surface, two EUVR datasets (Fig. 1) were collected in the experimental endstation of the SX700 beamline by measuring the specular reflectance in the wavelength range 10 nm – 20 nm (in steps of 0.25 nm), for the angular range between 4.5° – 87° (grazing incidence, in steps of 1.5°).

The two datasets were collected to test the sample’s homogeneity. A total of 41 wavelength scans were taken in each. The small difference between the two reflectivity maps regarding the different positions demonstrates the sample’s high homogeneity (Fig. 1(c)). The relative uncertainties of these two measurements peak at the steep angles of incidence and short wavelengths, i.e., for low reflectance values.

To enlarge the spectral range and to further examine the sensitivity of the modelling formalism against contamination and oxidation layers, an additional dataset was collected from the sample after cleaning it with H*. This dataset consists of 64 wavelength scans in the spectral range 8 nm – 23.75 nm (in steps of 0.25 nm) for the angular range between 3° – 85.5° (in steps of 1.5°). All EUVR data were collected using a nominally pure S-polarized radiation. The large angular range from near normal to grazing incidence supports the uniqueness in determining the optical constants of the thin film [53].

#### 5.2 Model parametrization

For a direct parametrization of the sample’s structure the WVD-based TFR (the WV spectrum, see Prieto et al. [54]) of an XRR measurement at 1750eV from the sample was calculated (Fig. 2) using multitaper spectral estimates followed with applying a reduced interference distribution (RID) kernel [54,55]. The direct variable of the transform has been translated from the grazing angle to the *corrected* momentum transfer vector ${Q_{z^{\prime}\; }}$as given by [56]:

The conversion shown in (9) is to correct for the refraction. To obtain well-resolved features, the signal was renormalized prior to the transformation by taking $I({Q_{z^{\prime}\; }}) \cdot \; Q_{z^{\prime}}^4$ [44].

With an argument from Smigiel et al. [43], given the direct relation between the penetration depth of the incident X-ray beam in the probed sample and the grazing angle, one could argue that the features resolved at lower grazing angles (smaller momentum vectors) can be attributed to layers closer to the surface, while those features occurring at higher grazing angles are to indicate layers closer to the substrate. Three features were resolved in the WV spectrum shown in Fig. 2(b), designated A, B and C. With the settings of the DC-magnetron sputtering tool adjusted at depositing a ruthenium film with a thickness of 50 nm, we could determine that feature A represents the targeted ruthenium deposition. With the sample being handled in the ambient conditions before the measurement, oxidation of the film and contamination from volatile molecules are expected. That is assumed to explain the presence of a top surface layer, feature B, which is resolved at the smallest momentum transfer vector, so it is expected to be the uppermost. Feature C could represent the substrate’s native oxidation layer known previously from the manufacturer. We assume feature C to represent the substrate’s oxide. The discontinuity between the center positions of the two features B and C can be explained with the occurrence of the latter below another layer which is the ruthenium in this stratification.

Feature C presumably represents the substrate’s native oxide (known to exist from the manufacturer) instead of a diffusion layer with the substrate, since SiO_{2} has a lower enthalpy of formation than ruthenium silicides [57].

The gain from using time-frequency analysis here is the direct verification of the presence of a surface layer on top of the ruthenium film, this verification is critical since the determination of optical constants is sensitive to the presence of such surface layers yielding different results depending on the assumed model as it will be demonstrated afterwards.

The surface layer resolved from the XRR data developed in an uncontrolled environment (ambient conditions). It is expected that this surface layer is comprised of the ruthenium film’s oxide and additional physiosorbed volatile molecules from the ambience. Surface ultra-thin layers could have a gradient owing to factors like the moisture in air, and it is assumed here to be the case [58]. The former assumption will be corroborated with the outcomes of the EUVR datasets analysis. Optimally, modelling such a surface layer requires a transition-layer approach where a stack of transition layers whose density profile trails a smooth transition gradient that represents the examined layer [59]. This approach can model any roughness distribution, but it is very computationally intensive. The alternative here is to *effectively* model the surface layer as only two separate layers, a carbonaceous contamination layer (with reduced density) above an oxide layer. This approximation is deemed necessary to reduce computational intensiveness. The EUVR data will be further treated using this four-layer model.

A preliminary four-layer model (on a substrate) describing the sample structure can be presupposed, ascendingly, the substrate native oxide, ruthenium deposition, oxidation layer and a carbonaceous contamination. To model the interfacial imperfections, the XRR data was optimized using the DE algorithm utilizing a logarithmic objective function. The forward calculation in the optimization process was simulated using Parratt’s reflectivity formalism for layered systems with Névot-Croce damping factors introduced to compensate the effect of interfacial imperfections [60,61]. The optical parameters of ruthenium and its oxide were left variables in the optimization, given the anticipated reliable high sensitivity from a 50 nm layer for the former and to compensate for the uncertainties in the estimated parameters of the latter, since the exact oxide’s stoichiometry (RuO_{x}) and its density are undetermined while the optical parameters of other materials were fixed from the Center for X-Ray Optics (CXRO) database [62,63].

The density of the film was approximated from the critical angle of the Cu-Kα XRR profile measured immediately after coating. This is important since the critical angle region is known to have a shallow penetration depth of the sample, and is highly sensitive to the surface [64]. Hence, oxidation concurrent with additional contamination could yield too low density for the film

#### 5.3 Curse of dimensionality mitigation

The optimization of the discussed four-layer on a substrate model regarding the experimental settings of the two measurements corresponding to the spectral range 10 nm – 20 nm requires 424 parameters. With the parameters representing the optical constants, the structural characteristics, and other parameters relevant for the calculation (Table 3). An even larger number for the optimization of the dataset collected after H*-cleaning. Ideally, all parameters must be variables in the MCMC-Bayesian framework to reveal the correlations and to examine the propagation of uncertainties. With such a tremendous number of variables, applying MCMC- Bayesian inferences framework is impractical here as it requires an enormous computational effort for the chains to reach the equilibrium state for all parameters [65]. With the aim here to determine the optical constants of pure ruthenium, the optical constants of the ultra-thin layers and the substrate were taken from previously published data (Table 3 merely details for the case EUVR data collected before the H*-cleaning) [63,66–68].

Moreover, the optical constants of ruthenium oxide are to be *partially* fixed using the atomic scattering factors (ASF) representation since the complex index of refraction for a compound with N different atoms can be expressed as [62]:

where r_{e} is the classical electron radius, λ is the wavelength, n is the number of atoms of type j per unit volume. f is the complex ASF for an atom with a real part related to the dispersion and an imaginary part related to the absorption:

With an assumed oxide stoichiometry of RuO_{2} but variable density, the ASFs of ruthenium are part of array **p** and the ASFs of oxygen are fixed as taken from the CXRO database [63,68]. To further exclude additional variables in the forward calculation, the roughness parameters of the stratification are to be replaced with the RMS values obtained from the XRR fitting, considering that XRR measurements are better suited for indicating such values than EUVR measurements given their low penetration depth to probe the buried interfaces under the ruthenium deposition. In such a problem, the trustworthiness of the determined optical constants could be ameliorated with the minimization of the simulations’ variables, via incorporating input from various analysis methods [69].

#### 5.4 Pre-cleaning collected EUVR data MCMC-Bayesian inferences

With 91 parameters left in array **p** 2.5×10^{5} forward calculations based on Parratt’s formalism coupled with Névot-Croce factors were iterated [60,61]. The stationary target distribution of the MCMC chains was attained after ca. 3.6×10^{4} iterations considering all the parameters. The calculations were done using emcee, a Python implementation of Goodman & Weare’s Affine Invariant MCMC Ensemble sampler to generate the proposals regarding 500 walkers for each parameter [70].

The first 4×10^{4} samples of the MCMC chains were discarded for the subsequent analysis (burn-in), since most of them have a low probability to describe the target parameters. To reduce autocorrelations and to tackle computer memory limitations, only every 15^{th} step was saved from the MCMC chains (thinning of 15, (Fig. 3(a)).

For MCMC sampling, convergence can be indicated using several empirical convergence diagnostics. Since a single diagnostic might fail, three diagnostics were checked for the convergence in this study: Gelman-Rubin criterion, autocorrelation analysis and Geweke-plots. All of which gave a positive *indication* for convergence [71].

The calculated Gelman-Rubin factors for all the chains are below 1.1 which is the recommended cutoff value to stop MCMC sampling. Geweke Z-scores were all found to be in the interval (-2,2). Another very important check is to examine the serial correlations of the MCMC chains; the autocorrelation of the samples decreased and was oscillating around zero right before terminating the algorithm [71].

To emphasize the significance of the model used here, an additional MCMC simulation was applied on an oversimplified single-layer (on a substrate) model regarding the EUVR data collected form Position 2. The determined optical constants from each model will be compared eventually in a plot, albeit the discussion will mainly focus on the four-layer model.

#### 5.5 Post-cleaning collected EUVR data MCMC-Bayesian inferences

With 137 variable parameters 2.5×10^{5} forward calculations were iterated. The number of variables here is larger since the data collected after cleaning the samples has a larger spectral coverage than the data collected before from Position 1 and Position 2. The stationary target distribution was attained after ca. 4.5×10^{4} steps. The geometrical model retrieved from the modes of the chains (excluding the burn-in) is akin to that shown in Fig. 3(b)) except for the thickness of the carbon contamination layer and the density of the oxide layer (Fig. 4). This verifies the sensitivity of the modelling scheme to surface layers and shows the anticipated slight increase in the base film’s thickness upon the reduction of its oxide. Hydrogen reduces a ruthenium oxide layer, as stated by I. Nishiyama, ”leaving only the base element” ref. [72].

Since the XRR measurements were taken before cleaning the sample, the parameterization of the RMS values of the roughness parameters considering the two uppermost interfaces in the optimization of the data collected after cleaning the sample was slightly different than that of datasets collected before. The roughness parameters of the two uppermost interfaces were scaled proportionally to the thickness of the carbon contamination layer obtained from the measurements prior to cleaning using a ratio. This approximation prevents unrealistic values of roughness versus thickness. The RMS roughness of ruthenium-ruthenium oxide interface was assumed to be the same as before H*-cleaning, since H*-cleaning has as stated by I. Nishiyama, “a negligible influence on the surface roughness of ruthenium” ref. [72].

## 6. Determined optical constants discussion

For the two datasets collected from Position 1 and Position 2, the first 36 thousand iterations (burn-in period) were discarded, the chains were integrated and yielded Gaussian probability distributions with numerous parameters exhibiting correlations. For the determined ASFs of ruthenium, the projected uncertainty defined as the confidence interval with 3-σ of the imaginary parts is generally higher than that of the real parts. Figure 5 shows the determined ASFs of ruthenium regarding a wavelength of 13.5 nm, with their covariances and correlations to selected parameters.

The measurements’ experimental uncertainty includes the uncertainty of (i) intensity measurement, (ii) stability of the synchrotron radiation and the beamline, (iii) and the photodiode noise. The MCMC-Bayesian framework includes the polarization as a variable, to account for minimal unpolarized fraction of the incident radiation. The uncertainty induced from the alignment procedure was also compensated in the simulation by introducing a single angular off-set. Moreover, the effect of the beam’s divergence was simulated on the EUVR data at a wavelength of 8.0 nm, and the difference was within the experimental uncertainty for the experimental data presented here. Merely the shortest wavelength was selected for the simulation of divergence effects, because Kiessig fringes are most susceptible to blurring and those fringes are the most prominent in our data at the wavelength of 8.0 nm [75].

The two variables of [Eq. (5)] *a* and *b*, are expected to partially compensate for the uncertainties that have not propagated due to fixing some parameters in the calculations like the optical constants of the substrate and ultra-thin layers materials taken from literature (Table 3). Additionally, with an argument adapted from scatterometry [76], this uncertainty modelling is more effective than identifying and calculating all potential error sources individually. Comparing the simulated error parameters (*a* and *b* in [Eq. (5)]) for the case of the single-layer model and for the case of the four-layer model used in treating the EUVR data collected from Position 2, while *b* in the two cases is in the same order (10^{−5}), *a* is two-orders higher in the case of the over simplified single-layer model (10^{−2}). This indicates that the two error parameters function as anticipated. The mean of total simulated uncertainty in the case of the single-layer model is more than four-fold the mean of the simulated uncertainty in the case of the four-layer model.

The reflectivity maps simulated using the modes of the posterior probability distributions show a remarkable resemblance with the measured datasets. For example, the goodness-of-fit is demonstrated for the data collected after the cleaning process in Fig. 6.

The calculated optical constants from the determined ASFs [Eq. (10)] were compared with the optical constants retrieved from listed ASFs from two known published datasets (Fig. 7) [63,68,77]. The determined optical constants from the three EUVR datasets are very close, although the determination from each dataset was conducted separately. Generally, our data agree with those from CXRO [63,68], below 15 nm and some deviations above. Both disagree with the optical constants retrieved from Chantler’s database (Fig. 7) [77].

Also, the inferences from an oversimplified (single layer) model leads to significant variations above 15 nm (approximately). MCMC chains of parameters relevant to the experimental settings like the polarization were compared regarding the two models. The mode value in the case of the single-layer model is ca. 8.6% (linear polarization) where it does not exceed 2.7% for the case of a four-layer model, considering the 3 collected datasets. Presumably, this can be considered as another verification for accuracy of our model since the latter model yields more realistic results known of our experimental settings (Table 1).

The significant drop in both the contamination layer’s thickness and the oxide’s density upon H*-cleaning demonstrates the sensitivity of the adopted modelling scheme to surface contamination and oxidation. It also further supports the parametrized model from the WVD-based formalism.

Although the optical constants determined from the three datasets show agreement (in the overlapping range), presumably, the optical constants determined from Post H*-cleaning dataset (Table 4) are those with the highest accuracy, because the sample then had the least contaminates and surface oxide. This presumption is based on a result of a study showing that contamination is known to affect the uncertainties of the calculated parameters from XRR data [78]. Although contamination was considered in the model, it is still simplified by using the double-layer approach instead of allowing for a continuous gradient profile.

Although we used XRR data to independently determine the number of surface layers, there are still presumptions made. For example, taking uniform discrete layers in modelling the stratification of the sample, implicitly stipulates that the density profile of the ruthenium layer is flawless, and this among other presumptions have not been checked directly. This, however, is supported by extensive investigations of ruthenium thin films elsewhere [57,79], as an important materiel in EUVL technology.

Furthermore, the comparison between the optical constants and the presumed structural characteristics of the sample, determined before and after H*-cleaning to remove the surface contamination and oxidation corroborates the used strategy.

## 4. Conclusions and potential outcomes

From reflectivity data collected in a synchrotron facility, we have presented the determination of optical constants of ruthenium in the EUV sub-spectral range from 8 nm to 23.75 nm, using MCMC-Bayesian inferences framework addressing the confidence interval for each parameter. Ruthenium was chosen in this work as an example because it is particularly interesting for EUVL as it is used in many optical components there. For wavelengths below 15 nm, our optical data resemble the CXRO data [63,68], while there is a clear discrepancy to the values from Chantler [77].

We successfully applied the WVD based TFR on experimental XRR data. The significance of resolving surface layers directly from the experimental data without prior knowledge is emphasized in the model selection. The outcome of this work encourages further applications of time-frequency analysis given its ability to resolve ultra-thin layers directly, especially surface layers. This is of crucial value since the optical characteristics of thin films in the EUV spectral range are highly sensitive to the presence of such surface layers that are not easily resolved with a trial-and-error analysis scheme of the XRR data. Summarizing, time-frequency analysis is shown to be a valuable metrology tool for quasi-model independent investigation of thin film stacks.

The reliability of applying MCMC-Bayesian inferences was tested and verified, with emphasis on the sensitivity for surface layers, where three different datasets yield very close values for the ASFs, but significant differences are observed in the surface layers characteristics upon cleaning the sample.

This surface modification at the same layer thickness clearly shows the robustness of the evaluation procedure and modeling with respect to surface contamination and oxidation. Investigating ruthenium thin films of different thicknesses would be a valuable further proof but was skipped here because we primarily aim to show how to extract the maximum information content from simple reflectance measurements.

Although the determined optical constants in our previous work are generally close those of this work, this approach conveys the following advantages: (i) the uncertainties of the results are determined, (ii) a reliable stopping criterion of the simulation process is set, (iii) multi-modalities and correlations between the model’s variable parameters are examined, (iv) Presumably, the model’s deficiencies have been partially accounted for by introducing the two error parameters *a* and *b* [Eq. (5)].

The presented analysis scheme is based on open-source and highly accessible software packages, promoting inverse-problem solving strategies of EUVR data. Based on the presented scheme, this work will pave the way for novel approaches utilizing reflectometry for the determination of optical constants of other materials, where the presented strategy circumvents the requirements of retrieving the optical constants using Kramers-Kronig relations, which is often hindered by missing optical constants in parts of the electromagnetic spectrum, this is a particular problem for alloys which are prominent candidates for alternative absorber materials on EUV masks [6,80].

## Appendix

- • Material relevant for section 4 (Hybrid model parametrization)

To illustrate the potential of a WVD-based approach in analyzing signals, the following example clarifies for the case of a time-dependent superposition of different components (Fig. 8(a)). Where the superposition denoted as M(t) is defined as:

The calculated TFR from the simulated signal (Fig. 8) represents how the modulation vary over both time and frequency, simultaneously [54,55]. The different frequencies present in the signal can be analyzed then.

- • Material relevant for section 6 (Determined optical constants discussion)

The determined optical constants obtained from data collected after H*-cleaning, with their corresponding uncertainties as retrieved from the MCMC evaluation are tabulated below.

## Funding

Electronic Components and Systems for European Leadership.

## Acknowledgements

This project has received funding from the Electronic Component Systems for European Leadership Joint Undertaking under grant agreement No 783247 – TAPES3. This Joint Undertaking receives support from the European Union’s Horizon 2020 research and innovation program and from Netherlands, France, Belgium, Germany, Czech Republic, Austria, Hungary and Israel.

The authors thank Dr. Sebastian Heidenreich and Analía Fernández Herrero (both from the PTB) for the fruitful discussions that improved this work. The authors also thank Vivien Tritschler (PTB) for proofreading the article.

## Disclosures

The authors declare no conflicts of interest.

## Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon request.

## References

**1. **A. Erdmann, H. S. Mesilhy, P. Evanschitzky, V. Philipsen, F. J. Timmermans, and M. Bauer, “Perspectives and tradeoffs of absorber materials for high NA EUV lithography,” J. Micro/Nanolithogr., MEMS, MOEMS **19**(04), 041001 (2020). [CrossRef]

**2. **N. Brimhall, N. Herrick, D. D. Allred, R. S. Turley, M. Ware, and J. Peatross, “Characterization of optical constants for uranium from 10 to 47 nm,” Appl. Opt. **49**(9), 1581–1585 (2010). [CrossRef]

**3. **N. Brimhall, N. Herrick, D. D. Allred, R. S. Turley, M. Ware, and J. Peatross, “Measured optical constants of copper from 10 nm to 35 nm,” Opt. Express **17**(26), 23873–23879 (2009). [CrossRef]

**4. **V. Bakshi, * EUV Lithography*, 2nd ed. (SPIE Monographs, 2018), Chap. 5-7.

**5. **S. Bajt, Z. R. Dai, E. J. Nelson, M. A. Wall, J. Alameda, N. Nguyen, S. Baker, J. C. Robinson, J. S. Taylor, M. Clift, A. Aquila, E. M. Gullikson, and N. V. G. Edwards, “Oxidation resistance of Ru-capped EUV multilayers,” Proc. SPIE **5751**, 118–127 (2005). [CrossRef]

**6. **V. Philipsen, K. V. Luong, K. Opsomer, L. Souriau, J. Rip, C. Detavernier, A. Erdmann, P. Evanschitzky, C. Laubis, P. Hönicke, V. Soltwisch, and E. Hendrickx, “Mask absorber development to enable next-generation EUVL,” Proc. SPIE **11178**, 24 (2019). [CrossRef]

**7. **P. Rochus, J.-P. Halain, E. Renotte, D. Berghmans, A. Zhukov, J-F. Hochedez, T. Appourchaux, F. Auchère, L. K. Harra, U. Schühle, R. Mercier and the EUI consortium, “The Extreme Ultraviolet Imager (EUI) onboard the Solar Orbiter mission,” in the 60th International Astronautical Congress (2009).

**8. **S. Bowyer, J. J. Drake, and S. Vennes, “Extreme ultraviolet astronomy,” Annu. Rev. Astron. Astrophys. **38**(1), 231–288 (2000). [CrossRef]

**9. **R. Soufli and E. M. Gullikson, “Reflectance measurements on clean surfaces for the determination of optical constants of silicon in the extreme ultraviolet–soft-x-ray region,” Appl. Opt. **36**(22), 5499–5507 (1997). [CrossRef]

**10. **R. Soufli and E. M. Gullikson, “Absolute photoabsorption measurements of molybdenum in the range 60–930 eV for optical constant determination,” Appl. Opt. **37**(10), 1713–1719 (1998). [CrossRef]

**11. **R. Soufli, A. L. Aquila, F. Salmassi, M. Fernández-Perea, and E. M. Gullikson, “Optical constants of magnetron-sputtered boron carbide thin films from photoabsorption data in the range 30 to 770 eV,” Appl. Opt. **47**(25), 4633–4639 (2008). [CrossRef]

**12. **A. Gottwald, K. Wiese, U. Kroth, and M. Richter, “Uncertainty analysis for the determination of B_{4}C optical constants by angle-dependent reflectance measurement for 40 nm to 80 nm wavelength,” Appl. Opt. **56**(20), 5768–5774 (2017). [CrossRef]

**13. **G. D. Scott, “Optical constants of thin-film materials*,” J. Opt. Soc. Am. **45**(3), 176–179 (1955). [CrossRef]

**14. **J. T. Cox, G. Hass, J. B. Ramsey, and W. R. Hunter, “Reflectance and optical constants of evaporated ruthenium in the vacuum ultraviolet from 300 to 2000Å*,” J. Opt. Soc. Am. **64**(4), 423–428 (1974). [CrossRef]

**15. **I. A. Artioukov, B. R. Benware, J. J. Rocca, M. Forsythe, Y. A. Uspenskii, and A. V. Vinogradov, “Determination of XUV optical constants by reflectometry using a high-repetition rate 46.9-nm laser,” IEEE J. Sel. Top. Quan. Elec. **5**(6), 1495–1501 (1999). [CrossRef]

**16. **N. Brimhall, “Extreme Ultraviolet Polarimetry with Laser-Generated High-Order Harmonics: Characterization of Uranium,” PhD thesis (Brigham Young University, Department of Physics and Astronomy, 2009).

**17. **N. Koster, B. Mertens, R. Jansen, A. van de Runstraat, F. Stietz, M. Wedowski, H. Meiling, R. Klein, A. Gottwald, F. Scholze, M. Visser, R. Kurt, P. Zalm, E. Louis, and A. Yakshin, “Molecular contamination mitigation in EUVL by environmental control,” Microelectron. Eng. **61-62**, 65–76 (2002). [CrossRef]

**18. **M. M. J. W. van Herpen, D. J. W. Klunder, W. A. Soer, R. Moors, and V. Banine, “Sn etching with hydrogen radicals to clean EUV optics,” Chem. Phys. Lett. **484**(4-6), 197–199 (2010). [CrossRef]

**19. **W. M. Haynes (Editor), * CRC Handbook of Chemistry and Physics*, 95th ed. (CRC, 2014), Sect. 4.

**20. **K. Sakurai, M. Mizusawa, and M. Ishii, “Significance of frequency analysis in X-ray reflectivity: towards analysis which does not depend too much on models,” Trans. Mater. Res. Soc. Jpn. **33**(3), 523–528 (2008). [CrossRef]

**21. **A. Schalchli, J. J. Benattar, and C. Licoppe, “Accurate measurements of the density of thin Silica films,” Europhys. Lett. **26**(4), 271–276 (1994). [CrossRef]

**22. **F. Scholze, B. Beckhoff, G. Brandt, R. Fliegauf, A. Gottwald, R. Klein, B. Meyer, U. D. Schwarz, R. Thornagel, J. Tuemmler, K. Vogel, J. Weser, and G. Ulm, “High-accuracy EUV metrology of PTB using synchrotron radiation,” Proc. SPIE **4344**, 402–413 (2001). [CrossRef]

**23. **C. Laubis, A. Kampe, C. Buchholz, A. Fischer, J. Puls, C. Stadelhoff, and F. Scholze, “Characterization of the polarization properties of PTB's EUV reflectometry system,” Proc. SPIE **7636**, 76362R (2010). [CrossRef]

**24. **A. Haase, “Multimethod Metrology of Multilayer Mirrors Using EUV and X-Ray Radiation,” PhD thesis (Technische Universität Berlin, Inst. Optik und Atomare Physik, 2017).

**25. **C. Laubis, A. Barboutis, C. Buchholz, A. Fischer, A. Haase, F. Knorr, H. Mentzel, J. Puls, A. Schönstedt, M. Sintschuk, V. Soltwisch, C. Stadelhoff, and F. Scholze, “Update on EUV radiometry at PTB,” Proc. SPIE **9776**, 977627 (2016). [CrossRef]

**26. **F. Senf, U. Flechsig, F. Eggenstein, W. Gudat, R. Klein, H. Rabus, and G. Ulm, “A plane-grating monochromator beamline for the PTB undulators at BESSY II,” J. Synch. Rad. **5**(3), 780–782 (1998). [CrossRef]

**27. **J. Lubeck, B. Beckhoff, R. Fliegauf, I. Holfelder, P. Hönicke, M. Müller, B. Pollakowski, F. Reinhardt, and J. Weser, “A novel instrument for quantitative nanoanalytics involving complementary X-ray methodologies,” Rev. Sci. Instrum. **84**(4), 045106 (2013). [CrossRef]

**28. **Q. Saadeh, V. Soltwisch, P. Naujok, and F. Scholze, “Validation of optical constants in the EUV spectral range,” Proc. SPIE **11147**, 49 (2019). [CrossRef]

**29. **D. S. Sivia, W. I. F. David, K. S. Knight, and S. F. Gull, “An introduction to Bayesian model selection,” Phys. D **66**(1-2), 234–242 (1993). [CrossRef]

**30. **M. P. Hobson and J. E. Baldwin, “Markov-chain Monte Carlo approach to the design of multilayer thin-film optical coatings,” Appl. Opt. **43**(13), 2651–2660 (2004). [CrossRef]

**31. **D. Windover, D. L. Gil, J. P. Cline, A. Henins, N. Armstrong, P. Y. Hung, S. C. Song, R. Jammy, and A. Diebold, “NIST method for determining model-independent structural information by X-ray reflectometry,” AIP Conf. Proc. **931**(1), 287–291 (2007). [CrossRef]

**32. **B. W. Treece, P. A. Kienzle, D. P. Hoogerheide, C. F. Majkrzak, M. Lösche, and F. Heinrich, “Optimization of reflectometry experiments using information theory,” J. Appl. Crystallogr. **52**(1), 47–59 (2019). [CrossRef]

**33. **A. R. J. Nelson and S. W. Prescott, “refnx: neutron and X-ray reflectometry analysis in Python,” J. Appl. Crystallogr. **52**(1), 193–200 (2019). [CrossRef]

**34. **D. S. Sivia and J. R. P. Webster, “The Bayesian approach to reflectivity data,” Phys. B **248**(1-4), 327–337 (1998). [CrossRef]

**35. **S. Sharma, “Markov Chain Monte Carlo methods for Bayesian data analysis in astronomy,” Annu. Rev. Astron. Astrophys. **55**(1), 213–259 (2017). [CrossRef]

**36. **J. Goodman and J. Weare, “Ensemble samplers with affine invariance,” Commun. Appl. Math. Comput. Sci. **5**(1), 65–80 (2010). [CrossRef]

**37. **S. Heidenreich, H. Gross, and M. Bar, “Bayesian approach to the statistical inverse problem of scatterometry: comparison of three surrogate models,” Int. J. Uncertain. Quan. **5**(6), 511–526 (2015). [CrossRef]

**38. **H. Kiessig, “Untersuchungen zur totalreflection von röntgenstrahlen,” Ann. Phys. **402**(6), 715–768 (1931). [CrossRef]

**39. **A. van der Lee, “Grazing incidence specular reflectivity: theory, experiment, and applications,” Solid State Sci. **2**(2), 257–278 (2000). [CrossRef]

**40. **E. Chason and T. M. Mayer, “Thin film and surface characterization by specular X-ray reflectivity,” Crit. Rev. Solid State Mater. Sci. **22**(1), 1–67 (1997). [CrossRef]

**41. **Y. O. Volkov, I. V. Kozhevnikov, B. S. Roshchin, E. O. Filatova, and V. E. Asadchikov, “Model approach to solving the inverse problem of X-ray reflectometry and its application to the study of the internal structure of hafnium oxide films,” Crystallogr. Rep. **58**(1), 160–167 (2013). [CrossRef]

**42. **E. Smigiel, A. Knoll, N. Broll, and A. Cornet, “Determination of layer ordering using sliding-window Fourier transform of x-ray reflectivity data,” Modelling Simul. Mater. Sci. Eng. **6**(1), 29–34 (1998). [CrossRef]

**43. **E. Smigiel and A. Cornet, “Characterization of a layer stack by wavelet analysis on x-ray reflectivity data,” J. Phys. D: Appl. Phys. **33**(15), 1757–1763 (2000). [CrossRef]

**44. **I. R. Prudnikov, R. J. Matyi, and R. D. Deslattes, “Wavelet transform approach to the analysis of specular x-ray reflectivity curves,” J. Appl. Phys. **90**(7), 3338–3346 (2001). [CrossRef]

**45. **R. Clinciu, “The application of the Wigner distribution function to X-ray reflectivity,” M.S. thesis (University of Warwick, Department of Engineering, 1992).

**46. **B. Tanner and D. K. Bowen, * X-Ray Metrology in Semiconductor Manufacturing*, (CRC, 2006), Chap. 13.

**47. **R. Storn and K. Price, “Differential Evolution – a simple and efficient Heuristic for global optimization over Continuous Spaces,” J. Global Opt. **11**(4), 341–359 (1997). [CrossRef]

**48. **L. Choen, “Time-frequency distributions-a review,” Proc. IEEE **77**(7), 941–981 (1989). [CrossRef]

**49. **X. Wang, Y. Xue, W. Zhou, and J. Luo, “Spectral decomposition of seismic data with variational mode decomposition-based Wigner–Ville Distribution,” IEEE J. Se. Top. App. E. Obs. Re. Sen. **12**(11), 4672–4683 (2019). [CrossRef]

**50. **S. C. Bradford, “Time-Frequency Analysis of Systems with Changing Dynamic Properties,” Ph.D. thesis (California Institute of Technology, 2007).

**51. **O. Rioul and M. Vetterli, “Wavelets and signal processing,” IEEE SP Mag. **8**(4), 14–38 (1991). [CrossRef]

**52. **J. P. Chi and C. T. Russell, “Use of the Wigner-Ville distribution in interpreting and identifying ULF waves in triaxial magnetic records,” J. Geophys. Res. **113**(A1), A01218 (2008). [CrossRef]

**53. **K. Lamprecht, W. Papousek, and G. Leising, “Problem of ambiguity in the determination of optical constants of thin absorbing films from spectroscopic reflectance and transmittance measurements,” Appl. Opt. **36**(25), 6364–6371 (1997). [CrossRef]

**54. **G. Prieto, R. Parker, and F. L. Vernon, “A Fortran 90 library for multi taper spectrum analysis,” Com. Geo. **35**(8), 1701–1710 (2009). [CrossRef]

**55. **L. Krischer, “mtspec Python wrappers 0.3.2,” Zenodo. (2016) https://doi.org/10.5281/zenodo.321789

**56. **O. Durand and N. Morizet, “Fourier-inversion and wavelet-transform methods applied to X-ray reflectometry and HRXRD profiles from complex thin-layered heterostructures,” Appl. Surf. Sci. **253**(1), 133–137 (2006). [CrossRef]

**57. **R. Coloma Ribera, R. W. E. van de Kruijs, S. Kokke, E. Zoethout, A. E. Yakshin, and F. Bijkerk, “Surface and sub-surface thermal oxidation of thin ruthenium films,” Appl. Phys. Lett. **105**(13), 131601 (2014). [CrossRef]

**58. **M. I. Mazuritskiy and A. A. Novakovich, “Surface transition-layer model used to study the fine structure of X-ray reflection spectra,” J. Sur. Inv. X-ray. Syn. Neut. Tech. **8**(6), 1291–1296 (2014). [CrossRef]

**59. **K. Stoev and K. Sakurai, “Recent theoretical models in grazing incidence x-ray reflectometry,” Rigaku J. **14**(2), 22–37 (1997).

**60. **L. G. Parratt, “Surface studies of solids by total reflection of X-Rays,” Phys. Rev. **95**(2), 359–369 (1954). [CrossRef]

**61. **L. Névot and P. Croce, “Caractérisation des surfaces par réflexion rasante de rayons X. Application à l’étude du polissage de quelques verres silicates,” Rev. Phys. Appl. **15**(3), 761–779 (1980). [CrossRef]

**62. **B. L. Henke, E. M. Gullikson, and J. C. Davis, “X-Ray interactions: photoabsorption, scattering, transmission, and reflection at E = 50-30,000 eV, Z = 1-92,” At. Data Nucl. Data Tables **54**(2), 181–342 (1993). [CrossRef]

**63. **An updated online version of [62] is available, https://henke.lbl.gov/optical_constants/asf.html.

**64. **R. S. Becker, J. A. Golovchenko, and J. R. Patel, “X-Ray evanescent-wave absorption and emission,” Phys. Rev. Lett. **50**(3), 153–156 (1983). [CrossRef]

**65. **L. Chen, “Curse of Dimensionality,” in L. Liu and M. T. Özsu, *Encyclopedia of Database Systems* (Springer, 2009). https://doi.org/10.1007/978-0-387-39940-9_133.

**66. **E. Filatova, V. Lukyanov, R. Barchewitz, J.-M. Andre, M. Idir, and P. Stemmler, “Optical constants of amorphous SiO2 for photons in the range of 60–3000 eV,” J. Phys.: Condens. Matter **11**(16), 3355–3370 (1999). [CrossRef]

**67. **A. Andrle, P. Hönicke, J. Vinson, R. Quintanilh, Q. Saadeh, S. Heidenreich, F. Scholze, and V. Soltwisch, “The anisotropy in the optical constants of quartz crystals for soft X-rays,” J. Appl. Crystallogr. **54**(2), 402–408 (2021). [CrossRef]

**68. **U. Schlegel, “Determination of the optical constants of Ruthenium in the EUV and soft x-ray region using synchrotron radiation,” Diploma thesis (Technische Fachhochschule Berlin, 2000).

**69. **L. V. Rodríguez-de Marcos, S. M. P. Kalaiselvi, O. B. Leong, P. K. Das, M. B. H. Breese, and A. Rusydi, “Optical constants and absorption properties of Te and TeO thin films in the 13-14 nm spectral range,” Opt. Express **28**(9), 12922–12935 (2020). [CrossRef]

**70. **D. Foreman-Mackey, D. W. Hogg, D. Lang, and J. Goodman, “emcee: The MCMC hammer,” Publ. Astron. Soc. Pac. **125**(925), 306–312 (2013). [CrossRef]

**71. **V. Roy, “Convergence diagnostics for Markov Chain Monte Carlo,” Annu. Rev. Stat. Appl. **7**(1), 387–412 (2020). [CrossRef]

**72. **I. Nishiyama, H. Oizumi, K. Motai, A. Izumi, T. Ueno, H. Akiyama, and A. Namiki, “Reduction of oxide layer on Ru surface by atomic-hydrogen treatment,” J. Vac. Sci. Technol., B: Microelectron. Process. Phenom. **23**(6), 3129 (2005). [CrossRef]

**73. **D. Foreman-Mackey, “corner.py: Scatterplot matrices in Python,” J. Open Source Software **1**(2), 24 (2016). [CrossRef]

**74. **J. D. Hunter, “Matplotlib: A 2D graphics environment,” Com. Sci. Eng. **9**(3), 90–95 (2007). [CrossRef]

**75. **K. Stoev and K. Sakurai, “Aberration effects in quick X-Ray reflectivity of curved samples,” IOP Conf. Ser.: Mater. Sci. Eng. **24**, 012014 (2011). [CrossRef]

**76. **A. Fernández Herrero, M. Pflüger, J. Probst, F. Scholze, and V. Soltwisch, “Applicability of the Debye-Waller damping factor for the determination of the line-edge roughness of lamellar gratings,” Opt. Express **27**(22), 32490–32507 (2019). [CrossRef]

**77. **C. T. Chantler, “Theoretical form factor, attenuation, and scattering tabulation for Z=1–92 from E=1–10 eV to E=0.4–1.0 MeV,” J. Phys. Chem. Ref. Data **24**(1), 71–643 (1995). [CrossRef]

**78. **D. L. Gil and D. Windover, “Limitations of x-ray reflectometry in the presence of surface contamination,” J. Phys. D: Appl. Phys. **45**(23), 235301 (2012). [CrossRef]

**79. **R. C. Ribera, R. W. E. van de Kruijs, A. E. Yakshin, and F. Bijkerk, “Determination of oxygen diffusion kinetics during thin film ruthenium oxidation,” J. Appl. Phys. **118**(5), 055303 (2015). [CrossRef]

**80. **V. Luong, V. Philipsen, E. Hendrickx, K. Opsomer, C. Detavernier, C. Laubis, F. Scholze, and M. Heyns, “Ni-Al alloys as alternative EUV mask absorber,” Appl. Sci. **8**(4), 521 (2018). [CrossRef]