-\section{Présentation de l'algorithme}
+\section{Introduction}
La principale difficulté soulevée par l'emploi d'algorithmes de type \textit{snake} orientés contour est le choix de la fonction d'énergie externe et la détermination de la nature des images auxquelles elle convient.
Dans l'approche orientée régions, les deux régions que sont l'extérieur et l'intérieur du contour (cas mono cible) sont prises en compte dans l'estimation de la forme du contour ; cela permet d'extraire des formes dans des images où les contours de la cible sont mal définis, en raison d'un fort niveau de bruit par exemple.
-Les algorithmes découlant de cette approche n'ont fait l'objet, à notre connaissance, d'aucune parallèlisation sur GPU, malgré le grand intérêt qu'elles revêtent dans l'interprétation d'images fortement bruitées ( RADAR, médicales,\dots ) et le besoin d'en réduire les temps d'exécution suffisamment pour permettre l'interactivité.
-Nous proposons dans la suite de ce chapitre de détailler tout d'abord l'algorithme séquentiel que nous avons pris comme référence, puis d'en présenter la version parallèle pour GPU que nous en avons conçu.
-L'algorithme a été décrit et proposé initialement en 1999 par Chesnaud \textit{et al.} dans \cite{ChesnaudRB99}. L'implémentation que les auteurs ont développé a continué d'être améliorée jusqu'à aujourd'hui et est employée comme brique élémentaire dans des algorithmes plus complexes. La version qui sert de référence ici est une implémentation séquentielle optimisée qui met aussi à profit les capacités de parallélisme des CPU actuels en employant le jeu d'instruction SSE2 des microprocesseurs. La description que nous en faisons dans les lignes qui suivent est très largement inspirée de \cite{ChesnaudRB99} à la différence que nous n'implémentons pas les critères de régularisation du contour ni de minimisation de la longueur de description pour nous focaliser sur la déformation du contour et sa convergence.
+Les algorithmes découlant de cette approche n'ont fait l'objet, à notre connaissance, d'aucune parallèlisation sur GPU, malgré le grand intérêt qu'elles revêtent dans l'interprétation d'images fortement bruitées ( imagerie ultrasonore, RADAR ) et le besoin de réduire suffisamment les temps de traitement pour permettre l'interactivité.
+
+Nous proposons dans la suite de ce chapitre de commencer par détailler l'algorithme séquentiel que nous avons pris comme référence, puis d'en présenter la version parallèle pour GPU que nous avons conçue.
+L'algorithme a été décrit et proposé initialement en 1999 par Chesnaud \textit{et al.} dans \cite{ChesnaudRB99}. L'implémentation que les auteurs ont développée a continué d'être améliorée jusqu'à aujourd'hui et est employée comme brique élémentaire dans des algorithmes plus complexes. La version qui sert de référence ici est une implémentation séquentielle optimisée et dont nous faisons ici une description très largement inspirée de \cite{ChesnaudRB99} à la différence que nous n'implémentons pas les critères de régularisation du contour ni de minimisation de la longueur de description pour nous focaliser sur la déformation du contour et sa convergence.
+\section{Présentation de l'algorithme}
\subsection{Formulation}
-À l'intérieur de l'image observée $\bar{v}$, soient $\bar{t}$ le vecteur composé par les niveaux de gris des $N_t$ pixels de la région cible $\Omega_t$ et $\bar{b}$ celui des $N_b$ pixels du fond $\Omega_b$. Les vecteurs $\bar{t}$ et $\bar{b}$ sont supposés non corrélés et sont caractérisés par leurs densités de probabilité (PDF) respectives $p^{\Theta_t}$ et $p^{\Theta_b}$ ; $\Theta_t$ et $\Theta_b$ étant les vecteurs des paramètres de leurs PDF. Dans le cas gaussien que nous supposerons ici, $\Theta = (\mu, \sigma)$ où $\mu$ est la moyenne et $\sigma^2$ est la variance.
+À l'intérieur de l'image observée $\bar{v}$, soient $\bar{t}$ le vecteur composé par les niveaux de gris des $N_t$ pixels de la région cible $\Omega_t$ et $\bar{b}$ celui des $N_b$ pixels du fond $\Omega_b$. Les vecteurs $\bar{t}$ et $\bar{b}$ sont supposés non corrélés et sont caractérisés par leurs densités de probabilité (PDF) respectives $p^{\Theta_t}$ et $p^{\Theta_b}$ ; $\Theta_t$ et $\Theta_b$ étant les vecteurs des paramètres de leurs PDF. Dans le cas gaussien que nous supposerons ici parcequ'il est considéré comme prépondérant dans les images naturelles, $\Theta = (\mu, \sigma)$ où $\mu$ est la moyenne et $\sigma^2$ est la variance.
On note $\Gamma$ le contour de la région cible ($\Gamma \in \Omega_t$), que l'on suppose continu en connexité à 8 voisins.
Le but de la segmentation est alors de déterminer la géométrie de $\Gamma$ qui maximise un critère de vraisemblance généralisée (GL).
-La vraisemblance sur l'ensemble de l'image, \textit{ie.} la région $\Omega$ est donnée par
+La vraisemblance sur l'ensemble de l'image, \textit{i.e} la région $\Omega$, est donnée par
\begin{equation}
P\left(\bar{v} | \Omega_t, \Omega_b, \Theta_t, \Theta_b\right) = P\left(\bar{v} | \Omega_t, \Theta_t\right)P\left(\bar{v}|\Omega_b, \Theta_b\right)
\right.
\end{equation}
-En intégrant \eqref{eq-teta} dans \eqref{eq-gl1}, il reste, à une constante près, le critère de vraisemblance généralisée suivant, noté GL, que l'on cherche à optimiser en déterminant la géométrie de $\Gamma$ qui en maximise la valeur et épousera alors au mieux la forme du contour de la cible.
+En insérant les expressions de \eqref{eq-teta} dans \eqref{eq-gl1}, il reste, à une constante près, le critère de vraisemblance généralisée suivant, noté GL, que l'on cherche à optimiser en déterminant la géométrie de $\Gamma$ qui en maximise la valeur et épousera alors au mieux la forme du contour de la cible.
\begin{equation}
\label{eq-gl}
GL=\frac{1}{2}\left(N_tln\left(\widehat{\sigma_t^2}\right)+N_bln\left(\widehat{\sigma_b^2}\right)\right)
\end{equation}
-\subsection{Optimisation des calculs}\label{snake-formulation}
+\subsection{Optimisation des calculs\label{snake-formulation}}
La maximisation de GL est effectuée en employant une technique itérative où sa valeur doit être calculée à chaque déformation du contour $\Gamma$.
Si l'on se reporte à l'équation \eqref{eq-teta}, on voit que l'obtention de la valeur de GL nécessite, à chaque évaluation d'une géométrie donnée de $\Gamma$, le calcul des sommes
\begin{align}
S_v(\Omega_b) &= \sum_{x_k\in \Omega_b} v_k & S_{v^2}(\Omega_b) &= \displaystyle\sum_{x_k\in \Omega_b} v_k^2 \nonumber\\
\end{align}
-Considèrons la région cible $\Omega_t$, les pixels de coordonées $(i,j)$ qui la composent, et généralisons l'écriture des sommes de \eqref{eq-sommes1} en
+Considérons la région cible $\Omega_t$, les pixels de coordonées $(i,j)$ qui la composent, et généralisons l'écriture des sommes de \eqref{eq-sommes1} en
\begin{equation}
\label{eq-sommes-gene}
- S_f(\Omega_t) = \sum_{i=i_{min}}^{i=i_{max}}\sum_{j=j_{min}(i)}^{j=j_{max}(i)}g\left(v(i,j)\right)
+ S_g(\Omega_t) = \sum_{i=i_{min}}^{i=i_{max}}\sum_{j=j_{min}(i)}^{j=j_{max}(i)}g\left(v(i,j)\right)
\end{equation}
-où $f$ représente la fonction de valeurs de niveaux de gris à sommer.
+où $g$ représente la fonction de valeurs de niveaux de gris à sommer.
En posant
\begin{equation}
L'équation \eqref{eq-sommes-gene} devient
\begin{equation}
\label{eq-somme-cumuls2}
- S_f(\Omega_t) = \sum_{i=i_{min}}^{i=i_{max}}\left[T_g(i,j_{max}(i))-T_g(i,j_{min}(i)-1)\right]
+ S_g(\Omega_t) = \sum_{i=i_{min}}^{i=i_{max}}\left[T_g(i,j_{max}(i))-T_g(i,j_{min}(i)-1)\right]
\end{equation}
qui représente une sommation sur le contour $\Gamma$ que l'on peut écrire
\begin{equation}
\end{equation}
La valeur de $C(i,j)$ est déterminée pour chaque pixel d'indice $l$ du contour en considérant les pixels d'indices $l-1$ et $l+1$ qui définissent les deux vecteurs $f_{in}$ et $f_{out}$ et leur code selon le codage de Freeman, comme l'illustre la figure \ref{fig-freeman}. La table \ref{tab-freeman} donne les valeurs de $C(i,j)$ selon les valeurs des codes de Freeman des vecteurs $f_{in}$ et $f_{out}$.
-Il faut noter que les valeurs dans la table \ref{tab-freeman} diffèrent de celles proposées initialement dans \cite{ChesnaudRB99}. Cette modification a été proposée plus tard pour permettre de s'adapter à la segmentation multi-cibles. Nous avons conservé la version la plus récente.
+Il faut noter que les valeurs dans la table \ref{tab-freeman} diffèrent de celles proposées initialement dans \cite{ChesnaudRB99}. Cette modification a été proposée pour permettre de s'adapter à la segmentation multi-cibles \cite{GallandBR03}.
\begin{figure}[htb]
\centering
\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/codage-freeman.png}
7 &0&0&0&0&0&-1&-1&-1\\
\bottomrule
\end{tabular}
- \caption{Valeur du coefficient $C(i,j)$ en fonction des valeurs des codes de Freeman des vecteurs $f_{in}$ et $f_{out}$.}
+ \caption{Valeur du coefficient $C(i,j)$ en fonction des valeurs des codes de Freeman des vecteurs $f_{in}$ et $f_{out}$ \cite{GallandBR03}.}
\label{tab-freeman}
\end{table}
L'intérêt de cette transformation est majeur :
\begin{itemize}
\item La sommation en deux dimensions sur la région $\Omega_t$ est ainsi réduite à une sommation à une dimension sur le contour $\Gamma$.
-\item Les valeurs $T_g(i,j)$ peuvent être calculées préalablement à la phase de segmentation proprement dite. Pour le cas gaussien qui nous intéresse, cela revient à pré-calculer les trois images \textit{cumulées} $S_1$, $S_x$ et $S_{x^2}$ définies par
+\item Les valeurs $T_g(i,j)$ peuvent être calculées préalablement à la phase de segmentation proprement dite. Pour le cas gaussien qui nous intéresse, cela revient à pré-calculer les trois images \textit{cumulées} $S_1$, $S_I$ et $S_{I^2}$ définies par
\begin{alignat}{4}
\label{eq-img-cumul}
- S_1(i,j) &= \sum_{x=0}^jx & \quad \text{,}\quad S_x(i,j) &= \sum_{x=0}^jv(i,x) & \quad \text{et}&\quad & S_{x^2}(i,j) &= \sum_{x=0}^jv(i,x)^2
+ S_1(i,j) &= \sum_{x=0}^jx & \quad \text{,}\quad S_I(i,j) &= \sum_{x=0}^jv(i,x) & \quad \text{et}&\quad & S_{I^2}(i,j) &= \sum_{x=0}^jv(i,x)^2
\end{alignat}
\item Les valeurs du coefficient $C(i,j)$ se calculent très facilement durant la génération du contour $\Gamma$.
\end{itemize}
-Par ailleurs, le choix d'un contour polygonal permet également d'améliorer l'efficacité de l'algorithme car lors de la phase de segmentation, le déplacement d'un sommet du polygone n'influe que sur les pixels des deux segments qui s'y rapportent, réduisant ainsi la quantité de calculs à effectuer à chaque nouvelle déformation du contour.
+Par ailleurs, le choix d'un contour polygonal permet d'améliorer l'efficacité de l'algorithme. Lors de la phase de segmentation, le déplacement d'un sommet du polygone n'influe ainsi que sur les pixels des deux segments qui s'y rapportent, réduisant ainsi la quantité de calculs à effectuer à chaque nouvelle déformation du contour.
\begin{upminfo}
- L'approche décrite dans ce chapitre n'est valide que si les segments formant le polygone du contour ne se croisent pas. Il est donc nécessaire, lors de la convergence de la segmentation, d'empêcher les croisements de segments. Une solution simple a été proposée dans \cite{ChesnaudRB99} et nous l'avons parallélisée dans le cadre des travaux présentés ici.
+L'approche décrite dans ce chapitre n'est valide que si les segments formant le polygone du contour ne se croisent pas. Il est donc nécessaire, lors de la convergence de la segmentation, d'empêcher ces croisements. Une solution simple a été proposée dans \cite{ChesnaudRB99}, que nous avons parallélisée et intégrée.
\end{upminfo}
-\subsection{Implémentation séquentielle}\label{snake-cpu-impl}
-Un des inconvénients des algorithmes de type \textit{snake} est l'influence du contour initial sur la convergence de la segmentation. Pour pallier simplement ce défaut, une technique progressive est adoptée, en initialisant le contour avec peu de sommets (4) puis en augmentant le nombre au fur et à mesure de la convergence. L'algorithme \ref{algo-snake-cpu1} décrit macroscopiquement la solution mise \oe uvre tandis que l'algorithme \ref{algo-snake-cpu2} en présente les détails.
+\subsection{Implémentation séquentielle\label{snake-cpu-impl}}
+Un des inconvénients des algorithmes de type \textit{snake} est l'influence du contour initial sur la convergence de la segmentation. Pour pallier simplement ce défaut, une technique progressive est adoptée, en initialisant le contour avec peu de sommets (par exemple 4) puis en augmentant le nombre au fur et à mesure de la convergence. L'algorithme \ref{algo-snake-cpu1} décrit macroscopiquement la solution mise \oe uvre tandis que l'algorithme \ref{algo-snake-cpu2} en présente les détails.
-\begin{algorithm}
-\label{algo-snake-cpu}
+\begin{algorithm}[h]
+\label{algo-snake-cpu1}
+\small
\caption{Principe mis en \oe uvre pour la convergence du snake polygonal}
Calculer les images cumulées\;
Initialiser le contour avec 4 sommets\;
}
\end{algorithm}
-\begin{algorithm}[h]
+\begin{algorithm}[H]
\caption{Détail de l'implémentation du snake polygonal}
-\label{cpualgo}
+\small
+\label{algo-snake-cpu2}
Lire l'image $\bar{v}$\;
- Calculer les images cumulées $S_1$, $S_x$ $S_{x^2}$ \nllabel{li-img-cumul}\tcc*[r]{en parallèle via SSE2}
+ Calculer les images cumulées $S_1$, $S_I$ $S_{I^2}$ \nllabel{li-img-cumul}\tcc*[r]{}
$n \leftarrow 0$ \tcc*[r]{indice de boucle niveau contour}
$N_n \leftarrow 4$ \tcc*[r]{nombre de n\oe uds}
$\Gamma \leftarrow \{\Gamma_0,\Gamma_1,\Gamma_2,\Gamma_3\} $\;
$d \leftarrow d_{max}$ \tcc*[r]{pas de déplacement des n\oe uds}
$l_{min} = 32$ \tcc*[r]{longueur mini des segments sécables}
$\Gamma_i \leftarrow \Gamma_0$ \tcc*[r]{sommet courant}
- $GL_{ref} \leftarrow GL(\Gamma, N_n, \bar{v}, S_y, S_{y^2})$ \tcc*[r]{voir à partir de ligne 18 pour le détail}
+ $GL_{ref} \leftarrow GL(\Gamma, N_n, \bar{v}, S_I, S_{I^2})$ \tcc*[r]{voir à partir de ligne 18 pour le détail}
\Repeter(\tcc*[f]{niveau contour}){$N_{add}=0$}{
$N_{add}\leftarrow 0$\;
\Repeter(\tcc*[f]{niveau n\oe ud}){$N_{move}=0$}{
$N_{move}\leftarrow 0$\;
\Pour{$i=0$ à $i=N_n-1$}{
- Calculer les positions $\{\Gamma_i^0, \dots, \Gamma_i^7\}$ \tcc*[r]{les 8 voisins de $\Gamma_i$ }
- \Pour{$w=0$ à $w=7$}{
- Soustraire à $GL_{ref}$ la contribution des segments $\Gamma_{i-1}\Gamma_i$ et $\Gamma_{i}\Gamma_{i+1}$\;
+ Calculer les positions $\{\Gamma_i^0, \dots, \Gamma_i^7\}$ \tcc*[r]{les 8 voisins de $\Gamma_i$ à distance $d$ }
+ $GL_w \leftarrow GL_{ref} - $ la contribution des segments $\Gamma_{i-1}\Gamma_i$ et $\Gamma_{i}\Gamma_{i+1}$\;
+ \Pour{$w=0$ à $w=7$}{
Discrétiser les segments $\Gamma_{i-1}\Gamma_i^w$ et $\Gamma_{i}^w\Gamma_{i+1}$\nllabel{li-bresen}\;
- Lire dans $S_1$, $S_x$ et $S_{x^2}$ les contributions des pixels de $\Gamma_{i-1}\Gamma_i^w$ et $\Gamma_{i}^w\Gamma_{i+1}$\nllabel{li-contrib-seg-deb}\;
+ Lire dans $S_1$, $S_I$ et $S_{I^2}$ les contributions des pixels de $\Gamma_{i-1}\Gamma_i^w$ et $\Gamma_{i}^w\Gamma_{i+1}$\nllabel{li-contrib-seg-deb}\;
Calculer les directions et lire les codes de Freeman \;
Calculer $GL_w$ incluant les contributions de $\Gamma_{i-1}\Gamma_i^w$ et $\Gamma_{i}^w\Gamma_{i+1}$ \nllabel{li-contrib-seg-fin}\;
- \lSi{$GL_w > GL_{ref}$}{
+ \Si{$GL_w > GL_{ref}$}{
$GL_{ref} \leftarrow GL_w$\;
$\Gamma_i \leftarrow \Gamma_i^w$\;
$N_{move} \leftarrow N_{move}+1$\;
}
$N_n \leftarrow N_n + N_{add}$\;
\lSi{$d > 1$}{ $d \leftarrow d/2$ } \lSinon{ $d \leftarrow 1$ }\;
- $GL_{ref} \leftarrow GL(\Gamma, N_n, \bar{v}, S_y, S_{y^2})$ \;
+ $GL_{ref} \leftarrow GL(\Gamma, N_n, \bar{v}, S_I, S_{I^2})$ \;
}
\end{algorithm}
-Les différentes sommations nécessaires au calcul de la valeur du critère $GL$ sont effectuées en parallèle à l'aide du jeu d'instructions SSE2. La taille des registres utilisables est de 128 bits et permet ainsi de traiter des images de 4096$\times$4096 pixels dont les niveaux de gris sont codés sur 16 bits. Cela ne laisse toutefois que 12 bits pour le codage des segments du contour et limite ainsi leur longueur à 4096 pixels. L'organisation d'un registre SSE 128 bits est donc la suivante :
-\begin{itemize}
-\item 24 bits pour les sommes à opérandes dans $S_1$
-\item 24+16 = 40 bits pour les sommes à opérandes dans $S_x$
-\item 24+32 = 60 bits pour les sommes à opérandes dans $S_x^2$
-\end{itemize}
-Soit un total de 124 bits, qui peuvent donc être représentés par un registre SSE2.
+% Les différentes sommations nécessaires au calcul de la valeur du critère $GL$ sont effectuées en parallèle à l'aide du jeu d'instructions SSE2, qui permet de travailler avec des registres de grande capacité (128 bits) et d'envisager d'y ranger côte à côte les opérandes des trois sommes pour les effectuer simultanément.
+% Si l'on cherche à traiter des images en niveaux de gris sont codés sur 16 bits, les sommes $S_1$, $S_X$ et $S_{X^2}$ vont utiliser :
+% \begin{itemize}
+% \item $N_c$ bits par opérande de chaque somme pour représenter les coordonnées des pixels.
+% \item $N_p$ bits pour traduire le nombre d'opérandes dans chaque somme.
+% \item 16 bits par valeur de niveau de gris dans $S_X$.
+% \item 32 bits par valeur de niveau de gris au carré dans $S_{X^2}$.
+% \end{itemize}
+% Les trois sommes utilisent donc, par opérande, un total de $\left(3\left(N_c+N_p\right)+16+32\right)$ bits devant être contenu dans un registre de 128 bits, ce qui nous donne un maximum de 26 bits pour $N_c+N_p$.
+% La longueur des segments pouvant être au maximum $\sqrt{2}$ fois supérieure au coté de l'image, on peut donc considérer qu'il est nécessaire d'avoir $N_p = N_c+1$ pour ne pas générer de restriction sur la longueur des segments. Cela nous conduit donc à $N_c = 12$ et $N_p=13$ ($12+13 = 25 < 26$).
+% La répartition retenue pour les données dans les registres SSE2 de 128 bits est alors la suivante :
+% \begin{itemize}
+% \item $N_c+N_p=25$ bits pour les opérandes des sommes de $S_1$.
+% \item $N_c+N_p+16=41$ bits pour les opérandes des sommes de $S_X$.
+% \item $N_c+N_p+32=57$ bits pour les opérandes des sommes de $S_{X^2}$.
+% \end{itemize}
\subsection{Performances}
-Les images de 1024$^2$ pixels de la figure \ref{fig-snakecpu-cochon512} montrent l'évolution du contour lors de la segmentation d'une image photographique prise en faible éclairement et bruitée artificiellement par un bruit gaussien d'écart type 25. Les paramètres de la séquence sont fixés empiriquement aux valeurs $d_{max}=16, l_{min}=8$.
-Les temps d'exécution indiqués sont mesurés sur Intel Xeon E5530-2.4GHz with 12Go RAM et sont les valeurs moyennes obtenues pour 10 exécutions.
+Les images de 1024$\times$1024 pixels de la figure \ref{fig-snakecpu-cochon512} montrent l'évolution du contour lors de la segmentation d'une image photographique prise en faible éclairement et bruitée artificiellement par un bruit gaussien d'écart type 25. Les paramètres de la séquence sont fixés empiriquement aux valeurs $d_{max}=16, l_{min}=8$.
+Les temps d'exécution indiqués sont mesurés sur Intel Xeon E5530-2.4GHz avec 12Go RAM.
\begin{figure}
\centering
\subfigure[Initialisation : 4 n\oe uds]{\includegraphics[height=3.5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/im000.png}}
\subfigure[Itération 10 : 244 n\oe uds 3~ms]{\includegraphics[height=3.5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/im010.png}}
\subfigure[Itération 13 : 256 n\oe uds 3~ms]{\includegraphics[height=3.5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/im013.png}}
\subfigure[Itération 14 : 256 n\oe uds 3~ms]{\includegraphics[height=3.5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/im014.png}}
- \caption{Évolution du contour lors de la segmentation d'une image de 512$^2$ pixels. La convergence est obtenue à l'itération 14 après 44~ms pour un total de 256 n\oe uds.}
+ \caption{Évolution du contour lors de la segmentation d'une image de 512$\times$512 pixels. La convergence est obtenue à l'itération 14 après 44~ms pour un total de 256 n\oe uds.}
\label{fig-snakecpu-cochon512}
\end{figure}
-La dépendance vis à vis du contour initial qui est un des principaux soucis liés au snake est ici fortement relativisée. La figure \ref{fig-snakecpu-compinit} montre le contour final segmentant l'image de test de la figure \ref{fig-snakecpu-cochon512} à partir d'un état initial très éloigné du précédent et \textit{a priori} très défavorable compte tenu du fait qu'il est loin de la cible et sans intersection avec elle. Toutefois, le contour final obtenu est très proche de celui obtenu à partir d'un état initial englobant la cible, malgré un n\oe ud qui s'est ``accroché'' au bord de l'image. La convergence est également plus longue à obtenir dans ce cas avec Un total de 17 itérations en 87~ms et 273 n\oe uds.
+La dépendance vis à vis du contour initial, qui est un des principaux soucis liés au snake est ici fortement relativisée. La figure \ref{fig-snakecpu-compinit} montre le contour final segmentant l'image de test de la figure \ref{fig-snakecpu-cochon512} à partir d'un état initial très éloigné du précédent et, \textit{a priori}, très défavorable compte tenu du fait qu'il est loin de la cible et sans intersection avec elle. Toutefois, le contour final est très proche de celui obtenu à partir d'un état initial englobant la cible, avec la particularité d'avoir identifié la peluche au fond et la zone sombre à la cible avec un n\oe ud à chaque coin de l'image. Cela est du à l'initialisation des zones avec comme cible une portion de la zone sombre, relativement homogène. La convergence est également plus longue à obtenir dans ce cas avec 87~ms pour 17 itérations et 273 n\oe uds.
\begin{figure}
\centering
\label{fig-snakecpu-compinit}
\end{figure}
-La dimension de l'image à traiter a également un effet sur le résultat et naturellement sur le temps de calcul. Si l'on conserve les mêmes paramètres d'optimisation que pour la segmentation de l'image 512$^2$ pixels et un contour initial dont les cotés sont à une distance des bords équivalente à 10\% des cotés de l'image, le résultat sur une image identique de 4000$^2$ pixels est obtenu en 1.3~s avec 1246 n\oe uds ; il est reproduit à la figure \ref{fig-snakecpu-cochon4ka}. Le nombre de pixels appartenant à la région cible est tel que l'amplitude des déplacements autorisés pour chaque n\oe ud ($d$) peut se révéler trop faible vis à vis du seuil d'acceptation des mouvements. On observe que les zones à gradient élevé ne posent pas de problème et sont détourées de la même manière, alors que dans le bas de l'image où figure une zone de gradient faible (ombre), la cible se trouve maintenant quelque peu surévaluée en surface là ou elle était plutôt sous évaluée dans l'image en 512$^2$ pixels.
-On parvient à un résultat très proche beaucoup plus rapidement en adaptant les paramètres à la taille de l'image, comme le montre par exemple la segmentation de la figure \ref{fig-snakecpu-cochon4kb}, effectuée avec $d_{max}=128$ et $l_{min}=32$ et qui converge vers un contour de 447 n\oe uds en moins de 0.7~s.
-Au delà des 16 millions de pixels (4000$^2$ pixels), l'implémentation séquentielle est toujours possible mais doit se priver des instructions SSE. Nous avons, avec leur accord, adapté le code des auteurs en ce sens et réalisé les mesures pour des tailles allant jusqu'à 150~MP. La table \ref{tab-snakecpu-speed-size} en synthétise les résultats en distinguant chaque fois le temps pris par les pré-calculs et celui nécessaire à la convergence de la segmentation.
+La dimension de l'image à traiter a également un effet sur le résultat et naturellement sur le temps de calcul. Si l'on conserve la même stratégie d'optimisation que pour la segmentation de l'image 512$\times$512 pixels et un contour initial dont les cotés sont à une distance des bords équivalente à 10\% des cotés de l'image, le résultat sur une image identique de 4000$^2$4000 pixels est obtenu en 1,3~s avec 1246 n\oe uds ; il est reproduit à la figure \ref{fig-snakecpu-cochon4ka}. Le nombre de pixels appartenant à la région cible est tel que l'amplitude des déplacements autorisés pour chaque n\oe ud ($d$) peut se révéler trop faible vis-à-vis du seuil d'acceptation des mouvements. On observe que les zones à fort contraste ne posent pas de problème et sont détourées de la même manière, alors que dans le bas de l'image où figure une zone de faible contraste (ombre), la cible se trouve maintenant quelque peu surévaluée en surface là ou elle était plutôt sous évaluée dans l'image en 512$\times$512 pixels.
+Ces deux contours correspondent chacun à un minimum local vers lequel l'algorithme du snake a convergé, mais les variances associées demeurent extrêmement proches.
+On parvient à un résultat très proche beaucoup plus rapidement en adaptant les paramètres à la taille de l'image, comme le montre par exemple la segmentation de la figure \ref{fig-snakecpu-cochon4kb}, effectuée avec $d_{max}=128$ et $l_{min}=32$ et qui converge vers un contour de 447 n\oe uds en moins de 0,7~s.
+Au delà des 16 millions de pixels (4000$\times$4000 pixels), l'implémentation séquentielle de référence ne permettait plus le traitement. Nous avons, avec l'accord des auteurs, adapté leur code en ce sens et réalisé les mesures pour des tailles allant jusqu'à 150~MP. La table \ref{tab-snakecpu-speed-size} en synthétise les résultats en distinguant chaque fois le temps pris par les pré-calculs et celui nécessaire à la convergence de la segmentation. On constate que les deux étapes et donc le temps total varient linéairement avec la taille de l'image.
\begin{figure}
\centering
- \subfigure[$d_{max}=16$ et $l_{min}=8$, 1246 n\oe uds en 1.3~s]{\label{fig-snakecpu-cochon4ka}\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/snakecpu-cochon4k.png}}
- \subfigure[$d_{max}=128$ et $l_{min}=32$, 447 n\oe uds en 0.7~s]{\label{fig-snakecpu-cochon4kb}\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/snakecpu-cochon4k-128-32.png}}
- \caption{Segmentation de l'image de test en 4000$^2$ pixels.}
+ \subfigure[$d_{max}=16$ et $l_{min}=8$, 1246 n\oe uds en 1.3~s]{\label{fig-snakecpu-cochon4ka}\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/snakecpu-cochon4k-t.png}}
+ \subfigure[$d_{max}=128$ et $l_{min}=32$, 447 n\oe uds en 0.7~s]{\label{fig-snakecpu-cochon4kb}\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/snakecpu-cochon4k-128-32-t.png}}
+ \caption{Segmentation de l'image de test en 4000$\times$4000 pixels. Le tracé du contour a été artificiellement épaissi pour le rendre visible à l'échelle de l'impression.}
\label{fig-snakecpu-cochon4k}
\end{figure}
+\begin{figure}[h]
+ \centering
+ \includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cochon4kc-128-8-t.png}
+ \caption{Segmentation de l'image de test en 4000$\times$4000 pixels avec une cible de petite taille. Le contour initial est la transcription de celui utilisé à la figure \ref{fig-snakecpu-cochon512}. Le tracé du contour a été artificiellement épaissi pour le rendre visible à l'échelle de l'impression.}
+ \label{fig-snakecpu-cochon4kc3}
+\end{figure}
+
+
\begin{table}[h]
\centering
{\bf Total} &0,51&4,08&5,7\\
\bottomrule
\end{tabular}
- \caption{Performances (en secondes) de la segmentation par snake polygonal sur CPU en fonction de la taille de l'image à traiter. Le temps sont obtenus avec la même image de test dilatée et bruitée et un contour initial carré dont la distance aux bords est proportionnelle à la taille de l'image. Seule l'image en 15~MP a pu être traitée par une implémentation utilisant SSE2.}
+ \caption{Performances (en secondes) de la segmentation par snake polygonal sur CPU en fonction de la taille de l'image à traiter. Les temps sont obtenus avec la même image de test dilatée et bruitée et un contour initial carré dont la distance aux bords est proportionnelle à la taille de l'image.}
\label{tab-snakecpu-speed-size}
\end{table}
-Enfin, il faut aussi considérer les tailles relatives de la cible et de l'image. Ainsi, si on fait l'hypothèse d'une cible de petite taille ``noyée'' dans une image de grandes dimensions, les résultats de la segmentation seront impactés en raison, cette fois, d'une moindre adaptation à la cible lors des toutes premières itérations, les plus grossières, où le nombre de n\oe uds et réduit et le pas de déplacement potentiellement grand vis à vis de la cible. Ce cas de figure est illustré par la segmentation reproduite à la figure \ref{fig-snakecpu-cochon4kc3} et qui met en évidence une qualité moindre par la confusion des zones les plus sombres de la cible avec le fond.
+Enfin, il faut aussi considérer les tailles relatives de la cible et de l'image. Ainsi, si on fait l'hypothèse d'une cible de petite taille \og noyée \fg{} dans une image de grandes dimensions, les résultats de la segmentation seront impactés en raison, cette fois, d'une moindre adaptation à la cible lors des toutes premières itérations, les plus grossières, où le nombre de n\oe uds est réduit et le pas de déplacement potentiellement grand vis à vis de la cible. Ce cas de figure est illustré par la segmentation reproduite à la figure \ref{fig-snakecpu-cochon4kc3} et qui met en évidence une qualité moindre par la confusion des zones les plus sombres de la cible avec le fond et confirme ainsi la nécessité d'adapter la stratégie d'optimisation au problème posé.
-\begin{figure}[h]
- \centering
- \includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cochon4kc-128-8.png}
- \caption{Segmentation de l'image de test en 4000$^2$ pixels avec une cible de petite taille. Le contour initial est celui utilisé à la figure \ref{fig-snakecpu-cochon4k}.}
- \label{fig-snakecpu-cochon4kc3}
-\end{figure}
-
\section{Implémentation parallèle GPU du snake polygonal}
-L'analyse de l'exécution du programme séquentiel révèle la prépondérance des blocs fonctionnels suivants, dans l'ordre d'importance, qui occupent à eux seuls plus de 80\% du temps total d'exécution :
+L'analyse de l'exécution du programme séquentiel révèle la prépondérance des blocs fonctionnels suivants, qui occupent à eux seuls 80\% du temps total d'exécution :
\begin{itemize}
- \item Le calcul de la contribution des segments (lignes \ref{li-contrib-seg-deb} à \ref{li-contrib-seg-fin} dans l'algorithme \ref{cpualgo})
- \item La génération des trois images cumulées, avant le début des itérations (ligne \ref{li-img-cumul}).
- \item La discrétisation des segments définis par les coordonnées de leurs extrémités (ligne \ref{li-bresen}).
+ \item Le calcul de la contribution des segments, pour environ 50\% (lignes \ref{li-contrib-seg-deb} à \ref{li-contrib-seg-fin} dans l'algorithme \ref{algo-snake-cpu2})
+ \item La génération des trois images cumulées, avant le début des itérations, pour environ 20\% (ligne \ref{li-img-cumul}).
+ \item La discrétisation des segments définis par les coordonnées de leurs extrémités, pour environ 7\% (ligne \ref{li-bresen}).
\end{itemize}
Cette proportion est globalement conservée lorsque la taille de l'image à traiter varie, comme le montre le graphique de la figure \ref{fig-snakecpu-chronos1}
\caption{Évolution du coût relatif des trois fonctions les plus consommatrices en temps de calcul en fonction de la taille de l'image à traiter.}
\end{figure}
-Si l'effort de parallélisation porte essentiellement sur ces fonctions coûteuses, l'ensemble du traitement est réalisé sur le GPU afin de réduire autant que possible les transferts entre le GPU et le système hôte qui, selon le volume concerné, sont susceptibles de grever considérablement la performance globale. L'hôte ne conserve que l'initiative du transfert initial et le contrôle de la boucle principale, ne nécessitant l'échange que d'un seul octet à chaque itération (le nombre de nouveau n\oe uds $N_{add}$).
+Si l'effort de parallélisation porte essentiellement sur ces fonctions coûteuses, l'ensemble du traitement est réalisé sur le GPU afin de réduire autant que possible les transferts entre le GPU et le système hôte qui, selon le volume concerné, sont susceptibles de grever considérablement la performance globale. L'hôte ne conserve que l'initiative du transfert initial et le contrôle de la boucle principale, ne nécessitant l'échange que d'un seul octet à chaque itération (représentant le nombre de nouveaux n\oe uds $N_{add}$).
-Les traitements étant totalement indépendants, nos traitons séparément la parallélisation des pré-calculs et celle de la segmentation.
+Les traitements étant totalement indépendants, nous traitons séparément la parallélisation des pré-calculs et celle de la segmentation.
-\subsection{pré-calculs des images cumulées}
-Pour réduire la quantité de mémoire requise, nous avons choisi de ne pas générer l'image $S_1$ mais plutôt d'en calculer les valeurs à la volée. L'expression en est simple et le temps pris par les opération élémentaires qu'elle met en jeu est largement compensé par le gain obtenu en économisant les accès mémoire qui auraient été nécessaires, ce qui n'est pas le cas des deux autres images $S_x$ et $S_x^2$ dont le calcul est quant à lui réalisé en appliquant la méthode des \textit{prefixsums} décrite dans \cite{BlellochTR90}.
+\subsection{Pré-calculs des images cumulées}
+Pour réduire la quantité de mémoire requise, nous avons choisi de ne pas générer l'image $S_1$ mais plutôt d'en calculer les valeurs à la volée. L'expression en est simple et le temps pris par les opérations élémentaires qu'elle met en jeu est largement compensé par le gain obtenu en économisant les accès mémoire qui auraient été nécessaires, ce qui n'est pas le cas des deux autres images $S_I$ et $S_I^2$ dont le calcul est quant à lui réalisé en appliquant une variante de la méthode des sommes préfixées (\textit{prefixsums}) décrite dans \cite{BlellochTR90} et qui permet d'évaluer les expressions de l'équation \eqref{eq-img-cumul}.
-Les sommations se font au niveau de chaque ligne de l'image, que l'on décompose en $n$ blocs de $bs$ pixels où $bs$ correspond aussi au nombre de threads exécutés par chaque bloc de la grille de calcul. La valeur $bs$ étant obligatoirement une puissance de 2 supérieure à 32, le bloc de pixels d'indice $n-1$ doit éventuellment être complété par des valeurs nulles. Chaque bloc de thread réalise son traitement indépendemment des autres, mais l'ensemble des sommes de bloc étant requise pour le calcul des sommes globales, une synchronisation est nécessaire à deux endroits du calcul. Nous avons choisi d'assurer ces synchronisations en découpant le traitement en trois \textit{kernels} distincts, rendant par la même occasion le code plus concis :
+Les sommations se font au niveau de chaque ligne de l'image, que l'on décompose en $n$ blocs de $bs$ pixels où $bs$ correspond aussi au nombre de threads exécutés par chaque bloc de la grille de calcul. La valeur $bs$ étant obligatoirement une puissance de 2 supérieure à 32, le bloc de pixels d'indice $n-1$ doit éventuellement être complété par des valeurs nulles. Chaque bloc de thread réalise son traitement indépendemment des autres, mais l'ensemble des sommes de bloc étant requise pour le calcul des sommes globales, une synchronisation est nécessaire à deux endroits du calcul. Nous avons choisi d'assurer ces synchronisations en découpant le traitement en trois \textit{kernels} distincts, rendant par la même occasion le code plus concis. Les diagrammes de la figure \ref{fig-calcul-cumuls} donnent le détail des opérations effectuées par ces trois \textit{kernels} pour l'image cumulée $S_I$ :
\begin{itemize}
-\item \texttt{compute\_block\_prefixes()} est le \textit{kernel} effectuant, en mémoire partagée, les prefixsums inclusif de chaque bloc, puis mémorisant les sommes des blocs, \textit{i.e} le dernier élément du prefixsum associé, dans deux vecteurs $V_x$ et $V_x^2$ en mémoire globale. L'ensemble des prefixsums est également mémorisé en mémoire globale. La largeur de l'image n'étant pas nécessairement une puissance de 2, il est nécessaire de faire du remplissage avec des valeurs nulles dans le dernier bloc (indice $n-1$).
+\item \texttt{compute\_block\_prefixes()} est le \textit{kernel} effectuant, en mémoire partagée, la \textit{prefixsum} inclusive de chaque bloc, puis qui en mémorise les sommes, \textit{i.e} le dernier élément, dans deux vecteurs $V_x$ et $V_x^2$ en mémoire globale. L'ensemble des prefixsums est également mémorisé en mémoire globale. La largeur de l'image n'étant pas nécessairement une puissance de 2, il est nécessaire de faire du remplissage avec des valeurs nulles dans le dernier bloc (indice $n-1$).
\item \texttt{scan\_blocksums()} est le \textit{kernel} effectuant les prefixsum exclusifs des vecteurs $V_x$ et $V_x^2$. Les résultat demeurent respectivement dans $V_x$ et $V_x^2$.
-\item \texttt{add\_sums2prefixes()} est le \textit{kernel} effectuant les additions de chaque élément d'indice $i$ des vecteurs $V_x$ (respectivement $V_x^2$ avec tous les éléments du prefixsum du bloc de même indice $i$ associé.
+\item \texttt{add\_sums2prefixes()} est le \textit{kernel} effectuant les additions de chaque élément d'indice $i$ des vecteurs $V_x$ (respectivement $V_x^2$) avec tous les éléments du prefixsum du bloc de même indice $i$.
\end{itemize}
-Les diagrammes de la figure \ref{fig-calcul-cumuls} donnent le détail des opérations effectuées par ces trois \textit{kernels} pour l'image cumulée $S_x$. La seconde image cumulée $S_x^2$ est obtenues exactement de la même manière en sommant non plus les valeurs $v_k$ mais $v^2_k$.
+La seconde image cumulée $S_{I^2}$ est obtenues exactement de la même manière en sommant non plus les valeurs $v_k$ mais $v^2_k$.
\begin{figure}
\centering
\subfigure[Détail des opérations effectuées par le \textit{kernel} \texttt{compute\_block\_prefixes()}. La valeur $bs$ correspond au nombre de pixels de chaque bloc, qui est aussi le nombre de threads exécuté par chaque bloc de la grille de calcul.]{\resizebox{0.9\linewidth}{!}{ \input{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/GPUcumuls.pdf_t}}}\vspace{1cm}
\subfigure[Détail des opérations effectuées par le \textit{kernel} \texttt{scan\_blocksums()}.]{\resizebox{0.9\linewidth}{!}{ \input{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/GPUscansomblocs.pdf_t}}}\vspace{1cm}
\subfigure[Détail des opérations effectuées par le \textit{kernel} \texttt{add\_sums2prefixes()}.]{\resizebox{0.9\linewidth}{!}{ \input{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/GPUaddsoms2cumuls.pdf_t}}}
+ \caption{Calcul des images cumulées $S_I$ et $S_{I^2}$ en trois étapes successives. a) cumul partiel bloc par bloc et mémorisation de la somme de chaque bloc. b) cumul sur le vecteur des sommes partielles. c) ajout des sommes partielles à chaque élément des blocs cumulés.}
\label{fig-calcul-cumuls}
- \caption{Calcul des images cumulées $S_x$ et $S_x^2$ en trois étapes successives. a) cumul partiel bloc par bloc et mémorisation de la somme de chaque bloc. b) cumul sur le vecteur des sommes partielles. c) ajout des sommes partielles à chaque élément des blocs cumulés.}
\end{figure}
-Les gains de performance de cette implémentation GPU comparée à l'implémentation CPU/SSE2 sont ceux de la table \ref{tab-speedup-cumuls}, soit un GPU environ 7 fois plus rapide pour des images de taille 15 à 150 millions de pixels. L'influence de la taille d'image sur le gain est faible, mais on peut toutefois noter que plus l'image est grande plus le gain est important.
+Les gains de performance de cette implémentation GPU comparée à l'implémentation CPU (mono thread) sont ceux de la table \ref{tab-speedup-cumuls}, soit un GPU environ 7 fois plus rapide pour des images de taille 15 à 150 millions de pixels. L'influence de la taille d'image sur le gain est faible, mais on peut toutefois noter que plus l'image est grande plus le gain est important.
+Les accélérations constatées peuvent sembler faibles en regard de ce que l'on attend d'un GPU, mais il faut rappeler que ce type d'opération (les réductions) n'est pas véritablement adapté à leur architecture en raison d'une grande inter-dépendance des données d'une étape de calcul à l'autre. Sans une implémentation optimisée, cette opération s'exécuterait même plus lentement sur GPU que sur un CPU.
On obtient des accélérations supérieures en rendant le calcul moins générique et en développant des versions spécifiques des trois \textit{kernels}, dédiées par exemple au traitement des images dont largeur est multiple de 256 pixels.
\begin{table}[h]
\centering
{\bf Accélération} &{\bf 6,5} &{\bf 6,9} &{\bf 7,0}\\
\bottomrule
\end{tabular}
- \caption{Accélération constatée, pour le calcul des images cumulées, de l'implémentation GPU par rapport à l'implémentation CPU de référence.}
+ \caption{Accélération constatée, pour le calcul des images cumulées, de l'implémentation GPU (C2070) par rapport à l'implémentation CPU de référence.}
\label{tab-speedup-cumuls}
\end{table}
-\subsection{Calcul des contribution des segments}
+\subsection{Calcul des contributions des segments}
-Le déplacement d'un des $N_n$ n\oe uds du contour $\Gamma$ vers l'une des 8 positions voisines permises, impose d'évaluer les contributions des 8 paires de segments associées, soit $16N_n$ segments pour la totalité du contour, que nous évaluons en parallèle au sein du \textit{kernel} \texttt{GPU\_compute\_segments\_contribs()}. Pour ce faire, chaque segment doit tout d'abord être discrêtisé en une suite de pixels puis, en conservant la règle \textit{1 pixel par thread} la contribution de chaque pixel est déterminée avant de toutes les additionner pour obtenir la contribution du segment.
+Le déplacement d'un des $N_n$ n\oe uds du contour $\Gamma$ vers l'une des 8 positions voisines à distance $d$\footnote{Sous réserve que la position considérée ne dépasse pas les limites de l'image.}, impose d'évaluer les contributions des 8 paires de segments associées, soit $16N_n$ segments pour la totalité du contour, que nous évaluons en parallèle au sein du \textit{kernel} \texttt{GPU\_compute\_segments\_contribs()}. Pour ce faire, chaque segment doit tout d'abord être discrétisé en une suite de pixels puis, en conservant la règle \textit{1 pixel par thread} la contribution de chaque pixel est déterminée avant de toutes les additionner pour obtenir la contribution du segment.
Les pixels représentant les n\oe uds font l'objet d'un traitement spécifique impliquant les codes de Freeman, pour ne pas fausser les contributions globales (voir paragraphe \ref{snake-cpu-impl}).
-Pour optimiser l'exécution de ce kernel et réduire l'effet de la disparité des longueurs des segments, nous créons un motif régulier en mémoire, basé sur la longueur $npix_{max}$ du plus grand segment et moyennant le recours au remplissage par des valeurs neutres là où c'est nécessaire pour les autres segments.
+Pour optimiser l'exécution de ce kernel et réduire l'effet de la disparité des longueurs des segments, nous avons crée un motif régulier en mémoire, basé sur la longueur $npix_{max}$ du plus grand segment et avons complété les blocs associés aux segments de longueur inférieure à $npix_{max}$ avec des valeurs neutres pour l'opération réalisée, c'est-à-dire des valeurs nulles.
Si $bs_{max}$ est la taille de bloc maximale admissible par le GPU, la taille $bs$ des blocs de threads/pixels employée pour le calcul des contributions des segments est alors déterminée de la façon suivante :
\[
\end{cases}
\]
-Dans notre implémentation, les calculs sont faits en mémoire partagée et la quantité nécessaire limite la taille de bloc admissible. Nous limitons celle-ci à 256 sur C1060 et 512 sur C2050. Toutefois, les tests ont montré que sur ces deux versions de l'architecture, La taille maximale conduisant aux meilleures performances est de 256 threads par bloc.
+Dans notre implémentation, les calculs sont faits en mémoire partagée et la quantité nécessaire limite la taille de bloc admissible. Nous limitons celle-ci à 256 sur C1060 et 512 sur C2050. Toutefois, les tests ont montré que sur ces deux versions de l'architecture, la taille maximale conduisant aux meilleures performances est de 256 threads par bloc.
-Le \textit{kernel} \texttt{GPU\_compute\_segments\_contribs()} calcule alors en paralèle pour tous les segments les coordonnées de tous les pixels qui les composent. Nous mettons pour cela en \oe uvre l'algorithme de Bresenham, \textit{i.e} la méthode du segment semi-ouvert, en distinguant les cas où
+Le \textit{kernel} \texttt{GPU\_compute\_segments\_contribs()} calcule alors en parallèle pour tous les segments les coordonnées de tous les pixels qui les composent. Nous mettons pour cela en \oe uvre l'algorithme de Bresenham, \textit{i.e} la méthode du segment semi-ouvert, en distinguant les cas où
\begin{itemize}
-\item la valeur absolue de la pente $k$ du segment à discrêtiser est supérieure à $1$; on applique alors la méthode au segment \textit{horizontal} semi-ouvert et on obtient un pixel par ligne.
-\item la valeur absolue de la pente $k$ du segment à discrêtiser est inférieure ou égale à $1$; on applique alors la méthode au segment \textit{vertical} semi-ouvert et on obtient un pixel par colonne.
+\item la valeur absolue de la pente $k$ du segment à discrétiser est supérieure à $1$; on applique alors la méthode au segment \textit{horizontal} semi-ouvert et on obtient un pixel par ligne.
+\item la valeur absolue de la pente $k$ du segment à discrétiser est inférieure ou égale à $1$; on applique alors la méthode au segment \textit{vertical} semi-ouvert et on obtient un pixel par colonne.
\end{itemize}
Cette distinction nous permet de conserver la règle \textit{1 pixel par thread} importante pour la régularité des motifs d'accès en mémoire et aussi pour \textit{charger} au maximum le GPU.
-La figure \ref{fig-structure-segment} représente la structure décrite ci-dessus pour la représentation en mémoire des segments à évaluer. La première ligne montre le détail du premier segment, avec la correspondance \textit{1 pixel par thread} et le découpage en un nombre de blocs suffisant pour permettre de décrire le plus long des segments.
-
-La seconde ligne présente l'ordre dans lequel sont concaténés les 16 groupes de blocs-segment associés au déplacement d'un n\oe ud particulier.
-
\begin{figure}
\centering
\resizebox{0.8\linewidth}{8cm}{ \input{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/contribs_segments.pdf_t}}
-\label{fig-structure-segment}
\caption{Structuration des données en mémoire du GPU pour l'évaluation en parallèle de l'ensemble des évolutions possibles du contour.}
+\label{fig-structure-segment}
\end{figure}
-Aux deux dernières lignes est décrite la concaténation des ensembles de 16 blocs-segment, avec la particularité de séparer la description des n\oe uds d'indices pairs et ceux d'indices impairs. Cela permet de moins s'écarter de l'heuristique d'optimisation en vigueur dans la version séquentielle : les statistiques globales comme la valeur de critère $GL$ y sont recalculées après chaque déplacement, alors que dans le cas d'un traitement parallèle, les statistiques de référence ne sont calculées qu'après le déplacement simultané des $N_n$ n\oe uds; les résultats et la convergence en sont potentiellement affectés, comme l'illustrent les situations présentées à la figure \ref{fig-cycle-contribs-segments}, car les segments composant véritablement le contour après déplacement des n\oe uds n'ont pas nécessairement été évalués (segments en rouge dans la figure \ref{fig-cycle-contribs-segments}).
+
+La figure \ref{fig-structure-segment} représente la structure décrite ci-dessus pour la représentation en mémoire des segments à évaluer. La première ligne montre le détail du premier segment, avec la correspondance \textit{1 pixel par thread} et le découpage en un nombre de blocs suffisant pour permettre de décrire le plus long des segments.
+
+La seconde ligne présente l'ordre dans lequel sont concaténés les 16 groupes de blocs-segment associés au déplacement d'un n\oe ud particulier.
+
+Les deux dernières lignes décrivent la concaténation des ensembles de 16 blocs-segment, avec la particularité de séparer la description des positions des n\oe uds d'indices pairs et ceux d'indices impairs. Cela permet de moins s'écarter de l'heuristique d'optimisation en vigueur dans la version séquentielle où les statistiques globales comme la valeur de critère $GL$ sont recalculées après chaque déplacement (figures \ref{fig-cycle-contribs-segments-a}, \ref{fig-cycle-contribs-segments-b} et \ref{fig-cycle-contribs-segments-c}) .
+
+En version parallèle, si les \og meilleures \fg{} positions de tous les n\oe uds sont calculées simultanément, le contour généré est constitué de segments qui n'ont pas été validés pendant la phase de déplacement des n\oe uds, comme l'illustre la figure \ref{fig-cycle-contribs-segments-e}. La valeur du critère $GL$ doit donc être calculée après coup sur les segments réels du nouveau contour. Dans l'absolu, nous ne sommes donc pas assurés d'améliorer réellement la valeur du critère par rapport au contour de l'itération précédente.
+Pour limiter ce phénomène, qui pourrait provoquer des oscillations et empêcher la convergence, nous avons effectué les déplacements en alternant ceux des n\oe uds d'indices pairs et ceux d'indices impairs. Cela permet de régler le problème lorsque le nombre de n\oe uds du contour est pair. Comme le montrent les figures \ref{fig-cycle-contribs-segments-e} et \ref{fig-cycle-contribs-segments-e}, un segment du contour demeure non validé lorsque le nombre de n\oe uds est impair et nous impose toujours de recalculer, \textit{a posteriori}, la valeur du critère $GL$ pour s'assurer de l'amélioration apporté par les déplacements des n\oe uds. La version parallèle reproduit, malgré cela, assez fidèlement la version séquentielle en effectuant une optimisation n\oe ud par n\oe ud et non une optimisation au niveau du contour complet.
+
\begin{figure}
\centering
- \subfigure[Contour de référence.]{\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-cpu1.png}}\quad
-\subfigure[Déplacement du n\oe ud $N_1$. Le critère est amélioré.]{\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-cpu2.png}}\quad
-\subfigure[Déplacement du n\oe ud $N_2$. Le critère est amélioré.]{\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-cpu3.png}}\\
-\subfigure[Déplacement en parallèle de tous les n\oe uds. Les segments du contour n'ont pas été évalués. On ne peut pas dire, a priori si le critère est amélioré.]{\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-gpu1.png}}\quad
-\subfigure[Déplacement en parallèle des n\oe uds impairs. Le critère est amélioré.]{\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-gpu2.png}}\quad
-\subfigure[Déplacement en parallèle des n\oe uds pairs. Un seul segment n'a pas été évalué.]{\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-gpu3.png}}
-\label{fig-cycle-contribs-segments}
+ \subfigure[Contour de référence.]{\label{fig-cycle-contribs-segments-a}\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-cpu1.png}}\quad
+\subfigure[Déplacement du n\oe ud $N_1$. Le critère est amélioré.]{\label{fig-cycle-contribs-segments-b}\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-cpu2.png}}\quad
+\subfigure[Déplacement du n\oe ud $N_2$. Le critère est amélioré.]{\label{fig-cycle-contribs-segments-c}\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-cpu3.png}}\\
+\subfigure[Déplacement en parallèle de tous les n\oe uds. Les segments du contour n'ont pas été validés. On doit recalculer le critère après les déplacements pour savoir s'il a été amélioré.]{\label{fig-cycle-contribs-segments-d}\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-gpu1.png}}\quad
+\subfigure[Déplacement en parallèle des n\oe uds impairs. Le critère est amélioré.]{\label{fig-cycle-contribs-segments-e}\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-gpu2.png}}\quad
+\subfigure[Déplacement en parallèle des n\oe uds pairs. Un seul segment n'a pas été évalué.]{\label{fig-cycle-contribs-segments-f}\includegraphics[height=3cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cycle-contribs-segments-gpu3.png}}
\caption{Comparaison des cycles de déplacement des n\oe uds. Ligne du haut : version séquentielle. Ligne du bas : version parallèle. Les segments en rouge sont des segments du contour non évalués, alors que ceux en pointillés sont les paires ayant reçu les meilleures évaluations parmi les 8 déplacements possibles des n\oe uds correspondant.}
+\label{fig-cycle-contribs-segments}
\end{figure}
-La représentation en mémoire des segments conduit à avoir un certain nombre non prévisible de threads inactifs dans la grille, sans que cela soit préjudiciable aux performances car cela n'engendre pas de branches divergentes.
+La représentation en mémoire des segments conduit à avoir un certain nombre non prévisible de threads inactifs dans la grille, sans que cela soit préjudiciable aux performances car cela n'engendre pas de branches d'exécution divergentes, qui sont à proscrire sur GPU.
-Les calculs liés à l'évaluation des contributions des pixels sont réalisés en mémoire partagée. Seule une très petite quantité de données doit être stockée en mémoire globale. Il s'agit, pour chaque segment (pas pixel) :
+Les calculs liés à l'évaluation des contributions des pixels sont réalisés en mémoire partagée. Seule une très petite quantité de données doit être stockée en mémoire globale. Il s'agit, pour chaque {\bf segment} :
\begin{itemize}
-\item des coordonnées de son milieu. Cela permet l'ajout efficace quand c'est nécessaire.
+\item des coordonnées de son milieu. Cela permet l'ajout de n\oe ud quand c'est possible, sans calcul supplémentaire.
\item les coordonnées des deux derniers pixels de chaque extrémité. Ils sont nécessaires pour calculer la dérivée aux extrémités et ainsi déterminer le code de Freeman des n\oe uds.
\end{itemize}
-Pour obtenir les contributions des segments, \textit{i.e} les sommes des contributions des leurs pixels, une première phase de réduction partielle est effectuée au niveau de chaque bloc.
+Pour obtenir les contributions des segments, c'est-à-dire les sommes des contributions des leurs pixels, une première phase de réduction partielle est effectuée au niveau de chaque bloc.
-Une synchronisation est alors nécessaire avant d'effectuer les sommes de l'ensemble des contributions partielles qui fournissent les contributions globales des segments. Le contour modifié est alors construit comme la suite des meilleures positions déterminées pour chaque n\oe ud, pour peu que ces nouvelles positions ne générent pas de croisement de segments.
+Une synchronisation est alors nécessaire avant d'effectuer les sommes de l'ensemble des contributions partielles qui fournissent les contributions globales des segments. Le contour modifié est alors construit comme la suite des meilleures positions déterminées pour chaque n\oe ud, pour peu que ces nouvelles positions ne génèrent pas de croisement de segments.
-La solution retenue pour vérifier l'absence de croisement est celle de l'implémentation séquentielle, parallélisée simplement par paire de segments. Cela n'apporte pas de véritable gain de performance par rapport à la version CPU, mais contraints de conserver les données en mémoire GPU pour limiter les transferts entre l'hôte et son périphérique, nous avons tâché de faire en sorte que cette fonctionnalité ne grêve pas les performances globales.
+La solution retenue pour vérifier l'absence de croisement est celle de l'implémentation séquentielle, parallélisée simplement par paire de segments. Cela n'apporte pas de véritable gain de performance par rapport à la version CPU, mais, contraints de conserver les données en mémoire GPU pour limiter les transferts entre l'hôte et son périphérique, nous avons fait en sorte que cette fonctionnalité ne grève pas les performances globales.
-Les calculs des statistiques globales du nouveau contour et du critère $GL$ sont effectués après l'obtention du nouveau contour. Les valeurs obtenues servent de référence pour les prochaines déformations du contour. Les techniques appliquées pour ces calculs sont de nouveau celles décrites au début ce paragraphe.
-Enfin l'ajout des nouveaux n\oe uds se fait simplement pour les segments suffisamment grands, en utilisant les coordonnées des pixels milieux mémorisées lors de la discrêtisation des segments.
+Les valeurs obtenues après détermination du nouveau contour, calcul des statistiques globales et évaluation du critère $GL$, servent de référence pour les prochaines déformations du contour. Les techniques appliquées pour ces calculs sont de nouveau celles décrites au début ce paragraphe.
+Enfin, l'ajout des nouveaux n\oe uds se fait simplement, pour les segments suffisamment grands, en utilisant les coordonnées des pixels milieux mémorisées lors de la discrétisation des segments.
\subsubsection{Cas particulier des segments dont la pente $k$ vérifie $|k|\leq 1$}
-Comme nous venons de le voir, les segments dont la pente $k$ vérifie $|k|\leq 1$ sont discrêtisés à raison de \textit{1 pixel par colonne} et comportent donc le plus souvent plusieurs pixels sur une ligne donnée, comme le montrent les schémas de la figure \ref{fig-segment-k<1}.
-D'après la formulation générale du snake faite au paragraphe \ref{snake-formulation}, le coefficient $C(i,j)$ est à appliquer en chaque point du contour. La technique de discrêtisation employée conduit à des coefficients $C(i,j)$ constants sur l'ensemble des pixels des segments dont la pente $k$ vérifie $|k|> 1$, mais ce n'est pas le cas pour ceux dont la pente $k$ est inférieure ou égale à $1$. Les quatre cas, un par quadrant, qui peuvent se présenter sont représentés à la figure \ref{fig-segment-k<1}.
+Comme nous venons de le voir, les segments dont la pente $k$ vérifie $|k|\leq 1$ sont discrétisés à raison de \textit{1 pixel par colonne} et comportent donc le plus souvent plusieurs pixels sur une ligne donnée, comme le montrent les schémas de la figure \ref{fig-segment-k<1}.
+
+D'après la formulation générale du snake faite au paragraphe \ref{snake-formulation}, le coefficient $C(i,j)$ est à appliquer en chaque point du contour. La technique de discrétisation employée conduit à des coefficients $C(i,j)$ constants sur l'ensemble des pixels des segments dont la pente $k$ vérifie $|k|> 1$, mais ce n'est pas le cas pour ceux dont la pente $k$ est inférieure ou égale à $1$. Les quatre cas, un par quadrant, qui peuvent se présenter sont représentés à la figure \ref{fig-segment-k<1}.
-D'un point de vue opérationnel, on constate en se reportant à la table \ref{tab-freeman}, que tout pixel dont les voisins immédiats sont sur la même ligne à un coefficient $C(i,j)=0$ ($F_{in}=f_{out}=0$). Les deux pixels des extrémités, n'ayant quant à eux qu'un voisin sur la même ligne, ont un coefficient qui dépend du quadrant :
+D'un point de vue opérationnel, on constate en se reportant à la table \ref{tab-freeman}, que tout pixel dont les voisins immédiats sont sur la même ligne, à un coefficient $C(i,j)=0$ ($F_{in}=f_{out}=0$). Les deux pixels des extrémités, n'ayant qu'un voisin sur la même ligne, ont un coefficient qui dépend du quadrant :
\begin{itemize}
\item dans les quandrant 1 et 2
\begin{itemize}
\subfigure[Quadrants 1 et 4]{\includegraphics[width=7cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/coeffs-pixels2.png}}\quad
\subfigure[Quadrants 2 et 3]{\includegraphics[width=7cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/coeffs-pixels1.png}}\\
\subfigure{\includegraphics[width=8cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/coeffs-pixels3.png}}
-\label{fig-segment-k<1}
\caption{Détermination des coefficients $C(i,j)$ des pixels du contour.}
+\label{fig-segment-k<1}
\end{figure}
-Les accès en mémoire, dans les images cumulées, aux contributions des pixels de coefficient $C(i,j)=0$ sont évités et une contribution nulle leur est automatiquement attribuée dès l'étape de discrêtisation au sein du kernel \texttt{GPU\_compute\_segments\_contribs()}.
+Les accès en mémoire aux contributions des pixels de coefficient $C(i,j)=0$, dans les images cumulées, sont évités et une contribution nulle leur est automatiquement attribuée dès l'étape de discrétisation au sein du kernel \texttt{GPU\_compute\_segments\_contribs()}.
\subsection{Performances}
Dans l'hypothèse la plus contraignante d'images en niveaux de gris codés sur 16 bits, l'implémentation parallèle que nous venons de décrire utilise de manière permanente 20 octets par pixel de l'image d'entrée, qui se détaillent en
\begin{itemize}
\item l'image d'entrée pour 4 octets par pixel (1 entier).
-\item l'image cumulée $S_x$ pour 8 octets par pixel (1 entier long)
-\item l'image cumulée $S_x^2$ pour 8 octets par pixel (1 entier long)
+\item l'image cumulée $S_I$ pour 8 octets par pixel (1 entier long)
+\item l'image cumulée $S_{I^2}$ pour 8 octets par pixel (1 entier long)
\end{itemize}
-auxquels il faut ajouter un maximum d'environ 50~Mo d'espace nécessaire à la mémorisation des variables temporaires des calculs et données diverses comme le contour lui-même (n\oe uds, milieus, Freemans, etc.).
+auxquels il faut ajouter un maximum d'environ 50~Mo d'espace nécessaire à la mémorisation des variables temporaires des calculs et des données du contour lui-même : coordonnées des n\oe uds, des milieux des segments, codes de Freeman.
-Sur un GPU de type C1060 disposant de 3~Go de mémoire, cela permet de traiter des image jusqu'à presque 150 millions de pixels.
-Il est possible de réduire cette empreinte jusqu'à 13 octets par pixel, mais cela soulève la question de l'alignement des données en mémoire qui est sans objet en employant les type entier et entier long (32 et 64 bits) pour la représentation des données et qui permet de préserver les performances maximales des opérations et accès aux données du GPU. On pourrait tout de même porter ainsi la limite de taille de l'image d'entrée à 230 millions de pixels.
+Sur un GPU de type C1060 disposant de 3~Go de mémoire, cela permet de traiter des images jusqu'à presque 150 millions de pixels.
+Il est possible de réduire cette empreinte jusqu'à 13 octets par pixel, mais cela soulève la question de l'alignement des données en mémoire, sans objet si on emploie les types entier et entier long (32 et 64 bits) pour la représentation des données et qui permet de préserver les performances maximales des opérations et des accès aux données du GPU. On pourrait tout de même porter ainsi la limite de taille de l'image d'entrée à 230 millions de pixels.
-La convergence de notre implémentation intervient en un nombre généralement plus réduit d'itérations vers un contour final qui diffère par essence de celui obtenu avec la solution de référence. Ces effets sont la conséquence déjà abordée de l'heuristique d'optimisation appliquée à l'implémentation parallèle qui conduit à l'adoption de certains segments non évalués au préalable (voir fig. \ref{fig-cycle-contribs-segments}).
+Comme l'heuristique d'optimisation appliquée à l'implémentation parallèle conduit à la création de certains segments non validés au préalable (voir fig. \ref{fig-cycle-contribs-segments}), notre implémentation peut converger plus tôt que la version de référence CPU, généralement une itération avant.
-Les comparaisons visuelle et de valeur du critère $GL$ qui peuvent être faites pour les images de taille inférieure à 4000$^2$ pixels nous renseignent toutefois sur la qualité de la segmentation obtenue. Pour les tailles au delà et jusqu'au maximum de 12000$^2$ pixels, le comportement est globalement conservé, mais on note qu'il n'est pas pertinent de permettre des tailles de segments trop petites vis à vis de la taille d'image, les déplacements des n\oe uds ne générant alors plus de variations significatives des contributions correspondantes.
+Les comparaisons visuelles et de valeur du critère final $GL$ après convergence nous renseignent toutefois sur la qualité de la segmentation obtenue avec des solutions très voisines de celles de la version de référence.
-La figure \ref{fig-snakegpu-result} présente deux segmentations effectuées sur des images de respectivement 100 et 150 millions de pixels alors que la table \ref{tab-snake-results} résume les performances obtenues sur l'image du \textit{cochon} en différentes tailles.
+La figure \ref{fig-snakegpu-result} présente une segmentation effectuée sur une image de 100 millions de pixels. La table \ref{tab-snake-results} résume les performances obtenues pour différentes tailles de la même image. Une implémentation CPU multi threads permettrait d'accélérer significativement le calcul des images cumulées. Nous n'avons pas toutefois évalué l'accélération réelle qu'une telle solution apporterait. En revanche, on peut affirmer que l'emploi du GPU reste pertinent car aucune solution multi threads classique n'est envisageable pour implémenter la partie segmentation de l'algorithme, la plus coûteuse en temps de calcul et donc celle qu'il est le plus judicieux de paralléliser.
+
\begin{table}[h]
\centering
\begin{tabular}{rrrrr}
\bottomrule
\end{tabular}
- \caption{Comparaison des temps d'exécution de l'implémentation GPU par rapport à l'implémentation CPU de référence, appliqués à une même image dilatée pour en adapter la taille.}
+ \caption{Comparaison des temps d'exécution de l'implémentation GPU (C2070) par rapport à l'implémentation CPU (mono thread) de référence, pour une même image dilatée (fig. \ref{fig-snakecpu-cochon512}) pour en adapter la taille.}
\label{tab-snake-results}
\end{table}
\begin{figure}
\centering
- \subfigure[5 itérations en 0,59~s]{\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cochon_it5_points.png}}\quad
-\subfigure[3 itérations en 0,35~s]{\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/Montserrat_3it.png}}
+ {\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/cochon_it5_points.png}}
+%\quad
+%\subfigure[3 itérations en 0,35~s]{\includegraphics[height=5cm]{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/Montserrat_3it.png}}
+ \caption{Segmentations d'une image de 100~MP en 0,59~s pour 5 itérations. Le contour initial conserve les proportions de celui de la figure \ref{fig-snakecpu-cochon512}. }
\label{fig-snakegpu-result}
- \caption{Segmentations de grandes images, avec le contour intial transposé de celui de la figure \ref{fig-snakecpu-cochon512}. a) image de 100~MP. b) image de 150~MP.}
\end{figure}
-\subsection{Discussion sur l'initialisation}
-Nous avons déjà discuté de l'influence du contour initial sur le résultat de la segmentation, mais il faut ajouter que la durée d'exécution est aussi impactée par le choix du contour initial, dans des proportions qui peuvent être importantes selon la distance, la taille et dans une moindre mesure la forme de la cible.
+\subsection{Détermination du contour initial au sens du maximum de vraisemblance}
+Nous avons déjà discuté de l'influence du contour initial sur le résultat de la segmentation, mais il faut ajouter que la durée d'exécution est aussi impactée par ce choix, dans des proportions qui peuvent être importantes selon la distance, la taille et dans une moindre mesure la forme de la cible.
-Ces effets se mesurent lors de la première itération, celle qui va cerner grossièrement la cible avec un polygone à quatre cotés. Si le contour initial se trouve très éloigné, comme dans la situation de la figure \ref{fig-snakecpu-cochon4kc3}, notre choix maintenant habituel d'un rectangle près des bords de l'image s'avère peu adapté et conduit à une première itération très longue. Dans un tel cas, pour une image de 10000$^2$ pixels, si la cible est un carré de 1000$^2$ pixels dont le sommet du bas à droite se confond avec celui du contour et que l'on approche par pas de 64 pixels, on devra dans le meilleur des cas déplacer les 4 n\oe uds du contour 110 fois de suite avant de pouvoir passer à la deuxième itération. Un pas de 128 permet de réduire ces valeurs, mais l'expérience montre qu'au delà, l'approche initiale de la cible est trop grossière et les itérations suivantes en pâtissent pour un résultat souvent dégradé.
+Ces effets se mesurent lors de la première itération, celle qui va cerner grossièrement la cible avec un contour à quatre n\oe uds. Si le contour initial se trouve très éloigné, comme dans la situation de la figure \ref{fig-snakecpu-cochon4kc3}, notre choix maintenant habituel d'un rectangle près des bords de l'image s'avère peu adapté et conduit à une première itération très longue. Dans un tel cas, pour une image de 10000$\times$10000 pixels, si la cible est un carré de 1000$\times$1000 pixels dont le sommet du bas à droite se confond avec celui du contour et que l'on approche par pas de 64 pixels, on devra dans le meilleur des cas déplacer les 4 n\oe uds du contour 110 fois de suite avant de pouvoir passer à la deuxième itération. Un pas de 128 permet de réduire ces valeurs, mais l'expérience montre qu'au delà, l'approche initiale de la cible est trop grossière et les itérations suivantes en pâtissent pour un résultat souvent dégradé.
En revanche, si les proportions sont celles de la figure \ref{fig-snakecpu-cochon512}, seules 31 passes de déplacement des 4 n\oe uds initiaux sont nécessaires.
-Pour optimiser l'initialisation, nous avons donc proposé de tirer parti du GPU pour évaluer une grande quantité de contours initiaux rectangulaires et réduire ainsi le coût de la première itération. Pour pouvoir employer la mémoire partagée comme tampon de données, il faut limiter le nombre de contours à évaluer. Nous avons donc effectué un échantillonnage spatial des images et déterminé le contour initial en deux temps, en mettant à profit la propriété qu'ont les segments horizontaux d'avoir une contribution nulle. Le principe mis en \oe uvre, illustré par la figure \ref{fig-smart-init} est le suivant :
+Pour optimiser l'initialisation, nous avons donc proposé de tirer parti du GPU pour évaluer une grande quantité de contours initiaux rectangulaires et réduire ainsi le coût de la première itération. Pour pouvoir employer la mémoire partagée comme tampon de données, il faut limiter le nombre de contours à évaluer. Nous avons donc effectué un échantillonnage spatial des images et déterminé le contour initial en deux temps, en mettant à profit la propriété qu'ont les segments horizontaux d'avoir une contribution nulle, comme on peut le vérifier en se reportant à la figure \ref{fig-freeman} et à la table \ref{tab-freeman}. Le principe mis en \oe uvre, illustré par la figure \ref{fig-smart-init} est le suivant :
\begin{enumerate}
\item on réalise un échantillonnage horizontal pour ne considérer que les colonnes d'indice $j=8k$.
\item on évalue alors tous les contours rectangulaires de diagonale $(0, j_L)-(J_H, H)$
\item on fait de même en échantillonnant verticalement : les lignes d'indice $i=8t$ permettent de décrire tous les contours de diagonale $(i_L, j_L)-(i_H, j_H)$. Le meilleur contour est celui retenu pour l'initialisation de la segmentatation.
\end{enumerate}
+Le gain de performance apporté par cette stratégie d'initialisation est variable selon la cible, mais dans des situations favorables comme celle de l'image de la figure \ref{fig-snakecpu-cochon4kc3}, on parvient à une accélération proche de 15 alors qu'elle n'est que d'environ 7 avec l'initialisation basique. Cette proportion est conservée pour les tailles supérieures et signifie que la phase de segmentation est tout de même effectuée 30 fois plus rapidement qu'avec l'implémentation CPU, grâce à une première itération optimisée.
+
\begin{figure}
\centering
\subfigure[Détermination de $j_L$ et $j_H$.]{\resizebox{6cm}{!}{\input{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/smart_init1.pdf_t}}}\quad
\subfigure[Détermination de $i_L$ et $i_H$.]{\resizebox{6cm}{!}{\input{/home/zulu/Documents/these_gilles/THESE/Chapters/chapter3/img/smart_init2.pdf_t}}}
\label{fig-smart-init}
- \caption{Détermination intelligente du contour initial en deux phases successives. a) La première étape repose sur un échantillonnage horizontal. b) La seconde étape repose sur un échantillonnage vertical. }
-\end{figure}
+ \caption{Détermination du contour initial au sens du maximum de vraisemblance, par deux phases successives. (a) La première étape repose sur un échantillonnage horizontal. (b) La seconde étape repose sur un échantillonnage vertical.}
+\end{figure}
\subsection{Conclusion}
-Nous avons conçu une implémentation parallèle de \textit{snake} polygonal orienté régions, ce qui à notre connaissance n'avait encore pas été réalisé, aucune publication n'étant parue à ce sujet.
-
+Nous avons conçu une implémentation parallèle de \textit{snake} polygonal orienté régions, ce qui n'avait pas encore été réalisé, n'ayant recensé à ce jour aucune publication y faisant référence. Elle a fait l'objet d'une publication et d'une communication à la conférence \textit{Computer and Information Technology} (voir \cite{6036776}).
Les objectifs étaient d'étendre les capacités de traitement de l'implémentation CPU de référence en terme de taille d'image en conservant des temps d'exécution acceptables ce qui, de l'avis des auteurs de la version CPU, impose de se situer \textit{a minima} sous la seconde pour pouvoir envisager l'intégration dans une application interactive.
-Sur ce point, les performances de notre version sont satisfaisantes, puisque nous avons repoussé la limite de taille de 16 à 150 millions de pixels et parvenons à segmenter ces grandes images en moins d'une seconde. Le temps de calcul dépend très fortement du contenu de l'image et la segmentation est le plus souvent obtenu en un temps plus court, mais il n'est pas impensable que certaines situations particulières puissent conduire à dépasser cette barre symoblique.
+Sur ce point, les performances de notre version sont satisfaisantes, puisque nous avons repoussé la limite de taille de 16 à 150 millions de pixels et parvenons à segmenter ces grandes images en moins d'une seconde. Le temps de calcul dépend très fortement du contenu de l'image et la segmentation est le plus souvent obtenue en un temps plus court.
+
+L'emploi du GPU dans notre implémentation ne parvient pas à être optimal car, par essence, la répartition des pixels d'intérêt est mouvante et ne permet pas de construire des accès mémoire coalescents. Les opérations de type réduction sont également nombreuses et ne sont pas les plus efficaces sur GPU. Dans notre situation, elles peuvent même représenter une perte de performances, car effectuées sur des vecteurs de tailles insuffisantes.
+
+S'il s'agit de parler d'accélération, notre implémentation divise les temps de traitement précédents par un facteur allant de 6 à 15 selon l'image et le contour initial adopté. Rappelons encore que l'implémentation CPU de référence n'est pas une implémentation naïve, mais une solution optimisée employant déjà les capacités de parallélisme des microprocesseurs modernes et affichant les performances les plus élevées dans ce domaine ; il n'était pas trivial d'en surpasser les performances, même avec un GPU.
-L'emploi du GPU dans notre implémentation ne parvient pas à être optimal car, par essence, la répartition des pixels d'intérêt est mouvante et ne permet pas de construire des accès coalescent à la mémoire. Les opérations de type réduction sont également nombreuses et ne sont pas les plus efficaces sur GPU. Dans notre situation, elles peuvent même représenter une perte de performances car effectuées sur des vecteurs de tailles insuffisantes pour que le GPU surclasse le CPU.
+Par nécessité, notre solution s'écarte cependant quelque peu de l'algorithme original pour permettre les déplacements simultanés des sommets d'indices pairs, puis impairs, du polygone. Ce faisant, on observe parfois la convergence de notre solution à l'avant dernière étape de la segmentation réalisée par la version séquentielle de référence itération, sans que cela n'influe significativement sur la qualité. En effet, seuls quelques n\oe uds voient leur position potentiellement modifiée d'un seul pixel (le pas de déplacement des dernières étapes) et le contour obtenu ne s'éloigne donc que très peu du contour obtenu par l'algorithme de référence.
-S'il s'agit de parler d'accélération, notre implémentation divise les temps de traitement précédents par un facteur allant de 6 à 15 selon l'image et le contour initial adopté. Rappelons encore que l'implémentation CPU de référence n'est pas une implémentation naïve, mais une solution optimisée employant déjà les capacités de parallélisme des microprocesseurs modernes et représentant l'\textit{l'état de l'art} du domaine ; il n'était pas trivial d'en surpasser les performances, même avec un GPU.
+La technique que nous avons proposée pour la détermination intelligente du contour initial permet d'augmenter encore les performances, surtout dans les grandes images lorsque la cible est petite vis-à-vis des dimensions de l'image. Il reste toutefois à concevoir une technique permettant de prévoir si cette recherche de contour initial serait génératrice de gain de performance.
-Par nécessité, notre solution s'écarte cependant quelque peu de l'algorithme original pour permettre les déplacements simultanés des l'ensemble des sommets du polygone. Ce faisant, la décroissance du critère n'est pas certaine à toutes les étapes de la segmantation et l'on observe cette conséquence en particulier lors des dernière itérations lorsque le pas de déplacement et aussi les variations du critère sont faibles. Ce comportement, lorsqu'il est observé, provoque parfois la convergence prématurée de la segmentation, mais n'influe toutefois que sur quelques n\oe uds et dans la mesure d'un pixel.
+L'analyse fine des séquences de segmentation montre enfin que les première étapes, qui mettent en \oe uvre les segments les plus longs, génèrent des grilles de calcul suffisamment chargées et homogènes pour présenter de bonnes performances. Les dernières étapes, en revanche, traitent un grand nombre de petits segments, générant beaucoup de trous dans la grille de calcul et induisant des performances moindres.
-La technique que nous avons proposée pour la détermination du contour initial permet d'augmenter encore les performances, surtout dans les grandes images lorsque la cible est petite vis à vis des dimensions de l'image. Nous ne sommes pas parvenu à concevoir une technique permettant de prévoir si la recherche intelligente du contour intial serait génératrice de gain de performance.
+Pour résumer, l'accélération globale obtenue est principalement déterminée par le calcul des images cumulées et des toutes premières étapes de déplacements. Une possibilité à explorer serait de construire une version hybride réalisant le début de la segmentation sur GPU, puis la terminant sur le CPU hôte. Ceci est envisageable en raison du très petit volume de données à transférer que constituent les paramètres du contour (2~ko pour 100 n\oe uds). Cette étude nous conforte également dans l'idée que la transposition pour GPU d'algorithmes séquentiels optimisés pour CPU, malgré des adaptations à l'architecture, ne semble pas être la démarche permettant d'atteindre les niveaux de performances attendus lorsqu'on met en \oe uvre ces processeurs graphiques.
-L'analyse fine des séquences de segmentation montre enfin que les première étapes, qui mettent en \oe uvre les segments les plus longs, générent des grilles de calcul suffisamment chargées et homogènes pour présenter de bonnes performances. Les dernière étapes, en revanche, traitent d'une plus grand nombre de petits segments, générant beaucoup de trous dans la grille de calcul et induisant des performances moindres. L'accéleration globale obtenue est ainsi généralement le fruit du calcul des images cumulées et des toutes premières étapes de déplacements. Une possiblité qui reste à explorer serait de construire une version hybride réalisant le début de la segmentation sur GPU, puis terminant sur le CPU hôte. Ceci est envisageable en raison du très petit volume de données à transférer que constituent les paramètres du contour (2 ko pour 100 n\oe uds).
