{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Analysing pore dimensions with HOLE\n", "\n", "Here we use `HOLE` to analyse pore dimensions in a membrane.\n", "\n", "**Last executed:** Feb 10, 2020 with MDAnalysis 0.20.2-dev0\n", "\n", "**Last updated:** January 2020\n", "\n", "**Minimum version of MDAnalysis:** 0.18.0\n", "\n", "**Packages required:**\n", " \n", "* MDAnalysis (Michaud-Agrawal *et al.*, 2011, Gowers *et al.*, 2016)\n", "* MDAnalysisTests\n", "* [HOLE](http://www.holeprogram.org)\n", "* matplotlib\n", "* numpy\n", "\n", "
\n", " \n", "**Note**\n", "\n", "The classes in `MDAnalysis.analysis.hole` are wrappers around the HOLE program. Please cite (Smart *et al.*, 1993, Smart *et al.*, 1996) when using this module in published work.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import MDAnalysis as mda\n", "from MDAnalysis.tests.datafiles import PDB_HOLE\n", "from MDAnalysis.analysis import hole\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Using HOLE with a PDB file" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`MDAnalysis.analysis.hole.HOLE` is a wrapper class that calls the [HOLE](http://www.holeprogram.org) program. This means you must have installed the program yourself before you can use the class. You then need to pass the path to your ``hole`` executable to the class; in this example, ``hole`` is installed at ``~/hole2/exe/hole``.\n", "\n", "HOLE defines a series of points throughout the pore from which a sphere can be generated that does not overlap any atom (as defined by its van der Waals radius). (Please see (Smart *et al.*, 1993, Smart *et al.*, 1996) for a complete explanation). By default, it ignores residues with the following names: \"SOL\", \"WAT\", \"TIP\", \"HOH\", \"K \", \"NA \", \"CL \". You can change these with the ``ignore_residues`` keyword. Note that the residue names must have 3 characters. Wildcards *do not* work.\n", "\n", "The PDB file here is the experimental structure of the Gramicidin A channel. Note that we pass `HOLE` a PDB file directly, without creating a `MDAnalysis.Universe`." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "h = hole.HOLE(PDB_HOLE, executable='~/hole2/exe/hole',\n", " logfile='hole1.out',\n", " sphpdb='hole1.sph',\n", " raseed=31415)\n", "h.run()\n", "h.collect()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This will create several outputs in your directory:\n", "\n", " - **hole1.out**: the log file for HOLE. \n", " - **hole1.sph**: a PDB-like file containing the coordinates of the pore centers.\n", " - **simple2.rad**: file of Van der Waals' radii\n", " - **run_n/radii_n_m.dat.gz**: the profile for each frame\n", " - **tmp/pdb_name.pdb**: the short name of a PDB file with your structure. As `hole` is a FORTRAN77 program, it is limited in how long of a filename that it can read. \n", " \n", " \n", " The pore profile itself is in a dictionary at `h.profiles`. There is only one frame in this PDB file, so it is at `h.profiles[0]`." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "scrolled": false }, "outputs": [ { "data": { "text/plain": [ "425" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "len(h.profiles[0])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Each profile is a ``numpy.recarray`` with the fields below as an entry for each `rxncoord`: \n", "\n", " - **frame**: the integer frame number\n", " - **rxncoord**: the distance along the pore axis in angstrom\n", " - **radius**: the pore radius in angstrom" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "('frame', 'rxncoord', 'radius')" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "h.profiles[0].dtype.names" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can then proceed with your own analysis of the profiles. " ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The pore is 42.4 angstroms long\n" ] } ], "source": [ "rxncoords = h.profiles[0].rxncoord\n", "pore_length = rxncoords[-1] - rxncoords[0]\n", "print('The pore is {} angstroms long'.format(pore_length))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Both `HOLE` and `HOLEtraj` (below) have the `min_radius()` function, which will return the minimum radius in angstrom for each frame. The resulting array has the shape (#n_frames, 2)." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([[0. , 1.19707]])" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "h.min_radius()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The class has a convenience `plot()` method to plot the coordinates of your pore." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAEMCAYAAAA4S+qsAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3hUZfbA8e9JMmkQpEUEAgSpAgpiaCsWbCAqllUXbCj6w4KLqNjdVWxrwbKKDRXbKuiuYFlAsYCNpYTepUMA6T19cn5/zA1GTMIEZnLvTM7neebJbXPvmUuYk7fc9xVVxRhjjAlGjNsBGGOMiRyWNIwxxgTNkoYxxpigWdIwxhgTNEsaxhhjghbndgDhVrduXU1PT3c7DGOMiRizZs3apqqppe2L+qSRnp5OZmam22EYY0zEEJG1Ze2z6iljjDFBs6RhjDEmaJY0jDHGBC3q2zRKU1BQQFZWFrm5uW6HUqbExETS0tLw+Xxuh2KMMQdUyaSRlZVFSkoK6enpiIjb4fyBqrJ9+3aysrJo2rSp2+EYY8wBVbJ6Kjc3lzp16ngyYQCICHXq1PF0ScgYUzVVyaQBeDZhFPN6fMaYqqnKJg1jjIlWUz+fyUdPfxaWc1vScMmXX35Jq1ataN68OU8++aTb4Rhjosj0/85i7D/Hh+XcljRc4Pf7GTRoEBMnTmTx4sWMHj2axYsXux2WMSZK+AuLiI0Lz9e7JQ0XzJgxg+bNm3PssccSHx9P3759+eyz8BQljTFVj9/vJzYuNiznrpJdbkt6ZcjbrJy3JqTnbNY+nVteuK7M/Rs2bKBRo0YH1tPS0pg+fXpIYzDGVF3+wvAlDStpGGNMlAln9VSVL2mUVyIIl4YNG7J+/foD61lZWTRs2LDS4zDGRCcraUSZTp06sXz5clavXk1+fj5jxoyhT58+bodljIkS4UwaVb6k4Ya4uDhGjBhBz5498fv9DBgwgLZt27odljEmSlj1VBTq3bs3vXv3djsMY0wUsuopY4wxQSuypGGMMSZY9nBfGKiq2yGUy+vxGWO8y1/oJyYaGsJFZBRwPrBFVds52z4CWjmH1AR2qWqHUt67BtgL+IFCVc043DgSExPZvn27Z4dHL55PIzEx0e1QjDERKJp6T70DjADeK96gqn8pXhaRZ4Hd5by/h6puO9Ig0tLSyMrKYuvWrUd6qrApnrnPGGMqKmp6T6nqDyKSXto+CfzJfzlwRrjj8Pl8NiOeMSZqVZXeU6cAm1V1eRn7FZgkIrNEZGB5JxKRgSKSKSKZXi5NGGNMOFSVpNEPGF3O/u6q2hE4FxgkIqeWdaCqjlTVDFXNSE1NDXWcxhjjaVHfe0pE4oBLgI/KOkZVNzg/twDjgM6VE50xxkQWf0Fh1Jc0zgKWqmpWaTtFpJqIpBQvA+cACysxPmOMiRj+wiJiY6MgaYjIaOB/QCsRyRKR651dfTmoakpEGojIBGe1HvCTiMwDZgDjVfXLyorbGGMiib/QT5wvCrrcqmq/MrZfW8q2jUBvZ3kV0D6swRljTJQI58N9XqmeMsYYEyJR3xBujDEmdKpKl1tjjDEhYEnDGGNM0Kx6yhhjTFBU1UoaxhhjglNUVARgScMYY8yh+QstaRhjjAmSv9APYG0alW3nlt3s373f7TCMMaZCig4kDStpVKqr0m/mw8fHuh2GMcZUiFVPucSX4KMgv9DtMIwxpkKsesolvvg4CvIK3A7DGGMqxG/VU+7wJfgoyLOShjEmshRXT9mAhZUsLj6OwgJLGsaYyGLVUy7xJcSRb9VTxpgIY9VTLvEl+Ci0hnBjTIQpLLCk4QprCDfGRCKrnnKJNYQbYyKRPafhkrj4OHtOwxgTcYpLGuGaI7xSk4aIjBKRLSKysMS2h0Vkg4jMdV69y3hvLxFZJiIrROTecMfqS7DqKWNM5Im2YUTeAXqVsv15Ve3gvCYcvFNEYoGXgXOBNkA/EWkTzkAD1VOWNIwxkSWqqqdU9Qdgx2G8tTOwQlVXqWo+MAa4MKTBHcQXH2e9p4wxEaeqNITfKiLzneqrWqXsbwisL7Ge5WwrlYgMFJFMEcncunXrYQVkDeHGmEhUnDSi+YnwV4FmQAdgE/DskZ5QVUeqaoaqZqSmph7WOeJ81hBujIk8UVU9VRpV3ayqflUtAt4gUBV1sA1AoxLrac62sLGGcGNMJDpQPRUbpdVTIlK/xOrFwMJSDpsJtBCRpiISD/QFPg9nXNYQboyJRAeeCA9Tl9u4sJy1DCIyGjgdqCsiWcBDwOki0gFQYA1wo3NsA+BNVe2tqoUicivwFRALjFLVReGM1RrCjTGRyO8MtBqu5zQqNWmoar9SNr9VxrEbgd4l1icAf+iOGy6+BB+FBX6KioqIiXG9QGaMMUEpbtOI84Xn692+DcsQFx+44cVFPWOMiQTFUzpEe5dbz/El+ACsXcMYE1F+a9Owkkal8iUEbrglDWNMJPEXRNHYU5HEV1w9ZY3hxpgIYpMwueS36ilLGsaYyBHuLreWNMpQ3BBuT4UbYyJJYZi73FrSKIM1hBtjIlGRV4cREZFqzpDlUSneGsKNMRGosKCQmBgJ2/NlQZ9VRGJE5AoRGS8iW4ClwCYRWSwiz4hI87BE6JI4awg3xkSgwgJ/2LrbQsVKGpMJjEZ7H3CMqjZS1aOB7sA04CkRuSoMMbqiuHoqP9dKGsaYyOEvKAxbewZUbBiRs1T1D9+gqroD+AT4RER8IYvMZQlJ8YAlDWNMZPEXFoWtPQMqUNIoLWEUE5Gahzom0sQ7SSMvJ9/lSIwxJniFbpc0ROQk4HzgRaAQaHvQqx1QDagZtihdcKCkYUnDGBNBAm0a7lZPvU5guPJ1wF5gEYFG8CUE5rXooKpbwhahS6ykYYyJRH6/P2wj3EJwSWMqcBcwG0gG3lDVjwFE5K5oTBhgJQ1jTGTyF/jDNsItBJE0VHWwiCSraraI1AYeFJHbgUcITJwUlRKspGGMiUCe6HKrqtnOzx2qegeBaqkrgHoi0iNs0bnIl+BDRMjLyXM7FGOMCVq4u9weVhlGVdeq6tXAycA9IvJ9aMNyn4gQn+iz6iljTEQJd5fboMswIiKq+rvqKFWdC/QqLm2Udkwki0+Kt+opY0xECXeX2wo9ES4ifxWRxiU3ikg8ECMi7wL9yzuBiIwSkS0isrDEtmdEZKmIzBeRccXPfJTy3jUiskBE5opIZgXiPmyJyQlW0jDGRJRwd7mtSNLoBfiB0SKy0RlzahWwnEAbxwuq+s4hzvGOc56SvgbaqeoJwC8EhikpSw9V7aCqGRWI+7DFJ8WTl2tJwxgTOfyF7ne5BUBVc4FXgFec4ULqAjmquqsC5/hBRNIP2japxOo04NJgzxduCUnx5GVb0jDGRA5/gR9fNW8MWHiAqhao6qaKJIwgDQAmlnVZYJKIzBKRgeWdREQGikimiGRu3br1sIOxNg1jTKTxRJfbyiAiDxAYpuSDMg7prqodgXOBQSJyalnnUtWRqpqhqhmpqamHHVNCUry1aRhjIoqXGsLDRkSuJTC+1ZVl9b5S1Q3Ozy3AOKBzuOOykoYxJtIUeWWU23ARkV7A3UCf4ocISzmmmoikFC8D5wALSzs2lKykYYyJNJ4paYjI1SKyVUSyRKS/s62riDwmIrOCPMdo4H9AK+c81wMjgBTga6c77WvOsQ1EZILz1nrATyIyD5gBjFfVL4P+lIcpwUoaxpgI44VRbov9HegNrAZuFZGvgdbAaGBIMCdQ1X6lbH6rjGM3OtdDVVcB7SsQa0jEJ1rSMMZEFn+hn7g4D3S5Bfap6kwAERkGbAZahqEHlWdY9ZQxJtK4PsptCcc4XV2XOa+saE4YAAnJVtIwxkSWcHe5rciZHwKOB650fqaIyDfAHGCOqn4YhvhclVgtkYK8AvyF/rD2RjDGmFBxfbrXYqo6suS6iKQRSB4nEHh2IuqSRnJKEgDZe3NIqVXd5WiMMebQwt3l9rDLMKqaBWRR9hPcES8pJRGAHEsaxpgI4Zkut1XRbyWNXJcjMcaYQ1NVT41yW+UklaieMsYYryvyFwGEdZTbCicNEbmsxNPZD4rIWBHpGPrQ3JdconrKGGO8rrCgEMBzw4j8TVX3ikh34CwCD+e9GtqwvCHJqqeMMRGkIC+QNBKS4sN2jcNJGn7n53nASFUdD4QvQhcVt2lYScMYEwnynUnjfIm+sF3jcJLGBhF5HfgLMEFEEg7zPJ5X3HvK2jSMMZEgP7cAgHiPJY3Lga+Ans4T4bWBu0IalUf8VtKw6iljjPf9ljTCV/lT4SZ2Z/jysSXWNwGbQhmUV/gSfMTGxVpJwxgTEYqrp8JZ0qhw0hCRv5e2XVUfOfJwvEVESE5JtDYNY0xEqIzqqcPpzLu/xHIigRn3loQmHO9JSkkie58lDWOM9xWPyu216qlnS66LyHACbRxRKTklydo0jDERobik4bXeUwdLBtJCcB5PSkpJtDYNY0xE8GqbxgJAndVYIBV4NJRBeUlSSpK1aRhjIoIne08RaMMoThqFwGZV9ZdzfERLTklix8adbodhjDGH5KnnNETkJ2dxYYnXUmCniOypwHlGicgWEVlYYlttEflaRJY7P2uV8d7+zjHLRaR/sNc8ElY9ZYyJFAWVUD0VdNJQ1e7OzxRVrXHwqwLXfAfoddC2e4FvVbUF8K2z/jsiUpvA7IFdgM7AQ2Ull1BKrm7VU8aYyFAZ1VOVPvyHqv4A7Dho84XAu87yu8BFpby1J/C1qu5Q1Z3A1/wx+YRcUkoS2XtzUdVDH2yMMS7y1HMaInJHeftV9bkjiKOe82Q5wK9AvVKOaQisL7Ge5Wz7AxEZCAwEaNy48RGEFWjT8Bf6KcgrCGv2NsaYI3VgwMIED1RPASnOKwO4mcAXdkPgJiBk82lo4E/6I/qzXlVHqmqGqmakpqYeUTw2aKExJlLk5xYQn+hDRMJ2jaBLGqo6DEBEfgA6qupeZ/1hYPwRxrFZROqr6iYRqQ9sKeWYDcDpJdbTgClHeN1DKjloYc3Uo8J9OWOMOWz5uflhrxE5nDaNekB+ifV8Sq9OqojPgeLeUP2Bz0o55ivgHBGp5TSAn0MlPIluU74aYyJFgVPSCKfDeU7jPWCGiIxz1i/it0bsQxKR0QRKDHVFJItAj6gngY9F5HpgLYHh1xGRDOAmVb1BVXeIyKPATOdUj6jqwQ3qIWdTvhpjIkV+ngeThqo+LiITgVOcTdep6pwKvL9fGbvOLOXYTOCGEuujgFEVCPeI2ZSvxphIkZ+bjy/M1VOHU9IAWO28NxFIEZFTna60UcemfDXGRIp8L1ZPicgNwG0EGqLnAl2B/wFnhDY0byiuntq/x5KGMcbbKiNpHE5D+G1AJ2CtqvYATgR2hTQqD6lZryZxvlg2rfzV7VCMMaZcedl5JCR5r/dUrqrmAohIgqouBVqFNizviE/w0fT4xizLXOl2KMYYU67sPTkk10gO6zUOJ2lkiUhN4FPgaxH5jECPp6jVMqM5v2SupKioyO1QjDGmTNl7skmukRTWa1QoaUjgMcPBqrpLVR8G/ga8ReljRUWNdt1bs393NqPu/9DtUIwxpkz79+Qc6LwTLhVKGs4QHxNKrH+vqp+ran45b4t4PfqdTO8bzuSjpz9j6mczD/0GY4ypZKrq2eqp2SLSKeSReFhsbCy3jrieZh3SeX7ga+zcHLXt/saYCJWfm4+/0E81L1VPOboA/xORlSIyX0QWiMj8UAfmNb54H/e+P5j9e3J4buBrNlS6McZTsp3HArxY0ugJNCPwXMYFBKZ/vSCUQXlVettG3PCPK5n2xSy+HPWd2+EYY8wB+w8kDY+VNFR1bWmvcATnRRcNPpcOZ7TjlSFvs9Ge3TDGeET2nmzAg0mjqouJieGutwcRGxfLU/1H4Pf73Q7JGGMOVE9V82D1VJV3dKO6/HXEDSyeuoyPn/7c7XCMMaZEm4aVNDzpjCu6c9rl3Xjv4Y9YMWe12+EYY6q44jl/krzynIaI3F1i+bKD9j0RyqAigYgw+JX/o0bdGjx59YsH5uY1xhg3/FY95ZGkAfQtsXzfQft6hSCWiFOjdgpDR93C2sVZ9rS4McZVxQ3hnilpAFLGcmnrVUannh3oc0tPPnlhPLO/XeB2OMaYKio3Ow8R8dQot1rGcmnrVcr/PX01aS3rM/y6l9m3a7/b4RhjqqDc/XkkJMcTGCIwfCqSNNqLyB4R2Quc4CwXrx8fpvgiQmJyAve+P5jtm3YyYvBbbodjjKmC8rLzSExOCPt1gk4aqhqrqjVUNUVV45zlGqqaAtx9yBOUQ0RaicjcEq89IjLkoGNOF5HdJY75+5FcM9RadWrOVX+7lG//9SPffzzV7XCMMVVMXk4+CV5KGodw+5G8WVWXqWoHVe0AnARkA+NKOfTH4uNU9ZEjuWY4XHH/JbTu3Jx/3jySbRu2ux2OMaYKyc0OVE+FW6iSRigr0c4EVkbi0CSxcbHc895fKcgrZPj1r9qghsaYSpOXnUditcSwXydUSSOU3459gdFl7OsmIvNEZKKItC3rBCIyUEQyRSRz69atIQzt0NJaNmDg8GuYNWken7/yVaVe2xhTdeVl53urpCEie0s0fu85qCG8QSiCEZF4oA/w71J2zwaaqGp74CUC082WSlVHqmqGqmakpqaGIrQKOf/Gs+l07om8cff7rFu6odKvb4ypenL353quITylRON3yVeKqsaFKJ5zgdmqurmU6+9R1X3O8gTAJyJ1Q3TdkBIR7nzzZhKSE3jqmpcoLCh0OyRjTJQLlDQ8lDSKiUiiiLRzXqGuQOtHGVVTInKMM0c5ItKZQOyebW2uU78WQ16/kV8yV/LBY5+4HY4xJsp5riFcROJE5GkgC3gXeA9YLyJPi4jvSAMRkWrA2cDYEttuEpGbnNVLgYUiMg94EeirHm9pPuWSLpzd/zQ+fGIsS6YvdzscY0wUCzyn4a2G8GeA2kBTVT1JVTsSmMGvJjD8SANR1f2qWkdVd5fY9pqqvuYsj1DVtqraXlW7qmpEPAwx6IXrSE2rw1PXvETO/ly3wzHGRKm87HwSvVTSIDCt6/+p6t7iDaq6B7gZ6B3qwKJFtaOqcdc7g9i44lfeuOt9t8MxxkQhVSV3f67n2jS0tOogVfVTxceeOpT2p7Xlsjsv4IvXJjFj4hy3wzHGRJmC/EKKitRzSWOxiFxz8EYRuQpYGrqQolP/R/vS9PjGPHv9K+zetsftcIwxUSQvOw/AW11ugUHAIBGZIiLPOq/vgcEEqqhMOeITfNz7/mD27tjHCzeNtKfFjTEhcyBpVPNQ0lDVDaraBXgEWAOsBoapamdVtSfYgnDsCU249tG+/DR2OpPeneJ2OMaYKJGbHZg5tDKqp4J+KE9EXuKPbRctReQiAFUdHMrAotWf7zif6RNm8/LgUbTr3pqGzeu7HZIxJsJtWL4JgOo1k8N+rYpUT2UCs5xXnxLLxS8ThNjYwKCGcb5YnrzqRXta3BhzRIqKivjw8U+o27A2J551QtivV5HqqXeLX8DOkuvONhOkoxvVZcjrN7J0xgreH1baMFvGGBOcT577L4v/9wvXPtqX+IQjfs76kA53lFtrxT1Cp17ajV4DzmD0P8Yx7/tFbodjjIlAK+etYdQDH3LyxZ05p//plXLNUA2Nbg7DLS9cS4Pmx/DU1S+xd+c+t8MxxkSQvJw8nrzqRWrUrcHtr98Y9rnBix3W0OgcNEe4s81UUFL1JO774DZ2/LqL52983brhGmOC9uY9H7Bm0XqGjrqFo+rWqLTrHu7Q6L+bI1xVKy/iKNMqoxnXPtqXH/8zja/enux2OMaYCDD1s5l8OmIiFw/uTaeeHSr12lY95QGX39WHDme04+XbRpHldJ0zxpjSbFm3leEDXqZFx6bc8NRVlX59SxoeEBMTwz3v3oovwccTV7xAfl6B2yEZYzyosKCQx6/4J/7CIh4Yc3ul9JY6mCUNj6jbsA53vnkzy2et4s17/uV2OMYYD3r3oY9ZPHUZQ16/0bUHgy1peMjJF3XmktvOY9yLE/hx7HS3wzHGeEjmpHmMeXIc515/Jj36nuxaHJY0POaGp66kVadmPHv9K2xa9Yep0o0xVdD2TTt56pqXSG/biFv+eZ2rsVjS8BhfvI8HxtyOiPBY3+etfcOYKs7v9/Pk1S+SszeHBz+6vVKGPy+PJQ0Pqt+0HkNH3cIvmStttj9jqrjRT4xj7ncLufWl62nSppHb4XgraYjIGhFZICJzRSSzlP0iIi+KyAoRmS8iHd2IszKcfFFn/jzkPD4dMZEfP5nmdjjGGBfM/2Ex7w/7mDOu6E7P63q4HQ7gsaTh6KGqHVQ1o5R95wItnNdA4NVKjaySXf/klbTu3Jzh1r5hTJWzc/MunrjiBY45th63vTqw0oYJORQvJo3yXAi8pwHTgJoiErUTUhS3b8TExPDoX56z9g1jqgh/oZ/H+73A3h37+NvHd5CckuR2SAd4LWkoMElEZonIwFL2NwTWl1jPcrb9jogMFJFMEcncunVrmEKtHMekH81dbw9i+axVjBz6ntvhGGMqwVv3fcC8KYsY8tqNNO/Q1O1wfsdrSaO7qnYkUA01SEROPZyTqOpIVc1Q1YzU1NTQRuiCP13YiT/ffj6fvfwlP/znf26HY4wJo+///T/+/ewXXHBzT86+5jS3w/kDTyWN4rnGVXULMA7ofNAhG4CS3QfSnG1R7/p/XEHrLi149oZX2bDCxqcyJhqtXZLF8AEvc1zXFtz8fH+3wymVZ5KGiFQTkZTiZeAcYOFBh30OXOP0ouoK7FbVKvEN6ov38eCY24mNjeGRS58lZ3+u2yEZY0Jo/55shl3yDInVEvn7v+/EF1/540oFwzNJA6gH/CQi84AZwHhV/VJEbhKRm5xjJgCrgBXAG8At7oTqjnpNUrnvg9tYvWAdzw98zebfMCZKqCrDB7zChhW/8uCY26nbsI7bIZUpzu0AiqnqKqB9KdtfK7GswKDKjMtrOvU6kese68eoBz6k5UnNuPSOC9wOyRhzhD5+5nN+Gjudgc9cQ/vT27odTrm8VNIwQep770Wc8ucuvHH3+8z5boHb4RhjjsDsbxcw6v4POO3yblx6x/luh3NIljQikIgwdNQgGrVuyGN/eZ7NayO7W7ExVdWW9dt4ot/zpLVqwJ1v3uyZB/jKY0kjQiWnJPHwuLvxF/p5+JJnyMvJczskY0wF5OcV8Ohlz1KQV8jDY+8iqbp3HuArjyWNCJbWoj73/WswK+eu4fkbX7eGcWMiyCu3vc3SGSu4651BNGr1h2eUPcuSRoTrct5JXPPw5Xz7rx/56OnP3A7HGBOEL179ivEjv+Yvd19I94u7uB1OhXim95Q5fFc8cAnrlmbx1n0fkJpWhzOvPMXtkIwxZZg+YTYj/voWXc7ryHWP9XM7nAqzpBEFYmJiGDpqEDs27WL4gJepdUxNOp55vNthGWMOsmLOah77y3M065DOA6OHEBsX63ZIFWbVU1EiPsHHw2PvIq1VA4b9+RlWzV/rdkjGmBK2rN/GA+f/g5Ta1Xn0i/sipuH7YJY0okj1mtV4fPz9JKckcX/vx9myzrriGuMF+3fv54HzniB3fy6Pj7+fOvVruR3SYbOkEWWOblSXJybcT86+XO7v/QR7d+5zOyRjqrT83Hwe/vNw1i/dyEP/GUrTdo3dDumIWNKIQk2Pb8KwcXezYfkmHr7kGZu8yRiX+Av9PHHlP5n73UKGjrqFjmed4HZIR8ySRpTq0KMdd71zK/O/X8zT/V+iqKjI7ZCMqVJUlRdufJ2fx83glheu46yrDmt6IM+x3lNR7Ix+3dmWtZ037vkXNWqn8NeXb4iIYQqMiXSqypv3/Isv357MVX+7lIsH93Y7pJCxpBHlLhvah93b9vLxM4EH/24dcT0xMVbANCZcVJW3HxzNx8M/p88tPbnm4cvdDimkLGlEORHhhievRAQ+evozioqUwa/cYInDmDBQVUbe9T7/ee4Let9wJoNeHBB1pXtLGlWAiHD9P65EYmIY8+Q4ivxFDHl9oCUOY0JIVXnltrf5dMRE+tzSk0EvDojK/2OWNKoIEWHA4/2IiRE+fGIsfr+fO0beFJFPpBrjNUVFRbx48xuMf+Mb/nz7+dw4/JqoK2EUs6RRhYgI1z7alzhfHO8N+5hdW3bz4JjbI/bJVGO8ID83n2cGvMKUMT/T996LGfB4v6hNGGBdbqscEeHqhy7jtlcHkvnlXIaeMYydW3a7HZYxEWnX1t3cddYjTBnzM9f/48qoTxjgkaQhIo1EZLKILBaRRSJyWynHnC4iu0VkrvP6uxuxRovzbzybh8bexdpF6xly8gNsWLHJ7ZCMiShZyzcxuNsDrJi9ir99fAd977ko6hMGeCRpAIXAnaraBugKDBKRNqUc96OqdnBej1RuiNHnT3068fS3D7FvVzZDTn6QFXNXux2SMRFh3dIN3Hn6Q+TszeHZKcM49dJubodUaTyRNFR1k6rOdpb3AkuAyJnKKoK16dqSf/78GL5EH3efOYxlM1e4HZIxnrZ28XqG9ngILSpi+OSHad25hdshVSpPJI2SRCQdOBGYXsrubiIyT0Qmikjbcs4xUEQyRSRz61Yb6fVQ0lo24LnvH6FazWrcffYjLJq6zO2QjPGk1QvXMbTHw0hMDMMnD6NJm0Zuh1TpPJU0RKQ68AkwRFX3HLR7NtBEVdsDLwGflnUeVR2pqhmqmpGamhq+gKPIMelH8+yUYdSqV5N7ez7KvO8XuR2SMZ6yesFa7j5zGLG+WJ6d/DCNW1fNyhDPJA0R8RFIGB+o6tiD96vqHlXd5yxPAHwiUreSw4xqRzeqy7NThnF047o80PsJZn09z+2QjPGEuZMXcvupfycuPo7hk4eR1rKB2yG5xhNJQwJdDt4Clqjqc2Ucc4xzHCLSmUDs2ysvyqqhTv1aDJ88jIYt6vO3C55k0rtT3A7JGFd99+GP3NfrMeo0qMULPz1GWov6bofkKk8kDeBk4GrgjBJdanuLyE0icpNzzKXAQhGZB7wI9FVVdSvgaFbr6KN45ruHaNu9Nc9c9zIv3zaKwoJCt8MyplL5/X7efegj/ilTIEsAAA/dSURBVHHVi7T5Uyte+Okx6jWx6m6J9u/djIwMzczMdDuMiOQv9PPG3e/zyQvjOf7U47jvX7eRmlbH7bCMCbttG3fw5FUvMm/KIs7ufxpDXruR+ASf22FVGhGZpaoZpe3zSknDeFBsXCw3PXct97z3V5bPWsXAE+7k+4+nuh2WMWE1/4fF3NRhKMtmrOCutwdx99u3VqmEcSiWNMwhnXXVqbw25xnSWjXgsb7P81T/l9i3a7/bYRkTUru37eG1O97hnrMfoUbdGoyY+STn9D/d7bA8x6qnTND8hX4+eOwTPnj8E+ITfZx11an0GdSLpu0aux2aMYdty/ptjPvnBP772iTyc/M5p//p/N/TV1OjTorbobmmvOopSxqmwlbOW8PYf45nypifyc8toF331nS7IIOOZ53Ase2beHIOAVVl7uSF/PjJdGrUrk6zDum079GWGrXD88WQn1dA5pdzSayWQKtOzah2VLWwXKc8qsreHftYMXcN+3dn06lXBxKTE8Jyrfy8An6ZuYIFPy5l+8YddOsT+H3w4lhMqsqaheuY+lkmP382g+WzVhETI/S4ojv97ruEJseluR2i6yxpWNIIiz3b9/LV25OZ9N4U1ixcD8BRdVNo3aUF9Y+tR2paHaodlUzzE5vS/MSmlT53x+5te1g+ezULfljM1M9msmbRehKTE8jPzaeoSImJERq3SeOY9KNJSI7HX1hEnfq1aPOnVpx6aVfifBWbOaCoqIil05fz3Yc/MeWjn9m9be+BfbXr16JJmzQyzmlPl/NPonHrhof8Qi0sKGTp9OX4/UW06dYSX/yh69V3/Loz8GX46XQW/LCEvJz8A/uSayRxRr/unHX1abTu3LzC/x6FBYX8PG4GC39aytYN2/El+Mjdn8uWddtYv2QDBfmBHna++DgK8gtpfFxDTr6oMyec1pY23VqSnFK5Q/D7/X42rdzMmkXr2bRyM1uztvPrmi2smLOaresDvfWP69qCP13YmdMu60b9Y+tVanxeZknDkkbYbd+0k9nfzGf2N/NZNW8tm1ZtJmdf7oH9dRrU4s+3X8D5N54Vlvk7fl2zhcyv5rFi9irWLd3AuiVZB760Y2JjaNWpGecNPJsefU9GYoRlM1cy++v5LJ+zii3rtlGQW4DECNuydpC9N4d6TVK59I4L6DmgB0nVEsu8btbyTfz0yTRmfTOf5bNWsX93NvGJPrpekEGvAWcgIiyftYqs5RtZMXs1q+avBaBBs3p0Oe8kOvXqQFrLBsTFx5GXk8/2jTtYPX8dc6csZM63C8jekwNAYrUEOvRox0nntKd5h3TqNKiNLyGOgrxCfl2zhfnfL2b2N/NZMm05qkqDZvXo3Lsj9Zqk0qRtI2LjYpn07mR+/M808nMLSE5Jom331rTt1oou53WkWYf0MpNYzv5cJr75Lf957gu2rt9OUvVE6jVJpSC/kMRqCdSuX4umbRtxXLeWHH/KcSSlJPH9x1MZP/JrlkxbTpG/iJjYGJoe35gmbdJIa9mARq0a0Pi4NNLbNQp5yTQvJ49v3v+BMU99yq+rtxzYXhx3ertGdOjRjq4XZFCnfq2QXjtaWNKwpFHpVJWcfbns2b6XJdOWM+HNb5j73UKSU5I4+eLO9OjXnY5nHl/hv3b9hX7ycvLZsWknK+asZvH/fiHzq7msX7YRgBp1Umh8XEMatWpI4+Ma0vT4xhzXNfi/clWVGRNm8+E/xrF46jJSalWjdZcW1Kx3FNVqJINCQV4Bm9ZsYd3iLLZt2AFA8xOb0rpzc9qe3JpufTICx5Ziy/ptTP/vLKaNn8WcbxdSkFdQ6nH1mqRy4pnH0/ncE4mNiyXzq7lkTprHplWbSz0+JkZomdGMzr070v3izqS3a1xqEti3az+zJs1jzrcLWDR1GWsWBUqIKbWrk96uEQ2OPYb4xECJZu/OfezauoeVc9ewd8c+jj/1OC4feiGde58Y9Bd9zr4clkxbzvwfFrN0xgqylm1ky7ptFH/v1K5fK1AaOfU4jm2fTmpabeKT4omNreDvhd/PvCmLmTz6J34aO519u/bTMqMZF9x0Dk2Pb0xay/quVBFGKksaljQ8Ycn05UwY+TU/jp3O/t3ZpNSqxknntKdTrxPp1KsDterVPHDsnu17WTFnNXMnL2TOdwvZsm4b+3buIz/391+y8Yk+2p/eloyeHeh87ok0bFE/ZPXoi6Yu44tXv2L9so3s3LyL7D05xMQIsb44UhvVoUmbNFp3bkHn3idSv2nFqzZy9ueydPpyNq/dRlGhH1+ij5pHH0WTNmkc3aj0EXI2r93K2sVZ7Ny8C3+Bn5i4WFLTatOqU3Oq16z4l+LOzbuYMXEOi6cuY83iLDav2YK/wE9RkZJSuzpHpdbgmPRU+tzSi7Z/alXh85cmLyePDct/ZeXcNUz9fCYzJ875XTUaBP5dk1OSSG1cl/antaXDGe1o+6dWv/uM+/dkM//7QNXjtC8y2bV1D0nVE+l+SRd6XteDE05t48k2lUhgScOShqfk5xUwc+Icpn4+k8wv57Lj110A1EytQbWa1cjZm3NgW0xsDK27tKDJcWlUr5lMUkoSickJ1KibQrP26TRuk2Z96CNcfl4B65ZksWreWnZt2U1edj45+3LI3ptL1i8bWTx12YH2kpTa1fEl+Cgq9LNra2BM0+QaSU4Jqwtdz+9IQlJ4GvurEksaljQ8q6ioiFXz1jLr6/lsWvkr+/fmkJgUT6PWDTm2fTrHdW1RZlWPqRrycvJY9PMyls9eza+rN1PkL0JEOKbp0bQ46VhOOK1NUJ0ETPAsaVjSMMaYoNkwIsYYY0LCkoYxxpigWdIwxhgTNEsaxhhjgmZJwxhjTNAsaRhjjAmaJQ1jjDFBs6RhjDEmaFH/cJ+IbAXWuh0HUBfY5nYQHmX3pnx2f8pn96d8h3N/mqhqamk7oj5peIWIZJb1hGVVZ/emfHZ/ymf3p3yhvj9WPWWMMSZoljSMMcYEzZJG5RnpdgAeZvemfHZ/ymf3p3whvT/WpmGMMSZoVtIwxhgTNEsaxhhjgmZJI4xE5BkRWSoi80VknIjULLHvPhFZISLLRKSnm3G6RUQuE5FFIlIkIhkH7avy9wdARHo592CFiNzrdjxuE5FRIrJFRBaW2FZbRL4WkeXOz1puxugWEWkkIpNFZLHz/+o2Z3tI748ljfD6GminqicAvwD3AYhIG6Av0BboBbwiIrGuRemehcAlwA8lN9r9CXA+88vAuUAboJ9zb6qydwj8TpR0L/CtqrYAvnXWq6JC4E5VbQN0BQY5vy8hvT+WNMJIVSepaqGzOg1Ic5YvBMaoap6qrgZWAJ3diNFNqrpEVZeVssvuT0BnYIWqrlLVfGAMgXtTZanqD8COgzZfCLzrLL8LXFSpQXmEqm5S1dnO8l5gCdCQEN8fSxqVZwAw0VluCKwvsS/L2WYC7P4E2H0ITj1V3eQs/wrUczMYLxCRdOBEYDohvj9xRxSZQUS+AY4pZdcDqvqZc8wDBIqOH1RmbF4QzP0xJlRUVUWkSj9HICLVgU+AIaq6R0QO7AvF/bGkcYRU9azy9ovItcD5wJn620MxG4BGJQ5Lc7ZFnUPdnzJUmftzCHYfgrNZROqr6iYRqQ9scTsgt4iIj0DC+EBVxzqbQ3p/rHoqjESkF3A30EdVs0vs+hzoKyIJItIUaAHMcCNGj7L7EzATaCEiTUUknkDngM9djsmLPgf6O8v9gSpZgpVAkeItYImqPldiV0jvjz0RHkYisgJIALY7m6ap6k3OvgcItHMUEihGTiz9LNFLRC4GXgJSgV3AXFXt6eyr8vcHQER6Ay8AscAoVX3c5ZBcJSKjgdMJDPe9GXgI+BT4GGhMYBqEy1X14MbyqCci3YEfgQVAkbP5fgLtGiG7P5Y0jDHGBM2qp4wxxgTNkoYxxpigWdIwxhgTNEsaxhhjgmZJwxhjTNAsaRhjjAmaJQ1jjDFBs6RhjAmaiLwkIrNFpJPbsRh3WNIwxgRFRKoBRwM3EhhPzVRBljSMqWQi8rCIDC2xPvUwz1NTRG45jPf1E5G5IrJARFRE5pRyTJKIfF9y8itV3Q/UB6YAL4pIvIj8ICI28GkVYknDRD0JcOV3PZhrq+qfDvP0NYEKJw1VHQ1cBuwlMMDdeaUcNgAYq6r+4g0iUgdIdt5X6EwM9S3wl4qHbiKVJQ0TEUQk3Zlv/QMRWSIi/xGRZGffHSKy0HkNKXH8MhF5j8C0so1E5CoRmeH8lf16aVPIisg1zpzu80Tk/RLb/3CNCl77ARH5RUR+AloddM19znuWiMgbzvzOk0QkqcQxn4rILGffQGfzk0Az5/M84xwXzGdsS2Dk00dU9QZV3VjKLb+SP46G+iAwHFhEYCpeCAwWeGUp7zfRSlXtZS/Pv4B0QIGTnfVRwFDgJAKjelYDqhP4QjvROb4I6OocfxzwBeBz1l8BrjnoGm0JzOVe11mv7fws6xrBXrv4uGSgBoHpa4eWuO4+5z2FQAdn28fAVSWOKY4liUAiquO8Z2GJY4L5jD5gFnBKOfc6Hvi1lPs/HRBgBPB/zvZYYKvbvx/2qryX1UWaSLJeVX92lv8FDAYKgHEaqG9HRMYCpxD4S3qtqk5zjj+TwJf3TGcmsyT+OBnNGcC/VXUbgP42fHT3Mq4hQV77FOe4bOe4subEWK2qc53lWQS+qIsNdoaSh8DETC0ITN1ZUjCfsRewSFV/LCMGCAw7vuugbY8RKJmoiCzBKWmoql9E8kUkRQPzUpsoZ0nDRJKDx/E/1Lj++0ssC/Cuqt4X2pCCunaw8kos+wl86SMipwNnAd1UNVtEpgCJpbw/mM/YGZh8iDhySp5fRDoAlwDdReRlZ9+CEscnALmHOKeJEtamYSJJYxHp5ixfAfxEYNKZi0Qk2ekSerGz7WDfApeKyNEAIlJbRJocdMx3wGVOgy8iUtvZXtY1gr32D85xSSKSAlxQwc99FLDTSRitga7O9r1ASgU/416gG+VQ1Z1ArIgUJ46nCMw+ma6q6UB7nJKGc6+2qWpBBT+TiVBW0jCRZBkwSERGAYuBV50v0nf4bTrYN1V1joikl3yjqi4WkQeBSU5vpgJgEIGZzIqPWSQijwPfi4gfmANcq6qzS7sGQJDXni0iHwHzCFQXzazg5/4SuMmpFloGTHPOu11EfhaRhcBEVb3rUJ8ReAMY5ZxrJXClqu4u5ZqTCJQsioBkVf2mxOfZLCLVnaTaAxhfwc9jIpjN3GcigvNF/F9VbedyKFFDRMYTaKeYXsq+jsDtqnr1Ic4xFrhXVX8JU5jGY6x6ypgqSESuJDB3famlHlWdDUwurctuiXPEA59awqharKRhjDEmaFbSMMYYEzRLGsYYY4JmScMYY0zQLGkYY4wJmiUNY4wxQbOkYYwxJmiWNIwxxgTt/wGhPKvBbmmJAAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "h.plot();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can also create a VMD surface from the `hole1.sph` output file, using the `create_vmd_surface` function." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "'hole1.vmd'" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "h.create_vmd_surface(filename='hole1.vmd')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To view this, open your PDB file in VMD.\n", "\n", "`vmd tmp*/*.pdb`\n", "\n", "Load the output file in Extensions > Tk Console:\n", "\n", "`source hole1.vmd`\n", "\n", "Your pore surface will be drawn as below.\n", "\n", "
\n", "
\n", "\n", "![sphpdb.png](sphpdb.png)\n", " \n", "
\n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "MDAnalysis supports many of the options that can be customised in HOLE. For example, you can specify a starting point for the pore search within the pore with `cpoint`, and a ``sample`` distance (default: 0.2 angstrom) for the distance between the planes used in HOLE. Please see the [MDAnalysis.analysis.hole](https://www.mdanalysis.org/docs/documentation_pages/analysis/hole.html) for more information." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Using HOLE with a trajectory\n", "\n", "One of the limitations of the ``hole`` program is that it can only accept PDB files. In order to use other formats with ``hole``, or to run ``hole`` on trajectories, we can use the ``hole.HOLEtraj`` class with an ``MDAnalysis.Universe``. While the example file below is a PDB, you can use any files to create your Universe." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "from MDAnalysis.tests.datafiles import MULTIPDB_HOLE\n", "\n", "u = mda.Universe(MULTIPDB_HOLE)\n", "\n", "ht = hole.HOLEtraj(u, executable='~/hole2/exe/hole',\n", " logfile='hole2.out',\n", " sphpdb='hole2.sph')\n", "ht.run()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that you do not need call `collect()` after calling `run()` with `HOLEtraj`. Again, the data is stored in `ht.profiles` as a dictionary of `numpy.recarray`s. The dictionary is indexed by frame; we can see the HOLE profile for the fourth frame below (truncated to 19.1126 angstrom from the pore axis). In this case, the `frame` field of each `recarray` is always 0." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "rec.array([(0, -23.6126, 15.98192), (0, -23.5126, 14.41088),\n", " (0, -23.4126, 12.88757), (0, -23.3126, 11.90096),\n", " (0, -23.2126, 11.27921), (0, -23.1126, 11.02236),\n", " (0, -23.0126, 10.96035), (0, -22.9126, 10.86575),\n", " (0, -22.8126, 10.77124), (0, -22.7126, 10.67909)],\n", " dtype=[('frame', '" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "ht.plot();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `plot3D()` function separates each frame onto its own axis in a 3D plot." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAV0AAADnCAYAAAC9roUQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOy9ebgcd3Xn/aml9+6736vlXm1X+2bLlmXLLI4HR8A4iXkdZgKYJGYMgWRCXh4CDA6DnXhCbJMMAeYlbyABgpiBGCdDBkIc28HBhMVosS3LkrGku+9r73vXMn90V93q7uq+3XeRdOX+Po8fRN/u6qrqqm+d3/ec8z2Crus00EADDTRweSBe6R1ooIEGGngtoUG6DTTQQAOXEQ3SbaCBBhq4jGiQbgMNNNDAZUSDdBtooIEGLiPkRf7eKG1ooIEGGqgfQqU/NCLdBhpooIHLiAbpNtBAAw1cRjRIt4EGGmjgMqJBug000EADlxEN0m2ggQYauIxokG4DDTTQwGVEg3QbaKCBBi4jGqTbQAMNNHAZ0SDdBhpooIHLiAbpNtBAAw1cRjRIt4EGGmjgMqJBug000EADlxEN0m2ggQYauIxYzGWsgQYqQtd1NE0jk8mgKAqyLCOKIpIkIYoioigiCBXNlhpo4DUJYZHBlA1rxwbKoOs6qqqiKErRv42/WYnWIGHjvwYZN/AaQcULvEG6DdSMUrIVBAFBEFAUBUVREEWx7P3W/xpk3MBrCBUv5Ia80MCi0HUdRVEYGxujqakJn89XRrB2MEjZbnsAiqKQy+WK/tYg4waudTRIt4GKMMjWkA4ikQhutxu/37+s7RoEWkqkpWRsRMejo6Ns3rwZSZJM3dgg5wYZN7DW0CDdBsqgaVqRTmtErKIoomla2ftXivgqkfHMzAxbtmwpkzZ0Xa8aGTcIuYGrEQ3SbcCEpmkoioKqqkC5PGAQ3eWGsR+VImPrA8L6GVEUkWW5QcYNXFVokO5rHEaSK5fLmVFsJWISRfGKkG4lLCZTGOVspZ8xomOrVNEg4wYuFxqk+xqFUWOrKMqiZGtAEARbeeFqQyUyhoXjVlWVbDZb9DerTGFExw0ybmCl0SDd1xhKybbS0t0OV0peWElUq6gwyHh2dhZFUdiwYQOArWbcqKhoYKlokO5rBJVqbOshjqtNXlhJWM+FIblIklRExo1a4wZWAg3SvcZhkG0wGCQWi9Hd3V1Tja0dKskLBkldS7BKFJUiY0MLz2azDTJuoGY0SPcahVFja0RoiqIQi8WWdeOXyguapjE5OcnQ0JBZvuX1evF6vfj9fnw+H263e82RTS0PkEbjRwNLRYN0rzGUNjQY2fqVkAaMOl1N0xgbG2N0dJTOzk4OHz5svieZTJJIJIhGo0xOTpJKpUwy9vl85n9XMxmXygj1oNbGDyuM30YURVwuV6Px4xpHg3SvEVRqaDBQqbGhHui6ztzcHKOjo6xfv56bb74Zh8OBpmlks1lEUcTv95d1rGmaRjKZJB6Pr2kyXg4WI+NgMEgoFGLbtm2Nxo9rHA3SXeNYrKHBwHJIV1EURkZGGBkZIRAIcMsttyDLtV86lchYVVVSqRTxeJxIJMLExATpdLqMjDVNW1b0WS8MwrscsB6TUTts7AM0Gj+uRTRIdw2inoYGA0sh3Ww2y/DwMDMzM2zatIldu3aRTqfrItxqkCSpIhkbMkUkEiGTyXDq1KlrOjIufag0Gj+uXTRIdw1hKQ0NBuoh3Uwmw9DQEPPz82zevJlbb70VURSZnp6+LM0RkiQRCAQIBAJAful95MiRMjKuFBkvl4wvZ1RtwKiZXgyNxo+1jwbprgEsp6HBQC2km0qlGBwcJBwOs3XrVnbu3Fm0zL7SzRGlZGxgNcj4cpPRSkgatTR+WM2CBEEgmUwSCARwOp2Nxo/LhAbpXsUwamxnZ2cBaG1tXXJ0Uo10k8kkAwMDxONxtm3bxt69e22/40qTbiWsNBlfiWPUNG3VdORqZNzf38/u3btRFKVRa3yZ0CDdqxCl3WPxeBxBEGhra1vyNu1INx6P09/fTzqdpre3l/379y9ZF74ab8alknEul8Pn85FKpS6bZnwlJA2j2cXhcBTp9I3Gj9VFg3SvIpQ2NFgTI6W1nfXCSpiRSISBgQEURWH79u1mBL0YrtZIt14sRsbDw8Mkk0kuXbpkkrE1Kvb5fLhcrhUlGk3TVixBWe/3lkbYjcaP1UWDdK8CVGpoMCBJUlmmul4IQn6W2fPPPw/A9u3baWlpqWsb1RosrkSkttIwyNjv9+P1eunq6gKKI+NwOMz4+PiKk/GVOn/1aMlLbfyoVN72WkWDdK8gFmtoMCCKolmHWy90XScYDDIwMEA2m+XQoUM0NTUtaVvVrB2vpZuolAAXi4zj8TihUGhZZLyamu5qYzEytpv4EYvFaGpqwu12lyXwrqVryQ4N0r0CqLWhwYAkSXWXaum6zuzsLIODg3g8Hvbu3cvZs2eXTLjGfpZGukYSLhwO4/F4TM8Fn8+Hx+NZs0RSC6qRcSKRIJFIEAqFGBsbI5PJVCXja2GlUIpqZDw2NkZvb6/tZ671xo8G6V4mLKWhwUA9ka6u60xPTzM4OEhTUxMHDx7E6/Uua9+t+2GQbjKZpL+/n0QiwbZt29i+fTuZTIZEIkE8HmdmZobZ7BgxeZZe6Qb8vmIytj1uPYc7/XvkHHehyneuyD4vBcslQEmSaGpqKnvAVSNjVVWLiHmlNeOrCYIgoKoqDocDSZLM118rjR8N0l1lLKehwUAtka6maUxNTTE0NERrays33HADbrd7WfteCkEQyGaznD17llQqRW9vLx0dHUC+e81wGOvs7ATgufAQJ0Lf5xc6f4V0IkMsFmNqaqrId8EaGbucMpJ2ClW7kaWJKVc3qpHxhQsXkGW5iIwlSSorbbtWyFhRlCLChZVp/DCCmuWs6FYbDdJdJRjJsWg0itfrXVJDg4Fqka6maYyPjzMyMkJHRweHDx/G5XItd/fLEI/HuXjxIpFIhOuvv5729vZF61olIX95uTwOmv3FSbvS5JRBNLfvE5mfmyKijuLz+fD7/Ze9YuJyL/WN6K2jo6MouVkpMl4pMr6SlSj1NoPU0vgB8NRTT/HSSy/x8MMPr9i+rjQapLvCsNbYKorCuXPnuOWWW5Z1E9vVx6qqyujoKOPj46xbt850/FppxGIx+vv7yWazbNq0CV3Xzeh2MchCfn8UvbzcrZIeSsJBIOAmHpOYn59nZGSEZDLJCy+8YBKMER2vxvEauBIdaaXfWSkyVhTFfGBVI2O/34/T6ax4LFc6ebcS57iUjCORCM3Nzcve7mqiQborBLtxOLIs19xTXw2SJJlPcsPxa3Jykg0bNtTt+FVrFBeLxejr6zNredva2shkMoyPj9e+3wXSVXVlkXda4cDlFNi4caP5yqlTpzh48KAZ9c3MzBCPx1EUBYfDUSRR+Hy+Zde7Xu0dabIsL4mMrefJ6XSiaVrZEn+tIxKJ1F0KebnRIN1lolJDw0pCkiQURaGvr4/p6Wl6eno4evRo3TeMETFX+1w0GqW/vx9VVc3GCQP1NkdUi3QrwwmUv9/hcNDS0lJ2Q2WzWZOMJycnSSQSqKqKy+Uqivi8Xm/N5+tKVBKsxHfWQsbBYJDR0VEymYyp0Y+OjhaR8VrWjCORCFu2bLnSu1EVDdJdIhZraFgpZDIZBgcHiUajdHd3m45fS0E10o1EIvT19QGVGyfqnT5hkK5qQ6KVoAsO7Ei3EpxOJ06ns+jhoOs62WyWeDxOIpFgbGyMRCKBpmm43e4iicLr9V4VZW2rudSvRMbGakaSpCIyliSprLRtJcl4JVZ/ldCIdK9B1NrQsFyk02kGBwcJhUJs3boVn8/Hpk2blrVNO9IMh8P09/cjCAI7duyoqodVa46wg5FIqy/SlRHqen85BEHA5XLhcrlob283X9d1nXQ6bZa1zc3NkUwmAfB4PCbBlGbHLweuVJ2u2+0uknJgITKOx+Omrp7NZleMjBVFWbWW52g02iDdawX1NjQsFUazQSwWY9u2bezZswdBEBgaGlr2tq2kGQqF6O/vR5Ikdu7cWVOJzZLlBW358sJKQBAEPB4PHo+nKBmoaRrpdNqMjKPRKLFYrCgxtdqDNq9EUqvSd1aTKQwpZzlkrKrqqmnJ0Wi0kUhby1hOQ4MBg+gWu6Hi8TgDAwNm/etijl9LgSiKBINBxsfHcTgc7N69u7x6oArqJV0zkUbtibR65YWVgHWKMeR14vXr1xMIBCoO2rQSzGJVArXgSkS69ZKfLMs0NzeXkVq9ZLyapBuJRIqkpqsJgiAIuq7rDdK1wUo0NBgwKg8qka6RuFIUhd7eXtra2lb85jP8F0KhEJqmsXfv3rIRObWg3v1aWiLt8pOuHQyNvto4IcNzwVolUFrW5nQ6a/q+qynSrQRFC5NWR/DIO5EEj/l6vWQM+WtyNRJ4V6umKwiCqOu6Bo1ItwhG2VcikWB6eprNmzcvW0YwKg9Ka0oNLRWgt7d3VZ7Ouq4zPz9Pf38/Ho+HlpYWduzYsSTCrQRFURgaGmJqasrURQ3CEQvRTL0lY8vVdJeLxaL5SjXGVpKZm5tjaGiIXC6HLMtlZW2l18OVGhFUD+lGc8/TH/84+5u/iVfetej7K5HxzMwMc3NzSJK04pqx0V58NcEgXEEQ/j3wxgbpUl5jq6oqwWBwRUpPrDW2VscvWZZr1lKt+1nLBWiMSh8YGMDr9XLgwAF8Ph/nz59fsRlniqIwPDzM1NQUPT093HjjjeRyOVMXnZmZYT4zCW0wPDqIx9NlEnL17ikHcPkTWVYslQArkUwulzPJeHp6mkQigaIoOJ1Ok1wMCety1s3Wu8xX9TgAkrC8h7au63i9XtsEnl1kLMuyWWliXEMOh6Oiq9lVCGPHVKDvNU26dg0NhsNR6djrpcLY1uzsLAMDA6bjV73RZi1OVIaz2MDAAH6/v8zsZjlj2K3fMTAwwOTkpFkvLIqieXN4PB7TeyGqbOC50f9Fa3sLLtVVNKXBiGqMCNC4kXTBiaAnl7WPVxvsaoyNXIHxkMrlcpw9e7aoxtha1rYaZFxvpLtSpKuqqm31Qr0yhUHGxkNLkqSr0q1NLzwNdF1/Gl6j8sJiNbayLC/Zv7b0ezKZDOfOnaOlpWVZjl/VtGFd15mZmWFgYICmpiauv/56PB5P2fuWQ7pGJ1wikUAUxaLmjEoXupFIk1wi65vWl23PuJFmZ2fNpfh1WxJ4XXEmQhMm6VxuXI4bVxAEnE4nbW1ttLW1MTU1xeHDh81rxjg3wWCQZDKJpmlFZW0rUWOsqmrNmjNYSXd5v0kl0q2EaisII9F54sQJ/uqv/orR0VF+4Rd+gX379nHPPffwxje+0Xabo6Oj/OZv/ibT09MIgsD73/9+PvShDxEMBnnHO97B0NAQW7du5fHHH7eV/o4fP86nPvUpAD75yU9y7733Vtx/QRBk4G7gKV3Xo68p0rUj20qm4cuJCK2OXwDbtm2ju7u76mceO3GWRCbHe287bPt3O6cxq41jc3Mzhw4dsiVbA0s5LoNsJyYm6OnpwefzsXXr1po+Wy2RVulGciT/HkGdR1VVxsfHTfI5e/ZsUVS82l69VypaEgQBt9uN2+22rTE2ImO7GmPj/FS0zizBUiJdES+CsLyoW1GUFXHAczgc5jX0zne+k9tuu40PfehDPPbYY7zyyitVZwrKssxnPvMZbrzxRmKxGIcPH+bYsWN87Wtf44477uD+++/n0Ucf5dFHH+XTn/500WeDwSAPPfQQp0+fRhAEDh8+zF133VUtL/M/gGeBPxAEIfmaIN16GxqWesNpmsbExATDw8Om49fY2FhNS8OTA+OMh6MVSdfqNKbrOlNTUwwODtZl41gP6SqKYhrq9PT0cOuttyJJEhMTEzV9HiwdaXUk0kTRhYhW1Ahy8uRJdu7cadvUsBp1tFejNmitMTbkG8hfc6lUyjw309PTpFIpBEFYdPz8UjRdSVx+Ena1SsaMyoWWlhZe97rXVX3vhg0b2LBhAwCBQIC9e/cyPj7Od77zHZ599lkA7r33Xm6//fYy0n3qqac4duyYSerHjh3jySef5F3velelr+sCbgLuA9quadK9XA0NqqoyNjbG2NgYXV1dHDlyxFy2WRNpBnRd58MffYw7/t1efuWXDwEQcDuJpSrPQTO2MzExwdDQEG1tbXXbONZCuqqqMjIywvj4uNl2XMsNYrckF8l/zi7SjedeRkcl4DhUvB2bOt1qTQ12U32XU7pV6ViuVlhrhq0wzk08Hi+beGy8P5lMEggEaj5eVY8vW8+F+uWFWrFUh7GhoSFefPFFbrnlFqanp00yXr9+PdPT02XvHx8fLwoKenp6FjOC+ihwHfBDYNM1R7or0dBg3Va1zxnR4MTEREXHLzvSFQSBkZEgk5MR87WAx0U0bU+6RjRz5swZurq6luyZW62N12oVuXHjxprJdrHvkwUHqg3pjiX/Ap0ce5u/UvKX2kvGrHW069atM1+vpBdbHcmM/73WXLYMVKsxNs5NKpViZGSEoaEhJEkqK2srfVCpemJFSNfOwHwlEA6H667RjcfjvP3tb+dzn/tcWSXRSgVpuq4PAUPG/79mSHclGxqgujlMLpdjeHiY6elpuru7qzp+VZrk6/U6SSQXXm9yu0hmcyiqhizldTarXCEIArt37zYn1C7nmKwwyHZsbIzu7u66rSKtsPWEFRy2ka4oOFE0uyqF5TdHVNKLDUeyeDzOxMSE6UhWaoKzmoYsVxpWj95QKMTmzZvx+/1VH1QGCWfcEZzy8uvJV1teqBW5XI63v/3tvPvd7+ZXf/VXAVi3bp1pmzo5OWl7v3V3d5sSBMDY2Bi33357xe8xOtEEQXACyponXaPsS1VVzpw5w3XXXbciTyij1Mt6cWSzWYaGhpidnWXz5s01OX7JskwikSh73edzkkhYSNeTj1xj6QzNHhdjY2OMjo6acsXg4OCyE0dWXdhKths3buTo0aM1k20puRrn204LrRTpCjjQbMhVF1bPe6GSI1mpCU44HC5L3K2m78KVgjWRVq1CwCxrUyNkEwFOXjhp1hhbo+NaiXS1DG/qMbvRdZ33vve97N27l9///d83X7/rrrs4fvw4999/P8ePH+dtb3tb2Wff8pa38IlPfIJQKATA008/zSOPPFLtu4ySsSys4UjXrsY2k8mgadqK/KBWWcDq+LVlyxZ27NhRMwHayQsAXo+LZGqhCcDvzpPuhYFB1GiobBrEStTYiqJIJpNheHiYsbExNmzYUBfZQm31wlZIyCg2iTRBcFK4BksgcznbgO304rNnz7Jjxw50XScejxONRov0YqshuOG7sBxcqcRdLRGnw+GgtbWV1tZWZoNZWts2sWXTkSIf4/HxcZLJZNGqwfpf6b2ympGuoccuhp/85Cf8z//5Pzl48CCHDuXzCg8//DD3338/v/Zrv8ZXvvIVtmzZwuOPPw7A6dOn+eIXv8iXv/xl2traeOCBBzhy5AgADz74YNVKCQBBEDqA3wO+seZIt1JDA+QvEKPtcrmQZZl4PG562Vodv+pBJdL1+ZzE4/lIV1VV0tEwAOFkil+sURuuB6qqMj8/z+zsLFu2bFmyjGCQf81TDkT7SFesEOmCAwENdBWWWZq0VOi6jiRJZpOCFYYmatgeDg8PFy3DrZHfUlcOlwtLKRmTBH+RdaaVbKw1xvF43Kwx1nW9SMKxGv6vJOqRF97whjdUfNg988wzZa/ddNNNfPnLXzb//3333cd9991Xz+4pQAa4Z82Rrqqq5HI5W9Pwleoki8fjhMNhIpEIu3btYt++fUu+QCqTrovp6SiDg4OMj4/jdeTra5vbO21v1qVGupqmmVJFIBCgu7ub7du3138gBZTKCIb3wsTERFFm3EjiSNhrupUjXSNqzAILNcdXS0VBpbllVr3YOsHCGvkZEyxKr9srdWz1RJyankUnWzWRVq3G2ChrSyQSZLNZnn/+eYCysrZaa4ztcLWa3QDouh4GHhYE4c41R7rGvHs7GJHuUmEdwtjU1MT69euXlbgCe9JVFAVVyRCO5Dt8jh49ykgwCpyoWMFQb6RrJdv169dzyy23EA6HmZ+fX/KxwALpGk0TP548w0nXEH9w6N20OZrMmysSiTA+Pk7akyGoz3MhfKGIjCtFurpgmJWsTBv2UrAUEqykF2cyGVMTnZ+fL2poMKLiKzVWvZ7jXE4LsFEz7PV66ezsZHZ2liNHjhTVGMdiMaampkin00U1xvWco6uZdA3ouv7EmiPdalhqpBsOhxkYGEDXddPxa2BgoCyy/Np/f5Ke3k5+8VftGxjsUDpUcnh4mMnJSTweB7mczrZt2wAIFDTdSrW6kiTVNNHASraroQsLgmBWbvT09LBp5zYe7z9JRldsS7henvgeAgLrvOtMw5f+/n7UliBCU5r+/v6iKDBfvQBX2t5xJUjQGvmV1henUini8TixWIzx8XFisRinT58u86NYrl5cyz7WAlXPJ4OX2xxhXSVZV0bW4MawzkwkEoTD4fzD20ZPL3Uhi0ajV62XrgFBEMTXNOkGg0H6+/uRZZnt27cXZW7ttnXiX3/O/Ey0btJVFIX+/n7TkevWW29lePQUmUwfiqIiyxKBQvVCpUjXWnlgB03TGB8fZ2RkpOJI9uWQrqZpjI6OEgqFaGpqMhNwk8H8zZit0HUmCw5yeqbM8GU08RxTaYXm5mZTH00mk6xvGWfPRhgfG8Tt1fH5fGbt9eWKBlc7sVXa0JDJZLhw4QL79++vaOxSWkmxWuNuKmElzW4WkzQqWWdaa4yt5yccDvPYY48RCoU4d+4cDoejSN6ww3333cf3vvc9urq6OHfuHADveMc7uHDhArBQ83vmzJmyz27dupVAIIAkSciyzOnTp2s6dkEQREBYc6Rb7carRV6w2h663W727NljOz1BluWy+tpAi5doqHYHrFwux9DQENFo1Gw4MKQRnzcfwaRSOQIBCbdDxilJxNL20ayd9wIUk21XV5ct2RpYCulat79+/Xra29vp6ekxb3pn4RLKafakKwkOUlq8fF8EB6DS3t5WFAWK2WnIgcspmgbhyWSSF1980YyiV5t4LrfGatQFV9KLrWVbpROPrWS8mkM2V9LsZqmVC5XOTzweR5IkPvnJT/Ld736XRx55hOuvv54///M/r7it97znPXzwgx/kN3/zN83XvvWtb5n//shHPlK1u+0HP/hB0XVbC65JE3M7ojRgOHENDg7i9/tNj9lq2yqtr21q8RKcjS26H6X1vD6fj82bNxe9x+fLR7aJRIZAIO+b0ORxEa0gL5RGutbGicXI1rqNWklX0zQmJycZGhoq2v5LL71UtA2HaAyftI/CK9bpCvmHjk4WgQXfCEHMv97Z2UzHul4gX65z4MABc0lul6gyyHi1jXBWA4tVEVjLtgyU6sXBYNC8XpdqgFMNqrYyka6o/BsHe76KoP0PdLFz8Q/UAL/fz5ve9CZkWebzn/98Tcd62223VZw7qOs6jz/+OP/6r/+6IvsHIAhCALgH2LXmSLfeSFfXdSYnJxkeHqa5ubmi7WEpDFnAikCrj6GL5b3YBrLZLIODg8zNzbFlyxYzsh0dHS17r7cQ6SaT1lpdJ7Eq8oKmaWVka/V5WAy1jFA3ztfQ0BDt7e1l2y+tXnCKeaKvLC/I9h1pBe1W0xXEop80/7qg57DuqSzLtp601saG2dlZUqkUgK02WivxXO5Idynft5hebE1Olc51M6LiemSUlZIXBG2MgPsiCWF1pjusxO/2ox/9iHXr1rFz586K3/HmN78ZQRD4wAc+wPvf//5q+yPpuq4CHwZuBYJrjnSh8oBEqw5rkNPIyAhtbW01O3FZt1WqoTa1eImGy+WFTCbD4OCgOW1i586di0ZbXm8h0i1pBa4W6UajUZ577jk6OzvrIlvrNipFulabyJaWlor+DqXE7RDzS8VchYm/ldqArZFuMcoTaZVupGpGOAYRh0IhRkdHi7TRpdTSriZWcj5ateSUVQ8dHh4mmUzywgsv4Pf7mRUV/tvwczx44E28ft22su0qBdKVl+u9oIXRBQGofSBqLUin00vyI7HD3/7t31ZzDOPHP/4x3d3dzMzMcOzYMfbs2cNtt9222GbfAPyJrus/vPJX3ApClmVyuRwjIyNmC+1NN920pAywXSKtqdVLNp0jk8ri8jiLOtW2bt3K7t27a37SGpFuaSvwfDxV9D5jmT8wMADALbfcsuSMtp3hjTFtor+/n38KzdG7YQP37d1b8zacwmKRrsPW2lEs1ONqJbW6umAc29KrF0RRtE3EGNqonURhJePL7b1wOSLrUj00l8tx7tw5Dhw4QDwep296iOlMgoH+ARzDs2XTKxQxL6uJy4109Siq7lvxxpelOoyVQlEUvv3tb5t1xHYwvLG7urq4++67OXnyZC2k2wpsFASh9ZohXUVRmJycZG5ujqamppo0zmqwJd2W/NSH2ekQ0eQ84XC4pk41g6is0Yyh6VrlhYDbxdBcvjPNINvh4WHa29u5/vrr6evrW1YJkTXSNYZW9vX14ff7OXToEH/2xHeZn5qiWp9N6SrD0HSrJdLsI938b6OXkWte+hH0FCuNStqoYQ4ej8eZmZkhkUjw/PPPF2nFq1m+daUmAUuSZJ4TTzYIwOH9B+l053g5+mNaOUImmcnrxe4B9ICTV86/uiy9WCCCpq/cYFQDK1Wj+/3vf589e/bQ09Nj+/dEIoGmaQQCARKJBE8//TQPPvhgtU0aEcrPgLcAr1uTpGu98Y3Idmpqio0bN+Lz+ZbVcWXArhnB5cufrlM/e5Fbbz/E3r17a7rg7Ebt+Gw03UBBXjA8c9vb281IPZvNroj3gqZpBINB+vr6cLvdHDx40EwoemSZ9CIld+Wabv6cZCuQ7mKRbmlXmm5mxxcqHirJSSsBO3PwU6dOccMNN5iJu9J235W2h7wSHWml12M8l/8d/LKTudwgJ+NPcLD79XR25slnMO4jkm1i27ZttnqxnR+F7VQWoqh67cNYa0U0Gq0r0n3Xu97Fs88+y9zcHD09PTz00EO8973v5bHHHiuTFkbJHGoAACAASURBVCYmJnjf+97HE088wfT0NHfffTeQD/Tuuece3vrWt1b8Hn3hwn0Y2A8cWJOkC8UVAps2bTLtFScnJ1dk+1bSTSaTDAwMMD2bNyreumk769evr/Zx221ZI29T0y3IC7quIyhZYukMkWi0TBZZrE63FkSjUWKxGCMjI+zbt6/Ma9UtyyQWacAo1YUN0k1mU2QymbLOIUmQ0VDRdBXRsqQ0It3yrrQ86Qp6uTPb5YQsy7YSRTabNSsGjHFCxvyy0nFCtRLplYp0i0hXWSDdYCa/ynCJC/P8jKkR1ZoZ7DT0Uj8KF1F0Vr5rLBwO10W6f/u3f2v7+te+9rWy1zZu3MgTTzwBQG9vLy+99FLd+6fr+iQwCXx/TZLu6Ogow8PDNdsrLgWGJPDyyy+TSCTo7e2l+fVdfIOf2SbTqsEuanY6JWRZJJFciGxlXUMHerZuK1vKLqexIRKJ0NfXZ2a8DVelUnhkB/PJ6sdW5r2QzpPm+NQEr86+WpawUqX8Piu6gtNCuotHuleWdCvBOkzSgNVbwDoyx0hqLdZhdqUiXWuEHsvlH/4Bh4tMKn8NuMSFKp9qUyMqNTMYo+eNc5JIJDi8NUws1c7UbF9RJcVyVwtroQVYyP/I0pok3a6uLjZu3LjsC/WLHz6OklP44BfeW/R6IpGgv7+fVCrF7t276ejoQBAEoo48EcTqaJCAyr4JbrfMwMAI0Wgnhw8fZvaVATjTRyyVMduCDdRS7lWKWCzGpUuX0HWdHTt20NzczE9/+tOK7/fIMqka5YV0Ot/CG41GAejauI7r11+PIAhFCatYPA4OOPX8CZo9bSYZ6+48GWtl1Qv56Ko00r0a55YZKPUWMFDNkcyqFVea8ryaKI90M7hFGYcoEVfiDCfaGEjMsDewJX8sehxJqK/iwG70vDOWwi2up9XdSiKRYHR0tGja8VJXC0uZGnG5UZAa1qaJucfjqdjuayzDa3lyzozOMTOyYAATj8fp7+8nnU6zfft2YrFY0U3ka8pfBMuJdK1DJV1OCZ+vmT179gAL/gvRdIaNdX1DMeLxOH19feRyOXbs2FFzP7pblkkp1asGjHbgbDZLb28v+/btQz7xPXLawkPFmrCai27i3DwcPHQAKesikUgQjUaJzk1DJ1y8+AoBOVBUPaDjKyLdq8FdbCmo5khmPJTGxsYIh8NomkYsFitK3q2maXpZpKtk8TvyUXhYjTOdCTCRmV8gXS2OU67Nq7YidAVZTCLSRntLu60TmSHdWFcLtejF0WiUHTt2LG//VhlCfozyujVJutVgN/GhElq6mrl4eqDIXWz79u20tbUhCAKXLl0qiggkScTf7CG2BNLN5XJMTU0xMDBAa2srhw8fpqV1kHR64eHRZJBulQGV1WBE6Ol0mh07dixqrFyKaok0o6V5YmKCdevWceONN5oXvkOUq5aMAWioNFn0wISS5pUIbOvtQcqsKyrjunm7g2h0lOnUQJH/6rWCUoliYmICVVVpb283TXCsSSor4fj9/mVV5RgojXRjuQwBOX/9xZX89e2TFuraV2YoZaGbUyzXXq2rBStK9eKxsTEymUzRTLfBwUHm5+e56aablrl/q46twG9fc6RrdKXVUijtbXYTmY3y81d+zo6d5SRlELhVh2tqrc9/wXiCT09P09nZyY033mg2aXg9zqLqBevInnqQSqXo7+8nHo+zY8cO2tvblxQheWRHmbygqirDw8N8b+gcPV3rObRlS1kE5hRkcppiq01KgtEmXDLht6DpCqJWvgRNttLkkEjKAeLxOKlUihdffLEoMWPccCtBQFcauq6bEZ3X67VtajA67gYHB81rsrSKoi5D8pJIN65kzEg3rhbG20sL99BKkK6g5wexClLtTmCV9GJFUcyo+Jvf/CYnTpzgiSee4Atf+AI33XQTf/zHf1xxm3ZmN3/0R3/EX//1X5sr24cffpg777yz7LNPPvkkH/rQh1BVlfe9733cf//9ix93YUYa0AvsXpOkW41QanEai0Qi+YhQT6FrOrt799DcVl7GYku6FbrSSmHt8BJFkc2bN7N169ai9/h8TmZmFrwcTHvHGknXqqtu376d/fv3L2s56nbIKJpGTlWRBKFoYOUPhTBdMYXDgXVlCT2nKJNbJNIt9V8QK9bp5pNpspims7OTzs5OIpEI+/btQxTFosRMPB639V/wer1rSpKoNmLKTqLQdb2oisLQRXVdN3XRxea6lWm6uSwBR6F2XMlffz4pn0jTdAWNzLK70UzSFZdvv2htC//Lv/xL7rvvPh588EE6OjoYHBys+lk7sxuAD3/4w3z0ox+t+DlVVfnd3/1d/uVf/oWenh6OHDnCXXfdxb59+xbbXQHQyTv062uSdKvB6Eqzg5HFB9i+fTvxw1me5IeEZ6M0d5aTrl0CLNDqZWYsVPH7DWOdgYEBmpubueGGG5idnS1aHp8fnOILj/+ELuSiNuBADfKCccMNDAwQDodNXXUlSMZTuPEHRkaYm5gwzc9lWcYz+hxpVbGtmXWIjop1uotFuqUdaXn4bEvG7IYnWo1fjGgwmUwWVQ6s5LJ8NVBv9YJ1XI5VF7X69JbOdStdIZQSfUzJsNGbvweSWhYQ8Mn5FZm2Ur4LJum2sNJiUSQSoa2tjQ0bNiw6J62a2U01nDx5kh07dtDbmzdieuc738l3vvOdRUnXcBcDzgEvrUnSXcz0pjTSDYVCDAwMIAiCmcUHaOnKX2ThmQhb9pV3oNh3pfnoOzde9l47sjVkhNIx7Kqqc25gijdt6SlqA/a5HEiCUNFTF+DChQsEg8Elz2wz9rX0c7quk03m6zODsWiZt4NHchDPZWxL15yiXLEjTS609ZaSbvVI14+ozdV0LJWMX6zL8rm5OXOkuLEst0bFV9qVbKXqdEt9eg0Y49WNjrt4PE4ymcTpdJJIJPD7/cSyGXxS/jdJq1nAZWq6ygqTrp2mu1xEIpFlG5h/4Qtf4Otf/zo33XQTn/nMZ8q2Nz4+zqZNm8z/39PTw4kTJ2radkFiGAYeXJOkWw1WogyFQvT19SHLMjt37izLIrd05X/88Ex00W0ZMDRdg7is3gVNTU0cOnSozMWsNGJu9hcSFKJAMpk1tyUIAgG3q2x6hJHESiQS9PT0sGvXriXfpAZpGnqe4S/c19eHks6T7sYtW8rqSV2izJyasPVvcAiVE2lGpFvalVYt0s3X6i6vI22xZbnVOB0oigSX2/lXL1a7TtduhdDf34/X68XpdOaTd0qGxHyIU6dOkfTmSTcZiuMIiKhSgXSXOTUC8qSrCytf2pXJZGpyD6yE3/md3+GBBx5AEAQeeOABPvKRj/DVr351xfZP13VdEITtQMc1R7oOh4P5+XmmpqZwOBwVTcqhONK1QyX/BSWnkkpkiCejVcnWQCXS1QBN08lkFNzufJQR8LhMI3Nj6OP09DSbN2+mtbWVrq6uZUVFVtINBoNcunQJr9fLoUOHiE2Ow8Wf21YweGSZlJqrIC9Ui3Tzx1UmL1R0GQMIIOiL+xbXi2rL8mQySSwWIxgMkk6nOXnyZFk97UoU8dvhSnWkud1uWltbCbQ0kzuvsXPTVm7cdiN/8/Nv4BAE4tEYUxOTZISXEbphYixI0jVuPqDqdWgT9BCaLgEr672wEpUtxngpgN/6rd/il3/5l8ve093dXWTTauQ7FoMgCGJBYrgJeOuaJF27qEDXdYLBoDnr7NChQ2VtrqXwt/oQJZHIrH2ka6vpFkxv/u0HP6VnW1dVsq20nYDXhSCAWlC2EonMAum6nURSaQYGBpiYmGDTpk1m193c3NyK+C+Ew2GGh4eRZZn9+/eb58kj5/chaaOJuyUHaVWxbdJwijIp1b59uGIireCfr9mY4ehCAIEk6DlYJd/Von2xzHaD/ArpyJEjRVHx6OioaRJuHZro9/uXPVjySnSkWYneaAEOyE4UIYuii7glh+lhEsqG6ItBe+sm1CQLc+4s0ysMrbiqXKOFyCl+WKVjXc45nJycNLXgf/iHf+DAgQNl7zly5AiXLl1icHCQ7u5uHnvsMb75zW8uum2LpvsiEFqTpAsLS07r+B2Px0Nvby+hUGhRwoX8zdbS2VQ10jWScsb3TM7kn3TdGzZz4EBtxdilEbMkijR53eQKBJpMZmlvz+uQDnQmZucRRZFbb721KLKqdyJwKWKxGLFYjMHBQfbs2VMmt7gLkYtdpOuWZDORVi4vSMyn8/4UgUCgqIKgYiJNkAHRNtLVBWO/4kDrqhreVINdy681WWVMPM5kMkWtz/Ua4Sw10lU0DakgS9ULa8lYzPBdcDjJaElUXcQjLjzsDAPzlqYNuFsXIrvSJObc3Jwp1xgNDVaHNkELktNW1kcXMLv8aoWd2c2zzz7LmTNnEASBrVu38qUvfQkoNruRZZkvfOELvOUtb0FVVe677z72799f03cWot2LwMU1S7pWLdXn85njd4xkQa1o7moiXCHSlWWZZDJpfo/f7+fgoX38PafJpmqPOO3IstnvJlt4LRZPMzIywsjICD6nTFB2lJWXwdL9F5LJJH19faTTaXw+H/v377cdVeRxVCZdj+QgXSIvGEbxsVCUnKzQ1NRkni+jgkDy6+CCWDKK4lGKlqQijoqRLoCgx9CFq2u6qzVZZV2SWluf7YxwqnWZLTXS/cS3n2E6Guf4fXfX/Vkr0Ru+C37ZRUZLoegiTbKlRrfCqJ5q0yvsGhpu3DaOonqZHx9f0Tl3kUikLICoBjuzm/e+97027yw2uwG48847bet3q0EQBIeu6zlBED4ANK1Z0n355ZcRBIHrrruuqIul3onALV3NhKfLSVfXdRKJBGNjY6TTafN7Jgttw/U0SNiTrodUPB9hnD59htfduoujR4/y48hPOTd1qebtVINRxxuLxcymiZdffrkicRvygp3/gkuSyWoqGjqappnddR0dHaxr7yKcGKWjo6MoYlNVldnoFIQhHA3y0thLRabhNDnI5ZJlpKML+chS0GfR2Vy2L1cjKnn1GlHx0NA0U1Ov0tXlLIuKFUVZUqSbyGTN1cl0Is59//yP/P6Ro9yxpXzyQymsydS4YjG70RKouoCvpDECah9KWSrXGHDF04TiPeZIqFITeavnQj3n42o2uylULVgji841S7oHDx60fb2WicBWtHQ1M35xwQ7SMPc2RrO3tbUVfZeh6dbTClxKlrqu45J1pgpmMZs3bzP7xtt8HqLpDDlVxVGyPK010jXqeIPBINu3by+q4622DeMGtvNfcBe8kYLRKNMTEyiKYnbXufrPkrUZ1yNJEm3N7RCG9s42Du88XGQaHtQlwtE5Tl08ZbZ1+v1+mv1NeGQQtXE06bB5ztYKBl+dxOGU6entNLvMnv7+AP/72y/zD3//n81uqng8ztTUFMFgkFgsVrYkX8zwJZnN0eLNJ2Wj2SxjsRhKjSshq8mO1Us3o82h6qLZGAGg6gkEXGaZ31IhEkEXDhYZhFeac2e0Bdfi0btSUyNWA4Wqhc8BHwE6WMuabqWor94leEtnE5HZqGnu3d/fj8fj4cCBA+i6Tn9/f9H7fQE3oiQSDdVuPWjsq9Xsxu0QUcmTaiazcBytvvzFHk6m6Qz4bLdTCUZp2czMDNu2bbMdH1Tt/BiRrq3/Qjb/2sTsNF3t7UValkNw1Fy9YDUNnwj58Lf62L755iIiGpvU6OoRmRw/wVRsG8lkkmAwSHt7+7KTVpcDf/3wE/gCHv7rX9xjvpbJKLhc+dutdMhmLpdj27ZtiKJongOr4YtVJ7Y2eSQyObpbCg0NhUDDV6O2WSQvKAvyQlhLougifnlh9ajoMeQ6HcbKoGcRhQSaXhyRLmfOnfFgWgMOYyfJd6R1ARvWLOlWQr03ZEtXE5lUlh//8Ce0tDcXjWZPp9NlUoUgCASaPXU7jeVyOX72s5+ZQx/PTr3AifOTuCieHtFWIN1QImVLunaEafgjTE5OFlU72KFaUmoh0l04ZsOxLBrJj3Pp3roZPRQv+pyziuGNKEgIiCg29bhOcR1ZbQooJaIe9OQGutdreDt7eeWVV8xEzXKTVpcDmXSOts5S8/MF0i2F4b1QOsECKnsvuFwuoskUopZvfEhk88TpddQ2UqioesHipTuVTKDqIgFp4dpT9QiyuLxpD4JeGENVo4F5LXPuJiYm+OxnP8tzzz1naucHDx7kzjvvtM1ZgL3vwsc+9jH+8R//EafTyfbt2/mbv/kbWxLfunUrgUAASZKQZZnTp0/XdCy6rhslDh8CuOZIt1YYJWbBeF6j7e7cxPYDW4veU0kfrtX0xqh46O/vR1EUbrnlFrNLrdmfr15wAfFE2vxMi3eBdEtROj1C0zTGxsYYHR1l48aN5vSMaqhJXsjlSKfT9PX1kUgk2LlzJ3uSbjj9KlldQ7LtSKscgTsFFzkb0nVLPUSyz9l+RhN3IKs/wu9J4Ha72bRpk6nd2yWtdF23z5hfgag4m1ZwuosjzkxGwemsTrp2qNTkkclkSD/1Ag4BBgcHOTuVl8imRkYYzeaqmqYbMM6NUb3gkx3E1AQg0CQvkJ2iRZGE5S3fBT3fOr/cxGipdv6lL32Jr371q8zNzeWDmbNn+cVf/MWKpGvnu3Ds2DEeeeQRZFnm4x//OI888gif/vSnbT//gx/8oCgiXwxCvkznzeQ7Q3JAZM2Sbr2DIK0wNFuXy8Xe6/fwHf6FbKycFCot55tafVU1Xasu7PV6ue666zhz5kzRCPhmvxsEAY/HSSK+0IFmRLrzFUhX0zR0XTenTaxbt870R6gF1UhXFATcksTE7AwvvPACO3bsoLOzM99UkMlHuhlNxWPTHKGhoeoqIuXnXBad5PTy1maXuImc/o+oegpJKK51zjo/jCf1DpzZ/44gFBuT2CWtSjPmxnLU2uAQCARWte3399/wh9zx628gm87Zkm6lSLfe6cOCIOBwOkkrKhs6Ozhw4AADLieMDtLb04OMUGSabtf6bEVcyeKTHUiCSEzJN6UESuQFl7g8L11Bz69oNLoWeWf9iMVi7Ny5k7e97W287W1vq/peO9+FN7/5zea/jx49yt///d+v5O75gbvJ+1oGgK1rlnSrwc4dDDAHMrpcLnNGWL8+BEDIpla30o0QaPEyPjhr+zfrd1ililI0+/Mk4/E6iMUXIt31zfmM73QkXvYZURSZn59nbGyM9vb2Mn+EWlCJdBVFyTdMAFqhRth6/G6jL19TcZc1RxTGsGsKDqlcU3QILnKaXaSb72NPKf34HcXF6Lq4iZzjfThzf0Gz53ryrnjVj8suY25tcBgeHrZt+12JSb+qqjH08ijxcJJMOofLjnSd9nrrUup0UwWN3VvYZqKQDOtsaaXd0qxTqfU5kUhw/vx5/H4/c/EofsmFruvEDC9d2ZpIiyCLu+vav1IIet6vRBMW7+CqF+FwmOuuu25FtvXVr36Vd7zjHbZ/EwSBN7/5zQiCwAc+8AHe//7317LJFPA5YDvwO8DJNUu6tdg7GjeSkSBzOBxlAxlb11X3X7BDU6uXn79QHOmGw2HT58Fu6CMU12O2FFqBXW4HsdhCFOhxOmjxuJkIR4s+Z8gUDoejyJO3XpSSrlWi6O7uJuDxINokqzyF6oUcmm1zBFDR3tEhumwj3YDjCJLgYyr9v9jheLTs7znHbyAr/8y29q8T1d+KMcqnHlRqcLAbo2N0VxmzveopXcoW3OJcXhfZjIKzJKrNZhWcLnvpZyl1uslsIXHmKnQRKvaJtEqtzydPnmTLli3E43HCqQQOVePUqVNM+vM17plIkpicr6jIywvL03RFbQxV8yDKK193HY1GVySR9id/8ifIssy73/1u27//+Mc/pru7m5mZGY4dO8aePXu47bbbqm5T1/UMcF4QBD+Q1XX9v65Z0q0Go2zMMLyp5sHQ3NmEIAiEpsI1b7+pxUsskjJHrBhDH3ft2lWxSNuQKgwZwPBfkJ0SsVi66L3rWwKMhfKka0TORrddMplcMuFCsUQxOTnJ4OBgkUThOeGo0JGWv5kzBbNyKxyFSLdSBYNDcJLTyknXIbbQ5LiZtDpgv7OCg4zrE3j09+Hlb4DKXqf1wC5JY0SEsVjMPC9G6VItFpGZVD7SdLodKDnVNtL1+eyN9ZcS6SYy+e8zIt1kLockCLhqSCgaJG/Oqxt20OFs5uabb+bHA89ACpyaxNjYGPFEBLalmJ9Noc0OLbn1WdBHyChdSKswlnEl6nS/9rWv8b3vfY9nnnmm4nEZPgtdXV3cfffdnDx5clHStRiYi4BfEIS3XZOkqygK58+fx+v1VjW8AZBkieauJkLT9qRrpw+3tPvRVI2f/ugkbl/ewWyxOkEj+jZItyXgKXy/WCQvANzS28PXf/Iif/cvz7Kzs8WMnOfn54nFlmcEIwgCkUiE0dFRWlpayiSKSnPSDK/djFY+OscYw57KZcDmeeAQ7CNdyHc5KVXGrWvSDczGbqMj8Dhp9S40adeix7gUGBGh0+k0pSGoXD1gFPQb/6UKFp1yIVlWpulmFdraKjcX1BvpJsxIN//bJXI5vA5HTdspjaxjuQydrvy+pQoyUO/GLWxwt5PT5jkTgo7WzTjTnqW1Pus6kvoqicwB5MDqkO5ybB2ffPJJ/vRP/5Qf/vCHZXq3AaPDMBAIkEgkePrpp3nwwQcX3ba+cLO8CHwL+I01S7p2F1coFDJnhPX09Ni20tqhdV0zwQqRrhGhmuU18TjhQsVDs7+dg4d31vQdpUk5n9uJxyWjChC3kG48HudIs8z/cTn4dt8U37jjFxBFwdzGcgxvQqGQaXRTyajHbmQP5K0dIa/pWvdB13WS0Txpnr/4c6bUEZOQDB8Gh+AkpZVr1JAnXaPjqRJGQ++gLXAOV+YBUp6vg7D4KKalopSQKlUPGA0eRk3t8CtjAEzP5pfnmpYreshWKxlbChLZ4kg3kcvhlWur0bUb1bPNn5df8qQrmQbmipZfcXndnbQ3r7NtfU4kEkxMTBCPx21bnz3OKAIh4unNNLWsfGlfNBqtmXTtfBceeeQRMpkMx44dA/LJtC9+8YtFvgvT09PcfXe+3VpRFO655x7e+ta31ryPuq6ngb8SBOGf1izpWmHoqZIksXv3boLBYF11m23rWwhN2pOuEaHmcjn6+vpIpVL07toEnEDL1r4kLCVdQRDoag2QiyjEYumioZJ7d+7kv7ib+a/f/j7fPfNz/p8b95nbWArpGqPYAbMbqJIzmkeWiWbLo1KjcSIvL+RJKRQKcfHiRVKuPOnu3LOL7d6NJiEZS/WoI0HcEeXixYsmGRszvSTBj6Yn0HUVQbD/zRQtQFT/OK18HGf2L8i6fr/uc1APFosWrQX9Rk2tK5OPkFrbC+SVSfLSS/m2Z4/HQzKRRtNypFKpFZnym8yUaLq53JIaIyDfkeaX8xFzWlUACa9pYJ4nXdmmZGyx1mfj93eLZ7lhG8xFmkkJ07S0tKyY7wLkvUVq9dJdqu9Cb28vL7300pL3URCEAHAnsGfNkq4gCEVka9VTo9FoXa3AreuaGT4/VvHvr776KplMhu3bt9PR0cHEUH6qQWiu9qW+XflZZ6uPybk5NE3n9OkX2b9/jzlU8q0tLfzdqXP8f98/wR17txPwuMrqdBeD1ehm165dtLS0mLPFKsEty0wny/9utAFnNZVcTuf5559HFEX279+PU5mEnz9LTlNsCWls9hQDyVk6mzrNseNGba2jIw5NMBccpzmw3raKQBAEshwlJ/9HHMo3UOQ3oEk313we6sFS240ziUKtayGBumXrJg4f3m+SUE55DtC4dOkS6XS6aGmuqmpZ9Dk8FcIhSWy0GSMFFnnBWSwv1ALrd+UrFjLmfLS0piALbuTCA1AtkK4k1taRZp3qawzYlHN9kIWssp4mSWJqaqpovl09rc+lMH6vKz39oxIsXrr3AnuBY2uWdI1aVbvklcPhMEuDakHr+hZC05GiCCCdTptzyLZt28bWrVsXKg868hdgeL7+VmAD2WwWl6gSTeWlhd27D9DRsZAMEASB/3LnG3n3l/6OL/3wFB996xtqbnHOZDLmvu/YsYOOjo6avBeg8hh2Q14YnZqkV2li3759ZvLCES20M9v4LwD4pCaSWozmluay2trJ2AQTCsxGTzM2sq2oisCQJ4wbK+v8EJJ6ElfmD0l5HoNlFuzbYamOX5nCtA+hcP043fnzZZBQNqvS1dVuljZZGzxyuRwvvvgimqaZDR6Pfus0fq+bP/vgL9nujxnpOg2Tohy+JXSjxZQsqq7T4nCj6RpZTcctLtCCUhixYxfp1gpBH0XHQU7rsPVdMM7DzMxMXcnLou+4elvDjR37d8BDwPyaJV1BENi/f79tZFKv01jruhZURSUWTOAOOBkYGCAUCtHb22tmuq0/qtfvwumSCS8h0rVOg1jX1kRaHS90pZXXse7Z0MmvHt7H4yfP8as37qOnxV810lUUhcHBQWZnZyvOULPzw7XC43CQyhWfu1wux8DAAE5EZI8br+otyhY7TH8F+31rktvRUEmoYQLyQumWKIpsaHors6G/hLbvcmjrVwHB9GiNxWJMT08TCoVIJBI0NzfT3vwhNvs/hjP9KbLuP101Q+x6kS78fkJBf7cm0lRVQ1G0oo4069J8amqKm266qcirN5HK4pbg1KlTZRMsfD7fgqbrWtB0W921LbGtDmPRXP6h3+x0k9VSqLpgVqoAKFr+GpeXUTImqT9HE3dASeOM3aoIqrc+l0bFpSuEqxg6+UYJ6eqMyZeJep3G2tbnCeTFn73I888/T3NzM7feeivr16+3HXQpCAItHX5Cc9UTQKWfmZqa4sSJEzidTm699Va2b9mAWuCM0rIxA//5TbfgcTr4syd/XJEwVVVlaGiIEydO4HK5OHr0KBs2bLB9+i8W6VqrF1RVZXBwkJMnT+Y9ZJ0uBFd5xGFUL9g5jUGedAHCq21pcAAAIABJREFUSnlDiSi46PF9kIRyjtnMP5gerR0dHWzbto2DBw/S2dnJrl27WLduHWllO+Ph/4hD+1fG+/+M8+fPMzw8zPz8PNms/fSKerDkSNeY6lyIIK0lY9lCI4OrQhuwAcOrd926dSiawPp1Hdx8883s37+f9vZ2crkco6OjPP/88/QPjyCLApPj4wSDQRLZbM2arjUxHM4WSNfhJm1rYB4FhCUPpRS0EUTtLJpo7wpoByN5uXHjRnbt2sWNN97IkSNH2L17N83NzaRSKYaGhvjsZz/L6173OkKhEJ/97Gd55plniETsBxJA3nehq6uraCpEMBjk2LFj7Ny5k2PHjhEK2U/6Pn78ODt37mTnzp0cP3685mPRdTMS+UtgFPjGmibdSjdHPZFuLpcjnss/zbNxlaNHj7Jx40Zz25Ik2W6rpd1PuAbS1TTNNCjXNI2jR4+yefNmRFGkq9VvPvzjcXvSbfV5+J1/d4QTA2P8qG/U1nvhZz/7Wdm2K6FWecHYrq7rHD16lJ6eHlySTMbGY8Et5Ze1leSFLme+82w6M2L793bnLxGQDzOW/AI5bd72PZIk0dLSQk9PD63dH0GRXs/ODd9mx9YMLpeLUCjE+fPnOXnyJGfOnKG/v5/p6WlTO15tGHW6xjdZmyMymQLp1lG9kMrkzCjWaPDYvHkz+/bt48iRIwRa2/E5neZMwGg6TTwY4syZM/T19Zm6qd1vbZUXzEjX4SZTcBgzkmiQlxckIYAgLIEq9Bju9McAN1n5N+r/vAXWh/HWrVvZv38/H/nIRzh+/Di9vb20tbXxxBNPcPbs2YrbeM973sOTTz5Z9Nqjjz7KHXfcwaVLl7jjjjt49NHyJp1gMMhDDz3EiRMnOHnyJA899FBFcrbZb4OkDgJRXddfWbPyQjXYRaelMNpep6am6OzJR2KSKpcRViUCb+nwMz1a+cSX+iPs2LEDRVGKlkIbOprQC8vRSpEuwH+46QD/+/QrfPapn/LRw1vNce/9/f20t7dz88031zyupBrp6rqOlslrfOFYzHYMe1q1qeEtGF6nbBogAPxSK07BTUiZtv27IAhs8f8B58PvZCTxObYH/rj6QQgiGdd/w5O6hxb+CNe6b8D69eYxGG2vsVisaIqF1YOhUvZ8qZFuulCna5CubaRrQ7qVHgipTA6vu/Jvmszm8Lud5lyvzE9/yNbujezdu7es5dfQlQ2NPJPJmNehEem2ON2ktBCKLuG3tAArWmxp0oKexZ3+GII+SMb1eRS9C1GcXPxzdSKTydDT08O999676HvtfBe+853v8OyzzwJw7733cvvtt5eZ3Tz11FMcO3bM7Go8duwYTz75JO9617vq2dVfAR5jLRveQOVIt1J0Cvml1cjIiDnD/ujRo2RTeSIJTZUvTWRZJpMpJ5PW9gAXzoyWva7rOtPT0wwMDBT5I0xPT5dtZ1NXC01NHtSpXFXSlSWRj/37N/DbX/8uTw1M4XKcIBAILKkduBLphsNhLl68SLaQgNzUu62sksCYk1YKt5h/X6XhlIIg0CS3EVXso1gAj7SVDZ7/xETqr2nPvpUW5+uL/l5GTkILGdcjuNO/hSvzx2Rcn4bCvDC7tldVVU0ysk4tsDqTGU00y5EXtMJuWjXdTIF07VzG7Eg+m1NRVA2PjZRjIJnNmTW6OVUlp2n4HM6Kx25MO56dnSUYDKIoCpFIhP5cXvJx5DRixFE0kRZ5QUpQ9Uj9LcC6hivzh0jaKTLO/4Yq34paaKZYaSzXwHx6etp8cK1fv57p6fLAwOAKAz09PYyPj9e0fUtzhAT8pSAI/7ymSbcS7G4aVVUZHR01xyZbhz56/BIev9u2FbhipNvuIxpKoioqkiwV+SM0NTWVEaJdyZgoCtzzlhv4+sUfcr6vehSwuz3ATRvb+P7QLPfe8Xp29yzN9amUdBOJBJcuXUJVVfbu3cvQ6DAM9ZNWFJpLehA8koOUmlvIxxbgEh2ICKRtTG0MBOR2Ykqw6r5t8PwngtlnGE48TED+FpKYv/krkaAmXU/W+UFc2c+jKo+jOOyNSiB//pubm4tuUF3XTWcyo9MqnU6TzWbp6+srcuVarCQpk8zi8jrJFqQEp4Uw04WHutsmcrUj3VShMqEa6SYyObMbrZLvgvXYrW3P4+Pj6LpOR0cHz1zI/ybzYxNc0l5GQ0SL5uUlv99Pjvq9dJ3ZzyGrT5N1/B6K45eA8oaMlcJKjuoRljjks0b8ELgO+PdrmnRrOUGG7jkyMsKGDRs4evSo7RO3ZV0zwclyuaCyvBBA13UioSRI+cYJt9tdNrPNQCWbyHccO8Q3vvgjXnhljGQ6i9ddHF0apKgoCvf/8u38+lf/D5//wWn+4td/eUkXiEG6mUyGvr4+0xbPiIw8jmpj2GWSaq7MhVkQBDySi5RqLy8ANEltzGTtNV1z3wQn23wP8PPofYwkP8c2/ycXPR5F/nUk9QWc2T9HEw+iSfsW/Yx1v0uHTKbTaV599VVaW1vLlug+n89covv9/qLrKJPM5M1u0gUpwUKw6QIR25Gune9CMl2oTFhEXjBsQBOF36rWOl1N05BlGbfbTU4W8ctOrj94HZNzF2F4hJ6WfE5jenqaZPMcetrD2aGzRWV8lRo85Nz/wqF8g5z8TnKOhSW/tTtvJbFc0l23bp05fn1yctKsLbaiu7vblCAAxsbGuP322+v6Hl3X/9D495om3WrQdZ3R0VFGRkbo6urilltuqap7btqzkYGzw2WvVyLL1o58FPbTfztJd29bRWexxbYjiSKd7QHGw1GO/9NpfuftrwMWhkrG43FzqCTAr+zYwN+9OsrT5/t4y4HaWpCt0DSNeDzO888/T29vb9H8NKg+ht0lyQSzqbKrRtM0HLrExNwUlzKXisawG4QSkNtIqBEUPWeO8LGD33GQ9e57mUp/jWbHUdpcv1h9BLsgknE9hCd1D67M/aQ834BljpYRRZH29vayJbpRxjQ9PU1/f7/Zbeb3+4mGYjg9DjIFwnRY9Nt0uhDpVtB0yyLdghxRNdLNZtnUlo9AE+aontrrdI2oM5hN0uosjIjK5hshNgS66C6MWn8xmKHVv4V163YWdZml0+miuXZ+v58Wzwlcuc+iSHeQdX6kqJxvtSLdcDhsGtEsBXfddRfHjx/n/vvv5/jx47Z+vG95y1v4xCc+YSbPnn76aR555JElf+c1R7qGe5bhGVqr5+y+W3dz8p9eJDwToaVrYQlqF+nG43Hmo3lTZo+jmRtuuGHR7Vebb9bZ7ietKnzr+y9x7OYdqMkQ8/PzZUMlAd64qYNzkSyfefInvG7HZgLu2nwINE1jdHSU0dFRBEHg6NGjtktmTxXS9coOkkoWCl+p6zqzs7P09fXhlGWcPg/t7e3EYjHm5+dJJBJmAkvz5m+4qdQQPd7qD4tu728TU04zlPgUPrmGyFVoJuN6FHf6vbgyD5Fx/dmS63crJdIqeTAYdbXxaAJBguHBUWSHSH9/nxkZpgqVDaXOY8Y2Sn+HVIGkq0a6mZzZjVbvfDRrydhcJkGXO292E1HyzT7NDp+5b4oeQxabbetprXPtovP/xvrOTxFJ93Jh5p34fCNFjmSlSeSVQjQaLZrXVw12vgv3338/v/Zrv8ZXvvIVtmzZwuOPPw7A6dOn+eIXv8iXv/xl2traeOCBBzhy5AgADz74YJFVaDVYXMbMf69p0i0a221JYLW1tdHa2sqWLVtqNqfe//q8SfP5n1zg9XcvtJhaSdfaVnv9TXuBHxALVl5SW1GNdAMBN01RFyFHmk99+Z/54/e9iV27dtkv3ySJP/ilN/Ker/wD//+/nuTjd76x6vdah2F2dXVx8803my28dnDLxf6sVvhllzlPKxKJcOHCBTweDzfeeCPNF18gJyhl3rVGAkuMqoi6xP8Z/yI3hf4DAV9T0XLV+juJgsx2/59wPnIP/fFPItZg6ahJB8k6/19c2c+iKo+hOOrKLJuop3rB2vIqCzLNbc20t3Xi9kzR3t5OPB5nZGSECxeGABgbHUIgUdRlZTc1ImnYNroqX7uJbK6oMQIWpKHFYJU0ZtIJ9jfnpZWomk+iNsl50tX0BKBWrF4w5tq1NiXwpD6PLqxHav8rtnqdRTp5JpPJr4YcDiYmJswJHitBwtFotOZEmp3vwv8l773DoyyzP+7PM72n94T0hBAgtIQOiooFFQR1FdvaG67ddfVnXcuqa3eta3cVRUVRsSKo1ECAUAMhpIckpE7vz/vHZIaZZCYk6L7v6vu9rlxXknmee5567nOf8z3fA7Bq1aoB/5s0aRL//ve/A39fdtllXHbZZcM+vqBEWuD337XRBQL0qYMHDxIVFRVIYO3YsWNYVWl5E7KJTY7mk6e+Yur8SYGHUiqV4nK52LNnD0ajMaC/IAgCWoOKtqbBk0N+RDK6Xq8XieChu9vEOeeM5u3vdrO/1U5GRmRmRkFiLOeWjubD8p2cUVLIqLTwLVA6Ozuprq7GYDAwceJElEpfd4DBeKuDebp6uRKTy4HVaqW6upqioqJAckYrU2F2DWwxFJzAmme+gi8Ov4I7r41MxZiwIuJ+I6zXx5KpuYuDlrtRaD4Eji5y45Zd0BfffQavJAevdPJR9/mt4LD5EmkOuwulOlQ4vabWBdSSlZ0Boiukykoul+N0Ojl8+HAgVmodJBwB4PF6sbvcgRLgY/F0pVJf8veww0JCn6drcfsYNIa+Vj2uvmaSMskgMVPRhMp+E+DGrnoWQRKDVktInBygoaEBh8OBx+MJ0d7o39dOqRyegtxvmUj7rSEIQhqQIYriRkEQtIBGFMXDv2uj63A4KC8vR6fTDZAqHG5Vmlwh45KH/sTTV7zCD+/+wtxLZuN0OqmtrcVisZCbm0tRUVGIV5KencDBvUPjHvY3usEeqEYrw2p1c8n8mazf28YLH69j6pjMAUk1OJIIu3ZOGT/sqeGhL3/i3SsXIQ3yXI1Gn6KXTCYbkNg7mhcXqQ270+nE1mPEJXoR5DImTZoU8nmsXE+TLXwLIz+KddM4YN3O+t4V5KdOIDl5BMlB/FqHw4HJZMJsNtPa2ordHoMkcTZO/Ve09BYilZ45OJNAEHAo70dtuxKV/RbsqpfwSofXxuWYK9IsDqISo3DYXQO6Rjj6kmsJ8bHI5Ue8O39z1Pr6+pBYaVWtz9iZersxGRQBRTY//F0jNMpjM7p+T9fsdmL3uInv09I1uZ1IkKKT+t4jl7ePTiZEaMQoulE67kQQ67GrXkCUZEX8TlEU0ev1gfvtP45f29fu12rp/jcQtJKYCkwBNgIzgELgud+10VUqlRHZAsPVXwCYs3gG376+mrfuXkrS6BjMdhNZWVlotdqQh8WP0ZOy+OztdbQ1d5OUNviNl0gkAS+zo6ODAwcOBNqxH+7Yw/c/HMBqcXLr4tlc89gnvPnlZq4/e/qAcfzGW6/Vcusp0/nbx9+zbPMuzps8NhD+cDgcFBQUHBN/sX8bdq/XS319PS0tLcSofS+nUzrQOMXI9fS4zEc1WnPjLqbBXsVXHa9xSep9SIUjwjAqlQqVShUSN3Q4R7Oz81LM6hepbUrAbooenEkgGLCr/oXKfgUq+w3YVc/hlZYM+zoMFw6bE5VGgdMxsCmlze5CIhGQyQZqDygUCjQaDTk5R3rANZh3APUopAKNjY1YrVZEUQyIwFhE3zhHFMb6whHDMLpSqZRWhy+G6xcwt3hcaKRHOkI4+4yuQpIwcBBRROF8HJlnIw7FPUdVfQvumuLHb9HX7n/R6AYhjSOy/uOADOivQPE7gz+mFg7HYnS9Xi9n3H4ipm4z3774E1OnTiUtLS2iETnlvMnI5FLeeOzrIY3vdrt9fagOHaKkpISioiKUSiXR0b5z6O62UJyTzLzpRSxbtYPaloGhi2Ce7dziPKbkZvDCqk2s37qdyspKUlNTmTRp0jETxtVy34thdTlpaWlhw4YNR0qBE/poVaJ3QIgiVmHALXoCyZiI40t1nBL3Z9qdDazvWXHU41EqNOgsNyGVqCD5FSZMGsX48eMD96WtrY3KykrKy8vZtWsXdXV1dHQLGCXPIAoxqOzXIXWvH/L5H3NFmtUZoIwNaNVjd6FShe/qEC6m63D5VkR5OVmMGjWKSZMmMXHiREaMGIFcLqepzSeUfri5gcrKShr7CP2Ca2ArpXDwJ9I67H1GV6XF5XVi94joZUeW9y6vT8JU3t/oiiJy1/PI3Z/glF+CW75gSN851BhucNlzcXExpaWlTJw4kYyMjEDZ8+rVqykpKaGhoYGHH36YpUuXUltbO6Tx9+3bx7hx4wI/BoOBZ555JmSbNWvWEBUVFdjmwQcfHNLY/aAFdvb9bgN2w+88pjvYyzGc8EJwc8b0kanMu/pEvnr5BxZefzp5E7KB8C9jQkoU5159HO8+8z1b11YzYUb4rLxfRNzpdDJ+/PgB7YNi/Ea3x0oWcPVZU/lp20GeWfozz9w8f0AnA3+Ywuv1snh0Jltqm/hwVz3PXHTmryZ3+8ML1QdrKZUpQ9gf/hfS3tecMnjJl672vZjru3dzWuLgsdR87QRG66azvucL8jXjSVZmD7q9VIwnSfg/Wjx/pdb8ALm6x8IyCUKXqmZEzxLGZb2ARryJDsdtCOrTh63XOlQ4LA5fTLd7YHjB1md0wyEcT9fmcCGTSlAEhSKCvcIGm+/+Txo7hsKEaL7e3I1UEGisraW6H5XLX/IcbPD833m4z9NNVGqxeHpxiVKiglqvu7yHEVAgDabgiSJy10soXG/jkp2NS37DkK7Pr+Xp9u9rl5+fz4knnsisWbOYPXs2O3fuxOl0kp09+LMEUFhYyPbt2wHfZJCWlhboChGMmTNn8uWXXw77WIPu6UhgviAI+cBCoFYQBOvv2ugCETmckcp3g9FfH8HfnPHCexNZ+8kmXrrpLZ5Yc9+AppLBOPOSafz4+TZee+RLnvvsBuRBpZ42my0gWp2fn4/dbg/bry0mxre86+nxLaNi9GquWjCFp97/iR82V3NS2ZG+YH4hcz8HOT0tjStmTeLlNZvZUNPItLwRQ7twYWA2m6mqqkIAohPiKSoqCvncoOjTWBBDPSpRFMlR+yrklrb8yFR9ETqZr5OuIAhh43Enxl5AnW03Xx5+jT+nPTAodxdAI5SQrl5Ck/VZWu3vkaIOFVAJV+ggiiIux3iczltJVD1GfUctu1qPQyqVBSXs9CEx02PX03Wi0ihxHupFqw8tzXY43BGTYuG+z2p3DcrR7enTYI7RqH1Vj3I5Grk8oNUbTOVqbm4ekLSy2+14PB7a+zzdeJWWHncDLq+UGPmR59PhbUUhSQo5PrnrNRSu13HJzsKp+OuQqXn/DZ6uRCJBIpEwb9485s2bd0xjrFq1itzcXDIzM3+z4wqyE/8ENgFxwAdAMnDe7zq8MBgG83T9SawNGzb4epKVlpKXlxe4WLpoLZc+cj77Ntfwwzs/DxqqkMtlXPm3eRxq6OKzN9cCvgTf3r172b59OykpKZSWlgbiTuEmiOgYn3fhN7oAZ84cRVFWIs9/tBaTxR7Y1263s3v3bmw2G2VlZWRlZfHnGRPIjIvmH1/9jN01vJAK+Aoxdu3axe7du8nLy0Mtk+MO8zLFKX3HaewzuqIoBjSCdVI1j468EpvHwRst3wRi2B6PB5fLhcvlwul04vH4eqwpBDWnxl9Gh6uZtd3Lh3ScyaoLiVGcQJP1eYyuLUfdXhAEFKo4PPpXcEvnkhn3ETNKfqBkbDEJCQmByauiooLNmzezZ88eWltbcTqdwwpNeb1enDYnCo0Cp901ILxgPwZPd1Cja+1TBtP4JkFLv1Y9fipXeno6RUVFgfBEZmYmSqUSl8tFVVUVu+prUAlSWusbOXi4BrcoJVFxJD5q99Sikmb5/hBF5M7nULhewSU7E6fiLhiG8th/oyItXGhmuFi6dGlE4ZoNGzZQUlLCqaeeyu7du4c9tiiKO0RR/Jcoig+Koni3KIqXi6J46u/e6A5H3tGfxNq0aROdnZ1MmDCBwsLCsFzeOYtnMGpaAW/d8yFOi2vQl3DctDymzS1m2Ws/Ub5uG1u2bCE6OpopU6aQkJAQIhMZjjam0yqRyaR0dR+Jh0olEm674Dh6zXZe+Wwj3d3dlJeX+3q05eRQUFAQqLBTyKT8bd4smrqNvPnL1qNftD643W6qq6vZunUrCQkJlJWVER0d3aepO/B8/VluY1/TRb8B9XuzowxZnJc2h3U9u9hoqgp011UoFMhksgBNyd+eZoR8FKM1M9jUu5J6a1VE9TP/akYQBLJ196KSZlJj+htOT3jVsoEDKHEoH8Ypvwy5+1N07luIjbKTkZERkEr0xwylUil2uz0QJ965cye1tbV0dHTgcDjCTpqmLp/Epz5Wh9PhDnSN8GMwoxtJe2GwwoheqwOJIAQKY6wuF5qjVKP5wxPJycmoVCrGjx+PS68iXRtFdHQ09XafeJO31UVFRQV7q3Zjc9chcaficTtROB/uCymcg1Nxz7AMLvx3PF2TyTRop++jwel0smLFCs4555wBn02YMIH6+noqKyu54YYbWLDg6HHr/hAEQS4Igqrfj+53b3Qjob+8Y3d3N5s3b6alpYWxY8dSXFw8qEKXIAhc8/QlmLvMrHpt3aAdG7xeLyedPxrRK/Ll21uZOnVqWBHxSEZXEARiYjQhni5AwYgE5s8s4vOfd/Pjhh0UFxeTkpISdrlelpPOaWMLeHPtVmoPH13rs6GhIUT0PCnpyDJSJZNhC7NK0EjlKCUyTKKbvXv3UldXR3d3d8iKYlHqLAq1Gbxc9zmdzt7AElAmkwXam6tUKhR9WrAnxC0mSpbAisMvYnR2BbziYKMecg0FLXn6J/CKDqpNtwY6GxwVggSX4nocinuQeLehtp2F0n4XMtdypO7VyMQK9Dpf+a+fVVJaWkpubi5arZbe3l6qqqrYvHkz27ZtC2jWWiwWetp95bNR8YY+ylg/T9fhQhnBcw1XkWY9iqfbbbVhUCmR9N2v/p7uUFFn6SFLF0tcXBw9Ct/EMWvUNMaNG0d8MiC4sfXqsB9egty9nOaeM9h/6E90dnVHnIAiIZxH/2vR09Pzqzi6X3/9NRMmTAjhE/thMBgCrIrTTjsNl8tFR0fHkMfuqz5ziaJo7/dj/t0b3cE8XZfLhdFopKKigrq6OkaNGhWRYhYOOWMzmXfNSWz6ZBvVWwdmRkVRpKmpiQ0bNhAVp+Hca45jx/o6tq+vCTveYFVp0dGhRte/5J+UqSQ+SsOKzS2o1OpB9XBvnjsNtULOo1/9HPaF8BeSWCwWbDYbkydPDit6rlMoAjSk4H1FUSROqQaditzcXFQqFR0dHezYsYMNGzawfft26g7WcXHMHNxeD88d/DTiiymRSHwJH6WBc1Juxina+ar7VSQyIeAReb3eAH3IH8bweDwohRHk6h/F5jnAftMNeCK0dw8Ht3wBNvVy3LL5vn5rzodQOW5Dbb8GjfVk9JLXkUl8Kw4/OyYxMZHc3FxKSkooKytj1KhRxMTEYLfbfZ011m4GoNfWjcVkRyon5D4fLbwQLqY7mKfbbbUHxG5geJ2A/ehyWGm29ZKji0UURfaaupAJkKaKRyqVIlP7GBKFiRtJNFTgkN+IPPp21GoNPT097N27d8AEFEk03Y/fOoH5awsjPvjgg4ihhdbW1sCzW15ejtfrDdHiGAx974ooCMIXgiB8IwjC24IgPC0Iwt2CIFz7u0+kRYLD4QgUCeTn5x8zheqCexax5sN1vHPXR0xYNzYQqwwWEfdn+DMuHcHqFdvDJtVgcKMbE62hs8uCy+UKLGfz8vIoLi5G0CVz10tf8+EPlcwuTow4RpxOww0nTOaRr35m5Y79zCspDHwWXLar1+vJycmJGGPTKxSY+1rfiKKI1+sNvExxSi1dTlsgk+zXIvXHm41GI6LJxAlCMSuN23ml4mNOiBoX2F6j0Qx4+RKVGcxLvJzP2l5kTfdHzE24KBB3r6urIyUlJVDd5T8ODaVkqh+iznY3e41Xkad7GoUkYUjelChJxqn8G07FHQhiG4JoRBC7kLm/RM9/KE7+FI/zYlzy8yBMCWx/zdqefb7JMi0rDY+7HLfoYNu2bYHklcVsRxB8se3+okuRYrpxhsiOQZfFRoz2yCrN5HSQqB2aIwHgFUXu3fEDAgKnpRayom0dhxxQGh2HSqrA7TXRZH0OGVL04m4cyntwyxegVQysNHM6nYGClsEU2f4b+DVauhaLhe+//55XXnkl8L+XX34ZgGuuuYaPP/6Yl156CZlMhlqtZunSpcMqD+/D2/i4uvFAFDAKSPjDGV2r1UpNTQ1WqxW5XD6gcmq40EVrOev2U3jnzk/44Z2fmXTmWKqrq9HpdAM0c31JtdN54Oq3+ezNtZxz9XEhYw2WkIuKUrNv3yHKy8vJzMwMEaSZOS6HGSXZvPnFZooz5hCljnzbFk4sZsX2fTz17XpmFmQiE31tv10uV6Bst6KiYlCPRKtQ0GG1hizv/VqjcUoNzVbjgH2CmwwmJSWRK+bSut/CKuMupmtLkFqtPqlAqxWZTIbBYAgYYq1WS7F+Gs32Gjb3fkuUJwlJYwwGg4FJkyaFNVQAcbI5SCX/pNZyN1XGy8jRPIlamhtYsvuvX+QKNimikIpIKgAe2TQ6ehahcr9GtOQV5K7X8Uin4JbNxyM9LmIc09jpC3FEx/smhuycEUyaNDFQceVwrgc87Ny5E7fbHWg7rtfrcTqdAyY/m8OFelBP18bI5CNVYkanI8AsGQo+tTSw1tbCXcXH0eFu553m74iRWzgl/nhE0U2N6VYc3kbGKiS4Vf/AIzsx4lgKhSKiIpu/sWhNTQ0Wi4WdO3eG6G0olcpf5f3+Gk9Xq9XS2Rkqqn/NNdcEfl+yZAlLliw55mMDEEXx43D//90bXf9N87dM7+3tDbQd37Bhw2/QZkoiAAAgAElEQVTyHVMXTuKXD8t57c730GVfTsmksWi12rDbjp+ex9STRrHstZ+ZfcY4ElOPPBThPF0/bc1q7cZkdlBWNhm5fOBtuem8mVx0//u8uXIHS+ZHLm2VSATuPn02F7y6jL9//A0L85LIz88nPj4+aJvBW/bo5HLqnI7A0jf4xYhTaljTdhCb2xXg9IaDIAgsyV7IjTuf4+3uVTw26iqkgi9s4HK5MJlMGI1G6urqAmpkSeoxxMh2s8r0Hmdl30RhfGHYsYONaYL8ODSKN6jqXUK19SqytPcSLTs+kLALTlT5zyUSjQ3ALeZQ130rhTFepJ6vkbm/ReW4HY+kEIfyUUTJQGqRP5Hmp09pDarA8el0OlwuD0lJ8UyYMCGk7bi/k4Pb7Q5oL+j1eqw2J+pB+ql1WWzEaI6EF4xOJ/ohCjt92bSXb20tnJ0xmsmJSfyt6jXSVNEkqRuIV6RRb74Lo3srBTIZavVTeGTThjRuMPorsnm9XioqKsjNzfUpkhmNASGc4ZT89sdwxG7+v4IgCBcBpwJGwATU/+5juk6nk3379rF161ZiY2MHMAZ+LSwWC83Nzcy+qgyH2cnWj/dGNLh+XHbHqQgCvPHYypD/Bxtdvyzixo0bMZlMFBXl4vWKWK3haW5JsXouO2MyW6tb2byvNeJ3e71eFHYTszPiWFXTgiYtK8TgQvg27MEUMK1MjtnpDPBsgzE7yVeu+uSenwe9BgBxCgPXZM2n2tLExy0/Bf4vl8uJjY0lKyuLMWPGBCh1ZqOF6ZyLVozmi64X+bH8G3bt2kV9fT1dXV0RKYBaeSGjY95DI83loPlOWhzPIVdIUKlUKJVK5HI5Uqk0cN7BNLb+CTu/kfZKC3EpbsKm/hKH4kEk3gNobAuRO19A8IZqTJi6zGii1Nj77p3OcMQgejxenE5PgKfrXxEkJCSQk5NDSkoKubm5gY6/DocDq8OFsbuTLVu2UFVVRVNTE729vXg8HhxuNxbHEQFzh9uN0+PBMAShmF09bfx9148UyA1cXTCRRw/8B7VUybRYFVJBxOl8hcPOH8mQqYnRvXlMBjcc/Bx3f3w8JycnEB8P7nTc0NBARUUFW7ZsYc+ePTQ2Ng5I1Aajp6fnf7YEWBAEpSAIfwEWAfuAWnw9Vxb87j1dm82GTqcLK4XoLyQ4FqqKX0TcZDKRkJBAbGwsp19j5osXv+PEC2cxcnJexH0TUqI556rZvPfsDyGVan6j29vby/79+1EqlQGhnrb2fYCvFDgqSh123LPnjGXlut188GMV80+YjFZ9xLvx6wjX1taSkpLCfeedwe6XPuKRr37mP1edg0x6ZH7t7+kGx20FQUCnUGKO8KAfn5zLJTkTeftgBbOSspmVlBN2Oz9mxI1hU89ePmxZzcToQvK0RwSn/cdcX19Peno6U6dORSKRkOFM5a2m+9mXtIox0bfiNHsCMXS32x2IF+r1egwGg48VIU1iVMzrNJifotX2Hyzu3WTp7kQrLwxcez+8Xm/IOfuvhcfjCWTlPR5Pn0csxS2fhyjEIHN/hNz1FnLXu7hlp+CWLcQrKcTcbUEfrcXS1+dOE1Qc4Rika4T/WCQSSaD0VaePwuMVyckawfjxY0P6upnNZjr6OLpeu4Wuri7sffdVf5TwwmG7hZsrviRWoea6qCKerl1Gj8vM4rRi9lhXkKUQ6XL/RIosmmTdf/BKB2b0jxWDaen6zztYDtTr9QbCE8GKbMFhGZvNFljV/o8iE7gAOE8UxZAs/O/e6EZHR0f0PP20seEYXZfLxcGDB0NExHt6ejh06BAX3ruIdcvLeeGGN3hm3YPIwoQB/Jh/yfQBlWoej4e6ujqUSiWFhYUhZaz+qjR/KXA4yKQSliyazG0vfMurn23k5vNnAUckHKOiokLKdm8/dQa3f/QtSzft4MJp4wLj+I1u/ySZf+mtUyhwejw4PR4UYa7dX0ZOY/3heu7Z/j2fHndRoGgiEq7OPIPdxlqerlnGU6OvRymR09PTE5Cd7B+3jVUkc3bKTbzf/A++Nb/B+al3kCb4jLW/3NdoNAYUupxOJyqVqi9OfCkjVEU0O55gZ/d5aKT5qGQZSAQtUkGDSjqCeNUpyKWxIc+F0+mkrq6Orq4uCgoKQgwxgItSJPLJSGTNKNwfIHevQO7+EhEB6+FC9DExWIw+1oMuyOj6BcnDCZj7zyfYWTDbfFWUOrUybF+33c3tsLaKeL2Orq4uqg77WAbdh1o5qDhilILb6Tg8bm7Z+iUmt4NXxp/B583fUWVp4NyUEvZYV5Amt5KtrCVVnkuq7g2QDL6SGy4iVXNGQv+SX2BAWObuu++msrKSlStX8ssvvzB58mQWL148pPGzsrLQ6/U+loZMxpYtoYU2oihy4403snLlSjQaDW+99RYTJkwY0thB9zMJsIiiWCsIgh6wAzJRFG2/e6M7GPy0saFodHo8Hurr6zl06BCZmZnk5+cHYkv+BJjGoOHqpy7mkfOe5fPnv2HRLadHHE+uOJJU++T1nxg7O4WOjg7i4uLCKt3H9lWlNTR2Mn5c5FLeoqwkZo1JZflPO5k+Oh2psxuJRMKYMWMGTD5zinKYkZ/JS6vLObE4l+Qo30MskUhwu92BpF7/uG1U3/XqsFpJDUM+V0hl/GPCqZz/y/vcX/k9z5UOrvmgk6n5S84i7tv3Jq8eXMFsaz4ej4dRo0ZFnDBHqEcyL+kKVrS9zFft/+b0xKuQCJKQct/+zAmTyeSLFTdn4HbejzRmE07tPlyyAyCx48WORzRSb34Sg3wSBsUEtNJiLN0JNNZ1kZ6ezuTJkwPnEmx4/eEXt5iCW3oLguRypN6tyDiIsaeCqKg2HF1vAIlo9UGdgPuMbqTEWH/2gqmvq7BeG/6Z9ZcA56SlkJeWiFmvg6pdFGZlotfpg2Qx7chkMrRaLa91V7Gzp40nxp3CLttetnvrOSkujXrHChKkZvJVDWSrTyNRe0yiLkfFb1EYEZyoTUhI4IMPPuDaa6/l8ssvRyKR0NLSMqzxVq9ePSDs5sfXX39NdXU11dXVbNq0iWuvvZZNmzYN95ClgFsQBLkoin4yuQv+AJ7uYC/7UJTGgsVu0tLSmDJlyoAHJDgWO21+KVNOn8j7D33KjEWTScoMI3vXhzGTsxgzNZOPX/2Z0VMvID8/PyBT1x+pqTGMLEzh/Q82MntWIdFR4b1HiUTCaZPS2FnXyWPvreal2+aTmBD+4REEgb+eNpNz/rWUJ75eyz//dAperxeNRsOePXvQ6XQBFoHBYAh4m2WpPq/y29oaLh07LuzYBYZ4biqaweO7f+KThl2cnTkm4nUAGK3N4njVGH7o2kpSfBTn5kTOiPsxRj+DHtdhfu76BIfXxoKk65BLBhqj4BcyuLGg0zkTo9EYMMZWqxWJqg15TCVWsQKj8xUQRJBISB51EWmGiSHPUzj2g98Qe71xiOKJOLxeentqScpKxWL0aSvHSq/A5boDr3Q89r7OvoMVRwR/p9/oGjThjW6XxScU74/pGvuofXE6PQkJCSGymC6Xizf2l/NDVz0Lo7JpaN7FV5LtjFdaMXvXYpA5KNY0k6+7lzjVGWG/77fAf7MpZWZmZogs5m+Bzz//nIsvvjjQ1sq/0vVP8oMh6F5uxqe78KUgCEuBesADjPrdJ9IGw9H0Fw4dOsTGjRtxOBxMnjyZrKyssDNyf+N99VMXI0gEXrrprbDkf6/XG6j4WnT1NKRSKV++uQWZTBaRYyuRCNx4w4lYLE5ef+OXsNu43W4aGhqwmo1cu6CU9h47KzaEL8TwIy3GwJWzJ7G6qpZvd1bj9XoD8dOcnByUSiUdHR1UVlayYcMGKisrEbp7yDFEsaGpcdCxL8gez5T4ETy+ew315vBVcH52xubNmzkreiql0SP5oGMN23qrBx3bj5mxZzE3/mL2W7byVtMDdDqH7tEoFAri4+PJzs5m7NixTJkyhfHFp5GmvgZv/d049j+Kt/kvYC6jzfE2Fe1n0Ny1ZtDqw+AKO7lcjlKpxNJjRRefT4/9LNQaCTKpGa3zKpSOO3H0fgeASlqJ4PwFryd00u1fkWY0+4yqTh2ejdBt8QmcJ8h/ROpeicnmWxrH8RpK+y0o7bcgd/4bwVvLhq4mXq6v4MSkXOaP0vKNZDuj1N2ota2opG7Gqg4ha1tCzY6kgCxmZ2fnUYWihov/VlPKY2UvCILA3LlzmThxIq+++uqAz5ubm8nIyAj8nZ6eTnNz85DH71sVWYAHgA/xiZefCJwNzPz/nacbTkT8aOGH/uMkZMRx0X3n8Nod77Hu03JmLJocGNvfp83fgVgmk3HO1b6k2pS5hSTlRI5/ZmXFc/bCiXy4bDMnnTCKsWN9N97r9dLc3BxoI28wGCgrK2Hz/nbe+3orx0/MIzctfLWMKIosnjKG73cf4PFv1jIxM4VYne8Y/FSd4GW61WrFZDIRK5PT2NXJhg0bAgUVfhqQ/3pJBIGHxp3Mwp/e4W/bvuHt6ecilxx5ubq7u6muriY6OjoQt73Fcy5/2/MqTxxYymOjriZDHb7VUDBKo+cSI09kRdvLvN54D6cmXMoYw4yj7hfuWrS2ttLU1ERWVhbJyaUIgoDHczHNvcs55HqWRvfNNDYVQPfpaKXjAuet1+vDdpP2eLxYeqwY4vR0ml1oo/Q49cuROV5D7lyG17oPmEus/BU0rsOILjUuyUxckjl4hExU0r1oJEYUjg4EsRNnbxtQRoL0JhSO0YhCEoLYhSAeQhB76eoyoFfkEc2D4IAeUwkwhWT5LwhiNILoQeb5idruD7hz+zxG6WwsyXiF/6srJkXhIEZlRiKITNT2UBL9GtrEkWFkMX0dHBQKRYga27HKYv63mlKaTKZjMrpr164lLS2N9vZ2TjrpJEaOHMmsWbN+k2Pyr1wEQZgGjAU6gSZADzQDP//ujS4MLu8YbCz9yRulUklJScmQy4H9VWjBOP26uaxeuo5XbnuHcSeMxuGxD+hH5sf8P0/npy8refepH7n2HycM+l3n/WkyP/2yn+f+tYoXn7+Anp5uDhw4QHx8fCDe2N7uS57ccO4MNu1u4LF3fuSlv4a27AlOkkmA+848jotf/5QnvlnHo2efFPa7g+Ol2YmJbKnaw+TJkwPVfb29vTQ2NuJwOFCpVAGDdOfIWdy183terd7E9YXTApKWXq+X4uLikLitRqrk7oILuW33Szy0/12eGHUtBvnR70OedhxXjHiYz1pfZEX7y9TZdnNywiUoJJH1M4LhnwBiYmIoLS0NWe5KpVJGxJ5NungGbbZPaJG8iUv9FKJkPDbXWIy9cmxtVrzWLDSyghDmhMPs8rWiidHRsL/TJ+soqHGr/oJbeR3d3q3ARuRxj+PUOBCc3yF3r0Lh9XnAhj6SgOhWImKg2zobAJ0mHpn7WwTMiGgRhWREIZomUyYp0Xqs6s8BL02u3Wjlh5BEf4O973y6bI1ctftzdDI39xWZeKhlHEqJjBE6IwhOJqrMTIh5DZU0fcB9D5bFDK42a29vx2azIZVKI8piRsJwE2lDRXAr+eHA37I9MTGRs846i/Ly8hCjm5aWRmPjkVVeU1PTsbR51wOjgXZgAj6Obgtw3R/C6EaCXC4PtMjev38/ACNHjvxVykR+SKUSljx/GTfPvJcnr/0X8/86N2wyy3ccMq6/fwF3XvQaP3y4k0mlkTOhKpWcJdfO4f/uW84zz37OqacUhlS++Y0pQLROzU1/mskDr3/Pxz/u4E8njovISMhLiuPKWRN5cfVmTthTw4mjcgc9P/8k8+yWcm4umxKoNPN/5jfEJpOJdKOLqcp4Xtm/iZjDVkagCnBQw3lGicoY7sq/kP+rep3HDrzP/YV/Ri45+qNokMVxYdpd/Nz1Keu6V1Br282kqBPJUBWikKhQSbSopFoUwpHMvd3umwzdbjejR48edKKVCEpSNItJUi+kzfYxrbaPsErf9un/a0GIk6CRL0LuOJWeHp8s5KEan9KZ1W2mu9OIWqs4EqcVZFhtvtirWp+PV6YD2RR6jDfQ0rCSaL2dxORckOYhCil4vF66nFsQhErQ/ZNeQQRciCgCYY0m46dkxUcjSnwG87BtH4maI2W2Do+bm7dvpMsp4ZnS+Tzd8hkO0coYQwsivZSpZIxwPhQwuJEgCEKg3Dk44eR2uwOGuLGxEYvFx9joX/YbvCrweDzDbjh5NAxHbCcYFosFr9eLXq/HYrHw3Xffce+994Zsc+aZZ/LCCy9w3nnnsWnTJqKiooYUz4WQlfdaYA/Qii/EsA/4D5D8hzC6kTxdr9dLS0sLHR0d5Ofn/6ZEapvNhlVmomxRCeUfV3LujWcNWjQxcvwITjp7Aj98spXqS5rIHxP+ofeVyfYyujietetbWLz4hJBS4/5G7ITSfL4r38+/P9/E9DGZJMfpg5c4IdteNG0cq6tqeezrtUzMTCVGG54PDHDJmBI+3reX13dsY25OLsXxRxI0wf3MEhMTEUWR2xpjuGrnV7xhPsAzOcfT3t5OfX19YInq9wz92gsj9SNYkn0WTx9cxtMHl3Fr7rmBirXBIBGkHBd3DjmaMazpXMbqzo8GbCNFRpoqD50jAY9JQkZiNiVJU1FJh7iyEVSkaC4kRXMhHq8Fj2jFi4tD1rdos30MkmVo40aTn3UbGpePJxqfEo/Z2EBUnJqNGzcil8vR6/W0tvqUqTRqH33Rz/0uLDwzZPIXAIlUSq/FQbRejUaj7ps8FYGJ1O500tRtZHpuOi6XC0EQOGQ2k6D2nZcoijy4cxWV3Yd4cOwJfND6FT1uE2P0DYiCg9lR2aTbr0T8FWpfMpmMmJiYkHcpEq9WrVaj0+kwm82BTtS/VdGS/30f7nhtbW2BLhFut5vFixdzyimnhOgunHbaaaxcuZK8vDw0Gg1vvvnmsRyfRRAEmyiKXkEQbEB5H1+3VjjKjHFs08n/y3C5XCFkf6fTSU1NDZ2dnajVaiZMmPCrb/b69euZOnVqgMfb3d1NXl4eWpWO6yfciUqn5NkND6EI08HXj54uE0vOfI74pGj+ufQaZPJQnqh/XB9dTcVV175DXm4ijz68KOT4169fz7RpR6qFWjuNXPLgUooyE3nihnmDLvcOtHdx8b8/YXZBVsQwgx+9DgcLP/kQnULBRwvORhlmidjV1cWBAweIiYmhO0rBlZuWMzsph6dLT0cqSAJLVL9XbLVaA2Wier2eX1x7eb/tR2bFlXBj9iJkkuEtF43uLg47mnCJDuxeCzaPmcOmQ9Rad2ORdSHiey7kgpKy6FOYEn0aKumx81Dtnha6Has4ZH0Pp7edpp+KefdSKTd/lcu7D4mMLDVw80OX43WrMJlM/OeDcn5YVc1tN4/D5XIRFRVFcnJyQDqw/72688WVHO6x8PpdAzVeH1yxhi8q9/HQgjnMGZnFvs4OzvvyMy4oKubmiWW8ebCCF2vKuSJnInXerdTb2ynSt2CQOTk5bi5jDBfR1NSERCIhNTX1mK/BUCCKYmCVWVdXF+CG+yej4LLfY3k3jUYj5513Hr/8Ej7p/P8l+iaXwEn1sRfuEUWxGv4AlLFguN1u6urqaGtrIzs7m/T0dA4ePPibzK4SiYSamhra2trIysqisLAwMO4NL13BvWc8xvsPL+fPf/9TxDEM0VpOOn8Uy1/ayop31rHw8lkh3Xb7j/vnS6bzrxd/ZPWaKuYcXzRgPL8HFGdQc+WZZTz70Tq+3bSfU6eOjHgMeYmxQw4zRCmVPDjrOK7+5iteqNjMrZOnBj6zWq1UV/sYCMHL9jtGz+Yfu9bw5J5fuKN4dlhBlGDthVGmRI7zjmRNZyVNPa1ck3AaqTFJQ4oVAhhksRhkvmomi8XC/pr95MizmJt3EQqlHKvHRI+rnc2937Gu+3Mqer+nLPpUpsfMRzJMIW4AlTSVFM1FJKoW0m7/lOqeHwCwq3/EYpqBS1VBRceHRCumoVUVYvVYUahEtFnL0ev1mF2VtHsNdLTPxF41EkFU+pbkeiUSbTWt3Q3otFa2dcxDKU1DJR2BWpbDhv0avqhsI0GvZFpeBr0uFzeu/h6AE7Jy+K69hhdrypmTlECt5xvq7FCgaydeIXBW3G2kqkcHSqDDJQR/a/hlMTUaDR0dHYwYMQKdTheYhE0mEx0dHYFJeLCebuHQ29sbUlz0vwI/71oQhI+BFKAGWABYBUGoAEx/CKPr9Xqpq6sLUD385aR2u33YHYHDjd3S0oLJZCIuLi4sj3fCiWOYe+lxfPrUl0ybP4mCSeENmUQiIX98ElNOKGLpi6vJGZuAxdlFSkpK2HFPO2UsP6zaw6v//onSSdnogyqd+iuAzZ81mh8ranjxkw1MLh5B7CDSgMMJM0xPH8E5I0fx1s7tHJeZSUl8IrW1tXR1dZGfnx9Svgk+GlmjpYd3D25lhDaa87IGtj/3ay/49x3DGMa2VfBiw+c81P4hi7rLSLbpA4Ixfq9Yp9OFfRndbje1tbV0d3dTUFAQojylk0Wjk0WTri5gquN0fulazs9dn7Cl93tOir+AAu3EISfjgiGVaEnRXESUMxF4n8lZn/Oa7WkyE84iUV1Nt+Mnup2r6bHMQqVOAc1uTG4PGsVI7O4G7IaXEAwyFJIkbF47JjrBA529FxIX34bbmoZX3o1Zuo+GTjlPfb2AnOQObjj9C3Ybl/HkhtPoscv4z5lz6fEe4r5taykyGJGo9nDQnkCutpN8bRxnJ96GVhrl6xfnctHT00NycjJOpzOgRwGDKLH9BghmL4SbhN1udyA80b+nW7BXHDxZ/BpZx/8mgq7jK8AIIBZfPDcLOAnQ/yGMbkNDAx6PJ0DR8mM4HYH7wy9Ic+DAgcBDkpqaGnEGvvzRxWz9bgfPXPUqz254CPkgyv9nXzudbeurefepVTzy1pURkwwSicBfrj+RG276D6+8toZbbz7ZJ2ijULBly5YAhctgMKBWq7ntgtlc8cgynvtoLfdfMTfi98skEu4983gu/vcnPP712qOGGW6fPI1NLU3c+v233JWVx6jsbMrKyiKuIG4vnk2z1cijO1eTqjYwK+noHVpPSJpIjj6VJ2s+4k3bGk5Pmco5ybPB5sFkMtHU1ITZ7FPzCjbEFouFhoYGMjIyyMvLG3RVk6zM4pyUm9nY/RXbjT/xedtLCEiQCXIEBNyiizRVHqXRcynQTkQqHP31MHVZ+pJmvpctOioNjXUGrQdnk5aejMa7gxiDhYnx3wX2EUUvRtcWep3lOD2tCIIMt2gC+2RM5i4mZJ1HUWyB77xbO/j3d+VoZB6unphLlPMq/rW3g70dEq6Z9C37LB/x+J6TSVaZyI7tpcqWSKamk7KYPBYm3xCYUDo6OqiuriY1NZWEhIRAhV3/cudgBbbfyhAfjb0gk8kGlDv7ZTFNJhOdnZ3U1dUF9Bf27t1LfX09MplsWHHixsZGLr74Ytra2hAEgauuuoobb7wxZJs1a9Ywf/78QFfhhQsXDki0DQWiKH4f8XyHPdr/IHJycsJ6tINJGA6G7u5u9u/fj1arDTAHdu3aNajXrI3ScMOLl3Pf/Cf44JHlXPzAuQO2MZvNvgdJ3cXiv5zAm499y5oVOzj5nNJBzi2Bc88pZemH5UybkktpaRbjxo0LdMUwGo0BnVqFQsG8siw+W3+QVeVVzCktjPhADifMYDeZuDQpjUcO7GOZqZuXU6cP+qBLBQmPTTiVP69bxm0VX/La1LMpiTl69jdbk8KTxdfxdsM3fNW2kVWHt3J68lROTSojLT0NieATMDKbzRw+fJiDBw8CPu+pp6cnkJXW6/WDvuRTYuYxOfpU6mx7qLftxS06A4mZautWPm19HrmgJEtTzKSok8hWj454vqYuM7oYDTarrzKss7udjg49EydORKFQYLVWoNWExvkFQUKUoowoRVngfx6vl9vf+gqACYVpPr1ZtZr7vtuM2SXyysXzyY7R8+q2Laxt2s3C1ARGeObw2F4nOqmE0oRUKm1e0tXdnJo4mZMTLkYiSHE4HOzfvx+v18v48ePDtqjys12CleYg1BAfTRJzMBwLTze45bwf/nLv1tZWvv32WyorKxk/fjxxcXG8//77YdvuBEMmk/Hkk08yYcIETCYTEydO5KSTTmLUqFEh2x1r6/VgCAMfGP/fQ5jKf8cYbiw3mFpWXFwccsOHUlI8cW4JJ148i4+f/JKpZ04if6KvPNFut3PgwAEsFgsKhYLx48dTUuKl4qdq3nziG8ZNyyMpLTyzQhRFzj17Ehs21vCvl1bzYvEF6HSqQLVVMJ3H6XSSlt7Dlup2nvtoHRJbBwadOqTUN5jg7gsz1EUMMwTHbeeVlkFSIg+s/YmnN2/ktsmDy/5pZApenLyAi9d9yHWblvPWtHPJN4QvVw6GUiLnqqwzOCWxjKUtq1nWsoZlLWsAUEkUKCVyRI8XQRQYEZ9EnCoKCRJUooxSby7WNmuAI+ynMYUrbhAECdma0WRrRod8/wniYmqsldRYd1BlLucDy1ZUEg06aTQIAlqpgWLdNEbqylBLtXQd6iY6KZp9e3zXKTM7PURbw2J1khB/dIriK8s3sXVfM9ctmkpRVhKiKPL412vZ1nCIvy+YQ3FaEt/XHuT1vbs5NSePG6dN45L1H4EgYVqyikpHDSmqXma4RxFXX8z+Lh9Nrqenh4KCgpDy6P4IFnz3I1gQqb9H7K+iG44h/i28Zn+599y5c+no6KC0tJQ77riDjo6OITGTUlJSAtQvvV5PUVERzc3NA4zubwFxIEMh8Lf0/vvvH2zfQT/8X0HwQ9EfjY2NISV94WC326mqqqK5uZm8vDxycnIGdAju7u5GoVAcVUt39MyRrHpvLdtW7WTOBdOpq6+jpqaG9PR0CgoKaGlpIT09HYlEwujSbL75cDPVO5s57sySkEnCf04ejweJRKCwIJkVX2yno9PCtKnhvVJfQkLLqJxklpDeJywAACAASURBVP+0B310HGccPxGJRBKQB6ytrQ30SfO43YzPSuXTbVU0dvVy4qgcBEHA5XJRU1NDQ0MD2dnZZGdnI5fLGRUXT4/Dznu7d5Kg0YbQyMJBI1MwOzmHL5v28kXTXmYlZROjiBw/DkaUXMf02NFMjxlNuiqBfE0acV4tSpuEDF0SqfoEDrt7OeTo5JCjiz2WetZYdmLSuBifU8zYrCJ0Oh1ut5vu7m7q6+tpbGykq6sLq9WK1+sNdCgOhiAIxClSyNOOY1L0XOIVaSgkKhQSFRqpnh7XYXaYfmZTz9fsNW+i4uVaPDEWnMVdNP8iEnVGC4o4N/EKn2f+yScVpKVFM21qqAShy+tEgs9wfbVuL699vokFs4u5/Ayf9/vR5t28sXYbf54+jsVTxrK9rZW/fP81xQmJPHzcHG7YvIJ6Szcz0qTsc9aTobZyc86lzM6cj0ajoampKRAXbW9vp6Wlhd7e3oB0pVwuH9QQ+o2pVCoNqHH5jXOwIFA4iUz//n60tLQcS3HBoFi7dm2gonQ4wud+1NXV8eijj/LII4+EhPfq6up48sknef/991mxYgXjx48fdMIaBA9E+uAP4ekO5tH6havD3ZTgfmR+IelIYw3F0wVfe5/rX7iUvy96iqf/8hIX3LsopPWOXzxHJpORkBLNZXecyr/u+4yV72/i9AunRixuyM9L4tyzS1n6UTkzpucxuSyyyMfIzETOOWEsH/5QyZxJeYwvSBvgEftDEx6TiVOzE1lRVcu/Vq5mTnYy3d3dAaW14OshCAJ/nTKdJqORh9b9TLJWy8yMgZ0UgpGuieKVKQu5bP3HXLB2Kf+cMI9piYPvE4wRmiR0djnV9dWUxWWRPTI77FK1x2Xmq7aNrGzbyMbuPYzR53Bd9vwQalSwLGRnZye1tbW4XC7U6tDVgP8llAlyRuunMVo/LWSMQ46D7OrZSG37fmzNXkaMicZitQJKelXNfNG+jTWdyxhrmIXRYsEq6+T7w+/RaN+HxWPE5jHjEh1opAakbUV8vlROSWEC1y/ysUM2HWzi6e/WM7sgi2uPL6Oht5cbvv+aJK2Ox+ecyC0VX7Cnt43xSSYOOOzkaiX8X/5fiZImcODAAbq6uiguLg7J7ns8ngBrwB8fF0UxkKTyn/9gYZngZ9iP/h5xsLPQ//PgMX4tjEYjWVlZx7Sv2Wxm0aJFPPPMMwMYEP7W6zqdjpUrV7JgwYLAau+3wh+Cp+v1eiMmzDZv3kxJSUmI5+oXpGlubmbEiBGkpaUd9WFobGxEFEVGjIgsuxicfPv+2XVUfLmDx364h6KpBYFttmzZwpgxYwIvtiiKPHTde+zcXMuTH15N8ghfRj9ccYPL5eHm2z6kt9fKi89fGMJm6A+708VlDy/D4/Hy77vOQR9BtQp8jQqXvPM5WxvbuWVyAfnxUdjtdhQKRdjQhNnp5LKvPudgTw+vnXYG45OSB712AE3WXm4o/5wDpk7OzxrHTUUz0AzS7gd8KxB/PLKgoGBIZds2j4PvD29hafOPuEQPI3UZpCjjkElk9LhMmN02CnQZHBc3jjSVbyKy2WwhXGJ/mXN/Q+yfwOvr62lra2PbR3v55tU13PvpLbS0WXnrn9/y1s930Cbbz6aer6m17mHjHZNIO/4QI05rJl2VR4w8CbVEh1qq52DbId59zY1c7aTs/ErUagkqez6ffJVAskHPm5cuxCl6uGDFcowOB8+eMo3H9q6hqtdCdlwHyATGGpK4v+A6erp7A4myjIyMIYXW/EUNwUpsHo8nwBoIFokfDvwG1mq1sm/fPrRa7QAD+WvjxHfddRdnnHEGc+dGThiHg8vl4vTTT+fkk0/mlltuOer2WVlZbNmyJaIM5CCIeAP+EEbXXyceDtu2baOwsBCNRhPSXSE5OTmiqlg4HDp0CLvdHshq9kdwt938/Hy8Ti83TL4bgOfLH0GjVw84Hv+xd7T2cPOiF0nJjOOhNy8LKZroj5qD7dx6+0eUTsrmrjtPG/Tl2lvXxpJ/fs70sZk8cOXcsNtaLBaqq6uxub08sm4vXkTeuXwhcTpNiEfsL2zwG2K3QsEtG9fS47Dz9ukLKIg9entqm9vFc1XreK92GyO00fy1+DhmJmYNOC6/tnF7e3ug391w0ens5YPmH2mwtnHI0YlH9BIl06KWKqm1tuLFS4IimpKoXCZF+TpaRMt1yCWyEH1e/7nb7fYADTE2NhZLo5NHz32OUy4/nmuf/TPvv7CK5W+s5cMt9wSMiNFiZvHiN1i4uJAFZ5UQqzgyOZltDq57fDndJhsP3zwdUd9KbU89z33Si9MlMntuJQWxeXy8JY7GHjfzJ7SxoQsOWQxkxpiRKRxckH4iZ8ZNC0xMhYWFYRNlw0HwasBviENF4vUDRNLDjdHU1ERzczMFBQUhnZzDxYn9GI4hXrJkCddffz1lZWWDbtf/uC655BJiY2N55plnwm7T2tpKUlISgiBQXl7O2WefTX19/bFw/SPu8IcILwwGP23MT5mJjo4O6a4wVEil0rDhBX+yKbjbLgAquPX1a7nzpL/z6q3vcNOrVwfGCe6T5vF4iI7XccVd83jmzk9Y9upPnH/9nIjHkZuTyMUXTeONN9ey8uudzDstcpPKoqwkrlowmZc+3cBnP+/mrNlHkkb+yrre3t5AifQTaRlc8dYKbl76DS9ffAaaCMk6vyG+NTuP+/fs5LIVy/nHhDLyEhMHJOuCoZbJ+evo45iTnMs9ld9xffln5OnjOC+rhHRNFBIE7CYzvYfaSU1MImFUHlKlFrPLgUSQoJRKkQ6xoCFOEcWS7LPCftbpNFLevZftxgNs6NrND4crAp9FybQU67MpixlJvj6dmOhEEtxJNNTU4XK5SE9Pp7W+neeu/jdxI2IYf+5I9u/fz6HGw+ijQ6urnH0KNClRqSEG1+3x8sDrP9B82Mg/b5jH2LQ07K4cnvm4E6fdyaPnl2HXJvL02kYaulzMKj7IDmMShywyRkQ70apEbsm5gHRbFFu3biU3N/dY444DcDSReH9TSf9KKDg0odFosFqt7N27F4PBQGlpaYhTc7TwRP/WSf7jCWeIe3t7h13Wv27dOt59913GjBnDuHE+nehHHnmEhoYG4Ne3Xh8q/hCeLhBRA7SyshKLxYJOpwvUUh8Lurq6aGtro6jIVxnWv2w3kjf27gPL+PAfn3Pnf/7CjIVl7N69m9TUVAwGw4C47Qv3fMaaL7Zz36uXMKYsMrfV6xV54KEV7NjRxFP//BPZWZE9Qa9X5K6Xv6aiqonnb1lAwYh4mpubaWpqIjMzc4AozS/767nto2+ZkpvOk+eejOwoK4ED3V1c/MVydHI5j44vReZwhnjE4VgTAC6vh6+b9/FWTQXVpo5Bv8OPKLmKP2WVMCspm+KoJGS/QXzQ7fWwz9xAs72DLpeRNkc3W3ur6XGZQ7ZLkEWRoU0kURJN1ZUVtFa1888195GUE4/RaOTxGz/G4XBy/m1lyGQy9Ho9Xd1uHn50FX+9/RRmzjgSYnr+o7V8smYXty2exekzRuHxevnbJz+wpqqWx86Zy6yCTP66ZhXfHDzAPdNnstlSz7eH9pOgt5IbreTm9IVY67rQ6/Xk5ub+VxS8hgL/BOw3xj09PXg8HuLj44mLi8NgMAy5sjAYR/OIAebNm8eKFStCRNv/x/DHDi/AQKNrtVrZv38/RqORjIyMiGGBocJoNFJfX09xcXFI2W5qauqgM6Hb5eb24x/k0ME2Xih/hA7zYZxOJwkJCRgMhpDloM3q4I7zX8VucfLksmswxERmSvT0WLnhpvfRaZU8/eR5ERsfAvSYbVz9j09wutxcc3IuuSNSyMrKiviyLt+6l0e++pkF40dy17xZR53pK9tbuWLlFyRrdTx70inkRMcMGpoINsRut5t1+3bS0tNFWkY6UrWKHqcNl9eDXCKlx2nDI4p4RC/bulpY0+bj5hYaErh37AmMHQL/d7jwil52tlWzuX4nSp0aWZSKBls7bY4uah+vwvupkfRHCjn/gjOYFF2I1+HlutOeJXVSHMfdUYLghU5LD9veP8yeig5uXDKOqCifJnF5dRdvfL2Ds48fw5JzpgPw9PcbeH/jDm6eO5Xzy8bwwNqf+HjfXm4pm8I+9yFWNu8jVmdhTko68yWTcPZYGTly5P9MGazRaKSqqor4+HjS09MDfcz8amSCIATCEv6f4fJ2/YbXbrfz1FNP8c4777Bnz57/mWsQBn98o+t0OgOSgzU1NRiNRvLz8zGbzUilUtLTB5eyOxosFguVlZWIokhKSgqZmZlDfnCa9rdw49T/Y+TkPO5adiO9vb0B78DhcASy5waDga5DNu659C1KpuTwt+cXD2rwtlc2cs99y5k5o4Dbbz054rYWi4WfNm7n2RV7SE+M4oXbzkI9SMUcwL9+3MRb67bzlxMmc9G08C17grH5UAu3rPoWh8fD7ZOncXZh0YDjCTbEfmPscrmIjo4OeP9HE8rucFj4pa2OF/at57DdzJT4EYyNScEgV9FuN9NiMyIAUxMyKYpKIFcfj0o6dE/Q6XRSXV2Nw+Fg5MiRISujnz/awJOXvUzxpaPpvEJCq6MLQRSQvCaBLQLemz2IRb5XxluhRvzGgGamixPPLmCSoZCuA04eeXsdRRnRXDBrBDKZlI1tRt7dXsNZJQXccdpMnqko562dlVxZMp6a/4e98w6Pqty6+G9Kei+kE1InJBBaCESKIAgoSrkgoDQVEbHB/SiX4lVB4VJEkaugYAF7AdErICAgTUoSOoH0Qnqvk0yf8/0RZ5xACAkkIUrW8+SBtHPemZyzzn733mttUQ6H8zJxtqnhWd8IAosd8PbyNrYc3m3odDrS0tIoLy8nNDS0Tl/79T9nIOLKykrkcrmxj9o0T3wrT4gLFy4wd+5cRo8ezeLFi1vFQ+IO8Pcn3ZqaGtLT0ykoKCAgIAAPDw9EIhG5ubmo1erbbi+B2mm7iYmJqNVq+vXr1+h8sGkLzf6tR/hgzjZmrH6CMS89VOdnlEplHTI6vT+Z375LYPSMSB6dGlVnftn12L7jDJ99cZIZTw9g3Ni6Pr2meVuZTEZCdgVLP9hH/25+LH92OGLxzclNLwi88sNBDsWnsfqxYQwJvfUcqjy5nFeOHiImL5c+nl4sGzgYX/sb9fGVlZUkJiZiZ2dHx44dUSgUt4yIryfiaq2aT1POcDAvmXR5KQJgIZbgYWVHjVZDkarW51UqEnO/uz/9OnSii6M7AbYu9XZNGIo/2dnZBAQE4ObmVuecmVezWTB4OYHd/Xjzl0WIpGIuVaby87aTXPo0nT7PhPDEzCFYSy2Jv5rLujcO0inMCY9pcL4yGUU51OzxwN7enDfnDqWrcycOXU3llR9/o5d3B16IDGZ7eio7C3IJ72BFvnUpxUoRPnbwkmMPvHGic+fOd1woay6UlZWRmJjYpG4JU5jaQRrIWKvV1umcMBTsVCoVa9eu5ciRI2zevJlu3W5ex2hD+PuT7tmzZ403sWkUUFhYaCwWNRVyuZzExETEYjFBQUFcuXKFqKioRv2uoUhmqg3/z+MbOHfgMm8fex2/rjdvPdPr9aya8zUXT6Uxe8Uw7FxrZ6sZIgPTnkpBEFi9di+nTqfyxutj6NHD1zjep7687fZDF9n4wymeHBnB04/eXH4MoNRoeeGLXSQVlLL+8YeI9L91g7sgCPyQGM+66FNo9Xomd+nKqCAZnrZ2VCuVFGZloVYqCQkJwdbWtrbzRKejWqPBUirF2sys0akJw2vS6HVUa9U4mFkavZUzq8tJrCziQlkee3MSKFbVziYTAZ1snOjp7IWfrROO5lbY68VYF1TSwdmFgICAG3YwNZUK5g9aRk2FgvUn38DZo9ZQ5/yJFP7z0lf0fTCU+WsnIBKJKCmR88/532Jpacb6dY9ja2tBYVUlz6/9kcpqBVaP5CGy02JZ6UxGtC2d3OyZO7YbPyYnsO9qIXYOCjS2GlRac6Js7XlM1xEbS2scHR2Nf3tbW9u7lsfVarW13S4KBaGhoVhZNU7s0hgY7CANf/fo6GhWrVqFWq1GJpMxe/ZsBg8efEu5bxvB3590r/fUNeD6AlhjYCrblclkODk5IQgCp06dquNjWx9uJm4AqCiqZG7Uv7G2t2bdsdeNbWT1oaq8hgWTPkQsFrP2m1nYOlgZeyoNF6VOp8PW1hZzC2veWf87FZVKXvv3cCorCnB1da23JU4QBNZ8cYR9pxN5bcaDDOkddJMV1KKsWsHzX+wmu6yCVeOHMVDWOGFDQbWcd2JO80tqcp2LSCoS4Wptg4BAjUZDjUaD7o9r0FIiZWRgECMCAunl4YmVSUTaVCK+/jXnKaq4WlFAcmUxVyoKuFSWT5laYfwZc7GEDpY2RLn6MlsWhYeV3Z/v15T3id5zjhV7FtFlQK1tZl5mCYsmf4SLhz2rPp+JpbU5Go2WJa/sJONaCW+/NZFOvi5odTr+9f4vXErJ4525o/DrZM9PKef44KeriCw1uEQWUF1tSWGWMw72erS2OjR6gacdghnWIYCgoCCjotC0jet2tud3CkMPen0F2OaGSqVi9erVHD9+nOXLl6NUKjl37hz9+/fnwQdvPUm6DeDvT7parbbeKa5VVVWkp6c3aktisAgsKioytuE0ZB5uiobI1hSXj8Xz2qNr6D+uD/O3Pt/ghZt8OZt/P72VLhGdeGXTVCSSG/XxBiJOTMpi46Yz2DuY8+LsSNzdXY035PX5P7VGx/z/7iIxs4h3/zmaMP+GI4fyGgUvf/ULifnFzB4cydMDejb6hovPzmJf3GWEP8a+lCgVFNfUIBaJsDYzw1pqho157b/xJcXsT0+lRqPB1sycl3v34fHQLnVmv9V5HbdJxIZ+7YT0VJy9PSm1gHOlueQpKjlckIYgCDzqE8q0gF5c2nKKL5bt4Ik3HqPrjEg6WjsiKHW8/tQ2yovlrPn6WTx8avtQ39/0G/v2x7H4Xw8zoH8wgiDwzjfH2fX7VRZNG8zD93Umq7SCmdv+h5lEwqbpIzmYl8C7py7iYW9NvkUpXhY2PGsbxINdezVoXWi6PTd9CFtbW9ch4qa2Rt7sfU5MTEQQBEJCQpp99M71OHv2LP/85z+ZMGECCxYsuGtR/R3i3iVdhUJBfHw8vXrdfC6Z6bRdHx+fG1IUBtRHuo0lW1Nsf2sXXy7fwez103n42YYHVR7ceY4Plv/M2Kf7M+2fN1owGnwSqqqqqK6x4u31v9HvvgCmTu5ujIqAOtN8bW1tqaxW8cJbP6JQa3hv3hh83BxvOLYplBoNK3YfY39cCvfLOrH0kftxsb15+51hOKUgCMhkskZvQ5VaLWfycvk87hInc7LwsLElwsMTewsLajQaCqqrSSsvI9TFlYEdfeni2oEQF1fM/4job0XEEomEjIwM7O3tCQgIuCE6zK2pZFvqWX7MjEMcW4L72+kIA93InO2FXgToBby+qcEmUYvL3EBmj3qAns7eHDhwhfc2/sZj4yN4anptV8LOw5f57/YTPDGsB8/9I4riqmpmfvYzcqWaj58aTXJVKfN/O4C1pQS5rZz7rFyZH9QPmZ//bRXKTIUNpoVKQ57U8PdvSk0iPz+fjIyMZu0FvhmUSiWrVq3i1KlTbN68uY5x0F8Q9y7parVazp07V69yxVS26+rqSkBAQINPVcPIHkPe0LSXsDFka4Ber2fFY+u5eOQqqw/+m+BeDbezbV6xm1+3n2HBugncN6yL8RgG1Y+fn5+xcPjd97F88dUpnn1mIGNG9wT+1N0bbka5XI5YLKZaK+GdnRcxk0rY8H9j6OjeMPEKgsA3MZfZeCgGK3Mpz97fm3ERoZiZpDBM1WTBwcF1zKqbAkEQOJiRzq6UJBJKipFr1NiYmeFsaUUnBwfO5uVRUPPHUEQzM4b6+dPXywd/B0dkzi5Ymvwd1Wq10fSmpqYGMzOzG2S+10fEcRdTeX3YKsw62uH13mCCXdwIsHPh8CexJO9IJuBJGSfDyqnUqHAvt8bqRzWhXTxZtWwcEomY6CuZLNm0l/vCO/HGrOHUqDU89/kusksr+GDaKK4pKlh89BBiMwGJo5onnYOZ2WvwbfeRN/Q+muZJKysr6yjMDO+BQeZsgFKpJD4+HnNzc2QyWYunLmJjY5k3bx6TJk1i3rx5f9Xo1hR/f9LV6XT1KsZulou9XrbbmKrw6dOniYyMRCwW1ymS3daMpxI5/9f/VcQSMet/fwPbBnpyNRotrz+zjWtJBaz68lmsHcWkpKTQoUOHG/K2er3AqjW/EB2TxquvPEpk7/oJ3TDV9WpqDqu+Pg3A0w8G0TXQ03gz3mx+VUZxGWv2/s6ZjFy8newZ0tkfV1trCkrLUFRVEuLtgb2TE1q9gIVUgoVUSpG8GqVGy7CwwAYj5MZCEATy5HKuFBdyPCuTAxlpVP0hBbeQSAhxdsHZygo7cwscRWI8NFrCfDvR+w8Tn4YiYpFWzJrHNqFRalh3bBmu3rXpgxP74nhn0Q6G/qMnz78+GoVOy//i4/h29UlUIi2FYwUGdQqkq9iDH766hKeLPe/NH4NYImbO179wObuA5eMeYEf2VU5eywEzPb6uUpZ1eYDITsEtmiO9/r0zdMyYti5aWFhgZ2dnnDAhk8luS4LdFCgUCv7zn/8QGxvL5s2bm1R7aeO4d0kX6qYFTGW7ISEhTRrHHhsbS9euXY1P4dslXAMSY1JYOuI/9BrejaXfzm3wWKWFlSx4fDNiCTy3Yijh3cNu+qCoqVHzyqs7uZZZwhvLxtK1S8NdB9fyyli86ReKK2p48R+96eHvRGVl5R+TiaU3TKgwRPonU7P4/OQFLmblo9M37lIxk4gZJPPj0e4h9AtqeqvRzaDR68iqrCSjopzo3BxSy8ooqammrKaaErXaeCEP6eTH/0VG4e94o4RUrVZTWlLGW5M3kX4hi8nrRuPfoyP29vaU5StZP/9/BIR5seyj6ZiZSVGrtSx55QeuZZay+I2RHNZmsDc1Ac1hPSLA5yFHBvsFcOhkBmnZ5TgFm5GrUoBGjI0NTPMP4tmIQViY3Xne9U4hCAJlZWUkJCQgkUgwNzdv0PSoORAdHc2CBQuYPHkyc+fO/TtEt6b4+5NuQ05jJ0+epHfv3o2S7dYHQxohIyODvLw8zM3NcXBwMBLRnfRO7tr4Kx8v+oonV0xi3D9H1vszBsnx1XPX+HZdND36B7F4w+MN5v0qKmpYtPQHSkurWf76GEI7N6zcKpcreH3Lr1xMyWP84K48O7YvluZmdab5GojYcCPa2tpSWlpKeWUlnQKC8HB1plKpJq+8CnOpBHOpBLVWh0qrxdHaCq1Oxw9nr7IvLoXyGiU9fT154YFIOnu6YvnH9lWn11Msr8HO0gJr89vb0hqa9svKyggJCUFvYU5cUSEXCwv47PJFFFot93n78HhoV+7z9sFSKiWptASFVsuRlb9w5NOjzN70NA9PH4xarSYrI4+Vs79DEASmLO6Lk6sdtrZ2fPNdHGfOZrN08Uj63ReEXKFiztv/I7ekkoHj/bmgzSPlSgXSUikiLz0aiRi9TmC0pwfz+w2uM8vtbsLguldQUEDnzp3rFPBUKtUN05wNE31vtSO6GRQKBStWrODcuXNs2bKFkJCQlnhZdxv3Lunq9XqOHj2KmZkZ/v7+t5TtmuJmRTKVSkVFRYWRiExVZQYybmwOTBAE1kx9n+jddVuSDGvPysoiNzcXf39/3N3d2fddLB+v+oVR0+7jyfn1O4cZUFRUxdJXd1JaWs3SxSOJ6OXX4Fo0Wh0f7jzFD0fi8O5gz+LpDxAeeCNZK5VKrl27Rn5+vrEoY2FhUScibuhBpNHp2HUxkU2/xVChqJVvW5ubodHp0Oj+8F0ViQh0c6anrwcPhwfT1fvWvZmGHL3BNN7Hx+eG96e4poYfEuP5PuEKBdW1OWErqRSFVovdyQLcvk+n/AFPKsb683hYF2aEdufdOd+RkVjAys9m4N/ZE7VazYdbDvPrgUQeHhFAr54uiCVSth1KIymnnDdmDqVPuD//2XOM3ReT6BHsQUx1Hq5m5rzaK5L7w7q0CUUZ1Hb3xMfHG1sMG7Mu0wexgYglEskNRFzfsU6dOsXChQuZNm0ac+bMabIc+C+Evz/pXm/vaGrjqNVqiYqKanSrS1OLZKbFCsOHQV1jGhHf7AIzNN9XlVSx+sC/8ZZ5UlxcTGpqKm5ubnUkx4Ig8Mnqvez9Nobxzw5k8ksNdz+UldewbPn/yLhWwpyXhjJ0yK1zZucTc1jz5RHyS6qI6OxNVBdf7GwsqVFqUCgU6BQVuLnYExIUgJuzPRKxiJSsQvQaFdZmeuODyNLSss7rv75qLlepOZmSSVZpJRUKJWYSCRZSCc42VhTLa7icXcCl7AKUGi3dO3owMLgTfq6OdPV2uyEvbPBuNTMzIzg4+JZ/a61ez4nsLBJKiilW1OCULOfA3O9xj+zEfW+NIbWynB8T4vHeX4lNqpJxrwxj8mP9EIlE/Pi/c3zy6e+MGdWDZ2fej06vZ8WnBzl8Lo1nHu5KqLcNH8YmcrGwAhdXS7KlNfR2dGLt8Idwc2gb0W1jJbyNhUajqaMuq66uRiwWY2dnR1xcHK6urvzyyy9cuXKFLVu2IJPJbn3QvzbuLdItKSkhOTkZe3t7goKCiIuLa7R65nol2e3mrwRBMPbQGrwWDE79hojY1tbWGA3kpxfyryFvIjWX8OR743DxciYoKOimgwQ3v7mbgzvPMen5wUx4blCD66yuVrFy9R4uXcpm7JiePP1k/xt6fq9HjVLNj0fi2HMygdziyka/7l4h3vxjUFfu6+qLRlPXa0GtVhv7SA0ft9oR1Kg1/Hwhge2xV8gsrTB+vYOdNWYSCQq1BkcLMzyszegT7M+E+3pg0cTcYMr5dF55aBWeAe78Z/9SrO2tEASB9W/8xImdFykf7EhBfaAsXgAAIABJREFUuAWBjk70Ubty6rurDOgfxL8WPIyAwNovjrA/Oonnxvblgb7BzP32F9IKyxAcQGcjMNXXnyGubjfYITZ3jrSxuFMJb2NhKNauX7+eX3/9lfLycjw9PenZsyebNm1qM9F+C+HeIN3S0lKjbFcmkxnnmV26dAl/f/8Gi2a302/bVOj1+jr50aqqKsRisXE7duV0Ip+++B3ufh1Y/eu/G+xo0Ov1bHz9fxz5+SKPTo3iyfnDG7yItVodn3z6O7v2XKRLmBcL54/AtREDE/V6PVcSUsjMzkUWFICzswt5JVXIFSqqalQUlVej1wt4utiRX1LFz79fpaBUTicPJ16e2J+IEG/j+1jfjkCj0dSRN9vb29+0oCJXqkgvLudMRi45ZZVU1dSgrJYj14vIKJdTqVDh6+zAjAG9GBoWgKXZrcm3IKOIfw15EzMLKWt+exUXz9oC20/bTvDF+gOMmnYfE+cOYV9aCt8ePEf1gRJUTiIcRrozOTycuNPZHI5JZdigEHQuYvbGJqHR6dA7Qf+O7rzcbwAyk1lyN+uaaA0ibkkJb32orq5m+fLlXL16lS1bthAUFGTckfTs2bNFz90G8PcnXb1eT3R0NP7+/jeYG8fHx+Ph4VGv6XFrkG1DUKvVpKamUlhYiJWVFSkxGWz/9168gt2Y98VzeHXyvKlLv16vZ+tb+/nl62juHxnOi2+MbXDqBMDhIwls/OAwEomY8eN68fCI8JuO/SkpKTH2MDd2yoZWp+f3i+ls/vE0eSVV+Hk6EebnjoujNU521lRWKykorSLUz51eIV54udrfQMTX+0xcn5pRKpUkJiYiEomQyWRYWloiCAKn/5gtll5cjo2FOUND/ene0QNLqZTE/GKqVGqcrC1xtrEm2N2ZIBs7Fg1bQWVxbVqnY+faeWpHd1/kv6/8SP8RXfjn6vGIxWIuXspi+Zs/4+ZpT++pYexKT6LofDnmlSJqPATUdiKkcpCYiRkQ5MLLAwfg10iPgNYg4taU8AqCwO+//86iRYuYOXMmzz//fIvmbmfMmMHu3btxc3MjLi4OqJX/T5o0iYyMDPz8/Pj++++bbHp+h/j7ky78ae94PZKTk3FwcKijqLnbZGso+qSlpeHu7o6vr6/xwjy1+wzrpm/CPaADU98ei8i8tlBlmh81nbG285PjfP3eb/QcEMT8tyZg1cA8NICcnDI+/vQ4sWcyEItF+Hg7YWNjgUgEiEQIej06nRofbzv63ReGi4s9NQo1nh4OuLk1zr9Upday73Qixy+kk5pbQkWVEr0gIBKBtaU51YraVFBwR1dmju5Dn7A/t7mm8mYDEen1euN0X4MnRn0G1oIgcCYjlz2XkjickE6Nura4aiYRY2dpQUWNEp0gIFLp8Ps5BV12Bct3/Yuu/Wsr6Ed2XWDja/8jtFcnXv1gKmbmUi7HZbP8zV24udmxauV4qpRqXvlwL1mFFUT09SW5qpyckkp6ezoyd1AEssCAO946NxcRGyS8er2ezp07t7iEVy6X8/rrr5OUlMSWLVsIDKx/cnVz4tixY9ja2jJ9+nQj6f7rX//C2dmZxYsXs3r1asrKylizZk2Lr8UE9zbpZmRkYG5ujpeX1x0pyZoLVVVVJCUlYWlpSVBQUL03wrkDl/jP4//FK8idN3cvwsLO3HgTVlRU3JAfjT2Uyier9uLl58L8tRPwDb51lJWaVsjJk6lkZpVSU6NGEPQoFEo0Gg06nZic3Ar01/Xf9o7wY+aMgfj4NC1q0On1VFarsDSXYmkuJbuwgtj4bLYfukheSRWdO3Wgd2hHgju6IvN1xcPZDp1eT2mlAp1ej1inIiU5CRsbG8zMzJDLayc7XC9vNiU7rU5HXoUcpUZLJxdHzKUSdHo9RSVVvDJmLUWXsykbG4J9hC/Duwbilafnu7cP0DXSn8UbnsDS2pzjvyfx9vpf8fBwYNWKcaQXlvP6R78CMGxwCD9ciUevF5ga7sf0of1bdMveFCIWBIGCggLS09NbRcIrCALHjx9n8eLFzJo1i9mzZ7dqzjYjI4NHH33USLohISEcOXIET09P8vLyGDx4MImJia22Hu4V0r2Z01h2djY6nQ5fX99mKZLdLtRqNSkpKdTU1CCTyW7pen/h8BVWTliPh78bb+5ehKP7n/2T1+vsKysrSbtSyJ5PL6FWaBk2oSfjnxmEo8ut87aCIFBYWEhaWhre3n+aZNfUqEhKKkCl0mJpZUZ8Qh4//nQetVrL5Mf78sjIcKxvEVXfChqtjj0n4vnlVCIpWcXo/7gexSKR8f8AZhIRnf3ceGnCAEJ8ayNcgzm2qbzZMKXAQELXj4tRK9WsnPguFw9f5cXNz6AN9+Dg1VQu7I3H8XwF5v72PDJ/CMO7B3PycBIff3qcsDAvXl36KIfOpfDf7Sews7NE5QQlKiXBTjYsHtGP8CD/Vr+eoH4ilkgkRqmvTCbDzs6uRddWVVXFa6+9RlpaGh999NEdeVffLq4nXUdHR8rLy4Ha69vJycn4eSvh3ibd/Px8SktLCQwMRCwWtzrhGprP8/Ly6jXIbgiXjsazYsI7OLk7sOSbufh17djgeXIyC/ny3YOcOZyCWCLCJ9gZN28HLCzMQZCgVevx6OhE/xFdCQj1pLq62iiHDgoKuqUZSllZNRs/OMzp6DQsLKQM7B9Mnz7+dA7xxM7OgooKBWqNDk8Phya/xyq1lvS8UhIyCikqr6amWg5aBe5ubpRWazkYm0JpZQ2RoT707+ZHzxBvfN0dEQQoKpeTW1yJRqNDjBZnazE11XLj5BA7OzusLa355J/fceFgHHM+nMnQqQMB2PXFKbat249bV3cK+tqTVV6FfbYOqwId/qFu3D+uC3uOJpCRXoLWAlSOEObmwAhZR8YP7INFMzh5NQcMRuxZWVl4eNQOwmzJYp0gCBw9epQlS5bw/PPPM2vWrLvWkdAQ6QI4OTlRVlbWmku6N0nXkEZQKpWkpqYaIyHDhefg4NBkNU1T0FDetilIjE1l1RP/RVGl4Nl1Uxk6deAt15yZXMCR3Re5HJ1OUV45giAgNRMjEkNFiQK9TsCxgzUunja4dnDGxs6a+0eGE9731hMiBEEgKamAXw9e4dixJBTKG0UpHTrYMXZMTx4d2e2W7WnXo7y8nKSkJJydnfH39ze+Z1U1Kn44fJm9pxIoKK1NLzjYWqLV6qlWquscw9fdkRcf60dkaEf0eh2lJWW8N+sTLh2K56F/3k/vMd2wt7fn9N5Uft4aQ9SDofxz9XgEAVas3cO52GvovMwp6SBgUSFCogEXL1v6d/ck0EbKoF7d24yiDGo7BeLj47GzsyMoKOiG66y5i3VVVVX8+9//JjMzky1bttCpU+N8llsK7emFuwSD01hDRTKdTmfMixqauA2yVkOhqjlGojQmb9sUlOSV8fZTH3DlRCJh/WRM+NdoegxpurJJEAQSrqZwZPc5shPLKSusRqXUoqzWoFJomPTSQB6d0u+mHRPXQ6PRkpJSRFJKASqlBlvb2vfu9xPJXLqcTecQD156YQh+DUwsNsB0PllISIix5a++15BbXMm5xBwSrhViJpEQ4O2CVwd7zKUScooq+GLfOXKLKvF0saNfuC/538UQt/8SM9dOYdQLw1Gr1Xyz8RA/b4umS5QPI6aFIiBl+w9JpKaV8cTjvZG4WPHZ3nMICMwaFYG3jdooVmkrPaYNSXhvhdshYkEQOHLkCEuXLuWll17imWeeaRPvxfWku3DhQlxcXIyFtNLSUtauXduaS7o3SFej0aDVaptcJFOr1UYSrqioQKVSGYtUBiJurBmHYTCmQqEgODi4WaeV6vV6Dmw7yjf/+Ymy/HKc3B0IivDH2d0RS1tLLG0ssLKzpNugMAJ7+N3w+xUVFSQlJeHg4HCDjaW8spq3/7WdS6cy6Dm4E31G+NHB07HOrqAphtiCIHDkaCJbPj5GVZWS4CA3wsK88PdzJSTEA3t7K8rKakhNLSQ7uxQrK7C1URMREdak9MvNoNboOHo+lQPRSVz55CjSlEK4L5D7nx7M/T38SDyYxI+f/M7gUd15YfkYCgoreXPFLnLzygnr48WlwjKqlVqCPO0YF+WDs50FISEhODo63pXcbX0wSHhdXFzw9789D97rcTMijo2NRa/XExsbS0VFBVu2bMHX9+Yjp1oTTzzxBEeOHKG4uBh3d3eWL1/O2LFjmThxIpmZmXTq1Invv/8eZ2fn1lzWvUG6CxcuxNbWlt69exMREXHbBQRDE7+pv4JhNI6BhK+fyKDT6YwRh7+/f7MQx82gUWs59VMssXsvcO1qNpXFVSjkSpTVf46hD+jeiR5DuuDu54Z7gCtmHWrtDBuKIHVaHV/+9xC7vzgFIug7tDN9h8tw9baiRlGNWq1utJDBgIqKGg4ciicmJo2U1ELU6hs9jw0QiWDoA6E89FBXgoPc66QldDo9qWlFyOVKOnSwo6PPrW8gQRDY9PJWft12lP4zH0DR1ZtTcRmIE8uwyq3Bv28npi98iIT4PL754hQ6vYDSWYLaTETfLr48HNkJM005HTp0QCwWU1VVhVKprONFa9q+11rQ6XSkp6dTVlbWLBLeW0GlUvHRRx/x448/olKp0Ov1ODs78/3337e49eNfGPcG6SYmJnL69Gmio6M5d+4carWarl27EhERQWRkJF26dLltM2a9Xm+slFdUVNTJD0NtM7ZhNPvd2m7p9XqqSuQc236a33dGkxSbhv4PA5lHXhrKzP9MbdTa8rNL2f9dLL/uOIuyRo2FpRnuHZ3QanQU51fSo78fDz7eFR0qY/+sqbT5ZnlrnU5PTk4ZySmFVFbWoFRW4uAgod99PaiSa/n1wBV277mERqPDzEyCr68zzs42KGrUpKYVoVD8mTvu1dOXcWN7ER7uU2/OWK/Xs2XeF+z9+DceWzCKqa+PR6vV8d6rP3FibxzWIa5kOYmxqNRjVqNHZybCJcSFgb0D6B/ug7qyEEtLS4KDg+tcM4IgoFKp6qSo1Gq10fCoqdMZmorWkvAaUFFRwdKlSyksLOTDDz+kY8faQm5RURHOzs5/Z8OaO8W9QbrXQ6lUcuHCBU6fPk1sbCxXrlzB2tqaiIgIevfuTe/eve+IJA3+o2Kx2Og/anDbMkTErR0FGVBSUkJiQiKWImtOfHmWg58fY9y8R5i+fEKjb9QauZLLMelcOZNBYU45UjMJEqmY0wfjcXC2Ye6qcYT29K3TtmU6HsjwHpi2bRmMiK5du1Zn4oUBcrmKs+cySE0rIj29mMpKBebmUgL8XQkL88bFxYbExHy+33GG6moVVlZmeHo44OPjzLAHw+jRvSN6nZ73XviEw1+fML7mqvIa3vnXDi7HpDNu5kAEZ2t2/XKJmmoVffsH8sxTA3B3sTNOvTAMJG0MbiZvNk1R3engSK1Wa2w3bA0JryAIHDhwgNdee4158+Yxffr0Vg0m1q9fz8cff4xIJCI8PJytW7e2mfHzjcS9SbrXw+DPEBsbayTia9eu4ePjQ2RkpJGMnZycGiQmlUpFSkoKSqXS2Adp+j3T/LBBxGCqJmtJs2aFQkFSUhIikYjg4GCsrKxqDXL+73P2fXKYkD6BPLNmCiGRt68USovP5e2F2ynILmPExEjGPNkPN+8/Ccp0PFBFRQXyKjmqGh2uHo5UVlZiZ2dHSEjIHUWDKpWW2DPpxF3JobCwioTEPCorlfh4OmCRnEVGdDJTXhvPhIWjSLyQxTuLdlBeIkc2MJgr2aVoNDr69vFn8hN9CQxwo7y8nMTExGYrlNXXR329vNnOzq5R10JrSnihtntkyZIllJaW8uGHH+Lt3bAJfnMjJyeHAQMGcPXqVaysrJg4cSIjR47kqaeeatV13CHaSfdmMJiTR0dHEx0dzZkzZ6iqqiI0NNRIwt27d8fS0pLq6moyMjKoqqoiICCADh063PIGML35DGSs1+vrRILXK6luB4Y8X0lJCcHBwTcUDQRB4LevfueLZTsoyy/nwSfv5x9zHsYnxOu2zqeoUfH1e7+x79sY9HqBTjJ3ukR0IrCLFx4+zkikYnIyirl0Oo3zJ1KoLKshINyVR6ZH4ORubSzQmObIryUUce5EMoGhXvS6PxizRhjWGKDRaNm/6wKfL/wCVUE54lBfug7vQXlyIbmJBUgszVA4WSG1tWDokFDGju6Jj48TGo2GlJQUFAoFnTt3bvYZZaYwdZ4zfJimZ0yHZ0JtUSspKQmdTtcqEl5BENi/fz/Lli1jwYIFTJ3auHRUcyMnJ4eoqCguXryIvb09Y8eOZc6cOQwfPrzV13IHaCfdpkCj0XD58mUjEV+8eJGqqirUajXTpk1j3LhxyGSy285nmbqNGfLDEomkTqdAY3smb6YmuxlqqhR8t+ondm06gE6rw6GDPd7BHtg522Jtb8XgJ/rT44HGT2Etzq/g2J5LXIpOI+liNqrrenZt7C3x7eyEl28HTu5NQFGtxtPXmU4yd+ycrFApVFRXK8iIL6QkX278PVsHS/oN70JgmBfu3k508HbE3fvmO5CkM7W9zDWVCsa9Mp7o8zlcu5AFgoDE1ZZOvXwZMFDGoIEyHB2t68hk60tztBYMPhOGB7JcLkcQBCQSCTU1NXTs2LFV6gRlZWUsXryYyspKPvjgA7y8bu9h3FzYsGEDr7zyClZWVgwfPpyvvvrqrq7nNtBOurcLQRAYM2YMrq6ujBw5kuTkZGJiYkhJScHNza1Oftjd3f22b1yDCbTh5qupqaljclNfy5ZcLm+SmswU5QUVnPgxhrTLmeSlFCAvq6asoIKqUjnT35zIP+Y+3OTXotPpyUkvprSwErm8BrmyFJ8AV2QyGebm5lSUVnN872XiYtLJTitCXqlAJBIhlojpGNiBAQ+H06O/Pxejkzm66zJJ53PRmHQ7dAzqwKTZg+nzQGck0toHnkat5Yd1u9j+1i6cvZwY8uxwDu66THmxnP4PdWXq3Adx86orYlAoFCQkJGBhYXFDoexuwzCFF2pVVNXV1Tf4TDg4ONwgb75dCILA3r17Wb58OYsWLWLy5Ml3ve+2rKyM8ePH89133+Ho6MiECRN47LHHmDp16l1dVxPRTrp3gtLS0nq367m5ucZoOCYmhuLiYoKDg40ta7169bojxZthYqtpldzGxgZbW1uqqqpQqVR07ty52XqBldUq/jv7I078GMv9E6KY/e6T2Dg0bbttSHOUlpYik8nuSLWl0+kpyi0nKy2flPhsDv94iZJ8OVa2ZoT29sHd1ZpT356m+Fox3l07oTC3prJCiaybD08tGEFI97qSaYOQID8/n5CQkNa2+msQgiCQk5NDdnZ2vaPrTfPkhojY4MVs6jPRlGuttLSURYsWoVAo2LhxI56eDc/Ray1s376dffv28cknnwDw+eefc/r0aTZt2nSXV9YktJNua0Cn0xEfH090dDSxsbGcO3cOnU5Ht27djNFwaGjobRfS9Ho96enp5OTkYGNjY1TfNWcEJAgCO9bt5usVO3F0s2f8vEe5f+J92LvU3wuq0+q48nsil48n4BPugbm7CG9v7xvamQRB4PDXJ9ixbjdOHg50GxRG+P2heAd7YO/auH5qnVZH7JFEju66QMyP0WhLy0EiQdzBFXMnW4LCPRg6vgeR99f2rpoe01AoM4ytv9vRnCmqq6tJSEjA1ta2XgnvzWCYzGAg4urq6jppqpsNjRQEgT179vDmm2+ydOlSHn/88TYj+IDaKcEzZswgNjYWKysrnnrqKXr37s3LL798t5fWFLST7t2AoZXo7NmzxMTEEB0dTUJCAg4ODsbe4d69e+Pt7X1LEjD4ETg5OeHv728kboPTlmlOUCKR1ElLNFbSa4qkM6l8uuQb4k8lI5FK6DIghI4hXnTwdaGDjwtarY644wmc/vkMVaXVxt/zlnkSMbwbwREBWNlakpdWwNWTScSfSqa8sILAnn4IeoH0S5lGG06puRTvYA8enT2MQY/3w8Kq/jSJQq7k121H+H7Nz9RUKuj/WBQPTB+Ei6cjbt6OqNQK487AMD7exsaGmpoadDodYWFhLS4kaAruRMJ7MxjSVKZEbGZmhr29PTExMQQGBrJlyxZ0Oh0bN27EvZFG662N119/ne+++w6pVErPnj35+OOP71r75W2inXTbCgRBoLi4uE5aIicnBz8/P2M03KtXLxwcal26ysvLyc7OvqWazBQajaZOWkKhUBhVVAYybkz+VxAEUi9kcOLHWM4fuEzBtSJqKhXG71vaWhJ2fzD+fb0Z/I8BJBxN5fedMVw9mYRG9WdBza2TK2H3yYgY3o0Bj/VFLBZTWSIn/nQSRZkllOSWcvHIVVLPZ2BpY0HPB8Nro2CZJxKJmKKsEi4fjyd69zmqy2vo/kAXnlk9mU5dfBpce25uLmlpadjb2xuNj9pKH3VLSHhvBoPMffHixcTExKDVapHJZAwZMoQlS5a02HnvcbSTbluGXq8nJSXFSMJnz55FLpdja2tLUVER7733Hn369Lnt5nCDisq0f1ir1Rp7Rg3N+43Z1srLqynKKqFKXkWlphwfX298fX1v8KzNSy1EpVDh6uOCs8et87qCIBB3PIETO2OI3nOO0ry63qc2jtZEPtSDkbOGEtInqMFjGQpl5ubmBAcH13nAGPLkNzODb8ywzDtBa0t4AYqLi5k/fz4ikYj3338fNzc38vPzycjIICoqqsXPf4+inXT/SkhJSWHixIlERUURFBTE+fPniYuLw9LSkp49exoj4oCA2x8LY+gZNRCxqZLMND98fVpCqVSSlJQEYJxP1twQBIHSvDLyUgvRaXV06OiCWydXpLfo2zUtlMlkskYZnNzMZ+N2Hki3gkHC6+npia+vb4vnUQVB4KeffmL16tW8+uqrTJjQeDVic6C8vJyZM2cSFxeHSCTi008/5b777mu1899ltJPuXwkqlYqioiJ8fP7cPguCQHl5ObGxscZCnaE3t1evXkZFnaur623fWIYKuantpVQqNRJPVVUVxcXFyGSyG6rrdxsVFRUkJiY2y3b9+hltlZW1I+hNDY+aImhpbQkvQGFhIfPnz8fMzIz33nuv3nlyLY0nn3ySgQMHMnPmTNRqNTU1NW3Kg7iF0U66f0fo9XqysrI4ffo0MTExxMbGUl5eTkhIiLFQ17179zuaEqBWq8nNzSUzMxOJRIJYLDaauxgI6G72uRoIrbq6ms6dOzcq5307qK9gadqydTND/OLiYpKTk/H19cXLy6tVotudO3eydu1ali1bxrhx4+5KZ0JFRQU9evQgLS2tTXVGtCLaSfdegVar5cqVK0ZviQsXLiASiejRo4dRyBESEtKo7bJKpSI5ORmNRkPnzp2xsrJCEASUSmWd/LDB9tJ0O97SLVmmSrzW8iS4Hlqttk40bOgUMBBwUVERAKGhoa1SsCsoKGD+/PlYWVmxYcOGu2q7eOHCBWbNmkVYWBgXL14kIiKCDRs2tNhDsQ2inXTvVQiCgFwu5+zZs8a0RFJSEi4uLkRERBAREUGfPn3qyGAFQSArK4vc3FzjTLeGcL2UtaqqqkXHIikUChITE5FKpUa1W1uBSqUiMzOT3NxcLCwsEAShjv9uU83gGwO9Xs8PP/zAunXreOONNxg7duxdjy7PnDlDVFQUJ06coG/fvsydOxd7e3vefPPNu7quVkQ76bbjTwiCQH5+PjExMcaIOD8/n6CgILy8vDh+/DjvvfcevXr1um0hR0uMRTKkU/Ly8hpdKGtNKJVKEhISMDMzQyaTYWZmZtwZmEbEph0Td2r7mJ+fz7x587Czs+Pdd99tM7n2/Px8oqKiyMjIAOD48eOsXr2aPXv23NbxDJNg/kJoJ912NIyysjJmzZpFQkICkZGRXL16FbVaTXh4uDE/HBYWdkf52zsZi2RaKPPz82tT5tm3kvDW9/OmznNVVVV1UjTXu43VB71ez/fff8/69etZsWIFo0ePbnOkNHDgQD7++GNCQkJYtmwZ1dXVvPXWW006hl6vN6aqtFotUqn0r0LA7aTbEBYuXMiuXbswNzcnMDCQrVu3Gqusq1at4pNPPkEikfDf//6XESNG3OXVtgxUKhX79u2rc/MqlUrOnz9fxwTe1ta2jsnP9T26TUFjxiJZWVmRlpaGXC4nNDS0zeUEb1fCez3qS9FA/SY3+fn5zJ07F2dnZ9avX9/mIn4DLly4YOxcCAgIYOvWrbftd/HNN9+wZcsWduzY0Wai+VugnXQbwq+//sqQIUOQSqUsWrQIgDVr1nD16lWeeOIJYmJiyM3N5cEHHyQpKalNRVmtCUEQKCkpqWMCn5mZia+vr9HkJyIi4pYm8A3BdCxSQUEB5eXlWFhY4OrqaiTi5swP3y5aQsJ7Pa43g1+7di1Xr16ltLSUadOm8eyzzxIcHNymfCSaA6aRrCAIzJs3j4sXL/Lmm2/Sv3//u7y6RqOddBuLH3/8kR07dvDVV1+xatUqAKNUcsSIESxbtuxeavC+JQwmPKYm8HK5nLCwMGNE3K1btyblbw25UUOhTCwW11GRmRqg3w05b2tKeA3Iy8tj7ty5uLi48I9//IPExERiYmJ455136NSpU4ufv7Wg0+nqBDVKpZIFCxbw0ksvodfrKSwspLKykhEjRrR1L4abkm7LzY35i+LTTz9l0qRJwJ8O9gb4+PiQk5Nzt5bWJiEWiwkMDCQwMJDJkycDtblbgwn81q1buXz5MmZmZvTs2dOYHw4KCrqBrEy7Jq7PjTo7O9fZRpsOh8zKyjLmh1tyLJKphLe1zHP0ej1ff/0177//PqtWrWLkyJGtHuXrdDqjMdPu3btb9FwSiYTi4mLee+89IiMjGTlyJHZ2dsyePRsfHx/Mzc2JiYlBoVAwadKkv0p+tw7uGdJ98MEHyc/Pv+HrK1euZMyYMcb/S6VSpkyZ0trL+1vB3NzcmGp44YUXEASByspKzpw5Q3R0NMvuuAJlAAAU+klEQVSWLSM1NRV3d3djNAwQGxvL1KlTiYyMvGUKx8LCgg4dOhiVVqbFqcLCQlJSUpp1LFJ5eTkJCQl4enrSu3fvVrnRc3NzmTNnDp6enhw7duyuqbk2bNhAaGioUZnXkjBYTk6fPp3t27ezbds2vv32WwoKCoyz2pYuXUp2djbAX45w4R4i3YMHDzb4/W3btrF7924OHTpk/EN6e3uTlZVl/Jns7OxWH9L3d4BIJMLBwYGhQ4cydOhQ4M+K/5EjR1izZg25ubl07NiR1NRUY364Z8+ejc7fikQibGxssLGxMZpxm45FyszMvK2xSKYSXoO6r6Wh1+v58ssv+eCDD1izZg0jRoy4a+SSnZ3Nnj17eOWVV3jnnXea9djXpxIM59u0aROWlpZ8/fXX9OnTB6lUire3N5cvX+bdd98lKSmJDz/8sFnX0ppoz+kC+/btY968eRw9erSORv3KlStMnjzZWEgbOnQoycnJ92whrSWwdetWdDodM2bMQBAErl69anRbO3/+PIIg1DGB79y58x2lDZoyFqm1JbxQSzpz5syhY8eOrFu3rkUKdE3BY489xpIlS6iqqmLdunUtkl64cOECXbp0QSqVMn36dPLy8tDr9cydO5cxY8ZQXFyMi4sLW7ZsobCwkFdffbXZ19ACaC+kNYSgoCBUKpUxhxgVFWV8kq5cuZJPP/0UqVTKu+++y8MPP9ws59y+fTvLli0jPj6emJgY4xYb7p02tVvBkDIwNYFPTEzEycnJmL6IjIzE29v7jgjx+rFIKpUKnU6HVColMDAQFxeXFn/Q6vV6Pv/8czZv3sxbb73FsGHD7vrWeffu3fzyyy9s2rSJI0eONAvpmka3SqWSmTNncvr0aR566CGmTZuGWq3moYceQi6XIxKJUKlUTJ8+nWeeeaZNvCdNQDvptjXEx8cjFot57rnnWLdunZF029vUGoYgCBQVFdUxgc/NzcXf37+OCby9vX2Tb1DDhGCDe5uZmZlRvNDcY5FMkZWVxcsvv0xAQABr165ttpl3d4olS5bwxRdfIJVKjQ+mcePG8eWXX97W8UyFDmfOnCE3N5ecnBxmzJjB119/zU8//cT27dt57rnnqK6uxsPDg2PHjvHII4+wcuXK5nxprYF20m2rGDx4cB3SbW9TazoMJvAGt7WzZ8+iVCrp0qWLkYi7du3aoOdBfRJeU7TEWCS9Xs+2bdv46KOPePvttxk6dGibjeSaK9ItKytj7ty5JCYmUlBQwKxZs1i6dCmFhYWsXbsWGxsbli9fTkxMDCdPnuT++++nV69ezfQqWhXtLWN/FbS3qTUdYrEYmUyGTCZj+vTpQG1L2YULF4iOjmbz5s1GE/hevXoZidjf3x9BEIiNjQVoUMJrIFjTHKvpWKT8/PwmjUXKzMzkpZdeQiaTceLEiTY1u62lkJ6eznPPPUefPn34/PPPWblyJfHx8RQWFuLm5sa0adOYP38+X375JVOnTqVPnz53e8ktgnbSbUE0pk2tHS0DCwsL+vbtS9++fYE/TeANueGdO3eSkJCAWq2mR48ePPXUUwiC0KS+TzMzM1xcXIxEbToWqaysjIyMjDpjkQzqte+++46tW7fy9ttvM2TIkDYb3Zpi8ODBDB48uNE/b5q7NbynVlZWSKVSrl69CtTu5saNG8dXX33FnDlzCAsLY/78+chkspZ4CW0G7aTbgrhVm1p9aG9TaxmIRCKcnJwYMWIEI0aMYN++fbz66qssWbIEtVrN8ePHeeedd6ioqKBz5843mMA39hyWlpZYWloap+yajkX68MMPOXXqFEqlklGjRpGZmYlGo2lT1pTNBQPhbty4keTkZLp27crMmTN58803WbZsGbt27WLUqFEsXryY559/noiICO6///5mK1S3ZbSTbhvD6NGjmTx5MvPmzSM3N5fk5ORm3Wbt27ePuXPnotPpmDlzJosXL262Y/+V0L9/f37//XejlPTxxx8HalMGBhP4r776ioULFyIWi41qut69eyOTyRpd2DREeF9//TUJCQl89tlnREZGcvHiRc6cOdPsqrm2AqVSyaRJk7C0tOS1115jzJgxaLVaZs+ezSOPPMJnn33GgAEDiIqK4oUXXvhbSZlvCcOW6iYf9xS0Wq2g1+tb5Vw7d+4UvL29BXNzc8HNzU0YPny48XsrVqwQAgICBJlMJvzyyy/Ndk6tVisEBAQIqampgkqlErp16yZcuXKl2Y7/d4RerxcqKyuFw4cPC6tWrRLGjRsndO3aVXjggQeEBQsWCN9++62QmpoqyOVyobq6+oaPuLg44YEHHhDmzJkjyOXyVl9/ZmamMHjwYCE0NFQICwsT3n333WY/x6VLl4Rff/3V+PmOHTuElJQUYe/evYJGoxEWLlwo9OrVSzA3NxeSk5MFpVIpTJ06VZg/f36zr6UN4aa8ek93L+j1euLj4/Hz86vXMlAQBGOby18h73YrnDp1imXLlrF//37gxk6JdjQOgiCQl5dXxwS+sLCQoKAgYzTcvXt3vvnmG7744gs2bNjAwIED78o1lJeXR15eHr169aKqqoqIiAh++uknwsLCmuX4hw4d4sUXX+TLL78kODiYlStXkpKSwoYNG/Dy8mLhwoWUl5fz6aef8uSTTxr70q9cuYKNjQ1+fn7Nso42iPbuhfqQn59Pv379GDFiBJmZmYwYMYLXXnvNuHUUiUQ3bCP1er3x/381S72cnBw6duxo/NzHx4fo6Oi7uKK/JkQiEV5eXowdO5axY8cCtYWjxMREoqOj+emnn5g9ezZ9+vThxIkTWFtb37W1enp6GmXRdnZ2hIaGkpOT02ykm5CQwMSJE8nPz2fChAl06NCBQ4cOYWdnh1qtRiKR8OijjwIgk8nYsWMHZ86cqSMGutdwT5NuQkICIpGId955h/LychYsWEBycjJ+fn5s3LiR3bt3ExoaysKFC/H39wfqJ1q9Xo8gCH+biLgdTYdEIiEsLIywsDCefvrpNul+lZGRwfnz540dHc2B7t27M2zYMPr06cOSJUv4+eef2bt3LxMnTsTc3ByRSMShQ4dYvXo1vXv3Jjk5GS8vr2Y7/18Rf61QrZlx7tw5HnnkEXx8fAgICKBPnz7s2LGDb7/9lt27d/PRRx/h4uLC+vXrAbh48SJTpkxh69at7N69G51OB9QSsUQiaXM32fVo74xoPbS1a0EulzN+/HjefffdZlW8SaVSHnnkEaytrZk1axZRUVFcvnyZhIQEAObNm8eUKVN48skn2bRp0z1PuHCPk+7Ro0cJDQ0FaseuWFpakpCQQFFRES+99BJBQUEMHz6c5ORkFAoFx44dY/fu3ZSUlHDgwAFKSkr48MMPGTt2LO+8847Rbq4+6PV6dDodt8ihtygiIyNJTk4mPT0dtVrNt99+y+jRo1vkXDNmzMDNzY2uXbsav1ZaWsqwYcMIDg5m2LBhlJWVtci521EXGo2G8ePHM2XKFMaNG9esx46KimLHjh3U1NSwfv16pkyZQk1NDb/++isqlQoPDw/69evHiy++2Kzn/SvjniVdlUpFcnIy1dXVAOzfv5+UlBQGDBhAUVERgYGBQO0FGx4ezrVr10hMTGTGjBksWLCADRs2sHnzZkpLS1m/fj35+fl8//33QN28rwH1RcOCILQqEUulUt5//31GjBhBaGgoEydOpEuXLi1yrqeeeop9+/bV+drq1auNTm1Dhw5l9erVLXLudvwJQRB45plnCA0NZd68eS12nk2bNrFp0yZUKhVRUVEkJia2KylvhoZaG1qzv6K1kZ6eLvTq1Ut44YUXhPDwcKF///7Cb7/9JpSWlgpRUVFCWlqaIAiCMHPmTGHFihVCQUGBMGLECOHgwYOCIAjC7t27hX79+gmhoaHC6tWrhdmzZwvTpk0TUlNTBUEQjK1nFy5cEObMmSP0799fWLp0qZCQkHDLten1+lZrXWtJpKenC126dDF+LpPJhNzcXEEQBCE3N1eQyWR3a2n3DI4fPy4AQnh4uNC9e3ehe/fuwp49e1rkXMuWLRMiIyOF6upqobi4uEXO8RfCTXn1ni2kxcbGIpFI2LhxI1VVVRQXFxuLZdOmTWPs2LH4+fmh1Wp55ZVX0Ov1FBQUEBERAdSOpImMjGTKlCnExcWRlZWFr6+vsfVMJBKRlZXF/PnzGTNmDE8//TQHDhygtLQUqPXq3b59O0qlksmTJ9OtWzfj2m4nHyj8ES23tVyiKQoKCoyVdA8PDwoKCu7yiv7+GDBgQKvtpF5//XWCg4Oxtra+qx0bbR33bJ/uqVOnuHz5MrNmzTJ+zdR6Ljc3l/j4eIKCgujUqRNHjhzhhRdeMOrG4+PjGTduHPHx8Tc9x6ZNm8jLy2PBggV1jFLy8/OZMmUKgwYNQq/Xc/HiRdatW0dgYCA5OTnExsbi5eV1UyWa8Edl3HS9bREZGRk8+uijxMXFAeDo6Eh5ebnx+05OTu153Xb8XdHep3s97rvvvhvsEk0JzMvLq06ldfDgwRw5csT4eXBwMFOmTGHAgAGEh4fTt29fHnjgAaOcUa/XI5FI0Gq1ODg4oFarMTMzQyQSsWPHDsLDw3nttdeMx7569Squrq689NJL2NjYkJqaSlRUFCtXrsTa2prTp0/j7++Pu7u7MZo1rPfYsWPExcUxatQoOnbsiEqlwsLCos2Rsru7O3l5eXh6epKXl4ebm9vdXlI72tHqaDt3ZCvD0O7VEK7fBZiShFQqZenSpaxbt46goCBKSkqM0azwR8+uu7u70TbQ0LNYXV1NeXk5nTt3BmrHeRvEGd9++y2WlpZ8+eWXnDp1ikOHDpGamkpZWRkjRozg5ZdfJjQ0lLi4ONLT0zl06BAqlYpr164RHx+PRqMB4P/+7//44YcfEIvFaLXaZnm/mgOjR4/ms88+A+Czzz5rd1prJPbt20dISAhBQUHtxce/Ae7ZSLcxhiW3yo+KxWKioqLq+N+a/t6QIUPYv38/3bp1IyIigtmzZ9O3b18EQeDy5ctA7RY8MzOTnj17UlxczKhRo4Da+VxDhgwhIyMDjUZjFHHY2dmRnZ3Nc889B8Dy5csJCgpCJpMZI/NNmzYZ13K3DFWeeOIJjhw5QnFxMT4+PixfvpzFixczceJEPvnkEzp16mTs9mguZGVlMX36dAoKChCJRMyaNYu5c+dSWlrKpEmTyMjIwM/Pj++//x4nJ6dmPXdLQafT8eKLL3LgwAF8fHyIjIxk9OjRzaYoa0fr454l3eaC8Ic/Q30kbm9vzwcffEBaWhpxcXFG39WhQ4fyxhtv8PLLL1NcXIyXlxdTpkxh0aJFxgJEeno6YrEYc3NzEhMTGTRoED4+PpSWlrJz5078/f3ZvHkziYmJPPTQQwwcOBBLS0tjQTApKYnNmzczaNAgvLy8CAkJqXf9Op2uRZR033zzTb1fP3ToULOexxRSqZS33367js/AsGHD2LZtG0OHDmXx4sWsXr2a1atXs2bNmhZbR3MiJiaGoKAgAgICgFo3tP/973/tpPsXxj2bXmgu1OfPcD0CAgIYPXo0QUFBAPTr14/169fj7+/PqFGjePXVV7GxsSE8PJzY2FhOnz7N+++/j16vZ9CgQZw5c8Yo3UxJSaGqqopJkyYBtUqjqKgoY1fAyZMn6dy5M8nJybzxxhucPHmSxx57jGef/f/2zjYkrzKM47+/tilslUVoQrO9kBS9YCGBLOhDAyuK+qKUUU8Yoy8xgigqoWxQrVaLiGCrNpoxstFsrhjblEESrNqbNDZpZSXswzanrTSX5OO/D+d+5NEey5fNZ7bzgwfOua/7vs91zgP/c537XPd9lo/4k0wm6e7uBhiVOzwwMEBjYyODg4NZncQxVYqLi0c+7ZK+zkBzczOJRAKARCLB1q1bs+nmpMi0Xkac/zq7iUV3hhgrYqlk9ZqaGgoKCoAoiikvL6e+vp7S0lLq6urIz8+npaVlRHRLSkro6OgYSU3r7OyksLBwJIpubW2loqKC9vZ2qqqqqKuro6Ghgb6+Prq6ukbaVFdXU1ZWRm1tLb29vdimra2NDRs2kJeXNyrynehsumQyydDQUMbJITNN+joDcapazPnEf6WMxZwjFKmabP+rQknKA1YA79oeCGXvA8eA9UADcBR4wfZJSXuB1UAlsMv2J5IqgbuBjbYPSMoH5tjuk7QS+B74BvgCKATeAzbZPjTBc5ln+4/JXoNzhaT5wJfAy7abJJ22XZBm/9X2rBjUlVQB1NuuDPvPAdh+NauOxUyZONLNEmHWyijBVUSupJy0eoO2V6cEN/AWsAhYC3QDx4GeYLse+A4oBb4NZUVAEkg9lz4J7JDUAiwHFtr+EdgHvAPsBv6UVCppraSvJb0o6cqx5yHpYmCFpAOS1kl6RdLNU78y00PSHGAL0U2jKRSfkFQc7MXAyWz5NwX2AtdIWiRpLvAAsC3LPsVMg1h0zyOCECcziHHOmHpHbD9q+x6ghigKTgZR/AkYApbY/jm0LQJO2z4h6S7gXuA+4EHgF6AziOdioNH2LuAv4A1gM3A7cAlwa/BHKb9s94WoawvwCFAFVAR7blrdf3xoTFKhpE2SGjPZJ0s41nqgw/aaNNM2IBG2E0DzdI81U9geAp4AdgIdwGbbh7PrVcx0iLMXZgHjiLCAYdtJoDsI4HHgBkkXAanlwy4D5gG/pZoDPbZPSVoAFBMJ9RXAfCIRJrRfArwJ7CeKnHMk7bbdLyk3CP0twOtAO3AnsN92f+hj2LaDGD4r6YztVam2ROJ8LLRdBnw+zUu1FHgYOCSpPZQ9D6wCNkt6DOgCqqd5nBnF9nZge7b9iDk7xKI7C8k0Dmx7OAjvcIiO9oX9HqBeUmoy/GFgrqQ9wED4dQHXAoO2z4R6w8DHwBqiSPc64HhKUNMi6w+Bp4G2tLYpRDSV/CmgnyhqBsglGu64KtieAT6SVAT8nqGfiV6Xrxh/+uUdU+kzJuZsE79Iu4CQJIc/XFIu0XDC5USR7I1EL+XagdeIRLbWdsZvYku6FHgcyLe9cmz/Y+oeAe63fTStbDHwEvAQsNT2HkmfAh/Y3jFeXzExs5040r2ACI/6OWE7CfyQZj4o6TbgJqKXck3AMkkHiYYAPgO22T4V6l9NJNobIRLx0OcoJJUT3dyPpoQ0jB+/DawDFgCpzIdmYOF4fcXE/B/4G1GaiYaZ2LiBAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "ht.plot3D();" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Ordering HOLE profiles with an order parameter\n", "\n", "If you are interested in the HOLE profiles over an order parameter, you can directly pass that into the analysis. Below, we use an order parameter of RMSD from a reference structure.\n", "\n", "
\n", " \n", "**Note**\n", "\n", "Please cite (Stelzl *et al.*, 2014) when using the ``orderparameter`` functionality.\n", "\n", "
" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "array([6.10501252e+00, 4.88398472e+00, 3.66303524e+00, 2.44202454e+00,\n", " 1.22100521e+00, 1.67285541e-07, 1.22100162e+00, 2.44202456e+00,\n", " 3.66303410e+00, 4.88398478e+00, 6.10502262e+00])" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from MDAnalysis.analysis import rms\n", "\n", "ref = mda.Universe(PDB_HOLE)\n", "rmsd = rms.RMSD(u, ref, select='protein', weights='mass').run()\n", "rmsd_values = rmsd.rmsd[:, 2]\n", "rmsd_values" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can pass this in as `orderparameter`. The result `profiles` dictionary will have your order parameters as keys. **You should be careful with this if your order parameter has repeated values, as duplicate keys are not possible; each duplicate key just overwrites the previous value.**" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "scrolled": true }, "outputs": [ { "data": { "text/plain": [ "odict_keys([6.105012519709198, 4.883984723991125, 3.663035235691455, 2.4420245432434005, 1.2210052104208637, 1.6728554140225013e-07, 1.2210016190719406, 2.4420245634673616, 3.6630340992950225, 4.883984778674993, 6.105022620520391])" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ht2 = hole.HOLEtraj(u, executable='~/hole2/exe/hole',\n", " logfile='hole3.out',\n", " sphpdb='hole3.sph',\n", " orderparameters=rmsd_values)\n", "ht2.run()\n", "ht2.profiles.keys()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "You can see here that the dictionary does not order the entries by the order parameter. If you iterate over the class, it will return each (key, value) pair in sorted key order." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.6728554140225013e-07 427\n", "1.2210016190719406 389\n", "1.2210052104208637 411\n", "2.4420245432434005 421\n", "2.4420245634673616 389\n", "3.6630340992950225 383\n", "3.663035235691455 423\n", "4.883984723991125 419\n", "4.883984778674993 369\n", "6.105012519709198 455\n", "6.105022620520391 401\n" ] } ], "source": [ "for order_parameter, profile in ht2:\n", " print(order_parameter, len(profile))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can use this to plot the minimum radius as a function of RMSD from the reference structure." ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEMCAYAAADeYiHoAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3hUddrG8e9DCC10QhMIoSNdqigWRF0Lu/a2ll1dX9S1u+q6ru+y6q6uva6FV10bq4iAvaHCIiq6lACh995CSwgkpDzvHzNoxAkkw0xmJrk/15UrM+ecOeceyjxzzvkVc3dERET2Vy3WAUREJD6pQIiISEgqECIiEpIKhIiIhKQCISIiIVWPdYBISU1N9fT09FjHEBFJKDNmzMhy96ah1lWaApGens706dNjHUNEJKGY2arS1ukSk4iIhKQCISIiIalAiIhISCoQIiISkgqEiIiEpAIhIiIhqUCIiEhIKhAiIgmqsKiYdzPWsWHnnqjsXwVCRCQBLdiQzZnPfM2Nb2YwetrqqByj0vSkFhGpCgqKinlm0jKenrSE+rWSMYPC4uhM/KYzCBGRBJG5bie/evprHvt8Maf2aMnEW44jOSl6H+MVfgZhZi8Bw4HN7t6jlG2OBx4HkoEsdz+u4hKKiMSX/MIinv5yKc9OXkajlBqMurQfJ3dvEfXjxuIS08vA08CroVaaWUPgGeAUd19tZs0qMJuISFyZs3YHt46dzeJNuzi7byv+MrwbDevUqJBjV3iBcPcpZpZ+gE1+DYx399XB7TdXRC4RkXiSV1DE458vYdSUZTSrV4uXftufE7o2r9AM8XiTujOQbGaTgXrAE+5e2tnGCGAEQFpaWoUFFBGJppmrt3Pb2Nks25LLBf3bcOfph9OgdnKF54jHAlEd6AcMA2oD35rZNHdfvP+G7j4KGAXQv3//6NzGFxGpIHkFRTzy2SJenLqCFvVr8coVAzmuc8i5fCpEPBaItcBWd88Fcs1sCtAb+FmBEBGpLP67chu3vz2HFVm5/HpQGn86tSv1alX8WUNJ8Vgg3gWeNrPqQA1gEPBYbCOJiETH7r2FPPjJIl75diWtGtbm31cO4qiOqbGOBcSmmesbwPFAqpmtBUYSaM6Kuz/n7gvM7BNgDlAMvODumRWdU0Qk2r5dtpU/jpvD6m27+c3gttx+SldSasbP9/ZYtGK6qAzbPAQ8VAFxREQq3K78Qh74eCGvTVtF2yZ1GDPiSAa1bxLrWD8TP6VKRKQKmLokiz+Om8P6nXv43ZB23HpyF2rXSIp1rJBUIEREKkB2XgH3f7SAN75fQ/vUFN6+ejD92jaOdawDUoEQEYmyyYs286fxc9mUncdVx7bn5pM6Uys5Ps8aSlKBEBGJkp27C7j3w/m8PWMtnZrV5ZlrjuKItEaxjlVmKhAiIlHw+fxN3DlhLltz93Lt0A7cMKwTNavH/1lDSSoQIiIRtGP3Xu5+fz4TZq2ja4t6vPibAfRs3SDWscKiAiEiEiGfZG7krncy2bF7LzcM68R1QztSo3riTrsTdoEwsxQgz92LIphHRCThbN2Vz8j35vHBnA10a1mfV64YQPfDEvOsoaQyFwgzqwZcCFwMDADygZpmlgV8CDzv7kujklJEpAyyduVTr1b1Cr3W/+GcDfzl3Uyy8wr4w0mdufr4DlGd5a0ilecMYhLwOfAnINPdiwHMrDEwFHjAzCa4++uRjyki8nO78gv5bvlWvlqSxdSlWSzdvIvLBrflnjNCTlYZUVty8vnLu5l8nLmRXq0bMPrcQXRtUT/qx61I5SkQJ7p7wf4L3X0bMA4YZ2axHXpQRCq1wqJiZq/dwdQlW5m6dAuzVu+gsNiplVyNge2asDk7j225e6Oawd15b/Z6Rr43j935Rdx+ShdGHNOe6pXkrKGkMheIUMVhHzNr6O47DrSNiEh5uTvLs3KZGjxDmLZsKzn5hZhBr1YNGHFse4Z0SqVf20bUrJ7ECY9Mjmqezdl53Dkhk88XbOKItIY8dG4vOjarF9VjxtJBC4SZ9QOGA08ChUD3/X56AClAw+jFFJGqImtXPl8vzWLqkiy+XprF+p15AKQ1rsPw3odxTKdUjurQpMLmZYZAoRo3cx33vD+P/MJi/nza4VwxpB1J1azCMsRCWc4gngeuAlYDOcA8YCGwgMBN6z6aN1pEwrVnbxHfr9zG10uz+GpJFgs2ZAPQoHYyR3dswnUdmzKkYyppTerEJN+GnXu4c/xcJi3awoD0RjxwTi/aN60bkywVrSwF4hvgNmAmUAf4P3d/C8DMblNxEJHyKCp25q3fyVfBM4TpK7ezt6iYGknV6Ne2Ebf9ogvHdEql+2ENYvoN3d15a/oa/vbBAgqKixn5y278ZnA61eLsrOG/d54Ytb4WBy0Q7n6DmdVx993BFkt3mdnNwD2A5oEWkYNas233DwXh62VZ7NgduF3ZtUU9fnNUW4Z0asrA9MZxM+z1uh17uGPcHL5aksWgdo158NxetG2SEutYITWoE722QWW6Se3uu4O/twG3mFlb4G9AczMb6u6TopZQRBLOjt17+XbZVr4K3ktYvW03AC3q1+LEw5sH7yOk0rRezRgn/aniYuff36/m/o8W4MC9Z3Tn4kFt4+6soaKE1ZPa3VcBl5rZI8A/zOyv7n5cZKOJSKLILyxixqrtP9xcnrtuJ8UOdWtW58j2Tbji6HSGdGpKh6YpmMXnh+2abbu5/e05fLt8K0M6pnL/2T1p0zg29z3iRXl6Upu7/+SSkrtnAKeY2dDSthGRysfdWbgx54fmp9+v2MaegiKSqhlHtGnIDcM6MaRjKr3bNIz7XsXFxc5r01bxwCcLqWbG/Wf35MIBbeK2kFWkcvWkNrNxwLvuvnrfQjOrAVQzs1cI9LZ+ObIRRSQebNi554eC8PXSLLJ2BTqkdWiawgUD2jCkYyqD2jemXq3E6S+7MiuX28fN4fsV2ziuc1PuP7snhzWsHetYcaM8BeIU4ArgDTNrB+wAagFJwGfA4+4+K/IRRSQWcvIK+G75NqYuzeKrJVtYtiUXgNS6NTi6YypDOqYypFMqLRsk3gdqUbHzr69X8PBni0hOqsZD5/bi3H6tddawn/L0pM4DngGeCQ6pkQrscfcd0QonIhWnoKiY2Wt2/NDaaNaaHRQFh7EY1K4JFw5IY0inVLq2qJfQH6TLtuzitrGzmbl6B8O6NuPvZ/WkRYNasY4Vl8K9SV0AbIhwFhGpQO7Osi25TF2yhalLtzJt+VZ2lRjG4urj2nN0xx+HsUh0hUXFvDB1BY9OXEzt5CQeu6A3Z/ZpldDFLto0YZBIFbIlJ59vlmX9cJawocQwFr/qcxjHdExlcAUPY1ERFm/K4baxs5m9die/6N6ce8/sQbN6Oms4mHIViOCcEHe4+31RyiMiEbRnbxHfrdj6wzAWCzfmAD8OY3F9jIexiLaComJGTVnOE58voW6t6jx10REM79VSZw1lVK4C4e7FZjYcUIEQiXPPTl7GYxMXx+UwFhVhwYZsbnt7Npnrsjm9Z0vuPqM7qXXjq2NevAvnEtMcMxsJ3Ltv0iARiT/frdhKat0a3H9Or7gaxiLaCoqKefzzxfxz0lIa1E7m2Yv7cmrPlrGOlZDCKRCNgeOAa8zsO2AOMMfdx0Y0mYgcsqb1anJc56axjlGhPp23iU/nbeKMPocx8pfdaZxSue6nVKRyFwh3Px/AzGoSmA+iJzAQUIEQkZhKTalJTl4h953Vk5O6NY91nIQXdismd88nMAT4zMjFEREJ34u/7U9yUjVqJVeNy2nRpmauIlJpJNIwH4kgvkfREhGRmFGBEBGRkMpdIMzsPDOrF3x8l5mNN7O+kY8mIiKxFM4ZxP+6e46ZDQFOBF4Eno1sLBERibVwCkRR8PfpwCh3/xBQQ2OROOLu7C1UP1Y5NOEUiHVm9jxwIfBRsD+E7mWIxInvlm/lguen8c2yrTSvrwHpJHzhNHM9n8DkQQ+7+w4zawncFtlYIlJeGWt28Mhni/hqSRZN69XknjO6c8GANrGOJQksnAJxa/B3j/1GRPysLC82s5eA4cBmd+8RYv3xwLvAiuCi8e5+Txg5RaqE+euzeXTiIj5fsJnGKTX482mHc8mRbavM2EsSPeEUiNwSj2sR+LBfUI7Xvww8Dbx6gG2+cvfh5Y8mUnUs3ZzDY58v4cM5G6hXqzq3ntyZ3x7djro11f9VIiOcsZgeKfnczB4GPi3H66eYWXp5jysiAau25vLEF0t4Z9Y6aicncf0JHblySHsa1FEvYomsSHzVqAO0jsB+ShpsZrOB9cCt7j4vwvsXSTjrd+zhqS+XMnb6GpKqGVce056rjm1PE81xIFFS7gJhZnMBDz5NApoCkbxHMBNo6+67zOw04B2gUylZRgAjANLS0iIYQSR+bM7J45lJy/j3d6txnIsHpXHt0I40UwslibJwziBK3hsoBDa5e2GE8uDu2SUef2Rmz5hZqrtnhdh2FDAKoH///r7/epFEtj13L89NWcYr36ykoMg5r19rrjuhI60bVc7pQSX+hHMPYlU0guxjZi0IFB03s4EE+lhsjeYxReJJdl4BL3y1gpemriB3byFn9mnFjcM6kZ6aEutoUsWUuUCY2VR3H2JmOQQuMZVs4+ruXr+M+3kDOB5INbO1wEggObiT54BzCcxWVwjsAS50d50dSKWXm1/Iy9+sZNSU5ezcU8BpPVtw04md6dy8XqyjSRVV5gLh7kOCvw/pX6u7X3SQ9U8TaAYrUiXkFRTx+rRVPDt5GVtz93JC12bcclJnerRqEOtoUsWV5wzilgOtd/dHDz2OSNWxt7CYMdPX8PSXS9iUnc/RHZtwy0ld6Ne2UayjiQDluwex78yhCzAAeC/4/JfA95EMJVKZFRYVM37WOp78Yglrt++hf9tGPH7BEQzu0CTW0UR+ojyXmO4GMLMpQF93zwk+/yvwYVTSiVQixcXO+3PW8/jnS1iRlUvPVg3425k9OK5zU/YbtkYkLoTTzLU5sLfE873BZSISgrvz6bxNPDZxMYs25dCleT2ev7QfJ3drrsIgcS2cAvEq8L2ZTQg+PxN4JXKRRCoHd2fy4i08+tli5q7bSfvUFJ686AiG92xJtWoqDBL/wukH8Xcz+xg4JrjocnefFdlYIontm2VZPPLZYmas2k7rRrV56NxenHVEK6onaeoUSRzhjsW0IvjaWkA9MzvW3adELpZIYpqxajuPTlzE10u30rx+Tf52Zg/O79+GGtVVGCTxhDMW05XAjQQG6MsAjgS+BU6IbDSRxJG5biePTlzMlws30ySlBv87vBsXD0qjVrLmZJDEFc4ZxI0EmrlOc/ehZtYVuC+ysUQSw+JNOTw2cTEfZ26kQe1kbj+lC78ZnE6K5mSQSiCcf8V57p5nZphZTXdfaGZdIp5MJI6tyMrlic8X8+7s9aTUqM6Nwzrxu2PaUb+W5mSQyiOcArHWzBoSGIZ7opltB6I6gJ9IvFi7fTdPfbGUt2euJTnJuOrYDlx1bHsapdSIdTSRiCtXgbBAo+0b3H0H8FczmwQ0AD6JRjiReLEpO49/TlrKG9+vxjAuG9yWa47vQLN6mpNBKq9yFYjgENwfAT2Dz/8TlVQicWLrrnye+88yXv12FUXFzvkD2nDd0I4c1rB2rKOJRF04l5hmmtkAd/9vxNOIxImduwv4v6+W89LXK8grKOLMI1px07DOpDXRZD1SdYRTIAYBF5vZKiCXwLwQ7u69IppMJAZ25Rfyr6krGPXVcnLyCjm9V0tuPrETHZtpTgapesIpEL+IeAqRGNuzt4jXpq3k2cnL2L67gBMPb84tJ3Wm22FlmgdLpFKKuylHRSpSfmERb36/hqcnLWVLTj7HdErlDyd3oU+bhrGOJhJz6s0jVVJBUTHjZqzlyS+WsH5nHgPTG/P0RUcwqL3mZBDZRwVCqpSiYue92et4/PMlrNq6m95tGvLAub0Y0jFVQ2+L7EcFQqqE4mLnk3kbeXTiYpZu3sXhLevzwmX9GXZ4MxUGkVKEM1ifARcD7d39HjNLA1q4u6Ydlbj1z0lLeWTiYjo0TeGfv+7LqT1aaE4GkYMI5wziGaCYwOit9wA5wDgCA/iJxKVNOXk0rJPMZzcfR5IKg0iZhNUPwt37mtksAHffbmYaiEbiXpKZioNIOYQzi0mBmSUBDmBmTQmcUYiISCUSToF4EpgANDOzvwNT0XwQIiKVTjijuU4BZgDDCAyzcaa7L4hCNhERiaGwRnN1957AwihlEhGROBDOJaaZZqYWSyIilZxGcxURkZA0mquIiIRU7ktMwdFcGwK/DP401AivEs+y8wqYtHALLRpoelCR8ih3gTCzG4HRQLPgz+tmdn2kg4lEyl/fm8fG7DzuPbNHrKOIJJRwLjH9jkBv6lwAM3sA+BZ4KpLBRCLho7kbGD9zHTcM60TftEaxjiOSUMJpxWRAUYnnRcFlInFlU3Yed06YS+/WDbj+hI6xjiOScMI5g/gX8J2ZTSBQGM4AXoxoKpFDVFzs3Dp2NnkFRTx6QR+Sk8L5LiRStYUz5eijZjYZGEJgPKbL3X1WpIOJHIrXpq3iqyVZ3HtmDzo0rRvrOCIJKZz5IGoBxwPHEBikr7qZLXD3vAhnEwnL0s053PfRAoZ2acolg9JiHUckYYVzielVAnNAPBl8/mvgNeC8SIUSCdfewmJuGpNBSs3qPHBuL80WJ3IIwikQPdy9W4nnk8xsfqQCiRyKJ75YTOa6bJ6/tB/N6qnfg8ihCHcspiP3PTGzQcD0sr7YzF4ys81mlnmQ7QaYWaGZnRtGRqmCpq/cxrOTl3F+/9b8onuLWMcRSXjhFIh+wDdmttLMVhLoAzHAzOaa2ZwyvP5l4JQDbRCckOgB4LMw8kkVlJNXwM1vZdCqUW3+8svusY4jUimEc4npgB/uB+PuU8ws/SCbXY/muZZyuOf9+azbvoe3rhpM3Zrh/LMWkf2F08w1quMumVkr4CxgKAcpEGY2AhgBkJam1ipV1SeZGxk7Yy3XDu1A//TGsY4jUmnEY++hx4E/uvtB57l291Hu3t/d+zdt2rQCokm82ZwT6C3do1V9bhzWOdZxRCqVeDwX7w+8GWyemAqcZmaF7v5ObGNJvHF3bn97Drn5hTx+QR9qVI/H7zsiiSvuCoS7t9v32MxeBj5QcZBQXv9uNZMXbeHuX3WnY7N6sY4jUumU+SuXmd1e4vF5+627rxz7eYNAy6cuZrbWzH5nZleb2dVl3YfIsi27+PuH8zmmUyqXHtk21nFEKqXynEFcCDwYfPwnYGyJdacAd5ZlJ+5+UVkP6O6/Leu2UnUUFBVz85gMaiUn8fB5valWTb2lRaKhPAXCSnkc6rlI1Dz1xRLmrN3JMxf3pXl99ZYWiZby3NXzUh6Hei4SFTNXb+fpSUs5u28rTuvZMtZxRCq18pxB9DazbAJnC7WDjwk+19c4ibrc/EJuHpNBywa1+euv1FtaJNrKXCDcPam0dWZ2c2TiiJTubx/OZ/W23bz5P0dSv1ZyrOOIVHqRajiuAiFRNXH+Jt74fg0jjm3PoPZNYh1HpEqIVIHQTWqJmi05+dwxbg6Ht6zPLSept7RIRYlURzndpJaocHfuGDeHnPxC3riwDzWrl3qlU0QirMwFwsxyCF0IDKgdsUQiJbz53zV8sXAz/zu8G52bq7e0SEUqz01q/e+UCrUyK5d7P5jP0R2bcPlR6bGOI1LllPsSk5nVAjoGny5197zIRhKBwqLA3NLVq5l6S4vESHnGYqpuZg8Ca4FXgFeBNWb2oJmpzaFE1D8nLSNjzQ7+flZPWjbQFUyRWChPK6aHgMZAO3fv5+59gQ5AQ+DhaISTqiljzQ6e/HIJZ/Q5jF/2PizWcUSqrPIUiOHA/7h7zr4F7p4NXAOcFulgUjXt3hvoLd28Xk3uOaNHrOOIVGnlGovJ3X/Wisndi1AzV4mQv3+4gJVbc3n4/N40qK0rlyKxVJ4CMd/MLtt/oZldAiyMXCSpqiYt3Mzo71Zz5ZB2HNUhNdZxRKq88rRiuhYYb2ZXADOCy/oT6ANxVqSDSdWydVc+t709h64t6nHrL7rEOo6IUL5+EOuAQWZ2ArBvKM2P3P2LqCSTKsPd+dP4uWTvKeC13w1Ub2mROFHufhDu/iXwZRSySBU1dvpaPpu/iTtP68rhLevHOo6IBJVnqI2nOMDNaHe/ISKJpEpZvXU3d78/jyPbN+bKIe1jHUdESijPGcT0Eo/vBkZGOItUMYVFxdz8VgbVqhmPnN9HvaVF4kx57kG8su+xmd1U8rlIOJ6fspwZq7bz+AV9aNVQvaVF4k2480Go34Mckrlrd/LYxMUM79WSM/qot7RIPIrUhEEiZbZnbxE3jZlFat2a/O3MHpjp0pJIPAp3Pog6Zpa9bxWBXtZqfiJl8o+PF7BsSy6v/24QDevUiHUcESmF5oOQCvWfxVt45dtVXH50OkM6qbe0SDzTJSapMNtz93Lb2Nl0alaXP57SNdZxROQgIjUntcgBuTt3TpjL9t17+dflA6iVrN7SIvFOZxBSIcbPXMfHmRu55aQudD+sQazjiEgZqEBI1K3ZtpuR781jYHpjRhyr3tIiiSKcOan7A38G2gZfv68VU68IZ5NKoKjY+cNbswF45PzeJKm3tEjCCOcexGjgNmAuUBzZOFLZjJqynO9XbuPh83rTpnGdWMcRkXIIp0Bscff3Ip5EKp3MdTt5dOIiTu3RgnP6top1HBEpp3AKxEgzewH4Asjft9Ddx0cslSS8vIIibh6TQaM6NbjvrJ7qLS2SgMIpEJcDXYFkfrzE5IAKhPzggU8WsmTzLl65YiCNUtRbWiQRhVMgBri75oSUUk1dksW/vl7Jbwa35bjOTWMdR0TCFE4z12/MrFvEk0ilsGP3Xm4dO5sOTVO449TDYx1HRA5BOGcQRwKzzWw5gXsQauYqQKC39F3vZJK1K5//u+xoatdQb2mRRBZOgTgl4imkUng3Yz0fzNnArSd3pmdr9ZYWSXThFIjflLL8nkMJIolt3Y49/O+7mfRr24irj+sQ6zgiEgHh3IPILfFTBJwKpJf1xWb2kpltNrPMUtafYWZzzCzDzKab2ZAwMkoFKi52/vBWBsXFzmPn96F6kkZwEakMyn0G4e6PlHxuZg8Dn5ZjFy8DTwOvlrL+C+A9d3cz6wW8RaBZrcSpF6euYNrybTx4Ti/Smqi3tEhlEYmvenWA1mXd2N2nANsOsH6Xu++buS4FzX8d1xZsyOahTxdxcrfmnNe/zP8MRCQBhDNY31x+/NBOApoS4fsPZnYWcD/QDDj9ANuNAEYApKWlRTKClMG+3tL1aydz/9nqLS1S2YRzk3p4iceFwCZ3L4xQHgDcfQIwwcyOBe4FTixlu1HAKID+/fvrTKOM1mzbzdx1OzmtZ8tD2s8jny1i4cYcXvptf5rUrRmhdCISL8K5B7EqGkFKOdYUM2tvZqnunlVRx63MPp+/iVveyiA7r5CF954S9sxu3yzL4oWpK7h4UBondG0e4ZQiEg/ibj4IM+sILAvepO4L1AS2RmLfVVlhUTEPf7aY5/6zjOSkwKUgD/Oca+eeAm59azbpTVL48+nqLS1SWVX4fBBm9gZwPJBqZmuBkQQG/sPdnwPOAS4zswJgD3BBiZvWEobN2Xlc98Ysvl+xjV8PSqN5vVo89vnisPf3l3cz2ZSTz7hrjqJODU1rLlJZVfh8EO5+0UHWPwA8EO7+5ae+XbaV69+Yxa78Ah49vzdn923Nc/9ZFvb+3pu9nncz1nPziZ3p06ZhBJOKSLzRfBCVVHGx8/yU5Tz06ULSU1MYfeUgurSod0j73LBzD3dNmEufNg25dqh6S4tUdpoPohLaubuAP4zN4PMFmzm9V0seOKcXdWse2qWg4mLn1rGzKShyHrtAvaVFqgLNB1HJzF27k2tGz2BTdh53/6o7lw1uG5H+Cf/6ZiVfL93KfWf1pF1qSgSSiki8C6dAfGNm3dx9fsTTSNjcnX9/v5q735tPat0ajLlqMH3TGkVk34s25vDAJwsZ1rUZFw1sE5F9ikj8C3c+iAwzW4Hmg4gLu/cW8ucJmUyYtY5jOzfl8Qv60DhC03zmFxZx05gM6tWszj/O6aXe0iJViOaDSHBLN+/i96NnsGTzLm45qTPXDe1ItWqR+xB/bOISFmzI5oXL+tO0nnpLi1Qlcd2TWg7s/dnruWPcHGomJ/HqFQM5plNk53/+bvlWnp+yjIsGtuHEbuotLVLVlLlAmNlUdx9iZjkEWi2V/Jrq7l4/4ukkpL2Fxdz30QJe/mYl/do24ulfH0HLBrUjeozsvAJueWs2aY3rcNfpmoJcpCoqc4Fw9yHB34fWmF4Oybode7h29Ewy1uzgyiHt+OOpXUmOQpPTv743jw079zD26qNIOcQmsiKSmMIZi6kmgeEw0ku+3t015WiUTV60mZvGZFBY5Dx7cV9OPcTRWEvz0dwNjJ+5jhtO6Ei/tpFpCSUiiSecr4bvAjuBGZToSS3RU1TsPPHFEp76cgldmtfj2Uv6Ra0vwqbsPO6cMJferRtw/bBOUTmGiCSGcApEa3dXS6YKsnVXPje+mcHUpVmc2681957Rg9o1whui+2D29ZbOKyji0Qv6ROXSlYgkjnA7yvV097kRTyM/MX3lNq779yy2797Lg+f04vwB0e2k9tq0VXy1JIt7z+xBh6Z1o3osEYl/4RSIIcBv1VEuetydF6eu4B8fL6RVo9qM//1RdD+sQVSPuXRzDvd9tIDjuzTlkkGavlVEwisQp0Y8hfwgO6+A28fO4ZN5Gzm5W3MeOq83DWonR/WYewuLuWlMBnVqJPGgekuLSJA6ysWR+euz+f3oGazZvoc/n3Y4Vx7TrkI+rJ/4YjGZ67J57pJ+NKtfK+rHE5HEoI5ycWLs9DXc9U4mDWon8+aIIxmQ3rhCjjt95TaenbyM8/q15pQeLSrkmCKSGNRRLsbyCooY+e48xkxfw1EdmvDEhUdU2JhHOfkF3PxWBq0a1Wbkr7pXyDFFJHGE01GuP3AnP+8op5vU5bQyK5drRs9kwYZsrj+hIzed2JmkCA60dzD3frCAddv38PzabxYAAA1ESURBVNZVgw95QiERqXzC+VQYDdwGzOXHGeWknD7J3MhtY2dTrZrxr98OYGjXZhWe4f3Z67l2aAf6V9DlLBFJLOEUiC3u/l7Ek1QRBUXFPPjJQv7vqxX0bt2Af17cl9aN6sQkS49W9blxWOeYHFtE4l84BWKkmb0AfEGJoTbcXXNSH8TGnXlc9++ZTF+1ncsGt+XPpx9OzerR6RV9IKl1a5JSI4nHL+hDjerqLS0ioYVTIC4HugLJ/HiJyQEViAP4emkWN745i917i3jiwj6c0adVzLKc07cVp/VsQZ0auu8gIqUL5xNigLt3iXiSSqq42Hlm8lIenbiY9k3r8uaIvnRsFtuGYGam4iAiBxXuWEzd3H1+xNNUMttz93LzWxlMXrSFM/ocxn1n9dTcCiKSMML5tDoSyNBYTAeWsWYH146eyZacfO49sweXDErTEBYiklDCKRAa6vsA3J3Xpq3i3g/m06xeLcZePZjebRrGOpaISLlpLKYIys0v5I7xc3l/9npO6NqMR8/vTcM6NWIdS0QkLLogHiFLNuVw9eszWJGVy22/6MI1x3WgWgX2ihYRiTQViAh4N2Mdd4ybS0rNJF6/chBHdUiNdSQRkUOmAnEI8guLuPeD+bw+bTUD0xvz1K+PoLmGyxaRSkIFIkxrtu3m2n/PZM7anVx1bHtu+0UXqmsOZxGpRFQgwvDlwk3cPGY2xe48f2k/ftFd8yiISOWjAlEOhUXFPDpxMc9MXka3lvV59pK+tG2SEutYIiJRoQJRRptz8rjhjVlMW76Niwa2YeQvu1MrueIH2hMRqSgqEGXw3fKtXP/GLLLzCnj4vN6c2691rCOJiESdCsQBuDujpiznwU8Xkda4Dq/+biBdW2jqbRGpGlQgSrFzTwG3jp3NxPmbOK1nCx44pxf1aiXHOpaISIVRgQghc91Ofj96Jut37OEvw7tx+dHpGmhPRKqcCm+4b2YvmdlmM8ssZf3FZjbHzOaa2Tdm1ruisrk7b36/mrOf/YaComLGXDWYK4a0U3EQkSopFmcQLwNPA6+Wsn4FcJy7bzezU4FRwKBoh9qzt4i73slk3My1HNMplccv6EOTujWjfVgRkbhV4QXC3aeYWfoB1n9T4uk0IOpNhpZv2cU1r89k8eYcbjqxE9ef0IkkDbQnIlVcvN+D+B3wcWkrzWwEMAIgLS0trAP8Z/EWrh09k+Qk45XLB3Js56Zh7UdEpLKJ2wJhZkMJFIghpW3j7qMIXIKif//+Hs5x0hrXoW/bRvzj7J4c1rB2WFlFRCqjuCwQZtYLeAE41d23RvNY7VJTePWKgdE8hIhIQoq74UfNLA0YD1zq7otjnUdEpKqq8DMIM3sDOB5INbO1wEggGcDdnwP+AjQBngk2Ly109/4VnVNEpKqLRSumiw6y/krgygqKIyIipYi7S0wiIhIfVCBERCQkFQgREQlJBUJEREJSgRARkZDMPawOyHHHzLYAq8J8eSqQFcE4saL3EV/0PuKL3kdobd095BhDlaZAHAozm14Z+lrofcQXvY/4ovdRfrrEJCIiIalAiIhISCoQAaNiHSBC9D7ii95HfNH7KCfdgxARkZB0BiEiIiGpQIiISEhVvkCY2SlmtsjMlprZHbHOEw4ze8nMNptZZqyzHAoza2Nmk8xsvpnNM7MbY50pHGZWy8y+N7PZwfdxd6wzhcvMksxslpl9EOss4TKzlWY218wyzGx6rPOEy8wamtnbZrbQzBaY2eCoH7Mq34MwsyRgMXASsBb4L3CRu8+PabByMrNjgV3Aq+7eI9Z5wmVmLYGW7j7TzOoBM4AzE/Dvw4AUd99lZsnAVOBGd58W42jlZma3AP2B+u4+PNZ5wmFmK4H+7p7QneTM7BXgK3d/wcxqAHXcfUc0j1nVzyAGAkvdfbm77wXeBM6IcaZyc/cpwLZY5zhU7r7B3WcGH+cAC4BWsU1Vfh6wK/g0OfiTcN/EzKw1cDqB6X8lhsysAXAs8CKAu++NdnEAFYhWwJoSz9eSgB9IlZGZpQNHAN/FNkl4gpdmMoDNwER3T8T38ThwO1Ac6yCHyIHPzGyGmY2IdZgwtQO2AP8KXvJ7wcxSon3Qql4gJA6ZWV1gHHCTu2fHOk843L3I3fsArYGBZpZQl/7MbDiw2d1nxDpLBAxx977AqcC1wUuyiaY60Bd41t2PAHKBqN8zreoFYh3QpsTz1sFlEiPBa/bjgNHuPj7WeQ5V8DLAJOCUWGcpp6OBXwWv378JnGBmr8c2UnjcfV3w92ZgAoFLy4lmLbC2xJno2wQKRlRV9QLxX6CTmbUL3vS5EHgvxpmqrODN3ReBBe7+aKzzhMvMmppZw+Dj2gQaQSyMbarycfc/uXtrd08n8P/iS3e/JMaxys3MUoINHghekjkZSLjWfu6+EVhjZl2Ci4YBUW+8UT3aB4hn7l5oZtcBnwJJwEvuPi/GscrNzN4AjgdSzWwtMNLdX4xtqrAcDVwKzA1evwe4090/imGmcLQEXgm2kqsGvOXuCdtMNME1ByYEvntQHfi3u38S20hhux4YHfwyuxy4PNoHrNLNXEVEpHRV/RKTiIiUQgVCRERCUoEQEZGQVCBERCQkFQgREQlJBUJEREJSgRARkZBUIEQk6szsKTObaWYDYp1Fyk4FQkSiKjjERTPgKiAh55SoqlQgJObM7K9mdmusc0RacAaw30dgP0XB2dAyzez9feM8Bdd5yUH0zKy6mW0pOQOcmf05OLPdnOB+Bu2333nB2e/+YGaH9JlgZrXN7D/BYUYAcPdcAsOPTAaeDG5Xw8ymmFmVHu4n3qlASIWygLD/3R3q66O9v/00BMpVIErJs8fd+wRnC9wGXFtiXS7QIzgoIAQGBvxhROLgtJTDgb7u3gs4kR/nQNm33+7B150KjCxP3hCuAMa7e1GJDE2AOkAOUAiBCW+AL4ALDvF4EkUqEBJRZnZL8JtuppndFFyWboF5v18lMJJmm+C32sVmNhXoUuL1l1hgPucMM3s+OPHOz16/3zHTLTBP72gLzNX7tpnVCa57JzhRzLx9k8WUkqe07Raa2cvBrKPN7EQz+9rMlpjZwAPlBv4BdAgueyjc97efb/n5pFYfEZj5DeAi4I0S61oCWe6eD+DuWe6+fv+dBofCHgFcZ8GR7fZnZr8P/r2uMrPrS8l3MfDufsvuAh4G5gHdSyx/J7i9xCt3149+IvID9APmAilAXQIfCEcA6QRmJTtyv+3qAPWBpcCtwOHA+0BycLtngMv2f32I46YTmDXs6ODzl4Bbg48bB3/XJvDh2yTU/g6wXSHQk8CXqRnBfRuBqWnfCb7mQLkzSxwj3Pe3K/g7CRgLnFJyHdCLwPwAtYAMAiP7fhBcXze4bHHweMftv9/9jrUDaB5i+TkE5lJIJlh0gOr7bVMD2Bji7+a74J/Z08D/lFiXBGyJ9b9b/ZT+ozMIiaQhwAR3z/XAnMzjgWOC61a5+7Tg42OC2+32wIxx++bgGEagePzXAsN9DwPah3h9KGvc/evg49eDWQBuMLPZwDQC38w7lbK/0rZb4e5z3b2YQMH7wgOfbnMJfPgdLHdJ4b6/2sHtNxIYvnpiyZXuPieY5SICZxMl1+0KHnMEgSkrx5jZb0s5zoHcAPzR3QvcfQNQwM+vQKQSKDAl/Q24J/hntoASZxAeuAy114LzNUj80Q0iqSi5ZdjGgFfc/U8/WRiYn/pgr99/3Ho3s+MJXHMf7O67zWwygW/ZP8lzkO3yS+yzuMTzYn78/3Og3D9ZdIDtDvT+9rh7n+Bls08J3IN4cr9t3iNwGed4Amc/Pwh+EE8GJpvZXOA3wMv7H8TM2gNFBObRLrk8Gejt7ouDz1sCWz1wH+EnOfnxzw0z6wOcDQwxs38G183d7zU1gbzS37rEks4gJJK+As40szoWaNp4VnDZ/qYEt6sd/Pb4y+DyL4BzzawZgJk1NrO2ZTx2WvCGLMCvgalAA2B78EO/K3BkKa8t63alKS13DlCvDNuVibvvJvBN/g8hWv+8BNzt7j/5ADazLmbWqcSiPsCq/fdtZk2B54Cng9/2S+oG1Dez9sEb6Pfz8wKFu28HksxsX5F4APiVu6d7YGa63pQ4gwjevM5y94KDvHWJEZ1BSMS4+0wzexn4PrjoBXeftf836eB2Y4DZBL6t/je4fL6Z3QV8FvwgKiDwbXljGQ6/iMCE9C8RmIrxWQLfhq82swXB9aVdwvmkjNuFVFpud58WvKGdCXzs7rcdwvvbd6xZZjaHwOWk10osX0uID20C9yCeskDT2EIC93tGBNftu3SVHFz3GhBqqtcjgNEEbn6nEGilNKqUiJ8ROGMoBuq4++clMm4ys7pm1tjdtwFDgQ/L+NYlBjSjnCS8YAH6wAPNQCXCzOxx4Ft3H1OGbfsCN7v7pWXYdjxwx75LVxJ/dIlJRA6mD4GWUAfl7jOBSVaio1woFphX+R0Vh/imMwgREQlJZxAiIhKSCoSIiISkAiEiIiGpQIiISEgqECIiEpIKhIiIhKQCISIiIf0/75vCwLveCh4AAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "\n", "min_radius = [[rmsd_i, p.radius.min()] for rmsd_i, p in ht2]\n", "arr = np.array(min_radius)\n", "\n", "plt.plot(arr[:, 0], arr[:, 1])\n", "plt.xlabel(r\"order parameter RMSD $\\rho$ ($\\AA$)\")\n", "plt.ylabel(r\"minimum HOLE pore radius $r$ ($\\AA$)\");" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "[1] Richard J. Gowers, Max Linke, Jonathan Barnoud, Tyler J. E. Reddy, Manuel N. Melo, Sean L. Seyler, Jan Domański, David L. Dotson, Sébastien Buchoux, Ian M. Kenney, and Oliver Beckstein.\n", "MDAnalysis: A Python Package for the Rapid Analysis of Molecular Dynamics Simulations.\n", "Proceedings of the 15th Python in Science Conference, pages 98–105, 2016.\n", "00152.\n", "URL: https://conference.scipy.org/proceedings/scipy2016/oliver_beckstein.html, doi:10.25080/Majora-629e541a-00e.\n", "\n", "[2] Naveen Michaud-Agrawal, Elizabeth J. Denning, Thomas B. Woolf, and Oliver Beckstein.\n", "MDAnalysis: A toolkit for the analysis of molecular dynamics simulations.\n", "Journal of Computational Chemistry, 32(10):2319–2327, July 2011.\n", "00778.\n", "URL: http://doi.wiley.com/10.1002/jcc.21787, doi:10.1002/jcc.21787.\n", "\n", "[3] O S Smart, J M Goodfellow, and B A Wallace.\n", "The pore dimensions of gramicidin A.\n", "Biophysical Journal, 65(6):2455–2460, December 1993.\n", "00522.\n", "URL: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1225986/, doi:10.1016/S0006-3495(93)81293-1.\n", "\n", "[4] O. S. Smart, J. G. Neduvelil, X. Wang, B. A. Wallace, and M. S. Sansom.\n", "HOLE: a program for the analysis of the pore dimensions of ion channel structural models.\n", "Journal of Molecular Graphics, 14(6):354–360, 376, December 1996.\n", "00935.\n", "doi:10.1016/s0263-7855(97)00009-x.\n", "\n", "[5] Lukas S. Stelzl, Philip W. Fowler, Mark S. P. Sansom, and Oliver Beckstein.\n", "Flexible gates generate occluded intermediates in the transport cycle of LacY.\n", "Journal of Molecular Biology, 426(3):735–751, February 2014.\n", "00000.\n", "URL: https://asu.pure.elsevier.com/en/publications/flexible-gates-generate-occluded-intermediates-in-the-transport-c, doi:10.1016/j.jmb.2013.10.024." ] } ], "metadata": { "kernelspec": { "display_name": "Python (mda0180)", "language": "python", "name": "mda0180" }, "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.6.7" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": false, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": true } }, "nbformat": 4, "nbformat_minor": 2 }