Organisation du cours

• 14 séances d’1h30 (CM/TP)
  • Implémentation en Python utilisant les librairies NumPy et SciPy
  • Utilisation des notebook jupyter
• Évaluation :
  • 1 ou 2 TP notés (note individuelle)
  • 1 projet en binôme
    • étude bibliographique, implémentation
    • rapport écrit, présentation orale
• Enseignant : Emmanuel.Caruyer@irisa.fr
  • CR CNRS, Équipe-projet Empenn
  • Traitement d’images médicales, neuroimagerie

Ressources

• Site web du cours : http://emmanuelcaruyer.com/imagerie-numerique.php
  • Transparents
  • Exemples de cours (Assurez-vous que vous pouvez les exécuter et les comprendre)
  • Sujets de TD et code source à compléter si besoin
• Bibliographie (complétée pendant le cours)
  • Stéphane Mallat. 2008. A wavelet tour of signal processing: the sparse way. Academic Press.

A wavelet tour of signal processing

Compression d’images

Une image compressée en JPEG et en JPEG2000 (même taille de fichier)

[source : wikibook https://fr.wikibooks.org/wiki/Le_format_JPEG_2000]

Progrès liés :

• à la représentation mathématique des images
• à la compréhension du système visuel humain

Débruitage (1)

Imagerie par résonance magnétique

Imagerie radar

Photographies numériques (différents réglages ISO)

Débruitage (2)

\[I[i,j] = I_0[i, j] + \epsilon, \qquad \epsilon \sim \mathcal{N}(0, \sigma^2)\]

Rapport signal sur bruit (Signal to Noise Ratio) \[SNR = \frac{1}{nm} \sum_{i=1}^n \sum_{j=1}^m I[i,j] / \sigma\]

Problèmes inverses (1) – Inpainting

Un masque de la grille à effacer

[source : G’MIC, David Tschumperlé]

Problèmes inverses (2) – super-résolution

800x600 (Taille originale de l’image: 381x231px)

source : G’MIC, David Tschumperlé

Problèmes inverses (3) – défloutage (deblurring)

Exemple de flou : bokeh – partie d’une photo hors zone de mise au point

Flou de mouvement

Problèmes inverses (4) – Tomodensitométrie (CT-scan)

Tomodensitométrie (CT scan Computerized Tomography)

Problèmes inverses (5) – Électro-encéphalographie (EEG)

\[z = A x + \epsilon, \qquad z \in \mathcal{R}^d, x \in \mathcal{R}^p, \qquad d \ll p\]

Localisation de source en EEG

[source : Jérémie Mattout]

Segmentation d’objets (1)

IRM cérébrale et tumeur segmentée

source: BRATS segmentation challenge

Segmentation d’objets (2)

source: DAVIS 2017 challenge

Segmentation d’objets (2)

source: DAVIS 2017 challenge

Représentation d’une image numérique

Image analogique (niveaux de gris) : \(f \in L^2([0, 1]^2)\)

Cas général : \(f: \Omega = [0, 1]^d \rightarrow [0, 1]^s\)

• \(d=2, s=3\) : images couleur
• \(d=3, s=3\) : video
• \(d=4, s=1\) : Images médicales (imagerie fonctionnelle), …

Représentation continue vs discrète, finie

En pratique, \(\Omega \subset \mathbb{Z}^d\) mais…

• Une représentation continue donne un sens à la notion de « régularité »
• Les transformées usuelles (Fourier, Ondelettes, …) sont plus facile à définir dans le domaine continu
• On peut généraliser certains des résultats du domaine continu au domaine discret asymptotiquement
• A contrario, certaines méthodes sont intrinsèquement discrètes (par exemple on peut voir \(\Omega\) comme un graphe)

De même, les valeurs de l’image sont discrètes :

• l’affichage se fait sur 255 niveaux de (rouge, vert, bleu)
• la quantification des niveaux de l’image peut être optimisée pour compresser l’image.

Espace de fonctions

• Espace des fonctions de carré intégrable \(\mathbf{L}^2(\mathbb{R})\)
• Produit scalaire : \[\langle f, g\rangle = \int_{-\infty}^{\infty} f(t) g^*(t)\mathrm{d}t\]
• base orthonormale (base de Hilbert): \(\{\psi_n(t), n\in \mathbb{N}\}\) telle que
  • \(i\neq j \implies \langle \psi_i, \psi_j \rangle = 0\)
  • \(\forall i \in \mathbb{N}, ||\psi_i||^2 = \langle\psi_i, \psi_i\rangle = 1\)
  • \(\forall f \in \mathbf{L}^2(\mathbb{R}), \exists (\alpha_n)_{n \in \mathbb{N}}\ \mathrm{ t.q. }\ \lim_{N} || f - \sum_{n=0}^N \alpha_n \psi_n|| = 0\)
• Décomposition dans une base orthonormale
  • La décomposition est unique \(\forall f \in \mathbf{L}^2(\mathbb{R}), f = \sum_{n=0}^\infty \langle f, \psi_n\rangle \psi_n\)
  • L’« énergie » est conservée (Parseval) \(\int_{-\infty}^\infty |f(t)|^2\mathrm{d}t = \sum _{n=0}^\infty |\langle f, \psi_n\rangle|^2\)

Filtre linéaire et convolution (1)

• Un filtre \(L\) est une transformation d’une fonction \(f\) vers une autre fonction, notée \(Lf\).
• On dit qu’un filtre \(L\) est linéaire si l’application \(f\longmapsto Lf\) est linéaire
• Soit \(f_\tau(t) = f(t - \tau)\) (\(f\) « décalée » de \(\tau\)). On dit qu’un filtre \(L\) est invariant par translation si l’application de \(L\) à \(f_\tau\) est telle que \(Lf\,(t-\tau) = Lf_\tau (t)\).

Filtre linéaire et convolution (2)

Exemple: moyennage uniforme sur une fenêtre \([t - \alpha, t + \alpha]\) \[Lf(t) = \frac{1}{2\alpha} \int_{t-\alpha}^{t+\alpha} f(t')\mathrm{d}t'\]

Filtre linéaire et convolution (3)

• La réponse impulsionnelle d’un filtre \(L\) est définie par \[h(t) = L\delta (t)\] (\(\delta\) est la distribution de Dirac).
  • NB : la distribution de Dirac correspond à une mesure de masse 1, dont le support est réduit à \(\{0\}\).
• Un filtre linéaire, invariant par translation est entièrement défini par sa réponse impulsionnelle : \[\forall f, Lf (t) = \int_{-\infty}^\infty h(u) f(t - u) \mathrm{d}u = h * f(t)\]
• L’opération \(h * f\) définie ci-dessus s’appelle le produit de convolution.

Filtre linéaire et convolution (4)

Exemple: moyennage uniforme sur une fenêtre \([t - \alpha, t + \alpha]\) \[Lf(t) = \frac{1}{2\alpha} \int_{t-\alpha}^{t+\alpha} f(t')\mathrm{d}t'\]

• La réponse impulsionnelle de \(L\) est \[h(t) = L\delta(t) = \frac{1}{2\alpha}\int_{t-\alpha}^{t+\alpha} \delta(t')\mathrm{d}t' = \begin{cases} 1/2\alpha & \text{si } t \in [-\alpha, \alpha] \\ 0 & \text{sinon.} \end{cases} \]

Polynôme trigonométrique et convolution

• Les polynômes trigonométriques \(e^{i \omega t}\) sont des valeurs propres de la convolution : \(f \mapsto Lf = h * f\). \[L e^{i\omega t} = \int h(u) e^{i\omega (t-u)} \mathrm{d}u = e^{i\omega t} \int h(u) e^{-i\omega u} \mathrm{d}u \]
• On reconnait la transformée de Fourier : \[ L e^{i\omega t} = \hat{h}(\omega) e^{i\omega t}\textrm{, avec } \hat{h}(\omega) = \int_{-\infty}^\infty h(u) e^{-i \omega u} \mathrm{d}u. \]
• On peut réécrire un filtre linéaire en passant par la transformée de Fourier \[Lf(t) = f * h(t) = \frac{1}{2\pi} \int_{-\infty}^{\infty} \hat{h}(\omega) \hat{f}(\omega) \exp(i\omega t) \mathrm{d}\omega\]

Le cas discret et fini (1)

• En pratique les signaux numériques sont finis et discrets : \(\bar{f}[n]\) est défini pour \(n\in \{0, 1, \ldots N - 1\}\).
• Pour définir une convolution il faut préciser les conditions aux bords :
  • par périodicité, on peut définir \(f[n] = \bar{f}[n \mod N]\)
• La convolution circulaire avec \(\bar{h}[n]\) définie pour \(n \in \{0, 1, \ldots N - 1\}\) est définie par \[f * h[n] = \sum_{p=0}^{N-1} f[p] h[n-p] \]

Le cas discret et fini (2)

• On définit la transformée de Fourier de \(f\) par \[\hat{f}[k] = \sum_{n=0}^{N-1} f[n] \exp\left(-\frac{2i\pi kn}{N}\right)\]
• La transformée de Fourier inverse est définie par \[f[n] = \frac{1}{N} \sum_{k=0}^{N-1} \hat{f}[k] \exp\left(\frac{2i\pi kn}{N}\right)\]
• Coût de calcul: \(\mathcal{O}(N \log N)\) en utilisant la transformée de Fourier rapide (ou fast Fourier transform, FFT).

Transformée de Fourier d’une image

• Soit \(f[m,n]\) une image définie pour \(m\in\{0, \ldots, M-1\}\) et \(n\in\{0,\ldots N-1\}\), la transformée de Fourier de \(f\) est définie par \[ \hat{f}[k, l] = \sum_{m=0}^{M-1} \sum_{n=0}^{N-1} f[m,n] \exp\left(\frac{-2i\pi km}{M}\right) \exp\left(\frac{-2i\pi ln}{N}\right) \] pour \(k\in \{0,\ldots, M-1\}\), \(l\in\{0,\ldots N-1\}\).

• On peut retrouver \(f\) à partir des coefficients \(\hat{f}\) en appliquant la transformée de Fourier inverse \[ f[m,n] = \frac{1}{MN} \sum_{k=0}^{M-1} \sum_{l=0}^{N-1} \hat{f}[k,l] \exp\left(\frac{2i\pi km}{M}\right) \exp\left(\frac{2i\pi ln}{N}\right) \]