diff --git a/README.md b/README.md index 7e3ff04..49a2b23 100644 --- a/README.md +++ b/README.md @@ -11,8 +11,30 @@ The main goal of this project is to standardize the handling of stellar spectra 3) Framework to mask telluric features from stellar spectra 4) Implement utility methods to interpolate, normalize and smooth stellar spectra -For information on installation, usage, and contribution check the official [documentation](https://kamuish.github.io/ASTRA/) +For information on installation, usage, and contribution check the official [documentation](https://kamuish.github.io/ASTRA/). +If you use it, please cite the paper: + +``` +@ARTICLE{2026JOSS...11.9413S, + author = {{Silva}, Andr{\'e} and {Faria}, J. and {Santos}, Nuno and {Sousa}, S{\'e}rgio and {Viana}, Pedro and {Martins}, J.}, + title = "{ASTRA: A Python Package for Cross-Instrument Stellar and Telluric Template Construction}", + journal = {The Journal of Open Source Software}, + keywords = {astronomy, Python, Cython, Instrumentation and Methods for Astrophysics, Earth and Planetary Astrophysics, Solar and Stellar Astrophysics}, + year = 2026, + month = jan, + volume = {11}, + number = {117}, + eid = {9413}, + pages = {9413}, + doi = {10.21105/joss.09413}, +archivePrefix = {arXiv}, + eprint = {2601.10439}, + primaryClass = {astro-ph.IM}, + adsurl = {https://ui.adsabs.harvard.edu/abs/2026JOSS...11.9413S}, + adsnote = {Provided by the SAO/NASA Astrophysics Data System} +} +``` ## Installation diff --git a/docs/index.md b/docs/index.md index c9a621b..d812284 100644 --- a/docs/index.md +++ b/docs/index.md @@ -66,6 +66,25 @@ Depending on the python version that is in use, the installation of TelFit might ## Changelog +### V1.2.7 (current) + +1) Added HARPS-N support +2) Updated plotting capabilities of DataClass + + +### V1.2.6 + + +1. Added preliminary support for PoET data +2. Updated spawn method of multiprocessing pools +3. Improved frame-rejection tools + + +### V1.2.4 + +1. Update init of ESO based pipelines + + ### V1.2.3 (5th December 2025) 1. Addition of new routines to oversample the stellar template diff --git a/docs/spectra/access_spectra.ipynb b/docs/spectra/access_spectra.ipynb index 908f607..e5780bc 100644 --- a/docs/spectra/access_spectra.ipynb +++ b/docs/spectra/access_spectra.ipynb @@ -11,48 +11,9 @@ }, { "cell_type": "code", - "execution_count": 27, + "execution_count": 5, "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "\u001b[32m2025-04-14 21:53:35.443\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.utils.UserConfigs\u001b[0m:\u001b[36mreceive_user_inputs\u001b[0m:\u001b[36m216\u001b[0m - \u001b[34m\u001b[1mGenerating internal configs of - \u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.444\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.utils.UserConfigs\u001b[0m:\u001b[36mreceive_user_inputs\u001b[0m:\u001b[36m221\u001b[0m - \u001b[1mChecking for any parameter that will take default value\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.445\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.utils.UserConfigs\u001b[0m:\u001b[36mreceive_user_inputs\u001b[0m:\u001b[36m228\u001b[0m - \u001b[34m\u001b[1mConfiguration using the default value: DISK_SAVE_MODE.DISABLED\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.446\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.utils.UserConfigs\u001b[0m:\u001b[36mreceive_user_inputs\u001b[0m:\u001b[36m228\u001b[0m - \u001b[34m\u001b[1mConfiguration using the default value: WORKING_MODE.ONE_SHOT\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.447\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m126\u001b[0m - \u001b[1mDataClass opening 3 files from a list/tuple\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.448\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m253\u001b[0m - \u001b[1mCreating frame from: /home/amiguel/spectra_collection/ESPRESSO/proxima/r.ESPRE.2019-07-03T01:43:39.634_S2D_A.fits\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.449\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[36mASTRA.Components.SpectrumComponent\u001b[0m:\u001b[36mregenerate_order_status\u001b[0m:\u001b[36m96\u001b[0m - \u001b[33m\u001b[1mResetting order status of Frame - ESPRESSO\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.466\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36massess_bad_orders\u001b[0m:\u001b[36m711\u001b[0m - \u001b[34m\u001b[1mRejecting spectral orders\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.468\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36massess_bad_orders\u001b[0m:\u001b[36m741\u001b[0m - \u001b[1mFrame -553975442344993652 rejected 37 orders for having SNR smaller than 5: 0-36\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.469\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m253\u001b[0m - \u001b[1mCreating frame from: /home/amiguel/spectra_collection/ESPRESSO/proxima/r.ESPRE.2019-07-14T02:07:49.063_S2D_A.fits\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.471\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[36mASTRA.Components.SpectrumComponent\u001b[0m:\u001b[36mregenerate_order_status\u001b[0m:\u001b[36m96\u001b[0m - \u001b[33m\u001b[1mResetting order status of Frame - ESPRESSO\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.488\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36massess_bad_orders\u001b[0m:\u001b[36m711\u001b[0m - \u001b[34m\u001b[1mRejecting spectral orders\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.490\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36massess_bad_orders\u001b[0m:\u001b[36m741\u001b[0m - \u001b[1mFrame -2975062627551510774 rejected 15 orders for having SNR smaller than 5: 0-12, 32-33\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.492\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m253\u001b[0m - \u001b[1mCreating frame from: /home/amiguel/spectra_collection/ESPRESSO/proxima/r.ESPRE.2019-07-20T01:43:40.032_S2D_A.fits\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.493\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[36mASTRA.Components.SpectrumComponent\u001b[0m:\u001b[36mregenerate_order_status\u001b[0m:\u001b[36m96\u001b[0m - \u001b[33m\u001b[1mResetting order status of Frame - ESPRESSO\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.508\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36massess_bad_orders\u001b[0m:\u001b[36m711\u001b[0m - \u001b[34m\u001b[1mRejecting spectral orders\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.511\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36massess_bad_orders\u001b[0m:\u001b[36m741\u001b[0m - \u001b[1mFrame 3535164837059905149 rejected 16 orders for having SNR smaller than 5: 0-13, 32-33\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.513\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m154\u001b[0m - \u001b[34m\u001b[1mSelected 3 observations from disk\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.514\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36m_collect_MetaData\u001b[0m:\u001b[36m369\u001b[0m - \u001b[1mCollecting MetaData from the observations\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.515\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[36mASTRA.data_objects.Target\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m73\u001b[0m - \u001b[33m\u001b[1mTarget dictionary not found in \u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.515\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.data_objects.Target\u001b[0m:\u001b[36mclean_targ_list\u001b[0m:\u001b[36m98\u001b[0m - \u001b[34m\u001b[1mParsing through loaded OBJECTs\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.517\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.Target\u001b[0m:\u001b[36m__init__\u001b[0m:\u001b[36m92\u001b[0m - \u001b[1mValidated target to be V V645 Cen\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.517\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36mshow_loadedData_table\u001b[0m:\u001b[36m885\u001b[0m - \u001b[1m--------------------------------------------------------------------\n", - "--------------------------------------------------------------------\n", - " subInstrument Total OBS Valid OBS [warnings] INVALID OBS \n", - "--------------------------------------------------------------------\n", - " ESPRESSO18 0 0 [0] 0 \n", - " ESPRESSO19 3 3 [0] 0 \n", - " Total 3 3 [0] 0 \n", - "--------------------------------------------------------------------\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.518\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36mload_instrument_extra_information\u001b[0m:\u001b[36m894\u001b[0m - \u001b[1mChecking if the instrument has extra data to load\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.519\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36mload_instrument_extra_information\u001b[0m:\u001b[36m901\u001b[0m - \u001b[1mCurrent instrument does not need to load anything from the outside\u001b[0m\n" - ] - } - ], + "outputs": [], "source": [ "from pathlib import Path\n", "\n", @@ -82,37 +43,20 @@ }, { "cell_type": "code", - "execution_count": 28, + "execution_count": 6, "metadata": {}, "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "\u001b[32m2025-04-14 21:53:35.531\u001b[0m | \u001b[34m\u001b[1mDEBUG \u001b[0m | \u001b[36mASTRA.base_models.Frame\u001b[0m:\u001b[36madd_to_status\u001b[0m:\u001b[36m574\u001b[0m - \u001b[34m\u001b[1mUpdating Frame (r.ESPRE.2019-07-03T01:43:39.634_S2D_A.fits) status to USER_BLOCKED: Filename rejected\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.532\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36mreject_observations\u001b[0m:\u001b[36m292\u001b[0m - \u001b[1mUser conditions removed 1 / 3 frames\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.533\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36mreject_observations\u001b[0m:\u001b[36m298\u001b[0m - \u001b[1mUpdated observation on disk!\u001b[0m\n", - "\u001b[32m2025-04-14 21:53:35.534\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[36mASTRA.data_objects.DataClass\u001b[0m:\u001b[36mshow_loadedData_table\u001b[0m:\u001b[36m885\u001b[0m - \u001b[1m--------------------------------------------------------------------\n", - "--------------------------------------------------------------------\n", - " subInstrument Total OBS Valid OBS [warnings] INVALID OBS \n", - "--------------------------------------------------------------------\n", - " ESPRESSO18 0 0 [0] 0 \n", - " ESPRESSO19 3 2 [0] 1 \n", - " Total 3 2 [0] 1 \n", - "--------------------------------------------------------------------\u001b[0m\n" - ] - }, { "name": "stdout", "output_type": "stream", "text": [ "Data Class from holding 0 OBS from ESPRESSO18, 3 OBS from ESPRESSO19\n", "Loaded sub-Instruments ['ESPRESSO19']\n", - "Current frameIDs [-553975442344993652, -2975062627551510774, 3535164837059905149]\n", + "Current frameIDs [565306508, 4035951396, 933504894]\n", "Current (rejected) frameIDs []\n", "---//---\n", - "Current frameIDs [-2975062627551510774, 3535164837059905149]\n", - "Current (rejected) frameIDs [-553975442344993652]\n" + "Current frameIDs [4035951396, 933504894]\n", + "Current (rejected) frameIDs [565306508]\n" ] } ], @@ -146,7 +90,7 @@ }, { "cell_type": "code", - "execution_count": 29, + "execution_count": 7, "metadata": {}, "outputs": [], "source": [ @@ -163,7 +107,7 @@ }, { "cell_type": "code", - "execution_count": 30, + "execution_count": 8, "metadata": {}, "outputs": [ { @@ -172,7 +116,7 @@ "[-21.3856633997537, -21.3845961288029]" ] }, - "execution_count": 30, + "execution_count": 8, "metadata": {}, "output_type": "execute_result" } @@ -192,7 +136,7 @@ }, { "cell_type": "code", - "execution_count": 31, + "execution_count": 9, "metadata": {}, "outputs": [ { @@ -201,7 +145,7 @@ "[-21.3856633997537, -21.3845961288029]" ] }, - "execution_count": 31, + "execution_count": 9, "metadata": {}, "output_type": "execute_result" } @@ -210,6 +154,33 @@ "data.collect_KW_observations(\"ESO QC CCF RV\", subInstruments=available_subInstruments, from_header=True)" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Or make a plot, for all available sub-Instruments:" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAn4AAAHVCAYAAABv4/bQAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjEsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvc2/+5QAAAAlwSFlzAAAPYQAAD2EBqD+naQAAZvlJREFUeJzt3XlcVOX+B/DPDMug4ICsMygILok7pbJYuUGBS25YaZiiXM1CMzWv2iJZqWmWmpqaJuhVy63rTVIMxV1cLmYqKaWBCzAoIQOorPP8/rhyfk7sCA44n/frdV63ec5znvM9j3X7dFaZEEKAiIiIiJ54ckMXQERERESPB4MfERERkZFg8CMiIiIyEgx+REREREaCwY+IiIjISDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIqMRGRkJmUxW7nLy5EkAQG5uLsLDw9GxY0dYWlrCzs4Onp6emDJlClJTU6XxPvroI73tGzdujPbt2+ODDz5AdnZ2ufs1NTVFs2bNEBISgpSUlFJ19u7du9waPTw89PpeuHABw4cPR4sWLWBhYYFmzZrhhRdewPLly/X6FRQUYNmyZXj66aehVCphY2ODDh06YMKECbh8+XKpGhISEjBq1Cg0a9YMCoUCzs7OCA4ORkJCQqm+JfMVGBgIW1tbyGQyREZGlvvnsG3bNvj4+MDGxgZ2dnbo1asXfvrpp0r//Ijo0ZkaugAiosft448/hru7e6n21q1bo7CwED179sTly5cxZswYTJ48Gbm5uUhISMCWLVswdOhQODs76223atUqWFlZITc3Fz///DPmzZuH2NhYHD9+HDKZrNR+8/LycPLkSURGRuLYsWO4ePEiLCws9MZs3rw5FixYUKpGa2tr6a9PnDiBPn36wNXVFePHj4dKpcKNGzdw8uRJLFu2DJMnT5b6BgUFYe/evRg5ciTGjx+PwsJCXL58GVFRUejRo4deoPzhhx8wcuRI2NraIjQ0FO7u7khOTsa3336LHTt24Pvvv8fQoUOl/hkZGfj444/h6uqKLl264NChQ+XO/fLly/H2229jwIAB+Oyzz5CXl4fIyEgMHDgQO3fuxLBhwyr50yOiRyKIiIxERESEACDOnDlTbp9t27YJAGLz5s2l1t2/f19otVrpd3h4uAAgbt++rddv2LBhAoA4ceJEhfudOXOmACC2bt2q196rVy/RoUOHSo+nf//+wsHBQdy5c6fUuvT0dOmvT58+LQCIefPmlepXVFQkMjIypN9XrlwRjRs3Fh4eHuLWrVt6fW/fvi08PDyEpaWluHr1qtSel5cn0tLShBBCnDlzRgAQERERZdbcpk0b0b17d6HT6aQ2rVYrrKysxKBBgyo9ZiJ6NLzUS0T0kKtXrwIAnn322VLrLCwsoFQqKx2jb9++AICkpKQK+z3//PN6+6xJrR06dICNjU2pdY6Ojnr9UM4xmZiYwM7OTvr9+eef4969e/jmm2/g4OCg19fe3h5r1qzB3bt3sWjRIqldoVBApVJVqebs7Gw4OjrqnQlVKpWwsrJCo0aNqjQGEdUcgx8RGR2tVouMjAy95a+//gIAtGjRAgCwceNGCCFqNH5J0Ho4UJUlOTkZANC0adNS64qLi0vVmJGRgbt370p9WrRogfj4eFy8eLHC/ZQc0+bNm1FUVFRh3927d8PNzU0KpX/Xs2dPuLm51fievN69eyM6OhrLly9HcnIyLl++jLCwMGi1WkyZMqVGYxJRNRj6lCMR0eNScsm1rEWhUAghhLh3755o27atACBatGghQkJCxLfffqt36bREyaXexMREcfv2bZGUlCTWrFkjFAqFcHJyEnfv3tXb7/79+8Xt27fFjRs3xI4dO4SDg4NQKBTixo0beuP26tWr3DrfeOMNqd/PP/8sTExMhImJifD19RX//Oc/xb59+0RBQYHeeDqdThrTyclJjBw5UqxcuVJcu3ZNr19WVpYAIAYPHlzhPA4aNEgAENnZ2aXWVXapNz09Xfj5+ekdk729vXRZnIjqFh/uICKjs3LlSjz11FN6bSYmJgCARo0a4dSpU5g3bx62bduGyMhIREZGQi6X46233sLixYuhUCj0tm3btq3e7w4dOmDDhg1o3LixXru/v7/ebzc3N2zatAnNmzcvVaObmxvWrl1bqv3hvi+88ALi4uKwYMEC7Nu3D3FxcVi0aBEcHBywbt06DBo0CAAgk8mwb98+LF68GJs2bcJ3332H7777DmFhYXjllVewZs0a2NjYICcnBwDQpEmTCuevZH12dnalff+ucePGaNu2LZo3b46BAwciJycHS5YswbBhw3D06FG0bt26WuMRUfUw+BGR0fHy8kK3bt3KXW9tbY1FixZh0aJFuHbtGg4cOIDFixdjxYoVsLa2xqeffqrXf+fOnVAqlTAzM0Pz5s3RqlWrMsctCZxarRbr16/HkSNHSoXIEpaWlqWCYlm6d++OH374AQUFBfj111/x73//G0uWLMHw4cNx7tw5tG/fHnhwH97777+P999/H2lpaTh8+DCWLVuGbdu2wczMDJs2bZJCXEkALE9VA2JZXn75ZZiammL37t1S2+DBg9GmTRu8//772Lp1a7XHJKKq4z1+REQVaNGiBcaNG4fjx4/DxsYGmzdvLtWnZ8+e8Pf3R69evcoNfXgQOP39/REUFIQff/wRHTt2xGuvvYbc3NxHrtPc3Bzdu3fH/PnzsWrVKhQWFmL79u1l9lWr1RgxYgSOHDmCNm3aYNu2bSgqKoK1tTXUajXOnz9f4b7Onz+PZs2aVelBl4f9+eefiI6Ols5ElrC1tcVzzz2H48ePV2s8Iqo+Bj8ioipo2rQpWrVqhbS0tFoZz8TEBAsWLEBqaipWrFhRK2OWKDmbWVmtZmZm6Ny5MwoLC5GRkQEAGDhwIJKSknDs2LEytzl69CiSk5MxcODAateVnp4OPHhw5e8KCwsrffCEiB4dgx8R0UN+/fVXKQQ97Nq1a/jtt99K3c/3KHr37g0vLy8sXboUeXl51d7+4MGDZT55vGfPHuChew//+OMPXL9+vVS/rKwsxMXFoWnTptKrW2bMmIFGjRrhjTfekJ50LpGZmYmJEyeicePGmDFjRrXrbd26NeRyObZu3apX982bN3H06FE8/fTT1R6TiKqH9/gRkdHZu3dvmZ8p69GjB2JiYhAeHo5BgwbBx8cHVlZW+PPPP7F+/Xrk5+fjo48+qtVaZsyYgZdffhmRkZGYOHGi1K7VarFp06Yytxk1ahQAYPLkybh37x6GDh0KDw8PFBQU4MSJE9i6dSvc3NwwduxY4EGYfe2119CvXz88//zzsLW1RUpKCjZs2IDU1FQsXbpUerilTZs22LBhA4KDg9GpU6dSX+7IyMjAd999V+qS9ooVK5CVlSV90m737t24efOmVKe1tTUcHBwwbtw4rFu3Dn5+fhg2bBhycnLw9ddf4/79+5g9e3atzi0RlcHQjxUTET0uFb3OpeQVJH/++aeYM2eO8PHxEY6OjsLU1FQ4ODiIAQMGiNjYWL3xyvtyR3n7LeuLIcXFxaJVq1aiVatWoqioSIhKXufy8P9t7927V4wbN054eHgIKysrYW5uLlq3bi0mT56s9/qZ9PR08dlnn4levXoJtVotTE1NRdOmTUXfvn3Fjh07yqz5/PnzYuTIkUKtVgszMzOhUqnEyJEjxYULF8rs36JFi3LrTUpKkvoVFhaK5cuXC09PT2FlZSWsrKxEnz59Ss0tEdUNmajpG0qJiIiIqEHhPX5ERERERoLBj4iIiMhIMPgRERERGQkGPyIiIiIjweBHREREZCT4Hr8GRKfTITU1FU2aNIFMJjN0OURERFRPCCGQk5MDZ2dnyOXln9dj8GtAUlNT4eLiYugyiIiIqJ66ceMGmjdvXu56Br8GpEmTJsCDP9TqfhydiIiInlzZ2dlwcXGRskJ5GPwakJLLu0qlksGPiIiISqnsVjA+3EFERERkJIw2+K1cuRJubm6wsLCAt7c3Tp8+XWH/7du3w8PDAxYWFujUqRP27Nmjtz4kJAQymUxvCQwM1OuTmZmJ4OBgKJVK2NjYIDQ0FLm5uXVyfERERER/Z5TBb+vWrZg2bRrCw8Nx9uxZdOnSBQEBAbh161aZ/U+cOIGRI0ciNDQUv/zyC4YMGYIhQ4bg4sWLev0CAwORlpYmLd99953e+uDgYCQkJCAmJgZRUVE4cuQIJkyYUKfHSkRERFRCJoQQhi7icfP29kb37t2xYsUK4MFrUlxcXDB58mTMmjWrVP9XX30Vd+/eRVRUlNTm4+MDT09PrF69Gnhwxi8rKwu7du0qc5+XLl1C+/btcebMGXTr1g0AEB0djf79++PmzZtwdnautO7s7GxYW1tDq9XyHj8iIiKSVDUjGN0Zv4KCAsTHx8Pf319qk8vl8Pf3R1xcXJnbxMXF6fUHgICAgFL9Dx06BEdHR7Rt2xZvvvkm/vrrL70xbGxspNAHAP7+/pDL5Th16lSZ+83Pz0d2drbeQkRERFRTRhf8MjIyUFxcDCcnJ712JycnaDSaMrfRaDSV9g8MDMTGjRtx4MABLFy4EIcPH0a/fv1QXFwsjeHo6Kg3hqmpKWxtbcvd74IFC2BtbS0tfIcfERERPQq+zqWWjBgxQvrrTp06oXPnzmjVqhUOHToEPz+/Go05e/ZsTJs2Tfpd8o4eIiIiopowujN+9vb2MDExQXp6ul57eno6VCpVmduoVKpq9QeAli1bwt7eHleuXJHG+PvDI0VFRcjMzCx3HIVCIb2zj+/uIyIiokdldMHP3NwcXbt2xYEDB6Q2nU6HAwcOwNfXt8xtfH199foDQExMTLn9AeDmzZv466+/oFarpTGysrIQHx8v9YmNjYVOp4O3t3ctHBkRERFRxYwu+AHAtGnTsHbtWmzYsAGXLl3Cm2++ibt372Ls2LEAgNGjR2P27NlS/ylTpiA6OhpffPEFLl++jI8++gj//e9/MWnSJABAbm4uZsyYgZMnTyI5ORkHDhzA4MGD0bp1awQEBAAA2rVrh8DAQIwfPx6nT5/G8ePHMWnSJIwYMaJKT/QSERERPSqjvMfv1Vdfxe3btzFnzhxoNBp4enoiOjpaeoDj+vXrkMv/PxP36NEDW7ZswQcffID33nsPbdq0wa5du9CxY0cAgImJCc6fP48NGzYgKysLzs7OePHFF/HJJ59AoVBI42zevBmTJk2Cn58f5HI5goKC8NVXXxlgBoiIiMgYGeV7/BoqvsePiIiIysL3+BERERGRHgY/IiIiIiPB4EdERERkJBj8iIiIiIwEgx8RERGRkWDwIyIiIjISDH5ERERERoLBj4iIiMhIMPgRERERGQkGPyIiIiIjweBHREREZCQY/IiIiIiMBIMfERERkZFg8CMiIiIyEgx+REREREaCwY+IiIjISDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIiIiMhIMfkRERERGgsGPiIiIyEgw+BEREREZCQY/IiIiIiNhtMFv5cqVcHNzg4WFBby9vXH69OkK+2/fvh0eHh6wsLBAp06dsGfPHmldYWEhZs6ciU6dOsHS0hLOzs4YPXo0UlNT9cZwc3ODTCbTWz777LM6O0YiIiKihxll8Nu6dSumTZuG8PBwnD17Fl26dEFAQABu3bpVZv8TJ05g5MiRCA0NxS+//IIhQ4ZgyJAhuHjxIgDg3r17OHv2LD788EOcPXsWP/zwAxITEzFo0KBSY3388cdIS0uTlsmTJ9f58RIREREBgEwIIQxdxOPm7e2N7t27Y8WKFQAAnU4HFxcXTJ48GbNmzSrV/9VXX8Xdu3cRFRUltfn4+MDT0xOrV68ucx9nzpyBl5cXrl27BldXV+DBGb933nkH77zzTpXqzM/PR35+vvQ7OzsbLi4u0Gq1UCqV1T5uIiIiejJlZ2fD2tq60oxgdGf8CgoKEB8fD39/f6lNLpfD398fcXFxZW4TFxen1x8AAgICyu0PAFqtFjKZDDY2Nnrtn332Gezs7PD000/j888/R1FRUbljLFiwANbW1tLi4uJSjSMlIiIi0mdq6AIet4yMDBQXF8PJyUmv3cnJCZcvXy5zG41GU2Z/jUZTZv+8vDzMnDkTI0eO1Evdb7/9Np555hnY2trixIkTmD17NtLS0vDll1+WOc7s2bMxbdo06XfJGT8iIiKimjC64FfXCgsL8corr0AIgVWrVumtezjEde7cGebm5njjjTewYMECKBSKUmMpFIoy24mIiIhqwugu9drb28PExATp6el67enp6VCpVGVuo1KpqtS/JPRdu3YNMTExld6H5+3tjaKiIiQnJ9f4eIiIiIiqyuiCn7m5Obp27YoDBw5IbTqdDgcOHICvr2+Z2/j6+ur1B4CYmBi9/iWh748//sD+/fthZ2dXaS3nzp2DXC6Ho6PjIx0TERERUVUY5aXeadOmYcyYMejWrRu8vLywdOlS3L17F2PHjgUAjB49Gs2aNcOCBQsAAFOmTEGvXr3wxRdfYMCAAfj+++/x3//+F9988w3wIPQNHz4cZ8+eRVRUFIqLi6X7/2xtbWFubo64uDicOnUKffr0QZMmTRAXF4epU6di1KhRaNq0qQFng4iIiIyFUQa/V199Fbdv38acOXOg0Wjg6emJ6Oho6QGO69evQy7//5OhPXr0wJYtW/DBBx/gvffeQ5s2bbBr1y507NgRAJCSkoIff/wRAODp6am3r4MHD6J3795QKBT4/vvv8dFHHyE/Px/u7u6YOnWq3n1/RERERHXJKN/j11BV9R09REREZFz4Hj8iIiIi0sPgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIiIiMhIGCX4XL140xG6JiIiIjJpBgl/nzp3h7e2NtWvXIicnxxAlEBERERkdgwS/w4cPo0OHDpg+fTrUajXGjBmDo0ePGqIUIiIiIqNhkOD3/PPPY/369UhLS8Py5cuRnJyMXr164amnnsLChQulz50RERERUe0x6MMdlpaWGDt2LA4fPozff/8dL7/8MlauXAlXV1cMGjTIkKURERERPXHq1Sfb7t69i82bN2P27NnIyspCcXGxoUuqV/jJNiIiIipLVTOC6WOtqhxHjhzB+vXrsXPnTsjlcrzyyisIDQ01dFlERERETxSDBb/U1FRERkYiMjISV65cQY8ePfDVV1/hlVdegaWlpaHKIiIiInpiGST49evXD/v374e9vT1Gjx6NcePGoW3btoYohYiIiMhoGCT4mZmZYceOHRg4cCBMTEzK7JOSkoJmzZo99tqIiIiInlQGear3xx9/xODBg8sMfRqNBpMnT0abNm0MURoRERHRE8sgwS8rKwsjR46Evb09nJ2d8dVXX0Gn02HOnDlo2bIlzpw5g4iICEOURkRERPTEMsil3pkzZ+LEiRMICQnBvn37MHXqVERHR0MulyM2NhY+Pj6GKIuIiIjoiWaQM3579+5FREQEFi9ejN27d0MIAU9PT0RFRTH0EREREdURgwS/1NRUtGvXDgDg5uYGCwsLjBo1yhClEBERERkNgwQ/IQRMTf//KrOJiQkaNWpkiFKIiIiIjIZB7vETQsDPz08Kf/fv38dLL70Ec3NzvX5nz541RHlERERETySDBL/w8HC934MHDzZEGURERERGRSaEEIYuojLHjx9Ht27doFAoDF2KQVX1A8xERERkXKqaEQxyj1919evXDykpKYYug4iIiKhBaxDBry5OSq5cuVJ6otjb2xunT5+usP/27dvh4eEBCwsLdOrUCXv27ClV45w5c6BWq9GoUSP4+/vjjz/+0OuTmZmJ4OBgKJVK2NjYIDQ0FLm5ubV+bERERERlaRDBr7Zt3boV06ZNQ3h4OM6ePYsuXbogICAAt27dKrP/iRMnMHLkSISGhuKXX37BkCFDMGTIEFy8eFHqs2jRInz11VdYvXo1Tp06BUtLSwQEBCAvL0/qExwcjISEBMTExCAqKgpHjhzBhAkTHssxExERETWIe/yaNGmCX3/9FS1btqyV8by9vdG9e3esWLECAKDT6eDi4oLJkydj1qxZpfq/+uqruHv3LqKioqQ2Hx8feHp6YvXq1RBCwNnZGdOnT8e7774LANBqtXByckJkZCRGjBiBS5cuoX379jhz5gy6desGAIiOjkb//v1x8+ZNODs7V1o37/EjIiKisjxR9/jVpoKCAsTHx8Pf319qk8vl8Pf3R1xcXJnbxMXF6fUHgICAAKl/UlISNBqNXh9ra2t4e3tLfeLi4mBjYyOFPgDw9/eHXC7HqVOnytxvfn4+srOz9RYiIiKimmoQwU8mk9XaWBkZGSguLoaTk5Neu5OTEzQaTZnbaDSaCvuX/G9lfRwdHfXWm5qawtbWttz9LliwANbW1tLi4uJS7eMlIiIiKtEggl8DuBpdJ2bPng2tVistN27cMHRJRERE1IAZ5AXO5Tl8+DDu3r0LX19fNG3aVGrPycmptX3Y29vDxMQE6enpeu3p6elQqVRlbqNSqSrsX/K/6enpUKvVen08PT2lPn9/eKSoqAiZmZnl7lehUBj9uwuJiIio9hjkjN/ChQvx4YcfSr+FEAgMDESfPn0wcOBAtGvXDgkJCXWyb3Nzc3Tt2hUHDhyQ2nQ6HQ4cOABfX98yt/H19dXrDwAxMTFSf3d3d6hUKr0+2dnZOHXqlNTH19cXWVlZiI+Pl/rExsZCp9PB29u71o+TiIiI6O8MEvy2bt2Kjh07Sr937NiBI0eO4OjRo8jIyEC3bt0wd+7cOtv/tGnTsHbtWmzYsAGXLl3Cm2++ibt372Ls2LEAgNGjR2P27NlS/ylTpiA6OhpffPEFLl++jI8++gj//e9/MWnSJODBPYjvvPMOPv30U/z444+4cOECRo8eDWdnZwwZMgQA0K5dOwQGBmL8+PE4ffo0jh8/jkmTJmHEiBFVeqKXiIiI6FEZ5FJvUlISOnfuLP3es2cPhg8fjmeffRYA8MEHH+Dll1+us/2/+uqruH37NubMmQONRgNPT09ER0dLD2dcv34dcvn/Z+IePXpgy5Yt+OCDD/Dee++hTZs22LVrl154/ec//4m7d+9iwoQJyMrKwnPPPYfo6GhYWFhIfTZv3oxJkybBz88PcrkcQUFB+Oqrr+rsOImIiIgeZpD3+P39vXweHh545513MHHiROBB8Grbti3u37//uEur1/gePyIiIipLvX6PX6tWrXDkyBHgQcj7/fff0bNnT2n9zZs3YWdnZ4jSiIiIiJ5YBrnUGxYWhkmTJuHo0aM4efIkfH190b59e2l9bGwsnn76aUOURkRERPTEMkjwGz9+PExMTLB792707NkT4eHheutTU1OlBy2IiIiIqHY0iG/10v/wHj8iIiIqS72+x68yZ8+excCBAw1dBhEREdETxWDBb9++fXj33Xfx3nvv4c8//wQAXL58GUOGDEH37t2h0+kMVRoRERHRE8kg9/h9++23GD9+PGxtbXHnzh2sW7cOX375JSZPnoxXX30VFy9eRLt27QxRGhEREdETyyBn/JYtW4aFCxciIyMD27ZtQ0ZGBr7++mtcuHABq1evZugjIiIiqgMGebjD0tISCQkJcHNzgxACCoUCBw8elL7cQWXjwx1ERERUlnr9cMf9+/fRuHFj4MF3bhUKBdRqtSFKISIiIjIaBrnHDwDWrVsHKysrAEBRUREiIyNhb2+v1+ftt982UHVERERETx6DXOp1c3ODTCarsI9MJpOe9qX/4aVeIiIiKktVM4JBzvglJycbYrdERERERs1gl3p1Oh0iIyPxww8/IDk5GTKZDC1btkRQUBBef/31Ss8IEhEREVH1GOThDiEEXnrpJfzjH/9ASkoKOnXqhA4dOiA5ORkhISEYOnSoIcoiIiIieqIZ5IxfZGQkjh49igMHDqBPnz5662JjYzFkyBBs3LgRo0ePNkR5RERERE8kg5zx++677/Dee++VCn0A0LdvX8yaNQubN282RGlERERETyyDBL/z588jMDCw3PX9+vXDr7/++lhrIiIiInrSGST4ZWZmwsnJqdz1Tk5OuHPnzmOtiYiIiOhJZ5DgV1xcDFPT8m8vNDExQVFR0WOtiYiIiOhJZ5CHO4QQCAkJgUKhKHN9fn7+Y6+JiIiI6ElnkOA3ZsyYSvvwiV4iIiKi2mWQ4BcREWGI3RIREREZNYPc40dEREREjx+DHxEREZGRYPAjIiIiMhJGF/wyMzMRHBwMpVIJGxsbhIaGIjc3t8Jt8vLyEBYWBjs7O1hZWSEoKAjp6enS+l9//RUjR46Ei4sLGjVqhHbt2mHZsmV6Yxw6dAgymazUotFo6uxYiYiIiB5mkIc7DCk4OBhpaWmIiYlBYWEhxo4diwkTJmDLli3lbjN16lT89NNP2L59O6ytrTFp0iQMGzYMx48fBwDEx8fD0dERmzZtgouLC06cOIEJEybAxMQEkyZN0hsrMTERSqVS+u3o6FiHR0tERET0/2RCCGHoIh6XS5cuoX379jhz5gy6desGAIiOjkb//v1x8+ZNODs7l9pGq9XCwcEBW7ZswfDhwwEAly9fRrt27RAXFwcfH58y9xUWFoZLly4hNjYWeHDGr0+fPrhz5w5sbGxqVH92djasra2h1Wr1wiMREREZt6pmBKO61BsXFwcbGxsp9AGAv78/5HI5Tp06VeY28fHxKCwshL+/v9Tm4eEBV1dXxMXFlbsvrVYLW1vbUu2enp5Qq9V44YUXpDOG5cnPz0d2drbeQkRERFRTRhX8NBpNqUurpqamsLW1LfdeO41GA3Nz81Jn6ZycnMrd5sSJE9i6dSsmTJggtanVaqxevRo7d+7Ezp074eLigt69e+Ps2bPl1rtgwQJYW1tLi4uLSzWPmIiIiOj/PRHBb9asWWU+OPHwcvny5cdSy8WLFzF48GCEh4fjxRdflNrbtm2LN954A127dkWPHj2wfv169OjRA0uWLCl3rNmzZ0Or1UrLjRs3HssxEBER0ZPpiXi4Y/r06QgJCamwT8uWLaFSqXDr1i299qKiImRmZkKlUpW5nUqlQkFBAbKysvTO+qWnp5fa5rfffoOfnx8mTJiADz74oNK6vby8cOzYsXLXKxSKcr9nTERERFRdT0Twc3BwgIODQ6X9fH19kZWVhfj4eHTt2hUAEBsbC51OB29v7zK36dq1K8zMzHDgwAEEBQUBD57MvX79Onx9faV+CQkJ6Nu3L8aMGYN58+ZVqe5z585BrVZX8SiJiIiIHs0TEfyqql27dggMDMT48eOxevVqFBYWYtKkSRgxYoT0RG9KSgr8/PywceNGeHl5wdraGqGhoZg2bRpsbW2hVCoxefJk+Pr6Sk/0Xrx4EX379kVAQACmTZsm3ftnYmIiBdKlS5fC3d0dHTp0QF5eHtatW4fY2Fj8/PPPBpwRIiIiMiZGFfwAYPPmzZg0aRL8/Pwgl8sRFBSEr776SlpfWFiIxMRE3Lt3T2pbsmSJ1Dc/Px8BAQH4+uuvpfU7duzA7du3sWnTJmzatElqb9GiBZKTkwEABQUFmD59OlJSUtC4cWN07twZ+/fvR58+fR7bsRMREZFxM6r3+DV0fI8fERERlYXv8SMiIiIiPQx+REREREaCwY+IiIjISDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIiIiMhIMfkRERERGgsGPiIiIyEgw+BEREREZCQY/IiIiIiPB4EdERERkJBj8iIiIiIwEgx8RERGRkWDwIyIiIjISDH5ERERERoLBj4iIiMhIMPgRERERGQkGPyIiIiIjweBHREREZCQY/IiIiIiMBIMfERERkZFg8CMiIiIyEgx+REREREbC6IJfZmYmgoODoVQqYWNjg9DQUOTm5la4TV5eHsLCwmBnZwcrKysEBQUhPT1dr49MJiu1fP/993p9Dh06hGeeeQYKhQKtW7dGZGRknRwjERERUVmMLvgFBwcjISEBMTExiIqKwpEjRzBhwoQKt5k6dSp2796N7du34/Dhw0hNTcWwYcNK9YuIiEBaWpq0DBkyRFqXlJSEAQMGoE+fPjh37hzeeecd/OMf/8C+ffvq5DiJiIiI/k4mhBCGLuJxuXTpEtq3b48zZ86gW7duAIDo6Gj0798fN2/ehLOzc6lttFotHBwcsGXLFgwfPhwAcPnyZbRr1w5xcXHw8fEBHpzx+/e//60X9h42c+ZM/PTTT7h48aLUNmLECGRlZSE6OrrMbfLz85Gfny/9zs7OhouLC7RaLZRK5SPOBhERET0psrOzYW1tXWlGMKozfnFxcbCxsZFCHwD4+/tDLpfj1KlTZW4THx+PwsJC+Pv7S20eHh5wdXVFXFycXt+wsDDY29vDy8sL69evx8OZOi4uTm8MAAgICCg1xsMWLFgAa2traXFxcanRcRMRERHB2IKfRqOBo6OjXpupqSlsbW2h0WjK3cbc3Bw2NjZ67U5OTnrbfPzxx9i2bRtiYmIQFBSEt956C8uXL9cbx8nJqdQY2dnZuH//fpn7nj17NrRarbTcuHGjRsdNREREBACmhi6gNsyaNQsLFy6ssM+lS5fqtIYPP/xQ+uunn34ad+/exeeff4633367xmMqFAooFIpaqpCIiIiM3RMR/KZPn46QkJAK+7Rs2RIqlQq3bt3Say8qKkJmZiZUKlWZ26lUKhQUFCArK0vvrF96enq52wCAt7c3PvnkE+Tn50OhUEClUpV6Ejg9PR1KpRKNGjWq4pESERER1dwTEfwcHBzg4OBQaT9fX19kZWUhPj4eXbt2BQDExsZCp9PB29u7zG26du0KMzMzHDhwAEFBQQCAxMREXL9+Hb6+vuXu69y5c2jatKl0xs7X1xd79uzR6xMTE1PhGERERES16YkIflXVrl07BAYGYvz48Vi9ejUKCwsxadIkjBgxQnqiNyUlBX5+fti4cSO8vLxgbW2N0NBQTJs2Dba2tlAqlZg8eTJ8fX2lJ3p3796N9PR0+Pj4wMLCAjExMZg/fz7effddad8TJ07EihUr8M9//hPjxo1DbGwstm3bhp9++slg80FERETGxaiCHwBs3rwZkyZNgp+fH+RyOYKCgvDVV19J6wsLC5GYmIh79+5JbUuWLJH65ufnIyAgAF9//bW03szMDCtXrsTUqVMhhEDr1q3x5ZdfYvz48VIfd3d3/PTTT5g6dSqWLVuG5s2bY926dQgICHiMR09ERETGzKje49fQabVa2NjY4MaNG3yPHxEREUlK3vWblZUFa2vrcvsZ3Rm/hiwnJwcA+D4/IiIiKlNOTk6FwY9n/BoQnU6H1NRUNGnSBDKZzNDl1KqS/1Lh2czq49w9Gs5fzXHuao5zV3Ocu7IJIZCTkwNnZ2fI5eW/ppln/BoQuVyO5s2bG7qMOqVUKvkPcg1x7h4N56/mOHc1x7mrOc5daRWd6SthVF/uICIiIjJmDH5ERERERoLBj+oFhUKB8PBwfqKuBjh3j4bzV3Ocu5rj3NUc5+7R8OEOIiIiIiPBM35ERERERoLBj4iIiMhIMPgRERERGQkGPyIiIiIjweBHREREZCQY/IiI6pGQkBDIZDJpsbOzQ2BgIM6fPy/1kclk2LVrl97vksXS0hJt2rRBSEgI4uPjDXQURFRfMfgREdUzgYGBSEtLQ1paGg4cOABTU1MMHDiwwm0iIiKQlpaGhIQErFy5Erm5ufD29sbGjRsfW91EVP8x+BER1TMKhQIqlQoqlQqenp6YNWsWbty4gdu3b5e7jY2NDVQqFdzc3PDiiy9ix44dCA4OxqRJk3Dnzp3HWj8R1V8MfkRE9Vhubi42bdqE1q1bw87OrlrbTp06FTk5OYiJiamz+oioYTE1dAFERKQvKioKVlZWAIC7d+9CrVYjKioKcnn1/lvdw8MDAJCcnFwndRJRw8MzfkRE9UyfPn1w7tw5nDt3DqdPn0ZAQAD69euHa9euVWucki9yymSyOqqUiBoaBj8ionrG0tISrVu3RuvWrdG9e3esW7cOd+/exdq1a6s1zqVLlwAA7u7udVQpETU0DH5ERPWcTCaDXC7H/fv3q7Xd0qVLoVQq4e/vX2e1EVHDwnv8iIjqmfz8fGg0GgDAnTt3sGLFCuTm5uKll14qd5usrCxoNBrk5+fj999/x5o1a7Br1y5s3LgRNjY2j7F6IqrPGPyIiOqZ6OhoqNVqAECTJk3g4eGB7du3o3fv3tDpdAAAU1P9//seO3YsAMDCwgLNmjXDc889h9OnT+OZZ54xwBEQUX0lEyV3/xIRUb2n0WigVqtx5swZdOvWzdDlEFEDwzN+REQNgBAC165dw+LFi+Hk5ISOHTsauiQiaoAY/IiIGgCtVou2bduiXbt2+P7772FhYWHokoioAeKlXiIiIiIjwde5EBERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIiIiMhIMfkRERERGgsGPiIiIyEgw+BEREREZCQY/IiIiIiPB4EdERERkJBj8iIiIiIwEgx8RERGRkWDwIyIiIjISDH5ERERERoLBj4iMRmRkJGQyWbnLyZMnAQC5ubkIDw9Hx44dYWlpCTs7O3h6emLKlClITU2Vxvvoo4/0tm/cuDHat2+PDz74ANnZ2eXu19TUFM2aNUNISAhSUlJK1dm7d+9ya/Tw8NDre+HCBQwfPhwtWrSAhYUFmjVrhhdeeAHLly/X61dQUIBly5bh6aefhlKphI2NDTp06IAJEybg8uXLpWpISEjAqFGj0KxZMygUCjg7OyM4OBgJCQml+pbMV2BgIGxtbSGTyRAZGVnun8OKFSvQrl07KBQKNGvWDNOmTcPdu3cr/fMjokdnaugCiIget48//hju7u6l2lu3bo3CwkL07NkTly9fxpgxYzB58mTk5uYiISEBW7ZswdChQ+Hs7Ky33apVq2BlZYXc3Fz8/PPPmDdvHmJjY3H8+HHIZLJS+83Ly8PJkycRGRmJY8eO4eLFi7CwsNAbs3nz5liwYEGpGq2traW/PnHiBPr06QNXV1eMHz8eKpUKN27cwMmTJ7Fs2TJMnjxZ6hsUFIS9e/di5MiRGD9+PAoLC3H58mVERUWhR48eeoHyhx9+wMiRI2Fra4vQ0FC4u7sjOTkZ3377LXbs2IHvv/8eQ4cOlfpnZGTg448/hqurK7p06YJDhw6VO/czZ87EokWLMHz4cEyZMgW//fYbli9fjoSEBOzbt6+SPzkiemSCiMhIRERECADizJkz5fbZtm2bACA2b95cat39+/eFVquVfoeHhwsA4vbt23r9hg0bJgCIEydOVLjfmTNnCgBi69ateu29evUSHTp0qPR4+vfvLxwcHMSdO3dKrUtPT5f++vTp0wKAmDdvXql+RUVFIiMjQ/p95coV0bhxY+Hh4SFu3bql1/f27dvCw8NDWFpaiqtXr0rteXl5Ii0tTQghxJkzZwQAERERUWpfqampwtTUVLz++ut67cuXLxcAxI8//ljpMRPRo+GlXiKih1y9ehUA8Oyzz5ZaZ2FhAaVSWekYffv2BQAkJSVV2O/555/X22dNau3QoQNsbGxKrXN0dNTrh3KOycTEBHZ2dtLvzz//HPfu3cM333wDBwcHvb729vZYs2YN7t69i0WLFkntCoUCKpWq0nrj4uJQVFSEESNG6LWX/P7+++8rHYOIHg2DHxEZHa1Wi4yMDL3lr7/+AgC0aNECALBx40YIIWo0fknQejhQlSU5ORkA0LRp01LriouLS9WYkZGhdy9cixYtEB8fj4sXL1a4n5Jj2rx5M4qKiirsu3v3bri5uUmh9O969uwJNzc3/PTTTxWOU5b8/HwAQKNGjfTaGzduDACIj4+v9phEVD0MfkRkdPz9/eHg4KC3NGvWDAAwZMgQtG3bFnPmzIG7uzvGjh2L9evX49atW+WOl5mZiYyMDCQnJ+Obb77B119/DScnp1LhqSRw3rx5Ezt37sTcuXOhUCgwcODAUmNevny5VI0ODg6YPn261Ofdd9/FvXv34OnpiR49emDmzJn4+eefUVhYqDeWj48PevXqhbVr16J58+Z47bXX8PXXX+P69eul6ktNTUWXLl0qnL/OnTvj5s2byMnJqWSm9bVt2xYAcPz4cb32o0ePAkCZD7oQUe3iwx1EZHRWrlyJp556Sq/NxMQEeHA26tSpU5g3bx62bduGyMhIREZGQi6X46233sLixYuhUCj0ti0JNCU6dOiADRs2SGeySvj7++v9dnNzw6ZNm9C8efNSNbq5uWHt2rWl2h/u+8ILLyAuLg4LFizAvn37EBcXh0WLFsHBwQHr1q3DoEGDAAAymQz79u3D4sWLsWnTJnz33Xf47rvvEBYWhldeeQVr1qyBjY2NFOSaNGlS4fyVrM/Ozq6078OeeeYZeHt7Y+HChWjWrBn69OmDS5cu4c0334SZmRnu379f5bGIqIYMfZMh1Z2kpCQxbtw44ebmJiwsLETLli3FnDlzRH5+foXbrVmzRvTq1Us0adJEACjzxvGXXnpJuLi4CIVCIVQqlRg1apRISUnR6xMdHS28vb2FlZWVsLe3F8OGDRNJSUnVOoaq1EJUVVV5uOPvkpOTxbfffivatWsnAIj3339fWlfycMfOnTtFTEyMOHTokLhy5Uq5+125cqWIiYkRO3bsEP379xdWVlbi0KFDpfpX9eGOh+Xn54vTp0+L2bNnCwsLC2FmZiYSEhLK7Juamiq+++474ePjIwCI4OBgIYQQWVlZAoAYPHhwhfsaNGiQAKD3oEuJih7uEEKImzdvimeffVYAEACEiYmJmDFjhvDy8hLW1tbVOmYiqj5e6n2CXb58GTqdDmvWrEFCQgKWLFmC1atX47333qtwu3v37iEwMLDCfn369MG2bduQmJiInTt34urVqxg+fLi0PikpCYMHD0bfvn1x7tw57Nu3DxkZGRg2bFi1jqEqtRDVpRYtWmDcuHE4fvw4bGxssHnz5lJ9evbsCX9/f/Tq1QutWrUqdywvLy/4+/sjKCgIP/74Izp27IjXXnsNubm5j1ynubk5unfvjvnz52PVqlUoLCzE9u3by+yrVqsxYsQIHDlyBG3atMG2bdtQVFQEa2trqNVqnD9/vsJ9nT9/Hs2aNavSgy5/16xZMxw7dgy///47jhw5gps3b2LRokW4ceNGqbOwRFQHDJ086fFatGiRcHd3r1LfgwcPVvks23/+8x8hk8lEQUGBEEKI7du3C1NTU1FcXCz1+fHHH/X6CCHErl27xNNPPy0UCoVwd3cXH330kSgsLHykWojKU5Mzfg/r2rWrUCgU0u/yXudS1f2W/H29YMECvfaanPF72IULFwQA8cYbb1TaNygoSACQXscyfvx4AUAcPXq0zP5HjhypcOzKzviVJSEhQQAQs2fPrvI2RFQzPONnZLRaLWxtbWt1zMzMTGzevBk9evSAmZkZAKBr166Qy+WIiIhAcXExtFot/vWvf8Hf31/qc/ToUYwePVp6ieuaNWsQGRmJefPm1Wp9RNXx66+/IiMjo1T7tWvX8Ntvv5W6n+9R9O7dG15eXli6dCny8vKqvf3BgwfLfPJ4z549wEP3Hv7xxx+lHuQAgKysLMTFxaFp06bSq1tmzJiBRo0a4Y033pCedC6RmZmJiRMnonHjxpgxY0a16y2LTqfDP//5TzRu3BgTJ06slTGJqHx8uMOIXLlyBcuXL8fixYtrZbyZM2dixYoVuHfvHnx8fBAVFSWtc3d3x88//4xXXnkFb7zxBoqLi+Hr6yv9CwkA5s6di1mzZmHMmDEAgJYtW+KTTz7BP//5T4SHh9dKjURl2bt3b5mfKevRowdiYmIQHh6OQYMGwcfHB1ZWVvjzzz+xfv165Ofn46OPPqrVWmbMmIGXX34ZkZGResFHq9Vi06ZNZW4zatQoAMDkyZNx7949DB06FB4eHigoKMCJEyewdetWuLm5YezYscCDMPvaa6+hX79+eP7552Fra4uUlBRs2LABqampWLp0qfRwS5s2bbBhwwYEBwejU6dOpb7ckZGRge+++67UJe0VK1YgKytL+qTd7t27cfPmTanOki+OTJkyBXl5efD09ERhYSG2bNmC06dPY8OGDXB1da3VuSWiMhj6lCNVX8nb/itaLl26pLfNzZs3RatWrURoaGiV91PZ5dXbt2+LxMRE8fPPP4tnn31W9O/fX+h0OiGEEGlpaaJNmzZixowZ4uzZs+Lw4cOiV69ews/PT+pjb28vLCwshKWlpbRYWFgIAOLu3bvVqoWoKkouuZa3REREiD///FPMmTNH+Pj4CEdHR2FqaiocHBzEgAEDRGxsrN54j3qpVwghiouLRatWrUSrVq1EUVGREA8u9VZUZ4m9e/eKcePGCQ8PD2FlZSXMzc1F69atxeTJk/W+3JGeni4+++wz0atXL6FWq4Wpqalo2rSp6Nu3r9ixY0eZNZ8/f16MHDlSqNVqYWZmJlQqlRg5cqS4cOFCmf1btGhRbr0PP9QVEREhunTpIiwtLUWTJk2En59fqXklorojEzV9QykZzO3bt0tdgvm7li1bwtzcHACQmpqK3r17w8fHR3otRVUcOnQIffr0wZ07d8r8MsDDbt68CRcXF5w4cQK+vr748MMPER0djTNnzpTqExcXBx8fHzRq1Ahz584t84GPli1b6tVZnVqIiIiobLzU2wCVvMi1KlJSUtCnTx907doVERERVQ591aXT6YCH3sx/7969UvsquZRU0veZZ55BYmIiWrduXSc1ERERkT4GvydYSkoKevfujRYtWmDx4sW4ffu2tK7ku5opKSnw8/PDxo0b4eXlBQDQaDTQaDS4cuUKAODChQto0qQJXF1dYWtri1OnTuHMmTN47rnn0LRpU1y9ehUffvghWrVqBV9fXwDAgAEDsGTJEnz88ccYOXIkcnJy8N5776FFixZ4+umnAQBz5szBwIED4erqiuHDh0Mul+PXX3/FxYsX8emnn1apFiIiIqoGQ19rprpT0f1MJZKSkgQAcfDgQamt5L6lsu5/Eg/u/enTp4+wtbUVCoVCuLm5iYkTJ4qbN2/q7f+7774TTz/9tLC0tBQODg5i0KBBpe49jI6OFj169BCNGjUSSqVSeHl5iW+++abKtRAREVHV8R4/IiIiIiPB9/gRERERGQne49eA6HQ6pKamokmTJpDJZIYuh4iIiOoJIQRycnLg7Oxc4YOcDH4NSGpqKlxcXAxdBhEREdVTN27cQPPmzctdz+DXgDRp0gR48Idak4+jExER0ZMpOzsbLi4uUlYoD4NfA1JyeVepVDL4ERERUSmV3QrG4EdERERUh4p1AqeTMnErJw+OTSzg5W4LE7lh7tVn8CMiIiKqI9EX0zB3929I0+ZJbWprC4S/1B6BHdWPvR6+zoWIiIioDkRfTMObm87qhT4A0Gjz8Oams4i+mPbYa2LwIyIiIqplxTqBubt/Q1lfyShpm7v7NxTrHu93NBj8iIiIiGrZ6aTMUmf6HiYApGnzcDop87HWxeBHREREVMtu5ZQf+mrSr7Yw+BERERHVMscmFrXar7Yw+BERERHVMi93W6itLVDeS1tkD57u9XK3fax1MfgRERER1TITuQzhL7UHHoS8h5X8Dn+p/WN/nx+DHxEREVEdCOyoxqpRz0BlrX85V2VtgVWjnjHIe/z4AmciIiKiOhLYUY0X2qv45Q4iIiIiY2Ail8G3lZ2hywB4qZeIiIjIeDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIiIiMhIMfkRERERGgsGPiIiIyEgw+BEREREZiQYT/DIzMxEcHAylUgkbGxuEhoYiNze3wm3y8vIQFhYGOzs7WFlZISgoCOnp6Xp9rl+/jgEDBqBx48ZwdHTEjBkzUFRUVOZ4x48fh6mpKTw9PfXa3dzcIJPJSi1hYWFSn969e5daP3HixEeaEyIiIqLqaDDBLzg4GAkJCYiJiUFUVBSOHDmCCRMmVLjN1KlTsXv3bmzfvh2HDx9Gamoqhg0bJq0vLi7GgAEDUFBQgBMnTmDDhg2IjIzEnDlzSo2VlZWF0aNHw8/Pr9S6M2fOIC0tTVpiYmIAAC+//LJev/Hjx+v1W7Ro0SPMCBEREVH1yIQQwtBFVObSpUto3749zpw5g27dugEAoqOj0b9/f9y8eRPOzs6lttFqtXBwcMCWLVswfPhwAMDly5fRrl07xMXFwcfHB3v37sXAgQORmpoKJycnAMDq1asxc+ZM3L59G+bm5tJ4I0aMQJs2bWBiYoJdu3bh3Llz5db7zjvvICoqCn/88QdkMhnw4Iyfp6cnli5dWuN5yM7OhrW1NbRaLZRKZY3HISIioidLVTNCgzjjFxcXBxsbGyn0AYC/vz/kcjlOnTpV5jbx8fEoLCyEv7+/1Obh4QFXV1fExcVJ43bq1EkKfQAQEBCA7OxsJCQkSG0RERH4888/ER4eXmmtBQUF2LRpE8aNGyeFvhKbN2+Gvb09OnbsiNmzZ+PevXsVjpWfn4/s7Gy9hYiIiKimTA1dQFVoNBo4OjrqtZmamsLW1hYajabcbczNzWFjY6PX7uTkJG2j0Wj0Ql/J+pJ1APDHH39g1qxZOHr0KExNK5+uXbt2ISsrCyEhIXrtr732Glq0aAFnZ2ecP38eM2fORGJiIn744Ydyx1qwYAHmzp1b6T6JiIiIqsKgwW/WrFlYuHBhhX0uXbr02Or5u+LiYrz22muYO3cunnrqqSpt8+2336Jfv36lLj8/fD9ip06doFar4efnh6tXr6JVq1ZljjV79mxMmzZN+p2dnQ0XF5caHw8REREZN4MGv+nTp5c6M/Z3LVu2hEqlwq1bt/Tai4qKkJmZCZVKVeZ2KpUKBQUFyMrK0jvrl56eLm2jUqlw+vRpve1KnvpVqVTIycnBf//7X/zyyy+YNGkSAECn00EIAVNTU/z888/o27evtO21a9ewf//+Cs/ilfD29gYAXLlypdzgp1AooFAoKh2LiIiIqCoMGvwcHBzg4OBQaT9fX19kZWUhPj4eXbt2BQDExsZCp9NJAervunbtCjMzMxw4cABBQUEAgMTERFy/fh2+vr7SuPPmzcOtW7ekS8kxMTFQKpVo3749zMzMcOHCBb1xv/76a8TGxmLHjh1wd3fXWxcREQFHR0cMGDCg0mMqeThErVZX2peIiIioNjSIe/zatWuHwMBAjB8/HqtXr0ZhYSEmTZqEESNGSJdUU1JS4Ofnh40bN8LLywvW1tYIDQ3FtGnTYGtrC6VSicmTJ8PX1xc+Pj4AgBdffBHt27fH66+/jkWLFkGj0eCDDz5AWFiYdKatY8eOerU4OjrCwsKiVLtOp0NERATGjBlT6l7Aq1evYsuWLejfvz/s7Oxw/vx5TJ06FT179kTnzp3rePaIiIiI/qdBBD88eCJ20qRJ8PPzg1wuR1BQEL766itpfWFhIRITE/WelF2yZInUNz8/HwEBAfj666+l9SYmJoiKisKbb74JX19fWFpaYsyYMfj444+rXd/+/ftx/fp1jBs3rtQ6c3Nz7N+/H0uXLsXdu3fh4uKCoKAgfPDBBzWaCyIiIqKaaBDv8aP/4Xv8iIiIqCxP1Hv8iIiIiOjRMfgRERERGQkGPyIiIiIjweBHREREZCQY/IiIiIiMBIMfERERkZFg8CMiIiIyEgx+REREREaCwY+IiIjISDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRYPAjIiIiMhIMfkRERERGgsGPiIiIyEgw+BEREREZCQY/IiIiIiPB4EdERERkJBj8iIiIiIxEtYLfxYsX664SIiIiIqpT1Qp+nTt3hre3N9auXYucnJy6q4qIiIiIal21gt/hw4fRoUMHTJ8+HWq1GmPGjMHRo0frrjoiIiIiqjXVCn7PP/881q9fj7S0NCxfvhzJycno1asXnnrqKSxcuBAajabOCs3MzERwcDCUSiVsbGwQGhqK3NzcCrfJy8tDWFgY7OzsYGVlhaCgIKSnp+v1uX79OgYMGIDGjRvD0dERM2bMQFFRkbT+0KFDkMlkpZa/H+vKlSvh5uYGCwsLeHt74/Tp09WuhYiIiKgu1ejhDktLS4wdOxaHDx/G77//jpdffhkrV66Eq6srBg0aVPtVAggODkZCQgJiYmIQFRWFI0eOYMKECRVuM3XqVOzevRvbt2/H4cOHkZqaimHDhknri4uLMWDAABQUFODEiRPYsGEDIiMjMWfOnFJjJSYmIi0tTVocHR2ldVu3bsW0adMQHh6Os2fPokuXLggICMCtW7eqXAsRERFRnRO1IDc3V6xZs0bY2toKuVxeG0Pq+e233wQAcebMGalt7969QiaTiZSUlDK3ycrKEmZmZmL79u1S26VLlwQAERcXJ4QQYs+ePUIulwuNRiP1WbVqlVAqlSI/P18IIcTBgwcFAHHnzp1y6/Py8hJhYWHS7+LiYuHs7CwWLFhQ5VrKkpeXJ7RarbTcuHFDABBarbbSOSMiIiLjodVqq5QRHul1LkeOHEFISAhUKhVmzJiBYcOG4fjx47WXSh+Ii4uDjY0NunXrJrX5+/tDLpfj1KlTZW4THx+PwsJC+Pv7S20eHh5wdXVFXFycNG6nTp3g5OQk9QkICEB2djYSEhL0xvP09IRarcYLL7ygd4wFBQWIj4/X249cLoe/v7+0n6rUUpYFCxbA2tpaWlxcXKo8Z0RERER/V+3gl5qaivnz5+Opp55C7969ceXKFXz11VdITU3F2rVr4ePjU+tFajQavUurAGBqagpbW9ty7yvUaDQwNzeHjY2NXruTk5O0jUaj0Qt9JetL1gGAWq3G6tWrsXPnTuzcuRMuLi7o3bs3zp49CwDIyMhAcXFxmeM8vJ/KainL7NmzodVqpeXGjRuVzBQRERFR+Uyr07lfv37Yv38/7O3tMXr0aIwbNw5t27at8c5nzZqFhQsXVtjn0qVLNR6/NrRt21bvGHv06IGrV69iyZIl+Ne//lWn+1YoFFAoFHW6DyIiIjIe1Qp+ZmZm2LFjBwYOHAgTE5My+6SkpKBZs2ZVGm/69OkICQmpsE/Lli2hUqn0HpQAgKKiImRmZkKlUpW5nUqlQkFBAbKysvTOtKWnp0vbqFSqUk/fljxpW964AODl5YVjx44BAOzt7WFiYlLqCd2/76eyWoiIiIjqWrUu9f74448YPHhwmaFPo9Fg8uTJaNOmTZXHc3BwgIeHR4WLubk5fH19kZWVhfj4eGnb2NhY6HQ6eHt7lzl2165dYWZmhgMHDkhtiYmJuH79Onx9fQEAvr6+uHDhgl6ojImJgVKpRPv27cut+9y5c1Cr1QAAc3NzdO3aVW8/Op0OBw4ckPZTlVqIiIiI6lx1nhi5c+eOGDFihLCzsxNqtVosW7ZMFBcXiw8//FA0atRIeHt7i++///5RH0wpU2BgoHj66afFqVOnxLFjx0SbNm3EyJEjpfU3b94Ubdu2FadOnZLaJk6cKFxdXUVsbKz473//K3x9fYWvr6+0vqioSHTs2FG8+OKL4ty5cyI6Olo4ODiI2bNnS32WLFkidu3aJf744w9x4cIFMWXKFCGXy8X+/fulPt9//71QKBQiMjJS/Pbbb2LChAnCxsZG72nhymqpiqo+sUNERETGpaoZoVrBb8KECcLV1VVMnz5ddOzYUcjlctGvXz8xYMCACl9LUhv++usvMXLkSGFlZSWUSqUYO3asyMnJkdYnJSUJAOLgwYNS2/3798Vbb70lmjZtKho3biyGDh0q0tLS9MZNTk4W/fr1E40aNRL29vZi+vTporCwUFq/cOFC0apVK2FhYSFsbW1F7969RWxsbKn6li9fLlxdXYW5ubnw8vISJ0+e1FtflVoqw+BHREREZalqRpAJIURVzw66uroiMjISffv2RXJyMlq2bIlZs2Zh/vz5dXtakgAA2dnZsLa2hlarhVKpNHQ5REREVE9UNSNU6x6/1NRUtGvXDgCkz5ONGjXq0aslIiIiojpXreAnhICp6f8/CGxiYoJGjRrVRV1EREREVMuq9ToXIQT8/Pyk8Hf//n289NJLMDc31+tX8nJjIiIiIqo/qhX8wsPD9X4PHjy4tushIiIiojpSrYc7quv48ePo1q0bvz5RS/hwBxEREZWlTh7uqK5+/fohJSWlLndBRERERFVUp8GvDk8mEhEREVE11WnwIyIiIqL6g8GPiIiIyEgw+BEREREZiToNfjKZrC6HJyIiIqJq4MMdREREREaiWi9wLs/hw4dx9+5d+Pr6omnTplJ7Tk5ObQxPRERERLWgWsFv4cKFyM3NxSeffAI8OKPXr18//PzzzwAAR0dHHDhwAB06dKibaomIiIioxqp1qXfr1q3o2LGj9HvHjh04cuQIjh49ioyMDHTr1g1z586tizqJiIiI6BFVK/glJSWhc+fO0u89e/Zg+PDhePbZZ2Fra4sPPvgAcXFxdVEnERERET2iagW/oqIive/uxsXFoUePHtJvZ2dnZGRk1G6FRERERFQrqhX8WrVqhSNHjgAArl+/jt9//x09e/aU1t+8eRN2dna1XyURERERPbJqPdwRFhaGSZMm4ejRozh58iR8fX3Rvn17aX1sbCyefvrpuqiTiIiIiB5RtYLf+PHjYWJigt27d6Nnz54IDw/XW5+amoqxY8fWdo1EREREVAtkgm9ZbjCys7NhbW0NrVYLpVJp6HKIiIionqhqRqjVL3ecPXsWAwcOrM0hiYiIiKiWVDv47du3D++++y7ee+89/PnnnwCAy5cvY8iQIejevTt0Ol1d1ElEREREj6hawe/bb79Fv379EBkZiYULF8LHxwebNm2Cr68vVCoVLl68iD179tRJoZmZmQgODoZSqYSNjQ1CQ0ORm5tb4TZ5eXkICwuDnZ0drKysEBQUhPT0dL0+169fx4ABA9C4cWM4OjpixowZKCoqktYfOnQIMpms1KLRaKQ+CxYsQPfu3dGkSRM4OjpiyJAhSExM1NtP7969S40xceLEWpsfIiIiospUK/gtW7YMCxcuREZGBrZt24aMjAx8/fXXuHDhAlavXo127drVWaHBwcFISEhATEwMoqKicOTIEUyYMKHCbaZOnYrdu3dj+/btOHz4MFJTUzFs2DBpfXFxMQYMGICCggKcOHECGzZsQGRkJObMmVNqrMTERKSlpUmLo6OjtO7w4cMICwvDyZMnERMTg8LCQrz44ou4e/eu3hjjx4/XG2PRokW1MjdEREREVSKqoXHjxiIpKUkIIYROpxNmZmbi2LFj1RmiRn777TcBQJw5c0Zq27t3r5DJZCIlJaXMbbKysoSZmZnYvn271Hbp0iUBQMTFxQkhhNizZ4+Qy+VCo9FIfVatWiWUSqXIz88XQghx8OBBAUDcuXOnyvXeunVLABCHDx+W2nr16iWmTJlSzSPXp9VqBQCh1WofaRwiIiJ6slQ1I1TrjN/9+/fRuHFjAIBMJoNCoYBara6rTCqJi4uDjY0NunXrJrX5+/tDLpfj1KlTZW4THx+PwsJC+Pv7S20eHh5wdXWVPisXFxeHTp06wcnJSeoTEBCA7OxsJCQk6I3n6ekJtVqNF154AcePH6+wXq1WCwCwtbXVa9+8eTPs7e3RsWNHzJ49G/fu3atwnPz8fGRnZ+stRERERDVVrff4AcC6detgZWUFPPiEW2RkJOzt7fX6vP3227VXIQCNRqN3aRUATE1NYWtrq3ev3d+3MTc3h42NjV67k5OTtI1Go9ELfSXrS9YBgFqtxurVq9GtWzfk5+dj3bp16N27N06dOoVnnnmm1H51Oh3eeecdPPvss+jYsaPU/tprr6FFixZwdnbG+fPnMXPmTCQmJuKHH34o97gXLFiAuXPnVmGGiIiIiCpXreDn6uqKtWvXSr9VKhX+9a9/6fWRyWRVDn6zZs3CwoULK+xz6dKl6pRY69q2bYu2bdtKv3v06IGrV69iyZIlpY4dD75ucvHiRRw7dkyv/eH7ETt16gS1Wg0/Pz9cvXoVrVq1KnPfs2fPxrRp06Tf2dnZcHFxqaUjIyIiImNTreCXnJxcqzufPn06QkJCKuzTsmVLqFQq3Lp1S6+9qKgImZmZUKlUZW6nUqlQUFCArKwsvbN+6enp0jYqlQqnT5/W267kqd/yxgUALy+vUsEOACZNmiQ9eNK8efMKj8vb2xsAcOXKlXKDn0KhgEKhqHAcIiIioqqq9qVenU6HyMhI/PDDD0hOToZMJkPLli0RFBSE119/HTKZrMpjOTg4wMHBodJ+vr6+yMrKQnx8PLp27Qo8+C6wTqeTAtTfde3aFWZmZjhw4ACCgoKAB0/mXr9+Hb6+vtK48+bNw61bt6RLyTExMVAqlXrfIP67c+fO6d3bKITA5MmT8e9//xuHDh2Cu7t7pcd07tw54MGlZCIiIqLHoVqfbBNCYODAgdi7dy+6dOkCDw8PCCFw6dIlXLhwAYMGDcKuXbvqpNB+/fohPT0dq1evRmFhIcaOHYtu3bphy5YtAICUlBT4+flh48aN8PLyAgC8+eab2LNnDyIjI6FUKjF58mQAwIkTJ4AHr3Px9PSEs7MzFi1aBI1Gg9dffx3/+Mc/MH/+fADA0qVL4e7ujg4dOiAvLw/r1q3D8uXL8fPPP8PPzw8A8NZbb2HLli34z3/+o3dZ2NraGo0aNcLVq1exZcsW9O/fH3Z2djh//jymTp2K5s2b4/Dhw1WeA36yjYiIiMpS5YxQnUeF169fL5o0aSJiY2NLrTtw4IBo0qSJ2LBhQ7UfQa6Kv/76S4wcOVJYWVkJpVIpxo4dK3JycqT1SUlJAoA4ePCg1Hb//n3x1ltviaZNm4rGjRuLoUOHirS0NL1xk5OTRb9+/USjRo2Evb29mD59uigsLJTWL1y4ULRq1UpYWFgIW1tb0bt371LHD6DMJSIiQgghxPXr10XPnj2Fra2tUCgUonXr1mLGjBnVfi0LX+dCREREZalqRqjWGb8XX3wRffv2xaxZs8pcP3/+fBw+fBj79u2rUVqlivGMHxEREZWlqhmhWu/xO3/+PAIDA8td369fP/z666/Vq5SIiIiIHotqBb/MzMxS7717mJOTE+7cuVMbdRERERFRLatW8CsuLoapafkPApuYmKCoqKg26iIiIiKiWlat17kIIRASElLuu+Xy8/Nrqy4iIiIiqmXVCn5jxoyptM/o0aMfpR4iIiIiqiPVCn4RERF1VwkRERER1alq3eNHRERERA0Xgx8RERGRkWDwIyIiIjISDH5ERERERoLBj4iIiMhIMPgRERERGQkGPyIiIiIjweBHREREZCQY/IiIiIiMBIMfERERkZFg8CMiIiIyEgx+REREREaCwY+IiIjISDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjESDCX6ZmZkIDg6GUqmEjY0NQkNDkZubW+E2eXl5CAsLg52dHaysrBAUFIT09HS9PtevX8eAAQPQuHFjODo6YsaMGSgqKpLWHzp0CDKZrNSi0WikPh999FGp9R4eHtWuhYiIiKguNZjgFxwcjISEBMTExCAqKgpHjhzBhAkTKtxm6tSp2L17N7Zv347Dhw8jNTUVw4YNk9YXFxdjwIABKCgowIkTJ7BhwwZERkZizpw5pcZKTExEWlqatDg6Ouqt79Chg976Y8eOVasWIiIiojonGoDffvtNABBnzpyR2vbu3StkMplISUkpc5usrCxhZmYmtm/fLrVdunRJABBxcXFCCCH27Nkj5HK50Gg0Up9Vq1YJpVIp8vPzhRBCHDx4UAAQd+7cKbe+8PBw0aVLl3LXV6WWqtBqtQKA0Gq1Vd6GiIiInnxVzQgN4oxfXFwcbGxs0K1bN6nN398fcrkcp06dKnOb+Ph4FBYWwt/fX2rz8PCAq6sr4uLipHE7deoEJycnqU9AQACys7ORkJCgN56npyfUajVeeOEFHD9+vNT+/vjjDzg7O6Nly5YIDg7G9evXq1VLWfLz85Gdna23EBEREdVUgwh+Go2m1KVVU1NT2Nra6t1r9/dtzM3NYWNjo9fu5OQkbaPRaPRCX8n6knUAoFarsXr1auzcuRM7d+6Ei4sLevfujbNnz0rbeHt7IzIyEtHR0Vi1ahWSkpLw/PPPIycnp8q1lGXBggWwtraWFhcXlyrNFxEREVFZTA2581mzZmHhwoUV9rl06dJjq6csbdu2Rdu2baXfPXr0wNWrV7FkyRL861//AgD069dPWt+5c2d4e3ujRYsW2LZtG0JDQ2u879mzZ2PatGnS7+zsbIY/IiIiqjGDBr/p06cjJCSkwj4tW7aESqXCrVu39NqLioqQmZkJlUpV5nYqlQoFBQXIysrSO9OWnp4ubaNSqXD69Gm97UqetC1vXADw8vIq9fDGw2xsbPDUU0/hypUrVa6lLAqFAgqFotz1RERERNVh0Eu9Dg4O8PDwqHAxNzeHr68vsrKyEB8fL20bGxsLnU4Hb2/vMsfu2rUrzMzMcODAAaktMTER169fh6+vLwDA19cXFy5c0AuVMTExUCqVaN++fbl1nzt3Dmq1utz1ubm5uHr1qtSnKrUQERER1TWDnvGrqnbt2iEwMBDjx4/H6tWrUVhYiEmTJmHEiBFwdnYGAKSkpMDPzw8bN26El5cXrK2tERoaimnTpsHW1hZKpRKTJ0+Gr68vfHx8AAAvvvgi2rdvj9dffx2LFi2CRqPBBx98gLCwMOlM29KlS+Hu7o4OHTogLy8P69atQ2xsLH7++WepvnfffRcvvfQSWrRogdTUVISHh8PExAQjR44EgCrVQkRERFTXGkTwA4DNmzdj0qRJ8PPzg1wuR1BQEL766itpfWFhIRITE3Hv3j2pbcmSJVLf/Px8BAQE4Ouvv5bWm5iYICoqCm+++SZ8fX1haWmJMWPG4OOPP5b6FBQUYPr06UhJSUHjxo3RuXNn7N+/H3369JH63Lx5EyNHjsRff/0FBwcHPPfcczh58iQcHByqXAsRERFRXZMJIYShi6Cqyc7OhrW1NbRaLZRKpaHLISIionqiqhmhQbzOhYiIiIgeHYMfERERkZFg8CMiIiIyEgx+REREREaCwY+IiIjISDD4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHxEREZGRaDDf6qW6U6wTOJ2UiVs5eXBsYgEvd1uYyGWGLouIiIhqGYOfkYu+mIa5u39DmjZPalNbWyD8pfYI7Kg2aG1ERERUu3ip14hFX0zDm5vO6oU+ANBo8/DmprOIvphmsNqIiIio9jH4GalincDc3b9BlLGupG3u7t9QrCurBxERETVEDH5G6nRSZqkzfQ8TANK0eTidlPlY6yIiIqK6w+BnpG7llB/6atKPiIiI6j8GPyPl2MSiVvsRERFR/cfgZ6S83G2htrZAeS9tkT14utfL3fYxV0ZERER1hcHPSJnIZQh/qT3wIOQ9rOR3+Evt+T4/IiKiJwiDnxEL7KjGqlHPQGWtfzlXZW2BVaOe4Xv8iIiInjB8gbORC+yoxgvtVfxyBxERkRFg8COYyGXwbWVn6DKIiIiojjH4NSBC/O9lytnZ2YYuhYiIiOqRkmxQkhXKw+DXgOTk5AAAXFxcDF0KERER1UM5OTmwtrYud71MVBYNqd7Q6XRITU1FkyZNIJM9WffgZWdnw8XFBTdu3IBSqTR0OQ0K5+7RcP5qjnNXc5y7muPclU0IgZycHDg7O0MuL//ZXZ7xa0DkcjmaN29u6DLqlFKp5D/INcS5ezScv5rj3NUc567mOHelVXSmrwRf50JERERkJBj8iIiIiIwEgx/VCwqFAuHh4VAoFIYupcHh3D0azl/Nce5qjnNXc5y7R8OHO4iIiIiMBM/4ERERERkJBj8iIiIiI8HgR0RERGQkGPyIiIiIjASDHz1Wbm5ukMlkpZawsDAAwNWrVzF06FA4ODhAqVTilVdeQXp6uqHLrheKi4vx4Ycfwt3dHY0aNUKrVq3wySef6H2XUQiBOXPmQK1Wo1GjRvD398cff/xh0Lrrg6rM3Q8//IAXX3wRdnZ2kMlkOHfunEFrri8qm7vCwkLMnDkTnTp1gqWlJZydnTF69GikpqYaunSDq8rfdx999BE8PDxgaWmJpk2bwt/fH6dOnTJo3fVBVebuYRMnToRMJsPSpUsfe60NjiB6jG7duiXS0tKkJSYmRgAQBw8eFLm5uaJly5Zi6NCh4vz58+L8+fNi8ODBonv37qK4uNjQpRvcvHnzhJ2dnYiKihJJSUli+/btwsrKSixbtkzq89lnnwlra2uxa9cu8euvv4pBgwYJd3d3cf/+fYPWbmhVmbuNGzeKuXPnirVr1woA4pdffjFozfVFZXOXlZUl/P39xdatW8Xly5dFXFyc8PLyEl27djV06QZXlb/vNm/eLGJiYsTVq1fFxYsXRWhoqFAqleLWrVsGrd3QqjJ3JX744QfRpUsX4ezsLJYsWWKQehsSBj8yqClTpohWrVoJnU4n9u3bJ+RyudBqtdL6rKwsIZPJRExMjEHrrA8GDBggxo0bp9c2bNgwERwcLIQQQqfTCZVKJT7//HNpfVZWllAoFOK777577PXWJ5XN3cOSkpIY/B5Snbkrcfr0aQFAXLt27TFUWH/VZO60Wq0AIPbv3/8YKqy/qjp3N2/eFM2aNRMXL14ULVq0YPCrAl7qJYMpKCjApk2bMG7cOMhkMuTn50Mmk+m9lNPCwgJyuRzHjh0zaK31QY8ePXDgwAH8/vvvAIBff/0Vx44dQ79+/QAASUlJ0Gg08Pf3l7axtraGt7c34uLiDFZ3fVDZ3FH5ajJ3Wq0WMpkMNjY2j7HS+qe6c1dQUIBvvvkG1tbW6NKly2Outn6pytzpdDq8/vrrmDFjBjp06GDAahsWU0MXQMZr165dyMrKQkhICADAx8cHlpaWmDlzJubPnw8hBGbNmoXi4mKkpaUZulyDmzVrFrKzs+Hh4QETExMUFxdj3rx5CA4OBgBoNBoAgJOTk952Tk5O0jpjVdncUfmqO3d5eXmYOXMmRo4cCaVS+djrrU+qOndRUVEYMWIE7t27B7VajZiYGNjb2xus7vqgKnO3cOFCmJqa4u233zZorQ0Nz/iRwXz77bfo168fnJ2dAQAODg7Yvn07du/eDSsrK1hbWyMrKwvPPPMM5HL+rbpt2zZs3rwZW7ZswdmzZ7FhwwYsXrwYGzZsMHRp9R7nruaqM3eFhYV45ZVXIITAqlWrDFJvfVLVuevTpw/OnTuHEydOIDAwEK+88gpu3bplsLrrg8rmLj4+HsuWLUNkZCRkMpmhy21YDH2tmYxTcnKykMvlYteuXWWuv337trhz544QQggnJyexaNGix1xh/dO8eXOxYsUKvbZPPvlEtG3bVgghxNWrV8u8N61nz57i7bfffqy11jeVzd3DeI+fvqrOXUFBgRgyZIjo3LmzyMjIeMxV1k/V+fvuYa1btxbz58+v4+rqt8rmbsmSJUImkwkTExNpASDkcrlo0aKFgapuGHgahQwiIiICjo6OGDBgQJnr7e3tYWNjg9jYWNy6dQuDBg167DXWN/fu3St15tPExAQ6nQ4A4O7uDpVKhQMHDkjrs7OzcerUKfj6+j72euuTyuaOyleVuSs50/fHH39g//79sLOzM0Cl9U9N/77T6XTIz8+v4+rqt8rm7vXXX8f58+dx7tw5aXF2dsaMGTOwb98+A1XdMPAeP3rsdDodIiIiMGbMGJia6v8tGBERgXbt2sHBwQFxcXGYMmUKpk6dirZt2xqs3vripZdewrx58+Dq6ooOHTrgl19+wZdffolx48YBAGQyGd555x18+umnaNOmDdzd3fHhhx/C2dkZQ4YMMXT5BlXZ3AFAZmYmrl+/Lr1/LjExEQCgUqmgUqkMVruhVTZ3hYWFGD58OM6ePYuoqCgUFxdL95Ta2trC3NzcwEdgOJXN3d27dzFv3jwMGjQIarUaGRkZWLlyJVJSUvDyyy8bunyDqmzu7OzsSv0HhpmZGVQqFf99URlDn3Ik47Nv3z4BQCQmJpZaN3PmTOHk5CTMzMxEmzZtxBdffCF0Op1B6qxvsrOzxZQpU4Srq6uwsLAQLVu2FO+//77Iz8+X+uh0OvHhhx8KJycnoVAohJ+fX5nzbGyqMncRERECQKklPDzcoLUbWmVzV3JpvKzl4MGDhi7foCqbu/v374uhQ4cKZ2dnYW5uLtRqtRg0aJA4ffq0oUs3uKr8M/t3fJ1L1chEea/BJiIiIqInCu/xIyIiIjISDH5ERERERoLBj4iIiMhIMPgRERERGQkGPyIiIiIjweBHREREZCQY/IiIiIiMBIMfERERkZFg8CMiIqIGLTk5GaGhoXB3d0ejRo3QqlUrhIeHo6CgoErbCyHQr18/yGQy7Nq1S2+dTCYrtXz//fd6fTZv3owuXbqgcePGUKvVGDduHP766y+9PllZWQgLC4NarYZCocBTTz2FPXv2VPtYIyMj0blzZ1hYWMDR0RFhYWHV2p7Bj4ioHgkJCdH7F4ydnR0CAwNx/vx5qc/f/+X0cH9LS0u0adMGISEhiI+PN9BRENWN3r17IzIyslT75cuXodPpsGbNGiQkJGDJkiVYvXo13nvvvSqNu3TpUshksnLXR0REIC0tTVoe/v758ePHMXr0aISGhiIhIQHbt2/H6dOnMX78eKlPQUEBXnjhBSQnJ2PHjh1ITEzE2rVr0axZs2od/5dffon3338fs2bNQkJCAvbv34+AgIBqjcFv9RIR1SNjxowRgYGBIi0tTaSlpYlffvlFDBgwQLi4uEh9AIh///vfer8jIiJEWlqaSEpKEvv27RNBQUHCxMREbNiwwUBHQlT7evXqJSIiIqrUd9GiRcLd3b3Sfr/88oto1qyZSEtLK/XPlijjn7e/+/zzz0XLli312r766ivRrFkz6feqVatEy5YtRUFBQbnjFBcXi/nz5ws3NzdhYWEhOnfuLLZv3y6tz8zMFI0aNRL79++v9JgqwjN+RET1jEKhgEqlgkqlgqenJ2bNmoUbN27g9u3b5W5jY2MDlUoFNzc3vPjii9ixYweCg4MxadIk3Llz57HWT1QfaLVa2NraVtjn3r17eO2117By5UqoVKpy+4WFhcHe3h5eXl5Yv349/pcH/8fX1xc3btzAnj17IIRAeno6duzYgf79+0t9fvzxR/j6+iIsLAxOTk7o2LEj5s+fj+LiYqnPggULsHHjRqxevRoJCQmYOnUqRo0ahcOHDwMAYmJioNPpkJKSgnbt2qF58+Z45ZVXcOPGjWrNC4MfEVE9lpubi02bNqF169aws7Or1rZTp05FTk4OYmJi6qw+ovroypUrWL58Od54440K+02dOhU9evTA4MGDy+3z8ccfY9u2bYiJiUFQUBDeeustLF++XFr/7LPPYvPmzXj11Vdhbm4OlUoFa2trrFy5Uurz559/YseOHSguLsaePXvw4Ycf4osvvsCnn34KAMjPz8f8+fOxfv16BAQEoGXLlggJCcGoUaOwZs0aaQydTof58+dj6dKl2LFjBzIzM/HCCy9U+V5GgJd6iYjqlTFjxggTExNhaWkpLC0tBQChVqtFfHy81KesS71lXYq6f/++ACAWLlz42Oonqk3z5s2T/lmwtLQUcrlcKBQKvbZr167pbXPz5k3RqlUrERoaWuHY//nPf0Tr1q1FTk6O1FbZZV0hhPjwww9F8+bNpd8JCQlCrVaLRYsWiV9//VVER0eLTp06iXHjxkl92rRpI1xcXERRUZHU9sUXXwiVSiWEEOLixYsCgN5xWVpaCjMzM+Hl5SXNBQCxb98+aYxbt24JuVwuoqOjqzCb/2Na9YhIRESPQ58+fbBq1SoAwJ07d/D111+jX79+OH36NFq0aFHlcUouR1V00zpRfTZx4kS88sor0u/g4GAEBQVh2LBhUpuzs7P016mpqejTpw969OiBb775psKxY2NjcfXqVdjY2Oi1BwUF4fnnn8ehQ4fK3M7b2xuffPIJ8vPzoVAosGDBAjz77LOYMWMGAKBz586wtLTE888/j08//RRqtRpqtRpmZmYwMTGRxmnXrh00Gg0KCgqQm5sLAPjpp59KPfChUCgAAGq1GgDQvn17aZ2DgwPs7e1x/fr1Co/1YQx+RET1jKWlJVq3bi39XrduHaytrbF27Vrp0lBVXLp0CQDg7u5eJ3US1TVbW1u9+/QaNWoER0dHvX8+SqSkpKBPnz7o2rUrIiIiIJdXfDfbrFmz8I9//EOvrVOnTliyZAleeumlcrc7d+4cmjZtKgWye/fuwdRUP06VBLyS//h69tlnsWXLFuh0Oqmu33//HWq1Gubm5mjfvj0UCgWuX7+OXr16lbnfZ599FgCQmJiI5s2bAwAyMzORkZFRrf8gZPAjIqrnZDIZ5HI57t+/X63tli5dCqVSCX9//zqrjag+SElJQe/evdGiRQssXrxY70Gokoc2UlJS4Ofnh40bN8LLy0t6gOrvXF1dpf9Y2r17N9LT0+Hj4wMLCwvExMRg/vz5ePfdd6X+L730EsaPH49Vq1YhICAAaWlpeOedd+Dl5SWdjXzzzTexYsUKTJkyBZMnT8Yff/yB+fPn4+233wYANGnSBO+++y6mTp0KnU6H5557DlqtFsePH4dSqcSYMWPw1FNPYfDgwZgyZQq++eYbKJVKzJ49Gx4eHujTp0+V54rBj4ionsnPz4dGowEeXOpdsWIFcnNzKzwLkZWVBY1Gg/z8fPz+++9Ys2YNdu3ahY0bN5a6lEX0pImJicGVK1dw5coV6WxYiZKzboWFhUhMTMS9e/eqPK6ZmRlWrlyJqVOnQgiB1q1b48svv9R7R19ISAhycnKwYsUKTJ8+HTY2Nujbty8WLlwo9XFxccG+ffswdepUdO7cGc2aNcOUKVMwc+ZMqc8nn3wCBwcHLFiwAH/++SdsbGzwzDPP6L2LcOPGjZg6dSoGDBgAuVyOXr16ITo6GmZmZlU+Jpl4+JlkIiIyqJCQEGzYsEH63aRJE3h4eGDmzJkICgqCTqeDiYkJdu/ejYEDBwJ/u4fPwsICzZo1w3PPPYe3334bzzzzjEGOg4jqJwY/IqIGRKPRQK1W48yZM+jWrZuhyyGiBoaXeomIGgAhBK5du4bFixdLL4AlIqouBj8iogZAq9Wibdu2aNeuHb7//ntYWFgYuiQiaoB4qZeIiIjISPCTbURERERGgsGPiIiIyEgw+BEREREZCQY/IiIiIiPB4EdERERkJBj8iIiIiIwEgx8RERGRkWDwIyIiIjIS/wcfwiqs7vETyAAAAABJRU5ErkJggg==", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "data.plot_quantity(\"BJD\", \"DRS_RV\", show=True)\n" + ] + }, { "cell_type": "markdown", "metadata": {}, @@ -221,7 +192,7 @@ }, { "cell_type": "code", - "execution_count": 37, + "execution_count": 11, "metadata": {}, "outputs": [ { @@ -230,7 +201,7 @@ "Text(0, 0.5, 'Flux')" ] }, - "execution_count": 37, + "execution_count": 11, "metadata": {}, "output_type": "execute_result" }, @@ -263,7 +234,7 @@ ], "metadata": { "kernelspec": { - "display_name": ".venv", + "display_name": "ASTRA", "language": "python", "name": "python3" }, diff --git a/docs/spectra/selecting_spectra.ipynb b/docs/spectra/selecting_spectra.ipynb index ccfdaad..ec0cac3 100644 --- a/docs/spectra/selecting_spectra.ipynb +++ b/docs/spectra/selecting_spectra.ipynb @@ -88,7 +88,7 @@ ], "metadata": { "kernelspec": { - "display_name": ".venv", + "display_name": "ASTRA", "language": "python", "name": "python3" }, diff --git a/pyproject.toml b/pyproject.toml index d254077..db17cab 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [project] name = "ASTRA-spectra" -version = "1.2.6" +version = "1.2.7" description = "Interface to open stellar spectra and construct stellar and telluric models" readme = "README.md" requires-python = ">=3.11,<=3.12.8" diff --git a/src/ASTRA/Instruments/ESO_PIPELINE.py b/src/ASTRA/Instruments/ESO_PIPELINE.py index c8c5bf1..fc42565 100644 --- a/src/ASTRA/Instruments/ESO_PIPELINE.py +++ b/src/ASTRA/Instruments/ESO_PIPELINE.py @@ -113,7 +113,7 @@ def __init__( "INS NAME": "INSTRUME", "INS MODE": f"HIERARCH {KW_identifier} INS MODE", "PROG ID": f"HIERARCH {KW_identifier} OBS PROG ID", - "OBS NAME": "HIERARCH ESO OBS NAME", + "OBS NAME": f"HIERARCH {KW_identifier} OBS NAME", } if override_KW_map is not None: for key, value in override_KW_map.items(): @@ -131,7 +131,11 @@ def __init__( file_path=file_path, frameID=frameID, KW_map=KW_map, - available_indicators=(available_indicators if override_indicators is None else override_indicators), + available_indicators=( + available_indicators + if override_indicators is None + else override_indicators + ), user_configs=user_configs, reject_subInstruments=reject_subInstruments, quiet_user_params=quiet_user_params, @@ -156,7 +160,9 @@ def __init__( if key in self.file_path.stem: break else: - raise custom_exceptions.InvalidConfiguration(f"{self.name} can't recognize {self.file_path}") + raise custom_exceptions.InvalidConfiguration( + f"{self.name} can't recognize {self.file_path}" + ) def load_instrument_specific_KWs(self, header) -> None: if self._internal_configs["use_old_pipeline"]: @@ -169,7 +175,9 @@ def load_instrument_specific_KWs(self, header) -> None: version = drs_version.split("/")[-1] # TODO: Check the number of the DRS version that has this fixed - version_sum = lambda a: sum(a * 10**b for a, b in zip(map(int, a.split(".")[::-1]), range(3))) + version_sum = lambda a: sum( + a * 10**b for a, b in zip(map(int, a.split(".")[::-1]), range(3)) + ) if version_sum(version) < version_sum("3.2.1"): # For older versions we always need to use the @@ -183,9 +191,13 @@ def load_instrument_specific_KWs(self, header) -> None: if berv_factor is not None: logger.warning("Recomputing the BERV from BERV_factor keyword") new_berv = (berv_factor - 1) * SPEED_OF_LIGHT * kilometer_second - diff_berv = (new_berv - self.observation_info["BERV"]).to(meter_second).value + diff_berv = ( + (new_berv - self.observation_info["BERV"]).to(meter_second).value + ) if abs(diff_berv) > 10: - logger.warning(f"Difference between BERV and BERV_factor is of {diff_berv} [m/s]") + logger.warning( + f"Difference between BERV and BERV_factor is of {diff_berv} [m/s]" + ) self.observation_info["BERV"] = new_berv @@ -250,7 +262,9 @@ def load_telemetry_info(self, header): } for name, endKW in ambi_KWs.items(): - self.observation_info[name] = header[f"HIERARCH {self.KW_identifier} METEO {endKW}"] + self.observation_info[name] = header[ + f"HIERARCH {self.KW_identifier} METEO {endKW}" + ] if "temperature" in name: # store temperature in KELVIN for TELFIT self.observation_info[name] = convert_temperature( self.observation_info[name], @@ -259,7 +273,9 @@ def load_telemetry_info(self, header): ) if self.observation_info["relative_humidity"] == 255: - logger.warning(f"{self.name} has an invalid value in the humidity sensor...") + logger.warning( + f"{self.name} has an invalid value in the humidity sensor..." + ) self.observation_info["relative_humidity"] = np.nan self.observation_info["airmass"] = header["AIRMASS"] @@ -271,13 +287,19 @@ def _load_ESO_DRS_KWs(self, header): ) # Load BERV info + previous RV - self.observation_info["MAX_BERV"] = header[f"HIERARCH {self.KW_identifier} QC BERVMAX"] * kilometer_second - self.observation_info["BERV"] = header[f"HIERARCH {self.KW_identifier} QC BERV"] * kilometer_second + self.observation_info["MAX_BERV"] = ( + header[f"HIERARCH {self.KW_identifier} QC BERVMAX"] * kilometer_second + ) + self.observation_info["BERV"] = ( + header[f"HIERARCH {self.KW_identifier} QC BERV"] * kilometer_second + ) berv_factor = header.get(f"HIERARCH {self.KW_identifier} QC BERV FACTOR", None) self.observation_info["BERV_FACTOR"] = berv_factor - self.observation_info["DRS_RV"] = header[f"HIERARCH {self.KW_identifier} QC CCF RV"] * kilometer_second + self.observation_info["DRS_RV"] = ( + header[f"HIERARCH {self.KW_identifier} QC CCF RV"] * kilometer_second + ) self.observation_info["DRS_RV_ERR"] = ( header[f"HIERARCH {self.KW_identifier} QC CCF RV ERROR"] * kilometer_second ) @@ -301,7 +323,9 @@ def _load_ESO_DRS_KWs(self, header): if obs_date.hour < 16: # noqa: PLR2004 # If before 4pm, the observation was from previous day obs_date = obs_date - datetime.timedelta(days=1) - self.observation_info["DATE_NIGHT"] = datetime.datetime.strftime(obs_date, r"%Y-%m-%d") + self.observation_info["DATE_NIGHT"] = datetime.datetime.strftime( + obs_date, r"%Y-%m-%d" + ) def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None): if self._internal_configs["use_old_pipeline"]: @@ -313,11 +337,15 @@ def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None): self.wavelengths = hdulist["WAVEDATA_VAC_BARY"].data self.qual_data = hdulist["QUALDATA"].data - SCIDATA_KEY = "SCIDATA" if overload_SCIDATA_key is None else overload_SCIDATA_key + SCIDATA_KEY = ( + "SCIDATA" if overload_SCIDATA_key is None else overload_SCIDATA_key + ) ERRDATA_KEY = "ERRDATA" if self._internal_configs["Telluric_Corrected"]: - logger.info("Loading S2D file from a non-DRS source (telluric corrected file)") + logger.info( + "Loading S2D file from a non-DRS source (telluric corrected file)" + ) SCIDATA_KEY += "_CORR" ERRDATA_KEY += "_CORR" @@ -327,7 +355,12 @@ def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None): if self._internal_configs["apply_FluxCorr"]: logger.debug("Starting chromatic flux correction") keyword = f"HIERARCH {self.KW_identifier} QC ORDER%d FLUX CORR" - flux_corr = np.array([hdulist[0].header[keyword % o] for o in range(1, self.N_orders + 1)]) + flux_corr = np.array( + [ + hdulist[0].header[keyword % o] + for o in range(1, self.N_orders + 1) + ] + ) fit_nb = (flux_corr != 1.0).sum() ignore = self.N_orders - fit_nb @@ -341,8 +374,12 @@ def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None): # corr_model = np.zeros_like(hdu[5].data, dtype=np.float32) corr_model = np.polyval(coeff, hdulist[5].data) - corr_model[flux_corr == 1] = 1 # orders where the CORR FACTOR are 1 do not have correction! - self.spectra = self.spectra / corr_model # correct from chromatic variations + corr_model[flux_corr == 1] = ( + 1 # orders where the CORR FACTOR are 1 do not have correction! + ) + self.spectra = ( + self.spectra / corr_model + ) # correct from chromatic variations self.flux_atmos_balance_corrected = True # TODO: understand if we want to include the factor in uncertainties or not! # self.uncertainties = self.uncertainties / corr_model # maintain the SNR in the corrected spectrum @@ -354,7 +391,9 @@ def load_ESO_DRS_S2D_data(self, overload_SCIDATA_key=None): # / corr_model if self._internal_configs["apply_FluxBalance_Norm"]: - logger.info("Normalizing the flux balance distribution due to dispersion") + logger.info( + "Normalizing the flux balance distribution due to dispersion" + ) # The physical sizes of the pixels (on the CCD) are the same # The flux that reaches eeach pixel is different, due to dispersion # The spectra will have a trend, even after removing the instrumental effect @@ -384,8 +423,12 @@ def load_ESO_DRS_S1D_data(self): logger.warning("SBART using air wavelengths!") self.wavelengths = full_data[wave_kw].reshape((1, self.array_size[1])) - self.spectra = full_data["flux"].reshape((1, self.array_size[1])).astype(np.float64) - self.uncertainties = full_data["error"].reshape((1, self.array_size[1])).astype(np.float64) + self.spectra = ( + full_data["flux"].reshape((1, self.array_size[1])).astype(np.float64) + ) + self.uncertainties = ( + full_data["error"].reshape((1, self.array_size[1])).astype(np.float64) + ) self.qual_data = full_data["quality"].reshape((1, self.array_size[1])) self.build_mask(bypass_QualCheck=False) @@ -420,7 +463,8 @@ def check_header_QC_ESO_DRS(self, header): if not self.is_skysub: extra_checks = { f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT CHECK": 0, - f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT FLUX_RATIO CHECK": 0, + f"HIERARCH {self.KW_identifier}" + + " QC SCIRED DRIFT FLUX_RATIO CHECK": 0, f"HIERARCH {self.KW_identifier}" + " QC SCIRED DRIFT CHI2 CHECK": 0, } nonfatal_QC_flags = {**nonfatal_QC_flags, **extra_checks} @@ -448,17 +492,23 @@ def check_header_QC_ESO_DRS(self, header): self._status.store_warning(KW_WARNING(msg)) if self._status.number_warnings > 0: - logger.warning("Found {} warning flags in the header KWs", self._status.number_warnings) + logger.warning( + "Found {} warning flags in the header KWs", self._status.number_warnings + ) espdr_to_num = lambda x: int("".join(x.split("."))) espdrversion = header["ESO PRO REC1 PIPE ID"].split("/")[-1] if self._internal_configs["USE_APPROX_BERV_CORRECTION"]: if espdr_to_num(espdrversion) >= espdr_to_num("3.2.0"): - logger.critical(f"Using approximated BERV correction in espdr/{espdrversion}") + logger.critical( + f"Using approximated BERV correction in espdr/{espdrversion}" + ) else: if espdr_to_num(espdrversion) < espdr_to_num("3.2.0"): - logger.critical(f"Not using approximated BERV correction in espdr/{espdrversion}") + logger.critical( + f"Not using approximated BERV correction in espdr/{espdrversion}" + ) @property def bare_fname(self) -> str: diff --git a/src/ASTRA/Instruments/HARPSN.py b/src/ASTRA/Instruments/HARPSN.py new file mode 100644 index 0000000..a07b1d8 --- /dev/null +++ b/src/ASTRA/Instruments/HARPSN.py @@ -0,0 +1,395 @@ +import datetime +from pathlib import Path +from typing import Any, Dict, Optional + +import numpy as np +from astropy.coordinates import EarthLocation +from astropy.io import fits +from loguru import logger +from ASTRA.Instruments.ESO_PIPELINE import ESO_PIPELINE +from ASTRA.utils import custom_exceptions +from ASTRA.utils.air_to_vac import airtovac +from ASTRA.status.flags import ( + ERROR_THRESHOLD, + FATAL_KW, + KW_WARNING, + MISSING_DATA, + NAN_DATA, + SATURATION, + SUCCESS, +) + +from ASTRA.utils.ASTRAtypes import UI_DICT, UI_PATH + +from ASTRA.utils.units import kilometer_second, meter_second +from scipy.constants import convert_temperature + + +class HARPSN(ESO_PIPELINE): + """ + Interface to handle HARPSN data; S1D **not** supported! + + This class also defines 2 sub-Instruments: + + * HARPSN - Until the ends of time (hopefully) + + + The steps to load the S2D data are described in the HARPS `DRS manual `_. The summary is: + + - Construct the wavelength solution & correct from BERV + - Load instrumental drift + - Construct flux noises: + + - Table 10 of `the user manual `_ gives max RON of 7.07 for red detector + - Noise = sqrt(obj + sky + n*dark*expTime + nBinY*ron^2) + + **User parameters:** + + Currently there are no HARPSN-specific parameters + + *Note:* Check the **User parameters** of the parent classes for further customization options of SBART + + """ + + sub_instruments = { + # "HARPSpre": datetime.datetime.strptime("2015-05-29", r"%Y-%m-%d"), + "HARPSN": datetime.datetime.max, + } + + _name = "HARPSN" + _default_params = ESO_PIPELINE._default_params + + def __init__( + self, + file_path, + user_configs: Optional[Dict[str, Any]] = None, + reject_subInstruments=None, + frameID=None, + quiet_user_params: bool = True, + ): + """ + + Parameters + ---------- + file_path + Path to the S2D (or S1D) file. + user_configs + Dictionary whose keys are the configurable options of ESPRESSO (check above) + reject_subInstruments + Iterable of subInstruments to fully reject + frameID + ID for this observation. Only used for organization purposes by :class:`~SBART.data_objects.DataClass` + """ + + logger.info("Creating frame from: {}".format(file_path)) + + coverage = [390, 700] + search_status = SUCCESS + if user_configs.get("use_old_pipeline", False): # For the old pipeline! + override_KW_map = { + "OBJECT": "HIERARCH TNG OBS TARG NAME", # "OBJECT", + "BJD": "HIERARCH TNG DRS BJD", + "MJD": "MJD-OBS", + "ISO-DATE": "DATE-OBS", + "DRS-VERSION": "HIERARCH TNG DRS VERSION", + "RA": "RA", + "DEC": "DEC", + "SPEC_TYPE": "HIERARCH TNG TEL TARG SPTYPE", + "MD5-CHECK": "HIERARCH TNG DRS BJD", # Missing the MD5 on the old pipe, so this is a stopgap + } + file_path, self.ccf_path, search_status = self.find_files(file_path) + override_indicators = ("CONTRAST", "FWHM") + + # The flag for the BERV correction is set at the time of loading the S2D arrays! + + else: # For the new pipeline + override_KW_map = {"OBJECT": "HIERARCH TNG OBS TARG NAME"} + override_indicators = None + + super().__init__( + inst_name="HARPSN", + array_size={"S2D": (69, 4096)}, + file_path=file_path, + KW_identifier="TNG", + frameID=frameID, + override_KW_map=override_KW_map, + user_configs=user_configs, + reject_subInstruments=reject_subInstruments, + quiet_user_params=quiet_user_params, + override_indicators=override_indicators, + ) + + if ( + user_configs.get("use_old_pipeline", False) + and not search_status.is_good_flag + ): + self.add_to_status(search_status) + + self.instrument_properties["wavelength_coverage"] = coverage + + self.instrument_properties["resolution"] = 115_000 + self.instrument_properties["EarthLocation"] = EarthLocation.of_site( + "Roque de los Muchachos" + ) + + # https://tngweb.tng.iac.es/weather/current + self.instrument_properties["site_pressure"] = 770 + + def check_header_QC_old_DRS(self, header): + logger.info("Currently missing QC checks for the old DRS") + + def _load_old_DRS_KWs(self, header): + if not self._internal_configs["use_old_pipeline"]: + raise custom_exceptions.InvalidConfiguration( + "Can't load data from old pipeline without the config" + ) + + self.observation_info["MAX_BERV"] = ( + header["HIERARCH TNG DRS BERVMX"] * kilometer_second + ) + self.observation_info["BERV"] = ( + header["HIERARCH TNG DRS BERV"] * kilometer_second + ) + + # Environmental KWs for telfit (also needs airmassm previously loaded) + + ambi_KWs = { + "relative_humidity": "HUMIDITY", + "ambient_temperature": "TEMP10M", + } + + for name, endKW in ambi_KWs.items(): + self.observation_info[name] = header[f"HIERARCH TNG METEO {endKW}"] + if "temperature" in name: # store temperature in KELVIN for TELFIT + self.observation_info[name] = convert_temperature( + self.observation_info[name], old_scale="Celsius", new_scale="Kelvin" + ) + + if self.observation_info["relative_humidity"] == 255: + logger.warning( + f"{self.name} has an invalid value in the humidity sensor..." + ) + self.observation_info["relative_humidity"] = np.nan + + for order in range(self.N_orders): + self.observation_info["orderwise_SNRs"].append( + header[f"HIERARCH TNG DRS SPE EXT SN{order}"] + ) + + self.observation_info["airmass"] = header["AIRMASS"] + + bad_drift = False + try: + flag = "HIERARCH TNG DRS DRIFT QC" + if header[flag].strip() != "PASSED": + bad_drift = True + msg = f"QC flag {flag} meets the bad value" + logger.warning(msg) + self._status.store_warning(KW_WARNING(msg)) + else: + # self.logger.info("DRIFT QC has passed") + drift = header["HIERARCH TNG DRS DRIFT RV USED"] * meter_second + drift_err = header["HIERARCH TNG DRS DRIFT NOISE"] * meter_second + except Exception as e: + bad_drift = True + logger.warning("DRIFT KW does not exist") + + if bad_drift: + logger.warning( + "Due to previous drift-related problems, setting it to zero [m/s]" + ) + drift = 0 * meter_second + drift_err = 0 * meter_second + + self.observation_info["drift"] = drift + self.observation_info["drift_ERR"] = drift_err + self._load_ccf_data() + + def _load_ccf_data(self) -> None: + """ + Load the necessarfy CCF data from the file! + """ + if not self._internal_configs["use_old_pipeline"]: + raise custom_exceptions.InvalidConfiguration( + "Can't load data from old pipeline without the config" + ) + + logger.debug("Loading data from the ccf file") + header = fits.getheader(self.ccf_path) + + for key in self.available_indicators: + full_key = "HIERARCH TNG DRS CCF " + key + self.observation_info[key] = header[full_key] + + self.observation_info["DRS_RV"] = ( + header["HIERARCH TNG DRS CCF RV"] * kilometer_second + ) + + RV_err = header.get("HIERARCH TNG DRS DVRMS", None) + if RV_err is None: + logger.critical( + "Couldn't find DRS error from the appropriate KW. Estimating it with photon noise!" + ) + RV_err = np.sqrt( + header["HIERARCH TNG DRS CAL TH ERROR"] ** 2 + + + # hdulist[0].header['HIERARCH ESO DRS DRIFT NOISE']**2 + + (1000 * header["HIERARCH TNG DRS CCF NOISE"]) ** 2 + ) + + self.observation_info["DRS_RV_ERR"] = RV_err * meter_second + + def build_HARPS_wavelengths(self, hdr): + """ + Compute the wavelength solution to this given spectra (EQ 4.1 of DRS manual) + Convert from air wavelenbgths to vacuum + """ + if not self._internal_configs["use_old_pipeline"]: + raise custom_exceptions.InvalidConfiguration( + "Can't construct wavelengths for new pipeline" + ) + + # degree of the polynomial + d = hdr["HIERARCH TNG DRS CAL TH DEG LL"] + # number of orders + omax = hdr.get("HIERARCH TNG DRS CAL LOC NBO", self.array_size[0]) + xmax = self.array_size[1] + + # matrix X: + # + # axis 0: the entry corresponding to each coefficient + # axis 1: each pixel number + + x = np.empty((d + 1, xmax), "int64") + x[0].fill(1) # x[0,*] = x^0 = 1,1,1,1,1,... + x[1] = np.arange(xmax) + + for i in range(1, d): + x[i + 1] = x[i] * x[1] + + # matrix A: + # + # axis 0: the different orders + # axis 1: all coefficients for the given order + + A = np.reshape( + [ + hdr["HIERARCH TNG DRS CAL TH COEFF LL" + str(i)] + for i in range(omax * (d + 1)) + ], + (omax, d + 1), + ) # slow 30 ms + + # the wavelengths for each order are a simple dot product between the coefficients and pixel-wise data (X) + wavelengths = np.dot(A, x) + + vacuum_wavelengths = airtovac(wavelengths) + return vacuum_wavelengths + + def load_old_DRS_S2D(self): + """ + load the data from the old HARPS-N pipeline. This will be mainly used for the comparison with the + HARPS-TERRA pipeline + """ + if not self._internal_configs["use_old_pipeline"]: + raise custom_exceptions.InvalidConfiguration( + "Can't load data from old pipeline without the config" + ) + + if self._internal_configs["use_old_pipeline"]: + self.is_BERV_corrected = False + + with fits.open(self.file_path) as hdulist: + # Compute the wavelength solution + BERV correction + wave_from_file = self.build_HARPS_wavelengths(hdulist[0].header) + + sci_data = hdulist[0].data # spetra from all orders + + # photon noise + estimate of max value for the rest + # from ETC calculator the readout noise should be the largest contribution + # assuming that it is of ~7e- (equal to manual) it should have a maximum contribution + # of 200 + flux_errors = np.sqrt(250 + np.abs(sci_data, dtype=float)) + + # Validate for overflows and missing data + quality_data = np.zeros(sci_data.shape) + quality_data[np.where(np.isnan(sci_data))] = NAN_DATA.code + quality_data[np.where(sci_data > 300000)] = SATURATION.code + quality_data[np.where(sci_data < -3 * flux_errors)] = ERROR_THRESHOLD.code + + self.spectra = sci_data.astype(np.float64) + self.wavelengths = wave_from_file + self.qual_data = quality_data + self.uncertainties = flux_errors.astype(np.float64) + + self.build_mask(bypass_QualCheck=False) + + def close_arrays(self): + """ + Reset the BERV correction flag if we are using the old pipeline version! + Returns + ------- + + """ + super().close_arrays() + if self._internal_configs["use_old_pipeline"]: + self.is_BERV_corrected = False + + def find_files(self, file_name: UI_PATH): + """ + Find the CCF and the S2D files, which should be stored inside the same folder + """ + logger.debug("Searching for the ccf and e2ds files") + + if not isinstance(file_name, Path): + file_name = Path(file_name) + + search_status = MISSING_DATA("Missing the ccf file") + ccf_path = None + + if file_name.is_dir(): + logger.debug("Received a folder, searching inside for necessary files") + # search for e2ds file + folder_name = file_name + + e2ds_files = file_name.glob("**/*e2ds_A.fits") + ccf_files = file_name.glob("**/*ccf_*_A.fits") + + for name, elems in [("e2ds_A", e2ds_files), ("ccf", ccf_files)]: + if len(elems) > 1: + msg = f"HARPS data only received folder name and it has more than 1 {name} file in it" + raise custom_exceptions.InvalidConfiguration(msg) + + if len(elems) < 1: + msg = f"HARPS data only received folder name and it has no {name} file in it" + raise custom_exceptions.InvalidConfiguration(msg) + + e2ds_path = e2ds_files[0] + ccf_path = ccf_files[0] + search_status = SUCCESS("Found all input files") + else: + logger.debug( + "Received path of E2DS file; searching for CCF with matching name" + ) + folder_name = file_name.parent + e2ds_path = file_name + file_start, *_ = file_name.stem.split("_") + + found_CCF = False + ccf_files = folder_name.glob("*ccf*_A.fits") + + for file in ccf_files: + if file_start in file.name: + ccf_path = file + found_CCF = True + + if found_CCF: + logger.info("Found CCF file: {}".format(ccf_path)) + search_status = SUCCESS("Found CCF file") + else: + logger.critical( + "Was not able to find CCF file. Marking frame as invalid" + ) + ccf_path = "" + + return e2ds_path, ccf_path, search_status diff --git a/src/ASTRA/Instruments/__init__.py b/src/ASTRA/Instruments/__init__.py index 762cc96..9ea20b2 100644 --- a/src/ASTRA/Instruments/__init__.py +++ b/src/ASTRA/Instruments/__init__.py @@ -19,10 +19,12 @@ "MAROONX", "CARMENES", "SimulatedSpirou", + "HARPSN", ] from .ESPRESSO import ESPRESSO from .HARPS import HARPS +from .HARPSN import HARPSN from .MAROONX import MAROONX from .CARMENES import CARMENES from .SimulatedSpirou import SimulatedSpirou @@ -31,6 +33,7 @@ instrument_dict = { "ESPRESSO": ESPRESSO, "HARPS": HARPS, + "HARPSN": HARPSN, "MAROONX": MAROONX, "SimulatedSpirou": SimulatedSpirou, "CARMENES": CARMENES, diff --git a/src/ASTRA/__init__.py b/src/ASTRA/__init__.py index 7edb00d..d95de3b 100644 --- a/src/ASTRA/__init__.py +++ b/src/ASTRA/__init__.py @@ -1,6 +1,6 @@ """ASTRA - interface for spectra.""" -version = "1.2.6" +version = "1.2.7" __version__ = version.replace(".", "-") __version_info__ = (int(i) for i in __version__.split("-")) diff --git a/src/ASTRA/data_objects/DataClass.py b/src/ASTRA/data_objects/DataClass.py index eeb3188..b6d176b 100644 --- a/src/ASTRA/data_objects/DataClass.py +++ b/src/ASTRA/data_objects/DataClass.py @@ -3,7 +3,6 @@ import hashlib from collections import defaultdict from collections.abc import Iterable -import hashlib from pathlib import Path from typing import ( TYPE_CHECKING, @@ -17,6 +16,8 @@ Union, ) +import astropy.units as u +import matplotlib.pyplot as plt import numpy as np import ujson as json from tabletexifier import Table @@ -37,7 +38,7 @@ from ASTRA.utils.BASE import BASE from ASTRA.utils.choices import DISK_SAVE_MODE from ASTRA.utils.custom_exceptions import FrameError, InvalidConfiguration, NoDataError -from ASTRA.utils.units import kilometer_second +from ASTRA.utils.units import convert_data, kilometer_second, meter_second if TYPE_CHECKING: from ASTRA.base_models.Frame import Frame @@ -807,7 +808,9 @@ def get_valid_frameIDS(self) -> list[int]: return out - def get_invalid_frameIDs(self, subinstrument: None | str = None) -> list[int]: # noqa: N802 + def get_invalid_frameIDs( + self, subinstrument: None | str = None + ) -> list[int]: # noqa: N802 """Get a list of the invalid frameIDs (by default, all of them). Args: @@ -982,6 +985,47 @@ def show_loadedData_table(self) -> Table: return tab + def plot_quantity( + self, + xx_var: str, + yy_var: str, + axis=None, + return_fig_and_axis: bool = False, + show: bool = False, + ): + """Plot any quantity from the frames against each other.""" + + if axis is None: + N = len(self.get_available_subInstruments()) + fig, axis = plt.subplots( + nrows=N, + sharex=True, + ) + if N == 1: + axis = [axis] + + for index, inst in enumerate(self.get_available_subInstruments()): + xx = self.collect_KW_observations(xx_var, [inst]) + yy = self.collect_KW_observations(yy_var, [inst]) + if len(xx) > 0: + if isinstance(xx[0], u.Quantity): + xx = convert_data(xx, as_value=True) + if len(yy) > 0: + if isinstance(yy[0], u.Quantity): + yy = convert_data(yy, as_value=True) + + axis[index].set_title(inst) + axis[index].scatter(xx, yy) + axis[index].set_xlabel(xx_var) + axis[index].set_ylabel(yy_var) + + fig.tight_layout() + if show: + plt.show() + + if return_fig_and_axis: + return fig, axis + def load_instrument_extra_information(self) -> None: """See if the given instrument is one of the ones that has extra information to load. diff --git a/src/ASTRA/spectral_modelling/GPmodel.py b/src/ASTRA/spectral_modelling/GPmodel.py index 0751086..ef51b43 100644 --- a/src/ASTRA/spectral_modelling/GPmodel.py +++ b/src/ASTRA/spectral_modelling/GPmodel.py @@ -3,9 +3,9 @@ from typing import NoReturn import numpy as np -from ASTRAModelParameters.Parameter import JaxComponent -from ASTRAutils import custom_exceptions -from ASTRAutils.paths_tools import build_filename +from ASTRA.ModelParameters.Parameter import JaxComponent +from ASTRA.utils import custom_exceptions +from ASTRA.utils.paths_tools import build_filename from ASTRA import astra_logger as logger @@ -23,9 +23,9 @@ # To avoid issues with possible needs of jax import numpy as jnp -from ASTRAspectral_modelling.modelling_base import ModellingBase -from ASTRAutils.status_codes import INTERNAL_ERROR, SUCCESS -from ASTRAutils.UserConfigs import ( +from ASTRA.spectral_modelling.modelling_base import ModellingBase +from ASTRA.status.flags import INTERNAL_ERROR, SUCCESS +from ASTRA.utils.UserConfigs import ( BooleanValue, DefaultValues, Positive_Value_Constraint, @@ -51,16 +51,24 @@ class GPSpecModel(ModellingBase): # TODO: confirm the kernels that we want to allow _default_params = ModellingBase._default_params + DefaultValues( - GP_KERNEL=UserParam("Matern-5_2", constraint=ValueFromIterable(["Matern-5_2", "Matern-3_2"])), + GP_KERNEL=UserParam( + "Matern-5_2", constraint=ValueFromIterable(["Matern-5_2", "Matern-3_2"]) + ), FORCE_MODEL_GENERATION=UserParam(False, constraint=BooleanValue), - POSTERIOR_CHARACTERIZATION=UserParam("minimize", constraint=ValueFromIterable(["minimize", "MCMC"])), + POSTERIOR_CHARACTERIZATION=UserParam( + "minimize", constraint=ValueFromIterable(["minimize", "MCMC"]) + ), OPTIMIZATION_MAX_ITER=UserParam(1000, constraint=Positive_Value_Constraint), ) def __init__(self, obj_info, user_configs): possible_folders = {} - for kern in self._default_params["GP_KERNEL"].existing_constraints.available_options: - possible_folders[f"modelling_information_{kern}"] = f"_SpecModelParams/{kern}" + for kern in self._default_params[ + "GP_KERNEL" + ].existing_constraints.available_options: + possible_folders[f"modelling_information_{kern}"] = ( + f"_SpecModelParams/{kern}" + ) super().__init__( obj_info=obj_info, @@ -98,7 +106,9 @@ def __init__(self, obj_info, user_configs): def _check_dependencies(self): if MISSING_TINYGP: - raise custom_exceptions.InternalError("Missing the custom tinygp installation for GP") + raise custom_exceptions.InternalError( + "Missing the custom tinygp installation for GP" + ) def _get_model_storage_filename(self) -> str: """Construct the storage filename for the model parameters for this observation! @@ -121,7 +131,9 @@ def load_previous_model_results_from_disk(self, model_component_in_use): # TODO: ensure that loaded information is in accordance with the current parameters! return super().load_previous_model_results_from_disk(model_component_in_use) - def generate_model_from_order(self, og_lambda, og_spectra, og_err, new_wavelengths, order) -> NoReturn: + def generate_model_from_order( + self, og_lambda, og_spectra, og_err, new_wavelengths, order + ) -> NoReturn: """Fit the stellar spectrum from a given order. If it has already been recomputed (or if it has previously failed) does nothing @@ -145,22 +157,32 @@ def generate_model_from_order(self, og_lambda, og_spectra, og_err, new_wavelengt return try: - solution_array, result_flag = self._launch_GP_fit(og_lambda, og_spectra, og_err, new_wavelengths, order) + solution_array, result_flag = self._launch_GP_fit( + og_lambda, og_spectra, og_err, new_wavelengths, order + ) except Exception as e: msg = f"Unknown error found when fitting GP to order {order}: {traceback.print_tb(e.__traceback__)}" logger.critical(msg) result_flag = INTERNAL_ERROR(msg) - solution_array = [np.nan for _ in self._modelling_parameters.get_enabled_params()] + solution_array = [ + np.nan for _ in self._modelling_parameters.get_enabled_params() + ] - raise custom_exceptions.StopComputationError("Unknown error encounterd") from e + raise custom_exceptions.StopComputationError( + "Unknown error encounterd" + ) from e - self._modelling_parameters.store_frameID_results(order, result_vector=solution_array, result_flag=result_flag) + self._modelling_parameters.store_frameID_results( + order, result_vector=solution_array, result_flag=result_flag + ) def _store_model_to_disk(self) -> NoReturn: # TODO: store information related with the GP parameters! return super()._store_model_to_disk() - def interpolate_spectrum_to_wavelength(self, og_lambda, og_spectra, og_err, new_wavelengths, order): + def interpolate_spectrum_to_wavelength( + self, og_lambda, og_spectra, og_err, new_wavelengths, order + ): """Interpolate the order of this spectrum to a given wavelength, using a GP. If the GP fit is yet to be done, then it is done beforehand. @@ -189,14 +211,18 @@ def interpolate_spectrum_to_wavelength(self, og_lambda, og_spectra, og_err, new_ t0 = time.time() t1 = time.time() - self.generate_model_from_order(og_lambda, og_spectra, og_err, new_wavelengths, order) + self.generate_model_from_order( + og_lambda, og_spectra, og_err, new_wavelengths, order + ) kern_type = self._internal_configs["GP_KERNEL"] try: fit_results = self._modelling_parameters.get_fit_results_from_frameID(order) except custom_exceptions.NoConvergenceError as exc: - logger.critical("Can't interpolate wavelengths from order that has not achieved convergence") + logger.critical( + "Can't interpolate wavelengths from order that has not achieved convergence" + ) raise exc param_names = self._modelling_parameters.get_enabled_params() @@ -217,7 +243,9 @@ def interpolate_spectrum_to_wavelength(self, og_lambda, og_spectra, og_err, new_ return mu, std def _launch_GP_fit(self, og_lambda, og_spectra, og_err, new_wavelengths, order): - initial_params, bounds = self._modelling_parameters.generate_optimizer_inputs(order, rv_units=None) + initial_params, bounds = self._modelling_parameters.generate_optimizer_inputs( + order, rv_units=None + ) param_names = self._modelling_parameters.get_enabled_params() result_flag = SUCCESS diff --git a/uv.lock b/uv.lock index 1857dad..c597444 100644 --- a/uv.lock +++ b/uv.lock @@ -77,7 +77,7 @@ wheels = [ [[package]] name = "astra-spectra" -version = "1.2.6" +version = "1.2.7" source = { editable = "." } dependencies = [ { name = "astropy" },