FR2846830A1

FR2846830A1 - Procede de filtrage spatio-temporel du bruit fluoroscopique

Info

Publication number: FR2846830A1
Application number: FR0213727A
Authority: FR
Inventors: Cyril Allouche
Original assignee: GE Medical Systems Global Technology Co LLC
Current assignee: GE Medical Systems Global Technology Co LLC
Priority date: 2002-10-31
Filing date: 2002-10-31
Publication date: 2004-05-07
Anticipated expiration: 2022-10-31
Also published as: US20040086194A1; US7379615B2; DE10350697A1; FR2846830B1

Abstract

Pour réduire le bruit fluoroscopique dans une image I acquise à une date t, on appaire les points de cette image avec les points d'une image l' acquise à une date t-1. Pour un point de coordonnées (x,y) de l'image I, on effectue une convolution avec un noyau U équivalent à un filtre passe bas dont les coefficients ont été modifiés en fonction du voisinage du point de coordonnées (x,y) dans l'image I. Pour le point appairé dans l'image l'on effectue une convolution avec le noyau U dont les coefficients ont été modifiés en fonction du voisinage du point de coordonnées (x,y) dans l'image l'. On associe linéairement le résultat des deux convolutions afin d'obtenir une valeur filtrée pour le point de coordonnées (x,y). On répète ces opérations pour chaque points de l'image I.

Description

I Procédé de filtrage spatio-temporel du bruit fluoroscopique La présente

invention a pour objet un procédé de filtrage spatiotemporel du bruit fluoroscopique. Le domaine de l'invention est celui de 5 l'imagerie et plus précisément celui de la réduction du bruit fluoroscopique dans les images. Le domaine de l'invention est plus particulièrement celui de

la réduction du bruit dans les images acquises en séquences temporelles.

Un but de l'invention est la réduction du bruit fluoroscopique dans des images. Un autre but de l'invention est la réduction du bruit fluoroscopique

dans des images issues d'un procédé d'acquisition fluoroscopique.

Un autre but de l'invention est la réduction en temps réelle d'un bruit fluoroscopique. Dans l'état de la technique on connaît des appareils d'acquisition 15 d'images de type scanner qui permettent d'obtenir des images de l'intérieur d'un organisme vivant, et en particulier de l'intérieur d'un corps humain. Cela permet d'obtenir des images d'organes internes et d'effectuer des diagnostics. Ces images sont obtenues par l'exposition d'un patient à un rayonnement qui est reçu sur un détecteur après avoir traversé le patient. Le 20 détecteur produit alors une images qui est interprétable par un praticien spécialisé. Ces images comporte du bruit dit fluoroscopique. Ce bruit fluoroscopique est la résultante du bruit quantique, du à la nature du

rayonnement, et du bruit électronique, du à la nature du détecteur.

Le bruit fluoroscopique est présent dans l'image et brouille donc le 25 signal, aussi appelé signal utile ou information, présent dans l'image ce qui

rend l'interprétation des images, difficile, hasardeuse, voir impossible.

Dans l'état de la technique on sait que pour améliorer le rapport signal sur bruit, ou SNR, c'est à dire pour augmenter le ratio d'information présent dans l'image, on peut augmenter l'intensité du rayonnement. Cependant ce 30 faisant on rend l'examen plus traumatisant pour le patient ce qui n'est pas

acceptable en terme de santé.

Dans l'état de la technique on connaît aussi une technique de filtrage du bruit fluoroscopique qui consiste à effectuer une moyenne temporelle entre la valeur de deux points, ou pixels, de même coordonnées dans deux 35 images. Les deux images appartiennent à une même séquence d'images représentant une même zone. En d'autres termes les deux images ont les même cadrage et paramètres d'exposition mais sont prises à des dates différentes. Si on filtre une image I à une date t, on effectue alors, pour chaque point de l'image I une moyenne avec le point correspondant de l'image l' obtenue à t-I. Cette méthode présente plusieurs inconvénients. Un premier inconvénient et la faible réduction du bruit dans l'image filtrée. Cette

réduction est de l'ordre de, au maximum, 42.

Un deuxième inconvénient est le problème de la rémanence, ou de 10 l'apparition de fantômes dans les images filtrées. En effet, plaçons nous dans le cas de l'acquisition d'images d'une artère, sachant que dans cette artère a été introduit un guide. Un guide est un objet métallique cylindrique qui est introduit dans les artères et qui est donc visible en radiographie. Sachant que ce guide se déplace rapidement, il est possible qu'il soit présent dans la zone 15 imagée à la date t-1, mais pas à la date t. Cependant comme l'image filtrée de la date t est obtenue en moyennant l'image acquise à la date t avec l'image acquise à la date t-I on obtient une image filtrée pour la date t faisant apparaître un guide alors que celui ci n'est pas présent à cette date. On a

donc ajouter à l'image de la date t une information erronée.

Dans la meilleur des solutions de l'état de la technique on a donc une réduction du bruit de 15% et un taux de rémanence de 10%, ce qui est très

insatisfaisant car source d'erreur ou d'impossibilité de diagnostic.

Dans l'invention on résout ces problèmes en utilisant un noyau de convolution, équivalent à un filtre passe bas, dont les coefficients s'adaptent 25 pour chaque point de l'image à filtrer. L'adaptation se fait en fonction d'un

écart entre la valeur du point à filtrer et la valeur des voisins de ce point. En effet un noyau de convolution est une matrice carré que l'on superpose à l'image à filtrer, à chaque coefficient de la matrice correspond alors un point de l'image. L'adaptation est aussi réalisée en fonction d'une statistique v(v) 30 de bruit connue pour la valeur v du point à filtrer.

Dans l'invention on prend aussi en considération le passé puisque que la valeur du point filtré est le résultat de la somme de deux convolutions pondérées par des facteurs y et (1- y). avec y dans l'intervalle [0, 1]. La première convolution est celle du paragraphe précédent, la deuxième 35 convolution est réalisé sur l'image précédant l'image à filtrer dans une séquence d'images. Toutes les convolutions sont réalisées à partir du même noyau pondéré selon son environnement. Pour la deuxième convolution, le noyau est pondéré en fonction de l'écart entre la valeur v du point à filtrer et de la valeur des points correspondant au noyau superposé à l'image, et en fonction de .o (v). On obtient ainsi un filtrage auto adaptatif et paramétrable, il s'agit donc d'un filtre de convolution spatio-temporel. Le paramètre y permet de régler la dépendance temporelle, et le paramètre A permet de régler la force

du filtrage.

L'invention a donc pour objet un procédé de filtrage spatio-temporel du bruit fluoroscopique dans une image caractérisé en ce que: - pour chaque point de coordonnées (x,y) d'une première image on pondère les coefficients U(k,l) d'un premier noyau de convolution de dimension D, équivalent à un filtre passe bas, en fonction d'un coefficient G 15 fonction de l'écart calculé entre l(x,y) et l(x+k, y+l), ou l(x,y) est l'intensité du point de coordonnées (x,y) de la première image, et k et I sont des indices permettant l'exploration des coefficients du noyau de convolution, on obtient ainsi un deuxième noyau de convolution de coefficients Up(k,l), - pour chaque point de coordonnées (x,y) de la première image on 20 pondère les coefficients U(k,l) du premier noyau de convolution en fonction du coefficient G fonction de l'écart calculé entre l(x,y) et l'(x+k, y+l) ou l'(x,y) est l'intensité du point de coordonnées (x,y) d'une deuxième image, on obtient ainsi un troisième noyau de convolution de coefficients Up'(k,l), - on calcul la valeur filtrée de l(x,y) par la formule: F(x,y)=( t (y*Up(k,l).I(x+k,y+l)+(1-r) *Up'(k,l).I'(x+k,y+l)) /N......(1) L = (Dk- 1)

=( 2......(2)

YG [0,11....(3)

L L N= % (y* Up(k,l)+(1- y)*Up'(k,l)).....(4) k=-Ll=-L

ou F(x,y) est la valeur de l(x,y) filtrée.

L'invention sera mieux comprise à la lecture de la description qui suit

et à l'examen des figures qui l'accompagnent. Celles ci sont présentées à titre indicatif et nullement limitatif de l'invention. Les figures montrent: Figure 1: un bloc diagramme illustrant les principes de fonctionnement

d'un filtre selon l'invention.

Figure 2: Une illustration d'une fonction de réponse fréquentielle d'un

filtre passe bas.

Figure 3: Une illustration de l'indexage d'une image pour la mise en oeuvre du procédé selon l'invention.

Figure 4: Une illustration d'étapes du procédé selon l'invention.

Figure 5: Une illustration d'une courbe gaussienne.

La figure illustre le fonctionnement d'un filtre selon l'invention. La figure 1 montre un bloc 101 correspondant à l'implémentation du filtre. Le 10 filtre est implémenté sous forme d'un programme ou d'un composant spécialisé. Le filtre, ou procédé selon l'invention, est donc mis en oeuvre via donc des codes instructions exécutés par une machine comportant des circuits logique de type microprocesseur. Ces codes instructions sont alors

enregistrés dans une mémoire de programme.

Ces composants logiques et mémoires font alors parti d'un appareil de traitement. Un tel appareil est soit l'appareil d'acquisition d'images lui-même, soit un ordinateur connecté à l'appareil d'acquisition d'images. Ici par ordinateur on entend ordinateur personnel, station de travail ou autre. Un

appareil d'acquisition d'image est, par exemple, un scanner.

La figure 1 montre qu'un filtre selon l'invention accepte en entrée au moins jusqu'à 5 paramètres ou entrées. On assimile ici le filtre 101 à un appareil mettant en oeuvre le procédé selon l'invention. Une première entrée est une image 102, prise à une date t, à filtrer. Une telle image a des caractéristiques illustrées par la figure 3. L'image 102 est assimilable à une 25 image bitmap. L'image 102 est donc une image numérique délivrée par l'appareil d'acquisition d'images ou obtenue à partir de signaux délivrés par cet appareil. Chaque point, ou pixel, de l'image 102 est repéré par ses coordonnées (x,y) et a une intensité l(x,y) . Ici l'intensité est équivalente à un niveau de gris. Les dynamiques de x,y et I dépendent du détecteur utilisé par 30 l'appareil d'acquisition d'images. De manière général x varie de 0 à Xmax, et

y varie de 0 à Ymax. Xmax et Ymax définissent alors la résolution de l'image.

1(x,y) représente la valeur du point de coordonnées (x,y), c'est à dire son

niveau de gris. I varie de 0 à Imax.

Le filtre 101 accepte aussi en entrée une image 103 prise à une date 35 tI. L'image 103 appartient à la même séquence d'images que l'image 102.

L'image 103 occupe, dans la séquence d'images, la place précédant l'image 102. Les images 102 et 103 ont la même définition et représentent la même zone de l'espace. Dans une variante préférée de l'invention l'image 103 est le résultat du filtrage, par le procédé selon l'invention, de l'image prise à la 5 date t-1. Une séquence d'images est réalisée au cours d'un examen en prenant plusieurs clichés successifs dune même zone, comme le ferait une caméra. On trouve de telles séquences d'images lorsque l'on cherche à imager un organe, ou une zone d'un organisme, au cours d'un cycle. Des exemples de cycles sont le cycle cardiaque, le cycle respiratoire, ou une 10 période de temps fixée arbitrairement.

Le filtre 101 accepte aussi en entrée des statistiques 104 de bruit. Le bruit traité est ici le bruit fluoroscopique dont on connaît les statistiques en fonction de la valeur du niveau de gris du point considéré. On connaît en particulier l'écart type de ce bruit car on sait que celui ci est en racine carrée 15 du nombre de photons arrivant sur le détecteur. Ce nombre de photons est lui même proportionnel au niveau de gris. On sait donc déterminer des statistiques, en particulier l'écart type, pour le bruit fluoroscopique. Une méthode simple pour modéliser un bruit fluoroscopique présent lors d'une radiographie, consiste à utiliser deux images successives d'une même zone 20 de manière à pouvoir appairer les points des deux images en fonction de la zone de l'espace qu'ils représentent. On regroupe les paires de points ainsi obtenues en sous groupes selon leur niveau de gris. Pour chaque sous groupe on calcul l'écart type a des lt(x, y) - lt-1 (x, y). On discrimine un sous groupe en en supprimant les points pour lesquels lt(x, y) - It-1(x, y) est 25 supérieur à la moyenne des lt(x, y) - It-1(x, y) plus k fois l'écart type. On répète ces calculs un certain nombre de fois et pour chaque sous groupe. On obtient alors des paires de point (v, a) à partir desquelles on effectue une régression itérative pour obtenir un modèle de bruit en c(v) = a.AIv + P.v + K, ou v est le niveau de gris et a, j3 et K des coefficients définissant le bruit. 30 D'autres méthodes sont possible comme la constitution de tables / abaques

permettant d'obtenir une statistique pour le bruit.

Le filtre 101 accepte d'autre part un paramètre y et un paramètre . Le paramètre y appartient à l'intervalle [0,1] et est d'autant plus proche de 1 que l'on accorde moins d'importance à l'aspect temporel du filtre. Le paramètre 35 permet de régler la force du filtre 101. La force permet de jouer sur la tolérance au bruit. Plus est faible, moins on filtre spatialement, plus est

fort, plus on filtre spatialement.

Le filtre 101 met encore en oeuvre un noyau U, ou noyau 105, de convolution de dimension D. Un noyau de convolution est une matrice carré 5 comportant D2 coefficients notés U(k,l) avec k et I appartenant à l'intervalle [(D-1)/2, ((D-1)/2]. Lorsque l'on réalise une convolution d'une image I par la noyau U, cela signifie que l'on produit une nouvelle image F pour laquelle: F(x, y) = 1 1 (U(k, 1).I(x + k, y + 1)) k=-Ll=- t Le noyau U est équivalent à un filtre passe bas tel que celui dont la courbe de réponse en fréquence est représentée à la figure 2. Un tel filtre a donc un 10 gain de 1 pour une fréquence nulle, et son gain décroît rapidement avec la fréquence. Son implémentation par un noyau de convolution, dans le cas particulier ou D vaut 5, obtenu en prenant les premiers coefficients de la transformée de Fourier inverse donne la matrice suivante: 1 \ k -2 -1 0 1 2

-2 0,0315715 0,0172795 -0,003931 0,0172795 0,0315715

-1 0,0172795 -0,00638812 0,0965786 -0,00638812 0,0172795

0 -0,003931 0,0965786 0,39044 0,0965786 -0,003931

1 0,0172795 -0,00638812 0,0965786 -0,00638812 0,0172795

2 0,0315715 0,0172795 -0,003931 0,0172795 0,0315715

On note que l'on a choisi D valant 5 pour des raisons de puissance de calcul et de compatibilité temps réel. Il faut en effet que le temps de filtrage reste compatible avec une utilisation temps réel et ne mette donc pas à l'épreuve la patiente humaine, ni ne constitue un goulet d'étranglement pour la 20 cadence de production des images. C'est à dire il faut que l'ont puisse traiter les images aussi vite que l'appareil d'acquisition les produit. De même le filtre passe bas employé n'est qu'un filtre passe bas parmi la multitude de filtres passe bas existant correspondant à autant de noyau de convolution. Le temps réel s'entent ici par une cadence de filtrage telle que l'on excède pas 25 1 00ms de traitement par image. Cette cadence est idéalement d'au maximum 30 ms par image. A l'heure actuelle ces performances sont atteintes pour un cot raisonnable avec D valant 5. Rien n'exclue à l'avenir que, conformément à la "loi de Moore", les contraintes de temps réelle et de

cots puissent être respectées avec D valant 7, 9 ou plus.

Lorsque l'on applique le noyau de convolution a des points situés sur le bord de l'image, c'est à dire des points tels que x < (D-1)/2, ou x > Xmax (D-1)/2, ou y < (D-1)/2, ou y > Ymax - (D-1)/2, alors le noyau de convolution déborde de l'image. Pour les coefficients du noyau de convolution 5 correspondant à des points situés en dehors de l'image à filtrer on considère qu'il correspondent à des points dont la valeur de niveau de gris est fixé

arbitrairement, en générale à 0.

Le filtre 101 produit en sortie une image 106 correspondant à l'image

102 filtrée. Cette image 106 est réutilisée comme image 103 pour le 10 traitement de l'image acquise à la date t+1.

La figure 4 montre une étape 401 préliminaire d'acquisition d'images.

Le procédé selon l'invention fonctionne de manière optimale avec deux images appartenant à une même séquence d'images. Cependant il fonctionne aussi avec une seule image, c'est le cas particulier ou y vaut 1. 15 Une séquence d'images est produite par un appareil d'acquisition d'images de type scanner. L'acquisition commence à une date origine que l'on note tO et valant O de manière arbitraire. Pour chaque date tn, l'appareil d'acquisition d'image produit une image. Toutes les images produites ont la même définition. Sinon on fait en sorte de ramener l'image ayant la plus grande 20 définition à la définition de celle qui a la définition la plus petite. Ces images sont disponibles sous forme numérique, chaque point d'une image ayant alors, comme précédemment décrit, des coordonnées (x,y) et une intensité

l(x,y) équivalente à un niveau de gris.

De l'étape 401 on passe à une étape 402 de traitement de l'image t-I. 25 Dans notre exemple on va décrire le traitement de l'image acquise à la date t > tO. L'étape 402 permet donc de produire l'image 103 utilisée par le filtre 101 pour le traitement de l'image 102 acquise à la date t. L'étape 402 est identique à ce qui va être décrit pour le traitement de l'image acquise à la

date t dans l'étape 403.

A l'initialisation du procédé de filtrage, c'est à dire au moment du

traitement de l'image acquise à la date tO, il n'y a pas d'image tO-1 pour jouer le rôle de l'image 103. Dans le cas de la première image filtrée, soit on fixe y à 1, soit on ne filtre pas l'image acquise à tO. Ne pas la filtrer revient à la considérer comme bonne et à l'utiliser telle quelle dans le rôle de l'image 103 35 pour le filtrage de l'image acquise à la date tO+1.

De l'étape 402, on passe à une étape 403 de traitement de l'image I acquise à la date t. L'étape 403 comporte plusieurs sous-étapes. Le traitement d'une image se fait point par point. Dans l'étape 404, initiant l'étape 403, on commence donc par initialiser le traitement en se positionnant 5 sur le premier point à traiter. Dans notre exemple on considère qu'il s'agit du point (0,0). D'une manière plus général on va décrire le traitement du point de coordonnées (x,y). De l'étape 404, on passe à une étape 405 de calcul

des Up(k,l) pour le point de coordonnées (x,y).

L'étape 404 est l'étape dans laquelle on pondère les coefficients U(k,l) 10 du noyau U en fonction du voisinage du point de coordonnées (x,y) dans l'image I acquise à la date t. Le voisinage d'un point de coordonnées (x, y) dans l'image I est formé par les D^2-1 points de l'image I les plus proches du point de coordonnées (x,y). Le point de coordonnées (x,y) à une intensité, ou niveau de gris, l(x,y) dans l'image I. Le résultat de cette pondération est un 15 noyau Up ayant des coefficients Up(k,l). Pour le point de coordonnées (x,y) les Up(k,l) sont obtenus par la formule suivante: Up(k,l) = U(k,l) x G(l(x+k, y+l) - I(x,y); a(l(x,y))). Dans une variante préférée, G(E; o) est la valeur de la courbe gaussienne centrée en 0, de déviation standard a et telle que G(0; a) = 1. La 20 figure 5 illustre une telle courbe gaussienne. a(l(x,y)) est l'estimation de

l'écart type du bruit fluoroscopique pour l(x,y).

Dans une autre variante G(E; a) est réduit à G(E) qui est une fonction

linéaire de type G(E) = -a. E + 1, avec a > 0.

De l'étape 405 on passe à une étape 406 dans laquelle on pondère 25 les coefficients U(k,l) du noyau U en fonction du voisinage du point de

coordonnées (x,y) dans l'image 1' acquise à la date t. Le point de coordonnées (x,y) a une intensité, ou niveau de gris, l'(x,y) dans l'image 1'. Le résultat de cette pondération est un noyau Up' ayant des coefficients Up'(k,l).

Pour le point de coordonnées (x,y) les U'p(k,l) sont obtenus par la formule 30 suivante: Up'(k,l) = U(k,l) x G(l'(x+k, y+l) - I(x,y); A.o(l(x, y))). Dans une variante préférée, on utilise la formule:

Up'(k,l) = U(k,l) x G(F'(x+k, y+l) - I(x,y); .a(l(x,y))) dans laquelle F'(x,y) est l'intensité du point de coordonnées (x,y) dans 35 I'image I' après filtrage par le procédé selon l'invention.

Le coefficient permet de régler la force du filtrage. Au même titre que y, est un paramètre du filtrage sur lesquels un utilisateur peut agir quand bon lui semble: au début du traitement d'une image, pendant le traitement;

au début du traitement d'une séquence, il n'y a pas de contrainte.

De l'étape 406 on passe à une étape 407 de calcul des F(x,y), ou

F(x,y) est l'intensité du point de coordonnées (x,y) de l'image I après filtrage.

F(x,y) est obtenu par la formule: L L F(x,y) = (r*Up(k,l).I(x+k,y+l)+(1-y) *Up'(k,l).I'(x+k,y+l)) IN Dans cette formule le coefficient y permet de régler la dépendance temporelle du filtre 101. y appartient à l'intervalle [0, 1]. Plus y est proche de 10 0, plus on attache d'importance à l'image 1'. A l'extrême, si y vaut 0, l'image I

n'a plus aucune importance. Plus y est proche de 1, moins on attache d'importance à l'image 1'. A l'extrême, si y vaut 1, l'image I' n'a plus aucune importance. Si y vaut 1 le procédé est équivalent à un filtre purement spatial.

En particulier on peut alors l'appliquer au filtrage d'images fixes 15 n'appartenant pas à une séquence d'images.

Le coefficient N est un coefficient de normalisation. Il est standard pour un noyau, il est obtenu par la formule: L L N= Z y (y* Up(k, I) + (1 - r) * Up'(k,)) k=-Li=-L Dans ces formules L vaut (D-1)/2, ou D est la dimension du noyau. D

est donc par nature impaire car D-1 doit être divisible par 2.

De l'étape 407 on passe à une étape 408 dans laquelle on vérifie si tous les points de l'image I ont bien été filtrés. A la fin du filtrage on est en possession d'une image F filtrée pour laquelle le point de coordonnées (x,y) a une intensité F(x,y), F(x,y) étant fonction de l(x,y) de l'image I à filtrer. Dans la mesure ou on connaît la définition de l'image 1, le nombre de points à 25 traiter est égal à sa définition vertical par sa définition horizontale. Un simple compteur permet de définir un critère de fin de traitement. On peut aussi ne traiter qu'une zone de l'image. Dans ce cas on connaît la définition de la

zone et donc on connaît le nombre de points à traiter.

S'il reste des points à traiter, on passe à une étape 409 d'itération de x 30 et y, sinon, on passe à une étape 410 de traitement d'une image suivante.

Dans l'étape 409 on augmente la valeur de x et ou de y en fonction du

mode de parcours de l'image. Dans la pratique on parcours l'image en lignes.

C'est à dire on augmente x d'une unité puis, quand x vaut Xmax, on augmente y d'une unité et on remet x à 0. On répète ce processus jusqu'à ce que x, respectivement y, valent Xmax, respectivement Ymax. D'autre 5 parcours sont possibles, en colonnes par exemple. D'une manière générale, tout parcours permettant de passer par tous les points de l'image I zone à

filtrer sont acceptables pour la mise en oeuvre du procédé selon l'invention.

De l'étape 409 on passe à l'étape 405 de calcul de Up(k,l).

Dans l'étape 410 on vérifie si la séquence d'images comporte une 10 image t+1 à traiter / filtrer. Si c'est la cas, on passe de l'étape 410 à l'étape 403 à la différence que l'image traitée n'est plus l'image acquise à la date t, mais l'image acquise à la date t+l, et ainsi de suite avec t+n, tant qu'il y a

des images à traiter.

L'étape 410 peut aussi comporter ou être suivie d'une étape de 15 stockage des images filtrées en vue de leur réutilisation pour une impression, un affichage ou un autre traitement. A ces fins l'appareil mettant en oeuvre le procédé selon l'invention comporte une mémoire de stockage, et ou des moyens de connexion à un écran, et ou des moyens de communication avec un tiers appareil traitant des images filtrées, et ou des moyens de traitement 20 d'images filtrées. Ces moyens de traitement sont, par exemple, des codes instructions mis en oeuvre de la même manière que le procédé selon l'invention. A.

Claims

REVENDICATIONS

1 - Procédé de filtrage spatio-temporel du bruit fluoroscopique dans une image caractérisé en ce que: - pour chaque point de coordonnées (x,y) d'une première image on pondère les coefficients U(k,l) d'un premier noyau de convolution de dimension D, équivalent à un filtre passe bas, en fonction d'un coefficient G fonction de l'écart calculé entre l(x,y) et l(x+k, y+l), ou l(x,y) est l'intensité du point de coordonnées (x,y) de la première image, et k et I sont des indices 10 permettant l'exploration des coefficients du noyau de convolution, on obtient ainsi un deuxième noyau de convolution de coefficients Up(k,l), - pour chaque point de coordonnées (x,y) de la première image on pondère les coefficients U(k,l) du premier noyau de convolution en fonction du coefficient G fonction de l'écart calculé entre l(x,y) et l'(x+k, y+l) ou '(x,y) 15 est l'intensité du point de coordonnées (x,y) d'une deuxième image, on obtient ainsi un troisième noyau de convolution de coefficients Up'(k,l), - on calcul la valeur filtrée de l(x,y) par la formule: F(x,y) = * Up(k,l).I(x +k,y+l)+ (1 -y)*Up'(k,l).I'(x+k,y+) N.....(1) k=-JI=-L

L = (D- 1..... (2)

r E [0,11....(3) L L N= y (y *Up(k,l)+(1-y)*Up'(k,1)).....(4) k=-Ll=-L

ou F(x,y) est la valeur de l(x,y) filtrée.

2 - Procédé selon la revendication 1 caractérisé en ce que: 20 Up(k,l) = U(k,l)xG(l(x+k,y+l)-l(x,y); o(l(x,y))), et U'p(k,l) = U(k,l)xG(l'(x+k,y+l) -l(x,y); A.ao(l(x,y))) avec G une fonction de pondération dépendant d'un écart entre la valeur du point à filtrer et de son voisinage, dépendant d'une statistique o de

bruit pour la valeur du point à filtrer, et dépendant d'une force A de filtrage.

3 - Procédé selon la revendication 2, caractérisé en ce que G est fonction de l'écart s calculé et d'une statistique o de bruit connu pour l(x,y), on note alors le coefficient G comme la fonction G(s, o), ou G est donc la

valeur en s d'une gaussienne centrée en 0 et de déviation standard v.

4 - Procédé selon la revendication 2, caractérisé en ce que G est t 9 fonction de l'écart E calculé de type: G(E) = -a. ú + 1, avec a > 0, et Up(k,I) = U(k,l)xG(l(x+k,y+l)-l(x,y)), et U'p(k,l) = U(k,l)xG(l'(x+k,y+l)I(x,y))

5 - Procédé selon l'une des revendications 2 à 4, caractérisé en ce

que est un nombre réel.

6 - Procédé selon l'une des revendications I à 5, caractérisé en ce

que la formule (1) devient: F(x,y)= * (yUp(k,l).I(x + k,y+l)+(1 -y)*Up'(k, l).F'(x + k,Y+l)) / N k=-LI=-I,

ou F'(x,y) est l'intensité filtrée du point de coordonnées (x,y) de la 10 deuxième image.

7 - Procédé selon l'une des revendications 1 à 6, une valeur de y

valant 0 implique une dépendante temporelle nulle.

8 - Procédé selon l'une des revendications 1 à 7, caractérisé en ce

que les première et deuxième images sont des images successives d'une 15 séquence d'images, la première image ayant une date t, et la deuxième

image ayant une date t-1.

9 - Procédé selon l'une des revendications 1 à 8, caractérisé en ce

que D est un nombre impaire.

- Procédé selon l'une des revendications 1 à 9, caractérisé en ce 20 que Dvaut 5.

11 - Procédé selon l'une des revendications 1 à 9, caractérisé en ce

que D est supérieur à 5.