{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Pluton emplacement test\n",
"\n",
"This Jupyter notebook employs the stochastic sampling algorithm of [Keller, Schoene, and Samperton (2018)](https://doi.org/10.7185/geochemlet.1826) as implemented in [Chron.jl](https://github.com/brenhinkeller/Chron.jl).\n",
"\n",
" \n",
"
If running this notebook as an online Binder notebook and the webpage times out, click the badge at left to relaunch (refreshing will not work). Note that any changes will be lost!
\n", "\n", "Hint: `shift`-`enter` to run a single cell, or from the `Cell` menu select `Run All` to run the whole file. Any code from this notebook can be copied and pasted into the Julia REPL or a `.jl` script.\n", "***\n", "\n", "## Load required Julia packages" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# Load the Chron.jl package\n", "using Chron\n", "using Plots" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "## Estimate the time of rheological lockup \n", "Given that the rheologically critical melt fraction ≈ 36% melt (see, e.g., [Caricchi 2007](https://doi.org/10.1016/j.epsl.2007.09.032) for experimental rheological data), which is also approximately the terminal porosity of a framework of close-packed spheres, a body of magma may be considered effectively \"emplaced\" by the time it cools to where F $\\leq$ 0.36. \n", "\n", "If we estimate the point within a MELTS zircon saturation distribution where melt fraction drops below 36%, we can use this to estimate the time by which a body of magma must have been emplaced.\n", "\n", "#### Input dataset (Try pasting in your own data here!)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "# Age and one-sigma uncertainty.\n", "ages = [101.333, 102.175, 102.074, 102.814, 102.335, 101.212, 101.398, 101.772, 102.456, 102.302, 102.039, 101.579, 102.135, 101.978, 102.325, 102.36, 101.53, 102.409, 101.277, 102.766,]\n", "uncert = [0.137, 0.141, 0.134, 0.149, 0.132, 0.141, 0.129, 0.129, 0.135, 0.142, 0.144, 0.128, 0.15, 0.144, 0.132, 0.13, 0.133, 0.133, 0.148, 0.14,]\n", "age_unit=\"Ma\"\n", "\n", "# Sort by age (just to make rank-order plots prettier)\n", "t = sortperm(ages)\n", "ages = ages[t];\n", "uncert = uncert[t];" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Input MELTS-derived $\\ \\mathcal{\\vec{f}}_{xtal}(t_r)$" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd0AUx+IH8Nndo/cO0kUUKYqCEsGCHU2woKIkscaoSUSNRmPyU58v5lmiJsauefo0mpigxoIVrMEuiNhRBEFReu93u/v74947EQ8BvX7fz193e3N7wwL7vZmdnaF4nicAAADailZ2BQAAAJQJQQgAAFoNQQgAAFoNQQgAAFoNQQgAAFoNQQgAAFoNQQgAAFoNQQgAAFoNQQgAAFoNQQgAAFpN04Jw3759ly5dUnYtVAjLssqugobDEZY3HGG54nme4zhl10LJNC0I//7776SkJGXXQoVUVVUpuwoaDkdY3nCE5YrjuNraWmXXQsk0LQgBAABaBEEIAABaDUEIAABaDUEIAABaDUEIAABaDUEIAABaDUEIAABaDUEIAABaTaDsCoAaYHmSmM+ffs6ffs6lFPKSWSg8TamhrvQQV8rXglJm/QAA3gGCUMMxDPOOe8ir5qdf5mhCutpQ/wxgfCwoyR5TS/mEHP7b65weQ/7Rifa1RBwCgPpBEGo4AwODd3k7xxNDAbU7lNGV1onexYbqYkPN9iMijlSzpI4jUosBAKgyBKGGo6h3aqXRFDHWabqYgCYmiEAAUE84ewEAgFZDEAIAgFZDEAIAgFZDEAIAgFZDEAIAgFZDEAIAgFbD7ROqpaKiQigUmpubS257qKqqqq2tNTU1ZRimqqqqsrJSUpiiKGtra0JIfn6+paUlwzBCobCkpOT13VpbW1MUVVFRkZaWxjBM69atjYyMXi92586dc+fO2draRkZGyuonSkhIoCiqe/fustohAIBsoUWoWiIiIiwtLS9fvizZ0rt3b0tLy3v37hFCVqxY4eLi0uF/goKCCCEikcjW1jYtLY0QcuPGDfFL3t7etra2kpJlZWXbtm1zdXWdPn36p59+6ujoeODAgQYf/fDhwz59+kjN0cb88ccfa9aseXOZ2NjYY8eONX+fAAAKhhahygkICNi5c2dwcDAhJDU1tbKysv7sMOHh4TExMY29Nygo6MWLF4SQ27dvd+jQQfyYEFJWVvbZZ58lJCSIs7O8vLy8vLzBe69du/bee+8tWLCg+VV98uRJRkZGg43l5eUmJibN3wkAgHKhRahyRo8effjw4aqqKkLItm3bxo8f/+77LCgoEAqF7u7u4qcmJiatWrWqX+D333+fPXv2mTNnPDw8/vOf/7y+h3Hjxs2fP1/8eMKECXPnzr169erKlSv/+OMPDw+PYcOGEULmzZvn5OQUEBDg6Oj4hrQGAFApaBFKkVtNbhfxCvggEx0SZNtwCjRTU9M+ffocPHgwMjLyjz/+uHr16j/+8Q/Jq7du3Zo7d674sbirszkf5ObmFhoaGhgYOHz48O7duw8cONDU1LR+gQ8//LC6uvro0aN//fWX1D38/PPPnTt37tatW15eXmJi4rVr1wwNDefOnZuRkbFlyxZxmU8++WTFihUURd27dy80NHTgwIFmZmbNPBQAAMqCIJTiQg63+QHXdLl3Zq5L7e0rZXWICRMmrFq1ytTUtGPHjg4ODvVfMjAwkDTmbGxsmvlBNE3Hx8fv27fvxIkT8+bNmzZt2oEDB3r27Nn8qlpYWOzatWvkyJEsy54/f97Q0PD1Mk5OTr/++uvjx4/r6uoIIQ8ePBD3xAIAqDIEoRQj3OkR7srsNO7bt++UKVOWL18+a9asBi95enp++eWXb7FPgUAwZsyYMWPGcBw3ZcqUb7/99sKFCy3ag5eXFyHExcWlbdu2r78qEomCg4ODg4NDQ0Npmt6zZ8/rlyEBAFQQrhGqIpqmp06dWl5eHh4eLo+dBwUFlZaWtuhdPM9/8sknY8eOtba2Xrx4sXijjo6OSCQSP05LS8vLy9u0adPo0aMHDRpUXFws22oDAMgJWoQqav78+ZLBKfWlpaWtW7dO8nTq1Kk0TRNCdu/ebWtrK944YMCABu9KTU397LPPoqKiPDw8nj59umTJkqlTp7aoPqtWrcrJydm7d29JSUnnzp179uw5YMAAX1/fjRs37tq1y8HBoUuXLjU1Ndu2bWvXrt2aNWvefUFgAADFYCTf7jXD8ePHrays1PfSVF1dXYcOHZycnOpvrKmpCQ0NNTY2rq2tFYlExfX07dtXIBCUlZXV1dVJNnp5ednZ2bEsGxYWJt6DsbGxnp7e5cuXT506lZeXN3369ClTpjRYqpBlWRsbm44dO75eq5qamlOnTi1dutTKysrIyOi9995LTk7u1q2bh4eHg4PD3bt3a2tre/bs2aNHD/HonsmTJ3fs2DEwMNDS0rKmpsbNzc3T07P+DnlCNGYx+7q6Oj09PWXXQpPhCMsVz/Msy+roNGPdUc1F8bwihkcqzIwZMzw9PaOjo5VdEXiT4lpioSlnNtw3KW84wnLFsmxdXV39m5W1ELpGoaGjR4+mpKTU39KtW7fevXvL8CMKa3maosx0ZbhLAIC3hCCEhqytrd3c3OpvsbCwkO1HGAmoH26x/wrEdUQAUD4EITQUFBQk74usNvpky33us/a0k5HGXCsEAHWF2ydACQQ0+cybXpSkiFkLAADeDEEIyvGVH3Mki3tQolFjtQBAHck3CMUzR79BbW0txzVsFohHMb1euK6uTnL7Nqg7M10yy5dZkoxGIQAombyC8MaNG15eXk5OTq6urmfPnn29QHl5+dChQ+3s7CwsLL7//nvxxkePHnl5eRkZGZmbmwcEBFy7dk28va6ubuzYsTY2NlZWVnPmzNGwWz601kxf+sxz7pZC5jcHAGiMXIKQ5/mxY8d+/vnnRUVFq1evjoqKer2Ft3Tp0tra2vz8/Dt37qxfvz4hIYEQYm1tfeDAgZqamoqKig8++ODjjz8WF16/fn1aWlpOTs7jx48PHTr0+qKy0BiWZV/feCmXL5XS5FY0IwGZ0wGNQgBQMrkEYWJiYnZ29meffUYIGTlypLGx8YkTJxqU2blz55w5c3R0dJydnT/66KOdO3cSQiwsLNq3b08IoWk6PDz8+fPn4sbfzp07o6OjDQwMrK2tP/nkE3FhaI7Xe6ev5PET/mZNVeMevune9JU8/mYhGoUAoDRyCcK0tLQ2bdpI5uxp3759Wlpa/QJVVVUvXrwQZ564wOPHjyWv7t+/f8OGDdOnT//Xv/4lngYsLS1NUtjb27t+4QY4jquqqpJMNlZWVibbH00D7HzETfCkVeSuBX2GzOtA/+MGGoUAoDRyuY+wtLTUyMhI8tTExKTBWgQlJSWEEEkZY2PjoqIiyatnzpwpKioqKioS38ddV1dXVVVlbGwstXADt2/f3r59+/Lly8VPjYyMrly5os3Lw1ZUVNR/WsuSvel6CQPqystVpREW5Uh+SNFPyKrxt1DLOGxwhEHmcITlSjw4UYPHIerr6zc5k6pcgtDa2rp+U6y4uFiyMIKkAEVRZWVl4qgrKSmxs7OTvLphwwZCyJ07dzp37jxw4EA7OztTU1PJskENCjfQsWPHkSNHYq7R+urP03gyg/O34trbGSuxPg2YEPK1P/fTQ/2/+qnrRDOYCVPecITlB3ONEjl1jXp7e6emplZWVhJCOI5LSUnx9vauX0BXV7dNmzY3btwQP01OTpb0fEq0bduW53lx48/Hxyc5OfkNhaGZ9jzmP26jcjePTm5HX8vn7xSrSiMVALSKvIIwICDg22+/zcnJWbZsmYmJiXjK5tjY2GnTponLTJkyZcmSJenp6WfPnv3zzz8nT55MCDlz5syxY8cyMzPv3LkzdepUT09P8WLoU6dOXbly5YMHD65du7Z169YpU6bIo9oar1xITj/nhriqXBDqM2SGD708RS27RgFA3clrrtE9e/bMnDkzODi4Xbt2hw8fFi8eW9+sWbOKiooGDRpkbGz8yy+/+Pn5EUIYhvn555/T09MNDAy6det28uRJ8fqu48aNy87OjoiI0NHRWbZsWWhoqJyqrdkOZ3I97WlLlVz/6Atv2uNP4aNS2tNMRcbxAIC2wHqEGq7+Wm5D4tjRramPVK9rVGzxDfZ5FdnaXc2uFGK1PHnDEZYrXCMkmGtUe5TUkb9zuHDV6xeVmOHD/JXBZVVo1DczAFB9qntaBNk6lMn1aUWbNjGKWJks9cgn7ejVt3GlEAAUCkGoLf5M50a3VvXLb7N8md1pXIkKTP8GANoDQagVimvJ5Vz+AxdV/3U7GJIPXOh/p6JRCACKo+pnRpCJfRncQCfaSF5jhGVpli+9/i4nQhQCgKIgCLVCTAYXqfL9omKdrCgXY3IwE0kIAAqCINR8RbXkej4/yEltftezfOmf7yIIAUBB1ObkCG/t2FOutwNtoA79omJDXennleRaPu6jAABFQBBqviNZ/Acu6tEvKsZQ5HNveh0ahQCgEAhCDSfkSHw2N9hZzX7Rk9vRx55yLxouKgwAIHtqdn6ElrpcQHuYUg6Gyq5HC5npkigPevN9VtkVAQDNhyDUcMez6XCVv31Qqi+86V9SuTr0jwKAnKnlKRKa78QLRr0uEEq0N6e8zan9GUhCAJAvBKEmSy3lq0SUv5VaBiEh5AtvesM9BCEAyBeCUJPFZvGDWrHqGoOEDHGls6tIUgHuowAAOUIQarKjWVxYKzVuUTEU+bQdvfm+Gv8IAKD6EIQaq7SO3Cjge9mpd4pM8aL3P+EKa5VdDwDQXAhCjXXsKdfLgTZg1Ltf0VqffOBM73io3nEOAKoMQaixzjznBzip7/XBl77wpjfe4zj1DnQAUF0IQo11IZfvYa8JQRhkS1nrk6NP0SgEALlAEGqmwlryvJL3tdCEICSERPtg6lEAkBcEoWa6lMsF2VKMhuQgiWxN3y7i7xSjexQAZA9BqJku5vIhdprzy9WlyRQv3EcBAHKhOedKqO9iLh9ipyntQUIIIdPaM3sec6V1yq4HAGgcBKEGquNISiEfZKtRQehgSAY64T4KAJA9BKEGSirg25hSJjrKroesRfvQ63EfBQDIGoJQA13O5YM1q19UrJstZaZLTjxDEgKALCEINdCFXL67JgYhIWSmL/3zHazWCwCyhCDUQJdzuRCNuJX+dWNa0/dKyO0iNAoBQGYQhJomrYxnaMrZSDODUIcmU7ywSCEAyBKCUNNc1JSZ1RozrT29NwPrUQCAzCAINY3m3UHYgI0+Ge5Gb32ARiEAyAaCUNNczNHwICSEzPKlN9zjhIhCAJAFBKFGKakjTyt5P02Za7sxvhZUOzOyLwNJCAAygCDUKNfy+QBrSqAFv9VZvvSaOwhCAJABLThlapOreZo2s1pj3nemC2vJlTzcRwEA7wpBqFGu5HHdtCMIaYp84U2vxSKFAPDOEISagyfkWj7fxUYrgpAQMrkdffIZ96wSjUIAeCcIQs2RVsobCahWhtoShCY6JMqD3oL7KADg3SAINcfVfG25QCgx04fe+oCrweSjAPAOBPLb9fPnz8+cOWNlZdWvXz8dHSlrAtXU1Jw8ebKqqqp///7W1tbijcXFxZcuXSopKfH39/fx8RFvLCgouHnzpuSN/v7+kvIgcTWPD9KaflExTzOqsxX1x2NuQlt8pQOAtySv08fly5d9fX1PnTq1ZMmS/v37i0SiBgUqKiqCgoLWrl175MgRb2/v1NRUQsi9e/dcXV3XrVsXFxfXq1evb7/9VrK3UaNGrfif9PR0OVVbrV3J49/TshYhISTah8GQGQB4F/JqES5evPjrr7/++uuv6+rqOnbsePjw4YiIiPoFdu7caWpqGh8fT9P07Nmzly1btmPHDkdHx7S0NFtbW0JISkpKp06dZs+eLW78+fj4xMfHy6m2GqCGJfdL+E5WWheEYU7Ul1fIhRy+u0bPsAoA8iOXFmFNTU18fPyoUaMIIbq6ukOHDj1y5EiDMrGxsRERETRNE0JGjhwZGxtLCDEzMxOnICGkVatWPM/X1dWJn1ZVVZ06dSopKUmyBepLLuS9zCkDOXZ1qyiaItE+9E+4uR4A3pZcTpzPnz/ned7R0VH81NHRMTk5uUGZ7OxsSQEnJ6eioqLq6moDAwNJgSVLlgwaNKhVq1bipzU1NWvWrHn48CHDMEePHm3durXUjy4sLMzNzRUI/vtzCQSCjz76SOoVSg1z8QUJsiZCobDBdqFQ+PpGDTPWnSxJJg+KOA8TJXy6Nhxh5cIRliuWZYVCoeScqXkYhhG3uN5ALj88y7KEEMlnMwzz+jVClmUZhpEUIITUL7N58+YjR44kJCSInw4ePDg8PJwQwnHc+PHj58yZc+DAAakfXVlZmZubm5SUJH5qYGAgaXdqtqt59GAnnmUb3lTHsqz416HB9CgywYNae5f6sYsS2oXacISVC0dYrtj/UXZF5KU553+5BKGDgwMhJD8/X9yey8vLE29pUCYvL0/8ODc318TExMTkv9/nd+zYsXTp0rNnz0qajJLIpGk6Kirq888/b+yjXVxc+vbtGx0dLdMfSA1cLxQt7cro6ze8TiYUCvX19ZVSJUWa3ZF47xMu6apvpafoj9aSI6xEOMJyxbIsTdNafoTl0lQyNjYODAyMi4sTP42LiwsNDSWEcBxXXFzM8zwhpFevXpLBL/Hx8b169RI/jomJ+b//+7+TJ096eHhI3XlycrKTk5M8qq2+8qpJhZBvY6a9o0XsDMgQV3rLfVwpBIAWk1e/8DfffDN16tTCwsI7d+7k5uZGRUURQlJTU729vYuLi83NzadMmdKxY8eZM2fa29uvWLFCPJrm9u3bH374YY8ePdasWSPez9y5c9u0afPll1/yPO/i4nL//v09e/Y01i+qtS7ncV1tKe2NQUIIIXP86P7HRLP9aH1G2VUBALUiryCMiIiwtbU9duyYt7f3ypUrjYyMCCEODg6bN282NDQkhNjb2yclJf36669VVVXnzp3z9/cnhFhZWW3cuLH+foyNjQkhEyZMiIuLy83Nbd++/e3bt93d3eVUbTV1LZ8PstH866Bv5mtBdbSi9jzmJuLmegBoCUrcUdnA7NmzR44cGRwcrPgKvaMZM2Z4enpq2zXCfsdEX3VgwpyktAnLy8slF181Xnw2P+cqmxIhUGTjWKuOsFLgCMsVy7J1dXX1R+xrIenfnU+cOBESEtKhQ4f169eXlpYquE7QIhxPEgu0aNGJN+jnSNGExD3DehQA0ALSg/D27dvx8fFeXl6zZs2yt7ePjIw8deqUgmsGzXS/hLc1oBQ/WlIFUYTM7UCvvKWxA8EBQB6kByHDMP369YuJicnMzFy8ePHVq1f79+/v7e29YsWKoqIiBVcR3uxqvtbNtf0Go1vTaWXkRgEahQDQXE0MK3B0dPz666/T0tJmz559//79+fPnu7q6zpw5Mzc3VzH1gyZdy+O7Igj/R0CTGT70qtu4jwIAmquJICwuLv7555/9/f1//PFHHx+ftWvXfvHFF7t37w4MDCwpKVFMFeHNtHAZwjf71IuOe8all6NRCADN0mgQXrhwYdy4ca1atZo3b16HDh3Onz9/586d6Ojo5cuX37lzp7i4GGtBqIIqEXlUyne0RBC+ZKJDJnvRWJsJAJpJ+n2EwcHBly9fdnV1Xbhw4SeffGJnZ1f/VQcHBxcXF4wmVQU3CnhfS0oPt5C/6ktfxnufcGEnBmOIAKBJ0oPQ09Pz22+/HTRokGSSzwauXLmip4dzjPJdy8cFQinsDMhQV3rrA+6bjri5HgCaIP008cknn/To0aNBCpaWlkpuojA1NUUQqoIreRgyKt1sP3r9Xa4O/aMA0BTpQRgZGXn37t0GG+/du9e/f3/5Vwla4BpGyjTC14LysSB/piMJAaAJLeg4qq2t1fKlOlRNXjUpF/IepghC6b70ZVbdQhACQBNeuUaYm5v77NkzQohQKExNTa3f+VlUVLRlyxbMdq1SruRxXW20fdGJNwhzpuZcJedf8L0ccJAAoFGvBOHvv/8+e/Zs8eNJkyY1KGpoaPjvf/9bQfWCZkC/6JtRhMz0pX+6w/VywLBaAGjUK0E4YsQIPz8/QkhkZOTixYu9vb0lL1lZWXl4eJiamiq6gtC4q/n8bF+c4t9kXBt6YaIwtZRup8WrFgPAm70ShC4uLi4uLoSQnTt3BgcHW1lZKalW0DSekCQsOtEUAwGZ4kWvv8utC8Y3BgCQTvpgmfDwcKSgintYylvoUtYYvdSUL7yZ3x9zxbXKrgcAqKqXLcL4+PjVq1dPmjQpMjIyKiqquLhY6htOnDihqLrBm9wo4DtboznYNAdDMtSV3oyb6wGgES+DkGXZ2tpakUhECKmrq6utxVdolZZcyHeyQhA2y7wOdO+joi99aX30jwLAa14GYVhYWFhYmPjx/v37lVQfaK7kQn62L5o4zeJlTgVYU7vTuMntcMQAoCGcF9TVzUK+E7pGm21uB2bVLY7D0kwA8BrpQbhv377Dhw+LHxcVFY0ePdrNzW3UqFH5+fkKrBs0KquCZyhib6DseqiPXg6UhR6JzcJEMwDQkPQgnDFjhiTz5s+ff/jw4e7du1+5cmXixIkKrBs0KrkQI2VabI4fvRIzrgHAa6QEYVlZ2YsXLwIDAwkhQqEwJiZm/vz5u3fv3rdv37FjxwoLCxVeSWgII2XeQoQbnVdDLuWiexQAXiElCCsrKwkhZmZmhJArV66UlpYOGTKEENKpUyee57OyshRcRXjdzULijyBsIZoiX/qiUQgADUkJQltbWx0dnVu3bhFCYmJibGxsOnbsSAgpKioihGABClVwo4APQNdoy01sS1/N5+6XoFEIAC9JCUKGYUaNGjV58uSPPvrol19++fDDD2maJoQkJSXp6OiI52ADJSqsJWVC3t0EQdhi+gyZ6sX8fAeNQgB4Sfpgmc2bN0dGRmZkZEyePHnJkiXijXFxce+//76RkZECqwdSJBfw/lZYfektRfvQezO4nGpl1wMAVIZA6lYTE5P169c32Pjzzz/Lvz7QNIyUeReWemR0a3rTPfafAZhmBgAIwQ316ghB+I6+6kBvus9VipRdDwBQDdJbhEKhcMeOHQcPHszMzKyufqUX6fHjxwqpGDQquZDH/NHvorUJ1dOB3vGQ+8IbhxEAGgnCTz/9dOfOnf7+/p07dzYwwPwlKqRKRJ5W8F7maBG+k6870KPPsNPa0wwOJIDWkxKEQqHwjz/+WLRo0T//+U/FVwjeLKWI9zKndNCSeTddbChHI7I3nRvjgUMJoO2knAWKiopqa2uHDRum+NpAk5KxDKGMzO/ILE/BLNwAIC0IbWxsXF1dHz58qPjaQJNuFvH+lghCGRjsTFEUOZqFKATQdlKCkKbpLVu2LF68OCkpSfEVgjdLQotQRihC5nek/3WTVXZFAEDJpA+WWbRoUXZ2dmBgoLW1tampaf2XMGpUieo48qCE74AWoYyMdKcXJXHnX/C9HHBIAbSX9CDs2bOnv7+/gqsCTbpdxHuYUobSf2nQYgxF5nWgl6WwvRxwTAG0l/T//5UrVyq4HtAcNwr4zriVXqbGedJLkrnEAj4QHc4A2gpjx9UJ5pSROR2azPajl93ENNwA2qvRIDx06FD37t0tLS2dnJzEW3744Yc1a9YoqmIgRRJWX5KDKV701Xz+dhGGjwJoKelBuHPnzmHDhhkYGAwdOlSy0d7eftmyZSzbrFF2lZWVCxcuDA8P/+qrr8QLGb7ut99+Gz58+Pjx4yXDU1+8eLFy5crRo0ePGTPml19+EYleTgd56NChkSNHfvjhh3///XdzfzjNIuLI3WK+I1qEsqbPkJk+9NIUNAoBtJSUIOR5/ttvv505c2Z8fPyECRMk20NCQvLy8rKzs5uz30mTJt24cWPGjBk5OTnDhw9/vcBvv/329ddfT5w4sXPnzn379n369Ckh5ODBg6mpqREREcOGDVu5cuU333wjLnzy5MlPPvlk9OjRffv2DQ8Pv3Pnzlv8qOruXgnvYkyZ6Ci7Hproc2/63HMutRSNQgCtxL/mxYsXhJCbN2/yPH/u3DlHR0fx9oqKCkLI9evXX39LA5mZmbq6uvn5+TzP19bWmpmZvf6ugICA//znP+LHo0ePXrhwYYMCBw4ccHV1FT8OCwtbuXKl+PH06dOnTZvW2EdHR0evXbu2yRqqo+2p7EdnRS19V1lZmTwqo3m+T2bHnWvx4eVxhOUPR1iuRCJRVVWVsmuhZFJahLq6uoSQBotOEEKePHlCCDEzM2syXJOTkz09Pa2trcV769Kly/Xr1+sXEIlEN2/e7N69u/hp9+7dGxQQf5yjo6P48fXr10NCQt5QWBtgyKhcRfvQx55yj9AoBNA+Um6fsLS09Pb23rhxY1BQEPW/hdB5nl+xYoWTk1ObNm2a3GlOTo6lpaXkqbW1dU5OTv0C+fn5LMtKylhZWTUokJ6e/t133+3bt48QIhQKi4qK3lC4vgcPHpw4ceLAgQP//fEEgh07djSYE0BNXc/TGWzPVlS07FJWZWUlhdXsm4Em5NM2zPdJ7IauLVuoEEdY3nCE5Ypl2bq6umYO/lBH+vr6AkETNwpLf3nFihVDhw7Nzs729vaurq5et27d3r17ExISdu3a1Zy/SCMjo9raWsnTqqoqIyOjBgUIIZIy1dXVxsbGklefPXs2YMCA7777rk+fPoQQgUCgp6fXWOEGXFxc2rZtK5kx3NDQ0MHBQQP+i1ie3CsVdnPSM9Zt2Rt5nn/D4YL65nYibfcKC4meq3EL/mBwhOUNR1iuxEGo5cvtSR81+sEHHxw6dOj58+fr168vKiqaMWNGRkbGrl27Pv744+bs1NnZOTMzk+f/28v05MkTFxeX+gVMTU3NzMwyMjIkBZydncWPc3Jy+vXrN23atOnTp4u3UBTl5OQk7phtUPh1hoaG7dq16/c/wcHBGpCChJCHpbydAWXewhSEFrHQI1O8cE8hgNZp9D7CDz744P79+0+fPk1MTHz48GFWVlYzU5AQEhISIhAIjhw5QghJTEzMyMgYPHgwIeTevXv79+8Xl4mMjNy+fcmS2GgAACAASURBVDshpLKy8s8//4yMjCSE5Ofn9+vXb+zYsV999VX9HYoL8zwvFAp37dolLqxVEnEHoULM9mP2ZnCZFbhSCKBFmphZxsnJKSAgwNPTs0XtKoFAsGHDhokTJ/bs2TMsLGzNmjXiq3Tnz59fsmSJuMzChQsTEhK6dOni4+PTqVOn8PBwQsjGjRvv3bu3evVqS0tLS0tLW1tbceHZs2dnZWX5+/v7+PhYWFiMHTv2bX5WdXYDi04ohJUe+dSL/uEWGoUAWoSSdGBKPH/+fOPGjefPn3/27BnP8/b29iEhIZ9//rmHh0eLdl1aWvro0SM3Nzfx8FFCSG1tbW1trWToCsuyt27dMjExkQzAqa6urqmpqb8TCwsL8QOO4+7cuaOrq+vl5fWGD50xY4anp2d0dHSLqqr6eh0RLezE9HNscRaWl5ebmJjIo0qaqqCGeO0V3owQOBk162jjCMsbjrBc4RoheX2wzN9//z1kyJDS0lJTU1M3Nzeaph89enT16tUNGzbs3r175MiRzd+1mZlZYGBg/S16enp6enqSpwzDdOrUqX4BAwODxn4fNE136NCh+Z+uSXhCUor4TmgRKoS1PpnYlv7hFre2G6PsugCAIrzSNVpaWhoZGWlqanr06NHi4uKUlJTk5OSCgoK///67TZs2EyZMyMzMVFZFtVlaKW+uS1npNV0SZOKrDsxvadzzKlwpBNAKrwTh77//XlhYePz48cGDB9P0f1+iKKpHjx7x8fECgWDr1q3KqKS2u1GIkTIKZWdAJrTFlUIAbfFKEJ47d27AgAE+Pj6vl3NwcBgzZszZs2cVVTF46UYBVl9StLkdmN2PuKeVaBQCaL5XgvDRo0cBAQGNFQ0ICHj48KH8qwQNYfUlxbM3IFPb09/dQKMQQPO9EoRlZWXm5uaNFbWwsCgtLZV/leAVvHg9XgShws3rwBzO4u6XoFEIoOFeCUKRSPSG+wVpmq6/QCAoxpNy3kBA2Wv12GblMNMlc/yYfyShUQig4RrePnHkyJHc3FypRdEvqhRYdEKJZvjQbWNEV/P4IFv8CgA0VsMgPHPmzJkzZ5RSFZAquZDvbK3sSmgrfYb8Xyd6URJ7clATs9cDgPp65d87JSWF49ARpFqSC/kpXk3MhAfyM6ktvfo2d+Y536cVGoUAmumVIGzOorugYJhlVLl0aPKvQHrOVTZpmIDG7wFAE6GpodKyK3kRT5ybN+klyMkod9pIQH5/jM4SAM2EIFRpmFNGRawKYr69zlVh0DSAJkIQqrQbBQRBqAres6W62lLr7qJRCKCBEIQqLbkQk6upih+60qtus7nVyq4HAMgaglClIQhVR2sTKsqD/tdNVtkVAQAZkx6EP//8c3Z2toKrAg3k15CyOr61KYJQVSzqxOx5zKWVYdI1AI0iPQhXrlzp6uo6bNiwY8eO4c5CZUkq4DtbNz7lHSictT6Z5cssSMR/BIBGkR6Ed+7c2bhx4+PHj99//30XF5f58+dnZWUpuGaAydVU0Gxf+lIufyUPjUIAzSE9CM3NzadMmXL79u3ExMT3339/7dq17u7u/fv337t3L8viGomC3Czk/RGEKsZAQBZ0ov8vEf8FAJqjicEyAQEBW7ZsycrKio6OPnXqVGRkpKen59q1a6urMXhO7nAToWqa1JZ+XkVOPEOjEEBDNBGELMvGxsZOnDhx/fr1FhYWM2fODAwMnDNnTnBwcE1NjWKqqJ2Ka0lBDd/WDEGocgQ0WRpIz7vGcohCAI3QaBA+ffp08eLFbm5uQ4YMyc3N/eWXX7Kzs9esWRMTE5OUlHT//v2TJ08qsqLa5kYh729FYXJL1TTcjTbRIbvSMGoGQBNIX1wmKipq7969enp6UVFR06ZNCwwMrP9qhw4d3N3d8/PzFVJDLXU9n++CflEVtrYbEx4nGuJCY30mAHUn/b/4xYsXq1evHj9+vLm5udQC+/bts7Ozk2fFtF1SAT/CDUGougKsqWGu9MIkdpmfsqsCAO9GehDu3LnTzs5OX1+//saamprnz5+3bt2aEOLj46OI2mmxpAJ+aSDm/VFpy7owPvtFI1rRvU2UXRUAeAfST7VBQUE3btxosDE5OdnDw0P+VQJSVEuKankPzCmj2sx0yb8C6dlJAhajZgDUWQvaHEKhUEdHR35VAYmkAr4zRsqog3GetKkO2foAo2YA1NgrXaM1NTXiGwQ5jisvLy8uLpa8VFJSEhMT4+TkpOgKaqWkAtxBqB4oQlZ2Foafp4e70fYGyq4NALyVV4Jw06ZNs2fPFj8OCwtrUJSiqFWrVimoXtotsYAf5Y4gVA/tzfjJ7eiZl9k/+zDKrgsAvI1XgrBfv35btmwhhMydO3fatGn1rwhaWlr6+vp6eXkpuoJaKamAX94FI2XUxqJOTKcDokOZ3FBX/NYA1M8rQejn5+fn50cIqaurGz58uKOjo5JqpdWKakkJRsqoFT2GbAphPjrH9nKgzXWVXRsAaCHpX2CnT5+OFFSW6/lYfUn99HKgBjlRCzAZN4AaetkivHLlyq5du4YOHTpgwIBvvvmmrKxM6hs2bNigqLppKcy1raZ+6Mr47hd96MEH2+HXB6BOXgbh06dPjxw54uPjM2DAgPj4+MZmUEMQylsSRsqoJws98nM3+tME9sZwgR7GzQCoj5dBOGrUqFGjRokfJyYmKqk+QBLz+RUYKaOeRrrTv6Xxy1LYxZ2RhABqAydc1VJQQ8qEfGuMlFFbG0OYTfe5u8WYbAZAbbxsEQqFwuYst2tqairP+mg78ZwyiEH15WBIvgtgPk1gL4QLMDcQgFp42SLcvXu3WTMosa7a4Fo+38UGp0/1NsWLNhCQjfcx7xqAenjZIgwODt64caMSqwKEkMQCflwbBKF6owjZGML0iBUNcaFcjPHbBFB1L4OwXbt27dq1k+3ea2pq9PT0qMa7+jiOE4lEurq4Cfm/rudzG4Kx1Kvaa2dGfenHTLvIHhuI3yaAqpPXYJmioqKwsDA7OzsrK6s1a9ZILbNgwQJzc3MbG5vIyEjJ5ckVK1aEhoba2tpu27ZNUvLEiROW9cTHx8up2sr1tJLneeJkhDaEJpjrR+dVk11p6CAFUHUvv64mJSXt27cvLCysV69eS5cuLS8vl/qGZcuWNWe/CxYsMDExKSwszMjICAoK6t27d8eOHesXOHny5M6dO1NTUy0tLcPCwlatWrVw4UJCiLm5+axZs5YtW1ZTUyMpLBQK27Zte/z4cfFTY2Pjlv6cauFaHt/VFuN4NYSAJtt6MmHHRQMcaTssTAGgwl4GYWpq6qZNm2xtbXv16rVjx468vDypb2hOELIs+9tvv8XHxwsEAk9Pz+HDh//666+rV6+uX2bnzp3jxo1zcHAghHz55Zdz5swRB+HUqVMJIevXr29YUYHAwsKihT+dmrlewHfBnDIapKMlNbEtPQMLUwCotpftjw8//LCkpOTLL78khDx8+LCkEc3ZaV5eXllZWfv27cVPvb29Hz9+3KBMWlpa/QIZGRks+6Z5GpOTky0tLVu3bj1//vza2trGinEcV1VVVfw/jU0Up5oS8/lADBnVLIs6MymF/MFMdJACqC65XMkvKSmhKMrIyEj81MTEpKio6PUykh5OY2NjlmXLysoaa/P5+/tfuXLFw8PjwYMHY8eOpSiqsYbprVu3tm/fvnz5cvFTIyOjq1evqsVdHxxPEvN1vQ1rKypkeS92ZWXlGwYrwbtr8giv70KPvSDwMxLa6eMu+7eBv2G5Ylm2rq7uze0Qtaavry8QNJF0jb4sEomOHTuWnJycnZ1tb2/v6+s7ZMgQfX395nywjY0Nz/NlZWXm5uaEkOLiYjs7uwZlrK2tS0tLxY9LSkp0dXXFhaVydnZ2dnYmhHTu3Hnx4sWLFi1qLAj9/f1HjRoVHR3dnHqqlAclvI0B62JpJNvd8jyvqZdUVUSTR7ifMZlaxH5+nT4xSIDT+VvA37BciYPQwECrr2NLD8Jnz56Fh4ffvHmTYRgrK6vi4mKhUOjh4REbGyvpz3wDKysrOzu7Gzdu9OnThxCSnJz8+ru8vb2Tk5PHjx8vKdDML31CoZBhNPCKC26l12AL/JnQo6KfbnOz/TAYCkDlSP+3HDdu3LNnz2JiYmpqanJzc2tqak6cOMFxXEREBMc1fbWDoqhPP/30H//4x9OnT48dO3b8+PGJEycSQrKysvr161dZWUkImTJlyq+//nrx4sXHjx8vXbp0ypQp4vfeuXPn1KlTRUVFqampp06dEo/ZiYmJuXjx4rNnz+Li4hYsWDB69GiZHQCVcR1BqLkENNnTh1l5i00uRO8ogMqR0iIsKio6e/bsvn37RowYId5C0/TAgQN///33bt263b9/38fHp8n9LliwoLy8vHfv3hYWFr/99pu7u7t4O8//90TQtWvXn3/++YsvvqiqqoqKipo2bZp4e2xs7JkzZ6ysrO7fv3///v1FixbZ2toWFRWtXLmyoKDAzs4uOjpaHXs+m5RYwEe2RnNBYzkbUT++x3x4lk0aJjDETfYAqoSSJJNEYWGhtbX1/fv3vby86m8vLy83NTW9deuWn5+fAmvYMjNmzPD09FS7pBRyxOJXYe7HOkayPkWWl5ebmJjIeKdQT4uO8PjzrJGAbAzRwL59+cHfsFzhGiGR2jVqZWUVGBgYGxvbYPvhw4ednJyac40QWiqliG9jSsk8BUHVbAhm4rP5I1noIAVQIS9PveXl5ZJV6b/77ruJEydmZmZGRETY29sXFBQcO3Zs+/bta9eubXIcKrwFXCDUEsY6ZEcvZuQp0c0IHUw3A6AiXqbavn37Jk2aVP+1DRs2bNiwof6WqKioMWPGKKhq2uR6Pv+eLYJQK4TYURPa0lMvsAf7o4MUQCW8DMLevXvHxMQosSraLLGAn+6NkTLa4p8BTLfDon+ncpPb4ZcOoHwvg9DNzc3NzU15NdFeFUKSUc77WaJFqC10abI7lOl1RBRiR7U3x+8dQMnwhVT5Egv4jpaUDn4V2qS9ObWsCxN5mq0WKbsqAFqv0ZEvKSkp+/btS09Pb7AMhaauBahE1/L5rhgpo30+aUeffcHPucribgoA5ZIehPv37x8zZoy1tTXHcQYGBrq6uunp6cbGxl26dFFw/bTBtXx+pBuCUBttDmECD4r2POaiPNAhAKA00v/9vvnmmw8++CArK2vw4MEff/zxw4cP79y54+LiMnDgQAXXTxvg3gmtZaxDdvdmZl1h08txZyGA0kgJwqqqqrS0tHnz5uno6BBC6urqCCFeXl6//PLLokWLKioqFF1HjZZTTapEfGtTBKGWCrSmFnZiIuLZKlwsBFASKUFYWVnJ87ylpSUhxNLSsrCwULy9Q4cO1dXVDx8+VGgFNd3VPK6rDRZb02rTvelAG2rsORatQgClkBKE1tbWRkZGWVlZhBBPT88zZ85UV1cTQhISEggh4oAEWUG/KBBC1gczTyv51bexkD2AEkgJQoqievfuffDgQUJIVFRUSUlJx44dhw4dOnz48ODgYFdXV4VXUpNdy+e72mCghLbTZ8j+fszqW+yJZ2gWAiia9FPw5s2bv/jiC0KImZnZqVOnAgMD8/PzJ0yYcPDgwWYunwvNwROSWIAWIRBCiLMRtaePYOJ5UQYGzgAolvTbJxwdHR0dHcWPAwICfv/9dwVWSYuklvCWepSNvrLrAaoh1IH61p8ZFs9eGiLAUiQACvOmTrmampq7d+/GxcWlpKSIl5UH2cKt9NBAtA/d1YYah4EzAAokPQg5jlu8eLGtra2vr+/AgQP9/f2tra1nzZpVW1ur4PppNoyUgdetD2ayq/gVKRg4A6Ag0vtf5s6d+9NPP40aNWrkyJF2dnbi9Qg3btyYl5eHblIZul7Aj26NkTLwCj2G7O/LBB1m/Syp953xPQlA7qQEYXV19caNGxcsWPDdd99JNkZERHTp0uWzzz5buXKl5PIhvIs6jtwt5jtb40wHDTkaUQf7Mx+cFMUNEnTAsiQAcialOVJSUlJTUxMVFdVg+5gxY3iez8nJUUjFNN/NQt7TlDLEmAiQJtCaWhfMDI1nc6uVXRUATSclCG1tbW1tbdPT0xtsT09P19PT8/DwUEjFNN/VPFwghDcZ5U5P8KSHxouwVBOAXEkJQoZhVq9ePX369HPnzkk2JiUljR8/fsmSJebm5oqrnUa7lMcH2yEI4U0WdaZbm1CTEjCIFECOXgZhbGxs4P+sWbOmpKSkd+/e1tbWPj4+dnZ2gYGB6enpf/75pxLrqmEu5vAhCEJ4I4qQ7T2ZzHL+uxsYRAogLy+vUJmYmLRu3VrytP5jkLnMCr6O49tg0Qloij5DDg8QvHdY5GlGPsSyhQBy8DIIQ0NDQ0NDlVcT7XIxlw+xw0kNmsVanxwewPQ5KmptQr1niy9PADKGc7FyXMrFBUJoAW9zaltPwcjTbFYFLhcCyFijQXj37t3x48eL55Tx8/OLioq6evWqImum2S7l4gIhtMz7ztS8DnTYCbYQ8zsByJT0IExISOjSpcvevXtbtWo1ZMgQNze3Y8eOhYSE/PXXXwqun0YqE5JHZbiVHlpshg89zJUadEJUIVR2VQA0iPTbuaOjo9u3b3/06FF7e3vxlpKSktGjR3/xxRdDhw5lGEaBNdRAV/L4zlaULrqloeX+1YXJT2CHxYuODhTo4R8RQBaknIyLiopSUlJWrVolSUFCiLm5+bp163Jych48eKDA6mmmS7kcLhDC26EI2RTCGOtQnySwHC4XAsiClCCsq6sjhJiamjbYLt6CBSje3SUMGYV3IKDJnt7M0wp+1hVW2XUB0ATSp1izt7dft25dg+3r16/X19dv166dQiqmsVieXMWcMvBuDAQkdqDgUi6/+AayEOBdSblGSNP0okWLPv/883v37o0ePdrBwSE/Pz82Nvb06dMLFy40MjJSfC01SUoh72REWeopux6g5kx1yIkwQY8jIlMdbrYfOhgA3p70wTKfffaZjo7Od99999VXX4m32Nra/vDDD3PmzFFg3TQTphgFWbHWJyfDmB5HWBsDMrYNshDgLUkJQpFIdO7cubCwsMmTJ2dkZBQVFZmZmbVu3Zqm8Z8mA5dy+QGOCEKQDRdj6kQY0/eYyFyXhLvgPxTgbUj5z8nLy+vfv39GRgYhxN3dPSAgoE2bNkhBWbmIW+lBptqbU4f6CyYnsH/nYBQpwNuQEm9WVlYGBgaVlZWKr43Gy6zga1ne0wxBCLLUxYb6o49g1GlRUgGyEKDFpAShnp7elClT1qxZI76PAmQoIYfvbo+2NchebwdqcwgTHie6W4wsBGgZ6YNl3NzcYmJi2rVrFxYW5uTkJBC8LPb1118rqm4aKCGH72mP5iDIxXA3WsiRAcfZE2GMnyX+zACaS3oQLl++PDc3lxCyefPmBi8hCN9FQg4/rT1ahCAvka1pmiL9jouODhQEYjJbgOaRHoQ5OTnvvutbt26lpKS0a9eua9euUgvk5eWdPn3a1NS0X79+enovb6xLT09PT0/v1q1b/XsWS0pK4uLidHV1BwwYYGho+O7VU7yCGvK8iu+Ar+ogTyPdaUJI+EnR8TCBvxX+2ACaJqV1UlBQcOnSpfT0dI7j3nq/69atGzhwYEJCwujRoxcsWPB6gVu3bnl7ex87duyHH37o0aNHdXU1IaSsrMzS0rJLly4DBgx48uSJpPCTJ0/at2+/b9++zZs3BwQElJSUvHXFlOhCLtfNlmJwagI5G+lObwxhBp0Q3cDYGYBmeCUIRSLR+PHjbW1tQ0JCPDw8vL29Hz58+BY7raqqWrRo0cGDB7du3XrmzJkff/zx9Sbm999/P2XKlF27dp05c4bn+T///JMQYmhomJiYWFhY2OBujZUrV4aHh8fExBw/ftzV1fWXX355i1opXUIO3wMjZUAhhrvR23oKBp8UXctHFgI04ZXz8saNG3/99dfevXsvXbp02rRpGRkZEydOfIudXrhwwdTUNCgoiBDi7u7eoUOHkydP1i/A8/yRI0dGjRpFCGEYZvjw4bGxsYQQgUDQunXr13cYGxs7cuRIQghFUSNGjBAXVjsJOXwPjJQBRRnsTP2np2BIHLIQoAmvXCM8fvx4z549T58+LX7q5+f3xRdflJaWmpmZtWin2dnZjo6OkqeOjo7Pnj2rX6CoqKi6ulpSxsnJ6fDhw43tjeO4Fy9e1C+cnZ3dWOHCwsLc3FwdHR3xU4FA8PHHH0ueKlGFkNwr5jtbcqxiJ0lmWZZV8EdqGVU+wgNaka0h9JA40YG+dFcbZdfmbanyEdYA7P8ouyLyQtM0RTXRAnklCJ88efLxxx9Lng4cOJAQkpGR4e/v36IPZlm2ft+mQCAQiUQNChBCJAv8MgzToEB9HMdxHNfMwhUVFXl5eYmJieKnBgYGw4cPV4WJwi++oDpaUgwnFL79hde3IRQKhUIsZy5HKn6EB9iTrd2o8Hh+Vw++t71aNg1V/AirO5ZlhUJh/XvkNIyOjk6Ti8m/8sNXV1fXH5Apzo+qqqqWfrC9vX1+fr7kaW5ubr9+/eoXsLa21tHRycvLs7KyIoTk5eU5ODg0WkWBwMbGJi8vz8vLq8nCrq6u/fr1i46Obmmd5e1qEdvLgejr6yr4c4VCob6+voI/VKuo/hEe0pr8ZchHnhb9GipQx3luVf8IqzVxu0XLj3DDbwHPnj1LSkoSPy4sLCSEpKam1r+3ISAgoMmdvvfee1lZWRkZGe7u7qWlpdevX9+6dSshpK6urra21sTEhKbpHj16xMXFtW/fnhASFxfXICkbCA0NjY+P79mzp7hwaGhoi35IVZCQw3/VoYlvJQBy0tOe+qufYPgp0dbuzFBXjNgCeEXDIPzxxx9//PHH+lsmTZpU/ynPN927Ym1tPXXq1IiIiEmTJu3duzc8PLxt27aEkG3btm3ZsuXmzZuEkPnz50dGRtbU1GRmZt66deu3334Tv/frr78uKSnhOO6f//ynhYXF6tWrjY2Nv/rqq759++rq6paVlR07duzGjRvv8jMrnpAjiQWYaxuUKdiOOjZQMDSefVxGsH4hQH2vBOHSpUsrKipkst+ffvrpjz/+uHnz5vjx48eNGyfe2KtXL1tbW/Hj/v37x8XFHTx40MXFJTEx0dLSUrzdz8+vqqpK0u4U91wHBgZevHgxJibGzMwsMTHRxcVFJpVUmKQC3sOUMlN0tyjAKwKsqatDmKHx7L0SflMIo4M0BCCEEEI1p4WnRmbMmOHp6alq1whX3uKeVvJruymha7S8vNzExETxn6s91O4IVwjJx+fYciG/r6/AQq/p8kqndkdYvbAsW1dXZ2BgoOyKKBO+EyrC+RdcqAP6RUElGOuQv/oxna2okFjRk3KN+h4M8HYQhHIn4siFXL4n5pQBlUFTZGUQM9uP7nGETcQ0bKD1cHaWuxuFvKsxZa3Vg5NBFU1uR28Mod8/KTqShSwErYYglLtzL3j0i4JqCnehjwwQTLvIrrur2IkeAFQJglDuzr3gerdCEIKK6mJDXR3CbH/ITb3AipCGoJUQhPIl4sglXCAE1eZoRJ17X/CknB9+SlSOucxA++AELV9JBby7CWWpDoPUQZuZ6ZKjAwVORlS3w6K0MlwyBO2CIJSvs7hACGpCQJNNIcwMH7rbYdHhTHSSghZBEMoX7iAE9TLFiz7YX/D5JW5FCoeGIWgJBKEcCXGBENRQiB11ZQjz1xNuzBm2ApcMQQvgHC1HiQV8GzNKLWaxAqjPyYj6+wOBmS4JPCi6X4KWIWg4BKEcnXvB97JHvyioJT2GbO3OzPKlex8VHXuKLARNhiCUo3PPcYEQ1Nu09vSB/oJpF9glybhkCBoLQSgvdRy5ksf3dMARBvXWzZa6NkwQn80NjWNL65RdGwA5wGlaXi7l8u0tKHOsQQjqz96AnB4saG1KuhwSXc9HyxA0DYJQXuKecQMd0S8KGkKHJmveY5YF0uFxoqU3ORZpCBoEQSgvJ7P5AU44vKBRRrjTScMEZ55zfY6KMisQhqAhcKaWi7xqklHOB9mgRQiaxtGIihskCHelux4S/fEYE9CAJkAQykV8NtfbgRbg6IImoinylR99erBgWQoXeZotrlV2hQDeDU7VchGXzQ9wQnMQNJmvBXV5iMBanwQcFP2dg25SUGMIQtnjMVIGtIOhgGwMYdYFMx+eZWdfYatFyq4QwFtBEMrerSLeVJdyM0EQglZ435lKiRC8qCYBB0XXcHMFqCEEoeydfMYPRL8oaBMrPbKnN7O4Mz00TjT/OlvDKrtCAC2BIJS9k8+4gbhxArRPZGs6JULneSXx3Y+rhqBOcL6WsUoRSczHXNugpWwNyK+hzLIu9Jgzoi+vsJW4agjqAEEoY+df8IE2lLGOsusBoDyj3Ol7I3WqRKTdXtFBLHYPKg9BKGMnnnEDHHFUQduZ65It3Zn/9GTmXOE+OsvmViu7QgCNwylbxg5n8uGu6BcFIISQ/o7U7RECF2PS4S/h5vscVnIC1YQglKWbhbwuQ7zNEYQA/2UoIMu6MGcGC357zIXEim4UIAxB5SAIZelQJj/EBSkI0JCPBfX3B4JPvej3T4o+v8gWYVY2UCUIQlk6lMkNdcUhBZCCImRSW/reSB2aIt77hL88QE8pqAqctWUmq4LPruKD7dAiBGiUhR5ZH8ycCBPsSuO6HhJdzkMYgvIhCGXmcBb/vjPNIAcBmuJvRZ3/QDDHjx59mh13js2uRByCMiEIZeZQJjcU40UBmociJMqDvj9K4GJM/A+Ivr3OltYpu06grRCEslFSR67n8/1xByFASxgJyPeBTEqEIL+GtNsrXHOHq8P996BwOHHLxrGnXC8H2lCg7HoAqKFWhtQvPZgz7wvOvuC99or+eIxhNKBQCELZOJTJo18U4F14m1OH+jM7ejFr7nJdD4pOPkMagoIgCGWgliXxkK+lQAAAH2hJREFU2dwHzjiYAO+qpz11eYhgXkd69hW222HRCcQhyB/68mTg1HPe14KyNVB2PQA0AkXIKHd6hBu9N4P76iprQOt+04kb5krT6HMB+UAjRgb2Z3Aj3HAkAWSJpsjo1vStCMEsL9GqW5zXPtGWB1w11nUCOZBXi5Dn+e3btx8+fNja2vrLL7/09fV9vcyVK1c2bNhQWVkZFRU1atQo8ca6urqffvrpwoULbm5u8+fPd3R0JITcunVrw4YNkjdGR0dL3aFSiDgSm8X9MwBtawDZoykyxIn7qL3gYi6/8hb33Q3R3A70VC/aAP9wIDvyasds2rRp2bJlU6ZM8fT0DA0NLSgoaFAgLS1t4MCB3bp1GzduXHR09IEDB8Tb582bd/jw4ejoaEJIv379WJYlhGRmZp4+fTrgf8zMzORU7bdw9gXvYUo5G6HXBkCOQuyog/2ZY2HM3zm8R4xwzR2uCq1DkBVePjw9Pf/66y/x47CwsFWrVjUoMGfOnAkTJogfb9q0qWfPnjzPl5aWGhsb3717l+d5juNat24dGxvL8/zhw4dDQkKa87nR0dFr166V1U/RHFMTRD+ksIr8xBYpKytTdhU0HI6wvL1+hG8WciNOiex21y29yZbWKaVSmkMkElVVVSm7FkomlxZhSUnJo0ePunfvLn7avXv369evNyhz/fr11ws8ePBAV1fX29ubEEJRVEhIiOSNT58+nTZt2vz58xMTE+VR57fD8uRgJhfhhuYggOJ0tKT29WXOvi+4X8x7/Cn89jomaYN3IpeO9pycHIqiLC0txU+trKxycnIalMnNza1foLq6uqSkJCcnR7Kx/httbW3HjRvn4eHx4MGD0NDQ3bt3Dxs2TOpHP3jw4MSJE5KOVoFAsGPHDlNTU9n+gBIX8mh7fYEdXVlRIadPeFeVlZUUhZyWIxxheWvsCDsLyMZAktme2pDKdNhP93PgP28rCrBCIrYMy7J1dXXii1AaSV9fXyBoIunkEoRGRkY8z9fW1hoaGhJCqqurjY2NG5QxNDSsqakRP66urqYoytDQ0MjIqLb25UplkjcGBQUFBQWJN9rY2CxfvryxIHRxcWnbtq3kVUNDQwcHB/mdp47dZkd5UMbGunLa/7vjef71gw8yhCMsb28+wj7GZKMdWd6NbEvlJl5hnIzIHD863AX3WjSXOAgNDLT69i+5BKG9vb2urm5GRoaPjw8h5MmTJ87Ozg3KODs7P3nyRPz4yZMn4rc4Ozvn5uZWVVWJE/TJkydhYWEN3ujl5ZWbm9vYRxsaGnp6evbr10+GP05jeEIOPOHjB+HGCQAlM9UhX/rSM3zo/Rnc98nc19e42X702DYYXArNIpeTuI6OzrBhw7Zv304IKSkpOXDggPjuiJKSkq1bt9bV1RFCIiMj9+zZI24Ubt++PTIykhDStm3b9u3b7969mxCSnp5+4cKFiIgIQsjTp095nieE1NbWbtu2TdI6VK4rebyZLvEyxzdPAJXAUCSyNX19mGBLd+boU97tT+GiJDanWtnVApUnr+9L33///YABAy5evPjs2bOwsLDevXsTQl68eDF16tTIyEhdXd3IyMiYmBgfHx9zc/Pa2tpTp06J37h27dpRo0bt2rXrwYMHCxYscHFxIYQsXLjw1KlTzs7O6enpbdq02bt3r5yq3SL7M7gRGCYDoHp6OVC9HJhHpfSau5z3PmGYE/2pFx3qgGu5IB0lbmnJg1AovH37toWFhbu7u3gLx3ElJSUWFhaSi3YPHz6srq729fVlGEbyxqqqqnv37jk6Ojo4OEg2ZmRkFBQU2Nvbv97LWt+MGTM8PT3FtyHKFccT1z9EcYOY9qrdIiwvLzcxMVF2LTQZjrC8veMRLqkju9O4rQ+4WpZMbkdPaEvb6MuwdmoP1wiJXOca1dHR6dy5c/0tNE3XHxRKCGnbtu3rbzQ0NAwMDGyw0d3dXRKoquB8Dm9nQFQ8BQHAXJdM96ane9NX8vgtD7i2McLBzvTU9nRPe/zzwn9hoMdb2v2I+6gNjh6A2njPlvpPTyZ9tE5XG+qzC2z7faLVt7k8XEEEBOHbqRaRg5lclAeOHoCasdAjM33puyMF/+7B3C3mvfYJR51m47N5rAWszTC4+G3EZnFdbCh7re5UB1BvIXZUiB2zRsj8nsbNu8aWC8mn7ejxbWn8X2shtGnexm+PefSLAmgAUx0yrT2dPFzwe2/mURnvvU/4/knRH4+x3pN2wdm8xQpqSEION9wVhw5Ac3S1of7dg8mO0vmoDb3jEee0Rzg5gT3/Qm6j6kGVoGu0xWLSucHOtLGOsusBALJmICAfetAfetDPq/g9j/mZl9mSOvJxG+rjNjSmztBgaNa02G+PuY8wTAZAo7UypOb40TcjBLEDmDqO9DvOdjkoWnOHe4ZlLjQRTugtk1rKp5fx/R3x3RBAK/hZUj90ZbLGCJZ3ZW4V8f5/iYIPi368zWVWIBE1B7pGW2Z5CveZNyPA9wcAbUJTpG8rqm8rRsgxZ57z+zK4wIOsuwk1wo0e6U55mOKbsXpDELZAejl/JIt7FInLgwBaSocmA52ogU7M5u7MuRf8/gyueyxrb0iNcKNHuFOYakpNIQhbYNlN7rP2tLnqLj4IAArC/K+NuD6YuZDL78/gBh7n9AWkvyPV35Hq7UCb4UShPhCEzZVZwR94wj1EcxAA6qEp0tOe6mnPrOlGbhXx8dn8pnvcuHOsvxU10Ike6ER1tqKwSrCKQxA21/IUbooXbamn7HoAgEqiCOloSXW0pL7yo2tYkpDDn3zGTTzP5Vbzg53p912ogU60Kb5IqyQEYbM8q+T3pnMPRuGvGACaps+I+0iZVUEks4I/msVvT+U++ZsNtKYGONEDHCl/NBNVCYKwWZancJPa0dZYxgwAWsjVmPrcm/rcm64SkfMv+Lhs7uNzXEENH+pA92lF9W1FeZohEpUMQdi0h6V8TDp3dySagwDw9gwFZJAzNciZIYQ8q+TPPOfPPOeX3uQYmgxwpAY4Un1a0Ra4+KIMCMKmzb3GzevIYFVrAJAVJyNqnCc1zpMQQh6U8HHZ/H8ecpP+ZtuaUX1aUX1a0SF2FOZxVBgEYRPOPOfvFvN7+zLKrggAaCYvc8rLnJrhQws5cjWPP/2cX3qTvVHIe5pSIXZUiB3V04FqZYjuUzlCEL4Jy5PZV9gfutK6mEoGAORMhybd7anu9tQ/OtN1HEkq4C/n8jEZ/IzLrKku1cOe6mFPdbOlvMwppKJsIQjf5D8POXM9EuGGGAQAhdKlSTdbqpstNduP8IS5X8In5PBnn/P/SuZK6vhgO6qHPd3Dngq0pnRwfnpnCMJGlQvJoiT2yAAcIgBQJooQb3PK25ya6kUIITnV5FIu93cOP/0S96iUD7CmQuyo92zpbnaUFcbavBWc5Rs17xr7gQvd2RqdEACgQuwNSIQbHeFGCCFlQnIpl7+Sx627y358jrc1oLraUEE2VFcbqjMai82GIJQuLps//pS/NQLHBwBUl6kOCXOiwpwYQgjHkwel/LU8/mo+vy2Ve1zO+1tR3Wyp92ypQGvKxRjf6RuFE70UJXXk0wR2e08G8yEBgLqgqf/2oE5oSwgh5UJyLZ+/nMvvfMRP///27jyqqSt/APj3JcgiJCEJCpFdEBQXFK1LqSBU60hdWZWZ0jp1nI72eOzxtODSqeOp2nZacWxrpeMydQqc6a89MIqCAWVTVkdRkGohshmWRAhLEkhI3v398TpvGFtpx4oxyffz13v33bx8c+Hk+27eve+WGWkCz0yg5k/gLJhIzZ9A4eIBI2Ei/BHby42rvKjnJ+EFFELIXPHGMetjMN9j3HsaUq0klUpysIa+dp9IxlPPTKCemUCFiGCKIzg4mDha08JE+KAzLfSVLlITjS2DELIcHo6UhyO1zgcAwEjgdi+pUpJqJfmyga7vtZngYJjuDLNE1BwXao6Y8uNb1wwN/Lr/L00D5LUrxv973sYRGwYhZKG4FEwXUtOF1MYAMBphSKfvMtrX9pAbPZDRSN6qonuGSLCYmiOm5oipYDEV5EzZWfQzRfD7/j8GhmG11LhrNjfU1aouhhBCVo1DwWQeNZlHrfH+vkSlg2vd5Ho3KWgnh+roxn7i7UTNFFLThdR0IcwUUX48ysaChqRiIvyekUBioWGxG/V6kAX9eRFC6H8ntBt5fxGGabjTR26pSG0PSW+EOhUt15KpAmqmiJohpKYIIEBA+fHMuNeIifB7KVXGQQMcWWS2f0mEEBob4zgwQ0jNEFIJk78v0Rqgvpfc7CH1KlLSSRr6oEVN3MZTAXwIEFCBzpQ/n/Lng7eTecxlxEQIAJBaR59pJRWrbSyps48QQmNkvA3Mc6HmjXjeiIGGFjX5rg/u9JF6FTnTQjf2Q7uWeDhSfjzw41N+fGoyDybzKV8e9bTNTLP2REgAdlcbs1tI/gourgSGEEKPxobDZDtY4fmf7KinoWWANPaDbIDI+klJJzQN0Hf7iYMN+DhRvjzKhwc+TpSXE+XLAy8nimeiBGnVidBAw2tXjHUqUrrKBp/RhxBCj5ctB6YIqCkCAPivEYiKQWhWk6YB0jQAtSpyro1uHoAWNbHngi+P8uFRXk7g5Uh5OoGHI+XpSEnGj22cVp0Ik4qNKh25GIWTJRBC6MmZ6ADMY1EfKGcSZPMAadPA3QFS3AltarpNQ3r14OlIuTuChyPl5gABAmrz1Md5H8uqM8Br0ziLJprHvVyEELJ4D0uQQ0ZoVZMOLbRpSOcgdA0+5ve16kQY5obzBRFC6Glnz4UAARXwg59YHxfsDSGEELJqmAgRQghZNUyEFi4vL8/UIVgyo9EolUpNHYUlGxgYKC0tNXUUlqy9vb2mpsbUUZgYJkILl5SUNDQ0ZOooLFZnZ+e2bdtMHYUlu379+oEDB0wdhSWTSqVpaWmmjsLEMBEihJD1IoSYOgTTw0SIEELIqmEiRAghZNUoC+sXR0VF3b1719PT09SBPC0KCwvDw8M5HLziGRM6na6qqmrx4sWmDsRi9fb2NjY2zps3z9SBWKz29vbe3t6goCBTBzJW1q1bt2XLltHrWFoirK6uvnfvHo/HM3UgT4umpiZfX19TR2GxCCHNzc3YwmPHYDB0dHTgpe3Y0Wg0arXa1dXV1IGMFV9fXz8/v9HrWFoiRAghhP4n+IsZQgghq4aJECGEkFXDRIgQQsiqYSJECCFk1bh79+41dQzol5LL5bm5ufX19SKRiB0xSwgpLS0tLi7u6enx9vamqO+XL9FoNGfPnr1x44a7u7uDg4PpojYnDQ0Nubm5Mpls0qRJ9vb2bLlCoWAa08HBQSwWM4Xd3d3Z2dkNDQ3e3t7jxo0zUchmpqamRiqV3rt3z8vLy8bmv5aHa2trq66unjRpElve1tb2z3/+Uy6X+/j4cLlcU8Rrfq5cuVJUVKRUKr29vZn5VDRNV1dXX7x4US6Xe3p6jmz227dvnz17tq+vb+RXhwXDUaNmLzc3NzExMTY21mAwZGdnZ2VlLVmyhKbp1atXy+XyiIiIkpISZ2fnvLw8GxsblUq1cOFCPz8/Pp9fXFxcVlaGQ/9/Ulpa2p49exISEhQKRVFRUWlpaWBgIACcP3/+pZdeCgsL4/F4nZ2dzNO3GxsbQ0NDn3/+eZVK1draWlZWJhAITP0JnnY7d+5MT0+PiYmpr6+XyWSVlZXsVYVer1+0aNG1a9eam5u9vb0B4PLly2vWrFm7dm19fb2Dg0N+fj7mwtERQhISEurr61944YXy8nIOh3Pp0iU7O7u1a9cyczQbGxvlcvnly5fd3d0BIDMzc9u2bdHR0ZcvX54/f/6pU6dM/QnGHkFmLioqat++fcx2cnJybGwsIaSurs7W1ra/v58QotFonJycqqqqCCEHDx5cvnw5U3nz5s1btmwxUdTmJDAwMD09ndmOjY1NSUkhhKhUKqFQeP78+Qcqb9q0iWlVmqaXLl360UcfPeFozc7w8LC9vX1FRQWzO2fOnL/+9a/s0XfeeSclJQUAmpubmZKIiIhDhw4RQnQ6XUBAQHZ29pOP2bw0NTVRFKVQKAghOp3OxcXl0qVLhJDGxka2zvLly5OTkwkhRqPR19eXadXu7m6BQFBXV2eiwJ8cvEdo9sRisUajYba1Wq2LiwsACAQCiqIGBwcBQKfTEUKEQiEA5OTkxMbGMpVjY2NzcnJMFLU5EYvFWq2W2dZoNExnJS8vz8vL69lnn7106dJ3333HVmZbmKKomJgYbOGfxOVyBQIB08JGo3FoaIjtDtbW1mZlZb355ptsZY1GU1RUFBMTAwC2trarVq3CFv5JTk5Otra2zLeBXq83GAxMC4+cZu7m5qbT6QCgrq6uq6srKioKAEQiUUREhDW0sM1PV0FPtw8++ODXv/71ihUrjEYjh8NJT08HAA8Pj5MnT4aHh8+aNevmzZsff/yxv78/AMjlcubXDwBwd3fv6OigaRofwDa6EydOvPTSSzk5OUql0t/f//XXXwcAmUxG03RoaOisWbNKSkpiY2MPHz5sMBgUCsXIFpbL5SaN3QxQFPXNN9+89tprgYGBMpls7dq1a9asAQCDwbBp06ajR4/a2dmxldvb2wFg0qRJzK67u3t9fb1JwjYjLi4umZmZy5cvnzlzZl1d3f79+2fNmjWywrfffpudnV1UVAQAcrnc1dWVvbdtJf/D+A1o9s6fP9/e3p6QkJCQkHD37t0LFy4AgFarPXz4cHh4eFxc3AsvvPCXv/ylv78fAJhkybyQy+XSNE3TtCmjNwcZGRlcLjc+Pn79+vXFxcVVVVUAMDQ0dOfOnZycnIyMjKtXrx4/fvzatWs0TRNC2MEFXC7XYDCYNHbzcPToUS8vr7i4uPj4+MzMzIaGBgD48MMPFy5cGBoaOrKm0WikKApb+H+i1+tTU1MXLFgQHx//4osvfvrpp/fv32ePdnR0rFmzZt++fbNnz4Z/tzB71Fpa2NS/zaJfys3Njb1T9fXXX/v7+xNC/v73v0+fPp2ts2DBgrS0NELI/PnzT58+zRQWFxdLJJInHq+Z6e/v53K5t2/fZnbff/995ibrp59+OmXKFLZaSEjIF198QQgRi8VXrlxhCk+ePPncc8898ZDNzM2bN21tbbVaLbPL3roWiUQbNmzYvHnzb3/7WwBITEysrKxUqVQA0NnZyVTetWtXUlKSyUI3E1lZWb6+vsxVGiEkMjKSvXWtUCiCgoL279/PVr569Sqfz2crJyQksEMQLBj2CM2ejY2NXq9ntnU6HTOCjsvlDg8Pk38PCdbr9UxHMCIigukyAoBUKl2yZIkJIjYrXC6XEDI8PMzssi0cGRnZ1dWlVquZwra2NubB0EuWLGGGjwK28M/D4XBommZbWK/XMy184sSJdevWLV26NCIiAgBCQ0Pd3NycnZ1nz57NtDAhJD8/nzmKRsH06tjffthvA5VK9atf/SouLm7Xrl1s5RkzZtjZ2ZWXlzM1i4qKrKGFcfqE2Xv33XfT0tK2b99uNBpTU1OTk5O3b9/e29sbHBy8cOHCyMjI0tLSS5cu3bhxY8KECa2trSEhIS+//DKfzz906FBhYWFISIipP8HTbuPGjdXV1b///e+7u7tTU1NPnz7N3MRKTExUKBTR0dE5OTmDg4MXL17kcDhVVVXLli3bsWOHSqVKT0+/fv06e8sQ/SiapiMjI2ma3rBhg0wmS0tLKykpmTNnDluBGfbMTp/4+uuvt2zZ8uabb964caO6urqmpganw45OrVaHhIQEBQWtWLGisrLy7NmzNTU17u7uL774YkVFBTt6Ljg4mFmu6IMPPvj888+3bt0qlUqHhoYKCwtNGv6TgInQEhQWFpaUlHA4nCVLlrBr4/X19WVmZra0tHh6eq5fv14kEjHlLS0tX375pcFgiIuLs+BFyB4jmqbPnDlz/fp1e3v7qKio4OBgptxoNGZmZtbX10+dOnXDhg3s+ILa2tpvvvnGzs7uN7/5Da4f9HPo9fqvvvrq9u3bzs7OMTExD8xtNRgMJ0+eTExMdHJyYkouX76cl5cnEoleeeUV9h8bjUKtVmdkZDQ3N0skkvXr10+YMAEAsrKylEolW8fb23v58uXMdk5OTllZmaen5yuvvGIN1xmYCBFCCFk1vEeIEELIqmEiRAghZNUwESKEELJqmAgRQghZNUyECCGErBomQoQQQlYNEyFCCCGrhqtPIPTompubIyMjH3Z04sSJFRUVTzKex6W2tra8vPzVV1/FNW+RNcBEiNCj4/P5zCOhGQcOHPDx8UlMTGR2eTyeieL6pfLz83fs2PHyyy9jIkTWAJ8sg9BjIxKJwsLCsrOzRxZqNJqBgQEXFxcbmx+/7uzp6bG3tx8/fjxbolQqxWLxw9aJ7Orq4vF4I+uzdDqdSqVydna2t7f/0dcODg729fW5uroyS+0QQpRKpVAoZJ8Pxzh06NCOHTuGhoZGrgX4M9E0rVAoRsag1Wr7+/tHaQGETAvvESI0VmprayMjI/l8vkQiEYvFf/zjH9kVAE6dOiUSicrLy+fOnSsWi3k83qZNm4aHh7Ozs728vCZOnOjs7Hzs2DH2VPv27fP19S0oKJg8eTKzCMPGjRuZNccZPT09SUlJQqFQIpHw+fz169czKxYBwJ07d0QiUUZGRnx8PI/Hk0gkAwMDdXV1S5cudXBwcHV1dXR0fOaZZ5gFBwDg3Xff3bNnDwBIJBKRSCQSibRabV1dnUgkunjxIvuOD5SsXLkyPj7+888/d3Nzk0gkn3zyCQDcunVr6dKlTAuIRKLdu3cbjcYxbHGEHgleoCE0JmQyWVhYWFBQUG5urru7u1Qq3blzJwDs27cP/t11S0pKeuONN9LS0vLy8t5++22apsvLyw8fPuzh4ZGamrp169bw8PBp06YBwODgYEdHx8aNGw8ePBgSEnLhwoXk5GQAOHXqFHO2ZcuWKZXK48ePz549u76+fvv27evWrWPWHDcajSqV6o033li5ciWznoCdnd39+/fnzp27e/duiUTS2tq6d+/eqKiohoYGFxeXDRs2yOXyY8eOpaenMz1FOzs75iTsYknsadmSgYGBqqqqGzduHDlyxNvbm1kvIiwsLCAg4Ny5cx4eHgUFBSkpKYSQAwcOPOG/BUI/wWQrISJkcYRC4Zo1a5jtpKQkiUTS19fHHt25c6eTkxOzTuRnn30GAJ999hl7dPbs2RRF3bx5k9nt7e21sbF57733mN2UlBQAOHHiBFs/OTmZw+G0tbURQo4fPw4A1dXV7NG8vDwAqKysJITcunULAFasWDFK5AqFgsPhsIs2f/TRRwAwNDTEVqipqQGA3Nzch5WEhYVxudw7d+6wFV599dWJEyeqVCq25O233x4/frxOpxslEoSePOwRIjQmpFLp1KlTq6qq2BIej6dWq5uamqZMmcKUsKveAEBAQIBCoZg5cyazKxAIXF1d7927N/KcMTExI7fff//9mpoaDw8PqVQqFot7e3sLCgqYowaDAQDq6urmz5/PlKxevfqBCDs7O7/66qumpiatVgsAdnZ2jY2Nv+Qjz5w5MyAgYGQLBAYGXr16lS1xcnLSarUymYzp5iL0lMBEiNDjRwhRKBS9vb3x8fEjy4VCoVKpZBOhUChkD9na2o7cZUr0ej27a29vLxAI2N1JkyYBAJMpu7q6+vr6fvhePT097K6bm9vIo2fOnElISPD19Q0NDRUKhRwOh8Ph9Pf3P+onBgBwdXUdudvV1aVUKn+0BTARoqcKJkKEHj+Kong83urVq0+fPv24zjk0NKRWq9nFaRUKBQBIJBIAEAgE7u7uzc3No4c0cnf//v2LFi0qKChgxqYaDIZDhw6N8nJmwCfT0WT8MGs+MMyVz+cvW7YsIyPjJz4YQqaGo0YRGhPh4eFSqbS3t/cxnvPChQvsNnMXcPr06cx7tbS0VFZW/vxTNTU1hYSEsKkrPz9/5EAYZgbkyFGpHh4eACCTydiS4uLi0d8iLCysoKBgZK8UoacTJkKExsTevXv7+vpWrVpVUVGh1Wrb29vPnTu3devWRz7huHHj3nrrrdLSUq1Wm5WV9d57761cudLf3x8Afve73/n5+cXHx2dnZ6tUqp6envLy8m3btrW2tj7sbMHBwf/4xz+uX7+u0+ny8/O3bt06csogk18//PDD8vLyf/3rXzRNCwSCBQsWpKamVlZWdnV1nTp16ujRo6MH/M4776jV6lWrVpWVlTEtkJub+4c//OGRWwChMYKJEKExMWfOHKlU2tfXt2jRIkdHR3d39/j4eLVa/cgndHR0/NOf/hQVFeXo6BgdHT137twvvviCOcTj8YqKimbMmBEdHS0SicRi8eLFi2tqah42rR4Ajhw54ujoGBISYm9vHxcXt3//fpFIxB599tln9+zZ87e//e25556bN28e0zVMS0ujKGrhwoVubm5//vOfDx8+PHrAs2bNKigo0Gg0oaGhTAvExMT8wtuQCI0FfLIMQo8NM1/+gVtlzc3NCoVCIBD4+Pg8woNaGDt37jx27JhKpdJoNN9++62zszPTF3xAd3e3TCYbP368l5cXn88f/ZwGg6GxsXFwcHDatGmjpMwHXtLQ0EDT9LRp0x724Jsfamlp6erq4vP5vr6+j9wCCI0dTIQImQE2EZo6EIQsEP40ihBCyKrh9AmEzMC6desCAwNNHQVClgl/GkUIIWTV8KdRhBBCVg0TIUIIIauGiRAhhJBV+38sdiKPlB9YUQAAAABJRU5ErkJggg==", "image/svg+xml": [ "\n", "\n" ], "text/html": [ "\n", "\n" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Distribution, AKA MZrn from MELTS zircon saturation (meltstzirc) calculations\n", "dist = [0.00102226,0.00105785,0.00109886,0.00113827,0.00118283,0.00122527,0.00127,0.00131584,0.00136402,0.00141291,0.00146517,0.00151593,0.0015771,0.00163667,0.00169765,0.00176478,0.00183148,0.00189186,0.00196526,0.00203994,0.00211244,0.00219109,0.00226588,0.00234415,0.00242908,0.00250853,0.00259262,0.00267285,0.00276366,0.00285942,0.00295,0.00304668,0.0031507,0.00324811,0.00334256,0.00345641,0.00357056,0.00368579,0.00381205,0.00393925,0.0040682,0.00420839,0.00434517,0.00448705,0.00462617,0.00478174,0.00494046,0.00509347,0.005262,0.00543599,0.00559582,0.00576977,0.00596128,0.00615107,0.00635583,0.00656449,0.00678063,0.00700385,0.00723641,0.00747424,0.00772279,0.00797097,0.00824662,0.00852537,0.00880685,0.00910631,0.00941175,0.00970515,0.0100419,0.0103883,0.0107337,0.0111116,0.0114879,0.0118724,0.0122876,0.0127015,0.0131353,0.0135732,0.0140491,0.0145477,0.0150343,0.0155732,0.0161175,0.0166661,0.0172398,0.0178859,0.018519,0.0191786,0.0198861,0.0206008,0.0213348,0.0221012,0.0228635,0.0236548,0.0244204,0.025197,0.0260081,0.0268075,0.0276271,0.0284577,0.0292575,0.0299837,0.0306664,0.0313003,0.0318339,0.0323055,0.03268,0.0329388,0.0330462,0.0330452,0.0329532,0.0326285,0.032175,0.0315716,0.0308023,0.0298608,0.0288332,0.0276613,0.0261425,0.0245284,0.0228927,0.0210703,0.0192917,0.0175039,0.015764,0.01409,0.0124369,0.0109181,0.00945112,0.00806163,0.00679755,0.00566667,0.0045914,0.0036195,0.00290947,0.00229281,0.00177751,0.0013579,0.0010393,]\n", "# Corresponding MELTS temperatures. \n", "T = [809.92,810.521,811.122,811.723,812.325,812.926,813.527,814.128,814.729,815.331,815.932,816.533,817.134,817.735,818.337,818.938,819.539,820.14,820.741,821.343,821.944,822.545,823.146,823.747,824.349,824.95,825.551,826.152,826.754,827.355,827.956,828.557,829.158,829.76,830.361,830.962,831.563,832.164,832.766,833.367,833.968,834.569,835.17,835.772,836.373,836.974,837.575,838.176,838.778,839.379,839.98,840.581,841.182,841.784,842.385,842.986,843.587,844.188,844.79,845.391,845.992,846.593,847.194,847.796,848.397,848.998,849.599,850.2,850.802,851.403,852.004,852.605,853.206,853.808,854.409,855.01,855.611,856.212,856.814,857.415,858.016,858.617,859.218,859.82,860.421,861.022,861.623,862.224,862.826,863.427,864.028,864.629,865.23,865.832,866.433,867.034,867.635,868.236,868.838,869.439,870.04,870.641,871.242,871.844,872.445,873.046,873.647,874.248,874.85,875.451,876.052,876.653,877.255,877.856,878.457,879.058,879.659,880.261,880.862,881.463,882.064,882.665,883.267,883.868,884.469,885.07,885.671,886.273,886.874,887.475,888.076,888.677,889.279,889.88,890.481,891.082,891.683,892.285,892.886,]\n", "# Both of the above must be in order of increasing temperature (=increasing age)\n", "# We are implicitly assuming a constant cooling rate by treating a \n", "# scaled temperature distribution as equivalent to a scaled time distribution\n", "\n", "# Plot distribution\n", "plot(T, dist, label=\"MELTS f_xtal\", ylabel=\"Probability Density\", xlabel=\"Temperature\", xflip=true, fg_color_legend=:white)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Run MCMC to estimate eruption age" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Configure model\n", "nsteps = 400000; # Length of Markov chain\n", "burnin = 150000; # Number of steps to discard at beginning of Markov chain\n", "\n", "# Run MCMC\n", "(tminDist, tmaxDist, llDist, acceptanceDist) = metropolis_minmax(nsteps,dist,ages,uncert; burnin=burnin);" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "Estimated lockup age:\n", " 102.16592593727708 +0.18860998126889683/-0.23459261427365163 Ma (95% CI)\n", "\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3dd1xT9/oH8G8WJOypBFAQHCjiwIkLLIpY3APrxtFaR911dnmts2ql+mu1anFUrVpvVUCtA0FUREUFHAgOQGXvkZB5fn+c3jQcEBchyPm8X/d1mzw5+eY5McmHszkURREAAAC24uq7AQAAAH1CEAIAAKshCAEAgNUQhAAAwGoIQgAAYDUEIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAq/H13YA+HTp0qGXLlp07d9Z3I5XQJ73jcDj6bqTeUalUPB5P313UOw3mbUl89Dj+RVHVehs7U0/3Vm87mlqt5nA4+B5VpVaruVwsAlXC6iCMiIiQSqX1LQjlcjmXyxUIBPpupN6RSCSmpqb67qLeaTBvy6y5C67w2pCmHStVMx96ZF9KuH75bUeTy+V8Pp/PZ/VPXLUkEomJiYm+u6hf8CkBgHqjuRdpF1CpknKFZF/SUzfAFlhABgAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgtYa81+i8efMuXbqUkJCg70YA4F8SiWTt2rVKpZJRT09LI1UPF8x7dv/OLSNzK0aZx+XEXo1u06aNztoEFmmwQRgZGSmRSIqLi/XdCABU8uzZs80/7674aC6jzs3OrWbqijJO4+bSxRcYZfOfBxUVVXP0PcA7aJhBKJFIvv/++7179547d07fvQAAk4G5dYX/l4wiN3q3utqpuTxiZMmocbgN87cL9KJhfphWrFixZMkSnD0BoCFpYc4x+98Jl0zd3Wr+gnO5XJxfrVo4v1pVDTAIHz16FBkZ2bNnzxMnTkgkkvDw8ICAgNc/DQBq26y5C37ZtrVqnWfr/A6j3RrG5xBCEUII4fTbZiQS1jCxgYHBO7wEG4hEIn23UO/o5E+DgoKCoKCgrl27dujQYcqUKampqVWnkcvly5cvb9++fb9+/aKiot5k2MLCwm3btk2bNi0wMFCt/nclikKhWLFiBT1UZGSkiYnJ2LFjnz59mpqaqlQq09LSamu+AOCtZOUVkKm/kV3ySv+bfZxSU+8wGp9DTAXETEDMBMTUxLhhnGq87mFBuSqdBCFFUb179/7111/379/P4XCqXSBbvXp1VFTU4cOHg4KChgwZkpWVRder7kumqWRkZNy4ccPW1vbYsWP0JRpoa9asuXTp0uHDh6dOnTp06FAOh7N06dKlS5fOnz/fzMxs1qxZuphHAABoGHSyatTa2nratGn07a+//trFxUUikRgZGWkmUKlUO3fuPHbsWJs2bdq0aXPkyJF9+/YtXbqUELJgwQKxWLxixQp6yk2bNiUnJ//666+EEHd39wMHDqSkpGzYsEEzlFqt3rFjx6FDh+ihjh49unfvXvrp5ubmNR87kZWVFR4erslggUAwe/ZsvZ+uXiaTcblc7UVeoMlkMqzvqqo+vy1qtUrfLUBDI5PJ3mp6gUDw2s2iOvzRT0tLy8/P3759+4gRI7RTkBCSlZWVm5vbpUsX+m7nzp0TExPp2998842vr69Kpfr666+3bNmye/fuiIiIGl4lJycnOzu72qE4HI65uXkNz1Wr1RKJRLMTtrGxsfaCJgAAsIEOg/Czzz579uxZWVnZ4cOHGQ/l5uYKBAJNOlpaWubk5NC3bW1tL1y40Ldv35iYmGfPnl26dEksFtfwKrm5uTweT7P/mPZQr2Vvb+/l5aVZeK0/cD3CasnlckNDQ313Ue/U57eFy8VmPKhluvi063A/2r///js5OfnAgQMff/xxdna29kPm5uYKhUIul9N3y8rKLC3/PU6oUaNGI0eOPHv27LBhw2pOQXoolUpVUVFR7VAAAAA10/n2MF9fX5FI9OjRo8aNG2uKYrHY0NAwOTm5bdu2hJDk5GRnZ2fNo1u2bDl69Ojdu3cnTJhgZGT09ddf1zC+nZ2dSCRKSUlp165d1aEAoKHKzs7et2+frjdnGBoafvrpp8bGxjp9FdAvnQRhSkqKiYmJWCymKOr333+XyWR04EVERKSmpk6dOlUoFI4ePTo4OHjXrl2pqamnTp3SHEERHBy8e/fuyMhIOzu7ixcvfvTRR4aGhkuWLCGEUBRVVFRUUlJCCCksLOTz+RYWFgYGBoGBgcHBwXv27ElLSztx4sTFixd1MVMAUK9ERkZ++9Nv6g5DdPoq3Jjfe/To0bVrV02lrKxMoVAQQkxMTOpgE8bhw4e9vLze8+/7pKSk+Pj4MWPGVFRU/PzzzwsXLnzbEaRS6c6dO+fPn08I2bNnz4ABAxwdHd92kNLS0r17937xxRdv+0Rd00kQJicnT5kyhcPhVFRUODo6Hj9+3MrKihBy48aN2NjYqVOnEkI2bNgwfPhwe3t7qVS6fPlyDw8P+rmdOnWKiIiws7MjhNja2kZERDx48IB+qKyszNXVlRBiaWnZsmVLa2vrlJQUQsi6detGjBhBD7V06dIOHTroYqYAoL4xcGpfMmyNTl/C/HEkozJ+/PgrV65YWFjk5uZ27Nhx3759b55SW7ZscXR0DAwMfPMGDh48aGdn955BmJiYuHv37jFjxkil0s2bN78qCH/44YdmzZqNGjWq6kMSiWTLli10EP74448tW7Z8wyCcN2/euHHjunXrRggpLS0NDg5mSxAGBATk5OTk5+eLRCLt/UWXLVumuW1vbx8bG5uXl2diYiIU/nuGiF69emkPZWtr6+3tTd82NTUtKCio+nJisTgmJiY/P9/Y2Fh7KAAAXVi8ePHy5ctlMllgYODixYv//PNPui6VShnnbVGr1aWlpZrd1x8+fFj1UGmpVCoUCqse515UVGRhYREWFqZdVCqVFEVVXRKVSqVcLrfqjiTl5eXaP8KWlpYvX758VXsPHjyoeqSBUqmsqKiwtrZOT0/XrqtUKrlc/trz1Ny8edPPz4++bW9v//jxY0bb1Y5QVlZWl+fI1OHOMtbW1oyjJqqysbGpreiytrZGCgJAnTE0NBwyZAi9yurUqVPNmjVr3bq1i4tLeHg4PcHy5cvFYrGXl5e9vf2ZM2eOHj165MiRDRs2uLq6zps3jxBy586dLl26eHh4ODo6/vDDD/Szxo8fv3jxYjc3t1atWmVmZg4YMIDe3FNWVjZhwoQmTZo4OTkFBgbSx31FRET4+PhMnDjR1dV1/fr12u1lZWX16dPHzc2tTZs2V69epYuFhYW2trb07SVLltDtOTg4nDt37tChQ8ePH1+3bp2rqyu9yGhlZbV27VpnZ+chQ4bk5eXRK+pooaGhrVu3bt68+eDBg+mL/Bw6dGjixIn0o0ql0srKqry8fM2aNXfu3Jk+fbqrq+v27dszMjKcnJzoae7evduhQ4fWrVuLxeIff/yRLs6bN2/58uVeXl5ubm4tWrR49OhRLf+bvUIDPNcoAEAdkMlkp06dat269cuXLydMmHDmzJmePXtGRUUNHjw4OTm5oqJiz549z549MzY2lslkZWVl1tbW58+fb9GiBb3TQ3l5+fDhw3/55ZeBAwcWFBT07Nmza9eu3t7eZWVlp0+fjo6OtrW1VavVxcXF9A72q1evzsvLS0tL43A4I0eO/Oqrr7Zv365QKC5fvvzHH38cOHCAcRaOBQsWtG3b9vLly8XFxb169aJjTK1W5+fnE0KePHly4MCBZ8+eGRkZabfXtm3bRYsW0SMUFhampqbSC4J5eXn0E2mxsbEJCQkCgWD06NH/+c9/Nm/eLJPJysvLNRMUFhZSFLVy5crw8PCVK1fS5xd7/vw5PYhSqRwzZszcuXNnz56dlpbWpUsXT09Pb29viUQSFhZ29epVOzu7xYsXr1mzZv/+/Tr+ZyQEV6gHAHhbmzZtcnV1bdSoUWFh4ebNmy9evNipU6eePXsSQry9vT08PCIiIszMzORy+Q8//JCQkGBoaGhtbc0Y5OrVqxwORyAQXLhw4fbt2x07drxw4Z/LLgYFBdHLbdorKkNDQ+fPn29gYCAQCBYtWnTq1Cm63rRpU3qjI2OtZnh4+IIFCwgh5ubm9J4Z2szNzaVS6Q8//JCYmFhte7RFixZxudyq60tnzZolFAp5PN7cuXNDQ0Pf4r0jhBCSlJSUnZ09c+ZMQoiTk9OYMWM0szN27Fg6swcMGPDw4cO3HfndIAgBAN7OnDlzbt26lZube/nyZWdn5/z8fO0gsbGxycvLs7KyioyMzMrKCggIaN68eVxcHGOQ7OxsmUx27H9MTU3d3NzohzRrL7UVFBTY2NhoJtAsnzVq1KjqxHK5nF7I07TEmMDGxiYiIiIjI2PgwIEtW7a8e/dutXNabSeEEHr/R0KItbV11V03XntMS0FBgZWVlSZfbWxsNLNjZmZG3zAwMNAca65rCEIAgLcjFAotLS01p3h1dXW9d+8e/euvVqvv3bvXvHlzQkiHDh127NiRlpY2fPjwTZs2EUL4fL5K9c/5V1u1aqVWq4ODg3f+z/jx42t4URcXF835IxMSEuiXeBUDAwNHR8f79+/Td+/du1d1Gk9Pz507d6anpwcEBGzevJnRXs00I9+/f5/emd/CwkKTiPT+/DSBQFB1TBcXl4yMDM30CQkJ9CD6gm2EAADvZeDAgcuXL1+wYMGoUaMOHTpkbGzs5+eXmJh48uRJb29vHo8XHx/fu3dvQkjbtm3379/v5OTk4uLSvXv3bt26BQYGzp8/XyAQxMbGdu7c2cfH51WvsnDhwgULFlhaWvJ4vGXLlv3nP/+puau5c+fOnTt38+bNGRkZf/zxR8uWLbUfjY+PDwsL69OnD4/HS0hI8PX1pdv7448/mjRp4urqqn3oZFU7duxwcnISCoUrVqz4/vvvCSHdu3efNGnSb7/95ujouG3bNs2U7u7uISEhUqnUw8PD1NSULjo6Og4bNmzy5MlffvnljRs3IiMjt2/fXvPs6BSCEAA+VJRSTiSFun0JFfNohxEjRmh2faQJBILo6OhNmzZt3ry5TZs2ly9f5vP5jRs3VigUW7du5XK5Q4YM+fzzzwkhn376qUgkSkxM5HA43bt3P378+K5du3bt2qVUKtu1a0dn1YgRI9zd3TWDjx07tlmzZoSQUaNGGRkZHTlyhF6OHDp0KCGkWbNmkyZNqrbzRYsWmZqabt++3cXF5eDBg0lJSYQQkUhEbzi0s7OTyWRbt27l8XgjRoyg25s5c6aJiUl8fDyXy+3atevSpUs1u+IbGRnRe7oSQqZNm+bl5fX7779nZWWtXbv2k08+IYSIxeJTp07t3LlTKBSuW7fO3d2dXmJet27d3r1779y5Y2tra29vP2fOHHqQAwcOBAcH//TTT2Kx+MaNG/TZNP39/TVrcZs2bTplypR3/Gd7Sxw2X29h+vTp9fCk2/RlmHDS7apKS0s1f1GCRn1+W0aMm/yXyId4TahUvRvK/WOhen0KY2L+ylbKMVtIu8qXL730C+9qiOqrG4SQ8iCBkdaf7pGRkYOGjdD1D5hQKIy9ernm9ZDwocMSIQB8kHx8fMqKqjnDBsDbws4yAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWw+ETAPBhOPhY3cz0fxftS4kmjh5EZFFpirJcS0lmp47tCSFqtZrD4VS9yB+oVCoej6fvLuoXBCEAfBg+u/LvKSt5Sz5RzfiDuHpVmiL+qvfjvZFnThJC5HI5n8/n8/ETxySVSuvymrcfBKwaBQAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNexSBQDvq6ysrOdHfmVSOaOenZNNBvnooyOAt4AgBID3VVJSkpzyuGJOKKPO2zWh2ukB6hUEIQDUAi5fQJw8GUWOwEAvzQC8FWwjBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWQxACAACrIQgBAIDVEIQAAMBqCEIAAGA1BCEAALAaghAAAFgN1yMEgIaiJOdhUtLYqTMIIWqVisPlcDhcQohIwA/etN7U1FTf/UE9hSAEgIYiPy1Pzv1D2eGfu6p//is6tfrL+XNat26tr76gnkMQAkDDwbW0V/eZzigKorbppRn4UGAbIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWQxACAACr4eoTAPAWfvl116ZtOxlFlVIhk8n00g/A+0MQAsBbSEi8/9TBm3QJrFTNesQ5NE9PHQG8LwQhALwlSwfi5FmpwsFGFviA4eMLAACshiAEAABWQxACAACrIQgBAIDVEIQAAMBqCEIAAGA1BCEAALAaghAAAFgNQQgAAKyGIAQAAFZDEAIAAKshCAEAgNUQhAAAwGoIQgAAYDUEIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWQxACAACrIQgBAIDVEIQAAMBqCEIAAGA1vr4bAADQrYrykmlfLDY2NmbUx48aEjRxgl5agnoFQQgADZxCJo2x+YhYOFSqPoqy+PsighAIghAAWMHdj4jdKlXk5UR+RU/dQP2CbYQAAMBqCEIAAGC1hhyEjx8/Dg8P13cXAABQrzXYbYQymezTTz99+vRpWlqavnsBAID6q8EG4erVq6dMmfL111/ruxGAD1V2dnZ5eTmjWFpaqpdmAHSnYQZhfHz88+fPFy9ejCAEeDcKhcKxSVMjW0dGvbwwlwx110tLADrSAINQqVR++eWXBw4c0HcjAB8wlUpFONySVUmMOm9tD730A6A7utpZRiaTxcfHx8XFSSSSV00jlUqvXbv2+PHjtxq5vLz86dOnjGJFRYVmqBs3bjx48GDChAnDhw/PyckZN27cO/QPAAAsoZMlwrNnz44fP97BwcHQ0DAtLe3w4cO+vr6Mae7fv+/n59e8efNnz575+fnt3r37tcMmJSWNGjUqKSlJpVIplUoej0fXHzx44Ofn5+LikpaW5uvru2fPnhcvXhBCioqK2rdvf+jQoVqfQQAAaDB0skTo7Ox8+/bthISEmzdvLlmyZObMmVWnWb58+ZQpU6Kiou7evXv69OkrV/45xcOtW7fy8/M1kxUVFd24cYO+bWNjs2PHjlu3bjGGWrly5aRJky5fvhwfH3/u3LnLly/TdYFAEBAQUPuzBwAADYhOlgjd3P49lVGXLl3WrFnDmKC8vPz06dM//vgjIcTKymro0KHHjh3r1asXIeT8+fNHjx69cOGCtbV1UVFR//79Bw0a1LVrV0KIjY1Nr169UlJStIeSSqVhYWEbN24khFhYWAwbNuzYsWPe3t6EEGNj459//rmGPhMSEk6ePPn999/Tdw0MDKKjo0Ui0fu/A+9DJpNxuVyBQKDfNuqhsrIyfbdQH+nobamoqHjFI9RbjPI2074VqjZGVioVLNwJtry8nKqVt+8DIRQKX/tzqvOdZX755ZeRI0cyihkZGWq12snJib7r7OwcGxtL316+fHlFRYWvr+/x48fHjRvn7e397bff1jB+ZmamUqnUHio6OvoNe3N3dx8xYkRgYCB918jIqFGjRm/4XN0xMDBAEL6Kqampvluoj3Txtrz6E8h5i1HeZtq3wqmNkfl8AQs/URwOx8TERN9d1C+6DcKNGzfGx8dfu3aNUZdIJHw+n8//59WFQqH24UqrVq2SSCRubm6zZs3atGlTzS8hkUh4PJ7mS8sYqmY8Hs/W1tbFxeUNpwcAgIZHh6dY27Zt286dO8+fP29pacl4yM7OTqFQlJSU0Hfz8/PFYrHm0cLCwsjIyG7dul2+fFl7e2G17OzsVCpVUVFRtUMBAADUTFdBuGfPns2bN1+4cMHRkXlALiHE1tZWex1mdHQ0vRWQEFJcXOzv7+/j43PlypWhQ4f6+vrm5eXV8ELW1taurq7VDgUAAPBaOlk1evr06c8++2zatGlHjx6lKwsWLDAwMNi4cWNsbOzx48e5XO6CBQsWLFigVCpv37798OHD8ePH01OuXLnS19d37dq1hJDvvvuOELJ8+fJdu3YRQmQy2datW+lc3Lhxo7Gx8dy5czkczoIFCxYuXEhR1N27dxMTE//8809dzBQAADRIOglCc3PzL7/8khBSWFhIV+idlDw9Pa2srOjKF198YWpqevDgQWtr6+joaHNzc7q+efNmQ0NDzVDfffedTCbT3C0sLOTxeEuXLi0uLlYqlXRx1qxZJiYmmqEsLCx0MVMAANAg6SQIe/bs2bNnz6r1fv36aW5zOJwpU6ZMmTKFMY12CjIqhoaG69evrzosh8OZPHny5MmT36tpAABgpYZ8PUIAAIDXQhACAACrIQgBAIDVEIQAAMBqCEIAAGA1BCEAALAaghAAAFgNQQgAAKyGIAQAAFZDEAIAAKshCAEAgNUQhAAAwGoIQgAAYDUEIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYja/vBgBAz65du/brb/sYRZVKpVar9dIPQB1DEAKw3fnz5/fdfEnaB1SqyiWEQhACKyAIAYAQJ0/SZ3qliqSI/LlMT90A1ClsIwQAAFZDEAIAAKth1SgAsNLzhL+uHDMKPckoCw0NUh8nm5mZ6aUp0AsEIQCwUkUpx91POmkHs/5tW6lUiiBkFQQhALAV34AYWTKLXGwwYh38kwMAAKshCAEAgNUQhAAAwGoIQgAAYLVKO8vMnj07Nze35iccPXpUl/0AAADUqUpBmJiY+PLly2qnU6vVqampddERAABAHaoUhJcvX652olOnTq1cuZIQ4uvrWxdNAQAA1JXXbCO8du2at7f30KFDRSLRqVOnLly4UDdtAQAA1I1XBuG9e/cCAwN79uyZnZ199OjR2NjYwYMH12VnAAAAdaCaIExNTZ0xY0b79u1jYmJ27tx579690aNHczicum8OAABA1yptI8zNzd28efPWrVtNTEzWrl07d+5ckUikr84AAADqQKUg9PLyevLkSY8ePRYtWmRmZnb16tWqT+jXr19d9QYAAKBzlYJQqVQSQq5du3bt2rVXPYGiKJ03BQAAUFcqBeHvv/9eUVGhr1YAAADqXqUg7NWrl776AAAA0AucaxQAAFgNQQgAAKyGIAQAAFZDEAIAAKshCAEAgNUQhAAAwGoIQgAAYDUEIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACsxtd3AwBQd/Ly8kpKShjFwsJCQsz00g9AfYAgBGARlxatiMicw+FoFyX5mcRvsb5aAtA7BCEAi1RUVCjWPCUGRpWqm/oTQumpIwD9wzZCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWQxACAACr4fAJAIB/UZQ6LS2tvLycUReLxSKRSC8tga4hCAEA/iUrL/cZPJrLq/TbqCgvnvv5pz9sWKevrkCnEIQAAP+iOBzp0qvErHGl6vmtFfIsPXUEOodthAAAwGoIQgAAYDUEIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWQxACAACrIQgBAIDVcGFegAYoPz8/IiKiap2i1HXfDEA9hyAEaIAOHz68ZF2w0LkDo65UqvTSD0B91pCDcOfOnZcvXz548KC+GwGoa2q1mmrrXzh6C/OBWGN9tANQrzXYIHz06FFYWFhCQoK+GwEAgHqtYe4so1arFy1atH79en03AgAA9V3DDMLg4OBhw4Y5ODjouxEAAKjvGmAQvnjx4tixY3379k1NTVUqlenp6fruCAAA6i9dbSO8c+dOXFzc48ePP//8c2dn52qnOXLkSHh4eKNGjb744gsnJ6fXjimTyW7evHnr1q2srKy1a9dyuf+m+NGjR8PDw21tbefMmaNSqTw8PDZu3CiXy0tKSn7++WesIwUAgFfR1RLhiBEjzp49+9NPP7148aLaCXbu3Ll06VJ/f39CSM+ePcvKyl47Znx8/MyZMyMiIjZs2EBRlKa+e/fuL7/8csCAAVwut0ePHra2tjt37ty5c+ePP/5oZWWFFAQAgBroaonw2bNnhBBbW9tqH6Uo6ocffti6deuwYcPGjRt3/fr1Q4cOffbZZ4SQVatWdejQYejQofSU4eHhMTEx33//PSGka9euiYmJKSkpoaGhjKF+/PHHESNGjBs3LjY29uDBgzNnziSEmJiY1HzsRGlp6c2bN01NTem7xsbG/v7+HA6nFub/PajVas3/gza1Wo23papq3xbtvxShVlAU1TA+fmz7HmmvO3wV/Rw+kZub++TJEx8fH/quj49PTEwMHYSBgYEDBgygKGrYsGFnz56dPn36qVOnahiqoKAgOTnZ29tbeyg6CPl8fq9evWp4bnZ29uPHj3Nycui7fD6/a9euIpHovefvvchkMi6XKxAI9NtGPSSVSnk8nr67qHeqfVvkcrnOsvBthtVZHOsy56sfW6lUSiQSHb5sXZFKpW+SDQ2GUCjk81+TdPoJwqysLD6fb25uTt+1tbW9desWfbt169ahoaEDBw68ceNGSEhIeHi4p6dnzUNxuVxLS0vNUNeuXXvDNpo3b+7l5TVt2rR3nQ+dEAgECMJqURRlYmKi7y7qnWrfFkNDQ52t23ibYXW2ekWX622qH1sgEDSYj1+DmZHaop+/C4RCoUqlUqn+OdtTRUWFkZGR5tH27dsvXrx43bp18+fPrzkF6aHUarVSqax2KAAAgJrpJwjt7e05HI5mP5rnz59rH/MXHh6+adOm/fv3//TTTydPnqx5KLFYzOPxnj9/Xu1QAAAANavTIHzw4MHly5cJISYmJn5+fgcOHCCElJSUnDx5csSIEfQ09HbBkydPTpw48b6O6j0AACAASURBVMKFC1988cWJEydqGNPIyMjf3//3338nhJSWlp44cUIzFAAAwGvpahuhv79/SkpKYWHhmDFjhELhxYsXnZ2djx07dvXq1XPnzhFC1q5dO3DgwGvXrj1+/LhPnz6aHWcyMjI02wVbt24dFhZ2/fp1+qGysrL27dsrFApCSMuWLS0tLekti2vWrPH3979+/fqTJ0+8vLx8fX11NFMAANDw6CoId+/eLZfLNXfp1ZXz5s2bMWMGXenYsWNKSsqtW7esrKzat2+vmXLq1Kna47Rr165du3b0bSMjo/Pnz2se0uwp1759++Tk5Fu3bllaWnbowLzuDAAAQA10FYSOjo5VixYWFtp3TU1N+/bt++ZjcrlcFxeXah9626EAAABoLDqaBAAAoCoEIQAAsBqCEAAAWA1BCAAArIYgBAAAVkMQAgAAq+nnpNsAAB+SwpcXki4HTZ/BKDeysd6wbo3er90G7wlBCADwOrlPkhUWSXLm+Tq4m+Z+/5/vDAwM9NIU1BYEIQDA63Ec3Emf6cziH/P10gzULgQhwIdNKpVqLkOm0TAuIQtQNxCEAB+2jl298vLyOZWvOa6oKOf0nvqqpwCANgQhwIetXFoh/y6eWNhXqv4yhqenfgA+ODh8AgAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWQxACAACrIQgBAIDVEIQAAMBqCEIAAGA1XI8Q4MNw4cKFYaMCKYpZl8oV+mgHoOFAEAJ8GLKysjht+pV/sp1R5yx11Us/AA0GghDgg8HhGxAjS313AdDQYBshAACwGoIQAABYDUEIAACshiAEAABWQxACAACrIQgBAIDVEIQAAMBqCEIAAGA1BCEAALAaziwDAPCO1DwDt05eHMJh1Nd/s3T06NF6aQneAYIQAOAdUXLps8HBhCfQLnKifk1KStJXS/AOEIQAAO+haUfCN6hUMRfrqRV4R9hGCAAArIYgBAAAVkMQAgAAq2EbIUC9c+7cudTUVEYxNjZWpVLpox2ABg5BCFDvfDZvcaZpC46JlXZR+SiGY99aXy0BNGAIQoB6h6KIPOAr4tC2UvX3OfzyPD11BNCQYRshAACwGoIQAABYDUEIAACshiAEAABWQxACAACrIQgBAIDVEIQAAMBqOI4QQG9kMll0dHQ19YqKum8GgLUQhAB6ExERMXzsJEPn9ox6aW6uXvoBYCcEIYDeqNVqUYuuRTNOMOrchfaUXhoCYCVsIwQAAFZDEAIAAKth1SgAQK3KuP/jrwl7j55klJs62F06G6aXjqBmCEIAgNpEleUXuvYr7D21UlUuebl9iJ46gtdAEAIA1DYLMXHyrFSRlempFXg9bCMEAABWQxACAACrIQgBAIDVsI0QoC5s3vpTxOUrjGJOVqZMLtRLPwCggSAEqAv7jxxPaORDxK0rVZ//yZfn6akjAPgHghCgrrTyJi17V6pk3Ccp1Zx0GwDqErYRAgAAqyEIAQCA1RCEAADAaghCAABgNQQhAACwGoIQAABYDUEIAACshiAEAABWwwH1ALUpLS1t0MgxMoWSUX/+MpP0rvYZwBaUWh0XF1e13qZNG5FIVPf9gAaCEKA2ZWRkpOWXl07Yzahztw/XSz9QXygq5DJp/3GfMcrSnOdb16+eMWOGXpoCGoIQoJbxRCbMi7ISwuFiMwS7qVWEZ1D45XVGWXRknlLJXH8AdQxfTgAAYDUEIQAAsBqCEAAAWA3bCAHe0eXLl2NiYhjF9PR0RZVdRgGgPkMQAryjVRt+jM4hPLsW2kVV+hN1RYW+WgKAd4AgBHhHFCGKrhMUHYdUql7ew730f3rqCADeBbYRAgAAq2GJEOA1Xr586eLaXC5jrvDkGFmQoAl6aQkaDkotkUgKCwsZZUNDQyMjI710xEIIQoDXKCsrEzZykn+TyKjzvmqNvWLgPSme3lr21b6vV6/TLlJqdWO7xumPH+mrK7ZBEAIA6A2lVqtHbZD1nVmpWvC8PLivnjpiI2wjBAAAVkMQAvxr5hfzOVW4ubmVlZfruzUA0BWsGgX4V35RCQnaRXpOrlS9d47aj4sDQJ2iKKrqHjQ8Hs/MzEwv/TRsCEIAgHqmvKAw+6W9kyujrKgo/+v4n4MHD9ZLUw0YghDYqLS01Mm1RYVMzqjL1YSM6qGXlgD+pajgGFlVbMlglE33TS4tLdVLRw0bghAauPT09NjYWEaxuLi4XKaQr05i1PkbfeqoLQCoNxCE9c7Dhw9FIlGrVq303Uj9IpFIoqKiBg0a9LZP3BK87df/nhPat9QuqqUlCoWCGFkyp/4AL59L6bsB+ICUlJRcu3bN399f343ULw0zCE+fPn3kyBGRSDRz5sz27dvru523c/DgQSsrq+XLl+u7kfolMTFx1apVNQdhSkpKSUkJo5idkyPtMl7qt6BSNfMBWdu71pvUDyQhayjzX2zfsy/0QhSj7unR5ssF895khLi4uHXr1iEIGRpgEJ4+fXrSpElbtmzJzs729va+c+dOs2bN9N3U26Eo/LbV5L///e/yVWuZ7xFFPb5/16J5B0a59OVT8rF7XbUGoEOK4twYa68YZeUPeeHL67tC3jAIoVoNMAg3bdr0zTffTJo0iRCSmJi4Y8eODRs26LspeEd+AUNTX7wkhFRUVOTm5LRs35kQUpiTld+8H+VT+WQcSjm551P45XXGCLx1veqqWQDdc+tLun1SqfIoKuvKr918+jMm5FLqJfPnNG3aVLuYnJwsk8l03eMHpwEGYWxs7JYtW+jb3t7eBw4c0G8/9V9qaqparWYUORyOs7Mzh8Nh1EtLS5VK5ik24+Pjc3NzGUWVSsXj8VxcXBj1tLS0qsdCXbt2bcPmHwlhvpxUIqEWniEG/5x9OIX+z8G5HDM74uRZaVIlrgIIrFSaKyP8G50WMcq836aOnbXY0MRCuygvLxYoyqseocjhcCwsLAhbcRrYWrjS0lIzM7P09PQmTZoQQk6ePLlkyZJHj6o/d22bNm3S09MNDQ3pu1wut23bttz3210i/UWGTKliFBVyOVHKNC+kUV5ebmxszCjm5eXxeDxLS+Z+HGVlZSYmJoxiRUUFn8/n8yv9QUNRVHmF3MjElDGxpLyMx+PxeJUmVioVRVIFVyBkTKyWS0VGRnweT7uoUsglUilfYMCYWK4iHJEJI8YouYSoVXy+gDGxQiblC41J5YhVV5RTNs5cGyfGxMrkqwKPAYTLpShKqVQJBHxCiPJxLMfcjmdbeWI1pYg/LegYwBhBlXyN2DjxrBwqVWVliuQYgQfzj2jVoyvErgXPvHGlGSnLVz1P5Lf2YfaWdJnr0JZralWpi4KX6vx0fgsv5sQPIrnNPLmiSn8BqLMfU5JiXrNOzInvX+Q19+IYVrr4gPrFfYpS85p4MCaWx/8taNOXU/kfRfXsNhGa8MQtGRMr7p4RtBvA2CdI+fg6x8KeZ1Np0YGoVYqEc4IOAxkjqJKvElsXnqW4UlVaqnxyg9/Wlzkjj65w7FryzBtpF6mSXFXGQ75bH+bISVEcx3Zck0qffHVeuroog9+8O3PkB5e4zTpzRZU+5OqsFKqijOfckRCiVKq4XA79dVbeu8hrwXw/Vc8TCYfLc2SuNlck/M1v8xGn8udW9SyOiMx5ds2ZE989LWg/kPFhVqXEECtHnnWTyh3LFfcjBO2ZG+dUyVeJrSvP0k67SEmKVc/i+O4fMef6UTRH7MYzs60018XZ6sxkvhtzm7fqYRSnaTuucaX3U5X/XP3slqGQeV0LpUwqFBpyGd9KlZKiKKGQ+eNQIa0wNjVl/JWslMsUCoXAgPnjUCGVGpuZM/7CVcrlPLVMVGVkPp/v5MT8EXhPI0aMmD17ds3TNLQgVCgUQqHw0aNHzZs3J4QcOXJk/fr1d+7cqXbikydPxsfHa/6Z+Xy+h4dH1WWgOlZQUMDj8czNzfXbRn2jUqlevnzJWM8DhJBnz55Vu+zOcrm5uUKh0NSU+ecgyymVyszMTHo5gSWaNWvm6so8NQFDQ1s1KhAIbG1t09PT6SBMT093dHR81cRDhw4dOnRoHXYHAAD1zod31NRrjRw5kt4uKJfLDx8+PHLkSH13BAAA9VdDWzVKCHnx4oWPj49YLC4oKLC3tw8LC6u6cQ4AAIDWAIOQECKXy2/duiUSiTp06IBtJwAAUIOGGYQAAABvqKHtLFOv5OXl3b9/39HRUXufJZVKdfHixYyMDG9vb+1T3hQVFSUmJtra2rq5uVUdKisr6969e5q7nTt31hz0k5ycfPXqVWdnZx8fnw9i8VehUNy/f7+wsLBv377a9WfPnkVFRdnb2/v6+vL+d9iGSqV6+PBhVlbWRx99VPXIluLi4ps3b2pXPDw8GjdunJGR8eDBA02xS5cu9X8vXIqiUlJS0tPTe/ToYWT0797t+fn5f//9t1Ao9Pf319TT09Nv3LihUqm6d+/+qt3N09PTL1261KhRo/79+2sOsFGr1RERES9evOjVqxe9Q1n9V+1Xg6Koy5cvP336tHv37q1bt6aLMpns6tWrL168aNq0qbe3d7Vfh/j4+ISEBDMzs759+2qOZ719+3ZBQQF928jIqEePD+AKJFKpNCEhgcPhdO3aVbseHx9/+/ZtNzc3L69/juFRq9W3bt1KSkqysbH56KOPqh4OoVAooqL+PW2b9m6WhYWFZ8+eFQgE/v7+VQ/fajgo0I1PPvnE0NDQ1NT0m2++0RTVavXgwYM7duw4ffp0a2vrM2fO0PU5c+YYGBiYm5vPnDmz2tEOHDhgY2PT73/u3btH148dO2ZjY/PZZ5+1bdt2/Pjxup6p9xcVFSUUCm1sbExMTLTrZ86csbKymj59uqen5+DBg9VqNUVRiYmJxsbGNjY2hBCJRFJ1tPv372vek169ehFCrl27RlFUSEiIra2t5qGHDx/Wzdy9s7y8PHNzc3pOk5KSNPWUlJRGjRqNHTt2wIABbdq0KSoqoijqyJEjNjY2I0eOHDNmjKmp6YEDB6oOGBkZaWlpOXXq1G7duvXr10+lUtH1UaNGtW/f/tNPP7W2tj558mTdzN37+OKLL+ivxueff65dnzp1auvWrWfMmGFra3vw4EG62KRJk549ewYFBbm5ufXp06eiooIx2uzZs5s3bz5+/Hg/Pz9bW9vExES63q9fPw8PD/rTMnnyZN3P1vvavXu3gYGBtbV1ly5dtOvbtm0Ti8UzZsxwcXFZsmQJXRwyZIi7u/vkyZN79uzp7Oz84sULxmi5ubkcDkfzfQkJCaHrqampdnZ2gYGBAQEBLVq0yM/P1/2c6QeCUFfS0tLkcvm4ceO0gzAyMtLe3r6srIyiqN9++61Tp050/fnz5xUVFfPmzashCD/++GNGUa1Wt2jR4tixYxRFFRYWWllZ3b59WyczU3uKiopycnJu3rzJCEJPT8/du3dTFFVWVubg4BAZGUnfzsjISE1NfVUQatu3b1/Lli3pBA0JCRkyZIjOZqL2yeXy1NRUlUrFCMLp06fPnj2boii1Wu3r67tlyxaKol6+fFleXk5PsG/fPkdHx6oD9u7dOzg4mKKoiooKV1fX8PBwiqJiYmJsbW2Li4spijp06JC7u7vu5+x90V+N+fPnawfh/fv3TU1N8/LyKIoKCwtzcnJSKpUURT1+/JieQCKRODo6Hj16lDHakydPNH8TTJ06VfO3Y79+/Y4cOaLrealF2dnZxcXF+/bt0w5CiURiZWVF/y2Ynp4uEolevnxJab0tFEX5+fktW7aMMVpubi6Xy636KnPmzJk2bRp9OyAgYO3atbU+I/VEAzx8op5o2rSpQMA8qUpYWNjAgQPps8mMHDkyLi4uIyODEOLo6PjaXVtLSkrOnj0bFxenOcNZUlJSWloafblqCwsLX1/fsLCw2p+TWmVubm5ra8sovnz58vbt26NGjSKEGBsbDxw4kJ4RY2NjsVhczSjV2bNnz9SpUzVrwzRvF50u9ZxAIKh2DWdYWBh9/A+Hwxk5ciT9ttjb22vWkYrF4qqnjiwuLo6OjqbfT0NDw0GDBtFPDAsL8/Pzo9cHDhs27NGjR0+fPtXlbNWCar8a4eHhPj4+1tbWhJABAwbk5+cnJiYSQjQr9EQikYWFRdV3xsXFRbOCXSwWy+X/Xpk5JSXl77//TktL09GM1K5GjRpVe55CkUhErxFt0qRJx44dz549S7TeFkKInZ3dq841GhUVFRUVpX39ltDQUPpTRAjRfPwaJARhnXr58qXmAH8zMzNTU9OXL1++4XNLS0t//vnnwMDATp060fGZkZFha2ur+ZlwcHB489HqlYyMDFNTU81mvHeYkadPn8bExEycOFFTKSoq+vnnn0ePHt25c+fMzMzabLeuKBSKnJwczQem6tuiUqnWrFkzffp0xhMzMjJ4PJ6dnR3jidofP5FIZGlp+YF+YLRnhM/nN27cmDEj//3vf/Pz8wMCmOfb08jMzNy1a9e0adPouyKRKCIiYuvWrR4eHgsWLHjVs+o57beFVPeBefDgwcmTJydPnlz1ufb29ps2bVqyZEmzZs3+/vtvQghFUZmZmQ4ODq8arSHBzjJ1SqlUam/A5/P5VU9gXa2xY8dOmDCBHmH06NErVqzYu3evSqXSHo3H473haPXN+8/I7t27Bw4caG9vT9+dOHFiUFAQIUSpVI4YMeKrr77as2dP7fVbR+jVvJp3hvG2UBRFrzX95ptvGE+k38+qT2S8z2/+8atvav4e3bx5c+bMmYcPH656wl5aSUnJsGHDgoKCBgwYQFf++usvev+sp0+fenp6DhkyhLEn1weh6r+v9uqQzMzMoUOHrl69uuolWq2trdPS0uhl5f/7v/8LCgrKyMigKEqtVr/q49fAYImwTtnb2+fk5NC3KyoqioqKNL/dNdPsRcnn80ePHk2fPVUsFufn52s+69nZ2W++IrFesbOzKy0traj45/IRbzsjKpXq999/nzp1qqbCeLvu3r1bi93WGUNDQysrK81lPbKzs7U/LQsXLoyPjw8NDa26E6CdnZ1SqdTsBql5P8Visebjp1Qq8/Pz3/DjV99of48oisrJydHMyN27dwcPHrx79+6PPmKesZpWVlb28ccfd+rUaf369Zqi5gPj4uLSrVu3V52duJ7T/vclhGRlZWm+Rzk5Of369Zs6deoXX3xR9YkcDkezxnjcuHFZWVmZmZlcLrdRo0av+vg1MAjCOuXj43PhwgU6us6fP+/i4lLD2W/VajW9mxajfvv2bfpZbm5u5ubm0dHRhBCFQnHp0qUP8c9YQkjTpk1dXFzOnTtHCFGr1RcuXKh5RkpLSyUSiebu2bNn5XL5xx9/XO3EcXFxH+4phn18fOj1VISQc+fO+fj40LdXrFgRGRkZHh6uvaFIKpWWlpYSQmxsbDw8POj3k6Ko8+fP0++nj49PREQE/Xf9pUuXGjdu/KEcQcHg4+MTFRVFb+u6ceMGl8v18PAghCQlJQUEBAQHB9MbzmkURWn+XpRIJEOGDGnZsuX27durPbhCKpU+fPjwAz23e7du3TIzM1NSUgghJSUlsbGx9AcmLy+vX79+gYGBy5cv156+qKhIeyspLS4uTigU0hvy+/btS3+KSOWPX8ODVaO6cuzYsQsXLty4cSMpKSkrK2vs2LE+Pj6DBg367rvvAgMDe/bsuXnz5tWrV9N/iJ05c+bEiRPXrl1Tq9UzZswYNGjQ4MGDnz592qJFi5ycHFtb208//dTU1NTe3j4+Pv7kyZPnz58nhAgEgiVLlkyZMmXu3LkXL150dnau/0FYWFi4bNmy3NxcmUxG7/v+/fffc7ncFStWfP75548fP6Y3+A8aNIgQolAo5syZU1ZWRgiZM2eOiYlJcHAwIWT8+PHu7u7r1q2jx9yzZ8+kSZO0d02aNm2ahYWFWCy+c+dOaGjoxYsX9TGvb2fx4sX0fgrffPONhYXF1q1bRSLRkiVL+vfvz+VyCwsLIyIi6AttHjx4cN26daNGjdL8rv3000+GhoYbN268du0aHZwrVqyYN28evReSVCqld3kYMGCAra3tyJEjvb29g4ODly9fzqt8ma166OzZs3/99VdMTIxSqdR8NXr06OHh4TFkyBB/f/+ff/558eLF9GKxj4+PhYVFREREREQEIWTIkCEBAQGlpaU2Njb37t1zd3dfsGDB9evXXV1dZ86cSQhp2rTpypUrs7Kyxo0b5+3tLRAI/vzzTwcHhyFDhuh5tl/nwYMHwcHBKSkpqampM2bM8PDwmDNnjqWl5axZs0aOHDlt2rTjx48PGjSoVatWhJBJkyZlZGRkZGTMmDGDENK+fftZs2YRQrp06fLNN99MnDhx9+7dV65cadOmTW5u7p49e1atWkV/mxYtWuTj42NoaFheXh4WFhYXF6ffudYdnFlGV2JiYug92Wg9e/Z0d3cnhJSUlOzduzc7O9vX11ez9ubOnTvaB4Z36tSpU6dOJSUlR44cmThxolAojIuLi4iIyM/Pd3BwGDlypPY6ijNnzkRHRzs6OgYFBWkfiF0/lZWVHTp0SHPXzMzsk0/+udz2xYsXIyIiGjduHBQURC/oqFQq7W17BgYG9Ja/c+fOWVpadunSha6HhIT0799fezeBW7duRUREFBQUVH276q39+/drVg4TQoKCggwMDAghiYmJx48fNzQ0nDBhAr1om5iYGBMTo/3cKVOmCASC27dv5+Tk+Pv/c8W76Ojov//+29raOigoSLO1rKysbO/evZmZmd7e3n5+fnU0b+/h7t27N27c0NylvxqEkIqKir1796anp3t5eWmW/3bt2qX9g9a5c2dPT0+FQrF3797Ro0dbWFicO3eOPhqHZmNjM2LECIVCceLECfrYXHd395EjRzIu8FkPvXjx4vTp05q7TZs2pf/dKYr6888/4+LiWrVqNWHCBDrP/vrrL+3rZjs5OdEbR48ePerp6dm8efOXL1+Ghoamp6fT5xno1q2bZuIHDx4cO3aMz+ePHz/e2dm5zmawjiEIAQCA1bCNEAAAWA1BCAAArIYgBAAAVkMQAgAAqyEIAQCA1RCEAADAaghCAP148uTJr7/+Sp8LpuGJiYlhXDP5/SUmJkZGRtbumAAEQQjwhtq1azd37txaHPDmzZszZszIz8+vxTFrtmHDBldX12XLlun6hbKzs/39/QsLC+m7x44dc3V1dXV1ffTokfZkBQUFLVu2dHV13b59+5sMq1AoBgwYkJycXPsdA7shCAHeyLNnz7Kzs/XdxbtTqVT/93//l5GRsXPnTqlUqtPXWrVqVevWrTVnriktLX369GlGRsbevXu1Jzt06NDz58+fPn2qicyaeXp69uvX76uvvqr1hoHlEIQAtaOsrCwrK+tV1wGWyWRZWVlVz3H8nqRSaXZ29pucH+r8+fPPnz/fsmVLUVHRiRMnqk6gUqmys7O1T/PGoFQqs7KytE93Xq3i4uL9+/dPmTKFUR82bNiBAwe035+9e/cOHz682kFycnIUCkXVelBQ0F9//dWAL4wHeoEgBHhfCQkJvXr1MjMzE4vFjRs3Xr9+vXYypaWlDR8+3NzcXCwWC4VCLy+vai8UnJqa6unp6evrW1BQsH79es0FUWnr16/XnEw1Ly/Pysrql19+GTdunKmpqZ2dnZubG30RkhqEhIQ4Ozt//vnnHTt2DAkJYTx69OjRJk2a2NnZWVhYzJo1a926dZrr+hJCpFLp/Pnzra2txWKxmZmZv79/enr6q17o+PHjEomkasJNmjQpKyvrwoUL9N179+7FxcXRJ4/VSEpK8vPzEwqFjRs3NjIy6tSp05UrV7QnGDRoEJfLPXjwYM0zC/BWEIQA7yUjI6Nv377Z2dmhoaG3b98eN27c8uXL16xZQz+al5fXu3fv69ev79y58969e5GRkX369KGvH6QtISGhd+/e1tbWf/31l5WVFX2tSu0JtCsURRUWFn799ddcLvfGjRuRkZEmJiYff/xxWlraq5osKCg4derU5MmTORzO5MmTL168qD3x1atXx40b16lTp+vXr0dHR2dmZm7ZskWzupKiqNGjRx84cGDjxo337t0LCwtLT0/38/N71bJjRESEi4tLo0aNGHUHB4d+/frt27ePvhsSEtK2bVtPT0/tafLz89u3b3/69OmHDx+eOXPGyMgoICBAe420SCTq0KHDpUuXXjWnAO+CAoA3YGJiEhgYWLW+YsUKLpd7//59TSUgIMDU1LS8vFzzaFxcXNUnHj58mBDy7Nmzc+fOmZmZTZ48WS6X0w99++23RkZG2hN/++23xsbG9G364qvu7u70Jewpinrx4oWBgcHChQtf1fxPP/1ECHn06BH9dIFAsGrVKs2jw4YNE4vFUqmUviuXy5s0aWJgYEDfpZfh/vzzT830Dx8+JIQcOXKk2tfy8PAYOHCgdoW+hEh8fPyhQ4eEQmFhYaFCobCzs9u8eTN9VYT//Oc/1Q5VUFDA5/N3796tXZw0aZJYLH7VnAK8AywRAryXu3fvenh4tGnTRlP55JNPSktL6aucnz9/vkOHDozlHm179uwJCAiYP39+SEiI9iUVX2v48OGaS8s6ODh07969huuqh4SE9OrVq2XLloQQW1vbgQMHhoSEqNVq+tHExMT+/ftrrnQvEAgGDhyoee65c+e4XK5QKLzwPy9evLCwsLh37161r5Wbm2tlZfWqnoVC4ZEjR8LDw/Py8saNG1d1mpycnG3bti1cuHDGjBnLli0zNDR8/Pix9gTW1ta5ubkULpsDtae+X3YLoJ5LT0/XvhQiIYS+/CF9XEReXl67du1qePq2bdvs7e3nzp1b7QXTa6C9DY8QIhaLb9++Xe2Ud+/evXPnzvDhw3/99Ve6YmZmlpqaGhUVRV/JOSsry8bGRvsp9AXKafSayYkTJ2pPwOFwiouLYl8ovgAABPxJREFUq305Pp//qj2GhELhmDFj9u3b17hx44CAADs7u7y8PO0JTp8+PWrUqKZNm/bq1cvKyorL5fJ4PPp6xRoKhYLP57/t2wVQAwQhwHsxNjZm/JrTq/vMzc0JIRYWFllZWTU8/eTJk9OmTevdu/eFCxc0FxDm8XhKpZKiKM3PPSMMCCFVX5QRjRohISFcLjcyMlL7aHQ+n//bb7/RQejg4MBoUnt3HnNzcz6fn5mZaWhoWMOMaDRu3JjRm7agoCAvLy8+n3/06NGqj65du7ZTp06RkZE8Ho8Qolar6ZW62vLz8xs3bvwmnQC8IawaBXgv3bt3j4+Pf/78uaYSGhoqFAo7duxICPH29r5z5w5j5Z42Jyen6OhoLpf70UcfvXjxgi46OjrK5XLNXYqioqKiGE88e/as5nZBQUFsbGzbtm2rji+Xyw8dOjR8+PCCyj799NPjx4/TO+B06tTpzJkzmqwtLy8PDQ3VjODt7S2Xy6s94qJaXbp0SUhIeNWj3bt3/+STTwYMGBAQEFD10WfPnnXs2JFOQULIpUuXqh7yGB8f37Vr1zdsBuBNIAgB3tTz58+PVJacnDx79mw+n//JJ588ePCgoKBgy5Ytf/zxx4wZM8zMzAghCxcuNDY2HjZsWFRUVHl5+fPnz3fs2ME4DE4sFkdERBgaGvbq1evJkyeEEB8fH4FAsGDBgvT09KdPn86cOTMlJYXRTEJCwldffVVQUJCamjphwgS5XD579uyqPZ84cSIvL2/8+PGM+vjx46VSKb1YtmzZMqlU2r9////+979//fWXn5+fSCTSTDl06FAvL69Zs2bt3bs3Nze3pKTk9u3bK1euvHHjRrXvUv/+/XNycugZqdbhw4fDwsIMDAyqPtShQ4djx47dunVLJpNdvHhxxowZmi2XtIKCgqSkpH79+r1qcIB3oeeddQA+ECYmJlW/Pps3b6Yo6vz5805OTnSFx+PNnDlTJpNpnnj79u0OHTpontKkSZP09HRKa69RerKCgoJu3bo1adIkOTmZoqhffvmFXhXJ4XDGjh27bNkyxl6ja9as6d69Oz2mhYXFsWPHqm3b39/fzMxMs0eohlqtbtasWffu3em7kZGRnTp14vP5Dg4Oa9as+fLLL+3t7TUTFxQUTJgwgc//d0tKp06dtHeU1SaTyRo3brxmzRpNRbPXaNWJGXuNJicnt27dmn4JU1PT/fv3Ozk5zZo1SzP9jh07jIyMCgsLq31pgHfDobDzFcAboI9VYBS5XC69GU+tViclJZWVlbVo0cLS0rLq058+fUofCO/q6vqGO3qUlZUlJyc3atSIsTNObm5uo0aNdu3aNX369JSUlKKiorZt22ovw70/Pz8/iUTCOJi9pKTk0aNHBgYGjo6O1tbWNTx99erV+/btS0pK0s7ON6RUKp88eVJeXt66deuqM9WlS5cePXoEBwe/7bAANdF3EgPA26GXCHft2lWLY549e5Y+ilGtVu/YsYMQsmnTpnceTSKRNG3aNCQkpNb6oyiKokJDQy0sLOhjJwBqEfYaBQASFBRUWlrapEmTzMzM4uLisWPHvs+lNkQi0a1btzTHKdaW7t27P3z4kHGkB8D7w6pRgA9MRUXF/v37+/Tp4+bmVltjlpSUXL9+PSsri8fjde7cuVWrVrU1MkD9hyAEAABWw+ETAADAaghCAABgNQQhAACw2v8DZHmlmA3ednQAAAAASUVORK5CYII=", "image/svg+xml": [ "\n", "\n" ], "text/html": [ "\n", "\n" ] }, "metadata": {}, "output_type": "display_data" }, { "ename": "LoadError", "evalue": "UndefVarError: `fit` not defined", "output_type": "error", "traceback": [ "UndefVarError: `fit` not defined", "", "Stacktrace:", " [1] top-level scope", " @ In[5]:24" ] } ], "source": [ "# The fraction of the way through the MELTS dist that lockup occurs\n", "lockup_fraction_mu = 0.2814311218465038 \n", "lockup_fraction_sigma = 0; #0.07475159714446021\n", "\n", "lockup_fraction = lockup_fraction_mu .+ lockup_fraction_sigma .* randn(size(tminDist))\n", "lockup_dist = tmaxDist .* (1 .- lockup_fraction) .+ tminDist .* lockup_fraction\n", "\n", "# Print results\n", "AgeEst = nanmean(lockup_dist)\n", "AgeEst_sigma = nanstd(lockup_dist)\n", "AgeEst_025CI = nanpctile(lockup_dist, 2.5)\n", "AgeEst_975CI = nanpctile(lockup_dist, 97.5)\n", "print(\"\\nEstimated lockup age:\\n $AgeEst +$(AgeEst_975CI-AgeEst)/-$(AgeEst-AgeEst_025CI) $age_unit (95% CI)\\n\\n\")\n", "\n", "# Plot results\n", "h = histogram(lockup_dist,nbins=100,label=\"Posterior distribution\",xlabel=\"Lockup Age ($age_unit)\",ylabel=\"N\",fg_color_legend=:white)\n", "# plot!(h,[wx,wx],collect(ylims()),line=(3),label=\"Weighted mean age\",legend=:topleft)\n", "savefig(h, \"histogram.pdf\")\n", "display(h)\n", "sleep(0.5)\n", "\n", "# Print out histogram\n", "edges = range(minimum(lockup_dist),maximum(lockup_dist),length=101) # Vector of bin edges\n", "hobj = fit(Histogram,lockup_dist,edges,closed=:left) # Fit histogram object\n", "print(\"Histogram N:\\n $(hobj.weights)\\n\\nHistogram bin centers:\\n $(cntr(edges))\")\n" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdZ1wURxsA8LlC70Wa0nuTolSlKUSjFAugIHZFfW1YEhRRIYq9lySWCJZoghVrFAVUVKJgUJEqKkVADjjqcVzb98MmFwQ0Uo4F7vn//DA77M49t5w8N7uzMyQMwxAAAAAgrMhEBwAAAAAQCRIhAAAAoQaJEAAAgFCDRAgAAECoQSIEAAAg1CARAgAAEGqQCAEAAAg1SIQAAACEGiRCAAAAQg0SIQAAAKE2oBIhk8mMiIjgb3K5XAKDESoYhvF4PKKjEBbwwe5N8MHuNQR+sAdUIqytrY2NjeVvMhgMAoMRKjwej8lkEh2FsIAPdq+BD3ZvIvCDPaASIQAAANBZkAgBAAAINUiEAAAAhBokQgAAAEINEiEAAAChBokQAACAUINECAAAQKhBIgQAACDUIBECAAAgxqVLlxT/oaSkJC8vj5fnzZvXm2FQe/PFwEBFIpHIZPhSBQDoHG9v78LCQrxsZ2d39uxZQ0NDhJCoqGhvhgGJEPQAMpksJiZGdBQAgH5GVFSUn/MoFIqsrKyCgkLvhwHf4kHPIJFIRIcAAABdAYkQAACAUINECAAAQKhBIgQAACDUIBECAAAQapAIAQAACDV4fKK76urqkpOTi4qKpKWlhw8fbmVlRXREAAAgWJWVlY2NjXi5oaFBRkYGL6uoqEhLSxMXVxdBIuw6Lpf7w+aY7Tt2UuTVkKohidnQ8m6ZucXQM7HHLCwsutNyQ0NDbm6upKSkqakp/0H1Z8+e0Wi0cePG9UTsXbdr165FixZJSUkRGwYAgEDbtm1LSEhACHE4nIqKiiFDhuD1e/bs8fPzIzS0riBhGEZ0DD2moqLCxsamvLwc32z9PUUQps2cfSnxIXPGcWQw4u+q5nrKzS1ij355/PB+17qGPB5v/fr1hw4dMjU1pdPpbDY7NjbWzc0NIbR79+7Xr1+fOHGiB99CF4iJiRUXF6uqqhIbhjAT9Acb8PF4PCaTKSkpSXQgfVdBQcH48ePz8/O735SxsXFCQoKJiUn3m+osAd4jpNFoHz58aF+ZkpLSpr6oqOjBgweZmZlsNvtzrZWWlqakpFRVVQkk1s67efPmhctXmSvv/ZsFEUISstzJ25gj54fMmd+1Znft2nXmzJlnz56lpaXl5eV9//333t7e79+/b70P/4oEH4fDaWhoaFPJZrM7PJ+NjY0cDucLMeD/+du31tLS8uWjmpqaOny5Dvdvbm5uU9PhuwAAAEETSCK8ePGipqamqqpqmz7ypUuXTE1Nt23bZm1tvX//frxy69atrq6uUVFR06ZNMzY2zsvLa9/g7t27bWxstm3bZmJigvfHCffj0V/YbouQvEb7H/HGR+Rmv37x4kUXmt2zZ09UVJSRkRG+uXDhQltb259++gnfrKqq8vDwsLKy0tHRuX//PkKIy+XOmTNn8ODBDg4OQ4YMyczMRAhVVlZOmDDBwMBAT09v1qxZeMo5ePBgSEiIu7u7kZHR0aNHtbS0+Kno9u3bI0aMQAgxmczFixfr6OiYmpp6eXnx+9Zr164dPHiwhYXFhg0bOgw7MDBQS0vLysrKyMjowYMHeGVeXp6VlZWlpaWVldWyZcsWLFiA11+8eNHQ0HDo0KF6eno3btzA38W8efP47yIjI6MLpw4AALoIE4A3b968fPkyLi5u2LBh/Eo2mz1kyJCrV69iGPbq1StJScmqqioMwxgMBn+fGTNmTJ8+vU1rHz9+lJCQyM7OxjDs8uXL2traHA6nw9ctLy9XU1Pjb9bX1/fce2pLTUsPhd1Ax1gd/pM2cY6Li+tsmzQaDSH05s2b1pXr1q1zd3fHMGzXrl2ioqLPnz/HMOzKlSvq6uoMBuP+/fsWFhZsNhvDsMbGxrq6OgzDfH19169fz+Px2Gx2YGBgTEwMhmE7duyQkJD466+/MAzjcrlOTk5nz57FXyIwMHD79u34awUEBDCZTAzDoqKipkyZgmHYtWvX9PT0aDQahmHh4eEIoYqKijaR478dDMP++OMPPT09vDxixIitW7fi70tXV3fatGkYhr1+/VpdXb2goADDsKysLFVVVRqNlpqaamZmxmKx8HdRW1vb2VMnVAT6wQatcbncpqYmoqPo0/Lz8w0NDXukKSMjo5ycnB5pqrME0iPU19e3tLSkUCitK9PS0phM5vjx4xFCFhYW5ubm169fRwhJSEjw91FWVm5zFELo2rVrVlZWpqamCCEfH5/6+vpnz54JIuxOYTYzkNhn7xxgYlIdXif8straWoSQrKxs60pZWVl+183d3d3GxgYh5OfnJykpmZ6erqKiUlpaum/fvoKCAikpKVlZ2cbGxmvXrtna2t67dy8lJcXKyioxMRE/3NPT09raGiFEJpNnzpx58uRJhFBdXd3NmzeDg4MRQufOnXNycnr48OHdu3cNDQ3xA69fvz5jxgxlZWWE0MqVKzuMfNCgQT///HNkZCQ+gLampqahoeHx48dLly5FCCkrK+PtI4TOnz9vbW39/v37u3fvlpeXq6urp6WlqaiolJWV7d27Nz8/X0pKSk5OrrOnDgAAuqz3Ro2WlJRoaWnxx0Dq6OgUFxfj5dzc3L179378+LG+vv7UqVPtD9TR0cHLFApFU1OzuLjY0dGxw1dhsVj8v/sMBmPo0KHa2tpfDqxr6wdpDNGspb375AZhK9jHN1paWp1tU1tbm0qlFhYWDho0iF9ZWFjIb0pJSYlfr6SkVF1d7eLicuXKldjY2O3bt2tqal68eBEhRCKRbt26xd/T3d0dL7RudurUqatXry4tLb1+/bqrqys+6KuysvLJkye5ubn4Pv7+/jwer6amxtzcHK9RVFRsP7l2VVWVjY3NzJkz8fFBZDK5oaGBSqVSqVT+KAN+bvv48WN5efn58+fxTXt7e0VFRUNDw6tXr544cWLXrl0aGhqXLl3S09Pr7NkTEhiG8Xg8Ho9HdCBCgfcPogPpu/CT01OnSBBn+2v+wvdeImQyma2XmBITE+N3dGRkZIYNG1ZaWhoXF5eRkcEfifufB7bBYrGampq2bt2Kb3K53ODg4GnTpn05MElJyS7kwgDf8dsvnGI6hXTws/yH3LpKfvr5eiIiIh4eHocPH+an+erq6t9///3IkSP45uvXr/ECk8l88+aNvr4+QsjNzc3NzY3D4UyfPv3IkSObN28WFxdfsmSJpaVlm/Zb5zA5OTkfH5+zZ89evnx5xYoVeKWxsbGfn1+bM6anp5ednc0PAGs3zPjx48e6urqbN29GCL158wYfoaOuri4pKfnixQu8D/rs2TP8l2hiYvL27Vv+O+JzcXFxcXHhcDizZs366aefdu7c2dmzJwwwDGtpaWlubm5/4QQIQocDx0BrTCaTx+MxGIzuN4Wf7R5pqjVxcXEq9T8yXe8lQjU1terqav5mVVWVvb09Xh48eHBoaChCyMjIaO3atW2G2KipqWVlZbU+UF1dvcOXEBUVVVBQSEpKwjcFOsp8+fJlB378mXVtE887ErXuJNHeiZ+cExm5rmtPlR46dMjV1XXGjBmBgYE1NTVbtmzx9PQMDAzEf1pSUhIREeHt7X3kyBFra2tLS8vk5OTMzEwHBwc2m52Xlzd27FgqlRoVFRUUFLRp0yZ1dfXc3Fwulzt37tz2rzVr1qxZs2YxmUxfX1+8Jjo6OjQ0lM1mm5qaFhcXFxQUREREhIaG2tnZDR8+3MDAYOvWre2/NBgYGLx8+fLKlSuysrJbt24VERFBCJHJ5MjIyKlTp65cuTInJ+f58+dOTk4IoTlz5uzfv3/VqlWTJk1qbm6+e/fu8uXL8/PzMzIyHBwcuFxubm7uqFGjunDqhAGJRBIXF5eSkuqPzyz3Rzwer/WFDdAe3pHokQ8kmUyWlJQk5LPde4nQxsamuLi4rKxMQ0ODxWKlpaVFR0e32afD77n29vabNm1is9kiIiIlJSVlZWX4fTJiKSgo3Lt985txPk15d5sdZiI1I9RcT8m/T354fMbMGRFr13StWSMjo8zMzMOHDx89elRaWnrdunXBwcF4T87e3v7IkSM0Gm3Hjh0mJiaHDh1CCGlra9++fRsfCBMWFjZjxgyE0HfffYcPr6XRaIaGhiEhIQghBweHNtcbR48ePWfOHBMTE3Fxcbxm/Pjx+IXWy5cvDx48eNKkSQghPT29P/744/Dhw0lJSZs2bRo2bFibp+nNzMxiY2PPnDkjJSW1devWa9eu4RdCV69ebWZm9vjxY0dHRzk5Ofw5Cmlp6YyMjP379+/evVtKSsre3l5eXl5bW/uPP/7YuXOnhITE0qVL8XcBAAC9RBAjcN69e7dt27YpU6YMHjx427Zt8fHxeP2sWbM8PDyuXLkyZcoUFxcXvPK77747cODAhQsXdu3apaqqeuDAAbze1taWP6zRyckpKCjoypUrbm5u8+bN+9zr9uaoURydTo+J2WI53FFBRX2wruHkqcEpKSmCftH+4sOHD/ivoKSkRE9P79atW0RHNEDAqNFeA6NG/9PAGDUqkB4hh8Oh0+k6Ojo6Ojp0Op3/SPXPP/984MCBs2fPGhkZ8e8SeXp63rhx4+HDh2pqar/++uvo0aPx+pCQEP4UAzdu3Ni9e/fZs2d9fHzwgYh9hLy8fETE2oiItUQH0hdlZ2fjF1qpVOqCBQvGjh1LdEQAANABmGINgH4GPti9BqZY+08wxRroYSwW6+0/enzoVBtZWVm///57d1pITk6+d+9epw5JTU3t7CF9WXl5+eHDh4mOAgDQXZAI+5Ds7GwDAwMvLy8vLy8NDQ0vLy/8EXtBKCws7GZOunPnzu3btzt1SGJiIj6LQv/l6+ubk5ODl6uqqi5dutSd1srKyvAZ1QEABIJE2LfgD9QXFhZWVFTU1dXxZxlFCDU1NdXU1HxlO/X19SwWq3UNh8MpLy/nz3fj5+d39OhR/k8ZDEbrh1v4WCzWV050jt8Ybl9Pp9O73LvFMKyysvILV+8bGhraf1fg8XiVlZXtd+ZyuV8+gTwer/1JqKurq6ys5HK5+GZWVhb/7VhaWrb5MtHc3IzPk4cQYrPZ7Vurqqpq/RQsi8XCp4dt7XNnEgAgIJAI+yhxcfHBgwfjCz6wWCxra+thw4aNGDHCzMwMf6ry5MmT/v7+/P0XL168bds2hFBubq6jo6ODg4O+vv7ChQvxv+BxcXFaWlo+Pj4WFhbr169HCP3222+TJ0/Gj7W1tbWxsXFxcTE2Nv7rr7/wSgcHh8jISEtLSxsbGycnp88tIoHbsGGDhoaGk5OTubk5v4X79++bmJg4OjpaWFiEhYW13r++vn7cuHHR0dHV1dWt57sJDAy8cOECQmjHjh3BwcEODg6jRo3S1NTkT+TNR6fTLSwsHBwc7O3tbW1t37x5g9cnJCRoaGh4enq6urr6+/sfO3YMIcTj8datW6etre3u7m5mZtY+9yCEtmzZoqWl5eHhYWRk9OTJE4RQbW2tm5ubnZ2dl5eXjo4OjUZbt25dWVnZ9OnThw8ffuXKlYyMDP7Ck5qampGRkdbW1sbGxqGhoXfv3h06dKiVlZWHhwf+jeTp06cGBgajRo0yNDT08fHBv5HMnj27qalp+PDhw4cPLykpYbFYS5cu1dXVHTlypK2tbUFBwRfOOQCgp8DCvB25cwfV1fXGCzk7o8GDW1dwudw1a9ZwOJyCgoLa2tqFCxcihCgUyrVr1zQ1NRFCR48eXbVq1e3btydNmhQWFlZUVKStrd3Q0PDrr79mZmZyudzAwMDw8PBp06a1tLR4e3ufPHlyzpw5ERERd+7cwf9q490UDofDX1bp0qVL+CR2p0+fXr58OZ516uvrc3NzX79+TSKRPD09T58+vWjRog7fxLVr1+Li4rKyslRUVA4ePBgcHJyVlUWn0ydPnhwbG+vj44MQ4o9gwsve3t4+Pj4bN26sqqpq3ftpbGzE00Zzc/PVq1dfvHihr69/4cKF6dOnFxQUtJ5gSFJSMjExEZ9aYceOHevWrfv999/pdPrs2bOvX7/u7Oyck5NjY2ODX3iMjY19+PBhXl6elJRUQkJCSEhI6ykaEEIXL168cOHC69ev5eTkkpKSpk2blp+ff+HCBWVlZXyVj4aGBnFx8ZiYmHPnzp0+fXrYsGEIoWfPnvGnHamtrW1sbMzNza2pqdHX16fRaJmZmRQKxcHB4cqVK4GBgXp6epmZmdLS0lwuNyQk5PDhw99//31sbKyVlVV6ejreyPbt24uKit68eSMmJnb8+PH58+enpKR05vMEAOgKSITtYBj69Vf0mVncepiYWJtESCKR8Mfe5eXl4+Li0tPTvb29KRRKU1PTli1bysvLa2pq8C6XjIzM1KlTT5w4ER0dfebMmREjRujo6GRnZ79//15VVfXu3bs8Hk9XVzcpKWnOnDna2trR0dFz5sxxd3dvPWEpjsPhbN26taysrK6ujt+fQwjNnz8fn5rIzc3tC6PCbt68GRISoqKighBatGjR999///bt25cvX+rp6eFZECHEnwwoJyfH1dU1JiaGP13O53z77bf4HHL+/v5Lly7NysqytbVtdebEaDTaL7/8UlFRUVFRgS969ezZM01NTWdnZ4SQqakp//bb+fPnHRwc8H4eQujNmzeVlZV4wK13wCdzZzAYtbW1796909LSevz48Z49e3x8fAwNDb8cLX66SCSSkpKSqanppEmTxMTEEEKOjo74qVNWVk5MTExNTf348WNZWVmHLVy4cMHPz+/hw4cIIWlp6UePHrFYrNbpHwAgCJAI2yGR0MmTRL04mUzGZ5tDCKmpqa1Zs8bb2zszM3PcuHERERG2trbV1dWXL1/Gd1i0aNHYsWPXr19//PjxjRs3IoSqqqooFMrdu3fxHRQVFfFpS69du/bzzz//8MMPgYGB+/btaz3jWk5OzqhRo9auXTt+/PjGxsbTp0/zf8Qfoy8iIvKFNZNra2v5c9ZQqVRZWVk6nV5dXY0vWNFGZWUlg8H4XF5pPd+uvLw8v6ygoNDmXmBqampQUNC6devs7Ozev3//+PFjhFBjY2Pr+Zn45aqqqsLCQv5paXOdFiFEo9EaGhr4O4SGhoqJiX3zzTdxcXHnzp2LiYkxNja+fv26oqLi504CanW6REVFW5fxU3fw4MG4uLiwsDAXF5dbt251eNmTRqPl5ubyr0KvWrUKEiEAvQASYd/F4XDwK2+JiYne3t5LlixBCLVel3jo0KG6urpRUVHl5eXjxo1DCBkZGTU3N4eFhampqbVuSllZOTIyMjIy8saNGwsXLmydCO/duzd69Ohly5YhhO7cudOFOA0MDPh33UpKSvBrgywWKyoqislk8udvw7m5ufn6+vr4+Jw/f97JyQlfc6q2tlZeXp7H47VODy9fvsQLdXV1RUVFBgYGrdu5devW1KlT8UvHJ//54mJmZpadnY2nQw6Hk56ejk99bmpqqqmpuWXLls+9BVNTUwzD8JusrY0ZM2bMmDEsFmv06NEXL16cP3++iIgIh8Ppwlm6fv16eHg43g++evUqXikqKtq6NTMzs2HDhvHnQAcA9A5IhH0Lj8fbvn07QqiqquqXX37BV4Q3Nzc/cODAjRs3mpqadu3a1Xr/RYsWzZgxY/369fg1TDU1taVLl/r4+KxZs0ZOTu7ly5eDBw8OCAiYP3++j4+PsrJyQkICvhwEn7m5eUxMzNWrV9ls9u7du7sQ88KFC62srLZt22ZlZbVly5a5c+cqKSmNHDnSysoqICBg8eLFTCazqqpq3rx5+P4eHh7nzp2bPHnymTNnRo0a5erqunz58ilTpiQkJLQeUfn27dvw8HAvL6+DBw+OHz++zbJW5ubmUVFRnp6eNBptz549eKWJiYm3t/e3334bEBCQmJgoKiqKTxEeGRnp4uIiLS3t6OhYVVV1//79Ns//hYeH29vbDxo0yMXFpba29tatW8ePHz9z5kx1dbWNjU1VVdW7d++GDh2KELKystq3b5+Hh4erq2unzpK5ufmxY8dUVFSeP39+5coV/LegqqoqIyMTFRWloaExderUqKgo/GKyjY1NWVnZ8+fP2/y6AQCCQImKiiI6hh7T2Nh45MiR1atX45ssFgu/T9Nf8Hg8CoXCZDKZTKaGhkZkZCQ+7bWhoSGew+h0+qZNm2RkZPjrM6iqqu7du/fUqVP8Bf+++eYbFRWVxMTEZ8+eycnJjR8/Xl5evqys7P79+2lpaSYmJlu2bBEXF8cwTEVFBe9TDhky5PLly1VVVZs2bZKWlvb09EQItbS0ODo64tcnuVzukCFD2sz4gGGYrq6ugYGBrKxsQEBAampqWlqaj4/P+vXr8VnCp0yZUldXd/v27Xfv3o0cOVJPT4/D4WhqahobG2tra48cOTI5Odne3t7X1zczMzMtLW3atGnDhw8fOnSompra/fv39fT0TExMrly5MmzYsB07drRZSMXCwkJcXPzq1atMJjM6OlpWVhbPTBMmTJCRkamoqJg7d+5ff/3l5uZmZmamrKwcFBSUlpZ27969qqqqMWPGGBsbt25NTk5u+vTpGRkZ+HLBXl5eZmZmZDL5yZMn9+7d+/Dhw4YNG0aOHIkQ8vT0bGhowGdP1dDQkJKSwlfVYDKZ7u7uePeXxWLZ2Njg9yDZbLaenp6+vr6Li0t5efmtW7dUVVVXrVqlqqpqbW1NJpN9fHyKioo+fPjg6OhoaGg4ceLE1NTUpKSkurq6cePGdbguY7/7YPdfGIZxOBx8QRXQoZqamrNnz/bIzJeHDh0KCgrq8JaKoMEUa/3bunXr8vLy8EcOBhL8yYoDBw509sCcnBxdXV1xcfHk5ORJkybl5OS0uUo8AAjDB7uPgCnW/tPAmGINLo32V2w2e9CgQYMHD+bfcAIIoXv37u3atYvNZqurq587d27gZUEAQI+DRNhfiYiICG4CNsLhg2C7YMmSJfioIgAA+EowswwAAAChBokQAACAUINECAAAQKhBIgQAACDUIBECAAAgHpfL5a931ssgEYIegGEYUZ9gAEC/xmKxIqJiVA0t3yMl14nTx0wK+vDhQy/HAI9PgB7A4/FaWlrguWMAQGdNmbPwNkuv+buniEytQSgxL9npG9/sP++3nkBf0KBHCAAAgBjv379/lF3c/E04Iv/dK8OMPSqtpx09EdebYUAiBAAAQIwXL140649sU9li4Jby518d7i8gkAgBAAAQQ0xMjMJutwo6q0lSQryj3QUFEiEAAABiODk5ieYmIg6rdaVs5vnAcZ69GQYkQgAAAMSQk5OLWLpA8dhEVPoKYTzUUCl7dZ0NqXTixAm9GQYkQgAAAJ3W1NS0LDzS1TvgbS3bwNb59NnfuraoX9jiBbd+2jwqc4fIGj2z30N2TjBLun4RX9O018B6hKAHcLlceHyi18AHu9fAeoSfw2KxrEeOfmsR0uI0C5HIiFErmxA+y27w/m2butwmgesRQo8QAABA55w+e65Yc1SL8xxEIiOEkKR8/dSfz12/+/HjR6JD6wpIhAAAADrnZvLjJtOxn1SRSM3GXk+fPiUoom6BRAgAAKAH9OptvR4FiRAAAEDnjPNwlsq59UkVhonnJdrb2xMUUbdAIgQAANA5M6YFa5ckiz36BWE8hBBi0GXPLQj28VJVVSU6tK6ARAgAAKBzREREnib/sVD1g9ouB0qEkeHJSYfnj927JZrouLoIVp8AAADQaVJSUvu2bVo8d8b48ePzM1KJDqdbIBECAICw2LhxY0tLC0KIRqO1tLQMGTIEr4+OjhYTEyM0NCJBIgQAAGGhoKCAJ8KsrKz6+npLS0uiI+oTIBECAICwCAsLwwsiIiIfPnwIDw8nNp4+AgbLAAAAEGqQCAEAAAg1SIQAAACEGiRCAAAAQg0SIQAAAKEGiRAAAPq0R48eKf5DTExMTk4OL/v4+BAd2gABj08AAECf5uDgUFhYiJf9/PxWr17t4uKCEKJS4Q94z4DzCAAAfRqVSlVQUOCXZWRk+JugR8ClUQAAAEINEiEAAAChBokQAACAUINECAAAQKhBIgQAACDUYNQoAAAIxL59+86cOYMQwjCsqKhIR0cHrw8LCwsJCSEyMvApSIQAACAQ/v7++AN/DAZj7NixR44cwev5y+GCPgISIQAACMSQIUPwnNfY2Egmk4cNG0Z0RKBjkAgBAOATYWFhZWVlCKG6urqamhpdXV28ft++fRoaGoSGBgQCEiEAYCBgMpmvX7/Gy01NTeLi4hQKBSEkLi5ubm7eqabGjh3b0NCAEEpNTX327FlAQABeLyMj06Mhg74CEiEAgEgXLlzg3zzLzc01NjYmkUgIIX9//wULFnx9O5WVlfz9379/Lycnh89DpqWldenSpU6FNHbsWLzAZrMrKir4iRAMVJAIAQBEsrOzk5eXx8t+fn7Lli2TlJRECGlra3eqHS0trfT0dLw8bdq08ePHBwcH92yoYKCCRAgAIJK2tjY/51EoFA8PD7gCCXoZPFAPAABAqEEiBAAAINQgEQIAABBqkAgBAAAINUiEAAAAhBokQgAAAEJNIImQTqdv2bLF39/fy8urubmZX89iscLCwvT09GxtbfmPuGZkZAQFBRkbG5ubm69cubKxsbF9gwsWLPD6R1hYmCBiBgAAIJwE8hwhnU4vLS21t7cPDw/ncrn8+q1btz59+vT+/fu5ubn+/v6WlpaGhoZ5eXmurq7R0dEsFis0NHTFihXHjh1r02BaWtrixYvxKWulpaUFETMAAIBeVlJSkpaWhpcbGhr++OOPV69eIYQ0NTUdHR17LQyBJEI9Pb0ff/yxtLQ0PDy8df3Ro0dPnDihqampqanp6+t74sSJrVu3tp79YeXKlevWreuwTSMjI5i7HQAABpKSkpLz58/jZXl5+fv374uIiCCEnJ2d+30i7FBtbW1ZWRk/mdna2qakpLTZJzU1dejQoR0evnLlSnFxcUtLy8jISE1NTYGGCq/EJScAACAASURBVAAAoBc4Ozs7Ozvj5YaGBqImFeq9RFhVVYVaTd8uLy9fWVnZeoc7d+6cPHmSP1tga6tXrzYwMCCRSD/++KObm9vLly87vEDKYrFoNBo+0y5CCMOwZcuWff/99z38TkA7XC63paWFx+MRHYhQaGpqwqelHngwDOuRd8fhcJhMZocDDjqFyWSy2ezut9PU1IRhWPfbQQhxudzm5ubuN9XS0tIjb43BYPB4vB55awL6YIuLi1Op/5Hpei8R4vmpsbFRTEwMIVRXV6ekpMT/aWpqakhISEJCgr6+fvtjp0+fjhfs7Ox0dHRSUlK8vb3b7yYqKqqsrJyTk4NvNjY2qqmp4R1tIFBcLldERASfKxkIGoZhA/VOOYlEkpKS6v67o1Kp4uLi3W9HXFxcRESkR842iUTqkXYoFIqEhET3mxITE+uRtyYpKUkmk3vkrRH4we69RKioqKigoJCTkzNy5EiEUG5urp6eHv6jtLS0yZMnnz171tXV9cuNUCgUGRmZ1iNR2yCRSPweIZVKhSwIAAA9buXKlQkJCQghDodTUVHB78Ds2bPHz8+P0NC6QlCJkE6n19XVIYRqa2s5HI68vDyJRJoxY8auXbscHByKiorOnz9/584dhFBGRoa3t/eBAweGDRtGp9NJJBK+Jsu5c+c4HM706dNpNNqbN2/s7e0xDPv5559LS0tHjBghoLABAGBgq66uXr1h87Wbt9kcTklN055NkUOGDOlsI2vWrFmyZAlerq+vl5WVxcsqKio9GWtvEUgi5HA4w4cPRwjp6em5ublJSEhkZWUhhKKioqZPn66srEylUiMjI/GBM3fv3pWTk1u/fv369esRQtLS0i9evEAIvXz5ks1mT58+vbGxcc6cOYWFhSIiIqampgkJCRoaGoIIGwAABrb379+P+HbSx9FruWtiEIl0ITf5/mifO/FxVlZWnWpHRUWln+a8DgkkEVKp1MLCwvb18vLy165dY7PZra9YhoeHt3nKArd161a8oKurm5OTw+PxeDzef97zBAAA8DmLwzeU++zEjP6+CYWZeVUqaM5Z/l1Gyh/EBkYsAvJK1+7bkclkMhkmhAOgrzh69Ojbt28RQs3NzUVFRSYmJnh9aGgo//Y/6Gv+epmFeX06FEPdpKSCxuPxhPkPrPC+cwBAd8jIyCgoKCgoKDCZzFevXin8Ay7b9GkdPp9AobaeAkwIwUcWANAVQUFBeCElJSU/P7/DGxygrxmirlZekY/UjP6taqDJipKFfIA9JEIAhAiXy01OTsbLNTU1oqKi+JNbFArFw8OD0NBAbziwZYP33PnV02LRID2EEKotVzg7d/emSKLjIhgkQgCECIvF2r59O14uKCgQExPT0tJCCElISEAi7OPS09NXRW37M6tw+uJVM6cGrP9+hbi4eGcbcXRwuP7L3vkrFr2vqOYh8hBF6UPbN3p5eQoi4H4EEiEAQkRCQiIxMREvr169Wl1dfdWqVcSGBL5G3Jlzq3Ydq5m0F00w+8Bu3vv4xKWRo58/vCshIdHZphwdHF49Ttq5c2dpaen+/fsFEW2/A4NlAACgT2Oz2Wt+2FoTehlpmCGEkIhEs9vi90YTDvx0tMttUigUGNbEB4kQAAD6tKysLJ6WLRL9ZC5fpvXkhDspBEU00EAiBACAPo3D4WBUsba1FBEWq4WIcAYgSIQAANCnmZqakt49RbxPHvWj5CW7O9kTFdIAA4kQAAD6NGlp6dlTJ8vGL0Utfy/7R8p/qPpw35oVS4gNbMCAm6UAACBAxcXFx+JOM6lSO/bsmz192qBBg7rQyLaodXon4mL2eFU0tEiLUYaZm/xy+4qysnKPR9vLWCxWU1MTXq6pqeFwOHhZVFRUSkqq18KAHiEAAAjKj8djh4+bElOiyVkYH/FSwsJ17OWr17rQDolEWjB3dvHrdCddxQtH9iRe+Q1/ALS/u379uv4/jI2NdXV18fLy5ct7MwzoEQIAgEC8fft2475jVUsSEVUMIcTVd6ocNnnhd17uLiP564d3FoVCIXU4X2j/NGnSpEmTJuFlY2PjhIQE/uztvQl6hAAAIBDxlxJqHOai1gM+JeQahk7iz2kA+ghIhAAA0FZDQ8Oy8MiVUduuP3g6PnBGfn5+FxopqaDxZNXbVDbLaJR/rOyJGEGPgUQIABhQcnNzi4qK/vrrLxqN1rUWqqurhzqPOlKj93HlE8YPuTf154/wC7mblNzZdiyN9MQqc9pUytGyjQz0uxYYEBC4RwhAP5CdnX3t2t+DLF6/fm1iYkKhUBBCZmZmPj4+hIbWhzAYjIBZC54V0WuGjM3IZ5xy/Xb1wtnfLV/c2XYiN+8oGRHGHR7497a+Y9W8i6Er/N+++LNT7UwJmBy9x6NiqB9S0v67qvSVwvsHo0dv7WxIQKAgEQLQD7BYLDqdjpevX78uLi6uqKiIEOIPPQcIoZn/C7sr48qaPwshxEWIyY3cEhdirK/j6z2+U+3cSX7AXbj+kyo5NYaofGVlpYqKyte3o6CgcPXM8anzQmpVhtKlhijVF6o0vr988VdRUdFOxQMEDRIhAAL07t27mpoavNz6z6iuri6eyb6StbW1tbU1Xr548eLq1auNjIy+fIiwaW5uvv/nc9aqw/9WUURqfbdtPfBdZxMhl8tFlHYL1VLFWCxWZ6OyGz48Lz314cOH48aNu/3okY2NzUAa8zlgQCIEQICOHz9++/ZthFBLS8ubN2/Mzc3x+piYmDFjxhAaWh/C4/EOHzm+9+cTjaLy5k4e80Kmrlm5rLPdpo8fP5L5VyD5BumWlZV1Nh4bK8ui/IfI2O3fKhaDXFOsoaHR2aYQQlQq1c7Ojkql2traduFw0AsgEQIgQDExMTExMQihvLw8Pz+/9PR0oiPqi/xnzE+sU2xccBuJSpZwWNuT9ib6Bjz8I6FTjSgpKWF1FW1r6z/Ky8t1Np4dG9c89p1aKXkUaQ5FCKHGavn4xZGrl5HJMLpwYILfKwCASK9evXr4htbot+XvZYaoooxvwl9zlO/evdupdmRkZIyGDCLl3/+k8t6uhbOCOxuSoaFh0oVT9qkbZLcOF988XDd2wpFV0/43f05n2wH9BfQIAQBESn30uMZobJtKusn4O/cfe3p6dqqp+BM/efgEfMh1atT3QC2NSs9//cZ88II5s7oQlbm5+Z/3bpw+fTohIeHChQtdaAH0I5AIAQBEaz9+hETCMKyzzairq2c/fXD5ypW10RstTI0jDq0fPnx4d+KiUCj4YypgYINLowAAIo1wdlLM/6NNpULuzTHuI7rQGplMnjxpkp2Fkb/vuG5mQSA8IBECAIg0dOjQEXqK0tfWI3YzQghx2RJ3d5mSKkaPHk10aEBYQCIEABDs4qnjP3hpa//oSVpjMOSA6/dDKcnXL8LzdqDXwD1CADrw/v17Ho+HEGKz2VwuV1xcHK/X0dGBMfQ9jkKhrFiyaMWSRbKystlpWTIyMkRHBITLlxJhY2NjVVUVlUpVUlKSkJDotZgAIJyPjw+DwUAI1dXVcblc/iwwr169kpSUJDQ0AEAP6yARpqamnjx5Mikp6e3bt3gNmUweOnToN998M2fOHGNj496NEAACvHr1Ci/s2rXr48ePO3fuJDYeAIDgfJIIExMTv//++8zMTD09PWdn59mzZysqKnI4nOrq6hcvXsTFxe3cuXPcuHE7duwwMzMjKmIAAACgB/2bCC9dujRjxoz58+fHxcVZWVm135XH4yUnJx87dszGxiY/P19bu920fgAAAEB/828itLS0LCwsVFVV/dyuZDJ59OjRo0ePfv36taysbK+EBwAAAAjWv4nQ0NDwK4/hz6APAAAA9HcwEBwAAIBQ++zjExkZGfHx8W/fvq2trW1dn5iYKPioAAAAgF7ScSKMj48PDg7W0NBgs9kSEhJSUlK5ubkSEhIODg69HB8AAAAgUB1fGl23bt3EiRPfvn07duzY4ODgV69e5eTk6Ovrjxo1qpfjAwD0cTwej8vlEh0FAF3XQSJkMBiFhYWrVq2iUqkIIRaLhRAyMDA4evToDz/80NDQ0NsxAgD6pPfv34/yCfAPXflnUa2m+fAfj/3ShbWTACBcB5dGm5qaMAxTUFBACCkpKVVXV+P1FhYWTCazoKDA1ta2V2MEAPQ9Hz9+dP52UsWEfZifE0KotKVxTXx4UWn59uhIokMDn5WSksLhcBBC+fn51dXVd+/exevd3d3xno9w6uCdKysrS0lJFRUVGRsbGxoabtmypampSUpKKiUlBSGkpKTU2zECAPqeLXsOVrqtwvSd/t4Wk24IOBC722nD9yukpKQIDQ181t69e/FJdKurq9ls9vbt2/F6Z2dnSISfIJFIo0ePvnz58jfffBMUFBQREWFpaWlsbJySkuLi4qKlpdX7UQIA+prUpxlc30WfVJEpmK59Tk4OrIjbZyUkJBAdQl/U8VeA2NhY/FuDrKxsUlLS3r173717t3DhwvXr1w+cRcIaGtCffxIdxABB4vEobDYSEyM6kJ6nXVAgS6ejf64gdZlkSYkzg9H9dhBCzgyG5OPHqLi4m+0YFhUpNTZ2OaQRDXWK2UlIdlDrSrmP7xSfP0efPnP19Ty4XEpyMur2+h4WFRVqWVndP9tqWVmWHz92vx1Kc7MHl9sjv/1hdLpCRgYaiKOTnBkMqSdPUGnpf++qpIRsbHrwpUkD6eZ2RUWFjY1NeXk5vtnQ0PClhc0ePkQ//NBLkQ10GIbxeDwKhUJ0ID2vqKiIxWJ9/bxLn8NgMDIzM52dnbsf0uPHj62trbu/GlRubq6YmJiurm7XDn9fXPy2HuMptLpExOOIljx3cbLv8tfl5ORkV1fX7n+QsrKylJWV1dTUutlORUUFjUaztLTsZjtcLvfBgwceHh7dbAchlJGRoaenh4/hGGAeP35sZWX1VdfV5eTQhQs9+drYAFJeXq6mpsbfrK+vJzAYocLhcPAxVgPPzp07V69e3f12cnNzjY2Nu98OhmEGBgZ5eXndaeH3Cxc1zYZJmjhLGzsa2Djdu5fUhUYYDIbp8BES03ajQ3R0jIXWP1W0dPn1t/PdCUxGRqZH/tsGBwf/+uuv3W/n119/DQwM7H47DQ0N0tLS3W8HwzAPD4+kpK78vvo+IyOjnJwcQl76k0ujERER/DGin3PkyJGezMMAgN51+uzvyw6cqw29gaQUEUJvassCwmZcPkx1dXHpVDsSEhJ/pd6L2bUv7ievysrKkSOcd506YG1tLZioARCgTxLh2bNny8vLv9wzhUQIQL+2fuuu2tBbSELu7215jZqgoys3hKUn3+psU2JiYj+sCx81wiE6Ovru5XM9HCj4x9OnT4OCgvByeXl5SEiIuLg4Qsje3v7cOTjtPeCTRKioqFhaWurk5DRz5kw/Pz+xgTj2AQBhxmKxGDzKv1kQN0ivrJJGUETgvw0dOpQ/yXN9fb2MjAx+FxZPh6D7PkmE6enpSUlJp06dmj17toiIyJQpU6ZPnz5y5EiiggMA9CwqlYq4rLa1GEbi8YgIB3wVcXFxPT09oqMYyD6ZYo1MJnt6ep46daqsrGzHjh2vX792cXExNTWNiooq7vZYbQAA4chkss5gdVT6qnUl6fUfzg52RIUEAOE6nnRbTk4uNDQ0NTU1Kytr7NixMTExCxYs6OXIAACCcPLQLrVzcyiZVxCLgZrrRdNODkmMPrA1iui4ACDMlxbmzcrKio2N/e233xBC3X+SBgDQF5iamr5MTZxDeaq401Ht0KilaqWv01LU1dWJjgsAwnQws0xtbW18fPypU6cePXpkamoaFhY2c+bM7j+aCgDoEVi3J8EYNGjQ0f07ZUUwdXX1VatW9UhUoL3k5OSnT58ihFgsFovF4k/s6eHhYW9vT2ho4BOf9Ahv374dEBCgpqa2du1aGxubZ8+eZWdnh4eHQxYEgHAcDmf73gODTW3fN4s4j/OfMntBVVUV0UGBL2lubqbT6XQ6vbGx0c3Njf4PJpNJdGjgE5/0CBcsWFBRUeHn5+fj4yMmJvbu3bt37961OSAgIKAXwwMA/G3mouUJlbJNyx8gqlg1QhczE9JHj89KS5GQkCA6NNCxcePGjRs3jugowH9re2m0paUlPj4+Pj7+cwd0/7IMAKCziouLE9Ozmxbf4ddwrf3Kafmxp878b8F8AgMDYAD4JBGeOXMG+uwA9EHPnz9nGLi1qWw29rz7+BdIhD2upKSEzWYjhCorK5uamt6+fYvXa2pqioiIEBoaEIhPEiE8Ow9A30SlUslcdttadouoiPAupio4c+bMwZMfk8lkMBheXl54fWJiIjzYPiDB/yIABIvL5R46cvzH2F/f0hrdfQK2RKx0dnL678M+5ejoKP5ddMPY9Yj87xJFsq8uT5rZAyv7DAzl5eXLly/Hy3/++WdhYeGVK1cQQurq6vv37+9UU/z5zHg8HpPJ7P6KV6CP+zcR3r59u7CwcN68eaKiol84oLq6evfu3QsWLNDW1hZ8eAD0bzwez2Ws3ytZm8aQ80hK8X5Zts/S8HUz/VYuXfTfB7eirKy8eOaU/bFB9Am7kJIWammSStlvwsjxn7xLQJH3O9LS0vyhfE5OTgoKCvj6AV9alBQAhFDrRKiqqrpw4cLNmzcHBwdPnjzZ1ta29aTb9fX1jx8//u233y5evGhhYfHdd98RES0A/cz5Cxdfixk1frvh720Ns5r5l7btHTl3RrCcnNwXD21r45rVzsOs125Z8uJVlubgwQtnT1u55BqZ/KU5MYSKjIwMjGkHXfNvIrS2ts7Jyfnxxx8PHz68e/duUVFRHR0dRUVFDodDo9FKSkp4PJ6Dg8Mvv/wSEBDQ5RWoARAqF24l1Q+d9kkVRYRl4pWWljZmzJjOtubl5enl5WloaHjj/CkjI6MeixIA4fbJPUJxcfGVK1eGhYWlpqYmJSW9fPmSRqOJi4sbGhra2dl5enrCRGugL+NwOA0NDe3rqVQqUdfHWlgsJNL2OT8uRaylpYWQeAAA7XUwWIZMJru6urq6una5US6Xm5eX99dff0lISEyaNIlfj2HYpUuX0tPTjY2NQ0JCqFQqQqimpubmzZt5eXlKSkoBAQGDBw9u3yCbzT59+vSbN2/s7OwmTpzY5cDAwJaWlubr64uXGxsbJSQkKBQKQmjkyJFXr14lJCR3B5s7f95v0RzaulLi7UMbG3jmAYC+QiCjRg8dOrRjxw4FBQVxcfHWiTAiIuLq1atz5849fvz4zZs38cf28dWWbWxsXr58uWHDhkePHrXvd06ZMqWqqmrChAnr1q17/vz5pk2bBBE26O9GjhxZU1ODl93d3aOjo93c2j5718tC58zaf9S9RM0UM/VECCEuW/L21tE2JpqamsQG1ncUFRUVFBTgZS6Xm5ycjI/S1NbWNjQ0JDQ0ICwEkggXL168fPnyM2fO7Nu3j19ZW1t78ODB9PR0ExOTefPmaWho5OTkmJqaxsfHS0tL4/s0NzcfP368zVjnrKysxMTEsrIyGRmZMWPGODk5fffdd7KysoKIHICeJS0t/eTO1flh4X9ei6hjk5XFeAtnBkes3vDfRwqNZ8+eHTlyBC8rKCgcOHAAH4Lg7+8PiRD0DoEkQvyaZxtPnz4dNGiQiYkJQkhWVtbR0TElJcXU1JSfBRFCPB6v/cSJycnJzs7O+D0ec3NzOTm5Z8+ejR49WhCRA9DjNDQ0bsSfzs7O9vPzK8gtIDqcPsff39/f35/oKIBQ670H6svLy1VUVPibqqqqZWVlrXd4+PDh7du3MzMz//PA8vLyDl+CzWbX19fPmzePvzlhwoRvv/22Z94A+Dwul9vS0tKnhvLzeDwWi9WdKQMZDMbOAz+ejr/EYrHYFLG1K5Z29oGH1jgcDplM7pEpDDEMa2lp6X5THA6Hw+F0vx0Wi4U/eN7Ndvog/H31qQ/2ANZTH+w2RERE8LECX9B7iZBKpfJ4PP4ml8tt3XHMysoKDAyMjY3V0dHp1IGtUSgUERGRYcOG4ZtMJlNDQ+M/TwHoERQKpU+dahKJRCaTuxxSRUWFy7cTK4bPbpl/DZEoP76+dd5tTFLC712eYQuPpKdOUY+cbTKZTCKR+k47fRD+vgbkW+ubBHG2v+Zhv95LhOrq6q27gGVlZfzLm3l5eWPHjt2zZ0+HI0LV1dUfPnyIlzEMKy8v19DQ6PAlyGSyhITEokV/z9nR0NAAk0r0DjKZzOPx+tR8xCQSiUqldjmklZGbij03YBZ/X05gO4SUqRjPX7H2wc1LXWuQSqWSSKQeOUXdfGt8+BeF7rfTg2+tr+HxeFwud0C+tT6opz7YXfDZLn9ZWdnmzZuDgoKCgoLwmmvXrt2/f7/Lr+Tk5MRkMp88eYI3np6ePnbsWITQmzdvPD09N27cyH8hXGZmJo1GQwh9++23f/75J3459NGjR1wu18HBocthAPA1Hj9Nx8zHflKla5dX+A6WIQNg4Om4R/j8+XMvLy8Oh6OlpUWn0/HKrKysM2fOvH79+j8bff78eXh4eEVFRXFxsZeXl729fUxMjISERHR0tL+//4QJExITExcuXIiPIJ87dy6DweAvgujo6Ig/HTFlypSIiIiZM2fq6OjMmzfP3d3d09Pz8uXLmzZtaj33GwCCwCORUPsrKlTRL1yZBwD0Ux3/l16wYIGFhcWVK1devHgREhKCV3p7e0dERFRWVrYeutIhHR2d8PBw/qaysjJeWLJkyciRIzMyMoKDg0eMGIFX7ty5s76+vv3OZ86c0dLSwssHDhx4+PBhfn7+ggULhg795NlkAARhkLxcRU0pUhzyb1VznQQZgywIwMDTwf9qOp2enp7+8OFDBQWF1rcZ8WEsZWVl/5kIFRUVPT09O/yRtbW1tbV16xp7e/sO97Szs2u96eLi4uLi8uXXBaCn7Ni4JnjtAvr0OCQzCCGEmuvkz4b+sGYl0XEBAHpeB4kQnwWx/Uhx/BopjCQGwmDsmG9i2eyla3zpZFmMRJFrqdq2IXx68FSi4wJgQCksLLx37x5erquri4+PV1NTQwjp6+v35sPiHSRCFRWVQYMG3bp1y9LSsnWPMD4+XlJS0tjYuNeCA4BAft7j/bzHb9iw4ePHj0eOJBIdDgADUE1NTUZGBl42MDB4//79hw8fUK/3uDqedHv58uUbN27EMExDQwPDsPz8/Pj4+M2bNy9btgwGqgChIisr29zcTHQUAAxMdnZ2/LtgBD7w1vGd/7Vr19JotIiICPxJdrwXGBISArNdAwAAGGA6ToRkMnnfvn3Lly+/d+8ejUaTkZFxd3e3sLDo5eAAAAAAQfvSWHBdXV3+vJ0AAADAgNRxIiwqKuJyue3r1dTU8KXCAAAAgIGh40To4ODw8ePHDn+kqqo6ceLEbdu2dWcmfgAAAKCP6DgR7tixY9WqVaampr6+vsrKyuXl5b///nt9ff2qVauys7N/+eWX/Px8/sMfAAAAQP/VcSK8cuVKUFDQgQMH+DXh4eETJkwoLCw8fPjwmDFj/Pz8MjMz28wRAwAQHlFRUdnZ2QghGo2WnZ0dGBiI12/cuNHc3JzQ0ADonA4SYUNDQ0JCwosXL1pXksnk+fPnz58/f8+ePb6+voqKiq9fv4ZECIDQcnd3xxMem83+8OEDfyXR/5yCEYC+poNE2NTUxOPxamtr29TX1tbyZ8dWUFCAudYAEGbu7u5EhwBAz+ggmamqqhoYGKxatQpfAhBXUFCwceNGfMmIpqamkpKSIUOGtD8WAAAA6F866BGSSKTjx4+PHz9eV1d32LBhysrKFRUV6enpysrK+F3D1NRUGxubNqtDAAAAAP1Rx5c33dzcsrKyli1bJiMjU1RUpKKismHDhqysLFNTU4TQmDFj0tLSxMXFezdUAAAAoOd9dmYZHR2dHTt2tKlks9kiIiICDgkAAADoPV874CU3N3fNmjX8JeMBAACAgeFLc40ihOrq6n777bfY2Ng///yTQqF4eHj0TlhAqDQ3NzOZzPb1EhIScAW+Z2EYxh8QzmQym5ub8QW3SSSSvLw8oaEBQJiOEyGPx0tJSYmNjb148WJzc7OIiMiePXumTp2qrq7ey/EBYbB58+affvoJIYRhWH19PX/2vpiYmEWLFhEa2kDT1NSkr6+Pl5lMJplM3rNnD0JIWlq6uLiY0NAAIEzbRFhUVBQXF3fy5Ml3796pqaktWrSISqX++OOPK1asICQ+IAxiYmJiYmIQQnV1ddra2jU1NURHNGBJS0vD6QWgjU8SYWBg4MWLF0VFRX19fQ8ePDhmzBgqlXrixAmiggMAAAAE7ZNE+PTpU1FR0a1bt86dO1dGRoaomAAAAIBe88mo0S1bttjZ2a1cuVJNTW369Ol37tzpcFVCAEAvS0pK8vpHWVnZ7Nmz8fLWrVuJDg2Afu+THmFwcHBwcHBBQUFsbOypU6fOnDkzePBg/q11ANp4/Pgxg8FACDU3N9fX16uqquL1zs7OxC7gjGFYS0tLhyNR+ykzM7Pw8HC8nJ+fb2hoSCKREEJqamqExgXAQNDBqFFDQ8MtW7Zs3rw5KSnp1KlTFy5caG5uNjc3DwgImDFjhp6eXu9HCfqmI0eOlJWVIYQqKiqqqqosLCzw+pMnTxKVCLlc7qbtu3+KPV1LlQ9aEaWnsvfUj3vMzMwICaYHqamp8XOeg4MD3LkAoAd99jlCMpns6enp6em5d+/eM2fOxMbGRkdHb9q0CS6WAr6TJ0/ihVOnTt26devcuXPExoMQWrgi/FyRSNPqPxFFhIVQRumrUZOmPU+6oaGhQXRoAIA+6r9nllFSUlq+fHlmZmZGRsb//ve/XogJgK6pra1NuPugyTsaUf6ZCHCIJc0jfNu+w4TGBQDo0/5jZpnWbG1tbW1tBRcKAN2Uk5PD07FDJFLrSp7BiCc3TxMVEgCg74PFdcHAISEhQWY1tq1taZKUIHLk/clYVgAAGchJREFUDgCgj4NECAYOCwsLaskLxKhtXSmVcXaK7xiiQgIA9H2duDQKQB9HpVIP7di0IMKnanwM0nNETdUyj44a1zybP3s9USGdO3fuxYsXCKGamhoajbZmzRq8PigoyMrKiqioAACtQY8QDCiT/HweXT7lXxkvsdHM5PyMLd9op927SeAimlJSUgoKCgoKCtra2t7e3gr/EBUVJSokAEAb0CMEA42RkdH5uCPu7u7R0dFubm7EBuPr6+vr60tsDACAL4MeIQAAAKEGiRAAAIBQg0QIAABAqEEiBAAAINQgEQIAABBqkAgBAAAINXh8AoAO7N27l8ViIYQePHjQ2Ni4fft2vH7FihXwCCAAAwwkQgA6QKfT8USopqbGYrHodDrREQEABAUSIQAd+OGHH4gOAQDQS+AeIQAAAKEGiRAAAIBQg0ujQqe2tjY+Ph4vv3//XllZWVpaGiEkLy8fGBhIaGgAAEAASIRCp7m5OSMjAy/fu3dPV1dXT08PIaSurk5oXAAAQAxIhEJHXV39yJEjeHnq1KkTJ06cMmUKsSEBAACB4B4hAAAAoQaJEAAAgFCDRAgAAECoQSIEAAAg1GCwTL9RW1tbU1ODlxsbG/FnHhBC8vLyioqKxMUFAAD9GyTCfiM+Pp4/9XNRUZGWlhaJREIIhYaGhoeHExpaz2hqarp16xaLw719+/aoUaNERESIjggAIBQgEfYboaGhoaGheFlKSiorK0tSUpLYkHrQH4l35y4Przce2xy4b8qh2/Kr1186dcTWxobouAAAAx8kQkC8ioqKGUu+oy26haSVEEJ1CNVVvfcNDix4/lhCQoLo6AAAAxwMlgHEOxt/ge44H8+Cf1PWqTcec+/ePeKCAgAIC+gRClxDQwOHw8HLLBaLv6yrjIwMlQrnHyGE8t+XcpTd21Q2KBi8KyomIhwAgHCBP8QCN3Xq1CdPniCEOBwOk8nkj/a8ceOGk5MToaH1FXpD1CmvirmfVkrXFWkOHkZMQAAAYQKXRgXuxo0bNTU1NTU1ly9fdnV1rfkHZEG+oIBJimnHUEvjv1V1FdLZ1728vIgLCgAgLKBHCHpAU1NTc3Nzlw/X1NQ8uClieeToGttpbCU98YosxazL544dlJKS6sEgAQCgQ9AjBN3y6PFjYzuXNYd+vV1Qq25kdejnY11rZ4r/xNePEveOlBQ7Pf/EZMP89IeuLiN7NlQAAOgQ9AhB12VmZk6Yv6pqxhmkpIUQqmhpWhe/opHRvGblsi60pqSkFDIteF3E2qCpU3s6UgAA+CzoEYKuC9+0s2riHjwLIoSQmFR94KF9Px/n8XiExgUAAJ3Q24mQx+Pl5eVVVlb28usCQcjNy0fanw7spIoiJZ3y8nKCIgIAgE4TSCIsLi52c3OTlZUlkUiNjf8OBSwtLbW0tJw0aZKFhcWSJUswDEMI0el0T09PBQUFEon07t27Dhu0srKSlZVVVFRUVFT89ttvBREz6AIRESriMNtUYsyGgTT3GwBgwBPIPUJJScmVK1eqqqq2eUJgw4YNzs7Ox44dq66utrKymjx5soeHh4iIyKJFiywsLExMTL7Q5tWrV93d3QURLegyn7FeP6b/znKa9W8V7Z2iCE9BQYGQeOh0ekZGRusym81GCCkoKAwbBo8kAgA6JpBEqKys7OfnV1pa2rqSx+P9/vvvDx48QAgpKSkFBAScO3fOw8NDWlp68uTJLS0tX26TxWLV19fLysoKImDQNdFrV9/y+Lakmc4YHozEpMi5ycp3fjh99hei4ikuLuYv0NHc3Hz58uVbt24hhKysrCARAgA+p/dGjdJoNAaDYWBggG/q6+tfv3796w+fOnUqj8eTl5c/ePCgj4/P53bjcDj8PgGDwTAxMRk0aFB3wgZfICsr+/JJ8q79h4+dnFxTU+M/wTfm/h9qampExWNlZZWYmEjUqwMA+qneS4T4zUJxcXF8U1JSsr6+/iuPvXTpkr6+PkLo5MmTQUFBOTk5mpqa7XdjsVh1dXXz58/HN7lc7ty5c+fMmdOFaNlsNj9P19TUiIiIyMjIIIRERES8vb270CBCqLm5mcvltr5p2h2NjY3dH5yJz/rWzZCWLZqvLCeVmJi4f2s0+ucX3TVNTU0YhvXUKRqompqa8KUogaDxeDwmkwmjoHuHgD7Y4uLi/zmrc+8lQlVVVYQQnU7Heww1NTVf33XAsyBCaObMmTt37nz8+PGUKVPa7yYqKqqkpPT8+XN8s6GhAc9eXdDU1HTt2jW8/OLFC0lJSUNDQ4SQlJTU1K4+5SYhIUGhUPhzjXaTtLR098ekUKlUcXHx7oeEf9S63w6XyyWRSD11igYqDMPgFPUOHo9HpVJh8FfvIPCD3XuJUFpa2tDQ8MmTJxMnTkQIPXnyxNbWtrONsNnsmpqaLqe3ryclJRUfH4+XV6xYoa2tHRYWJugXBQAA0PsEkgh5PN7x48fpdDpCKDY2VkpKCr8+uWzZsjVr1igoKOTm5t67d+/gwYP4/nFxcQwGAyF07tw5ZWXlOXPmUKnUyMhINpu9ffv2t2/fnjp1ytnZmUQi/fTTT7KysjB8FAAAQE8RSCLEMAwfsRIaGpqVlSUmJobXL168GMOwzZs3Kygo3LlzR0NDA69/8eIFg8EIDQ0tKioqKiqaNWsWQsjY2JjL5SKEZGVl6XT67t27SSSSra3t0aNH4UoFAACAniKQREihUI4cOdK+nkQiLV26dOnSpW3q9+7d237n6dOn4wVlZeX9+/f3eJAAAAAAgrlGAQAACDlIhAAAAIQaJEIAAABCDRIhAAAAoQaJEAAAgFCDRCikMAy7efPmy7zC67duZ2VlER0OAAAQpvdmlgF9x8ePHz39ppQoWNY5h+cy6/+YtWqCs+XR/TthBksAgBCCRCiMAmcvyh6xhmc6GiGEIVTlND0+fqnz6V9nzwghOjQAAOhtcGlU6NTX1+eWfMSz4L+V30Qcjv2VqJAAAIBAkAiFTlVVFUlBo22tvDqNVkVEOAAAQDBIhEJHVVUVqy5uW0t7x5/6FQAAhAokQqEjJSVlZ2Eskh7/bxWPK39j/feL5xIXFAAAEAYGywijs8cO+QbNep11pVrHTaSlTi772uLpARP9fImOCwAACACJUBjJysqm3Lj0119/hS5YNMLZMfzAZXV1daKDAgAAYsClUeFlY2Ojr6fj5OQEWRAAIMwgEfYz127cdPD0aRaRGTF2wrn480SHAwAA/R5cGu1PFq9aeza9pHbcXjRVO7O2fNHxmCu37v4e28EayAAAAL4S9Aj7jYKCgvjkZ7UhJ5CSNkIIyavXBR66l1/99OlTokMDAIB+DBJhv5GUnEI3n4A+nQ602tL/+u17RIUEAAADACTCfoPZ0sKjiretFRFvYjKJCAcAAAYISIT9ht0wW4WiB20qZd+muDkOJySeHhQTE6Ovr6////buPaipe8ED+C8kgQDK2xrKI8Rcr4wDC9pQmagUFdBBKaJVbGeLXrVa0bbWVgc7tt5tq3Ltbu2KndouBaxTR+v4KsYKDaBYrLI+KOiKFUMEY0yIFcsjhLzuH6fNUq94vULyk3O+n7/O+ck5fsUz+ea8ZbJx48Z1dnbKfldUVEQ7GgCwHy6WGTIUCsUY/ubzp77onbSU8DwIIYL/3Rv5S93M9P+iHW2gXn311RdffJGZNpvNXl5ezHRwcDC9UADAFSjCoeT7w/vWbfxw/9YEo0UQ4mmfMTVp+3eH+Xw+7VwDFRwcjM4DAFpwaNRNrl69euDwt9dabx86dMhmsz3eSnx9fT/9zy2Ga/UiQ2PzxZqvdm4PCAgY3JwAAFyDInSHte++Pyl7+c6u2Ja0D3OKT0fLJ2k0moGskMfj4W3yAACDAodGXe6o8ljhif9rX1nO3PnQGTezSXP++Zf+Un+6inY0AADAHqHrFRTvaZ+27g/3/0U9o7f5tLa20gsFAAC/QRG63C2djgRF3DdoD4rU6XRU8gAAQF84NOpyUonkkqGJRP3hbj+PtqaIiPvbcchpbW21WCyEkLa2ts7OTrVazYxHREQIhUKq0QAAHhWK0OXW5i6uWf0fvyw9QASezIjH5bI/BXuz4OVHCxcuvHHjBiHEZDKZTKbU1FRmvLKyUiKRUI0GAPCoUIQuN3nSpI1LX9j038/dG5vR6x0cpD07imc88s1u2rkGQWVlJTNhs9nMZrOPjw/dPAAAjwHnCN3h9VeXXj51PC/WMebCZ99tfb228ruQkBDaoQAAgBAUoduEhIRMnjw57OnQhIQE2lkAAOD/4dDow2i12jXvflimqhAKhI03dFveywsMDKQdCgAABhP2CPtVV1c3flrG/oCMext+Mq47W2T6t1jF1Js3b9LOBQAAgwlF2K+/vL7W8O+7HWNTiYeACLws8vm30je/vv6vtHMBAMBgQhE+mM1m07bdJeI/9x10RE87e/4irUgAAOAKKMIHs1qthP8Pt4TzeA4HjTQAAOAyKMIH8/Ly8vWwka5f/jCqa5REhFFKBAAALoEi7NffNq4P3L2I/Kr/bb5NHbJv+fZN79LMBAAAgw23T/Rr/twsb5HXmxvm6zrMHsQe9VTg/3y5DXcBAgCwDIrwYTJmpmfMTF+1alVUVNTbb79NOw4AAAw+HBr954RCoUCAbwwAAOyEIgQAAE5DEQIAAKehCAEAgNNQhAAAwGm4BoRztFrtm2++yUyfOXPmxo0bBw4cIISEhYVt27aNajQAAApQhJwzfPjwefPmMdMKhSIoKMjb25sQ4ufnRzUXAAAdKELO8fPzcxYhAADgHCEAAHAaihAAADgNRQgAAJyGIgQAAE5DEQIAAKehCAEAgNNQhAAAwGkoQgAA4DQUIQAAcBqeLDNklJWVHTx4kJnu7e197bXXmNcFT58+fc6cOVSjAQAMYSjCIUMsFj/zzDPMtJ+f3+jRo5npsLAweqEAAIY8FOGQERcXFxcXRzsFAADb4BwhAABwGooQAAA4zSVF2NzcvG3btoULF+bl5fUdv3PnzrJly+Ry+YIFCzQaDTN469atgoKCxYsX5+bm9rdCtVqdnZ0tl8tXrFhx9+5dV2QGAABuckkRnjt3rr6+3mKxqFSqvuOLFi3q6uoqKSmJjIycNWuWw+EghFy6dKm2tpbH4ymVygeuzW63p6enjxo1qqSkpL29ffHixa7IDAAA3OSSIpw3b15xcXF6enrfwebm5vLy8oKCgpiYmPz8/Dt37lRVVRFC0tLSdu/ePX/+/P7WplKpOjo6Nm/eHBMTU1BQcOzYsdbWVlfEBgAADnLfOcJLly7JZLKgoCBCiIeHh1wu/+mnnx5lwfr6+oSEBB6PRwgJCQmRSCQNDQ2uzQoAAJzhvtsnDAZDQECAczYwMNBgMDzGgkFBQf0t2NvbazQapVIpM+twOFasWLFy5coBpCaEEIvFYjabOzs7B7gek8lks9kGvp4nkM1mM5vNdruddhBO6OrqYr4XgqvZ7faenh5s2O7hog1bJBIxzx55CPcVob+/f3d3t3O2s7PT39//ERd0XllDCOno6OhvQU9Pz8DAwIqKCma2q6tLKpUOGzbs8UMTQggRCoVeXl4DX4+3tzefzx/4ep5ANptNKBT6+PjQDsIJDoeDlVvRE8hutwsEAmzY7kFxw3bfoVGJRKLRaCwWCzPb1NQUFRX1iAteu3aNme7t7W1paXnIgnw+f9TvoqKi8HkBAAAP574ilMvlYrH4q6++IoScOHGipaVl1qxZD/n5HTt2MCcRn3/++evXr586dYoQwlxxOm7cOPdkHhRJSUk8Ho/H46WkpFRUVPB+x/yLAACALpccGlWpVKmpqcw0j8dLT09XKpU8Hq+oqCg7Ozs/P7+9vb2oqIjZXaurq3MWG4/He/bZZ8+ePUsI+eKLLwICAuLi4vz8/AoLC7OysgIDA81m8/79+12R2XWqq6tpRwAAgH65pAhTUlKYewTvo1AoNBqNTqcbOXKkUChkBuPj4x/4w/X19c7p+fPnz54922AwhIaG8vl8V2QGAABucvdDt/l8fnh4+GMs6Onp+XgLAgAAPASeNQoAAJyGIgQAAE5DEQIAAKfhxbwPZjabP/nkE2b63LlzTU1NZrOZEOLl5bV69Wqq0QAAYDChCB/Mbrc73/cUGRkpEomYWZFIRDUXAAAMMhThg3l7e+fn59NOAQAALodzhAAAwGkoQgAA4DQ2F6FKpXI+4xtc6vbt2xcuXKCdgisqKip6e3tpp+AEvV5/7tw52im4orKyktaGzeYifOutt1paWmin4ITKyspPP/2UdgquWLt2bd8Xk4HrnDhxAhu226xbt06tVlP5q9lchOA2D3xaLMBQhw2bI1CEAADAaShCAADgNB6b9v1v3rw5evToSZMmMbM1NTXjx4/39vamm4oLdDrdnTt3YmJiaAfhhNOnT8fHx/v4+NAOwn63b99ua2uLjY2lHYQTfvzxx7i4uEHfsLOysnJzcx/+M6wqQkJIcXFxREQEM63RaCQSCY/HoxuJC7q7uzs6OkaOHEk7CCdgw3Ybk8l07949sVhMOwgnuGjDlkqlMpns4T/DtiIEAAD4l+AcIQAAcBqKEAAAOA1FCAAAnIYiBAAATmPta5iuXLly5syZUaNGJSUl4fo61zEajXV1dc7Z+Pj4kJAQinlY6datW1evXh0zZszTTz/tHDSbzWVlZb/++mtKSgouaxwsVqv1ypUrer1+2rRpzs+Nuro6o9HITItEIucNWjAQJpPphx9+0Ol0Mpls4sSJff+ourr6+vXrEyZMGDt2rHvCsPOq0T179qxevXrOnDnV1dUKhaKwsJB2ItY6duzYSy+9lJCQwMx+8MEHiYmJdCOxjEKhqK+vt9vtBQUFS5YsYQZ7enqSkpK8vb2lUunRo0dVKlV8fDzdnCxQV1c3ceJEHx8fo9FosVgEgt/2EzIyMpqbm0NDQwkhI0aM2LNnD9WYLOHv7z9u3DipVHrq1Kno6OgjR47w+XxCyPLly0+ePJmcnHzo0KGPPvooJyfHHWkcrGO1WiUSSWlpqcPhMBqNfn5+ly9fph2KtZRKZWJiIu0UbKZWq61WK/N9zjlYUlIyfvx4q9XqcDjee++9OXPm0AvIHh0dHTqd7ueffyaEWCwW5/isWbN27dpFMRgrNTU1MRPt7e2BgYHl5eUOh6OxsdHX19dgMDgcjuPHj4eHh/f9j3AdFp4jbGhoMBqNM2bMIIQEBwcnJycrlUraodisu7u7vLy8trYW7wZyBalUynxT7uvo0aNZWVnM+AsvvKBUKh1sPLTjZsOGDevvIPP169ePHz/e3Nzs5kgs5rzJ3d/f39fXl/n0UCqVSUlJI0aMIISkpKR0dHT0PfPiOiwsQq1WKxaLnYc1wsLCtFot3UjsZrVad+zYkZOTExMTc+3aNdpxOEGr1YaFhTHTYWFhZrPZeRILBp1IJKqurt6+fXt8fPzKlStpx2GboqIiT0/PKVOmEEK0Wm14eDgzzufzxWKxez69WXixjM1m63t1DJ/Pt1qtFPOw2/Tp09PT0wkhdrt9yZIla9asKS0tpR2K/Ww2m4fHb99imf1CbOSus3fvXuaX3NraGh8fn5GRwRxwgoGrqqrKy8srLS1lHjF636e3QCBwz4bNwj3C0NDQtrY255EivV7PnOUGV3AetfPw8FiwYIF7jmNAaGiowWBgpvV6vUAgeOqpp+hGYjHnRh4RETFx4sQLFy7QzcMaNTU1CxYs+OabbyZMmMCM9N2wCSF6vb7vldKuw8IijI2NFQgEZ8+eJYSYzeaTJ08yO93gaufPn3c+8RxcKjk5uaysjJkuLy+fPHnyP55HhEFnNpsvX74cGRlJOwgbnDlzZvbs2cXFxcnJyc7B5OTk6upqk8lECDl//rzNZnPP5dDsvH0iPz//yy+/zM3NLSsrs1gsFRUVtBOx1po1aywWi0QiaWxs3Lt378GDB9PS0miHYpWdO3devHjxyJEj0dHRY8aMWbVqVWxs7N27d+Pi4qZPny6TybZu3bpv377U1FTaSYe8np6eN9544969e/v27Vu6dOnw4cM//vjj9vb2zMzMKVOmeHl5HTp0yGaz1dTUiEQi2mGHNpPJFBoaKpFInHdbZWdnT506lRCSlpbmcDhmzpz52Wefvfzyyxs2bHBDHnYWISGktLT09OnTEolk0aJF2Gpdp6GhQaVS6fV6sVicmZkplUppJ2Kb77//vu/FijNmzGD2SPR6/a5duzo6OjIzM+VyOb2A7NHb21tSUuKcFYlEOTk5Npvt8OHDDQ0Ndrs9Ojp63rx5QqGQXkaWuO9XTQhRKBTMC017enp27dql0WgSExMzMzPdk4e1RQgAAPAoWHiOEAAA4NGhCAEAgNNQhAAAwGkoQgAA4DQUIQAAcBqKEAAAOA1FCAAAnIYiBGCJV155RSaTLVu2rO9gTU2NTCarqqqilQrgyYciBGAJvV6vVqsLCwv7PhXaZDKp1eru7m6KwQCecChCAPYIDw+PiIh45513aAcBGEpQhADsIRKJNm7cWFZWVllZSTsLwJCBIgRglYULF44dO3b9+vV4jDDAI0IRArAKn89///33a2trv/32W9pZAIYGFCEA28ydOzcxMTEvL89qtdLOAjAEoAgBWCg/P7+xsfHrr7+mHQRgCEARArDQc889l5KSsnHjRrPZTDsLwJMORQjATlu2bGlpafn8889pBwF40qEIAdhJLpfPnTu3tLSUdhCAJx2KEIC1Nm3aJBAIaKcAeNLxcLMRAABwGfYIAQCA01CEAADAaShCAADgNBQhAABwGooQAAA4DUUIAACchiIEAABO+zuDEMJDp0miOAAAAABJRU5ErkJggg==", "image/svg+xml": [ "\n", "\n" ], "text/html": [ "\n", "\n" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Plot eruption age estimate relative to rank-order plot of raw zircon ages\n", "h = plot(ylabel=\"Age (Ma)\", xlabel=\"N\", legend=:topleft, fg_color_legend=:white)\n", "plot!(h,1:length(ages),ages,yerror=uncert*2,seriestype=:scatter, label=\"Observed ages\")\n", "plot!(h,[length(ages)],[AgeEst],yerror=2*AgeEst_sigma, label=\"Bayesian lockup age estimate\",color=:red)\n", "plot!(h,collect(xlims()),[AgeEst,AgeEst],color=:red, label=\"\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "## Getting your data out\n", "Acccessing the command line with Julia's `;` command-line mode, we can use the unix command `ls` to see what's on the local filesystem and look for any files we have written" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Chron1.0Coupled.ipynb\n", "Chron1.0Coupled.jl\n", "Chron1.0CoupledConcordia.ipynb\n", "Chron1.0CoupledConcordia.jl\n", "Chron1.0CoupledSystematic.jl\n", "Chron1.0DistOnly.ipynb\n", "Chron1.0DistOnly.jl\n", "Chron1.0Radiocarbon.ipynb\n", "Chron1.0Radiocarbon.jl\n", "Chron1.0StratOnly.ipynb\n", "Chron1.0StratOnly.jl\n", "ConcordiaExampleData\n", "DenverUPbExampleData\n", "EruptionDepositionAgeDemonstration.ipynb\n", "EruptionDepositionAgeDemonstration.jl\n", "FCT\n", "Manifest.toml\n", "MyData\n", "PlutonEmplacement.ipynb\n", "Project.toml\n", "archive.tar.gz\n", "ffsend\n", "histogram.pdf\n" ] } ], "source": [ "; ls" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "If you're running in a Binder notebook, you can use a utility called `ffsend` to upload and transfer any files that you want. Edit the cells below to specify " ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "# Download prebuilt ffsend linux binary\n", "download(\"https://github.com/timvisee/ffsend/releases/download/v0.2.46/ffsend-v0.2.46-linux-x64-static\", \"ffsend\")\n", "\n", "# Make ffsend executable\n", "run(`chmod +x ffsend`);" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/bin/bash: ./ffsend: cannot execute binary file\n" ] } ], "source": [ "; ./ffsend upload histogram.pdf" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.9.2", "language": "julia", "name": "julia-1.9" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.9.2" } }, "nbformat": 4, "nbformat_minor": 4 }