{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# X-ray Tutorial Module\n",
"\n",
"**Lecturer:** Dipankar Bhattacharya
\n",
"**Jupyter Notebook Author:** Dipankar Bhattacharya & Cameron Hummels\n",
"\n",
"This is a Jupyter notebook lesson taken from the GROWTH Winter School 2018. For other lessons and their accompanying lectures, please see: http://growth.caltech.edu/growth-astro-school-2018-resources.html\n",
"\n",
"## Objective\n",
"Learn how to analyze x-ray data\n",
"\n",
"## Key steps\n",
"- Search for gamma-ray bursts in CZTI data\n",
"\n",
"## Required dependencies\n",
"\n",
"See GROWTH school webpage for detailed instructions on how to install these modules and packages. Nominally, you should be able to install the python modules with `pip install `. The external astromatic packages are easiest installed using package managers (e.g., `rpm`, `apt-get`).\n",
"\n",
"### Python modules\n",
"* python 3\n",
"* astropy\n",
"* numpy\n",
"* scipy\n",
"* matplotlib\n",
"\n",
"### External packages\n",
"None."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Load modules"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"from astropy.io import fits\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import math\n",
"import scipy as sc\n",
"import os"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Load event file with astropy"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Provide the path to the czti events file in the same format as given in the example below."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"filename = os.path.join('data', 'AS1A02_005T01_9000000948_06884cztM0_level2_common_clean.evt')\n",
"\n",
"dataHDU = fits.open(filename)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Combine data from all quadrants"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"CZTI contains four quadrants, here we obtain the time value for each quadrant and combine the data for them."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"times = np.concatenate((dataHDU[1].data['Time'], dataHDU[2].data['Time'],\n",
" dataHDU[3].data['Time'], dataHDU[4].data['Time']))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Make histogram and plot light curve"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Setting binsize for each interval and accumulate counts in each of them using histogram and plotting it. "
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"scrolled": true
},
"outputs": [
{
"data": {
"text/plain": [
"Text(0,0.5,'$counts\\\\ s^{-1}$')"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEKCAYAAAAB0GKPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3XeYlOX1//H3oS0sHVkQKaKAIigqrIqisaBRScHEEjUqscQUjYnGJFiSmG400cSv/uxGTIyKJdGoUbEFRYqrIkVUFkQBKUvv/fz+eO7dndk6s0zb2c/ruvbap9zzzGFnmDPPXc3dERERqU+zbAcgIiKNgxKGiIgkRAlDREQSooQhIiIJUcIQEZGEKGGIiEhClDBERCQhShgiIpIQJQwREUlIi2wHkEpdu3b1vn37ZjsMEZFG5Z133lnh7kX1lcurhNG3b19KSkqyHYaISKNiZp8mUk5VUiIikhAlDBERSYgShoiIJCSjCcPMOpnZE2b2oZnNMbMjzayLmU0ws7nhd+dQ1szsNjMrNbMZZjY0k7GKiEi8TN9h/BV4wd0HAgcDc4CxwCvuPgB4JewDnAoMCD+XAndmOFYREYmRsYRhZh2BLwD3A7j7NndfA4wGxoVi44DTwvZo4CGPTAE6mVmPTMUrIiLxMnmHsQ9QBvzNzN4zs/vMrC3Q3d2XhDJLge5huyewMObxi8IxERHJgkwmjBbAUOBOdz8U2Ehl9RMAHq0Xm9SasWZ2qZmVmFlJWVlZyoIVkcbphVlLmDxvZbbDyEuZTBiLgEXuPjXsP0GUQJaVVzWF38vD+cVA75jH9wrH4rj7Pe5e7O7FRUX1DlQUkTz33X+8yzn3Tsl2GHkpYwnD3ZcCC81s/3BoJPAB8AwwJhwbAzwdtp8BLgi9pYYDa2OqrkREJMMyPTXID4CHzawVMB+4kChpjTezi4FPgbNC2eeBUUApsCmUFRGRLMlownD36UBxDadG1lDWgcvSHpSIiCREI71FJC/t2pVU/xlJgBKGiOSluyfOz3YIeUcJQ0Ty0vSFq7MdQt5RwhCROv2/10u57ZW52Q4jaTt3weZtO1m4alO2Q8kbShgiUqebXviIWyZ8nO0wkvbynGVc+OA0jrnpNSWNFFHCEJG80rK5VWxPmb8KgLFPzchWOHlFCUNE8t6kUk0VkgpKGCKSVwyrdmzPDq2zEEn+UcIQkbxSWNC82rEB3dtlIZL8o4QhInnniH26ZDuEvKSEISJ5b9uOXXzw+TqiGYekoZQwRCQhf5/yabZDaJBBPTow9ZNVjLrtDd79bE22w2nUlDBEJCE///esbIfQIFt37KzYXrd5exYjafyUMEQkr/XpUlixfeX46WzfuSuL0TRuShgiTcyrHy5LeOTzlu076y+Uw/p0KeTIfntU7K/ZtJ2la7cA8XceVX3n7yVc+dj0tMfX2GR6ASURyaJ1W7Zz0YMlDN6rA89dcUy95bfuyL9v4/dMnF/RHvPghYfx9PTP2aNtK67/8qCKMi/OXgbArd84JCsx5iolDJEmYvzbC7n9tVIAZn++jm07drFg5UZe+3A5l35hX8yqD3jLB/t0jR+DEdt4/62/vV2xfd7wvenbtW2ddx5NnaqkRJqInz89i89iqqI2bN3BVeOn84f/fsii1ZsTusYz73+ervDS5qRB3RMqd+vLH/PxsvXsf/0L1c7t2uVc9dh0Xpi1JNXhNSpKGCJNRE3VSwtWRAnkmJte4+UPltV7jSnzG8ecTIWtotHew/dNfADf2s3b+eKtE+OOlY/bWLt5O0+9t5gfPPJe6oJshJQwRBqpnbu5BOnzM5ewYeuOiv1/TvuMNZu27W5YOaFPl0KeuXwEN3x1cMKPef2jsmrHrv3XzLj97Tub9sA/JQyRRmj252vZ7/r/8vT0xTWef2jyAr7yf2+yLdxVPD+zelXK9VXGVbz64XKG/mZCymPNliG9OlHYaveaaR+ZtjBF0eQHJQyRHLd1x07mLlsfd2xe2UZ27nIm1FKN9IunZzNz8VrWbYkGqn3/4XcTeq5djsYpSK2UMERy3O+fm8NJt07kkxUbWbxmc9x8SFu272J8ycKKxFC2fivXxVSjFP/2ZcY+mdziQVeNfz81gWfBui3bWbNpOzuqVNd9Yb+iamUvHNE3Q1HlDyUMkQyZV7aB4//0OpNKVyT1uMmhoXncWwsYceOrPPVuZTXUy3OW8dMnZvDAm58A8OyMz3l46mdxj3/07eSqVf4T0xPqxdlLk3pstq3cELXB7NGuIO74PecPq1a2bP3WhK45c9FaLvtn5R3aIb9+iVUb86OtJ1kZTRhmtsDMZprZdDMrCce6mNkEM5sbfncOx83MbjOzUjObYWZDMxmrSKrNWLSGT1Zs5PGS6h/g23bs4pS/TOS+N+YD8PIHy7hk3NtxYwIefGsBAB8vW8+cJeviHr9iQ/Th13Y36+zL/fXlucxctJa7/zcvJdfLtH27to3bb92y+hoZ5UYO7Fbntb5y+5u8Na+yd9iaTdsTTjb5Jht3GMe7+yHuXhz2xwKvuPsA4JWwD3AqMCD8XArcmfFIRdJgxYZtXPPUzLjpOTZs3cGHS9fz2+fmAHDdv2fy8pzlLFtb/YPp7onzufP1+A/yF2Yt44Q/vc7D0z6rVr4hbn35Y75xz2TmlW1MyfVyyePfPZJHvj28Yn/0oT2TXpHv0bc/49UPl3H5P9+tqA4EeHvBKlbn8d1HLlRJjQbGhe1xwGkxxx/yyBSgk5n1yEaAIomYvnANy9Ztqbfcm6UreGTaZ/xjyqes3riNp6cvrtbQvGHLjloeXbMVG7Yyf8VG3l+Yuum7N23LzxHPh/XtEje/FMCfzjw4qWv8bdICLnqwhGdnLOHFWUtxd1Zs2MqZd03myvG1z0E1vmRhRfVhY5TpqUEceMnMHLjb3e8Burt7eZ+/pUD5sMyeQOy9+6JwrGkPtZScsX7LdgpaNKdVi2Zs3raT0+6YxAE9OvDfH9Y8R9PMRfHVSHdPnM/dE6MqqKF9OlUc//3zc9iYpx/WueroAV3pu0chC1YmNiljrJ88MYOfPDGDg3tHr+Ebc2tvo/rpE5UdEC46eh/cnW07d1HQovYqs1yS6TuMo919KFF102Vm9oXYkx51/0hqZIyZXWpmJWZWUlZWfeCNSLoc8ftXuOCBqWzZvpOTbv0fQLW2hVg7dkV3EVXr14G4hX3uCUmk3G6Oz6vm8CrLl54yeM/UPkGOGnPk3ow6KH3/1vK7u8JWzfn1fz7gwr9N45Bfv8StEz5m6dotvPfZ6oqyd/1vHlu27+Taf83koBteqhgvk4iXZi/l75MXpDj6xGT0DsPdF4ffy83sX8DhwDIz6+HuS0KV0/JQfDHQO+bhvcKxqte8B7gHoLi4uGkPw5SM2rRtJ1Pmr2Lgz6vPPVS2fitTP1nJlw7qgZmxeM3miiqeru0LmL8isbaB+96cT+nyDSmL+ccn7cepB+3JibdUToFx4Yi+vNDIekM1xK9GH1jn+QtH7MMvn5ld6/lhe3fmnU9X13q+3PotO3hgUmW1019fmctfX5kbV2b5+q0M+sULFV8Gtu7YSasWdX9/X7tpO7e/Npd734iuff6RfeuNJdUydodhZm3NrH35NvBFYBbwDDAmFBsDPB22nwEuCL2lhgNrY6quRHLWzl3OL56exeX/fI+7J85n245djLjxVZ54Z1HS13pocmqXRf3ByAH079a+Yv/J7x3JYX0Tn28pn512SM9azx3WtzMPjDkspc+X7J3jy3OWVSSLbMlklVR34E0zex+YBjzn7i8ANwInmdlc4MSwD/A8MB8oBe4Fvp/BWEUarN+1z/PfWdE39hv/+2FcL5pcM2zvLjRrVjmteWzvoVjtC5reSghPfu+oiu2fnDyQjoUtG3SdqlWAiSpZsIqz7prM+wvXMO6tBWzJgWnXM/YucPf5QLWuCO6+EhhZw3EHLstAaHHmlW2gTcvm7NWpTaafWvJU8W9fjtuf9smqLEVSv32LqrevAJCfS2XUqUfH6l1tf3bKQP74wodJXaehr/d9b3zCtAWrOPfeKWzctpNvFPeOO79s3Ra6J9kdeHflQrfanDLyz//jtDsmpe367s7cZes1X08jN/HjzHawOPGAugeXpdP/+2bTGzP77s9PqvFLY33tDLvD3dm0LepO/cHn6yralcp7zFWd7mTd5szfuSphBFu27+Qnj0dz6CxP4yjOD5as46RbJ3LrhI/T9hySfqlY7/nYKvMb3XVe9ekryh0zoPpcSOkyaewJ/O5rlQ3Eow5qesOfurRtVePxFs3Sc6u1bssOrv3XTAb94kUemfYZo257Iy3Ps7uUMIL3PlvD4zGNkiULVrG+AXXP7s4dr5UyY1HNA6jWbY6+QSTS20Jy15oUfLtLdkXUbx3Vl757FNZ6/uj+XeP22zWw3aFnpzb07hw9zzEDutZTOr9VvaM49aA9Oefw3uzfPeo4cMnR+1ScG3vqwAY/z4gbX62YSv2ap2bWWObzNYmtiphOShhB7AygAGfcNZmDbniJf0z5lDteK+X/XpnL3f+bx0dL19dyhUjZ+q3c/OJHfO8f77J0bf2jfqVx+lKGv3V3aduKG746mME9OwJRMrj3gmKmXVvZ/Nerc3wVyo+/uF9GY2xsylflK79paNki2ohNvOXjI8oTR7f2rfnD14dwahjP0SmmIfy7x/bjwQur96QqTy4D92xf7VwyJldZ7fCkKqsDZoISRj2u//csbn7xI/484WP+8N8POfkvE3nn09Xc/OKHzCvbwJpN2/j5v2dV9JUvTzuL12xm+B9eabKTlOW7gphvnrHVF8neNcSqq7rj1AOjD6iNYYW8P54+pNpa1VefvD83fv2ghgfQxHz32H6ce0SfigRR2KoFL1/1Bf7fedXbbA7oEf9hf/rQXnxpSA9GV+mK26JZ9Y/UP5wevSa2O2+OWpSt38o590zhjDvfYsna9N+BNL2+cjVYt2U75943NeHyp9/5FgB3vDaPZhb1p169aRu3nzuU2Z+vjSu7YesOitoX1HQZyRNXnjiAnz8dDfi66fQh/OSJxNafKO/48KMTB7BP17YM2qsDAL/8yiB+9Z8P4sq2aB7/QVT+7ThW13YFnH14HwpaNmP24uojznt3Uc+/WPsWteP3X4tPsLFjVGI1q/Jh37tLIXecW39ngKcvG8HOUHtR0KIZN50+hB6dWnP+/dMaGHW8N0vLKu48Pl62gR4d0/sa6w4DePb9ho8HLO+4sGHrDrZs38lFD5bEnT/+T6/zzfum8OnKmkf2nvKXiVz9eMMWrBn31gIeTdHspNJwow+tfcDXecP7cMe5Q/nWUX3jjhe1L+D9hdGXi96dCxl9SE/26tSGub87lQtH7FPDlRL3tUN7cf2XB1U/XsfAtPoMDsmsVXN9ZNSl6k3ikF4d4/bPOqw3Q3pWzht2y1nJTXpY1ZR5lV12D9jNKq9E6NUnNUtSzv58HQf/6qUaz00qXVnrUpofLl3foBHAAL98ZjZja2kgk8zp0LqyHru8obl968qb9y8N6cENXx1csX/3+cN4a+wJNA+fLofGTDzYMkc/kH8dptXo161dliPJbUP37lyxfeS+e9RRMvL1ob24ejfamh6rYW2VdMrNd2eGLU1gSur6lK3fytY6JhBbvGYzqzdu45x7pwAwtcpgHnfn8ZKFLEhwjqGqj121cRsrN2zl2RmfsyvVs9VJwr6wXxH3nD+M7x/Xv9YyJw/eMy4x7E7ddoc2LenRsTU/OyW5HjrfPbYffz37kIr9r4e7pDY1VHVBZeNuExy/l5TWLZvTumX02t50xpAaX9vYLxMAl58wgOH7RqPBLzu+X4Ofu6CORaJSRW0YwIuz0j/x2t8mLeBvkxbEHes79rmK7eXrt1bUfb929XHs3LUrrj51zaZtTJy7gkN7d+InT7zP2FMPqDj3++fncO8bnzCoRwc+WLKO569oV1EfLunVsU38dBHNmxlfHLxnxmYTbd2yOZOvqTZRQr2qdgG96YwhXDPqgLi7JUmPZjV0brCQigfv1bHauZp0LmzJ6k2VXbv/fObB1d6L6aA7DKjoqphNV8UsunL8n17nxFsmUrq8sgvvn1/6mCseeY9LxpUwZf4qHpq8oOJc+YRkH4SptTWKPDO6dyhg6rXJf1jvjpoau1OhRfNm6pyRAzq1acn1Xzqg2iDOqm1HE396fNx+3xqmzE8HJQzg/845NNshMKl0ZbVj5VNQb962kwWh0bw8GZT37ZbsaW5W51rRqfC94/rFNZhfMXIAZwzrFVdXng69OrehZXOr1m1U0u+SY/bl4N6VX2K/f1w/pv/ypLgy7avcCRakccqSWKqSynHj317IT5+s3k0zdr59yT/jv3MkG7Zu54SB8WMtBu7ZIenlRBti36J2fPSbU2usPpHkdApVRcfv37D5wH5apX1q9q9Ojts/qGfHil5s6aaE0UCH9O7E9FrWT963a9uEF8ipT9VkUX7dZes0IDCfNXRK7FRSskhMYavmHFhH28O+Re0ouf5E9qhlfqpkta0y5csXB3VPy6DAmqhKqoH+8o1Daj1XWNA41ueV7CqfDVUfy7nr64f2pG2r5jSv4wN56rUj+fslh9d5na7tCpL6UC8fKDigjm7M3Ttkvs1JCaMGe8XMg/+b0YM5q7hXxb4ZfPKHUXGNTGNPHcjrVx/HCQOjW879u6uHUlNTPklf+fQeR/bryuC9OvD1ob1qfczNZwzh8uP707tL7RMK5pra7qrz1Z/POph3f3FSnXdb7Vu3pKBFcl8SH/n2cB6+5Ihaz3drX8DPvzyI/zu39vbVM4f1rvVcuihhVHF0/648dPERFY1Ih/bpzG9Oq5zquUthq2rfFL57bD/6dm1LUbv4jD+oR92Jo6b3YGzfeGk87r2gmCnXjKyYwqN/t3Y8d8UxDO1T2Tg9eK8OFeMdAA7s2ZGrT96/YgBfLutSGFWn7PKmNcbHzJJOBpccvS9m0RiZ2hzZbw9G9K99JmAz4+Kj92Hgnrn15VMJo4p/XHIE/bu1Y+86ppGuTae20Ruk/FbxqH6VIz3/+e3q3yZ+OHK/apPFHdK7U7Vykvtat2zOnjWs0BbruSuO4ZY6qjJzWYvmzejZqU3FeAGp3dUn789Hvzk17eMi2oUBgIUZXD5Xjd61uGjEPvzu+Tl0S6Ke8NvH7EuXwlacN3xvzj2iD50KW3Hfm1FvpqP6daWofUG9s9f2aUTVEyJSs2RX5uvXrS2T569MarT22Yf1xoAzhtVe7ZlqusOIEbvgzNmH92HmDSfTrX38t8a6Gq66tivgO8f2o21BC3p1LqxW5TR57Alx+7FzCNV0/Ux1lROR7PrJyQO567xhddYwVF0Qq1NhK75zbL+MjPAupzuM4JUfH0vLGuayh/g57n9wQu1zBNUndorq93/5RTq2aVnnbLPXjTqA+9/8hMnzV7IprOsrIvmnY5uWnBLWPKnN+O8cybIUzHu3O3SHEfQrakefWtotmjczPvj1ycz//SjGVJmmuqFqWiznipED4vaP6t+V+791GIWtorz+2KXDqz3mpSu/kJJ4RCS3DdqrA8cPbNjgv1TRHUaCyj+00+Ubxb256qRomuO7zhtG58Lqt5l7tKventK/SNNNS82S7d0jUh8ljAxItgGsvlvTcv27tdNoXKnVl4b0YNm6LazauI2/T/k02+FIHlDC2A0PXngYqzZuq/V8YasWXHXSfgyLmSjutnMOZeLHZWmbdVSkXMc2LbnypP3468tzsx2K5ImMt2GYWXMze8/Mng37+5jZVDMrNbPHzKxVOF4Q9kvD+b6ZjrU+x+3frc6RvBC1S8QO0PnqwXvxpzMPrugNVd4pqkXzht8pvHb1cQ1+rIhIorLR6P1DYE7M/h+BW929P7AauDgcvxhYHY7fGsrlnWP368bJg7vX2ZjeoU3ljeDNZwyJu2MB2CdDc+FL0+buFdPsS9OU0YRhZr2ALwH3hX0DTgCeCEXGAaeF7dFhn3B+pGVqSsYM2rNja+4+v5j96ljf4jejD+TLQ3rQu0sbzizuza9i1ocWyZTl67fyzqertQRwE5bpO4y/AD8FypeE2wNY4+47wv4ioHyynZ7AQoBwfm0oH8fMLjWzEjMrKSsrS2fsWTOif1duP3doRa+XrqG31JkZHOEpcnSYYHHDth31lJR8lbGEYWZfBpa7+zupvK673+Puxe5eXFRUlMpL56w9O7bmo9+ewneObfiC8SLJOrqOyfKkachkL6kRwFfNbBTQGugA/BXoZGYtwl1EL2BxKL8Y6A0sMrMWQEeg+jqmTZT62ItIpmXsDsPdr3H3Xu7eFzgbeNXdvwm8BpwRio0Bng7bz4R9wvlX3ZvY3MpJOPuw3hlb11dEmqZc+IT5GXCVmZUStVHcH47fD+wRjl8FjM1SfI3CjacP4c7zhmY7jCZj8vyV7NT3F2lisjJwz91fB14P2/OBausbuvsW4MyMBiaSoFYtmrFiQ+NYV71TmGamU2Fq1pSWpksjvUUaoJkZIwd2z3YYCTljWC8KWzVn1EE9sh2KNHJKGCJ5rm1BC84szvz6z5J/cqENQ0REGoGkE4aZ/SwdgYiISG6rt0rKzMbH7gKHkKfzOomISO0SacNY5+6XlO+Y2Z1pjEdERHJUIlVSv6uyf106AhERkdxWb8Jw908AzKxr2F+V7qBERCT3JNPo/UDaohARSYG35q0AYPvOXfWUlIZIJmHk3VoU+WrjVk0/LU3Tdf+aBcDna7ZkOZL8lEzC0MQ5Oa51y2gG21mfr81yJCKSj3SHkUcO2LMDAK2aazymiKReMp8s16QtChERyXkJJwx3n5XOQEREJLclVXdhZhPNrEPY/q6Z/cjMNGeyiEgTkGxld0d3X2dmw4BvA52Be1MfloiI5JpkpzffHtbXvgD4o7uPN7OSNMQlIiI5JtmEcRvwPtCayiVT26U0IhERyUlJJQx3f8jMngJ2uvtmM+sPTE5PaCIikkuSXnHP3TfEbJcCF6Y0IhERyUka4SUiIglRwhARkYRoHIaIiCRE4zBERCQhySaMquMwfgkMTuSBZtbazKaZ2ftmNtvMfhWO72NmU82s1MweK79jMbOCsF8azvdNMlYREUmhZBNG+TiMLwP/CccSHYexFTjB3Q8GDgFOMbPhwB+BW929P7AauDiUvxhYHY7fGsqJiEiWJJswJgJHAAcmOw7DI+VdcluGHwdOAJ4Ix8cBp4Xt0WGfcH6kmWmKdRGRLEk2YTzl7hvcfTNUjMO4O9EHm1lzM5sOLAcmAPOANe5evkTcIqBn2O4JLAzPswNYC+xRwzUvNbMSMyspKytL8p8jIiKJSihhmNlZZnYj0N7MDjCz2Mfdk+iTuftOdz8E6AUcDgxMKtqar3mPuxe7e3FRUdHuXk5ERGqR6B3GJOADol5RtwClZvaumT0LbE72Sd19DfAacCTQKTSkQ5RIFoftxUBvgHC+I7Ay2ecSEZHUSGhqEHdfDDxkZvPcfRKAme0B9AU+TOQaZlYEbHf3NWbWBjiJqCH7NeAM4FFgDPB0eMgzYX9yOP+qu2tdcRGRLEl28sFJMdsrSe4bfw9gnJk1J7qzGe/uz5rZB8CjZvZb4D3g/lD+fuDvZlYKrALOTiZWERFJraQShpmdAHwTWAPMAmYAs9x9a32PdfcZwKE1HJ9P1J5R9fgW4Mxk4hMRkfRJdrbaB4AfEXWJHULUBXYw0D/FcYmISI5JNmF86u7/DtuPpzoYERHJXUkP3DOzKzWATkSk6Un2DmMQcBDwMzN7B5gOTHd33W2IiOS5ZHtJnQ4QusWWJ4/hqHpKRCTvJdtLqoRo8sGZ4ec5d38wDXGJiEiOSbYN46tEdxOtgO8AC8xsQaqDEhGR3JNsldTnwOfACwBmdgDRKGwREclzyS7RunfsvrvPAfZLaUQiIpKTku0l9YiZ9QE+IWrDWAscmPKoREQawAzcoah9QbZDyUtJ3WG4+1FEM8heSLSeRSnwlTTEJSKStDOH9QLg9GE96ykpDZFsL6kuwJVAN6Lpzh9y99XpCExEJFkFLZrTpW0rClo0z3YoeSnZXlKPAuuJ1vMuBN40s2oTB4qISP5Jtg2jyN1vCtvPmtljwD+JBu+JiEgeS/YOY5WZHVS+E6YmL0xtSCIikouSvcP4PvCkmb1B1EtqMDAv5VGJiEjOSfYOoz1wFNGyqt2IJh88J9VBiYhI7kn2DuMhYKi7jwcws67AicCzqQ5MRERyS7J3GFvC0qkAuPsK4NepDUlEctHWHbsA+Hjp+ixHItmSbMKYb2anVjnWKlXBiEjuOqhnRwA2btuZ5UgkW5KtkvoB8F8zOx+Yghq9RZqMtgXJflxIvkl2apAlwDDgSaAImAGcm4a4REQkxyT9lcHddxIljCdTH46IiOSqZNswRESkicpYwjCz3mb2mpl9YGazzeyH4XgXM5tgZnPD787huJnZbWZWamYzzGxopmIVEZHqMnmHsQP4sbsPIpp76jIzGwSMBV5x9wHAK2Ef4FRgQPi5FLgzg7GKiEgVGUsY7r7E3d8N2+uBOUBPYDQwLhQbB5wWtkcTTZ/u7j4F6GRmPTIVr4iIxMtKG4aZ9QUOBaYC3UPvK4ClQPew3RNYGPOwReGYiIhkQcYThpm1I+ph9SN3Xxd7zt0d8CSvd6mZlZhZSVlZWQojFRGRWBlNGGbWkihZPOzuT4XDy8qrmsLv5eH4YqLlYMv1CsfiuPs97l7s7sVFRUXpC15EpInLZC8pA+4H5rj7LTGnngHGhO0xwNMxxy8IvaWGA2tjqq5ERCTDMjnWfwRwPjDTzKaHY9cCNwLjzexi4FPgrHDueWAUUApsAi7MYKwiIlJFxhKGu78JWC2nR9ZQ3oHL0hqUiIgkTCO9RUQkIUoYIiKSECUMERFJiBKGiIgkRAlDREQSooQhIiIJUcIQEZGEKGGIiEhClDBERCQhShgiIpIQJQwREUmIEoaIiCRECUNERBKihCEiIglRwhARkYQoYYiISEKUMEREJCFKGCIikhAlDBERSYgShoiIJEQJQ0REEqKEISIiCVHCEBGRhChhiIhIQjKWMMzsATNbbmazYo51MbMJZjY3/O4cjpuZ3WZmpWZrkyTtAAANMElEQVQ2w8yGZipOERGpWSbvMB4ETqlybCzwirsPAF4J+wCnAgPCz6XAnRmKUUREapGxhOHuE4FVVQ6PBsaF7XHAaTHHH/LIFKCTmfXITKQi0li98+lqduzcle0w8la22zC6u/uSsL0U6B62ewILY8otCsdERGrlwLotO7IdRt7KdsKo4O5O9HonxcwuNbMSMyspKytLQ2Qi0lg0bwYjB3bLdhh5K9sJY1l5VVP4vTwcXwz0jinXKxyrxt3vcfdidy8uKipKa7AiIk1ZthPGM8CYsD0GeDrm+AWht9RwYG1M1ZWIiGRBi0w9kZk9AhwHdDWzRcAvgRuB8WZ2MfApcFYo/jwwCigFNgEXZipOERGpWcYShrufU8upkTWUdeCy9EYkIiLJyHaVlIiINBJKGCIikhAlDBERSYgShoiIJEQJQ0REEqKEISIiCVHCEBGRhChhiIhIQpQwREQkIUoYIiKSECUMERFJiBKGiIgkRAlDREQSooQhIiIJUcIQEZGEKGGIiEhClDBERCQhShgiIpIQJQwREUmIEoaIiCRECUNERBKihCEiIglRwhARkYQoYYiISEJyOmGY2Slm9pGZlZrZ2GzHIyLSlOVswjCz5sAdwKnAIOAcMxuU3ahERJqunE0YwOFAqbvPd/dtwKPA6CzHJCLSZLXIdgB16AksjNlfBByRpVgaldtfm8fDUz/Ldhh57bOVm9i/e/tsh5EV1z41k8JWzbMdRo0+XbWJ7u1bZzuMvJXLCSMhZnYpcClAnz59shxNdnUqbMnFR+/DkrWbsx1K3hvQvR1nFvfKdhgZdUCP9nyjuDfrt27Pdii1GtC9HV87tGm9Lplk7p7tGGpkZkcCN7j7yWH/GgB3/0NtjykuLvaSkpIMRSgikh/M7B13L66vXC63YbwNDDCzfcysFXA28EyWYxIRabJytkrK3XeY2eXAi0Bz4AF3n53lsEREmqycTRgA7v488Hy24xARkdyukhIRkRyihCEiIglRwhARkYQoYYiISEKUMEREJCE5O3CvIcysDPg07HYFVmQxnGQ1pngbU6ygeNOtMcXbmGKFzMW7t7sX1VcorxJGLDMrSWTkYq5oTPE2plhB8aZbY4q3McUKuRevqqRERCQhShgiIpKQfE4Y92Q7gCQ1pngbU6ygeNOtMcXbmGKFHIs3b9swREQktfL5DkNERFIo5xOGmV1pZrPNbJaZPWJmrc1spJm9a2bTzexNM+sfyhaY2WNmVmpmU82sb8x1rgnHPzKzk2OOnxKOlZrZ2DTFe0KId5aZjTOzFqGsmdlt4blnmNnQmOuMMbO54WdMzPFhZjYzPOY2M7PdiPWHIabZZvajcKyLmU0IzzvBzDrnQqx1xHtm2N9lZsVVyif1moep9KeG44+FafVTHe/NZvZh+Bv+y8w65Xi8vwmxTjezl8xsr3A85967Med+bGZuZl1zIdba4jWzG8xscfjbTjezUTHls/peqJW75+wP0TKtnwBtwv544FvAx8AB4dj3gQdjtu8K22cDj4XtQcD7QAGwDzCPaMr05mF7X6BVKDMoxfFeRLTU7H7h2K+Bi8P2KOC/gAHDganheBdgfvjdOWx3DuemhbIWHntqA2M9EJgFFBLNWvwy0B+4CRgbyowF/pjtWOuJ9wBgf+B1oDimfNKveXi9zg7bdwHfS0O8XwRahDJ/jPn75mq8HWLKXEHl/6+ce++Gc72JlkT4FOia7Vjr+dveAFxdQ/msvhfq+sn5OwyiP3Abi76VFwKfAw50COc7hmMAo4FxYfsJYGT4ZjAaeNTdt7r7J0ApcHj4KXX3+e6+DXg0lE1lvBuBbe7+cTg/ATg9Jt6HPDIF6GRmPYCTgQnuvsrdV4fHnBLOdXD3KR69Mx4CTmtgnAcQ/cfZ5O47gP8BXyf+bzgu5vrZjLXWeN19jrt/VEP5pF7z8D45geh9U/Xfnsp4Xwr7AFOA8vVEczXedTFl2hL93yuPN9feuwC3Aj+NiTPbsdYXb02y/V6oVU4nDHdfDPwJ+AxYAqx195eAS4DnzWwRcD5wY3hIT6Jv84QXZi2wR+zxYFE4VtvxlMVLlPlbWGV1yRlE34Li4k0wrp5hOxXxzgKOMbM9zKyQ6FtYb6C7uy8JZZYC3XMg1rrirU2y8e4BrIn5MM9EvBcRfXvN6XjN7HdmthD4JvCLBsab9veumY0GFrv7+1XK5/J79/JQTfaAherfBsSb6vdCrXI6YYQ/4Gii27K9gLZmdh5wJTDK3XsBfwNuyV6UlWqKl+g/2dnArWY2DVgP7MxakIG7zyGqEnkJeAGYTpW4wrernOhGl0i8uaS+eM3sOmAH8HBWAqyirnjd/Tp3700U6+VZCzKoJdYC4FoqE1rOqONveyfQDziE6Avmn7MVY6JyOmEAJwKfuHuZu28HngJGAAe7+9RQ5jHgqLC9mMpvRS2IqqtWxh4PeoVjtR1PZbxHuftkdz/G3Q8HJhK1wcTFm2Bci6mswtjteN39fncf5u5fAFaHuJaFW3LC7+W5EGsd8dYm2XhXElVVtKhyPOXxmtm3gC8D3wxJOafjjfEwldWpufbenU30Re19M1sQrv+ume2Z7Vhrifdjd1/m7jvdfRdwL1GVEw2IN+Xvhbr+ITn7AxxB9EYoJGp8Ggf8gGgyrvJG5IuBJ8P2ZcQ3eo8P24OJb0SaT9SA1CJs70NlI9LgNMTbLZwvAF4BTgj7XyK+MW6aVzbGfULUENc5bHfxmhvjRu1GvOVx9QE+BDoBNxPf6H1TLsRaW7wx514nvtE76dcceJz4hsPvpzpe4BTgA6CoStlcjXdAzPkfAE/kwvuhrvdCOL6AykbvnHzvAj1izl9J1G6RE++FWv8d6bhoSgOEX4U/8Czg7+GP+DVgZviDvQ7sG8q2Dn+40vCC7xtzneuIehh8REyPB6L6xI/DuevSFO/NwJzw3D+KKWvAHeG5ZxL/gXdR+HeUAhfGHC8O154H3E4YfNnAWN8g+vB6HxgZju1BlNTmEvXm6JILsdYR79eI6my3AsuAFxv6mhP1PpkW/h2PAwVpiLeUqB56evi5K8fjfTK8hjOA/wA9c+H9UFOsVc4voDJh5Op79+8hnhnAM8QnkKy+F2r70UhvERFJSK63YYiISI5QwhARkYQoYYiISEKUMEREJCFKGCIiaWJ1TDYZU6a3mb1mZh+EyQl/GHOuxsk1zewkM3snTJD4jpmdkEAsNU7amgwlDBGRFDCz48zswSqHJwAHuvsQou6w19Tw0B3Aj919ENHYj8vMbFA4N4to3qmJVR6zAviKux8EjCHqolufO4kGix4C/BO4PoHHxFHCEKlBmPenfNrppVWmoX4rTc95qJndX8f5IjN7IR3PLenhtU82GVtmibu/G7bXE43Z6hn2a5xc093fc/fySVdnE014WgBgZl80s8nhbuJxM2tX/jBqnrQ1YS3qLyLS9Lj7SqI5fjCzG4AN7v6nND/ttcBv64ipzMyWmNkId5+U5lgk9S4imsqoVhat4XMoMLWuclWcDrzr7lstWgPkeuBEd99oZj8DriJaVqF80tbNwDqiu5mk6A5DJElmtiH8Ps7M/mdmT5vZfDO70cy+aWbTQt1yv1CuyMyeNLO3w8+IGq7ZHhjiYaZVMzs25o7mvXAe4N9EE1pKjrBo4aLpwH3AV2Net9iFj+qdbDLcCTxJNBvEutrKVXnMYKKJDb8TDg0nWk9jUohpDLB3OLfbk7bqDkNk9xxMtN7BKqJ5fu5z98NDw+UPgB8BfwVudfc3zawP0QI/B1S5TvlUFOWuBi5z90nhg2RLOF5CHXchknnufgREXyCAb7n7t2LPx0w2OdJrmVrDzFoSJYuH3f2pRJ7XzHoB/wIucPd55YeJ1vg4p0rZIqpP2pp09abuMER2z9uhDnor0fw+L4XjM4G+YftE4Pbwje8ZoENMvXK5HkBZzP4k4BYzu4JoYr3yevDlRFPnSyNgZqcQLej0VXffVEsZA+4H5rh7Qt/6Q2+r54gmCo2tnpwCjLDKZavbmtl+RDPkdgzbACcRtZUkRQlDZPdsjdneFbO/i8o7+GbAcHc/JPz0dPcNVa6zmWjyTADc/UaiOuc2RNULA8Op1qGsNA63A+2BCaGa6i4AM9vLzJ4PZUYQLQR3glVZ39vMvmbRQnFHAs+Z2YvhMZcTLfP6i5jHdHP3MqJlrB8xsxnAZGBg+MLxbeBJM3s/PN9Pkv3HqEpKJP1eIqqeuhnAzA5x9+lVyswBfly+Y2b93H0mMNPMDgMGEs2CvB/xVVeSI9z9daLZs2OP1TjWIfRwGhW23ySqSqqp3L+Iqp2qHv8ttVRNuvurwGGJXisZusMQSb8rgOIweOsD4LtVC7j7h0RVBuWN2z8ys1nhW+J2KpdyPZ6oKkIk4zS9uUiOMLMrgfXufl8dZSYCo919deYiE4noDkMkd9xJfJtInNDT5RYlC8kW3WGIiEhCdIchIiIJUcIQEZGEKGGIiEhClDBERCQhShgiIpKQ/w9AYs3SgesgFwAAAABJRU5ErkJggg==\n",
"text/plain": [
"