scripts/OAEaddition.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "661f0272-c0b9-404c-be85-91cd733d66ad",
   "metadata": {},
   "outputs": [],
   "source": [
    "import netCDF4 as nc\n",
    "from netCDF4 import Dataset\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.patches as mpatches\n",
    "from matplotlib.patches import Rectangle\n",
    "from matplotlib.legend_handler import HandlerPatch\n",
    "import xarray as xr\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import seaborn_image as isns\n",
    "import cmasher as cm\n",
    "import warnings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "0b0a6027-37b8-4dfa-8909-b04a2a6b2f0a",
   "metadata": {},
   "outputs": [],
   "source": [
    "warnings.filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "1a7451ac-2bad-4d38-b2e6-5d05326e3c8f",
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.style.use('seaborn-pastel')\n",
    "isns.set_context(mode=\"notebook\", fontfamily=\"MS PMincho\")\n",
    "sns.set_context(font_scale=0.8)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "f564f0ce-1d8f-43ae-af7b-306cf5e06534",
   "metadata": {},
   "outputs": [],
   "source": [
    "alk_mask = xr.open_dataset('/Volumes/ccHD/OAEdata/variables/masks/alkalinityMask.nc')\n",
    "mesh_mask = xr.open_dataset('/Volumes/ccHD/OAEdata/variables/masks/meshMask.nc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "d5bdb542-480a-43db-ae7e-195fce29e82e",
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh_mask = mesh_mask.tmask.isel(t=0, z=0)\n",
    "mask_alkalinity = mesh_mask + alk_mask.alk_flux\n",
    "mesh_mask_region = mask_alkalinity.isel(x=slice(520,595),y=slice(330,435))\n",
    "mesh_mask_study = mask_alkalinity.where((-40 < mask_alkalinity.nav_lon) & (mask_alkalinity.nav_lon < 40) \n",
    "                                                         & (30 < mask_alkalinity.nav_lat) & (mask_alkalinity.nav_lat < 90), drop=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "f6dedd9c-1c4c-4660-95e4-83157e1393ea",
   "metadata": {},
   "outputs": [],
   "source": [
    "mesh_mask_study = mesh_mask_study.where(mesh_mask_study)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "cf99d7a2-ebdd-4b19-ae10-e6d715b85b73",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAngAAAJaCAYAAAC88zoHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAA9hAAAPYQGoP6dpAABBCklEQVR4nO3de3RU9b3//9ckhBgICUYiiIBAJIFypKDcRMDWhEu1P4vFg1moBQ/xKPCNeOHaCgSUoKggtkA8R5NasUePrVIRxVRtQSkXlXpAQMDUQEQQBEkChJDL/v2hjpNAhsnMnsyezzwfa7FWZmdn789kM8Ob92s+n+2yLMsSAAAAjBEV6gEAAADAXhR4AAAAhqHAAwAAMAwFHgAAgGEo8AAAAAxDgQcAAGAYCjwAAADDUOABAAAYhgIPAADAMM1CPYBw98tf/lKdO3cO9TBgsy37v/Bpv/6dOgR5JAAA/KC4uFivvPLKefejwAtQ586dtXjx4lAPAzZ75O31Pu03M2NokEcCAMAP7r//fp/2o8ADAlC/EKTgAwA4AZ/BAwAAMAwFHgAAgGGIaIHv+Pq5O3+PQXwLAGgqdPAAAAAMQ4EHAABgGCJaoIl4xrfEtQCAYKKDBwAAYBgKPAAAAMMQ0QIhwGxbAEAw0cEDAAAwDB08RCw71r0LBiZjAAACRQcPAADAMBR4AAAAhiGiBRyMyRgAAH/QwQMAADAMBR4AAIBhiGgRUZw6c9YfzLYFADSEDh4AAIBhKPAAAAAMQ0QLGIDZtgAAT3TwAAAADEMHDzAckzEAIPLQwQMAADAMBR4AAIBhiGhhPJPWvgsUkzEAIDLQwQMAADAMBR4AAIBhiGgBSCK+BQCT0MEDAAAwDAUeAACAYYhoYRxmzdqPxZIBILzQwQMAADAMHTwAjcJkDABwPjp4AAAAhqHAAwAAMAwRLQDb1I9viWwBIDTo4AEAABiGAg8AAMAwRLQwAmvfOVND14XoFgCCiw4eAACAYejgAWhyrKUHAMFFBw8AAMAwFHgAAACGIaJFWGJShblYSw8AAkcHDwAAwDAUeAAAAIYhogXgaKylBwCNRwcPAADAMBR4AAAAhgnbiLa0tFSrV6+Wy+VSv379lJqaGuohAWhCzLYFgIaFbYH30ksvKSsrS1FRUcrLy6PAAwAA+I5jCry33npLu3btUps2bXTgwAGlp6erb9++mj59utav/+F/6r1791ZeXp4sy5JlWZKk2tpan89TVFSksWPHavTo0Zo+fbp7+7FjxzRlyhRNnjxZAwcOtO+JwTasfQdvmIwBAD9wRIF35swZLV++XH/5y18kfVuwjRw5UoWFherQoYM2bdp01s9kZmZq5cqVkqSMjAyfz5WSkqJJkyYpPz9fw4YNU58+fSRJSUlJGjFiBMUdAAAIe44o8E6dOqWSkhJVV1erWbNmioqKksvlkiRFRZ17HkhiYqLGjRvn1/mio6OVn5+vsWPH6u9//7vi4uL8HjsAAIDTOKLAa926ta6//noNGzZMixcv1meffaasrCxJUnl5uebMmaPq6mrt2bNHCxYsUFpaWsDnTElJ0cSJEzV16lQtW7Ys4OMBcCZv0T7xLQBTOaLAk6Sbb75ZFRUVuuuuu9SyZUu9+uqrkqT4+HjNmDFDUVFR2rZtm7Kzs1VYWGjLOcePH6833nhDa9as0Q033HDe/SsrK1VZWVlnW01NjS1jAQAAsIsj1sE7cuSIHnnkET3xxBPavHmzhg4dqn//93+XJGVnZ7tj2l69emnv3r22njsvL09z5szR4cOHz7vvwoULlZiYWOfPli1bbB0PAABAoBxR4K1bt05Dh34blbhcLs2bN08HDx7Uvn37NGHCBPdsWUmKjY219dxJSUlatGiROxL2ZtasWSotLa3zp3///raOB3U98vZ69x/Abp5/v/g7BsAkjijwevToof3797sfV1dXKz4+XnFxcdq8ebN7wsW2bds0aNAg28+fnp6u1NRUFRQUeN0vNjZWCQkJdf5ER0fbPh4AAIBAOOIzeD179tTgwYP19NNPKz4+Xvv27dPvfvc7XXzxxVq8eLEKCgoUHR2tHTt2aNGiRQGda//+/Vq7dq369u2r7t27u7fn5uZqzJgxgT4VAACAkHNEgSdJP//5z8+5ffjw4baep1OnTu718zw1b95cq1atsvVcAMILiyUDMIUjIloAAADYxzEdPABwKtbSAxBu6OABAAAYhgIPAADAMES0cBTWIkO4qf93lsgWgBPQwQMAADAMBR4AAIBhiGgBwEaspQfACejgAQAAGIYCDwAAwDBEtHAUO2IsZuLCiZhtC6Ap0cEDAAAwDB08GMdbZ4TuHpyC258BCCY6eAAAAIahwAMAADAMES0AOAxr6QEIFB08AAAAw1DgAQAAGIaIFhHFM+JiRm348TWiNPXaspYeAF/RwQMAADAMBR4AAIBhiGgRserHW6bGeuGOGLJhLJYMoCF08AAAAAxDBw9AyNFtsp9nd4/fLxB56OABAAAYhgIPAADAMES0wHdYIw+mYjIGEHno4AEAABiGAg8AAMAwRLTAObBGHiIFs20BM9HBAwAAMAwdPMAH/nQ26Po1LBidIn7fgWMyBmAOOngAAACGocADAAAwDBEtgCZBxBfemIwBhBc6eAAAAIahwAMAADAMES1go0ifyRnK6K6hc0f6NQkGZtsCzkcHDwAAwDAUeAAAAIYhogVgNH8jQ6Jd/zDbFnAGOngAAACGoYMHBIAuT131fx+h6uAwCcAZuA5A6NDBAwAAMAwFHgAAgGGIaAEEDR+4R0OcEucDpqKDBwAAYBgKPAAAAMMQ0QKNxMxZ/zglruX6OVND14XoFvAPHTwAAADDUOABAAAYhogWQJML9gK4nscgkg1vzLYF/EMHDwAAwDB08AAf0AVqOnRs4I1TJusATkcHDwAAwDAUeAAAAIYhogXgaKyPhoYEe7IOEM7o4AEAABiGAg8AAMAwRLQAwpKv8Vz9qI4Z0ZGB2baIdHTwAAAADEOBBwAAYBgiWgDGIYaFJxbPRiSigwcAAGAYOngAgIjC2oqIBHTwAAAADEOBBwAAYBgiWgAAxGQMmIUOHgAAgGEo8AAAAAxDRAsAwDlwuzOEMzp4AAAAhqHAAwAAMAwRLQAA5+Ht9nfEt3AiOngAAACGoYMHAEAAmIwBJ6KDBwAAYBgKPAAAAMMQ0QI+8IxdvH3YGkDkIZaFE9HBAwAAMAwFHgAAgGGIaAEAEYVIFZGADh4AAIBh6OABAM5rx2t7fdqv543dgjwS39ClQ6SjgwcAAGAYCjwAAADDENECACT5HsM+/9SEgM/l63qSRK2Af+jgAQAAGIYCDwAAwDBEtAAQwZoylvUVsSwQODp4AAAAhqHAAwAAMAwRLeADX2f8AeGMGNaZPN9/+L3BV3TwAAAADEMHDwAiiK+TKhA63hKD+t+jo4eG0MEDAAAwDAUeAACAYYhoAcBw3mLZppxYAaDphG2BV1paqtWrV8vlcqlfv35KTU0N9ZAAAAAcIWwLvJdeeklZWVmKiopSXl4eBR4AAMB3HFPgXX311SovL1d8fLwk6csvv9TmzZtlWZZyc3PVo0cPFRUVKScnRwkJCbIsS5ZlSZJqa2t9Pk9RUZHGjh2r0aNHa/r06e7tx44d05QpUzR58mQNHDjQ3icHo3jOWmN9PIQbIlnnsOP9o6FjMLsWjijwLMvSyJEjNXfuXEnSiRMnlJ+fr0suuUQ33HCDFi9erLS0NG3YsEGzZ8/W0qVLlZmZqZUrV0qSMjIyfD5XSkqKJk2apPz8fA0bNkx9+vSRJCUlJWnEiBEUdwAAIOw5YhZtbW2t7r33Xvfjp59+WnfccYdKSkpUUlKitLQ0Sd92+VatWqXKykolJiZq3LhxGjduXKPj2ejoaOXn5+vuu+9WRUWFnU8FAAAg5BzRwYuOjlZiYqIk6dChQ6qoqFCrVq20ceNGdenSxb1fVFSUWrdurZKSEl1++eUBnTMlJUUTJ07U1KlTtWzZMp9+prKyUpWVlXW21dTUBDQOhB9iWQC+4v0CoeKIDp6nxx9/XDfddJMk6ejRo2rVqlWd77dq1Upff/21LecaP368jhw5ojVr1vi0/8KFC5WYmFjnz5YtW2wZCwAAgF0c0cH73unTp7VmzRo9/vjjkqTk5GSVlZXV2ae0tFTJycm2nTMvL0/Dhg1Tv379zrvvrFmzdP/999fZNnv2bNvGAgB24ZZkTYcuHZzIUR28DRs2qHXr1u7HPXr0UHFxsftxVVWVysrK1LFjR9vOmZSUpEWLFikrK+u8+8bGxiohIaHOn+joaNvGAgAAYAdHFXjbtm1TbGys+/Gll16qbt26aefOnZKkd999V5mZmWrevLmt501PT1dqaqoKCgpsPS4AAEAoOCqiTUhIULt27epsW758uXJycpSWlqYvvvhC8+fPD+gc+/fv19q1a9W3b191797dvT03N1djxowJ6NgwF2vfwemIZCMD69vBV44q8CZMmKAJE+ouwtm2bVutWLHCtnN06tTJvX6ep+bNm2vVqlW2nQcAACBUHBXRAgAAIHCO6uABTkUsi3DDLcmajq+xqb/vI8Sy8AcdPAAAAMNQ4AEAABiGiBY4ByJZhANmzgJoCB08AAAAw9DBAwADMKnC+epPlvA1KfAnUWBiBujgAQAAGIYCDwAAwDBEtMB3mFgBoClxC0QEEx08AAAAw1DgAQAAGIaIFhGLSAThhnXvzMWsV9iNDh4AAIBhKPAAAAAMQ0QLAA7mLZZlcWMADaGDBwAAYBg6eAAQJujYAfAVHTwAAADDUOABAAAYhgIPAADAMBR4AAAAhqHAAwAAMAwFHgAAgGEo8AAAAAxDgQcAAGAYCjwAAADDUOABAAAYhluVIWLNzBha5/Ejb68P0UgA39x+z7N1HnPrMgANoYMHAABgGAo8AAAAwxDRAoCD9byxm/vrHa/trfM9z8jW1Li2MR+dqP+xCyCS0cEDAAAwDAUeAACAYYhoYYuK0hJVVRwN9TACEnfys1APAfCqb7qrzuN/rdvv/rrs0MdNPBp7xMRdpLjEjnW2MaMdCBwFHgJWUVqiDf99lWqqToV6KAHpHuoBAI3UPfmHrzcWPB66gQQgOqaFrrnzo7OKPACBocBDwKoqjqqm6pSu+P+eVXybtFAPx28Fm7eGeghAo3h28B6aPip0A/HTia93a/vqCaqqOKqlH3we6uEARqHAg23i26QpoV3vUA/DbxUty0I9BKBRjlVb7q/D+bUHwH5MsgAAADAMHTxELD7IDQAwFR08AAAAw1DgAQAAGIaIFgDCRP1blQFAQ+jgAQAAGIYCDwAAwDBEtIgozJyFSZ5/akKohwDAoejgAQAAGIYOHoxH1w4AEGno4AEAABiGAg8AAMAwRLQwDpEsEH5mZgx1f81rGAgcHTwAAADDUOABAAAYhogWRiDSASKPZ6wLoC46eAAAAIahwAMAADAMES3CEpEsAAANo4MHAABgGDp4CDm6cQAA2IsOHgAAgGEo8AAAAAxDRIuAPPL2esWd/EzdJRVs3qqKlmWhHhKAMMRHNQB70cEDAAAwDAUeAACAYWyJaHfu3KkPPvhArVq10i9/+UsdP35crVu3tuPQcCCiFAAAnC2gDt6pU6d0ww03qH///nrsscf04osvSpKOHz+u2bNnq7y83JZBAgAAwHcBdfBmzpypoUOH6uWXX1aLFi306KOPSpI6d+6sadOmKTc3VwsXLrRloAgdOnYA7FZ/gpZaXh7qIQFGCajAO3HihGbMmOF+7HK53F8nJCTQwQMAAAiBgCLamJiYOo8ty6rz+MCBA4EcHgAAAH4IqMA7efKkdu7c6X7s2cHbtGmTmjVjmT0AAICmFlAFNn/+fI0cOVI//elP9eMf/1gff/yxCgoKtHXrVr3++usqLCy0a5wAAADwUUAdvK5du2rjxo1KTk7WypUr9eGHH2rFihWKiYnRe++9p27dutk1TgAAAPgo4Az1oosu0sMPP2zHWOAgzJwNrh2v7W3wez1v5D9GgC/qv0/NzBgaopH4x9f32XB7XnCGoN7JYuXKlcE8PAAAAM7B5w7eO++8oxdeeMHnA585c0bFxcW67bbb/BoYAAAA/ONzgXfFFVdow4YNuuWWWxQdHS1J+vrrr7Vz504NHVq3fXz48GG9+eabuvPOO+0dLYImUiJZu6NRb8fzl+cxiWsRjL9jpvJ8H3NirBkp77NwBp8LvIsvvlizZs3S+PHj3dseeOABvfnmm2revPlZ+0+aNEl/+tOfbBkkAAAAfNeoSRaexZ0klZWVnbO4k6SePXtq2bJlfg8MwRGJ/4P0tQPia+fMjo7K809NaPB7t9/zbKPHhMjg7e8NnMmf91wndh8RfgKaZHHo0CGv3y8uLg7k8AAAAPBDQAXe4MGDdeutt2rXrl11tpeWlur+++/Xj370o4AGBwAAgMYLaB28GTNm6MEHH9RVV12l5s2bKzk5WWfOnNFXX32lrKwsLV261K5xwk9EsnX5E416E4zIzPOY3sZEZGsuJlaEr0h8z4UzBbzQ8cMPP6xp06Zp48aNOnr0qBITE9WvXz+1bdvWjvEBAACgkQIu8CQpMTFRI0eOtONQAAAACJAtBV5DlixZovvuuy+YpzAO7X3/+BvLNrSfZzTq7/Hs4G1MzLA1hx1/fxE6vG/DiQIq8K677roGv1dVVaUTJ05Q4AEAADSxgAq8M2fOKCsrSy6Xq872o0eP6oUXXtDkyZMDGhwAAAAaL6AC77HHHtPVV199zu9df/31+sMf/hDI4YGQcWIsVn9MzLA1lxP//oWbYN+2jAWM4XQBrYPXUHEnSd27d9e2bdsCOTwAAAD8ELRJFlVVVSoqKgrW4Y3Bh3Pt4dmxqt/N8ux00RmBU0T6WnczM4aq7FCCNu6U7hhwpZZ/Uha0c3l7n/W1q+brezVdOjhFQAXe/Pnzz7n99OnTevfdd9W/f/9ADg8AAAA/BFTg/fd//7fS09PPmmSRkJCgu+66S7feemtAgwMAAEDjBVTgLViwQL/61a/sGgsAGIu17pzJjo/JEMvCiQKaZEFxBwAA4DwBFXhffPFFg9979tlntXXr1kAODwAAAD8EFNH+8Y9/1PTp08/5vQkTJmjs2LH64x//GMgpgEarvwacZzQWzjNqvd0+jXXvwk+4/f0LNn9iTlYhABoWUAfPm7KyMu3evTtYhwcAAEADGt3B+8c//qHf//732rp1qz7//HPl5eWdtU9NTY0OHz6sX//617YMEgAAAL5rdIE3aNAgDRo0SOvXr9cLL7ygsWPHnrVPTEyMunTpoksuucSWQZ5LaWmpVq9eLZfLpX79+ik1NTVo50J4a2gR5PqRpxMjM2LZ8BbpixmbjJmzcDq/P4M3dOhQ7dmzR9dee62d4/HZSy+9pKysLEVFRSkvL48CDwAA4DsBTbLIysry+v1//etf6tq1q0/HOnPmjFasWKELL7xQ27dvV1pamrKysjR9+nStX//DB2l79+6tvLw8WZYly7IkSbW1tT6PuaioSGPHjtXo0aPrTBA5duyYpkyZosmTJ2vgwIE+Hw+wW0NdOzp24c2JHeJwV7+LFsxJF3TsEG4aVeDl5eUpKSlJY8aM0cmTJ3X06FGv+y9YsEDPPttwxORp1qxZmjp1qi655BLt3r1b77//viSpQ4cO2rRp01n7Z2ZmauXKlZKkjIwMn59DSkqKJk2apPz8fA0bNkx9+vSRJCUlJWnEiBEUdwAAIOw1qsCbN2+eEhMTNWbMGH344Yf66U9/KpfL5e6k1edyuXwq8EpLS/XPf/7T/Zm9tLQ0paWlSZKios490TcxMVHjxo1rzPDdoqOjlZ+fr7Fjx+rvf/+74uLi/DoOAACAEzWqwFu9erVatWolSerbt68mT56s3/72tw3uf//99/t03A8//FAXXXSRHn74YZ0+fVqff/65/uM//kPp6ekqLy/XnDlzVF1drT179mjBggXu4i8QKSkpmjhxoqZOnaply5b59DOVlZWqrKyss62mpibgsSC4GvqgeygjM2+TJzwRywIA/NGoAq9v377ur1u2bKkHH3zQ6/4zZ8706bhfffWVNm7cqKVLl6p9+/aqqKjQ4MGD9Y9//EPx8fGaMWOGoqKitG3bNmVnZ6uwsLAxw27Q+PHj9cYbb2jNmjW64YYbzrv/woULNW/evDrbiHQBAIDTBLTQcdu2bRv83qlTp3TxxRf7dJyWLVtq0KBBat++vSQpLi5Oqamp2rFjh7Kzs90xba9evbR3r73LDuTl5WnOnDk6fPjwefedNWuWSktL6/zp37+/reMBAAAIVECzaN9880397Gc/O+f3Vq9erQMHDvgU03bp0kUHDhyosy0qKkoxMTGaMGGCnnnmGblcLklSbGxsIEM+S1JSkhYtWqSsrCyNGTPG676xsbFnnT86OtrW8cB+vq6D5w9vMa+/xyeWBcKLt9m7zL5FqATUwdu+fXuD37vlllu0ceNGn47Tq1cvVVdX67PPPnNvq66u1sUXX6zNmze7i7tt27Zp0KBBgQz5nNLT05WamqqCggLbjw0AANDUGt3B++qrr/TVV19J+nadO8816r5XU1OjPXv21CnYzucPf/iDHnroIQ0ePFg1NTWaP3++2rZtq8WLF6ugoEDR0dHasWOHFi1a1Ngh17F//36tXbtWffv2Vffu3d3bc3Nzz9vBAwAACAeNLvBOnjyp999/XwUFBSoqKtJbb7111j4xMTHq2rWrnn76aZ+Pm5aWpueee+6s7cOHD2/sEL3q1KmTe/08T82bN9eqVatsPVdDgrkYJ7xrKK71FzFsXb7+Tk19/t6Ew63xYD/P93viWjSlRhd4Xbt21aRJkzR27Fg9/fTTmjFjRjDGBQAAAD/5PcmidevW+slPfmLjUICmZUcXyVvHKlK6VN5+B55dKs8OVv2fMfV3ZXfHGHU1ZRpC8oJwE9AkiwEDBjT4vRdeeOGckSsAAACCK6BlUqRvZ7Zu3779rDs6WJal5557zu/biQEAAMA/ARV4Tz/9tB577DF16NBBJSUl6tSpk/t7JSUlmjhxYsADBJzM1GixMbzFkA1NQuH3BgDBFVCBV1RU5F4KZf78+ZozZ477e0ePHj3nbFUAAAAEV0CfwUtOTnZ/HRsbq+PHj7sfX3TRRSopKQnk8AAAAPBDQB28Y8eOSfr283bDhw/Xk08+qZycHPf3y8rKAhocgPBSP3r1jGyJZevyjK9ZE893oZrN6usadsy2hVMEVOBdddVViouL0y9+8Qu9+OKLWrBggUaMGKGrr75an376qaqrq+0aJwAAAHwUUIF38803q2fPnmrRooUk6eWXX9by5cv17rvvqn379po3b54tgwQQnuja/cBbd5O7XJiDu1XAKQJeJqVHjx7ur10ul+6++24dPHhQCQkJ2rVrl6655ppATwEAAIBGCGiSxblER0fr4YcfVlxcnIYO5X8yAAAATS3gDl5DsrOz9dZbbwXr8ADQIF9vCxbKCJnbmAWXZ1TKxAdEIts7eJ6IZwEAAJpeUAu82NjYYB4eAAAA5xC0iFb6dn08AAil+jNSPWes1o9GI3HWb0PxpVNng/oTvdZ/LkS2iASN6uBt3bq1UQd3uVyN2h8AAACBa1QHb+HChZo+fbrPnblNmzb5NSgAAAD4r1EF3p///Gf9+c9/9nl/OngAQq3+IsKenBrJBvM2Zr7Gk/X3c2pk2xBiWES6RkW0U6ZMUW1trc9/7r333iANGwAAAA1pVAcvJyenUQefNm1ao/Y3Gf+bBJqOtzXmnNi183YbM6fwfA8Lt24eEIka1cFLTExs1MHbtWvXqP0BAAAQuKCugwcAAICmF9R18ACgqTgx1kRwNeX6duE+6QSRhw4eAACAYSjwAAAADENECyAseYtkPdeOq78OnufPOXFGbX3e1vHzdY08u9cCJK4EnI8OHgAAgGEo8AAAAAxDRAvACN5iWU/hEMt6W6jZk7fn6cuxvWnMDNVIXAQ5Ep8zwgsdPAAAAMPQwQPgaJG+vp23jpuvvxt/u3b+aMrOFreABBpGBw8AAMAwFHgAAACGIaIF0OT8jV19XfctUvgavfq6ZqDd7Fgvz6kxLBMr4HR08AAAAAxDgQcAAGAYIloATcLfW4v5qqGfC4d175qSv5GsZyQZ7NjUqbEsEE7o4AEAABiGAg8AAMAwRLQAmpy3mNDuWZ2NmbFLnGsvk2aamvRcEBno4AEAABiGDh6AJld/QoQdXTt/jlF/HJ7dPrp5/gmHTldDY2RyB0xCBw8AAMAwFHgAAACGIaIF0CQ8I8/6Ex88o1JuRxZcwYjHfVU/Ao07+Zm6SyrYvFUVLcuCdl5fY+NwiJcBX9HBAwAAMAwFHgAAgGGIaIOE2VhAw+rPUG3MWnWB8nYrNFNnznqLx/1RP8r0fL9zynsfcSsiHR08AAAAw1DgAQAAGIaIFgAiGDOYATPRwQMAADAMHTwAiCB2THCxYyLFzIyhKjuUoI07pTsGXKmEdr0bfXwmUgANo4MHAABgGAo8AAAAwxDRAnCUYN9Ky/N49c/lGVeauiaeN06ZcEH0CgSODh4AAIBhKPAAAAAMQ0QLIOS83UormLFh/eN5nitS4lp/bmNGhAo4Hx08AAAAw9DBA+Ao3tZpqz8pwpMd3b2GJmDU72yZ3NEDYAY6eAAAAIahwAMAADAMES2AoPHnNljhIJwnYJh6TQDURQcPAADAMBR4AAAAhiGiBdAk7Jjl6m0Wrd28jTfcZth6i2VDeUsyAMFDBw8AAMAwFHgAAACGIaIFYISG4ttgRJANLYgsNe0MW39mxBLJApGBDh4AAIBh6OABaBSnrKPmrTvWlLc3cwqTnguAwNHBAwAAMAwFHgAAgGGIaAEExInRoK/xrR3qP/9grtXnlHgcgPPRwQMAADAMBR4AAIBhiGgBGMGJ8aWvY/J3vTwnxuMAnIEOHgAAgGEo8AAgAlVXV6tDhw7avn27Lcf7/PPP9cQTT2jlypWaOXOmLrjgAp05c8aWYwNoPCJaAAHxnDXqlMgwlOPw9dzBnG3ri2bNmqmwsFDdu3c/776FhYVKTU1V586dG9xnxowZeu655xQXFydJ6ty5sw4dOqROnTrZNWQAjUCBBwAO0O7AB+q261W1PHFIJ+PbaW+Pm3To0n5BO19NTY0sy1JUlPcgp7a2VqtWrdL06dO97nfw4EHV1ta6H9999922jBOAfyjwbPTI2+tDPQQg6OpPCHDi5AYndhW9aXfgA121+SlZcsklS63KvtBVm5/SRwPuOW+R9/777+uWW27R448/rtOnT6u4uFhdunTR+PHjJUkff/yx3nnnHaWmpuqTTz7RLbfcoq5du+q9997TPffco23btulf//qXMjMzdc899+j06dM6cOCAOnXqpDvuuEMbNmzQp59+qsLCQnXt2lUZGRnnHMfEiRM1ZMgQ/eY3v9FNN9103sIRQHBR4AFAiHXb9aq7uJMklyxZcqnbp6vOW+ANHjxY3bp1U//+/ZWSkiJJysjIUHp6ulq3bq0ZM2Zo7dq1crlcSk9P10033aQ333xTP/nJT5SUlCRJ6tq1q370ox+pTZs2GjlypCRp2LBhuuOOOzRkyBB16tRJw4cP9xrRjh07VgMGDNDChQs1b948LVmyROnp6Tb8dgD4g/9iAUCItTxxyF3cfc8lSy3LD/p8jO+LO0nq06ePPvroI7333nvq3r27XC6XJKlFixZq2bKlPv3003Mew/PzeFVVVY15Cu4xPPPMM3rllVd033336ZNPPmn0MQDYgw4eANvUnzjQlPGoZ3TsxNjYm5Px7dSq7Is6RZ4ll060usSv51JVVaWYmBhZlqWYmJg634uOjnYXfHbasWOHevbsKUm6/PLLNW/ePL311lv6t3/7N9vPBeD86OABQIjt7XGTO5aV5I5r93a/yedjlJeXf/uzlqWtW7dqwIABGjx4sHbv3u3ep7KyUidOnFBqamqjxteyZUudOXNGZ86c0alTp865z2OPPVZnWZT9+/erd+/ejToPAPvQwQOAEDt0aT99NOAedft0lVqWH9TJVpdoT/eb9NWlfaWPfOvgrVq1SrGxsdq6dat+85vfqE2bNpKkBx98UEuWLNFll12mvXv3avny5YqOjta6detUVFSkdevWqXPnztq1a5feeust3XXXXfrggw9UVFSkzZs3a8CAARozZoxyc3PVqVMnTZ069ZzntyxLQ4cO1e23366WLVuqefPmfAYPCCEKPACN4pT40ynj8NX51r07dGk/94SKHa/tlT7SWcWdt8j79ttvlySNGTOmzvYBAwZowIABZ+1/7bXXqqSkxP148+bN7q/79etX53vXXnutrr32Wq/jf+6557x+H0DTIqIFgDC2bds27du3T88//3yohwLAQejgAUAY69Wrlz7//PNQDwOAw1DgAQiI3TNl/b2FVzgsaNwQb3FzOD8vAKFDRAsAAGAYOngAAmLHbcG8de3q3xrte/W7Xr52/kLVEat/Xm/jpWsHhBdvtyqdmTG0CUfyAzp4AAAAhqHAAwAAMAwRLYBGqR+Zekal3m5V5muE2lAk25j9nLhGHpEsgKYUtgVeaWmpVq9eLZfLpX79+jX61jsAAACmCtsC76WXXlJWVpaioqKUl5dHgQcAAPAdxxR406dP1/r1P8xC6d27t/Ly8vTll18qNzdXPXr0UFFRkXJycpSQkCDLsmRZliSptrbW5/MUFRVp7NixGj16tKZPn+7efuzYMU2ZMkWTJ0/WwIED7XtigOE8o1JfZ7b6GsM2RkOxbDDiT3/X6gNgjoZmzoZq1mx9jinwOnTooE2bNp21/c4779TixYuVlpamDRs2aPbs2Vq6dKkyMzO1cuVKSVJGRobP50lJSdGkSZOUn5+vYcOGqU+fPpKkpKQkjRgxguIOAACEPcfMoo2KOnsoJSUlKikpUVpamiTp6quv1qpVq1RZWanExESNGzdO48aNa3Q8Gx0drfz8fN19992qqKiwZfwAAABO4ZgOXnl5uebMmaPq6mrt2bNHCxYs0L59+9SlSxf3PlFRUWrdurVKSkp0+eWXB3S+lJQUTZw4UVOnTtWyZct8+pnKykpVVlbW2VZTUxPQOACTBCN6bYgdt/eyI2r19Tk7cWYvAN95W8zY234Rv9BxfHy8cnJylJubqzlz5ig7O1tHjx5Vq1at6uzXqlUrff3117acc/z48Tpy5IjWrFnj0/4LFy5UYmJinT9btmyxZSwAAAB2cUwHLzs72/11r169tHfvXiUnJ6usrKzOfqWlpUpOTrbtvHl5eRo2bJj69et33n1nzZql+++/v8622bNn2zYWAP5rqGvXmC6dHR1IXzt1rH0HOJ+vXTsnckQHr7y8XBMmTHDPipWk2NhY9ejRQ8XFxe5tVVVVKisrU8eOHW07d1JSkhYtWqSsrKzz7hsbG6uEhIQ6f6Kjo20bCwAAgB0cUeBVVlZq8+bNcrlckqRt27Zp0KBBuvTSS9WtWzft3LlTkvTuu+8qMzNTzZs3t/X86enpSk1NVUFBga3HBQAACAVHRLRt2rTR4sWLVVBQoOjoaO3YsUOLFi2SJC1fvlw5OTlKS0vTF198ofnz5wd0rv3792vt2rXq27evunfv7t6em5urMWPGBHRsAKFj963Q7JgUQQwLIFQcUeBJ0vDhw8+5vW3btlqxYoVt5+nUqZN7/TxPzZs316pVq2w7DwAAQKg4IqIFAACAfRzTwQOAxmjKNfe8Ra3ctgywn7fZq065FVhDnDI+OngAAACGoYMHAOdRv0vn2dFjIgXQtBrq7jmlc+YUdPAAAAAMQ4EHAABgGCJaAEazYz07AAg3dPAAAAAMQ4EHAABgGCJaAMbxFssy6xVwHm/r3oWSt5m5nmN24gxeOngAAACGocADAAAwDBEtAKMRyQKRoX7Ma0ds6mt07MS4lg4eAACAYejgATCat9uMAQidYE+s8KerZseYgtFJ9AcdPAAAAMNQ4AEAABiGiBaAcXre2M39df018epHtp6Ib4HI05RRsadgR7d08AAAAAxDgQcAAGAYIloARvOMa+vzNb4lugXM4oRbowV7ti0dPAAAAMNQ4AEAABiGiBZAxGpMfAsgcKGKRp0QyZ6P3bc7o4MHAABgGDp4AHAe3O7Mmbx1ZZxyw3cgVOjgAQAAGIYCDwAAwDBEtABwDr7e7oy4tmmFw4flgUDZsUYeHTwAAADDUOABAAAYhogWAM6j/np5rJHnfHavKQaEGzp4AAAAhqHAAwAAMAwRbYC27P+CWV0A4GB2zEgEwg0dPAAAAMPQwQMAAEFBwhU6dPAAAAAMQ4EHAABgGCJaAIBtvEVydkxu8DyGv/FfQz/H5AuYhA4eAACAYSjwAAAADENECwAIiK9Rqd23D/N2DH/iW9bLg0no4AEAABiGDh4A4LzsXs/skbfXK+7kZ+ouqWDzVk35RW9bj28HuzuOQFOigwcAAGAYCjwAAADDENECCLkdr+0N9RCME+hab+F+iyk71svz5OsxiHKdKRjXxemvETp4AAAAhqHAAwAAMAwRLYAm4WsM+/xTE4I8EjP5sxadk5gyY5W19EKnqX/Xdn8MwG508AAAAAxDgQcAAGAYIloATY4YNjLUj8ycckuzpozTTIme/RXs371Tfqd2zE63+7nQwQMAADAMHTwAgG3s7kIEY9JCqNYCZAKG7xMTTP3dNOXzooMHAABgGAo8AAAAwxDRAmhyt9/zbJ3HTLqAryJ90kI48Cfa5lrajw4eAACAYSjwAAAADENEC6BJ9Lyxm/trX29bBu+ceHukpoxQgz0r1dvxnPi791eonguxbHDRwQMAADAMBR4AAIBhiGgBIEyYFAsGQ0O/n2BEgb4u2OtEdoyXeNX56OABAAAYhg4eADhMuHWEPHnr7ISq6xXsW2I5cTJGOP8dgj3o4AEAABiGAg8AAMAwRLQAEGLEaeHHjgkd3HYNwUQHDwAAwDAUeAAAAIYhogWAECCWdQa7Z9hyXeEUdPAAAAAMQ4EHAABgGCJaAABCjBm1sBsdPAAAAMPQwQPQJHa8tjfUQwipYHz4PpxveB8OvHXVgvm7r388fzp69X+Gvx+Rhw4eAACAYSjwAAAADENEC6DJPf/UhFAPwQhOjN28xZpOHK+/wu25eIt5w+25wDd08AAAAAxDgQcAAGAYIloAaAImx5UNCffn6JT16Jy4Rp4Tx4S66OABAAAYhgIPAADAMES0ANAEwj2uNJWv8SLXr2F2LMwM+9HBAwAAMAwdPABN7vZ7nq3zmHXxYIdAO0fh0KXzt1sWDs8N9qKDBwAAYBgKPAAAAMMQ0QJoEj1v7Ob+esdre0M4kqZDLNa0WJvNGbgOzkAHDwAAwDAUeAAAAIYhooVtLjhdEuohIEwkNav7d6Xs0MehGUgQFGze6v46LoTjcKI7BlxZ5/GJr1tq+85vty//pMzWc3mLx02KDZ0ehwbjOjR0TDuev78fq3Di754CDwGrbpagmqhYdf7XY6EeCsJE9+S6jzcWPB6agQRB91APwME27jx7W3RMC8XEXSTJ3gIPiHQUeAhYVezF2vVvT6tZNW/Q8M2/1u2v8/ih6aNCMxAbeHbsnKR+t+x7vo63oZ+3W0zcRYpL7KiZGR3P+f1gTFQxdfJLuD0vu7uPrBFYV9gWeKWlpVq9erVcLpf69eun1NTUUA8polXFXqyq2ItDPQyEiWPVVp3HCe16h2YgNqho6cz/2DT0O/V1vOF8TQCEcYH30ksvKSsrS1FRUcrLy6PAAwAA+I5jCrzCwkIdOnRI1dXV2rJli3Jzc5WUlKTp06dr/fof2qe9e/dWXl6eLMuSZX3bBaitrfX5PEVFRRo7dqxGjx6t6dOnu7cfO3ZMU6ZM0eTJkzVw4ED7nhiA86p/6zJ/+HO7s3CIZuz+8LYTPwxeXzhcF9iLKN5+jinw7rvvPm3fvl1RUVFq27atcnJy9NRTT6lDhw7atGnTWftnZmZq5cqVkqSMjAyfz5OSkqJJkyYpPz9fw4YNU58+fSRJSUlJGjFiBMUdAAAIe44p8EaPHq2amhpFRUUpJSVFxcXFkqSoqHMv1ZeYmKhx48b5da7o6Gjl5+dr7Nix+vvf/664OBYzAAAA5nBMgTd//nz3188//7zuvPNOSVJ5ebnmzJmj6upq7dmzRwsWLFBaWlrA50tJSdHEiRM1depULVu2zKefqaysVGVlZZ1tVo3v8TCAb3netsxf9W935hnzesa1To1pwiEqBRC+HHUni5KSEk2ePFlHjhzRiBEjJEnx8fHKyclRbm6u5syZo+zsbNvON378eB05ckRr1qzxaf+FCxcqMTGxzp+S3btsGw8AAIAdHFXgdezYUcuWLdNtt92mzMxMSVJ2drY7pu3Vq5f27rX3JuV5eXmaM2eODh8+fN59Z82apdLS0jp/Oqb1sHU8AAAAgXJMRLtlyxb1799fkjR06FDdfPPN2rdvn+bPn69nnnlGLpdLkhQbG2vreZOSkrRo0SJlZWVpzJgxXveNjY096/yuaEfVyEDEqB/zeka2nnGtHXGwHYhkfef5u3JqxB7OInEB4Eh8zo6oTg4ePKjrrrtONTU1kr5dsiQxMVEtWrTQ5s2b3cXdtm3bNGjQINvPn56ertTUVBUUFNh+bAAAgKbmiA7eJZdcoieeeEL5+fm64IILtGvXLr322mtKTk7W4sWLVVBQoOjoaO3YsUOLFi0K6Fz79+/X2rVr1bdvX3Xv/sNdI3Nzc8/bwQPgXJ6dulB2gOjUwem8vSY8//56+7scDp2uSOzaeXJEgSdJd9111zm3Dx8+3NbzdOrUyb1+nqfmzZtr1apVtp4LAAAgFBwR0QIAAMA+jungAUAgGopjiEwB35kUVwb7uTj9vYUOHgAAgGEo8AAAAAxDRAsAcBSTYkIgVOjgAQAAGIYCDwAAwDBEtADCktNnsAFAKNHBAwAAMAwdPACAo4TyVnNAQ8ItNaCDBwAAYBgKPAAAAMMQ0QJwtHCLRWCv+tefyBbwDR08AAAAw1DgAQAAGIaIFkDIEcPCV/78XTEp1vX1+Zv0nOEfOngAAACGocADAAAwDBEtAMBo4RZr2vGRBWYf2yOcPz5CBw8AAMAwdPAANLlw/l8xzGVy1yvSJ6dEIjp4AAAAhqHAAwAAMAwRLQAAOItTJ6c4dVxOQwcPAADAMBR4AAAAhiGiBdAkmDmLcEMU6JvGvLb9+V35+94R6e85dPAAAAAMQ4EHAABgGCJaAAAC4C0KjPT4tr5Ij02bEh08AAAAw9DBAwAgSOhYIVTo4AEAABiGAg8AAMAwFHgAAACGocADAAAwDAUeAACAYVyWZVmhHkQ4GzVqlL766iv1799f0dHRoR4OvlNTU6MtW7ZwXRyG6+I8XBNn4ro4kxOuS3FxsV555ZXz7keBF6CysjIlJiaqtLRUCQkJoR4OvsN1cSaui/NwTZyJ6+JM4XRdiGgBAAAMQ4EHAABgGAo8AAAAw1DgBSg2NlZz585VbGxsqIcCD1wXZ+K6OA/XxJm4Ls4UTteFSRYAAACGoYMHAABgmGahHkA4KS0t1erVq+VyudSvXz+lpqaGekgQ18WJuCbOwzVxPq6R84TzNaHAa4SXXnpJWVlZioqKUl5eXlhdaJNxXZyHa+I8XBPn4xo5TzhfEwo8SVVVVerdu7d27NghSfryyy+Vm5urHj16qKioSDk5OUpISJBlWfr+I4u1tbU+H/+ee+7RZ599pr/85S+KiYmRJP3qV79SRUWFXn75Zfd+b7zxhtauXaunnnrKxmcXfv7xj3+ouLhYtbW12r59u0aPHq3+/ftr+vTpWr9+vXu/3r17Ky8vj+vShB566CHNnj3b/ZjXijOd67XSp08fv67JwYMHdeedd6pNmzYaNWqUJOnMmTN65plndN999+lnP/uZrWM3SUPvZXa/brhGvissLNShQ4dUXV2tLVu2KDc3V0lJSbb/++KIa2LB+q//+i8rKSnJ/fj666+3Pv30U8uyLOv999+37rnnHsuyLOv48ePW73//e+v3v/+9tXv3bp+Pf+DAAWvIkCF1thUWFlo33HBDnW2rV6+2vvnmGz+fhRmqq6utn/70p+7HJ06csK666irLsixr6dKl5/wZrkvwFRcXW8uXL7eio6PrbOe14kzneq34e00sy7Lmzp1rFRQU1NlWUVFhLV68OJBhGs3be5ndrxvL4hr56kc/+pFVU1NjWZZlvf7661Z2drZlWfb/+2JZob8mEd/BKy4uVmJiolq1aiVJKikpUUlJidLS0iRJV199tcaOHatFixYpMTFR48aNa/Q52rdvr6qqKh06dEjt2rWTJH3xxRfas2ePysvL3ec+duyYWrdubc8TC1NlZWXavXu3Tp8+rQsuuEAtW7ZUWVmZJCkq6txzgrguwXfZZZdp4sSJevTRR93beK0417leK/5ek4Y0b95c/+///T/bjmeaht7LgvG6aQjX6GyjR49WTU2NoqKilJKSouLiYkn2//vSkKa8JhE/i/bll1/WmDFj3I937dqlLl26uB9HRUWpdevWKikpCeg8N954o1atWiXp2zZvs2bN9LOf/Uxr1qyRJJ08eVLx8fEBncMEF154oQ4cOKALLrhAkvTBBx/o2muvlSSVl5drzpw5+vWvf62bb75Zu3fvDvh8XBf/8VpxrmC8VjxVVlbqxRdfdMfoOFtD72XBet3UxzU6t/nz57t/J88//7zuvPNOScF/zUhNf00iuoP3+uuv6+c//3mdbUePHnV3Cb7XqlUrff3117r88sv9PteoUaM0ZcoU3X333Vq/fr2GDBmijh07avny5crMzFRhYaEyMjL8Pr6JTp06pby8PD3xxBOSpPj4eM2YMUNRUVHatm2bsrOzVVhYGNA5uC7+47XiXMF4rWzdutXdNV27dq0GDhxow0gjg+d72ZtvvhmU143ENfJVSUmJHnnkEVVVVWnEiBGSgvOakUJ7TSK2g3fixAkdOHBAPXr0qLM9OTnZHQl+r7S0VMnJyQGdr0ePHjp48KC++eYbFRUVqXPnzhoyZIg+/PBDnT59WqWlpUpISAjoHCapqKjQE088ocWLF7tfHNnZ2e42eq9evbR3796Az8N18R+vFecKxmvlyiuv1KhRozRq1CgtWbKEzpCP6r+XBet1I3GNfNWxY0ctW7ZMt912mzIzMyUF5zUjhfaaRGyBt2nTJrVt21Zr167V2rVrVVFRobVr1youLs6dyUvfzrAtKytTx44dAz7n9ddfrz//+c/uCxwdHa3rrrtOf/rTn876H10kq6io0JIlS/TAAw8oMTFR5eXlKi8v14QJE9yzmSTZdqsYrot/evTowWvFgYL5WvleXFycbr31VluPaaJzvZd169YtaK8bT1yjc9uyZYv766FDh+r999/Xvn37gv6akZr+mkRsgZeRkaFRo0Zp5MiRGjlypOLi4jRy5Ehdc8016tatm3bu3ClJevfdd5WZmanmzZsHfM5Ro0bp17/+dZ146aabbtLMmTM1fPjwgI9viqVLl2rKlClq0aKFJOnZZ5/V6dOntXnzZrlcLknStm3bNGjQIFvOx3Xxz6WXXsprxYEqKyuD9lqpLy8vLyjHNcW53svat28ftNfNuXCNfnDw4EFdd911qqmpkfTtZK3ExES1aNGiyV4zUtNdk4j+DJ4kHT58WH/961/11VdfadmyZbr55pu1fPly5eTkKC0tTV988YXmz59vy7kGDhyojIwMtW/f3r0tIyNDGRkZdCW+89577+nZZ5/V66+/Lunb/91+8803uvfee7V48WIVFBQoOjpaO3bs0KJFi2w5J9fl/Pbt26fCwkIdPHhQ8+fP10033aQrrriC14oDtWnTxtbXysGDB/XRRx+ppKTE/XGJyspKvfrqq+rbt69NozaPt/cyu183XCPfXHLJJXriiSeUn5+vCy64QLt27dJrr72m5ORk2/99ccI1cVmePUkAAACEvYiNaAEAAExFgQcAAGAYCjwAAADDUOABAAAYhgIPAADAMBR4AAAAhqHAAwDD7dy5Ux9//HGohwGgCVHgATDS8ePHde+998rlcmn06NH68MMPQz2kOpYvX6527drpgw8+CNo5qqqqNHfuXPXv358CD4gwFHgAjNS6dWs9+eSTkqRp06Y5bkX/pKQkJSQknHVnjg0bNth2jpiYGM2bN099+vSx7ZgAwgMFHgCEQGZmpvbs2aPu3bu7t5WXl+uVV16x/VzR0dG2HxOAs1HgAUAIfP311/rrX//qflxdXa3f/e537huhA0AgmoV6AADgBKdOnVJOTo6io6MVHx+vzz77TNnZ2bryyislSRs3blR2drbi4+O1dOlSPfXUU6qqqtL+/fv15JNPqnfv3u5jlZWVadasWWrbtq0k6cCBA0pLS9OyZcv04IMPqmfPnsrJydHOnTtVXFyskydP6re//a3eeecdffPNN3rkkUckSTNnzpRlWZo9e7YWLlyod955Rz/5yU8kSS+88IJmzJihrKws5eTkuM9dWlqqBx54QMePH1dSUpIuu+wy1dbWnvV83377ba1YsUJt27bV119/rWuuuUZTpkwJzi8XQNOzAMBgkqyNGzd63ae2ttYaMWKE9T//8z/ubUePHrV69Ohh7d69272toKDA6tChgzVv3jyrpqbGsizLevHFF63BgwfXOV5WVpa1YsUK9+P//d//tTIyMqz/+7//sz777DPLsizr3XfftS677LI6Pzd37lxr3Lhx5xxjp06drL/97W91tv3qV7+y5s6d635cU1NjDRkyxFq+fLl7W0lJiXXxxRdbBQUF7m3r1q2z+vTpYx07dsz9c9dff7316KOPnvsXBCDsENECiHivvPKKvvzyS2VmZrq3JSUlacKECZoxY0adfS+88ELNnj1bUVHfvn1eeeWV+uc//1lnn1WrVumaa65xPx42bJjefvttdenSRSkpKZIkl8vVqDGea//62woLC7Vnzx7953/+p3tbhw4d1KNHjzr7zZo1Sw888IAuvPBCSVJUVJSmTZumxx9/XJZlNWpcAJyJAg9AxPvLX/6iQYMGnbV90KBBeuONN1RZWenelpSUVKewiomJ0cmTJ8/6Wc+fqa2tVfPmzRUTE2PzyOv629/+pv79+3udVHHq1Clt2rRJP/7xj+ts79mzp44cOaLi4uKgjhFA0+AzeAAi3sGDBzVkyJCztrdv315nzpzR0aNH1b59e5+P99hjj+nJJ5/U888/L5fLpRUrVmjatGm64IIL7Bz2WY4ePark5GSv+3zzzTeqra3VH//4RyUkJLi3W5al9PR0VVdXB3WMAJoGBR6AiFZdXa02bdro9OnTZ33vyJEjatasmTvK9FVUVJTatm2rRx99VM2aNVOHDh00btw4u4bcoMTERB0+fPi8+7hcLk2YMMEdF39v1qxZwRwegCZERAsgos2dO1cjRozQ1q1bz/repk2blJ6erri4uEYd8+GHH1bfvn11++2364EHHvC5uIuJiWlwmZSEhISzitCjR4/WeTxw4EB98sknZ/3skSNH3F/Hx8frxz/+8VmfG5Skjz76yKdxAnA+CjwAEeu3v/2tWrVqpdtuu00nT57UunXr3N8rKytTfn6+lixZ0ujjLly4ULfddps6dOigVq1aqW/fvsrOzta+ffu8/lyXLl20Z88eSd8ud/Lyyy+7v9e7d2+tXbvW/fiVV17R1q1bVVVV5d72i1/8QpZladWqVe5tf/jDH1RcXKwzZ864tz300EN66KGHdODAAfe2HTt2aNOmTY1+rgCcyWUxZQqAgQ4ePKjHHntMS5Ys0ejRozVw4ED3986cOaMNGzbojTfe0Pr16zVkyBCVlZVp3rx5io+PV4sWLbRv3z5NnDhRV1xxhSRp69atuvfee7V9+3YtXrxYd9xxh44dO6Zp06YpPz9f9913nxYuXKjy8nL98pe/1COPPKL4+HgVFxdr27ZtevXVVxUVFaUPPvhAH374oebPn6833nhDs2fP1syZMxUbG6uKigrdeuutOn78uDp06KClS5e64+Hi4mLdfvvtKi8vV6dOnfTzn/9cW7duVVFRkbKzs3XjjTdKkg4dOqTs7GzV1NSoXbt2uuqqq/Tqq68qJiZG999/v/uzhmvXrtWSJUuUmJio5ORkpaSk6L777mv07F4AzkSBBwA2Wrp0qfbu3avf/e53dbYfOnRI3bt31/Hjx0MzMAARhYgWAGzUr18/7dix46zZqHv27HF32QAg2OjgAYDN1q1bp1dffVXt27dXbW2tysvL3bcca+yEDQDwBwUeAACAYYhoAQAADEOBBwAAYBgKPAAAAMNQ4AEAABiGAg8AAMAwFHgAAACGocADAAAwDAUeAACAYSjwAAAADPP/AyeCk/a2F1ZZAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 700x700 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(ncols=1, figsize=(7, 7))\n",
    "\n",
    "im = mesh_mask_study.plot.pcolormesh(\n",
    "        ax=ax,\n",
    "        x=\"x\", \n",
    "        y=\"y\",\n",
    "        vmin = 0, \n",
    "        vmax = 4, \n",
    "        cmap=sns.cubehelix_palette(start=.5, rot=-.5, as_cmap=True),\n",
    "        add_colorbar=False, \n",
    "        )\n",
    "\n",
    "y = list(range(0,185))\n",
    "x = list(range(0,160))\n",
    "\n",
    "rect = mpatches.Rectangle((19, 18), 80, 102, linewidth=1, edgecolor='darkgoldenrod', facecolor='none', label=\"Cropped area\")\n",
    "rect = ax.add_patch(rect)\n",
    "\n",
    "ax.set_yticks(np.arange(min(y)+15, max(y), 40))\n",
    "labels_lat = [item.get_text() for item in ax.get_yticklabels()]\n",
    "lat_label = ['35°N', '50°N', '65°N', '70°N', '85°N']\n",
    "ax.set_yticklabels(lat_label, fontsize=10)\n",
    "ax.set_ylabel('Latitude', fontsize=14)\n",
    "\n",
    "ax.set_xticks(np.arange(min(x), max(x),30))\n",
    "labels_lon = [item.get_text() for item in ax.get_xticklabels()]\n",
    "lon_label = ['40°W', '25°W', '10°W', '5°E', '20°E', '35°E']\n",
    "ax.set_xticklabels(lon_label, fontsize=10)\n",
    "ax.set_xlabel('Longitude', fontsize=14)\n",
    "\n",
    "ax.plot(81.2, 62, \"o\", color='chocolate', markersize=4)\n",
    "ax.text(85, 60.5,'point S', fontsize=10, color='black')\n",
    "\n",
    "for location in ['left','bottom', 'right', 'top']:\n",
    "    ax.spines[location].set_linewidth(0.4)\n",
    "\n",
    "pass\n",
    "\n",
    "plt.savefig('/Users/chiaraciscato/Desktop/GEOMAR/OAEseasonality/out/alkalinityAddition.png')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7a06de06-1326-4a98-9ed6-e83d586107fb",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "geomar",
   "language": "python",
   "name": "geomar"
  },
  "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.18"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}