agoose77/Question_2_fit.ipynb

## Question_2_fit.ipynb
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from math import pi\n",
    "from numpy import sqrt, exp, array\n",
    "from units import m\n",
    "from matplotlib import pyplot as plt\n",
    "from scipy.optimize import curve_fit\n",
    "%matplotlib inline  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Define unit of cm\n",
    "cm = 1e-2 * m"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Load data\n",
    "data = [(5012, 1), (3981, 2), (2512, 3), (1413, 4), (525, 5)]\n",
    "counts, distances = zip(*data)\n",
    "\n",
    "counts = array(counts)\n",
    "distances = array(distances)  # Convert to meters\n",
    "\n",
    "t = 20 * 60 * 60"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Define count function\n",
    "def count(x, N, D):\n",
    "    return N / sqrt(4* pi * D * t) * exp(-x**2 / (4 * D * t))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\python34\\lib\\site-packages\\ipykernel\\__main__.py:3: RuntimeWarning: invalid value encountered in sqrt\n",
      "  app.launch_new_instance()\n"
     ]
    }
   ],
   "source": [
    "(N, D), pcov = curve_fit(c, distances, counts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAEACAYAAABLfPrqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmcjeX/x/HXZ8a+RtaQZG1sUdGGibIlfMuIfiRLylKo\nlG3sLdqVqK8U+pKMCiFRTPVVSdmNZSqEmJB9Z67fH+eeOt/CzGHMOTPzfj4e8+iey3Wf+3Ou0bxd\n13Wfc8w5h4iISCDCgl2AiIikPwoPEREJmMJDREQCpvAQEZGAKTxERCRgCg8REQlYisLDzPKbWYyZ\nrTezdWZW28wKmNkCM9toZp+ZWX6//v3NLN7r39CvvaaZrTazTWb26qV4QiIicumldOYxGpjnnLsG\nqA5sAPoBnzvnKgKLgP4AZhYBtAauAZoAY83MvMcZB3R2zlUAKphZo1R7JiIikmaSDQ8zywfUcc69\nC+CcO+2cOwC0ACZ53SYBLb3j5sA0r98WIB6oZWbFgLzOuWVev8l+54iISDqSkplHGWCPmb1rZsvN\n7N9mlgso6pxLAHDO7QKKeP1LANv8zt/htZUAtvu1b/faREQknUlJeGQBagJvOOdqAkfwLVn9/X1N\n9D4nIiKZRJYU9NkObHPO/eB9/yG+8Egws6LOuQRvSep37893AKX8zi/ptZ2r/R/MTEEkInIBnHOW\nfK+Ll+zMw1ua2mZmFbymBsA6YDbwgNfWAZjlHc8G2phZNjMrA5QDvveWtg6YWS1vA/1+v3POdt2Q\n/xoyZEjQa8gINapO1RnqX+mlzrSUkpkHwKPAFDPLCvwCdATCgelm1gnYiu8OK5xzcWY2HYgDTgHd\n3V/PqgcwEciB7+6t+an1REREJO2kKDycc6uAG87yR7efo/+zwLNnaf8RqBpIgSIiEnr0CvOLEBkZ\nGewSkpUeagTVmdpUZ+pKL3WmJUvrdbKUMDMXinWJiIQyM8OFyoa5iIjI3yk8REQkYAoPEREJmMJD\nREQCpvAQEZGAKTxERCRgCg8REQmYwkNERAKm8BARkYApPEREJGAKDxERCZjCQ0REAqbwEBGRgCk8\nREQkYAoPEREJmMJDREQCpvAQEZGAKTxERCRgCg8REQmYwkNERAKm8BARkYApPEREJGAKDxERCZjC\nQ0REAqbwEBGRgCk8REQkYCkKDzPbYmarzGyFmX3vtRUwswVmttHMPjOz/H79+5tZvJmtN7OGfu01\nzWy1mW0ys1dT/+mIiEhaSOnMIxGIdM7VcM7V8tr6AZ875yoCi4D+AGYWAbQGrgGaAGPNzLxzxgGd\nnXMVgApm1iiVnoeIiKShlIaHnaVvC2CSdzwJaOkdNwemOedOO+e2APFALTMrBuR1zi3z+k32O+cf\nPo3/FOdcCssTEZG0lNLwcMBCM1tmZl28tqLOuQQA59wuoIjXXgLY5nfuDq+tBLDdr32713ZWjy94\nnDrv1iF2S2wKSxQRkbSSJYX9bnHO7TSzwsACM9uIL1D8peo04Ylcawir/j5dZnehTIEyjLxtJLVL\n1k7NS4iIyAVKUXg453Z6/91tZjOBWkCCmRV1ziV4S1K/e913AKX8Ti/ptZ2r/axGjhjBnj1Qq3Zb\nKjU/RlRMFNcWu5bhtw3n2mLXpvgJiohkVLGxscTGxgbl2pbcvoKZ5QLCnHOHzSw3sAAYBjQA/nDO\njTKzp4ACzrl+3ob5FKA2vmWphUB555wzs++AR4FlwFzgNefc/LNc0znnWLECXngBPvsMOnQ+zmUN\n/s24tc9St3RdhkUOo1KhSqk3EiIi6ZyZ4Zyz5HtevJTseRQF/mtmK4DvgE+ccwuAUcAd3hJWA+A5\nAOdcHDAdiAPmAd3dXwnVA5gAbALizxYc/mrUgKlT4ccfIfFkDl5t+yj11/5EcatBnXfr8MDMB/hl\n3y+BP2sREbkoyc48giFp5vF3f/wB48bB669DtVoHuOLuV5jz+xiiIqIYVHcQJfKdc/9dRCTDS8uZ\nR7oKjyTHj8PkyfDSS5Dz8j2Uvu8F/nv0bTpU70C/W/tRJHeRc54rIpJRKTySCY8kiYkwe7ZvX2T7\ngZ2Uaf8sa5jCw9c/zBM3P0GBnAXSoFoRkdCg8EhhePhbssQXIl+v2Urp9iP5NedMet/Ui161e5E3\ne95LVKmISOhQeFxAeCTZuNG3nDVtYTzF2gzlj8s+p3+dJ+l+Q3dyZs2ZypWKiIQOhcdFhEeSXbt8\nG+tvxKwlT7MhnCzyHUMbDKRLzS5kC8+WSpWKiIQOhUcqhEeSQ4dgwgQYNflHTtw6iGzFN/B0o8F0\nuLY9WcJS+gJ7EZHQp/BIxfBIcuoUxMTA4An/ZWelQeS7Yicv3DmM+6q3Jsz0sSYikv4pPC5BeCRx\nDhYudPQf/wVrigzk8qLHeKnZCNrWaM5f7xwvIpL+KDwuYXj4W77c0efNOSzJNojCBbPzYtOR3Ff7\nDoWIiKRLCo80Co8km7ck0nPsDOafGEzhnEUZ1XgkHSLrpNn1RURSg8IjjcMjye69p+k2bgoz9w2l\nQGJFht82gofvugFNREQkPVB4BCk8khw8cpLu49/hg50jyXPwBgbeNJze91Uli27OEpEQpvAIcngk\nOXLiGI9OfpP3No8i2476PHbdUPp1qUCuXMGuTETknxQeIRIeSQ6fPMzjMa8xccMrhG1qTtdKgxnU\nszSFCwe7MhGRvyg8Qiw8kuw/vp8Bn7zEO2vHwpq23HvFAKJ7X0G5csGuTERE4RGy4ZFk95HdDF44\nikkr34UVnbg951MMeqwQtWoFuzIRycwUHiEeHkl2HNzB0EVPM3X1B2Rd0YPKBx+jf5/LaNoUwvSi\ndRFJYwqPdBIeSTbv28yw2BF8tO4T8qx+jPwbH+HJ3nm47z7Inj3Y1YlIZqHwSGfhkWTjno0MiR3C\ngk2xFNnUj4OLHqZ3zxw89BDkzx/s6kQko1N4pNPwSLI6YTXRi6NZunU5pbcOIj6mI506ZKN3byhZ\nMtjViUhGlZbhoZX5S6Ba0WrMajOL2e0+5LKbPyJPv0qsDZ9M1epn6NAB1qwJdoUiIhdHM4808NXW\nrxi4aCAJh/Zw/aHhLBpzDzVrhNG3L0RGorc/EZFUoWWrDBYeAM45Fvy8gEGLB3Hq9GnqnB7Bgjfu\nJG8eo29fuOce9PYnInJRFB4ZMDySOOeYtXEW0YujyZM1D01zjmTBmw3YsQMeeww6doTcuYNdpYik\nRwqPDBweSc4knmH6uukMiR1CyXwlaVNkJPPH38x//wvdukHPnujtT0QkINowzwTCw8JpW7UtcT3i\naFetHU9vbMuJVnfy1uzl7NoFFSv6QuSnn4JdqYjIPyk8gixLWBY61ejEpp6baFquKT2/uYu9DVrx\n8ZJ1FCoEN90ErVrB0qXBrlRE5C8KjxCRPUt2etTqQfwj8dxY8kZaz63P5hrt+Hz5T9StC/feC/Xq\nwZw5kJgY7GpFJLNTeISYXFlz8cTNTxD/SDwVL69Ig/dvZG2Zriz68VcefhgGD4YqVeCdd+DEiWBX\nKyKZVYrDw8zCzGy5mc32vi9gZgvMbKOZfWZm+f369jezeDNbb2YN/dprmtlqM9tkZq+m7lPJWPJl\nz0d0vWg2PbKJQrkKccOEGnx3WS/mfrmL11+H6dPh6qth1CjYvz/Y1YpIZhPIzKMXEOf3fT/gc+dc\nRWAR0B/AzCKA1sA1QBNgrNmfL4MbB3R2zlUAKphZo4usP8MrmLMgzzR4hrjucYRZGFXGVWah68eU\nj/Yybx6sXQtly8ITT8D27cGuVkQyixSFh5mVBJoCb/s1twAmeceTgJbecXNgmnPutHNuCxAP1DKz\nYkBe59wyr99kv3MkGUXzFOWVxq+w8qGV7D++n4pjKjJz3zDeePsgK1b49kGqVYP779fbn4jIpZfS\nmccrQF/A/8UXRZ1zCQDOuV1AEa+9BLDNr98Or60E4P9v4+1emwSgVP5SvNnsTZZ2WcrP+36m3Gvl\nmLbteUY8d4Sff4aICGjUCJo0gUWLIIO/XEZEgiTZN8QwszuBBOfcSjOLPE/XVP01NXTo0D+PIyMj\niYw836Uzn7IFyzL5X5OJ2x3H0NihlH+9PP1v7U+fJ7rSp092/vMf6NED8uWDN9+EGjWCXbGIpLbY\n2FhiY2ODcu1kX2FuZs8A7YDTQE4gL/AxcD0Q6ZxL8JakFjvnrjGzfoBzzo3yzp8PDAG2JvXx2tsA\n9Zxz3c5yzQz/CvPUtmLnCqIXR7M6YTWD6w2mQ/UOhFtW3nsP+vaFBx+E6GjIkSPYlYrIpRJSrzB3\nzg1wzl3pnLsaaAMscs61Bz4BHvC6dQBmecezgTZmls3MygDlgO+9pa0DZlbL20C/3+8cuUg1itdg\nzn1z+KDVB0xbO42IsRG8v3YK7dqfYfVq2LDBN/v49ttgVyoiGUFA721lZvWAx51zzc2sIDAdKIVv\nVtHaObff69cf6AycAno55xZ47dcBE4EcwDznXK9zXEczj4u0ePNiBi4ayKGThxjTZAz1rqrHjBnw\n6KPQujU8/bTegFEko9EbIyo8UoVzjo/Wf0Tvz3pTv0x9XrjjBcKPF6F3b1iyBMaPhwYNgl2liKSW\nkFq2kvTLzLgn4h7iusdROFdhqoytwozNbzFpciJjxvje/r1rVzhwINiVikh6o/DIBPJmz8uLDV/k\n8/s/Z/Lqydw04SaK11jB2rUQFuZ7u5M5c4JdpYikJ1q2ymQSXSITV06k/xf9aVO5DcNvG87yb/Pz\n4INQuzaMHg2FCgW7ShG5EFq2kksmzMLoVKMTcd3jOHrqKBFjI0goPI1VqxzFikHVqr73zVJ2i8j5\naOaRyX2z7Ru6ze1GkdxFeKPpG+z7qQKdOkGFCjB2LBQvHuwKRSSlNPOQNHNzqZv5seuPNC3XlJsn\n3Mzco4NZsvQYVapA9erw7ruahYjIP2nmIX/afnA7fT7rw/KdyxnTZAxXHG1Cp06+PZC33oKrrgp2\nhSJyPpp5SFCUzFeSmKgY3mj6Bj0/7cmIja34cOF2brsNrr8exozRpxiKiI/CQ/6hcbnGrO22lsqF\nK3P929eSte5LLP7qFO+/7/so3E2bgl2hiASblq3kvOL3xtNjXg92Hd7FmCbjWDn7FoYPhyefhMce\ngyzJvi+ziKQVvT2JwiOkOOeYvm46jy14jCblmtCt/HM89Ugh9u/3fZZ6tWrBrlBEQHseEmLMjHur\n3Mv6HuvJky0PTedWps3zE3jo4URuvx2GDIGTJ4NdpYikJc08JGArdq6g29xuhIeFM/T6cYyJrsbP\nP/tmIbVqBbs6kcxLy1YKj5CX6BJ5e/nbDFo0iHbV2lN171D6P5aX9u1h2DDIlSvYFYpkPlq2kpAX\nZmF0va4ra7uvZd/xPxicEMEzM2ewfYejenX48stgVygil5JmHpIqvt76Nd3mdqNU/lK0zDaGEX3K\n0rw5PPec73PUReTS08xD0p06peuw4qEV1L+qPgM31+b+CcM5dvIEVavC/PnBrk5EUptmHpLqfj3w\nK73m92Ld7+voUnws4/reTt268MorULBgsKsTybi0Ya7wyBDmbJrDI58+wnVFa5Pv25eZH3MFr78O\n99wT7MpEMiYtW0mG0KxCM9Z1X0fFwmWZXaIaUS+NZsCg00RFQUJCsKsTkYuh8JBLKlfWXDzd4Gm+\n7vg1a07OIlevWuQqv5Rq1eC99/R27yLplZatJM0455i6Zip9F/blpoLN2fDGM5QuUpC33oJSpYJd\nnUj6p2UryZDMjP+r9n/E9YijeNEs7G0TQc4bJ1GjpuOtt/R27yLpiWYeEjQ//PYDD895GE7l4uTH\n47j8TGXGj4dy5YJdmUj6pJmHZArXX3E9S7sspVOtNuxsFIk1fIratx7hpZfgzJlgVyci56OZh4SE\nXYd30XdhXxb99BUFvh9N7m0teGeCUblysCsTST/0Og+FR6a1ePNius/tTtbD5dj+79fp0/Eq+vWD\nrFmDXZlI6NOylWRat5W5jVXdVtHmlpug6/X8Z+sz1LzhJD/+GOzKRMRfsuFhZtnNbKmZrTCzNWY2\nxGsvYGYLzGyjmX1mZvn9zulvZvFmtt7MGvq11zSz1Wa2ycxevTRPSdK7bOHZGFBnAD8+tIwK9b9l\nb+vq3N5lMf37w/Hjwa5ORCCFy1Zmlss5d9TMwoElwKPAPcBe59zzZvYUUMA518/MIoApwA1ASeBz\noLxzzpnZUqCnc26Zmc0DRjvnPjvL9bRsJYDvtSGzNs6i59xehG+vS9bFLzLpjaLcckuwKxMJPSG3\nbOWcO+odZgeyAA5oAUzy2icBLb3j5sA059xp59wWIB6oZWbFgLzOuWVev8l+54iclZnRslJLNj4S\nR5umV7C7VRWaDhlLz0fPcPhwsKsTybxSFB5mFmZmK4BdwEIvAIo65xIAnHO7gCJe9xLANr/Td3ht\nJYDtfu3bvTaRZOXOlptRd4xiSdfFVL53GtNy30iFyB9YuDDYlYlkTllS0sk5lwjUMLN8wMdmVhnf\n7ON/uqVmYUOHDv3zODIyksjIyNR8eEmnqhSpwpIuXzJ51WT65G3G3W+3osUHIxnz4mVcdlmwqxNJ\nW7GxscTGxgbl2gHfqmtm0cBRoAsQ6ZxL8JakFjvnrjGzfoBzzo3y+s8HhgBbk/p47W2Aes65bme5\nhvY8JFl/HPuDvvMHMG3FbLJ/+QLv9LmPli3TZLlXJCSF1J6HmRVKupPKzHICdwDrgdnAA163DsAs\n73g20MbMsplZGaAc8L23tHXAzGqZmQH3+50jErCCOQsy4V9v8kXXj7i8+YvcN/927uywgd27g12Z\nSMaX7MzDzKri2xAP874+cM49bWYFgelAKXyzitbOuf3eOf2BzsApoJdzboHXfh0wEcgBzHPO9TrH\nNTXzkICcTjzNq0vGMviLEYQt78qYewfS4b5cmCYikonoFeYKD7lAvx36jQ7vP8aX8d9TY+frfDTq\nTkrotgzJJBQeCg+5SHM3LOT+aT04/Etlht88micfulKzEMnwQmrPQyQ9urPSHewYtJrOzWowcFtN\nKnR8gU0/nwp2WSIZhmYekuFt3P0zzd/syc+7t9Gj9Dhe6VOHMP2zSTIgLVspPCSVOed4/YsP6ftF\nH/LvbcDM7i9w87WFg12WSKrSspVIKjMzHr29Fb8PjqNK2cu59f3K/GvkvzlxUp99K3IhNPOQTOnT\n5atpO+VhTp12vHP3OO6td22wSxK5aFq2UnhIGjiTmEiXN95l0vYBXJ+9LXOfGE7hfPmCXZbIBdOy\nlUgaCA8L491HOrP64XUk7D/EFU9HMPLjD9A/XESSp5mHCOAcjJi0hBHLu1EiX3FmdxtDtRLlg12W\nSEC0bKXwkCDZmXCKZiNeZ2XeZ2hXoSdvte9Hjiw5gl2WSIooPBQeEmQTP9pO99m9yVpqJe9GvcHd\n1RoFuySRZCk8FB4SAg4cgDaDPmVhtp7cWLomHzzwKiXy6Y2yJHRpw1wkBOTPD5++3oRPGq9lXew1\nlH2xOoMXjOLoqaPJnyySwWnmIZICR45Aj8GbmJowkOzlltCvzlM8Xu8h7YdISNGylcJDQtTmzfDE\niyv55NAQspdZzqC6A+gT2Zls4dmCXZqIwkPhIaFu2zZ4/OVlfLx/MDlLrWdwvWh6Rd5P1vCswS5N\nMjGFh8JD0omdO6HPK0v48I/B5Cq+laG3DeHRyPsIDwsPdmmSCSk8FB6SzuzeDb1eXUzMnmhyF97L\n8NuG0bN+K8JM96RI2lF4KDwkndq71/Hoawv4YHc0eS87wYj6w+jRoAWmjzGUNKDwUHhIOrd/v+OR\n1+cwLWEwefOEM7LBcLrd3kQhIpeUwkPhIRnEwUOJ9BzzMVN3DSFf9rw83WAEDzdsoBCRS0LhofCQ\nDObQ4TP0GDedqb8NJX9YMUY2GEG3pnWDXZZkMAoPhYdkUIePnqbHm1OY8tsw8p8u6wuRu24MdlmS\nQSg8FB6SwR05dopu/57I1B0jyH+sKiMaDKdbi+vQapZcDIWHwkMyiSPHT9Dt7beZuu0Z8h2qxcj6\nw+l2T1WFiFwQhYfCQzKZw8eP0e2dcUzb9jx590YyvP5QureuRJheJiIBUHgoPCSTOnj8MN0njuGD\nX18mz67GDG8wmB5tyylEJEUUHgoPyeQOHD9I98mvMn3ra+Te3pJht0XTo11psmQJdmUSykLq8zzM\nrKSZLTKzdWa2xswe9doLmNkCM9toZp+ZWX6/c/qbWbyZrTezhn7tNc1stZltMrNXL81TEkn/8ufI\nx5Sug0kYGE/TOsV4Ir4ml9/fnZfH7+DUqWBXJ5KCmYeZFQOKOedWmlke4EegBdAR2Ouce97MngIK\nOOf6mVkEMAW4ASgJfA6Ud845M1sK9HTOLTOzecBo59xnZ7mmZh4ifnYf2cMj057nw80TyBnfnkH1\n+tGrczGyZw92ZRJKQmrm4Zzb5Zxb6R0fBtbjC4UWwCSv2ySgpXfcHJjmnDvtnNsCxAO1vBDK65xb\n5vWb7HeOiJxH4dyFmNb5ebb1W0fjRo5BOyModN+TPDt6D8eOBbs6yYwC2oYzs6uAa4HvgKLOuQTw\nBQxQxOtWAtjmd9oOr60EsN2vfbvXJiIpVCxPMaZ3Gs0vfVfTsNkhhuypSJE2gxjx4j6OHAl2dZKZ\npHj7zVuymgH0cs4dNrO/ryul6jrT0KFD/zyOjIwkMjIyNR9eJF0rma8kH3Ycx5b9T9H74xGM/KU8\no6J68fitvXi8Zz7y5Qt2hZIWYmNjiY2NDcq1U3S3lZllAeYAnzrnRntt64FI51yCtyS12Dl3jZn1\nA5xzbpTXbz4wBNia1MdrbwPUc851O8v1tOchEoD4vfE8Nns4C3/5jPDvH6f3zT154tHcFCgQ7Mok\nLYXUnofnHSAuKTg8s4EHvOMOwCy/9jZmls3MygDlgO+9pa0DZlbLfG8per/fOSJyEcpfXp5POr7H\nike/pF6b5bx8piwlW73CkwOPsWdPsKuTjCgld1vdAnwFrMG3NOWAAcD3wHSgFL5ZRWvn3H7vnP5A\nZ+AUvmWuBV77dcBEIAcwzznX6xzX1MxD5CKs2rWKvvOGsGTLMtxXA3johi70eyI7RYsGuzK5lPQi\nQYWHSKr44bcfePLTwSzbuo7ExYPoWPMB+j+ZlRK6VSVDUngoPERS1bfbvuXJz6JZ8+tmTn0xmPbV\n/4/+T2WhdOlgVyapSeGh8BC5JL7c8iX9FkSzacfvnFgwhHsr38uA/mGULRvsyiQ1KDwUHiKXjHOO\nz3/5nP4Lo/l11xGOfTqMlhX/xcABRqVKwa5OLobCQ+Ehcsk555gXP48Bn0ezezccmTOcRlffSfQg\no2rVYFcnF0LhofAQSTPOOWZumMmgLwZz6I/cHJo9nMgr7yB6kFGzZrCrk0AoPBQeImku0SUyfd10\nhiweyqkDhTk4cwQ3FoskOhpq1w52dZISCg+Fh0jQnE48zdQ1UxkaO4ysR67iwMcjqFbgZqKjoU6d\nYFcn56PwUHiIBN2pM6eYtGoSI74cQb4TEez7aATlcl1PdDTUr48+Zz0EKTwUHiIh48TpE0xYMYFn\nvn6GImeu448Zw7kivDrR0dC4sUIklCg8FB4iIefYqWO89eNbPPff57gqrC57YoZS4HQEgwZB8+YK\nkVCg8FB4iISsIyePMOb7Mbz07UtUytqQPTOGkPVQeQYNgnvugbCAPiVIUpPCQ+EhEvIOnjjI6O9G\nM3rpaK7N1ZzdMdGc2l2GgQPh3nshS4o/LUhSi8JD4SGSbuw7to+Xv32ZsT+M5ca8Ufw+YyD7fy3F\ngAHQrh1kzRrsCjMPhYfCQyTd2XN0Dy8seYHxy8cTeXk7Emb057eNxenXDx54ALJnD3aFGV8ofhiU\niMh5FcpViFF3jGJ9j/WULhnO+tsqc+PgJ/jgk92UKwdjxsCxY8GuUlKLwkNEUlXRPEV5pfErrOm2\nhoJFjrGyTkUaPDOAOV/8Qdmy8PLLcORIsKuUi6XwEJFLokS+Erxx5xuseGgFWfL/zg83VaDFy8P4\n8rsDXH01PPcc7NsX7CrlQik8ROSSKn1Zad5u/jbfdfmOI9l/5psbytF27LOsWHeYMmXg/vvh669B\n25zpizbMRSRNbdizgaGxQ4ndEkur8u3JsimK+e/cgGF06eILk8KFg11l+qS7rRQeIhne+t3rmbJm\nCjFxMZw4fYKb8rfiwLetWDK9No0bGQ8+6HsPLb3oMOUUHgoPkUzDOcfa39cSExdDTFwMh44focLp\nVvw6P4rEX2vTpXMYHTtC8eLBrjT0KTwUHiKZknOOdbvXEbPOFyR7Dx+iyO5WbJnXivoVbqLrg2E0\nbgzh4cGuNDQpPBQeIgLE7Y4jZl0MH6yNYee+/WT7+R7c2igevvNmunQO48org11haFF4KDxE5G/W\n715PTFwM7y2P4bd9ezmz9h6qZ4niiXtvoWXzcL0NCgoPhYeInNeGPRt4f9UM3v0+hoRDu8kSfzct\ny0cx+IFbqVgh865pKTwUHiKSQpv2bmJs7Aymrophz4mdFN93Nx1rRzHgvrrkypm5gkThofAQkQuw\nbmc8Iz+ewdzNMRwO20HV8Lvp3TCK9nXrkiUs479HvMJD4SEiF+mLlT8x8sMP+e8fMdhlv3JLwX/x\neJMoGleKzLBBElLhYWYTgGZAgnOumtdWAPgAKA1sAVo75w54f9Yf6AScBno55xZ47TWBiUAOYJ5z\nrvd5rqnwEJFUceoUTPjoF15dMIOfsseQrfBWmpRpycN1ooi8KpKs4Rlnpz3UwuNW4DAw2S88RgF7\nnXPPm9lTQAHnXD8ziwCmADcAJYHPgfLOOWdmS4GezrllZjYPGO2c++wc11R4iEiq27YNXnpnMxOX\nfsjpCjFQ8GfuqdyS+6pHUb9M/XQfJCEVHgBmVhr4xC88NgD1nHMJZlYMiHXOVTKzfoBzzo3y+n0K\nDAW2AouccxFeexvv/G7nuJ7CQ0QumTNnYMECGD1xK1/umUHe2jGczPMTd1duQVREKxpc3YBs4dmC\nXWbA0sOHQRVxziUAOOd2AUW89hLANr9+O7y2EsB2v/btXpuISJoLD4cmTWD+B6XZPOVxHsv3Hfmn\nLWf+pCr/VrjvAAAJiUlEQVR0nz6CYi8Wp+OsjsyLn8fJMyeDXW5ISq23HNM0QUTSpWLFoF8/2Lzy\nSv7Tsw83rv2G06+vZM3C6jw152mKvViMDjM7MGfTHE6cPhHsckPGhd5ykGBmRf2WrX732ncApfz6\nlfTaztV+TkOHDv3zODIyksjIyAssVUQkeWFhvnfxrV8f9uwpxeTJvRn/794UyLGDA60/ZOTvo7j/\n4/tpVqEZURFRNCzbkOxZgvvB7LGxscTGxgbl2ind87gK355HVe/7UcAfzrlR59gwr41vWWohf22Y\nfwc8CiwD5gKvOefmn+N62vMQkaBzDpYsgfHjYdYsqNP0N0o3/pDVZ2JY8/samlVoRqtrWtGoXCNy\nZMkR7HJDa8PczKYCkcDlQAIwBJgJxOCbTWzFd6vufq9/f6AzcIr/vVX3Ov73Vt1e57mmwkNEQsr+\n/TBlii9IDh6E1p13kq/WRyzYEcOqhFU0Ld+UqIgoGpVtRM6sOYNSY0iFRzAoPEQkVDkHP/zgC5GY\nGKhbF1o9sIsDV3zEhxtiWLFzBU3KNyEqIoom5ZqkaZAoPBQeIpIOHD4M06b5guS336BTJ7irbQI/\nHP6YmLgYfvztRxqXa0yriFY0Ld+UXFlzXdJ6FB4KDxFJZ1av9oXI1KlQqxY8+CDUrv87c37yBcmy\n35bRqGwjoiKiaFq+Kbmz5U71GhQeCg8RSaeOHYMZM3xBEh8PHTpAly6Qv/huZm6YSUxcDEt3LKVh\n2Ya0uqYVd1a4kzzZ8qTKtRUeCg8RyQA2bIC334bJk6FKFd9s5F//gsOJe/4Mkm+3fcsdZe8gKiKK\nZhWaXVSQKDwUHiKSgZw44bvVd/x4WLkS2rXzBUlEBOw9updZG2cRExfDN9u+oUGZBrSKaMVdFe4i\nb/a8AV1H4aHwEJEM6pdfYMIEePddKFPGFyKtW0OuXPDHsT+YtWEWM9bP4OutX1O/TH2iIqK4q+Jd\n5MueL9nHVngoPEQkgzt9GubO9c1GvvkG2rTxBUmNGr4/33dsH7M3ziYmLoavtn5F5FWRREVE0bxi\nc/LnyH/Wx1R4KDxEJBPZts03E5kwAQoX9oVI27aQz5ts7D++n9kbZzMjbgaxW2Kpd1W9P4PkshyX\n/fk4Cg+Fh4hkQmfOwMKFvtnIokVw993eLb+1wbxIOHD8AJ9s+oSYuBgWb15MndJ1iIqIokXFFhTM\nVVDhEYp1iYiklV27YOJE391aOXP6QqR9eyhQ4K8+B08c5JONnzBj/QyW/LqE3U/uVniEYl0iImkt\nMRG+/NI3G5k3D5o18wVJ3bp/zUYATpw+QY6sORQeoViXiEgw7dkD773nC5LERN+LDzt08O2TgPY8\nFB4iIufhnO8OrfHjYeZMaNjQNxtp2FDhofAQEUmB/ft976f14YewaJHCQ+EhIhKgtFy2Sq3PMBcR\nkUxE4SEiIgFTeIiISMAUHiIiEjCFh4iIBEzhISIiAVN4iIhIwBQeIiISMIWHiIgETOEhIiIBU3iI\niEjAFB4iIhIwhYeIiAQszcPDzBqb2QYz22RmT6X19UVE5OKlaXiYWRgwBmgEVAbamlmltKwhNcXG\nxga7hGSlhxpBdaY21Zm60kudaSmtZx61gHjn3Fbn3ClgGtAijWtINenhL1R6qBFUZ2pTnakrvdSZ\nltI6PEoA2/y+3+61iYhIOqINcxERCViafgytmd0IDHXONfa+7wc459yov/XTZ9CKiFyADPkZ5mYW\nDmwEGgA7ge+Bts659WlWhIiIXLQsaXkx59wZM+sJLMC3ZDZBwSEikv6k6cxDREQyCOdcqnwBjYEN\nwCbgqXP0eQ2IB1YC1yZ3LlAA3yxlI/AZkN9rLwgsAg4Br/3tGjWB1d5jvRqiNS72HmsFsBwoFMQ6\nbwd+AFYBy4DbUjKWIVZnKI3nDV4dSV8tQ3Q8z1dnyIyn359fie//pcdCcTyTqfOc45nGP/PSwFGv\nhuXA2JSO5VnrSkmnZB/EtwT1k1dcVu9JVvpbnybAXO+4NvBdcucCo4AnveOngOe841zAzUBX/vmL\neSlwg3c8D2gUgjUuBmqEyFhWB4p5x5WB7cmNZQjWGUrjmQMI846LAQl+34fSeJ6vzpAZT7/HjAE+\n4H9/KYfMeCZT51nHMwg/89LA6nP8XM85luf6Sq1bdVPy4r8WwGQA59xSIL+ZFU3m3BbAJO94EtDS\nO/+oc+4b4IT/BcysGJDXObfMa5qcdE6o1OjnXGOf1nWucs7t8o7XATnMLGsyYxkydYbgeB53ziV6\n7TmBREj272bI1OknJMYTwMxaAL8A6/zaQmo8z1Wnn7ONZ5rXCPzjTqwUjOVZpVZ4pOTFf+fqc75z\nizrnEgC8XxxFUlDH9nM8VqjUmGSimS03s0EprOGS12lmrYDl3l/G841lKNWZJGTG08xqmdlafEts\nD3u/pENuPM9RZ5Jgj2dRr8Y8wJPAMP73F1+ojGdydSY523gG4/+hq7w6FpvZrX7XON9YnlUwXyR4\nIfciu1Sv4vwuVY33OeeqAnWAOmbW7gKu4++i6zSzysCz+JbZLpVLVWdIjadz7nvnXBV8+woDzCzb\nRdZzLpeqzlAYz6QgGwK84pw7epE1pERq1un/WKk5nhfzM98JXOmcqwk8Dkz1Qu+CpFZ47MC3UZSk\npNf29z6lztLnfOfu8qZoSVOr31NQx9muEUo14pzb6f33CDAV3xTUv4Y0rdPMSgIfAe2dc1uSuUao\n1Rly4+lX10bgMFDlPNcItTpDbTxrA8+b2S9Ab3wh1/081wiVOvt7dZ5vPNO0RufcSefcPu94OfAz\nUOE81zi/5DZFUvIFhPPX5k02fJs31/ytT1P+2vi5kb82fs55Lr6Nn6fcuTenOgCv/63tO3w/HMO3\n8dM4lGr0Huty7zgrvg22rsEaS+Ayr19L/2ucbyxDqc4QHM+rgHD31wbldqBgCI7nWesMtfH82+MO\n4X83okNmPM9V5/nGMwg/80L8dVPE1fiWvS5LbizP+Xs/uQ4p/cJ329hGfLeU9fPaHvrbX7wx3hNe\nBdQ837lee0Hgc+/PFiQ9Ue/PNgN7gIPAr/x1p8F1wBrvsUaHWo347sL6wfthrwFewXu9TTDqBAbi\nu7VwOX+7lfB8YxkqdYbgeLYD1nr1/QDc5XdOKI3nWesMtfE81y/lUBvPc9WZ3Him8c/87r/9zJum\ndCzP9qUXCYqISMD0rroiIhIwhYeIiARM4SEiIgFTeIiISMAUHiIiEjCFh4iIBEzhISIiAVN4iIhI\nwP4fRawciCCNV+wAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0xdc1390>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(distances, counts)\n",
    "plt.plot(distances, c(distances, N, D))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "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.4.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
	{
	"cells": [
	{
	"cell_type": "code",
	"execution_count": 30,
	"metadata": {
	"collapsed": false
	},
	"outputs": [],
	"source": [
	"from math import pi\n",
	"from numpy import sqrt, exp, array\n",
	"from units import m\n",
	"from matplotlib import pyplot as plt\n",
	"from scipy.optimize import curve_fit\n",
	"%matplotlib inline "
	]
	},
	{
	"cell_type": "code",
	"execution_count": 31,
	"metadata": {
	"collapsed": false
	},
	"outputs": [],
	"source": [
	"# Define unit of cm\n",
	"cm = 1e-2 * m"
	]
	},
	{
	"cell_type": "code",
	"execution_count": 32,
	"metadata": {
	"collapsed": false
	},
	"outputs": [],
	"source": [
	"# Load data\n",
	"data = [(5012, 1), (3981, 2), (2512, 3), (1413, 4), (525, 5)]\n",
	"counts, distances = zip(*data)\n",
	"\n",
	"counts = array(counts)\n",
	"distances = array(distances) # Convert to meters\n",
	"\n",
	"t = 20 * 60 * 60"
	]
	},
	{
	"cell_type": "code",
	"execution_count": 33,
	"metadata": {
	"collapsed": false
	},
	"outputs": [],
	"source": [
	"# Define count function\n",
	"def count(x, N, D):\n",
	" return N / sqrt(4* pi * D * t) * exp(-x*2 / (4 D * t))"
	]
	},
	{
	"cell_type": "code",
	"execution_count": 34,
	"metadata": {
	"collapsed": false
	},
	"outputs": [
	{
	"name": "stderr",
	"output_type": "stream",
	"text": [
	"C:\\python34\\lib\\site-packages\\ipykernel\\__main__.py:3: RuntimeWarning: invalid value encountered in sqrt\n",
	" app.launch_new_instance()\n"
	]
	}
	],
	"source": [
	"(N, D), pcov = curve_fit(c, distances, counts)"
	]
	},
	{
	"cell_type": "code",
	"execution_count": 35,
	"metadata": {
	"collapsed": false
	},
	"outputs": [
	{
	"data": {
	"image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAEACAYAAABLfPrqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmcjeX/x/HXZ8a+RtaQZG1sUdGGibIlfMuIfiRLylKo\nlG3sLdqVqK8U+pKMCiFRTPVVSdmNZSqEmJB9Z67fH+eeOt/CzGHMOTPzfj4e8+iey3Wf+3Ou0bxd\n13Wfc8w5h4iISCDCgl2AiIikPwoPEREJmMJDREQCpvAQEZGAKTxERCRgCg8REQlYisLDzPKbWYyZ\nrTezdWZW28wKmNkCM9toZp+ZWX6//v3NLN7r39CvvaaZrTazTWb26qV4QiIicumldOYxGpjnnLsG\nqA5sAPoBnzvnKgKLgP4AZhYBtAauAZoAY83MvMcZB3R2zlUAKphZo1R7JiIikmaSDQ8zywfUcc69\nC+CcO+2cOwC0ACZ53SYBLb3j5sA0r98WIB6oZWbFgLzOuWVev8l+54iISDqSkplHGWCPmb1rZsvN\n7N9mlgso6pxLAHDO7QKKeP1LANv8zt/htZUAtvu1b/faREQknUlJeGQBagJvOOdqAkfwLVn9/X1N\n9D4nIiKZRJYU9NkObHPO/eB9/yG+8Egws6LOuQRvSep37893AKX8zi/ptZ2r/R/MTEEkInIBnHOW\nfK+Ll+zMw1ua2mZmFbymBsA6YDbwgNfWAZjlHc8G2phZNjMrA5QDvveWtg6YWS1vA/1+v3POdt2Q\n/xoyZEjQa8gINapO1RnqX+mlzrSUkpkHwKPAFDPLCvwCdATCgelm1gnYiu8OK5xzcWY2HYgDTgHd\n3V/PqgcwEciB7+6t+an1REREJO2kKDycc6uAG87yR7efo/+zwLNnaf8RqBpIgSIiEnr0CvOLEBkZ\nGewSkpUeagTVmdpUZ+pKL3WmJUvrdbKUMDMXinWJiIQyM8OFyoa5iIjI3yk8REQkYAoPEREJmMJD\nREQCpvAQEZGAKTxERCRgCg8REQmYwkNERAKm8BARkYApPEREJGAKDxERCZjCQ0REAqbwEBGRgCk8\nREQkYAoPEREJmMJDREQCpvAQEZGAKTxERCRgCg8REQmYwkNERAKm8BARkYApPEREJGAKDxERCZjC\nQ0REAqbwEBGRgCk8REQkYCkKDzPbYmarzGyFmX3vtRUwswVmttHMPjOz/H79+5tZvJmtN7OGfu01\nzWy1mW0ys1dT/+mIiEhaSOnMIxGIdM7VcM7V8tr6AZ875yoCi4D+AGYWAbQGrgGaAGPNzLxzxgGd\nnXMVgApm1iiVnoeIiKShlIaHnaVvC2CSdzwJaOkdNwemOedOO+e2APFALTMrBuR1zi3z+k32O+cf\nPo3/FOdcCssTEZG0lNLwcMBCM1tmZl28tqLOuQQA59wuoIjXXgLY5nfuDq+tBLDdr32713ZWjy94\nnDrv1iF2S2wKSxQRkbSSJYX9bnHO7TSzwsACM9uIL1D8peo04Ylcawir/j5dZnehTIEyjLxtJLVL\n1k7NS4iIyAVKUXg453Z6/91tZjOBWkCCmRV1ziV4S1K/e913AKX8Ti/ptZ2r/axGjhjBnj1Qq3Zb\nKjU/RlRMFNcWu5bhtw3n2mLXpvgJiohkVLGxscTGxgbl2pbcvoKZ5QLCnHOHzSw3sAAYBjQA/nDO\njTKzp4ACzrl+3ob5FKA2vmWphUB555wzs++AR4FlwFzgNefc/LNc0znnWLECXngBPvsMOnQ+zmUN\n/s24tc9St3RdhkUOo1KhSqk3EiIi6ZyZ4Zyz5HtevJTseRQF/mtmK4DvgE+ccwuAUcAd3hJWA+A5\nAOdcHDAdiAPmAd3dXwnVA5gAbALizxYc/mrUgKlT4ccfIfFkDl5t+yj11/5EcatBnXfr8MDMB/hl\n3y+BP2sREbkoyc48giFp5vF3f/wB48bB669DtVoHuOLuV5jz+xiiIqIYVHcQJfKdc/9dRCTDS8uZ\nR7oKjyTHj8PkyfDSS5Dz8j2Uvu8F/nv0bTpU70C/W/tRJHeRc54rIpJRKTySCY8kiYkwe7ZvX2T7\ngZ2Uaf8sa5jCw9c/zBM3P0GBnAXSoFoRkdCg8EhhePhbssQXIl+v2Urp9iP5NedMet/Ui161e5E3\ne95LVKmISOhQeFxAeCTZuNG3nDVtYTzF2gzlj8s+p3+dJ+l+Q3dyZs2ZypWKiIQOhcdFhEeSXbt8\nG+tvxKwlT7MhnCzyHUMbDKRLzS5kC8+WSpWKiIQOhUcqhEeSQ4dgwgQYNflHTtw6iGzFN/B0o8F0\nuLY9WcJS+gJ7EZHQp/BIxfBIcuoUxMTA4An/ZWelQeS7Yicv3DmM+6q3Jsz0sSYikv4pPC5BeCRx\nDhYudPQf/wVrigzk8qLHeKnZCNrWaM5f7xwvIpL+KDwuYXj4W77c0efNOSzJNojCBbPzYtOR3Ff7\nDoWIiKRLCo80Co8km7ck0nPsDOafGEzhnEUZ1XgkHSLrpNn1RURSg8IjjcMjye69p+k2bgoz9w2l\nQGJFht82gofvugFNREQkPVB4BCk8khw8cpLu49/hg50jyXPwBgbeNJze91Uli27OEpEQpvAIcngk\nOXLiGI9OfpP3No8i2476PHbdUPp1qUCuXMGuTETknxQeIRIeSQ6fPMzjMa8xccMrhG1qTtdKgxnU\nszSFCwe7MhGRvyg8Qiw8kuw/vp8Bn7zEO2vHwpq23HvFAKJ7X0G5csGuTERE4RGy4ZFk95HdDF44\nikkr34UVnbg951MMeqwQtWoFuzIRycwUHiEeHkl2HNzB0EVPM3X1B2Rd0YPKBx+jf5/LaNoUwvSi\ndRFJYwqPdBIeSTbv28yw2BF8tO4T8qx+jPwbH+HJ3nm47z7Inj3Y1YlIZqHwSGfhkWTjno0MiR3C\ngk2xFNnUj4OLHqZ3zxw89BDkzx/s6kQko1N4pNPwSLI6YTXRi6NZunU5pbcOIj6mI506ZKN3byhZ\nMtjViUhGlZbhoZX5S6Ba0WrMajOL2e0+5LKbPyJPv0qsDZ9M1epn6NAB1qwJdoUiIhdHM4808NXW\nrxi4aCAJh/Zw/aHhLBpzDzVrhNG3L0RGorc/EZFUoWWrDBYeAM45Fvy8gEGLB3Hq9GnqnB7Bgjfu\nJG8eo29fuOce9PYnInJRFB4ZMDySOOeYtXEW0YujyZM1D01zjmTBmw3YsQMeeww6doTcuYNdpYik\nRwqPDBweSc4knmH6uukMiR1CyXwlaVNkJPPH38x//wvdukHPnujtT0QkINowzwTCw8JpW7UtcT3i\naFetHU9vbMuJVnfy1uzl7NoFFSv6QuSnn4JdqYjIPyk8gixLWBY61ejEpp6baFquKT2/uYu9DVrx\n8ZJ1FCoEN90ErVrB0qXBrlRE5C8KjxCRPUt2etTqQfwj8dxY8kZaz63P5hrt+Hz5T9StC/feC/Xq\nwZw5kJgY7GpFJLNTeISYXFlz8cTNTxD/SDwVL69Ig/dvZG2Zriz68VcefhgGD4YqVeCdd+DEiWBX\nKyKZVYrDw8zCzGy5mc32vi9gZgvMbKOZfWZm+f369jezeDNbb2YN/dprmtlqM9tkZq+m7lPJWPJl\nz0d0vWg2PbKJQrkKccOEGnx3WS/mfrmL11+H6dPh6qth1CjYvz/Y1YpIZhPIzKMXEOf3fT/gc+dc\nRWAR0B/AzCKA1sA1QBNgrNmfL4MbB3R2zlUAKphZo4usP8MrmLMgzzR4hrjucYRZGFXGVWah68eU\nj/Yybx6sXQtly8ITT8D27cGuVkQyixSFh5mVBJoCb/s1twAmeceTgJbecXNgmnPutHNuCxAP1DKz\nYkBe59wyr99kv3MkGUXzFOWVxq+w8qGV7D++n4pjKjJz3zDeePsgK1b49kGqVYP779fbn4jIpZfS\nmccrQF/A/8UXRZ1zCQDOuV1AEa+9BLDNr98Or60E4P9v4+1emwSgVP5SvNnsTZZ2WcrP+36m3Gvl\nmLbteUY8d4Sff4aICGjUCJo0gUWLIIO/XEZEgiTZN8QwszuBBOfcSjOLPE/XVP01NXTo0D+PIyMj\niYw836Uzn7IFyzL5X5OJ2x3H0NihlH+9PP1v7U+fJ7rSp092/vMf6NED8uWDN9+EGjWCXbGIpLbY\n2FhiY2ODcu1kX2FuZs8A7YDTQE4gL/AxcD0Q6ZxL8JakFjvnrjGzfoBzzo3yzp8PDAG2JvXx2tsA\n9Zxz3c5yzQz/CvPUtmLnCqIXR7M6YTWD6w2mQ/UOhFtW3nsP+vaFBx+E6GjIkSPYlYrIpRJSrzB3\nzg1wzl3pnLsaaAMscs61Bz4BHvC6dQBmecezgTZmls3MygDlgO+9pa0DZlbL20C/3+8cuUg1itdg\nzn1z+KDVB0xbO42IsRG8v3YK7dqfYfVq2LDBN/v49ttgVyoiGUFA721lZvWAx51zzc2sIDAdKIVv\nVtHaObff69cf6AycAno55xZ47dcBE4EcwDznXK9zXEczj4u0ePNiBi4ayKGThxjTZAz1rqrHjBnw\n6KPQujU8/bTegFEko9EbIyo8UoVzjo/Wf0Tvz3pTv0x9XrjjBcKPF6F3b1iyBMaPhwYNgl2liKSW\nkFq2kvTLzLgn4h7iusdROFdhqoytwozNbzFpciJjxvje/r1rVzhwINiVikh6o/DIBPJmz8uLDV/k\n8/s/Z/Lqydw04SaK11jB2rUQFuZ7u5M5c4JdpYikJ1q2ymQSXSITV06k/xf9aVO5DcNvG87yb/Pz\n4INQuzaMHg2FCgW7ShG5EFq2kksmzMLoVKMTcd3jOHrqKBFjI0goPI1VqxzFikHVqr73zVJ2i8j5\naOaRyX2z7Ru6ze1GkdxFeKPpG+z7qQKdOkGFCjB2LBQvHuwKRSSlNPOQNHNzqZv5seuPNC3XlJsn\n3Mzco4NZsvQYVapA9erw7ruahYjIP2nmIX/afnA7fT7rw/KdyxnTZAxXHG1Cp06+PZC33oKrrgp2\nhSJyPpp5SFCUzFeSmKgY3mj6Bj0/7cmIja34cOF2brsNrr8exozRpxiKiI/CQ/6hcbnGrO22lsqF\nK3P929eSte5LLP7qFO+/7/so3E2bgl2hiASblq3kvOL3xtNjXg92Hd7FmCbjWDn7FoYPhyefhMce\ngyzJvi+ziKQVvT2JwiOkOOeYvm46jy14jCblmtCt/HM89Ugh9u/3fZZ6tWrBrlBEQHseEmLMjHur\n3Mv6HuvJky0PTedWps3zE3jo4URuvx2GDIGTJ4NdpYikJc08JGArdq6g29xuhIeFM/T6cYyJrsbP\nP/tmIbVqBbs6kcxLy1YKj5CX6BJ5e/nbDFo0iHbV2lN171D6P5aX9u1h2DDIlSvYFYpkPlq2kpAX\nZmF0va4ra7uvZd/xPxicEMEzM2ewfYejenX48stgVygil5JmHpIqvt76Nd3mdqNU/lK0zDaGEX3K\n0rw5PPec73PUReTS08xD0p06peuw4qEV1L+qPgM31+b+CcM5dvIEVavC/PnBrk5EUptmHpLqfj3w\nK73m92Ld7+voUnws4/reTt268MorULBgsKsTybi0Ya7wyBDmbJrDI58+wnVFa5Pv25eZH3MFr78O\n99wT7MpEMiYtW0mG0KxCM9Z1X0fFwmWZXaIaUS+NZsCg00RFQUJCsKsTkYuh8JBLKlfWXDzd4Gm+\n7vg1a07OIlevWuQqv5Rq1eC99/R27yLplZatJM0455i6Zip9F/blpoLN2fDGM5QuUpC33oJSpYJd\nnUj6p2UryZDMjP+r9n/E9YijeNEs7G0TQc4bJ1GjpuOtt/R27yLpiWYeEjQ//PYDD895GE7l4uTH\n47j8TGXGj4dy5YJdmUj6pJmHZArXX3E9S7sspVOtNuxsFIk1fIratx7hpZfgzJlgVyci56OZh4SE\nXYd30XdhXxb99BUFvh9N7m0teGeCUblysCsTST/0Og+FR6a1ePNius/tTtbD5dj+79fp0/Eq+vWD\nrFmDXZlI6NOylWRat5W5jVXdVtHmlpug6/X8Z+sz1LzhJD/+GOzKRMRfsuFhZtnNbKmZrTCzNWY2\nxGsvYGYLzGyjmX1mZvn9zulvZvFmtt7MGvq11zSz1Wa2ycxevTRPSdK7bOHZGFBnAD8+tIwK9b9l\nb+vq3N5lMf37w/Hjwa5ORCCFy1Zmlss5d9TMwoElwKPAPcBe59zzZvYUUMA518/MIoApwA1ASeBz\noLxzzpnZUqCnc26Zmc0DRjvnPjvL9bRsJYDvtSGzNs6i59xehG+vS9bFLzLpjaLcckuwKxMJPSG3\nbOWcO+odZgeyAA5oAUzy2icBLb3j5sA059xp59wWIB6oZWbFgLzOuWVev8l+54iclZnRslJLNj4S\nR5umV7C7VRWaDhlLz0fPcPhwsKsTybxSFB5mFmZmK4BdwEIvAIo65xIAnHO7gCJe9xLANr/Td3ht\nJYDtfu3bvTaRZOXOlptRd4xiSdfFVL53GtNy30iFyB9YuDDYlYlkTllS0sk5lwjUMLN8wMdmVhnf\n7ON/uqVmYUOHDv3zODIyksjIyNR8eEmnqhSpwpIuXzJ51WT65G3G3W+3osUHIxnz4mVcdlmwqxNJ\nW7GxscTGxgbl2gHfqmtm0cBRoAsQ6ZxL8JakFjvnrjGzfoBzzo3y+s8HhgBbk/p47W2Aes65bme5\nhvY8JFl/HPuDvvMHMG3FbLJ/+QLv9LmPli3TZLlXJCSF1J6HmRVKupPKzHICdwDrgdnAA163DsAs\n73g20MbMsplZGaAc8L23tHXAzGqZmQH3+50jErCCOQsy4V9v8kXXj7i8+YvcN/927uywgd27g12Z\nSMaX7MzDzKri2xAP874+cM49bWYFgelAKXyzitbOuf3eOf2BzsApoJdzboHXfh0wEcgBzHPO9TrH\nNTXzkICcTjzNq0vGMviLEYQt78qYewfS4b5cmCYikonoFeYKD7lAvx36jQ7vP8aX8d9TY+frfDTq\nTkrotgzJJBQeCg+5SHM3LOT+aT04/Etlht88micfulKzEMnwQmrPQyQ9urPSHewYtJrOzWowcFtN\nKnR8gU0/nwp2WSIZhmYekuFt3P0zzd/syc+7t9Gj9Dhe6VOHMP2zSTIgLVspPCSVOed4/YsP6ftF\nH/LvbcDM7i9w87WFg12WSKrSspVIKjMzHr29Fb8PjqNK2cu59f3K/GvkvzlxUp99K3IhNPOQTOnT\n5atpO+VhTp12vHP3OO6td22wSxK5aFq2UnhIGjiTmEiXN95l0vYBXJ+9LXOfGE7hfPmCXZbIBdOy\nlUgaCA8L491HOrP64XUk7D/EFU9HMPLjD9A/XESSp5mHCOAcjJi0hBHLu1EiX3FmdxtDtRLlg12W\nSEC0bKXwkCDZmXCKZiNeZ2XeZ2hXoSdvte9Hjiw5gl2WSIooPBQeEmQTP9pO99m9yVpqJe9GvcHd\n1RoFuySRZCk8FB4SAg4cgDaDPmVhtp7cWLomHzzwKiXy6Y2yJHRpw1wkBOTPD5++3oRPGq9lXew1\nlH2xOoMXjOLoqaPJnyySwWnmIZICR45Aj8GbmJowkOzlltCvzlM8Xu8h7YdISNGylcJDQtTmzfDE\niyv55NAQspdZzqC6A+gT2Zls4dmCXZqIwkPhIaFu2zZ4/OVlfLx/MDlLrWdwvWh6Rd5P1vCswS5N\nMjGFh8JD0omdO6HPK0v48I/B5Cq+laG3DeHRyPsIDwsPdmmSCSk8FB6SzuzeDb1eXUzMnmhyF97L\n8NuG0bN+K8JM96RI2lF4KDwkndq71/Hoawv4YHc0eS87wYj6w+jRoAWmjzGUNKDwUHhIOrd/v+OR\n1+cwLWEwefOEM7LBcLrd3kQhIpeUwkPhIRnEwUOJ9BzzMVN3DSFf9rw83WAEDzdsoBCRS0LhofCQ\nDObQ4TP0GDedqb8NJX9YMUY2GEG3pnWDXZZkMAoPhYdkUIePnqbHm1OY8tsw8p8u6wuRu24MdlmS\nQSg8FB6SwR05dopu/57I1B0jyH+sKiMaDKdbi+vQapZcDIWHwkMyiSPHT9Dt7beZuu0Z8h2qxcj6\nw+l2T1WFiFwQhYfCQzKZw8eP0e2dcUzb9jx590YyvP5QureuRJheJiIBUHgoPCSTOnj8MN0njuGD\nX18mz67GDG8wmB5tyylEJEUUHgoPyeQOHD9I98mvMn3ra+Te3pJht0XTo11psmQJdmUSykLq8zzM\nrKSZLTKzdWa2xswe9doLmNkCM9toZp+ZWX6/c/qbWbyZrTezhn7tNc1stZltMrNXL81TEkn/8ufI\nx5Sug0kYGE/TOsV4Ir4ml9/fnZfH7+DUqWBXJ5KCmYeZFQOKOedWmlke4EegBdAR2Ouce97MngIK\nOOf6mVkEMAW4ASgJfA6Ud845M1sK9HTOLTOzecBo59xnZ7mmZh4ifnYf2cMj057nw80TyBnfnkH1\n+tGrczGyZw92ZRJKQmrm4Zzb5Zxb6R0fBtbjC4UWwCSv2ySgpXfcHJjmnDvtnNsCxAO1vBDK65xb\n5vWb7HeOiJxH4dyFmNb5ebb1W0fjRo5BOyModN+TPDt6D8eOBbs6yYwC2oYzs6uAa4HvgKLOuQTw\nBQxQxOtWAtjmd9oOr60EsN2vfbvXJiIpVCxPMaZ3Gs0vfVfTsNkhhuypSJE2gxjx4j6OHAl2dZKZ\npHj7zVuymgH0cs4dNrO/ryul6jrT0KFD/zyOjIwkMjIyNR9eJF0rma8kH3Ycx5b9T9H74xGM/KU8\no6J68fitvXi8Zz7y5Qt2hZIWYmNjiY2NDcq1U3S3lZllAeYAnzrnRntt64FI51yCtyS12Dl3jZn1\nA5xzbpTXbz4wBNia1MdrbwPUc851O8v1tOchEoD4vfE8Nns4C3/5jPDvH6f3zT154tHcFCgQ7Mok\nLYXUnofnHSAuKTg8s4EHvOMOwCy/9jZmls3MygDlgO+9pa0DZlbLfG8per/fOSJyEcpfXp5POr7H\nike/pF6b5bx8piwlW73CkwOPsWdPsKuTjCgld1vdAnwFrMG3NOWAAcD3wHSgFL5ZRWvn3H7vnP5A\nZ+AUvmWuBV77dcBEIAcwzznX6xzX1MxD5CKs2rWKvvOGsGTLMtxXA3johi70eyI7RYsGuzK5lPQi\nQYWHSKr44bcfePLTwSzbuo7ExYPoWPMB+j+ZlRK6VSVDUngoPERS1bfbvuXJz6JZ8+tmTn0xmPbV\n/4/+T2WhdOlgVyapSeGh8BC5JL7c8iX9FkSzacfvnFgwhHsr38uA/mGULRvsyiQ1KDwUHiKXjHOO\nz3/5nP4Lo/l11xGOfTqMlhX/xcABRqVKwa5OLobCQ+Ehcsk555gXP48Bn0ezezccmTOcRlffSfQg\no2rVYFcnF0LhofAQSTPOOWZumMmgLwZz6I/cHJo9nMgr7yB6kFGzZrCrk0AoPBQeImku0SUyfd10\nhiweyqkDhTk4cwQ3FoskOhpq1w52dZISCg+Fh0jQnE48zdQ1UxkaO4ysR67iwMcjqFbgZqKjoU6d\nYFcn56PwUHiIBN2pM6eYtGoSI74cQb4TEez7aATlcl1PdDTUr48+Zz0EKTwUHiIh48TpE0xYMYFn\nvn6GImeu448Zw7kivDrR0dC4sUIklCg8FB4iIefYqWO89eNbPPff57gqrC57YoZS4HQEgwZB8+YK\nkVCg8FB4iISsIyePMOb7Mbz07UtUytqQPTOGkPVQeQYNgnvugbCAPiVIUpPCQ+EhEvIOnjjI6O9G\nM3rpaK7N1ZzdMdGc2l2GgQPh3nshS4o/LUhSi8JD4SGSbuw7to+Xv32ZsT+M5ca8Ufw+YyD7fy3F\ngAHQrh1kzRrsCjMPhYfCQyTd2XN0Dy8seYHxy8cTeXk7Emb057eNxenXDx54ALJnD3aFGV8ofhiU\niMh5FcpViFF3jGJ9j/WULhnO+tsqc+PgJ/jgk92UKwdjxsCxY8GuUlKLwkNEUlXRPEV5pfErrOm2\nhoJFjrGyTkUaPDOAOV/8Qdmy8PLLcORIsKuUi6XwEJFLokS+Erxx5xuseGgFWfL/zg83VaDFy8P4\n8rsDXH01PPcc7NsX7CrlQik8ROSSKn1Zad5u/jbfdfmOI9l/5psbytF27LOsWHeYMmXg/vvh669B\n25zpizbMRSRNbdizgaGxQ4ndEkur8u3JsimK+e/cgGF06eILk8KFg11l+qS7rRQeIhne+t3rmbJm\nCjFxMZw4fYKb8rfiwLetWDK9No0bGQ8+6HsPLb3oMOUUHgoPkUzDOcfa39cSExdDTFwMh44focLp\nVvw6P4rEX2vTpXMYHTtC8eLBrjT0KTwUHiKZknOOdbvXEbPOFyR7Dx+iyO5WbJnXivoVbqLrg2E0\nbgzh4cGuNDQpPBQeIgLE7Y4jZl0MH6yNYee+/WT7+R7c2igevvNmunQO48org11haFF4KDxE5G/W\n715PTFwM7y2P4bd9ezmz9h6qZ4niiXtvoWXzcL0NCgoPhYeInNeGPRt4f9UM3v0+hoRDu8kSfzct\ny0cx+IFbqVgh865pKTwUHiKSQpv2bmJs7Aymrophz4mdFN93Nx1rRzHgvrrkypm5gkThofAQkQuw\nbmc8Iz+ewdzNMRwO20HV8Lvp3TCK9nXrkiUs479HvMJD4SEiF+mLlT8x8sMP+e8fMdhlv3JLwX/x\neJMoGleKzLBBElLhYWYTgGZAgnOumtdWAPgAKA1sAVo75w54f9Yf6AScBno55xZ47TWBiUAOYJ5z\nrvd5rqnwEJFUceoUTPjoF15dMIOfsseQrfBWmpRpycN1ooi8KpKs4Rlnpz3UwuNW4DAw2S88RgF7\nnXPPm9lTQAHnXD8ziwCmADcAJYHPgfLOOWdmS4GezrllZjYPGO2c++wc11R4iEiq27YNXnpnMxOX\nfsjpCjFQ8GfuqdyS+6pHUb9M/XQfJCEVHgBmVhr4xC88NgD1nHMJZlYMiHXOVTKzfoBzzo3y+n0K\nDAW2AouccxFeexvv/G7nuJ7CQ0QumTNnYMECGD1xK1/umUHe2jGczPMTd1duQVREKxpc3YBs4dmC\nXWbA0sOHQRVxziUAOOd2AUW89hLANr9+O7y2EsB2v/btXpuISJoLD4cmTWD+B6XZPOVxHsv3Hfmn\nLWf+pCr/VrjvAAAJiUlEQVR0nz6CYi8Wp+OsjsyLn8fJMyeDXW5ISq23HNM0QUTSpWLFoF8/2Lzy\nSv7Tsw83rv2G06+vZM3C6jw152mKvViMDjM7MGfTHE6cPhHsckPGhd5ykGBmRf2WrX732ncApfz6\nlfTaztV+TkOHDv3zODIyksjIyAssVUQkeWFhvnfxrV8f9uwpxeTJvRn/794UyLGDA60/ZOTvo7j/\n4/tpVqEZURFRNCzbkOxZgvvB7LGxscTGxgbl2ind87gK355HVe/7UcAfzrlR59gwr41vWWohf22Y\nfwc8CiwD5gKvOefmn+N62vMQkaBzDpYsgfHjYdYsqNP0N0o3/pDVZ2JY8/samlVoRqtrWtGoXCNy\nZMkR7HJDa8PczKYCkcDlQAIwBJgJxOCbTWzFd6vufq9/f6AzcIr/vVX3Ov73Vt1e57mmwkNEQsr+\n/TBlii9IDh6E1p13kq/WRyzYEcOqhFU0Ld+UqIgoGpVtRM6sOYNSY0iFRzAoPEQkVDkHP/zgC5GY\nGKhbF1o9sIsDV3zEhxtiWLFzBU3KNyEqIoom5ZqkaZAoPBQeIpIOHD4M06b5guS336BTJ7irbQI/\nHP6YmLgYfvztRxqXa0yriFY0Ld+UXFlzXdJ6FB4KDxFJZ1av9oXI1KlQqxY8+CDUrv87c37yBcmy\n35bRqGwjoiKiaFq+Kbmz5U71GhQeCg8RSaeOHYMZM3xBEh8PHTpAly6Qv/huZm6YSUxcDEt3LKVh\n2Ya0uqYVd1a4kzzZ8qTKtRUeCg8RyQA2bIC334bJk6FKFd9s5F//gsOJe/4Mkm+3fcsdZe8gKiKK\nZhWaXVSQKDwUHiKSgZw44bvVd/x4WLkS2rXzBUlEBOw9updZG2cRExfDN9u+oUGZBrSKaMVdFe4i\nb/a8AV1H4aHwEJEM6pdfYMIEePddKFPGFyKtW0OuXPDHsT+YtWEWM9bP4OutX1O/TH2iIqK4q+Jd\n5MueL9nHVngoPEQkgzt9GubO9c1GvvkG2rTxBUmNGr4/33dsH7M3ziYmLoavtn5F5FWRREVE0bxi\nc/LnyH/Wx1R4KDxEJBPZts03E5kwAQoX9oVI27aQz5ts7D++n9kbZzMjbgaxW2Kpd1W9P4PkshyX\n/fk4Cg+Fh4hkQmfOwMKFvtnIokVw993eLb+1wbxIOHD8AJ9s+oSYuBgWb15MndJ1iIqIokXFFhTM\nVVDhEYp1iYiklV27YOJE391aOXP6QqR9eyhQ4K8+B08c5JONnzBj/QyW/LqE3U/uVniEYl0iImkt\nMRG+/NI3G5k3D5o18wVJ3bp/zUYATpw+QY6sORQeoViXiEgw7dkD773nC5LERN+LDzt08O2TgPY8\nFB4iIufhnO8OrfHjYeZMaNjQNxtp2FDhofAQEUmB/ft976f14YewaJHCQ+EhIhKgtFy2Sq3PMBcR\nkUxE4SEiIgFTeIiISMAUHiIiEjCFh4iIBEzhISIiAVN4iIhIwBQeIiISMIWHiIgETOEhIiIBU3iI\niEjAFB4iIhIwhYeIiAQszcPDzBqb2QYz22RmT6X19UVE5OKlaXiYWRgwBmgEVAbamlmltKwhNcXG\nxga7hGSlhxpBdaY21Zm60kudaSmtZx61gHjn3Fbn3ClgGtAijWtINenhL1R6qBFUZ2pTnakrvdSZ\nltI6PEoA2/y+3+61iYhIOqINcxERCViafgytmd0IDHXONfa+7wc459yov/XTZ9CKiFyADPkZ5mYW\nDmwEGgA7ge+Bts659WlWhIiIXLQsaXkx59wZM+sJLMC3ZDZBwSEikv6k6cxDREQyCOdcqnwBjYEN\nwCbgqXP0eQ2IB1YC1yZ3LlAA3yxlI/AZkN9rLwgsAg4Br/3tGjWB1d5jvRqiNS72HmsFsBwoFMQ6\nbwd+AFYBy4DbUjKWIVZnKI3nDV4dSV8tQ3Q8z1dnyIyn359fie//pcdCcTyTqfOc45nGP/PSwFGv\nhuXA2JSO5VnrSkmnZB/EtwT1k1dcVu9JVvpbnybAXO+4NvBdcucCo4AnveOngOe841zAzUBX/vmL\neSlwg3c8D2gUgjUuBmqEyFhWB4p5x5WB7cmNZQjWGUrjmQMI846LAQl+34fSeJ6vzpAZT7/HjAE+\n4H9/KYfMeCZT51nHMwg/89LA6nP8XM85luf6Sq1bdVPy4r8WwGQA59xSIL+ZFU3m3BbAJO94EtDS\nO/+oc+4b4IT/BcysGJDXObfMa5qcdE6o1OjnXGOf1nWucs7t8o7XATnMLGsyYxkydYbgeB53ziV6\n7TmBREj272bI1OknJMYTwMxaAL8A6/zaQmo8z1Wnn7ONZ5rXCPzjTqwUjOVZpVZ4pOTFf+fqc75z\nizrnEgC8XxxFUlDH9nM8VqjUmGSimS03s0EprOGS12lmrYDl3l/G841lKNWZJGTG08xqmdlafEts\nD3u/pENuPM9RZ5Jgj2dRr8Y8wJPAMP73F1+ojGdydSY523gG4/+hq7w6FpvZrX7XON9YnlUwXyR4\nIfciu1Sv4vwuVY33OeeqAnWAOmbW7gKu4++i6zSzysCz+JbZLpVLVWdIjadz7nvnXBV8+woDzCzb\nRdZzLpeqzlAYz6QgGwK84pw7epE1pERq1un/WKk5nhfzM98JXOmcqwk8Dkz1Qu+CpFZ47MC3UZSk\npNf29z6lztLnfOfu8qZoSVOr31NQx9muEUo14pzb6f33CDAV3xTUv4Y0rdPMSgIfAe2dc1uSuUao\n1Rly4+lX10bgMFDlPNcItTpDbTxrA8+b2S9Ab3wh1/081wiVOvt7dZ5vPNO0RufcSefcPu94OfAz\nUOE81zi/5DZFUvIFhPPX5k02fJs31/ytT1P+2vi5kb82fs55Lr6Nn6fcuTenOgCv/63tO3w/HMO3\n8dM4lGr0Huty7zgrvg22rsEaS+Ayr19L/2ucbyxDqc4QHM+rgHD31wbldqBgCI7nWesMtfH82+MO\n4X83okNmPM9V5/nGMwg/80L8dVPE1fiWvS5LbizP+Xs/uQ4p/cJ329hGfLeU9fPaHvrbX7wx3hNe\nBdQ837lee0Hgc+/PFiQ9Ue/PNgN7gIPAr/x1p8F1wBrvsUaHWo347sL6wfthrwFewXu9TTDqBAbi\nu7VwOX+7lfB8YxkqdYbgeLYD1nr1/QDc5XdOKI3nWesMtfE81y/lUBvPc9WZ3Him8c/87r/9zJum\ndCzP9qUXCYqISMD0rroiIhIwhYeIiARM4SEiIgFTeIiISMAUHiIiEjCFh4iIBEzhISIiAVN4iIhI\nwP4fRawciCCNV+wAAAAASUVORK5CYII=\n",
	"text/plain": [
	"<matplotlib.figure.Figure at 0xdc1390>"
	]
	},
	"metadata": {},
	"output_type": "display_data"
	}
	],
	"source": [
	"plt.plot(distances, counts)\n",
	"plt.plot(distances, c(distances, N, D))\n",
	"plt.show()"
	]
	},
	{
	"cell_type": "code",
	"execution_count": null,
	"metadata": {
	"collapsed": true
	},
	"outputs": [],
	"source": []
	}
	],
	"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.4.3"
	}
	},
	"nbformat": 4,
	"nbformat_minor": 0
	}