On peut dores et déjà avancer que la connaissance de l'origine d'une image et donc des propriétés des bruits associés qui en corrompent les informations, est un atout permettant de concevoir des techniques de filtrage adaptées à chaque situation. Toutefois, la recherche d'un filtre universel, bien qu'encore illusoire, n'est pas abandonnée, tant les besoins sont nombreux, divers et souvent complexes.
\section{Modèle d'image bruitée}
-On considère qu'une image observée est un ensemble de $N$ observations sur un domaine $\Omega$ à deux dimensions ($\Omega \subset \mathbb{Z}^2$). À chaque élément de $\Omega$, aussi appelé \textit{pixel}, est associé un indice unique $k \in [\![1;N]\!]$, une position $x_k=(i,j)_k \in\Omega$ et une valeur observée $v_k=v(i,j)_k$.
+On considère qu'une image observée, de largeur $L$ pixels et de hauteur $H$ pixels, est un ensemble de $N=HL$ observations sur un domaine $\Omega$ à deux dimensions ($\Omega \subset \mathbb{Z}^2$). À chaque élément de $\Omega$, aussi appelé \textit{pixel}, est associé un indice unique $k \in [\![1;N]\!]$, une position $x_k=(i,j)_k \in\Omega$ et une valeur observée $v_k=v(i,j)_k$.
La valeur observée peut, selon les cas, être de dimension $1$ pour les images représentées en niveaux de gris ou de dimension 3 pour les images couleur représentées au format RVB. Les dimensions supérieures, pour la représentation des images hyperspectrales n'est pas abordé.
L'image observée peut ainsi être considérée comme un vecteur à $N$ éléments $\bar{v}= (v_k)_{k\in [\![1;N]\!]}$.
-Les divers traitements appliqués aux images observées ont souvent pour but d'accéder aux informations contenues dans une image sous-jacente, débarrassée de toute perturbation, dont nous faisons l'hypothèse qu'elle partage le même support $\Omega$ et que nous notons $\bar{u}$.
+Les divers traitements appliqués aux images observées ont souvent pour but d'accéder aux informations contenues dans une image sous-jacente, débarrassée de toute perturbation, dont nous faisons l'hypothèse qu'elle partage le même support $\Omega$ et que nous notons $\bar{u}$. L'estimation de $\bar{u}$ réalisée par ces traitements est notée $\widehat{\bar{u}} = (\widehat{u}_k)_{k\in [\![1;N]\!]}$.
Le lien entre $\bar{u}$ et $\bar{v}$ peut être exprimé généralement par la relation $\bar{v}=\bar{u}+\sigma\epsilon$, où $\epsilon \in \mathbb{R}^N$ représente le modèle de perturbation appliquée à $\bar{u}$ et $\sigma$ représente la puissance de cette perturbation qui a mené à l'observation de $\bar{v}$.
Dans le cas général, $\epsilon$ dépend de $\bar{u}$ et est caractérisé par la densité de probabilité (PDF pour probability density function) $p(v|u)$.
\[ p(v \mid u)=\mathrm{e}\frac{u^v}{v!}\]
\section{Les techniques de réduction de bruit}
-La très grande majorité des algorithmes de réduction de bruit fait l'hypothèse que la perturbation est de type gaussien, même si le développement des systèmes d'imagerie radar et médicale a poussé les chercheurs vers l'étude des bruits multiplicatifs du type \textit{speckle} ou \textit{Poisson}.
-Un très grand nombre de travaux proposant des méthodes de réduction de ces bruits ont été menés, ainsi que beaucoup d'états de l'art et d'études comparatives de ces diverses techniques, que nous n'avons pas la prétention d'égaler.
+La très grande majorité des algorithmes de réduction de bruit fait l'hypothèse que la perturbation est de type gaussien, même si le développement des systèmes d'imagerie radar et médicale a favorisé l'étude des bruits multiplicatifs du type \textit{speckle} ou \textit{Poisson}.
+Un très grand nombre de travaux proposant des méthodes de réduction de ces bruits ont été menés, ainsi que beaucoup d'états de l'art et d'études comparatives de ces diverses techniques, que nous n'avons pas l'ambition d'égaler.
-Les techniques et implémentations que nous aborderons dans le chapitre suivant sont celles qui ont un lien direct avec les travaux que nous avons menés. Nous
-présenterons donc les principales classes d'algorithmes de réduction de bruit et les implémentations GPU qui leur ont été consacrées.
+Nous nous focaliserons sur les techniques en lien avec les travaux que nous avons menés et qui ont donné lieu à des implémentations efficaces susceptibles de fournir des éléments opérationnels rapides pour le prétraitement des images.
+
+La figure \ref{fig-ny-noises} montre une image de synthèse issue de la base de test COIL \ref{adresse}, supposée sans bruit et qui sera considérée comme référence, ainsi que deux versions bruitées, respectivement avec un bruit gaussien d'écart type 25 et un bruit impulsionnel affectant 25\% des pixels.
+\begin{figure}
+ \centering
+ \subfigure[Sans bruit]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/ny256.png}}
+ \subfigure[Bruit gaussien $\sigma=25$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/ny256_gauss25.png}}
+ \subfigure[Bruit impulsionnel 25\%]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/ny256_sap25.png}}
+ \caption{Images 256$\times$256 en niveau de gris 8 bits utilisées pour l'illustration des propriétés des filtres. a) l'image de référence non bruitée. b) l'image corrompue par un bruit gaussien d'écart type $\sigma=25$. c) l'image corrompue par un bruit impulsionnel à 25\%.}
+\label{fig-ny-noises}
+\end{figure}
+
+\subsection{Les opérateurs de base}
+\subsubsection{Les algorithmes de voisinage}
+L'opération la plus employée dans les procédés de traitement d'image est sans doute la convolution. Selon les valeurs affectées aux coefficients du masque, le filtrage par convolution permet de réaliser bon nombre de traitements comme la réduction de bruit par moyennage ou noyau gaussien ou encore la détection de contours.
+Si la fonction définissant le masque de convolution est notée $h$, l'expression générale de la valeur estimée de pixel de coordonnées $(i,j)$ est donnée par
+\begin{equation}
+\widehat{u}(x, y) = \left(\bar{v} * h\right) = \sum_{(i < H)} \sum_{(j < L)}v(x-j, y-i)h(j,i)
+\label{convoDef}
+\end{equation}
+Dans les applications les plus courantes, $h$ est à support borné et de forme carrée et l'on parle alors de la taille du masque pour évoquer la dimension du support.
+ La figure \ref{fig-ny-convo} présente les résultats de la convolution par deux masques \textit{moyenneurs} $h_3$ et $h_5$ de taille différentes, appliqués à l'image corrompue par un bruit gaussien : on voit la diminution des fluctuations mais aussi le flou apporté et qui rend les contours d'autant moins définis que la taille du masque est grande. La troisième image montre le résultat de la convolution de l'image de référence par un masque \textit{dérivateur} $h_{dx}$ selon l'axe horizontal. On y constate la mise en évidence, incomplète, des contours.
+Les matrices définissant les masques sont les suivantes :
+
+\[h_3=\frac{1}{9}\begin{bmatrix}1&1&1\\1&1&1\\1&1&1\end{bmatrix}, h_{25}=\frac{1}{25}\begin{bmatrix}1&1&1&1&1\\1&1&1&1&1\\1&1&1&1&1\\1&1&1&1&1\\1&1&1&1&1\end{bmatrix}, h_{dx}= \begin{bmatrix}0&0&0\\-1&1&0\\0&0&0\end{bmatrix}\]
+
+\begin{figure}
+ \centering
+ \subfigure[Moyenneur 3$\times$3]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/convo/ny256_gauss25_moy3.png}}
+ \subfigure[Moyenneur 5$\times$5]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/convo/ny256_gauss25_moy5.png}}
+ \subfigure[Dérivateur horizontal]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/convo/ny256_2dx.png}}
+\caption{.}
+\label{fig-ny-convo}
+\end{figure}
+
+Le filtrage médian \ref{médian\_tukey} est également une opération très employée en prétraitement pour sa simplicité et ses propriétés de préservation des contours alliées à une capacité de réduction de bruit gaussien importante.
+La valeur du niveau de gris de chaque pixel est remplacée par la médiane des niveaux de gris des pixels voisins. Un des intérêts de ce filtre réside dans le fait que la valeur filtrée est une des valeurs du voisinage, contrairement à ce qui se produit lors d'une convolution. Un autre est de bien filtrer les valeurs extrêmes et par conséquent de trouver naturellement son application dans la réduction du bruit impulsionnel.
+Toutefois, la non-linéraité de cette technique et sa complexité n'en ont pas fait un filtre très utilisé jusqu'à ce que des implémentation efficaces soient proposées, en particulier le filtre à temps de calcul ``constant'' décrit dans \ref{medianO(1)}. Il est à noter que le filtrage médian est souvent appliqué en plusieurs passes de voisinage restreint.
+La figure \ref{fig-ny-median} montre la réduction de bruit impulsionnel obtenu grâce au filtre médian, dans trois conditions distinctes : median 3$\times$3 en une ou deux passes, puis médian 5$\times$5.
+\begin{figure}
+ \centering
+ \subfigure[Médian 3$\times$3 une passe]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/median/ny256_sap25_med3.png}}
+ \subfigure[Médian 3$\times$3 deux passes]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/median/ny256_sap25_med3x2.png}}
+ \subfigure[Médian 5$\times$5 une passe]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/median/ny256_sap25_med5.png}}
+\caption{Réduction du bruit impulsionnel par filtre médian.}
+\label{fig-ny-median}
+\end{figure}
+
+Le filtre bilatéral \ref{bilatéral\_filter} est une composition d'opérations que l'on peut voir comme un filtre de convolution dont les coefficients ne dépendraient pas uniquement de la position du pixel courant par rapport au pixel central, mais également de la différence de leurs intensités (cas des images en niveaux de gris).
+Si l'on note $\Omega_k$ le voisinage du pixel d'indice $k$, l'expression générale du niveau de gris estimé est donnée par
+\[\widehat{u_k}=\displaystyle\frac{\sum_{p\in \Omega_k}\left(F_S(x_p, x_k)F_I(v_p, v_k)v_p\right)}{\sum_{p\in\Omega_k }\left(F_S(x_p, x_k)F_I(v_p, v_k)\right)} \]
+où $F_S$ et $F_I$ sont les fonctions de pondération spatiale et d'intensité. Classiquement, $F_S$ et $F_I$ sont des gaussiennes de moyennes nulles et d'écarts type $\sigma_S$ et $\sigma_I$.
+Ce filtre se prête également bien à une utilisation en plusieurs passes sans flouter les contours. Des approximations séparables du filtre bilatéral, comme celle proposée dans \ref{bilateral-sep}, permettent d'obtenir des vitesses d'exécution plus élevées que les versions standard. Une variante à temps de calcul constant à même été proposée en 2008 par Porikli \ref{dans bilateral-sep ref 1 porikli}.
+Ce filtre permet un bon niveau de réduction de bruit gaussien, mais au prix d'un nombre de paramètres plus élevé à régler, ce qu'illustre la figure \ref{fig-ny-bilat} où le filtrage de la même image a été réalisé avec 9 combinaisons de $\sigma_S$ et $\sigma_I$.
+\begin{figure}
+ \centering
+\subfigure[$\sigma_S=1.0$ et $\sigma_I=0.1$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_1_01.png}}
+\subfigure[$\sigma_S=1.0$ et $\sigma_I=0.5$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_1_05.png}}
+\subfigure[$\sigma_S=1.0$ et $\sigma_I=1.0$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_1_1.png}}\\
+\subfigure[$\sigma_S=2.0$ et $\sigma_I=0.1$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_2_01.png}}
+\subfigure[$\sigma_S=2.0$ et $\sigma_I=0.5$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_2_05.png}}
+\subfigure[$\sigma_S=2.0$ et $\sigma_I=1.0$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_2_1.png}}\\
+\subfigure[$\sigma_S=5.0$ et $\sigma_I=0.1$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_5_01.png}}
+\subfigure[$\sigma_S=5.0$ et $\sigma_I=0.5$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_5_05.png}}
+\subfigure[$\sigma_S=5.0$ et $\sigma_I=1.0$]{\includegraphics[width=4cm]{/home/zulu/Documents/these_gilles/THESE/codes/bilat/ny_gauss25_bilat_5_1.png}}\\
+\caption{Réduction de bruit gaussien par filtrage bilatéral de voisinage 5$\times$5. $\sigma_S$ et $\sigma_I$ sont les écarts type des fonctions gaussiennes de pondération spatiale et d'intensité.}
+\label{fig-ny-bilat}
+\end{figure}
+
+
+Beaucoup d'autres algorithmes basés sur des moyennes locales efféctuées sur des voisinages de formes diverses, variables et/ou adaptatives afin de sélectionner le plus finement possible les pixels pris en compte dans le calcul de la valeur filtrée.
+Le principal défaut de ces techniques dites de réduction de variance est de générer des aplats dans les zones homogènes et des marches d'escalier dans les zones de transition douce (staircase effect), ces dernières pouvant être considérablement atténuées comme il a été montré dans \ref{staircase-effect}.
+
+\subsubsection{Les algorithmes par dictionnaire}
+Il s'agit ici de décrire l'image à débruiter en utilisant une base de fonctions permettant de décomposer l'image en une combinaison linéaire des éléments de cette base. Les bases les plus employées sont les ondelettes \ref{mallat2009-deladallep15, daubechie} ainsi que les fonctions sinusoïdales (DCT \ref{irfu}). Les éléments de la base peuvent être prédéterminés ou bien calculés à partir des données de l'image, par exemple en s'appuyant sur une analyse en composantes principales \ref{Aharon-2006 deladallep67}.
+
+\subsection{Les techniques avancées}
+Les techniques de réduction de bruit les plus efficaces sont aujourd'hui celles qui reposent sur les propriétés d'auto-similarité ds images, on les appelles aussi les techniques par patchs. L'idée principale est, comme pour les techniques classiques à base de de voisinage, de rechercher un ensemble de pixels pertinents et comparables afin d'en faire une moyenne. Cependant, dans le cas des techniques à patchs, la recherche de cet ensemble ne se limite pas à un voisinage du pixel central, mais fait l'hypothèse qu'il existe des zones semblables au voisinage du pixel central, réparties dans l'image et pas nécessairement immédiatement contigues.
+Le moyennage s'effectue alors sur l'ensemble des ces zones identifiées.
+L'algorithme des Non-Local Means (NL-means \ref{nl-means}) fut le premier de cette lignée à être proposé, mais bien d'autres suivirent comme le BM3D et ses variantes qui représentent actuellement l'état de l'art en terme de qualité de débruitage \ref{bm3D}.
+ Les différences entre ces algorithmes résident essentiellement dans la méthode de recherche et d'identification des patchs similaires, incluant la possiblité de forme et taille variables. Une telle recherche est d'autant plus coûteuse en temps de calcul qu'elle est effectuée sur une zone étendue autour du pixel central et cela représente le principal inconvénient de ces techniques qui peuvent présenter des temps d'exécution prohibitifs dans l'optique d'un traitement en temps réel.
+
+\section{Les implémentations GPU des algorithmes de filtrage}
+Le fabricant de processeurs graphiques Nvidia, seul type d'équipements dont nous disposons, fournit des implémentations performantes de certains prétraitements et algorithmes de filtrage.
+%TODO
+%Ajouter qq mots sur FFT, DCT utilisés dans irfu et que nous n'avons pas cherché à améliorer.
+C'est en particulier le cas de la convolution qui a fait l'objet d'une étude et d'une optimisation poussées pour déterminer la combinaison de solutions apportant la plus grande vitesse d'exécution \ref{convolution-soup-gtc09}. L'étude a testé 16 versions distinctes, chacune présentant une optimisation particulière quant-à l'organisation de la grille de calcul, aux types de transferts entre l'hôte et le GPU ainsi qu'au types de mémoire employé pour le calcul sur le GPU. Les résultats montrent que l'emploi de texture comme mémoire principale pour le stockage des images à traiter apporte un gain d'environ 50\% par rapport à l'utilisation de la mémoire globale. Par ailleurs, les transactions par paquets de 128 bits apportent également une amélioration sensible, ainsi que l'emploi de la mémoire partagée comme zone de travail pour le calcul des valeurs de sortie. Le traitement de référence effectué pour les mesures est la convolution générique (non séparable) d'une image 8 bits de 2048$\times$2048 pixels par un masque de convolution de 5$\times$5 pixels, expression que l'on raccourcira déronavant en \textit{convolution 5$\times$5}.
+Le meilleur résultat obtenu dans les conditions détaillées précédemment, sur architecture GT200 (carte GTX280) est de 1.4~ms pour le calcul, ce qui réalise un débit global de 945~MP/s lorsque l'on prend en compte les temps de transfert aller et retour des images (1.5~ms d'après nos mesures).
+Nous continuerons d'utiliser cette mesure de débit en \textit{Pixels par seconde} pour toutes les évaluations à venir ; elle permet en particulier de fournir des valeurs de performance indépendantes de la taille des images soumises au traitement.
+
+Le filtre médian n'a pas fait l'objet d'autant de publications, peut-être en raison des implémentations CPU performantes et génériques que l'on a déjà évoquées \ref{median0(1)}. Néanmoins, une bibliothèque commerciale (LibJacket, ArrayFire) en propose une implémentation GPU dont nous avons pu mesurer les performances pour un masque de 3$\times$3 et qui est également prise comme référence par Sanchez \textit{et al.} pour évaluer les performances de leur propre implémentation appelée PCMF \ref{median sanchez x2}. Sur architecture GT200 (GTX260), les performances maximales de ces deux versions sont obtenues pour un masque de 3$\times$3 pixels avec respectivement 175~MP/s pour libJacket et 60~MP/s pour PCMF.
+La figure \ref{fig-compare-jacket-pcmf}, tirée de \ref{median sanchez}, montre que le débit permis par la libJacket décroit très vite avec la taille du masque pour passer à 30~MP/s dès la taille 5$\times$5, alors que le PCMF décroit linéairement jusqu'à la taille 11$\times$11 où il permet encore de traiter quelque 40~MP/s.
+Plus récemment, Sanchez \textit{et al.} ont actualisé leurs mesures sur architecture Fermi (GPU C2075) en comparant leur PCMF à la version ré-écrite en C de libJacket, nommée ArrayFire. Les courbes sont celles de la figure \ref{fig-compare-arrayfire-pcmf}, où l'on constate que les variations selon la taille du masque demeurent comparables, avec toutefois des valeurs de débit augmentées, avec près de 185~MP/s pour ArrayFire et 82~MP/s pour PCMF. Il faut aussi noter que certains codes sont plus performants sur l'ancienne architecture GT200 que sur la plus récente Fermi ; c'est le cas pour l'implémentation du médian incluse dans la bibliothèque ArrayFire et nous reviendrons sur les raisons de cette perte de performances constatée au passage à une architecture plus récente dans le chapitre consacré à notre implémentation du filtre médian.
+
+Le filtre bilatéral a été plus abordé et un certain nombre de publications font état d'implémentations véritablement rapides. Il est néanmoins parfois difficile de les comparer sans disposer des codes sources, en raison de conditions de test très variables, en particulier en ce qui concerne le modèle de GPU et la taille du masque . Ceci étant précisé, on peut prendre comme référence initiale la version proposée par Nvidia dans le SDK CUDA et nommée ``ImageDenoising''. Elle permet d'exécuter sur GPU GTX480 un filtre bilatéral 7$\times$7 sur une image, déjà en mémoire GPU, de 1~MPixels en 0.411~ms.
\section{Les techniques de segmentation}
\begin{figure}
\centering
- \subfigure[Image initiale comportant deux zones : le fond et le cochon (la cible)]{\label{fig-histo-cochon-a} \includegraphics[height=3cm]{/home/zulu/Documents/THESE/codes/cochon256.png}}\quad
- \subfigure[Histogramme des niveaux de gris]{\label{fig-histo-cochon-b} \includegraphics[height=3cm]{/home/zulu/Documents/THESE/codes/seg_histogramme/histo-cochon256.png}}\quad
- \subfigure[Image binaire représentant la segmentation. Seuil estimé à 101 après 4 itérations.]{\label{fig-histo-cochon-c} \includegraphics[width=3cm]{/home/zulu/Documents/THESE/codes/seg_histogramme/cochon256-seghisto-101-255.png}}\\
-\subfigure[Image initiale bruitée]{\label{fig-histo-cochon-d} \includegraphics[height=3cm]{/home/zulu/Documents/THESE/codes/cochon256-sig25.png}}\quad
- \subfigure[Histogramme des niveaux de gris]{\label{fig-histo-cochon-e} \includegraphics[height=3cm]{/home/zulu/Documents/THESE/codes/seg_histogramme/histo-cochon256-sig25.png}}\quad
- \subfigure[Image binaire représentant la segmentation. Seuil estimé à 99 après 5 itérations.]{\label{fig-histo-cochon-f} \includegraphics[height=3cm]{/home/zulu/Documents/THESE/codes/seg_histogramme/cochon256-sig25-seghisto-99-255.png}}
- \caption{Segmentation par analyse simple d'histogramme. Colonne de gauche : image d'entrée. Colonne centrale : histogramme des niveaux de gris. Colonne de droite : résultat de la segmentation.}
+ \subfigure[Image initiale comportant deux zones : le fond et le cochon (la cible)]{\label{fig-histo-cochon-a} \includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/cochon256.png}}\quad
+ \subfigure[Histogramme des niveaux de gris]{\label{fig-histo-cochon-b} \includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/seg_histogramme/histo-cochon256.png}}\quad
+ \subfigure[Image binaire représentant la segmentation. Seuil estimé à 101 après 4 itérations.]{\label{fig-histo-cochon-c} \includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/seg_histogramme/cochon256-seghisto-101-255.png}}\\
+\subfigure[Image initiale bruitée]{\label{fig-histo-cochon-d} \includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/cochon256-sig25.png}}\quad
+ \subfigure[Histogramme des niveaux de gris]{\label{fig-histo-cochon-e} \includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/seg_histogramme/histo-cochon256-sig25.png}}\quad
+ \subfigure[Image binaire représentant la segmentation. Seuil estimé à 99 après 5 itérations.]{\label{fig-histo-cochon-f} \includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/seg_histogramme/cochon256-sig25-seghisto-99-255.png}}
+ \caption{Segmentation d'une image en niveaux de gris de 128 $\times$ 128 pixels par analyse simple d'histogramme. Colonne de gauche : image d'entrée. Colonne centrale : histogramme des niveaux de gris. Colonne de droite : résultat de la segmentation.}
\label{fig-histo-cochon}
\end{figure}
Nous pouvons retenir que les premières d'entre elles, qui n'étaient pas spécifiquement dédiées à la segmentation d'images numériques mais au regroupement d'éléments répartis sur un domaine (1D ou 2D), ont été élaborées autour d'une mesure locale des liens basée sur la distance entre les éléments. La réduction du graphe est ensuite effectuée en utilisant un algorithme spécifique, comme le \textit{minimum spanning tree}, dont l'application a été décrite dès 1970 dans \ref{slac-pub-0672} et où il s'agit simplement de supprimer les liens \textit{inconsistants}, c'est à dire ceux dont le poids est significativement plus élevé que la moyenne des voisins se trouvant de chaque coté du lien en question.
L'extension a rapidement été faite aux images numériques en ajoutant l'intensité des pixels au vecteur des paramètres pris en compte dans l'évaluation du poids des liens.
D'autres critères de simplification ont aussi été élaborés, avec pour ambition de toujours mieux prendre en compte les caractéristiques structurelles globales des images pour prétendre à une segmentation qui conduise à une meilleure perception conceptuelle.
-Le principe général des solutions actuelles est proche de l'analyse en composantes principales appliquée à une matrice d'affinités qui traduit les liens entre les segments.
-On peut citer, par ordre chronologique, les méthodes reposant sur le \textit{graphe optimal} de Wu et Leahy \ref{wulealy_1993} et plus récemment \ref{cf_notes x5}. Le principal point faible de ces techniques réside essentiellement dans la difficulté à trouver un compromis acceptable entre identification de structures globales et préservation des éléments de détails. Cela se traduit dans la pratique par un ensemble de paramètres à régler pour chaque type de segmentation à effectuer.
+Le principe général des solutions actuelles est proche de l'analyse en composantes principales appliquée à une matrice de similarité qui traduit les liens entre les segments.
+Pour des images en niveaux de gris, l'expression générale des éléments $w_{ij}$ de la matrice de similarité $W$ est :
+\[w_{ij} =
+\begin{cases}
+\mathrm{e}^{\|v_i-v_j\|^2/\sigma_v^2}\mathrm{e}^{\|x_i-x_j\|^2/\sigma_x^2} & \text{si $\|x_i-x_j\|<r$}\\
+0 & \text{sinon}
+\end{cases}
+\]
+On construit ensuite la matrice de connectivité $D$, diagonale et dont les éléments sont :
+\[d_{i} = \displaystyle\sum_jw_{ij}\]
+
+Le système dont on cherche les valeurs propres $\lambda_k$ et les vecteurs propres associés $Y_k$ est alors le suivant :
+\[\left(D-W)Y=\lambda DY \]
+
+Parmi les méthodes reposant sur ce principe, on peut citer, par ordre chronologique, celles qui reposent sur le \textit{graphe optimal} de Wu et Leahy \ref{wulealy_1993} et plus récemment \ref{cf_notes x5}. Le principal point faible de ces techniques réside essentiellement dans la difficulté à trouver un compromis acceptable entre identification de structures globales et préservation des éléments de détails. Cela se traduit dans la pratique par un ensemble de paramètres à régler pour chaque type de segmentation à effectuer.
Cependant, elles sont employées dans les algorithmes de haut niveau les plus récents, comme nous le verrons plus loin.
-La figure \ref{fig-graph-cochon} montre un exemple de l'application de l'algorithme \textit{normalized cuts} décrit dans \ref{sm-ncuts}. Quatre paramètres sont à fixer pour réaliser la segmentation, ils correspondent respectivement aux facteurs de forme de la matrice d'affinités en termes de distance et de niveau de gris (2 paramètres), à la surface minimale admise pour un segment ainsi qu'au seuil de ségmentation des graphes. Le choix des valeurs pour les paramètres n'est pas immédiat et la figure illustre bien la difficulté de trouver un compromis acceptable. Enfin, les temps d'exécutions peuvent devenir très rapidement prohibitifs dès lors que l'on utilise des facteurs de forme étendus.
-
+La figure \ref{fig-graph-cochon} montre un exemple de l'application de l'algorithme \textit{normalized cuts} décrit dans \ref{sm-ncuts pami2000} et implémenté par Cour, Yu et Shi en 2004. Cette implémentation utilise des valeurs pré-établies des paramètres de calcul de la matrice de similarité produisant de bonnes segmentations d'objets et/ou personnes dans les images naturelles, mais requiert de prédéterminer le nombre de segments à obtenir. Les images de la figure représentent les résultats obtenus avec un nombre de segments variant de 2 à 5 et montrent qu'il difficile de trouver un compromis acceptable. Enfin, les temps d'exécutions peuvent devenir très rapidement prohibitifs, même avec des implémentations plus optimisées. Pour information, les résultats de la figure \ref{fig-graph-cochon} ont été obtenus en 1.5~s environ (Matlab R2010 sur CPU intel core i5-2520M @ 2.50GHz - linux 3.2.0)
+\begin{figure}
+ \centering
+ \subfigure[$s = 2$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/graphe/cochon128_ncuts_2seg.png}}
+ \subfigure[$s = 3$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/graphe/cochon128_ncuts_3seg.png}}
+ \subfigure[$s = 4$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/graphe/cochon128_ncuts_4seg.png}}
+ \subfigure[$s = 5$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/graphe/cochon128_ncuts_5seg.png}}
+ \caption{Segmentation d'une image en niveaux de gris de 128 $\times$ 128 pixels par simplification de graphe de type \textit{Normalized cut} pour un nombre $s$ de segments variant de 2 à 5.}
+\label{fig-graph-cochon}
+\end{figure}
+
\subsection{kernel-means, mean-shift et dérivés}
Parallèlement à la réduction de graphes, d'autres approches ont donné naissance à une multitude de variantes tournées vers la recherche des moindres carrés.
Toutefois, vraisemblablement du fait de sa simplicité d'implémentation et de temps d'exécution rapides, la communauté scientifique s'est beaucoup penchée sur cette méthode pour en compenser les défauts, jusqu'à en faire une des plus employées, en particulier par les statisticiens.
On compte aussi beaucoup de variantes telles les \textit{k-centers} \ref{k_centers} et les \textit{k-médians} \ref{k_medians} qui n'employent pas la moyenne arithmétique comme expression du ``centre'' d'un segment.
Des solutions ont aussi été apportées pour l'estimation de $k$ en employant, par exemple, un critère de vraisemblance pour choisir la meilleure valeur de $k$ dans un intervalle donné \ref{x-means}.
+À titre d'illustration et de comparaison, l'image du cochon a été traitée par une implémentation naïve de l'algorithme original des \textit{k-means} en donnant successivement au nombre de segments les valeurs $s=2$ à $s=5$. Les résultats sont reproduits à la figure \ref{fig-kmeans-cochon} et montrent encore une fois l'influence de $s$ sur la segmentation.
+\begin{figure}
+ \centering
+ \subfigure[$s = 2$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/kmeans/cochon128_kmeans_2seg.png}}
+ \subfigure[$s = 3$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/kmeans/cochon128_kmeans_3seg.png}}
+ \subfigure[$s = 4$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/kmeans/cochon128_kmeans_4seg.png}}
+ \subfigure[$s = 5$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/kmeans/cochon128_kmeans_5seg.png}}
+ \caption{Segmentation d'une image en niveaux de gris de 128 $\times$ 128 pixels par algorithme \textit{k-means} pour un nombre $s$ de segments variant de 2 à 5. Chaque couleur est associée à un segment. Les couleurs sont choisies pour une meilleure visualisation des différents segments.}
+\label{fig-kmeans-cochon}
+\end{figure}
Un algorithme initiallement proposé en 1975 par Fukunaga et Hostetler \ref{Lestimation_html} permet de manière plus générique de déterminer le nombre de segments, ou modes, ainsi que les points, ou pixels, qui les composent. Il cherche pour ce faire à localiser les $k$ positions ou le gradient de densité s'annule.
Il utilisé un voisinage pondére (ou \textit{kernel}) et détermine le centre de masse des segments en suivant itérativement le gradient de densité dans le voisinage autour de chaque élément du domaine. Lorsque l'algorithme à convergé, les $k$ segments sont identifiés et continennent chacun l'ensemble des points qui ont conduit à leur centre de masse respectif.
Étonnement, malgré ses qualités intrinsèques, cet algorithme du \textit{mean-shift} est resté longtemps sans susciter de grand intérêt, jusqu'à l'étude de Cheng \ref{meanshift_1995} qui en a demontré les propriétés et établi les lien avec d'autres techniques d'optimisation commme la descente/montée de gradient ou de filtrage commme le floutage.
Comaniciu et Peer ont alors étendu l'étude et proposé une application à la segmentation en utilisant l'espace colorimétrique CIELUV \ref{Computer Graphics by Foley, van Dam, Feiner, and Hughes, published by Addison-Wesley, 1990} et montré qu'elle permettait une meilleure identification des modes de l'image \ref{mean_shift 1999 2002}.
+Une implémentation de la variante proposée par Keselman et Micheli-Tzanakou dans \ref{yket1999} appliquée à notre image de test fournit les résultats reproduits à la figure \ref{fig-meanshift-cochon}. Pour se rapprocher des traitements précédents, nous avons identifié, par essais successifs, les tailles de voisinage conduisant à des nombre de segments identiques à ceux des figures précedentes (de 2 à 5). Le volume minimal admis pour un segment à été arbitrairement fixé à 100 pixels.
+\begin{figure}
+ \centering
+ \subfigure[$r=100 \Rightarrow s = 2$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/meanshift/cochon128_meanshift_r100m100.png}}
+ \subfigure[$r=50 \Rightarrow s = 3$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/meanshift/cochon128_meanshift_r50m100.png}}
+\subfigure[$r=35 \Rightarrow s = 4$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/meanshift/cochon128_meanshift_r35m100.png}}
+ \subfigure[$r=25 \Rightarrow s = 5$]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/meanshift/cochon128_meanshift_r25m100.png}}
+ \caption{Segmentation d'une image en niveaux de gris de 128 $\times$ 128 pixels par algorithme \textit{mean-shift} pour un rayon de voisinage $r$ de 100, 50, 35 et 25 pixels permettant d'obtenir un nombre $s$ de segments variant respectivement de 2 à 5. Le volume minimal admis pour un segment est fixé à 100 pixels. Chaque couleur est associée à un segment. Les couleurs sont choisies pour une meilleure visualisation des différents segments.}
+\label{fig-meanshift-cochon}
+\end{figure}
-Il est à noter que les segmentations basées sur des algorithmes de \textit{clustering} nécessitent une phase supplémentaire de génération des frontières inter-segments et d'affectation de la valeur de chaque segment aux éléments qui le composent.
+Il est à noter que les segmentations basées sur des algorithmes de \textit{clustering} comme ceux que l'on vient de présenter nécessitent le plus souvent une phase supplémentaire de génération des frontières inter-segments et d'affectation de la valeur de chaque segment aux éléments qui le composent.
Par ailleurs, dans les deux cas du \textit{k-means} et du \textit{mean-shift}, chaque itération génère une réduction de la variance (due au moyennage) et on peut donc rapprocher ces techniques de celles de réduction de bruit par minimisation de variance.
\subsection{Les contours actifs, ou \textit{snakes}}
Contrairement aux précédentes techniques et comme leur nom le laisse deviner, les éléments constitutifs de ces méthodes sont cette fois des \textit{contours} et non plus des \textit{régions}. De fait, ils définissent nativement une segmentation de l'image.
-Le principe général est de superposer une courbe paramétrique à l'image, le \textit{snake}, puis de lui appliquer des déformations successives destinées à rapprocher le \textit{snake} des contours de l'objet. Les déformations à appliquer sont guidées par l'évaluation d'une fonction d'énergie prenant en compte :
+Le principe général est de superposer une courbe paramétrique $S$ à l'image, le \textit{snake}, puis de lui appliquer des déformations successives destinées à rapprocher le \textit{snake} des contours de l'objet. Les déformations à appliquer sont guidées par l'évaluation d'une fonction d'énergie $E_{snake}$ prenant en compte :
\begin{itemize}
-\item l'énergie interne de la courbe, fonction de son allongement de sa courbure.
-\item l'énergie externe liée à l'image, fonction de la proximité de la courbe avec les zones de fort gradient et éventuellement une contrainte fixée par l'utilisateur comme des points imposés par exemple.
+\item l'énergie interne $E_{int}$ de la courbe, fonction de son allongement de sa courbure.
+\item l'énergie externe $E_{ext}$ liée à l'image, fonction de la proximité de la courbe avec les zones de fort gradient et éventuellement une contrainte fixée par l'utilisateur comme des points imposés par exemple.
\end{itemize}
-Ici encore, la résolution du problème revient à minimiser une fonction d'énergie sous contrainte et les diverses techniques de résolution numérique peuvent s'appliquer comme pour les autres classes d'algorithmes itératifs présentés précédemment.
+L'expression générique peut alors s'écrire
+\[E_{snake} = E_{int}+E_{ext}\]
+où
+\[E_{int} = \displaystyle\sum_{s\in S} \frac{1}{2}\left(\alpha\left|\frac{\partial x_s}{\partial s}\right|^2
++\beta \left|\frac{\partial^2x_s}{\partial s^2}\right|\right)ds\]
+et
+\[E_{ext} = \displaystyle\sum_{s\in S} -\left|\nabla\left[G_{\sigma}(x_s)\ast v_s\right]\right|^2ds\]
+
+L'idée générale de l'algorithme du \textit{snake} est de trouver une courbe $S$ qui minimise l'énergie totale $E_{snake}$.
+Ici encore, la résolution du problème revient donc à minimiser une fonction sous contrainte et les diverses techniques de résolution numérique peuvent s'appliquer comme pour les autres classes d'algorithmes itératifs présentés précédemment, avec ici encore, un nombre de paramètres à régler assez important. Notons également que dans le cas général, les paramètres notés $\alpha$ et $\beta$, que l'on qualifie aussi d'élasticité et de raideur, sont aussi des fonctions de l'abscisse curviligne $s$. La fonction $G_{\sigma}$ est la fonction d'attraction aux forts gradients de l'image.
Dans sa version originale proposée par Kass \textit{et al.} en 1988 \ref{snake_kass_1988}, l'algorithme dit du \textit{snake} présente l'intérêt de converger en un nombre d'itérations assez réduit et permet de suivre naturellement un \textit{cible} en mouvement après une convergence initiale à une position donnée, chaque position de convergence fournissant une position initiale pertinente pour la position suivante.
-Toutefois, il se montre sensible à l'état initial de la courbe et requiert souvent de celle-ci qu'elle soit assez proche de l'objet à ``entourer'', sous peine de se verrouiller dans un minimum local.
+Toutefois, il se montre particulièrement sensible à l'état initial de la courbe et requiert souvent de celle-ci qu'elle soit assez proche de l'objet à ``entourer'', sous peine de se verrouiller dans un minimum local.
La sensibilité au bruit n'est pas non plus très bonne du fait de la formulation locale de l'énergie.
Les ``concavités'' étroites ou présentant un goulot d'étranglement marqué sont par ailleurs mal délimitées.
Enfin, la fonction d'énergie étant calculée sur la longueur totale de la courbe, cela pénalise la bonne identification des structures de petite taille vis à vis de la longueur totale de la courbe.
+La figure \ref{fig-snake-tradi-cochon} illustre ces défauts en montrant quelques états intérmédiaires ainsi que le résultat final d'une segmentation réalisée à partir d'un contour initial circulaire et des paramètres à valeurs constantes et réglés empiriquement, en employant la méthode du snake original.
+On voit que la convergence est assez rapide mais que le contour ainsi détérminé ne ``colle'' pas bien à l'objet que l'on s'attend à isoler.
+\begin{figure}
+ \centering
+\subfigure[Les états initial et suivant chacune des trois premières itérations]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/snake/cochon128_tradi_snake_it3.png}}
+\subfigure[L'état du contour après la septième itération]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/snake/cochon128_tradi_snake_it7.png}}
+\subfigure[L'état du contour après la dixième itération]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/snake/cochon128_tradi_snake_it10.png}}
+\subfigure[L'état du contour après la centième itération. C'est le contour final.]{\includegraphics[width=3cm]{/home/zulu/Documents/these_gilles/THESE/codes/snake/cochon128_tradi_snake_result.png}}
+\caption{Segmentation d'une image en niveaux de gris de 128 $\times$ 128 pixels par algorithme dit du \textit{snake}, dans sa version originale. Les paramètres d'élastictié, de raideur et d'attraction ont été fixés respectivement aux valeurs 5, 0.1 et 5. }
+\label{fig-snake-tradi-cochon}
+\end{figure}
+
Il est cependant possible de contrôler la finesse de la segmentation mais au prix de temps de calculs qui peuvent devenir très longs.
-Les variantes les plus intéressantes sont :
+Parmi les variantes élaborées qui tentent de pallier ces défauts, les plus intéressantes sont :
\begin{itemize}
\item le \textit{balloon snake}, conçu pour remédier au mauvais suivi des concavités en introduisant une force supplémentaire de pression tendant à \textit{gonfler} le snake jusqu'à ce qu'il rencontre un contour suffisamment marqué. Cela suppose toutefois que l'état initial de la courbe la situe entièrement à l'intérieur de la zone à segmenter et est surtout employé dans des applications semi-automatiques où l'utilisateur définit au moins une position et une taille initiales pour la courbe.
\item le \textit{snake} GVF (pour Gradient Vector Flow), dont le but est de permettre qu'une initialisation lointaine de la courbe ne pénalise pas la segmentation. Une carte des lignes de gradient est établie sur tout le domaine de l'image et sert à intégrer une force supplémentaire dans l'énergie totale, qui attire la courbe vers la zone de fort gradient.
\item les \textit{snake} orientés régions, qui visent essentiellement à mieux caractériser les zones à segmenter et améliorer la robustesse vis à vis du bruit en employant une formulation de l'énergie calculée sur le domaine complet de l'image \ref{cohenSMIE93, ronfard}. Les premiers résultats confirment la qualité de cette méthode, mais la nécessité d'effectuer les calculs sur l'image entière générait des temps de traitement prohibitifs jusqu'à ce que Bertaux \textitat{et al.} proposent une amélioration algorithmique exacte permettant à nouveau un calcul en 1D, le long de la courbe, moyennant une simple étape initiale générant un certain nombre d'images intégrales \ref{snake_bertaux}. La section \ref{sec_contrib_snake} qui introduit notre contribution à cette technique en donnera une description détaillée.
\end{itemize}
-
+% ne faut-il pas mieux éluder le paragraphe ci-dessous
\subsection{Méthodes hybrides}
-Aujourd'hui, les algorithmes de segmentation les plus performants en terme de qualité emploient des techniques qui tentent de tirer le meilleur parti de plusieurs des méthodes ``historiques'' décrites précédemment.
-Le meilleur exemple, et le seul que nous citerons, est le détecteur de contour et l'algorithme de segmentation associé proposé par Arbelaez \textit{et al.} en 2010 \ref{amfm_2010}. Il compose avec la constructions d'histogrammes locaux pour générer une matrice de similitude (affinity matrix) et appliquer les techniques liées à la théorie des graphes pour réduire la dimension de l'espace de représentation (calcul des valeurs et vecteurs propres). Il utilise ensuite une technique adaptée de \textit{ligne de partage des eaux} (que l'on aurait rangée avec les mean-shift) pour regrouper les segments.
+Aujourd'hui, les algorithmes de segmentation les plus performants en terme de qualité emploient des techniques qui tentent de tirer le meilleur parti de plusieurs des méthodes ``historiques'' décrites précédemment.
+Le meilleur exemple, et le seul que nous citerons, est le détecteur de contour et l'algorithme de segmentation associé proposé par Arbelaez \textit{et al.} en 2010 \ref{amfm_2010}. Il compose avec la constructions d'histogrammes locaux pour générer une matrice de similitude (affinity matrix) et appliquer les techniques liées à la théorie des graphes pour réduire la dimension de l'espace de représentation (calcul des valeurs et vecteurs propres). Il utilise ensuite une technique adaptée de \textit{ligne de partage des eaux} \ref{watershed} (que l'on aurait rangée avec les mean-shift) pour regrouper les segments.
Les résultats sont très bons et des implémentations efficaces ont dores et déjà été écrites (voir section \ref{sec_ea_gpu}.
-
+%TODO
+%peut-être dire deux mots sur le partage des eaux (avec kmeans et meanshift) puisqu'il est employé dans gPb
