{ "cells": [ { "cell_type": "markdown", "id": "80bc7461-98bf-4a6a-a951-2c1666d43e7b", "metadata": {}, "source": [ "# Elastic scattering with `LocalOpticalPotential`\n", "\n", "`LocalOpticalPotential` provides a useful interface for defining optical potentials. " ] }, { "cell_type": "code", "execution_count": 1, "id": "709f54ea-cd6d-4a2a-9733-65f983b9562d", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "from matplotlib import pyplot as plt\n", "from pandas import DataFrame\n", "from tqdm import tqdm" ] }, { "cell_type": "code", "execution_count": 2, "id": "742d1f28-6eb3-4225-ba5d-8cc8b344287e", "metadata": {}, "outputs": [], "source": [ "from jitr.optical_potentials import LocalOpticalPotential\n", "from jitr.reactions.reaction import Reaction\n", "from jitr.rmatrix import Solver as SolverKernel\n", "from jitr.utils import utils\n", "from jitr.xs import elastic" ] }, { "cell_type": "code", "execution_count": 3, "id": "0bdc7e4a-c473-414b-92ea-a8e26683550a", "metadata": {}, "outputs": [], "source": [ "neutron = (1, 0)\n", "proton = (1, 1)" ] }, { "cell_type": "code", "execution_count": 4, "id": "0d708966-dbed-4ca3-a2e3-673aac93d2d8", "metadata": {}, "outputs": [], "source": [ "target = (208, 82)\n", "projectile = proton\n", "energy_lab = 24\n", "reaction = Reaction(target=target, projectile=projectile, process=\"El\")\n", "kinematics = reaction.kinematics(energy_lab)" ] }, { "cell_type": "code", "execution_count": 5, "id": "ff28c9ba-c458-45b7-9318-b807c3768f71", "metadata": {}, "outputs": [], "source": [ "# set the channel radius, number of nodes, and number of partial waves\n", "interaction_range_fm = 1.2 * (208 ** (1 / 3)) + 1\n", "channel_radius_dimensionless = utils.suggested_dimensionless_channel_radius(\n", " interaction_range_fm, kinematics.k\n", ")\n", "channel_radius = channel_radius_dimensionless / kinematics.k\n", "N = utils.suggested_basis_size(channel_radius_dimensionless)\n", "lmax = 50" ] }, { "cell_type": "code", "execution_count": 6, "id": "dd405939-51a3-4768-b5c2-7da5f4d57f52", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Compiling solver for 208-Pb(p,el) at 24 MeV\n", " - channel radius 13.94 fm\n", " - 25 nodes\n", " - 50 partial waves\n", "Done!\n" ] } ], "source": [ "# build a solver for the system and reaction of interest\n", "print(f\"Compiling solver for {reaction} at {energy_lab} MeV\")\n", "print(f\" - channel radius {channel_radius:1.2f} fm\")\n", "print(f\" - {N} nodes\")\n", "print(f\" - {lmax} partial waves\")\n", "\n", "solver = elastic.DifferentialWorkspace.build_from_system(\n", " reaction=reaction,\n", " kinematics=kinematics,\n", " channel_radius_fm=channel_radius,\n", " solver=SolverKernel(N),\n", " lmax=lmax,\n", " angles=np.linspace(0.1, np.pi, 180),\n", ")\n", "rgrid = solver.radial_grid()\n", "print(\"Done!\")" ] }, { "cell_type": "markdown", "id": "ee4c4c55-7a28-4f86-8683-03226733aced", "metadata": {}, "source": [ "## Use the built-in `LocalOpticalPotential` class\n", "\n", "Although, in general, one can define whatever interaction they want. This is just a handy tool." ] }, { "cell_type": "code", "execution_count": 7, "id": "d7c6b33c-d696-4751-94d9-8eb193006a8c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "['Vv',\n", " 'rv',\n", " 'av',\n", " 'Wv',\n", " 'rw',\n", " 'aw',\n", " 'Wd',\n", " 'Vd',\n", " 'rd',\n", " 'ad',\n", " 'Vso',\n", " 'Wso',\n", " 'rso',\n", " 'aso',\n", " 'rC']" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "omp = LocalOpticalPotential()\n", "omp.params" ] }, { "cell_type": "markdown", "id": "fbc1744e-1b30-44ab-8eae-c7a9e3ea051a", "metadata": {}, "source": [ "## Randomly generate a bunch of parameters" ] }, { "cell_type": "code", "execution_count": 8, "id": "721bd33e-c8c0-485d-8a2f-19e80a48b217", "metadata": {}, "outputs": [], "source": [ "means = np.array([56, 1.2, 0.7, 5, 1.2, 0.7, 13, 6, 1.3, 0.9, 8, 4, 1.1, 0.7, 1.2])\n", "std_devs = np.array(\n", " [3, 0.05, 0.03, 0.1, 0.05, 0.05, 1, 0.1, 0.1, 0.1, 0.2, 1, 0.1, 0.05, 0.1]\n", ")\n", "samples = np.random.multivariate_normal(means, np.diag(std_devs) ** 2, 1000)" ] }, { "cell_type": "code", "execution_count": 9, "id": "01571595-10ac-4d1d-8fdf-9b0967e60f3f", "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", "
VvrvavWvrwawWdVdrdadVsoWsorsoasorC
051.1004031.1932290.6940535.0011681.2048450.69862511.5935926.0023581.1805020.9871478.1704913.0768801.0856940.7341791.339886
152.8339761.2196700.6869454.8594371.2600520.71127811.8069695.9408251.2413720.9045508.0853284.3802841.1519380.7211911.342671
260.4906541.1116950.6903345.0177231.2665760.70978311.5676326.2078101.3800220.8020548.0889283.9474221.0423290.7715211.315199
355.8713651.1946290.7464195.1189071.1671640.68269013.3383356.0520001.3099500.5990377.9062944.4706661.1949990.6752121.157940
455.4605031.1175430.7332194.8964551.2100000.65896214.6175976.3075001.3415880.8905948.1724654.1205611.1020400.6200341.185585
\n", "
" ], "text/plain": [ " Vv rv av Wv rw aw Wd \\\n", "0 51.100403 1.193229 0.694053 5.001168 1.204845 0.698625 11.593592 \n", "1 52.833976 1.219670 0.686945 4.859437 1.260052 0.711278 11.806969 \n", "2 60.490654 1.111695 0.690334 5.017723 1.266576 0.709783 11.567632 \n", "3 55.871365 1.194629 0.746419 5.118907 1.167164 0.682690 13.338335 \n", "4 55.460503 1.117543 0.733219 4.896455 1.210000 0.658962 14.617597 \n", "\n", " Vd rd ad Vso Wso rso aso \\\n", "0 6.002358 1.180502 0.987147 8.170491 3.076880 1.085694 0.734179 \n", "1 5.940825 1.241372 0.904550 8.085328 4.380284 1.151938 0.721191 \n", "2 6.207810 1.380022 0.802054 8.088928 3.947422 1.042329 0.771521 \n", "3 6.052000 1.309950 0.599037 7.906294 4.470666 1.194999 0.675212 \n", "4 6.307500 1.341588 0.890594 8.172465 4.120561 1.102040 0.620034 \n", "\n", " rC \n", "0 1.339886 \n", "1 1.342671 \n", "2 1.315199 \n", "3 1.157940 \n", "4 1.185585 " ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "df = DataFrame(samples, columns=omp.params)\n", "df.head()" ] }, { "cell_type": "markdown", "id": "c9c0561d-13a8-4583-8810-eca698cc8a3a", "metadata": {}, "source": [ "## How do we calculate observables? \n", "We are using a `jitr.xs.elastic.DifferentialWorkspace`, which is set up specifically to interface with `omp` and other classes that have the same structure" ] }, { "cell_type": "code", "execution_count": 10, "id": "d637e670-caa5-4ecd-8530-0607460ca7d7", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Help on function xs in module jitr.xs.elastic:\n", "\n", "xs(self, central_potential: 'npt.ArrayLike', spin_orbit_potential: 'npt.ArrayLike | None' = None, coulomb_potential: 'npt.ArrayLike | None' = None) -> 'ElasticXS'\n", " Return differential and integral elastic observables.\n", "\n" ] } ], "source": [ "help(elastic.DifferentialWorkspace.xs)" ] }, { "cell_type": "markdown", "id": "771c868f-c329-47bd-869e-b2d3577dfd5d", "metadata": {}, "source": [ "Exactly the information that the solver workspace needs is what is provided by the `omp` class:" ] }, { "cell_type": "code", "execution_count": 11, "id": "d01c95e0-a51e-4070-81da-7f36d8d010ce", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Help on method evaluate in module jitr.optical_potentials.omp:\n", "\n", "evaluate(rgrid: 'ArrayOrScalar', reaction_model: 'reaction.Reaction', kinematics_model: 'kinematics.ChannelKinematics', Vv: 'float', rv: 'float', av: 'float', Wv: 'float', rw: 'float', aw: 'float', Wd: 'float', Vd: 'float', rd: 'float', ad: 'float', Vso: 'float', Wso: 'float', rso: 'float', aso: 'float', rC: 'float') -> 'tuple[PotentialArray, PotentialArray, PotentialArray | ArrayOrScalar]' method of jitr.optical_potentials.omp.LocalOpticalPotential instance\n", " Evaluate the local optical-potential terms on ``rgrid``.\n", "\n" ] } ], "source": [ "help(omp.evaluate)" ] }, { "cell_type": "markdown", "id": "f36b7420-4174-480b-af65-73853132f1e9", "metadata": {}, "source": [ "## Do you see the vision?" ] }, { "cell_type": "markdown", "id": "bf658a71-30bb-48b2-b1da-7e2fdc863b8e", "metadata": {}, "source": [ "Now running calculations is simple!\n", "Once it's been compiled, it's fast:" ] }, { "cell_type": "code", "execution_count": 12, "id": "e6aa31ae-92ad-4ecf-8c46-c40630237a0c", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "100%|██████████████████████████████████████████████████████████| 1000/1000 [00:09<00:00, 100.51it/s]" ] }, { "name": "stdout", "output_type": "stream", "text": [ "CPU times: user 9.68 s, sys: 265 ms, total: 9.94 s\n", "Wall time: 9.95 s\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\n" ] } ], "source": [ "%%time\n", "num_samples, num_params = df.shape\n", "xs_ratio = np.zeros((num_samples, solver.angles.size))\n", "Ay = np.zeros((num_samples, solver.angles.size))\n", "rgrid = solver.radial_grid()\n", "\n", "for i in tqdm(range(num_samples)):\n", " central_term, spin_orbit_term, coulomb_term = omp(\n", " rgrid,\n", " reaction,\n", " kinematics,\n", " *samples[i, :],\n", " )\n", " xs = solver.xs(central_term, spin_orbit_term, coulomb_term)\n", " xs_ratio[i, :] = xs.dsdo / solver.rutherford\n", " Ay[i, :] = xs.Ay" ] }, { "cell_type": "code", "execution_count": 13, "id": "2c5914f7-73a5-4493-9bcc-a8a93a2c4af1", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, '$\\\\frac{d \\\\sigma}{d\\\\Omega} / \\\\frac{d \\\\sigma_{\\\\text{R}}}{d\\\\Omega}$ [dimensionless]')" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAAG0CAYAAAAxRiOnAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAAVI9JREFUeJzt3Xl8VOWhPvDnzJ5JZib7HraAAgJhkxgVqpKCS3GvVFEoF+mF67VgilVqhetywbogt62Vyg8VWxeqVWvVskVQUWRJSAGBQFiSkD0kmck62zm/P2JGxgQyk5zJbM/38xklZ97znncmA3nynncRJEmSQERERBSmFP5uABEREZE/MQwRERFRWGMYIiIiorDGMERERERhjWGIiIiIwhrDEBEREYU1hiEiIiIKayp/NyDQiaKIyspKGAwGCILg7+YQERGRByRJQnNzM1JTU6FQXLzvh2GoF5WVlcjIyPB3M4iIiKgPysvLkZ6eftEyDEO9MBgMADrfTKPR6OfWEBERkScsFgsyMjJcP8cvhmGoF123xoxGI8MQERFRkPFkiAsHUBMREVFYYxgiIiKisMYwRERERGGNYYiIiIjCGsMQERERhTWGISIiIgprDENEREQU1hiGiIiIKKwxDBEREVFYYxgiIiKisMYwRERERGGNYYiIiIjCGsMQERERhTWGoSAgipK/m0BERBSyVP5uAPXM4RSxv7QRZefaUNvcAZ1aiRFJBoxJNSIuSuvv5hEREYUMhqEA1NBqw6eHqlDXbHUdszsdKCxtxOEKM26dkIa06Ag/tpCIiCh08DZZgGloteHtvWVuQeh8NoeIDw9UoLyhbYBbRkREFJoYhgLMsSoLbA7xomVsDhH/KKpAQ6ttgFpFREQUuhiGAszxmuZuxwShezm7U8LWb6s5uJqIiKifGIYCSK2lA41tdrdjBp0Kd0xMR2ykplv5KnMHCsoaB6p5REREIYlhKIAcr2lx+3pIvB5zsgcjI1aPrIzoHs/55uQ51Lf0PL6IiIiIescwFEB+eIts2ogERGiUAIDRKUZoVN2/XQ5RwvYjNZAk3i4jIiLqC4ahAFFt7oC5/ftbZHFRGrf1hDQqBUanGns8t8rcgYNnzT5vIxERUShiGAoQP+wVGp4Y1a3M+PToHgdTA8BXJ+vRanX4omlEREQhjWEoQJT+YN2gEYmGbmViIjUYHKfv8XyrXcTnx+t80jYiIqJQxjAUADrsTpw7bxB0jF6NBEPPW26Mz4i5YD3F1c04XMHbZURERN5gGAoAFU3tOH/88/AeeoW6DInT9zjNvstnx2pxtpGrUxMREXmKYSgAVDS2u309Iqn7eKEugiBgwqDoCz7vFCV8fLAK5h+sV0REREQ9YxgKABVN34chY4QaSUbdRcuPSjG6ptz3pN3mxHuFZ7ldBxERkQcYhvzM5hBRa/l+vFBmQmSv56iVCoxNM120jKXdjk37ynnLjIiIqBdBFYa++OILzJo1C6mpqRAEAR9++GGv5+zcuRMTJ06EVqvF8OHD8frrr/u8nd6oMrdDPG/AUGbChW+RnS8rIxpq5QXm2X+nw+7EB4UVKCxr5KKMREREFxBUYai1tRVZWVl46aWXPCp/+vRp3HTTTbj22mtRVFSEpUuX4v7778eWLVt83FLPnT9eSKdWIi06wqPzorQqTBoc22s5hyjh8+I6bNpX7jZjjYiIiDqp/N0Ab9xwww244YYbPC6/bt06DB06FC+88AIAYNSoUdi1axdefPFFzJw5s8dzrFYrrNbvQ4PFYulfo3tx9rzxQkPjI6FQXLy353yXD4nBkSoLLO29D5auMnfgnX3luH5Msse9T0REROEgqHqGvLV7927k5ua6HZs5cyZ27959wXNWr14Nk8nkemRkZPisfQ6niBpzh+vr4Ym9jxc6n0qpwI8uife4vM0h4p//rsQB7nRPRETkEtJhqLq6GklJSW7HkpKSYLFY0N7e3uM5y5cvh9lsdj3Ky8t91j5LhwMOsXMsj0ohYFCsd2EI6FyTaEh8z6tS90SSgJ3FdThS6dseLyIiomAR0mGoL7RaLYxGo9tjIAyK0/e4K70nZoxOhkHn3R3Pncdr0dzBtYiIiIhCOgwlJyejpqbG7VhNTQ2MRiMiIjwbqDxQetqY1VORWhVuGZ/mVZiy2kVsO1LTe0EiIqIQF9JhKCcnB/n5+W7Htm3bhpycHD+1qGdxURqMSu5fD1SCQYvrxyRfcFf7npSea8PBs039ui4REVGwC6ow1NLSgqKiIhQVFQHonDpfVFSEsrIyAJ3jfebOnesqv2jRIpw6dQq//vWvcezYMfzpT3/C3/72Nzz00EP+aP4FTRuR4NUssgvJTIjC1BGeD6gGgG9OnYPDKfb72kRERMEqqMLQ/v37MWHCBEyYMAEAkJeXhwkTJmDFihUAgKqqKlcwAoChQ4fik08+wbZt25CVlYUXXngB/+///b8LTqv3h6HxkRgS7/3A6QuZNDgW49Ivvjr1+VqtThytapbt+kRERMFGkLg08UVZLBaYTCaYzWbZB1M3tdkgSrjoLvR9IYoSPjhQgbIGz7biiNGrMe/KIRC8ucdGREQUwLz5+R1UPUOhJlqvkT0IAYBCIeDHlyX1ul1Hl8Y2O07UtsjeDiIiomDAMBSijDq1R9t1dNl/hgsxEhFReGIYCmGXD4mBMULtUdkaSwcqm3peiJKIiCiUMQyFMJVSgWlezC47WsVVqYmIKPwwDIW4EUkGxEd5Ni7peE0LnCLH0xMRUXhhGAoDo1I8mwXXYXfidH2rj1tDREQUWBiGwsClyQaPV6Y+Vs1bZUREFF4YhsKAQadGRoxnO9ufrmtFh93p4xYREREFDoahMDEyxeBROYco4UQN1xwiIqLwwTAUJkYkGjxehLGkjttzEBFR+GAYChMalQKZCVEela1s6oDIWWVERBQmGIbCSGaiZ2HI5hBR09zh49YQEREFBoahMJIWHeFx2bONXI2aiIjCA8NQGInUqjzeGPZso2c73hMREQU7hqEw42nvEMcNERFRuGAYCjPpsZ6FIZtDRLWF44aIiCj0MQyFGY4bIiIicscwFGYMOjWi9WqPynLcEBERhQOGoTDkae9QlZnjhoiIKPQxDIWhdA/3KbM5RNQ2W33cGiIiIv9iGApDaTGejxviIGoiIgp1DENhyBShhkGn8qhstZlhiIiIQhvDUJiKj9J6VK6W23IQEVGIYxgKUwkGz8JQQ6sNVofTx60hIiLyH4ahMOVpz5AkAbUWDqImIqLQxTAUpjztGQI4iJqIiEIbw1CYio5QQ60UPCpbwzBEREQhjGEoTCkUAuI8vFXGGWVERBTKGIbCmKfjhpo7HGi1OnzcGiIiIv9gGApj8VEaj8ty3BAREYUqhqEw5s0gao4bIiKiUMUwFMY8vU0GcHo9ERGFLoahMKZTK2GMUHtUto4bthIRUYhiGApzno4barE60G7jStRERBR6GIbCnDfjhupb2DtEREShh2EozHkzbqiOYYiIiEIQw1CYi9F7Pr2e44aIiCgUMQyFudhIDRSCZ9ty8DYZERGFIoahMKdUCDBFqDwq29BigyhKPm4RERHRwGIYIsREenarzCFKaGiz+bg1REREA4thiBAXyRllREQUvhiGCLEe9gwBQH0ze4aIiCi0MAyRV2GoroV7lBERUWhhGCL2DBERUVhjGCJoVAoYdJ7NKOO2HEREFGoYhgiAl71DHERNREQhhGGIAHgXhs618lYZERGFDoYhAuBlGGLPEBERhRCGIQLAniEiIgpfDEMEwLuFF8+1MAwREVHoYBgiAECERokIjdKjsh12J1qtDh+3iIiIaGAwDJFLrN6bcUPsHSIiotDAMEQu0Xq1x2XPtXIQNRERhYagC0MvvfQShgwZAp1Oh+zsbOzdu/ei5deuXYtLL70UERERyMjIwEMPPYSODm4p0RNPd68H2DNEREShI6jC0KZNm5CXl4eVK1eisLAQWVlZmDlzJmpra3ss/9Zbb+HRRx/FypUrcfToUWzYsAGbNm3Cb37zmwFueXCI8aJnqIEzyoiIKEQEVRhas2YNFi5ciPnz52P06NFYt24d9Ho9Xn311R7Lf/3117jqqqtwzz33YMiQIZgxYwbuvvvui/YmWa1WWCwWt0e4iPFizFA9b5MREVGICJowZLPZUFBQgNzcXNcxhUKB3Nxc7N69u8dzrrzyShQUFLjCz6lTp/Dpp5/ixhtvvOB1Vq9eDZPJ5HpkZGTI+0ICWLReA0HwrKzVLqKFM8qIiCgEBE0Yqq+vh9PpRFJSktvxpKQkVFdX93jOPffcgyeffBJXX3011Go1MjMzcc0111z0Ntny5cthNptdj/LycllfRyBTKgQYdV4MouZK1EREFAKCJgz1xc6dO7Fq1Sr86U9/QmFhId5//3188skneOqppy54jlarhdFodHuEk5hIb2aUcdwQEREFP5W/G+Cp+Ph4KJVK1NTUuB2vqalBcnJyj+c8/vjjuO+++3D//fcDAMaOHYvW1lb84he/wGOPPQaFIqSzYJ/E6DU4gzaPynJGGRERhYKgSQMajQaTJk1Cfn6+65goisjPz0dOTk6P57S1tXULPEpl5yrLkiT5rrFBzJtB1A0cRE1ERCEgaHqGACAvLw/z5s3D5MmTMWXKFKxduxatra2YP38+AGDu3LlIS0vD6tWrAQCzZs3CmjVrMGHCBGRnZ6OkpASPP/44Zs2a5QpF5M6bMMTbZEREFAqCKgzNnj0bdXV1WLFiBaqrqzF+/Hhs3rzZNai6rKzMrSfot7/9LQRBwG9/+1tUVFQgISEBs2bNwv/+7//66yUEvGgvxgx1zSiL0gbVx4iIiMiNIPF+0UVZLBaYTCaYzeawGEwtSRJe2lECu9Ozj8UdE9MxKE7v41YRERF5x5uf30EzZogGhiAIiPbqVhnHDRERUXBjGKJuvBtEzXFDREQU3BiGqBtv9ijjIGoiIgp2/Rr5+tFHH3l9zo9//GNERET057LkY97sXs+eISIiCnb9CkO33nqrV+UFQcCJEycwbNiw/lyWfCzai56hdpsTrVYHIjmjjIiIglS/b5NVV1dDFEWPHno9Zx0FA2/GDAHsHSIiouDWrzA0b948r2553XvvvWExPT3Y6dRKRGg8X5SS44aIiCiY9evexmuvveZV+Zdffrk/l6MBFB2hRrvN6VFZbstBRETBTLbZZO3t7Whr+36Dz9LSUqxduxZbt26V6xI0gLxaa4gbthIRURCTLQzdcssteOONNwAATU1NyM7OxgsvvIBbbrmFPUJByJvp9RwzREREwUy2MFRYWIipU6cCAN577z0kJSWhtLQUb7zxBn7/+9/LdRkaIN70DLXZnB7fUiMiIgo0soWhtrY2GAwGAMDWrVtx++23Q6FQ4IorrkBpaalcl6EB4k3PEMBtOYiIKHjJFoaGDx+ODz/8EOXl5diyZQtmzJgBAKitreUMsiDkTc8QwHFDREQUvGQLQytWrMCyZcswZMgQZGdnIycnB0BnL9GECRPkugwNEI1KgUitN9Pr2TNERETBSbZlg++8805cffXVqKqqQlZWluv49OnTcdttt8l1GRpA0XoNWq3tHpWtZ88QEREFKdnCUHt7O4xGI5KTkwF0Tq3/4IMPMGrUKEyZMkWuy9AAitFrUNHoWRjibTIiIgpWPp9af+utt3JqfZDyZhB1h92JFqvDh60hIiLyDU6tpwvyZsNWADjXwnFDREQUfDi1ni7I2xllHDdERETBiFPr6YKiI9QQBM/L17NniIiIghCn1tMFqZQKGHSe3yrjIGoiIgpGnFpPFxUbqYal3e5R2YZWKyRJguBNdxIREZGfyRaGACA5Odk1tb4Lp9UHtxi9BmfQ5lFZu1OCud3u9VgjIiIif5LtNhkAfPnll7j33nuRk5ODiooKAMBf/vIX7Nq1S87L0ACKi9R6VZ6DqImIKNjIFob+/ve/Y+bMmYiIiMCBAwdgtXYOpjWbzVi1apVcl6EBFhvl7R5lHERNRETBRbYw9PTTT2PdunVYv3491OrvB91eddVVKCwslOsyNMBivd2wtZU9Q0REFFxkC0PFxcWYNm1at+MmkwlNTU1yXYYGWIRGCb3G8w1bOb2eiIiCjWxhKDk5GSUlJd2O79q1C8OGDZPrMuQHsZGe9w41tNpgd4o+bA0REZG8ZAtDCxcuxJIlS7Bnzx4IgoDKykq8+eabWLZsGRYvXizXZcgPvAlDksTeISIiCi6yTa1/9NFHIYoipk+fjra2NkybNg1arRbLli3Dgw8+KNdlyA+8CUMAUNdsRYopwketISIikpdsYUgQBDz22GN4+OGHUVJSgpaWFowePRpRUVFyXYL8xNvp9bUW9gwREVHwkHXRRQDQaDQYPXq03NWSH8VEerd7fR1vkxERURDpVxjKy8vzuOyaNWv6cynyI4NODY1KAZvDs4HR51qsEEUJCgW35SAiosDXrzB04MABj8pxr6rgFxupQbW5w6OydqeEhjYb4qO8u71GRETkD/0KQzt27JCrHRTgvAlDQOcgaoYhIiIKBrLuTUahK87LGWW1zRw3REREwUHWAdT5+fnIz89HbW0tRNF9fMmrr74q56VogPVlej0REVEwkC0MPfHEE3jyyScxefJkpKSkcJxQiIk3eHfLi2GIiIiChWxhaN26dXj99ddx3333yVUlBRCjTg2tWgGr3bMZZR12J8ztdpgivJuWT0RENNBkGzNks9lw5ZVXylUdBaB4LxdfrGv2fMA1ERGRv8gWhu6//3689dZbclVHASje4N24ocomhiEiIgp8st0m6+jowCuvvILt27dj3LhxUKvdb49w0cXg5+22HN5MxSciIvIX2cLQwYMHMX78eADA4cOH3Z7jYOrQ4O0g6hpLB5yiBCVXoiYiogAmWxjiAoyhLz7Ku9tkDlFCXbMVySadj1pERETUf7KuM9TU1IQNGzbg6NGjAIDLLrsM//Ef/wGTySTnZchPtColjBFqWNrtHp9TaW5nGCIiooAm2wDq/fv3IzMzEy+++CIaGhrQ0NCANWvWIDMzE4WFhXJdhvzM294hjhsiIqJAJ1vP0EMPPYSbb74Z69evh0rVWa3D4cD999+PpUuX4osvvpDrUuRH8VFanKpr9bh8ZVO7D1tDRETUf7KFof3797sFIQBQqVT49a9/jcmTJ8t1GfIzbzdfbe5woMXqQJRW1juyREREspHtNpnRaERZWVm34+Xl5TAYDHJdhvwszsvbZABQxd4hIiIKYLKFodmzZ2PBggXYtGkTysvLUV5ejnfeeQf3338/7r77brkuQ34Wq9d4PVW+iuOGiIgogMl27+L555+HIAiYO3cuHA4HAECtVmPx4sV45pln5LoM+ZlCISAuSoNai+cbsXLcEBERBTJBkiRJzgrb2tpw8uRJAEBmZib0er2c1Q84i8UCk8kEs9kMo9Ho7+YEhB3HalFU3uRxeYUgYNE1w6BVKX3XKCIiovN48/Nb9lGter0eY8eOlbtaCiCp0RFehSFRklDR2I5hCVG+axQREVEf9SsM5eXl4amnnkJkZCTy8vIuWlauvcleeuklPPfcc6iurkZWVhb+8Ic/YMqUKRcs39TUhMceewzvv/8+GhoaMHjwYKxduxY33nijLO0JR6nR3i+iWM4wREREAapfYejAgQOw2+2uP1+IXHuTbdq0CXl5eVi3bh2ys7Oxdu1azJw5E8XFxUhMTOxW3maz4cc//jESExPx3nvvIS0tDaWlpYiOjpalPeHKoFPDoFOhucPh8TnlDW0+bBEREVHfyT5myJeys7Nx+eWX449//CMAQBRFZGRk4MEHH8Sjjz7arfy6devw3HPP4dixY1Cr1R5dw2q1wmr9fnCwxWJBRkYGxwz9wL8OVeFYdbPH5QUB+M9pmYjQcNwQERH5njdjhmSbWt/e3o62tu9/+y8tLcXatWuxdetWWeq32WwoKChAbm6u65hCoUBubi52797d4zkfffQRcnJy8MADDyApKQljxozBqlWr4HQ6L3id1atXw2QyuR4ZGRmytD/UpEZHeFVekoDyRvYOERFR4JEtDN1yyy144403AHSO05kyZQpeeOEF3HLLLXj55Zf7XX99fT2cTieSkpLcjiclJaG6urrHc06dOoX33nsPTqcTn376KR5//HG88MILePrppy94neXLl8NsNrse5eXl/W57KPI2DAG8VUZERIFJtjBUWFiIqVOnAgDee+89JCcno7S0FG+88QZ+//vfy3UZr4iiiMTERLzyyiuYNGkSZs+ejcceewzr1q274DlarRZGo9HtQd3FR2mgVXv38SljGCIiogAkWxhqa2tzbbuxdetW3H777VAoFLjiiitQWlra7/rj4+OhVCpRU1PjdrympgbJyck9npOSkoJLLrkESuX341RGjRqF6upq2Gy2frcpnAmCgFSTd71DTW12WDrsPmoRERFR38gWhoYPH44PP/wQ5eXl2LJlC2bMmAEAqK2tlaV3RaPRYNKkScjPz3cdE0UR+fn5yMnJ6fGcq666CiUlJRBF0XXs+PHjSElJgUbj/R5b5K4vt8pKalt80BIiIqK+ky0MrVixAsuWLcOQIUOQnZ3tCihbt27FhAkTZLlGXl4e1q9fj40bN+Lo0aNYvHgxWltbMX/+fADA3LlzsXz5clf5xYsXo6GhAUuWLMHx48fxySefYNWqVXjggQdkaU+4S4/xPgwVlTUhiCYwEhFRGJBtBeo777wTV199NaqqqpCVleU6Pn36dNx2222yXGP27Nmoq6vDihUrUF1djfHjx2Pz5s2uQdVlZWVQKL7PdxkZGdiyZQseeughjBs3DmlpaViyZAkeeeQRWdoT7lJMOug1SrTZLjw774fM7XacrGvB8ESDD1tGRETkuaBaZ8gfuDfZxW35thpHKi1enZMWE4G7JnPJAiIi8h2/7U2Wn5+P/Px81NbWuo3TAYBXX31VzktRgBgWH+l1GKpobEeNpQNJRu+39SAiIpKbbGOGnnjiCcyYMQP5+fmor69HY2Oj24NC06A4PZQK77dbKSzlZ4KIiAKDbD1D69atw+uvv4777rtPriopCGhVSqTHRKD0nHdrCBXXNCMrI7pPM9KIiIjkJFvPkM1mw5VXXilXdRREhsZHen2OJAHbjtTA4RR7L0xERORDsoWh+++/H2+99ZZc1VEQGZYQ1afzGlpt+OZUg8ytISIi8o5st8k6OjrwyiuvYPv27Rg3bly3XeLXrFkj16UowJgi1EgwaFHXbPX63ILSRmQmRiLFy9WsiYiI5CJbGDp48CDGjx8PADh8+LDbc4Lg/QBbCi6Th8TgX4d63jD3YkRJwt8LziJ3dBJGJnPpAiIiGniyhaEdO3bIVRUFoUuTDNh3ugH1Ld7v+WZ3SvjXoWpUmztwxbA46NTK3k8iIiKSiWxjhgDgyy+/xL333osrr7wSFRUVAIC//OUv2LVrl5yXoQAkCAKyh8X1q44DZU145YtT+EdRBUpqm7ltBxERDQjZwtDf//53zJw5ExERESgsLITV2jl+xGw2Y9WqVXJdhgLYiMQoxBu0/arDKUo4VdeKf/67Chu/PoPDFWaGIiIi8inZwtDTTz+NdevWYf369W6Dp6+66ioUFhbKdRkKYIIg4MrM/vUOna+xzY5tR2rwyaEqTsEnIiKfkS0MFRcXY9q0ad2Om0wmNDU1yXUZCnCZCVEYmSzvJqwnalrwfmEFOuyebwhLRETkKdnCUHJyMkpKSrod37VrF4YNGybXZSgIXDsyEQadrNveoaKpHe/sLUNTm/cDtImIiC5GtjC0cOFCLFmyBHv27IEgCKisrMSbb76JZcuWYfHixXJdhoKATq3Ej0cnQe4VFRrb7HhnXzkqmtrlrZiIiMKabL++P/rooxBFEdOnT0dbWxumTZsGrVaLZcuW4cEHH5TrMhQkBsdFYsrQWOyReYXpdpsT7xecxZ2T07lQIxERyUKQZJ6qY7PZUFJSgpaWFowePRpRUX3bqiFQWCwWmEwmmM1mGI1cFNBbO4prUVTWJHu9Bp0Kd08ZhEitvLfjiIgoNHjz81v2nyQajQajR4+Wu1oKUtdckgC7Q8S3lRZZ623ucOCTQ1W4c2I6FAqucE5ERH0naxjq6OjAwYMHUVtbC1F0nwp98803y3kpChKCIODHo5OgVAg4eNYsa90Vje3YUVyL6aOSZK2XiIjCi2xhaPPmzZg7dy7q6+u7PScIApxOTosOV4IgYPqoJOjUSuw9Le8YooNnzdCoFJg6IkHWeomIKHzINpvswQcfxE9/+lNUVVVBFEW3B4MQAcBVw+MxdUS87PXuP9OIPafOyV4vERGFB9nCUE1NDfLy8pCUxFsWdGGTh8TiupGJsk+7//rkORyvaZa3UiIiCguyhaE777wTO3fulKs6CmFZGdE+WYdo25EaLspIRERek21qfVtbG376058iISEBY8eOddufDAB++ctfynGZAcep9b5TUNqAL453H2PWH4lGLWZPzoBKKVvOJyKiIOSXqfVvv/02tm7dCp1Oh507d0I479d+QRCCNgyR70waHAtLuwNF5U2y1VlrsWJXST2uuTRRtjqJiCi0yRaGHnvsMTzxxBN49NFHoVDwt3LyzDWXJsDSYcepulbZ6iwqb8LwxCikx+hlq5OIiEKXbKnFZrNh9uzZDELkFUEQMPOyZETr1b0X9pAkdY4fsjvF3gsTEVHYky25zJs3D5s2bZKrOgojOrUSPxmXCrVSvhHVTW127CqRdzwSERGFJtlukzmdTjz77LPYsmULxo0b120A9Zo1a+S6FIWgBIMW141MwpZvq2Wr89/lTRidYkSSUSdbnUREFHpkC0OHDh3ChAkTAACHDx92e06Qew41haTRqUYcrjSjorFdlvokCfj8eB3umpwhS31ERBSaZAtDO3bskKsqCmNXD4/Hpn3lstVX0diOktpmDE80yFYnERGFFo52poCSGh2BYQmRstb55Yl6OEVZltMiIqIQ1K+eoby8PDz11FOIjIxEXl7eRctyzBB56qrh8Thd3wp5lgPtHEz977NNmDgoRp4KiYgopPQrDB04cAB2u9315wvhmCHyRnyUFiOTDThaJd9eYwVnGpGVHg2lgp9FIiJy168wdP44IY4ZIjmNS4+WNQy1WB04WmXBmDSTbHUSEVFo4JghCkip0RGIN2hlrbOgtBEybcVHREQhpN9jhjzFMUPkrbFpJuw4VitbfQ2tNpysa8XwxCjZ6iQiouDX7zFD5yssLITD4cCll14KADh+/DiUSiUmTZrUn8tQmBqZbMCuE3WwO+Xrzdl/poFhiIiI3Mg2ZmjNmjUwGAzYuHEjYmI6Z+00NjZi/vz5mDp1av9aSWFJp1ZiRJIBRyotstVZZe5AjaWDq1ITEZGLIMk0iCItLQ1bt27FZZdd5nb88OHDmDFjBiorK+W4zICzWCwwmUwwm80wGo3+bk7YqWxql3URRgAYnxGNa0cmylonEREFFm9+fss2gNpisaCurq7b8bq6OjQ3yzcriMJLanSErDvaA0BxTTMXYSQiIhfZwtBtt92G+fPn4/3338fZs2dx9uxZ/P3vf8eCBQtw++23y3UZCkOXJMm7lUa7zYnT9S2y1klERMFLtr3J1q1bh2XLluGee+5xLcSoUqmwYMECPPfcc3JdhsLQJUkG7D3dIGudR6q4XxkREXWSbcxQl9bWVpw8eRIAkJmZichIefeZGmgcMxQY/rL7DOpbbLLVp1QIWDh1GCI0StnqJCKiwDFgY4YOHjwIURTdjkVGRmLcuHEYN25ctyD07bffwuFw9OeSFKZGyHyrzClKOFYt3yw1IiIKXv0KQxMmTMC5c+c8Lp+Tk4OysrL+XJLC1KUyhyEAOFLFMERERP0cMyRJEh5//HHo9XqPytts8t3moPASE6lBolGLWotVtjprLVbUNVuRIPO2H0REFFz6FYamTZuG4uJij8vn5OQgIiKiP5ekMHZJkkHWMAR09g79yJAga51ERBRc+hWGdu7cKVMziHo3PCEKu07Uy1pncbUFU4fHQ6EQZK2XiIiCB3etp6ARE6lBfJRG1jpbrU6cOdcqa51ERBRcGIYoqGT6YJNVDqQmIgpvDEMUVEb4YKHEU3WtaLc5Za+XiIiCA8MQBZUEg1b2vcqcosTeISKiMOazMKRUcmVf8o3hPrhVdrjCLHudREQUHHwWhmTe5cPNSy+9hCFDhkCn0yE7Oxt79+716Lx33nkHgiDg1ltv9VnbyPd8EYYaWm0ob2iTvV4iIgp8soWhb7/91u1rQfDNVOVNmzYhLy8PK1euRGFhIbKysjBz5kzU1tZe9LwzZ85g2bJlmDp1qk/aRQMn2aiDQSfbHsMu7B0iIgpPfQ5Dr732Gn7+859j69atWLx4Md5++20523VBa9aswcKFCzF//nyMHj0a69atg16vx6uvvnrBc5xOJ+bMmYMnnngCw4YNu2j9VqsVFovF7UGBRRAEDEuQfwPgktoWDqQmIgpDfQ5DX3/9NV577TU8+eSTWLVqFZ5++ukey33zzTe46667kJWVhUmTJmHJkiWoqKhATk6O19e02WwoKChAbm7u9y9AoUBubi527959wfOefPJJJCYmYsGCBb1eY/Xq1TCZTK5HRkaG1+0k38tMkP9WmUOUcKSKvUNEROGmz2Govr4eO3bsQHJyMg4cOIDPPvvM7XlJkvDmm29i6dKlWLRoEbZv3453330XmZmZuPbaa3Hy5Mk+XdPpdCIpKcnteFJSEqqrq3s8Z9euXdiwYQPWr1/v0TWWL18Os9nsepSXl3vdTvK99Bg9tGr5h7wdqWRPIBFRuPH6p8l7770HALj11ltRXl6OWbNmoby8HGfPnu1WdvXq1di8eTOuu+46JCQkYNiwYfjlL3+JTz75BEOHDu1/63vR3NyM++67D+vXr0d8fLxH52i1WhiNRrcHBR6lQsDQOPlvldW32FDb3CF7vUREFLi8HoU6e/Zs/OlPf0JTUxMmTpyI66677oLT6EVRRHR0dLfjI0aMwIYNG7xubHx8PJRKJWpqatyO19TUIDk5uVv5kydP4syZM5g1a5ZbmwBApVKhuLgYmZmZXreDAkNmYhSOVTfLXu+xqmYkGnSy10tERIHJ654hSZKwZ88eKJVKbNy4EdnZ2Thz5kyPZdPT03vczPX5559HVlaWt5eGRqPBpEmTkJ+f7zomiiLy8/N7HIM0cuRIHDp0CEVFRa7HzTffjGuvvRZFRUUcDxTkhsRFQuWDDVaLq5t9ujQEEREFFq97hi655BK3mVv79u3Df/3Xf+HTTz91KycIAl5++WXceeedGDduHMaNG4fm5mZ8/PHHyMzMxIgRI/rU4Ly8PMybNw+TJ0/GlClTsHbtWrS2tmL+/PkAgLlz5yItLQ2rV6+GTqfDmDFj3M7v6qn64XEKPhqVAhmxepyul3ej1RarA2UNbRjsg9twREQUeLwOQ7GxsThy5AhGjx4NALj88stRUVHRY9nMzEwUFBRg8+bNOHbsGGJiYrB+/XpMmDABv/vd7/rU4NmzZ6Ourg4rVqxAdXU1xo8fj82bN7sGVZeVlUGh4C4j4WJ4YpTsYQgAjlY1MwwREYUJQfLyfsDevXtx3333Ydq0aRgzZgwOHjyIyspK/Otf/3Irp1AoUFJS0uu6PoHOYrHAZDLBbDZzMHUAarc58coXpyDKfFtLo1LgF9OGQa1ksCYiCkbe/Pz2umdoypQpKCoqwpYtW3DkyBFMnjwZc+bM6bHsokWLcPz4caSkpLhulXU9TCaTt5cm6iZCo0R6TATKZN5Kw+YQcaquFZcmG2Stl4iIAo/HPUMnT570auaVUqmE09m5mu+qVauwb98+jBo1CgUFBdi+fTuGDh2KkpKSvrV6ALFnKPAdPNuE/KMX346lL4YlROKW8Wmy10tERL7nk56hRYsWoaSkBMnJyV738vztb39DUVGR6+utW7fizTff9PTSRBc1PDEKnx2rhdwTwErPtaHD7oRO3fPSEUREFBo8HhCxbds2nD59GrNmzUJtbS0qKirw9NNPIzY2FsOHD7/ouTqdDkeOHHF9PWPGDBw+fLjvrSY6j16jQlp0hOz1OkUJJbUtstdLRESBxesxQ33p5dmwYQNmz56Na665BuPHj8ehQ4d8tqs9hacRSQacbWyXvd7i6maMSeP4NiKiUOb1VJm+9PJcdtllKCgowNSpU3HmzBkMHjy42+wzov4YnhgFX+Tr8sY2tFod8ldMREQBw+ueoR/28hw+fNijXh6NRoO77roLd911V58aSnQxUVoVMmL0ss8qkyTgeE0zJgyKkbVeIiIKHF73DP2wl2f9+vVYt25dt3LczoAG2qgU38z2+5Y72RMRhTSve4YA916es2fP4sYbb4RGo8GsWbNw8803Y/r06a4NUYkGyvDEKOwoVsDmkPezV9dsRUVTu08GaRMRkf/1e3nd1157DdXV1Xj77bdhMBiwdOlSxMfH44477sAbb7yBhoYGOdpJ1CuNSoHMBN9soVFU1uSTeomIyP9k2WtAoVBg6tSpePbZZ1FcXIw9e/ZgypQp+POf/4zU1FRMmzYNzz///AX3MCOSi69ulZXUtqC5w+6TuomIyL98svHSqFGj8Mgjj+Crr75CeXk55s2bhy+//BJvv/22Ly5H5DIoVo8obZ/u/l6UKEk4dNYse71EROR/Xm/U6qnzt+MIZtyOI/h8eaIO+880yl6vXqPEgquHQsXNW4mIAp43P7/79K/6e++912sZziYjf7ks1TeLJLbZnDhYwd4hIqJQ06cwNHv2bPz5z3/G7373O2zbtq3HHqCutYe++eYb3HXXXcjKysKkSZOwZMkSVFRUICcnp38tJ7qA2EgN0mJ8M/Nrz6kGtNuCv8eTiIi+16cwJEkS9uzZA6VSiY0bNyI7OxunT5/uVubNN9/E0qVLsWjRImzfvh3vvvsuMjMzce211+LkyZOyvACinozxUe9Qh92J3afqfVI3ERH5R5/GDI0cORLHjh1zfb1v3z6sXLkSn376qeuYQqHA6NGjsWvXLkRHR7udf+LECdx7773Ys2dP31s+QDhmKDjZnSLWf3kKVrv8610pBAFzrhiE+Cit7HUTEZE8fD5mKDY21m1/sssvv7zHafOiKHYLQgAwYsQIbNiwoS+XJvKIWqnAyGSDT+oWJQmfHauFKHJcHBFRKOhTGFq7di1uu+02LFy4EP/3f/+HBQsWIDU1tVu59PR07Ny5s9vx559/HllZWX25NJHHfHWrDAAqGtux+9Q5n9VPREQDx6sFWU6ePInMzExMmTIFRUVF2LJlC44cOYLJkydjzpw5bmUFQcDLL7+MO++8E+PGjcO4cePQ3NyMjz/+GJmZmRgxYoSsL4TohxKNOgyKlX/z1i77zjQg2aRDZkKUT+onIqKB4VXP0KJFizB06FDk5ORg6dKlqKiowNSpU3HPPfd0ux8nSRIEQUBBQQFmz54NQRAQExOD9evXY9OmTZg/f76sL4SoJ1ePiMd3ExtlJ0nAlm+rYW7nytRERMHMqzC0bds2nD59GrNmzUJtbS0qKirw9NNPIy4uDsOHD+9WftGiRcjMzMRTTz2F4uJiqFQqtLS0wGw245FHHpHtRRBdSJJRh0uSfDN2CACsdhH/OlTF8UNEREGsT/sW/O1vf0NRUZHr661bt+LNN990KyMIArZu3QoAWLVqFfbt24eKigp89NFH2L59O4YOHYqSkpK+t5zIQ1dmxqGktgVOHwWWKnMHdp86h6uGx/ukfiIi8q0+hSGdTocjR45g9OjRAIAZM2Zg+fLlFyzvSXgi8pVovQZj000+3Xl+35kGDIrVIyNW77NrEBGRb/RpNtmGDRswe/ZsPPjgg9iwYQOWLl3qWnG6J13hqcuMGTNw+PDhvlyaqE+uzIzzyQauXSQJ+PRQFWotHT67BhER+UafwtBll12GgoICTJ06FWfOnMHgwYPxr3/964LlvQ1PRHLTqpS4dmSiT6/RZnPi3YKzKPfR7DUiIvINn+1ar1AoIIrfr/5rs9nw4Ycf4tChQ4iNjcW9996LhIQEX1xaVlyBOrR8fLASJ2pafHoNlUJATmYcRqUYEenD3igiIrowb35++ywMhQqGodDSanVg4+4zPtmm44cUgoAh8XpcMSwOSUadz69HRETf8/l2HETBKlKrwpWZAzPrS5QknKprxdt7y7D5cBWa2mwDcl0iIvIO+/Ap7IxLM+HbSjNqLdYBuZ4kAUermnG0qvm7dY+icFmqCREa5YBcn4iILo49QxR2FAoB141M9NnK1BdTY+nAlyfq8epXp/H58TquXk1EFADYM0RhKcUUgTGpJhyqMPvl+jaHiMLSRhSWNiI2UoNBcXpkxOiRHhMBnZo9RkREA4lhiMLW1SPiceZcK5o7HH5tR0OrDQ2tNhSVNUEQgEGxeuSOToJRp/Zru4iIwgVvk1HY0qmVmD4qyd/NcCNJQOm5Nvz1m1IcrbKAkz2JiHyPPUMU1obGR+KyVCO+rbT4uylurHYRmw9XY9eJegxLiMSQ+EgkGrQwsLeIiEh2DEMU9qZdkoCKpnY0tQXeYOYWqwMHz5px8Gzn2Ca9RolLkgyYOCgGJj2DERGRHHibjMKeTq3E7RPTYdAF/u8GbTYnisqb8NrXp/HPf1dy6w8iIhkwDBEBMEWoceekdJ9u5ionSQJKalvwXsFZvLH7DOqaB2bNJCKiUMQwRPSdaL0Gt01Mg1YdXH8tzrXY8H7hWZxrYSAiIuqL4PpXn8jH4qO0uDkrFSqFH1Zk7Ic2mxN/LzyLhlZu+UFE5C2GIaIfSI/R4/oxyX5Zobo/Wq1OvLWnFF+V1MPqcPq7OUREQYNhiKgHI5IMuG1CGiK1wbUatN0pYe/pBrz21RkcrQqs5QKIiAIVwxDRBQyOi8S9VwzG4Di9v5vitXabE5sPV+OTg1XosLOXiIjoYhiGiC5Cr1Hh1vFpGJdu8ndT+uR4TTPe2H0GRwJsUUkiokDCMETUC4VCwPRRScjJjPN3U/qk1erElm+r8bf95TjbyHWJiIh+KDgWVSEKAFcMi0NGrB77TjfgdH2rv5vjtYrGdry7/yziDVqMT4/GJclR0KqCa0wUEZEvCBJ3grwoi8UCk8kEs9kMo9Ho7+ZQgKht7sCOY7WobOrwd1P6TK0UkJkQhbHpJqTHBN+4KCKii/Hm5zfDUC8YhuhCJEnCoQozvio5F/SDlJOMOkwaHINLkqIgBNuaAkREPWAYkhHDEPXG6nDi0FkzCssa0WoN7lCUFh2Ba0YmINGg83dTiIj6hWFIRgxD5CmHU0RReRP2nmmA1S76uzl9phAEjE41YsrQWJgi1P5uDhFRn3jz85sDqIlkolIqMHlILMakmbD71Dn8u7wJwfirhihJOFxhxtEqCy5LNSInMw56Df+pIKLQxan1RDLTqZW49tJE/OzyQUgwaP3dnD5zihIOnjXj9a/PoKC0EU4xCJMdEZEHeJusF7xNRv0hihIKyxrxzalzsDuD+6+aQafCuPRojEs3QafmlHwiCmze/PwOyp6hl156CUOGDIFOp0N2djb27t17wbLr16/H1KlTERMTg5iYGOTm5l60PJGcFAoBk4fE4r6cIciIDe7p680dDnxVUo8Nu05j35kG9hQRUcgIujC0adMm5OXlYeXKlSgsLERWVhZmzpyJ2traHsvv3LkTd999N3bs2IHdu3cjIyMDM2bMQEVFxQC3nMKZKUKN2yekYfygaH83pd9sDhG7TtTjr9+UoqS2xd/NISLqt6C7TZadnY3LL78cf/zjHwEAoigiIyMDDz74IB599NFez3c6nYiJicEf//hHzJ07t9vzVqsVVqvV9bXFYkFGRgZvk5FsDleYkX+0FmJw/dW7oASDFjmZcchMiPJ3U4iIXEL2NpnNZkNBQQFyc3NdxxQKBXJzc7F7926P6mhra4PdbkdsbGyPz69evRomk8n1yMjIkKXtRF3GpJlw8/hUaFRB9dfvguqarfioqBJfHK9DkP1uRUQEIMjCUH19PZxOJ5KSktyOJyUlobq62qM6HnnkEaSmproFqvMtX74cZrPZ9SgvL+93u4l+aGh8JO6YmI4ITegMRC4obcRH/66E1RHcC08SUfgJqjDUX8888wzeeecdfPDBB9Dpel5hV6vVwmg0uj2IfCHZpAu5QHSqrhVv7ylDfYu198JERAEiqMJQfHw8lEolampq3I7X1NQgOTn5ouc+//zzeOaZZ7B161aMGzfOl80k8liCQYvbJ6SF1FT1xjY73tlbhiOVFn83hYjII0EVhjQaDSZNmoT8/HzXMVEUkZ+fj5ycnAue9+yzz+Kpp57C5s2bMXny5IFoKpHHEo063D4xLWTGEAGA3Slhy7fV2HakBg5n8G5NQkThIej+9c3Ly8P69euxceNGHD16FIsXL0Zrayvmz58PAJg7dy6WL1/uKv+73/0Ojz/+OF599VUMGTIE1dXVqK6uRksLpwRT4Egy6nBzVipUitDaMf5whRnv7CtHU5vN300hIrqgoAtDs2fPxvPPP48VK1Zg/PjxKCoqwubNm12DqsvKylBVVeUq//LLL8Nms+HOO+9ESkqK6/H888/76yUQ9SgjVo+bxqVAIYRWIKprtuLNPWUoqW32d1OIiHoUdOsMDTRux0ED7URNMz49VB0y6xCdb8KgaFwxLC6kxkgRUWAK2XWGiMLBiCQDbhqXDGWI3TIDgANlTdiw6zR2FtfC0mH3d3OIiAAwDBEFpOGJBtw0LgVqZegFIptDxIGyJrz+1RnsOFaLVqvD300iojDH22S94G0y8qe6Zis+PliJprbQ7UVRKwWkRkcgPUaPQbF6JBm1EEJs3BQRDTxvfn4zDPWCYYj8rcPuxLYjNWGzKWqUVoVhCZEYm25CoqHnxVGJiHrDMCQjhiEKFN9WmvH58TpY7eGzbs+QeD0mDYrFoDi9v5tCREHGm5/fqgFqExH102WpJmTE6rH5cDUqGtv93ZwBcaa+DWfq2xAfpUFWRjRMEWrYHCKckoQItRIRGiVi9BqolRz+SER9x56hXrBniAKNKEr45vQ57D3dAP7tBZQKAYkGLdJiIjA0PhJp0REcc0RE7BkiCmUKhYArM+OREdPZS9QS5rOxnKKEKnMHqswd2H+mEXqNElE6FRxOCaIkIVKrgilCDUmS0NhmR0uHAya9GgkGLQxaFRQKASqFgKjvykXrNSG5rAERXRh7hnrBniEKZO02J7YeqcapulZ/NyVkaFQKpEbrMCg2EsMTo2CKUPu7SUTUBxxALSOGIQoGJ2qa8VVJPRpDeAq+v6SYdBgcF4n4KA3iorSI0at5G44oCPA2GVGYGZFkQGZCFA6UN+Grkno4Rf6OI5euW3BdNCoFEgxaJBi0iNVrEG/QItWkY0AiCmIMQ0QhQqEQMGlwDNJjIvDxwSpY2tlL5As2h4iKxna3GX3xBi2mDInFJUlRDEVEQYi3yXrB22QUjDrsTnx5oh7fVpo542wAGXQqjEw2YlSKAXFRWn83hyisccyQjBiGKJhVmdux41gdaiwdvRcmWSUata5gpNewE55ooDEMyYhhiIKdJEkoqW3B7lPncK7F5u/mhB2NSoEJGdGYNCQGWpXS380hChscQE1ELoIguAZYH6myYPfJc2G/NtFAsjlE7DndgH+fNWN4YhRGJEYhI1bPtYyIAgjDEFGYUCgEjEkz4dJkAwpKG7HvdAMcnHU2YDrsThyuMONwhRl6jRKjUowYm2ZCTKTG300jCnu8TdYL3iajUNXYakP+sVqUN7T5uylhLSNWj6x0EzIToqBgbxGRbHibjIh6FROpwZ2T0nGk0oIvT9Shzeb0d5PCUnlDG8ob2qDXKDEiKQqXJhuRFh3h72YRhRWGIaIwNzrViGEJkfiqpB6HKywQ2VnsF202J/5dbsa/y81INukwZWgsMhOi/N0sorDA22S94G0yCif1LVbsOlGP0/Xc6ywQxBu0mJARjVEpRg64JvISp9bLiGGIwlF5Qxu+Kql324aC/CdSq8TIZCNGphiQaND5uzlEQYFhSEYMQxTOTta1YPfJc6hrtvq7KfSd2EgNBsfpMSQuEnFRGkRqVBx4TdQDDqAmIllkJkRhWHwkSmpb8M2pc6jnoo1+19BqQ0OrDQfKmgAACkFApFaJ2EgNYiM1SIuOwKA4PRd4JPICe4Z6wZ4hok6SJOFUfSv2nW7g7bMAp1QISDRoodeqoFMpoBAEOCUJktR5y82oUyNGr0GCQYsIDUMThSb2DBGR7ARBQGZCFDITolB2rg3fnD7ntnM7BQ6nKHkcWPUaJTQqBVRKBbrutikEAXqNElFaFVRKBax2J+xOCYLQ+ZxWrUCUVgWjTo30mAhEavmjhIIbP8FE5LVBcXoMitOjoqkd+8804HR9K9jHHJzabM5+rTElCECSUYfU6AgYdCoYtCoIQmeq0qkVSDbqoFIq5GoukU8wDBFRn6VFRyBtfBrqmq0oKm9CcbUFdidTUTiRJKDa3IHqC/REKRUCkoxaJBp0iI/SIiW68/9EgYRjhnrBMUNEnuuwO3Gowox/lzehuYObwVLPTBFqZCZGYVCsHmnREdCo2HNE8uPUehkxDBF5TxQlFNc046uSeoYiuiiFICAuSoNovRrRERrERXUO7I7Va7hkAPULB1ATkV8pFAJGpRiRmRCFPafP4UBZE5wif++i7kRJQl2ztdtaVhpV53ij9JgIDI2PRKKRi02S77BnqBfsGSLqvzabAydrW3GithlnG9sZjMhrxgg1MhMiMTQ+EmnRERyUTb3ibTIZMQwRyavD7sSpus5gVHqujcGIvNbVa5QaHYFkkw4xejWMOjVvq5Eb3iYjooClUysxOtWI0alGdNidKKltwcGzZtRYuJAjecbmEFHW0IayhjbXMaVCgFqpgCAAKoWABIMWCQYtYvQaGHQqmCLUMOjUfmw1BTKGISLyG51aiTFpJoxJM6GiqR1HKi1oaLWiodWODnvf176h8OMUJTjF7z8zzR0OnKprdSuj1yiRZNQh0aBFolGL+CgtdGol1EoFlOxVCmsMQ0QUENKiI5AWHeH6us3mwLkWG+parKhsakdFY3u/FgckarM5cbq+FafrW7s9p1MrYYpQf9eDpIIxQo0orQqRWiX0GhUi1Eq3JQBEsXNF7q4FJql3kiTBIUqwOUR0fLequVOSoFYKSDT4d4A8wxARBSS9RgV9rAoZsXpMHBQDoDMgmdvtMLfbXRuW1jdb0dRu5wrY1C8ddic67M6L3q5VKgQoFQIcTgnieR84QQAEdIYihdA5m7Lrtp1a2fV/BTQqBXQqBXRqJbQqBbRqJSK+e+i1nduf6NSBtVec3SmipcOBFqsDVocIxXcB0Opwot3mhM0hwilKECXAIYoQpc6wY+lwoKXDAbtThFOSIIqdQainv6dpMRG4a3LGwL+48zAMEVHQ0GtU0GtUSDFFuB23Opyu6dn1LTaUNbTB0m73UyspVHXeiuv+01ySAAmdx0Wp6z9AO7zvydSoFIjUKKFWKaBWdN6+Uyg6w5aEzjBhc4hotzvhcErQaZSI1CihUiogoDOYOc8LHQqFAIUAWO0iOhyd4cXxXY+MVqVAxHe3CfHduSqlAmqFAJtTxLkWGywd4fGLBsMQEQU9rUqJ9Bg90mP0rmNV5naU1LagxmJFfYsV7bzFRkHA5hBhc4gel2+xOlDfx2u125xoAn9pABiGiChEpZgi3HqQzO121Fg699CqMrej1mKFg9P6iQgMQ0QUJroGx16SZAAAOJwiKps6UFLXjFN1rdw2hCiMMQwRUVhSKRUYFKfHoDg9rhsJ1LdYUXquDeUNbTjb2Aa7k71GROGCYYiICEB8VOe6M5MGx8ApSqhsasfZxnacbWxDjaWD4YgohDEMERH9gFIhICNWj4xYPYA4SJKEpjY76lqsONdiQ2ObDedabWhstXE7EaIQwDBERNQLQRAQE6lBTKQGSPr+uChKaGyzobbZihpLB+qarWhqs6PFyvFHRMGEYYiIqI8UCgFxUVrERWkxKuX7jSCtDifMbZ2LQza122Fu6/x/S4cdrd8tVEdEgYNhiIhIZlqVEolGJRKNPW8xYHeKaLV2rurbanV+938H2mwOtFidaLc70W5zwGoXOf2faAAwDBERDTC1UoFovQbRek2vZbv2c7I7RVjtIqzfLcpnczphc3QetztFdNg7VyW2Opywfvfn5g4HN7wl8gDDEBFRABMEwbW/lQfZqZsOu/O7Xiena3NMh9gZqtptXb1QTrR+1xPlFL/fTFMMh30YiMAwREQU0nRqJXRqJeL6cK7N0bmfVbvN6QpOHXYnrI7O3ijHd8sNCAIgShJarZ3l2mwOdDg6e7IYqCgYMAwREVGPNKrOndaNOnWf6+i6jSdKnbf8RBGwi53HbI6u/3f2VnUFrK5d0M/fAV4hCBAEwOGU4OjaCV3qDGHidxuoOsXODUg7z++8VtfXnTumf1+3JHVuZtpZT+efFULnzvNOUbzgDusUmhiGiIjIZ9RKhWtX9GDTFeTsDgkdjq7bjCKcYmcI69od3il1Bi1R+n5ne6cowS6KcDo7g5h4fkiTJFcdoti5370kwbUrPdDZ29YVzoTz2tRT2a7jTmdniOsKf111d9WhUHTW1BUEJcAVBMNdUIahl156Cc899xyqq6uRlZWFP/zhD5gyZcoFy7/77rt4/PHHcebMGYwYMQK/+93vcOONNw5gi4mIKNi4gpwGMKHvvWOB7vwg5/iuV8zp7Apu3wcvpyjB4ewMeTbH94P5HU4RNqfo1hvnCpLf9fQ5zgtqXbdOu0KYVuX/sBx0YWjTpk3Iy8vDunXrkJ2djbVr12LmzJkoLi5GYmJit/Jff/017r77bqxevRo/+clP8NZbb+HWW29FYWEhxowZ44dXQEREFDgEQYBS6Fx5XQP/BxN/ECQpuDrIsrOzcfnll+OPf/wjAEAURWRkZODBBx/Eo48+2q387Nmz0draio8//th17IorrsD48eOxbt26Xq9nsVhgMplgNpthNBp7LU9ERET+583P76CKgDabDQUFBcjNzXUdUygUyM3Nxe7du3s8Z/fu3W7lAWDmzJkXLG+1WmGxWNweREREFLqCKgzV19fD6XQiKSnJ7XhSUhKqq6t7PKe6utqr8qtXr4bJZHI9MjIy5Gk8ERERBaSgCkMDYfny5TCbza5HeXm5v5tEREREPhRUA6jj4+OhVCpRU1PjdrympgbJyck9npOcnOxVea1WC61WK0+DiYiIKOAFVc+QRqPBpEmTkJ+f7zomiiLy8/ORk5PT4zk5OTlu5QFg27ZtFyxPRERE4SWoeoYAIC8vD/PmzcPkyZMxZcoUrF27Fq2trZg/fz4AYO7cuUhLS8Pq1asBAEuWLMGPfvQjvPDCC7jpppvwzjvvYP/+/XjllVf8+TKIiIgoQARdGJo9ezbq6uqwYsUKVFdXY/z48di8ebNrkHRZWRkUiu87vK688kq89dZb+O1vf4vf/OY3GDFiBD788EOuMUREREQAgnCdoYHGdYaIiIiCT8iuM0REREQkN4YhIiIiCmsMQ0RERBTWGIaIiIgorDEMERERUVhjGCIiIqKwFnTrDA20rpUHuHs9ERFR8Oj6ue3JCkIMQ71obm4GAO5eT0REFISam5thMpkuWoaLLvZCFEVUVlbCYDBAEIRuz1ssFmRkZKC8vDxsF2UM9/cg3F8/wPcA4HsQ7q8f4HsABNZ7IEkSmpubkZqa6rYzRU/YM9QLhUKB9PT0XssZjUa/f+P9Ldzfg3B//QDfA4DvQbi/foDvARA470FvPUJdOICaiIiIwhrDEBEREYU1hqF+0mq1WLlyJbRarb+b4jfh/h6E++sH+B4AfA/C/fUDfA+A4H0POICaiIiIwhp7hoiIiCisMQwRERFRWGMYIiIiorDGMERERERhjWGon1566SUMGTIEOp0O2dnZ2Lt3r7+b5BOrV6/G5ZdfDoPBgMTERNx6660oLi52K3PNNddAEAS3x6JFi/zUYvn9z//8T7fXN3LkSNfzHR0deOCBBxAXF4eoqCjccccdqKmp8WOL5TVkyJBur18QBDzwwAMAQvP7/8UXX2DWrFlITU2FIAj48MMP3Z6XJAkrVqxASkoKIiIikJubixMnTriVaWhowJw5c2A0GhEdHY0FCxagpaVlAF9F/1zsPbDb7XjkkUcwduxYREZGIjU1FXPnzkVlZaVbHT19dp555pkBfiV919vn4Oc//3m313f99de7lQnmz0Fvr7+nfxcEQcBzzz3nKhPonwGGoX7YtGkT8vLysHLlShQWFiIrKwszZ85EbW2tv5smu88//xwPPPAAvvnmG2zbtg12ux0zZsxAa2urW7mFCxeiqqrK9Xj22Wf91GLfuOyyy9xe365du1zPPfTQQ/jnP/+Jd999F59//jkqKytx++23+7G18tq3b5/ba9+2bRsA4Kc//amrTKh9/1tbW5GVlYWXXnqpx+efffZZ/P73v8e6deuwZ88eREZGYubMmejo6HCVmTNnDr799lts27YNH3/8Mb744gv84he/GKiX0G8Xew/a2tpQWFiIxx9/HIWFhXj//fdRXFyMm2++uVvZJ5980u2z8eCDDw5E82XR2+cAAK6//nq31/f222+7PR/Mn4PeXv/5r7uqqgqvvvoqBEHAHXfc4VYuoD8DEvXZlClTpAceeMD1tdPplFJTU6XVq1f7sVUDo7a2VgIgff75565jP/rRj6QlS5b4r1E+tnLlSikrK6vH55qamiS1Wi29++67rmNHjx6VAEi7d+8eoBYOrCVLlkiZmZmSKIqSJIX+9x+A9MEHH7i+FkVRSk5Olp577jnXsaamJkmr1Upvv/22JEmSdOTIEQmAtG/fPleZf/3rX5IgCFJFRcWAtV0uP3wPerJ3714JgFRaWuo6NnjwYOnFF1/0beMGSE/vwbx586RbbrnlgueE0ufAk8/ALbfcIl133XVuxwL9M8CeoT6y2WwoKChAbm6u65hCoUBubi52797tx5YNDLPZDACIjY11O/7mm28iPj4eY8aMwfLly9HW1uaP5vnMiRMnkJqaimHDhmHOnDkoKysDABQUFMBut7t9HkaOHIlBgwaF5OfBZrPhr3/9K/7jP/7DbQPjUP/+n+/06dOorq52+56bTCZkZ2e7vue7d+9GdHQ0Jk+e7CqTm5sLhUKBPXv2DHibB4LZbIYgCIiOjnY7/swzzyAuLg4TJkzAc889B4fD4Z8G+sjOnTuRmJiISy+9FIsXL8a5c+dcz4XT56CmpgaffPIJFixY0O25QP4McKPWPqqvr4fT6URSUpLb8aSkJBw7dsxPrRoYoihi6dKluOqqqzBmzBjX8XvuuQeDBw9GamoqDh48iEceeQTFxcV4//33/dha+WRnZ+P111/HpZdeiqqqKjzxxBOYOnUqDh8+jOrqamg0mm4/AJKSklBdXe2fBvvQhx9+iKamJvz85z93HQv17/8PdX1fe/o3oOu56upqJCYmuj2vUqkQGxsbkp+Ljo4OPPLII7j77rvdNun85S9/iYkTJyI2NhZff/01li9fjqqqKqxZs8aPrZXP9ddfj9tvvx1Dhw7FyZMn8Zvf/AY33HADdu/eDaVSGVafg40bN8JgMHQbIhDonwGGIfLaAw88gMOHD7uNlwHgdv977NixSElJwfTp03Hy5ElkZmYOdDNld8MNN7j+PG7cOGRnZ2Pw4MH429/+hoiICD+2bOBt2LABN9xwA1JTU13HQv37Txdnt9tx1113QZIkvPzyy27P5eXluf48btw4aDQa/Od//idWr14ddNs29ORnP/uZ689jx47FuHHjkJmZiZ07d2L69Ol+bNnAe/XVVzFnzhzodDq344H+GeBtsj6Kj4+HUqnsNluopqYGycnJfmqV7/33f/83Pv74Y+zYsQPp6ekXLZudnQ0AKCkpGYimDbjo6GhccsklKCkpQXJyMmw2G5qamtzKhOLnobS0FNu3b8f9999/0XKh/v3v+r5e7N+A5OTkbhMqHA4HGhoaQupz0RWESktLsW3bNrdeoZ5kZ2fD4XDgzJkzA9PAATZs2DDEx8e7Pvvh8jn48ssvUVxc3Ou/DUDgfQYYhvpIo9Fg0qRJyM/Pdx0TRRH5+fnIycnxY8t8Q5Ik/Pd//zc++OADfPbZZxg6dGiv5xQVFQEAUlJSfNw6/2hpacHJkyeRkpKCSZMmQa1Wu30eiouLUVZWFnKfh9deew2JiYm46aabLlou1L//Q4cORXJystv33GKxYM+ePa7veU5ODpqamlBQUOAq89lnn0EURVdYDHZdQejEiRPYvn074uLiej2nqKgICoWi262jUHH27FmcO3fO9dkPh88B0NljPGnSJGRlZfVaNuA+A/4ewR3M3nnnHUmr1Uqvv/66dOTIEekXv/iFFB0dLVVXV/u7abJbvHixZDKZpJ07d0pVVVWuR1tbmyRJklRSUiI9+eST0v79+6XTp09L//jHP6Rhw4ZJ06ZN83PL5fOrX/1K2rlzp3T69Gnpq6++knJzc6X4+HiptrZWkiRJWrRokTRo0CDps88+k/bv3y/l5ORIOTk5fm61vJxOpzRo0CDpkUcecTseqt//5uZm6cCBA9KBAwckANKaNWukAwcOuGZKPfPMM1J0dLT0j3/8Qzp48KB0yy23SEOHDpXa29tddVx//fXShAkTpD179ki7du2SRowYId19993+ekleu9h7YLPZpJtvvllKT0+XioqK3P5tsFqtkiRJ0tdffy29+OKLUlFRkXTy5Enpr3/9q5SQkCDNnTvXz6/Mcxd7D5qbm6Vly5ZJu3fvlk6fPi1t375dmjhxojRixAipo6PDVUcwfw56+3sgSZJkNpslvV4vvfzyy93OD4bPAMNQP/3hD3+QBg0aJGk0GmnKlCnSN9984+8m+QSAHh+vvfaaJEmSVFZWJk2bNk2KjY2VtFqtNHz4cOnhhx+WzGazfxsuo9mzZ0spKSmSRqOR0tLSpNmzZ0slJSWu59vb26X/+q//kmJiYiS9Xi/ddtttUlVVlR9bLL8tW7ZIAKTi4mK346H6/d+xY0ePn/t58+ZJktQ5vf7xxx+XkpKSJK1WK02fPr3be3Pu3Dnp7rvvlqKioiSj0SjNnz9fam5u9sOr6ZuLvQenT5++4L8NO3bskCRJkgoKCqTs7GzJZDJJOp1OGjVqlLRq1Sq3oBDoLvYetLW1STNmzJASEhIktVotDR48WFq4cGG3X4qD+XPQ298DSZKkP//5z1JERITU1NTU7fxg+AwIkiRJPu16IiIiIgpgHDNEREREYY1hiIiIiMIawxARERGFNYYhIiIiCmsMQ0RERBTWGIaIiIgorDEMERERUVhjGCIiIqKwxjBEREREYY1hiIiIiMIawxARBQVJkrBmzRoMHToUer0et956K8xmc49lr7nmGgiCAEEQUFRUdNF6r7nmGixdulTWtv785z93Xf/DDz+UtW4ikh/DEBEFhYcffhgvv/wyNm7ciC+//BIFBQX4n//5nwuWX7hwIaqqqjBmzJiBa+R3/u///g9VVVUDfl0i6huGISIKeHv27MGaNWuwadMmTJs2DZMmTcLChQvx6aefXvAcvV6P5ORkqFSqAWxpJ5PJhOTk5AG/LhH1DcMQEQW8559/HtOnT8fEiRNdx5KSklBfX+9VPa2trZg7dy6ioqKQkpKCF154oVsZURSxevVqDB06FBEREcjKysJ7773ner65uRlz5sxBZGQkUlJS8OKLL/rkVhsRDRyGISIKaFarFZ988gluu+02t+MdHR0wmUxe1fXwww/j888/xz/+8Q9s3boVO3fuRGFhoVuZ1atX44033sC6devw7bff4qGHHsK9996Lzz//HACQl5eHr776Ch999BG2bduGL7/8slsdRBRcBr7/mIjIC4WFhWhvb8evfvUr/PrXv3Ydt9vtuPbaaz2up6WlBRs2bMBf//pXTJ8+HQCwceNGpKenu8pYrVasWrUK27dvR05ODgBg2LBh2LVrF/785z9j4sSJ2LhxI9566y1XHa+99hpSU1PleKlE5CcMQ0QU0I4fP47IyMhus8JuuukmXHXVVR7Xc/LkSdhsNmRnZ7uOxcbG4tJLL3V9XVJSgra2Nvz4xz92O9dms2HChAk4deoU7HY7pkyZ4nrOZDK51UFEwYdhiIgCmsViQXx8PIYPH+46VlpaihMnTuCOO+6Q9VotLS0AgE8++QRpaWluz2m1WjQ0NMh6PSIKDBwzREQBLT4+HmazGZIkuY797//+L2688UaMHj3a43oyMzOhVquxZ88e17HGxkYcP37c9fXo0aOh1WpRVlaG4cOHuz0yMjIwbNgwqNVq7Nu3z3WO2Wx2q4OIgg97hogooF133XXo6OjAM888g5/97Gd488038c9//hN79+71qp6oqCgsWLAADz/8MOLi4pCYmIjHHnsMCsX3vxMaDAYsW7YMDz30EERRxNVXXw2z2YyvvvoKRqMR8+bNw7x58/Dwww8jNjYWiYmJWLlyJRQKBQRBkPulE9EAYRgiooCWlJSE119/HQ8//DCeeuopXHfdddi1axcyMjK8ruu5555DS0sLZs2aBYPBgF/96lfdVrF+6qmnkJCQgNWrV+PUqVOIjo7GxIkT8Zvf/AYAsGbNGixatAg/+clPYDQa8etf/xrl5eXQ6XSyvF4iGniCdH7fMxFRCLjmmmswfvx4rF271ufXam1tRVpaGl544QUsWLDA7TlBEPDBBx/g1ltv9Xk7iKjvOGaIiELSn/70J0RFReHQoUOy1nvgwAG8/fbbOHnyJAoLCzFnzhwAwC233OIqs2jRIkRFRcl6XSLyHfYMEVHIqaioQHt7OwBg0KBB0Gg0stV94MAB3H///SguLoZGo8GkSZOwZs0ajB071lWmtrYWFosFAJCSkoLIyEjZrk9E8mMYIiIiorDG22REREQU1hiGiIiIKKwxDBEREVFYYxgiIiKisMYwRERERGGNYYiIiIjCGsMQERERhTWGISIiIgprDENEREQU1hiGiIiIKKz9f/nsiM3WHhyOAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.fill_between(\n", " np.rad2deg(solver.angles), *np.percentile(xs_ratio, [5, 95], axis=0), alpha=0.5\n", ")\n", "plt.xlabel(r\"$\\theta$ [deg]\")\n", "plt.ylabel(\n", " r\"$\\frac{d \\sigma}{d\\Omega} / \\frac{d \\sigma_{\\text{R}}}{d\\Omega}$ [dimensionless]\"\n", ")" ] }, { "cell_type": "code", "execution_count": 14, "id": "5a19ce84-5297-4b4f-a759-db23f8d1f597", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Text(0, 0.5, '$A_y$ [dimensionless]')" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkgAAAG0CAYAAADJpthQAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAAYmxJREFUeJzt3Xl8VOW9P/DPmZnMTJaZyb6HJQFZFAIEiFFBK1EQq7gWlRbkh3i1aquoVdpbl9rbeFtFblsr1evaYvW6W2sRRBHRyBYieyBsWSd7ZpKZzH5+f1CmJpAwmZyZOWfyeb9e83rJmXOe853JmPnmeb7P8wiiKIogIiIiIj9VpAMgIiIikhsmSERERER9MEEiIiIi6oMJEhEREVEfTJCIiIiI+mCCRERERNQHEyQiIiKiPjSRDkCJfD4fGhoaYDAYIAhCpMMhIiKiAIiiiK6uLmRnZ0OlGriPiAlSEBoaGpCXlxfpMIiIiCgItbW1yM3NHfAcJkhBMBgMAE6+wUajMcLREBERUSCsVivy8vL83+MDYYIUhFPDakajkQkSERGRwgRSHsMibSIiIqI+mCARERER9cEEiYiIiKgPJkhEREREfTBBIiIiIuqDCRIRERFRH0yQiIiIiPpggkRERETUBxMkIiIioj6YIBERERH1wQSJiIiIqA8mSERERER9MEEiIiIi6iMqEqRnn30Wo0aNgl6vR3FxMbZt2xbQdW+88QYEQcA111wT2gCJiIgUyOcTIx1CxCg+QXrzzTexYsUKPProo6ioqEBhYSHmzp2L5ubmAa87fvw4HnjgAcyaNStMkRIREclLc5cDDre31zG314eKmg68ub0G/7vlKOo7eyIUXWQpPkFatWoVli9fjqVLl2LixIlYs2YN4uLi8NJLL/V7jdfrxaJFi/D4448jPz8/jNESERHJx/ZjHSg/0tbr2JbDrfiiqgUNnQ7YnF68s7MOu+s6IxNgBCk6QXK5XNi5cydKS0v9x1QqFUpLS1FeXt7vdb/61a+Qnp6OZcuWBXQfp9MJq9Xa60FERKRkVocb1c3d2F1nQUuXEwBQ12HHt32SIa9PxMYDzahps0cgyshRdILU2toKr9eLjIyMXsczMjJgNpvPeM2WLVvw4osv4oUXXgj4PmVlZTCZTP5HXl7ekOImIiKKtMqaTvhEET5RxKaqZri9Pny6vwliP2VHXxxqHlY1SYpOkAarq6sLP/rRj/DCCy8gNTU14OtWrlwJi8Xif9TW1oYwSiIiotByeXzY22Dx/7uuowfv7KxDh93d7zWt3S7srrf0+3y00UQ6gKFITU2FWq1GU1NTr+NNTU3IzMw87fwjR47g+PHjuOqqq/zHfD4fAECj0aCqqgoFBQWnXafT6aDT6SSOnoiIKDL2N1rhdPt6HWu0OM56XfmRNozLMCBWqw5VaLKh6B4krVaLoqIibNy40X/M5/Nh48aNKCkpOe388ePHY8+ePaisrPQ/rr76anzve99DZWUlh86IiGhY2BNk0bXD7cWGA03wDoOhNkX3IAHAihUrsGTJEkyfPh0zZ87E6tWrYbPZsHTpUgDA4sWLkZOTg7KyMuj1epx33nm9rk9MTASA044TERFFo0ZLD1q7XUFff6S5Gx9+W4/vT85GjFrR/SwDUnyCtHDhQrS0tOCRRx6B2WzGlClTsG7dOn/hdk1NDVSq6P0BEhERDca++qHPxD7ease7FXW4flouNFGaJAmi2F+9OvXHarXCZDLBYrHAaDRGOhwiIqKAuL0+PL/5KFwe39lPDsCssamYPipZkrbCYTDf39GZ9hEREdFpDjd1S5YcAcC24+2nrcQdLZggERERDRP7GqSdpu90+/DN0bazn6hATJCIiIiGAYvdHZJ91XbXWWAZYP0kpWKCRERENAxUt3T3u0r2UHh9IiqjcK82JkhERETDwLFWW8jabrKefZFJpWGCREREFOUcbi/qO6QfXjulpcuJaJsUzwSJiIgoyh1vs8EXwgTG5fHB0hNddUhMkIiIiKLcsZbQDa+d0tzlDPk9wokJEhERURTz+UQcawtDgmRlgkREREQKUd/ZA6dbusUh+9PcFV2F2kyQiIiIotjREM5e+64WDrERERGRUlQ3d4flPnaXF12O6CnUZoJEREQUpZqsDljDOLssmgq1mSARERFFqUNNXWG9XzQVajNBIiIiilKHm8IzvHZKNBVqM0EiIiKKQs1WR9gXb4ymQm0mSERERFHoUJh7jwCgy+FBj8sb9vuGAhMkIiKiKHS4Obz1R6dESy8SEyQiIqIo09LlRKc9MlPuW7qZIBEREZEM1bTbI3Zv9iARERGRLNV39kTs3uxBIiIiItkRRRENEUyQOmwueH1ixO4vFSZIREREUaS12xXRmWRen4g2m/J7kZggERERRZFIDq+dEg11SEyQiIiIokhdR+QKtE9hgkRERESyUt8R+R6k1m5XpEMYMiZIREREUaKt2wm7DFayZg8SERERyYYc6o8AwOH2wuqIzEKVUmGCREREFCXkMLx2SqvCe5GYIBEREUWJJqsj0iH4KX2YLSoSpGeffRajRo2CXq9HcXExtm3b1u+57777LqZPn47ExETEx8djypQp+Mtf/hLGaImIiKTn8fpg6fFEOgy/NpuyC7UVnyC9+eabWLFiBR599FFUVFSgsLAQc+fORXNz8xnPT05Oxi9+8QuUl5dj9+7dWLp0KZYuXYpPPvkkzJETERFJp93mgk+UzwrWSk+QBFGU0bsZhOLiYsyYMQN//OMfAQA+nw95eXm455578PDDDwfUxrRp03DllVfiiSeeCOh8q9UKk8kEi8UCo9EYdOxERERS2d9gxSf7zJEOw0+jEnDX98ZApRIiHYrfYL6/Fd2D5HK5sHPnTpSWlvqPqVQqlJaWory8/KzXi6KIjRs3oqqqCrNnz+73PKfTCavV2utBRETK1mR14HBTF/Y1WKJi7zC5be/h8Ymw9Ch3Jpsm0gEMRWtrK7xeLzIyMnodz8jIwMGDB/u9zmKxICcnB06nE2q1Gn/6059w2WWX9Xt+WVkZHn/8ccniJiKiyNpxvB1fHm71/7u5y4nvjUuPYERD1ybDxRnbbC4kxWsjHUZQFN2DFCyDwYDKykps374d//Vf/4UVK1Zg06ZN/Z6/cuVKWCwW/6O2tjZ8wRIRkWREUcTmQy29kiMAqKzpRHVzd4SikkZrt7x6kICTdVFKpegepNTUVKjVajQ1NfU63tTUhMzMzH6vU6lUGDNmDABgypQpOHDgAMrKynDJJZec8XydTgedTidZ3EREFBnf1lmw80THGZ/bsL8J6UYdjPqYMEc1dE6PF10O+cxgO6VdZsN+g6HoHiStVouioiJs3LjRf8zn82Hjxo0oKSkJuB2fzwenU7k/RCIiOrtupwdfVbf2+7zD7cX7u+phc8ov0TgbOQ6vAcqeyaboHiQAWLFiBZYsWYLp06dj5syZWL16NWw2G5YuXQoAWLx4MXJyclBWVgbgZD3R9OnTUVBQAKfTiY8//hh/+ctf8Nxzz0XyZRARUYh9UdUCl8c34Dlt3S68vbMO1xflIkGnnK9IuSZIHTYXRFGEIMhnJluglPPT78fChQvR0tKCRx55BGazGVOmTMG6dev8hds1NTVQqf7dUWaz2fDjH/8YdXV1iI2Nxfjx4/HXv/4VCxcujNRLICKiEDvWasOhpq6Azm23ufD2jlpcOzUXpjhlDLe1ynQoy+09OZMtMU55hdqKXwcpErgOEhGRcvS4vFi79cSga3RitWp8f3IWcpPiQhSZdN7eWYfadnukwzijq6dkoyAtIdJhABhG6yARERENRBRF/HNvY1AFzD0uL96tqEd1c2A9T5HUJsMZbKcodSYbEyQiIopa3xxtx4m24HtWvD4Rm6pa4PEOXLsUSTanB3aXN9Jh9Euu9VFnwwSJiIiiUm27HVuPtQ25nS6HB7vrLRJEFBpy76GRe3z9YYJERERRx+P1YeOBJkhVZbvjePtZZ8BFitUh7+08OuwnZ7IpDRMkIiKKOtuOtaPDLl3iYHN68W1dp2TtSUmOC0R+l8vjQ7cC15ZigkRERFGltduJHf2slj0UO453yLK3xqqADWGtMk/izoQJEhERRQ2P14cN+5vg9Uk/pONwe/HuzjrZrbQt9x4kAOiSYWJ5NkyQiIgoanx6oBlmiyNk7XfY3Xi3og49Mpo1poTkw9oj/ySuLyZIREQUFXYcb8eBRmvI79Pa7cJbO2tlMztLCT1IShgG7IsJEhERKZ7Z4sCWATailVpbtwt/21aDKnNkF5G0uzzwhGA4UWpdTiZIREQUIi6PD3vqLIoYUgm3A2arZFP6A+Xy+PDxnkbsPNEe3ht/h1KGrpQS53cpfrNaIqLhwOH24oPKejR0OiAIwIjkOFw7NUeRu6SHwpHm7ojde/OhVqhVKkzJSwz7vZWSLCslzu9igkREJAMdNheau5xoszmRGKtFbnIsjPoYOD1edDk8+OdeM1q7Tu63JYrAiTY7Djd345wMQ4QjjzyzxRHxOpxNVc3QaVSYkBXeDcyVMn3e7RVhd3kQp1VO2qGcSImIooDXJ6LN5kS6Qe8/VtNmx7u76k4bItKohAHrS7Ydaw86QRJFEQ63D4IA6GPUQbUhF9UR7D06RRSBr4+0YXymIay9ekrqmbH2KCtBYg0SEVEYVdZ24P1d9bD8a1ZPj8uLT/aZz1g/c7bi25YuJ462DD452HyoBb/fWI01XxzBB5X1itwG4ruqmyNbKH2KtceNo622sN4z0j1ngyHHRTYHwgSJiChMHG4vth3rgM3pxYeV9XC4vVi/3zykbRi2Hx9cgfCWw63YeaIDvn8lRQ2dDuwPw9T4UGnrdkq6pchQ7Q7zdiRKSjqUNtVfOX1dREQKt/VYOxzukwsMtna7sHZrzZC/NBo6Hdh2rB15ybFIidfB6fGirduFLocHGrUArebk38Fen4j6zh5U1nSe1sZX1a0Yk54AnUYNURRh6XEjMU47pLjCRQ7Da991os2ODpsLSfHhef+U1IOkpFgBJkhERGFhsbvxbW1nr2NS/UX91RDX/7E5vSg/0oaxGQZsPtSC1i4nrirMxqjUeEniCxW7y4O9DfLq/RJF4Nu6TlwyLj3k93J7fbJa0ftslNTbBTBBIiIKi23H20OyP5hUKms7ses7vUt//7Yh7ElSl8ONdpsLKkGA0+NDc5cDrd0uxMWokZKgxYjkOKQk6ACcHK58p6JelsM2+xutuGhMKjTq0FaxKK1HRo4/q4EwQSIiCrFupycsW2AMRd86bY9PxN+/bcDCmXm9ZtyF0peHWwdcmVoQgIK0BEzJS8SXh1v9yx7IjdPtQ31nD0amhDa5VNIMNkA5SxKcwiJtIqIQ21XTIeveo/54fCI+2dcUlthtTs9Z64lE8WTN0ds769BkDd2GtFKoabeH/B5KW53a5fH5a/CUgAkSEVEIOT1e7K6zRDqMoLV2OVF+pC3k99lbb1FkEtmfcCRISutBApQ1zMYEiYgohPbUWeDy+CIdxpDsONGOhs6ekLXv84nYU6/cJPJMWrqcsLtC28OjtCErQFkxM0EiIgqhyj4z15RIFIFNVS0hW1DyaKtNcQXHZyOKQG176JJKAOiwu0LafigoaSYbEyQiohBp7XZGzRd/k9WBwyFacyjciyuGy4m20K2qfazVBrNF3nVYZ8IhNiIiQn1HaHsQwq38SBt8EtcJ1bTZcaIt9PU6kRCqOiSvT8QXVc0haTvUlPQHAxMkIqIQqYuyBKnd5pJ0WxKvT8SmQ8r8og9El8ODDpv0w2CVtR2y2l5lMDjERkREqOuIvp6Rb462DWnvuO+qrO1AW7fy6mgG44TEvUg2pwffHB3c/ntyoqSlCZggERGFQFu3E3YFbQMRqC6HB69+fRyVtZ1DKtruVvgXfaAOSrxAaGVtp6JnRTrcXsXEz5W0iYhCINqG177L5fHh84PNqDjRAX2MGgBw2cQMpBl0Abex8UCTYr4oh6LR4kBtux15yXFDbsvrE7GvQfnLIVgdbqQmBP5ZiZSo6EF69tlnMWrUKOj1ehQXF2Pbtm39nvvCCy9g1qxZSEpKQlJSEkpLSwc8n4goGNGcIJ1i6XGjyepAk9WBbccC7w2qqOnA0ZbQzfCSm+3HpekpO9rSDZtT+b2SSinUVnyC9Oabb2LFihV49NFHUVFRgcLCQsydOxfNzWcu/Nu0aRNuvvlmfP755ygvL0deXh4uv/xy1NfXhzlyIopm9Z3RV380kMPNXWgPoCC52erAlsOtYYhIPk602SXZGkXJK7J/l1Km+is+QVq1ahWWL1+OpUuXYuLEiVizZg3i4uLw0ksvnfH8tWvX4sc//jGmTJmC8ePH43//93/h8/mwcePGMEdORP3ZU2cJ6crNodZuc0XFX/qDIYr995SIooi2bif21lvwjz2NUbWlSKAG08N2Jh02F2qjpOhfKTPZFF2D5HK5sHPnTqxcudJ/TKVSobS0FOXl5QG1Ybfb4Xa7kZyc3O85TqcTTue/d422WuW9KzeRkjVZHfi8qhmiCBTnJ6NoZBI6bC602VzITYqFQR8T6RDPKhpnrwXiYGMXzs9PgSn23z8jm9ODd3fVo7XLOcCV0e9ISzdau51B197sqbcgRAuZh51ShtgUnSC1trbC6/UiIyOj1/GMjAwcPHgwoDYeeughZGdno7S0tN9zysrK8Pjjjw8pViI6O7fXh3V7zf4ehvIjbfjmaJv/i0GtEnBORgIm5SYiy6iHSiUAONlD4fL6oNOoIxV6L0ru/RoKnyji0/1NuHBMKjJNelgdbry7s06xa/ZISRSBr6pbsWBKzqCvtfS4o2qvOqUMsSk6QRqqJ598Em+88QY2bdoEvV7f73krV67EihUr/P+2Wq3Iy8sLR4hEw8qXh1tOq2P57l/NXp+IA41dONDYBV2MCrlJcXC6vWjucsLrEzE6NR7n5ZgwOjU+zJH3NhwKtPtT025HzbYapBp0cLq9iuktCIejLbZBz2jz+USs29sYVTP+lPKZUHSClJqaCrVajaampl7Hm5qakJmZOeC1Tz31FJ588kl8+umnmDx58oDn6nQ66HTyn5JIpFQ+n4jPq5oHVYTqdPtwpM/eYNXN3ahu7kZJQQrOz0+ROsyAWB1uxXwBhNJwH1Lrz5bqVtw8c0TA539zrA0Nncrbc20gNpcHHq8PGrW8y6DlHd1ZaLVaFBUV9SqwPlVwXVJS0u91v/3tb/HEE09g3bp1mD59ejhCJaJ+OD1efPBtvaQzdE4NzZ3icHvR7fTA7Q39X+HDdXiNAmO2OHCoqSugc4+12rD9WEeIIwo/UYRkq7GHkqJ7kABgxYoVWLJkCaZPn46ZM2di9erVsNlsWLp0KQBg8eLFyMnJQVlZGQDgv//7v/HII4/g9ddfx6hRo2A2mwEACQkJSEhIiNjrIBpuRFHEvgYryo9It3XFd5UfaUNtux1Wh6dXzUOGUY/rpuX4FziUGhMkOpt/7G7EP9AIlSCgMM+ES8al93re7fXhq+rWf61WHqEgQ8za40FinDbSYQxI8QnSwoUL0dLSgkceeQRmsxlTpkzBunXr/IXbNTU1UKn+3VH23HPPweVy4YYbbujVzqOPPorHHnssnKETDVvtNhf+sacx5MMwZ6oFarI68P6uelw7LSckRd31UTYcQqHjE0XsqumETqNGScHJIeFjrTZsPnR6LV60UcJUf0EcymY6w5TVaoXJZILFYoHRaIx0OESKcrzVho/3NsLpjmzRaU5SLG6YluufCScFh9uLNV8cidq/+il0ZoxKhtl6cluS4aA4PxkXFKSG/b6D+f5WfA8SESmD1ydix/F2fHO0HT4ZZBD1HT3YU29BYV6iZG02dPYwOaKgSLUdiVJYe1iDRESE2nY7Pq9qRlu3vIYNyo+2YVymQbJ6pGibbUQUKl0KGGJjgkREIWN3ebD5UCsONMpz9fkelxdbj7Xj4nPSJGmPBdpEgbEqYCkMJkhEFBJHWrrxyT5zxGuNzubb2k5MzjEhKX5oM2q6HG6YJdiQlGg46HZ44POJktYASk3R6yARkTw53F5s2N8k++QIOFkb9dWRoe8u/83R9mG5CStRMHyiiG6XvHuRmCARkeS+qm5Fj0s5u9lXN3ejeQi9P+02F/Y3yHMYkUiuumU+zMYEiYgkZbY4FLexpiieLNgO1tdHWmUxM49ISeS+mjYTJCKS1OdVzYqc6n60xYZGy+CLrJusDlT32ROOiM5O7nsWMkEiIslYetwwW5RbqPx19eB7kbYda1dkQkgUaexBIqJho65D2asA17TbsXcQw4NWhxtHW2whjIgoerEGiYiGjdp25a8DtKmqGW3dge0Rt6fOwtojoiB1O+W9WCQTJCKSjNJ7kADA7RXx8V4zPN6Blyjw+sRB9TYRUW+sQSKiYcFid8v+F16gWruc+Oxg82nHj7Xa4P5X4nSoqQt2BS1lQCQ3NqcXoox7YANeSfvDDz8cdOOXXXYZYmNjB30dESlPbRT0Hn3XvgYrUg06TBuRBADYW2/Bhv1NMMbGYPbYVHxb2xnZAIkUzieKsLm8SNDJc1OPgKO65pprBtWwIAg4fPgw8vPzBxsTESlQNAyv9fXloVYkx53cgmTjgZM9StYeNz7a3RjJsIiiRrfDo/wECQDMZjPS09MDOtdgMAQVEBEpU12H8gu0+/KJIj7e2whRBIuxiULgZKG2PtJhnFHACdKSJUsGNVz2wx/+EEajMaigiEhZOmyuqKk/6ksJ+8kRKZWcf28EnCC9/PLLg2r4ueeeG3QwRKRM0dh7REShJ+fFIoOaxdbT0wO7/d/1BidOnMDq1auxfv16yQIjIuU42sqtNoho8OS8WGRQCdKCBQvw2muvAQA6OztRXFyMp59+GgsWLGDPEdEw4/R4UdMWfQXaRBR6XdHWg1RRUYFZs2YBAN5++21kZGTgxIkTeO211/D73/9e0gCJSN6Ot9rh8bGAmYgGL+p6kOx2u3+W2vr163HddddBpVLh/PPPx4kTJyQNkIjkjTvZE1GwbE6PbBeLDCpBGjNmDN5//33U1tbik08+weWXXw4AaG5u5sw1omHE4/XheBs3ayWi4Hh8Inrc8lyRPqgE6ZFHHsEDDzyAUaNGobi4GCUlJQBO9iZNnTpV0gCJSL5q2u1weTgNnoiCJ9dhtqCWr7zhhhtw0UUXobGxEYWFhf7jc+bMwbXXXitZcEQkbxxeI6Kh6nJ6ENgS1OEVVILU09MDo9GIzMxMACen+b/33nuYMGECZs6cKWmARMOZw+1FZW0nGjp7cF6OCQVpCVCrhEiHBQAQRRHHWjm8RkRDE1U9SAsWLMB1112HO+64wz/NPyYmBq2trVi1ahXuvPNOqeMkGja8PhF1HXYcbbFhf6PVP4R1os2OBJ0GE7KMODfbiMS4GDR3OXG0xYbmLgesPW7oNGpcNjEDSfHakMfZ0u3kbvZENGRyXSwyqASpoqICzzzzDIB/T/PftWsX3nnnHTzyyCNMkIiC0O30oOJEB/Y2WPrd3qLb6cH24+3Yfrwd+hg1HGcobnx9Ww0uPicN5+WYQhpvPVfPJiIJyHW7kaASJE7zJ5KO1yfi6yOt2FXTCe8g1hM6U3IEAC6PDxv2N0GrUeGcjNBtGs3tRYhICnLtQeI0f6II6nK48fbOWuw43jGo5CgQG/Y3ocPmkrTNU0RRRH0nEyQiGjpbNCVIcpvm/+yzz2LUqFHQ6/UoLi7Gtm3b+j133759uP766zFq1CgIgoDVq1eHL1Ci72i3ufD61ho0dDpC0r7L48NHexrh9ko/Db/N5kIP64+ISAJR1YN0ww03oKamBjt27MC6dev8x+fMmeOvTQqXN998EytWrMCjjz6KiooKFBYWYu7cuWhubj7j+Xa7Hfn5+XjyySf9s/CIwq3L4ca7FXUhL3Ju7XLivV31kvckcXiNiKTi8vjg9MjvDy5BlOsa3wEqLi7GjBkz8Mc//hEA4PP5kJeXh3vuuQcPP/zwgNeOGjUK9957L+69995B3dNqtcJkMsFisXBIkQbN4fbi/3bUoq07NMNfZ6JRCZg5OhnF+SmStPeP3Y041NQlSVtERItLRiIlQRfy+wzm+zuoHiQA+PLLL/HDH/4QJSUlqK+vBwD85S9/wZYtW4JtctBcLhd27tyJ0tJS/zGVSoXS0lKUl5dLdh+n0wmr1drrQRQMt9eHDyrrw5ocASeX8//6SBsaLdL0/NR12CVph4gIkOcwW1AJ0jvvvIO5c+ciNjYWu3btgtPpBABYLBb85je/kTTAgbS2tsLr9SIjI6PX8YyMDJjNZsnuU1ZWBpPJ5H/k5eVJ1jYNH6Io4p97zSGrOQpExYnOIbfRxvWPiEhiUZMg/frXv8aaNWvwwgsvICYmxn/8wgsvREVFhWTBycXKlSthsVj8j9ra2kiHRAr0eVUzjkR4a47q5m5YetxDaoP1R0QkNTmuph3UOkhVVVWYPXv2acdNJhM6OzuHGlPAUlNToVar0dTU1Ot4U1OTpAXYOp0OOl3ox0YpOvl8IjYcaML+hsgPzfpEEd/WdmL2OWlBt1FlZu0REUnL5pJfghRUD1JmZiaqq6tPO75lyxbk5+cPOahAabVaFBUVYePGjf5jPp8PGzdu9C89QBRJLo8P71fWyyI5OmVvg8W/fclgtXY7uf4REUlOjqtpB5UgLV++HD/96U+xdetWCIKAhoYGrF27Fg888EDYtxlZsWIFXnjhBbz66qs4cOAA7rzzTthsNixduhQAsHjxYqxcudJ/vsvlQmVlJSorK+FyuVBfX4/KysozJnxEQ+Hy+PDerjqcaJNXQbPT7cO2Y+1BXbunziJxNEREgM0pv7rGoIbYHn74Yfh8PsyZMwd2ux2zZ8+GTqfDAw88gHvuuUfqGAe0cOFCtLS04JFHHoHZbMaUKVOwbt06f+F2TU0NVKp/54ENDQ29FrN86qmn8NRTT+Hiiy/Gpk2bwho7RS+P14e/f9sQ0YLsgWw/3g4RImaNDXyoze314YBZPj1hRBQ9up1Dq40MhSGtg+RyuVBdXY3u7m5MnDgRCQkJUsYmW1wHiQbi84n4aE9jxAuyAzEpx4TSiRlnPxHA3noLNuxvOvuJRESDJAjATy4dC5VKCOl9BvP9HXAP0ooVKwIOYNWqVQGfSxRt5DBbLVB76i3IMOoxKdcU0LlERKEgikC3ywOjPubsJ4dJwAnSrl27AjpPEEKb/RHJ2Y7j7ditsDqdzYdbMDI1rt9fTG3dTmypboXZIs/hQiKKDt0OhSZIn3/+eSjjIFK86uYubKlujXQYg+by+PDp/iZcNy2313GfT8Tmwy34ttYCn7J3JCIiBbDJbLHIoLcaIaJ/c7i92HigGUrNI0602VF+pM0//d/l8eHDbxuwq6aTyRERhUWXzBKkoGaxAcDGjRuxceNGNDc3w+frvabKSy+9NOTAiJTky8Otit9+45ujbdhV24EJmUY0WHrQbHVGOiQiGkbk1oMUVIL0+OOP41e/+hWmT5+OrKws1h3RsFbXYce+BmXVHfXH6fahsrYz0mEQ0TAkt+1GgkqQ1qxZg1deeQU/+tGPpI6HSFF8PhGfHVTu0BoRkVzIbcPaoGqQXC4XLrjgAqljIVKcw83daOt2RToMIiLFi4oE6bbbbsPrr78udSxEirP9eHBbdhARUW9RUYPkcDjw/PPP49NPP8XkyZMRE9N73QIuFEnDwfFWG1q6WMhMRCQFt1eEw+2FPkYd6VAABJkg7d69G1OmTAEA7N27t9dzLNim4YK9R0RE0upyeJSdIHHRSBruGi09qOvoiXQYRERRxeb0IM2gi3QYAIawDhJRJHTYXKjv7IHb64NOo4ZWI0AlCNCoVMhLjg1bD+a3nApPRCQ5ORVqB50gdXZ24sUXX8SBAwcAABMnTsSyZctgMp1900uiwbA63NhTZ8GBRiu6BlgnY1ymAXPPzYQ6xLtBOz1eVCtkM1oiIiWRU4IU1Cy2HTt2oKCgAM888wza29vR3t6OZ555BgUFBaioqJA6RhqmHG4v1u8z4+Utx7HtWPuAyREAVJm78OG39XB7fQOeN1SHzN1we7nwERGR1OQ0k00QxcEvcTdr1iyMGTMGL7zwAjSak51QHo8Ht912G44ePYrNmzdLHqicWK1WmEwmWCwWGI3GSIcTlaqbu/HZwSbYnIPfviPVoMP88zKRkhCacew3ttWgkTvbExFJLj8tHgum5ISs/cF8fwfdg/TQQw/5kyMA0Gg0+NnPfoYdO3YE0ySR3956Cz7a3RBUcgQArV1O/G1bDXbXdUobGIC2bieTIyKiEFH8EJvRaERNTc1px2tra2EwGIYcFA1fVeYufHqgachbd7i9IjYeaMa6vWZ4JBxy29dglawtIiLqTU5DbEElSAsXLsSyZcvw5ptvora2FrW1tXjjjTdw22234eabb5Y6Rhomqpu7sG6vWdJ9zQ40WvF/O+pgdbiH3JbXJ+KgmQkSEVGo2F1eeH3yqPEMahbbU089BUEQsHjxYng8J7O9mJgY3HnnnXjyySclDZCGh/0NVmzY3wRfCHZ9bbI68Oa2WlwzNWdI62scbu4KetiPiIjOThQBm8sDoz7m7CeHWFBF2qfY7XYcOXIEAFBQUIC4uDjJApMzFmlLq7K2E5uqmiXtOToTXYwKVxdmIzcpuM8pi7OJiEJv4Yw8ZCfGhqTtwXx/D2mhyLi4OEyaNGkoTdAw5vL48NnBZhxoDM+wldPtw3sV9SidmIEJWYNLbM0WB5MjIqIwkEsdUsAJ0ooVK/DEE08gPj4eK1asGPBcblZLZ9NkdeCfexrRYR96bdBgeHwi1u01o7XbiYvGpAa88nZlbUeIIyMiIkA+M9kCTpB27doFt9vt/+/+cLNaGojPJ2L78XZ8c7Q9JPVGgdpxvAOddjfmT8o668rbNqcHh5q4cjYRUTgoLkH67ga13KyWgtHW7cSnB5rQ0CmPoarq5m58vKcRV07KgmqAJGnniQ7ZzKogIop2chliC2qaf09PD+x2u//fJ06cwOrVq7F+/XrJAqPo4fb68FV1K9ZurZFNcnRKdXM3Pt7b2O9aSdXNXaio4fAaEVG4nG1bqXAJKkFasGABXnvtNQAnN62dOXMmnn76aSxYsADPPfecpAGSstW22/HXb05g27F22fbCHG7qxt+21aC5q3fy1trtxCf7hr5oJRERBU7RPUgVFRWYNWsWAODtt99GZmYmTpw4gddeew2///3vJQ2QlMnm9GDD/ia8vbMOnWEuxA5Ga7cLb2yrxedVzdh+vB27ajrwYWUDXJ7QbnxLRES92VzyWG8uqGn+drvdv6XI+vXrcd1110GlUuH888/HiRMnJA2QlKXH5cX24+3YXdepuB3vvT4RlTWdkQ6DiGhYc3l8cLi90MeoIxpHUAnSmDFj8P777+Paa6/FJ598gvvuuw8A0NzczIUTFaLR0oMDjVYcbbHBJ4pQCQI0KgEatQr6GDWyTHpkJ8YiNykWMeqzdzRa7G5U1HRgX4NFcYkRERHJi83pUWaC9Mgjj+CWW27Bfffdhzlz5qCkpATAyd6kqVOnShogScft9eFgYxcq6zrR2uUc8Nza9pNF+LoYFcZnGjA+04g0g65XsuTx+nC01Yb9DVacaLNHdNo+ERFFj26nBykJwW8NJYWgtxoxm81obGxEYWEhVKqTX5rbtm2D0WjE+PHjJQ3ybJ599ln87ne/g9lsRmFhIf7whz9g5syZ/Z7/1ltv4Ze//CWOHz+OsWPH4r//+78xf/78gO+ntK1G2rqd2FNvwYHGLjjcwY/tCgJg0MdAq1HB4fKixy2fTQWJiCh6XH5uBs7NNknebli2GsnMzERmZmavYwMlJaHy5ptvYsWKFVizZg2Ki4uxevVqzJ07F1VVVUhPTz/t/K+//ho333wzysrK8P3vfx+vv/46rrnmGlRUVOC8884Le/yhIIoiWrqcONpqw9EWG5qs0kytF0XA2iP/gmsiIlI2OWwMHnQP0saNG7Fx40Y0NzfD5+s90+ell16SJLhAFBcXY8aMGfjjH/8IAPD5fMjLy8M999yDhx9++LTzFy5cCJvNho8++sh/7Pzzz8eUKVOwZs2aM97D6XTC6fz3kJTVakVeXl5IepAcbi+sPW7E6zSI06r7XZnc6xNhd3lgc3rR7fTA6nDD2uNGa7cLTVaHomdfCQIgQOCQHRHRMFWYZ8Kl4zMkbzfkPUiPP/44fvWrX2H69OnIysqK2PYiLpcLO3fuxMqVK/3HVCoVSktLUV5efsZrysvLT9tLbu7cuXj//ff7vU9ZWRkef/xxSWIOhN3lRbvdhR7XyeTH5vSgpcsJl1eE1+eD0+2DZ4ChreR4LcbmJSBWq4Y+Rv2vZMoLm8uDbocHXQ4Puhxu2GUylVKtEjAmPQGTckxISdAiNkYNS48bG/Y3oa6jJ9LhERFRmHlkMNknqARpzZo1eOWVV/CjH/1I6ngGpbW1FV6vFxkZvbPMjIwMHDx48IzXmM3mM55vNpv7vc/KlSt7JVWnepBCQR+jxqjU+NOOOz1e1LbbUdvRg5YuJ1q7nXC6T/YSqVUC0gw6ZJr0GJuegNykuIDu5fb6YO1xo7PHjQ6bC81dTjR09oRtFVN9jBqFuSZMzktEgq73RzExTosbinKxu86CTVUt7E0iIhpGDPqYSIcQXILkcrlwwQUXSB2LbOl0Ouh0ka2m12nUGJNuwJh0g/+YKIrw+kQIgnDWDVfPJEatQkqC7uRMgbR/H7f0uP1JWLvNBUuPG51295AKvL8ry6THudkmjM8yDLiEgCAIKMxLRLxOjY/3mFkQTkQ0TBj0QZdISyaoCG677Ta8/vrr+OUvfyl1PIOSmpoKtVqNpqamXsebmppOKyA/JTMzc1Dny5kgCNCopR/eNMXGwBQbgzHpCb2O210etNtc6LC50WpzosPmQqfdjS6Hp98eHl2MCgadBmkGPTJNeoxIjkNyvHZQ8YxJN+DKyQL+sbuRSRIR0TBgVGoPksPhwPPPP49PP/0UkydPRkxM7xeyatUqSYI7G61Wi6KiImzcuBHXXHMNgJNF2hs3bsTdd999xmtKSkqwceNG3Hvvvf5jGzZs8K/lRP2L02oQp9UgN6n3cZ9PhN3thdvjg9vng0alglajgk6jCmiRyUAUpCXgsokZWLe3/6FQIiKKDortQdq9ezemTJkCANi7d2+v58JdsL1ixQosWbIE06dPx8yZM7F69WrYbDYsXboUALB48WLk5OSgrKwMAPDTn/4UF198MZ5++mlceeWVeOONN7Bjxw48//zzYY07mqhUwskaohCPQk7IMqLJ6sAubgdCRBTVFJsgff7551LHEbSFCxeipaUFjzzyCMxmM6ZMmYJ169b5C7Framr8C1kCwAUXXIDXX38d//mf/4mf//znGDt2LN5///2oWQMp2s0em4bWbpd/pW8iIooucVo1NBKNPgxF0Osgffnll/jzn/+Mo0eP4q233kJOTg7+8pe/YPTo0bjoooukjlNWlLaSdrTpcXnx+rYaLlpJRBSFMox63FI8IiRtD+b7O6gU7Z133sHcuXMRGxuLiooK/yKKFosFv/nNb4JpkihgsVo1rirMglYT+b8wiIhIWnIYXgOCTJB+/etfY82aNXjhhRd6FWhfeOGFqKiokCw4ov6kG/S4bKL0q6wSEVFkKTpBqqqqwuzZs087bjKZ0NnZOdSYiAJyToYBs89JjXQYREQkIUUnSJmZmaiurj7t+JYtW5Cfnz/koIgCVTQyGaUTMhCh3W6IiEhiclhFGwgyQVq+fDl++tOfYuvWrRAEAQ0NDVi7di0eeOAB3HnnnVLHSDSgSbkmzJ+UFdRq4kREJC9y6UEKKoqHH34YPp8Pc+bMgd1ux+zZs6HT6fDAAw/gnnvukTpGorM6J8MAnUaFj3Y3wuXxRTocIiIKklx6kIKe5g+c3JOturoa3d3dmDhxIhISEs5+URTgNH/5arT04IPKBvS4pNk3joiIwketEnDPpWNCtuj0YL6/h9SPpdVqMXHixKE0QSSpLFMsfjA9D+/tqlfUOklj0hNwbrYRyfFatNtc+OdeM3vCiGjYSdBpwr4jR3+CTpAcDgd2796N5uZm+Hy9f5FfffXVQw6MKFjJ8VrcNCMPH1Q2oMnqiHQ4Z5Vp0uOK8zL9K8cmxmlxY1EuPqhsQLfTE+HoiIjCRy71R0CQCdK6deuwePFitLa2nvacIAjwejm8QZEVr9PghqJcrN9vxuGm7kiH068EnQZXFWaftqx+ulGPH0zPw+vbauBw8/8nIhoe5FJ/BAQ5i+2ee+7BjTfeiMbGRvh8vl4PJkckF1qNCldOysIFBSmyXAYgRi3g+4VZJzf6PQNTXAwuP5dLGBDR8GGUUQ9SUAlSU1MTVqxY4d8QlkiuBEFAcX4KrpmSg3idOtLh+KkEAVdMykKWKXbA8wrSElA0MilMURERRZbie5BuuOEGbNq0SeJQiEJnVGo8fnj+SOSnxUc6FADAnAnpKEgLbNbnhQWpyDTpQxwREVHkyakGKahp/na7HTfeeCPS0tIwadKkXvuxAcBPfvITyQKUI07zV7avqlux7Vh72O8rCEBeUhwK80wYk24Y1LXHW214b1d9iCIjIpKHxSUjkZKgC1n7IZ/m/7e//Q3r16+HXq/Hpk2bek3JEwQh6hMkUrYLx6QiVqvG5kMtCH4VsMHJT4vHJePSYYoNrvt4VGo8Mox6RczKIyIKljHI35GhEFSC9Itf/AKPP/44Hn74YahUQY3SEUXUtBFJiI1RY/2+JvhCmCVpNSrMHpuGSbmmIbc1c3Qy/v5tgwRRERHJT6xWjRi1fHKKoBIkl8uFhQsXMjkiRZuQZYRWo8LHuxvh8UmfJKUZdLhyUhaS4rWStFeQFo/UBC1au12StEdEJCdGGRVoA0EWaS9ZsgRvvvmm1LEQhV1BWgKumZrT71T7YBXmmXDTjDzJkiPg5PD1jNHJkrVHRCQnxlj5FGgDQfYgeb1e/Pa3v8Unn3yCyZMnn1akvWrVKkmCIwqHvOQ4LL1wFHbXW7D9WDvsQ9jHTatRoXRCBsZlDq4IO1Bj0hKg1ai4DQkRRR259SAFlSDt2bMHU6dOBQDs3btX0oCIIkGjVmHaiCRMzjHhoLkLu2o70drlHFQbqf8aUkuWsNeoL41ahRHJcahulu/q4EREwZDTFH8gyATp888/lzoOIlnQqFU4L8eE83JM6LC50GDpQWOnA8fbbOhy9N4XTSUIUAlAhkmPopFJyE+ND8smiwVpCUyQiCjqyGkGGzCIBGnFihV44oknEB8fjxUrVvR7niAIePrppyUJjiiSkuK1SIrX4tzskzPQ2rqdsLu8MMbGwKDTQKWKzB4g+WnxUAlCSGffERGFm2KH2Hbt2gW32+3/7/6E4y9ookhISdAhJdJBANDHqJGTFIvadnukQyEikoxii7S/O6zGITaiyCpIi2eCRERRQxejgk4jn/0ygSCn+RNRZBWkB7aPGxGREshteA0YZA1SoDjNnyi0jPoYbj1CRFFDbgXawCBrkL6roqICHo8H48aNAwAcOnQIarUaRUVF0kZIRGc0PsvABImIooJRZlP8gSBrkFatWgWDwYBXX30VSUlJAICOjg4sXboUs2bNkj5KIjrN+EwDthxuhTcE26QQEYWTQYZDbIIoDn6ucE5ODtavX49zzz231/G9e/fi8ssvR0NDdG+oabVaYTKZYLFYYDQaIx0ODWMfftuAI1wTiYgU7qrCLIxJD80OBN81mO/voIq0rVYrWlpaTjve0tKCrq6uYJokoiBMzGKCTkTKJ8ci7aASpGuvvRZLly7Fu+++i7q6OtTV1eGdd97BsmXLcN1110kdY7/a29uxaNEiGI1GJCYmYtmyZejuHviv6eeffx6XXHIJjEYjBEFAZ2dneIIlCoH81HjEaeU1NZaIaLDkWKQdVIK0Zs0aXHHFFbjlllswcuRIjBw5ErfccgvmzZuHP/3pT1LH2K9FixZh37592LBhAz766CNs3rwZt99++4DX2O12zJs3Dz//+c/DFCVR6KhUAsazF4mIFEyrUUEfI78/9IKqQTrFZrPhyJEjAICCggLEx8dLFtjZHDhwABMnTsT27dsxffp0AMC6deswf/581NXVITs7e8DrN23ahO9973vo6OhAYmLigOc6nU44nf/euNRqtSIvL481SCQLrd1O/KX8RKTDICIKSmqCFj8qGRWWe4W8BumU+Ph4TJ48GZMnTw5rcgQA5eXlSExM9CdHAFBaWgqVSoWtW7dKeq+ysjKYTCb/Iy8vT9L2iYYiNUGHEclxkQ6DiCgochxeAwaRIO3evRs+ny/ghvft2wePx3P2E4NkNpuRnp7e65hGo0FycjLMZrOk91q5ciUsFov/UVtbK2n7REM1dURipEMgIgqKHAu0gUEkSFOnTkVbW1vADZeUlKCmpmbQAT388MMQBGHAx8GDBwfd7lDodDoYjcZeDyI5GZ0aj+R4baTDICIaNLltUntKwFGJoohf/vKXiIsLrCvf5XIFFdD999+PW2+9dcBz8vPzkZmZiebm5l7HPR4P2tvbkZmZGdS9iZRKEARMyUvEZwebz34yEZGMyLUHKeAEafbs2aiqqgq44ZKSEsTGxg46oLS0NKSlpQXUfmdnJ3bu3Onf3uSzzz6Dz+dDcXHxoO9LpHQTs434+kgbHG5vpEMhIgqYSaY1SAEnSJs2bQphGIM3YcIEzJs3D8uXL8eaNWvgdrtx991346abbvLPYKuvr8ecOXPw2muvYebMmQBO1i6ZzWZUV1cDAPbs2QODwYARI0YgOTk5Yq+HaKhi1CoU5pmw9Wh7pEMhIgqY4ou05Wjt2rUYP3485syZg/nz5+Oiiy7C888/73/e7XajqqoKdrvdf2zNmjWYOnUqli9fDuBkz9jUqVPx4Ycfhj1+IqlNG5EErUbR/1sT0TCii5HnGkjAENdBGq64FxvJ2dfVrdh6TLm9SAk6DaaNTEJLlxOHmrq4GS9RFEsz6PDD80eG7X6D+f6WZ+k4EQVt2sgk7KrthMsT+LIcchCjFnB+fgqm5CVCoz7ZCzZrbCrqOnoQr1NDq1HhrR11intdRNQ/uQ6vAQofYiOi0+lj1JialxjpMAYl3ajDLcUjMX1Usj85AoB4nQbjMg3ITYpDukHP9Z6IooxcC7QBiRKkO+6447Qp90QUOdNGJiFWAZvYCgJQNDIJN80YEdA6TkUKeV1EFBijXr4DWZIkSFdccQXmz5+Pxx57DDabTYomiWgI9DFqlOSnRDqMARn0Glw/LRezz0mDWiUEdI1Oo8aMUUkhjoyIwiXqh9gWLFiArVu3IiMjAxdccAHWrFkzqG1JiEh6k3JMSDXoIh3GGaUadFhUPBJ5QewhV5ibCIOM/+okosBF/RAbAKjValx55ZW477778J//+Z+YOHEi/v73v0vVPBENkkol4OKxZ190NdyMsTG4dmpO0ENlGrUKF45JlTgqIooEua6iDUiUIM2bNw8jR47ELbfcgt27d+MPf/gD1q5di/fffx/33nuvFLcgoiCMSInDeTmmSIfhF6tV49qpOUjQDa0HaEKWETlJg1+pn4jkI1arlvW6bZKsg1RZWYlJkyZBrT79L8Lx48eHfXPZUOM6SKQ01c1d+PxgC7qdnojFoNWocENRLjKMeknaa+ly4vWtNfBxKTciRco06XHzzBFhvedgvr8lSd2mTJlyxuQIAD7++GMpbkFEQzAm3YAflYxEdqI0yclgaVQCri7Mliw5Ak4uMDc5Tz69Y0Q0OHIeXgMkrEGyWq3YvHkzfv/73/c6np+fL9UtiGgI9DFqLJiSI2mSEgi1SsCVk7OCKsg+m5L8FFl30RNR/4yx8p5sEVR0NTU1qKys7PU4ceIERFFEfHw8fvKTn0gdJxFJQB+jxnXTcvBORR2arc6Q30+jEvD9wmyMTo0PSfv6GDXOyzGh4kRHSNonotCR8ww2YJA9SJdeeilSUlIwatQoLFmyBJ988gnS0tJQU1ODF198ESdOnEBXV1eoYiUiCehj1Lhpxgicn58S8PpDwdBqVLhmak7IkqNTikYmhfR1EFFoRNUQ25YtW3DHHXegtrYWHR0d+Oqrr/DnP/8ZgiBg5syZyMvLC1WcRCQhtUpASUEKFhWPCMmaQnHakz1VoRhW6ytBp8HELE6WIFIaOS8SCQwyQdq6dSu+/PJL3HXXXTh06FCoYiKiMElJ0OGaqTnQxUhXx5MUF4OFM/KQZQrfNPzpo5KgEtiLRKQUgiDvbUaAQSZIU6dOxebNm/GDH/wAc+fOxV133cU92IgULjVBh6smZ0syTDUmPQELZ4xAYtzZ91WTUmKcFqPTQjuUR0TSSdBpem1MLUdBRXfLLbdg3759SEpKwrnnngufzwev1yt1bEQUJnnJcbj83AwE2wkTq1Vj/qQsXFWYHbHNZMekJUTkvkQ0eHIfXgOGMM0/Li4Ov/71r7F161Z8//vfx5w5c/DUU0+hp6dHyviIKEzGZxoD3sJDq1GhID0B3xufjluKR+D2WfkYl2kIcYQDy0+L5zAbkULIfQYbEOQ0/+/Kz8/HBx98gPXr1+O+++7D008/jcbGRiliI6IwmzEqGd0ODyprO3sdN+g1mJBlhEGvQWKsFtmJetl1j+tj1MhO1KOug3+kEcndsEiQTrn88svx7bff4o9//KNUTRJRBFwyLg1pBh22HmuHtceNczIMmDMhHfqYyAydDUZ+WgITJCIFGFYJEgBoNBpuTkukcIIg4LwcEyZkGWG2OpCTqJxNYQvS4rH5UEukwyCis0iMk3+CJK8+ciKSDbVKUFRyBJyczZaSEN4ZdEQ0eEroQWKCRERRJT+Vs9mI5EyrUSFOK+81kAAmSEQUZQrSuR4SkZwpYYo/wASJiKJMlilWEfUNUtHFqMKypQuRVJQwvAZIXKRNRCQHE7KMKD/SFukwQkqrUWFKXiKKRiZBH6NGk9WB8iNtONZqi3RoRANigkREFCETsoz45mgbRDHSkUhLrRIwJj0B52QkYGRKPGK+sxZVhlGPBVOy8detNWjtckYwSqKBMUEiIooQU2wMchJjo2ZNJK1GhcLcREwZkYgEXf+/tgVBwMVj0/BORV0YoyMplBSk4ECjFZ12d6RDCTmlJEisQSKiqDQhyxjpECRhjI3BzTNH4KKxqQMmR6eMSIlDPjfuVZQJWUacn5+CH0zPQ5pBF+lwQo4JEhFRBI3NSECMWtl7s6UadFg4Iw/J8YNb22nW2DTuS6cQyfFaXDo+HQAQr9PghqLcQf+8lUQQAKNeGYNXTJCIKCrpNGqMzYjsBrpDkW7U4cai3IB6jfpKjtdicp4pBFGRlDQqAfMnZUGr+fdXsT5GjZmjkyMYVWgl6DSy28exP8qIsh/t7e1YtGgRjEYjEhMTsWzZMnR3dw94/j333INx48YhNjYWI0aMwE9+8hNYLJYwRk1E4TJ9ZBKU2JGSmqDFdVNzh7T/XUl+CmK18t8/bzibOTr5jENq4zIMihmGGiylrIEEKDxBWrRoEfbt24cNGzbgo48+wubNm3H77bf3e35DQwMaGhrw1FNPYe/evXjllVewbt06LFu2LIxRE1G4pCToMCZdWStrJ8bF4NppuUNObvQxapTkp0gUFUktzaDDjFFn7ilSqQQUjUwKc0ThkaigBEkQRWVOhD1w4AAmTpyI7du3Y/r06QCAdevWYf78+airq0N2dnZA7bz11lv44Q9/CJvNBo3mzF3ZTqcTTue/p81arVbk5eXBYrHAaIyOQlCiaNXc5cDab2oiHUZAEuNicENRLgx6ab5EfD4Ra7dx2r/cqAQBN8/MQ7pR3+85Hq8PL311DDanN4yRhd4FBSkojmDibrVaYTKZAvr+VmwPUnl5ORITE/3JEQCUlpZCpVJh69atAbdz6k3qLzkCgLKyMphMJv8jLy9vSLETUfikG/SKmNUldXIEnOyJuHhsmmTtkTQm5RoHTI4AQKNWYUpe9PUiJSmoAF2xCZLZbEZ6enqvYxqNBsnJyTCbzQG10draiieeeGLAYTkAWLlyJSwWi/9RW1sbdNxEFH7Fo1OgVsmzGEmjEjApx4Qbp+dJmhydkpccG1ShN4XOuMzARh7GKmx4OBBKGmKTXYL08MMPQxCEAR8HDx4c8n2sViuuvPJKTJw4EY899tiA5+p0OhiNxl4PIlKOTJMeV07OgkYGSZJKEJCTFIvzckyYfU4q/t9Fo1E6MSNkSYwgCBidKv8etOHCoNcg2zRw79EpSfHaqEtuTQraJ1F27/z999+PW2+9dcBz8vPzkZmZiebm5l7HPR4P2tvbkZmZOeD1XV1dmDdvHgwGA9577z3ExCjnB0ZEwSlIS8CCKTn4++4GuDy+sN9fEE7OTjo/PyXswwz5afHYU8/ZunJQkJ4AYRBTK3OTYnHQ3BXCiMInXqeGTqOcmZWyS5DS0tKQlnb2MfOSkhJ0dnZi586dKCoqAgB89tln8Pl8KC4u7vc6q9WKuXPnQqfT4cMPP4ReH1gmT0TKNyIlDrdeMAoNnT2o7+zB3noL3N7wzFOZPykL50RoXaYRyXGIUQthe63Uv8F+BvKS46ImQUqMVU79ESDDIbZATZgwAfPmzcPy5cuxbds2fPXVV7j77rtx0003+Wew1dfXY/z48di2bRuAk8nR5ZdfDpvNhhdffBFWqxVmsxlmsxleb3TNFCCiM4vXaTA2w4BLxqXj2mlDW2soUIV5poglR8DJgt+85LiI3Z9OStAFPrx2Sm5SbIiiCT8lDa8BCk6QAGDt2rUYP3485syZg/nz5+Oiiy7C888/73/e7XajqqoKdrsdAFBRUYGtW7diz549GDNmDLKysvwPFl4TDT85ibG4IcjVqvsSBJxxa5NUgw6zZTCTrCAt+gp+lWZMxuCG1wAgMU4Lg0K25jgbJRVoAzIcYhuM5ORkvP766/0+P2rUKHx3madLLrkECl32iYhCJM2gww9m5OGDynq0dbuCbqcwLxFFI5PwRVULqpu7oVYJyEuOxcXnpMtia4XRqfEQBIC/AiMn2FlpuUlxONBolTia8FPSFH9A4QkSEZEUTLEx+MH0PPz92wY0WR0YkRIPnUaF/Q2BfSmlJmgxa0wqNGoVrirMRpPVAVNsTFiG7wIVr9Mgw6iH2eKIdCjDUqpBh5zE4IbLcpNioyJBYg8SEZEC6WPUuG5aLnyiiBi1CqIoosvhQW27fcDr1CoBc8/L7NVLlHGWRQAjZUpeItZZAlsnjqQ1c1TyoIfXTomW+jHWIBERKZRaJSDmX4mOIAiYe27GWXuBLhmXhnSDPBOivsZnGpBuPH1z1EgShJMJQIZRD1NsDBJ0GsRp1Wdcs6rvscS4GEVsfpoUF4NzMoKvATPFKuN1DiROq6wp/gB7kIiI+mXQx6B0Qjo+2t14xudnn5OKybmJ4Q1qCARBwKwxaXinoi7SoQA4mRxdPjETE7NPX3zX6xPR1u1Ec5cTBv3J4UGdRoV2mwttNhdS4rVISdChtt2Ot3fK4/X0Z8bo4HuPTslPi0dlTac0AUVAUpyy6o8A9iAREQ1obIYBcyako+/32wUFKSgaeebd2OVsREocRqZEfshGEIDSCRlnTI6Ak7156UY9zssxYWRKPPQxagiCgJQEHc7JMCAl4WRPWF5yHApkvCWHMTYGEwLcWmQgSt92RGnDawATJCKis5qcm4j5k7L+NTMtDjdOz43ojuRDNWtsGrSayP76nzEqGeflmCRpa9aY1DPutZdm0KEgPSGixfLFo5OhkmCLm5xEZe+pp7QCbYBDbEREATknw4DcpFjEaZX/azPNoMONRbn4oLIB3U5P2O+v1ahQNFK6neqT4rWYlGvqNQSVGBeD66blIE6rgSiKaOl2orbdjpp2O463Dlx4L5WUBC0mZkmzd6cgCBiTnoDK2k5J2gs3pU3xB9iDREQUsGhIjk5JN+qxcGYekiPwxTUlL1HyXp2LxqTiorGpiNOqEadV49qpOf6flyAISDfoUTQyGddOzcXM0eEZGr2gIFWS3qNTxih4mI09SEREpBhGfQwuHZ8e1iJnrUaFaSOk6z06JUatwoxRyZiSlwib04PEAYqCLxyTii6HJ6RrC2Un6iVPaE72YKphdylvayzWIBERkaLkJcedtt9XKOuTCnMTEasNXU1QjFo1YHJ0ymUTMzAihOsLXTgmVfI2Tw2zKY1Br1HcFH+ACRIR0bBXUvDvgvPsRD1umTkiJEmSPkYtae3RUKhVAuadl4l4nfRf3OfnpyA3KTTJ19j0yG16HKxIDONKgQkSEdEwl5sUh7zkOBj0Gnx/cjaS4rU4P1/6Op0Lx6SEtPdosOJ1Glw+MfO0JRyGYkKWsVfCKbWcpFjoYpT11a3EAm2ANUhERIST6zqpVQLi/zWVfGpeEg40dqGlyylJ+5kmPSZJNK1fSqNS4zFtRBKOtnRjUq4JmaZYbD7UEvCedQk6DUamxKHH7YVWrcJlEzNCGq9aJWBUSjyqzF0hvY+UUpggERGRUmX32UhVpRJQOiEDb26vhU8Uh9S2IABzxqcPeTXpULloTCpmjU31x7dweh521XagytyNbqcbdpcXZ3oLUg06LJiSDaM+vAXI+WnKSpCUuIo2wASJiIj6kWnS47KJGVi/33zGBCEQWo0KJQUpSJfpBr4ATpuKr1IJKBqZ7F8pvaXLif/bUQuXx+c/Z3RqPK6YlBmR4uNRKfFQqwR4fUNLXMMlJYEJEhERRZmJ2UbYXB5sOdx6xue1GhWmj0xClikWzV0OtHa74Pb6IAJIN+hCsuZRuKUZdLhyUhY+qGyATxQxOdeE741Ll3SNo8HQx6iRkxiLmvbwLHg5FLFatWLXD1Nm1EREFDYzRiUjNkYNnUaFdKMeoiiivrMHFrsbk/MS/VtgjJDBHm+hMio1Ht8bnwanx4cZoyK/B19+WrwiEqRkhQ6vAUyQiIgoAH33TQtkraFoMzk3MdIh+OWnJWBTVUukwzgrpU7xBzjNn4iISHFMsTFINegiHcZZJSu0/ghggkRERKRIeX1WQJcjJQ+xMUEiIiJSoFCt1i0l9iARERFRWOUmxUq6CrjUtBpV2NeIkhITJCIiIgXSx6iRmiDfOiSlLhB5ChMkIiIihcqVcR2SkmewAUyQiIiIFEvOdUipCq4/ApggERERKZac65DSDfLdXiYQTJCIiIgUSs51SOlGecYVKCZIRERECibHOiSDXqP4PfiYIBERESlYXrL86pDSjcoeXgOYIBERESlalkl+yUi6ArZBORtFJ0jt7e1YtGgRjEYjEhMTsWzZMnR3dw94zX/8x3+goKAAsbGxSEtLw4IFC3Dw4MEwRUxERCStOK0GBr289p5PY4IUWYsWLcK+ffuwYcMGfPTRR9i8eTNuv/32Aa8pKirCyy+/jAMHDuCTTz6BKIq4/PLL4fV6wxQ1ERGRtDJl1osUDT1IgiiKYqSDCMaBAwcwceJEbN++HdOnTwcArFu3DvPnz0ddXR2ys7MDamf37t0oLCxEdXU1CgoKArrGarXCZDLBYrHAaDQG/RqIiIiksP14O7Ycbo10GACAOK0a/3FxYN+n4TaY72/F9iCVl5cjMTHRnxwBQGlpKVQqFbZu3RpQGzabDS+//DJGjx6NvLy8fs9zOp2wWq29HkRERHKRKaOi6GgYXgMUnCCZzWakp6f3OqbRaJCcnAyz2TzgtX/605+QkJCAhIQE/POf/8SGDRug1fa/4mdZWRlMJpP/MVAyRUREFG5pBp1sFoxU+gKRp8guQXr44YchCMKAj6EWVS9atAi7du3CF198gXPOOQc/+MEP4HA4+j1/5cqVsFgs/kdtbe2Q7k9ERCQlfYxaNpvDKn2ByFPkVfYO4P7778ett9464Dn5+fnIzMxEc3Nzr+Mejwft7e3IzMwc8PpTPUFjx47F+eefj6SkJLz33nu4+eabz3i+TqeDThcdP3AiIopOGUYd2m2uSIcRFQXagAwTpLS0NKSlpZ31vJKSEnR2dmLnzp0oKioCAHz22Wfw+XwoLi4O+H6iKEIURTidzqBjJiIiirQMox4HGrsiGoM+Rg1TbExEY5CK7IbYAjVhwgTMmzcPy5cvx7Zt2/DVV1/h7rvvxk033eSfwVZfX4/x48dj27ZtAICjR4+irKwMO3fuRE1NDb7++mvceOONiI2Nxfz58yP5coiIiIZEDlP9s0x6CHIphhoixSZIALB27VqMHz8ec+bMwfz583HRRRfh+eef9z/vdrtRVVUFu90OANDr9fjyyy8xf/58jBkzBgsXLoTBYMDXX399WsE3ERGRkqQl6KBWRTY5kUOSJhXFroMUSVwHiYiI5Gjt1hNotkauZOT6abkYkSK/veFOGRbrIBEREVFvOYmxEbu3ShCiqgeJCRIREVGUGJEcud6blAQttJroSSui55UQERENc7lJcRGrQ8pOjJ7eI4AJEhERUdTQalQRG+bKMkVueC8UmCARERFFkUgNs2UzQSIiIiK5ikSCFK9TwxQXHQtEnsIEiYiIKIpkGvVhL5bOjLLeI4AJEhERUVRRqQTkhbkXKSfKCrQBJkhERERRJ9zDbLlJ8l0cMlhMkIiIiKLMqDCuZq2LUSEtQRe2+4ULEyQiIqIokxinRUqCNiz3ykmMhSrCe8CFAhMkIiKiKFSQlhCW++QmRV+BNsAEiYiIKCqFL0GKvvojgAkSERFRVMow6pCg04T0HroYFdIN0Vd/BDBBIiIiikqCICA/LT6k98hJjIUgRF/9EcAEiYiIKGqFepgtWofXACZIREREUSsvOS6kq2rnRWmBNsAEiYiIKGqpVQLyU0MzzBanVSMtSuuPACZIREREUW1shiFE7SZEbf0RwASJiIgoqo1KCc0w2zkhSrzkggkSERFRFNOoVZIXaxv0GuQkRm/9EcAEiYiIKOqNy5S2t2dshiGqh9cAJkhERERRb0RyHPQxasnaGy9xwiVHTJCIiIiinFolYEy6NMNsiXExyDDqJWlLzpggERERDQPjJCqqlnq4Tq6YIBEREQ0DuUmxiNMObZhNoxJQmJsoTUAyxwSJiIhoGFCpBIzNGNow24QsI+JDvAGuXDBBIiIiGiaGsnaRIABFI5MkjEbemCARERENEzmJsTDog+sBGpOegKR4rcQRyRcTJCIiomFCEISgtx6ZPjJZ4mjkTdEJUnt7OxYtWgSj0YjExEQsW7YM3d3dAV0riiKuuOIKCIKA999/P7SBEhERyUQws9lGpsQh0xT9U/u/S9EJ0qJFi7Bv3z5s2LABH330ETZv3ozbb789oGtXr14d9auAEhER9ZVp0iMxLmZQ15yfnxKiaORLsQnSgQMHsG7dOvzv//4viouLcdFFF+EPf/gD3njjDTQ0NAx4bWVlJZ5++mm89NJLYYqWiIhIPkonZAS8svbIlDhkR/m+a2ei2ASpvLwciYmJmD59uv9YaWkpVCoVtm7d2u91drsdt9xyC5599llkZmYGdC+n0wmr1drrQUREpFR5yXG4ZeYIpCacvei6pGD49R4BCk6QzGYz0tPTex3TaDRITk6G2Wzu97r77rsPF1xwARYsWBDwvcrKymAymfyPvLy8oOMmIiKSA1NcDH4wIw+pBl2/54xKjUOWafj1HgEyTJAefvhhCIIw4OPgwYNBtf3hhx/is88+w+rVqwd13cqVK2GxWPyP2traoO5PREQkJzqNGlcXZiP2DCtsx2nV+N649DNcNTzIbjnM+++/H7feeuuA5+Tn5yMzMxPNzc29jns8HrS3t/c7dPbZZ5/hyJEjSExM7HX8+uuvx6xZs7Bp06YzXqfT6aDT9Z9hExERKZUpNgZXTsrCuxX18IkiACBGLWDBlBwkxg2fdY/6EkTxX++Gwhw4cAATJ07Ejh07UFRUBABYv3495s2bh7q6OmRnZ592jdlsRmtra69jkyZNwv/8z//gqquuwujRowO6t9VqhclkgsVigdFoHPqLISIiirD6zh4cae5GbYcdJfkpyE8b2rYkcjSY72/Z9SAFasKECZg3bx6WL1+ONWvWwO124+6778ZNN93kT47q6+sxZ84cvPbaa5g5cyYyMzPP2Ls0YsSIgJMjIiKiaJSTGIucYThbrT+yq0EajLVr12L8+PGYM2cO5s+fj4suugjPP/+8/3m3242qqirY7fYIRklERERKo9ghtkjiEBsREZHyDOb7W9E9SEREREShwASJiIiIqA8mSERERER9MEEiIiIi6oMJEhEREVEfTJCIiIiI+mCCRERERNQHEyQiIiKiPpggEREREfXBBImIiIioDyZIRERERH0wQSIiIiLqgwkSERERUR+aSAegRKIoAji5KzAREREpw6nv7VPf4wNhghSErq4uAEBeXl6EIyEiIqLB6urqgslkGvAcQQwkjaJefD4fGhoaYDAYIAjCac9brVbk5eWhtrYWRqMxAhFG3nB/D4b76wf4HgB8D4b76wf4HgDyeg9EUURXVxeys7OhUg1cZcQepCCoVCrk5uae9Tyj0RjxD0OkDff3YLi/foDvAcD3YLi/foDvASCf9+BsPUensEibiIiIqA8mSERERER9MEEKAZ1Oh0cffRQ6nS7SoUTMcH8PhvvrB/geAHwPhvvrB/geAMp9D1ikTURERNQHe5CIiIiI+mCCRERERNQHEyQiIiKiPpggEREREfXBBElizz77LEaNGgW9Xo/i4mJs27Yt0iGFTFlZGWbMmAGDwYD09HRcc801qKqq6nXOJZdcAkEQej3uuOOOCEUsvccee+y01zd+/Hj/8w6HA3fddRdSUlKQkJCA66+/Hk1NTRGMWFqjRo067fULgoC77roLQHT+/Ddv3oyrrroK2dnZEAQB77//fq/nRVHEI488gqysLMTGxqK0tBSHDx/udU57ezsWLVoEo9GIxMRELFu2DN3d3WF8FUMz0Hvgdrvx0EMPYdKkSYiPj0d2djYWL16MhoaGXm2c6bPz5JNPhvmVBOdsn4Fbb731tNc2b968XudE82cAwBl/LwiCgN/97nf+c+T+GWCCJKE333wTK1aswKOPPoqKigoUFhZi7ty5aG5ujnRoIfHFF1/grrvuwjfffIMNGzbA7Xbj8ssvh81m63Xe8uXL0djY6H/89re/jVDEoXHuuef2en1btmzxP3fffffh73//O9566y188cUXaGhowHXXXRfBaKW1ffv2Xq99w4YNAIAbb7zRf060/fxtNhsKCwvx7LPPnvH53/72t/j973+PNWvWYOvWrYiPj8fcuXPhcDj85yxatAj79u3Dhg0b8NFHH2Hz5s24/fbbw/UShmyg98But6OiogK//OUvUVFRgXfffRdVVVW4+uqrTzv3V7/6Va/Pxj333BOO8IfsbJ8BAJg3b16v1/a3v/2t1/PR/BkA0Ou1NzY24qWXXoIgCLj++ut7nSfrz4BIkpk5c6Z41113+f/t9XrF7OxssaysLIJRhU9zc7MIQPziiy/8xy6++GLxpz/9aeSCCrFHH31ULCwsPONznZ2dYkxMjPjWW2/5jx04cEAEIJaXl4cpwvD66U9/KhYUFIg+n08Uxej/+QMQ33vvPf+/fT6fmJmZKf7ud7/zH+vs7BR1Op34t7/9TRRFUdy/f78IQNy+fbv/nH/+85+iIAhifX192GKXSt/34Ey2bdsmAhBPnDjhPzZy5EjxmWeeCW1wYXCm179kyRJxwYIF/V4zHD8DCxYsEC+99NJex+T+GWAPkkRcLhd27tyJ0tJS/zGVSoXS0lKUl5dHMLLwsVgsAIDk5ORex9euXYvU1FScd955WLlyJex2eyTCC5nDhw8jOzsb+fn5WLRoEWpqagAAO3fuhNvt7vWZGD9+PEaMGBGVnwmXy4W//vWv+H//7//12sQ52n/+33Xs2DGYzeZeP3OTyYTi4mL/z7y8vByJiYmYPn26/5zS0lKoVCps3bo17DGHg8VigSAISExM7HX8ySefREpKCqZOnYrf/e538Hg8kQkwBDZt2oT09HSMGzcOd955J9ra2vzPDbfPQFNTE/7xj39g2bJlpz0n588AN6uVSGtrK7xeLzIyMnodz8jIwMGDByMUVfj4fD7ce++9uPDCC3Heeef5j99yyy0YOXIksrOzsXv3bjz00EOoqqrCu+++G8FopVNcXIxXXnkF48aNQ2NjIx5//HHMmjULe/fuhdlshlarPe1LISMjA2azOTIBh9D777+Pzs5O3Hrrrf5j0f7z7+vUz/VMvwdOPWc2m5Gent7reY1Gg+Tk5Kj8XDgcDjz00EO4+eabe21U+pOf/ATTpk1DcnIyvv76a6xcuRKNjY1YtWpVBKOVxrx583Dddddh9OjROHLkCH7+85/jiiuuQHl5OdRq9bD7DLz66qswGAynlRfI/TPABIkkcdddd2Hv3r296m8A9BpTnzRpErKysjBnzhwcOXIEBQUF4Q5TcldccYX/vydPnozi4mKMHDkS//d//4fY2NgIRhZ+L774Iq644gpkZ2f7j0X7z58G5na78YMf/ACiKOK5557r9dyKFSv8/z158mRotVr8x3/8B8rKyhS3JUVfN910k/+/J02ahMmTJ6OgoACbNm3CnDlzIhhZZLz00ktYtGgR9Hp9r+Ny/wxwiE0iqampUKvVp81QampqQmZmZoSiCo+7774bH330ET7//HPk5uYOeG5xcTEAoLq6OhyhhV1iYiLOOeccVFdXIzMzEy6XC52dnb3OicbPxIkTJ/Dpp5/itttuG/C8aP/5n/q5DvR7IDMz87SJGx6PB+3t7VH1uTiVHJ04cQIbNmzo1Xt0JsXFxfB4PDh+/Hh4Agyj/Px8pKam+j/3w+UzAABffvklqqqqzvq7AZDfZ4AJkkS0Wi2KioqwceNG/zGfz4eNGzeipKQkgpGFjiiKuPvuu/Hee+/hs88+w+jRo896TWVlJQAgKysrxNFFRnd3N44cOYKsrCwUFRUhJiam12eiqqoKNTU1UfeZePnll5Geno4rr7xywPOi/ec/evRoZGZm9vqZW61WbN261f8zLykpQWdnJ3bu3Ok/57PPPoPP5/MnkEp3Kjk6fPgwPv30U6SkpJz1msrKSqhUqtOGnqJBXV0d2tra/J/74fAZOOXFF19EUVERCgsLz3qu7D4Dka4SjyZvvPGGqNPpxFdeeUXcv3+/ePvtt4uJiYmi2WyOdGghceedd4omk0nctGmT2NjY6H/Y7XZRFEWxurpa/NWvfiXu2LFDPHbsmPjBBx+I+fn54uzZsyMcuXTuv/9+cdOmTeKxY8fEr776SiwtLRVTU1PF5uZmURRF8Y477hBHjBghfvbZZ+KOHTvEkpISsaSkJMJRS8vr9YojRowQH3rooV7Ho/Xn39XVJe7atUvctWuXCEBctWqVuGvXLv8MrSeffFJMTEwUP/jgA3H37t3iggULxNGjR4s9PT3+NubNmydOnTpV3Lp1q7hlyxZx7Nix4s033xyplzRoA70HLpdLvPrqq8Xc3FyxsrKy1+8Gp9MpiqIofv311+IzzzwjVlZWikeOHBH/+te/imlpaeLixYsj/MoCM9Dr7+rqEh944AGxvLxcPHbsmPjpp5+K06ZNE8eOHSs6HA5/G9H8GTjFYrGIcXFx4nPPPXfa9Ur4DDBBktgf/vAHccSIEaJWqxVnzpwpfvPNN5EOKWQAnPHx8ssvi6IoijU1NeLs2bPF5ORkUafTiWPGjBEffPBB0WKxRDZwCS1cuFDMysoStVqtmJOTIy5cuFCsrq72P9/T0yP++Mc/FpOSksS4uDjx2muvFRsbGyMYsfQ++eQTEYBYVVXV63i0/vw///zzM37ulyxZIoriyan+v/zlL8WMjAxRp9OJc+bMOe29aWtrE2+++WYxISFBNBqN4tKlS8Wurq4IvJrgDPQeHDt2rN/fDZ9//rkoiqK4c+dOsbi4WDSZTKJerxcnTJgg/uY3v+mVQMjZQK/fbreLl19+uZiWlibGxMSII0eOFJcvX37aH8rR/Bk45c9//rMYGxsrdnZ2nna9Ej4DgiiKYki7qIiIiIgUhjVIRERERH0wQSIiIiLqgwkSERERUR9MkIiIiIj6YIJERERE1AcTJCIiIqI+mCARERER9cEEiYiIiKgPJkhEREREfTBBIiIiIuqDCRIRKZYoili1ahVGjx6NuLg4XHPNNbBYLGc895JLLoEgCBAEAZWVlQO2e8kll+Dee++VNNZbb73Vf//3339f0raJSHpMkIhIsR588EE899xzePXVV/Hll19i586deOyxx/o9f/ny5WhsbMR5550XviD/5X/+53/Q2NgY9vsSUXCYIBGRIm3duhWrVq3Cm2++idmzZ6OoqAjLly/Hxx9/3O81cXFxyMzMhEajCWOkJ5lMJmRmZob9vkQUHCZIRKRITz31FObMmYNp06b5j2VkZKC1tXVQ7dhsNixevBgJCQnIysrC008/fdo5Pp8PZWVlGD16NGJjY1FYWIi3337b/3xXVxcWLVqE+Ph4ZGVl4ZlnngnJMB0RhQ8TJCJSHKfTiX/84x+49tprex13OBwwmUyDauvBBx/EF198gQ8++ADr16/Hpk2bUFFR0eucsrIyvPbaa1izZg327duH++67Dz/84Q/xxRdfAABWrFiBr776Ch9++CE2bNiAL7/88rQ2iEhZwt/PTEQ0RBUVFejp6cH999+Pn/3sZ/7jbrcb3/ve9wJup7u7Gy+++CL++te/Ys6cOQCAV199Fbm5uf5znE4nfvOb3+DTTz9FSUkJACA/Px9btmzBn//8Z0ybNg2vvvoqXn/9dX8bL7/8MrKzs6V4qUQUIUyQiEhxDh06hPj4+NNmo1155ZW48MILA27nyJEjcLlcKC4u9h9LTk7GuHHj/P+urq6G3W7HZZdd1utal8uFqVOn4ujRo3C73Zg5c6b/OZPJ1KsNIlIeJkhEpDhWqxWpqakYM2aM/9iJEydw+PBhXH/99ZLeq7u7GwDwj3/8Azk5Ob2e0+l0aG9vl/R+RCQPrEEiIsVJTU2FxWKBKIr+Y//1X/+F+fPnY+LEiQG3U1BQgJiYGGzdutV/rKOjA4cOHfL/e+LEidDpdKipqcGYMWN6PfLy8pCfn4+YmBhs377df43FYunVBhEpD3uQiEhxLr30UjgcDjz55JO46aabsHbtWvz973/Htm3bBtVOQkICli1bhgcffBApKSlIT0/HL37xC6hU//7b0WAw4IEHHsB9990Hn8+Hiy66CBaLBV999RWMRiOWLFmCJUuW4MEHH0RycjLS09Px6KOPQqVSQRAEqV86EYUJEyQiUpyMjAy88sorePDBB/HEE0/g0ksvxZYtW5CXlzfotn73u9+hu7sbV111FQwGA+6///7TVuN+4oknkJaWhrKyMhw9ehSJiYmYNm0afv7znwMAVq1ahTvuuAPf//73YTQa8bOf/Qy1tbXQ6/WSvF4iCj9B/G4fNRFRlLrkkkswZcoUrF69OuT3stlsyMnJwdNPP41ly5b1ek4QBLz33nu45pprQh4HEQWPNUhENGz86U9/QkJCAvbs2SNpu7t27cLf/vY3HDlyBBUVFVi0aBEAYMGCBf5z7rjjDiQkJEh6XyIKHfYgEdGwUF9fj56eHgDAiBEjoNVqJWt7165duO2221BVVQWtVouioiKsWrUKkyZN8p/T3NwMq9UKAMjKykJ8fLxk9yci6TFBIiIiIuqDQ2xEREREfTBBIiIiIuqDCRIRERFRH0yQiIiIiPpggkRERETUBxMkIiIioj6YIBERERH1wQSJiIiIqA8mSERERER9MEEiIiIi6uP/A/569BY/sHCEAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.fill_between(\n", " np.rad2deg(solver.angles), *np.percentile(Ay, [5, 95], axis=0), alpha=0.5\n", ")\n", "plt.xlabel(r\"$\\theta$ [deg]\")\n", "plt.ylabel(r\"$A_y$ [dimensionless]\")" ] }, { "cell_type": "code", "execution_count": null, "id": "621b5eb4-fcef-45d1-bd8b-ed1e36df19ab", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "id": "84abe41a-3939-4439-bebf-a8a55dfc598f", "metadata": {}, "outputs": [], "source": [] } ], "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 }