{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "d8ac8f62",
   "metadata": {},
   "source": [
    "# Traitement du Signal - TDM3 : Corrélation de signaux\n",
    "\n",
    "---\n",
    "\n",
    "Aujourd'hui, nous allons implémenter la fonction de corrélation entre deux signaux, la tester sur plusieurs signaux, analyser les résultats en fonction des entrées, et créer notre propre Shazam. Sympa comme programme je trouve, donc on y va dès maintenant !"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "215668ef",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1608e1e0",
   "metadata": {},
   "source": [
    "## Exercice 1 : Commençons par comprendre la corrélation\n",
    "\n",
    "\n",
    "L'intercorrélation compare un signal $f(t)$ et un signal $g(t)$ retardé. Elle mesure la similitude entre ces deux signaux. La formule de l'intercorrélation des signaux $f$ et $g$ est la suivante :\n",
    "\n",
    "\\begin{equation*}\n",
    "C_{fg}(\\tau) = \\int_{-\\infty}^{\\infty} f(t) g^{*}(t - \\tau) \\, dt\n",
    "\\end{equation*}\n",
    "\n",
    "Avec $g^{*}$ le conjugué complexe de $g$.\n",
    "\n",
    "Dans notre cas, nous avons des signaux discrets (non continus, il y'a un échantillonage), donc le calcul de l'intercorrélation devient :\n",
    "\n",
    "\\begin{equation*}\n",
    "C[k] = \\sum_n f[n+k] \\cdot g^*[n] \\cdot T_e\n",
    "\\end{equation*}\n",
    "\n",
    "$T_e$ représente ici le temps d'échantillonage.\n",
    "\n",
    "Avant de partir sur le développement de la fonction, on va d'abord tester avec NumPy. En effet, la librairie propose sa fonction *np.correlate*. Super !\n",
    "\n",
    "*\"Mais alors, pourquoi on la code la fonction du coup ?!\"*\n",
    "\n",
    "On va la développer ensuite pour 2 raisons :\n",
    "- Parce que ça vous permettra de mieux comprendre le fonctionnement\n",
    "- Parce qu'on va la customiser pour créer une variante de la corrélation, disponible nul part..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c46c812a",
   "metadata": {},
   "source": [
    "Bref, pour calculer la corrélation entre deux signaux... et bah, il nous faut 2 signaux. C'est pas plus compliqué... Créez donc :\n",
    "\n",
    "- $A$, un signal **périodique** de période T = 2, dont la période $a$ est définie par l'équation suivante : \n",
    "\n",
    "\\begin{equation*}\n",
    "a(t) = \n",
    "\\begin{cases}\n",
    "1 & \\text{si } |t| \\gt \\frac{T}{2}, \\\\\n",
    "0, & \\text{sinon}.\n",
    "\\end{cases}\n",
    "\\end{equation*}\n",
    "\n",
    "- $B$, un signal défini par l'équation suivante :\n",
    "\n",
    "\\begin{equation*}\n",
    "B(t) = \n",
    "\\begin{cases} \n",
    "1, & \\text{si } 0 \\leq t \\leq 1, \\\\\n",
    "0, & \\text{sinon}.\n",
    "\\end{cases}\n",
    "\\end{equation*}\n",
    "\n",
    "Les deux signaux sont échantillonés à 0.01 secondes. Le signal A s'étend de 0 à 10 secondes, tandis que le signal B s'étend de 0 à 2 secondes.\n",
    "\n",
    "*Note : La construction d'un signal périodique n'étant pas toujours simple, partez du principe que $A(t)$ est égale à 1 lorsque la partie entière de $t$ est impaire, et est sinon à 0*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "60679ffa",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Création des deux signaux\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "26337721",
   "metadata": {},
   "source": [
    "Calculez la corrélation entre les deux signaux via la fonction *np.correlate* (avec l'option *'full'*), et affichez le résultat sans axe temporel (on va gérer ça ensuite)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4a05a203",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Calcul et affichage de la corrélation entre A et B via la fonction correlate de NumPy\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c6f4d3e3",
   "metadata": {},
   "source": [
    "Bon, pour l'instant, le signal ne ressemble pas vraiment à grand chose pour plusieurs raisons :\n",
    "- l'axe temporel n'est pas correct (il n'y en a pas)\n",
    "- l'amplitude n'est pas très pertinente ici (on obtient des pics à 100 alors que les signaux originaux varient entre 0 et 1)\n",
    "\n",
    "On va donc corriger cela !"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b1629f7f",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "**Petit aparté sur les signaux numériques**\n",
    "\n",
    "Il faut bien faire la distinction ici entre un signal continu (ce que vous manipulez en TD), et un signal discrétisé. En numérique (donc sur ordi), les signaux sont obligatoirement discrétisés, c'est à dire:\n",
    "- Stockés via des tableaux\n",
    "- Echantillonnés avec un temps d'échantillonnage ($T_e$) ou une fréquence d'échantillonnage ($f_e$)\n",
    "- Inscrits dans un intervalle de temps\n",
    "\n",
    "De ce fait, lorsqu'on décrit un signal temporel numériquement, on a:\n",
    "- Un tableau représentant l'axe temporel, construit via l'intervalle de temps, et le temps/la fréquence d'échantillonnage\n",
    "- Un tableau représentant le signal associé à l'axe temporel\n",
    "\n",
    "Cela veut dire que savoir faire une corrélation, ce n'est pas juste appeler NumPy et afficher le résultat. Il faut 2-3 manipulations supplémentaires pour que le résultat obtenu ne soit pas juste un tableau NumPy, mais un signal temporel. Ici, il nous manque principalement l'axe temporel, qu'il va falloir retrouver.\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "443fc2ed",
   "metadata": {},
   "source": [
    "Concentrons nous d'abord sur l'axe temporel. On considère que les signaux $A$ et $B$ sont respectivement de taille $N$ et $M$\n",
    "\n",
    "**_QUESTION :_** Quelle est la taille du signal corrélé en sortie ? A quoi cela correspond par rapport aux signaux $A$ et $B$ ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7dd2a7d0",
   "metadata": {},
   "source": [
    "**_REPONSE :_** "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "142cbf9a",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Quel est le temps d'échantillonnage du signal corrélé par rapport au signaux $A$ et $B$ ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8e7d93bd",
   "metadata": {},
   "source": [
    "**_REPONSE :_**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7c6fe771",
   "metadata": {},
   "source": [
    "On a donc déterminé le temps d'échantillonnage, mais il va falloir retrouver l'intervalle de temps. Et pour déterminer, il faut revoir la nature même de la corrélation. La formule nous donne un résultat lorsque **les deux signaux se recoupent temporellement**. Si on reprend la formule suivante dans notre contexte avec $A$ et $B$:\n",
    "\n",
    "\\begin{equation*}\n",
    "C[k] = \\sum_n A[n] \\cdot B^*[n - k] \\cdot T_e\n",
    "\\end{equation*}\n",
    "\n",
    "Dans cette formule, le signal $B$ est retardé de $k$ dans le calcul, contrairement au signal $A$. La somme sur $n$ se fait en réalité uniquement sur les zones où les deux signaux se recoupent temporellement.\n",
    "\n",
    "On sait que $A$ évolue de $0$ à $N-1$, et $B$ de $0$ à $M-1$. Si $B$ a un retard de $k$ cases, cela signifie qu'il évolue entre $k$ et $M+k-1$. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "643335a4",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Quel est le plus grand retard $k$ que $B$ peut prendre, tout en s'assurant que les signaux $A$ et $B$ se recoupent temporellement ? \n",
    "\n",
    "*Note : N'hésitez pas à faire un schéma, avec le tableau $A$ fixe, et le tableau $B$ qui subit un retard $k$, ou une avance $-k$*"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e0ecf986",
   "metadata": {},
   "source": [
    "**_REPONSE :_**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c145b1d9",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Quelle est maintenant la plus grande avance (ou le plus petit retard) $k$ que $B$ peut prendre, toujours en s'assurant que les signaux $A$ et $B$ se recoupent temporellement ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18891785",
   "metadata": {},
   "source": [
    "**_REPONSE :_**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7e8c36c1",
   "metadata": {},
   "source": [
    "Vous avez donc normalement l'intervalle temporel où évolue votre signal corrélé, enfin juste les indices. Pour avoir les vraies valeurs temporelles, il suffit de multiplier par le temps d'échantillonnage. Construisez l'axe temporel de votre signal corrélé"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b8bc7d0a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Construction de l'axe temporel du signal corrélé\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3801e02a",
   "metadata": {},
   "source": [
    "Affichez maintenant le signal obtenu après corrélation, avec l'axe temporel correct."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fce3758a",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Affichage du signal corrélé avec le bon axe temporel\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1060c452",
   "metadata": {},
   "source": [
    "Normalement, vous devriez voir le signal évoluer entre -2 et 10 secondes. Si c'est le cas, on a un bon axe temporel. Mais l'amplitude n'est toujours pas bon.\n",
    "\n",
    "Afin d'illustrer le problème d'amplitude, refaites les mêmes opérations que précédemment, mais avec un temps d'échantillonnage de 0.0001 secondes, et trouvez la différence obtenue.\n",
    "\n",
    "*Pas besoin de recopier le code, modifiez juste Te au début du notebook et relancez les cellules*"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5777b499",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Que constatez-vous ? Comment $T_e$ influence sur le résultat de la corrélation ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cce33478",
   "metadata": {},
   "source": [
    "**_REPONSE :_** "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8663952f",
   "metadata": {},
   "source": [
    "En effet ici, le souci est lié avec le temps d'échantillonnage $T_e$, car le calcul réalisé par NumPy est une simple somme. Hors, on essaye de s'approcher de la formule originale qui est une intégrale. Il faut donc multiplier tout cela par le temps d'échantillonnage."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5519eba7",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "**RECAPITULATIF DE LA CORRELATION AVEC NUMPY**\n",
    "\n",
    "Normalement, vous avez ici **tout ce qu'il vous faut** pour effectuer une corrélation correcte entre deux signaux. A partir des signaux A et B (et leurs axes temporels), écrivez ci-dessous le code pour effectuer la corrélation avec NumPy (avec ce qu'on a vu précédemment pour corriger l'amplitude et l'axe temporel). Je vous invite ensuite à recopier le code de la cellule dans vos notes persos, afin de pouvoir faire une corrélation correcte à tout moment, en TDM ou en examen..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9c6dd210",
   "metadata": {},
   "outputs": [],
   "source": [
    "# GOLDEN CELLULE A COMPLETER\n",
    "# CORRELATION ENTRE DEUX SIGNAUX (avec NumPy)\n",
    "\n",
    "# Paramètres (MODIFIEZ AVEC VOS NOMS DE VARIABLES)\n",
    "# ... :  deux signaux temporels NumPy\n",
    "# ... : deux axes temporels associés à A et B\n",
    "# ... : le temps d'échantillonnage\n",
    "\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "04a2147f",
   "metadata": {},
   "source": [
    "**_QUESTION :_** A quels instants la corrélation est maximum ? Qu'est-ce que cela signifie ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7694d07d",
   "metadata": {},
   "source": [
    "**_REPONSE :_**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6b7ee3b4",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Exercice 2 : Bon, on va pouvoir coder maintenant !\n",
    "\n",
    "Maintenant qu'on sait le faire avec NumPy, on va coder la fonction à la main, pas à pas.\n",
    "\n",
    "Tout d'abord, créez une fonction de padding. Le padding consiste à augmenter la taille d'un tableau/vecteur aux dimensions souhaitées en rajoutant une valeur constante plusieurs fois pour compléter. Ici, nous voulons simplement rallonger la taille d'un signal, donc ajouter des valeurs à la fin pour avoir la même longueur. \n",
    "\n",
    "Votre fonction prendra en entrée un signal de taille $N$, une taille $M$ (entier), et une valeur constante à insérer $cst\\_value$. Si $M$>$N$, alors la fonction renverra le signal avec des $cst\\_value$ ajoutés à la fin de ce dernier pour qu'il ait une taille $M$. Dans le cas inverse, la fonction renverra le signal original.\n",
    "\n",
    "*Note : il existe une fonction np.pad, mais je vous demande de ne pas l'utiliser car :*\n",
    "- *Elle n'est pas très pratique dans notre cas (plus complexe que nécessaire)*\n",
    "- *On va être embêté si on l'utilise pour la suite du TDM*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e7fe92f2",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Développement de la fonction de padding\n",
    "def padding(signal,M, cst_value):\n",
    "    \"\"\"\n",
    "    Fonction de padding du signal pour l'agrandir à la taille M (si M supérieur à len(signal))\n",
    "    signal : tableau 1D Numpy, représentant un signal\n",
    "    M : entier, représentat la taille de padding attendue pour M\n",
    "    cst_value : une valeur constante à insérer pour pad le signal\n",
    "    \"\"\"\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5d25c1ae",
   "metadata": {},
   "source": [
    "Afin de savoir si votre fonction est correcte, voici quelques tests unitaires. Si un test ne passe pas, ça ne sert à rien d'aller à la suite du TP, validez d'abord votre fonction de padding !"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3ba22685",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Tests unitaires de la fonction de padding\n",
    "x = np.array([1,3,7,2,5,6,6,8,1])\n",
    "print(f\"Padding de x avec M = 15, cst_value = 0 : {'OK' if ((len(padding(x,15,0))==15) & (padding(x,15,0)==np.array([1,3,7,2,5,6,6,8,1,0,0,0,0,0,0]))).all() else 'KO'}\")\n",
    "print(f\"Padding de x avec M = 15, cst_value = 1.5 : {'OK' if ((len(padding(x,15,1.5))==15) & (padding(x,15,1.5)==np.array([1,3,7,2,5,6,6,8,1,1.5,1.5,1.5,1.5,1.5,1.5]))).all() else 'KO'}\")\n",
    "print(f\"Padding de x avec M = 1 (aucun padding attendu) : {'OK' if ((len(padding(x,1,0))==len(x)) & (padding(x,1,0)==x).all()) else 'KO'}\")\n",
    "print(f\"Padding de x avec M = len(x) (aucun padding attendu) : {'OK' if ((len(padding(x,len(x),0))==len(x)) & (padding(x,len(x),0)==x).all()) else 'KO'}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "63cf9a13",
   "metadata": {},
   "source": [
    "Développez maintenant une fonction qui effectue un décalage temporel de $k$ sur un tableau NumPy uni-dimensionnel. La fonction retournera un tableau NumPy de même taille, avec les valeurs du tableau d'origine décalées de $k$ cases, et le reste des valeurs du tableau remplacées par une valeur constante $cst\\_value$, comme précédemment.\n",
    "\n",
    "*Exemple : Soit le signal $x$ :*\n",
    "\n",
    "\\begin{equation*}\n",
    "x = \\begin{bmatrix}1 & 3 & 7 & 2 & 5 & 6 & 6 & 8 & 1 \\end{bmatrix}\n",
    "\\end{equation*}\n",
    "\n",
    "*Un décalage de $x$ avec $k = 2$ et $cst\\_value = 0$ donnera le signal $x'$ suivant :*\n",
    "\n",
    "\\begin{equation*}\n",
    "x' = \\begin{bmatrix}0 & 0 & 1 & 3 & 7 & 2 & 5 & 6 & 6 \\end{bmatrix}\n",
    "\\end{equation*}\n",
    "\n",
    "*Un décalage de $x$ avec $k$ = -5 et $cst\\_value = -1$ donnera le signal $x''$ suivant :*\n",
    "\n",
    "\\begin{equation*}\n",
    "x'' = \\begin{bmatrix}6 & 6 & 8 & 1 & -1 & -1 & -1 & -1 & -1 \\end{bmatrix}\n",
    "\\end{equation*}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "95135697",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Développement d'une fonction qui décale un signal d'un pas de décalage k.\n",
    "def shift(x, k, cst_value):\n",
    "    \"\"\"\n",
    "    Effectue un décalage de k sur le signal x\n",
    "    x : tableau 1D NumPy, représentant un signal\n",
    "    k : entier, représentant le pas de décalage à effectuer\n",
    "    cst_value : la valeur constante à insérer là où le signal n'existe pas\n",
    "    \"\"\"\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f536fdf8",
   "metadata": {},
   "source": [
    "Afin de savoir si votre fonction est correcte, voici quelques tests unitaires. Si un test ne passe pas, ça ne sert à rien d'aller à la suite du TP, validez d'abord votre fonction de décalage !"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f5887c05",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Tests unitaires de la fonction de décalage\n",
    "x = np.array([1,3,7,2,5,6,6,8,1])\n",
    "print(f\"Décalage de x avec S = 2, cst_value = 0 : {'OK' if (shift(x,2,0)==np.array([0,0,1,3,7,2,5,6,6])).all() else 'KO'}\")\n",
    "print(f\"Décalage de x avec S = -5, cst_value = -2 : {'OK' if (shift(x,-5,-2)==np.array([6,6,8,1,-2,-2,-2,-2,-2])).all() else 'KO'}\")\n",
    "print(f\"Décalage de x avec S = 0 : {'OK' if (shift(x,0,0)==x).all() else 'KO'}\")\n",
    "print(f\"Décalage de x avec S = 20, cst_value = 1.5 : {'OK' if (shift(x,20,1.5)==np.full((len(x)),1.5)).all() else 'KO'}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f8f8e339",
   "metadata": {},
   "source": [
    "Développez maintenant une fonction qui calcule le score d'intercorrélation entre deux signaux $f$ et $g$ avec un retard $k$ fixé (entier). Ici, $k$ désigne le décalage de cases du tableau à réaliser pour le signal $g$.\n",
    "\n",
    "Pour cela, il faudra, dans l'ordre :\n",
    "- Récupérer les tailles des deux signaux $f$ et $g$\n",
    "- Effectuer un padding sur les deux signaux pour qu'ils soient de même taille, en remplissant avec des 0\n",
    "- Décaler le signal $g$ de $k$ cases (toujours en ajoutant des 0)\n",
    "- Calculer le score d'intercorrélation entre les deux signaux, en suivant la fonction écrite préalablement\n",
    "\n",
    "A vous de jouer !"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1002a51e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Développement de la fonction de calcul du score d'intercorrélation entre deux signaux à un décalage k fixé\n",
    "def correlation_at_k(f, g, k, Te):\n",
    "    \"\"\"\n",
    "    Calcule le score d'intercorrélation entre le signal f et le signal g retardé de k\n",
    "    f, g : tableaux NumPy 1D, signaux à corréler\n",
    "    k : entier, retard du signal g\n",
    "    Te : le temps d'échantillonage (en secondes)\n",
    "    \"\"\"\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "92472cbd",
   "metadata": {},
   "source": [
    "On va maintenant pouvoir développer la fonction principale de corrélation. Cette fonction prendra en entrée deux signaux $f$ et $g$, ainsi que le pas d'échantillonage des deux signaux $T_e$. \n",
    "\n",
    "Nous avons utilisé depuis le début de l'exercice la variable $k$ pour représenter le décalage à effectuer sur le tableau (valeur entière). Cependant, $k$ n'a pas de valeur physique. De ce fait, en le multipliant par le pas d'échantillonage $T_e$, vous obtenez la valeur $\\tau$ qui est le retard dans la formule de corrélation.\n",
    "\n",
    "\\begin{equation*}\n",
    "    \\tau = k \\times T_e\n",
    "\\end{equation*}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cd351a38",
   "metadata": {},
   "source": [
    "Développez maintenant la fonction de correlation, prenant en entrée les signaux $f$ et $g$ et le pas d'échantillonage $T_e$. La fonction renverra :\n",
    "- L'axe temporel de la corrélation des deux signaux, à savoir les retards $\\tau$ possibles où $f$ et $g$ retardés se recoupent\n",
    "- Le score d'intercorrélation entre $f$ et $g$ pour chaque retard $\\tau$ (sous forme de tableau NumPy 1D de même taille)\n",
    "\n",
    "*Rappel : Pour l'axe temporel, on a déjà effectué le calcul dans l'exercice 1*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7cbccba2",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Fonction de corrélation entre les signaux f et g avec un pas d'échantillonage de Te\n",
    "def correlation(f,g,Te):\n",
    "    \"\"\"\n",
    "    Calcule la corrélation entre le signal f et le signal g\n",
    "    f, g : tableaux NumPy 1D, signaux à corréler\n",
    "    Te : le temps d'échantillonage (en secondes)\n",
    "    \"\"\"\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5b4a4a8f",
   "metadata": {},
   "source": [
    "Testez maintenant votre fonction de corrélation sur vos deux signaux $A$ et $B$, et tracez le résultat obtenu en fonction des valeurs de $\\tau$ possibles."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "edc4da4b",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Calcul de la corrélation entre A et B\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "64e73385",
   "metadata": {},
   "source": [
    "Vous pouvez maintenant comparer votre résultat de corrélation avec celui effectué par la fonction correlate de NumPy. Affichez les deux signaux calculés, et vérifiez que vous obtenez exactement le même résultat"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b942f70d",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Affichage des deux résultats de la corrélation (NumPy et la votre)\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e867ceff",
   "metadata": {},
   "source": [
    "Super, on a notre fonction de corrélation !\n",
    "\n",
    "...\n",
    "\n",
    "Cependant...\n",
    "\n",
    "Notre fonction n'est en réalité pas très robuste. En effet, en ajoutant une seule valeur abberante, on peut complètement la casser ! (emoji choqué)\n",
    "\n",
    "Reprenez votre signal $A$, et ajoutez une seule valeur aberrante (ex: à $t = 4.5$, $A(t)=100$). Recalculez la corrélation et affichez le résultat."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "90694121",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Ajout d'une valeur aberrante à A(t)\n",
    "...\n",
    "\n",
    "# Calcul de la corrélation entre A et B avec la valeur aberrante\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "64e575c8",
   "metadata": {},
   "source": [
    "Si vous avez bien codé, vous devriez avoir une courbe de corrélation toute chamboulée par ce désagrément. Cela est du au fait que notre fonction de corrélation est ultra dépendante de l'amplitude. De ce fait, une amplitude grande donne de meilleurs scores sur le calcul de corrélation. Le soucis ici est qu'on cherche avant tout à détecter un motif, indépendamment de l'amplitude. Pour cela, nous allons modifier notre fonction afin d'éviter ce problème."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8295b194",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Exercice 3 : ZNCC is the new correlation\n",
    "\n",
    "Dans ce nouvel exercice, nous allons reprendre la fonction de corrélation pour faire une corrélation normalisée, appelée ZNCC (Zero-Mean Normalized Cross Correlation). La formule est quasi la même que celle que nous avons vu précédemment, à quelques détails prêts:\n",
    "\n",
    "\\begin{equation*}\n",
    "ZNCC_{fg}(\\tau) =\\;\n",
    "\\frac{\\displaystyle \\sum_n \\bigl(f(n) - \\bar{f}\\bigr)\\,\\bigl(g^{*}(n - \\tau) - \\bar{g}^{*}\\bigr)}\n",
    "     {\\sqrt{\\displaystyle \\sum_n \\bigl(f(n) - \\bar{f}\\bigr)^2}\\;\n",
    "      \\sqrt{\\displaystyle \\sum_n \\bigl(g^{*}(n - \\tau) - \\bar{g}^{*}\\bigr)^2}}\n",
    "\\end{equation*}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4331834a",
   "metadata": {},
   "source": [
    "Typiquement, ça reprend la formule de corrélation utilisée préalablement, avec les signaux $f$ et $g$ normalisés, c-à-d en soustrayant leur moyenne $\\bar{f}$ et $\\bar{g}^{*}$, puis en divisant le score de corrélation par la multiplication de leur norme L2. On remarque cependant que contrairement à la formule de la corrélation originale, ZNCC est indépendante du temps d'échantillonage $T_e$.\n",
    "\n",
    "La difficulté ici est que le calcul des moyennes et des normes ne se font pas sur les signaux entiers, mais là où $f$ et $g$ se recoupent (leur intersection)."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8e68f23d",
   "metadata": {},
   "source": [
    "*\"Mais comment on va faire ?! On va devoir tout recommencer !!\"*\n",
    "\n",
    "Non, rassure-toi, on va pouvoir reprendre nos précédentes fonctions et on va faire les quelques modifications pour obtenir la ZNCC.\n",
    "\n",
    "Pour calculer l'intersection entre deux signaux, c'est difficile à faire avec nos paddings et décalages précédents. En effet, nous avons déterminé qu'on allait remplir tout ça par des 0, hors un 0 ne signifie pas qu'un signal n'existe pas à l'instant t, mais qu'il est nul. Du coup, il va falloir modifier cette constante afin de signifier que le signal n'existe pas à l'endroit où on calcule. On va utiliser donc une valeur spéciale de NumPy : *np.nan* (pour \"Not a Number\").\n",
    "\n",
    "Grâce à cette valeur spéciale (et la fonction *np.isnan()*, créeant un masque indiquant si chaque valeur sont des *np.nan* ou non), on va pouvoir déterminer facilement l'intersection entre les deux signaux, à savoir les indices temporels où les deux signaux ont des valeurs réelles (c'est-à-dire pas des *np.nan*)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "95a55247",
   "metadata": {},
   "source": [
    "Reprenez la fonction de calcul de score d'intercorrélation entre deux signaux à un retard $k$ (entier) fixé et modifiez-la pour en calculer le score de ZNCC. Pour cela, il faut, dans l'ordre :\n",
    "\n",
    "- Récupérer les tailles des deux signaux $f$ et $g$\n",
    "- Effectuer un padding sur les deux signaux pour qu'ils soient de même taille (avec comme constante *np.nan*)\n",
    "- Décaler le signal $g$ de $k$ cases (avec comme constante *np.nan*)\n",
    "- Coupez les signaux $f$ et $g^{*}$ à l'endroit de leur intersection\n",
    "- Si la taille de l'intersection est trop petite (inférieure au quart du plus petit signal), renvoyer 0\n",
    "- Calculer et retirer la moyenne des signaux $f$ et $g^{*}$\n",
    "- Calculer le score d'intercorrélation entre les deux signaux, en suivant la fonction écrite préalablement\n",
    "- Diviser ce score par la norme de $f$ et $g^{*}$\n",
    "\n",
    "*Note : La norme de $f$ et de $g$ pouvant être nulles, dans le cas où l'une ou l'autre vaut zéro, il faut lui rajouter $\\epsilon = 1e-12$ pour éviter une division par zéro* "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c719b657",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Développement de la fonction de calcul du score d'intercorrélation entre deux signaux à décalage k fixé\n",
    "def ZNCC_at_k(f,g,k):\n",
    "    \"\"\"\n",
    "    Calcule le score de ZNCC entre le signal f et le signal g retardé de k\n",
    "    f, g : tableaux NumPy 1D, signaux à corréler\n",
    "    k : entier, retard du signal f\n",
    "    \"\"\"\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b48db955",
   "metadata": {},
   "source": [
    "Reprenez votre fonction de calcul de corrélation qui prend en entrée $f$, $g$ et $T_e$. Modifiez-la pour qu'elle renvoit les valeurs de retards calculés en secondes et le score de ZNCC.\n",
    "\n",
    "*Note : Il n'y a quasimment rien à changer, il faut juste changer le nom de la fonction et l'appel à une fonction...*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bb846a21",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Fonction de ZNCC entre les signaux f et g avec un pas d'échantillonage de Te\n",
    "def ZNCC(f,g,Te):\n",
    "    \"\"\"\n",
    "    Calcule la ZNCC entre le signal f et le signal g\n",
    "    f, g : tableaux NumPy 1D, signaux à corréler\n",
    "    Te : le temps d'échantillonage (en secondes)\n",
    "    \"\"\"\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "52967dd4",
   "metadata": {},
   "source": [
    "Reprenez les signaux $A$ et $B$, calculez le score de ZNCC entre ces signaux et affichez le résultat."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "25adb627",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Calcul et affichage de la ZNCC entre A et B\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3d9e3e64",
   "metadata": {},
   "source": [
    "**_QUESTIONS :_** \n",
    "- Dans quel intervalle est inclut le score de ZNCC ?\n",
    "- A quels instants ZNCC est maximum ? \n",
    "- A quels instants ZNCC est minimum ? \n",
    "- Qu'est-ce que cela signifie ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f0768279",
   "metadata": {},
   "source": [
    "**_REPONSES :_** "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fd7b3c43",
   "metadata": {},
   "source": [
    "Reprenez maintenant le signal bruité $A$, comme fait précédemment, calculez le ZNCC entre $A$ et $B$ bruité, et affichez le résultat."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "60a77746",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A compléter\n",
    "# Calcul et affichage de la ZNCC entre A et B avec la valeur aberrante\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "11a4a206",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Que remarquez-vous ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d342fd26",
   "metadata": {},
   "source": [
    "**_REPONSE :_** "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "668d5a74",
   "metadata": {},
   "source": [
    "En conclusion, ZNCC est plus robuste car insensible à l'amplitude du signal. Avec ZNCC, on cherche avant tout un motif, qu'il soit fort, faible, ou même inversé. Malheureusement, ZNCC n'est pas disponible dans une librairie Python (je suis surpris comme vous), du moins pas sous la forme qu'on souhaite pour le traitement du signal."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "da577114",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## Exercice 4 : Où est caché DJ Khaled ?\n",
    "\n",
    "Grâce à votre superbe fonction ZNCC, vous pouvez trouver n'importe quel pattern dans un signal. Dans la série des patterns reconnaissables, on peut citer les tags audios des producteurs de musique, qui est l'équivalent d'une signature généralement au début du morceau. Parmi les tags audios célèbres, DJ Khaled en a plusieurs, dont son fameux \"Another One\".\n",
    "\n",
    "Dans cet exercice, l'objectif est simple: à partir du tag audio de DJ Khaled, à vous de retrouver dans la bibliothèque audio quel son a été produit par lui.\n",
    "\n",
    "*Note : En écoutant les WAV, vous pourrez remarquer que la qualité n'est pas terrible. C'est normal : la fréquence d'échantillonage a été fortement diminuée, afin de réduire drastiquement le temps de calcul...*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "91527427",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Pour lire un audio WAV directement depuis un notebook\n",
    "\n",
    "import IPython\n",
    "\n",
    "IPython.display.Audio('audios/tag.wav')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1dce2332",
   "metadata": {},
   "source": [
    "Dans un premier temps, chargez le fichier audio tag.wav. Pour cela, vous avez la fonction *read* du module *wavfile* du package *scipy.io*."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e341a3dc",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Chargement de l'extrait audio\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "72f59166",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Que renvoit la fonction *read* à la lecture du fichier WAV ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0a191780",
   "metadata": {},
   "source": [
    "**_REPONSE :_**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ad024fc4",
   "metadata": {},
   "source": [
    "Tracez le signal du tag audio, en reconstruisant l'axe temporel préalablement."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7a7b1728",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Reconstruction de l'axe temporel\n",
    "...\n",
    "\n",
    "# Tracé du signal audio\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9b4c6557",
   "metadata": {},
   "source": [
    "Faites de même avec le fichier *music_1.wav*. Chargez-le, reconstruisez l'axe temporel de la musique et tracez le signal."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "809b05b6",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Chargement de l'extrait audio music_1\n",
    "...\n",
    "\n",
    "# Reconstruction de l'axe temporel\n",
    "...\n",
    "\n",
    "# Affichage du signal\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7e35e736",
   "metadata": {},
   "source": [
    "Appliquez votre fameux algorithme ZNCC pour le tag audio et l'extrait musical, et affichez le résultat.\n",
    "\n",
    "*Note : le temps de calcul peut être long, vu la taille du signal... Et encore, ce n'est qu'un extrait de musique...*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "61a796a6",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Calcul de ZNCC pour le tag audio et l'extrait musical\n",
    "...\n",
    "\n",
    "# Affichage du score ZNCC en fonction du retard tau\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "812ede54",
   "metadata": {},
   "source": [
    "Si tout est ok, il devrait y avoir un pic qui ressort ! Récupérez la valeur de $\\tau$ où ZNCC est maximum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6747ca18",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Récupération du tau optimal\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee98dccd",
   "metadata": {},
   "source": [
    "Coupez music_1 de manière à avoir l'extrait où vous avez le meilleur score ZNCC (à partir du $\\tau$ optimal)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fc910b28",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Extraction de l'audio avec le meilleur score ZNCC\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "296d83eb",
   "metadata": {},
   "source": [
    "Exportez dans un nouveau fichier WAV l'audio du meilleur extrait trouvé via la fonction *write* du module *wavfile* du package *scipy.io*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "38e17ef9",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Export du meilleur extrait en fichier WAV\n",
    "..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "43869959",
   "metadata": {},
   "source": [
    "**_QUESTION :_** Alors, est-ce que l'extrait trouvé est pertinent ?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b82bd61d",
   "metadata": {},
   "source": [
    "**_REPONSE :_**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22a840eb",
   "metadata": {},
   "source": [
    "Vous pouvez maintenant continuer votre recherche avec les autres audios. Ecrivez pour finir les extraits audios où vous avez trouvé le tag."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "443148b9",
   "metadata": {},
   "outputs": [],
   "source": [
    "# A COMPLETER\n",
    "# Même processus que précédemment avec les autres extraits audios\n",
    "..."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "deep",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
