{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Tutorial for Strong castle model\n", "\n", "This tutorial shows how to do simulations using strong castle model and how to fit the experimental data using this model." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from cdsaxs.fitter import Fitter\n", "from cdsaxs.simulations.strong_castle import StrongCastleSimulation\n", "import numpy as np" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Simulation" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Prepare the data" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "pitch = 100 #nm distance between two trapezoidal bars\n", "qzs = np.linspace(-0.5, 0.5, 201)\n", "qxs = 2 * np.pi / pitch * np.ones_like(qzs)\n", "\n", "#Initial parameters\n", "dwx = 0.1\n", "dwz = 0.1\n", "i0 = 10\n", "bkg = 0.1\n", "y1 = 0.\n", "height = 10.\n", "bot_cd = 40.\n", "top_cd = 20.\n", "swa = [90., 90.0, 90.0]\n", "overlay = 10\n", "#fixed parameters not be fitted\n", "n1 = 2\n", "n2 = 1\n", "\n", "langle = np.deg2rad(np.asarray(swa))\n", "rangle = np.deg2rad(np.asarray(swa))\n", "\n", "\n", "overlay_params = {'heights': height,\n", " 'langles': langle,\n", " 'rangles': rangle,\n", " 'y_start': y1,\n", " 'bot_cd': bot_cd,\n", " 'dwx': dwx,\n", " 'dwz': dwz,\n", " 'i0': i0,\n", " 'bkg_cste': bkg,\n", " 'overlay': overlay,\n", " 'top_cd': top_cd,\n", " 'n1': n1,\n", " 'n2': n2,\n", " }" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "There are only 4 parameters that are different from the stacked trapezoid model, namely `overlay`, `top_cd`, `n1` and `n2.\n", "\n", "So the above parameters give us the following cross section:\n", "\n", "![overlay](../../../assets/images/overlay_ex.png)\n", "\n", "\n", "\n", "\n", "`n1` and `n2` correspond to the number of trapezoids in the bottom and top part of the model respectively. \n", "\n", " Discussion about other parameters can be found in the [stacked trapezoid model tutorial.](https://github.com/CEA-MetroCarac/cdsaxs/blob/main/Tutorials/stacked_trapezoid_simulation_and_fitting.ipynb)\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we will calculate the intensity exactly with strong castle model and plot the data." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Text(0.5, 1.0, 'Stacked Trapezoid diffraction simulation')" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAioAAAHICAYAAABgVMGnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/SrBM8AAAACXBIWXMAAA9hAAAPYQGoP6dpAABfVklEQVR4nO3dd3iUVfo+8HtKZlImM5PeSQg91FBFpCmINGFVXFmWquIqWFdd2f19EdaCuu4uu+piWVfUXRsgFkREFEGlSCeUhJJCSK+TPsnMnN8fyQwMSSBlMu+U+3Ndcylvpjzzpt055znnlQkhBIiIiIhckFzqAoiIiIhaw6BCRERELotBhYiIiFwWgwoRERG5LAYVIiIiclkMKkREROSyGFSIiIjIZTGoEBERkctiUCEiIiKXxaBCLuuHH36ATCbDxo0bu/R1EhISsGjRoi59DU+waNEiJCQkXPN+mZmZkMlkWL9+fYdeZ/369ZDJZMjMzLQdmzBhAiZMmGB3v4KCAtxxxx0ICQmBTCbD2rVrAQBnz57FzTffDJ1OB5lMhs8++6xDdThaS+/BVaxatQoymUyS1+7s10tr+H3tORhUyE5KSgruuOMOxMfHw9fXFzExMZg8eTJeeeUVu/s9//zzLvMLoKstWrQIMpnsmjf+UHSuRx99FN988w1WrFiB999/H7fccgsAYOHChUhJScFzzz2H999/H8OHD3daTadOncKqVavsQhZ1nT179mDVqlUoLy+XuhTqQkqpCyDXsWfPHkycOBHdunXDvffei8jISGRnZ2Pfvn34xz/+gQcffNB23+effx533HEHZs+eLV3BTnLfffdh0qRJtn9nZGRg5cqVWLp0KcaOHWs73qNHDynKc5q33noLFotFktfevn17s2Pff/89Zs2ahccff9x2rLa2Fnv37sWf/vQnLF++3JklAmgMKqtXr8aECROajT619B5cxf/7f/8PTz31lNRltNuePXuwevVqLFq0CHq93u5jaWlpkMv5t7gnYFAhm+eeew46nQ4HDhxo9k1fWFgoTVEuYPTo0Rg9erTt3wcPHsTKlSsxevRo/Pa3v231cdXV1QgICHBGiU7h4+Mj2WurVKpmxwoLC5t9nRYVFQFAs+Mtcfbnp6X34CqUSiWUSs/6daBWq6UugRyEcZNszp8/j/79+7f4Qz48PNz2/zKZDNXV1Xj33XebTXtkZWXhgQceQJ8+feDn54eQkBDMmTOnxaHw8vJyPProo0hISIBarUZsbCwWLFiA4uLiVms0Go2YMWMGdDod9uzZAwCwWCxYu3Yt+vfvD19fX0REROC+++5DWVmZ3WOFEHj22WcRGxsLf39/TJw4ESdPnmz/iWqBta9i165deOCBBxAeHo7Y2Nh2nRPrc+zevRv33XcfQkJCoNVqsWDBgmbvBQC+/vprjB07FgEBAQgMDMT06dPt3o+1x6el25V/7f/rX/9C//79oVarER0djWXLljUbTm+pR6W8vByLFi2CTqeDXq/HwoUL2zUMf/LkSdx4443w8/NDbGwsnn322RZHbS7v77CeJyEEXnvtNdt7WrVqFeLj4wEATzzxhN37tPZgnDp1Cr/5zW8QFBSEG264AQBw/PhxLFq0CImJifD19UVkZCSWLFmCkpKSZnXk5OTg7rvvRnR0NNRqNbp37477778f9fX1WL9+PebMmQMAmDhxoq2uH374odl7sCosLMTdd9+NiIgI+Pr6YvDgwXj33Xft7mPt4Xj55Zfx5ptvokePHlCr1RgxYgQOHDhwzXPc0NCA1atXo1evXvD19UVISAhuuOEGfPvtt7b7tNSjIpPJsHz5cmzYsAFJSUnw8/PD6NGjkZKSAgB444030LNnT/j6+mLChAnNvp5b6xFpS69OWz4nq1atwhNPPAEA6N69u+18W+to6fXT09MxZ84cBAcHw9/fH9dddx2++uoru/tYv28++eQTPPfcc4iNjYWvry9uuukmnDt37qp1U9fwrAhNnRIfH4+9e/fixIkTGDBgQKv3e//993HPPfdg5MiRWLp0KYBL0x4HDhzAnj17cNdddyE2NhaZmZlYt24dJkyYgFOnTsHf3x8AUFVVhbFjx+L06dNYsmQJhg4diuLiYnzxxRe4ePEiQkNDm71ubW0tZs2ahYMHD2LHjh0YMWIEgMapmfXr12Px4sV46KGHkJGRgVdffRVHjhzBzz//bBsJWLlyJZ599llMmzYN06ZNw+HDh3HzzTejvr7eYefwgQceQFhYGFauXInq6up2nROr5cuXQ6/XY9WqVUhLS8O6deuQlZVl+wFq/RwsXLgQU6ZMwYsvvoiamhqsW7cON9xwA44cOYKEhAT069cP77//vt1zl5eX47HHHrMLnqtWrcLq1asxadIk3H///bbXPHDggN35u5IQArNmzcJPP/2E3/3ud+jXrx82b96MhQsXtulc5efnY+LEiTCZTHjqqacQEBCAN998E35+fld93Lhx4/D+++9j/vz5mDx5MhYsWAAAGDRoEPR6PR599FHMnTsX06ZNg0ajsXvsnDlz0KtXLzz//PMQQgAAvv32W6Snp2Px4sWIjIzEyZMn8eabb+LkyZPYt2+f7Zzn5uZi5MiRKC8vx9KlS9G3b1/k5ORg48aNqKmpwbhx4/DQQw/hn//8J/74xz+iX79+AGD775Vqa2sxYcIEnDt3DsuXL0f37t2xYcMGLFq0COXl5Xj44Yft7v/BBx+gsrIS9913H2QyGV566SXcdtttSE9Pv+po16pVq7BmzRrb92xFRQUOHjyIw4cPY/LkyVc91z/++CO++OILLFu2DACwZs0azJgxA08++ST+9a9/4YEHHkBZWRleeuklLFmyBN9///1Vn6+t2vI5ue2223DmzBl8+OGH+Pvf/277mREWFtbicxYUFOD6669HTU0NHnroIYSEhODdd9/Frbfeio0bN+JXv/qV3f1feOEFyOVyPP744zAYDHjppZcwb9487N+/3yHvkdpBEDXZvn27UCgUQqFQiNGjR4snn3xSfPPNN6K+vr7ZfQMCAsTChQubHa+pqWl2bO/evQKAeO+992zHVq5cKQCITz/9tNn9LRaLEEKInTt3CgBiw4YNorKyUowfP16EhoaKI0eO2O77448/CgDif//7n91zbNu2ze54YWGhUKlUYvr06bbnF0KIP/7xjwJAi++lNQcOHBAAxDvvvGM79s477wgA4oYbbhAmk6lD58T6HMOGDbM75y+99JIAID7//HMhhBCVlZVCr9eLe++91+458/PzhU6na3bcymKxiBkzZgiNRiNOnjwphLh0Xm6++WZhNptt93311VcFAPGf//zHdmzhwoUiPj7e9u/PPvtMABAvvfSS7ZjJZBJjx45tdn5a8sgjjwgAYv/+/bZjhYWFQqfTCQAiIyPDdnz8+PFi/Pjxdo8HIJYtW2Z3LCMjQwAQf/nLX+yOP/300wKAmDt3brM6Wvr8fPjhhwKA2L17t+3YggULhFwuFwcOHGh2f+vX1IYNGwQAsXPnzmb3ufI9rF27VgAQ//3vf23H6uvrxejRo4VGoxEVFRV27ykkJESUlpba7vv5558LAOLLL79s9lqXGzx4sJg+ffpV72M9P5cDINRqtd3n4Y033hAARGRkpK0+IYRYsWJFs89ZfHx8i99XV54H6/u7/OulrZ+Tv/zlL81et7XXt369/fjjj7ZjlZWVonv37iIhIcH29W/9udOvXz9hNBpt9/3HP/4hAIiUlJRmr0Vdi1M/ZDN58mTs3bsXt956K44dO4aXXnoJU6ZMQUxMDL744os2Pcflfw03NDSgpKQEPXv2hF6vx+HDh20f27RpEwYPHtzsrxgAzYagDQYDbr75ZqSmpuKHH37AkCFDbB/bsGEDdDodJk+ejOLiYttt2LBh0Gg02LlzJwBgx44dqK+vx4MPPmj3/I888kib3ldb3XvvvVAoFHbH2npOrJYuXWr3F/L9998PpVKJrVu3Amj8a7O8vBxz5861e88KhQKjRo2yvecrPfPMM9iyZQvWr1+PpKQkAJfOyyOPPGLXeHjvvfdCq9U2Gxa/3NatW6FUKnH//ffbjikUCrum66vZunUrrrvuOowcOdJ2LCwsDPPmzWvT4zvid7/7XbNjl39+6urqUFxcjOuuuw4AbJ8fi8WCzz77DDNnzmxxFVFHlvZu3boVkZGRmDt3ru2Yj48PHnroIVRVVWHXrl129//1r3+NoKAg27+tjdzp6elXfR29Xo+TJ0/i7Nmz7a7xpptuspvuGzVqFADg9ttvR2BgYLPj16qlrdryOWmvrVu3YuTIkbYpPwDQaDRYunQpMjMzcerUKbv7L1682K6vqK3nmxzPY4LK7t27MXPmTERHR3d47wQhBF5++WX07t0barUaMTExeO655xxfrAsbMWIEPv30U5SVleGXX37BihUrUFlZiTvuuKPZN3JLamtrsXLlSsTFxUGtViM0NBRhYWEoLy+HwWCw3e/8+fNXnV663COPPIIDBw5gx44d6N+/v93Hzp49C4PBgPDwcISFhdndqqqqbE3AWVlZAIBevXrZPT4sLMzuh39nde/evdmxtp4Tqytr1Gg0iIqKss29W3/h3Hjjjc3e8/bt21tsfN62bRtWr16NFStW4Pbbb7cdt56XPn362N1fpVIhMTHR9vGWZGVlISoqqtn0ypXPdbXHX/le2/P4jmjp81NaWoqHH34YERER8PPzQ1hYmO1+1s9PUVERKioq2vw12xbW93/lyhTrVNGV575bt252/7Z+3bbUv3S5P//5zygvL0fv3r0xcOBAPPHEEzh+/HibarzyNXU6HQAgLi6uxePXqqWt2vI5aa+srKwWv7Ycfb7J8TymR6W6uhqDBw/GkiVLcNttt3XoOR5++GFs374dL7/8MgYOHIjS0lKUlpY6uFL3oFKpMGLECIwYMQK9e/fG4sWLsWHDBjz99NNXfdyDDz6Id955B4888ghGjx5t23Trrrvu6vDS1lmzZuGjjz7CCy+8gPfee8/uB7vFYkF4eDj+97//tfjY1uaru0pL/RWOPifWx7z//vuIjIxs9vErV29kZGRg3rx5mDx5Mp599tl2v54naenzc+edd2LPnj144oknMGTIEGg0GlgsFtxyyy2SLcduyZUjdVaiqdemNePGjcP58+fx+eefY/v27fj3v/+Nv//973j99ddxzz33dOg121JLa6NMZrO51cdbucLnpKPnmxzPY4LK1KlTMXXq1FY/bjQa8ac//QkffvghysvLMWDAALz44ou27vPTp09j3bp1OHHihC11t/TXlzeyDnXn5eXZjrX2Q2jjxo1YuHAh/vrXv9qO1dXVNVsJ0qNHD5w4caJNrz979mzcfPPNWLRoEQIDA7Fu3Tq759mxYwfGjBlz1SZM62qQs2fPIjEx0Xa8qKioy/9Caus5sTp79iwmTpxo+3dVVRXy8vIwbdo0AJcal8PDw+32d2lJbW0tbrvtNuj1enz44YfN/nq3npe0tDS781JfX4+MjIyrPn98fDy+++47VFVV2Y2qpKWlXbWmyx/f0nREWx/vCGVlZfjuu++wevVqrFy50nb8yrrCwsKg1Wqv+TXbnimg+Ph4HD9+HBaLxe7zkpqaavu4owQHB2Px4sVYvHgxqqqqMG7cOKxateqaQaUzgoKCWvwaz8rKsvtau1JbPydA+893S19bXXG+ybE8ZurnWpYvX469e/fio48+wvHjxzFnzhzccsstti/+L7/8EomJidiyZQu6d++OhIQE3HPPPV41orJz584W/1qw9kZcPmwaEBDQ4g8hhULR7DleeeUVmM1mu2O33347jh07hs2bNzd7jpZqWLBgAf75z3/i9ddfxx/+8Afb8TvvvBNmsxnPPPNMs8eYTCZbjZMmTYKPjw9eeeUVu+e3brveldp6TqzefPNNNDQ02P69bt06mEwmWxCfMmUKtFotnn/+ebv7WVn3EgEaezLOnDmDzZs3tzjFNWnSJKhUKvzzn/+0q/Htt9+GwWDA9OnTW31f06ZNg8lksguOZrO52S7GV3v8vn378Msvv9jV3troWFew/tV85efnyq8LuVyO2bNn48svv8TBgwebPY/18dZ9WdqyRHvatGnIz8/Hxx9/bDtmMpnwyiuvQKPRYPz48e15K626cpm1RqNBz549YTQaHfL8renRowf27dtnt6puy5YtyM7Ovurj2vo5Adp/vn/55Rfs3bvXdqy6uhpvvvkmEhISbH1b5Ho8ZkTlai5cuIB33nkHFy5cQHR0NADg8ccfx7Zt2/DOO+/g+eefR3p6OrKysrBhwwa89957MJvNePTRR3HHHXc4bMmdq3vwwQdRU1ODX/3qV+jbty/q6+uxZ88efPzxx0hISMDixYtt9x02bBh27NiBv/3tb4iOjkb37t0xatQozJgxA++//z50Oh2SkpKwd+9e7NixAyEhIXav9cQTT2Djxo2YM2cOlixZgmHDhqG0tBRffPEFXn/9dQwePLhZfcuXL0dFRQX+9Kc/QafT4Y9//CPGjx+P++67D2vWrMHRo0dx8803w8fHB2fPnsWGDRvwj3/8A3fccQfCwsLw+OOP25ZXTps2DUeOHMHXX3/d4lJoR2rrObGqr6/HTTfdhDvvvBNpaWn417/+hRtuuAG33norAECr1WLdunWYP38+hg4dirvuugthYWG4cOECvvrqK4wZMwavvvoqvvrqK7z33nu4/fbbcfz4cbu+BI1Gg9mzZyMsLAwrVqzA6tWrccstt+DWW2+1veaIESOuuqHdzJkzMWbMGDz11FPIzMxEUlISPv300zb3EDz55JO2re8ffvhh2/Jk60iDM2i1WowbNw4vvfQSGhoaEBMTg+3btyMjI6PZfZ9//nls374d48ePx9KlS9GvXz/k5eVhw4YN+Omnn6DX6zFkyBAoFAq8+OKLMBgMUKvVuPHGG+2Wg1stXboUb7zxBhYtWoRDhw4hISEBGzduxM8//4y1a9faNat2RlJSEiZMmIBhw4YhODgYBw8exMaNG7t859577rkHGzduxC233II777wT58+fx3//+99r7uDcns/JsGHDAAB/+tOfcNddd8HHxwczZ85scSO/p556Ch9++CGmTp2Khx56CMHBwXj33XeRkZGBTZs2cRdbVybJWqMuBkBs3rzZ9u8tW7YIACIgIMDuplQqxZ133imEEOLee+8VAERaWprtcYcOHRIARGpqqrPfgiS+/vprsWTJEtG3b1+h0WiESqUSPXv2FA8++KAoKCiwu29qaqoYN26c8PPzs1veW1ZWJhYvXixCQ0OFRqMRU6ZMEampqS0uVSwpKRHLly8XMTExQqVSidjYWLFw4UJRXFwshLBfnny5J598UgAQr776qu3Ym2++KYYNGyb8/PxEYGCgGDhwoHjyySdFbm6u7T5ms1msXr1aREVFCT8/PzFhwgRx4sSJVpdRtuZqy5NbWrra1nNifY5du3aJpUuXiqCgIKHRaMS8efNESUlJs+fduXOnmDJlitDpdMLX11f06NFDLFq0SBw8eNDu+Vq6Xb7MWIjG5ch9+/YVPj4+IiIiQtx///2irKzM7j5XLk8WovFzOH/+fKHVaoVOpxPz588XR44cadPyZCGEOH78uBg/frzw9fUVMTEx4plnnhFvv/12ly1PLioqalbDxYsXxa9+9Suh1+uFTqcTc+bMEbm5uQKAePrpp+3um5WVJRYsWCDCwsKEWq0WiYmJYtmyZXbLWN966y2RmJgoFAqF3VLllt5DQUGB7WtDpVKJgQMHNjtvrb0n6zm4ssYrPfvss2LkyJFCr9cLPz8/0bdvX/Hcc8/ZLYFvbXlyW89va9+rf/3rX0VMTIxQq9VizJgx4uDBg21antyez8kzzzwjYmJihFwut/u6aen7+vz58+KOO+4Qer1e+Pr6ipEjR4otW7a06b20VCc5h0wIz+sMkslk2Lx5s+06NB9//DHmzZuHkydPNmuQ0mg0iIyMxNNPP91sKL22thb+/v7Yvn37NTdGIuos66Z1Bw4ccOqF9IiIXJlXTP0kJyfDbDajsLDQ7iJylxszZgxMJhPOnz9vG5o8c+YMADZZERERScVjgkpVVZXddRgyMjJw9OhRBAcHo3fv3pg3bx4WLFiAv/71r0hOTkZRURG+++47DBo0CNOnT8ekSZMwdOhQLFmyBGvXroXFYsGyZcswefJk9O7dW8J3RkRE5L08pnvo4MGDSE5ORnJyMgDgscceQ3Jysm152zvvvIMFCxbg97//Pfr06YPZs2fjwIEDtk195HI5vvzyS4SGhmLcuHGYPn06+vXrh48++kiy90REROTtPLJHhYiIiDyDx4yoEBERkedhUCEiIiKX5dbNtBaLBbm5uQgMDOzQ1UuJiIjI+YQQqKysRHR09DU323ProJKbm9vsKp5ERETkHrKzsxEbG3vV+7h1ULFuMZ2dnQ2tVitxNURERNQWFRUViIuLa9OlItw6qFine7RaLYMKERGRm2lL2wabaYmIiMhlMagQERGRy2JQISIiIpfFoEJEREQui0GFiIiIXBaDChEREbksBhUiIiJyWQwqRERE5LIYVIiIiMhlMagQERGRy2JQISIiIpfFoEJEREQui0GFiJymymiC0WSWugwiciNuffVkInJ9ZdX1eH33efx0thin8irgo5BjaDc9JvWLwMLrE+Cj4N9LRNQ6BhUi6jKHssrw4AeHkWuosx2rN1mwL70U+9JL8fWJfLz6m2RE6fwkrJKIXJnkf8rk5OTgt7/9LUJCQuDn54eBAwfi4MGDUpdFRJ306eGL+PUbe5FrqENiaAD+cdcQ/PLHm/Dd78dj5YwkBKqVOJRVhun//AkncgxSl0tELkrSEZWysjKMGTMGEydOxNdff42wsDCcPXsWQUFBUpZFRJ10+EIZ/rDpOEwWgemDovDCbQMR6OsDAAgH0CNMg5v6heOB/x3GydwKLH3vIL588AaEaNTSFk5ELkcmhBBSvfhTTz2Fn3/+GT/++GOHHl9RUQGdTgeDwQCtVuvg6oioIwor6zDzlZ9QUGHE1AGR+Ne8oZDJZC3et6KuAbNf/RnpxdW4LjEY/717FJTsWSHyeO35/S3pT4QvvvgCw4cPx5w5cxAeHo7k5GS89dZbUpZERJ0ghMDDHx5FQYURvcI1+Mucwa2GFADQ+vrgjfnDEKBSYF96Kf727RknVktE7kDSoJKeno5169ahV69e+Oabb3D//ffjoYcewrvvvtvi/Y1GIyoqKuxuROQ6vjyeh73pJfDzUeCN+cOgUV97drlXRCD+MmcwAOCtH9ORUVzd1WUSkRuRNKhYLBYMHToUzz//PJKTk7F06VLce++9eP3111u8/5o1a6DT6Wy3uLg4J1dMRK2pazDjha2nAQAPTOiBxDBNmx87bWAUxvcOQ4NZ4LmvTndViUTkhiQNKlFRUUhKSrI71q9fP1y4cKHF+69YsQIGg8F2y87OdkaZRNQGb+1OR66hDtE6X9w7LrHdj/+/Gf2gkMuw43QBfjpb3AUVEpE7kjSojBkzBmlpaXbHzpw5g/j4+Bbvr1arodVq7W5EJL2iSiP+9cN5AMAfpvaFr4+i3c/RMzwQ869r/N5/9qtTsFgk6/MnIhciaVB59NFHsW/fPjz//PM4d+4cPvjgA7z55ptYtmyZlGURUTu9tzcTtQ1mDIrV4dbB0R1+nkcn9UagWonU/ErsOlPkwAqJyF1JGlRGjBiBzZs348MPP8SAAQPwzDPPYO3atZg3b56UZRFRO9TWm/HffVkAgN+N73HVVT7XovP3wdxR3QAAb+5Od0h9ROTeJN9Cf8aMGZgxY4bUZRBRB206fBFlNQ2IC/bDlP6RnX6+Rdcn4D8/ZWBveglSLhowMFbngCqJyF1xZyUi6jCLReA/P2UAAJaM6Q6FvOOjKVbRej/MGBQFoHG5MhF5NwYVIuqwnWmFSC+uhtZXiTuHO267gHvGNq4a+iolD7nltQ57XiJyPwwqRNRhHx1o3CLgrpHdENCGzd3aakCMDqO6B8NsEfj08EWHPS8RuR8GFSLqkJIqI3amFgIA5gyLdfjzz2kaodl46CIkvCQZEUmMQYWIOuTzo7kwWQQGx+rQKyLQ4c8/dUAk/FUKZJbU4FBWmcOfn4jcA4MKEXXIpqYpmdu7YDQFAALUSkwb2NhUu/EQp3+IvBWDChG12+m8CpzMrYCPQoaZgzq+wdu13NEUgrYcz0NtvbnLXoeIXBeDChG126amEY5J/SIQFKDqstcZmRCMuGA/VBlN+OZkfpe9DhG5LgYVImoXIQS2puQBAGYnx3Tpa8nlMvxqSONrfNX0mkTkXRhUiKhdUnIMyDXUwV+lwPjeYV3+ercMaOxT2X2mCNVGU5e/HhG5FgYVImqXbScap2Am9gnv0FWS26tfVCDiQ/xhNFmwM62wy1+PiFwLgwoRtZkQwhZUpgzo/HV92kImk+GWptf6+gT7VIi8DYMKEbXZucIqpBdXQ6WQ48a+4U573alN0z87UwtR18DVP0TehEGFiNrMOpoytlcoNA7cMv9aBsXoEKXzRU29GbvPFDntdYlIegwqRNRm2046d9rHSi6XYUr/SLsaiMg7MKgQUZvkG+pwMrcCMlnj/inOZg0qu9KKYLHw2j9E3oJBhYjaxDrlMjhWj+Au3OStNcMTgqBRK1FSXY+UHIPTX5+IpMGgQkRt8sOZxqXBE/p0/d4pLfFRyDGmZ0hjLWnsUyHyFgwqRHRNJrMFP54tBgCnbPLWmol9GlcacT8VIu/BoEJE13QkuxyVdSYE+ftgUKxesjrGN43mHLtYjtLqesnqICLnYVAhomv6oWkEY2yvMCjkMsnqiNL5oW9kIIQAfjzL6R8ib8CgQkTXtKupkVaq/pTLTbBO/6Ry+ofIGzCoENFVFVUacSKnAkDjiIrUrGFp99liLlMm8gIMKkR0VXvONzbR9o/WIixQLXE1wLD4IPirFCitrkdqfqXU5RBRF2NQIaKr2nu+BABwfY8QiStp5KOQY2T3YACXQhQReS4GFSK6qr3pjUFltIsEFeBSaNrTFKKIyHMxqBBRq3LKa5FVUgOFXIYRCcFSl2NzfY9QAMD+9BI0mC0SV0NEXYlBhYhaZZ32GRijQ6Cvj8TVXJIUpYXOzwfV9WZup0/k4RhUiKhV1qDiStM+QOPVlEcnNk3/nGOfCpEnY1AhohYJIbC3qVnVVRppL3d9T/apEHkDBhUiatGF0hrkGurgo5BheLzr9KdYWftUDmaVoa7BLHE1RNRVGFSIqEXWaZ/kuCD4qRQSV9Ncj7AAhAeqUW+y4PCFMqnLIaIuwqBCRC36JbMUAGx7lrgamUxmq+1gJoMKkadiUCGiFll/+Y9w0aACXApRB5pCFRF5HgYVImqmoKIOF0prIJcBQ7vppS6nVda9XQ5nlcHE/VSIPBKDChE1Yx1N6Rupdan9U67UJyIQWl8lquvNOJVXIXU5RNQFGFSIqBnrVMqIhCCJK7k6uVyG4U2jKr9kcPqHyBMxqBBRMwezmoKKC/enWFmnf9inQuSZGFSIyE6V0YRTuY3TKK64f8qVrKM+BzPLIISQuBoicjQGFSKyc+RCGSwCiAv2Q6TOV+pyrmlgrA4qpRwl1fU4X1QtdTlE5GAMKkRk54B1WbIbjKYAgFqpwJA4PQBO/xB5IgYVIrJzqKk/ZZiLN9Jezjr9cziLG78ReRoGFSKyMVsEjmUbAADD4t0nqAzt1hRUuJU+kcdhUCEim3OFVagymhCgUqBXeKDU5bRZclNQOV9UjfKaeomrISJHkjSorFq1CjKZzO7Wt29fKUsi8mpHmkYkBsfpoZDLJK6m7YIDVEgI8QcAHMkul7YYInIopdQF9O/fHzt27LD9W6mUvCQir3XkQjkAINmFt81vzdBuQcgsqcGRrDJM7BMudTlE5CCSpwKlUonIyEipyyAiAEeyG0dUkuPcpz/FKjk+CJ8eycHhprBFRJ5B8h6Vs2fPIjo6GomJiZg3bx4uXLggdUlEXqmirgFnC6sAAEPcckRFDwA4ml0Os4UbvxF5CkmDyqhRo7B+/Xps27YN69atQ0ZGBsaOHYvKysoW7280GlFRUWF3IyLHOJ5tgBBAt2B/hGrUUpfTbn0iAuGvUqDKaMK5psBFRO5P0qAydepUzJkzB4MGDcKUKVOwdetWlJeX45NPPmnx/mvWrIFOp7Pd4uLinFwxkeeyNtK6Y38KACgVcgyO1QPgMmUiTyL51M/l9Ho9evfujXPnzrX48RUrVsBgMNhu2dnZTq6QyHNZV8skN+3y6o6sIYsbvxF5DpcKKlVVVTh//jyioqJa/LharYZWq7W7EVHnCSEuG1Fxv0ZaK2vtxy6WS1sIETmMpEHl8ccfx65du5CZmYk9e/bgV7/6FRQKBebOnStlWURe52JZLcpqGuCjkKFvlPts9HalQbE6AI0b11UbTRJXQ0SOIGlQuXjxIubOnYs+ffrgzjvvREhICPbt24ewsDApyyLyOtYRiH5RWqiVCmmL6YQIrS8itGpYBHAqj832RJ5A0n1UPvroIylfnoiaHL/YeH0f64iEOxsUq8e3pwpwLLscIxLc4wrQRNQ6l+pRISJpHG8aURnUtGrGnQ2KaQxbKTkGiSshIkdgUCHychaLwImcxmkSjxhRaVq1ZB0lIiL3xqBC5OXSixuvmOzno0DPMI3U5XTawKYRlYziahhqGySuhog6i0GFyMsdy24ceRgQo4VS4f4/EoIDVIgN8gMAnOT0D5Hbc/+fSkTUKdZejoExemkLcSDrDrXHOP1D5PYYVIi8nHVp8uA49+9PsRoYa22oLZe2ECLqNAYVIi/WYLbgVK61kVYvbTEOZG0Ktk5rEZH7YlAh8mJnCiphNFkQ6KtEfLC/1OU4zICmhtqc8lqUVBklroaIOoNBhciLnWxaljwgWge5XCZxNY6j9fVBYlgAAO6nQuTuGFSIvNjJ3MZf4v2jPe8Cn9aN37ifCpF7Y1Ah8mInm/pTrFMlnmRgU88NgwqRe2NQIfJSFouwXbjPE0dUBnPlD5FHYFAh8lIZJdWoqTfD10eORA/YkfZKSdFayGVAQYURBRV1UpdDRB3EoELkpazTPn0jtVB4UCOtlb9Kid4RgQA4/UPkzhhUiLyUJzfSWg20NdSWS1sIEXUYgwqRlzrlwY20VtaN3ziiQuS+GFSIvJAQAidyPH9ExbrbbkqOAUIIaYshog5hUCHyQnmGOpTVNEAhl9n6ODxR36hA+ChkKK2ux8WyWqnLIaIOYFAh8kLWRtpe4Rr4+igkrqbrqJUK9IlsDGLWnhwici8MKkRe6FIjref2p1j1j2p8j9ZwRkTuhUGFyAtZf2l7cn+KVVLTezzFoELklhhUiLzQSS9opLWyvkeOqBC5JwYVIi9TVl2PXEPjTq1JXhBU+kZpIZMB+RV1KKkySl0OEbUTgwqRl7GOLCSE+CPQ10fiarqeRq1EQkgAANiubURE7oNBhcjLeFMjrVUSp3+I3BaDCpGXOdH0y9obpn2skqIYVIjcFYMKkZfxhmv8XKm/beUP91IhcjcMKkRepNpoQkZxNQDvmvqxvtf04mrU1JskroaI2oNBhciLpOZXQAggQqtGWKBa6nKcJiyw8f0KAaTmV0pdDhG1A4MKkRe5tNGb94ymWHE/FSL3xKBC5EW84YrJrbE21LJPhci9MKgQeZHTeY3THtZf2t7EOorErfSJ3AuDCpGXMFsEzhQ0BpV+XhlUGt9zan4lTGaLxNUQUVsxqBB5icySahhNFvj5KNAt2F/qcpyuW7A/NGoljCYLzhdVS10OEbURgwqRl0hrWu3SO0IDuVwmcTXOJ5fL0C8qEABwKo99KkTugkGFyEukNl3npm+k9037WFn7VE7msE+FyF0wqBB5Cev+IX0iAyWuRDq2lT+8OCGR22BQIfISaU2NtH2jvDioXLaXihBC4mqIqC0YVIi8QLXRhKySGgDePfXTOyIQPgoZDLUNyCmvlbocImoDBhUiL2BdlhwWqEZwgEriaqSjUsrRM7ypoZb7qRC5BQYVIi9g7U/p68X9KVbcSp/IvTCoEHmBNAYVG2tDLYMKkXtgUCHyAqn5XJpslWTboZZBhcgdMKgQeTghBJcmX6ZfU1i7WFaLiroGiashomtxmaDywgsvQCaT4ZFHHpG6FCKPUlhpRHlNAxRyGXqGa6QuR3I6fx9E63wBXJoSIyLX5RJB5cCBA3jjjTcwaNAgqUsh8jinmzY36x4aAF8fhcTVuIa+TX0qp7nxG5HLkzyoVFVVYd68eXjrrbcQFBQkdTlEHieN0z7NWK/5czqPIypErk7yoLJs2TJMnz4dkyZNkroUIo9kDSr9GFRsrE3FHFEhcn1KKV/8o48+wuHDh3HgwIE23d9oNMJoNNr+XVHBHzJE13LaNqLCFT9W/ZqmftLyK2GxCK+8mjSRu5BsRCU7OxsPP/ww/ve//8HX17dNj1mzZg10Op3tFhcX18VVErm3BrMF5wurAHAPlcslhPhDrZSjtsGMrNIaqcshoquQLKgcOnQIhYWFGDp0KJRKJZRKJXbt2oV//vOfUCqVMJvNzR6zYsUKGAwG2y07O1uCyoncR0ZxNerNFmjUSsTo/aQux2UoFXJbz04qp3+IXJpkUz833XQTUlJS7I4tXrwYffv2xR/+8AcoFM1XJ6jVaqjVameVSOT2rPun9I7QcHrjCn0jA3H8ogGn8ysxdWCU1OUQUSskCyqBgYEYMGCA3bGAgACEhIQ0O05EHZNm3ZE2iv0pV+rHJcpEbkHyVT9E1HVS83iNn9ZYV/5wK30i1ybpqp8r/fDDD1KXQORRbFvnRzCoXMm6l0p2aS0q6xoQ6OsjcUVE1BKOqBB5qIq6BuSU1wLgxQhbovdXIYpb6RO5PAYVIg91pumXb5TOFzp/jha0xDolxj4VItfFoELkoazTPuxPaZ2toZYjKkQui0GFyENZm0S5I23reHFCItfHoELkodI4onJNSU0Ntdat9InI9TCoEHkgIcSlqZ8oBpXWJIQEQKWUo6bejAvcSp/IJTGoEHmgXEMdKutMUMplSAzVSF2Oy1Iq5Lal29xPhcg1MagQeSDr9Wt6hmugUvLb/GourfxhQy2RK+JPMCIPZNvojf0p18SGWiLXxqBC5IHSGFTazLpDbSqXKBO5JAYVIg9k7bfox6XJ12Q9RxdKa1BZ1yBxNUR0JQYVIg9jNJmRXlQNgCMqbREUoEKktnEr/TMFHFUhcjUMKkQe5nxhNUwWAa2v0nYtG7o66xLuU2yoJXI5DCpEHiatoHHap2+kFjKZTOJq3IN1K/1UNtQSuRwGFSIPwxU/7ceLExK5LgYVIg+TmscdadsrqWlEhVvpE7keBhUiD8Nr/LRf99AAqBRyVNebkV3GrfSJXAmDCpEHKa+pR35FHQCgdwSDSlspFXL0imi81AB3qCVyLQwqRB7E2p8SG+SHQF8fiatxL/24Qy2RS+pQUFm4cCF2797t6FqIqJM47dNx1nPGixMSuZYOBRWDwYBJkyahV69eeP7555GTk+PouoioA6y/ZPtyR9p2szbUcit9ItfSoaDy2WefIScnB/fffz8+/vhjJCQkYOrUqdi4cSMaGrgFNZFUuDS546wXJ8wqqUGV0SRxNURk1eEelbCwMDz22GM4duwY9u/fj549e2L+/PmIjo7Go48+irNnzzqyTiK6BotF2KZ++nFpcrsFB6gQoVUDuDSFRkTS63QzbV5eHr799lt8++23UCgUmDZtGlJSUpCUlIS///3vjqiRiNrgYlktaurNUCnlSAgJkLoct2SdMmNDLZHr6FBQaWhowKZNmzBjxgzEx8djw4YNeOSRR5Cbm4t3330XO3bswCeffII///nPjq6XiFph7U/pGaaBUsEFfR1h20qfDbVELkPZkQdFRUXBYrFg7ty5+OWXXzBkyJBm95k4cSL0en0nyyOitrL2p3BH2o6zTplxLxUi19GhoPL3v/8dc+bMga9v61dm1ev1yMjI6HBhRNQ+XJrcedapH+tW+nI5L+pIJLUOjQ/v3LmzxdU91dXVWLJkSaeLIqL2O82lyZ2WGNa4lX6V0YSLZbVSl0NE6GBQeffdd1Fb2/ybuLa2Fu+9916niyKi9qlrMCOzuBoAR1Q6w0chR8/wpq302adC5BLaFVQqKipgMBgghEBlZSUqKipst7KyMmzduhXh4eFdVSsRteJcYRUsAgjy90FYoFrqctwat9Inci3t6lHR6/WQyWSQyWTo3bt3s4/LZDKsXr3aYcURUdtYf6n2jdRCJmNfRWdYG2pT2VBL5BLaFVR27twJIQRuvPFGbNq0CcHBwbaPqVQqxMfHIzo62uFFEtHVpXFHWofhEmUi19KuoDJ+/HgAQEZGBrp168a/3IhcRCp3pHUYa49PVmkNqo0mBKg7tDiSiBykzd+Bx48fx4ABAyCXy2EwGJCSktLqfQcNGuSQ4oiobS5d44crfjorRKNGeKAahZVGpBVUYmi3IKlLIvJqbQ4qQ4YMQX5+PsLDwzFkyBDIZDIIIZrdTyaTwWw2O7RIImpdcZURxVVGyGRA7wiN1OV4hL5RWhRWFuF0XgWDCpHE2hxUMjIyEBYWZvt/InIN1v6U+GB/+Ks4TeEI/SIDsftMERtqiVxAm3+qxcfHt/j/RCStVDbSOhyXKBO5jg5v+PbVV1/Z/v3kk09Cr9fj+uuvR1ZWlsOKI6JrS83jjrSOZr1eUmrTVvpEJJ0OBZXnn38efn5+AIC9e/fi1VdfxUsvvYTQ0FA8+uijDi2QiK4urYDX+HG0HmEabqVP5CI6NKGdnZ2Nnj17AgA+++wz3HHHHVi6dCnGjBmDCRMmOLI+IroKs0XgTAGnfhzNRyFHrwgNTuZW4FReBbqF+EtdEpHX6tCIikajQUlJCQBg+/btmDx5MgDA19e3xWsAEVHXyCqpRl2DBb4+csSHBEhdjkdhnwqRa+jQiMrkyZNxzz33IDk5GWfOnMG0adMAACdPnkRCQoIj6yOiq7Cu+OkdEQiFnBswOhKDCpFr6NCIymuvvYbRo0ejqKgImzZtQkhICADg0KFDmDt3rkMLJKLWnc5nf0pXse7yy6soE0mrQyMqer0er776arPjvCAhkXOlNf0S5Y60jpfUNKKSXVqLyroGBPr6SFwRkXfq8O5Q5eXl+OWXX1BYWAiLxWI7LpPJMH/+/DY9x7p167Bu3TpkZmYCAPr374+VK1di6tSpHS2LyKukcUSly+j9VYjS+SLPUIfU/EqMSAi+9oOIyOE6FFS+/PJLzJs3D1VVVdBq7S8r356gEhsbixdeeAG9evWCEALvvvsuZs2ahSNHjqB///4dKY3Ia9TUm5BVWgOAQaWr9IvSIs9Qh9N5FQwqRBLpUI/K73//eyxZsgRVVVUoLy9HWVmZ7VZaWtrm55k5cyamTZuGXr16oXfv3njuueeg0Wiwb9++jpRF5FXOFFRBCCBUo0aIRi11OR7J1qfChloiyXRoRCUnJwcPPfQQ/P0dt7eA2WzGhg0bUF1djdGjRzvseYk8lXVHWusvU3I868qfU7zmD5FkOhRUpkyZgoMHDyIxMbHTBaSkpGD06NGoq6uDRqPB5s2bkZSU1OJ9jUYjjEaj7d8VFfwrh7yX7Ro/EQwqXcUaVNLyK2C2CC4BJ5JAh4LK9OnT8cQTT+DUqVMYOHAgfHzsu+FvvfXWNj9Xnz59cPToURgMBmzcuBELFy7Erl27Wgwra9as4coioiapthU/DCpdJSEkAL4+ctQ1WJBZUo0eYRqpSyLyOjIhRLuvuCWXt97aIpPJYDabO1zQpEmT0KNHD7zxxhvNPtbSiEpcXBwMBgO0Wi7PJO8hhMDQZ75FWU0Dtjx4AwbE6KQuyWPNeu1nHMsux6u/ScaMQdFSl0PkESoqKqDT6dr0+7tDzbQWi6XVW2dCivW5Lw8jl1Or1dBqtXY3Im9UVGlEWU0D5DKgZzj/yu9KSWyoJZJUh/dRsaqrq4Ovr2+HHrtixQpMnToV3bp1Q2VlJT744AP88MMP+OabbzpbFpFHs+5ImxAaAF8fhcTVeLZLW+mzoZZICh0aUTGbzXjmmWcQExMDjUaD9PR0AMD//d//4e23327z8xQWFmLBggXo06cPbrrpJhw4cADffPON7SKHRNSySyt+OKrY1XjNHyJpdSioPPfcc1i/fj1eeuklqFQq2/EBAwbg3//+d5uf5+2330ZmZiaMRiMKCwuxY8cOhhSiNrCu+OnHRtouZ91ML89Qh/KaeomrIfI+HQoq7733Ht58803MmzcPCsWlYefBgwcjNTXVYcURUcusf9335TV+ulygrw/igv0AAKc4qkLkdB0KKjk5OejZs2ez4xaLBQ0NDZ0uiohaV2+y4HxRFQCgLzd7c4p+kexTIZJKh4JKUlISfvzxx2bHN27ciOTk5E4XRUStSy+uQoNZIFCtRIzeT+pyvAL7VIik06FVPytXrsTChQuRk5MDi8WCTz/9FGlpaXjvvfewZcsWR9dIRJdJbfqrvm9UoN0FQanrMKgQSadDIyqzZs3Cl19+iR07diAgIAArV67E6dOn8eWXX7IZlqiLneaKH6dLajrXZwuq0GC2SFwNkXfp8D4qY8eOxbfffuvIWoioDax7qLCR1nlig/ygUStRZTQhvaialy0gcqIOjagkJiaipKSk2fHy8nKHXKiQiFpn3UOFjbTOI5fLbFepPpVnkLgaIu/SoaCSmZnZ4lb5RqMROTk5nS6KiFpWUmVEYWXjJSZ41WTn6h/deD2lkznsUyFypnZN/XzxxRe2///mm2+g0126EJrZbMZ3332HhIQEhxVHRPbSmqZ94kP8EaDu9BUwqB2Sohun2k7kckSFyJna9ZNu9uzZABqvkLxw4UK7j/n4+CAhIQF//etfHVYcEdm71J/C0RRn698UVE7lVkAIwRVXRE7SrqBisTR2u3fv3h0HDhxAaGholxRFRC1L5Y60kukVHggfhQwVdSZcLKtFXLC/1CUReYUO9ahkZGQwpBBJ4HS+dWkyR1ScTaWUo3dTX9BJTv8QOU2HJ7m/++47fPfddygsLLSNtFj95z//6XRhRGTPZLbgTEHT1vkcUZFE/2gtTuZW4GRuBW4ZECV1OUReoUNBZfXq1fjzn/+M4cOHIyoqinO1RE6QWVKNepMF/ioFunHaQRKNK38u4mQuV/4QOUuHgsrrr7+O9evXY/78+Y6uh4haYb0gXp/IQMjl/ONACtaGWk79EDlPh3pU6uvrcf311zu6FiK6itR8NtJKrV+UFjIZUFBhRHGVUepyiLxCh4LKPffcgw8++MDRtRDRVVgvRshGWukEqJXoHhIAAJz+IXKSDk391NXV4c0338SOHTswaNAg+Pj42H38b3/7m0OKI6JLUnmNH5eQFK1FenE1TuQYML53mNTlEHm8DgWV48ePY8iQIQCAEydOOLIeImqBobYBOeW1AMAL4kmsf7QOW47n4RRHVIicokNBZefOnY6ug4iuwrrRW4zeDzo/n2vcm7oSG2qJnKtdQeW222675n1kMhk2bdrU4YKIqLlUbp3vMqxBJbOkBpV1DQj0ZXAk6krtCiqXX4SQiJzHtuKHjbSSC9GoEaXzRZ6hDqfzKjGye7DUJRF5tHYFlXfeeaer6iCiqziVx0ZaV9I/Wos8Qx1O5hoYVIi6WIeWJxOR85jMFluPSlI0g4orSIpuHF3mEmWirsegQuTiMoqrYWzaOj+haQ8PktalhloGFaKuxqBC5OJO5Vl3pA2EglvnuwRrUDlbUAmjySxxNUSejUGFyMVZ/2rvH81mdldhXSZusgicya+Suhwij8agQuTirBuLsT/FdchkMu6nQuQkDCpELkwIYZv66c+g4lLYp0LkHAwqRC4sv6IOpdX1UMhl6B3BPVRciXUq7gRHVIi6FIMKkQuzTvv0DNPA10chcTV0uYGxjUHlVG4FGswWiash8lwMKkQu7CT7U1xW95AABKqVMJosOFvAhlqirsKgQuTCTuWyP8VVyeUyDIhpHFU5frFc2mKIPBiDCpELO5nX2P+QFMWg4ooGxTUGlWMX2adC1FUYVIhclKG2AdmltQA49eOqBsfqAQApOeWS1kHkyRhUiFzU6aZlyTF6P+j9VRJXQy0Z2DT1k5pXiboG7lBL1BUYVIhcFDd6c32xQX4IDlDBZBFIza+Uuhwij8SgQuSibCt+2J/ismQyGQbFsqGWqCsxqBC5KO5I6x4GNU3/HMtmQy1RV2BQIXJBRpMZZwsapxI49ePaBrGhlqhLMagQuaCzBVUwWQR0fj6I0ftJXQ5dhXXq51xhFaqNJomrIfI8DCpELujUZf0pMplM4mroasK1vojU+sIigBM5nP4hcjQGFSIXxP4U92IdVUlhUCFyOAYVIhd0sumKvOxPcQ+D4/QAuEMtUVdgUCFyMRaLwOm8xkba/tE6iauhtrBu/JbCJcpEDidpUFmzZg1GjBiBwMBAhIeHY/bs2UhLS5OyJCLJXSitQZXRBJVSjsSwAKnLoTawTv1kltTAUNMgcTVEnkXSoLJr1y4sW7YM+/btw7fffouGhgbcfPPNqK6ulrIsIkmdaJr26RsZCB8FBz3dgd5fhfgQfwDAcS5TJnIopZQvvm3bNrt/r1+/HuHh4Th06BDGjRsnUVVE0kpp6nOwTieQexgUq0dWSQ2OXzRgbK8wqcsh8hgu9eeawdD4Azo4OLjFjxuNRlRUVNjdiDzN8aagYr0yL7kH6w613EqfyLFcJqhYLBY88sgjGDNmDAYMGNDifdasWQOdTme7xcXFOblKoq5lsQjbXhwDYzmi4k4uXfOHK3+IHMllgsqyZctw4sQJfPTRR63eZ8WKFTAYDLZbdna2Eysk6noZJdWoNJrg6yNHr3CN1OVQOwyI0UEmA/IMdSisrJO6HCKP4RJBZfny5diyZQt27tyJ2NjYVu+nVquh1WrtbkSexNqf0j9aByUbad1KgFqJnmGN4fI4L1BI5DCS/iQUQmD58uXYvHkzvv/+e3Tv3l3Kcogkd6ypv4GNtO4puZseAHAku0zaQog8iKRBZdmyZfjvf/+LDz74AIGBgcjPz0d+fj5qa2ulLItIMtYRlUHsT3FLyd2CAABHLpRLWwiRB5E0qKxbtw4GgwETJkxAVFSU7fbxxx9LWRaRJExmi20PlUFc8eOWhjYFlWPZ5TBbhMTVEHkGSfdREYLfyERW54qqUNdgQYBKgcRQ7kjrjnqGa6BRK1FlNOFMQSX6RbGPjqiz2K1H5CKsy1oHxOggl8skroY6QiGXYXBc47Td4QvsUyFyBAYVIhdh3SjMeiVeck/JcexTIXIkBhUiF8Gt8z3D0Hg9AI6oEDkKgwqRC6g3WXA6rxIAV/y4uyFNIyrpRdUor6mXuBoi98egQuQC0vIrUW+2QOfng27B/lKXQ50QHKBCQtOVlI9kl0tbDJEHYFAhcgHHc8oBNI6myGRspHV31mXKR7I4/UPUWQwqRC6AG715lqHxjUHlEPtUiDqNQYXIBRyzNdLqpS2EHGJ4wqWVPw1mi8TVELk3BhUiidU1mHGmgI20nqR3eCACfZWoqTfjdF6F1OUQuTUGFSKJncytgNkiEKpRI0rnK3U55AByuQzDm6Z/DmZy+oeoMxhUiCSW0rTRGxtpPcvwhGAAwMGsUokrIXJvDCpEEjuew43ePNHlIyq8rhlRxzGoEEnsWNNeG9ZrxJBnGBynh49ChsJKI7JLa6Uuh8htMagQSchQ04DzRdUALu1oSp7B10eBAU2jZAcyOf1D1FEMKkQSOpLd2GiZEOKP4ACVxNWQo42w9amwoZaooxhUiCRkvcKudSdT8izDmvpUOKJC1HEMKkQSsl5hN7mbXtpCqEuMbBpROVdYheIqo8TVELknBhUiiVgsAkebGmmTOaLikYICVOgbGQgA2J/OURWijmBQIZJIenEVKutM8PWR236Zkee5LjEEALA/o0TiSojcE4MKkUQOZ5UDAAbF6qFU8FvRU43q3jj9wxEVoo7hT0ciiVhX/LA/xbONbAoqaQWVKK2ul7gaIvfDoEIkEa748Q4hGjV6R2gAAL9w+oeo3RhUiCRQWdeAtKYrJifH6aUthrqctU9lH6d/iNqNQYVIAocvlEMIoFuwP8K1vGKypxvV3RpUOKJC1F4MKkQSONi0AdjwBE77eINRiY19Kqn57FMhai8GFSIJWHcqtW6xTp4tVKO2LUHfe56jKkTtwaBC5GQNZotto7fh8RxR8RZjeoYCAH46VyxxJUTuhUGFyMlO5lagrsECvb8PeoRppC6HnOSGpqDyM4MKUbswqBA5ma0/JT4IcrlM4mrIWUZ2D4ZSLsOF0hpcKKmRuhwit8GgQuRkB2yNtOxP8SYBaqVtz5yfz3NUhaitGFSInEgIgYOZjTvSjuCKH6/DPhWi9mNQIXKijOJqlFTXQ6WUY0CMTupyyMnG9GzcT2XPuWJYLELiaojcA4MKkRNZp30Gx+qgViokroacbXCcHgEqBcpqGnAqr0LqcojcAoMKkRNZt1C3bqlO3sVHIcfoHo2f+11niiSuhsg9MKgQOYkQwraFOoOK95rQJxwA8ENaocSVELkHBhUiJ7lQWoM8Qx18FDJeMdmLTegTBqDxek+GmgaJqyFyfQwqRE5iHU0ZHKuHn4r9Kd4qNsgfPcM1MFsEfjzH6R+ia2FQIXIS9qeQ1cSmUZUf0hhUiK6FQYXICYQQ2M/+FGpyqU+liMuUia6BQYXICbJLa5Fr7U+J10tdDklseEIQAlQKFFcZuUyZ6BoYVIicwNqfMihWD3+VUuJqSGpqpQLXN+1S+30qV/8QXQ2DCpET7Gm6tst1iby+DzW6qW/j9M+O0wUSV0Lk2hhUiLqYEAI/nWscUbFe64Xopn4RkMmA4xcNyC2vlbocIpfFoELUxdIKKlFcZYSfjwLD4rl/CjUKC1RjWNN+OhxVIWqdpEFl9+7dmDlzJqKjoyGTyfDZZ59JWQ5Rl/jpbOO0z8juwby+D9m5uX8EAGD7SQYVotZIGlSqq6sxePBgvPbaa1KWQdSlfjrXGFTG9uK0D9mbnBQJoLHZmrvUErVM0uUHU6dOxdSpU6UsgahLGU1m7G/a6I39KXSl7qEB6B2hwZmCKuxMK8Ts5BipSyJyOW7Vo2I0GlFRUWF3I3Jlh7PKUdtgRqhGjb6RgVKXQy7o5qZRlW9O5ktcCZFrcqugsmbNGuh0OtstLi5O6pKIrurnpmmfG3qGQCaTSVwNuaIp/RuDyg9pRaipN0lcDZHrcaugsmLFChgMBtstOztb6pKIrurHs43XcuG0D7VmQIwW8SH+qG0wY8dpbv5GdCW3CipqtRpardbuRuSqiquMOJ5jAACM7x0mcTXkqmQyGWYOigYAfHksV+JqiFyPWwUVIneyK60IQjT+xRyu9ZW6HHJhMwc3BpVdaUUw1HL1D9HlJA0qVVVVOHr0KI4ePQoAyMjIwNGjR3HhwgUpyyJyiO/TGofxJzZdKZeoNX0iA9E7QoN6swXb2VRLZEfSoHLw4EEkJycjOTkZAPDYY48hOTkZK1eulLIsok4zmS3YfaaxP2ViXwYVujbb9M/xPIkrIXItku6jMmHCBAghpCyBqEscyipDZZ0JwQEqDI7VS10OuYEZg6Px12/P4OdzxSiqNCIsUC11SUQugT0qRF1gZ1rjaMr43mFQyLksma6te2gAhsTpYbYIfHYkR+pyiFwGgwpRF9iZ2tifMqEPV/tQ2905vHFvqE8OZnO0magJgwqRg2WX1iCtoBJyGZclU/vMGBwFXx85zhZW4dhFg9TlELkEBhUiB7NuhT6yezD0/iqJqyF3ovX1wS1NO9VuOMgNLYkABhUih7MGFesvHKL2sE7/fHE0F7X1ZomrIZIegwqRAxVVGnEwqwwAcDODCnXAdYkhiA3yQ6XRhK9SuFSZiEGFyIG+PVUAIYBBsTpE6/2kLofckFwuw9yR3QAA7+7JZFMteT0GFSIHsk77TOFoCnXC3JHdoFLKkZJjwOEL5VKXQyQpBhUiB6moa8Ce88UAGFSoc4IDVLi16fo/7+7JlLYYIokxqBA5yI5TBWgwC/QIC0DPcI3U5ZCbW3R9AgBga0oeCirqpC2GSEIMKkQO8sWxXADAjKZrthB1xoAYHYbHB8FkEfjvviypyyGSDIMKkQOUVBnx49nGaZ9bhzCokGPcfUN3AI3TPxV1DRJXQyQNBhUiB9iakgezRWBgjA49wjjtQ44xpX8keoZrUFFnwvt7OapC3olBhcgBPj/aOO0zi6Mp5EByuQzLJvYAALz9UwZq6k0SV0TkfAwqRJ10sawGB7PKIJOxP4Ucb+agaHQL9kdpdT0+2H9B6nKInI5BhaiTrE2013UPQaTOV+JqyNMoFXI8MKFxVOX1XedRZeSoCnkXBhWiThBCYMPBiwCA2ckcTaGucdvQWCSE+KO4qh5v7DovdTlETsWgQtQJ+zNKkVFcjQCVgtM+1GVUSjmemtoXAPDWj+nIM9RKXBGR8zCoEHXCxweyATQuSQ5QKyWuhjzZlP6RGB4fhLoGC/66/YzU5RA5DYMKUQcZahqwtenqtr8e0U3iasjTyWQy/Gl6PwDApsMXcSirVOKKiJyDQYWogz47mgOjyYK+kYEYHKuTuhzyAsndgnDHsFgIATy58TjqGsxSl0TU5RhUiDpACGFbKnrXiDjIZDKJKyJv8f+m90OoRo3zRdV4bec5qcsh6nIMKkQd8PO5EqQVVMJfpcCvkmOlLoe8iN5fhWdm9QcArPvhPFIuGiSuiKhrMagQdcDbP6UDAOYMi4XO30fiasjbTB0YhakDImGyCDzwwSEYangdIPJcDCpE7XSusAo704ogkwGLx3SXuhzyUi/cNgixQX7ILq3F7zccg8UipC6JqEswqBC10zs/ZwAAbuobgYTQAImrIW+l8/fBunnDoFLIseN0Af75/VmpSyLqEgwqRO1QXGXEpsONO9HefQNHU0haA2N1+HNTv8raHWfx/j5eYZk8D4MKUTu8uTsddQ0WDIrV4brEYKnLIcJdI7vhoZt6AQBWfn4Cnx3JkbgiIsdiUCFqo6JKI97bmwkAeGRSLy5JJpfx6KReWDg6HkIAj35yFP/+MR1CsGeFPAODClEbvbHrPOoaLBgcp8fEPuFSl0NkI5PJ8PTM/ljQFFae/eo0Vn5+khvCkUdgUCFqg8LKOvx3f+P8/6McTSEXJJfLsPrW/vjTtMZt9t/fl4VbX/0JJ3K4zwq5N15FjagNXv4mDXUNFiR302N87zCpyyFqkUwmw73jEpEYFoA/bDqOMwVVmP3az5gzPBbLb+yFGL1fl75+ldGEsup6GGobUF7TgPLaepgtAgq5DEq5DEq5HD5KOUICVAgPVCNUo4ZcztBPVycTbjyRWVFRAZ1OB4PBAK1WK3U55KGOZpdj9ms/AwA23X89hsUHSVwR0bWVVBnxf5+fwNaUfACASiHHlAGR+FVyNMb2CoOPov0D6haLQGGlEVkl1cgqqUFWadN/S2qQVVKNijpTu55PrZQjPsQffSK1SI7TI7mbHknRWqiVinbXRu6lPb+/GVSIrsJiEbht3R4czS7HbUNj8Lc7h0hdElG7HMgsxd+2n8He9BLbMT8fBYbE6TEoVofYID9EaH3hr1JCpZSj3mRBTb0JJdX1yDfUNd4q6pBbXosLpTUwmixXfT1fHzn0firo/X2g9fOBSiGHyWKBySxgsgjUNZhRUl2PkiojWtqjTqWQY0CMFhP7hGNy/wj0iQjkVKsHYlAhcpBPDmbjyY3HEaBSYOfjExCu9ZW6JKIOOX6xHJ8ezsGW47korqrv8PMo5DLEBvmhW7A/EkICEB/i3/j/oQGIDfKDv6ptHQUmswW55XVIL67CiRwDjlwox5HscpRW29cWF+yHyf0iMW1gJIbFBzG0eAgGFSIHyC2vxZS1u1FZZ8KKqX1x3/geUpdE1GkWi8D5oiocyCzD2cJK5JTVoqDSCGODGfUmC1RKOdQ+CgT7+yBS54conS8itb6I1PkiPsQf0Xq/Dk0btYUQAhdKa7D3fAm+PVWAn84V243gdA8NwB3DYnH70FhE6vhHgztjUCHqJItFYP5/9uPncyUYEqfHxt+NhrKLfjgTUctq6k3YfaYY20/mY9vJfNTUNy63lsuAcb3DMP+6eEzsE86G3A6oqTehpKoepdX1KK2ph7HBAkBACMAiGs+xv1oJjVqB8EBfxAX7O/T1GVSoy9Q1mFFvtsBsFjALAbNFwEchh87PBwoP+mHxzs8ZWP3lKfj6yLH1obFIDNNIXRKRV6s2mvBVSh42HMzGgcwy2/H4EH8sGJ2AOcNjofXllcwv12C2IC2/EmcLK3G+sBrni6qQXlSNC6U1qG3HHjvTB0bhtXlDHVobgwp1iBACBRVGpBdV4XxxNTKKqpFTXoPiqnoUVxlRVGm0/UVzJZkM0Pr6IMjfB8EBKiSEBKBXRCB6R2jQKzwQsUF+bvNXz/70Evz27f1oMAs8M6s/5o9OkLokIrpMRnE1PtifhY8PZNtWGvmrFLhjWCwWjE5Az3Dv/MOitLoev2SU4kh2GY5kleN4TjnqGlpvflYp5AgOUCE4QAVfHzlkMhnkMkAGGSxCoLrejCpjAyb3i8TKmUkOrZVBha7JYhHIKq3B8YvlSLlowPEcA07lVqDK2L7lhXIZWuzcv1KgrxKjugfjusQQXJcYgqQorUsGl+zSGtz66k8oq2nAjEFReGVuMpv3iFxUTb0Jm4/kYP3PmThbWGU7Pq53GBZdH48JvT17WshoMuNQZhl+PFeMn84W40SuAVf+Rg/0VaJflBY9wzVIDA1Aj3ANuocEIDRQjQCVQrKfbwwqZEcIgezSWhzPaQolFw04kWtAZQt7HijkMnQL9kdiaAASwwIQF+yPMI0aoYFqhGnUCNaooFbKoZDJoJDLIJPJ0GC2NG7uVFOPspoGFFc1jsqcKajCmYJKpBdVo95sn+r1/j6Y0DsMtwyIxLjeYW1eKdCVSqqMuOvNfThbWIVBsTp8ct9o+PpwPwciVyeEwJ7zJVi/JxM7ThfYflknhPhjvodNC5VUGfF9aiF2nC7Aj2eLm41y947QYFh8MJK76TG0mx6JoRqXDGsMKl7M2jV/IqcCJ3INOJHTGEwMtQ3N7qtWypEUrcWgGB0GxuoxMEaH7qEBUCkd2zRqMltwKq8Ce8+XYG96CQ5klKL6sm8uXx85xvUKw9SBkZicFAmN2vmhpbCyDvPe2o+zhVWI0Krx+bIbuKqAyA1dKKnB+/sy8dGBbNsfY/4qBW4fGouF18ejZ3igxBW23/miKuw4VYAdpwtwKKvMbhQ7LFCNsT1DcUOvUNzQM9RttlBgUPESNfUmpBdV40xBJU7mVuBEjgGn8ipaHClRKeToFxWIgbE6DIrRY2CsDj3DNV22zPBqGswWHM0ut3XyZ5fW2j6mVspxU79wzBwUjYl9w50yonGhpAaL1v+C9KJqRGjV+ODe69CDzbNEbq21aaHre4RgzvBY3NI/Cn4q1xwxNVsEDl8ow45TBfj2dAHSi6rtPt4/WotJ/SIwOSkC/aO1bjk9zaDSSd+nFuAPm1IQ6KtEoK8PtL7Kxv9X+0Dr13jM+rHApo9przjW2QAgmhqZyqrrkWdo3BUy11CLvPI6ZJZU43xhFXINdS0+VqWUo19kIPrH6NA/WovBsXr0jgh0+EiJIwghcCqvAttO5OOr43lIL770DalRK3Fz/wjMHByNG3qGdkmo2nYiH09sPIbKOhNi9H744N5RiA8JcPjrEJE0hBDYe74E71wxLaRRKzFjUBTmDI/F0G7SbyRXZTThxzNF+PZ0AXamFqKs5tIouI9ChusSQzA5KQI39Yvo8ms2OYPbBZXXXnsNf/nLX5Cfn4/BgwfjlVdewciRI6/5uK4KKpsOXcTvNxzr1HP4+shtoUWtVMBH0XRRLoUcKoUcyqZ/N5gF6k0WNJgtqDdbUG002S7oZWpDl2pIgAo9wjRIitaif7QWA2KkGynpLCEETuZW4MtjufjyWK5dEAsOUGHqgEjMHByNkQnBnZ5zzTfU4S/fpGHT4YsAgGHxQXj1N8mI0rn/DwAiatnFshpsPHQRGw9dxMWySyO5kVpfTEoKx+SkSFyXGOyUaw2ZLQIpOQb8fK4YP58rxsHMMrtePp2fDyb2CcOkpAiM6x3mMT02Vm4VVD7++GMsWLAAr7/+OkaNGoW1a9diw4YNSEtLQ3h4+FUf21VBxVDbgOzSGlTWmVBZ12D/X2Pj/1fUmVBRe+ljVUYTKutMrS7f7Si1Uo5InS+idX6I0vsiRu+H2CC/pg5uDYICVA59PVdhaRr6/OJYLram5Nlt+R2p9cW0gVGY2DcMIxKC2zU9lF5UhY8OZOP9vVm2fQTuHdsdT97S1y3DHRG1n8UisD+jFBsPXcTXJ/Lsfm5r1Epc3yMEQ+ODkBynx6BYvUOmiIqrjEjJMeDERQOOXSzH/ozSZtP08SH+mNwvApOSIjA8PsijN5l0q6AyatQojBgxAq+++ioAwGKxIC4uDg8++CCeeuqpqz7WFXtUTGaLLbRUNIWbepMFJosF9SZhuzhXvdli2yxNpWwcZVEpZfBXKaH397Fd1IurThrP6d70Enx5LBdfn8i3++ZWKeUYGKPDwBgd+kQ27tcSFqiGStG4J0BptRF5hjocvVCOA5mlOHbRYHvssPgg/N+MJAyJ00vwrojIFdQ1mLH3fAm2NzWrFlUa7T6ulMvQNyoQvcID0T00AFE6X4RrfRHk79N4IUeFHPVmM+oaGkfF6xrMKKww4mJZDXLKa3GxrLbVqXqtrxLXJYZgTM9QjOkZih5hAZJPQTmL2wSV+vp6+Pv7Y+PGjZg9e7bt+MKFC1FeXo7PP//c7v5GoxFG46UvooqKCsTFxblUUKGuZTSZsSutCNtPFeCns8XIr2i5T6c1chkwoU845o7shkn9wr3mhwIRXZvFInA8x4BfMkpwOKschy+UofCK4NJRMhmQGBrQ+IdVrB7D44MwIEbnUTt6t0d7goqkm1cUFxfDbDYjIiLC7nhERARSU1Ob3X/NmjVYvXq1s8ojF6RWKnBz/0jc3D8SQghkFFfjeNNQamZxNXLKa1FSVY96kwVmIRCiUSFMo0a/KC2GJwTh+h6hiHCT5XtE5FxyuQxD4vS2UVYhBHINdUi5WI70pt268yvqUFRphKG2ATX1jRdyVPvIoVY2jo6rlQqEalSIDfK3TdV3C/ZHUrQWgR7WZ+Is0u+y1Q4rVqzAY489Zvu3dUSFvJNMJkNimAaJYRrMTo6Ruhwi8jAymQwxej+PWGXjziQNKqGhoVAoFCgoKLA7XlBQgMjIyGb3V6vVUKvVziqPiIiIJCZpS7FKpcKwYcPw3Xff2Y5ZLBZ89913GD16tISVERERkSuQfOrnsccew8KFCzF8+HCMHDkSa9euRXV1NRYvXix1aURERCQxyYPKr3/9axQVFWHlypXIz8/HkCFDsG3btmYNtkREROR9JN9HpTNccR8VIiIiurr2/P723G3viIiIyO0xqBAREZHLYlAhIiIil8WgQkRERC6LQYWIiIhcFoMKERERuSwGFSIiInJZDCpERETkshhUiIiIyGVJvoV+Z1g31a2oqJC4EiIiImor6+/ttmyO79ZBpbKyEgAQFxcncSVERETUXpWVldDpdFe9j1tf68disSA3NxeBgYGQyWRSlyO5iooKxMXFITs7m9c+6kI8z87B8+wcPM/Ow3N9iRAClZWViI6Ohlx+9S4Utx5RkcvliI2NlboMl6PVar3+m8AZeJ6dg+fZOXienYfnutG1RlKs2ExLRERELotBhYiIiFwWg4oHUavVePrpp6FWq6UuxaPxPDsHz7Nz8Dw7D891x7h1My0RERF5No6oEBERkctiUCEiIiKXxaBCRERELotBhYiIiFwWg4qbKy0txbx586DVaqHX63H33XejqqqqTY8VQmDq1KmQyWT47LPPurZQN9fe81xaWooHH3wQffr0gZ+fH7p164aHHnoIBoPBiVW7vtdeew0JCQnw9fXFqFGj8Msvv1z1/hs2bEDfvn3h6+uLgQMHYuvWrU6q1L215zy/9dZbGDt2LIKCghAUFIRJkyZd8/NCjdr79Wz10UcfQSaTYfbs2V1boJtiUHFz8+bNw8mTJ/Htt99iy5Yt2L17N5YuXdqmx65du5aXHmij9p7n3Nxc5Obm4uWXX8aJEyewfv16bNu2DXfffbcTq3ZtH3/8MR577DE8/fTTOHz4MAYPHowpU6agsLCwxfvv2bMHc+fOxd13340jR45g9uzZmD17Nk6cOOHkyt1Le8/zDz/8gLlz52Lnzp3Yu3cv4uLicPPNNyMnJ8fJlbuX9p5nq8zMTDz++OMYO3askyp1Q4Lc1qlTpwQAceDAAduxr7/+WshkMpGTk3PVxx45ckTExMSIvLw8AUBs3ry5i6t1X505z5f75JNPhEqlEg0NDV1RptsZOXKkWLZsme3fZrNZREdHizVr1rR4/zvvvFNMnz7d7tioUaPEfffd16V1urv2nucrmUwmERgYKN59992uKtEjdOQ8m0wmcf3114t///vfYuHChWLWrFlOqNT9cETFje3duxd6vR7Dhw+3HZs0aRLkcjn279/f6uNqamrwm9/8Bq+99hoiIyOdUapb6+h5vpLBYIBWq4VS6daX2HKI+vp6HDp0CJMmTbIdk8vlmDRpEvbu3dviY/bu3Wt3fwCYMmVKq/enjp3nK9XU1KChoQHBwcFdVabb6+h5/vOf/4zw8HCOtF4Df2K6sfz8fISHh9sdUyqVCA4ORn5+fquPe/TRR3H99ddj1qxZXV2iR+joeb5ccXExnnnmmTZPy3m64uJimM1mRERE2B2PiIhAampqi4/Jz89v8f5t/Rx4o46c5yv94Q9/QHR0dLOQSJd05Dz/9NNPePvtt3H06FEnVOjeOKLigp566inIZLKr3tr6Q+ZKX3zxBb7//nusXbvWsUW7oa48z5erqKjA9OnTkZSUhFWrVnW+cCIneeGFF/DRRx9h8+bN8PX1lbocj1FZWYn58+fjrbfeQmhoqNTluDyOqLig3//+91i0aNFV75OYmIjIyMhmjVomkwmlpaWtTul8//33OH/+PPR6vd3x22+/HWPHjsUPP/zQicrdS1eeZ6vKykrccsstCAwMxObNm+Hj49PZsj1CaGgoFAoFCgoK7I4XFBS0ek4jIyPbdX/q2Hm2evnll/HCCy9gx44dGDRoUFeW6fbae57Pnz+PzMxMzJw503bMYrEAaBytTUtLQ48ePbq2aHcidZMMdZy1yfPgwYO2Y998881Vmzzz8vJESkqK3Q2A+Mc//iHS09OdVbpb6ch5FkIIg8EgrrvuOjF+/HhRXV3tjFLdysiRI8Xy5ctt/zabzSImJuaqzbQzZsywOzZ69Gg2015De8+zEEK8+OKLQqvVir179zqjRI/QnvNcW1vb7OfwrFmzxI033ihSUlKE0Wh0Zukuj0HFzd1yyy0iOTlZ7N+/X/z000+iV69eYu7cubaPX7x4UfTp00fs37+/1ecAV/1cU3vPs8FgEKNGjRIDBw4U586dE3l5ebabyWSS6m24lI8++kio1Wqxfv16cerUKbF06VKh1+tFfn6+EEKI+fPni6eeesp2/59//lkolUrx8ssvi9OnT4unn35a+Pj4iJSUFKneglto73l+4YUXhEqlEhs3brT7uq2srJTqLbiF9p7nK3HVT+sYVNxcSUmJmDt3rtBoNEKr1YrFixfb/UDJyMgQAMTOnTtbfQ4GlWtr73neuXOnANDiLSMjQ5o34YJeeeUV0a1bN6FSqcTIkSPFvn37bB8bP368WLhwod39P/nkE9G7d2+hUqlE//79xVdffeXkit1Te85zfHx8i1+3Tz/9tPMLdzPt/Xq+HINK62RCCOHs6SYiIiKituCqHyIiInJZDCpERETkshhUiIiIyGUxqBAREZHLYlAhIiIil8WgQkRERC6LQYWIiIhcFoMKERERuSwGFSIiInJZDCpE5BJSU1MxceJE+Pr6onfv3ti6dStkMhmOHj0qdWlEJCEGFSKSXGpqKkaNGoWxY8fi5MmTePHFF7FgwQL4+PggKSlJ6vKISEK81g8RSe6mm25CfHw8/vOf/9iO3X777Th79iyOHz8uYWVEJDWl1AUQkXfLysrC999/j2PHjtkdV6lUGDx4sERVEZGr4NQPEUnq6NGjLU7xnDhxwhZUxo8fjyFDhmDIkCFQKBQ4ePCgFKUSkQQ4okJEkpLL5TCbzTCbzVAqG38kbdu2zS6o7Nq1CwDw9NNPY/z48Rg+fLhk9RKRc3FEhYgkNWzYMPj4+OCPf/wj0tPTsWnTJixbtgwA7KZ+1q5di8zMTKxdu1aiSolICgwqRCSp6Oho/Pvf/8Ynn3yCwYMH4+OPP8a9996LyMhIhIeHAwDWr1+P3bt34z//+Q9kMpnEFRORM3HVDxG5nMceewynTp3Ctm3bsHnzZrzxxhv4/PPPoVarpS6NiJyMIypE5HKOHz9um/ZZsmQJ0tPTMWrUKAwZMgRbtmyRuDoiciY20xKRy0lJScHixYsBAGVlZRJXQ0RS4tQPERERuSxO/RAREZHLYlAhIiIil8WgQkRERC6LQYWIiIhcFoMKERERuSwGFSIiInJZDCpERETkshhUiIiIyGUxqBAREZHLYlAhIiIil8WgQkRERC7r/wNHkAmoNHnrswAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "Simulation1 = StrongCastleSimulation(qys=qxs, qzs=qzs)\n", "\n", "intensity = Simulation1.simulate_diffraction(params=overlay_params)\n", "\n", "#plot\n", "import matplotlib.pyplot as plt\n", "plt.plot(qzs, intensity)\n", "plt.xlabel(r'$q_{z}$')\n", "plt.ylabel('Intensity')\n", "plt.title('Stacked Trapezoid diffraction simulation')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Fitting" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will do exactly the same steps as in the [stacked trapezoid model tutorial.](../Tutorials/stacked_trapezoid_simulation_and_fitting.ipynb) We will use the calculated intensities above as experimental and fit them with the strong castle model." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " height1 langle1 langle2 langle3 rangle1 rangle2 rangle3 \\\n", "0 10.03779 1.651405 1.480153 1.592135 1.54748 1.708779 1.546373 \n", "\n", " y_start bot_cd dwx dwz i0 bkg_cste overlay \\\n", "0 0.017236 40.000198 0.098872 0.100022 10.000337 0.099607 9.962454 \n", "\n", " top_cd \n", "0 19.999823 \n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAioAAAHACAYAAACMB0PKAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/SrBM8AAAACXBIWXMAAA9hAAAPYQGoP6dpAABi1klEQVR4nO3dd3gc1aH+8e/sqlu9WV1y73K3MWBaAFNDC+QaQgktIZQAgQC/JATIDYYEAglwyYUQDDchBkIgdAzGhmBj4265N8lWsbrV22p3fn+sLFu4qZ+V9H6eZx+8o9nddxdZfnXmzBnLtm0bERERER/kMB1ARERE5GhUVERERMRnqaiIiIiIz1JREREREZ+loiIiIiI+S0VFREREfJaKioiIiPgsFRURERHxWSoqIiIi4rNUVERERMRn9Zui8uWXX3LhhReSlJSEZVm88847HX4O27Z54oknGDlyJIGBgSQnJ/Pb3/62+8OKiIhIu/iZDtBdamtrmThxItdffz2XXnppp57jpz/9KQsXLuSJJ55gwoQJlJeXU15e3s1JRUREpL2s/nhRQsuyePvtt7n44otbtzU2NvKLX/yCf/zjH1RUVDB+/Hgef/xxTjvtNAC2bNlCZmYmGzduZNSoUWaCi4iISBv95tDP8dx22218/fXXLFiwgA0bNnD55ZdzzjnnsGPHDgDee+89hg4dyvvvv8+QIUPIyMjgxhtv1IiKiIiIQQOiqOzdu5eXX36ZN998k9mzZzNs2DDuueceTj75ZF5++WUAdu/ezZ49e3jzzTd59dVXmT9/PqtXr+Z73/ue4fQiIiIDV7+Zo3IsWVlZuN1uRo4c2WZ7Y2MjMTExAHg8HhobG3n11Vdb93vppZeYOnUq27Zt0+EgERERAwZEUampqcHpdLJ69WqcTmebr4WGhgKQmJiIn59fmzIzZswYwDsio6IiIiLS+wZEUZk8eTJut5vi4mJmz559xH1OOukkmpub2bVrF8OGDQNg+/btAKSnp/daVhERETmo35z1U1NTw86dOwFvMfnDH/7A6aefTnR0NGlpafzgBz9g6dKlPPnkk0yePJmSkhIWLVpEZmYm559/Ph6Ph+nTpxMaGsrTTz+Nx+Ph1ltvJTw8nIULFxp+dyIiIgNTvykqS5Ys4fTTTz9s+7XXXsv8+fNxuVz893//N6+++ir5+fnExsZywgkn8PDDDzNhwgQACgoKuP3221m4cCGDBg3i3HPP5cknnyQ6Orq3346IiIjQj4qKiIiI9D8D4vRkERER6ZtUVERERMRn9emzfjweDwUFBYSFhWFZluk4IiIi0g62bVNdXU1SUhIOx7HHTPp0USkoKCA1NdV0DBEREemE3NxcUlJSjrlPny4qYWFhgPeNhoeHG04jIiIi7VFVVUVqamrrv+PH0qeLyoHDPeHh4SoqIiIifUx7pm1oMq2IiIj4LBUVERER8VkqKiIiIuKz+vQcFRER8R1utxuXy2U6hvgAf39/nE5ntzyXioqIiHSJbdsUFhZSUVFhOor4kMjISBISErq8zpmKioiIdMmBkhIfH09ISIgW4BzgbNumrq6O4uJiABITE7v0fCoqIiLSaW63u7WkxMTEmI4jPiI4OBiA4uJi4uPju3QYSJNpRUSk0w7MSQkJCTGcRHzNge+Jrs5bUlEREZEu0+Ee+bbu+p5QURERERGfpaIiIiLSAQ899BCTJk3q0GNOO+007rzzTuM5+iJNphUREemAe+65h9tvv71Dj/nXv/6Fv79/DyXq31RURMQ3uOrBP9h0CpGjsm0bt9tNaGgooaGhHXpsdHR0D6Xq/3ToR0TM27scHk2CD+4xnUQGmMbGRu644w7i4+MJCgri5JNPZuXKlQAsWbIEy7L46KOPmDp1KoGBgXz11VeHHXJpbm7mjjvuIDIykpiYGO677z6uvfZaLr744tZ9vn3oJyMjg0cffZTrr7+esLAw0tLSeOGFF9pku++++xg5ciQhISEMHTqUX/3qVwNy5V8VFRHpfXu+hvWvH7y/8S2wPbDyRdj8b3O5pFvYtk1dU7ORm23bHcr685//nLfeeotXXnmFNWvWMHz4cObMmUN5eXnrPvfffz+PPfYYW7ZsITMz87DnePzxx/n73//Oyy+/zNKlS6mqquKdd9457ms/+eSTTJs2jbVr1/KTn/yEW265hW3btrV+PSwsjPnz57N582b++Mc/8uKLL/LUU0916P31Bzr0IyK9y+OB138AdaUQPRRSp0PO0oNff++nkDIDwru2mqWYU+9yM/bBT4y89uZH5hAS0L5/2mpra3n++eeZP38+5557LgAvvvgin376KS+99BLTp08H4JFHHuGss8466vM888wzPPDAA1xyySUAPPvss3z44YfHff3zzjuPn/zkJ4B39OSpp55i8eLFjBo1CoBf/vKXrftmZGRwzz33sGDBAn7+85+36/31FxpREZHeVbbTW1IAtn0AtWVQvMl7P2401O+Hdzs2UVGkM3bt2oXL5eKkk05q3ebv78+MGTPYsmVL67Zp06Yd9TkqKyspKipixowZrducTidTp0497usfOjpjWRYJCQmty84DvP7665x00kkkJCQQGhrKL3/5S/bu3dvu99dfaERFRHpX/qqDf97+CSS3/ECPGw1X/B88fyLs/BTyVkPK8X/Yi+8J9ney+ZE5xl67uw0aNKjbnxM47Cwgy7LweDwAfP3111x11VU8/PDDzJkzh4iICBYsWMCTTz7ZI1l8mYqKiPS4ncU1LN1Zyqo9+7my5BNmHfhC8WZYv8D754yTIW4kTPgerP8HfP0MXD7fUGLpCsuy2n34xaRhw4YREBDA0qVLSU9PB7zLva9cubLda55EREQwePBgVq5cySmnnAJ4r3+0Zs2aLq1xsmzZMtLT0/nFL37Rum3Pnj2dfr6+zPe/k0Skz2pwuXni38tJWfc0b7pPY5Odwc0B68ABTbYfAVYzbH3fu3N6y/D7rNu8RWXzv2F/DkRlGEov/d2gQYO45ZZbuPfee4mOjiYtLY3f/e531NXVccMNN7B+/fp2Pc/tt9/OvHnzGD58OKNHj+aZZ55h//79XVpCfsSIEezdu5cFCxYwffp0PvjgA95+++1OP19fZnyOSn5+Pj/4wQ+IiYkhODiYCRMmsGrVquM/UER82t6yOi79n2XErnue6/wW8mzoX7nvjDTGOXIBeM19RtsHHCgqCeNh2Bnes4CW/7mXU8tA89hjj3HZZZdx9dVXM2XKFHbu3Mknn3xCVFRUu5/jvvvuY+7cuVxzzTXMmjWL0NBQ5syZQ1BQUKdzffe73+Wuu+7itttuY9KkSSxbtoxf/epXnX6+vsyyO3ouVzfav38/kydP5vTTT+eWW24hLi6OHTt2MGzYMIYNG3bcx1dVVREREUFlZSXh4eG9kFhE2qOyzsV3n/uK3LIalgfdQTwtp3rOeRQ++X/YoYN5a9xzfG/F9wAo9E8l/oEsHI6W30B3fQ7/dwn4D4KfbYGgCEPvRI6noaGB7OxshgwZ0qV/mPsTj8fDmDFjuOKKK/jNb35jOo4xx/re6Mi/30YP/Tz++OOkpqby8ssvt24bMmSIwUQi0lVuj83tC9ayp6yOi8N2EO86uB4Fn/8WACt5Gt8750waNqURVLOXRfUjKV60g7vOGundb+jpEJkGFXuhYC0MPa3334hIO+3Zs4eFCxdy6qmn0tjYyLPPPkt2djZXXnml6Wj9gtFDP++++y7Tpk3j8ssvJz4+nsmTJ/Piiy+ajCQiXfTEwm18ub2EIH8HD6Vv8G5MmuL9r6vW+9+UaWBZBJ10Cx7Lj3+5T+aPi3awcFOh9+uWBXFjvH8u29W7b0CkgxwOB/Pnz2f69OmcdNJJZGVl8dlnnzFmzBjT0foFo0Vl9+7dPP/884wYMYJPPvmEW265hTvuuINXXnnliPs3NjZSVVXV5iYivmNTQSV//sJbLJ68aDiRe1oW/Tr3dweLB3iLCsCsn+D4dRmZJ3pPZf3FOxupaWz2fi2m5fBv+e7eiC7SaampqSxdupTKykqqqqpYtmxZ6xlA0nVGi4rH42HKlCk8+uijTJ48mZtvvpmbbrqJP//5yBPo5s2bR0REROstNTW1lxOLyNHYts3D723GtuGCzETO91sJrjqIGe4tJtOub9nTgqTJbR57/7mjyYgJoaS6kWc/3+ndGD3U+1+NqIgMaEaLSmJiImPHjm2zbcyYMUddee+BBx6gsrKy9Zabm9sbMUWkHRau2kxWdgFB/g4eOG8MbH7X+4XM//Ieypn4X5A6E6b9EALD2jw20M/Jry7w/ix46avdZJfWHjKioqIiMpAZnUx70kkntbkAE8D27dtbF975tsDAQAIDA3sjmoh0QMP+Ak784ExeCUhl2Sl/IzkiCHJXeL84rOU05KBwuGHhUZ/jjNHxnDoyji+2l/DbDzbzl++2FJX9OeBxg6P7VxwVEd9ndETlrrvuYvny5Tz66KPs3LmT1157jRdeeIFbb73VZCwR6aCVS94jjDpmOLbx40mB3sM19eXgFwQJE9r1HJZl8eCFY3E6LD7bUsym2jBwBoC7CSo1eioyUBktKtOnT+ftt9/mH//4B+PHj+c3v/kNTz/9NFdddZXJWCLSAW6PTf7mZa33g/KWQt433jtJk8EvoN3PNSwulPMneK+a/MJXeyCqZbkCzVMRGbCML6F/wQUXcMEFF5iOISKd9OnmQtIatsOBIzPZX3pHUgBSpnf4+W4+ZSjvri/g/Q37eGxUBsGl21rO/PlOt2UWkb7D+BL6ItJ32bbNC1/sZLwj++DG3V9AbsuISuqMDj/n+OQITh4ei9tjs7Y22rtRpyiLj8jJycGyLNatW3fUfZYsWYJlWVRUVHTptTIyMnj66ae79BzHM3/+fCIjI3v0NbpKRUVEOm31nv2U520j3KrHdgZ655RUF0DxJu8OKR0vKgA/OtV7avInhYO8G3ToR6RHfP/732f79u2mYxyTioqIdNprK/aSaXlHU6yECW2LSWQ6hA3u1POePDyW0QlhbG9uebxOURbpEcHBwcTHx5uOcUwqKiLSKdUNLj7cuO/gYZ+kyTD01IM7pM7s9HNblsV/TU9lj6elqOzPAXdz58OKHMHHH3/MySefTGRkJDExMVxwwQXs2tW2FH/zzTdMnjyZoKAgpk2bxtq1aw97ng8//JCRI0cSHBzM6aefTk5OzmH7fPXVV8yePZvg4GBSU1O54447qK2tbf16cXExF154IcHBwQwZMoS///3vx81/3XXXcfHFF/PEE0+QmJhITEwMt956Ky6Xq3Wf/fv3c8011xAVFUVISAjnnnsuO3bsaP36tw/9rF+/ntNPP52wsDDCw8OZOnUqq1atavf76AkqKiLSKR9s2EeDy8OMwD3eDUmTYMghy4Z3Yn7KoS6alEyZM5YG2x88zVB55IUgxQfZNjTVmrnZdrtj1tbWcvfdd7Nq1SoWLVqEw+HgkksuwePxAFBTU8MFF1zA2LFjWb16NQ899BD33HNPm+fIzc3l0ksv5cILL2TdunXceOON3H///W322bVrF+eccw6XXXYZGzZs4PXXX+err77itttua93nuuuuIzc3l8WLF/PPf/6T//mf/6G4uPi472Hx4sXs2rWLxYsX88orrzB//nzmz5/f5nlXrVrFu+++y9dff41t25x33nltysyhrrrqKlJSUli5ciWrV6/m/vvvx9/fv93voycYP+tHRPqmN1fnYeFhLIeMqMSOhMBwaKyCtFldev6oQQGcOTaRPdsGM8rKg7LdB5fVF9/mqoNHk8y89v8rgIBB7dr1sssua3P/r3/9K3FxcWzevJnx48fz2muv4fF4eOmllwgKCmLcuHHk5eVxyy23tD7m+eefZ9iwYTz55JMAjBo1iqysLB5//PHWfebNm8dVV13FnXfeCcCIESP405/+xKmnnsrzzz/P3r17+eijj/jmm2+YPt17ptxLL73UrosaRkVF8eyzz+J0Ohk9ejTnn38+ixYt4qabbmLHjh28++67LF26lBNPPBGAv//976SmpvLOO+9w+eWXH/Z8e/fu5d5772X06NGtWdv7PoKCgo6btzM0oiIiHba7pIbVe/Yz1FFIgLsW/IIhdhQ4/eG//g6X/C8kjO/y61w+LYUcOwEAV+nOLj+fyKF27NjB3LlzGTp0KOHh4WRkZAC0XsZly5YtZGZmtvkHeNastgV8y5YtzJzZ9jDnt/dZv3498+fPJzQ0tPU2Z84cPB4P2dnZbNmyBT8/P6ZOndr6mNGjR7frbJxx48bhdB5ctTkxMbF1JObA8x6aLyYmhlGjRrFly5YjPt/dd9/NjTfeyJlnnsljjz3W5lDY8d5HT9GIioh02D9X5wFwRVIZlOJdfdbZ8uNkSPddNXb2iDje8o8HD+zJ3sXwrg3SSG/xD/GObJh67Xa68MILSU9P58UXXyQpKQmPx8P48eNpamrq1kg1NTX86Ec/4o477jjsa2lpaV066+bAYZkDLMtqPXTVGQ899BBXXnklH3zwAR999BG//vWvWbBgAZdccslx30dPUVERkQ6xbZv3N+wD4Ozwvd6ikjylR17L6bAYnJwOuVC8by/De+RVpNtZVrsPv5hSVlbGtm3bePHFF5k9ezbgnSh6qDFjxvB///d/NDQ0tI6qLF++/LB93n333Tbbvr3PlClT2Lx5M8OHH/k7ePTo0TQ3N7N69erWQz/btm3r8josY8aMobm5mRUrVrQe+jnwvr99QeBDjRw5kpEjR3LXXXcxd+5cXn75ZS655JLjvo+eokM/ItIhmwqq2FteR5C/g/TaLO/GLpzhczzDhnjnpbiqCqlr0pk/0j2ioqKIiYnhhRdeYOfOnXz++efcfffdbfa58sorsSyLm266ic2bN/Phhx/yxBNPtNnnxz/+MTt27ODee+9l27ZtvPbaa20mswLcd999LFu2jNtuu41169axY8cO/v3vf7dOQh01ahTnnHMOP/rRj1ixYgWrV6/mxhtvJDg4uEvvccSIEVx00UXcdNNNfPXVV6xfv54f/OAHJCcnc9FFFx22f319PbfddhtLlixhz549LF26lJUrV7bOlTne++gpKioi0iEfbywEYM7wUBzFG70b007osddLTvFeTT3GrmDJtpIeex0ZWBwOBwsWLGD16tWMHz+eu+66i9///vdt9gkNDeW9994jKyuLyZMn84tf/KLNJFnwHvJ46623eOedd5g4cSJ//vOfefTRR9vsk5mZyRdffMH27duZPXs2kydP5sEHHyQp6eCE45dffpmkpCROPfVULr30Um6++eZuWd/k5ZdfZurUqVxwwQXMmjUL27b58MMPDztkBOB0OikrK+Oaa65h5MiRXHHFFZx77rk8/PDD7X4fPcGy7Q6cy+VjqqqqiIiIoLKykvDwcNNxRPo9e+VfWfDpf/hF1aX87YwGTlx2I0SkwV1ZPfeiBevghVMpsiP5zai3efbKnjnMJJ3T0NBAdnY2Q4YM6bGzPqRvOtb3Rkf+/dYcFRFpn6Y6+Ohe5nqa+cJvKFMdLQOyaT132AeAUO+ibzFUsWRrIQ0uN0H+zuM8SET6Cx36EZH2KViD5fHOEbk17AsCC1Z6t/fg/BQABsViY+FneQhsquCL7Tr8IzKQqKiISPvsPXgmw/j6VbD3a++dHpyfAoDTHyskBoA4q5KPsvb17OuJiE9RURGRdqnf7S0mLtuJhQ3NDRAQBvFHP82x27Qc/omzKliyvQS3p89OrRORDlJREZHj83hw5HsP9bwX+r2D21OmgaMX5ouEes9+SAuooaLOxbrc/T3/mtIhffi8DOkh3fU9oaIiIsdXtoNAVyX1dgDFk++AsETv9p4+7HNAy4jK9FjvhdQ+33r8i7VJ7zhwmmtdXZ3hJOJrDnxPHOlU6I7QWT8iclyunK/xB9bbw5g9NhWSHocV/wuTr+6dAGHeojI+oh7y4POtJdw7ZzR4PNBc7/OroPZnTqeTyMjI1uvLhISEYFmW4VRikm3b1NXVUVxcTGRkZJtrEXWGioqIHFfZ1v+QAGzxG8N1ieGQdBGMPXxlyx7TMqKS6l+DZcGWfVUUVjaQ8OEPIfs/cMdaCI3rvTzSRkKC98KRB8qKCEBkZGTr90ZXqKiIyHE5W+ankDLDzG/LLUUlsKGESamRrN1bweKtRczd/QW4aqFwAwz/Tu/nEsB7IbzExETi4+NxuVym44gP8Pf37/JIygEqKiJybHXlxDXsASBt4qlmMrRMpqWmmDPGxLN2bwUrN+9krqu2ZXuRmVzShtPp7LZ/nEQO0GRaETmmop3rANhrxzNj3AgzIVpGVKgp4vTR3tKSm7314NerCw2EEpHeoKIiIse0Y/dOAGoCBhMW1LXZ+512YESloYKxcYFEDwogpvmQ+RA1mhsh0l+pqIjIMRXn7wbAGZFoLkRQJDgDAHDUlTBrWAwp1iFL6ddoREWkv1JREZGj8nhs6krzAAiPTzMXxLIOOfxTzMnDY0m2Sg9+vVpzVET6KxUVETmq7cXVhDd7C0FcUobZMK0Taos4aVisRlREBggVFRE5qqU7yxhseZer94tMMhvmkAm1aTEhDPErP/g1jaiI9FsqKiJyVF/vKmUwLdfVCTM4RwXanKKMbbcdUXHVQmONmVwi0qNUVETkiJrdHlbsLiPBahm5MF5UDo6o0FBBkMe7hkoT/ge3i0i/o6IiIkeUlV+J1VhJkNWy0qjxotIyolKZBxW5AJTa4eR5YrzbtZaKSL+koiIibbkaoDCLZTtLSWiZn0JwFPgHmc2VNMX73z1LoWwHAGV+gykh0rtdE2pF+iUVFRFp66N74c8n07B14SGHfQxPpAVInAShCdBUA+v+AYArLIUSO9L7dS36JtIvqaiISFuFWQDEFC1rPeOHsK5fAbXLHA4YdY73zzs/BSA4bgjFB4qKDv2I9EsqKiLSVlUBAKM8u0n1q/RuCzc8P+WAkee2uRufOqK1qDRV7jMQSER6moqKiBzU3NR6CGW8I5vxYS2n/JqeSHvA0FPBL7j1btjgoa1nA9WU5JlKJSI9SEVFRA6qKQRsAMKseia5N3m3+0pR8Q+GYWccvB+ZRmyCd2n/5iod+hHpj1RUROSgqraHT6Lrc7x/CPeBybQHHJinAhCZSnrGUAACG0qO8gAR6ctUVETkoKr8I2/3hcm0B4w6D0JiICETAsMYO3IEABF2FXX1dYbDiUh3U1ERkYNaJtI22n5tt/vC6ckHDIqFW1fCDz8EICkxiWacAGzavtNkMhHpASoqInJQS1FZ6hl/cJvl9JYDXzIoBgLDALAcTqr9ogHYsWuXyVQi0gNUVETkoJZDP0s943A7W86uCUsAh9NgqONzh3iX1y8s2GM4iYh0NxUVEWnVXOktKvl2HHbCBO9GX5qfchQBUd5DU9Ulebg9tuE0ItKdVFREpFXzfm9RcUQk4Zcy1bvRV05NPobQuHQABrsL2V5UbTiNiHQno0XloYcewrKsNrfRo0ebjCQycHnc+NcXATA4ZRhMuBwiUmHcJYaDHZ8j3vtzY4SVx5q9+w2nEZHu5Hf8XXrWuHHj+Oyzz1rv+/kZjyQyMNWW4LTdNNsORgwbCilD4a6NplO1T/wYAEZaeXywZz9XzUw3HEhEuovxVuDn50dCgu8fAxfp79yVBTiBEiKZnO5jZ/kcT5y3qKQ6StiSUwBMMhpHRLqP8TkqO3bsICkpiaFDh3LVVVexd+/eo+7b2NhIVVVVm5uIdI+Cvd41SIqJZuTgMMNpOmhQDJ6QOAD89++ktKbRcCAR6S5Gi8rMmTOZP38+H3/8Mc8//zzZ2dnMnj2b6uojT4abN28eERERrbfU1NReTizSfxXl7QagITgBp8MynKbjDsxTGenIY+3eCrNhRKTbGC0q5557LpdffjmZmZnMmTOHDz/8kIqKCt54440j7v/AAw9QWVnZesvNze3lxCL9V1Wxdw0Sv8gUw0k6KX4s4J2nsnqPJtSK9BfG56gcKjIykpEjR7Jz55GXwQ4MDCQwMLCXU4kMDO4K76nJ4Ql9dCLqgREVK4//0Zk/Iv2G8Tkqh6qpqWHXrl0kJvr+ug0i/UllvYvQRu/VhxNThhpO00ktE2pHOPLYmF+phd9E+gmjReWee+7hiy++ICcnh2XLlnHJJZfgdDqZO3euyVgiA05WXiUJVhkAoXFphtN0UsuISrJVhqOpmp3FNYYDiUh3MFpU8vLymDt3LqNGjeKKK64gJiaG5cuXExcXZzKWyICzIbeUJKvceyfch66U3BHBURDqXepghJXP+rwKs3lEpFsYnaOyYMECky8vIi2KdmcRaLlocg4iIKKPjqiAd1SlppCRjjw25FVwxTSdGSjS1/nUHBURMcPatx6Ahthx4OjDPxZa5qnMcGxlc26p4TAi0h368E8kEekOhZUNpDduByA4bYrhNF00eBwAlzn/w6ulc2n+9BHDgUSkq1RURAa49XkVjHPkAOCf0seLyrhLsKffRAmRhFr1+C19Eup1qrJIX6aiIjLArd9bzjgrx3sncaLRLF0WGIp1/hPcm/o6OZ7B3m15q81mEpEuUVERGeCKczYyyGqk2RkMsSNMx+kWmalRrLZb3kveN2bDiEiXqKiIDGAej42zKAuAptix4HAaTtQ9JqZGssYz0nsnd4XZMCLSJSoqIgNYdlktw5u9l6wI7OsTaQ+RmRLJGo93RMXOWwket+FEItJZKioiA9jG/ErGt8xPcSZNMpqlO8WFBVIdNpwaOwirqRaKt5iOJCKdpKIiMoBtyq9gnCPbeydxktEs3W18agxrPcO9d3T4R6TPUlERGcCK924j3KrH7QiAuFGm43SriamRrGmdULvSbBgR6TQVFZEByrZtmou2AtAUNQKc/oYTda+JKRGaUCvSD6ioiAxQueX1xLsKAAiMH244TfcbnxLBWs8w753y3VBTYjaQiHSKiorIALWxoJI0qwgAR8xQw2m6X3iQP7Fxg9nlSfRuKFxvNpCIdIqKisgAtTG/koyWokLUELNhesiklEh22UneO2W7zYYRkU5RUREZoDYWVJFuFXrvRPe/ERWAzJQIdtstIyrlu8yGEZFOUVERGYBs22ZrXjkpVql3Q38tKqmR5NgJANhlOw2nEZHOUFERGYAKqxoIrC/A33JjOwMhLNF0pB4xNjGcvXjfm7tkh+E0ItIZKioiA9DG/CrSrWIArOgh4OifPwqC/J044rxrqTgrc6G5yXAiEemo/vnTSUSObNfn8K+b2ZGzl4wD81P66UTaA5KTM7xL6eOB/Tmm44hIB6moiAwkXz4JG14nafvfSD9wxk8/nZ9ywLiUiNZ5Kmieikifo6IiMpBU5QEwtnLJIUWlf4+ojEuKIPtAUdGZPyJ9jp/pACLSS2wbqr2He0baOcQ5Dpzx07+LypjEMJa1FJX6wm0EG84jIh2jERWRgaJ+PzQ3tN6Nsmq8f+jnh35CAvyoDfWWsYbC7YbTiEhHqaiIDBTV+w7fZjkhIrX3s/SygMHeM3/8K7INJxGRjlJRERkoWopKhTMGj215t0Wm9burJh9JTNoYAEKbiqGxxnAaEekIFRWRgaLKW1S22mmstEd5t/Xzwz4HjEhPo9wO9d4p1zV/RPoSFRWRgaJlRCWnKYLXm0/zbks/0VyeXjQ2KZzslmv+1GqeikifoqIiMlC0FJUiolgVeQ7cvgZO+qnhUL0jMiSAIv8UACp3LDecRkQ6QkVFZKBoOfRTZEczPjkcYoYNiPkpB+yJPQ2AqB1vgqvh2DuLiM9QUREZKKoLACi0oxiXFGE4TO9zj5hDvh1DsKsCNr9jOo6ItJOKishA0bLYW5EdzdikcMNhet/YlGhea/6O9843L5oNIyLtpqIiMhC4Xdg13qslF9lRjE0ceEVlXFIEr7tPp8l2Qv4qKFhrOpKItIOKishAUFOEhU2T7cQTHE18WKDpRL0uPiwQQuP50DPTu2HVy2YDiUi7qKiIDAQtE2mLiWJUYiSWZRkO1Pssy2JcUjgfuVuKStFGs4FEpF1UVEQGgpZTk4vtSEYnhhkOY8745HCK7CjvnZZDYSLi21RURAaClqJSaEczJmHgzU85YFxSBKW0nPFUU+y9orSI+DQVFZEBwK7ynppcZEcN7BGVpAhK7Jai4m6ExiqzgUTkuFRURAaAhvJ8wDtHZeTggVtUUqODCQgKodoO9m6oKTEbSESOS0VFZABoKM8DwBOaSJC/03AacyzLYmxiOKV2y+GvWs1TEfF1KioiA0HLHJXg2FTDQcwbn/yteSoi4tNUVET6O9smuKEIgNiEdMNhzBuXFE7pgXkqtTr0I+LrVFRE+rvaUoI89Xhsi6S0EabTGDf2kKJi1xQZTiMix6OiItLPuYq3AZBvxzIyNc5wGvOGxYVSYXmLSm35PsNpROR4fKaoPPbYY1iWxZ133mk6iki/UrpnEwB7HEkkRwYbTmOev9OBFTYYUFER6Qt8oqisXLmS//3f/yUzM9N0FJF+pyZ/CwBVIRkDcun8IwmNTgTAU63JtCK+znhRqamp4aqrruLFF18kKirKdByR/qd0BwDu6OGGg/iOmIQUAPzqSw0nEZHjMV5Ubr31Vs4//3zOPPPM4+7b2NhIVVVVm5uIHFtoTQ4AIYmjzQbxIcnJaQCENpcbTiIix+Nn8sUXLFjAmjVrWLlyZbv2nzdvHg8//HAPpxLpR5qbiGv2Lp8fP3S84TC+Y2jGUACCaaS6aj9h4RrNFfFVxkZUcnNz+elPf8rf//53goKC2vWYBx54gMrKytZbbm5uD6cU6dv252/HDw+1diBDh+jU5AMiI6NoIACA7Jwcs2FE5JiMFZXVq1dTXFzMlClT8PPzw8/Pjy+++II//elP+Pn54Xa7D3tMYGAg4eHhbW4icnSFu7MAyHcmMyjI33AaH2JZVPtFA5Cft8dwGBE5FmOHfr7zne+QlZXVZtsPf/hDRo8ezX333YfTOXCvRyLSXapbzvipDNGKtN/mCoqFmkJKC/NMRxGRYzBWVMLCwhg/vu0x80GDBhETE3PYdhHpHLvljJ/maB32+TZHWDzUQE1ZgekoInIMxs/6EZGec+CMn2Cd8XOYkCjvWirN1UW4PbbhNCJyNEbP+vm2JUuWmI4g0m80uz0kunLBgvgh40zH8TmhMd6iEumpILu0luHxoYYTiciRaERFpJ/am59HtFUNQEKGDqd+myPUu4x+rFXJ5n1ak0nEV6moiPRTBTu9k9VLHbE4gjRacJhQ7wUaY61KtqioiPgsFRWRfqoudz0AZYOGGU7iowbFAxCLioqIL1NREemn/Eo3A+CKHWs4iY8K9RaVeKuCHQVlhsOIyNGoqIj0U7E12wEITp1oOImPikzDDo5mkNXI3Pp/UFbTaDqRiByBiopIP1RZ28gQz14ABo+YajiNj/ILxLrgKQB+4nyX/HWLDAcSkSNRURHph7J3bibMqqcJP0KTxpiO47vGXcyysHNwWDZD/3MXNNaYTiQi36KiItIPlWevBaAwIB2cusbPsawf/wAFdjShjYWQ/aXpOCLyLSoqIv2QZ5/31OSqCK1IezzDUxP5xtPyOZVsMRtGRA6joiLSDw3avxUAR4IWejueMYlhbPekAuAu3Gw4jYh8m4qKSD/j9tgkNe4GIHLIZMNpfF9yZDC5/t6rS7v2bTKcRkS+TUVFpJ/JLSwhlSIABo+YZjiN77MsCzvOu9aMf8VOcDcbTiQih1JREeln9m1fjcOyKXNE4wyLMx2nT4hLGU6tHYjT44Ly3abjiMghVFRE+pmalqXzS0OGG07Sd4xJimSHney9U6x5KiK+REVFpJ+xSrwr0rqiRxlO0neMSQxvnVBrq6iI+JROFZVrr72WL7/UegMivii4Zg8AIYkjDSfpO0YMDmUH3qLSkL/RcBoROVSnikplZSVnnnkmI0aM4NFHHyU/P7+7c4lIJ1Q3uIhv3gdAXJrWUGmvIH8nNeEjAPAUay0VEV/SqaLyzjvvkJ+fzy233MLrr79ORkYG5557Lv/85z9xuVzdnVFE2ml7YSVpVjEAYUkaUekIv8RxAARX7QFXg+E0InJAp+eoxMXFcffdd7N+/XpWrFjB8OHDufrqq0lKSuKuu+5ix44d3ZlTRNphT/YuAi0XzTghPMV0nD4lOXUIlXYIDtxQpp9fIr6iy5Np9+3bx6effsqnn36K0+nkvPPOIysri7Fjx/LUU091R0YRaaf9ed4VaauCksDpZzhN3zImKYJttneeCjr8I+IzOlVUXC4Xb731FhdccAHp6em8+eab3HnnnRQUFPDKK6/w2Wef8cYbb/DII490d14ROYam4l0AuMLTDSfpe8YmhrPD4x2FchWqqIj4ik79ypWYmIjH42Hu3Ll88803TJo06bB9Tj/9dCIjI7sYT0Tay7Zt/KtywIKAeK2h0lFxYYGUBySCB6qLsok2HUhEgE4WlaeeeorLL7+coKCgo+4TGRlJdnZ2p4OJSMfk7a8nybMPnBCepDVUOsMRlQpl0Fy+x3QUEWnRqUM/ixcvPuLZPbW1tVx//fVdDiUiHbdlXxXplvcaP86YoYbT9E2hg4cBEFCTZziJiBzQqaLyyiuvUF9ff9j2+vp6Xn311S6HEpGO23pIUSFaRaUz4lO9a6mEuUrBraUWRHxBhw79VFVVYds2tm1TXV3d5tCP2+3mww8/JD4+vttDisjx5eXvJdRqwMbCitJk2s4YNmQIjbYfgVYznsp8HNEZpiOJDHgdKiqRkZFYloVlWYwcefhiUpZl8fDDD3dbOBFpv/pC79ofjSGJBPkFGk7TNw2NCyOPWDIopDhvJwkqKiLGdaioLF68GNu2OeOMM3jrrbeIjj44Lz4gIID09HSSkpK6PaSIHFt9kxv/qj3gD47YYabj9Fl+TgcVAYPBVUjR3h0kZJ5pOpLIgNehonLqqacCkJ2dTVpaGpZl9UgoEemY7UXVrfNTAmI1P6UrGgelQMV6aot2m44iInSgqGzYsIHx48fjcDiorKwkKyvrqPtmZmZ2SzgRaZ8t+6pI00TabuEXnQ4V4N6/13QUEaEDRWXSpEkUFhYSHx/PpEmTsCwL27YP28+yLNxud7eGFJFj21pYzdVWywhA7AizYfq4iMShsBuC6gpMRxEROlBUsrOziYuLa/2ziPiOitxNDHPsw2P54ciYbTpOn5bQcopybHMRlfUuIoL9DScSGdjaXVTS09OP+GcRMcu2bdJKvgCgNvkkwoLCDSfq20IHDwEgySpjbUEFJwyLM5xIZGDr9IJvH3zwQev9n//850RGRnLiiSeyZ4+WnhbpTYVVDcz2fANA0PgLDKfpB8KTceMg0GomJ0cTakVM61RRefTRRwkODgbg66+/5tlnn+V3v/sdsbGx3HXXXd0aUESObVd2DlMt7xoq/mPOM5ymH3D6URvgXbiyNH+n4TAi0qmLEubm5jJ8uPfqrO+88w7f+973uPnmmznppJM47bTTujOfiBxH45aPcFg2uYEjSI1IMR2nX2gOT4HSQuqKNR9PxLROjaiEhoZSVlYGwMKFCznrrLMACAoKOuI1gESk58TmLwKgMPEMw0n6j4AY7zw8qzIPl9tjOI3IwNapEZWzzjqLG2+8kcmTJ7N9+3bOO8873Lxp0yYyMjK6M5+IHIu7mZE1K71/HnWu2Sz9yKC4IbANEuwSdpfUMiohzHQkkQGrUyMqzz33HLNmzaKkpIS33nqLmJgYAFavXs3cuXO7NaCIHF1j8XaCaaTWDiR1zEzTcfoNKyoNgBSrhM37Kg2nERnYOjWiEhkZybPPPnvYdl2QUKR3Fe1YQxqw00ojMyLYdJz+I34MAJmO3fxvQRWXTDacR2QA61RRAaioqOCbb76huLgYj+fgMVzLsrj66qu7JZyIHFtt7noASkOG69pb3SlxEs2OIGI81VTuzQLGmk4kMmB1qqi89957XHXVVdTU1BAeHt7mB6SKikjvcRZvBqAxZozhJP2MXwD1CdMIK/iKyOJvsO0rVARFDOnUHJWf/exnXH/99dTU1FBRUcH+/ftbb+Xl5d2dUUSOIrLGu35KULIuBNrdgkacAsCE5o0UVTUaTiMycHWqqOTn53PHHXcQEhLSpRd//vnnyczMJDw8nPDwcGbNmsVHH33UpecUGSjs+gri3d4rJg8eMcVwmv7Hf6j3mkkzHVvYUqAJtSKmdKqozJkzh1WrVnX5xVNSUnjsscdYvXo1q1at4owzzuCiiy5i06ZNXX5ukf6uYs8GAArsaIamaqG3bpc8lSYrgDirkoLdG0ynERmwOjVH5fzzz+fee+9l8+bNTJgwAX//tlcX/e53v9uu57nwwgvb3P/tb3/L888/z/Llyxk3blxnookMGKW71hIF7PUbQlKA03Sc/scvkJKITJIrVuHYsxQ4y3QikQGpU0XlpptuAuCRRx457GuWZeF2uzv8nG63mzfffJPa2lpmzZp1xH0aGxtpbDx4rLiqqqrDryPSXzQVZAFQFT7ScJL+qzl1FlSsIr58tekoIgNWpw79eDyeo946WlKysrIIDQ0lMDCQH//4x7z99tuMHXvkUwHnzZtHRERE6y01NbUz8UX6heDyrd4/DNboY0+JHOO9LME4Vxa1DS7DaUQGpk4VlUM1NDR06fGjRo1i3bp1rFixgltuuYVrr72WzZs3H3HfBx54gMrKytZbbm5ul15bpM+ybQbX7wIgLH2S2Sz9WMSIWTTjIMHaz67sXabjiAxInSoqbreb3/zmNyQnJxMaGsru3bsB+NWvfsVLL73UoecKCAhg+PDhTJ06lXnz5jFx4kT++Mc/HnHfwMDA1jOEDtxEBqKmsj0Moo4m20nqCJ2a3GP8gyny947cluzs+gkEItJxnSoqv/3tb5k/fz6/+93vCAgIaN0+fvx4/vKXv3QpkMfjaTMPRUQOV9zyj2Y2KSTHqLD3pMrwUQC48nXmj4gJnSoqr776Ki+88AJXXXUVTufBsw0mTpzI1q1b2/08DzzwAF9++SU5OTlkZWXxwAMPsGTJEq666qrOxBIZMKqzvZM79wWP0IqpPcxKGA/AoP1bDCcRGZg6ddZPfn4+w4cPP2y7x+PB5Wr/hLPi4mKuueYa9u3bR0REBJmZmXzyySecdZZOAxQ5FmeR94yf2pjxhpP0fxFDpsAmSGrYidtj43SoGIr0pk4VlbFjx/Kf//yH9PT0Ntv/+c9/Mnly+y8z2tH5LCLiFVPt/e0+MHWS2SADwOAR0wDIYB85haUMS4oznEhkYOlUUXnwwQe59tpryc/Px+Px8K9//Ytt27bx6quv8v7773d3RhE5VG0pMe5SAOJHTDccpv9zhidQYUUSSQX529cwLGmO6UgiA0qn5qhcdNFFvPfee3z22WcMGjSIBx98kC1btvDee+/psI1ID6vY3TKR1pPAiNREw2kGAMuiZNAIAGr3rDUcRmTg6dSICsDs2bP59NNPuzOLiLRD+a5VRAI5AcMYoqXze4UrbizUrMSvRNchE+ltnRpRGTp0KGVlZYdtr6ioYOjQoV0OJSJH5ylYD0BlhFak7S0hLXOBYmt3mA0iMgB1qqjk5OQccan8xsZG8vPzuxxKRI4ufL935WYrUQu99ZbBI70Taod7ciiuqjOcRmRg6dChn3fffbf1z5988gkRERGt991uN4sWLSIjI6PbwonItzRUEe/KAyBq2DTDYQaO4MQxNOFHmFXPph2biZ+qz16kt3SoqFx88cWA9wrJ1157bZuv+fv7k5GRwZNPPtlt4USkraaCDQQABXY0w4dkmI4zcDj92Rc0jPSGbdRtWwwqKiK9pkOHfg5cITktLY3i4uI2V01ubGxk27ZtXHDBBT2VVWTAK9u5EoDt1hASI4IMpxlYCpPOBiA5V0swiPSmTs1Ryc7OJjY2truziMhxuPZ8A0Bx6Fgtnd/LnJmXAzCifj1UFRhOIzJwdPr05EWLFrFo0aLWkZVD/fWvf+1yMBE5XHjJGgDqE7TQW28bPnIM33hGMcOxjfo1bxB82p2mI4kMCJ0aUXn44Yc5++yzWbRoEaWlpezfv7/NTUR6QGU+kU2FuG2LsGEnmE4z4ESGBPBl4GkAuDe8YTaMyADSqRGVP//5z8yfP5+rr766u/OIyFHYuSuwgC12OiNSE0zHGZD2JZ+DK/svhJZvgpLtEDfSdCSRfq9TIypNTU2ceOKJ3Z1FRI6hftdSAFbbIxkxONRwmoFpSFoqX3pa1q/Z/G+zYUQGiE4VlRtvvJHXXnutu7OIyDE071kBQG5oJkH+WjrfhHHJESzztKwIXLjebBiRAaJTh34aGhp44YUX+Oyzz8jMzMTf37/N1//whz90SzgRadFUS2i5d0XapkRNpDVlfFIEL9ppAHgKN3XuNz0R6ZBOFZUNGzYwadIkADZu3NideUTkSPJX48BNgR1NQtpw02kGrLiwQMoGDQcXWPuzoakWAgaZjiXSr3WqqCxevLi7c4jIseR6D/us8YxkTGK44TADW0pyGiXZEcRZlVC8FVKmmo4k0q91qKhceumlx93HsizeeuutTgcSkcO5936DE1jlGcktKipGjUuOYMuuNOKcWVC0UUVFpId1qKgcehFCEek97qItOIG8wOHEhwWajjOgjU8KZ6udxilkQdEm03FE+r0OFZWXX365p3KIyNE0N+FXnQ9AQPxwLZ1v2LjkCD7xpALgKdKEWpGepr9jIr6uYi8OPNTZgSQmZ5hOM+AlRQSRHzgUAE/hRrBtw4lE+jcVFRFfV74bgD32YMYk6fCraZZlEZI8lmbbgV9jBVTvMx1JpF9TURHxcZ6yXQDk2IMZl6SJtL5gZHIc2Xai907RZrNhRPo5FRURH1ezbzsAuVYCw+O1dL4vGJ8czlbbO0+FIq0lJdKTVFREfFxD0U4AGsMy8Hfqr6wvGJ8UwRbPgRVqVVREepJ+6on4OL/KbAAC47Uira9Iiw4hx887odaVu9pwGpH+TUVFxJe5mwlvKAAgJm2M4TBygMNh0ZAwFY9tEVi5G6qLTEcS6bdUVER8WWUufrhptP0ZMlQjKr5kSGoKW1suUMiepWbDiPRjKioiPqwibxsAe+x4xiRGmg0jbYxPDmeFZ7T3joqKSI9RURHxYUV7tgBQGpBMcIDTcBo51ITkSJZ7vIfj7OyvDKcR6b9UVER8WF2h99TkxrAMs0HkMENjB7HZfzwAVulWqC01nEikf1JREfFhjv3eM37844YZTiLf5nBYpKakss2T4t2gwz8iPUJFRcSHRdTnAhCVOtpwEjmSzJRIVrQc/iFHRUWkJ6ioiPioypoGEj3e015Th403nEaOZFJqxMGiohEVkR6hoiLio3ZvX0+g5aKBAMIThpiOI0fQZkSlaCPUlJgNJNIPqaiI+Kj9u9cAUBA4DBw648cXJUYEYYXFs9GT4d2wY6HRPCL9kYqKiI/yFGwAoDZKK9L6KsuymJgSwSLPFO+G7R+bDSTSD6moiPio8MqtAPinTDScRI5lYkokn7lbisquz6G50WwgkX5GRUXEB9U3uclo3g1A/PBphtPIsWSmRrLRzqDMioKmGsjR4m8i3UlFRcQH7di9k3irAg8WUUMmmY4jxzAxJQIbBwtdk7wbtn9iNI9If6OiIuKDinesBqDQLwUrMNRwGjmWyJAA0mNCDpmn8hHYttlQIv2IioqID3LlrwegMmKU4STSHhNTIlnqGUezFQAVe6Fkq+lIIv2GioqIDwop916M0ErMNJxE2iMzJYJ6gsgOGOndULzFbCCRfsRoUZk3bx7Tp08nLCyM+Ph4Lr74YrZt22YykohxLreH5MadAEQNmWI4jbTHpNRIAHKawr0baorMhRHpZ4wWlS+++IJbb72V5cuX8+mnn+JyuTj77LOpra01GUvEqF35JQyhAIA4nfHTJ4xLisDpsMh1tRSV6kKzgUT6ET+TL/7xx20XR5o/fz7x8fGsXr2aU045xVAqEbPytq5itGVT4YgkMjzBdBxph+AAJyMHh1FcHOndoBEVkW7jU3NUKisrAYiOjjacRMScxpzlAJSGjQXLMpxG2mtiSgTFdqT3jkZURLqN0RGVQ3k8Hu68805OOukkxo8/8pViGxsbaWw8uOpjVVVVb8UT6TVRpasAaE45wXAS6YiJqZF8sDrKe0cjKiLdxmdGVG699VY2btzIggULjrrPvHnziIiIaL2lpqb2YkKRntfQ1MzIxo0ARI3R4c++JPOQERVbIyoi3cYnisptt93G+++/z+LFi0lJSTnqfg888ACVlZWtt9zc3F5MKdLzdm9bT6xVRSP+xI+aZTqOdMDIwWFU+XkPW1v15brmj0g3MXrox7Ztbr/9dt5++22WLFnCkCFDjrl/YGAggYGBvZROpPft3/IFADmBoxnlH2Q4jXSEv9NBSmIyjUV+BFrN3sM/kWmmY4n0eUZHVG699Vb+9re/8dprrxEWFkZhYSGFhYXU19ebjCVijH/+CgAq46YaTiKdMSUjmhIivXeqNU9FpDsYLSrPP/88lZWVnHbaaSQmJrbeXn/9dZOxRIxJrloHQMDQk8wGkU6ZkhZJyYEzf2o0T0WkOxg/9CMiXjWluSTbhXhsi9TM00zHkU6YkhbF2pai0rC/AB28E+k6n5hMKyKQv2ExALsc6cTExhtOI50RHx5EXUAsAMUFewynEekfVFREfER99jcAFIRPNJxEuiIgMgmA6pI8w0lE+gcVFREf4SzbDoA1eJzhJNIVUYO96zs1V+0znESkf1BREfER0fU53v9mHHllZukbklIyAAioL8Hj0Tw8ka5SURHxAaX7K0j0FAOQMWqy4TTSFclp3vWgYuz97C6tMZxGpO9TURHxAbu3rsNh2VQRSmh0ouk40gX+Ed7/fzFUsjanzHAakb5PRUXEBxRne6/vUxacoSsm93WD4vDgwGnZbN21y3QakT5PRUXEBzQXbQXAFT3CcBLpMocTV1AMAHl7c8xmEekHVFREDLNtm5Aq72/eg5LGGE4j3cEvIgGApsp9lFTr4oQiXaGiImJY3v56Ut3eNTfihmQaTiPdwRnunacSb1WwKqfccBqRvk1FRcQk22ZDbjlDLe+aGwEJow0Hkm4ROhiABMpZmbPfcBiRvs3otX5EBiyPB165AJpqKYj5GYGWC5cVgH9kmulk0h0Ge9fCOdu5ivuzdeaPSFdoREXEhOp9sGcp7FvHd7Y+DEBtaAY4nGZzSffIvALbGcQ4xx4CC1dR09hsOpFIn6WiImJCZW7rH4c27wTAMViHffqNkGisCd8D4AfOhazdq8M/Ip2loiJiQsXewzaF6Yyf/mXGjQCc51jBpm07DIcR6btUVERMqNgDQEH0CdTYQQBYg1VU+pWkyZRGZBJguYnc9rrpNCJ9loqKiAktIyrrHaP5ieunrE38Pow6z3Ao6W7uqdcDML1qIQ0ut+E0In2TioqICS1FZU1VGF96JlJ9+m/BL9BwKOlu8dMvxo2DYVYBm7duMR1HpE9SURExoaWobKiOxLJgUlqk2TzSI6zgKPYGjQKgZP3HhtOI9E0qKiK9zeOBSu9KtHl2LCPjwwgP8jccSnpKbfJsAAblfWU4iUjfpKIi0ttqisDdhAcnhUQzJT3SdCLpQdETzgZgdP0aGppchtOI9D0qKiK9reWwT6kzFjdOJqdFGQ4kPSlx/CnUEUSsVcnW9ctNxxHpc1RURHpbS1HJdsUAMDVdRaU/s/wCyR40CYDKTQvNhhHpg1RURHpbyxoqez2xRIX4MzR2kOFA0tMa004BIKJA81REOkpFRaS3tYyo5BPLtIxoLMsyHEh62uBJ5wIwujGL+vp6w2lE+hYVFZHe1lJU8uw4ZmREGw4jvSFpxCSqCCXIcrFpwzem44j0KSoqIr3MbrkgYZ4dx/QhKioDgeVwUDxoBAAFW1RURDpCRUWkN3k82BXeolLiHMy4pHDDgaS3WAmZALgLNhhOItK3qKiI9KbaYhzuRpptB4mpw/B36q/gQBE/choAiQ07KK5uMJxGpO/QT0mR3tQyP6WQaKYOiTMcRnpTWPoUAMZae1i6o8RwGpG+Q0VFpDcVbQRgtyeRGZqfMrDEjqTZ8ifcqmPT5o2m04j0GSoqIr2odtfXAKxnBJN1IcKBxS+AhsiRAFRlr8G2bcOBRPoGFRWRXuTZ6z3jY39UJiEBfobTSG8LSp0IQHLjTrYX1RhOI9I3qKiI9Ja6csJqcwAIHz7LbBYxwi95EuCdp/LF9mKzYUT6CBUVkd6StwqAXZ5EJo0aajiMGJEwAYCxjhwWb9WEWpH2UFER6SWVO5cBsM4ewXStSDswDR4HQLJVxo6cPVQ1uAwHEvF9KioivaRht3cibUlkJoMCNT9lQAqKgKgMAMayi692lJrNI9IHqKiI9AaPm4hy74qkQUNmGg4jRg05FYALHV+zeKvmqYgcj4qKSC+wS7YS5Kmj1g5k+DgVlQFt4lwAznOuYPnWvXg8Ok1Z5FhUVER6QfnWpQBssIcxdUis4TRiVNoJ2FFDGGQ1Mr3+KzYVVJlOJOLTVFREesH+nSsAKAwbT3CA03AaMcqysCZdCcBlzv+waGuR4UAivk1FRaQXOEq2ABCQPNFwEvEJmd8H4CTnJtZtyDIcRsS3qaiI9DCP283ghl0ApI+dbjiN+ISodFxpJwMwpfxdcsvrDAcS8V1Gi8qXX37JhRdeSFJSEpZl8c4775iMI9IjdmzfzCAaaLL9GDV2kuk44iP8Z94IwPXOj/li7WbDaUR8l9GiUltby8SJE3nuuedMxhDpUbs2eeenFAWm4x8QaDiN+IwxF1EaNoZQq4Go1X8ynUbEZxlddercc8/l3HPPNRlBpMfV7PGun9IUO8ZwEvEpDgfWmQ/B29/nzNoPKMvbQUzKCNOpRHyO5qiI9KC6pmYGVW4DICpjktkw4nNiMuew3n8igVYzlR//t+k4Ij6pTxWVxsZGqqqq2txEfNmK3eWMYi8AUUMmG04jPsey2D32VgDiCxaBrcXfRL6tTxWVefPmERER0XpLTU01HUnkmJZuyyPDKgTAarkgncihMmedRYPtT6inmoq8LabjiPicPlVUHnjgASorK1tvubm5piOJHFPutnX4WR6aAiIgLMF0HPFBwxKi2ek/EoDNKz41nEbE9/SpohIYGEh4eHibm4ivyimtJbTCOz/FkTAeLMtwIvFVzUne9XXqdy0znETE9xgtKjU1Naxbt45169YBkJ2dzbp169i7d6/JWCLd4vOtxYxyeEf9/BLGG04jvixt0ukApNZmUVBRbziNiG8xWlRWrVrF5MmTmTzZO8nw7rvvZvLkyTz44IMmY4l0XXUhISuf4SKn92KEDB5rNo/4tOhR3lVqRzry+WTVVsNpRHyL0XVUTjvtNGzNcpf+xtWA/T8n8l/1ZWCBxz8Ux5BTTKcSXzYolqqQdMLr9pC9bgmcqTPERA7oU3NURPqEfeuw6suoskP4fcBPsO7aCNFDTacSHxc45AQAYivWsa2w2nAaEd+hoiLS3XK/AWC5Zwz1E67GCokyHEj6gsAhswCYau3g9ZU6o1HkABUVkW5m560EYI1nBGeMjjecRvqM1JkATHLs5OM1O2hsdhsOJOIbVFREupNt49rjvQjhFr/RzBgSbTiQ9Blxo7Ej0xhkNXKL61U+21xsOpGIT1BREelOlXkE1BXhsp3EjTqBAD/9FZN2cjiwvvssAFf7fcbm//zLcCAR36CfoiLdyG6Zn7LFTuPMzAyzYaTvGXoqVRNvBOCa4t9RULjPcCAR81RURLpR2TbvuilZjOCUkXGG00hfFH7Bf1Pgl8Jgq4Kc935nOo6IcSoqIt2oKcc7P6UxYRohAUaXKZK+yj+Yomk/B2B8/us01Ow3HEjELBUVke7iaiCuxruqaOrEUw2Hkb4s86wfkG2lEE4t2957ynQcEaNUVES6Se7mr/GnmVI7ghOmTDEdR/owp9NJzugfAZC+fT52U63hRCLmqKiIdJM9qz4CIGfQBMKCAwynkb5uyvk3kmfHEWlXsvvTF0zHETFGRUWkG9i2TWj+VwD4jzjDcBrpDyJCQ9iQciUAzev/aTiNiDkqKiLdYN2ufMa6vfNTRp14keE00l9MOvtqAEY0biJr2w7DaUTMUFER6QYbl35IgOWmzD+RoMHDTceRfiIpfQS5QaNwWDarFv7DdBwRI1RURLrI5fbgzP4cgKZ0ne0j3St0oneELrX4c7LyKg2nEel9KioiXfTl9hKmezYAED/xHMNppL+JmnIxALMdG3l+4XqzYUQMUFER6aLPVqxlhCMfDw6cwzSiIt0sfiyu8HQCLRfunYtYtqvUdCKRXqWiItIFJdWNOHcuBKAxPhNCdLVk6WaWhf+4CwH4nvNLHn1vI26PbTiUSO9RURHpgg0f/JlfOV8BIHjc+YbTSL817lIAznKu4dfl9/Lxf5YbDiTSe1RURDrJ89Wf+M7WBwm0mslL+A7Mus10JOmvUqbCJS/Q5BzEdMd2Zi++jJodS02nEukVKioinVFfgb34UQBetC8m+ocLICDEcCjp1yZ+H/uWpWxyjCKcWgJeuxR2LjKdSqTHqaiIdMa613C669niSWXvpHsICdSS+dLzAmOHUPdf/2KJeyIBdgOe174P+3QmkPRvKioiHeXx0LziRQD+5j6LuTPTDQeSgWT6yBS+nPJHvnBn4vC4cC19znQkkR6loiLSUbsX41exmyo7mIK07zI2Kdx0IhlgfnbeBP4WfJX3zqZ/YdeVmw0k0oNUVEQ6yLXceyXbt9yncM2p4wynkYFoUKAfP77ycjZ70vG3Xax4+1nTkUR6jIqKSEcUZrWum/JlxEWcOjLOcCAZqKZmxFAxznvRwvjtr7FcC8FJP6WiItJeto3ng5/hwMP77pnMOXU2DodlOpUMYLMu+hENjmCGWvt44dVX2JivawFJ/6OiItJe6xfgyF1BnR3I8wHXc/HkZNOJZICzgsLxm/hfAPzMfoWbXvoPO4urDacS6V4qKiLt0VCJ/emDADzTfAmXnT6TIH+n4VAi4Hf6fXhCYhnn2MPtTX/h+/+7nNV79puOJdJtVFRE2mPTO1i1xWR7BvN+yCVcOTPNdCIRr/BEHJf9BRuLK/0Wc0bDQua+uJz31heYzdVQBa/9F3z2sNkc0uepqIi0g3vXYgD+7TmJm78zRqMp4luGnY512gMAPO7/Itfb73DPP1bw2ou/p2HR41CZ1/2vWVcOzU1H//pH98H2j+CrP0D57rZf87hhx2dQW9b9uaTfUVEROR6Ph6bt3qKyNWQa35+WajiQyBGcci9MuwEHNvf7L2Bt4I+4Mv+/CfrPo7j+NA33f54Ct6tzz523Cj64Byr2eu/nfAV/GAv/MxNqir3b6sph64dQV45n4zuw/rXWh1cufYn8inqKqxrYX9uE/dF98PfL4PkTIX9119639HuWbdt99nrhVVVVREREUFlZSXi4Ft2SnlG0bQWD/3E2NXYQn1+0gu9OyTAdSeTovnnRO5phuym2YtjnjmCiwzuiURo6Ci76H2JHTDv8ceteg09+AUNPgzN+CTHDvNvr98Oz06G2BMJT4MKn8bx1E44G7zyY4rCxvBt7E5fu/S3R7lLcWLhsP4IsFys9I5nu2E6JHcGsxmdoxo9LHV/yh4A/t76sywpgTcZNhA+bQfqEkwiJ0Cn/A0FH/v1WURE5jn89cy+Xlr3A6sCZTLn/EyxLpySLjytYB3WlNKTO5m8rcslb8lfucL9CtFWDy3ayNGg2DXETCUmbzKD0yQzJf4/oL3/Z+nDb4Ufx2OtZO+J20lc+wpi8Nw97iY2eDJKsUqKtmtZtVXYw4VY9ABs8Q7jGfojPnHcQSyW3ue+m2BPGq36PEmS5+HPzBYy08jjDua718bV2II+G/T/sYd/hvMRqTsh9Cb/BY2DSlRCe1GMfl/Q+FRWRbvLl9hL4v0s4xZlF0UkPM/isO01HEumw+iY3Hy1fT/LSXzKzcelR93ut+XQSrXJOd3ovdJjlyWCctQeHZXNr0x3c5vc2Yxy55NsxXO/3GJPDq3ik4gEC7CZyk85h36m/I85RQ0z5GgLHzCEwPN47mfarP8CgeO+oDDYNQ85k+xkvklteR1jWfML2LWdw7TaSKKbR9uNF9/lc61xIWEvp8eCgadwVBF3yJ/ALhIZKyF0JQ08Fp39vfIQD0/4c7yjbzB/DkNnd+tQqKiLdoLLexUVPfcbHjVcTZLngJysgfrTpWCKdZ9uUb1pE8cbPsfdtILZmG3Fu7xyTp5sv5enmy7Asi/P81vC4838IpQ6A9dHn8M2keQwJdTG+5APCJl3EoMEth4YKN0JVAYw4C4402lieDX+adPD+hCvg/Cch6Fs/s5ubaHj9hwTteL910yrPSDxYzHBsA2BD8Azqpt/KzHW/wKrKgzEXwuWvgEOT27tF6Q4o3gL+IZC7Apb9CZobICETfvTlkf//dpKKikhX7c/hb28sYO/eHP6f/z/whCbg+NnWbv2LKuIT6sqhfj929FCaPTZ+Dst7eLNsF7x1IzRWwQ8/gtD4zr/GokcgbyWceh9knHz0/dzN8N4dsP4f2Cf8hO0T7uGzbWUUr/2Q+yt/Q7B1hLOMpv4QLnhKfzc7Km8V9uJHaRyURGV0JsHZnxG+55PDdmtKPQn3nMcJTpnQrS+voiLda+VLsOplmH4DTLnG+9tLY433v/7BUFsKm96Gkm0w++6Dx5Jtu2/+8HA3U/uHyQyq3Xtw28S5cMmfj/4Ykf7K4wFHL58g2lQHASFtNuWtW0Tcuz8g0FPH++6ZLPFM4nd+L+CwbKoiRhOaMg7H6PNg/GV98+dOd/B4oCgLsMAvCPwCKW+Ewh2r8eQspaLWxULrBKjI5YH6JwimbfHz2BZZ9hAceGjCn780n8dHnhmcOz6R538wtVujduTfb79ufWXp8xqb3ewtq2N3aS3ZpbU48ldyw/Z7cOKB9+9k9/u/x7ZhCAU4LJv9dihh1OFneQAoXf0OzyU/zuSAPZyZ/7/YESlYFzxFSGqm4XfWfnlf/Z2U2r1U28HUhaYzeJATpt9kOpaIGb1dUuCwkgKQMuk7kLYMd1k24Z7x1K7M5cGtLh5xvkx45Vao3Aqb3qJhyycEXfzHIz5Hv9ZUh+f/LsWR+3WbzdEttwNO5pXWP3/hzmSrncpUv90UWIn83f8i9jpSsQCPDbVNzTgamxkUaLYqaERlANtf20RWfiVZ+ZVszK9kU0EV5fvLmWFtpsIOZYedwnsBvyDDUcQaz3CGWQVEWHVHfK4NniGEUs9QRyHNtqO1uAC4bCevBVzG1hE3M314IicMjSEpMtg7tLzoYYhMgzMeBL+A3nrrR1VcWUf10zMZZu/ljfBrueSnT+Pv1HJDIr6ooKKe979Yzrb1yxjRtIWbnO/jtGwKAodQe9ojDJ95AZaJotVLPB6bbUXVLNteyJTltzO5fjkNtj+VDCIQV+utyBHHrpCJRAc0M6riK5y2i4qxP6BpzuNEhYYc82ecbds0e+xu/zmoQz/SVmMNNVs+Y4NzLGtLHWzMr6QgN5uU6vWMdOSSbJXhxE0EtZzo2OSdOAo0W/742S5qghL56sx3iA3xI7VoEVZYAu6kyTj9AvGrKcDlDKbUP5Ga/cWMXHQDUeXrabICeSf0+8TVbOV0+xsAdniSebj5GmrsYM4O3cmNza8TYDd6M6afDFe8CoNiTH1KVDW4+NNzf+SX1b+hlmCa79hIRHSssTwi0j6NzW4+zNrH6iXv8tOKecRZVQBs9BtH1cQbmXLWXIKCgg2n7B51Tc0s3ZpH1fJXGVvwNgGeemwshjsKaLD9udX5K6z0E5mcFsWk1Egyk8MJCz7kl8D6/VCRCwkTjB4iU1EZqGwbKvZQte1LCouL2d4YSX3hTs4s+ztRVFFtB/Oyew7JVinfdXyNv+U+8tNEpmM1VHhPAcSC694/9gS4QzXVwdb3IW0WRKaCbVOz9l/4L/w5gQ2lh+3+jWcUYx17CaWe2kFp+F33LoFxQzr/GXTS/tomfvjSMh4qvZtJjl1UTrmViO8+2us5RKRrNm3fSenHj3FC2TsEtvzSVU44q9JvZsx37yY1ZpDhhB20fSG1We9Rm5uFszqPhmabMOoI/9botgcHBXNeIPmE7/WJtZ5UVLpTdaH39KzI9IPt00cmiTY2u9lTVsfO4hoqti7h9G3/TWLzka/pUWsHMshqbLPNPTgTZ9JEiB7qXZvAGeAtGIPHgacZcr/xTpZNntL1sPX74dMHYdO/8QSGUe0Xw2dBZ/No4Qyi67P5q//vSXWUkGfH8dfhz3DStCnMHhFHgF/PD9vmV9Rz518Xcdf+33KiczMeZyCOO7MgbHCPv7aI9Iz9hXvI/vAp0va+QyzeVXTfc89iadqPuDTDxaS0KAJGnOGTpzZ7PDZZeftp+ORhZubPP+I+FQEJlGfeSOqYE/BvqoC40RA7oldzdkWfKyrPPfccv//97yksLGTixIk888wzzJgx47iP66mi8p+sHSz790ucZ3/JBPcmABocIVT7xxLqriCwuZaakGRqwkfQGDUCd8worPjRBCaMInRQKKFNxfhVF0B4onfJ6aMdI/W4j/qXxLZt6hqbqSrcSWllDSVVDbgLNxJaspbiRn/mN53B+oog/GwXNzo/5G6/N/GzPLhsJxvsodQFxJDhLCfYHyrGX0fcydcRsWchrPqrd+GlWbd2TwHpIrfHZvWe/Sxbu4HLsn5Mqr2PPDuWX7uuZU3AdM4an8QFmUmcOCwGvx6YK/L5+t18+fb/cr3nLdIcJXj8B+H43ksw6txufy0R6X3uZhc73nuS4et/hx9tR5ELg4ZRe/IDDJ11MZbJheNsm/rSHLat+ZKSnWspKisjzZXNKc4sABa4T2df1AxSh49h+pBY0mJCsQZPAGffPR+mTxWV119/nWuuuYY///nPzJw5k6effpo333yTbdu2ER9/7PP2e6qoZL3+MBO2/AHwnq7lwkmg1Xzcx7ltiwYC2oxcNOFPhRVJjSOUSmc0xX6J4HAypjGLFFc2zZY/NVYYu/yGszxgJjnuOEIbixju2s53HKtIssqP+FqNtj+rPCOZ7NhJSMvr7U66gKoz5jEqLZngAN/7LeF4PBX5NL10HkHVOQDkeuL4yjOenXYS2YGjSRo/m/MnpjFjSDRORxdGtLa8T8Oqv1Gcn010fQ6hVgMAzeGp+F31undESUT6l73Lcf3zJhzVBewhkVhPWevhkxpCyAufhP/Q2aROOZuA5Ek9WwJqy2j22GSVO9iz+mPGbn2GkU2bD9utGT/WTnqYIWfeRGxoYM/lMaBPFZWZM2cyffp0nn32WQA8Hg+pqancfvvt3H///cd8bE8VlcqiHAJen0th+oXkJJxLKRE4yndi15RQYoez3+VPcE0uMXW7iG/MIalpD2mevUTiveZFs+2giCjiqCDgKPNA2qvJ9qPRCsRhQXlAEsWRE0lt3El8xbqDO4UmeC8iNvkHPnFIqktqSmDZH7HX/J93nswhKuxBfOGZyMqAGYSNPYsZY0cwY2hMh06d27fkLyQu+VmbbWWBKUScfBN+066D4MiuvwcR8U22DR43bsvJN5t3Uvv5k0wve5cIq7bNbo1WIBWhw/AkTiFi6mWEjDi1S4eIbNumqKKWvPWLiNzwF4aW/wcHNnV2YOsvmi7bSbYjnZroccTHDyYxJgrn2PMhaXKX3rKv6jNFpampiZCQEP75z39y8cUXt26/9tprqaio4N///vcxH+9Tk2ltG2pLaKoppyYkhWqXRXVdI43leVBXiqO+HL/afQRV78Vqrqc4chIlEZn4O2xCXWUklC1n8L7FBDRX4w5PwYoeinP0OQSOOB3LP/jw19qzFIo2QdoJMHiCmbUOelJTHexYCMWb8RRtxr37S/ybKtvs0mw7qCWYfP80yiMzsQePIzB+GKHRCQRYHpxNNbhKdtJclk1etU1xeQVXNizAYdm80Xwqu2NO4+LTZjJ60ol9v+CJSKc0NDaRtXopxRs/I2zfciZ5Nh82UbXciqQkKAPXoETs8GQckan4RSbhDE/EERKFy+XC1dSA29WEq6mBxtIcKN2JVVOI1bCfsIZChts5rWdUHqoZP3amfY/A0+8lI2NYn5gI2x36TFEpKCggOTmZZcuWMWvWrNbtP//5z/niiy9YsWJFm/0bGxtpbDx4WKWqqorU1FTfKCrSs9zNkLcS97aPqd/0IaGV2zv9VEvCv4t97hOcNjp+wPxQEJHj83hsNuVXsHXLeqqz1xBbvJRTm5cddf2ojqolhPVRZ1M45jpGjhzLmJBKnIOiYdDAWwah365MO2/ePB5++GHTMcQEpx+kz8KZPovQsx/2LuHfWE1hYQFF27+hOXcVIVXZRDbmE+apogk/GgikyJnI/sBkYoItEv1rGTRiNqed/rP+NwIlIl3mcFhMSI1iQuppwGnA3ZRUVLMm6z9UFe6mqTyXgNoCBjUUEtFcRpSnnFBqacYPN364LT+aLT+q/GOpCsnAE55CYHgsYbHJJI6awaD44ZzY5mePzixsjz516EcjKiIiIn1fR0ZUjP5aGRAQwNSpU1m0aFHrNo/Hw6JFi9ocCjogMDCQ8PDwNjcRERHpv4wf+rn77ru59tprmTZtGjNmzODpp5+mtraWH/7wh6ajiYiIiGHGi8r3v/99SkpKePDBByksLGTSpEl8/PHHDB6sY3ciIiIDnfF1VLrCp05PFhERkXbpM3NURERERI5FRUVERER8loqKiIiI+CwVFREREfFZKioiIiLis1RURERExGepqIiIiIjPUlERERERn6WiIiIiIj5LRUVERER8lvFr/XTFgdX/q6qqDCcRERGR9jrw73Z7ruLTp4tKdXU1AKmpqYaTiIiISEdVV1cTERFxzH369EUJPR4PBQUFhIWFYVmW6TjGVVVVkZqaSm5uri7S2IP0OfcOfc69Q59z79FnfZBt21RXV5OUlITDcexZKH16RMXhcJCSkmI6hs8JDw8f8H8JeoM+596hz7l36HPuPfqsvY43knKAJtOKiIiIz1JREREREZ+lotKPBAYG8utf/5rAwEDTUfo1fc69Q59z79Dn3Hv0WXdOn55MKyIiIv2bRlRERETEZ6moiIiIiM9SURERERGfpaLSx5WXl3PVVVcRHh5OZGQkN9xwAzU1Ne16rG3bnHvuuViWxTvvvNOzQfu4jn7O5eXl3H777YwaNYrg4GDS0tK44447qKys7MXUvu+5554jIyODoKAgZs6cyTfffHPM/d98801Gjx5NUFAQEyZM4MMPP+ylpH1bRz7nF198kdmzZxMVFUVUVBRnnnnmcf+/iFdHv58PWLBgAZZlcfHFF/dswD5KRaWPu+qqq9i0aROffvop77//Pl9++SU333xzux779NNPa0Xfduro51xQUEBBQQFPPPEEGzduZP78+Xz88cfccMMNvZjat73++uvcfffd/PrXv2bNmjVMnDiROXPmUFxcfMT9ly1bxty5c7nhhhtYu3YtF198MRdffDEbN27s5eR9S0c/5yVLljB37lwWL17M119/TWpqKmeffTb5+fm9nLxv6ejnfEBOTg733HMPs2fP7qWkfZAtfdbmzZttwF65cmXrto8++si2LMvOz88/5mPXrl1rJycn2/v27bMB++233+7htH1XVz7nQ73xxht2QECA7XK5eiJmnzNjxgz71ltvbb3vdrvtpKQke968eUfc/4orrrDPP//8Nttmzpxp/+hHP+rRnH1dRz/nb2tubrbDwsLsV155paci9gud+Zybm5vtE0880f7LX/5iX3vttfZFF13UC0n7Ho2o9GFff/01kZGRTJs2rXXbmWeeicPhYMWKFUd9XF1dHVdeeSXPPfccCQkJvRG1T+vs5/xtlZWVhIeH4+fXp69c0S2amppYvXo1Z555Zus2h8PBmWeeyddff33Ex3z99ddt9geYM2fOUfeXzn3O31ZXV4fL5SI6OrqnYvZ5nf2cH3nkEeLj4zXSehz6idmHFRYWEh8f32abn58f0dHRFBYWHvVxd911FyeeeCIXXXRRT0fsFzr7OR+qtLSU3/zmN+0+LNfflZaW4na7GTx4cJvtgwcPZuvWrUd8TGFh4RH3b+//g4GoM5/zt913330kJSUdVhLloM58zl999RUvvfQS69at64WEfZtGVHzQ/fffj2VZx7y194fMt7377rt8/vnnPP30090bug/qyc/5UFVVVZx//vmMHTuWhx56qOvBRXrJY489xoIFC3j77bcJCgoyHaffqK6u5uqrr+bFF18kNjbWdByfpxEVH/Szn/2M66677pj7DB06lISEhMMmajU3N1NeXn7UQzqff/45u3btIjIyss32yy67jNmzZ7NkyZIuJO9bevJzPqC6uppzzjmHsLAw3n77bfz9/bsau1+IjY3F6XRSVFTUZntRUdFRP9OEhIQO7S+d+5wPeOKJJ3jsscf47LPPyMzM7MmYfV5HP+ddu3aRk5PDhRde2LrN4/EA3tHabdu2MWzYsJ4N3ZeYniQjnXdgkueqVatat33yySfHnOS5b98+Oysrq80NsP/4xz/au3fv7q3ofUpnPmfbtu3Kykr7hBNOsE899VS7tra2N6L2KTNmzLBvu+221vtut9tOTk4+5mTaCy64oM22WbNmaTLtcXT0c7Zt23788cft8PBw++uvv+6NiP1CRz7n+vr6w34OX3TRRfYZZ5xhZ2Vl2Y2Njb0Z3eepqPRx55xzjj158mR7xYoV9ldffWWPGDHCnjt3buvX8/Ly7FGjRtkrVqw46nOgs36Oq6Ofc2VlpT1z5kx7woQJ9s6dO+19+/a13pqbm029DZ+yYMECOzAw0J4/f769efNm++abb7YjIyPtwsJC27Zt++qrr7bvv//+1v2XLl1q+/n52U888YS9ZcsW+9e//rXt7+9vZ2VlmXoLfUJHP+fHHnvMDggIsP/5z3+2+b6trq429Rb6hI5+zt+ms36OTkWljysrK7Pnzp1rh4aG2uHh4fYPf/jDNj9QsrOzbcBevHjxUZ9DReX4Ovo5L1682AaOeMvOzjbzJnzQM888Y6elpdkBAQH2jBkz7OXLl7d+7dRTT7WvvfbaNvu/8cYb9siRI+2AgAB73Lhx9gcffNDLifumjnzO6enpR/y+/fWvf937wfuYjn4/H0pF5eh09WQRERHxWTrrR0RERHyWioqIiIj4LBUVERER8VkqKiIiIuKzVFRERETEZ6moiIiIiM9SURERERGfpaIiIiIiPktFRURERHyWioqIiIj4LBUVERER8VkqKiJiXG1tLddccw2hoaEkJiby5JNPctppp3HnnXeyZMkSLMs67HbdddeZji0ivUBFRUSMu/fee/niiy/497//zcKFC1myZAlr1qwB4MQTT2Tfvn2tt88//5ygoCBOOeUUw6lFpDfo6skiYlRNTQ0xMTH87W9/4/LLLwegvLyclJQUbr75Zp5++unWfcvKypgxYwbnnHMOzz33nKHEItKbNKIiIkbt2rWLpqYmZs6c2botOjqaUaNGtdnP5XJx2WWXkZ6ezh//+MfejikihviZDiAi0h633HILubm5fPPNN/j56UeXyEChERURMWrYsGH4+/uzYsWK1m379+9n+/btrff/8Ic/8MYbb/Dvf/+bmJgYEzFFxBD9WiIiRoWGhnLDDTdw7733EhMTQ3x8PL/4xS9wOLy/R3322Wf8/Oc/57nnniM2NpbCwkIAgoODiYiIMBldRHqBJtOKiHE1NTXccsst/Otf/yIsLIyf/exnfPDBB0yaNInIyEgefvjhwx5z7bXXMn/+/N4PKyK9SkVFRHzSaaedxqRJk9qc9SMiA4/mqIiIiIjPUlERERERn6VDPyIiIuKzNKIiIiIiPktFRURERHyWioqIiIj4LBUVERER8VkqKiIiIuKzVFRERETEZ6moiIiIiM9SURERERGfpaIiIiIiPuv/A4dXVBNRE0/hAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "intensity_noisy = intensity + np.sqrt(intensity) * np.random.normal(0, 50, intensity.shape)\n", "plt.plot(qzs, intensity, label='original')\n", "plt.plot(qzs, intensity_noisy, label='added noise')\n", "plt.ylabel('Intensity')\n", "plt.xlabel('qz')\n", "plt.legend()\n", "\n", "\n", "initial_params = {'heights': {'value': height, 'variation': 10E-3},\n", " 'langles': {'value': langle, 'variation': 10E-3},\n", " 'rangles': {'value': rangle, 'variation': 10E-3},\n", " 'y_start': {'value': y1, 'variation': 10E-3},\n", " 'bot_cd': {'value': bot_cd, 'variation': 10E-3},\n", " 'dwx': {'value': dwx, 'variation': 10E-5},\n", " 'dwz': {'value': dwz, 'variation': 10E-5},\n", " 'i0': {'value': i0, 'variation': 10E-5},\n", " 'bkg_cste': {'value': bkg, 'variation': 10E-5},\n", " 'overlay': {'value': overlay, 'variation': 10E-3},\n", " 'top_cd': {'value': top_cd, 'variation': 10E-3},\n", " 'n1': n1,\n", " 'n2': n2,\n", " }\n", "\n", "StrongCastle1 = StrongCastleSimulation(qys=qxs, qzs=qzs, initial_guess=initial_params)\n", "\n", "Fitter1 = Fitter(Simulation=StrongCastle1, exp_data=intensity_noisy)\n", "\n", "bestfit,fitness = Fitter1.cmaes(sigma=10, ngen=100, popsize=1000, mu=10, n_default=15, restarts=10, tolhistfun=10E-5, ftarget=10, restart_from_best=True, verbose=False)\n", "\n", "print(bestfit)\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "15 parameters\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "100%|██████████| 200/200 [00:19<00:00, 10.10it/s]\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ " mean std count min max lower_ci \\\n", "height1 10.032762 0.154237 99987 8.643961 14.159510 10.031506 \n", "langle1 1.623662 0.050105 99987 0.190842 3.043659 1.623254 \n", "langle2 1.519561 0.066288 99987 0.058702 2.991455 1.519021 \n", "langle3 1.565897 0.066249 99987 0.156458 3.109113 1.565357 \n", "rangle1 1.489218 0.065244 99987 0.503257 2.847647 1.488687 \n", "rangle2 1.699659 0.067377 99987 0.002289 2.783351 1.699110 \n", "rangle3 1.529950 0.067966 99987 0.387656 3.065067 1.529396 \n", "y_start 0.457200 0.089482 99987 0.001071 1.319225 0.456471 \n", "bot_cd 39.984972 0.140185 99987 37.437256 41.138759 39.983830 \n", "dwx 0.098333 0.000941 99987 0.070365 0.127836 0.098326 \n", "dwz 0.100564 0.000768 99987 0.066336 0.108642 0.100558 \n", "i0 10.000036 0.000861 99987 9.958191 10.012195 10.000029 \n", "bkg_cste 0.099234 0.001000 99987 0.083602 0.113223 0.099225 \n", "overlay 9.806060 0.116891 99987 8.054958 13.368229 9.805108 \n", "top_cd 20.008583 0.087218 99987 18.643542 21.149993 20.007873 \n", "\n", " upper_ci uncertainity \n", "height1 10.034019 0.001257 \n", "langle1 1.624070 0.000408 \n", "langle2 1.520101 0.000540 \n", "langle3 1.566436 0.000540 \n", "rangle1 1.489750 0.000532 \n", "rangle2 1.700208 0.000549 \n", "rangle3 1.530504 0.000554 \n", "y_start 0.457929 0.000729 \n", "bot_cd 39.986114 0.001142 \n", "dwx 0.098341 0.000008 \n", "dwz 0.100570 0.000006 \n", "i0 10.000043 0.000007 \n", "bkg_cste 0.099242 0.000008 \n", "overlay 9.807013 0.000952 \n", "top_cd 20.009294 0.000711 \n" ] } ], "source": [ "with np.errstate(divide='ignore', invalid='ignore'):\n", " stats = Fitter1.mcmc_bestfit_stats(N=15, sigma = 100., nsteps=200, nwalkers=1000)\n", "\n", "print(stats)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.8" } }, "nbformat": 4, "nbformat_minor": 2 }