Coursimage09 PDF

Quelques méthodes mathématiques pour le traitement
d’image
2 janvier 2009
Ce cours est une introduction à la théorie mathématique de traitement de l’image . Il est donc
incomplet car les méthodes dans ce domaine sont nombreuses et variées. On se polarisera sur
les méthodes variationnelles.
M. Bergounioux Master 2 - 2008-2009

2
Chapitre 1
Introduction
1.1 Qu’est-ce qu’une image numérique ?

Une image numérique est composée d’unités élémentaires (appelées pixels) qui représentent
chacun une portion de l’image. Une image est définie par :
– le nombre de pixels qui la compose en largeur et en hauteur (qui peut varier presque à
l’infini),
– l’étendue des teintes de gris ou des couleurs que peut prendre chaque pixel (on parle de
dynamique de l’image).
1. Les images binaires (noir ou blanc)
Exemple, images les plus simples, un pixel peut prendre uniquement les valeurs noir ou
blanc. C’est typiquement le type d’image que l’on utilise pour scanner du texte quand
celui ci est composé d’une seule couleur.
2. Les images en teintes de gris
En général, les images en niveaux de gris renferment 256 teintes de gris. Image à 256
couleurs, simplement chacune de ces 256 couleurs est définie dans la gamme des gris.
Par convention la valeur zéro représente le noir (intensité lumineuse nulle) et la valeur
255 le blanc (intensité lumineuse maximale).
3. Les images couleurs
S’il existe plusieurs modes de représentation de la couleur, le plus utilisé pour le manie-
ment des images numériques est l’espace couleur Rouge, Vert, Bleu (R,V,B).
Cet espace couleur est basé sur la synthèse additive des couleurs, c’est à dire que le
mélange des trois composantes (R, V, B) donne une couleur.
3
4 CHAPITRE 1. INTRODUCTION
1.1.1 Quelques définitions

Pixels et niveaux de gris
Echantillonnage et quantification
1.1. QU’EST-CE QU’UNE IMAGE NUMÉRIQUE ? 5
L’échantillonnage est le procédé de discrétisation spatiale d’une image consistant à asso-

cier à chaque zone rectangulaire R(x, y) d’une image continue une unique valeur I(x, y). On
parle de sous-échantillonnage lorsque l’image est déjà discrétisée et qu’on diminue le nombre
d’échantillons.
Une image numérique est une image échantillonnée et quantifiée. La quantification désigne
la limitation du nombre de valeurs différentes que peut prendre I(x, y).
L’échantillonnage est une étape fondamentale qui doit tenir compte du contenu informa-
tionnel pertinent de l’image à analyser. Sur l’exemple ci- dessous, en 1d, le signal échantillonné
« ressemble » à une sinusoı̈de de fréquence 8 fois plus faible :
Ce phénomène appelé aliasing est encore pire en 2D, car il affecte la fréquence et la direction
des structures périodiques. Imaginons par exemple qu’on souhaite échantillonner l’image cor-
respondant aux bandes noires ci-dessous :
Avec un échantillonnage adapté, l’image numérique fait apparaı̂tre des structures conformes à
l’information présente dans l’image :
Mais en considérant seulement 1 échantillon sur 2, une structure différente apparaı̂t, dont l’ana-
lyse (ici des bandes verticales, plus épaisses) ne sera pas conforme à la réalité de l’objet :
Image originale Image

sous-échantillonnée
1.1. QU’EST-CE QU’UNE IMAGE NUMÉRIQUE ? 7
La quantification peut également faire apparaı̂tre des distortions dans les images. Comme pour
l’échantillonnage, il existe des règles pour déterminer la bonne quantification (le bon nombre
de bits) pour coder les images numériques. L’une dépend du capteur, et de sa capacité effective
à observer des signaux de valeurs différentes : le rapport signal sur bruit.
Le rapport signal sur bruit est défini à partir du rapport entre l’amplitude des niveaux de
gris mesurables par le capteur nmax − nmin et le niveau du bruit, en gros l’écart-type σn de la
perturbation aléatoire qui affecte les niveaux de gris. En prenant le logarithme, on a le nombre
de bits utile au capteur pour coder les images.
Outre les capacités du capteur, le nombre de bits réellement nécessaires pour coder une
image varie d’une image à l’autre, en fonction de leur contenu informationnel. Ce nombre
dépend de l’entropie, définie à partir de la distribution des niveaux de gris de l’image (statis-
tique). !
E=− pi log2 (pi ) ,
i≤N
où N est le nombre de niveaux de gris présents, pi est la proportion (0 < pi < 1) de points de
l’image ayant pour niveau de gris i. Cette grandeur représente le nombre moyen de bits par
pixel nécessaires pour coder toute l’information présente. Elle est utilisée dans les techniques
de compression sans perte pour adapter le volume de donnée des images à leur contenu infor-
mationnel.
Profil - Histogramme
Il est possible de tracer un trait sur le flanc du zèbre de l’image en niveaux de gris ci-dessous
et obtenir le profil correspondant, c’est à dire le niveau de gris de chaque point ou pixel traversé
par la ligne :
L’histogramme de l’image en niveau de gris ci-dessous permet d’établir une stricte corrélation
entre les données numériques codant la nuance de gris et la position des pixels de l’image.
1.2 Qu’est-ce que le traitement d’image ?

Les pages qui suivent sont extraites des cours en ligne [BMT, M].
1.2. QU’EST-CE QUE LE TRAITEMENT D’IMAGE ? 9
1.2.1 Quelques aspects du Traitement d’Image
– Filtrage / déconvolution (ou filtrage inverse)

– Compression
– Segmentation
– Restauration / reconnaissance
– Reconstruction tomographique
Dans ce cours nous évoquerons successivement trois aspects fondamentaux :

– Filtrage
L’outil mathématique essentiel pour le filtrage est la transformation de Fourier (ou toute
autre transformation du même type comme la transformation en ondelettes)
– Segmentation
La segmentation fait intervenir des notions d’optimisation, des outils géométriques et des
équations aux dérivées partielles
– Restauration
Les modèles variationnels en restauration utilisent de l’optimisation et de l’analyse fonc-
tionnelle fine (Banach, théorie de la mesure).
1.2.2 Applications
– Robotique - Industrie
– Assemblage, reconnaissance de pièces
– Contrôle de qualité
– Véhicule autonome
– etc ...
– Télédétection
– Météo
– Cartographie
– Analyse des ressources terrestres
– Astronomie
– Restauration
– etc...
– Applications militaires
– Guidage de missile
– Reconnaissance (aérienne, sous-marine, etc ...)
– etc ...
– Imagerie médicale
– Tomographie
– Aide au diagnostic
– Comptage (nombre de cellules)
– Suivi de formes anatomiques
– Restauration
– etc...
– Sécurité
– Reconnaissance (d’empreintes, visages, signatures)
– Détection de mouvement
– etc...
Chapitre 2
Traitement ponctuel des images

numériques
On s’intéresse d’abord aux traitements ponctuels qui consistent à faire subir à chaque pixel
une correction ne dépendant que de sa valeur. On trouve dans cette catégorie, les fonctions de
recadrage ou d’égalisation de dynamique, de binarisation ...
Sauf mention particulière, nous supposerons dans ce qui suit des images comportant N 2
pixels codés sur 256 niveaux de gris différents.
2.1 Correction ponctuelle d’une image -Recadrage de dynamique
Il s’agit d’une transformation du type f " = t(f ) qui permet de modifier la dynamique des
niveaux de gris dans le but d’améliorer l’aspect visuel de l’image. À un niveau de gris f de
l’image originale correspond le niveau t(f ) dans l’image transformée. On fait subir à chaque
pixel un traitement ne dépendant que de sa valeur. La transformation t(f ) peut être réalisée
en temps réel sur l’image en cours d’acquisition à l’aide d’une table de transcodage dans la-
quelle les valeurs de la transformation sont mémorisées. Un adressage de cette mémoire par
une donnée f fournit directement la valeur t(f ).
2.1.1 Transformation de recadrage
On suppose une image de départ présentant un histogramme concentré dans l’intervalle

[a, b]. Les valeurs a, b correspondent aux niveaux de gris extrêmes présents dans cette image. Le
recadrage de dynamique consiste à étendre la dynamique de l’image transformée à l’étendue
totale [0, 255]. La transformation de recadrage est donc une application affine qui s’écrit :
15
16 CHAPITRE 2. TRAITEMENT PONCTUEL DES IMAGES NUMÉRIQUES
Original Recadrage : a = 30, b = 200
Variantes pour le rehaussement des contrastes

Les types de correction donnés ci-dessous permettent d’accentuer le contraste dans une
plage précise de niveau.
Dilatation de la dynamique des zones sombres Dilatation de la dynamique des zones claires
Fonction de rehaussement de contraste.
2.1. CORRECTION PONCTUELLE D’UNE IMAGE -RECADRAGE DE DYNAMIQUE 17


 b
 f pour 0 ≤ f ≤ a
a
t(f ) =
 (255 − b) f + 255(b − a)

 pour a ≤ f ≤ 255
255 − a
Original Histogramme
Dilatation de la dynamique des zones claires
Dilatation de la dynamique des zones sombres

2.1.2 Égalisation de l’histogramme
L’histogramme d’une image est rarement plat ce qui traduit une entropie non maximale. La
transformation d’égalisation est construite de telle façon que l’histogramme de l’image trans-
formée soit le plus plat possible. Cette technique améliore le contraste et permet d’augmenter
artificiellement la clarté d’une image grâce à une meilleure répartition des intensités relatives.
Fonction d’aplatissement continue
Considérons l’histogramme continu h(f ) donné ci-dessous. En notant f " = t(f ), l’histo-
gramme égalisé h(f " ) doit s’approcher de la forme idéale décrite ci-dessous.
Deux surfaces élémentaires en correspondance dans les histogrammes initiaux et égalisés,
présentent le même nombre de points ce qui permet d’écrire :
&
" 256 f
f = t(f ) = 2 h(s) ds .
N 0
Histogramme d’origine Histogramme plat idéal

Fonction idéale d’égalisation d’un histogramme
Fonction d’aplatissement discrète

En remplaçant l’intégration continue par une somme, on obtient la transformation d’égalisation
discrète suivante :
f
256 !
f " = t(f ) = h(i) .
N2
i=0
2.1. CORRECTION PONCTUELLE D’UNE IMAGE -RECADRAGE DE DYNAMIQUE 19
Image égalisée Histogramme égalisé
2.1.3 Binarisation
Le but de la binarisation d’une image est d’affecter un niveau uniforme au pixels pertinents
et d’éliminer les autres.
Seuillage
Le seuillage consiste à affecter le niveau 255 aux pixels dont la valeur est supérieure à un
seuil S et 0 le niveau aux autres. Le graphe de la transformation correspondante est le suivant
Fonction « seuillage »
Seuillage à 70
Extraction d’une fenêtre d’intensité

Avec la transformation décrite ci-dessous, la nouvelle image ne visualise que les pixels dont
le e niveau d’intensité appartient à l’intervalle [a, b]. Sous réserve d’une connaissance a priori
de la distribution des niveaux de gris des objets de l’image originale, cette technique permet
une segmentation d’objets particuliers de l’image.
Fonction « fenêtre d’intensité »
2.2. ANAYSE DE LA NETTETÉ D’UNE IMAGE NUMÉRIQUE 21
Seuillage avec fenêtre d’intensité entre 30 et 100
2.2 Anayse de la netteté d’une image numérique
La focalisation automatique, ou mise au point, des systèmes de prise de vue traditionnels

(appareils photo,, graphiques caméras analogiques ... ) exploite généralement un télémètre qui
mesure avec précision la distance appareil-plan objet. Les systèmes de vision numériques ac-
tuels (appareils et caméras numériques ... ) réalisent automatiquement leur mise au point à par-
tir d’une analyse de l’image restituée. Pour cela un critère de netteté de l’image est déterminé
à partir des valeurs f (i, .j) := fij de l’intensité (luminance) des pixels.
Un grand nombre de critères a été proposé. Parmi les plus utilisés figurent les critères basés
sur l’analyse de la distribution des niveaux d’intensité pixel et ceux mesurant le contenu spec-
tral de l’image.
La première catégorie repose sur le fait que la défocalisation engendre l’uniformisation des
niveaux de gris, ou ce qui revient au même, que l’image nette présente l’histogramme le plus
large.
La seconde catégorie repose sur le principe similaire qu’une image focalisée contient plus
de fréquences spatiales élevées qu’une image floue.
Dans toutes ces méthodes, on recherche la distance de mise au point qui assure la majora-
tion du critère
2.2.1 Exemples de critères de netteté
Critères de netteté basés sur l’analyse de l’histogramme de l’image
On note hk la distribution statistique (histogramme) des niveaux k d’intensité des pixels.

Critère mathématique
!
Méthode de F = k hk
k>T
Mendelsohn T est un paramètre choisi
et Mayall proche du niveau de gris moyen de l’image
' '
j fij ∆ij
! i
Méthode de F = (k − T ) hk et T = ' '
i j ∆ij
(
k>T )
Mason et Green ∆ij = 2 (fi,j+1 − fi,j−1 )2 + (fi+1,j − fi−1,j )2
+ (fi−1,j−1 − fi+1,j+1 )2 + (fi−1,j+1 − fi+1,j−1 )2
(' ' )2
' ' 2 fij
i j fij i j
Variance des F = −
N2 N4
niveaux de gris
Critères de netteté basés sur l’analyse spectrale de l’image
Critère mathématique
!
Norme "1 du gradient1 F = (|fi,j − fi,j−1 | + |fi,j+1 − fi,j |)
i,j
!
Norme "2 du gradient F = (fi+1,j − fi,j )2 + (fi,j+1 − fi,j )2
i,j
Remarques générales
– Les tests des critères de netteté réalisés sur un grand nombre d’images différentes montrent
que leurs performances dépendent du type et du contenu de l’image analysée ;
– En règle générale, l’extremum des critères est peu prononcé lorsque l’image contient peu
de détails
2.3. CORRÉLATION D’IMAGES NUMÉRIQUES 23
– L’algorithme de Brenner est souvent utilisé car il présente généralement une bonne sen-
sibilité et la charge de calculs qu’il nécessite est raisonnable.
2.3 Corrélation d’images numériques
Nous considérons deux fonctions bidimensionnelles discrètes f (i, j) et g(i, j) avec 1 ≤ i, j ≤

N.
Ces fonctions sont représentatives de deux images numériques monochromes comportant
N 2 pixels. Leurs versions centrées et réduites sont définies respectivement par
f (i, j) − µf g(i, j) − µg
fc (i, j) = et gc (i, j) =
σf σg
où µ et σ sont respectivement la moyenne et l’écart type de chaque fonction.

La fonction d’intercorrélation entre fc et gc est définie par la relation :
1 !
ϕf g (k, l) = fc (i, j)gc (i − k, j − l) . (2.3.1)
N2
i,j
Il est souvent plus pratique de travailler sur des images dont les valeurs sont centrées et
réduites. Ceci permet d’une part d’éliminer de ϕf g (k, l) le carré des moyennes des images,
d’autre part de normaliser les fonctions de corrélation.
Pour éviter les effets de bord qui introduisent un biais dans le calcul de ϕf g (k, l) il convient :
– soit de s’assurer que la double sommation est réalisée sur des valeurs dont les coor-
données ne sortentjainais de l’image. Cela revient à définir un cadre d’intégration dont
le format, inférieur à celui de l’image, est choisi en fonction du domaine de calcul de
ϕf g (k, l) ;
– soit de remplacer dans l’expression (2.3.1) le diviseur N 2 par la valeur variable M qui
dépend de k et l selon la relation : M = (N − |k|)(N − |l|).
2.3.1 Application à la reconnaissance d’empreintes digitales
Les fonctions d’intercorrélation bidimensionnelles peuvent être utilisées en reconnaissance

d’images. La comparaison de l’empreinte digitale d’un suspect avec celles contenues dans un
fichier en est un exemple. les deux images représentées ci-dessous représentent les empreintes
digitales de deux individus que nous appellerons F et G.
Les caractéristiques générales de ces deux images sont données dans le tableau suivant :
Individu F Individu G
Format 256 x 256 256 x 256
Codage en niveaux de gris 0 à 255 0 à 255
Moyenne µ 132,8 127
Écart type σ 33,8 36,4
Les valeurs des deux images sont centrées et réduites. Pour éviter les effets de bords, l’in-
tercorrélation est réalisée entre deux zones carrées de 151 x 151 pixels (voir figure ci-dessous)
La fonction est estimée dans la gamme −20 ≤ k, l ≤ 20 grâce à la formule suivante
200
1 !
ϕf g (k, l) = fc (i, j)gc (i − k, j − l) . (2.3.2)
1512
i,j=50
La représentation 3D à gauche dans la figure ci-dessous est celle de la fonction d’autocorrélation

ϕf g (k, l) effectuée sur l’empreinte F . Son maximum estégal à 1 au centre du graphique qui cor-
respond ici à la position k = l = 0. La représentation de droite est la fonction d’intercorrélation
entre les empreintes F et G.
2.3. CORRÉLATION D’IMAGES NUMÉRIQUES 25
Ses valeurs ne dépassent pas 0.15 ce qui prouve que les deux empreintes ne sont pas iden-
tiques.
Dans la réalité, les empreintes ne sont pas toutes obtenues dans les mêmes conditions de
positionnement. Il est nécessaire d’ajouter un paramètre de rotation d’image à la fonction d’in-
tercorrélation ce qui peut alourdir considérablement les calculs.
2.3.2 Application à l’analyse de la texture d’une image
Par définition la texture d’une image est la structure spatiale sur laquelle sont organisés les
pixels. Les relations structurelles peuvent être :
– déterministes : c’est le cas de la répétition quasi périodique d’un motif de base (brique,
carrelage, mailles de tissu ... ) ;
– aléatoires : il n’y a pas de motif de base (sable, mur crépi...)
La fonction d’autocorrélation bidimensionnelle est bien adaptée pour mettre en évidence les
propriétés de la texture d’une image. La figure suivante représente deux exemples de texture.
Nous donnons ci-dessous la fonction d’auto-corrélation des textures « sable » et « tissu ».

2.4 Transformation de Hough

[,a transformation de Hough est utilisée pour détecter de manière systématique la présence
de relations structurelles spécifiques entre des pixels dans une image Par exemple une image
représentant un site urbain est composée de nombreuses lignes droites (immeubles, fenêtres)
en revanche, une vue de campagne en est quasiment dépourvue .
Hough propose une méthode de détection basée sur une transformation d’image permet-
tant la reconnaissance de structures simples (droite, cercle, ... ) liant des pixels entre eux. Pour
limiter la charge de calcul, l’image originale est préalablement limitée aux contours des objets
puis binarisée (2 niveaux possibles pour coder l’intensité pixel).
2.4.1 Principe de la méthode pour la recherche de ligne droite

Supposons que l’on suspecte la présence d’une droite ∆ reliant un certain nombre de pixesl
Pi . Soit le pixel P1 de coordonnées (x1 , y1 ). Une infinité de droites d’équation : y1 = ax1 + b
peuvent passer par P1 . Cependant, dans le plan des paramètres ab l’équation qui s’écrit b =
−ax1 + y1 devient une droite unique D1 (voir figures ci-dessous).
Un second pixel P2 = (x2 , y2 ) permet de définir une seconde droite D2 du type b = −ax2 +y2
dans le plan ab.L’intersection de D2 avec D1 fournit le couple " a" , b" ) qui sont les paramètres
2.4. TRANSFORMATION DE HOUGH 27
de la droite recherchée ∆ dans le plan image. Ainsi, tous les pixels Pi qui sont alignés sur ∆
possède une droite Di dnas le plan ab qui coupe les autres au point particulier (a" , b" ).
2.4.2 Problème pour la recherche de ligne verticale Représentation normale d’une

droite
Dans le plan image, une droite verticale possède des paramètres a etb infinis, ce qui ne
permet pas d’exploiter la méthode précédente. Pour contourner ce problème, on utilise la
représentation normale des droites. Cette représentation, de paramètres ρ et θ, obéit à l’équation
x cos θ + y sin θ = ρ . (2.4.3)

* +
& = cos θ
Cette équation représente le produit scalaire (projection) entre les vecteurs V et
sin θ
* +
−−→ x
OM = . Ainsi l’ensemble des points M d’une droite se projettent sur un même vecteur
y
particulier de coordonnées polaires (ρ, θ) . Les cas de droites horizontales et verticales sont
illustrés ci-dessous :
2.4.3 Transformation de Hough pour la détection de droites dans une image binaire
Domaine de variation des paramètres θ et ρ
Il est à noter que si (ρ, θ) sont les paramètres d’une droite (−ρ, θ + π) le sont également. Par
conséquent l’intervalle [0, π[correspond au domaine de variation complet du paramètre θ. Les
coordonnées cartésiennes x et y des pixels d’une image numérique sont généralement positives,
l’origine étant placée à un sommet de l’image. En considérant une image carrée comportant √
N × N pixels, les valeurs de ρ calculées par l’équation (2.4.3) peuvent être majorées par 2N .
Quadrillage du plan θρ
On suspecte dans l’image la présence d’une ou plusieurs structures caractérisées par l’ali-
gnement d’un certain nombre de pixels. Soit les deux intervalles [θmin , θmax ] dans lesquels on
cherche à calculer les paramètres de la représentation normale d’une ou plusieurs droites.
Le plan θρ limité aux intervalles de recherche, est subdivisé en cellules par quantification
des paramètres θ et ρ avec les pas respectifs ∆θ et ∆ρ (indexation par les indices p et q). La
figure ci-dessous donne un exemple de quadrillage du plan θρ. Un tableau A(p, q) dont les
valeurs sont initialement mises à zéro, est associé à ce quadrillage
Procédure de transformation
– Pour chaque pixel non nul Pi (xi , yi ) de l’image, on balaye l’axe des θ de θmin à θmax
suivant la trame du tableau et pour chaque θp on résout l’équation (2.4.3) ;
– le résultat ρ obtenu est arrondi à la valeur ρq du tableau la plus proche ;
– si à la valeur θp correspond la solution ρq la valeur A(p, q) est incrémentée d’une unité.
À la fin de cette procédure, A(p, q) = M signifie que M points de l’image sont alignés sur
la droite de paramètres approximatifs (θp , ρq ). Après transformation complète, le tableau peut
être représenté sous la forme d’une image en niveaux de gris (exemple plus loin) ou celle d’un
graphique 3D. La décision sur la détection de droites peut être prise après recherche des coor-
données des valeurs significatives du tableau.
La charge de calcul nécessaire pour réaliser la transformation de Hough est importante. Elle
dépend du nombre de paramètres recherchés :
- recherche de droites : 2 paramètres ;
- recherche de cercles : 3 paramètres.
Exemple pour la recherche de droite : pour N pixels non nuls de l’image binaire et K sub-
divisions de l’axe θ. il y a N K déterminations de l’équation (2.4.3). Une réduction du temps
de calcul peut être obtenue par l’utilisation de tables préenregistrées des conversions sin θp et
cos θp
Exemple
Nous considérons dans cet exemple une image binaire comportant 6 x 6 pixels dont les
valeurs sont données dans le tableau suivant
2.4. TRANSFORMATION DE HOUGH 29
Deux droites ∆1 et ∆2 comportant chacune 6 pixels alignés apparaissent clans cette image.
La procédure de calcul décrite au paragraphe précédent est utilisée avec les paramètres sui-
vants : θp varie de 0 à 180o par pas de 1o .
,
Plus généralement
, si on considère une image de taille N 1 ×N 2 , ρmin = − (N1 − 1)2 + (N2 − 1)2
et ρmax = (N1 − 1) + (N2 − 1) . Si on discrétise l’espace de Hough avec une résolution de
2 2
hθ (par exemple 1o ) pour θ et une résolution de hρ pour ρ le tableau de référence est donné par
θp = −90 + phθ avec p ∈ {0, Nθ } et − 90 + Nθ hθ = 90 ,
ρq = ρmin + qhρ avec q ∈ {0, Nρ } et ρmin + Nρ hρ = ρmax .
Le programme itératif suivant permet le calcul de la transformation de Hough :
Mi,j sont les valeurs binaires de l’image originale : Mi,j ∈ {0, 1}

Pour p de 1 à Nθ
Pour q de 1 à Nρ
Ap,q = 0
Pour i de 1 à N1
Pour j de 1 à N2
Pour p de 1 -
à Nθ( . p / . p /) 0
q = Arrondi h1ρ i cos π + j sin π − ρmin
180 180
Ap,q = Ap,q + 1 si Mi,j '= 0
Le même type de programme appliqué à l’image binaire (64 × 64) donne le résultat ci-dessous :
Chapitre 3
Filtrage
3.1 Filtrage linéaire des signaux 1D

L’exemple le plus courant de signal 1D est le signal sonore.
Définition 3.1.1 (Filtre) Un filtre est un système linéaire continu et invariant. En d’autres termes,
on se donne deux espaces normés X et Y ainsi qu’un opérateur linéaire continu, invariant par translation
A : X → Y.
Supposons qu’ on puisse définir la transformée de Fourier des signaux d’entrée f et de

sortie g (par exemple si ces fonctions sont dans L1 (R) ou dans L2 (R) et que le signal de sortie
le filtrage de f par un filtre linéaire est tel que
ĝ = H fˆ
où ĝ désigne la transformée de Fourier de g. La fonction H est appelée fonction de transfert

du filtre. Si la fonction de transfert H est dans L2 (R) ∩ L∞ (R), elle admet une transformée de
Fourier inverse h = F −1 H, bornée et à décroissance rapide, continue (sauf peut-être à l’origine)
et la relation
ĝ = ĥfˆ
entraı̂ne
g =h∗f .
La réponse à l’entrée f est donc un produit de convolution de l’entrée avec une fonction fixe
h que l’on appelle réponse impulsionnelle.
Définition 3.1.2 On appelle réponse indicielle d’un filtre sa réponse à l’échelon unité (fonction de
Heaviside) :
& t
h1 (t) = h(s) ds .
−∞
31
32 CHAPITRE 3. FILTRAGE
Notons que la plupart des filtres courants sont des filtres de convolution. Il est courant de
définir un filtre par la façon dont il modifie les fréquences du signal d’entrée, c’est-à-dire par
sa fonction de transfert H(λ) puisque les fréquences de l’entrée f et de la sortie g sont liées par
ĝ(λ) = H(λ)fˆ(λ) .
Un filtre idéal est un filtre qui élimine totalement les bandes de fréquence indésirables sans
transition et sans déphasage dans les bandes conservées.
Selon la bande rejetée, on rencontre 4 grandes catégories de filtres.
La région où les fréquences sont coupées est en gris ci-dessus. Par exemple le passe-bas idéal
est le filtre qui ne modifie par les fréquences λ telles que λ ≤ λc (fréquence de coupure) et
supprime les autres. D’où 1
1 si |λ| ≤ λc
H(λ) =
0 sinon
On reconnaı̂t h ∈ L2 (R) telle que ĥ = H. C’est
sin 2πλc t
h(t) = .
πt
Si on se limite à des signaux d’entrée d’énergie finie, f, h et H sont dans L2 (R) et g = h ∗ f . On
sait que g est continue, bornée et nulle à l’infini. h ∗ f est dans L2 (R) comme fˆ.
Toutefois un tel filtre n’est pas réalisable. On se contente de filtres passe-bas approchés
(voir ci-dessous un dispositif physique de réalisation d’un tel filtre) qui atténuent les hautes
fréquences au lieu de les supprimer (voir exercices) et qui entraı̂nent un déphasage :
3.2. FILTRAGE 2D : CONVOLUTION /MÉDIAN 33
3.2 Filtrage 2D : convolution /médian

Le célèbre format bitmap, qui tire son nom de l’anglais « bitmap » pour « carte de bits »
montre qu’une image est avant tout un domaine spatial sur lequel on peut se promener avec
la souris de l’ordinateur : les distances en pixels dans l’image I sont dès lors liées aux distances
réelles en mètres dans la scène réelle S.
La fréquence spatiale est un concept délicat qui découle du fait que les images appar-
tiennent au domaine spatial. Pour commencer on peut rappeler que la fréquence est une gran-
deur qui caractérise le nombre de phénomènes qui se déroulent au cours d’un temps donné :
en voiture le long d’une route vous voyez 2 bandes blanches PAR seconde : c’est une fréquence
temporelle. Il est ensuite facile de comprendre que ce concept de fréquence « temporelle » peut
aussi se traduire en disant qu’il y a 200 bandes blanches PAR kilomètre : c’est une fréquence
spatiale.
Dans une image, les détails se répètent fréquemment sur un petit nombre de pixels, on dit
qu’ils ont une fréquence élevée : c’est le cas pour les bords et les contours dans une image.
Au contraire, les fréquences basses correspondent à des variations qui se répètent peu car,
diluées sur de grandes parties de l’image, par exemple des variations de fond de ciel.
Nous verrons dans la suite que la plupart des filtres agissent sélectivement sur ces fréquences
pour les sélectionner, en vue de les amplifier ou de les réduire.
3.2.1 Filtrage spatial

Filtres de convolution
Le filtrage spatial est essentiellement une opération de convolution (2D). Si f est l’image à
filtrer (ou à rehausser) et g le filtre spatial (ou PSF ou masque) on a :
 

 

−1
f (x, y) ∗ g(x, y) = F F(f (x, y)) · F(g(x, y))} .

 2 34 5

G(u,v)
G est la fonction de transfert du filtre. On peut distinguer trois types de filtrage :
Une image numérique étant essentiellement discrète (pixels et niveaux de gris) nous allons
présenter les filtres dans le cas discret. Dans tout ce qui suit x et y sont des entiers (coordonnées
des pixels) et f est à valeurs entières (dans {0, · · · , 255}).
On ne fait pas en général une convolution globale mais une transformation basée sur le
voisinage d’un point (x, y) :
Le noyau de convolution du filtre κ est à support compact inclus dans [x1 , x2 ] × [y1 , y2 ] :
y2
x2 !
!
g(x, y) = (f ∗ κ)(x, y) = f (x − i, y − j)κ(i, j) .
i=x1 j=y1
Généralement le filtre est de dimension d impaire et est symétrique. Dans ce cas
[x1 , x2 ] = [y1 , y2 ] = [−d/2, d/2] ,
(d−1)/2 (d−1)/2
! !
(f ∗ κ)(x, y) = f (x + i, y + j)κ(i, j) . (3.2.1)
i=−(d−1)/2 j=−(d−1)/2
Filtre(i,j)
w1 w2 w3 ← y−1
w4 w5 w6 ← y
w7 w8 w9 ← y+1
↑ ↑ ↑
x−1 x x+1
Ici d = 3. On ne filtre pas les bords pour éviter des distorsions ; donc κ(0, 0) = w5 .
Sur cet exemple on a précisément
g(x, y) = w1 f (x − 1, y − 1) + w2 f (x, y − 1) + w3 f (x + 1, y − 1)
+w4 f (x − 1, y) + w5 f (x, y) + w6 f (x + 1, y)
+w7 f (x − 1, y + 1) + w8 f (x, y + 1) + w9 f (x + 1, y + 1) .
Afin
!de conserver la moyenne de l’image f , la somme des éléments du filtre est normalisée à
1: wi = 1.
i
Un filtre 2D est dit séparable s’il est possible de décomposer le noyau de convolution h2D en
deux filtres 1D appliqués successivement en horizontal puis en vertical (ou inversement) :
h2D = hV1D ⊗ hH
1D ,
où le symbole ⊗ désigne le produit de convolution. On peut alors traiter séparément les lignes
et les colonnes de l’image.
Pour qu’un filtre 2D soit séparable il faut et il suffit que les coefficients de ses lignes et de ses
colonnes soient proportionnels.
Exemple 3.2.1 (Filtres séparables )
α aα bα cα
a b c ⊗ β = aβ bβ cβ
γ aγ bγ cγ
Exemple 3.2.2 (Filtre de moyenne passe -bas )

1 1 1 1 1
1 1 1 1 1 1 1 1
1 1
· 1 1 1 · 1 1 1 1 1
9 25
1 1 1 1 1 1 1 1
1 1 1 1 1
Filtre 3 x 3 Filtre 5 x 5
Ce sont des filtres séparables.

Exemple 3.2.3 (Filtre gaussien)
1 2 1
1
· 2 4 2
16
1 2 1
Idéalement, on devrait prévoir un filtre de taille (6σ +1)×(6σ +1). En général un filtre gaussien
avec σ < 1 est utilisé pour réduire le bruit, et si σ > 1 c’est dans le but de fabriquer une image
qu’on va utiliser pour faire un « masque flou » personnalisé. Il faut noter que plus σ est grand,
plus le flou appliqué à l’image sera marqué.
Exemple 3.2.4 (Filtre binômial)
Les coefficients de ce filtre sont obtenus par le binôme de Newton. Un filtre 1D binômial d’ordre 4
1
est un filtre séparable donné par le vecteur v = (1 4 6 4 1). Un filtre 2D binômial d’ordre 4 est donné
16
par v " v :
1 4 6 4 1
4 16 24 16 4
1
· 6 24 36 24 6
256
4 16 24 16 4
1 4 6 4 1
Filtres médians
Ce ne sont pas des filtres de convolution, ni des filtres linéaires.
g(x, y) = médian{f (n, m) | (n, m) ∈ S(x, y) } ,
où S(x, y) est un voisinage de (x, y).

Exemple :
bruit
30 10 20 ↓
10 250 25 → 10 10 20 20 25 25 30 30 250
20 25 30 ↑
médiane
On remplace la valeur du pixel par la valeur médiane ou la valeur moyenne. Ce filtre est utile
pour contrer l’effet « Poivre et Sel » (P& S) c’est-à-dire des faux « 0 » et « 255 » dans l’image.
Image bruitée « Poivre et Sel »

Filtre de moyenne : rayon 3 Filtre de moyenne : rayon 5
Filtre de moyenne : rayon 7
Filtre médian : rayon 3 Filtre médian : rayon 5

Filtre médian : rayon 7
Si le bruit P& S est supérieur à la moitié de la dimension du filtre, le filtrage est inefficace.
Filtres passe-haut
L’image obtenue par un filtre passe-haut correspond en général à ce qui « reste » après un
filtrage passe-bas.
3.2.2 Filtrage fréquentiel

Transformation de Fourier 2D
Avant d’envisager le traitement numérique des signaux il est nécessaire de donner l’in-
terprétation des signaux bidimensionnels dans le domaine des fréquences.
La représentation fréquentielle des signaux 2D est l’extension directe de celle des signaux
monodimensionnels. La transformée F (u, v) de Fourier d’un signal 2D f (x, y) est
&&
F (u, v) = f (x, y)e−2iπ(xu+yv) dx dy . (3.2.2)
R2
La reconstitution du signal spatial se fait par transformation inverse :

&&
f (x, y) = F (u, v)e2iπ(xu+yv) du dv . (3.2.3)
R2
La transformée de Fourier est une fonction complexe, qui a pour chaque composante un
module et une phase.
Elle possède les mêmes propriétés que la transformation de Fourier 1D (linéarité, décalage,
dérivation, convolution, extension à L2 (R2 ).
On a vu qu’il est équivalent de multiplier la transformée de Fourier par une fonction de
transfert et d’agir sur la fonction dans la domaine spatial par convolution. Il est parfois plus
facile de mettre en œuvre le filtrage en agissant dans le domaine fréquentiel (après transformée
de Fourier ).
Filtre passe-bas
Nous en avons déjà parlé en 1D du filtre passe-bas idéal. Ici on définit une fréquence de
coupure δc au dessus de laquelle les fréquences sont annulées (filtre idéal) . La fonction de
transfert est alors 1 ,
1 si λ2 + µ2 ≤ δc
H(λ, µ) =
0 sinon
Le créneau centré H admet une transformée de Fourier inverse qui est le sinus cardnal qui
présente d’autant plus d’ondulations que la fréquence de coupure est petite. Cela entraı̂ne un
flou qui sera d’autant plus réduit que δc est grand.
Application d’un créneau « idéal » (δc . 15% de la taille de l’image) : on voit clairement les
ondulations
On a vu que le filtre passe-bas idéal n’est pas réalisable et on fait donc une approximation de
la fonction H précédente qui aura pour effet, non pas de couper les hautes fréquences mais de
les atténuer fortement. Le filtre suivant est le filtre passe-bas de Butterworth La fonction de
transfert est alors

1
H(λ, µ) = *√ +2n
λ2 +µ2
1+ δc
où δc est encore la fréquence de coupure.
En image un filtre passe-bas atténue les hautes fréquences : le résultat obtenu après un tel
filtrage est un adoucissement des détails et une réduction du bruit granuleux.
Filtres passe-haut
Le filtre passe-haut idéal est obtenu de manière symétrique au passe- bas par
1 ,
1 si λ2 + µ2 ≥ δc
H(λ, µ) =
0 sinon
Le filtre passe-haut de Butterworth est donné par
1
H(λ, µ) = * +2n
1+ √ δc
λ2 +µ2
Un filtre passe haut favorise les hautes fréquences spatiales, comme les détails, et de ce fait,
il améliore le contraste. Toutefois, il produit des effets secondaires
– Augmentation du bruit : dans les images avec un rapport Signal/ Bruit faible, le filtre
augmente le bruit granuleux dans l’image.
– Effet de bord : il est possible que sur les bords de l’image apparaisse un cadre. Mais cet
effet est souvent négligeable et peut s’éliminer en tronquant les bords de l’image.
Filtrage passe-bas et passe-haut avec un filtre de Butterworth

(n = 4 et δc . 0.15∗ taille de l’image)
3.3. FILTRAGE DIFFÉRENTIEL 45
Filtres passe-bande
Ils permettent de ne garder que les fréquences comprises dans un certain intervalle :
9 ε , ε
1 si δc − ≤ λ2 + µ2 ≤ δc +
H(λ, µ) = 2 2
0 sinon
ε est la largeur de bande et δc la fréquence de coupure.
3.3 Filtrage différentiel

Dans les modèles différentiels, on considère l’image comme une fonction continue f : I ×
I → [0, 255] dont on étudie le comportement local à l’aide de ses dérivées. Une telle étude n’a
de sens que si la fonction f est assez régulière. Ce n’est pas toujours le cas ! ! une image noir
et blanc sera discontinue (en fait continue par morceaux) les zones de discontinuité étant par
essence les contours.
Au premier ordre on peut calculer en chaque point M (x, y), le gradient de l’image :
∂f ∂f
∇f (x, y) = ( , ).
∂x ∂y
Grâce au plongement dans l’espace continu, un grand nombre d’opérations d’analyse peuvent
s’exprimer en termes d’équations aux dérivées partielles. Ceci permet de donner un fondement
mathématique satisfaisant aux traitements et aussi de fournir des méthodes pour les calculer,
par des schémas numériques de résolution.
Les filtres différentiels permettent de mettre en évidence certaines variations spatiales de
l’image . Ils sont utilisés comme traitements de base dans de nombreuses opérations comme le
réhaussement de contraste ou la détection de contours.
3.3.1 Calcul par convolution

En pratique, il faut approcher les gradients pour travailler avec des gradients discrets. Les
approximations les plus simples des dérivées directionnelles se font par différences finies cal-
∂f
culées par convolution avec des noyaux très simples : par exemple, l’approximation de se
∂x
fait par convolution avec [0 − 1 1]. En effet, dans ce cas, la formule générale de convolution
discrète (3.2.1) donne (avec dx = 3 et dy = 0 )
1 !
!
g(x, y) = f (x + i, y + j)κ(i, j) .
i=−1 j=0
κ(i, j)
0 −1 1 ←y
↑ ↑ ↑
x−1 x x+1
et
∂f
g(x, y) = −f (x, y) + f (x + 1, y) .
(x, y)
∂x
 
0
∂f
De même l’approximation de se fait par convolution avec  −1  :
∂y
1
κ(i, j)
0 ←y−1
−1 ←y
1 ←y+1
↑
x
et
∂f
g(x, y) = −f (x, y) + f (x, y + 1) . (x, y).
∂y
 
−1
On utilise plus souvent [−1 0 1] et  0  qui produisent des frontières plus épaisses mais qui
1
sont bien centrées. Ces opérations sont très sensibles au bruit et on les combine généralement
avec un filtre lisseur dans la direction orthogonale à celle de dérivation, par exemple par le
noyau suivant (ou sa transposée) : [1 2 1]. Le calcul des dérivées directionnelles en x et y revient
finalement à la convolution avec les noyaux suivants :
∂f ∂f
(x, y) . (f ∗ hx )(x, y) et (x, y) . (f ∗ hy )(x, y)
∂x ∂y
avec    
−1 0 1 −1 −1 −1
hx =  −2 0 2  et hy =  0 0 0 
−1 0 1 1 2 1
Ce sont les masques de Sobel.
De la même façon , l’approximation par différences finies la plus simple de la dérivée
 se-
1
∂2f  −2  pour
conde est la convolution par le noyau [1 − 2 1] pour l’approximation de et
∂x2
1
2
∂ f
l’approximation de . Le laplacien
∂y 2
∂2f ∂2f
∆f = +
∂x2 ∂y 2
peut donc être approché par l’un des deux opérateurs linéaires suivant
   
0 1 0 1 1 1
 1 −4 1  ou  1 −8 1 
0 1 0 1 1 1
Laplacien 4-connexe Laplacien 8-connexe
Norme du gradient
Gradient en x Gradient en y
Détection de contours par une convolution de type Sobel
T YPES DE MASQUE G RADIENTS PARTIELS A MPLITUDE D IRECTION
Masques de Roberts G1 , G2
, π
Substitution du pixel A= G21 + G22 θ=
*4 +
-1 0 0 -1 G2
supérieur gauche + arctan
0 1 1 0 G1
Masques de Sobel
1 0 -1 1 2 1 B * +
Gy
2 0 -2 0 0 0 Gx , Gy A = G2x + G2y θ = arctan
Gx
1 0 -1 -1 -2 -1
Masques de Prewitt
1 0 -1 1 1 1 B * +
Gy
1 0 -1 0 0 0 Gx , Gy A = G2x + G2y θ = arctan
Gx
1 0 -1 -1 -1 -1
Masques de Kirsh
Direction
5 5 5
-3 0 -3 Gi pour maximum des |Gi ] correspondant
-3 -3 -3
+ les 7 autres masques obtenus par i de 1 à 8 au Gi sélectionné
permutation circulaire des coefficients
Masques de Robinson
1 1 1
1 -2 1 Gi pour maximum des |Gi ] Idem
-1 -1 -1
+ les 7 autres masques obtenus par i de 1 à 8
permutation circulaire des coefficients
Laplacien discret Laplacien de Robinson
1 1 1 1 -2 1
1 -8 1 -2 4 -2
1 1 1 1 -2 1
3.3.2 Equation de la chaleur
Considérons un filtrage gaussien dans le cadre continu. On sait que si l’image de départ est
une fonction uo définie sur R2 (mais L∞ à support compact), l’image filtrée est la convolée de
uo avec un noyau gaussien
* 2 + * +
1 x1 + x22 1 1x12
Gσ (x) = Gσ (x1 , x2 ) = exp − = exp − 2 .
2πσ 2 2σ 2 2πσ 2 2σ
* +
1 1x12
On pose u(t, x) = (h(t, ·) ∗ uo )(x) où h(t, x) = G√2t (x) = 4πt exp − 4t . Comme h(t, ·) ∈
C ∞ (R2 ) a ses dérivées bornées et uo ∈ L1 (R2 ), la convolée u(t, ·) est aussi C ∞ et on peut calculer
∆u :
∂2u ∂2u
∀t > 0, ∀x ∈ R2 ∆u(t, x) = (t, x) + (t, x) = (∆h(t, ·) ∗ uo )(x) .
∂x21 ∂x22
Un rapide calcul montre que
* + * + * +
1 1x12 1x12 1 1x12
∆h(t, x) = − + exp − = − + h(t, x) ,
4πt2 16πt3 4t t 4t2
et on obtient
* +
1 1x12
∆u(t, x) = − + u(t, x) .
t 4t2
D’autre part, pour t > 0 on peut dériver directement u par rapport à t :

&&
∂u ∂h
(t, x) = (t, y)uo (x − y) dy
∂t R ∂t
et on obtient finalement
∂u
(t, x) − ∆u(t, x) = 0 sur ]0, t[×R2 .
∂t
D’autre part, avec

&& * +
1 1y12
u(t, x) = exp − uo (x − y) dy
R 4πt 4t
on obtient
lim u(t, x) =< δx , uo >= uo (x) ,
t→0
car la famille de Gaussiennes converge au sens des distributions vers la mesure de Dirac.
La fonction « filtrée » u vérifie l’équation aux dérivées partielles suivante (équation de la

chaleur) : 9
∂u
(t, x) − ∆u(t, x) = 0 dans ]0, T [×Ω
∂t (3.3.4)
u(0, x) = uo (x) ∀x ∈ Ω
où Ω est le « cadre » de l’image, i.e. l’ouvert de R2 où la fonction u est définie . On peut alors
imposer soit des conditions aux limites au bord de Ω (niveau de gris fixé) soit des conditions
aux limites périodiques en périodisant la fonction u si le cadre est rectangle par exemple.
On peut alors utiliser un schéma aux différences finies pour calculer la solution de l’EDP.
Suivant le temps d’évolution, on obtient une version plus ou moins lissée de l’image de départ.
Modèle de Peronna-Malik
Sur une image lissée, on peut plus facilement essayer de détecter les bords (ou contours).
On peut utiliser le détecteur de Hildrett-Marr : on cherche les zéros du Laplacien d’une image
u. Si en un point x, ∆u change de signe et si ∇u est non nul, on considère alors que l’image u
possède un bord en x.
Pour améliorer les résultats obtenus par l’EDP de la chaleur, Peronna et Malik ont proposé
de modifier l’équation en y intégrant le processus de détection des bords :

 ∂u
 ∂t (t, x) = div (c(|∇u|)∇u)(t, x) dans ]0, T [×Ω


∂u
= 0 sur ]0, T [×∂Ω (3.3.5)


 ∂n

u(0, x) = uo (x) ∀x ∈ Ω
où c est une fonction décroissante de R+ dans R+ .

Si c = 1, on retrouve l’équation de la chaleur. On impose souvent lim c(t) = 0 et c(0) = 1.

t→+∞
Ainsi, dans les régions de faible gradient, l’équation agit essentiellement comme l’EDP de la
chaleur, et dans les régions de fort gradient, la régularisation est stoppée ce qui permet de
préserver les bords.
Image originale Bords donnés par Hildrett-Marr

(après filtrage)
Filtrage (EDP Chaleur) Filtrage de Peronna-Malik
3.3.3 Mise en œuvre numérique

La mise en œuvre numérique se fait avec une discrétisation en différences finies. La condi-
tion de Neumann est assurée grâce à une réflexion de l’image par rapport à ses bords. En
image on considère souvent que la taille est donnée par le nombre de pixels de sorte que le pas
de discrétisation est h = 1. On peut discrétiser le gradient de différentes manières (centrée, à
droite, à gauche )
ui+1,j − ui−1,j ui,j+1 − ui,j−1
δx ui,j = , δy ui,j = , (3.3.6)
2 2
δx+ ui,j = ui+1,j − ui,j , δy+ ui,j = ui,j+1 − ui,j ,
δx− ui,j = ui,j − ui−1,j , δy− ui,j = ui,j − ui,j−1 .

La norme du gradient peut se calculer par un schéma ENO (Essentially Non Oscillatory).Deux
approximations possibles de |∇u| sont
B
∇+ ui,j = max(δx− ui,j , 0)2 + min(δx+ ui,j , 0)2 + max(δy− ui,j , 0)2 + min(δy+ ui,j 0)2 , (3.3.7)
ou
B
∇− ui,j = max(δx+ ui,j , 0)2 + min(δx− ui,j , 0)2 + max(δy+ ui,j , 0)2 + min(δy− ui,j 0)2 . (3.3.8)
Si l’opérateur gradient est discrétisé par différences finies à droite, alors une discrétisation pos-
sible de la divergence est donnée par
 1  2

 p − p 1 si 1 < i < N 
 p − p2i,j−1 si 1 < j < N
 i,j i−1,j  i,j
(div p)i,j = p1i,j si i = 1 + p2i,j si j = 1 (3.3.9)

 

 −p1 si i = N  −p 2 si j = N
i−1,j i,j−1
avec p = (p1 , p2 ) et N est la taille de l’image (carrée).

Chapitre 4
Quelques modèles de restauration

d’image
Etant donnée une image originale , on suppose qu’elle a été dégradée par un bruit additif v,
et éventuellement par un opérateur R . Un tel opérateur est souvent modélisé par un produit
de convolution , et n’est pas nécessairement inversible (et même lorsqu’il est inversible, son
inverse est souvent numériquement difficile à calculer). A partir de l’image observée f = Ru+v
(qui est donc une version dégradée de l’image originale u), on cherche à reconstruire u. Si on
suppose que le bruit additif v est gaussien , la méthode du Maximum de vraisemblance nous
conduit à chercher u comme solution du problème de minimisation
inf 1f − Ru122 ,
u
où 1 · 12 désigne la norme dans L2 . Il s’agit d’un problème inverse mal posé. Pour le résoudre
numériquement, on est amené à introduire un terme de régularisation, et à considérer le problème
inf 1f − Ru122 + L(u) .
u 2 34 5 2 34 5
ajustement aux données Régularisation
Dans toute la suite, nous ne considèrerons que le cas où est R est l’opérateur identité (Ru =
u). Commençons par un procédé de régularisation classique : celui de Tychonov.
4.1 Régularisation de Tychonov

C’est un procédé de régularisation très classique mais trop sommaire dans le cadre du trai-
tement d’image. Nous le présentons sur un exemple.
Soit V = Ho1 (Ω) et H = L2 (Ω) : on considère le problème de minimisation originel (ajuste-
ment aux données) :
(P) min 1u − ud 12H ,
u∈V
où ud est l’image observée (données) et le problème régularisé suivant : pour tout α > 0
(Pα ) min 1u − ud 12H + α1∇u12H .
u∈V
55
56 CHAPITRE 4. QUELQUES MODÈLES DE RESTAURATION D’IMAGE
Non seulement on veut ajuster u à la donnée ud , mais on impose également que le gradient soit
« assez petit » (cela dépend du paramètre α). Une image dont le gradient est petit est « lissée » ,
estompée. Les bords sont érodés et la restauration donnera une image floutée.
Il est facile de voir sur l’exemple suivant que la fonctionnnelle J(u) = 1u − ud 12H n’est pas
coercive sur V :
Ω =]0, 1[, un (x) = xn , ud = 0 .
On voit que
1 n
1un 12 = √ , 1u"n 12 = √ .
2n 2n − 1
On a donc
lim 1un 1V = +∞ et lim J(un ) = 0 .
n→+∞ n→+∞
Il n’est donc même pas clair (a priori) que le problème (P) ait une solution.
Proposition 4.1.1 Supposons que (P) admet au moins une solution ũ. Le problème (Pα ) admet une
solution unique uα . De la famille (uα ) on peut extraire une sous-suite qui converge (faiblement ) dans
V vers une solution u∗ de (P) lorsque α → 0.
Démonstration - Le problème (Pα ) admet une solution unique uα car la fonctionnelle
Jα = 1u − ud 12H + α1∇u12H ,
est coercive et strictement convexe (c’est en gros la norme de V à une partie affine près) . Mon-
trons maintenant que la famille (uα ) est bornée dans V . On pourra ainsi en extraire une sous-
suite faiblement convergente dans V vers u∗ ∈ V .
∀u ∈ V Jα (uα ) ≤ Jα (u) .
En particulier pour ũ :
J(ũ) ≤ J(uα ) ≤ Jα (uα ) = J(uα ) + α1∇uα 12H ≤ Jα (ũ) = J(ũ) + α1∇ũ12H . (4.1.1)
2 34 5 2 34 5
ũ est solution de (P) uα est solution de (Pα )
Par conséquent, pour α ≤ αo , Jα (uα ) est borné indépendamment de α. Ceci entraı̂ne que la
famille (uα )α≤αo est bornée dans H. D’autre part, avec (4.1.1) on a aussi
α1∇uα 12H ≤ J(ũ) + α1∇ũ12H − J(uα ) ≤ J(ũ) + α1∇ũ12H − J(ũ) = α1∇ũ12H ;
par conséquent la famille (uα )α≤αo est bornée dans V . On peut donc en extraire une sous-suite
qui converge (faiblement ) dans V vers u∗ . D’autre part l’équation (4.1.1) montre que
lim Jα (uα ) = J(ũ) = inf(P) .

α→0
Par semi-continuité inférieure de J il vient
J(u∗ ) ≤ lim inf J(uα ) = lim inf Jα (uα ) ≤ inf(P) .

α→0 α→0
4.2. L’ESPACE DES FONCTIONS À VARIATION BORNÉE BV (Ω) 57
Par conséquent u∗ est une solution de (P). !

Cherchons maintenant le moyen de calculer uα . Comme la fonctionnelle est strictement
convexe, une condition nécessaire et suffisante d’optimalité est
Jα" (uα ) = 0 .
Un calcul assez standard montre que

& &
∀u ∈ V Jα" (uα ) · u = (uα − ud )(x)u(x) dx + ∇uα (x)∇u(x) dx
&Ω Ω
= (uα − ud − ∆uα )(x) u(x) dx .

Ω
Par conséquent l’équation d’Euler qui fournit la solution uα est la suivante :
uα − ud − ∆uα = 0 , uα ∈ Ho1 (Ω) .
Comme pour la méthode des contours actifs décrite dans le chapitre suivant, on peut se conten-
ter d’approcher la solution uα en écrivant la formulation dynamique :
∂u
− ∆u + u = ud .
∂t
Le terme de régularisation classique L(u) = 1∇u122 ( régularisation de Tychonov) n’est pas
adapté au problème de restauration d’images : l’image restaurée u est alors beaucoup trop
lissée (en particulier, les bords sont érodés). Une approcheCbeaucoup plus efficace consiste à
considérer la variation totale , c’est à dire à prendre L(u) = |Du|.
Cela conduit à une minimisation de fonctionnelle dans un espace de Banach particulier,
mais bien adapté au problème : l’espace des fonctions à variation bornée que nous commençons
par présenter.
4.2 L’espace des fonctions à variation bornée BV (Ω)

4.2.1 Généralités
Dans ce qui suit Ω un ouvert borné de R2 de frontière Lipschitz et Cc1 (Ω, R2 ) est l’espace des
fonctions C 1 à support compact dans Ω à valeurs dans R2 .
Définition 4.2.1 Une fonction f de L1 (Ω) (à valeurs dans R) est à variation bornée dans Ω si
J(f ) < +∞ où
1& D
J(f ) = sup f (x) div ϕ(x) dx | ϕ ∈ Cc (Ω, R ) , 1ϕ1∞ ≤ 1 .
1 2
(4.2.2)
Ω
On note
BV (Ω) = {f ∈ L1 (Ω) | J(f ) < +∞ }
l’espace de telles fonctions.
Remarque 4.2.1 On rappelle que si ϕ = (ϕ1 , ϕ2 ) ∈ Cc1 (Ω, R2 ) alors

∂ϕ1 ∂ϕ2
div ϕ(x) = (x) + (x) .
∂x1 ∂x2
Donc, par intégration par parties
& & * +
∂ϕ1 ∂ϕ2
f (x) div ϕ(x) dx = f (x) (x) + (x) dx
Ω Ω ∂x1 ∂x2
& * +
∂f ∂f
=− (x)ϕ1 (x) + (x)ϕ2 (x) dx
∂x1 ∂x2
&Ω
= − ∇f (x) · ϕ(x) dx ,
Ω
où · désigne le produit scalaire de R2 .

Définition 4.2.2 (Périmètre) Un ensemble E mesurable (pour la mesure de Lebesgue) dans R2 est de
périmètre (ou de longueur) fini si sa fonction indicatrice χE est dans BV (Ω).
Commençons par donner une propriété structurelle des fonctions BV.
Théorème 4.2.1 Soit f ∈ BV (Ω). Alors il existe une mesure de Radon positive µ sur Ω et une fonction
µ-mesurable σ : Ω → R telle que
(i) |σ(x)|
& = 1, µ p.p. , et &
(ii) f (x) div ϕ(x) dx = − ϕ(x) σ(x) dµ pour toute fonction ϕ ∈ Cc1 (Ω, R2 )
Ω Ω
La relation (ii) est une formule d’intégration par parties « faible ». Ce théorème indique que la
dérivée faible (au sens des distributions ) d’une fonction BV est une mesure de Radon.
On rappelle qu’une mesure de Radon est une mesure finie sur tout compact et que grâce
au théorème de Riesz, toute forme L linéaire continue sur Cco (Ω) (fonctions continues à support
compact) est de la forme &
L(f ) = f (x) dµ ,
Ω
où µ est une (unique) mesure de Radon associée à L. Plus précisément
Théorème 4.2.2 ([Rudin] p. 126, [EG] p. 49 ) A toute forme linéaire bornée Φ sur Cco (Ω, R2 ), c’est-
à-dire
E F
∀K compact de Ω sup Φ(ϕ) | ϕ ∈ Cco (Ω, R2 ) , 1ϕ1∞ ≤ 1 , supp ϕ ⊂ K < +∞ ,
il correspond une unique mesure de Radon positive µ sur Ω et une fonction µ-mesurable σ (fonction
« signe ») telle que
(i) |σ(x)| = &1, µ p.p. , et
(ii) Φ(ϕ) = ϕ(x) σ(x) dµ pour toute fonction ϕ ∈ Cco (Ω, R2 ) .
Ω
(iii) De plus µ est la mesure de variation et vérifie
E F
µ(Ω) = sup Φ(ϕ) | ϕ ∈ Cco (Ω, R2 ) , 1ϕ1∞ ≤ 1 , supp ϕ ⊂ Ω . (4.2.3)
Démonstration du Théorème 4.2.1 - Soit f un élément de BV (Ω). On considère la forme linéaire

L définie sur Cc1 (Ω, R2 ) par &
L(ϕ) = f (x)div ϕ(x) dx .
Ω
Comme f ∈ BV (Ω),
E F
sup L(ϕ) | ϕ ∈ Cc1 (Ω, R2 ) , 1ϕ1∞ ≤ 1 := CL < +∞
pour toute fonction ϕ ∈ Cc1 (Ω, R2 ). Par conséquent
∀ϕ ∈ Cc1 (Ω, R2 ) L(ϕ) ≤ CL 1ϕ1∞ . (4.2.4)
Soit K un compact de Ω. Pour toute fonction ϕ ∈ Cco (Ω, R2 ) à support compact dans K, on peut
trouver (par densité) une suite de fonctions ϕk ∈ Cc1 (Ω, R2 ) qui converge uniformément vers ϕ.
Posons alors
L̄(ϕ) = lim L(ϕk ) .
k→+∞
Grâce à (4.2.4) cette limite existe et elle est indépendante de la suite (ϕk ) choisie. On peut donc
ainsi étendre L par densité en une forme linéaire L̄ sur Cco (Ω, R2 ) telle que
E F
sup L̄(ϕ) | ϕ ∈ Cco (Ω, R2 ) , 1ϕ1∞ ≤ 1 , supp ϕ ⊂ K < +∞ .
On conclut alors grâce au théorème de Riesz (que nous avons rappelé). !

D’après la propriété (4.2.3), J(u) = µ(Ω) ≥ 0 : c’est la variation totale de f .
Proposition 4.2.1 L’application
BV (Ω) → R+
u 5→ 1u1BV (Ω) = 1u1L1 + J(u) .
est une norme.
Démonstration - La démonstration est facile et laissée en exercice. !

On munit désormais l’espace BV (Ω) de cette norme.
Exemple 4.2.1 Supposons que
f ∈ W 1,1 (Ω) = { f ∈ L1 (Ω) | Df ∈ L1 (Ω) } ,
où Df est la dérivée de f au sens des distributions. Soit ϕ ∈ Cc1 (Ω, R2 ) telle que 1ϕ1∞ ≤ 1. Alors
& & &
f div ϕ dx = − Df · ϕ dx ≤ 1ϕ1∞ |Df | dx ≤ 1Df 1L1 < +∞ .
Ω Ω Ω
Donc f ∈ BV (Ω). De plus

&
J(f ) = sup{− Df · ϕ dx | 1ϕ1∞ ≤ 1 } = 1Df 1L1 ,
Ω
et 
 Df si Df '= 0 ,
σ= |Df |

0 sinon.
Donc
W 1,1 (Ω) ⊂ BV (Ω) .
En particulier, comme Ω est borné
∀1 ≤ p ≤ +∞ W 1,p (Ω) ⊂ BV (Ω) .
Toute fonction d’un espace de Sobolev est à variation bornée.
Remarque 4.2.2 D’après le théorème de Radon-Nikodym de décomposition des mesures, pour toute
fonction u ∈ BV (Ω), nous avons la décomposition suivante de Du (dérivée au sens des distributions) :
Du = ∇u dx + Ds u ,
où ∇u dx est la partie absolument continue de Du par rapport à la mesure de Lebesgue et Ds u est la
partie singulière.
4.2.2 Approximation et compacité

Théorème 4.2.3 (Semi-continuité inférieure de la variation totale) L’application u 5→ J(u) de
BV (Ω) dans R+ est semi-continue inférieurement (sci) pour la topologie séquentielle de L1 (Ω).
Plus précisément, si (uk ) est une suite de fonctions de BV (Ω) qui converge vers u fortement dans L1 (Ω)
alors
J(u) ≤ lim inf J(uk ) .
k→+∞
Démonstration - Soit ϕ ∈ Cc1 (Ω, R2 ) telle que 1ϕ1∞ ≤ 1. Alors

& &
u(x) div ϕ(x) dx = lim uk (x) div ϕ(x) dx .
Ω k→+∞ Ω
Donc, ∀ε > 0, ∃k(ϕ, ε) tel que

& & &
∀k ≥ k(ϕ, ε) u(x) div ϕ(x) dx − ε ≤ uk (x) div ϕ(x) dx ≤ u(x) div ϕ(x) dx + ε .
Ω Ω Ω
Comme &
uk (x) div ϕ(x) dx ≤ J(uk )
Ω
il vient &
∀k ≥ k(ϕ, ε) u(x) div ϕ(x) dx − ε ≤ J(uk ) ,
Ω
et donc &
u(x) div ϕ(x) dx ≤ lim inf J(uk ) .
Ω k→+∞
Comme c’est le cas pour tout ϕ, on obtient
J(u) ≤ lim inf J(uk ) .

k→+∞
!
Nous admettrons le résultat suivant
Théorème 4.2.4 (Approximation régulière) Pour toute fonction u ∈ BV (Ω), il existe une suite de
fonctions (uk )k∈N de BV (Ω) ∩ C ∞ (Ω) telle que
(i) uk → u dans L1 (Ω) et
(ii) J(uk ) → J(u) (dans R).
La démonstration est technique et utilise un procédé classique de régularisation par convolu-
tion. On peut se référer à [EG] p.172.
Remarquons que le résultat ci-dessus n’est pas un résultat de densité de BV (Ω) ∩ C ∞ (Ω) dans
BV (Ω) car on n’a pas J(uk − u) → 0.
Théorème 4.2.5 L’espace BV (Ω) muni de la norme
1u1BV (Ω) = 1u1L1 + J(u) ,
est un espace de Banach.
Démonstration - Soit (un )n∈N une suite de Cauchy dans BV (Ω). C’est donc une suite de Cauchy
dans L1 (Ω) : elle converge donc vers u ∈ L1 (Ω). D’autre part, elle est bornée dans BV (Ω) (toute
suite de Cauchy est bornée), donc
∃M > 0, ∀n J(un ) ≤ M .
D’après le théorème 4.2.3,
J(u) ≤ lim inf J(un ) ≤ M < +∞ .

n→∞
Par conséquent u ∈ BV (Ω). Soit ε > 0 et N tel que
∀n, k ≥ N 1un − uk 1BV (Ω) ≤ ε .
Donc
∀n, k ≥ N J(un − uk ) ≤ ε
et avec la semi-continuité inférieure de J en fixant n on obtient
∀n ≥ N J(un − u) ≤ lim inf J(un − uk ) ≤ ε .

k→∞
Ceci prouve que J(un − u) → 0. !

Terminons par un résultat important de compacité que nous admettrons.
Théorème 4.2.6 (Compacité) L’espace BV (Ω) s’injecte dans L1 (Ω) de manière compacte. Plus précisément :
si (un )n∈N une suite bornée de BV (Ω)
sup 1un 1BV (Ω) < +∞ ,

n∈N
alors, il existe une sous-suite (unk )k∈N et une fonction u ∈ BV (Ω) telle que unk converge fortement
vers u dans L1 (Ω).
Démonstration - [EG] p.176. !
Pour plus de détail sur les fonctions BV, on pourra se reporter à [EG]. Dans ce qui suit, nous
allons présenter plusieurs modèles qui se différencient par la nature du terme de régularisation.
En pratique, on pénalise par la norme de tout ou partie de la fonction u de l’image et le cadre
fonctionnel est très important.
4.3 Le modèle continu de Rudin-Osher-Fatemi

4.3.1 Présentation
Il faut trouver une alternative à la régularisation de Tychonov qui est trop violente. La
première idée consiste à remplacer le terme de régularisation 1∇u12H qui est en réalité une
pénalisation par une norme V , par une norme moins contraignante et régularisante. Rudin-
Osher et Fatemi ont proposé un modèle où l’image est décomposée en deux parties : ud = u + v
où v est le bruit et u la partie « régulière ». On va donc a priori chercher la solution du problème
sous la forme u + v avec u ∈ BV (Ω) et v ∈ L2 (Ω) et ne faire porter la régularisation que sur la
partie « bruit ». Cela conduit à :
1
(PROF ) min{ J(u) + 1v122 | u ∈ BV (Ω), v ∈ L2 (Ω) , u + v = ud } .
2ε
Ici J(u) désigne la variation totale de u et ε > 0.
Théorème 4.3.1 Le problème (PROF ) admet une solution unique.
Démonstration - Soit un ∈ BV (Ω), vn ∈ L2 (Ω) une suite minimisante. Comme vn est bornée
dans L2 (Ω) on peut en extraire une sous-suite (notée de la même façon) faiblement convergente
vers v ∗ dans L2 (Ω). Comme la norme de L2 (Ω) est convexe, sci, il vient
1v ∗ 122 ≤ lim inf 1vn 122 .

n→+∞
De la même façon un = ud − vn est bornée dans L2 (Ω) et donc dans L1 (Ω) puisque Ω est
borné. Comme J(un ) est borné , il s’ensuit que un est bornée dans BV (Ω). Grâce au résultat de
compacité du théorème 4.2.6 , cela entraı̂ne que un converge (à une sous-suite près) fortement
dans L1 (Ω) vers u∗ ∈ BV (Ω).
D’autre part J est sci (théorème (4.2.3) , donc
J(u∗ ) ≤ lim inf J(un ) ,

n→+∞
4.4. MODÈLE DISCRET 63
et finalement
1 ∗ 2 1
J(u∗ ) + 1v 12 ≤ lim inf J(un ) + 1vn 122 = inf(PROF ) .
2ε n→+∞ 2ε
Comme un + vn = ud pour tout n, on a u∗ + v ∗ = ud . Par conséquent u∗ est une solution du
problème (PROF ).
La fonctionnelle est strictement convexe par rapport au couple (u, v) et la contrainte est affine.
On a donc unicité. !
Nous aurons besoin d’établir des conditions d’optimalité pour la ou les solutions optimales
des modèles proposés. Toutefois les fonctionnelles considérées (en particulier J) ne sont en
général pas Gâteaux-différentiables et nous devons utiliser des notions d’analyse non lisse (voir
chapitre 3).
4.3.2 Condition d’optimalité du premier ordre

Le problème (PROF ) peut s’écrire de la manière (équivalente) suivante
1
min F(u) := J(u) + 1u − ud 122 . (4.3.5)
u∈BV (Ω) 2ε
La fonctionnelle F est convexe et
ū solution de (PROF ) ⇐⇒ 0 ∈ ∂F(ū) .
On aimerait utiliser la relation
0 ∈ ∂F(ū) ⇐⇒ ū ∈ ∂F ∗ (0) ,
mais elle n’est vraie que si l’espace est réflexif, ce qui n’est pas le cas de BV (Ω). Toutefois ce
sera vrai après discrétisation (ce que nous ferons dans la section suivante). On peut utiliser le
théorème A.3.4 pour calculer ∂F(u). L’application u 5→ 1u − ud 122 est continue sur L2 (Ω) est
J est finie sur BV (Ω) à valeurs dans R ∪ {+∞} e. D’autre part u 5→ 1u − ud 122 est Gâteaux
différentiable sur L2 (Ω). Le calcul du sous-différentiel pouvant se faire via la transformation
de Fenchel-Legendre dans le cas où l’espace est réflexif (voir Annexe). Nous allons donc faire
la discrétisation du problème.
4.4 Modèle discret

On va maintenant considérer des images discrètes (ce qui est le cas en pratique). Une image
discrète est une matrice N × N que nous identifierons à une vecteur de taille N 2 (par exemple
en la rangeant ligne par ligne). On note X l’espace euclidien RN ×N et Y = X × X. On munit X
du produit scalaire usuel
!
(u, v)X = uij vij ,
1≤i,j≤N
et de la norme associée : 1 ·1 X .
Nous allons donner une formulation discrète de ce qui a été fait auparavant et en particulier
définir une variation totale discrète que nous noterons de la même façon. Pour cela nous intro-
duisons une version discrète de l’opérateur gradient. Si u ∈ X, le gradient ∇u est un vecteur
de Y donné par
(∇u)i,j = ((∇u)1i,j , (∇u)2i,j ) ,
avec 1
ui+1,j − ui,j si i < N
(∇u)1i,j = (4.4.6)
0 si i = N
1
ui,j+1 − ui,j si j < N
(∇u)2i,j = (4.4.7)
0 si j = N
La variation totale discrète est alors donnée par
!
J(u) = |(∇u)i,j | , (4.4.8)
1≤i,j≤N
B
où |(∇u)i,j | = |(∇u)1i,j |2 + |(∇u)2i,j |2 . On introduit également une version discrète de l’opérateur
de divergence. On le définit par analogie avec le cadre continu en posant
div = −∇∗ ,
où ∇∗ est l’opérateur adjoint de ∇, c’est-à-dire
∀p ∈ Y, ∀u ∈ X (−div p, u)X = (p, ∇u)Y = (p1 , ∇1 u)X + (p2 , ∇2 u)X .
On peut alors vérifier que
 1  2

 p − p1i−1,j si 1 < i < N 
 p − p2i,j−1 si 1 < j < N
 i,j  i,j
(div p)i,j = p1i,j si i = 1 + p2i,j si j = 1 (4.4.9)

 

 −p1 si i = N  −p2 si j = N
i−1,j i,j−1
Nous utiliserons aussi une version discrète du laplacien définie par

∆u = div ∇u .
4.4.1 Le problème discret

On va remplacer le problème (PROF ) (dans la version (4.3.5)) par le problème obtenu après
discrétisation suivant
1
min J(u) + 1u − ud 12X . (4.4.10)
u∈X 2ε
Il est facile de voir que
B ce problème a une solution unique que nous allons caractériser. On
rappelle que |gi,j | = (gi,j
1 )2 + (g 2 )2 et que la version discrète de la variation totale est donnée
i,j
(de manière analogue au cas continu) par
J(u) = sup :u, ξ; .
ξ∈K
4.4. MODÈLE DISCRET 65
où
K := {ξ = div (g) | g ∈ Y, |gi,j | ≤ 1 , ∀i, j } ,
est la version « discrète » de
E F
ξ = div ϕ | ϕ ∈ Cc1 (Ω, R2 ) , 1ϕ1∞ ≤ 1 .
Donnons tout d’abord une caractérisation de la conjuguée de Fenchel de J.
Théorème 4.4.1 La transformée de Fenchel J ∗ de la fonctionnelle J « variation totale » approchée
définie sur X, est l’indicatrice de l’ensemble K̄, où
K = {ξ = div (g) | g ∈ Y, |gi,j | ≤ 1 , ∀i, j } ,
Démonstration - Comme J est positivement homogène la conjuguée J ∗ de J est l’indicatrice

d’un ensemble convexe fermé K̃ (proposition A.3.3).
Montrons d’abord que K ⊂ K̃ : soit u ∈ K. Par définition de J
J(u) = sup :ξ, u; ,
ξ∈K
où :·, ·; désigne le produit de dualité entre BV (Ω) et son dual. Par conséquent :ξ, u; − J(u) ≤ 0
pour tous ξ ∈ K et u ∈ BV (Ω). On déduit donc que pour tout u∗ ∈ K
J ∗ (u∗ ) = sup :u∗ , u; − J(u) ≤ 0 .
u∈K
Comme J ∗ ne prend qu’une seule valeur finie on a J ∗ (u∗ ) = 0, et donc u∗ ∈ K̃. Par conséquent
K ⊂ K̃ et comme K̃ est fermé :
K̄ ⊂ K̃ .
En particulier
J(u) = sup :u, ξ; ≤ sup :u, ξ; ≤ sup :u, ξ; = sup :u, ξ; − J ∗ (ξ) = J ∗∗ (u) .
ξ∈K ξ∈K̄ ξ∈K̃ ξ∈K̃
Comme J ∗∗ = J, il vient
sup :u, ξ; ≤ sup :u, ξ; ≤ sup :u, ξ; ,
ξ∈K ξ∈K̃ ξ∈K
et donc
sup :u, ξ; = sup :u, ξ; = sup :u, ξ; . (4.4.11)
ξ∈K ξ∈K̄ ξ∈K̃
Supposons maintenant qu’il existe u∗ ∈ K̃ tel que u∗ ∈ / K̄. On peut alors séparer strictement u∗
et le convexe fermé K̄. Il existe α ∈ R et uo tels que
:uo , u∗ ; > α ≥ sup :uo , v; .
v∈K̄
D’après (4.4.11) il vient

sup :uo , ξ; ≥ :uo , u∗ ; > α ≥ sup :uo , v; = sup :uo , v; .
ξ∈K̃ v∈K̄ v∈K̃
On a donc une contradiction : K̃ = K̄ . !

4.4.2 Algorithme de projection de Chambolle

Le résultat suivant donne la caractérisation attendue de la solution [Cham] :
Théorème 4.4.2 La solution de (4.4.10 ) est donnée par
u = ud − PεK (ud ) , (4.4.12)
où PK est le projecteur orthogonal sur K.
Démonstration - D’après la définition du sous-différentiel, u est une solution de (4.4.10 ) est
équivalent * +
1 u − ud
0 ∈ ∂ J(u) + 1u − ud 1X = 2
+ ∂J(u) .
2ε ε
Comme J est convexe, sci, propre on peut appliquer le corollaire A.3.3. Donc
ud − u ud − u ud − u
∈ ∂J(u) ⇐⇒ u ∈ ∂J ∗ ( ) ⇐⇒ 0 ∈ −u + ∂J ∗ ( ).
ε ε ε
Ceci est aussi équivalent à
ud − u ud 1 ∗ ud − u
0∈ − + ∂J ( ).
ε ε ε ε
ud − u
On en déduit que w = est un minimiseur de
ε
1w − uεd 12 1 ∗
+ J (w) .
2 ε
ud − u ud
Comme J ∗ est l’indicatrice de K, cela implique que est la projection orthogonale de
ε ε
ud 1
sur K. Comme PK ( ) = PεK (ud ), on peut conclure.
ε ε
!
Tout revient maintenant à calculer
PεK (ud ) = argmin { 1ε div (p) − ud 12X | |pi,j ≤ 1, i, j = 1, · · · , N } .
On peut le résoudre par une méthode de point fixe : po = 0.
pni,j + ρ (∇[ div pn − ud /ε])i,j
i,j =
pn+1 G
G
G .
G (4.4.13)
1 + ρ G(∇[ div pn − ud /ε])i,j G
Théorème 4.4.3 Si le paramètre ρ dans (4.4.13) vérifie ρ ≤ 1/8, alors

ε div pn → PεK (ud ) .
La solution du problème est donc donnée par

u = ud − ε div p∞ où p∞ = lim pn .
n→+∞
Chapitre 5
Méthode des contours actifs
5.1 Rappels sur la géométrie des courbes planes

Dans ce qui suit on considèrera des courbes planes paramétrées connexes, compactes.
5.1.1 Abscisse curviligne - longueur

On considère une courbe paramétrée (Γ, Φ) : Φ est une application d’un intervalle I de R2 .
C’est une paramétrisation de Φ :
Γ = {M (x, y) ∈ R2 | (x(t), y(t)) = Φ(t), t ∈ I }.
On suppose que I = [a, b] est un intervalle compact de R. Soit Σ l’ensemble des subdivisions σ
de [a, b] :
σ = {to = a, t1 , · · · , tn−1 , tn = b } .
On pose
n−1
!
"σ (Γ) = d(Φ(tn+1 ), Φ(tn )) ,
k=0
où d est la distance euclidienne dans R2 .
Définition 5.1.1 On dit que la courbe Γ est rectifiable si
sup "σ (Γ) < +∞ .

σ∈Σ
Ce nombre est la longueur de la courbe Γ : "(Γ). Elle est indépendante de la paramétrisation (régulière)
choisie.
Remarquons que la longueur vérifie la relation de Chasles.

On suppose maintenant que la courbe (c’est-à-dire Φ) est de classe C 1 .
67
68 CHAPITRE 5. MÉTHODE DES CONTOURS ACTIFS
Définition 5.1.2 Soit Γ une courbe rectifiable, connexe, compacte de classe C 1 . Soit to ∈ I. On pose
1
"(Γt,to ) si t ≥ to ,
S(t) =
−"(Γto ,t ) si t ≤ to ,
où Γa,b désigne la courbe obtenue lorsque t ∈ [a, b].
Théorème 5.1.1 La fonction réelle S est dérivable sur I et
dM
S " (t) = |Φ" (t)| = | (t)| ,
dt
où | · | désigne la norme euclidienne de R2 et M (t) = M (Φ(t)) = M (x(t), y(t)).
Corollaire 5.1.1 Sous les hypothèses précédentes la longueur de la courbe s’exprime par la formule :
& b
dM
"(Γ) = | (t)| dt . (5.1.1)
a dt
La fonction S est continue dans l’intervalle I. Son image J = S(I) ⊂ R est donc un intervalle
et la variable s = S(t) décrivant J est appelée abscisse curviligne de la courbe Γ. Le point
Mo = M (to ) est l’origine de l’abscisse curviligne.
On rappelle qu’un point régulier de Γ est un point où Φ" (t) '= 0. Si l’ensembles des valeurs du
paramètre t telles que M (t) soit un point régulier de Γ est dense dans I, la fonction S est un
homéomorphisme de I dans J = S(I). On peut donc redéfinir un paramétrage Ψ = Φ ◦ S −1 de
Γ appelé paramétrage de Γ par une abscisse curviligne.
Supposons maintenant que la courbe Γ soit paramétrée par une abscisse curviligne : Φ :
[a, b] → R2 . La longueur de la courbe est une quantité intrinsèque mais on peut l’écrire en
fonction de la paramétrisation : "(Γ) = "(Φ). On alors le résultat suivant qui nous servira par la
suite
Théorème 5.1.2 Soit X un ensemble de fonctions C 1 de R sur R2 représentant des courbes fermées
(valeurs de la fonction et de sa dérivée égales aux bornes) . Alors la fonctionnelle
" : X → R+
Φ 5→ "(Φ)
est Gâteaux différentiable en toute fonction Φ dont le gradient ∇Φ est non nul et
& b * +
d" ∇Φ
∀h ∈ X (Φ) · h = div (s) h(s) ds ,
dΦ a |∇Φ|
où | · | désigne la norme euclidienne de R2 .

5.1. RAPPELS SUR LA GÉOMÉTRIE DES COURBES PLANES 69
,
Démonstration - L’application N : R2 → R+ définie par N (x) = |x| = x21 + x22 est dérivable
en tout point x '= 0 et
x
∇N (x) = .
|x|
Si Ψ est une fonction C 1 non nulle de [a, b| dans R2 , le théorème des fonctions composées donne
d|Ψ| Ψ
(Ψ) = .
dΨ |Ψ|
Soit h ∈ X : calculons la Gâteaux-dérivée de " :
& b
d" ∇Φ
(Φ) · h = (s) · ∇h(s) ds .
dΦ a |∇Φ|
Une intégration par parties couplée à la condition aux limites donne (la courbe est fermée)
& b
d" ∇Φ
(Φ) · h = div( (s))h(s) ds .
dΦ a |∇Φ|
!
5.1.2 Étude géométrique locale d’une courbe paramétrée

On rappelle que si M (t) est un point régulier d’une courbe C 1 un vecteur tangent à la courbe
en M est défini par
dM
T = (t) = Φ" (t) .
dt
On appelle vecteur normal à la courbe un vecteur N orthogonal à T .
Supposons maintenant que la courbe est C 2 et que tous ses points sont réguliers. Le pa-
ramétrage par l’abscisse curviligne est alors aussi C 2 . On suppose donc que la courbe Γ est
paramétrée par une abscisse curviligne. Alors le vecteur tangent
dM
T (s) = (s)
ds
dT
est dérivable et le vecteur (s) est indépendant du choix de l’abscisse curviligne s.
ds
Définition 5.1.3 La courbure de Γ au point M (so ) est le nombre réel
dT
ρ(so ) = | (so )| ≥ 0 .
ds
Un point M de Γ est bi-régulier si et seulement si sa courbure est non nulle.

Précisons enfin tout cela dans un repère du plan. Soit Γ une courbe C 2 , connnexe, dont tous
les points sont réguliers. On la paramètre par une abscisse curviligne s. Étant donné un repère
(&ı, &j) on pose :
T (s) = cos(φ(s))&ı + sin(φ(s))&j.
dφ
La dérivée (s) est indépendante du choix du repère et de φ.
ds
Définition 5.1.4 La courbure algébrique de Γ au point M (s) est le nombre réel

dφ
ρa (s) = (s) .
ds
On peut alors montrer que

dT
ρ(s) = |ρa (s)| et (s) = ρ(s)N (s) .
ds
5.2 Méthodes des contours actifs

5.2.1 Introduction
La segmentation permet d’isoler certaines parties de l’image qui présentent une forte corrélation
avec les objets contenus dans cette image, généralement dans l’optique d’un post-traitement.
Les domaines d’application sont nombreux : médecine, géophysique, géologie, etc ... Dans
le domaine médical, la segmentation d’images est extrêmement compliquée. En effet, pour
chaque organe (cerveau, cœur, etc ...), l’approche est différente : l’outil de segmentation doit
donc pouvoir s’adapter à un organe particulier, suivant une modalité d’acquisition particulière
(scanners, radiographie, Imagerie par Résonance Magnétique, ...) et pour une séquence de
données particulière. L’objectif est la quantification de l’information, par exemple, la volumétrie :
volume d’une tumeur dans le cerveau, étude de la cavité ventriculaire cardiaque, etc ... C’est
à ce niveau que la segmentation de l’image est utilisée. En géophysique, la segmentation peut
permettre d’isoler des objets du sous-sol (failles, horizons ...) à partir de données sismiques
dans le but, par exemple, de modéliser ou d’exploiter un gisement.
En imagerie mathématique, deux types de segmentation sont exploités :
• la segmentation par régions qui permet de caractériser les régions d’une image présentant
une structure homogène
• la segmentation par contours qui permet de délimiter les différentes régions par leurs frontières.
C’est à ce dernier type de contourage que nous allons nous intéresser
Les points de contours sont les points de l’image pour lesquels la norme du gradient, dans
la direction de ce gradient, est maximale. Un seuillage est réalisé pour ne conserver que les
points de variation de niveau de gris significative. Une question légitime qui émerge alors, est
comment définir le seuil, autrement dit, pour quel critère de variation de niveau de gris un
point sera qualifié de point de contour ?
Un seuil choisi trop faible accordera l’étiquette « point de contour » à un nombre trop élevé
de points tandis qu’un seuil trop élevé ne permettra d’extraire que les points de fort contraste
et les contours détectés ne seront plus connexes. La représentation mathématique des frontières
d’un objet ne sera plus réalisée dans ce cas.
Pour pallier cette difficulté, une régularité sur la modélisation des contours doit être intro-
duite : les contours seront assimilés à des courbes possédant des propriétés de régularité et
satisfaisant le critère de détection énoncé précédemment.
Dans ce chapitre on s’intéresse à la méthode des contours actifs (encore appelés « snakes »)
introduits par Kass, Witkin et Terzopoulos [KWT], méthode qui intègre cette notion de régularité
5.2. MÉTHODES DES CONTOURS ACTIFS 71
des points de contour en introduisant une fonctionnelle interprétée en terme d’énergie pour les
propriétés mécaniques qu’elle revêt. En effet, cette méthode permet de faire évoluer en temps et
en espace la représentation du modèle vers la solution du problème de minimisation introduit
dans la modélisation. Ces méthodes de contours actifs font appel à la notion de corps élastique
subissant des contraintes extérieures. La forme prise par l’élastique est liée à une minimisation
d’énergie composée de deux termes :
• un terme d’énergie interne et
• un terme d’énergie externe.
Le minimum local obtenu par minimisation de cette fonctionnelle non-convexe est lié à la
condition initiale qui définit un voisinage de recherche du minimum.
De manière générale, les difficultés ma jeures rencontrées dans le processus de segmenta-
tion par contours actifs sont les suivantes :
– le modèle est non-intrinsèque du fait de la paramétrisation et n’est donc pas lié à la
géométrie de l’objet à segmenter.
– forte dépendance du modèle à la condition initiale, qui n’autorise pas à choisir une condi-
tion initiale éloignée de la solution.
– connaissance de la topologie de/ou/des objets à segmenter nécessaire, ce qui implique
(lorsqu’il y a plusieurs objets à segmenter) l’utilisation de procédures particulières du fait
de la paramétrisation.
– complexité des images / ambiguı̈té des données correspondantes, en particulier lorsque
les données des images manquent et/ou que deux régions de texture similaire sont adja-
centes.
– bruit sur les données.
5.2.2 Modélisation du problème

On peut définir un contour actif ou « snake »comme une courbe fermée qui minimise son
énergie, influencée par une contrainte interne et guidée par une force d’image qui pousse
la courbe vers les contours présents dans l’image. L’espace des formes est l’ensemble Φ des
courbes paramétrées (par l’abscisse curviligne) suivant :
1 D
v : [0, 1] → R2
Φ= v| , v(0) = v(1) .
s 5→ (x(s), y(s))
L’énergie du « snake » peut s’écrire sous la forme :
Es (v) = Ei (v) + Ee (v) .
Elle est constituée d’un terme de régularisation interne (Ei ) et d’un terme de potentiel d’at-
traction ou ajustement aux données (Ee ).
L’énergie interne peut se décomposer comme suit :
&
1 1H I
Ei (v) = α(s)|v " (s)|2 + β(s) |v”(s)|2 ds (5.2.2)
2 0
où |·| désigne la norme euclidienne dans R2 . α ∈ L∞ (0, 1) est le coefficient d’élasticité (résistance
à l’allongement) et β ∈ L∞ (0, 1) le coefficient de rigidité. Le premier terme de Ei (v) :
& 1
α(s)|v " (s)|2 ds
0
pénalise la longueur du « snake » ( augmenter α(s) tend à éliminer les boucles en réduisant la
longueur du contour).
Le second terme
& 1
β(s) |v”(s)|2 ds
0
pénalise la courbure ( augmenter β(s) tend à rendre le « snake » moins flexible).

On fait donc apparaı̂tre dans cette expression des propriétés mécaniques du comportement
d’un élastique (dérivée du premier ordre) et d’une poutre (dérivée du second ordre). Les modèles
déformables se comportent comme des corps élastiques qui répondent naturellement aux forces
et aux contraintes qui leur sont appliquées. Notons que le choix β = 0 autorise néanmoins les
discontinuités du second ordre.
En ce qui concerne l’expression de l’énergie externe, plusieurs expressions liées à une fonction
potentielle sont à notre disposition. Les contours que l’on souhaite déterminer sont :
– soit assimilés aux points de fort gradient de I , image donnée : l’expression de Ee fait
alors apparaı̂tre la norme du gradient :
& 1 & 1
Ee (v) = P (v(s)) ds := −λ |∇I(v(s))|2 ds , (5.2.3)
0 0
avec λ > 0.
– soit assimilés aux points de dérivée seconde nulle : si Gσ désigne un filtre gaussien, on
peut choisir
& 1 & 1
Ee (v) = P (v(s)) ds := λ |Gσ ∗ ∆I(v(s))|n ds , (5.2.4)
0 0
La modélisation du problème que l’on propose revient donc à trouver une fonction v̄ ∈
Φ∩X
min{Es (v) | v ∈ Φ ∩ X } . (5.2.5)
Ici X est un espace fonctionnel adapté à la formulation du problème de minimisation. En

effet, nous avons décrit les choses formellement jusqu’à présent. Il est clair que pour définir les
différentes énergies, la fonction v doit être dérivable au sens des distributions et de dérivées
intégrables. En pratique
X = H 1 (0, 1) ou H 2 (0, 1) .
Il s’agit donc d’un problème de minimisation de fonctionnelle (minimisation de l’énergie) que

l’on peut résoudre au cas par cas en utilisant les théorèmes et les méthodes classiques d’opti-
misation présentées dans le chapitre ??, section A.1.3 par exemple.
5.2.3 Conditions d’optimalité

Supposons avoir démontré que le problème (5.2.5) admet au moins une solution v̄. Nous
allons utiliser le théorème d’Euler-Lagrange (condition d’optimalité du premier ordre : cha-
pitre ??, Théorème A.1.5 ) pour déterminer l’équation aux dérivées partielles que satisfait v̄.
Supposons (pour simplifier) que X = H 2 (0, 1) de sorte que v " ∈ L2 (0, 1) et v” ∈ L2 (0, 1). On
supposera également que la fonction P : X → L1 (0, 1) est dérivable (au sens des distributions
) et que sa dérivée est dans L1 (0, 1).(Cet ensemble se note W 1,1 (0, 1)). La fonctionnelle coût
considérée est
& & 1
1 1H "
I
Es (v) = α(s)|v (s)| + β(s) |v”(s)| ds +
2 2
P (v(s)) ds ,
2 0 0
sur l’espace X . Elle n’est en général pas convexe de sorte que si l’on peut (éventuellement)
montrer l’existence d’une solution ū , celle ci n’est en général pas unique. D’après le Théorème
A.1.5, si ū est solution du problème alors
dEs
∀v ∈ X (ū) · v = 0 .
dv
Lemme 5.2.1 La fonctionnelle Es est Gâteaux-différentiable et sa Gâteaux-dérivée en ū, dans la direc-

tion v est
& 1
dEs H I
(ū) · v = α(s)ū" (s)v " (s) + β(s) ū”(s) v”(s) + ∇v P (ū(s))v(s) ds .
dv 0
Démonstration - La démonstration, facile, est laissée en exercice. !

Soit alors v ∈ D(0, 1) : on obtient
& 1J * + * 2 + K
d dū d2 d ū
α (s)v(s) + 2 β 2 (s) v(s) + ∇v P (ū)(s)v(s) ds = 0 ,
0 ds ds ds ds
c’est-à-dire * + * +
d dū d2 d2 ū
α + 2 β 2 + ∇v P (ū) = 0 ,
ds ds ds ds
au sens des distributions. Il s’agit ensuite de donner des conditions aux limites : ici elles pro-
viennent de la paramétrisation par abscisse curviligne :
dv dv
v(0) = v(1) = (0) = (1) .
ds ds
En définitive on est ramené à l’étude et à la résolution (numérique) d’une équation aux dérivées
partielles : * + * 2 +
 2
 d α dū + d

d ū
β 2 + ∇v P (ū) = 0 dans D" (0, 1)

 ds 2
ds ds ds
dv dv (5.2.6)

 v(0) = v(1) = (0) = (1)

 ds ds
v ∈ Φ ∩ H 2 (0, 1) .
5.2.4 Un modèle dynamique

L’équation écrite ci-dessus exprime un état statique ou état d’équilibre quand on est au mi-
nimum d’énergie. Mais bien souvent, alors même qu’on connaı̂t l’existence d’un infimum on
n’est pas capable de montrer l’existence d’un minimum et il faut approcher cet état d’équilibre
plutôt que le chercher exactement (il n’est parfois pas atteignable). Le principe des modèles
déformables et de considérer que le contour que l’on cherche est l’état d’équilibre d’un contour
qui va évoluer avec le temps. Le problème stationnaire est transformé en un problème dyna-
mique : on cherche la solution u(x, t) d’un problème d’évolution. La courbe cherchée ū est alors
donnée par
ū(x) = lim u(t, x) .
t→+∞
Une manière simple d’imposer un mouvement au contour est d’imposer sa vitesse d’évolution
∂u
(x, t) en posant
∂t
∂u
(x, t) = −∇u Es (u(x, t)) . (5.2.7)
∂t
En effet, la famille u(t, ·) de courbes (indexée par t) que l’on cherche doit être choisie de façon
à faire décroı̂tre l’énergie Es qui se rapprochera ainsi de son infimum. Si on fait (formellement)
un développement de Es au premier ordre, on a pour δt > 0
E(u(t + δt, ·)) − E(u(t, ·)) . :∇u E(u(t, ·)), u(t + δt, ·) − u(t, ·); . (5.2.8)
Le choix « u(t + δt, ·) − u(t, ·) = −δt∇u E(u(t, ·) »montre qu’on fait bien décroı̂tre l’énergie. En
faisant tendre δt vers 0 on obtient (5.2.7). Ceci conduit à une équation aux dérivées partielles
d’évolution (en général parabolique à laquelle il convient d’ajouter une condition initiale
u(t = 0, x) = uo (x) (on se donne le contour de départ)
et des conditions aux limites issues de l’analyse du problème stationnaire.

On calcule ensuite cette solution et on s’intéresse à son comportement symptotique (en
temps long).
En pratique, numériquement il suffit de s’arrêter lorsque t est assez grand (c’est-à-dire quand
deux valeurs consécutives sont assez voisines ).
5.2.5 Un exemple
Nous allons illustrer ce qui précède par un exemple. Supposons que la fonction α soit
constante, que β soit nulle (on ne contraint pas la rigidité du contour actif ). L’espace fonc-
tionnel des courbes est alors X = H 2 (0, 1) ∩ Ho1 (0, 1) (en prenant en compte la condition aux
limites v(0) = v(1) = 0. On supposera que l’image I est à dérivée dans L2 et que
v 5→ ∇I(v)
est continue de X dans L2 (0, 1). Choisissons dans un premier temps

1
P (v) = − λ|∇I(v)|2 .
2
La fonctionnelle d’énergie est

&
1 1H " I
Es (v) = α|v (s)|2 − λ|∇I(v(s))|2 ds , v ∈ H 2 (0, 1) ∩ Ho1 (0, 1).
2 0
On remarque immédiatement que la fonctionnelle Es peut ne pas être minorée ! ! On peut alors
prendre pour P (v) une fonction dite de détection de contours moins « sommaire » que celle que
nous venons de choisir, par exemple
1
P (v) = .
1 + |∇I(v)|2
Cette fonction a le mérite d’être positive et la minimiser revient à maximiser son dénominateur
donc le gradient de l’image. La fonctionnelle d’énergie est à présent
& J K
1 1 " λ
Es (v) = α|v (s)| +
2
ds , v ∈ H 2 (0, 1) ∩ Ho1 (0, 1).
2 0 1 + |∇I(v)|2
Cette fonctionnelle est minorée, donc l’infimum existe. Toutefois il est encore extrêmement
délicat de montrer qu’il est atteint, les propriétés de semi-continuité et surtout de coercivité de
Es n’étant pas évidentes. Dans ce cas, on adopte la démarche « dynamique » comme dans la
section précédente. Calculons :∇u E(u(t, ·)), u(t + δt, ·) − u(t, ·); pour δ > 0. Un calcul analogue
à la section précédente donne
:∇u E(u(t, ·)), u(t + δt, ·) − u(t, ·); =
& 1H I
αu" (t, s)[u" (t + δt, s) − u" (t, s)] + ∇u P (u)(t, s)[u(t + δt, s) − u(t, s)] ds .
0
Une intégration par parties couplée aux conditions aux limites de la forme (courbe fermée) :
dv dv
v(0) = v(1) = (0) = (1)
ds ds
donne
& 1H I
:∇u E(u(t, ·)), u(t + δt, ·) − u(t, ·); = −αu"" (t, s) + ∇u P (u)(t, s) [u(t + δt, s) − u(t, s)] ds .
0
Si on choisit
H I
u(t + δt, s) − u(t, s) = −δt −αu"" (t, s) + ∇u P (u)(t, s) p.p.
l’énergie va décroı̂tre (pour δt assez petit) grâce à l’approximation (5.2.8). On obtient
u(t + δt, s) − u(t, s)
= αu"" (t, s) − ∇u P (u)(t, s) p.p.
δt
et par passage à la limite lorsque δt → 0


 ∂u ∂2u

 − α 2 + ∇v P (u) = 0 dans D" (0, 1)
 ∂t ∂s
du du (5.2.9)

 ∀t > 0 u(t, 0) = u(t, 1) = (t, 0) = (t, 1) ,

 ds ds
 ∀s ∈ [0, 1] u(0, s) = u (s) donnée.
o
C’est une EDP parabolique non linéaire dont l’étude est en général classique.
5.3 La méthode des lignes de niveau (« Level set »)
Le principe des contours actifs est de faire évoluer une courbe. On a vu dans la section
∂u
précédente une formulation dynamique qui fait intervenir c’est-à-dire la vitesse d’évolution
∂t
du contour. On va donc s’intéresser à la façon de faire évoluer la courbe et plus généralement
à la notion de propagation de fronts.
5.3.1 Propagation de fronts
On se donne une courbe plane (on peut aussi se placer en 3D avec une surface) fermée que
l’on supposera régulière (on précisera cela plus tard). Cette courbe partage le plan en deux
régions, l’intérieur et l’extérieur. On oriente la courbe de façon à définir une normale extérieure
&n. On se donne la vitesse de propagation de la courbe : cette vitesse est portée par la normale :
&v = F &n ,
car on supposera qu’il n’y pas de déplacement tangentiel (pas de rotation ou de glissement par
exemple). La fonction F dépend de plusieurs facteurs :
1. locaux déterminés par l’information géométrique locale (comme la courbure et la nor-

male)
2. globaux : ce sont ceux qui dépendent de la forme (globale) ou de la position du front. La

vitesse peut par exemple inclure une intégrale sur tout le domaine.
Le problème le plus important est la modélisation du front de déplacement de la courbe, c’est-

à-dire la formulation de l’expression de la vitesse F . Nous allons supposer ici que la vitesse est
connue. On donnera ensuite quelques exemples.
On se donne donc une courbe régulière (par exemple C 2 dont tous les points sont réguliers)
Γo de R2 . La famille obtenue par déplacement de long de la normale à la vitesse F est notée Γt .
& = (x, y) est la position d’un point de la courbe et &n la normale extérieure (unitaire ) à la
Si φ
& On considère une paramétrisation de la courbe Γt par une abscisse
courbe, on aura F = &n · φ.
curviligne :
φ(s, & t) = φ(S,

& t) = (x(s, t), y(s, t)) où s ∈ [0, S] et φ(0, & t) .
On paramètre la courbe de façon que l’intérieur soit à gauche de la direction des s croissants.
&n(s, t) est une paramétrisation de la normale extérieure et κ(t, s) une paramétrisation de la
courbure.
5.3. LA MÉTHODE DES LIGNES DE NIVEAU (« LEVEL SET ») 77
On rappelle que dans ce cas la courbure est donnée par

J K
yss xs − xss ys
κ(t, s) = (t, s)
(x2s + ys2 )3/2
où xs désigne la dérivée par rapport à s et xss désigne la dérivée seconde par rapport à s.
On se restreint dans ce qui suit à des vitesses ne dépendant que de la courbure locale
F = F (κ) .
On obtient alors les équations donnant le mouvement de Γt :

 . /. /
yss xs −xss ys ys
 xt = F (x
 2 +y 2 )3/2
s s
2 2
(xs +ys ) 1/2
. / . / (5.3.10)

 yt = F yss xs −xss ys xs
2
(x +y )2 3/2 2
(x +y )2 1/2 .
s s s s
& est aussi sa lon-

La variation totale de la courbe Γt (ou de la fonction qui la paramètre φ)
gueur. Elle est donnée par
& S
"(t) = |κ(s, t)|(x2s + ys2 )1/2 ds .
0
La proposition suivante donne une idée de l’évolution du front lorsque la courbure s’annule
(courbe non convexe).
Proposition 5.3.1 On considère un front évoluant à la vitesse F (κ) le long du champ de vecteurs
normaux. Supposons que la courbe initiale Γo est simple, régulière et non-convexe, de sorte que κ(s, 0)
change de signe. Supposons que F est deux fois différentiable et que κ l’est aussi pour 0 ≤ s ≤ S et
0 ≤ t ≤ T . Alors, pour tout 0 ≤ t ≤ T
• Si F " (κ) ≤ 0 (resp. ≥ 0) dès que κ = 0 alors
d"
≤ 0 (resp. ≥ 0 )
dt
• Si F " (κ) < 0 (resp. > 0) et κ '= 0 dès que κ = 0 alors
d"
< 0 (resp. > 0 )
dt
5.3.2 La méthode « Level set »

La méthode des lignes de niveaux ( « Level set ») permet de s’affranchir de la paramétrisation
des courbes en terme d’abscisse curviligne. Le prix à payer est une augmentation de la dimen-
sion de l’espace dans lequel on travaille.
L’idée consiste à considérer une courbe plane comme une ligne de niveau d’une surface
3D d’équation z − Φ(x, y) = 0. On choisit le ligne de niveau 0, c’est-à-dire l’intersection de la
surface 3D, avec le plan z = 0.
5.3. LA MÉTHODE DES LIGNES DE NIVEAU (« LEVEL SET ») 79
Plus précisément on va chercher une fonction
Ψ : R3 → R
vérifiant
Γt = { X(t) = (x(t), y(t)) ∈ R2 | Ψ(t, X(t)) = 0 } . (5.3.11)
Pour cela, on va partir d’une courbe initiale Γo d’équation ψo , on pose
Ψ(0, X(0)) = ψo (X(0)) ,
et on va faire évoluer la surface en prenant en compte la vitesse d’évolution du front donnée

dX
par X " (t) = .
dt
L’équation (5.3.11) se traduit par
∀t ≥ 0 Ψ(t, X(t)) = 0 .
Comme nous cherchons une surface régulière (au moins C 1 ) nous pouvons dériver et la formule
de composition donne
∂Ψ
∀t ≥ 0 (t, X(t)) + ∇Ψ(t, X(t)) · X " (t) = 0 ,
∂t
où ∇Ψ désigne le gradient par rapport à X.
Or nous connaissons la vitesse de propagation du front :
X " (t) · &n = F .

Le vecteur normal est donné par

∇Ψ
&n = (X(t)) .
|∇Ψ|
On obtient donc une équation d’évolution pour Ψ :

 ∂Ψ (t, X) + F (t, X) |∇Ψ(t, X)| = 0 ,
∂t (5.3.12)

Ψ(0, X) donnée
Une fois l’équation résolue ( ! !) on obtient la courbe Γt en résolvant
Ψ(t, X) = 0 .
Pour certaines valeurs de F , l’équation (5.3.12) est une équation d’Hamilton-Jacobi :
∂Ψ
+ H(∇Ψ) = 0 ,
∂t
où H est le hamiltonien.
5.3.3 Application à l’image

Reprenons l’exemple de la section 5.2.5. L’équation (5.2.9) donne la vitesse d’évolution du
front :
∂u ∂2u
F = = α 2 − ∇v P (u) .
∂t ∂s
∂2u
∂s2
est à une constante près la courbure κ au point considéré et P (u) peut se réécrire en fonction
de Ψ : P̃ (Ψ). L’équation issue de la méthode des lignes de niveau est donc la suivante
* +
∂Ψ ∇Ψ
+ k div |∇Ψ| − P̃ (Ψ)|∇Ψ| = 0 .
∂t |∇Ψ|
5.4 Le modèle des ballons ( « Balloons »)

Le modèle des ballons ( « Balloons » ) a été introduit par Laurent Cohen [Coh]. En effet
un des principaux problèmes des « snakes » provient de la condition initiale. En effet, si la
condition initiale n’est pas assez proche de la solution, le contour n’évolue pas suffisamment et
tend à se localiser sur un minimum local non significatif. L’intérêt du modèle des ballons réside
dans la résolution de ce problème.
On ajoute une force supplémentaire que l’on peut appeler « force de gonflage ». La courbe est
assimilée à un ballon que l’on gonfle. Deux possibilités sont alors envisageables :
1. Soit la nuance (dans le cas des intensités) est assez forte et la courbe s’arrête.
2. Soit la nuance est trop faible et la courbe la surmonte pour aller chercher plus loin.
5.5. LE MODÈLE DE MUNFORD-SHAH 81
Grâce à ce modèle de ballons, ,on peut supprimer deux des inconvénients principaux des
« snakes » :
1. L’arrêt prématuré de la courbe sur un point non désiré.
2. Le choix d’une condition initiale très proche du contour à extraire.
Ceci est mis en évidence par la figure suivante.
On constate qu’avec le modèle classique (en haut), le contour s’arrête sur un point qui corres-
pond à un bruit, point non-significatif. Avec le modèle des ballons (en bas) , la courbe passe par
dessus et se colle parfaitement au contour même si la condition initiale est très loin du résultat
escompté.
Partant d’une courbe initiale orientée, on ajoute au modèle une force de gonflage définie
par : k1&n(v(s)), où &n est la normale unitaire extérieure à la courbe au point v(s). On aboutit
alors à l’expression suivante de la force extérieure F :
∇Ψ
F = k1&n(v(s)) − k .
|∇Ψ|
Si l’on change le signe de k1 , cela nous donne un dégonflage au lieu d’un gonflage. Notons
aussi que k1 et k sont approximativement du même ordre. Le paramètre k est cependant un
peu plus grand que k1 pour que le « snake » soit arrêté par les bons points.
5.5 Le modèle de Munford-Shah

La méthode des contours actifs précédente a le défaut majeur d’être dépendante de la
paramétrisation de la courbe (via l’abscisse curviligne). Une alternative est d’introduire des
contours actifs géométriques basés sur les propriétés géométriques intrinsèques de la courbe.
Une autre est de formuler le problème sans paramétrer la courbe mais en gardant une prin-
cipe de minimisation d’énergie. C’est le modèle de Munford-Shah.
On se donne un domaine Ω ouvert borné de R2 et une fonction uo : Ω → [0, 1] (on normalise)
qui représente les niveaux de gris de l’image à segmenter. On va chercher les contours sous
la forme d’un ensemble compact de Ω̄, K, reconstruit à partir des discontinuités de uo , ainsi
qu’une approximation régulière de uo en dehors de K qu’on appellera u. On cherche donc une

paire (K, u) qui va minimiser la fonctionnelle suivante :
& &
JM S (K, u) = |∇u|2 dx + α |u − uo |2 dx + β "(K) . (5.5.13)
Ω\K Ω\K
Remarque 5.5.1 (i) Le modèle original de Mumford-Shah fait intervenir la mesure de Hausdorff de K
et non sa longueur "(K) . &
(ii) On reconnaı̂t dans JM S (K, u) : un terme d’ajustement aux données |u − uo |2 dx, un terme
& Ω\K
2
de régularisation |∇u| dx pour u et un terme de régularisation pour les courbes constituant K :
Ω\K
"(K).
(iii) α et β sont des paramètres d’échelle et de contraste respectivement.
L’existence de solution(s) au problème
min JM S (K, u) ,
est un problème difficile. Il est encore largement ouvert.

Annexe A
Quelques outils mathématiques pour

l’image
A.1 Optimisation dans les espaces de Banach

Sauf mention du contraire, on considère dans toute la section un espace de Banach réflexif
V de dual (topologique) V " . On note 1 1V la norme de V et :·, ·; le crochet de dualité entre V et
V ".
A.1.1 Semi-continuité et convexité de fonctionnelles sur V

Définition A.1.1 Une fonction J de V dans R ∪ {+∞} est semi-continue inférieurement (sci) sur V
si elle satisfait aux conditions équivalentes :
– ∀a ∈ R, { u ∈ V | J(u) ≤ a } est fermé
– ∀ū ∈ V, lim inf J(u) ≥ J(ū) .
u→ū
Théorème A.1.1 Toute fonction convexe sci pour la topologie forte (celle de la norme) de V est encore
sci pour la topologie faible de V . !
En pratique ce résultat s’utilise sous la forme du corollaire suivant :
Corollaire A.1.1 Soit J une fonctionnelle convexe de V dans R ∪ {+∞} sci (par exemple continue)
pour la topologie forte. Si vn est une suite de V faiblement convergente vers v alors
J(v) ≤ lim inf J(vn ) .
n→+∞
A.1.2 Gâteaux-différentiabilité des fonctionnelles convexes

Définition A.1.2 Soit J une fonctionnelle de V dans R∪{+∞}. On dit que J est Gâteaux-différentiable
en u ∈ dom (J) si la dérivée directionnelle
J(u + tv) − J(u)
J " (u; v) = lim ,
t→0+ t
83
84 ANNEXE A. QUELQUES OUTILS MATHÉMATIQUES POUR L’IMAGE
existe dans toute direction v de V et si l’application
v 5→ J " (u; v)
est linéaire continue.
De manière générale on notera ∇J(u) la Gâteaux-différentielle de J en u. C’est un élément du

dual V " .
Si V est un espace de Hilbert, avec le théorème de représentation de Riesz (voir [Br] par
exemple) on identifie V et son dual ; on note alors
J " (u; v) = (∇J(u), v) ,
où ( , ) désigne le produit scalaire de V . L’élément ∇J(u) de V est le gradient de J en u.

Il est clair que si J est différentiable au sens classique en u (on dit alors Fréchet - différentiable),
alors J est Gâteaux-différentiable en u, et la dérivée classique et la dérivée au sens de Gâteaux
coı̈ncident.
La réciproque est fausse comme le montre le contre-exemple suivant : soit f de R2 dans R
définie par :
1
y si x = y 2 ,
f (x, y) =
0 sinon
La fonction f est continue en (0,0) et Gâteaux-différentiable en (0,0) mais pas Fréchet - différentia-
ble en (0,0).
Théorème A.1.2 Soit J : C ⊂ V → R , Gâteaux différentiable sur C, avec C convexe. J est convexe si
et seulement si
∀(u, v) ∈ C × C J(v) ≥ J(u) + :∇J(u), v − u; (A.1.1)
Démonstration - Supposons J convexe. Soient u et v dans C. Par convexité de J on a
∀t ∈]0, 1[ J(u + t(v − u)) − J(u) ≤ t(J(v) − J(u)) .
En divisant par t > 0 et en passant à la limite lorsque t → 0+ on obtient (A.1.1).

Réciproquement : on applique (A.1.1) à u + t(v − u) (t ∈ [0, 1]) et u, puis à u + t(v − u) et v pour
obtenir
J(u) ≥ J(u + t(v − u)) − t :∇J(u + t(v − u)), v − u; et
J(v) ≥ J(u + t(v − u)) + (1 − t) :∇J(u + t(v − u)), v − u; .

En faisant la combinaison convexe de ces deux inégalités on obtient
(1 − t)J(u) + tJ(v) ≥ (1 − t + t)J(u + t(v − u)) ,
et la convexité de J. !
A.1. OPTIMISATION DANS LES ESPACES DE BANACH 85
Théorème A.1.3 Soit J : C ⊂ V → R , Gâteaux différentiable sur C, avec C convexe. J est convexe si
et seulement si ∇J est un opérateur monotone , c’est-à-dire
∀(u, v) ∈ C × C :∇J(u) − ∇J(v), u − v; ≥ 0 . (A.1.2)
Démonstration - Soient (u, v) dans C × C. D’après le théorème précédent, si J est convexe, alors
J(v) ≥ J(u) + :∇J(u), v − u;
et
J(u) ≥ J(v) + :∇J(v), u − v; .
En sommant on obtient (A.1.2).

Réciproquement : soient (u, v) dans C × C, u '= v ; on définit ϕ de [0, 1] dans R de la manière
suivante :
ϕ : t 5→ ϕ(t) = (1 − t)J(u) + tJ(v) − J(u + t(v − u)) .
Il est facile de voir que ϕ est dérivable et que
∀t1 , t2 ∈ [0, 1] (ϕ" (t1 ) − ϕ" (t2 ))(t1 − t2 ) ≤ 0 ,
grâce à (A.1.2). Donc ϕ" est décroissante sur [0, 1]. De plus ϕ(0) = ϕ(1) = 0 et d’après le
théorème de Rolle, il existe a ∈]0, 1[ tel que ϕ" (a) = 0. On a donc le tableau suivant qui montre
que ϕ ≥ 0 sur [0, 1]. La convexité de J s’en déduit.
t 0 a 1
ϕ" = 0 =
ϕ" + 0 -
ϕ 0> =0 !
Remarque A.1.1 Supposons que D soit un opérateur strictement monotone :
∀(u, v) ∈ C × C, u '= v, :∇J(u) − ∇J(v), u − v; > 0 . (A.1.3)
alors J est strictement convexe. En effet, on peut reprendre la deuxième partie de la démonstration
du théorème précédent : ϕ" est strictement décroissante sur [0, 1] et donc ne s’annule qu’en un point
a au plus. ϕ est alors strictement croissante sur ] 0, a[ et strictement décroissante sur ] a, 1[ donc
strictement positive sur ] 0, 1[.
On définit de manière analogue la (Gâteaux) dérivée seconde de J en u, comme étant la

dérivée de la fonction (vectorielle) u → ∇J(u). On la note D2 J(u) et on l’appellera aussi Hes-
sien par abus de langage ; ce Hessien est identifiable à une matrice carrée n × n lorsque V = Rn .
A.1.3 Minimisation dans un Banach réflexif

Commençons par un résultat très général de minimisation d’une fonctionnelle semi-continue
sur un ensemble fermé de V .
Définition A.1.3 On dit que J : V → R est coercive si
lim J(x) = +∞ .
*x*V →+∞
Théorème A.1.4 On suppose que V est un Banach réflexif. Soit J une fonctionnelle de V dans R, semi-
continue inférieurement pour la topologie faible de V . Soit K un sous-ensemble non vide et faiblement
fermé de V . On suppose que J est propre ( c’est-à-dire qu’il existe un élément vo de K tel que J(vo ) <
+∞). Alors le problème de minimisation suivant :
1
Trouver u tel que
(P) (A.1.4)
J(u) = inf { J(v) | v ∈ K } ,
admet au moins une solution dans l’un des cas suivants :

– soit J est coercive i.e. lim J(v) = +∞,
*v*V →+∞
– soit K est borné.
Démonstration - On pose d = inf { J(v) | v ∈ K } ; d < +∞, sinon J serait identiquement égale
à +∞ sur K.
Soit une suite minimisante, c’est-à-dire une suite un de K telle que J(un ) → d. Montrons
que cette suite est bornée dans V . Si K est borné, c’est clair. Sinon, on suppose que J est
coercive. Si un n’est pas bornée, on peut en extraire une sous-suite encore notée un telle que
lim 1un 1V = +∞. La coercivité de J implique alors que J(un ) converge vers +∞ ce qui est
n→+∞
en contradiction avec le fait que d < +∞.
V est un Banach réflexif, donc sa boule unité est faiblement compacte ; on peut donc extraire
de la suite un une sous-suite encore notée un qui converge faiblement vers u dans V . On va
montrer que u est solution de (A.1.4).
Comme K est un ensemble faiblement fermé , la limite faible u de la suite un de K est bien un
élément de K.
D’autre part, on sait que J est sci faible ; on a donc
un 4 u dans V ⇒ J(u) ≤ lim inf J(un ) ,

n→+∞
ce qui donne
J(u) ≤ lim inf J(un ) = lim J(un ) = d ≤ J(u) .
n→+∞ n→+∞
On a donc J(u) = d et u ∈ K : u est bien solution. !

Un corollaire important concerne le cas convexe.
A.1. OPTIMISATION DANS LES ESPACES DE BANACH 87
Corollaire A.1.2 On suppose que V est un Banach réflexif. Soit J une fonctionnelle de V dans R,
convexe semi-continue inférieurement V propre et K un sous-ensemble convexe non vide et fermé de V .
Si J est coercive ou si K est borné , le problème de minimisation admet une solution. Si, de plus, J est
strictement convexe la solution est unique.
Démonstration - J étant convexe, sci elle est en particulier faiblement sci. De la même façon K
étant convexe fermé, il est faiblement fermé (voir [Br] par exemple) et le théorème précédent
s’applique.
Montrons à présent l’unicité lorsque J est strictement convexe ; soient u et v deux solutions
u+v
distinctes de (A.1.4). u et v sont dans le convexe K donc aussi. Par conséquent
2
u+v J(u) + J(v)
d ≤ J( )< =d.
2 2
Il y a contradiction, et de fait la solution est unique. !
Rappelons enfin le Théorème donnant une condition nécessaire d’optimalité du premier ordre.
Théorème A.1.5 Soient K un sous-ensemble convexe, non vide de V et J une fonctionnelle de K vers
R Gâteaux-différentiable sur K. Soit u dans V une solution du problème (P). Alors
∀v ∈ K, < ∇J(u), v − u >≥ 0. (A.1.5)
Démonstration - u est dans K ; donc
∀t ∈]0, 1[, ∀v ∈ K, u + t(v − u) ∈ K .
D’après (i),
∀t ∈]0, 1[, ∀v ∈ K J(u) ≤ J((1 − t)u + tv) ,
J((1 − t)u + tv) − J(u)
∀t ∈]0, 1[, ∀v ∈ K ≥0.
t
Par passage à la limite quand t → 0+ :
∀v ∈ K < ∇J(u), v − u >≥ 0 .
!
Corollaire A.1.3 Si J est convexe, sous les hypothèses du théorème précédent, la condition nécessaire
(A.1.5) est également suffisante si u ∈ K.
Démonstration - Soit u est dans K tel que
∀v ∈ K < ∇J(u), v − u >≥ 0 .
Or nous savons que si J est convexe ( Théorème A.1.2 ) < ∇J(u), v − u >≤ J(v) − J(u) , ce qui
permet de conclure que
∀v ∈ K J(v) − J(u) ≥ 0 .
!
A.1.4 Exemple : Projection sur un convexe fermé

Dans ce qui suit V est un espace de Hilbert muni d’un produit scalaire (·, ·) et de la norme
associée 1 ·1 et C est un sous-ensemble convexe et fermé de V .
Théorème A.1.6 Etant donnés C un sous-ensemble convexe, fermé et non vide de V et x un élément
quelconque de V . Alors le problème
min { 1x − y12 , y ∈ C }
a une solution unique x∗ ∈ C. De plus x∗ ∈ C est caractérisé par :
∀y ∈ C (x − x∗ , y − x∗ ) ≤ 0 . (A.1.6)
Démonstration - La démonstration est immédiate. Nous avons affaire à un problème de moindres

carrés. La fonction coût est continue, coercive et strictement convexe et l’ensemble C est convexe
fermé. On peut donc appliquer le corollaire A.1.2. La caractérisation de x∗ s’obtient par appli-
cation du théorème A.1.5 et de son corollaire. !
Nous avons une seconde caractérisation de x∗ de manière immédiate :
Corollaire A.1.4 Sous les hypothèses du théorème (A.1.6) on peut caractériser le projeté x∗ ∈ C de x
par :
∀y ∈ C (x∗ − y, y − x) ≤ 0 . (A.1.7)
Démonstration - Si x∗ est le projeté de x sur C, (A.1.6) donne :
∀y ∈ C (x − x∗ , y − x∗ ) ≤ 0 .
Donc
∀y ∈ C (x∗ − y, y − x) = (x∗ − y, x∗ − x) − 1y − x∗ 12 ≤ 0 .
Réciproquement : Soit y ∈ C et z = x∗ + t(y − x∗ ) ∈ C, pour t ∈]0, 1[. La relation (A.1.7) implique
∀t ∈] 0, 1[ (x∗ − z, z − x) = −t (y − x∗ , x∗ − x + t(y − x∗ )) ≤ 0
∀t ∈] 0, 1[ (y − x∗ , x∗ − x + t(y − x∗ )) ≥ 0 .
On fait ensuite tendre t vers 0+ pour obtenir (A.1.6). !
Le point x est le projeté de x sur C. L’application πC = V → C qui à x associe son projeté
∗
x∗ est la projection sur C. Le projeté πC (x) est donc le point de C qui est le « plus près » de x.
On définit de manière standard la fonction distance d’un point x à l’ensemble C par
d(x, C) = inf 1x − y1 . (A.1.8)

y∈C
Dans le cas où C est un convexe fermé, on vient donc de démontrer que
d(x, C) = 1x − πC (x)1 .
A.2. FORMULATION VARIATIONNELLE DES ÉQUATIONS AUX DÉRIVÉES PARTIELLES 89
Proposition A.1.1 La projection πC est continue de V dans C. Plus précisément on a
∀(x, y) ∈ V × V 1πC (x) − πC (y)1 ≤ 1x − y1 ,
c’est-à-dire πC est une contraction.
Démonstration - Soient x1 et x2 deux éléments quelconques de V . Appliquons la relation (A.1.6)

à x = x1 , x∗ = πC (x1 ) et y = πC (x2 ) ∈ C puis à x = x2 , x∗ = πC (x2 ) et y = πC (x1 ) ∈ C : on
obtient
(x1 − πC (x1 ), πC (x2 ) − πC (x1 )) ≤ 0
(x2 − πC (x2 ), πC (x1 ) − πC (x2 )) ≤ 0 ;
La somme des deux inégalités donne
(x1 − x2 , πC (x2 ) − πC (x1 )) + 1πC (x2 ) − πC (x1 )12 ≤ 0 ,
c’est-à-dire avec l’inégalité de Cauchy-Schwarzrtz
1πC (x2 ) − πC (x1 )12 ≤ (x2 − x1 , πC (x2 ) − πC (x1 )) ≤ 1x2 − x1 1 1πC (x2 ) − πC (x1 )1 .
Si πC (x2 ) = πC (x1 ) la relation que l’on cherche est évidente.

Sinon on divise par πC (x2 ) − πC (x1 ) et on obtient le résultat souhaité. !
Remarque A.1.2 1. Si x ∈ C alors πC (x) = x. Plus généralement si C = V alors πC = IdV .

2. Le théorème A.1.6 est faux si C n’est pas convexe ou si C n’est pas fermé.
3. La projection πC n’est pas différentiable en général, mais l’application x 5→ 1x − πC (x)12 l’est.
Exemple A.1.1 (Projection sur un sous-espace vectoriel) Dans le cas où C est un sous-espace vec-
toriel fermé de V , c’est bien sûr un convexe fermé non vide. L’opérateur de projection est dans ce cas
linéaire (c’est faux dans le cas général). Le projeté x∗ d’un élément x, sur C, est caractérisé par
∀y ∈ C (x − x∗ , y) = 0 .
Cela signifie que x − x∗ ∈ C ⊥ (l’orthogonal de C). On retrouve ainsi la classique projection orthogonale
sur un sous-espace vectoriel fermé.
A.2 Formulation variationnelle des équations aux dérivées partielles

A.2.1 Théorème de Lax-Milgram
Une application très importante des résultats de la section précédente est le théorème de
Lax-Milgram :
Théorème A.2.1 Soit V un espace de Hilbert et a une forme bilinéaire continue sur V ×V . On suppose
que a est V - elliptique (ou coercive), i.e. :
∃α > 0, ∀v ∈ V a(v, v) ≥ α1v12V . (A.2.9)
Soit L ∈ V " . Alors le problème : 1

Trouver u ∈ V tel que
(A.2.10)
∀v ∈ V a(u, v) = L(v)
a une solution unique.
Démonstration - On pose
1
J(v) = a(v, v) − L(v) .
∀v ∈ V,
2
Il est facile de voir que J est convexe à cause de la bilinéarité de a et de la linéarité de L.
J est strictement convexe grâce à la V -ellipticité de a. La V -ellipticité de a entraı̂ne aussi la
coercivité de J ; en effet, nous avons :
J(v) ≥ α1v12V − L(v) ≥ α1v12V − 1L1V ! 1v1V .
Enfin il est clair que J est continue pour la topologie forte de V , donc sci. De plus comme elle
est convexe, elle est aussi sci pour la topologie faible.
Par conséquent, d’après le théorème A.1.4 le problème
1
min { a(v, v) − L(v) | v ∈ V } , (A.2.11)
2
admet une solution unique u.
Grâce au Théorème A.1.5 on peut maintenant caractériser cette solution u dans le cas où a est
symétrique. Un calcul élémentaire montre que
∀v ∈ V ∇J(u).(v − u) = a(u, v − u) − L(v − u) .
Donc la solution u est caractérisée par :
∀v ∈ V a(u, v − u) − L(v − u) ≥ 0 ,
c’est-à-dire
∀v ∈ V a(u, v) = L(v) .
Lorsque a n’est pas symétrique le résultat est encore vrai mais la démonstration est moins
immédiate (elle est en tout cas différente). Nous renvoyons à [DL] Tome VI, p 1206. !
Nous allons interpréter ce résultat. Soit u dans V et Au ∈ V " défini par :
V →R
Au : v 5→ a(u, v) .
Soit maintenant l’opérateur A de V dans V " défini par
∀u ∈ V Au = Au .
Il est facile de voir que A est linéaire et que le problème (A.2.10) est équivalent à :
1
Trouver u ∈ V tel que
(A.2.12)
Au = L .
On a alors la proposition suivante, conséquence immédiate du théorème de Lax-Milgram :
Proposition A.2.1 L’opérateur A est un isomorphisme (algébrique et topologique) de V dans V " .
Démonstration - A est linéaire. Le théorème de Lax-Milgram prouve que A est un isomoprhisme

algébrique de V dans V " . La continuité de A provient de celle de a ; plus précisément, pour tout
u de V :
| :Au, v;V ! ,V | |a(u, v)|
1Au1V ! = sup = sup ≤ 1a1 1u1V .
v,=0 1v1 V v,=0 1v1V
:u, v;V,V ! désigne le crochet de dualité entre V et V " .
On conclut que A−1 est aussi continu grâce au théorème de l’application ouverte. !
Cas particulier important

Soient V et H deux espaces de Hilbert tels que V est inclus dans H avec injection continue.
On suppose également que V est dense dans H. On peut identifier H et H " par le théorème de
Riesz et on a grâce à la densité :
V ⊂H ⊂V" , (A.2.13)
avec injections continues et denses. H s’appelle l’espace pivot.
A.2.2 Généralités sur les espaces de Sobolev

Pour plus de précisions on pourra se référer à [DL].
Soit Ω un ouvert borné de Rn , (n ≤ 3 en pratique) de frontière régulière Γ. On appelle D(Ω) l’es-
pace des fonctions C ∞ à support compact dans Ω. Son dual D" (Ω) est l’espace des distributions
sur Ω.
∂u
Pour toute distribution u ∈ D" (Ω), la dérivée est définie de la manière suivante :
∂xi
L M L M
∂u def ∂ϕ
∀ϕ ∈ D(Ω) ,ϕ ≡ − u, .
∂xi D! (Ω),D(Ω) ∂xi D! (Ω),D(Ω)
∂u
On notera indifféremment la dérivée de u au sens des distributions Di u = = ∂i u.
∂xi
si α ∈ Nn , on note Dα u = ∂1α1 u · · · ∂nαn u et |α| = α1 + · · · + αn ; on obtient
∀ϕ ∈ D(Ω) :Dα u, ϕ;D! (Ω),D(Ω) = (−1)|α| :u, Dα ϕ;D! (Ω),D(Ω) .

Définition A.2.1 On définit les espaces de Sobolev H m (Ω) de la manière suivante :

∂u
H 1 (Ω) = { u ∈ L2 (Ω) | ∈ L2 (Ω), i = 1 · · · n } ,
∂xi
H m (Ω) = { u ∈ D" (Ω) | Dα u ∈ L2 (Ω), |α| ≤ m } .
Remarque A.2.1 H o (Ω) = L2 (Ω).
Nous allons énoncer une série de propriétés des espaces de Sobolev, sans démonstration. On
pourra consulter [DL] par exemple.
Proposition A.2.2 H m (Ω) muni du produit scalaire :

! &
(u, v)m = Dα u(x) Dα v(x) dx ,
|α|≤m Ω
est un espace de Hilbert.
Proposition A.2.3
!
H m (Ω) ⊂ H m (Ω)
et l’injection est continue, pour m ≥ m" .
Définition A.2.2
Ho1 (Ω) = { u ∈ H 1 (Ω) | u|Γ = 0 } .
C’est aussi l’adhérence de D(Ω) dans H 1 (Ω).
∂j u
Hom (Ω) = { u ∈ H 1 (Ω) | = 0, j = 1, · · · , m − 1} ,
∂nj |Γ
∂
où est la dérivée de u suivant la normale extérieure à Γ la frontière de Ω :
∂n
n
∂u ! ∂u
= cos(&n, e&i ) ,
∂n ∂xi
i=1
où &n est la normale extérieure à Γ et Ω est supposé « régulier » (de frontière C ∞ par exemple).
Définition A.2.3 (Dualité) Pour tout m ∈ N, on note H −m (Ω) le dual de Hom (Ω).
Théorème A.2.2 (Rellich) Si Ω est un ouvert borné de Rn , alors pour tout m ∈ N, l’injection de
Hom+1 (Ω) dans Hom (Ω) est compacte .
En particulier l’injection de Ho1 (Ω) dans L2 (Ω) est compacte. En pratique, cela signifie que
toute suite bornée en norme Ho1 (Ω) converge faiblement dans Ho1 (Ω) (après extraction d’une
sous-suite) et fortement dans L2 (Ω).
A.2.3 Problème homogène de Dirichlet : cas elliptique

Dans cet exemple on choisit V = Ho1 (Ω) muni du produit scalaire :
& !
n
(u, v)1 = [(Di u Di v) + uv] dx .
Ω i=1
On se donne la forme bilinéaire suivante

& ! n &
def
a(u, v) ≡ (Di u Di v) dx = ∇u(x) ∇v(x) dx .
Ω i=1 Ω
a est bilinéaire de V × V dans R ; elle est continue car :
|a(u, v)| ≤ 1u1V 1v1V ,
de manière évidente.
a est V -elliptique d’après l’inégalité de Poincaré
∀u ∈ Ho1 (Ω) 1u12L2 (Ω) ≤ C1∇u12L2 (Ω) . (A.2.14)
En effet on obtient alors
(1 + C)1∇u12L2 (Ω) ≥ 1∇u12L2 (Ω) + 1u12L2 (Ω) = 1u121 ,
c’est-à-dire
1
a(u, u) ≥ 1u121 .
1+C
&
Soit f ∈ H −1 (Ω) et L(v) ≡ f (x) v(x) dx pour tout v de Ho1 (Ω).
Ω
On sait que le problème variationnel
& &
∀v ∈ Ho (Ω)
1
∇u(x) ∇v(x) dx = f (x) v(x) dx ,
Ω Ω
admet une solution unique u ∈ Ho1 (Ω).
Quelle est l’interprétation de ce problème variationnel ?
Soit ϕ ∈ D(Ω) ⊂ Ho1 (Ω) ;

&
:−∆u, ϕ;D! (Ω),D(Ω) = ∇u(x) ∇ϕ(x) dx = :f, ϕ;D! (Ω),D(Ω) .
Ω
On obtient
−∆u = f dans D" (Ω) .
Comme u ∈ Ho1 (Ω), on a également : u ≡ 0 sur Γ. Ce problème est le problème de Dirichlet
homogène pour le laplacien :
1
−∆u = f dans Ω
(DH)
u = 0 sur Γ ,
et on a le résultat suivant :
Théorème A.2.3 L’opérateur −∆ (laplacien) est un isomorphisme de Ho1 (Ω) sur H −1 (Ω).
!
Par ailleurs Ho1 (Ω) ⊂ L2 (Ω) avec injection dense et continue. Donc
Ho1 (Ω) ⊂ L2 (Ω) ⊂ H −1 (Ω) .
On peut considérer (−∆) comme un opérateur non borné à valeurs dans H = L2 (Ω). Le do-
maine de (−∆) est alors
D(−∆) = { u ∈ Ho1 (Ω) | ∆u ∈ L2 (Ω) } .
On montre que si Ω est de frontière assez régulière
D(−∆) = H 2 (Ω) ∩ Ho1 (Ω) ,
et on même le résultat suivant :
Proposition A.2.4 L’opérateur −∆ est un isomorphisme de H 2 (Ω) ∩ Ho1 (Ω) muni de la norme du
graphe sur L2 (Ω).
!
A.2.4 Problème de Neumann homogène

On choisit maintenant V = H 1 (Ω) et la forme bilinéaire a est donnée par :
&
def
a(u, v) ≡ (u, v)1 = [∇u(x) ∇v(x) + u(x) v(x)] dx .
Ω
Soient f ∈ L2 (Ω) et g ∈ L2 (Γ) ;

on définit la forme linéaire L ∈ V " par
& &
L(v) = f (x) v(x) dx + g(x) v(x) dx .
Ω Γ
Le théorème de Lax-Milgram assure l’existence et l’unicité de la solution du problème varia-
tionnel
& & &
∀v ∈ H (Ω)
1
[∇u(x) ∇v(x) + u(x) v(x)] dx = f (x) v(x) dx + g(x) v(x) dx .
Ω Ω Γ
Ici aussi on peut interpréter cette équation variationnelle comme le problème de Neumann
pour l’opérateur (−∆ + I) :
9
−∆u + u = f dans D" (Ω)
(N ) ∂u
= g sur Γ ,
∂n
En conclusion, quand on parlera de solutions, ce sera toujours de solutions faibles , i.e. au
sens des distributions.
A.3. ANALYSE CONVEXE NON LISSE 95
A.3 Analyse convexe non lisse

A.3.1 Théorème de Hahn -Banach
Dans ce qui suit X est un espace de Banach réel de dual X " .
Le Théorème de Hahn-Banach, sous sa forme géométrique, permet de séparer des en-
sembles convexes. Il est très important en analyse convexe et sert en particulier à exhiber des
multiplicateurs en optimisation. Nous rappelons ici la forme géométrique de ce théorème qui
est la seule utile dans notre cas ainsi que des corollaires importants. Pour les démonstrations et
plus de précisions nous renvoyons au livre de H . BR ÉZIS [Br].
Définition A.3.1 (Hyperplan affine) Un hyperplan affine fermé est un ensemble de la forme
H = { x ∈ X | α(x) + β = 0 } ,
où α ∈ X " est une forme linéaire continue non nulle sur X et β ∈ R.
Dans le cas où X est un espace de Hilbert V (en particulier si V = Rn ), on peut identifier V
à son dual et tout hyperplan affine fermé est de la forme
H = { x ∈ V | (α, x)V + β = 0 } ,
où (·, ·)V désigne le produit scalaire de V , α ∈ V, α '= 0 et β ∈ R.
Définition A.3.2 (Séparation) Soient A et B deux sous-ensembles non vides de X . On dit que l’hy-
perplan affine H d’équation : α(x) + β = 0, sépare A et B au sens large si
∀x ∈ A α(x) + β ≤ 0 et ∀y ∈ B α(y) + β ≥ 0 .
On dit que H sépare A et B strictement s’il existe ε > 0 tel que
∀x ∈ A α(x) + β ≤ −ε et ∀y ∈ B α(y) + β ≥ ε .
Donnons à présent la première forme géométrique du théorème de Hahn-Banach :
Théorème A.3.1 Soient A et B deux sous-ensembles de X convexes, non vides et disjoints. On suppose
que A est ouvert. Alors, il existe un hyperplan affine fermé qui sépare A et B au sens large .
Corollaire A.3.1 Soit C un convexe non vide de Rn et fermé et x∗ ∈ C. Alors :

x∗ ∈ Int (C) si et seulement si il n’existe aucune forme linéaire séparant x∗ et C.
Citons ausi la deuxième forme géométrique du théorème de Hahn-Banach :
Théorème A.3.2 Soient A et B deux sous-ensembles de X convexes, non vides et disjoints. On suppose
que A est fermé et que B est compact. Alors, il existe un hyperplan affine fermé qui sépare A et B
strictement.
A.3.2 Sous-différentiel
Définition A.3.3 Soit f : X → R ∪ {+∞} et u ∈ dom f (i.e. f (u) < +∞). Le sous-différentiel de f
en u est l’ensemble ∂f (u) (éventuellement vide ) des u∗ ∈ X " tels que
∀v ∈ X f (v) ≥ f (u) + :u∗ , v − u; .
Les éléments u∗ sont appelés sous-gradients.
Remarque A.3.1 1. f : X → R ∪ {+∞} atteint son minimum en u ∈ dom f si et seulement si
0 ∈ ∂f (u).
2. Si f, g : X → R ∪ {+∞} et u ∈ dom f ∩ dom g, on a
∂f (u) + ∂g(u) ⊂ ∂(f + g)(u) .
3. Comme N
∂f (u) = { u∗ ∈ X " | :u∗ , v − u; ≤ f (v) − f (u) } ,
v∈X
∂f (u) est un sous-ensemble convexe, fermé pour la topologie faible *, comme intersection de convexes
fermés.
4. Pour tout λ > 0 on a ∂(λf )(u) = λ∂f (u) .
Lien avec la Gâteaux-différentiabilité
Théorème A.3.3 Soit f : X → R ∪ {+∞} convexe.

Si f est Gâteaux-différentiable en u ∈ dom f , elle est sous-différentiable et ∂f (u) = {f " (u)} .
Réciproquement, si f est finie, continue en u et ne posssède qu’un seul sous-gradient, alors f est
Gâteaux-différentiable en u et ∂f (u) = {f " (u)} .
Démonstration - Supposons que f est Gâteaux-différentiable en u ∈ dom f de Gâteaux -dérivée

f " (u). D’après une propriété standard des fonctions convexes (voir théorème A.1.2) on a, pour
tout v ∈ X
f (u + (v − u)) − f (u) ≥< f " (u), v − u > ,
f (v) ≥ f (u)+ < f " (u), v − u > .

Par conséquent f " (u) ∈ ∂f (u). Par ailleurs, si u∗ ∈ ∂f (u) on aura pour tous t > 0 et w ∈ X
f (u + tw) − f (u) ≥ t < u∗ , w > ,
soit en divisant par t et en passant à la limite
< f " (u), w >≥< u∗ , w > ,

c’est-à-dire f " (u) = u∗ .

Réciproquement : on suppose que f est finie, continue en u et ne posssède qu’un seul sous-
gradient u∗ . Soit v ∈ X . La fonction f est convexe, donc elle admet une dérivée directionnelle
f (u + tv) − f (u)
f " (u; v) = lim ,
t→0+ t
f (u + tv) − f (u)
éventuellement infinie. En effet, l’application t 5→ est une fonction croissante
t
de t sur ]0, +∞[ et donc elle possède une limite quand t → 0+ . Montrons la croissance de ϕ sur
s
]0, +∞[ : soit t ≥ s > 0. Posons λ = ∈ ]0, 1].
t
f (u + sv) − f (u) f (u + λtv) − f (u) f ((1 − λ)u + λ(u + tv)) − f (u)
ϕ(s) = = =
s s s
(1 − λ)f (u) + λf (u + tv) − f (u) f (u + tv) − f (u) f (u + tv) − f (u)
≤ =λ =
s s t
De plus, pour tous s > t > 0
f (u + tv) − f (u) f (u + sv) − f (u)
≤
t s
et par passage à la limite lorsque t → 0+
∀s > 0 sf " (u; v) ≤ f (u + sv) − f ((u) .
De même, pour tous s < 0 < t on a

f (u + tv) − f (u) f (u + sv) − f (u)
≥
t s
c’est-à-dire par passage à la limite
f (u + sv) − f (u)
∀s < 0 f " (u; v) ≥ .
s
En effet si s ≤ t < 0 on a encore
ϕ(s) ≤ ϕ(t) .
Il suffit de faire la démonstration précédente avec −s, −t et −v. Montrons maintenant
∀t > 0 ϕ(−t) ≤ ϕ(t) .
Avec la convexité, on a
* + * +
(u + tv) + (u − tv) f (u + tv) + f (u − tv)
f (u) = f ≤ .
2 2
Par conséquent
f (u + tv) − f (u) ≥ −(f (u − tv) − f (u)) ,
et donc
f (u + tv) − f (u) f (u − tv) − f (u)
≥ .
t −t
Supposons maintenant s < 0 < t avec par exemple s ≤ −t < 0 < t. D’après ce qui précède
ϕ(s) ≤ ϕ(−t) ≤ ϕ(t) .
De même si s < 0 < −s < t ϕ(s) ≤ ϕ(−s) ≤ ϕ(t) . Donc, pour tout t ∈ R
f (u + tv) ≥ f (u) + tf " (u; v) . (A.3.15)
Soit
C = {(v, t) ∈ X × R | f (v) ≤ t }) .
Comme f est convexe, il est facile de voir que C est convexe, non vide (car f est finie et conti-
nue). L’intérieur de C est un convexe ouvert et la relation (A.3.15) montre que la droite
∆ = {(u + tv, f (u) + tf " (u; v)) , t ∈ R }
ne rencontre pas int(C). D’après le théorème de Hahn-Banach, il existe u∗ ∈ X " et s∗ , α ∈ R,

avec (u∗ , s∗ ) '= (0, 0) tels que
∀t ∈ R s∗ (f (u) + tf " (u; v))− < u∗ , u + tv >= α , (A.3.16)
et
∀(w, t)) ∈ int C s∗ t− < u∗ , w >≥ α
ce qui entraı̂ne
∀(w, t) ∈ C s∗ t − < u∗ , w >≥ α . (A.3.17)
Prenons t = 0 dans (A.3.16) et (w, f (w)) ∈ C dans (A.3.17) : on obtient
α = s∗ f (u)− < u∗ , u > et s∗ f (w) ≥< u∗ , w > +α ; (A.3.18)
s∗ ne peut pas être nul : en effet si s∗ = 0, la relation (A.3.17) donne
∀w ∈ X 0 ≥< u∗ , w > +α =< u∗ , w > − < u∗ , u > .
Donc u∗ = 0 et de même α = − < u∗ , u >= 0.

D’autre part si on fait tendre t vers +∞ dans (A.3.17) on voit que s∗ est nécessairement positif.
On peut donc supposer s∗ = 1. La relation (A.3.18) donne alors
∀w ∈ X f (w) − f (u) ≥< u∗ , w − u > .
Cela signifie que u∗ est le sous-gradient de f en u. D’autre part (A.3.16) (avec t = 1) donne
f " (u; v) =< u∗ , v > ;
ceci prouve la Gâteaux-différentiabilité de f en u. !

Sous-différentiel d’une somme de fonctions

Théorème A.3.4 Soient f et g convexes, sci à valeurs dans R∪{+∞}. On suppose qu’il existe uo ∈dom
f ∩ dom g tel que f soit continue en uo . Alors
∀u ∈ X ∂(f + g)(u) = ∂f (u) + ∂g(u) .
Démonstration - La démonstration est laissée en exercice. Elle se fait suivant un schéma main-
tenant classique : si u∗ ∈ ∂(f + g)(u), on définit
C1 = {(v, a) | f (v)− < v − u, u∗ > −f (u) ≤ a } ,
et
C2 = {(v, a) | a ≤ g(u) − g(v) } .
Ce sont deux convexes et on sépare int C1 et C2 grâce au théorème de Hahn- Banach. !
Le résultat suivant s’en déduit mais on peut le montrer directement sans hypothèse de convexité
sur les deux fonctions :
Théorème A.3.5 Supposons que f est Gâteaux-différentiable, ϕ est convexe et que le problème suivant
min f (u) + ϕ(u) ,
u∈X
a (au moins) une solution ū. Alors

−f " (u) ∈ ∂ϕ(ū) .
Démonstration - Soit u ∈ X . On a pour tout t ∈ [0, 1]
f (ū) + ϕ(ū) ≤ f (ū + t(u − ū)) + ϕ(ū + t(u − ū)) ;
donc
− (f (ū + t(u − ū)) − f (ū)) ≤ ϕ(ū + t(u − ū)) − ϕ(ū) ≤ t(ϕ(u) − ϕ(ū))
grâce à la convexité de ϕ. Donc
f (ū + t(u − ū)) − f (ū)
− ≤ ϕ(u) − ϕ(ū) ,
t
et par passage à la limite quand t → 0+
−:f " (u), u − ū; ≤ ϕ(u) − ϕ(ū)
ce qui prouve le résultat. !
Nous admettrons enfin le résultat suivant :
Théorème A.3.6 Soit Λ linéaire continue de V dans Y (espaces de banach). Soit f convexe, sci de V
dans R ∪ {+∞} continue en au moins un point de son domaine (supposé non vide). Alors
∀u ∈ V ∂(f ◦ Λ)(u) = Λ∗ ∂f (Λu) ,
où Λ∗ est l’opérateur adjoint de Λ.
Pour plus de détails sur ces notions là on peut consulter [Azé, ET]. Nous terminons par un
exemple important.
Application à l’indicatrice d’un ensemble

Dans le cas où f est la fonction indicatrice d’un sous-ensemble non vide K de X :
1
def 0 si u ∈ K ,
f (u) = 1K (u) = ,
+∞ sinon
le sous-différentiel de f en u est le cône normal de K en u :
∂1K (u) = NK (u) = { u∗ ∈ X " | :u∗ , v − u; ≤ 0 pour tout v ∈ K } .
Dans le cas où X est un espace de Hilbert identifié à son dual, et K un sous-ensemble fermé,
convexe non vide de X , nous allons préciser le sous-différentiel de 1K en u (c’est-à-dire le cône
normal à K en u ) :
Proposition A.3.1 Soit u ∈ K, où K est un sous-ensemble fermé, convexe non vide de X espace de
Hilbert. Alors
λ λ
λ ∈ ∂1K (u) ⇐⇒ λ = c [u + − ΠK (u + )]
c c
pour tout c réel strictement positif, où ΠK est la projection de X sur le convexe K.
Démonstration - Remarquons tout d’abord que ∂1K (u) est un sous-ensemble de X . On rappelle
également que si ΠK est la projection de X sur le convexe fermé K, l’image ΠK (w) d’un élément
quelconque de X est caractérisée par
∀v ∈ K (w − ΠK (w), v − ΠK (w))X ≤ 0 ,
où (·, ·)X désigne le produit scalaire de X . Soit λ ∈ ∂1K (u) : λ est caractérisé par
∀v ∈ K (λ, v − u)X ≤ 0
c’est-à-dire, pour tout c > 0

* +
λ
∀v ∈ K u + − u,v − u ≤0.
c X
λ
D’après ce qui précède (en posant w = u + )
c
λ λ λ
λ ∈ ∂1K (u) ⇐⇒ u = ΠK (u + ) ⇐⇒ λ = c [u + − ΠK (u + )] .
c c c
!
A.3.3 Calcul « pratique » : transformation de Legendre-Fenchel

Définition A.3.4 Soit f : X → R̄ une fonction. La transformée de Legendre-Fenchel ou conjuguée
de f est la fonction f ∗ : X " → R̄ définie par
∀" ∈ X " f ∗ (") = sup { "(u) − f (u) } . (A.3.19)

u∈X
Remarque A.3.2 (a) Si f « prend » la valeur −∞, alors f ∗ ≡ +∞. Si f est propre (c’est-à-dire non
indentiquement égale à +∞) alors f ∗ est à valeurs dans R ∪ {+∞}.
(a) On notera désormais "(u) =< ", u >, où < ·, · > désigne le produit de dualité entre X et X " .
Ce produit de dualité coı̈ncide avec le produit scalaire dans le cas où X est un espace de Hilbert, après
identification de X et de son dual .
L’équation (A.3.19) s’écrit alors
∀u∗ ∈ X " f ∗ (u∗ ) = sup { :u∗ , u; − f (u) } .

u∈X
Définition A.3.5 Soit A ⊂ X un ensemble (non vide). La fonction d’appui de l’ensemble A est la
fonction σA : X " → R ∪ {+∞} définie par σA = (1A )∗ où 1A désigne l’indicatrice de A
1
0 si x ∈ A ,
1A (x) =
+∞ sinon.
Exemple A.3.1 Soit A un ensemble et f (x) = d(x, A). Alors f ∗ = σA + 1B ∗ .
Proposition A.3.2 Pour toute fonction f : X → R ∪ {+∞}, la fonction f ∗ est convexe et sci pour la
topologie faible *.
Démonstration - Par définition

f ∗ = sup ϕu ,
u∈domf
où dom f est le domaine de f ( i.e. l’ensemble des éléments u ∈ X tels que f (u) est finie) et
ϕu : X " → R est définie par
ϕu (u∗ ) =< u∗ , u > −f (u) .
Chacune des fonctions ϕu est affine et continue, donc convexe et sci pour la topologie faible de
X " . Il en est de même pour le sup. !
Exemple A.3.2 Soit f : u 5→ 1u1X . M f ∗ = 1B ∗ où B est la boule unité de X " .
Plus généralement
Proposition A.3.3 Soit f une fonction positivement homogène (propre) de X dans R ∪ {+∞}. Sa
conjuguée f ∗ est l’indicatrice d’un sous-ensemble K convexe et fermé de X " .
Démonstration - Soit f une fonction positivement homogène (propre) de X dans R ∪ {+∞}. Soit
u∗ ∈ X " . Deux cas se présentent :
• ∃uo ∈ X tel que :u∗ , uo ; − f (uo ) > 0 . Alors par homogénéité, pour tout λ > 0
:u∗ , λuo ; − f (λuo ) = λ[:u∗ , uo ; − f (uo )] ≤ f ∗ (u∗ ) .

Donc, en passant à la limite pour λ → +∞ on obtient f ∗ (u∗ ) = +∞.

• Dans le cas contraire
∀u ∈ X :u∗ , u; − f (u) ≤ 0 ,
et donc f ∗ (u∗ ) ≤ 0. Or par définition de f ∗ ,
:u∗ , 0; − f (0) ≤ f ∗ (u∗ ) ;
comme f est positivement homogène f (0) = f (n ∗ 0) = n ∗ f (0) pour tout n ∈ N et donc

f (0) = 0. On a donc finalement : f ∗ (u∗ ) = 0.
Posons K = {u∗ ∈ X ∗ | f ∗ (u∗ ) = 0 }. On vient de montrer que f ∗ = 1K . Comme f ∗ est convexe,
sci K est fermé , convexe.
!
Nous allons maintenant donner un résultat reliant f + g et f ∗ + g ∗ qui est le fondement de
la théorie de la dualité en analyse convexe :
Théorème A.3.7 Soient f, g : X → R ∪ {+∞} des fonctions convexes telles qu’il existe uo ∈domg
avec f continue en uo . Alors
inf (f (u) + g(u)) = max

∗ !
(−f ∗ (u∗ ) − g ∗ (−u∗ )) .
u∈X u ∈X
Démonstration - Posons
α = inf (f (u) + g(u)) et β = sup (−f ∗ (u∗ ) − g ∗ (−u∗ )) .

u∈X u∗ ∈X !
Soient u ∈ X et u∗ ∈ X " : par définition on a
−f ∗ (u∗ ) ≤ − < u∗ , u > +f (u) et − g ∗ (−u∗ ) ≤< u∗ , u > +g(u) ,
donc
−f ∗ (u∗ ) − g ∗ (−u∗ ) ≤ f (u) + g(u) ;
en passant au sup sur le terme de gauche et à l’inf sur celui de droite, on obtient
β ≤α.
Montrons l’inégalité inverse. Comme uo ∈dom f ∩ dom g, α ∈ R ∪ −∞. Si α = −∞, le théorème

est démontré. On peut donc supposer que α ∈ R. Soient
C = int({(u, t) ∈ X × R | f (u) ≤ t }) ,
et
D = {(u, t) ∈ X × R | t ≤ α − g(u) }) '= ∅ .
Comme f et g sont convexes, C et D sont convexes. Comme f est continue en uo , C est non
vide. De plus C ∩D = ∅. On peut donc appliquer le théorème de Hahn-Banach : on peut trouver
(u∗o , so ) ∈ X " × R\{0, 0} et c ∈ R tels que
∀(v, s) ∈ D < u∗o , v > +sso ≥ c ,

et
∀(w, σ) ∈ C c ≥< u∗o , w > +σso . (A.3.20)
Comme σ peut tendre vers +∞ d’après la définition de C, on obtient so ≤ 0.
Supposons que so '= 0. Dans ce cas, so < 0 et (quitte à diviser tout par |so |) on peut supposer
que so = −1. On obtient
∀(v, s) ∈ D − < u∗o , v > +s ≤ −c .
Soit u ∈ X et s = α − g(u) : le couple (u, s) est dans D. Donc
∀u ∈ X − < u∗o , u > +α − g(u) ≤ −c .
D’autre part l’inégalité (A.3.20) peut s’étendre à C̄ et par convexité
C̄ = {(u, t) ∈ X × R | f (u) ≤ t } ;
on peut donc l’appliquer au couple (u, f (u) pour tout u ∈ X ce qui donne
c ≥< u∗o , u > −f (u) .
Finalement,
g ∗ (−u∗o ) ≤ −c − α et f ∗ (u∗o ) ≤ c .
Donc
α ≤ −f ∗ (u∗o ) − g ∗ (−u∗o ) ≤ β ≤ α
ce qui finit la démonstration.
Si so = 0. Comme f est continue en uo on peut trouver une boule B(uo , R) avec R > 0 incluse
dans dom f . On a alors (uo , α − g(uo )) ∈ D et pour tout w ∈ B(uo , R), (uo + w, f (uo ) + εo ) ∈ C :
d’où
< u∗o , uo + w >≤ c ≤< u∗o , uo > .
Ceci entraı̂ne que u∗o = 0 et une contradiction puisque (u∗o , so ) '= (0, 0). !
Remarque A.3.3 Notons que dans le théorème le « sup » dans le terme de droite est toujours atteint
(c’est un « max » ) ce qui n’est pas le cas dans le terme de gauche. où l’inf n’est pas nécessairement
atteint.
Corollaire A.3.2 Soit f : X → R ∪ {+∞} une fonction convexe continue en u ∈ X . Alors
f (u) = max
∗ !
(< u∗ , u > −f ∗ (u∗ )) .
u ∈X
Démonstration - Posons g = 1{u} . On a g ∗ (u∗ ) =< u∗ , u > pour tout u∗ ∈ X " . Les fonctions f et
g sont convexes et f est continue en u ∈dom g. On a donc d’après le théorème précédent :
f (u) = inf (f + g)(u) = max

∗ !
(−f ∗ (u∗ ) − g ∗ (−u∗ )) = max
∗ !
(< u∗ , u > −f ∗ (u∗ )) .
u∈X u ∈X u ∈X
!
On peut généraliser ce résultat à des fonctions convexes sci.
Théorème A.3.8 Soit f : X → R ∪ {+∞} une fonction convexe semi-continue inférieurement. Alors,
pour tout u ∈ X
f (u) = max
∗ !
(< u∗ , u > −f ∗ (u∗ )) .
u ∈X
Démonstration - Voir [Azé] p. 89. C’est une démonstration similaire à celle du théorème (A.3.7).
!
Terminons par un résultat de bi-dualité que nous admettrons.
Théorème A.3.9 Soit f une fonction propre, convexe et sci de X dans R ∪ {+∞}. Alors f ∗∗ = f .
Lien avec le sous-différentiel

Théorème A.3.10 Soit f : X → R ∪ {+∞} et f ∗ sa conjuguée. Alors
u∗ ∈ ∂f (u) ⇐⇒ f (u) + f ∗ (u∗ ) = :u∗ , u; .
Démonstration - Soit u∗ ∈ ∂f (u) :
∀v ∈ X f (v) ≥ f (u) + :u∗ , v − u; .
Donc
f ∗ (u∗ ) ≥ :u∗ , u; − f (u) ≥ sup{:u∗ , v; − f (v) | v ∈ X } = f ∗ (u∗ ) .
On obtient : f (u) + f ∗ (u∗ ) = :u∗ , u;.
Réciproquement, si f (u) + f ∗ (u∗ ) = :u∗ , u; on a pour tout v ∈ X
:u∗ , u; − f (u) = f ∗ (u∗ ) ≥ :u∗ , v; − f (v) ,
:u∗ , v − u; ≤ f (v) − f (u) ,

c’est-à-dire u∗ ∈ ∂f (u). !
Corollaire A.3.3 Si f : X → R ∪ {+∞} est convexe, propre et sci, alors
u∗ ∈ ∂f (u) ⇐⇒ u ∈ ∂f ∗ (u∗ ) .
Démonstration - il suffit d’applique le théorème précédent à f ∗ et d’utiliser le fait que lorsque f

est convexe, propre et sci alors f = f ∗∗ . !
Les différentes photos et/ou diagrammes sont extraits des supports de cours en ligne sui-
vants
– Cours MASTER IAD, 2005-2006, Isabelle Bloch - ENST / Département Signal & Images,
Florence Tupin - ENST / Département Signal & Images, Antoine Manzanera- ENSTA /
Unité d’Électronique et d’Informatique
– Cours de Traitement d’Image, Max Mignotte, Département d’Informatique et de Re-
cherche Opérationnelle, Montréal,
http://www.iro.umontreal.ca/˜mignotte/ift6150
– Stéphane BRES – Département Informatique – INSA de LYON
– Valérie PERRIER - Laboratoire de Modélisation et Calcul de l’IMAG, INPG, Grenoble
Bibliographie
[AK] Aubert G. - Kornprobst P. , Mathematical Problems in Image Processing, Applied Mathema-

tical Science 147, Springer 2006
[Azé] Azé D. , Eléments d’analyse convexe et variationnelle, Ellipses, Paris, 1997.
[BP] Barbu V.- Precupanu Th., Convexity and Optimization in Banach Spaces, Sijthoff & Noord-
hoff, Bucarest 1978.
[BMT] Bloch I., Manzanera A. et Tupin F. , Cours MASTER IAD - 2005-2006, ENST /
Département Signal & Images, - ENSTA / Unité d’Électronique et d’Informatique
[Br] Brézis H., Analyse Fonctionnelle, Masson, Paris, 1987.
[Cham] Chambolle A., An algorithm for total variation minimization and applications, JMIV, 20,
89-97,2004.
[Coh] Cohen L., On active contours models and balloons, Computer Vision, graphics and image
processing : image understanding, 53 (2), 211-218, 1991.
[DL] Dautray R. - Lions J.L., Analyse mathématique et calcul numérique, 9 volumes, Masson, Paris,
1987.
[ET] Ekeland I. - Temam R., Analyse convexe et problèmes variationnels, Dunod, Paris, 1973.
[EG] Evans L. - Gariepy R., Measure theory and fine properties of functions, CRC Press, Boca Raton,
Floride, 1992.
[GW] Gasquet C. et Witomski P., Analyse de Fourier et applications. Masson, 1995.
[KWT] Kass M., Witkin A. et Terzopoulos D. , Snakes : Active contour models , International
Journal of Computer Vision, 1(4) :133-144, 1987.
[LG] Le Guyader C., Imagerie mathématique : segmentation sous contraintes géométriques, Thèse,
Décembre 2004.
[M] Mignotte M.,Cours de Traitement d’Image, Département d’Informatique et de Recherche
Opérationnelle, Montréal, http://www.iro.umontreal.ca/˜mignotte/ift6150
[PAF] Pallara D., Ambrosio L., Fusco N., Functions of bounded variations and free discontinuity
problems, Oxford Mathematical Monographs, 2000.
[Rudin] Rudin W. Analyse réelle et complexe. Masson, 1978.
[TPS] Tisserand É. - Pautex J.-F. et Schweitzer P, Analyse et traitement des signaux, Dunod, Paris,
2004.
107
108 BIBLIOGRAPHIE
Table des matières
1 Introduction 3
1.1 Qu’est-ce qu’une image numérique ? . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.1.1 Quelques définitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.2 Qu’est-ce que le traitement d’image ? . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.2.1 Quelques aspects du Traitement d’Image . . . . . . . . . . . . . . . . . . . 9
1.2.2 Applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2 Traitement ponctuel des images numériques 15

2.1 Correction ponctuelle d’une image -Recadrage de dynamique . . . . . . . . . . . 15
2.1.1 Transformation de recadrage . . . . . . . . . . . . . . . . . . . . . . . . . . 15
2.1.2 Égalisation de l’histogramme . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.1.3 Binarisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.2 Anayse de la netteté d’une image numérique . . . . . . . . . . . . . . . . . . . . . 21
2.2.1 Exemples de critères de netteté . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.3 Corrélation d’images numériques . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
2.3.1 Application à la reconnaissance d’empreintes digitales . . . . . . . . . . . 23
2.3.2 Application à l’analyse de la texture d’une image . . . . . . . . . . . . . . 25
2.4 Transformation de Hough . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.4.1 Principe de la méthode pour la recherche de ligne droite . . . . . . . . . . 26
2.4.2 Problème pour la recherche de ligne verticale Représentation normale
d’une droite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.4.3 Transformation de Hough pour la détection de droites dans une image
binaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
3 Filtrage 31
3.1 Filtrage linéaire des signaux 1D . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.2 Filtrage 2D : convolution /médian . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.2.1 Filtrage spatial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.2.2 Filtrage fréquentiel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.3 Filtrage différentiel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
3.3.1 Calcul par convolution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
3.3.2 Equation de la chaleur . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
3.3.3 Mise en œuvre numérique . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
109
110 TABLE DES MATIÈRES
4 Quelques modèles de restauration d’image 55

4.1 Régularisation de Tychonov . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
4.2 L’espace des fonctions à variation bornée BV (Ω) . . . . . . . . . . . . . . . . . . . 57
4.2.1 Généralités . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
4.2.2 Approximation et compacité . . . . . . . . . . . . . . . . . . . . . . . . . . 60
4.3 Le modèle continu de Rudin-Osher-Fatemi . . . . . . . . . . . . . . . . . . . . . . 62
4.3.1 Présentation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.3.2 Condition d’optimalité du premier ordre . . . . . . . . . . . . . . . . . . . 63
4.4 Modèle discret . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
4.4.1 Le problème discret . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
4.4.2 Algorithme de projection de Chambolle . . . . . . . . . . . . . . . . . . . . 66
5 Méthode des contours actifs 67

5.1 Rappels sur la géométrie des courbes planes . . . . . . . . . . . . . . . . . . . . . 67
5.1.1 Abscisse curviligne - longueur . . . . . . . . . . . . . . . . . . . . . . . . . 67
5.1.2 Étude géométrique locale d’une courbe paramétrée . . . . . . . . . . . . . 69
5.2 Méthodes des contours actifs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
5.2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
5.2.2 Modélisation du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
5.2.3 Conditions d’optimalité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
5.2.4 Un modèle dynamique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
5.2.5 Un exemple . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
5.3 La méthode des lignes de niveau (« Level set ») . . . . . . . . . . . . . . . . . . . . 76
5.3.1 Propagation de fronts . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
5.3.2 La méthode « Level set » . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78
5.3.3 Application à l’image . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
5.4 Le modèle des ballons ( « Balloons ») . . . . . . . . . . . . . . . . . . . . . . . . . 80
5.5 Le modèle de Munford-Shah . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
A Quelques outils mathématiques pour l’image 83

A.1 Optimisation dans les espaces de Banach . . . . . . . . . . . . . . . . . . . . . . . 83
A.1.1 Semi-continuité et convexité de fonctionnelles sur V . . . . . . . . . . . . 83
A.1.2 Gâteaux-différentiabilité des fonctionnelles convexes . . . . . . . . . . . . 83
A.1.3 Minimisation dans un Banach réflexif . . . . . . . . . . . . . . . . . . . . . 86
A.1.4 Exemple : Projection sur un convexe fermé . . . . . . . . . . . . . . . . . . 88
A.2 Formulation variationnelle des équations aux dérivées partielles . . . . . . . . . 89
A.2.1 Théorème de Lax-Milgram . . . . . . . . . . . . . . . . . . . . . . . . . . . 89
A.2.2 Généralités sur les espaces de Sobolev . . . . . . . . . . . . . . . . . . . . . 91
A.2.3 Problème homogène de Dirichlet : cas elliptique . . . . . . . . . . . . . . . 93
A.2.4 Problème de Neumann homogène . . . . . . . . . . . . . . . . . . . . . . . 94
A.3 Analyse convexe non lisse . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95
A.3.1 Théorème de Hahn -Banach . . . . . . . . . . . . . . . . . . . . . . . . . . 95
A.3.2 Sous-différentiel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96
A.3.3 Calcul « pratique » : transformation de Legendre-Fenchel . . . . . . . . . 100

Coursimage09 PDF

Transféré par

Droits d'auteur :

Formats disponibles

Coursimage09 PDF

Transféré par

Informations du document

Titre original

Copyright

Formats disponibles

Partager ce document

Partager ou intégrer le document

Options de partage

Avez-vous trouvé ce document utile ?

Ce contenu est-il inapproprié ?

Droits d'auteur :

Formats disponibles

Coursimage09 PDF

Transféré par

Droits d'auteur :

Formats disponibles

Quelques méthodes mathématiques pour le traitement

M. Bergounioux Master 2 - 2008-2009

1.1 Qu’est-ce qu’une image numérique ?

1.1.1 Quelques définitions

L’échantillonnage est le procédé de discrétisation spatiale d’une image consistant à asso-

Image originale Image

1.2 Qu’est-ce que le traitement d’image ?

1.2.1 Quelques aspects du Traitement d’Image

– Filtrage / déconvolution (ou filtrage inverse)

Dans ce cours nous évoquerons successivement trois aspects fondamentaux :

Traitement ponctuel des images

2.1 Correction ponctuelle d’une image -Recadrage de dynamique

2.1.1 Transformation de recadrage

On suppose une image de départ présentant un histogramme concentré dans l’intervalle

Original Recadrage : a = 30, b = 200

Variantes pour le rehaussement des contrastes

Dilatation de la dynamique des zones claires

Dilatation de la dynamique des zones sombres

2.1.2 Égalisation de l’histogramme

Histogramme d’origine Histogramme plat idéal

Fonction d’aplatissement discrète

Image égalisée Histogramme égalisé

Extraction d’une fenêtre d’intensité

Fonction « fenêtre d’intensité »

Seuillage avec fenêtre d’intensité entre 30 et 100

2.2 Anayse de la netteté d’une image numérique

La focalisation automatique, ou mise au point, des systèmes de prise de vue traditionnels

2.2.1 Exemples de critères de netteté

Critères de netteté basés sur l’analyse de l’histogramme de l’image

On note hk la distribution statistique (histogramme) des niveaux k d’intensité des pixels.

Critères de netteté basés sur l’analyse spectrale de l’image

2.3 Corrélation d’images numériques

Nous considérons deux fonctions bidimensionnelles discrètes f (i, j) et g(i, j) avec 1 ≤ i, j ≤

où µ et σ sont respectivement la moyenne et l’écart type de chaque fonction.

2.3.1 Application à la reconnaissance d’empreintes digitales

Les fonctions d’intercorrélation bidimensionnelles peuvent être utilisées en reconnaissance

La fonction est estimée dans la gamme −20 ≤ k, l ≤ 20 grâce à la formule suivante

La représentation 3D à gauche dans la figure ci-dessous est celle de la fonction d’autocorrélation

2.3.2 Application à l’analyse de la texture d’une image

Nous donnons ci-dessous la fonction d’auto-corrélation des textures « sable » et « tissu ».

2.4 Transformation de Hough

2.4.1 Principe de la méthode pour la recherche de ligne droite

2.4.2 Problème pour la recherche de ligne verticale Représentation normale d’une

x cos θ + y sin θ = ρ . (2.4.3)

θp = −90 + phθ avec p ∈ {0, Nθ } et − 90 + Nθ hθ = 90 ,

ρq = ρmin + qhρ avec q ∈ {0, Nρ } et ρmin + Nρ hρ = ρmax .

Le programme itératif suivant permet le calcul de la transformation de Hough :

Mi,j sont les valeurs binaires de l’image originale : Mi,j ∈ {0, 1}

3.1 Filtrage linéaire des signaux 1D

Supposons qu’ on puisse définir la transformée de Fourier des signaux d’entrée f et de

où ĝ désigne la transformée de Fourier de g. La fonction H est appelée fonction de transfert

Selon la bande rejetée, on rencontre 4 grandes catégories de filtres.

3.2 Filtrage 2D : convolution /médian

3.2.1 Filtrage spatial

G est la fonction de transfert du filtre. On peut distinguer trois types de filtrage :

Généralement le filtre est de dimension d impaire et est symétrique. Dans ce cas

[x1 , x2 ] = [y1 , y2 ] = [−d/2, d/2] ,

Exemple 3.2.2 (Filtre de moyenne passe -bas )

Ce sont des filtres séparables.