m-thesis-introduction/01_fourier.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Fourier Transforms"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "import numpy as np\n",
    "import scipy.fft as ft\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.gridspec as gridspec\n",
    "rng = np.random.default_rng()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Fourier transform (FT) of a quantity $t$ is defined as\n",
    "\n",
    "$$\n",
    "\\hat{f}(\\omega) = F(f(t)) = \\frac{1}{\\sqrt{2\\pi}} \\int_{-\\infty}^{\\infty} \\mathrm{d}t\\, f(t)\\, \\exp{ \\left(-i \\omega t \\right)}\n",
    ".\n",
    "$$\n",
    "\n",
    "In signal processing, taking the integral from $t = -\\infty$ to $t = \\infty$ is not practical. Together with a finite sample space $t$, the integral can be turned into a (discrete) sum, giving a Discrete Fourier Transform (DFT)\n",
    "\n",
    "$$\n",
    "F(t) = \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N-1} f(t_k)\\, \\exp{ \\left(-i \\omega t_k \\right)}\n",
    "\\label{eq:DFT} \\tag{1}\n",
    ",\n",
    "$$\n",
    "where the sum iterates over each sample $t_k$ up to the end ($k=N-1$) of the signal."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Interdependence $\\Delta t$, $T$, $\\Delta f$, $f_{min}$, $f_{max}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The highest frequency $f_\\mathrm{max}$ (also known as the Nyquist frequency) that can be (reliably) measured within a signal is dependent only on the sampling rate ($1/\\Delta t$) of the signal with\n",
    "$$\n",
    "    f_\\mathrm{max} = \\frac{1}{2} f_\\mathrm{sampling} = \\frac{1}{2} \\frac{1}{\\Delta t}\n",
    "    .\n",
    "$$\n",
    "Higher frequencies will fall in between two adjacent samples, causing them to go unnoticed in the data. Note that because of folding, such frequencies can actually appear as aliases.\n",
    "\n",
    "Likewise, the lowest frequency depends on the length of the signal $T$, such that over the whole signal one oscillation can be found\n",
    "$$\n",
    "    f_\\mathrm{min} \\sim \\frac{1}{T}\n",
    "    .\n",
    "$$\n",
    "Even lower frequencies will no longer be detected as a frequency. Instead, they will introduce a bias for the $0 Hz$ component in the spectrum."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# FTs of ideal signals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def ft_spectrum( signal, sample_rate, fft=None, freq=None, mask_bias=False):\n",
    "    \"\"\"Return a FT of $signal$, with corresponding frequencies\"\"\"\n",
    "    n_samples = len(signal)\n",
    "    real_signal = np.all(np.isreal(signal))\n",
    "    \n",
    "    if fft is None:\n",
    "        if real_signal:\n",
    "            fft = ft.rfft\n",
    "            freq = ft.rfftfreq\n",
    "        else:\n",
    "            fft = ft.fft\n",
    "            freq = ft.fftfreq\n",
    "\n",
    "    if freq is None:\n",
    "        freq = ft.fftfreq\n",
    "    \n",
    "    spectrum = fft(signal) / sample_rate\n",
    "    freqs = freq(n_samples, 1/sample_rate)\n",
    "    \n",
    "    if not mask_bias:\n",
    "        return spectrum, freqs\n",
    "    else:\n",
    "        return spectrum[1:], freqs[1:]\n",
    "\n",
    "    \n",
    "def plot_spectrum( ax, spectrum, freqs, plot_complex=False, plot_power=False, plot_amplitude=None):\n",
    "    \"\"\" Plot a signal's spectrum on an Axis object\"\"\"\n",
    "    plot_amplitude = plot_amplitude or (not plot_power and not plot_complex)\n",
    "    alpha = 1\n",
    "    \n",
    "    ax.set_title(\"Spectrum\")\n",
    "    ax.set_xlabel(\"f (Hz)\")\n",
    "    ylabel = \"\"\n",
    "    if plot_amplitude:\n",
    "        ylabel = \"Amplitude\"\n",
    "    if plot_power:\n",
    "        if ylabel:\n",
    "            ylabel += \"|\"\n",
    "        ylabel += \"Power\"\n",
    "    ax.set_ylabel(ylabel)\n",
    "\n",
    "    if plot_complex:\n",
    "        alpha = 0.5\n",
    "        ax.plot(freqs, np.real(spectrum), '.-', label='Real', alpha=alpha)\n",
    "        ax.plot(freqs, np.imag(spectrum), '.-', label='Imag', alpha=alpha)\n",
    "        ax.legend()\n",
    "\n",
    "    if plot_power:\n",
    "        ax.plot(freqs, np.abs(spectrum)**2, '.-', label='Power', alpha=alpha)\n",
    "    \n",
    "    if plot_amplitude:\n",
    "        ax.plot(freqs, np.abs(spectrum), '.-', label='Abs', alpha=alpha)\n",
    "\n",
    "    return ax\n",
    "\n",
    "\n",
    "def plot_phase( ax, spectrum, freqs):\n",
    "    ax.set_ylabel(\"Phase\")\n",
    "    ax.set_xlabel(\"f (Hz)\")\n",
    "\n",
    "    ax.plot(freqs, np.angle(spectrum), '.-')\n",
    "    ax.set_ylim(-1*np.pi, np.pi)\n",
    "    \n",
    "    return ax\n",
    "\n",
    "\n",
    "def plot_combined_spectrum(spectrum, freqs, \n",
    "                           spectrum_kwargs={}, fig=None, gs=None):\n",
    "    \"\"\"Plot both the frequencies and phase in one figure.\"\"\"\n",
    "    \n",
    "    # configure plotting layout\n",
    "    if fig is None:\n",
    "        fig = plt.figure(figsize=(8, 16))\n",
    "\n",
    "    if gs is None:\n",
    "        gs = gridspec.GridSpec(2, 1, figure=fig, height_ratios=[2,1])\n",
    "\n",
    "    ax1 = fig.add_subplot(gs[:-1, -1])\n",
    "    ax2 = fig.add_subplot(gs[-1, -1], sharex=ax1)\n",
    "\n",
    "    axes = np.array([ax1, ax2])\n",
    "    \n",
    "    # plot the spectrum \n",
    "    plot_spectrum(ax1, spectrum, freqs, **spectrum_kwargs)\n",
    "    \n",
    "    # plot the phase\n",
    "    plot_phase(ax2, spectrum, freqs)\n",
    "    \n",
    "    return fig, axes\n",
    "    \n",
    "\n",
    "def plot_signal( ax, signal, sample_rate = 1):\n",
    "    ax.set_title(\"Signal\")\n",
    "    ax.set_xlabel(\"t (s)\")\n",
    "    ax.set_ylabel(\"A(t)\")\n",
    "\n",
    "    ax.plot(np.arange(len(signal)) / sample_rate, signal)\n",
    "    \n",
    "    return ax\n",
    "\n",
    "\n",
    "def plot_signal_and_spectrum( signal, sample_rate, title=None,\n",
    "                             ft_kwargs = {}, spectrum_kwargs = {}):\n",
    "    \"\"\"Plot a signal and its spectrum in one go.\"\"\"\n",
    "    fig = plt.figure(figsize=(16, 4))\n",
    "    \n",
    "    if title:\n",
    "        fig.suptitle(title)\n",
    "\n",
    "    # setup plot layout\n",
    "    gs0 = gridspec.GridSpec(1, 2, figure=fig)\n",
    "    gs00 = gs0[0].subgridspec(1, 1)\n",
    "    gs01 = gs0[1].subgridspec(2, 1)\n",
    "      \n",
    "    # plot the signal\n",
    "    ax1 = fig.add_subplot(gs00[0, 0])\n",
    "    plot_signal(ax1, signal)\n",
    "    \n",
    "    # plot spectrum\n",
    "    signal_fft, freqs = ft_spectrum(signal, sample_rate, **ft_kwargs)\n",
    "    _, (ax2, ax3) = plot_combined_spectrum(signal_fft, freqs, spectrum_kwargs = spectrum_kwargs, fig=fig, gs=gs01)\n",
    "\n",
    "    # return the axes\n",
    "    axes = np.array([ax1, ax2, ax3])\n",
    "    \n",
    "    return (fig, axes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "sample_rate = 1/1e-3 # s\n",
    "n_samples = 2**10\n",
    "t = np.arange(n_samples) / sample_rate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### FTs: constant signal\n",
    "\n",
    "A constant signal (also termed offset or bias) results in only a single component in the Fourier transform, at the zero frequency.\n",
    "\n",
    "Since often such a component is included in a signal (for example when a signal is not on average zero), it can help to mask out the zero frequency when plotting the spectrum.\n",
    "Doing this for a constant signal, we get a totally flat spectrum."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA7wAAAEjCAYAAAAR7A2IAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3dedxdZX3v/c+XJIDIECCRKYQERZE6AEYEjYWKVXDCoQ6IRYaW8qhH26PHg62PejzHp7Z1PloVBCcsanHCihVqUQuKkggyozFhCAkkzIgyhPyeP/aK3t7eY+5h3fe6P+/Xa7+y11rXXvu7Fztc+e1rrWulqpAkSZIkqWu2aDuAJEmSJEkTwYJXkiRJktRJFrySJEmSpE6y4JUkSZIkdZIFryRJkiSpkyx4JUmSJEmdZMErSVIHJDkmyXmT8D6HJVk90e8jSdJ4sOCVJHVGklcnWZbkV0nWJvl2kqUT+H7jWvwlOS7JhZvz2qr6QlU9Z7yySJLUBRa8kqROSPLfgQ8B/x+wC7AQ+GfgqDZzSZKk9ljwSpKmvSQ7AO8GXl9VX62q+6rqoar6ZlX9j6bNVkk+lGRN8/hQkq2abYclWZ3kzUnWNaPDx/fZ//OSXJ3k3iQ3J3lLkkcC3wZ2b0aUf5Vk9yQHJflRkrua/Xw0yZZ99lVJTk7yiyR3JvlYeh4PfAI4pNnXXYN81uOSrGyyrEpyTJ/1F/Zp95wk1yW5O8k/J/l+kr/o2zbJ+5oMq5Ic2ee1xye5pnmPlUn+avz+a0mSNHkseCVJXXAIsDXwtSHa/B1wMLA/8GTgIODtfbbvCuwA7AGcCHwsyY7NttOBv6qq7YAnAP9ZVfcBRwJrqmrb5rEGeBj4G2Bek+tw4HX9srwAeGqT4xXAc6vqGuBk4EfNvub2/wBNkf0R4Mgmy9OBywZoNw84G3gbsDNwXdO2r6c16+cB/wicniTNtnVNxu2B44EPJjmw//tIkjTVWfBKkrpgZ+C2qtowRJtjgHdX1bqqWg/8L+DP+2x/qNn+UFWdC/wKeFyfbfsl2b6q7qyqnw72JlW1vKourqoNVXU98Eng0H7N3ltVd1XVjcAF9IrwkdoIPCHJI6pqbVVdNUCb5wFXNaPdG+gVybf0a3NDVZ1WVQ8DnwV2o3cqOFX1rar6ZfV8HzgPeOYoMkqSNCVY8EqSuuB2YF6S2UO02R24oc/yDc263+6jX8H8a2Db5vnL6BWRNzSnBh8y2JskeWySf0tyS5J76F1TPK9fs77FZ9/3GVIzqvxKeiPBa5N8K8m+AzTdHbipz+sK6D+51i19tv+6ebpt8xmOTHJxkjuaU6ufN8BnkCRpyrPglSR1wY+A+4EXD9FmDbBXn+WFzbphVdUlVXUU8Cjg68CXN20aoPnHgWuBfapqe+BvgQzQbsC3GkGW71TVn9Ibkb0WOG2AZmuBBZsWmlOVFwzQ7g801zV/BXgfsEtzavW5jPwzSJI0ZVjwSpKmvaq6G3gHvetuX5xkmyRzmpHKf2yanQW8Pcn85hrXdwBnDrfvJFs297jdoaoeAu6hd50uwK3Azs2kWZts17T5VTP6+v+M4qPcCizoO8lVvyy7JHlRcy3vA/ROu354gKbfAp7YHIvZwOvpXaM8ElsCWwHrgQ3NZFbe7kiSNC1Z8EqSOqGqPgD8d3oTUa2nd0rvG+iNyAL8H2AZcDlwBfDTZt1I/DlwfXOK8snAa5r3vJZeIb2ymZV5d+AtwKuBe+mNvn5pFB/jP4GrgFuS3DbA9i2AN9Mbmb6D3rXB/SfEoqpuA15ObzKq24H96H32B4YLUFX3Am+kN4p9Z/NZzhnFZ5AkacpI77IeSZLUVUm2oHcN7zFVdUHbeSRJmiyO8EqS1EFJnptkbnNN7qbriC9uOZYkSZPKgleSpG46BPglcBvwQuDFVfWbdiNJkjS5PKVZ6ogkxwCvraoJnVwmyWHAmVU1ohlfJUmSpLY4witNM0mWJvlhkrube2RelOSpVfWFiS52JUmaKQbrbyfw/a5P8uyJ2r80U81uO4CkkUuyPfBv9G5z8mV6tw95JiOYeVWSJI3MVOxvk8yuqg1tvb80XTnCK00vjwWoqrOq6uGq+k1VnVdVlyc5LsmFmxomeU6S65pfpv85yfeT/EWz7bgkFyZ5X5I7k6xq7rW56bXHJ7kmyb1JVib5q8n/qJIktWa4/vaiJP+36WOvTXL4phcm2SHJ6UnWJrk5yf9JMqvP9r/s08deneTAJJ8HFgLfTPKrJG9NsihJJTkxyY3AfyY5LMnqvkH7jgwneVeSf01yZrP/K5I8NsnbkqxLclMSzwbTjGLBK00vPwceTvLZJEcm2XGgRknmAWcDbwN2Bq4Dnt6v2dOa9fPo3avz9CRptq0DXgBsDxwPfDDJgeP9YSRJmqKG62+fBqyk14e+E/hqkp2abZ8FNgCPAQ4AngNs+sH55cC7gGPp9bEvAm6vqj8HbgReWFXbVtU/9nmvQ4HHA88dYfYXAp8HdgQuBb5D79/8ewDvBj45wv1InWDBK00jVXUPsBQo4DRgfZJzkuzSr+nzgKuq6qvN6U8fAW7p1+aGqjqtqh6m1znvBuzSvM+3quqX1fN94Dx6p3JJktR5I+hv1wEfqqqHqupL9H5Afn6z/Ujgr6vqvqpaB3wQeFXzur8A/rGqLmn62BVVdcMwcd7V7Guks6z/V1V9p+n//xWYD7y3qh4CvggsSjJ3hPuSpj0LXmmaqaprquq4ZpbkJwC7Ax/q12x34KY+rylgdb82t/TZ/uvm6bYAza/ZFzeTdNxFr4CeN76fRJKkqWuY/vbm+v1bndzQbN8LmAOsTXJX04d+EnhU025PercLG42bhm/ye27t8/w3wG3Nj9ublqHp76WZwIJXmsaq6lrgM/Q64r7WAr+9bVBzqvKIbiOUZCvgK8D7gF2qai5wLpAhXyhJUkcN0N/u0ecyIOhdf7uGXnH6ADCvquY2j+2r6o+adjcBjx7sbUaw/j5gm00LzbXB80fzWaSZxoJXmkaS7JvkzUkWNMt7AkcDF/dr+i3giUlenGQ28Hpg1xG+zZbAVsB6YEMzmZUTXEiSZowR9LePAt6YZE5zXe7jgXOrai29y4Den2T7JFskeXSSQ5vXfQp4S5KnpOcxSfZqtt0K7D1MtJ8DWyd5fpI5wNvp9dmSBmHBK00v99KbKOPHSe6j1/FeCby5b6Oqug14Ob3JqG4H9gOWMYLbKVTVvcAb6d2G4U7g1cA54/cRJEma8obrb38M7APcBrwH+LOqur3Zdiy9H4+vptePnk1vngyq6l+b9v/SvMfXgU2TXf098PbmVOi3DBSqqu4GXkevcL6Z3ohv/0uWJPWR37/8QFIXJdmCXod4TFVd0HYeSZKmqyTHAX9RVUvbziJpeI7wSh2V5LlJ5jbX5P4tvWtw+5/6LEmSJHWWBa/UXYfQmwnyNnr35HvxKG5pIEmSJE17ntIsSZIkSeokR3glSZIkSZ00u+0A42XevHm1aNGitmNIkjpi+fLlt1WV97ccA/tmSdJ42py+uTMF76JFi1i2bFnbMSRJHZHkhrYzTHf2zZKk8bQ5fbOnNEuSJEmSOsmCV5IkTYh1997P8hvubDuGJGkGs+CVJEkkOSPJuiRXDrI9ST6SZEWSy5McONw+b73nAV592sUWvZKk1ljwSpIkgM8ARwyx/Uhgn+ZxEvDxkez0wQ0buXjl7WMOJ0nS5rDglSRJVNUPgDuGaHIU8LnquRiYm2S3YfcL7L/n3HFKKUnS6FjwSpKkkdgDuKnP8upm3e9JclKSZUmWbbPFRgBWrPvV5CSUJKkfC15JkjQSGWBd/cGKqlOraklVLXn0bjuy/55z+fRFq9i48Q+aSpI04Sx4JUnSSKwG9uyzvABYM9yLTly6mOtv/zX/ee26CQsmSdJgLHglSdJInAMc28zWfDBwd1WtHe5FRz5hV3bfYWtOv3DVxCeUJKkfC15JkkSSs4AfAY9LsjrJiUlOTnJy0+RcYCWwAjgNeN1I9jt71hYc+/RF/Gjl7Vy15u4JyS5J0mBmtx1AkiS1r6q
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA7wAAAEjCAYAAAAR7A2IAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deZwlVX338c8XhkUj+yACAwwILrhEsANoXFAUQVGMiqJEgRiRqI9LNAaXiCHxiTHukaAIKCqKiKhjHAOEAMY8YpgBRGFAx1FkYJABhsWFZeD3/FHVcml6nV7u7duf9+t1X32rzqmq362umdO/e06dSlUhSZIkSVK/Wa/bAUiSJEmSNB1MeCVJkiRJfcmEV5IkSZLUl0x4JUmSJEl9yYRXkiRJktSXTHglSZIkSX3JhFeSpCGSHJbknBk4zr5JVo5R5ytJXjzdsQw55gVJ/nIK9rMwSSWZNxVxrWMMGyW5KsnDuxWDJKl7THglScNK8qokS5L8JsmqJN9N8rRpPN6Yyd8E93dEku+vy7ZVdVpV7T9VsayrJE8E/hj4Vrdjma2q6i7gFOBvux2LJGnmmfBKkh4kyV8DHwf+L7ANsCPwb8DB3YxrDno9cFpVVbcD6VXj7D3+MnB4ko2mOx5JUm8x4ZUkPUCSzYDjgDdW1VlV9duquqeqvl1Vf9PW2SjJx5Nc374+PphMDPbUJnl7khvb3uEjO/b//CRXJrkjyXVJ3pHkj4DvAtu1Pcq/SbJdkr2S/CDJre1+PpVkw459VZKjk/wsyZokx6fxWODTwFPafd06wmc9IsmKNpZfJDmsY/33O+rtn+TqJLcl+bckFw4O+R2sm+TDbQy/SHJgx7ZHJlnWHmNFktdP4NdxIHDhkHj/J8nH2nOyIslT2/XXtuf78I76L0hyaZLb2/L3d5RtnORLSW5u93Vxkm2GOUfbJrk8yTva5c2SnNz+Pq5L8o9J1m/L1m/Pw01JVgAvGOmDJXlFx+/6N0nuSnJBW7ZRu59fJfl1kk8neUhbNnh9/W2SG4DPtetfl2R5kluSLEqy3eCxqmolsAbYZwLnXpLUB0x4JUlDPQXYGPjGKHXeQ5M8PIlmyO1ewHs7yh8BbAZsD7wWOD7JFm3ZycDrq2oT4PHAf1XVb2mSu+ur6mHt63rgXuBtwPw2rv2ANwyJ5SDgT9o4Xg48r6qWAUcDP2j3tfnQD9Am2Z8EDmxjeSpw2TD15gNnAu8CtgKubut22rtdPx/4EHBykrRlN7YxbgocCXwsyZ5DjzNCfDu3+x16rMvbWL4MnN5+/l2BPwc+leRhbd3fAq8BNqdJPv8q998PfDjN72iHdl9HA78fEsNCmoT7U1X14Xb1qcDa9nh7APsDg/f7vq79rHsAA8DLRvp8VfXVwd81sB2wAvhKW/zPwKNorq9daa6j93Vs/ghgS2An4Kgkzwb+ieb3vy1wTXteOi2juUYkSXOICa8kaaitgJuqau0odQ4DjquqG6tqNfD3wKs7yu9py++pqsXAb4BHd5TtnmTTqlpTVZeMdJCqWlpVF1XV2qr6JfAZ4JlDqn2wqm6tql8B59MkSeN1H/D4JA+pqlVVdcUwdZ4PXNH2dq+lSZJvGFLnmqr6bFXdS5MQbkszFJyq+k5V/bwaFwLnAE8fR2yDSfodQ9b/oqo+1x7rqzQJ63FVdVdVnQPcTZMkUlUXVNWPq+q+qrqcJqEcPH/30Pyud62qe9tzfXvHcXYHLgCOraoTAdoe4AOBt7Y9/zcCHwMObbd5OfDxqrq2qm6hSUJHlWQ9msT9gqr6TPtFweuAt1XVLVV1B83Q+kM7Nruvjeuuqvo9zfV4SlVd0t6z+y6a3v2FHdvc0XFOJUlzhAmvJGmom4H5Gf3eyO1oetEGXdOu+8M+hiTMvwMGex1fSpNEXtMODX7KSAdJ8qgk/57khiS30yQ+84dU60w+O48zqrZX+RU0PZurknwnyWOGqbodcG3HdgUMnVzrho7y37VvH9Z+hgOTXNQOtb2V5rMP/QzDGRyGvcmQ9b/ueP/79phD1w0ee+8k5ydZneQ2ms86eOwvAmcDp6cZlv6hJBt07Ocw4Dqa3u1BOwEb0JyvW9vP8xlgcAbkB5wrHniNjOQD7Wd8c7u8NfBQYGnHMf6jXT9odVXd2bH8gOuxqn5Dcx1v31FnE+4/p5KkOcKEV5I01A+AO4HRHoVzPU3yM2jHdt2YquriqjqYJkn6JnDGYNEw1U8ArgJ2q6pNgXcDGabesIcaRyxnV9VzaXpkrwI+O0y1VcCCwYW2B3LBMPUeJM19zV8HPgxs0w6tXsw4PkObkP+cZmjvuvoysAjYoao2o7mvOe3+76mqv6+q3WmGaB9EM/x50PuBm4AvD96jS5PM3gXMr6rN29emVfW4tnwVTY/zoB1HCy7JocArgZdV1T3t6ptokvbHdRxjs3bo86Chv9sHXI/tcPCtaBL2QY8FfjRaPJKk/mPCK0l6gKq6jeZ+yeOTvDjJQ5Ns0PZUfqit9hXgvUm2bu9xfR/wpbH2nWTDNM+43axNcG6nuU8Xmp7LrdJMmjVok7bOb9re17+awEf5NbAgHZNcDYllmyQvapOju2iGXd87TNXvAE9oz8U84I0095COx4bARsBqYG2ayawm8rijxTx4CPdEbALcUlV3JtkLeNVgQZJnJXlCm8zeTjPEufPz3wMcAvwR8MUk61XVKpoh2R9JsmmS9ZI8MslgjGcAb06yoL1n+5iRAkuyB/CvwIvbYfEAVNV9NF88fCzts3OTbJ/keaN8zi8DRyZ5Uvslw/8FftgOgyfJ9jT3/F40xvmSJPUZE15J0oNU1UeBv6aZiGo1Tc/em2h6ZAH+EVhCM3nSj4FL2nXj8Wrgl+0Q5aNpJlqiqq6iSaRXtENZtwPeQZOk3UGTBH11Ah/jv4ArgBuS3DRM+XrA22l6B2+hSSyHTohFVd1Ek/h9iGaY7O40n/2usQJo7z99M00iuKb9LIsm8BlOBA7rmABrot4AHJfkDpovJc7oKHsEzXDl22kmdLqQIV9aVNXdwEtoeuNPae+3fQ1NIn8lzWc6k6aHHJrf0dk0PamXAGeNEtvBwBbA9ztmav5uW/a3wHLgovY6+U/uvwf8QarqPODvaHrTVwGP5IH3/L4KOLW9v1eSNIfER/tJkjR+bdK3Ejisqs6fgeN9GTijqr45ZmU9SNvj+yPgGe0kW5KkOcSEV5KkMbTDaX9Ic2/p39AMa96lnSFYkiT1KIc0S5I0tqfQTCB1E/BCmvtOTXYlSepx9vBKfSLJYcDhVTWRCXHW5Tj7Al+qqnHNUitJkiR1iz280iyT5GlJ/l+S29rnev5Pkj+pqtOmO9mVJGmuGKm9ncbj/TLJc6Zr/9JcNa/bAUgavySbAv9O82iWM2hmSn0645gtVpIkjU8vtrdJ5lXV2m4dX5qt7OGVZpdHAVTVV6rq3qr6fVWdU1WXJzkiyfcHKybZP8nV7TfT/5bkwiR/2ZYdkeT7ST6cZE2SX7TPBx3c9sgky5LckWRFktfP/EeVJKlrxmpv/yfJv7Zt7FVJ9hvcMMlmSU5OsirJdUn+sX3e9WD56zra2CuT7Jnki8COwLfbR3S9M8nCJJXktUl+BfxXkn2TrOwMtLNnOMn7k3wtyZfa/f84yaOSvCvJjUmuTeJoMM0pJrzS7PJT4N4kpyY5MMkWw1VKMp/m2ZjvArYCrgaeOqTa3u36+TTPFz2541mfNwIHAZsCRwIfS7LnVH8YSZJ61Fjt7d7ACpo29FjgrCRbtmWnAmuBXYE9gP2BwS+cDwHeT/M8602BFwE3V9WrgV8BL6yqh1XVhzqO9UzgscDzxhn7C4Ev0jzn+lKaZ2OvB2wPHAd8Zpz7kfqCCa80i1TV7cDTgAI+C6xOsijJNkOqPh+4oqrOaoc/fRK4YUida6rqs1V1L03jvC2wTXuc71TVz6txIXAOzVAuSZL63jja2xuBj1fVPVX1VZovkF/Qlh8IvLWqfts++/ljwKHtdn8JfKiqLm7b2OVVdc0Y4by/3dd4Z4b/76o6u23/vwZsDXywqu4
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "signal_constant = np.ones(n_samples)\n",
    "\n",
    "fig, axes = plot_signal_and_spectrum(signal_constant, sample_rate, \"Constant signal\")\n",
    "axes.flat[1].set_xlim(0,10);\n",
    "\n",
    "fig, axes = plot_signal_and_spectrum(signal_constant, sample_rate, \"Constant signal (masked zero)\", ft_kwargs={'mask_bias':True})\n",
    "axes.flat[1].set_xlim(0,10);\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### FTs: noise (normal distributed)\n",
    "\n",
    "A FT of normal distributed noise will not have easily identifiable peaks in the spectrum.\n",
    "Instead, the spectrum is on average a bit flat, with only low power in any individual peak."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA64AAAEjCAYAAAAyt1xdAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOydeZwcVbn3f0/3bJmsk4UsJJkQCInsMCBBRXGXK4oXZI0oXoLb9arX7VWvoqLe61UQVLhAAqJgwo7sGkDWkEyWgYQkZJ9kksmeSSeZySy91PP+UXWqT1VXdVfvPZPn+/lM0l1ddepU1amq8zvPcoiZIQiCIAiCIAiCIAiVSqjcFRAEQRAEQRAEQRCEdIhwFQRBEARBEARBECoaEa6CIAiCIAiCIAhCRSPCVRAEQRAEQRAEQahoRLgKgiAIgiAIgiAIFY0IV0EQBEEQBEEQBKGiEeEqCIIglAwi+jsRfaEI5X6ZiG4tdLnFgIiYiE4IuO7PiOiv1ufJRNRFROEC1eNOIvqJ9fkCImovRLk51uXPRPRL6/NpRLSoXHURBEEQKhMRroIgCEJgiGgrEe0hosHastlE9EqQ7Zn5Qmb+S4HrVAPgxwB+W8hyKw1m3sbMQ5g5kW49IrqWiBYGKO8rzPyLQtQtGzGeCWZ+G8BBIvpUIcoTBEEQBgYiXAVBEIRsqQLwzXJXQuNiAOuYeUehCyaiqkKXWQkUympbROYB+HK5KyEIgiBUDiJcBUEQhGz5LYDvEtEIrx+J6D1EtIyIDln/v0f77RUimm19PoGIXrXW209ED2nrzSCiF4joABGtJ6LL09TnQgCvattOsSyAXyCibVbZ/6X9XktEtxLRTuvvViKqtX67gIjaiej/EdFuAPdqy75PRHuJaBcRfYaI/oWINlh1/JFW/ruJaDERHbTWvc2yCmeEiI6zzkknEb0AYLTHcVVZ368lolZr3S1ENIuI3gXgTgDnWW7FB611/0xEdxDRc0R0BMAHdfdcbR8/ss7XViKa5XXdtH0vtD6/Zi1eae3zCmv5RUS0wjoPi4joNG37M4noTavuDwGoc52KVwB8WF0XQRAEQRDhKgiCIGTLcpjC4rvuH4hoJIBnAfwBwCgAvwPwLBGN8ijnFwCeB9AAYCKAP1plDAbwAoD5AI4BcBWA/yOik33qcyqA9R7L3wdgOoAPA7jBEnUA8F8AZgI4A8DpAN4N09VYMQ7ASACNAL6kLasDcCyAGwDMBfA5AE0AzrfKn2qtmwDwnzBF53nW/r/mU3c38wG0WNv+AoBnPLB1jv4A4EJmHgrgPQBWMPNaAF8BsNhyK9YHF64G8CsAQwF4uRKPs/Z7rLXfOUQ0PVOFmfn91sfTrX0+RERnAfgTTKvpKAB3AXjKGjSoAfAEgPthnudHAFzqKnMHgBjM6ycIgiAIIlwFQRCEnLgBwH8Q0RjX8k8C2MjM9zNznJkfALAOgFe8YgymOJzAzL3MrMTURQC2MvO9VhlvAngMwGd96jICQKfH8p8zcw8zrwSwEqZIBYBZAG5k5r3MvA/AzwFco21nAPgpM/cxc49W118xcwzAgzAF3u+ZuZOZ1wBYA+A0AGDmFmZutuq+FaZo+4BP3W2IaDKAcwD8xNr3awCeTrOJAeAUIhrEzLuseqTjSWZ+g5kNZu71WUft+1WYAxDpLN3puB7AXcy8hJkTVlxzH8wBg5kAqgHcyswxZn4UwDKPMjphXltBEARBEOEqCIIgZA8zrwbwDIAfuH6aAKDNtawNphXPzfcBEIClRLSGiP7NWt4I4FzLxfSg5e46C6ZF0IsITCuim93a524AQ3zq2GYtU+zzEHYdWlIkJWb3aL/3qPKJ6EQieoaIdhPRYQD/Dc3lNw0TAESY+YirbilY61wB07q6i4ieJaIZGcrfnuF3r31P8Fs5A40AvuO6hpOs8iYA2MHM7NqXm6EADua4f0EQBGGAIcJVEARByJWfwrSs6aJ0J0zRojMZQEriJGbezczXM/MEmC6l/2dlpt0O4FVmHqH9DWHmr/rU420AJ2ZRb3cdJ1vL7KplUZYXd8C0Mk9j5mEAfgRToGdiF4AG0jI2W3XzhJkXMPNHAYy39jdX/eS3SYb9e+1bnZcjAOq13/wGERTbYVqo9WtYb1ngdwE4loj0c+I4TiKaAKAG3i7ggiAIwlGICFdBEAQhJ5h5E4CHAHxDW/wcgBOJ6GoiqrIS9ZwE0zrrgIguI6KJ1tcITGGVsNY9kYiuIaJq6+8cLUbVzXMI4Iqr8QCAHxPRGCIaDdPt+a9ZbJ+JoQAOA+iyrKB+gtsBM7fBjB/+ORHVENH74O1iDSIaS0SftoRmH4AumOcOMC3BE4MmhHKh9n0+TJftR6zlKwBcQkT11uDCda7t9gCYqn2fC+ArRHQumQwmok8S0VAAiwHEAXzDaiOXwIwz1rkAwEvM3JfDMQiCIAgDEBGugiAIQj7cCMC20jFzB0zB8x0AHTDdgS9i5v0e254DYAkRdQF4CsA3mXkLM3cC+BiAK2Fa/HYD+F8AfhlmnwYww7LSBeGXMAXi2wBWAXjTWlYovgszEVInTAH3UPrVHVwN4FwAB2BatO/zWS8E8xzvtNb9AJIJoF6CGXO7m4i8zrsfu2EOIOyEOR3NV5h5nfXbLQCiMAXqX6zfdX4G4C+WW/DlzLwcpjX+NqvMTQCuBQBmjgK4xPoegeny/LirvFkwsyMLgiAIAgCAnCEmgiAIgtD/IKIvATiJmb9V7roI+UFEpwKYw8znlbsugiAIQuUgwlUQBEEQBEEQBEGoaMRVWBAEQRAEQRAEQahoRLgKQj+BiGYR0fMl2M8FRNRe7P0IgiAIgiAIQlBEuApChUFE7yOiRUR0iIgOENEbRHQOM89j5o+Vu36CIAiCMFDwe+cWcX9biegjxSpfEAYyVeWugCAISYhoGMypQL4K4GGY8xieD3O6C0EQBEEQCkQlvnOJqIqZ4+XavyBUMmJxFYTK4kQAYOYHmDnBzD3M/Dwzv01E1xLRQrUiEX2MiNZbo8T/R0SvEtFs67driWghEd1ERBEi2kJEF2rbfpGI1hJRJxG1EtGXS3+ogiAIglBWMr1z3yCiP1rv2XVE9GG1IRENJ6J7iGgXEe0gol8SUVj7/XrtPfsOEZ1FRPcDmAzgaSLqIqLvE9EUImIiuo6ItgF4yStkR7fUEtHPiOgRIvqrVf4qIjqRiH5IRHuJaDsRiYeWMOAQ4SoIlcUGAAki+gsRXUhEDV4rEdFoAI8C+CGAUQDWA3iPa7VzreWjAfwGwD1ERNZve2HOtTkMwBcB3EJEZxX6YARBEAShgsn0zj0XQCvM9+hPATxORCOt3/4CIA7gBABnwpx7Wg0eXwZzbuPPw3zPfhpABzNfA2AbgE8x8xBm/o22rw8AeBeAjwes+6cA3A+gAcBbABbA7NcfC3N+7bsCliMI/QYRroJQQTDzYQDvA8AA5gLYR0RPEdFY16r/AmANMz9uuRT9AcBu1zptzDyXmRMwX7DjAYy19vMsM29mk1cBPA/TPUoQBEEQjgoCvHP3AriVmWPM/BDMweBPWr9fCOBbzHyEmfcCuAXAldZ2swH8hpmXWe/ZTczclqE6P7PK6glY/deZeYHVB3gEwBgAv2bmGIAHAUwhohEByxKEfoEIV0GoMJh5LTNfy8wTAZwCYAKAW12rTQCwXduGAbgzAe/Wfu+2Pg4BAGtkudlKRHEQphAeXdgjEQRBEITKJsM7d4f1flW0Wb83AqgGsIuIDlrv0bsAHGOtNwnA5iyrsj3zKg72aJ97AOy3BqrVd8B65wvCQEGEqyBUMMy8DsCfYb5MdXYBmKi+WC7AExEAIqoF8BiAmwCMZeYRAJ4DQGk3FARBEIQBjMc791gtxAYw41N3whSZfQBGM/MI628YM59srbcdwPF+uwmw/AiAevXFip0dk82xCMJARISrIFQQRDSDiL5DRBOt75MAXAWg2bXqswBOJaLPEFEVgH8HMC7gbmoA1ALYByBuJW2SJA6CIAjCUUWAd+4xAL5BRNVW3Oq7ADzHzLtghtjcTETDiCh
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "signal_normal_noise = rng.normal(size=n_samples)\n",
    "\n",
    "plot_signal_and_spectrum(signal_normal_noise, sample_rate, \"Noise (normal distributed)\");"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### FTs: sine wave"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA74AAAEjCAYAAAAVGd21AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOx9eZxdRZX/93QSQFkkEBbZsij7EjYRfuIAjmJYBHGZAR0RR2XUQRzXwXFD1FEZHJSRQdlEAcMmICAEGCFsEkICISQhgezp7Pu+dL97fn/ce9+rOrXcezvd/Ui/8/18ku73qurU91TVvX3P/dZCzAyFQqFQKBQKhUKhUCj6KtqaTUChUCgUCoVCoVAoFIqehAa+CoVCoVAoFAqFQqHo09DAV6FQKBQKhUKhUCgUfRoa+CoUCoVCoVAoFAqFok9DA1+FQqFQKBQKhUKhUPRpaOCrUCgUCoVCoVAoFIo+DQ18FQqFQqFQKBQKhULRp6GBr0KhUCgUCoVCoVAo+jQ08FUoFApFIYhoMhGd2k22ZhPR+7tQ7mAiepmI1hLRpd3BpdVARHsQ0eNEtJKIbmo2H4VCoVAoegsa+CoUCoUCAEBEJxPR34hoNRGtIKLniOhdAMDMhzPz6CZT/BaA0cy8MzNf02QubzqUfKHwbQBvMPNAZv7sVtS1GxHdR0TriWgOEX2iq/mLbJVIv4SIxhHRZiK6RaTNJqItRDRIfD+BiJiIhoi87xf5LiKiZ0s2i0KhUCjexNDAV6FQKBQgol0APATgfwDsBmBfAD8EsLmZvAQGA5jsSyCi/r3MZVvF+wHc3Q12rgWwBcBeAD4J4DoiOryL+YtsFaUvAPBjADcH6p4F4IL8AxEdCeAtJXxUKBQKRR+CBr4KhUKhAICDAICZRzJzjZk3MvNjzDwRsNWw7PdvENHETB2+k4h2yA0R0bHGlOS7s/Qf+yolon2I6E9EtJSIZoWmMBPREwBOA/BrIlpHRAdlPP6diCYCWE9E/YvsEdExRPRSxu1OIroj55YpgO808t5i8o7ZLtEm+xPRvVnZ5UT06+z7bxLRnwTH/yGiX3ra4DIimpFxn0JE5xlptwI4AMCDWft8S5TdjohWAzgyy/Oqr53LgIh2BPBRAN9j5nXM/CyABwB8qmr+Iltl6mLme5n5fgDLA5RvBXCh8fnTAP7QBb//MWvb/N9mIhpd1Y5CoVAomgMNfBUKhUIBAK8DqBHR74noDCIaWJD/HwCMADAUwFEALgLSAAvAfQBuQaocjwRwns8AEbUBeBDAK0gV5r8H8G9E9EGZl5nfB+AZAJcw807M/HqWdAGAswDsCiCJ2cu43Y80ENoNqfL50QI/q3ANtUk/pGr6HABDsvJ3ZGVuAzCCiHbN8vYH8I8ZR4kZAN4L4G1I1fjbiOjtWft8CsBcAB/K2udK0X5bAJwEYEmWfqTw7yEiWhX495DgcRCAmtEHyNolpPjG8hfZqlqXD2MA7EJEh2Z98Y9I270SmPnOrO12ArAPgJlIx7dCoVAotgFo4KtQKBQKMPMaACcDYAA3AFhKRA8Q0V6BItcw8wJmXoE0IDw6+/5EAP2z9A5mvhfA2ICNdwHYg5mvYOYtzDwzq/v8CtSvYeZ5zLyxhL0TAQwA8MuM2z0AXixZTxmuoTY5AWmg9E1mXs/MmzLlEsy8EMDTAD6e5R0BYBkzj5cEmPnuzH7CzHcCeCOzXRZHIw0aHTDz2cy8a+Df2SL7TgBWi+9WA9g5UG8sf5GtqnWFkKu+HwAwFcD8QL77zaAfwP/KDNlLkD8iXW/+24o8FAqFQtEk6JoohUKhUAAAmPk1NFTKQ5CqYr+EsT7SwCLj9w1IAztkP+czMxvp8wJVDgawTxZg5OiHVNktC9N2kT0ftzkl6ynDNdQm+wOYw8ydAdu/B/BFpIH0P8Gv9oKILgTwNaSqMZAGhYN8eQMIBr4VsQ7ALuK7XQCs7UL+IltV6wrhVqQvGIYiPs35w8z8f/kHIroIwOdEnp8gDbx1Z3GFQqHYhqCKr0KhUCgcMPNUpNOVj6hYdCGAfYmIjO/2D+SdB2CWUBd3ZuYzq1CtYM/H7QDj9w0A3mp83rubuM4DcACFN+C6H8BRRHQEgLMB3C4zENFgpIHxJQB2Z+ZdAUwCYPrCspzAcAQCXyJ6RKxfNf89IrK/DqA/ER0obHs3HivIX2Sral1eMPMcpJtcnQng3iplTRDR+UhfBH2MmTu6akehUCgUvQ8NfBUKhUIBIjqEiL5ORPtln/dH+oA/pqKp5wHUAFySbTZ1LsLTcccCWJNtUPUWIupHREdQdoRSF1Bk73kAnQAuzbh9RHCbAOATWbkRAE7pJq5jkQbdPyOiHYloByJ6T57IzJsA3IN0+uxYZp7rsbEj0sB2KQAQ0WfgvpRYDGBYhEcw8GXmM/L1q55/Z4i865EGj1dk/rwHwLkIKNWx/EW2ytSV9eUOSBX4fln7+l4yfBbA+zKblUFExyDd9fzDzLy0KzYUCoVC0Txo4KtQKBQKIJ06+m4ALxDReqQB7yQAX69iJNtE6SNIg4xVSKfuPgTPsUjMXAPwIaRTcGcBWAbgRqSbN1VGkT2D20UAViLd5MhU/76SlV+F9Nic+7uDq1H2nUg3oGrP6jbxe6Q7LoeCxykAfoE0eF+c5X1OZPspgO9m61O/YSYQ0d4ABiJd39od+BLSI4GWIN3g6YvMXFdhMwX5P0rmj9oqkf5dABsBXIZ0vG3MvrPAzDOYeVyXPU4D7oEAno2o4QqFQqF4k4LspU4KhUKhUHQviOgFAL9h5t81m4sEEd0CoJ2ZnUCpl3kcgDQo3TvbaEyhUCgUCkU3QhVfhaJFQESfJKLHeqGeU4movafrUbx5QUSnENHe2RTUTyM92mdUs3m9WZHtEvw1AHdo0KtQKBQKRc9AA1+Foo+BiE4mor8R0WoiWkFEzxHRu5j5dmY+vdn8FC2Bg5GuJV2NdKr0x7JjexQCRLQjgDVIj9n5QZPpKBSKCgj9ve3B+mYT0ft7yr5C0dehxxkpFH0IRLQL0vWUXwRwF4DtALwXnvWVCkVPgZmvB3B9s3mUATNf1OT61yM9lkihUGxDeDP+vSWi/pFj0xSKlocqvgpF38JBAMDMI5m5xswbmfkxZp5IRBcR0bN5RiI6nYimZW+q/5eIniKiz2VpFxHRs0R0FRGtJKJZRHSGUfYzRPQaEa0loplE9C+976pCoVAoFE1D0d/b54jof7K/sVOJ6O/zgkT0NiK6iYgWEtF8IvoxEfUz0j9v/I2dQkTHEtGtSI9fezDbWO1bRDSEiJiIPktEcwE84VtuZCrFRHQ5Ed1NRLdl9l8looOI6NtEtISI5hGRzg5T9Elo4KtQ9C28DqBGRL8nojOIaKAvExENQnp8yrcB7A5gGoD/J7K9O/t+EIArAdxEVD//dAnS80Z3AfAZAFcT0bHd7YxCoVAoFG9SFP29fTeAmUj/hv4AwL1EtFuW9nukR6u9E8AxAE4HkL94/jiAywFciPRv7DkAljPzp5DuCv+h7JixK426TgFwKIAPluT+IaQ7yA8E8DKAR5HGBPsCuALAb0vaUSi2KWjgq1D0IWQb45yM9LzPGwAsJaIHiGgvkfVMAJOZ+d5sWtQ1ABaJPHOY+YbsKJbfA3g7gL2yev6SHQ3CzPwUgMeQTvFSKBQKhaLPo8Tf2yUAfsnMHcx8J9IXyWdl6WcA+DdmXs/MSwBcDeD8rNznAFzJzC9mf2OnM/OcAjqXZ7Y2lqT/DDM/mv39vxvAHgB+xswdAO4AMISIdi1pS6HYZqCBr0LRx8DMrzHzRcy8H4AjAOwD4Jci2z4A5hllGOnZoiYWGekbsl93AoDs7faYbDOPVUgD6UHd64lCoVAoFG9eFPy9nc/2maFzsvTBAAYAWJidub0KqcK6Z5ZvfwAzKlKZV5zFwmLj940AlmUvufPPgO49oOiD0MBXoejDYOapAG5B+gfZxEIA++UfsinM+6EEiGh7AH8CcBWAvZh5VwAPA6BoQYVCoVAo+ig8f2/3NZYHAen63AVIg9TNAAYx867Zv12
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "f = 1e2 # Hz\n",
    "if f > sample_rate/2:\n",
    "    print(\"Sampling a frequency above the sample_rate/2 gives problems\")\n",
    "\n",
    "signal_sine = np.sin( 2*np.pi*f*t )\n",
    "if False:\n",
    "    # aliased peak\n",
    "    signal_sine += 1/2*np.sin( 2*np.pi* (1.1* sample_rate/2) *t)\n",
    "\n",
    "fig, axes = plot_signal_and_spectrum(signal_sine, sample_rate, \"Single frequency at $f = {:g}$MHz\".format(f/10**6))\n",
    "axes.flat[1].axvline(sample_rate/2, color='r', label='Nyquist Frequency');\n",
    "axes.flat[1].legend();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A sine wave will give very clear peaks at the frequency of the sine wave."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### FTs: step function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA7sAAAEjCAYAAADzMBbxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3yV9fn/8deVAWETIGySMEWGgiDgwD1wYbXVuqrWqtXWrp/ja5el2tpdW7V1W7Wi1l20WLWgiANkKChDCCNsQtgbkly/P+478RCygJzc5yTv5+NxHjnn3J9zn+u+75Nzn+v+LHN3REREREREROqTlKgDEBEREREREaltSnZFRERERESk3lGyKyIiIiIiIvWOkl0RERERERGpd5TsioiIiIiISL2jZFdERERERETqHSW7IiIiIiIiUu8o2RUREREREZF6R8muiIjUOTM7zMw+MbOtZvb9OnzfOWZ2Ul29n4iIiERHya6ISANlZseb2YdmttnMNpjZB2Z2dLhsqZmdFse3vw14191buPu98XqT8tvh7v3d/d04vddKMxsUj3WH629sZo+ZWX54keATMzurivLvmtkuM9sW3r6ooMwlZjbPzLab2SIzGxmv+MP3a2Nmr4Tvl29ml5Vbnmtm481so5mtMbP7zSwtnjGJiEj9pWRXRKQBMrOWwOvAfUAboAvwS2B3HYWQA8ypo/eKOzNrB7QH5sXxbdKA5cCJQCvg58DzZpZbxWtucvfm4e2w2AVmdjrwO+CbQAvgBGDxwQZnZmPMbEw1xf4G7AE6AJcDD5hZ/5jlfwcKgE7AIIJt/c7BxiQiIg2bkl0RkYapD4C7P+vuxe6+093fcvfZZvZPIBt4LawRvA3AzDqb2Utmts7MlsQ2Pw5rUH9sZnPDWrl/mFlGRW9sZhOBk4H7w/X3MTM3s14xZZ4ws1/FrPsWM5sd1kL/K3bdZtbNzF4O41pvZveHz++3HeVres3s8LAGdFPYxHl0uW2q9H1jyvUiSEJTgPVhDLVeG+nu2919jLsvdfcSd38dWAIMOchV/hK4092nhOtb6e4rSxea2XXh8dxsZm+YWftDid/MmgFfBX7u7tvc/X1gHPCNmGLdgefdfZe7rwH+C/SPWccNZvYfM/ubmRWa2aowacfMepjZ6+Hzm83s7UOJV0REkp+SXRGRhmkBUGxmT5rZWWaWWbrA3b8BLAPOC2sEf29mKcBrwCyCWuBTgR+a2Zkx67wcOBPoSZBM/6yiN3b3U4DJfFnruKAG8V4MjCJIho4ArgYws1SCGup8IDeM7bnKtiN2hWaWHm7TWwS1st8DxppZbA1ohe9bbnvygFuAF8P3aevuRdVtUJiYbark9noNXt+BYD9XVUP+mzD5+8Bi+iqH+20okGVmeWa2Imwy3CRc/hPgBmA0kAWsBH5VXUzV6AMUlzves4hJZoG/ApeYWVMz6wKcRZDwljoCOIYgSW4PPAT8X7jsKeANglrjDsCYQ4xXRESSnJJdEZEGyN23AMcDDjwCrDOzcWECVZGjgSx3v9Pd97j74vB1l8SUud/dl7v7BuDXwKW1GPK97r4qXPdrBE1cAYYBnYFbw5rPXWGNYU2MAJoDvw23aSJB4hwbd2XvW96RwKexT5jZr81sspm9aGZNy7/A3c9199aV3M6tKvAwUR8LPOnu8ysp9n9AD4ILAA8T1HD3DJd1ANKBrwEjw+0aDPwsrMH9GXCpu+e5+x7gMYLPwKFoDmwu99xmgibUpSYRJL9bgBXAdODVmOVHEByvN929BJgbs6wnkAqkhp+DDw4xXhERSXJKdkVEGih3n+fuV7t7V2AAQdL4l0qK5wCdY2sfgZ8QJE2llsfczw/XV1vWxNzfQZA4AXQD8mtSk1qBzsDyMGkqlU+QHFb3vuUNIqilBMDMBgA93X0k8D/gmoOIr0JhLfs/Cfq+3lRZOXef6u5b3X23uz8JfACcHS7eGf69z91Xu3sh8Odw+alAI+DjmGP9X/ZPVPepnQZuB26vonZ6G9Cy3HMtga0x2/Um8DLQDGgHZBL0Ky41kOCiQ6kBfJnwXg6cD6yyYCCvNpXtGxERaRiU7IqICGHt4BMEyQMENb6xlgNLytU+tnD3s2PKdIu5nw2sOoAQdgCxtZ8da/i65UB2FX1ky29HrFVAtzDJKpVN0GS3xsLXD2Dfmt2RBE1qCf8eX8Hr3rAvR0ouf3ujfPnwNUZQy9oB+Kq77z2AUB0wAHffSFBzWtH+aQO8Uu5Yt3L3k/ZbYUztNPBbglrXymqnFwBpZtY75rkj+bIZdhuCz9D9YYK+HvgHYYJuZt0JBumKHVV6MOF+d/eJ7n4q0C9c79U12isiIlJvKdkVEWmAzKyvmd1sZl3Dx90Imu9OCYusJWgCW+pjYIuZ/Z+ZNTGzVDMbYOFURaHvmlnXsEbtJ8C/DiCkT4HLwvWOIhiFtyY+BlYDvzWzZmaWYWbHxSwvvx2xpgLbgdvMLD3s03oeYZ/fA9AkvMWeUzP5siZ0M0Eitw93PytmpOTyt8qmFHoAOJygH/LOSspgZq3N7Mxwf6SZ2eUEoy2/GVPsH8D3zKx92Gf7hwTNuGcCJ5vZUeG6WprZ+WGifdDcfTtBre2d4bE6jqAm9p/h8kKCAbduDGNuDVzFlzXmRwCflauJHwzMMrMLzax3GGMLgv2/T7NyERFpeJTsiog0TFuB4cBUM9tOkOR+DtwcLv8NQf/NTWZ2i7sXEySCgwgSkkLgUYIpcEo9QzDY0+LwdiADGv0gXP8mguaor1ZdPBATVy+CwahWAF+PKbLPdpR77R6CAZjOCrfn78CVVfSBrSyG7cCDwFwzWxE+vZEv900rYMOBrLMiZpYDfJvgGKyJqQW+PFz+RjiwFAT9cX8FrCPYtu8BX3H32FrRu4BpBDWu84BPgF+7+0fAncBLZraNoJnwKHevqpa8pr5DcGGgAHgWuNHdYwfYupBgQLB1QB5QBPwoXHYEMQmsmbUlaAHwOUHN+SSCz/V4ghrmibUQr4iIJDGrnXOXiIg0ZGa2FLjW3f8XdSyJwMwGAj9298vM7HqgsbvfF3VcIiIiDYlqdkWSiJldbmZv1cH7nBRTQyUiB8jdPwPyzWwywXRMj0cckoiISIOjZFckAZnZ8Wb2oZltNrMN4RyZR7v7WHc/I+r4RKR67v5jdx/p7l8NmzqLSAKq7Jwbx/dbamanxWv9IvKlykavFJGImFlLgkFibgSeJ5gCZCSwO8q4RKri7rlRxyAicqAS8ZxrZmkHOZ2aiJSjml2RxNMHwN2fdfdid9/p7m+5+2wzu9rM3i8taGZnmNkX4dXov5vZJDO7Nlx2tZm9b2Z/NLONZrbEzM6Kee03zWyemW01s8Vm9u2631QREZFIVXfO/cDM7gvPs/PN7NTSF5pZq3BO59VmttLMfmVmqTHLr4s5z841s6PM7J8EU5y9Fg4wd5uZ5ZqZm9m3zGwZMLGi7kSxNcJmNsbMXjCzp8P1f2Zmfczsx2ZWYGbLzUwtwaTBU7IrkngWAMVm9qSZnRVOCbIfM2sHvAj8GGhLMPfkseWKDQ+fbwf8HngsZvqQAuBcoCXwTeCe0qlGREREGojqzrnDCUaXbwf8AnjZgunVAJ4kGDG8F8E0WGcApRecLwLGAFcSnGdHA+vd/RsEI8efF04z9vuY9zqRYGqxM2sY+3kEU3dlEoym/ibBb/suBCOqP1TD9YjUW0p2RRKMu28hmEbDgUeAdWY2zsw6lCt6NjDH3V8OmzvdC6wpVybf3R8Jp2d5EugEdAjf5z/uvsgDkwimjBkZvy0TERFJLDU45xYAf3H3ve7+L4ILyOeEy88Cfuju2929ALgHuCR83bXA7919WniezXP3/GrCGROuq9I5tMuZ7O5vhr8BXgCyCKbd2kswX3huOF+1SIOlZFckAbn7PHe/2t27AgOAzsBfyhXrDCyPeY0TzDEaa03M8h3h3eYA4RXsKeFgHJsIkud2tbslIiIiia2ac+7KcnNM54fLcwjms14dzuO9iaAmtX1Yrhuw6ABDWV59kX2sjbm/EygML26XPob
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA7sAAAEjCAYAAADzMBbxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd3gc5dX38e9RccNFsuWCi9yNcW8YE0yA0GyKTWihBwIYCCQkD4EQHh5CSPImIYXQQg0lxKFDMMSEYnqxcTe4y0W23GS5yHK3pPP+MSOzltWt1a5Wv8917aXdnXtmzuysdvbs3czdEREREREREUkkSbEOQERERERERKS2KdkVERERERGRhKNkV0RERERERBKOkl0RERERERFJOEp2RUREREREJOEo2RUREREREZGEo2RXREREREREEo6SXREREREREUk4SnZFRKTOmdkRZjbbzArM7Md1uN/5ZnZCXe1PREREYkfJrohIA2Vmo83sczPLN7PNZvaZmR0VLltpZidHcfe3Ah+6ewt3vz9aOyl9HO7e390/jNK+1pjZkGhsO2IfH5rZbjPbHt4WV1D2RjObYWZ7zOzpUssam9nfzSw7/MFhtpmNjWbs4X5bm9lrZrYj3PfFpZZ3M7PJZrbFzNab2YNmlhLtuEREJDEp2RURaYDMrCXwJvAA0BroBPwK2FNHIXQF5tfRvqLOzDKAdsDCOtjdje7ePLwdUUG5tcBvgCfLWJYCrAaOB1oB/we8aGbdahqUmd1lZndVUuwhYC/QHrgEeNjM+kcs/xuQCxwODAnj+2FNYxIRkYZNya6ISMPUB8Ddn3P3Inff5e7vuPs8M3sWyATeCGsPbwUws45m9oqZbTSzFZHNj8Ma1F+Y2YKwVu4pM2tS1o7N7H3gRODBcPt9zMzNrFdEmafN7DcR2/6Zmc0La6FfiNy2mXUxs1fDuDaZ2YPh8wcdR+maXjM7Mqwt3Ro2cR5X6pjK3W9EuV4EiWMSsCmMIea1ke7+qrv/G9hUxrId7n6Xu69092J3fxNYAQwvKWNm14TnM9/M3jKzdocSj5kdBpwL/J+7b3f3T4FJwGURxboDL7r7bndfD/wX6B+xjevM7D9m9pCZ5ZnZWjM7JVzWw8zeDJ/PN7N3DyVeERGp/5Tsiog0TEuAIjN7xszGmll6yQJ3vwxYBZwV1h7eY2ZJwBvAXIJa4JOAn5jZaRHbvAQ4DehJkEzfUdaO3f07wCd8U0O5pArxXgCMIUiGBgFXAJhZMkENdTbQLYzt+fKOI3KDZpYaHtM7BLWyPwImmllkbWmZ+y11PFnAz4CXw/20cffCyg4oTMy2lnN7s4JVfxcmdJ9ZLfU/NrP2BOdsfvj4duA6YBzQFlhDUEt8KPoARaXO91wiklngPuBCM2tmZp2AsQQJb4lBwDEESXI74FHg5+GyfwBvEdQatwfuOsR4RUSknlOyKyLSALn7NmA04MDjwEYzmxQmPWU5Cmjr7ne7+153Xx6ud2FEmQfdfbW7bwZ+C1xUiyHf7+5rw22/QdDEFWAk0BG4Jayt3B3WGFbFKKA58PvwmN4nSJwj4y5vv6UNBuZEPmFmvzWzT8zsZTNrVnoFdz/T3dPKuZ1Zzn5+DvQgSOofI6i17lnF4y1TmPRPBJ5x90VhDe4dwEXunuXue4G/E7wHDkVzIL/Uc/lAi4jHHxEkv9uAHGAG8O+I5YMIztfb7l4MLIhY1hNIBpLD98FnhxiviIjUc0p2RUQaKHdf6O5XuHtnYABB0vjXcop3BTpG1j4CtxPUoJVYHXE/O9xebVkfcX8nQeIE0AXIrkpNahk6AqvDpKlENkEiWdl+SxtCUEsJgJkNAHq6+3HAe8APahDfQdx9mrsXuPsed38G+Aw4vabbC2vsnyXoR3tj+PRJQCPgy4hz/V8OTlQPqJ0GbgNuq6B2ejvQstRzLYGCiFjeBl4FDgMygHTgDxHlBxL86FBiAN8kvJcA44G1Fgy+1bqKL4OIiCQoJbsiIoK7LwKeJkgeIKjxjbQaWFGq9rGFu0cmWl0i7mcSDJBUVTuByNrPDlVcbzWQWUEf2dLHEWkt0CVMskpkEjTZrbJw/QEcWLN7HEGTWsK/o8tY7y37ZlTl0re3SpcvhwNWnXgj9m8ENbbtgXPdfV+4qDXwWqlz3crdTzho5xG108DvCWpdy6udXgKkmFnviOcG881AZa0J3kMPhsn8JuApwmTezLoTDKwVOQL1UMLX3d3fd/eTgH7hdq+owcsiIiIJRMmuiEgDZGZ9zexmM+scPu5C0Hx3alhkA0Fz2RJfAtvM7Odm1tTMks1sgIVTFYVuMLPOYY3a7cAL1QhpDnBxuN0xBKPwVsWXwDrg92Z2mJk1MbNjI5aXPo5I04AdwK1mlhr2fz2LsM9vNTQNb5HX1HS+qQnNJ0jkDuDuYyNGVS59O2gaIDNLM7PTwmNMMbNLgG8T1IYeJCzThLBpb8l6EUUeBo4k6NO8K+L5WcCJZjYs3E5LMxsfJsc15u47CGpt7w7P1bEENbHPhsvzCAbJuj6MPQ34Pt/UmA8CvipVEz8UmGtm55hZ7zDGFgSv/wHNykVEpOFRsisi0jAVAEcD08xsB0GS+zVwc7j8d8AdYXPUn7l7EUEiOIQgIckDniCYtqbEvwgGe1oe3qozoNFN4fa3EjRH/XfFxQMRcfUiGIwqB/heRJEDjqPUunsJBmAaGx7P34DLw1ruKguTuEeABWaWEz69hW9em1bA5upssxypBK/pxjDeHwFnu/ti2F9TfHtE+TuAXQTNiy8N798Rlu0KXEtwPtdH1Chf4u5fAHcDr5jZdoJmwmPcvaJa8qr6IcEPA7nAc8D17h45BdU5BAOCbQSygELgp+GyQUQksGbWhqAFwNcENecfEbyvJxPUML9fC/GKiEg9ZrVz7RIRkYbMzFYCV7v7e7GOJR6Y2UDgF+5+sZlNABq7+wOxjktERKQhUc2uSD1iZpeY2Tt1sJ8TImqoRKSa3P0rINvMPiGYjunJGIckIiLS4CjZFYlDZjbazD43s3wz2xzOp3mUu09091NjHZ+IVM7df+Hux7n7uWFTZxGJQ+Vdc6O4v5VmdnK0ti8i3yhv9EoRiREza0kw1+f1wIsEU4AcB+yJZVwiFXH3brGOQUSkuuLxmmtmKTWcTk1ESlHNrkj86QPg7s+5e5G773L3d9x9npldYWaflhQ0s1PNbHH4a/TfzOwjM7s6XHaFmX1qZn8ysy1mtsLMxkase6WZLTSzAjNbbmbX1v2hioiIxFRl19zPzOyB8Dq7yMxOKlnRzFqFczqvM7M1ZvYbM0uOWH5NxHV2gZkNM7NnCaY4eyMcFO5WM+tmZm5mV5nZKuD9sroTRdYIm9ldZvaSmf0z3P5XZtbHzH5hZrlmttrM1BJMGjwluyLxZwlQZGbPmNlYM0svq5CZZQAvA78A2hDMPfmtUsWODp/PAO4B/h4xfUgucCbQErgSuLdkqhEREZEGorJr7tEEo8tnAL8EXrVgejWAZwhGDO9FMA3WqUDJD87nA3cBlxNcZ8cBm9z9MoKR488Kpxm7J2JfxxNMB3ZaFWM/i2DqrnRgNsE0ZElAJ4IR1R+t4nZEEpaSXZE44+7bCKbRcOBxYKOZTTKz9qWKng7Md/dXw+ZO9wPrS5XJdvfHw+lZngEOB9qH+/mPuy/zwEcEU8YcF70jExERiS9VuObmAn91933u/gLBD8hnhMvHAj9x9x3ungvcC1wYrnc1cI+7Tw+vs1nunl1JOHeF29pVSbkSn7j72+F3gJeAtgTTbu0jmC+8WzhftUiDpWRXJA65+0J3v8LdOwMDgI7AX0sV6wisjljHCeYYjbQ+YvnO8G5zgPAX7KnhYBxbCZLnjNo9EhERkfhWyTV3Tak5prPD5V0J5r5eF87jvZWgJrVdWK4LsKyaoayuvMgBNkTc3wXkhT9ulzyG8Jov0lAp2RWJc+6+CHia4AIcaR3QueRB2Dy5M1VgZo2BV4A/Ae3dPQ2YDFiFK4qIiCSwMq65nSK6/0DQ33YtQWK6B8hw97Tw1tLd+4flVgM
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "t_0 = t[n_samples//4]\n",
    "signal_step = np.heaviside(t - t_0, 1)\n",
    "fig, axes = plot_signal_and_spectrum(signal_step, sample_rate, \"Step function $t_0 = {:g}ns$\".format(t_0 * 10**9), ft_kwargs={'mask_bias':True});\n",
    "axes.flat[1].set_xlim(0, 20);\n",
    "\n",
    "t_0 = t[n_samples//2]\n",
    "signal_step = np.heaviside(t - t_0, 1)\n",
    "fig, axes = plot_signal_and_spectrum(signal_step, sample_rate, \"Step function $t_0 = {:g}ns$\".format(t_0 * 10**9), ft_kwargs={'mask_bias':True});\n",
    "axes.flat[1].set_xlim(0, 20);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### FTs: delta function\n",
    "\n",
    "All frequencies are accounted for, but the phase alternates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA7sAAAEjCAYAAADzMBbxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOzdd5xcZdn/8c93d1NICOkJ6SEQSkLPCqEJSEeaNIOFIk0FFcGHBxWFh0d/oo+AoIBUQUQBkU5oooRQAuyGFgglhGwSEkjbFEL6Xr8/zlkclt3sJNnZKft9v17zysw597nPdTJn58w1dzmKCMzMzMzMzMxKSVm+AzAzMzMzMzNraU52zczMzMzMrOQ42TUzMzMzM7OS42TXzMzMzMzMSo6TXTMzMzMzMys5TnbNzMzMzMys5DjZNTMzMzMzs5LjZNfMzLImaZqk/Vuorq0kvSxpiaTvt0Sd67DvNyTt05r7NDMzs9blZNfMrA1Jk9VlaYK5UNJzkr4taZ2vBy2Q+J4PPBURXSLiqg2oZ60aizMiRkbEUzna3weSdsxF3Wn9Z0uqkrRC0i3NlO0g6SZJNel7/rKkQzLWD5U0VlKtpA8l/UFSRa5iz9hvD0n3Slqaxva1BuvzEpeZmZUWJ7tmZm3P4RHRBRgCXAr8N3BTHuIYAryRh/3mjKReQB9gcg53Mwv4BXBzFmUrgBnA3kBX4GfAXZKGpuuvAeYA/YAd03Lf3ZDgJF0s6eJmil0NrAT6Al8HrpU0MmN9i8dlZmZtj5NdM7M2KiIWRcQDwFeBkyRtCyCpv6R/SJor6f3GuhhLug0YDDwo6WNJ56fLL5D0XtqK+KakrzS2b0n/AvYF/pBuv2W6PCRtkVHuFkm/SJ9Pk/QjSa9JWiTpTkkdM8oOknRPGvf8tDWwqTg/be2VtI2kp9KW7jckHdEg1rXuN6PcFiSJZRkwP42hxVsjI+KeiLgPmJ9F2aURcXFETIuIuoh4CHgfGJUW2Qy4KyKWR8SHwKNAZtKJpNPT93KRpEck9dmQ+CV1Bo4BfhYRH0fEM8ADwDczijUZV9oT4WFJV0uaJ2mWpAMy6h8m6aF03SJJT2xIvGZmVryc7JqZtXER8SIwE9gr7c78IPAqMADYDzhH0kENtvkmMJ2klXjjiPhNuuo9YC+SVsT/Af4iqV8j+/wSMB44O93+nSzDPR44mCQZ2h44GUBSOfAQUAMMTWO/Yy1xkm7XLj3ex0laZL8H3C5pq2z22+CYpgA/Au5O99UzIlY3d0BpYrawicdD2fynZEtSX2BL/tOifiUwRlInSQOAQ0gSy/ryPwG+DRwB9AY+IGlV3hBbAmsavOev8tkke21xbQ/sRpIg9wGuI+mdUO/PwCMkrcZ9gYs3MF4zMytSTnbNzAySrrE9gC8AvSPikohYGRFTgRuAMdlUEhF/j4hZaSvincC7wC4tGOdVaf0LSJLU+rGxuwD9gf9KWzOXpy2GzRkNbAxcmh7vv0iS5hOy3G9DOwCvZC6Q9EtJ4yXdLalTww0i4rCI6NbE47AsjiEraWJ/O3BrRLyVLh5HkmQuJvnBowq4Ly3fB7gQOCEipkTESpLu7l/YwFA2BhY1WLYI6JLxusm4SJLdSyPisYioA95sUNfmQDlQnp4Hz25gvGZmVqSc7JqZGSQtoQtIxtH2z2xdBH5C0kLWLEknSnolY9ttgV4tGOeHGc8/IUmcAAYBNdm0pDbQH5iRJk31akj+P7LZb0M7krRSApB2Dd88IvYC/gl8ax3jaxFpi/1tJONkz85Y9hhwD9CZ5H3qDvw63Ww/oD3wYsb7+SifT1Q/0zoNXABcsJbW6Y+BTRos2wRYkmVc25H84FBvWz6b8H4dOBKYpWRyrh7N/PeYmVmJcrJrZtbGSfoCSXL3DMmY0/cbtC52iYhDG9k0GtQzhKQV+GygZ0R0AyYBWodwPgEyWz83zXK7GcDgJsbIRiPL6s0CBumzs1EPJumuu07SOrblsy27e5F0qSX9d89GtnskHU/c2OORhuXXIy6RtMj2BY6JiFXpqh4kPxL8ISJWRMR84E/AoRnr721wLnSNiH0a7iOzdZpk0rNL19I6/Q5QIWl4xrId+E/X6ibjkrQZyaRbb2dsuxMZ/+cR8a+I2A8YkdZ7crb/V2ZmVlqc7JqZtVGSNpF0GHAH8JeIeB14EVgs6b8lbSSpXNK2aULc0EfAsIzXnUkSy7lp/aeQJH/r4hXga+l+DyaZhTcbLwKzgUsldZbUUdIeTcSZ6QVgKXC+pHZK7r17OMn/ybraKH1kXlu785+W0EUkidxnRMQh6Rjfxh6HNCwvqSKdIKscKE+PdW0TYV0LbEMybnlZxn7nkUxW9Z20zm7ASfynZXoisK+kndP9biLpyDR5Xm8RsZSk1faS9L3ag6Ql9rYs4toeeL1BS/xO9TFLOlrS8DTGLiT//5/pVm5mZm2Hk10zs7bnQUlLSFpDfwpcDpwCEBFrSJK9HUkSjnnAjSQTTjX0K+DCtKvqjyLiTeAy4HmSBHM7YF3HS/4g3f9Cku6o9629eCIj7i1IJqSaSTLL9OfibLDdSpLJlw4hOdZrgBMzxrRmLU3i/gi8KWlmuriW//zfdSXpKr6hLgSWkXQX/kb6/ML6lWlL8U/S50OAM0nezw8zWoy/nhY/mmTirbnAFGA18MP0eJ4HLgH+Ieljkq7CB0fE2lrKs/Vdkh8G5gB/A74TEZm3oWoqru3JSF4l9SRp/Z+ULtqTZLzvEmAsSQvzv1ogXjMzK0JqmWuWmZmZNSRpO+DHEfE1SWcAHSLi9/mOy8zMrC1o8fv/mZmZWSIiXpdUI2k8SSvmifmOyczMrK1wy65ZEUm7Hp4UEQfmeD/7kIzhHJjL/ZiZmZmZ5YrH7JoVIEl7SnpO0iJJCyQ9K+kLEXF7rhNdMzOztqSpa24O9zdN0v65qt/M/sPdmM0KjKRNgIeA7wB3kdznci9gRT7jMjMzKzWFeM2VVLEe9ww3s0a4Zdes8GwJEBF/i4g1EbEsIh6PiNcknSzpmfqCkg6U9Hb6a/Q1ksZJOi1dd7KkZyT9VlKtpPclHZKx7SmSJktaImmqpDNb/1DNzMzyqrlr7rOSfp9eZ9+StF/9hpK6SrpJ0mxJH0j6haTyjPWnZ1xn35S0s6TbSO7l/WA6M/r5koZKCkmnSpoO/EvSPhmzutfX92mLsKSLJf1d0l/S+l+XtKWkH0uaI2mGJPcEszbPya5Z4XkHWCPpVkmHSOreWCFJvYC7gR8DPYG3gd0bFNs1Xd4L+A1wU8Y9MucAhwGbkNx25or6+2mamZm1Ec1dc3cFppJcRy8C7pFUf7/sW0lui7UFyf2eDwTqf3A+DriYZFK6TUhucTY/Ir5Jcnu0w9N7af8mY197k9wT+6AsYz+c5P7U3YGXgcdIvtsPILlt2HVZ1mNWspzsmhWYiFhMcq/IAG4A5kp6QFLfBkUPBd6IiHvS7k5XAR82KFMTETek9yC9FegH9E3383BEvBeJccDjJF23zMzM2oQsrrlzgN9FxKqIuJPkB+Qvp+sPAc6JiKURMQe4AhiTbnca8JuIeCm9zk6JiJpmwrk4rWtZluGPj4jH0u8Afwd6k9xbehVwBzBUUrcs6zIrSU52zQpQREyOiJPT2ZC3BfoDv2tQrD8wI2ObAGY2KPNhxvpP0qcbA6S/YE9IJ+NYSJI892rZIzEzMytszVxzP4jP3rqkJl0/BGgHzJa0ML2OXgf0ScsNAt5bx1BmNF/kMz7KeL4MmJf+uF3/GtJrvllb5WTXrMBFxFvALSQX4EyzgU9vDZR2T87qVkGSOgD/AH4L9I2IbsBYQGvd0MzMrIQ1cs0dkDH8B5LxtrNIEtMVQK+I6JY+NomIkWm5GcDmTe0mi+VLgU71L9KxwL3X5VjMzMmuWcGRtLWk8yQNTF8PAk4AJjQo+jCwnaSjJFUAZwGbZrmb9kAHYC6wOp24yhNZmJlZm5LFNbcP8H1J7dJxuNsAYyNiNsnwn8skbSKpTNLmkvZOt7sR+JGkUUpsIWlIuu4jYFgzob0DdJT0ZUntgAt
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "index_t0 = n_samples//8\n",
    "t_0 = t[index_t0]\n",
    "signal_step = np.zeros(n_samples)\n",
    "signal_step[index_t0] = 1\n",
    "\n",
    "fig, axes = plot_signal_and_spectrum(signal_step, sample_rate, \"Delta function $t_0 = {:g}ns$\".format(t_0 * 10**9), ft_kwargs={'mask_bias':True});\n",
    "axes.flat[1].set_xlim(0, 20);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Peculiarities $F^{-1}( F(A(t)) ) \\sim A(t)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# DFT\n",
    "\n",
    "A naive implementation of a DFT is simple to setup by implementing Eq. (1)\n",
    "\n",
    "$$\n",
    "\\hat{f}(\\omega) = F(t) = \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N-1} f(t_k)\\, \\exp{ \\left(-i \\omega t_k \\right)}\n",
    ".\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def naive_DFT( signal ):\n",
    "    \"\"\"\n",
    "    Calculate DFT in a naive way, returns a complex DFT.\n",
    "    \"\"\"\n",
    "\n",
    "    N = np.size(signal)\n",
    "\n",
    "    x_k = 2*np.pi/N*np.arange(N)\n",
    "    ft = np.zeros(N, dtype=np.complex128)\n",
    "    for j in range(N):\n",
    "        ft[j] = np.dot(signal, np.exp(- 1j * j * x_k ))\n",
    "\n",
    "    return ft"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### DFT: Cooley and Tukey (or Divide and Conquer)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A smarter DFT (Cooley and Tukey) can be implemented by splitting the sum in Eq. (1) into two parts: $k$ even and $k$ odd.\n",
    "\n",
    "With a uniform sample space $t$, the label $t_k$ can be rewritten to $t_k = k \\Delta t = $.\n",
    "\n",
    "This gives\n",
    "$$\n",
    "F(t) = \n",
    "    \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N/2-1} f(t_{2k})\\, \\exp{ \\left(-i \\omega \\, t_{2k} \\right)}\n",
    "    +\n",
    "    \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N/2-1} f(t_{2k+1})\\, \\exp{ \\left(-i \\omega \\, t_{2k+1} \\right)}\n",
    "    \\\\\n",
    "    =\n",
    "    \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N/2-1} f(t_{2k})\\, \\exp{ \\left(-i \\omega \\, 2k \\Delta t \\right)}\n",
    "    +\n",
    "    \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N/2-1} f(t_{2k+1})\\, \\exp{ \\left(-i \\omega \\, 2k \\Delta t \\right)}\\exp{ \\left(-i \\omega \\, \\Delta t \\right)}\n",
    "    \\\\\n",
    "    =\n",
    "    \\frac{1}{\\sqrt{2\\pi}} \\sum_{k=0}^{N/2-1} f(t_{2k})\\, \\exp{ \\left(-i \\omega \\, 2k \\Delta t \\right)}\n",
    "    +\n",
    "    \\frac{\\exp{ \\left(-i \\omega \\, \\Delta t \\right)} }{\\sqrt{2\\pi}} \\sum_{k=0}^{N/2-1} f(t_{2k+1})\\, \\exp{ \\left(-i \\omega \\, 2k \\Delta t \\right)}\n",
    "    \\\\\n",
    "    =\n",
    "    F(t_{even}) + F(t_{odd}) \\exp{ \\left(-i \\omega \\, \\Delta t \\right)}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cooley_tukey_DFT( signal ):\n",
    "    \"\"\"\n",
    "    Calculate a DFT according to the Cooley and Tukey algorithm (radix-2)\n",
    "    also known as Divide and Conquer.\n",
    "    \n",
    "    It divides the terms into 2 streams with odd and even index.\n",
    "    \"\"\"\n",
    "    radix = 2\n",
    "\n",
    "    N = np.size(signal)\n",
    "    \n",
    "    if N <= radix:\n",
    "        return naive_DFT(signal)\n",
    "    \n",
    "    if N % radix > 0:\n",
    "        raise ValueError(\"Provide a signal with  N%{radix} == 0 values.\".format(radix=radix))\n",
    "    \n",
    "    # divide and conquer\n",
    "    DFT_even = cooley_tukey_DFT(signal[::radix])\n",
    "    DFT_odd  = cooley_tukey_DFT(signal[1::radix])\n",
    "    \n",
    "    phase_terms = np.exp( -2j * np.pi * np.arange(0, N) / N )\n",
    "    \n",
    "    return np.tile(DFT_even, 2) + np.tile(DFT_odd,2) * phase_terms\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA74AAAEjCAYAAAAVGd21AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOy9ebxfVXU2/qyESSYJMzKjiBNOoOJb+zorOKC1g6BVsVpfbamt2vpqa9VS+zq0/rRWqwIiziiKCsggVUFAAyTMhASSQEgIkDmEJCS596zfH+ec79l77bX2OefekO+9uev5fCDfc/bea6+19j7n7nWePRAzw+FwOBwOh8PhcDgcju0V04atgMPhcDgcDofD4XA4HI8lPPB1OBwOh8PhcDgcDsd2DQ98HQ6Hw+FwOBwOh8OxXcMDX4fD4XA4HA6Hw+FwbNfwwNfhcDgcDofD4XA4HNs1PPB1OBwOh8PhcDgcDsd2DQ98HQ6Hw+FwOBwOh8OxXcMDX4fD4XA4HA6Hw+FwbNfwwNfhcDgcrSCiO4joJVtJ1r1E9IoxlDuGiG4ionVE9P6toctUAxHtR0RXENFqIvrGsPVxOBwOh2NbwQNfh8PhcAAAiOhFRPQ7IlpLRKuI6Foieh4AMPPTmfnKIav4YQBXMvMezPylIesy4dDxg8JHAdzNzDOY+V3jqGtvIvopEa0nokVE9Jax5m+T1SH9dCKaRUSbiOhckXYvEW0mon3F/ZuJiInoCJH3FSLfaUR0TUe3OBwOh2MCwwNfh8PhcICI9gRwMYD/ArA3gIMB/AuATcPUS+BwAHdoCUS0wzbWZbLiFQDO3wpyvgJgM4ADALwVwFeJ6OljzN8mqy19KYBPATjHqPseAKfWF0R0LIDHdbDR4XA4HNsRPPB1OBwOBwA8GQCY+QfMPMrMG5n5l8x8KxCzYdXvvyeiWyt2+IdEtEstiIieG0xJPr9K/5RWKRE9gYh+QkTLiegeawozEf0awEsBfJmIHiGiJ1d6/F8iuhXAeiLaoU0eET2HiG6sdPshEZ1X61YxgE8K8p4b6p2T3cEnhxLRBVXZlUT05er+PxDRT4SO/0VEX1R88BEiWlDpPoeI/ihI+w6AwwBcVPnnw6LsTkS0FsCxVZ7bND93ARHtBuCPAfwzMz/CzNcAuBDA2/rmb5PVpS5mvoCZfwZgpaHydwC8Pbh+B4Bvj8HuN1e+rf/bRERX9pXjcDgcjuHAA1+Hw+FwAMBdAEaJ6FtEdBIRzWjJ/2cATgRwJIBnAjgNKAMsAD8FcC5K5vgHAP5IE0BE0wBcBOAWlAzzywH8HRG9WuZl5pcBuBrA6cy8OzPfVSWdCuC1APYCUOTkVbr9DGUgtDdK5vOPW+zso6vlk+ko2fRFAI6oyp9XlfkugBOJaK8q7w4A3lzpKLEAwB8CeDxKNv67RHRQ5Z+3AbgPwOsr/3xO+G8zgBcCWFalHyvsu5iI1hj/XSz0eDKA0aANUPnFYnxz+dtk9a1Lw0wAexLRU6u2eDNKv/cCM/+w8t3uAJ4AYCHK/u1wOByOSQAPfB0Oh8MBZn4YwIsAMICzACwnoguJ6ACjyJeYeSkzr0IZED67un8CgB2q9C3MfAGA6w0ZzwOwHzOfwcybmXlhVfcpPVT/EjMvZuaNHeSdAGBHAF+sdPsxgBs61tNFV8snz0cZKP0DM69n5kcr5hLM/ACA3wL40yrviQBWMPNsqQAzn1/JL5j5hwDurmR3xbNRBo0JmPl1zLyX8d/rRPbdAawV99YC2MOoN5e/TVbfuizUrO8rAcwFcL+R72dh0A/gv2WG6iPI91GuN/96Tz0cDofDMST4miiHw+FwAACY+U40LOVTULJiX0SwPjLAg8HvDSgDO1T/3s/MHKQvNqo8HMATqgCjxnSUzG5XhLLb5Gm6LepYTxddLZ8cCmARM48Ysr8F4H0oA+k/h872gojeDuCDKFljoAwK99XyGjAD3554BMCe4t6eANaNIX+brL51WfgOyg8MRyI/zfmNzPw/9QURnQbg3SLPv6EMvH1ncYfD4ZhEcMbX4XA4HAmYeS7K6crP6Fn0AQAHExEF9w418i4GcI9gF/dg5tf0UbWHPE23w4LfGwDsGlwfuJV0XQzgMLI34PoZgGcS0TMAvA7A92QGIjocZWB8OoB9mHkvALcDCG1hWU7gWTACXyK6VKxfDf+7VGS/C8AORHS0kK1uPNaSv01W37pUMPMilJtcvQbABX3KhiCiU1B+CPoTZt4yVjkOh8Ph2PbwwNfhcDgcIKKnENGHiOiQ6vpQlAP8mT1F/R7AKIDTq82m3gB7Ou71AB6uNqh6HBFNJ6JnUHWE0hjQJu/3AEYAvL/S7U1Ct5sBvKUqdyKAF28lXa9HGXR/hoh2I6JdiOgP6kRmfhTAj1FOn72eme9TZOyGMrBdDgBE9E6kHyUeAnBURg8z8GXmk+r1q8p/J4m861EGj2dU9vwBgDfAYKpz+dtkdamrastdUDLw0yv/ah8Z3gXgZZXM3iCi56Dc9fyNzLx8LDIcDofDMTx44OtwOBwOoJw6+gIA1xHRepQB7+0APtRHSLWJ0ptQBhlrUE7dvRjKsUjMPArg9Sin4N4DYAWAs1Fu3tQbbfIC3U4DsBrlJkch+/e3Vfk1KI/N+dnW0DUo+ySUG1AtqeoO8S2UOy5bweMcAJ9HGbw/VOW9VmT7NICPVetT/z5MIKIDAcxAub51a+CvUB4JtAzlBk/vY+YBC1sxyP/YMX9WVof0jwHYCOAjKPvbxupeBGZewMyzxmxxGXDPAHBNhg13OBwOxwQFxUudHA6Hw+HYuiCi6wB8jZm/OWxdJIjoXABLmDkJlLaxHoehDEoPrDYaczgcDofDsRXhjK/DMUVARG8lol9ug3peQkRLHut6HBMXRPRiIjqwmoL6DpRH+1w2bL0mKqpdgj8I4DwPeh0Oh8PheGzgga/DsZ2BiF5ERL8jorVEtIqIriWi5zHz95j5VcPWzzElcAzKtaRrUU6V/pPq2B6HABHtBuBhlMfsfGLI6jgcjh6w/t4+hvXdS0SveKzkOxzbO/w4I4djOwIR7YlyPeX7APwIwE4A/hDK+kqH47ECM58J4Mxh69EFzHzakOtfj/JYIofDMYkwEf/eEtEOmWPTHI4pD2d8HY7tC08GAGb+ATOPMvNGZv4lM99KRKcR0TV1RiJ6FRHNq75U/zcRXUVE767STiOia4joP4hoNRHdQ0QnBWXfSUR3EtE6IlpIRP9n25vqcDgcDsfQ0Pb39loi+q/qb+xcInp5XZCIHk9E3yCiB4jofiL6FBFND9L/MvgbO4eInktE30F5/NpF1cZqHyaiI4iIiehdRHQfgF9ry41CppiIPklE5xPRdyv5txHRk4noo0S0jIgWE5HPDnNsl/DA1+HYvnAXgFEi+hYRnUREM7RMRLQvyuNTPgpgHwDzAPwvke0F1f19AXwOwDeIBuefLkN53uieAN4J4AtE9NytbYzD4XA4HBMUbX9vXwBgIcq/oZ8AcAER7V2lfQvl0WpPAvAcAK8CUH94/lMAnwTwdpR/Y08GsJKZ34ZyV/jXV8eMfS6o68UAngrg1R11fz3KHeRnALgJwOUoY4KDAZwB4Osd5Tgckwoe+Doc2xGqjXFehPK8z7MALCeiC4noAJH1NQDuYOYLqmlRXwLwoMiziJnPqo5i+RaAgwAcUNXzi+poEGbmqwD8EuUUL4fD4XA4tnt0+Hu7DMAXmXkLM/8Q5Yfk11bpJwH4O2Zez8zLAHwBwClVuXcD+Bwz31D9jZ3PzIta1PlkJWtjR/WvZubLq7//5wPYD8BnmHkLgPMAHEFEe3WU5XBMGnjg63BsZ2DmO5n5NGY+BMAzADwBwBdFticAWByUYZRni4Z4MEjfUP3cHQCqr9szq8081qAMpPfdupY4HA6HwzFx0fL39n6OzwxdVKUfDmBHAA9UZ26vQcmw7l/lOxTAgp6qLG7PEuGh4PdGACuqj9z1NeB7Dzi2Q3jg63Bsx2DmuQDORfkHOcQDAA6pL6opzIegA4hoZwA/AfAfAA5g5r0AXAK
      "text/plain": [
       "<Figure size 1152x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "f = 1e2 # Hz\n",
    "if f > sample_rate/2:\n",
    "    print(\"Sampling a frequency above the sample_rate/2 gives problems\")\n",
    "\n",
    "signal_sine = np.sin( 2*np.pi*f*t )\n",
    "\n",
    "fig, axes = plot_signal_and_spectrum(signal_sine, sample_rate, \"Single frequency at $f = {:g}$MHz\".format(f/10**6), ft_kwargs={'fft': cooley_tukey_DFT})\n",
    "axes.flat[1].axvline(sample_rate/2, color='r', label='Nyquist Frequency');"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}