diff --git a/Lightcurve/JWST_lightcurve.ipynb b/Lightcurve/JWST_lightcurve.ipynb new file mode 100644 index 0000000..37c0ec9 --- /dev/null +++ b/Lightcurve/JWST_lightcurve.ipynb @@ -0,0 +1,195 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "fec77a1f-d8c8-420e-aaeb-d38ef21a2455", + "metadata": {}, + "source": [ + "# Guide to import JWST Time-series Observations as stingray Lightcurve objects\n", + "\n", + "Author: Ping Hei Ng" + ] + }, + { + "cell_type": "markdown", + "id": "cb75e18b-bb5e-48f9-86ae-71377395d5dd", + "metadata": {}, + "source": [ + "## Downloads and imports\n", + "\n", + "For quick queries, we can use the astroquery package in python. To download, please follow the instructions on:\n", + "\n", + "https://astroquery.readthedocs.io/en/latest/#installation. " + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "821982e9-cfdf-4a8a-9dc7-8ac175354faa", + "metadata": {}, + "outputs": [], + "source": [ + "# Import packages\n", + "from astroquery.mast import Mast\n", + "from astroquery.mast import Observations\n", + "from astropy.io import fits\n", + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "from astropy.time import Time\n", + "from stingray import Lightcurve" + ] + }, + { + "cell_type": "markdown", + "id": "2e41958c-1d6e-4270-9feb-89134c898e17", + "metadata": {}, + "source": [ + "### Make a query and download data" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "5ef5bd2c-4990-44a0-bb5c-67cc86977295", + "metadata": {}, + "outputs": [], + "source": [ + "# Make a query with desired criteria\n", + "#e.g. instrument name (put * after instrument name to use wildcards), name of Principal Investigator, proposal ID\n", + "obs_list = Observations.query_criteria(instrument_name=\"NIRCam*\", proposal_pi='Pirzkal, Norbert', proposal_id='1076')\n", + "\n", + "# Examine the columns if needed\n", + "#print(obs_list.columns)\n", + "# Check the name of your desired target with its corresponding observation id if needed\n", + "#print(obs_list['obsid', 'target_name'])\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "3656987c-53c6-4322-96dd-3e2104dedc49", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "INFO: Found cached file ./mastDownload/JWST/jw01076-o014_t002_nircam_f444w-grismr-subgrism256/jw01076-o014_t002_nircam_f444w-grismr-subgrism256_whtlt.ecsv with expected size 826. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076-o014_t002_nircam_f444w-grismr-subgrism256/jw01076-o014_t002_nircam_f444w-grismr-subgrism256_x1dints.fits with expected size 3965760. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_02102_00001-seg001_nrcalong/jw01076014001_02102_00001-seg001_nrcalong_cal.fits with expected size 123840. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_02102_00001-seg001_nrcalong/jw01076014001_02102_00001-seg001_nrcalong_rate.fits with expected size 83520. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_02102_00001-seg001_nrcalong/jw01076014001_02102_00001-seg001_nrcalong_rateints.fits with expected size 89280. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_02102_00001-seg001_nrcalong/jw01076014001_02102_00001-seg001_nrcalong_uncal.fits with expected size 198720. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_03103_00001-seg001_nrcalong/jw01076014001_03103_00001-seg001_nrcalong_calints.fits with expected size 44688960. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_03103_00001-seg001_nrcalong/jw01076014001_03103_00001-seg001_nrcalong_o014_crfints.fits with expected size 44688960. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_03103_00001-seg001_nrcalong/jw01076014001_03103_00001-seg001_nrcalong_rate.fits with expected size 10549440. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_03103_00001-seg001_nrcalong/jw01076014001_03103_00001-seg001_nrcalong_rateints.fits with expected size 125899200. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_03103_00001-seg001_nrcalong/jw01076014001_03103_00001-seg001_nrcalong_uncal.fits with expected size 100713600. [astroquery.query]\n", + "INFO: Found cached file ./mastDownload/JWST/jw01076014001_03103_00001-seg001_nrcalong/jw01076014001_03103_00001-seg001_nrcalong_x1dints.fits with expected size 3965760. [astroquery.query]\n", + "CPU times: user 566 ms, sys: 54.8 ms, total: 621 ms\n", + "Wall time: 5.67 s\n" + ] + } + ], + "source": [ + "%%time\n", + "# Obtain the data by its observation id (use the first file as an example)\n", + "file_list = Observations.get_product_list(obs_list[0]['obsid'])\n", + "# Download the files to a new local folder (this takes a while)\n", + "file = Observations.download_products(file_list, productType=\"SCIENCE\")" + ] + }, + { + "cell_type": "markdown", + "id": "146f155e-de38-4906-ae08-8603cd7e3527", + "metadata": {}, + "source": [ + "### Open files and load light curves" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "a1c59b7a-c180-47ec-b8e0-34ea138550d8", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "Text(0.5, 1.0, 'Example JWST Grism light curves')" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAHFCAYAAAAaD0bAAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy81sbWrAAAACXBIWXMAAA9hAAAPYQGoP6dpAABmPklEQVR4nO3dd3xT1fsH8M9N2qa7ULqhLRvKKkv2VkZFREEFcbBkuNhSEBVEBMEfWicqKlOBr6KogExtBdmyKSJTZil00L2S8/ujNjR0Jk1yb5LP+0VeNDd3PPckuffJOeeeKwkhBIiIiIgciEruAIiIiIisjQkQERERORwmQERERORwmAARERGRw2ECRERERA6HCRARERE5HCZARERE5HCYABEREZHDYQJEREREDocJENmE5cuXQ5KkMh+xsbFyh2gWtWvXxogRI8y2vhEjRsDT01P//KGHHoKXlxcKCgoM5jty5AgkSUJwcHCJdezatQuSJOHDDz/UTzt9+jSeeeYZ1K1bF66urvDz80Pr1q3x0ksvIS0tDbGxseW+X8UfFbl58yZeffVVtGzZEt7e3nBxcUGtWrUwaNAg/Pzzz9BqtZUqi6LP0KVLlyo1vyWVFsuIESNQu3Ztk9Y3Z84cSJKE27dvVzjv/PnzsWHDBpO2Q2RPnOQOgMgYy5YtQ+PGjUtMb9KkiQzR2J6ePXti06ZNOHToEDp06KCfHhsbCw8PDyQkJODvv/82KOOi5LJnz54ACpOlzp07IyIiAm+88QZq166N27dv49ixY1i7di2mTZuG1q1bY+/evQbbfvTRR1GvXj383//9X6Xj3bdvHx5++GEIIfD888+jQ4cO8PT0xOXLl/HLL79g0KBB+PzzzzF69OgK19W/f3/s3bu31CRPCV5//XVMnDjR4tuZP38+HnvsMTzyyCMW3xaRkjEBIpvSrFkztG3bVu4wbFZREhMbG1siARo4cCB+//13/P777yUSID8/PzRr1gwAEBMTA5VKhdjYWHh5eenne+yxx/DWW29BCAFJkgzWDwAajQbVqlUrMb0sqampeOSRR+Dp6Yk///yzROLy9NNP4/jx40hKSip3PdnZ2XB1dYW/vz/8/f0rtW051KtXT+4QZJGVlQV3d3e5wyAHxCYwsitr166FJEn4+OOPDabPnj0barUa27dv109788030b59e/j6+sLb2xutW7fGV199hXvvD1y7dm089NBD2LhxI1q1agU3NzdERERg48aNAAqbMyIiIuDh4YF27drh0KFDBssXNUOdOnUK999/Pzw8PODv74+XXnoJWVlZFe5TWloapk2bhjp16sDFxQU1a9bEpEmTkJmZaXT5tGzZEtWrVzdoMtTpdNi1axd69OiB7t274/fff9e/lpeXh71796JHjx765qqkpCR4e3sbNK0VV5lmrcpYunQpbt68iUWLFpVZa9OiRQt9UgfcbVratm0bRo0aBX9/f7i7uyM3N7fUZqcjR47goYceQkBAADQaDUJCQtC/f39cvXrVYH9eeuklLFu2DI0aNYKbmxvatm2Lffv2QQiBd999F3Xq1IGnpyd69eqFc+fOmbS/pTWBpaamYvTo0fD19YWnpyf69++PCxcuQJIkzJkzp8Q6bt68iSeffBI+Pj4IDAzEqFGjcOfOHYN9yczMxIoVK/RNkD169Cg3rtzcXMydOxcRERFwdXVFjRo10LNnT+zZswcAcOnSJUiShOXLl5dY9t44i5rqDh8+jMceewzVq1dHvXr1EBMTA0mSSi276OhouLi4GDTv7dixA/fffz+8vb3h7u6Ozp07Y+fOnQbL3bp1C2PHjkVoaCg0Gg38/f3RuXNn7Nixo9z9JcfBGiCyKVqttkT/FUmSoFarAQBDhw5FXFwcpk6dig4dOqBt27b47bffMG/ePLz66qvo3bu3frlLly5h3LhxCAsLA1DY3PLyyy/j2rVreOONNwy2cezYMcycOROzZs2Cj48P3nzzTQwaNAgzZ87Ezp07MX/+fEiShOjoaDz00EO4ePEi3Nzc9Mvn5+fjwQcfxLhx4zBjxgzs2bMH8+bNw7///otffvmlzP3NyspC9+7dcfXqVbz66qto0aIFTp06hTfeeAMnTpzAjh07jEo4VCoVunXrhh07dqCgoABOTk44evQoUlJS0L17d2i1WsyePVs//759+5CdnW2QZHTs2BGbNm3CU089hXHjxqFdu3YG+2ou27dvh1qtxoMPPmj0sqNGjUL//v2xatUqZGZmwtnZucQ8mZmZ6N27N+rUqYNPPvkEgYGBSEhIwO+//4709HSDeTdu3IgjR47gnXfe0b/P/fv3x/Dhw3HhwgV8/PHHuHPnDqZMmYLBgwfj6NGjVU4EdTodBgwYgEOHDmHOnDn6ZsV+/fqVuczgwYMxZMgQjB49GidOnMDMmTMBAF9//TUAYO/evejVqxd69uyJ119/HQDg7e1d5voKCgoQFRWFXbt2YdKkSejVqxcKCgqwb98+XL58GZ06dTJp3wYNGoShQ4di/PjxyMzMROfOnREdHY3ly5dj3rx5+vm0Wi1Wr16NAQMGwM/PDwCwevVqPPvssxg4cCBWrFgBZ2dnfP755+jbty+2bt2K+++/HwDwzDPP4PDhw3j77bfRsGFDpKam4vDhwxXWGJIDEUQ2YNmyZQJAqQ+1Wm0wb05OjmjVqpWoU6eOiI+PF4GBgaJ79+6ioKCgzPVrtVqRn58v5s6dK2rUqCF0Op3+tfDwcOHm5iauXr2qn3b06FEBQAQHB4vMzEz99A0bNggA4ueff9ZPGz58uAAgPvjgA4Ntvv322wKA2L17t8G2hg8frn++YMECoVKpxMGDBw2W/f777wUAsXnz5nLLbfjw4cLDw8NgWkxMjAAg9uzZI4QQYvHixSI4OFgIIUR8fLwAIE6ePCmEEOLNN98UAER8fLx++ZycHPHII48YlH+rVq3ErFmzRGJiYpmxhIeHi/79+5cbb3GNGzcWQUFBJaYXvVdFD61Wq3+t6HPy7LPPlliu6LWLFy8KIYQ4dOiQACA2bNhQbhwARFBQkMjIyNBPK3qfW7ZsafBZKSrb48ePl7vOe2MRovC9Cg8P1z/ftGmTACCWLFlisOyCBQsEADF79mz9tNmzZwsAYtGiRQbzvvDCC8LV1dUgRg8PD4PPWHlWrlwpAIilS5eWOc/FixcFALFs2bISr5UV5xtvvFFi3kGDBolatWoZvJ+bN28WAMQvv/wihBAiMzNT+Pr6igEDBhgsq9VqRWRkpGjXrp1+mqenp5g0aVKl9pMcE5vAyKasXLkSBw8eNHjs37/fYB6NRoP//e9/SEpKQuvWrSGEwJo1a/S1REV+++03PPDAA/Dx8YFarYazszPeeOMNJCUlITEx0WDeli1bombNmvrnERERAIAePXoY9F8omv7vv/+WiP2pp54yeD5s2DAAMGhyutfGjRvRrFkztGzZEgUFBfpH3759Tb76rXg/oKL/u3fvro8/ICBAH1NsbCwCAwP1+wUUlu+PP/6I+Ph4vP/++xg6dChu3bqFt99+GxEREThz5ozRMRljypQpcHZ21j8efvjhEvMMHjy4wvXUr18f1atXR3R0ND777DPEx8eXOW/Pnj3h4eGhf15UHlFRUQY1PeW9/8aKi4sDADzxxBMG05988skyl7m3LFq0aIGcnJwSn+fK+vXXX+Hq6opRo0aZtHxZSnt/Ro4ciatXrxo0US1btgxBQUGIiooCAOzZswfJyckYPny4wfdBp9OhX79+OHjwoL5puF27dvoapX379iE/P9+s+0C2jwkQ2ZSIiAi0bdvW4NGmTZsS89WvXx9du3ZFTk4OnnrqqRJ9SA4cOIA+ffoAKOxr8ueff+LgwYOYNWsWgMKOs8X5+voaPHdxcSl3ek5OjsF0Jycn1KhRw2BaUFAQAJRbJX/z5k0cP37c4ITv7OwMLy8vCCEqddnzvZo3bw4/Pz/8/vvv+v4/RQkQAHTr1g2xsbHIzc3F3r17DZq/iouIiMCkSZOwevVqXL58Ge+99x6SkpL0TStVFRYWhlu3bpXoJzV16lR98ltW36DKXOnl4+ODuLg4tGzZEq+++iqaNm2KkJAQzJ49u8TJsqrvvymSkpLg5ORUYhuBgYFlLnPvZ0yj0QAo+XmurFu3biEkJAQqlXlPFaW9P1FRUQgODsayZcsAACkpKfj555/x7LPP6n+83Lx5E0Bhh/t7vxMLFy6EEALJyckAgHXr1mH48OH48ssv0bFjR/j6+uLZZ59FQkKCWfeFbBf7AJFd+vLLL7Fp0ya0a9cOH3/8MYYMGYL27dvrX1+7di2cnZ2xceNGuLq66qdbanyUgoICJCUlGZygig7E9560ivPz84Obm5u+D0dprxtLkiR0794dW7ZswYEDB5CammqQAHXv3h1z5szB3r17kZOTU2YCdO86J0+ejLlz5+LkyZNGx1Sa3r17Y9u2bdi8eTMee+wx/fTQ0FCEhoYCuJtwlBZPZTRv3hxr166FEALHjx/H8uXLMXfuXLi5uWHGjBlV34kqqFGjBgoKCpCcnGyQBFnzBO7v74/du3dDp9OVmQQVfX9yc3MNppeX2Jf2/qjVajzzzDP48MMPkZqaim+//Ra5ubkYOXKkfp6iz/tHH31U5tWERQmin58fYmJiEBMTg8uXL+Pnn3/GjBkzkJiYiC1btpSz1+QoWANEdufEiROYMGECnn32WezatQstWrTAkCFDkJKSop9HkiQ4OTkZNItlZ2dj1apVFovrm2++MXj+7bffAkC5V+E89NBDOH/+PGrUqFGi5qtt27YmD5zXs2dPZGZm4t1330VAQIBBE1f37t2RlJSEjz76SD9vcTdu3Ch1ndevX0daWhpCQkJMiulezz33HAIDAzF9+vQyt2kukiQhMjIS77//PqpVq4bDhw9bdHuVUZSUrlu3zmD62rVrq7RejUZT6RqhqKgo5OTklHqFV5HAwEC4urri+PHjBtN/+ukno2MbOXIkcnJysGbNGixfvhwdO3Y0GJKhc+fOqFatGuLj40v9PrRt27bUpDgsLAwvvfQSevfurYj3lpSBNUBkU06ePFniKjCgcAwVf39/ZGZm4oknnkCdOnXw6aefwsXFBf/73//QunVrjBw5Ul/D079/f7z33nsYNmwYxo4di6SkJPzf//2fvsnA3FxcXLB48WJkZGTgvvvu018FFhUVhS5dupS53KRJk7B+/Xp069YNkydPRosWLaDT6XD58mVs27YNU6dONajZKk1pv7aLkpoff/zRoHYFKBxrqUaNGvjxxx9Rs2ZNNGjQwOD1sWPHIjU1FYMHD0azZs2gVqvx999/4/3334dKpUJ0dHRli6Vc1apVw4YNGzBgwABERkYaDISYlJSEP/74AwkJCSZfibRx40Z8+umneOSRR1C3bl0IIfDDDz8gNTXV4GpBufTr1w+dO3fG1KlTkZaWhjZt2mDv3r1YuXIlAJjcLNW8eXPExsbil19+QXBwMLy8vNCoUaNS533yySexbNkyjB8/HmfOnEHPnj2h0+mwf/9+REREYOjQoZAkCU8//TS+/vpr1KtXD5GRkThw4IA+wTdG48aN0bFjRyxYsABXrlzBF198YfC6p6cnPvroIwwfPhzJycl47LHHEBAQgFu3buHYsWO4desWlixZgjt37qBnz54YNmwYGjduDC8vLxw8eBBbtmzBoEGDTCo3sj9MgMimFK8OL27p0qV47rnnMH78eFy+fBkHDx7Ud1qtW7cuvvzySzz++OOIiYnRX8779ddfY+HChRgwYABq1qyJMWPGICAgoFKjChurqLltwoQJmDdvHtzc3DBmzBi8++675S7n4eGBXbt24Z133sEXX3yhv7w+LCwMDzzwQIU1QFlZWaUmdU2aNEFQUBASEhIMmr+AwoSpa9eu2LBhQ6m1Uy+//DLWrVuHpUuX4tq1a8jMzIS/vz86duyIlStXVnqgw8ro0KEDTp48iQ8++AAbNmzA4sWLkZeXB39/f7Rp0wZLly4tt1NweRo0aIBq1aph0aJFuH79OlxcXNCoUSMsX74cw4cPN9s+mEqlUuGXX37B1KlT8c477yAvLw+dO3fG6tWr0aFDB1SrVs2k9X7wwQd48cUXMXToUP0wC2V1pndycsLmzZuxYMECrFmzBjExMfDy8kJkZKTB5fiLFy8GACxatAgZGRno1asXNm7caFIN5ciRIzF27Fi4ublhyJAhJV5/+umnERYWhkWLFmHcuHFIT09HQEAAWrZsqb+NjKurK9q3b49Vq1bh0qVLyM/PR1hYGKKjozF9+nSjYyL7JAlxz6hvRGRWI0aMwPfff4+MjAyrbzsyMhIajQYHDhyw+rbJMr799ls89dRT+PPPP02u/SIi1gAR2Z3c3Fzs27cPv/76K44fP46YmBi5QyITrVmzBteuXUPz5s2hUqmwb98+vPvuu+jWrRuTH6IqYgJEZGdu3LiBXr16ISQkBK+//jpefvlluUMiE3l5eWHt2rWYN28eMjMzERwcjBEjRhiMlkxEpmETGBERETkcXgZPREREDocJEBERETkcJkBERETkcByqE7ROp8P169fh5eVV6aHyiYiISF5CCKSnp5v13nQOlQBdv35dfw8hIiIisi1XrlxBrVq1zLIuh0qAvLy8ABQWoLe3t8zREBERUWWkpaUhNDRUfx43B4dKgIqavby9vZkAERER2Rhzdl9hJ2giIiJyOEyAiIiIyOEwASIiIiKHwwSIiIiIHA4TICIiInI4TICIiIjI4TABIiIiIofDBIiIiIgcDhMgIiIicjhMgIiIiMjhMAEiIiIih8MEiIiIiBwOEyCSjVanRZ42T+4wiIjIATEBItk8vvFxdPy2I7Lys+QOhYiIHAwTIJLN2ZSzyNPl4eTtk3KHQkREDoYJEBERETkcJkBERETkcJgAERERkcNhAkREREQOhwkQERERORwmQERERORwmAARERGRw2ECRERERA6HCRARERE5HCZARERE5HCYABEREZHDYQJEREREDocJEMlOQMgdAhERORgmQERERORwFJMA/fHHHxgwYABCQkIgSRI2bNhg8PqIESMgSZLBo0OHDvIES0RERDZNMQlQZmYmIiMj8fHHH5c5T79+/XDjxg39Y/PmzVaMkIiIiOyFk9wBFImKikJUVFS582g0GgQFBVkpIiIiIrJXiqkBqozY2FgEBASgYcOGGDNmDBITE8udPzc3F2lpaQYPUh4JktwhEBGRg7GZBCgqKgrffPMNfvvtNyxevBgHDx5Er169kJubW+YyCxYsgI+Pj/4RGhpqxYiJiIhIqRTTBFaRIUOG6P9u1qwZ2rZti/DwcGzatAmDBg0qdZmZM2diypQp+udpaWlMgoiIiMh2EqB7BQcHIzw8HGfPni1zHo1GA41GY8WoiIiIyBbYTBPYvZKSknDlyhUEBwfLHQoRERHZGMXUAGVkZODcuXP65xcvXsTRo0fh6+sLX19fzJkzB4MHD0ZwcDAuXbqEV199FX5+fnj00UdljJqIiIhskWISoEOHDqFnz57650V9d4YPH44lS5bgxIkTWLlyJVJTUxEcHIyePXti3bp18PLykitkIiIislGKSYB69OgBIcq+J9TWrVutGA0RERHZM5vtA0T2gzdDJSIia2MCRERERA6HCRARERE5HCZARERE5HCYABEREZHDYQJEREREDocJEBERETkcJkBERETkcJgAERERkcNhAkSykyDJHQIRETkYJkBERETkcJgAERERkcNhAkREREQOhwkQERERORwmQCQ73g2eiIisjQkQERERORwmQERERORwmAARERGRw2ECRERERA6HCRARERE5HCZARERE5HCYABEREZHDYQJEREREDocJEBERETkcJkBERETkcJgAkewSMhOQXZAtdxhERORAmACR7F778zX0W99P7jCIiMiBMAEiRUjOSZY7BCIiciBMgIiIiMjhMAEiMoEQQu4QiIioCpgAERkp+ZtvcLZDR+ScPi13KEREZCImQERGuvnWPGjv3MH1ma/KHQoREZmICRCRqSRJ7giIiMhETICIiIjI4TABIiIiIofDBIjIVLwSjIjIZjEBIiIiIofDBIiIiIgcDhMgIiIicjhMgIiIiMjhMAEiukfKd98hZe3aimdkJ2giIpvlJHcAREqiy85GwutvAAC8+/WDulo1eQMiIiKLYA0QKVKeNg+zds/CpgubrLpdUVCg/1uXm2vVbRMRkfUwATKjLSdv4JXvjiEnXyt3KDZv/dn1+Pn8z5ixa4bcoRARkR1iE5gZjV99GADQMNALY7rVlTka25aSkyJ3CEREZMdYA2QBtzLYdGKKl3a+BMGOxUREZAVMgEgx4q7G4U7uHbnDICIiB8AEiBRFkiR5A1BwDVRBnhZarU7uMIiI7AITIKIyyZyMFVOQp8UXE+Ow+rW9codCRGQXmAARmcqKtUVJ1zIhBJCRwv5lRETmwASIFElSUO0LERHZHyZAFsArmRyEYH8cIlIYnn8qjQkQWd2FOxew9dJWucMonREHj9yz5ywYCBGRkU6uBxY3Ai7vlzsSm8CBEMnqBm4YKHcIRGQD0gu0WH8zBf39feDv4ix3OMr3/ajC/9cMAaIvyRqKLWANEFFxcl+GT0R6U89cwYx/ruKxo+flDsW26Hg7pspgAkRkAwTYrk+OZ8utwoFRz2TmyBwJ2SMmQBYg+2B+9zqwFPh+NKAtqHhemT3606NYcmyJkobgISIiO8QEyBFsngac/B6I3yB3JBW6lX0Lnx79VO4wiIjIzjEBsgDFXgafmyZ3BMqn1PeOiIjMigkQUVkU3AwnhMDZgzeRfD1T7lCIiGySYhKgP/74AwMGDEBISAgkScKGDRsMXhdCYM6cOQgJCYGbmxt69OiBU6dOyRMskcwunUjCtq9OYc1c+cf7OLb9V1w+eUzuMIiIjKKYBCgzMxORkZH4+OOPS3190aJFeO+99/Dxxx/j4MGDCAoKQu/evZGenm7lSMkqirVEpeakyhaGYtzTMpf4rzKaM6/9HY8dX36C796aJXcoJtHlKuvCAK1Wi3/++Qc5Oda96kmny8X16/9DTs51q26XLIRN+ZWimIEQo6KiEBUVVeprQgjExMRg1qxZGDRoEABgxYoVCAwMxLfffotx48ZZM1SyguKXfSfnJKOaazX5gqEypd26KXcIJkvb8S/SdlyG77DGcG/hL3c4AAprwuPi4hASEoKxY8dabbsXL36MS/9+CrXaEz26szaPHINiaoDKc/HiRSQkJKBPnz76aRqNBt27d8eePXvKXC43NxdpaWkGDyIiAEjbcRkAkPLDWZkjuevYscLk4/p169bEJCXvAgBotRlW3W5FlDaiCNkXm0iAEhISAACBgYEG0wMDA/WvlWbBggXw8fHRP0JDQy0aJ9kBVh0TETkEm0iAitw7wKAQotxBB2fOnIk7d+7oH1euXLF0iP/FZZXNmIA/p4yhy+AVVkRyKn4sPZqWJV8gZJdsIgEKCgoCgBK1PYmJiSVqhYrTaDTw9vY2eDg2xWZmylEsob4VEyNfHBXhW0kOpt9f/8gdAtkZm0iA6tSpg6CgIGzfvl0/LS8vD3FxcejUqZOMkZWO7db2If/qVblDUDQhBLLYr86sFDuIKpEdUsxVYBkZGTh37pz++cWLF3H06FH4+voiLCwMkyZNwvz589GgQQM0aNAA8+fPh7u7O4YNGyZj1GQpkgzNdXd+/hmZe/ZafbuVocTz4m/LPsfRrRvlDoOIyCSKSYAOHTqEnj176p9PmTIFADB8+HAsX74c06dPR3Z2Nl544QWkpKSgffv22LZtG7y8vOQKmazEWndCvz492irbsRdMfsjaKur3SUUU+ItJgRSTAPXo0aPc6l9JkjBnzhzMmTPHekGRfJRwjOOBtkxCp5M7BHJATH7InGyiDxCRLITA7S+WIj02Vu5IFCcl4YbcIZADYL5jKhZcZSimBsieKLG/RiF+KYyREx+PnPh4AEDE36dljoaIiMyJNUBElaDL5JhAtkTka3F72Ulk7FPWva1ycnLw/fff48yZM3KHYkCOiw6q4uCdTPx6K1XuMMjGMQEiRVLaAfnKuPHyBqDcakVFytiXgJwzKUjdcF7uUAz88ccfOHnyJNasWYMjR47IHY6etS40MJcBh89i5MlLuJSdK3coZMOYADkU2zrIKUnWoUNyh2BA6OR9L5XeN0MYc5d3KxZlenq6/u+ffvoJ165ds97GbVBFef/VnDzrBEJ2iQkQUSVdHPwY8m/I0/n33vPApRO3ZYmDzCs1NbXM15KTk60XCJEDYgJEiqeU0XFzTp3CzfkL5Nl2Rr7B89wsI2o4yCb99NNPcodAZUjPycePR64iPSe/4plloYxjptIxASJF+uToJ3KHUCpdljw3ZDwRe/e2HLlZCjjoKr0NzEYVT/azs7NljITKM2ntUUxedwyT1h6VOxSqAiZAFqCocwMHrLMLeTla/d8n4thvpCLFKw0z9l6HKOD3wBYp6lhazM6/Ew3+J9vEBMhMdMU6pSqkxQb4dy8wt3qxCQo9mpBR9v90ARkpvPqlslJ/Oo/0XUwaicgQEyAz+eW4ssYbAQAs6yd3BGQu5WTVKQmZ0GmtW8OhtGEKStAallfepTsyBXJXQUEBTpw4YTDtypUrZc6faaWxp4QQyMhQ1rhElaWYH5tkk5gAmcmVZHn6hpB1FY0MrSTfztmPJS/GIiMlR+5QFCM9tuzEAgByL6RaJ5Bi1q1bV2La+fNlj1NkrQToypWvIQQvJzdVTr4WU9Ydxa8neHsYW8MEiMgI2pQUuUMo06HNl+QOwSYIrcCtL05UPKOZnT17tsQ0rVaLffv24datW1aPp8ilfz+XbdtVdUUB4wB9tfsifjhyDc9/c1juUMhITIAcCuuLiSDzIJLFJScnY8uWLfjkE2Ve9Si33Areqylnyq/ps4Zb6eyPZ6uYADkqLceRoSpQ6uU5ZHapd/7CnbRjcodBxmDnqEphAmQFeUq7BPenl4CF4UD6TbkjqRRbu08Rkb0oKEjHX389gUOHBkGnk7+5SYmUMlArGY8JkJkU/w4U/zpM//4YGr72q0I6Sf/3q/3IKiAvAzj4pbzhEJGi5effvXpOp1PAAJwKtGLvv3KHQCZiAmRh/ztUOILvV7svyhwJkfmwBcwyWJtQvtMZ2UjJZ/M9mYeT3AEQVUTxY86Q4gktEwt70PPgGah5OKgYf6FUCmuALKC0j97yPZdw8bZ1xvUoU0Yi8L9n5Y2BTFKZigGe4suWey4V+QlF3z+WlC2zlVz2VnouCqw8QCkZhwmQFT234qC8AcTOB+J5h2lbdOtyutwh3KNkmv/dvNeQl6PcG3im7bwsdwjkIE5eu4P73t6BoV/skzsUKgcTIDM5eb3iofbP38pkGz/ZrcsnjuLw5p/lDsNoIldb8UwOq3iiy2NXZa07WDg+0aF/lTtwKjEBMot8rQ5bT1XukvLflHT34PifgFT+Kibz+XPdKuTn8pYcRLLiD+1KYQJkBrn3jPNT3kdvx2kFJUC3zwAxzRX/ZeE4QLZl3w/rkHjpAs7/tV/uUMiG/XhTmbUnihvXjUzGBIjYL4iMJpVzlcmtSxewKnoCNix6C7f+5fAPQOE9vypy6NAhpKcrra+XfD74V5kDtf7xj3z3bSvXNd6LzFhMgCwot8BG+hYkKu8O53RXwoWK+5cpVcqNa3KHoAjnzp2rcJ6NGzdaIRKqKp1Sa8xvnpQ7ApvDBMgMyurYnJNXWlWpQr88pFjrF/1VuRkV8tGSOwxtpvJGLNbp7KHZRO53VhnKq/0k28IEyILYd8U8dMIeTh52RsHngBvzyr70WAih9C5vVlfeCZ0ne7JnTIBI8T4+8rHcIZARrp6WuSq+jAQn+1QSrs3cjbQd1r93k5KHv1BybEpUXkp4LjEDn8aeQ1aejXR/cHC8FYYFlXZc4bHGeLFXY+UOwSbE776O9KRs9H8xEmon+X7bFOTm3n2ipBoEXeGXL+MP9kuyBX9nKnM4hfI+0hPXHsGp62nWC6Y4nlyMxhogC/hq90Xk5PMXAFnfldMpOPeXNYZaUFBiQw7tnQs3EJesjKvnZEt+rh8BfplQbAKTocpgAmQhy/dcsp2PoE4L5Ml8nzIym8vxSXKHQPdgXxrLifn3JoYcO2+17Snyrfx+lNwR2CQmQBaSmpVf6oBZ1+8osFp31/8B80OAzNtyR0Jm8M9+5YyfIrGmCIB99LPJyDgjdwiKoMjPtK5A7ghsEhMgCzlxLRUdFuwsMV2xg2gBwNltckdANsKYGo2r8Sex57tvoKvEYID2Ki8vT+4Qqiwz86zcIRCZFTtBm8GO0yV/cf95js0QRACw7s0ZAACPatUR2ftBmaOhe7F5zkhKLC7br2CUBWuAzGDyumNyh0CkeCk3rssdAlGVbTx2Q+4QyEyYABGRxfz83ny5Q6AK3Ns/KSc3QaZIbMP6w1flDqEU91YBKbGaSnmYABGR8Xh8NYotdYLes6e73CEQWQUTICIi0hOi+BVFdzNd3tpHwUok2HyvKoMJEBHZDcEBSImokpgA0V0Zyhk/hpRNkWOhAMi/mSV3CDYnP98+rlg9l6XAMdashjU+pmACRHel8eoGIkeSnX1Z7hDM5kq27Y+1RNbFBEgG8zbGyx1CGfgrgiyJny+lSUk5UP4MHCPINtzbByg/C/ioDXByvTzx2AgmQDL4cvdFuUMgqhqeGO2aVpuN7Owr0GnlaVb6R6F3grcpSed4j7AKcCRoIiILs6XL4AEgNq5ZyYlW3If4jGyjl7GtEjY3x957U7EGiIrhr3oix8ITpzGy83iVoT1hAkR3sVmDLMoKny9+hu2CUtOybfEKHSU7N13uCGwSEyAiMppib6BpY01NRGaRlyF3BDaJCRAVo9CTGhE5lOPpHM+JLI8JEN11bI3cERCRVRlTY2a92rWbeQUVz3QP1v2RsUxKgK5cuYKrV+/eEffAgQOYNGkSvvjiC7MFRjLISZU7ArJrPEXZMt4LTMFNv+U5twPYHcPm4VKYlAANGzYMv//+OwAgISEBvXv3xoEDB/Dqq69i7ty5Zg2QiIgci1KHDbDB9AdYPRjYMRs4v1PuSBTHpATo5MmTaNeuHQDgf//7H5o1a4Y9e/bg22+/xfLly80ZHxFR5dniL3SyGTb98bpzteJ5HIxJCVB+fj40Gg0AYMeOHXj44YcBAI0bN8aNG7yfFJG9s8mmABkptUbDnpq1rLEnSr0JMJnGpASoadOm+Oyzz7Br1y5s374d/fr1AwBcv34dNWrUMGuARERESvDPTY63Y09MSoAWLlyIzz//HD169MCTTz6JyMhIAMDPP/+sbxojIrI6hda02AfrlW2eQt/HD3aelTsEMiOT7gXWo0cP3L59G2lpaahevbp++tixY+Hh4WG24IhIodgEZheEqPytHbKyrHcT52ytzmrbIsdlUg1Qr169kJ6ebpD8AICvry+GDBlilsCIiMiysrMvV3req1dXWjASQ6bU/yi1nxUpl0kJUGxsLPLy8kpMz8nJwa5du6ocFBEROS7WAJE1GNUEdvz4cf3f8fHxSEi4e2M4rVaLLVu2oGbNmuaLjojIDrB2wjj77mTKHYL94WewBKMSoJYtW0KSJEiShF69epV43c3NDR999JHZgituzpw5ePPNNw2mBQYGGiRhRERkDJ4UHc7+LwB3X6D5Y3JHIjujEqCLFy9CCIG6deviwIED8Pf317/m4uKCgIAAqNVqswdZpGnTptixY4f+uSW3RURk9+yoVsB+9sSCki8Av75S+DcTIOMSoPDwcACATidP+6yTkxOCgoJk2TYRVY0dnWvtB6/mszidTkClUkA5b5wEPMfbYRRn0mXwAPDPP/8gNjYWiYmJJRKiN954o8qBlebs2bMICQmBRqNB+/btMX/+fNStW7fM+XNzc5Gbm6t/npaWZpG4iByNKSNBH978E3oOH2OBaIiU61JSJur6e8odRiFdgdwRKIpJCdDSpUvx/PPPw8/PD0FBQQYHQ0mSLJIAtW/fHitXrkTDhg1x8+ZNzJs3D506dcKpU6fKHH16wYIFJfoNEZEdYy0TKYyOn0nFMikBmjdvHt5++21ER0ebO54yRUVF6f9u3rw5OnbsiHr16mHFihWYMmVKqcvMnDnT4LW0tDSEhoZaPFYiIpvAdkmLU9QVgMVjuXMV8KklXywKYFIClJKSgscff9zcsRjFw8MDzZs3x9mzZQ9NrtFo9DdtJSIisjZF1QDt/fju3xmJDp8AmTQQ4uOPP45t27aZOxaj5Obm4vTp0wgODpY1DiIiR5GSckDuEGyOTkk1QH9vlDsCRTGpBqh+/fp4/fXXsW/fPjRv3hzOzs4Gr0+YMMEswRU3bdo0DBgwAGFhYUhMTMS8efOQlpaG4cOHm31bRESOICfnulHzHz7yJHr2OAOVyuTrZxyOohIgMmDSp/iLL76Ap6cn4uLiEBcXZ/CaJEkWSYCuXr2KJ598Erdv34a/vz86dOiAffv26S/NJyJSKkX1AykmJXW/0csIkYcqXEBsMSn5lb+xqzUp9K0nmPgpvnjRencFLrJ27Vqrb5OIyBxMGTbAGvLzk+UOwWw+uZyIIcG+codRgnJrgJQal/WY1AeIiIgqT6kJkLt7HblDMJs8ocwbqCqqEzQZMKkGaNSoUeW+/vXXX5sUDBFRVSi1qUmp7Km8lLoryq0BIpMvgy8uPz8fJ0+eRGpqaqk3SSUiIvtgT0mTNbC8lMukBOjHH38sMU2n0+GFF14o99YURGQcHjyJKkep35TUrHy5Q6AymK0PkEqlwuTJk/H++++ba5VERHZBuYms8XEptT+TUks4PYf331Iqs3aCPn/+PAoK+GYTkTyUenImyxOKTYEUisVlWhPYvffeEkLgxo0b2LRpEwcmJCLZKLWm5caNG3KHUAbjy0upZaxUTMyUy6QE6MiRIwbPVSoV/P39sXjx4gqvECMicjSKTYDs6NzMvIyMZVIC9Pvvv5s7DiIiu+Xl5SV3CGaj02UD8JA7DKIqq1IfoFu3bmH37t34888/cevWLXPFRERkV1xdXeUOoXQmdJm6eu1b88dhBkqtAGLNlHKZlABlZmZi1KhRCA4ORrdu3dC1a1eEhIRg9OjRyMrKMneMRA6LB0+yKBM+X5kZ/5g/DpIBDy4mJUBTpkxBXFwcfvnlF6SmpiI1NRU//fQT4uLiMHXqVHPHSOS4mAHZBaV2HLanDrpK3ROFvvUEE/sArV+/Ht9//z169Oihn/bggw/Czc0NTzzxBJYsWWKu+IgcGi/rJkuSTGgDU2rSpNREQ6FhEUysAcrKykJgYGCJ6QEBAWwCIzIjScUEiCyIHy9yYCYlQB07dsTs2bORk5Ojn5adnY0333wTHTt2NFtwRESkLKbUGlmDcmumlBkXmdgEFhMTg6ioKNSqVQuRkZGQJAlHjx6FRqPBtm3bzB0jEVHl8FxjHJYXOTCTEqDmzZvj7NmzWL16Nf7++28IITB06FA89dRTcHNzM3eMREQ2zZ5qARRb0yJ3ALbGjj6TpjIpAVqwYAECAwMxZswYg+lff/01bt26hejoaLMER0RElqPUZMaesISVy6Q+QJ9//jkaN25cYnrTpk3x2WefVTkoIiKyBvs5PSt2TxQbGJmUACUkJCA4OLjEdH9/f+Xe84aISCYRERFyh2BGPKOTfTApAQoNDcWff/5ZYvqff/6JkJCQKgdFRGRPqlWrJncIZpObo8wfuezSQsYyqQ/Qc889h0mTJiE/Px+9evUCAOzcuRPTp0/nSNBERHYsvyBN7hBsinL7WSk1LusxKQGaPn06kpOT8cILLyAvLw9A4c3+oqOjMXPmTLMGSETKo9Nq5Q6BzMKUk6AyT5zKjIo1U0pmUgIkSRIWLlyI119/HadPn4abmxsaNGgAjUZj7viISIHysjniuzEUexm8SXEpdF+IjGRSAlTE09MT9913n7liISIb4erpJXcIJBeF5j9KbWpSZlQEmNgJmogcm5uXt9whkDmYcrNdZd4Jg8hoTICIiByVKU1gVqjS8HVWG71Mcr4y+6UptfWTmAARkT1R6NlGMqWmRaGs0dSk1Buu2hWFflesiQkQEZGDMi2Z4YmT7AMTIHtWvbbcERBZlx3VtFiDEKY0GzEBMoZSO2cTEyAisies1jeKSqrShcBUCfxIKhcTIHvGbx4RlUNSucgdApFsmAARETksU64Cs/wPqyg/H4tvw1qU+zNUuZFZCxMgIrIfPKZbnIDO4tto5e1u8W0QMQEiIiIih8MEiIjsh0IvAlPsvcCIHBgTICIih8XEzOKY/CoWEyAiIqo8ntCNolNqcfF9ZAJk3/gBJwfDj7xxTLoVmOULmW+jNbCUmQARETkopY5SrGPtBFkBEyAishuSEw9plue4NUA1q7nJHQKZEY8WRGQ3JFe13CHYFNPuuu64CZDG2fhTJq8AVC4mQERkkqD6DeUOgarIlCYwa5zQlZoyKHSUBTIREyAiIodlSqph+fTESaGZhiQZH5hSkzleBcYEiIhMZMrJwFGxGcQ4zvxskRUwAZKJzhqDQ/CYSxZ038OD5Q6BqkyZNUD2hLmvcjEBkkmBYkfHIqqc+vd1lDuEkvi1Mo5Cz84N3F3lDsFslFnCBDABko1Sx98gqixFNoEpMCQyXhsfD7lDKJV9fbx4DmICRER2Q5FJGQA/Pz+5QyiVUq8CA4D7vJWXBNX1Nz4m9v9SLiZAMuF3guxB2wGD5A7BJri6umLatGl49NFH5Q7FDBz34BXobT9NczwJMQEioiro/vQouUOwGZ6ennB1VdoJlJ2gjWFKzsA8Q7mYAMnEKl8Kd18rbISIHElBQZpVtmMv/STtZT/sERMgmVjlSzH4K8tvg4jIApg2kKUxAZKJVWqA/OpbYSNEysGTpv2wl/eSTWDKxQRIJvxOEJHcvDybyB1CmZSYOJh01ZwF4jCLgly5I5AdEyCZ8NJIIpJb8+afyB1CmezlCKnYQ/3R1XJHIDsmQGYwtbfxd8VW6neCiByHSuUidwhleiSgmtwhmIViO0HrCuSOQHZMgMygWS0fo5dR7K8CIiIFiPI3/rhqabwM3r4wAZILvxRERJXWqZqn3CGYdNjefOKG2eMwD2WOmm5NTIDMoEGA8V9MxVaLEpHFeHrKfxK3Vata1EF7hd4jrDynrltn3CQyns0lQJ9++inq1KkDV1dXtGnTBrt27ZI7JNSq7m70MqwWtQ2uTZR7lQzZnpo1a8odgs0o7RAp92HTro7b2jy5I5CdTSVA69atw6RJkzBr1iwcOXIEXbt2RVRUFC5fvix3aPbDCt/wz3t/bvFtmIt3/wflDoGMYVdnKMfGt9LC/tkidwSys6kE6L333sPo0aPx3HPPISIiAjExMQgNDcWSJUvkDs1ox66mlpj27ta/Mfunkyavc+Px61i45W9k52lNDyynZFzYOAXY8abp64x7F9jwov6I1imkE74b8J1Rq9Dq7u6TEAL5unxM/2M6vjph+mjXXxz/Al+e+LLceZzDwsp9vWg4A6HTQeTlQZuWhivjxiNl7VqTYsq/eROnG0fg4pAhZc7zwIgI+IaU3xRQFJdOJ6DN1yHzTi5+/uAIzuxPMCkubb4Omz49jr0/njNp+bJkpqbghwWzcf6v/SYtX5CXh/0bvsPJ37ebNa6C5BzcXnYSOedTARSWp9DqKr28Lk+LO9suIfvk7RKvDRw40OS4bt68iVWrVuHatWsmLZ+amorVq1fj1KlTJsdQmjt3juDI0RHIzDTt85GReRYnT01GYuLdk3KpNUCVSIqEECjQCeRodfj1ViqeOHoOV3JMq+04kpaFF+L/xV93Mk1avixbTyXguRUHkZxZGJdOJ1BgxOfrXGI6Fm87gwu3MswaF45+C6wZBuSZuL8XdwE/vQgkXzRvXBYiCRsZkCYvLw/u7u747rvvDO6oPHHiRBw9ehRxcXEllsnNzUVu7t3BntLS0hAaGoo7d+7A29vbrPHVnrHJ6GXcnNXQ6gTySvngO6kkSBIg4b//JUAlSZBQ+D8K/0GlkqDVCaTn3L2k8a/XHkANT03hkzkmXEnh5AYU5KD0Q1AZHeekMqaLYvs2O1U/X3xSPIZsLPskXxoXlQu0QgutMEzwVJLqv3KSIEEq8VySJKigAqTCee/k3jFY/sTwEwCA040jDKa7tW2D8FWr8HdEBc1gklTqkVlydjZq/wBA5Ofr/474+zQA4JPxvwEAguv7YOCkVlA7qZBw4Q7WL/qrgrhQ6luoUhnf+VGnK1yRs6saY2O6G7y2eMhDRq3Lxc0NednZJaZLKlXhe/bfB14qfMMK30OVBPz3v4TC13My0vXLBtdvhGFvL0b+rSzcXFxBuZRB5eEMCAGhExA5xT5jKgnQ3S1IyUVVolwL3/7/JgoA2sK/VZ7OCHmtQ6nbS0tLw3vvvVdhXC4uLsjLK3kC15dVsb/Le55drMyDgoIwfvx4AMDO3+pVGEPJbbtApXKBVlvaCVi653/o47h3OiBBiLuf+eCgQWjS5F0AwPmsHHTe/7f+tfNdm2PIsfM4lJZldLwA4FLWMaoMecW+05PDAxFdNxgzfziBNQeMb3Hw83SBEECBTkCnE0jPvXu8Lv7xclGr4KyWIFD4mRIQ//3/HwGD88XixyMxuE0t047zQOGxXu0C3HNMBABI6sJjm6QC/vvO3f1bBf2Jqeh5dvLdZZ9YCTQxPdEvTVpaGnx8fMx6/nYyy1qs4Pbt29BqtQgMDDSYHhgYiISE0n/VLliwAG++WYWaCwvLzi+7pqZAf8A1Pj/VOKtNjKho4yVPTneVEY+RebQpeXeervRfcbqiJMuMqXyNsWMRMGVy5WYuY1+KJzNV4e7tgqy0PDTuEAy1kxGVtmWUh05nekHl51ShdvE/pSU/wH+1aCauMz3plukB/UeXWcb7dU95ibzK/1JXuZV9iK3sQby05Af4r1bKxN+v1atXN2m5u9vOg7bMPiQlj12VDdPHp02JtRRX+ZIvKa8Kv/Uberj+95dp67idUXYNVPGPV55WB2Mq8D1dq3gKL8gu+3gvtKYfU71rmRySNdlMAlREuieLF0KUmFZk5syZmDJliv55UQ2QEux/9X7kFeiQmVeAAq1AdQ8XuDqpoP3v2yAA6ERh9l/0f9EvAp0o3G/df786M3O1GPjJnwAAT40TPDVVeFsnHi/8Pze9MLt39QGcXAGdtpxannK+JYv/GyQyuKXB8jph3KHs10G/wknlBJWkglpSQy0VJnkFouC/stD992tJ/FdGpT+HADLzMzFs8zAAwIN1yujjU8lfi/V/2wk4OUFydoakVkNSqyF0OugyTatCPtejJwDA74Xn9dMefKEFrv6djEYdgyq9nmfnd4JKLUGtVkFSS1CpJOh0wuQEZsXMws/XAyOr1il8zMdfQ6stQE5GOlw9POGscYVKrdafzIXQAQIG/wv9c3H3pC8E8rKz8e1rUwEAA6a8WqW4/EY2hdpHU/hzXAIklQRRoIPkrIaklgC1BJGvMzxbFf+MFKv00KXnI/GTowCAgOcjTY7pvvvuQ5cuXaDVapGdnQ1XV1e4uLhAkiSD5EcUL5dynp85cwY7d+4EAAwpp4m1Mjp1/ANC5CG/IA3OTj5Qqz0My8PgmFDy7+JXwV6//j9cvPgBAKBmzaGlrwKo8Krt2m4u2NCqAdQS4CRJcJIkqCUJmVot8k1I+t+5eAP/S0gBADwaWLWEcfOErlCrpLsPSUKeVgs3Fyc4/zctO18LrU4U1nDibnEWVrT8Nw3AhDVHcOjfwrj6Nq38MaFUk04ABXmFNUDufoCz+3+19kUnnKK/dfc8/+9R9JquAFjS6e56a7UpfXsKYzMJkJ+fH9RqdYnansTExBK1QkU0Gg00Go01wjNaoLdrxTMZYVy3uvj8jwuYPaCKVy1VDzdPQEWaPALEbwC6TjWYrDPyt1wtL/P+omgT2AZ/3fwLQxpV7UTgHBJS6nS1l5dJ63Nt2hQ5p07B+6EB+mmBtb0RWNu4Kl8v39I/X5pyaiTK4+bljOz0fNRqVLUTgbd/QJWWL0uNWuX31aqIayNfM0UCCM+7oytLVfiF7u3tDR8f8w0G6OXlhZ07d5rlmOjmZr6r2YKDHsXFix9A42J4HDe2D1CIxgVBmpLNzu5q07q6jq3lj/8lpKCll/FX+RbXvKYPmoSYr8vFxAca4JmvDqBnI/+qrah6HaBa1b43BgYtBX4YAzR/wnzrtDCbSYBcXFzQpk0bbN++3aAP0Pbt26vUqdBezIhqjNFd6yDAy7yJVZU99jWQ9laJL5rcXc+W9lmKpOwkBHlU/AvKpW5d5F24UGJ6rc/M3/m+9to10KalwalGjXLnqx5UtYOysZ55uxPyc7Rw91bWrRNe/GotdNoCaNytWx7lkZxUCHmjQ2H/JRP6W1mKm5sbpk+fDmcT+qZZkptbKLp03gsnJ8NkT1NK2ZV31DD32GrNvNxxrFNT+DrfPU2actgysutRhbo28Mfemb2qfqxvO8o8ARVp8QQQ2g7wMWNSZWE2dRXYlClT8OWXX+Lrr7/G6dOnMXnyZFy+fFnfmc+RSZKkvOQHAFTqUn9lGNsEZm7OKudKJT8AEPpp6TeMdG3c2JwhASjsOF1R8gMAGveyT2JG9ROqJGcXdZnJz4tfrUXT7g9gzCfLzL7dirh6esLdp5r+uTZVGXe4Vrk7Q1XF/hmW+JHg7u5e5QSoXt2pFc9kJI0mAGq1Yc1UuJsGPaob1qSWl+SoLDCycaDGGc7FEjGlXDIU7OMGdVWTa41ptdTlql4bUNlOWmEzNUBAYbt1UlIS5s6dixs3bqBZs2bYvHkzwsPN3GxDFmdLI2G71K4tdwhGMaavkDm4enqi3wuTrLrNsqg8lFW7YY/UauvVti1sVAvt953WPy/vqDEutIpNQhbSNER59zQDYP6qKRtkO6naf1544QVcunQJubm5+Ouvv9CtWze5Q1I212pyR1Cqxr7mrz0hcqrhCueavN2EJSUlyz/6fmmqOVXx6tdKMOWH28ORpfcTJPnZXAJERnr5sNwRlMrD+e5AfhG+EeXMqWBKqQ+/hyP/rlNpnODWtOImRDJdXl6SfBtX5leuXE5qpX4jlRqX9TABsnceyj8ZOKmU0RLr0alTma8FvjoTAOA7YoSVoiFSJj+/XrJtu7z8p6W3cjrCF2eV30n9Fhq/DJvAmABZwjuDmssdAplA7Vv2pdC+zz6LxieOw71DeytGZBob/JFMNkQlyfeDpUv1sps3XWyo863ZeZU+FEz5mAA58CfGvBYOLkx6Fj8eiaHtbOcyQKq8Ere2UGgTGJElyXkBw+NB5huvyRR29ZWXePpXRtuDHRhyXxgGtqwJ16rehoLIhlULDEbqzRulvubkbKUxhPjD1rJkzAJs8a2Ve8yzMjEBYg2QOTH5sXFsE6+yrsOGl/lataBgK0Zi+5R64rSlISwchimfFSdlDWoqB9YAEZHZaDzK7qOh1BM6Gct672MtVxc0dHeFq0qCu0ole+oVXM3N6GXkjpnKxhogolL4Pv1UqdOlMm/4aH2PTG4FtbOyvsKhTXkBgDGcnGzwN6gVR3FXSxJ+u68Rfm3bsMybXlvT+O515Q7BjOQvT7kp6+hJDkmJ4wA5+ZcxqqwCDsJFajaqjh7DGskdhgGVSgnNwMp5jyrSvLntJYwhIVW7gbCxnP67ezogf22Ku4sTBnBgQ7vBBIhks/7h9RjbYiymtJ0idyiVp6AEiGxfr15lj6njX1YSLjNXV8fuy8WmXPthg/WvZC8aVm+IhtUbyh0GOTiVp3z3D/Py8kKjRo1w5syZEq9FRCivZpSMx3xJuVgDRERmFVSvgdwhGEXtqcyrYZTQ50VpbLH2xSpXzXk5dq2cqZgAWYiXRkGVa0XDpA/+St44yCrkPm0G1pO5Vk/uAiCytvCOQO+35I7C5jABshQlHYQ7jAdevQE0f0zuSOyAkt7Y0sn+G9kGf6UT2bzOE4ybnzWMTIAsRXEfLRdl3ihQSao/8TgAwK1Nm0rNz/N8WVgwtiQocCC8PJvKHUal8JNF5qSgdhoiebnfdx/qx/4OJz+/smfir6YK2WI/DTnJXV5Nm74HIQR++72+rHHYLX4dFIsJEFExzkFB5c/Ak3vFWEQ2hx2uLUexX4ewTnJHIDs2gRHZMB9/44fmtzTeK8o4xZOPqKgoAEDfvn2tHkf7dputvk1j8ZNlRl6BckcgO9YAEdmw4PrV0POZxvh91d9yh6IndDZ2mlJQ5Uf79u3RokULuLlZP7H19FTWqOJKZWOfbioHa4CIjKHApoImnZU2NL+8pwhbb86RI/khckRMgCzE1g/CRKaSu1Ove6sAWbdPlmOLtS/sNqhcTIAsJKQaf8UZY0W/FXKHYAIe2Uolc7GovV2g9lbm6M5EpBxMgCxkzoAmcodgU0I8ldaMUwZW7FWC/Imh5KKEu9KTuSmtNmV3dE+5Q6AqYAJkIX5eGrlDsCm+rr5yh2A8pR2NFYNZIjkGVSW6OshyVWRQc+tv0wYxAbIQngKM46J2wa4hu7Cs7zK5QykX+3ZVrNMTw+QOgUgxZPmdFMAWiMpgAkSKUc21GgI9So5NMbfTXBmiIVN5+7ETMlmGtxObNsl8mABZiJOKRWuKUK9QPB/5vP75I/UfwaMNHpUxonuwBohINnXdNXi9XgjebxwqdyjK5ltX7ghsAs/SFhLqy6vATPVCyxfkDoHIITVsOFvuECr0YlgAngyuIXcYyvTKBWDyKaDVM3JHYhOYAFkI+4o4AHaCJjPw8fGROwS9miFD5Q5B+Yp97RV3mPeoAfjUAnxqAg/FyB2N4jEBIkWTe1C9EhR3xKMqk/ktdXd3lzcAsihfD5nGpNJ4ybNdG8IEiMgoTIBMpbhklsgKmtVUTg0fGWICRERERA6HCRApmiyDiJFFOGtc5Q6hdDJ/xHjzUyJ5OMkdAClIr9fljoDsGAdILF2bNm1w5coV1K9fX+5QiBwKEyC6q9s0uSMgO1anZVu5Q1AkJycnPPbYY3KHAQCQJGe5Q6i0SC83HEvPljsMsmFsAnM0NdvIHYFtK34VmEI79foGe8gdQqkUOzSEQsOSgyRJiGj8jtxh2KS6fsr83lHZmAA5mlHbgGln5Y7CdtnAybJZ95pyhyA7z84hcodgs4KDH0PTJu+X+ppGE2zlaMo2sqaf3CEYUGp+T2VjAuRo1E6Ap+3cq4mXThtPrebX2qNdEGo8yxtCmkKSJAQFPVzqa2Fho60cTdkaeSir87hiazipTDxSEpHdkdQquDXh7RLI/NQqG0l0mJBViAkQkREkZ9vpJEp3uUX66//26nXPjTR5oiAjTO/XCLWqu+G1/hEG01lbbXt4FRgpmtLGAXJt2lTuEMgEah+XYn9rZIyEzEUjU01Mreru2B3dCwCQcCdHlhgqhQlZhVgDZCVNQ7wxvGO43GFQFbGuwPa5NvaFZzd2FK9I+3ab5Q6hXBEeCh1YU4lG/ip3BIrEBMhKNk3oik71lXXVgi2o41NH7hAMqdX6P1Xe3jIGQqaSAFR7sK7cYSiep2cjuUMolyRJ+CQiTO4wbEN4J7kjUCQ2gZEirYpahbircRjedLjcoRiQ1GqELV8GXU4OnKpXlzscIoempCuvvN3YP9DWsAaIFKllQEtMbD0RGrXy+mt4dOgArx495A5D0XqPeQkqNX9fkeN4pW8jeLioMSOqsdyhlM6NP9juxQTIitgnjRxFiwf6YcLK7+QOo5i7NQX8GtoPJV151aymD47N7oPx3evJHUoZlFNbphRMgBxV8yeA4JZyR0F2TKVSVzyTtfDYT1bg9N8gpNH9FFALpKDmQaViHbWjGry08P85PoX/q13KnpeICIDETLJSlDkYu3Jqy5RCkW8TEZFZ8dhPjk5BzYVKwQTIqhT8AeSXg8yNVfBECsJj/L2YABHZgcA6HJOIrIFJbXHumrv93FyU2e5F5WAfICI70KBtIG5eTJM7DNvEc7rN0qjkTTq8XZ3xybDWUKsAV2cFdfovDSuASmACRERElSMpq5ajt583Wnu7o623h2wx9G8RLNu2jdJhPBC3UO4oFIUJkBWxmw2RAvB7aLLgoEFyh2BAo1Jhc5uGcodRgiKvlusezQToHspK58n6GvYr/P++0fLGQWRJCjwf2QJPD8P7gTk5yVfTQlWkpHG5FII1QI7usWXA5b1A7a5yR0J2rNfIcXKHQCZo3fob/LGrrdxhkEmY9VeENUCOzsUdqH8/4MSBEMlyQpu2kDsEMoGzc3U0b/YpAMDNrba8wRCZGWuALCjQW4Obabn650E+rjJGQ0RkPH//Priv7Qa4u9eROxQyCju7VYQJkAWtHt0e8zefxuTehZ30WoVVx5sPN0V4DXeZIyMiqhxJkuDt3VzuMIjMjgmQBTUI9MKyke0Mpg3vVFueYMiuKemu2EREtsBm+gDVrl0bkiQZPGbMmCF3WERERAbcXHjFlS2wqRqguXPnYsyYMfrnnp6eMkZDRERU0mNtamHrqQR0a+AvdyhUDptKgLy8vBAUFCR3GERk09hcSJbl6qzGqtHt5Q6jbPV7yx2BIthMExgALFy4EDVq1EDLli3x9ttvIy8vr9z5c3NzkZaWZvAgIsfDEVGIitGw9QSwoRqgiRMnonXr1qhevToOHDiAmTNn4uLFi/jyyy/LXGbBggV48803rRglERGREjDtr4isNUBz5swp0bH53sehQ4cAAJMnT0b37t3RokULPPfcc/jss8/w1VdfISkpqcz1z5w5E3fu3NE/rly5Yq1dIyIiIgWTtQbopZdewtChQ8udp3bt2qVO79ChAwDg3LlzqFGjRqnzaDQaaDSaKsVIZAt4FTwRVR5rhwCZEyA/Pz/4+fmZtOyRI0cAAMHBweYMiYgcjcSTAZEjsok+QHv37sW+ffvQs2dP+Pj44ODBg5g8eTIefvhhhIWFyR0eEdkwpj9EjskmEiCNRoN169bhzTffRG5uLsLDwzFmzBhMnz5d7tCIqBLcfarJun1N/WpIj7sqawxEpCw2kQC1bt0a+/btkzsMIuVSYB8gSZLw5FvvIj8nF+7ePrLEEPxae2hTcuEU4CbL9olIuWwiASIi2xTSMELW7as9XaD2dIHQCag8nQGdgNrDxWAeBeaORGQFTICIyO5JKglBU9pACEBytqnxX4lMo9PKHYHiMQEisgOC9RgVUrk7yx0CkfW4essdgeLxpxAREZG98a0rdwSKxxogInJovAye7JJfA+CxrwF308bacwRMgIiIiOxRs8GlT+fgnwDYBEZkH9gFiIjIKEyAiMix8ccwkUNiAkREREQOhwkQERERORwmQER2IKxpDQCAs0YtcyREpHxs9wV4FRiRXfCr5Ylhc9rDw0cjdyg2R+3rKncIRCQDJkBEdqJ6kIfcIdgkNZNGIofEJjAicmwcQoDIITEBIiIiIofDBIiIiMiRcCRoAEyAiIiIyAExASIix8Yfw0QOiQkQERERORwmQERERORwmAARERE5gqaDCv/v+KK8cSgEB0IkIiJyBI99DQz4AHD1ljsSRWANEBERkSOQJCY/xTABIiIiIofDBIiIiIgcDhMgInJokooDARE5IiZAROSQvPuEwynAHV5da8odChHJgFeBEZFD8u4VBu9eYXKHQUQyYQ0QERERORwmQERERORwmAARERGRw2ECRERERA6HCRARERE5HCZARERE5HCYABEREZHDYQJEREREDocJEBERETkcJkBERETkcJgAERERkcNhAkREREQOhwkQERERORwmQERERORwnOQOwJqEEACAtLQ0mSMhIiKiyio6bxedx83BoRKg9PR0AEBoaKjMkRAREZGx0tPT4ePjY5Z1ScKc6ZTC6XQ6XL9+HV5eXpAkyeC1tLQ0hIaG4sqVK/D29pYpQmVhmRhieZTEMjHE8iiJZWKI5VFSZcpECIH09HSEhIRApTJP7x2HqgFSqVSoVatWufN4e3vzQ3kPlokhlkdJLBNDLI+SWCaGWB4lVVQm5qr5KcJO0ERERORwmAARERGRw2EC9B+NRoPZs2dDo9HIHYpisEwMsTxKYpkYYnmUxDIxxPIoSa4ycahO0EREREQAa4CIiIjIATEBIiIiIofDBIiIiIgcDhMgIiIicjh2kwB9+umnqFOnDlxdXdGmTRvs2rWr3Pnj4uLQpk0buLq6om7duvjss88MXj916hQGDx6M2rVrQ5IkxMTElFjHnDlzIEmSwSMoKMicu1UlcpQJAFy7dg1PP/00atSoAXd3d7Rs2RJ//fWXuXbLZHKUR9Fr9z5efPFFc+6ayeQok4KCArz22muoU6cO3NzcULduXcydOxc6nc6cu2YSOcojPT0dkyZNQnh4ONzc3NCpUyccPHjQnLtVJeYuk6VLl6Jr166oXr06qlevjgceeAAHDhyo8natRY7y+OOPPzBgwACEhIRAkiRs2LDB3LtVJXKUyYIFC3DffffBy8sLAQEBeOSRR3DmzBnjAhd2YO3atcLZ2VksXbpUxMfHi4kTJwoPDw/x77//ljr/hQsXhLu7u5g4caKIj48XS5cuFc7OzuL777/Xz3PgwAExbdo0sWbNGhEUFCTef//9EuuZPXu2aNq0qbhx44b+kZiYaKndNIpcZZKcnCzCw8PFiBEjxP79+8XFixfFjh07xLlz5yy1q5UiV3kkJiYafD62b98uAIjff//dQntaeXKVybx580SNGjXExo0bxcWLF8V3330nPD09RUxMjKV2tVLkKo8nnnhCNGnSRMTFxYmzZ8+K2bNnC29vb3H16lVL7WqlWaJMhg0bJj755BNx5MgRcfr0aTFy5Ejh4+NjsL/Gbtda5CqPzZs3i1mzZon169cLAOLHH3+09K5Wmlxl0rdvX7Fs2TJx8uRJcfToUdG/f38RFhYmMjIyKh27XSRA7dq1E+PHjzeY1rhxYzFjxoxS558+fbpo3LixwbRx48aJDh06lDp/eHh4mQlQZGSkSTFbmlxlEh0dLbp06WJa0BYkV3nca+LEiaJevXpCp9NVLnALkqtM+vfvL0aNGmUwbdCgQeLpp582Inrzk6M8srKyhFqtFhs3bjSYHhkZKWbNmmXkHpifpctECCEKCgqEl5eXWLFihcnbtRa5yqM4pSVASigTIQp/bAIQcXFxlY7d5pvA8vLy8Ndff6FPnz4G0/v06YM9e/aUuszevXtLzN+3b18cOnQI+fn5Rm3/7NmzCAkJQZ06dTB06FBcuHDBuB2wADnL5Oeff0bbtm3x+OOPIyAgAK1atcLSpUuN3wkzkvszUjyO1atXY9SoUSVuxmttcpZJly5dsHPnTvzzzz8AgGPHjmH37t148MEHjdwL85GrPAoKCqDVauHq6mow3c3NDbt37zZiD8zPWmWSlZWF/Px8+Pr6mrxda5CrPJRMSWVy584dADCq3Gw+Abp9+za0Wi0CAwMNpgcGBiIhIaHUZRISEkqdv6CgALdv3670ttu3b4+VK1di69atWLp0KRISEtCpUyckJSUZvyNmJGeZXLhwAUuWLEGDBg2wdetWjB8/HhMmTMDKlSuN3xEzkbM8ituwYQNSU1MxYsQIk5Y3JznLJDo6Gk8++SQaN24MZ2dntGrVCpMmTcKTTz5p/I6YiVzl4eXlhY4dO+Ktt97C9evXodVqsXr1auzfvx83btwwbWfMxFplMmPGDNSsWRMPPPCAydu1BrnKQ8mUUiZCCEyZMgVdunRBs2bNKh2/3dwN/t5f1EKIcn9llzZ/adPLExUVpf+7efPm6NixI+rVq4cVK1ZgypQplV6PpchRJjqdDm3btsX8+fMBAK1atcKpU6ewZMkSPPvss5VejyXIUR7FffXVV4iKikJISIhJy1uCHGWybt06rF69Gt9++y2aNm2Ko0ePYtKkSQgJCcHw4cONiN785CiPVatWYdSoUahZsybUajVat26NYcOG4fDhw0ZEbjmWLJNFixZhzZo1iI2NLVELZux2rUWu8lAyucvkpZdewvHjx42uNbX5BMjPzw9qtbpEtpmYmFgiyywSFBRU6vxOTk6oUaOGybF4eHigefPmOHv2rMnrMAc5yyQ4OBhNmjQxmBYREYH169dXeh3mpoTPyL///osdO3bghx9+MHpZS5CzTF555RXMmDEDQ4cOBVD44+Hff//FggULZEuA5CyPevXqIS4uDpmZmUhLS0NwcDCGDBmCOnXqGL8jZmTpMvm///s/zJ8/Hzt27ECLFi2qtF1rkKs8lEwJZfLyyy/j559/xh9//IFatWoZFb/NN4G5uLigTZs22L59u8H07du3o1OnTqUu07FjxxLzb9u2DW3btoWzs7PJseTm5uL06dMIDg42eR3mIGeZdO7cucSliP/88w/Cw8MrvQ5zU8JnZNmyZQgICED//v2NXtYS5CyTrKwsqFSGhx61Wi3rZfBK+Ix4eHggODgYKSkp2Lp1KwYOHGj0OszJkmXy7rvv4q233sKWLVvQtm3bKm/XGuQqDyWTs0yEEHjppZfwww8/4LfffjPtB0Olu0srWNFleF999ZWIj48XkyZNEh4eHuLSpUtCCCFmzJghnnnmGf38RZfhTZ48WcTHx4uvvvqqxGV4ubm54siRI+LIkSMiODhYTJs2TRw5ckScPXtWP8/UqVNFbGysuHDhgti3b5946KGHhJeXl367cpKrTA4cOCCcnJzE22+/Lc6ePSu++eYb4e7uLlavXm29nS+FXOUhhBBarVaEhYWJ6Oho6+xsJclVJsOHDxc1a9bUXwb/ww8/CD8/PzF9+nTr7Xwp5CqPLVu2iF9//VVcuHBBbNu2TURGRop27dqJvLw86+18GSxRJgsXLhQuLi7i+++/NxgiIj09vdLblYtc5ZGenq7/HAEQ7733njhy5IjswwIIIV+ZPP/888LHx0fExsYazJOVlVXp2O0iARJCiE8++USEh4cLFxcX0bp1a4NL4YYPHy66d+9uMH9sbKxo1aqVcHFxEbVr1xZLliwxeP3ixYsCQIlH8fUMGTJEBAcHC2dnZxESEiIGDRokTp06ZcndNIocZSKEEL/88oto1qyZ0Gg0onHjxuKLL76w1C4aRa7y2Lp1qwAgzpw5Y6ldM5kcZZKWliYmTpwowsLChKurq6hbt66YNWuWyM3NteSuVooc5bFu3TpRt25d4eLiIoKCgsSLL74oUlNTLbmbRjF3mYSHh5daJrNnz670duUkR3n8/vvvpc4zfPhwC+5p5clRJqW9DkAsW7as0nFL/62IiIiIyGHYfB8gIiIiImMxASIiIiKHwwSIiIiIHA4TICIiInI4TICIiIjI4TABIiIiIofDBIiIiIgcDhMgItKLjY2FJElITU0FACxfvhzVqlUzmOeLL75AaGgoVCoVYmJiMGfOHLRs2bJK27106RIkScLRo0ertJ7y1nnvvhEpRY8ePSBJksGj6F55ZUlPT8ekSZMQHh4ONzc3dOrUCQcPHiwx3+nTp/Hwww/Dx8cHXl5e6NChAy5fvqx//fz583j00Ufh7+8Pb29vPPHEE7h586bBOg4fPozevXujWrVqqFGjBsaOHYuMjIwS21q+fDlatGgBV1dXBAUF4aWXXjKxRIBx48ZBkiTExMSYvI6KMAEishEjRoyAJEkYP358iddeeOEFSJKEESNGmHWbQ4YMwT///KN/npaWhpdeegnR0dG4du0axo4di2nTpmHnzp1m3W5pevTogUmTJpm8fKdOnXDjxg34+PiYL6hSzJkzR38Sc3Jygp+fH7p164aYmBjk5uYatS4mbfajR48eWL58eZmvjxkzBjdu3NA/Pv/883LX99xzz2H79u1YtWoVTpw4gT59+uCBBx7AtWvX9POcP38eXbp0QePGjREbG4tjx47h9ddf199VPTMzE3369IEkSfjtt9/w559/Ii8vDwMGDNDfm+/69et44IEHUL9+fezfvx9btmzBqVOnShxr3nvvPcyaNQszZszAqVOnsHPnTvTt29ekstqwYQP279+PkJAQk5avNOMHvSYiOQwfPlyEhoYKHx8fg/vdZGdni2rVqomwsLAqD41fNOR+SkpKqa+fOHFCABAXLlyo0nbuVXTLiCNHjpQ5T/fu3cXEiRPNuk5LmD17tmjatKm4ceOGuHbtmjh+/Lj48MMPRUBAgGjdurVIS0ur9Loqej/IdnTv3r3M2zQY+9nOysoSarVabNy40WB6ZGSkmDVrlv75kCFDxNNPP13merZu3SpUKpW4c+eOflpycrIAILZv3y6EEOLzzz8XAQEBQqvV6ucpuidZ0T3tkpOThZubm9ixY0e5cf/555+ia9euwtXVVdSqVUu8/PLLIiMjw2Ceq1evipo1a4qTJ0+K8PBw8f7775dfGFXAGiAiG9K6dWuEhYXhhx9+0E/74YcfEBoailatWhnMm5ubiwkTJiAgIACurq7o0qVLiSryzZs3o2HDhnBzc0PPnj1x6dIlg9eLN4EtX74czZs3BwDUrVsXkiTh0qVLpTaBLVu2DBEREXB1dUXjxo3x6aefGrx+4MABtGrVCq6urmjbti2OHDlidFlIkoQNGzYYTKtWrVqZv7LvrU1JSkrCk08+iVq1asHd3R3NmzfHmjVrDJbp0aMHJkyYgOnTp8PX1xdBQUGYM2dOhbE5OTkhKCgIISEhaN68OV5++WXExcXh5MmTWLhwoX6+1atXo23btvDy8kJQUBCGDRuGxMREAIVNeD179gQAVK9e3aCGTwiBRYsWoW7dunBzc0NkZCS+//57/XpTUlLw1FNPwd/fH25ubmjQoAGWLVtWYdwkn2+++QZ+fn5o2rQppk2bhvT09DLnLSgogFar1dfkFHFzc8Pu3bsBADqdDps2bULDhg3Rt29fBAQEoH379gbfmdzcXEiSBI1Go5/m6uoKlUqlX09ubi5cXFygUqkMtgNAP8/27duh0+lw7do1REREoFatWnjiiSdw5coV/TInTpxA3759MWjQIBw/fhzr1q3D7t27DZrJdDodnnnmGbzyyito2rSpsUVoNCZARDZm5MiRBiezr7/+GqNGjSox3/Tp07F+/XqsWLEChw8fRv369dG3b18kJycDAK5cuYJBgwbhwQcfxNGjR/Hcc89hxowZZW53yJAh2LFjB4DCBObGjRsIDQ0tMd/SpUsxa9YsvP322zh9+jTmz5+P119/HStWrABQWO3+0EMPoVGjRvjrr78wZ84cTJs2rUplYoqcnBy0adMGGzduxMmTJzF27Fg888wz2L9/v8F8K1asgIeHB/bv349FixZh7ty52L59u9Hba9y4MaKiogyS17y8PLz11ls4duwYNmzYgIsXL+qTnNDQUKxfvx4AcObMGdy4cQMffPABAOC1117DsmXLsGTJEpw6dQqTJ0/G008/jbi4OADA66+/jvj4ePz66684ffo0lixZAj8/P1OKiazgqaeewpo1axAbG4vXX38d69evx6BBg8qc38vLCx07dsRbb72F69evQ6vVYvXq1di/fz9u3LgBAEhMTERGRgbeeecd9OvXD9u2bcOjjz6KQYMG6T8nHTp0gIeHB6Kjo5GVlYXMzEy88sor0Ol0+vX06tULCQkJePfdd5GXl4eUlBS8+uqrAKCf58KFC9DpdJg/fz5iYmLw/fffIzk5Gb1790ZeXh4A4N1338WwYcMwadIkNGjQAJ06dcKHH36IlStXIicnBwCwcOFCODk5YcKECZYp6HtZrG6JiMxq+PDhYuDAgeLWrVtCo9GIixcvikuXLglXV1dx69YtMXDgQH0TWEZGhnB2dhbffPONfvm8vDwREhIiFi1aJIQQYubMmSIiIkLodDr9PNHR0QZNLsuWLRM+Pj7614uqvi9evKifNnv2bBEZGal/HhoaKr799luD2N966y3RsWNHIURhlbqvr6/IzMzUv75kyRKjm8AAiB9//NFgHh8fH30zw71NYJVpTnrwwQfF1KlTDbbZpUsXg3nuu+8+ER0dXeY67i2P4qKjo4Wbm1uZyx44cEAAEOnp6WXGnJGRIVxdXcWePXsMlh09erR48sknhRBCDBgwQIwcObLM7ZDlvf3228LDw0P/UKlUQqPRGEz7448/Sl320KFDAoD466+/ylz/uXPnRLdu3QQAoVarxX333SeeeuopERERIYQQ4tq1awKA/jNRZMCAAWLo0KH651u3bhV169YVkiQJtVotnn76adG6dWvx/PPP6+f55ptvRGBgoFCr1cLFxUVMmzZNBAYGioULF+r3FYDYunWrfpnExEShUqnEli1bhBBCNGnSRLi4uBjsv7u7uwAg4uPjxaFDh0RgYKC4du2afh2WbgJzsk6aRUTm4ufnh/79+2PFihUQQqB///4lft2fP38e+fn56Ny5s36as7Mz2rVrh9OnTwMovDqkQ4cOkCRJP0/Hjh2rFNutW7dw5coVjB49GmPGjNFPLygo0Hc+Pn36NCIjI+Hu7m627ZpCq9XinXfewbp163Dt2jXk5uYiNzcXHh4eBvO1aNHC4HlwcLC+mcpYQgiD8j5y5AjmzJmDo0ePIjk5Wd/x9PLly2jSpEmp64iPj0dOTg569+5tMD0vL0/fDPr8889j8ODBOHz4MPr06YNHHnkEnTp1MilmMs348ePxxBNP6J8/9dRTGDx4sEHNTs2aNUtdtnXr1nB2dsbZs2fRunXrUuepV68e4uLikJmZibS0NAQHB2PIkCGoU6cOgMLjhJOTU4nPUUREhL7pCgD69OmD8+fP4/bt23ByckK1atUQFBSkXw8ADBs2DMOGDcPNmzfh4eEBSZLw3nvv6ecJDg4GAINt+fv7w8/PT3/FmU6nw7hx40qt3QkLC8Onn36KxMREhIWF6adrtVpMnToVMTExJZrnzYEJEJENGjVqlL7t/JNPPinxuhACAAxOtkXTi6YVzWNORSfwpUuXon379gavqdVqs25XkqQS68rPz6/08osXL8b777+PmJgYNG/eHB4eHpg0aZK+yr6Is7Nzie0W7aexTp8+rT9pFF2B06dPH6xevRr+/v64fPky+vbtWyKG4oq2vWnTphIn0KK+HFFRUfj333+xadMm7NixA/fffz9efPFF/N///Z9JcZPxfH194evrq3/u5uaGgIAA1K9fv8JlT506hfz8fH1iUR4PDw94eHggJSUFW7duxaJFiwAALi4uuO+++3DmzBmD+f/55x+Eh4eXWE/Rj6jffvsNiYmJePjhh0vMExgYCKCw2d3V1VWfhBf90Dpz5gxq1aoFAEhOTsbt27f122rdujVOnTpV5v4/88wzeOCBBwym9e3bF8888wxGjhxZYTmYggkQkQ3q16+f/iRZ2qWm9evXh4uLC3bv3o1hw4YBKEwODh06pL+UvEmTJiU6Ee/bt69KcQUGBqJmzZq4cOECnnrqqVLnadKkCVatWoXs7Gx9Z0pTtuvv76/vgwAAZ8+eRVZWVqWX37VrFwYOHIinn34aQGFicfbsWURERBgdS2X8/fff2LJlC2bOnKl/fvv2bbzzzjv6vlSHDh0yWMbFxQVA4S/hIk2aNIFGo8Hly5fRvXv3Mrfn7++PESNGYMSIEejatSteeeUVJkAKdP78eXzzzTd48MEH4efnh/j4eEydOhWtWrUyqMG9//778eijj+p/+GzduhVCCDRq1Ajnzp3DK6+8gkaNGhkkC6+88gqGDBmCbt26oWfPntiyZQt++eUXxMbG6ucpumDB398fe/fuxcSJEzF58mQ0atRIP8/HH3+MTp06wdPTE9u3b8crr7yCd955R3+BRMOGDTFw4EBMnDgRX3zxBby9vTFz5kw0btxY35E/OjoaHTp0wIsvvogxY8bAw8MDp0+fxvbt2/HRRx+hRo0aqFGjhkHZODs7IygoyCAWc2ICRGSD1Gq1vimrqGalOA8PDzz//PN45ZVX4Ovri7CwMCxatAhZWVkYPXo0gMIq+sWLF2PKlCkYN24c/vrrr3LHKamsOXPmYMKECfD29kZUVBRyc3Nx6NAhpKSkYMqUKRg2bBhmzZqF0aNH47XXXsOlS5dMOjH36tULH3/8MTp06ACdTofo6OgStTXlqV+/PtavX489e/agevXqeO+995CQkGCWBKigoAAJCQnQ6XRISkpCbGws5s2bh5YtW+KVV14BUFjt7+Ligo8++gjjx4/HyZMn8dZbbxmsJzw8HJIkYePGjXjwwQfh5uYGLy8vTJs2DZMnT4ZOp0OXLl2QlpaGPXv2wNPTE8OHD8cbb7yBNm3aoGnTpsjNzcXGjRstlthR1bi4uGDnzp344IMPkJGRgdDQUPTv3x+zZ882+G4XNVMVuXPnDmbOnImrV6/C19cXgwcPxttvv23wHXj00Ufx2WefYcGCBZgwYQIaNWqE9evXo0uXLvp5zpw5g5kzZyI5ORm1a9fGrFmzMHnyZIMYDxw4gNmzZyMjIwONGzfG559/jmeeecZgnpUrV2Ly5Mno378/VCoVunfvji1btujjadGiBeLi4jBr1ix07doVQgjUq1cPQ4YMMWt5GsVivYuIyKyKOkGXpXgnaCEKxwd6+eWXhZ+fn9BoNKJz587iwIEDBsv88ssvon79+kKj0YiuXbuKr7/+usqdoIUo7DTZsmVL4eLiIqpXry66desmfvjhB/3re/fuFZGRkcLFxUW0bNlSrF+/vsJO0F27djXooHzt2jXRp08f4eHhIRo0aCA2b95sVCfopKQkMXDgQOHp6SkCAgLEa6+9Jp599lmDMi5tfJZ7y/les2fPFgD0nVN9fX1Fly5dxPvvvy9ycnIM5v32229F7dq1hUajER07dhQ///xziXKYO3euCAoKEpIk6ber0+nEBx98IBo1aiScnZ2Fv7+/6Nu3r4iLixNCFHY6j4iIEG5ubsLX11cMHDjQ7GM3Edk6SQgLdAQgIjKzxo0b47nnnpPlknkisj9sAiMiRUtMTMSvv/6KM2fO4P7775c7HCKyE0yAiEjR+vXrh5SUFHz44YclRrsmIjIVm8CIiIjI4fBWGERERORwmAARERGRw2ECRERERA6HCRARERE5HCZARERE5HCYABEREZHDYQJEREREDocJEBERETkcJkBERETkcP4fdPMah90OaKkAAAAASUVORK5CYII=", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "# Open one of the downloaded fits file as an example\n", + "filename = file[1]['Local Path'] \n", + "hdul = fits.open(filename)\n", + "\n", + "# Print the information of the file to check if needed\n", + "#print(hdul.info())\n", + "\n", + "# Load the corresponding fluxes for each exposure\n", + "data = [(hdul[2+i].data)['FLUX'] for i in range(len(time))]\n", + "\n", + "# Create a time array using the start and end time (in Modified Julian Date (MJD)) for each exposure: \n", + "time = hdul['INT_TIMES'].data\n", + "time_start = [time[i]['int_start_MJD_UTC'] for i in range(len(time))]\n", + "time_end = [time[i]['int_end_MJD_UTC'] for i in range(len(time))]\n", + "time_arr = [np.linspace(time_start[i], time_end[i], len(data[i]))for i in range(len(time))]\n", + "\n", + "# Create a list of stingray light curve objects for analysis\n", + "lc = [Lightcurve(time=time_arr[i], counts=data[i]) for i in range(len(time))]\n", + "\n", + "for i in range(len(lc)):\n", + " lc[i].plot()\n", + "\n", + "plt.xlabel('Modified Julian Dates')\n", + "plt.title('Example JWST Grism light curves')" + ] + }, + { + "cell_type": "markdown", + "id": "909fea8c-0714-47e6-93eb-961975d9b54d", + "metadata": {}, + "source": [ + "Note: This is only an example using an arbitrary NIRCam grism time-series observation. For specific light curves and quantites (e.g. uncertainties of the fluxes), please amend the cell above. " + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "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.11.7" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}