A NEW ALGORITHM FOR ANALYSIS OF EXPERIMENTAL MÖSSBAUER SPECTRA

A new approach to analyze the nuclear gamma resonance (NGR) spectra is presented and justified in the paper. The algorithm successively spots the Lorentz lines in the experimental spectrum by a certain optimization procedures. In Mössbauer spectroscopy, the primary analysis is based on the representation of the transmission integral of an experimental spectrum by the sum of Lorentzians. In the general case, a number of lines and values of parameters in Lorentzians are unknown. The problem is to find them. In practice, before the experimental data processing, one elaborates a model of the Mössbauer spectrum. Such a model is usually based on some additional information. Taking into account physical restrictions, one forms the shape of the lines which are close to the normalized experimental Mössbauer spectrum. This is done by choosing the remaining free parameters of the model. However, this approach does not guarantee a proper model. A reasonable way to construct a structural NGR spectrum decomposition should be based on its model-free analysis. Some modelfree methods of the NGR spectra analysis have been implemented in a number of known algorithms. Each of these methods is useful but has a limited range of application. In fact, the previously known algorithms did not react to hardly noticeable primary features of the experimental spectrum, but identify the dominant components only. In the proposed approach, the difference between the experimental spectrum and the known already determined part of the spectral structure defines the next Lorentzian. This method is effective for isolation of fine details of the spectrum, although it requires a well-elaborated algorithmic procedure presented in this paper.


Introduction
In Mössbauer spectroscopy, the primary analysis of the NGR spectrum structure is based on the representation of the transmission integral of an experimental spectrum f (x) by the following sum of Lorentzians: where x s is the position of the maximum of the sth Lorentz line on the velocity scale, and 2b s is its width.In the general case, the number n of lines and the values of parameters A s , x s , and b s are unknown.The problem is to find them for a given function f (x).The Mössbauer spectrum is always measured in some bounded velocity range.The maximal number of lines in the spectrum is limited by the number of probable positions of the Mössbauer atom in the crystal lattice and by the nature of changes of the nuclear energy levels (e.g., isomer shift, quadropole splitting, or hyperfine splitting).For most of metallic alloys and oxide compounds, the evaluation of the structural parameters of the Mössbauer spectrum is a mathematically difficult nonlinear problem.
In practice, before the experimental data processing, one elaborates a model of the Mössbauer spectrum.Such a model is usually based on some additional information on the concentration, the structural state, and the charge state of Mössbauer atoms in the material.Then, taking into account physical restrictions, one forms the shape of the lines which are close (in the sense of the minimum of χ 2 -criterion) to the normalized experimental Mössbauer spectrum.This is done by choosing the remaining free parameters of the model.However, this approach does not guarantee a proper, physically approved model.A reasonable way to construct a structural NGR spectrum decomposition should be based on its model-free analysis.This can be done directly by solving problem (1) and by further improving the model based on the results of other methods [1].
Some model-free methods of the NGR spectra analysis have been implemented in algorithms for the density distribution of hyperfine fields, for the density distribution of isomer shifts with lines in the Lorentz or Gauss forms [2], and for filtering and reducing noises [3].Each of these methods is useful but has a limited range of application.In fact, the previously known algorithms did not react to hardly noticeable primary features of the spectrum, but identify the dominant components only.This is caused by the least squares methods which are applied to the whole experimental spectrum without any analysis of its details.In the presented approach, the difference between the experimental spectrum and the known part of the spectral structure defines the next Lorentzian.This method is effective for isolation of fine details of the spectrum, although it requires (see [4]) a well-elaborated algorithmic procedure presented in this paper.

Results and discussion
Below we give a detailed description and a proof of the new algorithm for the Mössbauer spectra analysis.Let be the convolution of the functions ϕ s (x) and ϕ(x).Then and where the maximum is attained for t = x s only.Equality (4) follows easily from the following formula: Indeed, this equality implies the relations where p > 0, ρ = (p 2 + ξ 2 ) 1/2 , and ϕ = arg(p + iξ) (ϕ = 0 for ξ = 0), and the equality is attained only at the point ξ = 0.
Furthermore, we assume that and consider only derivatives of even order when the absolute value sign in (4) can be omitted, since where p = b + b s > 0, and ξ = t − x s .
Remarks.By the assertion, for large k, the left-hand side of ( 6) is close to the maximum of the first term (s = 1) of the sum in (7) attained by (4) for s = 1 at the point t = x 1 .Hence, if k is sufficiently large, we can approximate the value x 1 by a maximum point of the left-hand side of (7) which actually coincides with the left-hand side of (6).

Numerical algorithm
Suppose that the function f (x) has a unique local maximum, i.e., there is only one point x 1 such that f (x 1 ) = max x f (x).
We set and find the parameter b 1 and points and Alternatively, they can be found by keeping the value of the half-width of f (x) and using the following simple formulas: If the difference f (x)−ϕ 1 (x) is small enough, i.e., is comparable with the accuracy of the evaluation of the function f (x), the algorithm terminates by setting f (x) = ϕ 1 (x).Otherwise, or, if there are several points of local maxima of f (x), a different process is applied.This process is based on the asymptotic behavior of (6) and is described below.
For arbitrary a > 0, and b > 0, and for sufficiently large k, it is needed to find a point of maximum of the left-hand side of (6) and the maximum value.Further, let us consider this point as an approximate value for x 1 .After that values a 2 and a positive integer k (1) are found such that the left-hand side of (6) is greater than A = max x f (x) for k = k (1) and a = a (1) 2 and is smaller than ε for k = k (1) and a = a 1 , where ε is an admissible error for computing A 1 from the right-hand side of (6).Next, the segment [a 2 ] is divided into two equal parts.One of these parts, which satisfies assumptions analogous to those for [a 2 ] for some k (2) (≥ k (1) ), is taken as the next segment [a 2 ].After several iterations, we eventually obtain a number k (ν) and a segment [a 2 with error at most δ.
Taking the logarithm of both sides of equality (9) implies that Hence, Therefore, these values are small for large ν, and then the approximate equality holds.Once the parameters x 1 , b 1 , and A 1 are found, we can introduce the function repeat the same operations for this function as we did for f (x), and find x 2 , b 2 , and A 2 .Then these iterations must be continued until the number A n+1 becomes small enough.Let us remark that, although the algorithm uses the differentiation of the convolution of the experimental function f (x) with the function ϕ(x, b), first, one can differentiate the integrand and then compute the convolution ∂ r ϕ ∂b r * f .This order is preferable, because the numerical integration is a more regular operation than the incorrect numerical differentiation, and the derivatives of ϕ(x, b) are expressed analytically.
Let us now write explicitly the positive (see (9)) argument of the logarithm in the left-hand side of (10) whose maximum is to be found.Starting, as above, from the formula