{ "cells": [ { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Import supporting package" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import xarray as xr\n", "import numpy as np\n", "\n", "from uncertainties import ufloat\n", "from uncertainties import unumpy as unp\n", "from uncertainties import umath\n", "\n", "import matplotlib.pyplot as plt\n", "\n", "from DataContainer.ReadData import read_hdf5_file\n", "from Analyser.ImagingAnalyser import ImageAnalyser\n", "from Analyser.FitAnalyser import FitAnalyser\n", "from ToolFunction.ToolFunction import *\n", "\n", "from ToolFunction.HomeMadeXarrayFunction import errorbar, dataarray_plot_errorbar\n", "xr.plot.dataarray_plot.errorbar = errorbar\n", "xr.plot.accessor.DataArrayPlotAccessor.errorbar = dataarray_plot_errorbar\n", "\n", "imageAnalyser = ImageAnalyser()" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Start a client for parallel computing" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "D:\\Program Files\\Python\\Python38\\Lib\\site-packages\\distributed\\node.py:182: UserWarning: Port 8787 is already in use.\n", "Perhaps you already have a cluster running?\n", "Hosting the HTTP server on port 59390 instead\n", " warnings.warn(\n" ] }, { "data": { "text/html": [ "
\n", "
\n", "
\n", "

Client

\n", "

Client-9a0ca19d-ecb4-11ed-af18-9c7bef43b4fb

\n", " \n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", "
Connection method: Cluster objectCluster type: distributed.LocalCluster
\n", " Dashboard: http://127.0.0.1:59390/status\n", "
\n", "\n", " \n", " \n", " \n", "\n", " \n", "
\n", "

Cluster Info

\n", "
\n", "
\n", "
\n", "
\n", "

LocalCluster

\n", "

40e09024

\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", "\n", " \n", "
\n", " Dashboard: http://127.0.0.1:59390/status\n", " \n", " Workers: 6\n", "
\n", " Total threads: 60\n", " \n", " Total memory: 55.88 GiB\n", "
Status: runningUsing processes: True
\n", "\n", "
\n", " \n", "

Scheduler Info

\n", "
\n", "\n", "
\n", "
\n", "
\n", "
\n", "

Scheduler

\n", "

Scheduler-616be11d-c546-4d58-ad88-3bacd53b13c0

\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
\n", " Comm: tcp://127.0.0.1:59391\n", " \n", " Workers: 6\n", "
\n", " Dashboard: http://127.0.0.1:59390/status\n", " \n", " Total threads: 60\n", "
\n", " Started: Just now\n", " \n", " Total memory: 55.88 GiB\n", "
\n", "
\n", "
\n", "\n", "
\n", " \n", "

Workers

\n", "
\n", "\n", " \n", "
\n", "
\n", "
\n", "
\n", " \n", "

Worker: 0

\n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", "\n", " \n", "\n", "
\n", " Comm: tcp://127.0.0.1:59424\n", " \n", " Total threads: 10\n", "
\n", " Dashboard: http://127.0.0.1:59428/status\n", " \n", " Memory: 9.31 GiB\n", "
\n", " Nanny: tcp://127.0.0.1:59394\n", "
\n", " Local directory: C:\\Users\\Jianshun Gao\\AppData\\Local\\Temp\\dask-worker-space\\worker-ybrgy95q\n", "
\n", "
\n", "
\n", "
\n", " \n", "
\n", "
\n", "
\n", "
\n", " \n", "

Worker: 1

\n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", "\n", " \n", "\n", "
\n", " Comm: tcp://127.0.0.1:59422\n", " \n", " Total threads: 10\n", "
\n", " Dashboard: http://127.0.0.1:59427/status\n", " \n", " Memory: 9.31 GiB\n", "
\n", " Nanny: tcp://127.0.0.1:59395\n", "
\n", " Local directory: C:\\Users\\Jianshun Gao\\AppData\\Local\\Temp\\dask-worker-space\\worker-ms_3izxs\n", "
\n", "
\n", "
\n", "
\n", " \n", "
\n", "
\n", "
\n", "
\n", " \n", "

Worker: 2

\n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", "\n", " \n", "\n", "
\n", " Comm: tcp://127.0.0.1:59434\n", " \n", " Total threads: 10\n", "
\n", " Dashboard: http://127.0.0.1:59437/status\n", " \n", " Memory: 9.31 GiB\n", "
\n", " Nanny: tcp://127.0.0.1:59396\n", "
\n", " Local directory: C:\\Users\\Jianshun Gao\\AppData\\Local\\Temp\\dask-worker-space\\worker-qteiamcm\n", "
\n", "
\n", "
\n", "
\n", " \n", "
\n", "
\n", "
\n", "
\n", " \n", "

Worker: 3

\n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", "\n", " \n", "\n", "
\n", " Comm: tcp://127.0.0.1:59421\n", " \n", " Total threads: 10\n", "
\n", " Dashboard: http://127.0.0.1:59426/status\n", " \n", " Memory: 9.31 GiB\n", "
\n", " Nanny: tcp://127.0.0.1:59397\n", "
\n", " Local directory: C:\\Users\\Jianshun Gao\\AppData\\Local\\Temp\\dask-worker-space\\worker-r2y5vxh3\n", "
\n", "
\n", "
\n", "
\n", " \n", "
\n", "
\n", "
\n", "
\n", " \n", "

Worker: 4

\n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", "\n", " \n", "\n", "
\n", " Comm: tcp://127.0.0.1:59425\n", " \n", " Total threads: 10\n", "
\n", " Dashboard: http://127.0.0.1:59435/status\n", " \n", " Memory: 9.31 GiB\n", "
\n", " Nanny: tcp://127.0.0.1:59398\n", "
\n", " Local directory: C:\\Users\\Jianshun Gao\\AppData\\Local\\Temp\\dask-worker-space\\worker-3ykst_31\n", "
\n", "
\n", "
\n", "
\n", " \n", "
\n", "
\n", "
\n", "
\n", " \n", "

Worker: 5

\n", "
\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", " \n", "\n", " \n", "\n", "
\n", " Comm: tcp://127.0.0.1:59423\n", " \n", " Total threads: 10\n", "
\n", " Dashboard: http://127.0.0.1:59429/status\n", " \n", " Memory: 9.31 GiB\n", "
\n", " Nanny: tcp://127.0.0.1:59399\n", "
\n", " Local directory: C:\\Users\\Jianshun Gao\\AppData\\Local\\Temp\\dask-worker-space\\worker-egpqiwk3\n", "
\n", "
\n", "
\n", "
\n", " \n", "\n", "
\n", "
\n", "\n", "
\n", "
\n", "
\n", "
\n", " \n", "\n", "
\n", "
" ], "text/plain": [ "" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from dask.distributed import Client\n", "client = Client(n_workers=6, threads_per_worker=10, processes=True, memory_limit='10GB')\n", "client" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Read data" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Set file path" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# filepath = \"//DyLabNAS/Data/Evaporative_Cooling/2023/05/03/0043/*.h5\"\n", "# filepath = \"//DyLabNAS/Data/Evaporative_Cooling/2023/04/18/0003/2023-04-18_0003_Evaporative_Cooling_000.h5\"\n", "\n", "# filepath = \"//DyLabNAS/Data/Repetition_scan/2023/04/21/0002/*.h5\"\n", "\n", "filepath = r\"./testData/0002/*.h5\"\n", "\n", "# filepath = r\"./testData/0002/2023-04-21_0002_Evaporative_Cooling_0.h5\"\n", "\n", "# filepath = r'd:/Jianshun Gao/Simulations/analyseScripts/testData/0002/2023-04-21_0002_Evaporative_Cooling_0.h5'\n", "\n", "# filepath = \"//DyLabNAS/Data/Evaporative_Cooling/2023/04/18/0003/*.h5\"\n", "\n", "# filepath = \"//DyLabNAS/Data/Evaporative_Cooling/2023/05/04/0000/*.h5\"" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "groupList = [\n", " \"images/MOT_3D_Camera/in_situ_absorption\",\n", " # \"images/ODT_1_Axis_Camera/in_situ_absorption\",\n", "]\n", "\n", "dskey = {\n", " \"images/MOT_3D_Camera/in_situ_absorption\": \"camera_1\",\n", " # \"images/ODT_1_Axis_Camera/in_situ_absorption\": \"camera_2\",\n", "}\n", "\n", "dataSetDict = {\n", " dskey[groupList[i]]: read_hdf5_file(filepath, groupList[i])\n", " for i in range(len(groupList))\n", "}" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Rechunk the data for parallel computing" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "dataSet = dataSetDict[\"camera_1\"]\n", "\n", "scanAxis = dataSet.scanAxis\n", "\n", "# dataSet = dataSet.chunk(\n", "# {\n", "# # \"compZ_current_sg\": \"auto\",\n", "# \"sin_mod_freq\": \"auto\",\n", "# \"runs\": 2,\n", "# \"x\": \"auto\",\n", "# \"y\": \"auto\",\n", "# }\n", "# )" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Calculate absorption imaging" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## get OD images" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "
<xarray.Dataset>\n",
       "Dimensions:     (runs: 3, x: 1200, y: 1920)\n",
       "Coordinates:\n",
       "  * runs        (runs) float64 0.0 1.0 2.0\n",
       "Dimensions without coordinates: x, y\n",
       "Data variables:\n",
       "    atoms       (runs, x, y) uint16 dask.array<chunksize=(1, 1200, 1920), meta=np.ndarray>\n",
       "    background  (runs, x, y) uint16 dask.array<chunksize=(1, 1200, 1920), meta=np.ndarray>\n",
       "    dark        (runs, x, y) uint16 dask.array<chunksize=(1, 1200, 1920), meta=np.ndarray>\n",
       "    shotNum     (runs) int64 0 1 2\n",
       "    OD          (runs, x, y) float64 dask.array<chunksize=(1, 1200, 1920), meta=np.ndarray>\n",
       "Attributes: (12/96)\n",
       "    TOF_free:                          0.02\n",
       "    abs_img_freq:                      110.858\n",
       "    absorption_imaging_flag:           True\n",
       "    backup_data:                       True\n",
       "    blink_off_time:                    nan\n",
       "    blink_on_time:                     nan\n",
       "    ...                                ...\n",
       "    y_offset_img:                      0\n",
       "    z_offset:                          0.189\n",
       "    z_offset_img:                      0.189\n",
       "    runs:                              [0. 1. 2.]\n",
       "    scanAxis:                          ['runs']\n",
       "    scanAxisLength:                    [3.]
" ], "text/plain": [ "\n", "Dimensions: (runs: 3, x: 1200, y: 1920)\n", "Coordinates:\n", " * runs (runs) float64 0.0 1.0 2.0\n", "Dimensions without coordinates: x, y\n", "Data variables:\n", " atoms (runs, x, y) uint16 dask.array\n", " background (runs, x, y) uint16 dask.array\n", " dark (runs, x, y) uint16 dask.array\n", " shotNum (runs) int64 0 1 2\n", " OD (runs, x, y) float64 dask.array\n", "Attributes: (12/96)\n", " TOF_free: 0.02\n", " abs_img_freq: 110.858\n", " absorption_imaging_flag: True\n", " backup_data: True\n", " blink_off_time: nan\n", " blink_on_time: nan\n", " ... ...\n", " y_offset_img: 0\n", " z_offset: 0.189\n", " z_offset_img: 0.189\n", " runs: [0. 1. 2.]\n", " scanAxis: ['runs']\n", " scanAxisLength: [3.]" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dataSet = imageAnalyser.get_absorption_images(dataSet)\n", "\n", "dataSet" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Select region of interests" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "# imageAnalyser.center = (529, 962)\n", "# imageAnalyser.span = (100,100)\n", "# imageAnalyser.fraction = (0.1, 0.1)\n", "\n", "# imageAnalyser.center = (890, 1150)\n", "# imageAnalyser.span = (600,600)\n", "# imageAnalyser.fraction = (0.1, 0.1)\n", "\n", "imageAnalyser.center = (890, 950)\n", "imageAnalyser.span = (100,100)\n", "imageAnalyser.fraction = (0.1, 0.1)\n", "\n", "dataSet_crop = imageAnalyser.crop_image(dataSet)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAASIAAADQCAYAAACjk2wQAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABQ40lEQVR4nO29fbhdVX0u+v7ITsJGE0OkiRCQTQ01IFBFD6BYjEif4ieo+FE/qh699OljW7X2Vqyn9bSnnqter9ZzPdWHg4p6bJEiflDRXlHiFzUIEY1CkCgLScQgBCTIJslOxv1jjHfP337nGGutvXf2XnvBeJ9nP2PPueYcc8wxx/iN3/ewEAIqKioqBomDBt2AioqKikqIKioqBo5KiCoqKgaOSogqKioGjkqIKioqBo5KiCoqKgaOSogqDgjM7O1mttXMbjazPyhcc4yZbUzXfcbMlsx3OysWJiohGmJYxMC/oZkdD+DlAJ4A4GwA/2RmizKXvgfAB0IIawHcA+D189fKioWMgQ/iiunBzMYS1/FJAD8CcJSZ3e9+P8/MLk7/X2xm/8PMrjGzn5nZeen84Wb2TTO7wcx+ZGa/N8tmnQPgkhDC7hDCrQC2AjhF2m0AzgRwWTr1CQDnzvK5FQ8RjAy6ARUzwrEAXhNC+C4AxDlexOEAng5gHYAvIhKCVwD49xDCuxLncojeZGYfAPDMTH2XhBDeLefWAPiuO96Wznk8GsC9IYSJLtdUPExRCdFw4jYSoT7w+RDCfgA3mtnqdO57AD5mZovT7zfoTSGEtxyYplZU9EYVzYYTv5FjHzB4sPy22/1vABBC+CaAMwBsB3Cxmf2RPsDMPpBEN/27INOe7QCOcsdHpnMedwNYYWYjXa6peJiickQPDewws+MA3AzghQB2dbvYzI4GsC2E8L/MbCmAkwF80l8zTY7oiwD+2czeD+AIRNHxWqkvmNnVAM4DcAmA1wD4wjSeUfEQRuWIHhq4AMC/AbgGwB19XL8ewA/M7PsAXgbgg7N5eAjhxwAuBXAjgK8AeGMIYR8AmNmVZnZEuvRtAP7CzLYi6ow+OpvnVjx0YDUNSEVFxaBROaKKioqBoxKiioqKgaMSooqKioGjEqKKioqBYyjM92aHBGCFP4PGdYYhTfvcbwCwOJV7+nlCKpel8r5ptI5xm3w+n/sgGjrPbta28PyEO8f32Z9KvufB8vtv3DHPsf5RaUvufRbJNWwLo0UemcrxVLKP9rn6H5S2dgPvn4lxpN97c30xnWewD/bKc3Us8fwomm+nz1sk5T4pPTiGbrsrhPBbfbXW7FOrgFflAvp2A9gJvCOE8N/7qWshYCgIUSRCf+yOR9AMABKPXe43AKATcT8+c7xnfSqvmkbbGKXAyc7nbkUzYVcW2sLzO905vg8JAN9zbSqXp3KTu57nWP+6VB6Zytz78DmrpfxOKk9P5Q9T6Ynacen/rdLWbsgR3X7R7725vpjOM/g97pTn6lji+ePQfDt9HvuX7eH4yLl4cQz9H7dNo8GHvR7AIzI/3Argo8Bh06hr4BgSQkScnMq7AfCb6YflYN0h51ehGWCrUskJdFIqOWHZLX7F47VKPHSA+fZwwm5CHvz9blcnJzzrYVs56RW73LVHy7VbCveMuntK/bVRjsdcG1k/JxCPSSw7qeQk3Iv2QjFROM5hrVxDjoV17nTHpQVpXK5l2z0B4fjg8/he7BsuLKxrs7tXxwWv0e/pF1F+r+nQnwYHuxZ5LJ1RbYPFkBEiP6H54bny6MfkACShWIz2iqeimA7OXWg4g1vk3hF3Te5eoD2Aieen8lL0RlcnaWkL3/k2Oc+BzzbucP+zbUrolOP078VJxntWyjE5Mk8ItS3KyR4txztdPaxXiVW34asip3JtJDZ8r9Vo+o9t0+eOyXmgGR8kjioG62KwBs334bOPxkxQIkTDmORpyAhRRUUFUTmigcGzv2SB98o1utIpGw+0uRlCdQFHo9GZKLhKrpNjYhnKeopSnavQvJeKhGwzORm2fS2asK6bCm0k2BejaItiJVBs1bqBvI4LaLgAitKb0BZH+F58D+XigIb71PfIcZ/aJtZDzkjHCblrcqdXoK0DUuSU/vyWyvGtzVwLTOXcRzPn+scoKiGaZyxBHHy5gbdYjku6FC+bExy0nGT8/exUlnQ7HqpX8HL/C93/QDP5deJqXR4capxIrOvEVF6LMmEl+Luvn+3dKSVBAsuJRnHFD5lew8f3H0UYnbg6gfm+u1C2XnIcvCKVFHFH0RApflMuLgoSwitcO0p6KrZJdYkevLdkmKCubzUa3VI/YncZI2iPfp4fNgxjmysqKlA5ogFgP6YqG9ehYf9H3Tm481yZaJm6BW1Wf6dcy+OvuGepSOQtNEDDqVwr142gWTlLK5+y756bI1f2FeTBVXUlGm5FFcHaF14cUkWzQtucewflovjurNNblfhstTwpV5OzLpVAbsZzI+r2wDao4pvIWe2UuyFnluOEVATcKeefIe3YjIY74nhk26aHqqyed0wgfmAOjC1os8Cr5B410ec+mYoHBMWezWgGJUUi1qPWGBIVr3spESAdeN6yQr1KiQAp/LpYMtfnRD5CiQnrUp+knOWK/7Ovt0tJrETjoqC6ms9l6uWxJ+pA2TKqeiYPtm2LHGufnIi2T1BpodJ2+lLPU0T07z2WSra79N26oyqrKyoqBo5KiOYdSxGVwFwBb0J7dVI/Ea54foVXNlwtU4QXV9RJT61YxFY5v9y1kaLSWCq9b5N/3iiAB5CHWor4nB1ouDD1zSE8h8fnqe8PoX2hIjDQXv3Vt0pFRKDNrXTjMgh+721yD8G28/06aL8726Y+ZNxkhCJ17l6iWxu7WfCAtnVwGRol/ky8zBssRlVWzzMC4kfrpGM/YGil0sHSjY3OOS56cGKvQtvM3csxL0fc+BsJAdus3e8Hszq58X3U69ejNJH4XO/Aqf3D53kHP/9cb/XS4V8K8ch9p2NT2U8YDZ1ISQTZlhPl2OuiFCXv5WvleCXa/efdD3xd/h3U5YOiLNtG7/JxV5b0jdND1RFVVFQMHFU0m3fswVTRyltUejmD6aoG9A4+JbexHO2YL+UkeK+uaivRFqNKbVYxIXcNkYvZ4v03yTX6XuRqrkKbk9PnqTLeczClvlBO04seuXpy8Ep/tk05rhKX6mPo+O6PTmXOYdK30b8/71U/slwf8RpyaWzbnVJ6zp3nVHE/PYwuBUYzW9ot3Y/pJx8YMIaEEC1CnDg+LkkHBQeCOif+EGXw9Tmh1KIy7s6xProDUBxgO57h7mE7VDdDqN6JBChnstY2+UkHxD5RAkQooWA7lmd+I0oWKhKuNe5eNXOXLEhAO+BYrVYa3+XbotY51fd4kVstXixVRBvJnNfxUILvG95DvZ86uJIwcUGklRAo6/T6w8FLgNFMRrGlE6iEqKKiYn6weARYnCFEI0O4H8aQEKJ9iKuGV9RyNVHnOeVC/OpMLkJFmVJM04i7hwpa1qt5iDqZdnsrFdAW/XLdryu3cg5sI5WkY+iurAXauYf2oh2D1ysdBzm9LWj792hqDcVKdw3boO+l7fBtIbdJMYhcD7kMfoMtaPcp3125HI4l78NTMmYQ/I68dxMaDlmhxgRyQt7vS50fp4klaPKuDTmGhBAxe57Pg6MWL75Kp0s9nEzXyjGhSckWoyEIHHwlhz91SvOm65LToLoCTGTq5YRVVwD2xVlo3pmTnOIAJ4k3UbMOFQfYprFUqkcyCfF2d60GspZEjL1oW5MoZqlOx+d9Yl8wPk1FtH5iAdXj+axUqprXf7+SO4R60S9D70VAcRza7Z6hGf9g5AlRHxyRmR2FuKnm6nTHhSGED8o16xE3wbw1nbo8hPD3M2tsdwwJIaqoqGhhCfIzuD+6NgHgrSGETWa2DMD1ZvbVEMKNct23QgjPm11De2NICBH9iIhxNOKVZ8uBxm+jmwLQp8MA2uKQZ51V1BuVks/V1bPb8zU8ZSyVnp1XTohQK94GtBN3ER055u97kVfMA+0VPudrxT7XOKtStPwulMMYNDyEyCmMc2E0vg4v7qlIpmJkjptSEVq/ofpPdfvGJVF3U+Y5M5yGSwu39mGMCyHcgbQrcAhhl5ndhChqKCGaFwwJIToIceL7dJ8lSxEHdDfPXRVh1AnRm1o7cm8/znu9oJYUPylKy5kSOvUYBpp3PzWVTPfKe3OTr1fqVn2vZWgTHk0dq0Tl5My13upXgvaP3sv38gRILWrqALpJzntrF79tKYBaCccEygSnG1uilKIUeNwDi5CfwVFcO9HMrnNnLwwhXJirxszGADwJ7fzAAPBUM/sBgF8A+Mu0vfgBx5AQooqKihaWIh/jETdX2RxC+IteVZjZIwF8FsCbQwjK0m4CcHQI4X4zew6Az6NxKz+gGBJCFBBXEb/6qJKaKxzPa1T+uPtNla2aEsOvsCWlbik+LZeETMHn+p04WHcvkYzwogY5OLZNsylqAnnv+KdcB0Xbki+N73vey/i4kvjVcfWVFPX8fmzr3WiHVlB07hYbqCEVOrdKObr9Pb2i4T33qAr0UtZFz1VxLvP9uo2VLliKfDxHP5uqADCzxYhE6NMhhMv1d0+YQghXmtk/mdlhIYS7ZtbgMoaEEBmmkn7P9vJ8KcE6P/IatJO+EyUiknOcLH3lbk5wPqgVaCYHj30bvXXKt4nvRZ0VJ+7nkPdkzoH3HInGw7kj95beIyd2MRiV7Vei5u9RT22KP/R8/o78Poq2WKWWLoo4XmTS1B0sSaB00enlvMh6c/d4dUGJmOg9y9B2sl2GGWEJ8vEcOS5JYGYG4KMAbgohvL9wzWMA7AghBDM7BVFHcnfu2tliSAhRRUVFCwcjT4j6i3o9HcCrAWw2sxvSub8G8FgACCF8BMB5AP7EzCYQKe7LQwhz4i45p4TIzN4C4A2IstVmAK8DcDiASxCXwusBvDqE0MMhnSEePsxhufvfQ/eTym2Ap1kXyfUsk/I2lEMSuFrrljA5K1MpdIDvwLq9SMD6R+Q33uMTp2moRSlFyhYp/bUlRbNuJ+R9kHRbJirh2Rds61oA35BrNRZLFbbjaG/kqPeq75HnVskt8V3pP1Ta4smfK6WU4fO8MlsV2xqCw7o0VAjoz8LaBYuQ9yPqw8kxhPBtNFvWlq75EIAPzaRp08WcESIzWwPgzwEcH0IYN7NLAbwcwHMAfCCEcImZfQTA6wF8uEdtiITFE4VS9rxSalcvz4+lsiPXcvD6CVZiuTk56PXL5+T4YiU4uhGid6TUJPI3yTXb3LXA1PSobAMnR0mn4sHnlXb18Hl0/HHuGt3jjc/tZO6niMn39E6CvG67XKtikL6vJ7Dabg22zb1HaRNNQheF56NJV6vQ+jm2lqMv2akfLEWzE7nHEOYByUSqHFCMABg1sxEAhyD6LZwJ4LL0+ycAnDvHbaioeGiCOiL9O0B0bj4xZxxRCGG7mb0PwM8Rl5r/D1EUuzeEwGV4G8oBSg5LELkYvx9YKYaIVavlw+dA1mtL2fVWuntYn6Z7+Aby8IpMVY6TU1HO4Vg0il51LOQKu1zO70T7M6qTYj9Qaxm5grFUeuuTWn1Uoc730uRjQPPdNGYvx9mOyrWaYzz3rUuOpf1Eums/6nOUO86lU9GEeWOp5Dv0oxzvEwfjIcMRzaVodiiAcwAcA+BeAP+KZmuKfu4/H8D58ehRmEqERtBMBnWUIzTHzwiaiVmySCkWo60DKqX04PNLOqUcWDf1F1d1aVO3zIx67XQGe8kixImkqWvHUSZArEMJkNcrlRxQdSj6dCO9tsf26KVv6fbNNSulQq2CXnTXTQbWp7Ikunn0sRbnMAur2ULDXCqrzwJwawjhVwBgZpcjaupXmNlI4oqORIEdSV6gF8Z7jxjCxAYVFXOMEXTzrB4qzCUh+jmA08zsEMRl9FkArgNwNaJZ8BIAr0GM7p0mRtGsSuoDROTEAhU7dFdQgqtlB+0dShUaqZ1bacnxqK8MV1F61ucSo2kUfg4la0+vtgPNKs/+KvnusG7P3ajlqSQu70LbcZHilDor5pLJ5VKEeHgulJxwx53zdXQT0UoiOsEx55XjuoccUeKEjkbboluK0euByhH1Rghho5ldhmbj8+8jcjhfAnCJmf1DOvfR/mtlc3ehbNnotquCsvilHEY+3w6JVC5/TQ4+IyAnVyk9aq4unfi8RtNoECeieVclRCUC5D+7jtqS6d9PYN5PUWarXJsDdVD6zuow6glQSf+i8ESaxEKJpFoSc1DrW4mQ5zzvSQDZF5rvidiJtgioGyX0iYeQ1WxO/YhCCO8E8E45/TM0X62iomKmoJVMUQnRXGMslVvRXvH4Khp9T7Z3LcrcjFplCL96qmjUK0RgAmUFMzkudRbchTZXoQpS5Zg8B1EKB2EfHJnKje4aXe1LOZ49vOgKNKJnKcWGD71QdIvR6mX1U45pJRouRkWjkiXR95nuBst+077w2UEpmqsqQHca9hw1z/Hbdsur3gVVNJtvcIPFjjunplyFOjh6L1h1iMul4+iFkrNbN5cAJZaETwqv6Ss4GTjZN8g9Y2gmgYouJJYUL32MlDofaptymz8CU0UNnqPoyTbTOLpVyhz4HCV469COQ1PkApKVAGkmSIW3+HmPeqD5Fs+Q8yy7iVS6cJHI+eyiHG8zFBCqsrqiomLgqBzRfGMv4grzonT8HfTez4yv5h31SqKZssZcsU5HezVWnw8VBX08XClZVkkJ699JLXwb5VrCi2Ya7/YNOZ/rM43X0gT/OXi/Jw/eQ3ElNyNKW2crtqDNcZTyg/v30wyaypnwWoqpHdcebQvr+IYc+1CT6TiNlp4zHU7coTo0zjdGECcMCYb/mCelUs35RLeUo6XUISQUnghxAtHiocSrtNsG0Bb92GbW79vItpGA3iLHGnDqUTKBd5ss/e4A4k3oJSsgCQWfx8l+C8o6Gga2sn7fF7oHWcmFIReQq9k2VazKiVdaP8U6vge/uU9pwvv1/Ui0OB59NtBeOZ/6RIkjqoSooqJi3lA5ovnGfsSVxitUyZncItf2k6CsVzhGN2c3dZjT/MO8dx0aEYJhElw9yQnp81eiLRKRy1CuLZc3WbkptUjl3rdXmAvfx3ORuu1SKb6PfXYS2gYCtolcIsW9S1OZ2wutVzqVHGfZTcQEpoqr2k+aRD/nWFlSEXBc5IwWpayR08RByCum5zqUfQ4wJIRoH9ofr+QFq6/kdTeltBEqzvljHcifk+PSRN6CtnmW19CqRGsXCdQ42rotndycYB33fH12Lge6r9PXs1PKEnw/cEKy3SVLoQ+CJfHSIGESKPaRd4vg/+qxTSsWrYFjqbzWXUuwP5+fSn6/nPuFxiKy/ewrdVL0xJIoOUGyHblMo9PVMyXMgiPqc18zA/BBxNQ9DwB4bQhhhgqt7hgSQlRRUdHC7KxmE+i9r9mzEaPLj0XcGubDaLaIOaAYEkLEfc0895GLSeJvHl4MIcutXIauajnuwCe4BxruopO5FoirJVdQbROzK6oFaY27tpT4Tdn4XP5u9flh2z1XSQ6nNAS0baxjOZoVvCSW6H5xG9w9uUyMvi7PHeQSkQFlZ8URlKP8r3DXAN05QCrQO6ksiXe5OjRukZhwJftWt5WaJmYR4tHnvmbnAPhkSg/7XTNbYWaHp3sPKIaEEOm+ZkBDgHS749IA8/eXxCndDWML2hOHg361nFeCsbNLW0qOjV500i2gWVc3K5e6FpDwqK7FP3dx4TfeS4LPyfIdNO9aiqdSj27/e8mbWMWTE9Ge1CQmOUseEL9nSe9X0oF5XRLfS/tYdUb8Tl6U1vGn7fCEXVPLzFDaKYV4xE96IPY1WwPgdnfM/GEPV0JUUVHRQvcNFg/EvmbzhiEhRLSaeWhsUi/Lg/f5UJFCu8EnXKfPjHIb5CSenUqu1t5SVdpqprQ6+5WY93K1pBKUnITnBlTpuVeO1bfGo5RChPeo9ee4zHMYisG2qRLei6nsz5JTJDmUzWhbvvQ45+DoxTSgLKaqBRWujaVddbVOr1DXHOnKrbH/noHmm84waT4xy6DXXvuaITb6KHdczB82WwwJIVKciHbSdXUky6Wv0BStah1RAuUnSan/N8gxCdSj0SS6p5WHlhSy9iraeIsPRaWO3Et4E/oyuba0K4mfQNQ5chFUD3GKUOwDZsTc6dpCwqDHrJvteQbaBIfHqk8iUTkF7T7ncFX9oLcAkhDwm1Lfw3GiKXf9wlESpXPBtZD28Vv6PvZgG7+B9g4tM2REZqEj6mdfMwBfBPCnZnYJ4kf99Vzoh4ChJUQVFRUTS+OfYl9/VrN+9jW7EtF0vxXRfP+6WTa5iCEhRIapIkg37lA5Ik1QD/TPEo+4+hiWQQVqR+pSLus4lLcUVsdDcge70KzQatUqKYiXoR0nVorn4ghdjsZio86AXP3pKEoHToqe69CIVdRtkuNTHyE+9wo0oiWfq5ZLtUyNo+EU1HrGY91PbTHahgh+PxXrctyWKuxfmspL5VqOgZWujexbdXAl+NyVrv4Z+g8l7F5i2L20vTXZnsUB0dJcRp/7mgUAb5xFE/vGkBAiwptetek8HnPXAM0A2ZK5p5v1BYgDk4OFk660j5oObG8J0X3Tcsnl2Q7VdanuQc3ge9GeQKVgW9axFm1CR30L+0sniW6aCDQJ4n29cO1hn53knkdRaSyV7HuKd56osS90dxBOfj6fRPIsNERC8zWVYsCUiAONWKoEiND+9OcUOR2c9u3MtpzePzKCfSNtWrJ/0f4u7VmYGDJCVFFRQexZtAR7FrUJ0d5F+1AJ0ZyAO736FA6q/GTHq0+GF0/04ygnNJZKn4NZ/Yh65U8mdqItIim3o9yV51S4SjLyWyPRfbgDV31NAkbuQv1itqJR+JIDYr+xDloD+RzPyWh8H4+1r9hHD6DhkihGsU/IJfJeindr0I4x0+28aaVjv+5AmzPlmNFsmDlOiM/plTEx50yrU4nvrqJaLixkOWaCPViC3ZnAsr3Yi9mKffONISFENN9zYF6FNjurOyPohPWvukbOsV4VmXaibQUp5cTRPDh70Uy60v7ompbVW09YD8Ub1St1XB25dwQacYf6GW8Z897cOYylcq+UQDuXz3rXFv8cxtRNoLEglnJCqWVvL6am+QXa+6lpKg8/+TQ/FAksx4mKvEDTt6xXv1vOGbPXTinajq3oL5tnb+zBYuzJRL1O9NAPLUQMCSGqqKhQRI6oTYj2YP8AWjM7zCkhMrMVAC4CcAKiGv8/A7gZwGcQl9wOgJeGEO7pr8ZcpDRXtlz+YqBZxbajzQJzhVWxxzvB6WqlW0+rctT7Ly2WewgV69S/CWiHBnQK965CwzGohYhtU3H1ZDQcD3/jc2gRS6v+oclaeE/iTo4CcLs+x1sKgXa0/E73G/uE3A79kyjeeQdVcqrqL8S2E4zd80nONFLf7/jr6yRORnvfMrU6ao7pa9FwQqWYwFwyN/Y1HUG3Za7pjX0Ywb4MIdo/hEmr55oj+iCAr4QQzjOzJQAOQfRV+FoI4d1mdgGACwC8rXs1jDXzJl4VxdSaxPOlvek91LLhHdeUwKi+gu1Qy44nYLyHjn6ai8frOAjNQ8R6OVF9AC/r52TvlWy+455Jgqr6q/Re96ZDS0R8BMBh6dzudM3EWCwnY3xJQLxrgMawsa0kljnL2Ea5lsSD78G2a9J7oEVQJ/uLUCdSL3pq5k4SwNyOHZqeV/VHPGb9W1w72YaZOTTuxWLsyUzhiaILwcLFnBEiM3sUgDMAvBYAQgh7AOwxs3PQKBU+geia3IMQVVRUKHZjKXZnpvCeSoim4BgAvwLwcTP7XQDXA3gTgNXOTfyXaIexAwDM7HwA58ejR6GdHoIrOFcTJtZX349cpL0mH1MLg48f0ghytXjw3pwSW53oKH4o9+GjsXW1LIHddqerryP3crXWlCK70HBybOMhqaQIlrgQZVR2A7hLmnBvKn8rlTexbs9ZsJ/4m25xpBbMZ7v/6SdU2vYplzlRw046qVwpJZ+/193Pe/gtumV5VA5VMwaoCKhZJICZWs0iR9R2o54YQtXvXLZ4BHHk/FnafvqDiGLYJEIIwcyyKv6UsuBCADA7Qq5ZjanJyIHezme5c0pU1JqVM+2q/oNWk1xWRA7g0t5h6lznCVTJkqKiwEq0sy0ulmvoFU6iejJasWyWJnlI78e5cm4qD0+lzwpYyqH/iFTeQGK2He3Uvpx8PK8OiB20RSQSIBIkjoHcdtjjcg11UWp58yKu5qvSODi10k2g6WPqgnQ86PjzRKhk7ewPu7EEuzOEaO8Q7ic0l9lttwHYFkLgF7kMcQTtMLPDASCVpa0ZKioqumAfFiWFtf5VZfUkQgi/NLPbzezxIYSbATwLMfvbjQBeA+DdqfzC9Gv3Sl0N01B4S5RaQfwun/68TzGi4oA+RxXDXimqnAnrJ4dC7sBzYKpUVYdCzdO8FW3OSzk99hefv7W59rD0fnzluxZPfdzuVG6IxeHfvBV3/O0x8eDc9Bu78RJ57GQCibOAWz+T/tfVmorfdXJ+MRqFM0UzciJ8H82AeYqrTxOT5TYBAJq+olLZ10cOkpwQx8liOe/Btqk1dyyVXqzUtkwPe7EEezKh9lU0a+PPAHw6Wcx+hhi9exCAS83s9Yhf8qVd7i/AizCaD0bTin45lZ5FJvusugYOWm9R4jXKnutAL2UG9Ch5fxNnozFFa5wVRZmxVHLwHo222Z4T5e5UkpiybZsw2T93pcl8WLqHXfOfUknfxXNj8cCeUeC1SVLelsIL2GQ+/pGpJM38mm+3ikgkQJpi4240Zm0VtwiNi/P6GDa8tJ8awf69CXkx3rdJHWAdQZ8kQP67eOQWz9l5P0c/onb4/TCKZnNKiEIINwB4SuanZ83lcysqHg7YU+CIKiGaN+xCOy+yhk9sSCW5BZ+Wg75FXI3J5ZAd4Eq8Fe2dPNHjuJQ6Aq4uOrKp4vQr7lq+H8U6+sNwNfZWOlWy+3f29VKUWQ8gJeQ77EVTH8dbfiSPeUrkFtYu+Sk2r4gWtjf89kUAgItOeAMAYGxlBwDwk/ecNLWOpwP4Frmyi1JJsYvQsJEx9xvfnVyh+iKRK/H+PZqoTP2F1N9nFOXtfdQ5kqzeKNpifs6nSUFCwXtLYTbdUTmigSFnii/FBREcVH6QLitcS8LAe3I6h17sdMl6BzRtJTHppt/StBmdVDIY9cOpXIe2le/Zcg/b/PpUbseku8O+dIpdSrFqLJWkDw/Gvr8Lj8ZrV34cAPBN/B4A4IkrbwAAXLfjyfFaWtiensrLfNtYIRcDEv8N0vbNKG81zbpUVzSCdsoVQs31fGEfP6YuDfwG/PadTHs0j1OvDRqAhqDpHnnTQ1RWtxXTw6isHsI9ISsqKgAGvS5p/fWjrDazj5nZnWb2o8Lv683s12Z2Q/r72x71PdPMLjezH6e/y8xsfb/vMmQcEVeqLWizt8qOdnMILDn86T2ei+JKp9sXUUxYKedznJPGp6mC1qcq0XzTfHfmOKflbQeaFVqj0MkdsH46LZ6KSY7hiQzTSD91Uvm7qUw+QY888q5Uwzj+AP8OADgi7SpzGV4MANj/g3TxXVPvxWMAbElcYCCnQgdU3RzAx8dRrGJsFzkUcjHqGLoc5b4l90vuhuOGCv3T0XBpasXS5Pl8/k1ownZKTlWKVWhzQjNLjLYHS2cjml0M4EOIu72W8K0QwvN6VWRmz011/T2Av0PM23MygI+Z2Z+GEK7sVceQEKIliJOKA+IUTN2wDui9veU6tLcU1h0eFN4LtrRxn+555eV9NRlzkqhHt0+BoQ5y1CfxXk4CT+g0BQmJo8ayETsxaTX7Vjo1lsrnpvKpqUyi2/3bYoDZUb9zzWQtVJTevifZ6dfFbzF6Wnzv8e8eGs9/CU1GkC+T4JBAsG82pJLv6835JAxsZMkhdTGa/iFRoWhLgsN7VKTvoE3IS7uveKunWkA1blHjADejPWZmtpsHOSJFP4QohPDNtJ/ZgcD/CeDcEMIP3Lkb0r5q/y9i7uuuqKJZRcWQghyR/iVCdKKZXef+zp/BI55qZj8wsy+b2RO6XPcYIUIAgBDCD1EI4VIMCUe0B1NDHjajWb00RQOhykNvXZpw5zx0S+O9aHMbugrrqpmzllDEYx3q4Ogj+VUEJIdFS56KgjehsUCRa+LzaMFj31C5+7KmaZQWV0iTE8f/rFO/BADYlfrzKbgOI4lNWpTKC5a8BwDwpcc+BwDwEzweADD+y8QRvQHRfXVKGzW0Qy1gO91v5JL4zrqDrQ+VoJMo6+GxficN5/BpQJTrJMeS2/qo1xTK7R5bCi2aHvbhoEIakIOAPjdY7IJNAI4OIdxvZs8B8Hk0TmCK33Spp9tvkxgSQkT4gEW1fiiUPR1Debtmige6P/qEq5+6jZJ1Tk2x3hxMAqB6C91hZDnauWmukGONT1uFto6B4o94XC9b3zSDqTwYN0YLF61oqftIbJ6JqwEAn8MLcXzaHn1Jcrs+MfXrZThvyj04LdX1vwE8Mf1/f2rD7YkgHZYI0F0kUD69i9/5Amj6hsT+Fankt7gFbR2djg91QGTfjWR+04T7flEjVD/Xa6Fai+a9fDbK6aPkWX0gzPd+19cQwpVm9k9mdlgI4a7M5Y8zsy9mzhuA3+7neUNGiCoqKojdWIJFc+RHZGaPAbAjBaafgqjGubtw+TmpHEXkmgIitR0H8L5+nteTEJnZ8SGEG+Xc+hDChn4ecGDhk8x35DddtXRl8uKcKqnVStItYr+Uc9lvQYTUPrW2aH5j3YbnFrTDFs6WY92BdSeAf5bfyDlIWMiuT6fjVzYi2QnS/HTpU0+OHNC9SWa7PQWO/R3eiZd8+Yop977jqL8BADwn6SQ3JkvS8b8Th811/+XJ+PUFj4kXMwEMuY/J9ZV9QnFyO5p3/6j8Rl8k9gG/5xo034d9qgnQ9PuxvzejbenUTI2aIG052lkSSopn3a+O7dVz/WM2HJGZ/QuiZ+thZrYNwDuROiNtrngegD8xM+ZweXna5yyHawC8CzED68/TuaMQLXN/3c+79MMRXWpmnwLwXkRG/r2IYRtP7XrXAcUIInvOwaRJ7oFmZuV+U9CMXUqA7y0emhmxFzhIPSFkvdRXkCBywHTkOv8b9RZjqfQECJg68Pls1sM6JDfPMgA0ynIOMBAnEZcViNl7X4eLAQAnJu/mdVtua9HpsdR+imiH4AEAwNs3fSBecBcaHdTkGqAEl+IXYwPXovmW7C86i5LgKiHfhba5nhiRkkSaY+o+tOMJtY81M2THvZCm+9XMoVSveFFtZgSI2I0lOGiGHFEI4Q97/P4hRJN8P3gvoivsMSGEXQBgZssRuaH3IeYh64p+rGanIlK3awB8D8Av0CztFRUVA8L+5Fmtf/vn3xj+PADnkwgBkzqmP0Hcsron+uGIuEnSKCJHdGsIYZ63CdiHuGJ5xSLB5ZmrJ6+hwx/FlDH0zl+tvh/L0U6GpXFipTpWos2ZqPJ1mfw+irY4V8rmyJV3NdqJ1fa634DJ1XkkcRJjaNx2Egd09IsiB/HcJF5NJGvMS34UxbB/PSFaFB+/6jZY8sf70vNj7DKV19fgaemt4nh81snR4nbNfU/D+NZkQfvT9Nx/S22laHYHORRyNxvdu/I3za44lkqftkM3TeB3OtVdAzR97qPiCW9w8NeSE/IZEkohP/rdyJWuQfsbznQ7oSU4aGEEvYac2BZC2FdKfKjohxB9DzFn0H9CtLV8xMxeHEJ4ybSaOisETM2GdxKaAabxWiw1Q55niVXM4sDWRPmvBJJo0pxTAqRWNA74jWgn+CfUjJ8LtqRJvpNKWpVIYPnpvoP2/vUqfqTnTaTJuRSNqJSsZ3fsPAIAcNfKSGiPSc/dckIkwC/5eiRIPz9zFR778diXz33V1wAAPz7qcQCAv3jvPwEA3vJX/xcAYCvi+X0TizD6qijqjZ+WCBKl1AkSAn6DNGFHXglMXJzOjfHiVGqwML+t3+yAgbKc7BvlmGPAbx7QSf9rxks+V3VSfvpQB3ScXKvwRGf2aUBsYQS93mhmfxRCmOKlbWavQlmpOgX9EKLXhxCuS//fAeAcM3v19NpZUVFxoBEJ0YLgiN4I4HIz+8+IuemBqHUcBfDCfiroSYgcEfLnPjWNRh4APAKR1SZXsA1ltlbjgfx+VrrHlUZka3dcjLbDnYpOmt+YK+9KtHMuK/fElZfvMO7atF2uVV8a72jJVV3DTDSPssOKVCaOaP3KDalFPwUAbE51HrX79njBr2PRwRge+5JU73tj8YQd8Z6vXhqdkZ6C66bUcduGdU1KkJenckMqr2aDxlKZvtvENjTcpaYBYd/r70eiEbXYX52p9U72KzkY3a4baPpcuW01ZgDtxG46lnIWWOXAZ4ZoMRs8RxRC2A7gVDM7EwA9sK8MIXyt3zqGxI/oN4jsuN8qmgNKBws/fM6s6ie8v6YUdOjzHpEIapCtilU+TkkHHAc49U0dqQNoewuTWHESUqYZS+V6d26L/Mbna7qM1U1XPj2K8DTT35V0YG9IeYMe8cGoDgwx5RDO+LNrGwPQC2IR0vHvb/02AODWtTEPyI04HgCw60XLcO2aM+JFH0z3dvi+qY3LUn/uotVsBxodFx+oe8nxWhKI1e5/PoD9pilpqevzMYtq8dJAZhIrL9rrNy3pDj30287MelZKAzIAZTUAIITwdQBfn8m9Q0KIKioqFHuwGGFhiGazxpARIq6IJ6O8/5dm4vPRzxp1r+yzdsfRaFh3rpaqgO7H8qHWOFWUktvynBFXS4pVn5bzPoJfV2VySOQoKFYmLEbDJP0y5p1+4Ldjn3QQE+M/mk60KXm+rY/luzYD70icEM5MvyVj5PdfGJ/DUA+KZjfsfGJjNHowlSvYmPT9dj0gbT4VAKMG2Fi1oqnieQ2aFBujci1TiSjnSnhRWto2Cd1sYRRTRWStD2hzVzvd/6X81v1hD5YiLADR7EBgSAjRIkzdS97nqtHE6j6xuYePM9NtlnVfdGI5ymZ0QuOEvB6BrLya1/kenHzECJr2q9ftqXItdVFjaIuLuo8ZCVMScbYBoOYv0gysTvokmuJPvyPqHLf8TZwkq/8m/v6OjeOTDs53nRDTOd58wuOn3Lsi7bj44piaEb9YeQSufzAR1F9iajkJTY0y0rS3lQKF34t97yeybmLA/tM+4nPWyXXdwG/SSeV2tBczTUHL3/ktvPleF6bpYW+BI6q7eFRUVMwbdmMp9mU4olzYx0LHnBMiM1uEuP5uDyE8z8yOQdwB69GIpr5XhxD2dK9lHyKb7BOI6YrWawWEu+cm5KEKaC9uqRJZwwGIU9z/dKAsiYDkAnzSLN6v9XL158rqleTkolQc0JCPdHzC6slo+FVPi6FBV+9YDwA4b3XkYkIa31Q4r3tzas+twEUpmugNX74/1n7m9wEAdy6N7TgCvwAA/Dv+ID39AeDW1AQa/SZHHttM8Pvtcu33IRVAw5mMyb3fQNvwoPFhut8Y3HndKooxbZukZN970YxQnzG1mo12+W162IeDYAtIWT0bzAdH9CbEmc9Z9x4AHwghXGJmH0HM6P7h0s1TkdvjXCesD4AEprLcHJRqNSM42f2ELumA+Fz1Zvbe2+pBre1X4ubr7aSShEd1EV5XxPp9MKbHSVMPv43JFB13bnosAODJJ0/dKPJzK2MS+6VIawSlpNOAc78QRbK/s0iI3rk5WtauOSF6VlM0uzURim///JkAPc9uTiU/AYNgJ0ih/C4YXBh0o0rdmaOTytPRiKGEppelOKzB0V68Y/10C2Dwre45N4o24SF0USC2oz0uZiaa7cES7Msqq4dP0JlT0mlmRyImH70oHRuiivOydMkn0OwXWlFRMQ3EZPlLW397q2jWwj8C+Cs0pP/RAO4NIZAn9Z6JU5BSW6b0lo+SXyfQVkrqCpRLNsXVUbc39vX6OlaizQlplsCc9YXncxyPf74qr/eiWYUZ4qFKVIpuftfWP0v/d+R5Y6n0cU4AJsaAG9KplJ3jiJOjOHVzyq64Ow3mw5L17J43xrYe+qNxHLYlcUJ0xUncziuf81kAwFffFR0byU2teuwvcOd1kfOafO7jUzkZja/ZEHegEVk1GZxmHyVHcSfaRgxm3aTyWK1YOSW15rMu7d67FvldZv3z/bXAVI5oZgnRiL37l2Df/oyyev/wcURz1mIzex6AO0MI109nWxEihHAhgAtjXUdkAudKe12VzKlAOxBSkTO5an0Uf2jFokihSe7XoJ26o5QTJ9cGzeJHsC7W/VIAn5F6mbenk8r1sbBU1+kAvpt++odYMMiVgavEqclC1eHzTujg0D+Kky/8JFWbnK+5Nvz+xujY+IFTo+3/zhsf22SEvCCVl6RyFyejei+vRTM82W9jUvI7ci07Do1oxgWE/cXvxefxG2uOKKCdSoT1q+7oWjTe+godn57gsQ0zC3Yl9jy4BPZgW1m9b8/wme/nUjQ7HcALzKyDOOzORPSrXWFmHGFHYrZfo6LiYYp9+xZh30Tmb1/vDRb72NfMzOx/mNlWM/uhmZ2cu+5AYc44ohDC2wG8HYgZHQH8ZQjhlWb2r4jeK5cAeA1iZH8PcDshciPdrA0a2e7RKdyjMWC+DlVgawS4ihQ+0Ra5pTXuHNDudiqTN6EtznHFpRJct7Re7tpA9l8T1KfnhnT+3pObmK+kfz30+HsBAE9KshPzTtMU/Kkkf/3PP3lrFLbRKLRfmEItbl4ZFb5fOzyKlfv4nt9G0220nrU+Icc5+3sH2lsrEXx39jVFuCvQ9Bc5ECruNRqez/OckHJgan0cS+Ut7vpeqWXUrwgoZ3+YHnaPL4aNt0WzsLuvaX0xuu9r9mxEE8WxiOzkh9F2ZjtgGIQw+TYAl5jZPwD4Ppo8oF2gu3h4qEWqZMXI/VbK5kdCtBPN5B5LZSkVhFq3TsZUx7du8Gkt+BxOJA50EjPGV3GA+3w6bAufSwIl+XU6aHZ2jgYwPC4Fu9LixRimM26PotkZv4nljz/8ODzhjnjt0xD3OLtqZdQJHZKe9+QUgP2llBPried/Fze8Ipnp2J30sJ4kmhR3fP4oflMSJ84D9XAmYViLhrBQN0RPa7V+djAVJ7q2KOHTmD1iNdoqAfXcVs98n3my5OXdH8KepQi726IZ9vZWVvexr9k5AD6Z8gx918xWmNnhIYQ7utwzY8wLIUr5rTek/3+Gqc42FRUVM8HuEeDBzBTeswhI+5q5sxcmvWu/WAPgdndMw9LwEqIDBx+trNkNqUjstuUzobE9uXggoiMlofXn8kRrRkHWS26GGEvll9Hm2lTBzmOuzs9HO0SFbdN4u9SOY9ESja5OCm3uWfZafHzqBdE3EU/4x582evkYZI9fICZVY4J9RvJvTiLntgvXNspqDuNJjohi6cWpJNfTQcOh0ImJsWfknvh+fJk1aN6VnKRu3aTZEX1uabaFlkuCHJLWMYpmTJa2qtItinLc8QzTgTwI148O0Vg5233N5hVDQogOQZy8uY+tVqZ+9hEviUq5TIrj7n8PjVnipPDxZSUCRzGORDO3+R4/jWai5CThBLgK7bxGbOvWqedH0iT/j3HgtHQNHRuTPulx6R6KWV86KqaDfe6ZKbXMwZiUhDafEQnqqSnu7Y5EkN63Me7q8f/gv8QLJwDcn5ry+VTew/dlX1DPw4m7GQ0BYl9w0wMSYwYCex0LrWCU+Kn/I1GjyLYBU8FsyL4+9azWcTOB5tupSwj7nuPCL6IknGwr4waniX3Iq0v3Zc5NH9uBtLJEzKlhafh8wSsqKiIeLPz1CJjqE18E8EfJenYagF/PlX4IGBqOaBxxRfNiCq0hGltWwjLkc0MDbfaZq7NfADRzoi5Fq+X8OKau7kATKsAVVjmw09GsoLqianIuHp+CRnQgd8G2nDX12onETR16UrOB7FNiseSk6PNziPTN+t0pheI7YnHLscBhIXIOz7wv+gtdszyq/Kjw/odT3xov5m54vwSSDtyFeKn4Q7Foo7twlVzD8Az6Tb0ylRvcNQzDIIfC7zEmdbEhZ7vrJEvBpOK7V4gO0Bg8NFtDKZUI0M6rPk3sRjfRrCv62NfsSsQdOLYidvzrZtbI/jAkhIjJ870jok+rADSvwg+vQamL0RAv5TA1fWgnlcehIU5qclXP7n5kf04Ssum6XfF30LDrHOScBCQqJEwUDXehmUAUJdjmTirXpzIRqHt2Ak9PhDXRu5/eF9v0pOU3AGhSevz3pZECvWvzfwMAjITD8UX8HgDgNTdFM/qqU+PzPoOXAWj2N8PByQ91zDLJKFXU5QXL3fHdci37mITiPrlnB5r+0syMR7trgOb78ZsAzVgiAaJNhYSdFstuhMMHxALNd/QphnWzg2xwQW/MghD1sa9ZQMxFPS8YEkKkOqIRlDe942DlICW82VQHqxIM1rkd7Qh2EoKSgtGvnrrXPQeehhV4Isd293JL8MpPEjY+h5wD7/Xe3gCwvMkhnebY+La4u8bG46Me6U0ppyu3BvrEC6P/0gj24f7UtqtPfSoA4J6knH5Wcm34C7w/VrotJl3DDXC6fuXaOHHJHWxyv/ObchJrMCrrYn96ZbUSC44L9gnDYqhn2or2YqK5zbXOEbR1hQoNtN6FfCDsDPAgcimrD5RoNq8YEkJUUVHRwgTy4Wqzyy4yEAwJIdqLyAF4cUgdyDQj3vpU0gHRb0Wtu0EQuvKOotEXaTpX5XbG5XegLb7xeWo69nqEsVTyfcgV3CfHxFq0vbo18JcihkssP5aewwU9eZz8Yixavp52yDWpZbHO8WSx2odF+D18K93ylNTS2KarEC1sk3omn36J3XNW6vurqKSiLodmex9rRhGJ/gIqkvFaelP7PFOsj31OroN983epPNqV5NJ02/K1ct7r6fj9NdsivwnrpA7pRLSdIGeI3ZhM5TsFs4ulHQiGhBBpYjSgzApzkqtLvzenl7YgIrz7vyZY4+ziAFa/Je/Sr4SNI6Tk4n8fGtGEZmZO2LMxFRRPjkMzMSmSkdCqcp5tva8Jej04lSl16/2PiQ4/lz095pD9y0PeB6AJhj0Kt+OdaRKfl7K5/DRtpHj7nmjt/fWPUkg/fYfuQpOadhJ8P34Tvo8XOfkdNqRyfSpVJ0YC5N0t+P1V17ZGSp8sr6RfVMKU81XThYnfRF1OtqM9dmZIOWahI1poGBJCVFFR0cJuIJsnvxKiuYIhNjWXF0ZfoZTtrlsMmgrVPhBSd4MgN1V6DlfCde5avZesfi5fEf/nqs/nj6XSB7uyLraFViautBQr75x6fvTkuGcl0Dga8jitsPffENmZzz/tXHiswL14IIlpV+OZAIAv7YgxZcsO3TW1zhWp/BSAtNM07iEXwX4i50IRhr8/Gw0nxN94rI6H7Kv1aLgj5Xo1y6Omc/H1UflNro3cWm78aXZMHUsU47y7SbexOA2MA5lMsXlxbYFjSAgRoZYkoPE/0TzUM9HYqXVrFG3rnLLeunsD792CduQ85NpNcrwZbd8jEq2PyDEny1fQEBw1B5Ogym4i48cB46n996ZJfp408THR9L71gfi8sw6JBONmPH4yj/U1D0RxbX8KvPz1txM1Y9wqpUoDcA/FD/5IkPAoYbgUDSHgtx1LJftck9NdjoZoqf5FhzrrYh89A+0c4xtSqb5AnhDS25ttopsFCawsAlNQSifbJ/Yj70W9f2bVDRJDRogqKiom8SDyHFEVzeYKVFbnwNVSc8lw1fIilPqYqOVNrRmr3LXK6uvmfrqqHY223xLfgauwWuD2oq3cJOez1l3jcSLaaSrUa5h1kltwXursEiqTV6Qy+QDdP5FEtJTAaNXv/Bx3/iyFIH0v+Qk9Nb1HB1NLb1sYSaLRBC1gVKzTO5p9Qk5iDWJKcw++Jzm8XICw9jkx4a4B2jF8O9BwMxT1+G1Pd9cAU+MatR41khD8fovRfMNOKtUS2iceRD5Iq4pmcw0OolvQNlmrDiAHTsRerLDf45xdpAnqOJE08VUu9SjN55wkqmvwTpJ8HnUKfD/N08z3vNa1YSyVFH9UPOXvri5aXXgLHY3PSSXn1bmxeOCBQ4ANiQBRrzSS5IO1I1PP/4d79ETq85FETCY42SkOqRV0MZp+U7Fb9XNcFCbQJjT6fcZSyW/M525BY2JnHRwn/JY6XcbRdtEgSg6pflNGtqUUud8De5DniKr5vqKiYt5Q8mCpHNFcwydL15WGHFE3BWCJW2I3qPMiMDWNK9CslqUE+DlOqN98MxNopw3V/dPUXrsMjdL2JrmGHBE5IIZRjGPSOfBXiYMjI3d9KilJECui8vr+yw5rxDeuvJelOANayxg+ckMqw15MijsTfD/1rWJfMaOiv4b3aFoQcpgU9+5DIzbRUKCWSnJiY6nMcSOqUCfUedWH8fQKdvVjz+/aMgtUZfV842DEgdRJxz6+igNP04cS/jrdi57ERCP5/WTXDft8Wlf/3BxIGEr6LWR+13xEuvMHJ5gnMiREfC+KCSUHzpdisl/uSqcuSiWtZ5qT7PNJHFvrrk2R+5PxTqyLXUZvAuxE845s9w/db770yTs5qb8hx+wLvi8f5Am+7i+v3zhHBDSqXr+tOoaOZK4ZSyUJnC6M3s2ci0w/ObQyyDkzAlVZXVFRMY8oiWCVEM0VmDzfpwTVjIxUWKrI5CPdySmoqMRjdgfZgA7aIp5yYt1EQeWEeC9X59z+arqjyJg8z4czEBS9uMqfLteSw+M9V7i2JFHz8PScy9LplAZ2UgzjHpdXoYm0uCGV97pqATScqzcOkBslpzAq17Dv2WaKxEBj/VOLFK1c5K7G0XCMGhdG8L3Zvwy3GEczPvR5at300PSxmvVBx8UE2tbaXhxzAQ8iZshRVB3RXCFgKiu9HW2rBD+qTmTCp2wg8dABwGOy+jnioik8eI3GJ+XAezkQ+Q4c6HvR3tJIiQlL6qo2umup4+Jk5/PUu/do1+50zcFj6ThNpHXp2g3pdMoUi6ejUaF8O5WMVyPxIvG8w4u1pcnGCazxV3ejLcLopGd+IJ+mQ3a1bRFyFVN9zmriK+gPE5n6SrvCepcDXQhnOA13I0+I+uSIzOxsxL0GFwG4KITwbvn9tQD+bzQv+aEQwkWYA8xZqlgzO8rMrjazG83sx2b2pnR+pZl91cxuSeWhveqqqKjIYKLw10fOajNbBOB/IsbSHA/gD83s+MylnwkhPDH9zQkRAuaWI5oA8NYQwiYzWwbgejP7KoDXAvhaCOHdZnYB4ibEb+teFTM0ejFIV1i1npHb8CleddcHQlcoimaeu1F2Wv1HlBNa5q5V50cVBYnlyK/QQDsBl7fO6WaBvPaFqaQlyu/FJrtbcONDcmBbfQwWGvHrMjRMJ5XTtJaRMxrXNu5E06fKZbIfNVRmB5p31TS9KuJq+hN/Lftcd25h3/MZvu97QceCP1eyhOU4pX646C54EHmi0x9HdAqArWl7L5jZJYjeYzd2vWuOMGeEKCXaviP9v8vMbkIc9eeg0TJ8ApH570GIFiMOWD8wNdZLCZNm1VuMdm4aBQfGhsxvHOwcaLrBosYabUQzUUrsunp435epj9DdQvwg5oTns2ll8qZw347j0A7wZDBtqvd2r6NxbTwRzW5Xu8TCN873ZXvGUrkajc4mN4mBhrrxux6JtoWS78znqEVzl6uH9/CdmTJELW+E6na6wS8k+p14TF0V+8i7CWhqmRliD/Km+kiIeu1rltu3TFYfAMCLzewMAD8B8JYQwu2Za2aNedERpR0ln4Q4O1e73QB+ibbLMO85H8D58WhmW/JWVDyk0T3W7EDsa3YFgH8JIew2sz9GZBzOnGWdWcw5ITKzRwL4LIA3hxDuM7PJ30IIwcxy6jYk6n1hrOOIEFc2Lw5pYitaimi6GUulX1V5j4aHaIQ97/WKafWD4WqmidA3yfUeanEjciKjcg4aFtBJpRc5qGRdJr+R9WcdO9DETfFaruAk+hencn0qk8i2eSUaUVBj3MhV0BfIi4ra17yX5++U41G0h+eLUvkBOe+TnW2RcwSdOZXTI9ai7ZNG6EKo3DjbCzTfSzMuePA+HQfTxG7MJui1575lIYS73eFFAN47vQb2jzklRGa2GJEIfTqEcHk6vYN7aJvZ4SjLLRnkLuWAu0LOq2l8Au0Pz0mh1hivKyhtJaMWKcLfq3okTlgSpPvk92VoBr3qDUj41PJ2GxpxRJ83JnVwkoyiIaD0Sua7s34SdtbhiU5HzrEu6jKZvoP3+vQc2ucEn886v4M2PiLHuTQrmuuJ4OTn3NKFLOcpr6lE1NkUKLtvjErpM0WOpf85Hng8TUwgbzXrb4PF7wE41syOQez8lwN4hb9A9rp/AdpBfwcMc2k1M8TtNm8KIbzf/fRFAK9J/78GwBfmqg0VFQ9pMFXsDDZYDCFMAPhTAP+OSGAuDSH82Mz+3sxekC7782Tx/gGAP0c0NM0J5pIjOh3AqwFsNrMb0rm/BvBuAJea2esRl/OX9l+lZ3+5onXkGuVgcsnUFKUtfNah7dejqySvVW5kJdocnIqAuTaXFJglJ8xlaIc88LNSQaptOw4NF0GxkJzRmLSRmxmSoznOXcMFsuT75DlQ3sPnkYvSvcl4j+97FY26iT/K9Wr4hHKaOQ5KOdaSqDaBNifEvi9tO7UczZgdS6XmxO4TuxGTzin63HI6hHAl4kaK/tzfuv/fDuDtM2vc9DCXVrNvI99NANJ2D9OG/6gq1lAvsUOu5QRa7M5Rx8EJUzLbbkE7jYRucKiEidd7LrZkCVOxzu95palDFOoACLT3DNsuxzShb0LDhbP/6KXI96FOhfdy8o2hsSrqpKf+peOuBeI7sT46XVJEUmJztJwH2v3XT5yfitsjUirh89C8TvqNu6kI1CKruiPvjDvD9B/Eg8jPsKzWdWFjSDyrKyoqshhCopPDkBEiv5urKhc1tkh9dIBm1dIQDl3NvKikq9aElMp66x5m/toS6OOy012rynYNXfHpQVSkVJ8YvgOvGwPwz+l/ciLsW4o7yq2RA7wKbXFY04xQtGEbvZOfJpBTTsWLpiUrYwlenNOhTS641L8jaG+UUPrGniPjPaUMnjyvYT1A870058rDD0NCiOgJ7MWdkmOcmodJbFahnEOIrDi7I5fmgb9pPp3peMWyDg48tp2E0RMsJV4qwvjE8SQaFL3U0U8JlSeu7B8SoJKIwecdi0bk0vSoms+J1/lvpMGtJEh8vtexaTAvwbbl+kgJKN+ntGmiN9mXAqeJnEio468kUo9LCeSzeU4He/FQ2ep1SAhRRUVFGw8in3Rv+MLvh4QQ7UNckakMHUFbPOAHoYihK1I/2/t2U4by/5Lb0zr53SvAlTvjiqspTHw9XPXVwqaWnFFXLzmSUlYBHyKhTpRsY0eOeR3b+EN3z+ulDuVU9rqSbaTYQ5TCbVai4eym4Wo2eS3fVb87z+fEvZJlVbmcXIyYbiOuObHVDw1oj+HpYhyVEM0rGPTq5XolFhwcHCyqU+kWoKjpYH0AJu8r6SnYDk28DuQ3UPT38FqKJSeirZNSq50Gei53bWM9mtLUW2yAqZsCaMwc5LwSz6PRTHZtq05GH3A6Ib8RbBtzAG2VMge/y4qvwxMsFZlJTJYXfgfKTpYcF7pI5K7tBb+49XtPCZUjqqioGDgmkNcR9elItIAwZITIO+aVop4JVe76lVhXEVXucsXNRUyXxDce57Iu+m2QgEYJqwrizWiUqB1pg+ZThjt+oVyrPk47pFyN9jtqaIfuesr3XI2Gm1DLpHIxpXAVoOkTKuoZJ8f3X4M2x6CpO1gv7+kmwnVSeVyXa8hN67TQ8BTvLOudX4Fm3JUsfv69OJ6PxMxQRbMBwQ80te5sL5z38UElC5gO4NzH1XuJkue2f54SJ73WP08Hrnp06/N3opmQfJ9SwCVjwL6DtrcyJ46a79k2TuCtaCZdaSOBZVLuRHuikgApkaIObCfai8hYKvltyQ3kxLi1cq3uUab7nwFtiyRREqF8kK26DWhMoBePqUYojY9+UUWzioqKgaMSogHBWys0EnoslVyJlFNZi2bVUguY7hbr5W518+eHV9Zbf8/5cpTEuW4WNnJ4GtJCrELZD4X9RWU8HQ7vQ1u3UBJr2CdekV+KbCeoEPaR9pqsXrk41sl2eGugWvRUJMyl4yhdy+/DhPv+m5Q2SmBfKffmQ4D4nTTftWb0XI6mf36I2aFEiIZvG48hI0T8gN6buJTbRTe624KySFYysU6gvKWwpgHJDQhCCZCKRTsz1+pvpXi43HnVzWjyeSBvTgbaifapV1Kxy9eh0PiunWibwEs6L8Kn2tUNChigqyLNSpTzHvF9dA+0HAFW8aobqGOjgyZFvo7Udacr2celzRj7BZNUK6qyuqKiYt4wjny2htLOiwsXQ0KIGOLhrViaubDE4nt/mdIqWLLAjWBqyIGvn9C6uCKPu3q0ftZB5TGTup2Idu6pEsfgo8iVa1NlKzmJXGjEWCo7qdRVWh1FNyOv6AXaTpB+pVcH01ImAp+DXMW4XXJMeO6q5KdE7oZcoY6P5e5/zcZZiuT3UMWz+hzxntPRPXvjdFAiRFVHNEeg/O6JjhKAktMgB8hJaOtSeqWT8OdLBEGd+DgwdqFsUeOkIwHyei51EtRgUNVB7UXZC5vwIi0wNZZKRS/N46S5jdYBeKDwnDVyrT+vxFEdM3PxXSQ4fk8wDyXAuWT2vIZ9tEN+9xsXaB8TuouI3xxSRWedUioC54jQTKfhg5gNIepjX7OlAD4J4MmIMu3LQgidGTa2K+YsQ2NFRcVcg8pq/eutrO5zX7PXA7gnhLAWMVH4ew5Y0wVDwhFxOyGv8CNK2fN09fQWCg2bIOjMx1Uzp6zU6G71MaG/zU6UFdjKufAdNO820N6nrZRtwENFmtx2yL0+vUai8/rthWcCjViZU85TXKPzHp07Sz46ueyY5KLYx2p1Wp2pT/c18+lTeA+Q/9bKDRI+RIhjiSIZvxffvbTzsK9nppiVsrqffc3OAfBf0/+XAfiQmVkI4YBnQRoSQrQXZV1MaaNFxS40r1vaw0pFt5XIWz34G9CedMx02K2tatLl4N2ONiHNpY/wyHkgqw5F6/Qeweo1XAoK5fstxtT94n1ZSqPRQdtkvbZwLbEVbesfJ65uK+2dBUvidskD2o8FtWqyfvYFCWDHtatE6NTSR/j9+RbLtdPBsWfHtue28dgF4GlvAdBtO6F+9jWbvCaEMGFmv0Y0Wd41gwZ3RRXNKiqGEicB+G7ht40ArvlvZnad+zt/Hhs3bQwJR7QEU1f+ZWic9ErOfDmW2IcpAOUkYFzFFrtrVWnMVUyd6rjy+3i4UrI29ZPx59Sf5z75ned9yEBJuZvzsVJrUq/tj71jZemakt/NLrT7QPcfy4loqgD22yH5e7xzqXJCPisl0HAz3XKCaxwe+5PcnE9lorGBypHltoWaTlqTEj67CFizL+4FvsKd/xmAUSbB/9vcnQk99zVz12wzsxEAj0LjiHVAMSSEaD+miiaLUc7lopkTvWWFOgyNKVojJa/bhfZE0YBSFZnGUrkDZbFRJ/06d42mLeGgVcfJ1e68DnYd6Do57kSbgPba/tgvAiU9XD+ZATW/Eb9BLp+TOq2yVB0RCWBOlN4l17Lf+P45QsjfSKzUJSS3cIxJWwi+12zEsDZCCPvNXolohTuHZxE3Nvj57/RRRc99zdBs/fUfAM4D8PW50A8BAxLNzOxsM7vZzLaa2QWDaENFxfDjnw+KhI0qm5sBHIYQgqaTaKHPfc0+CuDRZrYVUd80Z3PV5ojAlR8YzYY/AfD7iAqy7wH4wxDCjeV7jgjAH7szI2jHhc0kT6/Gi+XyXPda9UvhG3vR5pY0AVtuNSYXwJWbdaj/C8UDH6WuKMVO+VW7pNzV7Yp8n5REi5K44y1gveCfW8qkoNYmcou3ufvVqECo35RXvCsHqw6a+t7+vcj5KPekbfTPUXH/bdeHEJ6CacDMzgCe8A3gxYg77f7iqBBCzmKyoDEI0awfs2EPeLOlWqRUZPK/6yBVNlon7nK0PXCJkhdxN9abBKi080MOpU3+PFSXsrpwXTePYJ30Jctbzpk0t/WzR7f39Hvd63P1vUqW0lz9Gpys6U00i6R/zhr5TTNC8v09cdX0v2elkoQvF7Rc8rzvHyGEb5odi7iBwdEIYfvQESFgMKJZzmyoOTphZudT41/25K2oqAC2Pjkucv+hTm5DgwWrrA4hXAjgQgAws18B//U3mAP/hVniMNQ29YOF2CZg4bVLnY76QghhE8q7Kg8FBkGI+jEbTkEI4bfM7Lrpys9zjdqm/rAQ2wQs3HY9HDEI0WzSbGhmSxDNhl8cQDsqKioWCOadI0qu4jQbLgLwsRDCj+e7HRUVFQsHA9ERhRCuBHDlNG+7cC7aMkvUNvWHhdgmYOG262GHefcjqqioqFDUoNeKioqBoxKiioqKgWPBE6KFEJdmZkeZ2dVmdqOZ/djM3pTOrzSzr5rZLak8dABtW2Rm3zezf0vHx5jZxtRfn0mWyflu0wozu8zMtpjZTWb21EH3lZm9JX27H5nZv5jZwQuhryoiFjQh6jOd5XxgAsBbQwjHAzgNwBtTOy4A8LUQwrEAvoY5DArsgjdhasb99wD4QErveQ9ius/5xgcBfCWEsA7A76b2DayvzGwNgD8H8JQQwgmI1tqXY2H0VQUWOCGCi0sLIewBwLi0eUUI4Y7kvYoQwi7EibUmteUT6bJPADh3PttlZkcCeC5itCPMzACciZjWc1BtehSAMxAjtxFC2BNCuBcD7itEC/FoyqtzCIA7MOC+qmiw0AlRX3Fp8wkzGwPwJMQ0eKtDCHekn36JcrTpXOEfAfwVYsImIKbxvDeleAAG01/HAPgVgI8nkfEiM3sEBthXIYTtAN4H4OeIBOjXAK7H4PuqImGhE6IFBTN7JIDPAnhzCGFK4uuUMGrefCHM7HkA7gwhXD9fz+wTI4gbn304hPAkAL+BiGED6KtDETmyYwAcAeARAM6er+dX9MZCJ0TTjkubK5jZYkQi9OkQwuXp9A4zOzz9fjgOTA7QfnE6gBeYWQdRZD0TUTezIokfwGD6axuAbSGEjen4MkTCNMi+OgvArSGEX4UQ9gK4HLH/Bt1XFQkLnRAtiLi0pHv5KICbQgjvdz8xlSZS+YX5alMI4e0hhCNDCGOI/fL1EMIrAVyNmNZz3tuU2vVLALeb2ePTqWch5poaWF8himSnmdkh6VuyTQPtq4oGC96z2syeg6gLYVzauwbQhqcD+BZiNnXqY/4aUU90KYDHImbRemkI4cAkJZ5e+9YD+MsQwvPM7LcROaSVAL4P4FUhhHndg9jMnoioQF+CmM39dYiL3sD6ysz+DsDLEC2g3wfwBkSd0ED7qiJiwROiioqKhz4WumhWUVHxMEAlRBUVFQNHJUQVFRUDRyVEFRUVA0clRBUVFQNHJUQVFRUDRyVEFRUVA0clRBWTSPuev9kdv4u5lyoq5hLVobFiEimzwOUhhJPN7CAAtwA4JYRw92BbVvFQx4Ld6bVi/hFC6JjZ3Wb2JMQ0Hd+vRKhiPlAJUYXiIgCvBfAYAB8bbFMqHi6oolnFFKQsB5sBLAZwbAhh34CbVPEwQOWIKqYghLDHzK5GzF5YiVDFvKASooopSErq0wC8ZNBtqXj4oJrvKyaRdibZirjbxi2Dbk/FwwdVR1RRUTFwVI6ooqJi4KiEqKKiYuCohKiiomLgqISooqJi4KiEqKKiYuD4/wGj3FhilK9OpgAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "dataSet_crop.OD.isel(runs=[0]).plot.pcolormesh(cmap='jet', vmin=0, col=scanAxis[0], row=scanAxis[1])\n", "# dataSet_crop.OD.plot.pcolormesh(cmap='jet', vmin=0, col=scanAxis[0])" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "## Remove the background" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "dataSet_crop['OD'] = dataSet_crop['OD'] + 500" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "dataSet_crop['OD'] = imageAnalyser.substract_offset(dataSet_crop['OD'])" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Test Fit" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fitAnalyser = FitAnalyser(\"Two Gaussian-2D\", fitDim=2)\n", "\n", "params = fitAnalyser.guess(dataSet_crop.OD, dask=\"parallelized\")" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fitResult = fitAnalyser.fit(dataSet_crop.OD, params).load()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fitCurve = fitAnalyser.eval(fitResult, x=np.arange(100), y=np.arange(100), dask=\"parallelized\").load()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# fitCurve.isel(**{scanAxis[0]:np.arange(30), 'runs':range(dataSet_crop.OD['runs'].size)}).plot.pcolormesh(cmap='jet', vmin=0, col=scanAxis[0], row=scanAxis[1])\n", "\n", "fitCurve.plot.pcolormesh(cmap='jet', vmin=0, col=scanAxis[0])" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "a = fitAnalyser.get_fit_value(fitResult)\n", "b = fitAnalyser.get_fit_std(fitResult)\n", "data = combine_uncertainty(a, b)" ] }, { "attachments": {}, "cell_type": "markdown", "metadata": {}, "source": [ "# Get the Ncount" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "Ncount = dataSet_crop.OD.sum(dim=(scanAxis[0], 'x', 'y'))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "Ncount.load()\n", "\n", "fig = plt.figure()\n", "ax = fig.gca()\n", "Ncount.plot(ax=ax)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fitAnalyser = FitAnalyser(\"Lorentzian With Offset\")\n", "params = fitAnalyser.guess(Ncount, x='sin_mod_freq', dask=\"parallelized\", guess_kwargs=dict(negative=True))" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fitResult = fitAnalyser.fit(Ncount, params, x='sin_mod_freq', dask=\"parallelized\")\n", "fitCurve = fitAnalyser.eval(fitResult, x=np.arange(40), dask=\"parallelized\").load()" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fig = plt.figure()\n", "ax = fig.gca()\n", "plt.errorbar([1], [1], yerr=[1])" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "fitCurve.plot.errorbar(yerr=fitCurve)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "np.ufunc(fitCurve)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "env", "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.8.10 (tags/v3.8.10:3d8993a, May 3 2021, 11:48:03) [MSC v.1928 64 bit (AMD64)]" }, "orig_nbformat": 4, "vscode": { "interpreter": { "hash": "c05913ad4f24fdc6b2418069394dc5835b1981849b107c9ba6df693aafd66650" } } }, "nbformat": 4, "nbformat_minor": 2 }