Download Convolution et Déconvolution en Traitement d'Images and more Study notes Physics in PDF only on Docsity!
Convolution / déconvolution.
1 Convolution discrète.
La convolution discrète est la base de nombreux algorithmes de modification d’images, et
particulièrement lorsqu’il s’agit de mimer, dans l’espace discret, des modifications localement linéaires
dans l’espace continu. Cette séance a pour but de vous familiariser avec la notion de convolution.
1.1 Principe
Les trois acteurs d’une convolution sont :
1 - Une image discrète représentant l’information que l’on souhaite modifier (image source),
2 - Une image discrète ne contenant que des zéros (image cible),
3 - Un noyau de convolution (image de dimension réduite).
Image source : Is Image cible : Ic Noyau : π
En supposant que les valeurs du noyau sont indexées en partant du centre du noyau tandis que les
images sont généralement indexées en ligne, colonne comme une matrice.
1.2 Algorithme
L’algorithme de la convolution est très simple. En supposant que le noyau ait une dimension 7x
comme représenté ci-dessus, la valeur du pixel (l,c) de l’image cible est obtenue par :
Quand vous programmez, faites attention aux effets de bord. Il existe plusieurs méthodes pour
calculer les valeurs sur les bords. La plus simple est de considérer que toutes les valeurs en
dehors du cadre de l’image sont nulles.
Vous pouvez aussi supposer un « effet miroir ». En supposant que l va de 0 à L, si (l+u < 0) ou
(l+u > N) on remplace (l+u) par (l-u) – et idem pour les colonnes.
1.3 Normalisation du noyau.
Dans la plupart des applications, il est important que la somme des éléments du noyau soient unitaire.
Quand ce n’est pas le cas, on remplace le noyau π par sa version normalisée π/ß avec :
0
1
2
3
-3 -2 -1 0 1 2 3
Ic(l,c) = ∑ ∑ Is(l+u,c+v).π(u,v)
u=-
3
v=-
3
ß = ∑ ∑ π(u,v)
u=-
3
u=-
3
1.4 Quel noyau?
On a vu en cours que la plupart des noyaux discrets sont obtenus par discrétisation de noyaux
continus. Un noyau continu est défini par une fonction de forme f(x,y) et un paramètre d’étendue μ.
On va s’intéresser ici aux noyaux moyenneurs, c’est-à-dire ceux provoquant un flou sur l’image.
Dans ce cas les fonctions de forme ƒ sont à support sur [-1,1]x[-1,1] – c’est-à-dire que f(x,y) = 0 si
(x,y)∉[ - 1,1]x[-1,1]. Enfin la majorité de ces fonction sont séparables, c’est-à-dire qu’elles
s’écrivent comme le produit de deux fonction : f(x,y) = fx(x). fy(y).
Voici un ensemble de fonction de forme pour (x,y)∈[-1,1]x[-1,1] :
- la fonction rectangle : ƒ(x,y) = 1,
- la fonction triangle : : ƒ(x,y) = (1-|x|).(1-|y|),
- la fonction d’Epanechnikov : ƒ(x,y) = (1-x
2 ).(1-y
2 ),
- la fonction circulaire : cos(
!
"
x).cos(
!
"
y),
vous pouvez en trouver beaucoup d’autres. Pour plus de simplicité, j’ai retiré les normalisation dans
l’espace continu. Toutes les fonctions ci-dessus sont positives et à support borné. Cependant elles ne
le sont pas toutes, comme par exemple la fonction sinus cardinal.
Pour passer d’une fonction ƒ à support [-1,1]x[-1,1] à une fonction ƒμ à support [-μ, μ]x[-μ, μ], il
suffit d’écrire :
Noyau d’Epanechnikov à support μ=4.
Pensez à toujours normaliser votre noyaux échantillonné (la somme des coefficient = 1).
1.5 Création de noyaux.
Créez des noyaux échantillonnés de support quelconque (entre 1 et 5). Vous n’êtes pas obligés de
choisir un nombre entier pour μ (vous pouvez choisir par exemple μ=1.2) par contre votre noyau
aura forcément une dimension entière de (2.m+1)x(2.m+1) avec m = [μ], où [.] est l’arrondi (e.g.
[1.2] = 1 ). Vous pouvez aussi créer des noyaux anisotropes (c’est-à-dire que l’étalement en x n’est
pas le même qu’en y) par exemple :.
Visualisez ces noyaux comme image et comme fonction bidimensionnelles. Pour l’image n’oubliez pas de
normaliser entre 0 et 255 et de passer en entier 8 bits.
Faites une fonction qui calcule directement le noyau (ici par exemple Epanechnikov). Dans un fichier
que vous appelerez Noyau.m saisissez :
function MonNoyau = FaireUnNoyau(delta_x,delta_y)
m = ceil(max(delta_x,delta_y)) ;
MonNoyau = zeros(m,m) ;
for u=-m:m
for v=-m:m
x = u/delta_x ;
y = v/delta_y ;
MonNoyau(u+m+1,v+m+1) = max(0,(1-xx))max(0,(1-y*y)) ;
end
end
MonNoyau = MonNoyau / sum(MonNoyau(:)) ;
ƒ (x,y) =
μ
ƒ ( , ) — —
x y
μ (^) μ
ƒ (x,y) = μ
ƒ ( , ) — —
x y
3.2 5.
function SignalCible = Convolution(SignalSource, Noyau)
[N] = length(SignalSource) ; % dimension signal source
SignalCible = zeros(size(SignalSource)) ; % creation signal cible
[n] = length(Noyau) ; % dimension du noyau
if( mod(n,2) )
% la dimension du noyau doit etre impaires
n = floor(n/2) ;
for k=1:N
for u=-n:n
l = k+u ;
if (l>=1) && (l<=N)
SignalCible(k) = SignalCible(k) + % votre code ici
end
end
end
end
Vous pouvez alors créer une fonction ayant pour argument le noyau permettant de filtrer les lignes
qui aurait pour prototype :
function ImageCible = Convolution1D(ImageSource, Noyau)
Pour filtrer les colonnes, il suffit de transposer la matrice ImageCible obtenue en filtrant les
lignes, rappeler la fonction avec le noyau permettant de filtrer les colonnes et retransposer le
résultat :
ImageCible = Convolution1D(ImageSource, Noyau_ligne) ;
ImageCible = ImageCible’ ;
ImageCible = Convolution1D(ImageSource, Noyau_colonne) ;
ImageCible = ImageCible’ ;
Comparez ce que vous obtenez avec différents noyaux 1D et 2D.
2 Deconvolution
2.1 Qu’est-ce que la déconvolution?
Dans de nombreux cas les images acquises en biologie – en microscopie optique par exemple ou en
imagerie par fluorescence, souffrent d’une mauvaise résolution optique. On peut généralement
modéliser cette dégradation par une convolution dans le domaine continu.
Sur le plan de l’image numérique, cette dégradation peut être vue comme une convolution avec un
noyau échantillonné qui est l’échantillonnage du noyau de convolution dans l’espace continu.
2.2 Algorithmes de déconvolution itératifs.
Ces algorithmes ont été vus en cours. On les rappelle ici.
Si M est l’image mesurée (de taille nxm), I l’image originale qu’on souhaite reconstruire qui a été
convoluée avec le noyau π, on écrit M = IÄπ.
L’inversion de cette équation, et donc l’estimation de I, supposerait qu’on possède un noyau inverse
(c’est-à-dire une noyau π
- 1 tel que πÄπ - 1 =∂ (∂ étant un noyau tel que ∂(0,0)=1, et ∂(u,v)=0 si
u.v≠0). Ce noyau n’existe généralement pas (dans le cas du noyau rectangulaire, le noyau inverse
est le noyau sinus cardinal qui est à support infini par exemple).
Une écriture matricielle de cette estimation induirait l’inversion d’une matrice trop volumineuse pour
être inversée avec sécurité, c’est pourquoi on se tourne vers des algorithmes itératifs minimisant
une distance entre l’image reconstruite convoluée avec le noyau et l’image mesurée. Ces algorithmes
consistent en une suite itérative Ik (k allant de 0 à une valeur K prédéfinie) convergeant vers I.
Comme nous utilisons des noyaux symétriques ( π(u,v) = π(-u,-v) ) ces algorithmes prennent une
forme assez simple.
o Algorithme additif (ou ART ou méthode de Schultz ou méthode de Hotteling).
L’objet de cet algorithme est de faire converger la suite d’image Ik de façon à ce que πÄIk se
rapproche de M au sens de la norme Euclidienne : || πÄIk+1 – M || ≤ || πÄIk – M ||.
Cet algorithme s’écrit :
Ik+1 = Ik + l.( π Ä ( M - πÄIk ) ), où l est une valeur assurant la convergence.
On prend généralement l = || π ||
- 2 , où ||. || est la norme euclidienne (racine de la somme des carrés
des valeurs).
o Algorithme multiplicatif (ou MLEM).
L’objet de cet algorithme est de faire converger la suite d’image Ik de façon à ce que πÄIk se
rapproche de M au sens de la différence des logarithmes. On essaye de rapprocher πÄIk de M de
façon à ce que la division terme à terme de πÄIk par M tende vers 1 (une image ne contenant que
des 1).
Cet algorithme s’écrit :
Ik+1 = Ik. ( π Ä ( M / ( πÄIk ) ) ), où / est la division terme à terme.
2.3 Programmation des algorithmes.
Vous disposez d’une image acquise par un microscope électronique « ImageMicroscope.tif » et d’une
coupe cytoplasmique d’ascidie « CoupeAscidie.tif ». Après avoir placé une mire sous le microscope, on
a pu identifier que la dégradation subie par l’image issue du microscope électronique était équivalente
à une convolution avec un noyau d’Epanechnikov avec μ = 5 , et, pour la coupe cytoplasmique, la
dégradation est équivalente à une convolution avec un noyau circulaire d’étendu μ = 8.
Essayez les deux algorithmes et comparer 1/ leur vitesse de convergence 2/ le résultat visuel.
Image obtenue par un microscope électronique.
Coupe cytoplasmique d’une ascidie.
Pour simplifier le travail, utilisez les fonctions de convolution que vous avez déjà créés.
Vous pouvez visualiser la convergence de votre algorithme en affichant l’erreur résiduelle qui est
la somme des valeurs absolues des écarts entre πÄIk et M.
