5.5 Filtrage.

5.5.1 Quand filtrer ?

Nous venons de constater (Fig.6.13) que la corrélation sur le bruit ne dépend pas uniquement de la géométrie d'acquisition: elle varie aussi avec la forme et l'intensité de l'objet lui-même (Par.6.3.2.4). Cette dépendance vis à vis de l'objet imagé rend difficile la modélisation du bruit sur les projections. Le noyau de filtrage à utiliser varie spatialement et sa détermination doit se faire pour chaque objet imagé. C'est pourquoi, si on conserve l'hypothèse que c'est principalement la projection qui induit une corrélation variant spatialement et si on suppose également que la corrélation pré-projection est invariante spatialement, alors nous pouvons envisager un filtrage du volume servant à la projection par un noyau fixe. On peut donc filtrer par simple multiplication dans l'espace de Fourier de l'objet. C'est le fait de filtrer le volume avant projection qui va introduire la corrélation faisant défaut au paragraphe Par.6.4.2.5.

D'autre part, on considère que la variance du bruit sur le volume dépend de la distribution radioactive initiale. Nous bruitons le volume avant de le filtrer. Ce filtrage, peut en quelque sorte matérialiser le parcours du $\beta$. Les répercussions sur le protocole pour construire un sinogramme de bruit sont illustrées Fig.6.14.

Figure: Protocole pour créer numériquement un sinogramme de bruit.

5.5.2 Quel noyau de filtrage ?

5.5.2.1 Forme générale du noyau.

Une deuxième question à laquelle il est nécessaire de répondre est la nature du noyau de filtrage que nous allons utiliser. Le choix du noyau de filtrage est un problème délicat, car le bruit en modifiant les valeurs du volume influe sur la corrélation des projections. De plus, nous avons accès uniquement aux projections bruitées et non au volume. A cette difficulté s'ajoutent des problèmes liés à la discrétisation notamment dans le calcul des FFT.

On va définir un noyau pour chaque direction $N_{x}$, $N_{y}$ et $N_{z}$. Les noyaux dans un plan transaxial seront les mêmes ( $N_{xy}=N_{x}=N_{y}$). Le noyau global de filtrage tridimensionnel sera construit sur ces noyaux monodimensionnels:

$$N(x,y,z)=N_{xy}(x)\times N_{xy}(y)\times N_{z}(z)$$

Par la suite $N$ désigne l'un de ces noyaux monodimensionnels. Chaque noyau $N(x)$ sera décrit par 3 paramètres $(p_{1},p_{2},p_{3})$. Comme on filtre avant projection, on fait l'hypothèse que les paramètres de ce noyau ne dépendent ni de l'objet, ni du nombre de coups dans le volume. On va chercher les paramètres à appliquer à tout volume émetteur, quel que soit son niveau de bruit et sa forme, pour que les projections calculées à partir du volume de bruit filtré présentent les mêmes caractéristiques en terme de corrélation que les projections réelles. Autrement dit, on fait l'hypothèse que le variation spatiale de la corrélation provient essentiellement de la projection.

Tous les noyaux de filtrage utilisés sont des courbes en ``cloche'' permettant de lisser le volume initial avant projection. Elles se caractérisent de la façon suivante:

$$\lim _{\vert x\vert\rightarrow \infty }N(x) \rightarrow$$ 0

$$\lim _{\vert x\vert\rightarrow 0}N(x) \rightarrow M$$

où $M$ est un nombre fini.

• Le premier paramètre $p_{1}$ est toujours calculé numériquement de manière à normaliser le filtre :
  $$\int ^{bm}_{-bm}N(x)dx\xrightarrow {bm\rightarrow \infty }1$$
  où $bm$ est une borne maximale d'intégration. On cherche à garder la même moyenne sur le volume après filtrage.
• Le deuxième paramètre $p_{2}$ permet de jouer sur la largeur à mi-hauteur du noyau.
• Le troisième paramètre $p_{3}$ va influencer la vitesse de convergence du noyau vers $0$ lorsque $x$ augmente.

Dans notre étude, seuls deux noyaux ont été envisagés.

• Le noyau exponentiel $N_{e}$ défini par:
  $$N_{e}(x)=p_{1}\times \exp -(\frac{\vert x\vert}{2p_{2}})^{p_{3}}$$

• Le noyau de type sigmoïde $N_{s}$ défini par:
  $$N_{s}(x)=\frac{p_{1}}{1+(\frac{x}{p_{2}})^{p_{3}}}$$

Si ce deuxième type de noyau a été envisagé en plus des noyaux exponentiels, c'est qu'il génère des corrélations qui ne décroissent pas aussi vite lorsque $x$ augmente. La forme générale de la corrélation ressemble au profil d'autocorrélation obtenu avec ce type de noyau. Pour réduire encore le champ d'investigation, nous avons considéré que $p_{3}$ est à valeur entière. Seules les valeurs $p_{3}=2,3$ pour le noyau sigmoïde et $p_{3}=2$ pour le noyau exponentiel ont été envisagées. D'autres cas ont été essayés mais ne seront pas exposés, car ils ne conduisaient pas à des résultats probants. Nous verrons pourquoi au Par.6.5.3.1.

5.5.2.2 Estimation des paramètres des noyaux.

On considère le sinogramme comme une série de plans de projection $\Pi _{\theta ,\phi }$ associés à la projection suivant $\vec{n}$. Afin de déterminer les paramètres et le type de noyau à utiliser, on va faire varier le triplet $(p_{1},p_{2},p_{3})$. Pour chaque triplet, on construit un sinogramme de bruit en accord avec le protocole défini Fig.6.14. Les profils d'autocorrélation suivant deux directions orthogonales (suivant $\vec{s}$ et $\vec{r}$) dans chaque plan sont calculés puis moyennés sur l'ensemble des directions $\theta$ considérées. $acCs$ et $acCr$ correspondent aux profils moyens déterminés sur le sinogramme de bruit calculé et $acRs$ et $acRr$ à ceux déterminés sur le sinogramme de bruit calculé à partir de données réelles. Le critère $\chi _{i}$ retenu pour choisir le meilleur filtre est calculé à partir de ces profils.

$$\chi _{i}=\frac{(acRi-acCi)^{2}}{acRi^{2}}$$ (5.10)

où $i$ est soit $s$ soit $r$ suivant que le profil est obtenu suivant une direction associée à $\vec{s}$ ou $\vec{r}$.

En pratique, pour déterminer l'ensemble des constantes $p_{1},p_{2}$ et $p_{3}$, nous fixons la vitesse de décroissance $p_{3}$. Puis, pour différentes valeurs du critère $p_{2}$ nous calculons par intégration la constante de normalisation $p_{1}$. Nous estimons alors la valeur du critère défini par Eq.6.10. La détermination du minimum est graphique. Le triplet optimal nous est fourni par la valeur de $p_{2}$ qui minimise ce critère.

5.5.2.2.0.1 Remarque:

Pour chaque critère $p_{2}$, il est nécessaire de recalculer une valeur optimale de $k_{2}$ . Autrement dit, pour chaque détermination de $\chi _{i}(p_{2})$, nous estimons la fonction $D_{KS}=f(k_{2})$, nous la minimisons et c'est pour cette valeur optimale qu'est estimé le critère $\chi _{i}$.

5.5.3 Variation du modèle en fonction du niveau de bruit.

Su les figures Fig.6.15.a,b,c

Figure 6.15: Variation de $\chi _i$ en fonction de $p_2$. On représente en vert la variation du critère suivant la direction $\vec{r}$ et en bleu celle suivant la direction $\vec{s}$.

[Noyau gaussien.]	[Noyau sigmoïde avec $p_3=2$.]	[Noyau sigmoïde avec $p_3=3$.]

sont représentées les variations du critère $\chi _{i}$ pour les cas envisagés. Globalement, nous observons sur chaque courbe un comportement analogue. En revanche, on n'observe pas réellement de valeur optimale. On constate tout d'abord, une zone plus ou moins constante du critère, puis, au delà d'une valeur seuil de $p_{2}$, le critère croit régulièrement. On a donc l'impression que dès l'introduction d'un filtrage ($p_{2}=0.1$ fois la taille du voxel!), le filtre conduit à une corrélation proche de celle observée sur le sinogramme réel. Ceci n'est pas vraiment le cas, car les profils d'autocorrélation montrent, pour de petits écarts de $\Delta r$ ou de $\Delta s$, une corrélation calculée nettement supérieure à celle observée dans le cas réel . Le critère n'est donc pas tout à fait optimal pour la sélection de la valeur $p_{2}$. On considère toutefois que la valeur optimale de $p_{2}$ correspond au moment où le critère décroche du plateau initial. Cette démarcation est très nette dans le cas du noyau gaussien (Fig.6.17.a noyau exponentiel avec $p_{3}=2$). Elle l'est moins pour l'autre type de noyau envisagé (Fig.6.17.b et c). Dans le tableau Tab.6.5.

Table 6.5: Valeurs optimales du triplet $(p_{1},p_{2},p_{3})$ et constantes $(k_{1},k_{2})$ correspondantes.

	$N_{s}$		$N_{e}$
$p_{1}$	0.54	1.03	0.57
$p_{2}$(voxels)	0.6	0.4	0.7
$p_{3}$	2	3	2
$k_{1}$	$1.17\times 10^{-4}$	$1.17\times 10^{-4}$	$1.17\times 10^{-4}$
$k_{2}$	$5.3\pm 0.5$	$2.9\pm 0.2$	$2.9\pm 0.2$

sont consignées les valeurs optimales pour chacun des cas considérés. Sur les figures Fig.6.16 et Fig.6.17 on donne, pour chaque cas et pour un plan $\Pi _{\theta ,\phi }$ particulier, les histogrammes

Figure: Histogramme des sinogrammes de bruits pour différents types de noyaux construits avec les triplets optimaux $(p_1,p_2,p_3)$ (Tab.6.5). Les distributions sont normalisées par l'écart type du bruit sur les sinogrammes (Eq.6.9). On représente en bleu le sinogramme de bruit calculé et en vert le sinogramme de bruit réel.

[Noyau gaussien.]	[Noyau sigmoïde avec $p_3=2$.]	[Noyau sigmoïde avec $p_3=3$.]

et les profils d'autocorrélation.

Figure: Profils d'autocorrélation sur sinogramme de bruit calculé ou réel pour différents types de noyaux construits avec les triplets optimaux $(p_1,p_2,p_3)$ (Tab.6.5). On représente les profils suivant la direction $\vec{r}$ (à droite) et suivant la direction $\vec{s}$ (à gauche) sur la même figure. Les profils vert (suivant $\vec{r}$) et rouge (suivant $\vec{s}$) correspondent aux profils d'autocorrélation calculés pour le sinogramme de bruit réel ($b^{p}_{v}$), les profils bleu (suivant $\vec{r}$) et violet (suivant $\vec{s}$) correspondent aux profils d'autocorrélation calculés pour le sinogramme de bruit calculé ($b^{p}_{b}$)

[Noyau gaussien.]	[Noyau sigmoïde avec $p_3=2$.]	[Noyau sigmoïde avec $p_3=3$.]

5.5.3.1 Type de filtre et Discrétisation.

Nous venons de constater que les filtrages ``optimaux'' s'obtiennent pour des valeurs très petites du paramètre $p_{2}.$ Par conséquent, les noyaux se traduisent par des fonctions représentées par des pics très étroits (entre $0.1$ et $1.5$ fois la taille du voxel). De ce fait, pour des petites valeurs du paramètre $p_{2}$ et une fois discrétisés, ces noyaux ont tous des profils similaires où seul le voxel central est non nul (Fig.6.18).

Figure 6.18: Noyaux gaussien $N_{e}(x)\vert _{p_{3}=2}$ et sigmoïde $N_{s}(x)\vert _{p_{3}=2}$ pour $p_{2}=0.1$.

Les corrélations induites sont, dans ce cas, semblables. Il faut faire croître le paramètre $p_{2}$ pour que les profils se différencient. D'autre part, les profils des courbes, pour chacun des deux noyaux, ne varient pas de la même manière en fonction du paramètre. Le noyau sigmoïde s'élargit plus rapidement que le noyau gaussien lorsque $p_{2}$ augmente. Sur les figures 6.19 et Fig.6.20, on représente les fonctions pour différentes valeurs du paramètre $p_{2}$.

Figure 6.19: Noyau $N_s(x)\vert _{p_3=2}$ pour $p_2=[0.1\,(vert),\,0.2\,(bleu),\,0.4\,(rouge),\,0.6\,(violet)]$.

Figure 6.20: Noyau $N_e(x)\vert _{p_3=2}$ pour $p_2=[0.1\,(vert),\,0.5\,(bleu),\,1.0\,(rouge),\,1.5\,(violet)]$.

Pour améliorer les choses, il faudrait construire des filtres plus précis en réduisant le pas d'échantillonnage spatial, c'est à dire diminuer la taille des voxels des volumes utilisés pour la projection.

5.5.3.2 Normalité de l'histogramme.

La distribution de bruit, une fois normalisée devrait suivre une loi normale (de moyenne nulle et d'écart type égal à 1). Nous allons donc comparer l'histogramme de la distribution de bruit $b^{p}_{b}$ à celui obtenu pour une distribution de bruit correspondant à une densité de probabilité normale.

5.5.3.2.1 Cas sans filtrage.

Dans un premier temps considérons la projection en l'absence de filtrage. Dans ce cas, nous constatons (Fig.6.21.a)

Figure 6.21: Comparaison des histogrammes avec la loi normale. On représente en vert l'histogramme normalisé (Eq.6.9) du sinogramme de bruit calculé $b^{p}_{b}$ et en bleu celui qui correspond à une distribution suivant une loi normale.

[Sans filtrage]	[Avec filrage gaussien]

que l'histogramme calculé s'étend plus qu'une loi normale ( $\sigma _{b_{b}^{p}}>1$). Il faut chercher l'explication dans le fait que pour construire les plans de détecteurs, nous intégrons l'information suivant plusieurs inclinaisons. Supposons donc que $b^{p}_{b}\vert _{l}$ représente la loi statistique suivie par le sinogramme de bruit pour une inclinaison particulière $\phi _{l}$. Pour simplifier, nous envisageons le cas où toutes ces lois sont de variance $\sigma _{b_{b}^{p}}$. Du fait, de l'écartement, nous intégrons dans un plan de détecteurs l'information provenant de $m$ inclinaisons différentes (où $m-1$ suivant la parité du plan et du segment considéré cf Ch.4). La loi statistique suivie par le sinogramme de bruit correspond donc à la somme des lois $b^{p}_{b}\vert _{l}$ pour les $m$ angles $\phi _{l}$ contribuant au sinogramme de bruit.

$$\overline{b^{p}_{b}}=\sum ^{m}_{l=1}b^{p}_{b}\vert _{l}$$

Donc, si maintenant nous cherchons la variance de la loi $\overline{b^{p}_{b}}$ effectivement suivie par le sinogramme de bruit, nous obtenons:

$$var(\overline{b^{p}_{b}})=m.\sigma ^{2}_{b_{b}^{p}}+\sum _{l_{1}\neq l_{2}}cov(b^{p}_{b}\vert _{l_{1}},b^{p}_{b}\vert _{l_{2}})$$

D'un autre côté, les variables $b^{p}_{b}\vert _{l}$ sont corrélées et vu les angles $\phi _{l}$ voisins considérés, les lois $b^{p}_{b}\vert _{l}$ ont tendance à fluctuer dans la même direction. Le terme de covariance sera donc positif. En l'absence de corrélation, nous aurions:

$$\sum _{l_{1}\neq l_{2}}cov(b^{p}_{b}\vert _{l_{1}},b^{p}_{b}\vert _{l_{2}})=0$$ (5.11)

Et si nous avions $m$ fois la même loi, nous aurions:

$$\sum _{l_{1}\neq l_{2}}cov(b^{p}_{b}\vert _{l_{1}},b^{p}_{b}\vert _{l_{2}})=m(m-1).\sigma _{b_{b}^{p}}$$ (5.12)

Les équations Eq.6.11 et Eq.6.12 nous donnent respectivement un minorant et un majorant pour ce terme de covariance. Ainsi, nous pouvons trouver un intervalle pour la variance dans un plan de détecteur qui vaut:

$$m.\sigma _{b_{b}^{p}}^{2}\leq var(\overline{b^{p}_{b}})<m^{2}.\sigma _{b_{b}^{p}}^{2}$$

La variance que nous observons réellement dans un plan est supérieure à celle que nous avons pour une seule inclinaison. Autrement dit, dans la normalisation effectuée (Eq.6.9), nous ne normalisons pas suffisamment (Fig.6.21.a). La normalisation devrait être:

$$\sigma _{\overline{b^{p}_{p}}}^{2}(r,s,\theta ,\phi )=k_{1}k^{2}_{2}k_{3}p(r,s,\theta ,\phi )$$

où $1<\sqrt{m}<k_{3}<m$. Nous n'avons pas besoin de calculer cette constante pour construire un sinogramme, car elle s'introduit de façon naturelle dans notre projection. Ce défaut sur la normalisation ne modifie pas la nature de la gaussienne obtenue. Il y a juste que nous ne pouvons pas comparer l'histogramme des valeurs avec la loi normale de façon directe.

5.5.3.2.2 Cas avec filtrage.

Dans tous les développements précédents, nous avons accepté l'hypothèse de non-corrélation sur notre bruit (Eq.6.8). Or, le filtrage effectué avant projection introduit une corrélation et infirme cette hypothèse. Toutefois, ce filtrage étant linéaire, nous n'altérons pas la nature du bruit sur le sinogramme (l'histogramme des valeurs reste gaussien). Cependant, le filtrage réduit la variance du sinogramme de bruit de sorte que la normalisation nous conduit à une constante $k_{3}<1$. (Fig.6.21.b). Pour estimer la variable $k_{3}$, nous cherchons a posteriori la largeur à mi-hauteur de l'histogramme de bruit que nous rapportons à la hauteur à mi-hauteur de l'histogramme relatif à la loi normale.

5.5.4 Conclusion.

Nous pouvons donc maintenant résumer le noyau de filtrage que nous utiliserons. Au vu des résultats précédents, nous n'envisageons pas le même type de noyau suivant la direction ($\vec{r}$ ou $\vec{s}$). Dans la direction portée par $\vec{s}$ nous utilisons un noyau gaussien alors que suivant la direction $\vec{r}$ nous envisageons un noyau sigmoïde. Le tableau Tab.6.6 résume, toujours pour le sinogramme $_{7}p_{v}$, l'intégralité des constantes.

Table 6.6: Constantes optimales pour le sinogramme $_{7}p_{v}$.

				$N_{e}(z)$			$N_{s}(x)$, $N_{s}(y)$
	$k_{1}$	$k_{2}$	$k_{3}$	$p_{1}$	$p_{2}$	$p_{3}$	$p_{1}$	$p_{2}$	$p_{3}$
$_{7}p_{v}$	$1.17\times 10^{-4}$	2.80	0.68	0.57	0.7	2	1.03	0.4	3