Farrow resampler.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Farrow resampler\n",
    "\n",
    "Based on:\n",
    "https://www.dsprelated.com/showarticle/149.php"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pylab as plt\n",
    "from scipy.signal import firwin\n",
    "from utils import *\n",
    "%matplotlib inline\n",
    "np.set_printoptions(precision=4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# test signal\n",
    "FREQ = 5\n",
    "N1 = 32\n",
    "N2 = 50\n",
    "x = np.arange(0, N1)\n",
    "y = np.sin(x*2*np.pi/(N1/FREQ))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "cMatrix = [\n",
    "-2.6766241133e-19, -7.8326517820e-05,  2.1673359859e-04, -1.1926223897e-18, -1.1626782202e-03,  3.4680189578e-03, -5.9158573370e-03,  5.8980610040e-03, -1.3010426070e-17, -1.3683354576e-02,  3.2542188101e-02, -4.7494867528e-02,  4.3379371285e-02, -4.1633363423e-17, -1.2008255824e-01,  6.0291326947e-01,  6.0291326947e-01, -1.2008255824e-01,  4.1633363423e-17,  4.3379371285e-02, -4.7494867528e-02,  3.2542188101e-02, -1.3683354576e-02, -2.6020852140e-18,  5.8980610040e-03, -5.9158573370e-03,  3.4680189578e-03, -1.1626782202e-03, -6.5052130349e-19,  2.1673359859e-04, -7.8326517820e-05, -2.7094466400e-19, \n",
    " 1.8250485395e-05, -3.0347364849e-04,  1.1023676037e-03, -2.3800654079e-03,  2.7421886180e-03,  6.7217025074e-04, -1.0687147953e-02,  2.6797569663e-02, -4.1593149400e-02,  3.9952935089e-02, -3.8338870620e-03, -7.7139156322e-02,  1.9218414293e-01, -2.9351804709e-01,  2.9105617226e-01,  7.8087126928e-01, -7.4306797975e-01, -3.5400012711e-01,  3.1972958168e-01, -1.8586046319e-01,  6.1597791470e-02,  1.7445879628e-02, -4.7695209593e-02,  4.3871248857e-02, -2.5864558193e-02,  8.7803206192e-03,  8.8423322827e-04, -3.5756073213e-03,  2.6663013673e-03, -1.1398700262e-03,  2.8506109980e-04,  1.2579585994e-06, \n",
    " 1.3580560084e-05,  3.5225747764e-04, -1.8863366779e-04, -2.3035136387e-03,  9.2919510111e-03, -2.0715650279e-02,  3.0951486679e-02, -2.7418054194e-02, -5.5590073439e-03,  7.6216672687e-02, -1.7084561350e-01,  2.4104051953e-01, -1.9477667773e-01, -1.2721196819e-01,  8.4380759156e-01, -7.8086984832e-01, -3.9179508427e-01,  7.7472965930e-01, -3.1713675930e-01,  2.1958985496e-02,  1.1308170522e-01, -1.3361542955e-01,  9.4847304983e-02, -4.3250724792e-02,  5.6002045697e-03,  1.1263158774e-02, -1.2918367164e-02,  8.2591818180e-03, -3.5800333740e-03,  1.0910858633e-03, -3.1689216527e-04, -3.0880465663e-06, \n",
    "-1.1015756330e-04,  2.4627628726e-04, -1.1304675345e-03,  3.5209008264e-03, -7.4034424511e-03,  1.0659603734e-02, -8.4504203851e-03, -5.2775764728e-03,  3.3468802168e-02, -6.9944065100e-02,  9.4642444937e-02, -7.3027124397e-02, -4.0786836481e-02,  3.0064745705e-01, -4.1186793611e-01, -1.4209542929e-06,  4.1186723632e-01, -3.0064697395e-01,  4.0786548904e-02,  7.3027238887e-02, -9.4642441064e-02,  6.9944007243e-02, -3.3468740814e-02,  5.2775369388e-03,  8.4504352827e-03, -1.0659603098e-02,  7.4034367573e-03, -3.5208962765e-03,  1.1304656052e-03, -2.4627595357e-04,  1.1015758329e-04,  1.8300879670e-06,  \n",
    "]\n",
    "NTAPS = 32\n",
    "ORDER = 3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def c(ixTap, ixPower):\n",
    "    return cMatrix[ixPower * (NTAPS) + ixTap]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "audio_delay_line = np.zeros(NTAPS)\n",
    "delay_line_sample = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def get_imp(tap, frac):\n",
    "    horner_sum = c(tap, ORDER)\n",
    "    for i in range(ORDER-1, -1, -1,):\n",
    "        horner_sum *= frac\n",
    "        horner_sum += c(tap, i)\n",
    "    return horner_sum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def get_cf(frac):\n",
    "    return [get_imp(i, frac) for i in range(NTAPS)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def process_sample(frac):\n",
    "    fir_c = get_cf(frac)\n",
    "    v = np.sum(audio_delay_line*fir_c)\n",
    "    return v"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x_out = np.arange(0, N1, N1/N2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "y_out = np.array([])\n",
    "for i in range(N2):\n",
    "    while(x[delay_line_sample] < np.floor(x_out[i])):\n",
    "        audio_delay_line = np.roll(audio_delay_line, 1)\n",
    "        audio_delay_line[0] = y[int(x_out[i])]\n",
    "        delay_line_sample += 1\n",
    "    v = process_sample(x_out[i] % 1)\n",
    "    y_out = np.append(y_out, v)\n",
    "# delay compensation\n",
    "x_out = x_out - (NTAPS-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x_ideal = np.arange(0, N1, 1/1000)\n",
    "y_ideal = np.sin(x_ideal*2*np.pi/(N1/FREQ))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "plt.plot(x_ideal, y_ideal, label=\"ideal\")\n",
    "plt.plot(x, y, 'o', label=\"sampled\")\n",
    "plt.plot(x_out, y_out, '.', label=\"Farrow\")\n",
    "plt.grid()\n",
    "plt.xlim(0, 15)\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  }
 ],
 "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.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}