FR 3034221 A1 20160930 FR 1552415 A 20150324 1 PROCEDE ET DISPOSITIF POUR DETECTER DES RADIOELEMENTS L'objet de l'invention concerne un procédé et un dispositif pour détecter des radioéléments contenus dans un objet, déterminer leur nature et leur activité. Elle permet notamment de tenir compte d'effets parasites provenant de l'environnement dans lequel se trouve l'objet à analyser. La spectrométrie gamma est une technique permettant de détecter les rayonnements gamma émis spontanément par des radioéléments. Cette technique permet d'identifier les radioéléments émetteurs gamma et de mesurer leur activité. Pour atteindre des objectifs qualitatif et quantitatif, le traitement d'un spectre gamma se base sur l'analyse des pics présents dans le spectre. Ces derniers correspondent au rayonnement gamma émis par un radioélément spécifique, par exemple, le 137Cs émettant un gamma à 661,5 keV ou le 60Co émettant un gamma à 1173,2 keV et un gamma à 1332,5 keV. La présence d'un pic dans un spectre gamma signifie que l'énergie du rayonnement gamma incident a été entièrement absorbée dans le détecteur, suite à une ou plusieurs interactions, effet photoélectrique, effet Compton ou production de paires, ayant eu lieu au sein de ce dernier. Les détecteurs utilisés pour réaliser des mesures de spectrométrie gamma peuvent être séparés en deux catégories : les scintillateurs et les semi-conducteurs. La qualité d'un spectre gamma peut être évaluée grâce à deux paramètres propres au détecteur utilisé: la résolution en énergie, c'est-à-dire la capacité à séparer deux pics proches en énergie définie comme étant la largeur à mi-hauteur d'un pic présent dans le spectre à une énergie donnée, et l'efficacité de détection, plus le matériau sera dense, plus la probabilité d'absorber l'intégralité du rayonnement gamma incident est importante. Ces deux paramètres varient de manière significative en fonction du détecteur utilisé. Les scintillateurs inorganiques (Nal, BGO) bénéficient d'une résolution dégradée par rapport aux détecteurs semi-conducteurs (CZT ou 3034221 2 HPGe). Néanmoins, ils présentent généralement une plus grande efficacité de détection. En effet, en raison de leur faible coût, il est possible de les fabriquer dans de grandes dimensions. La figure 1 montre, dans un diagramme Energie de canal, coup par canal, un exemple de spectre S(I) 5 obtenu à l'aide d'un détecteur HPGe et d'une source de 137Cs émettant un rayonnement gamma à 661,5 keV. La figure 2 compare des spectres S(HPGe), S(CZT), S(Nal) obtenus à l'aide d'un échantillon d'uranium et de différents types de détecteurs, présentant des performances différentes en termes de résolution en énergie. 1 0 Il est aussi connu d'utiliser des scintillateurs plastiques pour détecter des rayonnements gamma. Le coût de fabrication de tels scintillateurs étant assez faible, ces derniers peuvent être fabriqués dans de grandes dimensions pour différents types d'application. Du fait de leurs caractéristiques intrinsèques et en dépit de leurs dimensions très importantes 15 par rapport à d'autres détecteurs, la probabilité qu'un photon incident dépose toute son énergie au sein d'un détecteur par effet photoélectrique est extrêmement faible. La figure 3 montre un spectre S(II) obtenu à l'aide d'un scintillateur plastique de type EJ200, en présence de sources 60Co et 137Cs. Les formes caractéristiques de ces spectres sont dues aux fronts Compton 20 faisant suite aux interactions au sein du détecteur. On peut ainsi noter l'absence des pics photoélectriques caractéristiques de la présence de 137Cs et de 60Co. Pour cette raison, ces détecteurs sont généralement utilisés pour réaliser des mesures de comptage total, détection de tous les rayonnements gamma interagissant au sein du détecteur sans information sur leur énergie. 25 Une approche traditionnelle consistant à mesurer l'aire nette des différents pics d'intérêt n'est donc pas envisageable et une solution alternative doit être employée si l'on souhaite procéder à l'identification des radioéléments présents dans le spectre. Une solution pour traiter un spectre gamma, identifier les 30 radioéléments émetteurs et mesurer l'activité de ces derniers, consiste à analyser le spectre dans son ensemble et à ne pas se restreindre 3034221 3 uniquement à l'analyse des pics photoélectriques. L'approche suivie consiste à exprimer l'analyse du spectre sous la forme matricielle suivante : S=H.A, avec S : une matrice signal, matrice de dimensions (nbre_canaux, 1), où la 5 variable nbre_canaux est égale au nombre de canaux du spectre gamma, A : une matrice activité, matrice de dimensions (nbre_energies_incidentes, 1). La variable nbre_energies_incidentes correspond au nombre d'énergies incidentes défini par l'utilisateur et pris en compte dans la reconstruction. Ces énergies incidentes correspondent au nombre d'éléments de volume 10 élémentaires d'un problème classique de tomographie d'émission, H : un projecteur du problème, matrice de dimensions (nbre_canaux, nbre_energies_incidentes). Cette matrice comprend l'ensemble des efficacités de détection intervenant dans le problème. A titre d'exemple, un élément hij de cette matrice correspond à la probabilité qu'un photon 15 d'énergie incidente j soit détecté dans le canal i du spectre gamma. La matrice signal S correspond au résultat de la mesure, au spectre gamma à analyser, la matrice projecteur est généralement calculée par simulation, la matrice A correspond au résultat de la reconstruction. La figure 4A et la figure 4B illustrent un exemple de notion de 20 projecteur H, pour une énergie incidente de 660keV, figure 4A. Le spectre S(IV) de la figure 4B correspond à l'énergie incidente correspondant à une colonne entière du projecteur H. Le projecteur contient donc la réponse spectrale du détecteur pour chaque énergie incidente prise en compte dans le problème et définie par l'utilisateur. C'est sur cette grille d'énergies 25 incidentes que sera effectuée l'étape de reconstruction. Différentes méthodes permettent de prendre en compte l'ensemble du spectre gamma afin de procéder à une analyse qualitative et quantitative, permettant l'identification des radioéléments présents dans le spectre et l'activité du radioélément présent dans le spectre des énergies 30 incidentes. Cette méthodologie d'analyse est illustrée sur les figures 5 et 6. La figure 5 représente le spectre gamma à analyser, correspondant à une 3034221 4 simulation de la réponse d'un scintillateur plastique de type EJ200, en présence de trois sources gamma (241Am, 137C s, Rn S --Co). Le spectre reconstruit Sr à partir des résultats de l'analyse et des énergies incidentes est représenté en traits pointillés sur la figure. La figure 6 illustre le résultat 5 d'une reconstruction effectuée à l'aide d'un algorithme de type ML-EM, les énergies incidentes reconstruites à partir du spectre de la figure 5. Les méthodes d'analyse traditionnelles malgré tout leur intérêt et avantages présentent quelques limitations. Une des limitations intrinsèques est liée au caractère discret de la grille choisie pour effectuer la 10 reconstruction et obtenir les énergies incidentes. Dans toutes les approches classiques, comme celle de type ML-EM, cette grille est discrète, suivant un pas défini par l'utilisateur. Ce pas est l'une des limitations lors de l'étape de reconstruction puisqu'il ne sera possible de reconstruire que sur la grille définie initialement. En fonction de la complexité du problème traité, le pas 15 entre chaque énergie incidente peut être adapté. Un pas fin permettra une plus grande finesse lors de la reconstruction, mais aura un impact direct sur le calcul du projecteur. Ce dernier étant généralement calculé par simulation, un pas plus fin se traduira par une augmentation du temps de calcul, ce dernier pouvant devenir prohibitif pour certaines applications. 20 Une alternative consisterait à chercher les solutions sous forme d'un spectre discret dont le support est un espace continu. Une telle approche pour être intégrable dans le cadre d'un algorithme ML-EM demanderait à calculer le projecteur à la volée par simulation Monte-Carlo pour chaque énergie incidente au cours des itérations de l'algorithme. 25 Comme les valeurs d'énergies incidentes évoluent au gré des itérations, des calculs de réponses sont requis à chaque itération sans que les valeurs précédemment calculées puissent être réutilisées. L'ordre de grandeur de calcul d'une telle approche s'avère prohibitif pour les calculateurs actuels. Dans la reconstruction de type ML-EM, aucun a priori n'est introduit lors de 30 l'étape de reconstruction. Pour la représentativité du projecteur, une des principales limitations de ce type d'approche est liée à la représentativité du 3034221 5 projecteur utilisé lors de la reconstruction. Ce dernier est généralement obtenu suite à une étape de simulation, à l'aide d'un modèle décrivant le détecteur utilisé et son environnement. Néanmoins, il est généralement difficile d'aboutir à une modélisation fine, connaissance de la géométrie du 5 détecteur, limitations du code de calcul qui ne modélisent qu'imparfaitement les interactions au sein du détecteur, absence de prise en compte de certains phénomènes physiques, comme la collection de photons visibles par exemple dans le cas des scintillateurs plastiques. Si la modélisation du détecteur et/ou de l'environnement est imparfaite, un biais plus ou moins 1 0 important sera présent lors de la reconstruction. Une autre approche, plutôt que considérer une grille d'énergies incidentes sur laquelle on effectuerait la reconstruction est d'utiliser une base de données. Dans ce cas, on ne simule plus une grille d'énergies incidentes, mais directement la signature d'un élément donné, par exemple pour le 60Co 15 on simulera directement la réponse du détecteur pour deux énergies incidentes, émises respectivement à 1173,2 keV et 1332,5 keV. Avec cette approche la reconstruction est effectuée directement sur une grille de radioéléments. Elle présente aussi l'avantage de converger plus rapidement que l'approche avec grille, car le nombre de composantes à prendre en 20 compte est plus faible. L'un des inconvénients de cette méthode est sa limitation à ne prendre en compte que les radioéléments initialement présents dans la base de données. Si l'analyse d'un spectre gamma fait intervenir un radioélément absent de la base de données, la reconstruction sera erronée. 25 Dans la suite de la description, on désigne sous l'expression « spectre d'entrée » un ensemble d'énergies incidentes définies par un utilisateur et sous l'expression « spectre de sortie » un ensemble d'énergies mesurées. Le mot méta-donné et le mot projecteur désignent un même objet. L'invention concerne un procédé pour déterminer la nature des 30 radioéléments présents dans un objet et leur activité caractérisé en ce qu'il comporte au moins les étapes suivantes : 3034221 6 - une première phase de simulation numérique de réponses spectrométriques pour un ensemble d'énergie incidente E et un ensemble d'énergie de sortie mesurée E', afin d'obtenir un ensemble de données simulées, 5 - une deuxième phase de régression non paramétrique sur les données simulées, estimation non paramétrique de la quantité représentant la probabilité jointe des triplets (E,E',y) à partir de points simulés (Ei,Ei',yij) afin d'en déduire un méta-modèle S(E, E') pour tout couple d'énergie (E, E') sur une fonction continue, 10 - à partir du méta-modèle S(E, E'), la détermination de la nature et de l'activité des radioéléments présents dans l'objet. Selon une variante de réalisation, le procédé comporte les étapes suivantes : - on utilise n un nombre de points dans la grille d'entrée, énergies E, et 15 n' un nombre de points dans la grille de sortie, énergies E', - pour i = 1, ...,n et j = 1, ...,n', on dispose des données caractéristiques des intensités spectrales calculées À pour une énergie d'entrée EL et une énergie de sortie E'1, - on utilise une méthode d'estimation non paramétrique de la quantité 20 f(E,E',y) représentant la densité de probabilité jointe des triplets (E,E',y) à partir des points simulés (Ei,E'pyii), et on déduit un modèle S(E,E') pour tout (E, E') E rz2 où I': est un espace continu : S(E, = IE (yIE,E') fll, y . f(y1E,E)dy = fR y f(E,E',y) dy fR f (E, y) dy où le symbole IE(y) représente l'espérance mathématique de la variable aléatoire y, y est déduite des intensités spectrales calculées x,;. 25 On peut déterminer les valeurs À en utilisant un logiciel Monte- Carlo, les données À étant considérées comme des réalisations issues 3034221 7 d'une distribution de Poisson dont l'intensité ILS est estimée au moyen d'une procédure de régression non paramétrique et on introduit yii = + E) OÙ 0 E 1, puis on approxime la distribution de probabilité de yii pour des valeurs 5 suffisamment grandes de IL> > 10 par une loi gaussienne de moyenne log ILS et de variance pour la loi jointe f y), on choisit un mélange par processus de Dirichlet (DPM) comme distribution a priori et on exprime la distribution aléatoire en une somme sur l'infini de foi, composantes de f , y), f(E, y) =wk fok(E,E, y) k=1 10 paramétrée par 61k le paramètre associé à la kième composante de G mesure aléatoire définie par = 8ok(') k=1 avec w1 = V1, Wk = - V1) tel que Vk-Beta(1, a) et (Su(-) représente la fonction de Dirac localisée en u et Beta(a, b), pour 0 k,ik,ipk,flk) avec rik = = lek, k), - Xk(E) le vecteur de régresseurs centré avec É = (E,E'), et /3k le vecteur de coefficients de régression, 20 - la matrice Ek, dépendant du paramètre tpk E {0,1}, comme suit : Ek = k (01 (1 -h si R,h-1 (Tk T,k) .Rek avec Rip = °) Si tpk = 0 et Rek = k. 1 2 1 lPk = 1, permettant de choisir entre une composante alignée sur les axes E et E' (Ipk = 0) et une composante oblique orientée par la droite E = (IP k = 1). 3034221 8 chaque composante fek s'exprime alors : fOk(E, E, Y) = g I Tik, -7'r (YIK - g k(É), fiTc.4(É)) où, .7«. 62) représente la loi gaussienne de kt et de variance a2, , E) la loi gaussienne bivariée de moyenne fi E rz2 et de matrice de covariance E, où 5 la loi a priori du paramètre uk est une gaussienne, la variance -ck est distribuée suivant une loi gamma-inverse, tpk suit une loi a priori de type Bernoulli et l'a priori pour les coefficients de régression coefficients /3k est une loi normale (gaussienne) multivariée de dimension IX(É)I, et en appliquant une règle de Bayes sur l'expression de ffE,E',y) on obtient 10 l'expression f (yIE , E') Wk j'r2(EI> J^r (y lig - gk e-efek(É)) k=i L=1 W/ -7\r2n 1174,10 et le modèle probabiliste S(E, E)= IE (yIE , E) ov k N2 (É 111:' lk) le' - 'g (É>) JVZ(Lifi 1, El) k k ^ à partir de ce modèle probabiliste et des données observées par simulation , on estime la loi a posteriori f(E,E',yIE1,E'1,Y11. --- .En. , n' ynn,) et l'espérance conditionnelle .§.(E, E') = IE (yIE , , E1, E v --- En, E' n' ynn, ) 15 afin de déterminer les éléments présents dans l'objet et leur activité. Pour le calcul de l'a posteriori, on utilisera, par exemple, un schéma d'approximation Monte-Carlo par chaîne de Markov (MCMC), pour toute itération (t) de la procédure MCMC, on génère une réponse spectrale débruitée S(E,E')(0, pour T générations, la distribution a posteriori de la 20 réponse spectrale est approximée par l'ensemble des tirages S(E,E')(t) pour t = 1, T, et la réponse estimée s'exprime : T -,>-- T 1S(E , Er) t=1 3034221 9 Le schéma d'approximation peut comporter une étape d'échantillonnage par tranche utilisant un nombre aléatoire fini K de composantes pour chaque itération et en ce qu'il comporte les étapes suivantes : 5 on introduit des variables latentes de classification Kip définies pour i = 1, ...,n et j = 1, ...,n', tel que Ki = k si (Ei,ri,yii) est distribué suivant la kième composante du mélange f(E,E1,y), on définit un modèle pour les paramètres du mélange, pour tout i ri] ° (11Kii' 'Ki]) 1Kij ReKiJ. T ReKi; Kii,Ei, 01, 02, - N (yii I PLI - gKii(Ei, e -RKii équivalent à la vraisemblance f(K11, , E, E'n,, y , à partir de ces distributions de probabilités et par application de la règle de Bayes, on calcule la densité de probabilité conditionnelle f(w1, W2, --- 611, 02, --- I K11, --- Knn', E1, E 1, - - - , En, en', Y11, - - - , Y nn') en utilisant un échantillonneur de Gibbs qui à chaque itération (t) génère 15 successivement les échantillons suivants, pour tout k, W1, W2, --- I K11, --- Knn' kI Kll, --- Knn', E - -- , En, E T xi 11-k» Pk, K11, - - - Knn', E1, E ff, -- - , En, E --- Knn', E1, E -- - , En, en' µk'ki k, k, k, K11, - - - , Knn', --- En, en' le xi k, IN, K11, --- Knn', E1, E i, - - - , En, E in', Y11, - - - , Ynn' Ynn'i w1, w2, - - - , 01, 02, - - - ) 3034221 10 et au niveau des nombres de points i dans la grille d'entrée d'énergie et de points j dans la grille de sortie, pour tout i 0, fGammea,b)(X) - r(a) x e ba a-1 -bx 15 L'a priori gaussien sur les coefficients de régression /3k peut être remplacé par un a priori basé sur le tirage aléatoire de P points E' P) à partir de .7V'2(fik, Ek) et on prend pour (Yi, ,5"p) la plus proche valeur de yii correspondant à chaque point échantillonné où P est 3034221 11 plus grand que la taille du vecteur pk, on génère alors J^r(m,g, Fp) comme loi a priori (avec Ép = p)) r I= P \1 /3)el 10) par une loi gaussienne de moyenne log ILS et de variance L . 1i; La distribution de poisson sera exprimée de la manière suivante : Distribution PoissonM, pour n E N, Ân. fPoisson(Â) (ri) = e n 25 Le procédé selon l'invention repose sur l'estimation non paramétrique de la quantité f (E, E, y), représentant la densité de probabilité jointe des triplets (E,E1,y) à partir des points simulés (Ei,E'pyii), afin d'en déduire le modèle S(E, El) pour tout (E, E') E où rz est un espace continu : 3034221 16 S(E,É) = IE (yIE,É) fll, y. f(y1E,É)dy = fR y . ffE,É,y) dy fR f(E,E',y) dy où le symbole IE(y) représente l'espérance mathématique de la variable aléatoire y, ce symbole représente une espérance statistique utilisé plus loin dans la description pour caractériser l'intensité d'une énergie dans le spectre. Pour la loi jointe ffE,E',y), on choisit, par exemple, un mélange par 5 processus de Dirichlet (DPM) comme distribution a priori. Cette modélisation s'exprime à partir d'une distribution aléatoire G telle que G - DP(a, G0) où le symbole "-" signifie "est distribué suivant" et DP(a, G0) représente la distribution d'un processus de Dirichlet (Hjort N L, Holmes C, Müller P, Walker S G, Ghosal S, Lijoi A, Prünster I, The Y W, Jordan M I, Griffin J, 10 Dunson D B and Quintana F 2010 Bayesian Nonparametrics Cambridge Series in Statistical and Probabilistic Mathematics) avec une mesure moyenne Go et un paramètre a dit de concentration. Cette distribution de probabilité sur des distributions de probabilité aléatoires joue un rôle central en modélisation bayésienne non paramétrique. Une représentation courante 15 du DP est issue de Sethuraman (1994) qui exprime la mesure aléatoire G(-) comme : = Wk 8ok(-) k=1 avec w1 = V1, Wk = - V1) tel que Vk-Beta(1, a) et (Su(-) représente la fonction de Dirac localisée en u. Ici, Ok-GO représentent les paramètres associés à la kième composante de G. Le caractère non paramétrique est 20 issu du nombre de composante potentiellement infini intervenant dans la somme caractérisant l'a priori sur G. La distribution Dirichlet(a a2, --- CeK) est définie avec x = xl()' E I , xi > 0 pour 1 0 pour 1 Ifik,Ek) sN' (ylig - g k (É>), e- el 7 e..- k(É)) où, .7«. Iµ, 62) représente la loi gaussienne de kt et de variance 62, Me Irt, E) la loi gaussienne bivariée de moyenne fi E rz2 et de matrice de covariance E. Une distribution gaussienne multivariée .7V'p(161/) est définie pour 5 X E 1 e-1(x - ii)' -Sr' - (x- g) fjv- oisi)(X) = 2 (27-0 p /21n11 /2 où IAI représente ici le déterminant de la matrice A. Note : cette définition couvre également le cas univarié quand p = 1. La mesure a priori Go du processus de Dirichlet est choisie de la façon suivante. La loi a priori du paramètre iik est une gaussienne, la 1 0 variance -ck est distribuée suivant une loi gamma-inverse, tpk suit une loi a priori de type Bernoulli et l'a priori pour les coefficients de régression coefficients /3k est une loi normale (gaussienne) multivariée de dimension IX(É)I. En utilisant la règle de Bayes, nous déduisons de l'expression de f (E, E', y) que : . i(yiE,E) = y Wk j'r2 (É> 1111 k(-k) fl'k gk(É). Ei.iw1X2 E Yi, i Cette dernière expression permet de rendre compte de l'a priori 20 correspondant à l'approche bayésienne non paramétrique pour la réponse spectrale débruitée et interpolée. 3034221 19 Les caractères hétéroscédastiques et non gaussiens en tout (E, E') sont pris en compte dans le modèle de ffylE ,E') grâce au mélange de distributions gaussiennes .7^T (yliejl 1, générer Vk-Beta(1, a), poser wk = Vkrk_i, poser rk = rk_1(1 - Vk), tant que rk > u*. - Affecter à K* la valeur maximale de k. - Générer 61k pour 1 ui1) max(wk, 1/n .0 Ek) g k , k( où 15 1(A) est la fonction indicatrice : 1(A) = 1 si A est vrai et 1(A) = 0 sinon. 1 re x Générer K - - - ric,*-12-ic k=1 - Réordonner les étiquettes de composantes en suivant leur ordre d'apparition dans la génération des Kij. Affecter à Kn le nombre de 20 valeurs distinctes de Kij. Poser, pour tout k Kn, nk = #{1(ii = k}, le nombre d'observations yii assigné à la composante k. - Générer (w1, w2 - - - , u11, - - - , unn' I K11, - - - , Knn') - Générer (W1, W2, wkn,rkn) Dirichlet(nl, n2, , nkn, a) - Générer pour i n', ut] -Uniform(0,min(wk, 1/n. n')) 25 - Poser u* = minu({util) - Pour k > Kn, générer Vk-Beta(1, a) , poser Wk = Vk rk_i, poser rk = rk_1(1 - Vk), tant que rk > u*. - Affecter à K* la valeur maximale de k. - Générer (lbI ki K11, --- Knn:, , En, E'n,) pour k , 3034221 23 bat. bai' (hr-FEr )(Er/202 t,j: Kii=k )(E'rifk)2 2 {ij: Kii=k )a -r-E( - Calculer vo = nk nk - Calculer = btat. bt,at' nk nk )a.c+ )cit-F( 2 2 (bt-*{i,j: +(erifk) (Erilk)(erik) bt41E{ Kii=k)(Ei-iik)2+(erg'k) +2(E i- iik)(E'r g'k) ti - 5 - Générer v- Uniform(0,1) - Si O. Fî est la loi a priori uniforme discrète sur la collection d'effets de matrice, et a est le paramètre (positif) de concentration 5 du processus : F(-) = Wk 8X/Jk(*) k=1 avec w1 = V1, wk = Vk fi i (1- V1) tel que Vk-Beta(1, a) et (Su(-) représente la fonction de Dirac localisée en u. Avant de décrire les différents tirages aléatoires permettant d'explorer la loi a posteriori, il reste à introduire les variables Ki d'allocation du ième photon 10 observé à une composante k du mélange. L'échantillonneur de Gibbs consiste à alterner des tirages aléatoires suivant les lois conditionnelles suivantes, pour tout i 1, générer Vk-Beta(1, a), poser wk = Vk rk_i, et poser rk = rk_1(1 - Vk), tant que rk > u*. - Affecter à K* la valeur maximale de k. 5 - Générer Xk pour 1 ni) MaX(Wk, IP x* k,fic(PE;) où 1(A) est la fonction indicatrice : 1(A) = 1 si A est vrai et 1(A) = 0 sinon. - Générer Ki- -El Kn, générer Vk-Beta(1, a), poser Wk = V r Vk k-15 poser rk = rk_1(1 - Vk), tant que rk > u*, - Affecter à K* la valeur maximale de k. 25 - Générer (xk I Kn, p) pour k REVENDICATIONS1 - Procédé pour déterminer la nature des radioéléments présents dans un objet et leur activité caractérisé en ce qu'il comporte au moins les étapes suivantes : - une première phase (32) de simulation numérique de réponses spectrométriques pour un ensemble d'énergie incidente E et un ensemble d'énergie de sortie mesurée E', afin d'obtenir un ensemble de données simulées, - une deuxième phase de régression non paramétrique (33) sur les données simulées, estimation non paramétrique de la quantité représentant la probabilité jointe des triplets (E,E',y) à partir de points simulés (Ei,Ei',yij) afin d'en déduire un méta-modèle S(E, E') pour tout couple d'énergie (E, E') sur une fonction continue, - à partir du méta-modèle S(E, E'), (36) la détermination de la nature et de l'activité des radioéléments présents dans l'objet. 2 - Procédé selon la revendication 1 caractérisé en ce que : - on utilise n un nombre de points dans la grille d'entrée, énergies E, et n' un nombre de points dans la grille de sortie, énergies E', - pour i = 1, ...,n et j = 1, ...,n', on dispose des données caractéristiques des intensités spectrales calculées À pour une énergie d'entrée EL et une énergie de sortie E'1, - on utilise une méthode d'estimation non paramétrique de la quantité f(E,E',y), représentant la densité de probabilité jointe des triplets (E,E',y) à partir des points simulés (Ei,E'pyii), et on déduit un modèle S(E, E') pour tout (E, E') E où I': est un espace continu : rz 3034221 34 S(E,É) = IE (yIE,É) fR y - f , , y) dy = I y - f (yIE , dy R fR f , y) dy où le symbole lE(y) représente l'espérance mathématique de la variable aléatoire y, y est déduite des intensités spectrales calculées xv. 5 3 - Procédé selon la revendication 2 caractérisé en ce que l'on détermine les valeurs À en utilisant un logiciel Monte-Carlo, les données À étant considérées comme des réalisations issues d'une distribution de Poisson dont l'intensité ILS est estimée au moyen d'une procédure de régression non paramétrique et on introduit yii = + E) où 0 > 10 par une loi gaussienne de moyenne log ILS et de variance L pour la loi jointe f , , y) , on choisit un mélange par processus de Dirichlet (DPM) comme distribution a priori et on exprime la distribution aléatoire en 15 une somme sur l'infini de foi, composantes de f , , y), f(E, y) =wk fok(E, y) k=i paramétrée par 61k le paramètre associé à la kième composante de G mesure aléatoire définie par G(.) = Wk 8ok k=1 avec w1 = V1, Wk = - V1) tel que Vk -Beta(1, a) et (Su(-) représente la fonction de Dirac localisée en u et Beta(a, b), pour 0 4 - Procédé selon la revendication 3 caractérisé en ce que les composantes de la loi jointe foi, sont exprimées à partir des valeurs suivantes : - 0k = T k,ipk,f3k) avec fik =µlkl = lek, k), - Xk(E) le vecteur de régresseurs centré avec É = (E,E'), et /3k le vecteur 5 de coefficients de régression, - la matrice Ek, dépendant du paramètre tpk E {0,1}, comme suit : Ek = ek 2 1 = (1 1) si R-1 (Tk ° R avec R 1- ( °) si tpk = 0 et R ek. -ci). 0 1 lPk = 1, permettant de choisir entre une composante alignée sur les axes E et E' (Ipk= 0) et une composante oblique orientée par la droite E = 10 (Pk = 1). chaque composante foi, s'exprime alors : fek(E, E , y) = .7V'2(É> I fik, k) N (ylig - g k(É>), 0-2) représente la loi gaussienne de kt et de variance 62, .N2(- E) la loi gaussienne bivariée de moyenne fi E rz2 et de matrice de covariance E 15 où la loi a priori du paramètre µk est une gaussienne, la variance -ck est distribuée suivant une loi gamma-inverse, tpk suit une loi a priori de type Bernoulli et l'a priori pour les coefficients de régression coefficients /3k est une loi normale (gaussienne) multivariée de dimension IX(É)I, et 20 en appliquant une règle de Bayes sur l'expression de ffE,E',y) on obtient l'expression f fylE , E = w k j'r2 (É> fi 11:1' k) \ (ylig - g kg), e- el( k(É)) J2 (E irt1,zo et le modèle probabiliste S(E ,E') = E (yIE,E') wk (É> 111:i' \ le' - g (É) wi J\r2(r El) k k 3034221 36 A partir de ce modèle probabiliste et des données observées par simulation Épyii) , on estime la loi a posteriori f(E,e,y1E1_,E Y11 --- En, Ei Ynni) et l'espérance conditionnelle .§"(E, E') = IE (y1E,E',E1_,E v --- En, Ei Ynni) afin de déterminer les éléments présents dans l'objet et leur activité. 5 5 - Procédé selon la revendication 4 caractérisé en ce que pour le calcul de l'a posteriori, on utilise un schéma d'approximation Monte-Carlo par chaîne de Markov (MCMC), pour toute itération (t) de la procédure MCMC, on génère une réponse 10 spectrale débruitée S(E,E ft), pour T générations, la distribution a posteriori de la réponse spectrale est approximée par l'ensemble des tirages S(E,É)(t) pour t = 1, ...,T, et la réponse estimée s'exprime : T 1 T 1S(EEft) t=1 15 6 - Procédé selon la revendication 5 caractérisé en ce que le schéma d'approximation comporte une étape d'échantillonnage par tranche utilisant un nombre aléatoire fini K de composantes pour chaque itération et en ce qu'il comporte les étapes suivantes : on introduit des variables latentes de classification Kii, définies pour i = 20 1,...,n et j = tel que Ki= k si (Ei,ri,yii) est distribué suivant la kième composante du mélange f(E,E",y), on définit un modèle pour les paramètres du mélange, pour tout i 7 - Procédé selon la revendication 6 caractérisé en ce que l'on introduit une variable auxiliaire uij afin de ne générer qu'un nombre aléatoire fini K de composantes à l'itération (t) tout en évitant une troncature arbitraire du 15 modèle. 3034221 38 8 - Procédé selon l'une des revendications 5 à 7 caractérisé en ce qu'il comporte une étape de calcul de l'écart-type a posteriori et des intervalles crédibles à partir de l'ensemble {S(E,E')(t)} 1 7 65 ,E) 1 , E') - S , E')(t))2) . 5 9 - Procédé selon l'une des revendications 2 à 8 caractérisé en ce que l'on introduit un a priori Gamma4b, sur les paramètres d'échelle br et b'T dans la distribution des amplitudes de composantes conduisant à une distribution a posteriori Gamma : br- Gamma / Kn avec Gamma(a, b), pour x > 0, 1 (Pb Kn, +- k-Ck =1 1, b a xa-1 e-bx. fGamma(a,b) (X) r(a) 10 10 - Procédé selon la revendication 4 caractérisé en ce que l'on remplace l'a priori gaussien sur les coefficients de régression lek par un a priori basé sur le tirage aléatoire de P points (Ê''1, à partir de .7V'2(fik, Ek) et on prend pour (5,1, la plus proche valeur de yii correspondant à chaque 15 point échantillonné où P est plus grand que la taille du vecteur pk, on génère alors /3k - .7V'(Mig,Fig) comme loi a priori (avec Ép p)) / P \1 r= g k gk(Ep) - e -5719 \p=1 me = re - L gk( \P.1 t=1 3034221 39 11 - Procédé selon l'une des revendications précédentes caractérisé en ce que l'on génère un méta-modèle étendu noté S(E,E'',) où est un paramètre identifié par un indice entier et caractéristique d'un effet de matrice. 5 12 - Procédé selon l'une des revendications précédentes caractérisé en ce que l'on estime l'activité des radioéléments en exécutant les étapes suivantes : soit Nx le nombre d'émetteurs retenus pour l'élément x considéré et n-x,/ pour 10 1 = 1, ...,Nx, les probabilités d'émissions associées ainsi que vx,/ pour 1 = 1, ...,Nx, les énergies correspondantes, on définit la réponse du radionucléide tpx(E, 0, pour une énergie observée E' et un effet de matrice par Nx tp =71" S(v xj, E', 1=1 * X.on définit = fo tp xe') de et rcx/ = - pour tout 1 = 1, ...,Nx, x,e 15 on définit une réponse de radionucléide normalisée par Nx tpx =71-x* S(v xj, E, 1=1 on vérifie que fo- tpx* (El) dE' = 1 on exprime la densité de probabilité du ième photon observé dans le canal d'énergie E; pour i = 1, n J (Ei) = Wk E;) k=1 où p est un paramètre positif de conversion de l'indice de canal en énergie 20 (keV/canal). On introduit les variables Ki d'allocation du ième photon observé à une composante k du mélange, 3034221 40 On alterne des tirages aléatoires suivant des lois conditionnelles suivantes, pour tout i 13 - Procédé selon la revendication 12 caractérisé en ce que l'on calcule l'écart-type a posteriori des activités GAXk ''' T w(t)A(t) )2 1 k Xic,k )2 ( 1 Y ('4--- V - 1 L \ Xk _i t=1 et/ou le spectre d'entrée estimé, déconvolué de la réponse du système N (t) T Xk 1 ( S(E) R.--,' T Wkt) ITXk 1 (t) (t) (5v (t) (E) . ,,- /c Xk 't t=1 1=1 10 15