OH Benchmark For None and Simple.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9.28\n"
     ]
    }
   ],
   "source": [
    "import math\n",
    "from radis import get_version\n",
    "print(get_version())\n",
    "\n",
    "conditions = {\n",
    "    'wavenum_min':500,             # cm-1\n",
    "    'wavenum_max':10000,           # cm-1\n",
    "    'molecule':'OH',\n",
    "    'isotope':'1',\n",
    "    'pressure':1.01325,            # bar \n",
    "    'Tgas':300,\n",
    "    'mole_fraction':0.1,\n",
    "    'path_length':1,\n",
    "    'cutoff':1e-27,\n",
    "    'verbose':3,\n",
    "    'databank':'hitemp'\n",
    "    }\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Optimization = 'simple'\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using existing database HITEMP-OH\n",
      "\u001b[32mGenerated dataframe from /home/pipebomb/.radisdb/OH-13_HITEMP2020.h5 in 0.15s (24164 rows, 19 columns, 14.06 MB)\u001b[0m\n",
      "The following columns had the `object` format and were removed: ['globu', 'globl', 'locu', 'locl', 'ierr', 'iref', 'lmix']\n",
      "\u001b[32m... Fetching molecular parameters for all transitions\u001b[0m\n",
      "\u001b[32m... Fetched molecular params in 0.00s\u001b[0m\n",
      "HAPI version: 1.1.2.0\n",
      "To get the most up-to-date version please check http://hitran.org/hapi\n",
      "ATTENTION: Python versions of partition sums from TIPS-2017 are now available in HAPI code\n",
      "\n",
      "           It is free to use HAPI. If you use HAPI in your research or software development,\n",
      "           please cite it using the following reference:\n",
      "           R.V. Kochanov, I.E. Gordon, L.S. Rothman, P. Wcislo, C. Hill, J.S. Wilzewski,\n",
      "           HITRAN Application Programming Interface (HAPI): A comprehensive approach\n",
      "           to working with spectroscopic data, J. Quant. Spectrosc. Radiat. Transfer 177, 15-30 (2016)\n",
      "           DOI: 10.1016/j.jqsrt.2016.03.005\n",
      "\u001b[32mHAPI version: 1.1.2.0\u001b[0m\n",
      "Calculating Equilibrium Spectrum\n",
      "Physical Conditions\n",
      "----------------------------------------\n",
      "   Tgas                 300 K\n",
      "   Trot                 300 K\n",
      "   Tvib                 300 K\n",
      "   isotope              1\n",
      "   mole_fraction        0.1\n",
      "   molecule             OH\n",
      "   path_length          1 cm\n",
      "   pressure_mbar        1013.25 mbar\n",
      "   rot_distribution     boltzmann\n",
      "   self_absorption      True\n",
      "   state                X\n",
      "   vib_distribution     boltzmann\n",
      "   wavenum_max          10000.0000 cm-1\n",
      "   wavenum_min          500.0000 cm-1\n",
      "Computation Parameters\n",
      "----------------------------------------\n",
      "   Tref                 296 K\n",
      "   broadening_max_width  10 cm-1\n",
      "   broadening_method    fft\n",
      "   cutoff               1e-27 cm-1/(#.cm-2)\n",
      "   dbformat             hitemp-radisdb\n",
      "   dbpath               /home/pipebomb/.radisdb/OH-13_HITEMP2020.h5\n",
      "   dlm_log_pG           0.1375350788016573\n",
      "   dlm_log_pL           0.20180288881201608\n",
      "   folding_thresh       1e-06\n",
      "   include_neighbouring_lines  True\n",
      "   optimization         simple\n",
      "   parfuncfmt           hapi\n",
      "   pseudo_continuum_threshold  0\n",
      "   wavenum_max_calc     10005.0000 cm-1\n",
      "   wavenum_min_calc     495.0000 cm-1\n",
      "   waveunit             cm-1\n",
      "   wstep                0.01 cm-1\n",
      "   zero_padding         -1\n",
      "Information\n",
      "----------------------------------------\n",
      "   db_use_cached        True\n",
      "   export_lines         False\n",
      "   lvl_use_cached       True\n",
      "   total_lines          24164\n",
      "   warning_broadening_threshold  0.01\n",
      "   warning_linestrength_cutoff  0.01\n",
      "----------------------------------------\n",
      "\u001b[32mScaling equilibrium linestrength\u001b[0m\n",
      "\u001b[32mScaled equilibrium linestrength in 0.01s\u001b[0m\n",
      "Discarded 97.00% of lines (linestrength<1e-27cm-1/(#.cm-2)) Estimated error: 0.00%, Expected time saved: 2.3s\n",
      "\u001b[32mApplied linestrength cutoff in 0.0s (expected time saved ~ 2.3s)\u001b[0m\n",
      "\u001b[32mCalculated lineshift in 0.00s\u001b[0m\n",
      "\u001b[32mCalculated broadening HWHM in 0.00s\u001b[0m\n",
      "\u001b[32m> Calculating line broadening (725 lines: expect ~ 0.07s on 1 CPU)\u001b[0m\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/pipebomb/Desktop/radis/radis/misc/warning.py:302: MissingSelfBroadeningWarning: Self-broadening temperature coefficient Tdpsel not given in database: used Tdpair instead\n",
      "  warnings.warn(WarningType(message))\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[32m... Precomputed DLM lineshapes (92) in 2.9s\u001b[0m\n",
      "\u001b[32m... Initialized vectors in 0.0s\u001b[0m\n",
      "\u001b[32m... Get closest matching line & fraction in 0.0s\u001b[0m\n",
      "\u001b[32m... Distribute lines over DLM 0.0s\u001b[0m\n",
      "\u001b[32m... Convolve and sum on spectral range 68.8s\u001b[0m\n",
      "\u001b[32mCalculated line broadening in 71.72s\u001b[0m\n",
      "\u001b[32mCalculated other spectral quantities in 0.05s\u001b[0m\n",
      "\u001b[32mSpectrum calculated in 71.80s (before object generation)\u001b[0m\n",
      "\u001b[32mGenerated Spectrum object in 0.04s\u001b[0m\n",
      "Spectrum calculated in 71.84s\n"
     ]
    }
   ],
   "source": [
    "from radis import calc_spectrum\n",
    "import math\n",
    "s1 = calc_spectrum(**conditions,\n",
    "                  optimization='simple', \n",
    "                  )\n",
    "s1.name = 'LDM, simple ({0:.1f}s)'.format(s1.conditions['calculation_time'])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Optimization = None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using existing database HITEMP-OH\n",
      "\u001b[32mGenerated dataframe from /home/pipebomb/.radisdb/OH-13_HITEMP2020.h5 in 0.14s (24164 rows, 19 columns, 14.06 MB)\u001b[0m\n",
      "The following columns had the `object` format and were removed: ['globu', 'globl', 'locu', 'locl', 'ierr', 'iref', 'lmix']\n",
      "\u001b[32m... Fetching molecular parameters for all transitions\u001b[0m\n",
      "\u001b[32m... Fetched molecular params in 0.00s\u001b[0m\n",
      "\u001b[32mHAPI version: 1.1.2.0\u001b[0m\n",
      "Calculating Equilibrium Spectrum\n",
      "Physical Conditions\n",
      "----------------------------------------\n",
      "   Tgas                 300 K\n",
      "   Trot                 300 K\n",
      "   Tvib                 300 K\n",
      "   isotope              1\n",
      "   mole_fraction        0.1\n",
      "   molecule             OH\n",
      "   path_length          1 cm\n",
      "   pressure_mbar        1013.25 mbar\n",
      "   rot_distribution     boltzmann\n",
      "   self_absorption      True\n",
      "   state                X\n",
      "   vib_distribution     boltzmann\n",
      "   wavenum_max          10000.0000 cm-1\n",
      "   wavenum_min          500.0000 cm-1\n",
      "Computation Parameters\n",
      "----------------------------------------\n",
      "   Tref                 296 K\n",
      "   broadening_max_width  10 cm-1\n",
      "   broadening_method    voigt\n",
      "   cutoff               1e-27 cm-1/(#.cm-2)\n",
      "   dbformat             hitemp-radisdb\n",
      "   dbpath               /home/pipebomb/.radisdb/OH-13_HITEMP2020.h5\n",
      "   dlm_log_pG           0.1375350788016573\n",
      "   dlm_log_pL           0.20180288881201608\n",
      "   folding_thresh       1e-06\n",
      "   include_neighbouring_lines  True\n",
      "   parfuncfmt           hapi\n",
      "   pseudo_continuum_threshold  0\n",
      "   wavenum_max_calc     10005.0000 cm-1\n",
      "   wavenum_min_calc     495.0000 cm-1\n",
      "   waveunit             cm-1\n",
      "   wstep                0.01 cm-1\n",
      "   zero_padding         -1\n",
      "Information\n",
      "----------------------------------------\n",
      "   db_use_cached        True\n",
      "   export_lines         False\n",
      "   lvl_use_cached       True\n",
      "   total_lines          24164\n",
      "   warning_broadening_threshold  0.01\n",
      "   warning_linestrength_cutoff  0.01\n",
      "----------------------------------------\n",
      "\u001b[32mScaling equilibrium linestrength\u001b[0m\n",
      "\u001b[32mScaled equilibrium linestrength in 0.01s\u001b[0m\n",
      "Discarded 97.00% of lines (linestrength<1e-27cm-1/(#.cm-2)) Estimated error: 0.00%, Expected time saved: 2.3s\n",
      "\u001b[32mApplied linestrength cutoff in 0.0s (expected time saved ~ 2.3s)\u001b[0m\n",
      "\u001b[32mCalculated lineshift in 0.00s\u001b[0m\n",
      "\u001b[32mCalculated broadening HWHM in 0.00s\u001b[0m\n",
      "\u001b[32m> Calculating line broadening (725 lines: expect ~ 0.07s on 1 CPU)\u001b[0m\n",
      "\u001b[32m... Initialized vectors in 0.0s\u001b[0m\n",
      "\u001b[32m... Calculated Voigt profile (jit=True) in 0.0s\u001b[0m\n",
      "\u001b[32m... Initialized vectors in 0.0s\u001b[0m\n",
      "\u001b[32m... Get closest matching line & fraction in 0.0s\u001b[0m\n",
      "\u001b[32m... Aggregate center lines in 0.0s\u001b[0m\n",
      "\u001b[32m... Aggregate wing lines in 0.0s\u001b[0m\n",
      "\u001b[32mCalculated line broadening in 0.05s\u001b[0m\n",
      "\u001b[32mCalculated other spectral quantities in 0.05s\u001b[0m\n",
      "\u001b[32mSpectrum calculated in 0.13s (before object generation)\u001b[0m\n",
      "\u001b[32mGenerated Spectrum object in 0.02s\u001b[0m\n",
      "Spectrum calculated in 0.14s\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/pipebomb/Desktop/radis/radis/misc/warning.py:302: MissingSelfBroadeningWarning: Self-broadening temperature coefficient Tdpsel not given in database: used Tdpair instead\n",
      "  warnings.warn(WarningType(message))\n"
     ]
    }
   ],
   "source": [
    "s2 = calc_spectrum(**conditions,\n",
    "                   optimization=None,\n",
    "                  )\n",
    "s2.name = 'None ({0:.1f}s)'.format(s2.conditions['calculation_time'])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Comparison of performance and accuracy of both methods"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAmsAAAHACAYAAAD9fy9UAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAxOAAAMTgF/d4wjAACPZElEQVR4nO3de1xT9f8H8NcYMC4yQC6C4D3vaeYl85L3S2reSE0t09S8pHYR7ZdlIdlXvH7LsrS8VOYtI7VSK80bpnnXRE3xghcUEEW5CGxs+/z+oO3L2IBtbOy97f18PHg4zm2fs9eY733OOZ8jEUIIMMYYY4wxktzs3QDGGGOMMVY6LtYYY4wxxgjjYo0xxhhjjDAu1hhjjDHGCONijTHGGGOMMC7WGGOMMcYI42KNMcYYY4wwd3s3wJXIZDKEhITYuxmMMcYYs5OMjAwoFAqz1uFirRKFhIQgJSXF3s1weRqNBm5u3Klsb5wDHZwFDZwDHbbMIjIy0ux1+F3BXI5SqbR3Exg4B0o4Cxo4BzqoZcHFGnM5fIc1GjgHOjgLGjgHOqhlwcUaczkymczeTWDgHCjhLGjgHOiglgUXa8zlmHtiJ7MNzoEOzoIGzoEOallwscYYY4wxRhgXa8zleHl52bsJDJwDJZwFDZwDHdSy4GKNuZyCggJ7N4GBc6CEs6CBc6CDWhZOXaxpNBosWbIEXl5emDNnjtnrb9q0CS1atEBoaCjq1KmD2NhYqFQq6zeUMcYYY6wUTjso7vXr1zFmzBgkJydbdKLgypUrMXnyZGzevBlRUVE4f/48evTogaSkJKxfv94GLWaVhVr3tqviHOjgLGjgHOigloXT9qwNHDgQ/fv3x9dff232ug8ePMCMGTMwcuRIREVFAQCaNm2K2NhYbNiwAXv27LF2c1klota9Td733wP//GP1zXIOdHAWNHAOdFDLwmmLtR07diA6Otqi20Vs3rwZ2dnZGDRokN70wYMHAwBWr15tjSYy5hiGDwdeftnerWCMMZfltMWaJffe0kpISAAANGvWTG96SEgIwsLCcODAgQq1jdkXte5tR2CLsbw5Bzo4Cxo4BzqoZeG0xVpFXL58GQAQFhZmMC88PBx37txBXl5eZTeLWQm17m3KfvzxRwBAYmIibt68adVtcw50cBY0cA50UMuCizUjHj58CADw8fExmKedlpWVZXRdhUKB7Oxsoz+MOZohQ4YAKHpfR0dH27k1jDHmmpz2alBrkEgkJk0rLi4uDrGxsUbn+fn5IT8/H+7uRS+7SqWCu7s7hBBQq9Xw8PCARqPRPVar1dBoNPD09IRKpdI9LiwshBACMpkMSqVS91h71auXl5fuW0Fpj7XLSyQSeHp66h57eHhAqVTCzc0N7u7uusdSqRSFhYWQSqVwc3PTPZZIJLr9cJR90mg0UCgUTrVPtsqpuPj4eBQUFFhtn4QQuixc5b1HdZ80Gg3y8/Odap8cMSeNRqP7m3CWfXLUnAAgPz/fJvtkCYmgdmt5K9u/fz+6du2KmJgYk8daa9OmDU6cOIGcnBxUqVJFb17Lli1x+vRp5ObmwtfX12BdhUJR6lAhTZo0QUpKitn7wKwrPz8f3t7e9m6GQ5BIJBAAjgN4CoA1Py44Bzo4Cxo4BzpsmUVkZKTZtQAfBjWifv36AIDU1FSDeWlpaQgPDzdaqAFF1bNcLjf6w2iw9JsNsy7OgQ7OggbOgQ5qWXCxZkSnTp0AFJ1UXVxGRgZSU1PRuXNnezSLWYm2i5vZF+dAB2dBgz1zKHnag6spuf/U/ia4WAOQnp6OwsJC3e9Dhw6Fn58ftm3bprfc1q1bAQBjx46tzOYxK7PkjhbM+jgHOjgL+9m/fz8kEgkkEgmkUinGjBlTqc+fl5eHqVOnYvny5Vbf9oULF9C2bVtIJBLs37/f6tu3pgULFmDGjBm6vwVqfxMuX6wdOnQI1atXx4ABA3TTgoKCsHDhQmzYsAFbtmwBUPSmi4mJwfDhw9GzZ097NZdZAbXubUdgi1eMcyhbly5ddP+JF/8p7z/zhw8fGl1PIpHA3d0dwcHBaNasGcaPH49ffvkFGo3GaBbffPONwfqvvPKKye1/8803Dda/fv26ma9CxU2cOBHh4eE4f/58pT+3Kdq3b4/U1FR88sknlf7cmZmZ6NixIxISEjBs2DDd9OvXr5f6Hir58+abbxpsV61WY8GCBWjZsiWOHTtmlbbm5+fjyy+/RMeOHREREQFPT0+EhITg2Wef1XWkVMTLL7+Mn3/+GZ06dUJWVha9zyfhpObPny+qVasmAgMDBQDh6+srqlWrJgYMGKC33NmzZ0VQUJCYMGGCwTbWr18vmjdvLkJCQkStWrXEBx98IJRKpcVtioiIsHhdZj0qlcreTXAYAIQAxPGicXGtum3OoWz3798Xqampol27dgKAiI6OFqmpqeLhw4dlrqfRaERqaqpITU0V+De3H3/8UaSmpoqbN2+KU6dOiRUrVoj27dsLAKJx48biyJEjBtvJy8sTqamp4scff9Rtx93dXVy5cqXctqelpQlvb2/detr22CNzHx8fAUAsWrSo0p/bHF9//bUAIEaPHl0pz6fRaETnzp1FcHCwuHv3rt685ORkAUDUrVtXNGzY0OhPWFiYACBWrFiht+69e/fEU089JYKCgsT69etFrVq1BACxb98+i9uam5srnnjiCQFA9OnTRyQkJIhr166J33//XbRt21YAEG+++abF29e6ffu28Pf3Fz179hSFhYUV3l5pLKkFnLZYo4iLNRry8vLs3QSHYctijXMwTefOnQUAERMTY/a62mKptP8ov/rqK+Hu7i48PT3Frl27jC6zb98+AUBUrVpVABCvvPJKuc8bHR2tW97efQLz5s0TnTp1Ejdu3LBrO8pT2cXaV199JQCIBQsWGMzTFmvJycmlrv/iiy+KgIAAkZubqzf99OnT4vnnnxfp6elCCGGVYm3JkiUCgKhTp45QKBR68+7fvy/8/PwEAHHixAmLn0MrJibGaBFqTZbUAi5/GJS5HkvuF8usj3Owv1dffRX/+c9/oFQqMWTIEKSlpZW67PTp0wEA3333Ha5du1bqcvfu3cOKFSuMHh6zh1mzZuHAgQOoWbOmvZtChkajwfz58+Hp6YnXXnvNYL6npyeeeOKJUk+yT0tLww8//ICxY8cajIzQrFkzxMfHIzQ01Grt1d5VqFWrVgZtqlq1Kho0aAAAVjkv7o033oBUKsWSJUusOlRRRfGnJXM52gESmX1xDjRMnz4dkZGRyM7Oxrx580pdbsSIEWjUqBFUKhU++uijUpdbvHgx2rRpg2eeecYWzYVKpcKqVavw9NNPIzg4GHK5HE888QQmT56Mffv26a7qK3nOXfFxNtetW2cw79y5c+jTpw/8/PxQrVo1jB49Gnfv3gUAXL16FYMGDUJAQAACAgIwdOhQ3LlzR69dc+bM0dvmN998g23btuGZZ56Bv78/fH190bFjR/zyyy8W7XdGRgaio6NRv359yGQyBAYGomvXrti4caPZ2zp8+DCuXbuGNm3aGIwlCgDVq1fHmTNnUL16daPrr1ixAiqVClOmTDGYJ5VKzW5PeVq0aAEA+OeffwwKqIKCAiQnJwMwvOuQqe+V4gIDA9GqVStcvnwZR48etfq+WMzq/XusVHwYlAY+/GY68GFQu7PlYVCtKVOmCAAiPDzcYJ72MGhycrJYt26d7ty1q1evGix77949UaVKFbFnzx7detZ+34wYMUK4ubmJuLg4cf78eXHlyhWxYcMGUaNGDb1Dd9pz7qKjow1ev/z8fL15I0aMED179hSHDx8Wly5d0h0Ka9mypbhz544YPHiwOHz4sLh48aJ47733BADx5JNPCo1Go9tmTk6OSE1NFcOGDRMARJcuXUSLFi3E7t27xY0bN8Svv/4qGjZsKACIZcuWGexXWYdB//nnHxERESF8fX3F559/Li5evCgSEhJE//79BQAxadIks17DOXPmCABGz9Uuj1KpFGFhYaJ///4mLW+Nw6CFhYVi4MCBun29c+eO0Gg04tq1a2Lw4MECgPD09BRJSUl665n6Xilp3LhxAoD48MMPLW5zWSypBfirLXM5fPiNhrJyEELg+vXryMvLq8QWVZyPjw9q165N70qycrRs2RJA0UDgN2/eLPWQ4YgRIzB37lxcunQJ//nPf7B69Wq9+f/973/RvHlzdOvWzSZDNVy/fh0bN25Ev3798M477+im16tXDzVq1NDrzfP29oa3t7fRniMvLy+EhYXp5sXHx+PGjRsIDw8HUNRLtmfPHvz555947rnnsH37dt28jz76CHv27MGRI0dw+PBhdOjQAQBQpUoVVKlSRTfq/dGjR3H9+nXd4cCaNWvi119/RaNGjTB9+nT06NEDDRs2LHef1Wo1hgwZgtu3b+P777/XXbXZsGFDdOzYEa1bt8aKFSvQo0cPPP/88ya9jtoxRCMjI01avrjNmzcjLS0Nr7/+utnrWsrd3R3btm3Db7/9hrfffhvVq1eHu7s7VCoVAMDX1xfLly/XDWgPmPdeKUn7upQca9WeuFhjLscW3fTMfKXloNFoMGDAAOzYsaOSW2Qd/fr1w88//+xQXwqqVaume5yenl5qsebm5ob3338fL730EtauXYvZs2ejTp06AIAHDx5g2bJl+P77723Wznv37gEA7ty5A41Go/caP/300/jss88QFBRk9nZ79OihK8a0WrdujT///BPVqlUzmPfUU0/hyJEjOH36tK5YK2no0KEG523VqVMHffv2xbZt27B8+XKThuvYtm0bzp8/jxo1amDo0KF68yQSCSZOnIiJEydi2bJlJhdrN27cAFB0vpe5PvvsMzRu3Bg9evQwe92K+PzzzzFjxgy0bt0aO3fuRK1atXDt2jXs3r0bL774Ip566im95SvyXgkODgYAuww1UxrH+TRhzEqKD4DM7Ke0HO7du+ewhRoA7NixQ/cfhaMonkV5vYIjRoxAw4YNoVKp8J///Ec3/ZNPPkGDBg3w7LPP2qydTZo0QVhYGE6fPo1OnTohPj4e+fn5AIp6X6ZOnQo/Pz+zt1uvXj2DadpbBBqb5+/vD6CoQC1N48aNjU5v3bo1AODAgQMmtW3Xrl269YxlU7duXQDAkSNHoFarTdpmTk4OgKLbI5rj+PHjOHr0KKZNm2bWehX1+eefY+rUqYiIiMCePXvQp08fNGnSBM899xzi4uLw/PPPY/r06Xrv44q8V7Svi/Z1ooCLNeZyuGfNfLY4qFdaDsHBwejXr58NnrFy9OvXT/fN3FFkZmbqHpd3FZ+bmxtmz54NAPj222+RnJyMrKwsLF26FO+//75N2+nj44Pt27ejSZMmOHToEIYOHYrg4GBERUXpBvi1hLEeJm1hZGyetpemrOKotNcxJCQEAEy+kffNmzcBANu3b9cdai3+ox3QvaCgQC/HsmgPH5r7Wbh06VL4+/vj5ZdfNmu9ipo/fz4A4PXXXze4GtTHxweTJ0/Gxx9/rPf+q8h7Rfu6aF8nCvgwKHM5jnR4ypmVloObmxt++eUXhz1nTXtY0JGcOnUKQNFVgKYMcTFixAh89NFHunPXatasidq1a6N///62bipatWqFc+fOYf/+/di4cSO2bt2q++ncuTN+/vlnXa+YqcrqTbT0/MPy1jN3u0OGDMGHH35Y5jKBgYEmbUt71aRSqTT5+dPT0/HDDz9gypQpBsN12NLdu3d1hW1pvZWNGjUCUHSValxcnO61tfS9or3VVGXuZ3m4WGMup7CwkIeNIKCsHCQSiUMWPfYihIBCoYCnp6fZX0YKCwt1t9UreU5UaaRSKWbPno1Ro0Zh7dq18PHxwapVqyrtwgqJRIKuXbuia9eu+OKLL7Bt2za8/vrrOHDgABYsWKB3eNZecnNzjU7PyMgAAERERJi0nVq1agEoyumxxx6zStvCwsKQmJiI7Oxsk9cpa7iOylLWFzwAyMrKQm5urt7hTUveK9rDn8XP5bQ37mJgLocPg9LAOVjPjRs34O3tjYSEBLPXXbx4MVJTU+Hv7493333X5PW0564VFhYiIiICUVFRZj+3uS5dumQw2K67uzuGDBmChQsXAgBOnz5t83aYorQrCU+cOAEA6Ny5s0nb6d27NwCUeY/NESNG4IUXXjC5bdpBZNPT001avrCwEF9++SX69u1r9Bw+WwoODtb1GF68eNHoMtrpvr6+uit8K/Je0Q4ObcrVupWFizXmchxtWAVnxTnY35dffon3338fMpkMW7ZsMWvUealUikWLFmHgwIFYuHCh2T16KpVKN4Dsw4cPTVonNTUVS5cuxdWrVw3maU8up3Kngvj4eINzyJKTk7Fz5054enpi8uTJJm1nwIABaNasGW7evIn4+HiD+QcPHsSmTZvQqVMnk9umvYL1ypUrJi3/ww8/IDU11abDdZw4cQKNGzfG888/r3ehgJubG1566SUAwLJlywzOI1MoFFixYgUAICoqSve5UpH3yqVLlwAAHTt2rOBeWQ8fC2IuR6VSwcPDw97NcHmcQ9mSk5Px6NEjPHr0CEDRuTvnzp0zumzx0fSFEAY9JpmZmUhLS4NKpcK9e/dw7NgxfPPNN/jrr7/QtGlTrFy5Eu3atdNbR6lUIjMzU1dwZGRkwMvLC1WrVtWd5N2/f3+D89RKrgf8r6ciJCRE16OqvTuAXC5HQECAWa9N//798eGHH+KJJ56AEAJHjhzB//3f/6Fq1aqYMWMGACA/P193WAwoOiyZlpYGf39/SKVSZGZmGszTXkxgybySJ77369cPffr0wbx581C/fn38888/eOONN6BUKrFs2TJdr4329crKytK1Oy0tTXcBgVQqRXx8PHr06IExY8bg1q1b6NevHzQaDfbu3Yv3338fgwYNwqRJk0x+/Xr37g2ZTFZmb11xn376KRo1amTycB2ZmZm68+G0F2Fo34OA/vtA65tvvsHFixdx8eJFnD17Fq1atdLN++ijj3DkyBEcP34czz33HGbPnq0buiMmJgY3btxAo0aNsGjRIoO2mPJeKU6j0eDYsWPw9vZGr169TNrfSmHdcXlZWfgOBjQolUp7N8Fh4N87GJywwUj0nEPZtHcuMOdn37594sGDB6XOl0qlomrVquLxxx8X48aNE9u3bxcajcZoFsXvQFDyOcpS2nooMWK89m4IY8eONfk1KSwsFNu3bxejRo0SjRs3FnK5XPj5+YkmTZqI6dOni5SUFN2y2jsClPz5+uuvy9w3S+dpjR49Wvc8Bw4cEN26dRP+/v7C29tbtG/fXvz8888mvV4l71hx79498fbbb4sGDRoImUwmgoODRbt27cSqVass+lt65ZVXBACxd+/eMpc7duyYACA+//xzk7dd3nvX2J0D9u/fL4KDg0WnTp1Efn6+wXyFQiGWLl0qOnToIAICAoRUKhX+/v7i6aefFgsWLBCPHj3SW96c90pxO3bsEADE+PHjTd5fc1lSC0iEIHSnUicXGRlp8uXazHYKCwu5R8dEEokEAsBJAK0Bq97YmHOgo7KzUCgUePLJJ3H79m2cOHFCb+R5RzdmzBh8++23+PrrrzFmzBiz1q3MHO7cuYPGjRujRYsWJo/55gratGmDK1euIDEx0aI7PJjCklqAz1ljLofS2DmujHOgo7KzeOmll5CSkoLvv//eqQq1iqrMHKpXr44NGzbg8OHDiImJqbTnpWzGjBn4+++/sXnzZovuhGFLXKzZmEKhQHZ2tlmXSDPb4mE7aOAc6KjsLGbOnIkTJ07Y9G4Hjqiyc+jXrx927tyJL774Aps2barU56ZmxYoV2LBhA3bv3o2ePXuS+3yi1RonFBcXh9jYWACw6DYozPr4yL/5bHHdJudAR2VnUfI+js4gNzcXubm5ulsaZWVllXrxQWns8TfRs2dPXLlyBbdv367056akffv2SEpK0g39Qe3zic9ZszGFQqEbDblJkyZ8zhoB+fn58Pb2tnczHIL2nLVTAFrBuh9gnAMdnEXFzZkzR/fFvLh9+/ahS5cuJm2Dc6DDlllYcs4aF2uViC8woEGlUpHr4qbKlsUa50AHZ0ED50CHLbPgCwwYM4GlN3tm1sU50MFZ0MA50EEtCy7WmMvRDtLI7ItzoMNRs6D2H6qltPvhqDk4I2pZcLHGXA6P7UUD50CHo2WRl5eHqVOnYvny5ZX6vCkpKejTpw8kEgm++eYbq233zTffxPz58/l+uYRQ+5vgYo25HGrfmFwV51C22rVr6+6dKZFIyrxR9/Xr1/WW1f7MmTPHpOdypCwyMzPRsWNHJCQkYNiwYUaX+fXXX9GjRw8EBgaiSpUqaN26NVasWFGhcy7XrFmDxx9/HL/99pvF2yjNhAkT8PHHH2Pw4MF698Vk9kPtb4KLNeZynOXQiaPjHMp2/PhxpKam6u7ZuXnzZnz66adGl61RowZSU1ORmpoKAPjkk0+Qmppq9N6HxjhKFkIIREVF4datW9izZw9CQkIMlpk/fz769u0LHx8f/PHHHzh58iQ6d+6MyZMn44UXXjB7XwsLC9GnTx9MmzYNsbGxeOaZZ6y1OzqPP/44du/ejd9//x2vvfaa1bfPzEftb4KLNeZyTB3ziNkW51C2kJAQhIWFwdPTU3dV2owZM3DkyBGDZaVSKcLCwhAWFgYA8Pf3R1hYmG7MqPI4SharVq3CgQMHMHPmTKOF2sGDBzFr1iw0a9YMW7ZsQatWrdCwYUMsWbIEo0aNwg8//IDPPvvMrOd89OgRCgoKkJiYiDfeeANubrb5b7N58+Z48cUXsWrVKuzbt88mz8FMR+1vgos15nL4Nkfms8WguJyD6YYOHYq2bduisLAQw4YNw71796y6fUfIQqPRYP78+fD09Cy190l72Pf11183GHZh+vTpAID//Oc/Zu2vv78/9u7di7p161rWcDO8/vrrAIB58+bZ/LlY2aj9TXCxxlwOte5tV8U5mM7T0xPx8fEIDQ3FrVu38OKLL1r0+gkhsG7dOnTp0gUBAQHw9vZGo0aN8O677+L+/ft6y65bt87g/Lfz58+jf//+CAgIgK+vL7p06YKjR4+W+nwZGRmIjo5G/fr1IZPJEBgYiK5du2Ljxo1mt/3w4cO4du0a2rRpY7THMCMjA/v37wcAdO/e3WB+ixYtULVqVWRkZGDv3r0mP692/82Vl5eHJUuW4Mknn0RgYCACAwPRqlUrzJgxA8eOHTO6TpMmTRAaGoo//vhDd0ib2Qe1zycu1pjLoda97ao4B/NERkZi06ZNkEql2LVrFz788EOz1tdoNBgxYgRGjRqFunXrYv/+/bhw4QKio6OxevVqtGzZEpcvX9YtP2TIEKSmpiI6OhoAkJSUhMmTJyM6OhonT57E0qVLcezYMfTo0cPoAJ8XL17Ek08+iS+//BJvvfUWzp49i59//hl+fn4YOXIkJk+ebFb79+zZAwBo2rSp0fknTpyARqOBTCZDnTp1jC7TsGFDAEXnA9qSRqNB79698f7772Py5Mk4cuQIjh49iokTJ2LNmjVo27at0fU8PT3RrFkzAMAff/xh0zayslH7fOKhkpnLKSws5EvkCeAczNe1a1fMnz8fM2fOxNy5c9GuXTv07t3bpHUXLlyI77//HgMGDMCaNWt001999VX4+flhxIgRGDJkCE6dOgWpVAovLy+98942b96Ma9euoWbNmgCAevXq4cKFC/j444/xzTffYPbs2bptqtVqDBkyBLdv38b333+vu2qzYcOG6Nixo+7qzB49euD55583qf2JiYkAiopWY65evQoAqFatWqnbCA8P11vWVhISEvDnn39iypQpmDBhgm56gwYN4OPjg1GjRhldr7CwULd/2v1l9kHt84mLNeZy+A5rNJSbQ0EBoFRWTmOsxdMT8PKy6VNoD6P98MMPePHFF3H69GnUqFGjzHWUSiUWLVoEAJg2bZrB/IEDByIyMlLX+zV48GCDZXr06KEr1LS0PURnzpzRm75t2zacP38eNWrUwNChQ/XmSSQSTJw4ERMnTsSyZctMLtZu3LgBAKhatarR+dnZ2QAAHx+fUrehvddjVlaWSc9pKe05hcZ6HPv161fqRQ5CCAQHBwMoGo6F2Q+1/ye4WGMuRyaT2bsJDCbkULcu4Gjn7VSvDty+bfOnWbNmDc6fP48LFy5g6NChSEhIKPOwzcmTJ5GZmQkAaNOmjcF8mUyG1q1bIyUlBbt27TJarD322GMG07SF04MHD/Sm79q1CwDQunVro+d7aU/WP3LkCNRqtUk9GDk5Obq2Wkr7H7Al56CZo127dvD29sZPP/2E5557DlOmTEGPHj3g4eGBwMBATJ061eh6MplMt3/a/WX2Qe3/CS7WmMtRKpXwsnHvBytfuTlcu+aYPWuVoEqVKti6dSvatGmDo0ePIjo6uswhKW7evPlv8zzh7+9vMF+pVOqGwtD2YJUUFBRkME07ynvJAUS1z7d9+3ajFwNoT94uKChAZmam0WE4StJenVdaYSeXywEUndhfmoKCAr1lbSUiIgLbtm3D+PHjsWPHDuzYsQMBAQF47rnnMG7cOHTp0sXoekqlUrd/1K5GdDXU/p9w+mJt06ZNmD9/Pu7cuQNfX1+MGTMG7733nsFl3cY8evQIS5Yswffff4+MjAx4e3ujbdu2iI2NRePGjSuh9cwWqHVvu6pyc/DysvkhRUfWoEEDfPvtt4iKisKyZcvQsWPHMu9yUJbiWZTW62RJb9SQIUPKvRAiMDDQpG1pD28qSyng69WrBwBIT08vdRvaKyy1y9pSr169kJycjF9//RWbNm3CL7/8gnXr1mHdunUYOnQoNm7caFB4CiGgUCgAAL6+vjZvIysdtf8nnPpq0JUrV+Kll17CBx98gLt372L79u1YsWIFRo8eXe66Qgj06dMHixcvxrJly3D37l2cOHEC9+/fR5s2bXD+/PlK2ANmC9S6t10V51BxgwYNwjvvvAMAGD9+PC5evGh0uVq1agEoKnSMna8lk8mQkZEBAAbnpVlC+3yFhYV47LHHyvwx5YszAN2Av9pz00pq3bo13NzcoFAokJycbHSZpKQkAMYPBduCVCrFc889h3Xr1iE9PR2rVq2Cn58ffvjhB72LPLRkMpnu8GdZF0ow26P2+eS0xdqDBw8wY8YMjBw5ElFRUQCKLvmOjY3Fhg0bdJeBl+bQoUM4ePAgJk6ciK5duwIoGlF86dKlePToEZYuXWrzfWC2of3myuyLc7COjz76CL169UJubm6pJ+u3atVKd36ZsWErFAqFbnqvXr0q3CbtFaqljScGACNGjDCrJ7BBgwYASu85CwkJ0R1eNPb5fubMGdy/fx/BwcHo1q2byc9ricOHD2Pu3Ll607y8vDBu3Di8/fbbAIDTp08brKdQKJCWlgbgf8OMMPug9vnktMXa5s2bkZ2djUGDBulN1544u3r16jLX117FU3LUam33+a1bt6zUUsbos+3p2Kwi3NzcsGHDBtSqVQsXLlwwuoyHhwdmzpwJAEbPbdu2bRtu376N5s2bY8CAARVu04ABA9CsWTPcvHkT8fHxBvMPHjyITZs2oVOnTiZvs0OHDgCAK1eulLqM9g4Gn376qcE5X//9738BALNnzzbozdu5cyfq1atntftyJiUlYdGiRXj48KHBPO2N2kvrwbx06RIAoGPHjlZpC3MOTlusJSQkAIBugEEt7f32Dhw4UOb62oEXtd3mWtrf+VuP46J00qgr4xzKlpGRgbS0NCiVSuTn5yMtLU13qLKkoKAgbNmypczX9O2338awYcPw888/Y9y4cfj7779x/fp1rFq1CpMnT0aNGjUQHx+vO49KqVQiLS0Nubm5AIDc3FxdewAgLS1Nd4VpyWWlUini4+NRo0YNjBkzBh9//DGSkpJw8eJFfPHFFxg0aBAGDRqESZMmmfx69O7dGzKZrMzeumeeeQbz5s1DYmIioqKicOrUKSQlJSE6OhrfffcdhgwZYnToki+++ALXrl3D8uXLDe7kABQN9ZGWlqa3/8amFZeTk4O+ffti586duHbtGv755x8sW7YMixcvRu3atfXGX9PKz8/HpUuXEBERgdatW5v82jDrI/f5JJxUmzZtBACRnZ1tMO/JJ58UAMSjR4/K3Mb06dOFr6+v2LFjh1CpVOL69euiY8eOolatWiIlJcXsNkVERJi9DrO+vLw8ezfBYQAQAhCnAWHtjwvOoWy1atUS+Pd11/7UqlWrzHXWrFkjAIivv/7a6HyNRiPWrl0rnnnmGSGXy4VMJhMNGjQQb731lrh3757esvv27TN4fgBi3759QghhdF5MTIzeNu7duyfefvtt0aBBAyGTyURwcLBo166dWLVqlVAqlWa/Jq+88ooAIPbu3Vvmcjt27BDdunUT/v7+wsfHR7Rs2VJ88cUXQqPRGF1+48aNQi6Xi6FDhxqdP3r0aKP7W/I10crLyxPff/+9eP7550WDBg2Er6+vCAgIEM2bNxexsbHi/v37Rp/n448/FgDEvHnzyn8xmE3Z8vPJklpAIgSxSx6spEGDBrh8+TJUKpXBFTcdO3bEoUOHcOfOHd2I1qVZvHgxPvjgA6jVaiiVSjz77LNYvXo1qlevbnR5hUJR6rHuJk2aGB0kkVWu/Px83eCYrGwSiQQCwBkAT8K6V0hxDnQ4ShZ37txB48aN0aJFi3KPjjgahUKBunXrwsvLC+fOnXOIPJyZLf8mIiMjza4FnH7oDmOXm5tyCbpKpcILL7yAAwcOID4+Hj179kRmZiaio6Px9NNPY9u2bWjZsqXBenFxcYiNjTW6TT8/P+Tn5+vOl1CpVHB3d4cQAmq1Gh4eHtBoNLrHarUaGo0Gnp6eUKlUuseFhYUQQkAmk0GpVOoea4tELy8v3XhCpT3WLi+RSODp6al77OHhAaVSCTc3N7i7u+seS6VS3e033NzcdI8lEoluPxxln8S/l8c70z7ZKqeSCgoKrLZPQNF4W6703qO6T0II5Ofnk9+nwMBArFu3DlFRUXj33Xfx0UcfOUVOUqkUL730EnJzc7Fjxw4AMMjD0fbJ0f+eJBIJ8vPzbbJPlnDanrU2bdrgxIkTyMnJMRiUsWXLljh9+jRyc3NLHctmxYoVmDx5Mj799FO9cxwUCgXq1KkDPz8/XLhwwaDXjnvW6HOUXgQKuGfNNThaFrt378bIkSPx2WefYfjw4fZuToW988472Lp1K9avX8/nqhFBrWfNaS8wqF+/PoD/DYJYXFpaGsLDw8scdFB76XfJK3K0t2VJSkrC5cuXDdaTyWSQy+VGfxhjjFVcz549ceXKFTRv3tzeTbGK559/HomJiboL2xgryWmLNe0l4YmJiXrTMzIykJqais6dO5e5/qNHjwAYP2Tq5lb0spU2OCOjjdxVPi6Kc6DDEbPw9/dHkyZN7N0Mq2jTpg08PT0dMgdnRS0Lpy3Whg4dCj8/P2zbtk1v+tatWwEAY8eO1Zuenp6uG/8GgK4r+q+//tJbrrCwECdPnoSXlxd/C3JQ2nMImH1xDnRwFjRwDnRQy8Jpi7WgoCAsXLgQGzZswJYtWwAAFy5cQExMDIYPH46ePXvqlj106BCqV6+uNxjktGnTEBkZidjYWBw8eBBCCOTk5GDq1KlISUnB7Nmz+d5tjDHGGLM5p74adNKkSZDL5YiNjcWkSZPg4+ODCRMmYPbs2XrLyeVyBAYG6o0oHRISgmPHjmHOnDl46aWXkJWVBSEEHn/8caxbtw4vvvhiZe8OsxJq93xzBLa4gwHnQAdnQQPnQAe1LJz2alCKLLkChFmfo135Zk/aq0H/BtACfDWos+IsaOAc6OCrQRmzM0vHuWHWxTnQwVnQwDnQQS0LLtaYy/H09LR3Exg4B0o4Cxo4BzqoZcHFGnM5pQ1azCoX50AHZ0ED50AHtSy4WGMuh1r3tqviHOjgLGjgHOiglgUXa8zleHh42LsJDJwDJZwFDZwDHdSy4GKNuRylUmnvJjBwDpRwFjRwDnRQy4KLNeZytLcLY/bFOdDBWdDAOdBBLQtarWGsEri7O/VY0DZhi7M3OAc6OAsaOAc6qGXBxRpzOdS6t10V50AHZ0ED50AHtSy4WGMuh1r3tqviHOjgLGjgHOiglgWt1jBWCaRSqb2bwMA5UMJZ0MA50EEtCy7WmMspLCy0dxMYOAdKOAsaOAc6qGXBxRpzOdS+MbkqzoEOzoIGzoEOallwscZcDrVzEVwV50AHZ0ED50AHtSxotYaxSkCte9sRCBtsk3Ogg7OggXOgg1oWtAYScUIKhYLcDWFdHbXubVfFOdDBWdDAOdBBLQvuWbOxuLg4+Pv7w9/fH9nZ2fZuDgO9G/S6Ks6BDs6CBs6BDmpZcLFmY7NmzUJWVhaysrIgl8vt3RwGQKVS2bsJDscWH1ucAx2cBQ2cAx3UsuDDoDYmk8kgk8ns3QxWDLXbiLgqzoEOzoIGzoEOallwzxpjjDHGGGFcrDGXQ61721VxDnRwFjRwDnRQy4KLNeZyqHVvuyrOgQ7OggbOgQ5qWXCxxlyOELYYNYyZi3Ogg7OggXOgg1oWXKwxl6NWq+3dBAbOgRLOggbOgQ5qWXCxxlyOh4eHvZvAwDlQwlnQwDnQQS0LLtaYy9FoNPZuAgPnQAlnQQPnQAe1LLhYYy6HWve2I7DFoLicAx2cBQ2cAx3UsuBijbkcat3bropzoIOzoIFzoINaFlysMZdD7RuTq+Ic6OAsaOAc6KCWBRdrzOVQOxfBVXEOdHAWNHAOdFDLgos15nI8PT3t3QQGzoESzoIGzoEOallwscZcDrXbiLgqzoEOzoIGzoEOalk4dbG2adMmtGjRAqGhoahTpw5iY2PNCuDhw4eIjo5G/fr1ERoairCwMPTs2RM7duywYauZrVHr3nZVnAMdnAUNnAMd1LJw2mJt5cqVeOmll/DBBx/g7t272L59O1asWIHRo0ebtP7du3fRpk0b5Obm4vjx47h79y527tyJ06dPY/fu3TZuPbMlat3bropzoIOzoIFzoINaFk5ZrD148AAzZszAyJEjERUVBQBo2rQpYmNjsWHDBuzZs6fcbbz22mvw9/fH8uXLERAQAABo2bIl3n33XURERNiy+czGCgsL7d0EBs6BEs6CBs6BDmpZOGWxtnnzZmRnZ2PQoEF60wcPHgwAWL16dZnrX7t2DT/++CPGjBkDNzf9l2j69OmYOXOmVdvLKhe1G/RSoFarsWPZMvxZiYf4OQc6OAsaOAc6qGXhlMVaQkICAKBZs2Z600NCQhAWFoYDBw6Uuf4vv/wCAGjVqpVtGsjsSiaT2bsJ5LzzzjvoN20aCp97Dr/++qvBfFvcwYBzoIOzoIFzoINaFk5ZrF2+fBkAEBYWZjAvPDwcd+7cQV5eXqnrnz17FkBRZT1u3DjUqlULoaGh6NmzJ/bt22ebRrNKo1Qq7d0EchYvXgwAqA5g6NChlfKcnAMdnAUNnAMd1LJwymLt4cOHAAAfHx+DedppWVlZpa6fnp4OABg4cCAef/xxXLhwAYmJiQgICEDPnj3x008/lbquQqFAdna20R9GA7XubWoePXpUKc/DOdDBWdDAOdBBLQt3ezfAliQSw4M3xqaVlJ+fD6DogoK33noLAODr64s1a9Zg165dePPNNzFgwACj24qLi0NsbKzR7fr5+SE/Px/u7kUvu0qlgru7O4QQUKvV8PDwgEaj0T1Wq9XQaDTw9PSESqXSPS4sLIQQAjKZDEqlUvdYoVAAALy8vFBQUFDmY+3yEokEnp6eusceHh5QKpVwc3ODu7u77rFUKkVhYSGkUinc3Nx0jyUSiW4/HGWfNBoNFAqFU+1TRXMqSaFQGL10vaCgwGr7JITQZeEq7z2q+6TRaJCfn+9U++SIOWk0Gt3fhLPsk6PmBBTVArbYJ0tIBLXy0QratGmDEydOICcnB1WqVNGb17JlS5w+fRq5ubnw9fU1un7//v2xfft2zJs3D7NmzdKb17t3b+zatQuXLl1CgwYNDNZVKBS6AEtq0qQJUlJSLNwrZi35+fnw9va2dzNIkUgkEAAuAWiE/32r1E4/B6AZrPttk3Ogg7OggXOgw5ZZREZGml0LOOVh0Pr16wMAUlNTDealpaUhPDy81EINAGrWrAkACAoKMpgXGhoKAMjIyDC6rkwmg1wuN/rDGGOMMWYupyzWOnXqBABITEzUm56RkYHU1FR07ty5zPXbtm0LoGhg3JLu3bsH4H9FG3M8Xl5e9m4CA+dACWdBA+dAB7UsnLJYGzp0KPz8/LBt2za96Vu3bgUAjB07Vm96enq63gB4AwcOhL+/v8HguXl5eTh27Bjq16+v671jjkd7DgEzZIshOkrDOdDBWdDAOdBBLQunLNaCgoKwcOFCbNiwAVu2bAEAXLhwATExMRg+fDh69uypW/bQoUOoXr06BgwYoJvm7++PJUuWYP/+/Vi6dCkKCwuRk5ODyZMnIycnB1988UWl7xNjjDHGXJNTFmsAMGnSJKxduxaxsbEIDQ1F3759MWHCBKxdu1ZvOblcjsDAQN15alrjxo1DfHw81q9fj9DQUNSqVQsZGRk4ePAgevToUZm7wqyMWve2I7BFjxvnQAdnQQPnQAe1LJzyalCqLLkChFkfX3FlSHvVZxKAhjC8GvQ8gMfBV4M6K86CBs6BDr4alDFG2rMAkJtr72Ywxhj7FxdrzOVQ696m5lcAv3ToYPPbrXAOdHAWNHAOdFDLgos15nKoXeVDiXaY5wtnz2LdunU2fS7OgQ7OggbOgQ5qWXCxxhgz6ttvv7V3ExhjjIGLNeaCZDKZvZvAwDlQwlnQwDnQQS0LLtaYyynt3q2scnEOdHAWNHAOdFDLgos15nIkksocp98xVcYrxDnQwVnQwDnQQS0LLtaYy/H09LR3Exg4B0o4Cxo4BzqoZcHFGnM51Lq3qSo+AK4tvmNyDnRwFjRwDnRQy4KLNeZyqHVvuyrOgQ7OggbOgQ5qWXCxxlyOh4eHvZvAwDlQwlnQwDnQQS0LLtaYy7H1yPzOoDK+U3IOdHAWNHAOdFDLgos15nLc3PhtTwHnQAdnQQPnQAe1LNwrsnJBQQFycnKgUqng6+sLuVxurXYxZjPu7hV62zMr4Rzo4Cxo4BzooJaFWa1JSUnBhg0bsHv3biQmJiIjI0NvvpeXF+rXr4+OHTti4MCB6Nmzp1Uby5g1KJVKeHt727sZLo9zoIOzoIFzoINaFiYVa/fv38eMGTOwfv16qNVqvUv6i8vPz8fZs2dx9uxZLF++HI0aNcKnn36K7t27W7XRjkShUJC7BNjVUevepqq0v3Nr4Rzo4Cxo4BzooJZFucVaYmIi+vbti9u3b6Nly5bo0aMHmjVrhsjISAQHB0Mmk0EqlUKhUCAnJwdpaWm4dOkSDh8+jN27d6NXr16IjY3F7NmzK2N/yImLi0NsbCwAwM/Pz86tYQAglUrt3QTyKuMCA86BDs6CBs6BDmpZlFms3b59Gz179sRjjz2G9evXo1OnTiZttH///gCAe/fuYf78+Zg7dy6qVq2K1157reItdjCzZs3C9OnTAQBNmjSxc2sYABQWFpI7H4E6WxRvnAMdnAUNnAMd1LIosyVTpkxB+/btsXnzZosaHRwcjMWLF6N9+/YYPXo0+vfvjxo1aljcWEckk8kgk8ns3QxWDLVvTPaUnJyMz8eMQX07PDfnQAdnQQPnQAe1LEqtwIQQmD59Otq1a1fh6jIqKgoNGzYkVaUy10XtXAR7Gjt2LPYlJGCOHZ6bc6CDs6CBc6CDWhalVk8SicTkw56maNq0qdW2xVhFUOvetqf9+/cDAOwxVjfnQAdnQQPnQAe1LGiVjoxVAmrd21TZ+mpQzoEOzoIGzoEOallUuFjbsmWLNdrBWKWhdoNeiirjFeIc6OAsaOAc6KCWRYWLtUmTJlmjHYxVGpVKZe8mOCyNRgNIJEibOLHC2+Ic6OAsaOAc6KCWRbkHZOfNm1fm/Pz8fKs1hrHKQOk8BCpKfocs7QDojBkz8F8AV776CtdfeQVPP/20xc/JOdDBWdDAOdBBLYtyWxMXF4e8vLxSz1+h1lXIGLOdjz/+GP9FUTE3ZswYXLx40d5NYowxp1fuYdA333wT77zzDjQajdGfoKCgymgnY1ZDrXubqvK+iF29erVC2+cc6OAsaOAc6KCWRbnF2ttvv40NGzYY3LRdy9UGuWWOj1r3NkUS6F8Naov+c86BDs6CBs6BDmpZlFus+fn54cCBA/DwMD4S08mTJ63eKMZsydZDUjDTcA50cBY0cA50UMvCpNKxZs2atm4HY5VGrVbbuwlOoaIfZpwDHZwFDZwDHdSy4EFxmcsprZfYldnjMiHOgQ7OggbOgQ5qWZhdrHXr1s0W7WCs0mg0Gns3gYFzoISzoIFzoINaFmYXawcOHLBFO2xi06ZNaNGiBUJDQ1GnTh3ExsZadIXH7du34e/vz8OUOAlq3dsUlfdOt8ZfAudAB2dBA+dAB7UsaF3uYEUrV67E5MmTsXnzZkRFReH8+fPo0aMHkpKSsH79erO2NXXqVGRnZ9uopayyUevepsrWJ9hyDnRwFjRwDnRQy8Ipz1l78OABZsyYgZEjRyIqKgoA0LRpU8TGxmLDhg3Ys2ePydv66aefcObMGbRp08ZWzWWVjNo3JgrMLcsE+AIDZ8JZ0MA50EEtC6cs1jZv3ozs7GwMGjRIb/rgwYMBAKtXrzZpOzk5OZg6dSq++OIL+Pj4WLuZzE6onYvgqjgHOjgLGjgHOqhl4ZTFWkJCAgCgWbNmetNDQkIQFhZm8nl3s2bNQseOHdGnTx+rt5HZj6enp72bQE55r4gtztbkHOjgLGjgHOigloVTnrN2+fJlAEBYWJjBvPDwcJw+fRp5eXll9pYdPXoUmzZtwvnz523WTmYfKpUKUqnU3s0grTIuMOAc6OAsaOAc6KCWhVMWaw8fPgQAo8WYdlpWVlapxZpKpcKECRMQFxeHatWqmfXcCoUCCoXCvAazSkWte9sRSQBU02igVCot/gbKOdDBWdDAOdBBLQuzizVqt2Aoi7GhNkwZfmPRokWQy+UYP3682c8ZFxeH2NhYo/P8/PyQn5+vu+eYSqWCu7s7hBBQq9Xw8PCARqPRPVar1dBoNPD09IRKpdI9LiwshBACMpkMSqVS91hbJHp5eaGgoKDMx9rlJRIJPD09dY89PDygVCrh5uYGd3d33WOpVIrCwkJIpVK4ubnpHkskEt1+OMo+qdVqKBQKp9onS3MqS/ETbEv+1XcAcAfAwKFDsXXrVov2SaPRQKPRuNR7j+o+qdVq5OfnO9U+OWJO2uWcaZ8cNSchBPLz822yT5aQCEeqvkzUpk0bnDhxAjk5OahSpYrevJYtW+L06dPIzc2Fr6+vwbpXr15Fq1atcPjwYTRp0kQ3vUuXLjhw4EC5xWpZPWtNmjRBSkqKBXvErKmgoABeXl72bgYJEonE6JWgHwP4oV07/PXXXxAALgJo/O+84suXvOG7OTgHOjgLGjgHOmyZRWRkpNm1gFNeYFC/fn0AQGpqqsG8tLQ0hIeHGy3UAOD333+HEALdunVDWFiY7ufw4cMAoPt92bJlRteXyWSQy+VGfxgNTvj9xCFxDnRwFjRwDnRQy8Jmxdq9e/dQt25dW22+TJ06dQIAJCYm6k3PyMhAamoqOnfuXOq6r732GrKyspCWlqb30759ewDQ/T516lTb7QCzKZlMZu8mMHAOlHAWNHAOdFDLokLFWl5eHm7fvo2bN28a/Fy/fh03btywVjvNMnToUPj5+WHbtm1607du3QoAGDt2rN709PR0FBYWVlbzmJ0plUp7N4G8yrixGudAB2dBA+dAB7UsLCrWPv30UzRs2BB+fn6oWbMm6tSpY/DTtm1ba7fVZEFBQVi4cCE2bNiALVu2AAAuXLiAmJgYDB8+HD179tQte+jQIVSvXh0DBgywV3NZJaPWve2qOAc6OAsaOAc6qGVh9tWgCxYswKxZswAAcrkc/v7+cHMzrPnUajVu375d8RZaaNKkSZDL5YiNjcWkSZPg4+ODCRMmYPbs2XrLyeVyBAYGombNmka388wzz+Dy5cvIzMwE8L+x2y5cuICqVavadieYTVDr3qbK1h9WnAMdnAUNnAMd1LIw+2rQunXrIjQ0FGvXrkWDBg1KXS4tLQ3Vq1cnN1aJPVlyBQizvvz8fHh7e9u7GSSUdjXoJwC+f/ppHDlyBALAJQCN/p1nratBOQc6OAsaOAc6bJmFJbWA2T1rqamp2LRpU5mFGlA0+Ozo0aPN3TxjjABaBwAYY8y1mX3OWoMGDRAUFFTucnK5HF9//bVFjWLMlngco/JVxgUGnAMdnAUNnAMd1LIwu1h77733TCrC7Dl0B2Nl0Y4ozeyLc6CDs6CBc6CDWhZmHwYdNmwYlEolxo4di759+6JevXqoWrWqwS0U7t69a7ehOxhjFWfpbVEYY4xZl0U3cs/MzMTvv/+Ob7/91trtYczmqHVvuyrOgQ7OggbOgQ5qWZhdrH3++ed48803AQABAQGQy+VGv4Hbe+gOxkpTUFDAV1yZwNZDd3AOdHAWNHAOdFDLwuxi7bPPPkP79u2xdu3aMs9J0w7dwRhzPCW/fvHVoYwxZj9mF2s3b97Ed999V+7FAzx0B6OKWve2q+Ic6OAsaOAc6KCWhdlXg9apU8ekkft56A5GFbWrfByBLS414Bzo4Cxo4BzooJaF2cVadHQ01qxZU+5yPHQHY4wxxljFmX0YtEePHjh+/DgGDhyIUaNG8dAdzOFQu+eboxNCWDTMB+dAB2dBA+dAB7UszC7Wateurftg3r59u9UbxJitKRQKUlf52ENaWhoWzJwJ31LmV+Sen6biHOjgLGjgHOigloVF46xFRERAKpWWuQwP3cGo4sFegWnTpuGH+Hi8YMc2cA50cBY0cA50UMvComLtxIkTCA0NLXMZHrqDUeXp6WnvJthdfHw8AKCxHdvAOdDBWdDAOdBBLQuzLzDo3r27ScdyeegORpVCobB3Exg4B0o4Cxo4BzqoZWF2z9ru3btNWo6H7iiiUCjIhe7qqHVvU2TO2WqWXmDAOdDBWdDAOdBBLQuze9YAICcnB9nZ2cjJyTGYd/z4cRQWFla4Yc4iLi4O/v7+8Pf3R3Z2tr2bwwB4eHjYuwlklPZx1B1A2tGjNn1uzoEOzoIGzoEOalmYXaz9/fffCAgIQGBgICIiIgzmv/vuu6hevTp+++03qzTQ0c2aNQtZWVnIysqCXC63d3MYAKVSae8mkPc4gG3Ffi/tqlHA8qtGOQc6OAsaOAc6qGVhdrG2bt06CCEwfvx4o4dE582bh65du2LAgAE4deqUVRrpyGQyGeRyORdqhLi5WdSh7HKK97pFlrJMWwAdnngCmZmZZm+fc6CDs6CBc6CDWhYSYebX4tatW6Nz585YsmRJmcuNHz8e9+/fx9atWyvUQGcSGRmJlJQUezfD5anV6nKHnnF2EokEAkA2gNK+RpwF8AT+d/6atngr+YHxKYBbM2Zg0aJFZrWBc6CDs6CBc6DDlllYUguYXTpevnwZL730UrnLTZ48GSdPnjR384zZHLXubUfnC8sGyOYc6OAsaOAc6KCWhdnFmkKhQFBQULnLhYSEICMjw6JGMWZL1Lq3HZ2l10xxDnRwFjRwDnRQy8Ls1kRERODIkSPlLnfkyBEeFJeRxIcZTGNOEWbJRQacAx2cBQ2cAx3UsjC7WHv22Wfx1ltvlXmI8++//8aMGTPQp0+fCjWOMVvgoWVo4Bzo4Cxo4BzooJaF2YPizpo1Cxs3bkTbtm3Rrl07PPXUU7pbT2VkZOD48eM4dOgQ/P39MWvWLKs3mLGKovaNydFZetN3zoEOzoIGzoEOalmYXaxFRkZix44diIqKwqFDh3D48GG9+UIIhIWFYcuWLUbHYWPM3qidi0CVscOgflbcPudAB2dBA+dAB7UsLLqRe7t27XDp0iWsXLkSu3fvxs2bNwEAtWrVQu/evTF27FgeV4yRVVhYCHd3i976Ls/YPTgsvcCAc6CDs6CBc6CDWhYWt0QulyM6OhrR0dHWbA9jNkete9sZ8GFQx8ZZ0MA50EEtC1r9fIxVAmo36KXK1q8S50AHZ0ED50AHtSxKLdaEEFi+fDnu379vlSfavHkzLl68aJVtMVYRKpXK3k0gwxofR5ZeYMA50MFZ0MA50EEti1KLNYlEggsXLqB79+4W3fevuDlz5mDatGl8wQEjgdJ5CK6Mc6CDs6CBc6CDWhZlHgZduHAh3Nzc0LhxYyxZsgR5eXlmbfzAgQPo2rUr4uLisH79evj5WfNaMsaYLZna60brYAFjjDmfMktHb29v7NmzB4MHD8bMmTMxe/ZstGnTBk2bNkVkZCSqVq0KmUwGqVQKpVKJ7OxspKen48qVKzh27BhSU1Ph5+eHrVu3okePHpW1T3o2bdqE+fPn486dO/D19cWYMWPw3nvvlVs1P3r0CN999x2+++47XLx4ESqVCjVq1MCoUaMwffp0eHh4VNIeMGtTqVScn5VZehiUc6CBs6CBc6CDXBbCRF999ZWoW7eukEgkQiKRCDc3t1J/JBKJ8PLyEmPHjhV37twx9Sms7quvvhJSqVT8+OOPQgghzp07J8LCwsTIkSPLXXfEiBHC3d1drFu3TqhUKlFYWChWrVol3NzcRP/+/S1qT0REhEXrMetSKpX2boLdARACENn//mvs5zygW06UeFz851tA1KtXz+w2cA50cBY0cA502DILS2oBiRCmfyUWQuDw4cPYs2cPEhMTcePGDeTk5EClUsHHxwdhYWFo0KABOnTogN69eyMwMNBmRWZ5Hjx4gNq1a2PgwIFYu3atbvpXX32FiRMn4o8//kD37t1LXX/48OHw9/fHl19+qTd95MiR2LhxI/bu3YuuXbua1abIyEikpKSYtyPM6pRKJTw9Pe3dDLuSSCQQAHJQ+kC3FwA0RVGVBvx7IYGR5b4DMKduXVy9etWsNnAOdHAWNHAOdNgyC0tqAbPOoJNIJOjQoQM6dOhg1pPYw+bNm5GdnY1BgwbpTR88eDAmTpyI1atXl1msDRw4EPXr1zeY/vTTT2Pjxo04fvy42cUao0GtVtu7CQycAyWcBQ2cAx3UsqB1uYMVJSQkAACaNWumNz0kJARhYWE4cOBAmeuPGDHC6HTt5bwBAQEVbySzC1LnIdiZPS/54Rzo4Cxo4BzooJaF0w6Ke/nyZQBAWFiYwbzw8HDcuXPH7KtbAeD48ePw8PBAv379KtxGZh8ajcbeTXAI5lwNasbZFDqcAx2cBQ2cAx3UsnDanrWHDx8CAHx8fAzmaadlZWUZnV+alJQU/PTTT2WOGadQKKBQKMxvMKs01Lq3XRXnQAdnQQPnQAe1LJy2WNMydssIS24jIYTAlClT0KhRI3z00UelLhcXF4fY2Fij8/z8/JCfn68bNkSlUsHd3R1CCKjVanh4eECj0egeq9VqaDQaeHp6QqVS6R4XFhZCCAGZTAalUql7rC0Svby8UFBQUOZj7fISiQSenp66xx4eHlAqlXBzc4O7u7vusVQqRWFhIaRSKdzc3HSPJRKJbj8cZZ/UajUUCoVT7ZO5OZnCnJ41jUaD/Px8s/ZJ+7q50nuP6j6p1WqD/Bx9nxwxp+J/E86yT46ak/YzzRb7ZAmzrgZ1JG3atMGJEyeQk5ODKlWq6M1r2bIlTp8+jdzcXPj6+pq0vffeew/r1q3DX3/9herVq5e6XFk9a02aNOGrQQlQKBSQyWT2boZdaa8GLctFAI1R/tWg6wHMrl0bycnJZrWBc6CDs6CBc6DDlllYcjWo056zpr2SMzU11WBeWloawsPDTS7Uli5dijVr1uCPP/4os1ADiqpnuVxu9IfRQO1cBFfFOdDBWdDAOdBBLYsyi7X4+HiMHTsW3377bWW1x2o6deoEAEhMTNSbnpGRgdTUVHTu3Nmk7axZswZz587Frl27dAXg/fv3cf36dau2l1UeHsfINLa+wIBzoIOzoIFzoINaFmUWa59++in27t1rcg8UJUOHDoWfnx+2bdumN33r1q0AgLFjx+pNT09PR2Fhod60zZs3Izo6Gr/++qveECC//PIL5syZY5N2M9vTDr/C7ItzoIOzoIFzoINaFmWebfzPP/9g27ZteoPgSqVSpKamIjQ0tMwNFxQUYPPmzXj55Zet01IzBQUFYeHChZg6dSoGDRqEqKgoXLhwATExMRg+fDh69uypW/bQoUPo1KkTevXqhV9//RUAsGPHDrz00kvo1q0bduzYgR07duiWP3PmDI+z5sCodW87Oktv5M450MFZ0MA50EEtizKLtby8PISHh+tNM/VwR1ZWFl555RW7FWsAMGnSJMjlcsTGxmLSpEnw8fHBhAkTMHv2bL3l5HI5AgMDUbNmTd20RYsWobCwEL///jt+//13g22PHj3a5u1ntkGte5sqc4owPgzq2DgLGjgHOqhlUWaxVrt2bWzatAnvvvuu3nRLLz21h5EjR2LkyJFlLtOsWTPcu3dPb9r+/ftt2CpmT9rLupl9cQ50cBY0cA50UMuizGJt+PDhmD17Nr766ivUrFlTN1bJ4MGDy606lUql9VrJmBU56Wg1dmPpBQacAx2cBQ2cAx3UsiizWJs5cyauXLmC7777Djdv3tRNP3z4sEkbd6QeOOY6eBwj05jz13vr1i0olUqzDh1wDnRwFjRwDnRQy8KkQXEfPnyIpKQk5OXloXv37oiPj0fVqlXLXOf+/fsYOnQouVs22JMlA+Ex6ysoKICXl5e9m2FXpgyKexlAA5Q/KO5VAAsAeLz2Gj7//HOT28A50MFZ0MA50GHLLCypBcos1vbv34+EhAS0bNkSzz33HADAzc0NaWlp5V4Nmp6ejvDwcHJXVNgTF2s05Ofnw9vb297NsCtrFmsoPt+MQwecAx2OmsXNmzexdetW9OrVC40bN7Z3cyrMUXNwRrbMwpJaoMzDoPPmzcOJEycQFxenm5acnIyQkJByNxwSEmL27WcYqwzUurcr240bN+zdBACcAyWOmIVGo8EzHTvixq1bqOvjg7Pp6Qa3FnQ0jpiDs6KWRZmD4p49exYbN27ExIkTddO6du2KzMzMcjesVqvJ/KfAWHGl3bvVFSxduhT9a9dGgL0bAtfOgRpHzOLUqVO4eesWAKBaXh7WrVtn5xZVnCPm4KyoZVFmz1p2djYaNmyoN+369esmnYeWmZmJrl278jlrjBHy5ptvQgD4yd4NYayCit9xRoKicUEZc1Zl9qzVqFEDe/fuNZhu6lWe1C59ZQwAn8ALmNSzZutruTkHOhwxC2ccbcARc3BW1LIos2etb9++mDBhAr799lvUqlVLN87a66+/Xu6Jd/n5+U75x8QcX0FBAZ/EayJbflxxDnQ4ehaWjvVHjaPn4EyoZVFmsRYTE4MzZ87gwIEDOHjwoG76999/b9LGuVhjzLGVPUAPY/ZT8v8XZyjWGCtNmcVaQEAA9u3bhzNnzujGWRs3bhw+/vhj+Pv7l7nhhw8fYvr06VZtLGPWQK172x5M+RolMXE5S3EOdDhiFhKJxOaH6iubI+bgrKhlUWaxptWiRQu0aNECADB27FgMHz7cpHHW3nrrrQo3kDFro9a97ao4BzocPQtnKdocPQdnQi2LMi8wMGbfvn3l3r0AAKpWrYp9+/ZZ1CjGGGOMMVbEpJ614jp37mzSch4eHiYvy1hlota9TRVfDeo6HD0LZ7nAwNFzcCbUsjC7Z03r9u3bmDt3Lp599lk0b94ct/4dnPC///0vLly4YLUGMmZtBQUF9m4CA+dAiSNmUfwCA2cp1hwxB2dFLQuze9YAYOPGjZgwYQLy8vIghIBEItENULhu3Tq8/fbbWLJkCd544w2rNtYRKRQKciMhM2YKZzkPiDknHm2AuRKze9aOHz+Ol19+GRqNBsOHD8e7774LT09P3fyEhATExMRgxowZ2L9/vzXb6pDi4uLg7+8Pf39/ZGdn27s5DPTu+WYPnezdAHAOlDhDFs7Qs+YMOTgLalmYXawtWrQIjz32GJKTk7F+/Xp89NFHusFyAaBKlSp4//338eabb2LJkiVWbawjmjVrFrKyspCVlQW5XG7v5jDQu+cbZbbsu+Ac6HDELIwN3bF371688847uHTpkl3aVFGOmIOzopaF2cXan3/+icWLF5c7dMeLL76Is2fPWtwwZyGTySCXy7lQI4QPn5jG1q8S50CHo2chQdF51KO7d8e4BQvQqFEjTJ061eHuTe3oOTgTalmYXazdv38fjRs3Lne5qlWr4u7duxY1ijFbKn7YntkP50CHo2chAfDZZ5/haQD1AQgA5z7/HNu2bbNru8zl6Dk4E2pZmF2sBQUFmXS157FjxxAcHGxRoxizJWrd266Kc6DDEbMw1vNR/Ky1pwB89913ldYea3DEHJwVtSzMLta6dOmCt956C7dv3y51mVu3bmHWrFno3r17hRrHmC1Q6952JuZc7s450EE9i+PHj+Ojjz7CtWvXylyueLEmAXDy5EmbtsvaqOfgSqhlYXax9u677+LWrVuoX78+Ro8ejWXLlkGj0WDHjh1Yvnw5XnnlFTRp0gS3b9/GO++8Y4s2M1YhHh4e9m6CU3ocRV/mTL0qj3Ogg3IWmZmZ6NKhA9Lefx/16tWDXCLBzJkzy+1ZawkgPyUFf/75Z6W1taIo5+BqqGUhERZc7/zjjz/i5ZdfRn5+vsEfjBACvr6+WLduHQYOHGi1hjqDyMhIpKSk2Gz7BQUFyMzMRPXq1W32HM4gPz+f1D3fKpNEIoGpf/DXUTTEx03tukC56/oB2HP0KJ566qlyt+/KOVBDOYtFixZh19tvYzeAoQB+QNF7cefOnejfty9UALoC2A9gEICtxdbdDWBC7dpITk6u3EZbiHIOrsaWWVhSC1h0B4Pnn38e586dw5QpU1C/fn3IZDLIZDI0aNAAr7/+OhITE7lQq2Q5OTlo8cQTaB4RgejoaHs3hzQ3N4tv3MHK4Y6ii5BMwTnQQTkLhUIBzb+Pfyg2/c6dO7rH2i6Dkl8m/AFcv37dZm2zNso5uBpqWVh0BwMAqFOnDj777DNrtoVVwNKlS9E8KQmbAUj++18e464MxccFdBVZWVn48ccfbf485tz2xxVzoIpqFr/99hvmzZuHp01cvuQ7j9Z/t+WjmoMropYFrdYwix0/fhyN7N0IB6FUKl3uUEP//v3hc/AgzDkLQ8D8Ox24AdBoNOUuB7hmDlRRzOLvv//G8D59MABAmpH5xr4UlJxC6xTx8lHMwVVRy8LiYk0IgZ9//hm7du3CjRs3IJFIULNmTfTq1QsDBgwgdyWFs1Or1ZDauxEOglr3tq2p1WocPHgQAsBrZqwnATDBzOcyp1hztRwoo5jFO++8g1kA/g/AeyXmlfa/i6P3rFHMwVVRy8KiYu3q1asYMmSI7g4F2m84EokEK1aswBNPPIH4+HjUrVvXei1lZeJizXRSqWu9UsV7IMy9vsncV8qcw6CulgNlFLPIy8vTPf5PiXkSAHfv3tUVbRIAYQB+KbGcDEUF27lz5/D444/bqKXWQzEHV0UtC7OLtczMTHTv3h03b95EZGQkWrdujZCQEABARkYGjh8/jjNnzqB79+44deoUAgMDrd5oZoiLNdMVFhaSOx/BlooXT+b0dwuY3zNhTs+aq+VAGcUsyrpy2Q3Ae++9h6naZQE8Y2S5JgC+BTBu3DgcPXrUBq20Loo5uCpqWZjdkoULF+pOVh48eLDRZbZs2YJx48Zh4cKFiIuLq3AjWfk0Go3DdfnbC7VvTLZWvFibWsZyJUlgWc+aqcWaq+VAGaUskpKS8PLLL+Po0aPoXcoyhQBqoWhsP63S3nVdAIw6dsyKLbQdSjm4OmpZmP3/+08//YTFixeXWqgBQFRUFBYtWuRw92VzZNyzZjpq5yLYWvHi6TEz1zX325w5h0FdLQfKKGUxduxYjD56FJ+j7HH9WkC/p7i0Yo3OnpWPUg6ujloWZrfmxo0b6NGjR7nL9ezZEzdu3LCoUdayadMmtGjRAqGhoahTpw5iY2OhUqkqbf3KxMWa6QoLC+3dhEplwbjXRevB/J41NxSdKmEKV8uBMkpZHDp0CJNRdDFMWe9ctxLzSyvWtEOEW/p3UJko5eDqqGVhdrEmk8mQlZVV7nJZWVmQyWQWNcoaVq5ciZdeegkffPAB7t69i+3bt2PFihUYPXp0paxf2UoWa6+99hokEgmioqLs1iaqqHVv25JarUZ6erpF60oAPGHmOoEAJk6cqDdgaWlcKQfqqGZRVs9u8TtqSFB6sQYA4QB+//13azXLZqjm4IrIZSHM1KVLFzF+/Phylxs3bpzo3LmzuZu3iszMTCGXy8WoUaP0pn/55ZcCgPjjjz9sun5pIiIiLFrPFE899ZT4LyBE0eeXkP77OAQQhw8fttnzOiKlUmnvJlQKtVotnnnmGSH59z0hzPxJLvG7qduYBohx48aV2z5XycERUMhi9+7d4oUXXtB7n10o430WVexxT0D0K2PZXwERGhpq710sF4UcWBFbZmFJLWB2z9q4ceOwevVqREVFYe/evcjJydHNy8nJwd69ezF48GB8/fXXePXVV61UUppn8+bNyM7OxqBBg/Sma8+zW716tU3Xt4eSPWvN/v23F4Bp06bZoUV0UT2UbW07duzAnYMHy+xxsIW2AL799ttyl3OVHByBvbMoKChAz5490eH77+FVbHpZPWsl//Mq6yBnBIqG+qDO3jmw/6GWhdlXg7700kvYvn07Nm/ejJ9++gkA4OnpCaBoxF+g6NyAYcOG4cUXX7RiU02XkJAAAGjWrJne9JCQEISFheHAgQM2Xd8eFAqFXrGmfZt5Azh58iQkEgmeeOIJHD161K6HpymgdDm2rezduxcjR45Ei39/H1mJzx2Eog+6goICeHl5lbqcK+TgKOyVRXZ2Nt544w188803cAMwDUVXemqV1ZtQ8iBVWV9KmqHoUOn9+/cRFBRkWWMrAf9N0EEtC4tas3HjRrRo0QKLFy9GZmYmFAqFbl7VqlUxc+ZMvP3221ZrpLkuX74MAAgLCzOYFx4ejtOnTyMvLw8+Pj42Wb+yXbp0CefOndM7v6h4sQYAXwBY+/ffev95Nm/eHB07drRZu0SJE3qL/27uY2uur9Fo9O6wUdnPb6v1c3NzsW/fPgBFN7AufuVnV5gvvNjjh2as9yyAJwF4e3vj+eefR2hoqNHlNBqN0Suuytpvc6ZZuh6FNlT2eiX/JmzZBo1Gg19++d/wte0BDAPw07+/Ty+2XlnFWvH/YdQofwzBZwEEBwdj4sSJNr/Sz5LXSQih9zdhyvKWPEdFfnfW59BoNLh9+zYyMzNx/fp1AEDLli1x8uRJg/XsRSKM7YWJCgsLcfLkSdy8eRMAUKtWLbRq1cruFWmDBg1w+fJlqFQqg5MEO3bsiEOHDuHOnTsIDw+3+voKhUKveC2uSZMmSElJsXCvSvfYY4+hzdWr2Pjv7x8AqAPgFQB/A9iO/92uJQZFH1p1AXxi9ZboK/nhWdbv5ixrr3Ud4XkCAdQE0BnW91/o/0danq0Ayhvdyth/sKZM4/XotsGUZcIBPA9gCYDof6d9hqKeNUv8hqL3fttylvsYgC0OhlryGvA6dNbR/p6JovshN0fRl4j/W7MGr7zyipFnqJjIyEiza4EKFWtUaYsttVpt8A3qmWeewZ9//mlSsWbJ+nPmzEFsbKzR7bq5uZX6nI4gOzsbcrnc3s2oMGfYD94HGpxhHwDn2A/eBzqcYT9suQ8ZGRmlduqUxqIusJycHAghIJFI4Ofnpzfv+PHjaNGiBTw8zL0LofX4+/sDKLq3XJUqVfTmPXr0CADKDKEi68+aNQvTpxvvf5DJZA57vlh2djb8/f2RkpLi0H+EzrAfvA80OMM+AM6xH7wPdDjDflDcB7MP3P/9998ICAhAYGAgIiIiDOa/++67qF69On777TerNNAS9evXBwCkpqYazEtLS0N4eDh8fX1tsr5MJoNcLjf646iFGmOMMcbsx+xibd26dRBCYPz48di9e7fB/Hnz5qFr164YMGAATp06ZZVGmqtTp04AgMTERL3pGRkZSE1NRefOZZ/RU9H1nZFMJkNMTAwXnARwFjRwDnRwFnRwFrZh9jlrrVu3RufOnbFkyZIylxs/fjzu37+PrVu3VqiBlrh//z7q1KmDQYMGYe3atbrpX331FSZOnIhdu3ahZ8+euunp6emoWrWq7tCtueszx6Ht3s7KyiLTve2KOAc6OAsaOAc6KGZhds/a5cuX8dJLL5W73OTJk+122WtQUBAWLlyIDRs2YMuWLQCACxcuICYmBsOHD9crtA4dOoTq1atjwIABFq3PGGOMMWZLZl9goFAoTBpUMCQkBBkZGRY1yhomTZoEuVyO2NhYTJo0CT4+PpgwYQJmz56tt5xcLkdgYCBq1qxp0fqMMcYYY7ZkdrEWERGBI0eOGBQ3JR05cgTVq1e3uGHWMHLkSIwcWfbY7c2aNcO9e/csXp85Fj6fggbOgQ7OggbOgQ6KWZh9ztqUKVOwbds2/Pzzz2jVqpXRZf7++2/0798fAwYMwLJly6zSUGvZv38/unb933juo0ePxjfffGPy+vn5+Vi8eDE2btyI69evw8fHB40bN8a0adMwbNgwqFQqrFixAitXrsS1a9fg7u6Oxx57DK+88gpee+01G+wRY4wxxpyZ2eeszZo1C/n5+Wjbti2eeeYZREdHY8GCBViwYAFmzJiBzp07o1WrVnj06BFmzZplizZXSPv27ZGamopPPvnE6Py1a9fC398f3333ndH5kyZNwgcffIARI0bg7Nmz2Lp1K27evImdO3cCAObOnYtp06ahQ4cOOHHiBP744w+o1Wps3rzZVrvEGGOMMSdm0R0M/vrrL0RFRSE9Pd3o/eTCwsKwZcsWPP3001ZrqLV98803eOWVVwx61vr164edO3eiX79+2L59u946Dx48QNWqVVGjRg3dLbYAYPPmzXj06BFeeeUVBAQEIC8vD1lZWfD2Lroz5969e3Hy5EnMnDmzUvaNMcYYY87DojsYtGvXDpcuXcLKlSuxe/duvXuD9u7dG2PHjiVzuau5pk2bhoyMDLz++usG85KSkgDA4Hy9YcOGASgaAiQrKwsRERG6Qg0AunXrhm7dutmw1YwxxhhzVhbfcV0ulyM6OhrR0dHlL+xAnn32WTz77LNG5+Xn5wNAqTeqL28+Y4wxxpi5zD5nzdHs2rULnTt3RpUqVRAQEIAePXrg4MGDBstdv34dEolE99OlSxeDedoLEw4cOKC3LABIJBLUqVMHAHDjxg29+devX7f5fjLGGGPMOZldrKlUKiQkJCAhIUHv3pl37tzBqFGj8Pjjj6Nfv344evSoVRtqiW+++QbPPvssbt68iY0bN+Lvv//Gu+++i+joaOzdu1dv2Ro1aiA1NRU//vijwXZKzmvXrh1SU1N1P0DRfUSPHTsGAIiMjNSbX6NGDRvvKWOMMcacldnH637//Xf0798fbm5uWLVqFcaMGYPCwkL06dMH586dgxACFy5cwP79+3HmzBndTdErW0pKCqZMmQKJRIKff/4ZzZo1A1B0Xt327dtRu3ZtveWlUinCwsJQtWpVg22VnOfp6YmwsDC9ZcLCwlBQUKC3PGOMMcZYRZnds/bDDz/g8ccfx/Xr1zFmzBgAQHx8PBITE9GoUSMcPHgQe/bsQURERLn3D7Wl1atXIy8vDx07dtQValqhoaEYMmSInVrGGGOMMWY6s3vWTpw4gQULFiAyMlI37dtvv4VEIsGnn36KDh06AADi4uLsemumhIQEAEU3njemUaNGldkcxhhjjDGLmN2zdv36dTz++OO633NycrB//37UqFED3bt3101/8sknkZKSYp1WWuD27dsAgODgYKPz+TAlY4wxxhyB2cWan58fsrOzdb9v3boVSqUSzz//vN5ySqWSxH21Sg7aW950xhhjjDFKzC7WGjZsqLsVU0FBARYvXgyJRIIXX3xRb7ljx47pHSqtbBEREQCAu3fvGp2fk5NTmc1hjDHGGLOI2cXalClTsHjxYjRv3hyPPfYYzp07h/bt26Nly5YAinrU9u3bh9mzZ6N9+/ZWb7CpOnfuDAA4efKk0fmJiYmV2RzGGGOMMYuYXawNHToU77zzDq5cuYI7d+6gefPmWLt2rW7+lClT0L17d6SkpGD48OFWbaw5xo0bB19fX/z55584d+6c3rz79+8bHU+NMcYYY4wai+5gMG/ePDx48AAZGRk4c+aMbuR+AFiwYAGSk5ORnJyMTp06Wa2h5oqIiMDy5cshhED//v3xyy+/4ObNm9i/fz969+6tu0giPz8faWlpyMrKQlpaGjIzMwEU9RBqp6vVaqPz0tLSdLeYSktLQ0ZGBgDolk9LS0Nubq4d9p4xxhhjzkIihBD2boQt/fHHH5g7dy5OnDgBiUSCZs2a4a233kJeXh5eeeUV3XIvvPACvv/+e4P1R48ejTlz5ugVpMV9/fXXGDNmTKkXLMTExGDOnDlW2RfGGGOMuZ4KFWv379/H4cOHkZKSAolEgoiICLRv3x5BQUHWbCNjjDHGmMsye1BcAHjw4AGio6Oxfv16qFQq/Q26u+PFF1/EkiVLEBgYaJVGMsYYY4y5KrN71h4+fIiOHTviwoULAACZTKa7Z2ZmZiYUCgUAoEmTJvjzzz8REBBg3RYzxhhjjLkQsy8wiImJwYULF/Dqq6/izJkzyMvLw+3bt3H79m3k5eXhzJkzGD9+PC5cuMDnajHGGGOMVZDZPWs1atTAqFGjMG/evDKXe+edd7B+/XrcunWrQg1kjDHGGHNlZhdr3t7e+Pvvv9GgQYMyl0tKSkLz5s1RUFBQoQYyxhhjjLkysw+DRkREmHRfTTc3N4SGhlrUKMYYY4wxVsTsYm3EiBHYtGlTuctt2rQJzz33nEWNYowxxhhjRcw+DKpQKDBo0CA0a9YMU6ZMQa1atfTmp6am4quvvsKWLVuwd+9eHnONMcYYY6wCyizW6tata3S6RqPRXTggk8kQGhoKiUSC3Nxc3S2ZQkND4ePjg6tXr9qg2YwxxhhjrqHMYs3NzaJbh/5v4xIJ1Gp1hbbh6BQKhW7sueDgYD6PjzHGGHNhGRkZurrAVGXewUAikeDYsWMIDg42uzF3797F008/bfZ6ziYuLg6xsbEAAD8/P6SkpNi5RSw/Px/e3t72bobL4xzo4Cxo4BzosGUWkZGRZq9Tbs9aWlqaRb1B6enpCA8Ph0ajMXtdZ1K8Z61JkyZcrBEghDDpimZmW5wDHZwFDZwDHbbMIjIy0uxaoMzjnA8ePLD4sF21atXw4MEDi9Z1JjKZDHK5HHK53N5NYf/isf9o4Bzo4Cxo4BzooJZFmcWav7+/xRvOzs7GW2+9ZfH6jDHGGGPMgnHWTJWfn49vv/3WVptnzGJeXl72bgID50AJZ0ED50AHtSzKvMCgNAUFBfjiiy+wY8cO3L59G0ql0mAZV78KlNFVUFDAJ/ESwDnQwVnQwDnQQS0Ls4u1/Px8dO7cGSdOnCh3WT5RkjHGGGOsYsw+DLpgwQJcvnwZ3333HZKSkuDj44MDBw4gOTkZycnJ+PPPP/H666/Dx8cHf/zxhy3ajE2bNqFFixYIDQ1FnTp1EBsbC5VKZbP1NRoNlixZAi8vL8yZM8cKe8DsSSaT2bsJDJwDJZwFDZwDHdSyMLtn7ccff8TcuXPx4osvAijqPYuIiNDddqpWrVpo37495HI5vv32W3Tt2tWqDV65ciUmT56MzZs3IyoqCufPn0ePHj2QlJSE9evXW33969evY8yYMUhOTjZ7EDtGk0KhINW97ao4Bzo4Cxo4BzqoZWF2z9q1a9fQrVs33e8SiQTGhmobPnw49uzZU7HWlfDgwQPMmDEDI0eORFRUFACgadOmiI2NxYYNG8p9PkvWHzhwIPr374+vv/7aqvvC7IcPz9PAOdDBWdDAOdBBLQuzizUPDw/4+Pjofvf29jY6uJtSqcTdu3cr1roSNm/ejOzsbAwaNEhv+uDBgwEAq1evtvr6O3bsQHR0dIVvvcXo8PT0tHcTGDgHSjgLGjgHOqhlYXYFUrt2bZw5c0bv97Vr1xos98UXX1h9INiEhAQAQLNmzfSmh4SEICwsDAcOHLD6+pbcFoLRxoezaeAc6OAsaOAc6KCWhdnnrD3zzDOYOXMmatSogVatWqFPnz6YO3cu0tLS0LNnTwDAzp07sWfPHt2hRmu5fPkyACAsLMxgXnh4OE6fPo28vDy9nj9rrs+cA7XubVfFOdDBWdDAOdBBLQuzi7Vhw4Zh5cqVGDFiBJKSkjB9+nSsW7cOv/76K3777TcARffUksvlmDt3rlUb+/DhQwAwWkxpp2VlZZVabFV0fVMUvxcoo8nDw8PeTWDgHCjhLGjgHOigloVFPWvF75kll8tx+PBhfPjhh0hISIBarcZTTz2FWbNmoUGDBlZtrJaxitecKrii65clLi4OsbGxRuf5+fkhPz8f7u5FL7tKpYK7uzuEEFCr1fDw8IBGo9E9VqvV0Gg08PT0hEql0j0uLCyEEAIymQxKpVL3WFskenl56TIq7bF2eYlEAk9PT91jDw8PKJVKuLm5wd3dXfdYKpWisLAQUqkUbm5uuscSiUS3H46yT9nZ2fDx8XGqfXLEnBQKBeRyuVPtk6PmlJ2dDS8vL6faJ0fMKT8/H/7+/k61T46aU25uLjw9PW2yT5aQCGOXchLVpk0bnDhxAjk5OahSpYrevJYtW+L06dPIzc2Fr6+v1dffv38/unbtipiYmDLHWiurZ61JkyZGL8ZglUuhUJAbQ8cVcQ50cBY0cA502DKLyMhIs2sBh7rEsX79+gCA1NRUg3lpaWkIDw8vtVCzxvqmkMlkkMvlRn8YDdpvWcy+OAc6OAsaOAc6qGXhUMVap06dAACJiYl60zMyMpCamorOnTvbdH3mHIzdy5ZVPs6BDs6CBs6BDmpZOFSxNnToUPj5+WHbtm1607du3QoAGDt2rN709PR0FBYWWrw+c048Zh4NnAMdnAUNnAMd1LKg1ZpyBAUFYeHChdiwYQO2bNkCALhw4QJiYmIwfPhw3dAhAHDo0CFUr14dAwYMsGh95rykUqm9m8DAOVDCWdDAOdBBLQuHKtYAYNKkSVi7di1iY2MRGhqKvn37YsKECQYD88rlcgQGBqJmzZoWra+1YMEChIWF6caMW7x4McLCwjBw4EDb7CCzueK9rcx+OAc6OAsaOAc6qGXhUFeDOjpLrgBh1qdUKsndSsQVcQ50cBY0cA502DILp78alDFroHYugqviHOjgLGjgHOiglgWt1jBWCah1b7sqzoEOzoIGzoEOallwscZcDrUTR10V50AHZ0ED50AHtSy4WGMuh9oNel0V50AHZ0ED50AHtSy4WLMxhUKB7OxsZGdn27sp7F8qlcreTWDgHCjhLGjgHOiglgUXazYWFxcHf39/+Pv7c8FGBLXbiLgqzoEOzoIGzoEOallwsWZjs2bNQlZWFrKysvj+oIwxxhgzGxdrNlb8xu6MBmrd266Kc6CDs6CBc6CDWhZcrDGXQ61721VxDnRwFjRwDnRQy4KLNeZy+KYdNHAOdHAWNHAOdFDLgos15nLUarW9m8DAOVDCWdDAOdBBLQsu1pjL8fDwsHcTGDgHSjgLGjgHOqhl4XDF2qZNm9CiRQuEhoaiTp06iI2NNetEQHPXr+jzMXo0Go29m8DAOVDCWdDAOdBBLQtaZ9CVY+XKlZg8eTI2b96MqKgonD9/Hj169EBSUhLWr19v9fUr+nyMJmrd266Kc6CDs6CBc6CDXBbCQWRmZgq5XC5GjRqlN/3LL78UAMQff/xh1fUr+nzGREREmL0Os77CwkJ7N4EJzoESzoIGzoEOW2ZhSS3gMIdBN2/ejOzsbAwaNEhv+uDBgwEAq1evtur6FX0+Rhe5b0wuinOgg7OggXOgg1oWDlOsJSQkAACaNWumNz0kJARhYWE4cOCAVdev6PMxuqidi0DdhAkTsGLFCqtvl3Ogg7OggXOgg1oWDlOsXb58GQAQFhZmMC88PBx37txBXl6e1dav6PMxujw9Pe3dBDKEELhw4QLOnz+PuXPnYuTIkZBIJFi+fDny8/Nx9uxZ3bmbAJCdnY1u3bqhQYMG6N+/P86ePQulUmnRc3MOdHAWNHAOdFDLwmEuMHj48CEAwMfHx2CedlpWVpbR+Zasb+nzKRQKKBQKo20QQpQ6r6Lmz5+P/fv3Y//+/RgzZgykUilWr16NF154AQEBAdi1axeaNm2KgIAAmzy/I9FoNHBzs/33FIlEovtX+1PW77aYp1Kp8Mknn2Do0KH44Ycf4OnpaVJx9dprr+G1117Tm9a2bVscO3ZM9/vly5exfft2g3W7d++O4OBgvTZp21VcZeXAymfPLDIyMiCXyxEQEIBVq1ahRo0aqF69Oo4ePYoBAwbg559/Rs2aNXHz5k299bTvaQBo3bo1GjRogA0bNpT6PFFRUZDJZLrfjb0nS04ztowt8d8EDW5ubpg+fToaNWpk76boOEyxpmXKH5g11zd3+bi4OMTGxhqdJ5PJsGTJEl33qpubG9RqNSQSidHHbm5uuv9wpVJpmY/Pnj2LgwcPAgCuXbumKwovXLiAwMBABAQEIDs7GwUFBQCKjsdLpVIIISCEgJubW4Ufa18bjUaj+w+6Io+1bSzeXlMel9deDw8PqFQqm+6TNj9tu0pOF0LoPphLPi6+H8Ufa7cPwOhjbduLP9YWZtqeYkt7wQCgRo0aesVaaa5cuQKpVGp0X7XtdXNzg1QqRWFhocH0ko+L75M2D0seF3+9yntOWz+mtk8eHh5Qq9V22SeVSoXbt28jIyMDACCVSnH06FEAwM8//wwABoUaAPz222+6x6mpqcjKyjJYprjMzEzd+67k37mxfbLH5567uzuUSqVNPvec8bPcVvskhMAPP/wAHx+fUv/PLSws1I3FZs5jbVvM5TDFmr+/PwAgLy8PVapU0Zv36NEjACjzZunmrm/p882aNQvTp0832obGjRsjOjq61DZW1IULFzBx4kTs2rXLZs/hDIp/KLiSadOmYdq0aWjevLlJy9esWRPx8fF46qmnAADr16/Hjz/+WOY6mZmZpfZul+SqOVBELYvMzEzdf5pvvfUWqlWrhk8//VRvmWPHjsHd3R3169fHihUrEBQUhPbt2xtsKzU1FYGBgZXV9AqhloMrs2UWy5cvN3sdhynW6tevjxMnTiA1NRX169fXm5eWlobw8HD4+vpabX1Ln08mk+l1tRcnkUhKnWcNTz75pEk9H66uoKDApjlQ9dVXX5m1vFQqRUREhO738l4zc78xumoOFFHLIjw8XPd406ZNEEJg9OjRePToETp16qTrBdFesefl5YV27dohPj4eQ4YMAQAMHDgQX331FUJDQ+2yD5agloMro5aFwxwc79SpEwAgMTFRb3pGRgZSU1PRuXNnq65f0edjdFnaDe0s3nnnHZOuZpZIJKhevbpJ2yz5hcYUrp4DJdSzkEgkaNmyJTp27Ijff/9d79xMALoekAYNGujWGTNmjEMVagD9HFwJtSwcplgbOnQo/Pz8sG3bNr3pW7duBQCMHTtWb3p6ejoKCwstXt/c5ZnjoPRtyR7i4uJ0X0bKYuqJzi+88AKSkpLMboer50CJo2QhkUjQq1cvvd8BwN296CBRs2bNIIRAbm6uwRiZjsBRcnAF1LJwmGItKCgICxcuxIYNG7BlyxYARedoxcTEYPjw4ejZs6du2UOHDqF69eoYMGCARetbsjxzHBU5yd6VmHrhjqXndXAOdDh6FlWrVtX7vaxTYihz9BycCbUsHOacNQCYNGkS5HI5YmNjMWnSJPj4+GDChAmYPXu23nJyuRyBgYGoWbOmRetbujxzDNS6t+3Jz88POTk5RueZ2rNm6VADnAMdjpqFRCLB7du3TT5cT52j5uCMqGUhEdRa5MQiIyORkpJi72a4PB7LqMi7776Lhg0bYsyYMUbnN2rUCP/884+uh634UBHFvfHGG/jkk0/Mfn7OgQ7OggbOgQ5bZmFJLcDvCuZybDUwsaOZN28eRo8eXep8Uz6oFi1ahHnz5ln0/JwDHZwFDZwDHdSycKjDoIyxymPKOWszZsyohJYwxphr45415nK8vLzs3QSHULxYe+yxxwAAp06dstr2OQc6OAsaOAc6qGXBxRpzOdpbbrEipd2V4IknntA91hZuTz75pNWel3Ogg7OggXOgg1oWXKwx5uKioqJw8uRJvWmjRo3CmjVr7NQixhhjxfE5azamUCjInajo6qh1b1PQsmVLvd8DAgLg6emp+93UMdfMwTnQwVnQwDnQQS0L7lmzsbi4OPj7+8Pf3x/Z2dn2bg4Dve5tR2CLYo1zoIOzoIFzoINaFlys2disWbOQlZWFrKwsyOVyezeHsTJVq1YNgG2KM8YYY5bhw6A2JpPJyN1jzNVR696mJDg4GOnp6QbFWsnfvb29dVeIWopzoIOzoIFzoINaFtyzxlwOte5tSrp06QKg/J61Fi1a4OzZsxV6Ls6BDs6CBs6BDmpZcLHGGAMAqFQqfPbZZwAMi7XyfmeMMWY7fBiUuRw+LG2cVCrVPS6vGLNGscY50MFZ0MA50EEtC4frWdu0aRNatGiB0NBQ1KlTB7GxsVCpVDbdhkajwZIlS+Dl5YU5c+ZUcA+YvfFQKuYrXpy1bNkSr7zySoW3yTnQwVnQwDnQQS0Lh+pZW7lyJSZPnozNmzcjKioK58+fR48ePZCUlIT169fbZBvXr1/HmDFjkJycTC48Zhk+hFcxJQfQtRTnQAdnQQPnQAe1LBymZ+3BgweYMWMGRo4ciaioKABA06ZNERsbiw0bNmDPnj022cbAgQPRv39/fP3119bdIWY3xQd7ZcZVxgcV50AHZ0ED50AHtSwcpljbvHkzsrOzMWjQIL3pgwcPBgCsXr3aJtvYsWMHoqOj4ebmMC8VKwf3kJavMi4o4Bzo4Cxo4BzooJaFw1QgCQkJAIBmzZrpTQ8JCUFYWBgOHDhgk21ERkZa2mRGFLXubYoq4zXiHOjgLGjgHOigloXDFGuXL18GAISFhRnMCw8Px507d5CXl2fzbTDH5+HhYe8mkNa5c2cMHz5cb5otPrg4Bzo4Cxo4BzqoZeEwFxg8fPgQAODj42MwTzstKyvL6HxrbqM8fON2+pRKJby9ve3dDLL2799fKc/DOdDBWdDAOdBBLYtKL9by8/NNOmSp9eyzz+r9buwbvrnf+q2xjdLExcUhNjbW6Dw/Pz/k5+fD3b3oZVepVHB3d4cQAmq1Gh4eHtBoNLrHarUaGo0Gnp6eUKlUuseFhYUQQkAmk0GpVOoea4tELy8v3ejLpT3WLi+RSODp6al77OHhAaVSCTc3N7i7u+seS6VSFBYWQiqVws3NTfdYIpHo9oP3yTn3SUutVlt1nwoLC/X2j3Oy3z4plUqn2ydHzEmpVDrdPjlyTvn5+TbZJ0tIhBDCojUtdP36ddSpU8fk5bXNa9OmDU6cOIGcnBxUqVJFb5mWLVvi9OnTyM3Nha+vb6nbqsg29u/fj65duyImJqbMsdbK6llr0qQJUlJSSl2XVQ61Wq03ACwrm0QiQfPmzfH3339bdbucAx2cBQ2cAx22zCIyMtLsWqDSe9YiIyORnJxs9nr169fHiRMnkJqaivr16+vNS0tLQ3h4eJmFmrW2UR6+cTt91Lq3XRXnQAdnQQPnQAe1LCq9WHN3d0ft2rXNXq9Tp07YuHEjEhMT9QqtjIwMpKamGpwQbattMMfHw7CYzxYXGHAOdHAWNHAOdFDLglZryjB06FD4+flh27ZtetO3bt0KABg7dqzBOunp6SgsLKzQNpjz4cMMNHAOdHAWNHAOdFDLwmGKtaCgICxcuBAbNmzAli1bAAAXLlxATEwMhg8fjp49e+otf+jQIVSvXh0DBgyweBvMORUv4Jn9cA50cBY0cA50UMvCYYo1AJg0aRLWrl2L2NhYhIaGom/fvpgwYQLWrl1rsKxcLkdgYCBq1qxp8TYAYMGCBQgLC9Pdnmrx4sUICwvDwIEDrb+DrFJQ+8bkqjgHOjgLGjgHOqhlUelXg7oyS64AYdZX/FJuVj6JRIIWLVrg9OnTVt0u50AHZ0ED50CHLbOwpBZwqJ41xqyBWve2I7DFBQacAx2cBQ2cAx3UsuBijbkcat3bropzoIOzoIFzoINaFlysMZdD7Qa9jsAWrxnnQAdnQQPnQAe1LLhYYy5HpVLZuwkMnAMlnAUNnAMd1LLgYo25HD6B13y2+JbJOdDBWdDAOdBBLQtarXFCZd0rlDHGGGOsPNyzZmNxcXHw9/eHv78/srOz7d0cBnrd246gevXqVt8m50AHZ0ED50AHtSx4nDUbK96zFhwcjNDQUDu3yHLZ2dmQy+X2bkaFOcN+8D7Q4Az7ADjHfvA+0OEM+2HLfcjIyDD7iBsXa8wk2dnZ8Pf3R1ZWlkP/ETrDfvA+0OAM+wA4x37wPtDhDPtBcR/4MChjjDHGGGFcrDGTyGQyxMTEQCaT2bspLo+zoIFzoIOzoIOzsA0+DMpcCsXubVfEOdDBWdDAOdBBMQvuWWOMMcYYI4yLNcYYY4wxwrhYYy6Fz6eggXOgg7OggXOgg2IWfM4aY4wxxhhh3LPGGGOMMUYYF2uMMcYYY4RxscZIe/ToEVasWIEOHTogKCgI/v7+ePzxx7FgwQIUFhYaXWfTpk1o0aIFQkNDUadOHcTGxpZ6nzdbLevsbt++DX9/f0gkklKX4Rxs6+HDh4iOjkb9+vURGhqKsLAw9OzZEzt27DBYlrOwjUePHuHDDz9E06ZNERoailq1amHYsGH4559/jC7POViHRqPBkiVL4OXlhTlz5pS5LJXXvMIZCcYIGzFihHB3dxfr1q0TKpVKFBYWilWrVgk3NzfRv39/g+W/+uorIZVKxY8//iiEEOLcuXMiLCxMjBw5stKWdQWDBg0SAERpHyGcg22lp6eLxx57TEyYMEE8ePBACCHEyZMnRVBQkHjjjTf0luUsbEOj0YhnnnlG+Pn5ib179wohhLh7967o1q2b8PX1FefOndNbnnOwjuTkZNG5c2dRs2ZNAUDExMSUuiyV19waGXGxxkh74YUXxIQJEwymjxgxQgDQfUgKIURmZqaQy+Vi1KhRest++eWXAoD4448/bL6sK9i2bZuoXbu2aNOmjdFijXOwveeff160atVKqNVqvelLliwRCxcu1P3OWdjOwYMHBQAxY8YMvemJiYkCgHj11Vd10zgH62nevLlYvHix2LNnT5nFGpXX3FoZcbHGSNuwYYM4fvy4wfSlS5cKAGLBggW6aStWrBAAdN9etO7evSsAiBEjRth8WWeXnZ0tIiMjxc6dO0Xnzp2NFmucg21dvXpVABCfffZZuctyFrazceNGAUB88cUXetPz8vIEAPHss8/qpnEO1nPr1i0hhBD79u0rs1ij8ppbKyM+Z42RNmLECLRu3dpguvZYf0BAgG5aQkICAKBZs2Z6y4aEhCAsLAwHDhyw+bLObtasWejYsSP69OlT6jKcg2398ssvAIBWrVqVuyxnYTtNmzYFACQlJelN1/7esGFD3TTOwXoiIyNNWo7Ka26tjLhYYw7p+PHj8PDwQL9+/XTTLl++DAAICwszWD48PBx37txBXl6eTZd1ZkePHsWmTZvwySeflLkc52BbZ8+eBQAIITBu3DjUqlULoaGh6NmzJ/bt26e3LGdhO82aNcP06dOxcuVK7Ny5E2q1Gjdu3MDUqVNRq1YtzJw5U7cs51D5qLzm1sqIizXmcFJSUvDTTz9h2rRpiIiI0E1/+PAhAMDHx8dgHe20rKwsmy7rrFQqFSZMmIC4uDhUq1atzGU5B9tKT08HAAwcOBCPP/44Lly4gMTERAQEBKBnz5746aefdMtyFra1ZMkSzJkzB0OGDIGPjw9q166NKlWq4PDhw/zZZGdUXnNrZcTFGnMoQghMmTIFjRo1wkcffWR0GWPDSZQ2xIStlnU2ixYtglwux/jx401eh3Owjfz8fABAy5Yt8dZbb8HX1xfVqlXDmjVr4OvrizfffBOixI1pOAvrU6lUeP755zF//nzEx8cjNzcXaWlpCAoKwtNPP41Tp04ZrMM5VD4qr3lFM+JijTmU2bNn48yZM9i+fTu8vb315vn7+wOA0S7lR48eAQDkcrlNl3VGV69exYIFC/Dll1+a9OHCOdiW9tt4ly5d9Kb7+fnh6aefxvXr13WHXjgL21m1ahW2bNmCmJgY9O3bFx4eHqhWrRpWr14NlUqFESNGQK1WA+Ac7IHKa26tjLhYYw5j6dKlWLNmDf744w9Ur17dYH79+vUBAKmpqQbz0tLSEB4eDl9fX5su64x+//13CCHQrVs3hIWF6X4OHz4MALrfly1bBoBzsLWaNWsCAIKCggzmhYaGAgAyMjIAcBa2tGfPHgBAx44d9abLZDK0bt0aSUlJuqKZc6h8VF5za2XExRpzCGvWrMHcuXOxa9cu3Zv//v37uH79um6ZTp06AQASExP11s3IyEBqaio6d+5s82Wd0WuvvYasrCykpaXp/bRv3x4AdL9PnToVAOdga23btgUA3L1712DevXv3APyvaOMsbEfbK2Kst9nNrei/1uzsbACcgz1Qec2tlpFJA3wwZkfff/+9CAgIEMeOHdOb/vXXX4vRo0frfr93757w8/MrdfDBXbt22XxZV1LaOGucg209fPhQ+Pv7iy5duuhNf/TokahataqoX7++bhpnYTvvv/++0XHWlEqliIyMFF5eXiI3N1cIwTnYQnnjrFF5za2VERdrjLTt27cLDw8P0bt3bxETE6P3M3DgQL1iTQghli9frndbj/Pnz4uwsDAxfPhwg23ballXUVqxJgTnYGurVq0SAMQnn3wilEqlyM7OFi+//LLw8PAQu3fv1luWs7CNu3fvisjISFGtWjWRkJAgNBqNyM7OFhMmTBAAxEcffaS3POdgXeUVa0LQec2tkREXa4w0bUFQ2k/JYk0IIdavXy+aN28uQkJCRK1atcQHH3wglEql0e3balln1rFjR1GtWjXh4eEhAIhq1aqJatWqifv37+stxznYVnx8vGjTpo0ICAgQgYGBok+fPuLIkSNGl+UsbOPOnTtiwoQJombNmsLf31/I5XLRvn17sW7dOqPLcw4VN3/+fFGtWjURGBgoAAhfX19RrVo1MWDAAKPLU3nNK5qRRIgS13gzxhhjjDEy+AIDxhhjjDHCuFhjjDHGGCOMizXGGGOMMcK4WGOMMcYYI4yLNcYYY4wxwrhYY4wxxhgjjIs1xhhjjDHCuFhjjDHGGCOMizXGGGOMMcK4WGOMMcYYI4yLNcYYY4wxwrhYY4wxxhgjjIs1xhhjjDHCuFhjjDFmMo1GY+8mVApX2U/mGLhYY4wxG3vzzTchkUh0P9988429m2SRffv2oX379vZuRqUYMWIEVq9ebe9mMAaAizXGHErt2rX1/tPX/pR09+5duLm5QS6Xo7CwsNztdOnSpRJa77rmzp2L1NRUDBs2zN5NsdiqVavQq1cvjB492t5NsYqEhAS0a9cO7u7uRudPmjQJb7zxBiZNmlTJLWPMEBdrjDmQ48eP488//9T9fuzYMaSmphos99tvv0EIgZycHBw8eNDodhISEiCVSpGUlIQtW7bYtN2uzs/PD2FhYfD29rZ3Uyyye/duTJo0Ce+99x4mT55s7+ZU2MyZM/H2228jJCSk1GW6du2KzZs348svv8R//vOfSmwdY4a4WGPMgYSEhKBDhw6oWbMmAODMmTMICwszWO7XX381+rj4dk6fPo22bduifv36qFq1qu0azRyaQqHAxIkTUbVqVcycOdPezbGKcePG4ciRI2jZsmWZy/Xt2xddunTBhx9+iEuXLlVS6xgzxMUaYw6od+/eAIwXYhqNBrt27cLgwYMBADt37jS6jd9++w3PPvus7RrJnMLGjRuRnJyM8ePHw9fX197NsYpGjRqZvOz06dOhVCqxePFiG7aIsbJxscaYA+rTpw8AYM+ePQbnpB09ehQFBQWYP38+AODChQu4ceOG3jIKhQIHDhzQbQcoOodn6tSpeOKJJ+Dv7w8fHx88/vjjiImJQV5ent76Xbp0MThvbs6cObr5c+bM0ZtXu3ZtvfWvX7+OCRMmoFatWvD09ERISAj69u2L33//XW+5lJQUg3PrHj58iEmTJiEsLAwymQxNmzbFt99+q7debm5umefklWx/cR999JHBxQAHDx5Ep06d4Ovri4iICLzxxht49OgRAODkyZPo3r07qlSpgpCQEIwbNw7Z2dklIzNw6NAh9OrVC4GBgfD19UXHjh3xyy+/lLp8RkYGoqOjUb9+fchkMgQGBqJr167YuHGjwbLu7u56r71CocDs2bNRr149eHp6mnWe4rp16wAA3bp1K3UZIQQ2bNiArl27IjAwEJ6enqhRowYGDRqENWvW6N4/f/75p95rO2bMGNy6dQvDhg1DYGAggoKCMHjwYFy9ehVA0bmXL7/8MoKDg1GlShX07t0b//zzj8Hzl3y/GfvZv3+/SftbUteuXeHm5oZNmzYZPf+TsUohGGMOJzs7W3h4eAgAYv/+/XrzZs+eLfr16yeEEKJJkyYCgFi+fLneMr///rsICQkRGo1GN00mk4m6deuK+Ph4cfXqVZGYmCiWLl0q/P39RatWrcSjR490y96/f19s3LhRABD16tUTqampIicnRzc/JydHpKSkiODgYLFo0SJx9+5d3byDBw8Kf39/ERISItatWyeSkpLE7t27xdNPPy0AiLi4ON2yarVapKamih9//FEAEG3bthXdunUTGzZsEMnJyeLXX38VtWvXFgDEli1bdOtpNBq99Tp37qy3//fv3xepqakCgCj5MZibmytSU1PFsGHDBAAxduxYMWjQIHHq1Clx7tw5MXHiRAFADBo0SJw/f16MGDFCN2/cuHECgBg4cKDR3EaPHi0AiJdfflm0bt1a7N69W9y4cUP8+uuvomHDhgKAWLZsmcF6//zzj4iIiBC+vr7i888/FxcvXhQJCQmif//+AoCYNGmS3vJpaWni2LFjAoCoWbOmeO6558RHH30kkpKSxLFjx0TDhg0NXhNjCgoKhEwmEwDEnTt3jC6jVqvFyJEjdft16tQpceXKFREfHy8aNGggAIg333xTCCGEUqkUqamp4pNPPhEARL9+/US3bt3E7t27xZUrV8Rnn30m3N3dRY0aNcTdu3fF4MGDxe7du8Xly5fFxx9/LKRSqYiIiBC5ubnltt0UMTExQiqVlrtcvXr1BACRkJBgledlzFxcrDHmoDp16iQAiHfeeUdveqtWrXT/4c+YMUMAEP3799db5q233hIvvvii3rR69eqJo0ePGjzP2rVrBQCxaNEivekqlUpUr15dABB//fWXwXrbtm0T3t7e4uHDh7ppDx48ENWqVRMAxJEjR/SWz8/PF+Hh4cLNzU0cO3ZMb96+fft0hdX69ev15v30008CgOjevbtBG7TrlVaYGCvWtLSFVXBwsMjPz9dNV6vVombNmgKA6NKli948lUolIiIihEQiEbdu3Sp1m97e3iI9PV1v3rVr14Snp6fw9PQUFy9e1Ntm06ZNBQDx/fff662j0WhEy5YtBQARHx+vNy85OVm3f2+//bbevOXLl4sBAwYY3e/iEhMTBQAhlUr1Cvvi5s+fryu8Srp69aqQyWTijTfe0Jv+9ddf69p2/PhxvXkvvfSSACCaNWtmMG/48OECgFi3bl25bTeFqcVa586dBQDx+eefW+V5GTMXHwZlzEFpD2EWP2/t7t27OHXqlO5cNO0ye/fuhUKh0C3322+/6R0CBYArV67gqaeeMniedu3aAQB27NihN10qleLll18GUDSsQ0lfffUVhg0bBn9/f9201atXIz09HR06dEDbtm31lvfy8sKYMWOg0WjwxRdfGN1nPz8/DB06VG+adjtnzpwxuk5FDR06FF5eXrrf3dzcdCemN27cWG+eVCpFq1atIIQosz1Dhw5FaGio3rQ6deqgb9++UCqVWL58uW76tm3bcP78edSoUcNg3yUSCSZOnAgAWLZsWanPN2XKFL3fJ02ahJ9++qnU5bW0h88DAwONDhFTWFiIhQsXAgDeeOMNg/l169ZFr169IJPJjG6/cePGaN26td407e8qlcpgnvb9efr06XLbbk3BwcEAig7fM2YPXKwx5qC0Bdnff/+tG77jt99+w2OPPYZ69eoBADp27IgqVarg0aNHSEhIAADcunULly5dQq9evfS2l5WVhTlz5uCpp55CaGgo/Pz8UKVKFTzxxBMAgNu3bxu0YezYsQCATZs2IScnRzf91q1b+O233zBhwgS95Xft2gUARotCoOg/d6DofC5jatasCQ8PD71p2itZHzx4YHSditK+lsXJ5fJS52mL07La07hxY6PTtcXJgQMHdNO0r1nr1q2NFkza1+zIkSNQq9UG8729vXVXD5tLm2lpxdaJEyeQmZkJoPRMf/75ZyxYsMDoPFu8tqbYtGkTgoODsXDhQqjVagQHB5fafuB/+1/8Pc5YZeJijTEH1aJFC4SHhwMoKtKAol624j1mnp6e6N69O4D/XRX666+/olWrVnpjTKWnp6NVq1aIjY1FkyZN8P333+PUqVM4c+aMbj2lUmnQhvr166NTp0549OgRNm3apJu+evVqNGrUyGC0+5s3bwIAPv/8c1SpUsXgZ9q0aQCMF4YAEBQUZDBNW7zZ6vZAxoY10RZNxua5uRV9rBornLRK9qppaTNJSUnRTdO+Ztu3bzf6mg0YMAAAUFBQoCucjG3TEiqVCkBRj6Ex2rZ5enrq9aCayhavrSmGDx+Oe/fuIS8vD0II3Lt3D8eOHSt1ee3+a18Pxiqb8aGbGWMOoXfv3vjmm2/w22+/4eWXX8auXbuwfv16vWX69OmDn376Cb/++is+/vhjo0N2zJ07F1evXtVtr7jSRnjXGjt2LBISErBq1Sq8+uqr0Gg0WLNmDaZPn17qOq+//rru8J0xxnqQyppuCSGEScuV9ZyWtqe89YzNHzJkCD788MMy1wsMDDT7ucri4+MDwHihbg22eG1tQXsKgbMMXcIcDxdrjDmwPn364JtvvsHu3btx+PBh5Ofno3PnznrL9O3bFwBw6dIlJCUlYe/evQaDm2oPu5U8NGqKoUOH4vXXX8exY8eQmJiImzdv4u7duxg1apTBsrVq1cLFixcBAI899pjZz2UubY+IsZ6YjIwMmz9/aXJzc41O17YpIiJCN61WrVoAis4Pq4zXrDjtgMulDUWibZtSqUR2drbuEKaz0R7+rFatmp1bwlwVHwZlzIH17NkTUqkUDx48wJw5c9ClSxeDWxrVqFEDTZs2BQDMnj0bbm5uBufnlHVYqbRDklo+Pj4YPnw4gKILDb766is8//zzRg9ZagfzLe2QkxACXbt2LbNXzhxlneOkHcvLHhITE41OP3HiBADoFdzlvWZA0U3HX3jhBSu2sEiDBg0AAHl5eUYLzFatWulyPn78uNFtTJ061eACB0eTlpYGAGjYsKGdW8JcFRdrjDmwwMBA3dWQe/fuNbjCU0s7/YcffkCvXr0MzkHSFm/G7nbwww8/lNsO7YUG3377LXbs2IFXX3211OXCwsJw8OBBo8XHunXrsH//ft15dhXVoEEDeHt7IykpyaB3aM2aNVZ5DkvEx8cbnF+WnJyMnTt3wtPTU+/+mwMGDECzZs1w8+ZNxMfHG2zr4MGD2LRpEzp16mT1dgYHB+sKtitXrhjM9/DwwNtvvw3A+NWoiYmJWL58ucGgyI5ECIHLly9DIpGgQ4cO9m4Oc1FcrDHm4Iqff1ZesVZyea1Zs2ZBLpdjz549ePXVV3H69GmcP38e77//PlauXAmgqPctLS0NWVlZBuu3bdsWTZs2RVZWFurWrWtwKFbL398fW7Zsgb+/P5577jmsWbMG165dw7lz5zBv3jxMmDAB06ZNQ79+/XTrpKWl6QobpVKJtLQ0XS9PZmamrtej5LJA0XAgY8eORWFhIUaOHIkTJ07g6tWr+M9//gNPT0+99bSHIPPz85GWlob8/HwARVfJpqWlQa1WWzwvJydHb17z5s3Rs2dP7NmzBzdv3sTvv/+OPn36QKlU4r///a9eD45UKkV8fDxq1KiBMWPG4OOPP0ZSUhIuXryIL774AoMGDcKgQYMwadIk3ToZGRm6/dHmlpaWZtG5Z9oLGI4ePWp0/owZMzBixAhs27YN48ePx5kzZ3D16lVs3LgR/fr1Q8eOHXXDepR8DxV/zSydZ2uJiYnIzc1F27Zt+TAosx+7jvLGGKuw48eP6+4kUBqFQiH8/PzKHIn+/PnzIioqSlStWlW4u7uL6tWri1GjRok//vhDN4ApADF69Gij6//3v/8VAMSCBQvKbfPNmzfF5MmTRe3atYWnp6cICwsTXbt2FfHx8QaDrxZ/bu1PTEyMEOJ/g5UW/yk5AK5CoRAzZ84UERERwsPDQ9SrV08sWLBAqNVqvfUiIiKEEPoDthb/SU5OtnjeG2+8oTft66+/FvHx8aJ9+/aiSpUqwtvbW7Rv3178/PPPpb5m9+7dE2+//bZo0KCBkMlkIjg4WLRr106sWrVKKJVKvWVr1apltC379u0rN5uSrl69Ktzd3UWnTp1KXUaj0Yh169aJTp06CX9/f+Hr6yuaNWsm4uLi9O58UXyg3pKvh6XzbG3mzJkCgNiwYYPNn4ux0kiEMPGSKMYYYy5pxowZWLJkCfbs2VPmPUKdTXp6OurWrYsWLVro7mvKmD1wscYYY6xMCoUCzz77LJKTk7F//36HPgfNVLm5uXjuuedw9epV/PXXX4iMjLR3k5gL43PWGGOMlUkmk2HHjh3o0KGD0XMendH48eNRUFCAw4cPc6HG7I571hhjjJns8OHDBnemcEZ//fUX2rZtq7trAmP2xMUaY4wxxhhh/JWBMcYYY4wwLtYYY4wxxgjjYo0xxhhjjDAu1hhjjDHGCONijTHGGGOMMC7WGGOMMcYI42KNMcYYY4wwLtYYY4wxxgj7fz8XChx59ei0AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "(<Figure size 640x480 with 2 Axes>,\n",
       " [<AxesSubplot:>, <AxesSubplot:xlabel='Wavenumber (cm⁻¹)'>])"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from radis import plot_diff\n",
    "plot_diff(s1, s2, 'abscoeff')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### R = broadening_cutoff  / wstep * N_lines / log(spectral range / wstep)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Loading Scaling Parameter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Wavenumber Range\n",
    "spectral_range = s1.conditions['wavenum_max_calc'] - s1.conditions['wavenum_min_calc']\n",
    "# Grid resolution\n",
    "wstep = s1.conditions['wstep']\n",
    "# Number of points on the spectral grid\n",
    "points = spectral_range/wstep\n",
    "# Broadening cutoff to calculate convolution over\n",
    "broadening_cutoff = s2.conditions['broadening_max_width']\n",
    "# Number of lines\n",
    "n_lines = s2.conditions['total_lines']\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "DLM scales as for 'fft' broadening method - **spectral range / wstep * log(spectral range / wstep)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "5685249.671602007\n"
     ]
    }
   ],
   "source": [
    "DLM_Scale = points * math.log(points,10)\n",
    "print(DLM_Scale)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Legacy Method Scales as - **broadening_cutoff * spectral range / wstep^2 (convolution) * N_lines**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "22979963999979.375\n"
     ]
    }
   ],
   "source": [
    "Legacy_Scale = broadening_cutoff * spectral_range/math.pow(wstep,2) * n_lines\n",
    "print(Legacy_Scale)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## **R formula**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4042032.509981925\n"
     ]
    }
   ],
   "source": [
    "R = Legacy_Scale/DLM_Scale\n",
    "print(R)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Inference\n",
    "As we can see, R = 4e6 and already we are seeing that Legacy method is outperforming LDM by a big margin (**around 700 times faster for the given conditions**). We can infer from this for a value of R = 4e6, we are in Legacy territory, therefore we can safely assume that the critical value at which we want to switch between legacy and LDM will be larger than current R value "
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "radis",
   "language": "python",
   "name": "radis"
  },
  "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.6.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}