SFMExample_bal.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "30f616a3",
   "metadata": {},
   "source": [
    "<a href=\"https://colab.research.google.com/github/gtbook/gtsam-examples/blob/main/SFMExample_bal.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a6bc7c57",
   "metadata": {},
   "source": [
    "# SFM with BAL Files\n",
    "\n",
    "This is an SFM Example that reads from a BAL (Bundle Adjustment in the Large) file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "2ce38c72",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Note: you may need to restart the kernel to use updated packages.\n"
     ]
    }
   ],
   "source": [
    "%pip -q install gtbook # also installs latest gtsam pre-release"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "130b0885",
   "metadata": {},
   "outputs": [],
   "source": [
    "import argparse\n",
    "import logging\n",
    "import sys\n",
    "\n",
    "import gtsam\n",
    "from gtsam import (GeneralSFMFactorCal3Bundler,\n",
    "                   PriorFactorPinholeCameraCal3Bundler, PriorFactorPoint3)\n",
    "from gtsam.symbol_shorthand import C, P  # type: ignore\n",
    "from gtsam.utils import plot  # type: ignore\n",
    "from matplotlib import pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "ea6a378e",
   "metadata": {},
   "outputs": [],
   "source": [
    "logging.basicConfig(stream=sys.stdout, level=logging.INFO)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ec459c14",
   "metadata": {},
   "source": [
    "Run LM optimization with BAL input data and report resulting error:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "0586514e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "INFO:root:read 7 tracks on 3 cameras\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# Load the SfM data from file\n",
    "input_file = gtsam.findExampleDataFile(\"dubrovnik-3-7-pre\")\n",
    "scene_data = gtsam.readBal(input_file)\n",
    "logging.info(\"read %d tracks on %d cameras\\n\", scene_data.numberTracks(),\n",
    "                scene_data.numberCameras())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "ad00edf4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create a factor graph\n",
    "graph = gtsam.NonlinearFactorGraph()\n",
    "\n",
    "# We share *one* noiseModel between all projection factors\n",
    "noise = gtsam.noiseModel.Isotropic.Sigma(2, 1.0)  # one pixel in u and v\n",
    "\n",
    "# Add measurements to the factor graph\n",
    "for j in range(scene_data.numberTracks()):\n",
    "    track = scene_data.track(j)  # SfmTrack\n",
    "    # retrieve the SfmMeasurement objects\n",
    "    for m_idx in range(track.numberMeasurements()):\n",
    "        # i represents the camera index, and uv is the 2d measurement\n",
    "        i, uv = track.measurement(m_idx)\n",
    "        # note use of shorthand symbols C and P\n",
    "        graph.add(GeneralSFMFactorCal3Bundler(uv, noise, C(i), P(j)))\n",
    "\n",
    "# Add a prior on pose x1. This indirectly specifies where the origin is.\n",
    "graph.push_back(\n",
    "    PriorFactorPinholeCameraCal3Bundler(\n",
    "        C(0), scene_data.camera(0),\n",
    "        gtsam.noiseModel.Isotropic.Sigma(9, 0.1)))\n",
    "# Also add a prior on the position of the first landmark to fix the scale\n",
    "graph.push_back(\n",
    "    PriorFactorPoint3(P(0),\n",
    "                      scene_data.track(0).point3(),\n",
    "                      gtsam.noiseModel.Isotropic.Sigma(3, 0.1)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "b818e416",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create initial estimate\n",
    "initial = gtsam.Values()\n",
    "\n",
    "i = 0\n",
    "# add each PinholeCameraCal3Bundler\n",
    "for cam_idx in range(scene_data.numberCameras()):\n",
    "    camera = scene_data.camera(cam_idx)\n",
    "    initial.insert(C(i), camera)\n",
    "    i += 1\n",
    "\n",
    "# add each SfmTrack\n",
    "for j in range(scene_data.numberTracks()):\n",
    "    track = scene_data.track(j)\n",
    "    initial.insert(P(j), track.point3())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "adbd6d57",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "INFO:root:initial error: 2764.219293\n",
      "INFO:root:final error: 0.046138\n"
     ]
    }
   ],
   "source": [
    "# Optimize the graph and print results\n",
    "try:\n",
    "    params = gtsam.LevenbergMarquardtParams()\n",
    "    params.setVerbosityLM(\"ERROR\")\n",
    "    lm = gtsam.LevenbergMarquardtOptimizer(graph, initial, params)\n",
    "    result = lm.optimize()\n",
    "    # Error drops from ~2764.22 to ~0.046\n",
    "    logging.info(\"initial error: %f\", graph.error(initial))\n",
    "    logging.info(\"final error: %f\", graph.error(result))\n",
    "\n",
    "except RuntimeError:\n",
    "    logging.exception(\"LM Optimization failed\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee2fe836",
   "metadata": {},
   "source": [
    "Plot the SFM results."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "32cfb167",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQkAAAEPCAYAAABRMTF5AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABe+UlEQVR4nO2dd3hc1Zn/P0ejbvXmIsmSq9xt3DAJmNAhEGzT2STAkk1I3RA2C/sLKfAkLCTZkLIkeRISFkISCNgmBgxeygZICNgmYKvYli3ZKlYdadRH08/vD+lcZsYzoyl3pLF8P8+jx1aZM2dGut/7nve87/cIKSUGBgYGwUia6gkYGBgkNoZIGBgYhMQQCQMDg5AYImFgYBASQyQMDAxCYoiEgYFBSAyRMNAQQtwmhPjbVM/DILEwRCIKhBDnCiH+LoQYEEJYhBBvCyE2jH/vNiGEWwgx7PXxyPj3HhdCSCHE1X7j/WT867dNwcsJyvicFk71PAymluSpnsDphhAiB3gR+ALwDJAKnAfYvX7sHSnluUGGOArcCjw/Pl4ycD3QGMEckqWUrshnb2AQOUYkETmLAaSUT0kp3VLKUSnlK1LK6jAf/wLwUSFE/vjnlwPVQGewBwgh7hNCbBdC/F4IMQjcJoTIFUL8VgjRIYRoE0J8TwhhGv/5hUKIN8cjnR4hxJ/Gv145Hh0ke439hhDiXwI851vj/z04Hg3dKIQoEkK8KIToH4+g/iqEMP6GpjnGLzhyjgJuIcQTQogrvC72cLExFkXcNP75LcDvwnjcFmA7kAf8AXgCcAELgbOASwF1sX8XeAXIB8qA/45wjkgpN4//d7WUMktK+Sfg34CTQDEwE/gGYNT1T3MMkYgQKeUgcC5jF8ejgFkI8bwQYqbXj20av9uqj01+w/wOuEUIkQucD/w5jKd+R0r5ZymlB8gBrgDulFKOSCm7gR/zofA4gQpgjpTSJqXUKxnpBGYDFVJKp5Tyr9Jo/pn2GCIRBVLKw1LK26SUZcAKYA7wE68feVdKmef18a7f4//G2N34m8CLUsrRMJ621ev/FUAK0KGECPgVUDL+/bsBAewTQtQJIW6P/FUG5IdAA/CKEOK4EOI/dBrXIIExEpcxIqU8IoR4HLgjwof+Hvg2cEG4T+X1/1bGEqVFgRKYUspO4LMwthMDvDaeYxgY/5FMYHD8/7PCnbCUcoixJce/CSGWA38RQuyXUr4e7hgGpx9GJBEhQoglQoh/E0KUjX9eDtwMvBv6kafwM+AS4K2JftAfKWUHYzmHHwkhcoQQSUKIBUKI88fndL2aH9DHmMC4pZRmoA34lBDCNB5hLAjxVF3AfPWJEOKq8aSoYExk3OMfBtMYQyQiZwg4G9grhBhhTBxqGbvDho2U0iKlfD2GNf0tjG2/HmJMCLYzli8A2DA+v2HGkqRflVKeGP/eZ4F/B3qB5cDfQzzHfcAT40uaG4BFwGvAMPAO8Asp5RtRzt/gNEEYeScDA4NQGJGEgYFBSAyRMDAwCIkhEgYGBiExRMLAwCAkhkgYGBiExBAJAwODkBgiYWBgEBJDJAwMDEJiiISBgUFIDJEwMDAIiSESBgYGITFEwsDAICSGSBgYGITEEAkDA4OQGCJhYGAQEkMkDAwMQmKIhIGBQUgMkTAwMAjJRG7ZhredgUH8EVM9gVAYkYSBgUFIDJEwMDAIiSESBgYGITFEwsDAICSGSBgYGITEEAkDA4OQGCJhYGAQEkMkDAwMQmKIhIGBQUgMkTAwMAiJIRIGBgYhMUTCwMAgJIZIGBgYhMQQCQMDg5AYImFgYBCSifwkDOKA2+1mdHSUpKQkUlJSMJlMJCUZem2QmAgpQ/rKGKYzOiKlxOVy4XK5cDqdeDwe7XtJSUkkJycbonFmktCmM4ZITBIej0cTBiEETqdT+56UUvtQqCgjOTnZEI3pjyESZzJSStxutyYKQoz9PTgcDu3/gR5jiMYZhSESZypSSiwWCyaTifT0dE0UpJQhRSLQOOqjvb2dnJwccnNztaVJcnJy2GMZJCQJ/cszEpdxwuPx4HA4aGtrIz8/n4yMjKjHEkJoImCz2ZgxYwYejwebzab9jMlk8ok0DNEw0AtDJHTGe3khhAi4LJBSMjAwwIwZM0hJSYlofHXxe48rpTREwyBuGCKhI2oZoZKT6sN7Sed0OqmtrdWEREpJbm4u+fn55OXlkZwc+a/EO9JQ8zBEw0AvDJHQCbW8kFL6XLTeItHf309dXR0LFiygoKBAizoGBgbo6+ujqakJIQR5eXnk5+eTm5uLyWSKeC7BRGN0dBQhhBbFZGVlGaJhMCGGSMSId+1DoOWFEAKPx8OJEyfo6urirLPOIjMzE4fDAUBycjKFhYUUFhYCY5HGwMAAvb29NDY2YjKZfERDPWckKNFQc+vt7UVK6SMOycnJ2ochGgbeGCIRA/61D4EuLLfbzYkTJygsLGTjxo0Tbl2mpKRQVFREUVERMLZV2t/fT3d3Nw0NDTidThwOB8nJyWRnZ0e9FaqKt+DUPAqgFXYlJyeTlJRkiMYZjCESURCo9iHQRdTb28vJkyeZNWsWS5Ysieq5UlNTKSkpoaSkBID6+npMJhPt7e0MDQ2RlpZGfn4++fn5ZGVlRXUx+0ca6vW5XC5t+WSIxpmLIRIR4r+8CHSxSClpaGigr6+P8vJy0tPTdXv+5ORkcnNztUhjdHSUvr4+WltbGRoaIjMzU0uCzpgxIybR8H49SjTU91UkkpqaaojGNMcQiQgIZ3lhs9morq6moKCADRs20NraGnEOIRT+z5mRkUFGRgZz5sxBSqmJRlNTEyMjI8yYMUOLNDIyMrREaiQXdSDR6O7uZmhoiMrKSi3SUB+GaEwvDJEIg3BqHwC6u7s5duwYS5cupaCgAPDd3Yg3QggyMzPJzMyktLQUKSUjIyP09fXR0NCAzWYjKysLu92uJU6jfR71PphMJqSUOJ1On+WXf7OaIRqnL4ZITECg2gd/PB4PR48eZWRkhA0bNpCamqp9T+1uTAVCCLKyssjKyqK8vBwpJcPDwxw9epSWlhaamprIzs7WIo20tLSwx/aORoQQPlu1gUTDv+/EEI3TB0MkQhCs9sEbq9VKdXU1M2fOpKqq6pSf0TuSiGU8IQTZ2dlkZ2czc+ZMsrOzGRoaoq+vj0OHDuFyucjJydFEI9JqUO/n8RcNh8OB3W4HfJvVVN+JIRqJiyESAZBSYrPZ6OjoYM6cOUGXFx0dHZw4cYLly5drNQz+TOZyI1KSkpLIzc0lNzeXyspKPB6PVth18uRJ3G63VqPhXw0aSV4jlGi0tLRQUVFBamqqtjwxRCOxMETCD+8/4I6ODkpLS0/5GbfbzeHDh3G5XGzYsCHkHTeYSEyleAR73qSkJC2KAHyqQZubm5FSaqLhdrujvpC9RaOnp4fKykocDoeWJzEMeBILQyS8UI5RgJaQ82d4eJiamhrKysooKyub8EKJhxhMlriYTCYKCgq0JKzL5aK/vx+LxYLZbNa+lp+fT05OTtQl5P7NaoAhGgmEIRIErn0wmUw+CUcpJW1tbbS0tLBy5Uqys7PDGjseOYmpIjk5WasGnTFjBk6nk4yMDK0aNDk5WYtEoq0GVa9PCY7y0fAXDcOAZ/I440UiWO2D966Ey+Wirq6OpKQkNm7cGFGnZqLmJPQQG5PJ5FMNarfb6evro6Ojg/r6etLS0rTlSXZ2tm6FXVJK7Ha7lghV38/MzDREIw6csSLhX1odqDFL+T7U1dVRWVnJnDlzIn6eRFxu6DGfQInLtLQ0Zs2axaxZs4CxwjKVBB0eHiY9PV2LNPSsBnW5XBw8eJCzzjoL+LAt3nDt0oczUiTUPr5KvgX7I7Lb7Rw+fJjVq1czY8aMqJ5rOi03vJFSTnjHTk9PZ/bs2cyePTtoNajD4cBqtWrVoJHivTxReSTlpaGEzPDSiI0zTiTCqX1wOBzU1tYipQyrczMUibrciJVoSrv9q0GtVisHDhygsbERq9VKVlaWTwl5JHNRvyPDgEd/zhiRmMj3QaEKixYtWoTNZot5fZuUlJQwxVSJhBCCGTNmkJqaysqVK7Vq0L6+Po4ePYrdbg+7GlTlk4I9jyEasXFGiEQ4pdVSSo4fP05PTw9r164lIyODhoYG3Z4/kYg0CojXGN54V4POnTsXj8fjUw3qdDp9bP68S989Hk/YYj6Ra1dbWxtz5841RMOLaS8S3l6SoTo3a2pqyM3NZcOGDbpmx/WOJBIFvUXCn2DVoP39/Vo1qBKNlJSUqH9n3s1qHo+H7u5uysrKNNEAw7Vr2opEuMsLs9nM0aNHqaqq0jwa4jGXRB4v2jnoEY2Ei3c16Lx583yqQXt7e7Hb7TQ0NGg2f9EYCquIRH2oOQZz7VLJ0ukuGtNSJMLxffB4PBw7doyhoSHWr18fUQdkJASKJJxOJ4cOHUIIQUFBQUTNVInyBxmvbdRw8a4GLSwspLOzk7y8PCwWCydOnEAIoS1NwjUUdrvdp/xcKNcuxXR37ZpWIqGKbFwuV8hGodHRUaqrqykuLmbdunVx/aX6JxoHBwepra2loqKC5ORkLXz2eDzk5eVRUFAQtUt2pPOa6jH0WrJ4PB6falAYE+L+/n56eno0Q2EVieTk5ASMLMPJbQTKafi7diUnJ3P48GGqqqrCrsxNZKaNSKjahxMnTpCWlhawMQugs7OTxsZGli1bpjUyTTSuHheDKutubW1l9erVpKWl4Xa7tfBZ9UUol2xV4lxQUHBKtWKiFlNFSiQJx0jHSUlJobi4mOLiYmBsW7uvr4/Ozk6OHj1KamqqTwm5qrCNVJyDFXc98MADPPTQQ1RVVcX8+qaaaSES3rUP/j0XCrfbTX19PXa7nY0bN4YV3qulQiwXg0qI1dTUIIRg48aNmEwmn3AVOOVOqEqcT548qXlXFhQU4HQ6E6LsOJF2SMIZJzU1lZkzZzJz5kwgcDVoRkYGbrc7pnkp0RgZGSErKyuqMRKN01okAtnKmUwm3G63z8+pzs05c+awdOnSsP8A1AUey0VptVrp7e1lyZIllJWVhf047xJnVXjU19eH2WzG5XIxODio5TO8twMni0SJRiC6iCRQNWh7ezsjIyPs27dPMxTOz88nMzMz4nlardZpsdSA01gkgtU+JCUlaf0YAG1tbTQ3N7NixQpycnIieg7vJq9o6OzspKGhgezs7IACEe5FogqPVGm4EGO2dBaLhba2Ns0cpqCggLy8vLjnM7znFQvxXG5EgqoGzc/PRwjB/PnzNVE+fvz4KdWg3ifEB8NqtZKZmRn1nBKJ01IkQpVWq7u/y+Xi0KFDABF3bvqPFc386uvrsdlsrFmzhvr6+ojHCIZ6raqGQG0HKp+HEydOkJSUpEUZgVq2E2WpMJWRRKhxvEW5rKwsYDWov2gEGmuyxDrenFYiEU7tg8lkwmq1sm/fPioqKoImMMMhmkIom83GwYMHKSkpYcmSJVohVzwxmUw+RwU6HA4sFot2gE96eromGnrd3RJNJKK5CfgTaAsUAleDKtE4cuQIDodDF2/QQAghLgd+CpiA30gpH9Jt8DA5bUQinNoHKSU9PT2YzWY2btwYc+Io0uVGT08P9fX1U2apr0hNTfXJZ4yOjmKxWLTQ2e1209PTQ3Jyckz1IXosNxIxkpiIpKQkcnJyyMnJoaKiAo/Hw+DgoJYI/drXvobVamXXrl2cf/75Ye2iBUIIYQJ+DlwCnAT2CyGel1IeimrAKEl4kQj3SD2n06l1bs6cOVOXzHK4kYSUksbGRvr6+k4pzJpqPwnv7ksVOh88eBC73a45ZAczu9VrDqHGSISchMLtdke9LM3LyyMvL4958+bx4osvcskll/Dee+/xzDPP8Mc//jHaKW0EGqSUxwGEEE8DWwBDJBSBbOUCoZqAFixYQGZmJi0tLbo8fziRhMPhoLq6mtzcXNavX3/KHOPhJxHLeEKMnYFRWlrKjBkztPJmi8VCU1OTVv6s6jOCXXyJtNxINLFJTU0lOzub733ve7EOVQq0en1+Ejg71kEjJWFFItzlxYkTJzCbzZx11llkZmYyPDx8yhZotEyUuOzv76euro7FixdrRTv+BBOaqSzd9b44/c1uVdGRsqBTblIFBQU+W4F6FVPptdzQY5xgOYlIUYY6OhDoRU16407CiUSg2odA2O12ampqyM7O9uncjHZHIhDBlhtSSpqbm+ns7NTayoMRDzGIZ47Dv+hI5TNOnDjByMiI5vGgortYSLQIQK9xhoeH9RKJk0C51+dlQLseA0dCQolEOL4PAL29vRw5ciTgHVxPkQgUBbhcLmpra0lNTQ3btep0tq/LyMigtLRUc5MaHh7GYrEwMDBAbW2tFoVEks9QJOIWqF6RhE7VlvuBRUKIeUAbcBPwT3oMHAkJIxLh2Mp5PB4aGhoYGBhg3bp1Afeng5VlR4O/4AwNDVFTUxORKe506gj03gocHh5m7ty5OJ1OzbdSdV4WFBQEbaLyJtFEwu126zKOXssNKaVLCPFl4H8Z2wJ9TEpZF/PAETLlIqFKjoGQduijo6PU1NRQUFAQMEGoSEpK0i0n4Z0kVJWbq1atmjY1+bGg+mSys7O1fIYSjK6uLo4ePUpaWpomGoHcsRNtmaDXOHr2bUgpXwJe0mWwKJlSkVDLi0OHDlFZWRm0bLq7u5tjx4751B8EQ++chFpeuN3uqCs39USP3ZJ4GcakpKT4nMPh746dlZWlFXWlp6cnZCShx3JDx5xEQjBlf/H+R+oFuvur8mar1cqGDRvCamTSM7x3uVw0NDRQWVlJeXn5tFo66MFE70dGRgYZGRnMmTPHp7RZVSmmpKSQlpaG0+mMqUox0RKg06kDFKZAJALVPgTKI4yMjFBTU8OsWbNYsmTJpF+gXV1dtLe3U1ZWxty5cyf1uU8HorHU9y9tbmxsZGRkhOrqagAf051ILtZE3AI1RCJKgtU++EcS7e3tNDU1sXz5cnJzcydzipqt3fDwMJWVlVO+vAjEdPC4TEpK8jmHQzlJqaVlamqqtjTJysoK+VyJmJOI5rS3RGVSrgD/0mr/X4QSCZfLxeHDh6ds/W+z2aiurqawsJC1a9fS1tamW35DLybLpGWyxlB/C/5OUsoUpqWlRVvjK9Hwr0vR6+IGfd5fq9Vq5CQiIZwj9UwmEyMjI5w4cYLy8nLKysomfXmhai+WLFmidVPG6icRilheX6JEErESqi7B3xRmZGQEi8WitWorO33VPJVI+aLh4WFjuREu4dQ+SCm13oF169bp5uYT7p1OHcrT29t7Su1FvM7MkFLS0tJCf38/hYWFFBQUnHZu2aCPEW44EYAy2cnKytLyGYODg1gsFlpbWxkZGaGxsVGz059qHwdlUjNdiItIhHvmhdPppK6uDofDQVlZmW4CobYJJ/ojdjgc1NTUkJWVxfr16wOeLK53JOF2u6mtrcVkMjF79mzNLVtKqVUvhlOINNVM5ZLFu+sSYN++feTm5mrO2MnJyT6mO5MtrDr2biQEuotEuKXVqqx3/vz5mEwm+voGdNs3D8ebUj3/woULtT6FYOPohdVq5eDBg5SXlzNnzhwcDofmlq0KkZSbszKKKSwsjOjw3HBJhJyEXrsSQohTTIQtFotmIjxjxgytqCse76U/qsdluqCrSIRzpJ6UkqamJrq6urTOzfeO/4MbrzWzqPgdXvrFJpKWLo5pHqEubiklra2ttLW1ac8fzTiRogxpVqxYQW5u7injehciKaOY3t5enzV4QUEBHo8nISz19UCv+gZ/0tLSfPIZVqsVi8XCsWPHsNlsPi5SqvZGL8ECIycRFJWghFN3LxQqvM/MzGTjxo2YR81845Vv8Oj+JyDzl7S981muPutVfv+x75L/0L8h16yJai7BirNcLhd1dXWYTCbN2j4UelU32u12Tpw4EfZJYd5GMeXl5do5mBaLhe7ubk2EA53JMVkkyg7JRHj7Var3cmhoSDMRVoci5eTk6DaX6eSUDTqKhMo9BLuoLBYLhw8fZtGiRVrZ7nfe/A5P1jzJp1d9musvLqHl3V4e/M8NOA/dAzHcwQNFAMPDw1RXV0fkexlrJOFyuaipqUFKydq1a6NOqHmfg5mTk8PAwADp6ek+4bTKZwRqeosHiSISkYq490HE3oci9fT0MDw8zAcffBD0UKRwGR0dnZRlzWSh63Ij0J1XWbup3QvvP+Jvnfct/m3Tv1GaXkp9fT0XfTWLf/o8pJn+ioyhRsL/4lbFWStXroxI4WPZ3VCVhBUVFdhsNl3vmElJST4elmp78PDhw5od3WTY6+uRk4h1uRGr0KhDkWbMmIHT6WTx4sX09fXR1tbG0NAQGRkZPvmMcHfMEinxLIRIAoSU0i3GXkCalNIW7uPjugVqs9moqakhLy8v4O5BafbYHd1ut2vLg7FoPLZpqTJvj8ej9QlEU5wV7e6GqhpcuXIlOTk5tLS06JYHCGSP57096G2vf/z4cZKTk7VtVu9OzETYSk3EowKDHYrU0NCAzWbTTHcKCgoC9hJNxhIqXIQQQo794X0CaAYOAJcD24QQL0opnw9nnLiJhNls5ujRoz7FScEIlkOIlqSkJKxWK4cPH2bWrFkRndrlP04kIiGl1PwuvBvSJtMx299e32azaf6VKutus9lwOBxTcvKXN3okC+PZAeqdzygrK9PyGX19fVpncKCoTQ+hEEL8kLGL2wE0Av8spewf/97/Az4DuIF/lVL+b5BhksZ/5ivAd4UQacDXgDeBLwkheqSUf59oLrovNzweD0ePHmVkZCTsRJ3eImGz2Th69CirVq2K2s4cIltuOJ1Oqquryc7OPuWk8mAiEe0fUiSCk56ezpw5c7ROzKGhIerq6jhy5AiAdleMtKlKD/QIyyezb8M7n1FZWalFbX19fZw4cQKAXbt2aRYDMbYVvAr8v3Hjme8D/w+4RwixjDGHquXAHOA1IcRiKWWoC8gFOIEvA69LKb8vhNgAhHUH1FUkrFYrBw4coKSkhKqqqrAvAr3CMyklx44dY3BwkEWLFsUkEGpe4UQSQ0ND/OPAeyxZtJRZs2YFHCdey41IH5uTk0NGRgbLli0jKSmJvr4+bXmUlpbmU5sR77A5kZYb0XSA+kdtVquVf/zjH/T29rJ+/XrOP/98fvrTn0Y1HynlK16fvgtcN/7/LcDTUko7cEII0cCY9f47AYZRb+6+8cdfwFhUAZAHjIQzF11FwmKxUFVVpVXCTSZ2u53q6mry8/MpKyvT5Q8nnEiis7OT3z/3Cvfdfz7f2PAo3/zkIjwXXQTFxfzjf99h61fd/OCjL7B2WxesWwezZsEUr1nVxZmcnOzTVKXqCdT6OycnR6tc1PNUKv95TPUYoI/YZGZmcs0117Bnzx5effVV+vr6Yp7XOLcDfxr/fyljoqE4Of61U5BSusb//bYQ4iLgV1LKeiFEDvALoCmcJ9dVJMrLy3G5XHoOGRYqs19VVUVRURFNTU26LF9CRRJSSo4ePcqzR5/lB+2/w9X/WVzVVlJevQ0pBHLtWv6FhzG3bmLe858j8w9jdxRZUoJz505Yuzbm+emN9yE+3v0RLS0tWl2G2hrUAz0uzER1pVLvVyguvvhiOjs7qaurq/X71r1Syl0AQoh7GVsu/GH8e2Hb7AshtgK7GIsgLONfqwAGgT9LKR3hvKbEM0uIAFW92d3d7WNtr1elZLBx1IE81hQrPzjyAzYtOp+W+S4q7n4Ix8rrSHrlFboPHKX+7VIWb3mOjDt+yZDTSUZ9PeLgQWR5uTb/SO+Ck5UE9e6PmD9/Pk6nUztfdHBwELvdTltbW0ylzom03NDTSyLcvo3XXntN/XdFoO8LIW4FrgIukh/+0iOx2b+UMZH4DlDBmJluKmPXfY4QIk9KOeGSI6FEIpI/GqfTSU1NDRkZGT7nboB+ZriBlhuDg4PU1NRoRWGvFr3KutnrSL5l7Pkk63CvW4cYtfAv//cgW6q24BmagauqCvfFF384kI5lwJNBSkqKdh6HlJK9e/dq9oIOh8Mnyx/JUYHTUST0KMkWYwcF3wOcL6W0en3reeCPQoiHGUtcLmIs53AKUsovjv/3/FjmovvuRrSou3Y4Id/g4KDWHBYoUWgymbQS8VjwX260t7fT3NzMmjVrtLvF2aWBT10ryCjgZ1f+EIADBw6cIjank0D4oywHy8vLKS8v9zkq8MSJEyQnJ5Ofn09hYWFIV6lE2wLVI++io3XdI0Aa8Or4e/SulPLzUso6IcQzjJ0H6gK+FGxnY7yISjK2M/JHKWXT+NfzgBuklL8OZyIJE0mobdBQIiGlpK2tjdbWVlavXh00rNMrklB/wOqOabfb2bBhQ1RFWXouERLBLdsb/6MCVRemcpVSLtkFBQU+W+Kn2xZoOOh45sbCEN97AHggjDE8QohUYBUwRwjxUyAXeIixZcvpKRLBcLvdHDp0CCnlhM1ZenZvejwe3nvvPYqLi6M25E2ULdDJwr8Lc3h4mN7eXu0Uc1WbkWiRhF6Jy0TqAB1PTt4khPgRcBDoBj4jpXw93DESZrkRSiRUH0R5eTmlpaUTPo9ekUR/fz9Wq5UlS5ZoXgXRMJkVl4mG8HLJrqys9GmoUv4aqmzc+0DicEnESCKRREIIkQ58krFk55PAfGApMDUiEQvBRKKzs5Pjx4+zYsWKoIf3BBor1kiitbWVkydPkpmZGZNAQOItN6YS1VBVVFTEwMAAVVVVWp+J1WrVajPCtfTTa+mkVyRhtVon3PqcZMqB24AvSClrhRC5wK+EEOdKKW8KZ4CEFQmVB7DZbGzYsCGipFIsyw2Px8OhQ4fweDxs3LiRvXv3RjWON2faciMSvA8k9vZ6CNfSLxEjiUSyrpNSHgPOAxBCpEopB8a3VpeHO0ZCioTNZuPgwYOUlJRElQeIViTU886aNYu5c+fqdkEGa6F3uVxxqWSciKnwcQgHf6+HYJZ+amkCYxe3Hkcv6CUSw8PDCWc4I8ZOJV8JzBBCuBhr+joe7uMTLiehbN7COfczGNGIhKrajOV5g+EvEm63m7q6OgYHBzVDmcLCwrAMcKdbJBGKQJZ+yoZOWfq5XC5dfl96LjcSISeh2sSFECuArwPnMlZM1QR8BPgeY63jE6J7JBFtaJ2UlER7eztutzvs7tFgRNJVquztOzs7TzHF0Qvv90RFK7Nnz6aqqgq3243FYqGzs5P6+nrNZaqwsDDoezDVOYmp8EzwtvRTZeMDAwMcP36c5uZmOjo6NLGNxlFqGi43BGM1EucBLcCngE9KKb8ihPgXIOxEW0IsNxwOB21tbaSnp7N+/fqY/wDDjSTUHT0pKemUqk1vYr0olEgoh+6lS5eSn5+Pw+EgOTnZ5245MjLis12owuupaOVOZFQEpgyCs7KyTnHIjsTST89zQBNsuZEEjDImCnnjX8sDAlvEB2DKRaK/v5+6ujqKiopIS0vT5Q4Vzu7G6OgoBw8epLS0lPLy8qA/py7wWEWip6eHvr4+zaE7mL+EcpmqqKjA5XLR19dHV1cXx44dIyMjg8zMTF29N6IhkdyXVASQmpoas6WfHq8pUZYbfNj01QBkM1YjcY0Q4i3GKjX/J9yBpmy5IaWkubmZzs5O1q5dy+DgIENDQ7rNIZRIqLzH8uXLJ2xrD+cMj1BIKbFYLEgpI67W9G7lVlZqbW1tWCwW9u/fr4XXkUYZsV7kiSQSgao2o7X004NESVyO5yOEt2uVEOKrjHlPNEgpm8Mda0oiCZfLRW1tLampqWzcuJGkpCRGRkZ0u0MG+6WrrlGz2Rx23iMWM1yXy8XBgwcRQsR8QrmyUps1axYej4eFCxf6GMZkZGRouYx4O2brkRPRK68STtVmOJZ+TqdTF0s/m80WUz5NT7yEQo5/PkQERVSKSReJoaEhampqqKys9DmeXW8LO3+UMKWlpQU05Q1GuO5U/oyMjHDw4EHmz5/PyEhYBkBhI6UMGGX09vZq4bUqfc7Ly4tLLiNRtlCjifICWfpZLBZqa2vxeDwxW/olUu5I6vBGT6pItLW10dzczKpVq05Zt8VTJFRZ99y5c8M+c0MRzXaqWs4ot+zjx4/HtZjK27BVhdd9fX2YzWYaGhpIT0/X7qR6RBl61VkkQhGUsvRLS0tj7dq1Wh7I39IvnLLxRFqGCSEy/VrMvb8nIhGPuOQk/HG73Rw+fBi32x3U2l6vfgt/lGu3Ol4vGO+2vctbLW+xee5mNpVu8plXuO+nyrN0dXX5LGcmu3fDZDJppc+qvqC3t5cjR47gdDqx2Wz09fWRn58f1QU23WznvG8CwSz9Ghsbw7L0SyCh+JIQoh/4rZTSAz7isEUI8Z6U8mQ4A8U9klBNPGoXIdgbmJycrKtIyHF7+76+Ph97+0C82/YuVzx1BQ63g1RTKi/f/LImFOEuN1Q5N8CGDRvY17FPE51ZYtaU1TZ41xco74d9+/ZpJ3BHczDxdBSJYDsewSz9WltbAXxqM/QUCCHE14EfAsVSyp7xr4VrpQ9wM2NbnflCiMellN18WDvxReBbjLWLT0hcRaKrq4uGhoYJ7+Kgb3u30+lkdHQUl8vFunXrJvwjeqvlLRxuB27pxuF28FbLW5pIhBNJ2O12Dh48yMyZM5k7dy572/f6iM7/XPA/rJ+5XpfXFmtUYjKZSE5OZvHixQghtFyGOpjYO5cRz9O/EqnnItwxvC39AB9Lv+7ubh5++GEAmpubqaioiHo+Qohy4BLGiqDU1yK10u8A7mVMVKqEED+UUh4Z/54JCHsrMS7LDY/Hw7FjxxgeHp7wLq5NJDlZFxNddeZnSkoKixcvDuuXv3nuZlJNqdpFvXnuZu17E0USys5OmfDCqaKzr3sf60rWxfza4oF/lNHf309vby+NjY2kpaX5tHErEi2SiHWcaAupvC39qqqq8Hg83HvvvXzhC19g/vz5PPLII9FO6cfA3Yz5UyoisdKHsSjiH1LKl4UQvwAeEUJ8V0r5JmOOV2Fn03UXCZvNxvvvv09hYSFr166N6OyNWENy1Va+cuVKLQcSjkhsKt3Eyze/HDQnEUwk1PN529nBqaKzaeamKS+l9ifQ78V/q1DlMlSvRF5enlYurodZjF6heSK4WyUlJVFRUcHChQt57rnnov59CyGuBtqklAf93p+wrfTHyWa8oEpK+UUhxD8Bjwoh7gYygeFw56S7SIyMjLBgwYIJj/bzJ9binmPHjjE0NKS1lUe6fNlUuslHHBSBxlH5jsHBwYBt7P6iUybLAuZbor2bTpbgZGRkUFZWRllZmeZj2dvbi8Viwel00traSmFhoU+UES567W7oQTxcqUL9XpWVvsLLUv9e4BuMuVz7E7aV/jg1jFnnI4QwSSn/KIT4O7CdMXfu0RCP9UF3kSgqKtLFhDZclL19bm6uT+SiV47DP8JxuVzU1NSQmZkZMlLyFp22tjafMaSU2oHGSjySkpIQQiRsF6i3j+XIyAgNDQ0IIbSDfFTZc35+flgXXALtAuiWH7FarWE1d3lZ6StWAAghVgLzABVFlAHvCyE2EpmVPsBd40sT5Php4uNGuOuFEJ8Ltj0aiCnv3YgF5Zq9cOFCSkpKfL6np62+Ehu1U1NRUeFTCDYR3kLjLRApKSlIKXG73dpzuFwukpKStI9EJTk5WYsyPB6Plss4ceIEKSkpE1rSJZJI6H0wT7RIKWsA7Q9ZCNEErJdS9gghwrbSHx+ry+9z6fX/sAxwFaetSLS3t9PU1BSwMAv0sbCDD3c3VLNQODs1wcbwFgghhHaRqD9QFVmof4NFGVOd3/C/wJOSknzcspXvg6otyM3NpbCw0CfK0OvurQeng7+ljMBKX28mpZgqksdO9AvztrULVpgF+i43zGYzIyMjMflNqAtfXWCB3ifv6EGJifpQj1P/n0omen5/SzqVy1BRRrj+lZNFojplSykr/T4Py0pfbxIqklCl2cFEQh0KXFhYOKGtnR4i4fF46O7uBsYKpKL9QxJCMDAwEPaaFU4VDLUs6enpISUlBYfDoR2SM9l35EiWCsr3QZ3wbrPZ6O3tpb29Hat1bFmschl62NBFg545iQRpE9eVhBSJQHcZZdiyePFirWQ2nLGixeFwcPDgQVJTUykuLo5KINTyIjc3F6vVqhUsFRQUUFRUFHbzlRK848eP43K5WLJkCYC2JIk0+akH0UaM6enplJaWkp6ejsVioaioiN7eXpqamiJu4dbzbBW9lhvh/G2ebiTUciPYhX3y5ElaW1s1w5ZwiCWSGBoaorq6msWLF2Oz2aIK773zD8nJycydO1drvrJYLHR1dWl2darPIljRmdpRyc3NpaqqSnuPk5OTtSjDe6fE7XZrYhEPwdCrmMpkMvlEGXa7XVuWjI6OkpOTo+UyAkUZem2jJuARfwlFQkYSCo/Ho7U+T3Rqlz/RikRXVxeNjY2sXr2arKws2traIt7SDZaghLHX6N3iPTw8TE9PDwcPHgSgsLCQoqIizafRZrNpBxPNnj074OtU46akpPjkPlSUoS5IvaKMeFVcpqWlaS3cqk+it7eX5uZmbQu2sLBQizISqbQbDJGIiGirJwPZ6s+cOZOKioq42+pLKTl+/LjWEKbuLJG+FnVxhkpQKoTX6VbKQr6np4fm5maGh4fJyMhgaGiIZcuWhV2cFiyXod4LNbdY3bZiZSKh8e+TUOeLKqOYnJyciHeZgqGnv6UhEnFGiURfXx+HDh1iyZIlEVdu+o8VDm63m5qaGs1TwPviiURs5PhZGupxkZKSkqKdoan8DIqKimhoaKCpqYmioiKKi4sjWnKB7xZrS0sLKSkpMecy9CjLjuT5vM8XVVFGV1cXQ0NDvP/++1qUEeoUc73mEowEcsrWlYQSiaSkJDo7OxkeHmbt2rVhty4HGyucZYIyxFWFQYHGCefOqe7WE0UP4dDS0oLZbGbjxo1aRGOz2TQzG9WtWVRUFLYnhJSSxsZGRkdHWbt2rfY1/ygjHMGY6gYvFWWkpKTgcrlYuHAhFouF5uZmLcpQuYxwcg3T2ClbFxJmueF2uzGbzZhMppi2GxXhRAAqYlm2bJmWPPNnoi7QUPmHSJFSUl9fj8vl4qyzzvK5UNPT0336KJTz1NGjR7XzSpXjuD8ej4e6ujrS0tJYuXLlKTkS9TNKMLyToOr73nOZapFQqAjA/xRzlctoaWnRCr1CRRnGFmhoEiKSUHfzzMxMCgoKdFH1iSouT548ycmTJyeMWEKJjfeugroDR0uwHYxA+DtPjYyM0NPTQ21tLW63W0t+5uTk4HK5qK6upri4mLlz54Z8nYEKubx3TNRr1CMnoceFGWh3QwihHRU4f/58HA4HFouFlpYWzclabbOqKEPPI/4MkYgDylZt2bJljIyM6OIpAcF7N1TFpt1uDytiCbbc8E5QxioQE+1g+LO3bS9vtb7F5vLNnF16tmYfX1lZqRmhtLa2MjAwgNPpZM6cOWGNqwiU/PQWDqfTqUvyczLazf3P41BRxsmTY6ZMBQUF2O12XUTC4XAkjFO2nsRtuTERctzevru7Wyt3ttvt2O12XeYQKAJwOp0cPHiQ/Pz8sA8iDrTc8BeIWBgaGqK2tpYlS5YEXfJ4s7dtL1c+c6XmVbH7ht2cXXq29n1lhJKZmUlNTQ0LFy7EZrPxwQcfaH4RKvkZzuv3T34ODg7S0tKiHVEYbfJTz+VGuHhHGYAWZYyOjnLgwIGAUUakJErTmp5M+bkb3sfr6Wlh5z+WcqxasGABM2eGfcLZKZGEd4IyVoEwm81aTUa4OxZvtfpZ7bW+5SMSgGYUs3r1ap9su91up6enh4aGBkZHR8nLy9OSn+Es8fr6+qivr9fGDbbFGk4h12TazgVDRRmtra2sW7eO4eHhU6KMcM8WTaSuVr2ZdJFQ7dbl5eWn7CaYTCbdlhveW6Aqwbdq1aqIs8/eTVV6JShhbAeju7ubtWvXRnQgzOZyP6u98s0+3+/o6NByLf7jpqWl+TRe9fX1aaKRnp6u5TkCNbF1d3dz4sQJ1qxZo30/0BZruIVcUxFJhCIpKYmcnBxycnK0mhXvs0WzsrK0KCPY72u6CsWkisRE9vZ6tXfDhzmJEydO0NPTE7bXZqBxvNfjeu5g+NdkhMPZpWez+4bdPjkJNW5TUxP9/f2sXbs2rFyLt1WdSn7W1dXhcrm05Gdubi7t7e10dHSwdu3akGH4RIVc3l4ZiSYS/nj7V8rxA3x6e3upqakBPnTJzsnJ0bX6UwjxFeDLjLWD75ZS3j3+9UicsnVlUnISqprRYrGEvFj1PqBnYGCAtLS0sByzgyGEYHR0VDu+TY8djJycnAl3MEJxdunZPksMlYyVUrJ69eqoXqs63EcdVGyxWGhra9OOKVy4cGFE4wWLMtS/drtd+3+0S7fJ8qQQ4wf4+EcZbW1tHDlyhIyMDP7617/GvLMhhLiAMcPbVVJKuxCiZPzrkTpl60rcIwl1YWRkZEx4seolEmq3wGQysWzZsqjHkVKSkpJCSUmJdgdRVY+RHjQb6Q5GuKhq0dzcXCorK3UJd9UBNapLs7S0lN7eXj744AMtAon0PfCOMlpaWjRj3UgLubyZKuMa/yjDbDbT1NREfX09H/nIR/jMZz7DZz7zmWiG/gLwkJftXPf41yN1ytaVuIqEShbOmzcvrAtDD5Hwbik/fvx41ON4JygrKyuprKzE4XDQ09PD8ePHGRkZIT8/n+Li4gmrHiPdwQgX1c5eWloakZ3eRHg8Hmpra8nMzNQinry8PBYsWHDKe6CSn+HUt6gl0eDgoFYs5i0SgQq5Qi3vEiEHIISgpKSEO+64A7PZzG9/+1s6OjqiHW4xcJ4Q4gHABnxdSrmfyJ2ydSVuyw3VTbly5cqwk4WxikR7ezvNzc2cddZZpKen09DQEPEYoRKUqampPl2K3lWPwVq+1Q7GqlWrdK3rHxkZoaamhkWLFkXd3xIIVXxVVFQUsPjK/z3o7+/XTgNLS0vT3gP/AjXlaO50Olm5cqXPjpb3v4EKudT3/aMMvXZI9EAVUnk3pQXC3ykbNLfsexm7HvOBTcAG4BkhxHwid8rWlbiIRF9fHydPnmT9+vURJQujFQkpJUePHsVqtbJhwwaSk5M1T8lIxwk3Qemd+FNVj2azWVvDFxUVaYfPRrqDMRH9/f2a36aevQIOh4MDBw4wd+5cZs2aNeHP+3tbWq1Wenp6OHz4ME6nUzPXycnJob6+Xlv+TfS+BirkCuSVoZdI6NW3EaVTNnzolv0FYOe4ae0+IYQHKCJyp2xdiYtI5OXlRXQwjyKa0NHlcnHw4EFycnJYs2aNNkakY0XS4u2PEEKrepw3bx52u53a2lpGRkZITU3lxIkTFBcXh+1EFQq1FamiJb1QpfGxRCaZmZk+5jrKpu7AgQOkp6dTUVGB0+kMWzAn8soYHR0lMzMzIZKfOnWA/hm4EHhDCLEYSAV6gIicsvUmLiIRa5lyuIyMjHDw4EHmz58f1p0vGLEIhD8ul4vDhw9rQqmWJcqJKisrSwvJI63qa21t1Wor9DSSVTmT5cuXk5OTo8uYqrrz5MmTzJ8/n4KCglPMdYqLiyNq7faOMpqamnC5XBQUFERcyOVNgnWAPgY8JoSoBRzAreNRxZQ5ZUMC9G5Ei2qbXrlyZUx/2MoDQo8KykA7GP7NWMPDw5jNZq1MWn0v1F1Ijp8YZrPZTukOjRWLxaIVmumZM1El8Cp/AWjmOg6HQ/O1HB4eJjc3l6KiIgoLC8O6YE+cOMHg4KDPdm+0jlyJFElIKR3Ap4J8b0qcsuE0FAkpJc3NzXR1dbF+/fqoG2r0rqAMZwfD24lq/vz52O12LfGpDHKLi4vJzc31+eOvq6sjPT2dFStW6BqhdXV1aYlePRuTVG6jsrLylEOTYCz56W0gMzAwQE9Pj2a5r4TTv1Rd1duMjo76JD8hskIuvZOfMH29JGAKG7xCEWxry+PxcOjQIQCfno9oxtdTIFRpc6R347S0NB+PCIvFQkdHB0eOHCE7O5v8/Hw6OjqYOXMm5eXlEw8YASdPnqSrq4uzzjpL16VLpLmNpCRfy/3R0VEfcx2V/MzNzaWxsRGXy8Xy5csnTH5CeIceOZ1O3ZYbeta/JBIJF0moHQ5/d2S73a55Xs6dOzfsC9tfcPQWiNbWVrq6umLewfA3yO3p6eHQoUMkJyfT3d2NlJLi4uKY3Lrgw7vx8PAwa9as0eUCUYyMjFBdXc2yZcui9p/MyMigvLyc8vJyH2fx6upqkpOTteVKJJFPMK8Mj8eDzWYDiCn5CdPX3xLiKBKxmuF6i8Tg4CA1NTVUVVVRVFQU9lhqm0xdCHomKNW2q8PhiKoHIxTDw8M0NDSwevVq8vLysNlsmM1mjhw5gsPh8FmWRPIapJQcOXIEKSWrVq3SdekyODhIXV0dK1eu1O1iUTmb7u5uzRND9U94PB4t+RlOl6bCWzD6+/s5efIky5cv124eKj8V6aFH09WVChI4klB0dnZy/Phx1qxZE3FiyFskYjWp9UaVQmdnZ7N48WJdLzbV5u29dElPT/e5u/b29tLW1sbhw4fJycnRkn6hTsDyeDzU1NSQlZXF/PnzdZ2zdwt5uC3v4aDyMZmZmdqcs7OzNXMdZVE3NDQU9vugGBgY4MiRI6xZs0aLzryLuCL1yoj1sOBEJmFFQmX0BwcHfSzuoxlLiYUeywu17CkrK9O1FBrGKkbb2tpCLl1MJhMlJSWUlJRoTkuqdyAlJUVbsnjXUKhakpKSEt1zG6ra0ruFXA+UqKmeFH9SUlJ8HKdU8lOdBBZq12hgYIDDhw+zevVqn+WbijKiOfTIWG5EQbTLDZVMOnDgAJmZmVEVZXnPweVyTdgDEC5qB6OqqkqrMtQD1dMwMDAQVpu3wttpaeHChYyOjmI2mzl06BBOp5PCwkLy8vJobGyksrIyIrOdcOjs7NROVtOzotTtdmtnvoby5VSo3pK8vDzNiaunp4djx45hs9l8nMWHhoYCCoQ/+zv2+7TjT7TFauxuTCJSSurq6liwYEFMd2opJcnJyTQ1NTF79uyYD3KJdgdjIjweD0eOHEEIwapVq2JaCmVkZGgVjy6Xi/b2dmpqakhJSaG3txchRNi1CBPhvTsyUXjv78kZCrfbrUU9gY44CAdvZ3HvHhuV06msrAz5PgezCAy1xdrY2DglHamTQUKJhMViwWw2M3/+/JgFwu12s3jxYs1d6PDhw+Tm5lJSUkJBQUFEv1C9djD8UXfMvLw83dq8FVarlba2NtavX09WVhYDAwOYzWaOHz9OWloaxcXFQV2oJkKZ24SzO7K3bS+X/+lynG4nKaYU9ty4J6hQuFwuDhw4QGlpqW7biarHxmQy0dfXx+rVq7WI0O12a0lgZR4DE1sEem+xSil58MEHmTt3bsCakOlAXJcbkdDS0kJ7e7t24nS0eCcovfv+Vcei2Wzm2LFjzJgxQ7tQguU7VOei3W7nrLPO0nW70G63U11dTVlZme7766qK0juRqMLxRYsWYbVaMZvN1NXVaRb84ewSqDyR3W4PO+r5Q90fcLgdADjcDv5Q94eAIqGWmHPnztV9WdTf38+RI0e0fpeCggKtj0TdRAYHB8nOzqaoqIiPzP5ISItAhZSShx9+mIaGBl5++eWwEqanI1P+qlS47XK52LBhA62trVG3i4c6Rcu7YzFQeXRxcTElJSWaQLndbmpra5kxY4bulY7xavOGsTxBS0tLyCrKzMxMKioqtAtFHco7PDwc1B9CbZ8KISYsZvJG+HU5+38OH1Zozps3j+Li4ghe7cSonZdAiVV/85jBwUF6enpI6U3hv1b9F8ccx7h40cVsnLPxlHGllDzyyCN88MEH/OlPf5q2AgEgJkguRt2z7nK5JrzYlWlKUVGRFm63trYipQwrYaVNMsYCKVWH0N3djdvtJi8vj97eXubOnUtpqb7eHvFq84YPG8BWr14d1R+ttz+ExWIhPT2d4uJiCgoKaGhoICMjgwULFkT0/u5t28sVf7pCuyu/fOPLPpGE3W7nwIEDLFy4UHfBDCUQE6GcxXt6erBaraccq/jrX/+a//u//2PHjh16LEET2j03biLhdrtDOl8PDQ1RXV3N4sWLfe4e7e3t2O125s2bF9bz6F1B2dfXR01NDZmZmZongt5t3qtXr9Z1u1Cd82m1WlmxYoVuCbSRkRG6u7tpbm7GZDJRWloacecmBE9c2mw2Dhw4QFVVla6OXRCbQPjj7Sz++OOP87e//Q2bzcauXbsi9v4MgiES/ijXqlWrVp2yt6xOig7nzdezghI+3MFYuXKldq6ESqb29/eTnZ2t5TEizU+oQ4BXrVqla6+EWq4lJSXFZK4bCFVfMWvWLIqLi+nt7cVsNkdk3RcMq9VKdXU1S5YsCenkFA0qJ6N37QbAk08+yTPPPMMVV1zBq6++yk9+8hOWLl0a67CGSGiDjfcNqCxzoIulp6eH3t5eqqqqQo6vZwUljG3pdXR0sHr16oDho3fhUk9Pj7ZDUFxcHLKPwDv5uXz5cl23yVTeRLVg6ykQKk9QUVFxSiLR+85qsVjIzMzUxDOc0Fv1eOjpX6FQAqF3ZyvAn/70J373u9/x4osv6l1deWaKhDozUqFKmdPS0qiqqgp6sfT19dHR0RHS5TpUgjJS1EVss9lYvnx52BGCsqvr6enRmq+Ug7RCGcpmZGSwcOFCXS9i5dcwa9asqOsJgmGz2Th48GBYeQJv676enh4gsKO4WnKsL1pPujld1x4PhcVi4dixY6xZs0Z3gdi5cye//vWv2b17dzyKpgyRUGctBjq1y5/BwUGam5tZuXLlqZPROf/gvYMRaULOG4fDgdlsxmw2Y7PZKCwsJD8/n6ampri0eauLeN68ebrvzce6DFBu2mazGavVSkFBAS2eFm7eczMOt4NkkczOrTv52IKP6TrveArECy+8wM9+9jN2796t+9JonIQWibjv2/T19XHo0CGWLVsWVnIqmBmudx29HvZ4qk5BDzv61NRU7fg8t9tNR0cHtbW1mEwmhoaGMJvNYVnOh4PaPo1Hsk8VGcWy8+Lvpm2xWPjtO7/F7rLjwQMS9nbu1VUkent7aWhoiItA7Nmzhx//+Me89NJL8RKIhCeuInHy5ElaW1tZu3Zt2D4IgUTC/xTvWAVieHhYO5tDzx4MGLsTt7a2smbNGnJychgYGKC7u5uGhgYyMzMpKSmJyt8SxhqTDh06FJftU1VwpGfZeVJSEkVFRVyy8BJ+c/Q3uKSLlKQUSqwlvPfeez7LkmhRAqF3/wjA66+/zkMPPcTu3bt1/zs5nYjbcmN4eFj7g47kDup0Ovnggw/YuHGsgEXvHQzViq12MPTEu807kPWa2lLs6enxMZkJR0C9TwqP1Xgm2Njx2A1QF7FjpoN3O9/VtkFVHYJaogWy7gt37HgIxFtvvcU3v/lNdu/erXsFaAASerkRN5GQUuJwOCJ+nMfjYe/evZxzzjm6Jihh4h2MWFBt3uGOrQq4zGYzTqeToqIiSkpKAtYgqJPC4zFv5XO5Zs0a3cc2m83aSeShxlYOVGazmYGBAbKysiguLqawsDBoxKVa1OMhEG+//Tb33HMPL774ou52AEEwRCLSx/39739n06ZNuu5gNDQ0MDo6GtEORrhjK/fmlStXRjW2Ko3u7u5mZGTEp4Dr5MmT9PT0sGrVKt1Lf9va2jTR1LN2A8bEp6WlhTVr1kQ0tjrB22w209vbq0Vc3sa46pjBeAjbvn37uPPOO3nhhRd0TziH4MwUCRhLDkaKlJK//e1vmmNQouxgBMK7zTvUtm6kY6q7aldXFyaTiUWLFlFUVKSrSDQ3N2OxWFi1apWuogljkU9bWxtr1qyJec7KG8JsNmO328nIyGBkZIR169bpnqR8//33+dKXvsSuXbsCGt3EEUMkwkXlH9rb2+ns7MTtdmuNV9HkD9QOxpw5c3TvwVCnpcejzdvj8XD48GFMJhOzZs3SegjCLeAKhSrhVlGV3h4IJ0+e1PpH9Bafzs5OGhsbycnJYXh4WKuADdeyLhTV1dXccccd7NixQ69S60g4c0XC4XCE7U4VKEGp6g+6u7ux2+0UFRUxc+bMsHoH4rmDoSzsvA/h0YtQHhNWq1VLfAYr4AqFlJL6+nqklCxZskT3U9ZaWlro7e2NS3Si8hvqCADvCtje3l7tvI5oHMUPHTrE7bffzjPPPMOSJUt0nXeYGCIxEd4OP8HubC6Xi56eHp91+8yZMwM6Rqts/YoVK3Sv6lN1CvEQH1VFOXv27Akjn0AFXKEctNWZJWlpabpXfwKa/Z7/oTl64C8QgVDndahEsHo/vM1kAlFfX8+tt97KU089xfLly2Oeq81mY/PmzdjtdlwuF9dddx33338/FouFG2+8kaamJiorK3nmmWe861wMkQg6eJQVlCob3t3dzeDgIHl5eZSUlJCfn097e3vcdjDi2eatqijnz58fsaeCctA2m80MDg6Sm5urtXiruhPvpZGeqH4cq9Ual+VLd3c3TU1NER0i5HK5tPdDOWmrZYl3hNPY2MgnP/lJnnzySVavXq3LfNVWd1ZWFk6nk3PPPZef/vSn7Ny5k4KCAv7jP/6Dhx56iL6+Pr7//e+rh525IuF0OrUI4ZSBdSqxVh4I3d3ddHZ2aom+4uJiXUPerq4umpqadG/zhrGlUU1NDUuXLo25qk9Kqb0fFouFjIwMrFYrpaWlVFRU6DNhr+dqaGjA4XCwbNky3aMT1aYe6Q6J/xyVdV9vby9paWm0t7eTl5fH17/+dR577DHWrVun67wVVquVc889l1/+8pfccsstvPHGG8yePZuOjg4+9rGPUV9fr37UEIlTBo1TD0ZGRgYlJSXaH4T6PNoKR0W82rzhwyrKeDQ8ORwO3n//fTIzM7HZbBEXcIVC5TcA3VvUQR+BCITVauWXv/wljz76KDk5OVx33XV8/etf1/W9d7vdrFu3joaGBr70pS/x/e9/n7y8PPr7+7Wfyc/Pp6+vT32a0CIx6Z5beldQBtrBUNbqqsLx/fffJyUlRTuvItxliL/Hpd6htPKv8D4gRi+8OznVqWeqgOvw4cMTFnCFQkrJ4cOHSU5OZtGiRboLRLQ1FuEwMDDAzp07+cMf/sCqVat45ZVXdI8MTSYTBw4coL+/n23btlFbW6vr+JNNXCMJfws7vT0g1A5GOF6RamfAbDYjhNC2VoNdnG63m7q6uri0eUPkFZqRoDo5QzWB+Rdw5efnU1JSMqEDl0qAZmRk6H4SGMRXIDo7O7nuuuv40Y9+xAUXXKDr2MG4//77mTFjBo8++qix3AiEEgm9lxcQ2w6G3W7XBMPlcvncUeHDXYZ4tHnD2E5AX19fXLYKVX4jEkMXZSLT3d3t48DlX3+gTtXKyckJ214wEpRAhHOWR6SYzWauvfZa/vM//5NLL71U17H9nyclJYW8vDxGR0e59NJLueeee3jzzTcpLCzUEpcWi4Uf/OAH6mFntki4XC7dBaKtrY329nZWrVoVc9Wd0+nUajFsNhu5ubn09fWxcOFC3b0a1PJFJfr0Xr6o4+tiaV7zrz9ITU3Vdkrq6+vDPlUrUrxPA9NbICwWC9dccw3f+c53uPLKK3Ud25/q6mpuvfVWbVv/hhtu4Nvf/ja9vb3ccMMNtLS0MHfuXJ599lnvLfQzWyQcDsdp0YMBY1ucNTU1zJgxA7vdTkFBgRaCxzp3FaanpqbGZR3vfdaGnvkNq9Wq7eykpKRQWloadQVsMDo7Ozl58qQuZdz+9Pf3c+2113LPPfewdetWXcfWkTNXJP75n/+ZgoICtm7dGnPiL945Av82b9VD0d3dzcDAQNSnf6m5V1dXk5+fH5eeAFVLEI+GJ2WGO3v2bIqKirSCtnAKuMJBzz4PfwYHB7nuuuv46le/yvXXX6/r2Dpz5orE8PAwu3fvZseOHRw5coQLLriArVu3smHDhoguNHU+Rzx6MGDiJKJ/7UFWVpa2tTpRNKPmHo+TutTc29vb49LJGepUrYkKuMIhngIxPDzM9ddfzx133ME//dM/6Tp2HDhzRcKb0dFR9uzZw44dOzhw4ACbN29m69atnHPOOSH/qCLZwYiUaNq8VSuz6qFIT0+npKSE4uLiUy7S0dHRU7Yh9aSlpYWenp64NFMpt+zKysoJczNKRM1ms1bApeoxgglXR0eHJm56C4TVauWGG27glltu4bbbbtN17DhhiIQ/drudV199le3bt7N//34+8pGPsG3bNj760Y/6/FGpdXY8ejC827xjaXZSxwWazWZMJpNWi+F0OqmpqWHZsmUxn2jujxK34eFhXQ/jUcRyqpa3e7Z6T/wLuFTpfDgHDkfK6OgoN998M9dddx2f/exndV+WxomEnuSUiIQ3TqeTv/zlL2zfvp23336bjRs3smXLFo4cOcLcuXO54oordPcNUG3e+fn5VFRU6PaHNDo6itlspr29nZGREc0d3N/KLhaklBw9ehS3283SpUt1vwjUqVp6NbD5O3ClpaVht9tZt26d7hGE3W7nk5/8JFdeeSVf/OIXTxeBAEMkwsflcvHmm29y991309fXxznnnMPWrVu56KKLdKuKi2ebN4ztkx8/fpylS5cyODhId3e3Vt04c+ZMn7MoIkVKyaFDh0hJSYnLDolaHsXjVC0YO6u0tbWVrKwsXU4B88bhcHDrrbfysY99jDvvvPN0EggwRCIy1BLke9/7Hu+++y47duzg9ddfZ+nSpWzdupVLL7006jtzPNu8IbgdnNPp1HYFRkdHKSwspKSkZMI2Zm+8C5n0NrmBD9+bZcuW6X6qFoy9N52dndoSQxVwmc1m+vr6tGRwNAYyTqeT22+/nY0bN3L33XefbgIBhkhEhurp8Mbj8bB//362b9/OK6+8woIFC9iyZQuXX3552C3bqs07Ho1UUkofP4VQ62y1K9Dd3c3Q0FBY5dBut5uDBw9SXFwclwpQVaUZjxZ4GBOIrq6uoAlW72SwMpBRyeCJlpoul4s77riDpUuX8q1vfUsXgWhtbeWWW26hs7OTpKQkPve5z/HVr351Ik+IWDBEQk88Hg8HDx7k2WefZc+ePZSWlrJlyxY+/vGPBw2R49nmrXIELpeLpUuXRhQ2+5dD5+TkaLUY6mJS25Dx2kIdHBzUulD1PmIAorOzs1qtWh5DShnw2EAYE88vf/nLlJeX893vfle3CKKjo4OOjg7Wrl3L0NAQ69at489//jOPP/54KE+IWDBEIl5IKamtrWX79u289NJLFBYWsnXrVq688kotK3/ixAnN8FXvOgKPx0NdXR3p6ekxF3gp3wN1N50xYwYFBQW0tbVFZUQTDqqMO9A5IXqgh9+l97GBo6OjFBQUkJ6ezqxZs/j6179Ofn4+3//+93Xf4fFmy5YtfPnLX+bLX/5yqCatWDBEYjJQ/gbbt2/XTn1OTk5m1apVfPvb39Z9q83lclFdXU1RUZHuvQxSSnp6ejh06BDJycnayV/FxcW6VVT29fVRX18fl8N+ID6GuMqR7LHHHuO3v/0tubm5PPjgg1x66aW6R4iKpqYmNm/eTG1tLXPnzg3lCRELhkhMNqOjo1x77bU4HA5GRkZIS0vjE5/4BFu2bGH27Nkxh6Wq0CheOySqgGzp0qXk5ub61B0kJSVpbe7RXhjxPDsTxtb06qwQvcXZ4/HwzW9+E6vVyqc//Wmef/55br75ZtauXavr88DY7+H888/n3nvv5ZprrpnIOCYWDJGYbHp6etizZw+f+tSnkFLS0tLCzp07ee655/B4PFx11VVs27aNsrKyiAVDbRPGowIUxnIEdXV1QROsNptNa3N3u91ahBFuPiHcU7WiRQnE6tWrdV8CSCm5//776enp4dFHH9VdgLxxOp1cddVVXHbZZdx1113AmAOXsdw4ldNSJIIhpaSjo4OdO3eyc+dOrFYrV111FVu2bAnLQEWduh2PKkqIfAngfeSAw+GY0GkqnoYu8KGlfrwE4sEHH6S5uZnHH388rgIhpeTWW2+loKCAn/zkJ9rX//3f/z2UJ0QsGCKRqHR3d/Pcc8+xY8cO+vr6+PjHP86WLVsCejaqCzheuwCqCCvaJYD/kQOqFkN1aMazmQrGBEIliOMhED/60Y84dOgQv//97+Myf2/+9re/cd555/kcD/Cf//mfnH322aE8IWLBEInTgd7eXnbt2sXOnTvp6OjgsssuY+vWrSxbtox9+/YhhIjLFip8eCCwXnd4/yMHUlJScLlcrFu3Lm4RRDwF4r//+7/Zv38/Tz/9dFzmnwAYInG60d/fzwsvvMDOnTupqakhNTWVn//85xG3uIdDa2ur5sQdjztkc3MznZ2d5OTkaNZ0qrJRj5C9ublZs+KLh0D86le/4o033mD79u1xyaEkCIZInK48/PDDvPHGG1x//fU8//zz1NfXc+GFF7J161bWr18f80Xh3aYej31+/1O1lDWdqsWI9ciB5uZm+vv74zJ/KSWPPfYYL730Es8991zctjgTBEMkTleOHTvG/PnztTuu8sTYvn071dXVmifGpk2bIrorex9qE2mVZrjjT3Sqlmrp7urqoqenJ+IjB+J5rB/A7373O3bs2MHzzz8flzqOBMMQiemItyfGe++95+OJEWrZoM6sMJlMLF68WPdmpGhP1fI/ckBtrQa6QOMtEE8//TRPPvmkVhR3BmCIxHTH4XDwl7/8hR07dvD3v/9d88Q4//zzfe7KHo+H2tpaZsyYEZczK/Q6NVwdOdDd3Y3b7daKt2bMmEFTUxODg4NxMbsB2LlzJ48++igvvvhiXJrNEhRDJLx59tlnue+++zh8+DD79u1j/fr12vcefPBBfvvb32IymfjZz37GZZddpvfTxx2Xy8Vbb73F9u3beeuttzjrrLPYsmULGzZs4I9//CPbtm2LiyV9vE7V8j5yYHBwkOTkZK1ORG+Re+GFF/jZz37G7t274+JnEW8CdTCHSUKLxKQf87dixQp27tzJHXfc4fP1Q4cO8fTTT1NXV0d7ezsXX3wxR48ejWvRTDxITk7mwgsv5MILL8TtdvP222/zxz/+kc9//vNs2rSJyspKioqKdG2oUnb96enpLFiwQNeLNyUlhTlz5mC32xFCMHPmTFpbWzl8+LCuRw7s2bOHH//4x7z00kunpUC43W5MJhN2u52hoaG4eJpOFfFrnQvC0qVLqaqqOuXru3bt4qabbiItLY158+axcOFC9u3bN9nT0xWTycR5551Hc3Mzv/71r/nWt77FP/7xDy666CI+/elPs2PHDoaHh2N6Du8lTDyOGgA4fvw4w8PDrFy5klmzZrFy5UrOPvtsCgsL6ejo4N133+XQoUP09PQEPUU+FK+//joPPfQQL774YlzMgOKNEoihoSGuu+463nvvvamekq5MeiQRjLa2NjZt2qR9XlZWRltb2xTOSB+EEDz77LNaH8bZZ5+Nx+PhwIEDPPvsszz88MOUl5drnhiRlHu73W4fr854cPz4cUZGRk7ZJUlKSqKoqIiioiKfIweOHTsW0ZEDb731Fvfffz+7d+/W9e57++238+KLL1JSUqId2Bsv0xiTyYTD4eC6667jsssu4/LLLwfGKnr1PgVuKohLJHHxxRezYsWKUz527doV9DGBciOnoQ1ZQPwbtZKSkli7di0PPvgg77//Pg888AAtLS184hOf4Nprr+V3v/sdFosl5JjKraqoqChuAtHY2MjIyMiESUohBPn5+VRVVbFp0yYqKioYGhpi//79HDx4kPb2dpxO5ymPe/vtt/nGN77B888/f8q5HrFy2223sWfPHp+vPfTQQ1x00UUcO3aMiy66iIceeiim53jggQe0rtD6+nrS09O59NJLef7557nlllv47Gc/y5/+9KeYniMRiEsk8dprr0X8mLKyMlpbW7XPT548yZw5c/ScVkIihGDlypWsXLmS++67T/PEuO6668jOzubqq6/mE5/4BMXFxZpoep+qFa/3qLGxkdHRUVasWBGRWAshyMnJIScnh4ULFzI8PEx3dzcffPABycnJmnlOe3s7d999N88//3xcXsPmzZtpamry+dquXbt44403ADTT3FicpRYtWkReXh6dnZ2sXLmSNWvWcMUVV3Dbbbdx0UUXab4gpz1SylAfceP888+X+/fv1z6vra2Vq1atkjabTR4/flzOmzdPulyusMb6zne+I+fMmSNXr14tV69eLXfv3h2vaU8aHo9HHjt2TD700EPynHPOkR/72Mfkj370I7l//375r//6r7KxsVGOjIzo/jE8PCyrq6vl3r175fDwsK5j9/T0yHfffVeuXLlSFhYWym984xuysbExbu/hiRMn5PLly7XPc3Nzfb6fl5cX1bhOp1P7/zPPPCNnzpwpjx8/LqWU8ujRo9r3Pv7xj8sHH3wwnCEnug6n9GPSE5fPPfccZWVlvPPOO1x55ZXaNufy5cu54YYbWLZsGZdffjk///nPI9rZ+NrXvsaBAwc4cOAAH//4x+M1/UlDCMHChQu55557ePvtt3n88ccZGRnhsssuo6amhmeffZbW1taAy7RokVLS2NiIzWZj+fLlui/3MjIytG3CF198kcrKSh5//HFdn2MyUMVyr7zyCtdffz2f+cxnuOaaazh8+DCLFi2ivb2dT3ziE1RWVvIf//EfUzzb2Jn0xOW2bdvYtm1bwO/de++93HvvvZM8o8RHCEFFRQXV1dX8/ve/Z/Xq1ezYsYMvfOELjI6Oap4Y8+bNi+lMj8bGRux2e0SVmpFw6NAhPve5z/Hss89q+YvJZObMmXR0dGimMZEmFZUgCyF47rnn+MUvfsGll17KAw88QGZmJjfddBN/+MMfmDdvHp///Oe58sor4/EyJp1JjyTixSOPPMKqVau4/fbb9bIUSzieeuoprrjiCubMmcNXvvIVXn/9dXbt2kVRURF33XUXF1xwAT/4wQ+0qstwkeOl3PEUiCNHjnD77bfz1FNPBdwCnwyuvvpqnnjiCQCeeOIJtmzZEtHjhRDae3PhhRdSVlZGf38/UkruvfdebrzxRs4//3x6e3unjUDAaVSWffHFF9PZ2XnK1x944AE2bdpEUVERQgi+9a1v0dHRwWOPPTYFs5xalCfGjh076Ozs5PLLL2fr1q0hm8iUQDidzrgcGwjQ0NDApz71KZ588klWr16t+/iBuPnmm3njjTfo6elh5syZ3H///WzdujVm05g77riDvr4+BgYGaGlp4Zvf/CZr1qxh+fLlALz55pucf/75kU43obfxThuRCJempiauuuoqbW/8TEV5YuzYsYPm5mYuvvhitm3b5uP7MBkC0dTUxM0338xjjz3GunXrdB8/3ng8Hh+BbWxs1HaX7rvvPjIzM6mqqqKtrY3bbrst2lPMDZGIN2qdCfDjH/+YvXv38vTTT4f12D179vDVr34Vt9vNv/zLv0yLRJM/Q0ND7N69mx07dnD06FEuvPBCrr76al588UWuuuoqNm7cGBeBaG1t5cYbb+RXv/oVZ599tu7jxxtvgXjllVcoKSlhxYoVWuLyZz/7GcePH+cnP/kJ7733nk8fUoQktEgkTMVlLNx9990cOHAAIQSVlZX86le/CutxbrebL33pS7z66quUlZWxYcMGrr76apYtWxbnGU8u2dnZ3HTTTdx0001YrVZeeuklvvSlLwFj7ttut5uzzz5b1z6Z9vZ2brrpJh555JHTUiAATSBuvfVW7HY7nZ2dnHfeeVxzzTWcddZZLF26lP379wPEIhAJz7QQiSeffDKqx+3bt4+FCxcyf/58AG666SZ27do17UTCm8zMTGw2GxdccAH/9V//xWuvvcaTTz7JnXfeyUc/+lG2bdvGRz7ykZis9Do7O7nxxht5+OGHOffcc3WcffxRW7SKH/7wh6SkpPDEE09w4YUX8vrrr+NyucjPz2fevHlxOfks0ZgWy41o2b59O3v27OE3v/kNMCY2e/fu5ZFHHpnimcUXl8uFyWTyuRiUJ8b27dt55513OPvss9myZQubN2+OyFvSbDZzzTXXaCdrnW54L11hrHR86dKl3HfffSQlJXHXXXdx4YUXsnnzZu655x69dmoSerkxbbZAoyGQQE6XfpFQJCcnn/I6U1NTueyyy3j00Uc5cOAAn/zkJ3n55Zc599xz+fznP8/LL7+M3W4POW5vby/XX3893/3ud09Lgaivr+f222/H4XDgcrkANKexjo4Ovv71rzN37lzOPfdcKioqpmwrd7KZFsuNaDlT+0UmIpAnxo4dO7jvvvtYvnw5W7du5eKLL/YJtfv7+7nhhhu49957T9uKV6fTycDAACaTySc/k5OTQ3l5Oddffz2rV6+mr69Piz7PBM7o5YbL5WLx4sW8/vrrlJaWau5Ras87EiorK8nOzsZkMpGcnDztPAVgLNu/b98+tm/fzquvvsrChQu1HMYtt9zCnXfeyXXXXTfV04yJf/7nf+aWW27hggsu8MlP9Pf388wzz3Ds2DG++93v6u3endDh6xkdSSQnJ/PII49w2WWX4Xa7uf3226MSCMVf/vKXaeVI5E9SUhKbNm1i06ZNeDwePvjgA7Zv385dd93FN7/5zdNSIF577TXeffddhBBcf/312O12jh49ygUXXIAQQhOKvLw8Pve5z031dKeEMzqS0JPKykree++9aS0SwfB4PD4ly/EgXvUsb7zxBs3Nzbz88stUVFTw5JNPsmHDBu6//37WrFkDxORdGS4JHUkYIqET8+bNIz8/HyEEd9xxxxl714kHbrebxYsX+9SzPPXUU3HZqn7nnXd4/fXXsdlsXH311WzcuFH35whAQovEGb3c0JO3336bOXPm0N3dzSWXXMKSJUvYvHnzVE9rWhDvehbvSOGcc84hJyeHp59+WjuceO3atbo8z+nKGb0FqidqV6SkpIRt27ad9ia+iURbWxvl5eXa53r7n/ovJZYvX86NN95IWVkZCxcu1O15TlcMkdCBkZERhoaGtP+/8sorrFixIqzH3n777VpPgMJisXDJJZewaNEiLrnkkmnb+h4uU1HPsmLFCu666y5ycnLi+jynA9NaJDo7O3G73bq6NwWiq6uLc889l9WrV7Nx40auvPJKzTF5IibDsPV0Z6rqWeJxyvtpyQT+dqc1d955p/ziF7+ofe7xeKZwNsHx92JcvHixbG9vl1JK2d7eLhcvXjxVU0sInE6nnDdvnjx+/Li02+1y1apVsra2dqqnpSdT7mMZ6mNaRxI33XQTg4ODwFi3480338yvf/3rKZ7VxHR1dWn9A7Nnz6a7u3uKZzS1eNezLF26lBtuuCGmehaDyJjW8dTSpUt57733+OCDD/jGN77BqlWrwl4GGCQWH//4x0/bcu/TnWkZScjxHERGRgYbN27krrvuYvHixXz/+98Peliv2+2O6oi6eKAMW4GoDFsNDPRkWovEz3/+c20//Sc/+QmAJgTqZ9SuhMlkCnlK1WQSjWFroF2S++67j9LSUtasWcOaNWt46aWX4jZng2nMBEmL0xK73S5/+MMfynXr1skXXnhBbtmyRTocjoA/+8QTT8j169fLG264Qe7bt2+SZyrlTTfdJGfNmiWTk5NlaWmp/M1vfiN7enrkhRdeKBcuXCgvvPBC2dvbO+E4b775pvzHP/7hkwD9zne+I3/4wx/Gc/oG+jDlyclQH9MuJ9HT08Ndd92Fx+Phf//3fxkZGaGzs5OkpKRTTE0Btm7dyrXXXstf//pXHn30Ue3otr/+9a+cd955OJ1OUlJStJ8/fPgwOTk5lJaW6jLfp556KuDXX3/99YjGCXSsnYGBHkzUu3HaIYRIB84HaqWUbeNfewb4o5Tyz0KIJCmlR4xV43wM+AJQCLwN3AIsZWwZtg+4Tkp5WAhxNrB//HE3Mfa+PSWEMEkp3ZP+IoMghKgEXpRSrhj//D7gNmAQeA/4NynlmV2ZZRAxibEI1xEppU1K+b9KIMbpAmaMf19lJ5cAXwP2ADcAa4D3gRQp5QjwI+CrQoitwG+Bc8Yf1yqlfGp8LDeAECJJJKal1S+BBYy9tg7GXpOBQURMu+VGIKSUXxFC+AtiL5AOvCal7BVCuIF3+LDz1QlsBSzAv0op3xZCfAr4CVAkhCgFNgGvSikHJ+FlRIyUskv9XwjxKPDiFE7H4DRl2kUSwfCKIBRDQCPwrBDiD8BlwLtSyiEhxHrg00AJ8D9Syv8bf8xNwP1CiDnAV4Fbgb1CiJ8KIbIn5YVEgBBitten24Az+8Qig6iYdjmJSBFCpDEWESwA3gIuBG4GHgPmAUellE+P/2w/sAzIBH4DfFlKWTues6iRUlon/xWMIYR4irEcSxFjy6vvjH++hrHoqAm4Q0rZMSUTNDhtOeNFwp/xCz5bSvmaEOLzwGVSym1CiC3A96SUK4UQmcD/Y0wsfiGlbJzKORsYxJMzZrkRLlLKvVLK18Y/fR1QFUhfAXaP/98tpfwW0A88IoSYNbmzNDCYPIxIIkyEEOcD1UAV8A3gBcaSn9uA/5ZSGi4zBtMSQySiQAhxIXAtsBh4EPiblNIxtbMyMIgPhkgYGBiExMhJGBgYhMQQCQMDg5AYImFgYBASQyQMDAxCYoiEgYFBSAyRMDAwCIkhEgYGBiExRMLAwCAk/x/2bzM9sE1p9QAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_vals = gtsam.Values()\n",
    "for cam_idx in range(scene_data.numberCameras()):\n",
    "    plot_vals.insert(C(cam_idx),\n",
    "                        result.atPinholeCameraCal3Bundler(C(cam_idx)).pose())\n",
    "for j in range(scene_data.numberTracks()):\n",
    "    plot_vals.insert(P(j), result.atPoint3(P(j)))\n",
    "\n",
    "plot.plot_3d_points(0, plot_vals, linespec=\"g.\")\n",
    "plot.plot_trajectory(0, plot_vals, title=\"SFM results\")\n",
    "\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "341996cd3f3db7b5e0d1eaea072c5502d80452314e72e6b77c40445f6e9ba101"
  },
  "kernelspec": {
   "display_name": "Python 3.8.12 ('nbdev')",
   "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.8.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}