2 - CUDA Basics.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# GTC 2018 Numba Tutorial Notebook 2: CUDA Basics\n",
    "\n",
    "There are two basic approaches to GPU programming in Numba:\n",
    "\n",
    " 1. ufuncs/gufuncs (subject of this section)\n",
    " 2. CUDA Python kernels (subject of next section)\n",
    " \n",
    "We will not go into the CUDA hardware too much in this tutorial, but the most important thing to remember is that the hardware is designed for *data parallelism*.  Maximum throughput is achieved when you are computing the same operations on many different elements at once.  \n",
    "\n",
    "Universal functions are naturally data parallel, so we will begin with them."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Universal Functions\n",
    "\n",
    "NumPy has the concept of universal functions (\"ufuncs\"), which are functions that can take NumPy arrays of varying dimensions (or scalars) and operate on them element-by-element.\n",
    "\n",
    "It is probably easiest to show what happens by example.  We'll use the NumPy `add` ufunc to demonstrate what happens:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([11, 22, 33, 44])"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "a = np.array([1, 2, 3, 4])\n",
    "b = np.array([10, 20, 30, 40])\n",
    "\n",
    "np.add(a, b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ufuncs also can combine scalars with arrays:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([101, 102, 103, 104])"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.add(a, 100)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Arrays of different, but compatible dimensions can also be combined.  The lower dimensional array will be replicated to match the dimensionality of the higher dimensional array"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "c: [[ 0  1  2  3]\n",
      " [ 4  5  6  7]\n",
      " [ 8  9 10 11]\n",
      " [12 13 14 15]]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array([[10, 21, 32, 43],\n",
       "       [14, 25, 36, 47],\n",
       "       [18, 29, 40, 51],\n",
       "       [22, 33, 44, 55]])"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c = np.arange(4*4).reshape((4,4))\n",
    "print('c:', c)\n",
    "\n",
    "np.add(b, c)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the above situation, the `b` array is added to each row of `c`.  If we want to add `b` to each column, we need to transpose it.  There are several ways to do this, but one way is to insert a new axis using `np.newaxis`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[10],\n",
       "       [20],\n",
       "       [30],\n",
       "       [40]])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "b_col = b[:, np.newaxis]\n",
    "b_col"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[10, 11, 12, 13],\n",
       "       [24, 25, 26, 27],\n",
       "       [38, 39, 40, 41],\n",
       "       [52, 53, 54, 55]])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.add(b_col, c)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The NumPy documentation has a much more extensive discussion of ufuncs:\n",
    "\n",
    "https://docs.scipy.org/doc/numpy/reference/ufuncs.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Making ufuncs for the GPU\n",
    "\n",
    "Numba has the ability to create compiled ufuncs.  You implement a scalar function of all the inputs, and Numba will figure out the broadcast rules for you.  Generating a ufunc that uses CUDA requires giving an explicit type signature and setting the `target` attribute:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numba import vectorize\n",
    "\n",
    "@vectorize(['int64(int64, int64)'], target='cuda')\n",
    "def add_ufunc(x, y):\n",
    "    return x + y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a+b:\n",
      " [11 22 33 44]\n",
      "\n",
      "b_col + c:\n",
      " [[10 11 12 13]\n",
      " [24 25 26 27]\n",
      " [38 39 40 41]\n",
      " [52 53 54 55]]\n"
     ]
    }
   ],
   "source": [
    "print('a+b:\\n', add_ufunc(a, b))\n",
    "print()\n",
    "print('b_col + c:\\n', add_ufunc(b_col, c))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A lot of things just happened!  Numba automatically:\n",
    "\n",
    " * Compiled a CUDA kernel to execute the ufunc operation in parallel over all the input elements.\n",
    " * Allocated GPU memory for the inputs and the output.\n",
    " * Copied the input data to the GPU.\n",
    " * Executed the CUDA kernel with the correct kernel dimensions given the input sizes.\n",
    " * Copied the result back from the GPU to the CPU.\n",
    " * Returned the result as a NumPy array on the host.\n",
    "\n",
    "This is very convenient for testing, but copying data back and forth between the CPU and GPU can be slow and hurt performance.  In the next tutorial notebook, you'll learn about device management and memory allocation.\n",
    "\n",
    "You might be wondering how fast our simple example is on the GPU?  Let's see:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.1 µs ± 25.7 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)\n"
     ]
    }
   ],
   "source": [
    "%timeit np.add(b_col, c)   # NumPy on CPU"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "618 µs ± 5.9 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)\n"
     ]
    }
   ],
   "source": [
    "%timeit add_ufunc(b_col, c) # Numba on GPU"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Wow, the GPU is *a lot slower* than the CPU??  This is to be expected because we have (deliberately) misused the GPU in several ways in this example:\n",
    "\n",
    "  * **Our inputs are too small**: the GPU achieves performance through parallelism, operating on thousands of values at once.  Our test inputs have only 4 and 16 integers, respectively.  We need a much larger array to even keep the GPU busy.\n",
    "  * **Our calculation is too simple**: Sending a calculation to the GPU involves quite a bit of overhead compared to calling a function on the CPU.  If our calculation does not involve enough math operations (often called \"arithmetic intensity\"), then the GPU will spend most of its time waiting for data to move around.\n",
    "  * **We copy the data to and from the GPU**: While including the copy time can be realistic for a single function, often we want to run several GPU operations in sequence.  In those cases, it makes sense to send data to the GPU and keep it there until all of our processing is complete.\n",
    "  * **Our data types are larger than necessary**: Our example uses `int64` when we probably don't need it.  Scalar code using data types that are 32 and 64-bit run basically the same speed on the CPU, but 64-bit data types have a significant performance cost on the GPU.  Basic arithmetic on 64-bit floats can be anywhere from 2x (Pascal-architecture Tesla) to 24x (Maxwell-architecture GeForce) slower than 32-bit floats.  NumPy defaults to 64-bit data types when creating arrays, so it is important to set the `dtype` attribute or use the `ndarray.astype()` method to pick 32-bit types when you need them.\n",
    "  \n",
    "  \n",
    "Given the above, let's try an example that is faster on the GPU:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math  # Note that for the CUDA target, we need to use the scalar functions from the math module, not NumPy\n",
    "\n",
    "SQRT_2PI = np.float32((2*math.pi)**0.5)  # Precompute this constant as a float32.  Numba will inline it at compile time.\n",
    "\n",
    "@vectorize(['float32(float32, float32, float32)'], target='cuda')\n",
    "def gaussian_pdf(x, mean, sigma):\n",
    "    '''Compute the value of a Gaussian probability density function at x with given mean and sigma.'''\n",
    "    return math.exp(-0.5 * ((x - mean) / sigma)**2) / (sigma * SQRT_2PI)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.24195918], dtype=float32)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Evaluate the Gaussian a million times!\n",
    "x = np.random.uniform(-3, 3, size=1000000).astype(np.float32)\n",
    "mean = np.float32(0.0)\n",
    "sigma = np.float32(1.0)\n",
    "\n",
    "# Quick test\n",
    "gaussian_pdf(x[0], 0.0, 1.0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "46.7 ms ± 1.1 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)\n"
     ]
    }
   ],
   "source": [
    "import scipy.stats # for definition of gaussian distribution\n",
    "norm_pdf = scipy.stats.norm\n",
    "%timeit norm_pdf.pdf(x, loc=mean, scale=sigma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "13.8 ms ± 165 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n"
     ]
    }
   ],
   "source": [
    "%timeit gaussian_pdf(x, mean, sigma)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That's a pretty large improvement, even including the overhead of copying all the data to and from the GPU.  Ufuncs that use special functions (`exp`, `sin`, `cos`, etc) on large data sets run especially well on the GPU."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## CUDA Device Functions\n",
    "\n",
    "Ufuncs are great, but you should not have to cram all of your logic into a single function body. You can also create normal functions that are only called from other functions running on the GPU.  (These are similar to CUDA C functions defined with `__device__`.)\n",
    "\n",
    "Device functions are created with the `numba.cuda.jit` decorator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numba import cuda\n",
    "\n",
    "@cuda.jit(device=True)\n",
    "def polar_to_cartesian(rho, theta):\n",
    "    x = rho * math.cos(theta)\n",
    "    y = rho * math.sin(theta)\n",
    "    return x, y  # This is Python, so let's return a tuple\n",
    "\n",
    "@vectorize(['float32(float32, float32, float32, float32)'], target='cuda')\n",
    "def polar_distance(rho1, theta1, rho2, theta2):\n",
    "    x1, y1 = polar_to_cartesian(rho1, theta1)\n",
    "    x2, y2 = polar_to_cartesian(rho2, theta2)\n",
    "    \n",
    "    return ((x1 - x2)**2 + (y1 - y2)**2)**0.5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 1000000\n",
    "rho1 = np.random.uniform(0.5, 1.5, size=n).astype(np.float32)\n",
    "theta1 = np.random.uniform(-np.pi, np.pi, size=n).astype(np.float32)\n",
    "rho2 = np.random.uniform(0.5, 1.5, size=n).astype(np.float32)\n",
    "theta2 = np.random.uniform(-np.pi, np.pi, size=n).astype(np.float32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([1.0527985 , 2.3175654 , 0.37269294, ..., 2.4109268 , 1.1318382 ,\n",
       "       2.1343126 ], dtype=float32)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "polar_distance(rho1, theta1, rho2, theta2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the CUDA compiler aggressively inlines device functions, so there is generally no overhead for function calls.  Similarly, the \"tuple\" returned by `polar_to_cartesian` is not actually created as a Python object, but represented temporarily as a struct, which is then optimized away by the compiler."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Allowed Python on the GPU\n",
    "\n",
    "Compared to Numba on the CPU (which is already limited), Numba on the GPU has more limitations.  Supported Python includes:\n",
    "\n",
    "* `if`/`elif`/`else`\n",
    "* `while` and `for` loops\n",
    "* Basic math operators\n",
    "* Selected functions from the `math` and `cmath` modules\n",
    "* Tuples\n",
    "\n",
    "See [the Numba manual](http://numba.pydata.org/numba-doc/latest/cuda/cudapysupported.html) for more details."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercise\n",
    "\n",
    "Let's build a \"zero suppression\" function.  A common operation when working with waveforms is to force all samples values below a certain absolute magnitude to be zero, as a way to eliminate low amplitude noise.  Let's make some sample data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x10e634c50>]"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3Xl8FPX9P/DXOzcJhIQkQIBAuG+5wiWinIJixWptwdbiSevxtWr9VazWSq0ttfWoVetRb6tixQMVUUQEuQ3KfYZwhTMhhiuBXJ/fH5ndnd3sZq+Z+ex85v18PHgwOzu7894j7535zOfz/pAQAowxxtQVJzsAxhhj5uJEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjikuQHQAAZGdni/z8fNlhMMaYraxbt65MCJETbLuYSPT5+fkoLCyUHQZjjNkKEe0LZTtuumGMMcVxomeMMcVxomeMMcVxomeMMcVxomeMMcVxomeMMcVxomeMMcVxovfjm12l2H+8UnYYjDFmCE70OkIIPP3VLlz70lpc+PclssNhzO3VFXuw48gp2WEwm4qJkbGxoL5e4EdPL8eWQydlh6KUV1fsQVVNPW4Z01V2KLb20MdbAQDvzByBEV2yJEdjbyeqajBg9hcAgKJHLkFCvPrHu+q/whCdPFvTKMn3eOAznK2pkxSR/eXP+hQPfbwVf1u4XXYoypj2wmr85p3vZYdhawP/9IV7+cAPVRIjsQ4neg0RNVpXXVuPYyfPSYjG/u542zsZ3fLmOkmR2N/u0tNetz9af0hSJPZXfqYaQnhuL9p6RF4wFuJEr4mPa5zoAUBA+F3PmjZ/g3cy+myzM/6gzDD+saWN1lVUVkuIxN4qKqsx+OFFXuv+smA76uvV/xvnRI+GI/dJTyyTHYYyTp+r9bv+zdUhFdpjOoGaDic9yd/XcE19ZoXf9f/4YofFkViPEz2AZ78uwsEK/211Qv0fe0PV1wv0++Pnfu974MPNFkdjf48v2ul3/VFuUgzbvgBdpp/9erfFkViPEz2AisqagPeN+cfX1gWiAO61ZKwXlhXLDoEpgBM9gFdX7pUdgjJq6+ubvH/XUe4LHqqtQX40z9VyjzAWGk70zFA/fnZlk/cv2nbUokjs77UgByA9H1hoTSAO8MqKPbJDMBUn+hBsOXRCdgjKIPjv3cQam1t4QHYIyli6s7TJ+2drA9JUxYk+BFOeWi47BGXw4Ckmw4yX18oOQSrHJ/rKav9dAVn4BHdRslzRsdPBN2KO5/hE/8qKvbJDUMYb3E/eMNW1TV/Udik7zd0sjbKxpEJ2CKYJmuiJKI+IlhDRNiLaQkS/0da3IqJFRLRL+z9TW09E9BQRFRHRRiIabPaLiMaCTYdlh6CMBz/aEtJ2Rce4500w987bGNJ2y4K0PbPQXfPiGtkhmCaUI/paAL8VQvQGMALAbUTUB8AsAIuFEN0BLNZuA8AlALpr/2YC+LfhURso1H7f3JXNOFf9e5XsEGLeB98fDGm7QKU7WPgCjehWQdBEL4Q4LIT4Tls+BWAbgPYApgJ4TdvsNQBXaMtTAbwuGqwGkEFEuYZHbrFb3/xOdggxLZx6ISeqAg9QY+Hh6qrBzVtXIjsE6cJqoyeifACDAKwB0EYIcRho+DEA0FrbrD0Afb+wEm2d73PNJKJCIiosLZVz+nkmwC/41j9NarRu8fZjZodja29/u192CMoI56L2i9+o3f/bCL/93wbZIUgXcqInouYA5gG4UwjRVHuHv3PJRt9cIcQLQogCIURBTk5OqGEY6rEvGtcR+eZ3Y5GalIDmyTwnSzhKTzW+KLj29+Ox6r5xEqKxt00HG4/b6JDZDF2y0yREo6abLujsd72qlSxDSvRElIiGJP9fIcT72uqjriYZ7X/XIW8JgDzdwzsAiMkC2st2NT6TyGuVCgBYdPeFVoejnNbpKcht2Qxt01Ma3cfXPAKr9ZNsXr9hGL66Z4zf7blba/geuKyP3/XFZWp2Vw2l1w0BeAnANiHE47q75gOYoS3PAPCRbv0vtd43IwCccDXxxJqm+iDntmxmYST2V34mcH30j24f1Wjd4YqzZoZja/5OibvkNAcArH9wYqP7th7mQnLh+MfVAwAA7/16ZKP7nl2iZiXLUI7oRwG4FsA4Ilqv/bsUwBwAE4loF4CJ2m0AWACgGEARgBcB3Gp82NHbd/xMo3U3BjidY8G9vsq7D/2Xd1/kXm7j54i+pi60fuJO5Dvb2e8v7eVezkhNarT92j3lpsekkp8M6QAAKMhv1ei+90Ps7WQ3ofS6WS6EICHEeUKIgdq/BUKI40KI8UKI7tr/5dr2QghxmxCiqxCivxCi0PyXEb4SP3NFntehpdftB31O7/gUOXQdtSYwl7E9va/D+Ls+whpc+5J3f+5hnb0nA19+71iv2/y1DCzYyOF5tzQ+qleR40fG6g31+YXPau599PTEl7usDMfWEnz6d/vmIn8XHFmDU2e9e4N18vnR7JDpfZtHJAfmO5H6sM7ef+OtWzQ+21SRYxO9bzvohT1y0C7Du11+bK/WXrefWsyJPhQL7xyNOJ9E/8AU77OjQDN6scYyUhObvH9PWeNmSNbAd0DkmzcO97qd25ITvdKW7Srzuv3CtUMabZOe0vQfGPOvV9v0Ruu6tW4uIRI1+LbZA0B+VqqfLVkwSQneKS8hPg5T+tt+PGdQjk30zy31vrqekhjvd7ubR/MFWmadH3x6L825sr/f7X42tKMV4TiCqxeOyhyb6EOV5jNwiqfCa0zlGiFWe9dnspHmKf4H7iXGc40boyT7HOU31VXYrjjRB+E7I9L+cv8zyTtZoFISvgbkZZgcif35trf3zm3cDAYAI7pk+V3PAnv/1vP9rve9nvT79zdZEY6lONEH0To92et2RSUX5PL1r69Cu0jtO+xc1eHm0fCtQ981x/+1jb7tvH8ADvABSCPHTnkPyuse4nWihVuOmBGOVJzog5jQu43XbZ6RqrE3V4dW0OxSn4teCzbH5IBpqUIdsON7gfayf/F0l77+6DM/QgsHd67gRA9gZBOnwWnJ3hdp4+P4LYuUb+107ktvHC793Nhnm0M/Mm+foXbJE85aAF69YWjA+1KTvC+GicaFOJnOn6/oF/K2r6/kgT5NmTY0L/hGzBAL7hgtOwRTcaIHkJzgv2ulP/d/sNnESOzv58ND7/ZXxZNmNGlwx0zZIThGyyCD0uzOkYmea9aYx9/gHj0uHBe6YD1rpg/jvvSheu4X4U1drVpTmCMTPdcGMU5FZXh9jgs68VFqqIL8ZuKvAQZTscYm9wtv9OuA2V+YFIkcjkz0D/pcjQ/m2hGdTIrE/sLtbjqxT5vgGzEAjft3MxYpRyZ6vVD61v7xR/5no2HAquLjYW2fEO/9lTt+uvEUhKxBuD1BuEmSBeL4RD/lvOCndL7JaeXusgBbOk9VdXQXVHnSdY/v9v8Q1eMf/mSbQZEw1Tg+0XeKoArgur3R/UGq5E+fbI3q8TuOcO0glxW6iqr/nDYw7Me/vGKPkeHYWiRnN6/fMMyESGKD4xO970xIoTh5Vq0r8jK9tJyTk8sPuusd43uHdi3jN+O7mxWOrc2aF369mj7t/NcVUoHjE72/2unBnK3h+U79uahHTvCNAAzqyMXN/NEfkTcLUDbb101cRtuvuT5VQEOR3Tw5+EY25fhE71uGOBSRfImcINTufs9cE16fZifyLRcRiJPrt4QqkuZZ1Tg+0Yfqrgk93Mu+FQZZA9+pGKPdjjEjtG6h7pF6qDjRh2icz/yxjDF7uO58bt5yXKIvOhZZL4+2DplEmDHVDO4U2TWhTzeqU0bbcYl+wuPLInpcDp/+NVKnmzhkQm8+44nGuVou8GaW3JaRNRXe9tZ3Bkcij+MSPTPOf74p1t2KfLg+zzQFvLpiryHPU8YjjZkfjk703UKcWoz599fPtruXgxXgaspxBSdjDteGkgpDnudwxdngG7GAWqUlyQ7BFI5O9Jf2axvxY7lGi7dOEQw8c+HJXIAvt3lKQfjOBxuOF73OspxpTZj1l/RW3DvOwEhih6MT/R1RjCpcsqPUwEjs75YxXcPavp3u4vZL3/DoWH2X3fun9I74eeZvOGREOLb25prQ5jD2p1lS6JMQ2YmjE71vsbJwbD100sBI7C8rzFGF+oFq760rMTocWwt3ZqlHfhz69I1O8LHux+6T/7tAYiSxw9GJPlz6ujhcQCo6D13e173MbfTeUkIsf+ByDc80FVC/9i1lhxATHJXo66Ls3cHT4BlnVLds2SEoI9j0jYw5KtFHewE11PojTnCKK3gyZhuOSvTRap/JNVpc9pZVyg6BMRaioImeiF4momNEtFm37iEiOkhE67V/l+ruu4+IiohoBxFNMivwSOjryLdsFn7Vv7E9efSni9FdIp08Dd6uozz5CjNXKEf0rwKY7Gf9E0KIgdq/BQBARH0ATAPQV3vMs0QUM/2V9OUPwu0OyLx9ZfAUgNFeP7GzQyeMHeRUU8fVVY2yRJGpLoMmeiHEMgDlIT7fVADvCCHOCSH2ACgCEJPzc/GF1eg8+eUu9/Ir1w+N+vnqHHxEb7Q4vjhrmNJTagyMjKaN/nYi2qg17bg6/rYHoJ+Vo0RbF3MSo+hDz7zlGlDZ08lH9HoPRDhYqmebFu7lQxVVRoXjSPoJdOoVOQCJNNv9G0BXAAMBHAbwmLbe36GE33eKiGYSUSERFZaW8ihTO0tLCn+WLsD7B8LJiX5P6Wn38i9GdIroOcbq5ksY/eiSqGOyq6Jjp4NvFMRPhnRwL6/cHXk5hVgSUaIXQhwVQtQJIeoBvAhP80wJgDzdph0A+B2TLYR4QQhRIIQoyMkJba7RWMOjYxvkRVjnZnAnzwjQegc3K//98x3u5eSEyI69uF5Qg7nfesof/O2q0Ka29JWg60atSkmJiL5VRJSru/ljAK4eOfMBTCOiZCLqDKA7gLXRhRi7jKo46FT6tuTnlu2WGIlcZ6o9tegjHfx0flcegAYAL+rqJk0dGFmrsYoD0ELpXvk2gFUAehJRCRHdCOBRItpERBsBjAVwFwAIIbYAeBfAVgALAdwmhFB2RgVV2u9k6Zyd5l5evO2oxEjs76Ie9jwrNlO4pSRUFrRxVQgx3c/ql5rY/hEAj0QTVCwbmJeB9QcajuRPna2VHI29Zelqf+88Gn3bKmPMP+56EqbrR+W7l+foJt5g4ftpQV7wjRhjUXNMoq+sNuboO9IeJqwxVWt/MxZrHJPo9QMfbhsb+ahYBa/TMMaaoMK4BMck+t26vsqpfFTOGAvR4ROc6G3jhlcL3csjumRF/Dz9eSIDr8FN/H4w9dn/NN4xiV5vYF5GxI9tnR79cH+7e2GZZwLquyZGPu8uM09VtbK9mi1xk64WVtEx+1cXdWSi5wlEovPMkiL3cs+26RIjsT99Qp47c0RUzzU03zPS+P3vnTcPr77UdbR/4/dd6qk5dO+8TVE9VyxwZKJn0Tl9ztODKSXCIfsu+no3TqxJv7/cM4FLtBPbTOrb1r1cW+e893LZrjL38vldI2+eBdQ7GOREz6KS1Tw5qsdP101s7cTCZvoiqs2To+skME5X2OzfXzuvpMSirUfcy7++iOeb0ONEHyUnHoUaSf/2bTts/7bQcOkncMlITWpiy+D0R6FHTho7mYndGFE6WyWc6KP0xVau0RKNywe2cy//6OnlEiOR4y8LjBtd3cbhHQWW7vSUO4/2R1M1nOgjMKqbp/3vRGVNE1uyYPSFzVh0nF7E60C5p7+7E5sBm8KJPgIPXtbXvbxyd1kTWzLGZMhuzkf0epzoI9Alx3MU+uF6NSYmYEwlKtaUjwYn+gjwfLOMMTtxRMbiKf8YY07miET/Q2W1e5ln4jHO8M6tZIfAGAuBIxL9t3vL3cv/nDZQYiT2d7bGM2R/9tS+TWzJmL09MKV38I1swhGJ/tWVe93L3L82OjuPegY18SQsxumUlSo7BObjhlGdg29kE45I9DW19bJDUMblT69wL3PHhujoR1XPGJkvLxDF/Ga8MRVV4xSqd+OIRH+GS7aaIr1ZouwQbO2kbnJ5o6ZVzEj1fCb6ZjYnMePsaN46e1cDdUSiN8O0oTyxdXqKMYleP9K4ps45Z1+lpzz1aPq2M6bc89s3e0odf7b5sCHPaTdmVJ787f82GP6cVuJEH6GJfdrIDkEZF3Tz9IR6a81+iZFY67J/eWr7tMuIrkSxS+9czw/GXXPtnZzCoT97Gd45uhLFKuJEH6HWLZxdQMpIl52X615esMk5R6FnazxnLzwILzonqjw1p1pyk2Ij/O2KUJzuneNSxdHRt0+v2VPexJbq4uQUnam6TgJGXe9QCSf6CLVK83TTLC47IzES+0vlP0wWJafX3w+GE32Eclt62lR/+twqiZHYXyr3x2fMVI5K9GN7mlP+4PiZ6uAbMcaYJI5K9KpN+MsYY6FQPtEf0/VVjuOhnIZJT+HmFsbsQvlEv7LouHvZqL7KDGjt8PlJmTP8fHhH2SEYQvlE/+U2z+TdMy/sIjEStSRwMxhzgKzmye7lwr327fqrfKJfXcxH9Ga4flS+7BCU0SfXmPIHzHgX9ch2L+88elpiJNFRPtGXneYeMUY5UF7pXr56iLG1fvJaOfdHeHzv1oY+35s3Dncv82C+6Azp5Jlc5/jpcxIjiU7QRE9ELxPRMSLarFvXiogWEdEu7f9MbT0R0VNEVEREG4losJnBM2u9uXqfe9no69qzL3fWJCb19Z4EPLlfW0OfOzHe8+F8u/cHQ5871pnZSeCxRTtNe26zhXJE/yqAyT7rZgFYLIToDmCxdhsALgHQXfs3E8C/jQkz9un/cFW1StcMRgZn+kyHTQizvKjMvawfZW2EAXkZ7uWNJRWGPncsKtMdaXfI5Alc/Ama6IUQywD4XoWYCuA1bfk1AFfo1r8uGqwGkEFEuVDU9GGeK/J1DjhF3lhywrTndlrX13cLD7iXUxONPQpNSfSUlNDPrqaq5bs8P5o8U5d/kbbRtxFCHAYA7X9XI2N7AAd025Vo6xohoplEVEhEhaWlpRGGIdeU/p7fsHoHJHoz6fO8E9qV9cXbWqaaV9CsorIm+EY2t/2IZ3rLeyf3khhJ7DL6Yqy/wzK/f7VCiBeEEAVCiIKcHHNKE5hNP9K23jnzZZhCf0T/ns1n8wlF6SlrLuxVO2AaTX2J4vzsNImRxK5IE/1RV5OM9v8xbX0JAH13jA4ADkUeXmxLSfS8fedqnTltm1E66k659UdoLDrVDpixq7jUvt0erRJpop8PYIa2PAPAR7r1v9R634wAcMLVxKOi5ARPW+iCTUckRmJ/zXTtyi8t3yMxEmY3Tp3DIByhdK98G8AqAD2JqISIbgQwB8BEItoFYKJ2GwAWACgGUATgRQC3mhJ1jOid28K9vHALJ/po8AxLjJkn6OV+IcT0AHeN97OtAHBbtEHZhb6L4bKd9rygzBhTn9KHUfreG7Mu4avxjDFnUjrR19R5Ev2ILjwzPGMsfL++qKvsEKKmdKKv1fV57JDp3FoqRruwhz27wzIWiebJ9p/TWOlEf6jCM+kIX+wzziDdEHvGVDd1oN8xn7aidPZ7bulu93LLZuaNPnQCfS2fiX3aSIxELTktkoNvxKRqnuzps3K2xp7jZZRO9Gv2HA++EQtJcdkZ97LT6tKY6doRnUx53vsv7W3K88aygSadaTZL8jTdHD9jz7LnSid6LktgnD26RC/8V7WI2h3jupnyvLGsV9sWwTeKQP8OLU153liW3dycCqj6InFnztWasg+zKZ3os0z64J2oTveraVbNMf1F3opKex45hWLnUU+Jh26tm5uyD/1UjwcrqkzZR6yZ3M/8QrkXP7HM9H2YQelEb2ZZXZe5M0eYvo9YcOSE58K2WYlef+pdbtNT5FDoJ3AxaxqDfu09R/R2PQoNxW5dnZu0JPv3jjGL0oneCkkJzngLH/p4q3u5bztz5jhN0PWMmlt4oIkt7e31VZ5E3zXHnGqL+uYGlStY/uqNde7lgvxWTWzpbM7IUiZqkeLpzeOEOuoAEBdn/sVY/RmEyoyeqcufk2fVrUm/XzePcXIip7NA+J2Jkn5GG/1IXBadnUe59KxRrnlxjewQTKM/W0nisTIB8TsTJX1XQyfU/rbK8dPWTMzB1MGJPjB+Z6Kkn2XqznfWS4xELccsmoGJqcOKJkW7ckSit2rC4JIfKoNvxBhjFnNEok9JsKbbFTfdMMZikSMS/ZD8TEv2U8sXYxljMcgRif7hqf0s2Y++q5eqzOr3zVgs0xc2syNlE72+T3s8X6QxzJyrzpMdgjKsunbEovcLk4rPWUXZRP/V9mOyQ1ASJyfj3DvZ3Okts5tzCWSjXNAtW3YIUVE20XP3PHMQ+OwoGrW6C/Y9Tapc6cLzJBsnNyNFdghRUTbR/+ebYtkhKKOq2jPZAk+UEZ0jJz2lHcz+yZw6sJ3Je4gdrdLMrVTbqZXnTNaOReKUTfS7S88E34iFRJ+czDakkzU9pGTRj6Q2+wKfvlTxiSp1690A3tU6zXbAhuNllE30VvpZQZ7sEEz1/f4fLNuXWbMExQp9om+dbm5zgL5gmsr1/QHgpwUdTH1+/ee2bGepqfsyAyd6A1wzvKPsEEy148ip4BsZZFhntUvNFu4rl7LfWrMK38cIs68d6csrbDp40tR9mYETvQHys9TuW/78Muuud+ing1Ox7PO+43JO+z/47qCU/VrFynpmH284ZN3ODMKJ3gBOmXzECknxnnIVVTV1TWxpT09/VSRlvyq+l/pOAmN7tZYYSezjDGWAegWPPGXRD25buPmIxEjMISvhrtp9XMp+zbRDN/duskX1rOyKE70B9Gn+qIU9VFSUmeaZseveeRslRqKWrYft164cjH7CetY0TvQG0Lclb7fwwqWK2up6ovCMXawpj32xU3YItqF8ojdrIms9/UTMM15ea/r+VGbFHKpMDSsVbI4yi/KJ/uoh5vavBYBEnsKMMRbDlM9Q7TKayQ6BMcakimoMNhHtBXAKQB2AWiFEARG1AjAXQD6AvQB+KoSwbmilj5bNEoNvxBhjCjPiiH6sEGKgEKJAuz0LwGIhRHcAi7Xb0gxQfEi9lW4Z01V2CIyxCJjRdDMVwGva8msArjBhHyHTXyhl0bm4TxvZISgjJdGaVtMBHawr9qW6uyf2kB1CxKL9tgkAXxDROiKaqa1rI4Q4DADa/5YPWdtdetrqXTqC6vVSrDSxT1tL9vPgj/pash8nsHOni2jrpI4SQhwiotYAFhHR9lAfqP0wzASAjh2NLQo2/rGlhj4fa9DfwlKwqkuMt6YbKU+jaRz9Z/bDmWpkmlwD30hR/UQJIQ5p/x8D8AGAYQCOElEuAGj/+53TTwjxghCiQAhRkJOTE00YzESnznrqmFvVDPa3q/pbsh+ZsixKEk7I891bN7dkP78cme9ePnrKXiPgI070RJRGRC1cywAuBrAZwHwAM7TNZgD4KNogmTxvrdlv+T5Tk8ydkCMW3DG+uyX7yWjm+UHRT2OokvaZ1nSh1hcvtFvtoGj+otoA+EAbyZgA4C0hxEIi+hbAu0R0I4D9AK6OPkwmyzNL5FRbVNFJ3dlRixRruv121E3mXlMnoGLtr37trG9SrLHZj2bEiV4IUQxggJ/1xwGMjyYoOxrVLQsriuz1Kx+Kk2ftNz9mrCqWPL2lgDoX00/r5m392VDrZ3h7bmkxZl5on+7G9r2MHGOSbHxFPtboyz6r1Nygr58uQ02tOole/72QMfq9/Iy9pmbk7GSQ0d35grJR9Bd9l9pwfs5AZJewfvZrdZrh9FNAcM+i4DjRG+QXIzrJDkEZE3p7Bmadq1XniP7Oueul7t/KKSHNdrCiSnYItsKJ3iA8naBx9Edox212isysUceD98KidHZqL6lyZT1/CQ1zTsG5Tln0pj6zQnYItqJ0or/Kglr0/hTuk1asUzlnOdEzFjWlE72sazTvFh6Qs2MFLS8qkx0CY7andKIf3T1byn6/3qFOTxHZVheXyw6BMTf9nMZ2onSiT4iT8/LKTp+Tsl/GmLnOs2nZZ6UTvX7gDYvOhT14nABjwzq3kh1CRJRO9JzmjTPcpl9wxqO2jZRg08FZyn0Ddh495V7unJUmMRL7O3LCM5JzQAeektEo43tZOxfPJf2tmeTECfRjPOxUnkO5RD9//SH3sp0mBohFVbqujWTPA5mYZHUzWGqSgiUrJSnI95zZVlTVNLFlbFEu0W88eEJ2CMo4o6sQaHXTzf+b1NPS/VnJ+lHU6v5K92hjzaQjLq1bJLuXq21UnkO5RL9MoSJYsm09dNK9nGBxO6+sUc1m0Q/Z75Obbum+s5ure2b7K4tLBWc19yT6kh/sU29HuUQvU5v05OAb2YjMwUpDOmVK27cZ9BNVDMiz9nqHvrLq7tLTlu7bbIkSa0xtOWSf1gNO9AbKbq5Wop+/4VDwjUyi2vWVxdv8Tp1siY6tPLNMHbTRUWgo8iyaRtCf3Jb2GTzFid5Alw9oJzsEZTRP9kx+dsJGF70C2Vcub3aptrqEpNpF9ZwW8g6u3i0skbbvcHGiN9D04R1lh6AkFaqBxsprmPPZdtkhRE3frVFmSYKvtss7SwsXJ3oDpSV5jkLtdEU+1pVX2r8mfW2MJPotugvsdrX3uOfsyOpOAnbF75KB9GfFFVX2T06x4tONh2WHELUaGw2uiX2KtT9ZgBO9gfTtn1xmxziPL9opO4SoPbNkt+wQlME1rMLHid5ApMv0/GVkzBw3vPqt7BBshxO9SaqqeWYkxsxgp4FKsYITvUkeXbhDdgiMMQZA4USfkij3pS3cckTq/hlj5pjQu43sEMKmbKL/44/6yg5BGekpCcE3Yswh7FieQ9lE3621tVXtzLbl0AkcKK+Usu/WNp0nMxJf7ziGszXqXV+57LzcsB+zt+wMdhw5FXzDEBWXnsauo8Y9nyyT+9mvvr9SiV5fIVBf30MFU55ajtGPLgnrMW+t2Y9jp84G3zCIiX3kn6ruO34GH35/sNH6unqBZ78uwr8W73IPUhNC4JUVewKWTli1+zhWFx9vtH7roZNwdfkJAAAPUklEQVS47pVv8dD8Le51Gw5UYImNRkAGMqhj+EehY/7xNSY9uazR+rV7yrGyqAxCCLy6Yg9OVIZWomLcY0sx8YnGz2c3WSFUA12y4xjWH6iwIJrQ2Pqc/Pjpc3hs0U788Ud9kBAXhx4PfOa+Ly4Ginr89t0N6JSVijvGdw+4zX++KcafP92G60flezU3bT54Ap9vOYLW6SnI1R1RHyivRF6rVJyrrcPsj7fitxN7uEunHqyowvNLd+PBy/rg6Klz+P0HmzC3MAMf3TYKt7/1HT7ZeBj/vWk4RnXLdj/fiaoazPlsGx68rC8K95Wj6NhpXD+qs1eMN4/uYtRbErGL/v41AKBDZjM8unAHqmrqcH63LFRV1+H1VfsAAC8sK8am2ZOwZk85Zn+8Fd/tr8C/pg/Coq1HUXrqHK4Z3hEnqmow/cXVAIAxPXOQGB+Hh6f2wzNLityTun+84RCymyfjnkk9MfWZFQCAvXOmoLK6Fje9VogOmc0w58rzEBdH+PD7gyACpg5sb/2bEoZztZ6zlLfW7Mc1WrmOrYdO4tNNh3DPxT1RdroaQx/5EgvvHI1ebT2llM/W1OEPH27G/9aV4MmfDcSdc9cDAIblt8LaveVYvP0YOmQ2wwNT+iBNV6OosroWw/+yGPdc3BNf7/D+sRRC4G8Ld+DKwe3Ro00Lv99nAHht5V50ykrFmJ6tcbCiCrPmbTTl/QmHfmrGtXvKsWr3cSQnxuHdwgMoLj2DB6b0xp8/3QYAePqaQbjsvIYaWK+t3Iv87DRcJGH+ZRIx0N+7oKBAFBYWhv24e/63Ae+taygsdP2ofLyyYq/7vj1/vdSrX7tV8md92mjd3jlTAAAVldW4/OkVeP7aIeidm46zNXXo9YeFjbbvnJ2GPWWBi2B9cdeFuDiMI6MJvVvjS131xOSEOJzzKdHw24k98Jg2MOma4R3xlx/3d7+WXY9cgkQJQ839vZehyGmRjNJT5xqt3ztnSljP2bddesCSAVMHtkMcET7QzjJcn/Gd73yPD3WznL18XQEOnziL+z/YDAB45Mf98PPhnUKOwSivrNiD2R9v9VqX16oZDpSb01Vx40MXY+brhVhdXN7ovr1zpuCOt793V0fdPHsSfv3GukZlsbvkpKG4tOHv4MIeOSguPe3VtbJFcgI2zZ5kSvxNEUKg830LQt7+ioHtvL4TLqlJ8Vh7/wSvAn7hIqJ1QoiCoNvZOdH/7r0NASvIuf7wrBZpcoplO/98iYRZkdR8L9+/9XwMjqAZJVpfbDmCmW+ss3y/oZg2NA/vfHsg7MfdO7kXbhlj7cQjLkZ+N6PJVaEmelu30W8ssU/hfzuTkeQBNcs+J0kqwhXLXQIjSfJAw5mvCp5fan55DFsn+u0G9ghgsWd4F2vnqbVCxyy1OgnINKijtTN1meWvFpSONi3RE9FkItpBREVENMuMfSTEyb/g6muEgslJlvO7ZgffyGbSUxKl7Fd+A63xYqC/hW2YkuiJKB7AMwAuAdAHwHQi6mP0fmKlxrfeRT1ayw5BGfH8l2yYGDwmilos9KyzC7OO6IcBKBJCFAshqgG8A2CqSfuKKT8eFNvd7OwkztYNi7FFRg80s8XiGX2sMutPqT0A/RWWEm2d8lpwuQDDxPMfMmtCRqpaE8ibyaxE7+8v1KudhYhmElEhERWWlpaaFIb10qLoE8u88XvJmDHMSvQlAPJ0tzsA8BoxIIR4QQhRIIQoyMmxfqQYi32yLlwyphqzEv23ALoTUWciSgIwDcB8o3eS16qZ3/Wy+iozxlgsMuXcWAhRS0S3A/gcQDyAl4UQW4I8LGwpCfF+11852BGXAxhjLCSmHfoKIRYIIXoIIboKIR4xYx9tApTPzW3p/0ifhW/Olf1lh6CMfB4spYx+7dODbxSiJfeMMey5ArF1G0egOhddctQYGh0LruDuooaZeaGcuiwu7996vtT9G+nuiT2k7n/W5N6GPZcVpRxsnejzMv0fIeVnyU30H4TxB9U1Jw1PTR+EtgHOTvz1FZ4+LM/Plo3dMKozJvRujewA9bPn3z7K7/p/ThvoXk5J9N88pqqstCQsv3es+3Z6SoIhTYFXD+kQ8udmlsEdM3H/pcYkqCn9A09kMnfmCPRq2wKA95SeM0YaV7VTdsfbC7p7j9r+5chOKOiUGfbI+HByRTRs3X+tY1Yq3r55BJIT45AYF4f0ZgkgkPR6IvpJHkZ3z8ZtY7thQIcMrNxdhjbpKejbLh2d71uAgXkZ+PC2hmR7Sb+22He8EhMeX+p+7MguWXj5uqE4dbYGVTV1IBD2l1diVLcsvL3WM0zhmWsGY1yv1rj59UIsLypDckIcvvvDRHf3xHO1ddh++BSaa33826anoLK6DjktPHW/X71+qLvJq2fbFvjNO+vNe4PC8OXdF7nfk4F5GXhq2iD8UFntrhO/4cGLMeBPXwR8/KNXnYffzduIQR0z8Mp1QzHxiWUoPXUOn995IbrmpOH5ZcX4++c7cPvYbnh6SRH+dtV56JCZiqSEOFTX1mPt/ROQkhiPnBbJeH5pMQBg1X3jcLamHj9UVuPKZ1e695WaFI+aunosvnsM4uMJo+Z8hQem9MbVQ/LQMjU2ehBN7tcWjyzYhsd/OgA92rRAu4xmOHOuFrf+9ztsOti4SOC2P03GXxZswxur9+HhqX3xh4+24OGpfXHtyHx8OutTDOqYge/3eybY+L9x3TC8SxYW3DEaByuq0C6jGbr+vqGk7+yp/dCnXTrunbcJn95xAaY8tRwA8IfL+uDhTxpKKI/uno1vdpVhynm5+HTjYQDAh7eNwhXa5+1yfrcsU96fcMydOQI/e2E1/nplf0wf1tG9Xl/ZckBeBjYcqMDz1w5BWlICWqUlYd53JXhp+R4UPjAB2bra+2aydZniWHbR35dg6oB2uGtiD7+jEl3vu+99+i/JWzcPD1jvpb5egAioF56BRUII1IuG4e6hjoSs18pIxPmcOYz4y2L8YkRH3D4u8KQpVqmrFyA01DZxvS7Xurg4Qn29wNJdpbj+lW/dj7l7Yg/cNrYb4rX7XY9duPkIfv3mOmyZPcn9Q1hfL9zP43of3lqzHw9/shVbZk9CXBwFfG9dn1duyxSsnDUOQnjey7p6EZODvvSv08X1fVy09ShmvrEOzRLjceeE7vjVRQ3NTa7Xon+s6zHjH1+KyX3b4p6Le3p9Ri43vVaIhDjCc9cO8dq/673bPHsS+v3xczz3i8GY1LcthGj4rD9cfxB3v7sB2x+ejIS4ONTVCyTGk9d3XrZA72VNXcP75fueue7Xf0+i4Yh69Cpyffll1dO3u6cW78Lji3bid5N74tYx3Uzf323//Q6fbjqMl68rwLhesVsKOBa9sWovSk+dw90X95Qdim2Fmuht3XSjon9cPUC5+W6tdPPoLjhzrhY3+EyHaJaHr+iHTlmpXMwuAteOzJcdgmPwET1jjNmUI2aYYowxFhwnesYYUxwnesYYUxwnesYYUxwnesYYUxwnesYYUxwnesYYUxwnesYYU1xMDJgiolIA+yJ8eDaAMgPDsQN+zc7Ar9kZonnNnYQQQedijYlEHw0iKgxlZJhK+DU7A79mZ7DiNXPTDWOMKY4TPWOMKU6FRP+C7AAk4NfsDPyancH012z7NnrGGGNNU+GInjHGWBNsneiJaDIR7SCiIiKaJTuecBBRHhEtIaJtRLSFiH6jrW9FRIuIaJf2f6a2nojoKe21biSiwbrnmqFtv4uIZujWDyGiTdpjnqJQ5xc0GRHFE9H3RPSJdrszEa3R4p9LREna+mTtdpF2f77uOe7T1u8gokm69TH3nSCiDCJ6j4i2a5/3SNU/ZyK6S/tebyait4koRbXPmYheJqJjRLRZt870zzXQPprUMH+h/f4BiAewG0AXAEkANgDoIzuuMOLPBTBYW24BYCeAPgAeBTBLWz8LwN+05UsBfAaAAIwAsEZb3wpAsfZ/pracqd23FsBI7TGfAbhE9uvW4robwFsAPtFuvwtgmrb8HIBbtOVbATynLU8DMFdb7qN93skAOmvfg/hY/U4AeA3ATdpyEoAMlT9nAO0B7AHQTPf5Xqfa5wzgQgCDAWzWrTP9cw20jyZjlf1HEMWbPBLA57rb9wG4T3ZcUbyejwBMBLADQK62LhfADm35eQDTddvv0O6fDuB53frntXW5ALbr1nttJ/F1dgCwGMA4AJ9oX+IyAAm+nyuAzwGM1JYTtO3I97N2bReL3wkA6VrSI5/1yn7OaEj0B7TklaB9zpNU/JwB5MM70Zv+uQbaR1P/7Nx04/oyuZRo62xHO1UdBGANgDZCiMMAoP3vmow00Ottan2Jn/WyPQngdwDqtdtZACqEELXabX2c7tem3X9C2z7c90KmLgBKAbyiNVf9h4jSoPDnLIQ4COAfAPYDOIyGz20d1P6cXaz4XAPtIyA7J3p/7ZC260JERM0BzANwpxDiZFOb+lknIlgvDRFdBuCYEGKdfrWfTUWQ+2zzmtFwhDoYwL+FEIMAnEHD6XYgtn/NWpvxVDQ0t7QDkAbgEj+bqvQ5ByP1Ndo50ZcAyNPd7gDgkKRYIkJEiWhI8v8VQryvrT5KRLna/bkAjmnrA73eptZ38LNeplEALieivQDeQUPzzZMAMogoQdtGH6f7tWn3twRQjvDfC5lKAJQIIdZot99DQ+JX+XOeAGCPEKJUCFED4H0A50Ptz9nFis810D4CsnOi/xZAd+1KfhIaLuLMlxxTyLQr6C8B2CaEeFx313wArivvM9DQdu9a/0vt6v0IACe007bPAVxMRJnakdTFaGi/PAzgFBGN0Pb1S91zSSGEuE8I0UEIkY+Gz+srIcTPASwB8BNtM9/X7HovfqJtL7T107TeGp0BdEfDhauY+04IIY4AOEBEPbVV4wFshcKfMxqabEYQUaoWk+s1K/s561jxuQbaR2AyL9oYcCHkUjT0VtkN4H7Z8YQZ+wVoOBXbCGC99u9SNLRNLgawS/u/lbY9AXhGe62bABTonusGAEXav+t16wsAbNYe8zR8LghKfv1j4Ol10wUNf8BFAP4HIFlbn6LdLtLu76J7/P3a69oBXS+TWPxOABgIoFD7rD9EQ+8KpT9nALMBbNfiegMNPWeU+pwBvI2GaxA1aDgCv9GKzzXQPpr6xyNjGWNMcXZuumGMMRYCTvSMMaY4TvSMMaY4TvSMMaY4TvSMMaY4TvSMMaY4TvSMMaY4TvSMMaa4/w+UKVQKh+vmXQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x114e1c470>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Hacking up a noisy pulse train\n",
    "%matplotlib inline\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "n = 100000\n",
    "noise = np.random.normal(size=n) * 3\n",
    "pulses = np.maximum(np.sin(np.arange(n) / (n / 23)) - 0.3, 0.0)\n",
    "waveform = ((pulses * 300) + noise).astype(np.int16)\n",
    "plt.plot(waveform)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now try filling in body of this ufunc:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "@vectorize(['int16(int16, int16)'], target='cuda')\n",
    "def zero_suppress(waveform_value, threshold):\n",
    "    ### Replace this implementation with yours\n",
    "    result = waveform_value\n",
    "    ###\n",
    "    return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x11e88a400>]"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3Xl8FPX9P/DXOzcJhIQkQIBAuG+5wiWinIJixWptwdbiSevxtWr9VazWSq0ttfWoVetRb6tixQMVUUQEuQ3KfYZwhTMhhiuBXJ/fH5ndnd3sZq+Z+ex85v18PHgwOzu7894j7535zOfz/pAQAowxxtQVJzsAxhhj5uJEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjiuNEzxhjikuQHQAAZGdni/z8fNlhMMaYraxbt65MCJETbLuYSPT5+fkoLCyUHQZjjNkKEe0LZTtuumGMMcVxomeMMcVxomeMMcVxomeMMcVxomeMMcVxomeMMcVxomeMMcVxovfjm12l2H+8UnYYjDFmCE70OkIIPP3VLlz70lpc+PclssNhzO3VFXuw48gp2WEwm4qJkbGxoL5e4EdPL8eWQydlh6KUV1fsQVVNPW4Z01V2KLb20MdbAQDvzByBEV2yJEdjbyeqajBg9hcAgKJHLkFCvPrHu+q/whCdPFvTKMn3eOAznK2pkxSR/eXP+hQPfbwVf1u4XXYoypj2wmr85p3vZYdhawP/9IV7+cAPVRIjsQ4neg0RNVpXXVuPYyfPSYjG/u542zsZ3fLmOkmR2N/u0tNetz9af0hSJPZXfqYaQnhuL9p6RF4wFuJEr4mPa5zoAUBA+F3PmjZ/g3cy+myzM/6gzDD+saWN1lVUVkuIxN4qKqsx+OFFXuv+smA76uvV/xvnRI+GI/dJTyyTHYYyTp+r9bv+zdUhFdpjOoGaDic9yd/XcE19ZoXf9f/4YofFkViPEz2AZ78uwsEK/211Qv0fe0PV1wv0++Pnfu974MPNFkdjf48v2ul3/VFuUgzbvgBdpp/9erfFkViPEz2AisqagPeN+cfX1gWiAO61ZKwXlhXLDoEpgBM9gFdX7pUdgjJq6+ubvH/XUe4LHqqtQX40z9VyjzAWGk70zFA/fnZlk/cv2nbUokjs77UgByA9H1hoTSAO8MqKPbJDMBUn+hBsOXRCdgjKIPjv3cQam1t4QHYIyli6s7TJ+2drA9JUxYk+BFOeWi47BGXw4Ckmw4yX18oOQSrHJ/rKav9dAVn4BHdRslzRsdPBN2KO5/hE/8qKvbJDUMYb3E/eMNW1TV/Udik7zd0sjbKxpEJ2CKYJmuiJKI+IlhDRNiLaQkS/0da3IqJFRLRL+z9TW09E9BQRFRHRRiIabPaLiMaCTYdlh6CMBz/aEtJ2Rce4500w987bGNJ2y4K0PbPQXfPiGtkhmCaUI/paAL8VQvQGMALAbUTUB8AsAIuFEN0BLNZuA8AlALpr/2YC+LfhURso1H7f3JXNOFf9e5XsEGLeB98fDGm7QKU7WPgCjehWQdBEL4Q4LIT4Tls+BWAbgPYApgJ4TdvsNQBXaMtTAbwuGqwGkEFEuYZHbrFb3/xOdggxLZx6ISeqAg9QY+Hh6qrBzVtXIjsE6cJqoyeifACDAKwB0EYIcRho+DEA0FrbrD0Afb+wEm2d73PNJKJCIiosLZVz+nkmwC/41j9NarRu8fZjZodja29/u192CMoI56L2i9+o3f/bCL/93wbZIUgXcqInouYA5gG4UwjRVHuHv3PJRt9cIcQLQogCIURBTk5OqGEY6rEvGtcR+eZ3Y5GalIDmyTwnSzhKTzW+KLj29+Ox6r5xEqKxt00HG4/b6JDZDF2y0yREo6abLujsd72qlSxDSvRElIiGJP9fIcT72uqjriYZ7X/XIW8JgDzdwzsAiMkC2st2NT6TyGuVCgBYdPeFVoejnNbpKcht2Qxt01Ma3cfXPAKr9ZNsXr9hGL66Z4zf7blba/geuKyP3/XFZWp2Vw2l1w0BeAnANiHE47q75gOYoS3PAPCRbv0vtd43IwCccDXxxJqm+iDntmxmYST2V34mcH30j24f1Wjd4YqzZoZja/5OibvkNAcArH9wYqP7th7mQnLh+MfVAwAA7/16ZKP7nl2iZiXLUI7oRwG4FsA4Ilqv/bsUwBwAE4loF4CJ2m0AWACgGEARgBcB3Gp82NHbd/xMo3U3BjidY8G9vsq7D/2Xd1/kXm7j54i+pi60fuJO5Dvb2e8v7eVezkhNarT92j3lpsekkp8M6QAAKMhv1ei+90Ps7WQ3ofS6WS6EICHEeUKIgdq/BUKI40KI8UKI7tr/5dr2QghxmxCiqxCivxCi0PyXEb4SP3NFntehpdftB31O7/gUOXQdtSYwl7E9va/D+Ls+whpc+5J3f+5hnb0nA19+71iv2/y1DCzYyOF5tzQ+qleR40fG6g31+YXPau599PTEl7usDMfWEnz6d/vmIn8XHFmDU2e9e4N18vnR7JDpfZtHJAfmO5H6sM7ef+OtWzQ+21SRYxO9bzvohT1y0C7Du11+bK/WXrefWsyJPhQL7xyNOJ9E/8AU77OjQDN6scYyUhObvH9PWeNmSNbAd0DkmzcO97qd25ITvdKW7Srzuv3CtUMabZOe0vQfGPOvV9v0Ruu6tW4uIRI1+LbZA0B+VqqfLVkwSQneKS8hPg5T+tt+PGdQjk30zy31vrqekhjvd7ubR/MFWmadH3x6L825sr/f7X42tKMV4TiCqxeOyhyb6EOV5jNwiqfCa0zlGiFWe9dnspHmKf4H7iXGc40boyT7HOU31VXYrjjRB+E7I9L+cv8zyTtZoFISvgbkZZgcif35trf3zm3cDAYAI7pk+V3PAnv/1vP9rve9nvT79zdZEY6lONEH0To92et2RSUX5PL1r69Cu0jtO+xc1eHm0fCtQ981x/+1jb7tvH8ADvABSCPHTnkPyuse4nWihVuOmBGOVJzog5jQu43XbZ6RqrE3V4dW0OxSn4teCzbH5IBpqUIdsON7gfayf/F0l77+6DM/QgsHd67gRA9gZBOnwWnJ3hdp4+P4LYuUb+107ktvHC793Nhnm0M/Mm+foXbJE85aAF69YWjA+1KTvC+GicaFOJnOn6/oF/K2r6/kgT5NmTY0L/hGzBAL7hgtOwRTcaIHkJzgv2ulP/d/sNnESOzv58ND7/ZXxZNmNGlwx0zZIThGyyCD0uzOkYmea9aYx9/gHj0uHBe6YD1rpg/jvvSheu4X4U1drVpTmCMTPdcGMU5FZXh9jgs68VFqqIL8ZuKvAQZTscYm9wtv9OuA2V+YFIkcjkz0D/pcjQ/m2hGdTIrE/sLtbjqxT5vgGzEAjft3MxYpRyZ6vVD61v7xR/5no2HAquLjYW2fEO/9lTt+uvEUhKxBuD1BuEmSBeL4RD/lvOCndL7JaeXusgBbOk9VdXQXVHnSdY/v9v8Q1eMf/mSbQZEw1Tg+0XeKoArgur3R/UGq5E+fbI3q8TuOcO0glxW6iqr/nDYw7Me/vGKPkeHYWiRnN6/fMMyESGKD4xO970xIoTh5Vq0r8jK9tJyTk8sPuusd43uHdi3jN+O7mxWOrc2aF369mj7t/NcVUoHjE72/2unBnK3h+U79uahHTvCNAAzqyMXN/NEfkTcLUDbb101cRtuvuT5VQEOR3Tw5+EY25fhE71uGOBSRfImcINTufs9cE16fZifyLRcRiJPrt4QqkuZZ1Tg+0Yfqrgk93Mu+FQZZA9+pGKPdjjEjtG6h7pF6qDjRh2icz/yxjDF7uO58bt5yXKIvOhZZL4+2DplEmDHVDO4U2TWhTzeqU0bbcYl+wuPLInpcDp/+NVKnmzhkQm8+44nGuVou8GaW3JaRNRXe9tZ3Bkcij+MSPTPOf74p1t2KfLg+zzQFvLpiryHPU8YjjZkfjk703UKcWoz599fPtruXgxXgaspxBSdjDteGkgpDnudwxdngG7GAWqUlyQ7BFI5O9Jf2axvxY7lGi7dOEQw8c+HJXIAvt3lKQfjOBxuOF73OspxpTZj1l/RW3DvOwEhih6MT/R1RjCpcsqPUwEjs75YxXcPavp3u4vZL3/DoWH2X3fun9I74eeZvOGREOLb25prQ5jD2p1lS6JMQ2YmjE71vsbJwbD100sBI7C8rzFGF+oFq760rMTocWwt3ZqlHfhz69I1O8LHux+6T/7tAYiSxw9GJPlz6ujhcQCo6D13e173MbfTeUkIsf+ByDc80FVC/9i1lhxATHJXo66Ls3cHT4BlnVLds2SEoI9j0jYw5KtFHewE11PojTnCKK3gyZhuOSvTRap/JNVpc9pZVyg6BMRaioImeiF4momNEtFm37iEiOkhE67V/l+ruu4+IiohoBxFNMivwSOjryLdsFn7Vv7E9efSni9FdIp08Dd6uozz5CjNXKEf0rwKY7Gf9E0KIgdq/BQBARH0ATAPQV3vMs0QUM/2V9OUPwu0OyLx9ZfAUgNFeP7GzQyeMHeRUU8fVVY2yRJGpLoMmeiHEMgDlIT7fVADvCCHOCSH2ACgCEJPzc/GF1eg8+eUu9/Ir1w+N+vnqHHxEb7Q4vjhrmNJTagyMjKaN/nYi2qg17bg6/rYHoJ+Vo0RbF3MSo+hDz7zlGlDZ08lH9HoPRDhYqmebFu7lQxVVRoXjSPoJdOoVOQCJNNv9G0BXAAMBHAbwmLbe36GE33eKiGYSUSERFZaW8ihTO0tLCn+WLsD7B8LJiX5P6Wn38i9GdIroOcbq5ksY/eiSqGOyq6Jjp4NvFMRPhnRwL6/cHXk5hVgSUaIXQhwVQtQJIeoBvAhP80wJgDzdph0A+B2TLYR4QQhRIIQoyMkJba7RWMOjYxvkRVjnZnAnzwjQegc3K//98x3u5eSEyI69uF5Qg7nfesof/O2q0Ka29JWg60atSkmJiL5VRJSru/ljAK4eOfMBTCOiZCLqDKA7gLXRhRi7jKo46FT6tuTnlu2WGIlcZ6o9tegjHfx0flcegAYAL+rqJk0dGFmrsYoD0ELpXvk2gFUAehJRCRHdCOBRItpERBsBjAVwFwAIIbYAeBfAVgALAdwmhFB2RgVV2u9k6Zyd5l5evO2oxEjs76Ie9jwrNlO4pSRUFrRxVQgx3c/ql5rY/hEAj0QTVCwbmJeB9QcajuRPna2VHI29Zelqf+88Gn3bKmPMP+56EqbrR+W7l+foJt5g4ftpQV7wjRhjUXNMoq+sNuboO9IeJqwxVWt/MxZrHJPo9QMfbhsb+ahYBa/TMMaaoMK4BMck+t26vsqpfFTOGAvR4ROc6G3jhlcL3csjumRF/Dz9eSIDr8FN/H4w9dn/NN4xiV5vYF5GxI9tnR79cH+7e2GZZwLquyZGPu8uM09VtbK9mi1xk64WVtEx+1cXdWSi5wlEovPMkiL3cs+26RIjsT99Qp47c0RUzzU03zPS+P3vnTcPr77UdbR/4/dd6qk5dO+8TVE9VyxwZKJn0Tl9ztODKSXCIfsu+no3TqxJv7/cM4FLtBPbTOrb1r1cW+e893LZrjL38vldI2+eBdQ7GOREz6KS1Tw5qsdP101s7cTCZvoiqs2To+skME5X2OzfXzuvpMSirUfcy7++iOeb0ONEHyUnHoUaSf/2bTts/7bQcOkncMlITWpiy+D0R6FHTho7mYndGFE6WyWc6KP0xVau0RKNywe2cy//6OnlEiOR4y8LjBtd3cbhHQWW7vSUO4/2R1M1nOgjMKqbp/3vRGVNE1uyYPSFzVh0nF7E60C5p7+7E5sBm8KJPgIPXtbXvbxyd1kTWzLGZMhuzkf0epzoI9Alx3MU+uF6NSYmYEwlKtaUjwYn+gjwfLOMMTtxRMbiKf8YY07miET/Q2W1e5ln4jHO8M6tZIfAGAuBIxL9t3vL3cv/nDZQYiT2d7bGM2R/9tS+TWzJmL09MKV38I1swhGJ/tWVe93L3L82OjuPegY18SQsxumUlSo7BObjhlGdg29kE45I9DW19bJDUMblT69wL3PHhujoR1XPGJkvLxDF/Ga8MRVV4xSqd+OIRH+GS7aaIr1ZouwQbO2kbnJ5o6ZVzEj1fCb6ZjYnMePsaN46e1cDdUSiN8O0oTyxdXqKMYleP9K4ps45Z1+lpzz1aPq2M6bc89s3e0odf7b5sCHPaTdmVJ787f82GP6cVuJEH6GJfdrIDkEZF3Tz9IR6a81+iZFY67J/eWr7tMuIrkSxS+9czw/GXXPtnZzCoT97Gd45uhLFKuJEH6HWLZxdQMpIl52X615esMk5R6FnazxnLzwILzonqjw1p1pyk2Ij/O2KUJzuneNSxdHRt0+v2VPexJbq4uQUnam6TgJGXe9QCSf6CLVK83TTLC47IzES+0vlP0wWJafX3w+GE32Eclt62lR/+twqiZHYXyr3x2fMVI5K9GN7mlP+4PiZ6uAbMcaYJI5K9KpN+MsYY6FQPtEf0/VVjuOhnIZJT+HmFsbsQvlEv7LouHvZqL7KDGjt8PlJmTP8fHhH2SEYQvlE/+U2z+TdMy/sIjEStSRwMxhzgKzmye7lwr327fqrfKJfXcxH9Ga4flS+7BCU0SfXmPIHzHgX9ch2L+88elpiJNFRPtGXneYeMUY5UF7pXr56iLG1fvJaOfdHeHzv1oY+35s3Dncv82C+6Azp5Jlc5/jpcxIjiU7QRE9ELxPRMSLarFvXiogWEdEu7f9MbT0R0VNEVEREG4losJnBM2u9uXqfe9no69qzL3fWJCb19Z4EPLlfW0OfOzHe8+F8u/cHQ5871pnZSeCxRTtNe26zhXJE/yqAyT7rZgFYLIToDmCxdhsALgHQXfs3E8C/jQkz9un/cFW1StcMRgZn+kyHTQizvKjMvawfZW2EAXkZ7uWNJRWGPncsKtMdaXfI5Alc/Ama6IUQywD4XoWYCuA1bfk1AFfo1r8uGqwGkEFEuVDU9GGeK/J1DjhF3lhywrTndlrX13cLD7iXUxONPQpNSfSUlNDPrqaq5bs8P5o8U5d/kbbRtxFCHAYA7X9XI2N7AAd025Vo6xohoplEVEhEhaWlpRGGIdeU/p7fsHoHJHoz6fO8E9qV9cXbWqaaV9CsorIm+EY2t/2IZ3rLeyf3khhJ7DL6Yqy/wzK/f7VCiBeEEAVCiIKcHHNKE5hNP9K23jnzZZhCf0T/ns1n8wlF6SlrLuxVO2AaTX2J4vzsNImRxK5IE/1RV5OM9v8xbX0JAH13jA4ADkUeXmxLSfS8fedqnTltm1E66k659UdoLDrVDpixq7jUvt0erRJpop8PYIa2PAPAR7r1v9R634wAcMLVxKOi5ARPW+iCTUckRmJ/zXTtyi8t3yMxEmY3Tp3DIByhdK98G8AqAD2JqISIbgQwB8BEItoFYKJ2GwAWACgGUATgRQC3mhJ1jOid28K9vHALJ/po8AxLjJkn6OV+IcT0AHeN97OtAHBbtEHZhb6L4bKd9rygzBhTn9KHUfreG7Mu4avxjDFnUjrR19R5Ev2ILjwzPGMsfL++qKvsEKKmdKKv1fV57JDp3FoqRruwhz27wzIWiebJ9p/TWOlEf6jCM+kIX+wzziDdEHvGVDd1oN8xn7aidPZ7bulu93LLZuaNPnQCfS2fiX3aSIxELTktkoNvxKRqnuzps3K2xp7jZZRO9Gv2HA++EQtJcdkZ97LT6tKY6doRnUx53vsv7W3K88aygSadaTZL8jTdHD9jz7LnSid6LktgnD26RC/8V7WI2h3jupnyvLGsV9sWwTeKQP8OLU153liW3dycCqj6InFnztWasg+zKZ3os0z64J2oTveraVbNMf1F3opKex45hWLnUU+Jh26tm5uyD/1UjwcrqkzZR6yZ3M/8QrkXP7HM9H2YQelEb2ZZXZe5M0eYvo9YcOSE58K2WYlef+pdbtNT5FDoJ3AxaxqDfu09R/R2PQoNxW5dnZu0JPv3jjGL0oneCkkJzngLH/p4q3u5bztz5jhN0PWMmlt4oIkt7e31VZ5E3zXHnGqL+uYGlStY/uqNde7lgvxWTWzpbM7IUiZqkeLpzeOEOuoAEBdn/sVY/RmEyoyeqcufk2fVrUm/XzePcXIip7NA+J2Jkn5GG/1IXBadnUe59KxRrnlxjewQTKM/W0nisTIB8TsTJX1XQyfU/rbK8dPWTMzB1MGJPjB+Z6Kkn2XqznfWS4xELccsmoGJqcOKJkW7ckSit2rC4JIfKoNvxBhjFnNEok9JsKbbFTfdMMZikSMS/ZD8TEv2U8sXYxljMcgRif7hqf0s2Y++q5eqzOr3zVgs0xc2syNlE72+T3s8X6QxzJyrzpMdgjKsunbEovcLk4rPWUXZRP/V9mOyQ1ASJyfj3DvZ3Okts5tzCWSjXNAtW3YIUVE20XP3PHMQ+OwoGrW6C/Y9Tapc6cLzJBsnNyNFdghRUTbR/+ebYtkhKKOq2jPZAk+UEZ0jJz2lHcz+yZw6sJ3Je4gdrdLMrVTbqZXnTNaOReKUTfS7S88E34iFRJ+czDakkzU9pGTRj6Q2+wKfvlTxiSp1690A3tU6zXbAhuNllE30VvpZQZ7sEEz1/f4fLNuXWbMExQp9om+dbm5zgL5gmsr1/QHgpwUdTH1+/ee2bGepqfsyAyd6A1wzvKPsEEy148ip4BsZZFhntUvNFu4rl7LfWrMK38cIs68d6csrbDp40tR9mYETvQHys9TuW/78Muuud+ing1Ox7PO+43JO+z/47qCU/VrFynpmH284ZN3ODMKJ3gBOmXzECknxnnIVVTV1TWxpT09/VSRlvyq+l/pOAmN7tZYYSezjDGWAegWPPGXRD25buPmIxEjMISvhrtp9XMp+zbRDN/duskX1rOyKE70B9Gn+qIU9VFSUmeaZseveeRslRqKWrYft164cjH7CetY0TvQG0Lclb7fwwqWK2up6ovCMXawpj32xU3YItqF8ojdrIms9/UTMM15ea/r+VGbFHKpMDSsVbI4yi/KJ/uoh5vavBYBEnsKMMRbDlM9Q7TKayQ6BMcakimoMNhHtBXAKQB2AWiFEARG1AjAXQD6AvQB+KoSwbmilj5bNEoNvxBhjCjPiiH6sEGKgEKJAuz0LwGIhRHcAi7Xb0gxQfEi9lW4Z01V2CIyxCJjRdDMVwGva8msArjBhHyHTXyhl0bm4TxvZISgjJdGaVtMBHawr9qW6uyf2kB1CxKL9tgkAXxDROiKaqa1rI4Q4DADa/5YPWdtdetrqXTqC6vVSrDSxT1tL9vPgj/pash8nsHOni2jrpI4SQhwiotYAFhHR9lAfqP0wzASAjh2NLQo2/rGlhj4fa9DfwlKwqkuMt6YbKU+jaRz9Z/bDmWpkmlwD30hR/UQJIQ5p/x8D8AGAYQCOElEuAGj/+53TTwjxghCiQAhRkJOTE00YzESnznrqmFvVDPa3q/pbsh+ZsixKEk7I891bN7dkP78cme9ePnrKXiPgI070RJRGRC1cywAuBrAZwHwAM7TNZgD4KNogmTxvrdlv+T5Tk8ydkCMW3DG+uyX7yWjm+UHRT2OokvaZ1nSh1hcvtFvtoGj+otoA+EAbyZgA4C0hxEIi+hbAu0R0I4D9AK6OPkwmyzNL5FRbVNFJ3dlRixRruv121E3mXlMnoGLtr37trG9SrLHZj2bEiV4IUQxggJ/1xwGMjyYoOxrVLQsriuz1Kx+Kk2ftNz9mrCqWPL2lgDoX00/r5m392VDrZ3h7bmkxZl5on+7G9r2MHGOSbHxFPtboyz6r1Nygr58uQ02tOole/72QMfq9/Iy9pmbk7GSQ0d35grJR9Bd9l9pwfs5AZJewfvZrdZrh9FNAcM+i4DjRG+QXIzrJDkEZE3p7Bmadq1XniP7Oueul7t/KKSHNdrCiSnYItsKJ3iA8naBx9Edox212isysUceD98KidHZqL6lyZT1/CQ1zTsG5Tln0pj6zQnYItqJ0or/Kglr0/hTuk1asUzlnOdEzFjWlE72sazTvFh6Qs2MFLS8qkx0CY7andKIf3T1byn6/3qFOTxHZVheXyw6BMTf9nMZ2onSiT4iT8/LKTp+Tsl/GmLnOs2nZZ6UTvX7gDYvOhT14nABjwzq3kh1CRJRO9JzmjTPcpl9wxqO2jZRg08FZyn0Ddh495V7unJUmMRL7O3LCM5JzQAeektEo43tZOxfPJf2tmeTECfRjPOxUnkO5RD9//SH3sp0mBohFVbqujWTPA5mYZHUzWGqSgiUrJSnI95zZVlTVNLFlbFEu0W88eEJ2CMo4o6sQaHXTzf+b1NPS/VnJ+lHU6v5K92hjzaQjLq1bJLuXq21UnkO5RL9MoSJYsm09dNK9nGBxO6+sUc1m0Q/Z75Obbum+s5ure2b7K4tLBWc19yT6kh/sU29HuUQvU5v05OAb2YjMwUpDOmVK27cZ9BNVDMiz9nqHvrLq7tLTlu7bbIkSa0xtOWSf1gNO9AbKbq5Wop+/4VDwjUyi2vWVxdv8Tp1siY6tPLNMHbTRUWgo8iyaRtCf3Jb2GTzFid5Alw9oJzsEZTRP9kx+dsJGF70C2Vcub3aptrqEpNpF9ZwW8g6u3i0skbbvcHGiN9D04R1lh6AkFaqBxsprmPPZdtkhRE3frVFmSYKvtss7SwsXJ3oDpSV5jkLtdEU+1pVX2r8mfW2MJPotugvsdrX3uOfsyOpOAnbF75KB9GfFFVX2T06x4tONh2WHELUaGw2uiX2KtT9ZgBO9gfTtn1xmxziPL9opO4SoPbNkt+wQlME1rMLHid5ApMv0/GVkzBw3vPqt7BBshxO9SaqqeWYkxsxgp4FKsYITvUkeXbhDdgiMMQZA4USfkij3pS3cckTq/hlj5pjQu43sEMKmbKL/44/6yg5BGekpCcE3Yswh7FieQ9lE3621tVXtzLbl0AkcKK+Usu/WNp0nMxJf7ziGszXqXV+57LzcsB+zt+wMdhw5FXzDEBWXnsauo8Y9nyyT+9mvvr9SiV5fIVBf30MFU55ajtGPLgnrMW+t2Y9jp84G3zCIiX3kn6ruO34GH35/sNH6unqBZ78uwr8W73IPUhNC4JUVewKWTli1+zhWFx9vtH7roZNwdfkJAAAPUklEQVS47pVv8dD8Le51Gw5UYImNRkAGMqhj+EehY/7xNSY9uazR+rV7yrGyqAxCCLy6Yg9OVIZWomLcY0sx8YnGz2c3WSFUA12y4xjWH6iwIJrQ2Pqc/Pjpc3hs0U788Ud9kBAXhx4PfOa+Ly4Ginr89t0N6JSVijvGdw+4zX++KcafP92G60flezU3bT54Ap9vOYLW6SnI1R1RHyivRF6rVJyrrcPsj7fitxN7uEunHqyowvNLd+PBy/rg6Klz+P0HmzC3MAMf3TYKt7/1HT7ZeBj/vWk4RnXLdj/fiaoazPlsGx68rC8K95Wj6NhpXD+qs1eMN4/uYtRbErGL/v41AKBDZjM8unAHqmrqcH63LFRV1+H1VfsAAC8sK8am2ZOwZk85Zn+8Fd/tr8C/pg/Coq1HUXrqHK4Z3hEnqmow/cXVAIAxPXOQGB+Hh6f2wzNLityTun+84RCymyfjnkk9MfWZFQCAvXOmoLK6Fje9VogOmc0w58rzEBdH+PD7gyACpg5sb/2bEoZztZ6zlLfW7Mc1WrmOrYdO4tNNh3DPxT1RdroaQx/5EgvvHI1ebT2llM/W1OEPH27G/9aV4MmfDcSdc9cDAIblt8LaveVYvP0YOmQ2wwNT+iBNV6OosroWw/+yGPdc3BNf7/D+sRRC4G8Ld+DKwe3Ro00Lv99nAHht5V50ykrFmJ6tcbCiCrPmbTTl/QmHfmrGtXvKsWr3cSQnxuHdwgMoLj2DB6b0xp8/3QYAePqaQbjsvIYaWK+t3Iv87DRcJGH+ZRIx0N+7oKBAFBYWhv24e/63Ae+taygsdP2ofLyyYq/7vj1/vdSrX7tV8md92mjd3jlTAAAVldW4/OkVeP7aIeidm46zNXXo9YeFjbbvnJ2GPWWBi2B9cdeFuDiMI6MJvVvjS131xOSEOJzzKdHw24k98Jg2MOma4R3xlx/3d7+WXY9cgkQJQ839vZehyGmRjNJT5xqt3ztnSljP2bddesCSAVMHtkMcET7QzjJcn/Gd73yPD3WznL18XQEOnziL+z/YDAB45Mf98PPhnUKOwSivrNiD2R9v9VqX16oZDpSb01Vx40MXY+brhVhdXN7ovr1zpuCOt793V0fdPHsSfv3GukZlsbvkpKG4tOHv4MIeOSguPe3VtbJFcgI2zZ5kSvxNEUKg830LQt7+ioHtvL4TLqlJ8Vh7/wSvAn7hIqJ1QoiCoNvZOdH/7r0NASvIuf7wrBZpcoplO/98iYRZkdR8L9+/9XwMjqAZJVpfbDmCmW+ss3y/oZg2NA/vfHsg7MfdO7kXbhlj7cQjLkZ+N6PJVaEmelu30W8ssU/hfzuTkeQBNcs+J0kqwhXLXQIjSfJAw5mvCp5fan55DFsn+u0G9ghgsWd4F2vnqbVCxyy1OgnINKijtTN1meWvFpSONi3RE9FkItpBREVENMuMfSTEyb/g6muEgslJlvO7ZgffyGbSUxKl7Fd+A63xYqC/hW2YkuiJKB7AMwAuAdAHwHQi6mP0fmKlxrfeRT1ayw5BGfH8l2yYGDwmilos9KyzC7OO6IcBKBJCFAshqgG8A2CqSfuKKT8eFNvd7OwkztYNi7FFRg80s8XiGX2sMutPqT0A/RWWEm2d8lpwuQDDxPMfMmtCRqpaE8ibyaxE7+8v1KudhYhmElEhERWWlpaaFIb10qLoE8u88XvJmDHMSvQlAPJ0tzsA8BoxIIR4QQhRIIQoyMmxfqQYi32yLlwyphqzEv23ALoTUWciSgIwDcB8o3eS16qZ3/Wy+iozxlgsMuXcWAhRS0S3A/gcQDyAl4UQW4I8LGwpCfF+11852BGXAxhjLCSmHfoKIRYIIXoIIboKIR4xYx9tApTPzW3p/0ifhW/Olf1lh6CMfB4spYx+7dODbxSiJfeMMey5ArF1G0egOhddctQYGh0LruDuooaZeaGcuiwu7996vtT9G+nuiT2k7n/W5N6GPZcVpRxsnejzMv0fIeVnyU30H4TxB9U1Jw1PTR+EtgHOTvz1FZ4+LM/Plo3dMKozJvRujewA9bPn3z7K7/p/ThvoXk5J9N88pqqstCQsv3es+3Z6SoIhTYFXD+kQ8udmlsEdM3H/pcYkqCn9A09kMnfmCPRq2wKA95SeM0YaV7VTdsfbC7p7j9r+5chOKOiUGfbI+HByRTRs3X+tY1Yq3r55BJIT45AYF4f0ZgkgkPR6IvpJHkZ3z8ZtY7thQIcMrNxdhjbpKejbLh2d71uAgXkZ+PC2hmR7Sb+22He8EhMeX+p+7MguWXj5uqE4dbYGVTV1IBD2l1diVLcsvL3WM0zhmWsGY1yv1rj59UIsLypDckIcvvvDRHf3xHO1ddh++BSaa33826anoLK6DjktPHW/X71+qLvJq2fbFvjNO+vNe4PC8OXdF7nfk4F5GXhq2iD8UFntrhO/4cGLMeBPXwR8/KNXnYffzduIQR0z8Mp1QzHxiWUoPXUOn995IbrmpOH5ZcX4++c7cPvYbnh6SRH+dtV56JCZiqSEOFTX1mPt/ROQkhiPnBbJeH5pMQBg1X3jcLamHj9UVuPKZ1e695WaFI+aunosvnsM4uMJo+Z8hQem9MbVQ/LQMjU2ehBN7tcWjyzYhsd/OgA92rRAu4xmOHOuFrf+9ztsOti4SOC2P03GXxZswxur9+HhqX3xh4+24OGpfXHtyHx8OutTDOqYge/3eybY+L9x3TC8SxYW3DEaByuq0C6jGbr+vqGk7+yp/dCnXTrunbcJn95xAaY8tRwA8IfL+uDhTxpKKI/uno1vdpVhynm5+HTjYQDAh7eNwhXa5+1yfrcsU96fcMydOQI/e2E1/nplf0wf1tG9Xl/ZckBeBjYcqMDz1w5BWlICWqUlYd53JXhp+R4UPjAB2bra+2aydZniWHbR35dg6oB2uGtiD7+jEl3vu+99+i/JWzcPD1jvpb5egAioF56BRUII1IuG4e6hjoSs18pIxPmcOYz4y2L8YkRH3D4u8KQpVqmrFyA01DZxvS7Xurg4Qn29wNJdpbj+lW/dj7l7Yg/cNrYb4rX7XY9duPkIfv3mOmyZPcn9Q1hfL9zP43of3lqzHw9/shVbZk9CXBwFfG9dn1duyxSsnDUOQnjey7p6EZODvvSv08X1fVy09ShmvrEOzRLjceeE7vjVRQ3NTa7Xon+s6zHjH1+KyX3b4p6Le3p9Ri43vVaIhDjCc9cO8dq/673bPHsS+v3xczz3i8GY1LcthGj4rD9cfxB3v7sB2x+ejIS4ONTVCyTGk9d3XrZA72VNXcP75fueue7Xf0+i4Yh69Cpyffll1dO3u6cW78Lji3bid5N74tYx3Uzf323//Q6fbjqMl68rwLhesVsKOBa9sWovSk+dw90X95Qdim2Fmuht3XSjon9cPUC5+W6tdPPoLjhzrhY3+EyHaJaHr+iHTlmpXMwuAteOzJcdgmPwET1jjNmUI2aYYowxFhwnesYYUxwnesYYUxwnesYYUxwnesYYUxwnesYYUxwnesYYUxwnesYYU1xMDJgiolIA+yJ8eDaAMgPDsQN+zc7Ar9kZonnNnYQQQedijYlEHw0iKgxlZJhK+DU7A79mZ7DiNXPTDWOMKY4TPWOMKU6FRP+C7AAk4NfsDPyancH012z7NnrGGGNNU+GInjHGWBNsneiJaDIR7SCiIiKaJTuecBBRHhEtIaJtRLSFiH6jrW9FRIuIaJf2f6a2nojoKe21biSiwbrnmqFtv4uIZujWDyGiTdpjnqJQ5xc0GRHFE9H3RPSJdrszEa3R4p9LREna+mTtdpF2f77uOe7T1u8gokm69TH3nSCiDCJ6j4i2a5/3SNU/ZyK6S/tebyait4koRbXPmYheJqJjRLRZt870zzXQPprUMH+h/f4BiAewG0AXAEkANgDoIzuuMOLPBTBYW24BYCeAPgAeBTBLWz8LwN+05UsBfAaAAIwAsEZb3wpAsfZ/pracqd23FsBI7TGfAbhE9uvW4robwFsAPtFuvwtgmrb8HIBbtOVbATynLU8DMFdb7qN93skAOmvfg/hY/U4AeA3ATdpyEoAMlT9nAO0B7AHQTPf5Xqfa5wzgQgCDAWzWrTP9cw20jyZjlf1HEMWbPBLA57rb9wG4T3ZcUbyejwBMBLADQK62LhfADm35eQDTddvv0O6fDuB53frntXW5ALbr1nttJ/F1dgCwGMA4AJ9oX+IyAAm+nyuAzwGM1JYTtO3I97N2bReL3wkA6VrSI5/1yn7OaEj0B7TklaB9zpNU/JwB5MM70Zv+uQbaR1P/7Nx04/oyuZRo62xHO1UdBGANgDZCiMMAoP3vmow00Ottan2Jn/WyPQngdwDqtdtZACqEELXabX2c7tem3X9C2z7c90KmLgBKAbyiNVf9h4jSoPDnLIQ4COAfAPYDOIyGz20d1P6cXaz4XAPtIyA7J3p/7ZC260JERM0BzANwpxDiZFOb+lknIlgvDRFdBuCYEGKdfrWfTUWQ+2zzmtFwhDoYwL+FEIMAnEHD6XYgtn/NWpvxVDQ0t7QDkAbgEj+bqvQ5ByP1Ndo50ZcAyNPd7gDgkKRYIkJEiWhI8v8VQryvrT5KRLna/bkAjmnrA73eptZ38LNeplEALieivQDeQUPzzZMAMogoQdtGH6f7tWn3twRQjvDfC5lKAJQIIdZot99DQ+JX+XOeAGCPEKJUCFED4H0A50Ptz9nFis810D4CsnOi/xZAd+1KfhIaLuLMlxxTyLQr6C8B2CaEeFx313wArivvM9DQdu9a/0vt6v0IACe007bPAVxMRJnakdTFaGi/PAzgFBGN0Pb1S91zSSGEuE8I0UEIkY+Gz+srIcTPASwB8BNtM9/X7HovfqJtL7T107TeGp0BdEfDhauY+04IIY4AOEBEPbVV4wFshcKfMxqabEYQUaoWk+s1K/s561jxuQbaR2AyL9oYcCHkUjT0VtkN4H7Z8YQZ+wVoOBXbCGC99u9SNLRNLgawS/u/lbY9AXhGe62bABTonusGAEXav+t16wsAbNYe8zR8LghKfv1j4Ol10wUNf8BFAP4HIFlbn6LdLtLu76J7/P3a69oBXS+TWPxOABgIoFD7rD9EQ+8KpT9nALMBbNfiegMNPWeU+pwBvI2GaxA1aDgCv9GKzzXQPpr6xyNjGWNMcXZuumGMMRYCTvSMMaY4TvSMMaY4TvSMMaY4TvSMMaY4TvSMMaY4TvSMMaY4TvSMMaa4/w+UKVQKh+vmXQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10e674048>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# the noise on the baseline should disappear when zero_suppress is implemented\n",
    "plt.plot(zero_suppress(waveform, 15.0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "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.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}