pixelmap_gen.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Generating Pixel Maps from .edf data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 190,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "from sklearn import preprocessing\n",
    "\n",
    "from mne import Epochs, pick_types, find_events\n",
    "from mne.io import concatenate_raws, read_raw_edf\n",
    "\n",
    "import dhedfreader\n",
    "from PIL import Image"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Loading the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 222,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Extracting EDF parameters from data/SC4001E0-PSG.edf...\n",
      "EDF file detected\n",
      "Setting channel info structure...\n",
      "Creating raw.info structure...\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-222-915717e81f69>:2: DeprecationWarning: stim_channel will default to \"auto\" in 0.17 but change to False in 0.18, and will be removed in 0.19\n",
      "  raw = read_raw_edf('data/SC4001E0-PSG.edf')\n",
      "<ipython-input-222-915717e81f69>:2: RuntimeWarning: Interpolating stim channel. Events may jitter.\n",
      "  raw = read_raw_edf('data/SC4001E0-PSG.edf')\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Converting time column to int64...\n",
      "Converting time column to int64...\n",
      "Converting time column to int64...\n",
      "Converting time column to int64...\n"
     ]
    }
   ],
   "source": [
    "first = pd.DataFrame()\n",
    "raw = read_raw_edf('data/SC4001E0-PSG.edf')\n",
    "sampling_rate = raw.info['sfreq']\n",
    "\n",
    "raw_eeg_fpz = raw.to_data_frame(scaling_time=100.0)['EEG Fpz-Cz']\n",
    "raw_eeg_fpz = raw_eeg_fpz.to_frame()\n",
    "raw_eeg_fpz.set_index(np.arange(len(raw_eeg_fpz)))\n",
    "raw_eeg_fpz = raw_eeg_fpz.values.tolist()\n",
    "\n",
    "raw_eeg_pz = raw.to_data_frame(scaling_time=100.0)['EEG Pz-Oz']\n",
    "raw_eeg_pz = raw_eeg_pz.to_frame()\n",
    "raw_eeg_pz.set_index(np.arange(len(raw_eeg_pz)))\n",
    "raw_eeg_pz = raw_eeg_pz.values.tolist()\n",
    "\n",
    "raw_eog = raw.to_data_frame(scaling_time=100.0)['EOG horizontal']\n",
    "raw_eog = raw_eog.to_frame()\n",
    "raw_eog.set_index(np.arange(len(raw_eog)))\n",
    "raw_eog = raw_eog.values.tolist()\n",
    "\n",
    "raw_emg = raw.to_data_frame(scaling_time=100.0)['EMG submental']\n",
    "raw_emg = raw_emg.to_frame()\n",
    "raw_emg.set_index(np.arange(len(raw_emg)))\n",
    "raw_emg = raw_emg.values.tolist()\n",
    "\n",
    "first['EEG Fpz-Cz'] = raw_eeg_fpz\n",
    "first['EEG Pz-Oz'] = raw_eeg_pz\n",
    "first['EOG horizontal'] = raw_eog\n",
    "first['EMG submental'] = raw_emg"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 192,
   "metadata": {},
   "outputs": [],
   "source": [
    "def epoch_gen(data, epoch_length, fs):\n",
    "    print(\"Generating Epochs\")\n",
    "    eeg_data = data\n",
    "    L = epoch_length * fs # signal length\n",
    "    epochs = np.reshape(eeg_data, (-1, L))\n",
    "    return epochs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 193,
   "metadata": {},
   "outputs": [],
   "source": [
    "def extract_features_from_epochs(epochs, fs, epoch_length):\n",
    "    print(\"Extracting Features\")\n",
    "    N = epochs.shape[0]\n",
    "    L = epoch_length * fs # signal length\n",
    "    f = np.linspace(0, L-1, L) * fs/L\n",
    "    delta1, delta2, theta1, theta2, alpha1, alpha2, beta1, beta2 = 0, 4, 4, 8, 8, 13, 13, 30\n",
    "    all_indices = np.where((f <= beta2))    \n",
    "    delta_indices = np.where((f >= delta1) & (f <= delta2))\n",
    "    theta_indices = np.where((f >= theta1) & (f <= theta2))\n",
    "    alpha_indices = np.where((f >= alpha1) & (f <= alpha2))\n",
    "    beta_indices = np.where((f >= beta1) & (f <= beta2))\n",
    "    nr_features = 6 # number of features to be calculated\n",
    "    features = np.zeros((N, nr_features))\n",
    "    # calculation of delta, theta, alpha and beta band power ratios\n",
    "    for index in range(N):    \n",
    "        epoch = epochs[index, :]\n",
    "        Y = abs(np.fft.fft(epoch))\n",
    "        mean_total_power = np.mean(Y[all_indices])\n",
    "        features[index,:] = (mean_total_power, np.mean(f[all_indices] * Y[all_indices]) / mean_total_power, \n",
    "        np.mean(Y[delta_indices]) / mean_total_power, np.mean(Y[theta_indices]) / mean_total_power, \n",
    "        np.mean(Y[alpha_indices]) / mean_total_power, np.mean(Y[beta_indices]) / mean_total_power)\n",
    "    return preprocessing.scale(features)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 194,
   "metadata": {},
   "outputs": [],
   "source": [
    "def rescale_linear(array, new_min=0, new_max=255):\n",
    "    \"\"\"Rescale an arrary linearly.\"\"\"\n",
    "    minimum, maximum = np.min(array), np.max(array)\n",
    "    m = (new_max - new_min) / (maximum - minimum)\n",
    "    b = new_min - m * minimum\n",
    "    return m * array + b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 195,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating Epochs\n",
      "Extracting Features\n"
     ]
    }
   ],
   "source": [
    "eeg_data_1 = epoch_gen(raw_eeg_fpz, 30, 100)\n",
    "eeg_features_1 = extract_features_from_epochs(eeg_data_1, 30, 100)\n",
    "eeg_features_1 = rescale_linear(eeg_features_1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 196,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating Epochs\n",
      "Extracting Features\n"
     ]
    }
   ],
   "source": [
    "eeg_data_2 = epoch_gen(raw_eeg_pz, 30, 100)\n",
    "eeg_features_2 = extract_features_from_epochs(eeg_data_2, 30, 100)\n",
    "eeg_features_2 = rescale_linear(eeg_features_2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 207,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating Epochs\n",
      "Extracting Features\n"
     ]
    }
   ],
   "source": [
    "eog_data = epoch_gen(raw_eog, 30, 100)\n",
    "eog_features = extract_features_from_epochs(eog_data, 30, 100)\n",
    "eog_features = rescale_linear(eog_features)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 208,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Generating Epochs\n",
      "Extracting Features\n"
     ]
    }
   ],
   "source": [
    "emg_data = epoch_gen(raw_emg, 30, 100)\n",
    "emg_features = extract_features_from_epochs(emg_data, 30, 100)\n",
    "emg_features = rescale_linear(emg_features)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 210,
   "metadata": {},
   "outputs": [],
   "source": [
    "img = Image.new('RGB', (6,4), \"white\") # create a new black image\n",
    "pixelMap = img.load() # create the pixel map\n",
    "\n",
    "for j in range(img.size[0]):\n",
    "    pixelMap[j, 0] = (int(eeg_features_1[0][j]), 0, 0)\n",
    "\n",
    "for j in range(img.size[0]):\n",
    "    pixelMap[j, 1] = (int(eeg_features_2[0][j]), 0, 0)\n",
    "\n",
    "for j in range(img.size[0]):\n",
    "    pixelMap[j, 2] = (0, int(eog_features[0][j]), 0)\n",
    "\n",
    "for j in range(img.size[0]):\n",
    "    pixelMap[j, 3] = (0, 0, int(emg_features[0][j]))\n",
    "    \n",
    "\n",
    "img.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# EMG Power Spectrum Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 244,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x1200c9e48>]"
      ]
     },
     "execution_count": 244,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEKCAYAAAAB0GKPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzsvXeYG9d57/990XfRtmELm9hJkVSjqF4tyxIlFyWKnVi2E9fYuXHNdZzYcfJzEt9Up1wnceL4xr6WE8WOY19ZsiRbktWoQhVKYu9ctl1uAbahLDrO74+ZgwWxKGcGAwwAns/z4CF3FuVgMZj3vO37EmMMEolEIpFUw2L2AiQSiUTSGkiDIZFIJBIhpMGQSCQSiRDSYEgkEolECGkwJBKJRCKENBgSiUQiEUIaDIlEIpEIIQ2GRCKRSISQBkMikUgkQtjMXoCR9PX1sZUrV5q9DIlEImkZXn/99RBjLCBy37YyGCtXrsSuXbvMXoZEIpG0DER0WvS+MiQlkUgkEiGkwZBIJBKJENJgSCQSiUQIaTAkEolEIoQ0GBKJRCIRQhoMiUQikQghDYZEIpFIhJAGQyJpM944M4O/e+IIXjweMnspkjajrRr3JBIJ8KkH3sC5uQScO4Zx+KvbQURmL0nSJkgPQyJpI87NxnFuLgEASGZyOHAubPKKJO2ENBgSSRtxMhQDAHzxro0AgBPBqJnLkbQZ0mBIJG3E88eUvMVtG/sBAM8eCZq5HEmbIQ2GRNJGTEWTAIB1/R44bBaE42mTVyRpJ6TBkEjaiJdOTGHrii4QEW5ZH8CLJ2SllMQ4pMGQSNoExhgmwgk4bMrX2u2wIpHOIZbMmLwySbsgDYZE0ibEUllkcgy3blDyF7dsUGbijIcTZi5L0kZIgyGRtAmhiJK/6PM4AQADPhcAYHxOGgyJMUiDIZG0CaEoNxgOAMCgNBgSg5EGQyJpExYMhuJhDPoVg/Hy8JRpa9LDfCqDz/9wDz7z/TcxKcNpTYWUBpFI2oRgNAUACHgVg9HpUL7e86msaWvSw/0vncaP3xgBAHhcNvz5L19i8ookHOlhSCRtAs9h9Lgd+WNXrezGVCxp1pJ08R8vn87//ydvjpq4Ekkx0mBIJG1CKJpEd6cdduvC17q704HpWMrEVWljfC6B0dk43nHpED5xy2rMp7I4Oz1v9rIkKtJgSCRtQiiazOcvOL2e1jIYwyFF++rerUtx+8UDAIDHD4ybuSRJATKHIZG0CaFoapHB6HE7MDOfRi7HYLE0v8z5C6oW1tqAF8u6OwAAwUhrhdTaGelhSCRtQiiaRJ+32GA4kc0xhBOtoSnFK72W93TAYiGs6nPj6cOTJq9KwpEGQyJpE0KRZL4Hg9OrJsCnWiQs9eLxKVyhamEBgNdlw5wUUGwapMGQSNqAeCqLWCpbMiQFoGXyGMFIEo6CpP1bNw5gMpJEIt1apcHtijQYEkkbwEM5gTIGYyra/AYjlswglc3ltbCAhZ6SN07PmLUsSQHSYEgkbUCQd3l7zw9JtZKHUSxtAgDXr+kFAIzMxk1Zk+R8pMGQSNqAYuFBzoLBaP5Ko5DqBRUm7peqlVJPH5KJ72ZAGgyJpA3IX2yLDIbLboXbYcV0rPkTx3kPw73wHuxWC7xOG2IpOdOjGZAGQyJpA/jFtreoSgoAejyOFvEwSofVblzXhzGpuNsUSIMhkbQBoWgSPpcNTpt10e963M6WKKsNRZQ19rrP95KG/B04PhkFY8yMZUkKkAZDImkDSjXtcXrdrSEPEoom4e+w50fMcnwdiiDFcChmxrIkBdTVYBDRd4hokoj2Fxz7KhHtJaLdRPQEES0p89gPEtEx9fbBeq5TIml1QtHUebH/QnpaxGBMxRY3HgLADWv7AECKEDYB9fYwvgtge9GxrzHGLmWMXQ7gEQD/X/GDiKgHwFcAXAPgagBfIaLuOq9VImlZwvE0ujrtJX/X43ZgKpZq+pBOKLJYCwsAlnQplVIyj2E+dTUYjLEdAKaLjoULfnQDKHUW3wngScbYNGNsBsCTWGx4JBKJSjiehq+jvMFIZXKINfkgpVJquwAwoIbauDChxDxMyWEQ0Z8R0VkA70cJDwPAUgBnC34eUY+Veq6PE9EuItoVDAaNX6xEN08dmsCfPXpQzpRuAHPxNHyu8gYDAGaaPCwVjJYOSdmsFjhsFkSTsrTWbEwxGIyxLzPGlgN4AMCnanyubzHGtjHGtgUCAWMWKKmZdDaHj96/C//n+ZP4w5/sr/4AiW4yWcV74MnhYlpBgDCRziKSyJT0MADg5nUBTMj53qZjdpXUAwB+pcTxUQDLC35eph6TtAgPFozW/MWhCak4WkfCCWXnXc3DaOZeDG7MylV6DficmJRzMUyn4QaDiNYV/HgPgMMl7vY4gDuIqFtNdt+hHpO0CN945jgA4K/ffSkA4NkjUtqhXnBjXC7pzfsamlmAcCpaWtqEM+BzYTqWQjLT3HmYdqfeZbXfB7ATwAYiGiGijwL4SyLaT0R7oRiCz6r33UZE/wYAjLFpAF8F8Jp6+1P1mKRFmI6lcMPaXrzj0iEAwBMHJkxeUfsyO68YgrJVUp7aBQh/9PoIVn7xUXzygTd0P0clSgkPFjLocwEAJsPSyzCTuo5oZYzdV+Lwt8vcdxeAjxX8/B0A36nT0nTx8J5z+OazJ2CzEv70ni24fHmX2UtqSkZn44gkMrh1fT86HTZ1TGjz7m5bHe5h+MtUSbkdVjislpoMxtceVwIBj+4bw1djqXyYyyh4l3c5D6PfpxzfOzKH5T2dhr62RByzcxgtw9npeXzm+2/i3Fwc+0fn8LH7dyGTzZm9rKbkyLhSOc2/2Des7cNLJ6aavg+gVVkwGKUv4kSU78XQw/hcAhMFO/uHdhufTgxWCUmtCXgAAJMRmfg2E2kwBPmTnx4AAPzmTauxZakfoWgSP9s/bvKqmpNdp5RhN1svUjwwXqXDLwoSY6nmYQBK4ltvWe0Lx5X+h2+8bysA4PSU8R3XoWgSbocVHY7FWlgAsKy7Aw6rBeOyUspUpMEQIJPN4ReHJuF12vDJt6zFNz9wJQDg2y+cNHllzcmpKUXzh09/u3m9Iu1wMii1gOrB3Hx1g9Hr0e9hHDg3BwC4fVM/1vV7cK4Ow4ymoqmyFVKA4iX1+5wyh2Ey0mAI8JjqSfz+XRsBKFIF163uxe6zs4gkZLloMSMzcVy3uhdEBABY1aeEE14elnUL9SCSzMBpsywS7SukFj2pN8/Mwqsq4Q51dWBkxniDUa7Lu5ABn0v2YpiMNBgCPLz7HADgV7Yuyx9795XK/589IrvLC0lnc9g7MpfX/wGAlb1KLiMsjWtdiCQy8JbpweD4O+y6//4z86l8gUe/14mDY2GkDc7fhcp0eRcy4HPKkJTJSIMhwPPHgljZ23lefPWOzQMAgJ8fkHmMQnhz1UW9C5UsRITVfW4pEVInIok0vK7KBY8+lx3heBq5nLbCg1gyg9NT89i8xA8AWKEWMhjdRBeKptBbxcPo97owHIzJXgwTkQajCscmIkhmcrhuTd95x70uOwJepxREK4IbhUuW+s87Puh3yd0hlAa1T/7nG/i7J44Y9pzRZKaqwfB32JFj0DzqdFjNOw35lT6IS5Ypn+v4nHFhqUw2h5n50kq1hQTUHEc9QmISMaTBqMLLw1MAgHu3LtY+vGltH+bi6XyXqgSYVI3CgNpoxRn0ufD66RkzltRU/MNTx/Do3jH8w9PHcXrKmCKASCIDj7OKh6HqTHEZEVHGVMNwxQolJMUb6MbnjDvnp+dTYAwIVAlJbV3Rrb623HiYhTQYVXjjzCwcVguuXLF4HMftmwby95EoTOQNxvm7RR7Ou5CNazbHcP/O0/mfv/ncsCHPG02IeRjAQle4KPzz5IaCexoHx+a0LrMs1Zr2OIN+bqykwTALaTCqsOfsLPyddlgstOh3V6/qAQCckqMj85yZjsNmIXR3nr9bvFGdmjbRhGWR2RzDcDCKE8FovkS1Huw4phRIfOKW1ehxO/DaKWOqxpQcRuWkd5f6eWh9f6en5mEh5PML3PCcmDTunOfVW9W6x7nROiMn75mGNBgViCUzGA7FcMOa3pK/73U7YLUQnj0qhfU4O4en4O9YbGD71S/7RBN26v7NE0dw298+h7f+7XO4/e+fq1tHOs8H/OZNq3Hj2j4cn4xiXmNOoRSz8TS6KvRgAAs6U7MaVYN3Dk+hq1M5zwGlgGHjoFdzLqQSvHqr3AAoTofDCpfdghePy7yhWUiDUQHe0Xr1qtIGg4iwNuBBVGNcuJ3J5nJYVkLrp19NWE42WeL78HgY//LsCawOuPH+a1YgGEkaFioq5oFXTsNps6DX7ch7p7XullOZHOZT2YpNewDQpcqGaNX0yuYYlnd3nHdsaVeHoTPCeS9TtbAaAKzsdRte0isRRxqMCvCE95alvrL32XpRN/aMzMlpYCqhaAqblyz+e3HxuGbr1P3d/96j/HvHBvzenUpj5l/9/HBdLkrRRAarAx4QEbaoVWQvn5iq6TmrSZtz8h6GxpBUKJrEpiXnV7x119AEWIqIuuGqFlYDlOq7ZgxrXihIg1EBXga6YdBb9j5rAm4AwMFz4bL3uVBIq+WRgRLJS6fNiu5Oe1OFpIKRJPaPhjHoc+GuLYPwd9rx3quUuV1G5Rc44UQak5Ek7toyCADYqJ5TtfYzzMWVC7e/s3L832VXwjlaBlllcwzTsdSi6qVe1WAYFbrjBqNapRegVN8Fo0lkNfaTVOMXBydkibwA0mBU4OXhKawOuOG0lRZEA4BbNyhjYeuhr9NqKBeR8lPT+r0uvH66eSrKvrfzFADgT+7ZnJcx+R+3rgEA/HTPOUNf68CosqHgiVuX3Ypl3R01nzciwoOcrg5tAoTTsRRyJT7PHrcDSTUUZgSRRAZuhzWfJ6nEgN+FbI4ZVm2XyzH8cNdZfOx7u/CBb7+Cn7w5KlWVKyANRgVGZ+LosJc3FgAw5Ffiu4fHI41YUlMTVHfLpTwMTqiJymofVo3C7RcP5I9d1Kt4jG8YbNh4P8NVau4CAJb4O3BwrDbPlBsMn0D8v6vTrinpHSojOd7trn0gUyEiVV6cAd68Z9AG7QevncXv/Whv/ufP/ddu2S9UAWkwypBIZzEVS+HOzYMV7+d22mCzEPaONM/O2Sy4fHnAWzo8cuuGAGbnjQtl1MK52ThOT83j7ZcMLdrZvv+aFTg9HTN0naNqdzLvYwCUJO/RiWhNrxNNZvPPVQ1fhyIPIko5g9FruMHIwCOwfgB5jbLXT9V+UZ9PZfAHD+4DADz7u7fikU/fCAD4+lPHan7udkUajDK8cUY5IQu/4OW4fHmXoUnAViUUqTwEJ+B1Ip1liDRBgcD9O08BAN53zYpFv1vR04lEOmfoTvOpw5PwumxwFXisXNBPayK6kFiSx/+r79D9HXZNOYxyY1ON9jDCibRQSA0ALh5SCiq0vI9y/MuzJwAAX7hzA1b2ubFlqR8bBrx4/lhIquKWQRqMMnC9mssExrBuWerH4fEIEukLWxQtFK3csdurXnimouYbVz5j/PoSPTZv2dgPwNgwYzqbW7T5WKUWTNRSCMBLut3OyqFTQBEgjGgoAc93YHvr62FoMRhWC2HIAF2yTDaHf3z6OADg4zevzh/nOaxH947V9PztijQYZRibVU7IQtXVcgzIDlQASg6j02GFu0y1CzckZucxkpksToZiuGldXz7ZXciagAdWC+Enbxo3inQinMC2lT3nHVvQZarBYKgehttRPaSjx8Nw2CzwFn2eRnsYc/G0UA6GM+Bz1SwPwvNX9129AnbrwmXwHZcOAQB+cWiipudvV6TBKMOLJ0Lo8zgqVkhxtq1UdKYu9EqpUDSZVxQtRa9b+Z3ZelI/VwdiFSa7C7FaCL1uh2ElwMlMFqFoKm8gOHyjUUtvSiyZQafDWlK6phh/hx3RZEa4xyQYTSLgcS4yql6nDXYr6Z7gV0w4nhH2MAAlTPzSidpKYL/2uKIW/Du3rzvvuM1qwdp+D3YOyxn0pZAGowyz8yk4rGJ/Hp6I2z9qnCBbKxKMVJ6axmPhQZNDUg+pA7Heru4mS/GebctwdjqOjAENfDy8OVgUkuLNjHtH9RdMRJPVlWo53W7loizqZYSiqZJDjYiophnhhTDGEEmkq8qCFNLpsCHHoHva5dGJCMbmErhmVU9esqaQG9f2gTFIOf4SSINRhrHZBN62qfQOtBhe6nfiAp9ZHVJ3pOXg4nJmexivnZzG+gFPFeOm/G6PAdVv+0aUjUR/kffFvdez0/o9Uy0GY0GxVtBgVNgAdHfqnxFeSCyVRY6JVXlx+Iz4MZ1hqQdeVhSD/+gdm0r+nuewjBRYbBekwShBOJFGJJk5b8xoJWxWCzYOeg1z0VuVUDSJvjIltYDyd+rutJua9D47PY9IMoPbNlbeDNy0rk+9f+1hRn5eXLpscQHF1St7apoLH0tmyuaMiuEKwqIS55XmbPd6HJp1qUrB37tIlRdnqfq91BsCfmjPOfg77Hl5lmJW9ynFCA8amMNqF6TBKMGes8qusngIUCVW9HQaOoWs1VBkQdJVZxr0epymJr0f26dUv9ym7iLLsaxbKXZ40oDk52QkAbuV0F1C72nQ76ppoxFLZoUqpABtelK5HMNULFV2A9DdaYyeVDSvIyXuYSxVxRB3HNWexzh4LozZ+TTuvqR8f9Xynk50ddoNMYjthrDBICIfEfXwWz0XZTZ851JKRK8cQ36Xbhe5HeBeQ6WkN6DkMcz0MHj1y1UrFw/EKsRlt8LrtCFugPxFMFI6eQwoYbrpGv4ekWRGeHeuRbF2Np5GNsfKexhuhyGhxbAOg8HVFWI6+nm+tUPpvfjAtRdVvN/VK3vw9OH6jS14aPco3vI3z+Kurz/fUsUyVQ0GEX2CiMYB7AXwunrbVe+FmcnobAJEwErVNRVhwO9CJJExPT5vFsEqTXucXo8ToZg5f6NcjuG1UzPYvMRX8uJdzDWrewz5MgcjSQTKeKu9bgciyQySGX2GKZbMwCPqYWhIepfr8ub0uJ0IJ8QrrsqxIG0uHpIClP6ocxo9+myO4Se7z8HtsGLTUOXNIO89qcdArWyO4Qv/vRdjc3EcGgvj7548avhr1AsRD+N3AWxhjK1kjK1Sb6urPqqFeXl4CgGP87z67GoMeJULwpGJC1NTKpSXBalsMAIeZ74jvNHsP6ckn7lgZDWWdHXg8Hik5vLKyXByUcKbwyfZ6Q3vRBJpYVkNr9MGC4mFpKp17fPwlhapkVLwRkItfRgAsLTLhZc0SsPzwUvv2ba86obhJnVCJD9njOSNMzNIZXP5su6HdrdOrkTkingCwAXVkTYVTWoyFgBwxQoloXmhzhsWER4ElB11OJFBKtP4ITi/OKSEGN5+yRKh+/OqouEaR/BORhJlDcZC5Zh2g8EYQzgh3sNARPAJNu9V0wXjISQtneOlyEubazQYHqcN2RzT5AH85ytnAAAfvXFV1fuu6fcAAI7VYQP4lHoe/s7b1uPDN6xEOstaZsyzyFXxSwBeIqJ/JaJ/4Ld6L8xMJiNJ3H5x5aRoMfkB9Rdo7XawSgiDw3fUUyaEpZ5XZ2pfPFR+vkkhN6i7zFrCUqmMUgzQ7y0dkuJ9DnoS37FUFtkc09T05u+w50eiVqKazAsv5a3VYGiRZy+EFy2IqiswxvDzA+NY2tWB5SUmQhazNuCB3UoYq8P3+eXhKdithDUBD+6+ROkFerHGRsRGIWLW/xXA0wD2AWj72YjzqQwiiQwGBEQHC+l02OBz2fKSIhcaoWgSbocVHY7K8fS+Aj0pnrxsBHPxNN48M4ub1weE8hfAQvlmLZ8pN6S8Sa+YnrzMhnYDqudiKyoPEoomYbdS2efmOYdaSoIBpYTdbqWqYwSK4TL0O44Fccmy0uWxhTyvDkf61W3LhZ7fYiEM+TtwzuDvczbHsPvsLN6ihkW5AGWrXDdEPAw7Y+x/Msb+L2Psfn6r+8pMgo9/HCizI6zEgM91wWrQzMXT6Koy9Q1Y8DAaXVr79GHlc3nXZWLhKGChrLqWapm9aok2F+wrJu9x6QhJhfOzMMQNhs8laDAiSfS6S1d2AQUhqRqVh+fiivCgqBHn8ImFon+3f1Wro0qpE5dj0OfCEwfGNa2rGnwuytWrFNFLu9WCJX4Xfm7w69QLEYPxMyL6OBENaSmrJaLvENEkEe0vOPY1IjpMRHuJ6EEiKikFS0SniGgfEe0mooZWZPEcRLGMgwjdbgdyF6j+TDguJu8QyBuMxpbWcv2oSvX3xXDhvYTOCiZgIdR0cZmqHJ9Lvy6TllnYnK5Ou1Dcv1oTpi/vYdRuMLQYPA4RYf2AB08eqn6hTaSzePH4FC7q7axalFFIwOdEJseQM3AcLB/9e2mBVzTgdxlSvt0IRAzGfVDzGNBWVvtdANuLjj0JpeLqUgBH1ectx1sYY5czxrYJvJZh8DkYWpr2ONes6kEwYvy84VZA2SlWj3AuSJw31sN44VgIQ34XOgVUXQu5ZnVvTYUMwUgSROXnqnBdJj1/j2hS7ZLWkDDu7nRgWqAPQ9GRKn9xXUh61xaSiiYymnowChnwuXB2Ol61gIKP2xUNR3GuXdWDbI4Z6g2PTC/u8bppbR9GZ+O6+koaTVWDUVBKu0pLWS1jbAeA6aJjTzDG+F/lZQDLdK26jvDyxuU92uPrAz4Xcsx8+W4zEFUc7XRY4bJbGvo3GpuLI5bK4uZ1YuW0hdTakBmMJtHT6YCtQtVdj9upq6yWT9sT1ZICFC94Tm3Kq0QlWRBgwWDUOshIi7RJMbesVz7Pw+OVx9z+645hAMokRS1waaADNY7RLeTFEyEEvM7zwrdD6uvUOq63EYg07v1GqZsBr/0RAD8r8zsG4Akiep2IPm7AawkzEU5gVZ9bSNa8mAED5hu0KrPxlJDBICL0eZwN7fZ+ZI8iB3KXhnAUZ9Dvwlw8jfmUvt1fKFJZ8h1QCgH0hOj0yGp0d9rBWOULPWMMU1U8DJtVCdfVajCiNRiM69coVWxvVJiMeG42juOTUWwa8gnl2ArhjbvDBoqKToaTcNrOv+xeo8555yKVzYxISOqqgttNAP4YwLtqeVEi+jKADIAHytzlRsbYVgB3AfgkEd1c4bk+TkS7iGhXMBisZVkAKjdZVYPPOxi7wDSlGGOYjqXQ4xb7uynd3g00GHuVkMR1JabrVYM3qO0+q0+1NlhlRgigyoPo8jC4cJ+2kBRQWR4kHM8glc2VlDYvxC+YD6lELCWutlsMH27GK6BK8R+qMu1X3llambYSq3rdsFsp32NUK4wxjIcXq2Bzbax6NAkajUhI6tMFt98EsBWAR+8LEtGHALwDwPtZmRZaxtio+u8kgAcBXF1hfd9ijG1jjG0LBLSHHIqZiCR05S+AhQ/+xePaOlBbHWUoD0OPWyx52ed2NKzbO5XJYc/IHG7dENDlNfLdn94hR9VmhADQn8NIZECkhPlE4dPyKinWBgW79rs6tU3wK4UW8cRi3E4btq7owmSFc+nHb4zAbiVcs1r7ZsFiIfR7XYbN944kM5hPZRfls5w2K9YE3IZNMKwnetRqYwCqt0qWgIi2A/g9AO9ijJXsuCEiNxF5+f8B3AFgf6n7Gg1jDBPhBAbK1MxXo8ftABEuuNne/EQX9TD6PM6GNe79bL8SjuLxbq0Mqr0ievIYjDFFR6pqSMqJWCqr+byJqLMwtJSkcsXc6Vj5C301HSmOv8OOWRNDUoBSfbZvdK5kTmzXqWlMhJO4a0v5QVnVGPTXPg6WM5GvwFycH10T8LSECKFIDuOnRPSwensEwBEou/5qj/s+gJ0ANhDRCBF9FMA/AfACeFItmf2met8lRPSY+tABAC8Q0R4ArwJ4lDH2c13vTiPhRAaJdE63hwEAW5b4K+542hFuMMr1GhTTqyrWNmIEJp/L/e4r9dVXeJw2eF02XbvMaDKDZCZXVS6lR+eM7Ggis2jedjVEQlL84ttbJSTV1VHbTIx0NodUJgePxsq1QrjnwENPhXzl4QMAgM8VjWHVwqDPOA/jdTXXUjyqF1AS7EcnoqZI5mhB5JP6m4L/ZwCcZoyNVHsQY+y+Eoe/Xea+5wDcrf5/GMBlAusynCPjim5MqbGNogz4nPmRnBcK/ELXLWgw+jxKfbtos59eMtkcnjkSxMreTs1qqIUM+lzYcUx7fiyv4FuhnwFYMLRT0ZTw0C5AnbansSSVf0aVxqtWEx7keF22fOJdD7yMtBYPY/tmpZDhP185g8/dvj5//PB4GAfOhbG6z43VAd0RdAz4XHh03xiyOQarwNz0SnBDvGFwsTQN90LPTM9jbb/+9dabih4GEVkB/DFj7Dn19qKIsWhVzqq6NEs1fGmLCXhdhiXJWgU9HgZQ/+Y9LjZ41yX6QxIA4LRbdNXILwgyVt6A5HtTNIbpIgntCWO3wwq7lSqGkt5UE/zdVYy5x2lDtIbeAf5YvUlvQGmuvG1jPyYjyXxhQi7H8LH7lVaxv7j3Et3PDSzIrpw0QBxwIpxEV6e9ZDXhlRcp81mavcKyosFgjGUB5IioulhLGzARUT4sUXG6Ugz4nJiKpZretTQSPR4GUP9+le+/qqiTfuj6lTU9zy3rAwhFU5o7fkWTx71uffIgkURas+dERFXlQWbUyqdqO2q304Z5VQBRDzG1j6QWDwMAPn3bWgDAb//H64glM/jrx49gZCaOjYNeXcnuQq7OFz3UfiGfCCfKSg7xRHizV1iKfFJRAPuI6EkoCW8AAGPsM3VblUlMhpPwumyau4EL4fmPUDSpKbzQykzHUnDYLHALVuv01aCfJEoux/Dc0SDWBNw15aQA5TPNqiNLtUhL8NBO1bJaj74cRiSREVJeLcbfYa84x2I6lsyL41WC93/EUhld8h7RfEhKX5UU54oV3bj94n784tAkfus/Xs+X2d7/kbLFlcLwAhi+mayFiUiyrAglP0dPTTW3zLlIldT/A/BHAHagzSfuKRVStV1ceA+HUYmyVmA6lkJPp0O4WkdvCEYLPKzya1dpk4MoBZcm1/qZBqNJWC2ErioNjV6nDQ6rRfMkwkhSn6xGtZkYk+HqlV3Agme/HPCSAAAgAElEQVSgV9IiZkBIivNP79sKYKEn47kv3FrzdxlYyGdO6CyrLmSywvXFZbeiw27FK8PTJX/fLIh8Ul2Msa8XHiCiz9ZpPaZSS0kth19cXj05jStWVJ4b3S4oTXviyevuTgcshLr2YnznxZMAgDs3a+/uLoafE5ORBADx6KzSg+GApUpoh+tJaZ3tregwad/Z+zvsZaubtHhSHoMMRq0hKUC54D7xOzdjOBjDki5XXv68VjxOGzxOfVVyheRyDJORZMXry+qAG8kmD2WLeBgfLHHsQwavoymYCCd1yZoXwvMftSQDW43p+VTVEsxCrBblAlmvbu9MNodH946h1+0w5MIxoHOXKdKDwen1ODQp1mayOcTTWV27c17WXIqZ+RSyOVa1FBiofYiSEUnvQtYPeLF9yyAuXVZSBFs3/T6n7sZNzsmpGLI5VnaQFqCU5Dd7ZKKswSCi+4jopwBWFfRhPExEz6BIVLAdSGVyGJ2N11RSCygaOwM+Z9N/8EYyHUtVragpptddv9nezx5RSmB/RWfvRTH8oq91hnQomhK68AJqt7cGgxHTITzICXidCEaSJftg8pVdAhsnXtLL16KViA4tLDMY8Oorqy7kzJRSgVlOtRhQPNlQNIlMtnm9jEqf1EsAxgD0AfjbguMRAHvruSgzOBGMAljohK2FAZ/LkJhnI/j3nacQjKbwW7es1p3s1xqSApTeBD0zIET47kunAIjNbhbBbrXAaiHNoZdgJJkf9FONXrdDU8IzokPanBPwOJHK5hSF4aLzPSiYqAcAt3q+cE0rreTneRvkYdQLp92CeCoLxpjmQU8c/nctNxcFUPIlOabMUDEi/1IPyn5SjLHTAE4T0fsBnGOMJQCAiDqgyJKfasgKGwTvzub10LXQ73VitAVGLh48F8YfPaR0w1qJ8FkdHbGpTA6RREazweh1O7F3RJ+gXyWyOYYXjoewxO8y9Et3y/qAphr5XE4RmhMPSWlT8OXhHK2d3sCCMQhGEzUZDH6hj+r2MNLosFsrSr83Azes6cOzR4IIJ8Qk/EshUmJdqHbdrAZD5JP6Ic6f5Z0F8N/1WY558BBSpRijKP0+lyF12/Xmn589nv//A68sllYQgYvYifZgcHp1SnpXg4/U/PANxngXnAGNEhHHJhWPVXRWdY/bgflUVnjyGu+w1uVhqBetUnF50d6RwteO6hyiFNVZ5dVo8ga2hhBqMKKU7LsqnA9cMqSZlSJEDIaNMZb/Zqv/r5+eg0nwk6FcnbQW+r2t0bz34vEQVvZ24le3LcNkJKmrkY6XZ2rdefV5nIgmM4YLNT6uGoxfu7r2ctpCBn0uTMVSSAqOa+XG5VpBSfU+jaXGkRoSxnxTFCzxeQcjSXTYrUI9NQsehv6kd7OHo4CFUvnJGnoxRAoguNo1H+PajIgYjCAR5edfENE9AMoL0Lcok+EEfFV2AKIUNu81K4fHw5iZT+PGdX350tMnD05ofh79BoNfII3zMnI5hp/sPocVPZ26GskqMegvvysvRX4DIhiS4kq/os17eoYncSrtmPmFTSRW77BZ4LJbENZZJTWfyhpSUltvjPIwREQo7VZq6vneIgbjtwD8ARGdJaIzAH4fwCfqu6zGMxFO1lwhxeEXiT06h+40ggOjyjjIe7cuw43rlMlleiS8a/EwAGN7MZ5TK1m2b6m996IYfm68KfiZTgoK+HEWZp2LGQyeMNZzwfW5bHDYLBUNhij+Dv1DlKLJjKZZHmZhiMEQGKQFKMKERnSV1wuRAUonGGPXArgYwCbG2PWMsePVHtdqTEYSuiftFcPVJku5/M3CM0cUYb7NS3xw2qxY2tWBH+06q/l59BqMXi4PYmC3949fV3QxP2Jw/gIANgwo1U6iBi4YScLtsApf0Hk3+GxczGCE1byBHk+KiBDwOEvK8Aej2iZOVtOlqkSsRUJS/g47HNbSBlYUkVG9gBIurLXno56IzMMYIKJvA/hvxliUiDapsy3aCqUL0xgPY1l3J6wWauoPfmY+BbfDmp9C1+d1IhgtXZtfibzwoMZyZK5sG4oYF5J6/lgIAa8TgxVq3fUy5HfBYbUI7/4mIwlNHis3uKK79bl4GjYL6d6h816MYnR5GDoNRquEpIio7N9LhHgqi0gyI+Rt9ntLG/JmQSQk9V0AjwNYov58FMDn6rUgM2CMYWQmbpiHYbUoO7hmbt576cQU3nrxwmzhX7p8CdJZplkAbyqWgt1K+kNSBnkYIzPzmIundU/WqwYRod/nzE9Nq4ZIzLoQHzcYcbF8wFw8DX+HXXdfQH+JC2Ayk8VcPK153RGdfRjKtL3mD0kBqoHVGTEIaag86/e5MBVr3uY9EYPRxxjLl9YyxjJQSmvbBl4CWYtKbTEDPicmmnSnMBVNgrHzSzK5su6pqZKTc8sSiiTR6xZLkhbS4VAqcYxSrOVS5r+y1Zju7lIMamjIDEaTCGiouLNbFbVf0d06Nxh6KXUB5B6xFg/D7bTp7vSOJTP55r9mJ+DVLw8yqaG3pd/rBGP1nxWjFxGDESOiXgAMAIjoWgBzdV1Vg+ENWdevrU07v5CAt3l7MXhy++Z1C7vxFapM9i6NJX2haFKTjlQhvR6nYZVkLx6fgt1KuE6wjFUPAz6XcEgqGNbmYQDawjvheBreGg3GdCyFdMFOdo/aSKllCqLHadVVVpvLsZYJSQG1eRgLg7TEDAZQWwlvPRExGP8TwMMA1hDRiwC+B+DTdV1Vg+Gho1qFBwsZ8DVvLHI8P4x+4f1yCQutcdqpWEq4EqiYvgoieFpIZ3PYfXYWN67tq/m5KiEqQsdj1lp26oDS/Cg6IzuazMBXQ9MbX1uhweafhRa1A7fDpkutdj7Nhye1Rkiqv4SBFYUbGpGQt5Fy6vVApErqDQC3ALgeSjntZsZYW2lJ5bu8DWja4wz4XJjW0OjVSMbU91s4jJ6IsCbgxs/2j2t6rlAkqdtgGOVh8BkCN66rT/6CE/AqzYbVpqJp7cHg9LgdmvowaumS5lP+9o4sBAvG5hJwWC3Co3YBoFOduqd1GqGR0uaNoJSBFSUYSYIIQvI5/Du5U6PQZaMQqZJyAfgMgK8C+BMAn1SPtQ181q4RTXscrnt/bCJq2HMaxcRcAhZaaJ7j9KqidKIwxhCKphY9jyh9HqchsdpH940BAO7cPFDlnrVxUY8ilT4crCwSGIwqBlmzh9GpzcOopSSVy/AXesFjc3EM+J1V53cU4lE9hFhKm5dhtLR5veHhpP1q/5IWgpEket0OIc0s7vXHDVZAMAqRkNT3AGwG8I8A/kn9/7/Xc1GNptKsXb0s7VJyAmentSWRG8ELx0Po8zgXncDXre5FKJoUdrsjyQxS2VxNIanpWFLz7rSYHUeD6HRYsaxb+7hSLWxeoiiNVmtwXPAwtJ1TWj0Mj1N/DmNZdydsFsJ4gbf0/LEQhnzaxgpzD2FeY3fyPJ/n3SJJ71V9ymZBT2ltUKMXvmnI17T5TxGDsYUx9lHG2DPq7TehGI224YmDE4aGowBgY4kdXLMQT2XhsC3+6If8LjAmvmbexNbn1Zn0djuQYxDeVZdiOpbC6GwcN62rb/4CWNj9jVcJSfFucK0eRo/bgUgiU9Vg53IM0VRGl/Agx2ohDPhceePHGEMkkUanxpyCXj0pfn+tr2cWfBCXHoPx6skpTeeCUmHZugbjDbUyCgBARNegzWZ6e502XLva2Oqank4HbBZqyl6MYDRZsl+BXxCPT4qF0Xg4icfDtdLn5d3e+g3GM4eVjvV7Ll+q+zlEcdmt6O60480zleVBeDOiVsl3rvhbzYDOp7NgTJ+0eSGDflf+vUzFUkhnGW7V2MfCPQStiW8j53k3AofNgu5Oez7cKApjDOFEBlYNYT6t83R2nZrGT/ec07QuvYh8WlcCeEnVkQKAFQCOENE+AIwxdmndVtcg9v3JnYY/p8VC6Pc6m67aIZXJYTqWKhkuWanuog6eCws1wE1FteklFcMNTSiaxPoBsUFDxfCZGm+9uF/X47XisFlwpkqYcTKSwGXL/JouEoCyyQCAmVi6YjgrktA/PKkQl115L4wxnAopeRnejyOKW6eHwXMerZL0Bsp3x1diVu3cv1lDQUa/z5WfvCeS9/iv185ix7Eg3nnZkqr3rRWRT2t73VfRpvT7XE1XT82rPEqF4C7qVSRNpgW7r/lz6Q1J8WR5LYnvJw5OwOey5SVO6s1N6wJ4WvVqyjEZTmJFr/Z8SrdbyUlUy2NEDZpUd9vGAbx4fAqhaAqP7FUKB3isXhS+Bq3Ne7EWy2EA+gyGlvkinAHfQvOeiMzNeDhxXsVjPRGahwFgXJ3AtwrAPQDmGGOn1WOSMjTjbO/JCiWfRITl3R0YF+5mTinlghrneXO4ZzKls7SWMUXK5JYNjfEuAKXBcTqWqjjHYyKi7wvcIxiSys/CqNHDWK0ah4d2jyIUTcJlt2CdRk+P5yD0hqRapQ8DUCqltDbvhTR0eXO4dym62ZwIN25Cn4jB+DGALBGtBfAtAMsB/GddV9UmNONs78kqkwUHfC5hvaSpaBLdnWLlgqXwd9hhtZDuXoxXTk4jmclhmwFjdUUZ8vPGqtJ/o0Q6i9n5dL6sWgvc8FbL6XCBwq4aOr0B5GXtH903hp/vH8dVK3s0P4fepDcPSRkpx1NvuIehRaBTj4fBN3OvnhRTXRifS9RFcLMUIt/0nKofdS+Af2SMfQHAUH2X1R70e52Yi6fzY0ybgbyHUeaCNuh3YVzQKwpFk7p7MAAlz9Pr1t/tffCcUhNfTzmQYob8Soy/XOKbd4Lr2fHx+drhKvIgXAJdi4RHKexWCy5b5sehsTAyOZaXh9GC26kv6R1JZOB2WDXnecwk4HUikc5pMo48hKUlz8crLEVkYuKpLMKJTFN5GGkiug/AbwB4RD1m7DizNoXv4o82UfPepNp1Wq6bd9DnwpnpeSG1TKVpr7Zy5N4amveeOTKJDrsV69T5I41g/aDyWuU8jP3nlM5pPV9gp80Kl91SdYMxa5CHAQBfvOtiXLu6F7dfPIDfvGm15sd32vWFpCKJNLwGT0WsN/nRthryGMFIEg6bRZOMi9NmFQ5nH5mIAEDDchgi7+LDUKbu/Rlj7CQRrUKbNe7Vi8tXdAGA8I69EQQjiYpdp3zXemwyiouHfBWfayqaxCXLumpaT5/HoTsk9frpGfS4HbolvvXQ73XB67SVbd47N6v0aGhNHnNEBAi5wfAZYDCuW9Nbk4dmsRDcDiuiGpPekRqlTcwgkBcGTGJ1QGyTwmXutZ6jAz6XUC6RV7fxeeD1RkRL6iBj7DOMse+rP59kjP1V/ZfW+vBdpmhOoBEoA3LK70auWqnkA0R2N7XIgnD6PE5dU/emoknMp7K4eX39G/aKGfS7yupJjc8l0GG3YpnOL7CIwZiLp+Fz2ZomnKNInGtv3Ks1ad9o9IxqFR3NWsyAT0ztmm9Gtyz1a34NPejLVgpARN8hokki2l9w7GtEdJiI9hLRg0RUcntKRNuJ6AgRHSeiL9ZrjfXG57Khw25tKg9jMlJ5BCc3ctVUWRPpLKKCU8QqoTeHcVLdWd1+cX31o0ox6HflZ6gUMzaXwFCXS7fX09XhyHsQ5ZidT+XzHc2A12XTPESpVi0sM+B6UlpDUvoMhlPoujE+l4DHaWvY37JuBgPKpL7iHo4noUiNXAplct+Xih9ERFYA3wBwF4BNAO4jok11XGfdICIM+l35MEUzUO0E5snwaidrvgejVg/D68R8Kot5jeJ1zx4JAoBwaMBIetwODAdjJSU8Rmbj+UoqPfR6HFWrpCKJjK5Z3vXC67IjktBeVttKPRiA4v3ZraSptDak08MY9LkwO5+uWL4N8JJaY2WNKlHRYBCRlYj+Rs8TM8Z2AJguOvaEWnEFAC8DKDUe7WoAxxljw4yxFIAfQOn9aEkCHmfVRq9GMjufrjh/22lT5C+qhaR4otoIDwOAZi/j6cOTINKfK6iFS1T3f2z2/L9RNsew5+ysZtHBQvoEJN+bLf7vddmqVnYVE0u2zvAkjsVC6POIN+9lsjnd82L4Yw6cq6yOO9bAklqgisFgjGUB3Fin1/4IgJ+VOL4UwNmCn0fUYy3JUJcLmRzTpMj6tccP48qvPokb/vJpHFOrIIwglckhns5WHe2p9I9UMRg6ygVLwR+/b1TbEMe5eBpbVzSu/6KQzUsUg/HSidB5x7lXdpGOLm9Or0cJSVUSIAw3WYWR12XTJT7oaaGmPY6Wbu/pWAqMaRehBIAN6kCz0QrRiVyOYffZ2YaV1AJiIak3iehhIvp1IrqX32p5USL6MoAMgAdqeR71uT5ORLuIaFcwGKz16Qxn28oeZHNM2I09FYrhG8+cgNtpw+hsHF999JBha+HJVBGD8dzRyn9LnqjWO56Vs0YNKWnp9o4k0hidjeOGBvZfFHLpMsVgnCsqZnj1pDL0Rk8DHGeh+728x9VsHoYydU+8Sooxhlgyg84W8zAAtdtb0GDkZ3nr2FStUUvFKykjh2LGbNq0IGIwXACmANwG4J3q7R16X5CIPqQ+/v2sdMvkKJRucs4y9VhJGGPfYoxtY4xtCwTqO3FND8tUMbc9ZysrnHK+/tQxAMAX79qI7k47dhwNIq5x1kA5uMGoVo7pcdqQzrKKvRhGhaSWdnfAQtpGUvKGvUaVEhbjdtqwxO/Cs0fODzWenVa+3JuqlCNXgv89K4Wlwol0U+UwtFZJpbI5ZHKs5ZLegLbZ3iEdXd4cn8sOj9OGc7PlPX0+armWDYpWqn5ijLEPG/ViRLQdwO8BuIUxVk7y8zUA69R+j1EA7wXwPqPW0GjWqjuFagqngLLzevDNUThsFtx9yRCCkSS+8vAB/HTvOfzqtuVVH18N3hBWrUP4+rW9eHTfWEXxs2AkCa/TVvOUQquFEPA6NYk0Pqt6P9esMsfDAJREb7EX8MLxEPq9zrxMuR4CXi7IWPqilMuxmud5G43HaUMslQFjTKg6bEF4sHVDUiJKsnpH9XIG/S68UkEehBuTWoostCIyonU9ET3Fy2OJ6FIi+kOBx30fwE4AG4hohIg+CmVinxfAk0S0m4i+qd53CRE9BgBqUvxTAB4HcAjADxljB3S+P9NZ1t0Bl90i1NfAtWM+cM1FAIBfuVKpCTg8Zkweg1ffVJvZzKcPVqqUmoqlag5Hcfq9Lk2Dpg6Y7GEAwB2bBzA6Gz/P+xubi6OzxovggodROiQVSWbAmDFNe0bhdtqQY+JjRVttnnch3Cs6PF79OxmsUf7fZqGK4a+DqqpAUxkMAP8HSvlrGgAYY3uh7Porwhi7jzE2xBizM8aWMca+zRhbyxhbzhi7XL39lnrfc4yxuwse+xhjbD1jbA1j7M/0vbXmgIiwxN+xKN5disfU2dQfuXElAOXkXNvvwX6NCeFycNnsaoN9BqsI7AFK0tuo2KnWuSFHxsPYvnkQdp2ih0bAdZd+sluJlmZzDGen49i+pTaZtWohKf4ZGmWsjYAnr0UT3602z7uQq1Yp4R8RjzgYScLjtKFD5ybibZsGMB1Llg0N814grYO6akHkG9fJGHu16Ji2kogLnEG/C6Mz1XsxfrL7HHwu23mzqfu9Trx6arpqPbYIogaD92JUNBhRAw2Gz4WgYEhqLp7GRDhp+EhdrfAJf08dmgAA7FZzVEu6atvtdToUPalQmZ0lLw7o0TnlsB64Nc7EiOXHs7aeweDhpWqNrYD+pj3OkL8DuQojk0PRJC5b3tVQaRwRgxEiojUAGAAQ0bsBjNV1VW1Gv9eJ3Wdnka1QWnvg3Bzm4mls3zJ43vHbNiqzHkYEDE41pqIpuB3WqnmHPrcTNgvlk2olnyuW0j04qZh+r1MdEVpd8JB3eF/ZQEnzUjhsFlw85MNTao/No+oAIpFJhZUgooq9GKJhxUaiVbF2wcNozRwGIDb3PhRN6qqQ4vBQU3ndsgTWNLgPScRgfBLAvwLYSESjAD4HRYxQIggfIF8pJ8Dj8h+49qLzjl++XFFP2VlU86+H6VhSKCHLG5SeOVK6tDaTzWFmPqV7lncx/fkJY9W/hC8cU9a0cVB/JZJR3Li2F4wBr52axvdfPQOX3ZL/rGtB0dcqncPgifbmCklpm4kxr+Z9WjGHwRtbRUNSNXkYqrf68vDUot+F1dJyrSN1a0VEfHCYMXY7gACAjYyxG+WkPW3wsjdep1+Kpw8pO9VihVguKqYlKVyOqVhKeGfq67CV3fHzhqS+Gr4MhfAku4ibz/8OjZQ0L8f71OKE93xzJ+LpLO7YNFjlEWJU6ibm43MbGbeuhl4Po9WkQTj9XldDQlLr+5XmvVIbqaNq0n15T5MZDCI6QUQPAPh1ACvqv6T2Y8tSxQjwOv1SvHxyCgGvc1Ei12W3YmlXh1AOpBrTsZTwheb6NX1lVXZ5BU/AqCopNR+xZ6R6r8qDb45iy1IfLE2g1Lqqz533AAHg83esN+R5A15H2SqpqVgKHmfjZpiLoDXp3cpVUoByvlbbwCXSymCjWgeMbRz04pkS0kIvnVA2n9sa2IMBiIWkNkEJSfUC+JpqQB6s77Lai65OB/q9Tjx/rHSIZzKcwOx8umz8e0mXK3+C1IJiMMR2PAM+FyLJTMmLAN/x9BqU9Obd3tV2bdFkBpFEpiatJqP5lw9sxeYlPnz1ns2GhKMAxcOYjiVL5rymosaVMxuF3qR3K83zLiTgdVaVHh+Zmc/ftxZ8HYv7fQDk5fX1TEmsBRGDkYVSUpsFkAMwqd4kGnA7bWVzGD/bPw4AePslpUsyXXYrJiMJTbOEi2GMaeqd4Am3UonvUI315cW4nTYM+V1lk3ucI+NKnufdV5bSrDSHIX8HHv3MTfj161Ya9px9HidyDJgpMXlPi5fYKLSHpLKwW6mpvCQt9HtdCEYrz/Y+pk7ZrPU78raLBxBJZhZJ57xychqXL+9qeGm5yKuFAfxvACcBfJAxdh1j7BP1XVb7sX3LIM5Ox0vuGrlu0/VrS3cu37axX72AaFMELWQ+lUUqkxO+2HBBs7MzizvUT08px2qVNi+k0lAiDp+jvTrQeIXaRsKNeqnYdSiaNKzYwCh4LiIinPTOtGw4ClCq+tJZVvH7yJv2Lq1xIiVPar92aiZ/LJdjOBWKmdLHImIw7gOwA8BvA/gBEf0JEb21vstqP/iO/c0zM4t+9+LxEC5d5i+74xryKydNLXM1RHswOCtUxdXXSkgT7FVzDUaesEv8HSWNUyFnp+dhsxA2DHgNe91mpJIA4bSGwoVGYbUQOh1WTUnvVk14Aws5t0qVUhPhBKwWqvmz4uNzRwq+G8FoEjmmNPY1GpEqqYcYY18A8AkAjwH4EIBH6ryutuMtG5R+ikf2nt/CciIYRTKTww1ry48arVaPLYLW+v2l6s5mtsScg1gyi/UDHkMbhgJeJ85OxytedJ48OIGLejsb2qhkBuW6vRljSkiqyXIYgDYBwnC8udR2tcJ7K0ptpjgTYWWyZa3FGd2ddnTYrfi5GrYGgBeOKSX2K2qQ0deLSJXUj4noOICvA3AD+A0A5nZNtSDLezrhdljxC7UzmPPtF04CAO6osFvgbukLZZLmIugpx9y8xFfSqxkLxxeV/9bKmioijdFkBufmEljf5t4FUH4UaDieQSbHms7DANQxrYJT9yJNprarFR5mqjR7vdooZFGICOsHvdh1eiY/U+fxA4rxuNYE8U2RkNRfANjAGLuTMfa/GGPPMcaaZ0h1C3Hdmj6MzJy/i37iwDgcVst55ZnF8EqLWA0y53xOdHcVpdpClnZ15EehcnI5hvE6TPnaskQxQOXCbsdV3ZxKnli74OuwwWG1LCqtNWoGST3wuuwIJ8RybM02z0MrHQ6lea+Sxz8ZTqDfoMFGN69Tzvkd6obxuaNBLO/p0K1RVQsiBmMPgE8S0Y/U26eJqHW3BybCvQjuZRwaCyMUTWH7lsGqYZbLlncJKd6WQ3R4UiHcUBXupE5PzyOdZVjiN7ZhiCvPvni8dPkwd8kbqf1vFkSEXo9jUUhqIazYXElvAPBp8DCiydY2GAAw6O+oKJ1j5Kzt912jtL/92/Mn8fLwFJKZXM0SNHoRMRj/AuBKAP+s3raqxyQa2X6J0gn8f188BQD4ykOKavvHblpV9bGDPmfFE7Qa3MPQIot90zrlpCwc/sRHxi41WJIgUEWldTioeBhmzPA2g5IGI6qtcKGRKCEpMQ8jmmztKilA/T6W2cAlM1nMzKfzCga1MuTvwNKuDuwZmcVfPKZM4PzwDdWvGfVAxGBcxRj7IGPsafX2YQBX1Xth7YjPZceKnk7sPjuLj92/C6+emsaAzylUejfoc1XUoqrGXDwNr9MGq4Yk3PoBJa/Aw0HAwozhy1fUVi5YDBHhhrW9ZRsUd56Ywls39sNhM0/SvJH0eZyLqqSaOSTlcYqPaVXmebe4wfB35PXfiuGKDkbO2v7y2y+G3WrBqal53Lt1ab7ZtdGIfGpZIlrDGDsBAES0GkoTn0QHf/7Ll+AD334lH5b6y3svFXrcgN+FSCKD+VQGnTpKEufiafg7tUUSV/W50emw4s0CD2PfyBwcNktdEq9epx2haBKJdPY8Rd3JcAKRZKamSXatRp/HiSNFQ3pOTCpKvc3oYXicdiFpkFQmh1Qm1/IGg08LPFdCAPDgmGJIau3yLuTuS4Zwd5nG3kYisl37AoBniOhZInoOwNMAPl/fZbUvN67rw3c+tA0A8I33bcVbVPnyanD3VsugoUJm5rV3CBMRujrs+clegDJpzOey1aW0ldeVF0u5P/imMqTol69YavhrNit9HifG5s7v7j+uhuWasUPa47QimszkK3nK0eo6Upyb1BzCaIkiDS4bcoXBXngzIDLT+ykiWgdgg3roCGOsdunUC5jbNg5g1x/erkk2gDcLvXZyWlccfyaW0lQhxbl2TS8e3z+en9c8EU7kmxeZWogAABSSSURBVImMZqX6vp49MpmfhQ4slNpeX6fXbUZ4UvhkKIbVavghHE837UXIo653Pp2t6D3kZ2G0eNJ7SYXeqPG5BDrsVk0FJq2CSB+GC8pMjD8G8BUA/0M9JqkBrRozly5VLhSipYvFTM+n0K0xJAUAGwe9iKWy2HV6BmNzcUzFUtg4WJ9eiCvU0uLiXpUHXjmDi4d8bd+wV8ilyxRZ+8K+lPG5hGmx62p4nMq5Fa1SKdXK41kL4WXlY6X6lOYSGPK72vJ8FQlJfQ/AZgD/COCf1P//ez0XJVmMr8MGp82iey7GbCytKwfAZzz8aNcIHtmjdKlfUqM+TjksFsLafg/m4gsXHZ5wb9addb3gXiTfwWayOQSjyXzXf7PBPYZosvKGpl1CUl6XHW6HNT9xsZBfHJowvE+pWRD51LYwxjYV/PwMER2s14IkpSEiBLzlB+tUIpXJIZLMoEdHSGplnxtelw3/tetsPtF3Yx2b525c24fvvnQKZ6fnsbynE/++8xQA4N4LKH8BKBU2RAs72FA0hWyOGVp5YyRe1QBU68WItImHASjy/vOp899vJptDUoPIZ6sh4mG8QUTX8h+I6BoAu+q3JEk59BqMWVUmu0vnSfyr25YDUDrNb1rXp6k0Vys3r1eM0Z89egipTA7371SGO25dcWGp0ditFgQ8TuxQdYP2jyqFB83uYVQrrY21kcG4ZX0A+0fDyBRMppxQv5/tqkgg8qldCeAlIjqj/rwCwBEi2geAMcbE6kIlNdPvdeJkKKb5cdOqwdDjYQDAp96yFm+cmUEqk8NX3rlZ13OI8pYN/XDaLPj5gXHc/NfPAAA+dP3Kppiw12g8Lls+Z3VazWU0aw6Dq89WC0nxHEerJ72BhTzG8WA0P2N+n6rk3OhZ241CxMPYDmAVgFvU2yr12DsAvLN+S5MUo9fDmIlxHSl9VRvdbgce/O0b8OhnbjqveqkeEBH+5j2XAQDGwwkM+lz47FvX1fU1m5W3bRrAcDCGbI7hleEpuB1WXGSCQqkIvKqrWkgqn/RuYXlzDg/NFiohDKsbumaYOV8PRMpqTzdiIZLq9HtdmJlPI5XJaep45lpQWmRBzOSdly3Bqyen8e8vn8bOL93WltUmIgyp+YqXToRwemoebmd9+l+MwCOYw+Ahq1Ydz1oIH+R1aGyhwfLVk9PocTsuaA9D0iTwztFyekvlCKsGo0unh2EGX/2lLRj+87ub9gLZCO5SO3ufPDiBIxMRUwbmiMI3I5UkvwFF2txlt8DW4NGi9cDrsmNVnxvHJhWDwRjDC8dCuj35VqD1P7ULiHJzEqoxG1dyGK3WSHQh5i0KGfC50Ot24Htq4r9Zw1GAMnXP67JVNRhz8XTLnYeVWN7TiRePTyGRzmJ2Po1MjuHtTSDhUS+kwWghFkZDajMYc/E0rBZqi8qUC40b1y1U27zrsuYuLfZ32PPebDnm4ml0dbRPyemtqkTIjqNBvDysCGduXuo3c0l1RV5BWggektLqYczMp9Hdab+gwzutypfuuhibhnxY1t3Z9M1g/g57yZG+hbSbh/FLVyzFnz5yEN/beTr/vby6jWe2SIPRQnA5kUrD50uhV0dKYj6Dfhc+ccsas5chhM9lrzoTYy6exrLu5g2taaXH7cCKnk68fnoGqWwOGwe9ba2qLA1GC2G3WtDjdmj2MKZjqbY+iSXNgddlKzuTnRNJZODraK/Lzl/eewm+9fwwrET4nbetN3s5daW9PrkLgH6vU3MOY2Y+hdV97VkXLmkePAJjWiOJNHyu9glJAcD1a/twfZt2dhdTt6Q3EX2HiCaJaH/BsfcQ0QEiyhHRtgqPPUVE+4hoNxFJGZIC9DTvTcfS6Ha315dU0nxUC0kxxtpi2t6FTD2rpL4LpSO8kP0A7gWwQ+Dxb2GMXc4YK2tYLkQCHm0GgzGG2XmZw5DUH4/Thmgyc97Qp0LmU1nk2EJXuKT1qJvBYIztADBddOwQY+xIvV7zQiDgUwxGuS9lMZFkBpkca1v1TEnz4HXZkGOKYShFpI10pC5UmrUPgwF4goheJ6KPm72YZiLgcSKVzSEcrz4/GVAqpACgS3oYkjqzMBOj9LnJhQm9bZbDuJBoVoNxI2NsK4C7AHySiG4ud0ci+jgR7SKiXcFgsHErNIl+VV9ItLR2dl6VBWmj2ndJc8INQbk8BvcwvDKH0bI0pcFgjI2q/04CeBDA1RXu+y3G2DbG2LZAINCoJZpGwKOt23uuBXWkJK0JNwThMpVSeYMhQ1ItS9MZDCJyE5GX/x/AHVCS5RJo7/bmnbft1F0raU6qCRDmpc2lwWhZ6llW+30AOwFsIKIRIvooEf0yEY0AuA7Ao0T0uHrfJUT0mPrQAQAvENEeAK8CeJQx9vN6rbPV4HpSogaDf3n90sOQ1Bmu0jo3X9pg8PCo3Ly0LnUz9Yyx+8r86sES9z0H4G71/8MALqvXulodr9MGl90inMOYjqpJ7zYSfJM0J7ywYkad8FgMPy5LvFuXpgtJSSpDRAho6PYORZPo6rRrGrgkkejB32EHkSJ2WYqZWAodditc9tYfnnShIq8iLUi/1yUckgpGknnRQomknlgtBJ/LjrmyHkZa9gO1ONJgtCBa9KRC0WS+skoiqTddnfbyHsZ8SlbrtTjSYLQgAa8Tk2GxHEYwmkSfVxoMSWPo6rCXrZKakRI1LY80GC1Iv9eJcCKDRLq0BEMhoYj0MCSNw1dhiNJMTHoYrY40GC2IaC/GfCqDWCqLPq/c1UkaQ1eno+yY1qlYSubTWhxpMFqQfi+XB6lsMEIRJfkoPQxJo/B32EqW1SbSWUQSGfR55OallZEGowUR9TCCUeX3MochaRS9bidm59NIZ3PnHZ9WRTB75ealpZEGowXpzxuMyolvblCkhyFpFFyJYCp6vpfBf+6VZbUtjTQYLUivxwkLCYSkVA8jID0MSYNYEMc8fzMTiinnovQwWhtpMFoQq4XQ43ZiMlwlJBVJggiyWUrSMPLy+0XnJvcwZA6jtZEGo0Xp9zrzOYpyhKJJdHc6YLfKj1nSGHi4tNj75d6urJJqbeSVpEXp9zmrChAGZQ+GpMH0lQlJTUWTcNkt6HRIHalWRhqMFiXgqR6SCkWTsgdD0lAcNgt63I5FFXxT0RR63U4QkUkrkxiBNBgtSr/PialYCtkcK3ufUFQ2SkkaT8CzWOssFEvJ/EUbIA1Gi9LvdSGbY/n69lLIkJTEDJRwabGHkZQVUm2ANBgtSsBbOlbMiSUziKezsmlP0nACXieC4eIcRkr2YLQB0mC0KANlyhc5smlPYhb9XheC0SQYU8KljDFMxaRqcjsgDUaLMuhXDMbYXGkPIyRlQSQm0e91Ip1l+Rne4XgG6SyTHkYbIA1Gi9LvdYIIGC8zF0N6GBKzCBT1YvAub1mA0fpIg9Gi2K0W9HmcmKjqYchdnaSx9Bfl1/I6UrJKquWRBqOFGfS5MFbBw7CQoh4qkTQSLg/CvdwpdfMiz8XWRxqMFmbQ78LBc3MlfxeMptDjdsBqkY1SksZSLA8SikkdqXZBGowWptNhRSiaKjmqNRRNypixxBTcThs6HdZ8Bd/IzDwAKYLZDkiD0cJcs6oXADBeIo8RjCSlrLnENPq9C1pnr5+aQYfdCpsUwWx55CfYwqzo6QRQulJKehgSM+n3uvI5jPlUFusHPCavSGIE0mC0MLwXo9jDYIxJD0NiKgGvE6+cnAZjDCMz87hseZfZS5IYgDQYLUzeYBR5GNFkBslMTiYZJaYxpJ6bO09MIZzIYGlXh8krkhiBNBgtjMdpg9dpW+RhHJ2IApCjWSXmcd0aJb/2148fAQBceVG3mcuRGIQ0GC3OgN+FHUeD5x07O61UpSzr7jRjSRIJrlrVAwDYfXYWdivhihXSYLQD0mC0OG6HFZFk5rxjo7NxAMDmJT4zliSRwOey49O3rcWmIR/uu3qF7AdqE2xmL0BSG7dtHMDf/+IokpksnDZl/OXYXBxdnXZ0OuTHKzGPz9+xAZ+/Y4PZy5AYiPQwWpxl3Uoy8WQolj+269QMlvhlklEikRhL3QwGEX2HiCaJaH/BsfcQ0QEiyhHRtgqP3U5ER4joOBF9sV5rbAdW9Cp5it1nZvPHRmbicNrlXkAikRhLPa8q3wWwvejYfgD3AthR7kFEZAXwDQB3AdgE4D4i2lSnNbY8ly1T6tt5ZdTcfBrRZAZ3bh40c1kSiaQNqZvBYIztADBddOwQY+xIlYdeDeA4Y2yYMZYC8AMA99RpmS2Pw2bBoM+FZ49OAgB2DocALISqJBKJxCiaMW6xFMDZgp9H1GMlIaKPE9EuItoVDAbL3a2t2bLUj+FgDPFUFiMzSoXUdat7TV6VRCJpN5rRYGiCMfYtxtg2xti2QCBg9nJM4Zb1fQCAU1Mx/Gz/ODxOm1QGlUgkhtOMBmMUwPKCn5epxyRl4E1Rf//kUbx+egZrAm4Qybp3iURiLM1oMF4DsI6IVhGRA8B7ATxs8pqamouHfLBbCU8dVvIY925dZvKKJBJJO1LPstrvA9gJYAMRjRDRR4nol4loBMB1AB4losfV+y4hoscAgDGWAfApAI8DOATgh4yxA/VaZztgtRD+7YNX4c7NA3jvVcvxa1ctr/4giUQi0Qgxxsxeg2Fs27aN7dq1y+xlSCQSSctARK8zxsr2xRXSjCEpiUQikTQh0mBIJBKJRAhpMCQSiUQihDQYEolEIhFCGgyJRCKRCCENhkQikUiEkAZDIpFIJEJIgyGRSCQSIdqqcY+IggBOm70OjfQBCJm9iAYj3/OFgXzPrcFFjDEh5da2MhitCBHtEu2ybBfke74wkO+5/ZAhKYlEIpEIIQ2GRCKRSISQBsN8vmX2AkxAvucLA/me2wyZw5BIJBKJENLDkEgkEokQ0mCYDBF9nogYEfWpPxMR/QMRHSeivUS01ew1GgURfY2IDqvv60Ei6ir43ZfU93yEiO40c51GQ0Tb1fd1nIi+aPZ6jIaIlhPRM0R0kIgOENFn1eM9RPQkER1T/+02e61GQ0RWInqTiB5Rf15FRK+on/V/qVND2wZpMEyEiJYDuAPAmYLDdwFYp94+DuBfTFhavXgSwBbG2KUAjgL4EgAQ0SYoo3g3A9gO4J+JyGraKg1EfR/fgPK5bgJwn/p+24kMgM8zxjYBuBbAJ9X3+EUATzHG1gF4Sv253fgslMmgnL8C8PeMsbUAZgB81JRV1QlpMMzl7wH8HoDCRNI9AL7HFF4G0EVEQ6aszmAYY0+oI3gB4GUAfPj4PQB+wBhLMsZOAjgO4Goz1lgHrgZwnDE2zBhLAfgBlPfbNjDGxhhjb6j/j0C5gC6F8j7vV+92P4BfMmeF9YGIlgF4O4B/U38mALcB+JF6l7Z7z9JgmAQR3QNglDG2p+hXSwGcLfh5RD3WbnwEwM/U/7fze27n97YIIloJ4AoArwAYYIyNqb8aBzBg0rLqxf+GsuHLqT/3Apgt2BS13WdtM3sB7QwR/QLAYIlffRnAH0AJR7UVld4zY+wh9T5fhhLGeKCRa5PUFyLyAPgxgM8xxsLKhluBMcaIqG1KMonoHQAmGWOvE9GtZq+nUUiDUUcYY7eXOk5ElwBYBWCP+qVaBuANIroawCiA5QV3X6YeawnKvWcOEX0IwDsAvJUt1HS39HuuQju/tzxEZIdiLB5gjP0/9fAEEQ0xxsbUsOqkeSs0nBsAvIuI7gbgAuAD8HUoIWSb6mW03WctQ1ImwBjbxxjrZ4ytZIythOK6bmWMjQN4GMBvqNVS1wKYK3DrWxoi2g7FhX8XY2y+4FcPA3gvETmJaBWUhP+rZqyxDrwGYJ1aPeOAktx/2OQ1GYoau/82gEOMsb8r+NXDAD6o/v+DAB5q9NrqBWPsS4yxZer3970AnmaMvR/AMwDerd6trd4zID2MZuQxAHdDSfzOA/iwucsxlH8C4ATwpOpZvcwY+y3G2AEi+iGAg1BCVZ9kjGVNXKdhMMYyRPQpAI8DsAL4DmPsgMnLMpobAPw6gH1EtFs99gcA/hLAD4noo1BUpH/VpPU1kt8H8AMi+l8A3oRiSNsG2ektkUgkEiFkSEoikUgkQkiDIZFIJBIhpMH4/9u7nxCbwjCO49/flEyNkqlJNpQFO5GRDVM0TTayGGUxG1kMU8bKQlnJjgViQaYhMyXRqEkYGWlqoiE18m+yMEVKSs1CQ+GxOO/NaaZ759Dthvl9Vu+97znvuffU7bnvec95HjMzK8QBw8zMCnHAMDOzQhwwzABJByS9lOSnz83K8G21ZoCkV0BrRLzLvVd6YtfM8AzDDElngZXALUlTkvokjQJ9qd7BcUmPUh2PvWkfSTqT6lzclXRT0s7UN5mrb9Is6X5qN0jqlTSWaijsSO/vljQg6XaqHXEs99m2SXoiaVzSsKS6tE1T6q9LtReaannObH7yk94270XEvpS2ZAuwH9gObIqIaUmdZOlZNkhaCIxKukOWkXU1WY2LpWRPqffOcajDZCkk9qTiUWMpWSPA2jTmV2BC0mngC3AeaImIN5IaI+KHpH6ggyxbaiswHhEfq3U+zMpxwDCbbTAiplO7DVhTmj0Ai8lyXbUAl1MKk/eS7hUYt40sYd3B9LoeWJ7awxExBSDpBbACWAKMpBohRMSntG0vWY6ik2Rp4i/82dc0+z0OGGazfc61BXRHxFB+g5SltJxv/LrcWz9jrPaImJgx1kaymUXJdyr8NiPiraQPkraSFWjqqPBZzKrGaxhmlQ0BXSl9N5JWSWoARoBdaY1jGdnlrJJJYH1qt88Yqztld0XSujmO/RBoSRl8kdSY6+sB+oGr/0uiRvv7OWCYVdZDtj7xRNIz4BzZv//rwOvUdwl4kNvnCHBK0mOy2ULJUWAB8FTS8/S6rLQu0QkMSBoHruS6B4FF+HKU1ZBvqzWrAkkXgRsRcW2ubat0vGbgRERsrsXxzMBrGGb/HEmHgC68dmE15hmGmZkV4jUMMzMrxAHDzMwKccAwM7NCHDDMzKwQBwwzMyvEAcPMzAr5CUy9SekkxA8QAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from __future__ import division\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "data = np.array(raw_emg[:3000])\n",
    "\n",
    "ps = np.abs(np.fft.fft(data))**2\n",
    "\n",
    "time_step = 1 / 100\n",
    "freqs = np.fft.fftfreq(data.size, time_step)\n",
    "idx = np.argsort(freqs)\n",
    "\n",
    "\n",
    "plt.xlabel(\"frequency\")\n",
    "plt.ylabel(\"power spectrum\")\n",
    "plt.plot(freqs[idx], ps[idx])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}