create_daily.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "e3a39119",
   "metadata": {},
   "source": [
    "### CREATE DAILY INPUT FOR THE ML MODEL"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d36f5f34",
   "metadata": {},
   "source": [
    "Input:\n",
    "- the measured runoff at the gauging station,\n",
    "- the catchment area of the gauging station,\n",
    "- the ERA5 datasets for Potential Evapotranspiration, Temperature and Precipitation.\n",
    "\n",
    "Output:\n",
    "    - a .csv file is saved in conda/daily_input folder."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "889e3048",
   "metadata": {},
   "outputs": [],
   "source": [
    "from create_daily_csv import interpolate_df,readnetcdf_in_shp,xarray2df,check_data_gap\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import os"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "27466b2b",
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "shp_catchment_path = r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\Documents\\qgis\\final_delineation'\n",
    "shp_catchment_fileName =shp_catchment_path + r'\\ADIGE_BRONZOLO.shp'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "1a9d9e92",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Runoff\n",
    "catchment_name = 'CANDOGLIA_TOCE_2000_2019.csv'\n",
    "fileName = r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\data_discharge\\\\' + catchment_name\n",
    "#runoff = pd.read_csv(fileName, index_col=2, parse_dates=True)\n",
    "runoff = pd.read_csv(fileName, delimiter=',', index_col=1\n",
    "                     , parse_dates=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eb258f6c",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "f742defa",
   "metadata": {},
   "source": [
    "runoff['discharge_m3_s']=runoff['pretok (m3/s)']\n",
    "runoff=runoff.drop(columns='pretok (m3/s)')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4f649620",
   "metadata": {},
   "source": [
    "runoff=runoff.drop(columns='temp. vode (°C)')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "292a335e",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Date start: 2000-01-01, date end: 2019-12-31\n",
      "No missing dates\n"
     ]
    }
   ],
   "source": [
    "check_data_gap(runoff)\n",
    "runoff = runoff.reindex(pd.date_range(runoff.index.min() - pd.DateOffset(years=1), runoff.index.max()), fill_value=np.nan)\n",
    "runoff = runoff.rename(columns={'discharge_m3_s': 'Q'})\n",
    "runoff.drop(['id_station','quality'],axis=1,inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0af1fa25",
   "metadata": {},
   "outputs": [],
   "source": [
    "# ERA5 temperature\n",
    "\n",
    "era5_fileName=  r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\era5\\2m_temperature-19790101_20181231-eusalp-qm_era5.nc'\n",
    "\n",
    "t2m = readnetcdf_in_shp(era5_fileName,shp_catchment_fileName,plot=True,res=5500)['t2m']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "0e6f6ceb",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n of pixels287\n",
      "Date start: 1979-01-01, date end: 2018-12-31\n",
      "No missing dates\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEQCAYAAAC+z7+sAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAABG2klEQVR4nO2dd3xUVfr/32dKJplkJlQLJYlKkbKi4mLBgmKvu3YFZX+LhbUtqyj7da2srgVFdK2IuiJFXcUuimtvqCDggiiKLAgWQMD0TCbz/P64N2ESZpKZZCZT8rxfr/OamXuf+znPvXPnuXfOPc85RkRQFEVRshdHqh1QFEVRkosGekVRlCxHA72iKEqWo4FeURQly9FAryiKkuVooFcURcly0jbQG2MeNcZsMMYsi9H+dGPMl8aY5caY2cn2T1EUJVMw6dqP3hhzMFAOzBCRwS3Y9gWeBg4TkS3GmB1EZEN7+KkoipLupO0dvYi8B2wOX2aM2c0Y85oxZpEx5n1jzO72qvOB+0Rki72tBnlFURSbtA30UZgGXCoiQ4EJwP328n5AP2PMh8aYBcaYo1PmoaIoSprhSrUDsWKMKQAOAP5tjKlf7LFfXUBfYATQC3jfGDNYRLa2s5uKoihpR8YEeqx/H1tFZM8I69YBC0SkFlhtjPkaK/B/1o7+KYqipCUZ03QjIqVYQfw0AGMxxF79PHCovbwbVlPOd6nwU1EUJd1I20BvjJkDfAz0N8asM8aMBUYBY40xS4HlwEm2+evAL8aYL4G3gStF5JdU+K0oipJupG33SkVRFCUxpO0dvaIoipIY0vJhbLdu3aSkpCTVbiiKomQMixYt2iQi3SOtS8tAX1JSwsKFC1PthqIoSsZgjFkTbZ023SiKomQ5GugVRVGyHA30iqIoWY4GekVRlCxHA72iKEqWo4FeURQly9FAryiKkuVooFcURclyNNAriqJkORroFUVRshwN9IqiKFmOBnpFUZQsRwO9oihKlqOBXlEUJcvRQK8oipLlaKBXFEXJcjTQK4qiZDlpOcNUa5k2DWbPbtnu3Xet10MOSZytaqqmaqpmIjRFWtaMl6y6o589G5YsSbUXiqIo6YWRZFw+2sg+++wjrZkzdsQI6/Wdd5q369TJet26tWXNWG1VUzVVUzXbYhtr/IqGMWaRiOwTaV1W3dEriqK0yKxZUFICDof1OmtWYmzTmKxqo1cURWmWWbPgggugstL6vGaN9Rlg1KjYbWlim+boHb2iKB2Hv/1tW+Cup7LSWt4W2zRHA72iKB2HtWtjXx6PbZqjgV5RlI5DUVGLywOBp8jPz+eOwsL4NNIYDfSKonQcbr4ZvN7Gy7xeuPlmVqxYQXn5CKqqzqSyspIrt26lrzGcAlwCvBdmm2looFcUpeMwapSVWVlcDMZYr9Om8WFJCUOHDqWu7l0cjj788MMP3HjjjQS6dGEucB9wCOANBik94YQU70T8tNiP3hjTG5gB7ASEgGkicncEuxHAVMANbBKRQ2Ldtimt7Uffvz+sWgUFBc3b/fqr9Rrtn1lrbFVTNVUzczRra5+lpmYq4KSu7iOg1l4zHzhiO826uqVUVV1MXd2HuN3nUlv7eNL8bG1qU3P96GPpXhkErhCRz40xPmCRMeYNEfkyrIJOwP3A0SKy1hizQ6zbJpI0zP1SFCVFiFRSUXEE8BEApaW98HguJxB4nFBoqW3lxOHYDZfrRHJyzqO8vH9ELadzCAUFH1BWthu1tTOAU4AT22M3EkKLgV5EfgR+tN+XGWNWAD2B8GB9NjBXRNbadhvi2DZh9OhhFc2MVU3VVM3x46/m7rs/avjs85VSWno5AIMHD2bOnDkMHjw4Ls11696ld+/eOBwT8flOTEpmbDKIK2HKGFMC7AV80mRVP8BtjHkH8AF3i8iMGLetX38BcAFAUQY+1VYUJX3YtGkT06ZNY+DAgaxbtxBj8vj55xo++OADevbsye677x63pohwyy232J8CiXU4ycT8MNYYUwA8C4wXkdImq13AUOA44CjgWmNMvxi3BUBEponIPiKyT/fu3ePcDUVRMopYhxaIYwgCmTmzwfb6XXahqqqKJ554AmPyAPB4PIwcOZLdFy2KewiEVcbgcDi4//77KSkpIT//3VbtdqqI6Y7eGOPGCtSzRGRuBJN1WA9gK4AKY8x7wBBgZQzbKorSkYh1aIEYhiuoq/ucHXY4mo0bNzZsdiQwv7ycTsDeK1YAe8elGcn2KXvRpS4XUydNosulvVq//ymgxTt6Y4wBHgFWiMiUKGYvAAcZY1zGGC+wL7Aixm0VRelIxDq0QAt2IkJ5+dCGIN8Vq934M2AkMK1eozV1N7GtthfdHQziuPbalvcxzYil6WY4cA5wmDFmiV2ONcaMM8aMAxCRFcBrwBfAp8B0EVkWbdvk7IqiKBlBrEMLNPlcC0jY8traRwG48MILEWPYBJQCm4H/AKfFoNns8rBln9uvdc1ppDGx9Lr5ADAx2E0GJrdmW0VROhBFRVaTSaTlW7d9XLnzzpz8ww+sBsLvwX2Av1cvqqrWA3DdddfBa6/FpNls3c34ucJe5IykmQFoZqyiKO1LM8MQANTWvsAbb7zB2IIClgM5QF/b7ACHg3y/n1122QW3+3x8vnX06NGjRc1Y645kGwRWA8cDRodAUBRFiYEowxAwahTB4FtUVv6OI488kg9WruSSI45gS3ExK41Biov5cMYMfty6lffffx+vdxoOR88WNWOtO5qfL3TrhgBndesW3TbNyaqpBHUIBNVUzczSFAkQCNxFIPAoodDKhuUu19E4HP3xeK7A4ejdzn5W4fW+R03NTUAtdXWfAPn4/Rsbumom83imagiEjCENr1mKokShqmocgcBDEdZcQ37+39vFh0DgX9TU3IUxeeTkXAqsByY26pjjdA4nJ+eyhiCfiWRVoNchEFRTNTND8+6772b8eCvIP/DAA1x44YUEAgG6dv0Ip3NEu/g5b948jj32/wFWMlVV1WgAjOnJlCkTOOWUU+jdO/K/iWRODp4MtI1eUdKdTJnMOoa6Q6EQo0aNYvz48TiBpcC4W2/FzJ6Nx+PB5ToUK/0mds24bWfN4v8KCzn2WKun9+sTJ1JaWsqf/vQncnIm4PN9w/jx460gnyWTgyMiaVeGDh0qreGQQ6zSEoWFVomFWG1VUzWTojlzpojXK2K1TFrF67WWx2Gbaj/9/iopKPhCTjjhBAGkvzFS0ZKfydj3mTNlvscjgHQDWZBGxzPW+BUNYKFEial6R68o6UymTGbdTN01NTWUluZRXr4HL730EiNzc/lKBG8E21g146m/qd2smhoAPsBK4U/L45lgNNArSjqTKZNZR6ijHDh0zRpyc3MBMKaE559/nueqq7ezjaiRjH1fu5Zh9tvpidLMADTQK0o6E8Nk1uHLVmBNezcaWNmcbaKJUMcxwDtA165dycm5GJ/vK0466SR8xcWxacS577Fq/sl+O705u3jrT3M00CtKOhNjJueWLVv4rcPBQKyJrGdhjRv+VW5u+2Ryhvl5J9bcoR8Ax++1F5s2bSIv716M8Wxn20CCslhj0TReLxdijWKwNBGaGYAGekVJZ1rI5AwEplNRcSolJSUsXL2aE/bai8922onpWE0nA6qraZd+IrafH++0ExOAn4FdundnVqS+zknMYo1V8289rYza2X5/2zUzAM2MbYFMzzxUzezVrKw8i9raJwFwOg8jJ2csOTlnN6yvqXmI6upx9qcLgWsoLLTGURcpJxh8H5drf4zphIg0dGtsrv5QaDV1dd/gch1Maen7QC5+/3ACgckEg28RDM4HHPh8P+Bw7Ji0fU+EZmlpF0S2UFCwFKdzj7TxUzNjWyANr1mKkhTq6pbYQT4H2EpBwfZZmx7PhbhcB1FePgR4CHjInizbRTA4L8zSD5RiTA/8/vXN1ltW1h9rwOBtlIbNGedyHY3b/ftGQT5d8Xj+SnX1ROrqVjQE+mwlqwK9ZsaqZiZo/vJLHQ6Ho1Fi0Nq1a1mwYAFer5ezz/4NkMvGjV1wu90NNps2beLVV1/l1ltvZcUKa+DcvLz7yMnJa6b+gYgE8HqnUl19OcHgGwAMHTqU/v3707VrVxYtWsRHH32EyA9Mm/Y0F1xwesR9uvjii7n/fivIn3/++Tz22A8YsyOHHLKWAw88kCuuuIKCKH+n0/E7+uab39Ov30SmTatl9OjU+5k2k4MrSkdERKiuriYvr/Fdc3V1NSKexpmcLRAMfoDLdRDdunXjwQcf5JRTTmHVqlUMGTKEioqKRrb5+W6OOeYYnn76aa666iruueeehnUnn3wyF110EaecMrLFOo0xeDx/weU6jjVrutGlS5ftbCorK8nPz2fatGnA6dut/+9//8v999/PoEGDWLJkCS6Xi6eftta98UbMu5+WLFy4kNH1kT5L0YexSvaRwImnKyoq6NGjB16vl4Nzc/m3MRyWm4sxhry8PEpLHYhsiUnTumBcAUBVVRWnnnoqxhj69OlDRUUFlx55JP/85z/Jzb2LnJwJ1NbW8uKLL5Kbm8s999yDA5jWpQubH3qIZ599lpEjIwT5Zup3OvttC/JN7LzPPce4ceN48803CQQe207zimFW7/Pnt2zB9dRTRCRJwxUkS/Oefv0AyF+9uu2a6U60lNlUFh0CQTVbrRlHKnw0u+uvv15yci6W3Nz75bDDDhOsGeyiFofjNxIKhVpM2fd65wkgt956q/zy4INyttMpgAwDmR/Bz1AoJHeefbYMdzhkLGwbMiABwwBEsltwww0N+5ST82dZuHChXPf730tPYwSQPVIwXEEyNa+397WHMbL6rrtarxnP+dkMyRwCIeVBPVLRQK+ardYsLm78w6wvxcVR7TaBTAOZDOKMEMgP9XhEQFaDDLdLta3pdp/eYLd8552j1i0i4nafJ4CsX7++VX5G0mzrvje1mz17toBvu2PwZ5Cv0sjPRGiGQHYP28eioiIpKPgifs0WSIdAr003SnYR58TTFcDuwAXAlVhzgu633374fD+Sm/tP5syZw3x7bJQSrCSgDwCPrZGXNxO3+48AnPDjj9bk1U1YYg8DUFs7HYejDzvvvHOrJ8hudnkCNM866ywKC0spKFjRsPgzYCrQP438TISmAb4EZgK7Yj0QLy/fg4amOB0CQVHSlBjT1qV3b64GCoBNwETgPWBz7958/PHHOBw74fFcwplnnomrmZR9Y9x4vY9w3nnn8R3WD+oqoAz4BDgH2AuoqanB6RxKQcFn1sPbOFL2Y97PBGo6nbtbd4PFxUTsmJ0mfrZV0wCjgFXFxdxxxx0AVFaOib/+NEcDvZJd2Gnr/wOOA34HTHK7CU6aRE3NFKqqLqRnz5441q7lFmAfYDJwC3CQ10v+LbdE1WxEk1T4O+64g1EHHEB/Y5iM1TN9P6y7RbfTyVNPPUVBwUKM6RSzZlx2qtlmzcsvvxyHY3eCwZdYvnx5Vg2BkPL2+EhF2+hVszWaixYtkhkzZsiSm2+WPi5XCw9QHfL3U0+VUFGRiDFWe6z9kC9i/TNnWjZNbCPZze3eXUaBTCoslKcuuUQCgUCbNSPZqWbiNfPypgsgzz77bHyazZAObfQ6BEILpENKtGrGZutwDCAU+speasjNvZucnNGUltb3Gz8WeIHCwpbTRzJt31Wz7XYiAUpLi4BNFBR8idPZLyV+tjYk6xAISoegvlnE7f4jHs/FOJ17A1BQsJJQ6DsqK49KoXdKumPlOPwM3NNskM9EsirQ6xAIHVszN3c1NTWd2LLlkSYWfYG+aeOnaqafpojgcNwLgN//x6wbAkEfxiqZQQwZin369GHr1q0EHn885VmXqpk5mqEnnmDHHa1B2HJyLsSY/OT6mQqiNd6nsujDWNWMNUMx3Hbs2LECyJt2glNT29ZotmSrmpmvWZqX1/CQ3u0+VxyO3WXAgAFy/qGHyqycHLkI5MlE+dkMmhkbIxros1QzhgxJn2+TdOvWTZxYGaypzrpUzczSvLtz54ZgH62EEuFnM6Q0M9YY09sY87YxZoUxZrkx5s9R7EYYY5bYNu+GLT/aGPO1MeZbY8xf2/gHROmItJChGAqtpaJiGJs2beJGrAzWFjUyMJNTNZOnednWrYgIPt9aCgq+ZMuWLbwFPAwcbdu8nQg/U0QsbfRB4AoRGYCVA3KxMWZguIGxujvcD5woIoOA0+zlTqy5io8BBgJnNd1WUVokLBPxfGAIcDAwwefjpJNOoqysmFDoO2644Qb+luKJp1UzszUdjt44nQPo1KkThxYXcx5wl23ybQSNlStXsr5Hj9j9TBEtBnoR+VFEPrfflwErgJ5NzM4G5orIWttug718GPCtiHwnIgHgSeCkRDmvdBDsDMX3genAF8D7wJ2lpbz44osYswN5eU9y/fXXp23WpWpmruZuWN0TFzexram5nf79+9Nr/XqcWMMpGKxhLyQvL70yaKO16UQqWP+K1wL+JsunYt25vwMsAs61l58KTA+zOwe4N4r2BcBCYGFRUVGr2qi0jT6LNWfOlAe6dBFAPt5pJ9k6bZq888478u6774rfH9jONh2zLlUzczVH2O30vzz4YIOdMUUCyKhRo+SIwYPFE9ae//dTT42uGYW0yIw1xhQA7wI3i8jcJuvuxRo2ZCSQB3yMNdTIEOAoETnPtjsHGCYilzZXl2bGqmZTu7q6ZVRWnkMotIT8/P/gco2MaptKP1UzOzVra5+jsvJknM6DcLtPpbraelTpdp+O17ttIhaRckpLfYCPwsLSVvkZY0jejuYyY2PqR2+McQPPArOaBnmbdcBrIlIhIpuwBgIcYi/vHWbXC/ghHufjobUHSElvQqFNlJfvQyi0BACHo29qHVI6HG737/F4/kpd3fsNQR4OJi/voUZ2xhTgch0DlFFX92W7+xmVaLf6sq1JxQAzgKnN2AwA3sRqyvICy4DB9ufvgF2wpqtfCgxqqU7tXqma4XYffvihgDUzUygUSls/VTP7Na+66irp0aOH5OXNjmr3wAMPNDTh3HHHHWnRdBPLHf1wrLb1w+zuk0uMMccaY8YZY8bZF4sVwGtYz8k+xWqXXyYiQeAS4HWsh7hPi8jytl6clI5F586dAQiFQnFNxK0oiea2225j/fr15OScFdVm3LhxPPHEEwBMmDCBUGhDVNv2IpZeNx+IiBGRPURkT7u8KiIPisiDYXaTRWSgiAwWkalhy18VkX4ispuIpNFjaCVT2H333QG4+uqrMy69XjU7puZoY7jtttsAECltff2JItqtfiqLNt2oZiO7mTNlmMMhgKxqKcU8zdLrVbPjaj55ySUCiNf7XEyaOgRCjGigz1LN4mK50m7zfCkD0+tVs2NqVhdZ3S8B8fvrWtRMdRu9oqSUmjVrmGy/3zN8RYak16tmx9T0fP89newxisvL9+Dzzz9P2XAJGuiVtGeo2w3AFVj9cxvIsPR61ex4mu+88w6QRyi0nBkzZqRswnEN9Epas2LFCpbX1nKIw8Ed4SuyJL1eNbNbc8iQIfj9GwG791iqJhyP1qaTyqJt9KpZb3faaacJIJ/eeGNapMKrpmrGq+nxTBRAZtYvi2KbFkMgtCc6BEJku7q6VQSD88nJ+QPG5KWtn4nUhE4Y0xm/f3XCNDNl31Uz8zUDgaeoqjoT6I7f/xPGRG9ESfkQCJlCGl6zEkYotI7y8j5UV19EWdlAAoFZiPza8oYZzUfAr+Tk/L9UO6IocRMKrbWD/M7A4maDfLLRycFbIJWTGhcW1hEKrWbKlPe4+OKLAdh7771ZsmQJVVWjqaqC/fbbjw8//JAuXRwp8zMZmn7/WsrKhgPwyScnMWRI2zUzZd9VMzs0DznkHN57D/LzZ+Jy9dTJwbOWVmTV/dy7Nwf064cxhtJSF+XlfRk7dizV1dVMnjyZRZdfzuaePXkC6OtysWDBAkaOHNmsZlpOVtwCodBPAFx44YUMGTIkYzMkVbPjam7YYA194HAMan39iSJa430qS1Y8jI0jq646L08eAPk9jeeo9HgmisczSV5++WVZt27ddpobbbvT9923dRl9ydr3BNj6fD8KIFdffXXGZ0iqZsfUfOayywQQt/uMmDQ1MzZG0irQx5FVd36TAH9XJLsImh+FbeP3V8Sf0ZesfU+QJiD77LNPxmdIqmbH1bR+n+6YbDUzNhOJIQNuzZo1/G7NGh7Gmoy3DCtqj4+2fZNlQ8Lel5bmU1a2N48++mhGTFbcEiJ1ACxcuDDjMyRVs2NqVqxZA1hj1Leq/gSigT5ZtJABFwg8RklJCS/Yi6cBBRHsmtP0Yl0Y/uzzAbmEQosZO3Ysn+24Y3w+pSXVDe9e6949skkGZUiqZsfTHGZndHs8t7Su/gSigT5JXL/33gwwhhLgSOAF4PvcXK4aPJiysj5UVf2R3Nxc3rvmGsTr5TfhG8eZqTf1gQcoLKwiP/9TAIb99BMlxnAbUNGSZppiTD4+33oArvf5Mj5DUjU7lmZdXh5f1tYCkJMzBpFyFi1axDujRlGa1zgHRjNj4yRd2ujnzp3b0HbezR5et2kxprusXLnS2iDW7LsYbB977DE55JBDGtV1sMcjX9xyS7vse6I1d911VwEkOGNGxmdIqmbH0jziiCMi/vbryxf1zwY0MzY+0iEzNj//TSoqDsfh2I2Cgv9iTB6h0Dqqq6/BmG64XIdTWXkUYJKa0VdX9xWBwD2EQlY2rTE98flWUlrqbbVmMvyMblcHHAG8DYDfX44x+Wnop2qqZmTq6r6mpuZ2ams7AW/gcu2IwzGQQOAeAPLyHiIn54LtNFsbkpvLjM2qhKnUX7PKqag4CfDg9b7WMEyBw9ELr/df7eqJ07k7eXn3A1BdfSM1NTdQVXUJ8Gi7+tEaQqF1wK5ALU7ncHJz74wa5BUlXXE6++P1PhJ2EwjB4H/sQO/B7T6z3XzJqkDfXpmxdXV1TJ48mY8//hiv18uaNWtYsGA1IlaSz6hRo5g5s0+b609ERp+IcOutOVx9NQwfvprFi9uumQw/AYLBIE8++SSXXnopUIvbPZaamodbnCc23TMkVVM162293hPJyclhzpzZnHyyv5FdMjNjsyrQtxfz58/n//7v/xo+9+7dG3DgdB7KvvvWcOedd6bOuSZce+213HzzzfTo0YOnn36avn1T7VFkNm/eTNeuXRs+5+RcSl7ePehc4Eo2UVVVxemnn87JJ5/crvVqr5t4mTWLmjFjAPh3t26UTZ/O2rVr8fvXU1DwFh9++CE71ndvTHFKthQX84z9NH/VTTfRPVo3xUT72Yr9WRcW5BcuXEhu7t1t1kyXVHjVVE2ADz/8EIAepaU6BIKkc68bO315ov3U/COIPrRAAlKy26r5ru3nFe3pZytTzK+0fZ3sdmddKrxqdkzN2scfl40bN4rfH5TCQpE999xTAPk6Nzeipg6BECNJD/TFxRKyA1JXkFC0lGjbttGXGc02nuEKYrR9b8cd5RAQJ0ghyOomdn5/nbz++uuydu3axPvZihTzEEh3+7jWJUgzGX6qpmo2pxkCeQrkYPt3R4RyVNMgH6apgT5Gkh3oS0GG2l/YpeFfkjHbaxoT+QttahurXRy29SfVjiDfNLHz++vE5TpRAPF4PFZgTaSfzdhtd+xt23W2vzckUDMZfqqmajanOZNtAb07yOUgZ599trhcx4vb/Qc5//zz5fNmfm861k2acLPfzyKsoQoaPW5Ns5TsPPsJ5gNYAwlsBrYAF+fnU1rqJBh8EYCamhpe32GHxPrZiv15zv44PIGaLS5XTdVMoGYVcIP98RtgA3BncTGzZs0iP/8lvN7HmDZtGnsVF8defyKJdgVIZUnHO/pXX3214WpdEX4lTsM2+onHH9/ga6SSm3urbNiwoeHzQzk5DXpBkFc8Hvnwuuva3F4ZsvXe83ik9OGHtz/2tu2uILvZtqluV1VN1YxV0+f7RebNmyf3jhnT8Fs6uQ2a2nQTI8kK9B7PDQJIly5d5Ns770zYcAVx28Vhu/y22+ThLl3kXpA/FRTIpUceKfPnzxe/v7bB7t577204QX/q1UvEGLnE5wu7KPgkP/+zZusOhULi8dwkHs91ctNNN8nYsWOlZ+fO211cjjrqqEZ+vvTSS3L77bfLh9ddJ4D8FdIibV01VTNWTadzn0bn+Gleb5s0dQiEGEneEAj7Ap8CnSgs3JIgzfRI3Q4EZlJVdQ5O5yHU1b0LgDE743KNoLZ2DgAFBYtwOvdupBMK/Ygx+VRWnkYwOH+7eozphchWHI6+hEKLAXA4diUUGgB8ADSe79brfRa3e1vf4kw9nqrZkTT3AP5Lfv77OJ17Y4y3GdvY/WxtSNbJwdvMR4AH2Ept7fZBLZMxxuq/Xh/kXa6TyM9/Ga93NrAQyKWi4nCqq28iEJhFKLSWQOBRysp6UFpaaAf5y4DN5Oe/jc/3PX5/CL//ewoLy/D5PsfnW4PLdSyh0HfAK9QHeY/nxgY/nM792nO3FSUB9ALA4egZNcinDdFu9esL0BtrZKkVwHLgzxFsRmD9epfY5bqwdX+xt1sGzAFyW6oz3ZpuCgtFvv76awHE5/NJIBBIiGYy/IzXbuPGjXLcccfJZZddJsFgcDvb/Pz3xePxRGzrP/roo+Wee+6Jue5gMCg+3ybJz39TNmzYICIiy5cvl3nz5rV6f+KxVU3VTKSmx3OjALLXXnslRDPVvW6CwBUiMgBrIqSLjTEDI9i9LyJ72mUSgDGmJ9bt3j4iMhhwAu03kk8C6devH7fccgtlZWV80rt3VmTqAXTr1o2XzzqLu194AafbvZ2ty3Ugy5YtY+7cuTx2wQVc5fczGfipVy/mjR5tj0sTW91OpxOHoysu12ENWboDFy/m6HHj0voYqaZqRrLNzb2O888/n8WLF3PC3nvzY69e7ZvtGg/RrgDRCtYcGkfI9nf0L0ew7Ql8D3TBGlfnZeDIlupIxzt6EZFv77xTALk5zXvdpFQzzXtKqKZqJlJz06ZNMuqAAwSQv8Si2Qxp0+sGKAHWAv4my0cAvwBLgXnAoLB1fwbKgY3ArGa0L8BqFF5YVFTUqh1tj8zYfJCzwr9Qu7fIdpopzIxNqWaGZTOqpmomQrMvyGGxaDZDWgR6rDyhRcDJEdb5gQL7/bHAN/b7zsBbQHfADTwPjG6prnS9oxdjpCvIvk2/0DTLjE2pZoZlM6qmaiZC81y2Pbs6BmR6NM1mSHUbPcYYN/CsfUc+N0LzT6mIlNvvXwXcxphuwOHAahHZKCK1wFzggFjqTEe29OxJBdbVqxFpmKmXMs1s2x/VVM0YNK+z3xrgHeA8YG63bpG3TwEtBnpjzfrwCLBCRKZEsdnJtsMYM8zW/QWrmWc/Y4zXXj8Sq/dOxrF69Wp6bthANdZflgYyZLLidtPMtv1RTdWMQXM3rxcBQsA6e/EpGzfy4osvbq+bCqLd6tcX4ECsvyRfsK375LHAOGCcbXMJVhfKpcAC4ICw7W8EvsLqXvkE4GmpznRsuhkzZow4HA6ZMW6chIqK0j4zNqWaaZzNqJqq2R6an02aJICceuqpadF0o5mxLVBv63DsjTGFFBS8nTDN7MoSVE3VVM16ysr2IBT6Lx7PjdTUXBeXZmtDsmbGJoBQaBnG+Fs2VBSlQ1NTcyeh0H8BFx7Ptal2B8iyOWOTNTm4iFBaWsvw4b9G1U6nCYhVUzU7qmYwGGTx4sV4PB4GDhzIokWLOPLIYcCvbN7sx+Fo/t42EX4eeeTrvPEG/Otf0xkzxsSsmczJwbPqjj5ZBINWR6Ndd901xZ4oitIct912G8OGDWPIkCG43W72228/SksdlJZ2xul0tsvD0aOPPhqAvfbaK+l1xYoG+nqipDqHQmuorDwVgNGjR6d1SnbaaGbb/qhmxmj+eNttDR/tjoA4HLvhcFijr5500klcfswx1BYXJ83PskmTAFi1alXb9ylRRHtKm8rS7r1umkl1zst7RAB5/PHHs2+4gmRoZkDaumpmr+aHNB5475UJExrsnn766YblB4NMxp73OcF+brbrGNizp1RWVsasmRaZse1Z2j3QR0h1Xg2yZOedxek8VABrZMdsG64gGZoZlLaumtmpGQKZZAfbQW73dnZXsu1CsGOS/HzI1r/66qtj1tTulTFS/zAj7oexDod1uG3eAQ4Ns3e7/0Ag8Nh2dg0YQyd/qFnNqLbZptmMHaFQ42OvmqqZRM2hwI9AZaFsZ/crMAhYD0wAJifBzwKgAnA6D8aYHaiqmkOZ200B1ngwTTVHHGxpthS/otFhule2miYpzWvs17O8XrzeZ8nLmxbRLtr2cdlmm2a27Y9qZqzm77ECfTD45nZ2hViTbADcAUzt1InFixcjUpMwPx/v1o2DDjqIurr3CAafwe120wXIwRpPJhaNhBHtVj+VJdVt9P+y/3Z9N2VK+rR9Z4pmmrTVqqZq/pKXZzfROCQ/f0FEu+UgeYS36+dLfv6HCfXT7w9Jbu4dMn78eBk5aJAAMhCktomtttHHSJuGQLBTmEMgh+fmWoH+u+8iB8ZsGq4gGZrtkGKumqoZi+0zzzzTEMRramoi2m164AF5/fXX5cYbbxQwYkzJttnWkuDnrWecIYDMtNvm6221jT5GEjEEQnX19dTUTAIMfn+Q0lJHVNtYNVtrq5qqqZqJ0LwSuAOv9yXc7uNbsL0CmEJ+/tu4XCOS4qdIBaWlBTgce+DzLd3OrrUhubk2+qzKjE3ENau29iUA/P4qjNFHGIqS+fwZuAORH2KwvRSYQkXFoeTm3oYxPQiFlgH5OJ2Dcbt/3yZPREJUV1v97I3xtUkrzopT31TTtKRy9Mqzzz5bAJkzZ06LtrFqttZWNVVTNdtu6/fXCCBnnHFGDLYhcTp/G9Zm37h89dVXrfZz+vTpDTp77bWXfP/9943sUj7xSMbSiky9qbNnA/DYP/6RMM1MyCjUzFjVzEbN//3vfw2jRy5YsKBFTWMMBQWf8tNPP/Hiiy/y5CWX8MoOO2BvyWMTJrTaz48//hiAiw4/nA82baJXUZFmxraGRlfENvQ8wS5p1ZslUzTTsPeFanZMzcqGXjcI5MjSpUtbrfmjrXOS09lqPydPniyALPB4Itpqr5sYaXSgmmSgTQN5xH7KvV2wa2LbB6QzSKioqLFttmWxJkMzjTIkVbNja36PFZxdruPF769tk+YIW2t+G/x8/fXXBZCbothqoI+RRgeqycS+hJXc3KlSWCiybNky+eabb7aznWzbvQuNv9Bsm8g7GZoZNKGzamaP5myQkSAD7NI97Pfu9T7TJj8rw7Ti8XMzVpfOwkKR/Py3GjTmRalf2+hbg51p9h0wvcmqmpop1NY+y+DBg+nbty/3dW483ffv7Nevu3SJqBmtrlbZZptmtu2PamaE5pVOJ28Cq4AqYCPgA4Z7PLhcJ7bJz/rhCg6Pw8/FQBfA4/EQCMygpuZuAF7dYQeOjrX+RBLtCpDKkqg2+lvd7oarKCAveTxy9YknNlpWX25yu+UWkM9Bfg1b7veH4m9bjMc22zTTsK1WNbNbc9WqVQLIeU5n0vw8AaQAJJSXF5PmFS5XxDgTnDEjav3adBMj4QeqsrJS3E6n9HG55B6skShl5kxZv369OBx9xOHYXVasWCFvvfWW5OfnN/oyPt1pJ+nVpYsAUlDw5faBMZuyWJOhGUeWoGqqZms1Cwq+bfS7vefcc5Pm56Vs66DRnN1EEJdt6/F4ZNmyZeL1Pi1O5yFy9913N1u/ZsbGSHhmbCj0C2Vl3cjNnYLH85dGdttnqtUSCn1BXd2XVFWdi8t1MlBJMPga8BZwaBZlCaqmamaL5iuAlenqdp9FXt6shslGEu1ndfU11NTcDHjJyTmHQOB2wN/Irrb2OSorTwbA5TqKvLwncDi6x73vrQ3JHWb0yvoDJFJJeflAAIzp1uJ2xrhxOofidp8BOAgG5xIMvoYxRcCQ5DmsKEobOAaw2t9ra58BgkmryeO5CpfrGKCSQOAh4BTgSUQqAQgGF1BZeQ4Afv8W8vNfw+HonjR/4iWrhkConxz8kUd+pE+fDXTu3JlffhlN04t89Ml6c/jgg3f57rvvGDp0KAMHDqRzZxPFNlbN1tuqpmqqZnO2DuAFpk9/htNOO41gsBOPPfYYp59+ehL89CPyCitXrmTgwIGEQv8B/kNpKey2226sWrUKYwyffvoZ++zTqVX16+TgcVJSUoLH46G2tjbiX7nmOPDAAzn33HMZNGhQ3NsqitL+nHLKKZx33nlUVlZyww03JK0eYwz9+/dn69at5OX9C2O64Ha7WbVqFTvvvDMffPAB++wTseUk5WRloHc6nZx77rmUl5ezplevtErJ7hCa2bY/qpnWmsYYHn74YU4dNowVK1YkdtLtCPX7fD5ycsbg9/9CIBDg+++/Z91tt3HA2We374Tf8RDtKW0qSyJ63bx21VUCyB3Z0m0xUzSzrOueamaO5ni7S+NjLZ3zzZB0P5tBu1fGSPiBqisqEkB+H37wIXOHFsgUzTRMhVfNjqFZhtWtsUdL53wzJN3PZtDulTHSaHJwh4M9RPgSOAE4BNgDOCxTJ93OFM00nCRaNTuG5lIR9gSGAZ9E02yBpPvZDI3iVyvoMN0rG1FUxAygDnge+AswEri8oACR0Ha20TRaZdeRNbNtf1QzYzSvsN9OjkUjFpLhZ6qIdqtfX4DeWBOmrwCWA3+OYDMC+BVYYpfrwtZ1Ap4BvrI19m+pzkQOUxwAWQGyBOQwh0MAyc29PfPavjNFMw3aalWzY2ruRNjAY82d883QYdvogZ2Bve33PmAlMLCJzQjg5SjbPw6cZ7/PATq1VGfCZphqkmoceuIJAcTp3F/8/lDjjduYkh3xy+yomlHsVFM1k6l56rBhAsiG+vbxOANtu+17FNKqjd4Y8wJwr4i8EbZsBDBBRI5vYusHlgK7ShwVJXNy8PLykdTVvQV0xu0+ltzce3A4ukS1z5x0cNVUzY6tGQhMo6rqQvLz38LlOjRt/WzJLs6Q3EDC2uiNMSXAXoQ96whjf2PMUmPMPGPMIHvZrlgjhj5mjFlsjJlujMmPon2BMWahMWbhxo0b43GrgVgOUH7+S8DZwBZqa2dRVtaVUOh/rapPUZRUM5Xy8sOpqXmQUKh+8u9m7vQ6KtFu9ZsWrKO3CDg5wjo/UGC/Pxb4xn6/D9YAFPvan+8G/t5SXe0xObjfXycHHnigAPLss882axvL3642/T1UTdVUTQkGg+LxXCdu95lyyCGHyPLlyyUUCsmKFSvkl19+aWQ7b948gcgTeAcCgaT6mSzNlE88YoxxA88Cs0RkboSLRamIlNvvXwXcxhpNbB2wTkTq/wE8A+wdz4WoTTQ7CbCDyZOt5/NlZWVpn/2XUZrZtj+qGbNmqLiYRcawtagoLs1n//xnXC4XNTWTqK19knfffZdBgwbhcDgYMGAA+3btSqi4GGbNoq6ujj/96U8A5ORcwAMPPEDfvn3p5vPxfPfuuD2e9MxOTSXRrgCy7W7dADOAqc3Y7AQN7f3DgLVhn98H+tvvbwAmt1Rne00OXlpaKoCcNHRo+vdmyRTNNOh9oZqp0xzX5O56/sSJLWpOC7PPzb1LCgtF3nvvPRlz0EFyttMpfex1HpAD7Z5zgHg8/9e6c74ZsvWOPpZAf6B9YL9gW/fJY4FxwDjb5hKsrpdLgQXAAWHb7wkstLd/HujcUp3JmBy8oTTJagPklKYnSBTbrMtiTYZmGmRIqmZqNG8KC9i72K8D3e4GTb+/TqqrqyVUVNSg85Vt5wT5vmfPiHWHQMaHaRc5nTJx4sTWn/PN0GEDfSpKMiYHbyhNJvbdd999ZSDIphhss24i72RoJnmS6ABIXYI1k+FnR9OstoNwd5AfbJvd7GWXX365wLZZ3EpAZtH4zv+zGPzcBLK6rfveAtka6LM6MzaW5YFAgC+BblhtVBMIm77Atg2F1vPCCy+wcuedY68rgzIKE6qZxLoXYCViOIHTEqSZDD+zRfOdd95hSqdOkafzaGL7ry5WF+XLsBJvAK63X6dMmQJU4HBYkwH9DxgVtu05WL02WvKzK1ASzf9MyE5NJdGuAKks7dVGLyKycOFCuW/MGPmtwyGesDuMC10uKZs+XXJyJjQs6+T1ymceT4ua8dSvbfSxac7NyWl0BwjIwbvvLjfccIP4/QHx+0Py9ttvy9J//EOCeXkp8zNbNL3eVxqOsw9kDMjNYbahJ56QvLw54nT+tsHOC1IdpvlDXp64nE7xer3i861vqP+n++6Tv7pcMg2kFKS2Pfe9BbL1jj7lQT1SSVZmbCzZoU+D7N0kqDgcQ+SWW25p+Dy7W7fsyWJNhmaCMyS//PJLAasd9xmQ93bcsUnQd4gxXRs+711SIn/z+2U9pHUmZ7pq+v1lAn5xOp1y1FFHbXeBjVTGjh0rq6ZM2U6zrq5OQqFQeu17M2RroM+q0StjyYyF2DLVqqsnEQp9SW3tIOBaCgshGPyEysrjEPkFv7+c8NyvTMu+yxTNgoLllJcPBsDrfQW3+1gARIKIbKSsbADWMEvgcAxApAKRtQ06fn8txriS7memaPr9Qij0HQ7Hro1mUBMJEgqtpLp6IsHg60AtHs/fyM29iVDoR4zpQlXVH6itfdLeIherAW0qhYXRM8tb62emHM9MyYzNqjljE3nNys29Dth28AFcrn3xeP5OdfVFlJZaV5OcnMvIy7s7cRV3UESEQOBBHI4iXK7hWLl5I6ioOAGAvLzpDUEewBgXxuyM3/8/SkvnAmfh8+UhUkpNze0EAjMRWUNZWQm5uf8gJ+fclOxXulFZeQrB4HM4HH1xu0/H4diFQOBB6urCb6x6AKfg8fwfAA6H1eru9c4B5jRYhf82lDQn2q1+Kkt7ZMa25W/XxIkTxe12CyC77LJLQjST4We6a1ZXV8ukSZNkzJgxsuuuu0ZtFvjd734nwWAw7vrvu+++Bo0ePXpIQcGitNn3VGn6fL6Ix7hfv35y+eWXyyeffJIWfnZEzWQ23WTVHX17ceutt3LLLbcwbNgwVq9enWp3MhIRYf/992fx4sUA5OTksP/++9OvXz9+/vlnXn/9E0S20Lt3b5577rlW1XHRRRex3377MX78eN5//31gKHl5jwMd9+6+rKwMgJNPPpmLL76YkpISdtxxR/LzIw5BpWQJ2du9MlnYqdvG6WT411/zyy+/8NFHHzVrm00p7omwE9lC3759Wbx4MaOHDydUVERNbS0f/fAD/zriCObNm4ffvxm/P8jatWvbtD97n3wy7114IePHjwecVFWNYezYsak/RinSXLBgAQBz587lv6eeyq59+pA/aJAOF5DtRLvVT2VJ26abJl24vrT/9p65336RewCke1fIZGjG0CXO7T5PADlk991lSzt2hfT51orDMdhqcuveXVbn5rZZs7BQxO8PyV133SXHDhkiZzqdsjXNu1f+FNakFWpJsxkypUkkUzS1e2WMJD3QN0mzDoW1cfp8qxprZsJwBcnQjCG93uU6RVwuV0pS9n2+DQ3f2ZVt1AwVFYnbfU6DXn3pQ1j2biv9TMa+FxZaXVW9xgggw2PRbIZMCaCZoqmZsenC2rWNPgpQ37JZXr6/deWMYlvPgjVrqKn5Z4t2EZfHaptKzRjs6ureo6SkJCX743B059ZbbwXgjWbsWlq2FXCsXUtt7ROA9TxgC/Bb4FuaPAVIk++ytvYVBg8eTKV9nv4n1u2VjEcDfTw0Sad2AOXAmV4vIhsoKyuhtLSUSZMmMdHn43ngS6AWeB04GNgfqK6+jEDg6Yia0eqKyzaNh0AQqUZkM8FgMGX7c8YZZwDW6Hyt0fwV6Gy/d7mOwu+v5b777qNTcTEfA4cDs4AtbfQzLrsYbKuqziUUCvF89+5UYfWEj2l7JePRQB8PN98MXm/jZV4vsx96CKfzUETWUlhYyPXXX8/tpaX8HhiENT7L0VjjNddTVXUGP//8c1RNbr455vq3s41gV56Xx/vnnEMg8DSBwCNUVFS0WTMeu0ePPJLi4mJKS/OAOm644YbE1x2jbUlJCWMOOgiAf7VCc6T9undJCfn5rzUkZHHzzTi9XvrY6ze10c+47OKwPXHKFHJj1VSyg2htOqksadtGLxI1zdrvrxO3+3yr7XP4cPn888/l1SuvlBsKC+VakBsLC2X1XXdJKBQSl+uEhvbcurq6pA9X8N+dd5bO+dtGDwTE6XS2STMeu/pJ2QFxOPaU3Ny7ElZ3a/3cvHmzANLP5YpL87uwYxgIBLazq3z0UXHXf7dFRWkzBILP970AcvDBB8evGYVMafvOFE0dAiFGEjkEQry221LMSwHvtru8qLbHA6/g8VyNx3NTo3T0RPopUkdZWTEi63G5TiQYPAO4EVhJXt5T5OScHrdmvHbB4HtUVBxCTs6VBAK3J0QzEX6Wle1OKPQ1fv8WjOkUk2ZZ2d6EQospKFiE07n3dnYiQUpL/UAVDsdACgo+xhh/m/xMhO2vv04GrsLrfRa3++QEaSbDT9VsbUhO2OTg6U46XLOM8Tcb5LfxNMbsQk3NPygvH0QwGGm+9bZTU3MTIuvxeG4hP/8FrInRPwc8VFWNJhCY04JC6xCpIxT6joqKI6moOAQAj+fSpNTVWqzZL31YUx7HSrm97eaIa41xkZ//H4zZgVDoy6R9r/GzEbCeKSgdkGi3+qksad10k0DNmpoamTBh2zDI55xzjrzyyitWc04C/LzpppsEkOLi4oYhBOrt3njjjYZ6y8vL27RPkexGjx7doO/1euXll19us2ZbbcPt6qeRPP7440VEJBQKyWWXXSaAjBgxQvz+UETNv//97wLIoYce2mzdq1atatj/N998s9V+Jso2J+cvAkhtbW3CNNPld5Qtmtq9MkvJyclh8uTJvH/ttQA88cQTHHfccTidTma1MUvyyUsu4ZprrmEosDQUwvnkk43WH3744TzwwAMAvPLXvyY0Q/PD669n5syZALy2ww5UTJvGcccdF/U4pII5c6x/Mqd//DGVxnBQXh733HMPYE24UVl5NMFgWMazve/XXHcdA91uvlm6NLKwbbdrnz7sm5MDwMiRI6mtrY1q2x6ZsSJWu4DT6YxfU8l8ol0BUlk6yh29iDRkM1aDLMOaG7MTiNPhEL+/evsHbc1maAbl9ddflxP22ksA6QJSFcGuXnPLli0Nd50VUTSb+hmt7nC7x+3x/Oc0sWuX4xmj3aEDBwogU0AK7WMwwuGQDfffL8OHD284Ltdff73UzZjRaN+Ptdd9NmlSs8eoLuzB7bx581p9PBORGetw/EYA659dPJoJOO6qmfo7+pQH9UilQwX6CNmMz9rBwZgu4veH/dVuJvOxoOC/AoUNgWUIyMyWsi5F5I7OnQWQx8KC04ZItjFmci7beecGH1ancdblfV26bOsJFH5Rsv30+VY1Op49QbbYNovsZbu5XC0eo8m27R577CE+34a4j2eiMmNdrpPEuq+LUzMBx101Ux/otekm1UTIRvw9cCEgsplAYEpE2++Bm4C/ANPWrKG8/DfAr/zhD39gLVYyUPi8nNHq+tOWLQBcC4wACoEdgGfXrIm67TfACVhz7J62Zg21tW8QDH7GW2+9xZAffwRgIlDcQt2p5KItW3gXOBR4AjizfkVDBu2u+P2rGTRoEADrsZKkjgJKbdOcYJPZVCPs4wTgOOCLL76grGwHKiqOYt26de2eGVvf8+eTTz6JT1PJCjTQp5oI2YgGeLC4GIdjCNXVE5k6dSqbN2/mlW7dOAvIA4qwgvNUrIsCQF7eNB577DF6FxdvpxmtLm9xMfcD64B3qe9TApc5nYjUWH/7wrb9DugHvGzbPQNUVh5JRcUwRo4cSR3wb+BWez+aqzulFBVxMPAWVj+k8OX1GNOZZcuWUVtURFd72XysiwOAx+3eTjMSE3bcEb/fCrTB4Hxuuummds+Mzc21urXecccdOpF2B0QDfappJpsxN9fKVPzLX/5C165dOX7jRp7Emv/nj8AC4H+5uUwdPZr8/AXk5Jzfomak+v/k9TIbKzjXAZe4XPxQV0dpaS6lpQ6MMZg1azjc4WA3e7PnAfF6WXH77Xi9L+H1PseMGTP46PrrOTUTsi7jOEauf/yD/+XlsQoIzzo4dOTIxoZRNEfceSe//vorhYUC5DNr1iw+++MfeSwnh4qW6k9YZqz176Oqqio+TSU7iNamk8rSodroRZrNkMzPf0PGjx8vBx10kBx22GHy7Z13Jn1y8LoZM+yumfnicAyQrl27Sk7YhOn71LfnJiKTM9ZjlADbNvlp2/5I45Eqfb5f4tLMy3u40fbdHQ5rqOAkZsZu3bpVHI7dBJCXXnopfs0opN3vKMM1NTM2RtIhMzabNUOhDYhswOEY1CiTN938TLZmVdVlBAL1I5B+B+wSl6bX+xp1de9SU3Or/Xk+bvcRCfeznurqf1BT8zfgEQoL/5gQzWT4qZrWa2tDsmbGKgnB4dgBp3NwxOEaOhI5OReGfdol7u3d7qPJzb2FgoJvAaiuvjhBnm2PiBAITMV6xN5ykFeyk6yaM7ZHD6u8807zdp06Wa9bt7asGautanYkzUF8+ukn+Hw+9t+/LZq78ZvfDGbZsmXce+9MRo8enWA/4b777ueSSzbi8fyV3Nx0PZ6qCTBiRMtarSWr7ugVpb0YNmwYAwYMaLPOlClW99kxY8aQjGbUl1+2+kd5PH9PuLaSOWigV5RYSMaQAbNmccT551MEdAHM7Nltr7uJ7eKPP2b//fePcaA9JVtp8ds3xvQGZgA7ASFgmojc3cRmBPACsNpeNFdEJoWtdwILgfUicnxCPFeU9mLWLLjgAqistD6vWWN9BiKkpcWtuRk4LhSKrNlc3aOa1D1rFqHzz2djVRV/B75Zs4afgT3Ly1E6NrHc0QeBK0RkALAfcLExZmAEu/dFZE+7TGqy7s/Aijb6qiip4W9/2xZo66mstJYnQDMAPAWsiaQZY92hUIiu556Ls6qKnYD7sJK7BgBTNm1qvZ9KVtBioBeRH0Xkc/t9GVbA7hlrBcaYXlhZ4NNb66SipJRkDBkQtm192tWbkTRjrPuNN95gcygEwOVYmct1WHMWD/zpp9b7qWQFcbXRG2NKgL2ASLMp7G+MWWqMmWeMGRS2fCpwFVazT3PaFxhjFhpjFm7cuDEetxQluSRjyICwbesHPN4/kmaMdV9rD3W9HLgT687KEcVW6XjEHOiNMQXAs8B4ESltsvpzoFhEhgD/xMqQxxhzPLBBRBa1pC8i00RkHxHZp3v37rG6pSjJJxlDBtiaAvyAFZgHtGEIhKA9wNpAHdpAiUBMgd4Y48YK8rNEZG7T9SJSKta8bIjIq4DbGNMNGA6caIz5H/AkcJgxZmainFeUdmHUKJg2DYqLwRjrddq07R+GtkbTvtvew++PrBlj3T6fj9122y3xfipZQYtDIBgrDfJxYLOIjI9isxPws4iIMWYY1qCGxRImbvfMmRBLrxsdAkE1O5Lmr7+6cLkOIz9/fis0y4GzgJdxuY4iP/+1pPmpmu2jmaohEIYD52DdjS+xy7HGmHHGmHG2zanAMmPMUuAe4ExJRvZHC+gQCEom4nD0Jxh8g1Do5+3WBYNvUla2G7/+6qa6+paG5SJCZeU5WJObvwy4yct7uN18VjKLrBrUrD6FWIdAUM1M0pw3bx7HHnssRUVFfPvtt3Tvbo1z/+mnK+nfv38j2xdffJGioiJuvvlm/v3vf2PMLuTl3UZFxWlJ91M1k6sZa/yKRocZ1ExRMgo7i/WY445jUmEha9eu5ZprrgGgru5L9thjD1wuFy9efjnf9OgBwIknnsiee+7Jv//9b3bYYQd8vq9xu0/bTlMn/VbC0bxoRUkFTTJer/31V64D7rnrLoLyJcGgNUbNbWeeyQkPPgiVlfwXeBDYJyeHThddxMhJk+jd2x1VMyEZvEpWoHf0ipIKImS89gGqa2sJBl/G4RjI/PnzuerjjxvsBgP3An8IBPjdc8/h8/la1GxzBq+SFWigV5RUECHj9Wj7NT//U3y+5RxxxBEJmxxc6dhooFeUVBAhW/WfQKCoCJfrt83aRV2uk34rUdBAryipIErGq/sf/4jJLv7JwZWOjAZ6RUkFsWbbxpOVm4wMXiUryLpeN0uWtDwlV/3w3LFM3RWrrWqqZvyao6BkFJTYHx+2Sqx2bdOMx8+22almbLZLlsCee7as1xqyKtCffXZsdrvuat3wJNJWNVVTNVWzLbZ77hl7DIuXrAr0F1wQ1m1YURRFAbSNXlEUJevRQK8oipLlaKBXFEXJcjTQK4qiZDka6BVFUbIcDfSKoihZjgZ6RVGULEcDvaIoSpajgV5RFCXL0UCvKIqS5WigVxRFyXI00CuKomQ5GugVRVGyHA30iqIoWY4GekVRlCxHA72iKEqWo4FeURQly9FAryiKkuUYEUm1D9thjNkIrEmxG92ATSn2IRrp6pv6FR/qV3ykq1+QHr4Vi0j3SCvSMtCnA8aYhSKyT6r9iES6+qZ+xYf6FR/p6hekt2+gTTeKoihZjwZ6RVGULEcDfXSmpdqBZkhX39Sv+FC/4iNd/YL09k3b6BVFUbIdvaNXFEXJcjTQK4qiZDkdLtAbY3KNMZ8aY5YaY5YbY26MYLO7MeZjY0yNMWZC2PLexpi3jTEr7G3/nA5+ha13GmMWG2NeThe/jDGdjDHPGGO+so/b/mni11/s7ZYZY+YYY3IT4Vccvo0yxnxhl4+MMUPC1h1tjPnaGPOtMeav6eBXGpz7UY+XvT5V535z32NSzv1WISIdqgAGKLDfu4FPgP2a2OwA/Ba4GZgQtnxnYG/7vQ9YCQxMtV9h6y8HZgMvp8Pxstc9Dpxnv88BOqXaL6AnsBrIsz8/DfyhnY/ZAUBn+/0xwCf2eyewCtjVPl5L2/kci+ZXqs/9iH6lwbkf1a9knfutKR3ujl4syu2PbrtIE5sNIvIZUNtk+Y8i8rn9vgxYgRU0UuoXgDGmF3AcMD0R/iTCL2OMHzgYeMS2C4jI1lT7ZeMC8owxLsAL/JAIv+Lw7SMR2WJ/XAD0st8PA74Vke9EJAA8CZyUar/S4NyPdrxSfe5H9CuZ535r6HCBHhr+5i0BNgBviMgnrdAoAfbCusqng19TgauAUKL8SYBfuwIbgcfsv9XTjTH5qfZLRNYDdwBrgR+BX0VkfqL8aoVvY4F59vuewPdh69aRoIDaRr/CNUpI7bnf1K+ppMe5H+5XUs/9eOmQgV5E6kRkT6yr7zBjzOB4tjfGFADPAuNFpDTVfhljjgc2iMiiRPmSCL+w7pr3Bh4Qkb2ACiBhbc5tOF6dse6SdwF6APnGmNGJ8ise34wxh2IFiIn1iyLJpYFf9ctTeu439Stdzv0Ixyup5368dMhAX4/9V+od4OhYtzHGuLFO9FkiMjdN/BoOnGiM+R/WX/3DjDEz08CvdcC6sLugZ7BO/lT7dTiwWkQ2ikgtMBerrTXhNOebMWYPrOaGk0TkF3vxOqB3mFkvEtis1Aa/Un7uR/Er5ed+M99j0s/9WOlwgd4Y090Y08l+n4f1o/8qxm0NVpvbChGZki5+icj/iUgvESkBzgTeEpGE3KG20a+fgO+NMf3tRSOBL1PtF1aTzX7GGK/9nY7EanNOCLH4ZowpwrrAnCMiK8NWfQb0NcbsYozJwfo+X0y1X6k+96P5lepzvxm/knbutwpJ0VPgVBVgD2Ax8AWwDLjOXj4OGGe/3wnrilwKbLXf+4EDsf5GfwEsscuxqfaric4IEtvzoE1+AXsCC+3tn8fuoZAGft2I9aNdBjwBeNr5mE0HtoSdRwvDtj8Wq1fLKuBv6eBXGpz7UY9Xis/95r7HpJz7rSk6BIKiKEqW0+GabhRFUToaGugVRVGyHA30iqIoWY4GekVRlCxHA72iKEqKMcY8aozZYIxZFqP96caYL+3B1ma3aK+9bhRFUVKLMeZgoByYISLNZngbY/piDcR3mIhsMcbsICIbmttG7+gVRVFSjIi8B2wOX2aM2c0Y85oxZpEx5n1jzO72qvOB+8QeTK2lIA8a6BVFUdKVacClIjIUmADcby/vB/QzxnxojFlgjGlx6A9XEp1UFEVRWoE9eNwBwL+t0ScA8NivLqAvViZwL+B9Y8xgaWYYZA30iqIo6YcD2CrWyJlNWQcsEGtAvtXGmK+xAv9nzYkpiqIoaYRYQ0CvNsacBtagcmbbNIXPA4fay7thNeV815yeBnpFUZQUY4yZA3wM9DfGrDPGjAVGAWONMUuB5Wybaex14BdjzJfA28CVEjacdER97V6pKIqS3egdvaIoSpajgV5RFCXL0UCvKIqS5WigVxRFyXI00CuKomQ5GugVRVGyHA30iqIoWc7/BxxsZGKRfcBXAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "t2m = xarray2df(t2m.resample(time='1d').mean(), 'T','t2m')\n",
    "\n",
    "check_data_gap(t2m)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "f250450e",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(14610, 75)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t2m.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "90144c35",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n of pixels75\n",
      "Date start: 1979-01-01, date end: 2018-12-31\n",
      "No missing dates\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEQCAYAAAC+z7+sAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA7hklEQVR4nO3deXwV5fX48c/JviObIlsAEYorKCCCxaUuVEHRWkvF3YpQF1BatUq1rlSlihbFxuJX/VWsVRQUAbGIFhSogKxSUXYVZZElISHbPb8/7iWGkOUuc+/MvTnv12teuZmZ58mZm+Rk8szMeURVMcYYk7iS3A7AGGNMdFmiN8aYBGeJ3hhjEpwlemOMSXCW6I0xJsFZojfGmATn2UQvIi+IyDYRWRXk/peJyOcislpEJkc7PmOMiRfi1fvoRaQ/UAS8rKrHNbDv0cC/gLNUdZeIHK6q22IRpzHGeJ1nz+hV9T/AD9XXichRIjJLRJaIyDwR+Ulg0w3AM6q6K9DWkrwxxgR4NtHXoQC4RVVPBn4HPBtY3wXoIiIfi8hCERngWoTGGOMxKW4HECwRyQH6Aq+LyIHV6YGPKcDRwBlAW2CeiBynqrtjHKYxxnhO3CR6/P997FbV7rVs+xpYqKrlwAYR+QJ/4v80hvEZY4wnxc3QjaruxZ/EfwkgficGNk8Fzgysb4F/KGe9G3EaY4zXeDbRi8irwAKgq4h8LSLXA0OB60VkObAauCiw+3vAThH5HJgL/F5Vd7oRtzHGeI1nb680xhjjDM+e0RtjjHGGJy/GtmjRQjt06OB2GMYYEzeWLFmyQ1Vb1rbNk4m+Q4cOLF682O0wjDEmbojIprq22dCNMcYkOEv0xhiT4CzRG2NMgrNEb4wxCc4SvTHGJDhL9MYYk+As0RtjTILz5H30xlkFBTA5gskVP/rI//H0093rIxFiiKdjKCn5irKy78nN7YXPt5+UlDzHY7DqK7Fjib4RmDwZli2D7t3djsR4mc9XxrffPsu6dbcdsi03tzfHHz+T1NRmLkRmIuXJomY9e/ZUezLWOWec4f/44YfhtT/sMP/H3bvDjyHSPhIhBq8dQ2lpKc8//zzz589n8eLFrFu3rmq/7t2706FDBzp37sycOXP47LPPGDRoEG+//XbEMUT682hqJyJLVLVnbdvsjN6YRka1jAULlvCrX/2KLVu2HLTtxhtv5KmnniI9Pf2g9YMGDeLdd9+ltLSUHyd2M/HCEr0xCUhV+frrr/noo49o3rw5/fr147XXXqO8vA3FxRfQt69/vx49ejB37lxSU1PJzMyk2jSdBxk0aBDTp0/ngQceAB6O3YEYR1iiNybB7Nq1ixNPPPGQs/Xqfv3rXzNq1Ci6d+9OWlpag31efPHF3HjjjTzyyCOkpwsZGQ85GbKJMru90hiveOUV6NABkpL8H195JeT2mp/Ppc2asWXLFi486SQmT57MXXfdxe23386f/vQnMjMnk509l8mTJ9O7d+9Dk3wdMbRs2ZLCwkIASksfps5re5Eeg4kOVfXccvLJJ6txzumn+5dwNWniXyIRaR+JEEO97f/xD9WsLFX/XYf+JSvLvz6YPgLtnwMF9K5Q2wcZw6hRoxRQkbZhH0OkP4+mdsBirSOn2l03jYDddeONGOpt36EDbKqlnHh+Pmzc2GAfO9q1Y8jXXzMHyAG+AfJCaB9sDKWlpWRkZJCU1I3Kys/DOga76yY66rvrpsGhGxFpJyJzRWSNiKwWkZG17HOGiOwRkWWB5d4a25NF5DMRmR7+YRiTwDZvDm19NcuWLaNlIMn3AD4jkOSDbB9KDOeddx4AqamXhtXeuCOYMfoKYLSqdgP6ADeJyDG17DdPVbsHlgdqbBsJrIkwVmMSV/v2oa0PWLBgAT169ADgQWAp0DmE9qHEMGfOHD4KPNaann5PyO2NexpM9Kq6VVWXBl4X4k/YbYL9AiLSFrgA+Hu4QRqT8B5+GLKyDl6XleVfX4cpU6bQN3Cf5Mzf/54xIbYPNYYPA2MtOTmrEKnlXvowjsHERkh33YhIB/z/HS6qZfOpIrJcRGaKyLHV1o8H7gB8DfQ9TEQWi8ji7du3hxKWMfFv6FB/UaL8fBDxfywo8K+vxT333MOll/qHT15++WUGPPZYSO3DiaF7oIZGZeV/HTkGEztB30cvIjnAFGCUqu6tsXkpkK+qRSJyPjAVOFpEBgLbVHWJiJxRX/+qWgAUgP9ibNBHYEyiGDq01qQ4ffp0jjnmGDp16oTqLgoLO/PIIz/QqVMn3nrrLU444YR62zsRA0BmZmbgVT3nbE7EYBwXVKIXkVT8Sf4VVX2z5vbqiV9VZ4jIsyLSAugHXBhI/hlAnoj8Q1WvcCZ8YxLbHXfcweOPP37I+mbNmrFq1apqyTf6HnnkEQCSk3vH7GsaZwRz140Ak4A1qvpEHfu0CuyHiPQO9LtTVf+gqm1VtQMwBPjAkrwxwdm4cSOPP/44WVlZdOvWjcMOO4yUlHPJzCxg+/btMU3y4H9oyq88pl/XRC6YM/p+wJXAShFZFlh3N9AeQFWfAy4FRohIBVACDFEv3qBvTBy56KKLAJg9ezb9+vUDfrwPPsmFZ9qvvvpqpk6dSlHRyeTkLMV/uc7EA3tgqhHo2hXWrYOcnPDa79nj/9ikSfgxRNpHIsQQSnuf7wcKC5uTknIp2dmvuxJDbcrL36a4+CLgCPLyvkEkOewYPJh64lpED0yZ+Ge/UPFFtZzi4l8CkJZ2tcvRHCw19ULgZuB7Kivnuh2OCZJVr2wEWrf2L1YCwd0Ygm3/+OPjueOODzjvvPOYOfMCqlcOdvsYVJXU1AVUVsLKlUfRsWPofRwogWBix87ojfGYA+WFU1NTWbFiBfv27eOjjz6qu2JkDM2cOZPKyiWkp99Jx3CyvHGFJXpjPGbgwIGkpaUxffp0unfvTk5ODmeccQYvvPCC26Fx4GHG8vJ3XI7EhMISvTEec+6551JaWsr8+fMZNGgQvXv771u/8cYbUS1zNbaLL74YSMHn+5yvvvrK1VhM8CzRG+MVNSbt6LdxI2+//TaLFi1i7NixVFZWsn9/LcXE6mgf1qQfDfSRl5dHdrb/IuyvfvWr6MRgHGeJ3hgveOUVGDbMX89d1f9x2LCqRHnnnXcGdiwJq70TMRwgcgRA1YxTjsZgosISvTFecM89UFx88LriYv964MfnSuqY37WB9k7EAFBSUsK+fecA8PTTTzsfg4kKS/TGeEEDk3bsOfCUUZjtnYgBoKCgANVNpKXdyoABA5yPwUSFJXpjvKCBSTsOXJCts0SwE5N+BNHHihUrAEhNHRydGExUWKI3xgsamLTjwK2VSUm1Te7WcHsnYgD4wx/+AMD+/TUnkXMoBhMVluiN8YIGJu2YMWMG2dnZZGUVhNXeiRgAOnfuTErKOVRWfkhJSUnI7Y07LNEb4xVDh8LGjeDzwcaNbP7pTxk2bBgiwvvvv88VVzRQ4btG+7ASbBB9JCefBsD3338fnRiM4yzRG+NBq1atIj8/n+effx6AQYMGMWHCBJej8hPx16X/+uuvXY7EBMsSvTEeU1BQwEknnQTAX//6V1SVt99+m5QUr9QgrACgsrLS5ThMsLzyk2OiaOtWfz36A5ULQ3Xgzr5w2zvRRyLEUFv7ioqPqKj4NyLNSE29mtLSBygrewqAzMxJjBlzHWPGRDeGUO3f7y95PmDAZ6Snnx52DCZ2LNE3Ah4oemhqUC2jsLA1qjur1u3ffzsASUmdyMn5DJE8t8JrwLEA+Hwb3Q3DBM0SfSNg9ei9EUP19jfffDvPPLOT008/nSeffJKFCxdy//3307p1a+bMmUPTprUnebePASAvr5zCQrj66mIK6rgJqD5Wjz72LNEb44IlS5YAcPPNN9OjRw969OjBiBEjXI4qOCUl1wJw5plnuhyJCZZdjDXGBcOHDwdg06ZNLkcSmrVr11JR8Q6pqZfx61//2u1wTJAs0RvjgnXr1gHQt29flyMJzaRJkwBITraz+XjSYKIXkXYiMldE1ojIahEZWcs+Z4jIHhFZFljuDazPEJH/isjyQNv7o3EQxrguxDrsH330EfBjDZtw+og0hnD6OPXUUwOvSqMXg3Geqta7AEcCJwVe5wJrgWNq7HMGML2WtgLkBF6nAouAPg19zZNPPlmNc04/3b+Eq0kT/xKJSPvwdAz/+IdqVpaq/wYn/5KV5V9fo3129qf605/+VAH1//qF3kekMdT5HgTRx44dOxTQ9PSHwo4h0p9HUztgsdaRUxs8o1fVraq6NPC6EFgDtAnyj4iqalHg09TAYjf7mcQSZB121X3s29efefPmAXDdddeF3EekMUTaR7NmzQKvajmjt3r0nhXSGL2IdAB64D8zr+nUwBDNTBE5tlqbZBFZBmwD3lfV2toiIsNEZLGILD4wAbExcaGeOux//OMfGThwICJCYWEXoITJkyejqlXj3Q31EWkMQQuiDxEB0oFa5q61evSeFXSiF5EcYAowSlX31ti8FMhX1ROBvwJTD2xQ1UpV7Q60BXqLyHG19a+qBaraU1V7tmzZMrSjMMZNtdRb3wQkqfLQQw/x7rvvAiCSQUrKJQwZMiSoPupd73T7kPpIorLy8+jEYKIiqEQvIqn4k/wrqvpmze2quvfAEI2qzgBSRaRFjX12Ax8CA2q2Nyau1ajDrsC5Iihw/PHHU1pais/nIzd3HdnZUwJnxfX3AYRWyz1G9ej9Sqjx6+1cDCYqgrnrRoBJwBpVfaKOfVoF9kNEegf63SkiLUXksMD6TOBs4H8OxW6MN9Sowz61ZUvWBs7mV6xYQVpaWu3JvZ4+Qq7lHqN69P/+978BSEqq5TKd1aP3LNEGCqGIyGnAPGAl4AusvhtoD6Cqz4nIzcAI/GXtSoDbVfUTETkBeAlIxp/8/6WqtUxNc7CePXvqj5Mhm0gdeOTcSiDEJobf/va3TJw4keLiYjIzMx37+k70EUn7ESNG8NxzzwGQk/MlhYWdw4oh0p9HUzsRWaKqPWvb1mAJBFWdj/82yfr2mQAcUixbVVfgv3hrTKPxxRdfALB///6DEn0827dvX1WSz8wsIDk5vCRv3GFPxhrjoOeff54PPvgAgIb+W44Hqsrs2bMZE6iVPHbsWNLSbnA5KhMqK2pmjEM+/vhjhg0bBsDChQur3XMevwYOHMiMGTOqPu/Tp4+L0ZhwWaJvBGzikdjEUFq6HICsrJmcd94pjn99J/oItf3evf4qmzk5q0lKasfgwbmOxWBix4ZuGoEEGEHwvP37H2T//puAVJKTe7kdjmNUd5KWNpLk5GMQyXU7HBMmO6NvBGzikejHEKjjx3//+zG9ejWPytd3oo9Q2q9cuZITTqjgwQdbc8cdzsVgE4/Enp3RGxOhF154AYBzzjmHXr0S52z+1VdfBeD000OfF9Z4i53RGwPs3r2bmTNn0rlzZ3r27NngA07dunUjKyuLM888k7/85S/k5uYyfvz42AQbA3PnzmXs2LF07tyZU0459HqDiS+W6E2j9/HHH3PaaadVfd6+fXuWLVtG06ZND9nX59vEUUedxfr16wFYunQprVq1YvXq1XF9l42qsnHjRvbs2cOf//xnXnvtNQBGjRrlbmDGETZ0Y+JfBJNdvP/++5x22mkkJyVxZno6AJs3b6ZZs2bMnz+fefPmMXHiRM4++2yKis6ksLAD69evJzU1lTfeeINx48bxySef0GzmzKhP+hHN9o8++iidOnWiR48evPbaa6QBn7duzU2R3CJkvKOuQvVuLjbxiLMSeuKRICe7qNnHnj179L777quaAGRhenpV+0cD62pbUlIu0rfeeiusGCKd9MOJ96G29lOmTFFAjzzsMH0zLU3ngJaGexxBsIlHooN6Jh5xPanXtliid1ZCJ/r8/IOT24ElP7/OPkpLSzUnJ0cBbd68uc46/PBD2i8FvS03V//4xz/qkiVLdM+ePRHHUO97EGkfYba/4447qv6ILWrVKvLjCIIl+uioL9HbGL2Jb2FMdnHTTTdRVFTE6NGjGTdunH+oo4YeQI+iInigwRp8MZv0w8n227Zto3Xr1lRWVtKiRQs+/PBDjj3++MhiMJ5lY/QmvoUx2cWiRf5JzsaNGxd2H5HG4HgfIbZ/7LHHqKyspGPHjmzYsIFjjz3WJg5JYJboTXwLcbKLysovWblyJWdUf2onrib9cKZ9amoqACNHjiQnJ8eZGIx31TWm4+ZiY/TOSugxelX/xcL8fFUR/8daLsQuX778oIuqK1asCLmPSGNo8D2IYQzdunVTQL/77jvnj6MBNkYfHdQzRt/gxCNusIlHnGUTj0B2djbFxcUkJf2EJ54YzsiRI2Meg1feR9Uy9u5Np2/fvnz88ccxj8EmHomOiCYeMSbeXXnllRQXF5Oefi8ZGfcTRo5PKEVF/iddr732WpcjMbFiid4krKKiIk499VRWrVoFQHr6fS5H5L7y8nfw+ZaRk5PDdddd53Y4JkYs0TcCjbEefVnZa5SUDAd2k5R0ItnZs9m7NymmMTjd3ok+iotnASCynmbNwrsXw+rRxx+766YR8OBlmKgqKRlNSckQYDdpabeRm7uMpKTD3Q7LI7YAIFJ7KWWTmOyMvhHwYj369957j6KiIn7xi184HkPHjm/y3XcZbNmyhRYtWoTVR6QxRKO9E32IfAjAnj3hn+NZPfr40+B3W0TaichcEVkjIqtF5JBLWSJyhojsEZFlgeXeYNuaxmnAgAFceumlB40TT506lcGDB3PllVeyYMGCsPvesWMHvXr1OijJG/D5fEAhAA8//DDbtm1zNyATM8Gc0VcAo1V1qfjnElsiIu+r6uc19punqgPDbGsakbKysqrX//d//8dLL71Eu3bt2LRpU9X6L774gg8++ODHh3lC7L+iosKRWBNJRUUFIu1Q3cKYMWMYM2YMAwcO5J133nE7NBNlDZ7Rq+pWVV0aeF0IrAHaBNN5JG1N4ios9J9Vjhw5kubNm+Pz+di0aRMXX3wxc+fOZcSIEXz66aeccMIJLFq0KOSk3bZtW5YvX87WrVujEX7cSktLIy9vM3l5Pl4JlDCePn06IoKIkJ+fz/HHH8+AAQMYPHgwEydOZPfu3VXfLxPH6nqSqrYF6ABsBvJqrD8D2AksB2YCxwbbttr2YcBiYHH79u2j9fBYo+S1J2PLysq0SZMmetxxx+n27dt19aOPHvQ0ZumLL+of/vCHg55kzcn5X9BPhC5+8EEFtH///rpr167gjqORPRmrqrpr1y4966yzDnqf8zIzay3P/Pzzz+uXX37pSAz2ZGx04ESZYiAHWAJcUsu2PCAn8Pp84Mtg29a2WAkEZ3kt0auqdu3aVVu0aFFvHfUXXnhBU1NTqyWcw7Rv37761Vdf/dhRHe2H9OmjgPbo0aPh43CxFrxXY/CBjgE9RuSghN+9e3fNynpD8/L2hx2DJfroiDjRA6nAe8DtQe6/EWgRTlu1RO84Lyb6Dh06aO/evRuso+7z+XTatGmamjpERToooCeeeOKPHdXRflGrVgpofrVa6lGrJx/H9eiD6eObNm10ypQpOnjw4IOSfrgxWKKPjvoSfTB33QgwCVijqk/UsU+rwH6ISG/8Y/87g2lrGp/XXnuNjRs3cuGFFzZYR11EuPDCC8nKepW8vA384he/YPny5Xz66acH7QdQCZQD+4FTvvsOgCeffLLhgGJcCz7eYmj97bdccsklvPXWW2zfvr1q/b59+5yPwURFMDfT9gOuBM6qdvvk+SIyXESGB/a5FFglIsuBp4Ehgb8wtbaNxoGY+HHgIunw4cNDroE+evRoAObPn8/EiRNpm5REMpCP/xayNCAzsO9DDz3ExRdf3HBAjbAefbh9tGjRgszMAoCD7pJyLAYTHXWd6ru52NCNs7w2dPPkk08q4L9QGuLY8nfffVfrxcIs0FNArwEdkJSkvz71VC0sLAzuODw2Pu71GLKzP1FA77vvvrDa29BNdGBzxjZunk70qiHfLfLBBx/oxIkT9Z///Kd+//33+t0zz2hl+/b1tm/wOBrhXTfh9pGXV6aAnnXWWWG1t0QfHfUleqtH3wh4rR79+PHjue2229i1axeHBVkZywvlA9xu75UYMjMfZ//+O7jqqqt46aWXQm5v9eijo7569FbUzLjGiycZpn7vvPMO+/ffQVJSJx4IZuJ04wmW6E3MZWb6L5fOnDnT5UhMKN566y2GDh0KQHb2p+Tn57sckQmWJXoTc9dccw0ZGRn861//cjsUE6S5c+dyySWXUFhYSGbmqyQlNXM7JBMCK1PcCHhv4pF0SktbMX3690H36YVJO9xu72YMxcWTAcjNXU9hYUdKSmzikXhiZ/SNgBeHwkWa4vNtdDsMEwTVCsrL/05SUg+Skjq6HY4Jg53RNwJenHikRw9F9QiWLXMvhnhr71YMW7dup3Vr6NMnk48/tolH4pGd0RtXtGjRguXLlzNnzhy3QzENaNWqFQA9evRwORITLkv0xhUvvfQSqampnHvuuTZJiMeJCO3bt2ez1ayJW5bojStat27N6NGj8fl8TJs2ze1wTAOaNGlCcnKy22GYMFmiN64ZMWIEANOfego6dICkJP/HwOxHQXvllcjaO9FHAsdQUlLCypUrSU9PD70/4w111UZwc7FaN87yWq2b6gDtIhL1Yl5O9OGFgmLRjsGXmamfgz4Fehxoq2qF48aOHdtwH0GwWjfRgdW6ady8VuumupbJyezx+Sir2SA/HzZubDiGDh2gZrncWto70UfU2nskhr81b87wH344ZH0aUAbs2bOHvLw8x+66sVo3zrJaN8azfD4f3Wvb4IEJNxpbDM8HkvyNwAxgB/5T+VIRVJW8vLyg+zLeYoneuOqcrCw21LbBYxNuJHoMO3bsYAlwEvAc8HOgeagxGM+yRG9c1eSUU9gBrKi+MisLHn44uA4efti/f7jtnegjzmMoLy+nT58+AAxOTQ0/BuNddQ3eu7nYxVhnefli7Jw5cxTQic2aRXXCDSf68MKkH9GI4YILLlBAH3nkEWeOowF2MTY6sIuxjZuXL8YWFhaSl5fHjTfeyHPPPedKDPHSPhoxrF+/nqOOOoo+ffqwYMGCmMRgF2Ojwy7GGs/Kzc0FICcnx+VIGqcXX3wRgDvvvNPdQExUWaI3rmvZsiVr1qxxO4xGaenSpQBccMEFLkdioqnB6pUi0g54GWgF+IACVX2qxj5nANOg6gaKN1X1gcC2F4CBwDZVPc6xyE3QvFeP/mD79vVnxowpiOSRkzOf5OQTYh5DPLR3OgbVQvbufReR5rRsmVpvu2jEYGInmDP6CmC0qnYD+gA3icgxtew3T1W7B5bqk0m+CAyIPFQTLg9ehjlIZmYBkAsUsm/f+W6H04gcuEsn2ebvTXANntGr6lZga+B1oYisAdoAnwfzBVT1PyLSIZIgTWS8WI/+YM1Q3UNmZiZNm1aydasbMXi/vfMxJHPiiSewYsUKrrrqVm699VYqKiooKSlhw4YN7N27l8svv/yQGjdWjz7+hDTxSCBh9wAW1bL5VBFZDnwL/E5VV4fY9zBgGEB7e0Cj0RERSktLOeqoo9wOpVFZunQpKSkpTJgwgQkTJhyy/S9/+QsrVqwgKcku58WzoBO9iOQAU4BRqrq3xualQL6qFonI+cBU4OhQAlHVAqAA/LdXhtLWxL9NgTotmZmZLkfSuCQnJ7Nq1SpmzJhBSkoK69at4+ijj+bYY49lwoQJTJs2jUmTJnHDDTe4HaqJQFCJXkRS8Sf5V1T1zZrbqyd+VZ0hIs+KSAtV3eFcqCaRffvttwA0bdrU5Ugan2OPPZZjjz32kPX9+/cnPT2dtWvXuhCVcVKD/4+JiACTgDWq+kQd+7QK7IeI9A70u9PJQE0Ce+UV+gwZAsC3b7/tmTrsjT2Gss6dAdi7t+Y/8Cbu1PXI7IEFOA1/EbsVwLLAcj4wHBge2OdmYDWwHFgI9K3W/lX8F3PLga+B6xv6mlYCwVleLoFQvY76MYG651aP3hsxlAa+Hw/98pfB9xEEK4EQHdRTAsH1uja1LZboneXpRJ+fX5VYfh1ILF+Cf32I7Q9aarR3oo+otfdoDL7A9+Oq7Ozg+wiCJfroqC/R26V0465q9dIvCXzcWGN9sO2DWh+NPhI0Bgl8TNm3L/g+jCdZojfuqnYrbevAxwU11gfbPqj10egjgWPoCJRnZwffh/EkS/TGXdXqqHcKrLoXKL733pDbV4mjWvBejyFNhHdF6mhg4oUleuOuoUOhoADy82klwvBAFcuZTZqE3B4R/8eCAv/6MGIIq48EjmFXbi5JGRnB92E8yerRNwJerkdf0+LFi+nVqxdPP/00t9xyiysxeLW9GzHcfffdjB07lh07dtC8eXNHYrB69NFh9ehN3KioqADgpZdeskJbHtA5cC99UVGRy5GYSFiiN55y3HH+StZLlixh3LhxLkdj0tLSAP+8siZ+WaI3npKTk0NlZSWdOnXiiSeewOfzuR1So5YamCy8rKzM5UhMJEKqXmnik9cnHjlUEvv3X01p6X00bVqESJ7nJu1wo70bMZSX++/C6dVrNampxzgag4kdO6NvBOJxqFukJQCVlXZR3k0pKacDUF4+xeVITCTsjL4R8P7EI4dau/ZndO0KZ575JO+8c1Zc3rHidHt3YsijRYvmdOmymU8+cSYGm3gk9uyM3nhSly5duOyyy5g+fTqvv/662+E0Wl988QU7d+6kV69ebodiImCJ3njW2LFjAbjsssvw+Wxqg1hSVd5//31+8pOfAHDCCYdO2G7ihyV641mdOnWqusWysLAlqqUuR5T41q5dy7hx48jIyODcc88lJSWFN954g+uvv97t0EwELNEb99UzYcbo0aMZMGAAAGVlz4Tc3okYYtLeAzHceeeddO3ald///veUlZXRNz2dNX/+M7/4xS9Ci8N4T131i91crB69szxdjz6ICTNWrFihgKan3xVWeydiiGp7l2P45ptvtEmTJgpoK9DZoD+EexxBsHr00YFNPNK4eTrRBzFhxqhRoxTQrKx3wmrvRAxRbe9CDIWFhfrSSy/pxRdfrAQmGBmQkaHFkR5HECzRR0d9id5urzTuCmLCjHnz5gGQknJOWO2diCGq7WMcQ2XlMnJze1R9fuSRR/Kb3/yGBx56KLIYjGfZGL1xVwMTZlRWVrJkyRJSUs5BJD3k9k7EEPX2MYpBVSkvn0ZRkT/Jn3feeWzYsIFvv/2WBx54wJnjMJ5kid64q4EJM0pKSgBISjoxrPZOxBD19lGM4ZvRo7npppsQEZKSkiguHgyksnDhQmbNmkWHDh2cPQ7jTXWN6bi52Bi9szw9Rq/qv9iXn68q4v9Y7eLf2rVrFdCMjGfCau9EDDFpH6UYnn76aQW0U6dO2rdvX01PH6s5OaujF0MQbIw+OojkYizQDpgLrAFWAyNr2ecMYA+wLLDcW23bAOAL4Cvgroa+nlqid5znE309xowZo4BmZ3/oWgxeaR9OH88++6wCunXrVtdiqMkSfXTUl+iDuRhbAYxW1aUikgssEZH3VfXzGvvNU9WB1VeISDLwDHAO8DXwqYi8XUtbYw6xfPlyHnroITp37sy2bT91O5y49MEHHwBQWmoPmzVmDY7Rq+pWVV0aeF2I/8y+TZD99wa+UtX1qloG/BO4KNxgTeMycuRIAB577DFE7HJSqO655x7eeOMNAJo1a+ZyNMZNId1eKSIdgB7Aolo2nyoiy4Fvgd+p6mr8fxC2VNvna+CUOvoeBgwDaG9X+R0Vf/XoD7T5FEji2msvjsta7k63D6WP8vLpFBc/AkBu7ve0a5cb8xgaam9iJ+jTJBHJAaYAo1R1b43NS4F8VT0R+Csw9UCzWrqqtTq6qhaoak9V7dmyZctgwzJB0DisRw+QnNwN8KFqs0wFQ1UpKRnBnj3JFBcPAiAvby9JSYe7HJlxW1Bn9CKSij/Jv6Kqb9bcXj3xq+oMEXlWRFrgP4NvV23XtvjP+E0MxWM9eoDbb+/Pk08uYdGitfTp8xNXYvBS+4b6+OqrdRx99HMAXHTRRVx++eVcdlluTGMIhtWjj70Gz+hFRIBJwBpVfaKOfVoF9kNEegf63Ql8ChwtIh1FJA0YArztVPAmsV1yySUATJ061d1A4kRurj+pDxw4kKlTp3LZZZe5HJHximDO6PsBVwIrRWRZYN3dQHsAVX0OuBQYISIVQAkwJHC7T4WI3Ay8ByQDLwTG7o1pUL9+/cjOzubxxx8H7nI7HM874ogjAFi4cKHLkRivaTDRq+p8ah9rr77PBGBCHdtmADPCis40aiJC//79mTlzJjk5/yM5+SduhxQXdu7ciaoS+CfbGCuBYDygnjrq48ePB2D//j+G1d6JGGLS3qEY2iQn+5N8x47uxGC8qa4nqdxc7MlYZ3n6ydgG6qj7fD4FVKR9RHXYI4kh6u0djAHQM9yKIUj2ZGx0YPXoGzdPJ/og6qhfffXVgTIIC8Jq70QMUW3vYAyAXuJWDEGyRB8d9SV6G7ox7gqijvrFF18MQGnpI2G1dyKGqLZ3KIbvAy/z3IrBeJYleuOuemqgT58+naZNmzJ48GAAUlMvDKm9EzHEpL1DMRx4mKWTWzEYz7JEb9xVRw30LbffzqBBg9i9ezfHH388WVnvkZb2m6Dbe6EWfKxjyM7IAOBe/LfJNQNu69aNr776KjYxGO+qa0zHzcXG6J3l6TF61VproF9xxRUK6Pz588Nq70QMMW3vUAyftGqld4GemJpaNRcsoG+++WZsYgiCjdFHB/WM0Yt/u7f07NlTFy9e7HYYCePAI+fxUgJh2rRpDB48mMsuu4zXXnvNlRi82D7UPiorK5k9ezZXXXUVqsq2bdto1iwppjHUJtKfR1M7EVmiqj1r22ZDN8ZzHnvsMQAKCgpcjiS+JScn8/Of/5xx48axc+dOLrjgAny+7W6HZVxgid54zmmnnQbAjBn2QLUTrrrqKrKyspg1axaFhYezZ08Wl19+OTfccAO33347hYWFbodoosyGbhqBrl399ehzcsJrf6B+eJMm4ccQSh/79g2momIa2dlzSEk5y5UYvNg+kj58vq1UVLxPScnDwNqDtiUnn0xOTvC/b069Dx5MPXHNhm4auXj6hfL5fqCiYhqpqUOrkryJXFLSkaSlXYV/+ubZZGcvIC/PR3JyPyorl1Be/obbIZooCmmGKROf4qke/fr1uznqKLjttjY8+qg7MXi1vXMxnFPV/rvv3qB79+7s2DGE2bM/ol+/flGPwerRx56d0RtPKC0tZcuWH2edLC4udjGaxqNVq1YsXryYJk2a8POf/5ydO3e6HZKJAkv0xhOuueYa2rdvz1FHHQXYZNax1LZtW5555hkKCwuZNm2a2+GYKLChG+O6f//73/zzn/8EYMiQIfTp04drrrnG3aAamTMC4ynXX389nTp1qvrcJAZL9MZ1119/PQCbN2+mXbt2DextouGII45g8ODBTJ06lTPPPJMNGzbQoUMHt8MyDrGhm8bO5Qk3tm7dyubNm+mWmkq7/Py4nvQjnmMQEd566y3m3n03AL/t2NEmHkkkddVGcHOxWjfOqrO2SJgTTezevVu/+OILXf/EEzo7PV0ngO4A3Qv6bnq6PjZkiE6YMEFfeOEFXbhwoc6aNUszM/+hKSmDtG3btnrMMcdo3759dfHixXrd6acroK9HedIPJ/qIh4lHnIjhiEB9nPE28UhcwSYeadzq/MUKY6KJTz755KBiWU4sp4JWRnnSDyf6iIeJR5yIYVm1780Mm3gkbtSX6G2MvjELY6KJ+++/H4BzzjmHy95/nybAKvxjgNlAMnA+kL1lC1999RXLli2jWbNm3HTTcYi0Zu/eVgB8+umnzJo1i+b33stQahlDjLNJPxIphhPxPzvbBbgf+LlNPBL3Gkz0ItIOeBloBfiAAlV9qo59ewELgV+p6huBdSOBG/CXyH5eVcc7E7qJWPv2sGlT7etr8dlnn/Hee+/Rv39/Zs+e7R/D3bSJX9bcMT8f2ralbdu2VXdv3Hrrwbv06tWLXr16waRJIcUQ6TFEpY8EjOHAbFXHhtqH8aRgLsZWAKNVtRvQB7hJRI6puZOIJAOPAu9VW3cc/iTfG/+JwkAROdqJwI0DQphoQtXH+eefD8B9990XcnsnYohKe4vhkD62ABcHVv0+I8MmHkkEdY3p1LUA04Bzalk/CrgJeBG4NLDul8Dfq+3zR+COhr6GjdE7q94x0SAmmsjL82lS0nEK6OGHHx5ye1X3J/3wQgyxmHgk0hiWPvRQ1fj87MMPt4lH4ghOXYwFOgCbgbwa69sAH+Efoq2e6LvhH+5rDmQBC4C/1tH3MGAxsLh9+/axem8ahUh/sTIynlNAe/bsqZWVlWH1EWlyiPosV3HQPhYx3H777QroW2+9FbUYLNFHR32JPuj76EUkB5gCjFLVvTU2jwfuVNXKGv8trME/nPM+MAtYjn8o6BCqWqCqPVW1Z8uWLYMNy0TRl19+ycCBA9m/fzjQhPnz55OUZI9eJLKFCxfSpk2bqgnZTWII6q4bEUnFn+RfUdU3a9mlJ/BPEQFoAZwvIhWqOlVVJwGTAv08AnztSOQmaubNm0f//v0PWpeT8wHp6ekuRWSibdu2bZx11lmsXr2a008/3e1wjMOCuetG8CfqNar6RG37qGrHavu/CExX1amBzw9X1W0i0h64BDjVgbhNCLZu9U88cqC8bEP2759d9To7+wP27TuToqLg29fmwGQT4fYRaXsvxODlYyguvpny8tUALF58f739OxWDiZ1gzuj7AVcCK0VkWWDd3UB7AFV9roH2U0SkOVAO3KSqu8KM1YRJQ5x4JCVlIKWlD5Gaeh0pKWdGJyjjGRUViygvf520tGFkZv7N7XBMFDSY6FV1Pv574IOiqtfU+PynoYdlnBT6xCOn0LZtGw4//DOWLvXShBnxHYNXj+H++2fxpz/Bli0Pcvjh0Y/BCmPGnj0Zaw5RUVHBN998g8/nczsUEwNTpkyhadOm2E0QictuoTCHOPBA1IEHpEzi8vl8rFy5klNOOYXAzRQmAdkZvTlIeXk548ePB+Cpp2qtdGGC5PP5qm5H9f93dOh51TfffMPLL7/Mc889R8uWLZk4cSI9e/ZERNi7dy+TJ0+moqKC0tJSLr30UiDf0Rj/97//AfDTn9oIayKzM/rGrkYN872TJlFcXMzIkSPJzs4OuX281WF3OgYV4eeZmWRnZJCcnIyIICIkJyezZ08q5eUzufbaa2nXrh1HHnkkbdu25e6772bz5s0sWbKE3r17k5uUxGNNm9KkSRNGjBjBLbfcwu9+9zs6dOhAZeVXjh7DX//6VwBOOeUUZ94H4011PUnl5mIlEJwVSj16X2am5rdooW3btq3aLRFqwTsVQ07O5/rMM8/o8OHD9YorrtBbb71V3xo1qqr984HyAScnJemFJ52k1113nZ599tmakZFxSHnmLl266AknnKCvvvqq7i4o0AXp6dq1xj5X9OunmzZt0hdffFEBFTlcs7MX6X//+19dsGCBFhcXq6pqxcsva2lmZsjv4xFHHKGA+nw+Z2raB8GejI0OrB594xZKPfoVgQTTpk0b/y+/JkYteCf6yMh4ts6a+r8PtPl54PMva2mfm/u9ZmQ8rmPHjj20lES1GJaC/j/Q5TX6+O1vf1vn1xfQzMD3L9j3ce/evSoi2q1bt5DeB0v03mSJvpGr8xdL5JBf6m6BxDFv3ryq3er8xa6lvYJ/fQ2R9uF2DPfff39VUn3nnXd07969quqfbevkWhJvaSjHEMJxZGd/pOnpD+vrr7+uY8aM0aysrEO+9mjQkgbex7KyMj3ssMMU0HfffTekGCzRe1N9id7G6BuzGnXG1wSWe5s04bTTTgu5fYPro9FHlGPw+XwMHz6c++67D5GW5OZuY+DAgeTm5gLQpEkTZrVty/HAYfgnXVkHpDkYQ3UpKf3JyLibSy+9lAcffJB9+/bh8/nQ/HymBfb5C9ASeAfw1THZeo8ePdi9eze33nrrj3dXOfFeGk+yRN+Y1ahhviDw8ZQbbgirPeC5Ouw+n49Vq1b5/30NsY9Vw4fTrFkz/va3v9GzZ09yczeTlHToveYt/vxnVmRlsQt4F+gU6jEEcRz1ERF4+GEuzMqiEvgTUARcCCRv3sz111/P/v37q/YvL3+D1atX0759+6o7rCKNwXhcXaf6bi42dOOsYOvRP9G0qQI6f/78g3bxUi34+++/Xzt27Kjdu3fX/v3769FHH60ntG+vzzVrdlD7PXv26HHHHVc1nNG9e/eQYnjvzjur2vbp00fLy8vjqh791rZtdWjfvgcN6YwZM0bT0+9TQDt27KhFRUXROY4G2NBNdGBj9I1bsL9YK1euVEC7du1adSFW1f066pWVlZqd/YGKNKvzYiSggwYN0u3bt+t1111XtS4lJaXqdV7eLs3N3awffPCB7t69Wzds2KB9+vTRdu3aaU5Ojnbr1k27dOmiLVq0qGrzn//8x7H3wY330efz6ahRo2q8Vzn62WefxSyGmizRR4cl+kYulF+sjh07ViWEA2d8sU5Qa9eu1aeeekovv/xy7dev30FJqmPHjlpYWKglJSVaXl6uqqqlpaV65plnHrRfRkaGTpgwQVVVJ02aVO8fCECPPPJITUpKUkAHDBigI0aMCO0/G4ffA6f7qKio0Dlz5mh6+p2am7vFlRgOsEQfHfUlensy1hxk3bp1VU9z5uTkkJmZSUXF5Yik87Of/Y+8vDyaN2/OxIkTSU1Ndfzrl5eX06VLl4PWdevWjXXr+pOWdjvr13c5pE1aWhpz5szhd7/7HUuWLOH666/nyiuvrNp+1VVXsXjxYgoKlpCcfArXXltGy5YtKSoq4uyzz+aCCy5w/Di8Jjk5mbPOOouMjLPcDsW4QPx/CLylZ8+eunjxYrfDSBhdu/rr0efkBLd/Wdlkysqex+fbBCiqZUAp/kt8pVX7pab+irS0G0hOPhWRrDp68ztQg7xJk7r3UVX27Tudysp5pKf/gfT0PyGSFnT7hkTah9vtEy0GD6aeuCYiS1S1Z23b7K6bRiDUX6i0tMvJyZlLXt568vI2AN8AO2jSZD9NmihpaaMAKC9/jX37zmbv3sMoK3s5whh3UVTUk8rKeUAO6en3VSV5Y0xkbOimEQi9Hv3BDq0//iTwJKtXr+b555/nqaeeoqTkan772xWMGzcuyD5+VFpaSk7O4fh8FQwcOJC33377kEqKXq3lHsv2iRKD1aOPPTujN2E79thjGT9+PPv27aN379785S9/4ZZbbgmpj7Vr1zJ48GAqKioYPXo077zzjpXLNcZhluhNxLKyspg1axZt2rRhwoQJjBkzJqh2N9xwA127dmXWrFkA3HPPPdEM05hGyxK9cUTTpk1ZsMD/bO3DDz9MSUlJnfvu3buX5s2b8/e//52kpCSmTZtGRUUFTZs2jVW4xjQqluiNY9q1a0dBQQHgP8t/4403KC8vP2ifDRs20KVLF3744QdOO+009u3bx4UXXkhycrIbIRvTKFiiN5GpMVHFtenpnHzyyQD88pe/JC0tjaFDh1JZuYZ9+86hU6dOfP/999x1113MmzePjIwMm+zCmChrMNGLSDsRmSsia0RktYiMrGffXiJSKSKXVlt3W6DdKhF5VUQynAreuOyVV2DYMNi0yX8P56ZNpIwYweLbbuObb75hxIgRdOnShcmTJ1NUdAwVFf+mU6dOTJ8+nbFjx9bZB8OGWbI3xkHBnNFXAKNVtRvQB7hJRI6puZOIJAOPAu9VW9cGuBXoqarHAcnAECcCNx5wzz1QXHzwuuJiuOceWrduzbPPPssXX3zBHXfcQUrKuWRnz2HdunUHP4laTx/GGGc0mOhVdauqLg28LsRfsrxNLbveAkwBttVYnwJkikgKkAV8G1HExjs2bw5q/aOPPkp29nukpNTy+H2QfRhjwhfSGL2IdAB6AItqrG8DXAw8V329qn4DjAM2A1uBPao6u46+h4nIYhFZvH379lDCMm7xwsQjxpgGBZ3oRSQH/xn7KFXdW2PzeOBOVa2s0aYpcBHQEWgNZIvIFbX1r6oFqtpTVXu2bHno5A7Gg7ww8YgxpkFBlUAQkVT8Sf4VVX2zll16Av8MPNHYAjhfRCqAVGCDqm4P9PMm0Bf4hwOxG7cNHer/eM89/qGW9u39CfrA+lj1YYypV4OJXvzZexKwRlWfqG0fVe1Ybf8XgemqOlVETgH6iL+0YQnwM8DKUiaSoUMjT8pO9GGMqVMwZ/T9gCuBlSKyLLDubqA9gKo+V0c7VHWRiLwBLMV/985nQEEkARtjjAlNg4leVecDQVeZUtVranx+H3BfyJEZY4xxhJUpbiSWLQu/PGxRkf9jJOVlI+0jEWJIhGNwoo9ly6B79/C/vgmdJfpG4PLLI2vfqRNEWjk40j4SIYZEOAYn+ujePfKfSRMaS/SNwLBh/sUY0zhZUTNjjElwluiNMSbBWaI3xpgEZ4neGGMSnCV6Y4xJcJbojTEmwVmiN8aYBGeJ3hhjEpyoqtsxHEJEtgObwmzeAtjhYDheZ8eb2Ox4E5fTx5qvqrVO5uHJRB8JEVmsqj3djiNW7HgTmx1v4orlsdrQjTHGJDhL9MYYk+ASMdE3tolN7HgTmx1v4orZsSbcGL0xxpiDJeIZvTHGmGos0RtjTIKLm0QvIu1EZK6IrBGR1SIyspZ9morIWyKyQkT+KyLHVds2QES+EJGvROSu2EYfmkiONZi2XhPp9zawPVlEPhOR6bGLPDwO/CwfJiJviMj/An2cGtsjCI0Dx3tboN0qEXlVRDJiewShEZGMwDEsD8R9fy37iIg8HchHK0TkpGrbnM9VqhoXC3AkcFLgdS6wFjimxj6PA/cFXv8EmBN4nQysAzoBacDymm29tER4rA229doSyfFW2347MBmY7vbxRPt4gZeA3wRepwGHuX1M0TpeoA2wAcgMfP4v4Bq3j6mB4xUgJ/A6FVgE9Kmxz/nAzMC+fYBFgfVRyVVxc0avqltVdWngdSGwBv8PQXXHAHMC+/wP6CAiRwC9ga9Udb2qlgH/BC6KWfAhiuRYg2zrKRF+bxGRtsAFwN9jFnQEIjleEckD+gOTAtvKVHV3rGIPR6TfX/xTnmaKSAqQBXwbk8DDpH6BKdRJDSw173q5CHg5sO9C4DAROZIo5aq4SfTViUgHoAf+v5TVLQcuCezTG8gH2uL/odpSbb+v8XjyOyCMYw2mrWeFebzjgTsAX0yCdFAYx9sJ2A78X2Co6u8ikh27iCMT6vGq6jfAOGAzsBXYo6qzYxZwmAJDicuAbcD7qlrzeOvKSVHJVXGX6EUkB5gCjFLVvTU2/xloGniDbwE+Ayrw/3tUk+fvKw3zWINp60nhHK+IDAS2qeqSmAbrgDC/vynAScBEVe0B7AM8fc3pgDC/v03xn9F2BFoD2SJyReyiDo+qVqpqd/x/nHvXvKZE3TkpKrkqJdIOYklEUvH/oLyiqm/W3B744bk2sK/gH9vbgP/fvXbVdm2Lx//9i+BYG2zrRREc7xDgQhE5H8gA8kTkH6rq6WQQ4c/y19XOEN8gDhJ9BMd7HrBBVbcHtr0J9AX+EaPQI6Kqu0XkQ2AAsKrapq+pPSel1bE+InFzRh/45k8C1qjqE3Xsc5iIpAU+/Q3wn8AP0KfA0SLSMbB9CPB2LOIORyTHGkxbr4nkeFX1D6raVlU74P++fhAHST6S4/0O2CIiXQPbfgZ8HvWgIxDh7+5moI+IZAX6+Rn+MX7PEpGWInJY4HUmcDbwvxq7vQ1cFbj7pg/+IamtRClXxdMZfT/gSmBl4N87gLuB9gCq+hzQDXhZRCrx//BfH9hWISI3A+/hv6r9gqqujm34IQn7WOtqq6ozYhN6WCI53ngU6fHeArwSSATrCZwJe1gkv7uLROQNYCn+oavP8H6ZhCOBl0QkGf/J9L9UdbqIDIeq452B/86br4BiAt/DaOUqK4FgjDEJLm6GbowxxoTHEr0xxiQ4S/TGGJPgLNEbY0yCs0RvjDEuE5EXRGSbiKxqeG8QkctE5HPxF02b3OD+dteNMca4S0T6A0X469/UfIq25r5H4y/udpaq7hKRw1V1W31t7IzeGGNcpqr/AX6ovk5EjhKRWSKyRETmichPAptuAJ5R1V2BtvUmebBEb4wxXlUA3KKqJwO/A54NrO8CdBGRj0VkoYgMaKijeHoy1hhjGoVAAbi+wOv+yg8ApAc+pgBHA2fgr4UzT0SOq69ctSV6Y4zxniRgd6ACZk1fAwtVtRzYICJf4E/8n9bXmTHGGA8JFHTbICK/hKqpB08MbJ4KnBlY3wL/UM76+vqzRG+MMS4TkVeBBUBXEflaRK4HhgLXi8hyYDU/zjT1HrBTRD4H5gK/V9Wd9fZvt1caY0xiszN6Y4xJcJbojTEmwVmiN8aYBGeJ3hhjEpwlemOMSXCW6I0xJsFZojfGmAT3/wFuDT9lTJ30JwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# ERA5 evapotranspiration\n",
    "era5_fileName = r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\era5\\potential_evapotranspiration-19790101_20181231-eusalp-qm_era5.nc'\n",
    "\n",
    "pet = readnetcdf_in_shp(era5_fileName, shp_catchment_fileName,plot=True)['pet']\n",
    "pet = xarray2df(pet.resample(time='1d').sum(skipna=False), 'E','pet')\n",
    "\n",
    "check_data_gap(pet)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "cfe2f5a2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n of pixels75\n",
      "Date start: 1979-01-01, date end: 2018-12-31\n",
      "No missing dates\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEQCAYAAAC+z7+sAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA7hklEQVR4nO3deXwV5fX48c/JviObIlsAEYorKCCCxaUuVEHRWkvF3YpQF1BatUq1rlSlihbFxuJX/VWsVRQUAbGIFhSogKxSUXYVZZElISHbPb8/7iWGkOUuc+/MvTnv12teuZmZ58mZm+Rk8szMeURVMcYYk7iS3A7AGGNMdFmiN8aYBGeJ3hhjEpwlemOMSXCW6I0xJsFZojfGmATn2UQvIi+IyDYRWRXk/peJyOcislpEJkc7PmOMiRfi1fvoRaQ/UAS8rKrHNbDv0cC/gLNUdZeIHK6q22IRpzHGeJ1nz+hV9T/AD9XXichRIjJLRJaIyDwR+Ulg0w3AM6q6K9DWkrwxxgR4NtHXoQC4RVVPBn4HPBtY3wXoIiIfi8hCERngWoTGGOMxKW4HECwRyQH6Aq+LyIHV6YGPKcDRwBlAW2CeiBynqrtjHKYxxnhO3CR6/P997FbV7rVs+xpYqKrlwAYR+QJ/4v80hvEZY4wnxc3QjaruxZ/EfwkgficGNk8Fzgysb4F/KGe9G3EaY4zXeDbRi8irwAKgq4h8LSLXA0OB60VkObAauCiw+3vAThH5HJgL/F5Vd7oRtzHGeI1nb680xhjjDM+e0RtjjHGGJy/GtmjRQjt06OB2GMYYEzeWLFmyQ1Vb1rbNk4m+Q4cOLF682O0wjDEmbojIprq22dCNMcYkOEv0xhiT4CzRG2NMgrNEb4wxCc4SvTHGJDhL9MYYk+As0RtjTILz5H30xlkFBTA5gskVP/rI//H0093rIxFiiKdjKCn5irKy78nN7YXPt5+UlDzHY7DqK7Fjib4RmDwZli2D7t3djsR4mc9XxrffPsu6dbcdsi03tzfHHz+T1NRmLkRmIuXJomY9e/ZUezLWOWec4f/44YfhtT/sMP/H3bvDjyHSPhIhBq8dQ2lpKc8//zzz589n8eLFrFu3rmq/7t2706FDBzp37sycOXP47LPPGDRoEG+//XbEMUT682hqJyJLVLVnbdvsjN6YRka1jAULlvCrX/2KLVu2HLTtxhtv5KmnniI9Pf2g9YMGDeLdd9+ltLSUHyd2M/HCEr0xCUhV+frrr/noo49o3rw5/fr147XXXqO8vA3FxRfQt69/vx49ejB37lxSU1PJzMyk2jSdBxk0aBDTp0/ngQceAB6O3YEYR1iiNybB7Nq1ixNPPPGQs/Xqfv3rXzNq1Ci6d+9OWlpag31efPHF3HjjjTzyyCOkpwsZGQ85GbKJMru90hiveOUV6NABkpL8H195JeT2mp/Ppc2asWXLFi486SQmT57MXXfdxe23386f/vQnMjMnk509l8mTJ9O7d+9Dk3wdMbRs2ZLCwkIASksfps5re5Eeg4kOVfXccvLJJ6txzumn+5dwNWniXyIRaR+JEEO97f/xD9WsLFX/XYf+JSvLvz6YPgLtnwMF9K5Q2wcZw6hRoxRQkbZhH0OkP4+mdsBirSOn2l03jYDddeONGOpt36EDbKqlnHh+Pmzc2GAfO9q1Y8jXXzMHyAG+AfJCaB9sDKWlpWRkZJCU1I3Kys/DOga76yY66rvrpsGhGxFpJyJzRWSNiKwWkZG17HOGiOwRkWWB5d4a25NF5DMRmR7+YRiTwDZvDm19NcuWLaNlIMn3AD4jkOSDbB9KDOeddx4AqamXhtXeuCOYMfoKYLSqdgP6ADeJyDG17DdPVbsHlgdqbBsJrIkwVmMSV/v2oa0PWLBgAT169ADgQWAp0DmE9qHEMGfOHD4KPNaann5PyO2NexpM9Kq6VVWXBl4X4k/YbYL9AiLSFrgA+Hu4QRqT8B5+GLKyDl6XleVfX4cpU6bQN3Cf5Mzf/54xIbYPNYYPA2MtOTmrEKnlXvowjsHERkh33YhIB/z/HS6qZfOpIrJcRGaKyLHV1o8H7gB8DfQ9TEQWi8ji7du3hxKWMfFv6FB/UaL8fBDxfywo8K+vxT333MOll/qHT15++WUGPPZYSO3DiaF7oIZGZeV/HTkGEztB30cvIjnAFGCUqu6tsXkpkK+qRSJyPjAVOFpEBgLbVHWJiJxRX/+qWgAUgP9ibNBHYEyiGDq01qQ4ffp0jjnmGDp16oTqLgoLO/PIIz/QqVMn3nrrLU444YR62zsRA0BmZmbgVT3nbE7EYBwXVKIXkVT8Sf4VVX2z5vbqiV9VZ4jIsyLSAugHXBhI/hlAnoj8Q1WvcCZ8YxLbHXfcweOPP37I+mbNmrFq1apqyTf6HnnkEQCSk3vH7GsaZwRz140Ak4A1qvpEHfu0CuyHiPQO9LtTVf+gqm1VtQMwBPjAkrwxwdm4cSOPP/44WVlZdOvWjcMOO4yUlHPJzCxg+/btMU3y4H9oyq88pl/XRC6YM/p+wJXAShFZFlh3N9AeQFWfAy4FRohIBVACDFEv3qBvTBy56KKLAJg9ezb9+vUDfrwPPsmFZ9qvvvpqpk6dSlHRyeTkLMV/uc7EA3tgqhHo2hXWrYOcnPDa79nj/9ikSfgxRNpHIsQQSnuf7wcKC5uTknIp2dmvuxJDbcrL36a4+CLgCPLyvkEkOewYPJh64lpED0yZ+Ge/UPFFtZzi4l8CkJZ2tcvRHCw19ULgZuB7Kivnuh2OCZJVr2wEWrf2L1YCwd0Ygm3/+OPjueOODzjvvPOYOfMCqlcOdvsYVJXU1AVUVsLKlUfRsWPofRwogWBix87ojfGYA+WFU1NTWbFiBfv27eOjjz6qu2JkDM2cOZPKyiWkp99Jx3CyvHGFJXpjPGbgwIGkpaUxffp0unfvTk5ODmeccQYvvPCC26Fx4GHG8vJ3XI7EhMISvTEec+6551JaWsr8+fMZNGgQvXv771u/8cYbUS1zNbaLL74YSMHn+5yvvvrK1VhM8CzRG+MVNSbt6LdxI2+//TaLFi1i7NixVFZWsn9/LcXE6mgf1qQfDfSRl5dHdrb/IuyvfvWr6MRgHGeJ3hgveOUVGDbMX89d1f9x2LCqRHnnnXcGdiwJq70TMRwgcgRA1YxTjsZgosISvTFecM89UFx88LriYv964MfnSuqY37WB9k7EAFBSUsK+fecA8PTTTzsfg4kKS/TGeEEDk3bsOfCUUZjtnYgBoKCgANVNpKXdyoABA5yPwUSFJXpjvKCBSTsOXJCts0SwE5N+BNHHihUrAEhNHRydGExUWKI3xgsamLTjwK2VSUm1Te7WcHsnYgD4wx/+AMD+/TUnkXMoBhMVluiN8YIGJu2YMWMG2dnZZGUVhNXeiRgAOnfuTErKOVRWfkhJSUnI7Y07LNEb4xVDh8LGjeDzwcaNbP7pTxk2bBgiwvvvv88VVzRQ4btG+7ASbBB9JCefBsD3338fnRiM4yzRG+NBq1atIj8/n+effx6AQYMGMWHCBJej8hPx16X/+uuvXY7EBMsSvTEeU1BQwEknnQTAX//6V1SVt99+m5QUr9QgrACgsrLS5ThMsLzyk2OiaOtWfz36A5ULQ3Xgzr5w2zvRRyLEUFv7ioqPqKj4NyLNSE29mtLSBygrewqAzMxJjBlzHWPGRDeGUO3f7y95PmDAZ6Snnx52DCZ2LNE3Ah4oemhqUC2jsLA1qjur1u3ffzsASUmdyMn5DJE8t8JrwLEA+Hwb3Q3DBM0SfSNg9ei9EUP19jfffDvPPLOT008/nSeffJKFCxdy//3307p1a+bMmUPTprUnebePASAvr5zCQrj66mIK6rgJqD5Wjz72LNEb44IlS5YAcPPNN9OjRw969OjBiBEjXI4qOCUl1wJw5plnuhyJCZZdjDXGBcOHDwdg06ZNLkcSmrVr11JR8Q6pqZfx61//2u1wTJAs0RvjgnXr1gHQt29flyMJzaRJkwBITraz+XjSYKIXkXYiMldE1ojIahEZWcs+Z4jIHhFZFljuDazPEJH/isjyQNv7o3EQxrguxDrsH330EfBjDZtw+og0hnD6OPXUUwOvSqMXg3Geqta7AEcCJwVe5wJrgWNq7HMGML2WtgLkBF6nAouAPg19zZNPPlmNc04/3b+Eq0kT/xKJSPvwdAz/+IdqVpaq/wYn/5KV5V9fo3129qf605/+VAH1//qF3kekMdT5HgTRx44dOxTQ9PSHwo4h0p9HUztgsdaRUxs8o1fVraq6NPC6EFgDtAnyj4iqalHg09TAYjf7mcQSZB121X3s29efefPmAXDdddeF3EekMUTaR7NmzQKvajmjt3r0nhXSGL2IdAB64D8zr+nUwBDNTBE5tlqbZBFZBmwD3lfV2toiIsNEZLGILD4wAbExcaGeOux//OMfGThwICJCYWEXoITJkyejqlXj3Q31EWkMQQuiDxEB0oFa5q61evSeFXSiF5EcYAowSlX31ti8FMhX1ROBvwJTD2xQ1UpV7Q60BXqLyHG19a+qBaraU1V7tmzZMrSjMMZNtdRb3wQkqfLQQw/x7rvvAiCSQUrKJQwZMiSoPupd73T7kPpIorLy8+jEYKIiqEQvIqn4k/wrqvpmze2quvfAEI2qzgBSRaRFjX12Ax8CA2q2Nyau1ajDrsC5Iihw/PHHU1pais/nIzd3HdnZUwJnxfX3AYRWyz1G9ej9Sqjx6+1cDCYqgrnrRoBJwBpVfaKOfVoF9kNEegf63SkiLUXksMD6TOBs4H8OxW6MN9Sowz61ZUvWBs7mV6xYQVpaWu3JvZ4+Qq7lHqN69P/+978BSEqq5TKd1aP3LNEGCqGIyGnAPGAl4AusvhtoD6Cqz4nIzcAI/GXtSoDbVfUTETkBeAlIxp/8/6WqtUxNc7CePXvqj5Mhm0gdeOTcSiDEJobf/va3TJw4keLiYjIzMx37+k70EUn7ESNG8NxzzwGQk/MlhYWdw4oh0p9HUzsRWaKqPWvb1mAJBFWdj/82yfr2mQAcUixbVVfgv3hrTKPxxRdfALB///6DEn0827dvX1WSz8wsIDk5vCRv3GFPxhrjoOeff54PPvgAgIb+W44Hqsrs2bMZE6iVPHbsWNLSbnA5KhMqK2pmjEM+/vhjhg0bBsDChQur3XMevwYOHMiMGTOqPu/Tp4+L0ZhwWaJvBGzikdjEUFq6HICsrJmcd94pjn99J/oItf3evf4qmzk5q0lKasfgwbmOxWBix4ZuGoEEGEHwvP37H2T//puAVJKTe7kdjmNUd5KWNpLk5GMQyXU7HBMmO6NvBGzikejHEKjjx3//+zG9ejWPytd3oo9Q2q9cuZITTqjgwQdbc8cdzsVgE4/Enp3RGxOhF154AYBzzjmHXr0S52z+1VdfBeD000OfF9Z4i53RGwPs3r2bmTNn0rlzZ3r27NngA07dunUjKyuLM888k7/85S/k5uYyfvz42AQbA3PnzmXs2LF07tyZU0459HqDiS+W6E2j9/HHH3PaaadVfd6+fXuWLVtG06ZND9nX59vEUUedxfr16wFYunQprVq1YvXq1XF9l42qsnHjRvbs2cOf//xnXnvtNQBGjRrlbmDGETZ0Y+JfBJNdvP/++5x22mkkJyVxZno6AJs3b6ZZs2bMnz+fefPmMXHiRM4++2yKis6ksLAD69evJzU1lTfeeINx48bxySef0GzmzKhP+hHN9o8++iidOnWiR48evPbaa6QBn7duzU2R3CJkvKOuQvVuLjbxiLMSeuKRICe7qNnHnj179L777quaAGRhenpV+0cD62pbUlIu0rfeeiusGCKd9MOJ96G29lOmTFFAjzzsMH0zLU3ngJaGexxBsIlHooN6Jh5xPanXtliid1ZCJ/r8/IOT24ElP7/OPkpLSzUnJ0cBbd68uc46/PBD2i8FvS03V//4xz/qkiVLdM+ePRHHUO97EGkfYba/4447qv6ILWrVKvLjCIIl+uioL9HbGL2Jb2FMdnHTTTdRVFTE6NGjGTdunH+oo4YeQI+iInigwRp8MZv0w8n227Zto3Xr1lRWVtKiRQs+/PBDjj3++MhiMJ5lY/QmvoUx2cWiRf5JzsaNGxd2H5HG4HgfIbZ/7LHHqKyspGPHjmzYsIFjjz3WJg5JYJboTXwLcbKLysovWblyJWdUf2onrib9cKZ9amoqACNHjiQnJ8eZGIx31TWm4+ZiY/TOSugxelX/xcL8fFUR/8daLsQuX778oIuqK1asCLmPSGNo8D2IYQzdunVTQL/77jvnj6MBNkYfHdQzRt/gxCNusIlHnGUTj0B2djbFxcUkJf2EJ54YzsiRI2Meg1feR9Uy9u5Np2/fvnz88ccxj8EmHomOiCYeMSbeXXnllRQXF5Oefi8ZGfcTRo5PKEVF/iddr732WpcjMbFiid4krKKiIk499VRWrVoFQHr6fS5H5L7y8nfw+ZaRk5PDdddd53Y4JkYs0TcCjbEefVnZa5SUDAd2k5R0ItnZs9m7NymmMTjd3ok+iotnASCynmbNwrsXw+rRxx+766YR8OBlmKgqKRlNSckQYDdpabeRm7uMpKTD3Q7LI7YAIFJ7KWWTmOyMvhHwYj369957j6KiIn7xi184HkPHjm/y3XcZbNmyhRYtWoTVR6QxRKO9E32IfAjAnj3hn+NZPfr40+B3W0TaichcEVkjIqtF5JBLWSJyhojsEZFlgeXeYNuaxmnAgAFceumlB40TT506lcGDB3PllVeyYMGCsPvesWMHvXr1OijJG/D5fEAhAA8//DDbtm1zNyATM8Gc0VcAo1V1qfjnElsiIu+r6uc19punqgPDbGsakbKysqrX//d//8dLL71Eu3bt2LRpU9X6L774gg8++ODHh3lC7L+iosKRWBNJRUUFIu1Q3cKYMWMYM2YMAwcO5J133nE7NBNlDZ7Rq+pWVV0aeF0IrAHaBNN5JG1N4ios9J9Vjhw5kubNm+Pz+di0aRMXX3wxc+fOZcSIEXz66aeccMIJLFq0KOSk3bZtW5YvX87WrVujEX7cSktLIy9vM3l5Pl4JlDCePn06IoKIkJ+fz/HHH8+AAQMYPHgwEydOZPfu3VXfLxPH6nqSqrYF6ABsBvJqrD8D2AksB2YCxwbbttr2YcBiYHH79u2j9fBYo+S1J2PLysq0SZMmetxxx+n27dt19aOPHvQ0ZumLL+of/vCHg55kzcn5X9BPhC5+8EEFtH///rpr167gjqORPRmrqrpr1y4966yzDnqf8zIzay3P/Pzzz+uXX37pSAz2ZGx04ESZYiAHWAJcUsu2PCAn8Pp84Mtg29a2WAkEZ3kt0auqdu3aVVu0aFFvHfUXXnhBU1NTqyWcw7Rv37761Vdf/dhRHe2H9OmjgPbo0aPh43CxFrxXY/CBjgE9RuSghN+9e3fNynpD8/L2hx2DJfroiDjRA6nAe8DtQe6/EWgRTlu1RO84Lyb6Dh06aO/evRuso+7z+XTatGmamjpERToooCeeeOKPHdXRflGrVgpofrVa6lGrJx/H9eiD6eObNm10ypQpOnjw4IOSfrgxWKKPjvoSfTB33QgwCVijqk/UsU+rwH6ISG/8Y/87g2lrGp/XXnuNjRs3cuGFFzZYR11EuPDCC8nKepW8vA384he/YPny5Xz66acH7QdQCZQD+4FTvvsOgCeffLLhgGJcCz7eYmj97bdccsklvPXWW2zfvr1q/b59+5yPwURFMDfT9gOuBM6qdvvk+SIyXESGB/a5FFglIsuBp4Ehgb8wtbaNxoGY+HHgIunw4cNDroE+evRoAObPn8/EiRNpm5REMpCP/xayNCAzsO9DDz3ExRdf3HBAjbAefbh9tGjRgszMAoCD7pJyLAYTHXWd6ru52NCNs7w2dPPkk08q4L9QGuLY8nfffVfrxcIs0FNArwEdkJSkvz71VC0sLAzuODw2Pu71GLKzP1FA77vvvrDa29BNdGBzxjZunk70qiHfLfLBBx/oxIkT9Z///Kd+//33+t0zz2hl+/b1tm/wOBrhXTfh9pGXV6aAnnXWWWG1t0QfHfUleqtH3wh4rR79+PHjue2229i1axeHBVkZywvlA9xu75UYMjMfZ//+O7jqqqt46aWXQm5v9eijo7569FbUzLjGiycZpn7vvPMO+/ffQVJSJx4IZuJ04wmW6E3MZWb6L5fOnDnT5UhMKN566y2GDh0KQHb2p+Tn57sckQmWJXoTc9dccw0ZGRn861//cjsUE6S5c+dyySWXUFhYSGbmqyQlNXM7JBMCK1PcCHhv4pF0SktbMX3690H36YVJO9xu72YMxcWTAcjNXU9hYUdKSmzikXhiZ/SNgBeHwkWa4vNtdDsMEwTVCsrL/05SUg+Skjq6HY4Jg53RNwJenHikRw9F9QiWLXMvhnhr71YMW7dup3Vr6NMnk48/tolH4pGd0RtXtGjRguXLlzNnzhy3QzENaNWqFQA9evRwORITLkv0xhUvvfQSqampnHvuuTZJiMeJCO3bt2ez1ayJW5bojStat27N6NGj8fl8TJs2ze1wTAOaNGlCcnKy22GYMFmiN64ZMWIEANOfego6dICkJP/HwOxHQXvllcjaO9FHAsdQUlLCypUrSU9PD70/4w111UZwc7FaN87yWq2b6gDtIhL1Yl5O9OGFgmLRjsGXmamfgz4Fehxoq2qF48aOHdtwH0GwWjfRgdW6ady8VuumupbJyezx+Sir2SA/HzZubDiGDh2gZrncWto70UfU2nskhr81b87wH344ZH0aUAbs2bOHvLw8x+66sVo3zrJaN8azfD4f3Wvb4IEJNxpbDM8HkvyNwAxgB/5T+VIRVJW8vLyg+zLeYoneuOqcrCw21LbBYxNuJHoMO3bsYAlwEvAc8HOgeagxGM+yRG9c1eSUU9gBrKi+MisLHn44uA4efti/f7jtnegjzmMoLy+nT58+AAxOTQ0/BuNddQ3eu7nYxVhnefli7Jw5cxTQic2aRXXCDSf68MKkH9GI4YILLlBAH3nkEWeOowF2MTY6sIuxjZuXL8YWFhaSl5fHjTfeyHPPPedKDPHSPhoxrF+/nqOOOoo+ffqwYMGCmMRgF2Ojwy7GGs/Kzc0FICcnx+VIGqcXX3wRgDvvvNPdQExUWaI3rmvZsiVr1qxxO4xGaenSpQBccMEFLkdioqnB6pUi0g54GWgF+IACVX2qxj5nANOg6gaKN1X1gcC2F4CBwDZVPc6xyE3QvFeP/mD79vVnxowpiOSRkzOf5OQTYh5DPLR3OgbVQvbufReR5rRsmVpvu2jEYGInmDP6CmC0qnYD+gA3icgxtew3T1W7B5bqk0m+CAyIPFQTLg9ehjlIZmYBkAsUsm/f+W6H04gcuEsn2ebvTXANntGr6lZga+B1oYisAdoAnwfzBVT1PyLSIZIgTWS8WI/+YM1Q3UNmZiZNm1aydasbMXi/vfMxJHPiiSewYsUKrrrqVm699VYqKiooKSlhw4YN7N27l8svv/yQGjdWjz7+hDTxSCBh9wAW1bL5VBFZDnwL/E5VV4fY9zBgGEB7e0Cj0RERSktLOeqoo9wOpVFZunQpKSkpTJgwgQkTJhyy/S9/+QsrVqwgKcku58WzoBO9iOQAU4BRqrq3xualQL6qFonI+cBU4OhQAlHVAqAA/LdXhtLWxL9NgTotmZmZLkfSuCQnJ7Nq1SpmzJhBSkoK69at4+ijj+bYY49lwoQJTJs2jUmTJnHDDTe4HaqJQFCJXkRS8Sf5V1T1zZrbqyd+VZ0hIs+KSAtV3eFcqCaRffvttwA0bdrU5Ugan2OPPZZjjz32kPX9+/cnPT2dtWvXuhCVcVKD/4+JiACTgDWq+kQd+7QK7IeI9A70u9PJQE0Ce+UV+gwZAsC3b7/tmTrsjT2Gss6dAdi7t+Y/8Cbu1PXI7IEFOA1/EbsVwLLAcj4wHBge2OdmYDWwHFgI9K3W/lX8F3PLga+B6xv6mlYCwVleLoFQvY76MYG651aP3hsxlAa+Hw/98pfB9xEEK4EQHdRTAsH1uja1LZboneXpRJ+fX5VYfh1ILF+Cf32I7Q9aarR3oo+otfdoDL7A9+Oq7Ozg+wiCJfroqC/R26V0465q9dIvCXzcWGN9sO2DWh+NPhI0Bgl8TNm3L/g+jCdZojfuqnYrbevAxwU11gfbPqj10egjgWPoCJRnZwffh/EkS/TGXdXqqHcKrLoXKL733pDbV4mjWvBejyFNhHdF6mhg4oUleuOuoUOhoADy82klwvBAFcuZTZqE3B4R/8eCAv/6MGIIq48EjmFXbi5JGRnB92E8yerRNwJerkdf0+LFi+nVqxdPP/00t9xyiysxeLW9GzHcfffdjB07lh07dtC8eXNHYrB69NFh9ehN3KioqADgpZdeskJbHtA5cC99UVGRy5GYSFiiN55y3HH+StZLlixh3LhxLkdj0tLSAP+8siZ+WaI3npKTk0NlZSWdOnXiiSeewOfzuR1So5YamCy8rKzM5UhMJEKqXmnik9cnHjlUEvv3X01p6X00bVqESJ7nJu1wo70bMZSX++/C6dVrNampxzgag4kdO6NvBOJxqFukJQCVlXZR3k0pKacDUF4+xeVITCTsjL4R8P7EI4dau/ZndO0KZ575JO+8c1Zc3rHidHt3YsijRYvmdOmymU8+cSYGm3gk9uyM3nhSly5duOyyy5g+fTqvv/662+E0Wl988QU7d+6kV69ebodiImCJ3njW2LFjAbjsssvw+Wxqg1hSVd5//31+8pOfAHDCCYdO2G7ihyV641mdOnWqusWysLAlqqUuR5T41q5dy7hx48jIyODcc88lJSWFN954g+uvv97t0EwELNEb99UzYcbo0aMZMGAAAGVlz4Tc3okYYtLeAzHceeeddO3ald///veUlZXRNz2dNX/+M7/4xS9Ci8N4T131i91crB69szxdjz6ICTNWrFihgKan3xVWeydiiGp7l2P45ptvtEmTJgpoK9DZoD+EexxBsHr00YFNPNK4eTrRBzFhxqhRoxTQrKx3wmrvRAxRbe9CDIWFhfrSSy/pxRdfrAQmGBmQkaHFkR5HECzRR0d9id5urzTuCmLCjHnz5gGQknJOWO2diCGq7WMcQ2XlMnJze1R9fuSRR/Kb3/yGBx56KLIYjGfZGL1xVwMTZlRWVrJkyRJSUs5BJD3k9k7EEPX2MYpBVSkvn0ZRkT/Jn3feeWzYsIFvv/2WBx54wJnjMJ5kid64q4EJM0pKSgBISjoxrPZOxBD19lGM4ZvRo7npppsQEZKSkiguHgyksnDhQmbNmkWHDh2cPQ7jTXWN6bi52Bi9szw9Rq/qv9iXn68q4v9Y7eLf2rVrFdCMjGfCau9EDDFpH6UYnn76aQW0U6dO2rdvX01PH6s5OaujF0MQbIw+OojkYizQDpgLrAFWAyNr2ecMYA+wLLDcW23bAOAL4Cvgroa+nlqid5znE309xowZo4BmZ3/oWgxeaR9OH88++6wCunXrVtdiqMkSfXTUl+iDuRhbAYxW1aUikgssEZH3VfXzGvvNU9WB1VeISDLwDHAO8DXwqYi8XUtbYw6xfPlyHnroITp37sy2bT91O5y49MEHHwBQWmoPmzVmDY7Rq+pWVV0aeF2I/8y+TZD99wa+UtX1qloG/BO4KNxgTeMycuRIAB577DFE7HJSqO655x7eeOMNAJo1a+ZyNMZNId1eKSIdgB7Aolo2nyoiy4Fvgd+p6mr8fxC2VNvna+CUOvoeBgwDaG9X+R0Vf/XoD7T5FEji2msvjsta7k63D6WP8vLpFBc/AkBu7ve0a5cb8xgaam9iJ+jTJBHJAaYAo1R1b43NS4F8VT0R+Csw9UCzWrqqtTq6qhaoak9V7dmyZctgwzJB0DisRw+QnNwN8KFqs0wFQ1UpKRnBnj3JFBcPAiAvby9JSYe7HJlxW1Bn9CKSij/Jv6Kqb9bcXj3xq+oMEXlWRFrgP4NvV23XtvjP+E0MxWM9eoDbb+/Pk08uYdGitfTp8xNXYvBS+4b6+OqrdRx99HMAXHTRRVx++eVcdlluTGMIhtWjj70Gz+hFRIBJwBpVfaKOfVoF9kNEegf63Ql8ChwtIh1FJA0YArztVPAmsV1yySUATJ061d1A4kRurj+pDxw4kKlTp3LZZZe5HJHximDO6PsBVwIrRWRZYN3dQHsAVX0OuBQYISIVQAkwJHC7T4WI3Ay8ByQDLwTG7o1pUL9+/cjOzubxxx8H7nI7HM874ogjAFi4cKHLkRivaTDRq+p8ah9rr77PBGBCHdtmADPCis40aiJC//79mTlzJjk5/yM5+SduhxQXdu7ciaoS+CfbGCuBYDygnjrq48ePB2D//j+G1d6JGGLS3qEY2iQn+5N8x47uxGC8qa4nqdxc7MlYZ3n6ydgG6qj7fD4FVKR9RHXYI4kh6u0djAHQM9yKIUj2ZGx0YPXoGzdPJ/og6qhfffXVgTIIC8Jq70QMUW3vYAyAXuJWDEGyRB8d9SV6G7ox7gqijvrFF18MQGnpI2G1dyKGqLZ3KIbvAy/z3IrBeJYleuOuemqgT58+naZNmzJ48GAAUlMvDKm9EzHEpL1DMRx4mKWTWzEYz7JEb9xVRw30LbffzqBBg9i9ezfHH388WVnvkZb2m6Dbe6EWfKxjyM7IAOBe/LfJNQNu69aNr776KjYxGO+qa0zHzcXG6J3l6TF61VproF9xxRUK6Pz588Nq70QMMW3vUAyftGqld4GemJpaNRcsoG+++WZsYgiCjdFHB/WM0Yt/u7f07NlTFy9e7HYYCePAI+fxUgJh2rRpDB48mMsuu4zXXnvNlRi82D7UPiorK5k9ezZXXXUVqsq2bdto1iwppjHUJtKfR1M7EVmiqj1r22ZDN8ZzHnvsMQAKCgpcjiS+JScn8/Of/5xx48axc+dOLrjgAny+7W6HZVxgid54zmmnnQbAjBn2QLUTrrrqKrKyspg1axaFhYezZ08Wl19+OTfccAO33347hYWFbodoosyGbhqBrl399ehzcsJrf6B+eJMm4ccQSh/79g2momIa2dlzSEk5y5UYvNg+kj58vq1UVLxPScnDwNqDtiUnn0xOTvC/b069Dx5MPXHNhm4auXj6hfL5fqCiYhqpqUOrkryJXFLSkaSlXYV/+ubZZGcvIC/PR3JyPyorl1Be/obbIZooCmmGKROf4qke/fr1uznqKLjttjY8+qg7MXi1vXMxnFPV/rvv3qB79+7s2DGE2bM/ol+/flGPwerRx56d0RtPKC0tZcuWH2edLC4udjGaxqNVq1YsXryYJk2a8POf/5ydO3e6HZKJAkv0xhOuueYa2rdvz1FHHQXYZNax1LZtW5555hkKCwuZNm2a2+GYKLChG+O6f//73/zzn/8EYMiQIfTp04drrrnG3aAamTMC4ynXX389nTp1qvrcJAZL9MZ1119/PQCbN2+mXbt2DextouGII45g8ODBTJ06lTPPPJMNGzbQoUMHt8MyDrGhm8bO5Qk3tm7dyubNm+mWmkq7/Py4nvQjnmMQEd566y3m3n03AL/t2NEmHkkkddVGcHOxWjfOqrO2SJgTTezevVu/+OILXf/EEzo7PV0ngO4A3Qv6bnq6PjZkiE6YMEFfeOEFXbhwoc6aNUszM/+hKSmDtG3btnrMMcdo3759dfHixXrd6acroK9HedIPJ/qIh4lHnIjhiEB9nPE28UhcwSYeadzq/MUKY6KJTz755KBiWU4sp4JWRnnSDyf6iIeJR5yIYVm1780Mm3gkbtSX6G2MvjELY6KJ+++/H4BzzjmHy95/nybAKvxjgNlAMnA+kL1lC1999RXLli2jWbNm3HTTcYi0Zu/eVgB8+umnzJo1i+b33stQahlDjLNJPxIphhPxPzvbBbgf+LlNPBL3Gkz0ItIOeBloBfiAAlV9qo59ewELgV+p6huBdSOBG/CXyH5eVcc7E7qJWPv2sGlT7etr8dlnn/Hee+/Rv39/Zs+e7R/D3bSJX9bcMT8f2ralbdu2VXdv3Hrrwbv06tWLXr16waRJIcUQ6TFEpY8EjOHAbFXHhtqH8aRgLsZWAKNVtRvQB7hJRI6puZOIJAOPAu9VW3cc/iTfG/+JwkAROdqJwI0DQphoQtXH+eefD8B9990XcnsnYohKe4vhkD62ABcHVv0+I8MmHkkEdY3p1LUA04Bzalk/CrgJeBG4NLDul8Dfq+3zR+COhr6GjdE7q94x0SAmmsjL82lS0nEK6OGHHx5ye1X3J/3wQgyxmHgk0hiWPvRQ1fj87MMPt4lH4ghOXYwFOgCbgbwa69sAH+Efoq2e6LvhH+5rDmQBC4C/1tH3MGAxsLh9+/axem8ahUh/sTIynlNAe/bsqZWVlWH1EWlyiPosV3HQPhYx3H777QroW2+9FbUYLNFHR32JPuj76EUkB5gCjFLVvTU2jwfuVNXKGv8trME/nPM+MAtYjn8o6BCqWqCqPVW1Z8uWLYMNy0TRl19+ycCBA9m/fzjQhPnz55OUZI9eJLKFCxfSpk2bqgnZTWII6q4bEUnFn+RfUdU3a9mlJ/BPEQFoAZwvIhWqOlVVJwGTAv08AnztSOQmaubNm0f//v0PWpeT8wHp6ekuRWSibdu2bZx11lmsXr2a008/3e1wjMOCuetG8CfqNar6RG37qGrHavu/CExX1amBzw9X1W0i0h64BDjVgbhNCLZu9U88cqC8bEP2759d9To7+wP27TuToqLg29fmwGQT4fYRaXsvxODlYyguvpny8tUALF58f739OxWDiZ1gzuj7AVcCK0VkWWDd3UB7AFV9roH2U0SkOVAO3KSqu8KM1YRJQ5x4JCVlIKWlD5Gaeh0pKWdGJyjjGRUViygvf520tGFkZv7N7XBMFDSY6FV1Pv574IOiqtfU+PynoYdlnBT6xCOn0LZtGw4//DOWLvXShBnxHYNXj+H++2fxpz/Bli0Pcvjh0Y/BCmPGnj0Zaw5RUVHBN998g8/nczsUEwNTpkyhadOm2E0QictuoTCHOPBA1IEHpEzi8vl8rFy5klNOOYXAzRQmAdkZvTlIeXk548ePB+Cpp2qtdGGC5PP5qm5H9f93dOh51TfffMPLL7/Mc889R8uWLZk4cSI9e/ZERNi7dy+TJ0+moqKC0tJSLr30UiDf0Rj/97//AfDTn9oIayKzM/rGrkYN872TJlFcXMzIkSPJzs4OuX281WF3OgYV4eeZmWRnZJCcnIyIICIkJyezZ08q5eUzufbaa2nXrh1HHnkkbdu25e6772bz5s0sWbKE3r17k5uUxGNNm9KkSRNGjBjBLbfcwu9+9zs6dOhAZeVXjh7DX//6VwBOOeUUZ94H4011PUnl5mIlEJwVSj16X2am5rdooW3btq3aLRFqwTsVQ07O5/rMM8/o8OHD9YorrtBbb71V3xo1qqr984HyAScnJemFJ52k1113nZ599tmakZFxSHnmLl266AknnKCvvvqq7i4o0AXp6dq1xj5X9OunmzZt0hdffFEBFTlcs7MX6X//+19dsGCBFhcXq6pqxcsva2lmZsjv4xFHHKGA+nw+Z2raB8GejI0OrB594xZKPfoVgQTTpk0b/y+/JkYteCf6yMh4ts6a+r8PtPl54PMva2mfm/u9ZmQ8rmPHjj20lES1GJaC/j/Q5TX6+O1vf1vn1xfQzMD3L9j3ce/evSoi2q1bt5DeB0v03mSJvpGr8xdL5JBf6m6BxDFv3ryq3er8xa6lvYJ/fQ2R9uF2DPfff39VUn3nnXd07969quqfbevkWhJvaSjHEMJxZGd/pOnpD+vrr7+uY8aM0aysrEO+9mjQkgbex7KyMj3ssMMU0HfffTekGCzRe1N9id7G6BuzGnXG1wSWe5s04bTTTgu5fYPro9FHlGPw+XwMHz6c++67D5GW5OZuY+DAgeTm5gLQpEkTZrVty/HAYfgnXVkHpDkYQ3UpKf3JyLibSy+9lAcffJB9+/bh8/nQ/HymBfb5C9ASeAfw1THZeo8ePdi9eze33nrrj3dXOfFeGk+yRN+Y1ahhviDw8ZQbbgirPeC5Ouw+n49Vq1b5/30NsY9Vw4fTrFkz/va3v9GzZ09yczeTlHToveYt/vxnVmRlsQt4F+gU6jEEcRz1ERF4+GEuzMqiEvgTUARcCCRv3sz111/P/v37q/YvL3+D1atX0759+6o7rCKNwXhcXaf6bi42dOOsYOvRP9G0qQI6f/78g3bxUi34+++/Xzt27Kjdu3fX/v3769FHH60ntG+vzzVrdlD7PXv26HHHHVc1nNG9e/eQYnjvzjur2vbp00fLy8vjqh791rZtdWjfvgcN6YwZM0bT0+9TQDt27KhFRUXROY4G2NBNdGBj9I1bsL9YK1euVEC7du1adSFW1f066pWVlZqd/YGKNKvzYiSggwYN0u3bt+t1111XtS4lJaXqdV7eLs3N3awffPCB7t69Wzds2KB9+vTRdu3aaU5Ojnbr1k27dOmiLVq0qGrzn//8x7H3wY330efz6ahRo2q8Vzn62WefxSyGmizRR4cl+kYulF+sjh07ViWEA2d8sU5Qa9eu1aeeekovv/xy7dev30FJqmPHjlpYWKglJSVaXl6uqqqlpaV65plnHrRfRkaGTpgwQVVVJ02aVO8fCECPPPJITUpKUkAHDBigI0aMCO0/G4ffA6f7qKio0Dlz5mh6+p2am7vFlRgOsEQfHfUlensy1hxk3bp1VU9z5uTkkJmZSUXF5Yik87Of/Y+8vDyaN2/OxIkTSU1Ndfzrl5eX06VLl4PWdevWjXXr+pOWdjvr13c5pE1aWhpz5szhd7/7HUuWLOH666/nyiuvrNp+1VVXsXjxYgoKlpCcfArXXltGy5YtKSoq4uyzz+aCCy5w/Di8Jjk5mbPOOouMjLPcDsW4QPx/CLylZ8+eunjxYrfDSBhdu/rr0efkBLd/Wdlkysqex+fbBCiqZUAp/kt8pVX7pab+irS0G0hOPhWRrDp68ztQg7xJk7r3UVX27Tudysp5pKf/gfT0PyGSFnT7hkTah9vtEy0GD6aeuCYiS1S1Z23b7K6bRiDUX6i0tMvJyZlLXt568vI2AN8AO2jSZD9NmihpaaMAKC9/jX37zmbv3sMoK3s5whh3UVTUk8rKeUAO6en3VSV5Y0xkbOimEQi9Hv3BDq0//iTwJKtXr+b555/nqaeeoqTkan772xWMGzcuyD5+VFpaSk7O4fh8FQwcOJC33377kEqKXq3lHsv2iRKD1aOPPTujN2E79thjGT9+PPv27aN379785S9/4ZZbbgmpj7Vr1zJ48GAqKioYPXo077zzjpXLNcZhluhNxLKyspg1axZt2rRhwoQJjBkzJqh2N9xwA127dmXWrFkA3HPPPdEM05hGyxK9cUTTpk1ZsMD/bO3DDz9MSUlJnfvu3buX5s2b8/e//52kpCSmTZtGRUUFTZs2jVW4xjQqluiNY9q1a0dBQQHgP8t/4403KC8vP2ifDRs20KVLF3744QdOO+009u3bx4UXXkhycrIbIRvTKFiiN5GpMVHFtenpnHzyyQD88pe/JC0tjaFDh1JZuYZ9+86hU6dOfP/999x1113MmzePjIwMm+zCmChrMNGLSDsRmSsia0RktYiMrGffXiJSKSKXVlt3W6DdKhF5VUQynAreuOyVV2DYMNi0yX8P56ZNpIwYweLbbuObb75hxIgRdOnShcmTJ1NUdAwVFf+mU6dOTJ8+nbFjx9bZB8OGWbI3xkHBnNFXAKNVtRvQB7hJRI6puZOIJAOPAu9VW9cGuBXoqarHAcnAECcCNx5wzz1QXHzwuuJiuOceWrduzbPPPssXX3zBHXfcQUrKuWRnz2HdunUHP4laTx/GGGc0mOhVdauqLg28LsRfsrxNLbveAkwBttVYnwJkikgKkAV8G1HExjs2bw5q/aOPPkp29nukpNTy+H2QfRhjwhfSGL2IdAB6AItqrG8DXAw8V329qn4DjAM2A1uBPao6u46+h4nIYhFZvH379lDCMm7xwsQjxpgGBZ3oRSQH/xn7KFXdW2PzeOBOVa2s0aYpcBHQEWgNZIvIFbX1r6oFqtpTVXu2bHno5A7Gg7ww8YgxpkFBlUAQkVT8Sf4VVX2zll16Av8MPNHYAjhfRCqAVGCDqm4P9PMm0Bf4hwOxG7cNHer/eM89/qGW9u39CfrA+lj1YYypV4OJXvzZexKwRlWfqG0fVe1Ybf8XgemqOlVETgH6iL+0YQnwM8DKUiaSoUMjT8pO9GGMqVMwZ/T9gCuBlSKyLLDubqA9gKo+V0c7VHWRiLwBLMV/985nQEEkARtjjAlNg4leVecDQVeZUtVranx+H3BfyJEZY4xxhJUpbiSWLQu/PGxRkf9jJOVlI+0jEWJIhGNwoo9ly6B79/C/vgmdJfpG4PLLI2vfqRNEWjk40j4SIYZEOAYn+ujePfKfSRMaS/SNwLBh/sUY0zhZUTNjjElwluiNMSbBWaI3xpgEZ4neGGMSnCV6Y4xJcJbojTEmwVmiN8aYBGeJ3hhjEpyoqtsxHEJEtgObwmzeAtjhYDheZ8eb2Ox4E5fTx5qvqrVO5uHJRB8JEVmsqj3djiNW7HgTmx1v4orlsdrQjTHGJDhL9MYYk+ASMdE3tolN7HgTmx1v4orZsSbcGL0xxpiDJeIZvTHGmGos0RtjTIKLm0QvIu1EZK6IrBGR1SIyspZ9morIWyKyQkT+KyLHVds2QES+EJGvROSu2EYfmkiONZi2XhPp9zawPVlEPhOR6bGLPDwO/CwfJiJviMj/An2cGtsjCI0Dx3tboN0qEXlVRDJiewShEZGMwDEsD8R9fy37iIg8HchHK0TkpGrbnM9VqhoXC3AkcFLgdS6wFjimxj6PA/cFXv8EmBN4nQysAzoBacDymm29tER4rA229doSyfFW2347MBmY7vbxRPt4gZeA3wRepwGHuX1M0TpeoA2wAcgMfP4v4Bq3j6mB4xUgJ/A6FVgE9Kmxz/nAzMC+fYBFgfVRyVVxc0avqltVdWngdSGwBv8PQXXHAHMC+/wP6CAiRwC9ga9Udb2qlgH/BC6KWfAhiuRYg2zrKRF+bxGRtsAFwN9jFnQEIjleEckD+gOTAtvKVHV3rGIPR6TfX/xTnmaKSAqQBXwbk8DDpH6BKdRJDSw173q5CHg5sO9C4DAROZIo5aq4SfTViUgHoAf+v5TVLQcuCezTG8gH2uL/odpSbb+v8XjyOyCMYw2mrWeFebzjgTsAX0yCdFAYx9sJ2A78X2Co6u8ikh27iCMT6vGq6jfAOGAzsBXYo6qzYxZwmAJDicuAbcD7qlrzeOvKSVHJVXGX6EUkB5gCjFLVvTU2/xloGniDbwE+Ayrw/3tUk+fvKw3zWINp60nhHK+IDAS2qeqSmAbrgDC/vynAScBEVe0B7AM8fc3pgDC/v03xn9F2BFoD2SJyReyiDo+qVqpqd/x/nHvXvKZE3TkpKrkqJdIOYklEUvH/oLyiqm/W3B744bk2sK/gH9vbgP/fvXbVdm2Lx//9i+BYG2zrRREc7xDgQhE5H8gA8kTkH6rq6WQQ4c/y19XOEN8gDhJ9BMd7HrBBVbcHtr0J9AX+EaPQI6Kqu0XkQ2AAsKrapq+pPSel1bE+InFzRh/45k8C1qjqE3Xsc5iIpAU+/Q3wn8AP0KfA0SLSMbB9CPB2LOIORyTHGkxbr4nkeFX1D6raVlU74P++fhAHST6S4/0O2CIiXQPbfgZ8HvWgIxDh7+5moI+IZAX6+Rn+MX7PEpGWInJY4HUmcDbwvxq7vQ1cFbj7pg/+IamtRClXxdMZfT/gSmBl4N87gLuB9gCq+hzQDXhZRCrx//BfH9hWISI3A+/hv6r9gqqujm34IQn7WOtqq6ozYhN6WCI53ngU6fHeArwSSATrCZwJe1gkv7uLROQNYCn+oavP8H6ZhCOBl0QkGf/J9L9UdbqIDIeq452B/86br4BiAt/DaOUqK4FgjDEJLm6GbowxxoTHEr0xxiQ4S/TGGJPgLNEbY0yCs0RvjDEuE5EXRGSbiKxqeG8QkctE5HPxF02b3OD+dteNMca4S0T6A0X469/UfIq25r5H4y/udpaq7hKRw1V1W31t7IzeGGNcpqr/AX6ovk5EjhKRWSKyRETmichPAptuAJ5R1V2BtvUmebBEb4wxXlUA3KKqJwO/A54NrO8CdBGRj0VkoYgMaKijeHoy1hhjGoVAAbi+wOv+yg8ApAc+pgBHA2fgr4UzT0SOq69ctSV6Y4zxniRgd6ACZk1fAwtVtRzYICJf4E/8n9bXmTHGGA8JFHTbICK/hKqpB08MbJ4KnBlY3wL/UM76+vqzRG+MMS4TkVeBBUBXEflaRK4HhgLXi8hyYDU/zjT1HrBTRD4H5gK/V9Wd9fZvt1caY0xiszN6Y4xJcJbojTEmwVmiN8aYBGeJ3hhjEpwlemOMSXCW6I0xJsFZojfGmAT3/wFuDT9lTJ30JwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# ERA5 total precipitation\n",
    "era5_fileName = r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\era5\\total_precipitation-19790101_20181231-eusalp-qm_era5.nc'\n",
    "\n",
    "tp = readnetcdf_in_shp(era5_fileName, shp_catchment_fileName,plot=True)['tp']\n",
    "tp = xarray2df(tp.resample(time='1d').sum(skipna=False), 'P','tp')\n",
    "\n",
    "check_data_gap(tp)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "8242be04",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Concatenate and save\n",
    "daily_input = pd.concat([runoff, t2m, tp, pet], axis=1, join='inner')\n",
    "path=r'C:\\Users\\mmazzolini\\OneDrive - Scientific Network South Tyrol\\Documents\\conda\\daily_input\\\\'+catchment_name\n",
    "#save into a csv file.\n",
    "# if file does not exist write header \n",
    "if not os.path.isfile(path):\n",
    "   daily_input.to_csv(path)\n",
    "else: # else print the problem\n",
    "   print('file already exists')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1acfe16e",
   "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.9.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}