{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Reinterpreting by patching an existing HistFactory pdf spec\n", "\n", "An important pattern in High-Energy physics in the reinterpretation of analyses with respect to new signal models.\n", "\n", "The main idea is that a given phase space selection (an \"analysis\") designed for some original BSM physics signal may not only be efficient for that signal (indeed, likely it was *optimized* for that signal) but also be reasonably efficient for other signals (albeit not optimal). Thus, upon generating the new signal, one can pass the new signal sample through the analysis pipeline to obtain a new estimate of its distribution with the channels defined by the analysis.\n", "\n", "The final step is then to construct a new statistical model by swapping out the old for the new signal and evaluate new limits based on this new, modified models.\n", "\n", "In `pyhf` this final step is demonstrated here is very easy to perform as demonstrated in this notebook.\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First some basic import and plotting code we will use later" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] } ], "source": [ "import jsonpatch\n", "import pyhf\n", "from pyhf.contrib.viz import brazil\n", "\n", "%pylab inline" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "def invert_interval(test_mus, cls_obs, cls_exp, test_size=0.05):\n", " point05cross = {\"exp\": [], \"obs\": None}\n", " for cls_exp_sigma in cls_exp:\n", " yvals = cls_exp_sigma\n", " point05cross[\"exp\"].append(\n", " np.interp(test_size, list(reversed(yvals)), list(reversed(test_mus)))\n", " )\n", "\n", " yvals = cls_obs\n", " point05cross[\"obs\"] = np.interp(\n", " test_size, list(reversed(yvals)), list(reversed(test_mus))\n", " )\n", " return point05cross" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### The original statistical Model\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "data = [51, 62.0]\n", "original = pyhf.simplemodels.uncorrelated_background(\n", " signal=[5.0, 6.0], bkg=[50.0, 65.0], bkg_uncertainty=[5.0, 3.0]\n", ")" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "test_mus = np.linspace(0, 5)\n", "results = [\n", " pyhf.infer.hypotest(\n", " mu,\n", " data + original.config.auxdata,\n", " original,\n", " original.config.suggested_init(),\n", " original.config.suggested_bounds(),\n", " return_expected_set=True,\n", " )\n", " for mu in test_mus\n", "]" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nOzdd3hU1dbA4d9KDyQBQpUeikAo0hWkKgIi2MCCIla8oCDY9doAuyIiNrBwpYgCiogo2KUpXWpCCRA6hFASQkid/f2xJx8hpGdmzpT9Pk8eknPOnLMS75015+y91xKlFIZhGIbv8rM6AMMwDMNaJhEYhmH4OJMIDMMwfJxJBIZhGD7OJALDMAwfZxKBYRiGjzOJwDByEZH6IqJEJMDB510sInc78pyG4SgmERiWE5F4EemVZ9s9IrLCims7g1LqWqXU9OIeLyJdRSTF/nXWnpxScn3VLU0cIhJiP1ft0rze8E4O/dRjGIZjKKWWA2Gg71KAvUBFpVSWhWEZXsrcERhuT0SeFJFv82ybLCLv2b//S0ReF5E1IpIsIt+LSGSuY68XkW0ictp+bDP79plAXeAH+6fsp3Jd4k4R2S8iiSLyXK5z+YnIMyKyW0ROiMjcnGvZP23Psm8/LSJrRaR6rhgfsH/fSESWikiS/fxzSvl3iRSRGSJyVEQOiMhLIuJn39dURFbYr3FcRGbYX7bM/u8O++98o4jUEJEl9phPiMgfpYnH8FwmERieYBbQV0QqAtif398OzMh1zFDgPuASIAuYbD/2UuArYAxQFfgJ/cYfpJS6C9gPDFBKhSml3sp1vi5AE+Bq4MWc5AGMAm4EugM1gVPAh/Z9dwMVgDpAZWA4cC6f3+dl4BegElAbeL/kfxIAvgSSgAZAR3tcd9n3vQ4sACqik91U+/Zu9n+b2H/nBcDTwA6gCvrvN7aU8RgeyiQCw10ssH8iPS0ip4GPcnYopY6gP8neYt/UF0hUSq3P9fqZSqmtSqmzwAvArSLiD9wG/KiU+lUplQlMAEKBzkXEM04pdU4ptQnYBFxm3z4ceE4pdVAplY5+0xxkT06Z6ATQSCmVrZRar5RKzufcmUA9oKZSKk0pVeKxEBGph35Tf0wplWr/G01GJ8ica9QHath/j5WFnC4TndTqKqUylFLLCjnW8EImERju4kalVMWcL+ChPPunA0Ps3w8BZubZfyDX9/uAQPQn3Jr2nwFQStnsx9YqIp6jub5Pxf68Hv0G/l2uhBULZAPV7TH9DHwtIodF5C0RCczn3E8BAqyxP7K6r4hY8lMPCAGO54rlPXscAI8C5YB/RWSziAwp4DwArwKHgT9FJE5EHitFPIYHM4nA8BQLgFYi0gLoj34skludXN/XRX/KTUS/wdXL2SEiYj/2kH1TScvvHgCuzZ20lFIhSqlDSqlMpdQ4pVQ0+o6jP/qR1QWUUkeVUsOUUjWB/wAfiUijUsSRAlTKFUeEUqqt/RqHlFI5j8oeAabZZxpd9PsqpZKUUqOVUvWAgcDzInJlCeMxPJhJBIZHUEqlAd8As4E1Sqn9eQ4ZIiLRIlIOGA98o5TKBuYC14nI1fZP548D6cDf9tcdQz9jL64pwKv2RzOISFURucH+fU8RaWl/JJWMTka2vCcQkVtyTd88hX5zvui4wiil9gKrgLdEJNw+iN1YRLrYr3GbiNRUus78afvLsu2Ps3LGFXLiuV5EGtiTZBL6DqdE8RiezSQCw5NMB1py8WMh7Nu+QD/SCUF/CkYptQP9KOl99B3CAPTgcIb9da+jPwGfFpEnihHDe8BC4BcROYN+M77cvq8GOlklox8ZLS0g1g7AahFJsZ9rtFJqTzGunddg9GDwduAkMIfzj4Y6Aevt15gHPKiUyrkLehGYZ/+drweaAX8CZ9BjMROUUv+UIh7DQ4lpTGN4Cvujje3oAdDkXNv/AmYppT6zKjbD8GTmjsDwCPb58Y8BXxcwE8cwjFJyWiIQkWkikiAiWwvYL/ZFQXH2WQ1tnRWL4dlEpDz6ccs1wEsWh2MYXsdpj4ZEpBt6VsMMpVSLfPb3Qy/O6Yd+xvqeUuryvMcZhmEYzuW0OwL7opSThRxyAzpJKKXUKqCiiFzirHgMwzCM/FlZdK4WFy4COmjfdiTvgSLyIPAgQPny5ds1bdq0xBfbuXMnZ86cwc9PaNWqCf7+5Th+PJHAwEAqVqxYut/AMAzDQ6xfvz5RKVU1v30eUX1UKfUJ8AlA+/bt1bp160p8jtGjRzN58mQCyynefXc7IgHcd58/TZvWZMqUe6hQoTOPPz6Fbt26MWRIYYswDcMwPI+I7Cton5Wzhg5x4WrQ2pxf7elw7733HuEVQ0g/C9N2wd69WXzySTojRuwlLm4Uq1e3Y9myz1m1aizx8eNISPiZzp078c033zgrJMMwDLdgZSJYCAy1zx66AkiyF85ymoEDe4GCmV/AAw/C119D+fJ6X0AATJliY+DA3cTHj2XFir6kpa1m//4XiI8fR0zM9/TufQ1r1651ZoiGYRgu57RHQyLyFdADqCIiB9HT/gIBlFJT0OWA+wFx6KJe9zorlhy7dybqb7ZC9tWQ0Dy/uPW/kZEwcaICthMfP5bYWIiL82P//mepVetW9uypyZw5P/P8889TvXr1i09kGIbhIZyWCJRSg4vYr4CHnXX9/Pj5BSP+oM4B7WFhEsgOaLQLrr0W/P0Lfm2zZjBtmg34nZ07f+fHH2H6dH9GjKhKePhg1q07QnJyMv369cPPz6zTMwzDc3hciYnSDhbnuOnJSiyYcFpXn2mEvh+ZBS+8AFddVbJzZWZCoL3I8PjxEcTF+RETs5ywsBYcOnSImjVrIjm3GIZhGBYSkfVKqfb57fO5j66D+rTUD6iWoesrNgLuhdhStAIPzFVp/tlnk3nlldOsW9eStWvb0KVLa+6885aCX2wYhuEmfCoRbNiwgdG325ta7UffDQDUg28OwQfr9AByaW6SAgOhrj2ZJCdvZPDgRNq2/Y4tWwawd+9M+vTpzcqVhTWJMgzDsIZHrCNwlIiICCpViiAtMpWzu4DlwKXn93+7CILWwtVXQ9V8l10Uj78/9O4NYOPEiUWsX7+Ibdv8OHiwGmlptTl7Noz09HRq1qxZtl/IMAzDAXzqjqBRo0asWfMl1z6Jfjx0ED1fKUd3yBgG/2Q69rr168OMGTaqVfuSVasa8txzVxIVVZ/jx4879kKGYRil4FN3BAChoY1pFgE0AbYC64Gu9p1+QCWYtAt2LYcm5aF/f8dc9/xEomy6dt1BWBjs3dudjIyRzJlzlqZNm9OvXz/HXMwwDKMEfC4RPP30BD75FL2OGSCG84nATilY9DtsLwf9+uV+E3eMWrX0V2pqLNu3P8x77/nTpUs7evfuSkBAOEopM9vIMAyX8alHQwBRUVE0jy5HubZAZSC/tQMC3Arx18NmJ7dA8feHzz/PZvDgNaxaVY/ff3+YZs2amBXMhmG4jM8lgjFjxvDFF7257EqgNXqc4EQ+BwZClj88/y88Oxb27nVeTAEBEBYGWVmn2LnzIwIDd2Oz/Y+MjESSkpLwtLUehmF4Fp9LBAChoY2oL+g1BKDrmmbkf+zZM7BmI2yIcU1szZrB++/bOHfuY1avbsjAgZdzzTVXm2RgGIbT+Fwi2L9/P507f8LXjwA7gGpAOrCtgBdUANvD8FN1OJvlsjAByMpKpkOHHbRrt4EjRz7BZsti9+7drg3CMAyv53OJoHr16lx9dSduHAg0BTrYd6wq5EVBsOcsPP4jvDQWMgq4e3A0EejTB669NomdO4fz8cdRNG7cmIULF7omAMMwfILPJYLg4GCmTfuURx6GOlFAc/Tg8DEgofDX7tgHa2Pg9Gnnx5mfevUOcv/9iho1JnL2bAyHDh0iK8vFtymGYXgdn0sEAMHBtUlPD6aGvSo1Dez/bijihZfBuQfh+xQnBleIsDC4805ITV3KqlWXcfXVrbn55husCcYwDK/hk4ngueeeY9CgTNa+DewG2th3VCvGiwNh9j54dhL88ovzYiyKn18Wt92WSPfuq0lImEN2djaJiYlFv9AwDCMPn1tQBtCtWzdOnpyPrcIOPgX9VwgCDgBti3ECBas3Q2BmTk0h1xOBnj0BThATczu//96Ujz46ypo162jYsKE1QRmG4ZF88o6gb9++PPnk9dzSBwLLo5NAM3TJiaXFOIE/qDtgU2c4cs6poRZbgwbb6dUrhcDABShlw2azWR2SYRgewicTAYC/f30SjkCV7UA20BLIBP4EThbjBIGQnAVPrYCXX4VzFieEqCgYMSKLPXueYOnSzrRqFc1PP/1kbVCGYXgEn0wEmZmZXHrpaKZOhSNfoWcMRQGh9gP+Lf65Du6HpX/DHieuPC6po0dX4+e3C/gdpcydgWEYhfPJRBAYGMi4cU/Trx888DZQA11zqJX9gH/RdwnF0QiyH4Hlwc6ItHRq1ID33rNRrtxENm7sziuvPMUPP/xgdViGYbgpn0wEAM888zKdO4fSPZrzf4WW9n9TgF0lOFkIzDkIE+fDb785NMxSyyleevLkCr788h1mznzT2oAMw3BbPjlrCCA7O5uTJ+uQfHgnwash/XKgFlAR3awmrYQntMGiH2BnBd3hzF2qSPv767uD7OyVxMYOJSTkv2Rm+nHppZcW/WLDMHyCz94RvP/++9x0005WrICMX9BF5wS4zP59SWdg+oG6DY4PhFMO7nBWVkFBEBoKx47N5J572tO9+5WkpZU00xmG4a189o6gb9++ZGT8QlTUEvy6w5dH7TtaoqeQbkFPKa1UgpOGwkkbvLgBOu2B22/Tn8jdyciRZ9m3L53ExKnUrj2arKwsAgJ89n8GhmHgw3cEzZo14667bqZaNWidu1F9FeASYBkwHSjFpJtt6+Gzz2DzZoeE6lBVq0L79lnExY3h3Xc70KbNZRw6dMjqsAzDsJDPJgKAEyfCOXIEjq9Hrx/I0Qo9RnAa2FOKE7cAhkNiTQcE6UQ22zoiIvYREnLc6lAMw7CQTz8TGDDgKRo2hAoVIGgHZPRAjxM0B35G/3XWc76BTUlUh3d2gt8RaBwJdes6Lm5HadcO2rU7S0xMN+rX/5Sffz7D/fffb/olG4aP8elE8MEHH3H48CAaNEjHrx98m/OEJAK9wOwIunnNGSC85OdPz4DXx0GrhjBxgqOidrzs7DNMmnQ7EydC8+bRdOrU2eqQDMNwIZ9+NNS/f3/atGlMQAC0rJBnZ87jIRt64Lg0AiD7FuBmsLl5p8n+/eH996F8+VfJzDxNhqu67xiGYTmfTgRJSUls3VqRc+cg9gdgSa6dzdCrjZsBl5fhIrXg30yYHg87d5bhPE4mAi1awMmTP/H1161o3LgBa9assToswzBcwKcTwfLly7nvvhXs2QMZKRCanmtnCNAE2OeYa838CkY8BPscdD5n8vM7QLVqxwgPP2h1KIZhuIBPJ4LOnTvz1VdjqF8fHnkEuv8nzwGXoVcZ/wAsKtu1VFsI6Q/lq5ftPK5Qqxa8/XYWCQmDOXp0FosWLUIpN3+2ZRhGqfl0IoiMjKRPnxsoX17/3DIizwGNgHLohjXrgeQyXCwUUtvAq9vhdDJ4wvuqUhl8/vldDBgwgAULFlgdjmEYTuLTiQAgNjaVjRv198s+5sJP/v7oNQGnAEWJylMXZOMBuONumDev7Odyhc6dYexYaNr0e2w2N6udYRiGQ/h8Ihg7dhJTpuh581G1ICQyzwGXoWcOVUE3ty9ref9wONcUghqU8TwuIgLdu8OxY9P566+eDBhwrVmJbBhexqmJQET6isgOEYkTkWfy2V9XRP4UkX9FZLOI9HNmPPmZNGkSb7yhK3H+5z/Q8fo8B9QEKtu/TwLiynhBAfrC9AxITAdP6igZG7uSf/75jb173bB2hmEYpea0RCAi/sCHwLVANDBYRKLzHPY8MFcp1Qa4HfjIWfEUJDo6mkaNzofVIgLdsjJHTkXSRPRvUYqFZfk5nQmjJsPzz3tOMmjeHL78MovAwEc4dy6epKQkq0MyDMMBnHlH0BGIU0rtUUplAF8DN+Q5RqHX8QJUAA47MZ58nTp1ip9/zuDoUT2A++2zwK95DsppWFMDXZDOQY4CRwWyshx3TmcLDoZz5+J4//12NGhQj3//dcDAiWEYlnJmIqiFnm+T46B9W25jgSEichD4CRiV34lE5EERWSci644fd2yBtFOnTvHEEz+ybp1+Hn5d33ye31cC6gGbgOOU/fFQjg6wrzdsPeug87lQvXon6djxHNWrn7Y6FMMwysjqweLBwBdKqdpAP2CmiFwUk1LqE6VUe6VU+6pVq150krKoX78+a9bMo29f/fNdQ6BVp3wOvAw4ASwAFlL2QWMA0ad5ZRWMexVSUx1wThepXh2efjqDuLgBHD36I7/+mvc2yjAMT+HMRHAIqJPr59r2bbndD8wFUEr9g17PW8WJMV3Ez8+Ptm37ExBw/k9xaQCQ91N6NLpEXyh6PYGj7gqAU8dh+UrYvdtx53QVm+0sL754A3369CE2NtbqcAzDKAVnJoK1QGMRiRKRIPRg8MI8x+wHrgYQkWboRODy4vgbN25j3jw9bzQjA+aMBlblOSin5MQhIAxY58AA6kH2aNhZkm5obmTQoGzGjxcqVzaziQzDEzktESilsoCR6Mr+sejZQdtEZLyI5EzSfBwYJiKbgK+Ae5QFtQyWL1/Oxx8nkpSk+/s+NBL8885vAl2R9Bx6vGAneqGZo4TA1D0w71fwtGn6QUHQpYuNmJg7+euvN/jkk0+sDskwjBIQT6sh0759e7VunSM/jkNKSgr7979OQsJr/7/toQ0QeybPgdnAO+gh733oOVDNHRjIOfCbDNf0gGeecuB5XWjCBFi/viLbt8dToULe2t6GYVhFRNYrpdrnt8+nG9PkCAsLo2rVtiQk6J+zs6HWKYhN5fzkVtAlJ1qiHwuNzrPPEULBdjf4tXDweV1ozBhITEwiNfUbKlS43+pwDMMoBqtnDbmN2bO38OOP+vvkZPjtFSC/R96t0HcGu+w/pzk4kBqwOBF+PQR79zr43C4QEAA1aih27BjG88/fzOuvv251SIZhFMEkAruFC5excqX+c1SqBP8dh54ymldN9LymTcA3wJfOieetN+DxJyE9vehj3ZFSijVrvmP16m+xecrSacPwUSYR2C1evITJk8+/81/TDerm1ztAgNbo+U4V0Uvmjjo+nqwrodJACAxy/LldQQSefhpGjVrPsWP/I8uTlk8bho8xicAuKCiI8uVb/v/PKSkQuQNdaC6vy9AJIQs9yrLWCQFdAntqwtyDnntX4O+vv9avH8YVVzQzs4kMw02ZRGCXkJDA+PHb2WwfF0hKgo2fc34sILdwdNOabeh+BZtx/FiB3WdL4NbbPW9KaW7BwYpy5eIIDvbAFXOG4QNMIrALDQ3lt992cNT+mKdmTXjrQ6BtAS9oA5xBF6LLJP+BZQfIrgZZdUA8eH5XYCCMGwf1679LYuIiTp1y5AIMwzDKyiQCu/DwcA4e3Ebv3vpnEegQDXXKF/CCS9HlJg4At6DHDZyhIqQOhG/K0ibTDYiAUpnMnXsz9erV5ueff7Y6JMMw7EwiyCU4uBYBAefrPJw8CaHLyL/oRQB6Kul2IApw8qDud7vhsZcgJsa513G2+vUz6dIlnQYNzlkdimEYdiYR5LJo0SIee8xGRob+WSmIW8zFpfJytEavKdiK7mfszAKcApu2wdYdTryGC5QrB089lc3Ro0NJTl5HfHy81SEZhs8ziSCXzMxMMjKCyWm8VbkyfPUd+BX02OcS9BjBv0AC8A963MAZQsA2AjY3dtL5XSw7+wyPPNKdtm1bmx7IhmExkwhyuemmm/jxx7HkbnlQLQIuLaw9ZRvgCLoQnQ3d4N5ZAmDlCfhsOaxf78TruEifPqnccQdERmZbHYph+DSTCPLIvZYAICEBUr9GLyDLT0t0DaJ4oCG6DpEz39cUfDUFPvpEP7ryZJdcAoMGJbFlSx/27t3CsWPHrA7JMHySSQR5PPfcDCZMOP9z+fJwJp6CH/mUQ/cp2Ay0sx+33YkBCtgGQeg9uuGzNzhzZju9enXk1lsH4WnVcA3DG5hEkEeFClUIDw/7/5/Ll4evv4agloW8qDWQin5nbo5ODs4UCdsyYO5+8IaxVn9/GD48jf/8Jx2lMqwOxzB8julHkI/Nm6/l5MklF2x7fBNsKKhPezbwLrog3R1ODe0Cfr9A0L/w5SyIjHTddZ2pSpWbOXnyEa68sgv+/v5Wh2MYXqOwfgTmjiAfeccJjhyB+IlAQVM3/dH1h3ahHw2l4NzHQ3a2DhA+AMK8qP/L6tXz6dGjB5MnT7Y6FMPwGSYR5JGYmEjv3jP+vzcB6GmkNSqj3/AL0gb9aGgzsBSYB5x1YqAAleB4NEzfB95S3LNhQ3j2WejXL8XqUAzDZ5hEkEdkZCRt2rSmcuXz24KC4P13ILxpIS+sAtRBrynogH5c5KIpnl/9Dbfe4ZmNbPJzzTVw5MhLxMdP46+//rI6HMPweiYR5OHn58ecOQu44ooLP/77CbQKR5eeLkgbIBHd4L4Bujy1C6bIq0qQWhEyPGu4pwiKZ58dRp8+vTlw4IDVwRiGVzOJIB/+/iEEBja8YJ7+kSOw5ml0OYmCtACC0XcCl6PHC2KdGGiOMEi/E370ksdDOe6808a4cQFUrHjC6lAMw6uZRJCPWbNm0bNnHLmrJVevDlf1QT8CKkgQuhDdNqAWEIlTupcV5If9MHYSbHfBQLUrhIXBFVecY8uWfqxe/ZOpS2QYTmISQT5atGjBffd1umCbnx88MwaqNyzixe04X4huONDLOTHmKwuW/wUr17jwmi5w9uwRbrzxBh544B6rQzEMr2QSQT5at27Nyy8/lu/c/KY29OKxgtRA3w2sBwLt2wo73pFCwTYcjnZw0fVcJDAQXnghiyefTMdmMwvODMPRTCIoQGhoNMl5msEcOwZLnwO2FPHidugeBvvRA8YT0WsLXCEUfkuA+Vtgh4eXrM4tOhqCg1exffv9LFq0yJSiMAwHMomgALff/jhPPy0XbKteHYaPBooqBd0CPV6wHqiPnmnkzKqkeSn48DV4c4LnF6bLa8GCWQwYMID58+dbHYpheA2TCApw3333c9tt9S7aftuNEFWniBfnHjQuj65K6qKppIAuTHcDVB6qW0R6k06dYOxY6NTJLDgzDEcxiaAAN998M7fd1v2i7RkZUOsQUFT/9fboN/7NnJ9Kus3RURbiElhngyVHwZt6xYtA9+6wa9eDxMV9z9q1a60OyTA8nkkEBVBKkZpaj9N5Cs2lpsKKSRQ9TpB70LgRetqpBbN5Jk6B+4fBGWd1TrOIUhncddcgbrxxAOnp6VaHYxgezSSCAqSmptKhw8t8//2F2ytWhLfehYBO+b/uAjmDxgeBm4DbHR5mkTKbQGhHCAlx/bWdbdSoLMaPD0AkueiDDcMokEkEBShfvjxTp75Ht24X7+vQGi6revH2i+QeNK4FhBV+uFPUhMOXw48JFlzbyWrWhIYND7F1600sW/YHGRlmaqlhlIZJBIUYNmwU0dGXXrQ9IwPKbwD2FHGC3IPG59AN7qeh7xJc7OOVMHykd40X5IiJWclVV/XijTfesDoUw/BIJhEUIi0tjZ0765OUdOH2gADYuJCC+xPk1g49fXQTunPZIWCVgwMthgyBvQlwzAvvDGrXhv/+V3HLLTarQzEMj2QSQSG2b9/O0KG/kLchmp8fzPgC6g0sxkku4fygcXl0A5tNOL9XQV7VIGMEbHZ2G02LXHUVHDs2nkOHvmb9ehfV/zYML+HURCAifUVkh4jEicgzBRxzq4jEiMg2EZntzHhKqkWLFnz99Tt07HjxvgoVoFPli7fnK2fQeB/QCX2HYMWsRz/4bA98/jUct+DxlPMpxoy5i+7du3L0qAur/RmGh3NaIhARf+BD4FogGhgsItF5jmkMPAtcqZRqDoxxVjylERAQwC23jCQiIjjf/ccWAiuLcaIWQAh6+mhV9MrkNUCmgwItgczT8OU0+Gmx66/tCoMHZ/HUUyFUquSq1XuG4fmceUfQEYhTSu1RSmUAXwM35DlmGPChUuoUgFLK7Z5gHzlynN9+q0lmPm/aGYkQWJz5+UFAW3RvgiSgG3DxWjXXqAjqQXsMXigyErp1O8WWLQPYvn0TqamuqvhnGJ7LmYmgFpC7tdRB+7bcLgUuFZGVIrJKRPrmdyIReVBE1onIuuMufqaxcuVKXnllL3vymSH08svQfWgxT9QB3dN4Hbql5eWcr07qalVh1n7495C3PiKCI0f+5YorOvDoo49aHYphuD2rB4sD0A9KegCDgU9FpGLeg5RSnyil2iul2letWpwJ/I7Tt29fVqx4m8b5FJoTsY8TFKewWyWgCXrQOBM9TrAOiHdYqCWSlQVPjYE337Lm+s4WHg4jRmQyZEiQ1aEYhttzZiI4hP7sm6O2fVtuB4GFSqlMpdReYCdF1/Z0qYiICNq0uRG/Av5S62YBxR3i7ojuTbANEGAZsNQBQZaGP2T1gmr9Lbq+C/TpA9nZH3Ds2Ffsye+WzjAMwLmJYC3QWESiRCQIXWBhYZ5jFqDvBhCRKuhHRW73/9iNGxOYMyf/eZcN60G1uhTvrqAB52sO+aEfD+0Fjjgo0JJqCj8DMcn6DsFbvfLKUFq0iGaHNzVoMAwHcloiUEplASPR7zWxwFyl1DYRGS8i19sP+xk4ISIxwJ/Ak0opt+tU/ueffzJtWhrnzl28b+BAGHg3+hN+UQR9V3AYfS/UDt3sfoXjYi0pG/Dsh/DIaMj20ok2PXtmcdddQdSpE2p1KIbhlsTTOj21b99ercu7wsvJkpOTOXz4LY4efTXf/fvPwt1/ABHFOFk68A56vGAg8Bt6CupIoLjrEhxtK9RPgCn/heD8Z8p6hfDw9tSvv5Dy5SsR4o1V+AyjECKyXinVPr99Vg8We4SIiAiqVr2ywP3fTwO/KeiP10UJBtqgxwnOAFegR0/SHBBoabWAfVdBjJfPtExMXEeHDk146KGHrA7FMNyKSQTFNH/+LmbNyn9fjx7QeiDFSwSgp5La0PA2O8YAACAASURBVDOIwoD7uXhirYsp4NW/4Y23yHfNhDcIDoZ+/c7Qt6+XtW0zjDIyiaCYVq7cwNq1+T9OaNkShtyEngxbHFXQzWrWoaeRgp5NFFfWKMvmxCH4/S/Yu9faOJxp0CCoVu1/HD8+n6S81QQNw0eZRFBMn3/+OV99dXOB+xsFQXBJ5jt1BFLQw+gAvwJz0AnBKk0g6xE4YdVYhcsopky5g6ioumzb5sr+oYbhnkwiKCZ/f38iIi4vcP/iHyF9BkX3Ms7RCL3ILKd9ZSf0QrPVZYnSAULhnZ3w6zJIs3LcwsmaNk2nU6dsqlUzg8aGYRJBCYwbt5S5c/Pfd9VVcOeLFG/mEOi/fEd0EY7DQDX0TKLV6JlFFjpxEF57Cb77zto4nKlqVXjyybMcPHgPmZnnyPTWgRHDKAaTCErg2LEMTp/O/09WpQrc0gX8/EtwwtbognT/2H/uip49ZHU5/WrAEIj00sJ0uZ08uYKrr76UMWPcqvCtYbiUSQQlsHDhDzz+eNsC92cmQbU16LaUxRGKrkq6FTiNnkYaBSSWMVBHaAQfxsPBJDjr6iY6LuTvD3XrHqRePbdbx2gYLmMSQQmFhxc8TpCYCEd/QjegKa4r7P/mtK+8E7i+gGNd7EwaPDgCJk2yOhLneuAB6NjxG06e/JVsb11ebRiFMImgBE6fPs2QIb+wuICmLk2awCezwK9pCU5aEd24ZgP6TiJnCuopir8uwVkC4FxrCM93LaK3yWbOnIFERzc2BeoMn1PiRCAifiJS3CFRr1KhQgUiI2sWWIZBBBrXgtYXFdIuQmcgg/NjAweByejVx1a7HBaXgwOp4GHVSEosIuIMISFHycgw6wsM31KsRCAis0UkQkTKo59ox4jIk84Nzf2ICIsX/0nv3gW/06elQfp8YHMJTnwJujLpKvQCs5roRWfLsP6uAEizwVPTYOxY704GtWrBpEnnyMh4DqVseFodLsMoreLeEUQrpZKBG4HF6CHNu5wWlRsTEcLCOmAr4A06OBgyjoGklPDEndELzLag/6t0Rze8jylDsA509CzsTfLutQWg7+oSExczbFhHnnvuOavDMQyXKG4iCBSRQHQiWKiU8tlJ1xs2bKBXr5Vs3Jj/fhGY+jG0v66EJ24IVAf+Rhf+iUY3uv8Lt7gr4HI4fBMc8OK+BTlE4OjR9Rw7ttbcFRg+obiJYCq6qWJ5YJmI1EO3Yfc5UVFRXHNNJ8LDCz5GBHpW43wdoeIQ9F3BcWAX5+8KTgFHSx2u4whkA+PXw9RPvbd3Aej/fo89BkOHLic5eVXRLzAMD1esRKCUmqyUqqWU6qf0R6T9gJdPKsxfpUqVmD79Kxo3LryC5bppIF+W8OQtgHD0XQHou4JH0GMGbuJQLMyZAzFu8sjKWfz8QKl0fvvtevr06cmRI1a1kTMM5yvV9FF7MnjXwbF4jKCgqmRkNC30U3GbVlC7OcVrYZnDH72uIB5ddsIPqGDfV9xFas4WDWoUpLpRcnKm06cT2bhxJbt2xRZ9sGF4qLKsI/DZou5LliyhT59YCmuB278/DBlCyf9KOe0rV+ba9iPwGfrZjDuoCG/thNWbyLd9pzdp2BBmzsykcuWpVodiGE5TlkTgs6NoHTp04LnnhlKlSuHHXVERAkr6RCEEnQxiOF/JNAo4gZ646yZOJsCzj1Fgsx5vEhQEx4/P5Z13bmbChAlWh2MYDldoIhCRLSKyOZ+vLeg5Lj6pcuXKjB37IdWrF96JZsn3kDWVkg+rX46+k8gZK2iK/msvxa3uCtQgqHK11YG4hlLw00/fsXDhDFOGwvA6RfXUuhn9FnQgz/Y6uMdcFssoFczu3c2oWXMLBfVB794djgbBd6ElPHkFdGXSDeiKpBFAD3Tjmi32fe4gGj49DO2r6dYKYWFWB+Q8IvD00xAQEEdq6mbCw9tYHZJhOExRj4beBZKUUvtyf6E/4/rsYDHA77//zv33b2FzISuIq1eHYQMguDS9T7qg1w/kviuogW5v6UbOZcGIx+HlV7x71THoR0R+fudYs2YAw4ffy/Hjx60OyTAcoqhEUF0ptSXvRvu2+k6JyEN07dqVGTNeJjq68ONsaVBvOyV/PBQJXIZ+4z+DflQ0CPdbzy1wtgkEtdKfmn3B3r2HmDlzBn/++ZvVoRiGQxSVCAorn1bSBx5epXz58tx559NERJQv9LjkZNg5GyhkhlGBuqLHBHIa11RBzyjKxn3GCgDawd+1YEuS998VADRqBLNn22jZsoAytIbhYYpKBOtEZFjejSLyANb30bLcyZNJ/PlnFMnJBR9zySUw5XMILriNQcEqAy2BtUBOc5hU4GP7NjdiA178DkaN9v56RAAVKsCxYzNZsOAhpk41U0sNz1ZUIhgD3Csif4nIO/avpcD9wGjnh+fedu3axdixW9m0qfDjmjSALlVLeZGu6Kb2OXcFoejVx8vQbS3dyOls2J8MZ85YHYnrfPTRx7z55jjOefuCCsOrFTprSCl1DOgsIj3RBRAAflRK/eH0yDxAhw4dWLPmG1JSBhV5rN9fwEn0IHBJVEX/5degaxGVA3oBn6IHkq8q4fmcqRGcaQBrs6Cf1bG4yGOPQUZGEllZO3Cf6VyGUTLFrTX0p1LqffuXSQJ2AQEBtG9/M0FBRX/cTz8KoaeKPCx/3dCNa3Lqn9UCmqPvEtzt07cfTI6Ftz8AX5hUExICERGpbN7cn9dff4Fjx45ZHZJhlJhpVVlG+/fvZ9q0yiQkFH7ciy/C7aNKeZFq6AJ0qzlfc+gq9ICxm00nBUg/DYsXwcp/ij7WW+zde4hx417jiy8+tzoUwygxkwjKKCUlhZkzdxEXV/hx/v7QtzpIRikv1A1IRycD0APJ96JLVbubyqBGwtHmVgfiOrVrw6ef2rjuuo2mh4HhcUwiKKPo6GgOH95M585FH7vtH5AJQCGzjApUA72obBXnB4nroP8LumOboHCYexAWbabQ4nzepE4dSEycx6pVjzFt2jSrwzGMYjOJoIxEhCpVogkJaVDksZdeCu16Uvpyfd3QSWB1rm0H0Gu8D5bynE6kbPDuazBhom+sL8gxceIkxowZSUJRzwsNw02YROAAsbGxPP98NvHxhR9Xqxa89hRULKJqaYFqApeiB4lT7duqoVcd/4r71YP1A9tNEHKn+4XmTMOGwQcfZBEcvMvqUAyjWEwicICwsDDi4tI4caLoYwP84HIbkFjKi12NvitYYf85GD1OsA8oYpzCEjVgq4IZ8bB3r9XBuEZQENStm8nWrTcye/YH7Nmzx+qQDKNQJhE4QJ06ddi5cyvt2hVdbCczE5a/gS4pXRrVgVbodQU59Yvaost//op7NLrPx4y5MOxB30kGAKdOJfLQQ2MYP/4Fq0MxjEI5NRGISF8R2SEicSLyTCHHDRQRJSLtnRmPMwUFVSEs7LIijwsMhPHjoEnRa9AKljPO8Jf95wD0IrME3POuAFCtILgPRFxidSSuExYGEydmM3z4AWy20k4XMwznc1oiEBF/4EPgWvQs+MEiclGtThEJR5erWJ13nyfZtGkTt9++ny0X1Wq9WLt2cH3jMlysEtAB2Ih+8wf9F34APYbgjspBagd4YyekpvrO4HGDBpCWtpxNm+5l6tSpZmqp4ZaceUfQEYhTSu1RSmUAXwM35HPcy8CbuF3lnJKpW7cuUVEN8SvmX7RaAgT8WoYLdgUCgZx13gLUtn+fmu8r3MK6fXD7PbBggdWRuNasWbMZMWIEq1d79Ocdw0s5MxHU4sLOZgft2/6fiLQF6iilfizsRCLyoIisE5F17toMpFKlSixZ8hctWgQV6/gDe8BvE6UvEVEeuBLYDuzPtT0WmAi4a6WDMDhTF/xrWh2Ia/XvDx98oIiK8qFBEsNjWDZYLCJ+6Lesx4s6Vin1iVKqvVKqfdWqpS3j6Xz+/uUIDu5SrDLM110Hb3+BriRaWp3QCeE3zs/PrIceM1iCe87Z9AMGwEwbnMrwnUdEIhAdDdu338vSpZ+xbNkyq0MyjP/nzERwCL32NUdt+7Yc4ei6mn+JSDxwBbDQkweM9+3bR8+ey/mtGI2rgoKgVWWoV47Sv2EHoaeO7gdypqyXQw8m70XfLbipxAwY9TE8/wLY3HSmkzPYbOmMGDGckSMfxOZLv7jh1pyZCNYCjUUkSkSCgNuBhTk7lVJJSqkqSqn6Sqn66OIJ1yul3LCMWvHUrVuXMWNG0qRJ8f6sp07BuY/Rg76l1Q7d1vI3zk8dbY9eaPYz7ll+wu5QBhxI01NqfYUIvPRSNq++mkJmprs+vzN8jdMSgVIqCxiJfjuKBeYqpbaJyHgRud5Z17WSiPD66xPp0KFHsY6vWBEa1YbgsjT99EdXIk0ANufa1hdd02h/Aa9zBx3gQH9Y426ltJ2senUIDz/E5s3X8uGHk0gurMWdYbiAU8cIlFI/KaUuVUo1VEq9at/2olJqYT7H9vDku4Hczpzpxu7dRR8nAq+Oh+t7lfGC0cAlwJ/ovgUADdCTchuW8dzOJPrrjXXw1PNw9KjVAbnWli2bGD36MT791LS6NKxlVhY7mFKKu+76jE8/Lf5rbqwBfvFluKgf0Bu90vjvXNsr2P9189pnqamwYTNs97HSPA0awIcfKq69dqvVoRg+ziQCBxMRZs6czUsvdSz2a/5eDLYvKNsbdhS6a9kKIHcntC3AR+jBY3dVGbJHwYrSFuPzYE2aQELCDNaseYypU82dgWENkwicoGvXrrRsOaTYx/fuDcP/C5T1jbA3+nHLL7m2NQEqAj/g1gPHBMHvCfDuT/BrWRbaeajJk9/l0UcfYf9+dx7UMbyVSQROEhNzCd9+W7xjIyLgtmugVaUyXrQCesVxLJAzRhEEDABOAsvLeH5nU/DDfPjyG9+aUgpw773w4YeZhISssToUwweZROAkixevZPbswBJNjay7lbJXXOqErkW0GN3TGPSAcSv0YyN3nrEooAbB6cGQ4GM12vz9ISpKERs7hJkzx/PLL78U/SLDcBCTCJzkhRdeYPXqVwkMLP5rTm2HkH1lvHAg0Afd7yD3h8s+6DuGU/m9yI2EQpLAc5tg5mxIT7c6INfKykrnlVfG8cYbL5oCdYbLiKf9j619+/Zq3TrPmGWalnaAVavqUdylw2lpsOQEvFfWUtIK+BJd6WkUEGbfno1eY+AJ9gIz4MUXoGdPq4NxrZMnISKiKp07r6RcubKUqTWM80RkvVIq38oN5o7AiVat2s3TT4dz7lzxjg8Jgb41INxG2RrMCHpBWSZ6xXEOf/t516MXm7mzKGAEHCy6FbTXiYyEgIDjrF9/DaNGDePgQTdsSG14FZMInEhESE4uR0l6mB/ZD+nvADvLePEq6OpNG7mwsX0yevzgpzKe3xWqwf/i4btN8O+/Vgfjenv37uOLL6axeLGP1ew2XM4kAifq1q0bGzb8Q716xX9N3brQvTsEVHRAAN3Rj4V+4vwdRkWgB7ogXawDruFkCvhgErw5AbKzizzcq9StCzNm2GjXbjZZWSlWh2N4MZMInEhECA2tT7lybYv9JubvD/99EnoV3fWyaMHotQWH0SUAc3RC9z7+CY9oB2S7EdRdcNbHppQCVKoEycn/sGDBVQwceBNnz561OiTDC5lE4GS7du1iwIBdrFxZstf1jQDZ4IAAWqKnj/4GnLZv8weuB1LQJQHdXUVICIUXtsFvf0KGj00tBdi+fS1//72EQ4fMgjPD8UwicLKoqCh69epFZGTJXrflD2ARZZ/uKegFZYJeXZwzgakWcA26I4SH2Byri/QtvKhkoffr0gW++CKNzMyXUCrbTC01HMokAicLCAjgyy/nc8UVrUv0uptvhjc/Bv8SJpB8VQR6oVcb5+590Jnz1Uk94bFLTeAuON3G6kCsERwMx4/P48knL2fw4MFk+9qgieE0JhG4SFDQQPaWoPBbuXLQ4VI9ndQhb9LtgbroR0F56/8vBebinq0t82oIXx6E73bDGh+txpCaup6UlA2mw5nhMCYRuMjdd89jwgQp8etCl4HMoOxv0n7ocYFM4Mc85wtCzyLyoCma70+BF16CpCSrI3G9W2+Fxx/fxf79/+Xs2bMmIRhlZhKBi0ycOInXXit5Y7b6l0CD+kCWA4Kogu5nvB2IybX9cvQCrsXo4nQeQF0NfnfBqQCrI7GGCOzePYEuXRozatQoq8MxPJxJBC7Ss2dPevZ8rsSvu+46ePtZCA1xUCCd0N3MfgJS7dv8gBvt/36HZ4wXhEJaLXhmCyxdA4mJVgfkeoGB0KrVEVq2PGx1KIaHM4nAhU6frs5nn1UnpYRrgyoFwdXgmP7D/sANwDlgSa7tFYDrgEP2Lw9xLBnGj4MPPrY6EtcT0eWrmzZdwO7dTxETE2MeExmlYhKBCyUmJjJ3biLbtpXsdUrB1mng56g5/zWALuhm9ztybW+JLlJXx0HXcYVgsA2GE70gw4ffA1evfpt27S7jtddeszoUwwOZROBCbdu25eDBeK68smQdaETguf/CHc86MJhu6NXF33N+FpGgexmAThAesOoYgDqwNQPGbYVv50OWI8ZTPMwll8CwYVn07n3E6lAMD2QSgYtVq1abGjXuLnEHrkaN4I6mUCkQx0zzDAAGARnAfC4cFzgFzEEnCU+YUmr391r44H1KvIrbG4jotSepqR+xa9cTTJ8+3awzMIrNJAILvPbaYV59teSvU+kQPocL6waVRVXgWnTt/79zba+EXoAWC6xy0LVcoREwDGJrWx2Itb777h3uuece5s2bZ3UohocwicACDRq0JCqqPiWtEhAaCpeUgwqhDgymLRAN/MGF5ao7AU2BX3HMILWr1II5B2HKWvjuO6uDsUbHjvDuu9Cq1Y8oZe4KjKKZDmUWOX78W7ZtG1Ti1ykFa0/B01scGMw5YAp6jGA4EJJr+yfoNQwj0dVMPcUSCN0GX82EChWsDsY6WVl9mT69Mh999DHh4eFWh2NYyHQoc0OVK99AfHzVEtfYF4GOkdAiAcdN8wwFBgJJ6EJ3Ktf229BzVz0pCQBcA2n3wqpidofzVqtXL+GHH+YSFxdT9MGGzzKJwCJ//rmUe+89zj//lPy1GRlweC74/130scVWF92wZiuwKdf2GkBOvTxPKoXvDyoS3toBb8yEBT7a5OvKK2HWrEzgCbKykkhPT7c6JMMNmURgke7du/P++2/QunXJu8kHBcGEt+D+xxwcVFegHroWUd6VuvuASZS9haaL2RT8vBIWr6TEYzLeIiwMkpJW8MYbrWnevBn793vSoI/hCiYRWCQgIICRI5+mfv0bSvX6qCi4NQqiQnDcfH8/4Gb01NK5QO4PjzWBSOBb4LiDrucKAtwCu/rCn8cp8bRdb1K5cjx16yYQFpa3/Kzh60wisFhcXHdmzSrli22Q/Rn6ub6jVECvLzjOhXWHAoHB6CTxJZ71mCgAVAC8uhkeGA2LF1sdkDWaNYMXXzxLbGxPjh1bxtKlS60OyXATJhFYbMWK/fzySznSSvGp3t8fbuwHbbs6OKiGQB90ldLc7xUV0ckgBfgKXdLag9iAvemw21NWTDtJZuZxHn/8aq65phfx8fFWh2O4ATN91GIpKSlkZGxn8+aOlHYZb1ImDF0DyY4sraDQK4s3Arei1xrkiEGvLeiN532UUODvB880gcv8oGpVqwOyxrlzsH69H/fdN5WaNR+wOhzDBcz0UTcWFhZGZGR7IiMHlbrJSkQAdNiNY1cBC9AfqI1+RHQ0175ooC/6fz2eNglFIFvBaytgyFDfnU0UGgpdutjYuXMYCxc+yMMPP0xGRobVYRkWMYnADSilePjheF5/vXSvF4GsPRB+AMfWBgpAryMIAb7m4nGBk8AHXNgH2UOoSpDRHo7XtzoS6y1Z8ikLFkzn1KkEq0MxLOLURCAifUVkh4jEicgz+ex/TERiRGSziPwuIvWcGY+7EhEeeOBhbr+9Z6mnOD71FEx5C8o7umNXOHA7ukLpPCD3ArgIoDKwEIh38HWdzR/oBbOTYOpu+PVXSry4z1vceit88slZDh26lfT0o+zbt8/qkAwXc1oiEBF/4EN0WbNoYLCIROc57F+gvVKqFfAN8Jaz4nF3d999N8OH/w8/v6BSvb5cOahZDh6qheN7D9dC9zuOR3c2y0lWOXcMkeg7hqP5vdj9fb0cXnsNfvnV6kisU748JCf/w5NPNqNFi2h27dpldUiGCznzjqAjEKeU2qOUykC/VVwwaV4p9adSKqdh4ir0E2mf5e9fk2XLOrOxDI9aEpeD/IDjew9fBlwJrOfCmUShwJ3o6aUzgRMOvq4rRAF3wdq6kOXD6wwAunY9zS23ZFKhgiOLWRnuzpmJoBZwINfPB+3bCnI/un36RUTkQRFZJyLrjh/3pNVMJWOz2fj001388Ufpn+8MHgyTP4Q6hf2lS6sXutzEX8DqXNsrAXcD9YEwJ1zXFRrqBWeP/g1PPgNHPfTupqyqVoWhQzOJiRnE2rXP8MQTT3DunI8XbPIBbjFYLCJDgPbA2/ntV0p9opRqr5RqX9WL5/sFBwezatUaJk58tNTnCAyEFk3huWYQcAzHNqIXYADQBJ2yN+faVwW4BV2cLp3zXc88zNZ98O922OWjieA8xbffvslHH01i164S9lY1PI4zE8EhLux+W5t86mWKSC/gOeB6pZSnTUZ0uJo1a1Kv3jOkp0eUafAyIAGypwJrHBaa5o9eeVwPWED+tYfmAV+gF555mtqQPRLeTYWYZDjiw50f+/aFWbOySUu7i7NnY9ixY0fRLzI8kjMTwVqgsYhEiUgQeu7JwtwHiEgbYCo6CZi5a3b79p3kzjuz+f330p+jQQMY+TC07uGwsM7LKTdRHV2TKG8Ns65AMjADSMXzBMKpTBizAIbcBStWWB2QdSIjITV1O59+2o5mzZoxf/58q0MynMBpiUApldPO5Gd008O5SqltIjJeRK63H/Y2+qnyPBHZKCILCzidT2nQoAE33ngL0dGlHzvP6WH7Qhuo4IyFXyHoQeIIYDYXzhiqh077J4BZ6AY3HiizJtg6w57qVkdivejoNO69V9Go0S/YbGbhmbcxJSbc2OnTS9m4sSdlWSWmFIx6GmJOgBqC41P/aeBz9FjEUPRdQo4dwBx0L+E7HHxdF7umMoT+DkPu8N2yFDmCgzvw+uuRPPbYM/To0cPqcIxiMiUmPFRIyOXMnHkZy5aV/hwiMOBq6NsL5/zXrohOAIIeFzica18T9F1Dbydc18V+3QI/LIYVjl6j4YH27VvLpk2/Eh//k9WhGA5iEoEb8/f3Z+NGYe/esjXd7dMHnhoMQ+riuN4FuVUF7gOCgOnoJjY5GqJnFCngd+CYE67vCjVBjYIvysOak3DokO82uqlRAz77zEb9+m+zfft9fP/9PLZv3251WEYZmETgxgIDA1m58m8mTvwO/ZG7bHrYIPB9dPVQR4tEJ4Mw9LjA7jz7U9A1if7HhatLPEmYrvD6zAq4dxh88YXVAVknMFD/e+jQ/3j44TsYOfJeawMyysQkAjcXEhJCpUo9yci4rUwrjgHq1IEeXaBxY8fEdpEKwL3oBWaz0f0McoSjE0U59GyiOCfF4AIqAjK7wtaGkOLI0t8eyN8f3nsvi//8Zy3x8a+QnHyahAQzAdDTmMFiD6CU4sorO3H48HqmTcvCr4zpOykTHtoAh5PQC8AcLRV9V3AE3fqyZa59KehSFMeBm/Ls80A1QqDxMmjTCG66yeporDd5cnVWrMhi9+69hIeHWx2OkYsZLPZwIsK0aV+waNGMMicBgAqB0HEz+H2Gc/oJlEMPINdF9zj+m/MTn8KAe9BTTEOccG0XO5oCy3fBn3t0nwNfd8MNxxg6NIlTpz5GqWzOnvWknqa+y9wReJgdO4axbdtnREaW7TwbN8J3v8GayyHNWYXWMtBNbWLRNYr6oyuWgk4MOcMeMejCb6FOisPZFGCD5pXgziAgGTp1sjoo6x05Es3DDx/im2/mc9VVV1kdjs8zdwReZPr0itx7rx9lrb3XujWMewLebgXlzuGcchBB6PpD3dADxdNzXScnCZwB5gOfoh8XeSIB/GFbMrwwBV57m1L1oPY2IjG0aZNMpUqLsdky8LQPnb7EJAIPc//9w3n00dFUr17JIedrHgHVvgO/r3BsgbocfsBV6PpER9Bv+LlXIYejHyOl2/d5eDmb7BshZTC8uQdOpkGMM2ZoeYgaNeCFFxRJSRNYu7Y11157JW+88YbVYRn5MInAwzRs2JCxYyfSqtV3JCcHkplZtvOJwOgRMHo4RJSuJ07xtEDPKLKhVyLH5tpXF3gQ3e3sK6AMC+gsFwBUg7+Ow52fw8MPU+bZXt4gKSkWpf7hzJlvSE/Xqw7NHYL7MGMEHiolJYWWLRvSqlUCj5a+avUF9qTA6HmQEoleMewMyegWRYfRxel6oCuaAmQCi4Bq6CY4ni4T2ATtesETTSAgBSpX1snXVykFAQFhxMUNYubMeGbP/ooaNWpYHZZPMGMEXigsLIyHH36SoUOHOuycl/iD30IIKkPV0yJFoO8MWgPLgWmc76YWiJ5S2tn+8w5gE2UptWStQKA9rD8N9yyHu4fBx1OsDspaIpCdncLevV+QmLgK3fLO3B1YzdwReDilFDExtxIX9w0VylaJAoA9eyA7Al7eAwdSccSC5oJtBX5Av9Ffh26HmdvX6EVpLez7PXVWEUA2sBrqRcNTPaFBkP50HOrJv1MZKaUTQ3j4tTz44H4eeuhR7r//fqvD8lrmjsCLiQhbt97InXf6EeeA1boNGkDjKvBBa4j8Ht1J2llaACOAGuhppt9yYS2kW9EDzTHAFCDeibE4mz/QGfZVhJH/woj34K6hkOKJzXscJOcR2cGDiylffhtpabM4dy6ezMxMMss6+GWUiEkEXqBHj17ccccQGjSo67BzBiuIDoO26JBiaQAAFJlJREFUzhoryFER3e+4J/oOYQrni9b5oaee3od+I/2CfHrceR4FxNeEpFbwbSKkZ8OJE1ZHZZ0KFWD8eGje/C/WrGnC+PHX0LTppXhzf3J3Yx4NeZH09MOsWXMN27bF0LRp2c9ns+lPbUsT4fU/ICMUXU/IWQ6g7wpOA22AXkB5+74MdKJog35cdRQ9qOwFH2Uqn4PkSfDwCLjhBqujsd769bBiRSDvvPMMdeo8SmzsAZo1a0ZgTqU7o1TMoyEfERxck3nzujJmjB8nTxZ9fFH8/HQi6BoJlX6AwPk4d+C2DvAQesbQJuAD4F/7NYOAtugkkAR8hr5D8NSy1rmc8IPMy2FBKPydCAkJusy1r2rXDkaPzmT//pf544+6dOt2OcOH32d1WF7N3BF4mRMnTvD77z/TrNlXnDixyGHn3bcPkjNgbjqsOIZeD+DMdQfH0FNJD6DXGfRH3wFgv/ZG4Ff0QrQr0FNRvWTgtcJPkLYF5n8D5cpZHY21lIJVq6BGjWC6dHmQ4OB7+fbbvxg2bBhhYWFWh+dRCrsjMInAS9lsWcyffxO//baI2277v/bOPTyq6lrgvzWTEAjhIRAICZEAKhAQKPIs1NtyC0WLWlGBAr4gtZhS5VZakaK0qEWrV2mxj6vYVipFLFqLUC2PilJEQXmIvB+BGITwhgRISGbW/WOdEKCAGGYYkrN/33e+nD3nzDn7zGT22ns9I+u7/vBz8P5CCGcR3cH3TAN+TyypHcARYD62akgE7ic62VQvNgXADmjTBe7OgM/ehY4d4fLImYAqLbNmBZk0KcSKFbNp3/7bhMNhApHIxOgDziUI4s70oqPyEwjE8d57zZgzpw79+h0ikpOn/j2gQRD2pcOivZyaQC6SBDB1UEtMGLyPuZ33ALpi9oObgC5ALuVCYDOWxC5I5aQW0MpyF41eDPIb+Hp/GDei/GP2a1DajTeG6NABDhzox8qV3+CZZ+IpLq7DK6/MQPz6oUQAtyKowoTDYXbv3k0o9DobN/6I3buLadToi9/3ZfjbWnjuMQjfAKRG9tr/QT7wLyzQrCamDurEqdOZfOB3WMW0XkAmld8SVggEIfUy6FYAn7wKPxsPaWmx7ljsmTYNSkrq8OijD5OSMpw33phHr169qF+/fqy7dsnhjMU+JRAIkJKSQlpaNosX/4BhwwLs3BnZe2QGIb069Ez3ZqvRnFc0Ar4LDMfqJL8NTAaWA2WVwhp658QBMzGD80dYuofKShJQAz4vgte3Q85xmFUAe4th40bIrwIG84oyZAjcddchtmwZzd//nsqAAQN49tkHAZsIhUKhGPewcuBWBD4hNzeXl1+eyoABh8jLe8bTrUbm2mURopsLYcyjsC8IfCsy1z77TYGtwAIsb1ESpi66BrMXhLFAtPeBvcCPsEI4YarM9CcokPgnqFEKf3nJ2n5nyxaoVw9SU69g69Zr+fGP/8GcOW/Tvv3pYev+wxmLHaewdu1MbrhhMKNGlXB1BEtFhsPw3HNQEAfbu8KmQmymHk1LlGI2gSWYYIjH8hh1w7KZKhaXcJm3/yIWydyFci+kyswBoACSr4Q+ybDkSbitP/TtG+uOxZ61a2HGDPjlL/uQkXEHixYJmzZtZ/To0b6MSXDGYscpBAJtqV+/HS1apAFvoqoRMT4GAnDffbavCn9eAi89DeFBQJMLv/4ZEeBKb8vHBMJyYBlmZL4GaOGdW4IN/iswdVEaZoxuQ+Utm3mZbXuKYdoGIAAz9wK7oF01WDAHrr/esp76jcxM+PnP4dixuaxbN5fp04MsX16DYcOa0aDBDSxatIymTZvSrFmzWHc15rgVgU+xwV84fPgjhg//DiI7yM6O7D22boU//gm+lgWz9sGa9ZgnT4QN1v9BASYIPgKOYmqjdthKoSHmdvoJJjD2AP2946Ve/6qIiiW4HkKvQNZE6N8ZjhyA48chNdpG/UuYY8cs0Z9IDQYPFq6+ujVvvjmf+Pi6fPrpp7Rq1Yq4uKo5P3aqIYBRo1yFkDOgwOZNm0BKSGl0jFCokLBCIAqD4cYtVsIxfLG8XRQoK8N5zHutGiYYErFBv9h7rSxiuQDzSEqkasQkhICgfZ/xh6D4AGS2gWpxUFoKwaB/XVGLj5s6s0Z1EKnNylWHaZLWmBYtrkKBkpISql1qKqQOHWDSpAq91amGHGdFgCuvvPKEs8/Bg5+xdu1WmjWDmhGOam3e1H588dXN42XPdtAkbGCOBoIN6InYgHgEEwr7vS2BcoEQhwkEr/g8h73XEzFX1MqKF0sRVihOBOJg3RFIDEIoH7QEWnt5qUJhCFYRQ/r5kHBSZHwofJiml0NC9Z0cPnyAkpLarP50N5mtW9GwYSNC4TAaDlfZ1ULVfKozUUEp6hfKJoU7167lyQcf5Omn+7Pn2DQ++GABR45YZGskZ46FhfD449C8KxxsDQvz4OgSTH1TK3L3OSO7sFKZ6zG7wgFMXdUauArTu2+ivJzmIO/vAu9YC6KbfO9isQE4BqndoHM9eH8CtGsD48bY4cJCIhqIWHkoYu/eIubOhd69N9K8eS0+/DCNH/5wFsuWvU/Hjl3ZvXs3JSUlpFWRYA7/qIYcFeK2277Ne++9y4wZ1YADlJZCNCZFi5fAuLHQbhRsrAdFBzBVTROi6+65HxMI67HoZDDDcQYWnZyB2RVCwK+xlQJAA0wgtMOMzpUdxWpP1IHmnaFdTZjzPzBoMAy7w4z/OTnQtKmpk/xGbi688w4MHRpP3bodmTYtyLPPvs/nn6+mceO2rF69moKCArp3737JRjg7G4GjwhQVFbFhwwbatm3Jnj1/5dprs+nQ4SjZ2eGI32vXLmjY0Gy2z06Bt6dDk0cgL4DN2oNYqctoUQjkYG6oOZjbKZjqKgMb8GtidoQcrFBObyx+4TDwDpZBNR1zXa3MapZizOCeDumtIKMYFj0CQ7PhrlvhaKGli/7KV4hIZbzKRm6uuaf27QsJCU144okgS5bsY/36v1GrVkemT5/NkSNHuPfee2Pd1RM4QeCICKWlpUyYMIGWLTPo3TuJHTteZ/DgmQwZEqJnz8jeq6AA1q2DLl1gVxFMeBw2r4LaY2BfCVa8phrQOLL3PYUDlAuFXMpXAwEsFiHV+9sEsz/MpNwoXd3rW19M7VTiva+yzqaLgI3A5VC9HjTOg5znYdDD8PXOEN4Jb822SN9GjSAUKk9j7gf27rUI7zZtrP3IIzU4ciSR1167n6SkdmRnT6FBg3R++9vfArBmzRpSUlIuaioMJwgcUSE3N5dhw+4mO7sv7dvvYdWq+Ywfv5KRI5WWLcsjjiNBTg58/jn06GGCYfQPoUig/vdhyxEIeWoNIlCQ56wcBvKwKml5WERzWeqKAGZUvgyzvB3HVg4DMDXSh8BcTM3UEFsxNMDiHy4xx5TzohRzva0PVIPAOtC/Q49HoGUq7PkAFkyDKVMgpaF9f3l50K0bXGqOONGiuBgSPM+z3/8ekpLiGDGiLYmJmfTt+xadOrVh6tTfUKPGFTzwwFi6d+/OwIEDAdi1axfJyckEI6iHc4LAcVFYunQpo0bdz+TJ95Gc/BlvvPEaTz31MU8+GSI11coxFhZCkyYXrmfOz7dVwxVXwPEwDBkMjVtDq6GQcwRWToLSlliEMdiMvj7lFc8iQQjYhxmc84Hd3t9DJ50jWDnOmtjgWYLFNpStHB7CvJcWYyuPuqdt0QrEiwZlQ4lgK7bVINdDak1gAXw+H+79AzRJglX/gE+WWCR6IGBqloMH4atfjV33LyZLl1qtibZtbcJ0111x9OmTygMPfJOEhOa0avUzsrMHMnHiE8THp3DPPd9n4MCB9OnTp8L3dILAERMWLlzI5MmTefHF/0V1K089NYmJE99k4cIuBAIbePfdQ6xaBSNGmGA4dMhmixUpxhIOW4xCYqKpJR59FDLbQ2YvyC2Ap+6AZtdBwn/DjgIoeBnojEUfh7BBOIXIuLIWYUJhPyYo9p+0f/y0c6thdo8wJhxKKE+glwDcgwmRBZiQqeWdXwtbgZStgI5hK4tL1Q+wCLO5pHjtFcBWqP9dSKkOe2fA4U1w1yRokAALpsC+HfDryVAtAHPmwNGjcNtt9vacHPufqYo1GkpK4K23bJKTmQkFBXFkZSkjR36Lhx+eU+HrxkwQiEhf4FeYZnSKqj5x2vEEYCqWCGAfMFBVt53rmk4QVF62bt3KsmXLTix/x48fw5Qpf2TFit9QVLSNsWOnMnv2eubNy6SkJJ9p0/LZvl0ZPdrev3KlDfbdvFn+oUPmwVTzC2b5paWwapUZotPTbeb5459A75uhRTfY9Bn83/3Q4Q5I6gz5u2DLryH+Bii+ClMJLcbSUTTCZvS5mFG4bKZfjBXpOZeBWDFbwkHvmodO2g5ixupCzp7BNYDNtsPeOTWwGs4JWJT0IUwY1KDcRvFNTNgs995X3Tu/GqZKK4vyPooJkbgveIZocRT7bJK9dtnzfANqBiH8KkghdB4N9eJhyTMgYRj+KNSOg5d+AUmJMO4RO//5582IPchz/Z03z9pdulh7/XqoVas8lXdhoalxLmW1VXLyQNq0eaXC74+JIBCRIGZe6o1pVJcB31XVtSedkw20U9URIjIIuFlVB57ruk4QVC3KUl0ALFq0iA0bNpCVlQXAuHHjWLFiKTNmPM3x47sYOnQsOTk7eeutEZSW7icrayb5+QX85S/tCYUK+OlPN1FaWsr48UHC4aNMmaLEx8Odd9q9Zs40oXHdddaeP98Ggw4dytM5p6VZXp6pU2021jIT9hTAhDHQ73Zo2gF27ITpE6BHFtTJtMC4Zb+C5sOg2hVwaDvkz4DEW6A0DYrzQN/FfgkNsJn9cqyeQm1s9bARq62QgP1aylYoRZgufh82iBdhQqQQWz2cvIL4siRgA28QKwla5ggW8LYyF9kAlpIjeNqWgpURFSzT6+nHG2KCJoypioKUG8xPtqmEvM8gcNpWm/JgwDJDfZx3vz2Y0GvmXX+xd91r7HjgTahWBxrfCkkCG38JtdOh091QIw7++TCktYQbRppMfG4ktOwAQ38ACUF47D7o/FW4/W5ICMDDD8LXesLNN1s3Hn8cuneHXr1MtfO735mQ6dTJVqQzZtj/VWamzfDnzrX9Zs0szccHH9iMPzXV2p98Yq65ycnW3rzZ/hfr1LF2Xh5kZt5Mt26vV/DLPrcgQFWjsgHdgX+e1H4IeOi0c/4JdPf247CEwXKu615zzTXq8Cf5+fm6ZcuWE+3Zs2frq6++eqI9ceJEfeyxx1RVNRwO6eDBA3XYsKF67Ng2LSxcp127ttd+/a7V/fvn6969c7RVq3Tt16+L7tz5ku7Y8YJmZCTrTTd11O3bn9Jt236hqal19JZb2umWLWN08+bRmpycqAMGtNHly7+nM2f219q1q+mQIa11/vz+OnZsJ61ePah33nmV/vWv39Qbb2yqwSCalXWVvvBCd7366joK6O3DM3Ts4621dt04BfS2rDS984Gmis3xtd/wFL3p3pQT7T5ZDbTX3fWtLejXv1dPuw62a0kA7fn9utru5qQT7ba31NT07gnWDqKX/1eC1msZd+L9dVsFtUYDUcSun5Qe0PjacuJ+gQSUACfap+z7cZNzH5fguY8H4s59PC7+yx1/4YUXKvz7AT7Ss4yr0VwR3Ar0VdUsr3070FVVR550zqfeOXlee4t3zt7TrnUPpi0F0+puqGC3GmDCxk+4Z/YH7pn9wYU8c1NVTT7TgUvVtHQKqvo88PyFXkdEPtKzLY2qKO6Z/YF7Zn8QrWeOplloB2ZOK6OJ99oZzxGROMx8tS+KfXI4HA7HaURTECwDrhSRZiJSDUvdNeu0c2YBnimPW4F/abR0VQ6Hw+E4I1FTDalqqYiMxAzCQeAPqrpGRCZgRotZWOHAP4vIZszTetDZrxgRLli9VAlxz+wP3DP7g6g8c6ULKHM4HA5HZKnM+REdDofDEQGcIHA4HA6f4xtBICJ9RWSDiGwWkTGx7k+0EZE/iMhuL1bDF4hIuoi8IyJrRWSNiNwf6z5FGxGpLiJLRWSV98w/j3WfLgYiEhSRFSIyO9Z9uRiIyDYRWS0iK0Uk4qkVfGEjOJ90F1UNEbkWS0QwVVXbxro/FwMRaQw0VtXlIlIL+Bj4ThX/ngWoqaqFIhIP/Bu4X1U/iHHXooqI/AhLylFbVfvFuj/RRkS2AZ1OD7aNFH5ZEXQBNqvqVlU9DrwC3BTjPkUVVX0P88TyDaq6U1WXe/sFWNXhqlBI8qx42QMKvWa8t1Xp2Z2INAG+DUyJdV+qCn4RBGlYWq0y8qjiA4TfEZEMLDfnh7HtSfTx1CQrsdRt81S1qj/zJOAnlKfJ8wMKzBWRj72UOxHFL4LA4SNEJAl4DRilqoe/6PzKjqqGVLUDFr3fRUSqrCpQRPoBu1X141j35SLTU1U7AtcBP/BUvxHDL4LgfNJdOKoAnp78NWCaqlY8Z28lRFUPAu9glZKrKj2AGz2d+StALxF5ObZdij6qusP7uxv4G6bujhh+EQTnk+7CUcnxDKcvAutU9ZlY9+diICLJIlLX26+BOUSsj22vooeqPqSqTVQ1A/sd/0tVh8a4W1FFRGp6zg+ISE2gDxBRb0BfCAJVLQXK0l2sA15V1TWx7VV0EZHpwBKgpYjkicjwWPfpItADuB2bJa70tutj3ako0xh4R0Q+wSY881TVFy6VPqIR8G8RWQUsBeao6tuRvIEv3EcdDofDcXZ8sSJwOBwOx9lxgsDhcDh8jhMEDofD4XOcIHA4HA6f4wSBw+Fw+BwnCBwOh8PnOEHgcDgcPscJAofjAhGRjJPrPojIaBH5WQy75HB8KZwgcDgcDp/jBIHD4XD4HCcIHI7IICftx8esFw5HBXCCwOGIDE29TKAB4FogGOsOORznixMEDkdk2AdMBT7CUgTfISItYtslh+P8cNlHHY4LxCuLOVtVq2xlMEfVxq0IHA6Hw+e4FYHD4XD4HLcicDgcDp/jBIHD4XD4HCcIHA6Hw+c4QeBwOBw+xwkCh8Ph8Dn/D1ssM8imhS2OAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(1, 1)\n", "ax.set_title(\"Hypothesis Tests\")\n", "ax.set_ylabel(\"CLs\")\n", "ax.set_xlabel(\"µ\")\n", "artists = brazil.plot_results(test_mus, results, test_size=0.05, ax=ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Patching the likelihood to replace the BSM components\n", "\n", "A nice thing about being able to specify the entire statistical model using the ubiquitous JSON format is that we can leverage a wide ecosystem of tools to manipulate JSON documents.\n", "\n", "In particular we can use the [JSON-Patch](https://tools.ietf.org/html/rfc6902]) format (a proposed IETF standard) to replace the signal component of the statistical model with a new signal.\n", "\n", "This new signal distribution could for example be the result of a third-party analysis implementation such as Rivet." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "new_signal = [20.0, 10.0]\n", "patch = jsonpatch.JsonPatch(\n", " [\n", " {\"op\": \"replace\", \"path\": \"/channels/0/samples/0/data\", \"value\": new_signal},\n", " ]\n", ")" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'channels': [{'name': 'singlechannel',\n", " 'samples': [{'name': 'signal',\n", " 'data': [20.0, 10.0],\n", " 'modifiers': [{'name': 'mu', 'type': 'normfactor', 'data': None}]},\n", " {'name': 'background',\n", " 'data': [50.0, 65.0],\n", " 'modifiers': [{'name': 'uncorr_bkguncrt',\n", " 'type': 'shapesys',\n", " 'data': [5.0, 3.0]}]}]}]}" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "recast = pyhf.Model(patch.apply(original.spec))\n", "recast.spec" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### The Recasted Result\n" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "test_mus = np.linspace(0, 5)\n", "results = [\n", " pyhf.infer.hypotest(\n", " mu,\n", " data + recast.config.auxdata,\n", " recast,\n", " recast.config.suggested_init(),\n", " recast.config.suggested_bounds(),\n", " return_expected_set=True,\n", " )\n", " for mu in test_mus\n", "]" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd3RU1fbA8e/OpFdaaEGK0ouAIIoVrIigIIpir9h4PgWfT9GH5ffsggX1CVhRFBAbil2KIiICAlIEQ+81CQQSSDm/P+4dTZmWmTuJyezPWqzM3HvuuSdxOXvuKfuIMQallFKRK6qqG6CUUqpqaSBQSqkIp4FAKaUinAYCpZSKcBoIlFIqwmkgUEqpCKeBQKkSRKS5iBgRiXa43i9E5Bon61TKKRoIVJUTkQ0iclaZY9eKyNyquHc4GGPOM8a8FWh5ETlVRHLtfwft4JRb4l/TYNohIvF2XU2CuV7VTI5+61FKOcMY8wOQDNZTCrAeqGWMKazCZqkaSp8I1N+eiPxLRD4oc+wFEXnefj1bRB4XkQUisl9EPhGROiXKXiAiK0Qk2y7bzj7+NtAU+NT+ln1PiVtcISKbRGSPiNxfoq4oEblXRNaKyF4Rmeq+l/1t+x37eLaI/CIiDUq08Ub7dUsRmSMiOXb9U4L8u9QRkYkiskNENovIgyISZZ9rKyJz7XvsFpGJ9mXf2z9X27/zABFpKCJf2m3eKyIzg2mPqr40EKjq4B2gj4jUArD77y8DJpYoczVwPdAIKAResMu2Bt4D7gTSgc+xPvhjjTFXAZuA/saYZGPMUyXqOwVoA5wJjHIHD+AfwADgdKAxkAW8ZJ+7BkgDjgLqArcAeR5+n/8DvgZqA02AsRX/kwAwCcgBjgZ62O26yj73OPAxUAsr2I2zj59m/2xj/84fA/8GVgP1sP5+DwXZHlVNaSBQfxcf299Is0UkG3jZfcIYsx3rm+wl9qE+wB5jzKIS179tjFlujDkI/AcYLCIu4FJghjHmG2NMAfAMkACc5Kc9Dxtj8owxS4GlQGf7+C3A/caYLcaYw1gfmhfbwakAKwC0NMYUGWMWGWP2e6i7AGgGNDbG5BtjKjwWIiLNsD7UhxtjDtl/oxewAqT7Hs2Bhvbv8aOP6gqwglpTY8wRY8z3PsqqGkgDgfq7GGCMqeX+B9xW5vxbwJX26yuBt8uc31zi9UYgBusbbmP7PQDGmGK7bIaf9uwo8foQdn891gf4RyUC1iqgCGhgt+krYLKIbBORp0QkxkPd9wACLLC7rK730xZPmgHxwO4SbXnebgfAXUAi8KuILBORK73UA/AosA2YJSKZIjI8iPaoakwDgaouPgaOFZGOQD+sbpGSjirxuinWt9w9WB9wzdwnRETsslvtQxVNv7sZOK9k0DLGxBtjthpjCowxDxtj2mM9cfTD6rIqxRizwxhzkzGmMXAz8LKItAyiHblA7RLtSDXGHGffY6sxxt1Vdgfwuj3TqNzva4zJMcb80xjTDBgEPCAiJ1ewPaoa00CgqgVjTD4wDXgXWGCM2VSmyJUi0l5EEoFHgGnGmCJgKnC+iJxpfzsfARwG5tnX7cTqYw/UK8CjdtcMIpIuIhfar3uLSCe7S2o/VjAqLluBiFxSYvpmFtaHc7lyvhhj1gPzgadEJMUexG4lIqfY97hURBobK898tn1Zkd2d5R5XcLfnAhE52g6SOVhPOBVqj6reNBCo6uQtoBPlu4Wwj72J1aUTj/UtGGPMaqyupLFYTwj9sQaHj9jXPY71DThbRO4OoA3PA9OBr0XkANaH8Qn2uYZYwWo/VpfRHC9tPR74WURy7br+aYxZF8C9yxqCNRj8O7APmMJfXUM9gUX2Pd4Hhhpj3E9Bo4D37d/5AqAdMAs4gDUW84wx5qcg2qOqKdGNaVR1YXdt/I41ALq/xPHZwDvGmFerqm1KVWf6RKCqBXt+/HBgspeZOEqpIIUtEIjI6yKyS0SWezkv9qKgTHtWw3Hhaouq3kQkCau75WzgwSpujlI1Tti6hkTkNKxZDRONMR09nO+LtTinL1Yf6/PGmBPKllNKKRVeYXsisBel7PNR5EKsIGGMMfOBWiLSKFztUUop5VlVJp3LoPQioC32se1lC4rIUGAoQFJSUre2bdtW+GYrV64kLy+PpKQkgrleKaWqs0WLFu0xxqR7Olctso8aY8YD4wG6d+9uFi5cWOE6OnbsyIoVK2jXzvDdd6+QktINYwxRUTperpSq+URko7dzVfkpuJXSq0Gb8NdqT8c9/PDDAGQcdYjFi09kzJgLOO20U9izZ0+4bqmUUtVCVQaC6cDV9uyhE4EcO3FWWJxyyikALM2EwsIicnI+A5YRFfV7uG6plFLVQti6hkTkPaAXUE9EtmBN+4sBMMa8gpUOuC+QiZXU67pwtQVg+3YrxmxYDlOmwBVXwBlnHGTZstOpU2cE33zTgLvuuku7ipRSESdsgcAYM8TPeQPcHq77l+VyuawXAg1LdkhRzJtvPs3YsVH06tWLbt26VVaTlFLqbyFivv526tSJ3n07QyqYNqXP9e8P48cXU6fOJ1XTOKWUqkIREwgAjjm6KeTANxthU4nclSLQogVs3Ph/TJ58PaNGjUJzMCmlIkVEBYJpk2aDgQXPwkMPey7z4YdvMHnyOHJzcyu1bUopVVWqxToCp5hi+1t+feh5IRhjPQ2UdMstcOjQLrKzXyU5+U6kbAGllKphIuqJYOZMewygADbULx8EAKKiIDkZ/vhjONdd14tx48aVL6SUUjVIRAWCTp2Ot15kwYK18P5H3ssaAytX/sCqVfO8F1JKqRogorqGLr30GuvFfij8HV7+Ek7tCQ0bli/rcsHjjxtSUhZRXHyYqKi4ym2sUkpVkoh6IgBAsHaI7QzHPeQ5CLi5XHDo0Aq+/XYYU6dOraQGKqVU5YqoQPDhhx/S6oxYa1vuePhNILfQ/3WPPfYqd901jPz8/LC3USmlKltEBQKABo3j4TCQBwX7YOR/YaPXnHyWO+6ACRNSiYnRtQVKqZonogLB2LFjmfeevd3tMiAaVv7iPxDUqQOJiWtZt24kv/+uSeqUUjVLRAWClJQUYmPsX3kHkAJyNxx3UmDXP/fc83TufCxr1qwJWxuVUqqyRdSsoWuvvRbT6F2u7/MN7LWOFUbB3D3Qx8egsdtZZxmio1Np3jyAwkopVU1E1BMBQNP0RtbMIbuHiAL43/3w/vv+r61dGwYO3MuGDSPD2USllKpUERUIFi9ezEW9p1iB4JB9MAYOJkFySuD1fPXVS5x4Ylf27dsXjmYqpVSliqhAkJiYSJ06yZBAqd+8aBA06Rl4PUlJsGfPOjZv3ux4G5VSqrJFVCBo27Yt8+Y9SUIP4AhQ9Ne5n/fAoUPerizt6KNhwoT9HHOMKxzNVEqpShVRgQDA5UojpT7W6mL3vvXFMPUeePXVwOsRgczMh1i8eHEYWqmUUpUnogJBUVERnTvfwu4P7QO/2j+joKAbtOtSsfpGjfqA3r1PJycnx8lmKqVUpYqo6aMul4sWLY7iYPxeDmwAdpc4eTJIGy8XenHRRdC3bxdSU1MdbKVSSlWuiHoiAPj22ymcfZ/9Jrv0ubmbYMuWwOtq3hzatv2RQ4dWOtU8pZSqdBEXCKKj06gdB8QAB0uf++FJePnlitVnjOHhh6/g+eefd6qJSilVqSIuEAwceC1v3oQVCA5jDRrbis+FMwZXrD4RmD9/KbNmzXCwlUopVXkiaowA4LjjurN09ZdsLQC2YqWkds8CbQ0761S8zgcegIyM2s41UimlKlHEPRH83//9HyOfrwXtsNYRHC59fvavsHBhxeqMjYXdu6exZcs8Dhw44FRTlVKqUkRcIACoF58GafabraXPZX4MrwSxX31OTjHt2/dizJgxIbdPKaUqU8QFgieeeIIrz9sEC+wDv5YpcD4MDiKnXFoaXHFFIX36dA21iUopVakiLhC0bduW9u1Tob19oGzeuHRYURBc3Zdeaqhb9/tQmqeUUpUu4gLBgAEDeOeds6hzCtZv76FL//vv4cMPyx8PxK+/vs4rr1RwDqpSSlWhiAsEYK0lSIvCykLqYT/67OXw0XQwQWxR/PXXWQwb9g+2bdsWajOVUqpSRNz00a+//pr+/d/iSAFQH2tR2WEgrkShvnBhO2uNQEX17w+XXHISjRs3dqS9SikVbhH3RNC0aVPOPLM1TU4HMuyDxWUKxcHCMuknApWaCrGxP5KXty6EViqlVOWJuEDQtm1bxo+/meMvB5raBz3sQ7DoG3jiyeDuceCAYciQC/j444+DbaZSSlWaiAsEYO1JkATg3p5yffkyhbmwdisUFZU/509iIqxYsZqtWyuQwU4ppapIxI0R5OTkcPTRN5OfD5xvH1wGdC9T8HQ4rgm4gtiEzOWC8eML6dTpqNAaq5RSlSCsTwQi0kdEVotIpojc6+F8UxGZJSK/isgyEekbzvYAJCcnc+WV59G1N+D+nPayr8yirODvIwLbtk0gOzvIwQallKokYQsEIuICXgLOw1q+NURE2pcp9gAw1RjTFbgMCPsEfJfLxZgxo7j6DqAh1mwhL3sVr/0cbrk9+Hs9/vjndOjQjsLCwuArUUqpMAtn11APINMYsw5ARCYDFwIld3ExgHt7rzSgUibfR0enEV8A5AFJWKuLCyn/10iF5ChrnCCYLqLu3Q3NmnWkoKCA6OiI64VTSlUT4fx0ygA2l3i/BTihTJmHgK9F5B9YH8lneapIRIYCQ8Ga/hmqnj0HsGYN0AKohRUIsoF6ZQp2heYZwQUBgG7dIC5uDfHxcf4LK6VUFanqWUNDgDeNMU2AvsDbIlKuTcaY8caY7saY7unp6SHf9KabhnL22SBdgGb2QS8hcUk25OUFf69DhzYxbdpjbN++PfhKlFIqjMIZCLby13AsQBPKJX3mBmAqgDHmJyCe8t/LHXf77f9gxIhEUjtgjRMA5Houu/ZdGHpL8PfauRMGD/4Pb731VvCVKKVUGIUzEPwCtBKRFiISizUYPL1MmU3AmQAi0g4rEOwOY5sAa5/h3NxkknMB945kK7wUbgUdewWXdwigcWMYMyaK224bElwFSikVZmELBMaYQmAY8BWwCmt20AoReURELrCLjQBuEpGlwHvAtcYE+5EbuNtvv51Bg3az9Qn+2qBmlZfC7SH51ODyDrl17VpMTs604CtQSqkwCutUFmPM58DnZY6NKvF6JXByONvgyYABAzhw4CMWmB2scWFtZO+lawhg8R7YmgAZGd7L+DN+/Es0aBDHsGHDgq9EKaXCoKoHi6vEOeecw913d6btaVgb16dgTR894rn8ujfgviB2LStpzpz1fPDBe6FVopRSYRCRk9uLiorIzo6DPVj7EdTGmkK6F2jk4YKecHoja5wg2C6ikSOhQ4eLgm2yUkqFTUQ+EUyePJlevaYz/X5gHdDAPlF2TpNbC8g/JrRxgvh42LVravAVKKVUmERkIDj++OO5++5T6XMd0Ji/1hL4mLi6IBN++y20+77//gK6detMUTApTZVSKkwiMhC0bt2aYcPO5az+WCuL3QHAR5K5TR/A06NDu29iItSrV0RWVgjZ7JRSymEROUZQXFxMVlYUB7ZgjQ3UAgRYAnT1ctFZMKhdaPc97TQ4//wk6tUL+5o5pZQKWEQGgq1bt9K160hS07DWO18CxGJlRjJYQaGshrAlMfR7HziwgB07lpGe3gFXsEmMlFLKQRHZNZSens4TT9zM9dcDJ9kHU7H2LvaSkhpg3s8wZ05o916yBJo2PY65c+eGVpFSSjkkIgNBfHw8N998ORdeAMnugWJ3qgkfCS62zYaJb4d27zZt4LLL0mnSpEloFSmllEMiMhAA7N1bQGYmJLinjLrXD/jaEeECuOqR0O6bkADXX7+DjIyI7JVTSv0NRWwg6N37ah57DHaPwxoXcH9B3+PjolRY4aPrKFDGwOzZz7N27drQK1NKqRBFbCB48slHGTIEOtyGFQjcE3ka+75uzucwvWwO1QrKz4cBA57n2WefDa0ipZRyQMQGgssuu5qzzxaad8T6K6RizaHa6/u63ctg3vzQ7p2QAI88Usw991wXWkVKKeWAiA0Eu3btZu3aBPIygYNYf4kkYBnW7CFvhsAFw0O/f48eALNCr0gppUIUsYHg7rvv5t5785n5DNb2OAAJWEEhx8eF0bDE1/kK+OCD8bzxxhvOVKaUUkGK2EBw2223MXx4Y674D+VzDfkaMAa+eQveeSf0Nnz88R88++xToVeklFIhiNg5jD179iQurhmuDVuY5E4m1xhYjpWFtJX3a7OzYE926G0YPhw6dbo29IqUUioEEftEkJWVxZo1sH0Ff60dcK8l2O7n4ovhpMtDb0NaGmRnfxp6RUopFYKIDQTvv/8+Q4b8yNtjgZ/sg+6uIS87lZW0zKFxgo8/nsfFFw9wpjKllApCxAaCc845h5df7sMD/wHOtA8mYyWfq+/nYgMfPwZOjPMeOGBYv345OTkORRallKqgiA0EzZs359xzu9ClEyTUtQ8KVs4hP2sJEMhLg7Q6fsoF4IILYOLE40hLSwu9MqWUCkLEBoK8vDyWLs1h6VKI+73EiXis7St9JJ8DKO4PrXqF3g4R2LfvKwoLD4demVJKBSFiA8GGDRu46KL/MWUK5EzBSjMB1kb2xcBO/3UszYYjAYwn+DNz5n4yMhqxb9++0CtTSqkKithA0KxZMyZNuofrr4dj7y1xwp1ryNtG9m4F8NY/YPLk0NvSqBH07NmA3Nzc0CtTSqkKithAkJiYyLnnnk7LltAog792JWto/9zhp4IYkK7Qpm3obWndGu699zBNmzYNvTKllKqgiA0EAIsXb2b5cti3ANhvH3QPHAewv3zBGVC3vTNtyc9fz7p1cygoKHCmQqWUClBEB4JLL/0306bBgjf4qysoEWu9dWxgdfy8BQ4cCL0tixfDMcf00i0slVKVLqIDwdSpr3HttXDHC5ROKdEAKxOpPwfg1Vvhq69Cb0u7dnDrrUdxzDHHhF6ZUkpVQEQHgl69zqF5c2jbjNJZl+pirSXwlY4aIAUS+0G3bqG3JSEBBg/eSsOGiaFXppRSFRDRgeC33zJZtkzYsAhYWeJEAtaYwRL/dRzqDrEN/ZcLRFFRMZ9//izbtvnaOFkppZwV0YHg4Ycf5oUXhG8+BZlX4kSgU0gBiuCLXyArgMFlf/bsgYEDH2PSpEmhV6aUUgGK2DTUAI8//jiLF/9CUtI2bl0OO4rsE+n2z10BVJIFk16ERnfD+eeH1p4GDeCppxK46SbdwlIpVXki+omgXbt2tG6dTq1a0CC5xAl3DqFAvuXXhTpXwymnONOm44/Po7j4V2cqU0qpAER0IMjMzGTWrCNs2gS532FtUwlWvqFo4BD+B4wF9h0NRQnOtKmgAF555QnmzJnjTIVKKeVHWAOBiPQRkdUikiki93opM1hEVorIChF5N5ztKWv69On8+9+r2LAB1n4GlEz1UxsrLXWRx0tLOwQTP3FmnMDlgtGjZ/Pee++FXplSSgUgbGMEIuICXgLOBrYAv4jIdGPMyhJlWgH3AScbY7JExN9OAI664oorOOaYWcTHf8Y/3oCxG0uczAAygZgAKsqCTyZA53rQu3dobYqKgnHjijn33NtDq0gppQIUzieCHkCmMWadMeYIMBm4sEyZm4CXjDFZAMaYQIZnHdOgQQPat29OXBw0Silzsi6QS2DjBA2h6Qg4/XRn2lWvHuzb95kzlSmllB/hDAQZwOYS77fYx0pqDbQWkR9FZL6I9PFUkYgMFZGFIrJw924/GwVUwJ49e5gxYwt79sDSL4DfSpx0Dxh/F0BFLtiSAnn+xhMq4NFHX2bMmDHOVaiUUl5U9WBxNFZyh17AEGCCiNQqW8gYM94Y090Y0z09Pb3s6aCtW7eOu+76mD/+gAWzgD9KnHQnnwtgXwKA4r3wzMvO5B0CWL58C8uXL3amMqWU8iGcgWArcFSJ900ov0RrCzDdGFNgjFkPrKF01p+w6tSpE3PmPEzXrjBhPCReUuKk+4kgO8DK9sOc6bB+vTNte/RRePzxEBcmKKVUAMIZCH4BWolICxGJBS4Dppcp8zHW0wAiUg+rq2hdGNtUSkJCAm3btiU+3pqtkx5X4mQsVqqJAqyxAn+aQodH4dhjnWmbtYXl585UppRSPoRt1pAxplBEhgFfAS7gdWPMChF5BFhojJlunztHRFZiTdT8lzHG39bxTraRadN+weWytpw85A5L7hTUdbGeWXZhTSX1xQWrD8ORYoh1KLyOGvUBUVGHmTJlqjMVKqWUB2EdIzDGfG6MaW2MOcYY86h9bJQdBDCW4caY9saYTsYYBzZ+DJyIMHz4C8yeDTt2QM5vWIvI3Jpg7VyWFlh9BRvgngcgP9+Z9qWl5VG/fiDzV5VSKngRnWsIYNGiGWzYcDaJibC9NbyxocTJxlib2hcGWNlhWJcJu3aBE7tOXnklNGtWaUMmSqkIVdWzhqpcy5YdSEqy+uTrx5U56V7etjzAylpBmwecCQJu+/Z9Tk5OjnMVKqVUGREfCL744kd++AGMge/fBJaVOFnP/jkP68nAH4EV+6EokLIBeuihX+jWratzFSqlVBkRHwhefPEVpkwRRGDLakrnG4rG2rKyCAjwS3neMrjuJit5nBN69oRrrjmJoqJAkh4ppVTFRfwYweTJk1m8uA2wj1dfg3N/KFMgHSsr6S6g3FI3D2JAkmD/fqhb139xf046CerXL8blcoVemVJKeRDxTwT16tUjNbU2YE37rFV2ko47Kcb2ACtsDU1uciYIuG3f/iULFvzsXIVKKVVCxAeCn376iSlTrPmes2dDwbuUHg9wB4LNBOy3HHCyJ2f8+CxOPfVUDh486L+wUkpVUMQHgm+++YbnnttKURHk5oLsB46UKOBObXR04HUe+B4GXeJcMOjbF8aNu5jY2Fj/hZVSqoIiPhCMGDGCn3/uT1QU9OsH54wCSk4jrYO1LjqQNBNu9aHFCXD4sDNtbNYMOnZcTUyMLi5TSjkv4gNBUlISqal1EbHel1tL4MLarewPYH+AlbaEtP6QmOhYM1m7djEvvfQ0xjg4N1UppQgiEIhIlIikhqMxVWHdunWMH7+a3bshLw++HUP5BWS1gN1YeYcCtCzbma0r3X7+GYYNu4fMzEznKlVKKQIMBCLyroikikgS1sfkShH5V3ibVjk2bdrEs8/+xNatEBcHBQcpn1Kiif1zW+D1Zn0C11xrLVRzQu/e8OWXfWnVSlNOKKWcFegTQXtjzH5gAPAF0AK4KmytqkQnn3wya9e+SJcu1n7BT48FupQp1Mj+WYEnAtrBCYOgMNA8RX6kpEBS0jyKix2qUCmlbIEGghgRicEKBNONMQ6tm616MTExpKW1+PN93VhwSZlC7plDeypQcQuQHuDk+O4ff2QzdOglHHBqGzSllCLwQDAO2ICVcOF7EWlGwEkX/v7Gj5/JzJnW62nvg+udMgVqYf2lPHUb+bBkBzjZpZ+TA++9N4MVK1Y4V6lSKuIFFAiMMS8YYzKMMX2NNW1lE/BceJtWeSZN+pz5863XcXEQn0jpRWVRWJlIm1GhpBy7p8G/73VunKBTJ/j66zaceOKJzlSolFIEmWvIGGNE5FngA4fbUyWWLFnGvHkpFBfnc+GFsKw1zNxVplBDoKLf7k+CPg0caiTWdpoFBcvJz99CfHwT/xcopVQAQllHULYnvdqKjo4mNjbjz/fl1hKA9USQCyyoQMUZkHMUf65RcML69XDWWWexbNky/4WVUioAoQSCGrOyac6cOTzzzCEKC60tK795EPi9TCH3JjUV/PxduAp++cWBRtpSU2Hnzm3s2VORkWullPLOZ9eQiPyG5w98ARzs9Khaa9asYebMfVx+uTVNs34j2Fs2rY87EOytWN07v4DR+2Hyu0601Mpq+tprxZxyysnOVKiUinj+xgguwvrAL5t78yhgR1haVAVuvPFGzjgjk82bnwJg+Ci4aVGZQilY6SbygHwgPsDKz4YrOjrWVACKiw+SlTWbtLQziY6O+C0llFIh8tc19CyQY4zZWPIf1tTRZ8PfvMohIsTF/TX4mu5pjECwcg6BlW4iUPVgvcN7yuzYAR07XsyUKVOcrVgpFZH8BYIGxpjfyh60jzUPS4uqwOHDhxk58lPmzrXev/4iRL3loWBD+2dFMpEC836A774LpYWl1a8PXbpARkaG/8JKKeWHv0Dga3PGBCcbUpViY2OZNWsZ2+xcQi1aQFJzDwWblPkZoJ3zYOq0EBpYRlQU/PvfufTo0ch/YaWU8sNfIFgoIjeVPSgiNwJle9GrLREhM3Mxgwdb7wcMgNYXeijoHh4vu8bAn4Fw1UPBt8+b9eunsXt3RfqplFKqPH+B4E7gOhGZLSKj7X9zgBuAf4a/eZUnNrYhIn8NvNaPp/x8KffMoRkVrDwJljqcHig/H7p3H8WYMWOcrVgpFXF8BgJjzE5jzEnAw1i5hjYADxtjehpjasysIYA33niTp56y5oyuXAnf3gWsL1MoCYgB9lHhcYKZn8DUqaG30y0+Hm67DQYOPM+5SpVSESmguYfGmFnArDC3pUrt2LGDTZuiMAbS06HjqfBrkoeCdYCdWDOHkgOvf18m/JrNn91PTujfv5gWLbRrSCkVmojfqtJt5MiRTJ7cBxErEFx+M56XzLnHZ3dW8AaD4ex/hNbGsoyBuXMnsnjxYmcrVkpFFA0EJZRcS1A/HijyUMg9Y7Mim9QARMEiB7eudLv99s946KGHnK9YKRUxNBDYNm7cyE03fcmSJdb70SOB9zwUDCGxxuyJ8MILwV9flgg88EAxY8bc7lylSqmIo4HAFh8fz759h8nPt96ffRbEd/JQ0L1bWUMP5/w4VAx5xcG20LP27SE6er6zlSqlIooGAluDBg2YM2ci7j1f+vWDJp7yuiUAqVhjBJ66jnzpA+0dHCx2mzbtHd5++23nK1ZKRQQNBCXExZVO2VDX4PnDPh1YDgSRXnpxNhRVNID48eGHmTz/vK4nUEoFRwNBCQ888AKjR1uv58yBn/+F5w3rm2EtNttU8XvMfRlGjQq+jbSSbIoAAB0XSURBVJ7ccw988skdzlaqlIoYYQ0EItJHRFaLSKaI3Ouj3CARMSLSPZzt8ScmJo6YGCuFUqtW0GMwnjMquR8cyibnDkBhE8hoF2wLPatVC7Kzv3C2UqVUxAhbIBARF/AScB7QHhgiIu09lEvBSlfxc7jaEqgnnniCkSPbAtC4MZx9EdZ4QFmN7Z8HgEMVvElPqHd68G30ZuLEzxg2TGcPKaUqLpxPBD2ATGPMOmPMEWAy4CmV2/8BT2Jt91LlSq0liAL2eyjkHjAG2FrxeyzaC3sruNOZP9u25bFo0fcUOT0AoZSq8cIZCDIo3Xmyhb86VQAQkeOAo4wxPtO4ichQEVkoIgvDmW1zzpw5XHzxj2yy+/5ffAD4xEvhplgJOjw9Mfjxy3PwyP8F10Zvrr8eJk7sjcvl8C44Sqkar8oGi0UkChgDjPBX1hgz3hjT3RjTPT093V/xoKWlpdGsWQOMnXX06qug9mleCjcFCgFPu5n5YU6A7n2CbKQXIrBnz8cY42mLaaWU8i6cgWAr1t7Gbk0o3ZGSAnQEZovIBuBEYHpVDhh36dKFiRP/TbNm1vtTToF2x3kp7H62+Y3y6ar9aQ9FbYNroy/Tp2+madNG5OXlOV+5UqrGCmcg+AVoJSItRCQWuAyY7j5pjMkxxtQzxjQ3xjQH5gMXGGMWhrFNfpUcIygqgtrZeE453QBrH+PvgOyK32f+eli9Org2etOwIXTvnk5OTo6zFSularSwBQJjTCEwDPgKWAVMNcasEJFHROSCcN03VAMHPsDLL1uv9+6FGf/Ban1Z0fyVbiKIAePVE+HpZ4JrozedO8N99xXRsGEQ+S+UUhEroP0IgmWM+Rz4vMwxj8upjDG9wtmWQHXt2h1jrNw96elw0wiY4G0iTjOsbSu3YHVyVcSZcL7D6wkADh1axdq1P9CsWU+io8P6n1cpVUPoyuIynn12LIMGpQHWAOwlfcFV20th9wjIhiBudBRsrxXEdX4sWgQtW57GDz/84HzlSqkaSQOBByXHCQ5kQ531eB4Qdg8Y76LiCeiAuYtg7twgGuhD+/Zw661NaNmypbMVK6VqLA0EZYwfP56zz179Zzrq2bNh95t4HjCugzV99Ojg7rX9O3hlfHDXepOQAIMHb6V+fV1PoJQKjAaCMlq1akW/fq0oLLTen3YanP8AnnMOCVb30H4gmM/dvnCFwwvLAIqKDDNmPMP69eudr1wpVeNoICijd+/e/Pe/l5Bsb0xfrx4cfyzeh9UzsLqGfgviZrXgtzAk1sjNhcGDn+O1115zvnKlVI2jgcCDuLgMikvsJHZkPZDppbB7OGFmcPf6fjZMnBjctd6kpcHTT0cxfPjNzlaslKqRNBCUkZeXR/v2w5k8+a9jX0wGme3lAncm0izgcMXvd3A9fPEtpQKPE7p2LeLIkTnOVqqUqpE0EJSRkJDA1VcPolWrv46NGAHNb/RyQRJgdyOxPYgbngW9H4Qoh/9LGAMTJjzLhx9+6GzFSqkaRwOBB6NHj+H44/96n5EBLev7uMDdPRTECmNiYH5WENf5IQKTJi1h4sQ3na9cKVWjaCDwICamLgUF8X++P3QI8uYDO7xc0Nz+GcTWlQDrf4Kbb3e+e+jJJ4uZMOF6ZytVStU4Ggg8GDZsGEOGFPz5vrgY5k4E1nq5wL2w7Nggb+iCfBccOBDk9V6kpcHevd42VFBKKYsGAg/69u3LlVc2+/N9cjK89DZwkpcLGmL9JYMZIwDoAPVusD64nfa//01l6NChzleslKoxNBB40LdvX4YO7VnqWPsmkOhtLUEMUA9YQnB5h4BlObCvovsfByAr6xAbN/5GsdP9TkqpGkMDgReFhfU5cuSv96tXQ/x3eM8pdBRWGgpv6w383W8lXDoQtgf7VOHF9dfDCy90JsrpaUlKqRpDPx08WL58OZ07P8u8eX8d27ABsn/E82b2EFomUoAG0PAka7aP03bvfp+DB3WzGqWUZxoIPGjWrBkPPHAlLVr8deyMM+DmVwFvKandA8bbsfYyrqg6kHMWpDcI4lo/Zs7cR/36Ddm0KchpTUqpGk0DgQcpKSn86193/Ll3MUBMDByT4uOiulj5iIqwNqoJwoFC+G4lOL3lcKtW0L9/Y/8FlVIRSQOBFwUFtdmzp/Sxld8A3vYPiMJaWOYCgh303QKPD4Offgryei8aNYJbb91Ko0apzlaslKoRNBB4cfHFN/Hww6WPrV0J0Zt9XNQc64mgmY8yvjSGWgOtvYedZsxhfvppLJs3+/oFlFKRSAOBFyNGjOCqq+qVOvbgg9DhFh8XufMTZRLUjmVEQXZnyI33X7SicnPhzDMf5LnnnnO+cqVUtaaBwIt+/fpx1lmtSh0TgRZJPi5qBMQDnwIrgrxxIbz3LWzcGOT1XiQnw4MPGm6//TJnK1ZKVXsaCLzIz89ny5aUUmsJcnJg1RvAGi8XRQEtsWYNBbmegCL4eix88UWQ1/twyingcn3rfMVKqWpNA4EXX375JRde+DUld3tMSoJ9G4GDPi50P0T8AQSzmDcOuAEuvjqIawMwY8Y4xo0bF57KlVLVkgYCL3r06MELL1xOw4Z/HYuOhtffArr6uLCl/TMP2BncvU1jWOht4VqIPv98I48++hBFRcEMYiilaiINBF40btyYwYP7l0sElxwN9WJ9XJgEuBeFectW6o+BqVPhu++CvN6Hm26CL74YjMvlcr5ypVS1pIHAh71701hTZjxg3jzIewXf21K2sX8Gu0pYYMNCmOvwegKwMpzm5LyPMfpEoJSyaCDwYdiwpxgzpvSfKCbGGivA1+pfd/dQEHsYu5krodsNwV/vy6+/bufEE49l584g+66UUjWKBgIfHnvscUaN6lTq2PHHw7WjgFo+LszAGvRdBuzxUc6XOPhqh7X3sNMSEyEnZztbtgSZC0MpVaNoIPDhxBNPpEeP3uWOt0j0c6ELOBprmunPwd9/+fdw5TVQUOC/bEW0aAETJhTQpUs7ZytWSlVLGgh8MMawaFE8q1eXPj7jNZBJfi5ubf9c7bOUb8kQ09D5LSwBiopy2bnzQw4e9DUXVikVCTQQ+HHHHa/y4Yeljx3dHNKbAL66bdzjBPuBfUHevCXkXwS1vaW+DsHhw9Ct2408+OCDzleulKpWNBD4ICJ8/fW33Hln6QGBgQNhwLWAr01kUrC2r4Tgp5ECOw/D3I3g9Bf3uDi44ILDnHZaa/+FlVI1mgYCPzp37kzjxj3KHe9WC//ppt3TSP8IoQHZ8OAN8NVXIdThxZAh0KbNIucrVkpVKxoI/MjOzubjj6VUqgmAj1+BqFf9XOxON9HJZynfakH0uXBstxDq8CEz823+979ndXN7pSJYWAOBiPQRkdUikiki93o4P1xEVorIMhH5TkSCzeQfNkVFRTz00FcsKvPF+bTToOU5+M4ndBTWNNL1PsoEoOAEWBeG1NQA8+blcdttw5k/f354bqCU+tsLWyAQERfwEnAe0B4YIiLtyxT7FehujDkWmAY8Fa72BKtu3bqsXbuYQYNKHz/xRBh4Ib7/gu5ppCuAX0Jrx0e/wFxvu6OFoFcveOONhvTseYLzlSulqoVwPhH0ADKNMeuMMUeAycCFJQsYY2YZY9w97fOxNnv82zn66K7Ex5dvWod4YJufi1tirTBeEFobfv8Mxr7k/AKzmBho3nwHu3d/5GzFSqlqI5yBIAMouS/iFvuYNzcAHrPwi8hQEVkoIgt3797tYBMDs2rVKl58MZasrNLHp7wKMhHfu5G5p5HuxppKGqw+cM5/rM1xwmHkyOGMHDkyPJUrpf7W/haDxSJyJdAdeNrTeWPMeGNMd2NM9/T09MptHLB3714++mgzZTMyDBgAp96G72mkaUAd+/XSEBpRG+bkhnC9Hzt2bGbLlpXhu4FS6m8rnIFgK9ZwqVsT+1gpInIWcD9wgTEmhDRt4dOzZ082bvyMTmVm/7RsCX1Pwv9fsa39czG+F6H5sXkt/PNea/9hp40YAffdl+a/oFKqxglnIPgFaCUiLUQkFrgMmF6ygIh0BcZhBYFdYWxLSFwuF7Vrn4Cnr/51c8H1m58K3EPkcUB+CA0xsGY15Z5MnCACu3ZNZsOGpbppjVIRJmyBwBhTCAwDvgJWAVONMStE5BERucAu9jSQDLwvIktEZLqX6qrcjBmzeeSRpHKDtXO+haKP8J1yOgNrbwIBEkJoRGNwDYfmrfwXDcaqVUdo1eo4pk//2/5nUEqFQXQ4KzfGfA58XubYqBKvzwrn/Z20Z88etm2LYf9+Su1aNmAAFHeDd3N8XCxAN6y/xCqsKaVxQTRC4GAxfLoNzkml3O5poWrVCq66KoHOnTs4W7FS6m/tbzFYXB3ccMMNfPvtw+U+fOvWhV7HBFDBsVjrCqZg7VMQggmPwn9G+S9XUS4XXH31QRITQ8idrZSqdjQQVEBqavmcQwD5GyBupp+L44GOWE8HIab3KegIGSeFZ9MagLlzH+fZZ58NT+VKqb8dDQQVMHr0p/z3v+UHjNesgcIFgL8Mod2xZg3tsP8FqxMsa+47u0UoPvlkFQ8+eD9VsWZDKVX5NBBUQGxsAvHx5TcH6NcPhr8KJPmpoAngXgbxa2ht2XYQnnsfMjNDq8eTIUNg8uT61K2r00mVigQaCCrg/vvvZ8yYweWOx8VBj3oeLihLgOPt178R2lf6Avj8zfCkp05MhMTEjWzd+jLZ2dnO30Ap9beigaCCUlKO99g3/8diSHgP3+kmwEpJ7cJKUR3KXz8eim+EDpeEUIcfd911Lz17nkhhYWH4bqKUqnIaCCpoyJDxjB1b/viRIxCXB/hb9ZuAFQxW4XvtQSDqwrubIVzrv7p1O0z//g10rwKlajgNBBV0wgmn07x5+UUAp58O947Byi3kTzfgCPAKkOWnrB9/LIeLL4Ot5ZJ3hK5nTzj//J8pLt7ufOVKqb8NDQQV9OSTT3LVVZ5z93euBdGBfDtvAtTFCgJLQmxQXSisa21GHw7GHGbSpOt57bXXwnMDpVSV00AQhKSk7hw5Uv7477+BjAF2+qlAAHcsWUhog8YpkDsE9tfxXzRYb789k9GjH9McRErVUBoIKmjfvn106/Y/Pv20/Lmjj4Zju+A7LbVbJ6y//kFgTejtems1LFwYej2e3HUXvPnmUbhcrvDcQClVpTQQVFCdOnUYOvRqWrYsfy41FZ75L3RtHUBF7kFjgK/xP9vIjyXT4P4HIMdXzqMgpaXBoUNz2LVrOvv27XP+BkqpKqWBIAijR7/Caaed6vV8v1TgjwAqcq8pqE9gTxG+nA61boCElBDr8cIY6NNnCJdddikmXLktlFJVQgNBkFyuK1niZaD358kQ9QHWzCBfmgDtgExC28YSIAV21Yc3NsBBf6kugiACZ599iIED05Bw7ZeplKoSGgiCdN99H/PYYy6Pc/ivvhpuegyIDaCic+2fk4HvQm/XlC/gkkvDM530vPOgXbuPyMqaqWsLlKpBNBAE6amnnmbatDvxNH6akQGXdIEGgew5UAs4DSsJ3VwgxH3aTFNwtYe4xNDq8a6YceMu5rjjOpMTjgEJpVSl00AQpA4dOtCz531ERcV7PF94BOp9Q2DrBE7CCggAXxLSvsakQu75MD2MKYKSk7OIidnKwXD0QSmlKp0GghDs3VvA888fxerV5c/FxkLxDogNZKP5aOB8rACwjsAGmv2YtALufgB2hJLu2ouOHeGpp7IoLJzkfOVKqUqngSAEiYmJzJ27j/Xry58TgbEvwCWXB1hZK6CN/foLQp5OWlwES5bCagfWKHgiAitWjOS66wayYsWK8NxEKVUpNBCEIDU1lc2bt3HppZ6nkrpccFEGRAfald4HKzNpqv0zFLWg6J+woWmI9fiQl1fIZ599yuzZX4fvJkqpsNNAEKLY2FgyMv7BoUOez69eBIXPAZsCqKw21sDxRmAtoY0VAMTAO5vg07mwcWOIdXlQpw68+WYRvXr95HzlSqlKo4HAAWPHLuW661x4StvfpQsMvBwk0FxA7oHj9+1/Ic7SLDwEzz0F/5sQWj3eJCXB7t3vM2PGMB566CFdbKZUNaSBwAG9e5/BlVf2pqCg/LmEBLjjRjjjaAL7hh8D9AXygZVY4wWhfLYmQPE1sPlc2OdvgVsI3n33JSZMeE5TUChVDWkgcMAZZ5zB44+/R1KS56mkAJcnQ9xbwJ4AKmyNtdE9wC9Y6wtCUQ+2FcOIxfD8SxCOz+obb4SXXsohJ2c0xhhdcKZUNaKBwCExMXXZsKEXm7yMBSTGQJ1CSMwLsMK+QEf79XeEvm8BsGEjfDIdfghDl74I1KoFmzY9zh139OS6667TtNVKVRMaCBxy4MABbr/9e6ZN83y+YUN4+0347/ngCiRVTxQwEOvpAAIbbPanIZhh8F1jyA/jZ/TBgz9z8OCC8N1AKeUoqW6De927dzcLg0m8f+edeM0S55Cc/fuJitpKUaHvPBEbd0HWQaxZQv4YrLQT+UA6VvpqB3K+JRVD9H5o2hRcYfg6YIC42IbExBxNVFSU7mWglBO6dIHnngvqUhFZZIzp7umcPhE4KC01leTktkRF1SbXx4piVwHEFBDYILBgpamOA3YDWwAPs5Mq6mA+7M+FIw7U5YkAhw/v4NclP7Nkya8hz4RVSoVPdFU3oNIEGUUrSoCRw27ltdfGM2lSMXU8TBstKoKCIrh3JSwNdLFZPvAaVjDYA1wAdAixsflQOwVGtIKDv0Lv3hATE2KdZSxbVkReXhRNrh1K/fpDOHToEElJSc7eRCkVksjpGqpEWVlZzJs3iwYNHiU3d7HXcjty4Pr/Ql43oFEAFR8EXsXa9B6svQz6A6FmGt0EvA533AUDLwixLh/mzz+OF1/czpw5P3DMMceE70ZKqXK0a6iS1a5dm/PPv4hjj/2StWubesxFBHAkB2LXwjGBDtwmAbcAPez3q4BZhL4CuSlwDUxJh6XZVqK6cEz4SUtbTKdOe4iLmwVAfn6+8zdRSlWYPhGEUWFhIe3atSYlZRujRx/G08ZeeXnWorOvdsDzn0BePNA8gMq3Ah9hdRO1AjoDLQHvSxkCIkcgcTx0PxYe+k9odfmSlNSbgQNXcNdd9zBixIjw3UgpBegTQZWJjo5m+vQZfPDBTGJj65ObW/6bdkKC9fPs+tBwIdSeR2Df8DOAW4FzgPXANGA0MIPAFq15YWLg4Bkwryk8sxpW7YFXXoGdO4Ov05OdO2fRs+cuEhLGsnXr/9i9ezNvvPGGPiUoVQX0iaCSFBRkcdllZ7B06W+MG1fkcWezvDw4dAjWAKOXwd4FwLH4HwPIBj4ANpc4VhdrUVqoXfFrQd6F6x+CISfDgRwrq2pKSoj1ljFjRgLPPJPH/PmfcsIJ/di2bRvx8fHU8TTarpSqMF9PBBoIKtEnn3zChg1rGTQols2bn+GLLzbSowckJ5cvO3ESvPEq9HgYlgnk78Z6Uqjn4wa7gYXAUqxZRi6gLdZAdAxW11EdKr4OIRdIggbxUPdHWP0lvP8R1E6GXbusxHOhTgQyBtasgTZthOTkzrzwAnz44Wr27t1OXFwamzZtom7dujrjSKkgVVkgEJE+wPNYH0mvGmOeKHM+DpgIdAP2ApcaYzb4qrM6B4KS/vhjNa1bt+W22xpzySU7yM0tZsEC6NYN0tKsMjt3QoMGcLgI7n8clvwIcf+GQwJswPpAb+ahcoM1E2i5/a9kWosorP0OGgG9sQIDBD6ReDuwFVzHQ9sUOPAuHNgEk96FBBfMnm1NQT355Ir8NcpbswY2bIBzz40lNfVE7rhjK7m5Ln766Rvi4o5iwoQJ1K1bl0GDBgFQXFxMVJT2dCrlTZUEAhFxYfVynI21DOoXYIgxZmWJMrcBxxpjbhGRy4CBxphLfdVbUwIBwOLFi2nRogUpKTFMnz6eQYNG8MYb59Khww6WL1/Bu+8WcsMNVnqKFStg3jwYfJm10Pi/98CBg9D5X7D1CKybCkWxwJlYAWIp1irkY7D2N1iCdWE21tNCSVFYwSMaa+FaIlawaA3EAvvsOhOwBqOj7TJ1gB12nW2gdjzkvwOJyXDh3dA8DaY+AQ0bwL3/slJrvPWW9fuce65167lzoW5daNfOer92LaSmQnq69X73bkhMhN9/t7rNuneH2NhEbrmlmJYt6/PSS9eRlNSGU0/9Nz17dmfChPHExNTlzjvvpEePHlxxxRUAvP7667Rr146ePXsC8OOPP9KkSROaNWtm33ctderUoXbt2hhjyMrKIjExkfj4eIwxFBYW4nK5NNioastXIMAYE5Z/QE/gqxLv7wPuK1PmK6Cn/Toaa5hTfNXbrVs3UxMdPnzYLF682OTm5hpjjJk6dbJp3vwos3TpdLN37zfmueduMYD5/vs7TWbm3ebuu08ygJk162KzfPml5thjaxvATPz0BPPWzONNci2XAcw/JzY3t7zX1MQlRxkEc9GE+qbfK/VMVCwGwTTrF2/S2lhlASMu/nxd5f8kxPMh/hMRn+ejoqLCer2/8y6Xy+f56OjokK4P9by/+4favlD/fv7++/g7H85/nu4tIubVV18N+jMGWGi8fK6G84ngYqCPMeZG+/1VwAnGmGElyiy3y2yx36+1y+wpU9dQYKj9tg3gYbv4gNQjpDk11ZL+zpFBf+fIEMrv3MwYk+7pRLVIMWGMGQ+MD7UeEVlovD0a1VD6O0cG/Z0jQ7h+53B2eG4Fjirxvol9zGMZEYkG0rAGjZVSSlWScAaCX4BWItJCRGKBy4DpZcpMB66xX18MzDTh6qtSSinlUdi6howxhSIyDGtA2AW8boxZISKPYA1aTMfKp/m2iGRizU25LFztsYXcvVQN6e8cGfR3jgxh+Z2r3YIypZRSztJJ0UopFeE0ECilVISLmEAgIn1EZLWIZIrIvVXdnnATkddFZJe9ViMiiMhRIjJLRFaKyAoR+WdVtyncRCReRBaIyFL7d364qttUGUTEJSK/ishnVd2WyiAiG0TkNxFZIiKOp1aIiDGCQNJd1DQichpWuriJxpiOVd2eyiAijYBGxpjFIpICLAIG1PD/zgIkGWNyRSQGmAv80xgzv4qbFlYiMhzoDqQaY/pVdXvCTUQ2AN3LLrZ1SqQ8EfQAMo0x64wxR4DJwIVV3KawMsZ8jzUTK2IYY7YbYxbbrw9g7eGWUbWtCi87e0Cu/TbG/lejv92JSBPgfKyNW5UDIiUQZFA6W/8WavgHRKQTkeZAV+Dnqm1J+NndJO60gt8YY2r67/wccA9QXNUNqUQG+FpEFtkpdxwVKYFARRARScbaqudOY8z+qm5PuBljiowxXbBW7/cQkRrbFSgi/YBdxphFVd2WSnaKMeY44Dzgdrvr1zGREggCSXehagC7n/wDYJIx5sOqbk9lMsZkA7OAPlXdljA6GbjA7jOfDJwhIu9UbZPCzxiz1f65C2u38h5O1h8pgSCQdBeqmrMHTl8DVhljxlR1eyqDiKSLSC37dQLWhIjfq7ZV4WOMuc8Y08QY0xzr/+OZxpgrq7hZYSUiSfbkB0QkCWunckdnA0ZEIDDGFALudBergKnGmBVV26rwEpH3gJ+ANiKyRURuqOo2VYKTgauwviUusf/1repGhVkjYJaILMP6wvONMSYiplRGkAbAXBFZCiwAZhhjvnTyBhExfVQppZR3EfFEoJRSyjsNBEopFeE0ECilVITTQKCUUhFOA4FSSkU4DQRKKRXhNBAopVSE00CgVIhEpHnJfR9E5G4ReagKm6RUhWggUEqpCKeBQCmlIpwGAqWcISVex1RZK5QKggYCpZzRzM4EGgWcBriqukFKBUoDgVLO2AtMBBZipQi+WkSOqdomKRUYzT6qVIjsbTE/M8bU2J3BVM2mTwRKKRXh9IlAKaUinD4RKKVUhNNAoJRSEU4DgVJKRTgNBEopFeE0ECilVIT7fybuMNCOTpD3AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig, ax = plt.subplots(1, 1)\n", "ax.set_title(\"Hypothesis Tests\")\n", "ax.set_ylabel(\"CLs\")\n", "ax.set_xlabel(\"µ\")\n", "artists = brazil.plot_results(test_mus, results, test_size=0.05, ax=ax)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.5" } }, "nbformat": 4, "nbformat_minor": 2 }