{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Octahedral site of FCC Al\n", "Author: Jonathan Frassineti\n", "\n", "This example shows how to use the Celio approximation as implemented in UNDI.\n", "The simulations that follow reproduce the case of Al for different external magnetic fields along the [1.,1.,1.] axis.\n", "The comparison between the case with positive electric field gradient (EFG) and negative EFG is shown.\n", "For the physical aspects, see A.Yaouanc and P. de Réotier, Muon Spin Rotation, Relaxation and Resonance,\n", "Oxford University Press, 2010, chapter 7, Fig. 7.8." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "try:\n", " from undi import MuonNuclearInteraction\n", "except (ImportError, ModuleNotFoundError):\n", " import sys\n", " sys.path.append('../undi')\n", " from undi import MuonNuclearInteraction\n", "import matplotlib.pyplot as plt\n", "import numpy as np" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Define physical constants." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "angtom = 1.0e-10 # m, Angstrom.\n", "a = 4.0495*1.05 # Al lattice constant, in Angstrom, distorted by the muon presence.\n", "elementary_charge = 1.6021766e-19 # Elementary charge, in Coulomb = Ampere*second.\n", "epsilon0 = 8.8541878e-12 # Vacuum permittivity, in Ampere^2*kilogram^−1*meter^−3*second^4.\n", "h = 6.6260693e-34 # J*s, Planck's constant.\n", "hbar = h/(2*np.pi) # J*s, reduced Planck's constant.\n", "\n", "Quadrupole_moment_Al = 0.1466e-28 # m^2, quadrupole moment of Al.\n", "omega_E = 0.70e6 # Quadrupolar resonance frequency of Al, in us^-1." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Define the rotation matrix that brings the diagonal of the FCC lattice\n", "along z." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "R = np.array([[ 0.78867513, -0.21132487, -0.57735027],\n", " [-0.21132487, 0.78867513, -0.57735027],\n", " [ 0.57735027, 0.57735027, 0.57735027]])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Define the atomic positions for muon and Al nuclei." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "mu_pos = np.array([ 0 , 0., 0.])\n", "Al_pos = np.array([[ 0. , 0., -0.5],\n", " [ 0. , 0., 0.5],\n", " [ 0.5, 0., 0.0],\n", " [-0.5, 0., 0.0],\n", " [ 0. , 0.5, 0.0],\n", " [ 0. , - 0.5, 0.0]])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Define the atomic structure. Note that positions are rotated first and then\n", "transformed into Cartesian coordinates and SI units." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "atoms = [\n", " {'Position': R.dot(mu_pos)*angtom*a,\n", " 'Label': 'mu'\n", " }\n", " ]\n", "\n", "for i in range(6):\n", " atoms.append(\n", " {'Position': R.dot(Al_pos[i])*angtom*a,\n", " 'Label': '27Al',\n", " 'ElectricQuadrupoleMoment': Quadrupole_moment_Al,\n", " 'OmegaQmu': omega_E\n", " }\n", " )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Polarization function.\n", "\n", "The muon polarization is obtained with the method introduced by Celio." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "steps = 100\n", "tlist = np.linspace(0, 20e-6, steps) # Time scale, in microseconds." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Define the applied external magnetic fields, in Tesla." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Computing signal for positive EFGs...done!\n", "And now for negative EFG...done again!\n" ] } ], "source": [ "LongitudinalFields = np.array([0.0, 0.001, 0.00164, 0.003, 0.00329, 0.004, 0.005, 0.006, 0.0005, 0.002, 0.0045, 0.0055])\n", "\n", "#This is the Celio method for the polarization function with positive EFG.\n", "print(\"Computing signal for positive EFGs...\", end='', flush=True)\n", "\n", "signal_positive_EFG = np.zeros([len(LongitudinalFields), steps]) # Define the polarization signal with positive EFGs.\n", "\n", "for i, B in enumerate(LongitudinalFields):\n", "\n", " # Align the external field along z i.e. parallel to the muon spin.\n", " B_pos = B*np.array([0.,0.,1.])\n", "\n", " NS_positive = MuonNuclearInteraction(atoms, external_field=B_pos, log_level='warn')\n", "\n", " signal_positive_EFG[i] = NS_positive.celio(tlist, k=3)\n", "\n", " del NS_positive\n", "\n", "print('done!')\n", "\n", "#This is the Celio method for the polarization function with negative EFG.\n", "print(\"And now for negative EFG...\", end='', flush=True)\n", "\n", "# Align the external field along z i.e. parallel to the muon spin.\n", "B_neg = LongitudinalFields[2]*np.array([0., 0., 1.])\n", "\n", "# change sign\n", "for i, a in enumerate(atoms):\n", " if 'OmegaQmu' in a.keys():\n", " atoms[i]['OmegaQmu'] = -a['OmegaQmu']\n", "\n", "NS_negative = MuonNuclearInteraction(atoms, external_field=B_neg, log_level='warn')\n", "\n", "signal_negative_EFG = NS_negative.celio(tlist, k=3)\n", "\n", "del NS_negative\n", "\n", "print('done again!')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Integrate the different curves to reproduce Fig.7.8 (right)." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "area = np.trapz(signal_positive_EFG, x = 1e6*tlist, dx = 0.2e-6)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "...and the results are:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAeoAAAEPCAYAAACN5IcLAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAACB/ElEQVR4nOydd5wdVd3/32fm9rJ3e69pu+khCQkkEDrSEWkKIkVBBXwQGyo+PvpTH7vYEFQeRBRRkd5rqCFAQgmkt03ZbK+3tzm/P+buze5mN3Wzd3dz3nnNa2bOnJn5zuTufE79foWUEoVCoVAoFKMTLdMGKBQKhUKhGBol1AqFQqFQjGKUUCsUCoVCMYpRQq1QKBQKxShGCbVCoVAoFKMYS6YNOBzk5+fL6urqTJuhUCgUY4aVK1e2SSkLMm2HYk/GpVBXV1ezYsWKTJuhUCgUYwYhxLZM26AYHNX0rVAoFArFKEYJtUKhUCgUoxgl1AqFQqFQjGLGZR+1QqFQKDLHypUrCy0Wy13ADFSFcF8YwEeJROJz8+bNaxksgxJqhUKhUAwrFovlruLi4qkFBQWdmqapgBJ7wTAM0draOq2pqeku4LzB8qiSjkKhUCiGmxkFBQU9SqT3jaZpsqCgoBuz9WHwPCNoj0KhUCiODDQl0vtP6l0NqcdKqBUKhUKhGMUooVYoFArFuOM///lPVnV19YzKysoZ3/72t4sHy2MYBldddVVFZWXljClTpkx7/fXXXQdzrwcffDCrrq5uWl1d3TSXy3VUdXX1jLq6umkXXHBB9SE9RAol1AqFQqEYVyQSCW6++ebKp556asOGDRtWP/jgg7krV650DMz3wAMP+LZs2eKor6//6I477th2/fXXVx7M/S688MKedevWrVm3bt2aGTNmhO69994t69atW/Pwww/XH/LDoIRaoVAoFOOMl19+2V1VVRWdNm1azOFwyE984hMd//nPf7IH5nv00UezL7/88nZN0zjllFOCPT09lm3btlkH5nO5XEd98YtfLJs+ffrURYsWTVm6dKlrwYIFteXl5TPvu+8+3+F+HjU9S6FQKBSHjaW7Wis6IrGDalIeilyHLXRSacGOoY7v2LHDVlZWFuvdLy8vj7311luegfkaGxut1dXV6XwlJSWxbdu2WauqquJ984XDYe2kk07y33HHHQ2nnXbaxO985ztlr7322oZ3333XcfXVV9dcfvnl3cP1bIOhhFqhUCgU4wop9xxwLoTYI3GIfHukWa1WedFFF/UATJ8+PWy32w273S4XLFgQbmhosA2L0Xsho0IthLgbOAdokVLuMYdMmG/sN8BZQAi4Skr57shaqVAoFIqDZW8138NFZWVlrK+A7ty501ZaWhofmK+0tDReX1+fztfY2GirrKzcI5/FYpGaZvYUa5qG3W6XALquk0wm91T2YSbTfdT3AGfs5fiZwOTUch1wxwjYpFAoFIoxzAknnBCsr693rFu3zhaJRMRDDz2Ue+GFF3YNzHfeeed13XfffXmGYfDiiy+6vV5vcmCz92ggozVqKeWrQojqvWQ5H7hXmu0Ty4UQ2UKIEill496u6w9FhtNMhUKhUIwhrFYrv/zlL7efccYZU5LJJJdddlnb/PnzIwA/+9nPCgC+8Y1vtF5yySXdTz75pK+qqmqG0+k07rrrrvqMGj4EYrA2+hE1wBTqJ4Zo+n4C+ImU8vXU/ovALVLKFYPkvQ6z1o29eNK8f9//FTp3ZFFVO4sTF8w8rM+gUCgUYx0hxEop5fzhuNYHH3xQP3v27LbhuNaRwgcffJA/e/bs6sGOjfbBZIO1/Q9aspBS/gn4E0BB5URZVhBmxrEtPP9uJ5sfWcVV512Iru0xjW4PErKHkLGBuOwiKXtIyB504cYhKnBoldhEMUIcWI9BQvqJyVYSspOE7MQgjgUfFpGNVWRjFQVoYo8ZAUMipUGSIIYMkySEIaNowoaGE004sOBFCP2AbFQoFArF6GS0C/VOoKLPfjmwa18nVRXmMHPyjfzgjmd5rF1ydF0XxzTfRqHzJHJ9s9BF/5kChowQNDYQSK4iLLdglgU0LGQhkh4iiZ0E7WsgCUbUimyvIdszk7y8WjRtT4GV0iAmmwka6wkZ64nJpj5HNQRWJNF+aTZRiF2UYBV56MKFhgshrBgySEIGSOInLjuJyzbishNI7uUNaFhFHjZRgFXkYxW5WEUeVpGDwJqyQWAQIyF7zAIJPSSln4T0kySAlLuvL4QFHRe68KDjNu0TLnRciD4/IUkSgzBJGcYggpFaJ2UY2cdeAQhhR8OBjiO1bUcTdjSs9C2fSRIYMo4kjiSBJIkkAdIAIVLPoiHQEVgQWBHCkto3l93XE6n/W4lMrfvv9+bqza+l9nrvI3qtTy/CfEED7jHwOgPpew0GXLfvmWLAOYNdZ+DdBlx7kBGsCoVibDHahfox4EYhxD+BhUD3vvqne7HZbPzgpnPxPfYqkfeyCE1qp2fWs/REn8UminHpNSRlmKjcSVyaLTQWfGTrx+OSM2h4VfDeAy1sWdqBEZc4iiLkze2meEkHxSduwZ+1gc4unWR3Dt68PJyuXKSME5VNxGQzEnM8gl1UkKufil0rwypy0PEihEZSRkjILuKyg5hsImrsImiswyA8xBPpWEUOVpGHS5uCVWSj4UIXTgR2JDEMGSFJmITsJCZbickmgsZahmiEGOIubnThSQm6iZQBorKBJMEDupaJhoYjJaKmaEgMDKJIYvs4VzE8DFYI6Htsb+cMVmgYuC8YWIjpty/2p+AxWKFl4HUH5hvM5sH2B9qz53125+t/nslQv/n9fwcCLVWg65t3MJv3dp++1+xNH+paQ19XDJpPFehGM5mennU/cCKQL4TYCfwPmAohpbwTeApzatYmzOlZVx/oPb523hJemryKH905j+ADQT6e08Exp7aTPfMtdOHArpXh0WdgTVbR8qaXlc93sOn5zYTa4rjyrBz1mVKqj8/BV+HAW2pHaNC6rofWD9YT0TeQ1DuJRXfgLt+Epuk4LMV49bnYRSkubSK6cA9qly4cRNpy2bFMo229k7aNRXRsmoIUcVwlSVzFSZx54CvKJa86j8Kp+WSV2A+ohmQkJe2be2jZ3ER3WyORSAfxSIxEJIGUBomwTrjJQbzLgYx4cLqyyCpzk13lIKvUjqfIjqfYhsNnQegCNDAIEQ0ECAf8RCN+YoEYkZ4E0e4k0R6DSIdOpM1CpE0nGbKCYUW3aVhdOs4cK85cK84cC65cK85cHWeBxOpJoDkSaI44Qk9gJCXSkBhJiRHTMWI6yYiGEdOQCR0joYGhoVlAs0l0G+h2icVpoDsMdIdEsxkI3Ui1DPR+aCX9PpwA6W4M8/O1+5MszVp7urZt7KUWLun/MZeDfNrlINsDavayb9rAcwbb3522py29KUa/dDnIuYNfs+91B+STez6LZOB76bsMtG/w+8k+2+k02ff+A/MNRe8xI33v3ecYqVc9mH17s3Ff9+pra/93IDH6P1O/5zyA6/c7TwWmOpLI9KjvT+3juARuONT7nDx1FvqXNnDNnzawLDGL0Fei6K1BZFxgcei48qyEu3YRDyaxunVqluQw9YJCak7IQbf274/uDkmafV7q43PZ2TmHZEsUHmgh8WozMpjEXe5kzqVFTPt4IVqpvd+5ke44LWuC7Hy7m61LO2haFQBAswhyJjgpnO7F5tGJBZLEOpI0r42xemMPMtkDbMWVZ6VohoeimR5yJ7pSomfF5tGJdCcItcUItcVp2xCk+cMAzasDJMLmx8rqtlE0fSJZZQ48+VZceTa0fEE0O0HMnyTcGadre4T6VzoIth7o7AQNq9uKK8+0xebWsbl0LC5IxiXJuIF/V5SWNQHC7XGSsZH5yFjsGhanhm7T0K0C3aahWQRCgNBE/8pUH8yW7N35hLZ7relid5quDZqn71qziPRa00BYBJomELq5aNpuW3a3oIs97YHUPVPX01PXtAh0m0CzaOgWgWZLpevmInTRr2VeiP73Mm3d/Wy95wqLMK9n2f0MvfdP2z/YM6um9hFl92DgfRWEhtrqu/7esNqmGD5Ge9P3sHHCpClcceI7PHD/Wgq+dDrFm8NUGzEKSBLrimN1aEw4JY/KRdlY7P3FOZaQrN1l8P6OJI1d5g/bY4fKPA1rgZP45Coinyxn59J2/G+08sYvt/HGL7eh2zU8RTbcBTYCTVF6GlL90gJKZntZ/JUqak7IIb/OvUeBoJd4JEnbuiBNqwKm+H7kp/71TuReuqgtDo3C6R5mXlpM8UwPRbO85NY4TUHYD+KhJP6mKIHmGIHmGJHueN8KJvYsHWeuFUe2Kc7uAhs29/4NXpNSEg8ZhDvjhNvjhDrixIJJEhGDRCSJkZBpcdB0gW7XsDg0rA4dzZoSEQ0QApmUJGMGyZgkEUkSDxvEQ0liIfN68bBBIpw0Cwsxg2TMQCZBGhIpzfWe9qU2DNNWKUEmJdIgXcuXcbPGL43UMVLXTO6+bu99jGTqnITEMDDXqRYDmTTvYSTN2rQcLxWlVEGgbwFECLNFJl2w6V3rpAsevfu9BYze/+vewokQ5jFSBYJ0YaNPwaqfDaTS+xUwBhSgdHYXRFKFmt7r9X2O9PNopAsqmk4fO3fbDrvP6duy3LcQM7Bw07eA1LfAZNpG/8JSv4Jd6ppDzd5JdTuYBcOU7YPcQ7Nk2qWGYm8cMUIN8D9nXsbbb93Io0/CxKNncuqECTRZ4diJOnOrdaz67j+kSFxS32awpdVgfaNBNAH5HsFJU3UmFmjkecSA2oOVxKJSNjYX8+6KIA2vdaL7Y7hkAq0nTvEcL7MuK6FouofCGR5cufs3ytvq0CmZk0XJnKx0WjycpKchaopdZ5xYIInDZ8GVb8OVa8VbYkezHHzNxurSyZ3gInfCsLrnBcyPks1t1rp95fseha8wkXJ3wcBIibyRMDDikmRcYiTMlgtznSo0GGa6NPpfp2+hIF1oSfYvVCQT5jV6r50uhPTm6S3ApAsqpPNAb8EG6C3syNR2qsCDJP0c6QLQwMJNnwJS3+unC0XQ51lkWqvSzyZBylRTdzJVUEpd1zAkMmHaYL6j3c9onttnu89z9LWrt4ClGJ1cfPHF1S+++KIvLy8vsXHjxtWD5TEMg2uuuabipZde8jkcDuPuu++uP+6440IHeq8HH3ww69Zbby0H2L59u72wsDDucDiMqVOnhoYjgtYRJdRCCP7w1R9w2q3/YeOWDaz1J7mitpSl69y8vC6J0wYum9kc2dpjfgjsFphUpHFUpU5ZzkBx7k8gliRgRBFVCTSfG+KFbI9qVOUJFs2wkOcZnlKr1amTN2n4RVQxehHCrBGhm+PYTdQUvEyTLmT0FqJ6W0fo02LTpzDRd79focMAIyHNwleSPQpM6ULXgAJWXwZ+mtK361PQMVL3M5J97pcwr33gI4BGN9dcc03bTTfd1HL11VfXDJWnb5jLpUuXuq+//vrKVatWrTvQe1144YU9F1544RqABQsW1P7iF7/YsWTJkgMW/KE4ooQaoDorly9eMYF7r/8ZkaJL+Vl4DudV+bhxUSWJpCAUM5u6j52kU1OgUZot0PfSZLyjO8bf32/n7x90sLa1v0c0ASwoyuVjySLqXzaYXiE4Y7oN2yHUdhUKxehBaHv/PowpxplQn3nmmYH169fvNWDGUGEuB7oRdblcR1155ZUtr776apbP50v+6Ec/2nnLLbdU7Nq1y/bTn/50u4qedRj48tGn8PKVT7Pz9r9xlC543JjF+q4gv7ughk/MzNrn+dGEwSNru7hrRRsvbfUjJRxX5eHHp5VRm+9gcp4deyTKys2dvLsjwNubu/HkFCNEDqt2hDmpTueYCbbx8weuUCgUQ/CzwMMVW5Mtw9oEWKMXhr7hueCQg32oMJejGE0IvnnZf3H1R/8k4XKzZPkK3pk1ndPv2cjpk7L48WllzClx9mvmDsaSrGgI8cT6bv76XjttoQTV2Ta+e2IJn56Tx8RcOw2bWnn9kRXc+egq6lfvnu6tA+1eF4+cfTwLT5zPa+udvL4hwjETLMyptOBzDS3YUkoCUegKSSJxSTQO0QRoArKckOUUZDkFdlVLVygUigNChbkc5RybW8nUkz1s+L/nWbyolqmvrOAV3csbzGbeph48No2qbBtV2TaaAgk+aAqRNMCiwfl12Vx7dAGnTvASjyZ47eH3uf3uN9n4rlnAm7qgiqu/dzZF1bnkFHrJynPTVN/OmuX1vPHPh3ho6gymTa3BMLy8udnAbYdspyDbJdA1CMfNwWzBKHSHJUljHw8DlGULJhdrTC7Shq0vXKFQKA6V4aj5Hi7GSpjLI1aohRB8/4yLuW7Zizz5+MM88K/7mHbHSp6462GSp8+jcNFkgm4b27vj5Los3HJ8McdWeDi2wk22Q2PDyh385X9e4cV/rMDfGaJiSiGf/eG5HHf+LPLLsve4X/nkQuafNpXPAFs+auQLd7zNYwVlzM32cGalF92Tw45OA8MAp03gsEKBVzC5SCPbJfC5BC4b2C0CmwWSBvgjkp6wpC0g2dRs8PK6JC+vS1LiEyyYoFNbrKGp5nWFQqEYlPPOO6/rD3/4Q+G1117bsXTpUhXmcjQyK6uAhZ/9BC9sXMpldz7HW7ddx4zFE7jzG4/Q8/hysgs9fPGCOZTXFhLf1Ex4dYJ/1Lfz9jNr6Gz2o+kaC8+cxtmfW8ys4yfut7OHCTNKeO72s3nguU1c91wTS1s7OefND7j/f07A5bHv+wIpspyCshxz+/gpZu17Q1OSd7cZPPpegiwnzK/WmVOpqwFsCoXiiOLcc8+tWb58ubezs9NSVFQ065vf/Oaum2++uU2FuRwlzJ8/X65YsUckzEFZ093BZX+4A399AXd+aiannXgssUicFc+v4+UH3uWd59aSiO32LuL02Jl7ci3HnD2d+afV4ck+tDES7f4YZ/3iHd4xHFS0tvHANVNZsKDykK4ppWRTi8HbW5Ls6JA4rDCvWmd+tY7TpgRboVDsiQpzmVn2FubyiBdqKSXf3PwEz33hV5R7CnjkwfvR9d3zU8OBKCF/BJvDis1hxWrX6e2rGC6klHz3/rX8+KMglliMm6os/PiGo4flPg2dBss3J9nYbGDVYWqJxswKnfJ9zAlXKBRHFkqoM8vehPqIH3UkhODT+YvI+vw0tm9p5Ld/e6LfcafHTl6JD2+OC7vTOuwi3WvDDy6bxiufmUA+SX7eYmHGN19nfUPgkK9dlqNx4Xwrn11iZWqJxtpGg/vejPPHl2MsXZtgXWOS7rAcdPSjQqFQKDLPEd1H3cu0rFzmLDyFp5dM4ncfxrmiq4f87H3Ppx5uFtflsfUni/ncT97gH3EbM25fwycmefjOWTXMLHYe0rULvBpnzdY4dbpkfZPBhzuTvLM1Sa/zJIfV9Mpmt4LDkvIJnPKXrAmwWsCmg81iTgfLcYnUIDdzuptCoVAoDg9KqAFdE5zvWMjbZ7yL/yev869727nhvz6fEVusFp2/fmcJn3xmHTf9bQ0PyTIeuH0NJ9V4OH9qDgvL3cwpcWLfhxN9KSWxpMSQ4OwT8MNmEcws15lZrpNISlr9kl1dkja/QSQ1LSyS8oHc6/XQkBBPSuIJcw5331gWNh1KcwSl2RoVuRrluaKfz3SFQqFQHBpKqFMsyCli1uwZvDf3I+796z/4xEXnU1JanDF7zjyjjsULKvjVNx7l31vDvBupY+lWsyncpgtqcmy4bTpOi8CmawTjSbojSXqiBoFYkmDM2F1btgjyXRYK3BaOKnVxQrWXE6q9VGbbKMkWlGTD/vqN7nXA0hmUdIYkzd0GDZ2SNzclWUYSXYPKXEFNgUZNgUb+HsFLFAqFQnEgKKFO4bFaOFXMY93n36V1R4JP/uZZXvnplRm1KSvXzffuuoyTH3yfO7/xMK0xScl58/AcP5V2LIQTBqG4QSwpyXVaqM6247XreO0aLqu5CKAjnKQ9lKDRH+eh1V3cvbIdgJlFTq6dn88Vc/LwOfZPqIUQeB3gdQgq84AK87xoQtLQIdnSZrC11eCltUlYm8Rjh5oCjao8jfJcDZ9TxSxWKBSKA0EJdR+OzS7nKWbQMWcHuzY1s2LlKubPm5Vps1hy4RzmnjKFh29/lcfufI3uf73B4vNnsejcmcw9pRaXd//DRSYNyfu7gjy/rpMH1vbwX0/u4JvPNfDJWTl847hipuQfXOhJu0UwoVAwodBsZu8OS+pbDba2GWxsNvhwp+lezWOH0hyNAq8gzy3I85jOXOwWJeAKhWL42J8wlwBPPPGE92tf+1pFIpEQOTk5iXfeeWd977FEIsHMmTOnFRcXx5YuXbrpYOyYNWtWXSwW07q7u/VIJKIVFRXFAR599NFNtbW1sX2dD2p6Vj+klNy59UP+af0nyStWUFVUzt/uv+uwjPQ+WLpaAzz4m6W8+M8V+DtCWGw6MxdPpHRiPrnFWeSV+BCaIByIEglG8XeGaWvooq2hi/bGbvwdIYI9u6N8dRfksGPmZHZOrsLQNY7Ro9w4w8s5p07CmzM8fvQNKWntkezslDR0GuzqMugO9Q/lq2vgtoM75ZXNYTXXTpvpkc1lE3gcAp/TrNErj2sKxfAy3qZnPf300x6v12tcffXVNUMJdVtbm75w4cK6Z555ZuPkyZNjDQ0NlrKyskTv8e9973tFK1eudAUCAf1ghbqX3/72t3krVqxw33vvvdsHO7636VmqRt0HIQTHZ1fxarSStut6+PCX6/nFvU/zjavOzrRpabILPHz2h+dy1ffOYu3b23jrqdW8/8pGNr63g0BXeI/8FptOfqmPvNJsaudXkZXrxpvjxJ3tQksJesgfYVvLNh71W3i7pITLP4hR+tCznJXs4ZyTJ3DMOTMoLM85aJs1ISjyCYp8puMVgERS0hE0l56w6dc8EJWEYpJIHLrCZhCSSKy/oIMZdzfLAfkejTyPIM8rKPAK8j1CeWBTKBTA/oW5vOuuu3LPPvvszsmTJ8cA+or05s2brc8++6zvW9/6VuNtt91WNNj5F154YbXD4TA2bdrkaGhosP/xj3/ces899+SvXLnSfdRRRwUffPDB+uF4FiXUA5js81JXP52lpzcgl1/EH1cF+HxPAF+WJ9Om9UO36MxYNIEZiyak0yKhGJ3NPUhpzv92euzYndb9blL+LtDQHub/PbGVv4ly7kLjiVe2M/mXf+DY2lxOungui8+fecje2AAsuqAwS1C4j1lwhjSFOxQ1B7F1hyTdYXMgW3tAUt9u9Ataku0yp6LlewT5XrNpPdetBFyhyBT/850fVmzauHlYw1xOmjwx9P0ffueQg31s2LDBEY/HxYIFC2qDwaD2xS9+seXGG29sB7jhhhsqfvazn+3s7u7e6wCe7u5uy5tvvrnhH//4R/all146+aWXXlo3b9688KxZs6YuW7bMuWjRoj1rUAeIEuoB2HWNY+2Tec94m5xzutj07Yf56/+18183X59p0/aJw2WjpCb/kK5Rlufkj1dO439DCW5b1sxvrYLXJlWypb2Dl36zjNJbHmbB6VM54aKjOPr0qdgc1mGyfnA0sbvpO9+753FDSrpD0OI3aPWbTextAcnmFqPfNDK3HXJcqaZzp8DrMJvQPXaB2y5w283Cg0KhOHJIJBJi1apVrtdee21DMBjUjjnmmLolS5YEVq9e7cjPz08cf/zxoSeeeGKQL89uzj777C5N05g7d24oLy8vvmDBgjDAlClTwps3b7YroT5M1PmymNQ2lZXzl3PiaXO49577OPH0U5k1fUqmTRsx8lwWfnhqGV9dXMSf3mnj9resrMjLpUAmqF+9mRduepgC4wEWnjWNhWdOZ+7JtTgPIKDIcKEJQY4bctw6tX1m0yUNs1m9PSDTU8k6g5KdnQb+xv5zwXux6eC0mf3iThs4+/STO1P95s5UoaG379yqq0FwCsXeGI6a7+GivLw8lp+fn8jKyjKysrKMhQsX+lesWOFauXKl6/nnn88uKyvzRaNRLRgMaueff37No48+unXgNRwORzrkpc1mS39ZNE0jkUgMy8dBCfUglHucTG2cwkfyPbw3TiC+3ckn/28FK39cg9N+eGuQo40cp4VblhTzlcVFPLymkzvebuUVUQszaqlIRFj30RYe/voT5PvvZ/Zxk5h1/ERmLJ7IxFml6Jb9m/J1ONA1s9+6YJCysJRmn7g/IglGJcEYBKOScAxCMUk4Zm53BnudwAx9H4sGLjtmrdxm1szddpGqqae2Hea2cgSjOBwYhiSaYLfDorg5XTISh2jcPBZLQjwhiSchnjTHiMQNSKT248nxN6h4f7jooou6brjhhsp4PE4kEtHee+89z9e//vXma665pvP2229vAHNU+C9/+cuiwUR6pFBCPQiaEEzz5TAhOJm3PKv5+GlH8c/H3+aBf1v4zBWXZNq8jGDVBZfMzOWSmbns6I7xrw87+NeHnay0OGDONNzS4IOWNhwPbyLnzrcpCgWpm1nChJmlTJhVRvW0YkonFOBw73Vsx4gghMDjAI9j/4TTkJJoHFPA46aYh2LmfigKwZgkFJX0RCSN3WbaYJ89m4V0M7u7dzS7XeBKj3Q3a+12i+lBzpZy26pGuI9dDCmJJSCWMMWz73Y8CfEExJKSRBISKeFMGpKkAUlpxp1PGmAYkDAkCaNXaCGWMEU4nty3Hb2/JasusOhmS5BFM1uHrPr4LETuT5jLuXPnRk499dTuurq66ZqmccUVV7QeffTRkX1de6RR07OGoCMS4y/163iq4EEuchzLhzc/w/vvreKhx+6nuGTQAYBHJM2BOC9u9vPC5h5eqfeztdOcFqhJSW4ohKuhBW9jO9ktHWS1dVFQ5KV0Yj6FFTnppaAih4LybPJLfYe9z3skMORuAQ9GIdBbc4+aNXdzdLs5QC68HyHqNZH6sKY+rlZdoGt903Z/fK0pf+zm2ty2W8x57uZ+qhCQOq6JI6Pp3pBygOiZ3R+JpCl+yVRaMpkSxD7C2SuQpoj2pss+x1M11D6103jSFOT9EdFedM38/9VTiybMlqHefV2Q+n8WWC3m/7XdYvrnt1t2T2nsXfces+1n98x4m5411lDTsw6CXIeNKlsuNfFqntbe49ffvZmLr/g6n/z5I7zw82ux6KNnbnUmKfJYuWx2LpfNzgWgJRDnrZ1B3toZ5N1dIVYW+Ng+pQYAK5KyeIT8tk6y1u3C+uAHWGL9lSq70ENhuSneRZU5FFXlUlyVR3F1LgUVOVhto/8nqx1Ajd0wTLHubbKMxGWqxmXWmHZ//IcWh3BM7iEQg/XBD4ZgdwHAou8WB1PAzf3ebZE6QfQ9mT777M7Tu+otCIg+10lfk9R1B7wmKVML5nMY0hRXc9sUXCl31zh7jyUNmc7fm96bZ7joFVKr3iusuwtJvS0mtrSQ9raO7C4w2foVmkS6cKXGOij2xuj/6mWQWp+XNR21bLZtZW1+GydccAkPb0nyx3sf4oarL8q0eQC0tbaz/M232bFjJw07G2lsbMRqtVJcXERRcSE1NdUsOWExLvewzo4YkkKPlXPrsjm3Lhsw+4N3dMd5e2eQ5TuDLN8RYKXTRbSkFH3xfGbn25jvFcyUEXLaOmlv6KJlRydbP9rF28+sIR5NT2tE0wQFFTmU1ORRMiGf0gn5lE7Mp2xiAUVVuVismesTP1g0bXdf9nCSNPoLfr+m11SzaTxBn6ZUuVv4UuvegCxS7g7Q0tsAly4HDNg3endSQitT5xu9wpu6vnnMTN9zonwfgccUxF6B1wVoKbG0aWahSNN600WqoNG/Nto3TUuJq96n5mrR+otub6Gld79XmI+U1gfF6EMJ9V6Y5PNQ1FxEsZHPQ9Hl/PH6zxP6r1v5y6/fYMm8qcycNT0jdnV0dPLc0y/w/HMv8e7K95FSIoSgsKiAkpJiopEob2zcQltbO1JKHE4HJ528hLPPPYNFi48ZUU9rQggqs21UZtu4aIbpNCUcN3hzR4CXtvh5cbOfP28NIqWg0F3ImSdM5pxaH6dPysJtFXQ2+Wna1k5TfQeNW9tp3NpG49Z2Xn3wfYLdu2c9aLpGcXUuZZMKKJtYQPnkQsqnmGtf/uiaAz8S6Fpq5LoN+td5FQrFWCOjfdRCiDOA32CGbrpLSvmTAcd9wN+BSsxCxS+klH/Z13WHo4+6l6e3N/GGXMUb3tf5tfcaqoM5XHrxlXTl1XLfL7/GxLKCYbnP/tDc1MI9f/k7D//nUSKRKBMnTeC0j53MSSefQM2EKmy2/gO14rE4H320hqeeeJbnnnmB7u4epk2v4xvf+gpzjsq8D/Ne2oIJnt3UzVPru3lmUw+d4SQ2XXBSjZeza32cXeujJqf/1C8pJT0dIRo3t9KwuY2Gza3s2tRGw6YWdm1pIxbZXRP35rqomFJERW0hlbVFVNQWUVlXRG5xlqohKRQpVB91ZtlbH3XGhFoIoQMbgNOAncA7wKeklGv65Pk24JNS3iKEKADWA8VSyr06Mh9Ood7Q7ee5hkaeKf4P860T+Z73Ul588wM++9BWJsd38OyvbjjsNdTW1jbuvP0uHn34CUBy9rlncsVVlzFp0oR9nttLPBbn2Wde4Le//gMtza2cfe4Z3PSVGygsHLmCxv6QSEre2B7gsXVdPLG+m43tUQCmFTg4c4pZ0z6+yoPDOvQ7NwyD1p1d7NzQws6NrezY0MyO9S3sWN+MvzOUzuf2OamsLaSyrpjKuqLUUkxOkVcJuOKIQwl1Zhmtg8kWAJuklFsAhBD/BM4H1vTJIwGvML+aHqADSAy80OGk2uPGLqwclZjGa7xPS7KbU46dzU1rNnHPr+7jlz9L8rVbbjosH/ZIJMLf7vkHd9/1N+LxOJ+46Hyu+uynKS0tOeBrWW1WzjnvTE4+5QT+766/cu9f/sEbr73JD378Pxy/ZNGw236wWHTBCTVeTqjx8sszK9jYHuHJ9d08ub6b3y1v4ZdvNOOwCJZUezl5grkcVeJC7zOFSdM0iipzKarMZd6pdel0KSXdbUG2r2tix/pmtq1rZvu6ZpY98SHP3vtWOp8n25kW74raIiprCymfXEheqU8JuEKhGHEyKdRlQF+PNTuBhQPy/B54DNgFeIFLpZSDjuEUQlwHXAdQWVk5bEbadI0qj4to9yTIe49Ho29zres0vnzNJ4i1bONv/3qUDyI+7v3ulcNWs5ZS8sxTz/ObX91OU1Mzp5x6Il/+6o1UVJan83QGY3gcFqy6xupd3Ty4soEvnDiBQq+DtkCUpu4IdcXePUanu9wuvnTTFzn3vLP4xle/w5e++BWu/uwV3PBfn8diGX1DFibnOfjyIgdfXlREMJbk1foAz27q4YVNPXzzuQYAsh06S6o9nFDj5cRqL7OKnf2EuxchBNkFHrILJjHr+EnpdCklXS1+tq/rFe8mdqxv4fVHPugX6MThtlE6sYCSmjyKq3IpqsqlsCKHvFIfucVZeHNcSsgVilFAKBQSCxcurIvFYiKZTIpzzz2387bbbts1MJ9hGFxzzTUVL730ks/hcBh33313/XHHHRca7Jp748EHH8y69dZbywG2b99uLywsjDscDmPq1Kmhhx9+uP5QnyeTTd8XAx+TUn4utX8FsEBK+aU+eS4CFgNfASYCzwOzpZQ9e7v2cDZ9A2zuCfD8zhY2lb7BJtnAv7O/ik1YkVLyqe/fy/Kgj09mbePH3/mvQ/5Qr/5oLT//yW28/94q6qZO4Wu3fJn5R89NOUGQ2Cwaz3zUyBf+/i6P3LCYORXZPPNRI197YBWP3riYiQUe/r58G9955CNe+8ZJVOS66A7Hcdv0PUQ7Eonw85/8mgcfeISj5s7mF7f9L3n5eYdk/0jS5I+zdKufl7b4eWWrn00dZjN5tkPnuCqPKd7VZo37YPx4pwV8fQsNm1po2NjKzk2tNG/roHl7B4lY/0myVrsFX76HrDw3vjw3WXluvDkuPNlOvDkuvLkuvDlusnJdZBd6yS7wjIt544rxwXhq+jYMA7/fr/l8PiMajYqjjz669rbbbttxyimnBPvm+9e//uW7/fbbC19++eWNS5cudd98880Vq1atWnco916wYEHtL37xix1Lliw5IMEfrU3fO4GKPvvlmDXnvlwN/ESapYlNQoitQB3w9siYaFLpcWERgumR6ay0beTF2IecaZ+LEIL7/vsKvvvzP/Hk3+7HEuri69+8GZ9vHyGhBmH7th386c6/8OTjT5OTm833fnAr555/Frqus709xCfuWMatZ9dxwVHlHDMhj2+fVUdRljnA6owZJZwxY3dz+OnTi8h12yjPcQLws2fWsWxzO8/fvKSfWDscDv77e99k3vyj+H//8798+lOf5Xd/+CWTJk88xDc2MhR7rXxqVi6fmmXO4d7ZHePlrX5eqQ/wWr2fJ9Z3A+C1ayyu9LCk2suJNR7mlbr3yxOTEIKcoixyirKYvWRSv2PJpEFHYw+tOzvpaOqhvbGHjqZuutuCdLcH6WkL0FTfjr8zRLA7wlAFYneWg9ziLArKsykoNx2/lNTkUzIhj5KavGGJVKZQHGlomobP5zMAYrGYSCQSYrBK1KOPPpp9+eWXt2uaximnnBLs6emxbNu2zVpVVdXPwYPL5TrqyiuvbHn11VezfD5f8kc/+tHOW265pWLXrl22n/70p9svv/zy7sP5PJkU6neAyUKIGqAB+CRw2YA824FTgNeEEEVALbBlRK0ErJpGtdfNjm6oLirkochyzrAdhRACXdf44S2fpyLLwu8ffI0Hv3k/PzhzAheee+p+1a63b9vBn//4F5564lksFgufueoyPvf5q3l2QxePrWrkgqPKKc9xctq0IipzzY92tsvGdUuGFtNCr4OzZu4W7lOnFVGT706L9LqmHmr7DJg665yPUVVdyZdv/DpXXn4tP/vlj1h8/LGH8soyQrnPxqfn5PHpOWarQKM/zqv1pnC/utXPt583m8o9No3jqjycWOPllAle5gzo494fdF1LiWv2PvMmkwbB7jD+jhD+jiA9nSG6WwN0tvjpbPbT3thNa0MXWz7cRVdroN+5WXluyicXUjapgPLJBZRPKaRiSiGFlbnoyumOYgzwmy/9u2LbmqZhLXFWTSsO3fS7S/Ya7CORSDBjxoxp27dvt1955ZUtJ598cnBgnsbGRmt1dXV6cHJJSUlsMKEOh8PaSSed5L/jjjsaTjvttInf+c53yl577bUN7777ruPqq6+uGbdCLaVMCCFuBJ7FnJ51t5RytRDiC6njdwI/AO4RQnyIORn0FillRppTJvncbOoJcJKYy1+Sz/BhYhuzrNWAWfP6wvWfw6iczR3Pr+YHt36PFx9/grPO+RgLjplPUVFhv2s1NTbz4gsv88JzL/H+e6uw22186vJL+NhFFzNzYikA/1n5IQ6rzgVHlaNpgh9/YuZB235SbSEn1Zo2rN7VzTm/e50fnD+DTx9Tlc4zfcZU/v7Pu/mvG77Gl67/Krd8+ytc+qnR4dTlYCnxWrl0Zi6XztztNe2V+gAvb/Xz8lZ/vz7uE2u8nDTBFO6pBY5h7WvWdY2sXDdZuW5g76PsI6EYzfUdNNa3sWuzuezc2MLbz6zm+b/v/s5YbDrFVXmpmnc+RZVmX3l+abbZX57rwuHKvF91hSJTWCwW1q1bt6atrU0/++yzJ77zzjuOgX68B2vpGuxv32q1yosuuqgHYPr06WG73W7Y7Xa5YMGCcENDw2H/Q8vo6CEp5VPAUwPS7uyzvQs4faTtGowKtwubplEcqCLL5eTByPK0UPdy/TlHc+3H5vDvfzr585/v5TltOvbf/pspzjBut5uenh78PX66uszC14S6qXz++mu56OLz+fOKFj7511W8fWshbruFOz89D59z+PsvpxR5+e450zhvjlkgSCSNdE27qLiQv9x7J9/6xnf58Q9/wc4dDdz8tS+NqIOUw0mhx8rFM3K4OOV4pbeP+8UtPby02c8ja7sAKPZYOKnGy5LU4LQp+fYRGyTmcNmomlZM1bTiPY75O0M0bDSnnO3caM4Xb9zazqrXNhEN7ek03Oaw4M1x4/LacXrsODx2nG4bdpcNR2pxeh24vHZcXofZl57rJivPhS/Pgy/fndEIaIrxwb5qvoeb/Pz85HHHHed//PHHfQOFurS0NF5fX58W2sbGRltlZeUef0wWi0X2fgc1TcNut6dDWyaTycP+cRh9w3xHKbomqPG62OoPcWb2PB6IvkFzsosiPbtfPqvVyuVXXMoJZ5/DDfe+xayZPnrWvkOr4WJH9QxOcjYzuzIHWTWL7z+/g59ffAL5BR7OnWVhUoEn7fc4+zDVhqy6xtWLTd/b8aTBFf/3FkumFHD9iWYfrMvt4le//Sm/+Omv+dtf72dXQyM//Mn3cDodh8WeTDKwj3trZ5QXN/t5aUsPS7f6uf/DTgCKPBaOKXezoNzNwgo3c0tcZDtH/k/Hm+OibkE1dQuq+6X3On/paOymbVc37bu68XeF8HeECHSGCPmjhAMRwsEY/o4gkVCcaChGOBAlEowN2X+uaQJfgYfc4izyy8xm/sLyHAorTTeuxdV5uLzj73ehGPvs2rXLYrPZZH5+fjIQCIiXX34562tf+1rTwHznnXde1x/+8IfCa6+9tmPp0qVur9ebHNjsPRpQQn0ATMjysL47wMLkTB5gGY9E3+bzrsEr/OW5bh798smpvUv5qKGb257fwA1nncWkQg/b2oNEdDdeu/lfMKPMx4wy3wg9iUnSkJRluyjx9f/Y6rrOLd/+KuUV5fzip7/mc1ddz22/++moc44y3NTk2PncfDufm5+PlJKN7VFeqffzen2At3YGeXTd7m6oqmwbs4qczCx2MrXAwdQCJ3X5Dly2kW99EELgS400r5lRekDnGoZBJBgj1BMhkOpH7+kI0t0WoKPJT2dzDx1NPTRtbWfVq5sIB6L9zs8u8FBcnWqCr84zp6xVpqatlWSpGrkiI+zYscN61VVX1SSTSaSU4vzzz+/41Kc+1Q3QN8zlJZdc0v3kk0/6qqqqZjidTuOuu+6qz6jhQ6DCXB4ASUNyz4Z6JmZ5eNn7Eu8mtvDv7K/iEOOjL/DtrR1U5Dop8TnTaUtfepVvf+N/8Hjc/Oq3P82Yf/PRQEcowTsNQd5rDLOqKcQHTWE2tEf6RWeqyrYxtcBBbb6DaYUOphU4mV7oyEgN/HAQ6A7TVN+e6kdvp3FLG0317eza0k77ru5+tXMhBL58d3oqmtvnxJ3lwJXlwOm2Y3NasTms2B1W7C5z2+aw4vTa+zXHu7zDO2ZAMTjjaXrWWGS0Ts8ac+iaoMrjot4f5OP5C3klvprnoh9wnuPoTJt2yEQTSW7653vUFnu55+oF6fSTTl7Cvf/4Mzfd+HU+e+UX+Z//9y3OPvfMDFqaOXJdFj422cfHJu9u+YgmDDa1R1nTGmFta5j1rRHWtUV4pd5POL5btMqyrOka+OxiF0eVuJicZz/g0eaZxuNzMml2OZNml+9xLBaJ07rTjH7Wsr2T1oYuulr9dLX46WoJ0NbQTcgfIdgTHrRPfSgsVp2s1Lz03CIvucU+ckuyyCvJIq/ER16pj7wSH75897gZT6FQ9EUJ9QFSk+VmU0+QglgxU/RS/hNZxjn2eWhibH8g7Badu686mmzXngPYJk+ZxH3/+gtfv/nb3PrN7/P+e6v46jduwuFQ/ZN2i8b0IifTi5xATjrdMCTbu2OsbomwuiXM6uYwq5rDvLDFTzxpCrjbpjGn2MncUjcLy90sKHcxMXfkBq4NNzaH1YxeNmnfXSSGYRCLJIhF4kRDcWKReGo7RjgYI+SPEPZH8Xf2NsWnmuObe6hf20RXSwBjQKBp3aKRU5RFXnEW2UXelBc6L778XsczpvMZh9uG3WnF4bJhsVnQLRq6RUdLBclOx8/WtTH7f6EYXyihPkAqPS50IdgaCHFpzmJ+EHyA5fENLLLV7fvkUc7UEtNRi5SSHz65lnlVOen52Dk52dzx59/y+9/cwV//ch/vvbuKn/ziBwcUGORIQtME1Tl2qnPsnF27uwYeSxisa4vw7q4Q7zaGeHdXiLtWtvK75S0A5Dh1FpSZg9YWlrs5usxNvnv8/ZlqmpYeeU7ugZ+fTBp0tfhpb+yhfVdXyuFMD+2N3bQ3dtO0tZ11b9fT0x4acrDc/tkp0HQNi003m+mdVuxOW7oJ35VlNs97sl14c5x4fC48Oc60RzpPtguPz4nTO3YLYIrMM/6+AIcZq6ZR4XGytSfIpwqnUqT5+FfkjXEh1L1E4gYf7OjCqmv9HKdYrRZu/tqXWHjM0Xzn29/n05dezU1fuYFLPnkhuq4GDe0PNovGrGIXs4pdXJVKSyQla1rDvL0zxNsNQd7eEeSHLzdipPSlJsfG/FI3c0tdzCp2MqvYSanXekR/+HVdM5u9S3wwt2LIfMlEkkBXGH9niEBXmEBXmGgoRiS1JONJkgmDZNLASBpICUiJlBIjKTGS5rF4NFX7j8SJBuOE/BFCPRFad3YR6A4R6AyTiCeHtEPTtbS49/bVOz329OJw2bA5zYKAtbeWb9XRdQ2hCYQQe9T4hSbQdR3doqFZNGx2C1a7BavNgtVhwe6wmuMAUukWm7lWjnLGHkqoD4Iar5t6f4j2aIILHcfyh9AzrEs0UGcpy7Rpw4LTpvP3zy3EbjH/oLtCMXzO3cKw6Lhj+PeDf+e7t/6An/7vr3jisWf47+99k7qpUzJp9pjFoou0eH9ufj4AgWiSFbtCvNMQZMXOECt2BXlgdWf6nBynTm2+g0m5dibl2ZmQY6cy20aFz0aZ14rNoj7GALpFx5fvwZfvOaz3kVISDcXNAkFKuHcXDkIEusOEeiIEu8MEeyIEuyO07uwkHIgSDkSJhuPEwnEM4/AP7tU0gW7VsVj1VLO/pkbnj3KUUB8EVV4XGrC1J8jZhfP4a/hl/h15g+96Lsm0acOGw2r+4XaH43z89jf42PRivnXW1PTx/II8bv/jbTz95HP84me/5vJLr+aTl13MtV+4muzskZ1mNh7x2E1vaSfWeNNpneEEHzaHWdUU5sPmMJtS08f+/kFHv3OFgEK3hfIsG2VZVsqyTAGv8Fmp8NmoyrZR5rUdVKASxeAIIXC4bTjctv1yKzsYUkoSsSTxWIJkwiART2IkDLN2LyXSkOkaP5Cq7UuSiSSJeDJ9bjySIB5LEAubLQDxaMJsEUitE/Ekybh5Tm9rQjKe5N61w/c+FMOLEuqDwKHrlLqdbPEHWViYyzn2eTwQWUZTspNiPWffFxhDZDksnDenjBOm7DlASAjBWed8jOOOP5bf/voO/vH3f/HIQ49x+RWf5NNXfoqsLO8gV1QcLDlOC0uqvSyp7v9ew3GD7V0xtnfH2NFtrht64jT0xNjSEeXV+gBdkf7NsroGZV4bE3Jt1OSYNfJJeXam5DmYnGfHY1c1rJFGCGE2Xdsz81n+r99m5LaHhf0Nc3nHHXfk3nbbbcUAbrfb+MMf/rDt2GOPDQP84Ac/KLz33nsLpJR85jOfaf3ud7/bcjC2zJo1qy4Wi2nd3d16JBLRioqK4gCPPvroptra2ti+zgc1j/qgWd3Rw2tNbVwyoZykLcynum7j4/YF3Og+67DeN9M881ET86pyKPDa9zi2adMW7rz9Ll547iW8WV4uuvjjXHTpBZSVHZgTDsXwE4gmaeiJs707xrauGPVdUbZ1xtjaFWVzR5TmQKJf/lKvlboCh7nkO1JOXRyUHOF94+OZ8TSPen/DXD7//PPuOXPmRAoKCpL//ve/s374wx+Wrlq1at0777zjuOyyyya+++67ax0Oh3HCCSdM+eMf/7ht5syZ0aHuuS9++9vf5q1YscJ97733bh/suJpHfRioyXLxWhNs8QeZX5DDKbZZPBldyaedJ5CtuTNt3mGhIxjjaw98wLmzSwcNEjJp0gR+cdv/sm7tBv5859389S/3cc/df+f4JYu48OKPs2jxMVhtKv5yJvDYdWoLdGoLBp9SF4wl2dQeZWN7lA3tETa0mfPB//5+Oz3R3dOgfA6d2vzemrdZ+67JsVOVbaPIY1EirhgV7G+Yy9NOOy0t3CeddFLwxhtvtAF8+OGHzrlz5wa8Xq8BsHjxYv+//vWv7JkzZzb3Pf/CCy+sdjgcxqZNmxwNDQ32P/7xj1vvueee/JUrV7qPOuqo4IMPPlg/HM+jhPogcVksFDntbEsJ9eXO43k+9gEPRJZxreu0TJt3WMh12/j354+lMs+MWNc3oEdf6qZO4Ze/+QlNjc38598P89CDj/LqK2+QlZXFyaeewOlnnMrRR89Voj2KcNt0Zpe4mF3SPxqhlJJGf5y1rZH0srE9MmjfuMMiKPVaKfFaKfZYKfRYyXNZyHPp5Dkt5Los5Dh0cp0Wcl06OQ6L6ic/Avjsw/UVHzWHhzXM5YwiZ+j/Lqg+5DCXffnd736Xf9JJJ3UDzJkzJ/z//t//K2tqatLdbrd8/vnnfbNnzx70/O7ubsubb7654R//+Ef2pZdeOvmll15aN2/evPCsWbOmLlu2zLlo0aLwwT+piRLqQ6Da6+atlg4C8QSV1gJOss3g4chbXOJYjE8b1t/lqGFaqTnXOpE0uPqed5hflctNp04eNG9xSRE33vQFPv/Fz/Lmm2/z3DMv8PyzL/LIQ4/jdrs45tgFHH/icSw+7hgKCvJH8jEU+4kQgtIsG6VZNk6ZmNXvWChmsLkjyrbuKPWdMeq7YjT54zT646xpjfByvZ+OcJK99a5l2TXyXBby+y5uCwVuKwVuC4VuC8UeK0UeK0UeS3omgkKxL/YnzGUvjz/+uPfvf/97/rJly9YBzJ07N3LTTTc1nXzyyVNcLpcxbdq0kMUyuFyeffbZXZqmMXfu3FBeXl58wYIFYYApU6aEN2/ebFdCnWGqvC7eaulgmz/E9NwsPu04gZdiH/Jg5E2ucZ2SafMOK4aE8hwnRVl79lUPxGqzsuSExSw5YTGRSITly97mtVeX8dqrb/DiCy8DMHnKRI5dtJBjFy/kqLmzldezMYDLpjGz2HSLOhRJQ9IVSdIeStAZTtIZSdARStIRTtARTtAeStIRStAeTtAaSrCmNUJrMEEobgx6vRynTonHrLWXeq0U96vBWyhwmQKf49SVqI8S9lXzPdzsLcwlwFtvveW8/vrrq5588smNxcXF6VGXN998c9vNN9/cBnDjjTeWlZeXDzrwy+FwpENe2my2dLFU0zQSicSwNBkpoT4EcmxWsmwW6v1BpudmUWMpZIl1Gg9Fl3OJYxEebegP2FjHZtH48Sdmpfdf2dBKIJLg7FklezkLHA4HJ568hBNPXoKUkg3rN7HsjeUsX/YW99/3APfe8w/sdjtz583h2MULWHzcIiZMrFZ9n2MUXROp5u8D+9QEY0lagwlaggmaA3GaAnGa/InUOk5jIM6r2wI0+uPEkoNX2e0Wgc+u47XrOC0Cl03DadGw6QKbrmG3iNS2wGbRsOsCh0XDadVwWAROq5nfXKf2U2kOS/+8vdt2PeWYZBRjGJK4IYklJfFkaj0C87dHkv0Nc7lx40bbxRdfPPHuu+/eOmvWrH4DxRoaGixlZWWJjRs32p588snst99+e93IPUF/lFAfAkIIqj1uPursJpY0sOkaVzhP4NWeNTwYXc6VzpMybeKIce+yenZ1RzhjRvF+B5oQQlBbN5nauslc/dkrCIfCrFz5Hm++8TZvvrGcX/38d/zq57+jpLSY445fxPEnLGbBwnmqtn0E4LbpuG061Tl7b7GRUtIZTtIUiKeFvTUYpyuSpDuSpCuSJBgzCMcNwglzHYgZxJJJogmDWFISTUqiCYNoUhJJGP2CqRwMVl3gsAjsqcKAI1UIsKYWiybQBWi9HsYAidlKZc6ZNreTvduG7He81zop6bdtpM+VJCUkDJleekU5bkiSgzdWjCv2N8zld77znZKuri7Ll770pSoAi8UiP/roo7UA55133sSuri6LxWKRv/71r7cXFBQM7XruMKOmZx0iu4JhHtvWyGnlhUzMMr0f3er/B6sS9dzv+woe7cgQlXjSoC0QpcTnJBJP8sxHTZw/p/SQasKNu5p44/XlvP7aMt568x3C4TAOp4NFixdy4klLOOHE4/Ap5yqKYUZKU9hMcZeE4gaRPkIfTkgiccMU9YQp7pGEQSQuCScMon3SoklJLGEWBOIp0YwnTSHtFVcJmN5BTeHuFXGt39os2PaWgXv/qkTqPDDnxmvC3LPoAl0DXZgFA1tvQUHb3YrQW3Cw6RpWTXDt0QXjZnrWWERNzzqMFLscOHSNen8oLdRXOk/k8z13cn/ktXE7AnwgVl1Lx7F+6N0Gvv3wh9Tku5ldkX3Q1ywpLeaiSz7ORZd8nFgsxoq332XpS6/yysuv8dILr2Cx6MxfMI9TTzuJk089kdzc8eVsRpEZhBDYLQK7RSM708aMINdm2gDFkKga9TDwUkML2/whrqytQkvVIP838CAvx1Zzr++/KNazR8yW0YBhSN7c0s7iSeZI7qc/bKSuJIua/OGZXy6lZPVHa3nx+aW88NxSduzYia7rLFg4n4+deSonn3ICWb6sfV9IoVCkGU8OT8Yie6tRK6EeBrb0BHluZzPnVpVQ5jZrlS3Jbq7o/g1LbNO41XPRiNky0kgp6ZJBdibbaTa6CcoIIRklJKMIBDIp+PVvIkyqsnLTpSVkCze5modiLRubOPR51FJK1q/byHPPvsizTz9Pw85dWK1Wjj9hEWed/TGOP2Exdvu+R6YrFEc6Sqgzi2r6PsxUeJzoQrDNH0oLdaHu42LHIu6LvMqFiWPHTWStsIyxOrGd9+P1rErUsyXZQlDuOTVRQyABiST7M1baDMH3AzESfiudr1aQc1wDJTl2yrRcKvUCavQiJuiFTNCLDmi0vBCCuqlTqJs6hS/d9AU++nANzzz1PM889RwvvfAKHq+H0z92CueedyZz5s5Wo8cVCsWYQwn1MGDVNMrcDur9QY4tyk2LwWXO43kqupI7Qs/wa+81Y1YkYjLB8vgGXoh+wPL4BuIk0dCo1Us51TaLcj2PCi2PYj0Hr3DgEnbsmD6hk9IgmZMkKKN0GUGe2dnIbzfv4qJTKwhautge7uQF/UNCvJO+X4mWw2S9hMmWEmotZUzRS/fLgYwQgpmzpjNz1nRu/tqNvP3WSp564hmefvI5HvrPo5SVl3L2uWdyzrlnUFk1dAxjhUKhGE2opu9hojdIx6UTy8mx29Lpj0Xe4bbQ43zf80mW2KaNqE2HSqvRw7/Db/BM7D0CMkKu8HCybSZHWycxw1qJSxxck3IknkyH0fzKv95nU2uAu744g63JFjYnm9iQbGRjopFdxm4XlSVaDlMspdTqpen1/ta8Q8EQL734Co8/+iRvv7USKSWz58zk3PPP4rTTT1YjxxUKVNN3plFN3yNAldcM0rHNH+on1Gfb5/Jo9G1+F3ySuZaaMeEEZVeyg39EXuPZ6PsYSE60TecM+1HMtdSgi0MPf9gr0gAn1hUyvcxHoZ5NoZ7NmjesXDFlLhMLPASMMBuSjaxPNLA+sYsNiV28EludPrdEy2GSXswkSwkT9CJq9CJKtGw00d8jlcvt4pzzzuSc886kuamFp554hscfe5offv+n/ORHv+T4Exan+rMXqTnaCsU4IpFIMHPmzGnFxcWxpUuXbhp43DAMrrnmmoqXXnrJ53A4jLvvvrv+uOOOCx3ofR588MGsW2+9tRxg+/bt9sLCwrjD4TCmTp0aevjhh+sP9TmUUA8THquFPLuNbYEQc/Kz0+m60PmG+wKu7/kTt4ee4RbPBZkzch8EZYS/hV/hP5E30dA4yz6XTzmOO6wxts+bvTsEZktPhB8/vY5I3OCLJ3pw4WCmVs1c54R0nm4jxMbkLtYndrEp2cTmRCOvx9fR6wbCjpVKPZ8KPZ8KLY8KPZ9SPZdSLYcs4aKouJCrP/cZrvrsFaxbu56nnniWp598jqUvvoLL5eKEk47j9DNOZdHihWoQmkIxxvnhD39YNGnSpHAgEBi0hvHAAw/4tmzZ4qivr/9o6dKl7uuvv75y1apVB+yB7MILL+y58MIL1wAsWLCg9he/+MWOJUuWHLDgD4US6mGkyuvivbYuIskkDn3376LWUspljuP4e+RVToxNZ6FtSgat3BMpJc/HPuCPoefokAHOtM3lGtfJ5GsjO8WpMMvBsm+ejC3lo/mNzW187YEP+Os1C6grNm3xaS7ma5OYb52UPi8so9QnW9mabGFropltRitrEjtYanzEbj9O4BZ2irRsijQfhVo2RTU+Zt10Iif91znsWlnPiueW8/ILr/L0k8/hdDpZfNwxnHjyEo5bsohs1TyuUIwpNm/ebH322Wd93/rWtxpvu+22osHyPProo9mXX355u6ZpnHLKKcGenh7Ltm3brFVVVfG++Vwu11FXXnlly6uvvprl8/mSP/rRj3becsstFbt27bL99Kc/3X755Zd3H85nyahQCyHOAH4D6MBdUsqfDJLnRODXgBVok1KeMIImHhBVHhfvtnWxIxBmss/T79gVzhN5Pb6OXwQf4y+WG0eNx7KmZBe/CD7CysQW6vQyfui+jKmW8ozZk+/ZXYv12C3Mr8pNz79+c3M7Fl0wvyqn38A8p7Az1VJu2t2nEhyVcXYlO9hldNJodLAr2Umz0U2L0cXqxA56ZJ+gNlPNJfvGEyl6L0j41WaWvb6CF55fitAEE2dM4tjFCzl58fHMnDmdoSLpKBSK/jz5QbyizS+HNZxgvleEzp5t3WuwjxtuuKHiZz/72c7u7u4h++saGxut1dXV6WAbJSUlscGEOhwOayeddJL/jjvuaDjttNMmfuc73yl77bXXNrz77ruOq6++umbcCrUQQgduB04DdgLvCCEek1Ku6ZMnG/gDcIaUcrsQojAjxu4nhU47Dl1jmz+4h1DbhIVvui/g+p4/c3vo6Yw3gUspeSr6Ln8IPYNEcrPrXM6xz9ujfzeTHFWZw+2X7252/+2LG2kLRHnu5iWA+Qx7G0lvF1ZqLEXUMGhhmrCM0mr4aTW6aTG6aTV6zOXYbloXlpC8uRp9bQvxZa1sfbuZTXf+nb/d8Xc0l4WcOaVUzp3I1HlTmTltOmX2PIq0bPI177D04ysUioPn/vvv9+Xn5yeOP/740BNPPOEdKt9gg6kH+6ZYrVZ50UUX9QBMnz49bLfbDbvdLhcsWBBuaGiw7XHCMJPJasECYJOUcguAEOKfwPnAmj55LgMeklJuB5BStoy4lQeAEIIqj4ut/hCGlGkvZb3UWsq43HE8f4u8Ql2kjPMdCzJiZ4fh52fBR3grvpGjLDV8w/3xw9oPPVz831XzaegMI4QgnjQ47/dvcM3iai6ef3BTrZzCTqVup1IfOhZ2cHGElmO6aTK62Nqxi/feeo9NK9bR9O523lu2lPdYCjYNfUoWlhk+rFOzya8roaS8hALdR57mIU/zkqN5yBZusoWLbM2NRzjxCLsSdcW4Z18138PB66+/7nn++eezy8rKfNFoVAsGg9r5559f8+ijj27tm6+0tDReX1+fFtrGxkZbZWVlfOD1LBaL1DSzEqNpGna7PR3aMplMHvZ5t5kU6jKg73/gTmDhgDxTAKsQ4mXAC/xGSnnvYBcTQlwHXAdQWVk57MbuL1VeN+u7AzSFIpS69xzhfaXzJDYlG/lt6CnK9Nx+fa0jwbLYOn4WfISwjPFfrrM53370qKpF7w2XzcLkIrNw3BOOU5XrIs9j/o31ROKsbujhmAm5wzpf3S0c1Fgc1FDEscW1XHb+SWZxEmhrbeed997lrffe5cMPPmLbQ1uJxrYR4AMasuzYa3MQk9wkJ7mwTPaiVbgQA2IkO7HhFnYcwoZT2HAIKzYs2IQFG1YsQseKji40rOhY+m6njlnQsQoLdmHFTmotrDiw4hQ2nMKOW9hxCRtWLGN2Pr9Csb/cfvvtDbfffnsDwBNPPOH95S9/WTRQpAHOO++8rj/84Q+F1157bcfSpUvdXq83ObDZezRwQEIthDgGOAM4BigFnEAbsB54BXhEStm5v5cbJG1gO4QFmAeckrrXm0KI5VLKDXucKOWfgD+BOY96P20YdsrdTjRgWyA0qFDrQuM7nov5Us9dfC/wb27PupYqveCw2xWWMe4IPcPj0RVM0ou51XMR1fqo7knYK3keO3deMS+9/9DKnXzv8TU8fdPxTC0ZmUFw+QV5nHn6aZx5uhl4JR6Ls2nTZlZ/tI7VH61h/boNbHpwC7GY2QVmsVoorimlcHIJORML8UzIxVmTjVbsIEqCsIwSkXEixOkxwsRkgjhJkiSJyyQJkiSkYa5JYuzx57JvrOi4hQNP30Vz4hUOvMJJlnDi1VxkCSc+4SIrtZ0lXOhjpECnUOyNvmEuL7nkku4nn3zSV1VVNcPpdBp33XVXfYbNG5T9cngihLgS+BowHegBVgGtQBjIBWowa79R4N/A96WUe5ReBlzzWOB7UsqPpfa/BSCl/HGfPN8EHFLK76X2/w94Rkr5wN6unQmHJ315YlsjgXiCT04aukm2KdnF9T1/wims/C7rc+RqQ3ajHDIfxrfxk+DDNBqdXOJYxDXOU7CJ8TUYKhJP8uqGVk6fXgzAb17YSE8kznfOnprRGmQ8nmBb/TbWr9vAxg2b2bhxMxs3bKKluTWdx+VyMWFiNRMnTWDipBomTZrIxEkTKCwq2KvtyZRox1OCHpMJoiSIyThhGSOCuQ7LGCEZJSgjBPusAzJMwIjglxECMkyPDJNk8GDFAoFXOMkWLnI0Dzmah1zhITfVtJ+neckX5tornKrWPgZRDk8yyyE5PBFCfAAUAvcCnwHel4OouxDCB5wDXA6sFkJcLaX8114u/Q4wWQhRAzQAn8Tsk+7Lo8DvhRAWwIbZNH7bvmzONFUeF280t9Mdi+OzDR54oljP5ofeT/HVnnu4vufP/NjzaWosw1vDjckE94Rf4l+RNyjUfPzKexVzrDXDeo/RgsOqp0UaoDMUozscTwvGW1vamV2R3c/ZykhgtVqYNHkikyZP7Jfe0+Nn86YtfZatvP7qMh59+Il0Hm+Wl4kTa5g8ZSKTp0xi8uSJTK6djMdjjoLXhYaOhn0YgpuAObAmTIweI0y3DNEjQ/QYIbpliG4jRJcM0mUE6ZJBNiUa6ZQBgjK6x3VsWMjXssjXvP3WecKbFvYc4cEt7ErQFYr9YJ81aiHEl4E7pRwk8sLQ58wGiqWUz+4j31mYU6904G4p5Y+EEF8AkFLemcrzdeBqwMCcwvXrfd0/0zXqnlicf2zawaKiPGbl7X3+7frELr7t/ztREnzfcynzrBP3mn9/eTu2kd+FnmKn0c7Z9nlc7zrjoF1+jlV6R4W3+qMc++MXuXbJBG45oy7TZu2Vzs6utHhv2riFTRvNWnjAH0jnKa8oo7Z2shmMZFodU6fWkl+QlxF7IzJGhxGg3fDTZvhplz20GX5ajR7aUkur4SdOYo9zrej4hAuf5sYrHKkBdo50X71D2LBh6dcXL4RAQ6T+9WdgSm8eDS1dqNHRsPZeT+jYsGIXFmxYzPEBqX7+sTJuYzhRNerMMmxhLoUQ/wa+KqUc8VF8B0KmhRrgn5t24LFaOKeqZJ95m5NdfCvwd7Yn2/i863QusC/EcpCjgRuTndweepo34uuo0PK40XUWC2yTD+pa4wXDkCzb3E5VnouKXBcfNXTzrYc+5OcXz0o7UhnNSClpbmphw/qNrF+/kfXrNrJ+3QZ2bN+ZzlNQkM+06XVMnV7HtOlTmTatLmPiPRApJX4Zpt0I0CFNQe8yAnTJEF1GkB4ZIiAj6SUsY0RkfFBxHynSwo01JeAW7ClRtwsrVixYhY4NC1ZhwYqeKgBYsKD1GQS4u5Bh6VNA6M1vQccmLOm1FQu2VHrvPSxoaGgH3fpgSAOJJIlM/+uLhoaOwKJZlFBnkOH09X0R8Av6j9YGQAiRC9RJKZcdsIXjkCqviw87uoklDWz63kvnRXo2v/V+jh8GH+APoWd4PPIO17lOZ7G1br//ONcmdvKfyJu8HFuNDQvXOU/jIsexWMdZX/TBoGmC4ybvnoIViCbQNEFxlul05p36Dhq7I5w1oxjLPv6vMoEQguKSIopLilhy4nHpdL8/wIb1G1mzeh3r1q5n7ep1vPrKG+m5oYVFBWnRnjq9jmnTasnLH3nxFkKQlRqYVsP+d+8kZZIYCRLSIE6ShNwtNck9BtLJPfaMVF5DSgwMkqQG4slkuk8/TpKojBMjQVTG09uR1HaUOBEZT/X/m2ndRogYCeIykVqb1zPXiSH7+Q8FgUgLtoZAE1qqZWHA86af1VwbgwizYuyxP33Utal8a/eRdTLwGmYz9hFPpcfFB+3dNATD1GS595nfozn4sefTvBlfzx9Dz/HfgfuZppdznG0qc60TmKSX9Bt1G5MJ1ica+DCxjWXx9axO7MAt7FzoOIaLHYsoGGH3n2OJYybk8egNi9P7/3pnB69uaOWcmWbrx6YWP6XZTly20V3I8Xo9zJt/FPPmH5VOCwVDrFu3gTWr17Fm9VrWrF7HK0tfS4t3UXEhU6fWUjetlmnT6qidOoXCwr0PWssUutBxog8+P2QUI6U0B/mlCgUJDOIyQYLewX/JfoMAd69NoU8Lf6owkMBIXccc6Z+Ugwtwb5eAJkRa0PUBwq4N0mmQTAn6y/xgJF+T4gDYny/RJ4H/ASKYBbdvCyFeAt7FHFjW23nmS+VRAMUuBzZNY1sgtF9CDWbNY5GtjoXWyTwZfZeHosv5U/h5CINHOMgSLnprDm2Gn1iqabBGL+QG15mcZZ97xPVDDwc/u3AWOzpDaJr5Abvhvvco8Nr5++fMaf3+SByvY3gGbB1uXG4Xc+fNYe68Oem0QCCYrnGvWWPWvl95+fW0eOfm5VBXV0vd1CnU1k2mtm4KlVUV9Dp4UBwYQgiz2RrLmCpkXJlpAxRDsj9C/SvgVWAu8HNgBuZcahtgCCE2Y3oTm405bUsB6EJQ4XGyPRDap6vLPc/VOc9xNOc5jqbd8PNufAurEtsIyxgCsxnMJ1zMtlYz01KJT9u/goBicDRNUJVnvkMpJd87b3q6thJNJFn045e4bskEvnTKZKSURBPGiI8ePxQ8Hjfzj57L/KPnptNCwRDr121g3boNrFu7gbVr1nPvPfeRSCQBcDgdTJ48kSm1k5lSO4nJUyYxafJEsrIO3zRChWI4KSsrm+l2u5OapmGxWORHH320R6vwuAlzKaX0A0uBpUKIzwJXYArydEzxPgqYBrwHfPdQDRpPVHpcbO4J0haJUeA8uJpunublNPtsTrPPHmbrFIMhhODYibv7ceNJyRdOnMiCmlwAdnaGOeWXr/DrT87hrJklROJJonEDn2ts1Lh7cbldHDVvDkf1qXnHYjG2bK5n/boNrF+3kQ3rN/L8sy/x4AOPpPMUFRcyadIEJk6awISJNUyYWEPNhGq8Xs+eN1EoMswrr7yyoaSkZMhRieMyzKWUclqf3fdTi2IIKj1mwJhtgdBBC7Uis3jsFm44abebV10TXL24mrpis2a5fEs7V/3lHf7zhWOZX51LU3eE7R0hZlf4sFvGTq0bwGazmVO+pu4Ow9o74nzjxs1s2rCJjRs3s3nTVt55+920xzUwR51X11Sll6rqCqqrqygpLUbXx9Z7UBw5qDCXCpwWnSKnne3+EPMLRn/QC8W+Kc128q2zpqb3a/LdfP1jtdSl3JY+81Ej33t8DW99+xSKsnSWbWrjnfpOPn/ChDHVXN5L3xHnxy9ZlE5PJpPs3NnAls311G+pZ+vWbWzdUs8zTz+Pv8efzme1WimvKKWqqpKKynIqqyqoqCynorKc4uIiJeJHAM9+Y0NF24bg8Ia5nOIOfexnU/Y5TfiUU06ZLITg6quvbv3a1762x3SxcRPmUgjxKKarz/f254JCCAdwPRDqdVpyJFPpcfFOayehRAKXimE87qjKc/ercZ83p4zqfDeFXrMF5e36Du56bStfOtnM87sXN/Lapjb+dd0xCCHY1RXGbbOMuaZzXdepqqqkqqoSTl6STpdS0tnZxbat26iv3872bTvYvm0H9fXbeXPZ20Sjuz2ZWSw6JaUllJeXUVZeSllZCaVlJZSUllBcUkR+fp4a0KY4aN5444111dXV8YaGBsvJJ588Zfr06ZEzzzwz0DfPeApzuR1YLoR4H7gPeB1YJaVMt/sLIUoxw1aeC3wC0yXoNcNu7RikKiXU2wNh6rLVQJzxTq7bxom1u+cKf/nUKXzhhInpEeV5HjtVua70x+CHT65h9a4eXvn6SQA8sWoXdovOadMGj6E92hFCkJubQ25uTr/+bzAH7rQ0t7J92w527mxg544GduxooGFnA2vXrKOrq3+lxGLRKSgsoKiokMIic51fkE9+QR75+Xnk5eeRl5eLz5elBH0Usz8138NBdXV1HKCsrCxx9tlnd7355pvugUI9bsJcSim/JIT4NfBl4HuY07CkEKIHMwhHDmDFnIjwdirf36SUwz/rfwyS57Dhtuhs94eUUB+h9G3yvmxhJZct3B2G9apFNXQEd/f1/unVLeS4bGmh/sq/36cmz82XTjG9y7X4I+S77WnhH0tompZuRl/Ang6wgsEguxoaaWxspqmxiaamFpoam2hpbmXd2g28+vLrRCJ7+hbXdZ2cnGxycrPJyckhO8eHz+cjO9tHli8Lny+LrCwv3iyvufZ68WZ5cLlco3L+uOLQ6enp0ZLJJDk5OUZPT4+2dOnSrFtvvXXXwHzjKsyllHIz8CUhxFeBYzGDY5QCDqAdWAe8KqXcdrgMHasIIaj0uNjUEyApJbr6MCj60DuavJf/fGERPZHd34l4UpLs0zx35q9f4/Tpxfz4EzMBuPv1rcyvzmFWefaI2Hs4cbvdZvCRKYPHaJdSEggEaWtto621nfb2jvTS2dFBZ2c3nR2drFu7gZ7uHnp6/BjG0PUFTdNwe9x4vR48Hk9q7cbTZ9+b5cXjcePN8pr7KZH3ek3Bt1pVd9ZoZOfOnZYLLrhgEkAymRQXXnhhe2/T9bgNcznWGA2+vvtS7w/yzI5mzq0qoWyQGNUKxf5gGJJ/r9hBVZ6bYyfmEYwmmP4/z/L1j9Vyw0mTiMSTnPWb17jp1MmcP6eMeNJg+ZZ2ZpT6yHEf9m60UYdhGAT8AXp6/KnFFO+AP4DfH8Df4ycQCBAIBM20QICAP0ggEMDfEyAQCOxV6AGcTmdauLNSYu5JCbrH694t+Cnxd6cKAt5UutM1ekKCqqAcmWU4fX0rDoIytxNdCLb5Q0qoFQeNpgk+uWB3s7nbbuGj738MI1XYDkQT1JV48TnNgWnb2kNc8X9vc9uls7ngqHJ2dIS45cFVfPX0WuZV5dAdivP+zi7mlGePucFs+4OmaWT5ssjyHZw7XSkl4VDYFPeAKez+lPD7/X78PQF6enoIBILpYy0trWzZUp8qBARJJpN7vYeu63g8brKyskzBTzXN94p+li8rXcv3ejy43C7cbhculwuH04HDbsdmt2O32w64n94wDOLxOLFYnHh81LX2KvqghHoEsGoapS4H2wIhFjE6Ihopxgce++4/4XyPnT9cPi+9X5bt5P5rj2Fioel1LRBNEIkn6e3e/rChmyvvfpt/XncMx0zIY9mmNr72wAf86TPzmVHmY21jD4+838BnF9dQmOWgLRBlV1eY2mLvmJsjfjAIIXC5XbjcLuDAB/f1Cn3fmnrAb9bU/alafSAQoKfbLAj09Pjx9/hpbW5NtwL0nau+x/VNKxFIJALN6cGKgUVIhMWC4chCj0fQjRgJKYg7syHUhYyEiEtBIqsILdCOFgthWB0H+ZYUI4ES6hGiyuvi9aZ2uqIxsu1HXjOkYuRx2vR+XtamlmTx0PW7g5HMrvDxwBeOTTtvyXJaOXZiPnke8/e5pTXIX16v57JULf6ldS184z+reO0bJ1GR6+Khd3fy6xc28uAXF1HgtfPm5nZeXt/Cl0+dgtOmU98WpKErzDET8tA1QSiWwJDgtumjprn3QOnrDrgtEEUTgtxUt8JHDd247RZq8s2C0dL1LRR6HcwoK6SoCB5cuZPqyjIWV5njEv706mbmlfpYNCkfKSXffXQ1F0wp4NRpRcQSBp/96zucN6uIk2u8NLd385l/beDjk10cVwTN3UFuecvgnJIo8zxB2kNxft9QxMmeNqZZO+lMaNznn8iJjl3UWrrpwc79gQmc7mlmuitMp3RwT2sxFxb1MCdH0pyw8/VrX8nMS1XsEyXUI4TppaydbYGQEmrFqMDrsHJ09e7BbDPKfPzykt2uas+eVcJZM4vT+8dPzudPV8yjMMucI17gtXNUZTZeh/kZWdPYwz3L6vnK6aZns0feb+DXL2xk04/OBAS3L93EHS9vZvP/ngXAr57fwCPvNfDqN8ypaX94eROvbmjln9cdC5hC9u62Lu68Yl56f12jn19dOiedv74tyM8ump2+XlN3OL3/wyfW0BWO84uLzf1vP/whkVgyff7N/3ofIeBXl5j7X/jbStx2S/odXH7Xcgq9Dm5L5f/Yba8yqdDD7ZebPtM/8YdlzK3M5tefNKOXXXvvCo6blM/PU/f7xn9Wcfr0Yv73AnPg3/cfX80n5pYzLyXUtz2/kSuOrWLRpHyEEDyzuolin4NTpxVh1QXBaAKETl5+Hlk5OZw41c9xs0o4sbaQaCJJo3czJ0wp4KjKHMKxJIUrd7CwJo/aYi+hWIK5HzYxr+pEqvPdBKMJTtjYyqzybEqznQSjCU7b0UVtsZd8j51wLMnXrz2QX49iJFFCPUJk2azk2K1s94eYnZedaXMUiv2ib823xOekxLd7jMXxkws4fnJBev+zx9VwzeLq9DmfPLqSxZPy0zG+T6otJM9tTx+vK/ZyytTdc859Tms6RvhghGMGwViiz36SYHR3H7BhSJJ9xn65bDrxPgkFHjvRxO79qjxXv5CPdSXeflPpFk3MJ8ux+xN52cLKdO0Z4Osfq023PgDcdumcfsf/ce0xZPWJuvbczSfgsu++/nvfPQ1rn/jn79x6anpbCNGv9cOqa+kCB4DdovPlU3e7enXadD5zbHWfZ7dw4bzy9L7bbuGMGSX99hdPyu93vmL0csijvoUQQo6yoeOjbdR3L8ub21nV3s1VtdXYdOWgQaFQjB7UqO/MsrdR38OhFt1CiFP7JgghVPFsEKo8LgxgR3DYgqooFAqFYhDa2tr0M844Y0JNTc30CRMmTH/hhRf2iAdsGAZXXXVVRWVl5YwpU6ZMe/311w/KJ/mDDz6YVVdXN62urm6ay+U6qrq6ekZdXd20Cy64oPqQH4ThafpuB34uhPiilHJ5Km05cPQwXHtcUeRyYNc0tvtDTMxSYQEVCoXicHHddddVnH766T3PPPPMlkgkIgKBwB4V07ES5nI4atSdwMeBPwshZqbSVN/3IGhCUOFxsi0QSs99VSgUCsXw0tHRob311lveL3/5y20ADodD5ufn7zGpfagwlwPzuVyuo774xS+WTZ8+feqiRYumLF261LVgwYLa8vLymffdd5/vcD/PsAiqlHKbEOIS4N9CiAvoneKn2IMqr5tNPUFawlGKXWruokKhGN+0xh+tiMmWYQ1zaROFoQLr+UMG+1i3bp09Nzc3cfHFF1evWbPGNWvWrOCf//znHVlZWf1czY2VMJfDUaP+EEBKuRa4CngIlFePoaj0OBHANr/qp1YoFIrDQSKREGvXrnXdcMMNrWvXrl3jcrmM//7v/y4emO9gw1wed9xx/tEW5jKNECIfCEgpI71pUsor+2yvFEJ8CXhg+EwcX9h1nRKXg22BIAuLcvd9gkKhUIxh9lbzPVxUV1fHioqKYieffHIQ4NJLL+38yU9+sodQj5Uwl/usUQshdCHE94QQXUAz0COEeFAIkT1YfinlK1LKwsGOKUyqvC46onF6Ysq/rkKhUAw3lZWVieLi4tgHH3xgB3juueeyamtrIwPznXfeeV333XdfnmEYvPjii2M6zOUXgO8CLwPvABOAC4Ae4OrDZtk4ptrr5s3mDrb5Q8zMO+zjEBQKheKI43e/+932yy+/fEIsFhOVlZXR+++/vx7GaZhLIcT7wFtSys/3Sfs88HvALaUc2mt8hhitDk/68s9NO3BbLZxbVbLvzAqFQnGYUQ5PMsuhOjyZwJ59zv8CdKDqUAwTQpwhhFgvhNgkhPjmXvIdLYRICiEuOpT7jSaqvC4ag2Giyb3Hu1UoFArFkc3+CLUHs5m7L/7U2nuwN055L7sdOBOYBnxKCDFtiHw/BZ492HuNRqq9btNLWUCN/lYoFArF0OzvqO8yIcSEPvt6n/SuvhmllFv285oLgE29+YUQ/wTOB9YMyPcl4EHGmaezIqcdh66xLRBikk95KVMoFArF4OyvUP9niPRHBknbXz/fZUDfYfs7gYV9MwghyjAHrp3MPoRaCHEdcB1AZWXlfpqQOTQhqPS40l7KtDEan1ehUCgUh5f9EerDNbJ7MGUaOLLt18AtUsrkvgLNSyn/BPwJzMFkw2Hg4abK62JDd4CmUIRSt3PfJygUCoXiiGOfQi2l/OthuvdOoKLPfjmwa0Ce+cA/UyKdD5wlhEhIKR85TDaNKBVuF5qAen9ICbVCoVAoBiWTQZHfASYLIWqEEDbgk8BjfTNIKWuklNVSymrM5vfrx4tIA9h0jXK3k63+4KCu7BQKhUJxcFx88cXVubm5sydPnjx9b/meeOIJb11d3bRJkyZNP/roo2v7HkskEkydOnXaSSedNOlg7Zg1a1ZdXV3dtJKSkpk5OTmze8Nhrl+/fr9dj2YsypWUMiGEuBFzNLcO3C2lXC2E+ELq+J2Zsm0kqfa62R5ooyMaI89hz7Q5CoVCMS645ppr2m666aaWq6++umaoPG1tbfpNN91U+cwzz2ycPHlyrKGhoZ8m/vCHPyyaNGlSOBAI7O/Yqz3oDZv529/+Nm/FihXue++9d/uBXiOTNWqklE9JKadIKSdKKX+USrtzMJGWUl4lpRxqUNuYpdprBpXZqoJ0KBQKxbBx5plnBgoKChJ7y3PXXXflnn322Z2TJ0+OAZSVlaXzb9682frss8/6rr322iEdt1x44YXVl19+eeXChQunlJeXz3zyySc9F198cfWECROmX3jhhdXD9SwZFWoFuCwWip12tvqDmTZFoVAoDgvn//712r+8sTUPIJYwxPm/f732b29uywUIRhPa+b9/vfafb2/PAegMxfTzf/967X9W7swGaOmJWM7//eu1j77f4APY1RUetpbgDRs2ODo7Oy0LFiyonT59+tTf//736ciPN9xwQ8XPfvaznb3BOIaiu7vb8uabb274yU9+suPSSy+d/PWvf71548aNq9etW+dctmzZsAw+UkI9CqjJctMeiakgHQqFQjGCJBIJsWrVKtcLL7yw8YUXXtj485//vGTVqlX2+++/35efn584/vjj99nUefbZZ3dpmsbcuXNDeXl58QULFoR1XWfKlCnhzZs3D0t/Zsb6qBW76Q3SUe8PMUsF6VAoFOOMR288bn3vts2iyb77brvF6Luf47Il++4XZjkSffdLs517bc4+EMrLy2P5+fmJrKwsIysry1i4cKF/xYoVrpUrV7qef/757LKyMl80GtWCwaB2/vnn1zz66KNbB17D4XCkQ17abLb0qGBN00gkEsPiIEPVqEcBPpuVXLtNNX8rFArFCHLRRRd1vfnmm554PI7f79fee+89z8yZM8O33357Q3Nz86qGhoYP77nnni3HHHOMfzCRHilUjXqUUON18W5bF+FEEqfloAcYKg4j/lic7YEw7dEY3bE43bE44UQSKSUS04OPw6LjtljwWHV8Nit5Djt5DhvZNqvyPqdQjCDnnntuzfLly72dnZ2WoqKiWd/85jd33XzzzW19w1zOnTs3cuqpp3bX1dVN1zSNK664ovXoo4/eI251ptlnmMuxyFgIczmQ1nCUB7c2cGJJPnU5WZk2R5HCH0+wrrOH+kCI9ogZ0dWmaWTbrPjsVlwWHQ2BECAlhJIJgvEkwUSC7lgcI/XnZRGCIpeDktRS5LRj2ccgFYViJFFhLjPL3sJcqhr1KCHfYcNjtbDFH1JCPQrojMZ4v62Ljd0BJFDkdHBMYS7VXhc+m5V9ubQFSEpJVzROeyRKSzhKYzjCitZOAHQhKHY5KHM7KHc7yXfYVY1boVAMihLqUYIQggleNx91dhNNJrHrqvk7E4QSCd5s7mBjdwCLEEzLzWJ2XjZe64H/qehCkOewkeewMSXbjAgbTSZpCkVpCIZpCIZ5u6WTt+nErmuUuZyUe5yUuZ1kWS37VRhQKBTjHyXUo4iJWW5WdXRT7w9Rm33Qob4VB4GUkvVdft5s7iAuDY7Ky2ZWnm/YxwvYdZ0qr4uqlKObUCJBQzDCzkCIncEwW1IDCj1WC2VuB2UuJyUuB16bdVjtUCgOM4ZhGELTtPHXt3oYMAxDAMZQx5VQjyIKnXY8Fp0tPUEl1CNIMJ7gpYYWGkIRip0OTijNJ8e+3254DwmXxcJkn4fJPg9SSrpi8XRtu74nxPquAGAKd7HTQbHLTqHTQZ7Dhq5q3IrRy0etra3TCgoKupVY7x3DMERra6sP+GioPEqoRxFCCCZkeVTz9wiyKxjm+Z0txA2DJSX5TM32ZqzJWQhBjt1Gjt3GjFwfUko6ojEaQxF2hSLsCoXZ1GMKty4EeXazWT0/NbI8x25VvxnFqCCRSHyuqanprqamphmoacD7wgA+SiQSnxsqgxLqUYZq/h4ZpJSs6uhmeXMHWTYr51aXkDtCtej9RQiRmt5lTwt3IJGkORShJRylLRJlc0+QtV3+9DlOXSfbbiXLZsFrteK1WnBZdJwWHZdFx67r+1UTl1JiSIgbBnHDIGbI1NpIpUkShkHCkCSlxEjll0gEAoT5dbZoGroQWDSBTdPMRdewp5f9s0cxtpg3b14LcF6m7RgvKKEeZajm78NPUkpe2dXKhu4ANV4XJ5UWYtNHf6FfCIHXasHr8zDJ5wFMQfXHE3REY3RG43RFY3TF4uwIhAklAoNep1dALZpAYM7/RgikNEU3KSVJQw7dYTYEGqmLSZDmar+waQKHrpuLRRtk2xR0u65h18xtqybGxGA7Kc332DsNVhOpdz4GbFeMHpRQjzL6N38b2MeAgIwl4obB8ztb2B4IMS8/m/kFOWP6oymEIMtmJctmpXpAuS5hGATiScKJBOFkknDCIJJMkjAkcWnWhmVKVSWmiGhCoAuwCA2rJrDqGlbN3LZp5rYtJfLW1FoTAo09xadX+BNSmvc0DGJJs1YeTZq2RJMG4dQ6kkgSSiTpiMaJJJIk9uLjQYBpi65hS9nSa2evXRZhrnUh0rX63mfUUs878L9e9hYy+hZaUvabz2GknsXcjqfS0i0Mqby95+3dfoElZbdd13DoZsGkt/XDbbHgtui4rRbcVotqeTiCUUI9Cult/t7mD6an9SgOnWgyyVPbm2gORzm+OJ/pueN7vrpF08i2a2TbMzNiXAiBRQjzI3MQXecJwyCS7C/qsWSSaEroYynh79s0H0r0iqYkkS6MDNPzsFtcLZrAKnYXWJwWc91bIOgtHPTWoDUhkEikJN1N0FtYiqeep7eQEk4kBxV5t0XHY7XgtVpSayue1LbHasE2RloZFAeOEupRSG/z9+YeJdTDRTiR5PFtjXTFYpxWXsjELE+mTVLsA4um4dE0PIdQzuhtek72qeX2CqXR25TA7mZ6kepbF6K/2PaK70ggpSRmGAQTSYLxBIE+iz+eoCUcZUtPcI+uCYsQZu07VSN39i6proTeMQFmK4SmhH0MoYR6FKJGfw8vkWSSJ7Y30h2Lc1ZFMeUeV6ZNUowQQgh0QNfHjiAJIVJ98vqQAxwNKQmlhNyfEvFgIpFOawlHCSeTxI29tydY+3QbKEYvSqhHKZN8ZvP3lp4gU5VL0YMmmkzy5LZGOqMxzlQirRgnaEKkm7yL9pIvbph9/1HDIJIwiBrJ9DgBs9tApkf2K0YvSqhHKQUOOz6blY3dASXUB0ksafDk9ibaIzE+VlFMhRJpxRGGVdOw2jRUB9rYRrV3jFKEEEzxedgViuCPD1uc9COGpJQ8u7OZ1nCU08qL0i47FQqFYqyhhHoUMzk1V3ZT9+DzYRWDI6Xk5V2tNATDnFBaQE2WO9MmKRQKxUGjhHoUk2WzUuS0s7Hbv+/MijTLW8zoVwsKcqhTo+YVCsUYRwn1KGeKz0tHKqaxYt982NHNB+3dTM/J4qj87Eybo1AoFIeMEupRzsQsNxqwQTV/75Pt/hDLmtqp9rpYXJyn5ogqFIpxgRLqUY7DolPhcbGpO4CxF5eERzodkRjPNzST67BxSlkhmhJphUIxTlBCPQaY7PMQTCTZFYxk2pRRSTiR5OkdTVg1jTMripXzBoVCMa5QX7QxQLXXhU3TWK8Gle1BUkqe29lMKJHkYxVFeKzKNYBCoRhfZFSohRBnCCHWCyE2CSG+Ocjxy4UQq1LLMiHE7EzYmWksmsZkn4ctPUGiyWSmzRlVLGtqpzEU4cTSAoqcjkybo1AoFMNOxoRaCKEDtwNnAtOATwkhpg3IthU4QUo5C/gB8KeRtXL0MDXbS1JKNqpBZWnWdflZ3dnD7Fxfes65QqFQjDcyWaNeAGySUm6RUsaAfwLn980gpVwmpexM7S4HykfYxlFDvtNOvsPG2k5/Ogj9kUxLOMprjW2UuRwsLMrNtDkKhUJx2MikUJcBO/rs70ylDcVngaeHOiiEuE4IsUIIsaK1tXWYTBxdTM3Ooj0aozUSy7QpGSWcSPLczmacFp1Ty4vUCG+FQjGuyaRQD/Z1HbSqKIQ4CVOobxnqYlLKP0kp50sp5xcUFAyTiaOLST4PFiFY19WTaVMyhiElLza0EE4k+Vh5EU6LCgGqUCjGN5kU6p1ARZ/9cmDXwExCiFnAXcD5Usr2EbJtVGLXNSZkudnYHThiw9K929bFzmCYxcV5FDjtmTZHoVAoDjuZFOp3gMlCiBohhA34JPBY3wxCiErgIeAKKeWGDNg46piak0XckGzuCWbalBFnRyDEitZOpvg8TFU+vBUKxRFCxoRaSpkAbgSeBdYC/5ZSrhZCfEEI8YVUtu8CecAfhBDvCyFWZMjcUUOx0062zcraziOr+TsQT/BiQws5divHl+Qr96AKheKIIaPeIaSUTwFPDUi7s8/254DPjbRdoxkhBNNysljW3E5LOErhEdD8m5SS53c2k5SSj5UXKc9jCoXiiEJ98cYgddlerJrgw47uTJsyIrzT0kFzOMoJJQVk222ZNkehUChGFCXUYxCbrlGb7WVzd4BgPJFpcw4r2/wh3m/vZlqOl0nKqYlCoTgCUUI9RpmZ68MA1ozjvupAPMHSXS3k2W0sKsrLtDkKhUKREZRQj1F8NitVHherO3tIjMOpWkkpeSHVL31aeSEW1S+tUCiOUNTXbwwzM9dHJGmwaRxO1XqnpYOmcJQlql9aoVAc4SihHsOUuR3k2q182NE9rvx/b0/1S0/N9qpgGwqF4ohHCfUYRgjBjFwf7ZEYu0KRTJszLATiCV5K9UsvLlb90gqFQqGEeowzxefBZdFZ2dq578yjHENKXtjZQsJQ/dIKhULRi/oSjnEsmsacvGx2hSI0BsOZNueQeKelk6ZwhCWl+apfWqFQKFIooR4HTM3x4tR1VrZ1ZdqUg2abP8R77V1MzfYyxaf8eCsUCkUvSqjHAVZNY3aej53BME1jsK/aH4vzUkMLeQ7VL61QKBQDUUI9Tpiem4VD11jZNrb6qpOG5LmdLUgkp5cXqX5phUKhGID6Ko4TzFp1NjsCYVrC0Uybs9+82dxOayTKiaWF+GzWTJujUCgUow4l1OOI6TlZ2HWNt1s6xsS86g1dfj7q7GFWro8JWe5Mm6NQKBSjEiXU4wibrjEvP4edwTDbA6N7BHhrOMorjW2UuhwsLMrNtDkKhUIxalFCPc6YnptFts3KsuZ2kqO0Vh1OJHl2ZzNOXeO08iJ0ITJtkkKhUIxalFCPM3QhOLYoj+5YnNUdoy+yliElz+9sJpxIcnpFMU6LnmmTFAqFYlSjhHocUulxUuF2sqK1k3AimWlz0kgpeb2pnV2hCMeX5FPotGfaJIVCoRj1KKEehwghOLY4j7hh8M4oci26qqObNZ09zMnzUZetnJooFArF/qCEepySa7cxPTeLtZ09NI4CJyhbe4K82dzBBK+bhYVq8JhCoVDsL0qoxzELCnLxWC0sbWghljQyZkdLOMKLDS0UOu2cXFaAUIPHFAqFYr9RQj2OsekaJ5cV4o8nWNbcnhEbOiIxntrehNOic2ZFsfI8plAoFAeI+mqOc0pcDubkZbOuy0+9Pzii9+6Kxnh8WyOaEJxTVaJGeCsUCsVBoIT6CGB+YQ55Dhsv72ollEiMyD17YnEe39aIRHJuVYlyD6pQKBQHiRLqIwBdCE4pLSRhSJ7e3kzcOLz91V1RU6QT0hTpHBVbWqFQKA4aJdRHCLkOG6eVF9IWifLcjubD5rWsORThkfoG4obBOZUl5DnUXGmFQqE4FJRQH0FUed0sKclnRzDMy7tahz1wx1Z/kMe2NWLTNC6oKaNAOTRRKBSKQ8aSaQMUI8vUnCxCiSTvtHZi0zQWF+ehHeJ0KUNK3mvrYkVrJwUOO2dWKtegCoVCMVxktEYthDhDCLFeCLFJCPHNQY4LIcRvU8dXCSHmZsLO8cbc/Gxm5/lY3dnD49saCcYPfoBZVzTOo/W7eKe1k4lZbs6tVqO7FQqFYjjJWI1aCKEDtwOnATuBd4QQj0kp1/TJdiYwObUsBO5IrRWHgEgF7shz2Hh1Vxv/2dLAaeWFlLqd+32NuGGwprOHd1o60YXg1LJCJvk8h9FqhUKhODLJZNP3AmCTlHILgBDin8D5QF+hPh+4V5qdqcuFENlCiBIpZePImzv+mOLzkm+38+zOZh7b1kilx8nMXB/lbueQ3sPCiSQfdXSzurOHSNKgwu3khNICPFbVi6JQKBSHg0x+XcuAHX32d7JnbXmwPGXAHkIthLgOuA6gsrJyWA0dz+Q6bFxYU5YOmPHk9iaybVZKXA7cVgtui44hoS0SpS0Soz0axZBQ7XUxOy+bEpcj04+gUCgU45pMCvVgVbaBw5D3J4+ZKOWfgD8BzJ8///DMPRqn2HSN+QU5HJWfzZaeIGs6e9jqDxLp4x/cpmkUOGzMzDUjX6m50QqFQjEyZFKodwIVffbLgV0HkUcxTOhCMNnnYXKqrzkpJaF4Agl4rRYVTEOhUCgyQCZHfb8DTBZC1AghbMAngccG5HkM+Exq9PcxQLfqnx45dCHw2qxk2axKpBUKhSJDZKxGLaVMCCFuBJ4FdOBuKeVqIcQXUsfvBJ4CzgI2ASHg6kzZq1AoFApFJsjoUF0p5VOYYtw37c4+2xK4YaTtUigUCoVitKBciCoUCoVCMYpRQq1QKBQKxShGCbVCoVAoFKMYJdQKhUKhUIxilFArFAqFQjGKUUKtUCgUCsUoRgm1QqFQKBSjGCXUCoVCoVCMYpRQKxQKhUIxilFCrVAoFArFKEYJtUKhUCgUoxgl1AqFQqFQjGKEGfdifCGE8APrM23HPsgH2jJtxH6g7BxelJ3Di7Jz+KiVUnozbYRiTzIaPeswsl5KOT/TRuwNIcSK0W4jKDuHG2Xn8KLsHD6EECsybYNicFTTt0KhUCgUoxgl1AqFQqFQjGLGq1D/KdMG7AdjwUZQdg43ys7hRdk5fIwFG49IxuVgMoVCoVAoxgvjtUatUCgUCsW4QAm1QqFQKBSjmDEr1EKIM4QQ64UQm4QQ3xzkuBBC/DZ1fJUQYm4GbKwQQiwVQqwVQqwWQtw0SJ4ThRDdQoj3U8t3R9rOlB31QogPUzbsMU1jlLzP2j7v6X0hRI8Q4ssD8mTkfQoh7hZCtAghPuqTliuEeF4IsTG1zhni3L3+lkfAzp8LIdal/l8fFkJkD3HuXn8jI2Dn94QQDX3+b88a4twReZ9D2PivPvbVCyHeH+LckXyXg36HRuPvUzEEUsoxtwA6sBmYANiAD4BpA/KcBTwNCOAY4K0M2FkCzE1te4ENg9h5IvDEKHin9UD+Xo5n/H0O8htoAqpGw/sElgBzgY/6pP0M+GZq+5vAT4d4jr3+lkfAztMBS2r7p4PZuT+/kRGw83vA1/bjdzEi73MwGwcc/yXw3VHwLgf9Do3G36daBl/Gao16AbBJSrlFShkD/gmcPyDP+cC90mQ5kC2EKBlJI6WUjVLKd1PbfmAtUDaSNgwjGX+fAzgF2Cyl3JZBG9JIKV8FOgYknw/8NbX9V+Djg5y6P7/lw2qnlPI5KWUitbscKD9c999fhnif+8OIvc+92SiEEMAlwP2H494Hwl6+Q6Pu96kYnLEq1GXAjj77O9lTAPcnz4ghhKgGjgLeGuTwsUKID4QQTwshpo+sZWkk8JwQYqUQ4rpBjo+q9wl8kqE/gqPhfQIUSSkbwfxYAoWD5Blt7/UazJaTwdjXb2QkuDHVRH/3EE21o+V9Hg80Syk3DnE8I+9ywHdoLP4+j0jGqlCLQdIGzjPbnzwjghDCAzwIfFlK2TPg8LuYzbezgd8Bj4yweb0sllLOBc4EbhBCLBlwfDS9TxtwHvDAIIdHy/vcX0bTe70VSAD3DZFlX7+Rw80dwERgDtCI2bQ8kNHyPj/F3mvTI/4u9/EdGvK0QdLUnN4RZqwK9U6gos9+ObDrIPIcdoQQVsw/jvuklA8NPC6l7JFSBlLbTwFWIUT+CJuJlHJXat0CPIzZ5NWXUfE+U5wJvCulbB54YLS8zxTNvd0DqXXLIHlGxXsVQlwJnANcLqUc9EO8H7+Rw4qUsllKmZRSGsCfh7h/xt+nEMICfAL411B5RvpdDvEdGjO/zyOdsSrU7wCThRA1qdrVJ4HHBuR5DPhMarTyMUB3bzPPSJHqp/o/YK2U8ldD5ClO5UMIsQDz/6R95KwEIYRbCOHt3cYcXPTRgGwZf599GLK2MhreZx8eA65MbV8JPDpInv35LR9WhBBnALcA50kpQ0Pk2Z/fyGFlwJiIC4a4f8bfJ3AqsE5KuXOwgyP9LvfyHRoTv08FY3PUd6rAfxbm6MXNwK2ptC8AX0htC+D21PEPgfkZsPE4zGaiVcD7qeWsAXbeCKzGHE25HFiUATsnpO7/QcqWUfk+U3a4MIXX1yct4+8Ts+DQCMQxayGfBfKAF4GNqXVuKm8p8NTefssjbOcmzH7I3t/onQPtHOo3MsJ2/i3121uFKRYlmXyfg9mYSr+n9/fYJ28m3+VQ36FR9/tUy+CLciGqUCgUCsUoZqw2fSsUCoVCcUSghFqhUCgUilGMEmqFQqFQKEYxSqgVCoVCoRjFKKFWKBQKhWIUo4RaoVAoFIpRjBJqhUKhUChGMUqoFYq9IIT4uBDiKwd4zu+EEI8fJntuTgWlUH+7CsURgnJ4olDsBSHEPcCpUsr9Cv0ohJiIGUZwkZRyxWGwxwlsBb4lpfzLcF9foVCMPlSpXKEYXr4MfHA4RBpAShkG7gW+djiur1AoRh9KqBWKIUjVpq8EyoQQMrXU7yW/Hfg08I8B6ZoQwi+E+O6A9JzUNa/sk7ZQCPGYEKJRCBERQmwTQvx1wK3+CUwTQiw6tCdUKBRjAUumDVAoRjE/AAqAozHjXwNE95L/GCAbeG1A+hTAA7w3IP2o1Po9ACHE0alz7wM+B4SBycCMAee9D/QAZwDL9udBFArF2EUJtUIxBFLKzUKIViAmpVy+H6ccw+4oRX2Zm1q/OyD9KEzhX5va/zRQL6W8uk+elwaxyxBCrErdT6FQjHNU07dCMXyUAj1SytiA9HlAi5SyYUD6XGC1lDKe2m8BJgkhfiGEmLmPe7Wm7qdQKMY5SqgViuHDweBN43PZszYNZo26b3P4LzCb2z8BrBJCbBRC3DjEvcKA8xBsVSgUYwQl1ArF8NEO5PRNEEIIYA4D+qeFEIVAbd90KWVUSvk/UsoJwHTgA+B3QwwaywXahtV6hUIxKlFCrVDsnSj7X3NdB1iFEH3nXE/EHGCWHJD3S5h/f+8PdiEp5Rrg16ndwcaS1ADr99MuhUIxhlGDyRSKvfP/27l/V46iMI7j72c0KLEw+ReQ7AqLjDIpZOEPkDJZTBZm2/cPMFJMyiIsSgaDwiajZHkM5yopFnc48n4t90e3c7rTp+e595xroDciVoBz4DUzr7559qQ5jgEPzfnHj2TLEXFP+Q49RVn2BTAaEZeUtncXcAw8UgJ+AzgDTj9PEhE9lD/Jt3/1ZpL+BCtq6Wd7lHXLW5TQ/HZr0My8a56Z+XR7BHgG1oFNoAN0A7OUJVZzzSYmN5RW+C5wCKxRlmlNZObXanwaeAP2f/dqkv4CtxCVWhQRC8AOMJCZLxFxBJCZky3OcQA8ZeZ8W2NKqpcVtdSuDqV1vdpcDwMXbQ0eEUPAOKU6l/QPGNRSi5o29RLwEhGDQB8tBjXQDyxm5m2LY0qqmK1vSZIqZkUtSVLFDGpJkipmUEuSVDGDWpKkihnUkiRVzKCWJKliBrUkSRV7B4pfOIk7qTV1AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEPCAYAAACp/QjLAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAu2UlEQVR4nO3deXhU9dn/8fedjYQ17JAAgkoRNxYjglrr+uDjirRWrQtVW6wrttaqtba2z69qi7XaarVU3BW1FqltrbjjUkV2UAERFCRQCEtYE7Ldvz/OREMWModM5mT5vK5rrpn5zplz7njJ3Oe7m7sjIiJSVUrUAYiISNOj5CAiIjUoOYiISA1KDiIiUoOSg4iI1KDkICIiNaRFHUCidOvWzfv37x91GCIizcqcOXM2uHv36uUtJjn079+f2bNnRx2GiEizYmYraytXs5KIiNSg5CAiIjUoOYiISA1KDiIiUoOSg4iI1KDkICIiNSg5iIhIDZEmBzPra2ZvmNliM/vIzCbEyruY2Stmtiz23DnKOEVEWpuoaw5lwHXuPhgYCVxpZgcCNwKvuftA4LXYexERSZJIk4O7r3X3ubHX24DFQC5wJvBo7LBHgTGRBCgi0kpFXXP4kpn1B4YBM4Ge7r4WggQC9IgwNBGRVqdJJAczaw/8DbjW3beG+N54M5ttZrMLCgoaL0ARkVYm8oX3zCydIDE86e5TY8XrzKy3u681s97A+tq+6+6TgEkAeXl5npSARUSSYNq8fCZOX8qawiJysrO4fvQgxgzLTdr1ox6tZMBkYLG731XloxeAcbHX44C/Jzs2EZGoTJuXz01TF5FfWIQD+YVF3DR1EdPm5ScthqhrDkcBFwKLzGx+rOynwB3As2Z2KbAKODua8EREapeoO3t3Z2tRGRt27GLTjhI27Sjh1hc+oqi0fLfjikrLmTh9adJqD5EmB3d/B7A6Pj4hmbGIiMSr8s6+8ge88s4eYMywXIpLy9mwfRcbt5d89bwjeN64fRcbd5SwIfZ6044SyiriaxVfU1jUaH9TdVHXHEREGkUi7+x3lJSzOXZXv2lnCb/8R+139tf9dQE3P7+IHSXltZ6rbUYqXdpl0K19G3KzMzk0txNd22fQtX0burbLoEvscemjs1i3dVeN7+dkZ4WOf28pOYhIi1Pbnf2Nf1vIuq3FjNqvK9t3lbFjVznbikvZUlTK1qKy4Dn2fsvOUjbvLGHzzlK2FJVQWh7fnX15hXPuiH50aZdB11gS6NYh+OHv2j6Dthnx/eTe9L+Dd4sfICs9letHDwr/H2MvKTmISLNRVl7Big07+GzDDjZuL2HzzhI2bi9h046gqWZrUSlFpeUsL9hBebWmmuKyCm7/95I6z92+TRodM9PomJVOp6x09uvens7t0slum0HntsFz13YZdG6XweWPz2Hdtpp39rnZWdxy2oEN/jsrazhRjlZSchCRJqm0vIJ5qwr5aM0WFq/dyuK121i6bhslZRW7HdcuI5Uu7TPo0q4NndpmkJORyifrttd53gcvyqNdmzTat0mjQ2YanbLS6ZCZRlpq/IM3bzql8e/sxwzLTWoyqE7JQUQiUVufwJlDc5izcjN/n7+Gfy1ay6YdJQB0aZfBgb07Mm7UPgzu3ZGBPTrQrUMGndtmkJmeWuPcR93xOvm1dN7mZmdx4oE9Gxx7U7izb2zm3jLmjuXl5fns2bOjDkNE4lC9TwAgLcXomJXOph0lZKancOLgnpx2aG+G9etMjw5tCKZF7f35s9JTuX3sIS3qBzwRzGyOu+dVL1fNQUSSbuL0pTVG+5RVODt2lXHXt4fwPwf1on2bvf95ag139o0t1H99MxsJnEywvHYOkAVsAJYCM4Bp7r450UGKSMtSW5MPQElZBWOH90nINaJus2/u4uqBMbNxZrYI+A9wLdAWWEawgupm4AjgQSDfzB4xswGNE66INGfuzouL1pJaRxNRMsfxy57VW3MwswUES2Y/BlwEzPdaOirMrBNwGnA+8JGZXezuzyQ4XhFpptZtLeaHz8znP8s30rtTJhu3l1BS/tXIo2SP45c9i6dZ6WHgAXcv3tNB7r4FeBJ40syGAL0SEJ+ItACzPt/EFU/OZceuMv7vzIM4b0Q//rlwrfoEmjCNVhKRRuPuPP7+Sn71j4/p0zmLP1+Yx6BeHaIOS6rQaCURSari0nJufv5D/jZ3NScc0IO7zhlKp6z0qMOSOIUdrfQwsNXdJ1QpGxk7z7u19UWISOuzevNOfvDEHD7M38qEEwYy4YSBpKTEP09Bohe25nAc8OvKN2Z2JfCH2NuZZnaSu+9IVHAi0vwsXruVCyfPZFdpBZPH5XHC4IbPSJbkC7sTXE+CIayVJgBTgWMIRjT9OEFxiUgztOCLQs6d9D7pqSk8f+VRSgzNWNjksJlg4htmdjCwH3CLu78L3Aack9jwRKS5mPX5Js5/cCYds9J49rJR7N+jfdQhSQOETQ6zgG/GXo8DPnP3yjVwPwX2CXMyM3vIzNab2YdVyoaa2ftmNt/MZpvZiJAxikiSvb2sgAsnz6Rnxzb89bIj6dulbdQhSQOFTQ63ARea2RLgh8DjVT7LBcL2NzxCsBxHVb8FfunuQ4Gfx96LSBP16sfruPSR2fTv2o5nLhtFr06ZUYckCRCqQ9rdZ5rZ0QSzoNcDE6t8fCLBGkthzveWmfWvXgx0jL3uBKwJc04RSZ5/LVzLhKfncVBORx69ZATZbTOiDkkSJPQ8B3efRdC8VF0RkIjlMq4FppvZnQQ1myPrOtDMxgPjAfr165eAS4tIvF5fso4JT89jWL9sHvru4XTI1ByGliRUs5KZXWhmB5tZzd014Mfufm8CYroc+KG79yVouppc14HuPsnd89w9r3v37gm4tIjE4/0VG7n8ibkcmNNRiaGFCltzeJSg2WeXmX0EzIs95gN5ZvZtd/96A2MaRzBEFuCvBKu9ikgTsWj1Fr736Gz6dmnLIxePUGJoocImhy7AsNhjOHAUcAlQOfVxawJiWgN8A3gTOJ7d51WISIQ+Xb+dcQ9/QKesdB6/dARd2qmPoaUK2yFdCLwRewBgZt2BawiagM4Lcz4zmwIcC3Qzs9XAL4DvA/eYWRpQTKxPQUSitXrzTi6cPJMUM5743hH07qS9F1qyBi+85+4FwC1m1oUgObwU4rt1JZPDGhqXiCROwbZdXPDgTHbsKuPp8aMY0K1d1CFJI0vkqqwvE2wIJCItwLR5+UycvpT8wiLSUw3DmDL+CA7M6Vj/l6XZCzta6WdmdoqZ9a7l4x6EnwQnIk3QtHn53DR10Zd7PZeWBwsuf7Gp9r2fpeUJW3P4CdAecDNbD8wlGK1UDnwPuCGx4YlIFCZOX0pRafluZSXlFUycvlS7tbUSYZNDJ+BrBCOVDos9XwFkxz6/y8wuAOYAs919aoLiFJEkWlNYew2hrnJpecKOVnKCJTKWAlMqy81sP4JEUZk0vk9Qi6htspyINHGdstIpLCqtUZ6TrRFKrUVCOqTdfTmwnGDSGgBmpvUsRJqhGZ8UsLW4lBSDiip7O2alp3L96EHRBSZJVW+HtJn93cyGxXtCM8s0sx8BpzQoMhFJug/zt3DFE3MY1Ksjt511CLnZWRiQm53F7WMPUX9DKxJPzWEV8L6ZzQeeBN4BFrp7WeUBZpYDjABOB8YC+QQzp0WkmVi1cSfffXgWnbLSeeTiw+nZMZNzR6gBoLWqNzm4+9VmdjfBaqm3EnRKu5ltBXYBnYF0giU0Pogd97i7VzRKxCKScBu372Lcwx9QWl7B0+NH0bOj9mRo7eLqc4j1KVxtZtcBo4AjgBwgE9gILAHecveVjRWoiDSOHbvKuOSRWawpLOKp7x/B/j06RB2SNAFhRyuVADNiDxFp5krLK7jyqbksyt/Cny/M47B9ukQdkjQRiVw+Q0SaEXfnxr8t4s2lBdx21iGcdGDPqEOSJiTsHtIi0kLc+fJS/jZ3NdeeOJDvHKGOZ9mdkoNIK/TnGcu5743lnDeiHxNOGBh1ONIEqVlJpJW5/83l/OalJZx2aG/+78yDMLP6vyStjpKDSCty3xufMnH6Us4YksNd3x5CWqoaD6R2kSYHM3sIOA1Y7+4HVym/GrgKKAP+5e4/iShEkVoVlZQz87ONzF1VSEaq0TErnY6Z6XTMSqNTVjo9OmQyZ+VmJk5fyprCInKys7h+9KBIZxjf+/oy7nz5E8YMzeHOs5UYZM+irjk8AtxLlU2CzOw44EzgUHffZWY9IopN5EvlFc5Ha7bw9rINvLNsA3NWbqakvAIzcK/9O1XXJsovLOKmqYsAIkkQ97y6jN+/+gljh+Uy8ewhpKaoKUn2LFRyMLORwMnASIJJcFnABoJVWmcA09x9c7znc/e3zKx/teLLgTvcfVfsmPVhYhRJpNLyCp6Z9QV/eG0Z67ftAuCAXh0Yd+Q+HD2wOyP6dyElBbYVl7GlqJStRaVsKSrlminz2Fpcttu5ikrLk74fgrtz96vLuOe1ZYwdnsvEbykxSHziSg5mNg74MXAQsBVYCCwDioAuBDOmLwTuM7NngV+6+2d7GdPXgK+b2a+BYuDH7j5rL88lslfcnZc/XsdvXlrCioIdHN6/Mz89ZTBH7d+N7h3a1Di+TftUurX/qnxbtcRQKb+wiB27ymjXpvEr7UUl5dz8/CKmzsvnW4f14TffPFSJQeJW7/+hZraAYAvQx4CLgPmxfR2qH9eJoP/gfOAjM7vY3Z/Zy5g6E9RODgeeNbN967jmeGA8QL9+GqctiTFn5SZue3EJc1ZuZv8e7fnLRXmcOLhHqFE9OdlZX26xWd3ou9/i9rGH8PWB3RMVcg0rCrZz+RNz+WT9Nq49cSDXHD+QFCUGCSGeHqmHgQHufoO7z6vtRxrA3be4+5PufgrB+kuFexnTamCqBz4AKoBudVxzkrvnuXte9+6N9w9NWoctRaVcPWUe37z/PVZt2sntYw/hpQlf56QDe4Ye7nn96EFkpe++11VWeipXH78/GakpXDj5A37y3AK27Ky5oU5DvbhoLWfc+y7rtxXzyMUjuPbErykxSGjxrMp6d+XrWJPRde7+RT3fWQAs2MuYpgHHA2+a2deADIJ+DZFG82H+Fq54ci5rCouYcMJALvvGvrTN2Pumn8p+hdpGK1153P7c89oyJr21gjeWFvD/xhzM6IN6NfhvKC2v4I5/L2HyO58xtG82fzp/uHZuk71mdVQEaj/YrAIYGbujr/5ZF+AAd/9PiPNNAY4lqBmsA34BPA48BAwFSgj6HF6v71x5eXk+e/bseC8tAgR9C0/P+oJfvPARXdtlcO93hiVt8bkP87dw/XMLWbx2Kycd2JOrjtufIX2zmTYvP9QQ2IoKZ8ayAu5+dRkLvijku0f256enDCYjTUNVpX5mNsfd82qU15cczGwQQQ1jMcG8g7qSwxHAf9w9kn2jlRwkrJ0lZfxs2odMnZvP1wd24+5zhtK1fc3O5sZUWl7BpLdW8MCM5WwrLmNgj/as3LiTkvKvtkPJSk+tdRe24tJyps7NZ/I7K1hesIOeHdvws1MP5PQhOUn9G6R5qys5xFNvPpfgjr4YcOCnZvY6MJegc3p77LhOsWNEmrzlBdu5okqH7dXHD4xkJE96agpXHrc/F43ah2dmfcHtLy6hvNoNW1FpOXf8ewmHD+hCeblTVFrOPxeu4Yn3V7J5ZykH53bk7nOGcsohvVVbkISJp+bQAcgDhgMTgRVAH4K+gApgOfAxMARY5+6jGjPguqjmIPGau2ozFz88i9QU4+5zhnLM15rOYIb+N/4rruPM4MTBPfne0QMYMaCL1keSvbbXNQd33wa8AbxhZpcSzGdYSDDnYTgwDDgQmAf8PJFBiyTa28sKuOzxOXTv0IYnLj2Cvl3aRh3SbnLrGAKbnZXOT08ZTEqKkZZiDO2bTf9u7SKIUFqLsDvBHVjl7fzYQ6RZeHHRWiY8PY/9urfnsUtH0KND09sn+frRg7hp6iKKSsu/LMtKT+XWMw6KdF0maX2iXltJJCmmfLCKm59fxPB+nZn83cPplJUedUi12tMQWJFkimeG9N+BW919XjwnNLNM4Apgp7s/0MD4RBqscv+CYwd15/7zDyMrI5IBdXEbMyxXyUAiF0/NYRXwvpnNB54E3gEWuvuXi8eYWQ4wAjgdGAvkA5ckPFqRENyd3738Cfe+8SmnD8nhd2cP0WgekTjF0yF9tZndDVwL3EowZNXNbCuwi2AdpHTAgA9ixz3u7hW1nE4kaf74+qfc+8annHt4X3591iFadE4khLj6HNx9OXC1mV1HsG7SEQRLdmcCG4ElwFvuvrKxAhUJ44EZy7nrlU/45vA+3HbWIVpbSCSksKOVSgj2bZjROOGINNzD737GHf9ewulDcvjttw5VYhDZC2qAlRblqZmr+OU/Pmb0QT2569va2EZkbyk5SIvx3JzV3DxtEccN6s4fzxtOuvZIFtlr+tcjLcI/F67hJ88t4Kj9unH/BYdpVJJIA+lfkDR7Mz4p4IfPzOewfToz6aLDyExv2vMYRJqDBicH04pfEqE5Kzfzg8fnsH+PDjw47vAGbdAjIl9JRM1hi5mdWLXAzHTrJo1u6X+3cckjs+jZsQ2PXTKiyS6JIdIcJSI5bAQmmtnIKmXvJ+C8InVatXEnF06eSWZ6Co9fegTdOyR3kx6Rli4RyWEzMAb4i5kdEiuLq25vZg+Z2Xoz+7CWz35sZm5m3RIQo7Qg67cVc+FDM9lVVsFjlzS9ZbdFWoKEdEjHZkZ/G3jKzPYn2DEuHo8AJ1cvNLO+wEkE6zqJfGnLzlIumvwBBdt28fDFhzOoV4eoQxJpkRKRHBYBuPti4LvAVKBrPF9097eATbV89HvgJ8SfZKQV2FlSxiWPzmJ5wXYeuOAwhvfrHHVIIi1Wg4d2uPu4Kq/nmNnVwF/39nxmdgaQ7+4L6hsIZWbjgfEA/fr129tLShM2bV4+E6cvJb+wiDZpKZSUV3Dfd4Y3qa09RVqiUDUHM+sW26+hTu4+w9177E0wZtYWuJk4txt190nunufued2768eipZk2L5+bpi76ctvMXWUVpKUYJWVa8FeksdWbHMws1cxuNbNCYB2w1cz+ZmbZjRDPfsAAYIGZfQ70AeaaWa9GuJY0cROnL91tu0yA0nJn4vSlEUUk0nrE06z0A4I7+TeBWcC+wFnAVuDiRAbj7ouAL2sdsQSR5+4bEnkdaR4qawzVramjXEQSJ55mpe8Df3H34939Bnc/G7gSuMDMMhpycTObArwHDDKz1WZ2aUPOJy1Lx8za711ysrOSHIlI6xNPzWFf4MfVyp4B7gf2AZbt7cXd/bx6Pu+/t+eW5m3SW8vZWlxGqhnl/tWgtaz0VK4fPSjCyERah3hqDu0JmpCq2hZ71iBzSbj73viU215cwqmH9ua33zqU3OwsDMjNzuL2sYcwZlhu1CGKtHjxDmXNNbN9q7xPrVJeWPVAd1+RiMCk9XF37n51Gfe8towxQ3O48+whpKWm8M3D+kQdmkirE29yeK6O8mm1lGnRPQnNPRiF9Kc3l/Otw/rwm28eql3cRCIUT3JI6IgkkercndteXMxf3v6M80b05ddjDtG+zyIRqzc5uPujyQhEWqfyCudX//iIR99byUWj9uHW0w9SYhBpArQzikRm3dZirpkyj5mfbeJ7Rw/g5lMHo72jRJoGJQeJxFuxrT13lpTzu7OHqNNZpIlRcpCkKiuv4O5Xl3Hfm58ysEd7njl/OPv30IhokaZGyUEaReVqqmsKi8jJzuL60YMYuW9Xrpkyjw8+38Q5eX259YyDyMrQ4DaRpkjJQRKucjXVykXz8guLuO7ZBaSkQHpqCr8/ZwhnDVMzkkhTFndyiK2j9Azw+9gmPdJCFO4sYcHqLSz8opDlBdvJykilQ2Y6Hdqk0SEzjQ6Z6XRul07vTlnkZGfRKSt9j+erbTXVcnfapKbyr2u+zoBu7RrzzxGRBIg7Obh7iZmdCNzTiPFIEhSXlvPX2V8w6/PNLFxdyOcbd375WW52FiXlFWwtKmVXHfsmdGiTRk52FjnZmXTKSqe4tILisnKKS8spKq2oczXVopJyJQaRZiJss9K7wEiC5bulGaje9n/eiL68sGANn6zbTu9OmQzpk823D+/L0D7ZHNynEx0zv6oVlJRVsK24lG3FZWzaWcLawmLWFBaRH3usKSxiecEOMtNTyExPJTMtlY6ZaWSmpVBcS2LRaqoizUfY5HAdMM3MthMsnbGWavs8u7u26Woiamv7v/PlT+iYmcYjFx/OsYP2vGFfRloKXdu3oWv7NvSnHcS5E2v164JWUxVpbsImh0Wx53uovXnJ9+Kc0khqa/sHaJuRVm9iaIjKVVOrj1bSaqoizUfYH/JfUa2mIE1XXTumrdta3OjXHjMsV8lApBkLlRzc/dZEXtzMHgJOA9a7+8GxsonA6UAJsBy42N0LE3nd1iK7bTqbd5bWKFfbv4jUJ57NfhrTI8DJ1cpeAQ5290OBT4Cbkh1US7Dkv1vZsauc6mvYqe1fROIROjmY2TAzm2pmG8yszMyGx8pvM7PqP/R7FJsvsala2cvuXhZ7+z6g2VIhbS0u5QePz6FT23R+cfpB2klNREIL1axkZkcDrwIrgKeAq6p8XAH8AHgpYdHBJQQT7yROFRXOdc8uYPXmIqaMH8nh/bsw7sj+UYclIs1M2JrDHcB04CDgR9U+mwsMT0RQAGZ2M1AGPLmHY8ab2Wwzm11QUJCoSzdr989Yzisfr+Onpwzm8P5dog5HRJqpsMlhOHC/uzs1Ry1tALonIigzG0fQUX1+7Fq1cvdJ7p7n7nnduyfk0s3a28sK+N3LSzl9SA4XH9U/6nBEpBkLO5S1GGhbx2e9gS0NCwdi/RY3AN9w9531HS+B/MIirpkyj/17tOeOsYdo0xwRaZCwNYd3gGvNrOo6y5V39pcCr4c5mZlNAd4DBpnZajO7FLgX6AC8YmbzzeyBkDG2OhUVzlVPzaW03HnggsNo10bzEEWkYcL+itxCsL7SAuA5gsQwzszuAg4DDg9zMnc/r5biySFjavWem7uaeasKuevbQ9i3e/uowxGRFiBUzcHdFwDHAOuAmwHjqxFL33D3pYkNT+qzrbiU3760lOH9sjlLQ1RFJEFCtz+4+1zgBDPLBLoAheobiM59byxnw/ZdTB6Xp34GEUmYvW6cdvdiYE0CY5GQVm7cwUPvfMY3h/dhSN/sqMMRkRak3uRgZmE6md3dT2hAPBLCr/+1mPRU44aTtRyGiCRWPDWHFHaf0zAI6AV8TtD30BPoT7C3g/ockuTdTzfw8sfruH70IHp0zIw6HBFpYepNDu5+bOVrMxtDsI/DSHf/oEr5EQTLXGgL0SQoK6/gV//4mL5dsrj06AFRhyMiLVDYeQ7/B9xSNTEAuPtM4Fbg/yUoLtmDKR+sYum6bdx8ymAy01Pr/4KISEhhk8NAoK5FjNYD+zcsHKlP4c4S7nrlE0bt25XRB/WKOhwRaaHCJofPgMvq+Owygn4IaUT3vLaMLUWl/Pz0AzV0VUQaTdihrL8EnjSzDwlmSFd2SH8LOAA4P7HhSVVrtxTx5Pur+HZeXwb37hh1OCLSgoXdJvRpM9tAkCRuAtKBUmAWMNrdX0t8iFLpgTeXU+HOlcep9U5EGtfezJB+FXjVzFKAbsAGd69IeGSym3Vbi5ky6wu+dVgf+napa2FcEZHEaMgM6QqCTmhJgvvfXE5FhWoNIpIcoZNDbCOe84B+QPXZV+7u+yUiMPnK+q3FTPlgFWOH56rWICJJEXYP6VsI+hs+BOYDuxohJqnmgRkrKKtwrjpuYNShiEgrEbbmcClwj7v/sDGCkZrWbyvmyZkrOWtYLv26qtYgIskRdp5DV+AfjRGI1G7Sl7UG9TWISPKETQ4zgCGJuriZPWRm62PzJirLupjZK2a2LPbcOVHXa24Ktu3iiZkrOXNoDv27tYs6HBFpRcImh2uBi83sIjPrZmYp1R8hz/cIcHK1shuB19x9IPBa7H2r9Je3V1BSVsHVx6uvQUSSK+yP+SfAwcDDBLOjS6s9SsKczN3fAjZVKz4TeDT2+lFgTMgYW4QN23fx+HsrGTM0lwGqNYhIkoXtkP4Vu+/t0Bh6uvtaAHdfa2Y9Gvl6TdKDb3/GrrJyrjpefQ0iknxhl8+4tZHi2CtmNh4YD9CvX7+Io0mcLUWlPPH+Sk49NId9u7ePOhwRaYXCNislwzoz6w0Qe65zFra7T3L3PHfP6969e9ICbGyPv/c523eVcfk3NJ9QRKIRdhJcfftJJ2IP6ReAccAdsee/N/B8zUpRSTkPvfs5xw3qzoE5WnlVRKIRtuaQAli1RzfgKOBrsfdxM7MpwHvAIDNbbWaXEiSFk8xsGXBS7H2r8ezsL9i0o4TLj1Vfg4hEJ2yfw7G1lZvZfsA04LaQ5zuvjo8aWvtolkrLK5j01gry9unMiAFdog5HRFqxhPQ5uPtygjv8iYk4X2v1wvw15BcWccVx6msQkWglskO6gKBpSfZCRYVz/4zlHNCrA8cNapWjd0WkCUlIcjCzLsCPgOWJOF9r9OridXy6fjuXH7uf9oYWkciFHa30GTUnwWUQ7CMN8M1EBNWaTJuXz29fWsKaLcWkphhl5Y09x1BEpH5hZ0jPoGZyKAZWAn+N9T1InKbNy+emqYsoKi0HoLzC+dm0D0lNMcYMy404OhFpzcKOVvpuI8XRKk2cvvTLxFCpqLScidOXKjmISKSa4gzpVmNNYVGochGRZFFyiFBOdlaochGRZFFyiNAFR9RcLDArPZXrRw+KIBoRka8oOURo8X+30SYthV4dMzEgNzuL28ceov4GEYlc2NFKkiArCrbzz4Vr+P4x+3LT/w6OOhwRkd2o5hCR+95YTkZaCt87et+oQxERqWGvag5mNgQYBGRW/8zdH2toUC3dF5t2Mm1+PheN2ofuHdpEHY6ISA1hZ0hnA/8CRlYWxZ6rToxTcqjHn95cTqoZlx2jBfZEpGkK26x0G9AVOIYgMZwFHA88CawARiQ0uhZo7ZYinpvzBWfn9aFXpxoVLxGRJiFschhNkCDej71f7e5vuvtFwKvAhEQG1xL9ecYK3OHyY1VrEJGmK2xy6A2scPdygjWVOlT5bCpwaqICM7MfmtlHZvahmU0xs2Z/m71+WzFTPljF2OG59OncNupwRETqFDY5/BfIjr1eCYyq8lnC9rU0s1zgGiDP3Q8GUoFzE3X+qDz49meUlldwhbYAFZEmLuxopXcIEsI/gceBX5hZf6AMGAe8kODYssysFGgLrEnguZNu044Snnh/JWcMyaF/t3ZRhyMiskdhk8MvgZzY64kEndPnEPx4vwBcnYig3D3fzO4EVgFFwMvu/nIizh2V+974lOLScq48TrUGEWn6QjUruftyd3879rrU3a9z9z7u3sXdv+PuGxMRlJl1Bs4EBhAko3ZmdkEtx403s9lmNrugoCARl24Uqzbu5LH3Pufsw/oysGeH+r8gIhKxpjpD+kTgM3cvcPdSgs7uI6sf5O6T3D3P3fO6d++e9CDj9dvpS0hLSeFH/6MttkWkeQidHMxsmJlNNbMNZlZmZsNj5beZ2ckJimsVMNLM2lqwofIJwOIEnTup5q3azD8XruX7Xx9Az47NfsCViLQSoZKDmR0NvAccADxV7fsVwA8SEZS7zwSeA+YCi2LXmZSIcyeTu3Pbi4vp1j6D8d/QvAYRaT7C1hzuAKYDBwE/qvbZXGB4IoICcPdfuPsB7n6wu1/o7rsSde5kefnjdcz6fDM/POlrtG+jBXBFpPkI+4s1HBjr7m5mXu2zDUDTbfhPstLyCn7z7yXs36M95+T1jTocEZFQwtYcigmGrdamN7ClYeG0HE9/sIoVG3Zw48kHkJbaVPv9RURqF/ZX6x3gWjNLrVJWWYO4FHg9IVE1c9uKS7n71WUcMaALJwzuEXU4IiKhhW1WugV4F1hA0GHswDgzuws4DDg8seE1Tw/MWM7GHSU8fOpggsFWIiLNS9hJcAsIluteB9xMsGz3VbGPv+HuSxMbXvOzauNOHnz7M84cmsOhfbKjDkdEZK+EHkLj7nOBE2KrpHYBCt19Z8Ija4bKyiu49pl5ZKSl8JOTD4g6HBGRvRZ3zcHMMszseTM7BsDdi919jRLDV/7w+qfMXVXIr886hNzsrKjDERHZa3EnB3cvIVjWQkNvajHr803c+/oyxg7P5YwhOfV/QUSkCQv7Q/8uX+0fLTFbikq59un59Oncll+deXDU4YiINFjYPofrgGlmth2YBqzlq6GsALh7RWJCax7cnVumfch/txbz3A9GaSa0iLQIYWsOi4D9gHsIdoIrAUqrPVqV5+fl88KCNVx7wkCG9escdTgiIgkR9jb3V1SrKbRmqzbu5Od//4gR/btwhTbxEZEWJFRycPdb6/rMzI4FLmpYOM3H1uJSrn56Hmbw+3OHkpqiyW4i0nI0qIHczPYnSAgXAv0ItvS8JAFxNWmfrt/O+Mdms2rTTu79znANWxWRFid0cjCzTgT7Rl8EjIoVLyBYzntK4kJrml5fso4JU+aTkZbCU98fyYgBXaIOSUQk4eJKDmaWApxMkBDOADKBNcB9wJXAte7+VmMFGaVp8/KZOH0p+YVFdMxMY1txGQflduTPF+apxiAiLVa9ycHM7gTOB3oQLNn9PPAo8CrQka/WVmpxps3L56apiygqLQdga3EZqWZcNLK/EoOItGjxDGX9EUFieBHo5+7nu/vLsfkMjTZyycyyzew5M1tiZovNbFT930qsidOXfpkYKpW7c89ry5IdiohIUsWTHB4CtgGnAkvN7F4zG9G4YQHBXIqX3P0AYAiwOAnX3M2awqJQ5SIiLUW9ycHdvwf0Ai4A5gA/AN4zs8XADTRC7cHMOhIsDT45FkOJuxcm+jr1yamj6aiuchGRliKuGdKxFVifcvfRQF/gp0A5cCPBng53mNkFsWW8E2FfoAB42MzmmdmDZtYuQeeO2/WjB5GVnrpbWVZ6KtePHpTsUEREkir0Cqvuvtbdf+PuBwNHAH8CBgKPEay1lAhpwHDgfncfBuwgSES7MbPxZjbbzGYXFBQk6NJfGTMsl9vHBstvG5CbncXtYw9hzLDchF9LRKQpMfeGtwqZWTpwOnCRu49JwPl6Ae+7e//Y+68DN7r7qXV9Jy8vz2fPnt3QS4uItCpmNsfd86qXJ2RvBncvdfepiUgMsfP9F/jCzCrbb04APk7EuUVEpH5NeX3pq4EnzSwDWAFcHHE8IiKtRpNNDu4+H6hR1RERkcanLT9FRKQGJQcREalByUFERGpQchARkRqUHEREpAYlBxERqUHJQUREalByEBGRGpQcRESkBiUHERGpQclBRERqUHIQEZEalBxERKQGJQcREalByUFERGpQchARkRqUHEREpIYmnRzMLNXM5pnZP6OORUSkNWnSyQGYACyOOggRkdamySYHM+sDnAo8GHUsIiKtTVrUAezB3cBPgA51HWBm44HxsbfbzWxpI8bTDdjQiOdvbIo/Os05dlD8UWvs+PeprbBJJgczOw1Y7+5zzOzYuo5z90nApCTFNNvd85Jxrcag+KPTnGMHxR+1qOJvqs1KRwFnmNnnwNPA8Wb2RLQhiYi0Hk0yObj7Te7ex937A+cCr7v7BRGHJSLSajTJ5NBEJaX5qhEp/ug059hB8UctkvjN3aO4roiINGGqOYiISA1KDiIiUoOSQz3M7GQzW2pmn5rZjVHHE5aZPWRm683sw6hjCcvM+prZG2a22Mw+MrMJUccUhpllmtkHZrYgFv8vo44prOa+hI2ZfW5mi8xsvpnNjjqeMMws28yeM7MlsX8Do5J6ffU51M3MUoFPgJOA1cAs4Dx3/zjSwEIws2OA7cBj7n5w1PGEYWa9gd7uPtfMOgBzgDHN5b+/mRnQzt23m1k68A4wwd3fjzi0uJnZj4A8oKO7nxZ1PGHFhsPnuXuzmwRnZo8Cb7v7g2aWAbR198JkXV81hz0bAXzq7ivcvYRgzsWZEccUiru/BWyKOo694e5r3X1u7PU2gnW2cqONKn4e2B57mx57NJu7MS1hEx0z6wgcA0wGcPeSZCYGUHKoTy7wRZX3q2lGP04tiZn1B4YBMyMOJZRYs8x8YD3wirs3p/jvJljCpiLiOBrCgZfNbE5suZ3mYl+gAHg41qz3oJm1S2YASg57ZrWUNZs7v5bCzNoDfwOudfetUccThruXu/tQoA8wwsyaRdNe1SVsoo6lgY5y9+HA/wJXxppZm4M0YDhwv7sPA3YASe3zVHLYs9VA3yrv+wBrIoqlVYq11f8NeNLdp0Ydz96KNQm8CZwcbSRxaxFL2Lj7mtjzeuB5gqbi5mA1sLpKTfM5gmSRNEoOezYLGGhmA2IdQucCL0QcU6sR69CdDCx297uijicsM+tuZtmx11nAicCSSIOKU0tYwsbM2sUGMhBrkvkfoFmM2nP3/wJfmNmgWNEJQFIHYjTJVVmbCncvM7OrgOlAKvCQu38UcVihmNkU4Figm5mtBn7h7pOjjSpuRwEXAoti7fYAP3X3F6MLKZTewKOxUW8pwLPu3iyHhDZTPYHng3sM0oCn3P2laEMK5WrgydiN6Qrg4mReXENZRUSkBjUriYhIDUoOIiJSg5KDiIjUoOQgIiI1KDmIiEgNSg4iIlKDkoO0WGb2XTPzKo9yM8s3s2erTC5KZjx/NLN/JOhcB1T72+p6/MvMfmhmC81M/94lbpoEJ63B2QTLEaQC+wG3AK+Z2UHuviUZAZjZfsBlwJEJOuV/garr++cBfwRuBl6vdtw64AZgHPBwgq4vLZySg7QG893909jrd81sDfAKwQ/1v5MUw7XAAndPyIYzsbWavtwXwswOj72cVtt+F2b2GPBjlBwkTqpmSmtUubJrejIuZmZtgAuAp2r5bIGZPWJm3zezj82syMz+Y2b7mVmnWFPUOjPbbGb3xtabqs1QoBhYWsfnTwMHmlmiai7SwqnmIK1BqpmlETQr7QvcRrC/wptJuv5IIBt4u2phbM2cwUDn2Oc3xF7fH3t0JqjZnA+cQbDWzivA32u5xhDgQ3cvryOG+QRJ8WTgPw34W6SVUHKQ1qD6SqhrgNOSuDfESIJ9QBZWKz+YoPbynrufU1loZqcD3wLOdvfnYmVvAJcDB1ItOcQW9jsIeLKuANy9wswWxmIRqZealaQ1OAs4nGAt/zEESx+/aGaDG3JSM0sxs5/toamnUg6wNbbVbFXDYs8/r1beDlhYmRhisghu5jbWcv4DgEyC2sGeFMRiEamXkoO0Bh+6+2x3n+XufydoojHg1gae9xDgHK9/aeNMYFct5cOAVe5evZ9gGMEy8VUNiT0vqOU8Q/fwWVVFBElGpF5KDtLquHsRwfr4h1aWmVmHWIfvf8xssZk9UDkvwMz+bGZ3xl7nxDqOrwFeBHqa2Xwzu2MPl9xI0H9Q3TBgbtUCM+sF9KpeHju2nJpNU7DnxFFVF2BDPceIAEoO0gqZWVuC+Q4FVYqnAP929yMJ2vX3AU6JfXYLcJGZDSdICDe4+x+AqcBt7j7U3fe0v+8SIN3M+lSJIYUgOdWWBADm1VK+JJbYqhsKfBZHH8oA6h7NJLIbdUhLazDUzLoRNCX1Bq4iuIv+I4CZfQM4GuhjZr+Ofadj7Hjcfb2Z3Qu8C4xz98pZzocBz8Rx/bdizyMIJuMBDATaU3ty2A4sq6W8+rGVhsRiq1Nsu9KvAXfGEa+IkoO0Cn+t8rqAYB/hk929sl0/j2AL2B/V9mUz6wqMBbYAX8TKUglGG1W/w6/B3T83sw+A0wlqG/BVDaG25LDA3SuqXD+dYDTSY7XE1gvoQf1NSqcCJcDz9cUrAtomVAQzO4dg2Ykj3X17bNLaQHf/0Mw6Aa8C9wBlwAR3HxVrIprr7j3ivMZ3Y+fo7e47G+UP2fP1/w1scPcLk31taZ7U5yAS1CzeBOab2XyCyWr7m1k7gj6GB939CYImpDQzOxdYC8w2s6Vm9rs4rvE4kA9c0Qjx75GZDQWOA36Z7GtL86Wag0iSmNlIYLi7/ynJ1z0Z6OzuU5J5XWnelBxERKQGNSuJiEgNSg4iIlKDkoOIiNSg5CAiIjUoOYiISA1KDiIiUoOSg4iI1KDkICIiNfx/KGjZnMt9G6IAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Plot the polarization functions.\n", "\n", "fig1, ax1 = plt.subplots(1,1)\n", "\n", "for i in range(8):\n", " color = list(np.random.choice(range(256), size=3)/256)\n", " ax1.plot(1e6*tlist, signal_positive_EFG[i], label='{} mT'.format(LongitudinalFields[i]*1e3), linestyle='-', color=color)\n", "ax1.plot(1e6*tlist, signal_negative_EFG, label='{} mT'.format(LongitudinalFields[2]*1e3), linestyle=':')\n", "ax1.set_ylim([-0.1,1.1])\n", "ax1.set_xlim([0,20])\n", "ax1.set_xlabel(r't ($\\mu s$)', fontsize=16)\n", "ax1.set_ylabel(r'P$_{z}$ ($t$)', fontsize=16)\n", "\n", "plt.legend(bbox_to_anchor=(1.04,1), loc=\"upper left\")\n", "plt.show()\n", "\n", "# Plot the integrated areas.\n", "fig2, ax2 = plt.subplots(1,1)\n", "\n", "ax2.scatter(LongitudinalFields*1e3, area)\n", "try:\n", " # add interpolation if scipy is available\n", " from scipy import interpolate\n", " f = interpolate.interp1d(LongitudinalFields*1e3,area,kind='cubic')\n", " x=np.linspace(0,6)\n", " plt.plot(x, f(x))\n", "except:\n", " pass\n", "\n", "\n", "ax2.set_ylim([4,21])\n", "ax2.set_xlim([-0.5,6.5])\n", "ax2.set_xlabel(r'B$_{ext}$ ($mT$)', fontsize=16)\n", "ax2.set_ylabel(r'Area under P$_{z}$ ($t$) ($\\mu s$)', fontsize=16)\n", "\n", "plt.show()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.8.6" } }, "nbformat": 4, "nbformat_minor": 4 }