diff options
| author | Christian Kolset <christian.kolset@gmail.com> | 2025-11-16 15:27:32 -0700 |
|---|---|---|
| committer | Christian Kolset <christian.kolset@gmail.com> | 2025-11-16 15:27:32 -0700 |
| commit | 8f81b92b4942879edb32f9420514b7e402144a44 (patch) | |
| tree | 4a89503069455c9b570ebacca370f7aa67865e50 | |
| parent | 1873d28735e18866f7e0167045611440e4cbe20e (diff) | |
Completed Calibration of Spectrometer part
| -rw-r--r-- | tutorials/module_4/Spectroscopy problem.md | 83 | ||||
| -rw-r--r-- | tutorials/module_4/spectroscopy_problem/Spectroscopy.ipynb | 181 | ||||
| -rw-r--r-- | tutorials/module_4/spectroscopy_problem/spectroscopy.py | 110 |
3 files changed, 300 insertions, 74 deletions
diff --git a/tutorials/module_4/Spectroscopy problem.md b/tutorials/module_4/Spectroscopy problem.md index 18feeb9..5252d46 100644 --- a/tutorials/module_4/Spectroscopy problem.md +++ b/tutorials/module_4/Spectroscopy problem.md @@ -101,12 +101,34 @@ $$ $$ - Solve for C1 C2 and C3, use $I=195.54$ + +Well-known/defined Hg lines: + +Wavelengths [nm] +365.015 +435.833 +546.074 +576.960 +579.006 + + + +| Pixel peak | Wavelength (Literature) | +| ---------- | ----------------------- | +| | | + + + + + + Calibrate the y-axis +- Using Tungsten (W) - Find the response function of the spectrometer $$ R(\lambda)=\frac{I_{measured}}{I_{true}} $$ -- $I_{measured}$ => Intensity array from Mercury (Hg) +- $I_{measured}$ => Intensity array from Tungsten (W) - $I_{true}$ => Data from tungsten lamp using planks law of radiation $$ I_{\lambda,\Omega}(T)= \epsilon (\frac{2hc^2}{\lambda^5}\frac{1}{e^{hc/kT}-1}) @@ -124,7 +146,7 @@ I_{meas}^{plasma}(\lambda) = \frac{I_{meas}^{W}(\lambda)}{I_{true}^{W}(\lambda)} $$ - Measure the densities of an excited state of oxygen using $I(\lambda)$ +Measure the densities of an excited state of oxygen using $I(\lambda)$ $$ I(\lambda)=\frac{1}{4\pi}hvAnl\phi(\lambda-\lambda_0) $$ @@ -133,59 +155,4 @@ $$ n_{1,2}=\frac{\int_{\lambda_0-\Delta\lambda}^{\lambda_0+\Delta\lambda} I(\lambda)d\lambda}{\frac{1}{4\pi}hv_{1,2}A_{1,2}l}\tag{1} $$ - - -||Intensity||||Air <br>Wavelength (Å)||||Spectrum||||Reference|| -|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---| -|| | | | | | | | | | | | | | | | -|||||||||||||||| -|20|2026.860|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|400 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#2052.828)|2052.828|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|20|2224.711|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|10|2252.786|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|60|2260.294|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|400 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#2262.223)|2262.223|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|10|2263.634|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|1000 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#2536.517),c|2536.517|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|25|2652.039|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|40|2653.679|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|400 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#2847.675)|2847.675|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|30|2916.250|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|25|2947.074|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|250 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#2967.280)|2967.280|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|70|3021.498|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|90|3125.668|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|80|3131.548|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|80|3131.839|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|12|3208.169|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|10|3532.594|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|10|3605.762|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|600 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#3650.153)|3650.153|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|70|3654.836|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|50|3663.279|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|1000 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#3983.931),c|3983.931|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|400 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#4046.563)|4046.563|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|60|4339.223|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|100|4347.494|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|1000 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#4358.328)|4358.328|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|12 c|5128.442|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|15|5204.768|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|80 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#5425.253)|5425.253|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|500 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#5460.735)|5460.735|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|200 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#5677.105)|5677.105|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|50|5769.598|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|60|5790.663|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|12|5871.279|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|20 c|5888.939|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|15|6146.435|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|250 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#6149.475),c|6149.475|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|25|7081.90|Hg I|[F54](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#F54)| -|6|7346.508|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|250 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable4.htm#7944.555)|7944.555|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|6 h|9520.198|Hg II|[SR01](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#SR01)| -|200 [P](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable3.htm#10139.76)|10139.76|Hg I|[BAL50](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#BAL50)| -|50|13570.21|Hg I|[H53](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#H53)| -|40|13673.51|Hg I|[H53](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#H53)| -|50|15295.82|Hg I|[H53](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#H53)| -|50|17072.79|Hg I|[H53](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#H53)| -|25|23253.07|Hg I|[PBT55](https://www.physics.nist.gov/PhysRefData/Handbook/Tables/mercurytable7.htm#PBT55)|
\ No newline at end of file +|
\ No newline at end of file diff --git a/tutorials/module_4/spectroscopy_problem/Spectroscopy.ipynb b/tutorials/module_4/spectroscopy_problem/Spectroscopy.ipynb new file mode 100644 index 0000000..ca53b38 --- /dev/null +++ b/tutorials/module_4/spectroscopy_problem/Spectroscopy.ipynb @@ -0,0 +1,181 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "a8f4e56a-7e65-47d9-95a1-5ea952f52a5d", + "metadata": {}, + "source": [ + "Problem:\n", + "\n", + "- Import xls data into Python\n", + "- Plot the Intensity [a.u.] vs pixels\n", + "- Interpolate and convert x-axis from pixels to nm (true wavelength) using Hg lamp data (using data in file: **Lampa_Calibrare_Mercur.xlsx**)\n", + "- Find response function of the spectrometer using the tungsten lamp data from file: \"**Calibrare Intensitate Oxigen.xlsx**)\": $R=\\frac{I_{measured}}{I_{true}}$ (where True is computed by Planck's law of radiation (see notes in the pptx above)\n", + "- Convert y-axis from Intensity [a.u.] into Intensity in [W/(cm^2*sr*nm)] by dividing the measured Oxygen spectrum with the response function: $I_{oxygen, true}=\\frac{I_{oxygen, measured}}{R}$\n", + "- Once the spectra is in real units: compute the density of one of the oxygen lines by integrating underneath one of the peaks (see equation from Slide 39 - bottom). We will give all of the constants that are in this equation (see the \"I**ntensity_Calibration_Oxygen_Discharge_Solution.xlsx**\")\n" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "bec0d1d0-9534-462b-91ce-409b1b0a99f8", + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import pandas as pd\n", + "import matplotlib.pyplot as plt\n", + "\n", + "\n", + "Hg_data = 'Lampa_Calibrare_Mercur.xlsx'\n", + "#Ar_data = 'Spectru Descarcare Argon.xlsx'\n", + "Ox_data = 'Calibrare Intensitate Oxigen.xlsx'\n", + "\n", + "df_Hg = pd.read_excel(Hg_data, header=11, engine='openpyxl')\n", + "df_Ox = pd.read_excel(Ox_data, header=3, engine='openpyxl')\n", + "#df_Ar = pd.read_excel(Ar_data, header=14, engine='openpyxl')\n", + "\n", + "\n", + "# Trims mercury data to the first 2000 pixels\n", + "df_Hg = df_Hg[df_Hg['Pixels']<2000]" + ] + }, + { + "cell_type": "markdown", + "id": "a7909b57-25bb-4e96-9331-f4acbf6a5886", + "metadata": {}, + "source": [ + "## Plot Intensity vs Pixels for mercury" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "ec2a3865-26f7-4218-b003-11407f83452b", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAskAAAHUCAYAAADIlbU1AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlHJYcgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAd2lJREFUeJzt3Xl8E3X+P/DXJE3Sgza9aEulXIrIUQHR5XIXEGg5SmXdXXTh24VVQUVFBHS9La7CiquwP/BARFDAZd0V8MCtLcoh23JYqMohopSbUoTeRzJN5vdHmiEntJBmJunr+XjwoJl8MvOZdybJO5+85zOCJEkSiIiIiIhIplG6A0REREREasMkmYiIiIjIBZNkIiIiIiIXTJKJiIiIiFwwSSYiIiIicsEkmYiIiIjIBZNkIiIiIiIXTJKJiIiIiFwwSSYiIiIicsEkmYhahZUrV0IQBPlfSEgI2rdvjz//+c84deqU3G7KlCno1KlTi/Vjy5YtEAQBW7ZsaXKfv/nmm2Zvp7a2FtnZ2U3ajr906tQJU6ZMkW+fPn0a2dnZKCoqapHtvfDCC+jRowesVqu8TBAEPPTQQx7b/+c//2nyc+OqrKwM0dHR2LBhwxX2lojUJkTpDhAR+dOKFStwww03oK6uDtu2bcP8+fOxdetWfP/994iIiMCzzz6LRx55ROluXrXa2lrMnTsXADB06FBlO9No/fr1iIqKkm+fPn0ac+fORadOndCnTx+fbuv06dNYsGABVq5cCY2m5ceDYmJi8Oijj+Kxxx7DmDFjoNfrW3ybRNSyOJJMRK1Kr169MGDAAAwbNgzPP/88Hn/8cRQXF8sjgNdeey369u2rbCeDVN++fXHttdf6ZVv/+Mc/EB0djTvuuMMv2wOA+++/H0ePHsV//vMfv22TiFoOk2QiatUGDBgAADh27BgA93KLtWvXQhAELFmyxOlxzz//PLRaLfLy8uRl33zzDTIzMxEbG4vQ0FD07dsXH374oU/7O2XKFLRp0wY//fQTxowZgzZt2iAlJQWzZ8+GyWQCABw9ehRt27YFAMydO1cuMXEsdTh8+DAmTpyIhIQEGAwGdO/eHa+//rrTtuylIf/85z/x9NNPIzk5GVFRURgxYgQOHTrk1Hbv3r3IyMiQ15ecnIyxY8fi5MmTchvHcostW7bglltuAQD8+c9/lvuYnZ2NVatWQRAEFBQUuO3/Cy+8AJ1Oh9OnT3uNkdlsxvLlyzFx4kSfjCIvW7YM119/PQwGA3r06IEPPvjAY1lOYmIiRo4cibfeeuuqt0lEymOSTESt2k8//QQAclLp6q677sL999+P2bNny7XBX331FV588UU89dRTGDlyJABg8+bNGDx4MMrLy/HWW2/h448/Rp8+fXDnnXdi5cqVPu2zKIrIzMzE8OHD8fHHH+Puu+/GwoUL8fLLLwMA2rVrh5ycHADAPffcg4KCAhQUFODZZ58FABw4cAC33HIL9u3bh1dffRWfffYZxo4dixkzZsglGo6eeuopHDt2DO+88w7efvttHD58GOPGjYPFYgEA1NTUYOTIkTh79ixef/115OXlYdGiRejQoQOqqqo87sNNN92EFStWAACeeeYZuY/33nsv7rzzTiQlJbkl7Q0NDVi6dCl++9vfIjk52Wt8du7cifPnz2PYsGEe75ckCQ0NDW7/HGuX7d5++21MmzYNN954I9atW4dnnnkGc+fO9Vq3PHToUPzvf/9DeXm51/4RUYCQiIhagRUrVkgApB07dkiiKEpVVVXSZ599JrVt21aKjIyUSkpKJEmSpMmTJ0sdO3Z0emx9fb3Ut29fqXPnztKBAwekxMREaciQIVJDQ4Pc5oYbbpD69u0riaLo9NiMjAypXbt2ksVikSRJkjZv3iwBkDZv3tzkPu/evVteNnnyZAmA9OGHHzq1HTNmjNStWzf59rlz5yQA0vPPP++23vT0dKl9+/ZSRUWF0/KHHnpICg0NlS5cuODU1zFjxji1+/DDDyUAUkFBgSRJkvTNN99IAKQNGzZccn86duwoTZ48Wb69e/duCYC0YsUKt7bPP/+8pNfrpbNnz8rL/vWvf0kApK1bt15yOy+//LIEQH5OHQG47D/7c2OxWKSkpCSpf//+Tus4duyYpNPp3I4TSZKkvLw8CYD03//+95J9JCL140gyEbUqAwYMgE6nQ2RkJDIyMpCUlIT//ve/SExM9PoYg8GADz/8EOfPn8dNN90ESZLwz3/+E1qtFoBtNPqHH37ApEmTAMBpdHLMmDE4c+aMW3mCneRhVPNyBEHAuHHjnJbdeOONcsnIpdTX1+PLL7/Eb3/7W4SHh7v1tb6+Hjt27HB6TGZmptu2gIslKtdddx1iYmLwl7/8BW+99RYOHDhw2X5czgMPPADAVupgt2TJEqSmpuI3v/nNJR97+vRpCIKA+Ph4j/dPmDABu3fvdvtnH4m3O3ToEEpKSjBhwgSn5R06dMDgwYM9rjshIQEAnGZMIaLAxCSZiFqV999/H7t378bevXtx+vRpfPfdd14THkfXXXcdfv3rX6O+vh6TJk1Cu3bt5PvOnj0LAJgzZw50Op3Tv+nTpwMAfvnlF4/r3bp1q9tjjh49esm+hIeHIzQ01GmZwWBAfX39Zffj/PnzaGhowOLFi922O2bMGI99jYuLc9sWANTV1QEAjEYjtm7dij59+uCpp55Cz549kZycjOeffx6iKF62T54kJibizjvvxNKlS2GxWPDdd9/h66+/9jp9m6O6ujrodDr5S4yrtm3b4uabb3b716VLF6d258+fl/viqX+e2J8Xe2yIKHBxCjgialW6d++Om2++udmPe+edd7Bx40b86le/wpIlS3DnnXeif//+ACCPWD755JNeZ1Po1q2bx+X9+vXD7t27nZZdqt72asXExECr1SIrKwsPPvigxzadO3du9npTU1Oxdu1aSJKE7777DitXrsQLL7yAsLAwPPHEE1fU10ceeQSrVq3Cxx9/jJycHERHR8uj9ZcSHx8Ps9mMmpoaREREXNG2gYtfDuxfghyVlJR4fMyFCxfkPhBRYGOSTER0Gd9//z1mzJiBP/3pT1i2bBkGDRqEO++8E3v37kVMTAy6deuGrl274ttvv8W8efOate7IyMgrStovx3W01y48PBzDhg3D3r17ceONN/p8Pl9BENC7d28sXLgQK1euxJ49e5rdR7t+/fph0KBBePnll7Fv3z5MmzatSUnvDTfcAAD4+eef5dKQK9GtWzckJSXhww8/xKxZs+Tlx48fR35+vscvM0eOHAEA9OjR44q3S0TqwCSZiOgSampqMGHCBHTu3BlvvPEG9Ho9PvzwQ9x0003485//LM+vvHTpUowePRrp6emYMmUKrrnmGly4cAEHDx7Enj178O9//9uv/Y6MjETHjh3x8ccfY/jw4YiNjUV8fDw6deqEf/zjH7j11lvx61//Gg888AA6deqEqqoq/PTTT/j000/x1VdfNWtbn332Gd544w2MHz8eXbp0gSRJWLduHcrLy+XZPzy59tprERYWhjVr1qB79+5o06YNkpOTnZLPRx55BHfeeScEQZBLVy7HfvGUHTt2XFWSrNFoMHfuXNx33334/e9/j7vvvhvl5eWYO3cu2rVr53F6uR07diAuLg6pqalXvF0iUgfWJBMRXcL999+P48eP49///rc8itmlSxe88847+Pjjj7Fo0SIAwLBhw7Br1y5ER0dj5syZGDFiBB544AFs2rQJI0aMUKTvy5cvR3h4ODIzM3HLLbcgOzsbgG2Uc8+ePejVqxeeeeYZpKWl4Z577sF//vMfDB8+vNnb6dq1K6Kjo7FgwQJkZmbiD3/4A/bs2YOVK1di6tSpXh8XHh6Od999F+fPn0daWhpuueUWvP32205txo8fD4PBgPT0dHTt2rVJ/UlJScGvf/1rfPzxx83eF1fTpk3D22+/jW+//Ra//e1vMXfuXDzxxBPo27cvoqOjndpKkoRPPvkEEydOhCAIV71tIlKWIEmSpHQniIiIPPn000+RmZmJjRs3yicWNsVHH32EO++8E8eOHcM111zj0z6Vl5fj+uuvx/jx452S+i+//BJpaWnYv3+/XPJBRIGLSTIREanOgQMHcOzYMTzyyCOIiIjAnj17mjU6K0kSBg0ahH79+rldLbE5SkpK8NJLL2HYsGGIi4vDsWPHsHDhQvzwww/45ptv0LNnT7ntsGHDcN111zlNW0dEgYs1yUREpDrTp0/H//73P9x000147733ml2+IAgCli1bhk8++QRWq/WKL09tMBhw9OhRTJ8+HRcuXEB4eDgGDBiAt956yylBLisrw5AhQ5pcN01E6seRZCIiIiIiFzxxj4iIiIjIBZNkIiIiIiIXTJKJiIiIiFzwxD0fslqtOH36NCIjIzlHJhEREZEKSZKEqqoqJCcnX/KkXibJPnT69GmkpKQo3Q0iIiIiuowTJ06gffv2Xu9nkuxDkZGRAGxBj4qKavHtiaKI3NxcpKWlQafTtfj2AgXj4h1j4xnj4h1j4xnj4h1j4xnj4p2/Y1NZWYmUlBQ5b/OGSbIP2UssoqKi/JYkh4eHIyoqii84B4yLd4yNZ4yLd4yNZ4yLd4yNZ4yLd0rF5nKlsTxxj4iIiIjIBZNkIiIiIiIXTJKJiIiIiFwwSSYiIiIicsEkmYiIiIjIBZNkIiIiIiIXTJKJiIiIiFwwSSYiIiIicsEkmYiIiIjIBZNkIiIiIiIXTJKJiIiIiFwwSSYiIiIicsEkmYiIiIjIBZNkIiIiFXjlix/wxEffQZIkpbtCRGCSTEREpLiSinq8vvlnrN19AqVVJqW7Q0RgkkxERKS4w6VV8t+ixapgT4jIjkkyERGRwqwOFRYNFpZbEKkBk2QiIiKFOdYhN1g5kkykBkySiYiIFOY4dixyJJlIFRRNkrdt24Zx48YhOTkZgiBgw4YNbm0OHjyIzMxMGI1GREZGYsCAATh+/Lh8v8lkwsMPP4z4+HhEREQgMzMTJ0+edFpHWVkZsrKyYDQaYTQakZWVhfLycqc2x48fx7hx4xAREYH4+HjMmDEDZrO5JXabiIjIGcstiFRH0SS5pqYGvXv3xpIlSzze//PPP+PWW2/FDTfcgC1btuDbb7/Fs88+i9DQULnNzJkzsX79eqxduxbbt29HdXU1MjIyYLFY5DYTJ05EUVERcnJykJOTg6KiImRlZcn3WywWjB07FjU1Ndi+fTvWrl2Ljz76CLNnz265nSciImokOWTJIsstiFQhRMmNjx49GqNHj/Z6/9NPP40xY8ZgwYIF8rIuXbrIf1dUVGD58uVYtWoVRowYAQBYvXo1UlJSsGnTJqSnp+PgwYPIycnBjh070L9/fwDAsmXLMHDgQBw6dAjdunVDbm4uDhw4gBMnTiA5ORkA8Oqrr2LKlCl46aWXEBUV1RK7T0REBACQOJJMpDqKJsmXYrVasXHjRjz++ONIT0/H3r170blzZzz55JMYP348AKCwsBCiKCItLU1+XHJyMnr16oX8/Hykp6ejoKAARqNRTpABYMCAATAajcjPz0e3bt1QUFCAXr16yQkyAKSnp8NkMqGwsBDDhg3z2EeTyQST6eJ8lpWVlQAAURQhiqIvw+GRfRv+2FYgYVy8Y2w8Y1y8Y2w883VcxIYG+e96szmg481jxjPGxTt/x6ap21FtklxaWorq6mr87W9/w4svvoiXX34ZOTk5uOOOO7B582YMGTIEJSUl0Ov1iImJcXpsYmIiSkpKAAAlJSVISEhwW39CQoJTm8TERKf7Y2JioNfr5TaezJ8/H3PnznVbnpubi/Dw8Gbv85XKy8vz27YCCePiHWPjGePiHWPjma/i8v0FAYAWAJC/YxfKfgj80WQeM54xLt75Kza1tbVNaqfaJNnaWJN1++2349FHHwUA9OnTB/n5+XjrrbcwZMgQr4+VJAmCIMi3Hf++mjaunnzyScyaNUu+XVlZiZSUFKSlpfmlREMUReTl5WHkyJHQ6XQtvr1Awbh4x9h4xrh4x9h45uu46A+W4p1DRQCAPjf1w/Ab3Ad3AgWPGc8YF+/8HRv7L/+Xo9okOT4+HiEhIejRo4fT8u7du2P79u0AgKSkJJjNZpSVlTmNJpeWlmLQoEFym7Nnz7qt/9y5c/LocVJSEnbu3Ol0f1lZGURRdBthdmQwGGAwGNyW63Q6v74A/L29QMG4eMfYeMa4eMfYeOaruGi02os3BE1QxJrHjGeMi3f+ik1Tt6HaeZL1ej1uueUWHDp0yGn5jz/+iI4dOwIA+vXrB51O5zQ8f+bMGezbt09OkgcOHIiKigrs2rVLbrNz505UVFQ4tdm3bx/OnDkjt8nNzYXBYEC/fv1abB+JiIhccZ5kInVQdCS5uroaP/30k3y7uLgYRUVFiI2NRYcOHfDYY4/hzjvvxG9+8xsMGzYMOTk5+PTTT7FlyxYAgNFoxD333IPZs2cjLi4OsbGxmDNnDlJTU+XZLrp3745Ro0Zh6tSpWLp0KQBg2rRpyMjIQLdu3QAAaWlp6NGjB7KysvDKK6/gwoULmDNnDqZOncqZLYiIqMU5zW7BKeCIVEHRkeRvvvkGffv2Rd++fQEAs2bNQt++ffHcc88BAH7729/irbfewoIFC5Camop33nkHH330EW699VZ5HQsXLsT48eMxYcIEDB48GOHh4fj000+hdfjpas2aNUhNTUVaWhrS0tJw4403YtWqVfL9Wq0WGzduRGhoKAYPHowJEyZg/Pjx+Pvf/+6nSBARUevmME8yR5KJVEHRkeShQ4c6Xa/ek7vvvht333231/tDQ0OxePFiLF682Gub2NhYrF69+pLb6dChAz777LNLd5iIiKgFOH0UMkcmUgXV1iQTERG1FsyLidSHSTIREZHCHEeSJabMRKrAJJmIiEhhjonxZaoQichPmCQTEREpjIkxkfowSSYiIlIYz9sjUh8myURERApznOmJo8pE6sAkmYiIiIjIBZNkIiIihXF2CyL1YZJMRESkMM5uQaQ+TJKJiIgUxsSYSH2YJBMREakI82UidWCSTEREpDCnkWQOKxOpApNkIiIihTEtJlIfJslEREQKc5onWcF+ENFFTJKJiIgUxmoLIvVhkkxERKQ0JsZEqsMkmYiISGHO8yQzYyZSAybJRERECnO+4h4RqQGTZCIiIoUxMSZSHybJRERECnMaSWbGTKQKTJKJiIgU5lSTrGA/iOgiJslEREQK4+gxkfowSSYiIlKY8zzJzJiJ1IBJMhERkdKYGBOpDpNkIiIihTFFJlIfJslEREQqwkFlInVgkkxERKQw54uJMEsmUgMmyURERArjyXpE6sMkmYiISGHOs1so1g0icsAkmYiISGHO5RZEpAZMkomIiBTGxJhIfZgkExERKcyxJpnlFkTqwCSZiIhIRTi7BZE6MEkmIiJSGEePidRH0SR527ZtGDduHJKTkyEIAjZs2OC17X333QdBELBo0SKn5SaTCQ8//DDi4+MRERGBzMxMnDx50qlNWVkZsrKyYDQaYTQakZWVhfLycqc2x48fx7hx4xAREYH4+HjMmDEDZrPZR3tKRETknePoMRNmInVQNEmuqalB7969sWTJkku227BhA3bu3Ink5GS3+2bOnIn169dj7dq12L59O6qrq5GRkQGLxSK3mThxIoqKipCTk4OcnBwUFRUhKytLvt9isWDs2LGoqanB9u3bsXbtWnz00UeYPXu273aWiIjICybGROoTouTGR48ejdGjR1+yzalTp/DQQw/hiy++wNixY53uq6iowPLly7Fq1SqMGDECALB69WqkpKRg06ZNSE9Px8GDB5GTk4MdO3agf//+AIBly5Zh4MCBOHToELp164bc3FwcOHAAJ06ckBPxV199FVOmTMFLL72EqKioFth7IiIiG+bIROqjaJJ8OVarFVlZWXjsscfQs2dPt/sLCwshiiLS0tLkZcnJyejVqxfy8/ORnp6OgoICGI1GOUEGgAEDBsBoNCI/Px/dunVDQUEBevXq5TRSnZ6eDpPJhMLCQgwbNsxj/0wmE0wmk3y7srISACCKIkRRvOr9vxz7NvyxrUDCuHjH2HjGuHjH2Hjm67g0NFic/g7kePOY8Yxx8c7fsWnqdlSdJL/88ssICQnBjBkzPN5fUlICvV6PmJgYp+WJiYkoKSmR2yQkJLg9NiEhwalNYmKi0/0xMTHQ6/VyG0/mz5+PuXPnui3Pzc1FeHj4pXfOh/Ly8vy2rUDCuHjH2HjGuHjH2Hjmq7j8cEoAoAUAHPrxED6v/cEn61USjxnPGBfv/BWb2traJrVTbZJcWFiIf/zjH9izZw8EQWjWYyVJcnqMp8dfSRtXTz75JGbNmiXfrqysREpKCtLS0vxSoiGKIvLy8jBy5EjodLoW316gYFy8Y2w8Y1y8Y2w883VcTmwrxmfHDwMArr++G8YM7XLV61QKjxnPGBfv/B0b+y//l6PaJPnrr79GaWkpOnToIC+zWCyYPXs2Fi1ahKNHjyIpKQlmsxllZWVOo8mlpaUYNGgQACApKQlnz551W/+5c+fk0eOkpCTs3LnT6f6ysjKIoug2wuzIYDDAYDC4LdfpdH59Afh7e4GCcfGOsfGMcfGOsfHMV3ERNBfPo9doNEERax4znjEu3vkrNk3dhmrnSc7KysJ3332HoqIi+V9ycjIee+wxfPHFFwCAfv36QafTOQ3PnzlzBvv27ZOT5IEDB6KiogK7du2S2+zcuRMVFRVObfbt24czZ87IbXJzc2EwGNCvXz9/7C4REREAznRBpBaKjiRXV1fjp59+km8XFxejqKgIsbGx6NChA+Li4pza63Q6JCUloVu3bgAAo9GIe+65B7Nnz0ZcXBxiY2MxZ84cpKamyrNddO/eHaNGjcLUqVOxdOlSAMC0adOQkZEhryctLQ09evRAVlYWXnnlFVy4cAFz5szB1KlTObMFERG1OImZMZHqKDqS/M0336Bv377o27cvAGDWrFno27cvnnvuuSavY+HChRg/fjwmTJiAwYMHIzw8HJ9++im0Wq3cZs2aNUhNTUVaWhrS0tJw4403YtWqVfL9Wq0WGzduRGhoKAYPHowJEyZg/Pjx+Pvf/+67nSUiIvLCMUfmZamJ1EHRkeShQ4c269vz0aNH3ZaFhoZi8eLFWLx4sdfHxcbGYvXq1Zdcd4cOHfDZZ581uS9ERES+4vhJyEFlInVQbU0yERFRa8HEmEh9mCQTEREpzLHEgvkykTowSSYiIlKYxHoLItVhkkxERKQwpsVE6sMkmYiISGkSyy2I1IZJMhERkcJYbUGkPkySiYiIFMbEmEh9mCQTEREpzHl2C2bMRGrAJJmIiEhhTlfcY45MpApMkomIiIiIXDBJJiIiUpjk5W8iUg6TZCIiIoWx3IJIfZgkExERKYwn6xGpD5NkIiIipTmOJDNhJlIFJslEREQKk7zeICKlMEkmIiJSmMTLUhOpDpNkIiIihfFkPSL1YZJMRESkMKcp4JgxE6kCk2QiIiIVYY5MpA5MkomIiBTGxJhIfZgkExERqQjzZSJ1YJJMRESkMMe5kTmqTKQOTJKJiIgUxsSYSH2YJBMREakIr7hHpA5MkomIiFSEo8pE6sAkmYiIiIjIBZNkIiIihfECIkTqwySZiIhIRZgwE6kDk2QiIiKFMS0mUh8myURERApzHDxmwkykDkySiYiIVITVFkTqwCSZiIhIYZwbmUh9mCQTERGpCBNmInVgkkxERKQwp5pk5shEqsAkmYiIiIjIhaJJ8rZt2zBu3DgkJydDEARs2LBBvk8URfzlL39BamoqIiIikJycjD/96U84ffq00zpMJhMefvhhxMfHIyIiApmZmTh58qRTm7KyMmRlZcFoNMJoNCIrKwvl5eVObY4fP45x48YhIiIC8fHxmDFjBsxmc0vtOhERkUzy8jcRKUfRJLmmpga9e/fGkiVL3O6rra3Fnj178Oyzz2LPnj1Yt24dfvzxR2RmZjq1mzlzJtavX4+1a9di+/btqK6uRkZGBiwWi9xm4sSJKCoqQk5ODnJyclBUVISsrCz5fovFgrFjx6Kmpgbbt2/H2rVr8dFHH2H27Nktt/NERESNWG5BpD4hSm589OjRGD16tMf7jEYj8vLynJYtXrwYv/rVr3D8+HF06NABFRUVWL58OVatWoURI0YAAFavXo2UlBRs2rQJ6enpOHjwIHJycrBjxw70798fALBs2TIMHDgQhw4dQrdu3ZCbm4sDBw7gxIkTSE5OBgC8+uqrmDJlCl566SVERUW1YBSIiIiISG0UTZKbq6KiAoIgIDo6GgBQWFgIURSRlpYmt0lOTkavXr2Qn5+P9PR0FBQUwGg0ygkyAAwYMABGoxH5+fno1q0bCgoK0KtXLzlBBoD09HSYTCYUFhZi2LBhHvtjMplgMpnk25WVlQBspSKiKPpy1z2yb8Mf2wokjIt3jI1njIt3jI1nvo6L1Wpx+juQ481jxjPGxTt/x6ap2wmYJLm+vh5PPPEEJk6cKI/slpSUQK/XIyYmxqltYmIiSkpK5DYJCQlu60tISHBqk5iY6HR/TEwM9Hq93MaT+fPnY+7cuW7Lc3NzER4e3rwdvAquI+5kw7h4x9h4xrh4x9h45qu4HD+ugb0C8vjxE/j882M+Wa+SeMx4xrh456/Y1NbWNqldQCTJoijirrvugtVqxRtvvHHZ9pIkQRAE+bbj31fTxtWTTz6JWbNmybcrKyuRkpKCtLQ0v5RoiKKIvLw8jBw5EjqdrsW3FygYF+8YG88YF+8YG898HZf/bdiPgtJTAICUlBSMGdPzqtepFB4znjEu3vk7NvZf/i9H9UmyKIqYMGECiouL8dVXXzkln0lJSTCbzSgrK3MaTS4tLcWgQYPkNmfPnnVb77lz5+TR46SkJOzcudPp/rKyMoii6DbC7MhgMMBgMLgt1+l0fn0B+Ht7gYJx8Y6x8Yxx8Y6x8cxXcdFoLp5HLwiaoIg1jxnPGBfv/BWbpm5D1fMk2xPkw4cPY9OmTYiLi3O6v1+/ftDpdE7D82fOnMG+ffvkJHngwIGoqKjArl275DY7d+5ERUWFU5t9+/bhzJkzcpvc3FwYDAb069evJXeRiIjIeXYLTgJHpAqKjiRXV1fjp59+km8XFxejqKgIsbGxSE5Oxu9//3vs2bMHn332GSwWi1wfHBsbC71eD6PRiHvuuQezZ89GXFwcYmNjMWfOHKSmpsqzXXTv3h2jRo3C1KlTsXTpUgDAtGnTkJGRgW7dugEA0tLS0KNHD2RlZeGVV17BhQsXMGfOHEydOpUzWxARUYtjYkykPoomyd98843TzBH2+t7JkycjOzsbn3zyCQCgT58+To/bvHkzhg4dCgBYuHAhQkJCMGHCBNTV1WH48OFYuXIltFqt3H7NmjWYMWOGPAtGZmam09zMWq0WGzduxPTp0zF48GCEhYVh4sSJ+Pvf/94Su01EROQV50kmUgdFk+ShQ4dCusS7waXuswsNDcXixYuxePFir21iY2OxevXqS66nQ4cO+Oyzzy67PSIiIl9zLrcgIjVQdU0yEREREZESmCQTEREpzHH0mOUWROrAJJmIiEhFeBIfkTowSSYiIlIYR4+J1IdJMhERkcKcRo+ZMBOpApNkIiIiFWGOTKQOTJKJiIiUxsyYSHWYJBMREalIU64RQEQtj0kyERGRwiQvfxORcpgkExERERG5YJJMRESkMMcSC1ZbEKkDk2QiIiKFsdyCSH2YJBMRERERuWCSTEREpDDHEgvObkGkDkySiYiIVIQpMpE6MEkmIiJSGBNjIvVhkkxERKQmzJiJVIFJMhERkcJYh0ykPkySiYiIFMYUmUh9mCQTERGpiMSUmUgVQprSaNasWc1e8TPPPIPY2NhmP46IiKjVYV5MpDpNSpIXLVqEgQMHQq/XN2ml27dvx0MPPcQkmYiIqJlYnkykDk1KkgFg/fr1SEhIaFLbyMjIK+4QERFRa8MSCyL1aVJN8ooVK2A0Gpu80qVLlyIxMfGKO0VERNRacSSZSB2aNJI8efLkZq104sSJV9QZIiKi1oiJMZH6cHYLIiIihTkmySy9IFIHnyXJkydPxm233ear1RERERERKabJJ+5dzjXXXAONhgPTREREzcXRYyL18VmSPG/ePF+tioiIqNVifTKROnDol4iISGFMjInUp9kjyXffffcl73/33XevuDNEREStHfNlInVodpJcVlbmdFsURezbtw/l5eU8cY+IiOgKMDEmUp9mJ8nr1693W2a1WjF9+nR06dLFJ50iIiJqTZymgGPGTKQKPqlJ1mg0ePTRR7Fw4UJfrI6IiIiISFE+O3Hv559/RkNDQ7Mes23bNowbNw7JyckQBAEbNmxwul+SJGRnZyM5ORlhYWEYOnQo9u/f79TGZDLh4YcfRnx8PCIiIpCZmYmTJ086tSkrK0NWVhaMRiOMRiOysrJQXl7u1Ob48eMYN24cIiIiEB8fjxkzZsBsNjdrf4iIiK6M5OVvIlJKs8stZs2a5XRbkiScOXMGGzdubPblq2tqatC7d2/8+c9/xu9+9zu3+xcsWIDXXnsNK1euxPXXX48XX3wRI0eOxKFDhxAZGQkAmDlzJj799FOsXbsWcXFxmD17NjIyMlBYWAitVgvAdpnskydPIicnBwAwbdo0ZGVl4dNPPwUAWCwWjB07Fm3btsX27dtx/vx5TJ48GZIkYfHixc0NEREREREFuGYnyXv37nW6rdFo0LZtW7z66quXnfnC1ejRozF69GiP90mShEWLFuHpp5/GHXfcAQB47733kJiYiA8++AD33XcfKioqsHz5cqxatQojRowAAKxevRopKSnYtGkT0tPTcfDgQeTk5GDHjh3o378/AGDZsmUYOHAgDh06hG7duiE3NxcHDhzAiRMnkJycDAB49dVXMWXKFLz00kuIiopq1n4RERE1B+uQidSn2Uny5s2bW6IfboqLi1FSUoK0tDR5mcFgwJAhQ5Cfn4/77rsPhYWFEEXRqU1ycjJ69eqF/Px8pKeno6CgAEajUU6QAWDAgAEwGo3Iz89Ht27dUFBQgF69eskJMgCkp6fDZDKhsLAQw4YN89hHk8kEk8kk366srARgm/FDFEWfxcIb+zb8sa1Awrh4x9h4xrh4x9h45uu4WCWr/LfFag3oePOY8Yxx8c7fsWnqdnx2xT1fKykpAQAkJiY6LU9MTMSxY8fkNnq9HjExMW5t7I8vKSlBQkKC2/oTEhKc2rhuJyYmBnq9Xm7jyfz58zF37ly35bm5uQgPD7/cLvpMXl6e37YVSBgX7xgbzxgX7xgbz3wVl7NnNbCfJlRaWorPP//cJ+tVEo8ZzxgX7/wVm9ra2ia181mS/NRTT6GkpMTnFxMRBMHptiRJbstcubbx1P5K2rh68sknnWq0KysrkZKSgrS0NL+UaIiiiLy8PIwcORI6na7FtxcoGBfvGBvPGBfvGBvPfB2Xdef3AOW/AADatk3AmDE3XfU6lcJjxjPGxTt/x8b+y//l+CxJPnXqFE6cOOGr1SEpKQmAbZS3Xbt28vLS0lJ51DcpKQlmsxllZWVOo8mlpaUYNGiQ3Obs2bNu6z937pzTenbu3Ol0f1lZGURRdBthdmQwGGAwGNyW63Q6v74A/L29QMG4eMfYeMa4eMfYeOaruGgcBmQ0Gk1QxJrHjGeMi3f+ik1Tt+GzKeDee+89fPXVV75aHTp37oykpCSnoXez2YytW7fKCXC/fv2g0+mc2pw5cwb79u2T2wwcOBAVFRXYtWuX3Gbnzp2oqKhwarNv3z6cOXNGbpObmwuDwYB+/fr5bJ+IiIg8cZoAjmfxEamCojXJ1dXV+Omnn+TbxcXFKCoqQmxsLDp06ICZM2di3rx56Nq1K7p27Yp58+YhPDwcEydOBAAYjUbcc889mD17NuLi4hAbG4s5c+YgNTVVnu2ie/fuGDVqFKZOnYqlS5cCsE0Bl5GRgW7dugEA0tLS0KNHD2RlZeGVV17BhQsXMGfOHEydOpUzWxARERG1QleUJNfU1GDr1q04fvy42wU3ZsyY0eT1fPPNN04zR9jreydPnoyVK1fi8ccfR11dHaZPn46ysjL0798fubm58hzJALBw4UKEhIRgwoQJqKurw/Dhw7Fy5Up5jmQAWLNmDWbMmCHPgpGZmYklS5bI92u1WmzcuBHTp0/H4MGDERYWhokTJ+Lvf/978wJDRER0BZwuS61cN4jIwRXNkzxmzBjU1taipqYGsbGx+OWXXxAeHo6EhIRmJclDhw695M9KgiAgOzsb2dnZXtuEhoZi8eLFl7zoR2xsLFavXn3JvnTo0AGfffbZZftMRERERMGv2TXJjz76KMaNG4cLFy4gLCwMO3bswLFjx9CvXz+OvBIREV0Bjh4TqU+zk+SioiLMnj0bWq0WWq0WJpMJKSkpWLBgAZ566qmW6CMREVFQc/xVleftEalDs5NknU4nzx2cmJiI48ePA7CdRGf/m4iIiIgokDW7Jrlv37745ptvcP3112PYsGF47rnn8Msvv2DVqlVITU1tiT4SERG1GhxIJlKHZo8kz5s3T764x1//+lfExcXhgQceQGlpKd5++22fd5CIiIiIyN+aPZJ88803y3+3bds2KK4vT0REpCSnKeBYlEykCj674h4RERERUbBoUpJ80003oaysrMkrvfXWW3Hq1Kkr7hQREVFrIrESmUh1mlRuUVRUhG+//RaxsbFNWmlRURFMJtNVdYyIiKi1YIUFkfo0uSZ5+PDhTa6Tsk8RR0REREQUiJqUJBcXFzd7xe3bt2/2Y4iIiFoj5xP3lOsHEV3UpCS5Y8eOLd0PIiIiIiLV4OwWRERECnM8cY8n8RGpA5NkIiIihbHEgkh9mCQTERGpCBNmInVgkkxERKQw5sVE6tPsJHnKlCnYtm1bS/SFiIio1eNIMpE6NDtJrqqqQlpaGrp27Yp58+bxynpERERXi4kxkeo0O0n+6KOPcOrUKTz00EP497//jU6dOmH06NH4z3/+A1EUW6KPRERERER+dUU1yXFxcXjkkUewd+9e7Nq1C9dddx2ysrKQnJyMRx99FIcPH/Z1P4mIiIIWp4AjUp+rOnHvzJkzyM3NRW5uLrRaLcaMGYP9+/ejR48eWLhwoa/6SEREFNRYh0ykPs1OkkVRxEcffYSMjAx07NgR//73v/Hoo4/izJkzeO+995Cbm4tVq1bhhRdeaIn+EhERBTUmzETq0KTLUjtq164drFYr/vjHP2LXrl3o06ePW5v09HRER0f7oHtERETBj3kxkfo0O0leuHAh/vCHPyA0NNRrm5iYGBQXF19Vx4iIiFojJsxE6tDscovNmzd7nMWipqYGd999t086RURE1JpIrLEgUp1mJ8nvvfce6urq3JbX1dXh/fff90mniIiIiIiU1ORyi8rKSkiSBEmSUFVV5VRuYbFY8PnnnyMhIaFFOklERBTMJK83iEgpTU6So6OjIQgCBEHA9ddf73a/IAiYO3euTztHRETUGrDagkh9mpwkb968GZIk4bbbbsNHH32E2NhY+T69Xo+OHTsiOTm5RTpJRETUWvBiIkTq0OQkeciQIQCA4uJidOjQAYIgtFiniIiIWhOmxUTq06Qk+bvvvkOvXr2g0WhQUVGB77//3mvbG2+80WedIyIiam1YekGkDk1Kkvv06YOSkhIkJCSgT58+EATB43Q1giDAYrH4vJNERERBjZkxkeo0KUkuLi5G27Zt5b+JiIioZTBdJlKHJiXJHTt29Pg3ERERXT0mxkTqc0UXE9m4caN8+/HHH0d0dDQGDRqEY8eO+bRzDQ0NeOaZZ9C5c2eEhYWhS5cueOGFF2C1WuU2kiQhOzsbycnJCAsLw9ChQ7F//36n9ZhMJjz88MOIj49HREQEMjMzcfLkSac2ZWVlyMrKgtFohNFoRFZWFsrLy326P0RERJ6w2oJIfZqdJM+bNw9hYWEAgIKCAixZsgQLFixAfHw8Hn30UZ927uWXX8Zbb72FJUuW4ODBg1iwYAFeeeUVLF68WG6zYMECvPbaa1iyZAl2796NpKQkjBw5ElVVVXKbmTNnYv369Vi7di22b9+O6upqZGRkONVPT5w4EUVFRcjJyUFOTg6KioqQlZXl0/0hIiK6HF6imkgdmjwFnN2JEydw3XXXAQA2bNiA3//+95g2bRoGDx6MoUOH+rRzBQUFuP322zF27FgAQKdOnfDPf/4T33zzDQDbG8miRYvw9NNP44477gBgG+lOTEzEBx98gPvuuw8VFRVYvnw5Vq1ahREjRgAAVq9ejZSUFGzatAnp6ek4ePAgcnJysGPHDvTv3x8AsGzZMgwcOBCHDh1Ct27dfLpfREREjjg3MpH6NDtJbtOmDc6fP48OHTogNzdXHj0ODQ1FXV2dTzt366234q233sKPP/6I66+/Ht9++y22b9+ORYsWAbCdRFhSUoK0tDT5MQaDAUOGDEF+fj7uu+8+FBYWQhRFpzbJycno1asX8vPzkZ6ejoKCAhiNRjlBBoABAwbAaDQiPz/fa5JsMplgMpnk25WVlQAAURQhiqIvQ+GRfRv+2FYgYVy8Y2w8Y1y8Y2w883VcHEePrZIU0PHmMeMZ4+Kdv2PT1O00O0keOXIk7r33XvTt2xc//vijPMq7f/9+dOrUqbmru6S//OUvqKiowA033ACtVguLxYKXXnoJf/zjHwEAJSUlAIDExESnxyUmJsr10SUlJdDr9YiJiXFrY3+8fXo7VwkJCXIbT+bPn+/xUty5ubkIDw9vxp5enby8PL9tK5AwLt4xNp4xLt4xNp75Ki4VFVoAtot0lZeV4/PPP/fJepXEY8YzxsU7f8Wmtra2Se2anSS//vrreOaZZ3DixAl89NFHiIuLAwAUFhbKyauv/Otf/8Lq1avxwQcfoGfPnigqKsLMmTORnJyMyZMny+1cr/4nSdJlrwjo2sZT+8ut58knn8SsWbPk25WVlUhJSUFaWhqioqIuu39XSxRF5OXlYeTIkdDpdC2+vUDBuHjH2HjGuHjH2Hjm67gsPVoA1NjOpTFGR2PMmP6XeYR68ZjxjHHxzt+xsf/yfznNTpKjo6OxZMkSt+WeRlSv1mOPPYYnnngCd911FwAgNTUVx44dw/z58zF58mQkJSUBsI0Et2vXTn5caWmpPLqclJQEs9mMsrIyp9Hk0tJSDBo0SG5z9uxZt+2fO3fObZTakcFggMFgcFuu0+n8+gLw9/YCBePiHWPjGePiHWPjme/i4jxoEwyx5jHjGePinb9i09RtNDtJBoDy8nLs2rULpaWlTtOxCYLg0xkhamtrodE4T8Ch1WrlbXbu3BlJSUnIy8tD3759AQBmsxlbt27Fyy+/DADo168fdDod8vLyMGHCBADAmTNnsG/fPixYsAAAMHDgQFRUVGDXrl341a9+BQDYuXMnKioq5ESaiIiopUhe/iYi5TQ7Sf70008xadIk1NTUIDIy0q1kwZdJ8rhx4/DSSy+hQ4cO6NmzJ/bu3YvXXnsNd999t7y9mTNnYt68eejatSu6du2KefPmITw8HBMnTgQAGI1G3HPPPZg9ezbi4uIQGxuLOXPmIDU1VZ7tonv37hg1ahSmTp2KpUuXAgCmTZuGjIwMzmxBRERE1Ao1O0mePXs27r77bjkZbUmLFy/Gs88+i+nTp6O0tBTJycm477778Nxzz8ltHn/8cdTV1WH69OkoKytD//79kZubi8jISLnNwoULERISggkTJqCurg7Dhw/HypUrodVq5TZr1qzBjBkz5FkwMjMzPZaVEBER+RrnRiZSn2YnyadOncKMGTP8MntDZGQkFi1aJE/55okgCMjOzkZ2drbXNqGhoVi8eLHTRUhcxcbGYvXq1VfRWyIiIh9gwkykCs2+4l56erp8MQ8iIiIiomDU7JHksWPH4rHHHsOBAweQmprqdoZgZmamzzpHRETU2nAcmUgdmp0kT506FQDwwgsvuN0nCAIsFsvV94qIiKgVYYUFkfo0O0l2nPKNiIiIrp7kMH7MhJlIHZpdk+yovr7eV/0gIiIiIlKNZifJFosFf/3rX3HNNdegTZs2OHLkCADg2WefxfLly33eQSIiomDnOHossSqZSBWanSS/9NJLWLlyJRYsWAC9Xi8vT01NxTvvvOPTzhERERERKaHZSfL777+Pt99+G5MmTXK6GMeNN96IH374waedIyIiag04dkykPs1Okk+dOoXrrrvObbnVaoUoij7pFBEReZb9yX489/E+pbtBLYgn7hGpQ7OT5J49e+Lrr792W/7vf/8bffv29UmniIjIXZ3ZgpX5R/F+wTGcLKtVujvkQ7wsNZH6NHsKuOeffx5ZWVk4deoUrFYr1q1bh0OHDuH999/HZ5991hJ9JCIiAFaHRKqkoh7tY8IV7A35kmOKzHyZSB2aPZI8btw4/Otf/8Lnn38OQRDw3HPP4eDBg/j0008xcuTIlugjERHBOUkurTIp2BPyOSbGRKrT7JFkAEhPT0d6erqv+0JERJdgdUikTA28umkwkbz8TUTKafZIcpcuXXD+/Hm35eXl5ejSpYtPOkVERB4wewparEkmUp9mJ8lHjx6FxeI+gmEymXDq1CmfdIqIiNw5llsIEBTsCfmac00yE2YiNWhyucUnn3wi//3FF1/AaDTKty0WC7788kt06tTJp50jIqKLHJNkXpUtuDAvJlKfJifJ48ePBwAIgoDJkyc73afT6dCpUye8+uqrPu0cERFd5FiTzKQquPBLD5H6NDlJtlqtAIDOnTtj9+7diI+Pb7FOERGRO8ef4ZkkBxc+n0Tq0+zZLYqLi1uiH0REdBmOI8lWZlVBhU8nkfpc0RRwX375Jb788kuUlpbKI8x27777rk86RkREzhx/kmdSFbz43BKpQ7OT5Llz5+KFF17AzTffjHbt2kEQeIY1EZE/ONUks4aViKhFNTtJfuutt7By5UpkZWW1RH+IiMgLq0OWbGWOHFQkzlxCpDrNnifZbDZj0KBBLdEXIiK6BIk1yUGLzyaR+jQ7Sb733nvxwQcftERfiIjoEqyc3SJoSZzej0h1ml1uUV9fj7fffhubNm3CjTfeCJ1O53T/a6+95rPOERHRRc5JMjOpYMISCyL1aXaS/N1336FPnz4AgH379jndx5P4iIhajmMaxZrk4MLvPETq0+wkefPmzS3RDyIiugzH0WPWJAcXycvfRKScZtckExGRMnhZ6uDF55NIfZo8knzHHXc0qd26deuuuDNEROSdlSPJQYz15kRq0+Qk2Wg0tmQ/iIjoMhwvcMo8Krjw+SRSnyYnyStWrGjJfhAR0WVYecGJoMWaZCL1YU0yEVEA4uwWwYUlFkTqwySZiChAsCY5eElebxCRUpgkExEFCM5uQUTkP0ySiYgCBK+4F7z4dBKpj+qT5FOnTuH//u//EBcXh/DwcPTp0weFhYXy/ZIkITs7G8nJyQgLC8PQoUOxf/9+p3WYTCY8/PDDiI+PR0REBDIzM3Hy5EmnNmVlZcjKyoLRaITRaERWVhbKy8v9sYtERE3ifDERBTtCPic5nZRJRGqg6iS5rKwMgwcPhk6nw3//+18cOHAAr776KqKjo+U2CxYswGuvvYYlS5Zg9+7dSEpKwsiRI1FVVSW3mTlzJtavX4+1a9di+/btqK6uRkZGBiwWi9xm4sSJKCoqQk5ODnJyclBUVISsrCx/7i4R0SU5JsasSQ4ufDaJ1KfZl6X2p5dffhkpKSlO08916tRJ/luSJCxatAhPP/20fLGT9957D4mJifjggw9w3333oaKiAsuXL8eqVaswYsQIAMDq1auRkpKCTZs2IT09HQcPHkROTg527NiB/v37AwCWLVuGgQMH4tChQ+jWrZvH/plMJphMJvl2ZWUlAEAURYii6NNYeGLfhj+2FUgYF+8YG88CJS6i2CD/bbFY+D6jIF/HxfE7j9UqBXS8ecx4xrh45+/YNHU7gqTiwrYePXogPT0dJ0+exNatW3HNNddg+vTpmDp1KgDgyJEjuPbaa7Fnzx707dtXftztt9+O6OhovPfee/jqq68wfPhwXLhwATExMXKb3r17Y/z48Zg7dy7effddzJo1y628Ijo6GgsXLsSf//xnj/3Lzs7G3Llz3ZZ/8MEHCA8P90EEiIgu+qkCWHzANraRfo0VYzpYL/MIChSP79LCZBEAAPGhEp7ta7nMI4joStXW1mLixImoqKhAVFSU13aqHkk+cuQI3nzzTcyaNQtPPfUUdu3ahRkzZsBgMOBPf/oTSkpKAACJiYlOj0tMTMSxY8cAACUlJdDr9U4Jsr2N/fElJSVISEhw235CQoLcxpMnn3wSs2bNkm9XVlYiJSUFaWlplwy6r4iiiLy8PIwcORI6na7FtxcoGBfvGBvPAiUuO45cAA58AwDoct21GDOia4tvM1Bi42++jstThV/C1FgCGB4ejjFjfn3V61QKjxnPGBfv/B0b+y//l6PqJNlqteLmm2/GvHnzAAB9+/bF/v378eabb+JPf/qT3E4QBKfHSZLktsyVaxtP7S+3HoPBAIPB4LZcp9P59QXg7+0FCsbFO8bGM7XHRaPVXvxbo+H7jAr4Ki6OP+kKghAUseYx4xnj4p2/YtPUbaj6xL127dqhR48eTsu6d++O48ePAwCSkpIAwG20t7S0VB5dTkpKgtlsRllZ2SXbnD171m37586dcxulJiJSipWzWwQtiXNgE6mOqpPkwYMH49ChQ07LfvzxR3Ts2BEA0LlzZyQlJSEvL0++32w2Y+vWrRg0aBAAoF+/ftDpdE5tzpw5g3379sltBg4ciIqKCuzatUtus3PnTlRUVMhtiIiUxtktgpfE+S2IVEfV5RaPPvooBg0ahHnz5mHChAnYtWsX3n77bbz99tsAbD9JzZw5E/PmzUPXrl3RtWtXzJs3D+Hh4Zg4cSIAwGg04p577sHs2bMRFxeH2NhYzJkzB6mpqfJsF927d8eoUaMwdepULF26FAAwbdo0ZGRkeJ3ZgojI35zOs2ZOFVT4nYdIfVSdJN9yyy1Yv349nnzySbzwwgvo3LkzFi1ahEmTJsltHn/8cdTV1WH69OkoKytD//79kZubi8jISLnNwoULERISggkTJqCurg7Dhw/HypUroXWo71uzZg1mzJiBtLQ0AEBmZiaWLFniv50lIroMiSPJQUty+pvPLZEaqDpJBoCMjAxkZGR4vV8QBGRnZyM7O9trm9DQUCxevBiLFy/22iY2NharV6++mq4SEbUo1iQHMT6fRKqj6ppkIiK6yMqTu1oFPrdE6sAkmYgoQDiPJDOTCiYssSBSHybJREQBwvHEPRVfLJWuAKeAI1IfJslERAHC+cQ95fpBvsenk0h9mCQTEQUIp5pkplVBhb8MEKkPk2QiogDB2S2CF59OIvVhkkxEFCCsrEkOWnw6idSHSTIRUYDgyV2tA78AEakDk2QiogDBKeCCE5NiInVikkxEFCA4u0Vwcs2R+dQSqQOTZCKiAOFck6xgR8in+FQSqROTZCKiAOFck8zUKli4PpdKPrUNFiu2HCpFZb2oXCeIVIJJMhFRgGBNcnBS0zO5dNsRTFmxG//3zk6lu0KkOCbJREQBwsqa5KCkpu87HxWeBAB8d7JC4Z4QKY9JMhFRgHCqSVawH9SylLyaosXhGBv/+v/w3Mf7FOsLkdKYJBMRBQjH1InlFsFDTZcYtzj8RFF0ohzvFxxTsDdEymKSTEQUICRecS8ouU0Bp+BTy8OK6CImyUREAcJq5RRw1LIsLHYnkjFJJiIKEM4n7jGZCRZqupgIjyuii5gkExEFCOcp4BTsCPmUmmqSmSQTXcQkmYgoQDhfTES5fpBvqakm2VO5BUswqLVikkxEFCAcRxx54l7wUNMz6SkfbrBa/d8RIhVgkkxEFCBYkxyc1PSFx8qRZCIZk2QiogDBmuTg5P5UKvfkevry1cCDjVopJslERAHCqSZZuW6Qj6loINnjcdVgUVEHifyISTIRUYBwnieZiUvQUNGJe4KHZaxJptaKSTIRUYBgTXJwUtMUcBrBPU1mTTK1VkySiYgChPPsFgp2hFqUok+th6FklltQa8UkmYgoQHAkOTip6an0XG6hog4S+RGTZCKiACFxdoug5PpUKllvLngst2BNMrVOTJKJiAKEldNbBCU1nYTpIUfmSDK1WkySiYgCBMstgpOanklPJ+6xJplaKybJRERNIEkSDpyuRI2pQbE+OF9MhIlLsHB9KtX2zHJ2C2qtmCQTETXB1h/PYcz/+xrjX/+fcp1gtUVQUtMUcJ7wCxm1VgGVJM+fPx+CIGDmzJnyMkmSkJ2djeTkZISFhWHo0KHYv3+/0+NMJhMefvhhxMfHIyIiApmZmTh58qRTm7KyMmRlZcFoNMJoNCIrKwvl5eV+2CsiCgSffnsGAHC4tFqxPvCy1EFKRRcT8YRJMrVWAZMk7969G2+//TZuvPFGp+ULFizAa6+9hiVLlmD37t1ISkrCyJEjUVVVJbeZOXMm1q9fj7Vr12L79u2orq5GRkYGLBaL3GbixIkoKipCTk4OcnJyUFRUhKysLL/tHxGpW7heq3QXnBJjNZ3sRVdH7c8kv5BRaxUQSXJ1dTUmTZqEZcuWISYmRl4uSRIWLVqEp59+GnfccQd69eqF9957D7W1tfjggw8AABUVFVi+fDleffVVjBgxAn379sXq1avx/fffY9OmTQCAgwcPIicnB++88w4GDhyIgQMHYtmyZfjss89w6NAhRfaZiNRFHUkya5KDkVtNssqeW9YkU2sVonQHmuLBBx/E2LFjMWLECLz44ovy8uLiYpSUlCAtLU1eZjAYMGTIEOTn5+O+++5DYWEhRFF0apOcnIxevXohPz8f6enpKCgogNFoRP/+/eU2AwYMgNFoRH5+Prp16+axXyaTCSaTSb5dWVkJABBFEaIo+mz/vbFvwx/bCiSMi3eMjWdNiYtjjqxU/BosF+ertVolvs8oyJdxMXtYh1Lx9pSgi2JDs/rDY8YzxsU7f8emqdtRfZK8du1a7NmzB7t373a7r6SkBACQmJjotDwxMRHHjh2T2+j1eqcRaHsb++NLSkqQkJDgtv6EhAS5jSfz58/H3Llz3Zbn5uYiPDz8MnvmO3l5eX7bViBhXLxjbDy7VFyOnRIA2DLlTz77HCEK/A539KgG9h8AKyoq8fnnn/tt2zxmPPNFXMpMgOPHsSiKfn1uHZnNWrhed69g505c+KH5o8k8ZjxjXLzzV2xqa2ub1E7VSfKJEyfwyCOPIDc3F6GhoV7buV4hSJIkj1cNulQbT+0vt54nn3wSs2bNkm9XVlYiJSUFaWlpiIqKuuT2fUEUReTl5WHkyJHQ6XQtvr1Awbh4p8bYVNWLWJB7GBmpSejfOVaRPjQlLiX/O4pPjv8IALhtZBraGPz/9rnr04NAyQkAQJvISIwZM6jFt6nGY0YNfBmX0+V1yN7ztXw7RKfDmDHpV9vFK5L97WbUNDiPst18yy349XXxTV4HjxnPGBfv/B0b+y//l6PqJLmwsBClpaXo16+fvMxisWDbtm1YsmSJXC9cUlKCdu3ayW1KS0vl0eWkpCSYzWaUlZU5jSaXlpZi0KBBcpuzZ8+6bf/cuXNuo9SODAYDDAaD23KdTufXF4C/txcoGBfv1BSbN744jLW7T2Lt7pM4+rexivblUnHRaLQOf4coEz+nL+0C32dUwBdx0Ya4//SrVKw9DQxpNNor6g+PGc8YF+/8FZumbkPVJ+4NHz4c33//PYqKiuR/N998MyZNmoSioiJ06dIFSUlJTsPzZrMZW7dulRPgfv36QafTObU5c+YM9u3bJ7cZOHAgKioqsGvXLrnNzp07UVFRIbchopZxoqxpP3spzXEuW4tCJ1bxinvByf3EPWX64Q2PNWqtVD2SHBkZiV69ejkti4iIQFxcnLx85syZmDdvHrp27YquXbti3rx5CA8Px8SJEwEARqMR99xzD2bPno24uDjExsZizpw5SE1NxYgRIwAA3bt3x6hRozB16lQsXboUADBt2jRkZGR4PWmPiHwjTKf8rBFN4ZigNlit3hu2ZB+snN0iGLk+l2qb3UKhw51IcapOkpvi8ccfR11dHaZPn46ysjL0798fubm5iIyMlNssXLgQISEhmDBhAurq6jB8+HCsXLkSWu3FD+c1a9ZgxowZ8iwYmZmZWLJkid/3h6i1CVPB1GpN4TgNllJJg2Mypa40iq6G6wxrantulfrlhEhpAZckb9myxem2IAjIzs5Gdna218eEhoZi8eLFWLx4sdc2sbGxWL16tY96SURNZQgJjCS5wXIxUVBqJNkxWWHeEjxc5yFW268EahvZJvIXVdckE1HwM+guvg1ZVXzRAovVcY5iZfogsSY5KLmXWyjUEXhOiFX8siRqUUySiUhRgsOcrKKKix8brCoYSWZNclBSU5LsKSHmFfeotWKSTESKcpo1QsUfxmpIUK0stwhKrse9pGBVsqdfc/iFjForJslEpBqiRb0fxs4jyUySyXfsz6Wm8UcVJb8rekqImSRTa8UkmYhUo8Gi3nILx9E+pUa8Has8mLgED/vxFKKxfSQreaKcp5ksVFwFRdSimCQTkaIsFuWTz6ZwrENWqp+c3SI42Z9XbeNQsqIjyR4SYn4ho9aKSTIRKcqxdEFUcZKsjpFk5euiyffsI8chGsFtmb+x3ILoIibJRKQo0aHEQs3lFg0qGPF2TFZU/H2Cmsl+2IdoHZNkhfrCKeCIZEySiUhRzhfpUO+nsWPyoFy5heMt9caKmscql1tc/EhW4tmVJMljcq7mMiiilsQkmYgU5TRrhIpnt1BDuYXEkeSgZLW6l1soUeLg7ZjiFfeotWKSTESKcjwhTqmLdDSFYzLv6Sdpf1DDXM3ke64n7gHKlFt4O6Y4kkytFZNkIlKUU7mFmkeSVVAWwnmSg5P9cHKsSVbiS5C3ZJg5MrVWTJKJSFFOJ+4FyEiyp6uS+QPnSQ5O9uPJcSRZkX54OaZ4rFFrxSSZiBQVODXJjsk8R5LJd6wepoBTU00yk2RqrZgkE5GinEeS1fthrIaRZIvTiXvqjRU1j0UeSb74kazEIcZyCyJnTJKJSFGOH8yiiudJduynYiPJVo4kByNPI8lKzCjhbZs8cY9aKybJRKQoNVykoykaVDCzhGN4OJIcPOzPq9ap3ML//fD2+uMUcNRaMUkmIkWJDrW+oqprkpWvnbZwJDko2Z9XncPsFkpcTcTb1IYq/oGHqEUxSSYiRQXiSLJS8yQ7nbjHK+4FDauHeZKV+KXA2yb5qwW1VkySiUhRgTIFnOPsFkol81anE/cU6QK1gIs1ycpelprlFkTOmCQTkaKcT9xT74exGka8WZMcnOzfE5UeSfZ6xT0ea9RKMUkmIkU5lTGoeiRZBSfusSY5KNmPJ40ACI15siKXpfby8uOvFtRaMUkmIkU5lluoeSRZDSfuuSbn/Bk8ODhecc8+lqzEc+v1invMkqmVYpJMRIpSQxlDU6hhCjjXn71VHC5qBos8kixA0ziUrEhNMi9LTeSESTIRKarBaQq4wCi3UO5iIi63mbwEBbHB9sTqQjRyuYUiNcm84h6REybJRAEgd38JMhZ/jcNnq5Tuis85llio+bLUFqvyI97u5RaKdIN8zNz45dCg1UCwjyQrUZPsZZtq/oWHqCUxSSYKANNWFWLfqUo8/M+9SnfF59SQfDZFgwr66VqzzZHk4GBuHEnWh2jkmmQlnltOAUfkjEkyUQA5FJQjyYFSbqH8PMmu80gzdwkOpsYk2RCiuViTrMhIMqeAI3LEJJkogATjZ1WDCmaNaAo1jCS7xodX3QsOTiPJSk4B5/XEPT93hEglmCQTkWIkSVLFCXFN4VQWotC3FdeRZBWHi5rB5JAkX5zdQj3lFpwCjlorJslEAUQfElwvWdca2wYVl1uocSSZNcnBQU6StVqHmmT/98PbNnmcUWsVXJ+4REFOuHyTgOI6MhowI8kK9FOSJLf4MHcJDvZyC4POsdxCRRcT4XFGrRSTZCKVc/ywDLbPKtekz9RgUagnl+ZaFqJEkuzpCwRnHQgO9ing9A5TwCkyksxyCyInTJKJVM7sWIIQZJ9VruUDtWZ1JsmuSbESSbKnbTJ3CQ4m0Xbc22qS7UsVOMZ4xT0iJ6pOkufPn49bbrkFkZGRSEhIwPjx43Ho0CGnNpIkITs7G8nJyQgLC8PQoUOxf/9+pzYmkwkPP/ww4uPjERERgczMTJw8edKpTVlZGbKysmA0GmE0GpGVlYXy8vKW3kWiy6oX1Vune7Vca5DVmiS7juIqURbiaXo8jiQHB/txH2HQKjqS7O1wUvGkM0QtStVJ8tatW/Hggw9ix44dyMvLQ0NDA9LS0lBTUyO3WbBgAV577TUsWbIEu3fvRlJSEkaOHImqqovzyc6cORPr16/H2rVrsX37dlRXVyMjIwMWy8UP5IkTJ6KoqAg5OTnIyclBUVERsrKy/Lq/RJ44liAE23yloksmUKfSJNl1FFeJn589TY/HkeTgUGNuAABE6EPkkWQlXurevvxxJJlaqxClO3ApOTk5TrdXrFiBhIQEFBYW4je/+Q0kScKiRYvw9NNP44477gAAvPfee0hMTMQHH3yA++67DxUVFVi+fDlWrVqFESNGAABWr16NlJQUbNq0Cenp6Th48CBycnKwY8cO9O/fHwCwbNkyDBw4EIcOHUK3bt38u+NEDkyi80UsRIsVOq2qv982mftIcoNCPbk0VYwkO5zkqBFsCTJHkoNDjakxSTaEwH56rhKJqdjg+Vcr1iRTa6XqJNlVRUUFACA2NhYAUFxcjJKSEqSlpcltDAYDhgwZgvz8fNx3330oLCyEKIpObZKTk9GrVy/k5+cjPT0dBQUFMBqNcoIMAAMGDIDRaER+fr7XJNlkMsFkMsm3KysrAQCiKEIURd/tuBf2bfhjW4Ek2OJSXWdyul1VW4/IUN0VrUttsakzOfej1tSgSN8uF5d6k9npdoPF4vd+1jfGSqcVIEm2JMokihBFbYtuV23HjFr4Mi7V9bYk2aCFPJIsiv5/LdSbPW+vwWJtVl94zHjGuHjn79g0dTsBkyRLkoRZs2bh1ltvRa9evQAAJSUlAIDExESntomJiTh27JjcRq/XIyYmxq2N/fElJSVISEhw22ZCQoLcxpP58+dj7ty5bstzc3MRHh7ejL27Onl5eX7bViAJlricrAEcX6obc/IQpb+6daolNqdrAcd9+6WiCp9//rli/fEWl0oz4NjPEydP4fPPT/inU43O19v6IEjWxjILAV99+RWiDf7ZvlqOGbXxRVzKa7QABBTu+B9M9ba/v97+NYojrnrVzbLnFwGA+5eukpKSK3pd8pjxjHHxzl+xqa2tbVK7gEmSH3roIXz33XfYvn272332Ex3sJElyW+bKtY2n9pdbz5NPPolZs2bJtysrK5GSkoK0tDRERUVdcvu+IIoi8vLyMHLkSOh0VzayGIyCLS57jpcD3+2Sbw8eMhQpMVf2JUxtsdl/uhL4dod8W6MLxZgxQ/zej8vF5UxFPVC4Tb6dmNQOY8b09mcXUfxLDbD3fzDodTA3WGFpsGLosGFIjg5r0e2q7ZhRC1/GZc6uPAASRo+8DSuKd6JCNGHw4FvRM7nlP0cciUWngcP73JbHJyRgzJibmr4eHjMeMS7e+Ts29l/+LycgkuSHH34Yn3zyCbZt24b27dvLy5OSkgDYvuW2a9dOXl5aWiqPLiclJcFsNqOsrMxpNLm0tBSDBg2S25w9e9Ztu+fOnXMbpXZkMBhgMLgP4+h0Or++APy9vUARLHERrc5f1CyS5qr3SzWx0TiPWtWZLYr2y1tcBI3zT3PWxrZ+1RirEI0AUbAvCvFbP1RzzKjM1cbF1GCRrzxpjAiVL0ut0Wr9Hm+r4O1cB+GK+sJjxjPGxTt/xaap21D12T+SJOGhhx7CunXr8NVXX6Fz585O93fu3BlJSUlOw/Nmsxlbt26VE+B+/fpBp9M5tTlz5gz27dsntxk4cCAqKiqwa9fF0bqdO3eioqJCbkOklDrRcsnbgcx+pbFwvS0BrBUtqjwZzX2eZP/3ob7xeQ/TaeVEigJfreni6zlCf3EKOCVeBp6mGQQ4BRy1XqoeSX7wwQfxwQcf4OOPP0ZkZKRcH2w0GhEWFgZBEDBz5kzMmzcPXbt2RdeuXTFv3jyEh4dj4sSJctt77rkHs2fPRlxcHGJjYzFnzhykpqbKs110794do0aNwtSpU7F06VIAwLRp05CRkcGZLUhxrklxMM2bbE+So8N0qDVbYLFKMFusMIS07MlozeU69Z7F6v/nwD6Xbphei8rGE704NVfgq26c2SJUp0GI9uJlqZV4bj1NMwhwFhVqvVSdJL/55psAgKFDhzotX7FiBaZMmQIAePzxx1FXV4fp06ejrKwM/fv3R25uLiIjI+X2CxcuREhICCZMmIC6ujoMHz4cK1euhFZ78YN4zZo1mDFjhjwLRmZmJpYsWdKyO0jUBPVuSXIQjSQ3zlUeFabD6Yp6ALaSC7Ulya7JgxIja/YvS2F6LezjyJyZK/A5zpEMQP6VQInn1ttIMr+MUWul6iS5Kd9eBUFAdnY2srOzvbYJDQ3F4sWLsXjxYq9tYmNjsXr16ivpJlGLCuokuXEkOcIQghCNgAarhDrRgmhlu+XG8YIugDIjyfYLrYTptPJoI0f4Al+NyX61PdvHcYjW9uS6ziHuD97m/1biMuxEaqDqmmQicr8KXTDVJJsak2S9ViMnCfYLK6iJyeUiC0okDXKSrA+BRqPcaCP5lvOFRGyvBQAwK5Ekex1J9nNHiFSCSTKRyrkmxaYgqkm2n9WvD9EgMtSWJNjrbdXEdfReiSS5Vj5xTyOXW3AkOfDJSXLjyauGkMYk2cvV71qS6KWOiFfco9ZK1eUWROR+ol59Q/CMJNsTAZ1Wg6hQHYA6VNap72pUrl9MFEmSG5OpcH2IonWr5Fs1ZudyC3s9vuuvF/5gr0nOGtARMeE6mCxWLN16BCIPNGqlOJJMpHLBXZNs2xeDw0hylRpHkt1qkv2fNFyotV0aOzpcB21juUWDArXR5Fv2keQ29nILBUeS7TXJoToNZqV1w+Br4wEApiB6zyFqDibJRCrnVpNsDp7EyF53qQ/RICrMNrl7ZX0AjCQrUOZQVmNLkuMi9DDolEukyLeq5V8InMstXE8W9Ycah18rHPuiRH00kRowSSZSObd5koOw3EKvDYyR5DCdLZHxNp9sS7rQmCTHRhjkk7uU+EmefMvtxD0FR5Lt5wPYv7AaGo/3YDoPgqg5mCRTQDl4phIb9p5qVSeS2JPkmHDbB1cwlVvY6zHDDdrGmmSouiY5wmBLGpSYN/ZUuW0e6cQog6J1q+RbZbW24z0mXA/gYpKsxHNrf+1FhTqPJPM4o9aKJ+5RwLBaJYz5f19DkgCNRkBm72Slu3RV5n9+EIdLq/HGpJsQqvN+8Yx6OUnWo6xWDKor7l38UNbJMzWocSTZ/pN4VKgOv1Sbvc4n21IkScLx8zUAgI5xEXK5BWtFA1+Z/AtB4+itgolphf312DiSfHFUm8cZtU4cSaaAcaayHvYBvO9OlCval+b4cPcJ9J+3Cfk//SIvEy1WLN12BF/9UIqN35255OPtSXJ0EI4kVzn8vKvmmuSyxpPm4tsYAPh/SqzzNWbUmC0QBCAlNowjfEHEXkYTE2EbSbb/SqDE6/xclQkAECv3hccZtW5MkilgHGscSQOAYxdqFexJ8zz+0Xc4W2nCyvyj8rITDv2/3L5U1tkSyYTIUADBlSTbE+LI0BBV1ySXN/4kHh9pSx7sI8n+SpaPnbcdI8nGMBhCtHIixRP3Al9pla2MJi7C9gXMXlZl/2LmLzWmBpwqrwMAXNe2DQDn6eg4Jze1RkySKWDYEwUAOFNRp2BPms4xof2hpEr++6hDwn/yMkmyfeqvdtHBlyTb6zGNYTq5JrlKhSPJ9tE++0iyxSphyopd+PWCzSitrG/x7f9Uajt2OsdHAFC2bpV8x2KV5MQ0JTYMwMVRXPsx5y8/lVYDsB3j8qi27mKK4O1CI0TBjEkyBQzHJPl0ecsnJr5wsuxiMl9rvjhCWvzLxX05fokkWZIkuWYx2Wj7EPVlTbIkSVi35yS2/njOZ+tsDvuIekpMOCLlE/fUN5Jsf45uSIoCAJypqMeWQ+dwqrwO/9p9osW3b/+CdUNSJABlpwkj3zl+oRaiRYJeq0G7xtd3XOMXMXvpg78cajzGrk9sIy8LczhXQo1fXolaGpNkChjHL1wcfb1QY27SiOqbW37GoPlfouDn8y3ZNa8cyyp+qTbLcx4f/eXivlwqSS6rFeWf9uWRZB8mRl/sL8GsD7/F3St3yz/7+ktpVb08WtYxLhxRYfZyC3V9GFfWi/JzdEunGLf7f2wcgWtJP5yxJTDdGpPkCL16S1Oo6b48eBYAkNreKF8g5trGUod9pyr9+iVo99ELAIDeKdHyMp1WI/96cqYiMAYmiHyJSTIFjCPnapxu23+m9Ka0yoSXc37A6Yp6LNr04xVts/iXGmw5VHrJerzSynq88OkBFB4rc+/zL659rpXX69hP1wuG2P1wphKA7adYY+OJbd7aXom8A6UAbD/7bjnk39Hk/xSeBAD0uiYKEYaQiyPJKkv81uw4DotVwnUJbXBt2zZyMmO3/3RFi27f1GDBdyfLAQA9km0j2e1jbKOOJwKoNv9KfbG/BGt2HguqmtizlfV4v+AoXsuzvS/dcdM18n03JEUivo0BdaIF//2+xC/9kSQJO4ptAwm/6hzrdF9y45fzS32ZJwpWnAKOVOlUeR0Kj5Xhlk4xSIoKxc7iC/JPzu2MoThTUY/38o9iUv+OEC1WXJfQxmkaNdEKvLP9qHx7z/EyVJsa0MYQglPldSitrEfv9tHQNCY8pgYLvjxYCr1Wg6Hd2qK+wYpTZXW4/fXtqBeteGrMDZj2m2vd+ilJEmZ9+C22//QL1u4+ji2PDcW5KhPqRQuuS4h0K2P4+VwNwvQh2HO8zGV5Na5LaAOzxYo6swXnqkz48WwVlmz+CQBwU4cYRIfZ6gSPnq/BsfM1iArVITpcB0FwTtrs/ZIkQGr8WyMI0GgESJKE+gZb6covtQ34Yv/FD+FtP57DqF5JiDSEeFynN56SF4tVQlmtCKskQRBsI562fyKq6xtw6GwVFn9l27c/DewEAGgbaYAg2KZb23nkPJKjwxCm1yIuQt+s/lyKqcGCHUcuIDREg/5d4pzus0q2KenqaxpQXd+AalMDTpfXYfFXhwEADwy5FhqNgIRIg9Oo2pFzNagxNcgXg7hSFquE8zUmCBBQZ7ZAtNqOwTU7j6HGbEFSVCi6N5Z7dGqsTS44ch4HTleifWzYJZ83c4NVrmMOJF8fPof7VhUCsI1qTrg5BcDFY85Xx8XVkCQJ5gYrBAGoNVtgbrCi2tQA0WLFL1UmnCqvQ51oQa3ZglNldcj/+Rf87PCF/9dd4/GHfinybY1GwB03XYO3tx3BY//5FoXHytApPgKJUQYkRoUiMTIU+hANNBogRKNBiFZAiEaw/a0R5Pc0T/0EbDETLVaU1ZpRXivifLUZBT//ghMX6hCm0+Lmjs6/lvRNicZ3JyvwztdH0DbSgC7xEYj14WuyKUoq6lF4rAyixYrB18WjbaTBb9tuDQqPlWHF/4rxu37tMaxbgtLdURVBCqav5wqrrKyE0WhERUUFoqKiWnx7ZrMZz6zMwTldImrMVlitEvQhGpgbrHKyJZ9oJEkQLVaIFivCdFpIki0ZCdVpoQ/RwGqVIFolWK0SwvRalNeaoQ/RwBCitSVcgJx0ofE2YEsqLJKEuAgDJNh+KtdrNRAtVmgEAQadBlar7eILVklCvWiFVZLQYJEQqtMgJkKPkop6NFgl2N9y60SL08/Ieq1Gvizqr7vGY1i3BLzw2QGnWIRoBESGhqBetCJUp0FFrRlWuL+Jt4004JdqEyQJiNBrEarTQrRYUSda3E5McdwuAFwTHYZwva19jdkCk2hBjdkCSxNmOOgQG+42EhPfRo9rYsLx7WWms4uN0OOfUwegS9sI3PLSJnmmBQDyqKbjc+SJViNAAJo0v69GsH2QOr412P+63LtFSGN/mjqP8B19r8GrE3rLH7iZS7bju5POI7N6rQaGEA1C9VpYrbb9NIRooBEECAJs/yCgjSEEVkmCViPI/aw1N0Cn1UCrEWBusKKksh61jSPxycZQ6EI0qDVbUF0vou4Std4DusTig3sHQKMRMO/zg3h72xG0M4ai1myR55YN02kRYQhBuF4LCRKs1sa4C4Cmcf/srwNr46ZMDVbUixb5NXspL/8uFXfe0gEAUFErYtDfvpQvxoLGmBjDdCirNSO+jQGGEI287z+VViM6XIfocL0todIK0Go00GkEW7zsK5EACRIEQUCEXgvRIsHUYMGF8xcQExvjlBg5HgsXjw/3Y8ZbW3ufJUlCg1WCTqNBnWhBjakBotUKSbJ9WXZ8bO/2RtSJFhw7XwuNICAhyoAIfQgMOtu+huq0MIkWWylW4zGs12rk9yGLZHuPc/xfIwjQaTWwShIsVvuXSwligyT3o8FqhdVq2z9BEGBqsEKAhKo6E8xWAc2d5EQQgF7JRozr3Q733NrF7deJBosVj6wtwsbvLz09pCcaAQjRahBpCIFOq5GT9xpTAyyN8fB2wuejI67HIyO6Oi37+Vw1Mv7fdqcrf2oE28wX9k8CrSAgRKtBG0MIQrQCNADqamtgjIqEtjGR1wi2Y02rEaBt/NsQYnvOtRoBDRYJdaIFOq0AfYgG1SYLLFYrympEt18NEyINiArTQSsIiDBoEaLVIEynRbWpASEaAWF6LbSCgPM1ZpyrMqFjXDgEwfZFtMZkgUZjO75DtBpU14sI1WnlWGkE23tlqE4LQbC9ZgUIMFus8nL74IMk2dYpQUKYToswvRa/VJtharDKx51OI8D+0aKBhOMnT8EYl4CoMD1CdRo0WGyfvQ0WK0SLhAar7XNabJBw9HwNosN16JMSjYbGz+cGq+04DdVpIcD2HlMrWmBoPIYbrJJ8iXNL42evxSrJ7z0WqwRL47Fsaby988gF+f1nXO/kxtflxc8V22Mv7rPV8f/GWNjvb7BKOHGhFmH6ELSNNMjvMbYvc7YvcvUNFhjDdE7rEhssiKk7iaeyRkOn0zX7uG+upuZrTJJ9yJ9J8qGSKqQv2tai21BapCEEVY0XcQjVaTCieyKyM3siLkKP5duL8cm3p/FTaTX0IRqnxNEuJlyHu37VAe1jwvDcx/svm8waw3SoqhedPvDax4QhMlSHg41lD57oQzQY3SsJn39/BmJj8m+1AmaLFZGhIXho2HUY3j0Bf1q+C6cbRyD7pETj2YzuAARMX1OIs5UXT9LRCEBkqA5d2kagd/toTB92rTz929Yfz+Hl//6AH89WXdUFLSL0WsS1MaCdMRQP39YVq3YcxRf7z17x+i5FEIA2hhBEheoQGRpi+ztMhxHdE3HnLSlOCcJPpVV4dsN+7DleBkGwTz3VIt3ySqcVEBmqQ4RBizYGHboltkF2Zk9EN14RTbRY8eXBs+jbIQYb9p7CK18c8tnFRTSCbUQ7VKeBTqNBQpQBvVOicdctHdx+Bt9VfAGLvzqMohPlQV2bfG3bCJRWmuT3gkDQxhACQ4jtUusd4iLQxqBFmC4EMeE63NwpFgO7xMEYfulEwGqVkHfwLHYXX0BJZT1KK004W1WPs5X1aLDYEpyreW0Igu09LzZcj5gIPdJ6JGLqr7t4HInO/+kXvLXtCH4urcbpijq/vyY1AnB9YiREixVHfqnx+/bJP9KuseL1+0cxSQ5W/kySvz1Rjttf/x8A4Pbe7ZDWsx1Ei200ynZJUdsobmWdCEEQGr/p276h15stkABEGEJQa2pAg1VqHF26eMZ8ZKgOFqsV5gapcaTu4k+btr9t/0JDtLA0jhADtuTL0jiaZ2n8xmsf8dM2jixrBNs3ylpzA8pqRbSNNNhGtxtHJcJ0WkSF6RDfxoCzlfUQLVa0jbx4KV5XkiThbKUJFXUiwvValNfUY++Or3HX7aOh19sSm3NVJpTXmiFaJMRG6NEmNATHztfII222s8tDUVXfgB/PVkEfokGIRoOuiW2g12pwsKQSokVCjakB+hCNPHql12oQG6FHhCEEVfWiPHovwZZM2dcP2EYcKupERBi0TvvSYLGirFaU12cI0Vz2p0xJkmBqsP1kKkCARrA9MQIEp+fL/lzViRZYJUALK77+ahPGjxvj9kZkarBAkjxcFlpw/LPxGBDc726w2n52Fi1WRIXpEB2mg1XCJX8Cvpx60YKTZbVuybIkoTFJkBq/1Eg4X22GQadtHJmwHXcRhhCIDVY0WCUYQjSIDteja0IbnK8xyyeChuq0CNUCO77egtvHjkKbsOb9lFtea0ZVfQMkCagxN6DG1ABBsP0EbkuebSM3tlEq2/OiaRzlNITYRp/0IRqENo4E20afmxcve4lOeZ0ZMeF6/FJtQoNVgthgRX2DBVqNBuF6rS25stpGSC2No1a2X3EuHjeA7bmsEy22Y1GyYu/evbjpppsQotXCuWsXb9iXO97tuB/Oy23PYX2DBQIE6LS2WIXqbK8te2lIhCEE17Vtg2pzA4qOl8uvqQ6x4QCAc9UmVJsa5OPOJNrKSiTYfoHQaQV5lN7+PqjVwOFvQR59Q+Nj7J3VaWyPFxrbhTSOuFsl26ihWRSxK387xqQNR5tQg200Ua+FTqO54uO9uewji/bnsaHxObWXfDRYJLm/oY3/6kWL/EXVdQS7KepFC8prRbnEBLDFxL5NqyTBZG5AfsEO3PyrX0HQaJ1GQC1W+wim7fkK09v6ZLFCPnm3wSKhTWP5UqhOi94pRodzFkQc+6UWVSYRVivk0pY60YKo0BA0WG2fRw0WK6LDdQjXh+B8jclWciYIaLBaUS9aG38VtG2nTmyQR13R+N5ifx3rtLZfpvQhGlisUuOvCLZj2PG1ahItqDFZEBth26bZYpWPS/tnoElswMGDB3Fz71SYLLZ12T979Vrb/yEa2y8bIVoB+05VQq8VENo4Mm4/DrUawXZxIUAeOa4Xbc+HTqtBvWiRR8W1jaV2WsFWyiOP6Ns/lzUCYiP0GNAlDl8eLMWxCzXyZ8rFX+ts7/yaxnXIny+NbTSOnz2CrSTNVj5mhqXx2LQ0Hp8Wqy2nqDFZ5NciAAiwov7kQTx8l/tnU0toar7GmuQA1bltBFZO6Ydvv9mJB37XS04Gg01iVOhl2wiCgCRjKJKMtrZJkToc0Tt/QLeNNLjVsfVMNrqtKyZC71ar6q2tq8hQnfxGDgBajXNSb38zchWi1TS7xk4QBITqtPK0UZcT3fi/KIrQe7kCtj1xv9Qlsv0tVKfFdQmRPl+v6/EgiiL26S5OrdYc0eF6eZRZKWF6LTrEhaMDbMljSmMS6QuiKEI6LmFUz0S/fHh5EhWqw2+ub+u23F6brQRRFHE0zPbTv1Jx0WgE6DUC9H48Bz9Up0WS8dLvEaIo4twBCYOvjfN5bKJCdUhtf/n3YzUSRRGfVxzAmFvaNykuGTcm+6FXF429sZ1ft+dIFEV8/vlBxbbvTeCdyUEAbG8Ug6+NQ6dIdZy8QkRERBRMmCQTEREREblgkkxERERE5IJJMhERERGRCybJREREREQumCQTEREREblgkkxERERE5IJJMhERERGRCybJREREREQumCQTEREREblgkkxERERE5IJJMhERERGRCybJREREREQumCQTEREREblgkkxERERE5CJE6Q4EE0mSAACVlZV+2Z4oiqitrUVlZSV0Op1fthkIGBfvGBvPGBfvGBvPGBfvGBvPGBfv/B0be55mz9u8YZLsQ1VVVQCAlJQUhXtCRERERJdSVVUFo9Ho9X5BulwaTU1mtVpx+vRpREZGQhCEFt9eZWUlUlJScOLECURFRbX49gIF4+IdY+MZ4+IdY+MZ4+IdY+MZ4+Kdv2MjSRKqqqqQnJwMjcZ75TFHkn1Io9Ggffv2ft9uVFQUX3AeMC7eMTaeMS7eMTaeMS7eMTaeMS7e+TM2lxpBtuOJe0RERERELpgkExERERG5YJIcwAwGA55//nkYDAalu6IqjIt3jI1njIt3jI1njIt3jI1njIt3ao0NT9wjIiIiInLBkWQiIiIiIhdMkomIiIiIXDBJJiIiIiJywSSZiIiIiMgFk+QA9cYbb6Bz584IDQ1Fv3798PXXXyvdpRY1f/583HLLLYiMjERCQgLGjx+PQ4cOObWZMmUKBEFw+jdgwACnNiaTCQ8//DDi4+MRERGBzMxMnDx50p+74lPZ2dlu+5yUlCTfL0kSsrOzkZycjLCwMAwdOhT79+93WkewxcSuU6dObrERBAEPPvgggNZzvGzbtg3jxo1DcnIyBEHAhg0bnO731TFSVlaGrKwsGI1GGI1GZGVloby8vIX37upcKjaiKOIvf/kLUlNTERERgeTkZPzpT3/C6dOnndYxdOhQt+PorrvucmoTbLEBfPf6CbTYXC4unt5zBEHAK6+8IrcJxmOmKZ/RgfhewyQ5AP3rX//CzJkz8fTTT2Pv3r349a9/jdGjR+P48eNKd63FbN26FQ8++CB27NiBvLw8NDQ0IC0tDTU1NU7tRo0ahTNnzsj/Pv/8c6f7Z86cifXr12Pt2rXYvn07qqurkZGRAYvF4s/d8amePXs67fP3338v37dgwQK89tprWLJkCXbv3o2kpCSMHDkSVVVVcptgjAkA7N692ykueXl5AIA//OEPcpvWcLzU1NSgd+/eWLJkicf7fXWMTJw4EUVFRcjJyUFOTg6KioqQlZXV4vt3NS4Vm9raWuzZswfPPvss9uzZg3Xr1uHHH39EZmamW9upU6c6HUdLly51uj/YYmPni9dPoMXmcnFxjMeZM2fw7rvvQhAE/O53v3NqF2zHTFM+owPyvUaigPOrX/1Kuv/++52W3XDDDdITTzyhUI/8r7S0VAIgbd26VV42efJk6fbbb/f6mPLyckmn00lr166Vl506dUrSaDRSTk5OS3a3xTz//PNS7969Pd5ntVqlpKQk6W9/+5u8rL6+XjIajdJbb70lSVJwxsSbRx55RLr22mslq9UqSVLrPF4ASOvXr5dv++oYOXDggARA2rFjh9ymoKBAAiD98MMPLbxXvuEaG0927dolAZCOHTsmLxsyZIj0yCOPeH1MsMbGF6+fQI9NU46Z22+/XbrtttuclrWGY8b1MzpQ32s4khxgzGYzCgsLkZaW5rQ8LS0N+fn5CvXK/yoqKgAAsbGxTsu3bNmChIQEXH/99Zg6dSpKS0vl+woLCyGKolPskpOT0atXr4CO3eHDh5GcnIzOnTvjrrvuwpEjRwAAxcXFKCkpcdpfg8GAIUOGyPsbrDFxZTabsXr1atx9990QBEFe3hqPF0e+OkYKCgpgNBrRv39/uc2AAQNgNBqDJlaA7X1HEARER0c7LV+zZg3i4+PRs2dPzJkzx2lkLJhjc7Wvn2CODQCcPXsWGzduxD333ON2X7AfM66f0YH6XhPi8zVSi/rll19gsViQmJjotDwxMRElJSUK9cq/JEnCrFmzcOutt6JXr17y8tGjR+MPf/gDOnbsiOLiYjz77LO47bbbUFhYCIPBgJKSEuj1esTExDitL5Bj179/f7z//vu4/vrrcfbsWbz44osYNGgQ9u/fL++Tp2Pl2LFjABCUMfFkw4YNKC8vx5QpU+RlrfF4ceWrY6SkpAQJCQlu609ISAiaWNXX1+OJJ57AxIkTERUVJS+fNGkSOnfujKSkJOzbtw9PPvkkvv32W7m8J1hj44vXT7DGxu69995DZGQk7rjjDqflwX7MePqMDtT3GibJAcpxNAywHZSuy4LVQw89hO+++w7bt293Wn7nnXfKf/fq1Qs333wzOnbsiI0bN7q9STkK5NiNHj1a/js1NRUDBw7Etddei/fee08+ieZKjpVAjokny5cvx+jRo5GcnCwva43Hize+OEY8tQ+WWImiiLvuugtWqxVvvPGG031Tp06V/+7Vqxe6du2Km2++GXv27MFNN90EIDhj46vXTzDGxu7dd9/FpEmTEBoa6rQ82I8Zb5/RQOC917DcIsDEx8dDq9W6fWMqLS11+4YWjB5++GF88skn2Lx5M9q3b3/Jtu3atUPHjh1x+PBhAEBSUhLMZjPKysqc2gVT7CIiIpCamorDhw/Ls1xc6lhpDTE5duwYNm3ahHvvvfeS7Vrj8eKrYyQpKQlnz551W/+5c+cCPlaiKGLChAkoLi5GXl6e0yiyJzfddBN0Op3TcRSssXF0Ja+fYI7N119/jUOHDl32fQcIrmPG22d0oL7XMEkOMHq9Hv369ZN/lrHLy8vDoEGDFOpVy5MkCQ899BDWrVuHr776Cp07d77sY86fP48TJ06gXbt2AIB+/fpBp9M5xe7MmTPYt29f0MTOZDLh4MGDaNeunfxznuP+ms1mbN26Vd7f1hCTFStWICEhAWPHjr1ku9Z4vPjqGBk4cCAqKiqwa9cuuc3OnTtRUVER0LGyJ8iHDx/Gpk2bEBcXd9nH7N+/H6IoysdRsMbG1ZW8foI5NsuXL0e/fv3Qu3fvy7YNhmPmcp/RAfte4/NTAanFrV27VtLpdNLy5culAwcOSDNnzpQiIiKko0ePKt21FvPAAw9IRqNR2rJli3TmzBn5X21trSRJklRVVSXNnj1bys/Pl4qLi6XNmzdLAwcOlK655hqpsrJSXs/9998vtW/fXtq0aZO0Z88e6bbbbpN69+4tNTQ0KLVrV2X27NnSli1bpCNHjkg7duyQMjIypMjISPlY+Nvf/iYZjUZp3bp10vfffy/98Y9/lNq1axfUMXFksVikDh06SH/5y1+clrem46Wqqkrau3evtHfvXgmA9Nprr0l79+6VZ2jw1TEyatQo6cYbb5QKCgqkgoICKTU1VcrIyPD7/jbHpWIjiqKUmZkptW/fXioqKnJ63zGZTJIkSdJPP/0kzZ07V9q9e7dUXFwsbdy4Ubrhhhukvn37BnVsfPn6CbTYXO71JEmSVFFRIYWHh0tvvvmm2+OD9Zi53Ge0JAXmew2T5AD1+uuvSx07dpT0er100003OU2FFowAePy3YsUKSZIkqba2VkpLS5Patm0r6XQ6qUOHDtLkyZOl48ePO62nrq5Oeuihh6TY2FgpLCxMysjIcGsTSO68806pXbt2kk6nk5KTk6U77rhD2r9/v3y/1WqVnn/+eSkpKUkyGAzSb37zG+n77793WkewxcTRF198IQGQDh065LS8NR0vmzdv9vjamTx5siRJvjtGzp8/L02aNEmKjIyUIiMjpUmTJkllZWV+2ssrc6nYFBcXe33f2bx5syRJknT8+HHpN7/5jRQbGyvp9Xrp2muvlWbMmCGdP3/eaTvBFhtfvn4CLTaXez1JkiQtXbpUCgsLk8rLy90eH6zHzOU+oyUpMN9rhMadIyIiIiKiRqxJJiIiIiJywSSZiIiIiMgFk2QiIiIiIhdMkomIiIiIXDBJJiIiIiJywSSZiIiIiMgFk2QiIiIiIhdMkomIiIiIXDBJJiJqZbKzs9GnTx+frW/Lli0QBAHl5eU+WycRkdKYJBMRBaEpU6ZAEAQIggCdTocuXbpgzpw5qKmpwZw5c/Dll18q3UUiIlULUboDRETUMkaNGoUVK1ZAFEV8/fXXuPfee1FTU4M333wTbdq0Ubp7RESqxpFkIqIgZTAYkJSUhJSUFEycOBGTJk3Chg0bnMot6uvr0bNnT0ybNk1+XHFxMYxGI5YtWwYAkCQJCxYsQJcuXRAWFobevXvjP//5j9ftHjt2DOPGjUNMTAwiIiLQs2dPfP755y26r0REvsaRZCKiViIsLAyiKDotCw0NxZo1a9C/f3+MGTMG48aNQ1ZWFoYNG4apU6cCAJ555hmsW7cOb775Jrp27Ypt27bh//7v/9C2bVsMGTLEbTsPPvggzGYztm3bhoiICBw4cIAj10QUcJgkExG1Art27cIHH3yA4cOHu93Xp08fvPjii5g6dSr++Mc/4ueff8aGDRsAADU1NXjttdfw1VdfYeDAgQCALl26YPv27Vi6dKnHJPn48eP43e9+h9TUVLk9EVGgYZJMRBSkPvvsM7Rp0wYNDQ0QRRG33347Fi9ejDfeeMOt7ezZs/Hxxx9j8eLF+O9//4v4+HgAwIEDB1BfX4+RI0c6tTebzejbt6/H7c6YMQMPPPAAcnNzMWLECPzud7/DjTfe6PsdJCJqQaxJJiIKUsOGDUNRUREOHTqE+vp6rFu3DgkJCR7blpaW4tChQ9BqtTh8+LC83Gq1AgA2btyIoqIi+d+BAwe81iXfe++9OHLkCLKysvD999/j5ptvxuLFi32/g0RELYhJMhFRkIqIiMB1112Hjh07QqfTXbLt3XffjV69euH999/H448/jgMHDgAAevToAYPBgOPHj+O6665z+peSkuJ1fSkpKbj//vuxbt06zJ49Wz4JkIgoULDcgoiolXv99ddRUFCA7777DikpKfjvf/+LSZMmYefOnYiMjMScOXPw6KOPwmq14tZbb0VlZSXy8/PRpk0bTJ482W19M2fOxOjRo3H99dejrKwMX331Fbp3767AnhERXTmOJBMRtWI//PADHnvsMbzxxhvyyPDrr7+O8vJyPPvsswCAv/71r3juuecwf/58dO/eHenp6fj000/RuXNnj+u0WCx48MEH0b17d4waNQrdunXzWAdNRKRmgiRJktKdICIiIiJSE44kExERERG5YJJMREREROSCSTIRERERkQsmyURERERELpgkExERERG5YJJMREREROSCSTIRERERkQsmyURERERELpgkExERERG5YJJMREREROSCSTIRERERkYv/DwOqWGa5u3mXAAAAAElFTkSuQmCC", + "text/plain": [ + "<Figure size 800x500 with 1 Axes>" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "# Plot Hg intensity-pixel plot\n", + "plt.figure(figsize=(8,5))\n", + "plt.plot(df_Hg['Pixels'], df_Hg['Intensity'], linestyle='-')\n", + "plt.xlabel('Pixels')\n", + "plt.ylabel('Intensity [a.u.]')\n", + "plt.title('Pixel-Intensity (Hg)')\n", + "plt.grid(True)\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "57b60ccc-47c7-47d3-b18f-d4f9d4d80a37", + "metadata": {}, + "source": [ + "## Calibrate length dimension" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "80ee5a00-fc91-45e1-93e6-aa288451daed", + "metadata": {}, + "outputs": [ + { + "ename": "KeyError", + "evalue": "'Pixel'", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mKeyError\u001b[0m Traceback (most recent call last)", + "File \u001b[0;32m~/miniconda3/envs/main/lib/python3.13/site-packages/pandas/core/indexes/base.py:3805\u001b[0m, in \u001b[0;36mIndex.get_loc\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m 3804\u001b[0m \u001b[38;5;28;01mtry\u001b[39;00m:\n\u001b[0;32m-> 3805\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_engine\u001b[38;5;241m.\u001b[39mget_loc(casted_key)\n\u001b[1;32m 3806\u001b[0m \u001b[38;5;28;01mexcept\u001b[39;00m \u001b[38;5;167;01mKeyError\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m err:\n", + "File \u001b[0;32mindex.pyx:167\u001b[0m, in \u001b[0;36mpandas._libs.index.IndexEngine.get_loc\u001b[0;34m()\u001b[0m\n", + "File \u001b[0;32mindex.pyx:196\u001b[0m, in \u001b[0;36mpandas._libs.index.IndexEngine.get_loc\u001b[0;34m()\u001b[0m\n", + "File \u001b[0;32mpandas/_libs/hashtable_class_helper.pxi:7081\u001b[0m, in \u001b[0;36mpandas._libs.hashtable.PyObjectHashTable.get_item\u001b[0;34m()\u001b[0m\n", + "File \u001b[0;32mpandas/_libs/hashtable_class_helper.pxi:7089\u001b[0m, in \u001b[0;36mpandas._libs.hashtable.PyObjectHashTable.get_item\u001b[0;34m()\u001b[0m\n", + "\u001b[0;31mKeyError\u001b[0m: 'Pixel'", + "\nThe above exception was the direct cause of the following exception:\n", + "\u001b[0;31mKeyError\u001b[0m Traceback (most recent call last)", + "Cell \u001b[0;32mIn[4], line 5\u001b[0m\n\u001b[1;32m 1\u001b[0m C \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39mpolyfit(df_Hg[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mPixels\u001b[39m\u001b[38;5;124m'\u001b[39m],df_Hg[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mIntensity\u001b[39m\u001b[38;5;124m'\u001b[39m],\u001b[38;5;241m3\u001b[39m)\n\u001b[1;32m 3\u001b[0m lambda_Hg \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;01mlambda\u001b[39;00m p : C[\u001b[38;5;241m3\u001b[39m] \u001b[38;5;241m+\u001b[39m C[\u001b[38;5;241m2\u001b[39m]\u001b[38;5;241m*\u001b[39mp\u001b[38;5;241m+\u001b[39mC[\u001b[38;5;241m1\u001b[39m]\u001b[38;5;241m*\u001b[39mp\u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39m\u001b[38;5;241m2\u001b[39m\u001b[38;5;241m+\u001b[39mC[\u001b[38;5;241m0\u001b[39m]\u001b[38;5;241m*\u001b[39mp\u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39m\u001b[38;5;241m3\u001b[39m\n\u001b[0;32m----> 5\u001b[0m lambda_calibrated \u001b[38;5;241m=\u001b[39m lambda_Hg(df_Hg[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mPixel\u001b[39m\u001b[38;5;124m'\u001b[39m])\n\u001b[1;32m 8\u001b[0m \u001b[38;5;66;03m# Plot Hg intensity-wavelength plot\u001b[39;00m\n\u001b[1;32m 9\u001b[0m plt\u001b[38;5;241m.\u001b[39mfigure(figsize\u001b[38;5;241m=\u001b[39m(\u001b[38;5;241m8\u001b[39m,\u001b[38;5;241m5\u001b[39m))\n", + "File \u001b[0;32m~/miniconda3/envs/main/lib/python3.13/site-packages/pandas/core/frame.py:4102\u001b[0m, in \u001b[0;36mDataFrame.__getitem__\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m 4100\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcolumns\u001b[38;5;241m.\u001b[39mnlevels \u001b[38;5;241m>\u001b[39m \u001b[38;5;241m1\u001b[39m:\n\u001b[1;32m 4101\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_getitem_multilevel(key)\n\u001b[0;32m-> 4102\u001b[0m indexer \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mcolumns\u001b[38;5;241m.\u001b[39mget_loc(key)\n\u001b[1;32m 4103\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m is_integer(indexer):\n\u001b[1;32m 4104\u001b[0m indexer \u001b[38;5;241m=\u001b[39m [indexer]\n", + "File \u001b[0;32m~/miniconda3/envs/main/lib/python3.13/site-packages/pandas/core/indexes/base.py:3812\u001b[0m, in \u001b[0;36mIndex.get_loc\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m 3807\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28misinstance\u001b[39m(casted_key, \u001b[38;5;28mslice\u001b[39m) \u001b[38;5;129;01mor\u001b[39;00m (\n\u001b[1;32m 3808\u001b[0m \u001b[38;5;28misinstance\u001b[39m(casted_key, abc\u001b[38;5;241m.\u001b[39mIterable)\n\u001b[1;32m 3809\u001b[0m \u001b[38;5;129;01mand\u001b[39;00m \u001b[38;5;28many\u001b[39m(\u001b[38;5;28misinstance\u001b[39m(x, \u001b[38;5;28mslice\u001b[39m) \u001b[38;5;28;01mfor\u001b[39;00m x \u001b[38;5;129;01min\u001b[39;00m casted_key)\n\u001b[1;32m 3810\u001b[0m ):\n\u001b[1;32m 3811\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m InvalidIndexError(key)\n\u001b[0;32m-> 3812\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mKeyError\u001b[39;00m(key) \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01merr\u001b[39;00m\n\u001b[1;32m 3813\u001b[0m \u001b[38;5;28;01mexcept\u001b[39;00m \u001b[38;5;167;01mTypeError\u001b[39;00m:\n\u001b[1;32m 3814\u001b[0m \u001b[38;5;66;03m# If we have a listlike key, _check_indexing_error will raise\u001b[39;00m\n\u001b[1;32m 3815\u001b[0m \u001b[38;5;66;03m# InvalidIndexError. Otherwise we fall through and re-raise\u001b[39;00m\n\u001b[1;32m 3816\u001b[0m \u001b[38;5;66;03m# the TypeError.\u001b[39;00m\n\u001b[1;32m 3817\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_check_indexing_error(key)\n", + "\u001b[0;31mKeyError\u001b[0m: 'Pixel'" + ] + } + ], + "source": [ + "C = np.polyfit(df_Hg['Pixels'],df_Hg['Intensity'],3)\n", + "\n", + "lambda_Hg = lambda p : C[3] + C[2]*p+C[1]*p**2+C[0]*p**3\n", + "\n", + "lambda_calibrated = lambda_Hg(df_Hg['Pixel'])\n", + "\n", + "\n", + "# Plot Hg intensity-wavelength plot\n", + "plt.figure(figsize=(8,5))\n", + "plt.plot(lambda_Hg, df_Hg['Intensity'], linestyle='-')\n", + "plt.xlabel('Wavelength [nm]')\n", + "plt.ylabel('Intensity [a.u.]')\n", + "plt.title('Wavelength-Intensity (Hg)')\n", + "plt.grid(True)\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "4ac03b00-22c4-400d-9dca-cbf2549fd8f2", + "metadata": {}, + "source": [ + "## Calibrate intensity dimension" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "ca486dfb-3c7b-4a06-a426-266a9a9ce580", + "metadata": {}, + "outputs": [], + "source": [ + "\"\"\"\n", + "epsilon = # Surface emissivity\n", + "lambd =\n", + "\n", + "I_W_true = epsilon*((2hc**2)/(lambd**5)*1/(e**(hc/kT)-1))\n", + "I_W_meas= df_Ox['I_Tungsten [a.u.]']\n", + "R=I_W_meas/I_W_true\n", + "I_plasma_meas=R*I_plasma_true\n", + "\"\"\"" + ] + } + ], + "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.13.9" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/tutorials/module_4/spectroscopy_problem/spectroscopy.py b/tutorials/module_4/spectroscopy_problem/spectroscopy.py index f66734a..404e5ca 100644 --- a/tutorials/module_4/spectroscopy_problem/spectroscopy.py +++ b/tutorials/module_4/spectroscopy_problem/spectroscopy.py @@ -23,6 +23,7 @@ import pandas as pd import matplotlib.pyplot as plt + Hg_data = 'Lampa_Calibrare_Mercur.xlsx' #Ar_data = 'Spectru Descarcare Argon.xlsx' Ox_data = 'Calibrare Intensitate Oxigen.xlsx' @@ -35,45 +36,122 @@ df_Ox = pd.read_excel(Ox_data, header=3, engine='openpyxl') # Trims mercury data to the first 2000 pixels df_Hg = df_Hg[df_Hg['Pixels']<2000] +# %% Plot Intensity vs Pixels for mercury -"Plot Intensity vs Pixels for mercury" +# Add Vertical lines to help match pixel count +wavelength_true = np.array([365.015, 404.656, 407.783, 435.833, 546.074, 576.960, 579.006]) +bounds = np.array([635,1451]) +m = (bounds[-1]-bounds[0])/(wavelength_true[-1]-wavelength_true[0]) +c = bounds[0]-(m*wavelength_true[0]) +pixel_count = m * wavelength_true + c # Plot Hg intensity-pixel plot plt.figure(figsize=(8,5)) plt.plot(df_Hg['Pixels'], df_Hg['Intensity'], linestyle='-') plt.xlabel('Pixels') plt.ylabel('Intensity [a.u.]') +plt.vlines(pixel_count,0,17000,colors='r') plt.title('Pixel-Intensity (Hg)') plt.grid(True) plt.show() -"Calibrate length dimension" -C = np.polyfit(df_Hg['Pixels'],df_Hg['Intensity'],3) -lambda_Hg = lambda p : C[3] + C[2]*p+C[1]*p**2+C[0]*p**3 +# %% Calibrate length dimension using Mercury +# Arrays to be interpolated +wavelength_peak_nm = np.array([365.015, 404.656, 435.833, 546.074, 576.960, 579.006]) # From Literature +pixel_peak = np.array([635, 784, 902, 1324, 1443, 1451]) + +#Polynomial coefficients +C = np.polyfit(pixel_peak,wavelength_peak_nm,3) -lambda_calibrated = lambda_Hg(df_Hg['Pixel']) +# Calibration function +wavelength_Hg = lambda p : C[3] + C[2]*p+C[1]*p**2+C[0]*p**3 +# or +#def wavelength_Hg(p): +# return C[3] + C[2]*p+C[1]*p**2+C[0]*p**3 +wavelength_calibrated = wavelength_Hg(df_Hg['Pixels']) # Plot Hg intensity-wavelength plot plt.figure(figsize=(8,5)) -plt.plot(lambda_Hg, df_Hg['Intensity'], linestyle='-') +plt.plot(wavelength_calibrated, df_Hg['Intensity'], linestyle='-') plt.xlabel('Wavelength [nm]') plt.ylabel('Intensity [a.u.]') plt.title('Wavelength-Intensity (Hg)') plt.grid(True) plt.show() +# %% Calibrate intensity dimension using Tungsten +# Now that we have calibrated the wavelength dimension can can now use +# a tungsten source to calibrate intensity. +wavelength = df_Ox['Wavelength [nm]'] +I_W_measured = df_Ox['I_Tungsten [a.u.]']-df_Ox['I_Background [a.u.]'] -""" -"Calibrate intensity dimension" -epsilon = # Surface emissivity -lambd = - -I_W_true = epsilon*((2hc**2)/(lambd**5)*1/(e**(hc/kT)-1)) -I_W_meas= df_Ox['I_Tungsten [a.u.]'] -R=I_W_meas/I_W_true -I_plasma_meas=R*I_plasma_true -"""
\ No newline at end of file +# Calculate theoretical Intensity of Tungsten using Planks eqn. (TO BE CHECKED) +epsilon = 0.4 +h = 6.626e-34 # J/s +c = 2.99e8 # m/s +k = 1.381e-23 # J/K +T = 1800 # K + +I_planck = lambda wavelength0 : epsilon*((2*h*c**2)/(wavelength0**5)*1/(np.exp(h*c/(k*T))-1)) +I_W_true = I_planck(wavelength) + +#Response Function +R=I_W_measured/I_W_true +#I_plasma_meas=R*I_plasma_true + +# Plot Response function vs lambda +plt.figure(figsize=(8,5)) +plt.plot(wavelength, R, linestyle='-') +plt.xlabel('Wavelength [nm]') +plt.ylabel('R') +plt.title('Response function') +plt.grid(True) +plt.show() + +# Calculate Intensity of Oxygen +I_Ox_measured = df_Ox['I_Plasma Oxygen [a.u.]']-df_Ox['I_Background [a.u.]'] + +I_Ox_true = R * I_Ox_measured + + +# Plot Hg intensity-wavelength plot +plt.figure(figsize=(8,5)) +plt.plot(wavelength, I_Ox_measured, linestyle='-') +plt.xlabel('Wavelength [$nm$]') +plt.ylabel('Intensity [a.u,]') +plt.xlim(750,870) +plt.title('Wavelength-Intensity ($Ox_{measured}$)') +plt.grid(True) +plt.show() + +# Plot Hg intensity-wavelength plot +plt.figure(figsize=(8,5)) +plt.plot(wavelength, I_Ox_true, linestyle='-') +plt.xlabel('Wavelength [$nm$]') +plt.ylabel('Intensity [$W/m^2/sr/nm$]') +plt.xlim(750,870) +plt.title('Wavelength-Intensity ($Ox_{true}$)') +plt.grid(True) +plt.show() + + +# %% Calculate Dencity of Ox at upper state + +from scipy.interpolate import interp1d +from scipy.integrate import quad + +# Create an interpolated callable function +I_Ox_true = interp1d(wavelength, I_Ox_true, kind='linear', + fill_value=0, bounds_error=False) + +# Now integrate +delta = 1 +peak = 846.5 +bound_lower = peak - delta +bound_upper = peak + delta + +numerator, err = quad(I_Ox_true, bound_lower, bound_upper) |