FR 2689272 A1 19931001 FR 9203900 A 19920331 Dispositif de traitement d'informations d'ima Qes tridimension- nelles, avec extraction de lignes remarquables L'invention concerne le traitement d'images numériques tridimen- sionnelles Elle s'applique en particulier, mais non exclusive- ment, à des images dites médicales. En matière médicale, les appareils modernes tels que les scanners à rayons X ou les appareils de résonance magnétique nucléaire fournissent des informations numériques de mesure rapportées à des points de l'espace tridimensionnel Il y a donc trois coordonnées pour définir le point de l'espace, plus la grandeur mesurée, attachée à ce point On appelle ici "image tridimensionnelle" toute représentation numérique d'un jeu de grandeurs mesurées ou calculées ("valeur"), rapportées respecti- vement à un ensemble de points de l'espace. Le traitement des images tridimensionnelles est particulièrement délicat, en particulier dès que l'homme doit intervenir Par exemple, la visualisation d'une image tridimensionnelle sur un écran (bidimensionnel) ne peut se faire de façon immédiate: si l'on affiche sur cet écran tous les points de l'image en perspective, avec représentation de leur "valeur", le résultat sur l'écran sera parfaitement indéchiffrable. Le plus simple est alors de considérer seulement une coupe de l'image, c'est-à-dire un sous-ensemble bidimensionnel de celle- ci Mais il est clair qu'une coupe comporte beaucoup moins d'informations utiles que l'ensemble de l'image tridimensionnel- le Même en utilisant de nombreuses coupes, on ne retrouve pas toutes ces informations utiles, et encore moins leurs relations. Un moyen connu consiste à rechercher dans l'image tridimension- nelle des surfaces dites "d'iso-intensité"l (plus brièvement "isosurfaces") Ces iso-surfaces sont l'analogue dans le cas des images tridimensionnelles des lignes de niveau pour le cas des images bidimensionnelles (o l'altitude est considérée comme "valeur" en chaque point). Une iso-surface est définie par rapport à un seuil choisi, appliqué à ladite "valeur" mesurée Le seuil est appelé "iso- valeur" L'iso-surface est la surface qui sépare: les régions de l'espace pour lesquelles la valeur mesurée est inférieure ou égale au seuil, et les régions pour lesquelles la valeur mesurée est strictement supérieure au seuil. En théorie, la détermination d'une iso-surface supposerait que l'on connaisse une fonction continue capable d'associer à tout point de l'espace la valeur mesurée. En pratique, on ne dispose que d'une image numérique définie par échantillonnage pour un nombre limité de points de l'espace Une fonction d' interpolation permet de calculer la "valeur" pour des points situés entre deux points échantillonnés. L'homme du métier sait que les iso-surfaces possèdent des propriétés topologiques spécifiques: elles sont continues, sans trou, et ne s'intersectent pas elles-mêmes. Concrètement, l'article "Marching Cubes: A high résolution 3 D surface reconstruction algorithm", William E LORENSEN et Harvey E CLINE, Computer Graphics, volume 31, NO 4, Juillet 1987, décrit une technique d'extraction des iso-surfaces d'une image tridimensionnelle. L'usage de cette technique a permis des progrès significatifs, en particulier dans le traitement des images médicales Cepen- dant, cette technique souffre de différents inconvénients L'un d'entre eux est le fait que les surfaces obtenues présentent le plus souvent des trous, ce sur quoi on reviendra plus loin; un autre inconvénient est précisément le fait que cette technique ne permet d'extraire que des surfaces ("iso-surfaces"). La présente invention vient améliorer la situation dans ce domaine du traitement d'images, en particulier d'images médicales. Un premier but de l'invention est d'abord de permettre l'extrac- tion non seulement de surfaces (sans trous) mais aussi de lignes remarquables, telles que des "bi-iso-lignes" Une bi-iso-ligne est la ligne d'intersection de deux iso-surfaces, définies à partir de deux images différentes, ou plus généralement par deux fonctions continues différentes Fl et F 2 de l'espace tridimen- sionnel, rapportées à une ou plusieurs images. Un autre but de l'invention est d'obtenir une représentation des biiso-lignes, à partir d'un maillage polyédrique de l'espace tridimensionnel qui n'est pas nécessairement parallélipipédique, mais peut être quelconque (maillage de type "éléments finis"). Encore un autre but de l'invention est de permettre de détermi- ner des caractéristiques différentielles des iso-surfaces d'une image tridimensionnelle. L'invention a également pour but de permettre la détermination de lignes remarquables telles que les bi-iso-lignes avec un maillage polyédrique particulier, notamment le maillage constitué d'une grille régulière de cubes élémentaires identi- ques, ou celui constitué de tétraèdres complémentaires. L'invention a également pour but, à partir d'une image numérique tridimensionnelle donnée, de calculer des lignes remarquables de celle- ci qui dépendent des caractéristiques différentielles des iso-surfaces qu'elle contient. L'invention a aussi pour but remarquable de permettre la mesure dans une image tridimensionnelle d'un coefficient d'extrâmalité, qui caractérise les lignes de crête des iso-surfaces, et sera défini en détail plus loin. L'invention permet bien entendu d'extraire des iso-surfaces d'une image numérique tridimensionnelle Elle présente à cet égard deux avantages D'une part, elle permet de travailler à partir d'un maillage polyédrique quelconque de l'espace (et non nécessairement à partir d'une grille de cubes identiques, comme dans l'état de la technique) D'autre part, les moyens connus présentent l'inconvénient de fournir des iso- surfaces dans lesquelles existent des trous, inconvénient qui est surmonté par la technique de la présente invention, (hormis naturellement les bords de l'image). Selon une définition générale de l'invention, il est proposé un dispositif électronique pour le traitement d'informations d'image, qui comprend: des moyens de définition d'image propres à fournir, en correspondance d'un maillage polyédrique prédéterminé d'une portion de l'espace, une valeur numérique pour chaque sommet du maillage, des moyens de traitement de polygone, aptes à recevoir une représentation d'un polygone plan (F) à p sommets, p valeurs respectives pour ces p sommets, et un seuil (S), et à fournir, en correspondance du polygone F, une liste de segments orientés reliant des points des arêtes du polygone interpolés comme égaux au seuil, avec une convention de sens prédéterminée, des moyens de traitement de polyèdre aptes à recevoir une représentation d'un polyèdre avec les valeurs en ses sommets et un seuil, puis à solliciter les moyens de traitement de polygone pour toutes les faces du polyèdre, ce qui fournit une liste de segments orientés, puis à associer au polyèdre une liste de cycles fermés et orientés, constitués de segments de ladite liste de segments, et des moyens de traitement principaux aptes à présenter séquentiellement une pluralité de polyèdres dudit maillage aux moyens de traitement de polyèdre, avec les valeurs associées à leurs sommets et un seuil. Les cycles ainsi obtenus appartenant à une iso-surface de l'image pour ladite valeur de seuil. Des variantes intéressantes sont proposées dans le cas de mailles dites cubiques ou tétraédriques. Avantageusement, on y ajoute des moyens pour commander l'explo- ration maille par maille des polyèdres de ladite portion d'espace par les moyens de traitement principaux, soit à volonté, soit en suivant une ligne d'intersection de deux iso- surfaces. De préférence, les moyens de traitement principaux sont aptes à recevoir deux jeux de valeurs pour les points de ladite portion d'espace, ainsi que deux seuils ou iso-valeurs, et à déterminer des lignes caractéristiques par intersection des deux iso- surfaces correspondantes Il est particulièrement intéressant que le second jeu de valeurs soit constitué de grandeurs différentielles relatives au premier jeu de valeurs. Une application intéressante est alors le recalage d'une image par rapport à une autre, en vue de leur comparaison. D'autres caractéristiques et avantages de l'invention apparaî- tront à l'examen de la description détaillée qui va suivre, ainsi que des dessins annexés, sur lesquels: la figure 1 illustre le schéma de principe d'un dispositif de visualisation d'une image tridimensionnelle; la figure 2 illustre le schéma de principe d'un dispositif permettant le rapprochement de deux images tridimensionnelles relatives au même objet mais prises différemment; les figures 3 A à 3 C illustrent un exemple de maillage polyédrique de l'espace; la figure 4 illustre un automate, dit ARETE POINTO, et la figure 4 A est un dessin illustratif de celui-ci; la figure 5 illustre un automate, dit SEGMENTSCYCLE 0, et la figure 5 A est un dessin illustratif de celui-ci; la figure 6 illustre un automate, dit SEGMENTS_POLYGONE(), et la figure 6 A est un dessin illustratif de celui-ci; la figure 7 illustre un automate, dit CYCLES_POLYEDRE(), et la figure 7 A est un dessin illustratif de celui-ci; la figure réalisation la figure la figure la figure la figure la figure la figure la figure obtenu avec la figure la figure la figure la figure 8 est un schéma de principe d'un premier mode de de l'invention; 8 A est une variante de la figure 8; 9 illustre un automate, dit INTERDECALE(); illustre 11 illustre 12 illustre 12 A est une 12 B est un l'automate 13 illustre 14 illustre illustre 16 illustre un automate, dit INTERSEGMENTS(); un automate, dit INTERCYCLES(); un automate, dit SEGMENTSPOLYEDRES(); variante de l'automate de la figure 12; dessin comparable à la figure 7 A, mais de la figure 12 A; un automate, dit SCANMAILLES(); un automate, dit MARCHE LIGNES(); un automate, dit MARCHELIGNE(); un automate, dit MARCHE UNILIGNE(); la figure 17 est un schéma de principe d'un second mode de réalisation de l'invention, adapté au cas d'une grille régulière (un seul type de polyèdre, parallélipipédique); les figures 17 A, 18 A et 18 B sont des dessins utiles à la compréhension de la figure 17; B la figure 19 est un schéma de principe d'un troisième mode de réalisation de l'invention, tandis que la figure 19 A illustre un maillage tétraédrique; les figures 20 A et 20 B sont des dessins utiles à la compréhen- sion de la figure 17; la figure 21 illustre un automate, dit EXTRAITCOURBURES(>; la figure 22 illustre un exemple de mécanisme de recalage d'images, en correspondance de la figure 2; la figure 23 illustre un exemple de mécanisme de rééchantil- lonnage d'image utilisable dans le mécanisme de recalage d'images de la figure 22; les figures 24 et 25 sont deux exemples d'images obtenues à l'aide de la présente invention; et la figure 26 illustre la comparaison de ces deux images. Les dessins sont, pour l'essentiel, de caractère certain En conséquence, ils font partie intégrante de la présente descrip- tion Ils pourront donc non seulement servir à mieux faire comprendre l'invention, mais aussi contribuer à la définition de celle-ci. Par ailleurs, toutes les formules complexes sont reportées en annexe, à la fin de la description L'annexe 1 (plusieurs pages) donne les formules d'intérêt général L'annexe 2 définit les formules et la convention de notation des différentielles d'une fonction F définie dans l'espace tridimensionnel L'annexe 3 définit une structure COURBURES utile à la mise en oeuvre de l'invention. Sur la figure 1, une mémoire d'image MI 1 est reliée à une unité de traitement UT, par exemple du genre "station de travail", capable en outre de piloter un dispositif de visualisation VI, tel qu'un moniteur de haute définition. Sur la figure 2, deux mémoires d'images M Il et MI 2 ont reçu des informations d'images, échantillonnées suivant le même maillage de l'espace, mais relatives à des prises de vue différentes d'un même objet Une première unité de traitement UTRG va réaliser un recalage géométrique des deux images l'une par rapport à l'autre en utilisant leurs lignes remarquables: il en résulte une transformation géométrique qui permet de passer d'une image à l'autre Une seconde unité de traitement UTRI reprend la seconde image MI 2, pour la rééchantillonner dans le même maillage de l'espace que l'image MII, compte-tenu de la transformation géométrique donnée par l'unité UTRG L'image obtenue après rééchantillonnage, directement comparable à l'image MI 1, est stockée dans une mémoire MI 3 L'unité de visualisation VI peut alors réaliser toute visualisation combinée souhaitée des deux images M Il et MI 3, pour un médecin par exemple. L'invention sera tout d'abord décrite dans son application la plus générale, en référence à un maillage polyédrique quelconque M d'une portion de l'espace tridimensionnel On appelle "maille" chacun des polyèdres Il n'est pas nécessaire qu'il existe un motif répétitif, c'est à dire que le maillage soit régulier Les polyèdres sont supposés convexes, bien que des polyèdres concaves puissent opérer, au moins dans certains cas. Dans ce cas général, il est nécessaire de disposer d'un moyen, par exemple une liste en mémoire, définissant chacun des polyèdres convexes constituant le maillage polyédrique de l'espace; il faut y ajouter, pour chaque polyèdre, une défini- tion de chacune des faces du polyèdre. Par exemple, le maillage polyédrique peut être fourni sous la forme d'un fichier contenant la liste des polyèdres composant le maillage; chaque polyèdre est défini à son tour par une liste de polygones plans qui définissent chacune de ses faces; enfin chaque polygone plan est défini par la liste ordonnée de ses sommets, qui sont des points de coordonnées réelles de l'espace tridimensionnel Bien entendu, ce moyen de définition du maillage n'est donné qu'à titre d'exemple, d'autres formats pouvant être envisagés. A titre purement illustratif, la figure 3 A illustre un maillage très simple, constitué d'un prisme supérieur PS et d'un prisme inférieur PI, joints en un plan oblique OB La liste des polyèdres comprend deux mailles ces deux prismes Chaque prisme (le polyèdre) est défini par une liste ordonnée de ses cinq faces, qui sont des triangles ou des parallélogrammes. Chaque face est définie par la suite de ses points, parcourus par exemple dans le sens donné par la règle du tire-bouchon (bien connue en électro-magnétisme, encore appelée règle de la main droite, ou règle du bonhomme d'Ampère), par rapport à la normale N, orientée vers l'extérieur du polyèdre. La figure 3 B illustre les deux polyèdres séparés Elle illustre aussi, purement à titre d'exemple, une numérotation des faces du polyèdre supérieur PS Et la figure 3 C illustre un mode de parcours des sommets du parallélogramme supérieur PSO du polyèdre PS. Chaque noeud du maillage polyédrique est définissable par une représentation numérique de ses coordonnées (X, Y, Z) En ce noeud, on dispose d'une valeur numérique V, tirée d'une mesure, effectuée par exemple à l'aide d'un scanner à rayons X, ou calculée à partir d'une mesure, ou obtenue de toute autre manière sous forme numérique. Ces informations sont stockées dans la mémoire M Il de la figure 1, avec la définition complète du maillage polyédrique. On note généralement S la représentation numérique du seuil ou iso- valeur par rapport auquel seront définies des bi-iso-lignes ou des iso-surfaces. Dans la suite, on définira des listes ou suites, distinguées par des "noms" différents Cette notation fait référence à un ensemble chaîné de cellules de mémoire dont chacune contient une valeur et un pointeur vers la cellule suivante, comme cela est connu dans le langage de programmation dit LISP (acronyme de "LI St Processing") On distinguera ici les listes ouvertes, avec pointeur "NUL" en fin de liste, et les listes fermées ou "cycles" qui rebouclent sur leur point de départ Les cellules sont allouées dynamiquement en correspondance du "nom" de la liste (par exemple: CYCLE) Bien que la mémoire d'image puisse être utilisée pour cela, il est préférable d'utiliser une mémoire de travail. L'invention fait tout d'abord usage d'un premier automate de base 40, dénommé ARETE POINT 0, décrit en référence aux figures 4 et 4 A. Le bloc d'entrée 400 de cet automate reçoit les coordonnées (X, Y, Z) de deux points Pl et P 2 de l'espace tridimensionnel, les valeurs Vl, V 2 associées à ces deux points, ainsi que la valeur de seuil S. Le test 402 détermine si la valeur Vl est inférieure ou égale au seuil S Si oui, on crée ( 403) un signal Ll associé au point Pl et possédant le signe "-" Sinon, on lui donne ( 404) le signe "+" Ensuite, à partir d'une étape 412, on fait de même avec V 2 et le signal L 2 aux étapes 413 et 414. A l'étape 420, si les signaux Ll et L 2 sont de même signe, la sortie directe 495 ne donne rien (pas de point intermédiaire du segment Pl,P 2 qui possède la valeur S) S'ils sont de signes différents, l'étape 430 recherche si Ll possède le signe "-" Si oui, on recherche ( 432) si Vl est proche de S par valeurs inférieures, dans la limite d'une quantité epsilon, faible Si Ll a le signe "+", on fait de même à l'étape 434, mais avec la quantité V 2. il Si l'une ou l'autre de ces deux quantités est proche du seuil à epsilon près, on la diminue de epsilon, respectivement aux étapes 442 et 444. Après cela, l'étape 480 réalise une interpolation, sur la base des valeurs Vl et V 2, entre les points Pl et P 2, par exemple par la formule (XV), annexe 1 On détermine ainsi le point intermé- diaire P du segment P 1,P 2 o la valeur interpolée est égale au seuil S (figure 4 A) A la sortie 490, on associe à ce point P un signal L ayant le signe de L 2. L'homme du métier comprendra que la notation des signes par " 1 + 1 " et "-" est illustrative Elle peut être remplacée par toute autre notation à deux états, comme " O " et " 1 ", par exemple dans la suite, on parlera de signe ou de "label". Il est maintenant fait référence aux figures 5 et 5 A, qui illustrent un automate 70, dit "SEGMENTSCYCL Eoe". L'entrée 500 de l'automate 70 est un CYCLE, défini de la manière illustrée sur la figure 5 A par une suite ordonnée de points (non nécessairement coplanaires), rebouclée sur elle-même L'homme du métier comprendra que le nombre de points de cette suite est pair, mais inconnu dans le cas général d'un maillage polyédrique quelconque Ces points sont notés PCA à PCF sur le dessin (pour l'illustration, l'automate n'ayant pas à connaître cette identification des points) Chaque point est muni d'un signe, qui peut être celui précédemment décrit à propos de la figure 4. L'étape 510 de la figure 5 consiste à examiner le premier point du cycle (PCA) pour tester en 515 si son signe a l'un, prédéter- miné, de ses deux états possibles ("-" est choisi à titre d'exemple) Si le signe de ce premier point est "-", il va demeurer le premier point PCO examiné par cet automate Sinon, l'étape 520 prend le point suivant (PCB), comme premier point courant P. A l'étape 530, on part donc d'un premier point courant P (PCO). L'étape 540 examine le point suivant (PC 1) L'automate note ces deux points Pl et P 2 L'étape 550 consiste à ajouter le segment (Pl,P 2) à une liste de segments notée LISTE L'étape 560 considère le point suivant du cycle comme nouveau point Pi Si on est retourné au point de départ P (PCO), la LISTE complète est passée à l'étape de sortie 580 Sinon, le test 570 retourne à l'étape 540, pour compléter la LISTE. L'effet de ce mécanisme, illustré sur la figure 5 A, est de transformer un cycle fermé de points en une liste de segments ayant chacun pour origine un point dont le signe est "-" (par exemple) et pour extrémité un point possédant l'autre signe "+". Il est maintenant fait référence aux figures 6 et 6 A, pour la description d'un automate 30, dit "SEGMENTSPOLYGONE 0)". L'entrée 300 de cet automate est un polygone F, tel qu'illustré sur la figure 6 A Le nombre des sommets de F est pair mais inconnu, au départ (ici 6 points PFO à PF 5) A ces points sont associées des "valeurs" L'entrée comprend aussi le seuil S. Pour toutes les arêtes de ce polygone F ( 310), on appelle en 320 l'automate 40 ou ARETEPOINT() déjà décrit, avec comme paramè- tres les deux sommets de l'arête, leurs "valeurs", et le seuil S Si les valeurs des deux sommets se placent de part et d'autre du seuil S, l'automate 40 va fournir un point intermédiaire, engendré comme défini à propos de la figure 4 Ce point est stocké dans une liste notée CYCLE (fermée plus tard). Ceci se répète jusqu'à la dernière arête ( 335) De la sorte, on obtiendra par exemple les quatre points PI O à PI 4 illustrés sur la figure 6 A, avec, nécessairement, des signes alternés. L'étape 340 consiste à reboucler le CYCLE sur lui-même pour le fermer On ajoute ainsi un lien entre PI 4 et PI O Le CYCLE est illustré en trait tireté sur la figure 6 A. Après cela, l'étape 345 consiste à appeler l'automate 70 ou SEG- MENTSCYCLE () décrit à propos de la figure 5 Ceci fournit deux segments, de P Il à PI 2, et de PI 4 à PI 0, respectivement, différenciés par un tireté long sur la figure 6 A. A l'étape 350, la liste de (ici deux) segments ainsi engendrés est complétée par le fait qu'on lui associe une référence à la face polygonale F traitée Il en résulte en 360 une LISTE de segments orientés affectés à cette face polygonale, stockée en mémoire. Il est maintenant fait référence aux figures 7 et 7 A, pour la description de l'automate 20, dit CYCLESPOLYEDRE 0). L'entrée de cet automate en 200 est la définition d'un polyèdre tel que 'celui illustré sur la figure 7 A, avec le seuil S A chacun des sommets est associée une valeur mesurée. A l'étape 205, on démarre une boucle qui va traiter toutes les faces du polyèdre jusqu'à l'étape 218. Le corps de la boucle consiste en un appel 210 de l'automate SEGMENTSPOLYGONEO, qui vient d'être décrit Il en résulte pour la face considérée une LISTE de segments, comme rappelé par l'étape 215 Un exemple très simple est donné sur la figure 7 A pour la face supérieure du polygone. La seconde partie de l'automate consiste à analyser les liens entre tous les segments ainsi obtenus en LISTE, pour toutes les faces du polyèdre considéré Le but est de définir un CYCLE, ou plusieurs, deux dans l'exemple de la figure 7 A (l'ensemble sera noté LESCYCLES) Dans chaque CYCLE sont rangés des segments qui se joignent deux à deux et forment ensemble une boucle continue. L'homme du métier comprendra qu'il existe toujours des segments qui se rejoignent, puisque leurs extrémités sont des points interpolés sur des arêtes, correspondant à la valeur de seuil et aux valeurs associés aux points d'extrémité de l'arête concer- née, et ce par le même automate Quelle que soit la face considérée du polyèdre, un tel point sera toujours le même pour l'arête considérée. L'étape 220 détermine si la LISTE des segments est vide Tant qu'elle ne l'est pas, on prend un premier segment SEG 1 dans cette liste à l'étape 225 L'étape 230 consiste à introduire ce segment dans une liste provisoire dénommée CYCLE Après cela, une boucle allant des étapes 240 à 255 considère tous les autres segments de la liste Dès qu'on en trouve un, SEG 2, qui a une extrémité commune avec le segment SEG 1 en cours, à l'étape 245, alors l'étape 250 prend ce segment SEG 2 comme nouveau segment en cours SEG 1 et renvoie à l'étape 225 pour le retirer définitive- ment de la liste considérée, puis à l'étape 230 pour l'ajouter à la liste CYCLE, et ainsi de suite jusqu'à avoir traité tous les segments. Ainsi est obtenu un CYCLE L'étape 260 consiste à ajouter ce CYCLE à une liste globale notée LESCYCLES, après quoi la liste provisoire CYCLE est effacée. L'homme du métier comprendra en effet que d'une part, si l'on a un segment, il existe toujours une suite ou CYCLE de segments, qui va se refermer sur elle-même, à laquelle ce segment appartient; mais l'ensemble des segments que l'on va pouvoir constituer sur un polyèdre donné peut comporter plusieurs cycles de ce genre (fig 7 A). D'o la définition de la liste globale des cycles engendrés notée LESCYCLES et obtenue à l'étape 290. Il est maintenant fait référence à la figure 8 Ce qui est décrit plus haut correspond au rôle des automates 20, 30, 40 et de la figure 8, isolés dans le cadre SRF La figure 8 n'est pas le schéma séquentiel d'un mécanisme, comme c'était le cas pour les figures 4 à 7 Cette figure 8 est un schéma d'interac- tion des automates, les flèches étant dans le sens appelant vers appelé. L'homme de métier comprendra que la liste des segments contenue dans LESCYCLES, telle que finalement fournie à l'étape 290 de la figure 7, appartient à l'iso-surface associée au seuil S, pour la portion de l'espace considéré, et avec l'échantillonnage qui correspond au maillage polyédrique utilisé. Ainsi les blocs du cadre SRF constituent à eux seuls un moyen intéressant pour déterminer une iso-surface. Un avantage essentiel de ce moyen est qu'il évite (étapes optionnelles 432, 434, 442, 444 en particulier) l'existence de trous dans l'isosurface ainsi déterminée Un autre avantage est le fait que les moyens décrits conviennent pour tout maillage polyédrique de l'espace. La figure 8 ajoute d'autres fonctions intéressantes, touchant d'une part à différentes manières d'utiliser ledit moyen ( 10,60,70,80, 90) et d'autre part à différentes manières d'explorer la portion d'espace concernée (Ai; Bl à B 3). Pour obtenir des lignes remarquables, il convient de construire l'intersection de deux iso-surfaces dans l'espace. Sur la figure 8, le bloc MI correspond à au moins une mémoire d'images Plus précisément: ou bien le bloc MI contient une seule mémoire d'image, à laquelle on adjoint (de préférence) des automates convertisseurs dénommés VALEUR 1 () et VALEUR 2 (), capables d'associer deux valeurs différentes à chaque point de l'image considérée L'une de ces valeurs est par exemple la valeur mesurée pour le point de l'image, et l'autre représente une caractéristique différen- tielle de la valeur mesurée en chaque point, déterminée comme on le verra plus loin; le bloc MI peut également contenir deux mémoires d'images (individualisées alors en M Il et MI 2), auquel cas les deux automates convertisseurs VALEUR 1 () et VALEUR 2 () peuvent réaliser simplement une conversion de la valeur mesurée pour chacun des points des deux images dans la représentation numérique appropriée. Ce bloc MI fournit donc à un automate principal 10, également dénommé SEGMENTSPOLYEDRE(>, deux valeurs différentes pour chaque point de l'espace, les deux valeurs étant exprimées dans le même échantillonnage (maillage polyédrique). La figure 9 illustre un autre automate de base, référencé 90 sur la figure 8, et dénommé INTERDECALE(). L'entrée 900 de cet automate est la description d'une face polygonale F d'un polyèdre du maillage, les valeurs associées à chaque sommet de cette face, la valeur de seuil S, et la désignation d'un point B appartenant à l'une des arêtes de ladite face F. L'étape 910 détermine, par ses sommets Pl et P 2, l'arête de la face polygonale F qui contient ce point B Si le signe Ll associé au point Pl est -, comme déterminé par le test 920, l'étape 930 va diminuer la valeur Vi associée à ce point Pi, et ce de deux fois epsilon, epsilon étant le même que précédemment. Dans le cas contraire, l'étape 935 va faire la même diminution, mais à l'autre extrémité pour la valeur V 2 du point P 2. Après cela, il est prévu un appel en 940 à la procédure ARETEPOINTO, pour déterminer un point interpolé à l'intérieur du segment (Pi, P 2) Et ce sont les coordonnées dans l'espace de ce nouveau point, et le signe associé, qui sont retenus pour le point B recherché par cet automate. En correspondance avec la quantité epsilon décrite à propos de la figure 4, celle de la figure 9 a pour avantage essentiel d'assurer une bonne construction des iso-surfaces ainsi que de leurs intersections, en évitant des difficultés lorsque les points extrêmes d'un segment sont proches de la valeur de seuil. La figure 10 illustre un automate 80, dit INTER SEGMENT 50. L'entrée 800 de cet automate est la désignation de deux segments Les extrémités du premier sont Ai et A 2, les extrémi- tés du second sont Bl et B 2 Ils appartiennent à une face polygonale F du maillage polyédrique Et l'on désigne aussi le vecteur normal N à cette face polyédrique, orienté vers l'extérieur. La première étape 805 de cet automate détermine si le point Bl est confondu (ou très proche) avec l'un des points Al ou A 2 de l'autre segment Si oui, on appelle en 810 l'automate INTER DE- CALE() pour ce point Bl. Les étapes 815 et 820 font de même, mais avec le point B 2. Après cela, on procède à la recherche de l'intersection entre les droites-supports (ou axes) des deux segments (Ai, A 2) et (Bl, B 2). Une façon de faire est de calculer en 830 le déterminant DELTA associé au vecteur (A 1,A 2) au vecteur (Bl,B 2), et au vecteur normal N déjà cité. Le test 835 examine si ce déterminant est nul, auquel cas il n'y a aucune intersection et l'on sort, sans résultat, à l'étape 895 Si ce déterminant n'est pas nul, c'est que les droites- supports des deux segments se coupent L'étape 840 détermine alors quel est le point d'intersection P des droites supports des deux segments (Ai, A 2) et (Bl, B 2). Après cela, l'étape 845 détermine si ce point d'intersection est bien situé à l'intérieur des deux segments (Ai, A 2) et (Bl, B 2) eux-mêmes Si cette double condition n'est pas vérifiée, on sort à nouveau à l'étape 895 Si elle est vérifiée, on prend pour ce point P un signe qui est celui du déterminant DELTA, à l'étape 850 Et la sortie de l'automate est la définition du point P et de son signe L, à l'étape 890. Il est maintenant fait référence à la figure 11, pour la description de l'automate 60 de la figure 8, dénommé INTERCY- CLES(. L'entrée de cet automate et la définition de deux cycles notés CYCLE 1 et CYCLE 2 ainsi que du polyèdre POLYED, auquel appartien- nent ces deux cycles. Les deux cycles, CYCLE 1 et CYCLE 2, sont de la nature obtenue en sortie de l'automate 20 décrit à propos de la figure 7. Les deux premières étapes 610 et 620 de l'automate de la figure 11 correspondent au balayage de tous les segments SEG 1 de CYCLEI et SEG 2 de CYCLE 2 respectivement, et ce jusqu'aux fins de boucle 652 et 654. Le corps de ces deux boucles imbriquées comprend en 630 la détermination du fait que les deux segments SEG 1 et SEG 2 appartiennent ou non à la même face F du polyèdre POLYED. Si oui, on appelle en 640 l'automate INTER SEGMENT 50. en 650, le point fourni par cet automate est mis dans une liste que l'on appellera à nouveau CYCLE (L'attention est attirée sur le fait que cette notation est utilisée pour une bonne compréhension des interactions des différents automates. Mais il doit être clair que cette liste CYCLE n'est pas la même que la liste temporaire CYCLE de la figure 7) Après le parcours des deux boucles imbriquées, l'étape 660 effectue un rebouclage de cette liste CYCLE sur lui-même Et on appelle alors en 670 l'automate SEGMENTSCYCLE segments orientés appartient à l'une et à l'autre des iso- surfaces, de seuils respectifs Si et 52 des deux images M Il et MI 2 Ce seront donc des points d'une ligne remarquable, en l'espèce une bi-isoligne. Il reste alors à décrire en référence à la figure 12 l'automate principal 10 de la figure 8, qui gère l'interaction de tous les automates décrits jusqu'à présent. L'entrée de l'automate de la figure 12 à l'étape 100 est un polyèdre convexe appartenant au maillage polyédrique de l'espace, tel que considéré. La première étape 110 consiste à appeler la mémoire MI pour affecter par l'automate de conversion VALEUR 1 () une valeur mesurée ou calculée à tout sommet de ce polyèdre ("première image") Après cela, l'étape 120 appelle l'automate 20 de la figure 8, avec les sous-automates qui en dépendent, pour déterminer une liste LESCYCLE 51, qui correspond à l'iso-surface pour le seuil S, relativement à la première image. Bien entendu, si LES CYCLE 51 est vide, c'est-à-dire que cette isosurface n'existe pas, il est inutile d'aller plus loin, et le test 125 conduit directement à la sortie 190 avec une liste vide de segments orientés d'intersection. Dans le cas contraire, les étapes 130 et 140 font la même chose que les étapes 110 et 120, mais cette fois avec la "seconde image" (mesurée ou calculée) et son automate de conversion VALEUR 2 (). De même, si le test 145 révèle que la liste de cycles dite LESCYCLE 52 ainsi obtenue est vide, on saute à la sortie 190 avec une liste vide de segments orientés. La suite consiste à introduire une double boucle parcourant de à 184 la liste de LES CYCLE 51, et de 160 à 182 la liste LESCYCLE 52 Le corps de ces deux boucles imbriquées consiste en un appel de l'automate INTER_CYCLES() à l'étape 170, avec mémorisation en 180 des segments orientés fournis par cet automate INTER_CYCLES(). Ces deux boucles imbriquées ayant été parcourues, on obtient en une liste de segments orientés qui appartient à une ou plusieurs lignes remarquables communes aux deux images mesurées ou calculées, dans la portion de l'espace considérée, échantil- lonnée avec le maillage polyédrique choisi (On parle ici de "deux images"; on verra plus loin qu'il peut s'agir de deux valeurs ou quantités différentes relatives à une même image). La figure 8 A présente une version simplifiée de la figure 8. Elle se caractérise par le fait qu'elle fait l'économie des automates INTER_CYCLES() (bloc 60, figure 8), INTER_SEGMENTS() (bloc 80, figure 8) et INTER_DECALE() (bloc 90) Cette version simplifiée comporte le même bloc de calcul d'iso-surface SRF que la version de la figure 8, ainsi que les mêmes modes d'explora- * tion de l'espace (A 1; Bl à B 3) Elle diffère cependant en ce qu'elle comporte un seul bloc 1000 ou SEGS_POLYEDRES_BIS() à la place des blocs ( 10, 60, 80, 90) de la figure 8. L'homme de l'art comprendra que cette modification entraîne une grande simplification de la technique mise en oeuvre, avec en contrepartie une moins bonne approximation des lignes remarqua- bles reconstruites - La figure 12 A illustre plus en détail l'automate 1000 ou SEGS_POLYEDRES_BIS() de la figure 8 A Tout comme SEGMENTS_POLYE- DRES(), cet automate prend en entrée une maille polyédrique POLYED (bloc 1005) et retourne ( 1090) la liste des segments orientés correspondant aux bi-iso-lignes passant par ce polyèdre. Les blocs 1005 à 1030 de la figure 12 A sont exactement les mêmes que les blocs 100 à 130 de la figure 8 La différence est que maintenant l'automate CYCLES POLYEDRES() n'est pas relancé pour les valeurs calculées par VALEUR 2 () A ce stade, 1 'automate SEGSPOLYEDREBIS() passe successivement en revue tous les cycles CYCLE 1 de la liste LESCYCLE 51 entre les blocs 1040 à 1084; pour chaque cycle CYCLE 1, l'automate considère successi- vement (blocs 1050 à 1065) toutes les arêtes (Pl,P 2) de ce cycle, ainsi que les deux valeurs correspondantes (V 1,V 2), toutes deux calculées à l'aide de VALEUR 2 () (il est important de remarquer qu'à partir de l'étape 1020, on ne se sert plus des valeurs calculées par VALEUR 1 (>). Pl et P 2 appartiennent par construction aux arêtes de POLYED On calcule donc facilement par interpolation, par exemple par interpolation linéaire, la valeur Vl pour Pl et la valeur V 2 pour P 2, à partir des valeurs aux sommets de POLYED obtenues à l'aide de l'automate VALEUR 2 (). L'automate appelle alors (bloc 1055) l'automate ARETEPOINT() ou (figure 4), qui permet de calculer, le cas échéant, un nouveau point P interpolé entre Pl et P 2, ayant un label L Ce nouveau point, s'il existe, est placé dans une liste intermé- diaire de points appelée CYCLE, au bloc 1060. Après consultation de toutes les arêtes de CYCLE 1, il a donc été constitué une nouvelle liste intermédiaire de points CYCLE, qui est refermée sur elle-même au niveau du bloc 1070 On appelle ensuite l'automate SEGMENTSCYCLE() (bloc 1075) Il transforme CYCLE en une liste de segments orientés Ces segments sont soit envoyés directement pour affichage soit mémorisés dans une liste de segments au bloc 1080 La liste intermédiaire CYCLE est alors effacée, et l'on reboucle ( 1084) jusqu'à épuisement des cycles de la liste LESCYCLE 51 (). La sortie de cet automate est la liste des segments orientés ainsi constitués, si ceux-ci n'ont pas été envoyés directement à l'affichage. Cet exemple est basé sur celui de la figure 7 A, sur lequel on a déjà calculé les cycles correspondant aux valeurs de VALEUR 1 (). On calcule ensuite les valeurs aux sommets Po à PB à l'aide de VALEUR 2 () Puis les valeurs interpolées aux points (P 12, P 24, P 23) et aux points (PF 1, P 03, P 05) En comparant les valeurs obtenues par rapport au seuil 52, on obtient les labels indiqués sur la figure 12 B, et l'on obtient ainsi O segments générés pour le cycle Ci et le segment orienté P 2324, P 1224 pour le cycle C 2. on considérera maintenant l'application des automates décrits ci-dessus, en différentes variantes. Comme déjà indiqué, il existe plusieurs manières de se servir des résultats ainsi obtenus, en particulier soit les visualiser, soit les utiliser pour comparer deux images. Dans le cas d'une visualisation, l'ordre dans lequel les informations sont disponibles importe peu Il n'en est naturel- lement pas de même dans le cas d'une comparaison d'images. Lorsque l'ordre dans lequel on présente les informations graphiques pour affichage sur l'écran de visualisation n'a pas d'importance majeure, il est possible de faire usage de l'automate Ai de la figure 8, lequel est illustré plus en détail sur la figure 13, sous la dénomination SCAN_MAILLE 50. L'entrée A 10 de cet automate est la définition d'un maillage polyédrique d'une portion de l'espace On dispose par ailleurs des automates de conversion VALEUR 1 (), VALEUR 2 () pour convertir sous la forme numérique appropriée les informations de valeurs associées à chaque noeud du maillage polyédrique Il s'y ajoute deux valeurs de seuil Si et 52 Une boucle allant de A 20 à A 50 va pouvoir décrire toutes les mailles de la portion d'espace concernée Pour chaque maille, l'étape A 30 appelle l'automate principal 10 ou SEGMENTSPOLYEDRE 0. Celui-ci fournit une liste de segments orientés, qui est stockée en mémoire à l'étape A 40, ou bien (par exemple) directement transmise à un moniteur pour visualisation. 23 2689272 La sortie A 90 correspond à cette liste de segments orientés, qui peut être maintenant visualisée si elle ne l'a pas été en direct. Dans d'autres cas, un suivi individuel de bout en bout de lignes remarquables peut être requis Il est alors fait usage des automates BE à B 3 de la figure 8 qui seront maintenant décrits en référence aux figures 14 à 16, respectivement. L'entrée B 100 de l'automate BE (figure 14) contient la défini- tion d'une liste rangée de mailles M', correspondant aux mailles de départ des lignes recherchées, définies par rapport à un maillage polyédrique M, avec, à chaque noeud de la liste de mailles du maillage, une définition de deux valeurs par les automates de conversion précités VALEUR 1 () et VALEUR 2 () S'y ajoute la définition de deux seuils Sl et 52. Une boucle principale est établie entre les étapes B 110 et B 180. La première étape de cette boucle consiste (B 120) en l'appel de l'automate principal SEGMENTSPOLYEDRE( O ou automate 10 de la figure 8 Celui-ci fournit, pour la maille en cours, une liste de segments orientés. Une boucle interne s'étend alors entre les étapes B 130 et B 170. Pour chaque segment contenu dans la liste, il est procédé en B 140 à un appel d'une procédure MARCHELIGNE() que l'on décrira plus loin Cette procédure ajoute la ligne qu'elle calcule à une liste de lignes notée LLISTE, et ce à l'étape B 150. Le résultat final en B 190 est une liste de lignes remarquables LLISTE. Sur la figure 15, l'entrée B 200 de l'automate MARCHELIGNE(> reçoit un polyèdre POLYED, deux points Pl et P 2 appartenant à une arête de ce polyèdre, et leurs signes respectifs Li et L 2. L'étape B 210 appelle un autre automate MARCHEUNILIGNE() pour le point P 2 avec son sens ou signe L 2, et ce pour établir une première liste de points constituant une ligne remarquable partant de P 2, notée LISTE 2. L'étape B 220 fait de même, mais avec le point Pl et son sens Ll, pour remplir une autre liste de points, notée LISTE 1, partant de Pl. L'étape B 230 consiste dans un tri en ordre inverse de LISTE 1. L'étape B 240 réalise une concaténation ajoutant LISTE 2 à la fin de LISTE 1 (inversée) Le résultat concaténé est seul traité ci après, sous le nom de LISTE, qui représente la ligne remarquable passant par les points Pl et P 2. Ensuite, l'étape B 250 recherche la dernière cellule de LISTE (en fait, la dernière de LISTE 2), que l'on dénomme LAST L'étape B 260 recherche la première cellule de LISTE (en fait de LIST El inversée), que l'on nomme FIRST L'étape B 270 détermine si les deux points LAST et FIRST appartiennent au même polyèdre. Si oui, on attache le point LAST au point FIRST, et LISTE est fermée Sinon, LISTE est ouverte LISTE est le résultat à l'étape B 290; il représente une ligne qui est soit fermée sur elle-même, soit ouverte. L'homme du métier sait par ailleurs qu'on sait construire en traitement d'image ce qui s'appelle une "table de hachage" ("hashing table") C'est un outil de traitement numérique qui permet de stocker et de retrouver très facilement en mémoire des données associées à des points d'un espace de dimension quelconque. L'automate MARCHEUNILIGNE() sera maintenant décrit en référence à la figure 16 Il porte la référence B 3 sur la figure 8, et les références numériques de la figure 16 commencent par conséquent toutes par B 3. L'étape d'entrée B 300 consiste à recevoir la définition d'un polyèdre POLYED, ainsi que d'un point P avec son signe ou sens L. L'étape B 305 consiste à préparer sous la dénomination PTR une liste, initialement vide. On dispose d'une table de hachage H, supposée vide. Le test 310 va rechercher si le point P donné en entrée se trouve déjà dans cette table de hachage (Cette condition n'est naturellement pas vérifiée au départ) Si elle n'est pas vérifiée, on rajoute le point P à la fin de la liste PTR à l'étape B 320 A l'étape B 322, on le rajoute dans la table de hachage. A l'étape 325, on recherche quelle est la face F du polyèdre POLYED qui contient le point P. En vertu de la description connue du maillage polyédrique de l'espace, l'étape B 330 détermine une autre maille POLYED 2, qui satisfait la condition de partager cette même face F avec le polyèdre POLYED que l'on considère. Le test B 335 détermine si ce polyèdre POLYED 2 est extérieur au maillage M (extérieur à la portion de l'espace traitée). S'il est contenu dans le maillage, on prend ce polyèdre POLYED 2 à la place du polyèdre POLYED concerné, à l'étape B 340, et l'étape B 345 consiste en l'exécution de l'automate principal SEGMENTSPOLYEDRE() sur ce nouveau polyèdre. Ensuite, pour tous les points Pl et P 2 dans la liste engendrée par cet appel de l'automate SEGMENTSPOLYEDRE, on recherche en B 355 si le signe du point est "-" Suivant qu'il est "+" ou "-", l'étape B 362 ou B 364 va rechercher si P coïncide avec Pl ou P 2 respectivement, auxquels cas on prend comme point P, P 2 et Pl respectivement, aux étapes B 366 et B 368. Ceci se poursuit jusqu'à l'analyse du dernier segment de la liste établie par l'étape B 345, pour donner en B 390 une liste PTR finale. Entre-temps, un bouclage a été effectué entre les sorties des étapes B 366 et B 368 et l'étape B 310 Les sorties non traitées des étapes B 310 et B 335 se traduisent par une sortie directe en B 395, avec la liste PTR telle qu'elle existe jusqu'à présent. L'homme du métier comprendra que les automates BE à B 3 ainsi décrits en détail permettent, contrairement à l'automate Ai, d'effectuer une exploration le long des lignes remarquables de la portion d'espace que l'on considère. Bien évidemment, il existe différentes manières de définir l'ensemble des mailles de départ M', à l'entrée B 100 de l'auto- mate de la figure 14 Il est tout d'abord possible de passer toutes les mailles de M en revue séquentiellement On peut également tirer au sort de façon aléatoire les coordonnées d'un nombre donné de mailles de départ MI Il est encore possible de laisser à un utilisateur le soin de désigner M' par tout moyen approprié. Ce qui précède définit, pour un maillage polyédrique quelconque de l'espace: un moyen avantageux ( 20, 30, 40, 70) de déterminer une iso- surface pour une image; un moyen nouveau pour déterminer des lignes remarquables, comme des bi-iso-lignes ( 10, 20, 30, 40, 60, 70, 80, 90); deux modes d'exploration (Al; BE à B 3) d'une ou deux images, utilisables au choix selon l'application. Dans une première variante de réalisation (figure 17), le maillage de l'espace est constitué de cubes identiques, également appelés "voxels". Il en résulte alors des simplifications du cas général qui vient d'être décrit. Tout d'abord, le nombre de cycles et de segments que l'on va engendrer pour chaque voxel ont des valeurs maximum connues Il en résulte la possibilité de gérer la mémoire nécessaire pour la mise en oeuvre de l'invention sans mécanisme d'allocation dynamique de mémoire - L'automate SCANVOXEL() 1400 (Ai, figure 8) ne nécessite aucune gestion dynamique de la mémoire Au contraire, il travaille à l'aide d'un tableau mémorisé de taille fixe, puisqu'il consiste simplement à passer séquentiellement en revue tous les voxels contenus dans la portion d'espace considérée, maillée par une grille tridimensionnelle régulière, par exemple une maille cubique (Le mot cubique s'étend dans la suite à toute maille en forme de parallélépipède; l'homme de métier sait qu'il existe toujours un isomorphisme capable de faire passer d'un parallélé- pipède à un cube). Par ailleurs, dans l'automate 1430 (ou B 3, figure 8), la recherche du voxel voisin au voxel considéré est immédiate. De son côté, l'automate principal SEGMENT VOXEL() noté ici 1440 (correspondant à 10, figure 8), va maintenant engendrer au maximum six segments orientés, d'o là encore des simplifica- tions en ce qui concerne la gestion de la mémoire En bref, il exécute deux fois l'automate CYCLES VOXEL() ou 1445 pour les deux jeux de valeurs et de seuils définis, afin d'engendrer les deux listes de cycles orientés LESCYCLE 51 et LESCYCLE 52. Chaque liste de cycles ainsi obtenue est constituée d'au maximum quatre cycles, ayant chacun au maximum douze sommets De plus, le nombre total des sommets dans une liste de cycles ne dépasse pas douze. Chacun de ces cycles représente un morceau de l'iso-surface définie par interpolation trilinéaire à partir des valeurs des huit sommets du cube Plus précisément, les points de chacun de ces cycles correspondent à l'intersection de cette iso-surface avec les arêtes du cube élémentaire concerné. L'automate CYCLESVOXEL() ou 1445 de la figure 17 prend pour entrée les coordonnées du voxel, ainsi que les huit valeurs associées à ses huit sommets, et le seuil ou iso-valeur S considérée Cet automate va considérer successivement chacune des six faces du cube ou voxel, en constituant à chaque fois un tableau de quatre valeurs, correspondant aux quatre sommets du carré (parallélogramme) que constitue la face du cube (parallé- lépipède) considéré L'ordre des valeurs a de l'importance, et doit être tel que le sens de parcours du carré défini par la règle dite du tirebouchon donne une orientation qui pointe vers l'extérieur du cube L'automate SEGMENTCARRE () ou 1450 est appelé une fois pour chaque face du cube Il renvoie au maximum deux segments, formés chacun d'un couple ordonné de points de l'espace tridimensionnel. Les extrémités de ces deux segments correspondent à l'intersec- tion de l'iso-surface avec les quatre arêtes du cube appartenant à cette face. Le Demandeur a observé qu'en réalité l'intersection exacte de l'isosurface basée sur l'interpolation trilinéaire avec la face du cube est une hyperbole pouvant constituer au plus deux branches séparées En conséquence, il est estimé que les segments reconstruits par cet automate SEGMENT CARRE() soient considérés comme des approximations de ces branches d'hyperbo- les, approximations qui se confondent exactement avec les segments, mais seulement au niveau des arêtes du cube. Après les six appels de l'automate SEGMENTSCARRE() ou 1450, on obtient au maximum douze segments orientés qui permettent, comme pour l'automate CYCLESPOLYEDRE 0, de reconstituer des cycles s'enroulant autour du voxel, par comparaison des extrémités des segments Il existe au maximum quatre cycles distincts. L'automate CYCLESVOXEL() ou 1445 retourne tous les cycles ainsi constitués. On peut encore préciser la nature de l'automate SEGMENTCARRE(> 1450 En effet, elle a pour entrée quatre valeurs Vi à V 4, prises de façon ordonnée, et correspondant aux quatre sommets d'un carré. Il est alors possible de considérer directement tous les cas possibles de distribution des labels ou signes au sommet du carré Les quatre sommets reçoivent en effet chacun un label qui est "-" si la valeur associée au sommet est inférieure ou égale à l'iso-valeur ou seuil considéré S, et "+" dans le cas contrai- re Il existe donc au maximum seize cas possibles en ce qui concerne l'attribution des labels aux quatre sommets d'un carré. Il suffit alors de constituer un mot dit "code" de quatre bits, en utilisant les labels des sommets, dans l'ordre de parcours du carré, du poids le plus fort au poids le plus faible, comme illustré sur les différents dessins de la figure 18. La figure 18 A illustre le sens de parcours d'un carré en vertu de la règle du tire-bouchon déjà citée: lorsque le carré est parcouru dans le sens inverse des aiguilles d'une montre, la normale au carré est supposée orientée vers l'avant de la figure L'ordre des bits dans le sens des poids va de O à 3. Lorsqu'un sommet d'un carré est associé au label ou signe "+", on le représente par un cercle vide On le représente par un cercle plein dans le cas contraire (label "-" ou 0). Selon la figure 18 B, les codes O ou 15 correspondent respective- ment à quatre sommets de label "-" ou quatre sommets de label "+" Dans ce cas, l'automate SEGMENTCARRE(> s'arrête en ne retournant aucun segment. Comme précédemment, l'automate ARETEPOINT() est exécutée pour chaque arête du carré qui présente des sommets ayant des labels opposés Ceci permet d'engendrer au maximum quatre points intermédiaires, correspondant chacun à l'intersection de lliso- surface définie par interpolation trilinéaire avec les arêtes correspondantes du voxel. Il est alors inutile d'appeler l'automate SEGMENTSCYCLE 0, comme on le verra maintenant. En effet, les indications de lien entre les sommets suffisent pour définir la plus grande partie des codes binaires, ou de leur contre- valeur hexadécimale, qui est ici utilisée dans la description. Les cas particuliers concernent les codes 5 et 10, pour lesquels il existe une ambiguïté Cette ambiguïté est résolue en considérant le signe du label L que l'on peut associer à la moyenne des valeurs respectivement existantes aux quatre sommets du carré Ce label est par exemple "-" si cette moyenne V = (Vl+V 2 +V 3 +V 4)/4 est inférieure ou égale à l'iso-valeur S; il est "+" si V est supérieur à S. Le cas du code 5 se décompose donc en deux cas 5 a et 5 b De même, le cas du code 10 se décompose en deux cas 10 a et lob. On obtient donc en tout 18 cas de liaisons possibles à l'inté- rieur d'une face carrée du voxel. L'homme du métier aura remarqué que les différents liens ainsi engendrés sont les mêmes que ceux qui auraient été produits par l'automate SEGMENTSCYCLE(), sauf pour les cas 5 et 10 qui viennent d'être examinés Ces deux cas constituent une caracté- ristique remarquable de l'automate SEGMENTSCARRE() par rapport à SEGMENTSPOLYGONE( O La modification qu'ils impliquent a notamment pour but de réaliser un maillage polygonal de l'iso- surface correspondante, et ce avec la même topologie que l'iso- surface définie par l'interpolation trilinéaire En effet, si l'on change l'orientation de la face F, les liens constitués restent les mêmes en changeant simplement de sens. Enfin, il est intéressant d'associer à chaque arête une référence à la face F, cette fois-ci sous la forme d'un code noté ORI, compris entre 1 et 6, qui désigne l'orientation de la face du cube qui contient les segments concernés. L'automate SEGMENTCARRE() fournit alors comme résultat l'ensemble de ces segments orientés, qui sont au maximum au nombre de deux. Dans la partie droite de la figure 17, sous l'automate 1440, on pourrait là aussi utiliser les mêmes automates que sur la figure 8 Mais des simplifications sont intéressantes. Au lieu de INTERCYCLE(>, une variante intéressante dénommée INTERCYCLECUBEO) consiste à prévoir la fourniture d'au plus six segments orientés, renvoyés comme résultat direct à l'automate SEGMENTSVOXEL(>. En outre, il est possible de remplacer l'automate INTER SEG- MENT() par INTER_HYPERBOLE(O ou 1490. Ce nouvel automate prend pour entrée deux segments orientés SEGMENTA et SEGMENTB, formés de deux couples ordonnés de points (Al, A 2) et (Bt, B 2), ainsi que deux listes ordonnées des quatre valeurs (V 11, V 12, V 13, V 14) et (V 21, V 22, V 23, V 24), associées à la face F qui contient ces segments, et correspon- dant respectivement aux deux jeux de valeurs et aux deux seuils. INTER-HYPERBOLE() utilise encore l'indice "ORI" qui indique l'orientation de la face F considérée dans le cube ou voxel, parmi les six orientations possibles Cet indice ORI permet de retrouver immédiatement quelle est la direction N de la normale associée à la face F. Pour éviter les cas dégénérés d'intersection, un appel à l'automate INTERDECALE() permet le cas échéant de déplacer légèrement les points Bt et/ou B 2 en cas de coïncidence avec At ou A 2. Ensuite, on calcule exactement les intersections des deux hyperboles définies par interpolation bilinéaire à partir des valeurs (Vit, V 12, V 13, V 14) et (V 21, V 22, V 23, V 24) (Dans le cas de parallélépipèdes quelconques non rectangles -, on applique l'interpolation hyperbolique en se ramenant localement à une forme carrée, par changement de repère l'isomorphisme précité). Une façon de calculer ces intersections consiste à résoudre le système d'équations donné en Annexe I, formule (XVIII) Il est également fait référence à la figure 17 A. Les coordonnées u, v sont les coordonnées des points d'intersec- tion, dans le repère bidimensionnel associé à la face F Elles sont obtenues en éliminant le produit uv du système d'équation (XVIII) par multiplication de la première équation par V 23 +V 21- V 24-V 22, et de la seconde par V 13 +V 11-V 14-V 12, et en construi- sant la différence des deux équations ainsi multipliées Il en résulte la connaissance de v en fonction de u On remplace v par cette valeur dans l'une des deux équations, ce qui donne une équation du second degré en u, pouvant conduire à deux solutions d'intersection, ce qui fournit au plus deux points d'intersec- tion possibles P et P 2. Si ce calcul fournit zéro ou une solution réelle, alors l'automate INTERHYPERBOLE() s'arrête en ne retournant aucun point nouveau Si elle fournit deux solutions réelles, on affecte à chacun de ces points Pl et P 2 un label ou signe qui correspond au signe du déterminant de la matrice formé: par le vecteur tangent, en le point Pl ou P 2 considéré, sur l'hyperbole définie par Vll, V 12, V 13, V 14, par le vecteur tangent en ce même point Pl ou P 2 à l'hyperbole définie par V 21, V 22, V 23 et V 24, et par la normale N à la face considérée. L'automate INTERHYPERBOLE renvoie alors les points d'intersec- tion ainsi calculés, qui appartiennent bien aux parties des branches d'hyperboles définies d'un côté par (V 11, V 12, V 13, V 14) et de l'autre par (V 21, V 22, V 23, V 24), limitées respecti- vement par les segments Al, A 2 et BI, B 2, et ce dans leur ordre d'apparition le long de la branche d'hyperbole (AI, A 2). On observera que la face F est susceptible d'appartenir à deux voxels différents de la grille tridimensionnelle formant le maillage de base Toutefois, comme pour le cas de l'automate INTERSEGMENT(), le résultat est indépendant de l'orientation de la face F, à l'exception des labels associés aux nouveaux points créés, et de leur ordre d'apparition, qui s'inverse. La position de ces nouveaux points est la même dans les deux orientations possibles de la face considérée Ceci garantit la continuité des bi-iso-lignes reconstruites. L'explication de cela réside dans le fait que les intersections sont calculées à partir de valeurs (Vll, V 12, V 13, V 14) et (V 21, V 22, V 23, V 24) associées aux sommets de la face, et qui ne dépendent pas de l'orientation de cette face. Enfin, contrairement au cas rencontré en construisant l'inter- section de segments, les points calculés par l'automate INTERHYPERBOLE () appartiennent exactement aux bi-iso-lignes définies par interpolation trilinéaire sur les voxels. Finalement, la procédure INTER HYPERBOLE () retourne donc 0, 1 ou 2 points d'intersection, munis chacun d'un label "-" ou "+". Une autre variante de l'invention utilise un maillage tétraé- drique de l'espace. La figure 19 A montre comment une maille de base à cinq tétraè- dres peut représenter un parallélépipède, dont l'aptitude à remplir complètement l'espace sera immédiatement comprise. Le mécanisme global pour extraire les lignes remarquables est alors illustré sur la figure 19 Les blocs Ai, Bi, B 2, B 3 et MI, ainsi que 40 et 90 de celles-ci peuvent être les mêmes que dans le cas général d'un maillage polyédrique. L'automate SEGMENTTETRAEDRE() 1640 est par contre modifié par rapport à SEGMENTPOLYEDRE 0). En effet, avec un polyèdre de base tétraédrique, cet automate 1640 fournira au maximum deux segments orientés A cet effet, l'automate 1640 fait appel deux fois à un automate 1645 ou CYCLETETRA() Celui-ci est une version simplifiée de 1 'automate 20 de la figure 7: en effet, il prend en entrée un tétraèdre TETRA, une liste de valeurs associées chacune à l'un des sommets du tétraèdre, et un seuil ou iso-valeur S. Le sous-automate appelé, au lieu de SEGMENTPOLYGONE() de la figure 6, se simplifie en SEGMENTTRIANGLE() Cet automate a pour entrée une face triangulaire du tétraèdre (trois points et les trois valeurs associées) et d'un seuil ou iso-valeur S Il renvoie zéro ou un segment orienté On comprend dans ce cas qu'il n'est nullement nécessaire de faire appel à l'automate 70 pour parcourir ces segments. La première étape de cet automate SEGMENTTRIANGLE(> consiste à calculer un code de trois bits pour le triangle, de la même manière que cela était décrit pour le cas du carré précédemment, mais avec seulement trois sommets Il ne reste plus que les huit cas possibles de la figure 20 B Le parcours du tétraèdre est illustré sur la figure 20 A, avec les mêmes conventions de représentation que précédemment (dans la figure 18). Sur la figure 20 B, les cas O et 7 (en hexadécimal) sont ceux o aucun segment n'est généré Dans les autres cas, l'automate ARETE-POINT() 40 est appelé seulement pour deux des trois arêtes Cela suffit en effet pour déterminer le segment recherché. Une variante consiste à examiner directement les 16 cas possibles pour la distribution des labels aux quatre sommets du tétraèdre, auquel cas il n'est plus du tout nécessaire d'appeler la procédure 40 ou ARETEPOINT 0. Dans la partie droite de la figure 19, des simplifications interviennent également En effet, l'automate INTERCYCLES() référencé 60 sur la figure 8 peut être remplacé par une variante 1670 sur la figure 19, notée INTERBICYCLE Soi Cette variante tient compte du fait que le cycle engendré en CYCLE ne pourra posséder que 0, 2 ou 4 sommets (extrémités de segments) Elle fait l'économie de l'appel à l'automate SEGMENTSCYCLE 0. Dans le cas sans aucun sommet, il ne peut y avoir aucun segment. Si CYCLE a deux sommets, le segment correspondant, orienté du label "-" vers le label 11 + 11, est renvoyé comme résultat Si CYCLE possède quatre sommets, alors les segments sont constitués suivant les règles de liens du carré, telles que définies sur la figure 18 B, et en utilisant la moyenne des quatre valeurs associées pour lever l'ambiguïté des cas de code 5 et 10, comme déjà décrit. On s'intéressera maintenant à une variante particulièrement intéressante de l'invention. Dans ce qui précède, les automates de conversion VALEUR 1 () et VALEUR 2 () étaient supposés traiter des données relatives à deux images différentes. Dans ce qui va être maintenant décrit, l'un des deux automates continue à fournir une grandeur mesurée relative à une image. L'autre automate fournit une grandeur différentielle relative à ladite grandeur mesurée. Un automate utilisable à cet effet est celui dénommé EXTRAIT= COURBURES(), qui sera maintenant décrit en référence à la figure 21. L'étape d'entrée 1800 de cette figure 21 consiste en la réception d'une image, dont est désigné un point P, avec en outre une direction privilégiée de l'espace, représentée par exemple sous la forme d'un vecteur unitaire T dans l'espace tridimensionnel. La première étape de l'automate de la figure 21 est une étape dite DIFFERENTIELLES(, notée 1810 Elle consiste essentielle- ment en le calcul des dérivées successives jusqu'à l'ordre 3, pour l'image considérée, au point P. La notation et la définition de ces dérivées jusqu'à l'ordre 3 sont données dans l'annexe II à la présente description (la suite de l'annexe I sera discutée plus loin) On l'appellera "structure DIFFER" dans la suite. Pour calculer ces dérivées, un premier moyen utile consiste en une fonction d' interpolation INTERPOLE (), capable de définir les valeurs d'un point de l'image n'ayant pas nécessairement des coordonnées entières correspondant au maillage dans lequel est définie cette image Une telle fonction d'interpolation permet, à partir de l'image tridimensionnelle disponible, d'accéder à une fonction continue F de l'espace tridimensionnel, dont tous les points passent par les points de l'image telle qu'échantil- lonnée. Comme fonction d'interpolation, on peut utiliser une technique de filtrage linéaire Dans ce cas, la valeur V de la fonction F telle qu'interpolée en un point P correspond à la convolution de l'image tridimensionnelle IMAGE par une fonction tridimension- nelle de lissage G, définie à l'avance. Dans ces conditions, pour obtenir les différentielles au degré voulu de la fonction F au point P, il suffit d'effectuer la convolution de l'- image avec la différentielle correspondante de la fonction de lissage G choisie. Une description de cette technique se trouve dans l'ouvrage "COMPUTER VISION", Dana H BALLARD et Christopher M BROWN, Prentice-Hall INC, Englewood Cliffs N J 07632, 1982. En pratique, il est intéressant, pour les points P appartenant à l'image reçue IMAGE, de précalculer les 20 sous-images tridimensionnelles, de taille fixée par les utilisateurs, correspondant aux 20 différentielles possibles de la fonction G, de l'ordre O à l'ordre 3. De telles sous-images sont généralement appelées masques de convolution. Le calcul d'une différentielle donnée de la fonction F au point P, par exemple FXXY, avec les notations de l'annexe 2, s'effec- tue en calculant la somme des valeurs de l'image IMAGE autour du point P, pondérée par les valeurs de la sous-image correspondant à l'échantillonnage de la fonction GXXY au même point. * En règle générale, c'est une fonction gaussienne que l'on utilise pour fonction G, avec l'expression suivante: G = exp { ( x 2 + y 2 + z 2)} o exp{} désigne la fonction exponentielle. Encore une autre manière utile de réaliser l'étape 1810 de calcul des différentielles de la fonction F consiste à calculer une approximation de l'image tridimensionnelle IMAGE, sous la forme d'une fonction "SPLINE" tridimensionnelle notée F, constituée du produit tensoriel de trois fonctions polynomiales par morceau. Les différentielles de la fonction F sont alors des produits tensoriels de fonctions polynomiales, d'ordre inférieur, faciles à calculer à partir de la fonction F elle-même. Un exemple de réalisation de cette technique est décrit dans "Application des fonctions SPLINE au traitement d'images numériques", Philippe CINQUIN, Thèse de doctorat d'état, spécialité Mathématiques, Université scientifique, technologique et médicale de Grenoble, septembre 1987. "Tricubic spline functions, an efficient tool for 3 D images modelling", Philippe CINQUIN, Abstract congress ICIAM-1987, International conference on industrial and applied mathematics, Paris, juillet 1987. Après cette étape 1810 que l'on vient de décrire en différentes variantes, l'étape 1820 consiste en un automate dit CALCULE_- COURBURE() Cet automate reçoit en entrée la structure DIFFER, dont les constituants sont les différentielles de la fonction F jusqu'à l'ordre 3, telles que définies en annexe II. La première étape de l'automate CALCULE COURBURES() consiste à calculer le gradient GRAD de la fonction F, qui est le vecteur défini par la formule (I) de l'annexe I Sur ces formules, il est fait référence à l'ouvrage "Differential Geometry of Curves and Surfaces", Manfredo P DO CARMO, Prentice-Hall, 1976. On calcule ensuite sa norme ou module, élevé au carré, selon la formule (II) Les éléments de calcul de sa courbure gaussienne K sont définis par la formule (III) La courbure moyenne A est à son tour définie par la formule (IV). La formule (V) définit le mode de calcul du discriminant D associé à ces deux courbures, qui sont données par une équation du second degré (obtenue à partir de leur demi-somme A et de leur produit K déjà définis). Enfin, les formules (VI) et (VII) indiquent comment l'on peut calculer les courbures principales Kl et K 2 de la surface définie par l'équation implicite F(x,y,z) = S au voisinage du point P, avec F(P) = S. Quoique certaines des formules qui viennent d'être invoquées aient déjà été décrites, le Demandeur estime que, dans leur ensemble, ces formules constituent un moyen nouveau pour accéder à des caractéristiques différentielles d'une image tridimension- nelle Les formules (VIII) à (XIV) venant ci-après sont pour leur part considérées comme inédites. A partir de là, la courbure principale Tl est définie sous forme de direction par un vecteur, qui lui-même dépend de deux autres vecteurs notés ALPHA et BETA. Les trois composantes du vecteur ALPHA sont données par les formules (VIII) Les trois composantes du vecteur BETA sont définies par les formules (IX). Les trois composantes du vecteur Tl sont définies par la formule (X) Elles correspondent à la différence vectorielle entre le vecteur ALPHA et D fois le vecteur BETA, d étant le discriminant déterminé par la formule (V) S'il n'est pas nul, ce vecteur Tl définit la première direction principale de courbure; il est alors normé, c'est-à-dire divisé-par son module. De même, il est calculé un vecteur T 2 égal à la somme du vecteur ALPHA et de D fois le vecteur BETA, et ce vecteur T 2 est à son tour normé. Si l'un des deux vecteurs principaux Tl et T 2 est nul, alors il peut être remplacé par le produit vectoriel de l'autre vecteur principal (non nul) et du gradient GRAD. Si les deux vecteurs sont nuls, alors un drapeau ("flag") ERREUR est mis à 1. Un élément original des formules qui viennent d'être décrites pour le calcul des directions principales est que les axes principaux de l'espace (X, Y, Z) jouent un rôle symétrique dans les formules donnant les deux vecteurs principaux Tl et T 2 Plus précisément, les formules de leurs coordonnées suivant Y et suivant Z sont obtenues par permutation circulaire de (X, Y, Z) dans la formule décrivant la coordonnée suivant X. On s'aperçoit que jusqu'à présent seules les différentielles jusqu'à l'ordre 2 de la fonction F ont été utilisées. Après avoir ainsi déterminé les vecteurs de direction principale Tl et T 2, il est procédé au calcul des différentielles suivant x, y et z des courbures principales en Kl et K 2 Ceci s'effectue par différenciation des formules (VI) et (VII), d'ou un vecteur GRADK 1 = (Kl X, Kl Y, Kl Z) et un autre vecteur GRADK 2 = (K 2 X, K 2 Y, K 2 Z), comme indiqué par les formules (XII). Maintenant, on a fait intervenir les différentielles de la fonction F jusqu'à l'ordre 3 inclus. Une autre caractéristique particulièrement originale de l'inven- tion consiste à établir des coefficients d'extrêmalité, qui caractérisent les lignes de crête présentes dans une image tridimensionnelle. Un premier coefficient d'extrêmalité El est obtenu en réalisant le produit scalaire du vecteur GRADK 1 et du vecteur Tl, conformément à la formule (XIII). Un autre coefficient d'extrêmalité E 2 est obtenu de même par produit scalaire du vecteur GRADK 2 et du vecteur T 2, conformé- ment à la formule (XIV). Les signes des coefficients d'extrémalité El et de E 2 sont directement liés à l'orientation dans l'espace des vecteurs Tl et T 2. Cette technique de travail fonctionne dans la quasi-totalité des cas Il existe cependant des exceptions très rares o il n'est pas possible de calculer les vecteurs de direction principale Tl et T 2 C'est le cas par exemple quand le vecteur gradient GRAD est un vecteur nul, ou bien pour des points dits "ombiliques" (Voir l'ouvrage de DO CARMO) correspondant au cas o les deux courbures principales KI et K 2 sont égales Dans ces cas, cette impossibilité est matérialisée en mettant le drapeau ERREUR de la structure courbure à la valeur 1 Dans tous les autres cas, o les calculs se sont bien déroulés, ce drapeau est mis à zéro. Quoiqu'il arrive, l'automate CALCULECOURBURE So( fournit une structure COURBURE ainsi remplie, ou bien possédant un drapeau ERREUR à 1 La structure COURBURE est illustrée dans l'annexe 3. On s'intéressera maintenant à une variante particulièrement importante de l'invention Elle concerne l'extraction de lignes d'une image, dans le cas particulier o l'on s'intéresse à des lignes de crête. Les lignes de crête d'une iso-surface d'une image tridimension- nelle sont définies comme étant le passage par zéro du coeffi- cient d'extrêmalité pour cette surface Les lignes de crête sont donc les bi-iso-lignes dont l'une des iso-surfaces est définie à partir de l'image tridimensionnelle elle-même, et d'un seuil ou iso-valeur S fournie par l'utilisateur; l'autre iso-surface est calculée avec comme VALEUR 2 () le calcul du coefficient d'extrêmalité pour la même image, et un seuil ou iso-valeur associé qui est zéro. En conséquence, on utilise d'abord l'automate MARCHELIGNE() pour une grille tridimensionnelle régulière, avec comme automate de conversion VALEUR 1 (), la lecture en mémoire du voxel corres- pondant de 1 ' image tridimensionnelle (éventuellement convertie), et comme premier seuil Si, un seuil S choisi. L' automate VALEUR 2 () est alors en fait 1 ' automate servant au calcul de l'extrêmalité El par l'automate CALCULECOURBURE 0, avec un seuil 52 = 0. Cependant, cet automate CALCULECOURBURE() nécessite, pour un fonctionnement optimal, qu'il lui soit définit une direction privilégiée de l'espace T Cette direction ne peut être fixée une fois pour toutes pour tout l'ensemble de la mise en oeuvre de l'invention sur une image donnée; en effet, si on procède ainsi, cela entraîne des coupures indésirables le long des lignes de crête. En conséquence, il est estimé préférable d'orienter le calcul de l'extrêmalité à l'aide de la direction de la ligne de crête que l'on est en train de suivre. Ceci peut être réalisé en modifiant légèrement les automates MARCHEUNILIGN Eo) et SEGMENTSPOLYEDRE( 0, de sorte que ces automates acceptent en outre comme entrée l'indication d'une direction privilégiée T dans l'espace. En ce qui concerne l'automate MARCHEUNILIGNE 0, définir cette direction T correspond à passer la direction (P 1,P 2), en tant que paramètre, au niveau du bloc B 210 de la figure 15 De même, on passera la direction (P 2,P 1) au niveau du bloc B 220 de la même figure. Si aucun segment n'a encore été déterminé, on conserve la direction T qui a été passée comme argument d'entrée de l'automate MARCHE_UNI_LIGNE(>. Au moment o l'automate MARCHEUNILIGNE() appelle l'automate SEGMENTS_POLYEDRE 0, celui-ci va utiliser la direction T au niveau du bloc 130 de la figure 12, en passant la direction T à l'automate VALEUR 2 (), qui, dans ce cas, est l'automate EXTRAIT=- COURBURES() retournant l'extrémalité El. L' homme de métier aura également remarqué que 1 'automate SEGMENTSPOLYEDRE () est aussi utilisé dans 1 ' automate MARCHE LI- GNES(), au niveau du bloc B 120 de la figure 14 Dans ce cas, on peut se contenter d'une direction T arbitraire. Par ailleurs, dans ce mode de travail, on pourra considérer que l'exécution de l'automate MARCHE UNI LIGNE() est terminée au niveau du bloc B 345 de la figure 16, quand il n'est pas possible de calculer le coefficient d'extrêmalité pour l'un des sommets du voxel considéré Dans ce cas, concrétisé par le fait que le drapeau ERREUR du champ COURBURE vaut 1, on sort par le bloc B 395 de cette même figure 16, en retournant la liste PTR déjà constituée jusqu'à présent. L'invention s'applique également pour déterminer des lignes "iso-gaussiennes" Dans ce cas, l'entrée est une image tridimen- sionnelle, accompagnée de deux seuils Si et 52 L'automate VALEUR 1 () correspond simplement à la lecture en mémoire de l'image tridimensionnelle pour déterminer la valeur du voxel correspondant Si est le seuil d'iso-valeur définissant l'iso- surface dans l'image tridimensionnelle choisie pour l'utilisa- teur Le second automate d'entrée VALEUR 2 () correspond au calcul de la courbure gaussienne K par l'automate EXTRAITCOURBURE() déjà cité, et le second seuil 52 est la valeur de l'iso- gaussienne de l'iso-surface que désire l'utilisateur Les automates SCANMAILLE() ou MARCHELIGNE(> n'ont pas à être modifiés, car le calcul de la courbure gaussienne ne dépend pas d'une orientation dans l'espace, comme c'était le cas pour le coefficient d'extrêmalité des lignes de crête. Bien entendu, aucun caractère limitatif ne doit être attaché aux deux cas précités d'utilisation conjointe de la technique d'extraction de ligne et des moyens de calcul des caractéristi- ques différentielles des iso-surfaces, le tout sur une même image. On peut en envisager de nombreuses variantes Chaque variante correspond au cas o l'on modifie le mode d'élaboration des automates d' entrée VALEUR 1 () et VALEUR 2 () qui peuvent S ' attacher soit aux données d'images elles-mêmes, soit à leurs caractéris- tiques différentielles, soit à toute autre information pouvant avoir un rapport avec l'image considérée. Une autre application importante de l'invention sera maintenant décrite en référence aux figures 22 et 23. L'invention permet en effet le recalage et le rééchantillonnage de deux images tridimensionnelles IMAGE_ 1 et IMAGE_ 2, qui sont relatives à un même sujet, mais pris dans deux positions différentes Ce peuvent être des images médicales comme celles d'un scanner à rayons X ou celles que l'on obtient en résonance magnétique nucléaire, par exemple. Le schéma de principe le plus général en est donné sur la figure 2 déjà mentionnée. Les images IMAG El et IMAGE 2 sont stockées dans des mémoires MI 1 et MI 2, et l'on a prévu en attente une mémoire MI 3 pour stocker une image intermédiaire. Sur chacune des deux images M Il et MI 2, on a procédé à la mise en oeuvre de l'invention telle que décrite plus haut, afin d'extraire des lignes tridimensionnelles remarquables de l'image, comme par exemple ces lignes de crête Ces automates sont schématisés en E 100 et E 200 sur la figure 22 Il en résulte des données, dans des mémoires de lignes, respectivement stockées en M Ll et ML 2. A partir de là, une autre unité de traitement va procéder au recalage automatique des courbes dans l'espace tridimensionnel que représentent ces deux jeux de lignes Il s'agit de l'étape E 300. Son but est de fournir en sortie une transformation géométrique propre à permettre le passage de l'un des jeux de lignes à l'autre On suppose ici qu'il s'agit de la transformation qu'il faut appliquer aux informations LIGNE_ 2, pour les superposer (sensiblement) au jeu de données LIGNES-1. De façon générale, une telle transformation peut être décrite sous la forme d'une matrice de taille 4 x 4, qui donne la transformation en coordonnées homogènes existant de l'un à l'autre des jeux de données Plus précisément, si l'on note Pl: (Pi X, Pl Y, Pi Z, 1) le vecteur en coordonnées homogènes qui représente la position d'un point dans l'image IMAGE_ 1, et P 2 la position correspondante dans l'image P 2, la matrice M est définie alors par: P 2 = M x Pl Une technique utilisable pour le recalage de lignes tridimen- sionnelles est décrite par exemple dans le document "Smoothing and Matching of 3-D Space Curves", André GUEZIEC et Nicolas AYACHE, Rapport de Recherche INRIA No 1544, octobre 1991. La matrice ainsi déterminée est obtenue à l'étape 400. On peut alors, de différentes manières, s'en servir pour rappro- cher l'une de l'autre, automatiquement ou avec l'intervention d'un opérateur, les informations contenues dans les deux mémoires d'images M Il et MI 2. On décrira ci-après un exemple d'un tel rapprochement automati- que. Selon cet exemple, l'étape E 500 consiste à rééchantillonner l'image tridimensionnelle MI 2, compte tenu de la transformation définie par la matrice M. Sur la figure 23, l'étape d'entrée 2000 consiste à recevoir la matrice de transformation M, ainsi qu'en Tx, Ty et Tz des données fixant la taille de 1 'image tridimensionnelle désirée en sortie On a accès aux données dans la mémoire MI 2 qui définis- sent l'image IMAGE_ 2. L'étape 2010 de l'automate de rééchantillonnage d'image tridimensionnelle correspond à la préparation d'une image de sortie notée IMAGE_ 3, et qui est pour l'instant vide avec une taille Tx * Ty * Tz Cela signifie qu'on prépare un tableau tridimensionnel de valeurs, tableau dont les trois paramètres de taille sont Tx, Ty et Tz, et dont toutes les valeurs sont pour le moment considérées comme nulles. A l'étape 2020 de l'automate, pour tous les voxels possédant des coordonnées homogènes Pl:(Px,Py,Pz,1), de l'image en attente IMAGE_ 3, on effectue une boucle qui va jusqu'à avoir parcouru tout le tableau à l'étape 2080. Pour chacun de ces points à calculer de IMAGE_ 3, l'étape 2030 consiste à déterminer par la relation P 2 = M Pl, quel est le point correspondant dans l'image 2. A l'étape 2040, on détermine par interpolation linéaire dans l'image IMAGE_ 2, au point P 2, quelle est la valeur correspondant à ce point (qui a priori n'a aucune chance de tomber exactement sur un point déjà échantillonné dans IMAGE_ 2). L'étape 2050 affecte cette VALEUR qui vient d'être calculée à la case Pl de l'image IMAGE_ 3. En fin d'automate, on a à l'étape 2090 rempli la mémoire M 13 avec tous les points de IMAGE_ 3. Ceci étant réalisé, l'une des applications les plus simples (figure 22) consiste à effectuer sur un outil de visualisation approprié une présentation simultanée des deux images 3 D à l'attention d'un médecin, par exemple par couleurs différentes, ou mieux par un double curseur, les images étant disposées côte à côte. EXEMPLE PARTICULIER L'invention a été mise en oeuvre sur une station de travail modèle DEC 5000 (Digital Equipment) Les automates décrits ont été réalisés en langage C, qui permet l'allocation dynamique de mémoire, et l'utilisation de pointeurs. Une image tridimensionnelle de 64 * 64 * 64 pixels représentant une vue tomographique à rayons X d'un crâne a été prise comme IMAGE_ 1, en maillage cubique. Le résultat (projeté en noir et blanc sur un plan et simplifié par filtrage pour être lisible en plan) est donné sur la figure 24 Malgré cette simplification, on y distingue nettement les grands traits d'un crâne humain Il est à remarquer que la définition des lignes de crête est largement meilleure que celle de l'image initiale. Dans des conditions analogues, la figure 26 montre comment l'on peut superposer automatiquement (à l'aide des lignes de crête et après retraitement de l'une d'entre elles) deux images (figure 24 et 25) du même sujet prises dans des conditions différentes. Ceci était totalement inaccessible avec les moyens de l'état de la technique. Les moyens proposés selon l'invention sont extrêmement puis- sants, en ce sens que sans aucune intervention manuelle, le médecin peut par exemple: disposer des lignes de crête d'une image tridimensionnelle (l'opération était jusqu'à présent effectuée en dessin à la main, par exemple en chirurgie faciale); ou disposer d'une comparaison de deux images du même sujet obtenues dans des positions différentes, et ramenées à la même position. En d'autres termes, on obtient ainsi une chaîne de traitement automatique qui permet de prendre deux images médicales d'un même sujet, mais acquises dans des positions différentes, et de les transformer toujours en deux images du même sujet, mais ramenées à la même position géométrique Ceci facilite considé- rablement le diagnostic du médecin Par exemple, s'il s'agit de deux images du même sujet prises à des temps différents, il peut suivre l'évolution de symptômes d'une maladie On peut également comparer le même organe de deux patients différents On peut aussi suivre le déplacement d'un organe au cours du temps. Dans une autre application, la technique de recalage automatique est utilisée pour comparer et recaler une image médicale tridimensionnelle avec un atlas anatomique, afin d'attribuer un nom à chacune des parties de l'organe examiné, et ceci toujours de manière automatique. La comparaison d'images relatives à deux sujets différents a également un sens Après recalage selon la technique décrite, on peut observer les écarts entre lignes qui se correspondent, et calculer des déformations locales entre ces deux jeux de lignes (recalage avec déformation). Frmu Xe 1-1 Formules ( 1) GRAD: ( FX, FY, FZ); (II) H FX 2 + FY 2 + FZ 2; (III) K l FX 2 'FYY*zz E'Y 2 (FXOCFZZ FZ 2 FYZ 2) + 2-FY-FZ-(FXZ-FXY FXX-FYZ) + FXZ 2) + 2-FX-FZ-(FYZ-FXY FYY-FXZ) + FXY 2) + 2-FX-FY-(FXZ-FYZ FZZ'FXY) l / H 2; (IV) A l FX 2-(FYY + FZZ) 2-FY FZ-FYZ + FY 2-(FXX + FZZ) 2-FX FZ*FXZ + FZ (FXX + FYY) 2-FX-FY-FXY l / ( 2 H 3/2); (V) D A 2 K; (VI) Ki A + Dt/2; (VII) K 2 = A Dl/2; :Vs T ci +VH = (Zz I Z'x Izi) Zi(IX) vsse-aV Hralv7 7 'xl 1) T (X) (XS a S IZA XS'a -x zq): y 1 se (XI) ( /w H'z) / l AS'-X + ZZS-2;-z X - xx S a L* xm Z ZX 7 *z X Z R Z + 7 ai zz xm ZXS* z A Z& Z + XX.-zzz- 7 + XX z Z -x -z XXS-EXA + 7 X- z X - e-Z + ZZ -z XZ- X + Xxs-z L-Xzxa-XA -Z Ezx x 2 z 9 -ZM z E-ZZ - ZZA-úXA + x Xi-zza xa + xxa-zzú- X-Z + x M 9 EX xxa za- x-X E - zxs-xs-ZX&-z + Mxx-za-ZX 3-Z + ZX 3-xa 7 xa-z-A a-X -Z + ZXS'*XS- ?*z 3 *Z XX* Zz-ZX 9-*Z ZXS-*EL* Z + ZZS-ú* xg l -z Hd :(=/r H-) / l x XS-z XX -x-Z + ZZX-zx -Xa + xx S z 7ax xs X Zz ZX -zxs RZZ ZZX Sa-za ZX X - ZZS-z XJ + X Xa-z Zl x 3 + x X zz Z -X -Z + xza-úXg XXS*Z- -X -Z + ZX-* Tz X Z 7 Xx Zz-z X Z ZX 3 úXa Z + ZX&-xa-z Xa Z - zx J XS-zx EZ + M S E Za + XXS E Za Z Zaz Xi Z + Z& XS-z Z&Z az* za zxa + zza-zxa zú zxa-x -z Z Ez + xxa za zxa - zza- - x -z-Z zx-7 -zz -z + zx-x zz&-Z + XX 3 *-Zi-Z l =Ha TIV : (Z/E-Z) / l X-X - z XS + Z Za -X zx S - Xx*s S-z Xa-x zx&-zxa za-E + sxx zza xa zxa-zxa z& E + Xz S?X Z a-Z +Z LS *-Z RXZL+ *z Z - xx 2-Z Za-Xa + xxa zza xa z XXA zxa + zxa l X zx 3 ZX 3-xa-zxa z + xxa -úz + XXS-EZS zza-zxa zs + zx ezz - &.a-Z -zx& + zz- zx za ZX -X -ZZ -Z + Xxa -z -zx - aa-zs- a x î-zzza -x x- za Ez + zx + zxx r+zxx- zzaE - zx.-xa zz* z xx * zs * * z + zxa-& î*& E ZZ Hz zxx úz z l XV Hd UV Aw (za T a T':Ha Wl):yHaq (IIIA) Z-l HXHNNY ZZZ 689 Z 6 D ANNEXE 1-3 (XII) GRADK: (KIX, Kl Y, KIZ) Kl X a K 1/ax; K 2 X a K 2/ax; et GRADK 2:(K 2 X,K 2 Y,K 2 Z) avec K 1 Y a Kl/ay; Kl Z a Kl/az; K 2 Y a K 2/ay; K 2 Z a K 2/az; (XIII) E 1 = KIX-Tlx + Kl Y-Tly + KIZ-Tlz; (XIV) E 2 K 2 X-T 2 x + K 2 Y-T 2 y + Kl Z-Tlz; (s vi) * PIP 2 (V 2 Vl) (XVI) P A 1 + u A 1 A 2 (XVII) P Bl + v Bl B 2 avec: u v Al A 2 2 Al A 2 Bl B 2 A 2 Bl B 2 Bl B 2 2 \ AA 2 -Al BI B.B 2-Al Bli (k VIII) PI: (ul,vl) et 22: (u 2,v 2) soiutions de: Vll + (V 14-Vll) u + (V 12-Vll) v + (V 13 +Vll-V 14-V 12) uv = 51 V 21 + (V 24-V 21) u + (V 22-V 21) v + (V 23 +V 21-V 24-V 22) uv = 52 (XV) P Pl + ANNEXE 2 La structure DIFFER est constituée des 20 valeurs: F,FX, FY, FZ, FXX, FYY, FZZ, FXY, FXZ,FYZ,FXXX, FYYY, FZZZ, FXXY, FXYY, FXYZ,FYYZ, FYZZ,FXXZ et FXZZ. correspondant différentielles de F jusqu'à l'ordre 3, suivant les conventions suivantes: Fonction de départ: différentielles d'ordre 0: différentielles d'ordre 1: différentielles d'ordre 2: FXX= 82 F/ax 2; FXY= a 2 F/axay; F; FX a F/a 2 c; FYY= a 2 F/ay 2; FXZ= a 2 F/axaz; FY a F/ay; FZ a F/az; FZZ= a 2 F/az 2; FYZ= a 2 F/ayaz; différentielles d'ordre 3: FXXX= a 3 F/ax 3; FXXY= a 3 F/ax 2 ay; FYYZ= a 3 F/ay 2 az; FXZZ= a 3 F/axaz 2; FYYY= a 3 F/ay 3; FXYY= a 3 F/axay 2; FYZZ= a 3 F/ayaz 2; FZZZ= a 3 F/az 3; FXYZ= a 3 F/axayaz; FXXZ= a 3 F/ax 2 az; ANNEXE 3 La structure COURBURES peut contenir les valeurs GRAD,K,A,Kl,K 2,El,E 2, les vecteurs T 1 et T 2, ainsi que le drapeau ERREUR, avec: gradient de la fonction F Courbure gaussienne Courbure moyenne Courbures principales directions principales de courbures associées coefficients d'extrêmalité associés drapeau GRAD: K: A: K 1,K 2: Tl,T 2: El,E 2: ERREUR: REVENDICATIONS 1 Dispositif électronique pour le traitement d'informations d'image, caractérisé en ce qu'il comprend: des moyens de définition d'image (MI) propres à fournir, en correspondance d'un maillage polyédrique prédéterminé d'une portion de l'espace, une valeur numérique pour chaque sommet du maillage, des moyens de traitement de polygone ( 30; 1450; 1650), aptes à recevoir une représentation d'un polygone plan (F) à p sommets, p valeurs respectives pour ces p sommets, et un seuil (S), et à fournir, en correspondance du polygone F, une liste de segments orientés reliant des points des arêtes du polygone interpolés comme égaux au seuil, avec une convention de sens prédéterminée, des moyens de traitement de polyèdre ( 20; 1445; 1645) aptes à recevoir une représentation d'un polyèdre (POLYED) avec les valeurs en ses sommets et un seuil, puis à solliciter les moyens de traitement de polygone ( 30) pour toutes les faces du polyèdre, ce qui fournit une liste de segments orientés, puis à associer au polyèdre une liste de cycles fermés et orientés, constitués de segments de ladite liste de segments, et des moyens de traitement principaux ( 10; 1440; 1640) aptes à présenter séquentiellement une pluralité de polyèdres dudit maillage aux moyens de traitement de polyèdre, avec les valeurs associées à leurs sommets et un seuil, les cycles ainsi obtenus appartenant à une iso-surface de l'image pour ladite valeur de seuil. 2 Dispositif selon la revendication 1, caractérisé en ce que les moyens de traitement principaux ( 10; 1440; 1640) sont aptes à recevoir deux jeux de valeurs pour les points de ladite portion d'espace, ainsi que deux seuils ou iso-valeurs, et à déterminer des lignes caractéristiques par intersection des deux iso-surfaces correspondantes. 3 Dispositif selon la revendication 2, caractérisé en ce qu'il comprend en outre: des moyens d'intersection de cycles ( 60; 1470; 1670)) aptes à recevoir deux cycles du polyèdre, à rechercher, pour tout couple de segments constructible avec ces deux cycles, si ces segments sont dans la même face du polyèdre, auquel cas est recherchée une intersection de ces deux segments, et le point d'intersec- tion obtenu est retenu, pour constituer finalement une liste de segments orientés, et en ce que les moyens de traitement principaux ( 10; 1440; 1640) sont aptes à présenter séquentiellement deux pluralités de polyèdres dudit maillage aux moyens de traitement de polyèdre, avec les valeurs associées à leurs sommets et un seuil, puis les cycles obtenus aux moyens d'intersection de cycles, ce qui fournit des lignes remarquables par intersection de deux surfaces iso-valeur. 4 Dispositif selon la revendication 2, caractérisé en ce que les moyens de traitement principaux ( 10; 1440; 1640) sont aptes à déterminer des lignes remarquables en parcourant chaque cycle fourni, par lesdits moyens de traitement de polyèdre, pour une valeur du premier jeu en chaque point, à déterminer la valeur du second jeu en chaque sommet d'un cycle, à rechercher par interpolation, sur les segments du cycle, des points intermédiaires ayant pour valeur du second jeu le second seuil, et à constituer à partir de ces points des segments orientés, lesquels constituent une représentation approchée desdites lignes remarquables. 5 Dispositif selon l'une des revendications 1 à 4, caracté- risé en ce que les moyens de traitement de polygone compren- nent: un premier automate de base ( 40) apte à recevoir un couple de points, leurs deux valeurs respectives et un seuil, pour comparer les deux valeurs au seuil, et, si les sens de compa- raison sont différents, déterminer un point intermédiaire interpolé proche du seuil, en lui conférant un label lié au sens de la comparaison de l'un prédéterminé des deux points dans le couple, un second automate de base ( 70) propre à recevoir une suite de points formant un cycle, assortis chacun d'un label, à chercher le premier point possédant un label donné (-), et le prendre pour début, à établir un couple ou segment de ce point et du suivant, et à faire la liste de tels segments obtenus avec les autres points du cycle jusqu'à retour au point de début, et un premier automate intermédiaire ( 30), apte à recevoir p points associés à un polygone F, leur p valeurs respectives et un seuil, et, pour chaque arête du polygone, à mettre en oeuvre le premier automate de base ( 40), en mémorisant dans une zone dite "cycle" les points interpolés éventuels, à boucler ce "cycle" sur lui-même, puis à appeler le second automate de base ( 70) avec ce cycle bouclé, et à mémoriser la liste de segments orientés ainsi générés en correspondance du polygone F. 6 Dispositif selon la revendication 5, caractérisé en ce que le premier automate de base ( 40) comporte des moyens ( 432, 434, 442, 444) réagissant au fait que l'un des points dudit couple possède une valeur proche du seuil, en lui affectant une valeur écartée du seuil d'une manière prédéterminée. 7 Dispositif selon l'une des revendications 5 et 6, caracté- risé en ce que les moyens de traitement de polyèdre compren- nent: un premier automate supérieur ( 20) apte à recevoir un polyèdre (POLYED), avec les valeurs en ses sommets et le seuil, puis à appeler le premier automate intermédiaire ( 30) pour toutes les faces du polyèdre, ce qui fournit une liste de segments, puis, tant que cette liste est non vide, à en prendre le premier segment comme début d'un cycle, auquel est ajouté tout segment présentant une extrémité commune avec ledit premier segment, le cycle étant lui-même incorporé à une liste de cycles associée au polyèdre, et ce jusqu'à épuisement de ladite liste de segments. 8 Dispositif selon l'une des revendications 2 à 4, prises en combinaison avec la revendication 7, caractérisé en ce que les moyens d'intersection de cycles comprennent: un second automate intermédiaire ( 80) apte à recevoir deux segments ou couples de points d'une face polygonale F, et une normale orientée N, à établir par calcul de déterminant si les droites support des deux segments se coupent, puis déterminer si l'intersection fait partie des segments, auquel cas elle est retenue comme point P avec un label relié au signe du détermi- nant, et un second automate supérieur ( 60) apte à recevoir deux cycles du polyèdre, à rechercher, pour tout couple de segments constructible avec ces deux cycles, si ces segments sont dans la même face du polyèdre, auquel cas est appelé le second automate intermédiaire ( 80), et le point est ajouté dans un cycle, ledit cycle étant rebouclé, puis est appelé le second automate de base, pour donner une liste de segments orientés appartenant à la ligne d'intersection de deux surfaces iso-valeur. 9 Dispositif selon la revendication 8, prise en combinaison avec la revendication 6, caractérisé en ce que les moyens d'intersection de cycles comprennent en outre un troisième automate de base ( 90) réagissant au fait que l'un des points dans lesdits couples possède une valeur proche du seuil, en lui affectant ( 910-935) une valeur écartée du seuil d'une autre manière prédéterminée, et en reprenant ( 940-950) la détermination du point intermédiaire du couple, qui correspond au seuil concerné. 10 Dispositif selon l'une des revendications 1 à 9, caracté- risé en ce qu'il comprend des moyens (Al; Bl-B 3) pour commander l'exploration maille par maille des polyèdres de ladite portion d'espace par les moyens de traitement principaux. 11 Dispositif selon la revendication 10, caractérisé en ce qu'il comprend des moyens (B 1-B 3) pour commander l'exploration des polyèdres de ladite portion d'espace par les moyens de traitement principaux, en suivant une ligne d'intersection de deux iso-surfaces. 12 Dispositif selon l'une des revendications 1 à 11, dans lequel le maillage comprend une grille parallélipipédique régulière, caractérisé en ce que les moyens de traitement de polygone ( 1450) comparent directement, pour chaque face en parallélogramme, la configuration des points interpolés égaux au seuil à un jeu de configurations préétablies. 13 Dispositif selon la revendication 12, caractérisé en ce que l'interpolation entre deux sommets d'une face est réalisée en mode hyperbolique. 14 Dispositif selon la revendication 13, prise en combinaison avec la revendication 8, caractérisé en ce que les moyens d'intersection de cycles ( 1470) déterminent également les points d'intersection par interpolation hyperbolique. Dispositif selon l'une des revendications 1 à 11, dans lequel le maillage comprend une grille tétraédrique, caractérisé en ce que les moyens de traitement de polygone ( 1650) comparent directement, pour chaque face triangulaire d'un tétraèdre, la configuration des points interpolés égaux au seuil à un jeu de configurations préétablies. 16 Dispositif selon l'une des revendications 1 à 11, dans lequel le maillage comprend une grille tétraédrique, caractérisé en ce que les moyens de traitement de polyèdre ( 1645) comparent directement, pour chaque tétraèdre, la configuration des labels des sommets du tétraèdre à un jeu de configurations préétablies. 17 Dispositif selon l'une des revendications 1 à 16, prise en combinaison avec la revendication 2, caractérisé en ce que les deux jeux de valeurs appartiennent à deux images différentes. 18 Dispositif selon l'une des revendications 1 à 16, prise en combinaison avec la revendication 2, caractérisé en ce que l'un des jeux de valeurs appartient à une image, tandis que l'autre correspond à des valeurs différentielles relatives à la même image, ce qui fournit une ou plusieurs lignes remarquables relatives à ladite image. 19 Dispositif selon la revendication 18, caractérisé en ce que les valeurs différentielles appartiennent au groupe comprenant, pour un point échantillonné ou interpolé de l'image: le vecteur gradient en ce point, la courbure gaussienne en ce point, la courbure moyenne en ce point. 20 Dispositif selon l'une des revendications 18 et 19, caractérisé en ce que les valeurs différentielles appartiennent au groupe comprenant, pour un point échantillonné ou interpolé de l'image: les courbures principales en ce point, les directions principales de courbure en ce point. 21 Dispositif selon l'une des revendications 18 à 20, caractérisé en ce que les valeurs différentielles comprennent les coefficients d'extrémalité en un point échantillonné ou interpolé de l'image. 22 Dispositif selon l'une des revendications 18 à 21, caractérisé en ce qu'il comprend: des moyens (E 100, ML 1) pour extraire des lignes remarquables d'une première image, des moyens (ML 2) pour définir des lignes remarquables d'une seconde image, des moyens (E 300, E 400) pour déterminer, à partir desdites lignes remarquables, une transformation géométrique (M) faisant passer de l'une des images à l'autre, des moyens (E 500) pour rééchantillonner avec interpolation (MI 3) l'une des images, recalée en les mêmes points que l'autre, en vertu de ladite transformation géométrique, et des moyens (E 900) pour comparer l'image ainsi rééchantillonnée à ladite autre image. 23 Dispositif selon la revendication 22, caractérisé en ce que lesdits moyens (ML 2) utilisent une image de référence des lignes remarquables. 24 Dispositif selon la revendication 22, caractérisé en ce que lesdits moyens (ML 2) utilisent une seconde image (MI 2) relative au même objet que la première, et comprennent des moyens (E 200) pour extraire des lignes remarquables de cette seconde image.