{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "jLCfvrxgrNPX"
   },
   "source": [
    "# Exercise 18.2\n",
    "## Generation of air-shower footprints using WGAN\n",
    "In this tutorial we will learn how to implement Wasserstein GANs (WGAN) using tensorflow.keras.\n",
    "\n",
    "Recall the description of a cosmic-ray observatory in Example 11.2 and Fig. 11.2b.\n",
    "In response to a cosmic-ray-induced air shower, a small part of the detector stations is traversed by shower particles leading to characteristic arrival patterns (dubbed footprints, see Fig. 18.13).\n",
    "The number of triggered stations increases with the cosmic-ray energy.\n",
    "The signal response is largest close to the center of the shower.\n",
    "\n",
    "#### Tasks\n",
    " 1. Build a generator and a critic network which allows for the generation of $9 \\times 9$ air-shower footprints.\n",
    " 2. Set up the training by implementing the Wasserstein loss using the Kantorovich-Rubinstein duality.  \n",
    " 3. Implement the main loop and train the framework for 15 epochs. Check the plots of the critic loss and generated air shower footprints.\n",
    " 4. Name four general challenges of training adversarial frameworks.\n",
    " 5. Explain why approximating the Wasserstein distance in the discriminator/critic helps to reduce mode collapsing.\n",
    "\n",
    "This approach is a simplified version of: https://link.springer.com/article/10.1007/s41781-018-0008-x\n",
    "\n",
    "Training WGAN can be computationally demanding, thus, we recommend to use a GPU for this task."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "_Oc-7BCDrNPZ"
   },
   "source": [
    "### Software\n",
    "First we have to import our software. Used versions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "kCa9SQGrrNPZ",
    "outputId": "8d2f9bb6-453e-4310-86c8-f0475288f86f"
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "tensorflow version 2.4.0\n",
      "keras version 2.4.0\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from tensorflow import keras\n",
    "import matplotlib.pyplot as plt\n",
    "import tensorflow as tf\n",
    "tf.compat.v1.disable_eager_execution()  # gp loss won't work with eager\n",
    "layers = keras.layers\n",
    "\n",
    "print(\"tensorflow version\", tf.__version__)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "iin22Ul8rNPc"
   },
   "source": [
    "## Data\n",
    "To train our generative model we need some data. In this case we want to generate cosmic-ray induced air showers.  The showers were simulated using https://doi.org/10.1016/j.astropartphys.2017.10.006"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 36
    },
    "id": "PZwKf8QyrNPc",
    "outputId": "2f2401b9-c2e2-4cc2-e675-b38330c37482"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'airshowers.npz'"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import gdown\n",
    "\n",
    "url = \"https://drive.google.com/u/0/uc?export=download&confirm=HgGH&id=1JfuYem6sXQSE3SYecHnNk5drtC7YKoYz\"\n",
    "output = 'airshowers.npz'\n",
    "gdown.download(url, output, quiet=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "id": "657Ia0jArNPc"
   },
   "outputs": [],
   "source": [
    "file = np.load(output)\n",
    "shower_maps = file['shower_maps']\n",
    "nsamples = len(shower_maps)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "NSDISiI4rNPe"
   },
   "source": [
    "Now, we can have a look at some random footprints of the data set."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 297
    },
    "id": "S7quv3oprNPe",
    "outputId": "10149c08-48ea-443a-e8d2-cac6460f130f"
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ4AAAEYCAYAAABslZDKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABhqUlEQVR4nO2dd3gc1dX/P2eLiiW59wIGbMA2YAO26SCbZjokQCCBQIA4BPgBIYSXkITkTSHhJYWQhCQm1EDooReDAdFtsI0bCBsX3IvcZK3K1vP7Y1dC0vYyqxnrfp5nHmln5sy5e2e+e+7cufeMqCoGg8FgMBQLV1cXwGAwGAzdCxN4DAaDwVBUTOAxGAwGQ1ExgcdgMBgMRcUEHoPBYDAUFU9XOO3fv7+OHDmyKL4aGxupqKgoiq9MsVuZilmeefPmbVXVAcm2nzylQrdtD8fbLfLPVNVplhbOgXRnLXX38uSiJbvoqEsCz8iRI5k7d25RfNXU1FBdXV0UX5litzIVszwisjrV9q3bw8yZOTxuvXfIiv6WFcrBdGctdffy5KIlu+ioYIFHRNzAXGC9qp5eqOMauheKEtT4O57uhNGSoRDYWUuFfMZzHVCbq7HP52PTpk34fL6i2HWlz1AoVDSfTqsfRfFrKG7pZuSsJaed73x85qKjfH06pX4gsZbsQkHueERkOHAa8BvghmztfT4ftbW1iAiqypgxY6isrLTMrqt9BgIBamtrLffptPoBUCBIJOP9dzfy0ZLTzne+PrPVUSF8OqF+WrGzlgrV1XYncBNQlWwHEZkOTAcYNGgQNTU1bdtCoRCBQAC32004HKa+vh6PJ33RMrHz+XwdfBXDZzrbUChEXV2d5T4ztUtUR11RPxATS/dO43QnOWqpENdmtlqy2mcqu2x1VAifTqifVuyspbwDj4icDmxR1XkiUp1sP1WdAcwAmDhxorZ/CGdlayLZA7+ubMHU1dUxYMAA27S2EtVRl93xqBKwqVisJl8tWd26LuR1ko9trjoqhE8n1E8rdtZSIe54jgLOFJFTgTKgp4g8rKoXZXqAyspKxowZg8/no7KyMuPKzdWuq33W19dndRF1l/oBiCD4VbKy2Y3IS0tOO9/5+sxWR4Xw6YT6acXOWso78Kjqj4EfA8RaaTdmE3RayaVi87HrSp8ejydr++5SPwoEuum85kJoyWnnOx+fuegoX59OqR+wt5a6ZB6PwZCMaL+0PcViMDgJO2upoIFHVWuAmkIe09C9iCC0qGkPGS0Z8sXOWrJnqQzdFkUI2lQsBoOTsLOW7FkqQ7dFVQiou6uLYTA4HjtryQQeg61QhBb1dnUxDAbHY2ctmcBjsBXRB6L2bKUZCseqpZt4/B9vUrtwLX36V3L2JUdz3KkHIWLP4b9OxM5aMoHHYCvs3C9tKAwLZ6/g51c+QNAfIhJRtqzfwZ9/8hSLZq/g2l99rauLt9tgZy3ZZqydkxLwmSSh1vmMxLoHOi+pEJEyEflIRBaKyKci8r9ZO95NsPv5VlX++OMn8TcHiUS+mlXf0hzkzec/YdXSTRn7NElCU5NIS3bBFuHQSQn4TJJQ63xC9IFoDt0DfmCqqvpExAu8JyKvqOrsbA/kZJxwvjeu2Ub99saE24KBEO+/tpi99huckU+TJDQ1uWhJRMqAd4BSovHhKVX9eVYHyQBb3PH4fD5EhKqqKkQk4+ieq11X+3S73UXx6bT6gdhsa/XELSltorQ68sYWeyapshAnnO9IRCHJYxxVRTPILZarjnIpb1fZ5WsLibWUAa2NuPHABGCaiByeleMMsMUdT2VlJapKQ0MDqppVHqRc7LraZzgcLopPp9UPtPZLJ2yl9ReR9q/anBFLlgm0vTxtHjAK+JuqzsnK8W6AE873sJH9qexZjr85GLetpNTLEcePy9hntjrKpbxdZZevLaTUUnKbaOS3vBFnm8DjlAR8hfBpkoQmJ4LgjyTsi96qqhOT2alqGJggIr2BZ0TkAFVdkpVzh+OE8y0iXPurr3HbtY/gb/kq+JSWeTnihLGMGjcsY58mSWhqkmgpZQMOitOIs0XgAWcl4MvXp0kSmpwcn/G0s9edIvIWMA3oVoEHnHG+Jx+3P7976Ls8/JdZfLFkHb36VnL2t49i2vmTsvJpkoSmJomWUjbgonbWN+JsE3gMBsht7oGIDACCsaBTDpwI3G5B8QwFYv/xe/Drf13W1cXYrcl3Ho+VjTgTeAy2QhFaEne1pWII8GCsi8AFPKGqLxa8cAaDg8hFS8VqxJnAY7AVihDK/oHoIuBga0pkMDiTXLREkRpxJvAUka27GqmrbyQUjmRlp6r4mv0gQlV5qUWlsweqEIzYYpS/wULCkQhzlq1l7badVJaVcMyYvejZoyytXX1zC08vWMKby1YytdLLB6/VcNHECezRt7f1hXYYuWipWI04E3iKwLL1dfz2qbdYsnoTXo+biw8awIN/fIybz61m3B7JJ8upKi98+Bn3vDSbzTujIxyH9e/FlWccwckT9ytW8YtKdCSOuSx3Z2qWrODnj7+OPxgiHFHcIoQjEc45/ABuOrsajzvxj+XM2i+46dlXAWgJhThi1DD+8/FCHpu3iAsPHc/NJx5rcr21w85ayrtUIjICeAgYRPR51gxV/XO+x7UrgZYgH89aTKAlyKQTD6SyV4+U+9eu3cJlf36CpkB06GggFEYVFn25kcvuepJ/XvV1Juw9NKHt75+o4Zn3l9ASCLWtW715B//70Gus3LiN759xZErfG76s4+l7avhy6UYOmLw351x2HL37V2X5jYtNTt0DuwXdQUvvfLqSmx56mZZgKG7bs3M+ZXtDE7+/9PS4bXPXrOOmZ1+lJdTRLhiJQAQen7+InmWlXH1swec6Ohj7aqkQ4TAE/FBV54tIFTBPRF5X1c8KcGxbsWHlZn5w0m8J+KMXv0aUXz91PQccMTqpza3/ea0t6HSmJRDipw/P5IWfXRrXUvt8zRaeeW9JQoG2BEI89No8Tps8hj0G9Ul47BWfrufG8+4i6A8RDkdYtmgtMx+fw92v/Ii+A3tm+pWLTrR7oPBiEZG+GewWUdWdBXeeObu1llSVXz45K+E1DdASDPFO7So+W7uZsSMGddj2f7PejQs67WkOhrjng4/5zuGH0qPEPjnJuhI7aynvznRV3aiq82P/NwC1QPpZYJ1wQgK+u254mF07Gmn2tdDsa6Glyc//Tb8naZqPLzdvZ/WWHSmPuXVXI7Vrt8Stf+ytTwiEw0ntwpEIT76zKOn2e37zHC1NAcKx50mhQIjGXc08NePNlOVpT1clCQ1E3HFLAdgAzCU6MS7ZkrxCi0AhtGTnJKGfrNqAryWQcp9AMMxj7y/ssG5jfQOfb65LWw6XCG8sXVGw8na1Xb62ibRUIPLWUkE7AEVkJNEHU1nNdHVKAr7P565EIx2DzPbN9fjqm6jqXRG3/7pt9XjdLvyJb3gAcLmEddvqGbtHxxbeio3bOmTv7UwoHGHFhq1Jt3+xeG28TTDMog+XJy9MO7oqSSgqhCxopQG1qpryoamIfGKF41zIRUt2TxK6cceutMeLqMY11rb4fHjdbvyh5A0xAH8oxMZdDQUrb1fa5WsL2FpLkklSvkwQkUrgbeA3qvrfBNunA9MBBg0adOhjjz3Wti0UChEIBHC73YTDYUpKSvB40sfETOxa003kYtuZNZ9v6JDmA8DlFvY5cI+E+zcHQqzesoNIpzoeUOGlrjF6HJcIIwb0oqK0pMM+a+vqoyPZUtCrooyh/RJ3m635YhP+lo5dEwJU9i5n8Ih+HdYnqiOrzsmUKVPmpZo53Wf/gVp973lx6589+u6UdukQkTJVbcl3n2KQq5ZyPWeZ2uZ7nfhaAqzbVh+nh85UlZcyol+vts/+UJgVW7fR2WxgqZct7Vp1IjC4ZxV9e5Tn9T1ztSukjjKxzUVL+eoICqOlgtzxxFLRPw08kkgoALF8QDMAJk6cqNXV1W3brGxN1NTU0N5XPj7n8xm/+OZfCPqDqIK3xMOlPz8n4fEhmon3xFvvYeuujmngvzdpGP/8eD0APXuU8dZvvh43kue9xau4+V8v0ZTkdqm8xMvfrj2HCaMS98R84F/M7df9m0C7QFlS5uVP/72Ovcd2tElUR111xxNBCIQL30rLJKDYJOjkrCWrW9f5Xif+YIjjfvaPpNc0QI9SL3d8+zSOGbtX++/L1L/cy4b6jncz144axl3L17d9LnG7efPa0xhQGd/7kEt5s7UrpI7ytQV7a6kQo9oEuJfo7dcfczmGUxLwHVI9lrtm/ZRnZ7xOY0MTJ37zKCYfPz7p/i6XcNPXjuPWR15L+EC1zOvhB2cdk3D46BHj9mSfof1Zum4LgWDHLoZSr4fx+wxh/D6JR8MBHHnygdz4x29y/+0vsnXjTkaMHsT3f/H1uKCTjK5KEopCSIs7j0dEXlTV+KFURSZfLdk9SWip18MVJ0xmxutzOozUbMXtEgb3ruLI/ffssF5EuObYw/nVq2/RnGRgQqnHzdR990kZdLItb1fa5WsL2FpLhbjjOQq4GFgsIgti625R1ZezOYhTEvCNHDuM6++8NOP9Tz5kP/zBML976i0AguEwLhHKS7zccNYxfO2IAxLauV0u/nH91/n1w7N445Mv8HqiLZdQOMwpk/fn5gumpJ2zcMypEzjm1AkZl7UzXZIkFAgVfwLpd4vtMAl5a8nuSUIvP34SOxubefy9hUQ0qgeI3ukM7dOTe646F7cr/vx/bfw4lm7ZyhPzF8cFn3Kvl9ED+nHbGScVvLxdaZevrZ21lHfgUdX3SPpqJwPAmYeNZdoh+/LOZ6vYstNHZdMWam47h7KS1NVfXurlN5efwo0Nx7Fo1UYEYcI+Q+lZkX6Gt1NRhHCRxaKqG4vqMAndQUsiwo1nHceFx0zgqQ8Ws3zTNnr1KOOMSWOYPGpE0saUiHDLSdVM3Xcf7nn/Y+auXY8I7DuwH1ccMYlTxu1Liduec1a6CjtryZ7TWnczIpFGmvxPs//QF9h38C6WLv4Gzf4H8XrOx+1KP6GzT1UPjjtonyKUtOuxau5BKyKyigQvtlLVvS1zaohjYO9STj3WRVO4BI94GdajJKOsA4ePHMHhI0cA0WcqV11QbXFJnYudtWQCT4YEIyHeq1vE42vfZE3jZhRlYGkfzh1RzfGDDqWHJ/4uRFXZvuuP7PD9DXCh2gRARE9j26672Fb/G3pXXk6/Xj9GJHnLZM2unXy2NTrX58ABgxlWZd8JoPljeSut/YieMuA8IJMJcYYC0Bzaxbt1D7Fk5yxc4iKiEURcRDTE0PL9OXbgdxjWY0xS+1UNW/n3ytl8vHU1ZzUN4Ok5j/HtfQ5nYr89TbqcOOyrpW4ZeBqaWnhl7lK+3LydvlU9OGXi/gzr3yvp/lv99fzwk7+wI+ijOfzVEOcNLVuZseJ5HvjyFe4YfxV7V371sF9V2bLjRhqanyXRAI/WILSz8T6C4bUM7vv3+OwF2+r46buzWFy3Ca8r2nIJRsIcOngYvz7mBPbundnvpao6RpSqEI5YV1ZV3dZp1Z0iMg+41TKnBgAagtt4eNX1NIZ2EiH0VVs59ndt02IeX30zpw69kf17HdPBNhQJc+snz/PK+k8Ja4SQRvDTlzc3LuXDLSvZo6Iv9xx1Ef1Kc3sesjtiZy11u8Dz7AdL+N0TbyESTT3jdbu455U5nDJpf372zRPiHmw2h/xc/8ldbPXvJKzxWaVbIgFaIgFuWPBXZkz8EQPLoilsGpr+Gws6zSnLo9pMY8vr1PsepHfVpW3rF9dt4hvPPU5TKDr01N8ui8GH69dw5tMP88w532R03/4Jj9vcEuChx2fz3KsL8DX66d+3km+cPYlzzzgEd5IkjHZAEau7Bw5p99FFtNXW7XRQKJbMXcWrT35M/Y5GDj1mX04851AqKhPd/Ud4fPXNNIZ2ECH5RNCQBnh5wx/oUzqUQWVfdS/fMv9ZZm38HH+k48ACBZrCAZY3bOFb79zH01O+R4UndQZ3JzXE8sHOWrLvL5AFvLtkFbc/8Rb+YKhtOGcwHCEQCjNz3lL+8PTbcTavbprDzkBDwqDTnpaQn0dWvw60drH9IW3QaUW1me0Nd7Wl3omoMv3VZ9uCTtz+QGMwwPdfez7h9kAwxDX/8yhPvjAXX2P0Dm3rdh/3PvIuv7jjhaQpfuxCJCJxSwH5Q7vlt8AhwPmFdNBd+OdvX+Sn0+/nzRc+Ye67y3jgTzO5Ytof2LJhZ9y+q3zzaAhuSxl0WglrgA/rHm37/OnODcza+Dkt4eTzf0IaYXPzLp78cl7SfRa8tYTvjLmOkzznc+GI7/Hmo++mLYvTsVBHkIeWbBN4ipEH6a7n3k2eoDAQ4r/vL6ah6atuMVXlybVv4Y+kyHkTI0yENzbPpTnkxx9cSCgSn38tFRFtoNn/PgAfrF/DrkDqrAUKrPc1sHBL/CCS12tqWbthB4FAR5G3+EPMmbeKT5duyKhMXZGbKto94IpbCoWqTmm3nKiq01V1acEcdDHFygv22SereeWJj/A3B9syCvibg+za0cifb42f9/rRtqcIZtoQQ1nh+4imUD0ADyz/kEAksW7b44+EuP+LD4gkaCSu+Xw9Pz3jt6xbugEUtq7fzh+v+AdzX1uY4EiJcVqutkRaKiT5aMkWgad1hu66deuora3NuJKzsdu2q5HVW3amPJ7b5eLD2tVtnxtCTewIpM791MFe3Kxq3Ig/8CkJBnukRDWIP1gLwFurV9IYTB/s/OEQ765dHbf++ZkLaUkyO9wfCPLKG+lfn16Mc5IIRQiHXXGLlXTqMnAs+dR9trYzn/q4LUt7eyIRZdFHK2n0dXyuuaklsxyBrXikhK3+6LU9p25V2jQ7rewKtrDN3xi3/oW7ZxLsVF5/c4BHbnsqo+N2hR6s0JLVZKol2wQeEaGqqgoRyeqkZmrnD4Zwu1LfaqrSYUZ1KBKd7JkpAoQ0jBKANF1z8YRRjd7lNKdI/96eiCr+cPy+TU3J75ZUocGXPjNMMc5J4gJCRCVusZjvW+2gGORT99na7trZlLTLVkTwN3ds+Kim72LrsD8Q1ugxQlnYukQIROL337GlnkiCN//u3Fyf0XG7Qg9WaKkIZKQlWwSeyspKVJWGhgZUNat0FJnaDehViTtNEFFV9hs+oO1zlTf1S946E9QwA0p743YNBMnunSAiZXjcAwHYt28/yjNIJNjD62Vkr/j38Rx84B5J3+JYVubl0PF7JtzWnmKck2RoROIWK1FVu2QuyIt86j5b24nH7EtZeeJrvKpXOX36d7Qv9yQfNZoI1TCVnmgy28HlmduGNULfknjdHn3OZMoqOg468JS4Ofrrh2V03K7QgxVasppMtWSL0TzFyIPk9bj5+tEH8dg7C+Jyn0E0s+3w/r3Yb8TAr2xcHo4dMIE3N88nQvo7mJEVgxlS3o9IZApksH8HNERF+ckAnLPvWH47O36gQ5yJwql77xu3/vyzJvLqm0sIdWrhiUCJ182Jx41Ne+yuytWmChHru9b6AKOJzj2I+dV3LHVaBIqZF2zqGQfz6D/eIhAIEQl/dedTWubl8htPiRs1NqHPaXxQ9wghTf0+nlYqPH3pXxptIH17n8P59cKXaQqntnUhnDBkDOWekrhtx5x7OG8/+SEfv7oAjURwuV0M23cI3/zx1zMqjxNztdlZS7YIPFCcPEhXnn4EHy9by6pN2zsMMvC6XZSXevn9d8+Iszl/xFTerVuIP5I6kJS5SvjWntFcUS5XD6p6nMeuxv8A6Z/VgJuK8tNwu6Itu16lZVxx0ETuWzwvabdbucfL9ROPoNwb3+ocPrQPv/3p1/jZ754jElHC4Qgut1BVUcYdvziXHuXxwkxEV+Smghx6KbNARK4ArgOGAwuAw4EPganWeS0excoLVl5Ryp+fuJq7fv4M89//AoDefSu57MZpTDl9Qtz+B/U+mQ/qHsno2F4p47B+57YFr1OGjeP3S16nORxI+eS0xO3hu/senXCb2+3m50/dyNKPl7Ns7gqG7zeU8dXjcCXIC5cMp+VqA/tqyTaBpxiUl3i5/4ff4Ml3F/Gft+ZTV99Ij1IvZxw2lm+fMJGBveNP8F6VQ/h/o8/lL188lXR0W6nLy9nDjuHI/l8l/OzX84c0Nr9MOLKN1Hc/gsvVk/69bumw9sbJR9MUCvKfzxYSUY2+Wx7wuty4BK446FCmj5+U9KiHjt+T5x66mtnzVrJtu48Rw/py8IF74ErznKvrsbxL4DpgEjBbVaeIyP7AbVY63F3pN7An//v3S2hpDhBoCVLVu0fS+TE9PL04eci1zNz4F0Ka/BmkR0oY1mMMB/b5KuFnqdvL/UdfwsXv3kdTKECo06+pAKVuD78++Ez26zU4ZZn3mzSK/SaNyvxLOhr7aqlbBR6Ipma/aOohXDQ184FMJw+ZTP/SXsxY8Tzrm+twiYvWMaQ9vRVcutepnDC447uVPO7+jBj4AmvrziESqW/LVNAekR64pILhA57C6xnaaZvw86OmcskBB/Pgkk+Yu3E9CBw5dA8uPmACw6vS93t7vW6OOXx0xt/TFihWi6VFVVtEBBEpVdXPRWQ/Kx3u7pSVl1CWwV30uN7HI7iZufFOQAi2y+jhwoOIMLrqSE4ddgMu6TjxcXTPgTw79fvMWPouz65dgFsUdxhKxMXhA/bh6v2rObBPVm8J3/2xsZa6XeDJlUP77sc/+/6IVb6NLG1YQ0SVkRWDGdMzeY4or2cP9uz3b3bV/5Sd/o8IE6E1+bAHF728Y+nV+1e4vcmDw8hefbj1yCk0xSbP9XB7d/9Z19aKZZ2I9AaeBV4XkR1A/Jh0gyWM7V3NqJ6H83n92yzY8TLN4Xo8UsrIykM4tO9Z9C5JfMeiqgySOfx0+P38cMBKVvt7s/qLQbwx7kX6lJRDSQjVbyCy+2ZuzwmbaskEngxpDgV5dvVC7ln6IRua6kGgT0k53x51GBfucyi9S+Jft6v+d5Gd19BLA/R0lxNGiQBrcbGHuxyJLIPt30Z7/x4pOzHOvj7QzH9WzuOh5XPYGYi2DvuV9uA7ow/nG3sdQqU3dWoQR2JxK01Vz4n9+wsReQvoBbxqmUNDHCWuMg7qczIH9Tk5o/1Vw2j9/4D/ddBmyt2wf49GtrgD9PHshMhOaPgD2vwE9P034kqcw3DVth08vfhTNtTvol+PHpx1wBgOGDKocF/MbthYS90y8KxYt5UnXv+E5Wu30q93BV+bchCTx+2Z9PnHdn8TF771AJuad9HcmrZDYZu/ibtr3+XB5XN4tPpS9qrq12ajgYXojquBaMAQETyxux0RYnctEaAZ3flD6PsvpGRym/0a33YuqHmAhqC/Q36qLS0+/vxZDQ8t/4jHq7/D4B7JM1Vv2rqLl2s+Zcu2BkYO68cpx42lV1Xy99HbhizFIiIjgIeAQUSngMxQ1T+ns1PV9EMHDZbQ1ORn1lufsWz5Zvr2qeCk48cxfFjigKG7fgktr9GqpcS0QGgVuv0S6PcUIl81ypqDQa5/9iXe/3IN4YgSikRwAY8vXMyo/v2Ycd5Z9K9I/eZSx1KEIdSQvZa6XeB54IWPuPf52YRCYcKR6HOaOUtWM37fofzx+rPxeDr2Lasql737COsadxBMMETEHwkR8If4Vs2DzDrlGnrEhnLqrv8ltVDa04LW34oMiDYW/OEQ33r7Ibb7m4gkGMfTEg6xpaWBi955iJknX4U7wSsVHnpmDvc99SGqSjAUobTEwz8fe49brzmFKYfHD8G2DSpIOGuxhIAfqup8EakC5onI66r6WesOIjJfVVM+2MtkH0P+fLJwDbf84mlUlZaWIB6Pi8ee+ogzT53A1d+b2qErWUOrofm/QOoUUlFCEF4DzS9Cj+gw6YgqVzz+DAs3bOqQaDcCNAdD1G6u4/yHHueFyy+ioiSz0Z6OITctpaUQWirIIG8RmSYiS0VkuYjcXIhjWsHc2jXc9/xs/IFQW9ABaPYH+WTpemY880G8zdY1fNmwLWHQaUWBplCAl9Z+Gv0cWg6h7FKEENmIBhcB8Mq6z2gI+RMGnVbCqmzzN1Kz8Yu4be/PW8ED/51NIBgmGIqW2x8I4Q+E+OVfX+HLdZ2zmduMSIIlBaq6UVXnx/5vAGqBzk+ax4jIohTLYiBxqu8i4hQt5cr2HY3c8vOnaG4O0NISy0wQihAIhHjx1YW8PHNRh/218SGymhOnzWjjjLaP769aw5JNWzoEnfaEIhHqfI38d9FnCbc7nix0lAV5aynvwCMibuBvwCnAWOBCEUk/Q7ETxUjA98ALH3VIidMefyDEE7MWEgx1vEAfWv7RV91rKWgKB7lv2WwAtPkFMpu/0w71o03PAnD/F7NpCqWfaNcYCvDAF3Pi1t//9GxaEuTRAgiFwjz64tyMitQVSRFRot0DnRfoLyJz2y3TE5mLyEjgYKBzxewPnJFiOR04MvsCF45CaKkrElJmY/f8SwsIJ5kT19IS5KH/fNAxFU/Lc2StpfDG6J0ScP9H82hKk/ewJRTivo+SZ7VuxWlJQhNqqTDkraVCdLVNBpar6koAEXkMOAvIuAnRmgxPRFBVxowZk9GkqWztPv8ydcZoVWXztgaGD+rdtm5Vw7aM031uat4V/Se8ETJI/96RCESimaY3NGWWPwpgTeOOuHVfrK5Lun84oiyoXZ/2uMU6J4mQxL9LW1V1YsItrXYilcDTwPWquqv9NlV1wsi1vLSUT90X63wvWLQmLmt6e+q2NuD3hygri02M1syT9LYhXohsBfZk+dbtGZls3JXaT1fowUItJd8/g+elhdBSIQLPMGBtu8/rgLgESLEW6nSAQYMGUVNT07YtFAoRCARwu92Ew2Hq6+vxZJCrLBM7n8/X5uuCIwYSSnLLHSsjny+Zz/Lar24Ez/H3p9mdWa4ot7iivsITQPdIup+veRDvLLk2QQF6gbuGy3UPQu7MrpiSkLtDXQJcdOLQDl2JnSkr9XawaV9HrVh5TtIhObTMRMRLNOg8oqrxefmdQV5ayqfus9VSNnbtmXRwOWNHp8gVKMLs2e+3c3AtyTK9J9URbnBvAanhokF9CPRP/6p4EeK+W3uKVT+Fsm37XtlrKe3z0kJQtMEFqjoDmAEwceJEra6ubttmZWuipqaGVl+LH3uHx2rmtz336Mw+w/vzvUs6ZntY+GkN9yxN/z4QAU4YtB/fP7IabXoM3fUvIH7SKMA7S67l2APu6rS2DCqvw1VZzQsfPcNLaz9N+YwHwCMuvrHnIUw/uLrD+tn/fI2X3lqSMPiUlXq45qLjqK6e0LaufR210mV3PErWfdESfRp9L1Crqn/Mztp5JNOS1a3rQlwn776/jNvueJHmlvjuLxE4bNLeTL/iKx+Rut9BeGXCYyXWEUAJMvB9xNWL2W++w0NzFxBM1eAEjtl7JFd3+m7tKVb9FMoWyElLqroR2Bj7v0FEWp+X2i7wrAdGtPs8PLYuY4qVgO/CaYfw3DtLCIVb0E6/yaUlHn7wzePibfY+lHuWxg866EyZ28sV+8W6NctOh13ZZmFRpMfXAPjO6MN5fcPntCR45UF73OLiklGT49Zf+vXDeWv2MnxN/g7f0+N20bd3BaccNy5taboqSSiAZNtLCUcBFwOLRWRBbN0tqvpy1kfqWvLSUlckpMzW7ojDR7HHiH6s+rIuLllvWamX711W3WGdVHwX3fUrkjXi4nFB6VQklvfw4kMn8PC8BSktyrweph+RPP0UODNJKCTUUn8Raf+Qd0asIRNvm/x5ad4UYlTbx8BoEdlLREqAC4DE72ROQWVlJYMHD866crOx69+7kntvvZBRwwdQWuKhoryEHmVe+vbswW1XncbkcfFdAAPLq7hmzLGUu5O/5qDM7eWEYfsxoW90IJW4KqH8LCDTCZ4lUHYi4oq+4uCAPkM4Y8QBKX2Wu718c5+JjGw3d6iVwf178q/bvsWEMcMp8brpUV5CidfN0RP34V+/+SblZZm9sqEY56QzotHugc5LKlT1PVUVVT1IVSfElg5BR0QaRGRXgqVBRHYlO3aRyVtL+dR9Mc63x+3izjsu5ISp4ygp8dCjRwler5sx+w3hrt9/i71GDuhoUH4qJJgukJwSpPKKtk/DevXkFydNpSxJF1W518PFh07gsD2Gpz1yV+ih0Foi9qy03ZIs6CR9XloILeV9x6OqIRG5BpgJuIH7VPXTfI9rFXsO7sMjv76YFeu2snbzTnpXlXPQqKEpk2deOeZovG43dy6pwSXSNsqt1OVBUb4+cjw/nXByh/kHUvVjNDgfQquBVCPUvOAejPT8ZYe1vzrkdHp4Snh05TwEaZtEWub2EFHl0tGHcf3Y6qRHHTGkD3/7xTfYvrORHbuaGdivkqoKZ6QTyfaBaCaoalXhj1pYnKalXCkvK+GmH5zCNVcez+bN9VRVldG/X+LTI1IOff6Jbr+c9PPiyqDyKsR7UIe1544/gCE9q/jD2++zrG4rXrebUDjC4KpKrj3mCM4Yt39hvpgNyUVL6Z6XFkJLBXnGE2tdOqpbY5/h/dlneObTNi7f9wi+sdchPLt6ER/VrSZMhHG9h3D+XgfTvyy+NSKuHtD3cXTnVRBYQHRIaPv7XhdQCt79kT7/jN4ltd8qwk/Gn8z39juKJ1Z9wqId0R6XQ/vtwbkjJ9CnNLOX1PXtXUHf3g6ala05dbVljYgMpOM7RNZY7zU9TtRSrvQoL4m/w0mAlEyCvveiO64CDRLf7VYGKFTdgKvi0oTHOGqvPTlqrz3ZuKuBrY2N9CwrY4/evXbvvIc5aCmX56W5aKnbZS7Ih0pvKReNmsRFo1L3B7cirkqk70NosBZtfAD8bwPNgAfKTkEqLkO8B6Y8Rv+ySq4ac0zeZXcSVtzxtB1b5EzgD8BQYAuwJ9EJp+kffBm6DCmZBAPfh5ZZaNN9sZ4EN7j3hPILkR5fb3uuk4ohPasY0tP2N78FIwctZfy8NB8tmcCTJZ+u3cyHy1YTiSjjRw5h8qgRaVtN4h2D9L6dppYAO3Y1EXF/gqv3lUUqsQPJdOJUbvyK6AurZqnqwSIyBbjIUo+GAuFl5urRzJj9NZZs2sy1I/tx6xPnc9nkCZw/vgc9drOMNwUhSy2p6nu0ptBPT85a6raBJxQMs3VzPT37VNCjIv0ggLpdPq7+13Os2rydYDiMKpSVeOjdo4y/XnE2o4ck77ar2+Hjj4/U8O4nK3C7XXzzqEG88ufnuOFbUxiSwRyDboVae8cDBFV1m4i4RMSlqm+JyJ2WejTEsWnDTp598iNqF6+jb/9KzvjaRA6etFfSRlxElRuef5k3l6+iuV0mgo0NPv7w9vv855NFPHHxBfQuT/4cc8WyTTz16Gy+XLGFQUN6c/b5k5lw6MhCfzX7YGMtdcvA89zDH/DQXa8TDkeIhJUpp4/n6lvPoqQkcXW0BEJcfNfjbNrZ0GFuTJM/SLM/yCV/fYJnfvRtBiV4g+n2XU1cfOvD1Puao7ahMKrKuwtWsmDZeh7+5cUMSvJgtZVgIMgX81fhcgmjD9kbd6dEprsTguXPeHbGRuy8AzwiIluARks9Gjowb84KfvE/TxAORQjFUlTNm72CKScdwPU/Pj1h8PnXnLm8+cXKhK+CbwmFWFdfz/975gX+/c3zEvp86Zl5/OPPrxEMholElBVfbGb+Rys57ZxD+d618a8k2R2ws5YKkiS0EBQrD9K7ry7m/j/OpMnnx98cJBgIUfPSQv76v88mtXllwVK2twaOTijQEgjy0NuJcz3d+9yH7GpsibONRBRfk5+/P/1+QrtWfDsbufLgH/E/J/2KG4//X66a/D80N2aa9dr+ubviiLXSOi8F5CyiD9p+QPTdISuI5pfaLbD7+Q74Q/zqlifxtwTbgg5E87S99foSZr+3LM4mFIlwz+y5CYNOK8FwhE82bGLV9vgUUnVbdvH3O1/D7w8RaafDlpYgLz4zj08Xro2zSYQTc7VZqCPIQ0u2CDytM3TXrVtHbW1txpWci90jf38Tf6dZ0wF/NPg0NiT+QX/8/YU0B5InGgyGIzz3ceKJvS+99xmhcOIzHo4os+YsTZo0EeDBnz/OhhWbafG14G/0s+azdTz8myeT7t+eXOs1H9t8fLZiZeBR1UZVDatqSFUfVNW7VNXm6bozwwnne/Z7y+Imb7fS0hzk2Sc+ilu/aMMmgpH0TfdwJMLMpfHZ2me9vKhj4tF2BPxBnnv647TH7go9WKGlQpKPlmwTeESEqqoqRCSrk5qt3dZNiRNwutwudmxLbF/f2Jz2uE3++Lk6qkpTgtQg7QlHFH+SjNkAa5duINRueygQZnXturTlgdzrNR/bfHwCbUNAOy+FQkS+JiJfiEi9DSeQ5oUTzve2rQ2EgslP6Na6+GSdDX4/rgyGPYciEXY2x2t1y+Z6gkl8qsKWJL8J7ekKPVihpUKSj5ZsEXgqKytRVRoaGlDVrNJRZGu3z5ihCdeLCAOH9k64bXj/xOvb078qfq6MiNCvV+o5NOVlXspLk2cSGF89jtLyr4breMs8TJhyQNryQO71mo9tPj5bsbir7f+AM1W1l6r2VNUqVd0tRng44XzvMbI/Hm/iZ5Qiwt6j4l9FPaiqklCKpLetlHk8DO0Zfyr3Hj3oq2zXnXC7XYzef0jaY3eFHqzQUoHJWUu2GFxQzDxI3/nBydz8nX916G4rLfPyraumJh1ccPGxh7Doy400JeluK/V6+NaxByfcdsFJB/Ov52YnvKsp8bo5d+pBKYdjn3/jmaz9fD1vPPIuAMdfdCxfu+a0pPu3xwm5u+LIIbFhlmxW1VpLPXQRTjjfB0/am549y2lpDsR1uZWUeDjvoiPibPYb0J+BlT1YvSP1nYmqctrY/eLWH3/ygdx79xsJbTweF2efl35eniNztdlYS7YIPEDOSfCytdt//Ah+e//l3Pf7V1nx+Ub6DqjiwiuncPyZiQMHwNH7j2TCXkOYv3IDLcGOAcTrdjG0TxXnHXFQQtsLTz6EdxesZNnqLR1eQldW4mHPIX257MzDU5bX7XFz0wPXcMM9VyIiWY9oy7Ve87HNx6cALmvFMldEHgeepd37lB38KoUO2P18u1zC7/5yMT+6+kEafX6amwKUlnpRVa6+cRr77h/fIyEi/Kj6GG584VVakgwwKPN4OOeAsfTrEZ/Ro0dFKb/54zf5yQ2PEolEaGkOUloa/en70c/OYtiI+HyH+X7PQtjla2tnLdkm8BSTMeP34I5/J3yBZUJcLuGvl5/NH154l6dnL8bjjvZQBkJhph6wDz8773h6JOkuK/F6+PvN5/H824t59LVP2LrTh9fj5vvnHs05Uw6krCSzhJ0ebzc6VdaKpSfRnCsntVunwG4ReJzAsBF9+fd/r2POB1+w8otN9OzVg+NOGEev3snTQJ2832jqfI387q13iChtrzlwiVDqcTNln7259aQpSe3HHTSCR1+4nnfe+Iw1X25l4KBeTDlpHD17ZZZ6yrHYVEvd6NcsP7weNzefU821px7Fp2s3E45E2H/YAHpXlGdk+/XjJ/D14ycAre/tOMTaAjsViye9qep3rDu6IVPcHhdHHrsfRx4b3zWWjIsOncCUUXvz73kLeHvFKko9bk7dfzTfmXwoBw0ZnNa+vLyEk0+fkEepHYaNtWQCT5b0KPUyaVT6FOqG3HFZOOlNRBK9OawemKuqz1nn2VAIhvXqyc1Tj+XmqcdSU1PDd1O8vM1gXy3ZYlSbwdBG6wPRzkvhKAMmAF/EloOIvnDtcpM6x7BbkUhLhSVnLZk7HoOtECzPL3UQcJSqhgFE5O/Au8DRwGJLPRsMRcTOWjKBx2AvFCSDORt50AeoJNolAFAB9FXVsIj4k5sZnE4kEmHB27WsXLyWIXsN4LBp43fvQTs21lJetS4idxDNzRMgmqfnO6q6M59jGgxW9ksTnfS2QERqiDYKjwVuE5EKYJalnlNgtGQtLY1+bjr9dtYs3UgwEKKk1EtVnwr+NOsn9Bvcu6uLZxl21VK+z3heBw5Q1YOAZcCPcz2QkxLw5eszFAoVzafT6sfqJKGqei9wJNG5B88AR6vqv2J5p35UOE9ZUxAtOe185+MzGx39544XWPXpeloa/YSDYZp9LWzdsIO//OChopR1d0wSmo+W8rrjUdXX2n2cDZyby3Fak+GJCKrKmDFjMpo0latdV/sMBALU1tZa7tNp9QPW9UuLyP6q+rmItI5jb01JPFhEBqvq/MJ7zZxCaMlp5ztfn9no6K0nZhP0d8w8EglH+GhmNIFoupc5Oql+WrGzlgrZwXkZ8HiyjSIyHZgOMGjQIGpqatq2hUIhAoEAbrebcDhMfX09Hk/6omVi5/P5Ovgqhs90tqFQiLq6Ost9ZmqXqI66on6AWGJDS/qlbyB6/f0hsVemWuE0R3LSUiGuzWy1ZLXPVHbZ6Oj4ayYQDMT3O4nA22+/XZCy2qV+2rCxltJ+CxGZBSSanfWT1rHaIvITIAQ8kuw4qjoDmAEwceJErW43/t7K1kR0smZ1Tra5+kxnW1dXx4ABA2xzx5OojrrqjgeseXmVqk6P/U0+vd1irNaS1a3rQl4n+djmoqOH3n+GF+56lUC7HI0utzD55IO44rpLClJWu9RPe+yqpbSBR1VPSLVdRC4FTgeO12QvvUiDkxLwFcJnfX19VhdRd6kfwPKROCJyHvCqqjaIyE+BQ4BfqeonljmNYbWWnHa+8/WZjY6+ccNpLHpvKcsXriYSjuByC30G9uL6uzKbfO+k+mnDxlrKd1TbNOAm4DhVbcrnWE5KwJevT4/Hk7V9d6kfwfKROD9T1SdF5GjgBOAO4B/AYZZ6TUOhtOS0852Pz2x0VFpewh0v/w+fzVkeG049kIOnjMXtznx8lZPqB+ytpXyf8fwVKAVejz2cm62qV+Z5TEN3RtXquQetUjwNmKGqL4nIr610mCFGSxYjIow7fDTjDh/d1UUpDjbWUr6j2kblY28wJCKXkTgich/RbqotqprqTXnrReSfwInA7SJSig1SRxktGazA4swFOWupywVnMHRAQUIat2TAA8C0DPY7H5gJnByboNkX6Mr5OwaDNSTQUoHJWUu7cb4Ig1PJpXtAVd8RkZEZ7NdEu/eFqOpGYGPWDg0GB2BlV1s+WjKBx2ArJPlInP4iMrfd5xmxYcUGgyEBKbTU5ZjAY7AXse6BBGxV1YnFLo7B4FiSa6nLMYHHYDMsH4ljMHQT7Ksl2wwucFICPpMk1DqfrWk+Oi+GzHDa+S5WktBC+XRK/QAJtZQOEblPRLaIyJLsHWaOLe54nJSAzyQJtc5nGzm00kTkUaCa6LOgdcDPY9lzuw1OO9/FTBJaKJ9OqJ8OZK+lB4jOKcs8bXcO2OKOx+fzISJUVVUhIhlH91ztutqn2+0uik+n1Q8Q65eOxC1pzVQvVNUhqupV1eHdLeiA8853vj6z1VEhfDqhftpIoKW0JqrvANuzc5Q9trjjqaysRFVpaGhAVbPKg5SLXVf7DIfDRfHptPoBEBSJWDvrbXfFaec7X5/Z6qgQPp1QP63YWUu2CTxOScBXCJ8mSWgKFDDPdHLCaec7X5/Z6qgQPp1QP20k1pItpiXYIvCAsxLw5evTJAlNhSJhazMb7s447Xzn4zMXHeXr0yn1EyWhlmwxLcE2gcdgAMwdj8FQKGysJVsMLjAY2iORSNxiMBiyJ1sdxUaHfgjsJyLrRORyK8pl7ngM9kIVQqarzWDImxy0pKoXWlSaDpjAY7AXCoTNHY7BkDc21pIJPAaboRAxdzwGQ/7YV0sFecYjIj8UERWR/oU4nqEbo0S7Bzov3QSjJUPBSKQlm5B34BGREcBJwJp8juOkPEgmV5t1PlGFcDh+6QYUQktOO98mV5t1tgm1ZBMK0dX2J+Am4LlcD+CkPEgmV5t1PtvovqPY8tKS0863ydVmnW0bNtVSXoFHRM4C1qvqQhFJt+90YDrAoEGDqKmpadsWCoUIBAK43W7C4TD19fV4POmLlomdz+fr4KsYPtPZhkIh6urqLPeZqV2iOuqK+gFAFbVRy6xYFEJLhbg2s9WS1T5T2WWro0L4dEL9tGFjLaX9FiIyCxicYNNPgFuIdg2kJZaWYQbAxIkTtbq6um2bla2Jmpoa2vsqhs90tnV1dQwYMMA2ra1EddRldzyqEAxlvr+DsFpLVreuC3md5GObq44K4dMJ9dOGjbWUNvCo6gmJ1ovIgcBeQGsLbTgwX0Qmq+qmbArhpDxIJlebdT5bsWsrLV+s1pLTzrfJ1WadbSt21VLOXW2quhgY2PpZRL4EJqrq1lyO56Q8SCZXm3U+ow9E7dkvbRWF1JLTzrfJ1WadrZ21ZObxGGyFqqKhYFcXw2BwPHbWkqgWP4mciNQBq4vkrj+Q012YhditTMUsz56qOiDZRhF5NVaezmxV1WnWFcuZdHMtdffy5KIlW+ioSwJPMRGRuXZIA94eu5XJbuUx2BO7XSemPM7FZKc2GAwGQ1ExgcdgMBgMRaU7BJ6iv9Y1A+xWJruVx2BP7HadmPI4lN3+GY/BYDAY7EV3uOMxGAwGg40wgcdgMBgMRWW3CzwicoeIfC4ii0TkGRHpnWS/L0VksYgsEJG5FpRjmogsFZHlInJzgu2lIvJ4bPscERlZ6DK08zVCRN4Skc9E5FMRuS7BPtUiUh+rjwUicqtV5TE4A6OlhGUxWioEqrpbLUQTLXpi/98O3J5kvy+B/haVwQ2sAPYGSoCFwNhO+1wF/CP2/wXA4xbWyRDgkNj/VcCyBOWpBl7s6vNnFvssRksJy2O0VIBlt7vjUdXXVLU1JetsogkXi81kYLmqrlTVAPAYcFanfc4CHoz9/xRwvKTLh58jqrpRVefH/m8AaoFhVvgy7D4YLcVjtFQYdrvA04nLgFeSbFPgNRGZF3u/SSEZBqxt93kd8Rdn2z4xcdcD/Qpcjjhi3RAHA3MSbD5CRBaKyCsiMs7qshgchdFSJ4yWcseRSUJTvddEVZ+L7fMTIAQ8kuQwR6vqehEZCLwuIp+r6jvWlNgeiEgl8DRwvaru6rR5PtHcTz4RORV4Fhhd5CIaiozRUm4YLeWHIwOPJnmvSSsicilwOnC8xjpdExxjfezvFhF5hugtfaHEsh4Y0e7z8Ni6RPusExEP0AvYViD/cYiIl6hQHlHV/3be3l48qvqyiNwtIv01x9dcGJyB0VL2GC3lz27X1SYi04i+t/5MVW1Ksk+FiFS1/k/0IeqSAhbjY2C0iOwlIiVEH3g+32mf54FLYv+fC7yZTNj5EuvvvheoVdU/JtlncGu/uIhMJnptWCZeg/0xWorHaKkwOPKOJw1/BUqJ3vIDzFbVK0VkKPAvVT0VGAQ8E9vuAf6jqq8WqgCqGhKRa4CZREfl3Keqn4rIL4G5qvo80Yv33yKyHNhOVFBWcRRwMbBYRBbE1t0C7BEr7z+ICvb7IhICmoELrBKvwTEYLcVjtFQATMocg8FgMBSV3a6rzWAwGAz2xgQeg8FgMBQVE3gMBoPBUFRM4DEYDAZDUTGBx2AwGAxFxQQeg8FgMBQVE3gMBoPBUFRM4DEYDAZDUTGBx2AwGAxFxQQeg8FgMBQVE3gMBoPBUFS6JElo//79deTIkUXx1djYSEVFRVF8ZYrdylTM8sybN2+rqg5Itv3kKRW6bXs43m6Rf6aqTrO0cA6kO2upu5cnFy3ZRUddEnhGjhzJ3Llzi+KrpqaG6urqovjKFLuVqZjlEZHVqbZv3R7ig1fj3yRcNnRVf8sK5WC6s5a6e3ly0ZJddFSwrjYRcYvIJyLyYqGOaeh+KBAiHLd0J4yWDIUgkZbsQiGf8VwH1OZq7PP52LRpEz6fryh2XekzFAoVzafT6kdRghqJW7oZOWvJaec7H5+56Chfn06pH0isJbtQkK42ERkOnAb8BrghW3ufz0dtbS0igqoyZswYKisrLbPrap+BQIDa2lrLfTqtfgAiQIvap2VWbPLRktPOd74+s9VRIXw6oX5asbOWCvWM506ir8itSraDiEwHpgMMGjSImpqatm2hUIhAIIDb7SYcDlNfX4/Hk75omdj5fL4OvorhM51tKBSirq7Ocp+Z2iWqo66oH4i10ujWLye8kxy1VIhrM1stWe0zlV22OiqETyfUTyt21lLegUdETge2qOo8EalOtp+qzgBmAEycOFHbP4SzsjWR7IFfV7Zg6urqGDBggG1aW4nqqKvueFQhaE+tWE6+WrK6dV3I6yQf21x1VAifTqifVuyspULc8RwFnCkipwJlQE8ReVhVL8r0AJWVlYwZMwafz0dlZWXGlZurXVf7rK+vz+oi6i71A6AILdptp5flpSWnne98fWaro0L4dEL9tGJnLeUdeFT1x8CPAWKttBuzCTqt5FKx+dh1pU+Px5O1fXepHwWCNhWL1RRCS0473/n4zEVH+fp0Sv2AvbXUJfN4DIZkRBACuLu6GAaD47GzlgoaeFS1Bqgp5DEN3QtFaImY9pDRkiFf7Kwle5bK0G1RIGjTVprB4CTsrCUTeAy2QlUIqj3FYjA4CTtryQQeg62IILRoSVcXw2BwPHbWkgk8BluhCAGbttIMBidhZy2ZwGOwFdEhoOayNBjyxc5ass0gbycl4DNJQq3zqUT7pTsvhsxw2vk2SUKts02kJbtgi3DopAR8JkmodT4BIiq0RLwZ72/4Cqedb5Mk1DpbyE1LIlIGvAOUEo0PT6nqz7M6SAbY4o7H5/MhIlRVVSEiGUf3XO262qfb7S6KT6fVD5g7nnxw2vnO12e2OiqETyfUTys53vH4gamqOh6YAEwTkcOzcpwBtrjjqaysRFVpaGhAVbPKg5SLXVf7DIfDRfHptPoBe/dL2x2nne98fWaro0L4dEL9tJKLllRVgdYI540tBU81aguFOykBXyF8miShyYkmNrTFZek4nHa+8/VpkoSmJlctiYgbmAeMAv6mqnOyPkgabKNwJyXgy9enSRKaHFUhaNM0H07Aaec7H58mSWhqkmipv4jMbfd5Ruw1G+3sNAxMEJHewDMicoCqLsmpEEkwCjfYitZ+aYPBkB9JtLRVVSdmZK+6U0TeAqYBJvAYdl8iCH5zx2Mw5E0uWhKRAUAwFnTKgROB2wtdNqNwg62IvjXR3PEYDPmSo5aGAA/GnvO4gCdU9cVCl80EHoOtUIRQxAQegyFfctGSqi4CDramRF9hAo/BVmhu3QMjgIeAQUSHfs5Q1T9bUDyDwTHkoqVikXepjOgNhSTaPZD1vOYQ8ENVnS8iVcA8EXldVT9r3UFE+mZwnIiq7szWeaEwWjIUkhy1lJZCaKkQ4TCt6DPB5/PlNF49V7uu9NmaY6oYPp1WPzl2D2wENsb+bxCRWmAY0P4a3BBbJMWh3MAeWTkvLHlryWnnOx+fuegoX59OqR+wtNs6by3lHXgyFH1KnJQHyeRqs84nRJv5gcTdA2nnHwCIyEiifdSdJ73VqmrKvmsR+STjglpAvlpy2vk2udqss4WUWsqXvLVU0FKlED0iMh2YDjBo0CBqamratoVCIQKBAG63m3A4TH19PR5P+qJlYufz+Tr4KobPdLahUIi6ujrLfWZql6iOuqJ+IDrpLZS4eyDt/AMRqQSeBq5X1V2dNh+RgftM9ikKuWipENdmtlqy2mcqu2x1VAifTqifVlJoKV/y1lLBAk8a0RNrnc4AmDhxolZXV7dts7I1UVNTQ3tfxfCZzrauro4BAwbYprWVqI668o4nFMleLCLiJXr9PaKq/407rmpLWt8Z7FMMctWS1a3rQl4n+djmqqNC+HRC/bSSq5bSHrcAWipI4Ekn+nQ4KQ9SIXyaXG3JifZLZycWERHgXqJdAH/Myjhq/6Kqnp6tnRXkoyWnne98fZpcbanJRUv5kqmWCjGqLS/Rt+KkPEj5+jS52pKjCoHsH4geBVwMLBaRBbF1t6jqyxnafzdbh1ZQCC057Xzn49PkaktNjlrKl4y0VIg7nnxFbzC0kUsrTVXfI/UIm3T2G3O1LTBGS4aC0RV3PJlqqRCj2vISvcHQmbA1D0QBEJFVJHi/iKrubZnTDDFaMhQau2rJntNaDd0WVSEYtrR7oP3IuDLgPCCTCXEGg6Ows5ZM4MmCplAzb9fNYeHOz4hohH2r9uaEQUfTu6RnWlt/OMj7dcvY6m+gPNhMcyhAuaekCKV2FgqELeweUNVtnVbdKSLzgFstc2qIQ1VZuHETa+vrqSgp4Yg9RlDu9XZ1sXYr7KwlE3gy5N26j/jHiocRBH8kAMCS+mX8d90rnDNsGueOOJXos+GOqCr3r3yb+1a8jQsIaYRL/Ptzwpu3cdHIo/ne6Km4pLj9sLZGIazW9TaJyCHtPrqIttqMDorI26tW8bPXZrGzuQUREISwRrjo4AnceMzRuF1GDwXBxlrqloILhyN89MkqPv18A0MG92bKUfvRozz53cdH2xbwjxUPE4gEO6wPavTzcxtew+vycPbwk+Ns/7z0VZ5cM4eW8Fe2EZSWcJCHV73HtoCPnx5wdmG+2G6AIoTClv7w/KHd/yFgFXC+lQ4NX/HG8hVc98JLtIRCcdsenr+AtTvr+cuZpydsxAHs2NnI3x98m7feW8p504byzKzHufqyavbde5DVRXccdtZStws8gWCIG259gi9WbqG5JUh5mZcZD73DP39/EYMH9orbX1W5d9XjcUGnPf5IgCfXvcTJQ46j3F3Wtn5d03aeWD0bfyReZAAtkSAvr1/ABXsewagqIxyIDgGNWNs9MMWyg3dTFr2/jO1b6jnkuDH07Jt86G8oEuGmV15NGHQAmkMh3l71JbPXruWIPeLTfDU1B5h+48Ns3e4jHI6gwCeL13DNjx/l77d/i31GDijUV9otsLOWbHNP6/P52LRpEz6fz1K7V95YwtIVm2luiQaS5pYg9buauXPGGwn3r921nKZQc9rjCsL7W+d2WPfk6tmENW7QRwdCkRCPfvlB2uMXq3662idAOCJxi5V06jJwNMU+3/ff9iy3Xnw3f/7hI0w/9pfUb0tu+9aKlYQiqfXQEgxy38fzEm57reZT6nc1EQ5HOqz3+4Pc8/A7GZe52Nd1V2iwlWLqCDLXki3ueIqZnqXm/aX4/R1bXBFV5i78MuH+G1o2o/EjBuPwRwKsbdzQYV3trvWENJzSLozy+a4NKfdxUoLC/NN8iKUPRJPwfWwyiTQfuuJ8PzfjLfyxRlw4HOadFz7mjEsTN4S/2LaN5mDyngOIPhD/vG5rwm3vf7SCFn/83ZIqLFiyLm1ZobslCbWvlmxxx+Pz+RARqqqqEJGMo3sudr169iBR93FZaeIRNR5xk+nUihKXN+XnZJS6Usf/YtZPV/oEQEEjErdYiao6PuhA15zvHj2/6loWEbw9kg/fLXG7cSd5dtN5v0T0rCpLqF2A8rLMtFbs67orzkkbCbRkNZlqyRaBp7KyElWloaEBVc0qD1K2dueecQglJR1/6L1eF2edfGDC/cf12o9ImrsWgDJXKRP6jOuw7qQhB1LuTj1kuszt5cQhiX23Usz66UqfrUQiErcUEhHpIyKTReTY1qWgDrqIrjjfP7r7Eip6l+FyCwdV78NRpyTPln/cXnulHbFW4nZxyn6jE24746TxlJbEB5gSr5vTTzooo/IW+7ruinPSHit1BLlryRZdbcVMwHfA/sO48fsn8ed73iAYCqMRZdrUsVz2zcT1NaC0L2N6jmJJ/TIiRBLuA1DlrWRsz46COXHIgfyhNnW2ExfC6cNSd4s6KUFh3okNFdTC7gERuQK4DhgOLAAOBz4EplrmtEh0xfk++OixPDD3lxnZje7fj7EDB7Bo06akz3pc4uJbB09IuG38uOGcftJBvPDaQgKBaGOwrNTLyD36cdHXD8uovN0qSaiNtWSLwAPFTcB38pRxHH/sGOq2NtC7VznlZanvSq4ZfSk3LbyNhqCPcILgU+Yu5ab9r4wbAlruLuGvEy/h+x/fTyAc6mDrFhdel5s/HnoxVd6yzoeMw0kJCvPxCUIkbGmXwHXAJGC2qk4Rkf2B26x0WEzsfr7vPvtMvv7wo2xrauowus3tEkpcbu484zSGVFUltBURrr1iKiccO4bX3vqUXlVN/OyHp3HExH3wuDP/ge0uSULtrCVbdLV1BR63iyGDeqUNOgB9Snrxf+Nv4bB+B+MVLz3cZfRwl+MVD+N7jeF3B97MyIrhCW0P7LMHjx39/zhrxKGUujy4EAQ4Zeh4Hjnyaib16/IUYfbC+mc8La3vChGRUlX9HNivkA4MyelfUcGLl17MD44+kiFVVXhcLiq8Xs4eO5Znvv0tjh+1T9pjjN13CNd/7wQGD+zJMYeNzirodCtyeMYjIiNE5C0R+UxEPhWR61LsnrOWbHPHY3f6lPTiB/tdQUOwkS8b1xIhwvDyofQr7Z3WdniPvvzkgLO5ZdxZBCIhPnz3fa48qNryMjsWC2dbA+tEpDfwLPC6iOwAVlvp0NCRqtJSLp80kcsnpXyhbByqYfDXoI33Q3gVhC4msvX3UH4hUn4W4sr1zmA3JnsthYAfqup8EakC5onI66qa6PXrOWvJBJ4sqfJWcGDv/bOy0dAatOlhaHkJrzZC6Eoi9W8jFd9GPHtZVFKHooCF3QOqek7s31+IyFtAL+BVyxwaCoIGPkJ3/D8gANoYWxuC0DJouANtuB2tugFXxaVdWEqbkYOWYq812Bj7v0FEaoFhQFzgyUdLJvBkwc5AE0+vnk/NpmWENcKk/iO5YK9JDCmPz3gA0awH2ng3+P4BRIDWOQwhaH4cbX4K7XExUvWjpClCIDrj+8udOxCEkb177/a5rDT5GI7C+lF9uzieDJ1pDgSZ+dkXrN62k8qyEk4eO5rhfZLoyP8+uuP7QLK3KTdF/zT8iUhkJ66q660osiPJR0siMhI4GJiT1k+WWjKBJ0NqNi3lh3OfAlVaYilwPqvfyL9XzOaGcSdy0d7xo2q0cQb4ZgD+BEcMRZfmR1DcSM8fxu0RUWXG3I/559yPCUaiV1CZx8M1kw/jkgkHpwxWTkasGfY5X1VTDh/MZB9D/jy3oJZfvPgGItAUCOJ1u7jrzQ+Yut8+3P61kynxfPWzpJHt6M6rSR502tMMjfejJYcipcfEbZ295EseemUuazbtYEDvCi486VBOmLQvLtfuqSNIqKX+ItI+xcoMVZ0RZydSSfQV7Ner6q5O2/LWUkECj4hMA/4MuIF/qervCnFcu1Bbv5Efzn2qQ6JPgGAkOqTzT5/NYmh5L6YO+aoLTiPbwfcXIJD64NoMTQ+gFd9C3IM7bPrpG7N47vNamtuN/mkKBrnj/ffY6PPx42N2i+knHVGxqqttjIgsSrFdiHYVdCm7u5be/HwFP39xFi3Br67pYCwFzptLV3LDky/z1wvPbNumTY9l2WxvRn13xwWeu59+j0dfm09LIOp38/YGfnXfTN6cu4zbvn/67hl8Emtpq6qmfLgmIl6iQecRVf1vgl3y1lLegUdE3MDfgBOBdcDHIvJ8kodRSfH5fDmNV8/VLhvbvy99G384eaqPlnCQP302q2PgaXqSzAcNKtr0CFL11V3P8u3bePbz2oQJFZtDIR5a8AmXTDiYoUmGnrZSjPoppE+AFNOl8iGTB3PpZwpbSCG0ZOfzrar89tW3OwSd9vhDId5dvpoVddvYZ0A/VCPQ+CCZ3e20I7gEDa1BPNFEo8vW1PGf1+bjD3T02xII8cHiL6mZ/wVTJ+6b9rB2/o1KSpZakmg3yr1Arar+MclueWupEHc8k4HlqroSQEQeA84iwcOoZNg5D1JYI7y9aVnabG0bmnaytnE7IypiL+BreYXMBROAlpnQLvA89emnhMLJz52q8tzntXx/0uSk+zgxVxtqTVebqjph5FpeWrL7+V5Rt52tvsaE21oJRcI8t6CWG048GiJboj0C2SIeCC6CWOB5pmYRwWBiLTX7gzz2+idpA4+df6OSkpuWjgIuBhaLyILYultUtW0mfCG0VIjAMwxY2+7zOiDugYeITAemAwwaNIiampq2baFQiEAggNvtJhwOU19fj8eTvmiZ2Pl8vg6+svUZUeV7Ojpt4HFHhE9nz2OFO5bSI3QacELCfX3Ng3hnybWd1nrB81U5h+3axbVDhqb02Xvz5oTfrZVMv2OiOrLynKRDijS4wIbkpaV86j5XLWXjsykQZPqowUTSZKju3bg95icAoStJ1mxPrCMAF7iC4IqWdWhZA5cem1xLJV5PSh1BceonF5/pyFZLqvoemSanzIOiDS6IPcCaATBx4kStrq5u22Zla6Kmpob2vrL1GdEIP3rpt3HPdzpT4nLz2lGnMaAs2vUV2fZ3CH6ScN93llzLsQfc1XGlZz9c/b/X9vHuj+bwtzmz8Se56yn3eLj5mGOpHj8haZky/Y6J6qgr73goQjJDJ5NMS1a3rvO9Ttbv3MWP//IA/lDyO/kSj5urjzuc6mMno+FtaN0NJHtOmlBHAFKB9LoNKYuW9e6n3+Ph9+cSDMX/Cgtw7MH7cPm3quO2Zfs9C6mjfG0BW2upEIFnPTCi3efhsXUZY+c8SC5xcfaI8Ty1ej6hFA85D+ozvC3oAEj5+WhoKWhT+sJIOZSf12HVeeMO4C9zZic1UeCs/VJ3tToxVxt06zuevLRk9/M9rHdP9h88kIXrNiY/oMI5B4+N/u/qC+6hEP4y4+8RPUYISr66Ufxa9Xj+M3M+ie6cSks8XDQt/URWO/9GpcKuWirEhJCPgdEispeIlAAXAM9ne5DKykoGDx6cdeXmapeN7Xf3PYYenpKk95/lbi8/GndSp5WnZl4QVaT8nA6rBlRU8NPjqilPcGtd5vHw2xNOpGdZZjnerK6fQvoUBQnHL/kiIg0isivB0iAiu9IfoSjkrSW7n+9fnD6VHt7ErzAo93q4/KiJDKyKHkdEkIrvAj2yKI0LSqsRV9+2NYP7VfG/351GaYmHEk/0lQset4tSr4fvnnUEE/YdltGR7fwblYhEWioEhdBS3nc8qhoSkWuAmUSHgN6nqp/me1w7Mbi8Fw8fczlXzf4P2/2NNIcDKNDDXUKp28OfJp3PAX06Xrwi5dDrT+jO60g9yKAMev0WccWPTvvWQePZu08f7po9m/kboy+LO2z4cK49/AgmDs1MLI7EmsEFqYf/2YDuoKUxQwby0GXnccszr7Fmx048scnQLhGuqT6ciw/v9FqF8tPBdydEmiGDFzJCCVJ5Vdza4yfty0GjhvLM24tZsW4rQwf05OzjDmLPwX3y/k62xqZaKsgzntiIh9T5/x3OPlUDePWEa5m7bTWz61YSQTmw9zCOHTQajyvxi6ukbAr0uQvdeSMQbpfqA5AKQKDnb3GVn5zU7xEj9uCIEfHvn9+dKUb3gIgMBNpuGVV1jfVe09MdtHTA0EE8f/XFfLFlK2u211NZWsIhewzFm+AFcCLl0Pff6LZvgDaQenxwGfS6A/GOSbh1QJ9Kpp99RGG+hEOwq5ZM5oIsEBEm9R/JpP4jM7cprYaBH4L/DbT5GYjsAKlCev4Kyk4i2qNiaEOtFYuInAn8ARgKbAH2BGqBcansDIVn9MD+jB7YP+1+4tkb+j+H1v8YAvNRIkgs/ZTiQigF92Ck16+QkuTTC7odNtaSCTxFQMQLZdOQsmnRFe4apLy6C0tkbyxupf2K6AurZqnqwSIyBbjIUo+GOD7bsZlHv/iELxt2UOkt5YyRYzlx+Gi8yXoP3EP5wvtzXt78b/aUD9ijpJHmiJe3dg1jYWB/Dht0CUd5JxX5W9gfu2rJBB6DvbC4lQYEVXWbiLhExKWqb4nInZZ6NLSxK9DC9LefYuG2jQQjYcIafW7z7sZVeF1u/lV9LocOiH+31Ydb53H3igcIRIJ8SDSj+8nBcmZuGw2EWdj4MEsbVnDZXhfstjkMs8bGWtq90xwbnEkkwVI4dsYSIL4DPCIifwZST6c3FAR/OMQFsx7hk60baAmH2oIOQGMowM5AMxe/8Rifbt/UwW514zruXvEggUjyuXT+SIC36z5k1uZ3LSu/I7FOR5CHlkzgMdgKUXCF45cCchbQDPyA6LtDVgBnFNSDISEvrq5ldcMOApHkJ7Q5HOQXc1/vsO7Z9a8STBF0WvFHAjy17kUixXqvhs1JpKUCk7OWbNPV5qQEfPn6DIVCbfZ2L+tulCQUANX2Qwt50DpPXYOdz/c/P5tNUyh9AFm8fROrG3awZ1UfmkLNfLx9AZrRUGpoCfv5dNcyDuyVfHJ1sa9rJyUJzYZ8tGSLwOOkBHyF8BkIBKitrbXcp9Pqp5Vc+qVF5D7gdGCLqh6QYr+vAbcDA4lmTBFAVbVn9l7thd3P98pd2zIqi9flonbHFvas6sNW/3Y8Lg/BcOKM1p0Ja5j1TRuTBp5iX9ddmjIHy0e15awlW3S1+Xw+RISqqipEBJ/PZ6ldV/t0u91F8em0+gGiD0Rzy1zwADAtg/3+DzhTVXupak9Vrdodgg449HwnIrObmyRIysEFxf6eXXFO2rAoc0E7ctaSLQJPZWUlqkpDQwOqmlUepFzsutpnOBwuik+n1U8rEolf0qGq7wDbMzj8ZlWtzbpQDsDu53tUz/RzdiD6gsVxfQcBMKCsH6FIZnc7AG5xMax8cNLtxb6uu+KctCdbHWVJzlqyRVebkxLwFcJnfX19VrfN3aV+IJZfKnGLN6NX9mbAXBF5HHiWdu8kT/KmRUdh9/P9vXGH89OPXk37nGd8v6GMqOwNQLm7jMl9D+bDbfOIZPDAotxdxtieyd+tU+zruivOSSsptFQoctaSLQIPkPPDs5wfunWhT4/Hk7V9d6kfSNolkPaVvRnSE2gC2md1VcDxgQfsfb5P22MM93/+MUt31iUd2dbD4+UXkzom3D17+DQ+3rGAQCR14Cl1lXDeiNNxSeqOnGJf111xTlqxoHutPTlryTaBx2AALJ/0pqrfse7ohlSUuN08csI3+f47/2Ve3boOE0grPCWUuT3cU30uY/oM7GC3R49h/L9Rl/GX5fclnctT6irh+IFHM3Xg0ZZ/D8dgYy2ZwGOwHRaPxEnw5jDqgbmq+px1ng0AVd5SHj7+QpbtrOOx5Qv4smEHVd5STttzDFOHjWrLVt2Zyf0O5helN/L42uf4rH4ZXpcXlwglrhIGlfbn6yNO44h+hxb529gfu2rJBB6D7chxOPWjQDXRZ0HrgJ+r6r0Jdi0D9geejH3+OrAKGC8iU1T1+hyKbMiSfXsP4NaJJ2Zls0/lntwy5lq2B3aypmk9GxrW8NsDb2Z4j9SviO/OWJwyJ2ctmcBjsBWtL6/KFlW9MMNdDwKOUtUwgIj8HXgXOBpYnL1nQ7EJtnj4co2gAcG3y5vde+K6EblqKQty1pIJPAbbIRFLh+L0ASqJdgkAVAB9VTUsIv7kZoauxucPcPNzM3n7i1V43C6+O3Igv7j/Mfbs05s7zz2Nvfv3TWm/eOVG/vPGfFZt3M6gPlV8Y+oEjhi7526dVNSuWsor8IjIHURz8wSI5un5jqruzOeYhm6O9Rl1/w9YICI1RGdaHwvcJiIVwCxLPafAaCk1gXCYix98guV12wmEwwTCYVShJRhi2ZatnH/fYzw3/SKG9U48f3HGCx/y4My5+IPROXRfrNvKvGXrOG78Pvz68mm7Z/CxsZbynUD6OnCAqh4ELAN+nOfxDIZcMxdkROy5z5FE5x48Axytqv9S1UZV/VHhPGWN0VIKZn72BV9u20kgHH8xKNDoD3BXzQcJbRev3MiDM+fSEgih7TJiN/uDvL1gBa/M+dyqYnc52epIRO4TkS0isiTdvvloKa/Ao6qvqWrrtOLZQPyLNDLE5/OxadOmrNNC5GrXlT5bk4QWw6fT6qe1lVboGdcisn/s7yHAEGBtbBkcW9elFEpLTjvfmdo+OGc+TcHkE08jqrzy2TICofgsB4/Mmo8/mDj7QXMgyIMz5ybclmtZC2WXr20iLWXAA6RJPVUILRXyGc9lwOO5GDopAZ9JEmqdT4hlGrSmX/oGYDrRV/V2RoGpVjjNkZy05LTznY3tpl2Z/PAK9S1+BlR2/FlbuWEbmuKS2rhtV0HLWgi7fG0hNy2p6jsiMjLNbnlrKW3gEZFZQKLkRz9pHastIj8BQsAjKY4zPVZYBg0aRE1NTdu2UChEIBDA7XYTDoepr6/H40kfEzOx8/l8HXwVw2c621AoRF1dneU+M7VLVEddUT8AKLgyT82V+WFVp8f+Tin80TPDai0V4trMVktW+2zl0hH98Id6d1g3sNTLtaOHtX0WYMnHH8U9rzljbE98e5UlLYfX4074G5FtWbuyfhKSWEt5p54qhJbSfgtVPSHVdhG5lGg6+uNVk7crYl9uBsDEiRO1urq6bZuVrYmamhra+yqGz3S2dXV1DBgwwDZ3PInqyGmvRcj42CLnAa+qaoOI/BQ4BPiVqn5indcoVmvJ6tZ1Ia+TbG3XfrSAv77xLs3tusyuHT2Mu75YD0SDzjGjRnLVlPjfwvcWr+LmGS/R7I/vqiv1erj67KOork7dQ2T3+klGAi0VKvVUXlrKd1TbNOAm4DhVbcr1OE5KwGeShFrnE4j1S1s6BPRnqvqkiBwNnADcAfwDOMxKp+kohJacdr6zsT1n/FhmvPcR/lCYSIKYXOrx8IMpRyW0PeqAkRx9wF68t3glzYGvAlep182oYf04t/qggpa1EHb52gK21lK+z3j+CpQCr8dub2er6pW5HMhJCfhMklDrfIoqErZULK1je04j2s3wkoj82kqHGVIQLTntfGdqW1lawmOXXcDlj/yXzQ0+mgLRu5eKEi8ucfGX805n7JCBCW1FhNu+eyovzv6Mf8+cy6btDfSuLOf8KeM5f8oESr2Z/QzauX4SYWct5RV4VHVUPvYGQyIsnnuwXkT+CZwI3C4ipdjgvVRGS+kZ1rsnr1x1CR+vXs/by1fRr6meX55+AiftP4qSNM8+XC7hzCPHceaR44pUWnuQrZaySD0FeWjJZC4w2AvF6lba+USHi/5eVXeKyBCgK+fvGLJARJg8cjiTRw6PPlM5IPErrg3kpKUsUk9BHloygcdgO6wMPLHnJ/9t93kjsNEyhwZDF2JXLZnAY7AX1qf5MBi6BzbWkgk8BlshWN7VZjB0C+ysJRN4DPbC+iGgBkP3wMZa6vLRPK04KQ+SydVmnU9QJByJWwyZ4bTznY/PXHSUr0+n1E+UeC3ZBVvc8TgpD5LJ1WadT6AYo9p2W5x2vvP1ma2OCuHTCfXTho21ZIs7Hp/Ph4hQVVWFiGQc3XO162qfbre7KD6dVj+tSETjFkN6nHa+8/WZrY4K4dMJ9dMeu+rIFnc8lZWVqCoNDQ2oalbpKHKx62qf4XC4KD6dVj8Qm20dsk+XgJNw2vnO12e2OiqETyfUTyt21pJtAo9T8iAVwqfJ1ZaGiD3FYnecdr7z9Zmtjgrh0wn10wGbaskWgQeclQcpX58mV1sKbNwv7QScdr7z8ZmLjvL16ZT6AWytJdsEHoMBAFVI8Hpjg8GQJTbWkgk8Btth11aaweA07KolE3gM9kIBG803MBgci421ZAKPwWbYt3vAYHAW9tWSCTwGe2HjfmmDwVHYWEsm8BjshY27BwwGR2FjLRUkc4GI/FBEVET6F+J4hm5OJBK/dBOMlgwFxaY6yjvwiMgI4CRgTT7HcVICPpMk1DqfqEIoFL90AwqhJaedb5Mk1DrbhFqyCYXoavsTcBPwXK4HcFICPpMk1DqfUdS23QNFIC8tOe18mySh1tlGsa+W8go8InIWsF5VF4pIun2nA9MBBg0aRE1NTdu2UChEIBDA7XYTDoepr6/H40lftEzsfD5fB1/F8JnONhQKUVdXZ7nPTO0S1VFX1A8ACmrTB6JWUggtFeLazFZLVvtMZZetjgrh0wn104aNtZT2W4jILGBwgk0/AW4h2jWQFlWdAcwAmDhxolZXV7dts7I1UVNTQ3tfxfCZzrauro4BAwbYprWVqI667rUIaqsugUJitZasbl0X8jrJxzZXHRXCpxPqpw0baylt4FHVExKtF5EDgb2A1hbacGC+iExW1U3ZFMJJCfhMklDrfEbRnFppIjIN+DPgBv6lqr/L+iAWY7WWnHa+TZJQ62yj2FdLOXe1qepiYGDrZxH5EpioqltzOZ6TEvCZJKHW+YwOAc1OLCLiBv4GnAisAz4WkedV9bPcClFcCqklp51vkyTUOls7a8nM4zHYClUlEsy6e2AysFxVVwKIyGPAWYAjAo/BYAV21pKoFj+JnIjUAauL5K4/kNNdmIXYrUzFLM+eqjog2UYReTVWns6UAS3tPs+IPetARM4FpqnqFbHPFwOHqeo1hSu2PenmWuru5clFS0l1FLMpipa65I4nVWUVGhGZq6oTi+UvE+xWJjuVR1WndXUZnER31pIpT2rsrKWCZC4wGLqY9cCIdp+Hx9YZDIbsKIqWTOAx7A58DIwWkb1EpAS4AHi+i8tkMDiRomipOwwumJF+l6JjtzLZrTxZoaohEbkGmEl0COh9qvppFxdrd8Ru14kpT4Eplpa6ZHCBwWAwGLovpqvNYDAYDEXFBB6DwWAwFJXdLvCIyB0i8rmILBKRZ0Skd5L9vhSRxSKyQETmWlCOaSKyVESWi8jNCbaXisjjse1zRGRkocvQztcIEXlLRD4TkU9F5LoE+1SLSH2sPhaIyK1WlcfgDIyWEpbFaKkQqOputRBNtOiJ/X87cHuS/b4E+ltUBjewAtgbKAEWAmM77XMV8I/Y/xcAj1tYJ0OAQ2L/VwHLEpSnGnixq8+fWeyzGC0lLI/RUgGW3e6OR1VfU9XWPBGziY5DLzZtaSdUNQC0pp1oz1nAg7H/nwKOl3T58HNEVTeq6vzY/w1ALTDMCl+G3QejpXiMlgrDbhd4OnEZ8EqSbQq8JiLzYu83KSTDgLXtPq8j/uJs2ycm7nqgX4HLEUesG+JgYE6CzUeIyEIReUVExlldFoOjMFrqhNFS7jhyHk+q95qo6nOxfX4ChIBHkhzmaFVdLyIDgddF5HNVfceaEtsDEakEngauV9VdnTbPJ5r7yScipwLPAqOLXERDkTFayg2jpfxwZODRJO81aUVELgVOB47XWKdrgmOsj/3dIiLPEL2lL5RYMkk70brPOhHxAL2AbQXyH4eIeIkK5RFV/W/n7e3Fo6ovi8jdItJfc3zNhcEZGC1lj9FS/ux2XW0SfYnRTcCZqtqUZJ8KEalq/Z/oQ9QlBSxGJmknngcuif1/LvBmMmHnS6y/+16gVlX/mGSfwa394iIymei1YZl4DfbHaCkeo6XC4Mg7njT8FSglessPMFtVrxSRoUTfpncqMAh4JrbdA/xHVV8tVAE0SdoJEfklMFdVnyd68f5bRJYD24kKyiqOAi4GFovIgti6W4A9YuX9B1HBfl9EQkAzcIFV4jU4BqOleIyWCoBJmWMwGAyGorLbdbUZDAaDwd6YwGMwGAyGomICj8FgMBiKigk8BoPBYCgqJvAYDAaDoaiYwGMwGAyGomICj8FgMBiKyv8HNg96LSTbFa8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 8 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def rectangular_array(n=9):\n",
    "    \"\"\" Return x,y coordinates for rectangular array with n^2 stations. \"\"\"\n",
    "    n0 = (n - 1) / 2\n",
    "    return (np.mgrid[0:n, 0:n].astype(float) - n0)\n",
    "\n",
    "\n",
    "for i,j in enumerate(np.random.choice(nsamples, 4)):\n",
    "    plt.subplot(2,2,i+1)\n",
    "    footprint=shower_maps[j,...,0]\n",
    "    xd, yd = rectangular_array()\n",
    "    mask = footprint != 0\n",
    "    mask[5, 5] = True\n",
    "    marker_size = 50 * footprint[mask]\n",
    "    plot = plt.scatter(xd, yd, c='grey', s=10, alpha=0.3, label=\"silent\")\n",
    "    circles = plt.scatter(xd[mask], yd[mask], c=footprint[mask],\n",
    "                           s=marker_size, alpha=1, label=\"loud\")\n",
    "    cbar = plt.colorbar(circles)\n",
    "    cbar.set_label('signal [a.u.]')\n",
    "    plt.grid(True)\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "voaMgctlrNPe"
   },
   "source": [
    "## Wasserstein GANs\n",
    "To overcome the meaningless loss and vanishing gradients, [Arjovsky, Chintala and Bottou](https://arxiv.org/abs/1701.07875) proposed to use Wasserstein-1 as a metric in the discriminator.\n",
    "\n",
    "Using the Wasserstein distance as a metric has several advantages in comparison to the old min-max loss. The crucial feature of the Wasserstein distance is a meaningful distance measure even when distributions are disjunct. But before coming to the essential difference, let us try to understand the Wasserstein distance."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "DbpxbWzCrNPf"
   },
   "source": [
    "### Generator\n",
    "For training GANs we need to further define our generator and discriminator network. We start by defining our generator network, which should map from our noise + label space into the space of images (latent-vector size --> image size). Adding the label the input of to both the generator and discriminator should enforce the generator to produce samples from the according class."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "qwEBqLTJrNPg"
   },
   "source": [
    "#### Task\n",
    "Design a meaningful generator model!  \n",
    "Remember to check the latent and image dimensions. You can make use of the 'DCGAN guidelines' (See Sec. 18.2.3).\n",
    "Use a meaningful last activation function!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "XyXSTz58rNPg"
   },
   "outputs": [],
   "source": [
    "def generator_model(latent_size):\n",
    "    \"\"\" Generator network \"\"\"\n",
    "\n",
    "\n",
    "    return keras.models.Model(latent, z, name=\"generator\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "9aOBttl5rNPg"
   },
   "source": [
    "Now we can build and check the shapes of our generator."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "PmINju34rNPg",
    "outputId": "d40bcd7f-8b4f-4978-d741-d110bd40b185"
   },
   "outputs": [],
   "source": [
    "latent_size = 128\n",
    "g = generator_model(latent_size)\n",
    "g.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "lfASQfeorNPh"
   },
   "source": [
    "### Critic / Discriminator\n",
    "The task of the discriminator is to measure the similarity between the fake images (output of the generator) and the real images. So, the network maps from the image space into a 1D space where we can measure the 'distance' between the distributions of the real and generated images (image size --> scalar). Also, here we add the class label to the discriminator.\n",
    "\n",
    "#### Task\n",
    "Design a power- and meaningful critic model! Remember that you can make use of the DCGAN guidelines and check the image dimensions.\n",
    "Do we need a \"special\" last activation function in the critic?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "-RGjRFIErNPi"
   },
   "outputs": [],
   "source": [
    "def critic_model():\n",
    "    \"\"\" Critic network \"\"\"\n",
    "\n",
    "    return keras.models.Model(image, x, name=\"critic\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "bbSFKoCerNPi"
   },
   "source": [
    "Let us now build and inspect the critic."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "GwgSOjCkrNPi",
    "outputId": "9d37acba-df51-4da7-fe96-86a7deb2cb4b"
   },
   "outputs": [],
   "source": [
    "critic = critic_model()\n",
    "critic.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "WHH6tNHkrNPi"
   },
   "source": [
    "## Training piplines\n",
    "Below we have to design the pipelines for training the adversarial framework. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "9bc4CUC3rNPi"
   },
   "outputs": [],
   "source": [
    "def make_trainable(model, trainable):\n",
    "    ''' Freezes/unfreezes the weights in the given model '''\n",
    "    for layer in model.layers:\n",
    "        # print(type(layer))\n",
    "        if type(layer) is layers.BatchNormalization:\n",
    "            layer.trainable = True\n",
    "        else:\n",
    "            layer.trainable = trainable"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "9NsGMWYQrNPi"
   },
   "source": [
    "Note that after we compiled a model, calling `make_trainable` will have no effect until we compile the model again."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "hh1vHPBSrNPi"
   },
   "source": [
    "Freeze the critic during the generator training and unfreeze the generator during the generator training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "SabELgl3rNPj"
   },
   "outputs": [],
   "source": [
    "make_trainable(critic, False) \n",
    "make_trainable(g, True) # This is in principal not needed here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "dCWBEKcorNPj"
   },
   "source": [
    "Now, we stack the generator on top of the critic and finalize the generator-training step."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "M-zL_5tYrNPj",
    "outputId": "a82a240a-17c3-4306-a34e-0d15e639ac16"
   },
   "outputs": [],
   "source": [
    "gen_input = g.inputs\n",
    "generator_training = keras.models.Model(gen_input, critic(g(gen_input)))\n",
    "generator_training.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "d765RY6brNPk"
   },
   "source": [
    "We can further visualize this simple \"computational graph\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 312
    },
    "id": "msyfTsVJrNPk",
    "outputId": "10621acb-4271-4916-b456-408325e5f7e6"
   },
   "outputs": [],
   "source": [
    "keras.utils.plot_model(generator_training, show_shapes=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "b2thV6eOrNPk"
   },
   "source": [
    "#### Task\n",
    " - implement the Wasserstein loss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "2BbBke0wrNPk"
   },
   "outputs": [],
   "source": [
    "import tensorflow.keras.backend as K\n",
    "\n",
    "def wasserstein_loss(y_true, y_pred):\n",
    "    return \n",
    "\n",
    "generator_training.compile(keras.optimizers.Adam(0.0001, beta_1=0.5, beta_2=0.9, decay=0.0), loss=[wasserstein_loss])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "M_7mARHNrNPl"
   },
   "source": [
    "### Gradient penalty\n",
    "To obtain the final Wasserstein distance, we have to use the gradient penalty to enforce the Lipschitz constraint.\n",
    "\n",
    "Therefore, we need to design a layer that samples on straight lines between reals and fakes samples "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "VyjMTSgtrNPl"
   },
   "outputs": [],
   "source": [
    "BATCH_SIZE = 128\n",
    "\n",
    "class UniformLineSampler(tf.keras.layers.Layer):\n",
    "    def __init__(self, batch_size):\n",
    "        super().__init__()\n",
    "        self.batch_size = batch_size\n",
    "\n",
    "    def call(self, inputs, **kwargs):\n",
    "        weights = K.random_uniform((self.batch_size, 1, 1, 1))\n",
    "        return(weights * inputs[0]) + ((1 - weights) * inputs[1])\n",
    "\n",
    "    def compute_output_shape(self, input_shape):\n",
    "        return input_shape[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ARLKuVd0rNPl"
   },
   "source": [
    "We design the pipeline of the critic training by inserting generated (use generator + noise directly to circumvent expensive prediction step) and real samples into the sampling layer and additionally feeding generated and real samples into the critic."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "6D7zcsWErNPl"
   },
   "outputs": [],
   "source": [
    "make_trainable(critic, True)  # unfreeze the critic during the critic training\n",
    "make_trainable(g, False)  # freeze the generator during the critic training\n",
    "\n",
    "g_out = g(g.inputs)\n",
    "critic_out_fake_samples = critic(g_out)\n",
    "critic_out_data_samples = critic(critic.inputs)\n",
    "averaged_batch = UniformLineSampler(BATCH_SIZE)([g_out, critic.inputs[0]])\n",
    "averaged_batch_out = critic(averaged_batch)\n",
    "\n",
    "critic_training = keras.models.Model(inputs=[g.inputs, critic.inputs], outputs=[critic_out_fake_samples, critic_out_data_samples, averaged_batch_out])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "xXskjxTtrNPl"
   },
   "source": [
    "Let us visualize this \"computational graph\". The critic outputs will be used for the Wasserstein loss and the `UniformLineSampler` output for the gradient penalty."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 422
    },
    "id": "q8kzVv5-rNPm",
    "outputId": "aa5be2a6-359d-46e0-e2eb-adcb06cc53b4"
   },
   "outputs": [],
   "source": [
    "keras.utils.plot_model(critic_training, show_shapes=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "G9WcnuNLrNPm",
    "outputId": "d450187d-bbe4-4403-f688-118aadc070be"
   },
   "outputs": [],
   "source": [
    "critic_training.summary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "IULFk5ykrNPn"
   },
   "source": [
    "We now design the gradient penalty as proposed by in https://arxiv.org/abs/1704.00028"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "csNT3Jx7rNPn"
   },
   "outputs": [],
   "source": [
    "from functools import partial\n",
    "\n",
    "def gradient_penalty_loss(y_true, y_pred, averaged_batch, penalty_weight):\n",
    "    \"\"\"Calculates the gradient penalty.\n",
    "    The 1-Lipschitz constraint of improved WGANs is enforced by adding a term that penalizes a gradient norm in the critic unequal to 1.\"\"\"\n",
    "    gradients = K.gradients(y_pred, averaged_batch)\n",
    "    gradients_sqr_sum = K.sum(K.square(gradients)[0], axis=(1, 2, 3))\n",
    "    gradient_penalty = penalty_weight * K.square(1 - K.sqrt(gradients_sqr_sum))\n",
    "    return K.mean(gradient_penalty)\n",
    "\n",
    "\n",
    "gradient_penalty = partial(gradient_penalty_loss, averaged_batch=averaged_batch, penalty_weight=10)  # construct the gradient penalty\n",
    "gradient_penalty.__name__ = 'gradient_penalty'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "QCj2xws1rNPn"
   },
   "source": [
    "Let us compile the critic. The losses have to be given in the order as we designed the model outputs (gradient penalty connected to the output of the sampling layer)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "of1nxW_MrNPo"
   },
   "outputs": [],
   "source": [
    "critic_training.compile(keras.optimizers.Adam(0.0001, beta_1=0.5, beta_2=0.9, decay=0.0), loss=[wasserstein_loss, wasserstein_loss, gradient_penalty])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "z6bNJh9DrNPo"
   },
   "source": [
    "The labels for the training are (Remember we used as for the Wasserstein loss a simple multiplication to add a sign.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "cAp4NWkCrNPo"
   },
   "outputs": [],
   "source": [
    "positive_y = np.ones(BATCH_SIZE)\n",
    "negative_y = -positive_y\n",
    "dummy = np.zeros(BATCH_SIZE)  # keras throws an error when calculating a loss without having a label -> needed for using the gradient penalty loss"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "EM5pL3QvrNPo"
   },
   "source": [
    "### Training\n",
    "We can now start the training loop.\n",
    "In the WGAN setup, the critic is trained several times before the generator is updated.\n",
    "\n",
    "#### Task\n",
    " - Implement the main loop.\n",
    " - Choose a reasonable number of `EPOCHS` until you see a cood convergence.\n",
    " - Choose a meaningful number of `critic_iterations`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "1QjWP05WrNPp",
    "outputId": "a5d08d08-2eef-4bf6-fa12-9da898d51fa6"
   },
   "outputs": [],
   "source": [
    "EPOCHS = 1\n",
    "critic_iterations = 1\n",
    "\n",
    "generator_loss = []\n",
    "critic_loss = []\n",
    "\n",
    "iterations_per_epoch = nsamples // (critic_iterations * BATCH_SIZE)\n",
    "iters = 0\n",
    "\n",
    "for epoch in range(EPOCHS):\n",
    "    print(\"epoch: \", epoch)  \n",
    "    \n",
    "    for iteration in range(iterations_per_epoch):\n",
    "\n",
    "        for j in range(critic_iterations):\n",
    "\n",
    "            # complete the main loop here\n",
    "        \n",
    "        generated_maps = g.predict_on_batch(np.random.randn(BATCH_SIZE, latent_size))\n",
    "        \n",
    "        if iters % 300 == 1:\n",
    "            print(\"iteration\", iters)\n",
    "            print(\"critic loss:\", critic_loss[-1])\n",
    "            print(\"generator loss:\", generator_loss[-1])\n",
    "\n",
    "            for i in range(4):\n",
    "                plt.subplot(2,2,i+1)\n",
    "                footprint=generated_maps[i,...,0]\n",
    "                xd, yd = rectangular_array()\n",
    "                mask = footprint != 0\n",
    "                mask[5, 5] = True\n",
    "                marker_size = 50 * footprint[mask]\n",
    "                plot = plt.scatter(xd, yd, c='grey', s=10, alpha=0.3, label=\"silent\")\n",
    "                circles = plt.scatter(xd[mask], yd[mask], c=footprint[mask],\n",
    "                                       s=marker_size, alpha=1, label=\"loud\")\n",
    "                cbar = plt.colorbar(circles)\n",
    "                cbar.set_label('signal [a.u.]')\n",
    "                plt.grid(True)\n",
    "\n",
    "            plt.suptitle(\"iteration %i\" % iters)\n",
    "            plt.tight_layout(rect=[0, 0.03, 1, 0.95])\n",
    "            plt.savefig(\"./fake_showers_iteration_%.6i.png\" % iters)\n",
    "            plt.close(\"all\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "WU732pJQrNPp"
   },
   "source": [
    "### Results\n",
    "Let us plot the critic loss.\n",
    "We now expect to see a convergence of the total loss if the training was successful."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 338
    },
    "id": "7eFNdrXQrNPp",
    "outputId": "7b0a6504-45b1-40ac-a72a-8ae5e47165eb"
   },
   "outputs": [],
   "source": [
    "critic_loss = np.array(critic_loss)\n",
    "\n",
    "plt.subplots(1, figsize=(10, 5))\n",
    "plt.plot(np.arange(len(critic_loss)), critic_loss[:, 0], color='red', markersize=12, label=r'Total')\n",
    "plt.plot(np.arange(len(critic_loss)), critic_loss[:, 1] + critic_loss[:, 2], color='green', label=r'Wasserstein', linestyle='dashed')\n",
    "plt.plot(np.arange(len(critic_loss)), critic_loss[:, 3], color='royalblue', markersize=12, label=r'GradientPenalty', linestyle='dashed')\n",
    "plt.legend(loc='upper right')\n",
    "plt.xlabel(r'Iterations')\n",
    "plt.ylabel(r'Loss')\n",
    "plt.ylim(-6, 3)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ngwT2RTerNPp"
   },
   "source": [
    "In addition, we can visualize the generator loss."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 334
    },
    "id": "U17GdD5RrNPp",
    "outputId": "fcbbeb8a-ba42-42ee-9b5c-2fc7981f671d"
   },
   "outputs": [],
   "source": [
    "generator_loss = np.array(generator_loss)\n",
    "\n",
    "plt.subplots(1, figsize=(10, 5))\n",
    "plt.plot(np.arange(len(generator_loss)), generator_loss, color='red', markersize=12, label=r'Total')\n",
    "plt.legend(loc='upper right')\n",
    "plt.xlabel(r'Iterations')\n",
    "plt.ylabel(r'Loss')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "Ldj9t83XrNPq"
   },
   "source": [
    "#### Plot generated footprints\n",
    "To nicely see the training progress, we can create a gif of the generated samples during the training."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 305
    },
    "id": "Y_CVBVfFrNPq",
    "outputId": "43c7a5f5-d21b-4995-be1e-526735ce8f86"
   },
   "outputs": [],
   "source": [
    "import imageio\n",
    "import glob\n",
    "\n",
    "out_file = 'generated_shower_samples.gif'\n",
    "\n",
    "with imageio.get_writer(out_file, mode='I', duration=0.5) as writer:\n",
    "    file_names = glob.glob('fake_showers_iteration_*.png')\n",
    "    file_names = sorted(file_names)\n",
    "    last = -1\n",
    "\n",
    "    for i, file_name in enumerate(file_names):\n",
    "        animated_image = imageio.imread(file_name)\n",
    "        writer.append_data(animated_image)\n",
    "\n",
    "    animated_image = imageio.imread(file_name)\n",
    "    writer.append_data(animated_image)\n",
    "\n",
    "from IPython.display import Image\n",
    "Image(open('generated_shower_samples.gif','rb').read())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Name four general challenges of training generative adversarial frameworks."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####  Explain why approximating the wasserstein distance in the discriminator/critic helps to reduce mode collapsing."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": []
  }
 ],
 "metadata": {
  "accelerator": "GPU",
  "colab": {
   "name": "WGAN.ipynb",
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}