kinverarity1/fix.ipynb

## fix.ipynb
{
 "metadata": {
  "name": "Untitled1"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": "%matplotlib inline\n\nimport numpy as np",
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": "def fix(x, y, value=0):\n    '''If a line segment in array crosses below *value*, introduce\n    an interpolated point where it crosses and replace the one\n    that's below with nan.'''\n    x = np.asarray(x, dtype=float)\n    y = np.asarray(y, dtype=float)\n    xn = []\n    yn = []\n    assert x.ndim == 1\n    assert y.ndim == 1\n    for i in range(1, len(y)):\n        indices = np.array([i - 1, i])\n        ls_x = x[indices]\n        ls_y = y[indices]\n        if ls_y[0] < value and ls_y[1] > value:\n            func = interp1d(ls_y, ls_x)\n            xn += [ls_x[0], func(value), ls_x[1]]\n            yn += [np.nan, value, ls_y[1]]\n        elif ls_y[0] > value and ls_y[1] < value:\n            func = interp1d(ls_y[::-1], ls_x[::-1])\n            xn += [ls_x[0], func(value), ls_x[1]]\n            yn += [ls_y[0], value, np.nan]\n        elif ls_y[0] < value and ls_y[1] < value:\n            xn += list(ls_x)\n            yn += [np.nan, np.nan]\n        else:\n            xn += list(ls_x)\n            yn += list(ls_y)\n    return np.array(xn), np.array(yn)\n\nx = np.array([1, 2, 3, 4, 5])\ny = np.array([2, 1,-1,-1, 5])\nx2p, y2p = fix(x, y, 0)\nx2n, y2n = fix(x, y * -1, 0)\nprint 'x2p', x2p\nprint 'y2p', y2p\nprint 'x2n', x2n\nprint 'y2n', y2n\nplt.plot(x, y, marker='o', ms=10, mfc='none', mec='k', ls='', label='original array')\nplt.plot(x2p, y2p, color='g', marker='+', ms=6, label='array > 0')\nplt.plot(x2n, y2n * -1, color='r', marker='x', ms=6, label='array < 0')\nplt.xlim(0, 6)\nplt.ylim(-2, 6)\nplt.legend(loc='best')",
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": "x2p [ 1.          2.          2.          2.5         3.          3.          4.\n  4.          4.16666667  5.        ]\ny2p [  2.   1.   1.   0.  nan  nan  nan  nan   0.   5.]\nx2n [ 1.          2.          2.          2.5         3.          3.          4.\n  4.          4.16666667  5.        ]\ny2n [ nan  nan  nan   0.   1.   1.   1.   1.   0.  nan]\n"
      },
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 19,
       "text": "<matplotlib.legend.Legend at 0x90bde90>"
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAW0AAAD8CAYAAAC8TPVwAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XtcVHX+P/DXIIIzmkpqPty+3slMBW+rzojKQKm1SBKL\nZrYpai7efqYS6UYhEorpmvptv2WloSLboywxUVO6cMRkoNW0bDUVFXPVvKQmMjhMzvn9wTIxcpsZ\n5syZM7yejwcPOXAu7zPIiw9vzvkclSiKIoiISBF85C6AiIjsx9AmIlIQhjYRkYIwtImIFIShTUSk\nIAxtIiIFcSq009LSMHToUAwaNAibNm1ydU1ERFQLh0NbEAQYDAbk5+dDEAScOXNGirqIiKgGvo5u\nkJOTg6CgIERFReHWrVtYuXKlFHUREVENHA7tq1ev4vz589i5cyfOnDmDJ598Ej/++KPNOiqVymUF\nEhE1FvbcoO5we6Rt27YYNWoUfH190aNHDzRr1gzXrl2r8eDe+LZ48WLZa+D58fx4ft73Zi+HQ3vY\nsGHYs2cPAODixYsoLS1FmzZtHN0NERE5weH2SEREBPLy8jB48GBYLBa89dZbbIcQEbmJw6ENAK+/\n/rqr61AMvV4vdwmS4vkpG8/P+6lER5op9u5UpXKoR0NE1NjZm5tOjbSddf/99+PGjRvuPCS5UEBA\nAK5fvy53GUSNmltH2hyBKxu/fkTSsff7i3OPEBEpCEObiEhB3NrTronRaERBQQEMBgPKysqgVquh\n0+mg1Wqh0Wjctg8iIiWQtaedlJQEf39/m4CtGsAmkwkpKSl1HssV+3DE4sWLERgYiOeee67WdbKz\ns/HFF19g7dq1Th0jOTkZv/zyC958801ny5QEe9pE0vHIq0eqMhqN8Pf3R2Jios3HNRoNwsPDER4e\njtTUVBiNxlpHy67Yh6OWLFlS7zqRkZGIjIx0+hi8WYmIaiNbT7ugoAA6na7OdXQ6HQoKCiTdR6V3\n330XQUFB6NevH0aPHo1Tp04BAGJjY/Hkk0+iT58+WLRoEWJjY7Fq1SoAwO7duxEcHIz+/ftjypQp\n6NixI86dO4eNGzdaQ1uv1+Pll19GaGgounbtikmTJll/mi5btgxDhgxB3759ERgYiO3bt9dZ4+XL\nlxEVFYWhQ4eiW7duCAsLw9WrVwEAXbp0wYQJE9CrVy9s37692vLOnTsREhKCQYMGoXPnzkhKSgIA\nTJ8+3eaHXmZmJqKjo+t9vYhIHrKFtsFggFarrXMdnU4Hg8Eg6T4A4KuvvsLKlSshCAKOHDmCiRMn\nIioqyvr5O3fu4IcffsDy5cuhUqmgUqnwyy+/YNKkScjMzMThw4cRFhaGCxcuWEfJVUfLZ86cwb59\n+3D06FF89dVXyMvLw08//YQvv/wSeXl5+O6775CammoN0tp+Rfrwww8REhKC/Px8nDlzBhqNBhkZ\nGdbjBQUF4dixY4iKiqq2/MYbb2Dz5s3417/+BYPBgLS0NFy/fh1z5szBxo0bYbFYAADvvPMOZs6c\nWefrRUTyka09UlZWVm/LQqPRoKysTNJ9AMCePXswYcIE68RXkydPxgsvvIDi4mKoVCoMGzbMZn1R\nFJGXl4devXohKCgIADBp0iTMnTvXZh2gIkwrR90tWrRAYGAgrl+/jtDQUGzatAkZGRk4ffo0CgoK\nUFpaWmedc+fOxf79+/HGG2/g1KlT+OGHH2x+aA0fPtxm/arL2dnZyM7ORmZmJo4fPw5RFFFaWoq+\nffuia9eu2LlzJx566CFcunQJI0eOrLMOIpKPbCNttVoNo9FY5zpGoxFqtVrSfQCocWpEURRhNpsB\nAM2bN6+2TdOmTatt4+NT88tZ9fiVf2z49ttvodPpcPv2bYwePRoLFy60jnZrs3DhQixevBjt27dH\nXFwcRo0aZVNDixYtbNavXC4tLUW/fv1w5MgRDBw4ECtXrrSpf/bs2Xj//feRnp6OuLi4OmsgInnJ\nFtr29JoNBkOdPWtX7AMARo8ejQ8//NA6L3h6ejratm2LwMDAWq+CCQkJwcmTJ3H06FEAwCeffIKb\nN2/W+EfEmn4g7N+/H4MGDcK8efMwfPhwZGVl4e7du3XWmZOTg3nz5uHZZ59Fu3bt8Pnnn9e7DQCc\nOnUKJSUleO211xAREQFBEGAymazbxsTE4PDhw9i2bRumTp1a7/6ISD6yhbZWq62311xfz9oV+wCA\nxx57DPPnz0d4eDj69OmDjIwM7Ny509q/rimIAwIC8MEHH2DSpEkYOHAgcnJy4OvrC41GU22be7dX\nqVR45plncO3aNfTp0wfh4eHo168fbt68idu3b9d6zKSkJLz44ovQarWYOXMmYmJiUFRUVOe5AUDf\nvn0xZswYPPLIIxg+fDh++OEH/PGPf7Ru27RpU8TExECn0+H++++vd39EJB/Zr9P28/ODTqeDTqez\nXmNtMBhgMBhQXl5u13XaDd2HM0pKSpCamork5GSo1Wp8++23iIyMxIULF1x+LKmVlpYiNDQUb7/9\nNgYNGlTrerxOm0g69n5/yT5hlJLviHz11VeRlZWFpk2bomnTpli9ejVCQkIkO54U9u7di4kTJ2La\ntGlYsWJFnesytImko5jQJuXg149IOpzlj4jIC8k+YRQRUV04IZwttkfIbvz6kbu5e0I4OUk6YdSA\nAQPQqlUrAEC3bt2wYcMGZ3ZDRFSrmiaEE4oF6LvoJZsQTgkcDu07d+4AAHJzc11eDBFRpZomhKsM\n7UqVN9iFh4e7uTr5OPyHyO+++w5GoxGjR4/Go48+isLCQinqIqJGrqYb436z/GazbM+EcN7G4ZF2\n8+bNkZCQgGnTpuHUqVN44okncPLkyWrzbiQnJ1vf1+v10Ov19e773p+iznDFPjzVrl278PLLL8Nk\nMiE4OBgbNmzAfffdJ3dZRJKonBBOKBYgFAuwiBYs3b8UprsmNG/aHPoueui76OudEM5TCYIAQRAc\n31B0kMlkEsvKyqzLgwcPFv/zn//YrFPbbus73OLcxY6WI8k+3GnNmjXi119/Xe96V65cER944AGx\nqKhIFEVRXLhwoThr1iypy7PhxH8XIqelpqaKpaWl1uW9RXvFB1c9aLNOaWmpmJqa6u7SJGHv95fD\n7ZH09HTEx8cDAC5evIhbt26hQ4cOjv+08CAWiwUvvPACtFotevfujV69eiE/Px9A9YcgTJkyxbr8\nt7/9DSdPnsTIkSMxdOhQdOnSBVFRUTCZTMjMzLS5O/Knn37Cgw8+iN9+s/31LjAwEIsXL0ZQUBDW\nrl2LGzdu1FhjTk4OBg8ejO7duwMAZs6ciczMTIleESL53Tsh3NZjW9G7XW+bdeyZEM7bONwemTZt\nGqZMmYIRI0YAqAjx2qYktUflrz4AsGTfEizZV//jvOxV+etTfb755hv8/PPP1v8gy5cvx/Lly7Fj\nxw4Avz8EAQCmTJlis/zSSy9hypQpmDhxIn777TcMHDgQu3fvxvjx4xEfH4/jx4/jkUcewfr16xEb\nGwtfX9uXPCIiAhERESguLsaGDRswaNAg6PV6rF+/3ma98+fP43/+53+syw8++CBu3bqF27dvV5uS\nlcgbaLVarF69GuHh4TDfNSPreBb+70//Z7OOwWDAggULZKpQHg6Htq+vr/VpKa5wb7Am65MbtL9k\nIdnhfWi1Wrz22mt4++23cebMGQiCgJYtWwJAjQ9BqLr8+uuvIycnBytXrsSJEydw8eJF3L59G02b\nNsXzzz+P9957D3//+9+xadMm7N+/v846fHx80KRJk2rBDtT+NJsmTZo4dK5ESqHRaGAymZCamoqm\nPZuia+uueLrP09UmhGtMl/sBvCMSQMUf+ObNm4cXX3wRUVFR6NmzJ7Zs2WL9/L0PQai6PGHCBNy9\nexdPP/00IiIicP78eWvAxsXFYfDgwQgNDUVQUBA6depU7difffYZ1qxZg4sXL2L69On45ptvrNfA\nV9WpUyebK3UuXLiAgICAeh/wQKRkKSkpMBqNGJcxDm2vtMUrr7xivSNywYIFjS6wAQ8LbVdc9eHM\nPr744gtERkYiLi4Od+7cQVpamvUBAbWNcCvl5ORg3759CA4OxrFjx1BYWIinn34aANCxY0fodDrM\nnz8fb775Zo3bHz9+HK+88kq1R4Xda+TIkYiPj0dRURECAwOxbt06m+dYEnmrpv5NUfhrIQ79v0Po\n3Lqz3OXIjqENYMaMGZg4cSL69++PgIAAjB07FqtWrYIoijU+kKDq8rJly/DUU0+hffv26NSpE/78\n5z/j9OnT1s/HxsZi7ty5+NOf/lTjse3txz3wwANIT09HTEwMysvLERgYiM2bNzt8rkRKk1uci+73\nd2dg/xfnHpGQxWLBnDlz0LVrVyQkJMhdToM1tq8feYbp2dPRs01PxA+Nl7sUSXFqVpmVlJSgbdu2\nOHv2LObMmSN3OUSKVHnVSEyvGLlL8Rge1R7xJvfddx+uX78udxlEisbWSHUcaRORx9p6bCvG9xov\ndxkehSNtIvJIla2RQ389JHcpHoUjbSLySGyN1IyhTUQeia2RmrE9QkQeh62R2sk/0t61C7h50/Zj\nN29WfNyd+1CAXbt2oW/fvujZsyfGjx+PkpISuUsikgRbI7WTP7RDQoDExN9D9+bNiuUq05q6ZR8e\n4vjx40hKSqr28atXr2Lq1KnYtm0bfvzxR3Tr1g2LFi2SoUIi6bE1Ugd3TuZd6+Fu3BDFWbNE8ezZ\nin9v3HD8oA3Yx927d8W5c+eKQ4YMEXv16iU+8sgj4oEDB0RRFMXJkyeLkZGRYu/evcWFCxeKsbGx\n1uVFixaJJ06cEB977DFRp9OJnTt3FseOHSveuXNH3LJlizh06FDrMc6dOyf+4Q9/EM1mc7Xjl5WV\niZs3bxaHDx8uDhgwQHz33XerrbNlyxZxzJgx1uXi4mKxVatWDrxADSfRfxciG+W/lYttXm8jFt8o\nlrsUt7L3+8szQlsUK8IWcM3b2bMO1WswGMTx48dbl9PS0sTIyEhRFCtCe+TIkdbPxcbG2iwnJCSI\nmZmZoiiKotlsFoODg8Vt27aJ5eXlYvv27cVjx46JoiiKr776qvjyyy9XO/brr78udu3aVZw9e7Z4\n+PDhWmtMS0sTZ8yYYV02m82iSqUSS0pKHDrXhmBokzvsLdorDn5vsNxluJ2931/yt0eAinbGypXA\n2bPArFnAjRuOR/WNGxXbnj1bsa97e9x1qDqfdkJCAj755BOUlpYCsG8+7TZt2mDlypWYMWNGjfNp\nWywWbNq0CXFxcdWOXTl/dpMmTep8mITI+bSpkWBrpG7yh3Zl/3npUqBLl4p/q/an3bCPXbt2ISIi\nAj4+PoiKisKMGTNgsVisn69vPu333nsPXbp0wYIFCzBgwACb+bQ/+OADZGdn1zqfdnx8PI4ePYpB\ngwZh9uzZGDJkCN5///1q63Xq1AmXLl2yLnM+bfJGnGukfvKH9oEDFSHbunXFcuvWFcsHDrhtH1Xn\n0x44cCCysrIcmk87KSkJ48aNAwAUFhZanwNZdT7tmTNn1roPf39//OUvf8H+/fvx/vvvo6ioqNo6\nI0eOREFBgfVznE+bvBGvGqmf/NdpR0RU/1jr1jV/XKJ9yDmf9r169+6NZcuWVfs459OmxoCtkfpx\nPm0JcT5tIvuZ75rRYVUHHPpr43xCjeTzaV+5cgUdO3bEyZMnnd2FV+N82kSOYWvEPk61R8xmM+Li\n4qr9gY5+x/m0iRzD1oh9nArthIQEzJw5E2lpabWuk5ycbH1fr9dDr9c7cygiagQa41wjgiBAEASH\nt3O4p71x40ZcuHABiYmJCAsLw7p16/Dwww/b7pQ9ba/Erx9JJed0Dl7NfRWFzxfKXYps7P3+cji0\nQ0NDrVdUHDlyBA8//DA+/fRTtG/fvt6D85te2fj1I6k0lof31sXe7y+H2yP79u2zvh8WFoZ33nnH\nJrDrEhAQUO3yOVKOgIAAuUsgL9QYWyMN4dbrtPmHOSK6F68acUyDQjs3N9dVdRBRI8WrRhwj/x2R\nRNRosTXiOPnnHiGiRoutEccxtIlINmyNOI7tESKSBVsjzuFIm4hkwdaIcxjaRCQLtkacw/YIEbkd\nWyPO40ibiNyOrRHnMbSJyO3YGnEe2yNE5FZsjTQMR9pE5FZsjTQMQ5uI3IqtkYZhe4SI3IatkYbj\nSJuI3IatkYZjaBOR27A10nBsjxCRW7A14hocaRORW7A14hoMbSJyC7ZGXIPtESKSHFsjruPwSPvu\n3buYOnUqhg0bhuHDh+Pf//63FHURkRdha8R1HA7tnTt3wsfHB19//TVSU1ORmJgoRV1E5EXYGnEd\nh9sjY8eOxZgxYwAAxcXFCAgIcHlRROQ92BpxLad62k2aNEFsbCyysrLw8ccf17hOcnKy9X29Xg+9\nXu/MoYhI4dgaqZkgCBAEweHtVKIois4e9PLlyxgyZAiOHz8OtVr9+05VKjRgt0TkRaZnT0fPNj0R\nPzRe7lI8mr256XBPOyMjA2lpaQAAtVoNHx8f+PjwykEiqq6yNRLTK0buUryGw+2RmJgYxMbGIjQ0\nFGazGWvXroW/v78UtRGRwrE14noNao/UulO2R4gIbI04wt7c5M01RCQJXjUiDTajiUgSbI1Ig6FN\nRJLgDTXSYHuEiFyOrRHpcKRNRC7H1oh0GNpE5HJsjUiH7REicim2RqTFkTYRuRRbI9JiaBORS7E1\nIi22R4jIZdgakR5H2kTkMmyNSI+hTUQuw9aI9NgeISKXYGvEPTjSJiKXYGvEPRjaROQSbI24B9sj\nRNRgbI24D0faRNRgbI24D0ObiBqMrRH3YXvEDkajEQUFBTAYDCgrK4NarYZOp4NWq4VGo5G7PCJZ\nsTXiXgzteiQlJcHf3x86nQ7z58+HRqOxhvjq1athMpmQkpIid5lEsmFrxL0cbo+YzWY899xzGDFi\nBIYMGYLs7Gwp6vIIRqMR/v7+SExMRHh4ODQaDYRiARqNBuHh4UhMTISfnx+MRqPcpRLJhq0R93I4\ntDMzM9GuXTvk5eVhz549mDNnjhR1eYSCggLodDqbjwnFgs2yTqdDQUGBG6si8hyVrZGYXjFyl9Jo\nOBza48aNs7YDLBYLfH29t8NiMBig1WrrXEen08FgMLipIiLPwtaI+zmcuM2bNwcAlJSUYNy4cVi6\ndGmN6yUnJ1vf1+v10Ov1ThUop7KyMmtLpHKEvWTfEpy+cRrdA7pD30UPfRc9ysrK5C2USCZsjThP\nEAQIguDwdk4Nk8+fP4/o6GjMnj0bEyZMqHGdqqGtVGq1Gkaj0RrOAHC7/Da+OPMF/n3l3xjQYQBK\nS0uhVqvlLZRIBrxqpGHuHcwuWbLEru0cbo9cvnwZo0aNwooVKxAbG+vo5opSU7+6hV8LfBv3LZJC\nk5CUm4R+b/eD2EOEKIoyVUkkD7ZG5OFwaC9btgy//vorUlJSEBYWhrCwMNy5c0eK2mSn1Wqr9av1\nXfTwUfkgqmcUvo37FoPLBuOjqx9h4LsDsePEDoY3NRpsjchDJUqQMiqVymvCKykpCX5+ftDpdNDp\ndNbrtA0GAwwGA8rLy5G8JBk7TuxAspAMH5UPkvXJiOwRCZVKJXf5RJIw3zWjw6oOOPTXQxxpu4i9\nucnQtoO9d0RaRAvDmxqFnNM5eDX3VRQ+Xyh3KV6DoS0jhjd5u+nZ09GzTU/ED42XuxSvwdD2AAxv\n8kZsjUiDoe1BGN7kTdgakQZD2wMxvMkbsDUiDYa2B2N4k1KxNSIdhrYCMLxJadgakQ5DW0EY3qQU\nbI1Ih6GtQAxv8mRsjUiLoa1gDG/yRGyNSIuh7QUY3uRJ2BqRFkPbizC8SW5sjUiPoe2FGN4kF7ZG\npMfQ9mIMb3I3tkakx9BuBBje5A5sjbgHQ7sRYXiTlNgacQ+GdiPE8CYpsDXiHgztRozhTa7C1oj7\nMLTJ7vC298k81PiwNeI+DG2yqiu8k5KS4O/vbxPSVUPcZDIhJSVF7lMgmbA14j725qbDT2OvqrCw\nEGFhYQ3ZBblB1afHJ4UmISk3CQPfHYit32+Fn58fEhMTER4eDo1GgzUFa6DRaBAeHo7ExET4+fnB\naDTKfQokA/NdM7KOZyGmV4zcpVAVTof2ihUrMH36dJhMJlfWQxK6N7z/lvM3bGq2CTtO7LD+hN/+\n43abbXQ6HQoKCuQol2SWW5yL7vd3Zy/bwzgd2oGBgdi2bRvbIApUGd6T70xG6shU68h7x4kd1dbV\n6XQwGAwyVEly23psK8b3Gi93GXQPX2c3jI6ORnFxca2fT05Otr6v1+uh1+udPRRJxHTHhKf7Po0L\nxgt4//D7mPjJRJSaS9FvXT+0btYaUT2jME87D2VlZXKXSm5W2Ro59NdDcpfitQRBgCAIDm/ndGjX\np2pok2dSq9UwGo1YoFuABboFsIgWBL8djJ9v/4yk0CREPxINo9EItVotd6nkZmyNSO/eweySJUvs\n2k6y0CbPV9mvDg8PB1DRNmmraYuMpzLwROYTAIBWF1tBp9PJWSbJgK0Rz9Wgq0cA8IYNBdNqtdX6\n1VE9o9C/Q3989uxnmLVrFjbkb4BWq5WpQpIDrxrxbA0K7S5duiA/P99VtZCbaTQamEwmpKam4ssv\nv4TRaMQ87TwYjUZcP3Yd0WXR2HF3B/ac2yN3qeRGbI14Nt5cQ3XeEXni1xN4IvMJvBXxFqIfiZa7\nVHID3lAjD94RSS5z+NJhBncjwblG5GNvbvIPkVSvyh535R8nGdzei60Rz8fQJrswuBsHXjXi+Rja\nZDcGt3fjDTXKwNAmhzC4vRdbI8rA0CaHMbi9E1sjysDQJqcwuL0LWyPKwdAmpzG4vQdbI8rB0KYG\nYXB7B7ZGlIOhTQ3G4FY2tkaUhaFNLsHgVi62RpSFoU0uw+BWJrZGlIWhTS7F4FYWtkaUh6FNLsfg\nVg62RpSHoU2SYHArA1sjysPQJskwuD0bWyPKxNAmSTG4PRdbI8rE0CbJMbg9E1sjysTQJrdgcHsW\ntkaUy+EH+1osFsyYMQNDhw5FWFgYTp8+LUVd5IUqg/uTv09DdmEGAEAoFio+efMmsGuXfMV5u127\nKl5jVLzmucW56NusMzof+EHmwshRDof29u3bUV5ejvz8fCxfvhzx8Xz4J9mvf4f+WJjwKa4tiEN2\nYUZFaN+8CSQmAiEhcpfnvUJCKl7jmzchFAvY+c0W/H1fM77mCuRwe+TAgQN4/PHHAQBDhgzBwYMH\nXV4UebfgniOADXtwaNrjuPZsOLA1EVi6FGjdWu7SvFfr1hWvcWIiWg5rgo5vfoQ2H/yLr7kCORza\nt27dQsuWLa3LTZo0gcVigY+P7aA9OTnZ+r5er4der3e6SPIuQrEA4WcB5qlP4R/P/xPv/HUgLh1Z\nA30XPfRd9HKX55WEYgFCsYDWw5tiwTNr8eK4Vjh/9hPoxV/4mstEEAQIguDwdirRnme2VxEfHw+t\nVotx48YBADp27Ijz58/b7tTOR8FTI/bflsg/2p7BjBVfwffd9cBzz8ldlXf772v+zdE9CLzth/u/\nMnCk7UHszU2He9ohISHYvXs3AKCgoADBwcGOV0eNW2UPe+lSXAsbAt+9nwNxcUBGhtyVea/K1zw5\nGX2+u4T7N31k7XGTsjg80hZFEbNmzcL3338PAEhPT0ePHj1sd8qRNtVl166KP4C1bg2hWKj49Twv\nDxg7FtiwAYjm5YAuV/maf/MNbr00Dy2PHKsI7AMHgIgIuasj2J+bDoe2Kw9OZOPwYeCJJ4C33mJw\nS2X6dKBnT4BXfXkce3OTN9eQ5+jfH/jss4rgBhjcrmY2A1lZwCHeUKNkDG3yLAxu6eTmAt27A505\n14iSMbTJ8zC4pbF1KzCec40oHUObPBOD27XYGvEaDG3yXAxu12FrxGswtMmzMbhdg60Rr8HQJs/H\n4G4Ytka8CkOblIHB7Ty2RrwKQ5uUg8HtHLZGvApDm5SFwe0Ytka8DkOblIfBbT+2RrwOQ5uUicFt\nH7ZGvA5Dm5SLwV03tka8EkOblI3BXTu2RrwSQ5uUj8FdM7ZGvBJDm7wDg9sWWyNei6FN3oPB/Tu2\nRrwWQ5u8C4O7AlsjXouhTd7nnuA2Pv44CgoKYDAYUFZWBrVaDZ1OB61WC41GI3OxDWc0Gm3Or7mf\nH+I/+gh3Cwuhlrs4cjk+I5K81+HDKBk+HHuffBL3P/+8NaSrhpzJZEJKSorclTotKSkJ/v7+Nj+E\n7mRnw/TSS/jHX/6i+PNrTOzNTZ+GHCQrKwvPPvtsQ3ZBJBnjww/jw8mTEbNjB8JPnrSOqjUaDcK7\nd0ei2Qw/Pz8YjUaZK3WO0WiEv78/Es1mhHfvbj2/Zjt2oFV0tOLPj2rm9Ej7hRdeQE5ODvr3749/\n/vOftjvlSJs8wFdffQUACC8rA8aOBf7xD2DGDODcOWDMGGDnTnxZVASVSoXw8HCZq3Wc9fy6d7ee\nD/7wB6BdO6B9eyAnR9Hn19hIPtIOCQnB22+/zXAmj2UwGKDVaoGICODTT4E5c4ClS38PuM6dodPp\nYDAY5C7VKdbz69y54nzGjAH+93+B8nIgJ0fx50c1q/cPkRs2bMCaNWtsPrZx40aMHz8egiDUul1y\ncrL1fb1eD71e72yNRE4pKyv7/Q+NERHAu+8C06YBjz9eEd4ANABGHTxYMfpWmFEHD0JTte6uXYEX\nXwRmzbJe6qfRaFBWViZThVQXQRDqzNDa1Bva06ZNw7Rp0xzecdXQJpKDWq2G0WisCO5z54DVq4HM\nTCAhARgxAmjTBqbycly9cgX44x/lLtdhV69cgSk4GP5+fsAvvwC7dgELFwLZ2RXn27kzjEYj1Gpe\nQ+KJ7h3MLlmyxK7teMkfeS2dToeCggLbnm/nzkBIiHX566IiNOvVC1Bgz7dZYCAOoEpPOz+/4vxm\nzrSen6GoCDqdTu5SyYUadPWISqWCSqVyVS1ELqXVaiv6uenpvwc28HsPOD39976wAnn7+VHNGjTS\nDg0NRWhoqKtqIXIpjUYDk8mEVD8/6IqKoGvXznqdtqGoCAZfX5SXlyv2BhtvPz+qGW+uIa937x2D\n3n5HpLehXPakAAAE/0lEQVSdX2Nhb24ytImIPIBb7ogkIiL3YmgTESkIQ5uISEEY2kRECsLQJiJS\nEIY2EZGCMLSJiBSEoU1EpCAMbSIiBWFoExEpCEObiEhBGNpERArC0CYiUhCGNhGRgjC0iYgUhKFN\nRKQgDG0iIgVhaDtIEAS5S5AUz0/ZeH7ez+HQ/vXXXxEZGQm9Xo+hQ4eioKBAiro8lrf/p+H5KRvP\nz/s5HNqrV6/GyJEjIQgCNm7ciNmzZ0tRFxER1cDX0Q3mz58Pf39/AIDZbIZarXZ5UUREVAuxDuvX\nrxf79Olj83bw4EFRFEXx0qVLYv/+/cW8vLxq2wHgG9/4xje+OfhmD5VozzPb73H06FE888wzWLVq\nFUaPHu3o5kRE5CSHQ/vYsWOIjo7G1q1bERQUJFVdRERUA4dDOyoqCt9//z06d+4MAGjdujWysrIk\nKY6IiGw51R6pjcViwaxZs/D999/D398f69evR/fu3V21e49RWFiIRYsWITc3V+5SXMZsNmPq1Kk4\nd+4cTCYTXnnlFURGRspdlsvcvXsX06dPx8mTJ6FSqbBu3Tr07t1b7rJc7sqVKxg4cCC+/PJL9OjR\nQ+5yXGrAgAFo1aoVAKBbt27YsGGDzBW5VlpaGrKzs2E2mzFnzhxMnjy5xvUcvnqkLtu3b0d5eTny\n8/NRWFiI+Ph4bN++3ZWHkN2KFSuwZcsWtGjRQu5SXCozMxPt2rVDRkYGbty4gX79+nlVaO/cuRM+\nPj74+uuvsW/fPiQmJnrd/02z2Yy4uDg0b95c7lJc7s6dOwDgVQOlqgRBgMFgQH5+PkpLS7FixYpa\n13XpHZEHDhzA448/DgAYMmQIDh486Mrde4TAwEBs27YNLvwFxSOMGzcOKSkpACp+Y/L1denPc9mN\nHTsW77zzDgCguLgYAQEBMlfkegkJCZg5cyY6dOggdyku991338FoNGL06NF49NFHUVhYKHdJLpWT\nk4OgoCBERUUhMjISTz75ZK3rujS0b926hZYtW1qXmzRpAovF4spDyC46OtrrAg0AmjdvjhYtWqCk\npATjxo3D0qVL5S7J5Zo0aYLY2FjMnTsXEydOlLscl9q4cSPatWuHUaNGAYDXDSqaN2+OhIQE7N27\nF+vWrcOzzz7rVdly9epVHDp0CB9//LH1/Grj0tBu2bIlSkpKrMsWiwU+PpzeRCnOnz+P8PBwTJo0\nCRMmTJC7HEls3LgRJ0+exPTp01FWViZ3OS6Tnp6Ozz//HGFhYThy5AgmT56My5cvy12Wy/To0cMa\nZA899BDatGmDS5cuyVyV67Rt2xajRo2Cr68vevTogWbNmuHatWs1ruvSRA0JCcHu3bsBAAUFBQgO\nDnbl7klCly9fxqhRo7BixQrExsbKXY7LZWRkIC0tDQCgVqvh4+PjVQOKffv2QRAE5Obmol+/fti8\neTPat28vd1kuk56ejvj4eADAxYsXcevWLa9qAw0bNgx79uwBUHF+paWlaNOmTY3ruvT3/Keeegqf\nf/45QkJCAFS80N5KpVLJXYJLLVu2DL/++itSUlKsve3PPvsMzZo1k7ky14iJiUFsbCxCQ0NhNpux\ndu1a63QM5PmmTZuGKVOmYMSIEQAqssWbfuhGREQgLy8PgwcPhsViwVtvvVVrxrj0kj8iIpKW9/yo\nIiJqBBjaREQKwtAmIlIQhjYRkYIwtImIFIShTUSkIAxtIiIF+f+NrMkYNBuG2gAAAABJRU5ErkJg\ngg==\n",
       "text": "<matplotlib.figure.Figure at 0x8a062d0>"
      }
     ],
     "prompt_number": 19
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": "",
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}
	{
	"metadata": {
	"name": "Untitled1"
	},
	"nbformat": 3,
	"nbformat_minor": 0,
	"worksheets": [
	{
	"cells": [
	{
	"cell_type": "code",
	"collapsed": false,
	"input": "%matplotlib inline\n\nimport numpy as np",
	"language": "python",
	"metadata": {},
	"outputs": [],
	"prompt_number": 2
	},
	{
	"cell_type": "code",
	"collapsed": false,
	"input": "def fix(x, y, value=0):\n '''If a line segment in array crosses below value, introduce\n an interpolated point where it crosses and replace the one\n that's below with nan.'''\n x = np.asarray(x, dtype=float)\n y = np.asarray(y, dtype=float)\n xn = []\n yn = []\n assert x.ndim == 1\n assert y.ndim == 1\n for i in range(1, len(y)):\n indices = np.array([i - 1, i])\n ls_x = x[indices]\n ls_y = y[indices]\n if ls_y[0] < value and ls_y[1] > value:\n func = interp1d(ls_y, ls_x)\n xn += [ls_x[0], func(value), ls_x[1]]\n yn += [np.nan, value, ls_y[1]]\n elif ls_y[0] > value and ls_y[1] < value:\n func = interp1d(ls_y[::-1], ls_x[::-1])\n xn += [ls_x[0], func(value), ls_x[1]]\n yn += [ls_y[0], value, np.nan]\n elif ls_y[0] < value and ls_y[1] < value:\n xn += list(ls_x)\n yn += [np.nan, np.nan]\n else:\n xn += list(ls_x)\n yn += list(ls_y)\n return np.array(xn), np.array(yn)\n\nx = np.array([1, 2, 3, 4, 5])\ny = np.array([2, 1,-1,-1, 5])\nx2p, y2p = fix(x, y, 0)\nx2n, y2n = fix(x, y * -1, 0)\nprint 'x2p', x2p\nprint 'y2p', y2p\nprint 'x2n', x2n\nprint 'y2n', y2n\nplt.plot(x, y, marker='o', ms=10, mfc='none', mec='k', ls='', label='original array')\nplt.plot(x2p, y2p, color='g', marker='+', ms=6, label='array > 0')\nplt.plot(x2n, y2n * -1, color='r', marker='x', ms=6, label='array < 0')\nplt.xlim(0, 6)\nplt.ylim(-2, 6)\nplt.legend(loc='best')",
	"language": "python",
	"metadata": {},
	"outputs": [
	{
	"output_type": "stream",
	"stream": "stdout",
	"text": "x2p [ 1. 2. 2. 2.5 3. 3. 4.\n 4. 4.16666667 5. ]\ny2p [ 2. 1. 1. 0. nan nan nan nan 0. 5.]\nx2n [ 1. 2. 2. 2.5 3. 3. 4.\n 4. 4.16666667 5. ]\ny2n [ nan nan nan 0. 1. 1. 1. 1. 0. nan]\n"
	},
	{
	"metadata": {},
	"output_type": "pyout",
	"prompt_number": 19,
	"text": "<matplotlib.legend.Legend at 0x90bde90>"
	},
	{
	"metadata": {},
	"output_type": "display_data",
	"png": "iVBORw0KGgoAAAANSUhEUgAAAW0AAAD8CAYAAAC8TPVwAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XtcVHX+P/DXIIIzmkpqPty+3slMBW+rzojKQKm1SBKL\nZrYpai7efqYS6UYhEorpmvptv2WloSLboywxUVO6cMRkoNW0bDUVFXPVvKQmMjhMzvn9wTIxcpsZ\n5syZM7yejwcPOXAu7zPIiw9vzvkclSiKIoiISBF85C6AiIjsx9AmIlIQhjYRkYIwtImIFIShTUSk\nIAxtIiIFcSq009LSMHToUAwaNAibNm1ydU1ERFQLh0NbEAQYDAbk5+dDEAScOXNGirqIiKgGvo5u\nkJOTg6CgIERFReHWrVtYuXKlFHUREVENHA7tq1ev4vz589i5cyfOnDmDJ598Ej/++KPNOiqVymUF\nEhE1FvbcoO5we6Rt27YYNWoUfH190aNHDzRr1gzXrl2r8eDe+LZ48WLZa+D58fx4ft73Zi+HQ3vY\nsGHYs2cPAODixYsoLS1FmzZtHN0NERE5weH2SEREBPLy8jB48GBYLBa89dZbbIcQEbmJw6ENAK+/\n/rqr61AMvV4vdwmS4vkpG8/P+6lER5op9u5UpXKoR0NE1NjZm5tOjbSddf/99+PGjRvuPCS5UEBA\nAK5fvy53GUSNmltH2hyBKxu/fkTSsff7i3OPEBEpCEObiEhB3NrTronRaERBQQEMBgPKysqgVquh\n0+mg1Wqh0Wjctg8iIiWQtaedlJQEf39/m4CtGsAmkwkpKSl1HssV+3DE4sWLERgYiOeee67WdbKz\ns/HFF19g7dq1Th0jOTkZv/zyC958801ny5QEe9pE0vHIq0eqMhqN8Pf3R2Jios3HNRoNwsPDER4e\njtTUVBiNxlpHy67Yh6OWLFlS7zqRkZGIjIx0+hi8WYmIaiNbT7ugoAA6na7OdXQ6HQoKCiTdR6V3\n330XQUFB6NevH0aPHo1Tp04BAGJjY/Hkk0+iT58+WLRoEWJjY7Fq1SoAwO7duxEcHIz+/ftjypQp\n6NixI86dO4eNGzdaQ1uv1+Pll19GaGgounbtikmTJll/mi5btgxDhgxB3759ERgYiO3bt9dZ4+XL\nlxEVFYWhQ4eiW7duCAsLw9WrVwEAXbp0wYQJE9CrVy9s37692vLOnTsREhKCQYMGoXPnzkhKSgIA\nTJ8+3eaHXmZmJqKjo+t9vYhIHrKFtsFggFarrXMdnU4Hg8Eg6T4A4KuvvsLKlSshCAKOHDmCiRMn\nIioqyvr5O3fu4IcffsDy5cuhUqmgUqnwyy+/YNKkScjMzMThw4cRFhaGCxcuWEfJVUfLZ86cwb59\n+3D06FF89dVXyMvLw08//YQvv/wSeXl5+O6775CammoN0tp+Rfrwww8REhKC/Px8nDlzBhqNBhkZ\nGdbjBQUF4dixY4iKiqq2/MYbb2Dz5s3417/+BYPBgLS0NFy/fh1z5szBxo0bYbFYAADvvPMOZs6c\nWefrRUTyka09UlZWVm/LQqPRoKysTNJ9AMCePXswYcIE68RXkydPxgsvvIDi4mKoVCoMGzbMZn1R\nFJGXl4devXohKCgIADBp0iTMnTvXZh2gIkwrR90tWrRAYGAgrl+/jtDQUGzatAkZGRk4ffo0CgoK\nUFpaWmedc+fOxf79+/HGG2/g1KlT+OGHH2x+aA0fPtxm/arL2dnZyM7ORmZmJo4fPw5RFFFaWoq+\nffuia9eu2LlzJx566CFcunQJI0eOrLMOIpKPbCNttVoNo9FY5zpGoxFqtVrSfQCocWpEURRhNpsB\nAM2bN6+2TdOmTatt4+NT88tZ9fiVf2z49ttvodPpcPv2bYwePRoLFy60jnZrs3DhQixevBjt27dH\nXFwcRo0aZVNDixYtbNavXC4tLUW/fv1w5MgRDBw4ECtXrrSpf/bs2Xj//feRnp6OuLi4OmsgInnJ\nFtr29JoNBkOdPWtX7AMARo8ejQ8//NA6L3h6ejratm2LwMDAWq+CCQkJwcmTJ3H06FEAwCeffIKb\nN2/W+EfEmn4g7N+/H4MGDcK8efMwfPhwZGVl4e7du3XWmZOTg3nz5uHZZ59Fu3bt8Pnnn9e7DQCc\nOnUKJSUleO211xAREQFBEGAymazbxsTE4PDhw9i2bRumTp1a7/6ISD6yhbZWq62311xfz9oV+wCA\nxx57DPPnz0d4eDj69OmDjIwM7Ny509q/rimIAwIC8MEHH2DSpEkYOHAgcnJy4OvrC41GU22be7dX\nqVR45plncO3aNfTp0wfh4eHo168fbt68idu3b9d6zKSkJLz44ovQarWYOXMmYmJiUFRUVOe5AUDf\nvn0xZswYPPLIIxg+fDh++OEH/PGPf7Ru27RpU8TExECn0+H++++vd39EJB/Zr9P28/ODTqeDTqez\nXmNtMBhgMBhQXl5u13XaDd2HM0pKSpCamork5GSo1Wp8++23iIyMxIULF1x+LKmVlpYiNDQUb7/9\nNgYNGlTrerxOm0g69n5/yT5hlJLviHz11VeRlZWFpk2bomnTpli9ejVCQkIkO54U9u7di4kTJ2La\ntGlYsWJFnesytImko5jQJuXg149IOpzlj4jIC8k+YRQRUV04IZwttkfIbvz6kbu5e0I4OUk6YdSA\nAQPQqlUrAEC3bt2wYcMGZ3ZDRFSrmiaEE4oF6LvoJZsQTgkcDu07d+4AAHJzc11eDBFRpZomhKsM\n7UqVN9iFh4e7uTr5OPyHyO+++w5GoxGjR4/Go48+isLCQinqIqJGrqYb436z/GazbM+EcN7G4ZF2\n8+bNkZCQgGnTpuHUqVN44okncPLkyWrzbiQnJ1vf1+v10Ov19e773p+iznDFPjzVrl278PLLL8Nk\nMiE4OBgbNmzAfffdJ3dZRJKonBBOKBYgFAuwiBYs3b8UprsmNG/aHPoueui76OudEM5TCYIAQRAc\n31B0kMlkEsvKyqzLgwcPFv/zn//YrFPbbus73OLcxY6WI8k+3GnNmjXi119/Xe96V65cER944AGx\nqKhIFEVRXLhwoThr1iypy7PhxH8XIqelpqaKpaWl1uW9RXvFB1c9aLNOaWmpmJqa6u7SJGHv95fD\n7ZH09HTEx8cDAC5evIhbt26hQ4cOjv+08CAWiwUvvPACtFotevfujV69eiE/Px9A9YcgTJkyxbr8\nt7/9DSdPnsTIkSMxdOhQdOnSBVFRUTCZTMjMzLS5O/Knn37Cgw8+iN9+s/31LjAwEIsXL0ZQUBDW\nrl2LGzdu1FhjTk4OBg8ejO7duwMAZs6ciczMTIleESL53Tsh3NZjW9G7XW+bdeyZEM7bONwemTZt\nGqZMmYIRI0YAqAjx2qYktUflrz4AsGTfEizZV//jvOxV+etTfb755hv8/PPP1v8gy5cvx/Lly7Fj\nxw4Avz8EAQCmTJlis/zSSy9hypQpmDhxIn777TcMHDgQu3fvxvjx4xEfH4/jx4/jkUcewfr16xEb\nGwtfX9uXPCIiAhERESguLsaGDRswaNAg6PV6rF+/3ma98+fP43/+53+syw8++CBu3bqF27dvV5uS\nlcgbaLVarF69GuHh4TDfNSPreBb+70//Z7OOwWDAggULZKpQHg6Htq+vr/VpKa5wb7Am65MbtL9k\nIdnhfWi1Wrz22mt4++23cebMGQiCgJYtWwJAjQ9BqLr8+uuvIycnBytXrsSJEydw8eJF3L59G02b\nNsXzzz+P9957D3//+9+xadMm7N+/v846fHx80KRJk2rBDtT+NJsmTZo4dK5ESqHRaGAymZCamoqm\nPZuia+uueLrP09UmhGtMl/sBvCMSQMUf+ObNm4cXX3wRUVFR6NmzJ7Zs2WL9/L0PQai6PGHCBNy9\nexdPP/00IiIicP78eWvAxsXFYfDgwQgNDUVQUBA6depU7difffYZ1qxZg4sXL2L69On45ptvrNfA\nV9WpUyebK3UuXLiAgICAeh/wQKRkKSkpMBqNGJcxDm2vtMUrr7xivSNywYIFjS6wAQ8LbVdc9eHM\nPr744gtERkYiLi4Od+7cQVpamvUBAbWNcCvl5ORg3759CA4OxrFjx1BYWIinn34aANCxY0fodDrM\nnz8fb775Zo3bHz9+HK+88kq1R4Xda+TIkYiPj0dRURECAwOxbt06m+dYEnmrpv5NUfhrIQ79v0Po\n3Lqz3OXIjqENYMaMGZg4cSL69++PgIAAjB07FqtWrYIoijU+kKDq8rJly/DUU0+hffv26NSpE/78\n5z/j9OnT1s/HxsZi7ty5+NOf/lTjse3txz3wwANIT09HTEwMysvLERgYiM2bNzt8rkRKk1uci+73\nd2dg/xfnHpGQxWLBnDlz0LVrVyQkJMhdToM1tq8feYbp2dPRs01PxA+Nl7sUSXFqVpmVlJSgbdu2\nOHv2LObMmSN3OUSKVHnVSEyvGLlL8Rge1R7xJvfddx+uX78udxlEisbWSHUcaRORx9p6bCvG9xov\ndxkehSNtIvJIla2RQ389JHcpHoUjbSLySGyN1IyhTUQeia2RmrE9QkQeh62R2sk/0t61C7h50/Zj\nN29WfNyd+1CAXbt2oW/fvujZsyfGjx+PkpISuUsikgRbI7WTP7RDQoDExN9D9+bNiuUq05q6ZR8e\n4vjx40hKSqr28atXr2Lq1KnYtm0bfvzxR3Tr1g2LFi2SoUIi6bE1Ugd3TuZd6+Fu3BDFWbNE8ezZ\nin9v3HD8oA3Yx927d8W5c+eKQ4YMEXv16iU+8sgj4oEDB0RRFMXJkyeLkZGRYu/evcWFCxeKsbGx\n1uVFixaJJ06cEB977DFRp9OJnTt3FseOHSveuXNH3LJlizh06FDrMc6dOyf+4Q9/EM1mc7Xjl5WV\niZs3bxaHDx8uDhgwQHz33XerrbNlyxZxzJgx1uXi4mKxVatWDrxADSfRfxciG+W/lYttXm8jFt8o\nlrsUt7L3+8szQlsUK8IWcM3b2bMO1WswGMTx48dbl9PS0sTIyEhRFCtCe+TIkdbPxcbG2iwnJCSI\nmZmZoiiKotlsFoODg8Vt27aJ5eXlYvv27cVjx46JoiiKr776qvjyyy9XO/brr78udu3aVZw9e7Z4\n+PDhWmtMS0sTZ8yYYV02m82iSqUSS0pKHDrXhmBokzvsLdorDn5vsNxluJ2931/yt0eAinbGypXA\n2bPArFnAjRuOR/WNGxXbnj1bsa97e9x1qDqfdkJCAj755BOUlpYCsG8+7TZt2mDlypWYMWNGjfNp\nWywWbNq0CXFxcdWOXTl/dpMmTep8mITI+bSpkWBrpG7yh3Zl/3npUqBLl4p/q/an3bCPXbt2ISIi\nAj4+PoiKisKMGTNgsVisn69vPu333nsPXbp0wYIFCzBgwACb+bQ/+OADZGdn1zqfdnx8PI4ePYpB\ngwZh9uzZGDJkCN5///1q63Xq1AmXLl2yLnM+bfJGnGukfvKH9oEDFSHbunXFcuvWFcsHDrhtH1Xn\n0x44cCCysrIcmk87KSkJ48aNAwAUFhZanwNZdT7tmTNn1roPf39//OUvf8H+/fvx/vvvo6ioqNo6\nI0eOREFBgfVznE+bvBGvGqmf/NdpR0RU/1jr1jV/XKJ9yDmf9r169+6NZcuWVfs459OmxoCtkfpx\nPm0JcT5tIvuZ75rRYVUHHPpr43xCjeTzaV+5cgUdO3bEyZMnnd2FV+N82kSOYWvEPk61R8xmM+Li\n4qr9gY5+x/m0iRzD1oh9nArthIQEzJw5E2lpabWuk5ycbH1fr9dDr9c7cygiagQa41wjgiBAEASH\nt3O4p71x40ZcuHABiYmJCAsLw7p16/Dwww/b7pQ9ba/Erx9JJed0Dl7NfRWFzxfKXYps7P3+cji0\nQ0NDrVdUHDlyBA8//DA+/fRTtG/fvt6D85te2fj1I6k0lof31sXe7y+H2yP79u2zvh8WFoZ33nnH\nJrDrEhAQUO3yOVKOgIAAuUsgL9QYWyMN4dbrtPmHOSK6F68acUyDQjs3N9dVdRBRI8WrRhwj/x2R\nRNRosTXiOPnnHiGiRoutEccxtIlINmyNOI7tESKSBVsjzuFIm4hkwdaIcxjaRCQLtkacw/YIEbkd\nWyPO40ibiNyOrRHnMbSJyO3YGnEe2yNE5FZsjTQMR9pE5FZsjTQMQ5uI3IqtkYZhe4SI3IatkYbj\nSJuI3IatkYZjaBOR27A10nBsjxCRW7A14hocaRORW7A14hoMbSJyC7ZGXIPtESKSHFsjruPwSPvu\n3buYOnUqhg0bhuHDh+Pf//63FHURkRdha8R1HA7tnTt3wsfHB19//TVSU1ORmJgoRV1E5EXYGnEd\nh9sjY8eOxZgxYwAAxcXFCAgIcHlRROQ92BpxLad62k2aNEFsbCyysrLw8ccf17hOcnKy9X29Xg+9\nXu/MoYhI4dgaqZkgCBAEweHtVKIois4e9PLlyxgyZAiOHz8OtVr9+05VKjRgt0TkRaZnT0fPNj0R\nPzRe7lI8mr256XBPOyMjA2lpaQAAtVoNHx8f+PjwykEiqq6yNRLTK0buUryGw+2RmJgYxMbGIjQ0\nFGazGWvXroW/v78UtRGRwrE14noNao/UulO2R4gIbI04wt7c5M01RCQJXjUiDTajiUgSbI1Ig6FN\nRJLgDTXSYHuEiFyOrRHpcKRNRC7H1oh0GNpE5HJsjUiH7REicim2RqTFkTYRuRRbI9JiaBORS7E1\nIi22R4jIZdgakR5H2kTkMmyNSI+hTUQuw9aI9NgeISKXYGvEPTjSJiKXYGvEPRjaROQSbI24B9sj\nRNRgbI24D0faRNRgbI24D0ObiBqMrRH3YXvEDkajEQUFBTAYDCgrK4NarYZOp4NWq4VGo5G7PCJZ\nsTXiXgzteiQlJcHf3x86nQ7z58+HRqOxhvjq1athMpmQkpIid5lEsmFrxL0cbo+YzWY899xzGDFi\nBIYMGYLs7Gwp6vIIRqMR/v7+SExMRHh4ODQaDYRiARqNBuHh4UhMTISfnx+MRqPcpRLJhq0R93I4\ntDMzM9GuXTvk5eVhz549mDNnjhR1eYSCggLodDqbjwnFgs2yTqdDQUGBG6si8hyVrZGYXjFyl9Jo\nOBza48aNs7YDLBYLfH29t8NiMBig1WrrXEen08FgMLipIiLPwtaI+zmcuM2bNwcAlJSUYNy4cVi6\ndGmN6yUnJ1vf1+v10Ov1ThUop7KyMmtLpHKEvWTfEpy+cRrdA7pD30UPfRc9ysrK5C2USCZsjThP\nEAQIguDwdk4Nk8+fP4/o6GjMnj0bEyZMqHGdqqGtVGq1Gkaj0RrOAHC7/Da+OPMF/n3l3xjQYQBK\nS0uhVqvlLZRIBrxqpGHuHcwuWbLEru0cbo9cvnwZo0aNwooVKxAbG+vo5opSU7+6hV8LfBv3LZJC\nk5CUm4R+b/eD2EOEKIoyVUkkD7ZG5OFwaC9btgy//vorUlJSEBYWhrCwMNy5c0eK2mSn1Wqr9av1\nXfTwUfkgqmcUvo37FoPLBuOjqx9h4LsDsePEDoY3NRpsjchDJUqQMiqVymvCKykpCX5+ftDpdNDp\ndNbrtA0GAwwGA8rLy5G8JBk7TuxAspAMH5UPkvXJiOwRCZVKJXf5RJIw3zWjw6oOOPTXQxxpu4i9\nucnQtoO9d0RaRAvDmxqFnNM5eDX3VRQ+Xyh3KV6DoS0jhjd5u+nZ09GzTU/ED42XuxSvwdD2AAxv\n8kZsjUiDoe1BGN7kTdgakQZD2wMxvMkbsDUiDYa2B2N4k1KxNSIdhrYCMLxJadgakQ5DW0EY3qQU\nbI1Ih6GtQAxv8mRsjUiLoa1gDG/yRGyNSIuh7QUY3uRJ2BqRFkPbizC8SW5sjUiPoe2FGN4kF7ZG\npMfQ9mIMb3I3tkakx9BuBBje5A5sjbgHQ7sRYXiTlNgacQ+GdiPE8CYpsDXiHgztRozhTa7C1oj7\nMLTJ7vC298k81PiwNeI+DG2yqiu8k5KS4O/vbxPSVUPcZDIhJSVF7lMgmbA14j725qbDT2OvqrCw\nEGFhYQ3ZBblB1afHJ4UmISk3CQPfHYit32+Fn58fEhMTER4eDo1GgzUFa6DRaBAeHo7ExET4+fnB\naDTKfQokA/NdM7KOZyGmV4zcpVAVTof2ihUrMH36dJhMJlfWQxK6N7z/lvM3bGq2CTtO7LD+hN/+\n43abbXQ6HQoKCuQol2SWW5yL7vd3Zy/bwzgd2oGBgdi2bRvbIApUGd6T70xG6shU68h7x4kd1dbV\n6XQwGAwyVEly23psK8b3Gi93GXQPX2c3jI6ORnFxca2fT05Otr6v1+uh1+udPRRJxHTHhKf7Po0L\nxgt4//D7mPjJRJSaS9FvXT+0btYaUT2jME87D2VlZXKXSm5W2Ro59NdDcpfitQRBgCAIDm/ndGjX\np2pok2dSq9UwGo1YoFuABboFsIgWBL8djJ9v/4yk0CREPxINo9EItVotd6nkZmyNSO/eweySJUvs\n2k6y0CbPV9mvDg8PB1DRNmmraYuMpzLwROYTAIBWF1tBp9PJWSbJgK0Rz9Wgq0cA8IYNBdNqtdX6\n1VE9o9C/Q3989uxnmLVrFjbkb4BWq5WpQpIDrxrxbA0K7S5duiA/P99VtZCbaTQamEwmpKam4ssv\nv4TRaMQ87TwYjUZcP3Yd0WXR2HF3B/ac2yN3qeRGbI14Nt5cQ3XeEXni1xN4IvMJvBXxFqIfiZa7\nVHID3lAjD94RSS5z+NJhBncjwblG5GNvbvIPkVSvyh535R8nGdzei60Rz8fQJrswuBsHXjXi+Rja\nZDcGt3fjDTXKwNAmhzC4vRdbI8rA0CaHMbi9E1sjysDQJqcwuL0LWyPKwdAmpzG4vQdbI8rB0KYG\nYXB7B7ZGlIOhTQ3G4FY2tkaUhaFNLsHgVi62RpSFoU0uw+BWJrZGlIWhTS7F4FYWtkaUh6FNLsfg\nVg62RpSHoU2SYHArA1sjysPQJskwuD0bWyPKxNAmSTG4PRdbI8rE0CbJMbg9E1sjysTQJrdgcHsW\ntkaUy+EH+1osFsyYMQNDhw5FWFgYTp8+LUVd5IUqg/uTv09DdmEGAEAoFio+efMmsGuXfMV5u127\nKl5jVLzmucW56NusMzof+EHmwshRDof29u3bUV5ejvz8fCxfvhzx8Xz4J9mvf4f+WJjwKa4tiEN2\nYUZFaN+8CSQmAiEhcpfnvUJCKl7jmzchFAvY+c0W/H1fM77mCuRwe+TAgQN4/PHHAQBDhgzBwYMH\nXV4UebfgniOADXtwaNrjuPZsOLA1EVi6FGjdWu7SvFfr1hWvcWIiWg5rgo5vfoQ2H/yLr7kCORza\nt27dQsuWLa3LTZo0gcVigY+P7aA9OTnZ+r5er4der3e6SPIuQrEA4WcB5qlP4R/P/xPv/HUgLh1Z\nA30XPfRd9HKX55WEYgFCsYDWw5tiwTNr8eK4Vjh/9hPoxV/4mstEEAQIguDwdirRnme2VxEfHw+t\nVotx48YBADp27Ijz58/b7tTOR8FTI/bflsg/2p7BjBVfwffd9cBzz8ldlXf772v+zdE9CLzth/u/\nMnCk7UHszU2He9ohISHYvXs3AKCgoADBwcGOV0eNW2UPe+lSXAsbAt+9nwNxcUBGhtyVea/K1zw5\nGX2+u4T7N31k7XGTsjg80hZFEbNmzcL3338PAEhPT0ePHj1sd8qRNtVl166KP4C1bg2hWKj49Twv\nDxg7FtiwAYjm5YAuV/maf/MNbr00Dy2PHKsI7AMHgIgIuasj2J+bDoe2Kw9OZOPwYeCJJ4C33mJw\nS2X6dKBnT4BXfXkce3OTN9eQ5+jfH/jss4rgBhjcrmY2A1lZwCHeUKNkDG3yLAxu6eTmAt27A505\n14iSMbTJ8zC4pbF1KzCec40oHUObPBOD27XYGvEaDG3yXAxu12FrxGswtMmzMbhdg60Rr8HQJs/H\n4G4Ytka8CkOblIHB7Ty2RrwKQ5uUg8HtHLZGvApDm5SFwe0Ytka8DkOblIfBbT+2RrwOQ5uUicFt\nH7ZGvA5Dm5SLwV03tka8EkOblI3BXTu2RrwSQ5uUj8FdM7ZGvBJDm7wDg9sWWyNei6FN3oPB/Tu2\nRrwWQ5u8C4O7AlsjXouhTd7nnuA2Pv44CgoKYDAYUFZWBrVaDZ1OB61WC41GI3OxDWc0Gm3Or7mf\nH+I/+gh3Cwuhlrs4cjk+I5K81+HDKBk+HHuffBL3P/+8NaSrhpzJZEJKSorclTotKSkJ/v7+Nj+E\n7mRnw/TSS/jHX/6i+PNrTOzNTZ+GHCQrKwvPPvtsQ3ZBJBnjww/jw8mTEbNjB8JPnrSOqjUaDcK7\nd0ei2Qw/Pz8YjUaZK3WO0WiEv78/Es1mhHfvbj2/Zjt2oFV0tOLPj2rm9Ej7hRdeQE5ODvr3749/\n/vOftjvlSJs8wFdffQUACC8rA8aOBf7xD2DGDODcOWDMGGDnTnxZVASVSoXw8HCZq3Wc9fy6d7ee\nD/7wB6BdO6B9eyAnR9Hn19hIPtIOCQnB22+/zXAmj2UwGKDVaoGICODTT4E5c4ClS38PuM6dodPp\nYDAY5C7VKdbz69y54nzGjAH+93+B8nIgJ0fx50c1q/cPkRs2bMCaNWtsPrZx40aMHz8egiDUul1y\ncrL1fb1eD71e72yNRE4pKyv7/Q+NERHAu+8C06YBjz9eEd4ANABGHTxYMfpWmFEHD0JTte6uXYEX\nXwRmzbJe6qfRaFBWViZThVQXQRDqzNDa1Bva06ZNw7Rp0xzecdXQJpKDWq2G0WisCO5z54DVq4HM\nTCAhARgxAmjTBqbycly9cgX44x/lLtdhV69cgSk4GP5+fsAvvwC7dgELFwLZ2RXn27kzjEYj1Gpe\nQ+KJ7h3MLlmyxK7teMkfeS2dToeCggLbnm/nzkBIiHX566IiNOvVC1Bgz7dZYCAOoEpPOz+/4vxm\nzrSen6GoCDqdTu5SyYUadPWISqWCSqVyVS1ELqXVaiv6uenpvwc28HsPOD39976wAnn7+VHNGjTS\nDg0NRWhoqKtqIXIpjUYDk8mEVD8/6IqKoGvXznqdtqGoCAZfX5SXlyv2BhtvPz+qGW+uIa937x2D\n3n5HpLehXPakAAAE/0lEQVSdX2Nhb24ytImIPIBb7ogkIiL3YmgTESkIQ5uISEEY2kRECsLQJiJS\nEIY2EZGCMLSJiBSEoU1EpCAMbSIiBWFoExEpCEObiEhBGNpERArC0CYiUhCGNhGRgjC0iYgUhKFN\nRKQgDG0iIgVhaDtIEAS5S5AUz0/ZeH7ez+HQ/vXXXxEZGQm9Xo+hQ4eioKBAiro8lrf/p+H5KRvP\nz/s5HNqrV6/GyJEjIQgCNm7ciNmzZ0tRFxER1cDX0Q3mz58Pf39/AIDZbIZarXZ5UUREVAuxDuvX\nrxf79Olj83bw4EFRFEXx0qVLYv/+/cW8vLxq2wHgG9/4xje+OfhmD5VozzPb73H06FE888wzWLVq\nFUaPHu3o5kRE5CSHQ/vYsWOIjo7G1q1bERQUJFVdRERUA4dDOyoqCt9//z06d+4MAGjdujWysrIk\nKY6IiGw51R6pjcViwaxZs/D999/D398f69evR/fu3V21e49RWFiIRYsWITc3V+5SXMZsNmPq1Kk4\nd+4cTCYTXnnlFURGRspdlsvcvXsX06dPx8mTJ6FSqbBu3Tr07t1b7rJc7sqVKxg4cCC+/PJL9OjR\nQ+5yXGrAgAFo1aoVAKBbt27YsGGDzBW5VlpaGrKzs2E2mzFnzhxMnjy5xvUcvnqkLtu3b0d5eTny\n8/NRWFiI+Ph4bN++3ZWHkN2KFSuwZcsWtGjRQu5SXCozMxPt2rVDRkYGbty4gX79+nlVaO/cuRM+\nPj74+uuvsW/fPiQmJnrd/02z2Yy4uDg0b95c7lJc7s6dOwDgVQOlqgRBgMFgQH5+PkpLS7FixYpa\n13XpHZEHDhzA448/DgAYMmQIDh486Mrde4TAwEBs27YNLvwFxSOMGzcOKSkpACp+Y/L1denPc9mN\nHTsW77zzDgCguLgYAQEBMlfkegkJCZg5cyY6dOggdyku991338FoNGL06NF49NFHUVhYKHdJLpWT\nk4OgoCBERUUhMjISTz75ZK3rujS0b926hZYtW1qXmzRpAovF4spDyC46OtrrAg0AmjdvjhYtWqCk\npATjxo3D0qVL5S7J5Zo0aYLY2FjMnTsXEydOlLscl9q4cSPatWuHUaNGAYDXDSqaN2+OhIQE7N27\nF+vWrcOzzz7rVdly9epVHDp0CB9//LH1/Grj0tBu2bIlSkpKrMsWiwU+PpzeRCnOnz+P8PBwTJo0\nCRMmTJC7HEls3LgRJ0+exPTp01FWViZ3OS6Tnp6Ozz//HGFhYThy5AgmT56My5cvy12Wy/To0cMa\nZA899BDatGmDS5cuyVyV67Rt2xajRo2Cr68vevTogWbNmuHatWs1ruvSRA0JCcHu3bsBAAUFBQgO\nDnbl7klCly9fxqhRo7BixQrExsbKXY7LZWRkIC0tDQCgVqvh4+PjVQOKffv2QRAE5Obmol+/fti8\neTPat28vd1kuk56ejvj4eADAxYsXcevWLa9qAw0bNgx79uwBUHF+paWlaNOmTY3ruvT3/Keeegqf\nf/45QkJCAFS80N5KpVLJXYJLLVu2DL/++itSUlKsve3PPvsMzZo1k7ky14iJiUFsbCxCQ0NhNpux\ndu1a63QM5PmmTZuGKVOmYMSIEQAqssWbfuhGREQgLy8PgwcPhsViwVtvvVVrxrj0kj8iIpKW9/yo\nIiJqBBjaREQKwtAmIlIQhjYRkYIwtImIFIShTUSkIAxtIiIF+f+NrMkYNBuG2gAAAABJRU5ErkJg\ngg==\n",
	"text": "<matplotlib.figure.Figure at 0x8a062d0>"
	}
	],
	"prompt_number": 19
	},
	{
	"cell_type": "code",
	"collapsed": false,
	"input": "",
	"language": "python",
	"metadata": {},
	"outputs": []
	}
	],
	"metadata": {}
	}
	]
	}