Partie IV - Calcul scientifique · Chapitre 19

Introduction au traitement d’images

3 sections de cours · fiche d'exercices corrigee

1L'image comme tableau : types, pixels, transformations ponctuelles

Definition — Image numerique

Une image numerique est un tableau de pixels (picture elements). Chaque pixel porte une valeur codant une intensite lumineuse ou une couleur. Une image est donc une matrice, et le traitement d'images est de l'algebre lineaire appliquee.

TypeStructureValeur d'un pixel
Binaire\((h, w)\)0 (noir) ou 1 (blanc)
Niveaux de gris\((h, w)\)entier de 0 (noir) a 255 (blanc)
Couleur RVB\((h, w, 3)\)triplet \((R, V, B)\), chacun de 0 a 255
RVBA\((h, w, 4)\)+ canal alpha (transparence)
Convention d'indexation

img[i, j] designe la ligne \(i\) (de haut en bas) et la colonne \(j\) (de gauche a droite). L'origine est en haut a gauche, et l'axe des \(i\) descend : c'est l'inverse de la convention mathematique habituelle. Une confusion entre \(i\) et \(j\) produit une image transposee.

import numpy as np

img = np.zeros((100, 200), dtype=np.uint8)     # image noire, 100 lignes x 200 colonnes
h, w = img.shape
img[40:60, 80:120] = 255                       # carre blanc au centre
img[0, :] = 128                                # premiere ligne grise
Le type uint8 et le debordement
a = np.array([200, 100], dtype=np.uint8)
a + 100          # [ 44 200]  :  300 mod 256 = 44 !
a - 250          # [206 106]  :  -50 mod 256 = 206

Les entiers uint8 s'enroulent silencieusement. La regle absolue : convertir avant de calculer, saturer, reconvertir.

resultat = np.clip(img.astype(np.int16) + 100, 0, 255).astype(np.uint8)

Lecture, ecriture, affichage

import matplotlib.pyplot as plt
import matplotlib.image as mpimg

img = mpimg.imread("photo.png")        # (h, w, 3) flottants dans [0,1] pour un PNG
plt.imshow(img)
plt.axis("off")
plt.show()

gris = np.zeros((50, 50), dtype=np.uint8)
plt.imshow(gris, cmap="gray", vmin=0, vmax=255)   # vmin/vmax INDISPENSABLES
plt.colorbar()
plt.show()

plt.imsave("sortie.png", gris, cmap="gray")
Toujours preciser vmin et vmax

Sans eux, imshow normalise automatiquement entre le minimum et le maximum de l'image : une image dont les valeurs vont de 100 a 150 apparaitra en plein contraste, ce qui est trompeur. Avec vmin=0, vmax=255, l'affichage est fidele.

Conversion couleur vers niveaux de gris

def vers_gris_moyenne(rvb):
    """Moyenne arithmetique des trois canaux. Simple mais peu fidele."""
    return rvb.mean(axis=2)

def vers_gris(rvb):
    """Luminance perceptuelle (norme ITU-R BT.601)."""
    return 0.299*rvb[:,:,0] + 0.587*rvb[:,:,1] + 0.114*rvb[:,:,2]
Pourquoi ces coefficients ?

L'oeil humain n'est pas egalement sensible aux trois primaires : la retine contient environ 64 % de cones sensibles au rouge, 32 % au vert et 2 % au bleu, mais le vert domine la perception de luminance. Les coefficients \((0{,}299\,;\,0{,}587\,;\,0{,}114)\) refletent cette sensibilite. Une moyenne uniforme rend les zones vertes trop sombres et les zones bleues trop claires.

Transformations ponctuelles

Une transformation ponctuelle applique une fonction \(T\) a chaque pixel, independamment de ses voisins : \(g(i,j) = T\bigl(f(i,j)\bigr)\).

def negatif(img):
    """Inversion : T(v) = 255 - v."""
    return 255 - img

def seuillage(img, s=128):
    """Binarisation : T(v) = 255 si v > s, 0 sinon."""
    return np.where(img > s, 255, 0).astype(np.uint8)

def luminosite(img, k):
    """T(v) = v + k, avec saturation."""
    return np.clip(img.astype(np.int16) + k, 0, 255).astype(np.uint8)

def contraste(img, k):
    """T(v) = 128 + k(v - 128) : etirement autour du gris moyen."""
    return np.clip(128 + k*(img.astype(np.float64) - 128), 0, 255).astype(np.uint8)

def gamma(img, g):
    """Correction gamma : T(v) = 255 (v/255)^g. g < 1 eclaircit, g > 1 assombrit."""
    return (255 * (img/255.0)**g).astype(np.uint8)

def posterise(img, niveaux=4):
    """Reduit le nombre de niveaux de gris."""
    pas = 256 // niveaux
    return ((img // pas) * pas).astype(np.uint8)

Toutes ces fonctions sont vectorisees : aucune boucle. Sur une image \(4000\times3000\), soit 12 millions de pixels, une double boucle Python prendrait plusieurs dizaines de secondes contre quelques millisecondes en NumPy.

Histogramme et egalisation

def histogramme(img):
    """Effectif de chacun des 256 niveaux de gris."""
    return np.bincount(img.ravel(), minlength=256)

def egalise(img):
    """Egalisation d'histogramme : etale les niveaux sur toute la plage."""
    h = histogramme(img)
    cdf = h.cumsum()                          # fonction de repartition cumulee
    cdf_min = cdf[cdf > 0][0]
    table = np.round(255 * (cdf - cdf_min) / (img.size - cdf_min)).astype(np.uint8)
    return table[img]                         # application de la table, vectorisee
Principe de l'egalisation

Si \(F\) est la fonction de repartition empirique des niveaux de gris, alors la variable \(F(X)\) est approximativement uniforme sur \([0,1]\) (theoreme de la transformee inverse). En posant \(T = 255\,F\), on obtient donc une image dont l'histogramme est aussi plat que possible : le contraste est maximise, et l'entropie de l'image augmente.

2Geometrie, convolution et detection de contours

Transformations geometriques

def symetrie_horizontale(img):
    """Miroir gauche-droite."""
    return img[:, ::-1]

def symetrie_verticale(img):
    return img[::-1, :]

def rotation90(img):
    """Quart de tour dans le sens horaire."""
    return img.T[:, ::-1]

def rotation180(img):
    return img[::-1, ::-1]

def recadre(img, i0, j0, h, w):
    return img[i0:i0+h, j0:j0+w]

def sous_echantillonne(img, k):
    """Reduit la taille d'un facteur k (methode naive : on garde un pixel sur k)."""
    return img[::k, ::k]

def moyenne_bloc(img, k):
    """Reduction par moyenne de blocs k x k : bien meilleure qualite."""
    h, w = img.shape
    h, w = (h // k) * k, (w // k) * k
    return img[:h, :w].reshape(h//k, k, w//k, k).mean(axis=(1, 3)).astype(np.uint8)
Sous-echantillonnage et repliement

Garder « un pixel sur \(k\) » (img[::k, ::k]) provoque du repliement de spectre (aliasing) : les motifs fins produisent de faux motifs (effet de moire sur un tissu raye). Il faut d'abord lisser l'image (filtre passe-bas) puis echantillonner — c'est exactement ce que fait moyenne_bloc. Le theoreme de Shannon-Nyquist en donne la justification : il faut echantillonner a plus du double de la frequence maximale presente.

Filtrage par convolution

Definition — Convolution discrete

Le filtrage par un noyau (ou masque) \(K\) de taille \((2r+1)^2\) remplace chaque pixel par une combinaison lineaire de ses voisins :

\[ g(i,j) = \sum_{u=-r}^{r}\sum_{v=-r}^{r} K(u,v)\, f(i+u,\, j+v). \]

C'est une transformation locale : contrairement aux transformations ponctuelles, la valeur de sortie depend du voisinage.

NoyauMatriceEffet
Identite\(\begin{pmatrix}0&0&0\\0&1&0\\0&0&0\end{pmatrix}\) aucun
Moyenne\(\frac19\begin{pmatrix}1&1&1\\1&1&1\\1&1&1\end{pmatrix}\) flou uniforme
Gaussien\(\frac1{16}\begin{pmatrix}1&2&1\\2&4&2\\1&2&1\end{pmatrix}\) flou doux, sans artefacts
Sobel \(x\)\(\begin{pmatrix}-1&0&1\\-2&0&2\\-1&0&1\end{pmatrix}\) gradient horizontal
Sobel \(y\)\(\begin{pmatrix}-1&-2&-1\\0&0&0\\1&2&1\end{pmatrix}\) gradient vertical
Laplacien\(\begin{pmatrix}0&1&0\\1&-4&1\\0&1&0\end{pmatrix}\) contours (derivee seconde)
Nettete\(\begin{pmatrix}0&-1&0\\-1&5&-1\\0&-1&0\end{pmatrix}\) accentuation des details
def convolue(img, noyau):
    """Convolution 2D. Les bords sont laisses inchanges (traitement simple)."""
    img = img.astype(np.float64)
    kh, kw = noyau.shape
    rh, rw = kh // 2, kw // 2
    h, w = img.shape
    res = img.copy()
    for u in range(kh):
        for v in range(kw):
            # on decale l'image entiere : boucle sur le NOYAU (9 tours), pas sur les pixels
            res[rh:h-rh, rw:w-rw] += 0        # initialisation implicite ci-dessous
    # version vectorisee complete :
    sortie = np.zeros((h - 2*rh, w - 2*rw))
    for u in range(kh):
        for v in range(kw):
            sortie += noyau[u, v] * img[u:u+h-2*rh, v:v+w-2*rw]
    res[rh:h-rh, rw:w-rw] = sortie
    return np.clip(res, 0, 255).astype(np.uint8)
L'astuce de vectorisation

On ne boucle pas sur les pixels (des millions), mais sur les cases du noyau (9 pour un \(3\times3\)). A chaque tour, on ajoute une version decalee de l'image entiere, ponderee par le coefficient. Le cout passe de \(\Theta(hw \cdot k^2)\) tours d'interpreteur a seulement \(k^2\) operations vectorielles.

Detection de contours

SOBEL_X = np.array([[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]], dtype=float)
SOBEL_Y = np.array([[-1, -2, -1], [0, 0, 0], [1, 2, 1]], dtype=float)

def gradient(img):
    """Norme et direction du gradient (filtre de Sobel)."""
    f = img.astype(np.float64)
    h, w = f.shape
    gx = np.zeros((h-2, w-2)); gy = np.zeros((h-2, w-2))
    for u in range(3):
        for v in range(3):
            bloc = f[u:u+h-2, v:v+w-2]
            gx += SOBEL_X[u, v] * bloc
            gy += SOBEL_Y[u, v] * bloc
    return np.hypot(gx, gy), np.arctan2(gy, gx)

def contours(img, seuil=100):
    """Carte binaire des contours."""
    norme, _ = gradient(img)
    return (norme > seuil).astype(np.uint8) * 255
Pourquoi Sobel detecte les contours

Un contour est une variation brusque d'intensite, c'est-a-dire un endroit ou le gradient \(\nabla f = (\partial_x f, \partial_y f)\) est grand. Sobel approche ces derivees par differences finies, tout en lissant dans la direction perpendiculaire (les coefficients 1, 2, 1) pour reduire la sensibilite au bruit.

Le noyau se factorise : \(\begin{pmatrix}-1&0&1\\-2&0&2\\-1&0&1\end{pmatrix} = \begin{pmatrix}1\\2\\1\end{pmatrix}\begin{pmatrix}-1&0&1\end{pmatrix}\) : c'est un lissage vertical compose d'une derivation horizontale. Cette separabilite permet de calculer la convolution en \(2k\) operations par pixel au lieu de \(k^2\).

3Filtres interactifs et chaine de traitement

Filtres et transformations sur une image \(8\times8\)
Choisir une image source, puis appliquer un filtre.

Etude de cas complete

import numpy as np
import matplotlib.pyplot as plt

def image_test(n=200):
    """Image synthetique : disque, carre, degrade, plus du bruit."""
    img = np.zeros((n, n), dtype=np.float64)
    i, j = np.indices((n, n))
    img += 60 * j / n                                     # degrade de fond
    img[np.hypot(i - n//3, j - n//3) < n//6] = 200        # disque clair
    img[2*n//3:5*n//6, 2*n//3:5*n//6] = 120               # carre gris
    rng = np.random.default_rng(0)
    img += rng.normal(0, 12, (n, n))                      # bruit gaussien
    return np.clip(img, 0, 255).astype(np.uint8)

def convolue(img, K):
    f = img.astype(np.float64)
    h, w = f.shape
    kh, kw = K.shape
    out = np.zeros((h - kh + 1, w - kw + 1))
    for u in range(kh):
        for v in range(kw):
            out += K[u, v] * f[u:u+h-kh+1, v:v+w-kw+1]
    return out

MOYENNE  = np.ones((3, 3)) / 9
GAUSSIEN = np.array([[1,2,1],[2,4,2],[1,2,1]], dtype=float) / 16
SOBEL_X  = np.array([[-1,0,1],[-2,0,2],[-1,0,1]], dtype=float)
SOBEL_Y  = np.array([[-1,-2,-1],[0,0,0],[1,2,1]], dtype=float)

img = image_test()
lisse = convolue(img, GAUSSIEN)
gx = convolue(lisse.astype(np.uint8), SOBEL_X)
gy = convolue(lisse.astype(np.uint8), SOBEL_Y)
norme = np.hypot(gx, gy)

fig, ax = plt.subplots(2, 3, figsize=(13, 8))
for a, (m, t) in zip(ax.flat, [
        (img, "originale (bruitee)"),
        (255 - img, "negatif"),
        (np.where(img > 128, 255, 0), "seuillage a 128"),
        (lisse, "flou gaussien"),
        (norme, "norme du gradient"),
        (np.where(norme > 60, 255, 0), "contours seuilles")]):
    a.imshow(m, cmap="gray"); a.set_title(t); a.axis("off")
fig.tight_layout(); plt.show()
La chaine de traitement standard
  1. Lisser avant de deriver : le bruit a de fortes hautes frequences, que la derivation amplifie. Sans lissage prealable, la carte de contours est un champ de neige.
  2. Deriver (Sobel) et calculer la norme du gradient.
  3. Seuiller pour obtenir une carte binaire.
  4. (Algorithme de Canny : ajouter une suppression des non-maxima et un seuillage par hysteresis pour obtenir des contours fins et connexes.)

L'essentiel du chapitre

  • Une image est un tableau NumPy : \((h,w)\) en gris, \((h,w,3)\) en couleur ; img[i,j] = ligne \(i\), colonne \(j\), origine en haut a gauche.
  • Le type uint8 deborde silencieusement : convertir, calculer, np.clip, reconvertir.
  • Transformations ponctuelles (negatif, seuillage, gamma) contre locales (convolution) contre geometriques (rotation, symetrie).
  • Vectoriser en bouclant sur le noyau, jamais sur les pixels.
  • Sobel approche le gradient ; toujours lisser avant de deriver.
  • Sous-echantillonner sans lisser provoque du repliement (moire).
  • L'egalisation d'histogramme utilise la fonction de repartition cumulee.
Fiche d'exercices — chapitre 19Enonces gradues, du calcul d'application au niveau concours, avec corriges detailles.
Ouvrir la fiche →