{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "AirSeaFluxCode is developed to provide an easy and accessible way to calculate turbulent surface fluxes (TSFs) from a small number of bulk variables and for a viariety of bulk algorithms. \n",
    "\n",
    "By running AirSeaFluxCode you can compare different bulk algorithms and to also investigate the effect choices within the implementation of each parameterisation have on the TSFs estimates. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Getting started\n",
    "\n",
    "Let's first import the basic python packages we will need for reading in our input data, to perform basic statistics  and plotting"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# first import all packages you might need\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import netCDF4 as nc\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from tabulate import tabulate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### AirSeaFluxCode examples\n",
    "\n",
    "AirSeaFluxCode is set up to run in its default setting with a minimum number of input variables, namely wind speed; air temperature; and sea surface temperature. Let's load the code, import some real data composed for testing it (Research vessel data) and run AirSeaFluxCode with default settings (output height 10m, cool skin/warm layer corrections turned off, bulk algorithm Smith 1988, gustiness on, saturation vapour pressure following Buck (2012), tolerance limits set for both flux estimates and height adjusted variables (['all', 0.01, 0.01, 1e-05, 1e-3, 0.1, 0.1]), number of iterations are ten, non converged points are set to missing and Monin-Obukhov stability length is calculated following the ECMWF implementation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from AirSeaFluxCode import AirSeaFluxCode"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "inDt = pd.read_csv(\"data_all.csv\")\n",
    "date = np.asarray(inDt[\"Date\"])\n",
    "lon = np.asarray(inDt[\"Longitude\"])\n",
    "lat = np.asarray(inDt[\"Latitude\"])\n",
    "spd = np.asarray(inDt[\"Wind speed\"])\n",
    "t = np.asarray(inDt[\"Air temperature\"])\n",
    "sst = np.asarray(inDt[\"SST\"])\n",
    "rh = np.asarray(inDt[\"RH\"])\n",
    "p = np.asarray(inDt[\"P\"])\n",
    "sw = np.asarray(inDt[\"Rs\"])\n",
    "hu = np.asarray(inDt[\"zu\"])\n",
    "ht = np.asarray(inDt[\"zt\"])\n",
    "hin = np.array([hu, ht, ht])\n",
    "del hu, ht, inDt\n",
    "# run AirSeaFluxCode\n",
    "res = AirSeaFluxCode(spd, t, sst, lat=lat, hin=hin, hum=[\"rh\", rh], P=p, cskin=0, Rs=sw,\n",
    "                     tol=['all', 0.01, 0.01, 1e-05, 1e-3, 0.1, 0.1], L=\"tsrv\")\n",
    "flg = res[\"flag\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "res is the output of AirSeaFluxCode which is a dataFrame with keys: \"tau\", \"sensible\", \"latent\", \"monob\", \"cd\", \"cd10n\", \"ct\", \"ct10n\", \"cq\", \"cq10n\", \"tsrv\", \"tsr\", \"qsr\", \"usr\", \"psim\", \"psit\", \"psiq\", \"u10n\", \"t10n\", \"q10n\", \"zo\", \"zot\", \"zoq\", \"uref\", \"tref\", \"qref\", \"dter\", \"dqer\", \"dtwl\", \"qair\", \"qsea\", \"Rl\", \"Rs\", \"Rnl\", \"ug\", \"Rb\", \"rh\", \"rho\", \"cp\", \"tkt\", \"lv\", \"itera\". Let's plot the flux estimates."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAb8AAAEUCAYAAAC7yDhfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABzsElEQVR4nO29eXwU13bv+1samIQAgZgHCTBgwMyjjY1tJoEnwDa28cSxwXonNyc3uS/Jyzk39yXnJu++l3uTm3ySnCTHwsYjNjbGzIMBGzOPsphBIIFAAiGEBtAstbTeH927XF1dVV3VXdXdkvb389FH3dVVe6jatddea6+9NjEzJBKJRCJpT8RFuwASiUQikUQaKfwkEolE0u6Qwk8ikUgk7Q4p/CQSiUTS7pDCTyKRSCTtjoRoFyAapKamcnp6erSLIZFIJBIXyc7OvsfMvfV+a5fCLz09HadOnYp2MSQSiUTiIkR0w+g3afaUSCQSSbtDCj+JRBI2WVlZKC0tjXYxJBLLSOEnkUgc4cGDB9EugkRiGSn8JBKJRNLukMJPIpFIJO0OKfwkEolE0u5wTPgRURIRxTuVnkQikUgkbhGy8COiOCJ6nYi2E9FdAJcBFBPRBSL6eyIa4VwxJRKJRCJxjnA0v30AhgP4DYB+zDyYmfsAeALAMQB/R0RvOlBGiUQikUgcJZwIL/OYuUl7kJnLAWwAsIGIEsNIXyKRtCLkxtiS1kTImp+e4AvlHIlEIpFIIk1Iwo+I5hPRaiKa6Pue6WipJBKJRBJRmBlZWVnRLkbECFXz+08A/hzAm0Q0B8DEcApBRAuJKJeI8ojo1zq/P0VE94notO/vr6xeK5FIJJLgtDezdajCr5SZK5n5zwAsADAt1AL4lkf8G4BFAMYAWE5EY3ROPcjME31/f2PzWokDtPaRYXNzM3Jzc6NdjDaLG51nTU1Nu+uUJZEhVOG3XXxg5l8D+DSMMkwHkMfM15i5EcA6AIsjcK3EJq29E7p9+zb2798f7WJITMjKysK1a9eU72vXrsWVK1eiWKL2R2t/z60SkvBj5s0AQESpvu//GkYZBgIoVH0v8h3T8igRnSGinUQ01ua1IKJMIjpFRKdk9HmJJHbRvp8NDQ1RKonECerr61FXVxftYgQQboSXNQ6UgXSOaYcePwFIY+YJAP4VwCYb13oPMmcx81Rmntq7t+7GvpIgiBFhLI4M6+vrsWPHjmgXQyJptbj1Xm/atAlff/21K2mHQ7jCT0/42KUIwGDV90EAbqtPYOYHzFzt+7wDQKJP6wx6bVuGmWNSEEWDe/fuoaioyNY1e/bswYULF1wqkUQiAbzztrGovYcr/JzoeU8CGEFEQ4moA4DXAGxRn0BE/YiIfJ+nw1vuMivXtmV27tyJbdu2RSy/tiZor1+/jqtXr0a7GBIHycrKanPtVOIO4UR4ARzQ/JjZQ0S/AvAdgHgAa5j5AhH90vf77wG8DOAPiMgDoA7Aa+xt4brXhlum1sKtW7fki24D3/ipXVBaWoqePXsiPr71xZoPtU3LdyE81NMa7eFdCVf4/caJQvhMmTs0x36v+vw7AL+zeq3EHVpL5+LxeFBWVoa+ffsGPffu3btt8mXfuHEjZs2ahbFjxwY/2SFaS/uQSIAwzZ7MfN6pgkhaD7HeyZ0/fx6bN2+2fH5jY6OLpYkeLS0t0S6CRBKzhKv5gYimAvhLAGm+9AgAM/P4cNOWxBaxLvQEstNv3YTbztqiJh8J3PLmjtV+I2zhB2AtvKHOzgGQvU4EEY3q/Pnz6N+/P3r16hXlErVOYvXllEgk7uHETu6lzLyFma8z8w3x50C6EoscOXIEOTk5rufTWoSE0ahfagPuEu32Ee38Ja0LJ4TfXxPRB0S0nIheFH8OpCuJUdpaJ+NEfbKyslBVVeVAaZyjrT0nwFunYNFC2mK9I4Hb9622ttbV9O3ihPB7B95dHRYCeN7395wD6UpijLbSqWzdutWVdB88eOBKuu0F0TmatbNbt27hs88+0z3HavssLCxEfn5+aIWUhMzu3bujXQQ/nJjzm8DM4xxIRxLjtBXhV1xcHO0iSHRobm4Oeo6IFFJUVITBgwcHOVuf77//Ho2NjRg+fHhI14dLc3NzTK6/dNvhJdaivDih+R2T2whFDzmPFYi8J16CdWLMjIqKigiVxlnC8eiN5iCuqKgIH374ofI9KytLeidHCSeE3+MAzvg2lD1LROeI6KwD6bZZbt++jU2bNkW7GJJ2zu3bt7F+/fpoFyMAJ5Y6xCp688KxLPxOnz6N7du3Bz/RAmbP5cGDBzh37pwj+VjFCbPnQgfSaFcUFhbi7t270S6GbUTjLSkpwYABA6JcGvsYaYSx3Fm6icfjiXYR/FA/h5qaGuzbt8/SuWbHrF4bKdRtMFZ3SVGX59q1a7h3757reV64cAHnzp3DuHGRm0ELWfgR0V+Z/MwA/jbUtCWxiXgptm3bhszMzCiXxhhp9vRy7949HD58GLNmzYpIfk524iUlJbh9250NWszKWVdXh86dO7uSr1E5YkXzizUh7DbhmD1rdP4YwEoAfxF+0Vo/paWlaGlpwa5du3Dx4kVkZWVhz549ruSl3v1a0jY4efIkfvjhh5Cvz8/Pj7ktm+7fv+96HlbmOvVoaGjAZ5995kaRTMsRbeHX1NRkGOLPCYFoRcONhuANWfgx8/8WfwCyAHQG8C6AdQCGOVS+Vs3GjRtx9epV3Lx5E1euXAHg3UbHSSKp5bjdQJubmyPWEag9C2N1xJubm4u8vDzb11VWVjpfGIf46quv0NzcjIKCggBnGztmQCfNnvX19QAiI4T0zJ7RFn5bt27FV199FbPvgVuE5fBCRD2J6P8BcBZeE+pkZv4LZm59E1ouIRq2NMUF56uvvnJ1LZD6GbTlXd+3bAl/S8vc3FwcOHDA1jXHjx+3dB4zY/fu3Vi/fr2hJhip96W0tBSffvppRPLSEqrwy8/Pd9Rn4MGDB0EDBzhBrAWBCFn4EdHfw7uZbBWAccz8W2ZunX7TLtKWRlNu16W6uhplZWWu5iFQax3RfkY3b95ESUlJ2Oncv38fZWVliiYTDLN6X758GZcvX7aVfyjruERggI0bNyodsNvPQ51+JDr9YOUwEn61tbW6v33//fc4dOiQ37G6ujrHNmZ28v7HshNSOJrfnwIYAOC/AbhNRA98f1VEZCvUBREt9C2VyCOiX+v8/oZvGcVZIjpCRBNUvxX4llecJqJTYdTHFdz26FK/HNHuxJ0gGhH9o33fdu3ahb179yrf8/LyQtqRfOvWrdiwYUPQ85qbm3Hs2DHTrZzcWIStVx/RfktLS1FeXh5SGqGcI9rIrl27dK8rLS3Fd999FzSdcAgm/D7//HPs37/f9FrBxYsXsW/fPly6dMl2OdzUsn/66Se/79F+19SEM+cXx8ydmTmZmbup/pKZuZvVdIgoHsC/AVgEYAyA5TqL5q8DeNK3TdLfwjvHqOZpZp7IzFNDrY9bxNLDDoeGhgacPn062sUIitW1QrFmhlaXx472e/To0aBOVE1NTaiurla+V1dX4+zZs/jxxx8Nr3EzAon6ndCbe9V7Z27cuBE0AoyT79qNGzdw44a78fmF0DMze9rV5kpLS0MuTzjzqEZo23JZWRmOHj0aVppOEY7ZM2jvYeUcANMB5DHzNWZuhNdhZrH6BGY+ojKpHgMwyG55I0lBQYEyAo/kWh4387hx40ZIo0q7hFuH9hZf8+rVq4oTldHrtn79enzxxReGaeiZ/uLirHUNBQUFlj04QzWB3bp1C9999x3Ong2MndHU1GRJY9SivVfqPD0ej+G9rKystDy3GQzhWRqKw0t5ebmuedtKiLhguNmP5Obm6g5QW5vZcx8R/RERDVEfJKIORDSHiD4BsMJCOgMBFKq+F/mOGbESwE7Vdwawm4iyichw8RkRZRLRKSI6Fc7oyAo3b95URjyxovkxc9QXNWdlZTkytxWMUMye9fX1+PLLL6M6BxQuRgJLrfXp8f333wccs3oPd+/ejWPHjlk6Vw/1+6EeKKrzLyz0dg/q9ivOPX78OL755hvDNI0gooC8xffs7GzD665cuYIzZ84ETd8sXy2hCD9m1jWJxlrYt1jp//QIR/gtBNAM4Esiuk1EF4noGoCrAJYD+Cdm/thCOnpvme4dI6KnEbiOcBYzT4bXbPqHRDRb71pmzmLmqcw8tXfv3haKFTpGL3Q43Lp1C1lZWmuvcb5arl69ijVr1oRVBifQ62idQIx4rd5nrcv5p59+iqqqKly8eFH3/Lq6OlvP8KuvvkJ+fr7lkbi6PKGaZINdl5ubq3tcb1Dkhlm4qakJgPEzCvau6N1LPYEYCi0tLcr1ZpqfUblE3UJBb6sfK3XRe27q665eveqY8HHKmcaIVqX5MXM9M/87M88CkAZgLrxLHdKY+T1mPm0xqSIA6vDsgwAEhHYgovEAPgCwmJkVQzIz3/b9vwtgI7xm1JjBqYdqxbRkltfhw4cNf1uzZo2p8wMQfmcoTHPV1dWmwXxDvV9WnBPUdbC7t9hnn32mG0igoqJCd+3m/fv38f3334c0T6p3r+/duxe0DQR7RkbOE3pYNXsC1p/ZJ598Ynp+MOGnd1zUuaamxtYyGe29amho8HNAsRMKb9euXaZm5WDoCbE7d+6ElJa6fPv27Qt7JwWR3okTJ8JKJxZxIrA1mLmJmYuZuTKEy08CGEFEQ4moA4DXAPgtVPKZVr8F8BYzX1EdTyKiZPEZwAIA50OshiuITj7ac35mI1OPx+O6uU/rlOHE3IRTWH02evfo8OHDSt301l4FG1QIggmub7/9Ftu2bbOUlhO46RCUk5OjezyY8DNzZlq7dq3hVlXr1q3DrVu3/I5pzZ5Hjx5VrrfrMVxZWWlZyOila+YFaxWRbij9jJW6JiQ4EQY6tsygjgi/cGBmD4BfAfgOwCUAXzPzBSL6JRH90nfaXwHoBeDfNUsa+gI4RERnAJwAsJ2ZdyGGiKWHHQ6RqocT+YTTcYeT/8mTJ0O+Vo1R+WtqagI8EIXTQ1NTk+V6nz8ffHzohuYnMAqUrE7HyWf44MGDACcVrfCrr69XvF/NhJ+45siRI5buoxVCnafTu+9CixS/2bmPZoMPpwZDsdQfRl34AQAz72Dmkcw8nJn/h+/Y75n5977Pq5g5xbecQVnS4PMQneD7GyuujTbqBxwrml+4HDx4MORr3XChNkNoFllZWbh27RqamppMNVvhUBGMUMocTj31tEYj8+5HH31kuRPVxvt0o7MrKyvzM9NaEWxqs2MwrHTuHo9HsTBoBa72OrVLvpX8z58/79jSHyfDmxUVFQWYcAW5ubm2tieqqalRPocaMk/btrQe49nZ2WhsbGxdc34CIvorvT8nCtda0XN4cZqbN2+a5mvE9evXLZnPjh8/rkTVZ+YAM6WdetkRfi0tLUG9E4OhFnRlZWX48ccfTQMWh7Lu6Pbt25ZMXdp6FhcX67qoi864ublZWRhs1zRsZZmH1c7Pjuanx4YNG/zCrNnp4M06Wu39NPOe/Oyzz/wWsZuhnne7evWq4XIGK23Z7hRCqPN7RnzyySe6Wlx+fr6f+Xft2rW4fv26ocl07dq1SiSkHj16OFpGQXZ2tt/uHVu3brUcoShcnND81Ls6NMPrdZnuQLqtisrKSl2HBPEiqtciheMqLTBboGzGxYsX/Rqb0ULeM2fOICcnBzdu3MDq1att5VFYWIgvvvhCCVSt14kbCb+mpqawnAf08hEOKVYEttmaMfX127ZtQ05Oju3BzdatW3XNo0KQ2hESoQystHNfRlgRfiIto3LU1dVh586dyMnJMayXWrsQ6VRUVAQV/GZLZoRXa1NTk2F9tWbPUKitrQ1wnCorKzMdaAlBo56D15urtFI2s3PE/fZ4PIZe4jU1Nbh9+3ZAWDm9dPv27Ru0PGoOHTpkuh+jGlYtMykuLsb58+dRVFRkK79QCHsWk727OigQ0T9A47DSHvj6668RHx+PlStX+h13K2K73gS0lRdG2xkE85I0GoW3tLQYRgEpLi5GdXU1PvzwQ3Ts2BEdO3YMqaxOYDefgoICVFVVITk52TAtI4eNUBCduBjtRmORvp6mLTrp8vJydO3aFR06dAg4x6is6s68sLAQhYWFGDt2rO65x44dUzYwtfOunD9/HhMmTNDV/Kw4GVkVftnZ2ZgyZYoyD6hN+8iRI34WAKtBxdetW6d8Tk9P9/uNmQPCFpaXl6NXr16G6Wnvg7heCFkra3ydfCdzc3PR3NyMwYMHBz9Zg7B8uL1nqBtzfl3QTrY0am5u1g3PpF4T49SkuBarXmNO56HHxYsXDc1EDQ0Nup2kW4MCM5PT1q1bLaURbLG7GEBYuT/BnsnmzZv9vuvF5jRa4mBXIzdCT8sSdfvmm28MF7Eb1S0/Pz/gmPp5G0VXUZ+jjnVqxNq1a3WPNzU16Xb2VVVVtndwyM7ORl1dHVavXo2SkpKAdK9du+Z3/6yu91O3r5aWFng8HkWInj171s9cW1xcbClmqxpxT0XZ6uvrwwoPR0TweDwoKChQyuRUwAqPx6Nr+nV7gOzEnN85VdDpCwByAfxz+EWLbZqamvDhhx/6jfSIyPb6sVBx2onEqNxGQkp9/PTp07ZNudqYf6EGEW5qasJHH32kfNfOw4X6PMxMV+HOd+oJTbNYmlZNlU6iNnsadehGZjI9k6mZ8AvF5CvQG1hdvHgRmzZtCjheVVWlaNhazS8tLc0wD/GOP3jwIKx37MCBA7qbE7e0tGDjxo2KYNaa3oWQ0eZttu5T3Esh8BoaGgyXggjU6esNHvLy8rB7925cv34dW7du1fUd+O677/D555+b5qPl/PnzunVxezmUE5rfcwCe9/0tADCAmX/nQLoxjTB/qEOlNTc32x6hhYrTwk/PlEdEthaj2wnIrHVEMAsiXFRUhCNHjvh9Lykpwa5du1BfX2862nYyMgUzo6amRulEIhUcW+9eW10/aDd9ZkZBQYHupqtajEb+wUJ4acsuNuy10n61nbLR3JDe3K16rrW2ttYvP70ABnp5hhMi0GiLqJaWFlRUVBgOJkS5tVYks4Gd2DtPlNeOB21eXl5ANKi4uDjld7NA6jdu3EBtba0lwSXKaFS2cKLmWCFs4cfMN1R/t3zr9to8Rg8slMXiVjZ51L4Q4ZoNy8rK/CbC7UZc0XNdd1rwFxcX49y5c7hw4YLfi79jxw5s3rwZN2/ejOjuDC0tLbomPTPEfSorK1MEPBGhoqLCMNxYMDweDz7++OOQrtVDfQ/v3r2L3bt3Gwq/7du3o7i4GM3NzYrQ0nZ0wTQ/rdAU2puVd8euJ2BSUlLQfIMhtPL6+npdL2sjCgsLlXfMLDSh9t0zeueOHj2Kr776ylLewpwu7pd4VmaIQYneEg69MlkRcEbLiA4ePBg01KHbsYgdmfMjohQimk5Es8WfE+nGMk7NWd25cwdffvklAOCLL74I2kjNPLLsaH7akbGoT01NjdIJ2dX8QsXoJTh16hSOHj3qF7dT7akqyqjGzd2iT506ZXshtjj/xx9/VEy7RIRTp075hRtjZsO4olrUG/E6gVoTE+VVCzB1G7h16xa2bt3qNxjRztfoCT8zzV6sO7TSQV+5ciXoOWrU3qRa7JrE7a51s9oWtQMqs3dLax4M5nwlPC610XH0Ng4W5l29QYg69qkTXLp0KWg8XqeiyhjhxJzfKgAH4I3Q8t99/38bbrqxjhObyNbV1SkNLisrS9lnzUywihBaTps9xbVr165V5rvMhJ8TC1NF52OkTYkRtzBtVVRUBMwzaMsQ7hpBLXqed2YYmVm1glLPO0+7O7cRGzdutHReOAQLGG026tcTfk7t4WZ1+yQr2J1nDiYs1YJZHWfUbnxNO+/VyZMnbc+NlZeX+81pd+rUyfT8fv36+S1HcAoh3Iw0/s6dOzuanxYnNL8/BjANwA1mfhrAJADu7hkUA6g7xFDnlfRMKPfu3QsYpTU0NCgNfPPmzSguLnbcW7KlpUV3pGq07m3dunVKvUPVtg4cOGBoDrp7966ltT7C+8wtDh48iA8++ED5ru4AcnJyApwItM4pepFIzp8/H3Ob6WojEakdmJgZn3/+uakWpSYa0TrcRNQn2Oa2ai2soKBAuad23o9QFrzrOdGYoRbGFy9eDGpKTklJsST8rFouBMIMHa1txJzQK+uZuZ6IQEQdmfkyEY1yIN2YRj3astopaNHbnBPwRlcZN24c9u7dq3Tu3bt3V343ctkPp9Npbm5WzK+CkpIS0xejvLw8rJG42X3T0xT06qd2hHGa5ubmACeFYPMQRkJNezzcCCpO09zcjO3bt+uammpqalBbW+s332Q0+Lp3757rjgqRxqrmpn0XtB6XVjh//rzt9/j69etITU21fL7VJT+CuLg4Q7NnS0uL0pb1djeJZZwQfkVE1APAJgB7iKgCQOR9syOMukFr56GsYjZ309TU5KfVWBEy6g6ppqZGMVuqJ/0F2lGa3ugr2Ijw7NmzAQLc4/HA6mbBZtqrXnki5Ukr0NsGyoqXpbqTuHLliu48lR3NzymToRkej8fQFV60A7XgFwuRtXz77bfOFy7KhOpZK9ZH2rHShOref+rUqeAnhUhcXBzy8/N1pxQaGhoU86TZcp1YxAnhd4W9Wxn9loj2AegO4EkH0o1p1ILBjXVYYu8zO+zbtw9LliyxtPhZ6/nmlOnh/Pnzlu+H2YsejUgnWvTMVdpyJSYmBmg6VrzU7NxvtwICWM0jlA65V69etpa+tGXs3L9evXqFtHjcTVOzEGp65crPz8cjjzwCQH8aJ5ZxwvYyX3xg5v3MvEV9rC3y448/mq51iRalpaWOB8m1i53JfbddmcNFT4hrXbe1gq+ystJv0b0RkRBodjDroEPZEFWE4ho+fHjIZbLLhAkTIpaXHXbu3Gn5XGYOuhjdDG2oNCcw0+jcnHZwm5CFHxH9ARGdA/CwKsLLWSK6DkB/MquNEMu2baN5xEhhR6BFKhqOEf3793c8TadCPkUap6JpiBigZWVleO655yLq/DJjxgysWrUq4HhKSoqj+YwePdrW+XbugVkQcCuIOKnBjgVD7QGqNw88ZMgQ5XNDQ4OjWr7eNI0bhKP5fQFvVJfN+DnCy/MApjDzm3YSIqKFRJRLRHlE9Gud34mI/sX3+1kimmz1WjeI5Qn9YB5pbqPdKy6WCWeUnJmZiZ49e4as2UQjZJkW9YjeKH6nXcT8WFlZGQYMGGAaOcVJRBBkPUciO+sirTgiTZw40XJ6wRg/fjz69euHxx57LOQ01M+xX79+ShDx3r17AwAeffRR/OIXv8Bbb71lmk5aWhq6dOkSkKbe3P/ChQuVz0ePHnV0Pn7YsMiEhg5Z+DHzfWYuYOblmigvxnvC6EBE8QD+Dd6tkMYAWE5EYzSnLQIwwveXCeA/bFzrGsE6zxkzZlhOa/r06WGWJvp069bN9PdQ9gRbsGABJk2aFGKJzLHjIafHyy+/jLlz5zpUGnPmzJkTcGz69Ol+WsiiRYvw6quvYsaMGejatSuAnx1revbsGXB9586dld+tbuirZeTIkbrHRSeqJj4+3rADfvzxx5XPwe6pWusAgCVLluie9+KLL5qmI1CvJ9PbvUKL2PFDaCjhtM+ZM2fihRdewCOPPILZs3+ODTJo0CDLabz99tvKZyJS7v0LL7yAd955B4C3Xkbr5qZOnYopU6bgySefVKwxAwcOVH7v16+f3/na99xq4IFFixbhtddew9KlS/3KqyZS7xMQGzu5TweQx95d2RsBrAOwWHPOYgCfspdjAHoQUX+L17rG0KFDsWDBAgA/dwKjRv28ymP8+PEB23KotyUZMWIEAOCdd97B6NGjMWPGDLz00kt+5z/33HPK59mzZyMtLU23gQwdOjTg2OTJk/HMM88ACC6o33vvPbz77rt49913Dc/JzMxEZmamrlnilVde0XVSGT16NMaNG4eePXti2bJlynG9QMLz53unipcuXaoIprS0tKBmm9mzZ+O1114zPUcg7lNmZqZls6dYZmL0YuptvaLeymXFihVYsWKFpbzEcxKmuqeffhqAv+lJ/DZ69GjMmjULCxcuxKBBgzB48GB0794dEyZMUNrRe++9h3feeQcvv/yyoqUmJSUhIyMDL774It544w3DsvziF79Q2qgaoaU8//zzhhrLq6++CgDo06cPAG9HvHLlyoAOWNzTMWPGoEOHDkhISMDw4cOxfPlyZGRkYObMmVi5ciWWLVuGmTNnYtasWcjIyFA64AULFih5CFauXInXX3/db3Dz3nvv+Z0jOt2uXbv6vRsLFy7ESy+9hOXLlwPwCgY1asE7f/58ZGZmYsyYn8fb6nnHefPm+b2/APDmm28iMzMTzzzzDGbOnOn3m7rvmDt3LjIzMxVBIbYFE881PT0d/fr1w/jx45GYmIgZM2YoeT/88MN48sknER8fj8TERGjp3r07XnnlFb98p0yZgk6dOmHcuHEYMmQInnrqKSU9MZBKTU3F3LlzlTK8/vrrAWkLli1bhoyMDDz66KMAvGbUwYMHo1u3bopGCsBvMPTII49g+PDhmDBhQmSEoFi8GK0/AC8D+ED1/S0Av9Ocsw3A46rv3wOYauVa1W+ZAE4BODVkyBAOh/Lycn7w4IHyvb6+nltaWvju3bvc3NzMOTk5XF9fr/z+/vvv87fffsu1tbXKsZaWFsP0GxoauLS0lKurq7mlpYULCgrY4/H4nVNVVcXNzc1+xz0eDx85coQrKir44sWLSh7i/71797iiooIrKyv5zp077PF4lDzUNDc38+3bt7msrExJv7m5OaCcLS0t3NTUxNeuXWNm5vv37/Pu3bvZ4/Hwzp07+c6dOwHXvP/++/z+++9zc3MzNzY28s2bN/n999/nixcvGt4Pkdfdu3eVe8nMXFtby5cvX/Yr99WrV/2uq6io4Pv373N+fj6fOXMmIN28vDzeuXMn19XVcW1trVKnoqIirq+v54sXL7LH4+Hc3FxuaWnh8+fP8/vvvx+QTnNzs+4zVR+rrq5W6n/79m2uqanh6upq3rNnD7///vt84cIFfvDgAf/www/c1NTE9fX13NzczOfOnePm5ma+cuWKabvR5nvs2LGA4+Xl5QHPsqmpiS9dusSHDh3i/Px8rq+vV9p3fX09v//++3z+/Hluamryq4soi2gjemWrrq7msrIyv2MnTpzgK1eu8O3bt/2ue/Dggd97ZaWOwbhy5Qr/9NNPyvfy8nK+e/cuNzY2cm5uLpeXl3NZWRmfPHky4NqSkhJubGxU8hJ/ol5GlJeX85EjR5TvtbW13NTUZKm86nuspqGhIaAPCIXS0lIlnQcPHnBhYaHp+ffu3TP9vba2ls+ePcv37t3jQ4cOcVZWll/7amlp4YaGhoDrmpqalHvY2NjIpaWldqtiCQCn2ED2EEc5GgMRLQOQwcyrfN/fAjCdmf9Idc52AP8fMx/yff8ewP8F776BptfqMXXqVHZzXYxEIpFIog8RZTPzVL3f3I0cao0iAOrtfgcB0K4aNzqng4VrJRKJRCLxIxaE30kAI4hoKLyRYV4DoDUmbwHwKyJaB2AGgPvMXExEpRauDSA7O/seEYXrFpkK4F6YacQabbFOQNusV1usE9A269UW6wS0jnoZ7lIcdeHHzB4i+hW8u0HEA1jDzBeI6Je+338PYAeAZwDkAagF8I7ZtRby7B3snGAQ0Skjdbq10hbrBLTNerXFOgFts15tsU5A669X1IUfADDzDngFnPrY71WfGcAfWr1WIpFIJBIzYmGpg0QikUgkEUUKv9DR34iuddMW6wS0zXq1xToBbbNebbFOQCuvV9SXOkgkEolEEmmk5ieRSCSSdocUfhKJRCJpd0jhJ5FIJJJ2hxR+EolEIml3SOEnkUgkknaHFH4SiUQiaXdI4SeRSCSSdocUfhKJRCJpd0jhJ5FIJJJ2hxR+EolEIml3SOEnkUgkknaHFH4SiUQiaXdI4SeRSCSSdocUfhKJRCJpd0jhJ5FIJJJ2R0K0CxANUlNTOT09PdrFkEgkEomLZGdn32Pm3nq/tUvhl56ejlOnTkW7GBKJRCJxESK6YfSbNHtKJK2EhoYGNDQ0RLsYEkmbQAo/iaSVsHHjRnzzzTfRLkabpKmpCfv27Yt2MSQRRAo/iaSVUF1djZqammgXI6a5fPlySNpxeXk5rl696kKJJLGKFH4S25SWlka7CO0SZo52EWKeAwcO4MqVK9EuhqQVIIVfG4CZUVhYGNb1VjvWuro6bNy4MeS8JO2b7du3Iy8vL9rFCICIol0ESYRpNcKPiAYT0T4iukREF4joj33Hf0tEt4jotO/vmWiXNdI0NTVh586dIV+/evVqXL582dK5UvuQhMOtW7dw7dq1iOaZk5NjO8+amhrcuXPHpRJJYgFHhR8RJRFRvJNpqvAA+FNmHg1gJoA/JKIxvt/+iZkn+v52uJR/m6aiosLSeVL4uUtLSwuysrKiXYw2xcmTJ7F3717Tc7Sa34EDB7BlyxY3iyWJMmEJPyKKI6LXiWg7Ed0FcBlAsU8z+3siGuFMMQFmLmbmn3yfqwBcAjDQqfRbK0VFRdEuQqvn+vXrMTNP1NzcHO0ixAzMjKampqjknZDQLpdAu0pzczPKysp0fzt37lzEl/GEq/ntAzAcwG8A9GPmwczcB8ATAI4B+DsiejPMPAIgonQAkwAc9x36FRGdJaI1RJRicE0mEZ0iolNtyWFjx44dqKuri1h+Fy9ejFhekeLAgQP48ccfo10MiYacnBx89NFHYV0f6mBCCj/nuXTpEjZs2KD729GjR1FQUBDR8oQr/OYx898y81lmbhEHmbmcmTcw80sAvgozDz+IqCuADQD+hJkfAPgPeAXwRADFAP633nXMnMXMU5l5au/eutFuJBbIycmJdhEkJty9e7fNmE0fPHgQ1vUnT5401DTMuHfvHqqrq8PKWxKIx+OJdhH8CEv4MXNQm4SVc6xCRInwCr61zPytL/0SZm72Cd/VAKY7lZ9EEi1CnVs1Ehh1dXVyvtYE9Zzft99+i+Li4iiWJjKcOHEipMFBqMRa+wtZ+BHRfCJaTUQTfd8zHSuVfn4E4EMAl5j5H1XH+6tOWwrgvJvlkEgiweHDhx1N77PPPsPNmzcdTTMWqK6uVsxl3333HQCvdeKrr7wGJ6tezO1xqcPp06ct359IEGnhGI7m958A/DmAN4loDrxmRzeZBeAtAHM0yxr+FxGdI6KzAJ4G8F9cLkfM4USjaY8vvyAW6643txruc66vrw/r+khjpb4nT57E7t27AQA3bnhjGNfX1+P+/fsArAu/9kqsaWORJJxZ3VJmrgTwZ0T0dwCmOVMkfZj5EAC9XqrdLm2ora21fG5dXR06d+4cUj7Nzc2Ij3drBYtxnvfv30fPnj0jmm+0cbMzamlpCX5SK8Xtec7Tp09j9OjR6Nixo6v5RJr2LPzC0fy2iw/M/GsAn4ZfHIkdPv/8c7/vRg35/v37+Oyzz/yONTU1WWr4BQUF+PDDD0MvZIhcvHix3QRxzs3NVRyJwu2MzNZraoUfM2P//v2tIqYlM/st68nKyopox33ixIk2uawokvfQKC/18fLycuTm5kakPCELP2beDABElOr7/q9OFUriLHqN7qOPPsK5c+eCXiscKCI9Qow1z7BQaWlpCepuf/z4cZw8edKR/My8cbXPsLGxEbm5uZZ2M9i/f39Ut1MqKyvDjh3uGXmsmL7bopYUC3VSv+snT57E/v37I5KvExFe1jiQhiQMgjVgYbLUnldVVWU57Xv37vld7/RLU1tbGzUXfafn/BoaGpCfnw8A2LNnD9atW+f3u/beqYWjm51ROGnn5ubaCmgullzcunXLdl6HDh3Czp07/eYo9dpepDvuWBAUThMLml+0cEL4xZ63gESXcBuf2mymNaHdvn07rOgkrc0Zw4zLly/j+++/B+DVWLTbEK1evdpPizK6b2aC4/Tp04YLho2I5JyfWCqgt/QiWDvMz89HYWGh31o7cY36Xh0/fhz37t1zoriGAyC770xOTk6rcrKJJYGUnZ0d0fI4Ifxi5+7FEDU1NaisrIxIXsFGwuGMlI20vTNnzvidt23btrBChLnhcRktbVJtxjGql/ocIZTKysr87vHp06cN8yguLra8RksIDO0Aw00v13AErRgY6M1fHjhwQIn6cu7cOcsxaYNh9I6o62HV+9QpE3YkiCXhF+m9KqXm5xLbt2/H119/HZG8gjVgs7maYBqXOm11R/DTTz8FeJuqz9V25Ebs2rULZWVlAR2xWQeSlZVlKaSbWnM4e/asa1H6hVu9wEq99QSPHU3OTvgtIWi1AxYzKisrw2q/oXaqQmNMTEz0Oy604KtXr7o6H6xNe82aNcjOzgZgPrhcvXq1a2Vyk1gQfm5Op5jhhPD7jQNptDmiEaDYqOEcPXoUQOBovLq6GufPm8cEUAshrSDUepuq2bBhg+miajHKu3nzJgoLC21rIXbNpMeOHQsami0ULbGhoUFZUA34B2Nm5oB51WDzrOp7bPWeNDc3m2paIp2uXbtaSg8A7ty5EzHLhRqRp7Ytu6lNtbS04MKFCwCgmKsFzIy7d+/6lYmZ/bTulpYWW3tihkJ+fr5f4IOsrCxHgn7HgvCLFmFHb2Xm80Q0FcBfAkjzpUnen3h8uOnHOjU1NUhKSopqGYI1YCGItefduHFDGela6WjtmrLKy8vh8XgwfPhwv+P37t3Dt99+i8zMTKVcIv+qqiokJycHTZuZUVNTg44dOypaUFZWFlJSUrBs2TLlHCuEo0lo78nhw4d1F6h7PB6sWfOzb5iV+60+x8jT8d69e9i2bRvS0tIMlyyI+9CvXz/d43qo6yW0bHXnH6z84tycnByMHj3a9Fw1cXFxAfnbWc8aCvfv38elS5cAQHfA1tjYCODnOl2/fh179+7F2LFjMWzYMGzdutWvzOK8HTt2YMaMGejVq1fYZTx79ixKS0sxaNAg5ZgTc7jaNpCVlYUXXnghoK1YYe/evbh27Rqee+45DBgwwDAvj8fjZ7mIlgB2KnT5WnijvZwD0HZX0mqoqqrCl19+qXTi0cJq47l+/TpGjRrld8zO6NFuIxWjda3w0xM2ojPduHEjJk6cGDTtvLw8nD59GmPHjsWFCxcwc+ZMAP7zRCIfUW6jDttKJ1JWVobk5GR06NDB9Fq1R6R6x3ur7ttG99hojVldXR0aGxtN5/9EmnY6S3U5hEC4d+8eEhMTsWXLFrz++uummqS4PtQA0eqymlkYnCBYuy4pKVE+ezwe5Z25cOGCIjSBwGdUVFSE/v37hyX86urqsG7dOnTr1g3AzyHcAGfmbPXqXlpaGpLwExsGFxQU6Ao/8UzXr1+PqqoqpQ2pfQUKCwtt5xsqTm1mW8rMW5j5OjPfEH8OpR2zmGkMseRCrF6uEIwdO3YYrv9zqk5akzAz48CBAwC85sxjx44FTUM4g4j5TL1rtm/3xmEQpsZwOosNGzYo5mM1hw4dAvDzvRGaC+B/v40CTmvvqfoaK/dbT1MyysPsHO0z1zv35MmTihZ4+PBhU4civbLbqY9bCBOlGqua//Hjx7F9+3a/gYz6PgnhV19fr2irzIz6+nrlPpWUlODIkSOWy1tTU4OmpibXPHVDfU5mGE2l3L59G8DP72NNTQ3u3r2r+15FAqc0v78mog8AfA9A8a5g384LbkNECwH8M4B4AB8w899FIl+7MDM8Hk/AZL4eVVVVYGZlxBcsXfV/Iy5cuICePXuazuUUFRWhqakJ48aNC/jN6AW082IeP35ccbxQC0Ejt/7m5mbExcWhtLQU2dnZygsksKK5BtP8rM7P6nWSYu5SREtRawl6ZdB+P3jwoN9x9bykOMds7lQIC7M6iHQ6depkWKZTp075PXMjL0dhAgR+nndtbGwM0Ij1EPNqZrgdfPujjz7CuHHjMH26d/OXqqoqS4v8AW89y8vLDX/v0qWL8lloqy0tLX5ejJcuXcKVK1dQVVWFjIwM0/yuX7+uLJvQe8dCFVI5OTmKpUAvjWPHjmH8ePszVsOGDVO0Pz3E3KkgLi7OcEB+8uRJTJvmasRMx4TfOwAeBpCIn82eDMB14UdE8QD+DcB8AEUAThLRFmZ2fddVs8YnzD3FxcXo27cv4uLicPbsWRw/fhzTpk3DpEmTDK89cuSIMnoSJtX79+/j8OHDeOaZZwLOz8vLs1xObWcrUI+4mRl1dXUBLvt69W1pafEz/ejR0tKCpqYmJCQk+HkcfvnllwDMd6PPz89Ht27dsGXLFt3fRTBjLSdOnAg4JoTfrVu3UFVVhUGDBvl1WHpUVVUpWmmwOTKzpR56wk8vqr7abCoEmp6TjPBcFc9NLZSM8jbzEG1qasK9e/eQmpqq1Eev7ELruXHjhtI5ffzxx0o7PXr0KKZNm6Z7r4TGY+bSLoIDuEVzczNOnz6NsrIylJWV2Z5PNNMS9Z6B1klMPFOjdqvm+PHjisXAaIBZW1uLTp06GWrMFRUVWL9+vd/UTF5enjI94KSFyspAXY16WkBLTk5OqxF+E5g5UFWIDNMB5DHzNQAgonUAFgOI2pbjatf3rVu3Yv78+Rg6dKjS4G/cuOEn/CoqKnDz5k1MmDABgL/ZQASVvn79uqGQECPq/Px89O3bFz169LBcVjEa02pF2ligAPy8GgU1NTXKC33o0CGMGTMmYGnFsWPHcP78eTzyyCN+x0XHY7YEwePx2I71JybeBVoz06FDh5RntHDhQtO0iouLFa20vLwcFRUVSElJUX4X9y2YN6JeJ6snoNWd0Z07d3D9+nXd7Y2E+VF0jtp7npWVheeeew79+/cPcMQQaDvyb7/1jlWXL19uSZtXp9fU1ITExEScO3cOI0eO9PutpqbGr3MOd5NaJwh1bslMWOgtJVGf/8EHH9jKS23VMBK6n3/+OYYPH465c+eioqICPXr08HuXxSBcOJnU19f7zYs7Ocemzjc/Pz9grj/WcMrAfoyIxjiUll0GAlA/wSLfMdcxehG0Jqg9e/bg9OnTSievVf+FRqiHEJh6HaWW/fv34+uvv8bmzZuRlZWFyspKVFRUWNoZIdhci1EgbG1HVlJSgk8++cTvmBDmoUS+ICLbwk9rehEv+PXr1wPO3bVrl9/3jz76CFevXkVeXh6uXr3qp3VVVlZi/fr1uvNcweKkatclWh1x79mzJ+CYuh5mJrtt27Zh8+bNhmZxrVu/oLKy0k/4aU3Nepw9e1YZUGjb0g8//BCVpT9afvrpJ9fSNrIgmA0i6urqTKP4qO+j3rpW8Tzz8/Nx//59rF+/3q9t1NfXY+fOnQCAH3/8Efn5+fj00/D2H2Bmw6Di6n5NtK3q6uqohS0MhlOa3+MAVhDRdXjn/CK51EFvIifgyfg2280EgCFDhoSVoYjdKJwd1Bh5t2njIh46dAgzZsxAYmKiaUcYStQDMe9kZ5FyUVGRsphX70UT7txatm/frnhaAsDmzZsN8whlSYHbDhBampqaUFhYGNSUfPToUfTp0ydkN3wzk08wTp06Zfncu3fvKu1LmOXy8/MNBR8A7Ny5U7FC3L9/39AUqV5rmZ2drbSfHTt2+I36tQIgWvsnGt23pKQk5T1LSEhA3759/YTSSy+9FDQAQZ8+fZRNddWo131qOXHiBHJzczFz5kykpqaiT58+WLNmDZYtW4aUlJSgbV/9PgmrjPqYeuBWVVXliGOJeOZVVVXo0qULbt68iTNnzmDBggUBgpyZY0LLN8Ip4WduO3KXIgCDVd8HAQgYqjJzFoAsAJg6dWpYhu7Lly/7aWpCWHTq1AlffPEFhg0bFnCNVuu4ePEibt++jXnz5ilzRdeuXQu49ujRo34dSUVFBRISEiythbOD2lHDSsBrNW6OqN2I8B5sdwIrpiArO2KYYRahpnfv3qZBpO1qUUL45eXlobKy0pLXrxjFm2mWwptWS01Njd/8V21tre56soULF/oNRKO136C6bB6PR7m/Yp7bioOakSNMS0uLYXAFYdHIy8vDsWPH8NRTTwHwmukXLFgQNE9htVAL71OnTuHSpUuYMmWK3/6dRGQoTO/evYs+ffr4HTtx4gTOnDmDVatW+Q1WhHBdt24dMjIysHfvXgDA2rVrA9Jdu3atMjisq6vDoEGDkJiYqGuBiQZhCT8i+hMAhwHkMHO09qA5CWAEEQ0FcAvAawBedzND7Ut6+fJl1NTUKIubzTye1FRWVvrtWbd3716/RayCH3/8Ufm8fv16AMB7771nt9iuYeZsEYsEiw4Tza17gOCakd3RtLpztxoIWoQZ05roraI2cVdVVfmVWdzfXbt2+TliGDljuY1Wex87diz69OmDy5cvo7GxMahTFGD8TKx4uIpnIvqPW7duKfFLzRACV63tVVdXo7q6OiAogtlzvHnzJoqLi/3WI4qlRKtXr8aLL76I1NRU5Obm+nn1BhusqO9rQUEBmNmWJcdqwItQCVfzGwTvEoOHiegsgCPwCsOjzGzsE+wgzOwhol8B+A7epQ5rmDl4iwsD7UNXC75w0XNq0ZsXsGP6AoAOHTq0CiHVs2dPU3fycHFqFwA3cdosKIRfMI3STcw8YT0eT0iOTW4xfPhwDB8+HB6PBxcvXkRCQgJ69eplOZC4YMKECbbiqYY60LAb5kytKQLBLTfffvstJk6cGBBo3Y6naElJCTwej60+yO210mFNqDDznzHzYwD6AfivAMoBvAvgPBFFzNuSmXcw80hmHs7M/8Pt/MTchsApwWeEXiMIFqdSTUpKiunSCsETTzxhq1xu0LFjR1fTdzte5axZs8JOw07QaisI60I0NVqjjv3LL7/EmjVrwnbEcINZs2bZtrBMmjRJWU8ZqblNu+biuLg4vPTSS7au0dthxM6azCtXrqCkpEQx9wdb4wi4f/+c8iboDKAbgO6+v9sA9N0XJZbRhiKzgl7nu2zZMkuLkAcPHhz0HDVW0uzevbutNM14+umnw04jWCBvPYYOHWr5XO3ciR1SUlLwxBNP4LHHHsOMGTNCTseIaDofaHe+ENidXzYiOTkZr78e/myHur0SkdIBB/OYHjNmjHJe7969AdjXXMSm01qmT5+uG3QiVOLi4hyJNyq0+ccee8zyOyLqGBcXh6efftpvLZ92aYTbzm5hpU5EWUR0GMBXAB6F1+y5jJmnMvM7ThSwPTN79my/+IlTpkwJeo1Y4zdixAgAsDRxLjB6+YyYP39+0HPsNuCHH34YAPDss8/6racDgL59+9pKSw+7pqVnnnkGc+bM0f1t7NixAceISIkeYoWpU6cqac2dOxejR49GSkpKSBE22jNz585VHFPUgeZXrVrlSPpPPvkkHnvsMcPfH3/88YBjdpb29OrVy9CRqXv37o6umXNaqIwePdpSX6CmV69eGDFihJ9FauLEibY10nAI9y4MAdARwB14nU2KAFSGmWbMY2XdnB30OlEgMKrKlClTsHLlStO0xPkiUofwprNiGjEyM8ybN8/vuzrobf/+/U3TNPNq1DNxjhgxApmZmRg4cGBAxIhIewP26dMHgwYN8hsUqLVxbbgwwcSJEy0HO588eTIyMzMxa9Ysv3alfhbqpSSCaO8kEkmsWEDUGvfSpUuVz6KjT09Pt5SXkYYaFxdnydIBhB9IXEtSUlJYFgUt4p68+uqrpucZ9Uta7AyaRbvWcyISGql4d5w2/wfkF87FzLwQwDQA/+A79KfwhhfbTUT/PdzCxSpaYRAO48aNM42gLkwoIqyZWUPLyMhQBMbDDz+MzMxMpaGbzaWlpKSYjri0rt7qMmhf2nfffdfvu5ln5YoVKwAYm0+16wJDFX5WvPX00BsMzJgxQzG/MnPQDgSwZh42Q2+kPnLkyLDSdIvZs2c7nmaopnPx3N9++21LFpDOnTsrg0Y99J5Dr169FO1dEIqjhtEOGStXrrQt+IKZNEUdg3lSClOuk+i9U8OHD8e4ceP8IlO99dZbrs//h63/spfzAHYA2Amvt+dwAH8cbtqxijZ82MSJE0PujOLj4/3W9mlHqPPmzcPKlSt1l0DopdWtWzdkZmYGCKxgZpNgL4w6NJm6jKmpqX4aUHx8vLJA2gozZswwdBLR3uekpCTdNVdvvfWWaR525zIFei9qYmKiMtHPzOjevXtQZyJtx2bXvKwuh/jstBbslHPBww8/HPLcrNE7ZKYBqDtwUQciwssvv4zly5cDMNbQtbz11lt48cUXDX9Xr5sDvMLhpZdewuTJk3XPV9/T9PR0v3auzcdIwNltKwBMvVLT0tIUx7Zg5k/ttMOTTz4ZNG/tNVr02tncuXPx6KOP+v2mvdduEO6c338monVEVAjgAIDnAOQCeBGAs7bBGGbw4MG6Nn8riEn15cuXY/bs2UpjF/N7cXFxll8As5GSWecmOtIOHTr4mZiEGY6IAubbevfujZSUFDz66KN48803/fJRO2uIMhk15gkTJijzk1oee+wxrFy5UplM79ixo65TQ+fOnfH8888b1i/Ujl3vuri4OCWKj7hv6kl7vZdfe8zqSF7cO3U5RFsIZm62i5PrRkN1UTcScurNcBcvXuw392b0bHv27BmSOc4MMfDq2rUrFi1aZPjOi/qr28XIkSP9BrDi3Ro8eDCmTJmitKWhQ4fatlTovfdabVSQkJBgac5PL6i0lTV3RIRFixbh1VdfRWZmJh566CG/32PJXB+u5pcO4BsA05l5GDO/xcz/zsxnmLndbGqbmJiovDx2nTJEQ0xOTsbDDz+smEDN7O1GHZ+ZycYM0RnHxcX5je769++P5557DgMHDgzo0JYuXYouXbqYRo5Qp/3aa68F7WC0LzERIT4+3rQzFYJFnGMmWMLRuNTHxHG9zlqbR1paWkAnOWLECEOBr0ZvKybRSYYboi9S6G1qaoS6nur5XvU9TU1NVc5LTk62FH3Fbt5GiHb+2muvGVoTiAizZs3CggULMHbsWMUMnJyc7JdHXFwcnn32WWRkZPgJvw4dOmDp0qXo2LFjQF9iZJ155plnAoRVKJ7iavSsGURk6tW5dOlSZGRkYPDgwYqpWntfx48fj7fffjussjlFuMLvT5n5G2YuNjqBohXIL4IkJiYqL0awXQK0aG/P2LFjkZmZaWqqMWrYodzq7t27G67vi4uLw4ABAwIcb+x4iz388MNYuHAhEhMTg5pWrWyJIur40ksv4ZlnnsGyZcsAGDsZjB8/XrnGblgw9f1Um43FcSs7zvfo0UN3ZG5FOxJ1UZfDiYW/Zib0Dh06BMxj6s1Jm43gRRmtzrUJ1O1Kz8lHnCPOmzx5MhYtWoRXXnkFgL/Z0wpOmsMFzIyUlBRlamDUqFFYsmQJevXqFXDdwIEDAzYjZmYkJSVhxYoVWLx4sd/5L7zwgm6enTp1CrCsWH1H9bRMo2kEdfmFSVlN7969g2qHcXFxls3QbhOu8NtHRH9ERH7DUCLqQERziOgTACvCzCMmGTjw540j1BqAXTfiUASWk27P3bp1MzQ36ZVt8eLFtuY3O3bsqGgpZqZJo/wA6AreXr166XbiWuEQzthLfe3cuXMV06A4bkWTNNLGrQgxIayF9tShQwdbwk/vuWZkZChzT3oCJiUlJcDBpEuXLgHLN4TAMaNTp062nH3U7464x1oNV615jxo1Cl26dLG8hZdW21ZfZ6WdBNsUGQh8/4lIsUZYmXowe75GWq7eXptW24m6Dc+dOzfgmDYf0aaspu92lJZwCFf4LQTQDOBLIrpNRBeJ6BqAqwCWA/gnZv44zDxikkWLFimfO3ToEFInu2jRooA97qwQHx8f9uhJaEx6COGqV6e+ffuGNAkP/Pzydu7cGS+//HLA70b3cPz48YoZJiEhQfESVaMdQTuBtrOxq1lkZmYq91I9IrbTcQgnJuDnztrqAEtvsJGQkGCp/GprQFxcnJ8T09ChQ22ZG62aaPUce/QsKcHun1H9zEzNdpy0jFi8eLHpkgorwi8UjDaZFnTv3l1xyjEqw6pVqxQzq1k5H3vsMbz00kshC7VI79JiRrhLHep9c3yzAKQBmAtgMjOnMfN7zHzaiULGIuIhjhs3TukI0tLSkJCQYHm0O3jw4JDnLMK1JgvTjHZC2ggn3Y67deumu1bSqE79+/f3m9PQK0u/fv3w4osvOrIQXmDktj9ixAjL68YEWjOR1c5OjLTffPNNZcAV7rNPTk7Gs88+a3qO2slEaFti/ZXdHbvN5mETEhIUM76e8HMSszStBJAI9sz69u1rmkcwc6n6vx2YOeC6pKQkPProowCAYcOGKQ4wRoNmtTlZjdr0SkTo2LEjevXqpdyLYNacWMYxMczMTcxczMyVTqUJAET090R0mYjOEtFGIurhO55ORHVEdNr393sn87WKusFkZGQEnRR2I1+raDvsBQsW6I6G9bTRwYMHmy4p0I6cgy0/0DJ58mRL8UeNICKkpqZi9uzZfi+42kxmFyOBP3LkSFtzWeEgytClSxd07NgRgwcPVtpXsEXIRh0pEWHgwIGW74tW09emG2ywZ/Zc09LSFCcrdZt2Q0PQpqmuhx2Hl1Axy0MMKIIJPz0tWk/4ERHGjRuHzMxMZeD4+uuvB5iv+/btq1gl9Cwb6sGk+nhSUhI6d+4c1OtYW642o/lFiD0AHvFtjHsFwG9Uv+Uz80Tf3y+jUTirIzUnYunp5Ws1kogd9MwfRGS69karyYlztfdnypQpuuuipk6dGrIDQiwgRtlWsbKO6dVXXw0YWS9atEgJt6YX/1MtiPTm/IK1V/X1vXv3xiOPPOI3N5iRkRHURDh8+HA/c6UVjQfwvydiOYEaYSUxE+pm+Wk7XiHUrVpfUlNTdc31TmA1nN3ChQv9oqTMmTMH3bp1s9QPde3aNaBNPP3004qDk7g/VgYCHTp0sD3AtZp2pHA3fowDMPNu1ddjANxpfSESrQndUHZFt7v1iR1nhYceesiSCdWKeclJwtH87DBu3Ligmlj//v1x584d9OnTByNHjjRcHC0IFtlEr14vvfQSvvzySwDBw/Bp2+7y5cv9nrk6TJggLS0t4Fjnzp39tqpJSEiwPM8nyvCLX/wCFRUVynEi8hsMzZ07V+m4jaKhBEMt/MQ82JkzZyxrI0QUVmhDs77C7nwyAL/1taH2Q+r3I1gZnFh43maEHxF9xsxvEdEfM/M/O1UoE96FN4i2YCgR5QB4AOC/MbPhbphElAkgE3B2jVSwubA+ffrg7t27jj/0bt262d6b7pFHHtHdG1CPjIyMkNdotVeCdaJaLS7UsGtm+dl5DtoOM9SNQ1NSUgxjYgZDeA4HcxpTezgPGTLENMatUTrq46mpqYrGF6m2G8xRbOnSpSE/AycG4Uaa39KlS5GYmBjyoAMAlixZgk2bNsVUPxGu2XMKEaUBeJeIUoiop/rPaiJEtJeIzuv8LVad85cAPADW+g4VAxjCzJMA/J8AviAiw5l4Zs7y7TYxVcTLDJfly5frmoDUD1gsiQhnPkuPUEZhaWlplhe/pqWlObaAOJpY6Uy0Dhx2TZixgDBVG3WCkYyWr4eeFtm/f38/TdKOw4ueIAl2jXr5hNpKEakOuWPHjqaxT3v37m3Ji3vJkiUB99Pu/KUeRtf17t3b8nISIyK9x6EVwjV7/h7ALgDDAGi3A2bf8aAws2mkaCJaAW/otLnse8rM3ACgwfc5m4jyAYwEYG+L8zCwGu4H8F8X6AShNqJYXnfjNO+88w4SEhJw9OhR0/OWLFmCb775BrW1tQCsx4KMNuo2MGrUKDz22GOGz1fMObvx/K2kaWXAqa5PuFqxGdrlE5HskJ3IS+9eqp/BU089FVK6dub8QiWWhF+4Sx3+hZlHA1jDzEM1f5YEXzCIaCGAvwDwAjPXqo73JqJ43+dhAEYAuOZEnuGi57Xm9EMPtYOOpPBLTEx0dDPbUPInIkyaNCnoztFvvPEGAK9HbHp6esgdSLQgIvTu3Tuk2IlWAhY7idFuAeIdWbFiRchrSc2wYg5trajfayuh8/QgIrzzzjtOFQlAYGCAWLrXTjm8/Bcieh3eWJ9Kmsz8Nw6k/Tt49wzc47txx3yenbMB/A0ReeBdaP9LZi53IL+wmTZtGkaMGIHNmze79rBnzZoVcecRuzj9IoVK586ddR01BOpJ/5SUFCQmJsbslkHBiI+PR3JysqUd0kXHFG4cSLtMnjwZFy9eDDjuVAdpV8jFUoccbdye6oile+2U8NsE4D6AbPhMkU7BzLouhMy8AcAGJ/NyChGUdv78+ejTpw9OnjzpeB6JiYlKQ33vvfcsa3TDhg1DQ4Ojj6jVkJKSongUZmZmIisrK8oliiz9+vXzW3ITrhXghRdewJYtW2x1aPPnzzc0aYYr/KzO+dm9rjUgnqVezM1YIpbutVPCb5BvY1uJiqFDh9peXhAKdtz509LSTLWgtszChQuVZQBq1O79sfRyWmHVqlX44IMP/I4Z1cEoMHKo9OvXD0uWLEFycjIKCgosXWM29xctzS+WFl6HysCBA3Hz5s2QvUXdJhbNnk499SNENM6htNoUeg/diTiCEvvozSONGDHC9bBabhLtjrtPnz62PI/NNqZ1u4Ps1q1bQFi3+fPnt7r5XT369u2LJUuWRLsYATgZaN5pwl3ndw5er84EAO+QN6h1AwCCd5N3a2EL2jCx9LDbO3rPQr3MISUlxdbaxlihZ8+eIcU0jbTn7y9+8QvTwAmRMHtqva4jEYpQ8vOzifZgTU24Zs/nHClFGyYS7sMSf7p37657v7XH3njjDT+txWyni1gmWMgtp3d9D5VgEYPcNnvGAu1pqZEaK9tBRZqwhB8z3wAAIloGYBczVxHRfwMwGcDfArgRfhFbN+qHvWjRInTr1g0PHjywPEcisc+zzz5rqZMJZVlAa8Sow4m1jjgW54UkzuDkVmNO4ZQO+n/7BN/jADIAfALvAniJisGDB6N79+4YPHhwwG7ZEudISEhoE9FpQkUrPCIlTPQCbZvhRKxIPaTwjD1iUfg55e3Z7Pv/LID/YObNRPRbh9Ju9bz55puuLNqV2EN0ilb3MGwLTJ8+3dE9Ds2wE/vxjTfesBU4va0Qa9q22wwYMAClpaUxWW+nhN8tInofwDwA/5OIOqJ1bJcUEdwM1SSxTmJiItLT05Utgdoq/fv3VwT9xIkTDc9zMoxbRkYGBg0aZPl8N0zOUuOLPcaOHYuxY8eCmYNuoBxpnBJ+rwBYCOAfmLmSiPoD+HOH0pZIHCEuLi5im9BGE7PgyWrGjBnjmLdjrKwdXblypRSCMYiep220cUT4+WJufqv6XgzvrgsSiSRGIaI2Z5WI9ekFs3WOksgin4REIokJkpKSggYgb+0MHz4cKSkp0S6GBHJeTiKRxAhEFDPmU7eIi4tDampqtIshgRR+EolEImmHUCy6oLoNEZUi/AX4qQDuOVCcWKIt1glom/Vqi3UC2ma92mKdgNZRrzRm1o2m3i6FnxMQ0SlmnhrtcjhJW6wT0Dbr1RbrBLTNerXFOgGtv17S7CmRSCSSdocUfhKJRCJpd0jhFzptcRvwtlgnoG3Wqy3WCWib9WqLdQJaeb3knJ9EIpFI2h1S85NIJBJJu0MKP4lEIpG0O6Twk0gkEkm7Qwo/iUQikbQ7pPCTSCQSSbtDCj+JRCKRtDuk8JNIJBJJu0MKP4lEIpG0O6Twk0gkEkm7Qwo/iUQikbQ7pPCTSCQSSbtDCj+JRCKRtDuk8JNIJBJJu0MKP4lEIpG0O6Twk0gkEkm7IyHaBYgGqampnJ6eHu1iSCQSicRFsrOz7zFzb73f2qXwS09Px6lTp6JdDIlEIpG4CBHdMPpNmj0lEkmboba2FswMALh06RKam5ujXCJJrCKFn0QiaTN8/vnnyMvLAzPj4MGDKC0tDSmdlpYWNDY2Kt+FQJW0HdqE8COihUSUS0R5RPTraJcn0rS0tOCjjz4K+fr6+nrLLzcz4/r16yHnJQkOM6OhoSHgeH19Perq6kJKc82aNbh9+3bIZcrKykJFRYXl85kZJ06cCDk/uxQUFCht+MaNG1i9ejUA4OjRo7bKLThx4gQ+/vhj5fvq1atRUlLiSFklsUGrF35EFA/g3wAsAjAGwHIiGhPdUkUWj8eDpqamkK//9NNPcebMGd3fWlpa/IRdWVkZ9uzZEzTNwsLCiJmcmpubY868pSe8jLh58yZyc3OV73v27MEnn3wScN4333yDzz77zHZZmpub4fF4Qh601NfXAwAOHz5s+ZqWlhacPn3adl7Nzc24deuW4e/qdl5aWor79++joaEBu3fvVvK7du2a3znr16/Hhx9+aCn/kpIS5OTkKPkwM1paWgB4TapG1NTUWEq/PdPQ0OD3bKJNqxd+AKYDyGPma8zcCGAdgMVRLlPMc+jQIT+TUHV1te55X331Ffbs2YOysjIwM7799lsAXk3ArJPauXOnbmcrtJpNmzahpaUFBw4cQEVFhTJqb2pqQkNDA9asWYMdO3YA8HakzIz79++jpKQEeXl5qKiowIkTJ3D79m18+umn2L59e8j3QpSruLg4pGu19+7GjRu6wkuQlZXlZ1LbtWsX9u/fj6amJlRVVaGgoAAAlHMaGhpQWFjo1/l6PB5Fk1MLhKysrAANRfyuJ5AbGhqC1vvGDa/PgBXN0WwQ0tLSElRzzcvLw/bt2/HgwYOA33JycvDRRx/h8uXLqKiowMaNG/HVV1/hq6++AgCcPHnSVrmYGcyM5uZmpW7btm3DyZMncenSJQBeje+DDz4AAMNB36VLl7B27VpbA562iLifRly8eBF79+71O1ZcXIysrCy3i6ZLW/D2HAigUPW9CMAM7UlElAkgEwCGDBniSMZVVVW4e/cuhg8f7ne8sbERH3/8McaOHYtZs2bZTregoABVVVV46KGH0LlzZ8PztmzZggEDBij1KSkpQd++fS3lcfHiRRARevf2egHrNVpmRlVVFQDv6F/dYQPA9u3b8dxzz2HAgAF+x3fu3AkA+OGHH3Dw4EEsXrwYiYmJqKysVH4DoHQqly9f1i1jUVERCgoKsHv3bsN6iNH+nTt3wMxYvXo1VqxYgY4dO+Lq1avIzc1FRkYGzp49i/LycsyfPz8gjZKSEnz33Xeor69HZmamYV4AsHnzZpSUlCjnNTQ04IsvvkBmZqZyD+/duwfAqyVv2LABkyZNwrRp01BXV4eOHTsC8AqCiooKvwGE1nR9+vRpjBs3LkDbO3fuHLKzs9HY2IgVK1bgk08+wZIlS5RnuXnzZr96mFkFjhw5gqtXrwbUu7S0FBs3bkRmZqai+QBAfn4+UlJS0L17d8THx+PGjRs4ePAg3nzzTbS0tODDDz/E66+/Do/HA8DbhogIAHDw4EHk5ubivffeU47plQcA1q1bh8WLF6Nnz55ITExEY2OjItxyc3Nx4MAB5RqhmVrl4MGDmDlzJtauXYsxY8agZ8+e+OGHH5CZmYmePXuazhOq66NOD/CaSpOTk9GjRw/s3r0bPXr0wCuvvKKcV1FRgeTkZFRXV6NHjx6Wynr16lXs27cPAwcOxK1bt4K2Tzt4PB48ePAAPXv2dCS9HTt24NatW1i1ahXi4gL1KjEAKSsrw969e/H0008rA6tdu3Zh4cKFjpTDKiEJPyJKAlDPzLFga9J7iwJ6cmbOApAFAFOnTnVk9vrLL78EAAwbNszvhRCawIULF9ClSxfcuXMHhYWFGDhwIJ599tmg6YrO/ujRo6aN/c6dO7hz547SgW7evBkAMGvWLIwdOzbgfNEpp6amBvwmOu6dO3di0KBBuHHjhp/gNdKstm3bhvHjx2PGjBk4deoU0tPTUVj481ikqakJ33zzjWl9zTATfFrEPM8nn3yCN954A/v27QMQKFTq6upw69YtnD9/HkuWLFHumxW0WpX6hf7hhx/85pc2bNgAwKuxjB49Gl988YXy26effho0r5aWFl0z582bN5WBiBBMmzZtwuOPP66c4/F4sGbNGj/hJZ6nGNAkJycbCkbRVs6ePYtjx44px7///nsAwMyZMzF+/Hh89913ALxC8eLFiwC8mqKeQBKm3RMnTmDGjIDxKQB/Qb1582ZMnDgRU6dO9Zt/E4I1FI4ePYpLly5h5MiRaGxsRHFxMZKSkgDAkgayZ88eLFiwABcuXEB5ebmiIQLw+wwAlZWVqK2tRZcuXdDQ0ID169crv73zzjtITEw0zWvHjh0oKioCAOUdv3XrFkpLS5GcnBzQ7wjNS0/w6LFmzRoAcEygijJ+8MEHeOuttwIG7uJdEe/Fpk2bMGLECADeNg14+5/U1FTk5OQ4Kuj1sCT8iCgOwGsA3gAwDUADgI5EVApgB4AsZr7qWinNKQIwWPV9EIDQZ/atZFhUpJjkAG+nu3TpUmXkrTZ/qE0xt27dQnl5OZqamixraHV1dabaHxDYIR8+fBgDBgxAp06dsHbtWrz44ovo2bOnYrJUNyphYhPCr7Cw0E94WeHs2bPIz89HTU0NcnJybF3rFmvXrtU9XlBQgOPHj+P+/fuO5CNGruKFNuLo0aO20z579mzQc9TP/tChQ8pntYAQwk/8F53wihUrlPNaWlrQ0tKCwsJCDB06VDE7qgWfmuzsbF2hCHjb35QpUwDoa0p27n11dbUyiNGrmx2am5tx7tw5AF5NA/Dev06dOpleN23aNOU9LioqAjNbnv/8/PPPdY+fOXMGU6dONb1WCD416kHooEGD8MUXX+Dxxx/HiBEjcOTIEVy4cAETJ07E9OnTLZXPLQ4cOICMjAyUlpYiPj4ePXv2VN4VNVev/iw2Tpw4EVL/EypW5/z2ARgO4DcA+jHzYGbuA+AJAMcA/B0RvelSGYNxEsAIIhpKRB3gFdJb3MxQzyyyceNGAF5tbOvWrX6/qV/+b775xpbjwZ07d0Iq4/r16/HZZ5+hpaXF1IwjRvjBzgtGa5nw3717t1/nq9X67M7bCG0nGG55yBo5EGRnZwOAnyYqhJ/H44HH48Hu3buVDraurg4nTpxQ5rWCORAFc7ASbenWrVu6JnUxfxuMvLw85Ofn+x1Tm2HtsG7dOuWz+jlrO+WHHnoIy5Ytw9tvv43ly5dj0qRJWLlyJQDvvRMWhnCwqp2ZQURoampCSUkJGhoacOHCBQCw7GgUbFAdDqJ9bNy4UbH8VFZWml4TioNUOFh9AvOY+W+Z+SwzKy2PmcuZeQMzvwTgK3eKaA4zewD8CsB3AC4B+JqZL7iZZ9euXXWPNzc3645utC+/nTkKoxddCFsr7N+/X/f4hQsXFMEszIDRIJR5UafQdsBCO1bDzNi+fbuu+TY+Pt61shmhHkzl5eXpniM6wvXr1wdofoLy8nLlc3Z2tt/zD3ddm9qM9cEHH+Cnn35SfisoKMCmTZtsmZvVGDlnBcPqAK1fv35ISUlBp06dkJycDMAZYaXGiQhT4hkxs99culXcFH6tYV2kpSfKzEH96K2c4xbMvIOZRzLzcGb+H27npzXDPPLIIwC8c01GSwbESwQAvXr1spyXthExM27cuGFbSws2yq6trfUzQUSSYcOGRSVfPcR82Geffebn3ShM1lqc6hQnTZqEQYMGOZKWFvX8oLo9qR2Y1E5HjY2NivB0AmYO6OzD8Yx0u2MdPXp0wDEjBx3BvHnzws7XyHu6e/fuWL58ecBx9SBatFsA6NChg995OTk5uoPyYHUKB7tOSNEg6JtLRPOJaDURTfR9d3cWspUxc+ZMlJWVATCfi1A3zrKyMstOINoXvba2VnEysIN6lK2ehwjHeUAP4c0YDLX2nJDgvNOx8GCzM9BQU1dXh9LSUjQ3N5t2tk51IMOGDcMzzzxj+PvMmTNDzlOY4UtLS/1MdkamS7VzSXtCtEOj+6s3Ty/asXYAp/YqFtYBrZOPdtnH9u3bFXOzus1NnDjRb/AsEEs8tGlp2/zJkyd1NU0n2m5FRQWuX78e8I6oly8Z0a1bt7DzDwcrw9b/BODPAbxJRHMATHS1RK0Qu43oypUrflpEVVUVsrOzsWXLloCFtC0tLfjhhx/www8/ANCfBLeLOg8jTTVUVqxYgXfeecf0nCeeeAKvv/46Ro8ejeTkZFdMh5MmTQIADBw4MOQ0qqqq8OGHH5q+xE5pfqLDmjt3rq5HsBPaTrA5FyMi1UnZrWO3bt0wcuRIR/JOSEjA008/jSeffNLwHPUARGA02Bs6dKjyWXh1pqamIiUlRTmu9mQWdb99+zZqa2v95oiDtTHtfdN7n8Qg9/bt24rW7YTw27t3L/bs2RNgDQMQdG508ODBpr+7jZU3t5SZK5n5zwAsgNfbs10zZox/ABk9M4mWpUuXBhwTo7ycnBxkZ2fjzp07Aa799fX1yMvLU+Z27ETZAELXfMIhmAu3uF+PP/44XnvtNcfnU9SMGzfO9jVCIxKdg/AQ1CMc4arH8OHDddOM5hyKeD4LFixwPG21R6vaoiEGL2Z0795dWaYQLj179sTQoUMxatQow3P0hIWZ1UJ4fycnJ6Nv377o2bOnX3vUix26c+dOfP75536Lwe2+H/Hx8cri/bt37yrpMzO2bduGn376CQ0NDSE7DgnKy8sVhyqjuedYxspdVXxrmfnXAIIvUGrjaDsi7SJ3PfQml0XIJXV6WlPU8ePHlc9Hjx61baZ88cUXbZ2vZtiwYbrmFjtMmjTJcIRJRK7NO3Tv3h1AaCZVsS5QzHvpCb/Gxkbs378/pPiVXbt2VdY3ATDtcAXqeoQTozMUkpKSMGHCBKj3wDRy07f7PNXLJdT1suKMsWDBAsesBnPmzAl6jl7d1MskMjMz/cotfAEWLlyIxYsXB9RJCLXm5mZTQRSsDWvn/xMSEnDw4EF8+OGH2LRpEwDgwYMHStSh27dv45NPPlGma9SCtrCw0M+haP/+/Yaev3ZCCmr7tRdeeEG5n04tsrdLUOHHzJsBgIhSfd//1e1CxTqhjMKNOoVvv/3WT6CZBeE100C09OvXzzRfq2gnz62yZMkSAN69E9977z0sXuxuxLmXXnrJ73tqaioyMzMtz0GaoReS6+OPP/aLxyl44403gqanjjC0dOlSU1ObYPTo0Vi2bBkAe52OE3To0CFgvkpr/RBYXb+qRsxJqt8rK0ItPj7eMauBlXaul9ecOXPw5ps/r/JSv2/ifKN1hHFxcairqwtqWg8WkUorHPPz83WjJglNUwg9wbVr1/w0T3X0nNzcXNy/fx93794NcFKyM7DUTuf069dPEfgPP/yw5XScxE7LWeNaKVoZzKwIl3C5d+9ewDomu+iZ9pyYCyEiLFq0CK+++qrta/v06YO33nrLNHyak1gNF+Umq1atsmSGU8dAFPcnGPHx8X7zRZFETzAYDapCMbMXFxfj66+/9tN+rAo1p4SflY5cr86JiYno0qWL7vmibEb3qrS0VBEoZppfsAGsE9aTffv2KYM5MTcs1kU2NTVh06ZN+OSTT0JegC4GM+r5Y1FnJwaooWCn5bjnF9vKaGlpUUYrVhuem27Fep3iqFGjHAkP1KVLF8WEqCVYx6M28zgh/LTaneC9996zbP56+eWXwy6HEXY64tawDiozMxMvvviirqOHUXt+7LHHQsqrsrLS755EWvhZaT/aOuvN46vP0btH2o5eaPFGEYBee+21oOUKd+4O8M7ZiTIIs6eI8qN+LuppGDv5FhQUIDk52c/cL+r+0EMPhV7wMLDTcmL/bY0QeiGb3GDVqlWWzlO/uMOGDcPkyZMjUj4xp2GFUM2naoI50lgh2kJn5syZmDZtWtTLYZXU1FTdkbloc2qvRsC5QZ5VoeaUw4uVcms7ez2tfcGCBVi0aJFhmunp6X6B4IUAuHLlSsC5s2fPtuRpqy5XOIJE7YCjjnOqLpvH41GWTdlpw0eOHAm4H6Lckeir9JCaXwh4PB7bjhShPGCr16iF37x584LGDHQK4ZU4fvz4oOe66XWqvU9OmaRDxUzj7tevX9BYklrUc0qxwPjx412PbGNV+KmdcNxGdPaZmZmGz7h3796mLvxEhIULFyoBDczmb63OhYk0nn32WVeWD6jnDx88eICSkhI/Jx2rAxDtexrtPTjtCL/fuFaKVoiw87s5grci/OLj49GvXz9MmzbN0GvQ6c5z8ODBeOqppzBgwAB07949YEsjt7Byr0eMGIEXXnjBVhpG4erC4Y033sCKFSsAAJMnTw7I3067MZpTcgo72tOTTz7pN8dMREhLS3O8TFaFXyiDyoyMDNvXAPbNi0ZlS0hIULyonTBZejweEBEGDhwYMS3K4/EonqRWLUDaZ6oWfnp9lNvWEcvCj5mjE/gxBlmwYIElrzZ1zEorjfLVV1+15PmnLUuXLl0wadIkw2utdCQiuojaJm9U5iFDhmDkyJGIj4/Hq6++qnijqa91AysvQ7D7rJdGqFqM2X1NSkpCx44dsWrVKkyZMgVz5861fK2aUJyN7KLn5m+05dCoUaMChKVWmDgxr6q+PxMmTAg7PTVDhgxRtO+nnnrKsqXEronVrC2K+hltFSYCaRuhNrlWVVUFtGurjlShot6sWVjBgvkYmGl+bg/w9LA1W0xEU4loIxH9RERniegcEQXfd6UdoLcweezYsbY6gu7du2PUqFGIj49XtiR56623As5TRwCxIlSDnTN06FDFDNOrVy9LZkw99PYJtIqYzzPTIp0QfnqE+uKJbXsSExMNteu4uDgQkRL+StRh1qxZpmswu3Xrhri4OENnIyfR3rNOnTqFJXCccEJRp+F0EAQiUp7DyJEj/TRzM5KSkmw5kZm1RaPfUlNTg65fnDFjhqF1Q8zbObVllxWsDh6163qDrVl2W/OzuwJ4Lbyhzs4BCF9fb0N069ZNNzBtMHdngdq+rx71de7cGf3790dlZSXq6uowZ84cP7d+K529nc6DiDBz5kycPXvWloOJ1U7BqDPv0qUL7t+/b9rg9UxEVjYHDsa4ceOU+6tXLu0aJYEoa0JCQlABqidgzOb+IqHxOYFe23JiPlCdRijCLzU1VdmuSw8nzI3BCEX4JSUlGc5jzpkzBz/88AOIyPCeiHXCkZxPszqfr61ztIWf3VZVysxbmPk6M98Qf66UrJVg1AiFJhDsxU1MTET//v1NR5/PP/+8omV07twZSUlJSufghOYnfo+Pj1fMJcuWLdN1cQdC99ycP3++pUgaesydO9fPySg9PR2vvvqq7fBi8fHxAQJ48ODBhvN+Zh25eDnVXnLBsPpCuxn9RiCEb6jzds888wweffRRS+e+9957ttJW33e7mnn37t3x9NNPm54Ta9624lmbPXOxdMiszQuhF0kPSvV7adYetIMRp4Pq28Wu8PtrIvqAiJYT0Yviz5WStRL0QvPExcUpL2wwzW/SpEl4/vnnbTtdWHlZtOcG+33lypWKp2RKSoqu5jdr1qyQ3amHDh1qOBchyjBy5EgMGzYsQJMcNmwYunfvrkQ5mTlzpm2T4AsvvICePXsGmGdFLESzchkxZcqUkM3E0UZ0pgkJCSE5gQwaNEg3DJl6cBRqZB91hxosYINW0L366qtBB52xpvlZsRANHDgQq1atQq9evfzOUw8U7PQLTtGjRw+l/YwePVoZSAVbvB5MOxXrDN3CrvB7B95dHRYCeN7395zDZWo1rFq1SneOSr2rQbBGHaqbsJ3GHawjsLM0IC0tzZUXSwj/UaNG6e6NJvIU5xmtfzIrm7jGTrBrs7iDzIwpU6Zg2jTrsd4jqXE895z5qylClGk1ebtLMQTCUqDu9EJdm6leEB2s/Y4YMUJZbygCIYhrhCYiHNTEPGsknkO/fv38tjZSow3cLdptsEhFevdi4sSJymfhKxBJ1JaDhIQEpc56GrudOT+3dxOxK/wmMPNUZl7BzO/4/t51pWStAKP5DqvzFStWrLCtRWlHdnYEkd6GmIDXMccqbiwLALwmUT3nHi2JiYl4/vnnDX/Xux/aaDx6Ak3vRXv66acDvDTD4aGHHoroLhvqgZk2UkifPn2U+6EOXt65c2fT+2uGmYlY/VyCrUV7/vnnkZCQgD59+ljOWyzzEfdXLdwyMzOxePFiZGZmKlr/Y4895ueN7QZxcXEBQQAE6gHC448/rtwfO2t0n3rqKcTHx2PMmDGK0BPpGgl3MRixGkAjFOLi4pCZmRm0b+rTp4+f17xYGgR4hbjra0ltnn+MiPQj2koABAo7M+HXsWNHy8JL7BwRivAjImRmZvp1ctGKE2lEYmJigAnNyGzSv39/W2lrAxLoBSh46qmnAo6NGDEi4AVUB3S2qz3MmTMn7DiGTzzxhOFv6enphnN3WguD0W4dXbt2tbSjgh5WHVOEUDbSqrXPV93GjawU2mchtFejwdrYsWNtDfrcpG/fviENZkeOHImVK1eic+fOftpfRkaG4dy6MEPHxcVZGmwGw6ytaOuife8WLlzo57WqfjdCtT7Ywa7wexzAaSLKdXqpAxH9PRFd9qW7kYh6+I6nE1EdEZ32/f1edc0UXxnyiOhfKFpxclQYCb9wi9ahQwckJSUp81zh2vZDjcEYSYR50qzDt4PZvbI6ylTPFcWa00RKSgoyMjJ0BRsRKSa1IUOG6Ar7cLHaFsV903ZwVoIlGGni2vmjDh064N13341oBJhQEcthnCItLU13kDBo0CC/dm60xZIeeuXr1q0bXn/9dcNrtMJO+44ZOXXNmDHDsU2KzbAr/BYCGAHvprZivi80G0kgewA8wszjAVyBf0SZfGae6Pv7per4fwDI9JVphK98UWHq1KkYPnx4QKfilPADvFFDhB09VOE3dOhQTJ061baXpPA2jSRirVmomogW9csdqkCdMWOGpf333CLU3QeISBlMzJkzxxWTUrhtXDuY6Nevn5/mJjydgUAhqOf8lJCQELW4kXZwWviJNAH/6EJ62pQ6X/W6UzPtW6Cd4tGinfOzWse+ffu6usG1wNI6PyL6EwCHAeQwsyv+qcys3sL8GADT1eFE1B9AN2Y+6vv+KYAlAHa6Ub5gdO7cWXdU6vbLZzd9own4YEyZMiXiAjCUDtrsfqgdO4YOHYqDBw/aSnvatGno2LGj4rIdDc3voYceQmpqKn744Qe/fdnS0tIU4WZ0D0aOHIlevXo5EmRcD3W+/fv3R3Fxsa3rtR6YvXr18tMshBa3fPnyAO22V69ejuxiEg3i4+NdE35DhgzBTz/9hCVLlihTHerBZNeuXVFVVYUePXr4vW/du3cP+vyCtf9Jkybh+vXrynerdXR7rk9gVbwOAvDPAO4S0Y9E9P8S0bNE5NYWvO/CX4gNJaIcItpPRGLIPhBAkeqcIt+xmMIt12MnNMpIOl+Eg1WvrwULFgR40TnB448/7vc9lN3hnYKIkJKSErC90/Tp05WR/RNPPBFggSAixMfH23IisYt6tD5o0CAQka62oTbJqR0vBgwYYNomRWdrNF/ZWomLi3Nc0xHpibaanJyMxMREJCYm+s31CQvQ6NGj/a7XzpV27txZiQIlCCb8UlNTAxx+gi0xee2111wPzSaw9BYz858BABF1ADAVwGPwCqjVRFTJzJacYIhoLwC9Geu/VO0Y/5cAPPBGkwGAYgBDmLmMiKYA2EREY6G/y4Th0yCiTHhNpEF3RnYaN0ek4XTEo0aNwpEjRxwsjfPYuXduze+MGTMGhw4dUr4PHDgQJSUlMTPnp71HkQo0rkU9EJs4cSImTpyoOFupUQs/rcA0WzYSKY0g0rih+Yn0xD0z6idmz56Ny5cvo1u3bn5l0Ds/PT0dRUVefSMxMdG257d6E2cj3F7eoMbucKMzgG4Auvv+bgM4bnqFCmaex8yP6PwJwbcC3nnEN9h3l5i5gZnLfJ+zAeQDGAmvpqceigzylcco7yzfMo2pkRpZuIloqNHUQtoSdoIFCHf0SL6oRoSzn52TwnvkyJF+QtfImeGVV14JOGZ1gb1b5tpoIUy6bsz5aS1OZv3EypUr/dbvGgUJUJfx9ddf112PG4xYGTAC1uf8sgCMBVAFr7A7AuAfmbnCqYIQ0UIAfwHgSWauVR3vDaCcmZuJaBi8ji3XmLmciKqIaKavTG8D+FenyhPrSOFnnWHDhqGmpsbvWLgb465cuTIik/LBiIRXnBWseo/qdfJubInUGhBamVvt6LnnnkNycjLmzZtny9M5ISEBo0ePNrQKrVixIuQlO5GIrGMVqz3nEAAdAVwFcAteravS4bL8zpfHHt+DOubz7JwN4G+IyAOgGcAvmbncd80fAPgYXo10J6Lk7BINRGMOpxMfOnQoqqurnSpSzNKvX78A128zE1pycjKqqqpM04wFE9zkyZMd30aqNXhHthU6deqE559/3rW2JDRxsaNIMNTPXpSpa9euSh8hfrcj+LSa3sSJEwPmDqOF1Tm/hb41dGPhne/7UwCPEFE5gKPM/NfhFoSZdUOdMPMGABsMfjsFwNpOim0MMYIKZ9SYlJRkGLy6PTN16lTs27cv4HisCQY70UBiBaN7uGjRItMwe/PmzbMd3CDWISKlTrFgRdBzzpswYQIOHz4MwDuItOsk169fPxQUFCjfO3bsaHuZlVvY2cyWfRva7oBXwzoMYDiAP3apbBITIrlliUTiFEbCb/DgwaaDi2HDhjm23jMWiaWBlbosas/aHj16BHgZByOWg75bnfP7z/BqfLMANMEr+I4CWAPv3n4SG0yYMAFnzpwJK41obwciaf1Ee/4llpwfos2TTz5puG9kpNAK4JUrVyrmz3CeFRFZ8vSMNFbn/NIBfAPgvzCzvZWrEldIS0tDXl5etIvR5hg0aFDI2xu1NqLRGbW1e+gU6gDP0UJr9nRqLjIuLi4mLVVWzZ5/yszfmAm+WIir2Z6YPXs23njjjWgXo81hd4/A1kysjcQlsYHTXfnChVGLOmmKVeG3j4j+iIj8VocTUQcimkNEnwBYYXCtxAUSEhLCWuMlMcbo5W9r4zut2XPYsGGuL51Qr9WLtZ1F2jtute9YcXDRYtXsuRDeiC5fEtFQeJc5dAIQD2A3gH9i5tNuFFAiiTRtbTG1EVrhF8qiZTuoI728++67MeHhKPkZt0IxxipWlzrUA/h3AP9ORIkAUgHUMXOli2WTSKLCkCFD2oVJOT09HZWVlVHJWwZniF2k8DOAmZvgjbcpkbRJiKhdmJR79erl6E71ktZNexF6Aml3iAJy1Bu7mHUA7a1zkLQv3DZ7xpqDlSXhR0Sf+f7LBe0OMGHCBNuLRSWRwcy9Wwo/SVumvc35WdX8phBRGoB3iSiFiHqq/9wsYFtiyJAh6NOnDxISElrNXnrtDTOtXLszdWvgjTfesL31jKR90t6En1X72+8B7AIwDEA2/PfSY99xSRBidb2LJDgrVqxolV6gSUlJfruhSyShEGsmSyew6u35LwD+hYj+g5n/wOUySSQxR6hbuEgkrQU552eCFHyStk4sbFArkUjcx5bbIRF1BPASvLE+lWuZ+W+cLZZEEnlWrVrVbuY7JBItcs7PnM0A7sM779fgfHEkkughI45I2jNE5BeFp61jV/gNYmbptSGRSCTtCLd2m48mdoXfESIax8xyDz+JRCJpB7zyyivtV/gR0Tl4lzQkAHiHiK7Ba/YkeDd5j93teiUSiUQSMj169Ih2EVzBqub3HLzCD/AJPKcLQkS/BfAegFLfof/KzDt8v/0GwEoAzQD+MzN/5zs+BcDHADoD2AHgjznW/GklEolEEnNLHawKv/PQF3hCEDrlH/5PzPwPfhkQjQHwGoCxAAYA2EtEI5m5GcB/AMgEcAxe4bcQwE6HyiKRSCSSNorVRe7JbhfEhMUA1jFzA4DrRJQHYDoRFQDoxsxHAYCIPgWwBFL4SSySkZER7SJIJJIoEWu+3b8iorNEtIaIxDbPAwEUqs4p8h0b6PusPa4LEWUS0SkiOlVaWmp0mqQdkZaWhrS0tGgXQyKRRIGICj8i2ktE53X+FsNrwhwOYCK8+wX+b3GZTlJsclwXZs5i5qnMPLV3797hVUQikUgktmitc36OwMzzrJxHRKsBbPN9LQIwWPXzIAC3fccH6RyXSCQSSQzRrVu3mAsMHzO7qhJRf2YWO8QvhdfJBgC2APiCiP4RXoeXEQBOMHMzEVUR0UwAxwG8DeBfI11uiUQikZgTi/uXxozwA/C/iGgivKbLAgD/BwAw8wUi+hrARQAeAH/o8/QEgD/Az0sddkI6u0gkEknMkZiYGO0iBECxZoeNBFOnTuVTp05FuxgSiUQicREiymbmqbq/tUfhR0SlAG6EmUwqgHsOFCeWaIt1AtpmvdpinYC2Wa+2WCegddQrjZl1PRzbpfBzAiI6ZTSiaK20xToBbbNebbFOQNusV1usE9D66xVr6/wkEolEInEdKfwkEolE0u6Qwi90sqJdABdoi3UC2ma92mKdgLZZr7ZYJ6CV10vO+UkkEomk3SE1P4lEIpG0O6Twk0gkEkm7Qwo/mxDRQiLKJaI8Ivp1tMtjFyIqIKJzRHSaiE75jvUkoj1EdNX3P0V1/m98dc0lopjYA8i368ddIjqvOma7DkQ0xXcv8ojoX4hIL1h6xDCo12+J6JbveZ0momdUv8V8vYhoMBHtI6JLRHSBiP7Yd7zVPi+TOrX2Z9WJiE4Q0Rlfvf6773irfVamMLP8s/gHIB5APoBhADoAOANgTLTLZbMOBQBSNcf+F4Bf+z7/GsD/9H0e46tjRwBDfXWPj4E6zAYwGcD5cOoA4ASAR+HdIWQngEUxWK/fAvgznXNbRb0A9Acw2fc5GcAVX9lb7fMyqVNrf1YEoKvvcyK8MZNntuZnZfYnNT97TAeQx8zXmLkRwDp4N9tt7SwG8Inv8yfwbgosjq9j5gZmvg4gD957EFWY+QCAcs1hW3Ugov7wbYbM3rf1U9U1UcGgXka0inoxczEz/+T7XAXgErz7brba52VSJyNivk4AwF6qfV8TfX+MVvyszJDCzx5GG+u2JhjAbiLKJqJM37G+7NtRw/e/j+94a6qv3TrY2gw5yri2yXMkIaJ0AJPg1SjaxPPS1Alo5c+KiOKJ6DSAuwD2MHObeVZapPCzh60NdGOUWcw8GcAiAH9IRLNNzm0L9XVkM+Qo4uomz5GCiLoC2ADgT5j5gdmpOsdisl46dWr1z4qZm5l5Irz7o04nokdMTm819dJDCj97GG2s22pg5tu+/3cBbITXjFniM1XA9/+u7/TWVF+7dWgVmyEzc4mvQ2oBsBo/m51bTb2IKBFeIbGWmb/1HW7Vz0uvTm3hWQmYuRLAjwAWopU/KyOk8LPHSQAjiGgoEXUA8Bq8m+22CogoiYiSxWcAC+DdNHgLgBW+01YA2Oz7vAXAa0TUkYiGwreRcGRLbRlbdfCZb6qIaKbPE+1t1TUxg+h0fGg3eY75evnK8CGAS8z8j6qfWu3zMqpTG3hWvYmoh+9zZwDzAFxGK35WpkTb46a1/QF4Bl7vrnwAfxnt8tgs+zB4vbPOALggyg+gF4DvAVz1/e+puuYvfXXNRYx4bAH4El6zUhO8o8yVodQBwFR4O6h8AL+DL+JRjNXrMwDnAJyFt7Pp35rqBeBxeE1eZwGc9v0905qfl0mdWvuzGg8gx1f+8wD+yne81T4rsz8Z3kwikUgk7Q5p9pRIJBJJu0MKP4lEIpG0O6Twk0gkEkm7Qwo/iUQikbQ7pPCTSCQSSbtDCj+JRCKRtDuk8JNIJBJJu+P/By1twl0/tX+zAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 3 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(3, 1, sharex=True, sharey=False)\n",
    "fig.tight_layout()\n",
    "ax[0].plot(res[\"tau\"], \"-\", color=\"grey\", linewidth=1, alpha = 0.8)\n",
    "ax[1].plot(res[\"sensible\"], \"-\", color=\"grey\", linewidth=1, alpha = 0.8)\n",
    "ax[2].plot(res[\"latent\"], \"-\", color=\"grey\", linewidth=1, alpha = 0.8)\n",
    "ax[0].set_ylabel('tau (Nm$^{-2}$)')\n",
    "ax[1].set_ylabel('shf (Wm$^{-2}$)')\n",
    "ax[2].set_ylabel('lhf (Wm$^{-2}$)')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "You can save the output in a csv file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "res.to_csv(\"test_AirSeaFluxCode.csv\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and generate some statistics which you can save in a txt file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Input summary\", file=open('./stats.txt', 'a'))\n",
    "print('input file name: {}, \\n method: {}, \\n gustiness: {}, \\n cskin: {},'\n",
    "      ' \\n tolerance: {}, \\n qmethod: {}, \\n L: {}'.format(\"data_all.csv\", \"S88\", \"on\",\n",
    "                                                           0, \"all\", \"Buck2\", \"Rb\"),\n",
    "      file=open('./stats.txt', 'a'))\n",
    "ttl = np.asarray([\"tau  \", \"shf  \", \"lhf  \", \"L    \", \"cd   \", \"cdn  \",\n",
    "                  \"ct   \", \"ctn  \", \"cq   \", \"cqn  \", \"tsrv \", \"tsr  \",\n",
    "                  \"qsr  \", \"usr  \", \"psim \", \"psit \", \"psiq \", \"u10n \",\n",
    "                  \"t10n \", \"q10n \", \"zo   \", \"zot  \", \"zoq  \", \"urefs\",\n",
    "                  \"trefs\", \"qrefs\", \"dter \", \"dqer \", \"dtwl \", \"qair \",\n",
    "                  \"qsea \", \"Rl   \", \"Rs   \", \"Rnl  \", \"ug   \", \"Rib  \", \"rh   \"])\n",
    "header = [\"var\", \"mean\", \"median\", \"min\", \"max\", \"5%\", \"95%\"]\n",
    "n = np.shape(res)\n",
    "stats = np.copy(ttl)\n",
    "a = res.loc[:,\"tau\":\"rh\"].to_numpy(dtype=\"float64\").T\n",
    "stats = np.c_[stats, np.nanmean(a, axis=1)]\n",
    "stats = np.c_[stats, np.nanmedian(a, axis=1)]\n",
    "stats = np.c_[stats, np.nanmin(a, axis=1)]\n",
    "stats = np.c_[stats, np.nanmax(a, axis=1)]\n",
    "stats = np.c_[stats, np.nanpercentile(a, 5, axis=1)]\n",
    "stats = np.c_[stats, np.nanpercentile(a, 95, axis=1)]\n",
    "print(tabulate(stats, headers=header, tablefmt=\"github\", numalign=\"left\",\n",
    "               floatfmt=(\"s\", \"2.2e\", \"2.2e\", \"2.2e\", \"2.2e\", \"2.2e\", \"2.2e\")),\n",
    "      file=open('./stats.txt', 'a'))\n",
    "print('-'*79+'\\n', file=open('./stats.txt', 'a'))\n",
    "del a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}