{
 "metadata": {
  "name": "",
  "signature": "sha256:258f63577a368e82fccca83ad8461dfd62cc4c6ae6053b633251891965b7c5f3"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%matplotlib inline\n",
      "import matplotlib.pyplot as plt\n",
      "import numpy as np\n",
      "import math"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 9
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np\n",
      "# Should be 4 Hz, 10 Hz, and 20 Hz\n",
      "f = lambda t: 2*np.cos(4*2*np.pi*t) + 5*np.sin(10*2*np.pi*t) + 3*np.cos(20*2*np.pi*t)\n",
      "n = 256\n",
      "r = np.arange(0.0,1.0,1.0/n)\n",
      "s = f(r)\n",
      "F = np.fft.rfft(s)\n",
      "m = n/2\n",
      "p = lambda z: (abs(np.real(z))/m,abs(np.imag(z))/m)\n",
      "t = p(F)\n",
      "tol = 1.0e-8\n",
      "for i in range(0,len(t[0])):\n",
      "    if(t[0][i] > tol):\n",
      "        print str(t[0][i]) + '*cos(' + str(i) + '*2*pi*t)'\n",
      "    if(t[1][i] > tol):\n",
      "        print str(t[1][i]) + '*sin(' + str(i) + '*2*pi*t)'\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "2.0*cos(4*2*pi*t)\n",
        "5.0*sin(10*2*pi*t)\n",
        "3.0*cos(20*2*pi*t)\n"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "p2 = lambda z: math.sqrt(np.real(z)**2 + np.imag(z)**2)/m\n",
      "t = [p2(i) for i in F]\n",
      "freq = np.fft.fftfreq(n)[0:n/2]\n",
      "print len(t[0:-1]), len(freq)\n",
      "freq2 = [i * n for i in freq] # Normalize frequencies by multiplying by the sample number\n",
      "plt.plot(freq2, t[0:-1], 'r')\n",
      "print freq2\n",
      "d = dict(zip(freq2, t[0:-1]))\n",
      "for i in range(25):\n",
      "    print \"Frequency {0} = {1}\".format(i, d[i])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "128 128\n",
        "[0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 17.0, 18.0, 19.0, 20.0, 21.0, 22.0, 23.0, 24.0, 25.0, 26.0, 27.0, 28.0, 29.0, 30.0, 31.0, 32.0, 33.0, 34.0, 35.0, 36.0, 37.0, 38.0, 39.0, 40.0, 41.0, 42.0, 43.0, 44.0, 45.0, 46.0, 47.0, 48.0, 49.0, 50.0, 51.0, 52.0, 53.0, 54.0, 55.0, 56.0, 57.0, 58.0, 59.0, 60.0, 61.0, 62.0, 63.0, 64.0, 65.0, 66.0, 67.0, 68.0, 69.0, 70.0, 71.0, 72.0, 73.0, 74.0, 75.0, 76.0, 77.0, 78.0, 79.0, 80.0, 81.0, 82.0, 83.0, 84.0, 85.0, 86.0, 87.0, 88.0, 89.0, 90.0, 91.0, 92.0, 93.0, 94.0, 95.0, 96.0, 97.0, 98.0, 99.0, 100.0, 101.0, 102.0, 103.0, 104.0, 105.0, 106.0, 107.0, 108.0, 109.0, 110.0, 111.0, 112.0, 113.0, 114.0, 115.0, 116.0, 117.0, 118.0, 119.0, 120.0, 121.0, 122.0, 123.0, 124.0, 125.0, 126.0, 127.0]"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n",
        "Frequency 0 = 2.74086309204e-16\n",
        "Frequency 1 = 3.37372809122e-16\n",
        "Frequency 2 = 4.60802170076e-16\n",
        "Frequency 3 = 4.39674274228e-16\n",
        "Frequency 4 = 2.0\n",
        "Frequency 5 = 6.16910896269e-16\n",
        "Frequency 6 = 2.27035497305e-15\n",
        "Frequency 7 = 1.15771511635e-15\n",
        "Frequency 8 = 9.37044033561e-16\n",
        "Frequency 9 = 2.10382661632e-15\n",
        "Frequency 10 = 5.0\n",
        "Frequency 11 = 1.53818926081e-15\n",
        "Frequency 12 = 1.15742001833e-15\n",
        "Frequency 13 = 1.25699850103e-16\n",
        "Frequency 14 = 8.55994122917e-16\n",
        "Frequency 15 = 1.83961677303e-16\n",
        "Frequency 16 = 6.76533916228e-16\n",
        "Frequency 17 = 1.36323700262e-15\n",
        "Frequency 18 = 1.59761217397e-15\n",
        "Frequency 19 = 2.19614846786e-15\n",
        "Frequency 20 = 3.0\n",
        "Frequency 21 = 2.21257098187e-15\n",
        "Frequency 22 = 1.1456704224e-15\n",
        "Frequency 23 = 1.08055363968e-15\n",
        "Frequency 24 = 1.29150607432e-15\n"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAW8AAAEACAYAAAB8nvebAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAEdhJREFUeJzt3W+spGV5x/Hf75yzy7JSQdKGancTNkYSNS1KCKVV2rG1\ndmtU+sKkELRWE0NiG7am0Rb6grNv2pCmVZIGNBUIoYW+ACUQTWRVxkBMKNJdurLAsqmwUGUxBRRr\n3D175uqLZ2Z3zpxn/s/zzH07308yOXNmZp+5DjnzOxfXfc8zjggBAPKyNO8CAADjI7wBIEOENwBk\niPAGgAwR3gCQIcIbADK0MsqDbD8r6SeS1iWtRcQlVRYFABhspPCWFJIaEfFylcUAAEYzztjElVUB\nABjLqOEdkr5h+7u2P1llQQCA4UYdm7wrIn5o+1ck7bP9VEQ8VGVhAID+RgrviPhh++uPbH9F0iWS\nHpIk25wcBQAmEBETj6OHjk1sb7f9S+3rr5P0PkkHewqo5nL0qOINb6ju+BG6/vrrKz1+1Rfqp/5F\nrT/n2iOm73lH6bzPk/QV253H/1tEPDD1M49ifb24AAA2GBreEfF9Se+ooZbNWq3iAgDYIO13WNbQ\neTcajUqPXzXqny/qn5+ca58FTzt7sR2zmN+Ueuop6cILpePHqzk+AMyJbUWVC5ZzxcwbAEqlHd7M\nvAGgVNrhvb4uRRQXAMAp6Yd391cAgCTCGwCylHZ4d+bdzL0BYIO0w5vOGwBKEd4AkCHCGwAylHZ4\nM/MGgFJphzedNwCUIrwBIEOENwBkKO3wZuYNAKXSDm86bwAoRXgDQIYIbwDIUNrhzcwbAEqlHd50\n3gBQivAGgAwR3gCQobTDm5k3AJRKO7zpvAGgFOENABkivAEgQ2mHNzNvACiVdnjTeQNAKcIbADJE\neANAhtIOb2beAFAq7fCm8waAUoQ3AGSI8AaADKUd3sy8AaDUSOFte9n2ftv3V13QBnTeAFBq1M57\nj6RDkqLCWjYjvAGg1NDwtr1D0vslfUmSK6+oG2MTACg1Suf9OUmfkVR/gtJ5A0CplUF32v6ApJci\nYr/tRr/Hra6unrreaDTUaPR96HgIbwC/IJrNpprN5syO54j+Y2zbfyfpo5JOStom6fWS7omIP+16\nTAw6xlQ+9Snp5pulL3xBuvrqap4DAObAtiJi4lH0wLFJRFwXETsjYpekKyR9qzu4K8fMGwBKjbvP\nu/7dJlu3MjYBgB4DZ97dIuLbkr5dYS2bra9LW7YQ3gDQI+13WNJ5A0CptMO71So6b2beALBB2uHN\n2AQASqUf3oxNAGATwhsAMpR2eDPzBoBSaYc3M28AKJV+eDM2AYBNCG8AyFDa4c3MGwBKpR3ezLwB\noFT64c3YBAA2IbwBIENphzczbwAolXZ4M/MGgFLphzdjEwDYhPAGgAylHd7MvAGgVNrhzcwbAEql\nH96MTQBgE8IbADKUdngz8waAUmmHNzNvACiVfngzNgGATQhvAMhQ2uHNzBsASqUd3sy8AaBU+uHN\n2AQANkk/vOm8AWCTtMObmTcAlEo7vBmbAEApwhsAMpR+eDPzBoBN0g5vZt4AUCrt8GZsAgClhoa3\n7W22H7F9wPYh239fR2GKKC6MTQBgk5VhD4iIn9t+T0T8zPaKpIdtvzsiHq60slZLWlqSlpcZmwBA\nj5HGJhHxs/bVrZKWJb1cWUUd6+tFeC8t0XkDQI+Rwtv2ku0Dko5JejAiDlVblorAXl4uLoQ3AGww\naufdioh3SNoh6XdsNyqtSqonvF97TfriF6s5NgBUaOjMu1tE/Nj2VyVdLKnZuX11dfXUYxqNhhqN\nxvSVtVqnw7uqmfeRI9KNN0pXX13N8QGgrdlsqtlszux4jojBD7B/WdLJiHjV9pmSvi5pb0R8s31/\nDDvGRF55Rdq1S3r8cemyy6SjR2f/HI8+Kl11lXT48OyPDQAD2FZEeNJ/P0rn/UZJt9teUjFmuaMT\n3JWqY2yytlZcACAzo2wVPCjpohpq2aiO8D55srgAQGbSfYdlHTNvOm8AmUo3vOvY5014A8hU2uHN\n2AQASi12eNN5A8hUuuFdx8ybzhtAptIN7zpn3lXsUweACqUd3nXMvCXOWgggO4sd3p15N3NvAJlJ\nN7w7M++lpdMfzDBrhDeATKUb3p2Zt11cqhhtdMYmLFoCyEza4b28XFyvanRC5w0gU4sd3nTeADKV\nbnh3Zt5SdXu96bwBZCrd8O7MvKXq9nrTeQPIVNrhzcwbAEoR3t1fASAT6YZ3q7VxbMJWQQA4Jd3w\npvMGgL4WO7zpvAFkarHDm84bQKbSDW9m3gDQV7rhTecNAH0R3t1fASATix3ejE0AZCrd8K5j5r22\nJm3ZQucNIDvphnddnfeZZ9J5A8jOYof32loR3nTeADKz2OFN5w0gU+mGd10zbzpvABlKN7wZmwBA\nX4sd3oxNAGRqscN7bU3avp3OG0B20g3vus5tQucNIEPphjczbwDoa7HDm84bQKbSDe9Wa2N4s1UQ\nAE4ZGt62d9p+0PYTtr9n+5o6CtP6+saZN2MTADhlZYTHrEn6dEQcsH2WpMds74uIJyutrM6xyWuv\nzf7YAFChoZ13RLwYEQfa138q6UlJb6q6MLYKAkB/Y828bZ8v6Z2SHqmimA3qmHmzYAkgUyOHd3tk\ncrekPe0OfHR33indcMN4lTHzBoC+Rpl5y/YWSfdI+teIuLf3/tXV1VPXG42GGo3GxgccPSo999x4\nlbFVEMAvkGazqWazObPjDQ1v25Z0i6RDEfH5ssd0h3ep48eLyziqDu/O8bZupfMGULnexnbv3r1T\nHW+Uscm7JH1E0nts729fdo/1LCdOFJdxVD3z7nwEGh+DBiBDQzvviHhY076ZZ9LOu8qZ98mT0spK\nEd6MTQBkpp53WKY4Nul03isrdN4AsrO44X3y5OmxCZ03gMykG95VnxJ2ba3ouum8AWQo3fCm8waA\nvhY3vDudN7tNAGSonvCeZKsgC5YA0Fe6nXfVM2+2CgLIWLrhTecNAH0tbnizYAkgY4sb3mwVBJCx\ndMO7jpk3nTeATNUb3hGj/xu2CgJAX/VtFYwYL4BZsASAvurrvLu/jqKOBUu2CgLIVH3hvbw8XnjX\ncW4TOm8AmaovvF//+vQ6bxYsAWSq+vButYrO9qyz0grv3q2C4yymAsCcVR/eJ04UnxO5bVta4d3p\nvJeWqhnLAECFqg/v48elM84oLqnNvFfanwLHdkEAmamv8966dbwzC9a1VVBi0RJAdtLtvOvaKiix\naAkgO4sb3nTeADKWbnjXdW4Tic4bQHbSDe+6tgpKdN4AsrO44U3nDSBjixvebBUEkLF0twrWdW4T\nibEJgOwsbufNVkEAGVvc8KbzBpCxdMO71doY3lVsFaTzBpCpdMN7fX3jzLvKzpsFSwCZSTu869oq\nuLJC5w0gK4sb3mwVBJCxesK7s1UwpbfHs2AJIGP17PPudN4pnRKWBUsAGVvssQmdN4BMDQ1v27fa\nPmb74ETPkGp403kDyNgonfdtknZP/AypnhKWrYIAMjY0vCPiIUmvTPwMKXfebBUEkKk0Z94Rm99h\nyVZBADhlZRYHWV1dPXW90Wio0WicvnOSrYKtlmQXF4kFSwDZazabajabMzvezMN7k0m2CnbPu6Xq\nPgaNBUsANeltbPfu3TvV8dIcm3TPuyU6bwDoMcpWwbskfUfSBbaft/3xsZ4h1fCm8waQsaFjk4i4\ncqpnSDW82SoIIGNpjk3qmnmzVRBAptI8MVVdnTdbBQFkKs3OmwVLABgozbMKsmAJAAPV13lv2VIE\n8Ciz6zpm3ixYAshYfeFtjz73rrvzZsESQGbqC29p9Lk3WwUBYKDFDW+2CgLIWLXhHXF6q6A0+tik\n6pl35w9B5znovAFkptrwXlsruuZOF51K5909MpHYKgggO9WGd2ebYMeo2wWrDu/uxUqJrYIAslNt\neHfPu6XpOu9WqxjDzEJv583YBEBm0gzv3pl350MZZhXevZ03C5YAMpNmePd23tJsRyd03gAyt5jh\n3b1NUKLzBpCd6sO7s01QGu8dlks9pS0tzbbz7l2wpPMGkJE0O+/uT47v6CxazgJbBQFkLo+tgtLs\nxyZsFQSQsTQ7bxYsAWCgdMO7ypk3WwUBZC7N8K575k3nDSAzaYY3WwUBYKB0twpWPfNmqyCAjKXb\neVe9z5uxCYCMpblVsOqZNwuWADKXbufNVkEA6Gsxw5vOG0Dm0gzv3lPCSrP9KLR+b4+f1SlnAaBi\naYZ33VsFl5Zm/zmZAFAhtgp2MPcGkJHF7Lx7xyYS4Q0gK+luFaxy5t27YCmxaAkgK3TeHXTeADKy\nmOFN5w0gc4sZ3nTeADI3NLxt77b9lO1nbP/1WEdPdZ9371ZBic4bQFYGhrftZUn/LGm3pLdJutL2\nW0c+egZbBZvNZnFbpp33qfozRf3zlXP9Odc+C8M670skHYmIZyNiTdK/S7p85KNnMDYhvOeL+ucr\n5/pzrn0WhoX3r0l6vuv7F9q3jSaXDyCWGJsAyMrKkPtHO9nHBz9Yfvvhw9K2bae/P/NM6dVX+z++\n45lnpA9/eONtKyvS6qp0880jlTTQE09Ie/ZsvO2MM6RrrpHOOWf649fp6aelxx6bdxWTo/75yrn+\nWdR+8cXS9dfPpp6aOQacjMn2pZJWI2J3+/trJbUi4oaux3A2JwCYQER40n87LLxXJD0t6fcl/UDS\nf0i6MiKenPQJAQDTGzg2iYiTtv9C0tclLUu6heAGgPkb2HkDANI01Tssp3oDzxzY3mn7QdtP2P6e\n7Wvat59re5/tw7YfsJ3sqqXtZdv7bd/f/j6n2s+xfbftJ20fsv2bmdV/bft356DtO22fkXL9tm+1\nfcz2wa7b+tbb/vmeab+m3zefqk/rU/8/tH9/Hrf9Zdtnd92XfP1d9/2V7Zbtc7tuG6v+icN76jfw\nzMeapE9HxNslXSrpz9s1/42kfRFxgaRvtr9P1R5Jh3R6J1BOtd8o6WsR8VZJvyHpKWVSv+3zJX1S\n0kUR8esqxohXKO36b1Px+uxWWq/tt0n6ExWv5d2SbrJd7ekzhiur/wFJb4+ICyUdlnStlFX9sr1T\n0h9Ieq7rtrHrn+aHm+4NPHMQES9GxIH29Z9KelLFvvUPSbq9/bDbJf3xfCoczPYOSe+X9CVJnVXq\nXGo/W9JlEXGrVKynRMSPlUn9kn6i4o//9vZC/nYVi/jJ1h8RD0l6pefmfvVeLumuiFiLiGclHVHx\nGp+bsvojYl9EdM6T8YikHe3rWdTf9k+SPttz29j1TxPe072BZ87andQ7VfwCnBcRx9p3HZN03pzK\nGuZzkj4jqfskL7nUvkvSj2zfZvs/bf+L7dcpk/oj4mVJ/yjpqIrQfjUi9imT+rv0q/dNKl7DHTm8\nnj8h6Wvt61nUb/tySS9ExH/13DV2/dOEd7YrnbbPknSPpD0R8Vr3fVGs4Cb3s9n+gKSXImK/Tnfd\nG6Rae9uKpIsk3RQRF0n6P/WMGFKu3/abJf2lpPNVvNDOsv2R7sekXH+ZEepN9mex/beSTkTEnQMe\nllT9trdLuk5S97uCBu3zHlj/NOH9P5J2dn2/Uxv/ciTJ9hYVwX1HRNzbvvmY7V9t3/9GSS/Nq74B\nflvSh2x/X9Jdkn7P9h3Ko3ap+N14ISIebX9/t4owfzGT+i+W9J2I+N+IOCnpy5J+S/nU39Hv96X3\n9byjfVtybP+ZivHhVV0351D/m1X88X+8/TreIekx2+dpgvqnCe/vSnqL7fNtb1UxbL9viuNVzrYl\n3SLpUER8vuuu+yR9rH39Y5Lu7f238xYR10XEzojYpWKh7FsR8VFlULtUrDdIet72Be2b3ivpCUn3\nK4P6VSyuXmr7zPbv0XtVLBznUn9Hv9+X+yRdYXur7V2S3qLiTXlJsb1bxejw8oj4edddydcfEQcj\n4ryI2NV+Hb+gYgH8mCapPyImvkj6IxXvwDwi6dppjlXHRdK7VcyLD0ja377slnSupG+oWL1+QNI5\n8651yM/xu5Lua1/PpnZJF0p6VNLjKjrXszOr/7Mq/uAcVLHYtyXl+lX8H9oPJJ1QsT718UH1qvhf\n+iMq/lD9YYL1f0LSMyp2aXRevzdlUP/xzn//nvv/W9K5k9bPm3QAIEPz3gcJAJgA4Q0AGSK8ASBD\nhDcAZIjwBoAMEd4AkCHCGwAyRHgDQIb+H1TLpGi9CiuFAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x11229e7d0>"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 11
    }
   ],
   "metadata": {}
  }
 ]
}