3_3_scikit_image.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from __future__ import division, print_function\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "plt.rcParams['image.cmap'] = 'gray'\n",
    "plt.rcParams['image.interpolation'] = 'none'\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Scikit-Image\n",
    "Dans cette section, nous allons présenter les possibilités offertes par la librairie *scikit-image*. Cette dernière n'étant pas exhaustive, nous vous renvoyons donc sur la page principale de *scikit-image* pour plus d'informations:\n",
    "http://scikit-image.org/\n",
    "\n",
    "Cette section est librement basée sur les tutoriels proposés par `scikit-image`:\n",
    "https://github.com/scikit-image/skimage-tutorials"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Opérations de base\n",
    "### 1.1 Quelques généralités \n",
    "Dans `skimage`, les images sont gérées par des `numpy.ndarray`. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from skimage import data\n",
    "\n",
    "coins = data.coins()\n",
    "\n",
    "print(type(coins), coins.dtype, coins.shape)\n",
    "plt.imshow(coins, cmap='gray', interpolation='nearest');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "De base les images sont représentées par des entiers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "cat = data.chelsea()\n",
    "print(\"Shape:\", cat.shape)\n",
    "print(\"Values min/max:\", cat.min(), cat.max())\n",
    "\n",
    "plt.imshow(cat, interpolation='nearest');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Les fonctions `img_as_float` et `img_as_ubyte` permettent de les convertir en float ou en int."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from skimage import img_as_float, img_as_ubyte\n",
    "\n",
    "image = data.chelsea()\n",
    "\n",
    "image_float = img_as_float(image)\n",
    "image_ubyte = img_as_ubyte(image)\n",
    "\n",
    "print(\"type, min, max:\", image_float.dtype, image_float.min(), image_float.max())\n",
    "print(\"type, min, max:\", image_ubyte.dtype, image_ubyte.min(), image_ubyte.max())\n",
    "\n",
    "print(\"231/255 =\", 231/255.)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour éviter tout problème, skimage conseille de mettre son code sous cette forme :\n",
    "\n",
    "```python\n",
    "def my_function(any_image):\n",
    "   float_image = img_as_float(any_image)\n",
    "   # Proceed, knowing image is in [0, 1]\n",
    "```\n",
    "\n",
    "Il est recommandé d'utiliser la représentation flottante dans la mesure où car c'est celle qui est utilisée par défaut dans `scikit-image`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.2 Accesseurs\n",
    "Comme toute les images sont gérées par `numpy`, il est possible d'accéder aux pixels en utilisant les accesseurs donnés par `numpy`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#\n",
    "plt.imshow(cat);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous pouvons définir une version sous-échantillonnées de l'image de la manière suivante :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "image = data.camera()\n",
    "#\n",
    "plt.imshow(pixelated);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il est d'ailleurs possible d'utiliser `matplotlib` pour afficher les images avec la fonction `imshow`.\n",
    "L'attribut `cmap` permet d'utiliser différentes carte de couleurs pour les afficher."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "image = data.camera()\n",
    "face = image[80:160, 200:280]\n",
    "\n",
    "fig, (ax_jet, ax_viridis, ax_gray) = plt.subplots(ncols=3, figsize=(10, 10))\n",
    "ax_jet.imshow(face, cmap='jet')\n",
    "ax_viridis.imshow(face, cmap='viridis')\n",
    "ax_gray.imshow(face, cmap='gray');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Je vous renvoie vers cette [video](https://www.youtube.com/watch?v=xAoljeRJ3lU) pour une explication de `viridis` et de sa comparaison face à `jet`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercice\n",
    "À partir de l'image `data.camera()`, créer un masque de sorte que les niveaux de gris **inférieur à 87 soient blancs** et que **les pixels en dehors d'un cercle centré et du même diamètre que l'image soient noirs**.\n",
    "\n",
    "Vous aurez besoin de la fonction `np.ogrid` pour créer le cercle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "camera = data.camera()\n",
    "#\n",
    "\n",
    "plt.imshow(camera);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.3 Lecture / Ecriture\n",
    "`Skimage` propose différentes fonctions pour lire/écrire des images : "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "from skimage import io\n",
    "\n",
    "im = #\n",
    "plt.imshow(im);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "io.imsave('exo_img/skimage/input/salta_copy.jpg', im)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`ImageCollection` permet de lire un ensemble d'images."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "pano_imgs = io.ImageCollection('exo_img/skimage/input/pano/JDW_03*')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous introduisons ici une fonction pour comparer plusieurs images entre elles. Cette dernière nous sera utile par la suite de ce notebook."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def compare(*images, **kwargs):\n",
    "    \"\"\"\n",
    "    Utility function to display images side by side.\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    image0, image1, image2, ... : ndarrray\n",
    "        Images to display.\n",
    "    labels : list\n",
    "        Labels for the different images.\n",
    "    \"\"\"\n",
    "    f, axes = plt.subplots(1, len(images), **kwargs)\n",
    "    axes = np.array(axes, ndmin=1)\n",
    "    \n",
    "    labels = kwargs.pop('labels', None)\n",
    "    if labels is None:\n",
    "        labels = [''] * len(images)\n",
    "    \n",
    "    for n, (image, label) in enumerate(zip(images, labels)):\n",
    "        axes[n].imshow(image, interpolation='nearest', cmap='gray')\n",
    "        axes[n].set_title(label)\n",
    "        axes[n].axis('off')\n",
    "    \n",
    "    f.tight_layout()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "compare(#, figsize=(12, 10))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Filtrage\n",
    "Plusieurs filtres simples sont proposés dans le module `skimage.filters`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage import filters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Application d'un filtre TV:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "from skimage.restoration import denoise_tv_bregman\n",
    "\n",
    "denoised = denoise_tv_bregman(pixelated, 4)\n",
    "io.imshow(denoised);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Un filtre de détecteur de contour (Sobel):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "pixelated_gradient = #\n",
    "io.imshow(pixelated_gradient);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ou encore un filtre médian :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from skimage.morphology import disk\n",
    "selem = np.float64(disk(2))\n",
    "\n",
    "median = filters.rank.median(pixelated, selem)\n",
    "io.imshow(median);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Néanmoins ces filtres sont généralement définis uniquement pour des niveaux de gris. \n",
    "\n",
    "La librairie propose le décorateur `@adapt_rgb` pour généraliser son utilisation aux images couleur. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.color.adapt_rgb import adapt_rgb, each_channel, hsv_value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le décorateur se place avant la fonction. L'attribut définit la manière dont doit être appliqué le filtre sur les trois canaux. \n",
    "- `each_channel` applique indépendament le filtre sur chaque canal couleur. \n",
    "- `hsv_value` convertit l'image dans l'espace couleurs HSV et applique le filtre sur chaque canal indépendament."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "@adapt_rgb(#)\n",
    "def sobel_each(image):\n",
    "    return filters.sobel(image)\n",
    "\n",
    "\n",
    "#\n",
    "def sobel_hsv(image):\n",
    "    return filters.sobel(image)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from skimage.exposure import rescale_intensity\n",
    "\n",
    "image = data.astronaut()\n",
    "\n",
    "fig = plt.figure(figsize=(14, 7))\n",
    "ax_each = fig.add_subplot(121, adjustable='box-forced')\n",
    "ax_hsv = fig.add_subplot(122, sharex=ax_each, sharey=ax_each,\n",
    "                         adjustable='box-forced')\n",
    "\n",
    "ax_each.set_title(\"Sobel filter computed\\n on individual RGB channels\")\n",
    "ax_each.imshow(rescale_intensity(1 - sobel_each(image)));\n",
    "\n",
    "ax_hsv.set_title(\"Sobel filter computed\\n on (V)alue converted image (HSV)\")\n",
    "ax_hsv.imshow(rescale_intensity(1 - sobel_hsv(image)));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il est possible de définir sa propre convertion de couleurs en définissant une fonction de convertion. \n",
    "\n",
    "Cette dernière doit être de la forme suivante :\n",
    "``` python\n",
    "def as_personal_convertion(image_filter, image, *args, **kwargs):\n",
    "    # Tout un tas d'opération pour appliquer image_filter a l'image\n",
    "    ...\n",
    "    return filtered_image\n",
    "```\n",
    "\n",
    "Il est également possible de généraliser un traitement à un ensemble d'images en utilisant les possibilités offertes par Pythons."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from skimage.color import rgb2gray\n",
    "pano0, pano1, pano2 = [rgb2gray(im) for im in pano_imgs]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "compare(pano0, pano1, pano2, figsize=(12, 10))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il est tout à fait possible de définir ses propres filtres à partir des `numpy.array`. Considérer l'image de test suivante :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "bright_square = np.zeros((7, 7), dtype=float)\n",
    "bright_square[2:5, 2:5] = 1\n",
    "print(bright_square)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.imshow(bright_square, cmap='gray');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous définissons le filtre moyenneur de la manière suivante :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "mean_kernel = #\n",
    "\n",
    "%precision 2\n",
    "print(mean_kernel)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour l'appliquer, nous allons utiliser la fonction convolve fourni par `scipy.ndimage` :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from scipy.ndimage import convolve\n",
    "convolve?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le résultat obtenue après convolution de `bright square`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%precision 2\n",
    "print(convolve(bright_square, mean_kernel))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le résultat obtenu sur une image:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "new_image = #\n",
    "plt.imshow(new_image);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercice\n",
    "- Définissez un filtre qui calcule le gradient de l'image.\n",
    "- Ajouter un décorateur de sorte que si l'image est couleur, elle soit convertie en niveau de gris avant de calculer le gradient de l'image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.color import rgb2gray\n",
    "\n",
    "def as_gray(image_filter, image, *args, **kwargs):\n",
    "#"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "@adapt_rgb(as_gray)\n",
    "def compute_gradient(im):\n",
    "#"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous vous proposons de vérifier que les résultats sont corrects :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "v_img, h_img = compute_gradient(data.astronaut())\n",
    "\n",
    "plt.figure(figsize=(10, 10))\n",
    "plt.subplot(121)\n",
    "plt.imshow(v_img, cmap='viridis');\n",
    "\n",
    "plt.subplot(122)\n",
    "plt.imshow(h_img, cmap='viridis');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Retouche\n",
    "La librairie fournie également quelques fonctions de retouches des images contenues dans le module `skimage.exposure`. En voici quelques exemples :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage import exposure"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Une fonction pour ajuster le niveau gamma de l'image:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "exposure.exposure.adjust_gamma?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "im = data.moon()\n",
    "gamma_img = exposure.exposure.adjust_gamma(im, 1.)\n",
    "plt.imshow(gamma_img);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Une autre qui permet d'ajuster les variations d'histogrammes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "exposure.equalize_adapthist?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "log_img = exposure.equalize_adapthist(im, clip_limit=0.03)\n",
    "plt.imshow(log_img);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercice\n",
    "Nous vous proposons de programmer la version généralisée du `Midway Equilization Histogram` proposée par Julie Delon. Cet algorithme permet de limiter les variations de luminosité et de contraste d'une même série d'images.\n",
    "\n",
    "Il fonctionne de la manière suivante :\n",
    "- Calculer les histogrammes cumulés des deux images en utilisant `exposure.cumularive_distribution`\n",
    "- Trouver le niveau gris $n_2$ de l'image 2 qui possèdent le même distribution cumulées que le niveau de gris $n_1$  requête de l'image 1.\n",
    "- Définir le nouveaux niveaux gris par $n_m = (n_1 + n_2) / 2$\n",
    "- Mettre à jour $n_1 = n_m, n_2 = n_m$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def apply_midway_general(u_1, u_2):\n",
    "#"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vous pouvez contrôler les résultats avec les lignes suivantes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "im_1 = io.imread('exo_img/skimage/input/fitz_roy_1.png')\n",
    "im_2 = io.imread('exo_img/skimage/input/fitz_roy_2.png')\n",
    "\n",
    "midway_1, midway_2 = apply_midway_general(im_1, im_2)\n",
    "\n",
    "plt.figure(figsize=(17, 10))\n",
    "plt.subplot(221)\n",
    "plt.axis('off')\n",
    "plt.imshow(im_1)\n",
    "plt.title('Image 1 avant correction')\n",
    "\n",
    "plt.subplot(222)\n",
    "plt.axis('off')\n",
    "plt.imshow(im_2)\n",
    "plt.title('Image 2 avant correction')\n",
    "\n",
    "plt.subplot(223)\n",
    "plt.axis('off')\n",
    "plt.imshow(midway_1)\n",
    "plt.title('Image 1 apres correction')\n",
    "\n",
    "plt.subplot(224)\n",
    "plt.axis('off')\n",
    "plt.imshow(midway_2)\n",
    "plt.title('Image 2 apres correction')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Travailler avec des blocs\n",
    "Les blocs permettent de diviser l'image en sous imagettes **sans aucun recouvrement**. Attention ceci est différent de la notion de patches."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.util.shape import view_as_blocks"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "view_as_blocks(#)#"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "l = rgb2gray(data.astronaut())\n",
    "block_shape = (4, 4)\n",
    "view = view_as_blocks(l, block_shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "view.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "flatten_view = view.reshape(view.shape[0], view.shape[1], -1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "flatten_view.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A partir de ses derniers il est néanmoins possible de définir un certain nombre de traitement des images :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#\n",
    "\n",
    "plt.imshow(mean_view);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#\n",
    "\n",
    "plt.imshow(max_view);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#\n",
    "\n",
    "plt.imshow(median_view);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour extraire des patches, nous vous renvoyons vers la fonction `view_as_windows` qui permet d'extraire des sous parties d'une image avec un recouvrement.\n",
    "\n",
    "**Exercice Bonus :** Ecrivez une version simple du filtre non-local à partir de `view-as-windows`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## 5. Extraction de features\n",
    "La libraire offre également un certain nombre de fonction permettant l'extraction et l'appariement des features."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.1 Détection de features\n",
    "Nous allons ici utiliser l'`ORB` (d'autres sont disponibles dans la librairie). Je vous renvoie vers la documentation de `skimage` pour plus d'informations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.feature import ORB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ORB?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "orb = ORB(n_keypoints=800, fast_threshold=0.05)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "orb.detect_and_extract(pano0)\n",
    "keypoints0 = orb.keypoints\n",
    "descriptors0 = orb.descriptors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "orb.detect_and_extract(pano1)\n",
    "keypoints1 = orb.keypoints\n",
    "descriptors1 = orb.descriptors"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.2 Appariement des features\n",
    "La librairie possède également une fonction d'appariement des features :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.feature import match_descriptors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "match_descriptors?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "matches01 = match_descriptors(#)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from skimage.feature import plot_matches\n",
    "fig, ax = plt.subplots(1, 1, figsize=(15, 12))\n",
    "\n",
    "plot_matches(ax, pano0, pano1, keypoints0, keypoints1, matches01)\n",
    "ax.axis('off');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 5.3 Exercice : Appariement robuste\n",
    "Comme nous l'avons vu précédemment, la fonction appariement peut avoir quelques ratés. Par conséquent, nous vous proposant d'utiliser la méthode `RANSAC` pour faire un appariement plus sur."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.measure import ransac"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ransac?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.transform import ProjectiveTransform"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ProjectiveTransform?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Utiliser la méthode `RANSAC` pour obtenir un appariement robuste:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Select keypoints from \n",
    "#   * source (image to be registered): pano0\n",
    "#   * target (reference image): pano1, our middle frame registration target\n",
    "#"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(15, 12))\n",
    "plot_matches(ax, pano0, pano1, keypoints0, keypoints1, matches01[inliers01])\n",
    "ax.axis('off');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Warping\n",
    "Maintenant que certaines features sont mises en correspondances, on va pouvoir corriger la déformation des images et ainsi définir un panorama. Commençons par déterminer la taille de l'image résultante :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.transform import SimilarityTransform\n",
    "\n",
    "# Shape of middle image, our registration target\n",
    "r, c = pano1.shape[:2]\n",
    "\n",
    "# Note that transformations take coordinates in (x, y) format,\n",
    "# not (row, column), in order to be consistent with most literature\n",
    "corners = np.array([[0, 0],\n",
    "                    [0, r],\n",
    "                    [c, 0],\n",
    "                    [c, r]])\n",
    "\n",
    "# Warp the image corners to their new positions\n",
    "warped_corners01 = model_robust01(corners)\n",
    "\n",
    "# Find the extents of both the reference image and the warped\n",
    "# target image\n",
    "all_corners = np.vstack((warped_corners01, corners))\n",
    "\n",
    "# The overall output shape will be max - min\n",
    "corner_min = np.min(all_corners, axis=0)\n",
    "corner_max = np.max(all_corners, axis=0)\n",
    "output_shape = (corner_max - corner_min)\n",
    "\n",
    "# Ensure integer shape with np.ceil and dtype conversion\n",
    "output_shape = np.ceil(output_shape[::-1]).astype(int)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous pouvons appliquer une homomgraphie sur l'une des deux images de sorte à ce qu'elle puisse être associé ensemble :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from skimage.transform import warp\n",
    "\n",
    "# This in-plane offset is the only necessary transformation for the middle image\n",
    "offset1 = SimilarityTransform(translation= -corner_min)\n",
    "\n",
    "# Translate pano1 into place\n",
    "pano1_warped = warp(pano1, offset1.inverse, order=3,\n",
    "                    output_shape=output_shape, cval=-1)\n",
    "\n",
    "# Acquire the image mask for later use\n",
    "pano1_mask = (pano1_warped != -1)  # Mask == 1 inside image\n",
    "pano1_warped[~pano1_mask] = 0      # Return background values to 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le masque précédent est utilie pour ne pas faire une moyenne avec les pixels noirs :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Warp pano0 (left) to pano1\n",
    "transform01 = (model_robust01 + offset1).inverse\n",
    "pano0_warped = warp(pano0, transform01, order=3,\n",
    "                    output_shape=output_shape, cval=-1)\n",
    "\n",
    "pano0_mask = (pano0_warped != -1)  # Mask == 1 inside image\n",
    "pano0_warped[~pano0_mask] = 0      # Return background values to 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "compare(pano0_warped, pano1_warped, figsize=(12, 10));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###  7. Exercice : Panorama\n",
    "Les deux images sont désormais rectifiées, il ne reste plus qu'à les assembler sur une même image :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Add the three images together. This could create dtype overflows!\n",
    "# We know they are are floating point images after warping, so it's OK.\n",
    "#\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Track the overlap by adding the masks together\n",
    "#\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Normalize through division by `overlap` - but ensure the minimum is 1\n",
    "#\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le résultat précédent n'est pas parfait mais permet d'avoir une première approche du problème du panorama. Pour un résultat optimal, nous devrions appliquer un algorithme qui calcule automatiquement les bons masques à appliquer sur les deux images de sorte à ce que les divisions ne soient plus apparentes. \n",
    "\n",
    "Si tout ceci vous intéresse, nous vous renvoyons vers le tutoriel panorama de `Skimage` à `Scipy 2015` : \n",
    "https://github.com/scikit-image/skimage-tutorials"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}