{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Implementing Filtered Back Projection"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We'll use the classes of PyTomography to implement filtered back projection in SPECT."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "from pytomography.projections import SPECTSystemMatrix\n",
    "from pytomography.metadata import ObjectMeta, ImageMeta\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import torch"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The two foundational tools of image reconstruction are\n",
    "\n",
    "1. Forward projection $\\sum_{i} c_{ij} a_i$\n",
    "2. Back projection $\\sum_{j} c_{ij} b_j$\n",
    "\n",
    "Let's discuss what these operators actually mean. First, let's define our quantities. $c_{ij}$ is known as the system matrix, and may include information involving attenuation and PSF correction. $a_i$ is an arbtriary object and $b_j$ is an arbitrary image. \n",
    "\n",
    "It's worth now discussing what the indices $i$ and $j$ actually mean. You might think: objects are three dimensional, shouldn't there be at least 3 indices when we're doing linear operations? Consider the following: because we are in a discrete space, any 3 dimensional object can be converted to a single (albeit very long) one dimensional object: a 128x128x128 3D matrix can be converted into a single 1D vector of length 2097152. That's how many voxels there are in object space: you can think of index $i$ as indexing a single voxel.\n",
    "\n",
    "The same can be said for an image. If we have 64 projections of matrix size 128x128, then that can be thought of as a single vector of length 1048576. That's also how many individual detector elements there are.\n",
    "\n",
    "So in forward projection $\\sum_{i} c_{ij} a_i$, the system matrix $c_{ij}$ maps the contribution from voxel $i$ to a detector element $j$. In back projection $\\sum_{j} c_{ij} b_j$, the system matrix $c_{ij}$ maps the intensity in detector element back to every possible voxel $i$ that could have contributed to it. In reality, however, not every voxel that could have contributed to detector element $j$ does so with equal intensity; it is for this reason that forward projection followed by back projection does not yield the original image. \n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's experiment with these operators. First we'll make a 3D rectangle in object space:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "torch.Size([1, 128, 128, 132])"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = torch.linspace(-1,1,128)\n",
    "y = torch.linspace(-1,1,128)\n",
    "z = torch.linspace(-1,1,132)\n",
    "xv, yv, zv = torch.meshgrid([x,y,z], indexing='ij')\n",
    "object_truth = (xv>-0.2)*(xv<0.2)*(yv>-0.15)*(yv<0.15)*(zv>-0.1)*(zv<0.1)\n",
    "object_truth = object_truth.to(torch.float).unsqueeze(dim=0) # add batch dimension\n",
    "object_truth.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f5955a0d070>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAFWCAYAAAB6nKoDAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy88F64QAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAVWklEQVR4nO3df2zcdeHH8ddd2d2xbNfN1F235pKKv5CMrdrSy4FkITlp4lKzP4wVFto0MgWnjl2MG7D1psQViS6NodiwUPUPl1VNICZrSvBkGuWSapdGVAqZA9uQ3K0N4w46t4O7j3/w3fk91rH35/p27D48H8nnj73z+XX/9MXr/X7f4XMcxxEAAJfhf79fAABQHwgMAIARAgMAYITAAAAYITAAAEYIDACAEQIDAGCEwAAAGCEwAABGCAwAgBECAwDqzB/+8Ad1d3drw4YN8vl8euqppy57zfHjx/WZz3xGwWBQH/vYx/Szn/3M9XMJDACoM4uLi9q8ebOGh4eNzn/55Ze1detW3XbbbZqentZ9992nu+++W08//bSr5/r48UEAqF8+n09PPvmktm3bdslz9uzZo2PHjulvf/tbZezLX/6yXn/9dU1MTBg/i4YBAB6XyWSUSCSqxrq6upTJZFzd5xqbLwUAH1Tnzp1TsVis+XrHceTz+arGgsGggsHgcl9N2WxWkUikaiwSiahQKOjf//63rr32WqP7GAfGuz8IANQ7WzPy586d00c+8hFls9ma77Fq1Sq9+eabVWOpVEoHDhxY5tvZQ8MAgGUqFovKZrOam5tTOBx2fX2hUFA0Gr3oehvtQpKam5uVy+WqxnK5nMLhsHG7kAgMALBm9erVWr16tevrLjSdcDhcU+BcTjwe1/j4eNXYM888o3g87uo+LHoDQJ158803NT09renpaUnvbJudnp7W7OysJOn+++9Xb29v5fx77rlHp06d0ne+8x3NzMzoscce0y9/+Uvt3r3b3YMdQ5I4ODg4PHXYks/nHUnOmTNnnFKp5Po4c+aMI8nJ5/NGz3v22WeX/Dx9fX2O4zhOX1+fs2XLlouuaWtrcwKBgHPdddc5P/3pT11/TuPvYbDoDcBrDP/8XVahUFBjY6POnDlT8xrG2rVrlc/n/ydTUrYwJQUAMMKiNwBY4jhOTa3FVtP5X6NhAACM0DAAwBIaBgAAIjAAAIaYkgIAS5iSAgBANAwAsIaGAQCAaBgAYA0NAwAAERgAAENMSQGAJUxJAQAgAgMAYIgpKQCwhCkpAABEYAAADDElBQCWMCUFAIBoGABgDQ0DAAARGAAAQ0xJAYAlTEkBACAaBgBYQ8MAAEA0DACwhoYBAIAIDACAIaakAMASpqQAABANAwCs8XrDIDAAwBKvBwZTUgAAIwQGAMAIU1IAYInXp6QIDACwqF7++NeCKSkAgBEaBgBY4vUpKRoGAMAIDQMALPF6wyAwAMASrwcGU1IAACM0DACwxOsNg8AAAEu8HhhMSQEAjBAYAAAjTEkBgCVen5IiMADAEq8HBlNSAGDJhcCo5XBreHhYra2tCoVCisVimpycfM/zh4aG9MlPflLXXnutotGodu/erXPnzrl6JoEBAHVmbGxMyWRSqVRKJ06c0ObNm9XV1aXTp08vef6RI0e0d+9epVIpvfDCC3riiSc0NjamBx54wNVzCQwAqDOHDh3Sjh071N/frxtuuEEjIyNauXKlRkdHlzz/ueee0y233KI777xTra2tuv3223XHHXdctpW8G4EBAJYsd0qqUChUHefPn7/oGcViUVNTU0okEpUxv9+vRCKhTCaz5HvdfPPNmpqaqgTEqVOnND4+rs9//vOuPh+L3gBgyXIXvaPRaNV4KpXSgQMHqsYWFhZUKpUUiUSqxiORiGZmZpa8/5133qmFhQV99rOfleM4evvtt3XPPfe4npIiMADAkuUGxtzcnMLhcGU8GAxaea/jx4/r4MGDeuyxxxSLxXTy5Ent2rVLDz30kPbv3298HwIDAK4S4XC4KjCW0tTUpIaGBuVyuarxXC6n5ubmJa/Zv3+/7rrrLt19992SpBtvvFGLi4v66le/qgcffFB+v9nqBGsYAGDJldhWGwgE1N7ernQ6XRkrl8tKp9OKx+NLXnP27NmLQqGhoaHyzqZoGABQZ5LJpPr6+tTR0aHOzk4NDQ1pcXFR/f39kqTe3l61tLRocHBQktTd3a1Dhw7p05/+dGVKav/+/eru7q4EhwkCAwAsuVLf9O7p6dH8/LwGBgaUzWbV1tamiYmJykL47OxsVaPYt2+ffD6f9u3bp1dffVUf/vCH1d3dre9///uunutzDN/U5/O5ujEAXO1s/SRHoVBQY2Oj/vrXv2r16tWur3/jjTe0adMm5fP5y65hvJ9oGABgidd/S4rAAABLCAwAgBGvBwbbagEARggMAIARpqQAwBKvT0kRGABgUb388a8FU1IAACMEBgDACFNSAGAJaxgAACMEBgDACIEBADBCYAAAjHg9MNglBQAwQsMAAEtoGAAAiIYBANZ4vWEQGABgCYEBADDi9cBgDQMAYISGAQCWeL1hEBgAYAmBAQAwQmAAAIwQGAAAI14PDHZJAQCMEBgAACNMSQGAJV6fkiIwAMASAgMAYITAAAAYq5c//rVg0RsAYISGAQCWMCUFADBCYAAAjHg9MFjDAAAYoWEAgCU0DAAARMMAAGu83jAIDACwhMAAABghMAAARggMAIARrwcGu6QAAEZoGABgCQ0DAADRMADAGq83DAIDACwpl8sql8s1XVcPCAwAsMTrDYM1DACoQ8PDw2ptbVUoFFIsFtPk5OR7nv/6669r586dWr9+vYLBoD7xiU9ofHzc1TNpGABgyZVqGGNjY0omkxoZGVEsFtPQ0JC6urr04osvat26dRedXywW9bnPfU7r1q3Tr3/9a7W0tOhf//qX1qxZ4+q5BAYAWHKlAuPQoUPasWOH+vv7JUkjIyM6duyYRkdHtXfv3ovOHx0d1WuvvabnnntOK1askCS1tra6fk+mpADAEsdxKgvfbg43gVEsFjU1NaVEIlEZ8/v9SiQSymQyS17zm9/8RvF4XDt37lQkEtHGjRt18OBBlUolV5+PhgEAliy3YRQKharxYDCoYDBYNbawsKBSqaRIJFI1HolENDMzs+T9T506pd/97nfavn27xsfHdfLkSX3961/XW2+9pVQqZfyeNAwAsORCYNRySFI0GlVjY2PlGBwctPJe5XJZ69at0+OPP6729nb19PTowQcf1MjIiKv70DAA4CoxNzencDhc+fe724UkNTU1qaGhQblcrmo8l8upubl5yfuuX79eK1asUENDQ2XsU5/6lLLZrIrFogKBgNH70TAAwJLlNoxwOFx1LBUYgUBA7e3tSqfTlbFyuax0Oq14PL7ke91yyy06efJk1RcEX3rpJa1fv944LCQCAwCsWW5gmEomkzp8+LB+/vOf64UXXtC9996rxcXFyq6p3t5e3X///ZXz7733Xr322mvatWuXXnrpJR07dkwHDx7Uzp07XT2XKSkAsORK/TRIT0+P5ufnNTAwoGw2q7a2Nk1MTFQWwmdnZ+X3/7cPRKNRPf3009q9e7c2bdqklpYW7dq1S3v27HH1XJ9jGG0+n8/VjQHgamfrJzkKhYIaGxv1i1/8QitXrnR9/dmzZ7V9+3bl8/mqNYyrDQ0DS6qX37axjf8wwnLwW1IAAIiGAQDWeL1hEBgAYAn/PwwAgBEaBgDACIEBADDi9cBglxQAwAgNAwAs8XrDIDAAwBJ2SQEAjNAwAADG6uWPfy0IDACwxOsNg11SAAAjNAwAsIRFbwCAEa9PSREYAGAJgQEAMEJgAACMeD0w2CUFADBCwwAAS9glBQAw4vUpKQIDACwhMAAARrweGCx6AwCM0DAAwBIWvQEARrw+JUVgAIAlBAYAwAiBAQAwVi9//GvBLikAgBEaBgBYwpQUAMAIgQEAMEJgAACMEBgAACNeDwx2SQEAjNAwAMASfksKAGDE61NSBAYAWEJgAACMEBgAACNeDwx2SQEAjNAwAMASrzcMAgMALCEwAABGCAwAgBECAwBgxOuBwS4pAIARGgaW5PP53u9XAOqO4zg1/S5UvTQMAgMALGFKCgBg5EJg1HK4NTw8rNbWVoVCIcViMU1OThpdd/ToUfl8Pm3bts31MwkMALDkSgXG2NiYksmkUqmUTpw4oc2bN6urq0unT59+z+teeeUVffvb39att95a0+cjMADAkisVGIcOHdKOHTvU39+vG264QSMjI1q5cqVGR0cveU2pVNL27dv13e9+V9ddd11Nn4/AAICrRKFQqDrOnz9/0TnFYlFTU1NKJBKVMb/fr0QioUwmc8l7f+9739O6dev0la98peb3IzAAwJLlNoxoNKrGxsbKMTg4eNEzFhYWVCqVFIlEqsYjkYiy2eyS7/XHP/5RTzzxhA4fPrysz8cuKQCwZLm7pObm5hQOhyvjwWBw2e/0xhtv6K677tLhw4fV1NS0rHsRGABgyXIDIxwOVwXGUpqamtTQ0KBcLlc1nsvl1NzcfNH5//znP/XKK6+ou7u7MnbhuyLXXHONXnzxRX30ox81ek+mpADAkiux6B0IBNTe3q50Ol0ZK5fLSqfTisfjF51//fXX6/nnn9f09HTl+MIXvqDbbrtN09PTikajxs+mYQCAJVfqi3vJZFJ9fX3q6OhQZ2enhoaGtLi4qP7+fklSb2+vWlpaNDg4qFAopI0bN1Zdv2bNGkm6aPxyCAwAqDM9PT2an5/XwMCAstms2traNDExUVkIn52dld9vfwLJ5xhGG78tBMBrbP0kR6FQUGNjo771rW/VtFB9/vx5/fjHP1Y+n7/sGsb7iYYBAJaUy+WafnywlmveDwQGAFji9R8fJDAAwBICAwBgxOuBwfcwAABGaBgAYInXGwaBAQCWEBgAAGP18se/FgQGAFhCwwAAGPF6YLBLCgBghIYBAJZ4vWEQGABgCYEBADBCYAAAjBAYAAAjXg8MdkkBAIzQMADAEq83DAIDACwhMAAARggMAIARAgMAYMTrgcEuKQCAERoGAFji9YZBYACAJQQGAMAIgQEAMEJgAACMeD0w2CUFADBCwwAAi+qlLdSCwAAAS7w+JUVgAIAl5XJZ5XK5puvqAYEBAJbQMAAARrweGOySAgAYoWEAgCVebxgEBgBYQmAAAIwQGAAAI14PDBa9AQBGaBgAYInXGwaBAQCWEBgAACP8NAgAwAgNAwBgxOsNg11SAAAjNAwAsIQpKQCAEcdxappeIjAA4AOGhgEAMOL1wGDRGwAsubBLqpbDreHhYbW2tioUCikWi2lycvKS5x4+fFi33nqr1q5dq7Vr1yqRSLzn+ZdCYABAnRkbG1MymVQqldKJEye0efNmdXV16fTp00uef/z4cd1xxx169tlnlclkFI1Gdfvtt+vVV1919VyfY9iFfD6fqxsDwNXO1lRQoVBQY2OjtmzZomuucT/T//bbb+v3v/+98vm8wuHwZc+PxWK66aab9Oijj0p6p9lEo1F985vf1N69ey97falU0tq1a/Xoo4+qt7fX+D1pGABgyZWYkioWi5qamlIikaiM+f1+JRIJZTIZo3ucPXtWb731lj70oQ+5+nwsegOAJctd9C4UClXjwWBQwWCwamxhYUGlUkmRSKRqPBKJaGZmxuh5e/bs0YYNG6pCxwQNAwAsuRAYtRySFI1G1djYWDkGBwetv+PDDz+so0eP6sknn1QoFHJ1LQ0DACxZ7m9Jzc3NVa1hvLtdSFJTU5MaGhqUy+WqxnO5nJqbm9/zOT/84Q/18MMP67e//a02bdrk+j1pGABwlQiHw1XHUoERCATU3t6udDpdGSuXy0qn04rH45e89yOPPKKHHnpIExMT6ujoqOn9aBgAYMmV+uJeMplUX1+fOjo61NnZqaGhIS0uLqq/v1+S1Nvbq5aWlsqU1g9+8AMNDAzoyJEjam1tVTablSStWrVKq1atMn4ugQEAllypnzfv6enR/Py8BgYGlM1m1dbWpomJicpC+OzsrPz+/04g/eQnP1GxWNQXv/jFqvukUikdOHDA+Ll8DwPAB5bt72HcdNNNNX8P489//rPx9zDeLzQMALDE678lRWAAgCVe/3lzdkkBAIzQMADAEqakAABGyuVyTRuEapnGej8QGABgCQ0DAGCEhgEAMOL1hsEuKQCAERoGAFji9YZBYACAJaxhAACM0DAAAEZoGAAAI15vGOySAgAYoWEAgCVebxgEBgBYwhoGAMBYvbSFWhAYAGAJU1IAACNeDwx2SQEAjNAwAMASrzcMAgMALKl1txO7pADgA4aGAQAwQmAAAIx4PTDYJQUAMELDAABLvN4wCAwAsITAAAAYYVstAMAIDQMAYMTrgcEuKQCAERoGAFji9YZBYACAJQQGAMAIgQEAMOI4Tk1bZAkMAPiAqfUPf70EBrukAABGaBgAYInXGwaBAQCWEBgAACMEBgDACIEBADBSLpfl8/lcX1cvgcEuKQCAERoGAFjClBQAwAiBAQAwQmAAAIwQGAAAI14PDHZJAQCMEBgAYEm5XK75cGt4eFitra0KhUKKxWKanJx8z/N/9atf6frrr1coFNKNN96o8fFx9x/QMSSJg4ODw1OHLfl83pHk+Hw+x+/3uz58Pp8jycnn80bPO3r0qBMIBJzR0VHn73//u7Njxw5nzZo1Ti6XW/L8P/3pT05DQ4PzyCOPOP/4xz+cffv2OStWrHCef/55V5+TwODg4PjAHrZcCAzpndBwe1y41jQwOjs7nZ07d1b+XSqVnA0bNjiDg4NLnv+lL33J2bp1a9VYLBZzvva1r7n6nExJAYBFzv/9b1rdHG4Ui0VNTU0pkUhUxvx+vxKJhDKZzJLXZDKZqvMlqaur65LnXwq7pADgKlEoFKr+HQwGFQwGq8YWFhZUKpUUiUSqxiORiGZmZpa8bzabXfL8bDbr6v2MA8NtCgLAB0UgEFBzc7PrP8D/36pVqxSNRqvGUqmUDhw4sMy3s4eGAQDLFAqF9PLLL6tYLNZ8D8dxLvql23e3C0lqampSQ0ODcrlc1Xgul1Nzc/OS925ubnZ1/qUQGABgQSgUUigU+p8/JxAIqL29Xel0Wtu2bZP0znbedDqtb3zjG0teE4/HlU6ndd9991XGnnnmGcXjcVfPJjAAoM4kk0n19fWpo6NDnZ2dGhoa0uLiovr7+yVJvb29amlp0eDgoCRp165d2rJli370ox9p69atOnr0qP7yl7/o8ccfd/VcAgMA6kxPT4/m5+c1MDCgbDartrY2TUxMVBa2Z2dn5ff/dxPszTffrCNHjmjfvn164IEH9PGPf1xPPfWUNm7c6Oq5PofVbACAAb6HAQAwQmAAAIwQGAAAIwQGAMAIgQEAMEJgAACMEBgAACMEBgDACIEBADBCYAAAjBAYAAAjBAYAwMh/APRd6SAAzfvYAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 500x400 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(5,4))\n",
    "plt.pcolormesh(object_truth[0][:,:,64].T, cmap='Greys_r')\n",
    "plt.axis('off')\n",
    "plt.colorbar()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Before we do any projections, we need to get corresponding metadata for our object. In this case, we'll assume the voxel sizes are 1cm $^3$. For our image space, we'll assume 60 projections are taken at angular spacing of 6 degrees."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "angles = np.arange(0,360.,6.)\n",
    "object_meta = ObjectMeta(dr=(1,1,1), shape=object_truth[0].shape)\n",
    "image_meta = ImageMeta(object_meta, angles=angles)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With this metadata, we can create our forward and back projection networks. We'll model no phenomenon for now."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "system_matrix = SPECTSystemMatrix(\n",
    "    obj2obj_transforms=[],\n",
    "    im2im_transforms=[],\n",
    "    object_meta=object_meta,\n",
    "    image_meta=image_meta)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now use the `forward` method of `system_matrix` to model $g=Hf$ (convert the object $f$ into an image $g$)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "torch.Size([1, 60, 128, 132])"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "image = system_matrix.forward(object_truth)\n",
    "image.shape"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can look at a projection at 60 degrees for example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABJ4AAAFeCAYAAADEw9X5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy88F64QAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAfNElEQVR4nO3de5CVdf3A8c/hvrsgN5fLAkLeEzA0wUsqkIg5eClFCaM0mWi8lDmliVr6U1FTy5lUNNAYLxglYWppSmbljI2KqOHkXVG5CaSCoCmX7+8PZk+se2Fl+crCeb1mGNfnfJ9znrPDhz3nvc85p5BSSgEAAAAAW1iLrX0AAAAAAGyfhCcAAAAAshCeAAAAAMhCeAIAAAAgC+EJAAAAgCyEJwAAAACyEJ4AAAAAyEJ4AgAAACAL4QkAAACALISnZqZQKMTFF1+8tQ8D2ASzCtsu8wvbNjMM2wazSjXhqQGTJ0+OQqEQ+++//9Y+lKzWr18fV111VXzuc5+Ldu3axd577x2/+c1vGr3/e++9FxMmTIjKysqoqKiI4cOHx9y5czMeMdRUCrO6aNGiGDduXOyxxx7RoUOH6NSpUwwZMiRuvfXWSCnVWr9w4cI48cQTo1OnTrHDDjvEscceG6+99lqjb++xxx6Lgw8+OMrLy6NHjx7x/e9/P1atWrUl7xJERGnMb7VXX301TjrppOjWrVuUlZXFbrvtFhdccEGtdc8//3x85Stfifbt20eXLl3im9/8ZixbtqzRt3PvvffGvvvuG+3atYuddtopLrrooli7du2WvCtQVCoz/Morr8To0aOjc+fOUV5eHgcffHA88sgjda41wzRHpTKrkyZNimOOOSa6d+/eYPiaNWtWjBkzJnbeeecoLy+PPfbYI374wx/Ge++9V+f6psxlU59vbw9abe0DaM6mT58e/fr1iyeeeCJeeeWV2HXXXbf2IWVxwQUXxJVXXhnf+c53YvDgwXHPPffESSedFIVCIb7+9a83uO/69etj1KhR8eyzz8Y555wTO+64Y0yePDmGDRsWTz31VOy2226f0b2glJXCrC5fvjwWLFgQo0ePjp122inWrFkTs2fPjlNOOSVefPHFuPzyy4trV61aFcOHD48VK1bE+eefH61bt45rr702hg4dGs8880x07dq1wdt65pln4rDDDovPf/7z8Ytf/CIWLFgQ11xzTbz88svxwAMP5L6rlJhSmN+IDXM1bNiw6NWrV/zwhz+Mrl27xptvvhlvvfVWjXULFiyIQw89NDp27BiXX355rFq1Kq655pqYN29ePPHEE9GmTZsGb+eBBx6Ir371qzFs2LC47rrrYt68eXHZZZfF0qVL48Ybb8x5FylRpTDDb731Vhx44IHRsmXLOOecc6KioiKmTZsWI0eOjIcffjgOPfTQ4lozTHNVCrMaEXHhhRdGjx49Yp999okHH3yw3nUTJkyIqqqqGDduXOy0004xb968uP766+P++++PuXPnRllZWXFtU+eyKc+3txuJOr322mspItKsWbNSZWVluvjiiz+T242IdNFFF30mt5VSSgsWLEitW7dOZ5xxRnHb+vXr0yGHHJJ69+6d1q5d2+D+v/3tb1NEpLvuuqu4benSpalTp05p7Nix2Y4bqpXKrNbnqKOOShUVFTVm9Wc/+1mKiPTEE08Utz3//POpZcuWaeLEiZu8ziOPPDL17NkzrVixorht6tSpKSLSgw8+uGXvACWtVOZ33bp1acCAAWn//fdPH3zwQYNrTzvttFRWVpbeeOON4rbZs2eniEi/+tWvNnlbe+21V/rCF76Q1qxZU9x2wQUXpEKhkJ5//vnNvxNQh1KZ4dNPPz21atUqvfDCC8Vtq1evTn369En77rtvjbVmmOaoVGY1pZRef/31lFJKy5Yta/D2H3nkkVrbbr311hQRaerUqTW2N2Uum/p8e3vhpXb1mD59enTu3DlGjRoVo0ePjunTp9e4fP78+VEoFOKaa66JKVOmxC677BJt27aNwYMHx5NPPlnr+u66667Ya6+9ol27djFgwIC4++6745RTTol+/fpt8lgWLlwYp556anTv3j3atm0b/fv3j1//+tdb5H7ec889sWbNmjj99NOL2wqFQpx22mmxYMGC+Oc//9ng/jNnzozu3bvHcccdV9xWWVkZJ554Ytxzzz3x0UcfbZHjhPqUyqzWp1+/fvHBBx/Exx9/XNw2c+bMGDx4cAwePLi4bc8994zDDjssfve73zV4fStXrozZs2fHuHHjYocddihu/9a3vhXt27ff5P7waZTK/D700EPx3HPPxUUXXRRlZWXxwQcfxLp16+pc+/vf/z6OOuqo2GmnnYrbRowYEbvvvvsm5+/f//53/Pvf/44JEyZEq1b/O6n99NNPj5RSzJw5c4vcH6hWKjP86KOPxj777BN77LFHcVt5eXkcc8wxMXfu3Hj55ZeL280wzVGpzGpENOoYIiKGDRtWa9vXvva1iNjwctlqTZ3Lpj7f3l54qV09pk+fHscdd1y0adMmxo4dGzfeeGM8+eSTNZ7IRUTceeed8f7778d3v/vdKBQKcdVVV8Vxxx0Xr732WrRu3ToiIv70pz/FmDFjYuDAgXHFFVfEu+++G+PHj49evXpt8jjefvvtOOCAA6JQKMSZZ54ZlZWV8cADD8T48eNj5cqV8YMf/KC4dvny5Y26bx06dIi2bdtGRMTTTz8dFRUV8fnPf77GmiFDhhQvP/jgg+u9rqeffjr23XffaNGiZsMcMmRITJkyJV566aUYOHBgo44LNkepzGq1Dz/8MFavXh2rVq2Kv//97zFt2rQ48MADi6cDr1+/Pv71r3/FqaeeWuv6hgwZEg899FC8//770aFDhzpvc968ebF27drYb7/9amxv06ZNDBo0KJ5++ulGHTs0RqnM71/+8peIiGjbtm3st99+8dRTT0WbNm3ia1/7WkyePDm6dOkSERsejC9durTW/EVsmN/777+/wdusns9P7l9VVRW9e/c2v2xxpTLDH330UXTu3LnWmvLy8oiI4ttLmGGaq1KZ1aZasmRJRETsuOOOxW1NncumPt/ebmzlM66apTlz5qSISLNnz04pbTgVrnfv3umss84qrnn99ddTRKSuXbumd955p7j9nnvuSRGR7rvvvuK2gQMHpt69e6f333+/uO1vf/tbiojUt2/fGrcdnzgdcPz48alnz55p+fLlNdZ9/etfTx07dqxxyn5ENOrPtGnTivuMGjUq7bzzzrW+B6tXr04Rkc4777wGv1cVFRXp1FNPrbX9T3/6U4qI9Oc//7nB/aEpSmlWq11xxRU11hx22GHpzTffLF5efVrxJZdcUmvfG264IUVEjZcKfNJdd92VIiL94x//qHXZCSeckHr06FHvvvBplNL8HnPMMcX78Y1vfCPNnDkz/eQnP0mtWrVKBx10UFq/fn1KKaUnn3wyRUS67bbban2/zjnnnBQR6b///W+939Orr746RUSNfxOqDR48OB1wwAH17gufVinN8NFHH506deqUVq5cWeP6DzzwwBQR6ZprrkkpmWGap1Ka1Y1t6qV2dRk/fnxq2bJleumll4rbmjqXTX2+vb1wxlMdpk+fHt27d4/hw4dHxIZT4caMGRN33HFH/PznP4+WLVsW144ZM6bGb0AOOeSQiIjip0ctWrQo5s2bF+eff360b9++uG7o0KExcODAWLlyZb3HkVKK3//+93HiiSdGSqlG9T3iiCNixowZMXfu3PjSl74UERGzZ89u1P3r379/8esPP/ywzkLcrl274uUNaer+0BSlNKvVxo4dG/vtt18sW7Ys/vjHP8bbb79dY86qv97cudzU/maaLaWU5rf6EyEHDx4cd9xxR0REHH/88VFeXh4TJ06Mhx9+OEaMGNHo+a3vN7ub2r+h7wN8WqU0w6eddlrcd999MWbMmJg0aVJUVFTE5MmTY86cORHxv9kzwzRHpTSrTXHnnXfGLbfcEueee26ND8hq6lx6vryB8PQJ69atixkzZsTw4cPj9ddfL27ff//94+c//3k8/PDDMXLkyOL2jV+/HRHFQX333XcjIuKNN96IiKjzUwN23XXXmDt3br3HsmzZsnjvvfdiypQpMWXKlDrXLF26tPj1iBEjNnX3aikrK6vzfZj++9//Fi/PuT9srlKb1Wp9+/aNvn37RsSGCDVhwoQYMWJEvPjii1FWVlacuc2dy03tb6bZEkptfqvnZuzYsTW2n3TSSTFx4sR47LHHYsSIEeaXbUapzfCRRx4Z1113XZx33nmx7777Fo9r0qRJce655xafgJthmptSm9XN9eijj8b48ePjiCOOiEmTJtW4rKlz6fnyBsLTJ/z1r3+NxYsXx4wZM2LGjBm1Lp8+fXqN4dy4EG8spdTkY1m/fn1ERIwbNy5OPvnkOtfsvffexa+rX5O6KR07diz+Be/Zs2c88sgjkVKKQqFQXLN48eKI2PDa1Yb07NmzuHZjjd0fNlepzWp9Ro8eHVOnTo1//OMfccQRR0SXLl2ibdu2mz2XPXv2rLH2k/ubabaEUpvf6rnp3r17jTXdunWLiP89oN/U/FXPd3023r9Pnz619q9+PwloqlKb4YiIM888M7797W/Hv/71r+L7Ht5yyy0REbH77rtHhBmm+SnFWf20nn322TjmmGNiwIABMXPmzBpvIB7R9Lls6vPt7YXw9AnTp0+Pbt26xQ033FDrslmzZsXdd98dN910U6Ovr/rMhFdeeaXWZXVt21hlZWV06NAh1q1b16jiWz0UmzJt2rQ45ZRTIiJi0KBBcfPNN8fzzz8fe+21V3HN448/Xry8IYMGDYpHH3001q9fX+MNxh9//PEoLy8v/iCGLa3UZrU+1afnrlixIiIiWrRoEQMHDiye/r+xxx9/PHbeeed631g8ImLAgAHRqlWrmDNnTpx44onF7R9//HE888wzNbbB5iq1+f3iF78YU6dOjYULF9ZYs2jRouIxRET06tUrKisr65zfJ554olE/kyMi5syZU+OB8KJFi2LBggUxYcKERh07bEqpzXC1ioqKOPDAA4v//5e//CXKysqKLw0ywzQ3pTqrjfXqq6/GV77ylejWrVvcf//9NV4+WK2pc9nU59vbC+FpIx9++GHMmjUrTjjhhBg9enSty6uqquI3v/lN3HvvvbH//vs36jqrqqpiwIABcdttt8XEiROLf5n//ve/x7x584rDW5eWLVvG8ccfH3feeWc899xzMWDAgBqXL1u2rPhgNWLzXgd77LHHxtlnnx2TJ0+O66+/PiI2FO2bbropevXqFQcddFBx7eLFi2PFihWxyy67FD/VYPTo0TFz5syYNWtW8Xu2fPnyuOuuu+Loo4/eYp8wABsrxVn95HVUu+WWW6JQKBRP/Y/YMJfnnXdezJkzp/gJHC+++GL89a9/jR/96Ec19n/hhReivLy8eGp1x44dY8SIEXHHHXfET37yk2Kkuv3222PVqlVxwgknNOrYoT6lOL/HHntsnHXWWcUHx9W/qLn55psjIuLwww8vrj3++OPj1ltvjbfeeqv4m9WHH344XnrppTj77LOL69asWROvvvpqdOzYsfjgvH///rHnnnvGlClT4rvf/W7xN9c33nhjFAqFOr/f8GmV4gzX5bHHHotZs2bFaaedFh07dixuN8M0F2a1YUuWLImRI0dGixYt4sEHH6zzcXb19Td2LlesWBGLFy+Onj17Fv9d+DTPt7drn/GbmTdrM2bMSBGR/vCHP9R5+bp161JlZWU6+uiji+/8f/XVV9daF5949/x77703FQqFtPfee6drr702/fSnP01dunRJAwYMSP369Wtw3yVLlqS+ffum8vLydNZZZ6Vf/epX6YorrkgnnHBC6ty58xa539WfsjFhwoQ0derUNGrUqBQRafr06TXWnXzyySki0uuvv17ctnbt2nTAAQek9u3bp//7v/9LN9xwQ+rfv3/q0KFDg5+cBU1RirN61llnpf322y9deOGFacqUKenKK69MgwcPThGRvve979VYu3LlyrTLLrukbt26pauuuipde+21qU+fPqmqqiotXbq01v0YOnRojW1PPfVUatu2bdpnn33SjTfemC644ILUrl27NHLkyCbfDyjF+U0ppUsuuSRFRDr88MPTDTfckCZMmJAKhUIaO3ZsjXVvvvlm6tq1a9pll13SL3/5y3T55Zenzp07p4EDB9b4NKzq783JJ59cY//77rsvFQqF9OUvfzlNmTIlff/7308tWrRI3/nOd7bI/YBSnOH58+enIUOGpMsuuyzdfPPN6eyzz05lZWVpn332qfVJd2aY5qIUZzWllG677bZ06aWXpokTJ6aISMOHD0+XXnppuvTSS9P8+fOL677whS+kiEjnnntuuv3222v8eeihh2pcZ2Pnctq0aXV+0l5jn29vz4SnjRx99NGpXbt2afXq1fWuOeWUU1Lr1q2LH0vZmOFMacPg77nnnqlt27ZpwIAB6d57703HH3982nPPPTe579tvv53OOOOM1KdPn9S6devUo0ePdNhhh6UpU6Zs9n3d2Lp169Lll1+e+vbtm9q0aZP69++f7rjjjlrr6gpPKaX0zjvvpPHjx6euXbum8vLyNHTo0PTkk09ukWODupTirD700EPpqKOOSlVVVal169apQ4cO6Utf+lKaNm1a8aPYN/bWW2+l0aNHpx122CG1b98+HXXUUenll1+u83vwyfCUUkqPPvpoOuigg1K7du1SZWVlOuOMM2o9uIbNUYrzm9KGj6++7rrr0u67755at26d+vTpky688ML08ccf11r73HPPpZEjR6by8vLUqVOn9I1vfCMtWbKkxpr6nrSmlNLdd9+dBg0alNq2bZt69+5d7+3A5ijFGX7nnXfSsccem3r06JHatGmTPve5z6Uf//jH9f5cNMM0B6U4qymlNHTo0BQRdf555JFHahxbfX/qemzcmLmsLzw19vn29qyQ0hZ4pzA2y6BBg6KysrLRpxECW4dZhW2X+YVtmxmGbYNZpSEtNr2EplqzZk2sXbu2xra//e1v8eyzz8awYcO2zkEBtZhV2HaZX9i2mWHYNphVNocznj4D8+fPjxEjRsS4ceOiqqoqXnjhhbjpppuiY8eO8dxzz0XXrl239iECYVZhW2Z+YdtmhmHbYFbZHD7V7jPQuXPn+OIXvxg333xzLFu2LCoqKmLUqFFx5ZVXGkxoRswqbLvML2zbzDBsG8wqm8MZTwAAAABk4T2eAAAAAMhCeAIAAAAgC+EJAAAAgCwa/ebihUIh53HANm1beKs0Mwz1a+4zbH6hfs19fiPMMDSkuc+w+YWGNWaGnfEEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFq229gGUipTS1j6ETSoUClv7EKBZ6NatW+y8884REdGvX7/o1atXRETxvxtvq6ysjC5dukREREVFRbRqteX+WV27dm1ERKxevTreeeediIhYtmxZLFy4MCIi5s+fHxERCxcurLHttddei4iIpUuXbrFjgW3NDjvsED179oyIiKqqqqiqqoqIiK5du0bEhjnv1q1bRET06NEjdtxxx4iI6Ny5c3Tq1CkiItq3bx8REWVlZdGyZcs6b2fdunUREfHhhx9GRMSqVavivffei4iId999N5YvXx4REUuWLCnO5NKlS+M///lPREQsWrSo+N/FixdHRMTKlSubevdhu7Y1H1d7vAyb1tx/BkdE/Oc///Ez+DPkjCcAAAAAshCeAAAAAMhCeAIAAAAgC+EJAAAAgCyEJwAAAACyEJ4AAAAAyEJ4AgAAACAL4QkAAACALIQnAAAAALIQngAAAADIQngCAAAAIAvhCQAAAIAshCcAAAAAshCeAAAAAMhCeAIAAAAgC+EJAAAAgCyEJwAAAACyEJ4AAAAAyEJ4AgAAACAL4QkAAACALIQnAAAAALIQngAAAADIQngCAAAAIAvhCQAAAIAshCcAAAAAshCeAAAAAMhCeAIAAAAgC+EJAAAAgCyEJwAAAACyEJ4AAAAAyEJ4AgAAACCLQkopNWphoZD7WGCb1cgx2qrMMNSvuc+w+YX6Nff5jTDD0JDmPsPmFxrWmBl2xhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZCE8AAAAAZCE8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQhfAEAAAAQBbCEwAAAABZFFJKaWsfBAAAAADbH2c8AQAAAJCF8AQAAABAFsITAAAAAFkITwAAAABkITwBAAAAkIXwBAAAAEAWwhMAAAAAWQhPAAAAAGQhPAEAAACQxf8DfB1IZkmg40sAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 1500x400 with 5 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, axes = plt.subplots(1,5,figsize=(15,4))\n",
    "for i, proj in enumerate([0,5,10,15,20]):\n",
    "    axes[i].pcolormesh(image[0][proj].T, cmap='Greys_r')\n",
    "    axes[i].set_title(f'Angle={image_meta.angles[proj]}')\n",
    "    axes[i].axis('off')"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "At angles like 60 degrees, the cube is darkest in the center and lighter on the outside; this is like looking through a semi-transparent cube in real life, it's going to be darkest near the center if you're looking at it from an off angle."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also back project using the `backward` method of the `system_matrix` instance. Here we compute $H^T g$ and the normalization constant $H^T 1$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "tensor([[[[0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          ...,\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.]],\n",
       "\n",
       "         [[0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          ...,\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.]],\n",
       "\n",
       "         [[0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          ...,\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.]],\n",
       "\n",
       "         ...,\n",
       "\n",
       "         [[0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          ...,\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.]],\n",
       "\n",
       "         [[0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          ...,\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.]],\n",
       "\n",
       "         [[0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          ...,\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.],\n",
       "          [0., 0., 0.,  ..., 0., 0., 0.]]]])"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "object_bp, norm_bp = system_matrix.backward(image, return_norm_constant=True)\n",
    "object_bp"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But if we look at the new object:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "axial_slice = object_bp[0][:,:,64].T"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f5921114fa0>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZkAAAFICAYAAACCx4NEAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy88F64QAAAACXBIWXMAAA9hAAAPYQGoP6dpAABMwElEQVR4nO2dfZBX1X3/37vAAiK762rZdRtImdZRMD40kuDGhzF1B0RqNZJmqDuGRAZaw5qoHcVOhRhjpFKrCEmkZhofptCHTEda6YhuoZUkroBYiiIhZkqFRHfJDA8rKOzC8vvDud/c772f/ZzPuffc7y783q+ZHb6ce+7nnPvw/Z7zeTifU3Xy5MmTIIQQQgqgerA7QAgh5PSFgwwhhJDC4CBDCCGkMDjIEEIIKQwOMoQQQgqDgwwhhJDC4CBDCCGkMDjIEEIIKYzhg90BQgg5HTh69Ch6e3szn19TU4NRo0YF7NHQwDzIVFVVBWtUkhUvSx53JSWQjkcypLaqq4tR4FzXlUTrt3Z+/Lzoc7xedH1RmVTfVRa97FFZf39/qn68LNlf1z2W2tTQnmecqE+S3OR9kcqOHj2qtqmVafdKe7+l/lrejSwyQiA9d+15uvqtfd8leSGvK5Sso0ePYuLEiejq6soso6mpCbt37z7tBhpqMoQQkpPe3l50dXVh7969qK2t9T6/p6cH48ePR29vLweZwcA187NoC1ZNzDLLDEUI7dAiI14ny0zTR26e+xyiflKryNOOdm5R15n3vFBtWd/5ZD1freV0Y+zYsRg7dqz3eadzCkk6/gkhhBQGBxlCCCGFEcRcZnVIWlXlZP24EzlyOoY2a4UwjYTE2h/tmqP7FnfU+t43X7OQ5vjXAgWyoLUR0pzlMqFF/Thx4kTmdqQgjbzkubch+qO1L32nfYNArL8rlTRFnTx5MlN7NJcRQgghGcilyYScLVod+SFmr4OJbxBDCI0tfs+02bbUbtSWa+ap9bvoWZoWOi71UXqHfN+dPGHw2jPV+i3VH2ws9833XvlqNNZ6lbhn1GTSUJMhhBBSGBxkCCGEFEZhK/5Dqvhxh7GkeludzFkp2gznmymgqLUfvoTIpmCtb32eSSd8iHvlKrOsFfF1/Gd5TlmfbUhTTejnqVFUEEseaC5LQ02GEEJIYQRf8e87A9dGcG1W55IVIj9ZpVbka/XzaIRZw20HOleTkZzhu7QcSzizln/M1Zbk7JVyl1lzokntD9RvF76Of18Zvvg6xq3fLWsfLe0ORSe/BDWZNNRkCCGEFEZFFmMWca4U4jkYNtqirj2kRpjFR2AhxEJKq78oq2Y6GD4ZSV4In0wW7cZ6rkVm1u+Sb9uh300uxhxcqMkQQggpDA4yhBBCCqOiqf6zqtTxY5HDOI/zMUT9SuU6c12nb94nSW7WnFuaE9yaT8o3R5vVaW/tWxEr1qU2rc/RIsuFtp2D1eQVckuIOFnzlEkMRRMTzWVpqMkQQggpjIpoMpYZUOiFUr6hjCH6OBiLNqVjljBo1/Gs9y+PJqiFMPvKDdF/X8d/yGCNPGjX7uvIz5oZ2Vo/D4MdrixBTSYNNRlCCCGFwUGGEEJIYRRmLsuzgjrCmt/INzOAlazn5mkzqwPfeq/yZAHQZGgmOuvzCXG/Q74Lvmt5pDqVWt8TJ6vJKrSzP08wShLt/R6sNTHEBjUZQggJROSTyfLnw8aNG3HDDTegubkZVVVVWLNmzYB1/+zP/gxVVVVYtmxZWfn+/fvR1taG2tpa1NfXY+7cuTh8+HBZne3bt+Oqq67CqFGjMH78eCxdutSrn0DOQaaqqmrAmY5286qrq0t/audi9ZJ/RRFdU9FOfFf70XXG+6P95blH1jaS98VSJ17PVWbpo7Usa79d16IhPQPf5zcU3r0i28/6ndbquX6oB/OeFsWRI0dwySWX4Hvf+55a7/nnn8drr72G5ubm1LG2tjbs2LEDHR0dWLt2LTZu3Ij58+eXjvf09GDatGn45Cc/ia1bt+Kv//qv8cADD+Cpp57y6mtF18kQQsjpTKWiy2bMmIEZM2aodX71q1/hjjvuwEsvvYSZM2eWHdu5cyfWrVuHLVu2YMqUKQCAFStW4Prrr8ejjz6K5uZmrFq1Cr29vfjhD3+ImpoaXHjhhdi2bRsee+yxssHIxWlhLpNmKsmHnWe2nWfGq83cQmhqlhm79byQ/XDV0fph1UKy9kPCqoVobWrXlKcfWd+ZPO+v9fsQIf24no4aRNH09PSU/R07diyTnP7+ftx666245557cOGFF6aOd3Z2or6+vjTAAEBrayuqq6uxadOmUp2rr74aNTU1pTrTp0/Hrl27cODAAXNfTotBhhBCTgfGjx+Purq60t+SJUsyyXnkkUcwfPhwfP3rXxePd3V1Ydy4cWVlw4cPR0NDA7q6ukp1Ghsby+pE/4/qWKC5jBBCApHXXLZ3717U1taWykeOHOkta+vWrXjiiSfwxhtvDAlNsrDtl/MQoi2LDN92kuaBwSJPvyNc4arRZ2tbUb0QIdWWYwPVS7YvydCu09Wmdi8lGdZ7Xyl832HfZxAiDLroTAhDmdra2rJBJgs//vGPsW/fPkyYMKFUduLECfz5n/85li1bhv/7v/9DU1MT9u3bV3be8ePHsX//fjQ1NQEAmpqa0N3dXVYn+n9UxwLNZYQQEohKhTBr3Hrrrdi+fTu2bdtW+mtubsY999yDl156CQDQ0tKCgwcPYuvWraXzNmzYgP7+fkydOrVUZ+PGjejr6yvV6ejowPnnn4+zzjrL3J+Kbr9cxHmhZVjlhtAmQiDNIC2zcmt/tIWakgaRRftIyoiIa0PaTFmSG8mK13dt55zsVx7tJnksz4ZjGnk0B2vffAih7YfQtk73RZmHDx/GL37xi9L/d+/ejW3btqGhoQETJkzA2WefXVZ/xIgRaGpqwvnnnw8AmDRpEq677jrMmzcPK1euRF9fH9rb2zF79uxSuPMtt9yCb33rW5g7dy4WLlyIt956C0888QQef/xxr77SJ0MIIYGoVAjz66+/js9//vOl/999990AgDlz5uCZZ54xyVi1ahXa29tx7bXXorq6GrNmzcLy5ctLx+vq6vDyyy9jwYIFuOyyy3DOOedg8eLFXuHLAFB10nh1w4YNswk0zKK0vT7yEGLmLp2Tp2++M3zf/kSfJR+Itf/RKxB/FaKwxagsrslIqVui9yMqk+pLfYujaTJWWVE/sva7t7c3JTfL/YvQNKoQs21pr6WsMkL/OIa+vgjLdVrbtqZSctHT04O6ujr88pe/zORT6enpwSc+8QkcOnQot09mqBFEk/H9AY2/JJa1GHl+oEMPPBFSv7N+ySV8TXSaiSGLOUFaI2Ehaa4C5PtimWhY73G8nsUBHSLAIY7lfluPhzTz+L6j1sFUOieP2c5yzb7fraESpENoLiOEkGBwP5k03oNM6Nl/hK/j2kURGowlv5KLEAEDvmVZtK7kDFXarlkL1ZWCArR2XMekfmgypHfIIsOl5QyUO8van3h9V5mGJYhBwqVhali0liw/lslzrdeWJ0z+dP5RH2pQkyGEkEBQk0nDdTKEEEIKo7AV/77p+C0rqX3O9a1T9PqeotbaZJUhmb+ker6mucFYP+RaYZ885uv4zxNYECGZg6UyDd/6A52bxDcQIXSQiSZDOze0CY0UA81lhBASCJrL0gQZZPJsIpbVQZ9Hg/CVEfLcSmov2szdOiuW6ofQUi318szcI1x91QIWfLU4q7agOf6tffX9UcrqmHe1mfXarTJ8+1300gIXHGTS0CdDCCGkMHJpMpWa0WaRkVV+UVpZUfWzztyksFwpV5jWdp5Fja4+hZIRvz+SBpHst3RffN9b1zMJ4ZOxnOeaHWvvupZBIs+sO6vWkicM2TczeB6oyaShJkMIIaQw6PgnhJCAnM5aSRYKG2RChghWIoQ5xKr0rO37hm9LIbUDHU/WCRmwYF3xn5WiTG8hAifiaPdbqiNtqKY975AOd6t8zcTkm6nANwzaClfynxpQkyGEkEDQJ5Mm1yATzXLyOMsjBmMRn294q+uYr0blK19zRFtlWYMeLM5STW6I8GMrvm1J/Y6uV9ISfcOaXTN9S6CFJn8guQPVidfz1Yry5DpLyh+ojbz1XTCEeXChJkMIIYHgIJOG0WWEEEIKI4gmczpsQpZVfp4V+RFWc2PI9q2O/Khv8fxm2pqSCN93Ig9aWy5TSdIEmWXtj+WdkUxX1vvicr4nsW4lkDVnmG8/pPp55BexyRkpDprLCCEkEDSXpfEeZAZ7Rby2alsij5agHQuhwfj2I2vbrntluae+odSuMOuQWPuhaQSa1pKn/1k3KAu9qt5aph0LkSEgpHbjm6E5TlGaDgeZNPTJEEIIKQyzJlO0XT30bDeEb6gIrSJOiJBk7bgWPhvHMruUztNkWcN+Q85epfohFl66cpElrznEvi+ukGdfGZZFm3kWTYbID1bJxZVF/Z5Rk0lDTYYQQkhh0PFPCCGBoCaTJtcgUylHt7W+pAL7mkvy1HE5x/P2I6sJyFU2bNgwAO4wZZ++ZumHr4wQbWqO/+i++D53a1/zrL6XSOZEC711siVQwRpoUbSZVOtjkXCQSUNzGSGEkMLw1mSyaC9Fz1ojrAsAi5qdS20XHZJswdof34WIIYITfGfWrra08GBfTUa7B3lyelkyImf5DlgyJ/tqJr7OeOu9CKFx5NFyitIcqMmkoU+GEEICwUEmDc1lhBBCCiO4JhPSPFRJJ15RK/h9Zfk690MGFkSO7rhczXnsu5bDdxV2llXbybUtkglVyyMWvwfWZ2F5Bta0/hG+Zi0rVlm+pqiQ7bvMn74BC5WEmkwaajKEEEIKw1uTybOaOYRmEGJlcUTI8GYgbJaBEOHB1q15ffsRlcVn/RFZAy4GkqdhaSsuMx6inaxndfxbAgtcs9IQ34OQq+99t3D2JU/2Ak2eb38Ywjw4UJMhhJBARINMlj8fNm7ciBtuuAHNzc2oqqrCmjVrSsf6+vqwcOFCXHTRRRgzZgyam5vx5S9/Ge+9916ZjP3796OtrQ21tbWor6/H3Llzcfjw4bI627dvx1VXXYVRo0Zh/PjxWLp0qfc9qUh0Wd6ZmzVvVgh8F3u6ZPiG+SbL8rRpCdkdiGS7cS3AknNN0njz2M61c63atTSLjjQd6ZqsIcwaWu63OMmFlJqsONp98fVphdBarPcl5LbOp7MWMBBHjhzBJZdcgttuuw0333xz2bEPP/wQb7zxBhYtWoRLLrkEBw4cwDe+8Q380R/9EV5//fVSvba2Nrz//vvo6OhAX18fvvrVr2L+/PlYvXo1AKCnpwfTpk1Da2srVq5ciTfffBO33XYb6uvrMX/+fHNfq04an1BNTU365AAOfMt5obcXsPwI5+mHNsi4BiDfQUYz81jrSz/WI0eOBPCbHwBpkNH6NhQHmehatEHm2LFjar+lNgfql7V+sk8D1XelzLe2lWzT1ZZFfp4fet9+WOpYy3p7ey1ddNLT04O6ujps374dY8eO9T7/gw8+wMUXX4xDhw6htrbW69yqqio8//zzuOmmmwass2XLFnz2s5/Fu+++iwkTJmDnzp2YPHkytmzZgilTpgAA1q1bh+uvvx6//OUv0dzcjCeffBJ/+Zd/ia6urtLv/3333Yc1a9bgZz/7mbl/NJcRQkggKmUu8+XQoUOoqqpCfX09AKCzsxP19fWlAQYAWltbUV1djU2bNpXqXH311WUKxvTp07Fr1y4cOHDA3HYuc5k15FAjqyktdNCBrwaW1ww20LGkdpDFQW9py3czL9/Q3pBargvfwAzNJCYdC2HasfZNasea48w3NNoiy5r/zPoMsoZcW88LrWX5ktfx39PTU1Y+cuTIklUhK0ePHsXChQvxJ3/yJyUtqaurC+PGjSurN3z4cDQ0NKCrq6tUZ+LEiWV1GhsbS8fOOussU/vUZAghJBB5NZnx48ejrq6u9LdkyZJc/enr68OXvvQlnDx5Ek8++WSIS/TGrMnkmXH6znayklXzKKrN+PGs2k1oTUaazUtZmH0XYybPG2x8c2hp9yVOyBB6CatzXzs3j+9Lw/c7VLQGkcefN1TZu3dvmU8mjxYTDTDvvvsuNmzYUCa3qakJ+/btK6t//Phx7N+/H01NTaU63d3dZXWi/0d1LAyNXwRCCDkNyKvJ1NbWlv1lHWSiAeadd97Bf/zHf+Dss88uO97S0oKDBw9i69atpbINGzagv78fU6dOLdXZuHEj+vr6SnU6Ojpw/vnnm01lAAcZQgg55Th8+DC2bduGbdu2AQB2796Nbdu2Yc+ePejr68MXv/hFvP7661i1ahVOnDiBrq4udHV1laLpJk2ahOuuuw7z5s3D5s2b8dOf/hTt7e2YPXs2mpubAQC33HILampqMHfuXOzYsQP/9E//hCeeeAJ33323V1/NIcyjRo3yEiw2ltEZ7NqMrGinvda+bz80J78kN4u5TAsesIblRhElUZmUA6zo4A4rvk5h7VriIa3We5VsX3Lau/qYrCfVdwUDJMusodSS6c8aom055gq99g3ftvTDytGjR3PLAH4TwvzGG2/gzDPP9D7/8OHD+PSnP20OYf6v//ovfP7zn0+Vz5kzBw888EDKYR/xn//5n7jmmmsAfLwYs729HS+88AKqq6sxa9YsLF++vKz/27dvx4IFC7Blyxacc845uOOOO7Bw4UKva2Oqf0IIOcW45pprMk8EIhoaGkoLLwfi4osvxo9//GPv/sUZlO2Xfc9zbUZmCSywajvW3FshFnSGDAHW5Lr6HX2WFilGwQBx+dJiTC2DsuaIDu3U1vqhhWpLWy1nDWKwhjxrIcOud1NawGiRIbUvXaem3UianSZfwjdwIk/YtKU/ocgbwnw6Qk2GEEICwUEmTWGDTMhw4tBhsb5+F+086bgmP4//SNJatHpWzSdCmhVrPh8t3FfSZFxfpGR/rfWtCxItWqpLk7Fubewrw6KpSX2zakrJvrrqS/J9F0ZqaNpeiOUSgwUHmTSMLiOEEFIYNJcRQkggqMmkMQ8yReUpC7l9sC8WE5mrzBXC7Gui0Y5ZTW6aXMnBLaFtIGYJYfY17cTrRVjrW7c21t7X6HqlLajjROdKG6BZyZo1wGWO890uQAsUkIICQm4UKJnmNPKY7SqZz4yDTBpqMoQQEggOMmmG5CBj0VZCLAR0hUZbF0QmyRNOnDwmyXJpT5oWkgzZzVJPCx6Q+qVpPpqMeNvSLFqSawmbtd5v6dqlkO6B6sTlWrWQPFg0JKuWY+2Xdn1ZQ9dP5x/c/x8ZkoMMIYScilCTScPoMkIIIYVRmCaTZyVy8rws9SxrXEIEHfg65q3mMuue81azWlKua0W+tuI/q7nMZaLTsJgAAX1GKJl2InOQ5PiPPrtybyWPucx8EiGCAZLPxSXLMnt2BRb4mr9CmMm0HG2DDTWZNDSXEUJIIDjIpPEeZDQNJV5mJZodWbWKPCvmI3yzCFvlavV8tRDp/76r762aTPK8+GctICJPsETWLA5aAACgZxuW3tfkvbJuM605taXN36waiqQ9Wc/11eI0tHdHCye2fi+s/cgaEFHJsOVkGxxkyqEmQwghgeAgk6Yi2y/7krzhvgskJVwzZ4sfxeXXsfhd8iyk9A1r1vw6rtxi0ow6QrqXFu0wRA46q78j6rcUpm7Jxhwvc4XnJttyvSe++c+kNi34LtSUyLMA1Jq1WdM+Q1KJ3zAOMmkYXUYIIaQwaC4jhJBAUJNJU5i5zNecpeGbIj6Or3PatZo+KddaP2vfXKYxqR9J574r/Fgzl0mBGVkd/6G3bNBS5mvhxJK5Rwqz1gJbrI55yWxnkeFKu69hddpr9a2bxPluORDCNKYFG2Td0CwUHGTS0FxGCCGkMLzNZVlmCMlzsjjhLWih1C6ZljDoLHnEfDWTCE1TcmUKTmowroWR0sx6+PDhZWWu4AGfY8nPEZaAD9fCSGlbYu1Y8vqOHz+eat+V08uiocWvRcrgHMmIjmUJHrDg2sI52V9r0MNA//dBalPb9C3PEgpmYa4c9MkQQkggOMikobmMEEJIYQTRZKwqap5NwipFHke+JsM3F5n0f82Rr5nEpGwArjxiSce/ZKKT7o81x1iIFf/a7E8ySUlO+GRutshMOBDRNcTlW/KZSVkJpGvx/S5pwQ9WLKbL0ORZr5O8dpesSmoJ1GTS0FxGCCEBOZ0HjCwUFsJszdDrW2at79u+hlWWJcTY1cek9uHSWiRtJanJxGVYHf9JjUeSoV1naK3VmosqqXlJ2oLmtM+y4j+pTUiOcWkjM+v2zlGZK2+bNYw4ea5vaLRVVtZ8Zq6QZEk7jMgT5h0CajJpqMkQQkggOMikMQ8yIRY5hZg9aLPhEAsvNU0jj7/GupAyeS2urMBaWLMkQ/LrSH6UESNGAJCfe3RMur48z0DDmvMqqcFI5/X19ZU+J5+BNBOW/DtWLVsKSdbywkk+MG0BY56sykmNwLr/jEVWXJ5rMasvIWQwhLlyMLqMEEJIYXCQIYQQUhje5jLX6lrJ7JA3CCCLmcVizrCav7T6rpXwFnOQNVBAMo1JIcnJetoxQA7bjUxi1hBmS4CD7zNzmRC0UGBrmHAyhFkyjVn6aj02UBvJd8caqi3d5zxbOEtyfepbMwRoZb7mNVcYd6W2EojaoLmsHDr+CSEkEBxk0nhrMtZ6UohniLBiqV7I7L5ayLBvuLJVRpzkjFZyGLvCiZMajGtbZcmRbwlhljQkLZeadJ2+uPKIJcvi9SMNQtOypXsRx6LxSNsvS/22btOsBU5offTVaELkRnP1w9dpr2khWi46q4zQcJBJQ02GEEICwUEmDR3/hBByirFx40bccMMNaG5uRlVVFdasWVN2/OTJk1i8eDHOPfdcjB49Gq2trXjnnXfK6uzfvx9tbW2ora1FfX095s6di8OHD5fV2b59O6666iqMGjUK48ePx9KlS737GmTFv9VRJ60ZkEwpljU5IbIG+GJ18mvmMqm+7xoXV1nS1BV37EsmtOh4vN/JYABJvvQcrTnaksdcWM0gSYd/3FQT9TtupkqavVymFy0Vf1EkgxMAmznLVT/E7Dn5fZcc+lKbVse8FsBhaUfqa5FUSpM5cuQILrnkEtx22224+eabU8eXLl2K5cuX49lnn8XEiROxaNEiTJ8+HW+//TZGjRoFAGhra8P777+Pjo4O9PX14atf/Srmz5+P1atXAwB6enowbdo0tLa2YuXKlXjzzTdx2223ob6+HvPnzzf3leYyQggJRKUGmRkzZmDGjBkDylq2bBnuv/9+3HjjjQCA5557Do2NjVizZg1mz56NnTt3Yt26ddiyZQumTJkCAFixYgWuv/56PProo2hubsaqVavQ29uLH/7wh6ipqcGFF16Ibdu24bHHHitmkMnjXE/OIHxluR6AtkVwSKzh25rGIznQta2NpZX5UpmkaUTaiKS1SNpNnMgBbtG24p99NTxfXDPgpDM4PgOWwrEjLUGqL/W1qNlw8lqkrM0SvvnKirb9a88i/tm3H75aXCU1zTh5B5menp6y8pEjR2LkyJFesnbv3o2uri60traWyurq6jB16lR0dnZi9uzZ6OzsRH19fWmAAYDW1lZUV1dj06ZN+MIXvoDOzk5cffXVqKmpKdWZPn06HnnkERw4cABnnXWWqT/0yRBCyBBh/PjxqKurK/0tWbLEW0ZXVxcAoLGxsay8sbGxdKyrqwvjxo0rOz58+HA0NDSU1ZFkxNuwQHMZIYQMEfbu3Yva2trS/321mKGIeZBxbXIV4bvC3XI8hEM/tMlDM0XkCUrQ1tVoTn7NuR8/JpnQtLUzWpYBzfSXxVyWNKNmWcmddPhL5jJtHYtU3xfr+yWtH9FWvVvbKMpMpiXG9D3Pd8W/S16WfhVBXnNZbW1t2SCThaamJgBAd3c3zj333FJ5d3c3Lr300lKdffv2lZ13/Phx7N+/v3R+U1MTuru7y+pE/4/qWKC5jBBCAtHf35/5LxQTJ05EU1MT1q9fXyrr6enBpk2b0NLSAgBoaWnBwYMHsXXr1lKdDRs2oL+/H1OnTi3V2bhxY1nW8o6ODpx//vlmfwyQ01zmO8sIEZKcNZW8K/xYajPran2tv642NY1AkqEFA0hOfq1M0mQsIdJSPVdWAsv9c83INe1DWt0flUkhzJKTOio7fvy42m+tj1JafykkOamVSdcery/JSIZXuzQCjawbgknP2JWlwTcM2jdrgCYrNJWKLjt8+DB+8YtflP6/e/dubNu2DQ0NDZgwYQLuvPNOPPTQQzjvvPNKIczNzc246aabAACTJk3Cddddh3nz5mHlypXo6+tDe3s7Zs+ejebmZgDALbfcgm9961uYO3cuFi5ciLfeegtPPPEEHn/8ca++0idDCCGBqNQg8/rrr+Pzn/986f933303AGDOnDl45plncO+99+LIkSOYP38+Dh48iCuvvBLr1q0rrZEBgFWrVqG9vR3XXnstqqurMWvWLCxfvrx0vK6uDi+//DIWLFiAyy67DOeccw4WL17sFb4MAFUnjVd3xhlnfHyCY7ZumaH6+ies8q3HrAsGLQsppfraRmMujUDLGaYtpIxrJpEWYjkW/xyXGzkcs/p1XFmbffNxWRZeAmnNRNJa4ppJ9Dk6duzYMVP9uAkhKrMci8uN9y3ZX6nfUpl0P7RjWoixq762MNI3Z5gr1FmqN1CZVb50/MMPP0wdy0JPTw/q6urwL//yLxgzZoz3+UeOHMGsWbNw6NCh3D6ZoQZ9MoQQQgqD5jJCCAnEyZMnMznxT+cEmeZBJrpxrjBGybRkyUWmETcTaE57K9LKbw1f81qcpJnMN9eZK1w5+hxPUZ80k/ma14DfmMukY1JZsm/WzAZxfB3/Wphy0gwWL5PMWVFZXH7cxJUkZDRQvF3fQJX49fl+z3w3avNdwW8x0Un1rWVaQIErB11RP+qV8smcSlCTIYSQQHCQSRNk+2VXWfKYNKPVZEgzM5d2o82KteAFTTOJ/9+qBVnwzQHm2kAsWSZpOfGyKDdRvCwK9IiOxfMXRZ81jcq6yVkcX01GClOOPkeaieS07+3tLZVFn6N/4/I1bVmbRUu50VzPMVkmXacVzWkvaT7a7F+6t6G1lqyO/DzBBkXBQSYNNRlCCAkEB5k0jC4jhBBSGGZNRlulrOUhiqOZxqyBAllNaFJ9Kd5fWk8THZP66DIBafUj4o7oJNI1SSYMaXMuqV9SQEFk6oov1Ipi/UePHp06FgUFxMuSprm4KU0zl+VZJyOZy5KO/Pj9O3r0KIDytTBR2UcffZSqL5mYIvlav13PR5KbXCcTfye0NT/a/ZCOaWYn1/oXy9oZX9PYQHIH6qOrvrUsdOBGXG4W2UX1ZyhAcxkhhASC5rI05kFGWsUaMheZ1eGuzYC18GBXfyzOWC3HWPyzdTtlS94xSeOIO+G1kOTkv8BvtJB4CvHIyT927NhSWZRlNSqLP/9keLPUX9cmZxHac3fNgCUtIbk6Pu74j7SDuCZz5MgRAMAHH3yQ6qMl23T8s5aRId6PKMhAygIQlcXvrabJWPO2aRqVlrfNVwuxrrS3Ov4lLcvSHwk6/gcHajKEEBIIDjJpvEOYpTLrAk2pfjJzLJCe8Vo1Gc3X41o0aQnRlnwyki9EKvNd+CltKS3Z8qUsv0msPpm4dhNpMHV1dQB+o+0AclizlrVZ0wR9fTLaTBxIz+Ljvo1Ig4hfZ1wDBIADBw6UPkd+Gl9fkjV02OofyZrPzNpmVq1F6m8eGcl+uWRodXyOh4Q+mTSMLiOEEFIYNJcRQkggaC5LkyuEOSuaSSr5OV4n3r7Vka+FJrtMaBGSKivdB4t5QAs1jmPNYaVlAbBu1xyZveIhyZF5LPkvIDv+s+Yu83X8a+ae+GfJgS6FV0dEZrX4PYjuSzxQwBLYYF3db33/ktcG2EKMXQ70IlL9W81lmkkstEO/kiHMHGTSUJMhhJBAcJBJ4+34d80AtBl+1izM0uxICxRwITnVI6xh1lL9kBufSfUlLcSyFXL8mDSb1/KTRfW1Y3F5miZTiU3LkvnDpK2T4zIiDcaaoy265rjcqK1I24/fby08WNqC2KrBSt8DLfRf+95IudEswQlSPV/NRyLPIstkX33aCAEHmTTUZAghJBCMLkvD6DJCCCGFESTVfxxp3UuEZDaxmNKsGzhJ50jmB99NyDQZkvlLcrRrZi1Nrsss6Gtyiz5LJjfJ/KWZjLQV/4NtLtPkx81DyevTghni12LNLmH9Hlieu3WtlRakI82YfbNzaLNuyewT74d15b52zDcLgOVYKGguS0NzGSGEBOR0HjCy4K3JuGaemqZhmXVLbbic5dJ5WR3uA7VrabNSuPqhhdRqYc3SLF4LNpA2LdPadN1vLYRZKtOc35L8aOarBU5o9yB+LVqI9qnyniQzWUgZvyU0zUQKTpDalwIKpPO0QAFNUxusLADUZNJQkyGEkEDQ8Z/GPMhYQyuTxwY6VzuWnBFatZaBjmclGQ7rmllrC0ulrXkjNP+Sr98o/tk3e7QULq1litbCprP4wDSsi/2S8qxZr30zLlvfV2tZEpcfw5IrzKoRhPiB076DrueUrOfqT7J+XL7VD1UU1GTSMLqMEEJIYdBcRgghgaAmk8Y8yGTdeEoqc5lIkivyXWHTVrlZsQY9hHhRkmHQWUKBLSG18TLffFyWflQyhDmOFsJsMRVK5rXQ91szT0qr732xvodF/7BZc6JFaGYtKWDBaqKrpAmNg0waajKEEBIIDjJpggwyvjfI5SyPyiwbmrkIod1YAwvyBEfkxarhZXVO++ZcC+H8tuLbD9+FlL79Dq1RS/I0LV/7TlmtAr6EyOQ80P8HKstav0gYXZaGmgwhhASCmkwaRpcRQggpDO8V/1lWBfvmJ/Ppj5W4Ouqb4l+r43LoWtqR0NZGuPA1dfn2Leu5rgAR3za1dyBEH6Uy33sqkfXZuuQnTWiu9Snaiv9knfhn303L8szSNRm+wQPMXTY4UJMhhJBARINMlj8rJ06cwKJFizBx4kSMHj0av/u7v4tvf/vbqQF+8eLFOPfcczF69Gi0trbinXfeKZOzf/9+tLW1oba2FvX19Zg7dy4OHz4c7F5EcJAhhJBARI7/LH9WHnnkETz55JP47ne/i507d+KRRx7B0qVLsWLFilKdpUuXYvny5Vi5ciU2bdqEMWPGYPr06Th69GipTltbG3bs2IGOjg6sXbsWGzduxPz584PeD4COf0IICUYlzGWvvvoqbrzxRsycORMA8Du/8zv4h3/4B2zevLkka9myZbj//vtx4403AgCee+45NDY2Ys2aNZg9ezZ27tyJdevWYcuWLZgyZQoAYMWKFbj++uvx6KOPorm52fsaBsKsyVRVVaGqqgrDhg0r/VVXV6f+omNR/aqqqtKxeFnyzxdNlvQX76O1/eQMQ1Jt42WW2YlVPY766iqTSKrh8f5oarq1b1Idn/OymAh8zw3RR61t7Z66yPpsXW0k3z3rM9Zm0lm+S5X6nkt9sPw2xX+fQlMJc9nnPvc5rF+/Hj//+c8BAP/zP/+Dn/zkJ5gxYwYAYPfu3ejq6kJra2vpnLq6OkydOhWdnZ0AgM7OTtTX15cGGABobW1FdXU1Nm3aFOJWlKAmQwghQ4Senp6y/48cORIjR44sK7vvvvvQ09ODCy64AMOGDcOJEyfwne98B21tbQCArq4uAEBjY2PZeY2NjaVjXV1dGDduXNnx4cOHo6GhoVQnFPTJEEJIIPJqMuPHj0ddXV3pb8mSJak2/vmf/xmrVq3C6tWr8cYbb+DZZ5/Fo48+imeffbbSl2siiCbjCk1NllnDS0NufJaH6AWQwj+lelrYZzw3VsiwRU2WK6zUUqaFrUpl0rE8ZO239iys9V1lA/U1FJI8LdW/9m66yrLi+33T8olJ3zPXd08j9PMosr29e/eitra29P+kFgMA99xzD+677z7Mnj0bAHDRRRfh3XffxZIlSzBnzhw0NTUBALq7u3HuueeWzuvu7sall14KAGhqasK+ffvK5B4/fhz79+8vnR8KajKEEBKIvJpMbW1t2Z80yHz44YepQXrYsGGlQXrixIloamrC+vXrS8d7enqwadMmtLS0AABaWlpw8OBBbN26tVRnw4YN6O/vx9SpU4PeE7MmI82cpIVp0swjWV+SESfEpmUhybNI0Zfo3kQL46TNq0LM5uNlUVvxxXhRW8l/rf2QHMm+OcDihNBMtGuJ/o3fg+hz6Pst9UNqKyvS90yrVxTWxddRP6N3faCAmYHKXFqzZYO3UPg68ePnWbnhhhvwne98BxMmTMCFF16I//7v/8Zjjz2G2267DcDH9/bOO+/EQw89hPPOOw8TJ07EokWL0NzcjJtuugkAMGnSJFx33XWYN28eVq5cib6+PrS3t2P27NlBI8sAOv4JISQYlRhkVqxYgUWLFuFrX/sa9u3bh+bmZvzpn/4pFi9eXKpz77334siRI5g/fz4OHjyIK6+8EuvWrcOoUaNKdVatWoX29nZce+21qK6uxqxZs7B8+XLvvrvgIEMIIacQY8eOxbJly7Bs2bIB61RVVeHBBx/Egw8+OGCdhoYGrF69uoAeluOdu0xCymemmdBcanpSXZbyREl9y5OPS0LLXeZrytPWQMSvNx4YAORz2ktmAovJCEib0LRj8euMCLG5nctZLZnLLOYv63VaTYXJshDBA64V4JYN/Vz3O2mmyrNhmtV0pZkPpfN875/0fDS5oamEJnOqQU2GEEICwUEmjXmQ8XXQSjMxaXYraSHajFbTDOLnJY+7sjBLMpJhx1I/XEEMReCapVmc9tIMv6+vL1UW/Xv8+PHSseHDh6fKIqR7pT3j0I7/5PXF+6iVWe5BXL6m8RQdLmzF2g/N0e4baKFtRhY/rmkVVm3IGiAgUdR3lZuWpaEmQwghgaAmk8Zbk7HuaaHJcJUlcc2+NA1JIjrXep519hddi9W2Ld1TbVbsO+PUZvrSbD4+i4+O9/b2AgBGjBhROhavN1Cb8XubDEmPfy5Kk4n6H7/OqN/xsuj6ojJJk4nXlzQZSwhzHIuvTDomaZ9ayLVv2Lk19DqPf8lSP45mMbD6OZPn+bTvCweZNFyMSQghpDBoLiOEkEBQk0njbS6zbl1sPe6b60wyvUj1kuYYlwytj5J5TVLjLSGbvjnApHDbuPlGei6R810yqUTmoLj5KzIZRf/GP0f148eka0+2JZnLXFsWawEfmvlGM5dJ1x6/lmSZdA8kk5sWGi0dk8yTvia3PPnjknWS8pJomTh8ndRSvj5XgIAFLbOBK0dgUY52DjJpqMkQQkggOMikMQ8yyUWCcXw1GZcWUsQiyBAyrDNxrcyy6RjgP6P1XZAoOb/jW7N++OGHAMo1nojo3JqamlJZVC+6vvj7Ij0DTSOVsF57UkuIO/IlbSW6zujf+D2QNBnfPG++WoiGK0eg5XumZT+WsGZEtuYH02RI/bHUz4L2e5YHDjJpqMkQQkggOMikYXQZIYSQwsi14j+vQx+wOeGtpi7NnJUlx5hmcpNMaJqpSDoWfba2GX2Oy7CsN5FMGHFzT2QOOnbsWKnsgw8+KJMRNxlFe1zETWlRFoCob9o9iOPr+I9jybkmOe3j13nkyBEAv7ne+DEpeMCyzsP1zCJ50vukreFytZXMQaatv4oflxzoWvYCzYSWZ12XJl/qt9WUpgUIhIaaTBqaywghJBAcZNKYB5loxmcNYdbKrFmKfbUQ6bjV8W8JdZY0DkmrkMKJpfqSJhNpB9GxSEOIf5byiMWd8FkzC0vO7NGjRwNA2T4UkSYTL4v6G/U/3m9Jiwux4l+6vuTK/fg1RU79uLYSlX300UcAPt5BMCIqi9ePPksBBVGZdCz+zKRsBEmtKS7DuqFa8n5IxzRtxTdEWjrXFTZtCfOPU9Q200WGMGeRzUGGEEKIE2oyaXKFMFu1leQx1/4wFp+MVa61vqbxSMeSGkr8s+SPkLSW6HPct5Gs59KUNH9HhNUnEw/fjdqNZtTxY5HWJGlZyVBm17XE0XwyUpmmqUl5ypL52OKfpTDurD4Z60LU+P3QZr9SrrhIhpYnT8r3Je0BFLXt8oVo+cCsPhnNPyIdk7J6W7Qm6doluaHhIJOG0WWEEEIKg+YyQggJBDWZNEFS/Yd25FvKfE1jWUKBk/WswQO+aGGl0jGX014zpWh9lNqSwn4l534yUMFqMtL6Zs29Zc1dJqXzT5ZJIcxa/XgbWgCFNQuAVN8X7RlrOb0ks6okyxp+nNWE5soyoJnctECB0FkDJDjIpKEmQwghgeAgk8Y8yGi5pvI48i1lVvnWY9aFbhG+C+TiSI72ZN+08EzrBmhWpC2TIyRNJqof11qiWbwWju3KXRYihNmSu8yVty2p8cSDAiz142Wa5uPKf5a8Put1Wp3fUpmvtpQ1C3MczZHvi2+gQBxqMpWDmgwhhASCg0waRpcRQggpjMLWyUiEcIxnMVkNJM+aOl1Sy7U1K3EzSPK+Ses8fJH6LZlSfIMTJHNZ1Fbc3KPlJ9PWA1nXOyWvbaAybS2HtIZGMlNpWwMkgwjinyUZls3IrH1zre7X5Gr3RZPhWpmfdfV96BX3WqCAhGRCKwpqMmloLiOEkEBwkEmTa5DxvTHWsEhtpiJlmI1wBQgkkcJ+82hKUn+jGVs0m4rP4DQnpRawoGkt8TJLX+PytLBwV861ZJiy71bLyXPi/Rqo35oz25qjLTn7jzv+swYPaIEIA11fsr/W2b/VqW0JBvBdVW+VYe2bVsfnuLVOEXCQSUNNhhBCAtHf35/JFFhUws6hAB3/hBBCCsOsyeRZr5E0jbi2C0ianSRTmiZfQkvIlwXfflhXHfuaH7S2fIk/46SzXjIBac79LOYyC1Zzmeb81tabSMkotTUx8eOWoACpjwNdS/KYqyxC+97kMctofbSc55JhLUt+D4aSqYnmsjQ0lxFCSCA4yKQxDzLWWbRlhuobkuxaKa5lIwiJS772omj3T3PuF3VN0uwy3lZyRi85/qV++67u913xL5VpYc1SOLGkaUj/j+q7VusnNRgtRHqgMl8NzFcbiig6f5crUEXC0g/pOWrn0/E/dKAmQwghgeAgk8Y8yGg+jdCZmbVcZxa/R7yedcGWBZd92KKVSQseQ+cn05BmyhLJ+yXNxF0bqkVYt+y24Bva61oEmTxX8slYt0LW/C/WMumY9MwsIcmu+iHwDT+2LNaWZEjvkBYiLWlsp9NizF/96ldYuHAhXnzxRXz44Yf4vd/7PTz99NOYMmVKSd43v/lN/OAHP8DBgwdxxRVX4Mknn8R5551XkrF//37ccccdeOGFF1BdXY1Zs2bhiSeewJlnnundfw1GlxFCyCnEgQMHcMUVV2DEiBF48cUX8fbbb+Nv/uZvcNZZZ5XqLF26FMuXL8fKlSuxadMmjBkzBtOnTy/b+bWtrQ07duxAR0cH1q5di40bN2L+/PnB+0tzGSGEBKRojemRRx7B+PHj8fTTT5fKJk6cWNb+smXLcP/99+PGG28EADz33HNobGzEmjVrMHv2bOzcuRPr1q3Dli1bStrPihUrcP311+PRRx9Fc3NzsP56b1omlVlXcudx/EZoq+pD5jXLIs+3H5pqLzn+te0CrCv+kxkIBiKZoUDKU6aFPGcJzLDkLpOOa2Yh12ZhyXBsyVwmOf4lM5x0zBI2HT+edbW+dH3W+lYTqi8hUuxr74wUWBAiVDsPec1lPT09ZeUjR47EyJEjy8r+7d/+DdOnT8cf//Ef45VXXsFv//Zv42tf+xrmzZsHANi9eze6urrQ2tpaOqeurg5Tp05FZ2cnZs+ejc7OTtTX15cGGABobW1FdXU1Nm3ahC984Qve1zAQNJcRQkggokEmyx8AjB8/HnV1daW/JUuWpNr43//935J/5aWXXsLtt9+Or3/963j22WcBAF1dXQCAxsbGsvMaGxtLx7q6ujBu3Liy48OHD0dDQ0OpTiiCb78cEWI2b11UpuVCc+VJS+JbX0LShqwBC0ltSGrbFSig3Xvrc4lm9FLYb7T9chzf4I7kMRfWGX6EphnEtZXkeyJt6ubaTjkZUODK/OyrhWho2pBVk9Hk+ta3HstTliT+zmmh/0WHb0dy82gye/fuRW1tbak8qcUAH1/jlClT8PDDDwMAfv/3fx9vvfUWVq5ciTlz5mTseXFQkyGEkEDk1WRqa2vL/qRB5txzz8XkyZPLyiZNmoQ9e/YAAJqamgAA3d3dZXW6u7tLx5qamrBv376y48ePH8f+/ftLdULBQYYQQk4hrrjiCuzataus7Oc//zk++clPAvg4CKCpqQnr168vHe/p6cGmTZvQ0tICAGhpacHBgwexdevWUp0NGzagv78fU6dODdrfXI5/X0KsWbE4ywHd3GR1DlpMV3nMPZJjVNoSINlW3AkvmbOymi7ix5JmL9cGb1rusmT/XWW+/Y6jrZyXyrR1MtJ7oq218d20TDN1SSZRq+lKwte0qNXX+qGtzB+ozJcQMorMpJHHXGbhrrvuwuc+9zk8/PDD+NKXvoTNmzfjqaeewlNPPQXg42u788478dBDD+G8887DxIkTsWjRIjQ3N+Omm24C8LHmc91112HevHlYuXIl+vr60N7ejtmzZweNLAMYwkwIIcGoxCDzmc98Bs8//zz+4i/+Ag8++CAmTpyIZcuWoa2trVTn3nvvxZEjRzB//nwcPHgQV155JdatW4dRo0aV6qxatQrt7e249tprS4sxly9f7t13F1UnjVc3evToj0/IMRsNObu1buVrDZ+2hGNbMwtr2Q5860uZjqWV9pJWIeUYk7ZOlohswVoIs3RPNblSH31xzZQjNE1AC2E+duyY2r7myM+z/bK2aZklUECq55KhZTO2Ov4tAQu+Gpiv1ppH3kcffWSS4aKnpwd1dXX4yle+gpqaGu/ze3t78cwzz+DQoUNljv/TAWoyhBASiEpoMqcaQQYZ1w1K7kki2fetFJ1pOY/fRZORdcGoNPOUNAJrWLY0Q5V8D8OHD3fWl65F2mtG6oclVNxVxzd02aIlxEOYrT6Z5GzeutWy1T+maSF5/DSa/KLJGg4NcD+ZUw1GlxFCCCkMmssIISQQ1GTSeA8y1puhreCXzAmuEFnfNpNI5h7puHYsT16zCJcMzbQYITmdLTmeku0n24zfo8hsJPUj6+r+ELnl4mi5sazmsuT1xe+tdF80GdpGaa5gAIsT3prq31eGRNZjeVbV+zr3iwqRzgMHmTTUZAghJBAcZNLk2rRMwpJvLI+Gos2UJazBBtoCTaktTdPw7Zu0uDJCmln7an0uJ7/m4Na0OC0Ls0QlNJkIqd/aeckwZMCuySTbcuUps4b7av3VsGhKrvpSf7L2zfV9yKvxDCWNhoNMGmoyhBASCA4yaRhdRgghpDBymcuKWrNiMYm5zD3JYy6He9b1OnmCAbR+a2Yw67W4AiwGOg+Q82sl29dSrEv9lsxrvrjMM76O/4iov5K5TJJhDToI4bSXyizraVybhiXrW78DvpuRueT6blaXlUpoC9Rk0tBcRgghgeAgk8Z7kAmxhXKcPCvgNbRggxBOe2sQg6W/vhqHdG68fjQbt27AJs3ikzm6pBBmayirdC0ht/q1huVq4cdR8EX8HkRlkoZn1WTy5CKz1Ld+H4oKJy7aCS9tsWxBW0KRt08aHGTSUJMhhJCAnM4DRhYKG2RCaDwhFn1pIbu+/hTXy6MtXPR98Xxznblm4ElcPpOkT8YV2qvJz+ojiqOdm0eTiZCyK2sytJm1dS8Y37BmVzhx8j676luw+l+smlWIH2DNSjHYP/DUZNIwuowQQkhh0FxGCCGB6O/vz+RvDOmjHGoMyUHG1/yVPCbJypJLy9fcFBEix5lVlnR9muNfcvJLbUlmsmQ/QpjLtLBma7iybz4u322vrTKymrUGOtdXruWHytcx78p7lzV0+XQ1D9FclmZIDjKEEHIqwkEmjXmQCenIl8KJNUJkaJaQtAQpQMDXwWgNly46GEBC21RMC2GOY3kG0jPzdQpb60sOdK2+r9M+TggHvTWsWZORlJWlH5b7HdqME6IN3++NFvwTGg4yaajJEEJIIDjIpGF0GSGEkMIoTJPJMzInVdmQjnTAngXA1/FvMWtY19Box1zp/5P9ceU/k/qfNLlI/dCcwi7ThK/ZU6uvmcukPsbrJO+zy/RmMUX5rn9xlYWQKx3LaoYLge9mZEX1IzTUZNLQXEYIIYHgIJMm1yBTKQecS3uxyNAc+nG0leWu67VsCGbdftk32CCPtifNzrVNyyz9cWUzDom1H9ozkLJOJzNLx49rGoE1NFm7Bte76evI1/oROqtyhPY9c4VGW9ocij/MHGTSUJMhhJBAcJBJE3z75Yg8s9cQ9vpkf10agSXHWYjFmy6bfwifkzYz1LZJlmb90gxfez4hQrVD4OtPsWZhtoRBWzUIrW++OcMsbSTR/GhWmSE0nqz1i5KRB674T8PoMkIIIYVBcxkhhASC5rI0QcxlLrNW0Y7fiEqa9LQ2fVOP+4ZD50HbOlkLkbX2R1sxH/o9yGqO0ZzfUvCD615p/crq+Lf2u5IO8aHmfLduo1HJPtJcloaaDCGEBIKaTJpcg0xR+cx8znNhDZvWggfyaBpaSLIkS1tMGBJXmG2yTMqWrPXNGqbuex+lc61tSYsxk859a+hwaE3G4oR3HffVOi2yQtTPQ1brhLaMoEhOnjyZSSvhIEMIIcQJNZk0jC4jhJBTmL/6q79CVVUV7rzzzlLZ0aNHsWDBApx99tk488wzMWvWLHR3d5edt2fPHsycORNnnHEGxo0bh3vuuQfHjx8P3r+KaDJZTWIR0gptl0qa1cykmV6sbbqCAaS2BqrvalNay6HJ19axSGY760p4rc28z98l19qW5tzXZErHpefim0lCk2+d2WrrZIpOoS9h/Y6EcHQPxdl/pTWZLVu24G//9m9x8cUXl5Xfdddd+Pd//3f86Ec/Ql1dHdrb23HzzTfjpz/9KYCP3/uZM2eiqakJr776Kt5//318+ctfxogRI/Dwww9n6stAUJMhhJBARNFlWf58OXz4MNra2vCDH/wAZ511Vqn80KFD+Lu/+zs89thj+IM/+ANcdtllePrpp/Hqq6/itddeAwC8/PLLePvtt/H3f//3uPTSSzFjxgx8+9vfxve+9z309vYGux9AhkGmqqqq9OdL1lE+TvRAIll55FllxOtpf5aXxSpL+3PJk+5R8twQ/Thx4kTpL9m29OcjLy7Tdb+1tnzk+dz7rPe9iHcgifV++777WfH9jljl+ZLnt8tKUc9UYsGCBZg5cyZaW1vLyrdu3Yq+vr6y8gsuuAATJkxAZ2cnAKCzsxMXXXQRGhsbS3WmT5+Onp4e7NixI+PVy9DxTwghgci7Tqanp6esfOTIkRg5cmSq/j/+4z/ijTfewJYtW1LHurq6UFNTg/r6+rLyxsZGdHV1lerEB5joeHQsJMEHGWlETs4c4nUGO+Qw2acsWyz7yHfJslxziO2rpb7F/RPaYszoc1xu0nfjyo1m8edk2dNE6q8mL3lMyl3mIqTvQ3s+IbDKKuq7Z5Hrew8q+TvhIqtWEp0zfvz4svJvfvObeOCBB8rK9u7di2984xvo6OjAqFGjMve1UlCTIYSQIcLevXtRW1tb+r+kxWzduhX79u3Dpz/96VLZiRMnsHHjRnz3u9/FSy+9hN7eXhw8eLBMm+nu7kZTUxMAoKmpCZs3by6TG0WfRXVCQcc/IYQEIq9Ppra2tuxPGmSuvfZavPnmm9i2bVvpb8qUKWhrayt9HjFiBNavX186Z9euXdizZw9aWloAAC0tLXjzzTexb9++Up2Ojg7U1tZi8uTJQe+JtyYjmbWs5/jWl9BkWMOWNZOer1nLtbWxbyh11i2ZNaTthuNoJjHN/KSVWd8TqxlEu/YQJiDtHki4NhAL1R8f+b5bGmc1XWlYNyMLYeLylVEJs1olcpeNHTsWn/rUp8rKxowZg7PPPrtUPnfuXNx9991oaGhAbW0t7rjjDrS0tODyyy8HAEybNg2TJ0/GrbfeiqVLl6Krqwv3338/FixYIA5seaC5jBBCApHXJxOKxx9/HNXV1Zg1axaOHTuG6dOn4/vf/37p+LBhw7B27VrcfvvtaGlpwZgxYzBnzhw8+OCDQfsBAFUnjVdXU1PjJ9gzr5mvIz1e37p5lm+blv7myd+mbSCW514lZ/3xeyHdF2nBZTSb0bIOWzco067P+uXy1WRcGmayXnQtx44dSx2TAgDispKLU63BCdZjUZnvgk6pLI8D3fKsXMEaIbWPEEESodaE9PT0oK6uDp/97GcxfLj/3P348ePYvHkzDh06VOaTOR2gJkMIIYEYKprMUIKOf0IIIYVh1mSs+9xnxWrysGLpr8s5bQlYyLO+R+pjVue+ttYm/n/NnKGZNVymGm0tlERRZpPk/XNdu2+qf60/ecq0OpaN0lwystYLEXRgpZKz+aI2CaMmk4bmMkIICcTJk9xPJon3IJMn+3HWUGaXk98aTqz1VZMhHbOW+WIJoc6jPVnbtoQwa6HaecK4fbH2w6JVZAmbzhoKXNQPSwiNyqrR+fZDuqdFaGiDtZ0xNZk01GQIISQQ/f39mSaZgzUoVoIgg4zvTDWP/8WicUhkWQhaRM41X20rjhTyrGlUWpua1hL/bPVHaNsHF/0F0tqSrknK/SbtBWR9T7X7HSIkWaII7akSM/BKLaC0LgoNDTWZNIwuI4QQUhg0lxFCSCBoLkuTa5AJ6dB1hRPnrS/hck5nDQaQsOQkc8m3hjxrzmypTDLfWFL9W6mkKUC7dq0f0fW68rz5ypXqaCZIa245n7Zd9Sz3JQ+D/b6E3PLZBc1laajJEEJIIDjIpClskMnj3NdkZXX4u2RoizetoZiWma9v/yX5Lg0sOWOzOv41rPXzhFT7Lui0tJ/FkW9pX9NMLM9kIPkhNJMQ8ovI6JynflEyQkNzWRpqMoQQEghqMmkYXUYIIaQwcmkyvvnM8ji9LcesaDKsGQI0uRK+ZjLpPC0YQJJrNQVJubGS5ibXuprkvRlsc4h1LYoWKKAFWoRsc6AyC3nuUYgV/L5kDUCo5HXmgeayNDSXEUJIIGguSxNkkJFGYV/tJqtGM9DxkDJ85VrOs54r1fe9dtcL7Bsi6/u8rRpe1pmsdg+0zMtxub4bjuUJnPA5L0/9EOda70dRMnz7PdgaAQeZNNRkCCEkEDSXpaHjnxBCSGGYNRlf81Aek0rURpaklkkZmnyXDK2+5TyXrBCOf80Mo12n1TmtmYesG6UNRhCAdl98TVfWbQvyrHHJen1D5d5mlQXoa5s0ik4wOhRln4rQXEYIIYGgTyaN9yDju2o6jqbdhA5XLiIM2nelvVTP+jJZtxKw5C5zrTa3OLOtMqT+aPdF0zCzBB1YZsVxGZbtl13ta+daNz7TyqxkzTxgpagQ4+S5koY8VLUWqS0OMuVQkyGEkEBwkEkTZJDxDcsNscmZRJ4waKu8iEpuxJUkT/4z3/atPpnkMde2xxatL8uGab7+syLCiK3+v6J/WCqx2LIo7SaJ7yZkg/WjnfV3gdFlhBBCSAZoLiOEkEDQXJYm+CCTNTDA1+RWlAxNru8GZQP1zQfXdVhMQNZU/1azliQjax47Vx42qS2LXCkfm0RULxkAMBAW85dvmLpPvbyysvQtVNtFyRhKP9AcZNJQkyGEkEBwkEnjvRgzTkhtwdWWRtYwZV/NZ6i8CC7tRQsntsq1OP6lc7WMzi60EGbfBXq+znWpfiUCLCpNHi0hxAZvg0kl+sFBJg01GUIICQQHmTSMLiOEEFIYp4UmI5nGfPeLz+O0HahNIHv8u68DXTNhWR3/8b4WscZCMr34mlmKWoMSXa9rzY3FhJbHtBxyvUTIVfgueZbME/8/wHUyaU6LQYYQQoYCNJelyWUu025oVVVV6S9Jf39/6U8jXi/5VxTRNQ3WQ4/ajq4z3h/tL889sraRvC+WOvF6rjJLH61lWfvtuhYN6Rn4Pr+h8O4V2X7W77RWT/utAbL/8GfB993Kcs+XLFmCz3zmMxg7dizGjRuHm266Cbt27Sqrc/ToUSxYsABnn302zjzzTMyaNQvd3d1ldfbs2YOZM2fijDPOwLhx43DPPffg+PHjQe5DHPpkCCEkEJUYZF555RUsWLAAr732Gjo6OtDX14dp06bhyJEjpTp33XUXXnjhBfzoRz/CK6+8gvfeew8333xz6fiJEycwc+ZM9Pb24tVXX8Wzzz6LZ555BosXLw56PwCg6qTx6jQfgcveXESos5R5V2rT+vAsfazkdUayXNd54sQJALKvIlpgGJch+SwkGSNHjhywvlQm9VerL12Tbwiztu+QFEot1U+WHTt2LHUsvl2zJCNqy/IsBro+6+JRH1yyLG35fn+k+tJ1hiTPdYbqT09PD+rq6tDQ0OCdmzHqx/79+3Ho0CHU1tZ6nfvrX/8a48aNwyuvvIKrr74ahw4dwm/91m9h9erV+OIXvwgA+NnPfoZJkyahs7MTl19+OV588UX84R/+Id577z00NjYCAFauXImFCxfi17/+NWpqaryvYSCoyRBCyBChp6en7C8+6RmIQ4cOAQAaGhoAAFu3bkVfXx9aW1tLdS644AJMmDABnZ2dAIDOzk5cdNFFpQEGAKZPn46enh7s2LEj5CVxkCGEkFDkNZeNHz8edXV1pb8lS5ao7fX39+POO+/EFVdcgU996lMAgK6uLtTU1KC+vr6sbmNjI7q6ukp14gNMdDw6FpKKRJdFNzDENsZ5tmS21M+TIcAiIw/aNVmPWc1fWe+f9R5I9SzXYJWbt//SsSwyfM22Icj6nmSpV3T9kFufV4IoiCPLeQCwd+/eMnNZZLYeiAULFuCtt97CT37yE+82KwVDmAkhJBBZB77ovNraWrNPpr29HWvXrsXGjRvxiU98olTe1NSE3t5eHDx4sEyb6e7uRlNTU6nO5s2by+RF0WdRnVBU1FyWNbIiHqJYXV3t5VjzbdPaJ6uMvNEkWshnnpBXa+io1v88IbuSfK1tS3/y9C3EvdLun+s5+sjK8l75vNNZ3msr0ffXFXYcss1KYn1Gea7n5MmTaG9vx/PPP48NGzZg4sSJZccvu+wyjBgxAuvXry+V7dq1C3v27EFLSwsAoKWlBW+++Sb27dtXqtPR0YHa2lpMnjw5510oh5oMIYQEIq8mY2HBggVYvXo1/vVf/xVjx44t+VDq6uowevRo1NXVYe7cubj77rvR0NCA2tpa3HHHHWhpacHll18OAJg2bRomT56MW2+9FUuXLkVXVxfuv/9+LFiwwGmi84WDDCGEnEI8+eSTAIBrrrmmrPzpp5/GV77yFQDA448/jurqasyaNQvHjh3D9OnT8f3vf79Ud9iwYVi7di1uv/12tLS0YMyYMZgzZw4efPDB4P0Nsk6mTKCn01tqPikj/v+oHyE3BhtIXggHvq8M7dojkuajgWQMH56eQ0T1ozUd8bK4rFGjRg14TJKRXA8iOfTjMqRNzrR1Mlp9qS1pzYrUt6SMo0ePqvIlGRHRamnXO501/5mLvCakgUyWWdGuM24itFx7iN8VidDrZM4888xMvxsnT57E4cOHM62TGepQkyGEkEBUwlx2qlHY9sshNI084YvaDCXPFsG+JPvpG9qb5+WT+u/jZJb+HUiGVC95zBJk4CIuwze7gO8W0Vr78RX8vs9Iq2+VFfJHyVeW9btllZs1I3doa0YI+vv7M2sypyvUZAghJBDUZNIUtv1yyJvmmnkWvRdDaLu0Jj85q8vjk9Hqu8osSP2WZvoaoZ9dUp7UR40QmoRLhq9PpqgFyEUgPU/p+xsvy/oOhNQgQ8FBJg3TyhBCCCkMmssIISQQ1GTS5BpkrI5xS70Q29NqD6qofGK+WFZdA2FzXuXJpSTJ0Jz7QyVXlzUAITLr5XF++wYS+JrahvIPkOU7ZzWhJc9zyff9/akEHGTSUJMhhJBAcJBJE2SQ8XWW+m7+pYXMhiJE6HVIbSlEf3xDjbV+WOtZw5V95Wsbk2n9CKFlu4IrLNfqcuQX8SOT596G6I923wbaRtnSdoiFq0XBEOY01GQIISQQ1GTSMLqMEEJIYZwSmoyv+gzo+7/naSskIXM05VmzFGKdjE9/8rTp6+QNsU7GZU6y3PsicvqFIsQ6MO15uhz5QyUoJwTUZNKcEoMMIYScKpzOA0YWcq34J4QQAtTU1KCpqam0t0sWmpqaUFNTE7BXQwNzqn9CCCEDc/ToUfT29mY+v6amprTFxukEBxlCCCGFwegyQgghhcFBhhBCSGFwkCGEEFIYHGQIIYQUBgcZQgghhcFBhhBCSGFwkCGEEFIY/w+TGnbBS9TuSQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 500x400 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(5,4))\n",
    "plt.pcolormesh(axial_slice, cmap='Greys_r')\n",
    "plt.axis('off')\n",
    "plt.colorbar()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can see that it has been blurred, and that the values are much too large. How can we fix this?"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Example**: Filtered Back Projection. In this case the image estimate is given by\n",
    "\n",
    "$$\\hat{f}_i = \\frac{\\pi}{N_{\\text{proj}}} \\sum_i H_{ij} \\left( \\mathcal{F}^{-1}(|\\omega|\\mathcal{F}(g)) \\right)_i $$\n",
    "\n",
    "where the term in brackets involves applying a 1D convolution (in this case, multiplication in Fourier space with the Ramp filter) to the image along the $r$ axis."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "freq_fft = torch.fft.fftfreq(image.shape[-2]) \n",
    "filter = torch.abs(freq_fft).reshape((-1,1))\n",
    "image_fft = torch.fft.fft(image, axis=-2)\n",
    "image_fft = image_fft* filter\n",
    "image_filtered = torch.fft.ifft(image_fft, axis=-2).real"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can back project and normalize"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "object_fbp = system_matrix.backward(image_filtered) *np.pi / len(image_meta.angles)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7f592105af40>"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZcAAAFICAYAAACcDrP3AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy88F64QAAAACXBIWXMAAA9hAAAPYQGoP6dpAABi00lEQVR4nO2df5BeZXX4Dwkm4deCELKBzDpBrQIDBgwlE21HZro1KGXKtLUUqVBEHCyxjLGt0UJCy0BsgTQUaTMiGeoMFGydMk7DxKFLU2tJjYZ2prb8GARKBtyQiLjkl4Ek3z/8nr3nPft5z/O8u3dDiOczk+Hlvvd97nOfe+/e8/sctn///v2SJEmSJC0y5Y2eQJIkSXLokS+XJEmSpHXy5ZIkSZK0Tr5ckiRJktbJl0uSJEnSOvlySZIkSVonXy5JkiRJ6+TLJUmSJGmdw9/oCSRJkhwK7N69W/bs2TPu30+bNk1mzJjR4ozeWKpfLr/9278tIiKHHXbY6DZK7tdt9jv9jW7bt2/f6HdvectbRORnC6tMmTKlY397wWiMvXv3jtlPx9CLNXXq1PD89Hsdyx5Dt9lzonPXY+p/7bh6fj/96U9Hv9u5c+eYufkxjjjiiNHv9Pxef/310W2vvvqqiHRel1mzZomIyPTp0zv2sedi56Hn+d3vfnd027XXXisiIv/7v/875jx13rt27RrdtmPHDhER6evrExGRE044Ycwx7bx1/bZv3z66zV8jey2OPvroMed5+OGHj/ndj370IxERGRkZERGRo446avQ7XcsjjzxyzDmdfvrpIiJy++23j277xV/8RRHpvJ66pvaYxxxzjIg0a/rSSy+NOU/dx87b3vO6lnot6P6289Z52Hte96MxFLt++lnPhe5bew0I/X737t1jjqnnR/e3nYfup9vsGHp+r732Wtcx6G+NHd9/Z3/zta99LTy/Wnbv3i2nnHKKDA8Pj3uM2bNny7PPPnvIvGBSc0mSJJkge/bskeHhYdm8efOogNULIyMjMjAwIHv27Pn5fblYKcFrGP5ztE1RadZKFarN6Pj6/yKNBGOlLMVKSF5CtxoASTy6vx1Dv9djWclbt1kpi7QqL43ReZLk5efVbd66jTQunaM9Jz0HmrfFb6NrTJKxnhNJiXb9dE52P6+50PilMey94seo0aotdJ66jaRxr+na/UrXzH9H52S3kXainyOt2u7v14+en9Ja6fOl49oxSCPSMex18s+I1VLs9fZEf1dIQ6P7qW2OOeaYDi21lkOxxGM69JMkSZLWyZdLkiRJ0jo9m8XI3EGmDavy62dS+VVVJXMH/T+ZWMh8pvvpsa2q7U0ndn+rnnqTBqn3dhup3d58QWNEDkkysVhoDHuu/ndkOonMBTQPMjN4Exw5VO28yBTjTZ00bzuGOrXp/iMzEM1bnesEXTMy/yg6t8jMacejNSWiey0y+5XuV38Nas3b9hroOdAzSM57CsLwgQj22aZr4M+TTHE071IwUhvs379/XGOnWSxJkiR5w/nWt74lF154oZx88sly2GGHyYMPPlj8zfr16+W9732vTJ8+Xd75znfKPffcM6lznJBDPwr3I4lQJRKSnij8UCVTkjTtPGhclaBUyiJNh+Zh9/OO1FJ4I0nLXhqLJOUSpGGQ09lrLvaYkXMzklZprcghTWtLjm5aD7++9jvdPwpc8J/9/iS9R4EL/jh2DLqOuu6lYIk2HMtRYAZpDAo5uvW/pFmShhtZDUqBNrpuNG997mm9Lfq9jkXXotfAorY4UJrLjh07ZN68efLxj39cfuM3fqO4/7PPPisXXHCBXH311XLvvffK0NCQfOITn5CTTjpJFi1a1PN8a8hQ5CRJkjcZH/rQh+RDH/pQ9f6rV6+WU045RW677TYRETnttNPk29/+tvzlX/7lpL1c0iyWJElykDAyMtLxzyY7T4QNGzbI4OBgx7ZFixbJhg0bWhmfqNZcIocjqXSkCitWrY5MLKom20xmMktE5pko8572t6iZgExxlBdDZhrdzzs+7bnYufm1igIMuu2n6xs5PEtqOAU4RMf045f202trs879PWZ/pxngJaesP/5EzjMaPwoyKY1FeVNKFMhBWftR3grdh3QveDOqnW8pj8YHUJDZzd7TNK7eC/5Z6Xbu3hxG50nXIPob1hYTNYsNDAx0bF++fLnccMMNE57X8PCw9Pf3d2zr7++XkZER2bVrV0ceYFukWSxJkuQgwWf4q9/5zUjPmkvJoe9DgEUa6SNy8kfb6DhWAvPOe5GxjmWSwMjpTMfyodTdqAk7LY1BYZPdxqTfiTTroMe04aEkeavkGElcJSe0D9ElDcOOr8e02cz++PY8SYuNQuGjkGEKoaY56PrZY9L66XpTKHK343SjppqFnRtpwjVhzXaMWq2NnhH/nEXVEuwxKQufqnXoPOzfDq+x2L810fNT+tvVBhPVXPr6+sZVPqbE7NmzZcuWLR3btmzZIn19fZOitYikzyVJkuSQZ+HChTI0NNSx7eGHH5aFCxdO2jFbN4tFkj9pMBQC7J1YUfikHZek99oEsajGFNmrqa5apN2RlEXSrNc2Slob2eh9mCqtVelcaog0NQtpPTWaHEnZpXH98XuVUMk3Qn6KKBGYwp+jBFgLafLkp6gJN6ZzsUQVj6PnIUpatM82+Tu9lmI/637WFERJl/Qs+Xm/URyoUOTt27fL008/Pfr/zz77rPzXf/2XHH/88fK2t71NPv/5z8sLL7wgX/3qV0VE5Oqrr5YvfelL8sd//Mfy8Y9/XB555BH52te+JmvXru15rrWk5pIkSfIm43vf+56cffbZcvbZZ4uIyJIlS+Tss8+WZcuWiYjID3/4Q3n++edH9z/llFNk7dq18vDDD8u8efPktttuk6985SuTFoYskg79JEmSNx3nnXdeqO1Q9v15550n//mf/zmJs+pkQi8XcppGJ0zmJcJn7FpqqwL4UE0KV6QwTppH5KClMOLaGkYU2lljFrPnqc5mO4/IsUxr5UNB6bdkVozqgpXMkHp+kXOdmsrRMWh834CqtD8FotBaRY5iH0hh50HBLDQ3crLXZvTX1KWjkPwoEKBkevRmX3ueUQl9qk9Ga0Vl+0sNzPy51AbTtEHWFmtIs1iSJEnSOq2YxawkQRKpfytH0iphpRxK+FJsyKjOQ6VyCgooSeNeiixJTFSHKwpx9vvYMXz9JAs5VElrIwdz1IqZpL1Iu7LUNOmi6073jkKh5aVx/VgUvktaAUm3JDXTvP06k+RN95UdNwrHLdUI83OLwohJM+u1WRjNjZ57astN9df0udUxokZsBIVGR38nJpPUXBpSc0mSJElaJ18uSZIkSetUm8WijGAq925NHN60UnJWeucg7VNSe71zkExgpf7q3sxQKvNPpiyv1pM5w67Vzp07O/aP6o7Z8e13midEJgil1MApMjnRveAd6KWaVLp/VF6dHPo2yIMaZukYNH6UK0NmUQrQoLn5dbZ5GmT2pdwunysT1VyzRLXqyJRJJi0yp0VBLPaY/p6014ec9wqtc5QjFQWn2N9RMEuvOT5JO2QocpIkSUukz6Wh+uVCkmwk0VFbYZJCqBaQdzBayZ4kjqhGEoU3Ru2Io2x2kgSp7W6k/dD6kQNTtY9S2DatkXcOlzQXclx7DSRy9ouM1VxIKp9IbaeoeRVpdzof1QTtftH9Vwp08McWac6VrkV0/ewYer2pKoCfq92PGp/puLbSBQV++PMrVQ6n0Gw9vmos9nmgNAB/bJqHpUZjoWOWQtAPxT/mBxupuSRJkrREai4N1S+XyC4b1VmyRGGZkVRLPVMsJNVGrX6pR0TkDyIpkdZDe46QZEeSMUlZvsKzlSYpCZXs1V7TouTIUuXeKKExqjdFvUqiMN8o4bS0vw+5Fmk0lkjyJ+2RrjHNK1pT0gDJb0N9S6iCtR+/VAFZz0vvQ9KII+000hLseHab11gird0SaXSl/b11hPpDRa2s/ffJ5JArnCRJkrROmsWSJElaIs1iDT079EtON1VfKSuXoJpH0bHUzEBmN1KF6dhRG1gymdB8SK3XY6lZQqQJKaVMYzLnqVmETD4UWECmRiVyZEZtAexncvZGpkwKH49qO0Vl9clkVgos8GYoCsKg/SOzUckR7LeV2krTtdXj6vUvZcaTWUfvOz0/anIWObpr64hFz3bJHEVmwpry/nQ/6bFL7TdqAwWSdknNJUmSpCVSc2nouc1x7XfWqe2TqUqSdxRaXOsY93MiDYDGoOS/yElot1GiX01V3EhKtetCSWkzZszoOLb9LWkM5JCmxlNRMyo9PgURkOQY1Unr9b6ielzR/UTXk44RnW+pKZVfZwpFttqsD12mMUqSN2kRet/p8aPqxSJj7/lSqHNUn4yIAhFoP4XWhWq5UWh0zbxqvk8mTmouSZIkLZGaS0NGiyVJkiStMyGzWO3b1jvhKRPYjqX7R45jqoMVmRLIjFHbtIzMbhTDr/Om5khRFjKZraKe6+Rcj+qZWSiAguYRmYv0/KzDOFrvaPzaPBfav+ZYlPsUmf8ob4SabkX5UNT+gExJZM6LKlbQWlFgCznN6X6N8luIyMxK/0/PamT2o2AgMn1G5rCINIUdWNIsliRJ0hJpFmuY0MulV23GayR2fyutRNJ7lG1bqh7g501hwdTQikI1KShA96OK0KQtkUYU1WOKpHGSDkk7IU1EHcFReHLJ0R2Fk5IG1evD1Ou40RpEmos9jq4LaYWkPZKkTkEsdCx/r5WCRyis2rfGjiom23Gj+7tUS9CPFQXVlPajqhS1aQCe1FLeeFJzSZIkaYnUXBrSoZ8kSZK0Ts+FK6Ne1iJ1Wc2lQpTe0VgqMheZueg4UcazNYF4B3DJ7FZrFonm3W3+3b6LzIRRxjXlz5BzuJRr4iFTX+R4jYIwStIcOa6j+zMyNSpRwVU7BgWlREUvS9csOvdoXDsPbxaLsvFpDGvOjUzNUV5Wad50fj5wJwos8b/ttn+JLFw5+aRZLEmSpCXSLNbQ88uFJIgaSazb/jVSVq1TkcaodSpSEIEPNqCmZeTsjepPUcBA5JgmJysRaYOl9Y7Ck6NjUQM2alCmVQRKpf+jsNaoIgI17qL6WqSVeCjsuFQHy683ZZhbyLkeBU7QfRJprDQWhfTqtqjNtSX6A1jzvHXb5n9b0riUKDx9PPu1Qb5cGlI3TJIkSVqnWnMhu3Kv9ZuUkh3XJ5SVfB01oa4kDZEWYfGaAkkXJenTrxGFDEfaUq20FWkipDGU/A+RlE81tPy8qfJwryHJJY2YxvUSeqlFtv9dycdEvg4lSpQt4f0edJ9YovupJuxYhDWWbsexRJpFKaw+0rhqw8drNBGq8F36O9UGqbk0pOaSJEmStE469JMkSVrkUNRCxkP1y4XMBpFaHamdkZovMtZZXhsOSWNQ6G0UWhxlqZeOGZlUKAAgyj6mEOZaM4Nft/GYGfx60LlT6KpmtVuz2hFHHNExL/s9rYcfs9v+ula2xYEvO2/n+NOf/nTMcSKzTmR+IXMohVzT+kU1vWqaetnjR1UsSvdmFK5PRM9tbf09WlP9TEEHUaADQa0CqDlbvgAmn9RckiRJWiJ9Lg09v1yophI1R4qaeUXagR2XpFX/Ozu+HUN/o9KqJZIwo2S3klTrtSX7fVRhl8bVNbUVlkk6JKnW/7ak/UTnFzVus+vo9x9PEIGXnMfjXNfPUZM20gAiLYXuZVoPalEcaS7020iLpMTN6BrUJgpGocKRhiHS3GPTp0/vesza4JEoSKHbb7uNb58banU+WUmU+XJpSId+kiRJS+jLZTz/euXOO++UuXPnyowZM2TBggWycePGcP9Vq1bJu9/9bjniiCNkYGBAPvOZz3R0SG2bfLkkSZK8yXjggQdkyZIlsnz5cnnsscdk3rx5smjRInnppZdw//vuu0+WLl0qy5cvl8cff1zuvvtueeCBB+QLX/jCpM2x59piZMawpieK04+af0XHIkc0OfkjU4aq62QeK5X092a8kspfU9Or1OjLz62UG1JjbiutVYTPXxFhp6k/VsksRmsZZbrTMaNjkblVP9tzsSbdGqI1pUoHNc57/9lDpsbIrEP3U9Tegc6J0HWzZiZvDqvNlrfXoLZ5nx+DmovpWHaOZKo/ENn6k8nKlSvlqquukiuuuEJERFavXi1r166VNWvWyNKlS8fs/+ijj8r73/9++ehHPyoiInPnzpVLLrlEvvOd70zaHFNzSZIkaYmJmsVGRkY6/pFQvGfPHtm0aZMMDg6ObpsyZYoMDg7Khg0bcF7ve9/7ZNOmTaOms2eeeUYeeugh+fCHPzwJq/AzetZcSlnTpG14J3KpQqw/ppVCdH9yZEbSpK0xRe1aa6Q3khwpbJIkMKo8HGlh5NyMpM+o2gDNu9QmOoLm5seiAA2rJZDz22suVEeMHMFRaHapAne3+dvPparI/jvKBCepOWoIRuti7x2dU6T9UC0yCk+uDe+nem2RtkHagT7L9BxE0N8ManJH1gBa78nSXCbq0B8YGOjYvnz5crnhhhs6tm3btk327t0r/f39Hdv7+/vliSeewPE/+tGPyrZt2+SXfumXZP/+/fL666/L1VdffXCYxZIkSZLJZfPmzdLX1zf6/2p2nCjr16+Xm2++Wf76r/9aFixYIE8//bRce+21cuONN8r111/fyjE8Pb9corpZFivZecmk1zc7JdNZosSzbnPwRAliUaVXOnfyN+n+27dvH/1OpTiq1xYlaZb8FVFyHGl+/jv/WYS1n6jPCGkM1JLXaqVeGyBNh0Kzo2q+5HMhrSPS3ijRM/KFkTZbqiYeaZvU/pcSU/1vS1WU/TwiH5oIV8P2lEKoSdusIfIblu7NA9m7ZaKaS19fX8fLhZg5c6ZMnTpVtmzZ0rF9y5YtMnv2bPzN9ddfLx/72MfkE5/4hIiInHnmmbJjxw755Cc/KX/yJ38yKWuUPpckSZI3EdOmTZP58+fL0NDQ6LZ9+/bJ0NCQLFy4EH+zc+fOMS8QXyC4bdIsliRJ0hIHKolyyZIlcvnll8s555wj5557rqxatUp27NgxGj122WWXyZw5c2TFihUiInLhhRfKypUr5eyzzx41i11//fVy4YUXVkfp9Ur1y0XNEqX2rqU6TCLlWmQ+fJNU7dp2tGSCoDd2lKEdlZ8vmUx83Sk1Z4iwc9PPg8YqmcW8yYYCIuxaUTiuv2YU6hy1CoiCGkSadSAzF6Hrt2vXrjHfUctmOs+oKRpto7B6HZeyvqP6WhTgELUlsGul2+y9E4Ub0z1fcw+XQu2jbPmovh/V9IpM2CWTlg/NjszK3cadrD+oB+rlcvHFF8vWrVtl2bJlMjw8LGeddZasW7du1Mn//PPPd5zvddddJ4cddphcd9118sILL8iJJ54oF154odx00009z7WW1FySJEnehCxevFgWL16M361fv77j/w8//HBZvny5LF++/ADM7P8fs3ZHkhIiByk5dCmJkqRgn2hFUmLklLXjRcldJWe5l5rs/jQnjeygcamSMOETyuy89ZilUFAvwZLmQtI1SeiUBEh1rbzT3n7Xa0toPwe7P4W9R+HGNMeoVhhpvyRl0/Wne4iuj/6WtN5Ig7L3DjncveZEyYVUybpGC7L7R9YJShImhzslspIGRUnC/jqSdkWBRzTftjlQmsubgdRckiRJWiJfLg0ZLZYkSZK0Ts9msSgLXqQuY5wyai26PzlUKeZfM3WtuajbmN22kTPWx/VThjQ5K23ik28bUFLhvWOUsvFtcyx18kb5GdYJXlPPrNv5Kb02cKLxo3EJChSInNQKmf+icUuOXqpdpdeAzDQ6nr0GNXW+oux9u42CQbTSLeWaUVsKMm9HrRks3mxFzyWtB5nnoqCAKHeoFEARmXbbJjWXhtRckiRJktbpuc0xObxJmyHJgSQkCt/1kho5We02ldRskTeVykj7IYmnJvS31qFqz887V63G4MNP7f5RrTN7nqqtkTNeNSgrNVO76iiYQb+zkroek64ZhayT5E3SpNcaaG1LVZF9zTJaW6v5+XMnzZIyzC26zqTJ62+pAGEpI16JwrppjahygXea0zGj+9BCTvuoZtiMGTPGbIsy+aNKESJj780oXNmeVxSu3xapuTSkQz9JkqQl8uXSMCGfS5Q0WPItKCrpRL0qSklS+r3tqqZzUomKEtZKyYg6tyiZzvp5yB9Ekryfd9Q+2X6nY5EPxeJ/SxWhKQw2CjstrZUP5aV+GiV/SeTDIc2PNItIc1GiCsV0f5OP0K5pr4mE3kdTS6kHivcfUftx8o/6NsN2bqXqzP7+s+uiGkupqrTXLOiaURi7Emlv9vsDobkkDelzSZIkSVonzWJJkiQtkWaxhp6bhVFmMoWwRs5ha/6gMF+vspKpgkxltvaSOm1J/aXghCjTXfe3Kn9tVrMvcU/nSetBmfT6HWW/W/y5lNo5K1FoL7URoHOJGqZZ6LpEZrHIjEFmF2pnQM3ZPLQGtU52MtORGTJqd0BjkHmOgkH8c0bXzELmsG77dNvPP+/2GYlCv0tmKz9vCgCIzF10DSy9treuJV8uDam5JEmStES+XBqqXy5HHnmkiHRqLtS6N0qciqr6UphqTbMuGl+kCQ/V+ZYuHrU+pnBPhbQIkpoUHdfuT+tHUqf/jqQumjeNRRpRVHeN9o9qNUWJhJGmY49F8yJnchQsEVXMpcTXSJMiqTmqhUf3N4Wx07EoSZOadEVj0P60fj6YhjTi0nX392QpdDlKPVBKjdX834dS8qdqU/b8JktzSRpSc0mSJGmJ1FwaMlosSZIkaZ2eNReKnSd1mpxs5GSNzEBkkolMZTQuOXHJ8U7j6W81u5rUe7sedKzIXEROWX9+ZJKxY6hpw2ZB+3wbWj8LmXGiEvA6b3vN1BxKx9T9S+aLKO+D8hdo/aLGdP562v2jTHCaB60pmUo194ryiizeXETrF5mORcY+j9QiwlYn0HP3FQb8Z78tymGrzTmx+Oe8VFvMVyCgSgSlIKPJ5FDUQsZDai5JkiRJ6/SsuZCDj5ymJI2TU1Hf8lai8hIGZWBbIk2EMo171WZUKyBJ0G4jqcmPT+tnz0nXIQrlJc2oNiS1pqZXt2NF+BpkJQc2Odx9uCxlgpNT245Bx4+Izo/aYVPdLj8Pkrxr2+pS2Hv0TJGE7i0F3eYROeEjLYVqkNXWBVMopUGJ7mWRRtMia0BUIcQyWdpF+lwa0qGfJEnSEvlyaah+uaikUeqdEIUZ17SZtZ8pLJNa/UbJcXRMkuIi6ZDOjTQi8il5SpV+ve2axiq13VWJlFrK0nlGkh1pXFQfyv+uNnycqi0rFI5L15jWQ6GadeQvicYlSZrWlCpIR6Hftf1cojlaSItVKInSU3qOew3X9/uIdFooFK/xl7TT2pQGP++S76cN8uXSkD6XJEmSpHXSLJYkSdISqbk09PxyKdUEiupDRepsFKZaymrX30bNvEoXj4ITvFmETCxRGC8dn9aKnNqRA9uipiQyISnU7pbCx8lkF4UiW2pDxD3R/VJqoxyZDGtNQ95kQgEddj0pLN2H/pKpj+qf0Tx6rVlnqQnrp9BlMhf6efnfdpsHmYnttqhpGT0PvQabRIEFpdSDNsiXS0NqLkmSJC2RL5eG6pcLSaaRNB45WcmRTq1hSXqmcN+oqmuUDGihcWvaqJJEShpATSMsOwY1maJkOmqe5iU1K0lTlVmVrqNaXtSsyc5jx44dXc8zSrCj+ld+rnYepRBTv/9RRx3VdXy7H10naudMQSa6H9UdizTcqPkX3ROlJEe/Xyl02e9PIe5RyHC3uSk1FdJpP6qFF1U5ju79bnPLZmGTT2ouSZIkLZGaS0NGiyVJkiStU625UEZ1TX6EhfIAomZKOhaZR0r1hxQyu5G5SM06tkaXz2Ams0TUQMkSOZjtd75ZlZ0jmRRqttl10fHtuNpkLQqIIBNLTe6EHSMynfhj+O98GwE7bgTVYaOgFGpnoL+1tcgor8ObJsm0Sg5mMjnR2lJuDZ2fN6mWzNX+mHY+WhPNrkeU0R89B2SKixz61BaAzqWU86REZfvbJjWXhjSLJUmStES+XBp61lysMy3Kno3CJinD147r6yvZ8VV6ImmcapZRprseX2sU2XGtlOqlSNIO6JgWr8WQxmUz073kZcfXudn9I42Iwo5JOlSJ1e5H2dJ+TpG2adGsbHtuut41wRh2fzsGZXtHTm3SfnyFZ3tOui7UzI0aT0UZ8vb+03lTvTsK1Y3OkypPkCZPWrJ+T8+gNgcsVRKOqlhE9x8Fj0TVx6PK3iVNqrbCQRvky6UhNZckSZKWyJdLQ/XLhSRNSpIiO7KXfkp1p7wERlIiaTNkX49CEykc0vpcogtO35E9nsJT/XcUbhyF9EY+GpGxFYGtphNJgtSGmHwROoaVqL0WS5WKa2tXKSXflX6mHkM0R91m18NXt47WwM679noqJS3PJ1uWNGLSQP1zQ35GwmtNInFyMD1L5EuhMHldU3tdvOZUCiP251LSFOk+mazeLvlyachosSRJkjchd955p8ydO1dmzJghCxYskI0bN4b7v/LKK3LNNdfISSedJNOnT5d3vetd8tBDD03a/NIsliRJ0hIHSnN54IEHZMmSJbJ69WpZsGCBrFq1ShYtWiRPPvmkzJo1a8z+e/bskV/91V+VWbNmyT/8wz/InDlz5P/+7//kuOOO63mutVS/XChbWdVea84gVVihrPOakECrzupvqX1tVKMrKsdv50vqdGSionBZO65v/lXK1Pbzp7BWwpojvQO4VO8rcoiSGY3K9us5+EZOdoxS6GiUoU9jULCGblPncG01gyh0PmrrK9KcOznqowZyFM5MwSmEbyttx6CggMikFa1jyZzsnwMKMCDTHQWURNCzXVpTpba+WxscqJfLypUr5aqrrpIrrrhCRERWr14ta9eulTVr1sjSpUvH7L9mzRp5+eWX5dFHHx29P+bOndvzPHshzWJJkiQHCSMjIx3/bPSqsmfPHtm0aZMMDg6ObpsyZYoMDg7Khg0bcNxvfOMbsnDhQrnmmmukv79fzjjjDLn55puLZbEmQrXmEjUGK0nvup++MW0IsGJPUiUZlc527tw5+p1KStQsjBy71HrWjyUS10mLWttaKNzTh6mSozZyWpI0WSuxERRSTiGjUfMqkmqjSswUiKDXj2po+TnYuVkp17fYtVBzsSgsPUoMpbpWUQIpObWpdh4FPdTeJ3TP0/kpkSWBnpGorp7F/5Ya+1GzMAqT1nuCnmOL16rpb02pFt5kMVHNZWBgoGP78uXL5YYbbujYtm3bNtm7d6/09/d3bO/v75cnnngCx3/mmWfkkUcekUsvvVQeeughefrpp+X3f//35bXXXpPly5f3PN8a0ueSJEnSEhN9uWzevFn6+vpGt5MgPh727dsns2bNki9/+csydepUmT9/vrzwwgtyyy235MslSZLkUKevr6/j5ULMnDlTpk6dKlu2bOnYvmXLFpk9ezb+5qSTTpK3vOUtHZrgaaedJsPDw7Jnz54xLcbboOcM/ZJTjPJhvNmilMmuKq3WvLK/J0emd6jauUXOv5KZy5t/SvkokcOTnKyUJ+QdujTvUj/xKF+EzINkVvKmh1KtK2+6sfNRuzHNI6oVFvWYt/O20l1kMiHzX+TQp3WpyXMptWEgZ70/FpnuotpvlqieGR2bzG5RGX6ad2QOt+gzbf+Y+bWMqnsQlL9CuUk077Y5EA79adOmyfz582VoaEguuugiEfnZug8NDcnixYvxN+9///vlvvvuk3379o2ux1NPPSUnnXTSpLxYRNKhnyRJ0ir6gunlX68sWbJE7rrrLvnbv/1befzxx+VTn/qU7NixYzR67LLLLpPPf/7zo/t/6lOfkpdfflmuvfZaeeqpp2Tt2rVy8803yzXXXNPaeXt6ztAvSRVRvbGolazFj2slG82gL83DS6R2H1/LSISl2mheqnVEzZfssah2FWVZ+3BZOqadt8/GL0HaFdXtUslOj2VVdXKuR+2C6d6p0YRpPjRvqkVFDnoKuPD1rErj99oelxzdNZUFLKR9Ux0z/5xF1gORsfdYpEHbedCzrc9lFLZvP9dWt6bnIHLQR9eHrm3bHKhQ5Isvvli2bt0qy5Ytk+HhYTnrrLNk3bp1o07+559/vmN9BgYG5Jvf/KZ85jOfkfe85z0yZ84cufbaa+Vzn/tcz3OtJX0uSZIkLXEgy78sXry4qxls/fr1Y7YtXLhQ/uM//qPn44yX6pcL+QCoQmy0SCRtReG4NGaNlOPn5H+n86AwThpPv7PSJ1UGprDTGkkq8kWUNAw6jteqSomb+r21vXpfAWlLke+H+tDQeVILZn8epTHsufi5RQmWdlu0BqWERt8jh7Rfekai0FjyuVBINGknFG5O96tvb03hwSThR/c3aZGWUiJtt/1LKRDRmHTM0jUdL1lbrCF9LkmSJEnr5MslSZIkaZ1qs5hmyZMDtmQGihpPqToYhQaSWazU9lTVXppj1HCMaqf539n5UiMnqqWk41qzi4bQ1mboU6Y7OXu10RM51ykDPILMRTouObqjdroULh0Fd5Sy6yk82R/fzpGuQU3WtjWjkSlT0XUn052dRxT4QfjgCpHmvotqp9nvjj766DFzi0xaUTviknmu2/h2XLpmRFSrkNZWoXu+FODQBmkWa0iHfpIkSUvky6Wh+uWydetWEeHquFYaUmmcGnyRFBLVMKp1SJMU7MMxSSojLYwcwKRx0RgUzuq/s/PevXv3mDF0vpHTkhpmEbR+5ARXKOyZviN8FWpyspOWQnNTSGsrOWe9xF06T4USWv2x7bhUh420agoLrtEAoqRYkSYZsddkTtqfnP107Cggxx9HhLU8OkakgVAtMm/toEAbq+VRgMNk/THPl0tDai5JkiQtkS+XhnToJ0mSJK1TrbkcddRRIsLls+026j/gITU52kYmNutcp9o4PsO45OCLHLvkCCZTT2RaIROiQjH5lP9DwQ9Rzo7+Vk0oNB8RNgl6UyYFP9h5+GtAjuOoVpfdj76jcyezS1T6n+4F32CrNp/CbvOmG8o5ot+SiZlq1pEpKcroj2qRUSO7SGou5aR5M2TJtLZr1y4RYVNWbXOxXp3xUZuOtknNpSHNYkmSJC2RL5eG6pfLzJkzRYQlQWpZSg64KOSWWuaqFkSOetJ0SMOpdehHLW0pQ5+O6Z3a9reRJmKhFsKKSl5UWdnO21eJppa89tx1DJJSKRSZnOR+HhS0QedU69CPfks1yKJAh6gSgf0daR167hqMIdLcp9TSmDRt/Z4an/Xa3C5qOWyvU9Q+mYIlojB2Ctenygx0TGqU5udtz0mvFfU1idbArjc9q1F2/0TIl0tD+lySJEmS1smXS5IkSdI61WaxE044QUTYLEYmrciRWlIBo3L55OCLSpiT484XGizNLTLJUIY0OVLJUUvH9Kp+qeCmP46dBwVLRCYZCtZQByzlqJDpgUrpR478WrMYma3oWLofFfmkNdXzo+AUMuNSaXmfc0LX2K4Vmdui8vc0HzLV+m1UnaBkmvTHosKpkYmKIEd6ZJIutTXwcyNzFzW3S7PYgSUd+kmSJC2xb9++0N8X/e5Qo/rlQq1qKbw2asJD0lMkqZG0GgUAUJgqSUMqeZFkF82Dwo5Ltb+oLpT/ztYno3bIHis51gQn2POka6bH0nBzkbHBA7TeFt2mmoAdiwIdaCy/9qRFkIZmf6c18FSboGx5CsLQ86V5l8KqveRP199e4yhcn7RqnVspjN1fW9IKSauOWjiQZm7xz6NdW621RloEhTNH86D9a+rC0VxF6lImxkNqLg2puSRJkrREvlwaql8uIyMjIsKSTEm68VpBKUzZSyYkbZXqTvmwyVKrVQpZ9vOlUM2SdO33p7Ui/5Sfg92PkihtaKyuF9n2SdrzrX7tOajUbH9HkrRf51q/msVra7WNp+x3Xlou1ejyDbPsProuJOWSGYPWnZ4HWj9/z5B/gDRtup60LrQeUSVreraj5zfSoGtDnGubCfr50hrQWtG2tsmXS0NGiyVJkiStk2axJEmSlti/f/+4tKJDUXOpfrm8/PLLIlJ26FN4oDdzReWz/Wf7Ozt+qY98tB9lUkfnUjtvql3lxydnOIX00hjUWoDG9Q50a6bRUFp1fNvxrNNZzWw+219E5NVXXxWRzsCCqHYaheNG5en9PiLs9NXzomus8yUz1zHHHDO6zWfc2/Ok0vVqdiMzZM21sPONHPRksrXrEQXORCG6UeUEMl/Ze4fuJyUymZX289ednsto3hRsUmtObps0izWk5pIkSdIS+XJp6PnlYiW7yEkdSR8lycRLgNRKmAIACFJRyUkYJd35edljRq2b7X5R29jIKUvhpCQdRklj5Ni1WoeG35588smj27Qt7kknnSQindJ+X1+fiHCiJDmfqZ6UQveTn6sIV77Wc6D7KUrmVO1DpDl3PV+7Bscdd5yIcAVuW2naJ6HSNbbzjxqlRcEgdtwazc9qGBoiHrX6pfuKkjSj9sVE1FxMJK4IHWkdkdWjFPyQSZSTT2ouSZIkLZEvl4aMFkuSJElap1pz+dd//VcRKWfX15QCL2XnRg7P2lh4D82RpIWoeRWZu8hUZonMHWQu8lnnVP0gClKwv1VTjDU9Ucb9W9/6VhEReeKJJ0a3zZ07V0REfvCDH4iIyIknnjj6nY6npiSRxllOjcmi3BfKoyCsA92PS2ZTNQNZ89/27dtFpNNEtXXrVhFp1nnDhg2j3/34xz/u+K+dr3Xa+3UuNQsjk5PeC5RvQ20m/Fh2bv53djy7Lbq/o5L/kemYxqPnPXq26Vml5z3avxTc02t2fy1Z/qUhzWJJkiQtkWaxhuqXy2233TZmW6/1wJRazaVXSSmaW0lzobDGqOZRbaVkfyxy4lL4KzkrCToXXwfOfqfnYqVmamF93XXXdexnpf3IWV5TK85+jpqn1TbHot+Q9kvz1rVXze/tb3/76Heq5e3YsWN0G7W39tqG1ZaikHiL1+qpYVYpACAK0VWi4BG690tpAP68xvN8RpUw2vh70uvcJkK+XBpSc0mSJGmJfLk0VL9cNm7cWLVfJLGSjZTCBGtqB5EUUrLpR/tTb5pIAozs5iTVUn0j8kF5P0wpjJPWz28j7cduU0nbhuiqf4L655AGEM0xShas9RlE2+j8Ih8XaS66bdOmTaPfaaKprS1G5+JbU0fXotsc/fWmHjV075CWTFpb5K8jLSUKnSdKbZn9HOm3kUWhti1yzVwnk/S5NGS0WJIkSdI6aRZLkiRpiTSLNVS/XCJHbUmlixx2ijW76OdatTdyBKrZio5J49I8yOxCJcwJ78gnM0OteUnPpZRd7OtOURUBe05qDrO/0xpkinXo1zRWq91GJkT6/5pacd1+67fZa+Z/e/zxx49+1qoEtg4brb2vF0fhx5badtV+f3sda0y29js17VFtO2oNXKo84ecbVbiIgjbsNjpPquun+PDtbtSa1dvgQL5c7rzzTrnllltkeHhY5s2bJ3fccYece+65xd/df//9cskll8iv//qvy4MPPtjzcWtJs1iSJElL6MtlPP964YEHHpAlS5bI8uXL5bHHHpN58+bJokWL5KWXXgp/99xzz8kf/uEfyi//8i9P5DSrqNZcVHorvfEjp1xU94eqC6sESK1wrVRWU8+MpBeSYK3UpI5ucq5TAiQlwHmth+Ztx/AaETlZSVqlas7U0IxCNalybyR9R1pHVNW3FEYchZhSAACFd0daMs3Nn6ddA10Xqi5sr4FqfiRlk9ZL19FrIhRAUQot9s8BPYM2aMPPwx5T18WuTxRgQUEH9LxH9x89g1FDQvquNvR7ssxQB0pzWblypVx11VVyxRVXiIjI6tWrZe3atbJmzRpZunQp/mbv3r1y6aWXyp/+6Z/Kv/3bv8krr7zS8zx7ITWXJEmSltBosfH8E/lZx1/7j7qg7tmzRzZt2iSDg4Oj26ZMmSKDg4MdFSY8f/ZnfyazZs2SK6+8sv0TB/LlkiRJcpAwMDAgxx577Oi/FStWjNln27ZtsnfvXunv7+/Y3t/fL8PDwzjut7/9bbn77rvlrrvumpR5E9VmMa3VRLkh5CQkIpWVTENUglshc4c1DXnzUtSvvNs8vEmtZKKqORaZcMjsp0QO225j+JwdyruwkDkiOnciqojQa0sB+l1kWqs1zxHezEUmM5Ie7dy07hkFbUQZ4xbvQLdZ/vQcRPcFmdFqTHE2wCAK0iGzX2TOi4I2/Pee2sZkEWTePFgd+ps3bx5taSHSeR+Ml1dffVU+9rGPyV133SUzZ86c8Hi1ZChykiRJi0zEn9PX19fxciFmzpwpU6dOlS1btnRs37Jli8yePXvM/j/4wQ/kueeekwsvvHB0m+0R9OSTT8o73vGOcc+5G9UvF83YpqZepXaqfj9yvI/HKadEocjR+NFYpXmQtkROTV/LK9JS7G+pllYUwkpzozpY5OytqcfUa10wSxTCGlXppe9KD27NPGprY0Xh47SmJL0TFFIeZdfT3LwT3I4RVR6OtMiShuH3FxkbkFMbPh6NG41vP1OqQBRYQDX82uZAOPSnTZsm8+fPl6GhIbnoootE5GfrNDQ0JIsXLx6z/6mnnir//d//3bHtuuuuk1dffVVuv/12GRgY6Hm+NaTmkiRJ0hIHqvzLkiVL5PLLL5dzzjlHzj33XFm1apXs2LFjNHrssssukzlz5siKFStkxowZcsYZZ3T8Xrus+u1tUv1y0eQy8ifUVk4lrcZL2SJxOG4U6mrts95ebu3mtRJVlKimkJRFkNZGyYi+TS+FYJLGQHZ+mpeus+2PEtUxq7VvR2ta22q6Vlr2RNJ4aSx/fuSnsKgN3N5r0XqrdE1VjuneiZ4Rexwdl/rc0Lxq/Iy1fjU7b72fKBFXP0eh4iJxGDYlQXvtJLqedtxSiHgbHKhQ5Isvvli2bt0qy5Ytk+HhYTnrrLNk3bp1o07+559/ftJ61tSSmkuSJMmbkMWLF6MZTERk/fr14W/vueee9ifkyJdLkiRJSxzI8i8HO9UvF41goHBc2lbbxEjVU2tm0LDn2lDGmiZGFGZJjkBLTdZ0VBXAQuGs1ARKzVvaLng8Zef993RudC5kRonMVyVzhN+/NKeaB6zW7Eb7UwY9nZ9C7Qb0Po3W20LBF1G1CHJS6/hkAovq0lkHNt2veq/VXjsKjlH0/rZzjMaNUgkoQIjaiUfXs9Q6Y7JK3OfLpSE1lyRJkpbIl0tD9cuFQkdJs4gagkWaC0lqtQ6pSOKJtCtLNA8di6RQct7S3CgZis7TJ0VS4IKV4mg/76wkyS4KYbWfI4k+CrQoSYbjrS1WoqaKLp0nrQGFfut6R9o0BW3Ye5Pm4cele9TWBfO/s+Mq9p6j8UiD9+dAGgONEdUIrHWe0/NQo7lEjcTstvHcT72SzcIaUnNJkiRpidRcGrK2WJIkSdI61ZrLyy+/LCLlDPOokVSUvU+OxlpU7bZOc6+KlzKv9bNtihWpzpR7EJVcp3UhE1JkqtDxqac7rR/NkWpokUnQO4BrTQo0fuQ0nwg0Lh1fiZpX6fmS09zeJ1FNNHJgU3Z4TeuEUo6KbrOmL3+t6L6qzf+hPBpfs47GI3NUr89zaf0oC98fs9QgrNc51ZKaS0OaxZIkSVoiXy4N+XJJkiRpiXToN+TLJUmSpCVSc2noubZYKXHuzRCK7I9jj29tuz70tzQGhQXT3JQoBFjDX2tDkS1+G/l57DHVbm9Dbr39vjbxkPwJkxX2SXOLQpGj6r96vjNmzBj9jvq5RL4zhRJrrW+EQpGjhMPIj0V+xpp6c3bcqI6dXY/xhiLX9nuKQpGtL/RgDkXOl0tDRoslSZIkrZNmsSRJkpZIzaWh+uUSmXDaqC0WZfa2UVvMEtUWi7L3yQRBJhD6rTW3KVbV92OoScOulW9KZedrz8WfA5X0t+eipg8KXSUzGjUcqy3N7/fvteR+bZsEGpPm7c1cZM6zde/UvEUN4fxc7THpOlp0Wxu1xZTS/Rq1tojCfaO/AXZ8albnxxc5tGqL+fn8PJOaS5IkSUuk5tJQ/XIZGRkRkQPTLEwrArfZLMx+FyVaWfyxogq0dh6ESrylZmE6T5VSS+de47ylZmFWIo00kdpKwr02CyMiDSDaVtssjOYR1UIjSbqmyrEdS68nXTO6d+gZ0fFtbTGq2OzPqTYhuLZZGGn8uh4UFNJmszBK3v55bhb2ZiA1lyRJkpbIl0tDRoslSZIkrdNzbbFSXbCofhj1kdfP1uFd2+NeIVVYP1OeQW2v9mgeUcn4yClcqs3mVX0KfiDTTGT2K5nuKH/GnzuZdSy9OvQpR6XbmKX9a5uF0TZ/flRy36Jm1sisUzLn6bWy6+1NPfSdHcM3lbPfR+0Ses3/Iejcd+7cOeZ3bTzHZIrzJkEKzKH9S80B2yA1l4Y0iyVJkrREvlwaql8uRx99tIhw5VQrjUfZuJEkEzkTSxIYOdx9mCVJ5QQ50KOQawoLpmNRlr+eQ6TNlMKOScPR8VRCpgq+dIzIMd5rdnNJCo6qCyulas41xyrN2683Sc10T9PcqNo2aRHqmKdrS3P1Va79b/3cfSVuO34UolsbnEKO8agVeO29QETnWfsHWdev1noxEbK2WENqLkmSJC2RmktD9ctFbbulN36N3ZQSxCiJMvJTUFJa5PshqdxuUz8G9XghqS9q9Uv2XoLaJ0eJbf449rMdwydF2vUj7Ud/S306autgRT4a0jqiEGD/O7t/FC5Ncyv50KKW0D7MVoSTIv01G0/Iq++VQuG4pbH8cxCFOtv50jH9M2j3s3itp+TrIE1O14+eQQpn9udZCtcn0ucy+WS0WJIkSdI6aRZLkiRpidRcGqpfLuqE7DUU1P4mChkmhzSFNSuljGf9rGo9zZvGtaYBDTslc4Q3hZSoLf2vIdk0XwphJXyosjUzRGZIyq5WbLZ3bcitEgUK2P2jTHelFOpaEz4eXTO7trXh4xqG63/XjSikvDbYJKqXRcePghOiSgR2jtEzRC0ieg1j13UpBUREJliiZFpuk3y5NKTmkiRJ0hL5cmmofrnUhhHXSCb0HUlNkaRZG9IYSakUVl2qpup/a6XEKPSSpEpfCdeOQUEEdEwKI1ZJmpz3/jgijcPa1q7S2ma6H51nqdaaP09Lr2HBUfgpnR+NT/NWTVG3abKwSLOOtrpwVNGYrkUU5k2BFuR4p+cnCm2n0HzSLPyxKdTePpc1lcYttc+o14QpEKG2VmHNXCeTfLk0pOaSJEnSEvlyachosSRJkqR1qjWXc889d8y2yARG5hyFyonXOux6rVNE/x+p35TVXFtGPmq6FGWk27pqfv9SbgOdi5pxKN+BGmAdddRRItJpptGKDLqfdeiTecmbBEvOfsr8jxp90fqR+SdqkxA1zFKz0fz588d8t2PHjtFt1ADL5zdRJnjpvvUmHru/rn2pYoE3wUb3qIWe2dp2F1FFhIjoeSwF/PT696TXuU2E1Fwa0iyWJEnSEvlyaah+uXz2s58VEZYgShJp5ICLNJeobTDNo7aSayRdk6RGY1AdsSjMOHIwRy2Ea52yFDKqEq/VjKh98lvf+lYREXniiSdGtz3zzDMiIvLUU0+JiMiJJ544+p2Op9qNSKMt2Sq9SlTROKpCTI3VasfdtWuXiHRqEdu3bxeRTi1s69atIiLyrne9S0RE/u7v/m70u1NPPVVERH784x+Pma/V8vw62+tDYcFRMAhVeSAnvx/Lzs3/zo5HVQ9orGi9ay0J9FzSc+PnWKu5RPuXQqJ1Hp/+9Ke7zn887N+/f1x1wsbzcrnzzjvllltukeHhYZk3b57ccccdaGESEbnrrrvkq1/9qnz/+98XkZ9p6DfffHPX/dsgfS5JkiQtoZrLeP71wgMPPCBLliyR5cuXy2OPPSbz5s2TRYsWyUsvvYT7r1+/Xi655BL5l3/5F9mwYYMMDAzIBz/4QXnhhRfaOG2kWnP5wAc+ICIsUUV1f+xn0khIm9GFpn4NCkkhRBRGbMeIksD8XO0xbaimH9/uF2lokR2c7OwlLczPzX6n26w2oz4X1WBERAYGBkSkWftjjjlm9DsN0Y38GVRjis6T+vgQWi/NQhK9flYNyn6n52lDrvv6+kSkOd+FCxeOfnfmmWeKSKfPRbUUq7novCmcPZL2a6s+R9edQuej5GDSZvyY9rPVriILQe0zGIWqU/Vn8j16jaVUzyz6O9U2B8ostnLlSrnqqqvkiiuuEBGR1atXy9q1a2XNmjWydOnSMfvfe++9Hf//la98Rb7+9a/L0NCQXHbZZT3Pt4bUXJIkSQ4SRkZGOv7ZHCtlz549smnTJhkcHBzdNmXKFBkcHJQNGzZUHWfnzp3y2muvyfHHH9/a3D35ckmSJGmJiZrFBgYG5Nhjjx39t2LFijHH2LZtm+zdu1f6+/s7tvf398vw8HDVPD/3uc/JySef3PGCapueo8WsM5RMYKrSUn0gKkkeqbY6lj1mTe2o0n5qAiNzStQsrDRvqtvk97P7UH2ymkx0+x1l+eu46tS2pg01L9l6WPr9iy++OLpNnd8//OEPRaTTNPTqq6+KCJeij5zOFHQQNcyiEvDkLI+CJez1USnQmvhGRkZERGTmzJki0rkGairTtRBpTGpawcCeV3Qt7DyigJWoZQGNG/3W7qPXPQpdJoc+PZdRkAI9I6X9/LHouaRKARRyrWtkx6d7obYVQq9M1Cy2efPmUVOtSOcz1hZf/OIX5f7775f169ejubktMhQ5SZKkJSb6cunr6+t4uRAzZ86UqVOnypYtWzq2b9myRWbPnh3+9tZbb5UvfvGL8s///M/ynve8p+d59kL1y0Vtc+QUKznKotpBUbhi1DaYNAySVqJkM5IOqYVwJG2RYzdyypLj00qCXvuJzslCzb9IwqRtpCGqRKMSo/3uhBNOEBEOFSbtzUua9J3I2LUpBS5QQqOiEl8U6izSnJ+eL2nJ9neqtdlzUY2I1p2eh5rW21ZSp6AXamGt6HpQCHqUFGnHJ41Bx7XSrg88oXOi5ya6r0tpBn6+tAalhoST1Vb4QDj0p02bJvPnz5ehoSG56KKLRORn5zM0NCSLFy/u+ru/+Iu/kJtuukm++c1vyjnnnNPzHHslNZckSZKWOFDRYkuWLJHLL79czjnnHDn33HNl1apVsmPHjtHoscsuu0zmzJkz6rP58z//c1m2bJncd999Mnfu3FHfzNFHH92Rr9Ym+XJJkiR5k3HxxRfL1q1bZdmyZTI8PCxnnXWWrFu3btTJ//zzz3dodH/zN38je/bskd/6rd/qGGf58uVyww03TMocq18uagecSHx4pPZGzj+KtScnHjlNo1pDpCbXOiHJLBEdI3LQW/MF9RH3WCefHt+aKryZy64V5aHoeFTmXedG623ReWgQgeaU2PlYyGTixyVTiN2fxtXAA81zIfOHXVvfnM2uga4LBS5EfeGpDhuZ0QjqXa/jlVoX+O/JHGWDKmqCY0omb/88UrBJqfmXN6lFteLs/rWN+oiD1aHfC4sXL+5qBlu/fn3H/z/33HM9jz9RUnNJkiRpiQP5cjnYqX65qPRGobclLcKHB5YWMpKGKKSXpLKokRBJh9Hcaisg6zFpHpR9XFO1ljLv7fgUjuvXr9R4StfDhhurRK/f2XlRfTLdZrPf/fgWCpf1+9nvSAP07bDt8SksmK6Zfq/na9eA7mVyXPsQWtJm7TUjjdLfp6V7w1exsNv0niDNnLSmqBkfVbGIqlLYY5J2r7+NgnRK1YujkOtIW7LbJiPEVyRfLpbUXJIkSVpi375944pEm6zotTeSzNBPkiRJWqdac/nRj34kImxOsaq5qsKlXBb/XW1JbTL1RM5vinEn81Jt/ow/Jqnktbky/jj2M0kyOr4155HpwZv9rLNft5GpzJ6Lmm7UMU5mNzIr6jyo3zuZaSITCBVltETOdZ0blca3x9Tz0+8or8Oimfn23vHb7BjRvOmeJ9MQ3TtR9QoNpig9q/53tfkolMfj173bMamIaa8tOZQoB8Zen+hZbZs0izWkWSxJkqQl8uXSUP1y2bZtm4h0SiE+jLMW0jpomzrdrMShEgk5zS1eaqqtCkCQlEj1xkhr82OUbKteW7NzVGesDV1WSAvT9aOMasp0J2mSWjHTWtnj2/Podi7+Oxq3FHZM5+LDWUlSt+fiNSg7B39OImMz+mke9hrrPUGhxZF2UrtWURUDO3/9bOdW8/xS0AZpCqQdKOTkJw2HgnV8YIkfrwZ6fiNrx0TIl0tDai5JkiQtkS+XhuqXi4ZoUtIWJZ6RtkE21ag+EIV2KiTdRDb60rwjSZGkHJW8qPlSlExJkjTV6IrawZKWR2GqVBvLz8cekyRuOieq9Ott9LXhu1FVZArVpTGiitq0VtRemK4n+evo3H1TNtKIo4rZdjyqQRfVRIsSG0nTsfv7xmclbYmumddsKfTbHlN9XLYdNq2z/440NPIB6fFtyDUlK0eJrBMhXy4NGS2WJEmStE6axZIkSVrkUNRCxkP1y+XEE08UEe4jXzIvReYtMqMoURn8kmM8Ktuvnynbm+pwRQ2IqAd8ZEqwTkU1IUZZ+9Rn3W6jrHNvLio1U/LzFxnrzC5lTeu4FNJLY9B19L8hEwvNg0JMKeiB5qHo+ZZMQ1SxwIcR030eNaGj49L41pSjNdyoQoSOZdegpioArbddWzIJRm0gqIpAbQsMP++ojQX97SiZ8w6F2mIHO6m5JEmStES+XBqqXy6+ZpOFJFKCtIMo9JLGjKoWRxViyVnd7Rh+G4Xj6rmU6hr5MF/S8qImWnTupRpTKtVSa9uoHTFB505VeqP6cXSeJK36a1CbXBqF41KNKQpFjrB1qCjAQY+lraOts5q0lGiN6D7UdbbjkibirQD2O23VTImyFNYcaQVkvaDkWXouI82FoPX2z1Ktlkzjtk2+XBpSc0mSJGmJfLk0ZLRYkiRJ0jrVmouP5befSw41hWLWI3OA36c0PmXLk9mD+r1ThrE/BmU8k4mAzABkjqI8jZpzoixyyrin0vgK1b+KMtetaYhMVD5XhsanZm52Tf060Hc0Bh2LcndqzoXWgEr/W/R7uhY6fqlFhH8O6DxLpkzvhK91YNOzRHXByKHvj91rufxux/fUNhik49PcMkN/8kmzWJIkSUvky6Wh+uVC1XQJCjH1TjmSWiIHn5WUVLKMQpdFxkpxFIpMFWUj5zBJPiWnsw9LLVWQ9lpGad7dftcNkvpoHnoslfIpHJcywCkcV0NoKVs+OncKs6VsbNugzIdrl+btWzyX1qXXPwI6D6stUUUJfz9TCHqpircP+KitFEDPpc6XWjZbfEh+dO+L8D0ctWeOLCC+SZs9Zm1Iedvky6UhNZckSZKWyJdLQ/XLhTQWsmtTuLGXYKLaWyKNxEj1iqIKp1Gb45JkH4WHks27tkqvQjbsSCOiOmy0P4WdRppf1F/E4iVAugYkXZPPgLQCkqT9vEmTom0kUUdzLJ2LH5/8JVQ7LfJFkdYW1Wuj+4TqmVFfnsgXQf6jyO9Z8lPoc0hzpJ46GiavYdt2XB2rVPHcz9Nqs5T86b/rNm4b5MulIaPFkiRJktZJs1iSJElLpObSUP1yIRNYFEpL5hwyDUXmKh3DOmyj5mJkAqFWyX5e9lhRLa+S2YBMgt7pWAqB9CHRVLuMajtRyCgRZehTbSkK0CAHsM6Drg+FUFO7W+9kJRMOXZ8ovLtUEcE7wUstEWrMKZGpz2LnEdXQin4btRymsHoyTUaVHyJzq/0+qk5h10MDBSiwhRqJ+bnaY5DpmELK6Z6PmgNOhHy5NKTmkiRJ0hL5cmmofrlQqCFJpKSJeCmBJGvSImpqddnfWom+pvovSS9RshZpY6RFUEOmyPlMWliUNGihRk++MrGV4kgipTXyTmq7D52nT9yrrXpLknQ0r1Lot098JG2ppLH645fWyt/PdB8StKalNsGKXmO6/yiBuMYxTkE4UStm2laqtUcal7/n6frQvMlS4Odvj2UDBVJzmXxSc0mSJGmJfLk0ZLRYkiRJ0jo9Z+hbqAkY1Vfyb+VS7obP4i05GmtySCKzhx3PNmSqLemtRNUAorL6dnw9B8oRiHJ8ormRc5PMUVSeXvezY1D2tjfFRO0P7H7RvEv70Nr63A07x927d485F6VmDeyxrDlKA07oevqx7PFJWvV14ewx6Zmy+3nzY20NP73H7DnpNabnhhzoFPBD9fqovpueM1XfiEzS/jh2PwpEIBNi26Tm0pBmsSRJkpbIl0tDzxn6UWvWbkQSFWkHvu5UqQJA5AiMaoZZbUw/a4MtkbGtb8kBSy2HraROEp0/d6qXRY56aotcc36l9rUUMupDnEuh3zRfRbUHqxXq3KwEG2m4+lvrlNXfUiVjCq+lemY+cMKuAa2LdybbcyEtOarbZdfDO/7pGSmF8EcN8ihs2wfO2PtQM+ipEnh0/1HAT1Tzz46hz549z6jiRxTwUWprXluLr1fy5dKQPpckSZIW0RdML//Gw5133ilz586VGTNmyIIFC2Tjxo3h/n//938vp556qsyYMUPOPPNMeeihh8Z13Fp61lxK4aRRtdbRg0K4YpSMSDZYksrIPky1t/xcRUSOOuqojvnY76MQTErMorkRNXXYSuGTUW0x8n9QHTYKD418VnSeEZG9n+plKRTGWwpT9dDaRjXR7BpQNWKdo93Pa7FUd4wkdTuual+UzBmFzdaE/ts5Rc8qVSomrZqOFdULJP8e1Z4jzZ/2r6myHlkxus23DQ6U5vLAAw/IkiVLZPXq1bJgwQJZtWqVLFq0SJ588kmZNWvWmP0fffRRueSSS2TFihXya7/2a3LffffJRRddJI899picccYZPc+3htRckiRJ3mSsXLlSrrrqKrniiivk9NNPl9WrV8uRRx4pa9aswf1vv/12Of/88+WP/uiP5LTTTpMbb7xR3vve98qXvvSlSZtjvlySJElaYt++feP+V8uePXtk06ZNMjg4OLptypQpMjg4KBs2bMDfbNiwoWN/EZFFixZ13b8Nqs1i5ESLHJhEVJuIysKTikthk1G4cW14aJTBHJVZp3OJMphJNacaWlEYNAUz2G3e8UolzO0cKTQ2Kl0ehZ1GLa9Lzat8EzQyDdW2x9X5lJpC+eAEWhdyDlPjLgp5pcACWltvQqoNRaawYHrOyLzkvys1w6M/gNE1oCAJuq/9uZDJjOYWHTt6BkUO3gz9kZGRju3Tp0/vCHoREdm2bZvs3btX+vv7O7b39/fLE088geMPDw/j/sPDwz3PtZbUXJIkSVpiPM58+0IaGBiQY489dvTfihUr3uAzGj8957mQBFaq0eW3kcObJA1fK8lSSiRUogCDUnKfnxtJzSWHqpciS1qeSnlRMiq1fT7iiCPGHIPCT/1xLKRZUE00/S0lHNK4kcpfK30SkSRN92RU8TjSjEqVkhUdw2o12hyrlIzoHejk1KZQ2lLjPX8s0kQi7YrOMxqDQq7t+un9ap9pXa8oKbb27070t6hUI64NJqq5bN68Wfr6+ka3e61FRGTmzJkydepU2bJlS8f2LVu2yOzZs3H82bNn97R/G6TmkiRJ0hIT1Vz6+vo6/tHLZdq0aTJ//nwZGhoa3bZv3z4ZGhqShQsX4rwWLlzYsb+IyMMPP9x1/zbIDP0kSZI3GUuWLJHLL79czjnnHDn33HNl1apVsmPHDrniiitEROSyyy6TOXPmjJrVrr32WvnABz4gt912m1xwwQVy//33y/e+9z358pe/PGlzrH65UHauQuo6mZeohzmZl6Jcj2ge5CynnIwoYzdSaUsZz3aeijdRlJza3rxQckyS096vKZkxqBQ9zZPmSOfus9Sj7Gl7TMq3Uex3On5tuwbFzjEKXIjua6q4EJmG7Ly1nhldR0uU6e7bGXTDN8WiHBUyL/nzsPOgoJdobnSepTwhhYJTomCDWvMYNQs7WM1itVx88cWydetWWbZsmQwPD8tZZ50l69atG3XaP//88x1r8L73vU/uu+8+ue666+QLX/iC/MIv/II8+OCDk5bjIpKaS5IkSWscqJeLiMjixYtl8eLF+N369evHbPvIRz4iH/nIR3o+znjp+eVC0j5lz5KUT845ql7rj1UKGyQJhqqjKlEjMTq+1oCiJmCkrVD2NoV2UsBCqQ2yH4Nqvul4FABA7XSjCrG1mcw+hLbkPKV5+2vVa6i4PX7twxqdHzXw0jlRpWnK6I8qWxAk2UdVEijYQI9v7fU6np23fqaadbS2NcExJSuGn6P9vjbL3z+PpUAE2m+yOJAvl4Od1FySJElaIl8uDT2/XKIKviJs0/eJU7W9VaJksFK/EB9GWqulUI0z3za427lQK2g/X6rwTD6A2npSlCzoNb5Snaqo7loUxhwlstJ9QuG10TxK4ae0fl5zKrUcjjSLKMycrhmFLpMmR8f0PheqDEz1yUjz8/etCPsq/fMbJTbabaTN0Pi12mYUrh+lC9A1ptbr1FeGLA5t0Gu2vf3doUaGIidJkiStk2axJEmSlkizWEP1y0WbB0XtjkW4NXFNuGdtzaPaqgBkouo2vj0HMjn5tssWarBFobFk7qip/UVmCXsN1BlLpiFqGxuFOFv82tt96Fr5eVNAh3Uw14R+k+nJmjY0zJcgc0pN2Xk6Pu1P9cP0OtL6UNAGXVuqXRaZ22i+VNNL96d7mJ5LqswQPUv6nR0/ap5H845q1tHfgChUnKpY2GCGNItNPqm5JEmStERqLg3VLxcKs+xVs/BjicQNiGisqP6QlWS9JFBq0xtJ7yQlklM7kugjpzZpEVHdMbtmRx555Jhx/blESZKWKAyatB9yAOs1oMq2FqrS6+dN91VpbnrdtWUuObCj0Gtag9r6XVHNOnLQWynfa4ilsHq6n6KACKW2nTjtT0nQvk20XRe9BvTsRRoRnSeFg0dBDeTQJ226bfbv3z+usX+uXy5JkiRJTGouDRktliRJkrROteYSZcOWysh7522p0ZffRtnkVjVXUww5MH05b38shcwR3kRCJjDrJCSnpjeVlJpAdZuXCNfXouvi8z8om5xyLKJ6WaXWAt6RSmY6u35kHonWIdrfHsv3Vy8FfkR12Ghd9Fh2TdU0GZnF7DWj6+jNYhREUKqr5udmv6PKDH7cUhsGund8XhM9D2oes8eKTJ50H9bUv7NzpHuTzr1tUnNpSLNYkiRJS2S0WEP1y4Wc4JFDLZIOLVHdJJJ4IwnJhrr6SqjkvCdHupU6vSOVnLLkuCapPQoKiEK0SaKh86QMfa3BZPHh1Xb/qFEahcbSmur1sdpbTStmu5+fqz0W1dyiysdRHayogjTdt6VqurrOpF3pNnvNdP8oS90ehwIASNvwlgR6bihknSwQUdg2PSPRtdCGafYcKDCIgnXoOfBOe6qIQXOM6ti1RWouDam5JEmStERqLg0T8rkoZAePqqPW2khVKqPETSvpqTTUayXUUt8IPw8Kf7bSdZTQGCXfUbIgaST6XUkr8POmJE2CfEqkcdX4ikrjR+MSur8m89YSSbx2XAqlJSjJke4PRcez14y0CJ+0SGNRfTK7jvpb9XFYzZUqkvswX/KD0DEt0fOl50z3ld2maxnVHoxCl+lvB62tJX0uk09GiyVJkiStk2axJEmSlkjNpaFns1ip1hA54b1jNKr7ZD+TM5wyk8lJ7YMCarPgI/OShcwSUXY11QXTz1FzKTtv3S8yOdL5RSGv9rdRiHipPbN3Ckcl6e0xI/NEFBpdOpZCc6RzIbMlHZPqwfnQ71LYtv6WHPrk1CazcHQu+tuS6bEmo5/qn9H3UYi7hZ5pXzuNnmMy/1HQQVQdoVTKvw3S59KQmkuSJElLpObSUP1yIelJIWmIKpuSk1Ah6ZDGJ2mcpKZIEogqMJN0HTkESVIjCVClJwrLpJBbcj77efnPfgzFStk6RyvRlypd+/9XRzGFYVOILGl5UU05gmpM6bwpUIA0Ef1sw4L9HC00vq4laS7RPURQmHQkqdt7hwIL/PHtvRm1f45SBEr1xrxloNQsLGqQF2nVlkg7pX3o785kaQqpuTSk5pIkSdISqbk0ZLRYkiRJ0jo9ay6kstaakpSSk9BDJrZSnSUPlQmnzOvIOUymu6iygD2WZivbY1LdM38sMvWReh+Vy6cAg1IdKW+apNLukTmFTEl2HlGuh0J11Wi9oyZQtFa1Zk7FrhU1w4vGoHWOKhwodJ9QAywyQ0WmVLqHyRxFkDnZ33e1NdEo0CIiCvyg2nI0p1INxDZIzaUhzWJJkiQtkf1cGnpuFlYKYSVp3GdB1y4+NSgjxyRdGB8WStItSY4kvUfVbkkLo0x+P397fNJOIg2DtCvKvFZKbWAjqZ2g2lV+LMrsJk2EHNJ+LLu/zXSnsF1/DnaOvUrIUXgy3QukRUYN4ei5ISmbnPw0ht+fzoXuhei60zNin6UoZJ2eqaimXGTtoHte7wV6RmyQSmbovzGk5pIkSdIS+/btG9eL6+c6Wiyys1OVXgtJFjXHihIVLVF4pYbNlrQfhW4MskmThF4TPkySN0E1mKJ+GhQqTJJjrxIS2c1rEjej6tL+sz8Xgo4ZVcCNKj1TUmkt0ZpSqHP03JBWRf4JSj622oOHkgb9XP05+P+n+4TmpveannNtImaU4FlKcPRJl6UeUxRSXqPFjofUXBoyWixJkiRpnTSLJUmStESaxRp6frmUanpFZoAax7vIWBMBOUNLmca+3Hwp7LjGHFYyyUXNjsiEQ2NQDTIlctSS89vvY8coXR8fxFAK0Y3qx0XZ3lHod23ZeWq14E0nFgrfrb2v6Boouu6l6xM186LgBz9XO0Z0zShoozYUOQpPjkzfUbUJO67FV3Uorbc3TdL+FChC92vbpFmsITWXJEmSlsiXS8OE2hyPt2ZPqUKx1wBKjncKqYwk7mi8Umixp5SY5Z2UJNVGoctRBd9ux/QSYK22FEm1tC5WQ/IBC/Y7PSdKnIuaUZFESmNYadk3qLLziMLNSVKPtBnS1iMtqXQdfVBAaX/Skv09Vmqi5TWd2lbg0dzonErj+iAQqqhuicxO9Kz2WntwIqRZrCEd+kmSJC2hmst4/k0WL7/8slx66aXS19cnxx13nFx55ZWyffv2cP9Pf/rT8u53v1uOOOIIedvb3iZ/8Ad/ID/5yU96Om6+XJIkSQ5hLr30Uvmf//kfefjhh+Wf/umf5Fvf+pZ88pOf7Lr/iy++KC+++KLceuut8v3vf1/uueceWbdunVx55ZU9HbfaLBb1tS69dSPzAqnTXoUvmYEix2HkmCzFuntzDpWuJ/MI1RaLaqFFjdVq1WU6pjIe53pUq4zMYj4fhioG9GriKJk0yYmsn/X8yCxG1Jakp+uu1D4PhHfkkzPeoutrz8+fAwWxWHzgB13zkmmtNrvf70+VPsicV5PXFpmj7f5kzmubg80s9vjjj8u6devku9/9rpxzzjkiInLHHXfIhz/8Ybn11lvl5JNPHvObM844Q77+9a+P/v873vEOuemmm+R3f/d35fXXXw//NlhSc0mSJGmJg80stmHDBjnuuONGXywiIoODgzJlyhT5zne+Uz3OT37yE+nr66t+sYiMI1psPJJ0JGnUOO1LUkZtaHH0HdWR8pJxFMIqUlctmCTBaI4kfdL8o+xjcrKSplPb+pVqxHnnup3D7t27x4xPaxpljJPGpZJ3X1/fmHH12tlM9poHmCrsUoWISAqmel8W0iwiqPJDFKZN1gBylvu1Kt3fkRWg5rtu2/zfgJKVocYSYqH9D9YM/ZGRkY7t06dPx+Z2tQwPD8usWbM6th1++OFy/PHHy/DwcNUY27ZtkxtvvDE0pRGpuSRJkrSEdqIczz8RkYGBATn22GNH/61YsQKPs3TpUjnssMPCf0888cSEz2dkZEQuuOACOf300+WGG27o6beZ55IkSXKQsHnz5g5NvJvW8tnPflZ+7/d+Lxzr7W9/u8yePVteeumlju2vv/66vPzyyzJ79uzw96+++qqcf/75cswxx8g//uM/jumVVKLnwpUlalRVcmCTSk6qeeTgi7J4S2ow5Sr4+PiSmYacsTX5NuP9zn5P5rmoOGBN4yc/rh+DzDqUOa5rWtuygL4jon7sRE3VCFoDKrlvoRwcJTIX0TyifBEal8x4dA1qyvDb+Uf3fGRmjZ5Bmr/I2JYc9PfBnktNaf7SvTOZeSUT8Z/09fV1vFy6ceKJJ8qJJ55Y3G/hwoXyyiuvyKZNm2T+/PkiIvLII4/Ivn37ZMGCBV1/NzIyIosWLZLp06fLN77xDZkxY0b9Sfx/0iyWJEnSEgebQ/+0006T888/X6666irZuHGj/Pu//7ssXrxYfud3fmc0UuyFF16QU089VTZu3CgiP3uxfPCDH5QdO3bI3XffLSMjIzI8PCzDw8NFf6JlQmaxXhckkuwonJRCOyNJncYlyKFK5bu9pEZORdJcbIinvvEjBymFukatmO08VFXdtWtXeH4KBVBoVntUnYCudW2wQakMew29jktrQPdTFPbus/3t/hY9BpV2P+KII6rG8NomaR0kvVPghA/HtucVNcij72j9Ik2ErAz0TNH9RI3ESDulEHRPr2HvbTFRh/5kcO+998rixYvlV37lV2TKlCnym7/5m/JXf/VXo9+/9tpr8uSTT8rOnTtFROSxxx4bjSR75zvf2THWs88+K3Pnzq06bvpckiRJWuJgfLkcf/zxct9993X9fu7cuR3HP++881qZT/XLZSIH89KHlTiosm1UX4v8GiQheWpttnRcsleTFEyamQ+5JU2HtJlIc6EqwPbYUdMyGjeqsEuaC0mY/v6IfGP2cxQ+Xtq/5lh2TN9KWCQOg43WL6qxR3470pJpXJpjFJpL0jsl/ZKfzF/30jNO4cweCjGu9UFS6HeUIlCjwVgm8w+4Ml5fTtYWS5IkSZIK0iyWJEnSEgejWeyNomezWG1oLDXmocZd1DjJjxWZa2rmZI9dmmNUTpxMWpQxro5gkbhmVbSNzH8a827PJSotX9t4iswikSOfTDc+jJTWsbaBXPRdqZGUn4c1X0WtB2oDF2hufp2p3QBdM+0/b39ba+IjE7Ped1qVoFRTzpt9S88ZVVWIGtPRMaOQaLrno1YVemz7HNPzVvP3oS3y5dKQmkuSJElL5MulofrlUls7iDQQr7FYabJGYyk5zRWSyvTYVsKjdqpKFIpMTb2ophclHEXNmuz6qdRJiZAkxem5ROtYW9uJmr5RsETU1IlCuinowP/O7qfYumC0f00FZkrOjULKoxbVdryoJppdF3JS6zzsNj1X0n4jzdKumd53NfeEPa+o7h2tBznQdd4liwIFG/j9SSOOgnDsfRUFwtCx2iZfLg2puSRJkrREvlwaMlosSZIkaZ2ea4tFjji7zaqlqgJHda0iNbXUUzvKXKcs66ipE5kNyLlJZgM1S0S1l0pmAzUvUH4EOejJxOed6rVZ1nR+CpmGLL4MP5mByNFNTmeFzCNkRrPHqrnHas1iNC9aU5+LYcegumrRtSWzWOQYp2dP70MbMBCNQWMpJROVnaedvz+Wh0r/+7na41PzNDKtRQ3NyKzYNpnn0pBmsSRJkpZIs1hDzw59ylqmrGbKDq51rHlpxYc7eshprxIVZU1HGc8k1ZLUTBnSNZn8kdRv53vkkUeOmQ9Ba+SdyBRAQRJpyZmtUAtmdUjrvCncvFTpNyIKFSfnrc6HAgsi7Zs0wKj6s8jYWnh0LQg7hoaZb9++fczvaN40Xx+wYMu1Ry24awId7Daqnuw1GLsfPb9R6+iSdcTPzf6t0WOR5kJBLG2TL5eG1FySJElaIl8uDdUvl6i3CtVvIolHKSXYeW2DpDgrrajkGElPkT/BzqlUP8wfk7SZyLZLY1nJzmsbpZBK0sy89ljSNmmdvaZKErK16aum4BM47RiUfEcht/53dn+6d0gz8yHdfk7d5kY1rMi2T/496usSJY5GLZXtPRFVf671oXjtyn5PYc2RT45C4fU60jNoocrRvs0y3TtRTT77nV4LSiqmsPS2yZdLQ0aLJUmSJK2TZrEkSZKWSM2loZWXS1Qun6BMYPs7VWmpRlLkvLfjqsmGspujzOso7LlUE41MFX6sKJjAjkHmADqmframGF0PXSM7bwq5jcJ8o0ZVFu+0jc7Jzje6d2jeVI6dzGJkRiOHsTc/lgIdyITk19lWFqBKDn4sOw8yt0Z/eMhJTWY02uad/HQv0/g0N11vG0QQmcjsGFSJIZq/X8tSuHQv3RMnyv79+8cVLJAvlyRJkqQr431J5MvFESVfERTCSpK311jIeW+3RY5Ucvr533XDS28kSVP4Zu16kKPbawU0VrRWItzMy49HodmRQ5fCQ8kpG+1vUUmX1k8hhz5pLjS+rkvt/nSNo7Wi3+q623XReZQq9ypUty2qKm3x91+pbXENpXpjUaM0elb9fCy6VqW6XzVNy0qBNpNFvlwaUnNJkiRpiXy5NGS0WJIkSdI6rWsuUUw5qc6Uh+JrRlnTRpR1HvUmJzNDad7R/uTYjWozkUOaiBp8KRR0QA5mMgfQelD+Qg1R3o0lyvKvdVyXsuS7Hb9XibDUA76mWgPdE1RLi86dTD5ksqPWBj6bnXK2qGJBr03caG56TpS9T9UMoutO9zfl7EStHyarpH6J1Fwa0iyWJEnSEvlyaei5zXFpW+T8VqmGpJuOSTmNxToEozBYyoyn8MnaUE1/zKhRmT0GzYMq+NJ4NRJXqXKvb85GknepdpWnFPzgw6RLtZ2oWrU/vj1P3b9Ud4ok9F7OhTTu2vXTdd+9e/eYcWvDU6NrYO+dKHw4Gos0gOg5LoXC++eMgjBseHKpZXS3c6JnhTTFqKpGdMy22Ldv37i0pp/rl0uSJEkSk5pLw4Q0l6jlK/W0IDsx+Qz0MyVH6vhUO4g0BupfQjWpSApWSZT8QnSeFPbsk/RK5+5txuQDKLUtVomOpL6onwuNR+tCxyRJmvZTfNVq2s8eU/e3171UYbrbPvR9aU39+PY4PrG3NBZdR38P0z1B/hK65yl0me5XPQeqxxVVQY98OXTu9v6jend6baM/sJSAq1C4uaXXFIGJkC+XhowWS5IkSVonzWJJkiQtkZpLQ88vl5KZgcwjkZOVzADeeUsqcalWk6r4uq3kVFTIaUpqddQ6lUxwpRDkGiITBZWAVyjDnDKvI3NOKfzZO3Yj06BIsx62rpRfI6qXFmWH+8/+mGS2jFrs+uPYMcjUQyafaLyJ1LzSeUett0tmX39f27Ug8ykFpWgNv6iWWxS6bI+v19i2cojM2rR+0bM9kWZ1teTLpSE1lyRJkpbIl0tDzy+XUs2emtDfqBqs3UYSikqMdEzrkPQhwFaaJymOHLv+WCWpmcbw51KqU+W3kXZA+1MoskIJpDRHIqpTFVW3pnOjenBRcit9V6pT5SsZl+Zd40SulZqppluUaFobQhtpXNG5lO5Xfy7kGKdnxF4Db0koBQX4xmD+ez8Pshr461IKRY+qlbdNhiI3pOaSJEnSEqm5NGS0WJIkSdI6E3Lo95p1TG/nqIS+7k85LRSvT2YUdeSXcg+iOkXk8I6y5K1KHvU11/OK8lzIdEJml6jvODlUyaRF50L/H+VAULWEyIlLJg2/j6VkuqlpWlab00LnSQEO3gQXXU+7LcpWt+ekpiR6RsgsRqbp6B6mc4rK2tsx9PnSqgQUyGPR+drr5PPZ7HnS3wA/j5J5OboGbZOaS0OaxZIkSVrkUHxRjIfD9udKJEmSTIjdu3fLKaecIsPDw+MeY/bs2fLss8/KjBkzWpzZG0e+XJIkSVpg9+7dHXlbvTJt2rRD5sUiki+XJEmSZBLIaLEkSZKkdfLlkiRJkrROvlySJEmS1smXS5IkSdI6+XJJkiRJWidfLkmSJEnr5MslSZIkaZ3/B6izhntwnGSSAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 500x400 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(5,4))\n",
    "plt.pcolormesh(object_fbp[0][:,:,64].T, cmap='Greys_r')\n",
    "plt.axis('off')\n",
    "plt.colorbar()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The cube is no longer blurred, but artifacts are present. Such artifacts are not present when using algorithms like OSEM for reconstruction."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "torch",
   "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.9.16"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "324f5aa1485c6e81c0125c86ad1340c603dab08bc11babd97bd730cfafa2ec4a"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
