cde_net_ictals.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ad hoc analysis for ICTALS\n",
    "\n",
    "### Outline\n",
    "1. Show changes in single neuron behaviour\n",
    "2. Show change in avalanche size\n",
    "3. Quantify avalanche size in GABRG2\n",
    "4. Show network cycle differences"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Housekeeping"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import re \n",
    "import glob\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import cde_net_ictals as cn\n",
    "import networkx as nx\n",
    "from scipy import stats\n",
    "\n",
    "from IPython.display import set_matplotlib_formats\n",
    "set_matplotlib_formats('pdf', 'svg')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define folder locations and such\n",
    "#---------------------------------------------------------------------\n",
    "F        = {}\n",
    "F['dom'] = '/Users/roschkoenig/Dropbox/Research/1901 Fish2p/suite2p data/Project'\n",
    "F['gb2'] = F['dom'] + os.sep + 'LIGHTSTIM-GABRG2--'\n",
    "F['gct'] = F['dom'] + os.sep + 'LIGHTSTIM-WILDTYPE'\n",
    "F['ptz'] = F['dom'] + os.sep + 'PTZ' +os.sep+ 'filtered'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Locate and load relevant datasets "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### PTZ dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Load actual datasets - PTZ \n",
    "#---------------------------------------------------------------------\n",
    "conditions = {'bln':'Baseline', 'p05':'PTZ 5mM', 'p20':'PTZ 20mM'}\n",
    "spatial_conversion = [.5, .5, 15]\n",
    "Fish = []\n",
    "\n",
    "for key in conditions.keys():\n",
    "    folds = glob.glob(F['ptz'] + os.sep + '*'+key)\n",
    "    for f in folds:\n",
    "        tracefile = glob.glob(f +os.sep+ '*all_sig.npy')[0]\n",
    "        coordfile = glob.glob(f +os.sep+ '*all_cord.npy')[0]\n",
    "        fish                = {}\n",
    "        fish['Voxels']      = np.load(coordfile)\n",
    "        fish['Coordinates'] = np.multiply(fish['Voxels'], spatial_conversion)\n",
    "        fish['Data']        = np.load(tracefile)\n",
    "        fish['Condition']   = conditions[key]\n",
    "        Fish.append(fish)\n",
    "        \n",
    "        if fish['Coordinates'].shape[0] != fish['Data'].shape[0]: \n",
    "            raise Exception('The data and coordinates do not match in length for: \\n \\n' \n",
    "                            + coordfiles[k] + ' \\n and \\n' + tracefiles[k])\n",
    " \n",
    "        if fish['Data'].shape == (0,0):\n",
    "            raise Exception('Data are empty for : \\n \\n' \n",
    "                            + coordfile + ' \\n and \\n' + tracefile)\n",
    "    \n",
    "    \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:20: RuntimeWarning: divide by zero encountered in log\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlMAAAJCCAYAAADky0LWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzs3X10XPd93/nPbwYDzhAyAVJiChCkIrFxmbMWYYNBclzqYVurBZMdy6ZVmXK67XG7SdycdNeIkqUNqpI8YXQkOmwrQ91kU9fJlt22sWiWpmLPSYlW2qye6iQUIYPupowT2g4JArUocmALGJADzG//uLiDebh3noF5er/O8aF4cWfulc4x9dXv9/19vsZaKwAAAFQn0OgXAAAAaGUUUwAAADWgmAIAAKgBxRQAAEANKKYAAABqQDEFAABQA4opAACAGlBMAQAA1IBiCgAAoAZdG/mwO+64w951110b+UgAAICqvPnmm9estdtL3behxdRdd92lc+fObeQjAQAAqmKM+V4597HNBwAAUAOKKQAAgBpQTAEAANSAYgoAAKAGFFMAAAA1oJgCAACoAcUUAABADSimAAAAakAxBQAAUAOKKQAAgBpQTAEAANSAYgoAAKAGFFMAAAA1oJgCAACoAcUUAABADSimAAAAakAxBQAAUAOKKQAAgBpQTAEAANSAYgoAAKAGXY1+gfVwZmpGx89e1NVEUjv6Ijp8YI8ODg82+rUAAEAbarti6szUjI6cvqBkakWSNJNI6sjpC5JEQQUAAOqu7bb5jp+9mCmkXMnUio6fvdigNwIAAO2sZDFljPldY8z3jTHfyrp23Bjz34wx08aYrxpj+tb3Nct3NZGs6DoAAEAtylmZ+teSfjrv2n+SdI+1dkjSn0k6Uuf3qtqOvkhF1wEAAGpRspiy1r4i6XretUlr7fLqb78haec6vFtVDh/Yo0gomHMtEgrq8IE9DXojAADQzurRgP6/SHqhDt9TF26TOaf5AADARqipmDLG/BNJy5L+XZF7PiXpU5J055131vK4soV631LPj03oPQtz6unpV6h3TBLFFAAAqL+qT/MZYz4p6cOS/mdrrfW7z1r7RWvtiLV2ZPv27dU+rmzxS3HF3ohpdmFWVlazC7OKvRFT/FJ83Z8NAAA6T1XFlDHmpyV9VtJHrLWL9X2l2kycn9DSylLOtaWVJU2cn2jQGwEAgHZWTjTC70n6L5L2GGOuGGN+TtL/Iek9kv6TMeYtY8xvr/N7lm1uYbai6wAAALUo2TNlrf1Zj8u/sw7vUhf9K1azQeN5HQAAoN7aLgF97J3rCqfTOdfC6bTG3rnu8wkAAIDqtd1svmjXNunaO5rY2qe5rqD6l1c0diOhaNftjX41AADQhtqumNKDTyn6tU8reuXq2rVQRDrwVOPeCQAAtK222+bT0CHpoeel3l2SjPPrQ8871wEAAOqs/VamJKdwongCAAAboD2LqRZyZmqG0TcAALQwiqkGOjM1oyOnLyiZWpEkzSSSOnL6giRRUAEA0CLar2eqhRw/ezFTSLmSqRUdP3uxQW8EAAAqRTHVQFcTyYquAwCA5kMx1UA7+iIVXQcAAM2HYqqBDh/Yo0gomHMtEgrq8IE9DXojAABQKRrQG8htMuc0HwAArYtiqsEODg9SPAEA0MLY5gMAAKgBxRQAAEANKKYAAABqQDHVaNMnpefukWJ9zq/TJxv9RgAAoAI0oDfS9Enpa5+WUqshnfOXnd9LDGoGAKBFsDLVSC8dXSukXKmkcx0AALQEiqlGmr9S2XUAANB02OZrpN6dztae1/U6OzM1QzgoAADrgJWpRnrwKSmUN4cvFHGu19GZqRkdOX1BM4mkrKSZRFJHTl/QmamZuj4HAIBORDHVSEOHpIeel3p3STLOrw89X/fm8+NnLyqZWsm5lkyt6PjZi3V9DgAAnYhtvkYbOrTuJ/euJpIVXQcAAOVjZaoD7OiLVHQdAACUj2KqweKX4ho9NaqhE0MaPTWq+KV43Z9x+MAeRULBnGuRUFCHD+yp+7MAAOg0bPM1UPxSXLE3YlpaWZIkzS7MKvZGTJIU3R2t23PcU3uc5gMAoP6MtXbDHjYyMmLPnTu3Yc9rdqOnRjW7MFtwfaBnQJOPTDbgjQAAgMsY86a1dqTUfWzzNdDcwlxF1wEAQPOhmGqg/p7+iq4DAIDmQzHVQGP7xhQOhnOuhYNhje0ba9AbAQCAStGA3kBuk/nE+QnNLcypv6dfY/vG6tp8njF90hmgPH/FGVfz4FPrnm8FAEAnoJhqsOju6PoUT9mmT0pf+7SUWg3pnL/s/F6ioAIAoEZs83WCl46uFVKuVNK5DgAAasLKVCeYv1LZ9SLOTM2QVwUAQBZWpjpB787Krvs4MzWjI6cvaCaRlJU0k0jqyOkLOjM1U/s7AgDQoiimOsGDT0mhvDl8oYhzvQLHz15UMrWScy2ZWtHxsxdrfUMAAFoWxVQnGDokPfS81LtLknF+fej5ipvPryaSFV0HAKAT0DPVKYYO1Xxyb0dfRDMehdOOvojH3QAAdAZWplC2wwf2KBIK5lyLhII6fGBPg94IAIDGY2UKZXNP7XGaDwCANRRTqMjB4UGKJwAAsrDNBwAAUAOKKQAAgBpQTAEAANSAYqpDxC/FNXpqVEMnhjR6alTxS/Hqvmj6pPTcPVKsz/l1+mR9XxQAgBZDA3oHiF+KK/ZGTEsrS5Kk2YVZxd6ISZKiu6Plf9H0Selrn14bmjx/2fm9VHOGFQAArYqVqQ4wcX4iU0i5llaWNHF+orIveunoWiHlSiWd6yWcmZrRvcde1t3jcd177GXm+QEA2gYrUx1gbmGuouu+5q9Udn2VOyDZnevnDkiWRMwCAKDlsTLVAfp7+iu67qt3Z2XXVzEgGQDQziimOsDYvjGFg+Gca+FgWGP7xir7ogefkkJ5c/hCEed6EQxIBgC0M7b5OoDbZD5xfkJzC3Pq7+nX2L6xyprPpbUm85eOOlt7vTudQqpE8zkDkgEA7cxYazfsYSMjI/bcuXMb9jw0h/yeKckZkPzsw3vpmQIANC1jzJvW2pFS97EyhXXHgGQAQDujmMKGYEAyAKBd0YAOAABQA4opAACAGlBMYWMw0w8A0KbomcL6Y6YfAKCNUUxh/RWb6edRTJ2ZmuHkHwCgZVBMYf1VMNOPOX4AgFZDzxTWXwUz/ZjjBwBoNRRTqEj8Ulyjp0Y1dGJIo6dGFb8UL/2hCmb6MccPANBqKKZQtviluGJvxDS7MCsrq9mFWcXeiJUuqIYOSQ89L/XukmScXx963rNfym9eH3P8AADNip4plG3i/ISWVpZyri2tLGni/ETpoclDh8o6uXf4wB7POX6HD+yp6p0BAFhvFFMo29zCXEXXq8EcPwBAq6GYQtn6e/o1uzDreb2emOMHAGgl9EyhbGP7xhQOhnOuhYNhje0bq++DSEsHALQQVqZQNrcvauL8hOYW5tTf06+xfWOl+6UqQVo6AKDFGGvthj1sZGTEnjt3bsOehxb03D1OAZWvd5f02Lc2/n0AAB3LGPOmtXak1H2sTKG5VJiWTqM6AKDR6JlCcykzLd0dOzOTSMpqbezMmamZ9X9HAACylCymjDG/a4z5vjHmW1nXthlj/pMx5turv25d39dExygzLZ2xMwCAZlHOytS/lvTTedfGJb1krX2vpJdWfw/Ursy0dMbOAACaRcmeKWvtK8aYu/Iuf1TS31j96xOS/lDSZ+v4XuhkZaSl7+iLaMajcGLsDABgo1XbM/VXrLWzkrT664/U75XQjqoakFzE4QN7FAkFc64xdgYA0AjrfprPGPMpSZ+SpDvvvHO9H4cm5A5Iduf6uQOSJVWdUXVweFCDl7+uXeeP60fs2/q+2a7L+w7rJ4fzd6QBAFhf1a5M/XdjzIAkrf76fb8brbVftNaOWGtHtm/fXuXj0MqKDUiu2vRJ/eSFz6lfbytgpH69rZ+88DnS0gEAG67aYur3JX1y9a8/KenF+rwO2tG6DEh+6ehaSrorlXSuAwCwgcqJRvg9Sf9F0h5jzBVjzM9JOibpbxtjvi3pb6/+HvDkNwi5pgHJFYR7AgCwnso5zfezPj96sM7vgjY1tm8sp2dKqsOA5N6dPmNn1sI9SUgHAGwEEtCx7qK7o4rtj2mgZ0BGRgM9A4rtj9U2ILlEuCcJ6QCAjcJsPmyI6O5obcVTPjeH6qWjztZe706nkFq9XiwhndUpAEA9UUyhdRUJ9yQhHQCwUdjmQ1vyS0InIR0AUG+sTKEtHT6wR6999bf0y/qydphrumrv0Bf0Cd134Jca/WoAgDZDMYWmE78U18T5Cc0tzKm/p19j+8Yq7rc6GHxdHw59SV2rJwh3mms6FvySuoLvl1R87h8AAJVgmw9NxR09M7swKyubGT1T8Sy/l45mCilX18oSoZ4AgLqjmEJTqdvoGUI9AQAbhGIKTaVuo2eywjvLug4AQJUoptBU6jZ6pkSop+QEe9577GXdPR7XvcdeJtATAFAViik0lbF9YwoHwznXqho9M3RIeuh5qXeXJOP8+tDzmVwqEtIBAPXCaT40FffUXq2n+SQVDfUkIR0AUC8UU2g6dR894+FqIqmPBF7TZ7pOZnKofmP5kL6WuG9dnwsAaD8UU+hIn7ztj/WZ1Je02dyStJpDFfqStoW6Ja1vIQcAaC/0TKEjfSb0QqaQcm02t/SZ0AsNeiMAQKuimEJH2pz0jlrwuw4AgB+KKbSs+KW4Rk+NaujEkEZPjVaWkk4OFQCgTiim0JJqHjtTRg4VAADloJhCS6p57EyJHCoAAMrFaT60pLqMnfHJoTozNaPjZy/qaiKpHX0RHT6wh+wpAIAviim0pP6efs0uzHper8WZqRm99tXf0gv6snZsuqari3foC1/9hKRfoqACAHhimw8tqW5jZ/K8Ff+ijpovamfgmgJG2hm4pqPmi3or/sWavhcA0L4optCSorujiu2PaaBnQEZGAz0Diu2P1Zyc/vO3/q1n/tTP3/q3NX0vAKB9sc2HlrUeY2d2BN6p6DoAAKxMAVmWIt49V37XAQCgmELbqibUc/PPHNVyXi/WcjCszT9zdL1eEwDQ4tjmQ1tyQz3dLCo31FNS8a3BoUPO/yleOirNX5F6d6rrwacyEQp/8vv/UrvOH9eP2Lf1fbNdl/cd1k9+5B+t798MAKCpGWvthj1sZGTEnjt3bsOeh841emrUMzphoGdAk49MVvWdf/L7/1L3vPmEIlkN6knbrW/9xNMUVADQhowxb1prR0rdxzYf2lJdQj3z7Dp/PKeQkqSIuaVd549X/Z0AgNZHMYW25BfeWUuo54/Yt32uX6v6OwEArY9iCm1pPUI9v2+2+1y/o+rvBAC0PooptKX1CPW8vO+wkrY751rSduvyvsO1vi4AoIVxmg9tq96hnj/5kX+kP5FWT/Nd0/fNHbr8E5zmA4BOx2k+dKz4pbgmzk9obmFO/T39Gts3Vl3xNX0yJ0pBWVEKAIDWVe5pPlam0JGqzqHKN31S+tqnpVTS+f38Zef3EgUVAHQIeqbQkSbOT2QKKdfSypImzk9U9kUvHV0rpFyppHMdANARKKbQkeqWQzV/pbLrAIC2QzGFjlSvHKpFnwHIftcBAO2HYgodqV45VL+RelSLeXEJi7Zbv5F6tOZ3BAC0BhrQ0ZHcJvNaT/OdePendD1wS5/pOqkd5h1dtbfrN5YPyaRuSc/dwwk/AOgARCMAeSqJTLj32MuaSeQ2oH8k8Jo+3/07iujm2sVQRHroeQoqAGghDDoGquBGJswuzMrKZiIT4pfinvcfPrBHkVAw59pnQydzCymJE34A0MYopoAslUYmHBwe1LMP79VgX0RG0mBfRDvMO95fzgk/AGhL9EwBWaqJTDg4PKiDw4NrF57b6YR35uvdWevrAQCaECtTQJa6RCY8+JSW804KLgfD0ntHnab0WJ/z6/TJWl4VANAkKKaALPWITDizcq/GUz+vK+k7lLZGV9J36IXU/Vqe+nerK1Z2bewMBRUAtDy2+YAspSITyjnpd/zsRc3c2q9T2p+59lr3p9WV14uVaUrnhB8AtDSKKSBPdHfUMwqh3OHIV/OiEiRph7nm/TCa0gGg5bHNB5Sp3JN+O/oiBZ+9au/w/lKa0gGg5VFMAWUq96SfV/bUF/SJgqZ0BULSrQUa0gGgxVFMAWUq96SfV/bUfR/7JXV99F9IvbskGSmyTTJGSl4XDekA0NoYJwOUKb9nSnJO+sX2xyqe6efM7fPKotolPfatGt8UAFAP5Y6ToQEdKFOtw5HPTM3o+NmLuppI6i/CV7yXhWlIB4CWQzEFVMDrpF85cQlnpmZ05PQFJVMrkqSr6du1M+Bxwo+GdABoOfRMATUodzDy8bMXM4WUJP3G8iEt2u7cLwtFpAef2ojXBgDUEcUUUINy4xLys6d+P31fJiVdMk6v1EPPE+AJAC2IYgqoQblxCV7ZU7+fvk+Pbv5XUiyx1nTO7D4AaDkUU0ANyo1L8MqeioSCOnxgj/Ob6ZNONAKz+wCg5VBMATUodzCyV/bUsw/v1cHhQeeGl446s/qyubP7AABNjdN8QA0qiUs4ODy4VjzJOeF377GXiUoAgBZHMQXUqJq4hLKjEkzA6aHq3emc9KNBHQCaDtt8QJ2VE5dQVlSCJNkV0UMFAM2NYgqos3LiEkpGJZjcZnVJ9FABQJOimALqrJy4hJJRCTbt/eX0UAFA06GYAuqsnLiEklEJvmNlLBlUANBkKKaAOvOLS3hg5wMaPTWqoRND+q2/+If6xN982z8q4cGnnPEyXuifAoCmYqy1G/awkZERe+7cuQ17HtAo+af5Htj5gF788xdzeqnCwbBi+2OZU35npmZ0/OxFXU0ktaMvoi/8D9/WT/7Fv1gN8vTQu2stOR0AUHfGmDettSMl76OYAtbf6KlRzS7MFlwf6BnQ5COTBVEJkrPt9+zDe3XwxfdJ8vr/qXH6qwAA66LcYoptPmADlGpKz49KkKRkakXHz17075/y7asCAGwkiilgA5RqSs+PSnBdTSS9+6dCEec6AKDhKKaADVBqhp9XVELm+tAh6aHnnR4pGefX9/9dJ3Mq1sfpPgBoMHqmgA3i1ZT+ypVXNLcwpy2h7bp++UEt3nh/5v5Mz1TWPD9JTuH0tU/nDkYORZyCi3EzAFA3NKADTcwdOZN9ui9kNil4/ZCuzb1PO/oiOnxgT2EhJTkrUV4n/DjdBwB1VW4xVdOgY2PMY5J+Xs5RowuS/qG1dqn4pwB4jZxJ2Zu6Y+d/1p/88mckOVEJ9x57OROVkCmu/FLQSUcHgIaoumfKGDMo6dOSRqy190gKSvpEvV4MaGelTve5UQkziaSspJlEUkdOX9CZqRlO9wFAk6m1Ab1LUsQY0yVps6Srtb8S0P5Kne4rGpXgmY5unK0/mtEBYMNVXUxZa2ck/VNJfylpVtK8tXYy/z5jzKeMMeeMMefefvvt6t8UaCNep/u6TJeSy0kNnRhS4vbPqWvLVMHnriaSeaf7JMkoE+rJqBkA2HC1bPNtlfRRSXdL2iGpxxjz9/Lvs9Z+0Vo7Yq0d2b59e/VvCrSR6O6oYvtjGugZkJFRb3evjDFK3EzIyirQnVB44HRBQZWJUBg65DSb9+5SQTp6KunEJgAANkQt23x/S9J3rLVvW2tTkk5L2l+f1wLaX3R3VJOPTGr6k9PaHNqsVDqV83MTSGnT9rOZ30dCQR0+sCf3S2hGB4CGq6WY+ktJHzTGbDbGGEkPSvrT+rwW0Fn8GtIDoYRu+/Fxvee9n1eod0qPvfCW7j32stOILhVpOrf0TwHABqmlZ+qPJJ2SdF5OLEJA0hfr9F5AR/FrSJeRjJHUdUPp27+i4Jap3JN9ns3oq+ifAoANUdNpPmvt56y1P26tvcda+/ettTfr9WJAJ/FqSM+Xve2XOdlX0Iyeh/4pAFh3zOYDmkB+Q7ofE0pk/jozHNltRvf7HP1TALCuKKaAJpHdkD7QM+B5j031Zf66YDgy/VMA0BAUU0AT8tr2s+mQbr59QJLPyT76pwCgISimgCaUv+0XCbxHxnYrvOMFbXnv5/WJv/l24RBk+qcAoCGMtbb0XXUyMjJiz507t2HPA9pB/FJcsTdiOYORQ2aTgtcP6drc+3KHILtifSoI85QkGSmW8LgOAMhnjHnTWjtS6j5WpoAmN3F+IqeQkqSUvanFnq8VDkF2MQwZADYMxRTQ5PwCPbNP9mWiElwMQwaADUMxBTQ5v0DP7JN9UlZUgsQwZADYQBRTQJPzPtkXkAK3dNuPj6vnrx5T15apwqiEUsOQv/qLFFQAUAddjX4BAMVFd0clOb1TcwtzCgdu02J6UYGuRUmS6U4oPHBaoz/6o95f4BfaaVecFSrJKbwAAFVhZQpoAdmBnn2R22QCKzk/N4GUXr/+f3t/uFjTOXEJAFAziimgxfg1pF99d1b3Hns591SfVDzMU2LcDADUiGIKaDHFGtI9YxLcZnQT9P5C4hIAoCYUU0CL8WxIt05UQs9fPaZU5FxuTILkFFQf+23iEgBgHdCADrSY7Ib0q+/OSpKMcX7mNqN/f1aSPpT7QbfJ/KWjTgHlFZeQfR8AoCysTAEtyG1ID6xszRRSLhNIKfJXJr0/WCougWZ0AKgYxRTQwmyX95w9v+sZfk3nNKMDQMUopoAWNuDXjG6thn7nAf3ayxXGJZgAvVMAUCGKKaCFeTWjS04Ple26oa987znvgsovLsEN8qSgAoCyUUwBLSy6O6rY/pgGegYKWqAkp3/qP3znXxX+oFhcQiopnf4FTvgBQJkopoAW5zaje9RSkqR08Ib3D4YOSTbt/8UMRAaAslBMAW0isLK1ouuSSgd2csIPAEqimALaxN+5+xdk06Gca9Y6vVOjp0YVvxQv/FCpUTMSJ/wAoASKKaBNfO5Df18f/9HHZJa3ylpJdi3Mc3ZhVrE3YoUFlds71bvL/4sZNwMARVFMAW3kcx/6+5r+uVe047YBJ+A8y9LKkibOTxR+yA3yfPhfea9S3VqgbwoAiqCYAtrQ3MKc5/XZhVn/LT93lSqyLfd68jqN6ABQBMUU0Ia2hLb7/sx3y09yCqrunsLrNKIDgC+KKaAN3fz+gYJm9GxLK0t6/LXHvQsqRs0AQEUopoA2dG3ufVqafVjpW31OM7qHtE17r1AxagYAKkIxBbShHX0RLf9gWAt/MS6b6vO9z7MpnVEzAFARiimgDR0+sEeRkDMq5ubbxbf8CprVS42a+eovUlABQBaKKaANHRwe1LMP79Xg6grV0uzDstZ43tvf0194sdioGVaoACAHxRTQpg4OD+r18Q+tFVRXDxWsUIXMJo3tG/P+gmJhnZzuA4AMiimgzblbfu4KlduUnr7Vp6XZh5Wa/4D3B0uNmpm/zOoUAEgy1u+ozzoYGRmx586d27DnAXCcmZrRr578play/v/etWVKm7afVSCU0MBtAxrbN6bo7mjuB6dPOj1SdsX7i0MRp79q6NA6vj0ANIYx5k1r7Uip+1iZAjrAweFBpfMKqfDAaQW6E5IpMbvvY7/tv0JFQzoAUEwBnWJH31pBtGn7WZlAKufnvkGe7uk+PzSkA+hwFFNAh8iOSzChhOc9vkGeQ4ek3l3+X05DOoAORjEFdIjsuIRSQZ7jr44XDkQu2ZDOuBkAnYliCuggblzCbzw4rnAwXPTegj6qYmGeUvEoBQBoYxRTQAeK7o7qwzs+LdnifwQUjJsp1pB+a4G+KQAdiWIK6EBnpmb05f9nu5JXP1501IzkrFDlbPe5K1SRbbk3Jq/TiA6gI1FMAR3o+NmLSqZWCoI8/RQ0pQ8dkrp7Cm+kER1AB6KYAjrQ1UQy89fLPxjWwl+Ma+nqo76rVAXbfZJ/wzmN6AA6DMUU0IGyM6dcyz8Y1q3Zh30/M7cwl3vBr+HcBNjqA9BRKKaADpSdOZXt1g+GfWMT+nv6cy/4RSUQ4gmgw1BMAR3IzZwKGlPws6XvH5A8tvuSy0nvRnSvqATGzADoIBRTQIfKn9fncpvSe7t7c64nbib0xGtP6P4v36+hE0NOqOdtPZJNez+AFSoAHYJiCuhgXr1TkvQjgf3aHNpccH3ZLitxMyEruxbqub1IWCen+wB0AIopoIN59U5FQkEdPrCnsOHcw9LKkia29pUYM3OZ1SkAbY1iCuhg2fP6jKTBvoiefXivDg4PFjac+5hL/aD4mBmJ7T4Aba2r0S8AoLEODg/q4PBgwfWxfWOKvRHT0spS0c/39/Q7zeiSUzSlkoU3uQ3p0tq9ANAmKKYAeIrujkqSHn/tcaV9mszDwbDG9o05v3GLpNO/4P2FbkN69r0A0AbY5gOQcWZqRvcee1l3j8d177GXlZr/gJ657xmFg+GCe41MJhk9E5kwdEjq3eX/ABrSAbQhiikAkpxC6sjpC5pJJGUlzSSSOnL6glLzH1Bsf0wDPQMyMurt7lUoEJKVE6uQOdXnFlR+YZ4uGtIBtBmKKQCS1oYfZ0umVvSrJ7+p1PwHNPnIpKY/Oa3Noc1KpVM59y2tLOnx1x53CqpiYZ4uGtIBtBGKKQCScocfZ1uxVkdOX9CZqRlJHjP6VqVtem2FauiQ9LHf9l+hSiWlP/hsXd4bABqNYgqAJP8AT8lZoTp+9qIkjxl9WdweKklrK1S+X3qd1SkAbYFiCoAk/+HHLnflamzfmGdDumt2Ybb8hnSa0QG0AYopAJKKDz+W1lauorujiu2PKWD8//gYf3Vc93/5fqeoevAp/4fOX6npnQGgGVBMAcg4ODyof3bo/Z4rVIu3ljN9U9HdUd/IBFfiZsLpobqtR4ps876pt8hcPwBoERRTAHK4K1R9kVDO9RuLqZxGdHeFqphMD9XPfN67Gf3WAn1TAFoexRSAAgeHB9WzqXBAQnYjuuQUVAM9A0W/a25hbq0ZPX+FKnmdmAQALY9iCoAnv6iE/OulGtIzp/+GDkndPYU3kIoOoMVRTAHw5BeVEDAms9UnrW339Xb3et4/tzCnvSf2avTUqOLL170fRiM6gBZGMQXAk19UQn6Ip+QUVK/97Gs6dv+xgqIqZ+zM9m2K92wufFhka31fHgAZacMGAAAgAElEQVQ2EMUUAE/FohLye6dc0d1RbQ55FEurlozRxLa+wh8kr0ufv5veKQAtiWIKgK+Dw4NKW+v5M7+eKr9xM5mfB32CQWlGB9CiKKYAFOXXO+V3vdi4GUmyku6/c9B7u49mdAAtiGIKQFFevVORUFCHD+zxvL/U6T4Zo0QwqCe33+5dUNGMDqDFFAbJAECWg8ODkqTjZy/qaiKpHX0RHT6wJ3M9X3R3VJI0cX5Cswuzvt+bMkaPb7/d+czC4toPaEYH0GKM9emHWA8jIyP23LlzG/Y8AI03dGIoc6LPSzidVuza9bWCKtgtffQ3nVwqAGggY8yb1tqRUvfVtM1njOkzxpwyxvw3Y8yfGmP+ei3fB6C5nZma0b3HXtbd43Hde+zlnHgEP6V6qJYCAU1szTrht3JLOv0L0q9tk2K90nP30JQOoKnV2jM1Iek/Wmt/XNL7Jf1p7a8EoBmdmZrRkdMXNJNIykqaSSQL8qa8jO0bU5cp3lEw1+Vxws+uOL/OX+aUH4CmVnUxZYzZIukBSb8jSdbaW9baRL1eDEBzOX72opKplZxrfnlT2aK7o3r6vqd9E9IlqX95xfdnkjjlB6Cp1bIytVvS25L+L2PMlDHmS8aYgsFbxphPGWPOGWPOvf322zU8DkAjlTurz0t2Qnr+Sb9w2mrsRhn/HcYpPwBNqpZiqkvSPkn/p7V2WNKCpPH8m6y1X7TWjlhrR7Zv317D4wA0UqV5U1685viFQ5ul7ttKf7h3Z9nPAYCNVEsxdUXSFWvtH63+/pSc4gpAG6o0b6qYmys3M3+dWEkq1rvJO3MqW/I6fVMAmlLVxZS1dk7SZWOM+yfpg5L+v7q8FYCm487qG1xdiQoak+mZKudUn2vi/ISWVpZyri0Zo/Htt2t05w7/ourWgvTiP6agAtB0aj3N979J+nfGmGlJH5D0TO2vBKBZHRwezKxQraxm1JV7qs/lO7vPGM2GujS+/XbtvWuXd2G1ckv6g8/W8rcAAHVXUzFlrX1rtR9qyFp70Fp7o14vBqA5VXuqz1Uqd0rGZAqr2B3bCgsqtvsANBlm8wGoSC2n+qQyZvdlKQj0dBGTAKCJUEwBqEitp/rcE33l8gz0JCYBQBOhmAJQkXqc6ovujmqgZ6Cse61U2D9FTAKAJkIxBaAi9TrVV/Z2X37/VCgiPfhUta8PAHVHMQWgYn6n+n75hbc0fHSyrKLK3e4rd4VqKRDQxO23Sw89Lw0dqun9AaCeKKYAVMXrVJ8k3VhMlR2VEN0d1eQjk55jZrzMBo3TfB7rk567h1N9AJoCxRSAqhQ7vVdJVIJUwSqVtXo6+ENJVpq/LJ3+lPT1Xyn7OQCwHiimAFSl1Om9cqMSXO4qVdGCyhi9sOU9Wc3oVjr3u6xQAWgoiikAVfE61ZetkgHI2Uo2pheMnrHkTgFoKIopAFVxT/X1RUKeP1+8tVzR6T6Xu+UXMEX+eMo/4UfuFIAGopgCULWDw4N663Oj+sKjHygoqippRM8X3R3VM/eVHvWZSUg3hqZ0AA1DMQWgZgeHB9WzqavgeqWN6Nmiu6N6dM+jJe+b6wpKNq1MU/rXPk1BBWBDUUwBqItaZ/Z5eeKDT+jY/ceKbvn1L+fFM6SS9FAB2FAUUwDqotaZfX6iu6Oyq8GgBaxVMmA0dNeu3JEz85dreiYAVIJiCkBd1GNmn5/+nn7fnyWCQdn8hnQZtvoAbBiKKQB1kT2zz0ga7Ivo2Yf36uDwYM3f7RmXYK3TeJ4l05BOXAKADVTYMQoAVTo4PFiX4ilfdHdUkvTsHz2r+VvzzsW8Qso117W6OkZcAoANwsoUgHVzZmpG9x57WXePx3XvsZeriklwRXdHtTm0ueR9mYb03p1VPwsAKkExBWBdnJma0ZHTFzSTSMpKmkkkq86dcs0tzBW/wVrNdgU1umtQ8eGPVf0cAKgExRSAujszNaNfPflNJVO5sQXJ1Ip+9eQ3qy6oijWiZ3qojNFsV1Cxv4wrfile1XMAoBIUUwDqyl2RWvGJM1ixtuoVqrF9Y/4/zG9Gtyk9++oTFT8DACpFMQWgro6fvViwIpWv2mT06O6o+jb1lX3/vE0p/odPVvwcAKgExRSAuio38bzaZPTxnxovjEnwY4wmLn21qucAQLkopgDUVbmJ59Umo0d3RxXbH9NAz0BZ988GpPhvMgAZwPqhmAJQV15J6PlqTUaP7o5q8pFJHbv/WOmbjVFss1X87C9LX/+Vqp8JAH4opgDUlVcS+t/74J3rkozuhnmW4iSj90rnfpcVKgB1Z3wHiK6DkZERe+7cuQ17HoD2N3pqVLMLs6VvtFZ96bTGf5hSdOzb6/9iAFqeMeZNa+1IqftYmQKwIeqZhp7Nc26fF2OUCAb1ZO8mxSfeywoVgLphNh+AdedmT7mRCW4auqSat/vcrb6J8xOaffeq78w+V8oYTfR0Kfq1TzsXhg7V9HwAYGUKwLrzyp6qNmvKi9uQXryMWjPbFdTewTt0/7lfIyUdQM0opgCsO79MqWqzpvz0d5cZ6Lk6diYRDOjJ1/4JBRWAmlBMAVh3fplS1WZN+Rn74BGFTaiiz6TsisZfHdf7/837tffEXo2eGqW4AlARiikA684ve2rx1nLdGtGl1UDP+369MNCzjFPLaZuWJM0uzCr2RoyCCkDZiEYAsCHOTM0o9vv/VYlkKud6JBSsW+6Up+fu0d6tKtmYnm+gZ0CTj0yuzzsBaAlEIwBoKgeHB9WzqfAAcT0b0b3El69X9bm5hbk6vwmAdkUxBWDDbFQjeraJ27dVvCrlYqsPQDkopgBsGL+GcyvVNcgz21ywukLKyurJ15+koAJQEsUUgA1z+MAehQLexY0b5Fnvgqo/vxm9Aql0ShPnJ+r4NgDaEcUUgA1zcHhQt4X9By+sR/9U2eNmfMwuzGr01KiGTgwRmwDAE8UUgA2VWEwV/Xm9+6eiu6OK7Y9poGdAxlr1Lq8oUOEp5tmFWVlZzS7MavzVcfKoAOSgmAKwoUoFddY7yFNaGzcz/d0reu3yjJ55+x31Lq+Uzp8q8nPyqAC4KKYAbCi/AE/JyZw6fGDP+j28d6ckKbqwqNcuz+jY2+8onE7n3mPt2v9KWFpZoqcKAMUUgI11cHhQzz68V4OrK1DB1diCwb7I+oZ3StKDT0mhtZWv6MKiYteur61SWauItc4fjKvz+0ohjwqAfycoAKyTg8ODmaLpzNSMjp+9qKuJZKb5fN0KqqFDzq8vHZXmr0jG+e/Jm4G1wilZYSZVf09/XV8RQOthZQpAw5yZmtGR0xc0k0jKyolH+OUX3tLw0cl1yZyS5BRUj31LiiWkn/gHmtjap6VAdX8UhoNhje0bq/MLAmg1FFMAGub42YtKplYKrt9YTK1L5lSO6ZPSN/+95rq8+7dKCZiAYvtjiu6O1vnFALQaiikADVMsBmG9Z/bppaNSKqn+5cJiriRr9fHkiqLvLtT/vQC0HIopAA1TKgZhPWf2af6KJGnsRqLwRF8pxujFbqv4fz7srHAB6GgUUwAaplhMgrQ+mVMZWTEJsWvXNZBaLjsSQZKWAgFNbNnsrHAB6GgUUwAaxo1J6IuECn627plTWTEJ0YVFTV65qgvfvaxj5QZ6SprtCmpoq0hDBzqcsRWOVajFyMiIPXfu3IY9D0DryI5I2NEX0eEDe9Y3c0pytui++ouSLeybivds1sTWPs11BWUkpcuMTBjoGdDYvjEa04E2YIx501o7UvI+iikAzaIhBVWsT1LxPwfjPZsVu2NbxREKFFZAayu3mCK0E0BTcDOn3KiEmURSR05fkLSOIZ6S0zs1f7noLdGFRUnSxNY+zbpRCmWsVLnz+yRRUAFtjJ4pAE3BK3Nq3eMRpIIRM17cLb9KCilXPef3xS/FNXpqVEMnhujTApoIxRSApuAXgzCTSK5/IvpDz0u9uzx/7G7xzYa6yp7Xl68e8/vil+KKvRHT7MKsrGxm1YuCCmg8iikATaFYDMKNxdT6jplxR8x4FFS1jJtx1WN+38T5CS2tLOVcq+eqF4DqUUwBaAqlMqekDRgz47HlV+24mWzJ5WTNK0h+q1v1WPUCUBuKKQBNwc2cKmVd+6hytvyM1LurunEzeRI3EzVvyfmtbtVj1QtAbSimADSNg8ODGiwj9Xxdx8y4W36xhPTYtzR2M1j5uBkPtW7Jje0bUzgYzrkWDoY1tm+s1lcDUCOKKQBN5fCBPQoFijd5r+uYmTzR+59S7NqNtXEzNahlSy66O6rY/pgGegZkZDTQM6DY/pgkccIPaDBCOwE0nTNTM3rs5FuetYuR9NyjH1j/MM9sX/8V6dzvanTngHOqr0q93b3aHNqsuYU59ff01xzo6Z7wy25MDwfDiu2PkWsF1EG5oZ2sTAFoOgeHB31Dya3WOcTTy4f/ufTwFzW2sFzTlt/i8mJOtMH4q+N6+htPV/19nPADmgPFFICm5LeVV05P1boYOqRoepNi165rILUsY61MhSv7qXSq4NoLF1+oamsufimu2YVZz59xwg/YWBRTAJqSV1RCJBTU4QN7GvRGkuYvK7qwqMkrVzX93eIjaCpR6UqSu73nhxN+wMZiNh+ApuRu5W344ONiTFCya1EJ/csr3j1U7opVmWnpla4keW3vuTjhB2w8iikATevg8GBji6d8NjdzauxGQrE7tuUkpIfTaYWtVSJYfthnpStJxYovms+Bjcc2H4CWcGZqRvcee1l3j8d177GX1y8FvZi8cTPRhcWcHqqB1LJi165rvsLxM5WuJPkVXwM9AxRSQANQTAFoememZnTk9AXNJJKycoYfr+tYGT8e42aye6jGbiQ0sbXP7yCir0oLIAI8geZCMQWg6R0/e1HJVO4W27qOlfHjjpsxhVt48Z7Nit2xzemhKrNXSpICJlDxaT6/AE9WpYDGoGcKQNPzGx+zrmNl/Awdcn792qel1NrzJ7b15fROlStt05mTeZUUQ9Hd0Zz745fiGj01mgkEfWDnA3rlyiuaXZhVwASUtmkN9AzUHBQKoBDFFICmt6MvohmPwmkjx8rkcAuql45K81ek0GbNVdBwnm9pZUmPv/a4pLWCKn4pronzE2Wlpecnoc8uzOqFiy9kfp626cz1ago3AMUxTgZA03N7prK3+iKhoJ59eG/jT/tNn5ROf6rmUTOS1GW6dFv3bUrcTBT8rNiYmNFTo74Bnl4GegY0+chkTe8KdIINGydjjAkaY6aMMV+v9bsAwMvB4UE9+/BeDfZFZOSkoDdFISVJf/BZSVZjNxKFo2Yq/I/VZbvsWUhJxcfEVJpTRUI6UF/12OYbk/SnkrbU4bsAwFPTZU5JzqpU8rok51SfJE1s7dNcV1D9yyv60dQtfSMSqaghvRi/Iqi/p7+ilSljjIZODNVl2DKAGlemjDE7JUUlfak+rwMALWL6pPTVX8y5lB2TMHnlqr4X6q5bISU5RZPbaD50Ykijp0YVvxT3jEooJm3TmWHLsTdiVc0GBLCm1m2+L0j6jKTqx6gDQKuZPumc5stLRM8311V9U3q+cDCsB3Y+oNgbMc0uzOYUQ5IKohLK5W4fukXa3hN79f5/837tPbE3U6wBKK7qBnRjzIcl/U/W2l8yxvwNSf+7tfbDHvd9StKnJOnOO+/8ie9973s1vC4ANIHn7pHmSw86Ht25w392X4UrVn2b+qTlm0qsFJ5q9Goor7QpPRwMe877K9b4DrS7jWhAv1fSR4wx35X0ZUkfMsb82/ybrLVftNaOWGtHtm/fXsPjAHS6phgpIzlxCGXwakoPp9PaXMV/xCZuJpRYXvT8mVcvVSVbfwET8B2cXKzxHYCj6mLKWnvEWrvTWnuXpE9Ietla+/fq9mYAkKVpRspIUu/Osm7zmt330VtGtwJVnv3xWc3ymtWXn5Le292rgMcf+aFAKJND5YfTf0BxhHYCaAl+I2UeO/mWJG3sSb8HnypIQPcTXVjMnPRTKKLRu39My6n56p+dt0UYTqdzZvL5hX2OnhrV/K3C527u2qzNoc1FtwS3dHNYGyimLsWUtfYPJf1hPb4LALz4jY6xVjp86puSNrCgyk9A790pvXdUOvc7xT+XSmruVqKmE369K842oRu/MHYjoeiLn5UeXFD8tp6CJHS3Qd2vWHILrIACSvucJVpcXlT8Upy+KcAHCegAWsK9x172HCnjGuyL6PXxD23gG3koozHdtym9HFl/Xvel0xp/50bBqtesx6pXb3ev56pUJcpNTa9kDA7Q7DYsAR0ANsLhA3uK/rwhQ4/zPfiUFCo+L3DsRqLiZPQMYzL/SwSDenL77Yr3bHZ+5q56eai1kJLK65tyZwTmRzcQr4B2RzEFoCUcHB7U1s0h358HjGn8Kb+hQ9JDz0u9u3xvyawk1UHKGI1vv11779ql0Z071Js/zqaOyumbmjg/UXAqkNOA6AQUUwBaxuceep9CAe9+oxVrG3/KT3IKqse+JcXmpcg2z1sGlouHfVZkdaVqNtSldwMBhQK5BWc4GHYyqmrk9k0V47d6xWlAtDuKKQAt4+DwoI5//P3qi6wVDF61VTK1ouNnL27gm2WZPun0TsX6MnP78o3dSChs/FfZqrVsjFIrt9SXtjJy+pxi+2Ma/6nxisbNeEmlUyVXmLwiGopdB9oF0QgAWkr+wOO7x71XSxrSQ+WOmSkRmRBdWJSuXdfElohmu4IKSEobo4AJKG3TMjJy1tmqYIyWbFrPvvNDRYeflLKavyfOT1SUip7Pa4Upu+F8S/cWhQIhpdKpzM/DwXBOdAPQjliZAtDS+nz6qHb0FW8EXxcvHS0re0qSoj+44axQWav0alRC2qYVDoZ1aM8hdZnq/1t3KRDQxJbNzvu4z9sd1eQjk7rwyQs6dv+xiub3ufJXmPIbzudvzctam9lWdJPVH3/tcWb9oa1RTAFoWWemZvTu0nLB9VDQlDz9ty7KHDPjmtjap6VA7h/DSytLeuXKK3r6vqdlVH0e1VxX0Pd9orujFY2bcSWXkznFkFfD+bJdlrVW4WA4k6zu/srpPrQriikALevXvvZflUoXbof1dHdtbCK6q8wxM665rqD39XevKvqV/1XPfv9a1TEK/csrRd/HqxDKFglG1Nvdm3MtcTOhx199XPd/+X4NnRgqGgTKrD90EoopAC3pzNSMbiymPH82n/S+vu68cqZCEd9Tff0+p/r6l1ek5HVFFxb0wWSy4oIqnE5r7AeLzvv4KHbCrst06XP7PyfjkdSeVlqJm4nqe7pKPBtoRRRTAFpSsdN6DemXkvJypozz60PPSz/zec8wz7EbCYXzsqHC6bQT7Lnqe6Hu8sfPrA5Tjl27rujfOr429saD3wm7gAno6fueVnR3VImb3iGgxZQTxcDpPrQbiikALanYab2G9Eu5MjlTCefXoUNrRZbJ3daLLiwqdu26BlLLMtmFUFawp99WoJ+xGwlp8+0a/bMvaejEkG/Tt1fPVDgY1sf/2sc1cX5Ce0/srei5ro/+2EeLRjHUcrovfimu0VOjRf++gEZgNh+AllRsVt9gX0SHD+xpTN9UMdMnpTO/JKXL34asdJZfVzotY4xSWatZ4WBYsZ0/rejUV53ZgSYo2RXFt+/Ss1s2az7t/HOMBCNatss50QbVGOgZ0AM7H9ArV17R7MJsJvJhoGcgM6uv0hl+7snB7F6scDCs2P4Ys/+wbsqdzUcxBaAlnZma0ZHTF5RMefcdRUJBPfvw3uYsqP7gs76BnvniPZsVu+N2Lfkkv5erd3lFr13OTYWv13cX4xZSvd29MsZo/ua8tnRv0eLyYkEeVbHCaPTUqGfDe7kDmIFqMOgYQFs7ODyoZx/eq0Gf/qiGpqAXM3RI+ux3nHEzsfmic/wkKZoOK3b3x6rKhco2HwysDUVe5UQzrF8hJa3FIszfms80rs/fmi9Y/Sp1yo9RNWhmFFMAWtbB4UG9Pv4h3zSmhqSgV6pUNlV3j6J/49c1+chkbQWVMZrYmtsYXmk/1nqbW5jz7YtiVA2aGcUUgJbnd3ovYEzjBh57yZ7b9/m7nf+VihiYv+x8ZvpkVUGb2fKLJ79ohkbZ0r0lJ1E9O+TTr2GeUTVoBhRTAFre4QN7FAkVrrKsWKsjpy80R0Hlzu2bvyzJOj1TZfZNaf6y9OI/VvTdBcX2x6peodqyktbozh0aumuXRnfu0AOLiwp7hJ42ilfYp7v9F90dzfy9G5nMEGeaz9EMaEAH0BbOTM3oV09+Uysef6YN9kX0+viHGvBWWZ67Z7WQqkFkm9NvJWnoxFBlwZnWKiTlnvKzVh+9Y0SvLF2taQDyejMymv7kdKNfAx2IBnQAHeXg8KDSPv9x2BS9UxXO7fOUvO6scKm6XqFUXvjnkjF6ZemqJh+ZrGkO4HqjLwrNjmIKQNvo2xzyvN4UvVMVzu1zxXs252zNxf/z4bX+KeP991uJudUVKb+CZaBnoGSi+XoKB8N6YOcDhHWiqVFMAWgLZ6Zm9O7SsufPmqJ3ymtuXwlODtQ2zYa6ZI3RbKhLsa23Kf7qUad/6to7CpTTqmGt70ia/hXn88UavDeyHSRbwAT00R/7qF788xc9m9KBZkHPFIC2UCwRPVtD09GnT0ovHXW2/EKbpdRC0dv90s8HUsuanPnvkl3R0F27ZP1m91mrgKS05F1MWatjb19X9PBVSfJNJa+4P6uO3NBPr+vW2rLS04FqkYAOoKPcPR4v+1/3TZOO/vVfkc79rvziEfwKJWOtpr/rNLP7FVy93b26uXKz4HRczj1uKnrvLmflzGcwsl/6eLPoMl0KBUJKrjjFdN+mPo3/1DgFFmpGAzqAjuKXNeWladLRvz2pYjlTfjlQ2dfHbiQUTueu3ISDYRljihZS4XRaR67fcH4zf9mJbVhtbs/X7FlOy3Y5U0hJUuJmQk++/iRbgdgwFFMA2oJf1pSfVjjh51kopdMau5HI/D66sKjYtesaSC3LWKuBUK9i+2Oavznv+73GWsWuXVd0YXHtYirpbEF6KLXC04wnAVPpVNHxNEA9UUwBaAvurL6gX/9QnkpWstZNiRN+BYVSarmwCFq9b/KHQU3ve0qTf/c1RXdHi8YJ+K6FFSnu/IJCqwkQdYM3KxEwlf/rirl92CgUUwDaxsHhQf2zQ+8vuUJl5KxkNdyDT0kliorowqImr1zV9Hcva/LK1YJCSpIT5vnYt3J6nopuzRmj8e23O1EL2cOP84q77Dl5i6lFBTz+lTG7MCtTZgGb/W6VZkc9c98zFRdu5FNho9CADqDtnJma0fGzF4ue7vvusWjmvquJpHY06pTfiY9I3/l/a/+eyDbpZz6fU1A9/Y2n9cLFF4p+LJxO5612Gal7s+IhKbZ9m5YqLJRK+WD/B/W9H36v4oZ2I1PRicKAAtqyaYvmb85z4g9VowEdQMc6ODyo18c/pEGfrbzBvojOTM3oyOkLmkkkZSXNJJIbn0U1fVK68sd5F4109//oFEeVSF4vaCJ/4oNPlNxOWwoENLE1O5TTSrcWNLG1t+JCyt2Kc3/t7e5V36a+zCy9R/c8qrfefquqk4GVFFKRYETBQFCJmwmyqbAhCs/TAkCbOHxgj46cvqBkau30WyQU1OEDe3T87MWc69LaKb8NW5166ajT+J3DStcvSd095Q9CXhXvNpo4d1RzU09nVmPKKUJmu4IaumuX+pdXNHYjoejCoua6ym/mz7y5tbrwyQuF77WaX1Vqlawe3K3A/ILNHZgsyTNLK/s9vX4GFEMxBaBtuUWR11beYy+85fmZDT3l59fwPX9FxSITvLhp6UsBI2WtxvRt6lPiZqL4h42RlZyE9TucFbH+5RXP/KpivHqU4pfiir0RKxrTUE/Fms7dfybuu7i/d/n9jIIKpdAzBaBjZPdIBYzRiseff4N9Eb0+/qGNeaHn7nEynvL17pJ+cFWy3jlTXnzDO5dXdDMYqGjLbiC1rLEbidXiLKsbpMhYmnAwrNj+WEHhsdGBn34rU5J/mnqxnw30DGjykcn6viRaBj1TAJAlv0fKq5BytwA3jNe8vlDEuV5BISXJd1vuB8GAYm+/o4HlFRlJvYGI+lbSTmHk8x/Tc13BtViG5bVYhkd/8MNMTENvIJLTD+VVSEkbH08wtm/Md86gXyElyfdnxCugHGzzAWhb2StRMt61Q9AYpa1tzGk+9+SdO6+vd+faWJeXjnqvWvnw25brX15RdGHROa3Xu0t67I+l6ZMaPXdUs0HvVSYr6f47BzUfCKh/xerZt99ZO+13fXXLMLJN+ux3cj7n1XPU39O/oStT46+OZ/7aPQE40DOgsX1jmjg/UfG7EK+AcrAyBaAt5a9E+XU0pK3Vd45F9fr4hxozq2/okJMRFUvkZkV5rVoFQpLxXoEqJy0906M1dEhzwSJ//BujRDAoa4xmuwKK3bEtN49Kcprjs04Our1RswuzOSfoHtj5QMEqUb5SP6+WlVUoEMo0knutWJV6r2YfpYPmQDEFoC15ndbz0hRJ6F6GDkkPPe+sJsk4v256j+/2X1lp6VmhnJWsuBTGJ6zKGj8zcX6ioMl8aWVJr1x5RbH9sUwvU358grtFWE2Sejmyx8pEd0dLPitgAiW3LoF8bPMBaEvlnspriiR0P0OHnP9Nn3QKlxJRCZntPC9uL5acVaTFlM99Pjx7srJOI/r1Fs0tzCm6O1pWUbJep/5mF2Y1dGIos/U4+cik5ylDvyb6bMQnwAvFFIC2tKMvUjQBXZL6IqHGbO1VYvqkE8ZZkEdVghv6mbyR04tVbVRB/7LHiljeSpdXP1K5K2BuQVJuX5PbD1XshF627K3H/OeVWxjl/7MjPgEutvkAtKXDB/YUndEXCQUV+8j7NvCNquQZ7FmG7h5nvExeL5bXdpwk/6YySeG01dj8u3lXjfTe0czv/E7QuT1H2XP+Rk+N6ulvPJ3z+/iluKK7o5p8ZO1TndYAACAASURBVLKsLb/+nn5d+OQFffyvfbzkvdmWVpY0/uq4Rk857z75yKSmPzmtyUcmSxZEfluZ7jYiOhcrUwDaUn5gZ28kJGOkxGKqcXP4quEX7JkRlOSxajR/2VnRknLm9RU76j+wvKy5YFBbVtIyRs5pvkwqukdS+zf/vXTnB6WhQ0VXerxWdLLT0PNXeMqJI5hbmFP8Ulwv/vmLJe/1Us2qUrGtTHQ2QjsBoJkVC/Z87FvSr20rnknl3rfKL0RzILWsyStXJTkJ6hXp3bUW6eCh3ODO3u5evfazr5V1f7FwTql4QGf+fc/c90xZBZXvPzuCPdsWoZ0A0A6KBXtKpcM981a2PLfjciIUqvgPbHcVLCsqIVu5Kzfzt+YVvxQvGWHgbh8W+95yCin3vnKHIJfaykTnopgCgGbmFZHw0PNrq0A+uVMZWU3iUlY8wIr1j1CoRiqZE5WQrZIYhonzEzkRBkZGvd29nmnr9QrULLfvKf+9iE+Ai54pAB0nOxm9pfqnvPzEP5DO/Y7/z90VrCzR3VHpL7+hiT8/pbmuYCZDyq+givds1sTWPs11BbN6qDzu9envGts3VvYJwuwYgwd2PqBXrryS6cEa/6nxnMJlbN+YnnjtCS3b5ZLfW85zR0+Nam5hTlu6t8gYo/mb8wWn/MqNeUBnoZgC0FHcZHQ30HMmkdSR0xckqTkLqvxohPzG8g//c2n6y9KthcLPRrZ59jHFL8UVu/IftbQ6fmY21KXYHU6UQn6RFO/ZnDPwuNi9+atgLq/m9Pmb81pc9i7e3BiDYk3qLmNM0Z3Jcnun3GdIznZjqecC2djmA9BRvJLRk6kVHT97sUFvVIJXNEL+ltqHv+DdV/Uzn/f8Ss8j/j4p5xNb+zKFVM6927YWPs9jFczlxh64MQRP/fWnFK7wAFT+dtzE+Qml0inf+8PBsJ657xkdu/9YRc8p9VwgHytTADqKXzJ6uYnpG84vGiH7es7A5MtOH1V2wZW3OuV7xN8j5dwz+VxyZvuZoNMAX+I0n5fo7qj0H35BE1t7NdcVdBaXjPfgZb93L9aA3repL2dbMHsAcjWKbT9mbwOSkN6ZKKYAdBS/ZPSmndHXu9MnGiFvS80tZIptCa7yTSvPTznv3aX+23b432tX1lakKiikXNGubYpecf7eRnfu0Gyo9L+SspvO/f4+JGlpeW3lrZyTeuUotf0oiYT0DsU2H4CO4pWMHgkFm3dGX6lohGzlbAnKLx7BZsUjrD2jZJRCkVN8JWX9vY3dSCicLt3blB1DUCxCIXtrbr236NxnkZDeuViZAtBR/JLRH3vhLR0/e7H5TvblbOFdyZmzV8B3SzB3ZSvTEP6NZzV3K+Gc0FtYVjQdlpTMeUb0678izV3VRN8W/9N8Xitn0tqAZve93zsq/devrg1sjmyT3v93pW9PKjp/RdpsNLG1V7Opee/vy3r3+KW4nv2jZ4ueEJxdmNV9v3dfTkP5eim25UhCevsjAR1Ax8o/2Sc5q1TPPry3uQqqcvmlpctID38xtwDzGqAciuRmWH39V4rHLlT6/V6C3dJHfzPns6WSxuOX4nWLRKiXYonsJKS3LhLQAaCEljvZV8qDT8kZB5PPFm7FlbMl+Oa/LvPBZX6/l5Vb5W1DZiWNT5yfaKpCyn03EtI7F8UUgI7Vcif7Shk6JN/QpfnLzmrR9MkiK1jK3SosNarG63PTJ6XP3+3//aWeqdJJ46W2zXq7e8t/dh2Eu5wCyiu53cho/NVx7T2xV/d/+f66NcOjudAzBaBjtdzJvnL07vIvZM78khM/sHKryOdXTwn6zNkr+rnpk84zimQ/FX1mlmJJ48VO8fV29+rmys3Knl+jxM1Ezqm96O6o51Zk4mZCT77+ZOY+tA9WpgB0rJY72VcOr9N/rnSqeCGVfUqwohN6xvncS0crL6SC3UXDPr2M7RtTlylcCwgFQjLGFG1KD6zTv/a8AkW9tiJT6RSn+9oQxRSAjnVweFDPPrxXg30RGUmDfZHWbT53uYORK5U/QNnvZKAn63yu1GdGfs45weeKbCtoPi9HdHdUT9/3dM52Xt+mPv36vb+u+Zv+J/f6NvXpmfufqehZlSg3UJTTfe2HbT4AHcVryPHr4x9q9GvVn5tOXo7INumxb+Ve8wsL9XvW9Mnin4lsk749KSVvVJWYns9vG/DYHx9T4mai4Hpvd69e/cSril+KlzWvr9jpPD9burf4nkTMlh08ivbAyhSAjuFGIcwkkrJaG3J8Zmqm0a9WP24kQSXN47feLeyRKrZdmM+uOM9876gUCBX+3ASdZ8xflmTXktkr7csqIX4p7rsy9W7qXT39jacVeyNW1uDjB3Y+4Lud6CWggBaXF0sWUqFAiNN9bYhiCkDHaLsoBC/lRhJk84gnyGwX9u4q7ztSSWfl6eBvFW7lhXsLe7VqSU73MXF+QtbnNOOKXdFX/uwrRfupsr1y5RVFd0d1W/dtZd1vZYsOXZbWtiJpPm8/bPMB6BhtF4XgpaJepxKfGzrk/C/WJ9/IhfzvcD+TLdZX/jNrUKoXqZwVqfzvKtaDlc2viHMZGb36iVfLfj5aC8UUgI5RKgrBq5/q4PCg7/WmVEmvU/7nav5OK8V6ndWo933MWamavyKZgPe2Y2TrauZViTE5ZSoWmSCprF6p7O+KX4rLGKN6TAqxsho9NZoTPDq3MKf+nn796Ht+VH/83/9YaZtWwAT08b/2cT3xwSckOVuX2feO7Rsre2Wrls+iMoyTAdAxio2PkeT5s7/zE4P6D2/OtM7ImXLHuGTLHyNTj+8sJRAqzLwq9R4lFBszEwqE9PB7H9aLf/5iya2+Su6tVJfpkjGm5Jbgo3se1fCPDCv2RiznHcLBcE6AqZ/4pXjVn8UaxskAQJ5iUQh+/VS/90eXW6vPKqfXyaz2L3mNmFmVH4lQznf27pK6eyp/NxNc+/ym99S9j8orMkFa61V64oNPFCSrP7rnUc+IhVeuvOJZSAVMQI/ueTRz2q9Sy3a5ZCElSV/5s69o4vxEwTvk51n5qeWzqBwrUwAg6e7xeDldQRlG0neOtch/4fv2PBkpVhgjUNt3FpP1vPV4pzoaOjHk2QdlZDT9yemi99SLkSn5Dn7KeX+UVu7KFD1TADpSfh9U3+aQbiwWrhgYI3n9N2dLjZzx63kq1idV7XcWYwJOEdW70+mXSl4v+53q3Tskqej3+fVfbenekvnOevVT+fEr1NycqmL/TPzen4yr9cE2H4CO45U39e7SskLB3O2wUMB4bpCFgqa1Rs54ZUZlj46pxntHK/+MXVEmZ+rmD51RMmW8k9v/M7swKyur2YVZxd6IlTU02OuzT7z2hJ58/cmi3+eXMbW4vFhRXlW9hYNhje0bK/nPZGzfmMLBsOdnUX8UUwA6jld/VCpt1dPdldNPdVu4S2mPxYGe7q7mbD7349XzVEOjtyTnpF4t0imp+7ay3qnevUNefUv53+eXMZVKpyrKq8rnXZ6XJ2ACmQbyUv9MorujBf1hNJ+vH7b5AHQcv1yp+WRKb31ubcXl7nHvlY/5ZIXDfJuBV/5TLeqREZW8IX32OyVv88uPKmfGXSVz8PLv9cuYqmVFqpYeK2ttphgq55+J38gd1B8rUwA6jl+/U/71cu/rSJX0W5lgTd/h1+dTTv9PJT1C+ff6fTZgqv9XZy2fzX6fWv6ZoP4opgB0nMMH9igSyv0XvJHTO/WBX5vU8NFJ3T0e18JN7z6qxVvLuns8rnuPvdxec/3KNX1SurVQ3r2hiHTXfSqMZzDefVfTJ50gz1if8+v0yZr6f7w+22W6FMqbIej1fV6flaRNgU0Fn/f6Ti9eny3XAzsfKPlu7j3xS3GNnhrV0IkhjZ4aVfz/b+/+g+MozzuAf5/bW1krOZashCBZkkOYBKcxvwwefjloWruxIcYxMJgEJpmh7YzbIW0MdNxAh3GEQ8ekbkPETNoZBmiSSSAYQxTAJHbGtKWQMQVsbOMQE+oELNkKECMBRrKl09s/9lba29vd2733tHvSfT8zjKTb99fKWHp8+7zPe2hb5NcoPpZGIKKa5Ozm6x8chiB4k7+ZEcyuz2Lww1E0WSaOnxzDaG6ydVUX8JwKcQp4OpXQ9z7o395bpNNv7HybbbMbE93N5+676flNGDpZ+MgvK1nMrpuNoRNDBWP6tfVyAq/hXLwiqN6im3fuuhMPH3y4qM3qT60uKjjqVyzU7zUW9iwUtTRC2cGUiHQC+CGAVgDjAO5VSoVmAzKYIqJqs+Sup32PmHFrb7bw3K1LA9s612vC3WdGL4ngHJIc1r6pE7j5lfCx3W1SsHzrct8yA22NbdhxzY5Ibb3CjrZxCoKWmjNorjjH5gTN772vWpVEnakxAH+vlNotIh8B8JKI/FIp9WuNMYmIEhXlkGOnTU0clFxKnMTzKG3dbYLaV/hA5LjiJMBHTXgPC3bCxnBfC2qnW7IhTtI+2crOmVJKHVVK7c5//j6AVwHUyPvcRDRTREkmd9owIR3xEs+bOkq3d18PaqtTXLQC4iR7R00AD0tEb21sjTTnVCTIh41LwSpSGkFETgOwCMDzlRiPiCgp61cswPqtewvyoNws05go0Ll+xQLfw5CrvoDnvi32mXdDfXZg4hTG3LnRfqwmhl1Qs6nTvubOYfL2W7bBJ2cqAzvbw6PldOCjnwJefACBR8cMHbYf7wWOjcI2lSzvUIKTa+X3KK3eqEdXRxc+99DnSuZI+Ql698jMmBP5V96DigFgcGQQl/7kUgydGMKcujkwM2ZRzayo70xlJQsFhZya/P/ZPT9Fpx1MichsAI8CuEkp9Z7P9bUA1gLA/PnzdacjIqq4XEAgNbfBxDdXLZxILnc+uo+hWb9iQXUnn3uTuocOA7032ufkOAcNO79Mhw7bbR3efk983U4YX3XPZJBlzQVOflB8aDEA/O6/gTeeRWAg5bzuHvuc6/2DL/faEgionArjfsU52xrb0NXRhUdfexRjaqzkWBlkMO4XbHo0ZBuw4eINBcnf3oT24dzwROL60MmhsouAuu/BLclNaTOJ1m4+ETEBPAlgu1LqO6XaMwGdiKpNWAL6jEgsj5Mw7ghLHPcmg5czvvNOWJx5g+afIqWSzqMmmgPRE8J1EtrjKHUPTECfNOUJ6CIiAO4H8GqUQIqIqBqFJY/PiMTycpK3w/p4r5Uzvl8gFXWshJLRSyWdx0nSjvrYTSehPY5S98AE9Ph0stSWAPgqgKUi8nL+vy9UaF1ERIkISx6fEYnl5SRvhyWOe18vZ/ywiuhxEtanUKkE8DhJ2lETwnUS2uModQ9MQI9PZzffs0opUUqdrZQ6N//fU5VcHBHRVFu/YgHMTHHeiWlI9SeWR7Fsg1340i1jAkadf3vTmkwG9/ZzrpUav2AuT+BkWsD5NwSPHTae3/xTpFTV9XXnrUNWSj/cqTfqseaMNSXbBlV0D6p07girvF6q0rtOZXkqxAroRFTzevf0o/vxAxjMH2DcYGYwyzQw+OGob5K5Uz19WiWhx9nNF3bNm/y9bwvw828Aw8cCJnclmlstwOXftj9/4iZgNH8kjWSA0y4Fjh3K50v51KR3+p59bfGc7muagnbwNc9qxq0X3FqQHB5UHV0gUFBoa2ybqKwe1Baw37m64NQLsPftvQVV0RuyDfhw7MOy7sNZLwDf+3HW2FTXBBEpqOTO6ueTkijaSUQ0I1y5qH0iGOrd04/bHtuPdz+0A6v+wWHc9tj+iXbOdac8gvd6VTr7Wv9Aw+817+4/lZt8R8gvkPKWMsiYhTsF3UHR2DDw5i5g9w8B93Z+NW7v/Jt8oXhdY8OTc/beWNh/+Bjws68F31NEYTv4RsaKX3OCDm+fWcasoiNZgtoCdk7VroFdReOXG0i51+vMe/uztxfsPFT57/HQySHUG/XYdOkmBlEaeNAxEZHL5u0HC+pIAcDwaA6btx+MdH3a27mxuM7T6LD9epS246P+ZRKccV76fmEgFJWzhp0b/fvnTvqvMYae3T2+gRQAjORG0LO7+MQ0vz5x2k4V9xp6dveElnAIWi9Fx3emiIhcSh0ZM+OPlIlzpEsld/JFUWo+zZ1+pXaxxdltl9TOvDBxdh5yB58evjNFRORS6siYGX+kTJwjXSq5ky+KUrv9NHf6ldrFFme3XVI788LE2XnIHXx6GEwREbmsX7EAlln4C997pEzY9Wkv6i6+oLaldgqef4PdJi73bj+//kad9k6/sJ1zcXbblbszr5Li7DzkDj59fMxHRORS6siYaXmkTBxOArd3919YAnvcc//mX1S8G2/hVcBvdxT2CdtJOAW7+ZwEbGf3m1O53L0rL6zPwPGB0B1xYeN3dXTh8dcfL3s3n2VYmJWd5bsrz/no3k3ot+OQysfSCEREREQ+WBqBiEiDt/aU99Bjp43zDlWTZUIEgbWppgW/elROXaco71RFGct7zZoLjJ2YrDnlrSflfYcrsCZWfiwAGH43+jojcupPlXr3Kayv37tdQPC7WlHnDKqN5cevXpZ3rjl1cyAiGDwxGOndOeI7U0RERXr39GP9I3sxOl7489E0BJuvOce33pSXZRrYdPVZ0yeg8qsZZVrAOdcDex8sfn3VPcGBStBYq+6xP/de8zLqgEVfLZ7Xq6imlUepdUbkV3+q3qgvqiUVta8jK1mICEZdpR6ccYHimlR+c4aNH8TMmPjWkm8VBG1Rxoh6zzNJ1HemGEwREXksuetp9AeUOmhvtvDcrUtD23jbTgt3n5mvPu7hvBPk1dQJ3PxKvLGaOu2PfteizhtX2DojWr51ue+7Pm2NbdhxzY6y+oZpa2wDgEhzljO+d5w4Y0S555mEj/mIiMoUVjOqVL2pqONUnaAaTUEBTVhNpzi1qoJUIpCKO2eAOLWkymkTp4/3Wrn1odz94ozBelT+WBqBiMgjrGZUqXpTUcepOkE1moLqQpVT76lUnago88alWXsKiFdLqpw2fn2izllufSh3vzhjsB6VPwZTREQe61csgJmRotdNQ0LrTblNu9pTQfWlzr8het2pUmM5daK817yMOv95vcJqWkVZZ0RxaklF6evIShamp2aWM27UOcupXWVmzIJxoo7BelTB+JiPiMjDSRoP283nrTc17XfzhdWXmn9RvN18UWpVRdnNNzFvurv54tSSCusbdzdflDm945fit5vPe3/czRcfE9CJiIiIfERNQOdjPiIiIiINfMxHROTDXZBzXrOFP/vMKXhy79Gix35A+OPAmrdvS/SjX568BXjp+5OP886/wf9Rn/sRX9DYUQqNhrSJUjDTr9Cl33EucYTNG7dwaFj7O3fdiUdeewTjahwZyeCCUy/AG++/UXQv3vvq6ujCM33PxHrcqVPwdLrgYz4iIo9SBTkdTo66p7ZnQXHPmrZvC9B7I+AqSgnAThpf/b3C4ObJW4AX7y8eI2MA4z5/DhkTUOPFJRSCCn56C3iGFBbdNruxZMHMUoUuyylwGVYcFIhWxDPKWHve2oOHDz4ceV1hSt2nTsHTasCinUREZYpSkLOUaVWwc6oEFe8Eiotp3tFSudpSUQqNhhQWXd45r2TBzCiFLuMWuAwrDgpEK+IZZaw/fPgHjKvxyOsqJew+dQqeVgMW7SQiKlMlim1Oq4KdUyVOYc9KBVJhY7nnDCksOtDin04ct9Bl3AKX5RQHjdtn4PgAFCr7Jkol1zddMQGdiMijEsU2p1XBzqkSp7BnpYp0ho3lnjOksGiUgplTUbAzbN64hUPDXs9IZX/1h92nTsHT6YTBFBGRR6mCnI6MTOZNubmLe9a0ZRvs3CYvo664mOb5N/iPkQn4c8iY/kFTUMFPbwHPkMKiUQpmlip0WU6By7B54xYODWu/5ow1sdYVptR96hQ8nU6M7u7uxCa79957u9euXZvYfERE5fhM2xx0zLWwv38IH4yMob3Zwupz5+HNP36IkTE712Rug4lNV5+NFQtbsevQHwte/6erzmLyOQCcuhCYexrw++eAsfxjT6sFuOLu4p11Z6wAPngbOLoPgLIDpcV/CVz4N8CRl4ET7+WDJ2XnPn3hn4E/ucJ/7EtvAZrn5/u9b7e/7K7COU9dGNjmjLlnoH12Ow788QCOjx5HW2NbUaFLb5umuiZYWQsncid820cRNm+UNUUdq6ujC8dGjuHVY69CQSEjGVzYeiFyKld0L977Wnn6ShwbORZpDaXWMR3ccccdR7u7u+8t1Y4J6EREREQ+WLSTiIiIKAEMpoiIiIg0MJgiIiIi0sA6U0REJXiPlnF26nlfcyed3967Hw89fxg5V15qu087v7FTSV6PcvxKGvN4j6Nx8z0+5jAAAZxaSmYjkBsFxk/69194FfDbHWUfOxNVnCNV4h5lE2e8eqMeI7mRiVpTDdkGDI8NRzoqxn0EjR/LsDArO0v7SJ3piAnoREQh/I6WMTMCCDCam/z5aZkGNl1t7+K7vXc/frTrTd/x3O38xnZfT0zI0SoVDajizhN0HI1b0PEx5Ypx7EzU702cI1WitNUdLyr3mHfuujP2ETTT6diYIExAJyKqgM3bDxad0Tc6rgoCKQAYHs1h8/aDAICHng84QsXTzm9s9/XE7NxYHIiMDtuvpznPzo3hgRQA5E7ahyNXIpDyW08Fvjc9u3uKgpmR3Ah6dveU1VZ3vKjcYz7y2iNa/Wc6BlNERCHiHAvjtM2VeMffaRc0duJH0YQcrZLqPFHnr+RRNN55K/C9iXOkSpS2lRgvKqd/uWf5zbRjY4IwmCIiChHnWBinrSE+ZdF92gWNnfhRNCFHq6Q6T9T5K3kUjXfeCnxv4hyponOUTZzxonL6l3sEzUw7NiYIgykiohB+R8uYGYFpFAZMlmlMJKZfd2Fn4Hjudn5ju68nJuRolVTnCTqOxi3o+JhyxTh2Jqo4R6qUe5RNnPGico9ZzhE0M/HYmCA8ToaIKITf0TLdX1yI5Z9tLXhtw6rPTiSNL/3MqXjngxM40P8e3A/8vO38xnZfT0zI0SqpzuN3HI2b7/Ex78HezZdnNtpf+z0KtFqAc68Hjr9T1rEzUcU5UqWco2zijGcZFnKu70VDtgFj42Mlj4rxHkHjxzIszK6brXWkTrXhcTJEREREGribj4iIiCgBDKaIiIiINDCYIiIiItLAYIqIiIhIA4MpIiIiIg0MpoiIiIg0MJgiIiIi0sBgioiIiEgDgykiIiIiDdm0F0BENN317ulH9+MHMDg8CgCY22Dim6sWAgA2bz+I/sFhGCLIKYX2Zmvi7L3N2w/iyOAw5uVfK+cYmd49/RUZJxX7tgA7NwJDh+3DilXOPq7FOfdu50ZgqM8+VHjZBvsIl4k+fYA11243/K7r82OTY0EA5+gTqwW4/NuTx8C453a3c5MMoMYn1+R3hIx7Pe51ToFth7ahZ3cPBo4PoLWxdeLcO/drXR1deKbvmYI2K09fWdB3Tt0ciAgGTwwiIxmMq/GJj5ZhYSQ3UnBkjPu1jGSw5ow1WPTxRdj0/CYMnRwCAJhiYkyN+R4105BtwIaLN0ys9ejxo0XzNtU1QUQwdGJoYn1DJ4YK7qGa8TgZIiINvXv6sf6RvRgdL/xZmhHAyAhGc8U/Y82MAIKCa5ZpYNPVZ8UKhHr39OO2x/ZjeHTyrLVyxknFvi3AE18HRn3O3MuYgAiQOzn5mmkB51wP7H3Qv08URh2w+nv250FzBzEtYNU9hYGS3z34tauAbYe2oftX3RjJjUy8lpUsRASj46OB/eqNeqz+1Gr87PWfFfTVJZDAM/qC2mcz2dC1Bqk36tF9SXcqAVXU42QYTBERaVhy19PoHyzzl7tHe7OF525dqj133HFScfeZ+XeFYph4x0lDU6f9Me7cTt+bX5n8OugevO0qYPnW5Th6/GhZfZ13f6aztsY27LhmR+LzRg2m+JiPiEjDkQoFUuWMFdS+kmuaMkN98fvoBlLlzhvUN2gsnTkCDBwfKLvvdA+kAL37TwIT0ImINMxrtlIbK6h9Jdc0ZZo64vcRozLzljO30zfs61Kva2htbC27b0am/696nftPwvT/DhMRpWj9igV2DpRHRgDTKH4dsHOmvNcs05hITI8zt2UWBhjljJOKZRvs/CI/GdPOb3IzLeD8G4L7RGHU2fOGzR3EtCYT4x1+4/i1q4B1561DvVFf8FpWsjAzZmi/eqMea85YU9RXl8D//+2w9qXWGqTeqJ9Itq9WfMxHRKTBSfROYzef035a7uZzErTj7uabf1FldvO55y53N1/BPUztbj4n+brc3XyLPr6Iu/mmEBPQiYiIiHxETUDnYz4iIiIiDQymiIiIiDQwmCIiIiLSwGCKiIiISAODKSIiIiINDKaIiIiINDCYIiIiItLAYIqIiIhIA4MpIiIiIg0MpoiIiIg0MJgiIiIi0sBgioiIiEgDgykiIiIiDQymiIiIiDQwmCIiIiLSwGCKiIiISAODKSIiIiINDKaIiIiINDCYIiIiItLAYIqIiIhIQ1ans4hcBqAHgAHgPqXUXRVZFRFRlerd04/N2w+if3AYhghySqG92cL6FQsAALc9tg/Do+MT7QWAcn0MY2YAV9cCGQGuv3A+Fn+iBZu3H8SRwWHMy8975aL2wLV5P7Z7+lSlfVuAnRuBoT6gqQNYtgE4+9rC6z//BjB8zP7abASys+yvxQBUDmjqLOzn7eNm1Nl/OOMno6/RagEWXgX8dgcwdHhyXqvFvj78LmDNnfy8qQP49HLgwE8n12C1AJd/O/ze/Nr4fK+2jR1Dz0dbMGAIWhvbsO68dVh5+sro9wNg26Ft6Nndg4HjA2htbEVXRxd+8btfYOjkkG/7tsa2km0qSSBQUGgr8/6mkihV6q93QEcRA8BrAD4PoA/ACwCuU0r9OqjP4sWL1YsvvljWfEREaevd04/bHtuP4dFc0TUzIxgdL+/naRxGRpBzzWOZOZyW+AAACFZJREFUBjZdfRYABK7Ny+lTlQHVvi3AE18HRocnXzMtYNU9dkCxbwvQeyMwPlp6LKcfEL1P0ow6YPX3wu/N3cYt/73aVifo/lgLRjKTD5vqjXp0X9IdOeDYdmgbun/VjZHciO4dJSLu/ZVLRF5SSi0u2U4jmLoYQLdSakX+69sAQCm1KagPgykims6W3PU0+geHSzdMWHuzBQCx1tbebOG5W5dO1ZLKd/eZ9js9Xk2dwM2vBF8P0tRpf4zTJ2lR7s1p45Zvv7xjHo6axQ+a2hrbsOOaHZGWsHzrchw9fjTuylMV5/7KFTWY0nnM1w7A/afeB+BCn4WsBbAWAObPn68xHRFRuo5UYSAFlLeuar0XDPWFvx50Pe541STKvfldy782kDV8uwwcH4i8hDhtq0U1rVknAV18Xit6m0spda9SarFSavEpp5yiMR0RUbrm5d8Bqjbzmq3Ya6vWe0FTR/jrQdfDxovbJ2lR7s3vWv611jH/R7utja2RlxCnbbWopjXrBFN9ADpdX3cAOKK3HCKi6rV+xQJYpv+7AGbG79+XlWd45rFMA+tXLAhdm5fTpyot22DnOrmZlv26cz1jRhvL6RenT9KMutL35m7jlv9erXt3EPXjhTsX6o16rDtvXeRlrDtvHeqN+jgrT1Xc+5tqOsHUCwA+LSKfFJE6AF8G8HhllkVEVH2uXNSOTVefNZGjZIgd2LQ3W9i85hx890vnwjILf6yK52MYM+QnckaAr1w0H/+65hy0N1uQ/LxOInnQ2rwf3X2q0tnX2knjTZ0AxP7oJJ8716/8t8ldc4C9m8/5WvIBpbufXx83ow7I1MVbp9UCLP6ryZwsZ16rJT+PFH7e1Gm3d6/BailMLPdbp7eNW/57tTL7UXS/8y7acgoCO5cobnL2ytNXovuSbrQ1tkEgaGtsw5cWfAlNdU2BfaK0qSTJ/y0q5/6mWtkJ6AAgIl8A8F3YpREeUEr9U1h7JqATERHRdJFEAjqUUk8BeEpnDCIiIqLpjBXQiYiIiDQwmCIiIiLSwGCKiIiISAODKSIiIiINDKaIiIiINDCYIiIiItLAYIqIiIhIA4MpIiIiIg0MpoiIiIg0MJgiIiIi0sBgioiIiEgDgykiIiIiDQymiIiIiDQwmCIiIiLSwGCKiIiISAODKSIiIiINDKaIiIiINDCYIiIiItLAYIqIiIhIA4MpIiIiIg0MpoiIiIg0MJgiIiIi0iBKqeQmE3kbwBuJTQh8DMA7Cc5XTXjvtalW771W7xvgvfPea0+S9/4JpdQppRolGkwlTUReVEotTnsdaeC9895rSa3eN8B7573Xnmq8dz7mIyIiItLAYIqIiIhIw0wPpu5NewEp4r3Xplq991q9b4D3Xqt471VkRudMEREREU21mf7OFBEREdGUmpHBlIhcJiIHReR1Ebk17fUkSUQeEJG3ROSVtNeSJBHpFJH/FJFXReSAiKxLe01JEZF6EflfEdmbv/c70l5T0kTEEJE9IvJk2mtJkoj8XkT2i8jLIvJi2utJkog0i8hWEflN/u/9xWmvKQkisiD/5+38956I3JT2upIgIjfnf8a9IiIPiUh92mtyzLjHfCJiAHgNwOcB9AF4AcB1Sqlfp7qwhIhIF4APAPxQKXVm2utJioi0AWhTSu0WkY8AeAnAlbXw5y4iAqBRKfWBiJgAngWwTim1K+WlJUZEbgGwGMAcpdQVaa8nKSLyewCLlVI1V29IRH4A4H+UUveJSB2ABqXUYNrrSlL+910/gAuVUknWcEyciLTD/tn2WaXUsIhsAfCUUur76a7MNhPfmboAwOtKqUNKqZMAfgJgdcprSoxS6hkAx9JeR9KUUkeVUrvzn78P4FUA7emuKhnK9kH+SzP/38z6V1IIEekAsBLAfWmvhZIhInMAdAG4HwCUUidrLZDKWwbg/2Z6IOWSBWCJSBZAA4AjKa9nwkwMptoBHHZ93Yca+aVKNhE5DcAiAM+nu5Lk5B9zvQzgLQC/VErVzL0D+C6AfwAwnvZCUqAA7BCRl0RkbdqLSdDpAN4G8B/5x7v3iUhj2otKwZcBPJT2IpKglOoH8C8A3gRwFMCQUmpHuquaNBODKfF5rWb+lV7rRGQ2gEcB3KSUei/t9SRFKZVTSp0LoAPABSJSE494ReQKAG8ppV5Key0pWaKUOg/A5QC+ln/MXwuyAM4D8O9KqUUAjgOotfzYOgBfBPBI2mtJgojMhf2U6ZMA5gFoFJGvpLuqSTMxmOoD0On6ugNV9FYgTZ18vtCjAH6slHos7fWkIf+o478AXJbyUpKyBMAX87lDPwGwVER+lO6SkqOUOpL/+BaAn8JOc6gFfQD6XO/AboUdXNWSywHsVkr9Ie2FJOTPAfxOKfW2UmoUwGMALkl5TRNmYjD1AoBPi8gn85H7lwE8nvKaaIrlk7DvB/CqUuo7aa8nSSJyiog05z+3YP/Q+U26q0qGUuo2pVSHUuo02H/Xn1ZKVc2/VqeSiDTmN1sg/4hrOYCa2MWrlBoAcFhEFuRfWgZgxm828bgONfKIL+9NABeJSEP+5/0y2LmxVSGb9gIqTSk1JiJ/C2A7AAPAA0qpAykvKzEi8hCAPwXwMRHpA/BNpdT96a4qEUsAfBXA/nzuEAD8o1LqqRTXlJQ2AD/I7+zJANiilKqpEgE16lQAP7V/ryAL4EGl1C/SXVKi/g7Aj/P/aD4E4C9SXk9iRKQB9o71v057LUlRSj0vIlsB7AYwBmAPqqgS+owrjUBERESUpJn4mI+IiIgoMQymiIiIiDQwmCIiIiLSwGCKiIiISAODKSIiIiINDKaIiIiINDCYIiIiItLAYIqIiIhIw/8Dlq9bCMSByJ4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x720 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "pconds = {'bln':'Baseline', 'p05':'PTZ at 5mM', 'p20':'PTZ at 20mM'}\n",
    "lvals  = {}\n",
    "durdat = np.ndarray((0,0,0))\n",
    "[fig, ax] = plt.subplots(1, figsize = (10,10))\n",
    "pi     = 0 \n",
    "\n",
    "for p in pconds.keys():\n",
    "    files = glob.glob(F['dom'] +os.sep+ 'PTZ' +os.sep+ \n",
    "                      'criticality' +os.sep+ 'size' +os.sep+ '*' +p+ '.npy')\n",
    "    bins      = np.exp(np.linspace(0,8,1000))\n",
    "    allcounts = np.ndarray([0,len(bins)-1])\n",
    "    \n",
    "#     for f in files:\n",
    "    f = files[8]\n",
    "    x = np.load(f)\n",
    "    [cnts, vals] = np.histogram(x,bins)\n",
    "    allcounts = np.vstack((allcounts,cnts))\n",
    "    \n",
    "    \n",
    "    ax.scatter(np.log(bins[:-1]),np.log(np.mean(allcounts,0)))\n",
    "#     ax.axis('equal')\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 162,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAD8CAYAAACLrvgBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAE/hJREFUeJzt3X9sXfV5x/H3g+OuLl3rMDwETrqgLcrEhkTgCtiYpnZdSWirkrGqBWkjqlAzqWVqtykb2T9s7SQyZT86pA2JFUbYWhijKUUrrRsBUrc/YHEIa/iVkrWlxPyIt2BYh7WG8OyP+zW7SR3761+ca/v9kq7uvc8595zHiX0/vt/zPT6RmUiSVOOUphuQJC0ehoYkqZqhIUmqZmhIkqoZGpKkaoaGJKmaoSFJqmZoSJKqGRqSpGormm5gvp1++um5Zs2aptuQpEVl7969/5mZA9Ott+RCY82aNQwPDzfdhiQtKhHxTM16Dk9JkqoZGpKkaoaGJKmaoSFJqmZoSJKqLbnZU7N1z74Rdgwd4Lmxcc7q72PrhnVsWj/YdFuS1FUMDdqBsW3XfsaPHgNgZGycbbv2AxgcktTB4Slgx9CBNwJjwvjRY+wYOtBQR5LUnQwN4Lmx8RnVJWm5MjSAs/r7ZlSXpOXK0AC2blhHX2/PcbW+3h62bljXUEeS1J08EM7/H+x29pQkTc3QKDatHzQkJGkaDk9JkqoZGpKkaoaGJKmaoSFJqmZoSJKqGRqSpGqGhiSpmqEhSapmaEiSqhkakqRqhoYkqZqhIUmqZmhIkqoZGpKkatOGRkSsjogHI+KJiHg8Ij5V6qdFxO6IeLrcryz1iIgbI+JgRHwrIs7v2Nbmsv7TEbG5o35BROwvr7kxImKqfUiSmlHzSeM14Pcy8xzgYuCTEXEOcB1wf2auBe4vzwEuA9aW2xbgJmgHAHA9cBFwIXB9RwjcBHy843UbS/1k+5AkNWDa0MjM5zPzkfL4v4EngUHgcmBnWW0nsKk8vhy4PdseAvoj4kxgA7A7M49k5kvAbmBjWfaOzHwoMxO4/YRtTbYPSVIDZnRMIyLWAOuBh4EzMvP5sugF4IzyeBB4tuNlh0ptqvqhSepMsQ9JUgOqQyMi3g58Cfh0Zr7Suax8Qsh57u04U+0jIrZExHBEDI+Oji5kG5K0rFWFRkT00g6ML2TmrlJ+sQwtUe4Pl/oIsLrj5atKbar6qknqU+3jOJl5c2a2MrM1MDBQ8yVJkmahZvZUALcAT2bmX3QsuheYmAG1GfhKR/3qMovqYuDlMsQ0BFwaESvLAfBLgaGy7JWIuLjs6+oTtjXZPiRJDVhRsc4lwG8C+yPi0VL7Q2A7cFdEXAM8A3ykLLsPeD9wEHgV+BhAZh6JiM8Ce8p6n8nMI+XxJ4DbgD7ga+XGFPuQJDUg2ocKlo5Wq5XDw8NNtyFJi0pE7M3M1nTreUa4JKmaoSFJqmZoSJKqGRqSpGqGhiSpmqEhSapmaEiSqhkakqRqhoYkqZqhIUmqZmhIkqoZGpKkaoaGJKmaoSFJqmZoSJKqGRqSpGqGhiSpmqEhSapmaEiSqhkakqRqhoYkqZqhIUmqZmhIkqoZGpKkaoaGJKmaoSFJqmZoSJKqGRqSpGqGhiSpmqEhSapmaEiSqhkakqRqhoYkqZqhIUmqZmhIkqpNGxoRcWtEHI6IxzpqfxQRIxHxaLm9v2PZtog4GBEHImJDR31jqR2MiOs66mdHxMOl/o8R8ZZS/7Hy/GBZvma+vmhJ0uzUfNK4Ddg4Sf0vM/O8crsPICLOAa4Efq685m8ioicieoC/Bi4DzgGuKusC/GnZ1s8ALwHXlPo1wEul/pdlPUlSg6YNjcz8JnCkcnuXA3dm5v9m5neBg8CF5XYwM7+TmT8E7gQuj4gAfgW4u7x+J7CpY1s7y+O7gfeW9SVJDZnLMY1rI+JbZfhqZakNAs92rHOo1E5W/wlgLDNfO6F+3LbK8pfL+j8iIrZExHBEDI+Ojs7hS5IkTWW2oXET8NPAecDzwJ/PW0ezkJk3Z2YrM1sDAwNNtiJJS9qsQiMzX8zMY5n5OvC3tIefAEaA1R2rriq1k9X/C+iPiBUn1I/bVln+zrK+JKkhswqNiDiz4+mvARMzq+4Friwzn84G1gL/BuwB1paZUm+hfbD83sxM4EHgw+X1m4GvdGxrc3n8YeCBsr4kqSErplshIu4A3g2cHhGHgOuBd0fEeUAC3wN+CyAzH4+Iu4AngNeAT2bmsbKda4EhoAe4NTMfL7v4A+DOiPgTYB9wS6nfAvx9RBykfSD+yjl/tZKkOYml9st7q9XK4eHhptuQpEUlIvZmZmu69TwjXJJUzdCQJFUzNCRJ1QwNSVI1Q0OSVM3QkCRVMzQkSdUMDUlSNUNDklTN0JAkVTM0JEnVDA1JUjVDQ5JUzdCQJFUzNCRJ1QwNSVI1Q0OSVM3QkCRVMzQkSdUMDUlSNUNDklTN0JAkVTM0JEnVDA1JUjVDQ5JUzdCQJFUzNCRJ1QwNSVI1Q0OSVM3QkCRVMzQkSdUMDUlSNUNDklTN0JAkVTM0JEnVpg2NiLg1Ig5HxGMdtdMiYndEPF3uV5Z6RMSNEXEwIr4VEed3vGZzWf/piNjcUb8gIvaX19wYETHVPiRJzan5pHEbsPGE2nXA/Zm5Fri/PAe4DFhbbluAm6AdAMD1wEXAhcD1HSFwE/DxjtdtnGYfi9I9+0a4ZPsDnH3dV7lk+wPcs2+k6ZYkacamDY3M/CZw5ITy5cDO8ngnsKmjfnu2PQT0R8SZwAZgd2YeycyXgN3AxrLsHZn5UGYmcPsJ25psH4vOPftG2LZrPyNj4yQwMjbOtl37DQ5Ji85sj2mckZnPl8cvAGeUx4PAsx3rHSq1qeqHJqlPtY8fERFbImI4IoZHR0dn8eUsrB1DBxg/euy42vjRY+wYOtBQR5I0O3M+EF4+IeQ89DLrfWTmzZnZyszWwMDAQrYyK8+Njc+oLkndarah8WIZWqLcHy71EWB1x3qrSm2q+qpJ6lPtY9E5q79vRnVJ6lazDY17gYkZUJuBr3TUry6zqC4GXi5DTEPApRGxshwAvxQYKsteiYiLy6ypq0/Y1mT7WHS2blhHX2/PcbW+3h62bljXUEeSNDsrplshIu4A3g2cHhGHaM+C2g7cFRHXAM8AHymr3we8HzgIvAp8DCAzj0TEZ4E9Zb3PZObEwfVP0J6h1Qd8rdyYYh+Lzqb17cM0O4YO8NzYOGf197F1w7o36pK0WET7cMHS0Wq1cnh4uOk2JGlRiYi9mdmabj3PCJckVTM0JEnVDA1JUjVDQ5JUzdCQJFUzNCRJ1QwNSVI1Q0OSVM3QkCRVMzQkSdUMDUlSNUNDklTN0JAkVTM0JEnVDA1JUjVDQ5JUzdCQJFUzNCRJ1QwNSVI1Q0OSVM3QkCRVMzQkSdUMDUlSNUNDklTN0JAkVTM0JEnVDA1JUjVDQ5JUzdCQJFUzNCRJ1QwNSVK1FU03oHr37Bthx9ABnhsb56z+PrZuWMem9YNNtyVpGTE0Fol79o2wbdd+xo8eA2BkbJxtu/YDGByS3jQOTy0SO4YOvBEYE8aPHmPH0IGGOpK0HM0pNCLiexGxPyIejYjhUjstInZHxNPlfmWpR0TcGBEHI+JbEXF+x3Y2l/WfjojNHfULyvYPltfGXPpdzJ4bG59RXZIWwnx80nhPZp6Xma3y/Drg/sxcC9xfngNcBqwtty3ATdAOGeB64CLgQuD6iaAp63y843Ub56HfRems/r4Z1SVpISzE8NTlwM7yeCewqaN+e7Y9BPRHxJnABmB3Zh7JzJeA3cDGsuwdmflQZiZwe8e2lp2tG9bR19tzXK2vt4etG9Y11JGk5WiuoZHANyJib0RsKbUzMvP58vgF4IzyeBB4tuO1h0ptqvqhSerL0qb1g9xwxbkM9vcRwGB/Hzdcca4HwSW9qeY6e+qXMnMkIn4S2B0RT3UuzMyMiJzjPqZVAmsLwLve9a6F3l1jNq0fNCQkNWpOnzQyc6TcHwa+TPuYxItlaIlyf7isPgKs7nj5qlKbqr5qkvpkfdycma3MbA0MDMzlS5IkTWHWoRERp0bEj088Bi4FHgPuBSZmQG0GvlIe3wtcXWZRXQy8XIaxhoBLI2JlOQB+KTBUlr0SEReXWVNXd2xLktSAuQxPnQF8ucyCXQF8MTO/HhF7gLsi4hrgGeAjZf37gPcDB4FXgY8BZOaRiPgssKes95nMPFIefwK4DegDvlZukqSGRHti0tLRarVyeHi46TYkaVGJiL0dp06clGeES5KqGRqSpGqGhiSpmqEhSapmaEiSqhkakqRqhoYkqZpX7luCvCyspIViaCwxXhZW0kJyeGqJ8bKwkhaSobHEeFlYSQvJ0FhivCyspIVkaCwxXhZW0kLyQPgSM3Gw29lTkhaCobEEeVlYSQvF4SlJUjVDQ5JUzeGpZcgzxiXNlqGxzHjGuKS5cHhqmfGMcUlzYWgsM54xLmkuDI1lxjPGJc2FobHMeMa4pLnwQPgyU3PGuLOrJJ2MobEMTXXGuLOrJE3F4Skdx9lVkqbiJw0dZ6rZVQ5bSfKTho5zsllU7+zrZduu/YyMjZP8/7DVPftG3twGJTXK0NBxTja7KgKHrSQZGjrepvWD3HDFuQz29xHAYH8fN1xxLmOvHp10/ZGxcc6+7qtcsv0BP3VIy0BkZtM9zKtWq5XDw8NNt7HkXLL9AUamOWu8r7eHX79gkAefGvW4h7TIRMTezGxNu56hoRonTsU9mQA6v6N6Twne/tYVjL161BCRupihoXnXOXtqtt81hojUnQwNLaia4aoaEyHy0qtH6YngWCaDhon0pjM0tKAmG646cWhqLiYLE0NFWjhLJjQiYiPwV0AP8PnM3D7V+obGm+fEk/3e87MDfGnvyLTHPebLKQGvJz8SKife9/f1EsFJA+g9PzvAg0+NMjI2Pmndg/paDpZEaERED/Bt4H3AIWAPcFVmPnGy18w1NCbeCDvfQObr/mRvUDN5/cRfo52ux+n2NdvlJ3sznayn+fzk0S1qg6r2fqpA65b7bu7R3ubvU/hSCY1fAP4oMzeU59sAMvOGk71mLqFRO0OoSb2nBAQcPdY9/2/d2JO0XPX19nDDFefOODhqQ6PbT+4bBJ7teH6o1BbEZH+sr9scfT277s25G3uSlquF/ksN3R4aVSJiS0QMR8Tw6OjorLfjJU8lLQUL+V7W7aExAqzueL6q1I6TmTdnZiszWwMDA7PemZc8bVZPxJTLp14qacJCvpd1e2jsAdZGxNkR8RbgSuDehdrZZH+sr9v0nhL09nTX2+d89NTX28NVF62e9I8lfu6j5/G97R/gu9s/wOc+eh6D5QdiImSmu+/v62Xl23qPq3XXv6A0fxb68s1dfT2NzHwtIq4FhmhPub01Mx9fqP11XgrV2VOzmz3VOf22c53pnk/M+Gj91GlTXrNjqqsOztTJrg8y2VTi2f6fdeMMm6XUo73N3+ypWl09e2o2PE9DkmZuqcyekiR1EUNDklTN0JAkVTM0JEnVDA1JUrUlN3sqIkaBZ+awidOB/5ynduaTfc2Mfc1ct/ZmXzMz275+KjOnPTt6yYXGXEXEcM20szebfc2Mfc1ct/ZmXzOz0H05PCVJqmZoSJKqGRo/6uamGzgJ+5oZ+5q5bu3NvmZmQfvymIYkqZqfNCRJ1QyNIiI2RsSBiDgYEdc13c+EiLg1Ig5HxGNN99IpIlZHxIMR8UREPB4Rn2q6J4CIeGtE/FtE/Hvp64+b7qlTRPRExL6I+Oeme5kQEd+LiP0R8WhEdM1f+4yI/oi4OyKeiogny+WfGxcR68q/1cTtlYj4dNN9AUTE75Tv+8ci4o6IeOu878PhqfYPMvBt4H20Lym7B7gqM59otDEgIn4Z+AFwe2b+fNP9TIiIM4EzM/ORiPhxYC+wqel/s4gI4NTM/EFE9AL/CnwqMx9qsq8JEfG7QAt4R2Z+sOl+oB0aQCszu+qcg4jYCfxLZn6+XE/nbZk51nRfncp7xwhwUWbO5fyw+ehlkPb3+zmZOR4RdwH3ZeZt87kfP2m0XQgczMzvZOYPgTuByxvuCYDM/CZwpOk+TpSZz2fmI+XxfwNPsoDXb6+VbT8oT3vLrSt+M4qIVcAHgM833Uu3i4h3Ar8M3AKQmT/stsAo3gv8R9OB0WEF0BcRK4C3Ac/N9w4MjbZB4NmO54fogjfAxSIi1gDrgYeb7aStDAE9ChwGdmdmV/QFfA74feD1phs5QQLfiIi9EbGl6WaKs4FR4O/KcN7nI+LUppuaxJXAHU03AZCZI8CfAd8HngdezsxvzPd+DA3NSUS8HfgS8OnMfKXpfgAy81hmnkf7mvIXRkTjw3oR8UHgcGbubbqXSfxSZp4PXAZ8sgyJNm0FcD5wU2auB/4H6JpjjQBlyOxDwD813QtARKykPUJyNnAWcGpE/MZ878fQaBsBVnc8X1VqmkI5ZvAl4AuZuavpfk5UhjMeBDY23QtwCfChcvzgTuBXIuIfmm2prfyGSmYeBr5Me7i2aYeAQx2fEu+mHSLd5DLgkcx8selGil8FvpuZo5l5FNgF/OJ878TQaNsDrI2Is8tvD1cC9zbcU1crB5xvAZ7MzL9oup8JETEQEf3lcR/tyQ1PNdsVZOa2zFyVmWtof389kJnz/lvgTEXEqWUiA2X451Kg8Zl6mfkC8GxErCul9wKNT0w5wVV0ydBU8X3g4oh4W/n5fC/tY43zasV8b3AxyszXIuJaYAjoAW7NzMcbbguAiLgDeDdwekQcAq7PzFua7Qpo/+b8m8D+cvwA4A8z874GewI4E9hZZrWcAtyVmV0zvbULnQF8uf0ewwrgi5n59WZbesNvA18ov8h9B/hYw/28oQTs+4DfarqXCZn5cETcDTwCvAbsYwHODnfKrSSpmsNTkqRqhoYkqZqhIUmqZmhIkqoZGpKkaoaGJKmaoSFJqmZoSJKq/R/yIQLCtDmwQwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "for ai in range(allcounts.shape[0]):\n",
    "    plt.scatter(np.log(bins[:-1]),allcounts[ai,:])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Network cycle difference "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Reducing dimensionality of the data to regional agerages"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "for f in range(len(Fish)): Fish[f] = cn.spacek(Fish[f])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "arange = []\n",
    "for a in range(10):\n",
    "    cd = []\n",
    "    for f in range(len(Fish)):\n",
    "        print('Now working on fish ' + str(f+1) + ': ' + Fish[f][\"Condition\"]) \n",
    "        (dt,lb,cs) = cn.average(Fish[f][\"Data\"], Fish[f][\"KLabel\"], Fish[f][\"Coordinates\"])\n",
    "\n",
    "        # Network construction\n",
    "        #=====================================================================\n",
    "        # Calculate correlation coefficient\n",
    "        #---------------------------------------------------------------------\n",
    "        A  = np.corrcoef(dt)\n",
    "        np.fill_diagonal(A, 0)\n",
    "\n",
    "        # Z-score and binarise according to z cut off\n",
    "        #---------------------------------------------------------------------\n",
    "    #     A  = stats.zscore(np.ravel(A)).reshape(A.shape)\n",
    "        bA = A > (a / 10)\n",
    "\n",
    "        # Generate Network x object from this\n",
    "        #---------------------------------------------------------------------\n",
    "        G = nx.from_numpy_matrix(bA)\n",
    "\n",
    "        # Calculate cycle density\n",
    "        #=====================================================================\n",
    "        cycs = nx.algorithms.cycle_basis(G)\n",
    "        no_E = G.number_of_edges()\n",
    "        if no_E > 0: ncd = len(nx.algorithms.cycle_basis(G)) / no_E\n",
    "        else:        ncd = 0\n",
    "\n",
    "        cd.append(ncd)\n",
    "        \n",
    "    arange.append(cd)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "swp = np.zeros((3,10))\n",
    "ai  = 0\n",
    "for a in arange:\n",
    "    swp[0,ai] = np.mean(a[0:9])\n",
    "    swp[1,ai] = np.mean(a[10:19])\n",
    "    swp[2,ai] = np.mean(a[20:])\n",
    "        \n",
    "    ai = ai + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x3407f3eb8>,\n",
       " <matplotlib.lines.Line2D at 0x3407fc048>,\n",
       " <matplotlib.lines.Line2D at 0x3407fc198>]"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xl4k1Xax/HvnaR7Swu0INCy7/tSFkUWd0AFVwRHZ8Z9BnRU3EdxZmScd0aFcRlccN8QUVFQUVwAUQG1rEXWUpaWtWyFQrck5/0jaSm10AJJnjS5P9eVK8mT05y7AX45nDw5R4wxKKWUCi02qwtQSinlexruSikVgjTclVIqBGm4K6VUCNJwV0qpEKThrpRSIUjDXSmlQpCGu1JKhSANd6WUCkEOqzpOTk42zZs3t6p7pZSqlZYsWbLHGJNSXTvLwr158+ZkZGRY1b1SStVKIrKlJu10WkYppUKQhrtSSoUgDXellApBGu5KKRWCqg13EXlNRHaLyKrjPC4i8qyIZInIShHp6fsylVJKnYyajNzfAIac4PGhQBvv5VbghdMvSyml1OmoNtyNMQuAfSdoMgJ4y3gsBpJEpJGvClRKKXXyfHGeexMgp8L9XO+xHT547t/4bvEbLN30BTYRHGLDZrdhF8Fms+Gw2bDbbIgIgoAAYgPK7nuuRcRzTATAc+1t57lvK3/8aFvb0bbe5/IctyM2GyJ2z32bDcR7X2xIhcsxx8vbVW5rB5vdU6vNXulxO+J9jLLfqux38v6+lY9JVW3LXh/v71Ox7THPc5JtbWIrf+x4t8va2byvZ9ntKn8O2zE1K6VqLqBfYhKRW/FM3dC0adNTeo5567/kI7Pal2WpWqCqwLdJ1bftYkfk2OuyN5HyCzZsNpunDcdpU+lS5fN5n8fGb9vabXbPdYXbDpujyscc4jh6zFbpeIXbZTVXfB6HODx9VnqesjaR9kgibZE4bA59kwwjvgj3bUBahfup3mO/YYyZAkwBSE9PP6WduW8Z8TIX7D1IUYmTolInhaVOiopLKSr13C9xeq6LSpwUlzopLnVR7CyluNRFSanL+7iLEmep577LSUmpCzBIxYu4EQw2jPcxsOGG8mtDtEOoG2ujTpSNpOgIEqOFhCgbCVFCnWg78ZGe23GRQqxDwLgxxoUxbnB7ro3bVX68/Np9tJ1xu4AK7Y372MeNy/OYu9RzcTkxbifGXQoup/e495ir1PN8FY95L5TfdmME72/t/XOrfPEGxG+Oex8zgNsRg4mIxURE446IwTiiMRExuCOiMI4YjCMKtyMa44jy3o7COKJxOyIx9kiMzY7buDHGYDC4jRu38bz+buPGjRvM0duV21V5wY3b7cZlXBgMLrfLc+w47Z1uZ/ltl3H95nmq+tmK7VzGhcvtOvbauE7lr73PRNgiysM+wh5xzHWkPfLo42W3yx73/kzZ8Yo/U3btsDmOuV+xbYwjhviIeOIj4olxxOibTAD4ItxnAbeLyDSgL5BvjPHLlAxAk8QEmiQm+PQ5jTEUO90Ul7opLHV5LiWe66JKt4vKH3dzsKiUvQXF7CkoYVNBMXt2l7DvcDHuKt62HDahfnwkyfFR1I+PIjk+kpT4qPJjyd7bKfFR1IuLxGG36CxVtxvcpeAqAVep91JS4X7JcY4VQ3EBFB+EonzvxXu7+CAcPgBFO4/e9wb1cdmjIDoRout4rqPqVLqf6L1fqU1MEiQ0BltwnuVrjCkPeZfbhdM4cbvdOI3zmDcCp3GWv8FUfHMoe8Opqn3FN5Cy26XuUkpdpZS4SyhxlVDqLj3muuLxsnaFzkLyi/M9x8raeduWtXFX9+d3AjaxERcRR0JEAnGR3uuIOOIj48vfAOIj4z1tIhOOaVvx8Sh7lA//ZEJPteEuIu8Bg4FkEckF/gZEABhjXgRmA8OALOAIcIO/ivUXESE6wk50hJ1Ez692ylxuw/4jJewtKGFPQbH34r19qJi9hz23s3YdYk9BCSWuqv+R1I2NOCb0k+OjSEmIon6c983AezslIYroCPtp1XwMmw1sUeDw4z8cY6Ck4LdvAOVvChXvVzh+cNvRY87C4z+/IwZS2kJKe8+lQQfPdVIzy0NfRHCIAwcO8OEfW6C53K7fvGGU3y57E6jw5lHkLKKgtMBzKfFcHy49zKGSQxwuPcyewj1sObil/PESd0m1NUTYIsqDPj4ivso3iLLbZW8i3Rt0JzEqMQCvkPXEmFOaHTlt6enpJtwXDjPGcKjYeTT0D1V6MygorvAmUUJBsbPK52mUGE231CS6N02iW2oSXVMTiYuybE24wHCWVP2GcGQv7NkAeWth91o4tP3oz5SHfgdo0N5zndIuKEJfHavEVeJ5Ayg5zKFSzxtA2ZtCxTeIqt4oKrapPA2WFJXEven3MrzV8Fo7NSQiS4wx6dW203CvPYpKXeVBv9cb/nmHilm/q4AVuQfYsvcIADaBNg0S6J6WRLe0JLqnJdG2Ybx1Uz1WKjwAees8YZ+3Fnav8VwfqjBzGBELyW2PjvBT2nvCP7Gphn4tZoyhyFVUHvx5R/J4btlzLM9bTu8zejO+33haJLawusyTpuEehvYdLmFF7gGWbz3guc45wIEjpQBER9jo0iTxmMBvkhTGH2yVh/4azwg/b43n/jGhH3d0pJ/S7mj4J6Zp6NdSbuNmxoYZTFoyiSJnETd3uZmbutxUq+bvNdwVxhi27jvC8pwD5Zdftx+kxOmZ50+Oj/SEfaon8LulJpEYe3qfOdR6hfs9Ib97zbHhX7DzaJuIOE/Yl43wy6Z56qRq6NcSewr38OQvTzJ702ya12nO+H7j6dOoj9Vl1YiGu6pSidPNup2HWJ6zn+U5+SzP2c/GvMPlj7dMjjtmdN++UQJRjlr8yZ+vHBP6a4+Gf8XQj4w/Or3TqDv0vB4iYqyrWVVr4baFTFg8gdyCXIa3Gs696fdSN7qu1WWdkIa7qrGDRaWszMkvn8pZnnOAvEPFAETabXRoXIceaUl0S0uke1pdmtePDd/pnMqO7Ds6wq8Y/gW7oH5rGDEZmvazukp1AkXOIqasnMLrq14nLjKOe3rdw2WtLwvav+Ma7uqUGWPYkV/E8pwDrMg5wLKcA2Tm5lNY6jnzIDEmwjOyT00sP0OnfnztmbMMiI3zYNZfID8H+t4G5z0KkXFWV6VOIGt/Fo8tfoxlu5eR3jCd8WeOp2ViS6vL+g0Nd+VTTpebrLyC8g9rl209wPpdh8q/sJVWL4ZLuzbmxrNbkKxB71FcAN/+A36e4jndcviz0HKw1VWpE3AbNx9v+JiJSyZS6Czk5i43c3OXm4PqA1cNd+V3R0qcZOZ6pnN+yt7H3HW7ibTbGNU7jVsGtiS1bqzVJQaHzT/CrNthXzb0+iNc8JjnG7UqaO0p3MNTGU/xefbnNKvTjPH9xtO3UV+rywI03JUFNuYV8NJ3G5mx1LO00IjuTfjz4Ja0buDb5SJqpZIjMP9fsGgyJDSCS56GthdaXZWqxqLti/jn4n+y9dBWLm15Kff2vpd60fUsrUnDXVlm24FCXl6QzbRftlLsdHNRxzMYc04ruqYmWV2a9XIzYOZYz4eu3UbDRf+CWGvDQp1YkbOIlzNf5rVVrxEXEce4XuO4rPVl5ctWB5qGu7Lc3oJiXv9xM28u2syhIicD2iQzZnBr+rWsF7RnIgSEsxgWPAnfT4LY+nDJJOhwqdVVqWpkH8jmH4v+wdLdS+nZoCePnvkorZJaBbwODXcVNA4VlfLO4q28+sMm9hQU06NpEmMHt+bc9g2w2cI45Hes8Izid2ZCp8th6JMQn2J1VeoE3MbNzKyZPJXxFEecR7ix843c0uUWoh3RAatBw10FnaJSFx9k5PDSgmxy9xfSrmECY85pxcVdGoXnujfgWTL5x2fgu/94vgQ19AnoclX5blsqOO0t3MvEjIl8mv0pTROa8ki/Rziz8ZkB6VvDXQWtUpebT1ds54X5G9mwu4Cm9WK5bVBLruyZ6tvli2uT3Ws9o/htGdB2qGeqpk5jq6tS1Vi8YzETFk1g66GtXNzyYu5Lv4/6MfX92qeGuwp6brfh6zW7eH7+RlbkHCAlIYpbBrTg2r7NiA/1JYur4nbB4hdg7gTPZiUXPQ49rtNRfJArdhXzSuYrvJL5CrGOWMb1GsflbS732weuGu6q1jDGsHDjXp6fn8WPWXtJjIngD2c154azmlM3LtLq8gJv70aYdQds+RFanuP58lPSqe05rAInOz+bxxY9xpJdS/z6gauGu6qVlucc4Pl5WXy1ehexkXZG92nKLQNackZi4D6wCgpuN2S8Cl//zTNyP//vkH6TrjoZ5IwxzNzo+cD1cOlhbuh0A7d2vdWnH7hquKtabf2uQ7w4fyMzV2zHJnBlz1RuG9SKFslhtj7L/i3w6Z2QPQ+a9Yfhz0H9wJ9+p07OvqJ9TMyYyKyNs0iNT2V8v/Gc1eQsnzy3hrsKCTn7jjBlQTbvZ+TgdLkZ1qURYwa3pmPjOlaXFjjGwLJ3YM7Dns3Iz30E+v0ZbGH64XMt8vOOn5mweAKbD25mWIth3Nf7PpJjkk/rOTXcVUjZfaiI137YzDuLt1BQ7OScdimMPac16c3D6NudB7fDZ+Ng/RfQJN2znHCD9lZXpapR7Crm1cxXeSXzFaId0dzd626ubHPlKX/gquGuQlL+kVLeWrSZ1xduZt/hEvo0r8eYc1oxqG1KeHzr1RhY9RHMvg9KCmDQA9D/TrCH+Q5atcCm/E1MWDyBX3b+wt297ubGzjee0vNouKuQdqTEyfu/5DBlQTY78ovo1LgOYwa3ZkjnM7CHw7deC/Jg9r2w+hM4o6tnFN+oq9VVqWoYY/gs+zMGpg4kMerUVgbVcFdhocTp5pNl23jxu41k7zlMy+Q4nhnVgy6pYbKk7uqZ8Pm9ULgPzh4HA+8FR/CsPa58T8NdhRWX2/Dlqp38a/YaDhaW8voNvcNnPv7IPvjyIVg5zbNZ94jJkNrL6qqUn9Q03PWkWRUS7Dbh4q6N+OBPZ5KSEMX1r/7Mj1l7rC4rMGLrwRUvwbUfQPFBePV8+Go8lBZaXZmykIa7CimNk2KYdls/mtaL5YY3fmHu2l1WlxQ4bS+EMYug5+9h4bPwQn/YssjqqpRFNNxVyGmQEM20W/vRrmECt729hNmZO6wuKXCiE+HSZ+D3M8FdCm9eCuvnWF2VsoCGuwpJdeMiefeWvnRNTeL2qUuZsTTX6pICq+VguO17aNgJ3r8esudbXJAKNA13FbLqREfw9k196NeyPvd8sIJ3f9pidUmBFZME13/sWa7gvdE6RRNmNNxVSIuNdPDaH3tzTrsGPPzxKl75PtvqkgIrtp5niqZOY3j3ati2xOqKVIBouKuQFx1h58XrejGsyxn88/M1PPftBqw6BdgS8Q3g97M8Qf/2FbBzldUVqQCoUbiLyBARWSciWSLyYBWPNxWReSKyTERWisgw35eq1KmLdNh4dlQPrujRhIlfr+eJOevCK+ATm8AfZkFkHLw1AvLWW12R8rNqw11E7MBkYCjQERgtIh0rNXsEmG6M6QGMAp73daFKnS6H3cZTV3fj2r5NeWH+Rv7x6Wrc7jAK+LrNPSN4scFbw2FfmE1RhZmajNz7AFnGmGxjTAkwDRhRqY0BytZgTQS2+65EpXzHZhMev6wzN53dgjcWbuahGZm4wingk1t75uCdxfDmCDiQY3VFyk9qEu5NgIp/A3K9xyr6O3CdiOQCs4E7qnoiEblVRDJEJCMvL+8UylXq9IkIj1zcgb+c25r3M3IYN305pS631WUFTsOOnrNoivI9I/hDO62uSPmBrz5QHQ28YYxJBYYBb4v8drFiY8wUY0y6MSY9JSXFR10rdfJEhHEXtuOBIe2ZuXw7Y99dSrHTZXVZgdO4O1z3IRza5ZmDPxwmSzWEkZqE+zYgrcL9VO+xim4CpgMYYxYB0cDpbTeiVAD8eXAr/n5pR75avYtb31pCYUkYBXxaH7j2fdi/Gd6+DAr3W12R8qGahPsvQBsRaSEikXg+MJ1Vqc1W4DwAEemAJ9x13kXVCn/s34InruzKgg15/PH1nykodlpdUuC0GACj3oW8dfDOVVB8yOqKlI9UG+7GGCdwOzAHWIPnrJhfReQxERnubXYPcIuIrADeA/5owuo8M1XbjeydxtPXdCdjy36uf/Un8o+UWl1S4LQ+H65+A7Yvg6nXQMkRqytSPqDruStVwZxfd3LH1GW0bhDP2zf1oX58GG18kfkhfHSzZ12a0dMgItrqilQVdD13pU7BRZ3O4OU/pLMxr4BRUxaz62CR1SUFTperYMT/IHsefPBHcIXR/15CkIa7UpUMapvCmzf2YfuBQka+tIjc/WE0TdHjOhj2FKz/wjOKd4XR5w8hRsNdqSr0a1mft2/uy/7DJYx8cRGb9hy2uqTA6XMLXPhPz+bbs24Hdxh9ByCEaLgrdRw9m9blvVv7UeR0M/KlRazfFUZnkpx1B5zzMKx4Dz4fB3p+RK2j4a7UCXRqnMj7t/ZDgGteWsSqbflWlxQ4A++Ds++GJa/DnIc14GsZDXelqtGmYQLTbzuT2EgHo19ezJItYfJlHxE472/Q90+weDLMe9zqitRJ0HBXqgaaJ8cx/U9nUj8ukutf/YlFG/daXVJgiMCQf0PPP8CCJ2HBU1ZXpGpIw12pGmqSFMP0286kSVIMf3z9Z+av2211SYEhApf8F7qMhLkTYJGu6F0baLgrdRIa1Inm/dvOpHWDeG55K4MvV4XJioo2O1z2AnQYDnMegozXrK5IVUPDXamTVC8ukqm39KNLk0TGTl3KJ8sqr6MXouwOuPJVaHMhfDYOlr9ndUXqBDTclToFiTERvH1TX/o0r8fd05fz3s9brS4pMByRMPJtaDEQZo6BXz+2uiJ1HBruSp2iuCgHr9/Qm0FtU3hoRiav/bDJ6pICIyIaRr8HaX0932Jd94XVFakqaLgrdRqiI+y8dH0vhnQ6g8c+W83keVlWlxQYkXFw7XQ4oytM/z1snGt1RaoSDXelTlOUw87/ru3BZd0b8+ScdTw5Zy1hseJ1dB247iNIbgvvXQubf7S6IlWBhrtSPuCw25g4sjuj+6Qxed5GHvtsNc5w2Jc1th5c/wkkpcHUkZCry3gHCw13pXzEbhP+dXkXbuzfgtd/3MyVLyxk7c6DVpflf/Ep8PtZEJcC71wBO1ZaXZFCw10pnxIRxl/Sgf9d24Pc/YVc8uwPTPpqXehvvl2nEfxhFkQmePZj3b3W6orCnoa7Uj4mIlzStTHfjBvE8G6NeXZuFhc/+0Por0mT1NQT8DYHvDUc9m60uqKwpuGulJ/UjYtk0jXdef2G3hwpdnLViwv5x6e/cjiUN+Cu38ozReN2wpvD4UCYnP8fhDTclfKzc9o14Ktxg7i+XzNe/3EzFz29gO835Fldlv80aA/XfwwlhzwBf3CH1RWFJQ13pQIgPsrBYyM6M/22M4m027j+1Z+5/8MV5B8J0X1KG3WD62bA4TzPFE1BCL+ZBSkNd6UCqE+Lesy+cwBjBrfio6XbOP+/3/HlqhAd2aame77odCDHs+F2OJz7H0Q03JUKsOgIO/cPac/Msf1JiY/iT+8s5c/vLGH3oSKrS/O95v3hwgmw5QfInm91NWFFw10pi3RuksjM2/tz/5B2fLt2NxdMWsAHGTmh9+3Wnr+HOqkw//909B5AGu5KWSjCbmPM4NZ8cecA2jaM574PV/L7134mZ98Rq0vzHUcUDBgHOT9B9jyrqwkbGu5KBYFWKfG8f+uZTBjRiaVb9nPR0wt448dNuNwhMtLtcZ1n9D5PR++BouGuVJCw2YTrz2zOV+MG0bt5Pf7+6WqufnEhWbsPWV3a6XNEwcB7IPdn2Pit1dWEBQ13pYJMk6QY3rihN/+9phvZew4z7JkfeO7bDZTW9oXIul8HiWkw/986eg8ADXelgpCIcHmPVL4ZN4gLOzVk4tfrufS5H8jMzbe6tFPniIQB90DuL5Clo3d/03BXKoglx0fxv2t7MuX6Xuw/UsKIyT/wf1+soai0li5E1v13kNhUz5wJgBqFu4gMEZF1IpIlIg8ep81IEVktIr+KyFTflqlUeLuw0xl8dfcgrumdxkvfZTPk6QUszt5rdVknzxHpmXvflgFZ31hdTUirNtxFxA5MBoYCHYHRItKxUps2wENAf2NMJ+AuP9SqVFhLjIng/67oytSb++I2MGrKYv76cSYHi2rZEgbdrtXRewDUZOTeB8gyxmQbY0qAacCISm1uASYbY/YDGGN2+7ZMpVSZs1onM+eugdwyoAXTft7KhZMWMHftLqvLqjlHJAy8F7YtgQ1fW11NyKpJuDcBcircz/Ueq6gt0FZEfhSRxSIyxFcFKqV+KybSzsMXd2TGmP4kxkRw4xsZ3DltGXsLiq0urWa6X+tZ/11H737jqw9UHUAbYDAwGnhZRJIqNxKRW0UkQ0Qy8vJ0lTilTlf3tCQ+veNs7j6/LbMzd3DBfxcwc/m24F/CwB4BA++D7Uthw1dWVxOSahLu24C0CvdTvccqygVmGWNKjTGbgPV4wv4Yxpgpxph0Y0x6SkrKqdaslKog0mHjzvPb8PlfBtC0Xix3TlvOzW9msCO/0OrSTqzbaEhqpqN3P6lJuP8CtBGRFiISCYwCZlVq8wmeUTsikoxnmibbh3UqparRtmECH/35LMZf0pGFG/dywaQFTM/Iqf4HrVI+el8G6+dYXU3IqTbcjTFO4HZgDrAGmG6M+VVEHhOR4d5mc4C9IrIamAfcZ4yphedpKVW72W3CTWe3YM5dA+nSJJH7P1zJN6uD+MPWbqOgbnMdvfuBWDU3l56ebjIyMizpW6lwUFTq4soXFpK7v5DP7jibtHqxVpdUtWXvwMyxMHoatBtqdTVBT0SWGGPSq2un31BVKkRFR9h5/nc9cbsNt09dSokzSNem6ToK6rbQ0buPabgrFcKa1Y/jyau7siI3n3/NXmN1OVWzOzxz7ztWwLovrK4mZGi4KxXihnRuxI39W/DGws3MzgzS/Vq7XqOjdx/TcFcqDDw4tD3d05J44MOVbN5z2OpyfsvugEH3w86VsG621dWEBA13pcJApMPG5N/1xG4Xxry7NDhXlewyEuq11NG7j2i4KxUmmiTFMGlkN1bvOMhjn622upzfsjtg4P2wMxPWfm51NbWehrtSYeTc9g3506BWTP1pK58sq/xF8yDQ5Wqo18qzW5M7SM/uqSU03JUKM/de2JY+zevx148zg29/1rK5912ZsE5H76dDw12pMOOw23h2dA9iIuyMeXcpR0qcVpd0rM5XQf3WOno/TRruSoWhMxKjeXpUdzbsLuCRT1YF1yqSdgcMegB2rYK1n1ldTa2l4a5UmBrQJoW/nNuGGUu38UFGrtXlHKvzlVC/jY7eT4OGu1Jh7C/ntaF/6/qMn7mKNTsOWl3OUTa7Z/S++1dYU3kRWlUTGu5KhTG7TXj6mh4kxkQw9t2lFBQH0fx75ysguS189x8dvZ8CDXelwlxKQhTPje7B5r2HeWhGZvDMv5eP3lfDmplWV1PraLgrpejbsj73XNiOT1ds552ftlpdzlGdLofkdjBfR+8nS8NdKQXAnwe1YnC7FCZ8uprM3Hyry/Gw2T3nveetgdWfWF1NraLhrpQCwGYT/juyO8nxkYyZuoT8wlKrS/IoG73r3PtJ0XBXSpWrGxfJc9f2ZMeBIu77YEVwzL/b7DD4AchbC6s/trqaWkPDXSl1jF7N6vLg0PZ8tXoXr/6wyepyPDpeBint4bsnwB2EK1oGIQ13pdRv3HR2Cy7s2JB/f7GWpVv3W13O0TNn8tbCrzp6rwkNd6XUb4gIT17djUZJ0dz+7lL2Hy6xuiTv6L2Djt5rSMNdKVWlxJgInr+2F3sKShg3fTlut8Xz7zabZ+59zzodvdeAhrtS6ri6pCYy/pIOzFuXx4sLNlpdDnQYAQ06es+c0dH7iWi4K6VO6Lp+zbikayOemrOOxdl7rS3GZvPMve9ZD6tmWFtLkNNwV0qdkIjw7yu70rx+HH95bxl5h4qtLajDcGjQSUfv1dBwV0pVKz7KweTf9SS/sJS73l+Gy8r597K5970bYNVH1tUR5DTclVI10qFRHSaM6MyPWXt59tsN1hbT/lJo2NkzencF0UqWQUTDXSlVY1enp3Jlz1SenbuBHzbssa6Qsrn3vVk6ej8ODXelVI2JCBMu60SbBvHcOW0Zuw4WWVdM+0ugYRcdvR+HhrtS6qTERjp4/nc9KSx1ccfUZThdFi3mVTb3vm8jrPrQmhqCmIa7UuqktW6QwOOXd+bnzfuY+PV66wppfwmcoaP3qtQo3EVkiIisE5EsEXnwBO2uFBEjIum+K1EpFYwu75HK6D5NeWH+Ruau3WVNESIw6EHYlw2ZH1hTQ5CqNtxFxA5MBoYCHYHRItKxinYJwJ3AT74uUikVnP52aUc6NqrD3e+vIHf/EWuKaH8xnNEVFjyho/cKajJy7wNkGWOyjTElwDRgRBXtJgD/ASz8hEUpFUjREXae/11PXG7D7VOXUeK0YP5dBAY/5B29Tw98/0GqJuHeBMipcD/Xe6yciPQE0owxn/uwNqVULdA8OY4nrurK8pwD/PuLtdYU0W4oNOrmWTFSR++ADz5QFREbMAm4pwZtbxWRDBHJyMvLO92ulVJBYliXRvzxrOa89uMmvly1M/AFlI3e92+CldMC338Qqkm4bwPSKtxP9R4rkwB0BuaLyGagHzCrqg9VjTFTjDHpxpj0lJSUU69aKRV0/jqsA93SkrjvwxVs2Xs48AW0HQKNusOCJ8EVJPu/Wqgm4f4L0EZEWohIJDAKmFX2oDEm3xiTbIxpboxpDiwGhhtjMvxSsVIqKEU6bPxvdA8EGDt1KUWlAV7Uq3z0vhlW6Oi92nA3xjiB24E5wBpgujHmVxF5TESG+7tApVTtkVYvlokju7Nq20H++fnqwBfQ9iJo3ENH79Rwzt0YM9sY09YY08oY87j32KPGmFlVtB2so3alwtcFHRty28CWvLN4KzOXb6v+B3ypbPR+YAuseC+wfQcZ/YaqUsrn7r2oHenN6vLQjEyydhcL0XX+AAAOcklEQVQEtvM2F0LjnmE/etdwV0r5XITdxnPX9iA6ws7Yd5dSWBLA+ffy0ftWWD41cP0GGQ13pZRfNEqM4elrurN+9yH+PuvXwHbe5gJo0gu+fwqcJYHtO0houCul/GZg2xT+PKgV72fk8M3qAK4/U3H0viI8R+8a7kopv7rr/LZ0aFSHB2dksu9wAEfRrc+HJumwYGJYjt413JVSfhXpsDFpZDfyC0t45JNMjAnQ/qtlo/f8rbD83cD0GUQ03JVSftehUR3uvqAtszN3MmvF9sB13Po8SO0N34ff6F3DXSkVELcNbEXPpkmM/2QVO/MDtHisCAx+EPJzYPk7gekzSGi4K6UCwm4TJo7sTqnL8MBHKwM3PdPqPEjt4517Lw5Mn0FAw10pFTAtkuN4aFh7vlufx9Sftwam07LR+8FcWPZ2YPoMAhruSqmAuq5vM85unczjn68J3OqRrc6FtH6w4CkosWjHqADTcFdKBZTNJjxxVVfsNuHeD1bgcgdgekYEzv8bHNoBP0/xf39BQMNdKRVwjZNi+Pulnfhl835e/SE7MJ02O8uz7swP/4XCA4Hp00Ia7kopS1zRswkXdmzIU3PWs37XocB0eu54KDoAC58NTH8W0nBXSllCRPjXFV1IiHYwbvpySl0B2Fy7UVfofBUsfgEOBXA5BAtouCulLJMcH8Xjl3dh1baDPDc3KzCdnvNXcJV4lgQOYRruSilLDel8Blf0aMLkeVmsyAnAXHj9VtDz97Dkddi3yf/9WUTDXSllub8N70SDhCjGTV8emL1XB94PtgiY/3/+78siGu5KKcslxkTwxFVd2Zh3mCfnrPN/h3UaQd/bYOV02BXgteYDRMNdKRUUBrRJ4fp+zXjtx00szt7r/w7Pvgui6sC3E/zflwU03JVSQeOhYe1pVi+Wez9YQUGx07+dxdSFs++E9V/A1sX+7csCGu5KqaARG+lg4shubD9QyD8/W+3/Dvv+CeIawDf/gEAtZBYgGu5KqaDSq1k9bh3Yimm/5DB3rZ/PRY+Mg0H3w9aFkPWNf/sKMA13pVTQufuCNrQ/I4EHPspkv7+35uv5B0hqBt/+A9wB+CJVgGi4K6WCTpTDzsSR3ThwpITxM1f5tzNHJJz7COzMhF9n+LevANJwV0oFpU6NE7nzvDZ8tnIHn/p7a77OV0GDTjDvcXCV+revANFwV0oFrT8NakX3tCTGz1zF7oN+3JrPZoPzHoV92SGzoYeGu1IqaDnsNiaO7EZRqcv/W/O1vcizocd3T4TEhh4a7kqpoNYqJZ4HhrRn3ro83v8lx38dhdiGHhruSqmg94czm3Nmy/pM+Gw1Ofv8OKoOoQ09NNyVUkHPZhOevLorIp6t+dz+3JovRDb0qFG4i8gQEVknIlki8mAVj48TkdUislJEvhWRZr4vVSkVzlLrxvLopR35adM+XvvRj0v1hsiGHtWGu4jYgcnAUKAjMFpEOlZqtgxIN8Z0BT4EnvB1oUopdXWvVM7v0IAn5qwja7cft+YLgQ09ajJy7wNkGWOyjTElwDRgRMUGxph5xpiyibDFQKpvy1RKqaNb88VF2hk3fYX/tuYLgQ09ahLuTYCKH1Hneo8dz03AF1U9ICK3ikiGiGTk5eXVvEqllPJqkBDN45d3YWVuPs/P2+i/jmr5hh4+/UBVRK4D0oEq/y9jjJlijEk3xqSnpKT4smulVBgZ1qURI7o35rm5G8jMzfdPJ7V8Q4+ahPs2IK3C/VTvsWOIyPnAw8BwY0yxb8pTSqmqPTa8M/XjI/27NV8t3tCjJuH+C9BGRFqISCQwCphVsYGI9ABewhPsu31fplJKHSsxNoL/XNmVDbsLmPT1ev90Uos39Kg23I0xTuB2YA6wBphujPlVRB4TkeHeZk8C8cAHIrJcRGYd5+mUUspnBrdrwLV9m/Ly99n8vGmffzqppRt6iF/XajiB9PR0k5GRYUnfSqnQcbjYydBnvgfgizsHEBfl8H0nP78Ms++F330IbS7w/fOfBBFZYoxJr66dfkNVKVWrxUU5eOrqbuTsP8Ljs9f4p5NauKGHhrtSqtbr06IetwxoydSftjJ/nR8+9quFG3pouCulQsK4C9rStmE8D3y0kvwjfthwo2xDj7n/rBUbemi4K6VCQnSEnUkju7O3oIRHZ/lha76yDT32b6oVG3pouCulQkbnJonccW4bZi7fzuzMHb7voO1FkNYX5v8n6Df00HBXSoWUMee0omtqIg9/nMnuQz7emk8Ezv87FOwM+g09NNyVUiElwm5j0shuHC5x8dcZmb7fmq/ZWdD6gqDf0EPDXSkVclo3SOD+i9rxzZrdfLAk1/cdnPdo0G/ooeGulApJN/ZvQd8W9Xjs09Xk7vfx/HijrtD5yqDe0EPDXSkVkmw24amru2GM4b4PVvp+a75zHg7qDT003JVSISutXizjL+nIouy9vLlos2+fPMg39NBwV0qFtGt6p3Fu+wb8+4u1vt+aL4g39NBwV0qFNBHh31d0ITbSzuWTF/L24i2+m6IJ4g09NNyVUiGvQZ1oPhnbn65piYz/ZBUjX1rEhl0+GsUH6YYeGu5KqbDQrH4c79zUl6eu7kZWXgHDnv2e/369nmLnae7iFKQbemi4K6XChohwVa9Uvhk3iIu7NOKZbzcw7JnvT3+jjyDc0EPDXSkVdpLjo3h6VA/euKE3xU43I19axEMzMskvPMXVHiPjYND9sHUhZH3j22JPkYa7UipsDW7XgK/uHsgtA1rw/i9buWDSd3yRuePUliwo29Djm+DY0EPDXSkV1mIjHTx8cUdmjj2blIQo/vzuUm55awk78gtP7onKNvTYFRwbemi4K6UU0CU1kZlj+/PQ0Pb8kJXHBZMW8NaizSd32mQQbeih4a6UUl4Ou43bBrXiq7sG0aNpEo/O/JWrXlzIup01PG0yiDb00HBXSqlKmtaP5a0b+zBpZDc27TnMJc99z8Sv1lFUWoPTJoNkQw8Nd6WUqoKIcEXPVL69ZzCXdm3Mc3OzGPbM9yzO3lvdDwbFhh4a7kopdQL14iKZdE133r6pD6VuN6OmLObB6jbhDoINPTTclVKqBga0SeGruwZx28CWfLAkl/MmfcfnK09w2qTFG3pouCulVA3FRNp5aFgHZo7tzxmJUYydupSb38xg+4EqTpu0eEMPDXellDpJnZsk8smY/jxycQcWbtzLBZO+440fN+GqfNqkhRt6aLgrpdQpcNht3DygJV/dPZBezevx909Xc+ULC1m78+DRRvVbQY/rLdnQQ8NdKaVOQ1q9WN68oTfPjOpOzr4jXPLsDzw5Z+3R0yYHPQA2B8z7V0Dr0nBXSqnTJCKM6N6Eb8YNYkT3Jkyet5Ghz3zPwo17jm7okfkB7FwVsJpqFO4iMkRE1olIlog8WMXjUSLyvvfxn0Skua8LVUqpYFc3LpKJI7vxzk19cbkN1778E/d/uIL8nmM9G3rMDdyGHtWGu4jYgcnAUKAjMFpEOlZqdhOw3xjTGvgv8B9fF6qUUrXF2W2SmXPXQP40qBUfLd3GeS+sZE2rG2D9lwHb0KMmI/c+QJYxJtsYUwJMA0ZUajMCeNN7+0PgPBER35WplFK1S0yknQeHtufT28+mSVIMVyztSr6tLsVfPhqQDT1qEu5NgJwK93O9x6psY4xxAvlAfV8UqJRStVnHxnWYMaY/917Sk2eclxO1/Sd++up9v/fr8HsPFYjIrcCtAE2bNg1k10opZRm7Tbjp7Bbktv8bK19dRYO6CX7vsybhvg1Iq3A/1Xusqja5IuIAEoHfrK5jjJkCTAFIT08Pjo0GlVIqQFKTk0h94OuA9FWTaZlfgDYi0kJEIoFRwKxKbWYBf/DevgqYa05pnyqllFK+UO3I3RjjFJHbgTmAHXjNGPOriDwGZBhjZgGvAm+LSBawD88bgFJKKYvUaM7dGDMbmF3p2KMVbhcBV/u2NKWUUqdKv6GqlFIhSMNdKaVCkIa7UkqFIA13pZQKQRruSikVgsSq09FFJA/Ycoo/ngzs8WE5tZ2+HsfS1+MofS2OFQqvRzNjTEp1jSwL99MhIhnGmHSr6wgW+nocS1+Po/S1OFY4vR46LaOUUiFIw10ppUJQbQ33KVYXEGT09TiWvh5H6WtxrLB5PWrlnLtSSqkTq60jd6WUUidQ68K9us26w4WIpInIPBFZLSK/isidVtcUDETELiLLROQzq2uxmogkiciHIrJWRNaIyJlW12QVEbnb++9klYi8JyLRVtfkb7Uq3Gu4WXe4cAL3GGM6Av2AsWH8WlR0J7DG6iKCxDPAl8aY9kA3wvR1EZEmwF+AdGNMZzxLl4f8suS1Ktyp2WbdYcEYs8MYs9R7+xCef7iV97YNKyKSClwMvGJ1LVYTkURgIJ69FjDGlBhjDlhblaUcQIx3p7hYYLvF9fhdbQv3mmzWHXZEpDnQA/jJ2kos9zRwP+C2upAg0ALIA173TlO9IiJxVhdlBWPMNuApYCuwA8g3xnxlbVX+V9vCXVUiIvHAR8BdxpiDVtdjFRG5BNhtjFlidS1BwgH0BF4wxvQADgNh+RmViNTF8z/8FkBjIE5ErrO2Kv+rbeFek826w4aIROAJ9neNMTOsrsdi/YHhIrIZz3TduSLyjrUlWSoXyDXGlP1v7kM8YR+Ozgc2GWPyjDGlwAzgLItr8rvaFu412aw7LIiI4JlPXWOMmWR1PVYzxjxkjEk1xjTH8/dirjEm5Ednx2OM2QnkiEg776HzgNUWlmSlrUA/EYn1/rs5jzD4cLlGe6gGi+Nt1m1xWVbpD1wPZIrIcu+xv3r3u1UK4A7gXe9AKBu4weJ6LGGM+UlEPgSW4jnLbBlh8E1V/YaqUkqFoNo2LaOUUqoGNNyVUioEabgrpVQI0nBXSqkQpOGulFIhSMNdKaVCkIa7UkqFIA13pZQKQf8PUJZCjYoAXT8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(np.transpose(swp))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "f = 0\n",
    "a = 0.5 \n",
    "        \n",
    "(dt,lb,cs) = cn.average(Fish[f][\"Data\"], Fish[f][\"KLabel\"], Fish[f][\"Coordinates\"])\n",
    "\n",
    "# Network construction\n",
    "#=====================================================================\n",
    "# Calculate correlation coefficient\n",
    "#---------------------------------------------------------------------\n",
    "A  = np.corrcoef(dt)\n",
    "np.fill_diagonal(A, 0)\n",
    "\n",
    "# Z-score and binarise according to z cut off\n",
    "#---------------------------------------------------------------------\n",
    "bA = A > (a / 10)\n",
    "\n",
    "# Generate Network x object from this\n",
    "#---------------------------------------------------------------------\n",
    "G = nx.from_numpy_matrix(bA)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/roschkoenig/.local/lib/python3.6/site-packages/networkx/drawing/nx_pylab.py:611: MatplotlibDeprecationWarning: isinstance(..., numbers.Number)\n",
      "  if cb.is_numlike(alpha):\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.collections.LineCollection at 0x349d892e8>"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAADU9JREFUeJzt3GGI5Hd9x/H3xztTaYym9FaQu9Ok9NJ42ELSJU0Raoq2XPLg7oFF7iBYJXhgGylVhBRLlPjIhloQrtWTilXQGH0gC57cA40ExAu3ITV4FyLb03oXhawxzZOgMe23D2bSna53mX92Z3cv+32/4GD+//ntzJcfe++dndmZVBWSpO3vFVs9gCRpcxh8SWrC4EtSEwZfkpow+JLUhMGXpCamBj/JZ5M8meT7l7g+ST6ZZCnJo0lunP2YkqT1GvII/3PAgRe5/lZg3/jfUeBf1j+WJGnWpga/qh4Efv4iSw4Bn6+RU8DVSV4/qwElSbOxcwa3sRs4P3F8YXzup6sXJjnK6LcArrzyyj+8/vrrZ3D3ktTHww8//LOqmlvL184i+INV1XHgOMD8/HwtLi5u5t1L0stekv9c69fO4q90ngD2ThzvGZ+TJF1GZhH8BeBd47/WuRl4pqp+7ekcSdLWmvqUTpIvAbcAu5JcAD4CvBKgqj4FnABuA5aAZ4H3bNSwkqS1mxr8qjoy5foC/npmE0mSNoTvtJWkJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJamJQcFPciDJ40mWktx1kevfkOSBJI8keTTJbbMfVZK0HlODn2QHcAy4FdgPHEmyf9Wyvwfur6obgMPAP896UEnS+gx5hH8TsFRV56rqOeA+4NCqNQW8Znz5tcBPZjeiJGkWhgR/N3B+4vjC+NykjwK3J7kAnADef7EbSnI0yWKSxeXl5TWMK0laq1m9aHsE+FxV7QFuA76Q5Nduu6qOV9V8Vc3Pzc3N6K4lSUMMCf4TwN6J4z3jc5PuAO4HqKrvAq8Cds1iQEnSbAwJ/mlgX5Jrk1zB6EXZhVVrfgy8DSDJmxgF3+dsJOkyMjX4VfU8cCdwEniM0V/jnElyT5KD42UfBN6b5HvAl4B3V1Vt1NCSpJdu55BFVXWC0Yuxk+funrh8FnjLbEeTJM2S77SVpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDUxKPhJDiR5PMlSkrsuseadSc4mOZPki7MdU5K0XjunLUiyAzgG/BlwATidZKGqzk6s2Qf8HfCWqno6yes2amBJ0toMeYR/E7BUVeeq6jngPuDQqjXvBY5V1dMAVfXkbMeUJK3XkODvBs5PHF8Yn5t0HXBdku8kOZXkwMVuKMnRJItJFpeXl9c2sSRpTWb1ou1OYB9wC3AE+EySq1cvqqrjVTVfVfNzc3MzumtJ0hBDgv8EsHfieM/43KQLwEJV/aqqfgj8gNEPAEnSZWJI8E8D+5Jcm+QK4DCwsGrN1xg9uifJLkZP8Zyb4ZySpHWaGvyqeh64EzgJPAbcX1VnktyT5OB42UngqSRngQeAD1XVUxs1tCTppUtVbckdz8/P1+Li4pbctyS9XCV5uKrm1/K1vtNWkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgYFP8mBJI8nWUpy14use0eSSjI/uxElSbMwNfhJdgDHgFuB/cCRJPsvsu4q4G+Ah2Y9pCRp/YY8wr8JWKqqc1X1HHAfcOgi6z4GfBz4xQznkyTNyJDg7wbOTxxfGJ/7P0luBPZW1ddf7IaSHE2ymGRxeXn5JQ8rSVq7db9om+QVwCeAD05bW1XHq2q+qubn5ubWe9eSpJdgSPCfAPZOHO8Zn3vBVcCbgW8n+RFwM7DgC7eSdHkZEvzTwL4k1ya5AjgMLLxwZVU9U1W7quqaqroGOAUcrKrFDZlYkrQmU4NfVc8DdwIngceA+6vqTJJ7khzc6AElSbOxc8iiqjoBnFh17u5LrL1l/WNJkmbNd9pKUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpoYFPwkB5I8nmQpyV0Xuf4DSc4meTTJN5O8cfajSpLWY2rwk+wAjgG3AvuBI0n2r1r2CDBfVX8AfBX4h1kPKklanyGP8G8ClqrqXFU9B9wHHJpcUFUPVNWz48NTwJ7ZjilJWq8hwd8NnJ84vjA+dyl3AN+42BVJjiZZTLK4vLw8fEpJ0rrN9EXbJLcD88C9F7u+qo5X1XxVzc/Nzc3yriVJU+wcsOYJYO/E8Z7xuf8nyduBDwNvrapfzmY8SdKsDHmEfxrYl+TaJFcAh4GFyQVJbgA+DRysqidnP6Ykab2mBr+qngfuBE4CjwH3V9WZJPckOThedi/wauArSf49ycIlbk6StEWGPKVDVZ0ATqw6d/fE5bfPeC5J0oz5TltJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaGBT8JAeSPJ5kKcldF7n+N5J8eXz9Q0mumfWgkqT1mRr8JDuAY8CtwH7gSJL9q5bdATxdVb8L/BPw8VkPKklanyGP8G8ClqrqXFU9B9wHHFq15hDwb+PLXwXeliSzG1OStF47B6zZDZyfOL4A/NGl1lTV80meAX4b+NnkoiRHgaPjw18m+f5aht6GdrFqrxpzL1a4FyvcixW/t9YvHBL8mamq48BxgCSLVTW/mfd/uXIvVrgXK9yLFe7FiiSLa/3aIU/pPAHsnTjeMz530TVJdgKvBZ5a61CSpNkbEvzTwL4k1ya5AjgMLKxaswD85fjyXwDfqqqa3ZiSpPWa+pTO+Dn5O4GTwA7gs1V1Jsk9wGJVLQD/CnwhyRLwc0Y/FKY5vo65txv3YoV7scK9WOFerFjzXsQH4pLUg++0laQmDL4kNbHhwfdjGVYM2IsPJDmb5NEk30zyxq2YczNM24uJde9IUkm27Z/kDdmLJO8cf2+cSfLFzZ5xswz4P/KGJA8keWT8/+S2rZhzoyX5bJInL/VepYx8crxPjya5cdANV9WG/WP0Iu9/AL8DXAF8D9i/as1fAZ8aXz4MfHkjZ9qqfwP34k+B3xxffl/nvRivuwp4EDgFzG/13Fv4fbEPeAT4rfHx67Z67i3ci+PA+8aX9wM/2uq5N2gv/gS4Efj+Ja6/DfgGEOBm4KEht7vRj/D9WIYVU/eiqh6oqmfHh6cYvedhOxryfQHwMUafy/SLzRxukw3Zi/cCx6rqaYCqenKTZ9wsQ/aigNeML78W+MkmzrdpqupBRn/xeCmHgM/XyCng6iSvn3a7Gx38i30sw+5Lramq54EXPpZhuxmyF5PuYPQTfDuauhfjX1H3VtXXN3OwLTDk++I64Lok30lyKsmBTZtucw3Zi48Ctye5AJwA3r85o112XmpPgE3+aAUNk+R2YB5461bPshWSvAL4BPDuLR7lcrGT0dM6tzD6re/BJL9fVf+1pVNtjSPA56rqH5P8MaP3/7y5qv5nqwd7OdjoR/h+LMOKIXtBkrcDHwYOVtUvN2m2zTZtL64C3gx8O8mPGD1HubBNX7gd8n1xAVioql9V1Q+BHzD6AbDdDNmLO4D7Aarqu8CrGH2wWjeDerLaRgffj2VYMXUvktwAfJpR7Lfr87QwZS+q6pmq2lVV11TVNYxezzhYVWv+0KjL2JD/I19j9OieJLsYPcVzbjOH3CRD9uLHwNsAkryJUfCXN3XKy8MC8K7xX+vcDDxTVT+d9kUb+pRObdzHMrzsDNyLe4FXA18Zv27946o6uGVDb5CBe9HCwL04Cfx5krPAfwMfqqpt91vwwL34IPCZJH/L6AXcd2/HB4hJvsToh/yu8esVHwFeCVBVn2L0+sVtwBLwLPCeQbe7DfdKknQRvtNWkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJauJ/Acz2XLpusNoKAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "cs = Fish[f]['Coordinates'][:,0:2]\n",
    "ps = {}\n",
    "for c in range(cs.shape[0]): ps.update({c:cs[c,:]})\n",
    "nx.draw_networkx_edges(G, pos = ps)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Working on fish file 0\n",
      "Working on fish file 1\n",
      "Working on fish file 2\n",
      "Working on fish file 3\n"
     ]
    }
   ],
   "source": [
    "# Load actual datasets - GABRG2 \n",
    "#---------------------------------------------------------------------\n",
    "conditions = {'gb2':'GABRG2 Knockout'}\n",
    "spatial_conversion = [.5, .5, 15]\n",
    "Fish       = []\n",
    "\n",
    "for key in conditions.keys():\n",
    "    tracefiles = glob.glob(F[key] + os.sep + '*realtrace.npy')\n",
    "    coordfiles = glob.glob(F[key] + os.sep + '*realcoord.npy')\n",
    "    tracefiles.sort()\n",
    "    coordfiles.sort()\n",
    "    \n",
    "    for k in range(len(tracefiles)):\n",
    "        try: \n",
    "            print('Working on fish file ' +str(k))\n",
    "            fish = {}\n",
    "            fish['Voxels']      = np.load(coordfiles[k])\n",
    "            fish['Coordinates'] = np.multiply(fish['Voxels'] , spatial_conversion)\n",
    "            fish['Data']        = np.load(tracefiles[k])\n",
    "            fish['Condition']   = conditions[key]\n",
    "            Fish.append(fish)\n",
    "\n",
    "            if fish['Coordinates'].shape[0] != fish['Data'].shape[0]: \n",
    "                raise Exception('The data and coordinates do not match in length for: \\n \\n' \n",
    "                                + coordfiles[k] + ' \\n and \\n' + tracefiles[k])\n",
    "        except: \n",
    "            print('Skipping ' + tracefiles[k])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Spatial k-Means for Fish 0\n",
      "Spatial k-Means for Fish 1\n",
      "Spatial k-Means for Fish 2\n",
      "Spatial k-Means for Fish 3\n"
     ]
    }
   ],
   "source": [
    "for f in range(len(Fish)): \n",
    "    print('Spatial k-Means for Fish ' + str(f))\n",
    "    Fish[f] = cn.spacek(Fish[f])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n",
      "Now working on fish 1: GABRG2 Knockout\n",
      "Now working on fish 2: GABRG2 Knockout\n",
      "Now working on fish 3: GABRG2 Knockout\n",
      "Now working on fish 4: GABRG2 Knockout\n"
     ]
    }
   ],
   "source": [
    "arange = []\n",
    "for a in range(10):\n",
    "    cd = []\n",
    "    for f in range(len(Fish)):\n",
    "        print('Now working on fish ' + str(f+1) + ': ' + Fish[f][\"Condition\"]) \n",
    "        (dt,lb,cs) = cn.average(Fish[f][\"Data\"], Fish[f][\"KLabel\"], Fish[f][\"Coordinates\"])\n",
    "\n",
    "        # Network construction\n",
    "        #=====================================================================\n",
    "        # Calculate correlation coefficient\n",
    "        #---------------------------------------------------------------------\n",
    "        A  = np.corrcoef(dt)\n",
    "        np.fill_diagonal(A, 0)\n",
    "\n",
    "        # Z-score and binarise according to z cut off\n",
    "        #---------------------------------------------------------------------\n",
    "    #     A  = stats.zscore(np.ravel(A)).reshape(A.shape)\n",
    "        bA = A > (a / 10)\n",
    "\n",
    "        # Generate Network x object from this\n",
    "        #---------------------------------------------------------------------\n",
    "        G = nx.from_numpy_matrix(bA)\n",
    "\n",
    "        # Calculate cycle density\n",
    "        #=====================================================================\n",
    "        cycs = nx.algorithms.cycle_basis(G)\n",
    "        no_E = G.number_of_edges()\n",
    "        if no_E > 0: ncd = len(nx.algorithms.cycle_basis(G)) / no_E\n",
    "        else:        ncd = 0\n",
    "\n",
    "        cd.append(ncd)\n",
    "        \n",
    "    arange.append(cd)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 57,
   "metadata": {},
   "outputs": [],
   "source": [
    "swp = np.zeros((1,10))\n",
    "ai  = 0\n",
    "for a in arange:\n",
    "    swp[0,ai] = np.mean(a)\n",
    "    ai = ai + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'cd' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-5-91fec7bf6847>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mplt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mviolinplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcd\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;36m9\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mshowmeans\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;32mTrue\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      2\u001b[0m \u001b[0mplt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mviolinplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcd\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m10\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;36m19\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mshowmeans\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;32mTrue\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0mplt\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mviolinplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcd\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m20\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mshowmeans\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;32mTrue\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mp1\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mstats\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mmannwhitneyu\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcd\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;36m9\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mcd\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m10\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;36m19\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'cd' is not defined"
     ]
    }
   ],
   "source": [
    "plt.violinplot(cd[0:9],[1], showmeans = True)\n",
    "plt.violinplot(cd[10:19],[2], showmeans = True)\n",
    "plt.violinplot(cd[20:],[3], showmeans = True)\n",
    "\n",
    "[s, p1] = stats.mannwhitneyu(cd[0:9], cd[10:19])\n",
    "[s, p2] = stats.mannwhitneyu(cd[10:19], cd[20:])\n",
    "[s, p3] = stats.mannwhitneyu(cd[0:9], cd[20:])\n",
    "print(str(p1) +' '+str(p2)+' '+str(p3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x1ae0da2898>]"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAHadJREFUeJzt3Xl4VOXd//H3dyZ7CARI2BIggQQEVySioiiCbREVfo+tLVQrdaOu2O2xWrv97Pb0sa2VilpcqtgqrdZaqlhURMUVgriwiARkRwhLQoDsuZ8/EjHEQAaYcM/yeV1XLnLOuTPz4VzwmTtn5pxjzjlERCS2BHwHEBGR8FO5i4jEIJW7iEgMUrmLiMQglbuISAxSuYuIxCCVu4hIDFK5i4jEIJW7iEgMSvD1xFlZWS4vL8/X04uIRKVFixZtc85ltzXOW7nn5eVRXFzs6+lFRKKSma0NZZwOy4iIxCCVu4hIDFK5i4jEIJW7iEgMarPczewhM9tqZksOsN3MbKqZlZjZ+2Z2cvhjiojIoQhl5v4wMOYg288DCpu+JgP3HnksERE5Em2Wu3PuVWDHQYaMB2a4Rm8BmWbWM1wBRUTk0IXjc+45wPpmyxua1m0Ow2N/zsI1O5j/USmJwQAJwQCJQSMxGGhaNpKa/kxsti0hECApwUgIBPZfv2/8/o8TDFh7RBcROWrCUe6tNWGrN2Y1s8k0HrqhT58+h/Vk76zdydSXSg7rZ0NlRmPRB4zEhKYXh6B97kUgKSFAWlKQ9KQE0pKCpCUHSWv6Pj0pgdSkIOnJQVITE0hPDjaOSUrYb1tKQpCAXkxEJMzCUe4bgN7NlnOBTa0NdM5NB6YDFBUVHdadub91dn8mn9WPugZHXb2jpr6BuvoGausdtfUN1NY3UNfgqKlr/PPTdbX1rmncZ2M/9/MNDdTWOeoaGprWuxbj93+e6roGKqrq2LKrir019U1fdVTVNhzS3+nT0m/8s/ErPTmB1MTGPz9b/+mLSAJpicGmF4wEuqQnkZOZSmZaImZ6oRCR8JT7LOAGM5sJnAqUO+fa5ZDMp8ysaQYNqQTb86kOS32DY29NHZU19expKvx95V9dx56aeipr6pq2Na7bW9t8Wz27q+sorahmz6ePU11PZW39QZ83NTFIr8wUcjqnkZOZQq9OqfTKTCWncyo5mal075hCUoI+/SoSD9osdzN7HBgJZJnZBuCnQCKAc+4+YDYwFigB9gKXt1fYaBEMGBkpiWSkJIb1cRsaHJW19fsV/t6aOrbtrmZjWRWbyirZVFbJxrJKlm0qZ9vumv1+3gy6ZSQ3Fn7TV6+mr0+XO6YmaPYvEgPaLHfn3MQ2tjvg+rAlkgMKBIz05ATSk0P7hauqtp7N5Y2lv3FnY+lvKqtkU3klSzaW8/zSLdTU738IKT0puG+232vfC0DjbwE5nRtn/4lBzf5FIp23q0JK+0tJDJKflU5+Vnqr2xsaHNv31Ow349/3AlBWxQcbytm+Z//Zf8Cge8eUZsWfSk5mCjmdUzm9XxapSZF3mEwkHqnc41ggYGRnJJOdkcyJvTNbHVNZU8+m8srPXgB2Vu47BPTehjL+s+STfbP//tnpTLvkZI7p0fFo/jVEpBUqdzmo1KQg/bM70D+7Q6vbGxoc23ZX8+76Mm57egnj736dn407lgmn9NaxexGPdPBUjkggYHTrmMIXj+3BczeNYFh+F2596gOmzHyXiqpa3/FE4pbKXcImq0Myj1w+jJvHDGT2B5u58I+vsWRjue9YInFJ5S5hFQgY140sYObk06iua+Cie97gkTfW0PihKhE5WlTu0i5OyevC7CkjGFGYxU9nLeWavyyifK8O04gcLSp3aTed05N4YFIRPzp/EHOXb+X8P85n8bqdvmOJxAWVu7QrM+OqEf148trhAFx835vc/+pqHaYRaWcqdzkqTuqdybNTRnDuoO78cvZyrnqkmJ0tTpASkfBRuctR0yk1kXsvPZn/P+5Y5q/cxtip81m45mD3gRGRw6Vyl6PKzJg0PI+nrhtOckKACdPfYtq8EhoadJhGJJxU7uLFcTmd+PeNZzL2+J7cMWcFk/68gNKKat+xRGKGyl28yUhJZOqEk/j1Rcez4OMdjJ06nzdKtvmOJRITVO7ilZkxcVgf/nXDGXRMSeCSB9/mzhc+ol6HaUSOiMpdIsIxPTry7xvP5KIhudw1dyWXPPAWW3ZV+Y4lErVU7hIx0pIS+N1XT+S3F5/Ie+vLGXvXfF75qNR3LJGopHKXiPOVobn8+8YzyM5IZtJDC/jNfz6krv7QbjouEu9U7hKRCrpl8PT1ZzBxWB/ufXkVE6a/xaaySt+xRKKGyl0iVkpikF9fdDxTJw5h+eZdjJ06nxeXbfEdSyQqqNwl4o07sRfPTBlBTmYqV80o5hfPLKOmTodpRA5G5S5RIT8rnaeuG843h+fxwGsfc/Gf3mT9jr2+Y4lELJW7RI3khCA/G3cs9116MqtLdzN26nye+2Cz71giEUnlLlFnzHE9mT1lBP2yO3DtX9/hJ/9aQlVtve9YIhFF5S5RqXeXNJ741ulcPSKfGW+u5cv3vsHH2/b4jiUSMVTuErWSEgLcdv5gHpxUxMaySi6YOp9/vbvRdyyRiKByl6g3elB3Zk8ZwaCeHblp5rv8Z4mOw4uo3CUm9MpMZebk0yjo1oE7X1ip68NL3FO5S8xICAa4cVQBK7ZU8PyyT3zHEfFK5S4x5YITetEvK5275pboJtwS11TuElOCAeP6cwpYvnkXL+hSBRLHQip3MxtjZivMrMTMbmllex8zm2dmi83sfTMbG/6oIqEZf1Iv+nZNY+pLKzV7l7jVZrmbWRCYBpwHDAYmmtngFsN+BPzdOTcEmADcE+6gIqFKCAa4/pwClmzcxbwVW33HEfEilJn7MKDEObfaOVcDzATGtxjjgI5N33cCNoUvosih+68hOfTukspdL2r2LvEplHLPAdY3W97QtK65nwGXmtkGYDZwY2sPZGaTzazYzIpLS3WHHWk/icEA148s4L0N5bqbk8SlUMrdWlnXcio0EXjYOZcLjAUeNbPPPbZzbrpzrsg5V5SdnX3oaUUOwUUn55KTmcpdczV7l/gTSrlvAHo3W87l84ddrgT+DuCcexNIAbLCEVDkcCUlBLh2ZH8WryvjtZJtvuOIHFWhlPtCoNDM8s0sicY3TGe1GLMOGA1gZoNoLHf9LizeXVyUS89OKTr2LnGnzXJ3ztUBNwBzgOU0fipmqZndbmbjmoZ9D7jazN4DHge+6fQ/SSJAckKQa0f2p3jtTt5cvd13HJGjxnx1cFFRkSsuLvby3BJfqmrrOfuOeeR1Tedv3zrddxyRI2Jmi5xzRW2N0xmqEvNSEoNcc3Z/3v54B29p9i5xQuUucWHisD5kdUjmjy+t9B1F5KhQuUtcaJy99+P1ku0Ur9nhO45Iu1O5S9y45NS+ZHVI4q65mr1L7FO5S9xITQpy9Yh+zF+5jXfW7fQdR6Rdqdwlrlx6Wl+6pCcxVbN3iXEqd4kr6ckJXDUin5dXlPLe+jLfcUTajcpd4s5lp+eRmZao2bvENJW7xJ0OyQlceUY+cz/cypKN5b7jiLQLlbvEpUln5NExJUGzd4lZKneJSx1TErnizHyeX7aFZZt2+Y4jEnYqd4lblw/PJyM5QWetSkxSuUvc6pSWyOVn5PHckk9Y8UmF7zgiYaVyl7h2xZn5pCcFmarZu8QYlbvEtcy0JCYNz2P2B5tZuUWzd4kdKneJe1eN6EdqYpC755X4jiISNip3iXtd0pP4xul9+fd7m1hVutt3HJGwULmLAFeP6EdyQpBpL2n2LrFB5S4CZHVI5tLT+vD0uxtZs22P7zgiR0zlLtLk6rP6kRgM6Ni7xASVu0iTbhkpfP3UPvxz8UbWbd/rO47IEVG5izRzzdn9CQaMe17W7F2im8pdpJnuHVOYeEpvnly0gfU7NHuX6KVyF2nhmpH9CZhx7yurfEcROWwqd5EWenZK5aun5PJE8Xo2lVX6jiNyWFTuIq24dmQBAPe+rNm7RCeVu0grcjJT+crQXP62cD2flFf5jiNyyFTuIgdw3cgCGpzjPh17lyikchc5gN5d0rjo5BweX7COrbs0e5foonIXOYjrzymgrsHxp1dX+44ickhCKnczG2NmK8ysxMxuOcCYr5rZMjNbamaPhTemiB99u6bz/07K4a9vr6W0otp3HJGQtVnuZhYEpgHnAYOBiWY2uMWYQuBW4Azn3LHAt9shq4gXN4wqoKaugfvna/Yu0SOUmfswoMQ5t9o5VwPMBMa3GHM1MM05txPAObc1vDFF/MnPSmf8STk8+uZatu/W7F2iQyjlngOsb7a8oWldcwOAAWb2upm9ZWZjwhVQJBJcf04BVXX1PPDax76jiIQklHK3Vta5FssJQCEwEpgIPGBmmZ97ILPJZlZsZsWlpaWHmlXEm4JuHbjghF7MeGMNO/fU+I4j0qZQyn0D0LvZci6wqZUx/3LO1TrnPgZW0Fj2+3HOTXfOFTnnirKzsw83s4gXN44qYG9tPQ9q9i5RIJRyXwgUmlm+mSUBE4BZLcY8DZwDYGZZNB6m0btPElMGdM9g7HE9efiNNZTt1exdIlub5e6cqwNuAOYAy4G/O+eWmtntZjauadgcYLuZLQPmAf/tnNveXqFFfLlxdAG7q+t46PU1vqOIHJQ51/Lw+dFRVFTkiouLvTy3yJG45tFFvL5qG6/9YBSdUhN9x5E4Y2aLnHNFbY3TGaoih+jG0QVUVNXxsGbvEsFU7iKH6NhenTh3UHcefG01FVW1vuOItErlLnIYbhpdyK6qOma8udZ3FJFWqdxFDsPxuZ0YdUw37p+/mt3Vdb7jiHyOyl3kME0ZXUjZ3loe1exdIpDKXeQwndQ7k7MHZHP//NXsrdHsXSKLyl3kCEwZXciOPTX85S3N3iWyqNxFjsDQvp05syCL6a+uprKm3ncckX1U7iJH6KZzC9m2u4bHFqzzHUVkH5W7yBE6Ja8Lp/fryn2vrKKqVrN3iQwqd5EwmDK6kNKKamZq9i4RQuUuEgan9+/KsPwu3KvZu0QIlbtImNw0upAtu6p5onh924NF2pnKXSRMhvfvytC+nbnn5VVU12n2Ln6p3EXCxMy4aXQhm8ur+Meijb7jSJxTuYuE0YjCLE7qncm0eSXU1DX4jiNxTOUuEkZmxrfPLWRjWaXOWhWvVO4iYXb2gGxGFGZx19yVuteqeKNyFwkzM+O28wdRUVXLH18q8R1H4pTKXaQdHNOjI18t6s2MN9fw8bY9vuNIHFK5i7ST735xAInBAP/z3HLfUSQOqdxF2km3jBSuG9mfOUu38Pbq7b7jSJxRuYu0oyvP7EfPTin84tnlNDQ433EkjqjcRdpRalKQm8cM5ION5Tz9rk5skqNH5S7SzsafmMMJuZ24Y84K3dBDjhqVu0g7CwSM28YOYnN5FQ/MX+07jsQJlbvIUXBqv66MObYH976yiq27qnzHkTigchc5Sm457xhq6xv4/Qsf+Y4icUDlLnKU5GWlc9npefyteD3LNu3yHUdinMpd5CiaMqqQTqmJ/Gr2cpzTRyOl/ajcRY6iTmmJTBlVyGsl23h5RanvOBLDQip3MxtjZivMrMTMbjnIuK+YmTOzovBFFIktl57Wl/ysdH7x7DJq63XNd2kfbZa7mQWBacB5wGBgopkNbmVcBjAFeDvcIUViSVJCgFvPO4ZVpXuYuWCd7zgSo0KZuQ8DSpxzq51zNcBMYHwr434O/C+gz3mJtOELg7tzan4X7nxxJbuqan3HkRgUSrnnAM1v576had0+ZjYE6O2ceyaM2URilpnx4wsGs3NvDdPm6ZrvEn6hlLu1sm7f2/xmFgDuBL7X5gOZTTazYjMrLi3Vm0kS347L6cRFQ3L582trWL9jr+84EmNCKfcNQO9my7nApmbLGcBxwMtmtgY4DZjV2puqzrnpzrki51xRdnb24acWiRH//aWBBALwm/986DuKxJhQyn0hUGhm+WaWBEwAZn260TlX7pzLcs7lOefygLeAcc654nZJLBJDenRKYfJZ/Xnm/c0sWrvTdxyJIW2Wu3OuDrgBmAMsB/7unFtqZreb2bj2DigS6751Vj+6ZSTzi2eX6cQmCZuQPufunJvtnBvgnOvvnPtl07qfOOdmtTJ2pGbtIqFLT07g+18cyOJ1ZTzz/mbfcSRG6AxVkQjw5aG5DOrZkf957kOqanXNdzlyKneRCBAMGD86fxAbyyr58+trfMeRGKByF4kQZxRkMfqYbtwzr4Rtu6t9x5Eop3IXiSC3jh3E3tp6/vCirvkuR0blLhJBCrp14NJT+/DY2+tYuaXCdxyJYip3kQhz07kDSE9O4Fezl/uOIlFM5S4SYbqkJ3HjqALmrSjl1Y90mQ45PCp3kQg0aXgevbuk8qvZy6lv0IlNcuhU7iIRKDkhyC1jBvHhJxU8Uby+7R8QaUHlLhKhxh7fg6F9O/Pb5z9id3Wd7zgSZVTuIhHKrPHEpm27q7nv5VW+40iUUbmLRLAhfToz7sRe3D9/NZvKKn3HkSiicheJcDePGYgD7pizwncUiSIqd5EIl9s5javOzOefizfy3voy33EkSqjcRaLAtSP7k9UhiV8+u1zXfJeQqNxFokBGSiLf+cIAFqzZwZyln/iOI1FA5S4SJb5W1JsB3Tvw6+c+pKauwXcciXAqd5EokRAM8MOxg1i7fS8z3lzjO45EOJW7SBQZObAbZw3IZurclezcU+M7jkQwlbtIlLlt7CB2V9cx9aWVvqNIBFO5i0SZgT0y+NopfXj0zbWsLt3tO45EKJW7SBT67hcGkJwQ4NfPfeg7ikQolbtIFMrOSOa6cwp4YdkW3ly13XcciUAqd5EodeWZ+eRkpvKLZ5fRoGu+Swsqd5EolZIY5OYxA1m6aRdPLd7oO45EGJW7SBS78IRenNg7k9/OWcHeGl3zXT6jcheJYoGA8ePzB/HJriruf/Vj33EkgqjcRaJcUV4Xxh7fg/teWcWWXVW+40iEULmLxIAfjDmG+gbH757XNd+lkcpdJAb07ZrOpOF9eWLRBpZuKvcdRyKAyl0kRtwwqpDM1ERd812AEMvdzMaY2QozKzGzW1rZ/l0zW2Zm75vZXDPrG/6oInIwnVIT+fa5A3hj1XZe+nCr7zjiWZvlbmZBYBpwHjAYmGhmg1sMWwwUOedOAJ4E/jfcQUWkbV8/tQ/9stP55ezl1Nbrmu/xLJSZ+zCgxDm32jlXA8wExjcf4Jyb55zb27T4FpAb3pgiEorEYIAfnjeI1aV7eOztdb7jiEehlHsOsL7Z8oamdQdyJfBcaxvMbLKZFZtZcWlpaegpRSRkowd14/R+XfnDix9RXlnrO454Ekq5WyvrWn23xswuBYqAO1rb7pyb7pwrcs4VZWdnh55SREJmZtx2/iDKKmu5W9d8j1uhlPsGoHez5VxgU8tBZnYucBswzjlXHZ54InI4jsvpxMVDc3no9TXMXb7FdxzxIJRyXwgUmlm+mSUBE4BZzQeY2RDgTzQWu96mF4kAP73wWI7t1ZHrH3uHRWt3+o4jR1mb5e6cqwNuAOYAy4G/O+eWmtntZjauadgdQAfgCTN718xmHeDhROQoSU9O4KFvnkKPjilc+chCSrZW+I4kR5H5OtmhqKjIFRcXe3lukXiybvteLrr3DZKCxj+uG07PTqm+I8kRMLNFzrmitsbpDFWRGNenaxoPX34Ku6rqmPTQAsr36hM08UDlLhIHjsvpxPRvDOXjbXu4asZCqmrrfUeSdqZyF4kTwwuy+P1XT6J47U6mPL6Yet2aL6ap3EXiyIUn9uInFwzm+WVb+PG/lugCYzEswXcAETm6Lj8jn60V1dz78iq6ZSTz7XMH+I4k7UDlLhKHbv7SQEorqvnDiyvJzkjmklN1IddYo3IXiUNmxq8vOp7tu6v58dNLyOqQzJeO7eE7loSRjrmLxKnEYIBpl5zMCbmZ3Pj4YhZ8vMN3JAkjlbtIHEtLajyLNbdzKlc9spAVn+gs1lihcheJc13Sk5hxxTBSk4JMemgBG8sqfUeSMFC5iwi5ndN45Iph7Kmp47IH32bnnhrfkeQIqdxFBIBjenTkgcuKWL+zkiseWUhljc5ijWYqdxHZ59R+XZk64STeXV/GDY+9Q53uwxq1VO4isp8xx/Xk9vHHMffDrfzwnx/oLNYopc+5i8jnfOO0vpTuqmLqSyV0y0jh+18a6DuSHCKVu4i06jtfGEDp7mrunldCdkYyk4bn+Y4kh0DlLiKtMjN+Pv44Sitq+Nm/l5LVIZnzT+jpO5aESMfcReSAEoIB7v76EIb26cx3/vYub6za5juShEjlLiIHlZIY5IFJRfTtmsa3Zixi2aZdviNJCFTuItKmzLQkHrliGB1SEpj05wWs37HXdyRpg8pdRELSKzOVGVcMo6augcseWsD23dW+I8lBqNxFJGSF3TN46JtFbCqr5IqHF7Knus53JDkAlbuIHJKhfbtw99dP5oON5Vz713eo1VmsEUnlLiKH7AuDu/Or/zqeVz8q5QdPvk+DbrYdcfQ5dxE5LBOG9aG0oprfvfAR2R2TufW8Qb4jSTMqdxE5bDeMKqB0dzV/emU13TJSuPLMfN+RpInKXUQOm5nx0wuPpbSimp8/s4ysDkmMPynHdyxBx9xF5AgFA8adXzuJU/O78P0n3mP+ylLfkQSVu4iEQUpikOmXFdE/uwPXPLqIDzaU+44U91TuIhIWnVITeeSKYWSmJXH5wwtYu32P70hxLaRyN7MxZrbCzErM7JZWtieb2d+atr9tZnnhDioika97xxRmXDmM+gbHNx5cQGmFzmL1pc1yN7MgMA04DxgMTDSzwS2GXQnsdM4VAHcCvwl3UBGJDv2zO/DQN0+htKKayx9ewG6dxepFKDP3YUCJc261c64GmAmMbzFmPPBI0/dPAqPNzMIXU0SiyZA+nbnn0pNZvrmCax5dRE2dzmI92kL5KGQOsL7Z8gbg1AONcc7VmVk50BXQxZ9F4tQ5A7vxmy+fwPefeI+Rd8wjPVmfvP7UlNGFXHhir3Z9jlD2dmsz8JbnGocyBjObDEwG6NOnTwhPLSLR7CtDczFg7odbfEeJKJ1SE9v9OUIp9w1A72bLucCmA4zZYGYJQCdgR8sHcs5NB6YDFBUV6WIUInHgy0Nz+fLQXN8x4k4ox9wXAoVmlm9mScAEYFaLMbOASU3ffwV4yTmn8hYR8aTNmXvTMfQbgDlAEHjIObfUzG4Hip1zs4AHgUfNrITGGfuE9gwtIiIHF9I7HM652cDsFut+0uz7KuDi8EYTEZHDpTNURURikMpdRCQGqdxFRGKQyl1EJAap3EVEYpD5+ji6mZUCaw/zx7PQpQ2a0/7Yn/bHZ7Qv9hcL+6Ovcy67rUHeyv1ImFmxc67Id45Iof2xP+2Pz2hf7C+e9ocOy4iIxCCVu4hIDIrWcp/uO0CE0f7Yn/bHZ7Qv9hc3+yMqj7mLiMjBRevMXUREDiLqyr2tm3XHCzPrbWbzzGy5mS01s5t8Z4oEZhY0s8Vm9ozvLL6ZWaaZPWlmHzb9OznddyZfzOw7Tf9PlpjZ42aW4jtTe4uqcg/xZt3xog74nnNuEHAacH0c74vmbgKW+w4RIe4C/uOcOwY4kTjdL2aWA0wBipxzx9F46fKYvyx5VJU7od2sOy445zY7595p+r6Cxv+4OX5T+WVmucD5wAO+s/hmZh2Bs2i81wLOuRrnXJnfVF4lAKlNd4pL4/N3k4s50Vburd2sO64LDcDM8oAhwNt+k3j3B+BmoMF3kAjQDygF/tx0mOoBM0v3HcoH59xG4LfAOmAzUO6ce95vqvYXbeUe0o2444mZdQD+AXzbObfLdx5fzOwCYKtzbpHvLBEiATgZuNc5NwTYA8Tle1Rm1pnG3/DzgV5Aupld6jdV+4u2cg/lZt1xw8wSaSz2vzrnnvKdx7MzgHFmtobGw3WjzOwvfiN5tQHY4Jz79Le5J2ks+3h0LvCxc67UOVcLPAUM95yp3UVbuYdys+64YGZG4/HU5c653/vO45tz7lbnXK5zLo/GfxcvOedifnZ2IM65T4D1ZjawadVoYJnHSD6tA04zs7Sm/zejiYM3l0O6h2qkONDNuj3H8uUM4BvAB2b2btO6Hzbd71YE4Ebgr00TodXA5Z7zeOGce9vMngTeofFTZouJgzNVdYaqiEgMirbDMiIiEgKVu4hIDFK5i4jEIJW7iEgMUrmLiMQglbuISAxSuYuIxCCVu4hIDPo/zQofgwX6ZFUAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(np.transpose(swp))"
   ]
  },
  {
   "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.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}