EP3274863A1 - Method and device for detecting radioelements - Google Patents

Abstract

The invention relates to a method for determining the nature of the radioelements contained in an object and the activity thereof, characterised in that it comprises at least the following steps: a first phase (31) of digital simulation of spectrometric responses for an incident energy set E and a measured output energy set E', in order to obtain a simulated data set; a second phase of non-parametric regression (32) on the simulated data, non-parametric estimation of the quantity representing the joint probability of the triplets (E, E', y) from simulated points (Ei, Ei', yij) in order to deduce therefrom a meta-model S(E, E') for any energy couple (E, E') on a continuous function, from the meta-model S(E, E'), (34) determining the nature and the activity of the radioelements contained in the object.

Description

 METHOD AND DEVICE FOR DETECTING RADIOELEMENTS

The object of the invention relates to a method and a device for detecting radioelements contained in an object, to determine their nature and their activity. In particular, it makes it possible to take into account spurious effects originating from the environment in which the object to be analyzed is located.

Gamma spectrometry is a technique for detecting gamma radiation emitted spontaneously by radioelements. This technique makes it possible to identify gamma emitting radioelements and to measure their activity. To achieve qualitative and quantitative objectives, the treatment of a gamma spectrum is based on the analysis of the peaks present in the spectrum. These correspond to the gamma radiation emitted by a specific radioelement, for example ¹³⁷ Cs emitting a gamma at 661.5 keV or ⁶⁰ Co emitting a gamma at 1173.2 keV and a gamma at 1332.5 keV. The presence of a peak in a gamma spectrum means that the energy of the incident gamma radiation has been fully absorbed in the detector, following one or more interactions, photoelectric effect, Compton effect or pair production, occurring within this last.

 The detectors used to perform gamma spectrometry measurements can be separated into two categories: scintillators and semiconductors. The quality of a gamma spectrum can be evaluated by two parameters specific to the detector used: the energy resolution, that is to say the ability to separate two near peaks in energy defined as the half-height width of the detector. a peak present in the spectrum at a given energy, and the detection efficiency, the denser the material, the greater the probability of absorbing all of the incident gamma radiation. These two parameters vary significantly depending on the detector used.

Inorganic scintillators (Nal, BGO) have a degraded resolution compared to semiconductor detectors (CZT or HPGe). Nevertheless, they generally exhibit greater detection efficiency. Indeed, because of their low cost, it is possible to manufacture them in large dimensions. FIG. 1 shows, in a channel energy diagram, blow per channel, an example of spectrum S (1) obtained using an HPGe detector and a ¹³⁷ Cs source emitting gamma radiation at 661, 5 keV. Figure 2 compares S (HPGe), S (CZT), S (Nal) spectra obtained using a uranium sample and different types of detectors, with different performance in terms of energy resolution.

It is also known to use plastic scintillators for detecting gamma radiation. The cost of manufacturing such scintillators is quite low, they can be manufactured in large dimensions for different types of application. Because of their intrinsic characteristics and despite their very large dimensions compared to other detectors, the probability that an incident photon will deposit all its energy in a photoelectric detector is extremely low. FIG. 3 shows a spectrum S (II) obtained using a plastic scintillator of the EJ200 type, in the presence of ⁶⁰ Co and ¹³⁷ Cs sources. The characteristic forms of these spectra are due to the Compton fronts following the interactions within the detector. It can thus be noted the absence of the photoelectric peaks characteristic of the presence of ¹³⁷ Cs and ⁶⁰ Co. For this reason, these detectors are generally used to perform total count measurements, detection of all the gamma radiation interacting within the detector without information on their energy. A traditional approach consisting of measuring the net area of the different peaks of interest is therefore not possible and an alternative solution must be used if it is desired to identify the radioelements present in the spectrum.

A solution for treating a gamma spectrum, identifying the radioelements transmitters and measure the activity of the latter, is to analyze the spectrum as a whole and not to restrict only to the analysis of the photoelectric peaks. The approach followed is to express spectrum analysis in the following matrix form:

S = H.A, with

 S: a signal matrix, matrix of dimensions (number_channels, 1), where the variable number_channels is equal to the number of channels of the gamma spectrum, A: an activity matrix, matrix of dimensions (nbre_energies_incidents, 1). The variable nbre_energies_incidentes corresponds to the number of incident energies defined by the user and taken into account in the reconstruction. These incident energies correspond to the number of elementary volume elements of a classical problem of emission tomography,

 H: a projector of the problem, matrix of dimensions (number_channels, number_energies_incidentes). This matrix includes all the detection efficiencies involved in the problem. By way of example, an element hij of this matrix corresponds to the probability that a photon of incident energy j is detected in channel i of the gamma spectrum.

 The signal matrix S corresponds to the result of the measurement, to the gamma spectrum to be analyzed, the projector matrix is generally calculated by simulation, the matrix A corresponds to the result of the reconstruction.

 FIG. 4A and FIG. 4B illustrate an example of the H projector concept, for an incident energy of 660 keV, FIG. 4A. The spectrum S (IV) of FIG. 4B corresponds to the incident energy corresponding to an entire column of the projector H. The projector thus contains the spectral response of the detector for each incident energy taken into account in the problem and defined by the user. . It is on this grid of incident energies that the reconstruction step will be carried out.

Different methods make it possible to take into account the gamma spectrum in order to carry out a qualitative and quantitative analysis, allowing the identification of the radioelements present in the spectrum and the activity of the radioelement present in the spectrum of the incident energies. This analysis methodology is illustrated in FIGS. 5 and 6. FIG. 5 represents the gamma spectrum to be analyzed, corresponding to a simulation of the response of a plastic scintillator type EJ200, in the presence of three gamma sources ( ²⁴¹ Am, ¹³⁷ Cs, ⁶⁰ Co). The reconstructed spectrum Sr from the results of the analysis and the incident energies is shown in dashed lines in the figure. FIG. 6 illustrates the result of a reconstruction carried out using an ML-EM type algorithm, the incident energies reconstructed from the spectrum of FIG. 5.

 Traditional methods of analysis despite all their interest and advantages have some limitations. One of the intrinsic limitations is related to the discrete character of the grid chosen to carry out the reconstruction and obtain the incidental energies. In all conventional approaches, such as the ML-EM type, this grid is discrete, following a step defined by the user. This step is one of the limitations during the reconstruction step since it will only be possible to reconstruct on the grid defined initially. Depending on the complexity of the problem being addressed, the pitch between each incident energy can be adapted. A fine step will allow greater finesse during the reconstruction, but will have a direct impact on the calculation of the projector. Since the latter is generally calculated by simulation, a finer step will result in an increase in computation time, which may become prohibitive for certain applications.

An alternative would be to look for solutions in the form of a discrete spectrum whose support is a continuous space. Such an approach to be integrable within the framework of an ML-EM algorithm would require to calculate the projector on the fly by Monte-Carlo simulation for each incident energy during the iterations of the algorithm. As the values of incident energies evolve according to the iterations, computations of answers are required at each iteration without the previously calculated values being able to be reused. The order of magnitude of calculation of such an approach proves prohibitive for the current calculators. In ML-EM type reconstruction, no priori is introduced during the reconstruction step. For the representativity of the spotlight, one of the main limitations of this type of approach is related to the representativeness of the projector used during reconstruction. The latter is generally obtained following a simulation step, using a model describing the detector used and its environment. Nevertheless, it is generally difficult to arrive at a fine modeling, knowledge of the geometry of the detector, limitations of the code of computation which model only imperfectly the interactions within the detector, absence of taking into account of certain physical phenomena, as the collection of visible photons for example in the case of plastic scintillators. If the modeling of the detector and / or the environment is imperfect, a more or less important bias will be present during the reconstruction.

Another approach, rather than considering a grid of incident energies on which one would perform the reconstruction is to use a database. In this case, we no longer simulate a grid of incident energies, but directly the signature of a given element, for example for the ⁶⁰ Co we will directly simulate the response of the detector for two incident energies, emitted respectively at 1173.2 keV and 1332.5 keV. With this approach the reconstruction is performed directly on a grid of radioelements. It also has the advantage of converging faster than the grid approach, because the number of components to be taken into account is lower. One of the disadvantages of this method is its limitation to take into account only the radioelements initially present in the database. If the analysis of a gamma spectrum involves a radioelement missing from the database, the reconstruction will be erroneous.

 In the remainder of the description, reference is made to the expression

"Input spectrum" a set of incident energies defined by a user and under the expression "output spectrum" a set of measured energies. The meta-given word and the word projector designate the same object.

The invention relates to a method for determining the nature of the radioelements present in an object and their activity, characterized in that it comprises at least the following steps: A first phase of numerical simulation of spectrometric responses for a set of incident energy E and a set of measured output energy E ', in order to obtain a set of simulated data,

A second phase of non-parametric regression on the simulated data, nonparametric estimation of the quantity representing the joint probability of the triples (E, E ', y) from simulated points (Ei, Ei', yij) in order to deduce a meta-model S (E, E ') for any energy pair (E, E') on a continuous function,

· From the meta-model S (E, E '), the determination of the nature and the activity of the radioelements present in the object.

 According to an alternative embodiment, the method comprises the following steps:

 • we use n a number of points in the input grid, energies E, and n 'a number of points in the output grid, energies £ ",

• for i = l, ..., n and j = l, ..., n ', we have the characteristic data of the calculated spectral intensities À _{£ y} - for an input energy _£ and an output energy E ,

 • a nonparametric estimation method of the quantity f (E, E ', y) representing the joint probability density of the triplets is used

(E, E ', y) from the simulated points (Ει, Ε, γ ^), and we deduce a model S (E, £ ") for all (E, £") E KL ² where IL is a space continued :

where the symbol E (y) represents the expected value of the random variable y, y is deduced from the calculated spectral intensities λ ·.

The values λ _£ can be determined _; - using Monte-Monte software

Carlo, the X data being considered as outputs of a Poisson distribution whose intensity is estimated by means of a non-parametric regression procedure and introduce y _tj = log (À _{£ y} - + ε) where 0 <ε «1, then

we approximate the probability distribution of y _tj for sufficiently large values of> 10 by a Gaussian law of mean log and variance, for the joined law f (E, E ', y), we choose a mixture by Dirichlet process (DPM) as a priori distribution and we express the random distribution in a sum over the infinite of f _{Q] i} components of f {E, E, y), f (E, E ^' , y) = ^ w _k fe _k (E, E ^' , y)

 k = l

parameterized by 9 _k the parameter associated with the k ^th component of G random measure defined by

G (-) = w _fc <¾ (·)

 fc = l

with w ₁ = V _x , w _k = V _k a) and 5 _U (-) represents the Dirac function localized in u and Beta (a, è), for 0 <x <1 with The components of the attached law f _Qk are, for example, expressed from the following values:

• 0k = (fik>fk> ≠ kk with ¾ = (μ¾, μ ' _¾ ), ¾ = (τ _¾ , τ' _¾ ),

• ¾ (E) the vector of regressors centered with È = (Ε, Ε ^' ), and /? _fe the vector of regression coefficients,

· The matrix Σ _fe , depending on the parameter ifj _k G {0,1}, as follows: Σ _k =

= f (1 ^" 1 ¹ ) ⁵¹

= 1, to choose between a component aligned on the axes E and E ^' {ifj _k = 0) and an oblique component oriented by the line E = E Wk =!) ^■ each comaante f _9k then expresses itself:

where, Χ · μ, σ ² ) represents the Gaussian law of μ and variance σ ² , Κ ₂ {· \ μ, Σ) the bivariate Gaussian law of mean μ GR ² and covariance matrixΣ, where

the a priori law of the parameter _¾ is a Gaussian, the variance r _k is distributed according to a inverse-gamma law, ifj _k follows a Bernoulli-type prior distribution and the a priori for the coefficients of regression coefficients _k is a normal distribution (Gaussian) multivariate of dimension | ^ (E) |, and

by applying a Bayes rule on the expression of f {E, E, y) we obtain the expression

and the probabilistic model

from this probabilistic model and data observed by simulation

(Ei.e ' _j .yi _j ), we estimate the posterior distribution f E, E', y \ E _x , E, y _xx , ..., E _n , E ' _n y _nni ) and the expectation conditional SE, E ') = E (y \ E, E'', E _X , E', y _xl , ..., E _n , E ' _n ,, y _nni ) in order to determine the elements present in the object and their activity.

For the calculation of the posterior, we will use, for example, a Markov chain Monte Carlo approximation scheme (MCMC), for all iteration (t) of the MCMC procedure, we generate a denoised spectral response S { E, E ') ^, for T generations, the posterior distribution of the spectral response is approximated by the set of draws S {E, E') ^ for t = 1, ..., T, and the estimated response expresses himself:

t = i The approximation scheme may include a wafer sampling step using a finite random number κ of components for each iteration and including the following steps:

we introduce latent variables of classification K ^, defined for ί = Ι,.,., η and; = Ι,.,., Η ', such that K _t = k if (E _£ is distributed according to the fth component of the mixture f (E, E y), we define a model for the parameters of the mixture,

for all i≤ n, j <n,

equivalent to the likelihood

f _{(Kll> ■■■>} Knn _'> E _\, E _1} ..., E _n, E _n;

from these probability distributions and by applying the Bayes rule, we calculate the conditional probability density

f (w ₁ , w ₂ , ..., θ, θ ₂ , ... IK _11t ..., K _nn >, E _x , E, ..., E _n , E ' _n >, y _xx , ..., y _nn using a Gibbs sampler which at each iteration (t) successively generates the following samples, for every k,

and at the level of the number of points i in the energy input gate and of points j in the output gate, for all i≤ n, j≤ n,

ι, γιι _> μ ₁ , μ ₂ , ..., μ μ ₂ , ..., τ ₁ , τ ₂ , -, τ ₂ , -, βι, β2>■■■>Ψι>^'Ψ2> ■■ ■> ^w i> ^w 2> ■■ from T iterations of the Gibbs sampler, we obtain an estimate of the meta-model for spectrometry

 T t = i

 It is possible to introduce an auxiliary variable in order to generate only a finite random number κ of components at the iteration (t) while avoiding an arbitrary truncation of the model.

The method may include a step of calculating the posterior standard deviation and credible intervals from the set {S {E, E ') ^}

According to an alternative embodiment, a Gamma a priori (<p _fc , ¾) is introduced on the scale parameters b _T and b ' _T in the distribution of the amplitudes of components leading to a Gamma a posteriori distribution:

with Gamma (a, ô), for> 0,

b ^a _a _ ₁ _ _bx

Gamma (a, fc) 0Ό ^- | ^~ (a) ^ ^

The Gaussian a priori on the regression coefficients β _1ί can be replaced by a priori based on the random draw of P points (Ei, E'i, ..., Ë _P , E'p) from ^ \ ₂ ( _¾ , Σ _¾ ) and we take for ( _1; ..., y _P ) the nearest value of y _tj corresponding to each sampled point where P is larger than the size of the vector /? _fe , we then generate _k ~ Ν (Μ _β , Γ _β ) as a priori law (with Ë _p = (Ë _P , '' _P ))

According to an alternative embodiment, the method generates an extended meta-model denoted Ξ (Ε, Ε ^' , ξ) where ξ is a parameter identified by an integer index and characteristic of a matrix effect.

 The method can estimate the activity of radioelements by performing the following steps:

or Ν _χ the number of issuers selected for the element χ and π _χ for l, ..., N _x , the associated emission probabilities and v _x for 1 = 1, ..., Ν _χ , the corresponding energies,

the response of the radionuclide -ψ _χ {Ε ', ξ), for an observed energy E' and a matrix effect ξ, is defined by

1 = 1

we define ^ = J ₀ Ψχ _{, ξ} ίΕ ¹ ) dE'et π ^* _ιξ = ^ for all / = 1, ..., Ν _χ , we define a radionuclide response normalized by

we check that J ₀ ¾ (£ ") dE '= 1

one expresses the probability density of the ^ith observed in photon energy for i = Ι,.,., η f {E [) = _{Xk JWK r,? k} ^{E p i ^')

k = l where p is a positive parameter for converting the channel index into energy (keV / channel).

Introducing the K variables _t allocation of the i ^th measured at a photon k component of the mixture,

We alternate random draws according to the following conditional laws, for all i≤ n,

for all k

w _lt w ₂ , ... IK _lt ..., K _n

Xk \ Κ-L, ..., K _n , E _lt ..., Ε _η , ξ ₁ , ξ ₂ , ..., p

Κι _> -, K _n , E, -, Ε _η , χ ₁ , χ ₂ , ..., p

as well as

using a finite random number of components at the iteration (t), from T iterations of the Gibbs sampler, we obtain an estimation of the radioelements activities involved in the mixture for all k,

with p a Gaussian a priori centered on μ _ρ and of standard deviation σ _ρ .

For example, the posterior standard deviation of activities is calculated by performing the following steps:

and / or the estimated, deconvolved input spectrum of the system response Other features and advantages of the present invention will appear better on reading the description which follows given by way of illustration and in no way limiting attached to the figures which represent:

 FIG. 1, an example of a gamma spectrum obtained with an HPGe detector,

 • Figure 2, a comparison of uranium spectra obtained using different types of detectors,

 FIG. 3, an example of a spectrum obtained using a plastic scintillator,

 FIG. 4a and FIG. 4B, an example of energy incident at 660 keV and the spectrum measured by a plastic scintillator corresponding to this incident energy,

 FIG. 5, an example of a gamma spectrum to be analyzed and a spectrum reconstructed according to the prior art,

 FIG. 6, the energies reconstructed from the spectrum presented in FIG. 5,

 FIG. 7, an example of a system for implementing the method according to the invention,

 FIG. 8, a block diagram of the steps implemented by the method according to the invention,

 FIGS. 9 and 10, an example of a spectrum of the incident energies obtained by the process according to the invention,

 • Figures 1 1 and 12, a second example of spectrum of incident energies obtained by the implementation of the invention.

Figure 7 is an example of a system for implementing the method according to the invention. The system comprises an object 10 which is capable of spontaneously emitting gamma radiation due to the presence of radioactive elements. The gamma radiation is received by a device 20 comprising a radiation detection module, 21, connected to an acquisition electronics making it possible to obtain a spectrum, a processor 23 or computer unit adapted to perform the steps of FIG. method according to the invention, processing the data of the spectrum obtained, a display module 24 of the results, a storage means 25 for saving the results, for example. The storage means can also store incident energies defined by a user in the context of a given application.

 The method according to the invention is based in particular on a use of all the information present in the spectrum to be analyzed in order to carry out a succession of iterative steps of deconvolution of the incident spectrum. The method uses a continuous projector or meta-model. The response of the spectrometry system will be sought by a simulation technique by calculating whenever necessary the complete energy response of the system for a given configuration of the input spectrum. We will search for the system's answer in the form of a discrete spectrum whose support (incident energies) is a continuous space. This continuous model will make it possible to calculate the response of the spectrometry system for any incident energy E and any observed output energy E ', meta-model denoted S (E, E'). A Bayesian reconstruction step uses the meta-model to determine the radioelements contained in the object to be analyzed.

The method will comprise a first phase during which a numerical simulation of spectrometric responses S, 32 will be carried out for a grid of incident energies Ei, 31, followed by a second phase in which a non-parametric regression is applied. (statistical prediction), 33, on the simulated data. During the second phase, one will inject, for example of the physical knowledge, 34, in the form of a priori to guide the regression in regions of the incident energy where one does not have data simulated. In essence, this makes it possible to indicate that on restricted domains of energy, the signature of the spectral response has either characteristics proportional to the input energy, or constant regardless of the input energy E ,. From the meta-model 35, and the knowledge a priori, it will be possible to go back, by application From a deconvolution algorithm, 36, to the composition of the spectrum 37. The chosen regression model being a non-parametric model, the method does not assume any linear or polynomial model and can thus be adapted to any what non-linearity in the answer.

The method will seek to estimate a continuous function of E ² in

E and E 'from a potentially infinite number of parameters, or components of the model.

To implement the method is used, for example, n the number of points in the input grid (energy E) and n 'number of points in the output grid (£ energies ") ■ ^{Pou r} i = l , -, n and; = 1, ..., ^η ', the characteristic data of the spectral intensities calculated at _i; - for an input energy E _t (user-defined data) and an output energy are available; E (measured data) To illustrate the process according to the invention, it is assumed that the values λ _i have been obtained by means of Monte Carlo simulation software, therefore the data of the spectral intensities λ _i) - are considered as outputs from a Poisson distribution whose intensity will be estimated by means of a non-parametric regression procedure As Xi _j is then a positive or zero quantity, we introduce y _tj = log ( A _{i -} + ε) where 0 <ε «1. Under these conditions, the probability distribution y _{t j} can be approximated for sufficiently large values of (typically>

10) by a Gaussian mean log the _£ - and variance.

The distribution of fish will be expressed as follows: Poisson distribution (l), for n G N,

_ λ ^η _ _{λ η}

/ Fish (λ) () - _\ ^en

The method according to the invention is based on the non-parametric estimation of the quantity f {E, E y), representing the joint probability density of the triplets {E, E y) from the simulated points (Ει, Ε, γ ^ ), in order to deduce the model 5 (£ ^" , £") for all (£, £ ") E KL ² where IL is a continuous space:

where the symbol E (y) represents the expected value of the random variable y, this symbol represents a statistical expectation used later in the description to characterize the intensity of an energy in the spectrum. For the attached law f {E, E, y), one chooses, for example, a mixture by Dirichlet process (DPM) as a priori distribution. This modeling is expressed from a random distribution G such that

G ~ DPO, G ₀ )

where the symbol "~" means "is distributed next" and DP (a, G ₀ ) represents the distribution of a Dirichlet process (Hjort NL, Holmes C, Miiller P, Walker SG, Ghosal S, Lijoi A, Prunster I The YW, Jordan MI, Griffin J, DB Dunson and Quintana F 2010 Bayesian Nonparametrics Cambridge Series in Statistical and Probabilistic Mathematics) with a mean G ₀ and a parameter called concentration. This probability distribution on random probability distributions plays a central role in non parametric Bayesian modeling. A common representation of DP comes from Sethuraman (1994) who expresses the random measure G (-) as:

with w ₁ = V _lt w _k = V _k a) and 5 _U (-) represents the Dirac function located in u. Here, 9 _k ~ G ₀ represent the parameters associated with the k ^th component of G. The nonparametric character is derived from the number of potentially infinite component intervening in the sum characterizing the a priori on G.

The distribution D \ ùch \ and ( ₁ , ₂ , ..., _K ) is defined with x = (x ₁ , x ₂ , ..., x _K ) ^' e R ^K , xi> 0 for 1 <i < K, E -To ¹ <<1, x _K = 1 -Σϊ ^ ₀ ^{1 χ} ί ^and i> ⁰ P ^{or r 1} ί ί K K, fr _v - \ _ Γ (α ₁ + α ₂ + - · + α _κ ) ", a ₁ -l _v a ₂ -l _v Kl

/ Dirichlet ( _ai , a ₂ a _K ) W - _Γ ( _αι ) Γ (α ₂ ) ... Γ (α) ^{Xl 2} - ^Χκ '

 The distribution Beta (a, b), is defined for 0 <x <1:

The random distribution f (E, E, y) will then express: f (E, E ^' , y) = ^ w _k fo _k {E, E ^' , y)

 k = l

where fe _k denotes the k ^th component of f {E, E, y), parameterized by 9 _k .

The components f _Qk in the case of the application of the DPGLM to the characterization of the energy response meta-model, are written by introducing the following additional notations:

• 0 _k = G½ A, with ec _¾ = (μ¾, μ ' _¾ ), r _k = (τ _¾ , τ' _¾ ),

• X _K (É) the vector of regressors centered with Ë = (Ε, Ε ^' ), and _k the vector of regression coefficients. Note that the regression model for the outcome, knowing covariates (Ε, Ε ^' ) can be chosen linear (ie

 T

X _k {E) = (l, E - μ ^ Ε ^' - μ' _¾ ) where v ^T represents the transpose of the vector v) or polynomial in the DPGLM approach (eg X _K (Ë) = (l, E - μ ^ Ε - -μ _¾ ) (£ -μ ' _¾ ), (£ - μ _κ ² , (E - μ' _k )) for a quadratic regression),

• Moreover, we define the matrix Σ _fe , depending on the parameter xfj _k G {0,1}, as follows: Σ _k = R ^. ^ _Τ °) - ^Β Ψ * ^{WITH R} _k = (J °) ^if ^ = ^{0 and}

^R Pk ⁼ ~ r (i ^~ i) ^if ^ fe ^{= The} discrete binary variable thus makes it possible to choose between a component aligned on the axes E and E {ifj _k = 0) and an oblique component oriented by the line E = E ( 4> _k = 1).

Equipped with these notations, each component of the random distribution is expressed as: fe _k (E, E ^' , y) = N ₂

where, Χ · μ, σ ² ) represents the Gaussian law of μ and variance σ ² , Κ ₂ {· \ μ, Τ) the bivariate Gaussian law of mean β GE ² and covariance matrixΣ.

A multivariate Gaussian distribution Ν _ρ (μ, Ω) is defined for XGR ^p ,

 1

where | 4 | represents here the determinant of the matrix Note: this definition also covers the univariate case when p = 1.

The pr / or / G ₀ measurement of the Dirichlet process is chosen as follows. The a priori law of the parameter _¾ is a Gaussian, the variance r _k is distributed according to a inverse-gamma law, ifj _k follows a Bernoulli-type prior distribution and the a priori for the coefficients of regression coefficients β _ίι is a law normal (Gaussian) multivariate dimension

Using the Bayes rule, we deduce from the expression of f (E, E ^' , y) that: f (y ï) = Σ _{y Λ} "(r | ol ^■ * ■ (¾) ■" - *'">) ■ Finally, by calculating the conditional expectation E E ^' ) from the expression of E ^' ), the expression of the S (E, E ^' ) spectral response model is obtained:

This last expression makes it possible to account for the a priori corresponding to the nonparametric Bayesian approach for the denoised and interpolated spectral response. The heteroscedastic and non-Gaussian characters in all (Ε, Ε ^' ) are taken into account in the model of E ^' ) by mixing Gaussian distributions M (| /? [· X _K (E), _E

From this choice of probabilistic model and the data observed by physical simulation (E _£ the process will estimate the posterior distribution f {E, ..., Ε _η , Ε '_η; γ _ηη ·) and the conditional expectation of this posterior law: S (E, E ^' ) = E yE, E, E _x , E, y _n , ... , E _n , E '_n; y _nn <representing the intensity in the spectrum of an output energy E, for any input energy E.

 For the exact calculation of the posterior distribution the method uses, for example, a Gibbs sampler Markov Chain Monte Carlo (MCMC) scheme to generate samples of the target law. In order to allow the MCMC sampling procedure for objects of infinite dimension (infinity of component of the DP), the method adopts a so-called slice sampling approach according to Kalli et al. (Kalli M, Griffin JE and Walker SG 201 Slice Sampling Mixture Models, Statistics and Computing, 21, 93-105), where only a finite random number K of components is required at each iteration.

For any iteration (t) of the MCMC procedure, it is thus possible For T generations, the posterior distribution of the spectral response is approximated by the set of draws S (E, E ^' ) ^ for t = 1, ..., T, and the estimated response (posterior mean) s' Express :

 It is also possible to calculate in a similar way the posterior standard deviation and the credible intervals from the collection.

An example for performing the Gibbs sampling procedure will now be given. The generative model for DPM mixing parameters is given by:

V _X , V, ... ~ Beta (l, a)

be w ₁ = V _lt w ₂ = y ₂ (l - Vi), ■■■

τι ¹ > ^τ 2 ¹ > ■■■ ~ Gamma (a _T , b _T )

τ ^, τ ^ ¹ , Gamma (a _T , δ _T )

β ₁ , β ₂ , ... ~ Μ _ι ^ _{) 1} {Μ _β , Γ _β )

-φ ₁ , -φ ₂ > - ~ Bernoulli (p)

let θ ₁ = (β ^ τ ^ β ^ -ψ ^, θ ₂ = (μ ₂ , τ ₂ , β ₂ , ψ ₂ ), ...

This generative model is entirely equivalent to the probability density a priori f (w ₁ , w ₂ , ..., θ _χ , θ ₂ , ...), and corresponds to a priori of the Dirichlet process type for the measurement. random G (-) = Σ _{¾ = 1} w _fe ¾ (·) ~ ÛP a, G _Q ).

The description of the modeling of the DPM parameters is completed by a generative model of the observed data resulting from the physical simulation. For this, we introduce latent variables called K ^, defined for i = l, ..., n and j = l, ..., n ', such that Ki = k if (Ei.E'j.yij) is distributed according to the k ^th component of the mixture f (E, E', y).

 Loop on all n input energies E and n 'energies E' of the output grid

 For all i≤ n, j≤ n,

with μ _Κί . = (μ _Κί μ ' _κ .), Σ * _y =, τ' * _.. "^ y _ij \ Κ _φ Ε ₀ Ε), θ ,, θ ₂ , ... ~ M (y _y I / ¾ ^■ X _Ki] {E _it E) ^ ^K ^ ⁱ ). This generative model is entirely equivalent to the likelihood:

f (.Kn _> ■■■, Knn '> E ±, E i, ..., E _n , E _n ;

From these probability distributions and by applying the Bayes rule, we try to calculate the conditional probability density θ ₂ , ... IK ₁ , ..., K _nn ', E _x , E _x , ..., E _n , E _n ', yii, -, y _n n- which represents the distribution of weights and component parameters conditionally to observations and classification variables.

 This inference is carried out by means of a Gibbs sampler which at each iteration (t) of the algorithm successively generates the following samples, for all k,

w _1} w ₂ , ... | K _11} ..., Κ _ηη ·

I fej k>Kn>^■■■> Knn '> E, E -L, ..., E _n, E _η ·

^T k \ E _x , ..., E _n , E _η · τ k \ fej k '^ k'Kll> -, Κηη'Έί, Ε ..., E _n , E _n ^> fc fc E _x , ..., E _n , E _n ^> k \ fej k ' ^Tk ' ..., E _n , E _n ', yii, -, y _n n' and at the level of the number of points i in the energy input grid and points j in the output grid, for all i≤ n, j≤ n,

ΐ μ _χ , μ ₂ , ..., μ _χ , μ _v ..., τ _χ , τ ₂ , ..., τ _χ , τ ₂ , ..., β _χ , β ₂ , ..., ψ ₁ , ψ ₂ , ..., w ₁ , w ₂ , ...

The nonparametric behavior of the process is characterized by the potentially infinite collection of parameters w _k and 9 _k . In order to reduce the computational difficulty induced by the infinite number of parameters, it is possible to use truncation of the DP at a sufficiently high number of components (Ishwaran and James 2001).

An alternative solution proposed by (Kalli et al., 2011) is to introduce an auxiliary variable that makes it possible to generate a finite random number κ of components at the iteration (t) while avoiding an arbitrary truncation of the model. In this case, the different stages of random generation of the algorithm will be as follows: Random initialization.

- Generate for 1 <i≤ n, j≤ n ', u _i; ~ Uniform (0, 1 / n.n '), and ask u ^* = min _u ({y})

- Generate l ^ -BetaCl.a); Ask w ₁ = V ₁ and r _x = 1 - V.

- For k> 1, generate V _k ~ Beta (l, c), put w _k = V _k r _k _ ₁ , ask r _k =

r ^ Ul ^ as r _k ≥u ^* .

 - Assign to K * the maximum value of k.

- Generate 6 _k for 1 <k≤ K *, from the prior distributions.

At the iteration t = l, ..., T,

· Generate

w ₂ , ...) for ί≤ η, j≤ η '

- Calculate for 1 <k≤ κ ^* ,

- X _k = (w _k > u _y ) max (w _k ^ / n.nD: Ar ₂ (£ _i; | ^

 1 (A) is the indicator function: 1 (A) = 1 if A is true and 1 (A) = 0

 if not.

- Generate Κ _ίΓ -J-Σ i

Lk = \ ^A k

 • Reorder component labels in order

of appearance in the generation of Κ ^. Assign to κ _η the number of distinct values of Κ ^. For all k κ _η , let _k = # {κ ^ = k] be the number of observations y _i; - assigned to component k.

• Generate (w _1; w ₂ ..., u, ..., u _nn \ K _llt ..., K _nn

- Generate (w _lt w ₂ , ..., w _Kn , r _Kn ) ~ Dirichlet (n ₁ , n ₂ , ..., n _Kn , a)

- Generate for i <n, j≤ n ', it _i; ~ Uniform (0, min (w _{/ c} , 1 / n.n '))

- Put u ^* = min _u ({tt _i; })

- For k> κ _η , generate l ^ -BetaC ^ a), put w _k = V _k r _k _ ₁ , ask r _k = r _{fe - 1} (l - V _k ), as long as r _k ≥ u ^* .

 - Assign to K * the maximum value of k.

• Generate ..., E _n , E ' _n ) for k≤ K *,

 Calculate

 - Generate ν ~ Uniform (0,1)

- If <-, put ip _k = 0 else ask x / j _k = 1

 Vo + Vi

• Generate (τ _¾ | μ ^ μ '^ ψ ^ Κ ^, K _nn >, E _Xl E, ..., E _n , E' _n ) for / c </ ^* ,

- Si = 0 T _fc ¹ ~ Gamma

- = l -μ ',)

· Generate (τ ^ \ μ ^ μ ' _{] ί} , ψ _{1 {} , Κ ₁₁ , ..., Κ _ηη ; Ε ₁ , Ε' ₁ , ..., Ε _η , Ε ' _η ·) for / c < κ ^*

 - Si = 0

-

• Generate (ji _k \ T _k , r ' _k , ip _k , K _llt ..., K _nn ; E _x , E, ..., E _n , E' _n ) for k≤ κ

- Calculate = (Σ _{{ij K j = k]} E _t, Σ _{[ij: K j = k...}} E) - Install Σ _k = R _k ^ _T °) - _k with _k = (JJ). if tf _fc = 0 and - Calculate ν _μ = (Γ _μ ^_1 + n _k Σ _fe ^_1 ) ^~

- Calculate Μ _μ¾ = ν _μκ . (Γ _μ ^_1 .τ ^η _μ + Σ _fc ^_1 -¾

- Generate μ ^ - ₂ {Μ _μκ , ν _μκ )

• Generate μ ^ μ ' _ίι > ^τ ^ ^τ >' Φ ^^ ιι _> ..., Κ _ηη <, Ε _χ , Ε, ..., Ε _η , Ε ' _η <, γ _χχ , ..., y _nn ) for k≤ K *,

 - To pose

- To pose

 - Generate ^ - ^^ (Μ ^, Γ ^).

For this random draw, we approximate the observation variance of yt _j (itself approximated to a Gaussian random variable) σ? - = _e -fiI-Xk (È ⁱ j) to the value of the observation Poisson A _{i ;} - = e ^{~ yi)} corresponding. This approximation, which greatly simplifies the inference, is all the more valid if A _{£ y} - is high.

• Calculate at each iteration (t) for any choice of prediction grid Ë = (E, E '),

The denoised spectral response can, optionally, share or not the same grid (E _£ , E) as the initial physical Monte-Carlo simulation. Indeed, the user can choose any point (£, £ ") of interest, in particular, it is possible to interpolate between the points of the initial grid.

Finally, from T iterations of the Gibbs sampler, we obtain an estimate of the meta-model for spectrometry τ

 SOW E ') * -} S. E ^

 t = i

 It is then possible to calculate the posterior standard deviation of the spectral response estimator

t = i

 The hyper parameters can be considered fixed in the estimation of the response of the spectrometric system according to prior physical knowledge on the detector. On the other hand, these can be estimated by means of an additional level of hierarchy which allows them to be assigned a so-called vague a priori.

• As an example, it is possible to place a Gamma a priori (<p _fc , ¾) on the scale parameters b _T and b in the distribution of component amplitudes. This choice leads to a Gamma posterior distribution.

Example:

 • It may be relevant to estimate the Dirichlet process concentration parameter a following the approach of Escobar and West (1995).

"The parameter p of the Bernoulli law intervening in the mixture model of the component type ifj _k can itself be estimated via a hyperprior BetaC ^, ^) - The posterior distribution for p then follows expression:

p ~ Beta + # {k κ _η ·. i _k = 0}, _πι + # {k κ _η ·. i _k = 1})

• Note also that it may be effective to replace the Gaussian priori on the regression coefficients _k by a so-called empirical a priori based on the random draw of P points {^ Έ '^ ..., Ë _P , Ë' _P ) from · 7 \ Γ ₂ ( _¾ , Σ _¾ ) and take for ( _1; ..., y _P ) the nearest value of y _tj corresponding to each sampled point. P must be larger than the size of the vector k. We then generate β _1ί ~ Ν (Μ _β , Γ _β ) as a priori law (with

Ep (^ p E)) ■

The expression for the posterior distribution of _k follows directly from it.

 At the end of these different steps, the method has a continuous meta-model characterizing any type of detector for nuclear spectrometry, the meta-model taking into account all the physical interactions involved in the spectral response for an energy. given input.

 The meta-model thus generated could be injected into a nuclear spectrum convolution method using non-parametric probabilistic modeling such as that proposed by Barat et al 2007 (Barat, E., Dautremer, T., Montagu, T., "Nonparametric Bayesian Inference in Nuclear Spectrometry, "Nuclear Science Symposium Conference Record, 2007. NSS '07 .I EEE, vol.1, no., pp.880-887, Oct. 26, 2007-Nov. 3, 2007), to determine the radioelements present in the spectrum and their activity.

 Depending on the environment in which the object to be analyzed is located, spectral signature shape variations due to the presence of matrix effects at the input source may occur. The meta-model defined according to the invention can take into account several spectral responses estimated as has been explained above. For this, a solution is for example to use different shields of different thicknesses and different materials, to be simulated and intervene in the meta-model.

Another solution is to conserve a Polya tree to model the diffuse interactions leading to continuous bottoms that may not be modeled in the response. Assuming that the meta-model takes into account a collection of matrix effects as previously described, the desired input spectrum then consists only of discrete peaks. It is then possible to substitute a priori Dirichlet process on mono-energetic lines, a priori Dirichlet process on radioelements of a database of radionuclides. The spectrum sought is then constituted by a discrete sum of elements of the periodic table whose number of components is unlimited. It then becomes possible to assign to each radioelement a prior probability of presence in the sample analyzed. The meta-model is again used to characterize the spectral response of the detector of the considered element. This response is directly determined by the discrete sum of the individual spectral responses corresponding to each monoenergetic line considered. The use of a radionuclide database makes it possible, in particular, to consider the energy calibration of the spectrum observed as uncertain and to estimate it from the data. For each proposed value of keV / channel it is necessary to evaluate the spectral response of the system for each energy of the output grid. This operation uses the continuous meta-model, eliminating the need for interpolations when using a discrete model.

In the case where the meta-model takes into account a collection of matrix effects, to directly estimate the activity of radionuclides, we use an extended meta-model noted Ξ (Ε, Ε ^' , ξ) where ξ represents a parameter characteristic of the matrix effect associated with the source. The collection of matrix effects is discrete of size M and each element ξ is identified by an integer index.

 The activity estimation algorithm is based on a Markov Chain Monte Carlo (MCMC) approach and more particularly on a Gibbs sampler, an example of which is detailed below.

Before explaining the description of the algorithm, various parameters are set. The spectrum is considered a mixture by Dirichlet process of different standardized radionuclide responses. An example of construction of this mixture is given below. Note any radionuclide of a given table of radioelements. This table being discrete, we will find in the following each element by an integer index. We denote N _x the number of emitters selected for the element χ considered and n _x for l = 1, ..., N _X , the associated emission probabilities as well as v _xX for l = 1, ..., N _X , the corresponding energies. We define the response of the radionuclide ψ _χ (Ε ', ξ), for an observed energy "and a matrix effect ξ, by:

we define Α _{χ ξ} = / _ο ^° ° tl = 1, ..., N _X. We define

then, a radionuclide response normalized by:

 1 = 1

we check that J ₀ ^° ° i / ^* _ξ (£ ") dE '= 1.

The probability density of the ^1st photon observed in the energy channel then expresses for i = 1, ..., n:

k = l

where p is a positive parameter for converting the channel index into energy (keV / channel). We suggest for p a gaussian a priori centered on μ _ρ and of standard deviation σ _ρ :

 This Dirichlet process mixture is based on the probability measurement F generated by a Dirichlet process:

F ~ DP (a, F ° XF ^) Ρχ {χ) ^s j (x) represents the probability distribution a priori to observe an element χ where / is the number of radionuclides of the base selected for the analysis and p °> 0. F ° is the discrete uniform prior law on the matrix effects collection, and a is the (positive) process concentration parameter:

fc = l

with w ₁ = V _x , w _k = and _U (-) represents the Dirac function located in u.

Before describing the different random draws for exploring the posterior, it remains to introduce the K variables _t allocation of the ^ith photon observed at k component of the mixture.

 The Gibbs sampler consists in alternating random draws according to the following conditional laws, for all i≤ n,

..., χ _ί , χ ₂ , ..., ξ _ί , ξ ₂ , ..., ρ

and for all k

w ₁ , w ₂ , ... | K _lt ..., K _n

X _k \ Κ- _L , ..., K _n , E _lt ..., Ε _η , ξ ₁ , ξ ₂ , ..., p

Κι _> -, K _n , E, -, Ε _η , χ ₁ , χ ₂ , ..., p

as well as

p | K _lt ..., K _n , E _lt ■■ ·, Ε _η , χ ₁ , χ ₂ , ..., ξι, ξ ₂ , ■■■

As for the determination of the meta-model, in order to minimize the number of computations, the approach proposed by (Kalli et al., 2011) and the introduction of an auxiliary variable _£ which makes it possible to generate only a finite random number κ components at iteration (t) while avoiding arbitrary truncation of the model are used. Under these assumptions, the various stages of random generation of the algorithm are detailed below.

 Random initialization.

- Generate for 1 <i <n, ιι _έ ~ 11ηιτοπτι (0, 1 / n), and ask u ^* = min _u ({uj) - generate l ^ -BetaCl, a); ask w ₁ = V ₁ and r _x = l-V.

- For k> 1, generate V _k ~ Beta (l, a), put w _fe = V _fe r _fe _ ₁ , and ask r _k = r _fc _ _! (l - Vfc), as r _fe> u ^*.

 - Assign to K * the maximum value of k.

- Generate χ _κ for 1 </ c <K *, from the prior distribution F °.

- Generate ξ _{} ι} for 1 <k≤ K *, from the prior distribution F °.

- Generate p ~

At the iteration t = l, ..., T,

• Generate (K _I \ E _i ^' , w ₁ , w ₂ , ..., x _lt x ₂ , ..., ξ ₁ , ξ ₂ , ..., p) for i <n,

 - Calculate for 1 <k≤ K *,

- At _k = l (w _k> u _t) max (w _fe, l / n) ^ ^* _¾¾ (£ ^"-) or 1 (A) is the indicator function 1 (A) = 1 if A is true and 1 (A) = 0 otherwise.

- Generate K _T ~ -J-Σ i S _K (KÔ,

• Reorder the component labels in their order of appearance in the K _t generation. Assign to κ _η the number of distinct values of K _t . For all k <κ _η , let _k = # {K _t = k} be the number of observations assigned to component k.

• Generate (w ₁ , w ₂ ..., u _n \ K, ..., K _n )

- Generate (w _lt w ₂ , ..., w _Kn , r _K J ~ Dirichlet (rr _lJ n ₂ , ..., n _Kn , a),

- Generate for i < 1 / n)),

- Put u ^* = min _u ({uj),

- For k> K _n , generate l ^ -BetaC ^ a), put w _k = V _k r _k _ ₁ , ask r _k = r _{fe - 1} (l - V _k ), as long as r _k ≥ u ^* ,

 - Assign to K * the maximum value of k.

· Generate Κ ₁ , ..., Κ _η , Ε ₁ ^' , ..., Ε _η , ξ ₁ , ξ ₂ , ..., ρ) for k≤ K *,

- Calculate for everything; = 1, ...,],

- Generate ^ - -Σ ^J _{j = 1} Vk (x _k ).

Lj _{= 1} Vk, j Generate (ξ \ Κ _1} ..., Κ _η , Ε, ..., Ε _η , χ ₁ , χ ₂ , ..., p) for k

 - Calculate for all m =

• Generate p \ ..., K _n , E, ..., Ε _η , χ ₁ , χ ₂ , ... ,, ξ ₁ , ξ ₂ , ...) by a Metropolis-Hastings stage,

- Generate a random walk of standard deviation ε _ρ around p

- Calculate

 - Generate a uniform variable v ~ Uniform (0,1).

- If v≤ min l, g _MH ) assign p = p ^* , otherwise leave p unchanged.

Finally, from T iterations of the Gibbs sampler, we obtain an estimate of the activities of the radioelements involved in the mixture, for all k,

In addition, it is also immediate to calculate the posterior standard deviation of activities:

 It is also possible to obtain the estimated, deconvolved input spectrum of the system response

 N (t)

 T% k

^ ® ^¾ r ^wit) ¾ ^ ¾ ^(f)

t = i 1 = 1 ^K In the proposed version the algorithm allows the estimation of keV / channel conversion gain from only the data as well as the determination of a (potentially different) matrix effect for each element involved in the mixture. It also makes it possible to set a prior probability of the occurrence of radionuclides in the context under consideration and to obtain directly the activities of the constituent elements of the mixture as well as the associated uncertainties.

 Figures 9 and 10 illustrate the spectrum of the energies obtained from the spectrum analysis of Figure 10 using the method according to the invention. The dotted curve represents the simulated spectrum, the solid curve the reconstructed spectrum.

Figures 1 1 and 12 illustrate a second example of spectrum. Figure 11 shows the spectrum of incident energies, with an emission region of the first peak at ⁶⁰ Co, 1173.2 keV. Figure 12 the same spectrum with a zoom on the region of interest.

Claims

1 - Method for determining the nature of the radioelements present in an object and their activity, characterized in that it comprises at least the following steps:

 A first phase (32) of numerical simulation of spectrometric responses for a set of incident energy E and a set of measured output energy E ', in order to obtain a set of simulated data,

A second non-parametric regression phase (33) on the simulated data, non-parametric estimation of the quantity representing the joint probability of the triplets (E, E ', y) from simulated points (Ei, Ei', yij) in order to to deduce a meta-model S (E, E ') for any energy pair (E, E') on a continuous function,

 · We use n a number of points in the input grid, energies

E, and n 'a number of points in the output grid, energies £ ", for i = l, ..., n and j = l, ..., n', we have the characteristic data of the spectral intensities calculated at _{£ y} - for an input energy _£ and an output energy £ ", -,

A nonparametric estimation method of the quantity f {E, E ', y), representing the joint probability density of the triplets (E, E', y) from the simulated points (Ει, Ε, γ) is used ^), and we deduce a model S (E, E ') for all (E, E') ER ^{2 where} M is a continuous space:

where the symbol E (y) represents the expected value of the random variable y, y is deduced from the calculated spectral intensities λ ·, • from the meta-model S (E, E '), (36) the determination of the nature and activity of the radioelements present in the object.

2 - Process according to claim 1 characterized in that one determines the values À _{£ y} - using Monte-Carlo software, the data À _{£ y} - being considered as outputs from a Poisson distribution of which intensity is estimated by means of a non-parametric regression procedure and we introduce y _tj = log (A _i; - + ε) where 0 <ε «1, then

we approximate the probability distribution of y _tj for sufficiently large values of> 10 by a Gaussian law of mean log and variance, for the joined law f (E, E ', y), we choose a mixture by Dirichlet process (DPM) as a priori distribution and we express the random distribution in a sum over the infinite of f _Qk components of f (E, E, y), f {E, E ^' , y) = ^ w _k fo _k { E, E ^' , y)

 k = l

parameterized by 9 _k the parameter associated with the k ^th component of G random measure defined by

G (-) = Σw _fc ¾ (·)

 fc = l

with w ₁ = V _x , w _k = Ι ^ Π? = ι (1 - ½) such that V _k ~ Be \ .a (l, a) and 5 _U (-) represents the Dirac function localized in u and Beta (a, ô), for 0 <x <1 with

3 - Process according to claim 2 characterized in that the components of the attached law f _Qk are expressed from the following values:

• 0 _k = G½ A, with ec _¾ = (μ¾, μ ' _¾ ), r _k = (τ _¾ , τ' _¾ ), • X _K (Ê) the vector of regressors centered with = = (Ε, Ε ^' ), and _k the vector of regression coefficients,

• the matrix Σ _fe , depending on the parameter ifj _k G {0,1}, as follows: Σ _k = xfj _k = 1, to choose between a component aligned on the axes E and E (4> _k = 0) and an oblique component oriented by the line E = E W> _k =!) ^■

each component _{θ¾ is} expressed then:

where, ^ \ (· | μ, σ ² ) represents the Gaussian law of μ and variance σ ² , Κ ₂ {· \ μ, Τ) the bivariate Gaussian law of mean μ GE ² and covariance matrixΣ

or

the a priori law of the parameter _¾ is a Gaussian, the variance r _k is distributed according to a inverse-gamma law, ifj _k follows a Bernoulli-type prior distribution and the a priori for the coefficients of regression coefficients _k is a normal distribution (Gaussian) multivariate of dimension and

by applying a Bayes rule on the expression of f {E, E, y) we obtain the expression

and the probabilistic model From this probabilistic model and the data observed by simulation {EuE ' _j .yi _j ), we estimate the posterior distribution f E, E', y \ E ₁ , E ' ₁ , y ₁₁ , ..., E _n , E ' _n y _nni ) and the conditional expectation SE, E') = E y \ E, E '', E _X , E ', y _xl , ..., E _n , E' _n ,, y _nni ) to determine the elements present in the object and their activity. 4 - Process according to claim 3 characterized in that for the calculation of the posterior, using a Markov chain Monte Carlo approximation scheme (MCMC),

for any MCMC iterative procedure, a denoised spectral response is generated

for T generations, the posterior distribution of the spectral response is approximated by the set of draws S ^ E. E ") ^ for t = 1, ..., T, and the estimated response is:

 T

S {E, E ^' ) ~ - S (E, E ^' )

 t = i

5 - Process according to claim 4 characterized in that the approximation scheme comprises a step of sampling per slice using a finite random number κ of components for each iteration and in that it comprises the following steps:

we introduce latent variables of classification K ^ defined for ί =

Ι,.,. , η and; = Ι,.,. , η ', such that K _t = k si is distributed according to the fth component of the mixture f (E, E y), we define a model for the parameters of the mixture,

for all i≤ n, j <n,

equivalent to the likelihood (# 11, ■■■, Knn _'> E ±, E i, ..., E _n, E _{n; yn> ->} y _n n \ w _1, w _2, ..., _Θ χ e ₂ , ...) from these probability distributions and by applying the Bayes rule, we calculate the conditional probability density

f (w ₁ , w ₂ , ..., θ, θ ₂ , ... I end, ..., K _nn ; E _x , E, ..., E _n , E ' _n >, y _xx , ..., there _nn

using a Gibbs sampler which at each iteration (t) successively generates the following samples, for all k,

w _1} w ₂ , ... I # 11, ..., # _ηη ·

I fej k>#11> Knn '>Ei> E i, ^■■■ , E _n , E _n ^>

Tk \ _x , ..., E _n , E _η · τ k \ _x , ..., E _n , E _n ^> fej E _x , ..., E _n , E _n ^> k \ fej k ' ^Tk ..., E _n , E _n ', n, ..., y _n n' and at the level of the number of points i in the input grid of energy and points j in the output grid, for all i≤ n, j≤ n,

ΐ μ _χ , μ ₂ , ..., μ _χ , μ _v ..., τ _χ , τ ₂ , ..., τ _χ , τ ₂ , ..., β _χ , β ₂ , ..., ψ ₁ , ψ ₂ , ..., w ₁ , w ₂ , ... from T iterations of the Gibbs sampler, we obtain an estimate of the meta-model for spectrometry: 6 - Process according to claim 5 characterized in that one introduces an auxiliary variable in order to generate a finite random number κ of components at the iteration (t) while avoiding arbitrary truncation of the model. 7 - Method according to one of claims 4 to 6 characterized in that it comprises a step of calculating the standard deviation a posteriori and credible intervals from the set {S {E, E ') ^ }

8 - Method according to one of claims 1 to 7 characterized in that one introduces a prior Gamma (<p _fc , ¾) on the scale parameters b _T and b ' _T in the distribution of component amplitudes leading to a Gamma posterior distribution:

with Gamma (a, ô), for> 0,

/ Gamma (a, 6) (· * ·) ^- Γ (α) ^

9 - Process according to claim 3 characterized in that the Gaussian a priori is replaced on the regression coefficients _k by a priori based on the random draw of P points {^ Ε '^ ..., Ë _P , È ' _P ) from J \ f ₂ Q½, Σ _fe ) and we take for ( _1; ..., y _P ) the nearest value of y _tj corresponding to each sampled point where P is greater than the size of the vector _k , we then generate _k ~ Ν (Μ _β , Γ _β ) as a priori law (with Ë = (Ë _p , Ë ' _p ))

10 - Method according to one of the preceding claims characterized in that one generates an extended meta-model noted Ξ (Ε, Ε ^' , ξ) where ξ is a parameter identified by an integer index and characteristic of an effect of matrix. 1 1 - Method according to one of the preceding claims characterized in that the activity of the radioelements is estimated by performing the following steps:

Ν _χ the number of issuers selected for the element χ and π _χ for 1 = 1, ..., N _X , the associated emission probabilities and v _x for l = 1, ..., Ν _χ , the corresponding energies,

the response of the radionuclide -ψ _χ {Ε ', ξ), for an observed energy E' and a matrix effect ξ, is defined by

we define | ^ = J ₀ ^°° ^^ (£ ") d £" and π * _{ι ξ} = for all / = 1, ..., N _X , we define a radionuclide response normalized by

 1 = 1

we check that J ₀ ¾ (£ ") dE '= 1

one expresses the probability density of the ^ith observed in the photon energy channel for i = Ι,.,. , Η f {E [) = _{Xk JWK r,? K} ^{E p i ')

 k = l

where p is a positive parameter for converting the channel index into energy (keV / channel).

Introducing the K variables _t allocation of the i ^th measured at a photon k component of the mixture,

 We alternate random draws according to the following conditional laws, for all i≤ n,

for all k

w ₁ , w ₂ , ... | K _lt ..., K _n Xk \ KL, ..., K _n , E _lt ..., Ε _η , ξ ₁ , ξ ₂ , ..., p

Κι _> -, K _n , E, -, Ε _η , χ ₁ , χ ₂ , ..., p

 as well as

 using a finite random number of components at the iteration (t), from T iterations of the Gibbs sampler we obtain an estimation of the radioelements activities involved in the mixture for all

5k,

 t = l

with p a Gaussian a priori centered on μ _ρ and of standard deviation σ _ρ .

12 - Process according to claim 11 characterized in that one calculates o the standard deviation a posteriori activities

and / or the estimated, deconvolved input spectrum of the system response

5