{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Interpolation/Fractional Delay/Resampling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[back to overview page](index.ipynb)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "WARNING: this is work-in-progress!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Digital Audio Resampling Home Page: https://ccrma.stanford.edu/~jos/resample/\n",
    "\n",
    "Delay-Line and Signal Interpolation: https://ccrma.stanford.edu/~jos/pasp06/Delay_Line_Signal_Interpolation.html\n",
    "\n",
    "Time-Varying Delay Effects: https://ccrma.stanford.edu/~jos/pasp06/Time_Varying_Delay_Effects.html\n",
    "\n",
    "MUS420/EE367A Lecture 4A: Interpolated Delay Lines, Ideal Bandlimited Interpolation, and Fractional Delay Filter Design: https://ccrma.stanford.edu/~jos/Interpolation/Welcome.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "TODO: a lot of explanations\n",
    "\n",
    "For an introduction on NumPy and how to create simple signals, see [Creating Simple Audio Signals](http://nbviewer.ipython.org/urls/raw.github.com/mgeier/python-audio/master/simple-signals.ipynb)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First, let's set the sampling rate we'll be using:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fs = 44100  # Hz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A nice sine tone as test signal ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "duration = 2  # seconds\n",
    "frequency = 440  # Hz\n",
    "\n",
    "t = np.arange(int(duration * fs)) / fs\n",
    "sine = np.sin(2 * np.pi * frequency * t)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(sine[:100])\n",
    "plt.ylim(-1.1, 1.1);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO: fade in/out!\n",
    "# TODO: write WAV file"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "According to [JOS](https://ccrma.stanford.edu/~jos/pasp06/Linear_Interpolation.html):\n",
    "$\\hat{y}(n+\\eta) = (1-\\eta) \\cdot y(n) + \\eta \\cdot y(n+1)$\n",
    "\n",
    "That's a one-zero FIR filter.\n",
    "\n",
    "TODO: block diagram like in [JOS](https://ccrma.stanford.edu/~jos/pasp06/Fractional_Delay_Filtering_Linear.html)?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def dumb_linear_interpolation(x, fraction):\n",
    "    \"This is a very dumb implementation!\"\n",
    "    x = np.asarray(x)\n",
    "    assert x.ndim == 1, \"x must be one-dimensional!\"\n",
    "    \n",
    "    y_len = len(x) - 1\n",
    "    y = np.empty(y_len)\n",
    "    for i in range(y_len):\n",
    "        y[i] = (1 - fraction) * x[i] + fraction * x[i + 1]\n",
    "    return y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sine2 = dumb_linear_interpolation(sine, 0.5)\n",
    "\n",
    "len(sine2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note: length is one sample less"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(sine2[:100])\n",
    "plt.ylim(-1.1, 1.1);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks OK? Sounds OK?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO: same with noise, listen to both versions (WAV file) -> low-pass?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "OK, this is actually FIR filtering (TODO: explain), so we can do it with a convolution:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO: explain mode='valid'\n",
    "\n",
    "sine3 = np.convolve(sine, [0.5, 0.5], mode='valid')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "max(abs(sine2 - sine3))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When I tried it, there wasn't even a rounding error. YMMV."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check the frequency response of our interpolation filter to check if it's really a low-pass filter."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO: better function to display frequency response?\n",
    "# TODO: use scipy.signal.freqz()?\n",
    "# TODO: in dB, log frequency?\n",
    "# TODO: x-axis: frequencies (np.fftfreq()?)\n",
    "\n",
    "fftlength = 1000\n",
    "plt.plot(abs(np.fft.rfft([0.5, 0.5], fftlength)));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Yes, that looks like a low-pass filter!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(abs(np.fft.rfft(np.column_stack((\n",
    "  [1  , 0  ],\n",
    "  [0.9, 0.1],\n",
    "  [0.8, 0.2],\n",
    "  [0.7, 0.3],\n",
    "  [0.6, 0.4],\n",
    "  [0.5, 0.5])), fftlength, axis=0)))\n",
    "plt.ylim(0, 1.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "TODO: plot phase delay like in [JOS](https://ccrma.stanford.edu/~jos/pasp06/Fractional_Delay_Filtering_Linear.html)?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO: many more things!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "TODO: [One-Multiply Linear Interpolation](https://ccrma.stanford.edu/~jos/pasp06/One_Multiply_Linear_Interpolation.html)?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "TODO: [First-Order Allpass Interpolation](https://ccrma.stanford.edu/~jos/pasp06/First_Order_Allpass_Interpolation.html)?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<p xmlns:dct=\"http://purl.org/dc/terms/\">\n",
    "  <a rel=\"license\"\n",
    "     href=\"http://creativecommons.org/publicdomain/zero/1.0/\">\n",
    "    <img src=\"http://i.creativecommons.org/p/zero/1.0/88x31.png\" style=\"border-style: none;\" alt=\"CC0\" />\n",
    "  </a>\n",
    "  <br />\n",
    "  To the extent possible under law,\n",
    "  <span rel=\"dct:publisher\" resource=\"[_:publisher]\">the person who associated CC0</span>\n",
    "  with this work has waived all copyright and related or neighboring\n",
    "  rights to this work.\n",
    "</p>"
   ]
  }
 ],
 "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.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}