![\"Launch](\"https://static.mybinder.org/badge_logo.svg\"){kind=icon}
\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, DelimitedFiles" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "***\n", "## Enter sample information\n", "First, let's enter some basic information about your samples. How many are there, what are the sample names (needs to be a valid filename, BTW), and what are the stratigraphic heights and uncertainties? Then, we'll enter the ages as .csv files." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "nSamples = 5 # The number of samples you have data for\n", "smpl = ChronAgeData(nSamples)\n", "smpl.Name = (\"KJ08-157\", \"KJ04-75\", \"KJ09-66\", \"KJ04-72\", \"KJ04-70\",)\n", "smpl.Height .= [ -52.0, 44.0, 54.0, 82.0, 93.0,]\n", "smpl.Height_sigma .= [ 3.0, 1.0, 3.0, 3.0, 3.0,]\n", "smpl.Age_Sidedness .= zeros(nSamples) # Sidedness (zeros by default: geochron constraints are two-sided). Use -1 for a maximum age and +1 for a minimum age, 0 for two-sided\n", "smpl.Path = \"MyData/\" # Where are the data files? Must match where you put the csv files below.\n", "smpl.inputSigmaLevel = 2 # i.e., are the data files 1-sigma or 2-sigma. Integer.\n", "\n", "smpl.Age_Unit = \"Ma\" # Unit of measurement for ages and errors in the data files\n", "smpl.Height_Unit = \"cm\"; # Unit of measurement for Height and Height_sigma" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note that smpl.Height *must* increase with increasing stratigraphic height -- i.e., stratigraphically younger samples must be more positive. For this reason, it is convenient to represent depths below surface as negative numbers.\n", "***\n", "Now let's see what's in our current directory (we'll use `;` to activate Julia's command-line shell mode, followed by a unix command, in this case `ls`" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Chron1.0Coupled.ipynb\n", "Chron1.0Coupled.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", "DenverUPbExampleData\n", "DispersionExample.jl\n", "EruptionDepositionAgeDemonstration.ipynb\n", "EruptionDepositionAgeDemonstration.jl\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": [ "Equivalently, we can also run unix commands using the `run()` function" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Chron1.0Coupled.ipynb\n", "Chron1.0Coupled.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", "DenverUPbExampleData\n", "DispersionExample.jl\n", "EruptionDepositionAgeDemonstration.ipynb\n", "EruptionDepositionAgeDemonstration.jl\n", "Manifest.toml\n", "MyData\n", "PlutonEmplacement.ipynb\n", "Project.toml\n", "archive.tar.gz\n", "ffsend\n", "histogram.pdf\n" ] } ], "source": [ "run(`ls`);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now that we know how to access the command line, let's make a new folder for our example data. This can be called whatever you want, just make sure it matches `smpl.Path` above" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ ";mkdir -p MyData/" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, let's make some files and paste in csv data for each one. For now, I'm pasting in example data from Clyde et al. (2016), [10.1016/j.epsl.2016.07.041](https://doi.org/10.1016/j.epsl.2016.07.041)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "# You can just paste your data in here, in the following two-column format.\n", "# The first column is age and the second column is two-sigma analytical uncertainty.\n", "# You should generally exclude all systematic uncertainties here.\n", "data = [\n", "66.12 0.14\n", "66.115 0.048\n", "66.11 0.1\n", "66.11 0.17\n", "66.096 0.056\n", "66.088 0.081\n", "66.085 0.076\n", "66.073 0.084\n", "66.07 0.11\n", "66.055 0.043\n", "66.05 0.16\n", "65.97 0.12\n", "]\n", "\n", "# Now, let's write this data to a file, delimited by commas (',')\n", "# In this example the filename is KJ08-157.csv, in the folder called MyData\n", "writedlm(\"MyData/KJ08-157.csv\", data, ',')" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "data = [\n", "66.24 0.25\n", "66.232 0.046\n", "66.112 0.085\n", "66.09 0.1\n", "66.04 0.18\n", "66.03 0.12\n", "66.016 0.08\n", "66.003 0.038\n", "65.982 0.071\n", "65.98 0.19\n", "65.977 0.042\n", "65.975 0.066\n", "65.971 0.082\n", "65.963 0.074\n", "65.92 0.12\n", "65.916 0.088\n", "]\n", "writedlm(\"MyData/KJ04-75.csv\",data,',')" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "data = [\n", "66.13 0.15\n", "66.066 0.052\n", "65.999 0.045\n", "65.989 0.057\n", "65.98 0.11\n", "65.961 0.057\n", "65.957 0.091\n", "65.951 0.066\n", "65.95 0.11\n", "65.929 0.059\n", "]\n", "writedlm(\"MyData/KJ09-66.csv\",data,',')" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "data = [\n", "66.11 0.2\n", "66.003 0.063\n", "66.003 0.058\n", "65.98 0.06\n", "65.976 0.089\n", "65.973 0.084\n", "65.97 0.15\n", "65.963 0.055\n", "65.959 0.049\n", "65.94 0.18\n", "65.928 0.066\n", "65.92 0.057\n", "65.91 0.14\n", "]\n", "writedlm(\"MyData/KJ04-72.csv\",data,',')" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "data = [\n", "66.22 0.27\n", "66.06 0.11\n", "65.933 0.066\n", "65.918 0.087\n", "65.92 0.34\n", "65.916 0.067\n", "65.91 0.18\n", "65.892 0.09\n", "65.89 0.063\n", "65.89 0.15\n", "65.88 0.2\n", "65.812 0.069\n", "65.76 0.15\n", "]\n", "writedlm(\"MyData/KJ04-70.csv\",data,',')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Alternatively, you could download .csv files that you have posted somewhere online using the Julia `download()` function, or using the unix command `curl` throught the command-line interface" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For each sample in `smpl.Name`, we must have a `.csv` file in `smpl.Path` which contains each individual mineral age and uncertainty.\n", "\n", "***\n", "## Configure and run eruption/deposition age model\n", "To learn more about the eruption/deposition age estimation model, see also [Keller, Schoene, and Sameperton (2018)](https://doi.org/10.7185/geochemlet.1826) and the [BayeZirChron demo notebook](http://brenh.in/BayeZirChron). It is important to note that this model (like most if not all others) has no knowledge of open-system behaviour, so *e.g.*, Pb-loss will lead to erroneous results.\n", "\n", "#### Boostrap relative pre-eruptive (or pre-depositional) mineral age distribution" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdZ1wUV9sG8DMzuxSlN2nSmzQJgoooNhR7i2KILTZiNyYxJmp8TDWJRn0SNZbYjUkAEwVFEuwtsXdFEJAqoDQpS9mZeT/s866Ii4JuZa//Lx92Z2d37z2SvfacOXOG4nmeAAAAaCta1QUAAACoEoIQAAC0GoIQAAC0GoIQAAC0GoIQAAC0GoIQAAC0GoIQAAC0GoIQAAC0GoIQAAC0GoIQAAC0mqYGYXV19dKlS5u/P8uyiisGZEKbKx/aXPnQ5son9zbX1CAsLi7etWtX8/evrq5WXDEgE9pc+dDmyoc2Vz65t7mmBiEAAIBcIAgBAECrIQgBAECrIQgBAECrIQgBAECrIQgBAECrIQgBAECrIQgBAECrCVRdQGtzt4z//Cr3Vy5nKKQMBMRQhxgKiV0b6kN/2teUUnV1AADQGIJQblLL+c+vckfyuAW+zH+7CkUsX1lPKupJZT25XsKHJ4ojXejPAhlTXVUXCgDywzCMqkuA14UglIOsSn7ZZS4pl5vvw/wUKjQUSjY/7f+F21GTPehPL7PecfWfdWKmedI0OocArYK+vr6qS4DXhWOEr+tBBd8jgXUyIGmRwsUB9P+nYGNmumR9NyZpoGBvOhd8QJxWziu3TABQCIrCr1qNhyB8LQUi0u8wu6gj/VknxqiJCGyooxl1YrBgigcdkcTmVSELAQBUD0Ojr668jgxKEk/2oGd7t+z3xGxvmuNJn0T29BCBFYZVAABUCj3CV1QtJkP+FvewphYHvEobzvWhRzpRQ/8WV9bLvTQAAGgBBOGrqOPIm0fErobU2pBXnzC2IpgJMKeGJ4trcV1PAADVQRC+ioknWEMhtTWMeZ2j5BQhG0IZM11q/AmWxeFCAAAVQRC22L5M7nYpv6c381oxSAghhKHInl5MWR3/yUX0CgEAVANB2DIV9eS9f7kfuzE6cmo5XYb81kew5z53pgC9QgB5qq2tLf1/qq5FgRITEy9evPiCHVJSUn7//XdCSE1NzerVq1/hLUQi0dq1ayW3t27dmpub+wovUlFR8eOPP77CE5UAQdgyyy6zEfZULxt5njlkrks2hjLTTrM16BYCyM/WrVttbW2DgoI6duxobGy8aNEinm/xz83JkyffuXOnqUfFYnG/fv1EItHrVfpa4uLizpw584Idbt68+fPPPxNCRCLR999/39RuK1eujIuLk/lQdXW1NEHXrFmTmZnZzNrmz59//vx5ye2Kior//ve/zXyikiEIW+BGCf9rOrciWP4rKg1zpP3NqM+uIAkB5KlLly7p6enZ2dnXrl3bsGFDw8DgOK6mpub5pzRKtX/++ae8vLzRDtInchx35MgRlm38f25VVVV9fb30Nsdxz79RdXX189t5nq+qqnppVRKVlZXPb2xURsPsNzU1zcvLk97lOK7hR7tz505WVlajVxCLxZWVlebm5tnZ2Q23syzbnPi/ePHi48ePJbdtbW3v37/f8NGmXuGln0vuEITNxfFkxhl2RTBjqaeQ198QyuxM4y4+wgApgPw5OjoaGxtLvmFZll24cKGtra27u3uvXr2k3/6//vpr+/btO3To4OHhcfz4cULIwoULMzIyRo8e7erqumfPnuzs7KCgIB8fH29v78DAQELIm2++SQjx8/NzdXW9ePHijh07Ro8ePWDAAHd39927dycmJrq6ur7xxhs2NjazZs2SxN7Jkye7d+8+btw4f39/Gxub7du3S959+PDhCxcu9Pb29vLy6tWrV0FBgWR7bGysq6urv7+/q6vr4cOHJRtPnDjh5OTUsWPH3r17FxcXP/95CwoKwsLCvLy8vL29z549K9lYWlpqaWkpuf3RRx/Z2NiEhITY2dn9/fffe/fu3bdv34oVK1xdXd9//31CiJmZ2ddff+3k5DRs2LDHjx9bW1tLXzwhIaFDhw5ubm5Dhw6VROnevXsnTJggeVQsFpuZmVVVVX311VdXr16dNm2aq6vrunXr8vPzHR0dJftcu3YtICCgQ4cONjY2a9askWycP3/+J598EhIS4uXl5e7ufu/evdf9V28+XjNlZ2e3b9+++fs/efLkNd9xwx22e0I995qv8kJ777MdYutrxIp8DyV6/TaHllLDNr9TyiXnKuO/o3lcPfvMW69fv759+/aLFi1auHBhnz593n777bq6Op7nt23b1qFDh5KSEp7nFyxY0LdvX57n7927Z2BgcO3aNZ7n4+Pjzc3NS0tLeZ739PQ8d+6c5AUXLlz44YcfSm7n5eXxPF9bW0sIqaiokGzcuHGjUCg8e/Ysz/Msy2ZnZ0seqqysDAkJiY2N5Xk+OTmZEPLLL7/wPH/79m1DQ8O0tDSe5/v06ePv719eXs5xXHR09Lhx43iev3Hjhq2tbXp6Os/z169ft7a2fvz4cVVVlbW1dUJCAs/zly9fFgqFq1evbtTsb7311syZM3meLysr8/X1DQ8P53n+8ePHFEXxPH///n1ra2tJf7Gmpubx48c8z7/zzjurVq2SvgIhZPr06SzLsixbWFgoEAgk2318fMLCwkQikVgsHjly5Pvvvy9p0pEjR0p2kHSFJR88JCTk4MGDku3Z2dlt27aV7ODh4bFu3Tqe5x88eGBpaXnixAme56dNm+bi4vLw4UOe5z/44IMJEyY09Ucl979zrCzTLIUi8p/L7JFBAoWuKhjlSsdl8l9fYz/rhPXsoZXYncZdfKyMcQ6KEHdjpn3bZ/4f1dfXd3Fxqaurq6ioOHPmTG5urrOzc0JCQnR0tKmpKSFk0aJFNjY2VVVVSUlJ4eHhHTt2JIQMHTq0Xbt2586dGzRoUMNXk3TgAgICBgwYYGtrK7OMkJCQbt26EUJomra3tz9y5MiFCxcqKiooirp06dLo0aMJIdbW1lFRUYQQb2/v8PDwpKSkOXPmEEKmTp1qZGRECFmwYEHXrl0JIbGxsYGBgRkZGRkZGYQQS0vLCxcutG3b1tDQcMiQIYSQwMDAXr16PV/GoUOHLl++TAgxNjaeMmVKYmJiw0eNjY1FItHKlStHjRrl5+enqyv7mjgffPABTcsYNZw1a5aenh4hZN68edHR0S847ihTSkpKYWHhzJkzCSGOjo5jx46Nj4/v2bMnISQqKkrS9YyIiFi8eHGLXvZ1IAib5YPz7DRP2t9M4avrrg9lOv5RP9yRDrTASr7QGnytgGPqzWdjYxMdHS25PWnSpG+//Xbjxo0lJSUWFhaSjZIbxcXFxcXF5ubm0idaWFhID25JzZs3z8jI6Jdffpk+ffqbb765Y8eO599ROvZICFm6dOmJEyeio6PbtGmTnp4uPfRlZmYmXarbwsJCOrZpZmYmuWFubv7kyZP6+vrCwsK8vLzY2FjJ9pCQEFNT04KCgoalNrwtUVdXJzmw1/AzNmRhYXHs2LFNmzYNHDiwTZs2MTExAQEBL/4sDTWss6SkpNGj/MtmJJWUlJiZmUkj1sLCQjr7RvI7gBCio6NTV1f34teRIxwjfLkLj/jTBfzSN5Tx/7O1PvmuMzPjLM6wB5AzmqYlUeTi4nLz5k3Jxps3b+rp6dnZ2bm6uko31tXVpaSkuLm5EUKEQqF0LgzDMFOnTk1MTMzNzT127Ni5c+cYhqFpuuFkmYYXo/jzzz9XrVo1adKkMWPGNNwnKytLGoq3bt2SvBEh5Pbt29KNDg4OQqHQ09PT1tZ2UwNdu3Z1cXG5f/++ZFRWsnOjT6qjo2Nvb9/w1Z5vjcDAwE2bNmVnZw8ePFjSpRMIBM/P+pFJ+sq3b992dXUlhJiYmEgTMS0tTbpnw9aTcnFxyc/Pl+5/48YNyYuoEHqEL/fNde4jf7qNsppqkju98S4Xm8FFuuBnCsBrefjw4ebNmzmOu3v37u+//37w4EFCyNy5c/v27evl5eXg4PDxxx+/9957DMOMHj162bJlS5cuHTBgwJYtW9zc3EJCQgghPj4+P/30U35+fmBg4OHDh83Nzd3d3bOysurq6pydnRmG8fLy+uabbwIDA3v37t3o3d3d3Tdv3kzT9PHjx0+fPj127FjJdoZhoqOj58yZc+TIkdzc3FGjRkm2792719fX18bG5r333ps3bx4hZNq0aT/++OPChQtHjhxZVVV15MiR999/39/f38/Pb+bMmdOmTYuPj284EVRq3rx58+bN+/777/Pz83/77TcPD4+Gj16/fv3gwYNhYWEMw9y4caNv376EEF9f399++619+/aurq6dO3d+Qatu3LjR0dFRT09v8eLFX375JSGka9euEydO3LZtm729fcOTBX18fLZv3y4Sifz8/AwNDSUb7e3tR4wYMWnSpIULF164cOHEiRPr1q1r/r+pIjDLly9XbQWv5smTJ9u2bVuwYEEz96+rq2tqHPzF7pXzn11hd/QUyOsM+uZwN6IWnOdmdtDs6/e+cpvDK0ObN8SyLEVRpaWlZWVlktmJQUFBhBAbG5u+ffvGx8efPXs2MjLyww8/pChKKBRGRUWdPn06MTHRy8vrp59+krRkr169ioqK0tPT27dvb2JikpiYeOjQocLCwtWrV3fo0IEQEhERkZqampWV5e/vb2hoaGVl5e/vLymgT58+58+fT0hIcHZ2njVrloODg5eXV0ZGxpUrV+bMmbNhwwZCyLZt2yRDl7t27ZozZ86NGzeSk5PffvttyVFDXV3dSZMmXbly5c8//7x7966Pj09ISIhAIBg5cuTly5fj4+M7deo0duxYT09PBweHhp89JCSEoqhff/21trZ2yZIl1tbWkmmuIpEoIiKCoqh//vknISHhwoULAwcOfO+992iaDgwMFAgEd+7c0dPT8/X1ramp6du3r1AoJIRI5tT079+fEFJbWzt37twDBw6cPn167ty5koOdhoaGISEhMTExd+/eXbZsmZGRUXh4OMMw3bt3LysrS0tLa9eunZOTU11dXXh4OCFk+PDhBQUF+/btY1l269atkuIlk2icnJwk76ivry85UPo8uf+dUy8dz1VPOTk5oaGhjU5teYGKigrp75EWmXKKdTGklr6h7M5ZeKJ4rAs93UuDO4Wv3ObwytDm6u/IkSOLFi2SzGRpqG/fvjNnzpTMpoEXk/vfuQZ/zypBXhW/P4ub2cLLDcrFN52Zz65y1WLlvzMAgHbBMcIXWX2Lm+JBm6tiqCnIgupiSW24y33ohx8rAK1HeHj4891BQsjRo0eVXwxI4Eu2SaW1ZEcqN89HZU30dTD93XW2tFZV7w8AoBUU1SOsr6+/devWjRs33N3dJaeXNlJeXp6UlJSSkmJubj569Gjp+j0xMTFlZWWS2zY2NkOHDlVQhS+17g43wpF2MFDZfBVPY2qIA73mFvs5zq8HAFAYRQXh+PHjL1y4IBaLR40aJTMIx44dKxQKg4KCzp8/v2TJknPnzvn4+BBCli9f7u7uLsnFRlN+lalaTNbfYY8NVvHQ8Red6IA/xbO9mXb6qi0EAKDVUtQX/a5du3R1dWfPnt3UDjExMdJFBCIjI3/++Wfp0qsfffRRaGioggprpm2pXGg72ttExacv2LWlxrvRX19j/xuCTiEAgEIo6gDYS0/ykKYgIYSm6Yb7//3335s3b75w4YKCanspMUdW3+Q+8leLA6iLA5hf7nO5VRp5lgsAgPpT/azRs2fPJiUlXblyRXLXw8OjsLDw4cOHS5YsiYyMXL9+vcxn1dfXl5eXf/jhh9ItvXv3lqyPIFNNTY3kzNDmiH1AtW9DdTQSy7pambIZEjLehfr+Wt2KThqWhS1qc5ALlbc5wzAMo12jF2KxWCBQ/Replqivr5ec3d/8v3OhUPjSv0kV//vdvXt3zJgxW7dudXFxkWzZv3+/5MaSJUu8vLymTZv2xhtvPP9EiqJompau/UoIsbCweMGnbdH/n9vTycwORH3+f17gSzrFk8UBxERH1aW0hBZ+J6qcOrR5M9erbDUqKyuxiIHSSJbqbtHfecPVX5uiyiBMS0vr37//t99+K7m4ZSOOjo6urq6pqakyg1AgEBgaGjb/Oh1CobCZvyAeVPA3S8UjnYVCtfkOdzQmg9qz29MpNRmtbabmtznIC9pc+dDmyif3NlfqF2tNTc2lS5ckPxizsrIiIiKWLl0qva4xIUQsfrqSSmZmZnp6uqenpzIrJIRsT+XGudG6apOCEgv96R9uc3WcqusAAGh1FNUj3L17965du1JSUiiKunPnzuTJk99+++3MzMzg4ODS0lITE5NJkyaVlpbGxcXFxcURQkJCQj7//PNr165NmjSpa9euPM8fOHAgOjpa5lWyFIfjyc40fn8/NYtBQvzMKF9Tsvc+946HJnUKAQDUn6KCMDQ01MbGRnpXcrkpR0fH5ORkAwMDQsh333335MkT6Q6SK0AGBARs3rw5NTWVYZgFCxb4+fkpqLymHMnnLfRIgLk6XvRhoT8z5xw70V2zL0kBAKBuFBWELi4u0vkvUm3atJFcg4MQIvN6VwKBIDQ0VIUnEW67x01W1y5XX1uqrYAk5fKD2iMJAQDkRk2/9FWipJb8lcu97aq+bfKhH73yhnZNyQMAUDT1/dJXvl/ucwPb06ZqfFnTMS50ThU5X6RhJxQCAKgzBOFT21O5Keo6LirBUGSeD736FiaPAgDIjVp/7yvTtWK+pJb0sVX3w29TPOjj+VxGBTqFAADygSD8n22p3DselPpPyDQQkqme9Lrb6BQCAMgHgpAQQuo48ls6N9FdM1pjZgd6932uSvzyPQEA4KU046tf0f58wHU0p1wM1b4/SAghxMGA6taO/i0dnUIAADlAEBJCyO407h0N6Q5KzOpAb7iLIAQAkANN+vZXkLI6craQH+6oSU3R356qrCf/4jwKAIDXpknf/gqy/wHX14420Kjl4ylC3vWif0KnEADgtSEISVwmN8ZZM44ONjTFk07I5opEqq4DAEDDaXsQltWRM4X8oPaa1w4mOmSkI709FZ1CAIDXonkBIF/xWVxfW9pQo8ZFpeb60D/d5VgcKAQAeA3aHoSxmjkuKhFgTlm3IYdzkIQAAK9Oq4PwST05XcAPctDgRpjVgd5wF9ejAAB4dRqcAa/vQBbX24Y20sxxUYmxLvSVx3xaOTqFAACvSKuDMC6TH62x46ISugyZ5E5vuYcpMwAAr0h7g7Cinpx6yA3VqPPoZXq3A70jlavF+CgAwCvR+Bh4ZQeyuJ4aPi4q4WJI+ZlR8dnoFAIAvArtDcK4TF5z54s2Ms2T/jkFQQgA8Cq0NAgr6snJh9wQTZ4v2tCbzvS1Ej79CabMAAC0WCtJgpZKyObCbGhjHVXXISc6NIlypXemoVMIANBiWhqEsRmtZ1xU4l0vensqj1VmAABaShuDsEpMjuW3nnFRiQ4mVPu2WGUGAKDFWlUYNFNyHtfZijJpLeOiUlM96a04oRAAoIW0MQgPZvNDNPByEy/1lit9qoDLr0anEACgBVphHrwYT8jhHH6wQ6s6QCjRVkBGO9M70xCEAAAtoHVBeOkRb6xD3IxaYRCS/z+hkEMUAgA0m9YF4cFsbmhr7A5KBFtSxjrkxEMkIQBAc2ldECZk861svmgjUzzpnzFlBgCg2VpzJDwvr4rPruRDrFptj5AQMt6NPpzDFdequg4AAA2hXUF4MIcf0J4WtOoPbaJDBrWn995HpxAAoFladSY852A2N6R9a+4OSkz2wOgoAEBzaVEQisTkdAEfYd/6P3JfO6qynlwtxpQZAICXa/2pIHUknws0p0x1VV2H4lGETHCntqeiUwgA8HJaFIQHW/t80YYme9C/pnM1uGw9AMDLaEsw8IQk5vBDWu8ZhI04GlAB5tSBLHQKAQBeQluC8EYprccQD2NtCUJCyGQPGqOjAAAvpS1BmJhHD3fUohQkhIxyoi894rMrMWUGAOBFtCUIk/Jp7TlAKKHHkLGuWIMbAOAltCIbCkQks4oKbaddPUJCyBQPensqluAGAHgRrQhCAUVWBYqFWvFZn9HJgjIUkuP5iEIAgCZpRThY6JFIRy09k2CyB70DU2YAAJqmFUGozca70QnZXHmdqusAAFBXCMJWzkKPhNvRv6ajUwgAIBuCsPXDCYUAAC+AIGz9BthTD6vJjRJMmQEAkAFB2PrRFJnoTu1MQ6cQAEAGBQZhbW1tRkZGVVVVUzuIRKITJ05cunSJ55/prKSlpR05cqSoqEhxtWmbqZ707jSuVktnzgIAvIiigrBPnz5GRkYeHh5//fWXzB3S09M9PDy++OKLSZMmDRo0qL6+XrL9iy++CAsLW7t2rbe394EDBxRUnrZxNqS8TalDOegUAgA0pqgg/PLLL4uKijw9PZva4Ysvvhg+fPjRo0evXLmSk5Pzxx9/EELy8vJWrFhx5syZgwcPbt68+b333uM4fHfLB6bMAADIpKgg7Natm7Gx8Qt2+OOPPyZOnEgI0dXVHTNmjCQI4+Pjg4ODXV1dCSHDhg0rLS29cuWKgirUNpHO9LlCPqcKU2YAAJ4hkLn1/fffHz16dLdu3RT0rmVlZRUVFY6OjpK7jo6Of//9NyEkNzfXycnpf5UJBLa2tjk5OUFBQc+/As/zIpEoJiZGuiUoKEj63OdxHKflnUtdmoxyonancR/7K2nNVbS58qHNlQ9trnwtanOafnl/T3YQJiUlrVmzxs/PLzo6esKECS/u272CmpoaQohQKJTc1dXVFYlEku0CwdOSpNufV1tbW11d/fvvv0u3iEQia2vrpt6xtrZW+nZaa5wjNf1fwXyPeuUkIdpc+dDmyoc2V74WtbmOjk7DWJFJ9sM3b948fvy45CjdwoULhw4dGh0dHR4e3rJim2ZpackwTHFxsZmZGSHk0aNHNjY2hBBra+v09HTpbo8fP5Zsf56enp65ufm+ffua+Y4sy7Zp0+a1C9dsvRyIzgXxlSf6PayVEYVoc+VDmysf2lz55N7msvuMDMOEh4fHxMRkZWUtX778/Pnz/fr18/b2/vbbb0tKSl7/XRmG6dy584kTJyR3T548GRISQggJCQk5d+6cWCwmhNy/f7+4uDggIOD13w6kJntiygwAwLP4Zqirq3v//fcl+xsYGMybN6+goODFT9m7d+8333xjZWU1fvz4b775Jjc3l+f5TZs2BQcHS3bYt2+fpaXlrl27PvvsMzMzs/z8fMn2rl27jh07NjY2tmvXrnPnzm3q9bOzs9u3b9+c4iWePHnS/J1bscJq3mRn3ZM6ZbwX2lz50ObKhzZXPrm3+UtGTktLS3ft2rV58+Y7d+74+Pi8++67eXl5W7Zs+eOPP27evGliYtLUEysrK0tLSydPnix5EUknLzAwcMqUKZIdRo0apaurGxcXZ2RkdPbsWekQ6OHDh9euXXvgwIFx48bNmDHj9ZMeGrLSJ71s6JgMbqonFhUCACCEEIrnZc+nP3PmzObNm2NjYzmOGzVq1MyZM8PCwiQPPXz40N3dffv27WPGjFFiqc/IyckJDQ3Nzs5u5v4VFRWGhoYKLUlTHMzmv77Gnhv2kt9Arw9trnxoc+VDmyuf3Ntc9rdht27d/vnnH0dHx08//XTq1Knt2rVr+KiNjY2Dg0N5ebkc6wClGdiemn2OXC/hO5op6TwKAAB1JjsI3d3dFy9ePHDgQIZhZO7w77//6urqKrIwUBSGIpPcqe2p3Nqusv9xAQC0iuwDRVOnTu3Ro0ejFCwvLz9y5IjktpGREYJQc033ovfe52qwBjcAQFNBGBkZefv27UYb79y5069fP8WXBArXvi31hjn15wOcRwEA0JK1Rmtra/X09BRXCijTNC/653sIQgCAZ48RFhYW5ubmEkLq6+vv3bvXcPCzpKRk06ZNzs7Oyi4QFGO4Iz33HJtWzrsbY8oMAGi1Z4Jw79690hPnpSf8SbVp0+bnn39WUl2gYDo0meBGb0vlVgRjygwAaLVngvDNN9/08/MjhERGRi5fvtzb21v6kLm5uaurq5GRkbILBIWZ7kX3PCj+vBMjxLn1AKDFnglCBwcHBwcHQsjOnTu7detmbm6uoqpAGTyMKXdj6lAON8IRSQgA2kv2N+DQoUORgtpguhe9JQVTZgBAqz3tESYnJ3///fdTpkyJjIyMiooqLS2V+YSkpCRl1QYKN9qJXvAPm13JOxhgygwAaKmnPUKWZWtrayWrY9fV1dU2QXWlgvzpC8hbrvSONNnrzQIAaIMmF91Wc1h0W16ul/DD/mYzxgoYefcJ0ebKhzZXPrS58sm9zTFLQtt1NKPs25JD2ThSCABaSnYQxsXFxcfHS26XlJSMHTvWyclpzJgxjx49UmJtoCQzO9A/3UUQAoCWkh2E8+bNk2bexx9/HB8f371793///VdyoV1oZcY401eL+bRyjRwkBwB4TTKC8MmTJw8fPgwKCiKE1NfXx8TEfPzxx3v27ImLi0tMTCwuLlZ6kaBYugyZ5E5vwdKjAKCVZARhVVUVIcTY2JgQ8u+//5aXlw8bNowQ8sYbb/A83/z5KaBB3u1A70jFhZkAQBvJCEIrKyuhUHjjxg1CSExMjKWlZceOHQkhJSUlhBBcgKJVcjGkgiyp2Ex0CgFA68gIQoZhxowZM23atHHjxm3ZsuXtt9+maZoQcvnyZaFQKFmDDVqfmR3on+4gCAFA68ieLLNx48bIyMjMzMxp06Z98cUXko1///334MGD27Ztq8TyQHkGt6cfisiVx5gyAwDaBSfUw1NfXeNyKvmN3eVzYSa0ufKhzZUPba58OKEeFCjai47L5MrrVF0HAIASCWRura+v37Fjx/79+7OyskQiUcOH0tPTlVIYqIClHgm3o/fc52Z74xcSAGgL2UE4ffr0nTt3BgQEBAYG6uvrK7kmUKGZHejZ51gEIQBoDxlBWF9f/9tvvy1btuyzzz5TfkGgWj1tKIqQY/l8H1tcmAkAtIKMH/4lJSW1tbUjRoxQfjWgDub70mtv4TwKANAWMoLQ0tLS0dExNX7V5rMAACAASURBVDVV+dWAOpjgRl98xN3D0qMAoB1kBCFN05s2bVq+fPnly5eVXxConC5DpnnR63FyPQBoB9mTZZYtW5aXlxcUFGRhYWFkZNTwIcwa1QZzvZkOcfXLAxkzXVWXAgCgYLKDMCwsLCAgQMmlgPqw0ieD29PbUrkP/TB9FABaOdlBuHLlSiXXAermPV961BH2PR9agCgEgFYNX3IgW6AF5WhA9mfhSCEAtHJNBuGBAwe6d+9uZmZmb28v2fLdd9+tXbtWWYWB6uE8CgDQBrKDcOfOnSNGjNDX1x8+fLh0o7W19YoVK1gWF2/VFsMd6QIRufAI51EAQGsmIwh5nl+8ePH8+fOTk5Pfeecd6fbQ0NCioqK8vDzlVQcqxVBktjf9X3QKAaBVkxGEhYWF+fn5kydPbrTd2tqaEFJUVKSMukA9TPWk/8rl8qrQKQSAVktGEOro6BBCGl10ghDy4MEDQoixsbHiqwJ1YSQkE93pNegUAkDrJSMIzczMvL29N2zYwPM8Rf1v5WWe57/99lt7e3s3NzflVggq9r4fvS2Ve1yj6joAABRD9nmE33777fDhw/Py8ry9vUUi0Y8//hgbG3v69Ondu3dLoxG0hH1bapQTve4OuzxQPleuBwBQK7JnjQ4ZMuTAgQP5+fnr1q0rKSmZN29eZmbm7t27x48fr+T6QB180pHecIerqFd1HQAACiC7R0gIGTJkyJAhQ3JzcwsLC42MjNzc3NAX1FquRlQfW3pzCvcBVlwDgFanySCUsLe3l55QD9rskwB6YJJ4tjeth/FRAGhdZPzAz8/PX7p0aY8ePZydnZ2cnLp27frBBx/gohNarqMZFWhO7UzD9FEAaG0aB+GpU6e8vb2/+uqrGzduGBkZmZqapqWlrV692sfHJy4uTiUlgppY+gbz3XVOjCgEgNblmSAsLy+PjIw0MjI6dOhQaWnp9evXr169+vjx41OnTrm5ub3zzjtZWVmqKhRUrqsVZd+W/J6BJASAVuWZINy7d29xcfHhw4cHDRpE0/97iKKoHj16JCcnCwSCzZs3q6JIUBefBDBfXeM4rDMDAK3IM0F44sSJ/v37+/j4PL+fjY3NW2+9dfz4cWUVBupogD3VVkAOZqNTCACtxzNBmJaW1qlTp6Z27dSpU2pqaotevaamhufRfWhVFgfQn19FnxAAWo9ngvDJkycmJiZN7WpqalpeXt7M1y0uLo6IiLC2tjYzM1uzZs3zO4wfP96sgW7dukm2d+3aVbpx6NChzf4goCQjnWgBTfZlolMIAK3EM+cRisXiF5w1T9O0WCxu5ut++umnRkZGxcXFGRkZXbp06dOnT8eOHRvusHHjxvr6/y1VEhkZ2aVLF8ntJ0+e7N27V3JXIHjJaY6gEp93YuaeY0c40gKcXg8Amq9x0hw8eLCwsFDmrs0fF2VZds+ePUeOHGEYxt3dfdSoUbt27fr+++8b7mNgYCC5UVBQcOrUqY0bN0ofMjQ0NDU1be4nAKXrb0fZtSF77nPveCAJAUDjNQ7CY8eOHTt27DVftKioqKKiokOHDpK7HTp0OH36dFM7b9++PTQ01NXVVbpl5MiRYrE4ICBg5cqVLzhmyXFcaWmp9K6BgYFQKHzNyqGZvunMjD3GRrnSulhoBgA03DNBeP36dY6Tw7GfsrIyiqLatm0ruWtgYFBSUtLUzrt27VqyZIn07saNG/38/HieX7NmzcCBA1NSUszMzJ5/lkgkKigocHFxkW756KOP5syZ09S7VFZWvsongSZ00CMdDHXW3aiPdmtytBxtrnxoc+VDmytfi9pcT0/vpX2kZ4JQXhfdtbS05HleOvWmrKysXbt2Mvc8derUw4cPR40aJd0SFhYmufHFF1/s2bPnzJkzw4YNe/6J+vr6tra22dnZza/K0NCwBZ8BXubbrny/w+JoX33Dpv/G0ObKhzZXPrS58sm3zRVyjMfc3Lxdu3ZXrlyR3L169ap0mLSRbdu2vf32223atHn+IY7jxGIx5suoLT8zqo8tve4Opo8CgGZTSBBSFBUdHb1s2bIHDx4cPHgwMTFxypQphJDs7Oy+fftWVVVJdqusrNy3b5/kIYmsrKxdu3alpqbeu3dPMs7ZvXt3RVQIcvFZIL36JltSq+o6AABeg6L6W0uWLKmsrIyIiDAxMfn111+dnJwk2xuennHx4sURI0YEBQVJt9A0HRcX99VXXwkEgsDAwGPHjhkZGSmoQnh97sbUSCd69U32yyDMmQEATUVp6MovOTk5oaGhzT9GWFFRgXF8Rcit4gP+EF8ZKXAwaHwGKtpc+dDmyoc2Vz65tznOA4PXYt+WmuNDL7yAI4UAoKlkB+F///vfvLw8JZcCGmqRP3PxEX/ioUYOLQAAyA7ClStXOjo6jhgxIjExUS5nFkIrpi8g33Wm55xjcc1eANBEsoPw1q1bGzZsSE9PHzx4sIODw8cff9yiM/ZA24x2pu3akE0pSEIA0Dyyg9DExCQ6OvrmzZuXLl0aPHjwDz/84Ozs3K9fv9jYWJZllVwiaIQ1IcznV9nHNaquAwCghV4yWaZTp06bNm3Kzs6eO3fukSNHIiMj3d3df/jhB5FIpJz6QFN4m1CRzvR/ruB3EgBomJcEIcuyCQkJkydPXrdunamp6fz584OCgj744INu3brV1ODHPzzj807MH5nc9RLMmgEATdJkEObk5CxfvtzJyWnYsGGFhYVbtmzJy8tbu3ZtTEzM5cuX7969+9dffymzUFB/prrkP4HM7LMskhAANIjslWWioqJiY2N1dXWjoqJmzJjRcPEXQoi/v7+zs/OjR4+UUiFokule9JZ73K40bpI7TlEFAM0gOwgfPnz4/fffT5o0SXL5iOfFxcU1dUEJ0GYMRbaHMeGJ4nBbCovjAYBGkB2EO3fubNeunZ6eXsONNTU1+fn5kksA+vj4KKM60ED+ZtTMDvSMs+zeEFWXAgDQDLLHr7p06SK9iJLU1atXG15HHqApn77B5FWRuGysxA0AGqAFB3Lq6+tfep1fAEKIgCZbw5hPrgoKcZYNAKi9Z4ZGa2pqJCcIchxXUVFRWloqfaisrCwmJsbe3l7ZBYJmesOcGufMzj3HxvRFvxAA1NozQfjTTz+9//77ktsDBgxotCtFUatWrVJSXaD5PvER9zoqjMvkRjtjBikAqK9ngjA8PHzTpk2EkIULF86YMaPhEUEzMzNfX18vLy9lFwgaS5chW3sww5PFPW1oS72X7w8AoBLPBKGfn5+fnx8hpK6ubuTIkXZ2diqqClqJLlbUBDf63TPsvnCm8XV7AQDUg+wxqzlz5iAFQS6+DmYKqvkfbuHCFACgpp72CP/999/du3cPHz68f//+n3zyyZMnT2Q+Yf369cqqDVoDIU1+68N0PiDubEWFWKFbCABq52kQ5uTkHDx40MfHp3///snJyU2toIYghJZyMKB+7iGIOsZeGSkw01V1NQAAz6J4XiNXSM7JyQkNDW3+5YIrKioMDQ0VWhI00qjN5//DZlaQA/1xsFCB8HeufGhz5ZN7m2NeOyjJqi7Moxp+3W0cLAQA9fJ0aLS+vr45l9s1MsJayvAqJAcLuxwQB1tSXXGwEADUxtMg3LNnz5QpU176BA0dSgV14GhAbQhloo6z/wwTWOuruhoAAEJIwyDs1q3bhg0bVFgKaINRTvTNEjLkL/GpIYI2sq99AgCgVE+/ijw9PT09PVVYCmiJZYF0+hN+4kk2pg9DY4gUAFQNk2VA2ShCfg5jimv4Ty+zqq4FAKBBj/Dy5ctxcXEDBgzo2bPn119/XVFRIfMJK1asUFZt0Grp0CQuXBASL3Y04KK98GsMAFTpaRDeu3fvp59+srKy6tmz544dO4qKimQ+AUEIcmGuS+L7MT0Pid2MqD62GCEFAJXBCfWgKM1p86P5/Ljj4uSBAj8zZKEc4O9c+dDmyocT6qFV6WtLbQhlIpLEd8o08gcZALQCTU5gF4vFiYmJV69ezcvLs7a29vX1HTZsmJ4eLisHcjbKia4Sk/BE8fHBAk9j9AsBQNlkB2Fubu7QoUOvXbvGMIy5uXlpaWl9fb2rq2tCQkKHDh2UXCK0ehPcaDFH+h9mTwxmnA2RhQCgVLKHRidOnJibmxsTE1NTU1NYWFhTU5OUlMRx3KhRozgOa0WC/E32oBf60/0Ps/nVGCMFAKWSEYQlJSXHjx/fuHHjmDFjBAIBIYSm6YiIiL1796akpNy9e1fpRYJWmONNR3vRfRPZgpcveQsAIDcyglAyj9THx6fRdskW9AhBcRb605Pc6e4J4owK9AsBQElkBKG5uXlQUFBCQkKj7fHx8fb29jhGCAr1cUd6aQAddpC9VowsBABleDpZpqKiQnpV+s8//3zy5MlZWVmjRo2ytrZ+/PhxYmLitm3bfvjhB8lgKYDivONBm+iSgUni2L6C7taYOwMAivX0hPrt27dr0GWYcEK9+nvNNj/+kB97VLypOzPSCWe7Nhf+zpUPba58cm/zp9273r17x8TEyPGlAV5HbxvqYIRg+N/i8jryjgeyEAAU5WkQOjk5OTk5qa4SgMY6W1LHBwuG/s3eKeO/CcY1mwBAIfBDG9Salwl1frjgymN+8F/i8jpVVwMArVGTM1+uX78eFxeXkZHR6DIUycnJiq8K4CkzXZI0QPD+ebbzAXFCf8YDy7ABgFzJDsJ9+/a99dZbFhYWHMfp6+vr6OhkZGQYGBgEBwcruT4AQoiAJj+EMOvvcD0Pin/rI+hpgyxUO0/qSfoTnmswl06XIa6GlD6mmYPak/1H+sknnwwZMiQmJiY6OtrOzu7LL79MSUmJjIyMiIhQcn0AUrO9aU9j6q1j4g/8mQ/8cMRQxSrqyYmH3LVicr2Ev/qYL6rhXQ0pQYODLbUsyajg/cyoECuqWzuqmxVl1xb/aKCOZARhdXX1/fv3d+7cKRQKCSF1dXWEEC8vry1btvTu3XvGjBkGBgbKLhOAEEJIuB11YYRg7FH2TAG/PYwx1VV1QdpHzJHkPH7PfS4xhwu2pDpZUGOcqa+CaDcjinku5kRicukx/08Rv+c+P+cca6NPLfCjo1xpHUxOAHUiIwirqqp4njczMyOEmJmZFRcXS7b7+/uLRKLU1NTAwECl1gjQQPu21Mkhgo8usEH7xb/3ZYIs0MlQktul/JZ73O/pnLMhNd6N/qGb0PxlP0T0BaSHNdXDmiKE8IRJzuO/v8EuvsjN8aHf9aLN8DsG1IOMH2YWFhZt27aVnKvu7u5+7NgxkUhECDl9+jQhRBKQACokpMmarszKLvTgv8Qb7mDxW4W7W8a/dYwNTxSb6JDTQwXnhglmedMvTcFGKEL621F/DRQkDWDSynn3mPoF/7IV9YqpGKAlZAQhRVG9e/fev38/ISQqKqqsrKxjx47Dhw8fOXJkt27dHB0dlV4kgAyjnOhzQwXbU7mRyezjGlVX00qllvPjT7C9D4kDLaj7Y4XLAxk3o9ftgvuZUdvCmNujhVVi4v+H+EieWixWBdpM9lD9xo0bZ8+eTQgxNjY+cuRIUFDQo0eP3nnnnf3791NUs/43ePTo0dSpU4OCgqKiomQuhLZ+/frI/zdhwgTp9rS0tNGjRwcFBc2aNausrOyVPhRoC1cj6uwwgacJCfhT/Fcuvk/lqUhEppxiexwUe5tQ9yOFH/nTbeU6/9Nan2zuzmzqzkw7zc44g64hqJLsILSzs/P29pbc7tSp0969e8+dO7d+/XpLS8tmvu748ePFYvGuXbtsbW2HDRv2/AqlFy9ebNOmzZgxY8aMGTNy5EjJRpZlBw4c6O3tvXPnzsePH0dHR7/ShwItokOTb4KZPb2Yd8+w7/3L1rCqLkjzcTzZlML5/VFvqUfSIoWLA2gDoaLeq78ddeNNAUeI/x/io/n4KQOq8aLfeDU1Nenp6Xl5ee3atXNzc2vbtm0zXzQtLe3kyZNFRUVGRkbfffedtbX1mTNnevTo0Wi3jh07jhkzpuGWpKQksVj8+eefE0J++OEHBweH/Px8W1vblnwi0Ea9bKhrowQzz7LB+8V7ejMdzTCD5hVdL+FnnmEZmhwdJPA1VUYzGgnJ5u7M33n85JPsLG/6446YUQrKJvtvjuO45cuXW1lZ+fr6RkREBAQEWFhYvPfee7W1tc150Zs3b3p6ehoZGRFCGIbp1KnT9evXn9/tl19+GThw4Ny5c9PT0yVbbty4IT1n39ra2s7O7vbt26/ysUD7mOiQX3szizrS/Q+Lv7jK1WMOTQtVi8kH59mIw+KpnvSpIUpKQan+dtT54YLYTG7GGVaMfztQLtk9woULF65Zs2bMmDGjR49u166d5HqEGzZsKCoq2rt370tftKioyMTERHrX1NS00TpthJD+/fsPHjzYwMDg0KFDnTp1unHjhoODQ6MnmpmZFRYWynwLkUhUUFDg7OwsuUtR1OzZs999992mSqqqqmrm0U2QF5W0+Qhr0qU/mXtB0DmT3ti53sdEu0bbXrnNLxXT0/8VBJnz5weIzXX5qkq5l/ZyhoQk9qImnhMMOVy/qzvbhtGMfzt8tyhfi9pcT0/vpZfRlfGwSCTasGHD0qVLJUOUEqNGjQoODp45c+bKlSvt7Oxe/KLGxsbV1dXSu5WVlcbGxo32efvttyU3Bg4cePfu3T179ixevNjY2LhhZFZUVDTMxYb09fUtLS2PHj0q3WJlZfWCM/15nsc6AEqmqjZ3NyBJg0lsJjf8JPWOB/1FJ0aoNYNtr9DmYo58f5NbfYtd3YUZ50YTosqT+wwISRxI5p5j+x0VJEYwGrESDb5blE/ubS7jG6KsrKympiYqKqrR9rfeeovn+YKCgpe+qKOjY0ZGBsv+b97C/fv3X3yBJzs7u/LyckKIk5NTWlqaZKNIJMrLy3vBExmGcWkAf4vQ0Bhn+tIIwbViPiRefK1YM/oWyne7lO8SLz5byN8YJRznpha/FxiKbAhlIl3oHgfZ1HL8w4EyyPjTt7KysrKyysjIaLQ9IyNDV1fX1dX1pS/atWtXU1PTPXv2EEKSk5OLiooGDRpECLl06dLGjRsl+5w/f15y4/LlywcOHOjduzchZMSIEXfv3v3nn38IIVu3bnV3d/f19X31Dwfazb4tdXiAYI43HZEkXnwRE0qfwROy9hbX+5B4Vgc6vj/TTl/VBT1rSQC99A2632E2qxJZCIrHy7J7924nJ6fjx49Lt1y6dMnPz++7776Tuf/zTp06ZWdn5+rqamlpGR8fL9m4efPmLl26SG47OjoaGxvb2dmZmJisXLlS+sRff/3VzMzM1dXV0dHx4sWLTb1+dnZ2+/btm1kMz/NPnjxp/s4gF+rT5gXVfORRsUdM/Yl8TtW1KFYz27ygmh+YVN/1QH36E7VukB9vs+4x9Q+rVV3HC6nP37n2kHubU/z/n+GXkJDw2WefSQMyPT29rKzM3NxcMlmmqKiobdu2Xl5ely5dambEsixbUFBgZWUlWbz7eaWlpSKR6PmzI+rq6h49emRjY0PTTY7V5OTkhIaGyjxVX6aKigpDQ8Nm7gxyoW5tHp/FzTnHDWxPfRPcalfrbk6bJ+XyU0+xUzyp/7zBCNRiNPRFvrrG/Z7BnRgsUNuFSdXt71wbyL3Nn06WMTQ0dHFxkd5tePvVMAzz4mk1pqampqamz2/X0dF56XwcgJYa5kj3sqUXX2R99tV/15kZ56Z1F3KqZcnHF9k/HvB7ezOack3HJQF0eR0/6C/xkYECxZ3XD1ruaY9Qs6BHqP7Uts0vPuJnnmWNdciGUMazdV3v/gVtfqeMH3ecdTWitnTXsA4xT8iMM+z9J/yhCIEeo+pqnqO2f+etmNzbXO1HRgDkLdiSOj9cMNyR7pEgXnaZFYlVXZCC8YT8dJfrdVA8x5uO66thKUgIoQj5KZRpp0+NO85q5M92UHtNBuHt27cnTZokWVPGz88vKipKOs8TQNMxFJnnQ18bJUgtJ977xH8+aLVrmTyqIcP/ZrencmeHCqZ6auoPX5oiO3oyRTX8Z1cw9xfkT/b/GKdPnw4ODo6NjZUsme3k5JSYmBgaGvrHH38ouT4AxbFtQ/3Wh9kexvznCtf/sPhuWWvrb/yVy7/xp9jXlJwdKnDX8EFgHZrs6yvYkcr/0Xp/tYCqyD5GGBAQwDDMoUOHrK2tJVvKysrGjh1748aN3NxchlH9OD2OEao/DWpzMUd+ust9eY2d4EYvC2SMNHZShrTNK+rJwvNsUi6/s6fGzItpjiuP+YF/iZMHCvzVZl11Dfo7bzWUcYywpKTk+vXrq1atkqYgIcTExOTHH38sKChISUmR49sDqAMBTeb60LfeFJbXEc+Y+vV3NHvN7mP5fMc/xBwhN94UtKYUJIQEWlD/7crgUswgXzKCsK6ujhAiuXZEQ5ItzbwABYDGsdQjW3owxwYLTj7kffeJYzM5jRsqFbHUxxfZiSfZ/4bQm7trcNf2Bd5ypce6UCOPiOs0+ccKqBXZS6xZW1v/+OOPjbavW7dOT0/P09NTKYUBqEYHEyqmL/NjN+bLq1zPg+LzRRqThkm5fJfDOo9E5PabgqEOmjovpjm+DGJMdMgH/2LiDMiHjKtP0DS9bNmyWbNm3blzZ+zYsTY2No8ePUpISDh69Oinn37a/MvzAmiu/nZU+EjBjjQu8hjrb0aWBzKdLNR3jDGzgl/wL3e3jF8TJB7hrqPqchSOpsgvvQVB+8UxGVykS2uOfFAO2VdpmjlzplAo/Pzzzz/88EPJFisrq+++++6DDz5QYm0AqkRTZIoHPd6N3pHKjUhmvUzIimAmSM3iUCQmP9zmVt9iZ3agf+8jqKsWqboiJTESkri+THiiOMCc8tDwCbGgcjKCUCwWnzhxYsCAAdOmTcvMzCwpKTE2NnZxcXnByp8ArZUOTaK96Enu9OYUbkQyG2xBfehPh7ZT/Tcvx5OYDG7xJS7Ykro8QmDfliKE1Km6KmXyN6OWvsFMOMGeHirQwZcTvAYZfz5FRUX9+vXLzMwkhDg7O3fq1MnNzQ0pCNpMlyFzfej7kYJwO2ryKbbLAfFv6SqbWcryZG8657tP/MNtbmsY83sfxl4TLmCrCHN9aOs21OKLOFgIr0VGvJmbm+vr61dVVSm/GgB1pseQ2d50ymjBkgB6cwrn8rv42+vcIyXO42d5suc+57tPvOEOtzaEOTdM0Lt1nR3RUhQh28OY2Ez+cI7GzGkCNSQjCHV1daOjo9euXSs5jwIAGqIpMsyRPjZYkNCfuVfOe8bWD/1b/Fs6p9A1S7Mq+c+vcu4x4s0p3LpuzJmhgv52Wh2BUma65JdezNTT4vxqZCG8ItmTZZycnGJiYjw9PQcMGGBvby8QPN1t0aJFyqoNQK0FmFPbwpgfxcz+B9yu+9ysc+xwR3q0M93TmpLXBYOqxGRfJrcjlbtZyr/lQsf2VevJq6rS3Zqa1YEZf5xNHiRg0DzQcrKXWLO2ti4sLJT5BDW5bBOWWFN/2tbmhSLyWzp3IIu7+JgPMKf62lJ9bemuVpSwhYfX6zhy8RF/8iF/uoD7p4jvYU29404PcaB1m7Gyoba1uRTHk/BEcYQ9vaijsmczaG2bq5ACL8zbUEFBgRzfA0AbtNMn833p+b60SEzOFPJH8rgF/7J3ynhXQ8rDmPIwJh7GlLsRpS8gxjqEpghNiLEO9biGz68m2VV8XhXJq+JvlfKXH/NeJlQPa2pmB3pvb1rjrpqkEpLLUwQfEA9sT6nPMqSgKWQE4ePHj1NTU62trZ2cnDBZFKCl9AWknx3Vz44hhFSLSdoTPrWcTy0nx/P5LSlcHUfK6wjHE5YnT+p5M13Krg1xMKDs2hIPY2qoIx1iRRm2xqXRFM3BgFrZmZl4gr0wAmdTQMs8E4RisXjq1Km7d++WjH96enrGx8d7eHioqDYAjddGQDqaUR3RR1GKie50Qjb/n8vsimDVXyEHNMgzP5w2bNiwa9eu3r17f/311zNmzMjMzJw8ebKqKgMAaKmfQpldafypArWYygCa4pke4eHDh8PCwo4ePSq56+fnN3v27PLycmNjY1XUBgDQMhZ65KdQevJJ9tooAUaYoZme6RE+ePCgf//+0rsRERGEEMkSMwAAGmGYI93ThvroApabgeZ6JghFIlGbNm2kdyUXmqiurlZ2UQAAr2FtCJOUyyflYoAUmqXxrNHc3NzLly9LbhcXFxNC7t27p6v7dAZ3p06dlFYcAMArMBKSrT2Yd06yN98UGLf+y1LB63rmhHonJ6esrKwXPwEn1EMzoc2VD23e0IwzLMuTLT0UO4MUba58ij2h/uuvv66srJTjqwMAqMr3XRn/feK/cvkIe5y+Ai/yTBC+/fbbqqoDAEC+2grI5h7MtNPsDcwghRfCAgwA0Gr1taX62FCLMIMUXghBCACt2equzMFsPjlPLSY3gHpCEAJAa2asQzZ2Z949w1bWq7oUUFcIQgBo5Qa1p7q3o5ZcwgApyIYgBIDWb20IE5fJnyvEACnIgCAEgNbPTJds7E5POcXWoFsIz0EQAoBWGOpA+5pRX11DEkJjCEIA0BbrujGbU7grjzFACs9AEAKAtrDWJ98GM1NOsfWcqksBdYIgBAAt8o4HbdOGrL6JJISnEIQAoF02dWdW3WTvlmGAFP4HQQgA2sXBgFr2BvPuGZZDFAIhBEEIAFpotjfN8eSnuxggBUIQhACghWiKbAtjll9hMyrQKwQEIQBoJQ9j6kM/ZvppFkkICEIA0FIf+NEV9WRHKgZItR2CEAC0lIAmW3swiy6wuVXoFmo1BCEAaC8/M2qWNz3zLNZd02oIQgDQaosDmKxK8ls6Bki1F4IQALSaDk229mAW/Ms+qlF1KaAiCEIA0HbBltREd3rOOQyQaimB4l766NGjCQkJ5ubm06dPt7a2bvRoXl5eQkJCSkqKhYVFVFSUq6urZPvPP/9cXFws1V6eTgAAIABJREFUue3g4BAVFaW4CgEAJD7vxATtF8dkcJEu6B5oHUX9k8fGxkZFRbm5ueXm5nbt2rWysrLRDjNmzDh37pyLi0tRUZG/v//ly5cl21evXn3r1q3S0tLS0tKKigoFlQcA0JAuQ3b2ZOb/wxaKVF0KKJ2ieoQrVqxYtWrVxIkTCSE9evT45Zdf3n333YY7xMXF6erqSm6XlJTs3LmzU6dOkrszZswIDQ1VUGEAADIFWlBTPOnoM+yBfoyqawGlUkiPsKKi4urVq/369ZPcDQ8PP336dKN9pClICKmurjYyMpLe3bdv37fffvvXX38pojYAgKb8J5B5UMHvxQxSLaOQHuHDhw8pirKwsJDctbKyOnnyZFM7Jycnnzhx4scff5TcDQ4O1tHRKSkpmT59evfu3ffu3SvzWbW1taWlpVOnTpVuiYiIGDJkSFPvIhKJGAa/8pQKba58aPPXt6kLNfQY09lYbNumWfujzZWvRW2uo6MjELwk6RQShDo6OjzPsywrFAoJIfX19Q37fw1du3ZtwoQJe/bssbOzk2zZuXOn5MaCBQvc3NwuXLjQuXPn558oFAp1dHSCg4OlW7y9vZt6F0JIXV3dCx4FRUCbKx/a/PUFW5NZ3vx7l4UHwps1YIY2V74WtTlNv/zfUSFBaG1tTdN0Xl6eZC5obm6ura3t87vduHFj4MCBmzZtGjRokMwXcXJyyszMlBmENE23bdt2xowZzSyJYRj8alMytLnyoc3lYukbpPMB8Z4MapL7y79D0ebKJ/c2V8gxQj09vYiIiN9//50QUltbe+DAgWHDhhFCqqqqkpOTxWIxISQlJWXQoEFr1qwZPny49Im1tbUs+79TeW7fvn3//n1fX19FVAgA0BQhTXb2ZBaeZzNxkSbtoKhZo5999tnAgQNv3ryZlpbm6Og4ePBgQkh2dnb//v1LS0tNTEymTZtWUVGxatWqVatWEULCwsJWr1598+bNkSNHBgcH8zx/7NixRYsW+fj4KKhCAICm+JtRH3VkJpxgTw4RMJSqqwEFo3heUT95Hj16dPr0aVNT07CwMEk3tqam5vbt2wEBAQzDpKSkVFVVSXc2MTFxdXXlef7u3bupqak0TXfs2NHR0bGpF8/JyQkNDc3Ozm5mMRUVFYaGhq/5iaBF0ObKhzaXI44n/Q+Le9vSSwJeNHKGNlc+ube5AoNQoRCE6g9trnxoc/nKq+I77Rcf6CfoYtVkrxBtrnxyb3MsJgQAIJtdW2pDKDPuBFtZr+pSQJEQhAAATRrlRIe2oz48j/W4WzMEIQDAi6zvxhx7yMdnYbmZVgtBCADwIgZCsj2MmXGWza/WyBkV8FIIQgCAlwhtR832ZqKOsWJ0C1sjBCEAwMt90pFuIyD/uYKDha0QghAA4OVoivzSW7A3nU/IRq+wtUEQAgA0i5ku+bU3M/00m1WJg4WtCoIQAKC5ulpR7/kybx1j69EtbEUQhAAALbCoI22lTy27jIOFrQeCEACgBShCtvZg9qbzB3BmYWuBIAQAaBkLPbK/HxN9hr1ZgoOFrQGCEACgxd4wp9Z3Y4Yns49rcZUmjYcgBAB4FaOd6bEu1PizwjoMkWo4BCEAwCv6KogxFpIP/sXEGc2GIAQAeEU0RbaG1B9/yG9KQa9QgyEIAQBenYGA/zOc+c9l9lQBJs5oKgQhAMBrcTemfuktiDwqvlWKLNRICEIAgNfV15Za25UZmITV1zSSQNUFAAC0Bm+50iW1pN9h9swQgZW+qquBlkCPEABAPmZ506OcqKF/i6vEqi4FWgJBCAAgNyuCGV9TauxRMS7hq0EQhAAAckMRsrE7wxMy7TTL4XChhkAQAgDIk5AmsX0FOVX8lFMsiyzUBAhCAAA5ayMgB/sL8qv5ccdZjJGqPwQhAID86QtIQn9BlZgfdwJZqO4QhAAACqHLkH3hghqWRB3HFe3VGoIQAEBRdGgS04epYfmo42wtluZWVwhCAAAFkvQLaUL6HRYX16q6GpAFQQgAoFg6NPm9L9PHlupyQHyvHBNJ1Q6CEABA4ShClgcyiwPongfFuE6FukEQAgAoyRQPek8vwZij4r3pmDyjRhCEAADKE25HJQ8UfHKRW3IJp1WoCwQhAIBS+ZtRl0YIrhbz3Q+KH1RgmFT1EIQAAMpmqUcORQiiXOhuCeLEHGShiiEIAQBUgCJkvi/9Z7hgzjl2/j84416VEIQAACrTxYq6MFyQ/oTvFi++VoyuoWogCAEAVMlCjxyMEHzUkR6YJJ7/D1tZr+qCtA+CEABA9cY407dHC2tY0vEP8V+56BoqFYIQAEAtmOmSTd2Z9aHMzLPsuOMsJpQqDYIQAECNDLCnbr0pcDcmQfvF8/9hi0SqLkgLIAgBANRLGwFZHsjcixTqC4h3XP3HF9myOlXX1KohCAEA1JG5LvkmmLkyUlAkIp6x9Z9eZh9Wq7qmVgpBCACgvhwMqG1hzJkhgrJa4ruvfuIJ9spjHDuUMwQhAIC6czemfuzGpI8V+ptTo46wYQfFe+5zONFCXhCEAACawUSHfOhHp48VzPehf8/g7H+tH3ecPZTDY1Wa1yRQdQEAANACDEXedKbfdKYf15CYDG7FNXbySX6kEz2oPdXXljYQqro+DYQgBADQSBZ6ZJY3PcubflDB/5nFr7/DTTjBdrakBrSnB9hTPqYUpeoKNYUCgzAlJeXmzZvu7u4BAQEydygpKTl58qSBgUGvXr2Ewqc/Yy5dupSZmRkYGOjq6qq48gAAWgcnQ2qBL7XAl64Sk2P53OEcftQR7pGID2lHhVjR3dpRnS0pQ/QUm6aoINy0adOyZcsiIiKOHz8eHR396aefNtrhzp07vXr1CgsLy8vLo2n62LFjurq6hJAPP/xw37593bt3nz179po1a8aNG6egCgEAWpm2AjLUgR7qQAghRSLyTxF3rpD/7Ap3tZi3bUP5m1F+ZpSvKfE3o5wNKQYdxv9H8bz8Z+KKRCJ7e/v4+PjQ0ND79+/7+/s/ePDAysqq4T5RUVHt27f/7rvvxGJx586dFyxYMGHChAcPHnh7e6elpdnZ2SUlJU2bNu3BgwcCgYy0zsnJCQ0Nzc7ObmZJFRUVhoaGcvhs0Gxoc+VDmyufRrS5mCP3yvlbpfyNEv5WKblZwudX844GlIshcTGiXA0pZ0Ni04ayb0va6VNCtZ9DKfc2V0iP8MyZM23atAkNDSWEuLm5+fr6JiUlTZw4seE+8fHxp06dIoQIBII333zzwIEDEyZMOHjwYEhIiJ2dHSGkf//+1dXVly9f7tKliyKKBADQEgKa+JhSPqbUWJf/ballSWYFn15BMp7wGRX8qQJSIOJyq8gjEW+mS9rpU1b6xEKPMtcl5nrEXJcy1SVGQmKkQxkJ/6+9Ow9r4ugfAD45CYchQIhJICBvFUVURBRROQUl9aoYj1qw9QRE2/fRXqi1+uojXq0X+mpFW+WtitZqFVBQ8SjUI6jIISBCg3LJEQhEIPf+/pjfu2+eEGJQATHz+Wt3mN2dnSz5ZmdnZwDDDNApBDMSeG+aW7slEFZVVfF4PHyVx+NVVlZqZ2hsbGxra8PzODo6pqSk6GxIJBK5XG5lZaXeQKjRaFpbWw8dOoSn+Pr6urm5dVYktVqtVqvf4JyQLkN13vNQnfe8PlrnZAAG9QOD+gHAxdMIAAANRqhtB3UyUNsOxHJMLAONClAswSQK0KIALUqsRQFalIRmBabQAKkSWJCBGQkwqAQCAP0oGInw/ykAAGsqAd5bWlIA9b93mSQCoFN1C2NNAV1qpnWhgTmuxtY5kUgkvKrbULcEQqVSSSKR/ncMMlmpVOpkAADgeSgUikKhMGZD7T0oFIr79+/jKUwmc+DAgQaK1NmukG6C6rznoTrvee9fnTMpgEkBQ41remxTAYWG0KTAAABSJUGNgTYVkKsBAKBFCeD7jS+VAH/TUaUBL1W6YalJ1rVHdFYajfF1TqFQtMOKXt0SCDkcTn19Pb5aV1cXHBysncHe3p5MJtfX19vZ2cEMXC4XblhSUqK9IYfD0XsIMzMzGxubI0eOGFkkpVJJo9G6eiLIm0B13vNQnfc8E69zeObsnj2oVKp4u3XeLU9Fvb29y8vLYU8WqVQqFAr9/PwAAGq1Wi6XAwCIROKECROuXbsG81+7dg1m8PPzy8rKkslkAIDCwkKpVDpq1KjuKCGCIAiCQN1yR8hisZYsWSIQCCIjI5OSkvh8Pnx6d/jw4Z9++unRo0cAgNjY2E8++USj0YhEogcPHhw7dgwA4O3tPWrUqFmzZoWFhcXHx69cufLd746FIAiC9Gnd1U923759MTEx+fn5AoHg5MmTMNHX1/ebb76By3w+/+LFi+Xl5dbW1kKhkMlkwvSLFy+Ghobm5+d/++23cXFxb6UwL1++zMzMfCu7QozU2toKewUjPaatre3WrVu9XQrT0t7efvPmzd4uhWmRyWQ3btx4yzvF+qbnz5/zeDwjM2dmZo4bN65by4PouH37tre3d2+XwrTcu3fPy8urt0thWrKzsz09PXu7FKbl4cOHI0aMeLv7fOffnEQQBEGQ7oQCIYIgCGLSUCBEEARBTFq3jDXaA8rKyoYPHw5HcXul5ubmp0+fjh49urtLheCam5tLSkrGjBnT2wUxIS0tLcXFxd7e3r1dEBMilUqLiopQnfckqVRaWFho/NCbYWFhMTExhvP01UCo0WiOHz+uPZCbASqVqqamxsjMyFuB6rznqdXqqqoqJyen3i6ICUF13vM0Gk1FRYWzs7OR+V1cXF45o19fDYQIgiAI8lagZ4QIgiCISUOBEEEQBDFpKBAiCIIgJg0FQgRBEMSkdcug2z2vpqYmIyODTqeHhoaamZl1zCCXy9PS0qRSaUhICJv9vzlDSktLs7KyeDzexIkTXzl5I6JNIpGkpaWRSCQ+n693bPTCwsJHjx6Zm5v7+/vD+bYAAPn5+bW1tXCZQqEEBAT0XIn7vtra2qtXr/br14/P53e8zuvq6vLy8vBVT09PvNrLysoyMzMdHByCg4OJRPTztwuam5vT0tIIBAKfz6fT6Tp/FYlEZWVl2ikBAQEUCqWgoODFixcwhUwmBwYG9kxp3w/V1dVPnjxxdXV1cHDQm0EsFqenp5uZmfH5fEtLSzz93r17hYWFHh4eXZ226H3oNZqTkxMSEjJ9+nSRSCSTyW7duqUzVVV7e7u/v7+lpeWAAQNSUlIyMjI8PDwAABcuXFiyZMmsWbOEQuHgwYNPnz7dS2fQ91RWVvr4+Pj4+CiVyoKCgjt37rBYLO0M33//fWJi4vjx46VS6V9//ZWamjpu3DgAwLx58/Ly8hwdHQEAdDr9999/750T6IPy8vKCgoKmTZv2/PlzqVSamZlpbm6uneH8+fNLlizx8vKCq9u2bYPLqampn332WVhY2IMHD5ydnc+fP98Lpe+bampqxo4dO3r0aAzDcnJy7t69q/0zGgBw4sQJOHMOAKC6urq6urq2tpZKpS5YsCA7Oxu+PmRhYXHhwoUeL3tfNWHChNzcXI1Gs3fv3mXLlnXMUFZWNmHChKCgIIlEUl5efufOHQaDAf77ncPn81NSUlavXr169eouHPXtDl3aK2bOnLlhwwYMw1Qq1ahRoxITE3Uy/Pzzz6NHj1apVBiGrV+/XiAQwPShQ4eeOHECw7Dm5mZ7e/t79+71aLn7stWrVy9YsAAuh4WFwfrXJhKJYIVjGPbVV199+OGHcHnu3LkHDx7sqWK+V+bMmbN27VoMw9Rqtbe399GjR3UynDt3LiAgoOOGHh4ex44dwzBMKpVyOJzMzMzuL+x7IjY29uOPP4bLc+fOhfXfmQULFsTExMDliIiIvXv3dnv53kfwq8PX1/fw4cN6M0RGRkZHR2MYptFoJk2atHPnTgzDamtrzc3Nnz59imFYTk4OnU6XSqXGH7TPN5JoNJrU1FSBQAAAIJFIM2fOTElJ0cmTkpIyc+ZMEokEAJg9e3ZqaiqGYWVlZSUlJWFhYQAAOp0+efLkjhsinUlJSZk9ezZcFggEHatuwIABsMIBABwOR6FQ4H8qLy9PS0vTaVBCXgm/zolEYlhYmN7LtbW1NT09PTs7W6lUwpTnz5/n5eXBDa2srEJDQ9F1brzk5GRYdQCA2bNnG6i65ubm33//ffHixXjKs2fPLl++XFpa2u2lfL9of3XohX/5EAgE/EO5evWqm5vbwIEDAQAjR460t7fv0pRkfT4Q1tfXK5VK2NQGAHBwcKiqqtLJU1VVhbc1Ozg4yGQysVhcVVVla2uLNy7p3RDpTHV1tXaVGqg6iUQSHx+/dOlSuEqlUu/cuRMfH+/l5bV48WKs77fM94ympqa2tjbD1zkAQCaTHThwIDw83MPDo7y8HABQXV3NYDCsrKwMb4jopfPVYaDqTp069cEHH+Dt0hQKJTs7+8CBA2PGjPn00081Gk1PFNcEqNXq2trajh9KVVUV/t8Bun6d9/nOMmq1GgCA93MhkUgqlapjHryDAPytoVKp1Gq1du8YvRsindGuPQNV197eLhAIQkJCPv74Y5hy7Ngx+BHU1NR4enqeO3cO/8WNGGDMdT5jxgzYwqHRaCIiIr7++uvffvsNXedvwsjrHADw888/az/QSkhIgNd5XV2dp6fnmTNn8H8B5E1oNBoMwzp+KDrXOZlM7tJ13ufvCFksFpFIrK+vh6u1tbVcLlcnD4fDqaurwzOQyWQWi8XhcJqamvDKqq2t5XA4PVbsvo7NZhuucwCAXC6fNWsWm80+dOgQnqjdXhoUFJSTk9MDpX0P2NnZUalUw3WO1y2RSJw3b96jR48AAGw2WyKR4E3T6DrvEg6H88rrHABQUFCQm5s7f/58PAX/LFgsVnBwMLrO3xYKhWJnZ9fxQ9H+kgcGPyy9+nwgJJPJfn5+6enpcPXKlSuwp7JGoxGLxbBFIjAw8MqVK3gGf39/IpE4cOBAFot148YNAIBKpbp+/XpQUFDvnEMfFBgY2LHOAQCNjY3wt4VSqZw7d66FhcXx48f1tvirVKq8vDw0WrGRCASCv79/xzrHMAy/zrU9fPgQdll0cXFxcnK6du0aAECtVmdkZKDr3HhBQUHaXx0dr3PoyJEjYWFhTCaz4x7UanVubi66zt+QQqGQSCRwWfv7PD09HX4o/v7+OTk5YrEYAFBRUVFWVjZ+/PguHOBNu/i8A9LS0hgMxtatW5cuXerk5NTU1IRh2LNnzwAAlZWVGIY1NjY6OjouW7Zs69atDAbj6tWrcMP9+/fzeLwff/xxxowZY8eOVavVvXkafUp+fj6dTv/uu++++eYbGxub0tJSmE4mk//8808Mw2JjYykUyqJFiyIjIyMjI2F3O5lMNn78+O+//z4uLm7cuHHu7u4vX77szdPoUzIyMqytrePi4qKiohwcHMRiMYZh1dXVAACRSIRh2IoVK7744oudO3cuXLjQysrqxo0bcMOffvrJwcHhxx9/DAsLGzVqFN6bF3mloqIia2vrtWvXrlmzhsFgFBcXw3Rzc/OMjAy4LJfL7e3tr1y5gm+lVqt9fHzWr1+/detWX1/fIUOGtLS09ELp+6ZDhw5FRkay2Wx/f//IyMjc3FwMw06fPs3j8WAGoVBIp9M3bty4atUqJpNZUVEB0yMiInx8fPbs2ePp6fn555936aCkjRs3vsW43SsGDhwYHBxcUFDA5XL3798Pf5eRyeQBAwZ4e3tTKBRzc/NPPvmkoqKitbV1y5Ytvr6+cENvb++hQ4cWFRWNHDly9+7det/ER/RisVizZs0qKiqi0Wjx8fH4LCcODg4+Pj5WVlYUCmXs2LGOjo5cLpfL5fJ4vOHDh5NIJDs7u7q6Orlczufz9+3bp/MmHGKAi4vLpEmTCgoK+vfvf+DAAXt7ewAAmUx2dnYeO3YslUqF0bGhoWHQoEHx8fHwZVkAgJeX14gRIwoLC4cPH75nzx6dt2wRA5hMpkAgKCoqolAo+/btc3V1helcLtfHxweOIyEWix0dHQUCAf6MikAg2Nvbw+t80qRJ8fHx2i99I4Y1NTWZm5sHBQWNGDGCy+W6u7szGAwLCws3N7fhw4cDABwcHKZPn15YWNivX78DBw7gc73NmDGDRqOJRCKBQPDll192aYCU9+GFegRBEAR5bX3+GSGCIAiCvAkUCBEEQRCThgIhgiAIYtJQIEQQBEFMGgqECIIgiElDgRBBEAQxaSgQIibk+PHjBQUFPXa4urq6s2fPwmFCX0NRUdHhw4flcnlnGY4fP/748ePXLd07p6Gh4fjx4/i8zV3yxx9/VFRUvPUiISYCBULkPSGXy/9hkFqtXrhwYWpqao8V6csvvzx58qThOWUM+PPPP6OiolpbW/X+FcMw7dN58ODB0aNHX7Og74bS0tKFCxcWFxfD1cuXLycnJxu5bWpqanR0dLcVDXnP9fnZJxAEIpFI2rPB7d+/n0QiLV++HE8hEombN2/29/fvmfLk5uaePHny/v373bR/AoGwefNmPz8/uJqcnLxt27YlS5Z00+F6AI/H27x5s4uLC1yNj49vbW2dPn26MduuXbt24MCBN27cQCOpIq8BjSyDvJ9GjBhBpVKNjEP19fV0Oh0fYw+O2A7HMOtIIpHI5fL+/fsb3mdkZKRQKIRTQGiTSqUajcba2rqzDWUymUQiYbFYCQkJ0dHRYrHY1tZWLBZTKBQ6nd7ZVhs3bty2bZtMJtP71/b29ubm5v79++PjTtXV1dnY2FAolM522NTUBACwsbExcI5vor6+3sbGhkzu9Lf4lClTWltbO06vqlarGxoa+vXrZ2FhoZ0+ceJEW1vbs2fPdktxkfcaahpFTAiHwzlw4ABc3rRpk7Ozc1ZWlqurK4vFotPpcGTww4cP9+/fn8Vi2dvbnz9/XnvzlJSUYcOG2djYsNlsHo938uTJzg6kUChOnToFZwfERUdHs9lsOp3OYDC4XO62bdvwn6EvXrywtbU9cuTIp59+SqfTORwOPq3oo0ePPDw8mEymjY3N9OnTGxoatE/n3//+NwDg22+/3b59u1wut7W1tbW1hSFcJBLZ2tomJibOnz8f7rOxsbG4uHjy5Mk0Gq1///4WFhZeXl5ZWVn4DuPj421tbYVC4ciRI21tbT/66KPY2NgBAwboxNegoKAZM2Z0PGtYAO1Zt3RSVq1a5enpmZGR8cEHH7BYLCsrq0WLFimVSvjXhw8fcjicu3fvAgBCQkKuXbt2+/ZteEZ8Ph8AIBaL58yZQ6PR2Gy2paWli4tLRkYGfiyBQHDhwgUYvxGkS1DTKGJCXrx48fLlS7gsk8levHixbNmy9evXDxkyJDExcevWrfX19Tk5OYmJiQwGY/369QsWLHj27JmdnR0AIDk5eebMmeHh4QkJCTQa7ZdffomIiLC2tp46dWrHAwmFwpcvX+pMBNPe3r5///7BgwcrlcrTp0+vWbOGxWLB5lyNRtPU1LRu3bqJEydeunRJo9EwGAy4VURERGxsbEhIyP3791euXDlnzpzr168TCAQMw/DTWbx4cU1NTVJS0pkzZ8B/5+9Vq9VNTU1ff/01n89PS0tTKBQWFhZisdjDwyM2NpbL5VZWVv7rX/+aOnVqSUkJvMGVyWRNTU3z5s1bsWLFwYMH29vb2Wz2jh07zp49GxERAcuTl5d38+bNX3/9VW8NNzU16URN7ZS2traysrIVK1Zs2rRp8ODBycnJmzZt8vHxiYqKAgAoFIoXL17AzkFxcXErV66UyWS7du0CAMAb6H/+85937949f/68m5ubRCK5c+cOPuE2AMDHx0elUmVmZuoN0ghiyNuaOwNB3inDhw/38vLSSQQAwPswDMPWrFkDAEhNTYWrSqXS3t7ewsKipqYGpsAOmadOnYKrrq6uwcHBcIJsKCgoKDAwUO/R4Y1aeXm5gRJOmzZt4sSJcBne/40bN047A7yRWrdunU5KVlYWhmFwDsLt27fDP23YsMHMzEx786dPnwIAgoODDZShsbGRTCYfOXIEru7YsQMAEB8fr53H19fX398fX12+fLmdnV17e3vHvcGAt3v37s5SIiMjAQC3b9/GM4wePTo0NBQu37lzBwBw8+ZNuPrhhx9qHxfDsMGDB8fExHR2LlKpFACwadMmA+eLIHqhO0LEdJFIpJCQELhMJpPhZFJsNhumuLq6EggE2Cm/vLy8pKQkODhYuy3O0dHx0qVLevcMZ9C2tbXVTlSpVOfPn8/Pz4dvCFRUVOBzjUJ6b2UEAoH2cnR09MOHDydMmGDkOXbcZ11d3enTp0UiEeyPamZmVlpaqp3ho48+0l5dvnx5eHj448eP3d3dW1tbT5w4sXTp0teey4nBYIwbNw5fdXd3z87ONnLbUaNGJSYmWlhYzJkzx8vLS6c7rpWVlZmZmfY05QhiJBQIEdNlaWlJpVLxVSqVqj1vHJlMJhKJsKUOhq7//Oc/SUlJOjtRKBTaO8G3BQBoT2Le2Njo5+dXXV0dGhrKZrNpNJq5uTmcPhqHx+DOEplMJpVKraysNP4cdfZ56dKl2bNnOzk5+fn52djYEIlEEonU0tKinUenH9Ds2bNXr1599OjRXbt2nTp1qqWl5U36pur0vjEzM1MoFEZuCyf2O3LkyA8//MBkMiMiIjZv3mxlZQX/imGYWq3u+FkgyCuhQIggrwa7a27fvj0mJsaY/DCWiMVi/Hs/KSmpuLi4pKQEn8R48eLF+DtzkN6pRBsaGjgcDlyWSCQKhQJfNYbOPuPi4kaPHn3jxg14O6XRaPbt26ezifaDNwAAlUpdtGhRQkJCXFxcQkJCYGDg0KFD9R6LRCIRCATt8K8TYt+QnZ1dQkLCwYMH79+/f+7cuV27dsnlctgKDQCQSCQqleqVvXkRpCPUaxRBXm3w4MEcDue3337DjHvdaMyYMQCA/Px8PEUkEtnY2OBRUCaTabeyGpCtyWfUAAADZklEQVSWlqazPGzYsI7Z+vXrp1QqXzmKjUgk8vT0xBsVr1+/3t7e/soyLF++XCKRbNiwQSgUwud8epHJZDabXVZWhqd0fPnBeFZWVnrLRiaTfXx8duzYERoaCruYQnl5eQAAb2/v1z4iYrJQIESQVyMSiVu2bLl58ya8jWtvbxeJRL/++uvWrVv15h82bBibzf7rr7/wlJEjR4rF4r1798pksr///js8PLyxsdGYQ+/Zs+fixYttbW23bt366quvRowYofedcXd3d41Gs3379rt37z58+LCzvY0cOfLMmTP379+Xy+UZGRlRUVHGPO1zcnLi8/k7d+5kMpmzZs0ykHPy5MknT568fPlyfX19cnLyd999Z8w56uXu7p6Xl3fq1Kns7OwnT54AAJYtW3bt2rWGhgalUpmVlSUUCr28vPD8t2/ftrKyGjt27GsfETFZKBAiiFEWLVr0yy+/pKenu7m5WVhY/OMf/1i1alVnL6QTicRFixYlJSXht2jz58//7LPPVq1aZW5uPmjQIBqNZuDWStvu3buXLFliaWkZGBjIZDL/+OMPvWO2hYaGrl69ev/+/ePHjzdwV7R7924bG5sxY8bQaLSwsLCNGzca2ZYYFRWFYdjChQvxYQf02rJli6ur65QpU1gsVlRUVMd2V+N98cUX06ZN+/zzz729veFLJkVFRXw+397enkqlBgQE+Pr6/vDDD3j+pKSk8PBwc3Pz1z4iYrLQyDLI+wm+XaDzuEutVhOJRL2P4ozf7ZMnT6RSKYvF4vF4BsYRraysdHV1PXv27JQpU/DEmpqayspKJyenLj3Kksvljx8/ptFobm5u2oV/vdNRqVRlZWWtra1ubm7Gh43Dhw9HR0cXFxe7uroazolhmEgkam5udnd372rXFbVabXho1ubm5oqKCoVC4ezsDN/vhIRC4fjx43Nzc93d3bt0RAQBKBAiSPdZt27dlStXhELhm4Ted0Fra6unp+eQIUMuXrzY22XRb+rUqTweT3tQGwQxHgqECNJd5HJ5VVWVs7Pza09A8S4ICAgoKChQKpX37t1zc3Pr7eLo9/fff7PZbJ3RRxHESCgQIghiyLFjx8hk8sSJE7lcbm+XBUG6BQqECIIgiElDvUYRBEEQk4YCIYIgCGLSUCBEEARBTNr/AcKodw4xW8RRAAAAAElFTkSuQmCC", "image/svg+xml": [ "\n", "\n" ], "text/html": [ "\n", "\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Bootstrap a KDE of the pre-eruptive (or pre-deposition) zircon distribution\n", "# shape from individual sample datafiles using a KDE of stacked sample data\n", "BootstrappedDistribution = BootstrapCrystDistributionKDE(smpl)\n", "h = plot(range(0,1,length=length(BootstrappedDistribution)), BootstrappedDistribution,\n", " label=\"Bootstrapped distribution\", xlabel=\"Time (arbitrary units)\", ylabel=\"Probability Density\", \n", " fg_color_legend=:white)\n", "savefig(h, joinpath(smpl.Path,\"BootstrappedDistribution.pdf\"))\n", "display(h)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Run eruption/deposition age model" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Estimating eruption/deposition age distributions...\n", "1: KJ08-157\n", "2: KJ04-75\n", "3: KJ09-66\n", "4: KJ04-72\n", "5: KJ04-70\n" ] } ], "source": [ "# Number of steps to run in distribution MCMC\n", "distSteps = 10^6\n", "distBurnin = distSteps÷100\n", "\n", "# Choose the form of the prior closure/crystallization distribution to use\n", "dist = BootstrappedDistribution\n", "## You might alternatively consider:\n", "# dist = UniformDistribution # A reasonable default\n", "# dist = MeltsVolcanicZirconDistribution # A single magmatic pulse, truncated by eruption\n", "# dist = ExponentialDistribution # Applicable for survivorship processes, potentially including inheritance/dispersion in Ar-Ar dates\n", "\n", "# Run MCMC to estimate saturation and eruption/deposition age distributions\n", "smpl = tMinDistMetropolis(smpl,distSteps,distBurnin,dist)\n", "results = vcat([\"Sample\" \"Age\" \"2.5% CI\" \"97.5% CI\" \"sigma\"], hcat(collect(smpl.Name),smpl.Age,smpl.Age_025CI,smpl.Age_975CI,smpl.Age_sigma))\n", "writedlm(smpl.Path*\"results.csv\", results, ',');" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "AgeDepthModel.pdf\n", "AgeDepthModel.svg\n", "AgeDepthModelHiatus.pdf\n", "BootstrappedDistribution.pdf\n", "DepositionRateModelCI.pdf\n", "KJ04-70.csv\n", "KJ04-70_distribution.pdf\n", "KJ04-70_rankorder.pdf\n", "KJ04-70_rankorder.svg\n", "KJ04-72.csv\n", "KJ04-72_distribution.pdf\n", "KJ04-72_rankorder.pdf\n", "KJ04-72_rankorder.svg\n", "KJ04-75.csv\n", "KJ04-75_distribution.pdf\n", "KJ04-75_rankorder.pdf\n", "KJ04-75_rankorder.svg\n", "KJ08-157.csv\n", "KJ08-157_distribution.pdf\n", "KJ08-157_rankorder.pdf\n", "KJ08-157_rankorder.svg\n", "KJ09-66.csv\n", "KJ09-66_distribution.pdf\n", "KJ09-66_rankorder.pdf\n", "KJ09-66_rankorder.svg\n", "age.csv\n", "age_025CI.csv\n", "age_025CI_hiatus.csv\n", "age_975CI.csv\n", "age_975CI_hiatus.csv\n", "age_hiatus.csv\n", "agedist.csv\n", "agedist_hiatus.csv\n", "distresults.csv\n", "height.csv\n", "height_hiatus.csv\n", "lldist.csv\n", "lldist.pdf\n", "lldist_hiatus.csv\n", "results.csv\n" ] }, { "data": { "text/plain": [ "6×5 Matrix{Any}:\n", " \"Sample\" \"Age\" \"2.5% CI\" \"97.5% CI\" \"sigma\"\n", " \"KJ08-157\" 66.0703 66.034 66.0935 0.0151394\n", " \"KJ04-75\" 65.9319 65.8859 65.9706 0.021838\n", " \"KJ09-66\" 65.9368 65.8931 65.9775 0.0213176\n", " \"KJ04-72\" 65.9552 65.9224 65.9769 0.0137651\n", " \"KJ04-70\" 65.8309 65.7533 65.8947 0.0364563" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Let's see what that did\n", "run(`ls $(smpl.Path)`); sleep(0.5)\n", "results = readdlm(smpl.Path*\"results.csv\",',')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To see what one of these plots looks like, try pasting this into a markdown cell like the one below\n", "```\n", "![](\"MyData/KJ04-75_rankorder.svg\")
\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For each sample, the eruption/deposition age distribution is inherently asymmetric, because of the one-sided relationship between mineral closure and eruption/deposition. For example \n",
"(KJ04-70):\n",
"\n",
"![](\"MyData/KJ04-70_distribution.svg\")
\n",
"\n",
"(if no figure appears, double-click to enter this markdown cell and re-evaluate (`shift`-`enter`) after running the model above"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"***\n",
"## Getting your data out\n",
"As before, we can use the unix command `ls` to see all the files we have written. Actually getting them out of here may be harder though."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"AgeDepthModel.pdf\n",
"AgeDepthModel.svg\n",
"AgeDepthModelHiatus.pdf\n",
"BootstrappedDistribution.pdf\n",
"DepositionRateModelCI.pdf\n",
"KJ04-70.csv\n",
"KJ04-70_distribution.pdf\n",
"KJ04-70_rankorder.pdf\n",
"KJ04-70_rankorder.svg\n",
"KJ04-72.csv\n",
"KJ04-72_distribution.pdf\n",
"KJ04-72_rankorder.pdf\n",
"KJ04-72_rankorder.svg\n",
"KJ04-75.csv\n",
"KJ04-75_distribution.pdf\n",
"KJ04-75_rankorder.pdf\n",
"KJ04-75_rankorder.svg\n",
"KJ08-157.csv\n",
"KJ08-157_distribution.pdf\n",
"KJ08-157_rankorder.pdf\n",
"KJ08-157_rankorder.svg\n",
"KJ09-66.csv\n",
"KJ09-66_distribution.pdf\n",
"KJ09-66_rankorder.pdf\n",
"KJ09-66_rankorder.svg\n",
"age.csv\n",
"age_025CI.csv\n",
"age_025CI_hiatus.csv\n",
"age_975CI.csv\n",
"age_975CI_hiatus.csv\n",
"age_hiatus.csv\n",
"agedist.csv\n",
"agedist_hiatus.csv\n",
"distresults.csv\n",
"height.csv\n",
"height_hiatus.csv\n",
"lldist.csv\n",
"lldist.pdf\n",
"lldist_hiatus.csv\n",
"results.csv\n"
]
}
],
"source": [
";ls MyData"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We could use the trick we learned before to view the SVG files in markdown cells, which you should then be able to right click and download as real vector graphics. e.g. pasting something like\n",
"```\n",
"![](\"AgeDepthModel.svg\")
\n",
"```\n",
"in a markdown cell such as this one"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Meanwhile, for the csv files we could try something like `; cat agedist.csv`, but agedist is probably too big to print. Let's try using ffsend instead, which should give you a download link. In fact, while we're at it, we might as well archive and zip the entire directory!"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [],
"source": [
"# Make gzipped tar archive of the the whole MyData directory\n",
"run(`tar -zcf archive.tar.gz ./MyData`);"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [],
"source": [
"# Download prebuilt ffsend linux binary\n",
"isfile(\"ffsend\") || download(\"https://github.com/timvisee/ffsend/releases/download/v0.2.65/ffsend-v0.2.65-linux-x64-static\", \"ffsend\")\n",
"\n",
"# Make ffsend executable\n",
"run(`chmod +x ffsend`);"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/bin/bash: ./ffsend: cannot execute binary file\n"
]
}
],
"source": [
"; ./ffsend upload archive.tar.gz"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"You could alternatively use the ffsend command in this way to transfer individual files"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Keep in mind that any changes you make to this online notebook won't persist after you close the tab (or after it times out) even if you save your changes! You have to either copy-paste or `file`>`Download as` a copy.\n",
"***\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"***\n",
"[](https://doi.org/10.17605/osf.io/TQX3F) "
]
}
],
"metadata": {
"@webio": {
"lastCommId": null,
"lastKernelId": null
},
"kernelspec": {
"display_name": "Julia 1.9.0-beta4",
"language": "julia",
"name": "julia-1.9"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.9.0"
}
},
"nbformat": 4,
"nbformat_minor": 2
}