{ "cells": [ { "cell_type": "markdown", "id": "b342fc18-ac2f-4fb2-ae02-d7bc4cb12fa8", "metadata": {}, "source": [ "# Uncertainty quantification of partial wave transmission coefficients\n", "\n", "This notebook extends the uncertainty-propagation workflow to angular reaction observables. It shows how the same model choices that affect transmission coefficients flow through to angle-dependent reaction outputs.\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "38e6ed69-1acd-4a24-9ad4-1beb7be4769c", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from IPython.display import Math, display\n", "from periodictable import elements" ] }, { "cell_type": "code", "execution_count": 2, "id": "edfb7c75-cb6f-4a40-b272-2046c81038c7", "metadata": {}, "outputs": [], "source": [ "from matplotlib import pyplot as plt" ] }, { "cell_type": "markdown", "id": "e9650684-15b4-42d3-bbd7-9102360898af", "metadata": {}, "source": [ "## Compare mass-model inputs to Koning-Delaroche" ] }, { "cell_type": "code", "execution_count": 3, "id": "1ab2056f-8d65-4bd6-9bc0-6e963df78615", "metadata": {}, "outputs": [], "source": [ "A, Z = (24, 12)" ] }, { "cell_type": "code", "execution_count": 4, "id": "97b4f1d1-cf76-407a-9c7f-c1ecf678d4ed", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle ^{24} \\rm{Mg}$" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "name_core = str(elements[Z].symbol)\n", "display(Math(f\"^{{{A}}} \\\\rm{{{name_core}}}\"))" ] }, { "cell_type": "markdown", "id": "28e11d77-551e-4b45-a63f-3224d2483d0e", "metadata": {}, "source": [ "## Compute transmission coefficients\n", "Let's use `jitr` to calculate UQ'ed transmission coefficients using KDUQ Fermi energies and those from a variety of mass models" ] }, { "cell_type": "code", "execution_count": 5, "id": "08b94977-1a6c-4b0b-a14d-4e6485896117", "metadata": {}, "outputs": [], "source": [ "from tqdm import tqdm\n", "\n", "import jitr" ] }, { "cell_type": "code", "execution_count": 6, "id": "da0c0c5e-0447-49ef-b5ad-80e3476bebdc", "metadata": {}, "outputs": [], "source": [ "neutron = (1, 0)\n", "proton = (1, 1)\n", "projectile = proton\n", "target = (A, Z)" ] }, { "cell_type": "code", "execution_count": 7, "id": "a52e68e5-4102-446f-8b07-5459edbc7218", "metadata": {}, "outputs": [], "source": [ "# we have 416 samples from the KDUQ posterior\n", "kduq_omp_samples = jitr.optical_potentials.kduq.get_samples(proton)" ] }, { "cell_type": "code", "execution_count": 8, "id": "e1cab469-859c-4eaf-868b-23bde1896196", "metadata": {}, "outputs": [], "source": [ "# com_energy_grid = np.logspace(-1, 1.3, 100)\n", "lab_energy_grid = np.array([65, 200])\n", "range_fm = 15\n", "lmax = 20" ] }, { "cell_type": "code", "execution_count": 9, "id": "34e5d61b-bf31-4353-867b-fd91babc04f9", "metadata": {}, "outputs": [], "source": [ "reaction = jitr.reactions.Reaction(target=target, projectile=projectile, process=\"EL\")" ] }, { "cell_type": "code", "execution_count": 10, "id": "2252119f-5d0e-44c4-bcbc-adce1c344236", "metadata": {}, "outputs": [], "source": [ "def set_up_grid(core, lab_energy_grid):\n", " solvers = []\n", " for _i, Elab in enumerate(tqdm(lab_energy_grid)):\n", " kinematics = reaction.kinematics(Elab)\n", " a = range_fm * kinematics.k + np.pi / 2\n", " N = jitr.utils.suggested_basis_size(a)\n", " solvers.append(\n", " jitr.xs.elastic.IntegralWorkspace(\n", " reaction=reaction,\n", " kinematics=kinematics,\n", " channel_radius_fm=a / kinematics.k,\n", " solver=jitr.rmatrix.Solver(N),\n", " lmax=lmax,\n", " smatrix_abs_tol=0,\n", " )\n", " )\n", " return solvers" ] }, { "cell_type": "code", "execution_count": 11, "id": "8f3777be-61d4-4f80-b403-970e246cdb4a", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|███████████████████████████████████████████████████████████████████████████████████| 2/2 [00:14<00:00, 7.29s/it]\n" ] } ], "source": [ "solvers = set_up_grid(target, lab_energy_grid)" ] }, { "cell_type": "markdown", "id": "d2a16a52-51cf-4897-845e-f07e7bef205f", "metadata": {}, "source": [ "## Run the uncertainty propagation" ] }, { "cell_type": "markdown", "id": "355ff9ba-863e-4387-b2fc-54d16e18015c", "metadata": {}, "source": [ "### KDUQ" ] }, { "cell_type": "code", "execution_count": 12, "id": "7447cd48-7623-4d4c-91e3-b58ea3a9e9d5", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ " 97%|████████████████████████████████████████████████████████████████████████████▌ | 403/416 [08:09<00:23, 1.81s/it]/home/kyle/umich/jitr/src/jitr/optical_potentials/kduq.py:467: RuntimeWarning: overflow encountered in exp\n", " d2 = d2_0 + d2_A / (1 + np.exp((A - d2_A3) / d2_A2))\n", "100%|███████████████████████████████████████████████████████████████████████████████| 416/416 [08:26<00:00, 1.22s/it]\n" ] } ], "source": [ "tcoeff_kduq = np.zeros((lab_energy_grid.size, kduq_omp_samples.shape[0], 2, lmax + 1))\n", "for j, sample in enumerate(tqdm(kduq_omp_samples)):\n", " for i, _Elab in enumerate(lab_energy_grid):\n", " rgrid = solvers[i].radial_grid()\n", " central_params, spin_orbit_params, coulomb_params = (\n", " jitr.optical_potentials.kduq.calculate_params(\n", " projectile,\n", " target,\n", " solvers[i].kinematics.Elab,\n", " *sample,\n", " )\n", " )\n", "\n", " tplus, tminus = solvers[i].transmission_coefficients(\n", " jitr.optical_potentials.kduq.central(rgrid, *central_params),\n", " jitr.optical_potentials.kduq.spin_orbit(rgrid, *spin_orbit_params),\n", " jitr.optical_potentials.kduq.coulomb_charged_sphere(rgrid, *coulomb_params),\n", " )\n", " tcoeff_kduq[i, j, 0, :] = tplus\n", " tcoeff_kduq[i, j, 1, :] = tminus" ] }, { "cell_type": "code", "execution_count": 13, "id": "f155b4aa-bed3-423d-91be-cf3eaee5be39", "metadata": {}, "outputs": [], "source": [ "data = [np.zeros(lmax + 1)] * 8\n", "names = [\n", " \"E=65 MeV, j = l + 1/2\",\n", " \"E=65 MeV, j = l - 1/2\",\n", " \"E=65 MeV, j = l + 1/2, err\",\n", " \"E=65 MeV, j = l - 1/2, err\",\n", " \"E=200 MeV, j = l + 1/2\",\n", " \"E=200 MeV, j = l - 1/2\",\n", " \"E=200 MeV, j = l + 1/2, err\",\n", " \"E=200 MeV, j = l - 1/2, err\",\n", "]" ] }, { "cell_type": "code", "execution_count": 14, "id": "adaa45cb-e322-409f-b70d-90118647f355", "metadata": {}, "outputs": [], "source": [ "from pandas import DataFrame as df" ] }, { "cell_type": "code", "execution_count": 15, "id": "b7d119a9-4ec6-4b6a-99ad-848510c4e8e4", "metadata": {}, "outputs": [], "source": [ "data = df.from_dict(dict(zip(names, data, strict=False)))" ] }, { "cell_type": "code", "execution_count": 16, "id": "af00fea2-869c-4fb8-ba1d-bbfe700f3e27", "metadata": {}, "outputs": [], "source": [ "plus_color = \"tab:blue\"\n", "minus_color = \"tab:orange\"" ] }, { "cell_type": "code", "execution_count": 17, "id": "5dd6e94d-5dd3-4e1d-9c22-723c2325a163", "metadata": {}, "outputs": [], "source": [ "ci_plus = np.percentile(tcoeff_kduq[0, :, 0, :], 50, axis=0)\n", "ci_plus_errs = np.percentile(tcoeff_kduq[0, :, 0, :], 84, axis=0) - np.percentile(\n", " tcoeff_kduq[0, :, 0, :], 16, axis=0\n", ")" ] }, { "cell_type": "code", "execution_count": 18, "id": "2f308f0c-10a0-4094-8973-fac4573e331e", "metadata": {}, "outputs": [], "source": [ "ci_minus = np.percentile(tcoeff_kduq[0, :, 1, :], 50, axis=0)\n", "ci_minus_errs = np.percentile(tcoeff_kduq[0, :, 1, :], 84, axis=0) - np.percentile(\n", " tcoeff_kduq[0, :, 1, :], 16, axis=0\n", ")" ] }, { "cell_type": "code", "execution_count": 19, "id": "bde08fb2-72ec-4592-abf8-e244587fe381", "metadata": {}, "outputs": [], "source": [ "ci_minus[0] = None\n", "ci_minus_errs[0] = None" ] }, { "cell_type": "code", "execution_count": 20, "id": "cce3bfea-52bd-4ccd-a15d-76ec8090d3cb", "metadata": {}, "outputs": [], "source": [ "data[\"E=65 MeV, j = l + 1/2\"] = ci_plus\n", "data[\"E=65 MeV, j = l + 1/2, err\"] = ci_plus_errs\n", "data[\"E=65 MeV, j = l - 1/2\"] = ci_minus\n", "data[\"E=65 MeV, j = l - 1/2, err\"] = ci_minus_errs" ] }, { "cell_type": "code", "execution_count": 21, "id": "de35e3e0-dd65-457b-b9bc-b20fe3fc551f", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAk4AAAGGCAYAAACNCg6xAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAASR9JREFUeJzt3Xl4E+X6PvB7knRJt9SWFhqWsggilEKVRWWXTUU2F5CDyOJPBQqyieCC4kYBBTmiB4TvURAPuCMILiACZRHZbMtSSzmspUBP1aZ0bzPv74/QQGjaTtokk7b357pyXcnM5M1DQ5q7z7wzIwkhBIiIiIioUhq1CyAiIiKqKRiciIiIiBRicCIiIiJSiMGJiIiISCEGJyIiIiKFGJyIiIiIFGJwIiIiIlKIwYmIiIhIIQYnIiIiIoUYnIiIiIgUYnAiIiIiUojBichDxcXFoWPHjggMDET9+vUxfPhwnD171u62EydOhCRJeP/9991bJBFRHcPgROShdu3ahSlTpuC3337Djz/+iL/++gv3338/SkpKbLbbvHkzfv31VxiNRpUqJSKqOxiciDzUjz/+iDFjxqBNmzaIiYnBqlWr8Mcff+DEiRPWba5cuYKJEydi7dq18PLyUrFa+4QQuPPOO9G/f3+1SwEApKSkQKfT4V//+pfapRBRDcXgRFRDmEwmAEBISIh12bhx4/Dss8+iXbt25T7v7NmzkCQJkiShQYMGZTpWpZKTk63bNW3a1Ck1f/LJJzhy5Ahef/11p4xXXbfddhtGjhyJ1157DVevXq3SGBs2bEC/fv0QGhoKX19fNGvWDCNHjsSFCxfKbNu0aVPrz/TmW69evRS9nrvfv3/84x+QJAnr16+vcLvs7Gz4+fkhODgY+fn5VX49oppGp3YBRFQ5s9mM5557Dg888AAaNWoEAHj//feRm5uLmTNnKhpDp9PhypUr+P777zF48OAy6//9739Do3He31KyLGPevHno3r077rrrLqeNW13PP/88Pv30U7z33nt46aWXFD9PCIEJEyZg5cqVaNGiBR577DEEBgYiPT0du3btwrlz59C4ceMyzzMYDJg2bVqZ5Y6GG3e9f08++STWr1+Pjz76CCNHjix3u/Xr1yM/Px9jxoyBXq+v9usS1RiCiDyaLMvi//2//ydatmwpMjIyhBBCJCcni/DwcHHmzBnrdpGRkWLZsmVlnn/mzBkBQPTo0UMYDAYxZMiQMtsUFxeL+vXri/79+wsfHx8RGRlZ7bo3b94sAIhVq1ZVeyxni46OFpGRkcJsNit+ztKlSwUAMWnSJFFSUlJmfXFxcZllkZGR1f5Zuvv9k2VZNGvWTGg0GnHu3Llyt+vcubMAIA4ePFjl1yKqibirjsiDCSEwadIk/Pzzz9i+fTvCwsIAAPv378f//vc/3HrrrdDpdNDpdDh37hymTp2KDh062B1Lr9fjsccew5YtW5CRkWGzbvPmzbhy5QrGjx9v97klJSWIi4tDixYt4Ovri1tvvRVxcXE4ffo0JEnC2LFjyzzn448/hiRJePjhh22W79y5E5IkYd68edizZw969eqFwMBABAcH4+GHH8apU6cU/3yqOtbw4cNx7tw57NixQ9Hr5Ofn47XXXkPz5s3xz3/+E1qttsw2Op1rG/jVef9KxcfHY9CgQahXrx58fHzQsmVLvPzyy8jLy7NuI0kSxo0bB1mW8fHHH9sd5/jx4zhw4ACio6PRsWPH6v/jiGoQBiciDyWEQGxsLLZs2YJffvnFZjfQ0KFDkZSUhISEBOvNaDRizpw5+Oqrr8odc/z48SgpKcHatWttln/00UcICQnB0KFDy33eiy++CACIjY3Ffffdh3fffdfuLqjS2nfs2IHbbrsNt9xyi91t9u/fjz59+sBgMGDKlCno2bMnNmzYgHvuuQenT5+u4CdT/bHuvvtuAMD27dsVjb9161b8/fffGDp0KMxmM7755hssWLAAK1asqDToFRYWYvXq1Zg/fz7ef/99/Pbbbw79225U1fcPAJYvX45evXph7969GDhwIJ599lk0atQIb731Fvr164eioiLrtmPHjoVGo8Hq1ashhCgzVmmgevLJJ6v8byGqsdRteBFReSZOnCiCg4PF7t27xaVLl6y3wsJCu9tXtqtuwIABQgghoqKiRNu2ba3rL126JHQ6nZgyZYoQQpTZ1fPzzz8LAKJDhw4iNzfXujw9PV3Ur19fABBjxoyxec3jx48LAGLUqFFl6tmxY4cAIACIFStW2KxbsWKFACAefPDBin841RzLZDJZd38pMXfuXAFAzJo1S7Rq1cr6mgCERqMRM2fOtPu8yMhIm21Lb506dRKnTp1S9NrVff+EsLwfOp1OtG/fXmRmZtqsi4uLEwDEO++8Y7P8vvvuEwDEzz//bLO8dLegj4+P+PPPPxX9G4hqE3aciDzU8uXLkZWVhe7duyMiIsJ627dvX7XGHT9+PI4fP27tfKxZswYlJSXl7ub59NNPAQCvvPIK/Pz8rMsjIiIwdepUu89JS0sDANSvX7/cOlq1aoWnnnrKZtlTTz2Fli1bYsuWLfjf//6n+N/k6FhBQUHw9fW11lmZ0l1jS5YsgcFgwIEDB3D16lXEx8ejVatWWLx4MZYvX17meePGjcP27dtx5coV5Obm4vfff8fo0aNx8OBB9OnTp0pH9jn6/gHAhx9+iJKSEixbtgyhoaE2655//nmEhYWVOYqutJv00Ucf2Swv3S04ZMgQmyM8ieoKHlVH5KGEnV0kFSnvrOI3e/zxxzF79mx89NFH6NKlCz7++GPExMSUOzcqMTERANCtW7cy67p27Wr3OX/++ScAIDg4uNw6unbtWuYoMI1Gg65duyI1NRWJiYno27evgn9R1cYKCQlBZmamovFlWQYAeHt749tvv7WebLR79+748ssv0b59eyxevBgTJ060ed6rr75q87hDhw745JNPAABr167FqlWrMGPGDEU1lHL0/QMsuzIB4KeffrK7e9LLywt//PGHzbIhQ4YgLCwMGzZsgMlkgsFgAHA9SHE3HdVVDE5EdUxYWBgGDRqEzz77DI8++ihSUlKwbNmycrfPzs6GRqNBvXr1yqwrr6NUenh6QUFBueOW99zS5aXnrVKiKmPl5+fbdNAqUhoaOnbsWOYM7VFRUWjevDlOnTqFrKysCsNiqWeeeQZr167F3r17HQ5Ojr5/APDXX38BAN566y3Fr+Pl5YXRo0djyZIlWLduHSZOnIjLly/jhx9+QJMmTRSHWqLahrvqiOqgJ598EtnZ2Rg7dix8fX0xatSocrcNCgqCLMt2uzNXrlyx+5zSo/9Kv7DtKe+5pctLw4oSjo4lyzJMJpO1zsrcdtttAMrvoJUuV3oiyNIQmpubq2j7mzny/gGW9xCwhGAhRLk3e68DWM4RBVi6ZCUlJRg3bpxTz/lFVJPwfz5RHTRgwAA0bNgQFy9exNChQ8s98g0A2rdvDwDYu3dvmXXlzbdq27YtNBoNUlJSyh1379691l1gpWRZxr59+yBJkvV1lXB0rNTUVMiyXOEZ12/Uu3dvAJazc9+suLgYp06dgr+/v+IgVjo/qapn+Hbk/QOALl26ALi+y06pNm3a4K677sLhw4eRlJRkPcXEuHHjqlQ3UW3A4ERUB2m1Wnz77bfYsGED4uLiKty2tJvx+uuv23RULl++jH/+8592nxMcHIzo6GgcOnSoTKApdfLkSaxatcpm2apVq3Dy5EkMHDiwTAjp1asXJEnCzp07qz1WaXDp2bOn3dpu1qJFC/Tv3x+nTp3C//3f/9msW7BgAbKysjBs2DCbczn98ccfNudHunH57NmzAVgub1IVjrx/ADBp0iTodDpMmTIF58+fL7M+KysLv//+u93nlnadJk2ahOTkZPTt2xeRkZFVqpuoVlDvgD4icoebD2evjL3D2f/xj38IAKJFixZi5syZYvLkySI8PFwMGjRIABDjxo0rM85rr70mAIi9e/faLC89hcCAAQOEt7e3GDx4sHjhhRfE4MGDhSRJol69euK///1vmfG6d+8uAIjdu3dXe6zHH39c6HQ6cfHiRUU/EyGEOHXqlAgPDxcAxMCBA8XMmTPFvffeKwCIyMhIcenSJZvtX331VREYGCgGDhwoJk2aJGbNmiWGDBkivLy8BADxwgsvKHpdZ7x/QgixcuVKodVqha+vr3jooYfErFmzxIQJE6xnG3/mmWfsjpednS38/f2tp1L47LPPFNVBVFsxOBHVcs744i0uLhZvvPGGaNasmfD29hbNmzcX8+fPF7/99psAIKZOnVpmnIsXLwqdTicmTpxos7w07Lz66qti9+7domfPnsLf318EBQWJYcOGidTU1DJjybIsQkJCRNOmTW0ubVKVsXJzc0VAQIAYOnSoop/Hjc6fPy/Gjh0rGjRoILy8vETjxo1FbGysuHLlSpltd+7cKYYPHy5atmwpgoKChE6nEw0aNBBDhgwRP/30k+LXdFZwEkKIAwcOiMcee0wYjUbh5eUl6tWrJ+644w4xZ84ckZycXO6Y48aNEwBESEiIKCgoUFw7UW3E4EREVbZq1SoBQPzrX/+yu/7xxx8Xt9xyi8jOzrYuuzHsKHX06FEBQHzwwQc2y6syVmnNu3btUvwcIqJSnONERJW6fPlymaOuLl68iDfffBNarRYPPvig3ee9+eabyM/Pr/Rw+crs3r0b9evXr/RabJUpKSnB/PnzMXjwYPTo0aNaYxFR3cTzOBFRpRYsWIAtW7age/fuCA8Px/nz57F582ZcvXoV8+bNs7mO3o0iIyOxZs2ack8XoNTEiRPLnFyyKs6fP48nnngCo0ePrvZYRFQ3MTgRUaXuu+8+nDhxAlu2bMHff/8NX19fREdHY9KkSZUeGTZ8+HA3VVm55s2bY968eWqXQUQ1mCRu7r8TERERkV2c40RERESkEIMTERERkUIMTkREREQK1brJ4bIsIz09HYGBgZAkSe1yiIiIyMMJIXD16lUYjcZKL2Bd64JTenp6uYdGExEREZXnwoULaNSoUYXb1LrgFBgYCMDyjw8KClK5GiIiIvJ02dnZaNy4sTVDVKTWBafS3XNBQUEMTkRERKSYkik+nBxOREREpBCDExEREZFCDE5ERERECtW6OU5ERETOZjabUVxcrHYZVEVeXl7QarVOGYvBiYiIqBxCCFy+fBlZWVlql0LVFBwcjAYNGlT7HI8MTkREROUoDU3h4eHw8/PjiZVrICEE8vLykJGRAQCIiIio1ngMTkRERHaYzWZraAoNDVW7HKoGvV4PAMjIyEB4eHi1dttxcjgREZEdpXOa/Pz8VK6EnKH0fazuXDUGJyIiogpw91zt4Kz3kcGJiIiISCEGJyIiIiKFGJyI3CSvqARN52xB0zlbkFdUonY5RORG/PzXHgxOREREtdisWbMwaNAgp4753HPPYejQoU4ds6ZgcCIiIqrFXnrpJaxbt86pYyYkJKBDhw5OGy8+Ph6DBg2C0WiEJEn49ttvy9123LhxePnllxEXF4dOnTohMDAQ4eHhGDp0KFJSUpxWU3kYnIiIiNxo33//dOvrBQcHIzAw0KljJiYmKgpOvXr1wurVqyvdLjc3F+3bt8cHH3xQ4XZmsxmbN2/G4MGDsWvXLsTGxmL//v3Ytm0biouL0b9/f+Tm5ir8V1QNT4BJRETkYkII6/2l206iT+twt5zmIDMzE2FhYTh69CiioqKcMmZaWhoyMzPRvn17p4wHAPfffz/uv//+Srfbt28fvLy80KlTJ/z4448261avXo3w8HAcPnwYPXr0cFptN2PHiYiIyMX2nrreZTqWno341Ey3vG5iYiJ8fHzQunXrMuvmz5+PgICACm/nz58v87yEhAQYDAY0a9bMHf8EG5s2bcKgQYPshk6TyQQACAkJcWkN7DgRERG5kBAC721PtT7WSMDirSno0bKey7tOiYmJaNu2LXS6sl/3EyZMwPDhwyt8vtFoLLMsISHBqd0mR2zcuBHvvvtumeWyLGPatGno2rWr0zpr5WFwIiIicqH41EwcS8+2PpYFkJRmQnxqJnq2CnPpa1c0iTskJKRK3ZmKgtP8+fMxf/586+P8/Hzs378fkydPti47ceIEmjRp4vDrJicnIz09HX369CmzLjY2FseOHcOePXscHtdR3FVHVKooF5hnsNyKXDu5kIjqBiEEFm9NgeamxlJp1+nGuU+ukJiYWGHIqequuvLC2IQJE5CQkGC9dezYEa+//rrNMntdLCU2bdqEfv36wdfX12b55MmTsXnzZuzYsQONGjWq0tiOYMeJiIjIReJTM5GUZiqz3B1dp6KiIiQnJ5cbnKqyq+7q1as4ffq04i6WXq9HeHg4br31VseKt2Pjxo14+umnrY+FEJgyZQo2bNiAnTt3um3OFYMTERGRC5R2myQJsNdYklw81yk5ORnFxcXlBqeq7KpLTEyEVqtF27ZtnVGiVU5ODk6dOmV9fObMGSQkJCAkJARNmjRBRkYGDh06hE2bNlm3iY2Nxbp167Bx40YEBgbi8uXLAACDwQC9Xu/U+m7EXXVEREQuUGSWkZ6Vbzc0AZYwdSmrAEVm2SWvn5CQgMjISAQHBzt1zNatW8PHx8dpYwLAoUOHEBMTg5iYGADAjBkzEBMTg1deeQUA8N1336Fz586oV6+e9TnLly+HyWRCr169EBERYb19/vnnTq3tZuw4ERERuYCPTotNk7vhr9wiFBSb8ciKXwEAX024G75eWgBAaIA3fHRal7z+wYMH0blzZ6eOOXnyZJuJ3pXZuXOnou169epV4XyvjRs3YvDgwTbLXD0/rDwMTkTX5BWVwO/G+96qlkNEtYAxWA9jsN7mwr5tjEHw83bd129BQQGOHz+Or7/+Gm+++abLXsedunXrhpEjR6pdBgAGJyIiIpfz89bh7IKBbnmtpUuXYuHChRgxYgSeeOIJt7ymqz3//PNql2DF4ETkLkW5OOv7DwBAXtF5wNugckFEVBvNmTMHc+bMUbuMWouTw4mIiIgUYnAiIiIiUsjjgpPZbMbcuXPRrFkz6PV6tGjRAm+88YZqs+eJiIiISnncHKeFCxdi+fLlWLNmDdq2bYtDhw5h3LhxMBgMePbZZ9Uuj4iIiOowjwtO+/btw5AhQzBwoOXog6ZNm2L9+vU4cOCAypURERFRXedxu+ruuecebN++HSdPngRgOb37nj17cP/996tcGREREdV1HtdxmjNnDrKzs9G6dWtotVqYzWa89dZbGDVqlN3tCwsLUVhYaH2cnZ3trlKJiIiojvG4jtMXX3yB//znP1i3bh2OHDmCNWvW4J133sGaNWvsbh8XFweDwWC9NW7c2M0VExERUV3hccFp1qxZmDNnDh577DG0a9cOo0ePxvTp0xEXF2d3+xdeeAEmk8l6u3DhgpsrJiIiorrC44JTXl4eNBrbsrRaLWTZ/tWjfXx8EBQUZHMjIiLyKEW5wDyD5VaU69aXnjVrFgYNGuTUMZ977jkMHTrUqWPWFB43x2nQoEF466230KRJE7Rt2xa///47lixZgvHjx6tdGhERUY3z0ksvQavVOnXMhIQEdOvWzWnjxcfH4+2338bhw4dx6dIlbNiwodxgNm7cODRs2BD+/v745ptv8Mcff0Cv1+Oee+7BwoULcdtttzmtLns8ruO0bNkyPPLII5g0aRJuv/12PPfcc3jmmWfwxhtvqF0aERGRY0xpQHoCcDnp+rLLSZZl6QmA6aLLSwgODkZgYKBTx0xMTESHDh2cNl5ubi7at2+PDz74oMLtzGYzNm/ejMGDB2PXrl2IjY3F/v37sW3bNhQXF6N///7IzXVtR8/jOk6BgYFYunQpli5dqnYpREREVVdSCKzsDeRm2C7/6L7r9wPCgWnHAJ2PS0rIzMxEWFgYjh49iqioKKeMmZaWhszMTLRv394p4wHA/fffr+i0Q/v27YOXlxc6deqEH3/80Wbd6tWrER4ejsOHD6NHjx5Oq+1mHtdxIqrN9pij0LdwEfadzlK7FCJyNa03YGiI8r9qNUBQQ8t2LpKYmAgfHx+0bt26zLr58+cjICCgwtv58+fLPC8hIQEGgwHNmjVzWd3l2bRpEwYNGgRJksqsM5lMAICQkBCX1uBxHSei2koIgUUlI3BKNMLSHWfRp10Tux9+IqolJAm492Xg04fL2UC2rHfh74HExES0bdsWOl3Zr/sJEyZg+PDhFT7faDSWWZaQkODUbpMjNm7ciHfffbfMclmWMW3aNHTt2tVpnbXyMDgRmdKA3ExI+TnWRdKVo0B2gOWBf9i1vxqrZ+/pLCSJFgCAY5dyEJ+aiZ6twqo9LhF5sBZ9AGMMcCkJEObryyUtEBFtWe9CCQkJ5c5FCgkJqVJ3pqLgNGfOHCxcuLDC5ycnJ9vtgFUmOTkZ6enp6NOn7M8sNjYWx44dw549exwe11EMTlS33TAHQX/DYv3agdcfOGEOghAC7+08Bw3MkKGFRgIWb01Bj5b12HUiqs3K6zoJs8u7TYCl4/Tkk0/aXTd//nzMnz+/wuefOHECTZo0sVmWkJCABx54wO72M2fOxNixYyscs3nz5hWuL8+mTZvQr18/+Pr62iyfPHkyNm/ejPj4eDRq1KhKYzuCwYnqttI5CLmZAOydK8w5cxDiUzNx7FIOAMshwbIAktJM7DoR1QUt+gAR7YFLiZbHbuo2FRUVITk5udzuUFV21V29ehWnT58ut4sVFhaGsDDX/E7buHEjnn76aetjIQSmTJmCDRs2YOfOnW6bc8XgRHWbG+YgCCGweGsKNJIlMJVi14mojpAkoOds4LN/WB67qduUnJyM4uLicoNTVXbVJSYmQqvVom3bts4o0SonJwenTp2yPj5z5gwSEhIQEhKCJk2aICMjA4cOHcKmTZus28TGxmLdunXYuHEjAgMDcfnyZQCAwWCAXq8v8xrOwqPqiK7NQRDSTSeIk7SWuQnV/KswPjUTSWkmm9AE2HadiKiWa9bz+v2I9i7vNgGWXWqRkZEIDg526pitW7eGj49zT59w6NAhxMTEICYmBgAwY8YMxMTE4JVXXgEAfPfdd+jcuTPq1atnfc7y5cthMpnQq1cvREREWG+ff/65U2u7GTtORNe6TpIL5iCUdpskCRCi7HqJXSeiuuHGz3evF1zebQKAgwcPonPnzk4dc/LkyZg8ebJTxwSAXr16Qdj7JXnNxo0bMXjwYJtlFW3vSuw4EQFAiz4wBbdFibB8JEqEBtkh7ar9V2GRWUZ6Vr7d0ARYwtSlrAIUme1fi5GIaqFmrjs5IwAUFBTg8OHD+PrrrzFgwACXvpa7dOvWDSNHjlS7DADsOBEBAASAxSWP4nVpHgBAJ8lYXPIo5gGozt+FPjotNk3uhr9yi1CQlwPftZYjUQpGfw9fP8vpDkIDvOGjc+51pIio7lq6dCkWLlyIESNG4IknnlC7HKd4/vnn1S7BisGJCJZ5SJ9ktsTD3s3RXnMaiXJzrMlogXudcNSbMVgPY7AeeTmAn+YsACAvIgB+AQYnVE5ENYK3PzDP5JaXmjNnDubMmeOW16qLuKuO6rzrR71JWFQyAqlyQywqGQGNJGHx1hTV9qMTEZHnYXCiOu/Go972yu3Qr+ht7JXb8ag3IiIqg8GJ6rQbj3qzp/SoN3adiIgIYHCiOo5HvRERkSM4OZzqNB71RkSVYce5dnDW+8jgRHUej3ojInu8vLwAAHl5eS69hAe5R15eHoDr72tVMTgRERHZodVqERwcjIyMDACAn58fz/BfAwkhkJeXh4yMDAQHB0Orrd4eBAYnqjmKcoH5167U/WK65bwoREQu1KBBAwCwhiequYKDg63vZ3UwOBEREZVDkiREREQgPDwcxcXFapdDVeTl5VXtTlMpBiciVzOlAbmZkPJzrIukK0eBbMvkc/iHAYaGKhVHREpotVqnffFSzcbgRORKJYXAyt5AbgZunFqqXzvw+oOAcGDaMUDn4/byiIjIMTyPE5Erab2vdZPK+6hpgKCGlu2IiMjjMTgRuZIkAfe+DKC8E2jKlvU8UoeIqEZgcCJytRZ9AGMMhHTT/AhJCxhjLOuJiKhGYHAicrVrXSdJmG2XCzO7TURENQyDE5E7tOgDc/1o60PBbhMRUY3E4ETkDpKE4m6zrj9kt4mIqEZicCJyE7lpD+t9c/1odpuIiGogBicid7mhu1TcYza7TURENRCDE5EK5MjuapdARERVwOBEREREpBCDE9UYeUUldu8TERG5C4MTERERkUK8yC95PlMakJsJKT/Huki6chTIDrA88A+7dj04IiIi12JwIs9WUgis7A3kZkB/w2L92oHXHwSEA9OOATqfar2UX4ABmGey3K/WSEREVFtxVx15Nq33tW5Sef9VNUBQQ8t2RERELsbgRJ7t2nXeALmcDWSegZuIiNyGwYk8X4s+gDHGcn23G9Ww6735eevs3iciopqDwYk837WukyTMtst5vTciInIzBieqGVr0sVzf7RpRw7pNRERUOzA4Uc0gSSjuNuv6Q3abiIhIBQxOVGPITXtY75vrR7PbREREbsfgRDXHDd2l4h6z2W26SV5RCZrO2YKmc7bwkjRERC7C4EQ1khzZXe0SiIioDmJwIiIiIlKIwYmIiIhIIQYnIiIiIoUYnIiIiIgUYnAiIiIiUojBiYiIiEghjwxOFy9exOOPP47Q0FDo9Xq0a9cOhw4dUrssIiIiquM87hLtf//9N7p27YrevXvjhx9+QFhYGFJTU3HLLbeoXRoRERHVcR4XnBYuXIjGjRvj448/ti5r1qyZihURERERWXjcrrpNmzahY8eOePTRRxEeHo6YmBisWrVK7bKIiIiIPC84nT59GsuXL0fLli3x008/YeLEiXj22WexZs0au9sXFhYiOzvb5ka11x5zFPoWLsK+01lql0JERHWQxwUnWZZxxx13YP78+YiJicHTTz+Np556CitWrLC7fVxcHAwGg/XWuHFjN1dM7iKEwKKSETglGmHpjrMQQqhdEhER1TEeF5wiIiLQpk0bm2W33347zp8/b3f7F154ASaTyXq7cOGCO8okFew9nYUk0QIAcOxSDuJTM1WuiIiI6hqPmxzetWtXpKSk2Cw7efIkIiMj7W7v4+MDHx8fd5RGKhJC4L2d56CBGTK00EjA4q0p6NGyHiRJUrs8z1CUi7O+/wAA5BWdB7wNKhdERFT7eFzHafr06di/fz/mz5+PU6dOYd26dVi5ciViY2PVLo1UFJ+aiWOXciBDCwCQBZCUZmLXiYiI3MrjglOnTp2wYcMGrF+/HlFRUXjjjTewdOlSjBo1Su3SSCVCCCzemgLNTY2l0q4T5zoREZG7eNyuOgB48MEH8eCDD6pdBnmI+NRMJKWZyiy/sevUs1WYCpUREVFd43EdJ6IblXabypvGJLHrREREbuSRHSeiUkVmGelZ+SgvFwkBXMoqQJFZho9O697iHOXtD8wr2zkjIqKag8GJnKcoF5hvtNx/Md0SFKrJR6fFpsnd8FduEQqKzXhkxa8AgK8m3A1fL0tQCg3w9vzQREREtQKDE3k8Y7AexmA98opKrMvaGIPg583/vkRE5F6c40RERESkEIMTERERkUIMTkREREQKMTgRERERKcTgRERERKQQgxNRLbLHHIW+hYuw73SW2qUQEdVKDE5EtYQQAotKRuCUaISlO87ybOpERC7A4EROc+N5lm68T+6x93QWkkQLAMCxSzmIT81UuSIiotqHwYmoFhBC4L2d56CBGQCg4TX8iIhcolrBKS4uDgBw5MgRFBYWOqUgInJcfGomjl3KgQzLpWdkASSlmdh1IiJysmpds6JXr14AgAULFuD48ePQaDRo27YtoqOjER0djU6dOqF+/frOqJOIyiGEwOKtKdBIlsBUqrTr1KNlPUiSpF6BRES1iEMdp5vb/nfffTcA4IsvvsDx48dx4MABzJw5E/Xr18e2bdvwwAMPYO7cuc6rlojKiE/NRFKaySY0Aew6ERG5gkMdp7fffhuxsbHw97dc9X769OnW7lJUVBT0ej06deqETp06WZ9z55134o033nBu1eRZTGlAbiak/BzrIunKUSA7wPLAPwwwNFSpuNqttNskSYC96UwSu05ERE7lUHDq1q0b4uLi8OabbwIAevfujaSkJPzwww84fvw4JEmy2VU3aNAg7N+/3yWFk4coKQRW9gZyM6C/YbF+7cDrDwLCgWnHAJ2P28ur7YrMMtKz8u2GJsASpi5lFaDILMNHp3VvcUREtZBDwalRo0b47LPPrMFp8ODBGDx4sHV9QUEBjh07hqSkJPzyyy8YNGgQvLy8nFsxeRatt6WblJsJQLazgQYIamjZjpzOR6fFpsnd8FduEQrycuC79gEAQMHo7+HrZ+n4hQZ4MzQRETmJQ8GpsLAQZ86cQXFxsd1A5Ovri44dO6Jjx45OK5A8nCQB974MfPpwORvIlvXcTeQyxmA9jMF65OUAfpqzAIC8iAD4BRjULYyIqBZyaHK4LMsQQiAmJgbPP/88duzYgZISnuiwzmvRBzDGQEg3dTUkLWCMsawnIiKqBRwKTiEhIWjUqBHatGmD1atXo0+fPggNDcWwYcPw73//G3///ber6iRPdq3rJAmz7XJhZrfJHUxpQHqCZUL+NdKVo0B6guVmuqhaaUREtY1DwSksLAwPP/wwvvjiC1y5cgW//vorpk6dirS0NDz99NNo3LgxvvjiC1fVSp6sRR+Y60dbHwp2m9yjdHL+yp42E/L1awcCK3tabqt6WbYjIqJqc/jM4aGhoSgsLIQkSejSpQtef/11HDx4EOnp6XjzzTd56oG6SpJQ3G3W9Ycu6Db5eetwdsFAnF0wEH7e1Tp3a+1ROjm/3I8yJ+cTETmTw8FpypQp+PDDD8ssr1+/PqZNm4YlS5Y4pTCqeeSmPaz3zfWj2W1yh9LJ+XaPaAQ4OZ+IyLkcDk4GgwEjR47E1atX7a7v169ftYuiGuqGL+fiHrP5Ze0unJxPROQ2VbrIb1hYGAIDA51dC9UicmR3tUuoOzg5n4jIbaoUnACguLgYFy5cQEpKCv766y9n1kREjuLkfCIit3AoOF29ehXLly9Hz549ERQUhKZNm+L2229HWFgYIiMj8dRTT+HgwYOuqpVqgD3mKPQtXIR9p7PULqVuccPkfCIiciA4LVmyBE2bNsXHH3+Mvn374ttvv0VCQgJOnjyJX3/9Fa+++ipKSkrQv39/3HfffUhNTXVl3eSBhBBYVDICp0QjLN1xFqK8C6iRS3ByPhGR6yk+pvvgwYOIj49H27Zt7a7v3Lkzxo8fj+XLl2P16tXYvXs3WrZs6bRCyfPtPZ2FJNECAHDsUg7iUzPRs1WYylXVITdNztey20RE5HSKg9P69esVbefr64sJEyZUuSCqmYQQeG/nOWhghgwtNBKweGsKerSsB4lf4G7HyflERK5R5cnh9ly4cAHjx4935pBUQ8SnZuLYpRzIsBwSLwsgKc2E+NRMlSsjIiJyHqcGp7/++gtr1qxx5pBUAwghsHhrCjQ3NZZKu06c60RERLWFQ9et2LRpU4XrT58+Xa1iqGaKT81EUpqpzPIbu06c60RERLWBQ8Fp6NChkCSpwg4C57PULaXdJkkC7P23kDjXiYiIahGHdtVFRETgm2++gSzLdm9HjhxxVZ3koYrMMtKz8u2GJsASpi5lFaDIXN611IiIiGoOhzpOd955Jw4fPowhQ4bYXV9ZN4pqHx+dFpsmd8NfuUUoyMuB79oHAAAFo7+Hr18AACA0wBs+Om1FwxAREdUIDgWnWbNmITc3t9z1t956K3bs2FHtoqhmMQbrYQzWIy8H8NOcBQDkRQTAL8CgbmFERERO5lBw6t694nPD+Pv7o2fPntUqiIiIiMhTORSciCrk7Y+mBesAACe8/VUuhoiIyPmceh4nIiIiotqMwYmIiIhIIQYnIiIiIoWqHZzS0tIgyzxHD5Ha/Lx1du8TEZHzVDs4tWnTBmfPnnVCKURERESerdrBiSe8JCIiorqCc5yIiIiIFHJ4IsQnn3xi87ikpATffPMNwsPDrcueeOKJ6ldGRERE5GEcDk4ff/yxzePi4mJ89dVX0Ov1ACzXq2NwIiIiotrI4eB087XoAgMDsW7dOjRv3txpRZGLFOUC842W+y+mAzy7NxERkUM8fo7TggULIEkSpk2bpnYpREREVMd5dHA6ePAgPvzwQ0RHR6tdChEREVH1g9OLL76IkJAQZ9RiIycnB6NGjcKqVatwyy23OH18IiIiIkdVOzi98MILCA4OdkIptmJjYzFw4ED07dvX6WPXVXlFJXbvO4uftw5nFwzE2QUDeeZqIiKqlTzy2+2zzz7DkSNHcPDgwUq3LSwsRGFhofVxdna2K0sjIiKiOszj5jhduHABU6dOxX/+8x/4+vpWun1cXBwMBoP11rhxYzdUSURERHWRJDzsminffvsthg0bBq1Wa11mNpshSRI0Gg0KCwtt1tnrODVu3BgmkwlBQUFurd3T5eWY4PdOE8v9587DL8CgckVERETqy87OhsFgUJQdFO+qe+WVV/DVV1/BaDSia9euGDZsGDp06GCzzffff48uXbogNDS0SoUDQJ8+fXD06FGbZePGjUPr1q0xe/Zsm9AEAD4+PvDx8any6xEREREppXhXXXJyMoYNG4aZM2eioKAAo0aNQseOHfH1119btwkNDa32WcMDAwMRFRVlc/P390doaCiioqKqNTYRERFRdSgOTkuXLkVqairGjh2L//73vxg7diz69++PV155BXfccQd27NiB4uJinDlzxpX1UlWY0oD0BEhXrnfypCtHgfQEy810UbXSiIiIahKH5zhlZ2fjl19+we7du3Ho0CGcPXsWV65cgZeXF/Ly8rBmzRo8/vjjrqpXUX1K91PWCSWFwLtRQG5G+dsEhAPTjgE67vIkIqK6xyVznEoFBQVh6NChGDp0aJl1a9aswddff61qcKKbaL0BQ0MgNxOAbGcDDRDU0LIdERERVcippyMYM2YMPvvsM2cOSdUlScC9L8N+aIJl+b0vW7YjIiKiCjkcnPLz87Fnzx6cOHGizLqCggJ8/vnnTimMnKhFH8AYAyHZHpEISQsYYyzriYiIqFIOBaeTJ0/i9ttvR48ePdCuXTv07NkTly5dsq43mUwYN26c04ukarrWdZKE2Xa5MLPbRERE5ACHgtPs2bMRFRWFjIwMpKSkIDAwEF27dsX58+ddVR85S4s+MNePRumhALKkYbeJiIjIQQ4Fp3379iEuLg716tXDrbfeiu+++w4DBgxA9+7dcfr0aVfVSM4gSSju9py1uaQRMkRvdpuIiIgc4VBwys/Ph053/UA8SZKwfPlyDBo0CD179sTJkyedXiA5T7w5GolycwBAotwc8XK0yhURERHVLA6djqB169Y4dOgQbr/9dpvl77//PgBg8ODBzquMnEoIgfd2ncctJY/iFd2neKdkBEzbTqJHqzBI7DoREREp4lDHadiwYVi/fr3dde+//z5GjhwJD7tmMF0Tn5qJY5dysFtuj35Fb2O33A5JaSbEp2aqXRoREVGN4fCZwz0dzxxelhACQz7Yi2MXTZBveLc1EhDV0ICNsV3ZdaJK5RWVoM0rPwEATrw+AH7eDp8/l4jIIzmSHRR3nBw9cu7iRV7/zFPEp2YiKc02NAGALMCuExERkQMUB6dOnTrhmWeewcGDB8vdxmQyYdWqVYiKisLXX3/tlAKpeoQQWLw1pdyD5yQJWLw1hbtYiYiIFFDcaz9x4gTeeust9OvXD76+vrjzzjthNBrh6+uLv//+GydOnMDx48dxxx13YNGiRXjggQdcWTcpVGSWkZ6Vj/JykRDApawCFJll+Oi09jciAoCiXJz1/QcAIK/oPOBtULkgIiL3c3iOU35+PrZs2YI9e/bg3LlzyM/PR7169RATE4MBAwYgKirKVbUqwjlOZaVn5eOv3CIU5OXAd60l0BaM/h6+fgEAgNAAb0QY9GqWSDVAXo4Jfu80sdx/7jz8AhiciKh2cCQ7ODy7U6/X45FHHsEjjzxS5QLJvYzBehiD9cjLAfw0ZwEAeREB/OIjIiJykMMX+SUiIiKqqxiciIiIiBTiiVjqEm9/NC1YBwA44e2vcjFEREQ1DztORERERAoxOBEREREpxOBEREREpBCDExEREZFCDE5ERERECvGoujrEz1uHswsGql0GERFRjcWOExEptscchb6Fi7DvdJbapRARqYIdJyKqmCkNyM2ElHcVX5m7wwfF+G7bdvSpFwMJEuAfBhgaql0lEZFbMDgRUflKCoGVvYHcDOgBLPW+tvwqgJXX7geEA9OOATofdWokInIj7qojovJpvQFDQwhIdlcLSEBQQ8t2RER1AIMTEZVPkoB7X4YEYX81BHDvy5btiIjqAAYnIqqQaH4vTmpbokTY/rooERqc1LaEaH6vSpUREbkfgxMRVSj+1J94M/8h6CTZZrlOkvFm/kOIP/WnSpUREbkfgxMRlUsIgcVbU7BbRCNRbm7tOpUIDRLl5tgtorF4awqEsL8rj4iotmFwIqJyFZllpGflQwgJi0setXaddJKMxSWPQggJl7IKUGSWKxmJiKh2kEQt+1MxOzsbBoMBJpMJQUFBapdDVOOlZ+Xjr9wiFOReRZv/xMBPKkJuSBTOPLQFkCSEBngjwqBXu0wioipzJDvwPE5EVCFjsB7GYD3ycgA/qQgAoO3zAqIaBatbGBGRCrirjogcJkd2V7sEIiJVMDgRERERKcTgRERERKQQgxMRERGRQgxORERERAoxOHmSolxgnsFyK8pVuxoiIiK6CYOTB8krKrF7n4iIiDwDgxMRERGRQgxORERERArxzOGewJQG5GZCys+xLpKuHAWyAywP/MMAQ0OViiMiIqJSDE5qKykEVvYGcjNw49W+9GsHXn8QEA5MOwbofNxeHhEREV3HXXVq03pf6yaV91ZogKCGlu2IiIhIVQxOapMk4N6XAcjlbCBb1kuSO6siKsPPW2f3PhFRXcLffp6gRR/AGANxKQmSMFsXC0kLKSLasp5Ibd7+wDyT2lUQEamKHSdPcK3rdGNoAmB5zG4TERGRx2Bw8hCi+b04qb0VJcLylpQIDU5qW0I0v1flyoiIiKiUxwWnuLg4dOrUCYGBgQgPD8fQoUORkpKidlkuF3/qT7yZ/zB0kmWuk06S8Wb+Q4g/9afKlREREVEpjwtOu3btQmxsLPbv349t27ahuLgY/fv3R25u7b12mxACi7emYI+IRqLcHACQKDfHHhGNxVtTIIRQuUIiIiICPHBy+I8//mjzePXq1QgPD8fhw4fRo0cPlapyrfjUTCSlmQBIWFQyAvN0n2BRyQjIQkJSmgnxqZno2SpM7TKJiIjqPI8LTjczmSxH8YSEhNhdX1hYiMLCQuvj7Oxst9TlLKXdJkkChAD2yu3Qr+ht63pJAhZvTUGPlvUgcZI4ERGRqjxuV92NZFnGtGnT0LVrV0RFRdndJi4uDgaDwXpr3Lixm6usniKzjPSsfJS3N04I4FJWAYrM5Z3niYiIiNxFEh48gWbixIn44YcfsGfPHjRq1MjuNvY6To0bN4bJZEJQUJC7Sq2W9Kx8/JVbhIK8HPiufQAAUDD6e/j6Wa5VFxrgjQiDvqIhiIiIqIqys7NhMBgUZQeP3VU3efJkbN68GfHx8eWGJgDw8fGBj0/NvoabMVgPY7AeeTmAn+YsACAvIgB+AQZ1CyMiIiIbHhechBCYMmUKNmzYgJ07d6JZs2Zql0REREQEwAODU2xsLNatW4eNGzciMDAQly9fBgAYDAbo9bV8d5W3P5oWrAMAnPD2V7kYIiIiupnHTQ5fvnw5TCYTevXqhYiICOvt888/V7s0IiIiquM8ruPkwXPVXc7PW4ezCwaqXQYRERGVw+M6TkRERESeisGJiIiISCEGJyIiIiKFGJyIyCPk5ZiAeQZgnsFyn4jIAzE4ERERESnE4ERERESkEIMTERERkUIedx4nIqpjTGlAbiak/BzrIunKUSDbcpFr+IcBhoYqFUdEZIvBiYjUU1IIrOwN5Gbgxgsq6dfecCLYgHBg2jFAV7Mv5k1EtQN31RGRerTe17pJ5f0q0gBBDS3bERF5AAYnIlKPJAH3vgxALmcD2bJektxZFRFRuRiciEhdLfoAxhjIku2vIyFpAGOMZT0RkYdgcCIidUkSRO+XoRG2XSdJyBC92W0iIs/C4EREqouXo5EoN0eJsPxKKhEaJMrNES9Hq1wZEZEtBiciUpUQAou3ncSSkkehkyxdJ50kY0nJo1i87SSEECpXSER0HYMTEakqPjUTSWkm7LrWdQKARLk5dsnRSEozIT41U+UKiYiuY3AiItUIIbB4a8q1aUwSFpWMQKrcEItKRgCQIEnA4q0p7DoRkcfgCTCJSDVFZhnpWfkozUV75XboV/S2db0QwKWsAhSZZfjotCpVSUR0HYMTEanGR6fFpsnd8FduEQrycuC79gEAQMHo7+HrZ7nkSmiAN0MTEXkMBiciUpUxWA9jsB55OYCf5iwAIC8iAH4BBnULIyKyg3OciIiIiBRicCIiIiJSiLvqiMgzePujacE6AMAJb3+ViyEiso/BiYg8gp+3DmcXDFS7DCKiCnFXHREREZFCDE5ERERECjE4ERERESnE4ERERESkEIMTERERkUIMTg7IKypB0zlb0HTOFvycfEXtcoiIiMjNeDoCB9x4hfal206iT+twSJbLuhORpzKlAbmZQEk+8NF9lmXjfwR0est9/zDA0FC9+oioRmFwcsDeU39a7x9Lz0Z8aiZ6tgpTsSIiqlBJIbCyN5CbYbu8NEABQEA4MO0YoPNxb21EVCNxV51CQgi8tz3V+lgjAYu3pth0oYjIw2i9r3WTyvtVpwGCGlq2IyJSgMFJofjUTBxLz7Y+lgWQlGZCfGqmilURUYUkCbj3ZQByORvIlvXc5U5ECjE4KSCEwOKtKdDc9LuVXSeiGqBFH8AYA1my/XUnJA1gjLGsJyJSiMFJgfjUTCSlmSDflI/YdSKqASQJovfL0AjbrpMkZIje7DYRkWMYnCpR2m0q73erxK4TkceLl6ORKDdHibD8yisRGiTKzREvR6tcGRHVNAxOlSgyy0jPykd5uUgI4FJWAYrM5c2hICI1CSGweNtJLCl5FDrJ8jnVSTKWlDyKxdtO8o8eInIIT0dQCR+dFpsmd8NfuUXlbhMa4A0fndaNVRGRUqW72gFL16m95jQS5ebYJUcD13a187QiRKQUg5MCxmA9jMF6tcsgIgfduKtdCAmLSkZgnu4TLCoZAUCy7mrv0bIeT2ZLRIowOBFRrXXzrva9cjv0K3rbuv7GXe3sGhOREgxORFRruWVXe+klXcrDS7oQ1SoMTkRUq7l0V3t5l3S5ES/pQlSr8Kg6IqKq4iVdiOocBicioqriJV2I6hwGJyKi6uAlXYjqFM5xIiKqjtJLuvznYdvF1y7pUu3THHDyOZFHYXAiIqqmeDkawXJztJXOQifJKBEaHBdNkSVHo2d1BubkcyKPw+BERFQNpZd0uaXkUazxXgjg2iVdih/F39tOokersKp3nUonn+dmwv48KidMPmdHi8ghDE5ERNXg0ku6lE4+//Thcjao5uRzV3e0GMqoFmJwIiKqIrdc0qVFHwhjDOT0RGhv6DqZoYHG2B5SdSafu7Kj5Y7djK4OZhy/do9fRR4bnD744AO8/fbbuHz5Mtq3b49ly5ahc+fOapdFRGTllku6SBKOtpqC6PTxNou1kJHUagqiqzP53JUdLVfvZnR1MOP4tXv8avDI4PT5559jxowZWLFiBbp06YKlS5diwIABSElJQXh4uNrlEREBcM8lXYQQePloON6wM/l87tFwbOwpqnfknqs6Wq7ezejqYMbxa/f41eCR53FasmQJnnrqKYwbNw5t2rTBihUr4Ofnh48++kjt0oiIbBiD9YhqaCj3FmGo3uVe4lMzkXQxG4tLHoVOsnyB6CQZi0seRdLFbMSnVrArQ4lrHS3tTV9OWsg42mpK9U7eeS2UmW/6qjFDA1Hdc1y5+uSjHL92j18NHhecioqKcPjwYfTt29e6TKPRoG/fvvj111/LbF9YWIjs7GybGxFRbXDjHKp42TL5HAAS5eaIl6Otc6hE6b7CKr7Gy0fDkSg3R4mwfCWUCA0S5eZ4+Wh4tcZ2aSgDXBvMOH7tH7+KPC44ZWZmwmw2o379+jbL69evj8uXL5fZPi4uDgaDwXpr3Lixu0olInIp2zlUlsnnqXJD6+TzG+dQVZUrO1ouDWWA64MZx6/d41eRxwUnR73wwgswmUzW24ULF9QuiYjIKUrnUG2eYrm9EDsRhc/8ihdiJ1qXbZrStcpzqFzd0XL1bkZXBzOOX7vHryqPC0716tWDVqvFlStXbJZfuXIFDRo0KLO9j48PgoKCbG5ERLWFK+dQubKj5Y7djK4OZhy/do9fVR4XnLy9vXHnnXdi+/bt1mWyLGP79u24++67VayMiKh2cWVHy9W7GV0dzDh+7R6/OjzydAQzZszAmDFj0LFjR3Tu3BlLly5Fbm4uxo0bp3ZpRES1ijFYD2Nw9Y78s6fsqRq6oRAT8cIN21TnVA32gtmNJx+t7jm0OH7tHr86JKHWTsJKvP/++9YTYHbo0AHvvfceunTpUunzsrOzYTAYYDKZuNuOiKgWS8/Kr/QcWtXZlcnxa/f4N3IkO3hscKoqBiciIiJyhCPZwePmOBERERF5KgYnIiIiIoUYnIiIiIgUYnAiIiIiUojBiYiIiEghBiciIiIihTzyBJjVUXp2hezsbJUrISIiopqgNDMoOUNTrQtOV69eBQA0btxY5UqIiIioJrl69SoMBkOF29S6E2DKsoz09HQEBgZCkiSnj5+dnY3GjRvjwoULNfIEm6xfXaxfXaxfPTW5doD1q83V9QshcPXqVRiNRmg0Fc9iqnUdJ41Gg0aNGrn8dYKCgmrkf75SrF9drF9drF89Nbl2gPWrzZX1V9ZpKsXJ4UREREQKMTgRERERKcTg5CAfHx+8+uqr8PHxUbuUKmH96mL96mL96qnJtQOsX22eVH+tmxxORERE5CrsOBEREREpxOBEREREpBCDExEREZFCDE4O+uCDD9C0aVP4+vqiS5cuOHDggNolKRIXF4dOnTohMDAQ4eHhGDp0KFJSUtQuq0oWLFgASZIwbdo0tUtxyMWLF/H4448jNDQUer0e7dq1w6FDh9Quq1Jmsxlz585Fs2bNoNfr0aJFC7zxxhuKLk2ghvj4eAwaNAhGoxGSJOHbb7+1WS+EwCuvvIKIiAjo9Xr07dsXqamp6hRrR0X1FxcXY/bs2WjXrh38/f1hNBrxxBNPID09Xb2Cb1LZz/9GEyZMgCRJWLp0qdvqq4yS+pOTkzF48GAYDAb4+/ujU6dOOH/+vPuLtaOy+nNycjB58mQ0atQIer0ebdq0wYoVK9Qp9iZKvqcKCgoQGxuL0NBQBAQE4OGHH8aVK1fcWieDkwM+//xzzJgxA6+++iqOHDmC9u3bY8CAAcjIyFC7tErt2rULsbGx2L9/P7Zt24bi4mL0798fubm5apfmkIMHD+LDDz9EdHS02qU45O+//0bXrl3h5eWFH374ASdOnMDixYtxyy23qF1apRYuXIjly5fj/fffR3JyMhYuXIhFixZh2bJlapdmV25uLtq3b48PPvjA7vpFixbhvffew4oVK/Dbb7/B398fAwYMQEFBgZsrta+i+vPy8nDkyBHMnTsXR44cwTfffIOUlBQMHjxYhUrtq+znX2rDhg3Yv38/jEajmypTprL6//vf/6Jbt25o3bo1du7ciaSkJMydOxe+vr5urtS+yuqfMWMGfvzxR3z66adITk7GtGnTMHnyZGzatMnNlZal5Htq+vTp+O677/Dll19i165dSE9Px0MPPeTeQgUp1rlzZxEbG2t9bDabhdFoFHFxcSpWVTUZGRkCgNi1a5fapSh29epV0bJlS7Ft2zbRs2dPMXXqVLVLUmz27NmiW7duapdRJQMHDhTjx4+3WfbQQw+JUaNGqVSRcgDEhg0brI9lWRYNGjQQb7/9tnVZVlaW8PHxEevXr1ehwordXL89Bw4cEADEuXPn3FOUA8qrPy0tTTRs2FAcO3ZMREZGinfffdfttSlhr/4RI0aIxx9/XJ2CHGSv/rZt24rXX3/dZtkdd9whXnrpJTdWpszN31NZWVnCy8tLfPnll9ZtkpOTBQDx66+/uq0udpwUKioqwuHDh9G3b1/rMo1Gg759++LXX39VsbKqMZlMAICQkBCVK1EuNjYWAwcOtHkPaopNmzahY8eOePTRRxEeHo6YmBisWrVK7bIUueeee7B9+3acPHkSAJCYmIg9e/bg/vvvV7kyx505cwaXL1+2+T9kMBjQpUuXGvk5BiyfZUmSEBwcrHYpisiyjNGjR2PWrFlo27at2uU4RJZlbNmyBa1atcKAAQMQHh6OLl26VLg70tPcc8892LRpEy5evAghBHbs2IGTJ0+if//+apdWxs3fU4cPH0ZxcbHN57d169Zo0qSJWz+/DE4KZWZmwmw2o379+jbL69evj8uXL6tUVdXIsoxp06aha9euiIqKUrscRT777DMcOXIEcXFxapdSJadPn8by5cvRsmVL/PTTT5g4cSKeffZZrFmzRu3SKjVnzhw89thjaN26Nby8vBATE4Np06Zh1KhRapfmsNLPam34HAOW+R6zZ8/GyJEja8z1xxYuXAidTodnn31W7VIclpGRgZycHCxYsAD33Xcftm7dimHDhuGhhx7Crl271C5PkWXLlqFNmzZo1KgRvL29cd999+GDDz5Ajx491C7Nhr3vqcuXL8Pb27vMHwnu/vzWuov8UuViY2Nx7Ngx7NmzR+1SFLlw4QKmTp2Kbdu2ecw8AkfJsoyOHTti/vz5AICYmBgcO3YMK1aswJgxY1SurmJffPEF/vOf/2DdunVo27YtEhISMG3aNBiNRo+vvTYrLi7G8OHDIYTA8uXL1S5HkcOHD+Of//wnjhw5AkmS1C7HYbIsAwCGDBmC6dOnAwA6dOiAffv2YcWKFejZs6ea5SmybNky7N+/H5s2bUJkZCTi4+MRGxsLo9HoUd18T/6eYsdJoXr16kGr1ZaZvX/lyhU0aNBApaocN3nyZGzevBk7duxAo0aN1C5HkcOHDyMjIwN33HEHdDoddDoddu3ahffeew86nQ5ms1ntEisVERGBNm3a2Cy7/fbbPeZInIrMmjXL2nVq164dRo8ejenTp9fI7l/pZ7Wmf45LQ9O5c+ewbdu2GtNt2r17NzIyMtCkSRPrZ/ncuXOYOXMmmjZtqnZ5lapXrx50Ol2N/Szn5+fjxRdfxJIlSzBo0CBER0dj8uTJGDFiBN555x21y7Mq73uqQYMGKCoqQlZWls327v78Mjgp5O3tjTvvvBPbt2+3LpNlGdu3b8fdd9+tYmXKCCEwefJkbNiwAb/88guaNWumdkmK9enTB0ePHkVCQoL11rFjR4waNQoJCQnQarVql1iprl27ljms9uTJk4iMjFSpIuXy8vKg0dj+qtBqtda/vmuSZs2aoUGDBjaf4+zsbPz222814nMMXA9Nqamp+PnnnxEaGqp2SYqNHj0aSUlJNp9lo9GIWbNm4aefflK7vEp5e3ujU6dONfazXFxcjOLiYo/9PFf2PXXnnXfCy8vL5vObkpKC8+fPu/Xzy111DpgxYwbGjBmDjh07onPnzli6dClyc3Mxbtw4tUurVGxsLNatW4eNGzciMDDQuj/YYDBAr9erXF3FAgMDy8zF8vf3R2hoaI2ZozV9+nTcc889mD9/PoYPH44DBw5g5cqVWLlypdqlVWrQoEF466230KRJE7Rt2xa///47lixZgvHjx6tdml05OTk4deqU9fGZM2eQkJCAkJAQNGnSBNOmTcObb76Jli1bolmzZpg7dy6MRiOGDh2qXtE3qKj+iIgIPPLIIzhy5Ag2b94Ms9ls/SyHhITA29tbrbKtKvv53xz0vLy80KBBA9x2223uLtWuyuqfNWsWRowYgR49eqB379748ccf8d1332Hnzp3qFX2Dyurv2bMnZs2aBb1ej8jISOzatQuffPIJlixZomLVFpV9TxkMBjz55JOYMWMGQkJCEBQUhClTpuDuu+/GXXfd5b5C3Xb8Xi2xbNky0aRJE+Ht7S06d+4s9u/fr3ZJigCwe/v444/VLq1KatrpCIQQ4rvvvhNRUVHCx8dHtG7dWqxcuVLtkhTJzs4WU6dOFU2aNBG+vr6iefPm4qWXXhKFhYVql2bXjh077P5fHzNmjBDCckqCuXPnivr16wsfHx/Rp08fkZKSom7RN6io/jNnzpT7Wd6xY4fapQshKv/538zTTkegpP5///vf4tZbbxW+vr6iffv24ttvv1Wv4JtUVv+lS5fE2LFjhdFoFL6+vuK2224TixcvFrIsq1u4UPY9lZ+fLyZNmiRuueUW4efnJ4YNGyYuXbrk1jqla8USERERUSU4x4mIiIhIIQYnIiIiIoUYnIiIiIgUYnAiIiIiUojBiYiIiEghBiciIiIihRiciIiIiBRicCIiIiJSiMGJiIiISCEGJyKqM9LS0vDII49g7ty5apdCRDUUgxMR1RnTp09Hy5Yt8eWXX6pdChHVUAxORFQnmEwm7Ny5E926dYPRaFS7HCKqoRiciKhO+OWXX9C9e3fs2rULXbt2VbscIqqhGJyIqE7YvXs37rrrLmzevBnDhg1TuxwiqqEYnIioTjh06BD+/vtv+Pr64o477lC7HCKqoXRqF0BE5A5nz57F//73P7z22mtql0JENZgkhBBqF0FE5Go+Pj5o3bo1fv/9d2g0bLYTUdXwtwcR1QleXl5YvHgxQxMRVQt/gxBRrbdmzRrk5ubCx8cH+/fvx+bNm9UuiYhqKM5xIqJaraCgAN988w3WrFmDp59+Gk2bNsUnn3yidllEVENxjhMRERGRQtxVR0RERKQQgxMRERGRQgxORERERAoxOBEREREpxOBEREREpBCDExEREZFCDE5ERERECjE4ERERESnE4ERERESkEIMTERERkUIMTkREREQKMTgRERERKfT/AYeMHoU/Y/80AAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig = plt.figure(figsize=(6, 4))\n", "ls = np.arange(lmax + 1)\n", "plt.errorbar(\n", " ls,\n", " ci_plus * (2 * ls + 1),\n", " ci_plus_errs * (2 * ls + 1),\n", " linestyle=\"none\",\n", " marker=\"^\",\n", " label=\"$j = l + 1/2$\",\n", ")\n", "\n", "plt.errorbar(\n", " ls[1:],\n", " ci_minus[1:] * (2 * ls[1:] + 1),\n", " ci_minus_errs[1:] * (2 * ls[1:] + 1),\n", " linestyle=\"none\",\n", " marker=\"v\",\n", " label=\"$j = l - 1/2$\",\n", ")\n", "\n", "# plt.yscale(\"log\")\n", "plt.xticks([0, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20])\n", "plt.xlabel(r\"$\\ell$\")\n", "plt.ylabel(r\"$(2 \\ell +1) \\mathcal{T}_{lj}$\")\n", "plt.title(r\"$^{24}$Mg(p,p) 65 MeV\", fontsize=14)\n", "plt.legend()\n", "plt.tight_layout()" ] }, { "cell_type": "code", "execution_count": 22, "id": "ee3f8958-cb7c-4da6-90a9-f69a7f235f16", "metadata": {}, "outputs": [], "source": [ "ci_plus = np.percentile(tcoeff_kduq[1, :, 0, :], 50, axis=0)\n", "ci_plus_errs = np.percentile(tcoeff_kduq[1, :, 0, :], 84, axis=0) - np.percentile(\n", " tcoeff_kduq[1, :, 0, :], 16, axis=0\n", ")" ] }, { "cell_type": "code", "execution_count": 23, "id": "3134f66d-8c3c-4e88-aefc-204f4fcdea1d", "metadata": {}, "outputs": [], "source": [ "ci_minus = np.percentile(tcoeff_kduq[1, :, 1, :], 50, axis=0)\n", "ci_minus_errs = np.percentile(tcoeff_kduq[1, :, 1, :], 84, axis=0) - np.percentile(\n", " tcoeff_kduq[1, :, 1, :], 16, axis=0\n", ")\n", "ci_minus[0] = None\n", "ci_minus_errs[0] = None" ] }, { "cell_type": "code", "execution_count": 24, "id": "4e4ce1ce-cf10-4016-88f2-3ff083f359c3", "metadata": {}, "outputs": [], "source": [ "data[\"E=200 MeV, j = l + 1/2\"] = ci_plus\n", "data[\"E=200 MeV, j = l + 1/2, err\"] = ci_plus_errs\n", "data[\"E=200 MeV, j = l - 1/2\"] = ci_minus\n", "data[\"E=200 MeV, j = l - 1/2, err\"] = ci_minus_errs" ] }, { "cell_type": "code", "execution_count": 25, "id": "81e608ad-d1c3-4c77-85c8-b4d1ead13306", "metadata": { "tags": [ "nbval-check-output" ] }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAk4AAAGGCAYAAACNCg6xAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAATlVJREFUeJzt3XlcVPX+P/DXgWETWQIBGRVFskxQ1FxKTfRCLpmGWZqZufSrq+GCdk2tXCoTtZt5NdPse11abLm3zKWuXTMU9xQvIGqKuQtIlIKOss7n98c4AyMzMANn5szA6/l4zOMx53POnHkLwrx5n895fyQhhAARERER1chF6QCIiIiInAUTJyIiIiILMXEiIiIishATJyIiIiILMXEiIiIishATJyIiIiILMXEiIiIishATJyIiIiILMXEiIiIishATJyIiIiILMXEiIiIishATJyIHlJSUhC5dusDHxwchISEYPnw4zp8/b/LYiRMnQpIkfPDBB/YNkoioAWLiROSAdu/ejcmTJ+PQoUPYvn07/vzzTwwcOBBlZWVGx23btg0HDhyAWq1WKFIiooaFiRORA9q+fTvGjBmDdu3aoVOnTvj444/x66+/4sSJE4Zjrl69iokTJ+LTTz+Fm5ubgtGaJoTAgw8+iH79+ikdCgDg1KlTUKlU+PDDD5UOhYicGBMnIidQUFAAAAgICDCMjRs3DlOmTEH79u3Nvu78+fOQJAmSJKFp06ZVKlZ6J0+eNBzXqlUrWWL+5JNPcPToUbz11luynK+u7r//fowcORJvvvkmbty4YfHrrly5gmXLlqFfv34ICwuDu7s7mjZtimHDhuHQoUPVvvbw4cN47LHH4O/vD29vbzz00EP4+uuvZTveFHt+z5999llIkoQvvvii2uMKCwvRqFEj+Pv74/bt27V6LyJHwcSJyMGVl5fjb3/7Gx577DE0b94cAPDBBx9Ao9HglVdesegcKpUKV69exQ8//GBy/z//+U+4uLjAxUWeXwlarRbz58/HI488goceekiWc8rh1VdfRV5eHpYvX27xa1asWIFp06bh7Nmz6NevH1555RX06tULmzdvRo8ePfDVV1+ZfF1ycjJ69uyJvXv3Yvjw4ZgwYQJyc3MxYsQIvPfee3U+vib2+J6/8MILAIC1a9dWe9wXX3yB27dvY+TIkfDy8qr1+xE5BEFEDkur1Yr/9//+n2jTpo3Iy8sTQghx8uRJERwcLM6dO2c4rmXLlmLFihVVXn/u3DkBQPTu3Vv4+fmJJ554osoxpaWlIiQkRPTr1094eHiIli1b1jnubdu2CQDi448/rvO55NahQwfRsmVLUV5ebtHx33zzjdi1a1eV8ZSUFOHm5ibuueceUVRUZLSvtLRURERECA8PD/G///3PMH79+nVx3333CXd3d3H+/PlaH18de37PtVqtCA8PFy4uLuLChQtmj+vWrZsAIA4fPlyr9yFyJKw4ETkoIQRefvll/PTTT9i5cyeCgoIAAAcPHsTvv/+Oe++9FyqVCiqVChcuXMDUqVPRsWNHk+fy8vLCM888g++//x55eXlG+7Zt24arV69i/PjxJl9bVlaGpKQkREREwNPTE/feey+SkpJw9uxZSJKEsWPHVnnNunXrIEkShg0bZjS+a9cuSJKE+fPnY+/evejTpw98fHzg7++PYcOG4cyZMxZ/fWp7ruHDh+PChQtITk626H2efPJJxMTEVBl/5JFH0LdvX1y7dg3Hjh0z2vfzzz/jt99+w7PPPmv0PfHz88Nrr72GkpISbNiwodbHW6Iu33MASElJweDBg9GkSRN4eHigTZs2eOONN3Dr1i3DMZIkYdy4cdBqtVi3bp3J8xw/fhy//PILOnTogC5dulj1byByREyciByQEAIJCQn4/vvv8fPPP6NFixaGffHx8cjIyEBaWprhoVarMWvWLPz73/82e87x48ejrKwMn376qdH42rVrERAQgPj4eLOve+211wAACQkJGDBgAN5//30kJiaajT05ORn3338/7rnnHpPHHDx4ELGxsfDz88PkyZMRExODTZs2oUePHjh79mw1X5m6n+vhhx8GAOzcudOq9zFFPylfpVIZje/atQsATE6M79+/PwDdnZO1Pd5Stf2er1q1Cn369MG+ffswaNAgTJkyBc2bN8c777yDRx99FCUlJYZjx44dCxcXF6xfvx5CiCrn0idU+st6RE5P2YIXEZkyceJE4e/vL/bs2SNycnIMj+LiYpPH13Sprn///kIIIaKiokRkZKRhf05OjlCpVGLy5MlCCFHlss1PP/0kAIiOHTsKjUZjGM/OzhYhISECgBgzZozRex4/flwAEKNGjaoST3JysgAgAIjVq1cb7Vu9erUAIB5//PHqvzh1PFdBQYHhUlZdXLhwQXh4eIjQ0FBRVlZmtO+pp54SAMSRI0dMvrZx48aiRYsWtT6+OnX9nh8/flyoVCoRHR0t8vPzjc6dlJQkAIi///3vRuMDBgwQAMRPP/1kNK6/JOjh4SH++OMPi+IncnSsOBE5oFWrVuH69et45JFHEBoaanjs37+/TucdP348jh8/brgbbMOGDSgrKzN7yeazzz4DAMydOxeNGjUyjIeGhmLq1KkmX3P58mUAQEhIiNk47rvvPrz44otGYy+++CLatGmD77//Hr///rvF/yZrz+Xr6wtPT09DnLVRWlqK0aNHo7i4GIsXL4arq6vRfv1dkH5+fiZf7+vrazimNsdbw9rv+UcffYSysjKsWLECgYGBRvteffVVBAUFVbmLztwkcf0lwSeeeMLojlAiZ6aq+RAisjdh4pJHdcx1Fb/bc889h5kzZ2Lt2rXo3r071q1bh06dOpmdG5Weng4A6NWrV5V9PXv2NPmaP/74AwDg7+9vNo6ePXtWuZvLxcUFPXv2RFZWFtLT0xEXF2fBv6h25woICEB+fr5F57+bVqvF2LFjkZKSghdffBGjR4+u1Xnsxdrv+cGDBwEAP/74o8nLmW5ubvj111+Nxp544gkEBQVh06ZNKCgoMCSA+kSKl+moPmHiRNSABAUFYfDgwfjyyy/x9NNP49SpU1ixYoXZ4wsLC+Hi4oImTZpU2WeuoqS/3byoqMjsec29Vj9uTXWlNue6ffu2UQXNUlqtFuPHj8fGjRvx3HPPYfXq1SaP0ycO5v4dhYWFRvO/rD3eGtZ+z//8808AwDvvvGPxe7i5uWH06NFYunQpNm7ciIkTJyI3Nxf/+c9/EBYWZnESTOQMeKmOqIF54YUXUFhYiLFjx8LT0xOjRo0ye6yvry+0Wq3J6szVq1dNvkZ/95/+A9gUc6/Vj5u7ZCXHubRaLQoKCgxxWkqr1WLcuHHYsGEDRo4cifXr15vtgdSmTRsAQFZWVpV9ubm5uHnzpuGY2hxvLWu/54AuWRNCmH2Yeg9A1x8KAD799FOUlZVh3LhxsvUHI3IE/N9M1MD0798fzZo1w5UrVxAfH19tJSM6OhoAsG/fvir7zM23ioyMhIuLC06dOmX2vPv27YNWqzUa02q12L9/PyRJMryvJaw9V1ZWFrRabbUd1++mT5o++eQTjBgxAp9++mmVeU2V6dsX/Pe//62y78cffzQ6pjbHW8ua73n37t0BVFyys1S7du3w0EMPITU1FRkZGYaWFOPGjat13EQOSbl56URka3ffYaV3+PBhsWnTJqMmmkJUvcNqx44dAoDo1KmTuHXrlmE8JydHNG3a1ORddUII0bFjR+Hj41OlyWRt74SLiYkRAERycnKdz7VhwwYBQHz00UdV9plSXl4uxowZIwCIp59+WpSWltb4mtLSUtG6detqG1pW/tpbe3x16vo9P3bsmFCpVOL+++832dTy2rVr4ujRoybf++OPPxYARM+ePQUA8eijj1oUM5Ez4RwnogaoS5cuFjUjjIuLw7PPPouNGzeiffv2iI+PR3FxMb7++mt0794dW7duNXkZZujQoZg3bx4OHjyIHj16VNnfv39/TJkyBT/88AMiIyNx/PhxbN26FU2aNME//vGPKsfrK0p390uqzbl27NgBlUqFxx9/vMZ/PwC89dZb2LBhAxo3boz77rsPCxYsqHJMfHy80WRrlUqF//u//0P//v3Ru3dvPPPMM/Dx8cE333yDCxcu4O9//7vR+nDWHl8bln7Po6Ki8OGHH2LixIm4//778dhjjyEiIgI3btzA2bNnsXv3bowdO9bk/K4RI0YgMTHRUKHkpHCql5TO3IjIdsxVH8wxtfxGaWmpePvtt0V4eLhwd3cXrVu3FgsXLhSHDh0SAMTUqVOrnOfKlStCpVKJiRMnGo3rq0Tz5s0Te/bsETExMcLb21v4+vqKoUOHiqysrCrn0mq1IiAgQLRq1cqo2lObc2k0GtG4cWMRHx9v0ddDCGGoNlX3WLduncnXHjp0SAwYMED4+voKLy8v0a1bN/Hll1+afS9rjzdFju+5EEL88ssv4plnnhFqtVq4ubmJJk2aiM6dO4tZs2aJkydPmj3fuHHjBAAREBBQZSkaovqAiRMR1Yr+ssyHH35ocv9zzz0n7rnnHlFYWGgYq5zsWOrYsWMCgFi5cqXReG3OpY959+7dFr+GiKgyTg4nomrl5uZWuYvqypUrWLBgAVxdXc1e8lqwYAFu375d7a3vltizZw9CQkKqXVfNEmVlZVi4cCGGDBmC3r171+lcRNRwcY4TEVVr0aJF+P777/HII48gODgYFy9exLZt23Djxg3Mnz/faB29ylq2bIkNGzaYbRdgqYkTJ2LixIl1OgcAXLx4Ec8//7zDN6wkIsfGxImIqjVgwACcOHEC33//Pa5duwZPT0906NABL7/8Mp599tlqXzt8+HA7RVmz1q1bY/78+UqHQUROThJ31+CJiIiIyCTOcSIiIiKyEBMnIiIiIgsxcSIiIiKyUIOcHK7VapGdnQ0fHx9IkqR0OERERKQgIQRu3LgBtVpd46LUDTJxys7ONnsLNRERETVMly5dQvPmzas9pkEmTj4+PgB0XyBfX1+FoyEiIiIlFRYWokWLFob8oDoNMnHSX57z9fVl4kREREQAYNH0HU4OJyIiIrIQEyciIiIiCzFxIiIiIrJQg5zjREREZI3y8nKUlpYqHQbVkpubG1xdXWU5FxMnIiIiM4QQyM3NxfXr15UOherI398fTZs2rXP/RiZOREREZuiTpuDgYDRq1IhNk52QEAK3bt1CXl4eACA0NLRO52PiREREZEJ5ebkhaQoMDFQ6HKoDLy8vAEBeXh6Cg4PrdNmOk8OJiIhM0M9patSokcKRkBz038e6zlVj4kRERFQNXp6rH+T6PjJxIiIiIrIQEyciIiIiCzFxIiLHVaIB5vvpHiUapaMhqrVbJWVoNet7tJr1PW6VlCkdDtUBEyciIqJ6bMaMGRg8eLCs5/zb3/6G+Ph4Wc/pLJg4ERER1WOvv/46Nm7cKOs509LS0LFjR9nOl5KSgsGDB0OtVkOSJHz33Xdmjx03bhzeeOMNJCUloWvXrvDx8UFwcDDi4+Nx6tQp2WIyh4kTERGRHe3/7Q+7vp+/vz98fHxkPWd6erpFiVOfPn2wfv36Go/TaDSIjo7GypUrqz2uvLwc27Ztw5AhQ7B7924kJCTg4MGD2LFjB0pLS9GvXz9oNLa9rM8GmERERDYmhDA8X7bjNGLbBtulzUF+fj6CgoJw7NgxREVFyXLOy5cvIz8/H9HR0bKcDwAGDhyIgQMH1njc/v374ebmhq5du2L79u1G+9avX4/g4GCkpqaid+/essV2N1aciIiIbGzfmYoqU2Z2IVKy8u3yvunp6fDw8EDbtm2r7Fu4cCEaN25c7ePixYtVXpeWlgY/Pz+Eh4fb459gZMuWLRg8eLDJpLOgoAAAEBAQYNMYWHEiIiKyISEElu/MMmy7SMB7/z2F3m2a2LzqlJ6ejsjISKhUVT/uJ0yYgOHDh1f7erVaXWUsLS1N1mqTNTZv3oz333+/yrhWq0ViYiJ69uwpW2XNHCZORFR7JRpg4Z1frK9lA+7eysZD5IBSsvKRmV1o2NYKIONyAVKy8hFzX5BN37u6SdwBAQG1qs5UlzgtXLgQCxcuNGzfvn0bBw8exKRJkwxjJ06cQFhYmNXve/LkSWRnZyM2NrbKvoSEBGRmZmLv3r1Wn9davFRHRERkI0IIvPffU3C5q7CkrzpVnvtkC+np6dUmObW9VGcuGZswYQLS0tIMjy5duuCtt94yGjNVxbLEli1b8Oijj8LT09NofNKkSdi2bRuSk5PRvHnzWp3bGg5XcUpJScG7776L1NRU5OTkYNOmTWZ7RUyYMAEfffQR3n//fSQmJto1TiIiopqkZOUj43JBlXF7VJ1KSkpw8uRJs4lTbS7V3bhxA2fPnrW4iuXl5YXg4GDce++91gVvwubNm/HSSy8ZtoUQmDx5MjZt2oRdu3bZbc6VwyVO+lsSx48fjyeffNLscZs2bcLBgwdrnbkSUQPHy4xkY/pqkyQBpgpLko3nOp08eRKlpaVmE6faXKpLT0+Hq6srIiMj5QjR4ObNmzhz5oxh+9y5c0hLS0NAQADCwsKQl5eHI0eOYMuWLYZjEhISsHHjRmzevBk+Pj7Izc0FAPj5+cHLy0vW+CpzuMTJklsSr1y5gsmTJ+PHH3/EoEGD7BQZERGR5UrKtci+fttk0gTokqmc60UoKdfCQ+Uq+/unpaWhZcuW8Pf3l/Wcbdu2hYeHh2znBIAjR46gb9++hu3p06cDAMaMGYP169dj69at6NatG5o0aWI4ZtWqVQB0vaIqW7duHcaOHStrfJU5XOJUE61Wi9GjR2PGjBkWZ7zFxcUoLi42bBcWFlZzNBERUd15qFyxZVIv/KkpQVFpOZ5afQAA8O8JD8PTTZcoBTZ2t0nSBACHDx9Gt27dZD3npEmTjCZ612TXrl0WHdenT59q53tt3rwZQ4YMMRqz9fwwc5wucVq8eDFUKhWmTJli8WuSkpLw5ptv2jAqIiKiqtT+XlD7exkt7NtO7YtG7rb7+C0qKsLx48fxzTffYMGCBTZ7H3vq1asXRo4cqXQYAJwscUpNTcU//vEPHD161KrrwbNnzzaU/QBdxalFixa2CJGIiKiKRu4qnF9kn6kly5Ytw+LFizFixAg8//zzdnlPW3v11VeVDsHAqdoR7NmzB3l5eQgLC4NKpYJKpcKFCxfwyiuvoFWrVmZf5+HhAV9fX6MHERFRfTRr1ixcu3YNq1evhpubm9Lh1DtOVXEaPXo04uLijMb69++P0aNHY9y4cQpFRURERA2FwyVONd2SGBgYaHS8m5sbmjZtivvvv9/eoRIREVED43CJU023JBIREREpxeESp5puSbzb+fPnbRcMERERUSVONTmciIiISElMnIiIiIgsxMSJiIiIyEJMnIiIiIgsxMSJiIiIyEJMnIiIiGytRAPM99M9SjR2fesZM2Zg8ODBsp7zb3/7G+Lj42U9p7Ng4kRERFSPvf7669i4caOs50xLS0PHjh1lO19KSgoGDx4MtVoNSZLw3XffmT123LhxeOONN5CUlISuXbvCx8cHwcHBiI+Px6lTp2SLyRwmTkRERLZScBnITgNyMyrGcjN0Y9lpQMEVm4fg7+8PHx8fWc+Znp4ua+Kk0WgQHR2NlStXVntceXk5tm3bhiFDhmD37t1ISEjAwYMHsWPHDpSWlqJfv37QaGxb0XO4BphEJKMSDbBQrXv+Wjbg7q1sPEQNSVkxsKYvoMkzHl87oOJ542AgMRNQedgkhPz8fAQFBeHYsWOIioqS5ZyXL19Gfn4+oqOjZTkfAAwcOBADBw6s8bj9+/fDzc0NXbt2xfbt2432rV+/HsHBwUhNTUXv3r1li+1urDgRERHZgqs74NcM5j9qXQDfZrrjbCQ9PR0eHh5o27ZtlX0LFy5E48aNq31cvHixyuvS0tLg5+eH8PBwm8VtzpYtWzB48GBIklRlX0FBAQAgICDApjGw4kRERGQLkgT85Q3gs2FmDtDq9ptIAuSSnp6OyMhIqFRVP+4nTJiA4cOHV/t6tVpdZSwtLU3WapM1Nm/ejPfff7/KuFarRWJiInr27ClbZc0cJk5ERES2EhELqDsBORmAKK8Yl1yB0A66/TZU3STugICAWlVnqkucZs2ahcWLF1f7+pMnT5qsgNXk5MmTyM7ORmxs1a9ZQkICMjMzsXfvXqvPay0mTkRERLZiruokym1ebQJ0FacXXnjB5L6FCxdi4cKF1b7+xIkTCAsLMxpLS0vDY489ZvL4V155BWPHjq32nK1bt652vzlbtmzBo48+Ck9PT6PxSZMmYdu2bUhJSUHz5s1rdW5rMHEiIiKypYhYIDQayEnXbdup2lRSUoKTJ0+arQ7V5lLdjRs3cPbsWbNVrKCgIAQFBdUq3pps3rwZL730kmFbCIHJkydj06ZN2LVrl93mXDFxIiIisiVJAmJmAl8+q9u2U7Xp5MmTKC0tNZs41eZSXXp6OlxdXREZGSlHiAY3b97EmTNnDNvnzp1DWloaAgICEBYWhry8PBw5cgRbtmwxHJOQkICNGzdi8+bN8PHxQW5uLgDAz88PXl5essZXGRMnIiIiWwuPqXgeGm3zahOgu6TWsmVL+Pv7y3rOtm3bwsND3vYJR44cQd++fQ3b06dPBwCMGTMG69evx9atW9GtWzc0adLEcMyqVasAAH369DE617p162q8XFgXTJyIiIhsrXJ1qc9sm1ebAODw4cPo1q2brOecNGkSJk2aJOs5AV3yI4Qwu3/z5s0YMmSI0Vh1x9sS+zgREdmCgmuTkYMLt11zRgAoKipCamoqvvnmG/Tv39+m72UvvXr1wsiRI5UOAwATJyIionpl2bJliIuLwxNPPIHnn39e6XBk8eqrr6JFixZKhwGAl+qIiIhsz90bmF9gl7eaNWsWZs2aZZf3aohYcSIiIiKyEBMnIiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIqqGUv2CSF5yfR+ZOBEREZng5uYGALh165bCkZAc9N9H/fe1thyuHUFKSgreffddpKamIicnB5s2bUJ8fDwAoLS0FG+88QZ++OEHnD17Fn5+foiLi8OiRYuqLERIRERUF66urvD390deXh4AoFGjRpDs0PGb5CWEwK1bt5CXlwd/f3+4urrW6XwOlzhpNBpER0dj/PjxePLJJ4323bp1C0ePHsWcOXMQHR2Na9euYerUqRgyZAiOHDmiUMRERFRfNW3aFAAMyRM5L39/f8P3sy4cLnEaOHAgBg4caHKfn58fduzYYTT2wQcfoFu3brh48SLCwsLsESIRFVwGNPlA2e2KsdwMQHVnRXLvIMCvmTKxEclIkiSEhoYiODgYpaWlSodDteTm5lbnSpOewyVO1iooKIAkSdWu/lxcXIzi4mLDdmFhoR0iI6qnyoqBNX0BzV1/ga8dUPG8cTCQmAmoarmCOhMzcjCurq6yffCSc3PqxKmoqAgzZ87EyJEj4evra/a4pKQkvPnmm3aMjKgec3XXJS2afABaEwe4AL7NdMfVhj0SMyKiWnLau+pKS0sxfPhwCCGwatWqao+dPXs2CgoKDI9Lly7ZKUqiekiSgL+8AdNJE3Tjf3lDd1xt6BMzs7+e6piYERHVgVMmTvqk6cKFC9ixY0e11SYA8PDwgK+vr9GDiOogIhZQdwKkuy5dSK668YjY2p/b1okZEVEdOF3ipE+asrKy8NNPPyEwMFDpkIgaHn1yI8qNx0W5PEmNLRMzIqI6cLjE6ebNm0hLS0NaWhoA4Ny5c0hLS8PFixdRWlqKp556CkeOHMHnn3+O8vJy5ObmIjc3FyUlJcoGTtTQRMQCodEV23ImNbZOzIiIasnhEqcjR46gU6dO6NSpEwBg+vTp6NSpE+bOnYsrV65gy5YtuHz5Mjp27IjQ0FDDY//+/QpHTtTASBIQM7NiW+6kxpaJGRFRLTncXXV9+vSpdj0ZrhlE5EDCYyqeh0bLm9ToE7Mvn9Vts9pERA7A4SpORA1OiQaY76d7lGiUjsY6lZOYPrPlT2psmZgREdUCEycikkd4b/nPaevEjIjISg53qY6IZFAfO2/bIjEjIrISEyei+oadt4mIbIaX6ojqG3beJiKyGSZORPUNO28TEdkML9UR1Uf6zts5GcZNJCVXILRDw747rT7O/yIiu2HiRFQf6atOnw0zHm/ovZA4/4uI6oiX6ojqK3berorzv4iojpg4EdVXtl4SxRlx/hcR1RETJ6L6jJ23q9LP/5JcjcdZkSMiCzBxIqrP2Hm7Kn3VqfKkeYAVOSKyCBMnooaCnbcrcP4XEdUSEyciqrVbJWUmnzs8zv8iolpi4kREDRPnfxFRLTBxIiKHZdOKFud/EVEtMHEiIuL8LyKyEBMnIiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIiKyEBMnIiJnU6IB5vvpHiUapaMhalCYOBERERFZSKV0AEQNVsFlQJMPlN2uGMvNAFReuufeQYBfM2ViIyIik5g4ESmhrBhY0xfQ5BmPrx1Q8bxxMJCYCag8av02t0rK0Kjyc/dan4qIiMBLdUTKcHW/U00y9yPoAvg20x1HREQOw+ESp5SUFAwePBhqtRqSJOG7774z2i+EwNy5cxEaGgovLy/ExcUhKytLmWCJakuSdIvKQmvmAC0XnSUickAOlzhpNBpER0dj5cqVJvcvWbIEy5cvx+rVq3Ho0CF4e3ujf//+KCoqsnOkRHUUEQuoOwGSq/G45Kob56KzREQOx+HmOA0cOBADBw40uU8IgWXLluGNN97AE088AQD45JNPEBISgu+++w7PPPOMPUMlqht91emzYcbjopzVJjvg/C8iqg2HqzhV59y5c8jNzUVcXJxhzM/PD927d8eBAwfMvq64uBiFhYVGDyKHEBELhEZXbLPaRETk0JwqccrNzQUAhISEGI2HhIQY9pmSlJQEPz8/w6NFixY2jZPIYpIExMys2Ga1iYjIoTlV4lRbs2fPRkFBgeFx6dIlpUMiZ2LrLs3hMRXPQ6NZbXJ2BZeB7DRdTy693AzdWHYaUHBFocCISA4ON8epOk2bNgUAXL16FaGhoYbxq1evomPHjmZf5+HhAQ+P2vfCIbKpytWlPrNZbXJmdurPRUTKcaqKU3h4OJo2bYqdO3caxgoLC3Ho0CE8/PDDCkZGJJPw3kpHQHXB/lxE9Z7DVZxu3ryJM2fOGLbPnTuHtLQ0BAQEICwsDImJiViwYAHatGmD8PBwzJkzB2q1GvHx8coFTUQEmL9T0oD9uYicncMlTkeOHEHfvn0N29OnTwcAjBkzBuvXr8err74KjUaDl156CdevX0evXr2wfft2eHp6KhUyEVEFfX+unAzdZH89yRUI7cA5bEROzuESpz59+kAIYXa/JEl466238NZbb9kxKiIiC7E/F1G95lRznIiInAL7cxHVW0yciIjkxv5cRPUWEycihd0qKTP5nJwc+3MR1UtMnIiIbIH9uYjqJSZORES2xv5cRPWGw91VR0ROxN0brYo2AgBOuHsrHAwRke0xcSKqz5jYEBHJipfqiIiIiCzExImIiIjIQkyciIiIiCxUp8QpKSkJAHD06FEUFxfLEhARERGRo6rT5PA+ffoAABYtWoTjx4/DxcUFkZGR6NChAzp06ICuXbsiJCREjjiJ7K/gMqDJB8puV4zlZgAqL91z7yDAr5kysRERkSKsSpyEEJAqNXF7+OGHAQBff/01AOD27dvIzMxERkYGduzYgXnz5uGxxx7D22+/LWPIRHZQVgys6Qto8ozH1w6oeN44GEjMBFQe9o2NiIgUY1Xi9O677yIhIQHe3rrbmqdNm2aoLkVFRcHLywtdu3ZF165dDa958MEHmTiR83F111WTNPkAtCYOcAF8m+mOIyKiBsOqOU69evUyzGsCgL59++LKlStYvHgxOnfujKioKIwYMQLvvPMOtm7dCgA4ePCgvBET2YMk6RZlNZk0QTfORVtt704fqlZFGwH2oSIiB2BVxal58+b48ssvsWDBAgDAkCFDMGTIEMP+oqIiw6W6n3/+GYMHD4abm5u8ERPZS0QsoO4E5GToVrfXk1yB0A5ctJWIqAGyKnEqLi7GuXPnUFpaajIh8vT0RJcuXdClSxfZAiRSjL7q9Nkw43FRzmoTEVEDZdWlOq1WCyEEOnXqhFdffRXJyckoKyuzVWxEyouIBUKjK7YlV10VitUmIqIGyarEKSAgAM2bN0e7du2wfv16xMbGIjAwEEOHDsU///lPXLt2zVZxEilDkoCYmRXbtqg2cR4PEZHTsCpxCgoKwrBhw/D111/j6tWrOHDgAKZOnYrLly/jpZdeQosWLQytCYjqjfCYiueh0aw2kUVulZSZfE5Ezs3qzuGBgYEoLi6GJEno3r073nrrLRw+fBjZ2dlYsGABWw9Q/VO5utRnNuc2ERE1YFYnTpMnT8ZHH31UZTwkJASJiYlYunSpLIEROQqjykGzHgpGQkRESrM6cfLz88PIkSNx48YNk/sfffTROgdFRERE5IhqtVZdUFCQ3HEQERERObxaL/JbWlqK3Nxc3Lp1C0FBQQgICJAzLiIiIiKHY9Wluhs3bmDVqlWIiYmBr68vWrVqhQceeABBQUFo2bIlXnzxRRw+fNhWsRIRkT2UaID5frpHiUbpaIgcisWJ09KlS9GqVSusW7cOcXFx+O6775CWlobTp0/jwIEDmDdvHsrKytCvXz8MGDAAWVlZNgu6vLwcc+bMQXh4OLy8vBAREYG3334bQgibvSeRs9v/2x9Kh0B1VXAZyE4DcjMqxnIzdGPZaUDBFYUCI2o4LL5Ud/jwYaSkpCAyMtLk/m7dumH8+PFYtWoV1q9fjz179qBNmzayBVrZ4sWLsWrVKmzYsAGRkZE4cuQIxo0bBz8/P0yZMsUm70nkjCr/MbFsx2nEtg2GxHYKzqmsGFjTF9DkGY+vHVDxvHEwkJgJqDzsGxtRA2Jx4vTFF19YdJynpycmTJhQ64AssX//fjzxxBMYNGgQAKBVq1b44osv8Msvv9j0fcmBlWiAhWrd89ey2YH7jn1nKqpMmdmFSMnKR8x9trm5Y/9vfyDugRCbnJsAuLoDfs0ATT4ArYkDXADfZrrjiMhmrG5HUJ1Lly5h/Pjxcp7SpB49emDnzp04ffo0ACA9PR179+7FwIEDbf7eRLYk5+U0IQSW76y4ZO4iAe/995Ssl7TvrmjxcrkN6RedNpk0QTfOxaeJbE7WxOnPP//Ehg0b5DylSbNmzcIzzzyDtm3bws3NDZ06dUJiYiJGjRpl8vji4mIUFhYaPYgcha2Sj5SsfGRmV/xf1wog43IBUrLyZTk/YLqiRTYUEatbZFpyNR7n4tNEdmNVO4ItW7ZUu//s2bN1CsZSX3/9NT7//HNs3LgRkZGRSEtLQ2JiItRqNcaMGVPl+KSkJLz55pt2iY3IWra4nCaEwHv/PQUXSZcw6emrTr3bNKnzXCdzFS05zm0XdxZXBoATznJpV191+myY8bgtFp8mIpOsSpzi4+MhSVK1fxHb4xfmjBkzDFUnAGjfvj0uXLiApKQkk4nT7NmzMX36dMN2YWEhWrRoYfM4iWpiq+QjJSsfGZcLqoxXrjrVNTmrrqJlq3lUBF1VKTQayEnXbUuuQGgHVpuI7MSqS3WhoaH49ttvodVqTT6OHj1qqziN3Lp1Cy4uxqG7urpCqzV97d/DwwO+vr5GDyJHYIvLafpqk7m8S5JhrlPlilZltphHRXeRJCBmZsU2q01EdmVV4vTggw8iNTXV7P6aqlFyGTx4MN555x18//33OH/+PDZt2oSlS5di6NChNn9vIrnYKvkoKdci+/ptmHu5EEDO9SKUlJubZFwzfUVLe9d72GIeFZkQHlPxPDSa1SYiO7LqUt2MGTOg0ZjvInvvvfciOTm5zkHVZMWKFZgzZw5efvll5OXlQa1W469//Svmzp1r8/cmkoutLqd5qFyxZVIv/KkpQVFpOZ5afQAA8O8JD8PTTTepOLCxOzxUrtWdxqzKFS1TyZnkbHOdnFHlr2uf2aw2EdmRVYnTI488Uu1+b29vxMTEVHuMHHx8fLBs2TIsW7bM5u9FZAu2Tj7U/l5Q+3vhVkmZYayd2heN3Gu9PKWBNRWt2iZnprBPlBnhvZWOgKhBqftvUSKymlLJhxxsXdGqjJ3PicjRMHEiqokNblu3Z/JhC7asaFVmz87nRESWYOJEpBB7JR/Oyun7RBFRvSRr53Aiu7PzavFyLolC1bNH53MiImvV+U/by5cvQ61WV+mrRGRzdlotnvNs7M8enc+JiGqjztlOu3btcP78eRlCIbKSfrV4s/+N5Vktnuux2R/7RBGRo6pz4sQOwaQYO6wWb26eDf/f2449Op8TEdUWr6+Rc7uzWryw0WrxnGdjf/bofE5EVFtWz3H65JNPjLbLysrw7bffIjg42DD2/PPP1z0yIkvcqTpJNlgtnvNslOHsrRqIqH6zOnFat26d0XZpaSn+/e9/w8vLC4BuvTomTmRXEbEoD+kA16u6O+uE5ApJhtXibbUkCtWMrRqIyFFZ/Vvo7rXofHx8sHHjRrRu3Vq2oIisIkko7TUDrt+M1m3KWG3iemxERFQZ5zhRvaBtVbFeV3lI3atNnGdDRESmsO5N9UOlqk9p75lwrWMViPNsiIjIlDonTq+99hoCAgLkiIVIFtqWj8hyHs6zaTj2//YH4h4IUToMi90qKUOjys/r1qqMiKxQ50t1s2fPhr+/vwyhEBHZz90d4dkXiogswTlORNQgsSM8EdUGEycianDYEZ6IaouJExE1OHbpCO/ujVZFG9GqaCPg7i3feYlIURYnTnPnzkW7du0QFxeHefPmIS0trcoxP/zwA/7444+qLyYichCVO8JXxqoTEVnC4sTp5MmTGDp0KF555RUUFRVh1KhR6NKlC7755hvDMYGBgewaTkQOTd8RXntXfsR1CInIEhYnTsuWLUNWVhbGjh2L3377DWPHjkW/fv0wd+5cdO7cGcnJySgtLcW5c+dsGS+RWXvLoxBXvAT7z15XOhRyUJU7wpsisepERDWwuCFNs2bN8PXXX6OwsBA///wz9uzZgyNHjuDmzZv47bffMGTIENy6dQsbNmywZbxEJgkhsKRsBM6I5liWfB6x7cO4FApVYU1HeDY3JSJTrO7k5+vri/j4eMTHx1fZt2HDBnzzzTd47rnn5IiNyGL7zl5HhogAAGTm3OQCvGQSO8LXoOAyoMkHym5XjOVmACrdIu7wDgL8mikTG5GDkLUF8pgxY/DMM8/IeUqiGgkhsHzXBbigHFq4Gib5cgFeMoUd4c0oKwbW9AU0ecbjawdUPG8cDCRmAioP+8ZG5ECsbkdw+/Zt7N27FydOnKiyr6ioCF999ZUsgRFZKiUrH5k5N6GFrkrASb5EteDqfqeaZO5jwQXwbaY7jqgBsypxOn36NB544AH07t0b7du3R0xMDHJycgz7CwoKMG7cONmDpHqgRAPM99M9SjSynZa3lhPJRJKAv7wBQGvmAK1uP6u41MBZlTjNnDkTUVFRyMvLw6lTp+Dj44OePXvi4sWLtorPrCtXruC5555DYGAgvLy80L59exw5csTucZCyeGs5kYwiYgF1J0C6a46X5Kobj4hVJi4iB2LVRf39+/fjp59+QpMmTdCkSRNs3boVL7/8Mh555BEkJyfD29s+3XGvXbuGnj17om/fvvjPf/6DoKAgZGVl4Z577rHL+5NjqHxruanCkuQkc50auatwftEgpcMgqqg6fTbMeFyUy1dtKtEAC9W6569ls6s6OR2rEqfbt29Dpap4iSRJWLVqFSZNmoSYmBhs3LhR9gBNWbx4MVq0aIF169YZxsLDw+3y3uQ4eGs5kQ1ExAKh0UBOum5bcgVCO7DaRHSHVYlT27ZtceTIETzwwANG4x988AEAYMiQIfJFVo0tW7agf//+ePrpp7F79240a9YML7/8Ml588UW7vD85BnvdWs6KEDUokgTEzAS+fFa3LWe1iagesGqO09ChQ/HFF1+Y3PfBBx9g5MiRdpmMe/bsWaxatQpt2rTBjz/+iIkTJ2LKlClmm28WFxejsLDQ6EH1g9rfC1HN/NBO7WsYa6f2RVQzP0Q180Oon5eC0RE5qfCYiueh0aw2EVViVeI0e/Zs/PDDD2b3f/jhh9Bqzd2RIR+tVovOnTtj4cKF6NSpE1566SW8+OKLWL16tcnjk5KS4OfnZ3i0aNHC5jESETmtytWlPrNZbSKqxOLEydo7565cuWJ1MJYKDQ1Fu3btjMYeeOABszHOnj0bBQUFhselS5dsFhvdpeAykJ2m6z6sl5uhG8tOAwps9/+EyFHs/+0PpUOovfDeSkdA5FAsTpy6du2Kv/71rzh8+LDZYwoKCvDxxx8jKioK33zzjSwBmtKzZ0+cOnXKaOz06dNo2bKlyeM9PDzg6+tr9CA70HciXhNj3H147QDd2JoY4OM+uuOI6pnK0xaW7TjNnmJE9YTFk8NPnDiBd955B48++ig8PT3x4IMPQq1Ww9PTE9euXcOJEydw/PhxdO7cGUuWLMFjjz1ms6CnTZuGHj16YOHChRg+fDh++eUXrFmzBmvWrLHZe1It6DsRa/JhuqkeOxFT/bXvTEWVKTO7kOsnEtUTFlecAgMDsXTpUuTk5OCDDz5AmzZtkJ+fj6ysLADAqFGjkJqaigMHDtg0aQJ01a9Nmzbhiy++QFRUFN5++20sW7YMo0aNsun7kpXYiVhx+jsCzy8axLXY7EgIgeU7swzb7GRPVH9Y/ZvUy8sLTz31FJ566ilbxGOxxx9/HI8//riiMZAF7nQiFjkZkER5xTh7w1A9lpKVj8zsirt3K3eyZ9WJyLlZvcgvkVXuVJ2MkiaAvWGo3uL6iUT1GxMnsr2IWJSHdDBsCq57RfUY108kqt+YOJHtSRJKe82o2GS1ieqpyusnmiKx6kTk9Jg4kV1oW1X0gikP4dwmqp+sWT+RiJwTb7Mh+6j0J3hp75lwlbnaxPXkyBHYa/1EuHujVZFuUfUT7t51OxcRWYWJE9nN3vIozC97HrPK2yNO6WCIbETt7wW1vxdulZQZxtqpfdkOgqie4KU6sgshBJaUjcAZ0RzLks9zjgeRA6uc9FV+TkRMnMhO9p29jgwRAQDIzLnJO4uIiMgpMXEimxNCYPmuC3CBrpcT+9kQEZGzYuJENpeSlY/MnJvQQjchlv1siIjIWTFxIptiF2UiIqpPmDiRTbGLMhER1SdMnMhm2EWZiIjqGyZOZDPsokxERPUNO7KRzditizIphh3blbf/tz8Q90CI0mEQNRhMnMim2EWZSH6VL28v23EasW2DITn6otkFlwFNPlB2u2IsNwNQeemeewcBfs2UiY3ICvz0IiJyMvvO/GF4npldiJSsfMTcF6RgRDUoKwbW9AU0ecbjawdUPG8cDCRmAioP+8ZGZCXOcSIiciJCCCzfmWXYdorWHq7ud6pJ5j5yXADfZrrjiBwcEyciIieSkpWPzOxCw7ZTtPaQJOAvbwAwdyOIVrff0S83EoGJExGR03DqhrIRsYC6EyDddTOI5Kobj4hVJi4iKzFxIrvQ3311ftEgTgwnqiWnbiirrzqJcuNxUc5qEzkVJk5ERE6gXjSUjYgFQqMrtlltIifExIkqlGiA+X66R4lG6WiIqJJ60VBWkoCYmRXbtqg28fcY2RivmRAROYF601A2PKbieWg0q03kdJg4ERE5iXrRULZydanPbM5tIqfjRD9tZDPs6EtESgjvrXQERFZz+jlOixYtgiRJSExMVDoU56Tv6LsmxriL79oBurE1McDHfXTHERERNXBOnTgdPnwYH330ETp06KB0KM6LHX2JiIgs5rSJ082bNzFq1Ch8/PHHuOeee5QOx3mxoy8REZHFnDZxSkhIwKBBgxAXF6d0KM7vTkdfwY6+5GDYOFUh7t5oVbQRrYo2Au7eSkdD5FCc8jfRl19+iaNHj+Lw4cMWHV9cXIzi4oo5OoWFhdUc3QDdqTpJnw0zHmdHXyIiIiNOV3G6dOkSpk6dis8//xyenp4WvSYpKQl+fn6GR4sWLWwcpROKiEV5SAeUC12SVA4XCFabiEhmlVspVH5O5CycLnFKTU1FXl4eOnfuDJVKBZVKhd27d2P58uVQqVQoLy+v8prZs2ejoKDA8Lh06ZICkTs4SUJaxAS4Srq2xK7Q4th9k1ltIiIiqsTpEqfY2FgcO3YMaWlphkeXLl0watQopKWlwdW1atdcDw8P+Pr6Gj3ImBAC8040Q7o2HACQrm2NN44FO/a6V0Rkc/t/+0PpEIgcitMlTj4+PoiKijJ6eHt7IzAwEFFRUUqH57RSsvKRmavBkrJnkKVthiVlI5BxpdCxV1snIpuo/AfTsh2n+QcUUSVOlziR/PSrrrtIwD5tezxa8i72advDxRlWWyci2e07U1FlyszmH1BElTnlXXV327Vrl9IhOLWUrHxkXC6oMq4VQMblAqRk5SPmviAFIiMiexNCYPnOLMO2/g+o3m2aQOKcRyJWnBo6fbXJ3O9DiVUnogYlJSsfmdkVLVsq/wFFREycGrySci2yr9+GubxICCDnehFKys11Fiei+qLyZfvKeNmeqEK9uFRHteehcsWWSb3wp6YERaXleGr1AQDAvyc8DE833R2KgY3d4aGqerciEdUvTn3ZvuAyoMkHym5XjOVmACov3XPvoDvrchLVDRMngtrfC2p/L6NmdO3UvlzigqgBqXzZ3lRhSXLkuU5lxcCavoAmz3h87YCK542DgcRMQOVh39io3uGlOiIicu7L9q7ud6pJ5j7SXADfZrrjiOqIJQUy0C+oSkQNj1Nftr+z3ibuXm/TQMt1N0k2TJyIiAiAk1+2j4gF1J2AnAzdAuV6kisQ2oHrbpJsnOCngYhIfqyw1jPmqk6inNUmkhXnOBERUf0QEQuERldsS666KhSrTSQjJk5ERFQ/SBIQM7Nim9UmsgEmTkREVH+Ex1Q8D42Wt9pUogHm++keJRr5zktOhYkTERHVH5WrS31ms9pEsmPiREREdrW3PApxxUuw/+x1275ReG/bnp8aJCZORERkN8KtEZ4rfQ1nRHMs23WR69+R02HiREREdrPvzB+G55nZhUjJypf1/JV7UFV+TiQX9nEiIrIB9omqSgiB5TuzDNsujrz+HZEZrDgREZFdpGTlIzO70LCtFUDG5QLZq05EtsTEiYiIbE4Igff+ewoudxWW9FUnznUiZ8FLdc6g4DKgyQfKbgNrB+jGxm8HVF66595Bd1YGJyJyTClZ+ci4XFBlvHLVKea+IAUiI7IOEydHV1YMrOkLaPKMx/UJFAA0DgYSMwGVh31jIyKygL7aJEmAqcKSxLlO5ESYODk6V3ddNUmTD0Br4gAXwLeZ7jgiIgdUUq5F9vXbJpMmQJdM5VwvQkm5Fh4qV/sGZ4nKVX+93AxW/RsoJk6OztyK3wZarsVERA7NQ+WKLZN64U9NCYpKy/HU6gMAgH9PeBiebrpEKbCxu2MmTaz6012YODmDiFhA3QkiJwOSKDcMC8kVUmgHrvxNRA5P7e8Ftb+XUW+ldmpfNHJ38I8hVv3pLryrzhncqTpVTpoA6LZZbSJqcPQ9os4vGuT4iYez01f9TSZNAKv+DQ8TJychWv8Fp13vRZnQfcvKhAtOu7aBaP0XhSMjIqrn7lT9Id11KVFy1Y2z6t+gMHFyEiln/sCC28OgknR/9agkLRbcfhIplZYvICIiG9BXne6q+oNV/wbJKROnpKQkdO3aFT4+PggODkZ8fDxOnTqldFg2o7+Vd6/ogHRtawBAurY19ooObBxHRHSXveVRiCtegv1nr8t30ohYIDS6YpvVpgbLKROn3bt3IyEhAQcPHsSOHTtQWlqKfv36QaPRKB2aTegbx2mFhCVlI5ClbYYlZSOgFRKXKyAiqkQIgSVlI3BGNMey5PPy/WEpSUDMzEpvxGpTQ+WUswq3b99utL1+/XoEBwcjNTUVvXv3Vigq27i7cdw+bXs8WvKuYT8bxxERVdh39joyRAQAIDPnprwdycNjKp6HRstfbSrRAAvVuuevZQPu3vKen2ThlBWnuxUU6Nr4BwQEKByJ/KxpHEdE1JAJIbB81wW4QDcXSfZ18Cr/cdpnNqtNDZRTVpwq02q1SExMRM+ePREVFWXymOLiYhQXFxu2CwsLTR7niCo3jjPHYRvHEZFT0rc7cDYpWfnIzLkJQPf7UO518G6VlKGR/nmzHobn1LA4feKUkJCAzMxM7N271+wxSUlJePPNN+0Ylbz0jeOIiMg0/bQGF0mXMOm5cDoDycypL9VNmjQJ27ZtQ3JyMpo3b272uNmzZ6OgoMDwuHTpkh2jJCIiW6u4icZ4vHLViUgOTllxEkJg8uTJ2LRpE3bt2oXw8PBqj/fw8ICHB9cQIiKqj+6+ieZuvImG5OSUiVNCQgI2btyIzZs3w8fHB7m5uQAAPz8/eHnxkhYRUUNizU00DjkftOCybi28stsVY7kZgOrO55l3kG69PHIITpk4rVq1CgDQp08fo/F169Zh7Nix9g+IiIgUU/kmmqLScjy1+gAA4N8THoanmy5RctibaMqKgTV9AU2e8fjaARXPGwcDiZmAildOHIFTJk7slE1ERJXpb6K5VVJmGGun9nX8RZBd3XXVJE0+TC8k7AL4NtMdRw7BqSeHExER2ZPsy7no18EzmTRBN84O5Q6FiRMREZEFbLacS0Ssbt076a5LiVwPzyExcSIiIrvRN9c8v2iQ419Gu4up5Vxkoa86iXLjca6H55CYOBEREdXA5su5RMTq1r/Tk7vaVKIB5vvpHiUaec7ZQDFxIiIiqoF+ORetieVcZCFJKOrxt4ptVpscFhMnIiKialRezqUyuatO2la9Dc/LQzpwbpODcq4LzERERHamX87lbnIvIly5ulTaeyZc5ag2sbmm7Jg4ERERmaHUci7alo/U/SRsrmkTvFRHRERkhjXLuchB1j5R+uaaZj/q2VyzNlhxIiKiekPf7kAu9lzO5e4+UbHtw+pWxdK3OfhsmJkD2FyzNpg4ERERVcNey7mY6hNV57lT+uaaORnGfaIkVyCUE9Brg5fqiIiIFGazPlFsrik7Jk42sDcrH3FLd2OvXP09iIioXrNpn6iIWF17gzuELZZyaUANNpk4yUwIgSU//oozeTex5Mdf5esqS0RE9ZLN+0RJEkp7zajYlLPaVHAZyE7TtTjQy83QjWWnAQVX6v4eDoZznGRWud+HrP09iIioXrJHnyhtq964IbzgI91GgX8k/OSoNjXQdgesOMlICIF3t/9q2JZ9LSMiIqpXKveJMkWS6XNEAFhSOgJZ2mZYUjYCsnwqNdB2B0ycZJSSlY/M7ELDtuxrGRERkWL0rQ7OLxok2x119uoTte/sdXyq7YdHS97F5/n3yvO5pJ94DnOx1c92B7xUJ5PK16i1lX4AXGzUVZaIiJyfPfpEVb5jTwtXeT+XbNnuQL9cjDkKLRfDxEkmdlvLiIiI6hVb94nS37EHE3fs1flzyVyTzbpOQDc3f6oyheZP8VKdDOx1jZqIiMgaNr9jDzC0OygXujcphwtEXdsdOPD8KSZOMrD3WkZERESW0F8N0d71+STrHFxJQlrEBLhKujdxhRbH7ptct7lNDjx/ipfqZFD5GrU5cq1lRERE9ZPc6+xVvhpi6g97Saa5TkIIzDvRDO9owxHtcg7p2taYcywYm2NE3eZQOehyMUycZKK/Rk1EROQIrLkaUpc/7FOy8pGZq8ESl2cwX/UJlpSNQMaVwrrPobLV/Kk6YuJERERUD9nrjj39HKp92vZ4tORdADLeUR4RiwL/SHhfOwmVpEWZcMGtwEj4Krg4MRMnIiKieqry1RA5LwPq2fqOcgHgvbKn8ZY0HwCgkrR4r+xpzAegVIMfTg4nIiIiq9njjvKUrHx8kt8G6drWAIB0bWtsyItQtLE0EyciIiKymq3vKK+4DChhSVnFcjEukqRoix+nvlS3cuVKvPvuu8jNzUV0dDRWrFiBbt26KR0WERFRvWfrO8orXwbcJyrmTwHKNpZ22sTpq6++wvTp07F69Wp0794dy5YtQ//+/XHq1CkEBwcrHR4REVG9Z6s7yu3VSqE2nPZS3dKlS/Hiiy9i3LhxaNeuHVavXo1GjRph7dq1SodGREREdeDIjaWdsuJUUlKC1NRUzJ492zDm4uKCuLg4HDhwoMrxxcXFKC4uNmwXFhbaJU4iIiKyniM3lnbKxCk/Px/l5eUICQkxGg8JCcGvv/5a5fikpCS8+eab9gqPiIiI6shRG0s77aU6a8yePRsFBQWGx6VLl5QOiYiIiJyQU1acmjRpAldXV1y9etVo/OrVq2jatGmV4z08PODh4WGv8IiIiKiecsqKk7u7Ox588EHs3LnTMKbVarFz5048/PDDCkZGRERE9ZlTVpwAYPr06RgzZgy6dOmCbt26YdmyZdBoNBg3bpzSoREREVE95bSJ04gRI/D7779j7ty5yM3NRceOHbF9+/YqE8aJiIiI5CIJpXqWK6iwsBB+fn4oKCiAr6+v0uEQERGRgqzJC5xyjhMRERGREpg4EREREVmIiRMRERGRhZx2cnhd6Kd1cekVIiIi0ucDlkz7bpCJ040bNwAALVq0UDgSIiIichQ3btyAn59ftcc0yLvqtFotsrOz4ePjA0mSZD9/YWEhWrRogUuXLjndXXuMXRmMXRmMXRmMXRmM3TwhBG7cuAG1Wg0Xl+pnMTXIipOLiwuaN29u8/fx9fV1uv+ceoxdGYxdGYxdGYxdGYzdtJoqTXqcHE5ERERkISZORERERBZi4mQDHh4emDdvHjw8PJQOxWqMXRmMXRmMXRmMXRmMXR4NcnI4ERERUW2w4kRERERkISZORERERBZi4kRERERkISZOMlu5ciVatWoFT09PdO/eHb/88ovSIdUoKSkJXbt2hY+PD4KDgxEfH49Tp04pHVatLFq0CJIkITExUelQLHblyhU899xzCAwMhJeXF9q3b48jR44oHVaNysvLMWfOHISHh8PLywsRERF4++23LVqywN5SUlIwePBgqNVqSJKE7777zmi/EAJz585FaGgovLy8EBcXh6ysLGWCvUt1sZeWlmLmzJlo3749vL29oVar8fzzzyM7O1u5gCup6ete2YQJEyBJEpYtW2a3+KpjSewnT57EkCFD4OfnB29vb3Tt2hUXL160f7B3qSn2mzdvYtKkSWjevDm8vLzQrl07rF69Wplg72LJ51FRURESEhIQGBiIxo0bY9iwYbh69ardYmTiJKOvvvoK06dPx7x583D06FFER0ejf//+yMvLUzq0au3evRsJCQk4ePAgduzYgdLSUvTr1w8ajUbp0Kxy+PBhfPTRR+jQoYPSoVjs2rVr6NmzJ9zc3PCf//wHJ06cwHvvvYd77rlH6dBqtHjxYqxatQoffPABTp48icWLF2PJkiVYsWKF0qFVodFoEB0djZUrV5rcv2TJEixfvhyrV6/GoUOH4O3tjf79+6OoqMjOkVZVXey3bt3C0aNHMWfOHBw9ehTffvstTp06hSFDhigQaVU1fd31Nm3ahIMHD0KtVtspsprVFPtvv/2GXr16oW3btti1axcyMjIwZ84ceHp62jnSqmqKffr06di+fTs+++wznDx5EomJiZg0aRK2bNli50irsuTzaNq0adi6dSv+9a9/Yffu3cjOzsaTTz5pvyAFyaZbt24iISHBsF1eXi7UarVISkpSMCrr5eXlCQBi9+7dSodisRs3bog2bdqIHTt2iJiYGDF16lSlQ7LIzJkzRa9evZQOo1YGDRokxo8fbzT25JNPilGjRikUkWUAiE2bNhm2tVqtaNq0qXj33XcNY9evXxceHh7iiy++UCBC8+6O3ZRffvlFABAXLlywT1AWMhf75cuXRbNmzURmZqZo2bKleP/99+0eW01MxT5ixAjx3HPPKROQFUzFHhkZKd566y2jsc6dO4vXX3/djpFZ5u7Po+vXrws3Nzfxr3/9y3DMyZMnBQBx4MABu8TEipNMSkpKkJqairi4OMOYi4sL4uLicODAAQUjs15BQQEAICAgQOFILJeQkIBBgwYZff2dwZYtW9ClSxc8/fTTCA4ORqdOnfDxxx8rHZZFevTogZ07d+L06dMAgPT0dOzduxcDBw5UODLrnDt3Drm5uUb/d/z8/NC9e3en+9kFdD+/kiTB399f6VBqpNVqMXr0aMyYMQORkZFKh2MxrVaL77//Hvfddx/69++P4OBgdO/evdpLkY6kR48e2LJlC65cuQIhBJKTk3H69Gn069dP6dCquPvzKDU1FaWlpUY/r23btkVYWJjdfl6ZOMkkPz8f5eXlCAkJMRoPCQlBbm6uQlFZT6vVIjExET179kRUVJTS4Vjkyy+/xNGjR5GUlKR0KFY7e/YsVq1ahTZt2uDHH3/ExIkTMWXKFGzYsEHp0Go0a9YsPPPMM2jbti3c3NzQqVMnJCYmYtSoUUqHZhX9z6ez/+wCurkfM2fOxMiRI51iLbLFixdDpVJhypQpSodilby8PNy8eROLFi3CgAED8N///hdDhw7Fk08+id27dysdXo1WrFiBdu3aoXnz5nB3d8eAAQOwcuVK9O7dW+nQjJj6PMrNzYW7u3uVPwzs+fPaIBf5JfMSEhKQmZmJvXv3Kh2KRS5duoSpU6dix44dDjG3wFparRZdunTBwoULAQCdOnVCZmYmVq9ejTFjxigcXfW+/vprfP7559i4cSMiIyORlpaGxMREqNVqh4+9PiotLcXw4cMhhMCqVauUDqdGqamp+Mc//oGjR49CkiSlw7GKVqsFADzxxBOYNm0aAKBjx47Yv38/Vq9ejZiYGCXDq9GKFStw8OBBbNmyBS1btkRKSgoSEhKgVqsdqmrvqJ9HrDjJpEmTJnB1da0ys//q1ato2rSpQlFZZ9KkSdi2bRuSk5PRvHlzpcOxSGpqKvLy8tC5c2eoVCqoVCrs3r0by5cvh0qlQnl5udIhVis0NBTt2rUzGnvggQcc4s6cmsyYMcNQdWrfvj1Gjx6NadOmOV3lT//z6cw/u/qk6cKFC9ixY4dTVJv27NmDvLw8hIWFGX52L1y4gFdeeQWtWrVSOrxqNWnSBCqVyil/dm/fvo3XXnsNS5cuxeDBg9GhQwdMmjQJI0aMwN///nelwzMw93nUtGlTlJSU4Pr160bH2/PnlYmTTNzd3fHggw9i586dhjGtVoudO3fi4YcfVjCymgkhMGnSJGzatAk///wzwsPDlQ7JYrGxsTh27BjS0tIMjy5dumDUqFFIS0uDq6ur0iFWq2fPnlVutT19+jRatmypUESWu3XrFlxcjH+FuLq6Gv4adxbh4eFo2rSp0c9uYWEhDh065PA/u0BF0pSVlYWffvoJgYGBSodkkdGjRyMjI8PoZ1etVmPGjBn48ccflQ6vWu7u7ujatatT/uyWlpaitLTUYX92a/o8evDBB+Hm5mb083rq1ClcvHjRbj+vvFQno+nTp2PMmDHo0qULunXrhmXLlkGj0WDcuHFKh1athIQEbNy4EZs3b4aPj4/hOrGfnx+8vLwUjq56Pj4+VeZieXt7IzAw0CnmaE2bNg09evTAwoULMXz4cPzyyy9Ys2YN1qxZo3RoNRo8eDDeeecdhIWFITIyEv/73/+wdOlSjB8/XunQqrh58ybOnDlj2D537hzS0tIQEBCAsLAwJCYmYsGCBWjTpg3Cw8MxZ84cqNVqxMfHKxf0HdXFHhoaiqeeegpHjx7Ftm3bUF5ebvj5DQgIgLu7u1JhA6j56353kufm5oamTZvi/vvvt3eoVdQU+4wZMzBixAj07t0bffv2xfbt27F161bs2rVLuaDvqCn2mJgYzJgxA15eXmjZsiV2796NTz75BEuXLlUwap2aPo/8/PzwwgsvYPr06QgICICvry8mT56Mhx9+GA899JB9grTLvXsNyIoVK0RYWJhwd3cX3bp1EwcPHlQ6pBoBMPlYt26d0qHVijO1IxBCiK1bt4qoqCjh4eEh2rZtK9asWaN0SBYpLCwUU6dOFWFhYcLT01O0bt1avP7666K4uFjp0KpITk42+X98zJgxQghdS4I5c+aIkJAQ4eHhIWJjY8WpU6eUDfqO6mI/d+6c2Z/f5ORkpUOv8et+N0dqR2BJ7P/85z/FvffeKzw9PUV0dLT47rvvlAu4kppiz8nJEWPHjhVqtVp4enqK+++/X7z33ntCq9UqG7iw7PPo9u3b4uWXXxb33HOPaNSokRg6dKjIycmxW4zSnUCJiIiIqAac40RERERkISZORERERBZi4kRERERkISZORERERBZi4kRERERkISZORERERBZi4kRERERkISZORERERBZi4kRERERkISZORNSgXb58GU899RTmzJmjdChE5ASYOBFRgzZt2jS0adMG//rXv5QOhYicABMnImqwCgoKsGvXLvTq1QtqtVrpcIjICTBxIqIG6+eff8YjjzyC3bt3o2fPnkqHQ0ROgIkTETVYe/bswUMPPYRt27Zh6NChSodDRE6AiRMRNVhHjhzBtWvX4Onpic6dOysdDhE5AZXSARARKeX8+fP4/fff8eabbyodChE5CUkIIZQOgohICR4eHmjbti3+97//wcWFBXgiqhl/UxBRg+Xm5ob33nuPSRMRWYy/LYioQdqwYQM0Gg08PDxw8OBBbNu2TemQiMgJcI4TETU4RUVF+Pbbb7Fhwwa89NJLaNWqFT755BOlwyIiJ8A5TkREREQW4qU6IiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIiKyEBMnIiIiIgsxcSIiIiKy0P8Hdu1tRkL7ZMUAAAAASUVORK5CYII=", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig = plt.figure(figsize=(6, 4))\n", "plt.errorbar(\n", " ls,\n", " ci_plus * (2 * ls + 1),\n", " ci_plus_errs * (2 * ls + 1),\n", " linestyle=\"none\",\n", " marker=\"^\",\n", " label=\"$j = l + 1/2$\",\n", ")\n", "\n", "plt.errorbar(\n", " ls[1:],\n", " ci_minus[1:] * (2 * ls[1:] + 1),\n", " ci_minus_errs[1:] * (2 * ls[1:] + 1),\n", " linestyle=\"none\",\n", " marker=\"v\",\n", " label=\"$j = l - 1/2$\",\n", ")\n", "\n", "# plt.yscale(\"log\")\n", "plt.xticks([0, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20])\n", "plt.xlabel(r\"$\\ell$\")\n", "plt.ylabel(r\"$(2 \\ell +1) \\mathcal{T}_{lj}$\")\n", "plt.title(r\"$^{24}$Mg(p,p) 200 MeV\", fontsize=14)\n", "plt.legend()\n", "plt.tight_layout()" ] }, { "cell_type": "code", "execution_count": 26, "id": "a368a569-03be-4eec-9dfc-29c01784f571", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
E=65 MeV, j = l + 1/2E=65 MeV, j = l - 1/2E=65 MeV, j = l + 1/2, errE=65 MeV, j = l - 1/2, errE=200 MeV, j = l + 1/2E=200 MeV, j = l - 1/2E=200 MeV, j = l + 1/2, errE=200 MeV, j = l - 1/2, err
00.732306NaN0.103491NaN0.793631NaN0.083665NaN
10.7143970.7339120.1063250.1016220.7817210.8069500.0867710.084428
20.7074900.7402070.1010530.0939850.7653260.8096280.0890570.087703
30.6789680.7301610.1003710.0884270.7415240.8089220.0904510.094959
40.6563830.7252340.0959360.0753560.7094590.8044220.1016300.100298
50.6813550.6898160.0740870.0820020.6640130.7973200.1112310.107043
60.5899530.5389680.0835850.0959310.6051660.7814230.1274920.114318
70.3589190.3293510.0674990.0656630.5371470.7514750.1540210.130638
80.1697560.1654480.0358430.0343350.4657120.6981400.1783900.146219
90.0722020.0732350.0172260.0164290.3931390.6216290.1780450.154005
100.0294730.0308090.0079920.0079030.3197780.5271260.1628890.147182
110.0120780.0126680.0038800.0036370.2491480.4229260.1256350.130248
120.0049440.0052160.0018460.0017900.1875340.3192090.0896970.102278
130.0020360.0021630.0008900.0008560.1367520.2298550.0612290.077375
140.0008410.0008960.0004240.0004250.0962770.1581280.0416660.055007
150.0003560.0003770.0002110.0002070.0655680.1058520.0277580.039233
160.0001520.0001590.0001000.0000950.0437960.0695920.0189470.025934
170.0000650.0000680.0000470.0000460.0290160.0445600.0130110.017157
180.0000280.0000290.0000230.0000220.0190500.0285280.0088260.011763
190.0000120.0000120.0000110.0000110.0124000.0181590.0062080.007886
200.0000050.0000050.0000050.0000050.0080330.0115650.0042650.005293
\n", "
" ], "text/plain": [ " E=65 MeV, j = l + 1/2 E=65 MeV, j = l - 1/2 E=65 MeV, j = l + 1/2, err \\\n", "0 0.732306 NaN 0.103491 \n", "1 0.714397 0.733912 0.106325 \n", "2 0.707490 0.740207 0.101053 \n", "3 0.678968 0.730161 0.100371 \n", "4 0.656383 0.725234 0.095936 \n", "5 0.681355 0.689816 0.074087 \n", "6 0.589953 0.538968 0.083585 \n", "7 0.358919 0.329351 0.067499 \n", "8 0.169756 0.165448 0.035843 \n", "9 0.072202 0.073235 0.017226 \n", "10 0.029473 0.030809 0.007992 \n", "11 0.012078 0.012668 0.003880 \n", "12 0.004944 0.005216 0.001846 \n", "13 0.002036 0.002163 0.000890 \n", "14 0.000841 0.000896 0.000424 \n", "15 0.000356 0.000377 0.000211 \n", "16 0.000152 0.000159 0.000100 \n", "17 0.000065 0.000068 0.000047 \n", "18 0.000028 0.000029 0.000023 \n", "19 0.000012 0.000012 0.000011 \n", "20 0.000005 0.000005 0.000005 \n", "\n", " E=65 MeV, j = l - 1/2, err E=200 MeV, j = l + 1/2 \\\n", "0 NaN 0.793631 \n", "1 0.101622 0.781721 \n", "2 0.093985 0.765326 \n", "3 0.088427 0.741524 \n", "4 0.075356 0.709459 \n", "5 0.082002 0.664013 \n", "6 0.095931 0.605166 \n", "7 0.065663 0.537147 \n", "8 0.034335 0.465712 \n", "9 0.016429 0.393139 \n", "10 0.007903 0.319778 \n", "11 0.003637 0.249148 \n", "12 0.001790 0.187534 \n", "13 0.000856 0.136752 \n", "14 0.000425 0.096277 \n", "15 0.000207 0.065568 \n", "16 0.000095 0.043796 \n", "17 0.000046 0.029016 \n", "18 0.000022 0.019050 \n", "19 0.000011 0.012400 \n", "20 0.000005 0.008033 \n", "\n", " E=200 MeV, j = l - 1/2 E=200 MeV, j = l + 1/2, err \\\n", "0 NaN 0.083665 \n", "1 0.806950 0.086771 \n", "2 0.809628 0.089057 \n", "3 0.808922 0.090451 \n", "4 0.804422 0.101630 \n", "5 0.797320 0.111231 \n", "6 0.781423 0.127492 \n", "7 0.751475 0.154021 \n", "8 0.698140 0.178390 \n", "9 0.621629 0.178045 \n", "10 0.527126 0.162889 \n", "11 0.422926 0.125635 \n", "12 0.319209 0.089697 \n", "13 0.229855 0.061229 \n", "14 0.158128 0.041666 \n", "15 0.105852 0.027758 \n", "16 0.069592 0.018947 \n", "17 0.044560 0.013011 \n", "18 0.028528 0.008826 \n", "19 0.018159 0.006208 \n", "20 0.011565 0.004265 \n", "\n", " E=200 MeV, j = l - 1/2, err \n", "0 NaN \n", "1 0.084428 \n", "2 0.087703 \n", "3 0.094959 \n", "4 0.100298 \n", "5 0.107043 \n", "6 0.114318 \n", "7 0.130638 \n", "8 0.146219 \n", "9 0.154005 \n", "10 0.147182 \n", "11 0.130248 \n", "12 0.102278 \n", "13 0.077375 \n", "14 0.055007 \n", "15 0.039233 \n", "16 0.025934 \n", "17 0.017157 \n", "18 0.011763 \n", "19 0.007886 \n", "20 0.005293 " ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# NBVAL_CHECK_OUTPUT\n", "data" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.12.3" } }, "nbformat": 4, "nbformat_minor": 5 }