Morphologie mathématique

Une forme (en bleu), sa dilatation morphologique (en vert), et son érosion morphologique (en jaune) par un élément structurant en forme de diamant.

La morphologie mathématique est une théorie et technique mathématique et informatique d'analyse de structures qui est liée avec l'algèbre, la théorie des treillis, la topologie et les probabilités.

Le développement de la morphologie mathématique est inspiré des problèmes de traitement d'images, domaine qui constitue son principal champ d'application. Elle fournit en particulier des outils de filtrage, segmentation, quantification et modélisation d'images. Elle est également utilisable en traitement du signal, par exemple pour filtrer les variations d'une mesure (physique, biologique) au cours du temps.

Aperçu général

Une des idées de base de la morphologie mathématique est d'étudier ou de traiter un ensemble à l'aide d'un autre ensemble, appelé élément structurant, qui sert de sonde. À chaque position de l'élément structurant, on regarde s'il touche ou s'il est inclus dans l'ensemble initial. En fonction de la réponse, on construit un ensemble de sortie. On obtient ainsi des opérateurs de base qui sont relativement intuitifs.

Des propriétés que l'on retrouve souvent dans les opérateurs morphologiques sont :

• la non-linéarité ;
• la non-inversibilité ;
• l'idempotence.

Ceci implique en particulier une perte d'information ; bien utilisés, ces opérateurs permettent d'éliminer des structures ne respectant pas certains critères, comme de largeur ou de volume.

La morphologie mathématique s'intéresse aussi aux ensembles et aux fonctions aléatoires.

Le principal domaine d'application de la morphologie mathématique est le traitement d'images. Elle fournit, en particulier, des outils de filtrage, de segmentation et de quantification. Depuis son apparition, en 1964, elle connaît un succès grandissant et désormais contribue à garnir la boite à outils de tout traiteur d'images.

Bref historique

La morphologie mathématique a été inventée en 1964 par Georges Matheron et Jean Serra dans les laboratoires de l'école des Mines de Paris. Son développement a toujours été fortement motivé par des applications industrielles. Dans un premier temps, il s'est agi de répondre à des problèmes dans le domaine de l'exploitation minière, mais très vite ses champs d'applications se sont diversifiés : biologie, imagerie médicale, sciences des matériaux, vision industrielle, multimédia, télédétection et géophysique constituent quelques exemples de domaines dans lesquels la morphologie mathématique a apporté une contribution importante.

La morphologie mathématique reste un domaine actif de recherche. En témoignent les nombreuses publications scientifiques sur le sujet, ainsi que les symposiums internationaux sur la morphologie mathématique qui ont lieu tous les deux ou trois ans.

Quelques exemples de thèmes de recherche actuels :

• ligne de partage des eaux : parallélisation, approche topologique, hiérarchisation ;
• extension à des espaces non-euclidiens, à des fonctions vectorielles (images en couleurs, images multi-spectrales, etc.) et à des images non conventionnelles;
• liens avec l'apprentissage profond (réseaux morphologiques; pré-traitement et post-traitement);
• approches stochastiques pour la segmentation et le filtrage;
• architectures électroniques efficaces pour l'implantation d'opérateurs morphologiques;
• développement d'algorithmes rapides;
• traitement de nuages de points;
• modélisation de structures physiques, en particulier dans le domaine des sciences des matériaux.

Opérateurs de base

La morphologie mathématique peut être développée dans le cadre abstrait de la théorie des treillis. Cependant, une présentation plus pratique, visant un utilisateur potentiel d'outils de traitement d'images, plutôt qu'un mathématicien, est ici adoptée.

Cas ensembliste

Plaçons-nous dans $E=\mathbb {Z} ^{2}$, souvent utilisé comme modélisation du support des images binaires à deux dimensions, même si tout ce qui est présenté dans cette section reste valable dans $\mathbb {R} ^{d}$, où $d$ est un entier strictement positif. Soit $B$ un sous-ensemble de $E$, appelé élément structurant. Si $x$ est un élément de $E$, alors nous noterons $B_{x}$ l'ensemble $B$ translaté de $x$ : $B_{x}=\{b+x\mid b\in B\}$

L'élément structurant joue en quelque sorte le rôle de modèle local, ou de sonde. Il est promené partout sur l'image à traiter, et à chaque position on étudie sa relation avec l'image binaire, considérée comme un ensemble. Ces relations peuvent être du type « est inclus dans l'ensemble », ou « touche l'ensemble », par exemple.

En pavage carré, les éléments structurants les plus classiquement utilisés sont la croix, constituée de l'origine et des quatre points les plus proches, et le carré, constitué de l'origine et des huit points les plus proches. Ces deux éléments structurants correspondent respectivement à deux définitions possibles du voisinage d'un pixel ou du type de connexité de l'image. En pavage hexagonal, l'élément de base est l'hexagone centré.

On introduit aussi le symétrique d'un ensemble, noté ${\breve {B}}$ : ${\breve {B}}=\{-b\mid b\in B\}$

Si $B$ est symétrique, on a ${\breve {B}}=B$.

Dilatation et érosion

Soit $X$ un sous-ensemble de $E$. La dilatation morphologique avec l'élément structurant $B$ est définie comme la somme de Minkowski: $\delta _{B}(X)=X\oplus B=\{x+b\mid b\in B,x\in X\}=\cup _{x\in X}B_{x}$

Une autre formulation plus intuitive est : $\delta _{B}(X)=\{x\mid {\breve {B}}_{x}\cap {X}\neq \emptyset \}$

La dilatation morphologique n'est, en général, pas inversible. L'opération qui en quelque sorte tente de produire l'inverse de la dilatation est l'érosion morphologique: $\epsilon _{B}(X)=X\ominus B=\{x\mid B_{x}\subset X\}$

La dilatation et l'érosion sont les opérateurs de base de la morphologie mathématique. Pratiquement tous les autres peuvent être définis à l'aide de ceux-ci, en utilisant des compositions de fonctions et des opérations ensemblistes.

• 
  Image originale (en noir : l'objet ; en blanc : le fond).
• 
  Dilatation par un carré 3x3 : les pixels noirs et gris font partie de l'ensemble résultant.
• 
  Érosion par un carré 3x3 : seuls les pixels noirs font partie de l'ensemble résultant.

Propriétés algébriques de la dilatation et de l’érosion

La dilatation est une transformation extensive si B contient l’origine : $X\subseteq \delta _{B}(X)$L’érosion est anti-extensive si B contient l'origine : $\epsilon _{B}(X)\subseteq X$La dilatation et l’érosion sont des transformations croissantes (comme l’union et l’intersection) : $X\subseteq Y\Rightarrow \delta _{B}(X)\subseteq \delta _{B}(Y)$$X\subseteq Y\Rightarrow \epsilon _{B}(X)\subseteq \epsilon _{B}(Y)$La dilatation et l’érosion ne sont pas des transformations idempotentes : $\delta _{B}(\delta _{B}(X))\neq \delta _{B}(X)$$\epsilon _{B}(\epsilon _{B}(X))\neq \epsilon _{B}(X)$Par contre la dilatation et l’érosion vérifient la propriété d’itérativité qui permet de construire des dilatations ou des érosions avec des éléments structurants homothétiques : $\delta _{3B}(X)=\delta _{1B}(\delta _{1B}(\delta _{1B}(X)))$$\epsilon _{3B}(X)=\epsilon _{1B}(\epsilon _{1B}(\epsilon _{1B}(X)))$La dilatation est une transformation continue et l’érosion est une transformation semi-continue supérieure. Cette propriété découle directement de la propriété de l’intersection dans la soustraction de Minkowski.

La dilatation est distributive par rapport à l’union et l’érosion par rapport à l’intersection : $\delta _{B\cup B'}(X)=\delta _{B}(X)\cup \delta _{B'}(X)$$\epsilon _{B\cap B'}(X)=\epsilon _{B}(X)\cap \delta _{B'}(X)$Soit une famille d’éléments structurants homothétiques et une homothétie de rapport $\lambda$. La compatibilité de la dilatation et de l’érosion avec les homothéties s’écrivent : $\delta _{1B}(X)={\frac {1}{\lambda }}\delta _{\lambda B}(\lambda X)$$\epsilon _{1B}(X)={\frac {1}{\lambda }}\epsilon _{\lambda B}(\lambda X)$

Propriétés topologiques de la dilatation et de l'érosion

La dilatation, comme l’union, préserve la connexité. La dilatation n’est pas une transformation qui préserve l’homotopie. En effet, elle connecte des éléments disjoints et bouche des trous. L’érosion n’est pas une transformation qui préserve l’homotopie. En effet, elle sépare les parties connexes et supprime des éléments lorsqu’ils sont petits. L’érosion, comme l’intersection, ne préserve pas la connexité.

Ouverture et fermeture

La composition d'une dilatation morphologique avec l'érosion par le même élément structurant ne produit pas, en général, l'identité, mais deux autres opérateurs morphologiques, l'ouverture morphologique : $\gamma _{B}(X)=X\circ B=\delta _{B}(\epsilon _{B}(X))$ et la fermeture morphologique : $\phi _{B}(X)=X\bullet B=\epsilon _{B}(\delta _{B}(X))$L'ouverture peut être caractérisée géométriquement: elle donne l'union de tous les $B_{x}$ inclus dans $X$. Ainsi, la forme de l'élément structurant permet de choisir les structures qui peuvent le contenir.

La fermeture est le dual de l'ouverture : la fermeture du complémentaire d'un ensemble est égale au complémentaire de l'ouverture de cet ensemble.

On notera que si l'élément structurant $B$ n'est pas symétrique, on devra utiliser l'élément symétrique ${\breve {B}}$ pour le second opérateur (dilatation dans le cas de l'ouverture et érosion dans le cas de la fermeture).

Propriétés de l’ouverture et de la fermeture

• La fermeture et l'ouverture sont des opérations croissantes et idempotentes, deux propriétés qui définissent les filtres morphologiques.
• La fermeture est extensive ($X\subset \phi (X)$), et l'ouverture est anti-extensive ($\gamma (X)\subset X$).
• L’ouverture et la fermeture sont des opérateurs semi-continu supérieur. C’est l’érosion qui impose cette propriété, la dilatation étant continue.
• Du point de vue topologique, l’ouverture ne préserve pas la connexité alors que la fermeture la préserve. L’ouverture et la fermeture ne sont pas des transformations homotopiques.

• 
  Image originale
• 
  Fermeture par un carré 3×3 : les pixels noirs et gris font partie de l'ensemble résultant
• 
  Ouverture par un carré 3×3 : seuls les pixels noirs font partie de l'ensemble résultant

Transformation en tout ou rien

On peut aussi prendre deux éléments structurants $A$ et $B$ pour définir des transformations. Si on demande en chaque point $x$ à $A$ d'être à l'extérieur de l'ensemble et à $B$ à l'intérieur on obtient la transformation en tout ou rien (hit or miss transform en anglais) : $\eta (X)=\{x\mid A_{x}\subset X^{c};B_{x}\subset X\}$ où $X^{c}$ désigne le complémentaire de l'ensemble $X$. Cette transformation permet de détecter certaines configurations précises de pixels. Parmi  les configurations les plus utilisées on aura :

• Les points extrêmes et les points triples du squelette obtenu par amincissement.
• Les configurations correspondant aux différents nombres de connexités.

Épaississement et amincissement

En ajoutant le résultat de la transformation à l'ensemble initial on obtient un épaississement : $\operatorname {ep} (X)=X\cup \eta (X)$en enlevant le résultat de l'ensemble initial on obtient un amincissement : $\operatorname {aminc} (X)=X-\eta (X)$

Applications : squelette, skiz, enveloppe convexe

• En prenant des suites d'amincissements, on peut réduire progressivement l'ensemble initial (comme si on l'épluchait). De cette façon on peut calculer différents types de squelettes, dont des squelettes homotopiques.
• Des épaississements homotopiques permettent aussi de construire des Skiz (squelettes par zone d’influence).
• Avec une autre configuration on pourra obtenir par épaississement les enveloppes convexes.

Extension aux fonctions

Une image à niveaux de gris peut être modélisée comme une fonction de $\mathbb {Z} ^{2}$ dans $\mathbb {Z}$. Soit $f$ une fonction appartenant à cet ensemble. On a alors :

$\delta _{B}(f)=\sup\{f_{b}\mid b\in B\}$

$\epsilon _{B}(f)=\inf\{f_{b}\mid b\in B\}$

L'ouverture et la fermeture de fonctions s'obtiennent comme dans le cas ensembliste :

$\gamma _{B}(f)=\delta _{B}\,\epsilon _{B}(f)$

$\phi _{B}(f)=\epsilon _{B}\,\delta _{B}(f)$

L'ouverture et la fermeture morphologiques constituent déjà des outils intéressants de filtrage d'images. Cependant, ils peuvent modifier le contour des objets, propriété qui peut être malvenue. Les opérateurs par reconstruction et plus généralement les nivellements, introduits plus loin, permettent de pallier cet inconvénient.

Épaississements et amincissements ne sont pas, en général, des opérateurs croissants. Par conséquent, leur application aux fonctions (en pratique, aux images à niveaux de gris) n'est pas triviale. Plusieurs extensions ont été proposées dans la littérature.

Exemple d'utilisation : détection de contours

La détection de contours représente une tâche importante en traitement d'images. La morphologie mathématique propose des outils non-linéaires de détection de contours, comme le gradient et le laplacien morphologiques.

Le gradient morphologique, aussi appelé gradient de Beucher du nom de son inventeur, est défini par :

$\operatorname {grad} _{B}(X)=\delta _{B}(X)-\epsilon _{B}(X)$

Il correspond, en quelque sorte, à la version morphologique du module du gradient euclidien.

Le laplacien morphologique est construit de façon analogue :

$\operatorname {lap} _{B}(X)=\delta _{B}(X)+\epsilon _{B}(X)-2\operatorname {id} (X)$

où $\operatorname {id}$ correspond à l'opérateur identité.

Opérateurs connexes, nivellements

Tous les opérateurs définis dans les sections précédentes l’ont été dans un cadre euclidien à savoir que l’espace de définition de l’image sert de référence pour les opérateurs. Dans cette section, on va reprendre les deux opérateurs de base que sont l’érosion et la dilatation mais en restant dans un sous espace $X$ de $\mathbb {E}$, nommé espace de référence. Les transformations euclidiennes deviendront ainsi des transformations géodésiques (appelées aussi transformations conditionnelles).

La géodésie est la science de la mesure de la forme et des dimensions de la terre. Ainsi une distance géodésique correspond au chemin le plus court pour aller d’un point à un autre en restant à la surface du globe. La longueur de ce chemin, à la différence de la distance euclidienne, ne correspond pas à un segment droit mais à celle d’un arc géodésique. La définition d’arc géodésique implique la notion de connexité par arc. Un espace topologique $\mathbb {E}$ est dit connexe par arcs si tout couple de points de $\mathbb {E}$ est relié par un chemin dont le support est inclus dans $\mathbb {E}$.

De l’euclidien au géodésique pour les ensembles

Distance euclidienne et distance géodésique

La distance géodésique utilise les mêmes axiomes que la distance euclidienne, seul, le chemin est différent.

Soit quatre points $\left\{x_{1},x_{2},x_{3},x_{4}\right\}\in X$.

Sur la partie gauche de l’image, on a représenté les segments de droite associés aux différentes distances euclidiennes. Sur la partie droite, on a représenté les trois arcs géodésiques reliant $x_{1},x_{2}$ et $x_{3}$. On remarque que le point $x_{4}$ ne possède aucun chemin géodésique car il appartient à une composante disjointe de celle contenant les trois autres points.

Comparaison entre distance euclidienne et géodésique

La distance géodésique vérifie bien les axiomes d’une distance. On a, en effet : $d_{X}(x_{1},x_{2})=d_{X}(x_{2},x_{1})$$x_{1}=x_{2}\Rightarrow d_{X}(x_{1},x_{2})=0$$d_{X}(x_{1},x_{3})\leq d_{X}(x_{1},x_{2})+d_{X}(x_{2},x_{3})$

À ces axiomes, il faut ajouter un quatrième. Lorsqu'il n’y a aucun chemin géodésique, on a :$d_{X}(x_{1},x_{4})=\infty$Pour les mêmes points, il est possible de comparer les distances géodésiques par rapport à $X$ et les distances euclidiennes. On aura toujours une distance géodésique plus grande ou égale à la distance euclidienne avec la possibilité d’avoir une distance euclidienne finie et une distance géodésique infinie :$d_{X}(x_{1},x_{2})\geq d(x_{1},x_{2})$On notera au passage que, pour relier deux points, il peut y avoir plusieurs arcs géodésiques équivalents, alors que le chemin euclidien est unique.

Élément structurant géodésique

Dans $\mathbb {R} ^{2}$ l’élément structurant isotrope centré en $x$ et de taille $\lambda$ est le disque fermé $\lambda B$, défini par :$\lambda B_{x}=\left\{y:d(x,y)\leq \lambda \right\}$Le disque géodésique, $\lambda B_{X}$ se définit par rapport à l’ensemble de référence $X$ en remplaçant la distance euclidienne par la distance géodésique par rapport à $X$. On a alors :$(\lambda B_{X})_{x}=\left\{y:d_{X}(x,y)\leq \lambda \right\}$

Boule géodésique (cyan) et euclidienne (disque entier)

La figure ci-contre illustre la différence qui existe entre un élément structurant euclidien et un élément structurant géodésique. Sur cette figure, on a placé le disque $B$ sur différentes positions $x_{1},x_{2},x_{3},x_{4},x_{5}$. Le disque euclidien reste le même indépendamment de la position x. Par contre, le disque géodésique change de forme ou peut disparaître en fonction de sa position. Du disque initial, il ne reste que ce qui vérifie la distance géodésique. Ainsi les points du disque qui ne sont pas commun avec $X$ sont ignorés (positions ${\textstyle x_{3},x_{4}}$). Si le centre du disque n’appartient pas à $X$, alors le disque géodésique n’existe pas ($x_{5}$).

Les transformations géodésiques de base pour les ensembles

En plus de l’ensemble à transformer que l’on appellera $M$ (pour marqueur), il faut introduire l’ensemble de référence géodésique $X$.

La dilatation géodésique

On part de la définition euclidienne de la dilatation en remplaçant la boule euclidienne par la boule géodésique. On peut alors écrire : $\delta _{X,B}(M)=\left\{z:(B_{X})_{z}\Uparrow M\right\}=\left\{z:(B_{X})_{z}\cap M\neq \emptyset \right\}$Ensuite, on répète cette opération élémentaire $\lambda$ fois : $\delta _{X,nB}(M)=\delta _{X,1B}(...\delta _{X,1B}(\delta _{X,1B}(M))...)$ Les figures suivantes illustrent l’effet d’une dilatation géodésique par un hexagone. Une composante connexe de l’ensemble $M$ ne peut être dilaté que si elle intersecte $X$, sinon elle disparaît. La dilatation s’arrête lorsque les frontières de $X$ sont atteintes.

• 
  Ensemble de référence X (jaune et rouge) et marqueur M (rouge et bleu)
• 
  Dilatation géodésique (hexagone taille 15) de M dans X
• 
  Dilatation géodésique (hexagone taille 40) de M dans X

L’ensemble $X$ se comporte comme un masque dans lequel on peut modifier $M$ et $M$ se comporte comme un marqueur permettant d’envahir une composante connexe de $X$.

L’érosion géodésique

De la même manière, on introduit l’érosion géodésique en remplaçant la boule euclidienne par la boule géodésique dans l’expression définissant l’érosion. On peut donc écrire :$\epsilon _{X,B}(M)=\left\{z:(B_{X})_{z}\subseteq M\right\}$Comme l’élément structurant utilisé est symétrique, on a dualité entre l’érosion et la dilatation géodésique. Cette dualité s’exprime de manière légèrement différente de celle du cas euclidien car c’est le complémentaire par rapport à $X$ qu’il faut prendre en compte. Cette relation de dualité s’écrit alors avec l’opérateur de différence symétrique :$\epsilon _{X,B}(M)=X/\delta _{X,B}(X/M)$Les figures illustrent l’effet de l’érosion géodésique sur un ensemble $M$ par rapport à $X$ en utilisant un élément structurant hexagonal.

• 
  Ensemble de référence X (jaune et rouge) et marqueur M (rouge)
• 
  Érosion géodésique (hexagone taille 15) de M dans X

On remarque que les parties de $M$, totalement incluses dans $X$ et n’ayant pas de frontières communes, s’érodent comme dans le cas euclidien. Lorsqu’il y a des frontières communes, celles-ci ne sont pas affectées par l’érosion.

La reconstruction pour les ensembles

Définition de la reconstruction

La reconstruction d’un ensemble à partir d’un autre est une des principales applications de la dilatation géodésique. On part donc de deux ensembles ; le premier est appelé l’ensemble des marqueurs noté $M$, le second est l’ensemble de référence ou masque $X$. Par définition, la reconstruction est une dilatation géodésique infinie, par rapport à $X$, des marqueurs $M$. Cela s’écrit :$\gamma ^{Rec}(M:X)=\delta _{X,\infty B}(M)$Lorsque toutes les composantes connexes contenant des marqueurs seront envahies, on ne pourra plus modifier l’image. Ceci constitue le test d’arrêt de la procédure. Les figures suivantes illustrent la reconstruction de $X$ à partir des marqueurs $M$.

• 
  Ensemble de référence X (jaune et rouge) et marqueur M (rouge)
• 
  Reconstruction de X à partir des marqueurs M

Propriétés de la reconstruction

Propriétés par rapport à l’ensemble de référence

• La reconstruction est idempotente : la reconstruction s’arrêtant lorsque l’on ne modifie plus rien, elle vérifie cette propriété.

• La reconstruction est anti-extensive par rapport à $X$ : ne sont reconstruites que les composantes connexes contenant un marqueur.$\gamma ^{Rec}(M:X)\subseteq X$
• La reconstruction est croissante :$X_{1}\subset X_{2}\Rightarrow \gamma ^{Rec}(M:X_{1})\subset \gamma ^{Rec}(M:X_{2})$La reconstruction correspond donc à une ouverture algébrique de l’ensemble de référence.

Propriétés par rapport aux marqueurs

• La reconstruction est idempotente

• La reconstruction est extensive par rapport à $M$ :$\gamma ^{Rec}(M:X)\supseteq M$
• La reconstruction est croissante :$M_{1}\subset M_{2}\Rightarrow \gamma ^{Rec}(M_{1}:X)\subset \gamma ^{Rec}(M_{2}:X)$

La reconstruction correspond donc à une fermeture algébrique des marqueurs.

Applications de la reconstruction pour les ensembles

On ne citera ici que les plus importantes.

L’érosion reconstruction

Dans l’érosion-reconstruction, les marqueurs $M$ seront les érodés euclidiens de $X$ et la dilatation euclidienne sera remplacée par la dilatation géodésique des marqueurs par rapport à l’ensemble $X$. Cela s’écrit :$\gamma _{\lambda B}^{Rec}(X)=\gamma ^{Rec}(\epsilon _{\lambda B}(X):X)$Cette érosion-reconstruction est illustrée par les figures suivantes et comparée avec l’ouverture géodésique utilisant le même élément structurant.

• 
  Ensemble X1 (jaune et rouge) et érosion hexagonale de taille 11 (M1 = rouge)
• 
  Reconstruction de X1 à partir de M1
• 
  Ouverture géodésique de X1 (rouge : hexagone de taille 11)

Élimination des composantes connexes coupant le bord du champ d’image

Pour certaines applications analytiques, il est nécessaire d’éliminer les composantes connexes de $X$ qui coupent le bord du champ d’observation $Z$. Pour cela, il faut les isoler en les reconstruisant à partir d’un marqueur constitué de tous les pixels du bord du champ $Z$, noté $M=\partial Z$ et les soustraire de $X$. La procédure va donc être la suivante :$Y=X/\gamma ^{Rec}(\partial Z:X)$L’ensemble $Y$ ne contient que des composantes connexes totalement incluses dans $Z$. Ceci est illustré par les figures suivantes.

• 
  Ensemble X2 (jaune) et bord du champ dZ (rouge)
• 
  Composantes connexes de X2 reconstruites, coupant le bord du champ (rouge)
• 
  Composantes connexes de X2 totalement incluses dans le champ Z

Bouchage de trous

La reconstruction est utilisée dans une opération qui sert souvent en traitement d’image : le bouchage de trous. Comme dans l’application précédente, le marqueur est $M=\partial Z$. L’algorithme est le suivant :$Y=\left(\gamma ^{Rec}(\partial Z:X^{C})\right)^{C}$Les figures suivantes illustrent la séquence d’opérations.

• 
  Ensemble X3 (jaune) et bord du champ dZ (rouge)
• 
  Complémentaire de l'ensemble X3 (jaune) et bord du champ dZ (rouge)
• 
  Reconstruction du complémentaire de X3 à partir de dZ
• 
  Bouchage des trous de X3

Les érodés ultimes

• Définition

Soit un ensemble $X$ présenté sur la figure suivante. Les érodés ultimes apparaissent au cours d’une succession d’érosions de $X$ par un élément structurant $B$ convexe. Ils sont constitués par l’union des composantes connexes qui disparaissent lors de l’érosion de taille immédiatement supérieure.

Érodés ultimes (rouge) de l'ensemble X4

Soit $\epsilon _{1B}$ une érosion digitale élémentaire, et $\epsilon _{iB}$ son itérée d'ordre i. L'érosion ultime dérivée de $\epsilon$ et notée $U_{i}(X)$ se définit alors comme le ‘’résidu’’ entre $\epsilon _{iB}$ et les ouvertures par reconstruction de chaque érosion dans la précédente :$U_{i}(X)=\epsilon _{iB}(X)/\gamma ^{Rec}(\epsilon _{(i+1)B}(X):\epsilon _{iB(X)})$Les érodés ultimes correspondent à l’union de toutes ces reconstructions en faisant varier de 1 à imax lorsqu’il n’y a plus rien par érosion.$U(X)=\bigcup _{i}(U_{i}(X))$

• Propriétés des érodés ultimes
  • L’érosion ultime est anti extensive.
  • Les éléments constituants l’ensemble des érodés ultimes sont disjoints.
  • Une érosion unitaire sur les érodés ultime donne un ensemble vide.
  • L’érosion ultime est idempotente.

Ces érodés ultimes peuvent également être obtenus à partir des maxima de la fonction distance en utilisant une méthode géodésique pour les fonctions que l’on voit ci après. En restant dans le domaine ensembliste, les érodés ultimes servent à marquer les parties convexes des objets et peuvent être utilisés pour segmenter des agrégats de particules convexes.

Transformations géodésiques pour les fonctions

Comme pour le cas ensembliste, on aura une fonction de référence $f(x)$ et une fonction qui sera transformée $m(x)$. Les deux opérateurs de base sont toujours la dilatation géodésique et l’érosion géodésique.

Dilatation géodésique pour les fonctions

La dilatation géodésique élémentaire de $m(x)$ sous $f(x)$ s’exprime de manière analogue à celle utilisée pour les ensembles.

Dans le cas des fonctions, ces éléments seront des éléments plats et convexes isotropes. On a, en effet :$\delta _{f,1B}(m(x))=\delta _{1B}(m(x))\wedge f(x)$Pour une dilatation géodésique de taille quelconque, on aura également :$\delta _{f,\lambda B}(m(x))=\delta _{f,1B}(...\delta _{f,1B}(m(x))...)$

Interprétation sur une fonction f(x) définie dans R1-R

Les images en niveaux de gris sont des fonctions $f(x)$ définies dans $\mathbb {R} ^{2}\times \mathbb {R}$. Pour faciliter la vision du comportement des opérateurs, on va utiliser une fonction $f(x)$ définie dans $\mathbb {R} ^{1}\times \mathbb {R}$. Soit également une fonction $m(x)$ définie dans le même espace. Pour notre exemple, la seule possibilité est un ‘’élément structurant plat’’ constitué d’un ‘’segment droit’’ de longueur $\lambda$ centré.

Les figures suivantes illustrent le cas d’une dilatation géodésique de taille $\lambda =2n+1$.

• 
  Fonction f(x) (jaune et rouge) et marqueurs m1(x) (rouge)
• 
  Dilatation de taille 20 des marqueurs m1(x) (rouge) sous la fonction f(x) (jaune et rouge)

On remarque que les marqueurs dilatés restent toujours sous la fonction $f(x)$. Avec un élément structurant plat, il reste des parties de $f(x)$ qui sont inaccessibles à la dilatation. Ce sont des parties convexes incluant les maxima.

On remarque que les marqueurs sont généralement choisis de telle sorte que l’on ait :$m(x)<f(x)$

Érosion géodésique pour les fonctions

Par définition, l’érosion géodésique élémentaire est donnée par :$\epsilon _{f,1B}(m(x))=\epsilon _{1B}(m(x))\vee f(x)$Comme pour la dilatation géodésique, on aura par itération :$\epsilon _{f,\lambda B}(m(x))=\epsilon _{f,1B}(...\epsilon _{f,1B}(m(x))...)$L’érosion géodésique pour les fonctions se déduit également de la dilatation géodésique pour les fonctions par dualité. Appelons $M$ ‘’le niveau de gris maximum supporté par l’image’’. On aura alors :

$\epsilon _{f,B}(m(x))=M-\delta _{f,B}(M-m(x))$

Cette expression est d’ailleurs celle utilisée pour construire n’importe quelle érosion géodésique.

Interprétation sur une fonction définie dans R1-R

Les figures illustrent le comportement de l’érosion géodésique pour une fonction $\mathbb {R} ^{1}\times \mathbb {R}$.

• 
  Fonction f(x) (jaune) et marqueurs m2(x) (jaune et rouge)
• 
  Érosion de taille 15 des marqueurs m2(x) (jaune et rouge) sur la fonction f(x) (jaune)

On remarque que la fonction érodée conserve des résidus dans les parties concaves de la fonction $f(x)$. Ce résultat est symétrique de celui obtenu par dilatation géodésique.

La reconstruction pour les fonctions

Comme nous l’avons fait pour les ensembles, il est possible de reconstruire des fonctions géodésiquement par rapport à une autre fonction. Deux cas sont à envisager.

• La reconstruction de $m(x)sousf(x)$ par dilatation géodésique infinie (ouverture)

$\gamma ^{Rec}(m(x):f(x))=\delta _{f,\infty B}(m(x))$

• La reconstruction de $m(x)surfx)$ par érosion géodésique infinie (fermeture)

$\varphi ^{Rec}(m(x):f(x))=\epsilon _{f,\infty B}(m(x))$

La première figure illustre le résultat obtenu dans le cas de la reconstruction de $m1$ sous la fonction $f$ et la seconde figure, la reconstruction de $m2$ sur la fonction $f$.

• 
  Reconstruction des marqueurs m1(x) (rouge) sous la fonction f(x) (jaune et rouge)
• 
  Reconstruction des marqueurs m2(x) (jaune et rouge) sur la fonction f(x) (jaune)

Application de la reconstruction pour les fonctions

Maxima d’une fonction

Les maxima régionaux d’une image sont les points de l’image à partir desquels il n’existe que des chemins descendants. Soit $f(x)$ une image. À partir de cette image, on va construire une image des marqueurs $m(x)$ en soustrayant un niveau de gris à $f(x)$. on aura donc :$m(x)=f(x)-1$Ensuite, on effectue une reconstruction de $m(x)$ sous $f(x)$ et par différence avec $f(x)$, on obtient les maxima régionaux. On a donc :$Max(f(x))=f(x)-\gamma ^{rec}(f(x)-1:f(x))$

Maxima (rouge) d'une image en niveaux de gris

Minima d’une fonction

Pour la recherche des minima régionaux d’une fonction, on applique le même principe. Tout d’abord, on forme l’image des marqueurs :$m(x)=f(x)+1$Ensuite, on effectue une reconstruction de $m(x)$ sur $f(x)$ et par différence avec $f(x)$, on obtient les minima régionaux. On a donc :$Max(f(x))=f(x)-\varphi ^{rec}(f(x)+1:f(x))$

Minima (rouge) d'une image en niveaux de gris

Maxima et minima étendus

La recherche des maxima et minima d’une fonction donnent de très bons résultats si l’image n’est pas bruitée. En présence de bruit, la notion de maxima et de minima étendus, appelés aussi Hmax et Hmin, permet de n’extraire d’une image que les extrema significatifs. L’algorithme est semblable à celui des maxima et des minima. Seule, la construction des marqueurs est légèrement différente. En effet, au lieu de translater l’image d’un niveau de gris (en moins ou en plus), on effectue une translation de h niveaux de gris en plus ou en moins. Les Hmax et Hmin s’écrivent alors de la manière suivante.$Hmax(f(x))=f(x)-\gamma ^{rec}(f(x)-h:f(x))$$Hmin(f(x))=f(x)-\varphi ^{rec}(f(x)+h:f(x))$

Les figures suivantes illustrent la construction des Hmax dans le cas d’une fonction $\mathbb {R} ^{1}\times \mathbb {R}$.

• 
  Fonction f(x) (jaune et rouge) et marqueurs m(-h) (rouge)
• 
  Fonction (fx) (jaune et rouge) et Hmax de f(x) (rouge)

A titre d’exemple, on montre le résultat dans le cas de l’image d’un pavage en niveaux de gris bruité. Les maxima régionaux de ce pavage sont inutilisables à cause du bruit. Par contre, les Hmax permettent de visualiser les plateaux de chaque pavé.

• 
  Image (en niveaux de gris) d'un pavage bruité
• 
  Maxima du pavage bruité
• 
  Hmax (h=30) du pavage bruité
• 
  Filtrage morphologique des Hmax par une ouverture de taille 1 suivie d'une fermeture de taille 10.

La recherche des maxima significatifs peut être améliorée par un filtrage de l’image binaire résultante, comme le montre la figure suivante. Les ‘’Hmin’’ sont obtenus et traités de manière analogue.

Segmentation

Segmenter une image à niveaux de gris consiste à produire une partition du support de l'image, de façon que les régions de la partition correspondent aux objets présents dans l'image.

Les filtres morphologiques constituent une aide précieuse dans un processus de segmentation. En particulier, les nivellements permettent de filtrer les images tout en préservant les contours importants, ce qui simplifie l'opération de segmentation proprement dite. Dans certains cas, un filtrage important peut de lui-même produire une partition pertinente. Mais l'outil morphologique le plus connu en segmentation d'images est la ligne de partage des eaux.

Il existe plusieurs algorithmes de segmentation par ligne de partage des eaux. L'idée de base consiste à simuler une inondation de l'image, vue comme un relief topographique où le niveau de gris correspond à l'altitude. Les frontières entre régions de la partition ont alors tendance à se placer sur les lignes de crête. Typiquement, on applique cet opérateur au gradient de l'image (norme du gradient euclidien, ou gradient morphologique) que l'on cherche à segmenter, et par conséquent les frontières se placent de façon privilégiée sur les lignes de gradient élevé.

Plusieurs algorithmes de calcul de ligne de partage des eaux ont une complexité linéaire en fonction du nombre de pixels de l'image, ce qui les place parmi les méthodes de segmentation les plus rapides.

Quantification : Paramètres de base

À l’origine la morphologie mathématique a été conçue pour traiter et analyser des images de matériaux ou images biologiques afin d’en extraire des informations quantifiées sous forme de paramètres ou de fonctions:. Ici, on se limitera aux images 2D définies dans l’espace $E=\mathbb {Z} ^{2}$ et aux sous-espaces. Dans ce cas l’espace $\mathbb {Z} ^{2}$ est représenté par une grille de points. Deux cas sont envisagés: la grille carrée (pavage carré) et la grille triangulaire (pavage hexagonal). En ce qui concerne les paramètres, on sait qu’ils peuvent être obtenus à partir de la caractéristique d’Euler-Poincaré ou nombre de connexité des différents espaces, notés $N_{n}(X)$ pour l’espace $\mathbb {R} ^{n}$.

Les nombres de connexité dans l’espace discret

• 
  Pavage carré d'un ensemble X
• 
  Ensemble X et configuration de voisinage pour obtenir N1 (pavage carré)
• 
  Ensemble X et configurations de voisinage pour obtenir N2 (pavage carré)

Espace R⁰ : N₀ (X)

Cet espace correspond au réseau de points associés aux pixels.Sur l’image binaire, $N_{0}(X)$ est égal au nombre de pixel à 1.

Espace R¹ : N₁ (X)

Les droites $\Delta$ utilisables dans $\mathbb {Z} ^{2}$ correspondent à des pixels alignés. Les extrémités de segments de ces droites coupant $X$ correspondent (en sortie) à des transitions de pixel de type 1 0. L’image binaire associée $Y$ est obtenue par une transformation en tout ou rien. D'un point de vue pratique cela revient à vérifier, pour chaque pixel $x$, la configuration de voisinage $T(x)$. Les éléments 1 de la configuration sont relatifs à l’ensemble $X$ et ceux à 0 au complémentaire. On aura donc :Pour l'ensemble : $Y=\eta _{T}(X)$Pour la mesure : $N_{1}(X)=N_{0}(Y)$

Les éléments structurants dans les différents pavages sont :

• En pavage hexagonal dans la direction 0 : $T(x)={\begin{bmatrix}&.&&.&\\.&&1&&0\\&.&&.&\end{bmatrix}}$.

(Les autres orientations à 60° et 120° sont obtenues par rotation de la configuration.)

• En pavage carré dans la direction 0 : $T(x)={\begin{bmatrix}.&&.&&.\\.&&1&&0\\.&&.&&.\end{bmatrix}}$.

(Les autres orientations à 45°, 90°, 135° sont obtenues par rotation de la configuration.)

Espace R² : N₂ (X)

Rappelons que $N_{2}(X)$ correspond au nombre de composantes connexes diminué du nombre de trous qu’elles contiennent.

• Pavage hexagonal

Pour déterminer ce nombre avec le pavage triangulaire, on utilise la relation d’Euler : $\chi =N_{2}(X)=s-c-f$En pavage hexagonal s représentant le nombre de sommets (pixels à 1), c le nombre de côtés de type 1-1 (à une rotation près) et f le nombre de triangles ayant les 3 sommets à 1. Un calcul élémentaire sur toutes les combinaisons donne le résultat suivant :Pour les ensembles : $Y1=\eta _{T1}(X)$ et $Y2=\eta _{T2}(X)$Pour la mesure : $N_{2}(X)=N_{0}(Y1)-N_{0}(Y2)$

Les éléments structurants dans les différents pavages sont :

$T1(x)={\begin{bmatrix}&.&&0&\\.&&1&&0\\&.&&.&\end{bmatrix}}$ et $T2(x)={\begin{bmatrix}&1&&0&\\.&&1&&.\\&.&&.&\end{bmatrix}}$.

• En pavage carré (8 connexité) on fait le même raisonnement ce qui donne les éléments structurants :

$T1(x)={\begin{bmatrix}.&&.&&.\\.&&1&&0\\.&&0&&0\end{bmatrix}}$ et $T2(x)={\begin{bmatrix}.&&.&&.\\.&&.&&1\\.&&1&&0\end{bmatrix}}$.

Les paramètres métriques de base associés

Comme pour les nombres de connexité, les paramètres métriques de base doivent vérifier les conditions de Hugo Hadwiger. L’ensemble ${X}$ doit être un ensemble aléatoire stationnaire et constitué d’une union finie de convexes. La mesure $W$ doit avoir les propriétés suivantes :

• Invariante par translation de ${\vec {t}}$ :${W(X)=W({X_{\vec {t}})}}$
• Compatible avec homothéties (agrandissement $\lambda$) :${W(X)={\frac {W({\lambda }X)}{\lambda ^{k}}}}$
• C-additive :${W(X)+W(Y)=W(X\cup Y)+W(X\cap Y)}$
• Continue ou semi continue

Dans R¹

Le paramètre métrique est la longueur totale de l’ensemble $X$ notée $L_{1}(X)$.Elle se calcule à partir de $N_{0}(X)$ et la taille du pixel $\mathrm {d} D$. On a en effet :${L_{1}}(X)=N_{0}(X)\times \mathrm {d} D$

Dans R²

Ces paramètres métriques sont:

• L’aire de l’ensemble $X$ notée $A(X)$

Elle se calcule à partir de $N_{0}(X)$ et l’aire du pixel $\mathrm {d} A$. On a en effet : ${A}(X)=N_{0}(X)\times \mathrm {d} A$

• Le périmètre de l’ensemble $X$ noté $L_{2}(X)$

Pour obtenir ce périmètre, on va utiliser la relation de Cauchy (géométrie intégrale) qui relie la variation diamétrale d'un ensemble à son périmètre :$L_{2}(X)=\int _{0}^{\pi }D_{2}(X,\alpha )\times \mathrm {d} \alpha =\pi \times \mathbb {E} [{D(X)}]=\pi \times \mathbb {E} [{N_{1}(X)}]\times \mathrm {d} D$ avec $\mathrm {d} D$ taille du pixel. On notera que l'estimation de ce périmètre à un aspect statistique. Le nombre de connexité $N_{1}(X,\alpha )$ doit être estimé dans plusieurs directions.

• 
  Illustration de la relation de Cauchy
• 
  Illustration de la relation de Crofton
• 
  Illustration de la relation de Meunier

Vers l'espace R³

La géométrie intégrale permet également d’accéder à des paramètres de $\mathbb {R} ^{3}$ en utilisant les nombres de connexité des espaces inférieures.

• Ainsi, la surface de la frontière $\partial X$ notée $S(X)$ est obtenue à partir du nombre de connexité $N_{1}(X)$ par la relation de Crofton :$S(X)={\frac {1}{\pi }}\int _{0}^{4\pi }A(X',\omega )\times \mathrm {d} \omega =4\times \mathbb {E} [A(X')]=4\times \mathbb {E} [{N_{1}(X)}]\times \mathrm {d} A$$X'$ est la projection de $X$ sur le plan $\Pi$ perpendiculaire à la direction $\omega$ et $\mathrm {d} A$ l'aire associée à la ligne $\Delta$.

• L'intégrale de courbure moyenne est estimée à partir de la variation diamétrale dans ${\mathbb {R} ^{3}}$ notée $D_{3}(X,\omega )$en utilisant la relation de Meusnier:$M(X)=2\pi \times {\mathbb {E} [D_{3}(X,\omega )]}$Cette variation diamétrale est reliée au nombre de connexité dans ${\mathbb {R} ^{2}}$ par la relation :$D_{3}(X,\omega )=N_{2}(X,\omega )\times \mathrm {d} D_{3}$Ce qui donne finalement :$M(X)=2\pi \times {\mathbb {E} [N_{2}(X,\omega )]\times \mathrm {d} D_{3}}$$\mathrm {d} D_{3}$ représente la distance entre 2 plans.

De l’analyse globale à l’analyse locale (stéréologie )

Les images destinées aux études scientifiques sont souvent obtenues à partir d’un microscope dont le champ est plus petit que l’ensemble à analyser. Dans ce cas, on dit que l’analyse est locale par opposition à l’analyse globale où l’ensemble est totalement visible.

Les paramètres globaux précédemment définis doivent être transformés en paramètres locaux ramenés à l’unité d’espace.

Paramètres locaux de l'espace R⁰

• L’unique paramètre est la fraction de points :${N_{P}(X)={\frac {N_{0}(X\cap Z_{0})}{N_{0}(Z_{0})}}}$

Paramètres locaux de l'espace R¹

• Fraction linéique :${L_{L}(X)={\frac {L_{1}(X\cap Z_{1})}{L(Z_{1})}}}$Dans ce contexte statistique, on démontre aisément que l'estimation de la fraction linéique est égale à la fraction de points :${L_{L}(X)=N_{P}(X)}$

• Nombre de connexité par unité de longueur:${N_{L}(X)={\frac {N_{1}(X\cap Z_{1})}{L(Z_{1})}}}$

Paramètres locaux de l'espace R²

• Fraction surfacique :${A_{A}(X)={\frac {A(X\cap Z_{2})}{A(Z_{2})}}}$On a également :${A_{A}(X)=N_{P}(X)}$

• Périmètre spécifique :${L_{A}(X)={\frac {L_{2}(X\cap Z_{2})}{A(Z_{2})}}}$En utilisant la relation de Cauchy, on obtient :${L_{A}(X)=\pi \times N_{L}(X)}$

• Nombre de connexité par unité de surface :${N_{A}(X)={\frac {N_{2}(X\cap Z_{2})}{A(Z_{2})}}}$

Paramètres locaux de l'espace R³

• Fraction volumique :${V_{V}(X)={\frac {V(X\cap Z_{3})}{V(Z_{3})}}}$Ici encore, on a :${V_{V}(X)=N_{P}(X)}$

• Surface spécifique :${S_{V}(X)={\frac {S(X\cap Z_{3})}{V(Z_{3})}}}$En utilisant la relation de Crofton, on a :${S_{V}(X)=4\times N_{L}(X)}$

• Intégrale de courbure moyenne par unité de volume :${M_{V}(X)={\frac {M_{3}(X\cap Z_{3})}{V(Z_{3})}}}$L'intégrale de courbure moyenne peut être estimée à partir du nombre de connexité dans ${\mathbb {R} ^{2}}$ selon la relation de Meusnier :${M_{V}(X)=2\pi \times N_{A}(X)}$

• Nombre de connexité par unité de volume :${N_{V}(X)={\frac {N_{3}(X\cap Z_{3})}{V(Z_{3})}}}$Ce paramètre de nature topologique n'est pas accessible à partir des espaces de dimension inférieure.

Quantification 2 : Granulométrie et dispersion

Les paramètres stéréologiques sont des paramètres moyens. De plus, ils ne sont pas nombreux. On conçoit aisément qu’ils soient insuffisants pour donner une description assez complète de la structure. Si l’on accepte de perdre l’aspect stéréologique, la morphologie mathématique permet d’obtenir de nombreuses informations quantitatives supplémentaires. Cette quantification dépend souvent d’un paramètre de taille $\lambda$ associé aux transformations d’images. La quantification conduira à une opération de tri $\Psi _{\lambda }(X)$ dont le dénombrement ou la mesure conduira à une fonction granulométrique. La dispersion d’un ensemble dans un autre est également quelque chose d’important à connaître. La stéréologie ne fournit qu’un paramètre dérivé qui est loin de répondre à la question.

L’analyse granulométrique

Axiomes de l’analyse granulométrique (Georges Matheron)

La méthode de tri doit vérifier les règles suivantes :

• La méthode de tri doit vérifier la propriété de croissance$X\subset Y\Rightarrow \Psi _{\lambda }(X)\subset \Psi _{\lambda }(Y)$
• La méthode de tri est anti-extensive$\Psi _{\lambda }(X)\subset X$

• Dans le cas de deux tris successifs, c’est le paramètre le plus grand qui l’emporte. Le résultat est indépendant de l’ordre des opérateurs. Cela s’écrit :$\Psi _{\lambda 2}(\Psi _{\lambda 1}(X))=\Psi _{\lambda 1}(\Psi _{\lambda 2}(X))=\Psi _{sup(\lambda 2,\lambda 1)}(X)$ Cette dernière condition va réduire considérablement le nombre des opérateurs candidats.

Les fonctions granulométriques

On distingue les granulométries en nombre et les granulométries en mesure.

• Si le tri granulométrique $\Psi _{\lambda }(X)$ consiste à ne conserver que tout ce qui est plus petit que $\lambda$, alors la fonction de distribution granulométrique en nombre sera définie par :$F(X,\lambda )={\frac {Nb(X)-Nb(\Psi _{\lambda }(X))}{Nb(X)}}$

• Au lieu de compter les objets, on peut reporter une mesure qui leur est associée (la masse dans un tamis).On obtient alors une fonction de distribution granulométrique en mesure définie par :$G(X,\lambda )={\frac {Mes(X)-Mes(\Psi _{\lambda }(X))}{Mes(X)}}$ En analyse d’image la mesure est le volume dans $\mathbb {R} ^{3}$, l’aire pour $\mathbb {R} ^{2}$ et la longueur pour $\mathbb {R} ^{1}$.

Analyse granulométrique individuelle

Ensemble de particules individuelles (jaune)

Ce type d’analyse n’est possible que si l’ensemble à analyser est constitué d’une collection d’objets totalement disjoints. Chaque objet est isolé puis mesuré selon un critère de taille (aire, périmètre diamètre de Féret…). Le résultat de la mesure permet de mettre cet objet dans une classe de taille.

• Les problèmes du masque de mesure

Pour faire les mesures qu’on vient de citer, il est nécessaire que l’objet soit totalement inclus dans le champ de mesures. On doit donc éliminer ceux qui coupent le bord du champ. On a vu que cela est facilement réalisé par morphologie mathématique . Cependant, on a d’autant plus de chance d’éliminer un objet que sa taille est grande. Cela va introduire un biais dans l’analyse granulométrique. Pour résoudre ce problème, il faut pouvoir connaître la probabilité qu’un objet $X_{i}$ a d’être inclus dans le champ $Z$. Or, on a vu que $\epsilon _{B}(X)$ donne l’ensemble des points de $X$ où $B$ est totalement inclus dans $X$.

Particule individuelle et rectangle circonscrit associé

Ce raisonnement peut être transcrit pour résoudre notre problème en cherchant à éroder le masque rectangulaire $Z$ par $X_{i}$. Il est facile de voir qu’on obtiendra rigoureusement le même résultat si on remplace $X_{i}$ par le rectangle circonscrit minimum $R_{i}$ de même orientation que $Z$. La probabilité d’inclusion $\Pr(X)$ se calcule alors facilement : $\Pr(X_{i})={\frac {A(\epsilon _{R_{i}}(Z))}{A(Z)}}={\frac {(a_{Z}-a_{R})\times (b_{Z}-b_{R})}{a_{Z}\times b_{Z}}}$ Dans cette expression, $a$ représente le côté horizontal et $b$ le côté vertical du champ (indice Z) ou du rectangle (indice R). Le biais sera alors corrigé en incrémentant la classe de taille non pas de 1 mais de $1/\Pr(X_{i})$. Cette méthode corrective a été proposée par Lantuéjoul.

Analyse granulométrique par ouverture avec un élément structurant à 2 dimensions

Le milieu complémentaire des objets de la figure 1 ne peut pas être traité par cette méthode car la notion d’objet individuel n’a plus de sens. Cependant les axiomes de Matheron se vérifient lorsque l’ouverture est faite avec un élément structurant convexe. En effet, un élément structurant convexe permet de construire une famille de même nature dont tous les membres se déduisent de l’élément de taille 1 par un rapport d’homothétie de taille $\lambda$. Ce type de granulométrie est une granulométrie en mesure car l’ouverture ne possède pas de bonnes propriétés topologiques (un objet peut être scindé en deux par ouverture). Pour une image définie dans $\mathbb {R} ^{2}$, la seule mesure utilisée est l’aire de l’ensemble ouvert.

• 
  Ensemble Booléen de disques circulaires (X= BD)
• 
  Ouverture hexagonale de taille 5 pixels sur BD et masque érodé (cyan)
• 
  Ouverture hexagonale de taille 20 pixels sur BD et masque érodé (cyan)

En analyse locale, il faut tenir compte du masque de mesure et donc travailler dans un masque érodé, on a donc : $G(X,\lambda B)={\frac {A_{A}(X\cap \epsilon _{2\lambda B}(Z))-A_{A}(\gamma _{\lambda B}(X)\cap \epsilon _{2\lambda B}(Z))}{A_{A}(X\cap \epsilon _{2\lambda B}(Z)}}$ Il existe un cas particulier où la granulométrie peut être établie en nombre. C’est lorsque l’ensemble ${\hat {X}}$ est constitué d’objets convexes disjoints. Dans ce cas, on aura : $F({\hat {X}},\lambda B)={\frac {N_{A}({\hat {X}}\cap \epsilon _{2\lambda B}(Z))-N_{A}(\gamma _{\lambda B}({\hat {X}})\cap \epsilon _{2\lambda B}(Z))}{N_{A}({\hat {X}}\cap \epsilon _{2\lambda B}(Z)}}$

Analyse granulométrique par ouverture avec un élément structurant linéaire

L’élément structurant linéaire est traditionnellement noté $\ell$. Les règles s’appliquent de la même manière que pour l’ouverture bidimensionnelle mais ici les granulométries en mesure et en nombre peuvent toujours se calculer puisque l’intersection d’un ensemble $X$ par une droite $\Delta$ donne toujours des segments de droites convexes par définition. Les granulométries correspondantes sont données par les expressions suivantes pour les granulométries en mesure et les granulométries en nombre :

• Granulométries en mesure$G(X,\ell )={\frac {L_{L}(X\cap \epsilon _{2\ell }(Z))-L_{L}(\gamma _{\ell }(X)\cap \epsilon _{2\ell }(Z))}{L_{L}(X\cap \epsilon _{2\ell }(Z)}}$

• Granulométries en nombre

$F(X,\ell )={\frac {N_{L}(X\cap \epsilon _{2\ell }(Z))-N_{L}(\gamma _{\ell }(X)\cap \epsilon _{2\ell }(Z))}{N_{L}(X\cap \epsilon _{2\ell }(Z)}}$

La fonction P(l)

En fait, il n’est pas nécessaire de passer par l’ouverture pour obtenir ces granulométries, mais que l’on peut s’arrêter à l’érosion qui donne la fonction $P(\ell )$. Cette fonction est définie par : $P(\ell )={\frac {A(\epsilon _{\ell }(X)\cap \epsilon _{\ell }(Z))}{A(\epsilon _{\ell }(Z))}}$

Cette fonction possède un certain nombre de propriétés remarquables :

• La valeur à l’origine n’est autre que la teneur de la phase.$P(0)=A_{A}(X)$

• La fonction $P(\ell )$ décroît de manière monotone.
• La dérivée à l’origine, dans la direction $\alpha$, est égale, au signe près, au nombre de connexité spécifique pour cette même direction.$P'(0,\alpha )={\frac {d(P(0,\alpha ))}{d\ell }}=-N_{L}(X,\alpha )$

• La fonction $P(\ell )$ possède des propriétés granulométriques.

• 
  Ensemble Booléen à grains Poissonniens (X = BP)
• 
  Érosion linéaire de taille 10 pixels sur BP et masque érodé (cyan)]
• 
  Érosion linéaire de taille 20 pixels sur BP et masque érodé (cyan)

D’après les relations précédentes, on a immédiatement : $F(X,\ell ,\alpha )={\frac {P'(0)-P'(\ell ,\alpha )}{P'(0)}}$ et $G(X,\ell ,\alpha )={\frac {P(0)-P(\ell ,\alpha )+\ell \times P'(\ell ,\alpha )}{P(0)}}$

• La fonction $P(\ell )$ a une propriété stéréologique avec la ‘’’notion d’étoile’’’.

Étoile bi et tridimensionnelle

Supposons que l’ensemble $X$ soit transparent et que le complémentaire soit opaque. À partir d’un point $x$ appartenant à $X$, on peut définir un domaine $Y_{x}$, constitué par tous les points y visibles de x. $Y_{x}$ sera appelé ’’’ l’étoile de dimension 2’’’ associée au point x.

Illustration de la fonction étoile de X dans R2 (X : jaune et cyan, St(X) : jaune)

En recommençant la même opération pour tous les points de $X$, on peut définir dans $/R^{2}$ une étoile moyenne caractérisée par son aire. Cela s’écrit : $St_{2}(X)={\frac {\int _{X}^{}{A(Y_{x})dx}}{A(X)}}$ Considérons l’élément de surface $\ell \times d\ell$ orienté selon $\alpha$. Cet élément appartient à l’étoile $Y_{x}$ et aura pour probabilité conditionnelle le rapport : $P(\ell ,\alpha )/P(0)$ En utilisant la définition et en supposant le milieu isotrope, on a : $St_{2}(X)={\frac {2\pi }{P(0)}}\int _{0}^{\infty }{\ell \times P(\ell )d\ell }$ Le même raisonnement peut être fait pour $\mathbb {R} ^{3}$. L’étoile dans $\mathbb {R} ^{3}$ est définie par : $St_{3}(X)={\frac {\int _{X}^{}{V(Y_{x})dx}}{V(X)}}$ Ce qui donne dans le cas isotrope : $St_{3}(X)={\frac {4}{P(0)}}\int _{0}^{\infty }{\ell ^{2}\times P(\ell )d\ell }$ L’étoile dans $\mathbb {R} ^{3}$ définit un volume moyen en mesure et dans $\mathbb {R} ^{2}$, une aire moyenne en mesure. Si $X$ est une union de convexes disjoints, l’étoile représente un ensemble convexe moyen. Puisqu’elle est mesurable à partir de $\mathbb {R} ^{1}$, l’étoile a des propriétés stéréologiques.

Analyse de la dispersion

L’étude de la dispersion suppose au moins un ensemble et son complémentaire, tous les deux non vide. Les paramètres stéréologiques qui ont été définis ne portent que sur un seul ensemble, l’analyse granulométrique également. En morphologie mathématique, il existe une fonction qui permet effectivement de tester l’état de dispersion d’un ensemble dans un autre. Cette fonction s’appelle la ‘’’fonction covariance’’’. Elle correspond à la mesure de l’érodé par un élément structurant constitué de deux points distants de h. Comme l’érosion est construite à partir de la soustraction de Minkowski, il est facile d’obtenir le résultat de l’érosion par h puisque cet élément structurant ne contient que 2 points distants de h. Pour cela, il suffit de translater l’image et de prendre l’intersection avec le translaté.

Illustration de l'érosion par un bipoint h obtenue par translation et intersection

On aura donc : $\epsilon _{h}(X)=X\cap X_{\vec {h}}$

La covariance simple

La covariance est surtout utilisée dans le cas local. Dans ce cas, on définit la fonction covariance dans un masque de mesure Z par : $C(X,h)={\frac {A(\epsilon _{h}(X)\cap \epsilon _{h}(Z))}{A(\epsilon _{h}(Z))}}$

Propriétés de la covariance

Comme la fonction $P(\ell )$, la fonction $C(h)$ possède un certain nombre de propriétés. Ainsi, on a :

• La valeur à l’origine n’est autre que la teneur de la phase.$C(0)=P(0)=A_{A}(X)$
• La dérivée à l’origine.$C'(0)=P'(0)=-N_{L}(X)$
• Entre l’origine et l’infini, la covariance est toujours inférieure à la valeur à l’origine. Les variations dans ce domaine dépendent de l’état de dispersion. En particulier, la présence de périodicité se traduit par des oscillations de la covariance et des paliers intermédiaires. Des paliers intermédiaires apparaissent pour des superpositions de structures à des échelles différentes. Des exemples sont donnés dans la section suivante.
• Dans le cas d’un ensemble aléatoire stationnaire, la valeur à l’infini correspond à la probabilité pour que 2 points indépendants tombent dans $X$. C’est donc égal au carré de la teneur :$C(\infty )=(A_{A}(X))^{2}$ Cette valeur est atteinte de manière asymptotique.
• Comme la covariance a un caractère directionnel, elle renseigne sur l’anisotropie de l’ensemble $X$.

Exemples

À titre d’exemple, nous prenons des cas limites :

• Un ensemble 2D isotrope (schéma booléen).

Dans l’exemple choisi la décroissance de la covariance se poursuit jusqu’à la valeur asymptotique. Comme l’analyse n’a été conduite que sur un champ, la concordance entre l’asymptote théorique et l’asymptote expérimentale n’est pas parfaite. Le léger passage par un minimum montre un petit effet de répulsion entre les disques.

• Un ensemble anisotrope composé de lamelles irrégulières (image d’un eutectique métallique).

L’aspect périodique de la structure se traduit par des oscillations de la courbe de covariance. Le premier minima correspond à une épaisseur moyenne d’une lamelle et le premier maximum à l’épaisseur moyenne du couple lamelle-complémentaire.

• Un ensemble constitué d’amas.

La covariance est plus complexe à interpréter mais on peut estimer ma distance moyenne entre amas.

• 
  Érosion par h(10 pixels) sur un ensemble Booléen de disques (X :jaune et rouge, érosion par h : rouge, érosion du masque : cyan)
• 
  Covariance sur l'ensemble Booléen de disques
• 
  Érosion par h(24 pixels) sur un ensemble lamellaire (X :jaune et rouge, érosion par h : rouge, érosion du masque : cyan)
• 
  Covariance sur l'ensemble lamellaire
• 
  Érosion par h(24 pixels) sur un ensemble d'amas (X :jaune et rouge, érosion par h : rouge, érosion du masque : cyan)
• 
  Covariance sur l'ensemble d'amas

Les modèles probabilistes (Ensembles aléatoires)

Les paramètres stéréologiques sont peu nombreux et la caractérisation de la taille par des fonctions granulométriques a un caractère stéréologique uniquement pour les granulométries linéaires. Le même raisonnement peut être tenu pour la fonction covariance.

Si on accepte de rester dans l’espace, $\mathbb {R} ^{2}$ les autres granulométries, l’étude de la forme et de l’anisotropie permettent de compléter les informations morphologiques. Avec ces méthodes, on se trouve en possession d’un grand nombre de paramètres et de plusieurs distributions, mais en contrepartie, la lisibilité de la caractérisation n’est plus évidente.

Les modèles probabilistes constituent une autre approche beaucoup plus synthétique. Ils sont conçus pour caractériser des ensembles aléatoires. Bien entendu, tous les ensembles réels ne peuvent pas être décrits par de tels modèles. En premier lieu, ces ensembles doivent être, spatialement stationnaires , pour être facilement modélisables.

Généralités sur les modèles probabilistes

Pour que l’on puisse observer une morphologie, il faut qu’un ensemble $X$ ne remplisse pas tout l’espace. On aura donc au minimum un milieu à deux composantes: l’ensemble $X$ et son complémentaire $X^{C}$. Les objets constituant cet ensemble peuvent être aussi bien des points que des droites ou des sous-ensembles quelconques. L’ensemble $X$, ainsi obtenu, sera donc un ensemble aléatoire topologiquement fermé dénommé RACS (RAndom Closed Set) et bien décrit dans les livres de Matheron, Serra, Stoyan et Jeulin. Cette restriction est importante pour conserver de bonnes propriétés, mais peu gênante dans la pratique. Ainsi, un RACS restera un RACS après une opération d’érosion de dilatation… Le choix d’un modèle dépend d’un certain nombre de connaissances a priori. Dans le cas des matériaux ou systèmes, le nombre de phases morphologiquement discernables est le premier élément. On aura donc deux grandes catégories :

• Les matériaux ou systèmes monophasés décrits à partir d’une partition de l’espace.
• Les matériaux ou systèmes polyphasés décrits à partir des RACS polyphasés.

La capacité de Choquet

Un RACS peut être caractérisé par une probabilité d’événements correspondant à des mesures morphologiques comme la probabilité d’inclusion d’un compact $K$ dans un ensemble ou son complémentaire. C’est le rôle attribué à la capacité de Choquet $T(K)$ définie par :

$T(K)=Pr(X\cap K\neq \emptyset )$

On peut aussi définir $T(K)$ à partir de la probabilité $Q(K)$ que l’intersection entre $X$ et $K$ soit vide :

$Q(K)=1-T(K)=Pr(X\cap K=\emptyset )=Pr(K\subset X^{C})$

Illustration de la Capacité de Choquet

De la même façon qu’une fonction de distribution définit une variable aléatoire, la connaissance de la capacité de Choquet pour tout compact $K$ permet de définir complètement un modèle probabiliste. Il est évident que tous les compacts possibles ne peuvent pas être testés. On se contentera des plus simples.

Propriétés des RACS

Les RACS peuvent vérifier ou non un certain nombre de propriétés.

L’infini divisibilité

Un RACS $X$ est infiniment divisible si $X$ est équivalent à l’union de n RACS indépendants et de même nature. L’intersection d’un modèle infiniment divisible par un sous espace conserve la nature du modèle engendré. Ceci constitue une caractéristique de nature stéréologique. Pour un tel RACS et un compact $K$ donné, la capacité de Choquet est une expression de la forme :

$Q(K)=1-T(K)=\exp(-\Psi (K))$ avec : $\Psi (\emptyset )=0$

La stabilité

Un RACS $X$ infiniment divisible est stable par union si la fonction $\Psi (\lambda K)$ vérifie l’équation :

$\Psi (\lambda K)=\lambda ^{\alpha }\times \Psi (K)$ avec $(\lambda >0,\alpha >0)$

La calculabilité

Un RACS possède la propriété de calculabilité si pour certains compacts $K$, les capacités de Choquet possèdent des formules explicites. Ceci permettra de vérifier si une structure réelle peut correspondre à la réalisation d’un RACS sans effectuer de simulation. Lorsque la capacité de Choquet n’est pas calculable, on peut utiliser des grandeurs caractéristiques liées aux paramètres du modèle.

Le processus ponctuel de Poisson

Le point de départ de tous les modèles probabilistes est un processus ponctuel aléatoire. Il faut un processus où le nombre de points tombant dans un sous-ensemble $Z1$ soit indépendant de celui tombant dans $Z2$. Le processus suivant la loi de Poisson vérifie cette condition. La probabilité pour que n points d’un processus de Poisson de densité $\theta$ appartiennent à un ensemble $Z$ est donné par :

$\Pr(Z)={\frac {(\theta \times mes(Z))^{n}}{n!}}\times \exp(-\theta \times mes(Z))$

Les modèles ensemblistes de base

Les partitions de l’espace

Les modèles de partition aléatoire sont des ensembles qui découpent l’espace $\mathbb {R} ^{n}$ en plusieurs sous-ensembles fermés et bornés appelés classes. L’union de tous ses sous-ensembles remplit tout l’espace $\mathbb {R} ^{n}$. Les principaux modèles de partition sont : la partition de Voronoï, la partition de Johnson Mehl, la mosaïque de Poisson et les feuilles mortes. Si tous ces modèles ont pour origine un processus ponctuel de Poisson, leur construction et leurs propriétés sont très différentes. Le dernier sera présenté après le modèle de Boole Matheron qui est un modèle polyphasé.

La partition de Voronoï

Pour construire une partition ou diagramme de Voronoï, on implante des points $p_{i}$ selon un processus obéissant à la loi de Poisson de densité $\theta$. À chaque point $p_{i}$, correspond une zone d’influence $X_{i}$ qui est définie par :

$X_{i}=\left\{x:d(x,p_{i})<d(x,p_{i\neq j})\right\}$ avec $d(x,p)$ la distance de $x$ à $p$

Partition de Voronoï (frontière bleue) avec 100 germes poissoniens (rouge)

Cette zone d’influence est un polygone convexe dans $\mathbb {R} ^{2}$ et un polyèdre convexe dans $\mathbb {R} ^{3}$.

La simplicité du modèle a conduit de nombreux auteurs à l’utiliser pour décrire des structures cellulaires ou granulaires. Mais, la partition de Voronoï n’est pas infiniment divisible puisqu’une partition de Voronoï dans $\mathbb {R} ^{3}$ n’engendre pas une partition de Voronoï dans $\mathbb {R} ^{2}$. De plus, on ne connaît pas d’expression analytique de la capacité de Choquet de cette partition pour les compacts usuels. Pour comparer une section plane de la structure réelle à un modèle de Voronoï dans $\mathbb {R} ^{2}$, on dispose de quelques caractéristiques résumées par Miles, fonction de la densité de Poisson $\theta$.

• Paramètres stéréologiques dans $\mathbb {R} ^{2}$
  • Périmètre spécifique : $L_{A}(X)=2\times {\sqrt {\theta _{2}}}$
  • Nombre spécifique de cellules : $N_{A}(X)=\theta _{2}$
  • Nombre spécifique de sommets : $P_{A}(X)=2\times \theta _{2}$

• Paramètres moyens associés à une cellule $X'$
  • Aire moyenne : ${\bar {A}}(X')=\theta _{2}^{-1}$
  • Nombre moyen de sommets : ${\bar {P}}(X')=6$

Dans le cas d’une partition de Voronoï dans $\mathbb {R} ^{3}$, il existe des relations similaires utilisant la fonction $\Gamma$ et la densité de Poisson $\theta _{3}$.

• Paramètres stéréologiques dans $\mathbb {R} ^{3}$
  • Surface spécifique : $S_{V}(X)=4\times (\pi /6)^{1/3}\times \Gamma (5/3)\times \theta _{3}^{1/3}$
  • Longueur spécifique des joints triples $J_{3}$ : $L_{V}(J_{3})=(16/15)\times (3/4)^{1/3}\times \pi ^{5/3}\times \Gamma (4/3)\times \theta _{3}^{2/3}$
  • Nombre spécifique de cellules : $N_{V}(X)=\theta _{3}$

De plus, les densités de la partition et les fonctions "étoile" , $St_{3}$ ,$St_{2}$ et $St_{1}$ sont reliées à la densité du processus de Poisson $\theta$ par les relations suivantes :

• Densité de la partition dans $\mathbb {R} ^{1}$ : $\theta _{1}=1,45\times \theta _{3}^{1/3}$
• Densité de la partition dans $\mathbb {R} ^{2}$ : $\theta _{2}=1,46\times \theta _{3}^{2/3}$
• Étoile dans $\mathbb {R} ^{1}$ : $St_{1}(X)=0,95\times \theta _{3}^{-1/3}$
• Étoile dans $\mathbb {R} ^{2}$ : $St_{2}(X)=1,04\times \theta _{3}^{-2/3}$
• Étoile dans $\mathbb {R} ^{3}$ : $St_{3}(X)=1,24\times \theta _{3}^{-1}$

Ces fonctions étoile sont calculées à partir des moments de la fonction $P(\ell )$, (probabilité d’inclusion d’un segment $\ell$ dans un grain).

On a :

• Densité de la partition dans $\mathbb {R} ^{1}$ : $\theta _{1}=1,45\times \theta _{3}^{1/3}$,
• Densité de la partition dans $\mathbb {R} ^{2}$ : $\theta _{2}=1,46\times \theta _{3}^{2/3}$,
• Étoile dans $\mathbb {R} ^{1}$ : $St_{1}(X)=0,95\times \theta _{3}^{-1/3}$,
• Étoile dans $\mathbb {R} ^{2}$ : $St_{2}(X)=1,04\times \theta _{3}^{-2/3}$,
• Étoile dans $\mathbb {R} ^{3}$ : $St_{3}(X)=1,24\times \theta _{3}^{-1}$.

Le modèle de Johnson Mehl

Le modèle de Johnson-Mehl est également basé sur un processus ponctuel de Poisson. Cependant, le modèle est séquentiel (fonction du temps). Chaque séquence $dt$ est composée de deux processus élémentaires :

• Création de germes ponctuels selon un processus de Poisson de densité $\theta =\mu \times dt$
• Epaississement conservant l'homotopie des germes, d’épaisseur $v\times dt$, selon une vitesse $v$.

Cependant, tous les germes créés ne vont pas nécessairement donner naissance à un « grain ». Si le germe apparaît dans un grain déjà formé, il disparaît. La construction s’arrête lorsque le complémentaire des grains a complètement disparu. Les grains constituant la partition ont des frontières hyperboliques dans $\mathbb {R} ^{2}$ et hyperboloïdes dans $\mathbb {R} ^{3}$. Elles ne sont donc pas toujours convexes, mais on connaît leurs caractéristiques, en particulier la fonction de distribution du nombre de voisins.

• 
  Construction du modèle de Johnson Mehl (quelques étapes)
• 
  Construction du modèle de Johnson Mehl (quelques étapes)
• 
  Construction du modèle de Johnson Mehl (quelques étapes)
• 
  Construction du modèle de Johnson Mehl (résultat final)

Comme dans le cas de la partition de Voronoï, le modèle de Johnson-Mehl ne possède pas de propriétés stéréologiques. Dans le cas ou $\mu$ est constant, il existe des relations liées à la fonction $P(\ell )$ :

• Densité de la partition dans $\mathbb {R} ^{1}$ : $\theta _{1}=1,28\times \theta _{3}^{1/3}$,
• Densité de la partition dans $\mathbb {R} ^{2}$ : $\theta _{2}=1,22\times \theta _{3}^{2/3}$,
• Étoile dans $\mathbb {R} ^{1}$ : $St_{1}(X)=1,10\times \theta _{3}^{-1/3}$,
• Étoile dans $\mathbb {R} ^{2}$ : $St_{2}(X)=1,50\times \theta _{3}^{-2/3}$,
• Étoile dans $\mathbb {R} ^{3}$ : $St_{3}(X)=2,20\times \theta _{3}^{-1}$.

La partition de Poisson

La partition de l’espace $\mathbb {R} ^{2}$ selon un processus de Poisson est faite par des droites de Poisson. Les droites de Poisson sont construites de la manière suivante. Soit une droite $\Delta$ d’orientation comprise entre $\omega$ et $\omega +d\omega$ et passant par l’origine du plan. Sur cette droite, on réalise un processus ponctuel de Poisson de densité $\lambda ^{2}\times d\omega$. En chacun de ces points, on implante une droite de Poisson perpendiculaire à $\Delta$. Dans le cas d’une mosaïque isotrope, $\lambda$ est constant et la valeur de $\omega$ est choisie selon une loi de probabilité uniforme.

Construction des droites de Poisson

L’espace $\mathbb {R} ^{2}$ est alors partagé en une infinité de polygones aléatoires appelés polygones de Poisson.

Mosaïque de Poisson générée par 20 droites de Poisson

Une construction similaire dans l’espace $\mathbb {R} ^{3}$ conduira à un espace partagé en une infinité de polyèdres de Poisson. L’angle $\omega$ est alors compris entre 0 et $2\pi$ stéradians. Les droites de Poisson sont remplacées par des plans de Poisson perpendiculaires à $\Delta$ selon une densité $\lambda _{2}\times d\lambda$.

À la différence du modèle de Voronoï, la mosaïque de Poisson possède des propriétés stéréologiques. D’une part, le paramètre $\lambda _{3}$ du modèle permet de caractériser le polyèdre moyen $X'$ par son volume moyen, ${\bar {V}}(X')$, sa surface moyenne, ${\bar {S}}(X')$ et son intégrale de courbure moyenne, ${\bar {M}}(X')$ . On a en effet les relations :

• ${\bar {V}}(X')=6/(\pi ^{4}\times \lambda ^{3})$
• ${\bar {S}}(X')=24/(\pi ^{3}\times \lambda ^{2})$
• ${\bar {M}}(X')=3/\lambda$

D’autre part une mosaïque de Poisson dans $\mathbb {R} ^{3}$, de paramètre $\lambda _{3}$, intersectée par un plan, génère une mosaïque de Poisson dans $\mathbb {R} ^{2}$ de paramètre $\lambda _{2}$ avec :

$\lambda _{2}=\pi \times \lambda _{3}/2$

La mosaïque de Poisson est rarement utilisée pour modéliser une partition de l’espace. Par contre, elle permet de générer des grains aléatoires pour des modèles polyphasés.

Un dernier mode de partition ne sera pas abordé dans cette section, il s’agit du modèle de feuilles mortes que nous verrons plus en détail dans la section suivante.

Les RACS multiphasés

Les ensembles fermés aléatoires multiphasés constituent un groupe très important. On peut les classer en trois groupes.

• Le premier groupe est le ‘‘ modèle Boole Matheron ‘‘. C’est le plus connu de tous au même titre que la partition de Voronoï dans la catégorie précédente.
• Le deuxième groupe rassemble les ‘‘modèles de feuilles mortes‘‘ qui, comme le ‘‘modèle de Johnson Mehl’’ sont des modèles séquentiels.
• Le troisième groupe est une extrapolation, dans le cadre polyphasé, des partitions de l’espace. Ce sont les partitions multiphasées.

Le modèle de Boole Matheron

Ce modèle, également appelé schéma Booléen, est construit de la manière suivante. En chaque point d’un processus de Poisson, de densité $\theta$, on place un grain primaire. Le schéma Booléen est l’union de tous ces grains primaire (figure de gauche).

• 
  Ensemble de Boole Matheron (jaune), (30 germes, disques taille 20)
• 
  Ensemble de Boole Matheron (jaune), (disques de taille variable)
• 
  Ensemble de Boole Matheron à grains Poissoniens (jaune)

Pour cet ensemble, on a généré un processus ponctuel et remplacé chaque point par un dodécagone de taille unique (grain primaire). La seconde figure représente un modèle de Boole Matheron construit avec des disques de taille variable. Pour la dernière figure, les grains primaires du modèle sont des polygones de Poisson, obtenus par tirage au sort à partir d’une partition comme nous l’avons dans la section précédente.

Propriétés

Le modèle de Boole Matheron possède de très bonnes propriétés. Il est infiniment divisible, stable et calculable. En effet, si $X'$ est le grain primaire du schéma Booléen, on a la relation :

$Q(K)=\exp(-\theta \times E(Mes(\delta _{K}(X'))))$

$E(Mes(\delta _{K}(X'))$ est l’espérance de la mesure de Lebesgue de l’ensemble $X$ dilaté par le compact $K$.

La capacité de Choquet peut encore s’écrire :

$Q(K)=\Pr(x\in \epsilon _{K}(X^{C}))$

$\epsilon _{K}(X^{C})$ est l’érosion du complémentaire de $X$ par le compact $K$. Pour tester un modèle de Boole Matheron par la capacité de Choquet (ou la fonctionnelle complémentaire), il suffit d’estimer les teneurs des érodés par une ou plusieurs familles de compacts $K$. Chaque famille étant définie par l’ensemble de compacts homothétiques. Un modèle de Boole Matheron sera défini par la densité de Poisson $\theta$ et le grain primaire caractérisé par une forme et une distribution de taille.

Si le grain primaire est convexe et de géométrie simple, le modèle de Boole Matheron sera calculable par les paramètres stéréologiques du grain moyen $X'$ : ${\overline {V}}(X'),{\overline {S}}(X'),{\overline {M}}(X')$. Dans le cas de grains sphériques, les paramètres stéréologiques sont calculables à partir des moments 3, 2 et 1 de la distribution granulométrique $f(a)$. On a, en effet, les relations :

• ${\overline {V}}(X')={\frac {\pi }{6}}\times E(a^{3})$
• ${\overline {S}}(X')=\pi \times E(a^{2})$
• ${\overline {M}}(X')=2\pi \times E(a)$

Enfin, comme une union de modèles de Boole Matheron est toujours un modèle de Boole Matheron, on dispose de beaucoup de nombreuses solutions pour aborder la modélisation d'une structure réelle.

Certains compacts sont particulièrement intéressants pour tester un modèle de Boole Matheron. Ce sont les points $(K=P)$, les segments $(K=\ell )$, les hexagones de taille r $(K=H(r))$, le bipoint $(K=h)$ et, pour certains modèles, le triplet de points défini par les sommets d’un triangle équilatéral .$(K=T(h))$ Pour les compacts convexes et en appelant $q$ la teneur de l’ensemble complémentaire, on a les relations suivantes :

• $Q(P)=q=\exp(-\theta \times {\overline {V}}(X'))$
• $Q(\ell )=q\times \exp(-\theta \times \ell \times {\frac {{\overline {S}}(X')}{4}})$
• $Q(H(r))=q\times \exp(-\theta [-{\frac {3r}{4}}\times {\overline {S}}(X')+{\frac {3{\sqrt {3}}}{8}}\times r^{2}\times {\overline {M}}(X')])$

Pour le bipoint $(K=h)$, on utilise le covariogramme géométrique moyen $K(h)$ à la place de la covariance $Q(h)$. Le covariogramme est lié à la covariance par l’expression :

• ${\frac {K(h)}{K(0)}}=2-{\frac {\log Q(h)}{\log(q)}}$

Dans le cas de grains de Poisson, on a :

• ${\frac {K(h)}{K(0)}}={\frac {K(h)}{V(X')}}=\exp(-\pi \times \lambda \times \ell )$

Lorsqu’on a des grains sphériques, le covariogramme géométrique est fonction de la distribution. En appelant $a$, la taille maximale du grain primaire, on a :

• $K(h)={\frac {\pi }{6}}\times \left[\int _{h}^{\infty }{a^{3}\times f(a)\times da}-{\frac {3h}{2}}{\int _{h}^{\infty }{a^{2}\times f(a)}\times da}+{\frac {h^{3}}{2}}{\int _{h}^{\infty }{f(a)}\times da}\right]$

Le modèle de feuilles mortes

Le modèle de feuilles mortes est un schéma booléen séquentiel. La version monophasée du modèle est due à Jeulin. La construction est la suivante. Des grains primaires sont générés par un processus de Poisson de densité $\theta (t)$. À la différence du schéma booléen, les grains peuvent se recouvrir. Les plus anciens peuvent disparaître sous les plus récents. Le processus peut être arrêté au bout d’un temps t. Si le support n’est pas entièrement recouvert, le processus ressemble un peu à un schéma booléen. On peut aussi poursuivre jusqu’à stationnarité. La partition recouvre alors intégralement le support.

Dans le cas d’un modèle de ‘‘feuilles mortes biphasé‘‘, les grains primaires de la phase 1 et de la phase 2 apparaissent successivement avec les densités respectives $\theta _{1}(t)$ et $\theta _{2}(t)$. Le processus est répété jusqu’à stationnarité. Il reste une structure biphasée, où les deux phases sont imbriquées les unes dans les autres.

• 
  Construction d'un modèle de feuilles mortes biphasé avec des disques bleus et jaunes
• 
  Fin de la construction d'un modèle de feuilles mortes biphasé avec des disques bleus et jaunes

Ces modèles sont infiniment divisibles et donc ils génèrent des modèles équivalents dans des sous espaces. La calculabilité n’est pas aussi forte que dans le cas du schéma booléen. Les capacités de Choquet ne sont calculables que pour le bi point distant de h $(K=h)$ et pour le triplet de points distant de h $(K=T(h))$ . Lorsque le bi point est testé, on utilise la fonction $\phi (h)$ définie par :

$\Phi (h)={\frac {2-r_{1}(h)-r_{2}(h)}{2-p\times r_{1}(h)-q\times r_{2}(h)}}$

avec $p$ teneur de la phase 1 et $q$ celle de la phase 2 et $r_{i}={K_{i}(h)}/{K_{i}(0)}$

À titre d’illustration, on présente quelques réalisations avec des disques circulaires ou des polygones de Poisson comme grains primaires.

• 
  Modèle de feuilles mortes biphasé avec des disques bleus et jaunes
• 
  Modèle de feuilles mortes biphasé avec des grains Poissonniens bleus et jaunes

Les partitions de l’espace multiphasées : la partition de Poisson

Pour construire une partition de Poisson polyphasée, il faut construire une partition monophasée et attribuer les classes à une phase donnée de manière aléatoire. Les paramètres du modèle (Miles) sont la teneur de chaque phase et $\lambda$ qui caractérise la partition de Poisson. Dans le cas de la partition de Poisson biphasée, les propriétés du système monophasé sont conservées.

Mosaique de Poisson biphasée (bleue et jaune)

Les lois analytiques sont connues pour $(K=h)$ , $(K=\ell )$ et $(K=T(h))$

• $P(\ell )=p\times \exp(-q\times \pi \times \lambda \times \ell )$
• $Q(\ell )=q\times \exp(-p\times \pi \times \lambda \times \ell )$
• $C(h)=p^{2}+p\times q\times \exp(-\pi \times \lambda \times h)$
• $Q(T(h))=q^{3}-q\times (1-3q+2q^{2})\times \exp(-{\frac {3}{2}}\times \pi \times \lambda \times h)$

Le modèle de Stienen

Dans un modèle de Boole Matheron ou un modèle de feuilles mortes, les grains primaires peuvent se recouvrir. Pour construire le modèle de Stienen, on part toujours d’un processus ponctuel de Poisson de densité $\theta$. Mais chaque point sera remplacé par la plus grande sphère contenue dans la cellule de Voronoï correspondante. Dans ces conditions, les sphères ne se recouvrent pas mais peuvent se toucher (figure de gauche). Il a été généralisé en réduisant la taille des sphères par un facteur $\alpha [0,1]$ (figure de droite).

• 
  Modèle de Stienen avec alpha = 1 (rouge) et partition de Voronoï associée (frontière bleue)
• 
  Modèle de Stienen avec alpha < 1 et partition de Voronoï associée (frontière bleue)

Dans le cas du modèle initial ($\alpha =1$), la teneur des sphères est constante puisque l’on a :

$\alpha =2\times (V_{V}(X))^{{1}/{3}}$

De plus, la distribution des sphères est connue puisqu’elle est directement liée à la distribution des distances des plus proches voisins d’un processus de Poisson. On a en effet :

$F(D)=1-\exp(-{\frac {\pi \times \theta }{6\times V_{V}(X)}}\times D^{3})$

Pour $\alpha \neq 1$, les sphères ne sont plus en contact (figure 14). On peut calculer par intégration numérique la fonction de corrélation d’une paire de points en fonction de la distance r. Enfin, pour ce modèle, on ne dispose que d’une expression complexe de la covariance.

Bibliographie

En français

• Georges Matheron, Éléments pour une théorie des milieux poreux, Masson, Paris, 1967.
• Michel Schmitt et Juliette Mattioli, Morphologie mathématique, Masson, Paris, 1993, (ISBN 2-225-84385-6) .
• Laurent Najman et Hugues Talbot (dir.), Morphologie mathématique, 1 : Approches déterministes, Hermès - Lavoisier, Paris, 2008, (ISBN 978-2746218413) .
• Laurent Najman et Hugues Talbot (dir.), Morphologie mathématique, 2 : Estimation, choix et mise en œuvre, Hermès - Lavoisier, 2010.
• Jean Marc Chassery, Annick Montanvert, Géométrie discrète en analyse d'images, Hermès, 1991, (ISBN 2-86601-271-2).
• H. Poincaré Calcul des probabilités. Carré, Paris, 1912.
• Michel Coster, Jean Louis Chermant, Précis d'analyse d'images, Presse du CNRS, 1989, (ISBN 2-87682-020-X) .
• Robert T. DeHoff, Frederick N. Rhines, Microscopie Quantitative, (traduit par Jean Montuelle), Masson 1972.

En anglais

• Georges Matheron, Random Sets and Integral Geometry, Wiley, New York, 1975.
• Georges Matheron, Estimating and Choosing, Springer–Verlag Berlin, Heidelberg, 1989.
• Jean Serra, Image Analysis and Mathematical Morphology, vol. 1, Academic Press, Londres, 1982, (ISBN 0-12-637242-X) .
• Jean Serra (dir.), Image Analysis and Mathematical Morphology, vol. 2 : Theoretical Advances, Academic Press, Londres, 1988, (ISBN 0-12-637241-1) .
• Charles R. Giardina et Edward R. Dougherty, Morphological Methods in Image and Signal Procesing, Prentice-Hall, New Jersey, 1988.
• H.J.A.M Heijmans, Morphological image operators, Academic Press, coll. "Advances in Electronics and Electron Physics", Boston 1994.
• Pierre Soille, Morphological image analysis, Springer-Verlag Berlin, Heidelberg, 1999 (2^e édition 2003).
• Gonzalez R.C., Woods R.E., Digital Image Processing, 3^e éd., Prentice Hall, 2008.
• Goutsias J., Batman S., Morphological Methods for Biomedical Image Analysis, Handbook of Medical Imaging, volume 2 : Medical Image Processing and Analysis, M. Sonka & J.M. Fitzpatrick (eds.), SPIE Optical Engineering Press, 2000, p. 175-272.
• Laurent Najman and Hugues Talbot (Eds).Mathematical morphology: from theory to applications, ISTE-Wiley, (520 pp.) June 2010, (ISBN 978-1-84821-215-2) .
• Henri Poincaré, Papers on Topology: Analysis Situs and its Five Supplements, (traduit par John Stillwell), American Mathematical Society, Providence, R. I., 2010, (ISBN 978-0821852347) .
• E.E. Underwood, Quantitative Stereology, Addison Wesley, 1970.
• Luis A. Santaló, Integral Geometry and Geometric Probability 2nd édition, Cambridge Mathematical Library, 2004, (ISBN 978-0521523448) .
• Ewald R. Weibel, Stereological Methods. Vol. 1: Practical methods for biological morphometry , Academic Press, 1979, (ISBN 978-0127422015) .
• C. Lantuéjoul. On the estimation of mean values in individual analysis of particles, Microscopica Acta, 1980, 5, 266-273.
• Dominique Jeulin (editor), Advances in Theory and Application of Random Sets, World Scientific, 1997, (ISBN 978-9810230012).
• Sung Nok Chiu, Dietrich Stoyan, Wilfrid S. Kendall, Joseph Mecke, Stochastic Geometry and its Applications, Wiley Series in Probability and Statistics, 2013, (Online), (ISBN 978-1118658222).
• D. Jeulin, Multi-component random models for description of complex microstructures, Mikroscopie, 1980, 37, 130-137.
• G. Choquet, Theory of capacities, Ann. Inst. Fourier, 1953, 5, 131-295.
• R.E. Miles, The random division of space, Suppl. Adv. Prob., 1972, 243-266.
• R.E. Miles, Poisson flats in euclidean spaces. Part 2: homogeneous Poisson flats and the complementary theorem, Adv. Prob., 1971, 3, 1-43.
• H. Stienen, The sectionning of randomly dispersed particles, a computer simulation, Mikroskopie (Wien) 1980, S37, 74-78.
• M. Schlater, D. Stoyan, The covariance of the Stienen model. Advances in Theory and Application of Random Sets. World Scientific 1997, 157-174.
• D. Stoyan, K. Wiencek, Spatial correlation in metal structures and their analysis, Materials characterization, Elsevier Sc Pub. 1991, 26, 167-176.

En Allemand

• Hugo Hadwiger, Vorlesungen über Inhalt, Oberfläche und Isoperimetrie, Springer-Verlag, Berlin u. a. 1957, (ISBN 3-540-02151-5) .

Articles connexes

• Centre de morphologie mathématique
• Somme de Minkowski
• Ligne de partage des eaux (segmentation)
• Automate cellulaire

Liens externes

• Centre de morphologie mathématique de MINES ParisTech
• Cours de morphologie mathématique de Jean Serra
• Smil - A Simple (but efficient) Morphological Image Library
• Introduction à la morphologie mathématique par Isabelle Bloch
• Cours de traitement d'image, dont morphologie, par Christian RONSE