Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

Cours MATLAB - 24.11.2009 - Michele Volpi

[email protected]

FGSE Cours MATLAB

Traitement d’Image
Opérateurs Morphologiques
1 Introduction
Le domaine du traitement d’image (en anglais Image Processing) est composé de toute technique visant à élaborer
et analyser quantitativement des images. Aujourd’hui, ces techniques sont entièrement dédiées aux analyses et au
traitement d’images digitales. Due à ce binôme très strict entre image et ordinateur, les techniques d’Image Processing
sont regroupées dans le domaine des sciences informatiques, plus précisément dans la famille du traitement du signal.
Des exemples classiques de traitement d’image sont la mesure d’indicateurs statistiques de texture (distribution
des valeurs des pixels; corrélation, variance, moments angulaires,...), le filtrage d’image (filtres morphologiques, con-
volution, détecteurs d’angles, wavelets,...) et la transformation d’image (composantes principales, fusion d’images,
extraction de caractéristiques,...).
Comme d’autres logiciels scientifiques, MATLAB offre des fonctions spécifiques pour les tâches de traitement
d’image dans la Toolbox d’Image Processing (regroupement de fonctions spécifiques). Pratiquement, les techniques
de traitement d’images consistent en opérations sur des matrices deux dimensionnelles (MxN), qui correspondent à
des images en échelle de gris; matrices cubiques (MxNx3) qui correspondent a des images en couleurs (RGB) et a des
hypercubes (MxNxT) ou d’autres informations sur les couleurs (réflectance ou luminance, quantifié en fourchettes de
longueurs d’ondes réfléchies par les objets) sont enregistrées en chaque élément t ∈ T de l’image I (images multispec-
trales, hyperspectrales, ultraspectrales).
Dans ce TP nous allons nous familiariser avec les technique de morphologie mathématique et les opérateurs de
filtrage qui leur correspondent (érosion, dilatation, ouverture, fermeture et ouverture / fermeture géodésiques (filtres
de reconstruction)).

2 Les Filtres Morphologiques

2.1 Érosion et Dilatation
L’érosion (en MATLAB, imerode) et la dilatation (imdilate) sont deux techniques de base appartenant au domaine
de la morphologie mathématique issus de la théorie des sets visant à filtrer une image. Ce groupe d’opérateurs se base
sur le filtrage des hautes et basses valeurs de l’image I en appliquant un élément structurant B (filtre). L’image érodée
/ dilatée par l’élément structurant B est l’ensemble des pixels x ∈ I (souvent noté IB ) tel que B soit entièrement
contenu (érosion) ou tel que le centre de l’élément soit contenu (dilatation) dans l’objet considéré. Le pixel filtré x̂
prend la valeur booléenne 1 si l’élément est entièrement contenu (ou partialement pour la dilatation) ou 0 sinon, pour
le cas de base ou on considère des images en 1 bit (binaires).
Élement
Structurant

Objet original

Figure 1.1: érosion (gauche) et dilatation (droite) de l’objet binaire en utilisant un élément structurant B circulaire.

Plus précisément érosion ǫB (I) et dilatation δB (I) sont définis (pour une image binaire en noir et blanc) comme:

ǫB (I) = I  B, (1.1)
δB (I) = I  B. (1.2)

1
Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

Figure 1.2: Érosion (gauche) et dilatation (droite) des objets originaux (centre) en utilisant un élément structurant B
circulaire de rayon 2 (lx) et 5 (dx) pixels.

L’extension de ces techniques aux cas de l’image à échelle de gris, nécessite une adaptation par rapport à l’assignation
des nouvelles valeurs de l’image filtrée. Si avant il s’agissait d’opérations booléennes, ces opérations deviennent
arithmétiques assignant la valeur maximale (dilatation) ou minimale (érosion) contenue dans la surface de l’élément
structurant au pixel filtré x̂ sur lequel B est centré. On définit l’érosion ǫB (I) et dilatation δB (I), pour une image à
échelle de gris, comme suit:

ǫB (I) = inf(IB ), (1.3)

δB (I) = sup(IB ). (1.4)

Érosions Image Originale Dilatations

Figure 1.3: Érosions (vers gauche) et dilatations (vers droite) de l’image originale (centre) en utilisant un élément
structurant B circulaire de rayon [2, 4, 6, 8] pixels. (Image example de MATLAB, ’cameraman.tif ’)

D’un point de vue purement pratique, l’érosion (en prenant le minimum de valeurs de gris contenu dans l’élément
B) efface les éléments clairs contenus dans la taille de l’élément. Au contraire, la dilatation efface les éléments sombres
de l’image. (Rappel: plus une valeur est grande, plus elle tend vers le blanc!)
Érosion et Dilatation sont les deux opérations de base dans le domaine de la morphologie mathématique. Par
la suite, on verra comment la combinaison ou l’itération de ces éléments nous amènent à obtenir d’autres filtres et
opérateurs très utilisés lors de l’extraction d’information structurale d’une image ou bien lors d’un lissage, ou encore
pour éliminer des objets dans l’image.
Dans MATLAB ces deux opérations sont contenues dans la toolbox de traitement d’image:

IM2 = imdilate(I, SE): I l’image orginale; SE l’élement structurant B (SE = strel(type,size))

IM2 = imdilate(I, NHOOD): I l’image orginale; NHOOD matrice binaire spécifiant le voisinage de B
IM2 = imerode(I, SE || NHOOD) fonctionne de même manière.

Les deux fonctions donnent en sortie l’image filtrée IM2, avec la même classe de l’image en entrée. Les classes
acceptées sont toutes classes acceptées pour une image valable (logical, uint8,uint16,...).

3 Ouverture et Fermeture
Les opérateurs ouverture et fermeture sont une combinaison des opérateurs élémentaires qu’on a vu au chapitre
précédent appliqués en série en respectant un certain ordre. La définition de ouverture (γB (I)) et fermeture (φB (I)),
par un élément structurant B, est donnée par la combinaison de dilatation et érosion, comme suit:

2
Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

γB (I) = δB (ǫB (I)) = (I  B)  B = I  B, (1.5)

φB (I) = ǫB (δB (I)) = (I  B)  B = I  B. (1.6)

Pratiquement, ces opérateurs appliquent un filtrage externe ou interne à l’objet. Plus l’élément structurant B est
grand, plus la géométrie de l’objet (en termes d’angles et formes) sera uniformée. Un objet qui présente une géométrie
régulière sera peu filtré par rapport aux autres. Plus B est grand, plus les objets clairs ou sombres seront filtrés sur
la base de leur tonalité et dimension. Avec ces types d’opérateurs, la forme de B (le type de voisinage à considérer)
est aussi très influente. Par la suite nous verrons des exemples. Ces deux opérateurs exploitent la combinaison de
l’érosion et de la dilatation (dans un ordre défini), pour filtrer (et non éliminer comme on a vu auparavant) des parties
sombres ou claires de l’image.

Ouvertures (openings) Image Originale Fermetures (Closings)

Figure 1.4: Ouvertures (vers gauche) et fermetures (vers droite) de l’image originale (centre) en utilisant un élément
structurant B circulaire de rayon [2, 4, 6, 8] pixels.

Ces deux opérations sont présentes en tant que fonctions indépendantes en MATLAB (I.P. Toolbox), mais lors de
ce TP on ne les utilisera pas. À titre informatif, les fonctions sont:

IM2 = imopen(IM,SE): IM image, SE él. struct. B, IM2 image filtrée

IM2 = imopen(IM,NHOOD): IM image, NHOOD voisinage de B, IM2 image filtrée.
IM2 = imclose(IM, SE || NHOOD) fonctionne comme la précédente.

4 Ouverture et Fermeture Géodésique (Opening et Closing “by Recon-

struction”)
Les filtres mentionnés auparavant ont une faiblesse remarquable: ils ne préservent pas la forme des objets présents
dans l’image qu’on traite avec la croissance de B. Dans les tâches de filtrage il est souvent utile d’avoir les contours
des formes qui sont l’objet du filtrage. Différemment des opérateurs présentés, qui distribuent ou effacent des valeurs
précises, les filtres géodésiques produisent une image filtrée avec les critères présentés mais contraignant l’opération
dans les limites de l’objet lui-même et non B. Ces filtres reproduisent itérativement l’image originale sur la base d’une
image filtrée (IM ; Marker Image), dilatée ou érodée au départ, et ses subséquents filtrages.
Ouverture géodésique ρδ (IM ) et fermeture géodésique ρǫ (IM ) sont définis, pour une image I et son “marker” IM :

k−1
ρδ (ǫB (I)) k
= ρδ (IM ) = min xIM , δB k

(IM )| δB (IM ) = δB (IM ) , (1.7)
k−1
ρǫ (δB (I)) = ρǫ (IM ) = max xIM , ǫkB (IM )| ǫkB (IM ) = ǫB

(IM ) . (1.8)

En deux mots, l’ouverture géodésique d’une image est donnée par la reconstruction d’une version dilatée de l’image
(marqueur IM ) en assignant, pour chaque pixel sur lequel l’élément structurant est centré (x), le minimum (équivalent
k
à l’infimum) entre xIM (la valeur du pixel-centre de l’image originale) et δB (IM ), la valeur des dilatations pour des
pixels contenus dans B à l’itération k. La fermeture géodésique, de façon opposée, assigne le maximum entre la valeur
originale et les subséquentes itérations d’érosions (donc le supremum dans B) de l’image à l’iteration k. Ces filtres
sont aussi connus comme filtres de reconstruction parce-que ce processus se fait jusqu’à quand l’image filtrée produite
à l’itération k est égale a celle produite à l’itération k − 1.
Visuellement, ces opérateurs filtrent les objets claires et sombres selon la grandeur de l’élément structurant qu’on
utilise au départ (qui nous définit la grandeur des objets à filtrer). Dans l’itération de reconstruction, on utilise
toujours une connectivité B d’ordre 8 (Figure 1.5). Ce type de connectivité est standard lors de l’application des
filtres géodésiques, ainsi il est possible de reconstruire n’importe quel type de forme.
Pour que soit plus clair ce concept de reconstruction itérative, l’algorithme est donné (voir Algorithme 1).
Voyons maintenant des exemples. Ce type de filtrage donne, pour une image en noir et blanc, l’élimination des
objets plus petits de l’élément structurant utilisé au départ (pour créer l’ IM ) en gardant la forme des autres objets

3
Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

_@@ O ? 1 1 1
@@ ~~~
@@ ~
@ ~~~
o x /
~~ @@@ 1 1 1
~ @@
~~~ @@
~  
1 1 1

Figure 1.5: Connectivité de type 8: visuelle (gauche) et matricielle (droite).

Algorithme 1 Calculer ρδ (IM ) et ρǫ (IM )

1: Input: Image I 1: Input: Image I
2: IM = ǫB (I) 2: IM = δB (I)
k k−1 k−1
3: while δB (IM ) 6= δB (IM ) do 3: while ǫkB (IM ) 6= ǫB (IM ) do
4: G = δB (IM ) 4: G = ǫB (IM )
5: IF = min(G, IM ) 5: IF = max(G, IM )
6: if IF 6= IM then 6: if IF 6= IM then
7: IM ← IF 7: IM ← IF
8: end if 8: end if
9: end while 9: end while

composant l’image. La “puissance” du lissage est donc contrôlée par la grandeur de l’élément structurant utilisé pour
la première érosion ou dilatation. Comme dit auparavant, l’opérateur dans la boucle d’itération, prend en paramètre
un élément structurant B de connectivité 8.

Ouverture Géodésique Image Originale Fermeture Géodésique

Figure 1.6: Ouvertures (vers gauche) et fermetures (vers droite) géodésiques de l’image originale (centre) en utilisant
un élément structurant B circulaire de rayon [2, 4, 6, 10] pixels.

Dans l’exemple (Figure 1.6) on a un exemple de filtrage géodésique en noir et blanc. La fermeture efface les
éléments qui ne contiennent pas entièrement B. Au contraire, l’ouverture remplie les zones blanches (assigne la valeur
1) qui sont contenues par un élément structurant, en ne modifiant aucune autre forme.
Pour une image à échelle de gris, ces opérateurs lissent les éléments qui composent le pixel selon le critère choisi.
L’exemple à la Figure 1.7 montre comment les éléments plus ou moins uniformes de l’image (éléments du fond, visage,
main) sont filtrés par rapport à leur dimension et à leur tonalité de gris.

Ouvertures géodésiques Image Originale Fermetures géodésiques

Figure 1.7: Ouvertures (vers gauche) et fermetures (vers droite) géodésiques de l’image originale (centre) en utilisant
un élément structurant B circulaire de rayon [2, 4, 6, 15] pixels.

4
Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

5 Morphologie Mathématique dans le Domaine de la Télédétection

Selon la résolution de l’image (taille réelle du pixel au sol) la variance dans un objet et celle entre les objets peut
être plus ou moins élevée. Normalement, lors de l’analyse qualitative et quantitative d’une image, la variance est de
préférence réduite de façon d’avoir des objets homogènes dans sa surface et une grande variance entre différent objets.
Ces opérateurs ont la capacité de lisser et de rendre uniforme plusieurs objets dans l’image et donc d’en changer les
statistiques locales et globales selon différents besoins. Du point de vue pratique, ces filtres ont été crées pour traiter
des images à échelle des gris, mais souvent en travaillant avec des images satellites, on dispose de plusieurs bandes
spectrales. Or, l’assomption de base, pour des images ayant 3 ou 4 bandes (photo aériennes, QuickBird, Spot,...)
le spectre est très contigu (visible et proche infrarouge), et, malgré des différentes distributions de valeurs selon la
réflectance relative à chaque objet, les bandes nous montrent les mêmes objets pour chaque élément t (sous image à
échelle de gris) de cette matrice MxNxT. Ce discours ne peut pas être fait lors de l’analyse d’images hyperspectrales
/ ultraspectrales captées avec un senseur qui enregistre plusieurs portions du spectre. Dans ces cas, il y a plusieurs
aspects qui doivent être considérés pour garder l’intégrité sémantique de l’image (infrarouge thermique, microondes)
qui peuvent fausser un ultérieur traitement (e.g. classification, régression pour estimation de paramètres biologiques,
NDVI,...).
Les filtrages de base érosion et dilatation (Figure 1.8), agissent sur les éléments claires ou sombres.

Érosions Image Originale Dilatations

Figure 1.8: Érosions (vers gauche) et dilatations (vers droite) de l’image originale (centre) en utilisant un élément
structurant B circulaire de rayon [2, 6, 10, 14] pixels.

Les ouvertures et les fermetures, par contre nous donnent des images filtrées qui gardent les structures principales
avec un niveau de lissage proportionnel soit à la taille de B soit aux valeurs de gris plus ou moins grands pour chaque
bande.
Dans le cas des opérateurs ouverture et fermeture, la forme et la taille de l’élément structurant jouent un rôle
important sur l’effet du filtrage. Jusqu’à maintenant, nous avons utilisé qu’un disque de taille variable pour filtrer les
images. Mais que se passe-t-il lors de l’utilisation de différentes formes?
Avec la fonction SE = strel(type, parameters) on peut initialiser B avec différentes formes. Ici se trouve une
courte liste du paramètre type qui caractérise la forme de l’élément structurant:

‘arbitrary’ → défini par l’utilisateur dans un ma- ‘disk’ → disque de rayon R

trice NHOOD
‘rectangle’ → rectangle de taille [M N]
‘pair’ → 2 cellules décalés, [X Y] indique le saut
‘line’ → ligne de longueur R et de direction DEG
‘diamond’ → diamant de rayon R ‘square’ → carré de taille R
‘periodicline’ → ligne tratillé (périodique) ‘octagon’ → octogone de rayon R

Avec des éléments structurants B de forme différente, le filtrage se fait avec des critères qui suivent la forme de
l’élément. Dans la Figure 1.9 des exemples sont montrés. Dans 1.9(a), le filtrage est cohérent avec les angles d’un
diamant (forme rhomboı̈dale), en ressortant les structures qui sont reliées avec les pistes de l’aéroport. En 1.9(b),
un objet linéaire vertical, nous ressorte de l’information par rapport aux routes qui suivent cette direction. Dans le
dernier cas (1.9(c)), un élément B nous préserve la forme de champs circulaires, mais fourni de l’information “lissée”
par rapport à l’image originale.
Si les opérateurs d’ouverture et de fermeture sont beaucoup utilisés, les plus répandus en absolus sont les opérateurs
géodésiques. Le fait de préserver formes et contours en fait des opérateurs de filtrage très puissants. À titre d’exemple,
ces opérateurs peuvent effacer des objets contenus dans l’image considérée, mais en ne modifiant pas la structure des
éléments plus grands. Un exemple est montré dans la Figure 1.6 et 1.10.

5
Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

(a)

(b)

Ouverture Image Originale Fermeture

(c)

Figure 1.9: Ouvertures (vers gauche) et fermetures (vers droite) de l’image originale (centre) en utilisant un élément
structurant B circulaire de rayon [2, 6, 10, 15] pixels. Ici, les formes des éléments structurants utilisés sont, dans l’ordre,
le diamant (a), la ligne verticale (b) et le cercle (c).

Ouvertures géodésiques Image Originale Fermetures géodésiques

Figure 1.10: Ouvertures (vers gauche) et fermetures (vers droite) géodésiques de l’image originale (centre) en utilisant
un élément structurant B circulaire de rayon [2, 4, 6, 10] pixels.

6 TP: Morphologie Mathématique

Dans cette partie, on va se familiariser avec les filtrages présentés auparavant. Dans cette fiche, uniquement les points
principaux (et les fonctions à utiliser) sont marqués. A vous de créer les fonctions de filtrage!
Suivez aussi les instructions dans le fichier nommé Cours matlab MM.m.

1. Lisez et visualisez l’image

2. Maintenant, deux possibilités s’offrent à vous pour appliquer les opérateurs morphologiques: l’application sur
une version à échelle de gris de l’image et sur une version ‘couleur’.
(a) VARIANTE 1: faite une ACP (Analyse par Composantes Principales) sur l’image, de façon à extraire (avec
la première composante) l’information structurale et spatiale. Cette ‘image’ ainsi obtenue doit être scalée
entre 0 et 1. Utilisez la fonction fournie ISTAN = standardimg(I,type) avec l’image à standardiser comme
paramètre I et 2 comme paramètre type (type de normalisation, ici un “stretch” dans [0; 1]).
(b) VARIANTE 2: filtrez chaque bande de l’image. D’abord standardisez-les (standardimg.m). Attention! Les
fonctions qu’on va construire prennent une bande à la fois. Pour traiter l’image entière, il faut donner en
paramètre une bande à la fois:

imfilter = zeros(size(image,1),size(image,2),size(image,3)); % pre−allocation

% matrice pour vitesse
for i = 1:no bandes
imfilter(:,:,i) = filtres morphologiques(image(:,:,i),SE,p switch);
end

6
 Cours MATLAB - Analyse morphologique d’image - UNIL-FGSE H 09/10

3. Créez des fonctions de filtrage (voir fichier Cours matlab MM.m). N’oubliez pas d’initialiser d’abord les éléments
structurants B soit par la fonction SE = strel(type,par) soit en le déclarant soi-même comme une matrice
binaire de voisinage.
4. Visualisez et sauvegardez les images filtrés. Après avoir vérifié le fonctionnement des trois fonctions, regroupez
les dans un seule fonction (nommez-la filtres morphologiques). Ainsi, vous avez la fonction complète.
5. Testez différents types de SE et observez les changements, aussi par rapport aux différents types d’opérateurs.
OPT.1 Lancez la fonction exemple ocr avec une image à échelle de gris, définissez un SE et le paramètre 'op recon'
ou 'cl recon' pour choisir le type de paramètre. Il s’agit d’une démonstration itérative de comment le filtres
de reconstruction marchent.
OPT.2 Ouvrez GoogleEarth et cherchez des images avec des structures particulières à traiter. Une fois qu’on a sauve-
gardé l’image en .jpg et on l’ouvre avec MATLAB, on a à disposition 3 bandes (RGB), comme si on travaillait
avec des photos aériennes. Recoupez l’image en MATLAB pour enlever la légende et des graphiques qui ne sont
pas parties de la photo aérienne. Filtrez-la avec des opérateurs et des éléments structurants de votre choix.

7
Références

[1] M. Fauvel, J. A. Benediktsson, J. Chanussot, and J. R. Sveinsson. Spectral and spatial classification of hyperspectral
data using svms and morphological profiles. IEEE Transactions on Geoscience and Remote Sensing, 46(11):3804
– 3814, 2008.
[2] R. C. Gonzalez, Richard E. Woods, and S. L. Eddins. Digital Image Processing Using MATLAB. Pearson Prentice
Hall, 2004.
[3] Kumar Navulur. Multispectral Image Analysis Using the Object-Oriented Paradigm. CRC Press, 2007.
[4] Maria Petrou and Pedro Garcia Sevilla. Dealing With Texture. John Wiley and Sons, 2006.

[5] J. Plaza, A. Plaza, and C. Barra. Multi-channel morphological profiles for classification of hyperspectral images
using support vector machines. Sensors, 9:196 – 218, 2009.
[6] D. Tuia, F. Pacifici, M. Kanevski, and W. J. Emery. Classification of very high spatial resolution imagery using
mathematical morphology and support vector machines. IEEE Transactions on Geoscience and Remote Sensing,
47(11):3866 – 3879, 2009.