Traitement d’images

Arcimboldo painting

Giuseppe Arcimboldo était un peintre du 16ème siècle qui a créé des portraits imaginaires à partir de fruits, légumes, fleurs, poissons et livres.

Objectif

Dans cette feuille, vous allez créer une séquence d’images transformant pas à pas le smiley ci-dessous en une peinture à la Arcimboldo, puis les assembler pour réaliser une petite animation! Vous pourrez ensuite épicer votre réalisation avec Deep Art.

Smiley

Voici les étapes du projet:

  1. Extraire la partie droite de l’image pour avoir une image carrée

  2. Extraire le smiley de son fond

  3. Recadrer l’image sans changement de résolution, de sorte à obtenir une image de taille 256x256 centrée sur le smiley,

  4. Réduire la résolution à 128x128, avec une technique anti-crénelage, par lissage puis sous-échantillonage.

  5. Extraire les contours du smiley pour en faire une image noir et blanc.

  6. Mettre en toile de fond un papier-peint représentant des fruits

  7. Choisir trois images de fruits et les extraires de leur fond.

  8. Remplacer les yeux et la bouche par ces fruits

Lors de ces étapes, guidées pas-à-pas, vous apprendrez à:

  • Appliquer un filtre pixel-à-pixel pour extraire le premier plan d’une image

  • Appliquer un filtre de convolution, construire une pyramide Gaussienne et choisir le bon niveau de lissage avant de sous-échantillonner.

  • Appliquer un filtre de convolution par différence pour extraire un contour.

# Automatically reload code when changes are made
%load_ext autoreload
%autoreload 2
from PIL import Image, ImageDraw, ImageFont
import matplotlib.pyplot as plt
from scipy import signal

from intro_science_donnees import *
from utilities import *

Lecture de l’image et extraction d’une sous image carrée

Exercice: Chargez l’image smiley.jpg situé dans le dossier media et affectez-la à la variable img.

Indication: consultez la feuille image de la Semaine4 si besoin.

# VOTRE CODE ICI
raise NotImplementedError()
img
assert img.width == 640
assert img.height == 427

La bibliothèque PIL fournit une méthode crop pour recadrer des images. Elle prend en paramètre une liste (left, upper, right, lower) donnant successivement les coordonnées des points en haut à gauche et en bas à droite de la zone rectangulaire (box) à extraire. Voici, par exemple, l’image recadrée pour les points de coordonnées (0,0) et (200,100):

img.crop(box=(0,0,200,100) )

Exercice: Recadrer l’image pour ne garder que sa partie droite et de sorte à ce qu’elle soit carrée. Affectez le résultat à la variable img_carree.

Indication : Faites un dessin sur papier. Quel coin garder ? Comment trouver les autres coordonnées ? Pour cela, commencez par identifier la plus petite dimension entre la largeur (width) et la hauteur (height).

# VOTRE CODE ICI
raise NotImplementedError()
img_carree
# Vérifications: on obtient une image de la bonne taille
assert isinstance(img_carree, Image.Image)
assert img_carree.size == (427, 427)

Ce sera la première image de votre animation!

Pour produire cette dernière, vous stockerez les images successives dans une liste images, puis les assemblerez à la fin. Commençons par celle-ci:

# Initialise une liste d'images vide
images = []
# Rajoute une image à la fin
images.append(img_carree)

La méthode append permet d’ajouter des élements à la fin d’une liste. Nous la réutiliserons pour mettre à jour images au fur et à mesure du TP.

image_grid(images)

Le générique de l’animation

Lorsque l’on produit une œuvre, ici une animation, l’attribution de l’œuvre à ses auteurs est un élément indispensable. Pour cela, vous allez maintenant écrire votre nom et prénom sur l’image.

Exercice :

  1. Affectez votre nom à la variable name.

# VOTRE CODE ICI
raise NotImplementedError()
assert isinstance(name, str)
  1. Étudiez en détail les quatre commandes suivantes qui respectivement:

    • font une copie de l’image

    • construisent un canevas pour dessiner sur cette image

    • choisissent une fonte (si “LiberationSans-Regular” n’est pas disponible, essayez “Arial” à la place)

    • dessinent un texte

img_title = img_carree.copy()
canvas = ImageDraw.Draw(img_title)
font = ImageFont.truetype("LiberationSans-Regular.ttf", size=40)
canvas.text(xy=(100,10), text=name, font=font)
  1. Affichez l’image

# VOTRE CODE ICI
raise NotImplementedError()
  1. Quel argument permet de définir où le texte sera écrit ?

VOTRE RÉPONSE ICI

  1. Ajoutez l’image produite à la liste images

# VOTRE CODE ICI
raise NotImplementedError()
image_grid(images)
assert len(images) >= 2
assert images[-1] is img_title
  1. ♣ Modifiez ce qui précède pour personnaliser votre générique. Voici quelques idées:

    • explorez les options de canvas.text pour dessiner le texte selon vos préférences (couleur, fonte, position, …);

    • ajoutez non pas une seule image à images, mais une séquence d’images de sorte à réaliser un générique qui défile;

    • explorez les fonctions de canvas pour réaliser tout dessin de votre goût sur l’image.

    Indications:

    • rajoutez autant de cellules ci-dessus que vous en aurez besoin;

    • rassemblez les commandes pour rajouter du texte à une image en une fonction réutilisable dans utilities.py

Détection de l’avant-plan

L’objectif est maintenant de déterminer la position du smiley sur l’image. Pour cela nous allons commencer par détecter quels pixels sont dans l’avant-plan (le smiley), et non dans l’arrière plan.

Un premier essai

Reprenons la fonction foreground_filter du Projet1, qui est fournie dans utilities.py. Prenez le temps d’en consulter le code pour vous remémorer son fonctionnement.

show_source(foreground_filter)

Exercice : Appliquez cette fonction à img_carree, affectez le résultat à foreground0, et visualisez le résultat:

# VOTRE CODE ICI
raise NotImplementedError()
plt.imshow(foreground0, cmap='gray');
# Vérifications:
assert isinstance(foreground0, np.ndarray)    # foreground est un tableau ...
assert foreground0.dtype == np.dtype('bool')  # de booléens ...
assert foreground0.shape == img_carree.size   # de même taille que l'image

Exercice :

  1. Le résultat est-il probant? Réessayez avec différentes valeurs de seuil theta. Laquelle donne le meilleur résultat?

    VOTRE RÉPONSE ICI

  2. Est-ce satisfaisant? Quel est le problème?

    VOTRE RÉPONSE ICI

Le smiley a la jaunisse

Pour améliorer cela, nous allons utiliser un filtre pixel-à-pixel mesurant non pas le niveau de gris mais le niveau de jaune de chaque pixel.

Exercice:

  1. Implémentez la fonction yellowness_filter dans le fichier utilities.py en utilisant la formule \(y = r + g - b\)\(y\) est l’intensité du jaune (yellow).
    Indication: si vous n’avez plus en tête comment extraire et calculer avec des différentes couches de couleur de l’image, consulter la feuille d’extraction d’attributs.

show_source(yellowness_filter)
  1. Appliquez ce filtre à img_carree et gardez le résultat dans l’objet img_yellowness

# VOTRE CODE ICI
raise NotImplementedError()

Affichons le résultat:

plt.imshow(img_yellowness, cmap='gray')
plt.colorbar();
# Vérifications: 
assert isinstance(img_yellowness, np.ndarray)      # img_yellowness est un tableau numpy ...
assert img_yellowness.dtype == np.dtype('float64') # de nombre flottants ...
assert img_yellowness.shape == img_carree.size     # de même taille que l'image
  1. Quelle est la valeur maximale théorique que peut prendre un pixel d” img_yellowness ? Et minimale ?

    VOTRE RÉPONSE ICI

Ajoutez img_yellowness à la liste d’images:

# VOTRE CODE ICI
raise NotImplementedError()
image_grid(images)

Généralisons!

Notez que la formule \(r+g-b\) est un produit scalaire! Nommément, c’est le produit scalaire \(v.w\) de la couleur du pixel – représentée par le vecteur \(v = (r,g,b)\) – avec le vecteur \(w=(1,1,-1)\). Le rôle du coefficient \(-1\) est de discriminer les couleurs du blanc – représentée par le vecteur \((255,255,255)\) – qui contiennent du bleu en plus du jaune.

Une autre façon de discriminer ces couleurs, est de diviser le produit scalaire par la norme de \(v\). Cela nous donne alors \(\frac{v . w}{|v|}\), où \(w=(1,1,0)\). À un coefficient multiplicateur et centrage près, on retrouve la formule \(\frac{\overline v.\overline c}{|\overline v|,|\overline c|}\) de la corrélation de \(v\) avec la couleur jaune \(c=(255,255,0)\). La corrélation est maximale lorsque \(v\) est proportionnel à \(c\), c’est-à-dire que \(v\) représente une couleur jaune quelconque, qu’elle soit sombre ou claire.

Exercice ♣:

  1. Implémenter dans utilities.py la fonction color_correlation_filter calculant, pour chaque pixel d’une image, sa corrélation avec une couleur donnée.

    Indication: Transformez l’image en array. Puis pour chaque pixel, calculez la corrélation avec la couleur à l’aide de la fonction np.corrcoef(). Utilisez une compréhension. Usage : la fonction de Numpy np.corrcoef(u, v)[0, 1] calcule la corrélation entre deux vecteurs u et v.

show_source(color_correlation_filter)
  1. Appliquer cette fonction à img_carree avec la couleur jaune \((255,255,0)\) et affecter le résultat à img_yellowness.

yellow = (255,255,0)
# VOTRE CODE ICI
raise NotImplementedError()
  1. Affichez le résultat:

# VOTRE CODE ICI
raise NotImplementedError()
  1. Le smiley n’est en fait pas en jaune pur (la corrélation n’est pas exactement de 1). Extraire une zone du smiley et en calculer la couleur moyenne. Puis utiliser cette couleur pour appliquer le filtre.

# VOTRE CODE ICI
raise NotImplementedError()
plt.imshow(img_yellowness, cmap="gray")
plt.colorbar();
  1. Quels sont les avantages et inconvénients de cette méthode?

    VOTRE RÉPONSE ICI

  1. Ajoutez img_yellowness à la liste d’images:

# VOTRE CODE ICI
raise NotImplementedError()
image_grid(images)

Avant-plan par jaunisse

Exercice:

  1. Utilisez un seuil sur img_yellowness pour calculer les pixels du smiley; affectez le résultat à la variable foreground.

# VOTRE CODE ICI
raise NotImplementedError()
plt.imshow(foreground, cmap='gray');
# Vérifications:
assert isinstance(foreground, np.ndarray)    # foreground est un tableau ...
assert foreground.dtype == np.dtype('bool')  # de booléens ...
assert foreground.shape == img_carree.size   # de même taille que l'image
  1. Ajoutez foreground à la liste d’images:

# VOTRE CODE ICI
raise NotImplementedError()
image_grid(images)

♣ En boîte

Exercice:

  1. Calculez les coordonnées de la plus petite boîte rectangulaire (bounding box) englobant le smiley. Pour cela, il faut ignorer les points épars. On pourra par exemple fixer un seuil \(s\) et ne considérer que les lignes et colonnes contenant au moins \(s\) pixels.

  2. Dessinez la boite sur l’image d’origine, et ajoutez le résultat à la liste d’images. Pour cela utilisez la classe ImageDraw.Draw de PIL et sa méthode .rectangle()

TODO: rajouter solution

Recadrage (crop)

Calcul du centre du smiley

Exercice: Calculez les coordonnées du centre du smiley et affectez le résultat à la variable center.

Indications:

  • Vous pouvez par exemple calculer le barycentre des coordonnées \((i,j)\) des pixels dans le smiley.
    Consultez le code de la fonction elongation de la semaine 2.

  • Si vous l’avez calculée ci-dessus, vous pouvez aussi vous baser sur la boîte englobante.

# VOTRE CODE ICI
raise NotImplementedError()
assert np.linalg.norm(center - np.array([275, 230])) < 5

Exercice (\(\clubsuit\)):

  1. Dessinez le centre sur l’image a l’aide d”ImageDraw.Draw() et de sa méthode .ellipse().

  2. Affichez le résultat et ajoutez le résultat à la liste d’images.

Indication: Pour dessiner sur l’image, il faut au préalable la reconvertir depuis notre tableau de booléen vers une image couleur. C’est ce que fait la première ligne ci-dessous.

img_barycentre = Image.fromarray(foreground).convert("RGB")

# VOTRE CODE ICI
raise NotImplementedError()
assert isinstance(img_barycentre, Image.Image)
assert img_barycentre.size == img_carree.size
assert list(np.array(img_barycentre)[int(center[0]), int(center[1]), :]) == [255,0,0], "Le centre du smiley n'est pas rouge"

Recadrage du smiley

Exercice:

  1. Recadrez votre image pour qu’elle soit centrée sur le smiley et de taille 256 x 256.

# VOTRE CODE ICI
raise NotImplementedError()
plt.imshow(img_cropped);
  1. Procédez de même pour foreground en affectant le résultat à foreground_cropped.

    Indication: Utilisez Image.fromarray pour reconvertir foreground en Image.

# VOTRE CODE ICI
raise NotImplementedError()
plt.imshow(foreground_cropped);
  1. Ajoutez foreground_cropped à la liste d’image! Vous pouvez aussi ajouter l’image recadrée maintenant, ou bien plus tard selon votre goût.

# VOTRE CODE ICI
raise NotImplementedError()
image_grid(images)

Réduction de la résolution

Nous allons voir comment réduire la résolution d’une image. Cela peut être typiquement utile pour réduire le volume des données et accélérer les calculs.

Essayons de réduire la résolution de notre image à 128x128 pixels.

img_downsampled = foreground_cropped.resize((128, 128))
plt.imshow(img_downsampled, interpolation="none"); # imshow applique un anti-crénelage par défaut

Comme on le voit ci-dessus, l’opération de sous-échantillonnage permettant de réduire la résolution d’une image est sujette au problème de crénelage (aliasing en anglais).

Pour mitiger cela, il existe plusieurs techniques d”anti-crénelage; elles consistent à lisser l’image obtenue, en donnant à chaque pixel de l’image obtenue une valeur interpolée à partir des pixels qu’il «représente» dans l’image d’origine.

Cependant, pour interpoler entre du noir et du blanc, il faut avoir du gris à disposition! C’est pourquoi la méthode resize n’a pas pu la lisser notre image foreground_cropped qui est en noir et blanc.

Dans tous les autres cas, resize mets automatiquement en œuvre de l’anti-crénelage. Il suffit donc de convertir au préalable notre image en niveaux de gris pour en bénéficier:

# VOTRE CODE ICI
raise NotImplementedError()

Consultez la documentation de convert! Cette méthode vous resservira.

Suppression du bruit par lissage

À ce stade, notre avant-plan contient de nombreux points isolés qui viennent du bruit dans l’image (des pixels du fond qui sont par hasard de la même couleur que l’objet, et réciproquement):

plt.imshow(foreground_cropped, cmap='gray');

Pour réduire le bruit, une technique classique est de lisser (ou flouter) l’image, en remplaçant chaque pixel par une moyenne avec ses pixels voisins. Ainsi, des pixels blancs isolés au milieu de pixels noirs deviendront gris sombre. Réciproquement pour des pixels noirs isolés au milieu de pixels blancs deviendront gris clair. Il restera alors à appliquer un nouveau seuillage pour les recataloguer en pixels blancs ou noirs.

Lissage par filtre de Gauß

Une première façon de lisser une image est d’utiliser un filtre de Gauss, comme celui implanté dans la bibliothèque SciPy.

sigma=1

from scipy.ndimage import gaussian_filter
img_filtered = gaussian_filter(foreground_cropped.convert("L"), sigma=sigma)
plt.imshow(img_filtered, cmap='gray');

Après seuillage, on obtient:

seuil = 100

foreground_cropped_clean = Image.fromarray(img_filtered > seuil)
plt.imshow(foreground_cropped_clean, cmap='gray');

Exercice

  1. Faites varier ci-dessous la valeur du paramètre sigma et celle du seuil et observez le résultat. Que se passe-t’il si sigma est trop faible ou trop élevé?

    VOTRE RÉPONSE ICI

  1. Proposez ci-dessous des valeurs pour extraire au mieux l’avant-plan.

# VOTRE CODE ICI
raise NotImplementedError()

img_filtered = gaussian_filter(foreground_cropped.convert("L"), sigma=sigma)
foreground_cropped_clean = Image.fromarray(img_filtered > seuil)

fig = plt.figure(figsize=(15,6))
ax = fig.add_subplot(1,3,1)
ax.set_title("avant-plan original")
ax.imshow(foreground_cropped, cmap="gray")
ax = fig.add_subplot(1,3,2)
ax.imshow(img_filtered, cmap="gray")
ax.set_title("lissé")
ax = fig.add_subplot(1,3,3)
ax.imshow(foreground_cropped_clean, cmap="gray")
ax.set_title("avant-plan nettoyé");
  1. ♣ Comme beaucoup de traitements d’images, lisser est un traitement local: un pixel est transformé en fonction de sa valeur et des valeurs de ses voisins plus ou moins proche. Cela peut s’exprimer par un produit de convolution. Dans le filtre de Gauss, on convolue avec une fonction Gaussienne (courbe en cloche).

    À quelle propriété de la fonction Gaussienne correspond le paramètre sigma (noté 𝜎). Expliquer avec vos mots pourquoi augmenter 𝜎 floute davantage l’image.

    VOTRE RÉPONSE ICI

♣ Lissage par pyramide Gaussienne

Dans cette section, nous allons explorer un autre procédé de lissage. Au lieu d’effectuer une seule convolution avec une fonction compliquée, nous allons itérer plusieurs fois une convolution très simple.

À chaque itération chaque pixel de la nouvelle image est obtenu à partir du pixel d’origine et de ses huit voisins immédiats en y appliquant des coefficients. Ces coefficients sont donnés par une matrice \(3\times3\) appelée le noyau (kernel en anglais). Nous utiliserons le noyau suivant. Le \(.25\) au centre spécifie par exemple que le pixel d’origine contribuera avec un coefficient de \(1/4\), tandis que celui du dessus contribuera avec un coefficient de \(1/8\):

ker = np.outer([1, 2, 1], [1, 2, 1])
ker = ker/np.sum(ker)
ker

Nous appliquons seize fois de suite la convolution à l’image foreground_cropped en affichant les résultats intermédiaire pour visualiser l’évolution. La convolution est calculée en utilisant signal.convolve2d; les options précisent le comportement au bord.

import scipy.signal

M = foreground_cropped

fig = plt.figure(figsize=(20,20))
for k in range(16):
    M = scipy.signal.convolve2d(M, ker, boundary='wrap', mode='same')

    fig.add_subplot(4,4,k+1)
    plt.imshow(M, cmap='gray')

On note que, au fur et à mesure, les pixels blancs isolés deviennent de plus en plus sombres: ils sont en effet influencés par leurs voisins noirs, puis par les voisins des voisins et ainsi de suite. On note aussi que l’image devient de plus en plus floue; il ne faudrait donc pas aller trop loin.

Il ne reste plus qu’à appliquer à nouveau un seuil, pour obtenir une image propre:

foreground_cropped_clean = M > 0.7
plt.imshow(foreground_cropped_clean, cmap='gray');
smiley = transparent_background_filter(img_cropped, foreground_cropped_clean)
smiley

Rajoutez les images que vous souhaitez à votre liste d’images

# VOTRE CODE ICI
raise NotImplementedError()
image_grid(images)

Exercice : Rappelez ci-dessous les étapes effectuées jusque là :

VOTRE RÉPONSE ICI

Extraction des contours

Nous souhaitons maintenant extraire les contours de l’image. Remettons l’image en réels:

M = foreground_cropped_clean * 1.0

Le principe est de calculer la valeur absolue de la différence entre chaque pixel et son voisin de gauche. S’ils sont identiques, on obtient zéro. S’ils sont différents – on est sur un bord – on obtient une valeur positive. De manière équivalente, on calcule la différence entre l’image et l’image décallée de un pixel.

Note: lors de cette opération de différence, on perd une bande de 1 pixel en largeur; pour conserver l’image carrée, on a aussi supprimé une bande de 1 pixel en hauteur.

contour_horizontal = np.abs(M[:-1, 1:] - M[:-1, 0:-1])
plt.imshow(contour_horizontal, cmap='gray');

Exercice :

  1. Notez que c’est déjà presque parfait, sauf lorsque le contour est horizontal. Pourquoi ?

    VOTRE RÉPONSE ICI

  1. Pour améliorer cela, procédez de même verticalement et affectez le résultat dans contour_vertical:

# VOTRE CODE ICI
raise NotImplementedError()
plt.imshow(contour_vertical, cmap='gray');

Maintenant, c’est au tour des contours verticaux d’être peu détectés. Qu’à cela ne tienne, il suffit d’additionner les deux résultats:

contour = contour_horizontal + contour_vertical

Et voilà le travail!

plt.imshow(contour, cmap='gray');

Superposition d’images

Nous allons maintenant vous montrer comment superposer deux images. Récupérons nos images de fruits de la semaine 3:

from intro_science_donnees import data
fruit_dir = os.path.join(data.dir, 'ApplesAndBananasSimple')
fruits = load_images(fruit_dir, "*.png")
image_grid(fruits)

Et choisissons une banane:

banana = fruits[18]

On convertit les deux images en tableaux NumPy :

M = np.array(smiley)
B = np.array(banana)

On choisit les coordonnés où l’on veut superposer la banane :

i = 80
j = 130

On extrait dans P la zone de l’image qui contiendra l’image, on calcule l’avant-plan de la banane et, pour tous ces pixels de l’avant-plan, on affecte les couleurs de F à P :

P = M[i:i+32, j:j+32]
F = foreground_filter(banana)
P[F] = B[F]

Enfin on réinsère P dans l’image d’origine :

M[i:i+32, j:j+32] = P

Et voilà :

plt.imshow(M);
images.append(M)

Exercice : Superposer de même un autre fruit là où vous le souhaitez!

# VOTRE CODE ICI
raise NotImplementedError()

Exercice : Rajoutez l’image à la liste.

# VOTRE CODE ICI
raise NotImplementedError()

VOTRE RÉPONSE ICI

VOTRE RÉPONSE ICI

image_grid(images)
assert len(images) >= 7

Production de l’animation

Nous allons maintenant assembler toutes les images pour produire l’animation illustrant toutes les étapes du traitement d’images ci-dessus.

La cellule suivante utilise les widgets de Jupyter pour construire une mini application interactive:

import ipywidgets

# Vue
output = ipywidgets.Output()
play = ipywidgets.Play(min=0, max=len(images)-1, value=0, interval=500)
slider = ipywidgets.IntSlider(min=0, max=len(images)-1, value=0)
controls = ipywidgets.HBox([play, slider])
application = ipywidgets.VBox([controls, output])

# Controleur
def update(event):
    with output:
        output.clear_output(wait=True)
        fig = Figure()
        ax = fig.add_subplot()
        ax.imshow(images[slider.value], cmap='gray')
        display(fig)

slider.observe(update)
output.on_displayed(update)
ipywidgets.jslink((play, "value"), (slider, "value"))

application

Pour pouvoir exporter l’animation, nous allons produire une vidéo à l’aide de la bibliothèque OpenCV que l’on peut appeler directement depuis Python. Consultez le code de l’utilitaire video fourni pour en savoir plus. Un autre outil classique est ffmpeg.

video(images, "video.mp4")

Selon votre navigateur, vous pourrez visualiser la vidéo directement ici:

from IPython.display import Video
Video('video.mp4')

Sinon, téléchargez la.

Conclusion

Vous êtes arrivé à la fin de ce TP où, vous avez appris à créer des diaporamas avec Jupyter et à effectuer de petits traitements d’images tels que ceux que vous mettrez en œuvre pour prétraiter vos données lors de votre projet final.

Il ne vous reste plus qu’à revenir à la feuille d’accueil mettre à jour votre rapport de TP et vérifier la qualité de votre code.