La transformée de Fourier en traitement d'images

Ce cours introduira l'utilisation des transformées de Fourier en traitement numérique d'images. Nous verrons en quoi ces opérations permettent un gain de temps intéressant pour l'application de certains types de filtres.

Ce cours est divisé en deux parties distinctes. Une partie sera assez théorique, nous verrons d'un point de vue purement théorique l'utilisation des transformées de Fourier sans détailler une implémentation possible.
Ensuite, un autre cours présentera la bibliothèque FFTW qui permet de calculer des transformées de Fourier de signaux (à une ou plusieurs dimensions) et nous l'utiliserons pour nos images.

  • Cours théorique sur l'utilisation des transformées de Fourier : ici.
  • Les sources C++ (ici)
  • Le cours sur l'utilisation de la bibliothèque FFTW : ici.

Commentez Donner une note à l'article (5)

Article lu   fois.

L'auteur

Site personnel

Liens sociaux

Viadeo Twitter Facebook Share on Google+   

I. Chapitre 1 - Introduction

Le traitement numérique d'images est un domaine vaste, qui admet de plus en plus d'applications, que ce soit en imagerie médicale ou pour les particuliers. Nous avons pu constater dans le cours « Introduction au traitement numérique d'images » qu'une grosse partie permettant l'analyse des images utilise la notion de filtre. Malheureusement, les méthodes simples qui permettent l'application de filtre à tout type d'images sont parfois très gourmandes en temps d'exécution.

La transformée de Fourier est un outil mathématique, qui permet de remplacer des opérations gourmandes par d'autres, plus rapides. Cette technique ne pouvait être fortement employée avant, car à l'époque, il n'existait d'algorithme permettant de calculer rapidement la transformée d'une image. Heureusement pour nous, en 1960, on a découvert un algorithme rapide permettant le calcul de la transformée de Fourier, ce qui va créer une petite révolution dans le domaine. Dans la littérature, on trouve souvent la dénomination fft, comme fast fourier transform.

On peut également constater que la transformée ne sert pas uniquement à l'application de filtres, mais permet par exemple, de traiter directement les fréquences des images (à la manière des signaux) ou encore de compresser des images.

Ce cours ne traitera pas du vaste domaine d'application de la transformée de Fourier, mais de l'utilisation de celle-ci pour l'application de filtre. Ceci se fera surtout d'un point de vue théorique et nous n'entrerons pas dans l'implémentation effective de cette transformée. Néanmoins, il existe pour la plupart des langages des bibliothèques de fonctions permettant de réaliser ces transformations (notamment la bibliothèque FFTW disponible sur plusieurs plates-formes pour le langage C).

II. Chapitre 2 - Transformée de Fourier

De nombreux traitements utilisent la notion de filtre. Ces filtres peuvent s'écrire sous la forme d'un produit de convolution entre une matrice (ou noyau de convolution) et d'une image.

Malheureusement, en appliquant de gros masques de convolution, on peut constater que le temps de calcul devient très important. En effet, pour une image de taille kitxmlcodeinlinelatexdvpn^2finkitxmlcodeinlinelatexdvp et un masque de taille kitxmlcodeinlinelatexdvpm^2finkitxmlcodeinlinelatexdvp. Le nombre de calculs à réaliser pour appliquer le masque est de kitxmlcodeinlinelatexdvpn^2 \cdot m^2finkitxmlcodeinlinelatexdvp. La taille du masque peut être au maximum de kitxmlcodeinlinelatexdvpn^2finkitxmlcodeinlinelatexdvp (s'il est plus gros, certains coefficients seront inutiles pour le produit de convolution). Le temps de calcul devient alors très important. Pour cela, on peut appliquer les masques dans l'espace des fonctions de Fourier qui rendront les calculs beaucoup plus rapides, car il existe un algorithme de calcul de la transformée très rapide. Dans le cas où le masque est très petit, il vaut mieux éviter de passer dans l'espace de Fourier.

II-A. Définitions

II-A-1. Espace continu

Dans l'espace continu kitxmlcodeinlinelatexdvp\mathbb{C}^2finkitxmlcodeinlinelatexdvp, on définit la transformée de Fourier d'une fonction kitxmlcodeinlinelatexdvpf \in \mathcal{L}^1(\mathbb{R}^2)finkitxmlcodeinlinelatexdvp par :

kitxmlcodelatexdvp\forall(u, v) \in \mathbb{C}^2, \mathcal{F}(f)(u, v) = \int \int_{\mathbb{R}^2} f(x,y) \cdot \exp(-2i\pi(ux+uv))dx\,dyfinkitxmlcodelatexdvp

La transformée de Fourier inverse étant définie pour toute application kitxmlcodeinlinelatexdvpFfinkitxmlcodeinlinelatexdvp de l'espace de fonction kitxmlcodeinlinelatexdvp\mathcal{L}^1(\mathbb{R}^2)finkitxmlcodeinlinelatexdvp par :

kitxmlcodelatexdvp\forall(u, v) \in \mathbb{C}^2, \mathcal{F}^{-1}(F)(u, v) = \frac{1}{4\pi^2}\int \int_{\mathbb{R}^2} F(x,y) \cdot \exp(+2i\pi(ux+uv))dx\,dyfinkitxmlcodelatexdvp

à noter que kitxmlcodeinlinelatexdvp\mathcal{F}^{-1}(\mathcal{F}(f)) = ffinkitxmlcodeinlinelatexdvp sur kitxmlcodeinlinelatexdvp\mathbb{R}^2finkitxmlcodeinlinelatexdvp.

En traitement d'images, la transformée Fourier ne nous servira que dans l'espace des fonctions kitxmlcodeinlinelatexdvp\mathcal{L}^1(\mathbb{C}^2)finkitxmlcodeinlinelatexdvp.

II-A-2. Espace discret

Dans l'espace discret kitxmlcodeinlinelatexdvp\mathbb{Z}^2finkitxmlcodeinlinelatexdvp, on peut définir la transformée de Fourier pour deux types de suites. Soit une suite à support bornée, soit une suite périodique de période kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp.

Pour une suite à support bornée, nous réduirons les possibilités en choisissant une suite appartenant à l'espace des fonctions de kitxmlcodeinlinelatexdvp[[-n/2,n/2]]^2finkitxmlcodeinlinelatexdvp dans kitxmlcodeinlinelatexdvp\mathbb{R}finkitxmlcodeinlinelatexdvp (bien que kitxmlcodeinlinelatexdvp\mathbb{C}finkitxmlcodeinlinelatexdvp soit possible). Pour une suite périodique, la définition n'est pas différente, car nous allons accéder seulement à la fonction sur kitxmlcodeinlinelatexdvp[[-n/2,n/2]]^2finkitxmlcodeinlinelatexdvp. On définit la transformée de Fourier discrète par :

kitxmlcodelatexdvp\forall(a, b) \in [[-n/2,n/2]]^2, \mathcal{F}(f)(a, b) = \sum_{i=-n/2}^{n/2}\sum_{j=-n/2}^{n/2} u(i,j)\cdot\exp\left(-2i\pi\left (\frac{ai}{n} + \frac{bi}{n} \right )\right )finkitxmlcodelatexdvp

Et la transformée inverse étant définie par :

kitxmlcodelatexdvp\forall(a, b) \in [[-n/2,n/2]]^2, \mathcal{F}^{-1}(f)(a, b) = \frac{1}{n^2}\sum_{i=-n/2}^{n/2}\sum_{j=-n/2}^{n/2} u(i,j)\cdot\exp\left(+2i\pi\left (\frac{ai}{n} + \frac{bi}{n} \right )\right )finkitxmlcodelatexdvp

à noter que kitxmlcodeinlinelatexdvp\mathcal{F}^{-1}(\mathcal{F}(f)) = ffinkitxmlcodeinlinelatexdvp sur kitxmlcodeinlinelatexdvp[[-n/2,n/2]]^2finkitxmlcodeinlinelatexdvp.

Le fait que la fonction soit centrée autour du point (0, 0) n'est pas obligatoire, mais cela permet une application des filtres dans l'espace de Fourier de manière plus simple.

II-A-3. Transformée de Fourier d'une image

Dans la définition précédente, les fonctions sont à valeurs dans un espace de dimension 1. Mais les images sont en général à valeurs dans un espace vectoriel de dimension supérieur. Par exemple 3 pour des images en couleurs RGB. Ainsi, dans le cas où la dimension est supérieure à 1 (par exemple kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp). Pour étendre la formule de la transformée de Fourier, on définit une base kitxmlcodeinlinelatexdvp(e_k)_{k\in[[1,n]]}finkitxmlcodeinlinelatexdvp de l'espace vectoriel. Puis on détermine la projection de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp dans chaque dimension. On définit ainsi une sous-image kitxmlcodeinlinelatexdvpf_kfinkitxmlcodeinlinelatexdvp à valeur dans kitxmlcodeinlinelatexdvp\mathbb{R}finkitxmlcodeinlinelatexdvp définie par :

kitxmlcodelatexdvp\forall k \in[[1,n]], f_k =\,<f,e_k>finkitxmlcodelatexdvp

On reconstruit ensuite la transformée de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp par la formule :

kitxmlcodelatexdvp\mathcal{F}(f) = \sum_{k=1}^n \mathcal{F}(f_k)\cdot e_kfinkitxmlcodelatexdvp

Avec des mots, cela revient à déterminer la transformée de Fourier pour chaque canal de couleurs, à calculer leurs transformées de Fourier, puis à regrouper ensuite les données.

II-B. Propriétés

Nous allons uniquement énoncer et démontrer certaines propriétés par rapport au produit de convolution.

II-B-1. Espace continu

Théorème 1

Soit kitxmlcodeinlinelatexdvpf\in\mathcal{L}^1(\mathbb{R}^2)finkitxmlcodeinlinelatexdvp et kitxmlcodeinlinelatexdvpg\in\mathcal{L}^1(\mathbb{R}^2)finkitxmlcodeinlinelatexdvp bornée, alors :

kitxmlcodelatexdvp\mathcal{F}(f*g) = \mathcal{F}(f)\cdot\mathcal{F}(g)finkitxmlcodelatexdvp

La multiplication est la multiplication usuelle dans l'espace des fonctions kitxmlcodeinlinelatexdvp\mathcal{L}^1(\mathbb{C}^2)finkitxmlcodeinlinelatexdvp. La condition que kitxmlcodeinlinelatexdvpgfinkitxmlcodeinlinelatexdvp soit bornée n'est pas forcement nécessaire, en pratique cela suffit et cela facilite la preuve. Celle-ci est disponible en annexePreuve dans l'espace continu.

II-B-2. Espace discret

Le théorème sur un espace discret admet quelques restrictions sur la nature des suites.

Théorème 2

  • Soit kitxmlcodeinlinelatexdvphfinkitxmlcodeinlinelatexdvp une suite de kitxmlcodeinlinelatexdvp[[−n/2,n/2]]^2finkitxmlcodeinlinelatexdvp à valeurs dans kitxmlcodeinlinelatexdvp\mathbb{R}finkitxmlcodeinlinelatexdvp.
  • Soit kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp une suite de kitxmlcodeinlinelatexdvp\mathbb{Z}^2finkitxmlcodeinlinelatexdvp de période kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp.

Alors :

kitxmlcodelatexdvp\forall(a,b) \in [[-n/2,n/2]]^2, \mathcal{F}(h*f)(a,b) = \mathcal{F}(h)(a,b)\cdot\mathcal{F}(f)(a,b)finkitxmlcodelatexdvp

Remarquons quelque chose d'important, dans l'expression kitxmlcodeinlinelatexdvp\mathcal{F}(h)\cdot\mathcal{F}(f)finkitxmlcodeinlinelatexdvp, la périodicité de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp n'est pas importante. Elle l'est dans le produit de convolution, car le produit avec une fonction bornée n'est pas bien défini comme on a pu le voir. Il y a plusieurs manières de remplir les zones manquantes. Nous verrons dans la partie application comment passer outre le problème.

La preuve est accessible en annexePreuve dans l'espace discret.

II-B-3. Applications successives de filtres

Supposons que l'on dispose de deux filtres kitxmlcodeinlinelatexdvph_1finkitxmlcodeinlinelatexdvp et kitxmlcodeinlinelatexdvph_2finkitxmlcodeinlinelatexdvp de taille kitxmlcodeinlinelatexdvpmfinkitxmlcodeinlinelatexdvp, d'une suite kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp de période n et que l'on souhaite calculer de manière rapide :

kitxmlcodelatexdvph_2*(h_1*f)finkitxmlcodelatexdvp

A priori, kitxmlcodeinlinelatexdvph_1finkitxmlcodeinlinelatexdvp et kitxmlcodeinlinelatexdvph_2finkitxmlcodeinlinelatexdvp sont finis. On dispose de deux possibilités.

On peut calculer directement kitxmlcodeinlinelatexdvph_2*h_1finkitxmlcodeinlinelatexdvp. On tronque ou on ajoute des 0 au produit pour que sa taille soit kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp. On se ramène ainsi directement au cas précédent.

On tronque ou on ajoute des 0 à kitxmlcodeinlinelatexdvph_1finkitxmlcodeinlinelatexdvp et à kitxmlcodeinlinelatexdvph_2finkitxmlcodeinlinelatexdvp pour que leur taille soit kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp. Comme kitxmlcodeinlinelatexdvp(h_1*f)finkitxmlcodeinlinelatexdvp est périodique de période kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp, on peut écrire ceci :

kitxmlcodelatexdvp\mathcal{F}(h_2*(h_1*f)) = \mathcal{F}(h_2) \cdot\mathcal{F}(h_1*f)) = \mathcal{F}(h_2) \cdot \mathcal{F}(h_1) \cdot \mathcal{F}(f)finkitxmlcodelatexdvp

Il suffit donc de multiplier les transformées des filtres (que l'on a complétées ou tronquées).

II-C. Module et phase

Les images que l'on obtient en appliquant une transformée de Fourier donnent une image complexe. On pourrait alors se demander comment les représenter. En général, on calcule le module kitxmlcodeinlinelatexdvpF_mfinkitxmlcodeinlinelatexdvp et la phase kitxmlcodeinlinelatexdvpF_pfinkitxmlcodeinlinelatexdvp de l'image, et l'on représente son module.

Ces deux images réelles peuvent être définies ainsi :

kitxmlcodelatexdvpRe(\mathcal{F}(f)(x,y)) = F_m(x,y)\cos(F_p(x,y)) \\ Im(\mathcal{F}(f)(x,y)) = F_m(x,y)\sin(F_p(x,y))finkitxmlcodelatexdvp

Où kitxmlcodeinlinelatexdvpRefinkitxmlcodeinlinelatexdvp et kitxmlcodeinlinelatexdvpImfinkitxmlcodeinlinelatexdvp désignent les parties réelles et imaginaires. On constate alors que kitxmlcodeinlinelatexdvpF_pfinkitxmlcodeinlinelatexdvp n'est pas unique.

En général, pour représenter la transformée, on représente uniquement le module. Voir la figure 2.1.

À noter que l'image est carrée, on a complété l'image du papillon par du noir pour la rendre carrée. En effet, la bibliothèque que j'ai utilisée pour calculer la transformée de Fourier semble avoir certains problèmes lorsque les images sont rectangulaires. De plus, on travaille plus fréquemment avec des images carrées (que ce soit du point de vue théorique que pratique).

Image non disponible
Fig. 2.1 - Module du canal bleu du papillon

Mais en fait, la phase de l'image complexe dispose de beaucoup d'informations. Voici deux images dont ont a permuté les modules (figure 2.2).

Image non disponible
Fig. 2.2 - Inversion des modules de deux images

II-D. Phénomène de repliement de spectre

Lorsque l'on réalise l'acquisition d'une image (une image réelle), on effectue un échantillonnage (ce qui correspond à une discrétisation de l'image réelle). En effet, l'image réelle étant définie de manière continue, on ne garde que les valeurs espacées d'une certaine distance kitxmlcodeinlinelatexdvpD_efinkitxmlcodeinlinelatexdvp (on considère que le pas d'échantillonnage est identique à la verticale et à l'horizontale pour simplifier les explications).

Par exemple, lorsque l'on échantillonne l'image kitxmlcodeinlinelatexdvpf :(\mathbb{R}, \mathbb{R}) -> f(x, y)finkitxmlcodeinlinelatexdvp, on ne considère que les points kitxmlcodeinlinelatexdvpf(n \cdot D_e, m \cdot D_e)finkitxmlcodeinlinelatexdvp où kitxmlcodeinlinelatexdvp(n, m) \in\mathbb{Z}^2finkitxmlcodeinlinelatexdvp. En unités normalisées, on écrit plus simplement kitxmlcodeinlinelatexdvpf(n,m)finkitxmlcodeinlinelatexdvp. La fréquence d'échantillonnage est définie par kitxmlcodeinlinelatexdvpf_e = \frac{1}{D_e}finkitxmlcodeinlinelatexdvp.

À noter que ceci est un modèle théorique, pour des acquisitions réelles, il n'est pas forcément possible de prendre une valeur ponctuelle, on réalise plutôt une moyenne sur un petit intervalle.

Intuitivement, on remarque que si le pas d'échantillonnage est trop important (c'est-à-dire si la fréquence d'échantillonnage est trop petite), alors l'image sera de mauvaise qualité (elle ne correspondra pas vraiment à l'image réelle). Il est possible théoriquement de quantifier la valeur que doit prendre la fréquence, le théorème de Nyquist-Shannon nous apprend que l'on doit choisir une fréquence d'échantillonnage au moins deux fois supérieure à la fréquence maximale de l'image. Dans le cas contraire, on obtient le phénomène de repliement de spectre (ou aliasing).

Considérons la figure 2.3. Les hautes fréquences dans cette image sont très présentes, à tel point que l'on voit apparaître des motifs particuliers au centre (la dénomination est : effet de Moiré). On peut donc supposer que la fréquence d'échantillonnage est légèrement trop faible. À cause de la périodisation de l'image et du fait que la fréquence maximale est supérieure à la fréquence de coupure, on voit apparaître sur les bords de la figure 2.4 (qui correspond au module de la transformée) un repliement de spectre (les hautes fréquences ont tendance à revenir vers le centre à cause de la périodisation). Les hautes fréquences réelles viennent à se superposer avec des fréquences plus basses, ce qui rend impossible de différencier une basse d'une haute fréquence.

Maintenant, considérons cette même image (en figure 2.5), mais avec une fréquence d'échantillonnage réduite par 4. On peut constater que dans l'image réelle, les motifs circulaires ont totalement disparu. On peut également voir que le module de la transformée (figure 2.6) laisse apparaître un repliement encore plus important.

Effet dans l'espace réel : dans l'espace réel, le repliement de spectre se caractérise principalement par :

  • l'apparition de nouvelles formes ;
  • la disparition de motifs particuliers (dans notre exemple, on voit que les motifs circulaires ont disparu dans l'image sous-échantillonnée).
Image non disponible
Fig. 2.3 - Image légèrement sous-échantillonnée
Image non disponible
Fig. 2.4 - Module de la transformée de l'image, on voit apparaître un repliement du spectre sur les bords
Image non disponible
Fig. 2.5 - Image fortement sous-échantillonnée, les motifs circulaires ont disparu
Image non disponible
Fig. 2.6 - Module de la transformée de l'image fortement sous-échantillonnée

III. Chapitre 3 - Application

III-A. Produit de convolution

Nous allons voir des exemples d'applications concernant le produit de convolution. Comme nous l'avons vu précédemment, lorsque l'on passe dans le domaine de Fourier, le produit de convolution devient un produit simple. Mais des problèmes surviennent dus au fait que l'image est considérée comme périodique. Nous allons voir plusieurs stratégies pour régler nos soucis.

Pour l'instant, on suppose que l'on dispose d'une image kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp définie sur kitxmlcodeinlinelatexdvp[[-1,1]]^2finkitxmlcodeinlinelatexdvp. On souhaite lui appliquer une matrice de convolution kitxmlcodeinlinelatexdvphfinkitxmlcodeinlinelatexdvp. Sa taille peut être plus petite, plus grande ou égale à celle de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp.

Nous rappelons que le produit de convolution entre kitxmlcodeinlinelatexdvphfinkitxmlcodeinlinelatexdvp et kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp n'est pas bien défini sur le bord. En effet, si :

kitxmlcodelatexdvph = \begin{pmatrix} 1 & 2 & 1\\ 2 & 1 & 2\\ 1 & 2 & 1 \end{pmatrix}finkitxmlcodelatexdvp

Et

kitxmlcodelatexdvpf = \begin{pmatrix} 2 & -1 & 0\\ 4 & -2 & -1\\ 1 & 2 & 3 \end{pmatrix}finkitxmlcodelatexdvp

On remarque qu'il manque des informations pour calculer la convolution en kitxmlcodeinlinelatexdvpf(-1,-1)finkitxmlcodeinlinelatexdvp (on a considéré que la fonction kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp était centrée en kitxmlcodeinlinelatexdvp(0, 0)finkitxmlcodeinlinelatexdvp). La convolution donnant en kitxmlcodeinlinelatexdvp(-1,-1)finkitxmlcodeinlinelatexdvp : kitxmlcodeinlinelatexdvp1*2 + 2*(-1) + 2*4 + 1*(-2)+ \dotsfinkitxmlcodeinlinelatexdvp (il manque des valeurs autour de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp pour la calculer correctement).

Il y a alors plusieurs stratégies qui s'offrent à nous, soit on agrandit l'image kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp, soit on ne l'agrandit pas.

III-A-1. Agrandissement

Zero-padding : si on l'agrandit, on peut compléter kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp par des valeurs nulles (zero-padding). Dans notre exemple, on aura kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp qui vaudra :

kitxmlcodelatexdvp f = \begin{pmatrix} 0 & 0 & 0 & 0 & 0\\ 0 & 2 & -1 & 0 & 0\\ 0 & 4 & -2 & -1 & 0\\ 0 & 1 & 2 & 3 & 0\\ 0 & 0 & 0 & 0 & 0 \end{pmatrix}finkitxmlcodelatexdvp

Dans ce cas, la convolution en kitxmlcodeinlinelatexdvp(-1,-1)finkitxmlcodeinlinelatexdvp : kitxmlcodeinlinelatexdvp1*2 + 2*(-1) + 2*4 + 1*(-2) + 0*(1 + 2 + 1 + 2 + 1)finkitxmlcodeinlinelatexdvp. Ce qui peut fausser la convolution (par exemple, rendre plus sombre l'image si on applique un flou).

Extrapolation : on peut également extrapoler les valeurs de l'image :

kitxmlcodelatexdvpf = \begin{pmatrix} 2 & 2 & -1 & 0 & 0\\ 2 & 2 & -1 & 0 & 0\\ 4 & 4 & -2 & -1 &-1\\ 1 & 1 & 2 & 3 & 3\\ 1 & 1 & 2 & 3 & 3 \end{pmatrix}finkitxmlcodelatexdvp

Dans ce cas, la convolution est : kitxmlcodeinlinelatexdvp1*2 + 2*(-1) + 2*4 + 1*(-2) + 1*(-1) + 2*2 + 1*2 + 2*2 + 1*4finkitxmlcodeinlinelatexdvp.

Afin d'extrapoler les valeurs, j'ai utilisé une méthode assez simple, il est évidemment possible d'utiliser d'autres techniques que je ne détaillerai pas.

Une fois l'image agrandie de cette manière, on peut donc calculer sa transformée de Fourier.

Mais pour utiliser le théorème précédent, a priori, il faudrait compléter le filtre kitxmlcodeinlinelatexdvphfinkitxmlcodeinlinelatexdvp par des 0 (pour que sa taille soit la même que celle de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp). Dans ces deux cas, kitxmlcodeinlinelatexdvphfinkitxmlcodeinlinelatexdvp deviendra :

kitxmlcodelatexdvph = \begin{pmatrix} 0 & 0 & 0 & 0 & 0\\ 0 & 1 & 2 & 1 & 0\\ 0 & 1 & 1 & 1 & 0\\ 0 & 1 & 2 & 1 & 0\\ 0 & 0 & 0 & 0 & 0 \end{pmatrix}finkitxmlcodelatexdvp

III-A-2. Non-agrandissement

Si on décide de ne pas agrandir l'image, elle sera considérée comme périodique (à cause de la transformation de Fourier qui périodise l'image). Ainsi, l'image deviendra automatiquement (on a borné l'image pour l'exemple) :

kitxmlcodelatexdvpf = \begin{pmatrix} 3 & 1 & 2 & 3 & 1\\ 0 & 2 &-1 & 0 & 2\\ -1 & 4 &-2 &-1 & 4\\ 3 & 1 & 2 & 3 & 1\\ 0 & 2 &-1 & 0 & 2 \end{pmatrix}finkitxmlcodelatexdvp

Mais il est important de bien voir que l'on n'aura pas à modifier kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp, l'image sera automatiquement considérée comme périodique (kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp restera donc dans kitxmlcodeinlinelatexdvp[[-1,1]]^2finkitxmlcodeinlinelatexdvp).

Ici, il n'y a pas à augmenter la taille du filtre kitxmlcodeinlinelatexdvphfinkitxmlcodeinlinelatexdvp étant donné que le filtre est de même taille que l'image.

On peut ainsi appliquer le théorème, calculer les transformées de Fourier et calculer le produit. Une fois réalisé, il faut appliquer la transformée de Fourier inverse et tronquer l'image kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp dans kitxmlcodeinlinelatexdvp[[-1,1]]^2finkitxmlcodeinlinelatexdvp.

III-A-3. Quand utiliser la transformée de Fourier ?

Nous avons vu que l'on devait compléter le filtre par des 0 dans le cas où l'image était plus grosse que le filtre. Pour un filtre assez petit, il vaut souvent mieux appliquer directement le produit de convolution plutôt que de passer dans l'espace de Fourier.

Les algorithmes rapides permettant de calculer la transformée d'une image sont en kitxmlcodeinlinelatexdvpt\log tfinkitxmlcodeinlinelatexdvp où kitxmlcodeinlinelatexdvptfinkitxmlcodeinlinelatexdvp désigne la taille de l'image. Avec kitxmlcodeinlinelatexdvpt=n^2finkitxmlcodeinlinelatexdvp, on trouve que la complexité est en kitxmlcodeinlinelatexdvpn^2\log nfinkitxmlcodeinlinelatexdvp. Or, si le filtre est de taille kitxmlcodeinlinelatexdvpm^2finkitxmlcodeinlinelatexdvp, le calcul du produit de convolution s'effectuera en kitxmlcodeinlinelatexdvpn^2\cdot m^2finkitxmlcodeinlinelatexdvp. On en déduit donc que si kitxmlcodeinlinelatexdvpm^2finkitxmlcodeinlinelatexdvp est inférieur (en ordre de grandeur) à kitxmlcodeinlinelatexdvp\log nfinkitxmlcodeinlinelatexdvp, le passage dans le domaine de Fourier ne sera plus bénéfique.

III-B. Application pour le flou gaussien

La transformée de Fourier peut être très utile pour appliquer des filtres gaussiens dont la matrice de convolution est assez importante. En effet, peu importe la taille de la matrice de convolution, le temps sera toujours identique. On rappelle que la fonction gaussienne est définie comme suit :

kitxmlcodelatexdvpg_\sigma(x,y) = \frac{1}{2\pi\sigma^2}\exp \left ( -\frac{x^2+y^2}{2\sigma^2}\right )finkitxmlcodelatexdvp

On pourrait calculer le masque, le passage dans le domaine de Fourier et le produit. Mais cette fonction admet une propriété importante, on peut facilement calculer sa transformée de Fourier continue. Comme kitxmlcodeinlinelatexdvpg_\sigmafinkitxmlcodeinlinelatexdvp est intégrable sur kitxmlcodeinlinelatexdvp\mathbb{R}^2finkitxmlcodeinlinelatexdvp :

kitxmlcodelatexdvp\begin{align} \mathcal{F}(g_\sigma)(u,v) & = \frac{1}{2\pi\sigma^2}\iint_{\mathbb{R}^2}\exp \left ( -\frac{i\pi}{\sigma^2}(x^2+y^2)(ux+uv)\right )dx\,dy\\ & = \exp \left( \frac{-\sigma^2(u^2+v^2)}{2}\right ) \end{align}finkitxmlcodelatexdvp

On peut donc appliquer une matrice de convolution sans se préoccuper de la taille du noyau, car la fonction est bien définie partout dans l'espace de Fourier. Cependant, il faut éventuellement compléter kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp comme on l'a vu précédemment. On obtient directement le produit de convolution avec la formule suivante :

kitxmlcodelatexdvp(g_\sigma*f)(x,y) = \mathcal{F}^{-1}(\mathcal{F}(f)\cdot\mathcal{F}(g_\sigma))finkitxmlcodelatexdvp

Exemple d'application : voir la figure 3.1. Ici, nous avons tout d'abord rendu l'image carrée en la complétant par du noir (en bas) ce qui explique le noir qui apparaît en bas. Celui qui apparaît en haut s'explique par le fait que l'image est considérée comme périodique (celui-ci provient donc du bas).

Image non disponible
Fig. 3.1 - Flou gaussien obtenu en passant dans le domaine de Fourier

Une image où l'on constate mieux ce phénomène est l'image 3.2. Le rouge s'est diffusé à droite à cause de la périodicité et le rouge ne s'est pas diffusé en bas, car on a complété pour avoir l'image carrée.

Image non disponible
Fig. 3.2 - Flou gaussien obtenu en passant dans le domaine de Fourier

III-C. Filtre sur les fréquences

Une première application directe de la transformée de Fourier est de pouvoir permettre de retirer des fréquences de l'image de manière très simple et rapide.

III-C-1. Passe-bas, passe-haut

Passe-bas

L définition de base d'un filtre passe-bas est de ne laisser que les basses fréquences d'une image (par exemple en dessous d'une fréquence de coupure kitxmlcodeinlinelatexdvpf_cfinkitxmlcodeinlinelatexdvp).

Soit kitxmlcodeinlinelatexdvp f\in\mathbb{I}finkitxmlcodeinlinelatexdvp une image et kitxmlcodeinlinelatexdvp\mathcal{F}(f)finkitxmlcodeinlinelatexdvp sa transformée.

On définit la fonction kitxmlcodeinlinelatexdvpH\in\mathcal{F}(\mathbb{I},\mathbb{I})finkitxmlcodeinlinelatexdvp par :

kitxmlcodelatexdvp\begin{align} H(f)(u,v) & = 0\ \mathrm{si}\ \left \| u,v \right \| > f_c \\ & = \mathcal{F}(f)(u,v)\ \mathrm{sinon} \end{align}finkitxmlcodelatexdvp

La fonction passe-bas kitxmlcodeinlinelatexdvpG\in\mathcal{F}(\mathbb{I},\mathbb{I})finkitxmlcodeinlinelatexdvp peut être définie par :

kitxmlcodelatexdvpG(f) = \mathcal{F}^{-1}(H(f))finkitxmlcodelatexdvp

Son effet sera similaire à un lissage d'image.

Passe-haut

De la même manière, un filtre passe-haut ne doit laisser passer que les hautes fréquences. On peut le définir comme suit :

On définit la fonction kitxmlcodeinlinelatexdvpH\in\mathcal{F}(\mathbb{I},\mathbb{I})finkitxmlcodeinlinelatexdvp par :

kitxmlcodelatexdvp\begin{align} H(f)(u,v) & = 0\ \mathrm{si}\ \left \| u,v \right \| < f_c \\ & = \mathcal{F}(f)(u,v)\ \mathrm{sinon} \end{align}finkitxmlcodelatexdvp

Puis la fonction passe-haut est définie par :

kitxmlcodelatexdvpG(f) = \mathcal{F}^{-1}(H(f))finkitxmlcodelatexdvp

Son effet sera similaire à une détection de bords.

III-C-2. Filtre de Butterworth

Filtre passe-haut

Le filtre passe-haut de type Butterworth peut être défini comme suit.

On choisit un entier kitxmlcodeinlinelatexdvpnfinkitxmlcodeinlinelatexdvp, on définit la fonction kitxmlcodeinlinelatexdvpHfinkitxmlcodeinlinelatexdvp à valeur dans kitxmlcodeinlinelatexdvp\mathbb{C}finkitxmlcodeinlinelatexdvp :

kitxmlcodelatexdvp\begin{align} H_n(u,v) & = 0\ \mathrm{pour}\ u=v=0\\ & = \frac{1}{1+\left ( \frac{f_c}{\left \| u,v \right \|}\right )^{2n}}\ \mathrm{sinon} \end{align}finkitxmlcodelatexdvp

La fonction filtre passe-haut de Butterworth sera définie par :

kitxmlcodelatexdvpG(f) = \mathcal{F}^{-1}(\mathcal{F}(f)\cdot H_n)finkitxmlcodelatexdvp

Exemple d'application : vous pouvez voir la figure 3.3

On peut remarquer que le filtre passe-haut ressemble à une détection de bords.

Image non disponible
Fig. 3.3 - Butterworth passe-haut

Filtre passe-bas

De la même manière que précédemment, on peut créer un filtre passe-bas de type Butterworth en définissant une fonction kitxmlcodeinlinelatexdvpH_nfinkitxmlcodeinlinelatexdvp par :

kitxmlcodelatexdvpH_n(u,v) = \frac{1}{1+\left( \frac{\left \| u,v \right \|}{f_c}\right )^{2n}}finkitxmlcodelatexdvp

Ensuite, la fonction filtre passe-bas de type Butterworth sera définie par :

kitxmlcodelatexdvpG(f) = \mathcal{F}^{-1}(\mathcal{F}(f)\cdot H_n)finkitxmlcodelatexdvp

Exemple d'application : vous pouvez voir la figure 3.4.

Image non disponible
Fig. 3.4 - Butterworth passe-bas

On peut remarquer que cette fonction ressemble à une fonction de lissage.

IV. Chapitre 4 - Annexe

IV-A. Preuves

IV-A-1. Preuve dans l'espace continu

On utilise le fait que g soit bornée et que f soit intégrable pour montrer que le produit de convolution existe. En effet, on a :

kitxmlcodelatexdvp\forall(u,v)\in\mathbb{R}^2, |f(x,y)g(x-u,y-v)|\le\left \| g \right \|_\infty|f(t)|finkitxmlcodelatexdvp

Le produit de convolution existe donc. En utilisant le théorème de continuité sous le signe, on montre que le produit de convolution est continu.

Ensuite, on peut montrer que le produit de convolution est intégrable sur kitxmlcodeinlinelatexdvp\mathbb{R}^2finkitxmlcodeinlinelatexdvp.

Soit kitxmlcodeinlinelatexdvp(A,B)\in\mathbb{R}_+^2finkitxmlcodeinlinelatexdvp,

kitxmlcodelatexdvp\begin{align} \int_{-A}^A\int_{-B}^B|f*g|dx\,dy & \le \int_{-A}^A\int_{-B}^B \left ( \iint_{\mathbb{R}^2} |f(u,v)g(x-u,y-v)|\,du\,dv\right )\,dx\,dy \\ \\ & \le \iint_{\mathbb{R}^2} \int_{-A}^A\int_{-B}^B |g(x-u,y-v)|\,dx\,dy\ f(u,v)\,du\,dv \\ \\ & \le \iint_{\mathbb{R}^2} \iint_{\mathbb{R}^2} |g(x,y)\,dx\,dy|\ |f(u,v)\,du\,dv \\ \\ & = \iint_{\mathbb{R}^2}|f|\cdot\iint_{\mathbb{R}^2}|g| \end{align}finkitxmlcodelatexdvp

Donc kitxmlcodeinlinelatexdvpf*gfinkitxmlcodeinlinelatexdvp est intégrable sur kitxmlcodeinlinelatexdvp\mathbb{R}^2finkitxmlcodeinlinelatexdvp.

On peut démontrer le gros de la preuve comme suit :

kitxmlcodelatexdvp\begin{align} \mathcal{F}(f*g)(x,y) & = \iint_{\mathbb{R}^2}(f*g)(u,v)\exp(-2i\pi(ux+vy))\,du\,dv \\ & = \iint_{\mathbb{R}^2}\left( \iint_{\mathbb{R}^2} f(a,b)g(u-a,v-b)\,da\,db \right )\exp(-2i\pi(ux+vy))\,du\,dv \\ & = \iint_{\mathbb{R}^2} \left( \iint_{\mathbb{R}^2} g(u-a,v-b)\exp(-2i\pi(ux+vy))\,du\,dv \right ) f(a,b)\,da\,db \\ & = \iint_{\mathbb{R}^2}\left ( \iint_{\mathbb{R}^2} g(u,v)\exp(-2i\pi(u+a)x+(v+b)y)\,du\,dv \right ) f(a,b)\,da\,db \\ & = \left( \iint_{\mathbb{R}^2}g(u,v)\exp(-2i\pi(uv+vy))\,du\,dv \right ) \left( \iint_{\mathbb{R}^2}f(a,b)\exp(-2i\pi(ax+by))\,da\,db \right ) \\ & = \mathcal{F}(f)(x,y)\cdot\mathcal{F}(g)(x,y) \end{align}finkitxmlcodelatexdvp

IV-A-2. Preuve dans l'espace discret

Soit kitxmlcodeinlinelatexdvp(a, b) \in [[-n/2,n/2]]^2finkitxmlcodeinlinelatexdvp.

kitxmlcodelatexdvp\begin{align} \mathcal{F}(h*f)(a,b) & = \sum_{i=-n/2}^{n/2}\sum_{j=-n/2}^{n/2}(h*f)(i,j)\cdot\exp\left( -2i\pi\left(\frac{ai}{n} + \frac{bj}{n} \right )\right ) \\ & = \sum_{i=-n/2}^{n/2}\sum_{j=-n/2}^{n/2} \left( \sum_{i'=-n/2}^{n/2}\sum_{j'=-n/2}^{n/2} h(i',j')f(i-i',j-j')\right ) \exp\left( -2i\pi\left(\frac{ai}{n} + \frac{bj}{n} \right )\right ) \\ & = \sum_{i'=-n/2}^{n/2}\sum_{j'=-n/2}^{n/2} h(i',j')\left( \sum_{i=-n/2}^{n/2}\sum_{j=-n/2}^{n/2} f(i-i',j-j') \exp\left(-2i\pi\left(\frac{ai}{n} + \frac{bj}{n} \right )\right )\right ) \\ & = \sum_{i'=-n/2}^{n/2}\sum_{j'=-n/2}^{n/2} h(i',j')\left( \sum_{v_1=-n/2+i'}^{n/2+i'}\sum_{v_2=-n/2+j'}^{n/2+j'} f(v_1,v_2)\exp\left(-2i\pi\left(\frac{a(v_1+i')}{n} + \frac{b(v_2+j')}{n}\right ) \right )\right) \\ & = \sum_{i'=-n/2}^{n/2}\sum_{j'=-n/2}^{n/2} \left( h(i',j')\exp\left(-2i\pi\left(\frac{ai'}{n}+\frac{bj'}{n} \right ) \right ) \right ) \\ & \cdot \sum_{v_1=-n/2+i'}^{n/2+i'}\sum_{v_2=-n/2+j'}^{n/2+j'} \left(f(v_1, v_2)\exp(-2i\pi\left(\frac{av_1}{n}+\frac{bv_2}{n} \right )) \right ) \end{align}finkitxmlcodelatexdvp

Grâce à la périodicité de kitxmlcodeinlinelatexdvpffinkitxmlcodeinlinelatexdvp, on peut réorganiser la deuxième somme de sorte que l'on ait :

kitxmlcodelatexdvp\begin{align} \mathcal{F}(h*f)(a,b) & = \sum_{i'=-n/2}^{n/2}\sum_{j'=-n/2}^{n/2}\left( h(i',j')\exp\left(-2i\pi\left(\frac{ai'}{n}+\frac{bj'}{n} \right )\right ) \right ) \\ & \cdot \sum_{v_1=-n/2}^{n/2}\sum_{v_2=-n/2}^{n/2}\left( f(v_1,v_2)\exp\left(-2i\pi\left(\frac{av_1}{n}+\frac{bv_2}{n} \right )\right ) \right ) \\ & = \mathcal{F}(h)(a,b)\cdot \mathcal {F}(f)(a,b) \end{align}finkitxmlcodelatexdvp

Vous avez aimé ce tutoriel ? Alors partagez-le en cliquant sur les boutons suivants : Viadeo Twitter Facebook Share on Google+   

  

Les sources présentées sur cette page sont libres de droits et vous pouvez les utiliser à votre convenance. Par contre, la page de présentation constitue une œuvre intellectuelle protégée par les droits d'auteur. Copyright © 2017 Humbert Florent. Aucune reproduction, même partielle, ne peut être faite de ce site et de l'ensemble de son contenu : textes, documents, images, etc. sans l'autorisation expresse de l'auteur. Sinon vous encourez selon la loi jusqu'à trois ans de prison et jusqu'à 300 000 € de dommages et intérêts.