test pour g2g.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/anaconda3/envs/G2G2/lib/python3.7/site-packages/unagi/plotting.py:30: MatplotlibDeprecationWarning: You are modifying the state of a globally registered colormap. In future versions, you will not be able to modify a registered colormap in-place. To remove this warning, you can make a copy of the colormap first. cmap = copy.copy(mpl.cm.get_cmap(\"viridis\"))\n",
      "  IMG_CMAP.set_bad(color='black')\n"
     ]
    }
   ],
   "source": [
    "\n",
    "\n",
    "%pylab inline\n",
    "from galaxy2galaxy import problems\n",
    "from astropy.visualization import make_lupton_rgb\n",
    "import tensorflow as tf\n",
    "Modes = tf.estimator.ModeKeys"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['attrs2img_cosmos',\n",
       " 'attrs2img_cosmos128',\n",
       " 'attrs2img_cosmos128_euclid',\n",
       " 'attrs2img_cosmos32',\n",
       " 'attrs2img_cosmos_psf_euclide',\n",
       " 'img2img_cosmos',\n",
       " 'img2img_cosmos128',\n",
       " 'img2img_cosmos32',\n",
       " 'img2img_cosmos_hsc',\n",
       " 'img2img_hsc',\n",
       " 'img2img_hsc_anomaly',\n",
       " 'img2img_hsc_anomaly_small',\n",
       " 'img2photoz_hsc']"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# To list the problems currently available in G2G\n",
    "problems.available()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\n",
    "# Let's create an instance of the hsc_problem\n",
    "problem128 = problems.problem('attrs2img_cosmos_psf_euclide')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def wiener_tf(image, psf, balance):\n",
    "    trans_func = tf.signal.rfft2d(tf.signal.ifftshift(psf))\n",
    "    wiener_filter = tf.math.conj(trans_func) / (tf.dtypes.cast(tf.math.abs(trans_func),'complex64') ** 2 + balance)\n",
    "    deconv = tf.signal.irfft2d(wiener_filter * tf.signal.rfft2d(image),fft_length=image.shape)\n",
    "    deconv = tf.keras.backend.clip(deconv, -1, 1)\n",
    "\n",
    "    return deconv\n",
    "\n",
    "def tikhonov_tf(img, psf, tau):\n",
    "    if len(img.shape) == 2:\n",
    "        return(wiener_tf(img, psf, tau))\n",
    "    tikho_list = []\n",
    "    for i in range(img.shape[0]):\n",
    "        deconvolved = wiener_tf(img[i], psf[i], tau[i])\n",
    "        tikho_list += [deconvolved]\n",
    "    return tikho_list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From /opt/anaconda3/envs/G2G2/lib/python3.7/site-packages/tensor2tensor/data_generators/problem.py:651: The name tf.logging.info is deprecated. Please use tf.compat.v1.logging.info instead.\n",
      "\n",
      "INFO:tensorflow:Reading data files from /Users/hk263532/Documents/galaxy2galaxy/data/attrs2img_cosmos_psf_euclide/attrs2img_cosmos_psf_euclide-train*\n",
      "INFO:tensorflow:partition: 0 num_data_files: 1\n",
      "WARNING:tensorflow:From /opt/anaconda3/envs/G2G2/lib/python3.7/site-packages/tensor2tensor/data_generators/problem.py:687: parallel_interleave (from tensorflow.python.data.experimental.ops.interleave_ops) is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "Use `tf.data.Dataset.interleave(map_func, cycle_length, block_length, num_parallel_calls=tf.data.experimental.AUTOTUNE)` instead. If sloppy execution is desired, use `tf.data.Options.experimental_determinstic`.\n",
      "WARNING:tensorflow:From /Users/hk263532/Documents/galaxy2galaxy/galaxy2galaxy/data_generators/galsim_utils.py:75: The name tf.FixedLenFeature is deprecated. Please use tf.io.FixedLenFeature instead.\n",
      "\n",
      "WARNING:tensorflow:\n",
      "The TensorFlow contrib module will not be included in TensorFlow 2.0.\n",
      "For more information, please see:\n",
      "  * https://github.com/tensorflow/community/blob/master/rfcs/20180907-contrib-sunset.md\n",
      "  * https://github.com/tensorflow/addons\n",
      "  * https://github.com/tensorflow/io (for I/O related ops)\n",
      "If you depend on functionality not listed there, please file an issue.\n",
      "\n"
     ]
    }
   ],
   "source": [
    "dset = problem128.dataset(Modes.TRAIN, data_dir='/sparseastro/Euclid/ShapeDeconv/attrs2img_cosmos_psf_euclide')\n",
    "dset2 = problem128.dataset(Modes.TRAIN, data_dir='/sparseastro/Euclid/ShapeDeconv/attrs2img_cosmos_psf_euclide')\n",
    "#from tensorflow.compat.v1.spectral import irfft2d\n",
    "def pre_proc_unet(dico):\n",
    "    print(dico['inputs'].shape)\n",
    "    print(dico['psf'].shape)\n",
    "    dico['psf'] = tf.signal.irfft2d(tf.cast(dico['psf'], complex64)[...,0])\n",
    "    print(dico['psf'].shape)\n",
    "    dico['psf'] = tf.signal.fftshift(dico['psf'])\n",
    "    dico['inputs'] = wiener_tf(tf.compat.v1.image.resize_image_with_pad(dico['inputs'], 128, 128)[...,0], dico['psf'], 0.05)\n",
    "    dico['psf'] = tf.expand_dims(dico['psf'], axis=-1)\n",
    "    dico['inputs'] = tf.expand_dims(dico['inputs'], axis=-1)\n",
    "    dico['inputs'] = tf.image.crop_to_bounding_box(dico['inputs'],33,33, 64,64)\n",
    "    print(dico['psf'].shape)\n",
    "    print(dico['inputs'].shape)\n",
    "    return dico\n",
    "dset = dset.map(pre_proc_unet)\n",
    "\n",
    "dset2 = dset2.batch(1)\n",
    "dset = dset.batch(1)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<DatasetV1Adapter shapes: {inputs: (?, 64, 64, 1), inputs2: (?, 64, 64, 1), ps: (?, 64, 33), ps2: (?, 64, 33), psf: (?, 128, 65, 1), psf2: (?, 128, 65, 1), targets: (?, 64, 64, 1)}, types: {inputs: tf.float32, inputs2: tf.float32, ps: tf.float32, ps2: tf.float32, psf: tf.float32, psf2: tf.float32, targets: tf.float32}>"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From <ipython-input-6-91ef2c814408>:2: DatasetV1.make_one_shot_iterator (from tensorflow.python.data.ops.dataset_ops) is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "Use `for ... in dataset:` to iterate over a dataset. If using `tf.estimator`, return the `Dataset` object directly from your input function. As a last resort, you can use `tf.compat.v1.data.make_one_shot_iterator(dataset)`.\n"
     ]
    }
   ],
   "source": [
    "# Build an iterator over this dataset, and extract a batch\n",
    "it = dset.make_one_shot_iterator().get_next()\n",
    "it2 = dset2.make_one_shot_iterator().get_next()\n",
    "\n",
    "sess = tf.Session()\n",
    "batch = sess.run(it)\n",
    "batch2 = sess.run(it2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensorflow.python.framework.ops.Tensor"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "type(it['inputs'])\n",
    "#it.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "ename": "IndexError",
     "evalue": "index 1 is out of bounds for axis 0 with size 1",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mIndexError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-8-534d8a2fa3aa>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      2\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mi\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m25\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[1;32m      3\u001b[0m     \u001b[0msubplot\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m5\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m5\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mi\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 4\u001b[0;31m     \u001b[0mimshow\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mbatch\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'inputs'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mi\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[0m\u001b[1;32m      5\u001b[0m     \u001b[0maxis\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'off'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      6\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mIndexError\u001b[0m: index 1 is out of bounds for axis 0 with size 1"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAOgAAACBCAYAAADKUvbqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAABTgUlEQVR4nO29Z7Bl2XXf91t7n3jjuy93np7pScAAQ2SAoBktiaQtU8llumzJdkmlUJY/OJVV/iCX/cGhXGWVZcuBdrFsSSWzHEiatimJVFkUk4lIhBlgQs9MT+d+8b4bT9zbH9a59zVA4A0AApye4f1XdXX3jefee9bZa6/1//+XeO9ZYYUVHk2Yt/sAVlhhhW+OVYCusMIjjFWArrDCI4xVgK6wwiOMVYCusMIjjFWArrDCI4xVgL4NEJGfFZE9EXnhm9wvIvI3ReS6iHxJRD740H0/LiIvN/f9tT+4o17h7cAqQN8e/E/Aj59x/08ATzZ//iLw3wKIiAX+VnP/e4B/UUTe8z090hXeVqwC9G2A9/7XgaMzHvJTwN/2it8B1kTkHPBR4Lr3/nXvfQH8XPPYFd6lWAXoo4kLwK2H/n+7ue2b3b7CuxTBWXf+8B/9T/3kYsRsRzAVRCOPzaFsw9qrBTarOXwupU6EOoay48FAMBF6bzrwkA0Ms/Oe9m2Y7wjZxZIrvwAnV0PGjzvSewZTgxcou9C+4wmnntmOoUrBFrB2vSKY1gBUHUvRMUgN0djhA5htWaYXhPYtv3yt2a5gc7AZ2KKhM3oQB+LB5p7um3OKtYi8b5meN+x8akYxiLj9o4bODUOQecTB6Ano3ICtL0y59WMdAGwO6YHHBVB2hcmH5qRfTuncdoiH6TlDHUE4AR9APvC4J2f0f7VF+37FJz767/CFl/8OH/jL/4Uv+sLkiQozM8THBuCXv8HP4QH5Jrf/HojIX0TTY9rt9oeeeeaZb+F0WOF7hc997nMH3vutb/d5ZwZo2Q2oUvAWwmNPlQr5AKoU+m8YbJYjDkwJUoMthCqFcArxcc1sJ8Dmnv51cBaCKYT7AaMr+prpAwMGzAyCzIMI3kDRFSaXHcm+BuLJ4wE2s4RTaD0oKXYtZUeoE0P/dYfNIT6EOhXq5nQNMsABAsHcU7UEF0DvVkXZMtSxcPRsi6IvuEiD+uD7Wjirr1W1ACOEY0/rrhDOPMVajKn189ocpueEdM/TuV3jTYrNIV8zFD29sIQTCCf6vYVjIXuzRdkWRpcDsq0E95pQpYLNIN6zBDMhGnrQlfHSQz/FReAuEH2T238PvPc/A/wMwIc//GH/2c9+9ts9N1b4LkJE3vxOnndmilu2BBfqRdvmUMdQ9D1lz1GlBhdaEF2NgpknGnqCKQRTTzirKFuCOOjcLjC1Bk10Isy3BS+Q7vnmtT3RRF/DBULZFdxapcFQQLbhme0K803BG6GOhWLgmT+RU8WCKT3x0FOlULWhTsDOPeLBhWBLDUBvhXBSYUtPHQqTy8L0oiPbcgCMr3iybU90oqteHYMLoPXAYQvPbCvQ4Mz0mPN1XWFb9wvWv1oSzPQY8nWPVBCNPabSx4ZjT+u+UCeQben7eqvHakpPuick+574xAP8EvDnmmrux4ET7/094DPAkyJyVUQi4Kebx67wLsWZK6g4aN+r8VY4fK8QjaB1R/BWOHgf1B9pYSrY+LIjyDyjywHpoZ7sR88kTC96pIZsMyaYQ7YJ2W7FxucseKgSiE70JC9ToWoJ+boHB1d/DopujTjY+dSE6ZUOk13LrT9m2f6UIx4KJy4iyD3e6qoLGohYCOYw2YJsq8bmpln1PPc/kuIDfVw4hXhokAqCmWftuqOOhNmO0LnpyQfCyVOe5NOaso6erNj4vMVbqFKh+yYEuSPbijh62tK95Vh/qaa4a4hGNVXLsPchYf1FvUjMzkFyCHd//m8zufsaVTbl1b/1H7Lxg3+M1r2KKhE6P/hJ0BT3J4HrwAz41wC895WI/FXgH+qn5Ge99y9+b06NFR4FnBmgk0sGZ/XfndseLx5vBKmhfc/jBYqeMN80eBHKrqa5Nodw6hEv2FzTwHxNiA+hdc/qfQ68gdEVg801DW3f0xtdCN4IoyuWsgdHz/SJTnQ1DKaCOE+QedI9Q9kCPNjSUyHEx57k2DHfNCQH+pj0wDHdsVQt6N+oma8b6kSwuWe+rYEdZOACvfiAfi5TQu+6oWh7vIHo2FLHQr4O+bqj86Zhtm3wVldcZ0VX3kjI1i1FR6hbNUXPapYw0u3BxX/2zxGfeEaPgymE5NATreuq7oYerxrAf/0b/Sbe+1/mG+9RV3gX4swUN9v05JuOYs0TjTym1Gd4q+lbeuiIxp6ipyets1AlutcL5h5TNvurkRZOorGne7NCapBa0z9vNfBcIEQjRzDT96ljoWpBvl4zu1pSp2BqsJnuU02tgZgPNCXGQTCDaOKJRjV1IgRz6NytCeaOOoWqA627GUFGU2Ty1JGnajfFnpZQttE9aQDiNPV2EUgF8ZFeHOrI41KHi6BYg2xLL1YugrJlKNtCti4UfZBKP0cdCeHEUye6v3UWqrau/sFc7wdNn1dYYYEzV9CqUxMfWKKhcHJNV0ac7rHq1BCOPcmxro6LVXW+KYgX0kNPel+rqnUkxENPHcLoSoAPtOrrDWx/rqRqGcqWYbZlyTaFOvZUqeHcb+cgMHwiYv0rM6pOyL3vD5hvGuKhp7VXMvszU8bDFr3PxZz7jSHja10efCTCVDDfFuabARiYXnBgPSavmFwQijXP+ksl03MxZReKPhRrnjr2uKTm4q8Kk/OW45+Y0v+VFulBjS090UlF74al7BgOn4Oq5fHWE57oHnm2KxQ9jw88di5sfEGY72hQpntQrDlc4si2Da5TgxPyga78dayp8worLHBmgO7+hmF6XveF3Td0H4XXvZsXTUOnu7L8/2xHA9EUkPcsLhIcUIeaTroQ6kTwVve3Noc6MVSpBqzU+tqmErIB2DygSoTj9ztMlSJNC0VqMJXHlJ75y2tI7Jle8Iye7FJ0DLbQSqspQSotRLW9pqLDZ7sAREPBFFoB9kYLVq17GiDZRkDZcuDAv9mi/0ZO0Q8YXwxIUqMZw0AvJPGxEMyFKtHKrZlA+gCKnlaoXeRxzef1Bvqv6CbZlGCqQIOypQUpL3rMK6ywwJkBmu6XzLZj6hjtwHlNLeMTyHsabC4Am3mkqZKK01RUvFZkETD+NDC90cosoidkHQlVogWicOwJJx6TC2VHKNt6nxnklJ1U0+ZKMJUHgWItIBoKZQ/q1DPf0Iw9mOkJj9Ogwev7uhjmm7rHNSX4QLCZBycEmcfUHlsIptZAdSGEE8GUDm+1hVSlQtHVVVI8RCNIjhyTCwZTaLsoHjrEm2XaCvp5XQjxiRbOAMKZo+hqy6eOtBpts5UFzQqnODNAH3wkXp7Mhx8vscOA+Miwdt2Rrwt14tn4smN6zuINdG55Zuc00LZ+d8b0fAdTQPt+zeH7DKYUohPYeDHn5ImI8WVIjiBfF4q+J5xA53aFeLj7AwHB1CDOI7cTOndrTOkpWwHRWNPqvY9CMNFUMhwJZQfioSc6dhy+X4iOhSDT4Mw2oew5kgNDdrkA65nciohPNPjmW4a8D/EJbH5hxoOPtSj62ja5/7EWwdwTH3ttKc0FN9bX7tytSQ4LqiRpiltQpUYLV5knOXGUbUuVwnxLW0QuhLKnRaVorO2V2a6QHELnfvUH9duv8A7AmQG6/bslwydC5tvCY/8bjC9qpdZUnmgILtQVMzrxlB1hdFUD0BQwO5/Suucp2zC6aqkT3xSHYLYTkg+EquOY7RqioSfZ1/6g1AFB5omGQnzitAVzaMjWzJIUMXxaCQzd14Tk2OGsrpj9GyUnV0P2PiysvQTeesq2UKcQD6F1X0iOa3pvBLhAWzP5JdO0k3yzYsPo8VRX+lKwRplNNtPsYXTVgGjRaPBKic0cVWKZ7QjJkWYA4jXFrRIou5ayq9lCMId8Xb9bk4N1ejGLJo680PZNtma/97/6Cu8YnFnFNYUjnGrjPj7MMJVWKvOuWZLOsoH+e5F2mgoQyPrNieybk3OiVVWpvaZzBUQnSoczJcRjR9Va7Me0xbHYtxVdvQC4sGnfVCzTVlNBkHuCDEzpcIHgWtqXDRqanwv0ItK5V2taOqmJRw5vmrTbalXVVPqeVcIpTXCulWhxupeuU0/Z9lQtTzQsqSNDthHgQpotgP5dNoSJYKb0yMUfqfQCFk5o9tSa3taxpvverIpEK5zizBX06D0x6b5j48WS/Q91mZ4Xyr6mtyYHbzyjZ2u6rwZEQ0/vDU/R0f2ji1Au6hj6r9eEU0fZNuR9Q5B72i9VBNOKe9+f6v4sEMq2J0x1L+YtlG0txlz8xG3u/uZFohF07pV0b3tOrkYcfn9J9WpEfOyJR56jZ3S/nN625ANI9x3pQU3RD+ncKUjuTbj5xzfwxiIO4iNP55buDWe7QjjR9lA49TDRvmSd6F4x7xnm20IwgfmFGj8owHvGlwPGjyk9UOpmrxsI2ZbDFMLur9xh8r5d8p7BG+15mkovHEfPWPJ1KNYMsyslrTdDopvuD+inX+GdgDMDdHbOU7YN4TSijjVNDCeG+Y4jPlASgriAbEu5rmuvOOpNDU6pNThNpats2dLqZ7EG1VAYXzZUaYjU2g7JS2HwFY80Pr3pvrYdghns/cpFdr9UUifC7R8NSfa1LXHlfxdcVDFft4weM3RuaaU4H3jESbMXNNgM5tshZXeN5EBX+joSZucFaVZqU2mrpU41bQ8y3X/Od4V4KLQf1HTueoZPhsQvGGyesPchvZh0b0DRhWjqCeaO46cCNr6oKfGtP3WRsgNVx1NtFWz9WoQ4z+S8xVTKcW7vOTq3AvIBHDy/SnFXOMWZAWpKPVlxSvr2BiQQomOjzf7KE06EumEPmcqTHC6a8bIkNiixXXmt2v5QGl2+7ui+YZYVXRexTIOLrlZTpfY4K0jlcEFAuVUSH0Z6XA8XPD0EmaOolR8czE7Ta3Hoyt5UZoPZItWGsKnaIs3q16S82bpWe+1cn+utweb6hqaCcO6XBTRTs2RcLQkHpbaBpOZ0I5FbxC9aTB4XCj4Q5uvaGlphha/HmQHauqdFj2jqCCY140sRrgubL1bMNix1IkQjT2vfYXOHC4XNzx5Rd2MOnm/jRZlFRR8QLfx0bjlVhvQtBJ7em7VyaTuGo+c8g68I8dRz+LzQuw4IlM9PGO218QGE7Tmt+yHiPXsfDDVlLKF13xNOHbYwSC30blUEM00Xj69FlF2hjrTV0r/usHMQJyQHGkSzc3pBsYUG9eiaIz40bLxQc/icpY5AnNEWTqRtoe3PTMl2YrI1S3uvpuhoyyQaebKB9kd3Pj1l/wNtXCwkL1hM5bGlp/96yckTEZNLUF7NCN5MiIYqt1thhQXegosLTYOT+Dig6AHi6d2EsqdtjXAkHD9lcVarkPvPb2gBaKQroqmgfVcrpMFc08u8Z1RfemgoW55sw5CvQ92qmW8FVInV/d1juorHn+/gAz1xk093mG9pBdmHzYfIPcFcyfplT1ew2YalvBJQp9C57ZjHTf8zh9FVJcivvexJjyrKtmV8Fbo3mtW2L1z4xzV15Blf1M9lM92jdu7Wqi8VuPnjbeIjSI40yPO+VrmjE5Z846P3tsi2dHWWSnCRVorLnuHCr1ekB4bp3YT2A1W3VPGqSLTCKd5SzSKV0vhcoEUhBMrUEEy1we8DrUwGzV4uW/dIrcyhqiUEpQq9sw3BFkI4czjrKUrtneZr2jMMR+BCSzDX9C8+1rR4wbGtEsEWnvZ915AboKz1ZC7bommoEa2kViz7jS7QC4KugA0LyuoeGYG8b6lD1WvGI6cCgK7uX+tQiRaqdYVw5qkSg3ht/7jI45oKbFWL8pSNcnMXKev0vFasAfKBfhYAF3rGF4JGKKAECuCbyK9X+MOKMwM0HAnxieo055tG1R4B5H1oP3CI9xw9HZDuaS8vnDqOnw7BK2k921QSezysydfB5oZoWGFKR7aeUMcw3/F0bkLrVs3sxBJNVXsZTVUls6DAVS0QL7QeFATTkqodMrkQUbYhGwj5pmPtJWCu7CYXNi2PAmYXlJ0UZHqxSIZK5Tm5EmhRqNI+aDyscKGQZ8LwCYMtIRp6TEtbROHUM7lgCKaqyAkmStzX9FkvDlJDfi0nPkqw4skulgw+H+CNMHxfRf+rGpQ+MBx9sCI6tHTehOP3O+zEkBytVtAVTnFmgMZD7eGZWk/s/g09sW//mJBtWIKsERxXmnIePhdSh7pSTc4bptcKZBrgolAdDyLY+3CsPUwL8ZG+h6kg7+qedryuJ31y6LUnGoI1pxzV8cWI9NAyvmgZfn/O4LdighmUhfYr60aJkhwrcb9qwfqLnnBWU4fCvR8QeteVjVS1oHvLUYfCwQc9Bz4gHBm6b3q6tzy28ASZo2oHOKvFonjoyftCtiXEx7pCgzo1HL4noOp4+p9S5pPNHdHYcvKkUKWe6FDbP9EIBi87psNAK8oJSKm2Mgsa4AorwFu2WYTWfZVcZQNN+0zpSe9bXKjaxvSBVlwxmsqZvEmNa4gehJhCmgBuKp452vf0gNWmfdkSXKyvt7ATwWsF1c714jC+IlQdNH0utKfIMCSc6EXEhYJ4pe3VkVZobQGMms+ypTYp0fCUvJBvOFr3RSurNbTuqrA72xD6b2ik5D1LOFKvJS96vIs0NTnS/qiLlNSxsH5BYLZlkFr9lnCqY23f9WQbjb3K4cvc+ge/CM6xde1jtAf/NFI1gnNARP5d4F966Hd6Ftjy3h+JyA1gDNRA5b3/8HflbFjhkcOZATq/WhBMoyUNzosWegav1Bw9bXGxp3+jYu8DIS6C7g3fmHJ5bXvMF1QiqM+rmDs59kvfIU1bVaeZryk5oXVPiE50rxuNlfva2is5ei5ShtDMUlSqy+y+YYlPqka7afFG950+ZFmUCifq1jC6ClW/ZvBFbRGVLajXK8pOpI8bG7Y/l1F2A+78kGHwijKe5ltC701dfcu2YEuPTHTfmBw2t1sl8NtcBdhVAuW2BltypC2oaORZu55zFMaUsefmp36eC//KX6ZT9Xn97/4N2q8+R7C7i29+Ee/9fw785wAi8seBf9N7/7BV54947w++R+fFCo8IzgzQYC9k8ljN5Alh8CWvLZNYCGfSONoJs62A7k2nvcGpY7atgWJzKNZ0PygOtr7gMJU6/blQGr8jZfAEU+i+qQwccZ6qLcy3tM3jAuHeJ2LiI7D3LMmh5+QpiI6F3U/N2ftg2vQ9PfMdIRzB+gueooPyfVuNA98emLtWSfMdj82FK/+HsPcBmF2AaAijK7GaheE5fE+glihHntmW0aqxVaZRtiFkG54qDQln2mZygb6XCz0bL1YUfasqnUR9mqq28MY/FzH4CvjXbpC2Nxhk65gS1h/7PkbXX2D9wi5l+xv+FP8i8L9813/9FR55nBmg0UgonVnu6SbndHUQ5xsNpBaCnNU0NjmUpaYxnDlmVtky0dCTN2mws6hONNQUND7UlROgtVdTdJVk3rrvl84H3ZsaFM6qFUnYpK2HzyXk683+1DdWKzNtuRQd/b/4hkQggNPXne00xZ22JT7WzxmNPEGuDnz+oT1v2Wkq1Y2B2eTSaUW46qjEro5hvqlsK1sIdapSsyqFfE2WFd3WfXX8m5oxQW+N6Xm9OIXdNWbjm9QR1PHXlnFFpIW60P/Vh272wK+IiAf++8bBb4V3Ic4O0BNt5rsAolEF5wx1uvAlUjpd0VVXAWX9SKPZ1OIKWEwB6XHN8ZOqIEF0BQumakLWuV8x3wgoW5AcVVSpUuHaDxxlW5oAzSjbKXVfKGKlAZYd4eRpD05fK5jryb5oaZgaZO7xOcug8wF07lbUUUDZEbI1IT3QqrEtXNNC0dZKOPV4aWxLHMv2j/ZmtS9ax/pedSzkmyp7C+aQrWkQ14lQdTyu0GNr33GcXDOUbY+3Qr5d44ZGW1iB4GKPD39Pn+WPA7/1dentJ733d0VkG/hVEXmpcatf4mFf3MuXL3+Hp8cKbzfODND+jZLjJ0OKXeXkulirpOPLelKZZnVLjmtM5Sk6luk5Q9XSk3zRsxxfDCg7kBxoo//uj0I4tqQHjpOr4dKfaHQlokp0xRpfsmSberKW7RQvShTwDbE+mHnWviJMLwjpHmx8JePgfQknT0K5Blf+74rRlYD5trDxQs3wmiXb8LTv6R7U1H7p/1PHQjZQ0ru3kOwbem+WlB1L2TNUrUY1k3mu/p8Zx0+njB/3JPtKfo/HNbwKdWzIe8LR+z2DF4TWrZpoaCj6qm45fF6Ij6Dl15geHbP+eUOdCvc4wZ/v4oLGK/hr8dN8XXrrvb/b/L0nIr+AjoT49a97zNf44n5HZ8cKbzvOdlR4c0y2tkbZNRw/pcSCcKIEgvhQ92POwmzbqiYz0ZUKr5XZhSTNlJ72XV2FEOi+YjH1wjlPi0jOwvyCnsCm9IyvgM1l6WRPk3ZGIy1EuUBXxmikbYqD9yUkRw5xhjyzFN1mT1tD0cjjbCEgml7XsVpEFH11sM+2POkDpfupNYsB70kOPEVfby/bwv1PtAim0LvumVzSfbZUAe07und2Iax9Vb+L0WVL2dFU3hZg94XWA0/cv8jN2QFjOSIp+0y+9Ltc+ON/lvhQvoaoICJ94IeAf/mh29qA8d6Pm3//UeA/+q6eFSs8MjhbDzrLiEeOYArzizUY1TFWLd3nJSc6eiFfE7INVaqoX9DCle+0gZ8cOWUXxULvZo2da3tGGg2mt0K+pa2NYA5122GKxrcobdoPTgn5QX7qCGjnGqzT80q3i4eO1gNl/Cz0qUVXAzWYNV5Kiz9Ggz8feOqdnGjoNeXNlRtcR+ooH471tauWMHm8wkXQflBRJ55sp2J+pVSlTleF4507tQrC1zzFulvavMRDTzyqsc6y+Sf+FG/+Hz/DK3/vP2Pr/PMM/Daj3/4thp/5rYd/gj8J/Ir3fvrQbTvAb4rIF4FPA/+P9/4ffDdPihUeHZy5gt7+qfNEQ09rzzF+X0V0YmntOaYfyZkPU7xRCVmkbugUaxp8LhROHje677O6krmg6Z22YP2lWv1jU3BTGrd4dUBo7Tl61ydEkzZV7Cg7wvgx2PiSvsfJNaMV3UJ1m0VPU9/kAO7/oCO9G7D2qrKaTGUwhTA7J3Tf9KQHNdnAKrk/c5jSMT0XEU6FLE9ITtRWpYqF0WMCyLKSWycqR1v/XYu3nuG1kPYdmIql6qlAPD3U/fn9T1i6b3jSA50xE05UkTPfEsp2gKlgp3yW8//cs2SbOuZifN9wvvoE/a+eLL9/7/3/hI4qfPi214Hnv4vnwAqPMM6WmxVaBTUlbPxmhKn05G1/OlWHgWauSJ3oKhlMYLqri3I40n1eHavXbOderZzYRNSFwWtrI9tspGolJHcDxpdgutPTAUeZrly913X/J97TvSH6nAzWXs2484M6EyUeetI7AfGxrshHzwZqJzL29F7XVbR4PKBzV/ejdWzpv+4oepqWRicwvmDV67bniY/1O3ARjK4B3mMKFZCb6jTddg8Ed2TZ+c1DHvzAOtMLkD7QNhEI0VizCXFa3EoPHVUszLYNwVzdIRbFpTIVxtd639MffIV3Fs5McW3R+Ah56N6pMLWKmJNj1Tr6hoeAp9E5akAuHPzUNkQb+OG4Ipy5ZWHGlM0eNjjVUkZjTTenlx1lZ7FPVC2qakZFWyFJM3NlXi6nlYlTQbnN9DX1MUqgt6XHxZrOhuOaOlYW0WxHub40086qFpRdXS0XliemhLLrlsocb/XxptRjCsee1gMPB8fUkVB29bupUzXEDmZ+SdoPZuqrawt9j4VFinj9jC6E6bmVYHuFU5zdZhlrWlmHsP98qILmCE6enNP9VEpr3zG8Zth8oQIP+88HbH65JhsYDn8op/1CQjDXkzwfBNjSs/Z6Rfv1E6ZX+4wvWQYvObINQ9FttJiNYDpodl1lR82g4yN97+kFj7eOqiU8sD3WX64bnqxea4q+UKwJl351ysm1FifXhPLxjPSFlI2vVrReesD5Yovh4zEnPzaj/TstohMdAZHu+8ZT1zPfUkvO9g2HVPa0SHXiCeeqZnnwYUPvdUgPK+7/ySfUfOy2YfjemvSOJd3zDF6acf/jbcouzViLiLIj5ANP/3VHPrCU6xW9f+wZXQ4ZX11ZnqxwijMDdLpjaO070sOacOYpWzouL5yk2LwhA2QwfEJfJshgck7TxPB2vCQHJCNlGBU9XaHC5zbUXMzDyVM6O7R91zO5KHTfgHCqS/NsV3uw6Z4nX9c0uv8qTC4vWEieqqV727Lj2f10RRUbsoFw8ydaytvNoPuPE7J1OHw2AM4x27RUqdD7tRRb6MqpLRjtfY4e01ZM1RZOehZbgBNdgbMNoXVPvxclugvZmmV6AVr3lZQRnhjSPU889gyfbC3pe8WaXnDCCWx9watPkUB6K2S+4XTq285q9sMKpzgzxV34vIK2NEwzz6RzR4nh3jZDh0KW7nxlTwtBNpMmPdY2StGYS3ur6WUda3rnB6WaXxfaPglnEE3c11Zb7al1STRWQoA4fV/kdL6LN0oVNCUUA0edKoGi92aBqdQXaDE7xUUqiVs46y1SaBfozJSFa2CdNgqZ5nuoY68DiyNZOt17q1pXm3mt+o6EaOoRpyuzs41mdq48YQSikTpJ6F5b9a7egJufec1c4Q8Z3lJu5qwSDbKfHOE+12f95ZpoWOG3QiqEcOJI9zVAsnVlybgYXKDzMAHm24ay72nfEja+kjG8FuvetPBMz0cq3O4rvzdfE+2rNoUbH8DwvY7+S5agqdpGIwgaJs/mFwuqlmVy3rL3IUM01FS1/5JOM1s44odTZe90bmVUccJ8R3jwcc/gRSFqhvTma43qptA9KV7J77PzGnTpgXoweQvT82p/Esw90dhx4dcKXGSoE4O/0/gqdQyzc474yNB64Nn44oS9D3dxIUzOBQRzfb+qA/5YSR+9Gxb+/Pf+h1/hnYGzLU8uQ/cN6NypmLzQIxk1V/r1cGmyNXxaSA7UQUEqiIcalKYUgplTc+iBsm68hcP3Jsx2NXjb9x2912XJxd344oSj93Yod/U29bH1TC6rzC2ce+zQM7psl8Wa2z8aNu/n2f6sY75hGF+WpmCkK9T4crRUz9z9p1KVxRWQ3rO09mtcIEwvqrF1lUB2uSDdDwlnIM4htVFt6FzT7vgYkkPH3kebFHfTUnwsWX5v3kDrnrpOrL2s1jB1JGQ7qfozof3kOtLiUDz05GsGU3niyWoPusIpzkxxizUlIgTzmvhQ08oq1TTSWQ0ik0tj/KxN/UX1M5hqgck1lc1o5JcVVkQb+mVqmlklah1S9qLlYKVFD9WLYApZuvgtZoOavCFEJJ468UvhtLcop1W0NRON1L5TnAZlnai2k8ZFsGgbypbgbVNtnkBwFCqPNmH5ORf93IUrob6eXow0LXZLMoULdAjxQma2SMfLtupDg7mS+hezbYLMN/tzoUzP/ElW+EOGsz2JBgVVmgJ6pVf5lrD2akn2ZEwdC9ufLxlfDvACrQcVoyTAW92v1rGya/x6QXIcYirdA3qjBl7T88LgFUe2ZsjX4M4PR0Qn2het2p7xhcY7dgJ4bUNUogN561jI1gzdN0xj8+k5elbpiCZvHBtOPOlewfCJmPi4meTdoTHFhuyxnDqJsHMhGgmtPUc0qtn8suPOD0dkGzoFLdvWfW98BN2bDmdVLte9cUrNM5VdOhYOr1m6t2vCiQ7v9Y1XUZkK6QO1h2ndyzh5vK1718Iwu1SRz4zSEldYocHZ4wf/r5ii7Tl8LlmqQQDmW6HOwew7qlaolpczp/svo6tC2RHiY0dyBPXNGG804CaXPdufrZePE+dp7dfEI6GKhfRI09qDToDNtcjjjXByDYKZYfCyo2zpSRyfOIZPC3Yu6shXN73Mts5JmVwwHL0nQSo1+PJG+b3xkZIGWg8ixpd1Fey95rn3SRV+966HSKUFq/Z9R52qC6B4JTyYQkcxZptqzh1OPdFEV9i8r1ae012Lt1YphstsAk6uQXxkSR8IO58rmW+qwCC9F0DT011hhQXOFmzPHWVqcUHjVDcXgtlCTQJezHJ1cI3jQNVhOQ8zmDXDlUbCdEf3YuJ1+leQ6Z6ybBnyNW3ftPaUt+pMk9JCY5qtUjZT6h64Tpqe5KQhADQmZVUz0yWYa9HJhSClVporr5VbbzSQbaEtlQX5oI6buTKuSd2bgKzDUxMvbzVAo7HHzPQzYhaTtdV1cDmXtPYgspwf46xmE6bS4Ur5RqQKoUCpkm6u+tE6+Ua/xAp/WHF2H3S7STGnnuKTE+qXOrTu64nbf13bBNNd2wyhNTqpeq1xoK+b6WDNPuzo+ZpwZEnvCZML6nDQfuCYX7SMnq6xg5zkH6bMN2RJHfRWkMLTvqepZ50aRpcsPmi8iiaNxUoFtnDkfdsMAfYcPu9p3TGsveZUCtdW2Vi+Jksx9nyrSVE9TC4JnRvSsKeUEeUC7Xu60CMi1F5pgOJU7mazZpT9mjC5VtK6oUbaLoJk6MGrwZg32nqa7zh61w229IwvBYyecLTuCLufmTPbjpicP5XYrbACvEWAlh1ZEgrM57uIhdk5mCG07qqiZHrBEx+pI0H/NTXRqlNdDcKJJzl2tO/mZOst0gNP742cez8Qkw8EHzTBOTEE+y3Gj6kJWXSizgZlR48h21Tf2sV+L9331CGMHjNUbddcIAKSo6YvGWjqiVHzrulFrRqHY6XnTa6AVIbeG7qvVu8gz/iy4EMdDBVkEB972vcrDt6vlWKbqdzOBTpPpnvLcXLVUAw8F37FMLoC0wtKRIhGnmhUEw89R89GBDN44n+b8+CjbUDov1ZSx/q6R88kpIe6f41OVrabK5ziLfugiz1RkKl+s0qblK6jK1E01PurVMh7hmji8TM1GdP0UcgH0XLcXtkNlI9rVL4V79mlrMw2gWFLT+eOevFWLShTjyllSY5YUA6LnvrdLlLN2a4S2heazrLdjJY3OtzIG3UpnFwUfKDOD+rRwtLGpLaeagDhbV35soGOq69jLTAF0yaFzXxjpq2FpHBSI7X2i8KJTvsuukL7fk04aSrPWUVy5BsChA51MpWO16iarCEcf+9+7BXeeTizZNi5U5DuVyRHNa29muTIE40bq8yONth7bzoNvDZMLwpB7kgPK9JD1X8WbWF4LVgS4yfnLFEzBr5KYP0rNa0HSiofvJwTTpXnOnhhRDjWdknVVWH4olBVNcOZ6tQTjrUCazMoL+Xk2zV1oqylYs2Rr6u+s+w5sk3P2itzorFauaCKMlyoUjlbgKmFaqNctpRGj5mmdaIMqKqlrKfWgSPvC8HU076jMjNT6yqbHnrmu57RNXXOT48cyVFNtRbTvZ2TDGvuz17h5b/zn/DVv/sfc/A7v0rZ1uOJG+meiPywiJyIyBeaP3998buIyI+LyMsicl1E/tr39hRZ4e3E2fNBn1XTnYXawgvYuSfyTaGk8pSpEgC8UUOu0SWd4TK/WHH+/4WqbZh/fIL/dEcnVedKSK8j8ClEYyUXZBtCfhRy/KSh7HnKtK9i8Vw4bjWmXnNdQfO+EB95dj6bc/OPJMhM6N50JIcx8y1hvlsTbGWUowiZWeSpCf52m+hEOLmWquVm4jl6VgkNurLD+lcd8w2heqoiGoXNaAotGoVTCKdmOYMl21AyvAs0s5hvh3RuO/qv10x3As7/eoGLDfe+P6B1x2ALQx2rPWkVOG79Dz/PtX/mL8H2Gq/93N+A+nnird0laaPBb3jv/9mHbxARC/wt4I8At4HPiMgvee+/8l08L1Z4RHBmgNbpYiaLpofBvFn50sZOxMtyzsnCtNlU6vgX71uqtOk5nsQkzaQzgOi+w9TqPTTbCqgTTVPDSU26L9hSmG/r+HipoH1b97LeCPMtQ9kDjJBtRJQ9hwsN2bphvuMpNku6OxN+4spX+a0Hj3P3wRrWOvLUNQ4PNCtdU4haLKSV7lerFNx9LaV6w1LHWUdKJJAKjIU60HYOsKQT1pFQdqzahm6Hy4FJPoAq0CpuOIH81k3CwSZycRObw8blDzD/0gt0P/J7AvQb4aPA9Ua4jYj8HPBTwCpA34U4O0BDXWl8CKZq+n+1Bm7R940ESx3rFmPi1bzZE52oObQXSO5pMWQR8OZ240MUNHaaYcO6OS7oOZhPLUfPqaNeOBY2v1wrt7WnXkNl1zfkfItvF1SRZ15a7NNjnhoMeXbtPv/B1qf5r4Ip/2v5AUbjFoSOsmtwgVlWcR+GqXROjNRC+476EXmjwdx+UJINAvUfch6pBCOa5i/aKjbXC1eVKkd4cl4tV8LJKZm/6nhae+APhoTdNYqep30HBscph/5ecxH6msP6RGNtchf4d7z3LwIXgFsPPeY28LHfz0mwwqOLsx0VKmjdgGToqJKmz5kIGy/WDB+3lF11UYjGWuho3yu5+0/FuEi5t8FMqXbRiUq3FquDprTK/jn3m2OGz3SYbxqCoymTKynzLUPvdU/7vgbm5II67uGg94anTtWFID52dF+MKLu6xywOU+KtA3568ClmvuSg7DAat+j8dot8TdlJ2XYznqHps3bvaKo+PWcYvOQpUxg/Dutfrcj7lsllmO6Ey4LTYnJbMG+cJNCeZzhhSSnc/HLBbDvAWaF9v+TwvVHj2CeYRa+0gPhYyDbg+IMR+Q3DyRPm4dksnweueO8nIvKTwC8CT7K0Yvsa/J7ezMp2892Bt5iwrX+7QKjjJsBExdHxiRaMpFaRdJUKtgx1JJ/VKmslmlJWbU/6QFdZnAqz0z1PPDRML7XV+LqAw49ukTdu9F7g5LFGZzpTbqs4Aa+ECR/A+IrORakjXZEJPeMi5p9Mn+HvZuvcz7psrE043lG6olQaGItVL5yqY58X7dfWob6uyYXhtVDdCO8oJVG8EiA6d7U4lK8J7Xta0KJJ38VpgM43AvI1NeDOB5E630/1QlbHIBtrFG8MyQfqcC83T+jkXZJGzQPgvR8tfgfv/S+LyH8jIpvoinnpoZ/pIrrCfg1WtpvvDpxZxV20LxYOfYv5m0VPq6bJkVMGTwvydR2wtEjRvKEZVOspNhzimgJRrkGdDD3tezWTc7Zh8XiGT7FkG9WxyrxmF7TVY+qFvvTUeiTfcGpHEjZMo7BmOEv5R3vP8MtffS83RwO22xOKTTUps5k0BtRqjRIPHUVX53baXIPHhdqmmZ73lF1lN5Ud7bsuLEtAiQdBplKzcOqWYxK9QDZQB/lizTO96JbkCZ2zKiQXL1EM95lnB1BUjL/yu5xLn1ZBQWO8LSK7IiLNvz/a/FaHwGeAJ0XkqohEqG/uL30vTo4V3n6cuYJWKU1f0ZOt63j2eOQYPWbJ1wRxWsmUGozT/VkwEWypU61PrhrSB8K533acXLH4UINzclFwofYfg6n2Sr3o8KVo6nABHL5P6L+qesvxZUPZddi5HkuxBsFEuPA5x8FzgVL4rO4ZT273yR6s0z2B/WcDuATxvqV1Ty05vYHJRUO1Cdmm5tzBTNUlRacZ+hsvpm17pruGy7+iQ5Um5wNu/GkhuSN0bnr2P2hI9pWqp9+HZgfRiaeqhGgoDK5rPj3btOx/rObcrxuisePae/8Er/9X/wP7zrN99SOM/u0BR3/v/6PzTzzwbwL8GeCviEgFzIGf9t57oBKRvwr8Q8ACP9vsTVd4F+ItuLg0ShFlxwSZXzraLZrstmikXA6SY7U2wUN8UtO6r/zXk6s6JrCO9LWCObAwD1v43trGk6hxNggn0oxkUPvKqiWYuikozbTnmA0s8bFvmEuCu5VgvA4wqlpCOLQcFOu0R9oKqmNRdtA9NT8bPuNo3dF0XYtEmsZHQyViBHMl5N/5wZR4qK743ZdDncLdF7qve2a7MNvVi0sda4pcdvSC40O1g0kO9OIRnli8ePKeJXjiOT5w6T2NSkd4sOfZeO6T9FPdHHvv/2vgv/5Gv4v3/peBX/7unQYrPKo427i6MfGqU5ZzL71dVC39co4n0MzAdMvHeBHixvB5kUKaUlPfYOYJR55oqCLoYOqVfN+QEWgGMC1mZYYTTzhVCp4LVWUilfYoo4k+Xxwke0IwbUjyHsITIb3XDBqOtNfpQohHNcmRBkIw0+NRZ3yPC/SYFk58poT5par5DI7WA4ctNLtoHSofuVqrSY7V7mTBdgL9HPmab6ae6fbAhXoxKXrNBa7RlkaHZmmsvcIKC5ytB/U6mMgUMH6qJjsyxENDcqB0ubrlSfaEsu21r9gEZ9713P+YJd3TVbB1X2VjC6ZNOPPEJ45wUjHfDEnfKDFFzc0/1qZzS/eI822HKQ1upFKvzi0d8ZcPhHRfBytNLuno+ipWCmJ3zy/HBJ7/J1Nm5xOm5+xykhpFM4dlzRJknkv/qMaFjjI1zHYN6y/oe8x2Va2DwOSCYeMzgos8R8+GWtiZa3ZxciWgfdez9qqhc3NO2QupU0OVGOabagi29qojWzdML3gufOQu+6MLxMeedF+dIZJDT+9GxuixdLlHXWGFBc4M0N6bFdOdgKytnNn4WH1gxasHLSdaRKnaKjtbiJPDkTB4RQ2iS5qJX+u6RzNVY0BWGbwJmFywmnoOHYOXHLNdw2wX2ndOB/lOLgvxkZzOZEmkCdpG1pZ50j0hHulsFpsJ48dSTOVJ9x2TiwY71wvO5BL0r6uUDg+TXZ28vZjXYkqvs0QrTZOzdU1b1Z9IpXPxoY69OHrGqvNCKkC6lM1JDRhd5etYDcaCKdz64jl6Q209hTPH8FpAtilUqeptw6mntVf9AfzsK7xTcLZxdeZ0QHZjQh1MPdHULyVZUntsrlxXU0EVC97oahmN6iXR3lQNI8mi/Uev7gLzDTX2qlJZ9hTrSPm2rQeaLpcdoWqphMsFLFNPZR6pntSF2ptcOP+ZCmbbupLZ8tTo2htZ6ld1GrdQdoWqrc8p21qtFt+4DFanWlcdvqTDeIPML42rXQhVWz+DC06lZYuZpM6ezoVp3VVBu/r+1viF9UunWTX9SrC9wtfizBV0dCWiauveqez4ZiaKIz4uOX4mJd+sAUs09EsfIYCy7Tl6JqJ70+l4+WlN1TLUkaEOITmqOX4mZPxETfdVSzzUVsThn54Tf77N1hcq0rtzXv9THartgs5LkY5ZCCDfADlCg7ctjB/XlSoaCUXPkg88xUaFnanDXjG2zB8vyGaWYGxYe0UvAtNd04jPNe3uv1Gz90GDizWt7944HbFYttWrN5yqeDQfGCYXDa37fjk6Ij2sWX9hDsC9H+guZ4fa0hPtqU1K3heO36uSsmhitbc6hvRQK97ZhlB2wu/dr73COw5v0WbRlC8+BkT3ebMNS3RSER+CyS3pgc5ccZEGSute4+pXeY6fEWxmSPd0Upj4ZirZSOjerEkODePL2jcM5lAex7Tnug87frZD73WQ10KSY11tqtQQTvQ4TAhSeXqvKkup6MHOl0rmG5bZbsDmlwvGF9SaZe1z0XJKdvdmzvhSRD4QsnW1LpFa55HGw2aVFW2xmEpXvs6dmrIl5ANN5aORp7XnGT5hVPt54pltW+YbHfB6fz7QNlIVSzM8SplH7Vtqx5LsFwS7Fh8I2UCN1cKx8p1XWGGBs42ro6ZoIepggNGBSXWkPNPF/k/HzGuaakvVPtriobTWaxpYx3ryL8Y0hFPXVE4FU3uiQwtO3f7qRLWhSxVNqUSHxarnwibNnevtC4sUWyjt0M50L+cttPdqbKYnfjCrMA2dbuHMEMwedoFoRh7GpxVnWzhcqBcB3xDjxfmm76nV7LIjmo4v9qCNU4Opm7Td6v49mOlz9P2lEZtrWr2Yc7PCCgucGaBVS1sk2YY0A5M0MPJBwHzHL+05FkWR+bZ69rhQCQbdG9C/7hi8ki3tNoNMKXonVy2jKzomIsg8rQcl/etKEpjt6Pvla6rHPH7G4gIhnDnViDYDl0y58B7SPu3B+0PGlyx1JNz7ZIvZOWmCTkn5k0seH5qG1QTxkXrvdm/VdG/Wuh9uabFmvq3c49ZezehSwOiqMLtSEp9oD3V02dK54/QC0swZDeZeecupLGe5pIeVTjerYfMzx7io4SJvRcy3lIFVdiA+UnaWFpxWWEFxZorbuq9c1bIDw6ch2ddVYHLekBzoCrP/fQGdW1qZLNtKk1tUW/OB4CIhyJWPOtv1jK85utctrfuOeFhz66JlctFQtmPSQ6eUPQfd1yc8+GgHv5uTvJgSTJ1OV2tZum/qclnH2g6Zb8HsoqPzppLNvYVs25HsqaP7fMPSvenUouVczPScZgC7vz1j+FSK78P6V2bMN1oqSYvVeaGOhYP3B82KCcFxwMEHHeHYEJ2ok/7yi5zrxalKhOTI4a2+R9G1jaODsP/xgcrzKpgPtFkaHat1iniYbRqy7e/hr73COw5nrqALj9lworKyBRk8nKhxtDR+Ra6ZnWIzJYTrymYamxPPbNMsSeosnl95XKym1N7qylV0dM8bzD2ziy2tgs6DZUrqQk0jbdk8P1BGDwLxgRZ1Fit1eGIIp40Yuy2Ec096UFKm6gpYpVD2w2VAerv401R5FwuZVzaVVJAcyNeYVSvxoLnfNQSMxshM7VQ0MLUSrLNQF4U0F+r3anPVkVaxNFuKVY67winODNBsQ9sX3Zs13Tc8plCGzeaXpk2AQed2M6S3I4QTT+duSTjzTC6rhUl66Bk/pv5GnVvQe9XSvldTpYbhEwHRSKdYi4PxY7oSJUPP7R8T4kND70X1qK0TQ96z5Bte96gNgX/0VI04uPwPRpqitgEP25+v1OA6UQI+Huy8WrZpyq7j3scDioZdNL2QkK83zwfydQ2m7i0dwhRksPHVkvUXhOhYqBMdDoxXUYAX7Ruvf3W6tO0EdQ50oe53izW/3LuL09deiATGjwnFGqs96ApfgzNT3M4ddTHI+4bp+dNe3e0f7ahvbK1B3H9Dl7jjJy3RNNDG/EzY+3CsIxjmML6idL2omUNStpoDmLHUQJoK5puCG8IT/3vO/vtTyo6m2rNt9ZBN74M3Ophoet7T/4qlTuH1P9Wlc8vrVHAHkwtWLwp3K5wNmG8aprstTp72dF8T4teEoitLWt/Rewztu82q1xfad7XPu1DwLEZPjK/oXjHdE+absPFiSbI34/V/vsdBGGKLkNk5zTBspoSOfKCf8+L/WzK6EmIK6N4qqFOdb2FKoXOvZHSp0b2usEKDs6l+tboiVG0twlRpM9l61lDn/Ol8E1N7wrESEBbC5TrSolB05JVJVGvlcr6l9MFg6plc1qFF8VBnlZi6MfHqhgQzj3hZjioUpzNc6pilMVg0gUIaZUxjJu1CrYx6gbIdkA9kme7ahvNbJZrK1rGqaYL5IrVtCBOxVqYX5H5T6JwZbxrJW+khEaX2dWO9v7GHCaY6gW3BZV4MMc4GevHCeVxkmK8HS26zVCoCD8erItEKpzh7D2p1NcnWPb2b9XK169+oiY91r2hzFTTPtg3duxVVSyh7jaN6CfGRprrtu9piqFowvaCyr9a+Q56Y4gJo3y9JDtXdz1th//mAaKqaUQ0wfX6+5rWdkYAPNUWMxp74QPeZ4lXXWSUwvQDHz8LsUkUd6963e0MrrpNLeizZesPjfdM1xw54XcmnF2HyeE001AptmRrVszbpq81hct5y+FxCfKw2o6aC3uuejRcL+m+UlG30c809h42NS9VWFtX4MWG+rXvoqm2XvOUVVljgzBXUG2n8aoX955WL276nk82yTaFs6Uqqj/FUiaG1V5OvGbInNA0UB0U/UIMtr2bQVTPVrHU3Y/fvJNiiwIXCwUdr0tsB6b5n/SXHfNNQ9A35huPEGFp7nsd/YcKDj3YRD+kXLdNdaVY0OHyvIRppQJRdXZnDqceUBtBeZrYhtO95wqkqV2Y76lGU7pfk/YhsS5hfrtj6lHJ0R21d+V0IRSB0bnnKNsuKdTDVERamVKJ+mQqjJ4QqDdW6s60XDJtB7w3NQgAOj17mzi/9Ih7H4P0fZ/NjP7bcmwKIyL8E/HvNTzEB/or3/ovNfTeAMVADlff+w9/9U2OFRwFnBmjRU1J6fKQeuAt6ndqEaGD6ACiV0DDfkCUlbzniL1RiQnqozgdVohO+qgSOn20RTv2yitp+w+hJnTSrSkt1nFKpCBwg20pUjpbTjC5sCOlNEJhGU+qs9jPTg5rxJbvk4dpcPZbCSc1sOyAa6Yo12w0bsjzE+wGmUo5xNFT2kLNKpsgHjTA7a9LsWK1YotIvq7VVS4n2wUy9mTQ1188VjTx14Lj12z/P43/iLxN2+lz/X/8GvSvvJV3ffbhI9AbwQ977YxH5CdS+5GFzsB/x3h98906FFR5FnJni5gMwtU4fa91XRk++jtp4OA0Q/5D5c7ahbYU6aQTLjV1K0RPadzLiE3VL6N2qKDvC4QdVaTLdsRRtYftzOcGUxk1ePX1coKykcKQFnOG1sAk0Tzh1pIdO9ZyhrpzhVCulvpm7GR8XzM6pY1+xdmpbEkxKqkSIR8p8Gl1RF3tTqDGZbwY4JYcs54P6ppdbRzoXxkWarpZdVayUC+F4oPeJh+7turEqFfJ1iCaO7M5Ngq0N/NUNXDdg7ckPMHn5BeyChA9473/be3/c/BS/g3oPrfCHDGeuoMm+rgrzTXPqhxt6eq9q66BKm+JOJNhKR/Xhld0jzaBavK6Ed364RdXxVC1Hp+HrUmuTHrQH+ua/6pD7jvjQ4C1c+PUZVTvgxk9ZohOro+3fyJnuhpSNY/2i77kw74JG0VIqY6eOUvrXtd1RdqDsCQfvi/E2JtvQCWQIePHUqcPODYiQr2nQYzz9V5shSS1h8JIORUJ0TMPCiFs9kgRn4er/VXDyWEzZFcYXVXsKUMee6Y4le+WENFhj57Ml3giZ7zOa3SQfiNqZ/l78eeDvP/R/D/yKiHjgv28MwlZ4F+LsPahdTPvSYlGQQTCXpcuCC079cnWOiSGYe3wgCJrqmkJXLVMKNhPC0DLb1mZ+et+Q99Wxz1uwtxKdczLVi8FsN8bUnvXf1clpZVuPqWqdrtILd4dw2qzeMQ+Zaas0Tsch6kUjHDfBEmk0hyP9XOHIM72kRtnRUA3DlqlsxJLcUCVNVbbyVG2lGdaVvrctdCRiPlAChAuUSO+DBZFD1N6l2R64SKhiQz0GP9X9fnL4tVVcEfkRNEB/4KGbP+m9vysi28CvishL3vtf/7rnrWw33wV4yyouzWpUJRpo8bFvrC7VIqRqoSuQ0ZR44d6uxtT6OjaDdN/Rvu9I95QXKzW0Hui/y54G/drL0L7rCcdKG5xcsGR9y9ZnhsuVuI4MZUfZQFKrO1809iRHtRIQ2jooCdGKqi1OifRSK9PJ5k1gF5oBtO85Bi/PiY8gGkF67JRBNdW0WTWrOgRYU1j9fFos0vfLdnTPmh45pttmyS6KpqorNVWzZ7cQtvuU4yFly+jYiPqEoN/D5iowX0BE3g/8j8BPee8Pl7+L93ebv/eAX0Dd5r/2t/P+Z7z3H/bef3hra+s7PkFWeHtxZoC29hx1KkwuCumB7suKvhaOWg+0ZdG+7ZfTp1tN8aeOtRdoGhf6+Y4QD2um5wyznxgTHzdV1BCyDaUGth840sNaGUqR2nVWLR23MHm815hGQzQum72pWqf4QJjtGA7eHxAP1UokHsLsnGNyUZics8QnNc5C2dcpYtPzQrYF6y8qQ2J6zvDgYy3KdlPUams12Dbm1NmmSs/WXvZ07jmqtjB82iyd9KMRdN4wJEcOU2qbaP2rJYNXS46eNWy8OGXrCzOqFqx/pWKQXGRSHzDLjyhNzfCV36X9zHON4Lsxwxa5DPw88Ge9968sfhMRaYtId/Fv4I8CL3wvTo4V3n6cmeK6h+aEVIksp1DvfGrM/U92ydY9/ddgvq2Ss/UXtUeJgGma82UXsk1HPrBEI8/shR5lR1cSZ3WQb9GF46cs7TuNOqaC6ERTzIXP7GLk4GwnXg5fmm8tZrqoA8P0nO6VbQHJnmmcHnS/uNifRhNPdaIuCidPWDq3HOZIPX3jRghetb72e0gfNMbWcWOBYnX/uZCRVakSMupQGt8kOOiEy/bP/ve1MZV+1mxgmW8Jg3/5p3j1f/4ZvPOsvf+juPdt8+Dv/zbpoQP+LYC/DmwA/01jj7top+wAv9DcFgB/z3v/D36f58EKjyjOns0SyTJIXXiqwbRv3Mf9cJdyTW0UqrbONAlnDqkMBk09NcCEuuMoOpZgpvKzyQWjjgPNRLTxY0KxoQLuBQkgnDQGYOZ06BJoMcnmNER5dT+wmY6dOHk8Xuo5k6Nm+nYNdWxUm5kLNnfNiD9h/Lhj8LInHpaUrZh45KgjYdY1jWZTi2HJkW/SZ9FxF1MVZZed0+9naaUizeyYjoNKnRkml1Vs0Lmp0rKqDa1rT7LLv9+4NMBBWtH/xPez9srSdvMvAH/h63+TZmjS87+P33yFdxDODNDpBSWF+wDiN7WoUXbg3r/wJN2bjv5rMNuC/ivaPyxbRsXbosJutcuE5G6gA3WtjhMbvKoppwtUmTItLVIKndsV03NBU3DxS++h5NgRj3T1joc12bolmCmdrrWn3kfZesDWlyrKlqFsa9sjGmsL5fgpS/uuJ75eM1+3tPZq4hOYnbeId1SJJdsSsi27DLTeTXUOHF82mEonk+XrnsFXVCRuKiVApAc1pnTc+ZGQ7c84kqOS4dMx4ZEhPhZ2Pj3j8L0p04uekx+bs/V/JrQfOEaThFkznbz3Rs3GZwLKDpxcPXPXscIfMrylHtQ3gVRHWpHFay+wSszSb9bmUBuYtYTxVUcwE9Zf9EtiQ75ZkxwpWcAFMN0xynkVLdoEMwDD+JKmiKbSwpKaTav/7UJnWrYCem8WlB3LydWAvQ8ZwrGyg0abuu81pR6n1OqGILVhtiPMNwPiY0+2bnWidg3hWAc0lZ0AO9dKrM100LB4LYppNVsnaRc9mG9rkSoaAVi1B7XajqqSCBd5xCnJ4ujZFPGQ7At53WK+qdXlYOrZ/owao1WJXlDCie6hV1hhgbNd/fJmLMJEV1FTNm7zsTKLii5LcywEtQSJNKiX9h0C2FMiOaZhCHVYDv61+eJ1G4vLhyRZ3pxqPr1RaxFl+WgKW25UlH0lTGRbXiu4NHYtoSwHM9UxyrNtnPSUNaStGhdpyi1e97+mVCldHclSxSJNFZbGta9Yc1SJsq3yNbN0BczXDKZohAFG73dWW0HpA0++rp8nHnni4wqbqeVL1dbHhytPohUewtmufk+IThWrF0GkTfmyrftABPKBp3NHHz+9AOf+iSpg7n/S07ptCebQf9GSbTbKlqEWlKrk1A2vbAXUEWx/PuPwParL1EKRxxTCyTM1/dcgnDim5wKOnk6wueovvQ0JR6pDPXo+QMZCOPYcPu8p+oZoaEn3HeFUV6mqpcbXpnK40HL8VIyLAN+s5F4ldFtfLJnuBhx+yDH4omn6r7DxFd0rF129tpVdvfj0X4FiTfeX6y/qd2BqHaG4//4AFyvtL7uWEd6MWfviIS//pQ2CqbD2Mkyv1MxKYTJdpbgrnOLMAN36QkW2ppOyBy8X1ImhahmShuljSth4oeTkqjbmW/fUeQ9g7UVDfKJFl6IvpPvaqpBKpWPZhu7purcNLlYvonwQqvfsTLQ1c6RV2PZNiwsds+2AyWVRu5FNOPiIMPiSrjjDaxHxAUQn+pxoZAgmelGZbxtNaR30btTKFOoZZuc88XFDYJjq7E6t9gqjKwFVKoTHZikGsJlO4V4YXe9+KifbCMkGKh6Ih0rqGF8W4mMwmbZ1ooaIEE486UsJ4RSmTwxwiaMuLVUCyX11wPd2tYKucIozL9fBpG70lYLNa1wkFG0VTi8qpC40VG3dk5lKeaouaPSicuqV641Qh42WM9LbpIairSf7giQu9amAezFqMJyogsYFGuAL3yE6lR5HYz2yaKt4o6LqpY4zZJk2R6N6qQOVZhRFMKNxntfj0DZOQ66fN457i4pwstCj6gHbwi1XXm8a5/u23mcq3b/qflidCdMD3X/mawYzN5ha+bzBXI/F5is96AqnOHMFPXxvzPycp1yr8KbVONB5vPX0rhtcCLM/PWT+4gbhSLWNyb7S4qaXHVIK4VhI9zzD96ghWDBrWET3PGvXHUdP2yZwINnLma+nuFAJ66ZRiKgfrTrmDa4XzDcjqiMh3YuoY12Z+m9UHD2jsjYNLG31BJkabksNQe4JRyWT8yEuFHY+U5Pu5bjQMHwipmovSBhqpF10dWWsQ9375uuQ7utFo07g7g/ExMeaig+uVxy8LyDbrokPrBLfK+31Zpv6/umBKmycbSaVf1GN1WbnPPGhUhxtvlpBVzjF2Y4KDtq3BW4HVCl03lR63PGzWiixBUx/Z5ONG04nYl9TSZYdg9wwzV5PT7jua3bpRBDMdJ+3/7y6A0aNomT4VKpjBgtVfZSpnsjO6upTx8LB+1Kl3TVFp53PlbhQmG1Z2vcayZeFoPGsdYFw/B5PfKzV3my9jVQa/OOLlv0PpEgl9F93TC/qkKR83dB/zeEiYXaxIt23xCc65j6ca89XamHjRSW7Vy3D4bOBKmFetfRu6sVi/JgajYlfuFMI04v6+dp3PZNLerHqvKmGYVWrKYitsEKDt7Q8sYXKvPJBwxCqPMFcRc4mh3SmKeVij7cwfw69pnd1I9OKh0rtqyPB1Op15AMlAQS5KkSKnlnapeicl+ZkNdrq8YtZKA0v1lutxi5G+mlbRWVmVXL6fFtosQsaAsVilOLCANvoZwjmUIlQdpQL7AXCkaVKhLDyxONGANC0gUBfo25IHAsbGKlOzbURTaG9UbXLwt9IvYNZpsYmb15vFZ8rPIQzA9SWekLWsVCsqddrMNNJXem+J8gcZdtw/B6VoSV7je1HrcF5+LyunuHIsPnlmvmGFmbqVPde7TvC2qsTsq2E2abVVkhXWzriDMmhBlxhVMJmC09y4NUcO4I6dRw9EyzHxo+fcLTftKy/VDM9HxBO1K9359OevKeb4a1fu8P02R3yNUsydMTDRsvZgsErFVnfcvQ+z/SCru7nfqtm70MBiNDar5mcb6Zyzz2H7wmXF6XkUI278zWP+ABTQnKkAm1TqrXK6H0FW78RYmqd3L37mTnjizGHzwsbX9ShUw/TK1dY4ew96Adqei8H9G9UmDIgHjmCzJMNjGorO9r/iw8EpGlvvF9TwPZtoXW3GW3Q9pSpMnsGL2m/NB8I5ZZndj5ldNmSbXi6N9X5D5Sovxg1WCXQv6vvXaXChd/Q6unhe23j4KAMp/O/BtNduPdJS+9VNd2e72hKaTP9c/CDFzQrEBhdsmqdiV6ERpcDTOkZvKjzHcRBtmbJtiuqxAJBMyWt+fLmLG07O7c9g5d1da4SSA9UJzu+fFoA2vqNkGjimG8Y5ruO2z+UEsyhcxMe/FBFcBTQur9aQlc4xVukuFqdzXuWYK5pY9mW5YwVDBTRaYNfPMtBvnWiKZ84VLLWzCtZVDVVfiWULVmmgi5Q9o/e58j7thGF62vXsVB0hLCn7nim0CptHTXeP5EsObfe6n508V5VS4+p/cCpp26iKXaVqCXJohKM13RW54VqBdvkBpurqXbVkoc4vvo9mbIpHEWyTFFdcHpcUmvfOJo6io56+toMinWHHBrSfY89CQjmsvTTXWEFeIs2S3rXUrU9x8/q/8ePwdEHauZbDdMmhumVirLjl8ydtVege6MROCcLcbIQTtWTaLajARyOfONNpIEWjoV8TZZCah0VoZXTwasFVSzMtgzzHWH/+wJOnjDLWaAuhGK9ZviUjlsYvOSWxZZoqNO5y45nfqGm/+k7WqHtKX2wajfspGYUvbdavS36Os0724D0gaFz29O5py0aaVo/ZbdxLjxWD93JRWH4tH4v803DfEuo4wVdUmVs0/NKeOjcQvugkVq37P6ODjCOj903+ilW+EOKs7m4D/yyiDHbFdp3PGsvW7X3WFM6Xf+lQPWfc+i/XnL4XNjIyLyqQELI1wx7H9C3MiUMXikYXY6YXhRa9zzxUFPO+abQ2tdU9vhaQNHzuMRx60dD4mHjV7voP86UGNC5V5IMDcmhpbVfaTCEQueOugJWLejczZlvpcy3DQc/ckknbx8oGT451NfL1wwmVzrg+MmazmtWrUWzZkxESzh+KiB94Mm2hNE2bP2uI1/TrCLIaMyuWepas4FlflH9gn0A+z9Q0XklJD7ypMeO6uWA9l1H/wv73PqTu8qIul9/D3/uFd5peMs2izig1hQ1OVb+6GwnJJxoUJlCWwPewGw7UGe9ZoK1qdXq0oUNAaFZ8abnQoqeBlz7vhIgyrQRbXqWsznjY6FKDVVXpWdSKXHA5jR7PWF0SQfe6kVDJ1g7y3LkQ5146tDo6MKp4KynbiZmm4aQsJieDcr4ifdtY9epGUDRbarBmToCll2rNMc1rSCLP52SLbUnc0q0r2IhvWOIxrp3luw0YTGl8ofLllCe76u6JxSqZEX1W+EUbynYNgBevXzi44roOGf0WEhy7LG5BmAdG8quMHxKXQdAV9z6SIcClS0lIkgzaOn4Gd332bnQfn1Ifq5LmYaNO6DgSq1+Joe6eo2uyXIqWDTSlLJONF2eXlK+bjhWooSpdJUuelC1XWO+HTTDjbQFM9tWkkX7nmO2o2lx675+ljj3tPZPfYvKtpLwo6HQvqdysmzNYDZgfEWphcFM308NzFTxUjep++DlumlPGdK7VjOSAI4OXubV3/wlcI7Npz/GoP6ntRXTWzoqCPBfAj8JzIB/1Xv/+ea+H2/us8D/6L3/T78H58YKjwDODNB4rA198Z6T9YDxpRA5HzK6CuOrgp0H9N7wRJOFs50wfMZTdRxmvaB9N1H5VqFFGql1Be2+qY6AZRve+DOb1Kmm0jufrTl81pJvCK17spwYdukflUx3Q/I1Yb6lEplg7um94ZohvjXBpObo2Zj0wNG5lTHfiRk+acm2dGxh2VJX9+RI1LXAa1ZQ9Gn8kUr2PxgvfXjr1J9Wa2cqY8v7hvlWgtSazno5Vd+076o5mamge7vg5GpE2Wkc92/X2Lymji3JsSNvwatf/UWe/tG/REt6vPCr/yXp+5+jne5qb1TxE8CTzZ+PAf8t8DERscDfAv4IcBv4jIj8kvf+K9/902OFtxtnC7Z3TBOguueUWpYUtmgkhCNPclwzOWepE+2P1qFgM4vspeQDlvabeR/EaJq7GGgbjRvJmdWqa9FpfH6yplLajDfM1oMlESIaa0GnjoX5ptp21pGh2lCjLt/0EfOeIZipc+BsS5aEgaqlLRBbqAVmOAUczLdDTA5RqWQGTV8B0wTjwtmhIRp4EbJt9WRKDx3jS0pZ9IUnWw8o+kresBlMdy2mQD18pzWT/VvE3U3CwQaVh/6TH2D86gu0PrC7NFoDfgr42957D/yOiKyJyDngMeB646yAiPxc89hVgL4LceaGZ3bOM72gMrJiUFMnervNRIsbNyrig5xiTX14wqkadrVve7Y/W5Kva2DbrLGebGZvls2IwGjkiQ994y+rfj9aFZVG4qYXh8kFQ9VSNlB8qEOGqhRmj5eULSHvCdNd09hiCi5WhwSbeTq3nVaOvYrDy7YQzB3hpKLo6bS1aOSZnlPL0HTf0bmrs2eisdqbdG+WxEOtxgaZVmxdCNl5HbXYfWVIPtCeaNkRprtWectd3cdOzy94ygXhpKI+OKYlPepUK8Wyu0Y5PWlIIcuv/wJw66Gf43Zz2ze7fYV3Ic5cQa//tX/r0e+a//m3+f1/j2vQW0NE/nngjz34m5/+C83//yzw0eNP/+a/8fDDvsFT/Rm3f/17LH1xgVxE3unOf5vAO3nUxdPfyZPODNAVvme4DVx66P8Xgbvf4mOib+G5NG7zPwMgIp99pw9Yeqd/BhH57HfyvFVN/+3BZ4AnReSqiETATwO/9HWP+SXgz4ni48CJ9/7et/jcFd4lWK2gbwO895WI/FXgH6Ktkp/13r8oIn+5uf+/A34ZbbFcR9ss/9pZz30bPsYKfwAQLRKu8G6GiPzFd/qApXf6Z/hOj38VoCus8AhjtQddYYVHGKsAfRdBRH5cRF4Wkesi8te+wf0iIn+zuf9LIvLBt+M4vxm+heP/YRE5EZEvNH/++ttxnN8MIvKzIrL3zVpa39H3771f/XkX/EELRq8Bj6OtmC8C7/m6x/wkOghYgI8Dn3q7j/vbPP4fBv7vt/tYz/gMPwh8EHjhm9z/bX//qxX03YOP0lAAvfcFsKAAPowlfdB7/zvAgj74KOBbOf5HGl6HKB+d8ZBv+/tfBei7B98KBfBRpgl+q8f2CRH5ooj8fRF57x/MoX3X8G1//6s+6LsH3woF8FuiCb5N+FaO7fPAFe/9RER+EvhFVO3zTsG3/f2vVtB3D34/9MFHAW95bN77kfd+0vz7l4FQRDb/4A7x941v+/tfBei7B78f+uCjgLc8fhHZbYTsiMhH0fP38A/8SL9zfNvf/yrFfZfA/z7og48CvsXj/zPAXxGRCpgDP+2b8uijABH5X9BK86aI3Ab+AyCE7/z7XzGJVljhEcYqxV1hhUcYqwBdYYVHGKsAXWGFRxirAF1hhUcYqwBdYYVHGKsAXWGFRxirAF1hhUcYqwBdYYVHGP8/jfml5TtqaWgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 720x720 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "\n",
    "figure(figsize=(10,10))\n",
    "for i in range(1):\n",
    "    subplot(141)\n",
    "    imshow(batch2['inputs'][i])\n",
    "    subplot(142)\n",
    "    imshow(batch['inputs2'][i])\n",
    "    subplot(143)\n",
    "    imshow(batch['inputs'][i])\n",
    "    subplot(144)\n",
    "    imshow(batch['psf'][i])\n",
    "    axis('off')\n",
    "    \n",
    "figure(figsize=(10,10))\n",
    "imshow(batch['inputs'][0])\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Test de raccordement avec le UNet"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "WARNING:tensorflow:From /opt/anaconda3/envs/G2G2/lib/python3.7/site-packages/tensorflow_core/python/ops/resource_variable_ops.py:1630: calling BaseResourceVariable.__init__ (from tensorflow.python.ops.resource_variable_ops) with constraint is deprecated and will be removed in a future version.\n",
      "Instructions for updating:\n",
      "If using Keras pass *_constraint arguments to layers.\n",
      "chat\n",
      "chien\n",
      "ours\n",
      "Train on 1 steps\n"
     ]
    }
   ],
   "source": [
    "def DenseBlock(n_layers, n_kernels, input_layer, activation_function='swish',\n",
    "                                            axis_concat=3, concat_input=True,dilation_rate=1):\n",
    "    #concat_input: implies that we have a skip concatenation between input and output of block\n",
    "    connect_input = input_layer\n",
    "    for n in range(n_layers):\n",
    "        bn = tf.keras.layers.BatchNormalization()(connect_input)\n",
    "        act = tf.keras.layers.Activation(activation_function)(bn)\n",
    "        conv = tf.keras.layers.SeparableConv2D(n_kernels, 3, padding='same', use_bias=False,\n",
    "                                        kernel_initializer='he_uniform')(act)\n",
    "        if n == 0:\n",
    "            concat = conv\n",
    "        else:\n",
    "            concat = tf.keras.layers.Concatenate(axis = axis_concat)([concat, conv])\n",
    "        connect_input = concat\n",
    "    if concat_input:\n",
    "        return tf.keras.layers.Concatenate(axis = axis_concat)([input_layer, concat])\n",
    "    return concat\n",
    "\n",
    "# Define the modified U-net\n",
    "img_rows = img_cols = 64\n",
    "nb_scales = 4\n",
    "growth_rate = 12\n",
    "nb_layers_per_block = [4,5,6,7]\n",
    "activation_function= 'relu'\n",
    "gamma=1\n",
    "\n",
    "\n",
    "inputs = tf.keras.Input(shape=(64,64,1), name=\"inputs\")\n",
    "targets = tf.keras.Input(shape=(64,64,1), name=\"inputs2\")\n",
    "\n",
    "#INPUT CONV\n",
    "net = tf.keras.layers.Conv2D(32, 3, padding='same', use_bias=False, kernel_initializer='he_uniform')(inputs)\n",
    "\n",
    "#CONTRACTING PATH\n",
    "skip_connect = []\n",
    "\n",
    "for scale in range(nb_scales-1):\n",
    "    block = DenseBlock(n_layers=nb_layers_per_block[scale],\n",
    "                   n_kernels=growth_rate,\n",
    "                   input_layer=net,\n",
    "                   activation_function=activation_function,\n",
    "                   concat_input=True)\n",
    "    skip_connect.append(block)\n",
    "    batch_norm = tf.keras.layers.BatchNormalization()(block)\n",
    "    activation = tf.keras.layers.Activation(activation_function)(batch_norm)\n",
    "    conv_transi = tf.keras.layers.Conv2D(32+np.sum(nb_layers_per_block[:scale+1])*growth_rate, 1, \n",
    "                         padding='same', use_bias=False, kernel_initializer='he_uniform')(activation)\n",
    "    net = tf.keras.layers.AveragePooling2D(pool_size=(2, 2), strides=(2, 2), padding='same')(conv_transi)\n",
    "\n",
    "#BOTTLENECK\n",
    "block = DenseBlock(n_layers=nb_layers_per_block[-1],\n",
    "                   n_kernels=growth_rate,\n",
    "                   input_layer=net,\n",
    "                   activation_function=activation_function,\n",
    "                   concat_input=False)\n",
    "\n",
    "\n",
    "#EXPANSIVE PATH\n",
    "for scale in range(nb_scales-2, -1, -1):   \n",
    "    up = tf.keras.layers.Conv2D(nb_layers_per_block[scale+1]*growth_rate, 2, activation = activation_function,\n",
    "             padding = 'same', kernel_initializer = 'he_normal')(tf.keras.layers.UpSampling2D(size = (2,2))(block))\n",
    "    net = tf.keras.layers.Concatenate(axis = 3)([up, skip_connect[scale]])\n",
    "    block = DenseBlock(n_layers=nb_layers_per_block[scale],\n",
    "                   n_kernels=growth_rate,\n",
    "                   input_layer=net,\n",
    "                   activation_function=activation_function,\n",
    "                   concat_input=False)\n",
    "\n",
    "#FUSION AND SKIP CONNECT\n",
    "print(\"chat\")\n",
    "batch_norm = tf.keras.layers.BatchNormalization()(block)\n",
    "activation = tf.keras.layers.Activation(activation_function)(batch_norm)\n",
    "outputs = tf.keras.layers.Conv2D(1, 1, activation = 'linear')(activation)\n",
    "print(\"chien\")\n",
    "\n",
    "model = tf.keras.Model(inputs=[inputs, targets], outputs=outputs)\n",
    "model.compile(optimizer = tf.keras.optimizers.Adam(lr=1e-3), loss = 'mse', metrics=['accuracy'])\n",
    "print(\"ours\")\n",
    "history = model.fit(dset, steps_per_epoch=1,epochs=1, max_queue_size=1, shuffle= False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Object `tf.keras.fit` not found.\n"
     ]
    }
   ],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}